方法
联立方程组,使用matlab的solve求解符号方程组;
化简表达式;
编写C++代码。
方程组
{ ( x − a 1 ) 2 + ( y − b 1 ) 2 = R 1 2 ( x − a 2 ) 2 + ( y − b 2 ) 2 = R 2 2 \begin{cases}(x-a_1)^2+(y-b_1)^2=R_1^2\\(x-a_2)^2+(y-b_2)^2=R_2^2\end{cases} {(x−a1)2+(y−b1)2=R12(x−a2)2+(y−b2)2=R22
matlab求解代码:
clear all
syms a1 b1 R1 a2 b2 R2 x y;
[x,y]=solve((x-a1)*(x-a1)+(y-b1)*(y-b1)-R1*R1,(x-a2)*(x-a2)+(y-b2)*(y-b2)-R2*R2);
%simplify(x)
%simplify(y)
再用simplify()化简表达式,得到解析解。
C++实现
typedef struct
{
float x;
float y;
}Point2f;//类似opencv
typedef struct
{
Point2f c;
float r;
}CIRCLE;
void IntersectionOf2Circles(CIRCLE c1, CIRCLE c2, Point2f &P1, Point2f &P2)
{
float a1, b1, R1, a2, b2, R2;
a1 = c1.c.x;
b1 = c1.c.y;
R1 = c1.r;
a2 = c2.c.x;
b2 = c2.c.y;
R2 = c2.r;
//
float R1R1 = R1*R1;
float a1a1 = a1*a1;
float b1b1 = b1*b1;
float a2a2 = a2*a2;
float b2b2 = b2*b2;
float R2R2 = R2*R2;
float subs1 = a1a1 - 2 * a1*a2 + a2a2 + b1b1 - 2 * b1*b2 + b2b2;
float subs2 = -R1R1 * a1 + R1R1 * a2 + R2R2 * a1 - R2R2 * a2 + a1a1*a1 - a1a1 * a2 - a1*a2a2 + a1*b1b1 - 2 * a1*b1*b2 + a1*b2b2 + a2a2*a2 + a2*b1b1 - 2 * a2*b1*b2 + a2*b2b2;
float subs3 = -R1R1 * b1 + R1R1 * b2 + R2R2 * b1 - R2R2 * b2 + a1a1*b1 + a1a1 * b2 - 2 * a1*a2*b1 - 2 * a1*a2*b2 + a2a2 * b1 + a2a2 * b2 + b1b1*b1 - b1b1 * b2 - b1*b2b2 + b2b2*b2;
float sigma = sqrt((R1R1 + 2 * R1*R2 + R2R2 - a1a1 + 2 * a1*a2 - a2a2 - b1b1 + 2 * b1*b2 - b2b2)*(-R1R1 + 2 * R1*R2 - R2R2 + subs1));
if(abs(subs1)>0.0000001)//分母不为0
{
P1.x = (subs2 - sigma*b1 + sigma*b2) / (2 * subs1);
P2.x = (subs2 + sigma*b1 - sigma*b2) / (2 * subs1);
P1.y = (subs3 + sigma*a1 - sigma*a2) / (2 * subs1);
P2.y = (subs3 - sigma*a1 + sigma*a2) / (2 * subs1);
}
}
完成!第一次发csdn博客!!