Odpowiednik gdalbuildvrt w języku Python

12

Czy istnieje sposób wykonania tego samego zadania, co narzędzie gdalbuildvrt przy użyciu powiązań GDAL Python? Do tej pory nie znalazłem innego sposobu niż utworzenie vrt pojedynczego zestawu danych i ręczna edycja pliku XML. Chciałbym utworzyć vrt z wielu rastrów (zasadniczo wykonujących mozaikę). Czy jest to możliwe przy użyciu czystego Pythona? Inną opcją jest użycie podprocesu do wywołania gdalbuildvrt.

Brian
źródło

Odpowiedzi:

10

Szczerze mówiąc, łatwiej jest to zrobić, używając gdalbuildvrt w a subprocesslub os.system.

Jeśli chcesz to zrobić za pomocą Pythona, możesz to zrobić. Używając standardowych metod tworzenia zestawów danych w GDAL Python, możemy łatwo stworzyć podstawowy zestaw danych VRT .

from osgeo import gdal

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

Zauważ, że tworzymy zestaw danych początkowo bez pasm. Z dokumentacji VRT wynika, że zestawy danych VRT są jednym z niewielu typów zestawów danych, które mogą przyjmować AddBandargumenty.

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

Teraz dla każdego pasma musimy ręcznie ustawić elementy metadanych:

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)

SetMetadatItempobiera dwa argumenty, pierwszy ciąg elementu metadanych, drugi sam element. Oznacza to, że nie można podgrupować elementu metadanych, dlatego w przypadku źródeł danych należy ustawić całą zawartość jako ciąg.

Zauważ, że możemy użyć tej metody do tworzenia złożonych źródeł ( ComplexSource), które zawierają tabele przeglądowe wartości, źródła filtrów jądra ( KernelFilteredSource) o dowolnych rozmiarach i kształtach oraz maski pasków ( MaskBand).

om_henners
źródło
Dziękuję @om_henners - Skończyłem używać podprocesu do wywołania gdalbuildvrt. Mam większe doświadczenie z Pythonem niż z wierszem poleceń, więc miałem nadzieję, że mogę to zrobić bezpośrednio w Pythonie, ale nie warto kłopotać się tworzeniem ciągów XML zgodnie z opisem. Z pewnością dobrze jest wiedzieć, że mogę to zrobić, jeśli zajdzie taka potrzeba w przyszłości.
Brian
Właśnie znalazłem przypadek użycia równoważnego języka python: dodawanie nieobsługiwanych funkcji. Na przykład format pliku vrt obsługuje overviewselement, ale gdalbuildvrt go nie używa. Dziękujemy za podanie kodu pośredniczącego, w jaki sposób można go dodać w Pythonie.
matt wilkie
@om_henners jest jakiś sposób na drv.CreateCopy ('path / to / file.vrt', input_ds) z bezwzględną ścieżką do pliku input_ds w pythonie? istnieje opcja relativeToVRT = „1”, ale jak to zmienić lub ustawić podczas tworzenia VRT?
Dmitriy Litvinov
8

Od wersji GDAL 2.1 narzędzia CLI są dostępne jako funkcje biblioteczne, a tak naprawdę narzędzia CLI wywołują teraz wewnętrznie.

Na przykład:

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

Jest odpowiednikiem:

from osgeo import gdal

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

Te dostępne opcje CLI bezpośrednio mapować do parametrów BuildVRTOptions , plus tam kilka dodatków, takich jak callbacków postępu dostępny.

rcoup
źródło
Wydaje się, że dotyczy to tylko niektórych narzędzi CLI. Na przykład próbuję nakłonić gdaladdo do pracy, ale nie pojawia się. To samo z gdalwarp. Czy wiesz, czy oni również zamierzają je wspierać? Byłby bardzo pomocny.
fpolig01
@ fpolig01 większość z nich tam jest - patrz RegenerateOverviews()i Warp()w dokumentacji API . Argumenty zasadniczo pasują do poleceń CLI.
rcoup
@rccoup Dziękujemy za odpowiedź. Czy RegenerateOverviews () to to samo co gdaladdo? Czy masz przykład tego działania? Próbuję zrobić coś podobnego do gdaladdo -r przeciętnego „D: \ image.tif”
fpolig01
@ fpolig01 ten post sugerujeBuildOverviews() (czego właściwie szukałem, kiedy znalazłem RegenerateOverviews) - może spróbować?
rcoup
8

Odpowiedź @rcoup działała tylko dla mnie, jeśli zmodyfikuję ją w następujący sposób:

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

W przeciwnym razie plik nie zostanie zapisany na dysku.

JensL
źródło
JensL dzięki! czy potrafisz wyjaśnić intuicję my_vrt = None, aby zapisać na dysk? Po prostu wydaje się naprawdę dziwny
mmann1123
3
@ mmann1123 : W przeciwnym razie to nie zadziałało, a ja miałem na uwadze, że samouczek interfejsu API GDAL powiedział: „Zauważ, że metoda CreateCopy () zwraca zapisywalny zestaw danych i że musi zostać poprawnie zamknięty, aby zakończyć zapisywanie i opróżnianie zestawu danych na dysk W przypadku Pythona dzieje się to automatycznie, gdy „dst_ds” wykracza poza zakres. ” Gdy nie ma już closingPythona, musisz sprowadzić vrtgo z zakresu, przypisując go do None.
JensL
Właściwie właśnie naprawili ten problem (patrz osgeo-org.1560.x6.nabble.com/… )
umbe1987