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

Massimo Di Stefano massimodisasha@yahoo...
Tue Mar 2 17:29:42 CST 2010

```Hi,

from here :

http://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.cholesky.html

seems it is applicable only on a square matrix.
is that true ?

http://geomatica.como.polimi.it/corsi/misure_geodetiche/MG_03SistRif2.pdf

scroll to the page : 15

i can see it isn't a common system :  Y = AX
but it is :  Y = AX + B + V

the points are :

p1:
4402553.569 727053.737 4542823.332
4399375.518 703845.639 4549214.880
4412911.336 701094.214 4536517.955

p2 :
4402553.334 727053.937 4542823.474
4399375.347 703845.876 4549215.105
4412911.150 701094.435 4536518.139

tring to rewrite the code i did :

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[1::3, 3] = L[:,2]
A[2::3, 3] = -L[:,1]
A[ ::3, 4] = -L[:,2]
A[2::3, 4] = L[:,0]
A[ ::3, 5] = L[:,1]
A[1::3, 5] = -L[:,0]
A[ ::3, 6] = L[:,0]
A[1::3, 6] = L[:,1]
A[2::3, 6] = L[:,2]

V = np.zeros((3*G.shape[0],1),float)
V[ ::3, 0] = G[:,0] - L[:,0]
V[1::3, 0] = G[:,1] - L[:,1]
V[2::3, 0] = G[:,2] - L[:,2]

i used the code in the previouse mail :

N = np.dot(A.T.conj(), A)
T = np.dot(A.T.conj(), Y)
C = np.dot(linalg.inv(N), T)

to solve a system like :

Y = AX

how can i change the code to solve :

Y = AX + B + V

the results in the dpf are :

x0     -9.256 m
y0     -23.701 m
z0     16.792 m
μ    0.00000046

any suggestion, also on how to apply a decomposition can give me a great help!

Il giorno 02/mar/2010, alle ore 19.00, scipy-user-request@scipy.org ha scritto:

> Hi,
>
> I see your coordinates are very big numbers. These can lead to an
> ill-conditioned (AtA) matrix which can lead to a wrong inverse. Try
> using cholesky decomposition instead of using linalg.inv. The method
> is outlined here http://en.wikipedia.org/wiki/Cholesky_decomposition
> and numpy has routines for doing just this. The other option is to
> subtract the centre-of-mass from your coordinates i.e. P1 -> P1x -
> P1xave, P1y - P1yave etc. Do the same for P2. Calculate the
> transformation again. The shifts should now be zero but your rotations
> and scale should be statistically independent and correct. The
> difference between the two coordinate systems' centres of mass gives
> the shifts. Transformations far away from the origin rarely behaves
> well.
>
> Good luck.

-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/scipy-user/attachments/20100303/548a759d/attachment.html
```