[gdal-dev] Numerical instability with thin plate spline transform

[Jan Hartmann]

On 09/09/2011 08:11 AM, Even Rouault wrote:
> You might try latest GDAL trunk, in particular 
> http://trac.osgeo.org/gdal/changeset/22876 . It offers the option to 
> use libarmadillo to speed-up matrix inversion in thinplatespline.cpp 
> (speed-up when TPS are in the thousands), but perhaps as a side 
> effect, you'd get also more numerical stability.

I tried to compile the SVN version a few times the last couple of days, 
but it always broke down: Can someone do a "make clean" in their sandbox 
and repair the error?

  g++ -fPIC -L/home/a907hart/local/lib -L/usr/lib  gdaldem.o 
commonutils.o  -L/home/a907hart/src/gdal -lgdal  -larmadillo 
-L/home/a907hart/local/lib -lgeos_c -L/home/a907hart/local/pgsql/lib 
-lpq -lpthread -lm -lrt -ldl  -L/home/a907hart/local/lib   
-lcurl            -o gdaldem
gdaldem.o:(.data.rel.ro._ZTV21GDALGeneric3x3Dataset[vtable for 
GDALGeneric3x3Dataset]+0x30): undefined reference to 
`GDALDataset::CloseDependentDatasets()'
gdaldem.o:(.data.rel.ro._ZTV22GDALColorReliefDataset[vtable for 
GDALColorReliefDataset]+0x30): undefined reference to 
`GDALDataset::CloseDependentDatasets()'


Jan

>> On Thu, Sep 8, 2011 at 5:23 AM, Jan Hartmann<j.l.h.hartmann at uva.nl>  wrote:
>>> Not sure whether this can be considered a bug, so I give it for what it
>>> is worth.
>>>
>>> I'm doing thin plate spline transformation from one set of projected
>>> coordinates to another. Both sets have values between -600000 and 600000
>>> (meters). A typical set of gcps looks like:
>>>
>>>   -gcp 62402 -74383 181915 315371 \
>>>   -gcp -100169 -24213 19685 366661 \
>>>   -gcp 118323 233822 239994 623190 \
>>> ...
>>>
>>> When transforming the the first two columns of this list with
>>> gdaltransform -tps, using the same gcp-list, the results should be
>>> exactly equal to the last two columns:
>>>
>>>   (from gdal_tps.cpp:)
>>>   * The thin plate spline transformer produces exact transformation
>>>   * at all control points and smoothly varying transformations between
>>>   * control points with greatest influence from local control points.
>>>   * It is suitable for for many applications not well modelled by
>>> polynomial * transformations.
>>>   *
>>>
>>>
>>> However, they aren't. With a few control points the differences are only
>>> in millimeters, but with a few hundred gcps the differences become
>>> meters, and above thousand points the error can get to fifty meters. The
>>> problem gets worse when some control points are very near to other ones.
>>> I was happy to discover that dividing all numbers by 1000 solves the
>>> problem, but there certainly is a numerical instability in the
>>> algorithm.
>>>
>>> Of course, thin plate spline transformations normally use scan
>>> coordinates as input, and these are almost always small numbers. Even
>>> so, I can imagine problems with very large rasters and the rubber
>>> sheeting applications I am working with.
>>>
>>> Jan
>>>
>>>
>>>
>>> _______________________________________________
>>> gdal-dev mailing list
>>> gdal-dev at lists.osgeo.org
>>> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>> _______________________________________________
>> gdal-dev mailing list
>> gdal-dev at lists.osgeo.org
>> http://lists.osgeo.org/mailman/listinfo/gdal-dev
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.osgeo.org/pipermail/gdal-dev/attachments/20110917/58a68ce7/attachment.html