[gdal-dev] Does GDAL API for Python use <Scale> and <Offset> in VRT?

Even Rouault even.rouault at mines-paris.org
Wed Feb 8 06:01:41 EST 2012


Selon Anton Korosov <anton.korosov at nersc.no>:

> Hello everybody!
>
> I'm trying to read values from a VRT file in Python.  When I modify
> <Scale> or <Offset> it seem not to make any effect.
>
> I have a simple VRT of a Landsat7 file (below). It has neither scale nor
> offset. I launch the given script twice:
> 1st time to test the value of scale and the value from a band and to
> change the value of scale.
> 2nd time I test value of scale and value from a band. The output shows
> that the scale has changed, but the value is the same. I would expect
> that the second time the value should be 70.
>
> Does it mean that ReadAsArray() doesn't interpret <Scale> and <Offset>
> tags? Ot do I do (expect) something wrong?

Yes, ReadAsArray() doesn't interpret <Scale> or <Offset>. They are
interpretation hints for the user, but GDAL won't alter the values read from the
dataset.

Alternatively, you can use the -unscale option of gdal_translate so that scale
and offset are taken into account.

>
> Thank you!
> Anton
>
>
>
> Script:
> l7 = gdal.Open('L71211005_00520100824_B10.vrt')
> print l7.GetRasterBand(1).GetScale()
> print l7.GetRasterBand(1).ReadAsArray()[3000, 3000]
> l7.GetRasterBand(1).SetScale(2)
> l7.FlushCache()
> l7 = None
>
> Output:
>  > python test.py
> 1.0
> 35
>  > python test.py
> 2.0
> 35
>
> VRT-file before test:
> <VRTDataset rasterXSize="9211" rasterYSize="9021">
>    <SRS>PROJCS["WGS 84 / UTM zone 33N",GEOGCS["WGS
> 84",DATUM["WGS_1984",SPHEROID["WGS
>
84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",15],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","32633"]]</SRS>
>    <GeoTransform>  3.6178500000000000e+05,  3.0000000000000000e+01,
> 0.0000000000000000e+00,  8.6874150000000000e+06,
> 0.0000000000000000e+00, -3.0000000000000000e+01</GeoTransform>
>    <Metadata>
>      <MDI key="AREA_OR_POINT">Point</MDI>
>    </Metadata>
>    <VRTRasterBand dataType="Byte" band="1">
>      <Metadata />
>      <ColorInterp>Gray</ColorInterp>
>      <SimpleSource>
>        <SourceFilename
> relativeToVRT="1">L71211005_00520100824_B10.TIF</SourceFilename>
>        <SourceBand>1</SourceBand>
>        <SourceProperties RasterXSize="9211" RasterYSize="9021"
> DataType="Byte" BlockXSize="9211" BlockYSize="1" />
>        <SrcRect xOff="0" yOff="0" xSize="9211" ySize="9021" />
>        <DstRect xOff="0" yOff="0" xSize="9211" ySize="9021" />
>      </SimpleSource>
>    </VRTRasterBand>
> </VRTDataset>
>
> VRT-file after test:
> <VRTDataset rasterXSize="9211" rasterYSize="9021">
>    <SRS>PROJCS["WGS 84 / UTM zone 33N",GEOGCS["WGS
> 84",DATUM["WGS_1984",SPHEROID["WGS
>
84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",15],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","32633"]]</SRS>
>    <GeoTransform>  3.6178500000000000e+05,  3.0000000000000000e+01,
> 0.0000000000000000e+00,  8.6874150000000000e+06,
> 0.0000000000000000e+00, -3.0000000000000000e+01</GeoTransform>
>    <Metadata>
>      <MDI key="AREA_OR_POINT">Point</MDI>
>    </Metadata>
>    <VRTRasterBand dataType="Byte" band="1">
>      <Metadata />
>      <Scale>2</Scale>
>      <ColorInterp>Gray</ColorInterp>
>      <SimpleSource>
>        <SourceFilename
> relativeToVRT="1">L71211005_00520100824_B10.TIF</SourceFilename>
>        <SourceBand>1</SourceBand>
>        <SourceProperties RasterXSize="9211" RasterYSize="9021"
> DataType="Byte" BlockXSize="9211" BlockYSize="1" />
>        <SrcRect xOff="0" yOff="0" xSize="9211" ySize="9021" />
>        <DstRect xOff="0" yOff="0" xSize="9211" ySize="9021" />
>      </SimpleSource>
>    </VRTRasterBand>
> </VRTDataset>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>




More information about the gdal-dev mailing list