Замена значений в растрах с помощью GDAL/Numpy
по адресу http://gis-lab.info/qa/remap-gdal.html
Описание и скрипт.
Задача: заменить в растре одно значение на другое.
Для решения можно использовать возможности GDAL для чтения/записи растров и Numpy для операций замены. Данный способ использовался, например, для нормализации ошибочных значений NODATA в данных SRTM 5x5.
Реализация
Все основные моменты работы с растрами изложены в статье «Работа с растрами при помощи GDAL и Python». Остановимся на ключевых деталях реализации:
Так как мы хотим только поменять значение, то важно сохранить растр в том же числовом формате, что и исходный. Для этого нужно при чтении получить исходный числовой формат:
gd = gdal.Open( input )
intype = gd.GetRasterBand( 1 ).DataType
и использовать его же при сохранении результирующего растра:
outRaster = driver.Create( output, xsize, ysize, 1, intype )
Также полезным может оказаться установка значения NoData, если для его установки и заменяются значения:
outRaster.GetRasterBand( 1 ).SetNoDataValue(-32678)
Основная работа производится с помощью функций обработки массивов Numpy, значительно более быстрых чем попиксельная обработка. Цикл поиска-замены значения выглядит следующим образом:
temp1_bool = numpy.equal(raster,inval)
numpy.putmask(raster,temp1_bool,outval)
Скрипт
Запуск скрипта производится следующим образом:
remap.py input_raster output_raster input_value output_value
В качестве аргументов используется имя исходного и конечного растра, а также исходное и конечное значение. Например, для замены значения 255 на -32678:
remap.py input.tif output.tif 255 -32678
Скачать скрипт.
Предупреждение: на данный момент описанная реализация относительно медленная, операции на больших растрах могут занимать значительное время.