Пример использования утилит GDAL для вычисления яркостных характеристик снимков (на примере выявления рубок по зимним данным Landsat)
Пример решения технических задач, возникающих при выявлении изменений отражательных свойств поверхности по разновременным материалам ДЗЗ
Одним из способов выявления изменений между сериями разновременных снимков является проведение различных математических операций между значениями яркости пикселей. Настоящая статья создана с целью помочь начинающим пользователям освоить основы техники данного приема с использованием свободного ПО и является опытом её автора. Для примера выбрана достаточно простая практическая задача и соответствующие ей исходные данные.
Используемое ПО
Командная строка Windows 7, Python27, GDAL 1.11.1
Исходные данные
Зимние безоблачные снимки Landsat 8 OLI панхроматический(_B8) канал одной сцены формат TIF, уровень обработки L1T в проекции UTM(отмечу что для моего региона в средней полосе Европейской части, удалось найти порядка по 1 такому снимку за весь период работы сенсора(с апреля 2013 г.) в периоде февраль-март за 2015 и 2016 г.).
Краткое описание действий
Участки вырубленные за период между датами съемок на более позднем снимки будут иметь значительно большие значения яркости, чем на раннем (естественно при условии наличия снежного покрова и не слишком большого временного интервала). Расчет разности яркости и последующая классификация результата по пороговому значению осуществляется с использованием утилиты gdal_calc.py.
Создание общей области покрытия
Для корректной работы утилиты gdal_calc.py с несколькими растрами необходимо, чтобы входящие растры имели одинаковый экстент. Получить сведения об экстенте можно используя утилиту gdalinfo.exe:
gdalinfo rastr
- здесь и далее приведен пример синтаксиса командной строки
- rastr - полный путь к набору растровых данных
Результат работы утилиты, будет представлен в виде текста в командной строке:
Необходимо определить максимальные значения для левой и нижней, и минимальные для правой и верхней границ используемых изображений. Обрезка растров производится с помощью утилиты gdalwarp.exe:
gdalwarp.exe -te Xmin Ymin Xmax Ymax -tr 15 -15 inputRastr outputRastr
- -te экстент обрезки
- -tr размер ячейки по x,y
Для получения правильного результата необходимо также, чтобы валидные значения яркости одного растра не попадали на значения NODATA (0) другого, т.е. получить область пересечения растров. Данную операцию необходимо провести для каждого из растров приведенных к одинаковому экстенту используя gdal_calc.py:
gdal_calc.py -A inprastr1 -B inprastr2 --outfile=outrastr1 --calc="A*(A>0)*(B>0)" --NoDataValue=0 gdal_calc.py -A inprastr1 -B inprastr2 --outfile=outrastr2 --calc="B*(A>0)*(B>0)" --NoDataValue=0