将“二次”曲面拟合到 3D 数据点答案

【问题标题】：Fitting "quadratic" surface to data points in 3D将“二次”曲面拟合到 3D 数据点
【发布时间】：2018-03-02 02:13:04
【问题描述】：

我正在尝试将二次平面拟合到 python 中的数据点云。我的平面函数的形式是

   f(x,y,z) = a*x**2 + b*y**2 + c*x*y + d*x + e*y + f - z

目前，我的数据点没有与之相关的错误，但是，如有必要，可以假设一些错误。根据here 的建议，我计算出从点 p0=(x0,y0,z0) （对应于我的数据点）到平面上点 p=(x,y,z) 的垂直距离987654322@。然后我结束了

def vertical_distance(params,p0):
    *** snip ***
    nominator = f + a*x**2 + b*y**2 + c*x*y - x0*(2*a*x-c*y-d) - y0*(2*b*y-c*x-e) + z0
    denominator = sqrt((2*a*x+c*y+d)**2 + (2*b*y+c*x+e)**2 + (-1)**2)
    return nominator/denominator

最终，我认为我需要最小化的是 vertical_distance 函数。我可以很高兴地提供一个起始参数列表（参数）和二维数据指向它的数组，但是，我不确定如何在 3D 中实现这一点。 ODR 包似乎只允许包含 x、y 或二维的数据。此外，如何将平面 (p) 上的点实现到最小化例程中？我猜在拟合操作期间，点会根据参数优化而变化，因此此时平面的精确方程也会有所不同。

【问题讨论】：

这些链接用于线性平面。

标签： python python-2.7 math 3d data-fitting

【解决方案1】：

我猜想«quadratic surface» 会比«plane» 更正确。

问题是要适合 z = ax^2 + by^2 + cxy + dx + e y + f 到给定的一组点 P。

要通过优化来做到这一点，您需要制定残差函数（例如，垂直距离）。

对于P残差中的每个3D点p是

|p_2 – ap_0^2 + bp_1^2 + c*p_0*p_1 + dp_0 + ep_1 + f|

您需要最小化所有残差，即它们的平方和，改变参数 a…f。

从技术上讲，以下代码应该可以解决上述问题。但是拟合问题是多极值的，如果没有好的起点或搜索的全球化，这样的例程可能无法找到正确的参数集。

import numpy
import scipy.optimize

P = numpy.random.rand(3,10) # given point set

def quadratic(x,y, a, b, c, d, e, f): 
  #fit quadratic surface
  return a*x**2 + b*y**2 + c*x*y + d*x + e*y + f

def residual(params, points):
  #total residual
  residuals = [
    p[2] - quadratic(p[0], p[1],
    params[0], params[1], params[2], params[3], params[4], params[5]) for p in points]

  return numpy.linalg.norm(residuals)

result = scipy.optimize.minimize(residual, 
                                 (1, 1, 0, 0, 0, 0),#starting point
                                 args=P)

【讨论】：