题意: 求任意两多边形的交,可能是凹的
思路:把顶点数为n,m的多边形切割成n-2个三角形和m-2个三角形.每两两求面积交
如果两面积方向相同: 都为正,代表凸重叠；都为负,代表是两个凹重叠； 总面积减去交面积
面积方向不相同就相加
#include<cstdio>
#include<vector>
#include<cmath>
#include<string>
#include<string.h>
#include<iostream>
#include<algorithm>
#define PI acos(-1.0)
#define pb push_back
#define F first
#define S second
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
const int N=500;
const int MOD=1e9+7;
const double eps=1e-10;
template <class T>
bool sf(T &ret){ //Faster Input
char c; int sgn; T bit=0.1;
if(c=getchar(),c==EOF) return 0;
while(c!='-'&&c!='.'&&(c<'0'||c>'9')) c=getchar();
sgn=(c=='-')?-1:1;
ret=(c=='-')?0:(c-'0');
while(c=getchar(),c>='0'&&c<='9') ret=ret*10+(c-'0');
if(c==' '||c=='\n'){ ret*=sgn; return 1; }
while(c=getchar(),c>='0'&&c<='9') ret+=(c-'0')*bit,bit/=10;
ret*=sgn;
return 1;
}
int sign(double x) { //三态函数，减少精度问题
return abs(x) < 1e-10 ? 0 : x < 0 ? -1 : 1;
}
struct Point { //点的定义
double x, y;
Point(double x=0.0, double y=0.0) : x(x), y(y) {}
Point operator + (const Point &rhs) const { //向量加法
return Point(x + rhs.x, y + rhs.y);
}
Point operator - (const Point &rhs) const { //向量减法
return Point(x - rhs.x, y - rhs.y);
}
Point operator * (double p) const { //向量乘以标量
return Point(x * p, y * p);
}
Point operator / (double p) const { //向量除以标量
return Point(x / p, y / p);
}
bool operator < (const Point &rhs) const { //点的坐标排序
return x < rhs.x || (x == rhs.x && y < rhs.y);
}
bool operator == (const Point &rhs) const { //判断同一个点
return sign(x - rhs.x) == 0 && sign(y - rhs.y) == 0;
}
};
typedef Point Vector;
double polar_angle(Vector A) { //向量极角
return atan2(A.y, A.x);
}
struct Line { //有向直线的定义，左边就是对应的半平面
Point p,s; //直线上任意一点
Vector v; //方向向量
double ang; //极角，即从x正方向旋转到向量v所需的角（弧度）
Line() {}
Line(Point p, Point s, Vector v) : p(p), s(s) ,v(v) {
ang = polar_angle(v);
}
bool operator < (const Line &rhs) const {
return ang < rhs.ang;
}
Point point(double a) {
return p + v * a;
}
};
double cross(Vector A, Vector B) { //向量叉积，两向量组成的三角形的有向面积的两倍，可判断两向量的方向
return A.x * B.y - A.y * B.x;
}
bool point_on_left(Point p, Line L) { //判断点在直线左侧
return cross(L.v, p-L.p) > 0;
}
Point line_line_inter(Point p, Vector V, Point q, Vector W) { //两直线交点，参数方程
Vector U = p - q;
double t = cross(W, U) / cross(V, W); //两直线有唯一交点的前提：cross(V, W)非0
return p + V * t;
}
//Point p[N],p1[N];
vector<Point> p,p1;
Line Q[N];
typedef vector<Point> polygon;
bool HPIcmp(Line a,Line b){
if(fabs(a.ang - b.ang) > 1e-7)return a.ang < b.ang; //斜率排序
//斜率相同
return sign(cross((a.p - b.p),(b.s - b.p))) < 0;
}
polygon half_plane_inter(vector<Line> line){
int n=(int)line.size();
int tot=n;
sort(line.begin(),line.end(),HPIcmp);
polygon ans;
tot = 1;
for(int i = 1;i < n;i++)
if(fabs(line[i].ang - line[i-1].ang) > 1e-7) //去掉斜率重复的
line[tot++] = line[i];
int head = 0, tail = 1;
Q[0] = line[0];
Q[1] = line[1];
for(int i = 2; i < tot; i++){ //如果双端队列底端或顶端两个半平面的交点在当前半平面之外，则删除
if(fabs(cross((Q[tail].s-Q[tail].p),(Q[tail-1].s-Q[tail-1].p))) <1e-7 ||
fabs(cross((Q[head].s-Q[head].p),(Q[head+1].s-Q[head+1].p))) <1e-7)
return ans;
while(head < tail && sign(cross(line_line_inter(Q[tail].p,Q[tail].v,Q[tail-1].p,Q[tail-1].v)-line[i].p,line[i].s-line[i].p))>0)
tail--;
while(head < tail && sign(cross(line_line_inter(Q[head].p,Q[head].v,Q[head+1].p,Q[head+1].v)-line[i].p,line[i].s-line[i].p))>0)
head++;
Q[++tail] = line[i];
}
while(head < tail && sign(cross(line_line_inter(Q[tail].p,Q[tail].v,Q[tail-1].p,Q[tail-1].v)-Q[head].p,Q[head].s-Q[head].p))>0)
tail--;
if(tail <= head + 1) return ans;
//保存点到ans
for(int i = head; i < tail; i++)
ans.pb(line_line_inter(Q[i].p,Q[i].v,Q[i+1].p,Q[i+1].v));
if(head < tail - 1)
ans.pb(line_line_inter(Q[head].p,Q[head].v,Q[tail].p,Q[tail].v));
return ans;
}
int ks,n,m;
double CulArea(polygon &t){
double ans=0;
for(int i=1;i<(int)t.size()-1;i++) ans+=cross(t[i]-t[0],t[i+1]-t[0])/2.0;
return ans;
}
int main(void){
while(scanf("%d%d",&n,&m)==2){
p.resize(n);
p1.resize(m);
polygon t;
double area1=0,area2=0;
for(int i=0;i<n;i++) scanf("%lf%lf",&p[i].x,&p[i].y),t.pb(p[i]);
if(sign(area1=CulArea(t))<0) reverse(p.begin(),p.end()),area1=-area1;
t.clear();
for(int i=0;i<n;i++) scanf("%lf%lf",&p1[i].x,&p1[i].y),t.pb(p1[i]);
if(sign(area2=CulArea(t))<0) reverse(p1.begin(),p1.end()),area2=-area2;
vector<Line> L;
double res=area1+area2;
for(int i=0;i<n;i++){
Line t;
t={p[i],p[(i+1)%n],p[(i+1)%n]-p[i]};
L.pb(t);
}
for(int i=0;i<n;i++){
Line t;
t={p1[i],p1[(i+1)%n],p1[(i+1)%n]-p1[i]};
L.pb(t);
}
for(int i=1;i<n-1;i++){
for(int j=1;j<m-1;j++){
// cout <<i <<" "<<j<<endl;
vector<Point> s1;
vector<Point> s2;
s1.pb(p[0]);
s1.pb(p[i]);
s1.pb(p[i+1]);
s2.pb(p1[0]);
s2.pb(p1[j]);
s2.pb(p1[j+1]);
double ss1=CulArea(s1);
int f1;
if(sign(ss1)<=0) f1=0,reverse(s1.begin(),s1.end());
else f1=1;
double ss2=CulArea(s2);
int f2;
if(sign(ss2)<=0) f2=0,reverse(s2.begin(),s2.end());
else f2=1;
vector<Line> L;
for(int i=0;i<3;i++){
Line t;
t={s1[i],s1[(i+1)%3],s1[(i+1)%3]-s1[i]};
L.pb(t);
}
for(int i=0;i<3;i++){
Line t;
t=Line(s2[i],s2[(i+1)%3],s2[(i+1)%3]-s2[i]);
L.pb(t);
}
polygon t=half_plane_inter(L);
double nowarea=CulArea(t);
if(f1==f2) res-=nowarea;
else res+=nowarea;
}
}
printf("%.2f\n",res);
}
return 0;
}
京公网安备 11010502036488号