I found the culprit. If remove this section from each band definition in
VRT file then everything works fine:
<ComplexSource>
<SourceFilename relativeToVRT="1" shared="0">dummy.tif</SourceFilename>
<SourceBand>3</SourceBand>
<SourceProperties BlockXSize="128" BlockYSize="128"
RasterXSize="40961" RasterYSize="139265" dataType="Byte" />
<SrcRect xOff="0" xSize="1" yOff="0" ySize="1" />
<DstRect xOff="0" xSize="1" yOff="0" ySize="1" />
<ScaleRatio>0</ScaleRatio>
<ScaleOffset>0.0</ScaleOffset></ComplexSource
Rasterio's code that builds this code was altered 2 years ago after
conversation in OSGeo/gdal#1135 <https://github.com/OSGeo/gdal/issues/1135>.
But it has this undesirable consequence I'm experiencing. Any ideas how to
overcome that?
On Thu, Aug 27, 2020 at 12:41 PM Denis Rykov <[email protected]> wrote:
> I was able to reproduce this issue with pure GDAL. When you read data with
> boundless=True in rasterio it creates an intermediate VRT file. This is the
> example of file that being created in my case:
>
> <?xml version="1.0" encoding="UTF-8"?>
> <VRTDataset rasterXSize="40961" rasterYSize="139265">
> <SRS>PROJCS["WGS 84 / UTM zone 47N",GEOGCS["WGS
> 84",DATUM["WGS_1984",SPHEROID["WGS
> 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",99],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32647"]]</SRS>
> <GeoTransform>443070.0,1.0,0.0,4366312.0,0.0,-1.0</GeoTransform>
> <VRTRasterBand band="1" dataType="Byte">
> <NoDataValue>0.0</NoDataValue>
> <ColorInterp>Red</ColorInterp>
> <ComplexSource>
> <SourceFilename relativeToVRT="1"
> shared="0">dummy.tif</SourceFilename>
> <SourceBand>1</SourceBand>
> <SourceProperties BlockXSize="128"
> BlockYSize="128" RasterXSize="40961" RasterYSize="139265" dataType="Byte" />
> <SrcRect xOff="0" xSize="1" yOff="0" ySize="1"
> />
> <DstRect xOff="0" xSize="1" yOff="0" ySize="1"
> />
> <ScaleRatio>0</ScaleRatio>
> <ScaleOffset>0.0</ScaleOffset>
> </ComplexSource>
> <ComplexSource>
> <SourceFilename relativeToVRT="0"
> shared="0">/vsicurl/https://*.vrt</SourceFilename>
> <SourceBand>1</SourceBand>
> <SourceProperties BlockXSize="128"
> BlockYSize="128" RasterXSize="40961" RasterYSize="139265" dataType="Byte" />
> <SrcRect xOff="0" xSize="40960" yOff="0"
> ySize="139264" />
> <DstRect xOff="-6961.0" xSize="40960.0"
> yOff="-105176.0" ySize="139264.0" />
> <NODATA>0.0</NODATA>
> <OpenOptions /> </ComplexSource>
> </VRTRasterBand>
> <VRTRasterBand band="2" dataType="Byte">
> <NoDataValue>0.0</NoDataValue>
> <ColorInterp>Green</ColorInterp>
> <ComplexSource>
> <SourceFilename relativeToVRT="1"
> shared="0">dummy.tif</SourceFilename>
> <SourceBand>2</SourceBand>
> <SourceProperties BlockXSize="128"
> BlockYSize="128" RasterXSize="40961" RasterYSize="139265" dataType="Byte" />
> <SrcRect xOff="0" xSize="1" yOff="0" ySize="1"
> />
> <DstRect xOff="0" xSize="1" yOff="0" ySize="1"
> />
> <ScaleRatio>0</ScaleRatio>
> <ScaleOffset>0.0</ScaleOffset>
> </ComplexSource>
> <ComplexSource>
> <SourceFilename relativeToVRT="0"
> shared="0">/vsicurl/https://*.vrt</SourceFilename>
> <SourceBand>2</SourceBand>
> <SourceProperties BlockXSize="128"
> BlockYSize="128" RasterXSize="40961" RasterYSize="139265" dataType="Byte" />
> <SrcRect xOff="0" xSize="40960" yOff="0"
> ySize="139264" />
> <DstRect xOff="-6961.0" xSize="40960.0"
> yOff="-105176.0" ySize="139264.0" />
> <NODATA>0.0</NODATA>
> <OpenOptions /> </ComplexSource>
> </VRTRasterBand>
> <VRTRasterBand band="3" dataType="Byte">
> <NoDataValue>0.0</NoDataValue>
> <ColorInterp>Blue</ColorInterp>
> <ComplexSource>
> <SourceFilename relativeToVRT="1"
> shared="0">dummy.tif</SourceFilename>
> <SourceBand>3</SourceBand>
> <SourceProperties BlockXSize="128"
> BlockYSize="128" RasterXSize="40961" RasterYSize="139265" dataType="Byte" />
> <SrcRect xOff="0" xSize="1" yOff="0" ySize="1"
> />
> <DstRect xOff="0" xSize="1" yOff="0" ySize="1"
> />
> <ScaleRatio>0</ScaleRatio>
> <ScaleOffset>0.0</ScaleOffset>
> </ComplexSource>
> <ComplexSource>
> <SourceFilename relativeToVRT="0"
> shared="0">/vsicurl/https://*.vrt</SourceFilename>
> <SourceBand>3</SourceBand>
> <SourceProperties BlockXSize="128"
> BlockYSize="128" RasterXSize="40961" RasterYSize="139265" dataType="Byte" />
> <SrcRect xOff="0" xSize="40960" yOff="0"
> ySize="139264" />
> <DstRect xOff="-6961.0" xSize="40960.0"
> yOff="-105176.0" ySize="139264.0" />
> <NODATA>0.0</NODATA>
> <OpenOptions /> </ComplexSource>
> </VRTRasterBand>
> <MaskBand>
> <VRTRasterBand dataType="Byte">
> <SimpleSource>
> <SourceFilename relativeToVRT="0"
> shared="0">/vsicurl/https://*.vrt</SourceFilename>
> <SourceBand>mask,1</SourceBand>
> <SourceProperties BlockXSize="128"
> BlockYSize="128" RasterXSize="40961" RasterYSize="139265" dataType="Byte" />
> <SrcRect xOff="0" xSize="40960"
> yOff="0" ySize="139264" />
> <DstRect xOff="-6961.0" xSize="40960"
> yOff="-105176.0" ySize="139264" /> </SimpleSource>
> </VRTRasterBand>
> </MaskBand>
> </VRTDataset>
>
> If I read data from this file with gdal using
> resample_alg=gdalconst.GRIORA_NearestNeighbour then GDAL takes into account
> *.vrt.ovr file and sends very few HTTP requests to the server (~30):
>
> ds = gdal.OpenEx("/tmp/rasterio.vrt")
> image = ds.ReadAsArray(xoff=0, yoff=0, xsize=5671, ysize=5648, buf_xsize=383,
> buf_ysize=385, resample_alg=gdalconst.GRIORA_NearestNeighbour
>
> If I do the same but using resample_alg=gdalconst.GRIORA_Cubic then GDAL
> sends a huge amount of requests to the server (~1k) because overviews are
> not used:
>
> ds = gdal.OpenEx("/tmp/rasterio.vrt")
> image = ds.ReadAsArray(xoff=0, yoff=0, xsize=5671, ysize=5648, buf_xsize=383,
> buf_ysize=385, resample_alg=gdalconst.GRIORA_Cubic
>
> Is it expected or might there be something wrong with that VRT file?
> Thanks in advance for any help.
>
> On Wed, Aug 26, 2020 at 6:18 PM Denis Rykov <[email protected]> wrote:
>
>> I have remote *.vrt raster and *.vrt.ovr accessible through HTTP. When I
>> run the following script with rasterio:
>>
>> with rasterio.open("http://*.vrt"";) as src:
>> image = src.read(indexes=[1, 2, 3], **{
>> "window": Window(col_off=6961, row_off=105176, width=5671,
>> height=5648),
>> "resampling": Resampling.cubic,
>> "boundless": True,
>> "out_shape": (3, 383, 385),
>> "masked": True
>> })
>>
>> depending on "resampling" algorithm GDAL sends different amounts of
>> requests to the server. In the case of "cubic" it doesn't take into account
>> overviews and sends requests directly to *.tif files (900 in my case). In
>> case of "nearest" everything is ok (only 60 requests, *.vrt.ovr is taken
>> into account).
>>
>> Does GDAL check the resampling algorithm of overviews and in case it
>> differs from the option specified in read() method they are bypassed or it
>> works differently?
>>
>
_______________________________________________
gdal-dev mailing list
[email protected]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
