题意:逆时针给出\(n(n<=10)\)个凸多边形的顶点坐标,求它们交的面积。
其他模板题推荐
还有一些前置知识。两向量\((x_1,y_1),(x_2,y_2)\)的叉乘为\(x_1y_2-x_2y_1\),结果为正说明向量\((x_2,y_2)\)在向量\((x_1,y_1)\)逆时针方向,结果为负则在顺时针方向。求两直线交点的公式如下图所示:
代码与原博客稍有不同。
#include<bits/stdc++.h>
#define ll long long
using namespace std;
inline int read() {
char ch = getchar(); int x = 0, f = 1;
while (ch < '0' || ch>'9') { if (ch == '-') f = -1; ch = getchar(); }
while ('0' <= ch && ch <= '9') { x = x * 10 + ch - '0'; ch = getchar(); }
return x * f;
}
const int N = 1e3 + 10;
const double EPS = 1e-5;
int T, tot;
struct node {//一个点,两个坐标
double x, y;
};
node p[15][55];
node ppx[N];
node operator - (node a, node b) {//两个点相减,得到向量ab
node t;
t.x = a.x - b.x;
t.y = a.y - b.y;
return t;
}
double operator ^ (node a, node b) { return a.x * b.y - a.y * b.x; }
//这里的ab应该是向量,然后求叉乘
struct Line {//一个向量,起点s,终点e
node s, e;
};
Line L[N], que[N];
double getAngle(node a) { return atan2(a.y, a.x); }//这里的a应该是向量
double getAngle(Line a) { return atan2(a.e.y - a.s.y, a.e.x - a.s.x); }//求向量a的arctan值
bool cmp(Line a, Line b) {//对所有向量进行极角排序
double A = getAngle(a);
double B = getAngle(b);
return A < B;
}
node getIntersectPoint(Line a, Line b) {//求两直线交点
double a1 = a.s.y - a.e.y, b1 = a.s.x - a.e.x, c1 = 1.0 * a.s.x * a.e.y - 1.0 * a.e.x * a.s.y;
double a2 = b.s.y - b.e.y, b2 = b.s.x - b.e.x, c2 = 1.0 * b.s.x * b.e.y - 1.0 * b.e.x * b.s.y;
node t;
t.x = (1.0 * c1 * b2 - 1.0 * c2 * b1) / (1.0 * a2 * b1 - 1.0 * a1 * b2);
t.y = (1.0 * c1 * a2 - 1.0 * c2 * a1) / (1.0 * a2 * b1 - 1.0 * a1 * b2);
return t;
}
//判断向量a是否在向量bc交点的右侧
bool onRight(Line a, Line b, Line c) {
node o = getIntersectPoint(b, c);
if (((a.e - a.s) ^ (o - a.s)) < 0) return true;//可以自己画图a.s a.e o三个点
return false;
}
double HalfPlaneIntersection() {
sort(L + 1, L + tot + 1, cmp);
int head = 1, tail = 1;
que[1] = L[1];//构造双端队列
for (int i = 2; i <= tot; i++) {
while (head < tail && onRight(L[i], que[tail], que[tail - 1])) tail--;
while (head < tail && onRight(L[i], que[head], que[head + 1])) head++;
que[++tail] = L[i];
//极角相同的向量,保留靠左的那一个
if (fabs(getAngle(que[tail]) - getAngle(que[tail - 1])) < EPS) {
tail--;
if (((que[tail].e - que[tail].s) ^ (L[i].e - que[tail].s)) > EPS)que[tail] = L[i];
}
}
while (head < tail && onRight(que[head], que[tail], que[tail - 1])) tail--;
while (head < tail && onRight(que[tail], que[head], que[head + 1])) head++;
if (tail - head < 2) return 0;//剩下的直线无法构成多边形
double ans = 0;
int tot_jd = 0;
for (int i = head; i < tail; ++i) {
ppx[++tot_jd] = getIntersectPoint(que[i], que[i + 1]);
}
ppx[++tot_jd] = getIntersectPoint(que[tail], que[head]);
for (int i = 2; i < tot_jd; ++i) {
double x1 = ppx[1].x, y1 = ppx[1].y;
double x2 = ppx[i].x, y2 = ppx[i].y;
double x3 = ppx[i + 1].x, y3 = ppx[i + 1].y;
ans = ans + (x2 - x1) * (y3 - y1) - (y2 - y1) * (x3 - x1);
}
return ans / 2.0;
}
int main() {
int T; cin >> T;
for (int t = 1; t <= T; ++t) {
int n; cin >> n;
for (int i = 1; i <= n; ++i) {
cin >> p[t][i].x;
cin >> p[t][i].y;
}
for (int i = 1; i < n; i++) {
L[++tot].s.x = p[t][i].x;
L[tot].s.y = p[t][i].y;
L[tot].e.x = p[t][i + 1].x;
L[tot].e.y = p[t][i + 1].y;
}
L[++tot].s.x = p[t][n].x;
L[tot].s.y = p[t][n].y;
L[tot].e.x = p[t][1].x;
L[tot].e.y = p[t][1].y;
}
printf("%.3lf\n", HalfPlaneIntersection());
return 0;
}