Паншерпенинг Ландсата 8

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Паншерпенинг Ландсата 8

Сообщение Александр Мурый » 27 ноя 2013, 08:28

bolotoved писал(а):
rhot писал(а):В GRASS никак нельзя преобразовать 16-бит в 8-бит?
Можно: http://grass.osgeo.org/grass70/manuals/r.rescale.html
Но нужно разбираться как правильно настроить гистограмму.
r.rescale.eq ?
Редактор материалов, модератор форума

Аватара пользователя
rhot
Гуру
Сообщения: 1727
Зарегистрирован: 25 янв 2011, 17:50
Репутация: 194
Ваше звание: доктор
Откуда: Архангельск

Re: Паншерпенинг Ландсата 8

Сообщение rhot » 27 ноя 2013, 12:50

Получилось с помощью r.rescale.eq без атмосферной коррекции.

Привожу результаты паншерпенинга. Хочется добиться такой же насыщенности и контраста как и после атмосферной коррекции, пока не выходит. :(
Вложения
collage.png
collage.png (2.42 МБ) 9880 просмотров
___________(¯`·.¸(¯`·.¸ Scientia potentia est _/ {SILVA}:::{FOSS}:::{GIS} \_ Знание сила ¸.·´¯)¸.·´¯)___________

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Паншерпенинг Ландсата 8

Сообщение Александр Мурый » 27 ноя 2013, 13:04

rhot писал(а): Хочется добиться такой же насыщенности и контраста как и после атмосферной коррекции
i.landsat.rgb тут не поможет?
Редактор материалов, модератор форума

Аватара пользователя
rhot
Гуру
Сообщения: 1727
Зарегистрирован: 25 янв 2011, 17:50
Репутация: 194
Ваше звание: доктор
Откуда: Архангельск

Re: Паншерпенинг Ландсата 8

Сообщение rhot » 27 ноя 2013, 13:15

Пробовал и до, и после - не помогает. Я считаю, что корень проблемы кроется в том, что снимки не обработаны i.landsat.toar.
___________(¯`·.¸(¯`·.¸ Scientia potentia est _/ {SILVA}:::{FOSS}:::{GIS} \_ Знание сила ¸.·´¯)¸.·´¯)___________

Strix
Интересующийся
Сообщения: 31
Зарегистрирован: 23 сен 2013, 18:45
Репутация: 7

Re: Паншерпенинг Ландсата 8

Сообщение Strix » 27 ноя 2013, 15:17

rhot писал(а):Привожу результаты паншерпенинга. Хочется добиться такой же насыщенности и контраста как и после атмосферной коррекции, пока не выходит. :(
756, Метод Brovey -- странный результат...
У меня в той же комбинации каналов без атмосферной коррекции преобразование Бровея в GDAL (в связке gdal_calc.py + gdal_contrast_stretch) дает такое изображение:

Изображение

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Паншерпенинг Ландсата 8

Сообщение Александр Мурый » 27 ноя 2013, 15:24

Strix писал(а): У меня в той же комбинации каналов без атмосферной коррекции преобразование Бровея в GDAL (в связке gdal_calc.py + gdal_contrast_stretch) дает такое изображение:
Распишите подробнее, пож-та, как именно вы делали.
Редактор материалов, модератор форума

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Паншерпенинг Ландсата 8

Сообщение Александр Мурый » 27 ноя 2013, 15:55

Господа, есть предложение сделать совместную статью "Паншарпенинг данных Landsat 8 с помощью GDAL и GRASS". Заготовку на вики сделал.
Редактор материалов, модератор форума

Strix
Интересующийся
Сообщения: 31
Зарегистрирован: 23 сен 2013, 18:45
Репутация: 7

Re: Паншерпенинг Ландсата 8

Сообщение Strix » 27 ноя 2013, 15:56

Мне стыдно показывать этот bash-скрипт гуру-программистам, т.к. в силу незнания bash наверняка реализация получилась через одно место :), тем не менее, с помощью форума gis-lab и гугла я написал нечто, что дает приемлемый для моих целей результат. Скрипт запускаю в Ubuntu в каталоге с обрезанными по нужному региону каналами сцены Landsat 8 (на полной сцене не пробовал -- боюсь, что не хватит оперативки в моем компьютере). В системе должен быть установлен dans-gdal-scripts и, как писал выше, gdal_calc.py должен быть допилен по рецепту Voltrona. Да и еще, это в общем-то вариант скрипта для одной комбинации каналов; в расширенном его варианте я получаю все комбинации каналов сразу.

