#!/usr/bin/env pytest # -*- coding: utf-8 -*- ############################################################################### # $Id$ # # Project: GDAL/OGR Test Suite # Purpose: Misc tests of VRT driver # Author: Even Rouault # ############################################################################### # Copyright (c) 2013, Even Rouault # # Permission is hereby granted, free of charge, to any person obtaining a # copy of this software and associated documentation files (the "Software"), # to deal in the Software without restriction, including without limitation # the rights to use, copy, modify, merge, publish, distribute, sublicense, # and/or sell copies of the Software, and to permit persons to whom the # Software is furnished to do so, subject to the following conditions: # # The above copyright notice and this permission notice shall be included # in all copies or substantial portions of the Software. # # THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS # OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, # FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL # THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER # LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING # FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER # DEALINGS IN THE SOFTWARE. ############################################################################### import os import shutil import struct import sys import gdaltest from osgeo import gdal, osr ############################################################################### # Test linear scaling def test_vrtmisc_1(): ds = gdal.Translate("", "data/byte.tif", options="-of MEM -scale 74 255 0 255") cs = ds.GetRasterBand(1).Checksum() ds = None assert cs == 4323, "did not get expected checksum" ############################################################################### # Test power scaling def test_vrtmisc_2(): ds = gdal.Translate( "", "data/byte.tif", options="-of MEM -scale 74 255 0 255 -exponent 2.2" ) cs = ds.GetRasterBand(1).Checksum() ds = None assert cs == 4159, "did not get expected checksum" ############################################################################### # Test power scaling (not in VRT file) def test_vrtmisc_3(): ds = gdal.Open( """ data/byte.tif 1 2.2 0 255 """ ) cs = ds.GetRasterBand(1).Checksum() ds = None assert cs == 4159, "did not get expected checksum" ############################################################################### # Test multi-band linear scaling with a single -scale occurrence. def test_vrtmisc_4(): # -scale specified once applies to all bands ds = gdal.Translate( "", "data/byte.tif", options="-of MEM -scale 74 255 0 255 -b 1 -b 1" ) cs1 = ds.GetRasterBand(1).Checksum() cs2 = ds.GetRasterBand(2).Checksum() ds = None assert cs1 == 4323, "did not get expected checksum" assert cs2 == 4323, "did not get expected checksum" ############################################################################### # Test multi-band linear scaling with -scale_XX syntax def test_vrtmisc_5(): # -scale_2 applies to band 2 only ds = gdal.Translate( "", "data/byte.tif", options="-of MEM -scale_2 74 255 0 255 -b 1 -b 1" ) cs1 = ds.GetRasterBand(1).Checksum() cs2 = ds.GetRasterBand(2).Checksum() ds = None assert cs1 == 4672, "did not get expected checksum" assert cs2 == 4323, "did not get expected checksum" ############################################################################### # Test multi-band linear scaling with repeated -scale syntax def test_vrtmisc_6(): # -scale repeated as many times as output band number ds = gdal.Translate( "", "data/byte.tif", options="-of MEM -scale 0 255 0 255 -scale 74 255 0 255 -b 1 -b 1", ) cs1 = ds.GetRasterBand(1).Checksum() cs2 = ds.GetRasterBand(2).Checksum() ds = None assert cs1 == 4672, "did not get expected checksum" assert cs2 == 4323, "did not get expected checksum" ############################################################################### # Test multi-band power scaling with a single -scale and -exponent occurrence. def test_vrtmisc_7(): # -scale and -exponent, specified once, apply to all bands ds = gdal.Translate( "", "data/byte.tif", options="-of MEM -scale 74 255 0 255 -exponent 2.2 -b 1 -b 1", ) cs1 = ds.GetRasterBand(1).Checksum() cs2 = ds.GetRasterBand(2).Checksum() ds = None assert cs1 == 4159, "did not get expected checksum" assert cs2 == 4159, "did not get expected checksum" ############################################################################### # Test multi-band power scaling with -scale_XX and -exponent_XX syntax def test_vrtmisc_8(): # -scale_2 and -exponent_2 apply to band 2 only ds = gdal.Translate( "", "data/byte.tif", options="-of MEM -scale_2 74 255 0 255 -exponent_2 2.2 -b 1 -b 1", ) cs1 = ds.GetRasterBand(1).Checksum() cs2 = ds.GetRasterBand(2).Checksum() ds = None assert cs1 == 4672, "did not get expected checksum" assert cs2 == 4159, "did not get expected checksum" ############################################################################### # Test multi-band linear scaling with repeated -scale and -exponent syntax def test_vrtmisc_9(): # -scale and -exponent repeated as many times as output band number ds = gdal.Translate( "", "data/byte.tif", options="-of MEM -scale 0 255 0 255 -scale 74 255 0 255 -exponent 1 -exponent 2.2 -b 1 -b 1", ) cs1 = ds.GetRasterBand(1).Checksum() cs2 = ds.GetRasterBand(2).Checksum() ds = None assert cs1 == 4672, "did not get expected checksum" assert cs2 == 4159, "did not get expected checksum" ############################################################################### # Test metadata serialization (#5944) def test_vrtmisc_10(): gdal.FileFromMemBuffer( "/vsimem/vrtmisc_10.vrt", """ bar baz foo.tif 1 """, ) ds = gdal.Open("/vsimem/vrtmisc_10.vrt", gdal.GA_Update) # to trigger a flush ds.SetMetadata(ds.GetMetadata()) ds = None ds = gdal.Open("/vsimem/vrtmisc_10.vrt", gdal.GA_Update) assert ds.GetMetadata() == {"foo": "bar"} assert ds.GetMetadata("some_domain") == {"bar": "baz"} assert ds.GetMetadata_List("xml:a_xml_domain")[0] == "\n" # Empty default domain ds.SetMetadata({}) ds = None ds = gdal.Open("/vsimem/vrtmisc_10.vrt") assert ds.GetMetadata() == {} assert ds.GetMetadata("some_domain") == {"bar": "baz"} assert ds.GetMetadata_List("xml:a_xml_domain")[0] == "\n" ds = None gdal.Unlink("/vsimem/vrtmisc_10.vrt") ############################################################################### # Test relativeToVRT is preserved during re-serialization (#5985) def test_vrtmisc_11(): f = open("tmp/vrtmisc_11.vrt", "wt") f.write( """ ../data/byte.tif 1 """ ) f.close() ds = gdal.Open("tmp/vrtmisc_11.vrt", gdal.GA_Update) # to trigger a flush ds.SetMetadata(ds.GetMetadata()) ds = None data = open("tmp/vrtmisc_11.vrt", "rt").read() gdal.Unlink("tmp/vrtmisc_11.vrt") assert '../data/byte.tif' in data ############################################################################### # Test set/delete nodata def test_vrtmisc_12(): gdal.FileFromMemBuffer( "/vsimem/vrtmisc_12.vrt", """ foo.tif 1 """, ) ds = gdal.Open("/vsimem/vrtmisc_12.vrt", gdal.GA_Update) ds.GetRasterBand(1).SetNoDataValue(123) ds = None ds = gdal.Open("/vsimem/vrtmisc_12.vrt", gdal.GA_Update) assert ds.GetRasterBand(1).GetNoDataValue() == 123 assert ds.GetRasterBand(1).DeleteNoDataValue() == 0 ds = None ds = gdal.Open("/vsimem/vrtmisc_12.vrt") assert ds.GetRasterBand(1).GetNoDataValue() is None ds = None gdal.Unlink("/vsimem/vrtmisc_12.vrt") ############################################################################### # Test CreateCopy() preserve NBITS def test_vrtmisc_13(): ds = gdal.Open("data/oddsize1bit.tif") out_ds = gdal.GetDriverByName("VRT").CreateCopy("", ds) assert out_ds.GetRasterBand(1).GetMetadataItem("NBITS", "IMAGE_STRUCTURE") == "1" ############################################################################### # Test SrcRect/DstRect are serialized as integers def test_vrtmisc_14(): src_ds = gdal.GetDriverByName("GTiff").Create( "/vsimem/vrtmisc_14_src.tif", 123456789, 1, options=["SPARSE_OK=YES", "TILED=YES"], ) gdal.GetDriverByName("VRT").CreateCopy("/vsimem/vrtmisc_14.vrt", src_ds) src_ds = None fp = gdal.VSIFOpenL("/vsimem/vrtmisc_14.vrt", "rb") content = gdal.VSIFReadL(1, 10000, fp).decode("latin1") gdal.VSIFCloseL(fp) gdal.Unlink("/vsimem/vrtmisc_14_src.tif") gdal.Unlink("/vsimem/vrtmisc_14.vrt") assert ( ' 0 data/vrtmisc16_tile1.tif 1 0 data/vrtmisc16_tile2.tif 1 0 """, ) gdal.GetDriverByName("GTiff").CreateCopy( "/vsimem/vrtmisc_16.tif", gdal.Open("/vsimem/vrtmisc_16.vrt") ) ds = gdal.Open("/vsimem/vrtmisc_16.tif") cs = ds.GetRasterBand(1).Checksum() assert cs == 206 gdal.Unlink("/vsimem/vrtmisc_16.tif") gdal.Unlink("/vsimem/vrtmisc_16.vrt") ############################################################################### # Check that the serialized xml:VRT doesn't include itself (#6767) def test_vrtmisc_17(): ds = gdal.Open("data/byte.tif") vrt_ds = gdal.GetDriverByName("VRT").CreateCopy("/vsimem/vrtmisc_17.vrt", ds) xml_vrt = vrt_ds.GetMetadata("xml:VRT")[0] vrt_ds = None gdal.Unlink("/vsimem/vrtmisc_17.vrt") assert "xml:VRT" not in xml_vrt ############################################################################### # Check GetMetadata('xml:VRT') behaviour on a in-memory VRT copied from a VRT def test_vrtmisc_18(): ds = gdal.Open("data/byte.vrt") vrt_ds = gdal.GetDriverByName("VRT").CreateCopy("", ds) xml_vrt = vrt_ds.GetMetadata("xml:VRT")[0] assert gdal.GetLastErrorMsg() == "" vrt_ds = None assert ( 'data/byte.tif' in xml_vrt or 'data\\byte.tif' in xml_vrt ) ############################################################################### # Check RAT support def test_vrtmisc_rat(): ds = gdal.Translate("/vsimem/vrtmisc_rat.tif", "data/byte.tif", format="MEM") rat = gdal.RasterAttributeTable() rat.CreateColumn("Ints", gdal.GFT_Integer, gdal.GFU_Generic) ds.GetRasterBand(1).SetDefaultRAT(rat) vrt_ds = gdal.GetDriverByName("VRT").CreateCopy("/vsimem/vrtmisc_rat.vrt", ds) xml_vrt = vrt_ds.GetMetadata("xml:VRT")[0] assert gdal.GetLastErrorMsg() == "" vrt_ds = None assert '' in xml_vrt vrt_ds = gdal.Translate( "/vsimem/vrtmisc_rat.vrt", ds, format="VRT", srcWin=[0, 0, 1, 1] ) xml_vrt = vrt_ds.GetMetadata("xml:VRT")[0] assert gdal.GetLastErrorMsg() == "" vrt_ds = None assert '' in xml_vrt ds = None vrt_ds = gdal.Open("/vsimem/vrtmisc_rat.vrt", gdal.GA_Update) rat = vrt_ds.GetRasterBand(1).GetDefaultRAT() assert rat is not None and rat.GetColumnCount() == 1 vrt_ds.GetRasterBand(1).SetDefaultRAT(None) assert vrt_ds.GetRasterBand(1).GetDefaultRAT() is None vrt_ds = None ds = None gdal.Unlink("/vsimem/vrtmisc_rat.vrt") ############################################################################### # Check ColorTable support def test_vrtmisc_colortable(): ds = gdal.Translate("", "data/byte.tif", format="VRT") ct = gdal.ColorTable() ct.SetColorEntry(0, (255, 255, 255, 255)) ds.GetRasterBand(1).SetColorTable(ct) assert ds.GetRasterBand(1).GetColorTable().GetCount() == 1 ds.GetRasterBand(1).SetColorTable(None) assert ds.GetRasterBand(1).GetColorTable() is None ############################################################################### # Check histogram support def test_vrtmisc_histogram(): tmpfile = "/vsimem/vrtmisc_histogram.vrt" ds = gdal.Translate(tmpfile, "data/byte.tif", format="VRT") ds.GetRasterBand(1).SetDefaultHistogram(1, 2, [3000000000, 4]) ds = None ds = gdal.Open(tmpfile) hist = ds.GetRasterBand(1).GetDefaultHistogram(force=0) ds = None assert hist == (1.0, 2.0, 2, [3000000000, 4]) gdal.Unlink(tmpfile) ############################################################################### # write SRS with unusual data axis to SRS axis mapping def test_vrtmisc_write_srs(): tmpfile = "/vsimem/test_vrtmisc_write_srs.vrt" ds = gdal.Translate(tmpfile, "data/byte.tif", format="VRT") sr = osr.SpatialReference() sr.SetAxisMappingStrategy(osr.OAMS_AUTHORITY_COMPLIANT) sr.ImportFromEPSG(4326) ds.SetSpatialRef(sr) assert ds.FlushCache() == gdal.CE_None ds = None ds = gdal.Open(tmpfile) assert ds.GetSpatialRef().GetDataAxisToSRSAxisMapping() == [1, 2] ds = None gdal.Unlink(tmpfile) ############################################################################### # complex scenario involving masks and implicit overviews def test_vrtmisc_mask_implicit_overviews(): with gdaltest.config_option("GDAL_TIFF_INTERNAL_MASK", "YES"): ds = gdal.Translate( "/vsimem/cog.tif", "data/stefan_full_rgba.tif", options="-outsize 2048 0 -b 1 -b 2 -b 3 -mask 4", ) ds.BuildOverviews("NEAR", [2, 4]) ds = None gdal.Translate("/vsimem/cog.vrt", "/vsimem/cog.tif") ds = gdal.Open("/vsimem/cog.vrt") assert ds.GetRasterBand(1).GetOverview(0).GetMaskFlags() == gdal.GMF_PER_DATASET assert ds.GetRasterBand(1).GetMaskBand().IsMaskBand() assert ds.GetRasterBand(1).GetOverview(0).GetMaskBand().IsMaskBand() ds = None gdal.Translate( "/vsimem/out.tif", "/vsimem/cog.vrt", options="-b mask -outsize 10% 0" ) gdal.GetDriverByName("GTiff").Delete("/vsimem/cog.tif") gdal.Unlink("/vsimem/cog.vrt") ds = gdal.Open("/vsimem/out.tif") histo = ds.GetRasterBand(1).GetDefaultHistogram()[3] # Check that there are only 0 and 255 in the histogram assert histo[0] + histo[255] == ds.RasterXSize * ds.RasterYSize, histo assert ds.GetRasterBand(1).Checksum() == 46885 ds = None gdal.Unlink("/vsimem/out.tif") ############################################################################### # Test setting block size def test_vrtmisc_blocksize(): filename = "/vsimem/test_vrtmisc_blocksize.vrt" vrt_ds = gdal.GetDriverByName("VRT").Create(filename, 50, 50, 0) options = ["subClass=VRTSourcedRasterBand", "blockXSize=32", "blockYSize=48"] vrt_ds.AddBand(gdal.GDT_Byte, options) vrt_ds = None vrt_ds = gdal.Open(filename) blockxsize, blockysize = vrt_ds.GetRasterBand(1).GetBlockSize() assert blockxsize == 32 assert blockysize == 48 vrt_ds = None gdal.Unlink(filename) ############################################################################### # Test setting block size through creation options def test_vrtmisc_blocksize_creation_options(): filename = "/vsimem/test_vrtmisc_blocksize_creation_options.vrt" vrt_ds = gdal.GetDriverByName("VRT").Create( filename, 50, 50, 1, options=["BLOCKXSIZE=32", "BLOCKYSIZE=48"] ) vrt_ds = None vrt_ds = gdal.Open(filename) blockxsize, blockysize = vrt_ds.GetRasterBand(1).GetBlockSize() assert blockxsize == 32 assert blockysize == 48 vrt_ds = None gdal.Unlink(filename) ############################################################################### # Test setting block size through creation options def test_vrtmisc_blocksize_gdal_translate_direct(): filename = "/vsimem/test_vrtmisc_blocksize_gdal_translate_direct.vrt" vrt_ds = gdal.Translate( filename, "data/byte.tif", creationOptions=["BLOCKXSIZE=32", "BLOCKYSIZE=48"] ) vrt_ds = None vrt_ds = gdal.Open(filename) blockxsize, blockysize = vrt_ds.GetRasterBand(1).GetBlockSize() assert blockxsize == 32 assert blockysize == 48 vrt_ds = None gdal.Unlink(filename) ############################################################################### # Test setting block size through creation options def test_vrtmisc_blocksize_gdal_translate_indirect(): filename = "/vsimem/test_vrtmisc_blocksize_gdal_translate_indirect.vrt" vrt_ds = gdal.Translate( filename, "data/byte.tif", metadataOptions=["FOO=BAR"], creationOptions=["BLOCKXSIZE=32", "BLOCKYSIZE=48"], ) vrt_ds = None vrt_ds = gdal.Open(filename) blockxsize, blockysize = vrt_ds.GetRasterBand(1).GetBlockSize() assert blockxsize == 32 assert blockysize == 48 vrt_ds = None gdal.Unlink(filename) ############################################################################### # Test support for coordinate epoch def test_vrtmisc_coordinate_epoch(): filename = "/vsimem/temp.vrt" gdal.Translate(filename, "data/byte.tif", options="-a_coord_epoch 2021.3") ds = gdal.Open(filename) srs = ds.GetSpatialRef() assert srs.GetCoordinateEpoch() == 2021.3 ds = None gdal.Unlink(filename) ############################################################################### # Test the relativeToVRT attribute of SourceFilename def test_vrtmisc_sourcefilename_all_relatives(): shutil.copy("data/byte.tif", "tmp") try: src_ds = gdal.Open(os.path.join("tmp", "byte.tif")) ds = gdal.GetDriverByName("VRT").CreateCopy("", src_ds) ds.SetDescription(os.path.join("tmp", "byte.vrt")) ds = None src_ds = None assert ( 'byte.tif<' in open("tmp/byte.vrt", "rt").read() ) finally: gdal.Unlink("tmp/byte.tif") gdal.Unlink("tmp/byte.vrt") ############################################################################### # Test the relativeToVRT attribute of SourceFilename def test_vrtmisc_sourcefilename_source_relative_dest_absolute(): shutil.copy("data/byte.tif", "tmp") try: src_ds = gdal.Open(os.path.join("tmp", "byte.tif")) ds = gdal.GetDriverByName("VRT").CreateCopy("", src_ds) path = os.path.join(os.getcwd(), "tmp", "byte.vrt") if sys.platform == "win32": path = path.replace("/", "\\") ds.SetDescription(path) ds = None src_ds = None assert ( 'byte.tif<' in open("tmp/byte.vrt", "rt").read() ) finally: gdal.Unlink("tmp/byte.tif") gdal.Unlink("tmp/byte.vrt") ############################################################################### # Test the relativeToVRT attribute of SourceFilename def test_vrtmisc_sourcefilename_source_absolute_dest_absolute(): shutil.copy("data/byte.tif", "tmp") try: src_ds = gdal.Open(os.path.join(os.getcwd(), "tmp", "byte.tif")) ds = gdal.GetDriverByName("VRT").CreateCopy("", src_ds) ds.SetDescription(os.path.join(os.getcwd(), "tmp", "byte.vrt")) ds = None src_ds = None assert ( 'byte.tif<' in open("tmp/byte.vrt", "rt").read() ) finally: gdal.Unlink("tmp/byte.tif") gdal.Unlink("tmp/byte.vrt") ############################################################################### # Test the relativeToVRT attribute of SourceFilename def test_vrtmisc_sourcefilename_source_absolute_dest_relative(): shutil.copy("data/byte.tif", "tmp") try: path = os.path.join(os.getcwd(), "tmp", "byte.tif") if sys.platform == "win32": path = path.replace("/", "\\") src_ds = gdal.Open(path) ds = gdal.GetDriverByName("VRT").CreateCopy("", src_ds) ds.SetDescription(os.path.join("tmp", "byte.vrt")) ds = None src_ds = None assert ( 'byte.tif<' in open("tmp/byte.vrt", "rt").read() ) finally: gdal.Unlink("tmp/byte.tif") gdal.Unlink("tmp/byte.vrt") ############################################################################### # Test Int64 nodata def test_vrtmisc_nodata_int64(): filename = "/vsimem/temp.vrt" ds = gdal.Translate( filename, "data/byte.tif", format="VRT", outputType=gdal.GDT_Int64 ) val = -(1 << 63) assert ds.GetRasterBand(1).SetNoDataValue(val) == gdal.CE_None assert ds.GetRasterBand(1).GetNoDataValue() == val ds = None ds = gdal.Open(filename) assert ds.GetRasterBand(1).GetNoDataValue() == val ds = None gdal.Unlink(filename) ############################################################################### # Test UInt64 nodata def test_vrtmisc_nodata_uint64(): filename = "/vsimem/temp.vrt" ds = gdal.Translate( filename, "data/byte.tif", format="VRT", outputType=gdal.GDT_UInt64 ) val = (1 << 64) - 1 assert ds.GetRasterBand(1).SetNoDataValue(val) == gdal.CE_None assert ds.GetRasterBand(1).GetNoDataValue() == val ds = None ds = gdal.Open(filename) assert ds.GetRasterBand(1).GetNoDataValue() == val ds = None gdal.Unlink(filename) ############################################################################### # Check IsMaskBand() on an alpha band def test_vrtmisc_alpha_ismaskband(): src_ds = gdal.GetDriverByName("MEM").Create("", 1, 1, 2) src_ds.GetRasterBand(2).SetColorInterpretation(gdal.GCI_AlphaBand) ds = gdal.GetDriverByName("VRT").CreateCopy("", src_ds) assert not ds.GetRasterBand(1).IsMaskBand() assert ds.GetRasterBand(2).IsMaskBand() ############################################################################### # Check that serialization of a ComplexSource with NODATA (generally) doesn't # require opening the source band. def test_vrtmisc_serialize_complexsource_with_NODATA(): tmpfilename = "/vsimem/test_vrtmisc_serialize_complexsource_with_NODATA.vrt" gdal.FileFromMemBuffer( tmpfilename, """ 0 i_do_not_exist.tif 1 1 """, ) ds = gdal.Open(tmpfilename) ds.GetRasterBand(1).SetDescription("foo") gdal.ErrorReset() ds = None assert gdal.GetLastErrorMsg() == "" fp = gdal.VSIFOpenL(tmpfilename, "rb") content = gdal.VSIFReadL(1, 10000, fp).decode("latin1") gdal.VSIFCloseL(fp) gdal.Unlink(tmpfilename) print(content) assert "1" in content ############################################################################### # Test bugfix for https://github.com/OSGeo/gdal/issues/7486 def test_vrtmisc_nodata_float32(): tif_filename = "/vsimem/test_vrtmisc_nodata_float32.tif" ds = gdal.GetDriverByName("GTiff").Create(tif_filename, 1, 1, 1, gdal.GDT_Float32) nodata = -0.1 ds.GetRasterBand(1).SetNoDataValue(nodata) ds.GetRasterBand(1).Fill(nodata) ds = None # When re-opening the TIF file, the -0.1 double value will be exposed # with the float32 precision (~ -0.10000000149011612) vrt_filename = "/vsimem/test_vrtmisc_nodata_float32.vrt" ds = gdal.Translate(vrt_filename, tif_filename) nodata_vrt = ds.GetRasterBand(1).GetNoDataValue() assert nodata_vrt == struct.unpack("f", struct.pack("f", nodata))[0] ds = None # Check that this is still the case after above serialization to .vrt # and re-opening. That is check that we serialize the rounded value with # full double precision (%.18g) ds = gdal.Open(vrt_filename) nodata_vrt = ds.GetRasterBand(1).GetNoDataValue() assert nodata_vrt == struct.unpack("f", struct.pack("f", nodata))[0] ds = None gdal.Unlink(tif_filename) gdal.Unlink(vrt_filename)