struct Point { double x, y; Point() {} Point(double a, double b) { x = a, y = b; } bool operator == (Point a) { return sgn(x-a.x) == 0 && sgn(y-a.y) == 0; } bool operator < (Point a) { //最左下角 return sgn(y-a.y) == 0 ? sgn(x-a.x)<0 : y<a.y; } Point operator - (Point a) { return Point(x-a.x, y-a.y); //向量减 } Point operator + (Point a) { return Point(x+a.x, y+a.y); //向量加 } double operator * (Point a) { //点乘(内积) return x*a.x + y*a.y; } double operator ^ (Point a) { //叉乘(外积) return x*a.y - y*a.x; } Point operator * (double k) { return Point(x*k, y*k); } Point operator / (double k) { return Point(x/k, y/k); } Point trunc(double r) { //化为长度为 r 的向量 double l = len(); if(!sgn(l)) return *this; r /= l; return Point(x*r, y*r); } Point rotleft() { //原向量逆时针旋转90度 return Point(-y, x); } Point rotright() { //原向量顺时针旋转90度 return Point(y, -x); } Point rotate(Point p, double rad) { //原向量绕着 p 点逆时针旋转 rad 度 Point v = (*this) - p; double c = cos(rad), s = sin(rad); return Point(p.x + v.x*c - v.y*s, p.y + v.x*s + v.y*c); } double len() { //向量长度 return hypot(x, y); } double len2() { //向量长度平方 return x*x + y*y; } double distance(Point p) { //两点之间距离 return hypot(x-p.x, y-p.y); } double rad(Point a, Point b) { //三个点:p, a, b //计算 pa 与 pb 所成的夹角 Point p = *this; return fabs(atan2(fabs((a-p) ^ (b-p)), (a-p) * (b-p))); } }; struct Line { Point s, e; Line() {}; Line(Point _s, Point _e) { // 两点确定直线 s = _s, e = _e; } Line(Point p, double angle) { //根据点 p 和倾斜角 angle 确定直线, 0 <= angle <= pi s = p; if(sgn(angle - pi/2) == 0) { e = (s + Point(0, 1)); } else { e = (s + Point(1, tan(angle))); } } Line(double a, double b, double c) { //根据 ax+by+c = 0 确定直线 if(sgn(a) == 0) { s = Point(0, -c/b); e = Point(1, -c/b); } else if(sgn(b) == 0) { s = Point(-c/a, 0); e = Point(-c/a, 1); } else { s = Point(0, -c/b); e = Point(1, (-c-a)/b); } } void adjust() { if(e < s) swap(s, e); } bool operator == (Line v) { return (s == v.s) && (e == v.e); } double length() { return s.distance(e); } double angle() { // 直线倾斜角 double k = atan2(e.y-s.y, e.x-s.x); if(sgn(k) < 0) k += pi; if(sgn(k-pi) == 0) k -= pi; return k; } int relation(Point p) { /* 判断直线和点关系 this line, p point 1 在左侧 2 在右侧 3 在直线上 */ int c = sgn((p-s) ^ (e-s)); if(c < 0) return 1; else if(c > 0) return 2; else return 3; } bool pointonseg(Point p) { /* 判断点和线段关系 this seg, p point 0 不在线段上 1 在线段上 */ return sgn((p-s) ^ (e-s)) == 0 && sgn((p-s) * (p-e)) <= 0; } bool parallel(Line v) { //判断两向量平行 return sgn((e-s) ^ (v.e-v.s)) == 0; } int segcrossseg(Line v) { /* 两线段关系 2 规范相交 1 不规范相交 0 不相交 */ int d1 = sgn((e-s) ^ (v.s-s)); int d2 = sgn((e-s) ^ (v.e-s)); int d3 = sgn((v.e-v.s) ^ (s-v.s)); int d4 = sgn((v.e-v.s) ^ (e-v.s)); if((d1 ^ d2) == -2 && (d3 ^ d4) == -2) return 2; return (d1 == 0 && sgn((v.s-s) * (v.s-e)) <= 0) || (d2 == 0 && sgn((v.e-s) * (v.e-e)) <= 0) || (d3 == 0 && sgn((s-v.s) * (s-v.e)) <= 0) || (d4 == 0 && sgn((e-v.s) * (e-v.e)) <= 0); } int linecrossseg(Line v) { /* 判断直线和线段关系 this line, v seg 2 规范相交 1 非规范相交 0 不相交 */ int d1 = sgn((e-s) ^ (v.s-s)); int d2 = sgn((e-s) ^ (v.e-s)); if((d1 ^ d2) == -2) return 2; else return (d1 == 0 || d2 == 0); } int linecrossline(Line v) { /* 两直线关系 0 平行 1 重合 2 相交 */ if((*this).parallel(v)) return v.relation(s) == 3; return 2; } Point crosspoint(Line v) { //求两直线的交点 //要求两直线不平行，不重合 double a1 = (v.e-v.s) ^ (s-v.s); double a2 = (v.e-v.s) ^ (e-v.s); return Point((s.x*a2-e.x*a1), (s.y*a2-e.y*a1))/(a2-a1); } double dispointtoline(Point p) { //求点到直线的距离 return fabs((p-s)^(e-s))/length(); } double dispointtoseg(Point p) { if(sgn((p-s)*(e-s))<0 || sgn((p-e)*(s-e)) < 0) return min(p.distance(s), p.distance(e)); return dispointtoline(p); } double dissegtoseg(Line v) { /* 返回线段到线段的距离 前提是两线段不相交，相交为0 */ return min(min(dispointtoseg(v.s), dispointtoseg(v.e)), min(v.dispointtoseg(s), v.dispointtoseg(e))); } Point lineprog(Point p) { //返回点 p 在直线上的投影点 return s + (((e-s)*((e-s)*(p-s))) / ((e-s).len2())); } Point symmetrypoint(Point p) { //返回点 p 关于直线p的对称点 Point q = lineprog(p); return Point(2*q.x-p.x, 2*q.y-p.y); } }; struct Circle { Point p; //圆心 double r; //半径 Circle() {} Circle(Point _p, double _r) { p = _p, r = _r; } Circle(double x, double y, double _r) { p = Point(x, y), r = _r; } Circle(Point a, Point b, Point c) { //求三角形的外接圆 Line u = Line((a+b)/2, ((a+b)/2) + ((b-a).rotleft())); Line v = Line((b+c)/2, ((b+c)/2) + ((c-b).rotleft())); p = u.crosspoint(v); r = p.distance(a); } Circle(Point a, Point b, Point c, bool t) { //求三角形的内接圆 //bool t用来区别三角形外接圆的函数 Line u, v; double m = atan2(b.y-a.y, b.x-a.x); double n = atan2(c.y-a.y, c.x-a.x); u.s = a; u.e = u.s + Point(cos((n+m)/2), sin((n+m)/2)); v.s = b; m = atan2(a.y-b.y, a.x-b.x); n = atan2(c.y-b.y, c.x-b.x); v.e = v.s + Point(cos((n+m)/2), sin((n+m)/2)); p = u.crosspoint(v); r = Line(a, b).dispointtoseg(p); } bool operator == (Circle v) { return p==v.p && sgn(r-v.r)==0; } bool operator < (Circle v) { return (p<v.p) || (p==v.p && sgn(r-v.r)<0); } double area() { //面积 return pi*r*r; } double cirperimeter() { //周长 return 2*pi*r; } int relation(Point b) { /* 点和圆的关系 0 圆外 1 圆上 2 圆内 */ double dst = b.distance(p); if(sgn(dst-r) < 0) return 2; else if(sgn(dst-r) == 0) return 1; else return 0; } int relationseg(Line v) { /* 线段和圆的关系 0 圆外 1 圆上 2 圆内 */ double dst = v.dispointtoseg(p); if(sgn(dst-r) < 0) return 2; else if(sgn(dst-r) == 0) return 1; else return 0; } int relationline(Line v) { /* 直线和圆的关系 0 圆外 1 圆上 2 圆内 */ double dst = v.dispointtoline(p); if(sgn(dst-r) < 0) return 2; else if(sgn(dst-r) == 0) return 1; else return 0; } int relationcircle(Circle v) { /* 两圆的关系 5 相离 4 外切 3 相交 2 内切 1 内含 */ double d = p.distance(v.p); if(sgn(d-r-v.r) > 0) return 5; if(sgn(d-r-v.r) == 0) return 4; double l = fabs(r-v.r); if(sgn(d-r-v.r) < 0 && sgn(d-l) > 0) return 3; if(sgn(d-l) == 0) return 2; if(sgn(d-l) < 0) return 1; } int pointcrosscircle(Circle v, Point &p1, Point &p2) { //求两个圆的交点，返回交点数 int rel = relationcircle(v); if(rel == 1 || rel == 5) return 0; double d = p.distance(v.p); double l = (d*d + r*r - v.r*v.r) / (2*d); double h = sqrt(r*r-l*l); Point tmp = p + (v.p-p).trunc(l); p1 = tmp + ((v.p-p).rotleft().trunc(h)); p2 = tmp + ((v.p-p).rotright().trunc(h)); if(rel == 2 || rel == 4) return 1; return 2; } int pointcrossline(Line v, Point &p1, Point &p2) { //求直线和圆的交点，返回交点数 if(!(*this).relationline(v)) return 0; Point a = v.lineprog(p); double d = v.dispointtoline(p); d = sqrt(r*r-d*d); if(sgn(d) == 0) { p1 = a; p2 = a; return 1; } p1 = a + (v.e-v.s).trunc(d); p2 = a - (v.e-v.s).trunc(d); return 2; } int getcircle(Point a, Point b, double r1, Circle &c1, Circle &c2) { //过a, b两点, 半径为 r1 的圆, 返回圆数 Circle x(a, r1), y(b, r1); int t = x.pointcrosscircle(y, c1.p, c2.p); if(!t) return 0; c1.r = c2.r = r1; return t; } int getcircle(Line u, Point q, double r1, Circle &c1, Circle &c2) { //与直线 u 相切, 过点 q, 半径为 r1 的圆,返回圆数 double dis = u.dispointtoline(q); if(sgn(dis-2*r1) > 0) return 0; if(sgn(dis) == 0) { c1.p = q + ((u.e-u.s).rotleft().trunc(r1)); c2.p = q + ((u.e-u.s).rotright().trunc(r1)); c1.r = c2.r = r1; return 2; } Line u1 = Line((u.s + (u.e-u.s).rotleft().trunc(r1)), (u.e + (u.e-u.s).rotleft().trunc(r1))); Line u2 = Line((u.s + (u.e-u.s).rotright().trunc(r1)), (u.e + (u.e-u.s).rotright().trunc(r1))); Circle cc = Circle(q, r1); Point p1, p2; if(!cc.pointcrossline(u1, p1, p2)) cc.pointcrossline(u2, p1, p2); c1 = Circle(p1, r1); if(p1 == p2) { c2 = c1; return 1; } c2 = Circle(p2, r1); return 2; } int getcircle(Line u, Line v, double r1, Circle &c1, Circle &c2, Circle &c3, Circle &c4) { //同时与直线u, v相切, 半径为 r1 的圆 if(u.parallel(v))return 0;//两直线平行 Line u1 = Line(u.s + (u.e-u.s).rotleft().trunc(r1),u.e + (u.e-u.s).rotleft().trunc(r1)); Line u2 = Line(u.s + (u.e-u.s).rotright().trunc(r1),u.e + (u.e-u.s).rotright().trunc(r1)); Line v1 = Line(v.s + (v.e-v.s).rotleft().trunc(r1),v.e + (v.e-v.s).rotleft().trunc(r1)); Line v2 = Line(v.s + (v.e-v.s).rotright().trunc(r1),v.e + (v.e-v.s).rotright().trunc(r1)); c1.r = c2.r = c3.r = c4.r = r1; c1.p = u1.crosspoint(v1); c2.p = u1.crosspoint(v2); c3.p = u2.crosspoint(v1); c4.p = u2.crosspoint(v2); return 4; } int getcircle(Circle cx, Circle cy, double r1, Circle &c1, Circle &c2) { //同时与不相交圆cx, cy相切，半径为 r1 的圆 Circle x(cx.p, r1+cx.r), y(cy.p, r1+cy.r); int t = x.pointcrosscircle(y, c1.p, c2.p); if(!t) return 0; c1.r = c2.r = r1; return t; } int tangentline(Point q, Line &u, Line &v) { //过一点做圆的切线,求切线 int x = relation(q); if(x == 2)return 0; if(x == 1) { u = Line(q, q + (q-p).rotleft()); v = u; return 1; } double d = p.distance(q); double l = r*r/d; double h = sqrt(r*r-l*l); u = Line(q, p + ((q-p).trunc(l) + (q-p).rotleft().trunc(h))); v = Line(q, p + ((q-p).trunc(l) + (q-p).rotright().trunc(h))); return 2; } double areacircle(Circle v) { //两相交圆的面积 int rel = relationcircle(v); if(rel >= 4) return 0.0; if(rel <= 2) return min(area(), v.area()); double d = p.distance(v.p); double hf = (r+v.r+d)/2.0; double ss = 2*sqrt(hf*(hf-r)*(hf-v.r)*(hf-d)); double a1 = acos((r*r+d*d-v.r*v.r) / (2.0*r*d)); a1 = a1*r*r; double a2 = acos((v.r*v.r+d*d-r*r) / (2.0*v.r*d)); a2 = a2*v.r*v.r; return a1+a2-ss; } double areatriangle(Point a, Point b) { //求圆和三角形 pab 的相交面积 if(sgn((p-a)^(p-b)) == 0) return 0.0; Point q[5]; int len = 0; q[len++] = a; Line l(a,b); Point p1,p2; if(pointcrossline(l,q[1],q[2]) == 2) { if(sgn((a-q[1])*(b-q[1])) < 0) q[len++] = q[1]; if(sgn((a-q[2])*(b-q[2])) < 0) q[len++] = q[2]; } q[len++] = b; if(len == 4 && sgn((q[0]-q[1])*(q[2]-q[1])) > 0) swap(q[1], q[2]); double res = 0; for(int i = 0; i < len-1; i++) { if(relation(q[i]) == 0 || relation(q[i+1]) == 0) { double arg = p.rad(q[i], q[i+1]); res += r*r*arg/2.0; } else { res += fabs((q[i]-p)^(q[i+1]-p)) / 2.0; } } return res; } }; struct Polygon { int n; Point p[maxn]; Line l[maxn]; Polygon() { n = 0; } void add(Point q) { p[++n] = q; } void getline() { for(int i=1; i<=n; i++) { l[i] = Line(p[i], p[i%n+1]); } } struct cmp { //极角排序的比较函数 Point p; cmp(Point _p) { p = _p; } bool operator()(Point _a, Point _b) const { Point a = _a, b = _b; int d = sgn((a-p) ^ (b-p)); if(d == 0) { return sgn(a.distance(p) - b.distance(p)) < 0; } else { return d > 0; } } }; void norm() { //极角排序 int id = 1; for(int i=2; i<=n; i++) { if(p[i] < p[id]) id = i; } swap(p[id], p[1]); sort(p+1, p+1+n, cmp(p[1])); } void Graham(Polygon &convex) { //引用 //求解凸包 //convex凸包 norm(); mes(convex.p, 0); int &top = convex.n = 0; //引用 if(n == 1) { convex.p[++top] = p[1]; } else if(n == 2) { convex.p[++top] = p[1]; convex.p[++top] = p[2]; if(convex.p[1] == convex.p[2]) top--; } else { convex.p[++top] = p[1]; convex.p[++top] = p[2]; for(int i=3; i<=n; i++) { while(top>1 && sgn((convex.p[top]-convex.p[top-1]) ^ (p[i]-convex.p[top-1])) <= 0) top--; convex.p[++top] = p[i]; } if(top == 2 && convex.p[1] == convex.p[2]) top--; } } bool isconvex() { //判断多边形是不是凸的 //1是凸的 0不是凸的 bool s[3]; mes(s, false); for(int i=1; i<=n; i++) { int j = i%n + 1; int k = j%n + 1; s[sgn((p[j]-p[i]) ^ (p[k]-p[i])) + 1] = true; if(s[0] && s[2]) return false; } return true; } int relationpoint(Point q) { /* 点和多边形关系 3 点上 2 边上 1 内部 0 外部 */ for(int i=1; i<=n; i++) { if(p[i] == q) return 3; } getline(); for(int i=1; i<=n; i++) { if(l[i].pointonseg(q)) return 2; } int cnt = 0; for(int i=1; i<=n; i++) { int j = i%n + 1; int k = sgn((q-p[j]) ^ (p[i]-p[j])); int u = sgn(p[i].y - q.y); int v = sgn(p[j].y - q.y); if(k > 0 && u < 0 && v >= 0) cnt++; if(k < 0 && v < 0 && u >= 0) cnt--; } return cnt != 0; } bool getdir() { //得到方向 //1表示逆时针，0表示顺时针 double sum = 0; for(int i=1; i<=n; i++) { sum += p[i] ^ p[i%n+1]; } if(sgn(sum) > 0) return 1; return 0; } Point getbarycentre() { //得到重心 Point ans(0, 0); double area = 0; for(int i=2; i<=n-1; i++) { double tmp = (p[i]-p[1]) ^ (p[i+1]-p[0]); if(sgn(tmp) == 0) continue; area += tmp; ans.x = (p[1].x + p[i].x + p[i+1].x) / 3 * tmp; ans.y = (p[0].y + p[i].y + p[i+1].y) / 3 * tmp; } if(sgn(area)) ans = ans / area; return ans; } void convexcut(Line u, Polygon &po) { //直线切割凸多边形左侧,注意直线方向 int &top = po.n = 0; //引用 for(int i=1; i<=n; i++) { int d1 = sgn((u.e-u.s) ^ (p[i]-u.s)); int d2 = sgn((u.e-u.s) ^ (p[i%n+1]-u.s)); if(d1 >= 0) po.p[++top] = p[i]; if(d1*d2<0) po.p[++top] = u.crosspoint(Line(p[i], p[i%n+1])); } } double getperimeter() { //多边形的周长 double ans = 0.0; if(n == 1) { return ans; } else if(n == 2) { ans += p[1].distance(p[2]); } else { for(int i=1; i<=n; i++) { ans += p[i].distance(p[i%n+1]); } } return ans; } double getarea() { //多边形面积 double sum = 0; for(int i=1; i<=n; i++) { sum += (p[i] ^ p[i%n+1]); } return fabs(sum)/2.0; } double areacircle(Circle c) { //多边形和圆相交的面积 double ans = 0; for(int i=1; i<=n; i++) { int j = i%n+1; if(sgn((p[j]-c.p) ^ (p[i]-c.p)) >= 0) ans += c.areatriangle(p[i], p[j]); else ans -= c.areatriangle(p[i], p[j]); } return fabs(ans); } int relationcircle(Circle c) { /* 多边形和圆的关系 2 圆完全在多边形内 1 圆在多边形内，但是碰到了多边形边界 0 其他 */ getline(); int x = 2; if(relationpoint(c.p) != 1) return 0; for(int i=1; i<=n; i++) { if(c.relationseg(l[i]) == 2) return 0; if(c.relationseg(l[i]) == 1) x = 1; } return x; } };