Код: Выделить всё

#!/bin/bash

mkdir ./mult
for i in *B5.TIF *B6.TIF *B7.TIF; do gdalwarp -tr 15 15 -r bilinear $i ./mult/$i; done
find . -name "*B8.TIF" -exec cp {} ./mult \;
cd ./mult
gdal_calc.py -A *7.TIF  -B *5.TIF -C *6.TIF -D *8.TIF --outfile=7.tif --calc="(D*A)/(A+B+C)"
gdal_calc.py -A *7.TIF  -B *5.TIF -C *6.TIF -D *8.TIF --outfile=5.tif --calc="(D*B)/(A+B+C)"
gdal_calc.py -A *7.TIF  -B *5.TIF -C *6.TIF -D *8.TIF --outfile=6.tif --calc="(D*C)/(A+B+C)"
mkdir ./8bit; for i in *.tif; do gdal_contrast_stretch -histeq 100 $i ./8bit/$i; done
cd ./8bit
gdal_merge.py -separate 7.tif 5.tif 6.tif -o 756Lt.tif

KolesovDmitry
Гуру
Сообщения: 810
Зарегистрирован: 22 авг 2007, 14:58
Репутация: 123
Откуда: Казань

Re: Паншерпенинг Ландсата 8

Сообщение KolesovDmitry » 27 ноя 2013, 16:22

Немного смущает строка такого вида:

Код: Выделить всё

gdal_calc.py -A *7.TIF  -B *5.TIF -C *6.TIF -D *8.TIF --outfile=7.tif --calc="(D*A)/(A+B+C)"
у меня есть подозрение, что она будет работать, только если в каталоге лежит только один файл, оканчивающийся на ...7.TIF. Вы не проверяли, что случится, если таких файлов будет несколько?

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Паншерпенинг Ландсата 8

Сообщение Александр Мурый » 27 ноя 2013, 16:25

Strix писал(а):gdal_calc.py должен быть допилен по рецепту Voltrona.
По вашей ссылке выше по теме я рецептов не нашёл (ну и блог Voltrona закрыт для простых смертных). В чём заключается эта доделка? Выложите, пож-та, допиленную версию <gdal_calc.py>.
Редактор материалов, модератор форума

bolotoved
Гуру
Сообщения: 920
Зарегистрирован: 30 дек 2008, 14:11
Репутация: 236
Откуда: Ханты-Мансийск
Контактная информация:

Re: Паншерпенинг Ландсата 8

Сообщение bolotoved » 27 ноя 2013, 16:28

