gdalwarp: suboptimal clipping

[hamish]

Hi,

there are clipping issues with both -order 3 and -tps methods of gdalwarp. (with latest GDAL SVN build; but it has been like this for years)

the attached images demonstrate.

-tps: for heavily "bent" images the edges of the image are slightly clipped off. In the example image this is most obvious in the bottom left, but is also there in the middle-right.

-order 3: regular part of image is fine, but a reflection of some sort makes the image much larger than it should be and leaves behind an artifact.

the third image is the input paper roll scan. GCPs were applied using gdal_translate.

The geotiff input image is available for download here: (5mb)

​http://bambi.otago.ac.nz/hamish/grass/gdal/18_R_gcp.tif.bz2

thanks, Hamish

[hamish]

3rd order gdalwarp showing reflection artifact

[hamish]

heavily bent tps gdalwarp output with corners slightly clipped

[hamish]

reduced version of input image

[hamish]

Hi,

slightly OT to this bug, but some positive feedback: the latest TPS code seems to produce a much more stable image than the old. Compare the warped image attached to this ticket with the same produced from gdalwarp circa version 1.2.6:

​http://bambi.otago.ac.nz/hamish/grass/gdal/18_R_geo_Ord_tps.jpg

cheers, Hamish

[warmerdam]

I'm reviewing this. In particular I'm wondering if the "bounds computation" properly goes to the corner and if the failure to select the right output bounds implies that the TPS transformation is not invertable at points other than the GCPs (which would explain the problem).

[warmerdam]

I have confirmed that the thin plate spline is not invertable at locations other than the gcps. The extent for the output file is computed by forwarding transforming through the tps transform. But actually filling of it is accomplished by backwards transforming to sample the source image. The non-invertability means that the approach cannot accurately compute the ideal bounds. This can impact other non-invertable transforms.

pixel to geo:

warmerda@gdal64[230]% gdaltransform -tps 18_R_gcp.tif 
1200 4900
2262335.24467225 5410436.31834668 0

geo to pixel:

warmerda@gdal64[4]% gdaltransform -i -tps 18_R_gcp.tif 
2262335.24467225 5410436.31834668 0
1177.4299815143 4872.53662873868 0

In theory we could use some sort of iterative approach to get an essentially invertable transformation; however that could be much more computational expensive and it isn't clear that the results are more accurate - just more well behaved and predictable.

It would be interesting if the whole warper had a standard way of forcing an iteration based inverted transformation for transformations that are not otherwise properly invertable.

I'm going to treat this as an enhancement request, and have no plans to act on it. But I claim it is not technically a bug. It is open for anyone else who would like to try and address it.

[hamish]

I have confirmed that the thin plate spline is not invertable at locations other than the gcps.

Re-checking the GCPs I see that the lowest ones are in fact not at the bottom corners, just at the last cross-line. See the GCP marks over-plotted here:

​http://bambi.otago.ac.nz/hamish/grass/gdal/sidescan_warp.html

So that makes sense.

I wonder as a practical solution to minimize the effects if tps output images could have their edges expanded up/down/left/right by say 5% of the width/height, rather than tight clipping? ie round outwards a little bit.

any comment on the "-order 3" reflection artifact?

thanks, Hamish

[Even Rouault]

I've made a try with the issue and did an implementation based on Frank's suggestion. It is just a proof of concept, as far as the quality of coding is concerned, and for the speed too... It gives very nice results with the TPS algorithm. The image is not cropped at all. The drawback is that it took 84 minutes to do the wrapping...

Anyway, I'm attaching the patch to gdal_tps.cpp if someone wants to experiment...

As far as order 2 and 3 are concerned, I've never managed to make them work. I often get reflections too.

[Even Rouault]

Proof-of-concept not-production-ready patch to improve the computation of the invert of the TPS interpolation

[hamish]

Replying to rouault:

As far as order 2 and 3 are concerned, I've never managed to make them work. I often get reflections too.

FWIW, "-order 2" has always worked well for me in GRASS (we use identical [ie same parent/origin code] for imagery/i.rectify/crs.c for our i.rectify module as gdal/alg/gdal_crs.c does for gdalwarp), and order 3 has been a problem there too- see:

• "bad target region calculations with i.rectify order=3" ​http://intevation.de/rt/webrt?serial_num=3166

Hamish

[Even Rouault]

This ticket has been automatically closed because Trac is no longer used for GDAL bug tracking, since the project has migrated to GitHub. If you believe this ticket is still valid, you may file it to ​https://github.com/OSGeo/gdal/issues if it is not already reported there.