[gdal-dev] Forcing gdal translate to only rotate/scale

[Even Rouault]

Le samedi 20 septembre 2014 18:56:11, Jesse McGraw a écrit :
> Is it possible to force gdal_translate to only rotate/scale an image
> when given a list of more than two GCPs?
> 
> I'm manually georeferencing maps which are all North-up and I'd like
> to supply as many GCPs as possible to average out my errors but if I
> provide more than two it will potentially skew the image.

Jesse,

GDALGCPsToGeoTransform() would indeed return an affine transform matrix, 
computated with least-squares, but with potential shearing.

If you want only rotation and scaling (I assume uniform scaling), that would 
mean to solve

Xi = X0 + S * cos(theta) * xi - S * sin(theta) * yi
Yi = Y0  + S * sin(theta) * xi + S * cos(theta) * yi

where the unknown variables are X0, Y0, S and theta.

Unfortunately, I think this would be difficult to analytically solve it because 
of the non-linear terms in cos(theta), sin(theta).

Perhaps a pragmatic way would be to use GDALGCPsToGeoTransform()  that would 
return (X0, A, B, Y0, C, D) such that

Xi = X0 + A * xi + B * yi
Yi = Y0 + C * xi + D * yi

And then you can try the following approximation :
S ~= sqrt(A* D - B * C)
theta ~= 0.5 (atan2(-B,A) + atan2(C,D))

(that would be exact if the transform is only scaling and rotation)

You can then recompute X0 and Y0 as the average of
Xi - (  S * cos(theta) * xi - S * sin(theta) * yi ) and
Yi - (  S * sin(theta) * xi + S * cos(theta) * yi )

This will be not the mathematical optimal solution, but hopefully it would not 
be very far from the optimal provided that the amount of shearing is low. If 
you want something more exact, you can use this solution as the initial 
solution for an interative method as described in 
http://en.wikipedia.org/wiki/Non-linear_least_squares

Note: I didn't try, so there might be errors in the formulas. And they are 
based with indexing pixels from the bottom left corner of the image, so some 
sign tweaking will be necessary.

Even

> 
> Thanks!
> -Jesse
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/gdal-dev

-- 
Spatialys - Geospatial professional services
http://www.spatialys.com