Python-эквивалент gdalbuildvrt

12

Есть ли способ выполнить ту же задачу, что и утилита gdalbuildvrt, используя привязки GDAL Python? До сих пор я не нашел никакого другого способа сделать это, кроме как создать vrt одного набора данных и вручную отредактировать xml. Я хотел бы создать VRT из нескольких растров (по сути, выполняя мозаику). Возможно ли это с использованием чистого Python? Другой вариант - использовать подпроцесс для простого вызова gdalbuildvrt.

Брайан
источник

Ответы:

10

Честно говоря, это проще сделать с помощью gdalbuildvrt в subprocessили os.system.

Если вы хотите сделать это через Python, это можно сделать. Используя стандартные методы создания набора данных в GDAL Python, мы можем легко создать базовый набор данных VRT .

from osgeo import gdal

drv = gdal.GetDriverByName("VRT")
vrt = drv.Create("test.vrt", x_size, y_size, 0)

Обратите внимание, что изначально мы создаем набор данных без полос. Из документации по VRTS что VRT наборы данные один из немногих типов наборов данных , которые могут принимать AddBandаргументы.

vrt.AddBand(gdal.GDT_Float32)
band = vrt.GetRasterBand(1)

Теперь для каждой полосы мы должны установить элементы метаданных вручную:

simple_source = '<SourceFilename relativeToVRT="1">%s</SourceFilename>' % source_path + \
    '<SourceBand>%i</SourceBand>' % source_band + \
    '<SourceProperties RasterXSize="%i" RasterYSize="%i" DataType="Real" BlockXSize="%i" BlockYSize="%i"/>' % (x_size, y_size, x_block, y_block) + \
    '<SrcRect xOff="%i" yOff="%i" xSize="%i" ySize="%i"/>' % (x_offset, y_offset, x_source_size, y_source_size) + \
    '<DstRect xOff="%i" yOff="%i" xSize="%i" ySize="%i"/>' % (dest_x_offset, dest_y_offset, x_dest_size, y_dest_size)
band.SetMetadataItem("SimpleSource", simple_source)
band.SetMetadataItem("NoDataValue", -9999)

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

Обратите внимание, что мы можем использовать этот метод для создания сложных источников ( ComplexSource), которые содержат справочные таблицы значений, источники фильтра ядра ( KernelFilteredSource) произвольных размеров и форм и маскирующие полосы ( MaskBand).

om_henners
источник
Спасибо @om_henners - я использовал подпроцесс для вызова gdalbuildvrt. У меня больше опыта работы с Python, чем с командной строкой, поэтому я надеялся, что смогу сделать это непосредственно в Python, но это не стоит того, чтобы возиться с созданием строк XML, как вы описали. Конечно, приятно знать, что я могу сделать это, если понадобится в будущем.
Брайан
Только что нашел вариант использования для получения эквивалента Python: добавление неподдерживаемых функций. Например, формат файла vrt поддерживает overviewsэлемент, но gdalbuildvrt не использует его. Спасибо за предоставление заглушки, как это может быть добавлено в python.
Мэтт Уилки
@om_henners есть ли способ drv.CreateCopy ('путь / к / file.vrt', input_ds) с абсолютным путем к файлу input_ds в python? есть параметр absoluteToVRT = "1", но как его изменить или установить при создании VRT?
Дмитрий Литвинов
8

Начиная с GDAL 2.1 инструменты CLI доступны в виде библиотечных функций, и фактически это то, что инструменты CLI теперь вызывают внутри.

Например:

gdalbuildvrt -r cubic -addalpha my.vrt one.tif two.tif

Является эквивалентом:

from osgeo import gdal

vrt_options = gdal.BuildVRTOptions(resampleAlg='cubic', addAlpha=True)
gdal.BuildVRT('my.vrt', ['one.tif', 'two.tif'], options=vrt_options)

В Доступные опции CLI непосредственно сопоставить с параметрами BuildVRTOptions , плюс есть некоторые дополнительные , как прогресс обратных вызовов доступны.

rcoup
источник
Похоже, что это только для определенных инструментов CLI. Например, я пытаюсь заставить gdaladdo работать, но он не появляется. То же самое с gdalwarp. Вы знаете, планируют ли они поддержать их также? Было бы очень полезно.
fpolig01
@ fpolig01 большинство из них там - смотрите RegenerateOverviews()и Warp()в справочнике по API . Аргументы обычно соответствуют командам CLI.
встреча
@rccoup Спасибо за ответ. Является ли RegenerateOverviews () таким же, как gdaladdo? У вас есть пример того, как это работает? Я пытаюсь сделать что-то похожее на gdaladdo -r среднего "D: \ image.tif"
fpolig01
@ fpolig01 этот пост предлагаетBuildOverviews() (что на самом деле я искал, когда нашел RegenerateOverviews) - может быть, попробовать?
rcoup
8

Ответ @rcoup работал только для меня, если изменить его следующим образом:

from osgeo import gdal 

vrt_options = gdal.BuildVRTOptions(resampleAlg='cubic', addAlpha=True)
my_vrt = gdal.BuildVRT('my.vrt', ['one.tif', 'two.tif'], options=vrt_options)
my_vrt = None

В противном случае файл не записывается на диск.

JensL
источник
JensL спасибо! Можете ли вы объяснить интуицию my_vrt = None для записи на диск? Просто кажется действительно странным
mmann1123
3
@ mmann1123 : В противном случае это не сработало, и я имел в виду, что в учебнике по API GDAL сказано: «Обратите внимание, что метод CreateCopy () возвращает доступный для записи набор данных, и его необходимо правильно закрыть для завершения записи и сброса набора данных на диск. В случае с Python это происходит автоматически, когда «dst_ds» выходит из области видимости ». Так как closingдля python нет, вы должны вывести свой vrtиз области видимости, назначив его None.
JensL
На самом деле, они только что исправили эту проблему (см. Osgeo-org.1560.x6.nabble.com/… )
umbe1987