[gdal-dev] "wrong" overview selected when warping from ortho to webmercator

[thomas bonfort]

I'm trying to warp/extract a webmercator tile from a large +proj=ortho
image. The issue is that the further I move away from the center of the
image the quality of the output result degrades. This seems to be due to
which overview is being selected from the source image.

input image created with -t_srs "+proj=ortho +lon_0=99.8 +lat_0=61.52" .
gdalinfo in-ortho.tif:

Size is 270322, 145529
Coordinate System is: [snip]
Origin = (-4730365.833427881821990,2876021.392761887982488)
Pixel Size = (35.000000000000000,-35.000000000000000)
Corner Coordinates:
Upper Left  (-4730365.833, 2876021.393) (  2d20'24.47"W, 40d39'26.68"N)
Lower Left  (-4730365.833,-2217493.607) ( 47d47' 8.30"E, 19d46'55.15"N)
Upper Right ( 4730904.167, 2876021.393) (158d 3'24.38"W, 40d38'56.55"N)
Lower Right ( 4730904.167,-2217493.607) (151d49'12.07"E, 19d46'34.13"N)
Center      (     269.167,  329263.893) ( 99d48'20.20"E, 64d28'44.92"N)
Band 1 Block=512x512 Type=Byte, ColorInterp=Red
  Overviews: 135161x72765, 67581x36383, 33791x18192, 16896x9096, 8448x4548,
4224x2274, 2112x1137, 1056x569

extracting webmercator tile 7,96,12 with "gdalwarp -t_srs epsg:3857 -ts 256
256 -te 10018754.171394620091 15967389.460660178214 10331840.239250702783
16280475.528516260907 -r cubic -to
AREA_OF_INTEREST=90.000000,80.647035,92.812500,81.09321 in-ortho.tif
bad.tif" gives:

GDAL: GDALOpen(in-ortho.tif, this=0x55597052f230) succeeds as GTiff.
GDAL: Using GTiff driver
Using band 4 of source image as alpha.
Creating output file that is 256P x 256L.
GDAL: GDALDriver::Create(GTiff,bad.tif,256,256,4,Byte,0x5559705570a0)
Processing in-ortho.tif [1/1] : 0WARP: Copying metadata from first source
to destination dataset
GTiff: ScanDirectories()
GTiff: Opened 135161x72765 overview.
GTiff: Opened 67581x36383 overview.
GTiff: Opened 33791x18192 overview.
GTiff: Opened 16896x9096 overview.
GTiff: Opened 8448x4548 overview.
GTiff: Opened 4224x2274 overview.
GTiff: Opened 2112x1137 overview.
GTiff: Opened 1056x569 overview.
WARP: Selecting overview level 3 for in-ortho.tif
GDAL: GDALOverviewDataset(in-ortho.tif, this=0x55597051dce0) creation.
GDALWARP: Defining SKIP_NOSOURCE=YES
WARP: SetAlphaMax: AlphaMax not set.
WARP: SetAlphaMax: AlphaMax not set.
GDAL: GDAL_CACHEMAX = 796 MB
GDAL: GDALWarpKernel()::GWKCubicNoMasksOrDstDensityOnlyByte()
Src=8130,1295,105x99 Dst=0,0,256x256
...10...20...30...40...50...60...70...80...90...100 - done.
GDAL: Flushing dirty blocks:
0...10...20...30...40...50...60...70...80...90...100 - done.
GDAL: GDALClose(bad.tif, this=0x55597056bd20)
GDAL: GDALClose(in-ortho.tif, this=0x55597052f230)

Where we can see that the warper is using a window of 105x99 in the 3rd
overview. The output image is of course blurry.

Selecting a tile closer to the ortho center (7,104,38) /opt/bin/gdalwarp
-t_srs epsg:3857 -ts 256 256 -te 12523442.714243277907 7827151.696402047761
12836528.782099360600 8140237.764258129522 -r cubic -to
AREA_OF_INTEREST=112.500000,57.326521,115.312500,58.813742 in-ortho.tif
good.tif gives:

GDAL: GDALOpen(in-ortho.tif, this=0x562469aa9230) succeeds as GTiff.
GDAL: Using GTiff driver
Using band 4 of source image as alpha.
Creating output file that is 256P x 256L.
GDAL: GDALDriver::Create(GTiff,good.tif,256,256,4,Byte,0x562469ac2f90)
Processing in-ortho.tif [1/1] : 0WARP: Copying metadata from first source
to destination dataset
GTiff: ScanDirectories()
GTiff: Opened 135161x72765 overview.
GTiff: Opened 67581x36383 overview.
GTiff: Opened 33791x18192 overview.
GTiff: Opened 16896x9096 overview.
GTiff: Opened 8448x4548 overview.
GTiff: Opened 4224x2274 overview.
GTiff: Opened 2112x1137 overview.
GTiff: Opened 1056x569 overview.
WARP: Selecting overview level 3 for in-ortho.tif
GDAL: GDALOverviewDataset(in-ortho.tif, this=0x562469a97ce0) creation.
GDALWARP: Defining SKIP_NOSOURCE=YES
WARP: SetAlphaMax: AlphaMax not set.
WARP: SetAlphaMax: AlphaMax not set.
GDAL: GDAL_CACHEMAX = 796 MB
GDAL: GDALWarpKernel()::GWKCubicNoMasksOrDstDensityOnlyByte()
Src=9741,5481,355x359 Dst=0,0,256x256
...10...20...30...40...50...60...70...80...90...100 - done.
GDAL: Flushing dirty blocks:
0...10...20...30...40...50...60...70...80...90...100 - done.
GDAL: GDALClose(good.tif, this=0x562469ac4680)
GDAL: GDALClose(in-ortho.tif, this=0x562469aa9230)

Where we can see that we are still using the 3rd overview but with a source
window of 355x359 pixels which of course gives a much better result.

* I know that in the first case I can use "-ovr AUTO-2" to force the
selection of the "correct" overview, but that's going to be inefficient
when extracting close to the center.
* I added the -to AREA_OF_INTEREST= options while trying to debug the issue
but it makes no difference with or without it
* The issue is not specific to a +proj=ortho input projection, I'm having
similar issues with an +proj=laea projection

Is there an easy way to fix this? Or could there be an initial first step
that programatically sets the -ovr AUTO-N parameter depending on the size
of the source window that is going to be used?

best regards,
Thomas
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20200717/48cb9d8a/attachment.html>