计算几何摸黑(3)-圆及有关计算(上)
圆的表示方法
和直线一样，圆也可以使用参数式表示。显而易见，一个唯一的圆可以用一个点表示圆心，和一个实数表示半径。所以，我们可以这样定义一个圆。
class circle_base{public: double r; Point c; Circle(Point C = Point(0.0, 0.0), R = 0.0) : c(C), r(R) { }}C; 约定：这里用C.c表示圆C的圆心，C.r表示圆C的半径。
通过圆心角求点的坐标
我们作一条过圆心，平行于\(x\)轴的直线\(lx\)，再作一条过圆心，平行于\(y\)轴的直线\(ly\)，那么这又构成了一个坐标系。这样，在圆上且和\(lx\)正半轴有一个夹角的点唯一。我们把那个角称为圆心角。所以，我们可以用三角函数和圆的成分来实现这个功能。注意：在这里我用到了类继承，也可以直接将Circle类声明中添加point函数。
class Circle : public circle_base{public: Circle(circle_base &t) : circle_base(t) { } inline Point point(const double theta) { return Point(x + std::cos(theta) * r, y + std::sin(theta) * r); }} 直线和圆的交点
我们可以通过解方程法求出直线和圆的交点。 若直线为\(AB\)，我们可以设交点为\(P=A+t(B-A)\)。这里要用到高中课本上的圆方程公式：
\[ (x - x_0) ^ 2 + (y - y_0) ^ 2 = r ^ 2 \] 由上面的设出的式子，我们可以这样推导：\(P = (ax, ay) + t(bx - ax, by - ay)\)
$ = (ax - tax + tbx, ay - tay + tby)$
\(\because P \in \odot C\)
\(\therefore (ax - tax + tbx - C.c.x)^2 +(ay - tay + tby - C.c.y)^2=r^2\)
\(\because bx - ax = \vec {AB}.x, by - ay = \vec {AB}.y\)
\(整理得，(t(\vec {AB}.x) + (ax - C.c.x))^2 + (t(\vec {AB}.y) + (ay - C.c.y))^2 - r^2 = 0\)
\(令a = (\vec {AB}.x), b = (ax - C.c.x), c = (\vec{AB}.y), d = (ay - C.c.y)\)
\(\therefore (at+b)^2+(ct+d)^2 - r^2 = 0\)
\((a^2+c^2)t^2+(ab+cd)t+b^2+d^2-r^2=0\)
\(令e = (a^2 + c^2), f = (ab + cd), g = (b^2 + d^2 - r^2)\)
\(\therefore et^2 + ft + g = 0\)
\(令\Delta= f^2 - 4eg\)
易知， 如果\(\Delta > 0\)，那么直线将与圆有两个交点，\(\Delta < 0\)，直线将与圆有相离，\(\Delta = 0\)，直线将与圆有一个交点。而把这个一元二次方程的解带入，并用三角函数计算就可以求出角了。实现如下。
int getLineCircInter(const Line &L, const Circle &C, double &t1, double &t2, double &s1, double &s2) // 意义: L交C于圆心角为t1和t2的点，分别为s1, s2{ double a = L.v.x, b = L.p.x - C.c.x, c = L.v.y, d = L.p.y - C.c.y, e = a * a + c * c, f = 2 * (a * b + c * d), g = (b * b + d * d - C.r * C.r); // 详细推导过程见上文 double delta = f * f - 4 * a * c; // 判别式 if (dcmp(delta) < 0) return 0; if (dcmp(delta) == 0) { t1 = t2 = -(f / (2 * e)); s1 = s2 = C.point(t1); return 1; } t1 = (-f + std::sqrt(delta)) / (2 * e); t2 = (-f - std::sqrf(delta)) / (2 * e); s1 = C.point(t1); s2 = C.point(t2); return 2;} // 注意，这里的返回值的含义是交点个数。