圆是任意一点拥有唯一的圆心角。所以在定义圆的时候,可以加一个圆心角求坐标的函数。
首先是辅助宏:
const double eps = 1e-10;
const double PI = acos(-1); // π
const double TWO_PI = 2*PI;
const int maxn = 100 + 5;
const int INF = 10000;
圆的存储结构定义:
struct Circle
{
Point c;
double r;
Circle(Point c, double r):c(c), r(r){}
Point point(double a)
{
return Point(c.x + r*cos(a), c.y + r*sin(a));
}
};
直线与圆的交点。假定直线为AB,圆的圆心为C,半径为r。
我们可以通过解方程组的方式解决。
设交点为P=A+t(B-A),代入圆方程整理后先得到(at+b)^2 + (ct+d)^2 = r^2 进一步整理后得到一元二次方程et^2 + ft + g=0 根据判别式的值可以区分三种情况,即午无交点,一个交点和两个交点。
函数返回的是交点的个数,参数sol存放的是交点本身。注意上述代码并没有清空sol,这样可以反复调用这个函数,把所有交点放在一个sol里。
int getLineCircleIntersection(Line L, Circle C, double& t1, double& t2, vector<Point>& sol)
{
double a, b, c, d, e, f, g, delta;
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;
//判别式
delta = f*f - 4*e*g;
//相离
if(dcmp(delta) < 0)
return 0;
//相切
if(dcmp(delta) == 0)
{
t1 = -f / (2*e);
t2 = t1;
sol.push_back(L.point(t1));
return 1;
}
//相交
t1 = (-f - sqrt(delta)) / (2*e);
sol.push_back(L.point(t1));
t2 = (-f + sqrt(delta)) / (2*e);
sol.push_back(L.point(t2));
return 2;
}
线段与圆的交点(包括端点):类似处理,只要判断交点是否在线段上就可以。判断的时候我们可以通过判断t1、t2的值来判断,因为如果在线段上的话,t1,t2一定是大于等于0小于等于1的。
int getSegmentCircleIntersection(Point A, Point B, Circle C, double &t1, double &t2, vector<Point> &sol)
{
double a, b, c, d, e, f, g, delta;
Vector v;
int cnt;
v = B - A;
a = v.x;
b = A.x - C.c.x;
c = v.y;
d = A.y - C.c.y;
e = a*a + c*c;
f = 2* (a*b + c*d);
g = b*b + d*d - C.r*C.r;
delta = f*f - 4*e*g;
cnt = 0;
if(dcmp(delta) < 0)
return -1;
else if(dcmp(delta) == 0)
{
t1 = -f / (2*e);
t2 = t1;
if(dcmp(t1) >= 0 && dcmp(t1-1) <= 0)
{
sol.push_back(A + v*t1);
cnt++;
}
return cnt;
}
else
{
t1 = (-f - sqrt(delta)) / (2*e);
t2 = (-f + sqrt(delta)) / (2*e);
if(dcmp(t1) >= 0 && dcmp(t1-1) <= 0)
{
sol.push_back(A + v*t1);
cnt++;
}
if(dcmp(t2) >= 0 && dcmp(t2-1) <= 0)
{
sol.push_back(A + v*t2);
cnt++;
}
return cnt;
}
}
圆和线段是否相交,相切不算,线段不考虑端点.
bool CircleIntersectSegment(Point A, Point B, Circle C)
{
double t1, t2;
int c;
vector<Point> sol;
c = getLineCircleIntersection(Line(A, B-A), C, t1, t2, sol); //线段与圆的交点
if(c <= 1) //少于2个肯定不相交
return false;
if(dcmp(t1) > 0 && dcmp(t1-1) < 0) //端点不在圆上
return true;
if(dcmp(t2) > 0 && dcmp(t2-1) < 0)
return true;
return false;
}
两圆相交。
假定圆心分别为C1和C2,半径为r1和r2,圆心距为d,根据余弦定理可以算出C1C2到C1P1的角da,根据向量C1C2的极角a,加减da就可以得到C1P1和C1P2的极角。有了极角,就可以很方便地计算出P1和P2的坐标了。
int getCircleCircleIntersection(Circle C1, Circle C2, vector<Point>& sol)
{
double d, a, da;
Point p1, p2;
d = Length(C1.c - C2.c);
if(dcmp(d) == 0)
{
if(dcmp(C1.r - C2.r) == 0) //两圆重合
return -1;
return 0;
}
if(dcmp(C1.r + C2.r - d) < 0) //两圆相离
return 0;
if(dcmp(fabs(C1.r - C2.r) - d) > 0) //两圆内含
return 0;
a = angle(C2.c - C1.c); //向量c1c2的极角
da = acos((C1.r*C1.r + d*d - C2.r*C2.r) / (2*C1.r*d));
//c1c2到c1p1的角、
p1 = C1.point(a-da);
p2 = C1.point(a+da);
sol.push_back(p1);
if(p1 == p2)
return 1;
sol.push_back(p2);
return 2;
}
过定点做圆的切线。先求出pq的距离和pc的夹角ang,则向量pc的极角ang就是两条切线的极角。注意切线不存在和只有一条的情况。
过点p到圆C的切线,v[i]是第i条切线的向量。返回切线条数。
int getTangents(Point p, Circle C, Vector *v)
{
double dist, ang;
Vector u;
u = C.c - p;
dist = Length(u);
//没有切线
if(dist < C.r)
return 0;
//p在圆上,只有一条切线
else if(dcmp(dist - C.r) == 0)
{
v[0] = Rotate(u, PI/2);
return 1;
}
ang = asin(C.r / dist);
v[0] = Rotate(u, ang);
v[1] = Rotate(u, -ang);
return 2;
}
两圆的公切线。根据两圆的圆心距从小到大排列。
情况一:两圆完全重合,有无数条公切线。
情况二:两圆内含,没有公共点。没有公切线。
情况三:两圆内切。有一条公切线。
情况四:两圆相交。有两条外公切线。
情况五:两圆外切。有三条公切线,其中一条内公切线,两条外公切线。
情况六:两圆相离。有四条公切线,其中内公切线两条,外公切线两条。
情况三和情况五中的内公切线都对应于“过圆上一点求圆的切线”,只需连接圆心和切点。旋转90度即可直到切线的方向向量。
返回切线的条数,-1表示无数条切线。
int getTangents(Circle A, Circle B, Point *a, Point *b)
{
double d, base, ang;
int cnt;
Vector u;
if(A.r < B.r) //选择圆半径大的为A,便于后面计算
{
swap(A, B);
swap(a, b);
}
cnt = 0;
u = B.c - A.c;
d = Length(u); //圆心距
base = angle(u); //两圆心向量极角
if(dcmp(d) == 0)
{
if(dcmp(A.r - B.r) == 0) //重合
return -1;
return 0;
}
if(dcmp(A.r - B.r - d) > 0) //内含
return 0;
if(dcmp(A.r - B.r - d) == 0) //内切
{
a[cnt] = A.point(base); //根据极角求坐标
b[cnt] = B.point(base);
cnt++;
return 1;
}
ang = acos((A.r-B.r) / d); //向量u与圆心与交点连线之间的角
a[cnt] = A.point(base + ang);
b[cnt] = B.point(base + ang);
cnt++;
a[cnt] = A.point(base - ang);
b[cnt] = B.point(base - ang);
cnt++;
if(dcmp(d - A.r - B.r) == 0) //外切
{
a[cnt] = A.point(base);
b[cnt] = B.point(base + PI); //加PI:因为在圆B中,圆心之间的向量为-u,所以要加上180度
cnt++;
}
else if(dcmp(d - A.r - B.r) > 0) //相离
{
ang = acos((A.r+B.r) / d); //用三角函数公式可求
a[cnt] = A.point(base + ang);
b[cnt] = B.point(base + ang + PI);
cnt++;
a[cnt] = A.point(base - ang);
b[cnt] = B.point(base - ang + PI);
cnt++;
}
return cnt;
}
点在圆内,圆周上不算。
bool isInCircle(Point P, Circle C)
{
return dcmp(Length(P-C.c) - C.r) < 0;
}
判断点是否在圆内,包括边界
bool InCircle(Point P, Circle C)
{
return dcmp(Length(C.c - P) - C.r) <= 0;
}
求圆心为O(0,0),半径为r,点AO与BO之间组成的夹角组成的圆弧的长度
double arcArea(Point A, Point B, double r)
{
return r * r * fabs(Angle(A, B))/ 2;
}
球面相关问题。
经纬度转换为空间坐标。经线圈(-180度-180度)和纬线圈(0度-360度)的概念相信大家并不陌生,但如何把经纬度转化为对于球心的空间左边呢?可以先算出z坐标:用半径乘以sin(纬度)就可以了,然后乘以cos(纬度),投影到x,y平面,在按照平面上的方法乘以经度的正余弦。
lat为纬度角度, lng为经度角度
x=rcos(lat)cos(lng)
y=rcos(lat)sin(lng)
z=rsin(lat)
角度转换成弧度。 此为逆时针的角度,若为顺时针,加上负号。
double torad(double deg)
{
return deg / 180.0 * PI;
}
经纬度转换为空间坐标
void get_coord(double R, double lat, double lng, double &x, double &y, double &z)
{
lat = torad(lat);
lng = torad(lng);
z = R * sin(lat);
x = R * cos(lat) * cos(lng);
y = R * cos(lat) * sin(lng);
}
球面距离:已知球面两点,如何求出它们的最短路?注意,只能沿着球面走,不能穿过球的内部。从表面走的话,最近的路径是走圆弧,具体来说是走大圆圆弧。用一个穿过球心的平面截这个球,截面就是一个大圆。怎么求大圆弧长呢?你无需想象那个大圆的空间位置,而可以把它们想象成一个平面问题:求半径r,弦长为d的圆弧长度。
圆心角为2arcsin(d/2r),因此弧长为2arcsin(d/2r)r. (弧长=圆心角(弧度)* 半径)