matlab7.0多点三维坐标拟合圆求圆心坐标及拟合半径的实现

设该圆心坐标为r(x,y,z)，半径为R，则拟合最优函数及其导数满足：

所要做的就是求解偏导数为0的4个非线性方程，可以用matlab的fsolve求解

输入的离散点在
变量X中
X=[x1 y1 z1;
x2 y2 z2;
];
是一个n行，3列的矩阵，1，2，3列分别为x,y,z坐标
p0=mean(X);%以重心为初始圆心
r02=mean(sqrt(sum((X-repmat(p0,[size(X,1) 1]))^2,2)));%以到重心的平均距离为初始半径
f=@(p)(X(:,1)-p(1))^2+(X(:,2)-p(2))^2+(X(:,3)-p(3))^2-p(4)^2;
p=lsqnonlin(f,[p0 r02]);
x0=p(1);y0=p(2);z0=p(3);
radii=abs(p(4));
x0 y0 z0 为圆心坐标
radii为半径。

很久没搞立体几何运算，所以给不出答案。但可以给点建议。
空间不在一条线上的3点，可以确定1个圆。
空间圆的参数方程：
x = Cx + r (Uxcos(t) + Vxsin(t));
y = Cy + r (Uycos(t) + Vysin(t));
x = Cz + r (Uzcos(t) + Vzsin(t));
Cx,Cy,Cz 是圆心 x,y,z 座标。r 是半径。
Ux,Uy,Uz 是 与 圆所在平面的 法线 方向 有关 的 单位向量（规一化后的向量）
Vx,Vy,Vz 是 与 圆所在平面的 切线 方向 有关 的 单位向量 （规一化后的向量）
-------
空间不在一条线上的3点，可以组成2个向量，例如
( x3-x1) i + (y3-y1) j + (z3-z1) k
( x2-x1) i + (y2-y1) j + (z2-z1) k
2个向量 的 叉乘积 的 方向 就是 法线方向。
2个向量 的 点乘积 的 方向 就是 切线 方向。
叉乘积公式 用 行列式记：
白白白白| i j k |
a X b = | a1 a2 a3|
白白白白| b1 b2 b3|
点乘积 公式简单，略。
把 方向 系数 代入 空间圆的参数方程 求解便可。
------
还有个关系式 是 半径的 平方 等于 圆心到 圆上一点 空间距离。
r r = (x-Cx)^2 + (y-Cy)^2 + (z-Cz)^2

1、在平面内画一个圆，如下图所示；

2、在圆内任意做两条不平行的弦，如下图所示；

3、在园内的两条弦上分别作中垂线，下图为弦1的中垂线；

4、继续做弦2的中垂线，如下图所示；

5、两条中垂线的交点位置，即下图中的点A便是圆的中心点及圆的圆心。

---