如何使用python执行坐标仿射变换?

dai*_*ini 4 python coordinate-transformation

我想为此示例数据集执行转换.
在一个坐标[primary_system]系统中有四个已知点坐标x,y,z,并且接下来四个已知点具有属于另一个坐标系[secondary_system]的坐标x,y,h.那些点对应; 例如,primary_system1 point和secondary_system1点完全相同,但我们在两个不同的坐标系中有它的坐标.所以我在这里有四对调整点,并希望根据调整将另一个点坐标从主系统转换到二次系统.

primary_system1 = (3531820.440, 1174966.736, 5162268.086)
primary_system2 = (3531746.800, 1175275.159, 5162241.325)
primary_system3 = (3532510.182, 1174373.785, 5161954.920)
primary_system4 = (3532495.968, 1175507.195, 5161685.049)

secondary_system1 = (6089665.610, 3591595.470, 148.810)
secondary_system2 = (6089633.900, 3591912.090, 143.120)
secondary_system3 = (6089088.170, 3590826.470, 166.350)
secondary_system4 = (6088672.490, 3591914.630, 147.440)

#transform this point
x = 3532412.323 
y = 1175511.432
z = 5161677.111<br>
Run Code Online (Sandbox Code Playgroud)


目前我尝试使用四对点中的每一对来平均x,y和z轴的平移,如:

#x axis
xt1 =  secondary_system1[0] - primary_system1[0]           
xt2 =  secondary_system2[0] - primary_system2[0]
xt3 =  secondary_system3[0] - primary_system3[0]
xt4 =  secondary_system4[0] - primary_system4[0]

xt = (xt1+xt2+xt3+xt4)/4    #averaging
Run Code Online (Sandbox Code Playgroud)

...等等y轴和z轴

#y axis
yt1 =  secondary_system1[1] - primary_system1[1]           
yt2 =  secondary_system2[1] - primary_system2[1]
yt3 =  secondary_system3[1] - primary_system3[1]
yt4 =  secondary_system4[1] - primary_system4[1]

yt = (yt1+yt2+yt3+yt4)/4    #averaging

#z axis
zt1 =  secondary_system1[2] - primary_system1[2]           
zt2 =  secondary_system2[2] - primary_system2[2]
zt3 =  secondary_system3[2] - primary_system3[2]
zt4 =  secondary_system4[2] - primary_system4[2]

zt = (zt1+zt2+zt3+zt4)/4    #averaging
Run Code Online (Sandbox Code Playgroud)

所以上面我尝试计算每个轴的平均平移向量

mat*_*fee 9

如果它只是一个平移和旋转,那么这就是一个称为仿射变换的变换.

它基本上采取以下形式:

secondary_system = A * primary_system + b
Run Code Online (Sandbox Code Playgroud)

哪里A是3x3矩阵(因为你是3D),并且b是3x1平移.

这可以等效地写出来

secondary_system_coords2 = A2 * primary_system2,
Run Code Online (Sandbox Code Playgroud)

哪里

(有关详细信息,请参阅维基页面).

所以基本上,你想要解决这个等式:

y = A2 x
Run Code Online (Sandbox Code Playgroud)

for A2,其中y包括secondary_system1个卡在末端的x点,并且是primary_system1个卡在末端的点,并且A2是4x4矩阵.

现在,如果x是方阵,我们可以解决它:

A2 = y*x^(-1)
Run Code Online (Sandbox Code Playgroud)

但是x4x1.但是,你很幸运,有4x4对应的套装y,所以你可以x像这样构造一个4x4:

x = [ primary_system1 | primary_system2 | primary_system3 | primary_system4 ]
Run Code Online (Sandbox Code Playgroud)

其中每个primary_systemi都是4x1列向量.与...相同y.

一旦你有了A2,将一个点从system1转换为系统2,你只需:

transformed = A2 * point_to_transform
Run Code Online (Sandbox Code Playgroud)

你可以numpy像这样设置它(例如):

import numpy as np
def solve_affine( p1, p2, p3, p4, s1, s2, s3, s4 ):
    x = np.transpose(np.matrix([p1,p2,p3,p4]))
    y = np.transpose(np.matrix([s1,s2,s3,s4]))
    # add ones on the bottom of x and y
    x = np.vstack((x,[1,1,1,1]))
    y = np.vstack((y,[1,1,1,1]))
    # solve for A2
    A2 = y * x.I
    # return function that takes input x and transforms it
    # don't need to return the 4th row as it is 
    return lambda x: (A2*np.vstack((np.matrix(x).reshape(3,1),1)))[0:3,:]
Run Code Online (Sandbox Code Playgroud)

然后像这样使用它:

transformFn = solve_affine( primary_system1, primary_system2, 
                            primary_system3, primary_system4,
                            secondary_system1, secondary_system2,
                            secondary_system3, secondary_system4 )

# test: transform primary_system1 and we should get secondary_system1
np.matrix(secondary_system1).T - transformFn( primary_system1 )
# np.linalg.norm of above is 0.02555

# transform another point (x,y,z).
transformed = transformFn((x,y,z))
Run Code Online (Sandbox Code Playgroud)

注意:这里当然存在数值误差,这可能不是解决变换的最佳方法(您可能能够做某种最小二乘的事情).

此外,转换primary_systemx为的错误secondary_systemx(对于此示例)为10 ^( - 2).

你必须考虑这是否可以接受(它确实看起来很大,但与你的输入点相比可能是可以接受的,这些输入点都是10 ^ 6的顺序).