Strix писал(а):
rhot писал(а):Привожу результаты паншерпенинга. Хочется добиться такой же насыщенности и контраста как и после атмосферной коррекции, пока не выходит. :(
756, Метод Brovey -- странный результат...
У меня в той же комбинации каналов без атмосферной коррекции преобразование Бровея в GDAL (в связке gdal_calc.py + gdal_contrast_stretch) дает такое изображение:
Потому как в QGIS теперь мощные инструменты настройки гистограммы, rhot - хочет увидеть такую же картинку в GRASS, какую Strix автоматом получает в QGIS.

Кстати, чтобы в QGIS получить картинку похожую на паншарп, можно вообще его не делать, а применить фильтр (помоему, умножения), положив панхром под мультиспектр и настроив прозрачность.

Strix
Интересующийся
Сообщения: 31
Зарегистрирован: 23 сен 2013, 18:45
Репутация: 7

Re: Паншерпенинг Ландсата 8

Сообщение Strix » 27 ноя 2013, 17:14

bolotoved писал(а):Потому как в QGIS теперь мощные инструменты настройки гистограммы, rhot - хочет увидеть такую же картинку в GRASS, какую Strix автоматом получает в QGIS.
Я гистограмму в QGIS не настраивал, т.к. у меня изображение, полученное с помощью данного приема, выглядит одинаково где угодно: в любом графическом вьюере, редакторе etc.

Strix
Интересующийся
Сообщения: 31
Зарегистрирован: 23 сен 2013, 18:45
Репутация: 7

Re: Паншерпенинг Ландсата 8

Сообщение Strix » 27 ноя 2013, 17:21

KolesovDmitry писал(а):Немного смущает строка такого вида:

Код: Выделить всё

gdal_calc.py -A *7.TIF  -B *5.TIF -C *6.TIF -D *8.TIF --outfile=7.tif --calc="(D*A)/(A+B+C)"
у меня есть подозрение, что она будет работать, только если в каталоге лежит только один файл, оканчивающийся на ...7.TIF. Вы не проверяли, что случится, если таких файлов будет несколько?
Скрипт я писал изначально для обработки каналов одной сцены Landsat 8, находящихся в одном каталоге. Ситуацию с несколькими сценами в одном каталоге я не предусматривал.
Последний раз редактировалось Strix 27 ноя 2013, 17:34, всего редактировалось 1 раз.

Strix
Интересующийся
Сообщения: 31
Зарегистрирован: 23 сен 2013, 18:45
Репутация: 7

Re: Паншерпенинг Ландсата 8

Сообщение Strix » 27 ноя 2013, 17:29

Александр Мурый писал(а):
Strix писал(а):gdal_calc.py должен быть допилен по рецепту Voltrona.
По вашей ссылке выше по теме я рецептов не нашёл (ну и блог Voltrona закрыт для простых смертных). В чём заключается эта доделка? Выложите, пож-та, допиленную версию <gdal_calc.py>.
Там по ссылке нужно было пролистать вниз до сообщения "Подводные камни gdal_calc" (прошу прощения что не уточнил). Привожу сам рецепт:
...открываем gdal_calc.py в любимом текстовом редакторе, и строки 230-232

Код: Выделить всё

myval=BandReadAsArray(myFiles[i].GetRasterBand(myBands[i]),
                      xoff=myX, yoff=myY,
                      win_xsize=nXValid, win_ysize=nYValid)
приводим к виду

Код: Выделить всё

myval=BandReadAsArray(myFiles[i].GetRasterBand(myBands[i]),
                      xoff=myX, yoff=myY,
                      win_xsize=nXValid, win_ysize=nYValid).astype(numpy.float32)
Сам файл:
gdal_calc.py

Аватара пользователя
rhot
Гуру
Сообщения: 1727
Зарегистрирован: 25 янв 2011, 17:50
Репутация: 194
Ваше звание: доктор
Откуда: Архангельск

Re: Паншерпенинг Ландсата 8

Сообщение rhot » 27 ноя 2013, 18:38

Делал как здесь:

Код: Выделить всё

r.info -r map=band1
r.recode input=band1 output=band1_recoded rules=- << EOF
0.:1.65006412008:0:255
EOF
r.colors map=band1_recoded color=grey
Автор утверждает, что i.pansharpen путает цвета. Я проверил оба способа - результат один и тот же. По мне так, для Landsat 8 лучше делать паншерпенинг без коррекции методом IHS. Оцените сами...
collage.png
collage.png (4.39 МБ) 9780 просмотров
Для статьи предлагаю взять одну и ту же сцену. У Strix в примере, похоже, земли совсем другого назначения - поля и луга, у меня же - чисто лесфонд. :)
___________(¯`·.¸(¯`·.¸ Scientia potentia est _/ {SILVA}:::{FOSS}:::{GIS} \_ Знание сила ¸.·´¯)¸.·´¯)___________

Ответить

Вернуться в «GRASS»

Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей и 4 гостя