三点定圆原理与C++实现
文章目录
1 原理2 C++实现
1 原理
根据我们小学二年级就学过的三点定圆定理：
不
共
线
的
三
个
点
可
唯
一
确
定
一
个
圆
不共线的三个点可唯一确定一个圆
不共线的三个点可唯一确定一个圆
且，不共线的三点相互连接必然构成一个三角形，这个三角形称为圆的内接三角形，这个圆称为三角形的外接圆。三角形三边垂直平分线的交点即为三角形外接圆的圆心。
有了以上知识，由不共线的三个点确定一个圆就非常easy了
已知平面中不共线的三点
A
(
x
1
,
y
1
,
z
1
)
A(x_1,y_1,z_1)
A(x1,y1,z1)，
B
(
x
2
,
y
2
,
z
2
)
B(x_2,y_2,z_2)
B(x2,y2,z2)，
C
(
x
3
,
y
3
,
z
3
)
C(x_3,y_3,z_3)
C(x3,y3,z3)，互相连接构成三角形 ΔABC，
L
a
b
L_{ab}
Lab，
L
b
c
L_{bc}
Lbc，
L
a
c
L_{ac}
Lac 分别为三条边的垂直平分线，且相交于一点
O
O
O，该交点即为外接圆圆心。
三角形三边斜率：
{
k
a
b
=
y
2
−
y
1
x
2
−
x
1
k
b
c
=
y
3
−
y
2
x
3
−
x
2
k
a
c
=
y
3
−
y
1
x
3
−
x
1
\begin{cases} k_{ab}=\cfrac {y_2-y_1}{x_2-x_1}\\ k_{bc}=\cfrac {y_3-y_2}{x_3-x_2}\\ k_{ac}=\cfrac {y_3-y_1}{x_3-x_1}\\ \end{cases}
⎩⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎧kab=x2−x1y2−y1kbc=x3−x2y3−y2kac=x3−x1y3−y1 三边中点:
P
a
b
(
x
1
+
x
2
2
，
y
1
+
y
2
2
)
,
P
b
c
(
x
2
+
x
3
2
，
y
2
+
y
3
2
)
,
P
a
c
(
x
1
+
x
3
2
，
y
1
+
y
3
2
)
P_{ab}(\cfrac {x_1+x_2}{2}，\cfrac {y_1+y_2}{2}), P_{bc}(\cfrac {x_2+x_3}{2}，\cfrac {y_2+y_3}{2}), P_{ac}(\cfrac {x_1+x_3}{2}，\cfrac {y_1+y_3}{2})
Pab(2x1+x2，2y1+y2),Pbc(2x2+x3，2y2+y3),Pac(2x1+x3，2y1+y3)
根据直线的点斜式方程，可求三条垂直平分线方程如下：
{
L
a
b
=
y
−
y
1
+
y
2
2
+
1
k
a
b
(
x
−
x
1
+
x
2
2
)
=
0
L
b
c
=
y
−
y
2
+
y
3
2
+
1
k
b
c
(
x
−
x
2
+
x
3
2
)
=
0
L
a
c
=
y
−
y
1
+
y
3
2
+
1
k
a
c
(
x
−
x
1
+
x
3
2
)
=
0
\begin{cases} L_{ab}=y-\cfrac {y_1+y_2}{2}+\cfrac {1}{k_{ab}}(x-\cfrac {x_1+x_2}{2})=0\\ L_{bc}=y-\cfrac {y_2+y_3}{2}+\cfrac {1}{k_{bc}}(x-\cfrac {x_2+x_3}{2})=0\\ L_{ac}=y-\cfrac {y_1+y_3}{2}+\cfrac {1}{k_{ac}}(x-\cfrac {x_1+x_3}{2})=0\\ \end{cases}
⎩⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎧Lab=y−2y1+y2+kab1(x−2x1+x2)=0Lbc=y−2y2+y3+kbc1(x−2x2+x3)=0Lac=y−2y1+y3+kac1(x−2x1+x3)=0
三个方程，两个未知量，因此任选其中两个方程求解圆心坐标
O
(
x
0
,
y
0
)
O(x_0,y_0)
O(x0,y0)，则圆心
O
O
O到点
A
、
B
、
C
A、B、C
A、B、C任一点的距离就是外接圆的半径，即
r
=
∣
O
A
∣
=
∣
O
B
∣
=
∣
O
C
∣
=
(
x
0
−
x
1
)
2
+
(
y
0
−
y
1
)
2
r=|OA|=|OB|=|OC|=\sqrt {{(x_0-x_1)^2}+{(y_0-y_1)^2}}
r=∣OA∣=∣OB∣=∣OC∣=(x0−x1)2+(y0−y1)2
至此，可求得外接圆的方程为
(
x
−
x
0
)
2
+
(
y
−
y
0
)
2
=
r
2
(x-x_0)^2+(y-y_0)^2=r^2
(x−x0)2+(y−y0)2=r2
2 C++实现
以圆：
(
x
−
1
)
2
+
(
y
−
1
)
2
=
1
(x-1)^2+(y-1)^2=1
(x−1)2+(y−1)2=1为例，选择圆上三点
p
1
(
0
,
1
)
、
p
2
(
1
,
2
)
、
p
3
(
1.5
,
0.75
+
1
)
p_1(0,1)、p_2(1,2)、p_3(1.5,\sqrt{0.75}+1)
p1(0,1)、p2(1,2)、p3(1.5,0.75
+1) 进行验证。
大家也可以自行输入其他点进行验证，只需在main.cpp中修改三个点的坐标即可。
注意： 需要用到Eigen库求解二元一次方程组，PointXY类型的点；可自己定义一个结构体PointXY，也可直接使用PCL库里面定义好的结构体
mian.cpp
#include "deifen_circle_with_3points.h"
int main()
{
PointT p1, p2, p3; //定义三个点
p1.x = 0.0;
p1.y = 1.0;
p2.x = 1.0;
p2.y = 2.0;
p3.x = 1.5;
p3.y = sqrt(0.75) + 1.0;
DefineCircl3Points dc; //定义三点定圆对象dc
dc.setThreePoints(p1, p2, p3); //设置三点
dc.defineCircle(); //执行三点定圆
return 0;
}
deifen_circle_with_3points.h
#pragma once
#include
#include
using namespace std;
using namespace Eigen;
typedef pcl::PointXY PointT;
class DefineCircl3Points
{
public:
/**
* @brief : 输入三点
* @param[I]: p1 (x1,y1)
* @param[I]: p2 (x2,y2)
* @param[I]: p3 (x3,y3)
* @param[O]: none
* @return : none
* @note :
**/
void setThreePoints(PointT &p1, PointT &p2, PointT &p3);
/**
* @brief : 三点定圆
* @param[I]: none
* @param[O]: none
* @return : none
* @note : 输出圆心、半径、圆方程
**/
void defineCircle();
private:
PointT m_p1, m_p2, m_p3; //输入的三点
bool is_set3Points = false; //是否输入三点
};
deifen_circle_with_3points.cpp
#include "deifen_circle_with_3points.h"
/**
* @brief : 输入三点
* @param[I]: p1 (x1,y1)
* @param[I]: p2 (x2,y2)
* @param[I]: p3 (x3,y3)
* @param[O]: none
* @return : none
* @note :
**/
void DefineCircl3Points::setThreePoints(PointT & p1, PointT & p2, PointT & p3)
{
//判断是否三点共线
if ((p2.x - p1.x)*(p3.y - p1.y) - (p3.x - p1.x)*(p2.y - p1.y) == 0)
{
PCL_ERROR("->三点共线，无法确定一个圆！\a\n");
system("pause");
abort();
}
m_p1 = p1;
m_p2 = p2;
m_p3 = p3;
is_set3Points = true;
}
/**
* @brief : 三点定圆
* @param[I]: none
* @param[O]: none
* @return : none
* @note : 输出圆心、半径、圆方程
**/
void DefineCircl3Points::defineCircle()
{
if (!is_set3Points)
{
PCL_ERROR("->请输入三个点！\a\n");
system("pause");
abort();
}
float kab, kbc, kac; //三边斜率
kab = (m_p2.y - m_p1.y) / (m_p2.x - m_p1.x);
kbc = (m_p3.y - m_p2.y) / (m_p3.x - m_p2.x);
kac = (m_p3.y - m_p1.y) / (m_p3.x - m_p1.x);
PointT Pab, Pbc, Pac; //三边中点
Pab.x = (m_p1.x + m_p2.x) / 2;
Pab.y = (m_p1.y + m_p2.y) / 2;
Pbc.x = (m_p2.x + m_p3.x) / 2;
Pbc.y = (m_p2.y + m_p3.y) / 2;
Pac.x = (m_p1.x + m_p3.x) / 2;
Pac.y = (m_p1.y + m_p3.y) / 2;
//矩阵方程Ga=d
Matrix2f G;
G << 1 / kab, 1, 1 / kbc, 1;
Vector2f d;
d << (Pab.y + 1 / kab * Pab.x), (Pbc.y + 1 / kbc * Pbc.x);
Vector2f a;
a = G.colPivHouseholderQr().solve(d);
float r; //圆半径
r = sqrt(pow((a[0, 0] - m_p1.x), 2) + pow((a[1, 0] - m_p1.y), 2));
cout << "圆心：" << "(" << a[0, 0] << ", " << a[1, 0] << ")" << endl;
cout << "半径：" << r << endl;
cout << "->圆方程："
<< "(x - " << a[0, 0] << " )^2 + (y - " << a[1, 0] << " )^2 = " << r << "^2" << endl;
}
输出结果：
圆心：(1, 1)
半径：1
->圆方程：(x - 1 )^2 + (y - 1 )^2 = 1^2
当然，也可以直接调用PCL库RANSAC圆拟合的方式计算圆心和半径，只需将最大迭代次数设置为1即可
ransac.setMaxIterations(1);
相关链接
PCL随机采样一致性：RANSAC 圆拟合（二维圆 + 空间圆）
线性方程组/矩阵方程求解(方法汇总)
Markdown输入公式、符号等语法指令汇总