3D任意平面含噪圆周点的圆心优化算法求助

1. 先拟合点集所在的3D平面 #

所有含噪点都落在某一3D平面上，第一步要确定平面参数：

• 计算点集的质心 $\mathbf{C} = (\bar{x}, \bar{y}, \bar{z})$，其中 $\bar{x} = \frac{1}{n}\sum_{k=1}^n x_k$，同理可得$\bar{y}$、$\bar{z}$
• 构造协方差矩阵 $\mathbf{M}$：
  $$
  \mathbf{M} = \begin{bmatrix}
  \sum(x_k-\bar{x})^2 & \sum(x_k-\bar{x})(y_k-\bar{y}) & \sum(x_k-\bar{x})(z_k-\bar{z}) \
  \sum(y_k-\bar{y})(x_k-\bar{x}) & \sum(y_k-\bar{y})^2 & \sum(y_k-\bar{y})(z_k-\bar{z}) \
  \sum(z_k-\bar{z})(x_k-\bar{x}) & \sum(z_k-\bar{z})(y_k-\bar{y}) & \sum(z_k-\bar{z})^2 \
  \end{bmatrix}
  $$
• 对$\mathbf{M}$做特征值分解，最小特征值对应的特征向量就是平面的单位法向量 $\mathbf{n}=(a,b,c)$
• 代入质心求平面方程的常数项：$d = -(a\bar{x} + b\bar{y} + c\bar{z})$，平面方程为 $ax + by + cz + d = 0$

2. 将3D点投影到平面，转为2D圆拟合问题 #

在拟合平面上建立局部正交坐标系：

• 以质心$\mathbf{C}$为原点
• 取与$\mathbf{n}$正交的任意单位向量作为$\mathbf{u}$（比如若$\mathbf{n}$不是$(0,0,1)$，则取$\mathbf{u} = \frac{(0,1,0) \times \mathbf{n}}{|(0,1,0) \times \mathbf{n}|}$，否则取$\mathbf{u}=(1,0,0)$）
• 计算$\mathbf{v} = \mathbf{n} \times \mathbf{u}$，确保$\mathbf{u}$、$\mathbf{v}$正交且单位化
• 每个3D点$\mathbf{P}_k$转换为2D坐标：$\mathbf{Q}_k = \mathbf{P}_k - \mathbf{C}$，$u_k = \mathbf{Q}_k \cdot \mathbf{u}$，$v_k = \mathbf{Q}_k \cdot \mathbf{v}$

3. 2D含噪圆的最小二乘拟合 #

代数法（线性化求解） #

将圆的方程改写为 $u^2 + v^2 + Au + Bv + C = 0$，其中圆心为 $(-A/2, -B/2)$，半径 $r = \sqrt{(A/2)^2 + (B/2)^2 - C}$。
构造线性方程组求解$A,B,C$：

• 对每个点$(u_k, v_k)$，有 $u_k^2 + v_k^2 = -A u_k - B v_k - C$
• 写成矩阵形式 $\mathbf{y} = \mathbf{X}\boldsymbol{\beta}$，其中：
  $\mathbf{y} = [u_1^2+v_12, u_2^2+v_22, ..., u_n^2+v_n2]^T$
  $\mathbf{X} = \begin{bmatrix} -u_1 & -v_1 & -1 \ -u_2 & -v_2 & -1 \ ... & ... & ... \ -u_n & -v_n & -1 \end{bmatrix}$
  $\boldsymbol{\beta} = [A, B, C]^T$
• 通过最小二乘求解：$\boldsymbol{\beta} = (\mathbf{X}^T\mathbf{X}){-1}\mathbf{X}^T\mathbf{y}$，再转换为圆心坐标

非线性优化法（更鲁棒） #

若想进一步提升精度，可直接最小化点到圆的距离误差：
$$E = \sum_{k=1}^n \left( \sqrt{(u_k - u_0)^2 + (v_k - v_0)^2} - r \right)^2$$
用Levenberg-Marquardt算法求解这个非线性最小二乘问题，初始值可选用代数法的结果。

4. 将2D圆心转回3D空间 #

得到2D圆心$(u_0, v_0)$后，对应的3D圆心为：
$$\mathbf{O} = \mathbf{C} + u_0 \cdot \mathbf{u} + v_0 \cdot \mathbf{v}$$

直接3D空间优化的替代方案 #

如果不想做平面投影，可直接在3D空间中带约束优化：

• 约束条件：圆心$\mathbf{O}$必须在拟合平面上（即 $a x_0 + b y_0 + c z_0 + d = 0$）
• 误差函数：$E = \sum_{k=1}^n \left( |\mathbf{P}_k - \mathbf{O}| - r \right)^2$
• 利用约束消去一个变量（比如用平面方程把$z_0$表示为$x_0,y_0$的函数），将问题转化为二元变量的非线性优化，再用Levenberg-Marquardt或梯度下降求解。

方法对比 #

你之前用的三点取平均法对噪声和 outliers 抗性较差，因为单个三点组合可能受极端点影响。上述方法利用所有点的信息，通过最小二乘或非线性优化抑制噪声，拟合结果更稳定鲁棒。

内容的提问来源于stack exchange，提问作者alexismue