Растровая арифметика в GDAL

Растровый калькулятор это круто, но когда возникает необходимость обработать кучу растров начинаются проблемы. Ну нет в калькуляторе  (ни в одном, ни в другом) пакетного режима, нет. Вот в моем TODO он есть. Правда, от осознания этого факта легче не становится, когда у тебя 100500 растров для которых надо расcчитать, к примеру NDVI или выполнить замену пикселей по хитрому условию.

Тут на помощь приходит GDAL и Python. Буквально за 10 минут пишется скрипт, благо есть документация и на английском, и на русском. Но опять же, писать скрипты на каждый чих как-то не очень целесообразно.

Открою маленький секрет. В GDAL 1.8.0 есть замечательный инструмент с незатейливым названием gdal_calc.py. Этот небольшой (около 300 строк) скрипт на Python работает с растрами, имеющими одинаковые размеры (проверка на соответствие проекций не выполняется), и поддерживает базовые арифметические и логические действия. Пользоваться просто:

# сумма двух растров
gdal_calc.py -A input1.tif -B input2.tif --outfile=result.tif --calc="A+B"
# среднее значение двух растров
gdal_calc.py -A input1.tif -B input2.tif --outfile=result.tif --calc="(A+B)/2"
# разность каналов
gdal_calc.py -A input.tif --A_band=1 -B input.tif --B_band=2 --outfile=result.tif --calc="A-B"

На вход можно подавать до 26 растров, этого должно хватить для большинства применений. Как видно из примеров, есть поддержка скобок, можно манипулировать каналами, а консольная природа позволяет легко использовать скрипт для пакетной обработки. Его бы чуть-чуть допилить и получится отличный консольный калькулятор.

Мітки:
4 коментарі в “Растровая арифметика в GDAL
  1. Boris сказав:

    Спасибо вещь занятная. Но и тут, и в тексте скрипта, операции, которые можно использовать и их результаты приведены как очевидные. А можно пальцем ткнуть в место, где они описаны. Ну или их доступное подмножество. Потому как в документации по 'Numpy' как-то с очевидностью в операторах небогато.

  2. alex сказав:

    Можно, хотя тыкать особо и некуда. Но дело не в «небогатой очевидности» документации.

    Поддерживаются только базовые арифметические действия (+, -, *, /) и операции сравнения (>, < , ==, !=, >=, <=).

  3. Boris сказав:

    Я не знаток питона, поэтому подскажите, что есть результат сравнения? В бытность С – для "Да" это было -1. Здесь судя по примерам это с очевидностью не так.

    И второе для файлов с несколькими каналами как организован расчет выходных данных? Или подразумевается, что выходной канал всегда один?

  4. alex сказав:

    В числовом представлении истина (True) это 1, ложь (False) — 0. Так, выражение A*(A==100) занулит все пиксели растра значение которых не равно 100.

    Для файлов с несколькими каналами ситуация такая. Желаемый канал исходного растра указывается при помощи параметра –A_band, где A — идентификатор соответствующего растра. Если такого параметра нет и растр многоканальный — берется первый канал. Результат всегда пишется в 1й канал выходного растра. Это может быть как новый одноканальный растр, так и существующий одно- или многоканальный растр с теми же размерами.

    P.S.: опять комментарий в спам засунуло…

Залишити відповідь

Ваша e-mail адреса не оприлюднюватиметься. Обов’язкові поля позначені *

*