拟合点到平面算法，如何迭代结果?

编程入门行业动态更新时间:2024-10-22 18:32:00

本文介绍了拟合点到平面算法，如何迭代结果?的处理方法，对大家解决问题具有一定的参考价值，需要的朋友们下面随着小编来一起学习吧！问题描述

更新:我修改了Optimize和Eigen and Solve方法以反映更改.现在，所有函数都返回相同"向量，以实现机器精度. 我仍然对本征方法感到困惑.具体来说，如何/为什么选择特征向量的切片没有意义.直到正常情况与其他解决方案匹配之前，这只是反复试验.如果有人可以纠正/解释我真正应该做的事情，或者为什么我做过的事情行得通，那么我将不胜感激..

谢谢亚历山大·克雷默(Alexander Kramer)解释了我为什么要切片的原因，只允许选择一个正确的答案

我有一个深度图像.我想计算深度图像中像素的粗表面法线.我考虑周围的像素(在最简单的情况下为3x3矩阵)，然后将平面拟合到这些点，然后计算该平面的法向单位矢量.

I have a depth image. I want to calculate a crude surface normal for a pixel in the depth image. I consider the surrounding pixels, in the simplest case a 3x3 matrix, and fit a plane to these point, and calculate the normal unit vector to this plane.

听起来很简单，但是最好先验证平面拟合算法.搜索SO和其他各种站点，我看到了使用最小二乘法，奇异值分解，特征向量/值等的方法.

Sounds easy, but thought best to verify the plane fitting algorithms first. Searching SO and various other sites I see methods using least squares, singlualar value decomposition, eigenvectors/values etc.

尽管我不完全理解数学，但我已经能够使各种片段/示例正常工作.我遇到的问题是，每种方法的答案都不同.我原以为各种答案将是相似的(不完全是)，但是它们似乎有很大的不同.也许有些方法不适合我的数据，但是不确定为什么我得到不同的结果.有什么想法吗?

Although I don't fully understand the maths I have been able to get the various fragments/example to work. The problem I am having, is that I am getting different answers for each method. I was expecting the various answers would be similar (not exact), but they seem significantly different. Perhaps some methods are not suited to my data, but not sure why I am getting different results. Any ideas why?

这是代码的更新后的输出 :

Here is the Updated output of the code:

LTSQ: [ -8.10792259e-17 7.07106781e-01 -7.07106781e-01] SVD: [ 0. 0.70710678 -0.70710678] Eigen: [ 0. 0.70710678 -0.70710678] Solve: [ 0. 0.70710678 0.70710678] Optim: [ -1.56069661e-09 7.07106781e-01 7.07106782e-01]

以下代码实现了五种不同的方法来计算平面的表面法线.算法/代码来自互联网上的各个论坛.

The following code implements five different methods to calculate the surface normal of a plane. The algorithms/code were sourced from various forums on the internet.

import numpy as np import scipy.optimize def fitPLaneLTSQ(XYZ): # Fits a plane to a point cloud, # Where Z = aX + bY + c ----Eqn #1 # Rearanging Eqn1: aX + bY -Z +c =0 # Gives normal (a,b,-1) # Normal = (a,b,-1) [rows,cols] = XYZ.shape G = np.ones((rows,3)) G[:,0] = XYZ[:,0] #X G[:,1] = XYZ[:,1] #Y Z = XYZ[:,2] (a,b,c),resid,rank,s = np.linalg.lstsq(G,Z) normal = (a,b,-1) nn = np.linalg.norm(normal) normal = normal / nn return normal def fitPlaneSVD(XYZ): [rows,cols] = XYZ.shape # Set up constraint equations of the form AB = 0, # where B is a column vector of the plane coefficients # in the form b(1)*X + b(2)*Y +b(3)*Z + b(4) = 0. p = (np.ones((rows,1))) AB = np.hstack([XYZ,p]) [u, d, v] = np.linalg.svd(AB,0) B = v[3,:]; # Solution is last column of v. nn = np.linalg.norm(B[0:3]) B = B / nn return B[0:3] def fitPlaneEigen(XYZ): # Works, in this case but don't understand! average=sum(XYZ)/XYZ.shape[0] covariant=np.cov(XYZ - average) eigenvalues,eigenvectors = np.linalg.eig(covariant) want_max = eigenvectors[:,eigenvalues.argmax()] (c,a,b) = want_max[3:6] # Do not understand! Why 3:6? Why (c,a,b)? normal = np.array([a,b,c]) nn = np.linalg.norm(normal) return normal / nn def fitPlaneSolve(XYZ): X = XYZ[:,0] Y = XYZ[:,1] Z = XYZ[:,2] npts = len(X) A = np.array([ [sum(X*X), sum(X*Y), sum(X)], [sum(X*Y), sum(Y*Y), sum(Y)], [sum(X), sum(Y), npts] ]) B = np.array([ [sum(X*Z), sum(Y*Z), sum(Z)] ]) normal = np.linalg.solve(A,B.T) nn = np.linalg.norm(normal) normal = normal / nn return normal.ravel() def fitPlaneOptimize(XYZ): def residiuals(parameter,f,x,y): return [(f[i] - model(parameter,x[i],y[i])) for i in range(len(f))] def model(parameter, x, y): a, b, c = parameter return a*x + b*y + c X = XYZ[:,0] Y = XYZ[:,1] Z = XYZ[:,2] p0 = [1., 1.,1.] # initial guess result = scipy.optimize.leastsq(residiuals, p0, args=(Z,X,Y))[0] normal = result[0:3] nn = np.linalg.norm(normal) normal = normal / nn return normal if __name__=="__main__": XYZ = np.array([ [0,0,1], [0,1,2], [0,2,3], [1,0,1], [1,1,2], [1,2,3], [2,0,1], [2,1,2], [2,2,3] ]) print "Solve: ", fitPlaneSolve(XYZ) print "Optim: ",fitPlaneOptimize(XYZ) print "SVD: ",fitPlaneSVD(XYZ) print "LTSQ: ",fitPLaneLTSQ(XYZ) print "Eigen: ",fitPlaneEigen(XYZ)

推荐答案

优化

平面a * x + b * y + c * z = 0的法线向量等于(a，b，c)

Optimize

The normal vector of a plane a*x + b*y +c*z = 0, equals (a,b,c)

优化方法找到a和b的值，使得a * x + b * y〜z(〜表示近似值)在计算中完全不使用c的值.我没有在这台计算机上安装numpy，但我希望将模型更改为(a * x + b * y)/c应该可以解决此问题.对于所有数据集，它不会给出相同的结果.此方法将始终假定一个平面穿过原点.

The optimize method finds a values for a and b such that a*x+b*y~z (~ denotes approximates) It omits to use the value of c in the calculation at all. I don't have numpy installed on this machine but I expect that changing the model to (a*x+b*y)/c should fix this method. It will not give the same result for all data-sets. This method will always assume a plane that goes through the origin.

产生相同的结果. (区别在于机器精度的大小.)

produce the same results. (The difference is about the size of machine precision).

选择了错误的特征向量.就像在SVD和LTSQ中一样，对应于最大特征值(lambda = 1.50)的特征向量为x=[0, sqrt(2)/2, sqrt(2)/2].

The wrong eigenvector is chosen. The eigenvector corresponding to the greatest eigenvalue (lambda = 1.50) is x=[0, sqrt(2)/2, sqrt(2)/2] just as in the SVD and LTSQ.