[SciPy-User] helmert implementation (7 parameters - geometric transformation)

Charles R Harris charlesr.harris@gmail....
Mon Mar 1 22:08:20 CST 2010


On Mon, Mar 1, 2010 at 6:17 PM, Massimo Di Stefano
<massimodisasha@yahoo.it>wrote:

> Hi All,
>
> i'm tring to implement a function to perform a geometric trasformation
> between a set of points in 2 different reference systems
>
> it is a 7 parameters trasformation (helmert) , based on last square method
> it use as input 2 set of x,y,z coordinates in 2 different reference systems
> and give as output 7 parameters needed to performs transformation
> from one reference system to the other one.
>
> http://en.wikipedia.org/wiki/Helmert_transformation
>
> this function is used in geodesy to "reproject" coordinates between 2
> different Datum
>
> i wrote it as :
>
> import numpy as np
> from scipy import linalg
>
> def helmert(p1, p2):
>    L = np.loadtxt(p1)
>    A = np.zeros((3*L.shape[0],7),float)
>    A[ ::3, 0] = 1.0
>    A[1::3, 1] = 1.0
>    A[2::3, 2] = 1.0
>    A[ ::3, 3] = L[:,0]
>    A[1::3, 3] = L[:,1]
>    A[2::3, 3] = L[:,2]
>    A[1::3, 4] = L[:,2]
>    A[2::3, 4] = -L[:,1]
>    A[ ::3, 5] = -L[:,2]
>    A[2::3, 5] = L[:,0]
>    A[ ::3, 6] = L[:,1]
>    A[1::3, 6] = -L[:,0]
>    G = np.loadtxt(p2)
>    Y = np.zeros((3*G.shape[0],1),float)
>    Y[ ::3, 0] = G[:,0] - L[:,0]
>    Y[1::3, 0] = G[:,1] - L[:,1]
>    Y[2::3, 0] = G[:,2] - L[:,2]
>    N = np.dot(A.T.conj(), A)
>    T = np.dot(A.T.conj(), Y)
>    C = np.dot(linalg.inv(N), T)
>    print C
>
>
> from a pdf i find online
> it has numerical examples
> and i get different results :'(
>
>
> p1 :
>
> 4402553.334 727053.937 4542823.474
> 4399375.347 703845.876 4549215.105
> 4412911.150 701094.435 4536518.139
>
>
> p2 :
>
> 4402553.569 727053.737 4542823.332
> 4399375.518 703845.639 4549214.880
> 4412911.336 701094.214 4536517.955
>
>
> pdf results :
>
> x0    -9.256 m
> y0    -23.701 m
> z0    16.792 m
> Rx    -0.0001990982 rad
> Ry    0.0001778762 rad
> Rz    0.00015 rad
> μ      0.00000046
>
>
> my results :
>
> [[ -9.91564629e+00]
>  [ -2.36172028e+01]
>  [  1.57283853e+01]
>  [ -2.68063331e-07]
>  [  3.04330048e-06]
>  [ -2.76148185e-06]
>  [ -2.31598150e-06]]
>
>
>
> anyone here can give it a try ?
> the results seems wrong but
> i haveb't yet find a solution.
>
> thanks to ALL!
>
> regards,
>
> Massimo.
>
>
>
>
>
> _______________________________________________
> SciPy-User mailing list
> SciPy-User@scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/scipy-user/attachments/20100301/4ef0d873/attachment-0001.html 


More information about the SciPy-User mailing list