11

Is there a way to perform the same task as the gdalbuildvrt utility using GDAL Python bindings? So far I have not found any way to do this other than creating a vrt of a single dataset and manually editing the xml. I would like to create a vrt from multiple rasters (essentially performing a mosaic). Is this possible using pure Python? My other option is to use subprocess to simply call gdalbuildvrt.

10

Honestly it's easier to do this by using gdalbuildvrt in a subprocess or os.system.

Should you wish to do this through Python it can be done. Using the standard dataset creation methods within GDAL Python we can easily create the base dataset VRT.

from osgeo import gdal

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

Note that we are creating the dataset with no bands initially. From the documentation on VRTs that VRT datasets are one of the few dataset types that can accept AddBand arguments.

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

Now for each band we have to set the metadata items manually:

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 takes two arguments, the first a string of the metadata item, the second the item itself. This means that you can't subset a metadata item, so for data sources you have to set the entire contents as a string.

Note that we can use this method to create complex sources (ComplexSource) that contain look-up-tables of values, Kernel filter sources (KernelFilteredSource) of arbitrary sizes and shapes, and Mask Bands (MaskBand).

  • Thank you @om_henners - I ended up using subprocess to call gdalbuildvrt. I am more experienced with Python rather than command line so I was hoping I could do this directly in Python, but its not worth the trouble of messing with creating XML strings as you described. It's certainly good to know that I can do that if need in the future though. – Brian Dec 17 '12 at 3:19
  • Just found a use-case for having a python equivalent: adding unsupported features. For example the vrt file format supports an overviews element, but gdalbuildvrt doesn't use it. Thanks for providing a stub how this might be added in python. – matt wilkie Feb 4 '14 at 20:30
  • @om_henners is there any way to drv.CreateCopy('path/to/file.vrt',input_ds) with absolute path to input_ds file in python? there is relativeToVRT="1" option, but how change it or set while createing VRT? – Dmitriy Litvinov Mar 21 '17 at 7:03
2

Since GDAL 2.1 the CLI tools are available as library functions, and in fact that's what the CLI tools now call internally.

For example:

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

Is the equivalent of:

from osgeo import gdal

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

The available CLI options directly map to the parameters of BuildVRTOptions, plus there's some extras like progress callbacks available.

1

The answer of @rcoup only worked for me, if modify it as follows:

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

Otherwise, the file is not written to disk.

  • JensL thanks! can you explain the intuition of my_vrt = None to write to disk? Just seems really odd – mmann1123 Mar 28 at 14:25
  • 1
    @mmann1123: Otherwise it didn't work and I had in mind that the GDAL API Tutorial said: "Note that the CreateCopy() method returns a writable dataset, and that it must be closed properly to complete writing and flushing the dataset to disk. In the Python case this occurs automatically when "dst_ds" goes out of scope." SInce there is no closing for python, you have to bring your vrt out of scope, by assign it to None. – JensL Apr 2 at 6:52

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.