二维计算几何复习

 1 struct Point{
2     double x,y;
3     Point(){
4         x=y=0.0;
5     }
6     Point(double xx,double yy){
7         x=xx;y=yy;
8     }
9 };
10 typedef Point Vec;
11
12 Vec operator + (Vec v1,Vec v2){
13     return Vec(v1.x+v2.x,v1.y+v2.y);
14 }
15 Vec operator - (Point p1,Point p2){
16     return Vec(p1.x-p2.x,p1.y-p2.y);
17 }
18 Vec operator * (Vec v,double k){
19     return Vec(v.x*k,v.y*k);
20 }
21 Vec operator / (Vec v,double k){
22     return Vec(v.x/k,v.y/k);
23 }
24
25 bool operator < (const Point &a,const Point &b){
26     if(a.x==b.x)return a.y<b.y;
27     else return a.x<b.x;
28 }
29
30 double Myabs(double x){
31     if(x<0)return -x;
32     else return x;
33 }
34
35 const double eps=1e-10;
36 int dcmp(double x){
37     if(Myabs(x)<eps)return 0;
38     if(x>0)return 1;
39     else return -1;
40 }
41
42 bool operator == (const Point &a,const Point &b){
43     return dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)==0;
44 }

1 double Dt(Vec v1,Vec v2){
2     return v1.x*v2.x+v1.y*v2.y;
3 }
4 double Getlen(Vec v){
5     return sqrt(Dt(v,v));
6 }
7 double Getang(Vec v1,Vec v2){
8     return acos(Dt(v1,v2)/Getlen(v1)/Getlen(v2));
9 }

1 double Crs(Vec v1,Vec v2){
2     return v1.x*v2.y-v2.x*v1.y;
3 }
4 double GetS2(Point A,Point B,Point C){
5     return Crs(B-A,C-A);
6 }

1 Vec Rotate(Vec v,double rad){
3 }
4 Vec Getn(Vec v){
5     double L=Getlen(v);
6     return Vec(-v.y/L,v.x/L);
7 }

 1 Point GetLineIntersection(Point P,Vec v,Point Q,Vec w){
2     Vec u=P-Q;
3     double t=Crs(w,u)/Crs(v,w);
4     return P+v*t;
5 }
6 double DistanceToLine(Point P,Point A,Point B){
7     Vec v1=B-A,v2=P-A;
8     return Myabs(Crs(v1,v2))/Getlen(v1);
9 }
10 Point GetLineProjection(Point P,Point A,Point B){
11     Vec v=B-A;
12     double t=Dt(v,P-A)/Dt(v,v);
13     return A+v*t;
14 }
15
16
17 bool SegmentProperIntersection(Point a1,Point a2,Point b1,Point b2){
18     double c1=Crs(a2-a1,b1-a1),c2=Crs(a2-a1,b2-a1);
19     double c3=Crs(b2-b1,a1-b1),c4=Crs(b2-b1,a2-b1);
20     return dcmp(c1)*dcmp(c2)<0&&dcmp(c3)*dcmp(c4)<0;
21 }
22 bool OnSegment(Point p,Point a1,Point a2){
23     return dcmp(Crs(a1-p,a2-p))==0&&dcmp(Dt(a1-p,a2-p))<0;
24 }
25
26 double PolygonArea(Point* p,int n){
27     double area=0;
28     for(int i=1;i<n-1;++i){
29         area+=Crs(p[i]-p[0],p[i+1]-p[0]);
30     }
31     return Myabs(area)/2;
32 }

 1 int isPointInPolygon(Point p,vector<Point>poly){
2     int wn=0;
3     int n=poly.size();
4     for(int i=0;i<n;++i){
5         if(OnSegment(p,poly[i],poly[(i+1)%n]))return -1;
6         int k=dcmp(Crs(poly[(i+1)%n]-poly[i],p-poly[i]));
7         int d1=dcmp(poly[i].y-p.y);
8         int d2=dcmp(poly[(i+1)%n].y-p.y);
9         if(k>0 && d1<=0 &&d2>0)wn++;
10         if(k<0 && d2<=0 &&d1>0)wn--;
11     }
12     if(wn!=0)return 1;
13     else return 0;
14 }

 1 int ConvexHull(Point *p,int n,Point *ch){
2     sort(p,p+n);
3     int m=0;
4     for(int i=0;i<n;++i){
5         while(m>1&&Crs(ch[m-1]-ch[m-2],p[i]-ch[m-2])<=0)--m;
6         ch[m++]=p[i];
7     }
8     int k=m;
9     for(int i=n-2;i>=0;--i){
10         while(m>k&&Crs(ch[m-1]-ch[m-2],p[i]-ch[m-2])<=0)--m;
11         ch[m++]=p[i];
12     }
13     if(n>1)m--;
14     return m;
15 }

 1 double RotatingCaliper(Point* ch,int m){
2     int q=1;
3     ch[m]=ch[0];
4     double ans=0;
5     for(int p=0;p<m;++p){
6         while(Myabs(Crs(ch[q+1]-ch[p+1],ch[p]-ch[p+1]))>Myabs(Crs(ch[q]-ch[p+1],ch[p]-ch[p+1])))q=(q+1)%m;
7         ans=max(ans,max(Getlen(ch[p]-ch[q]),Getlen(ch[p+1]-ch[q+1])));
8         ans=max(ans,max(Getlen(ch[p]-ch[q+1]),Getlen(ch[p+1]-ch[q])));
9     }
10     return ans;
11 }

O(n2)

 1 typedef vector<Point> Polygon;
2 Polygon CutPolygon(Polygon poly,Point A,Point B){
3     Polygon newpoly;
4     int n=poly.size();
5     for(int i=0;i<n;++i){
6         Point C=poly[i];
7         Point D=poly[(i+1)%n];
8         if(dcmp(Crs(B-A,C-A))>=0)newpoly.push_back(C);
9         if(dcmp(Crs(B-A,C-D))!=0){
10             Point ip=GetLineIntersection(A,B-A,C,D-C);
11             if(OnSegment(ip,C,D))newpoly.push_back(ip);
12         }
13     }
14     return newpoly;
15 }

O(nlogn)

 1 bool OnLeft(Line L,Point p){
2     return Crs(L.v,p-L.P)>0;
3 }
4 Point GetLineIntersection(Line a,Line b){
5     Vec u=a.P-b.P;
6     double t=Crs(b.v,u)/Crs(a.v,b.v);
7     return a.P+a.v*t;
8 }
9
10 Point p[maxn];
11 Line q[maxn];
12 int HalfPlaneIntersection(Line *L,int n,Point *poly){
13     sort(L,L+n);
14
15     int h,t;
16     q[h=t=1]=L[0];
17     for(int i=1;i<n;++i){
18         while((h<t)&&(!OnLeft(L[i],p[t-1])))--t;
19         while((h<t)&&(!OnLeft(L[i],p[h])))++h;
20         q[++t]=L[i];
21         if(dcmp(Crs(q[t].v,q[t-1].v))==0){
22             --t;
23             if(OnLeft(q[t],L[i].P))q[t]=L[i];
24         }
25         if(h<t)p[t-1]=GetLineIntersection(q[t-1],q[t]);
26     }
27     while((h<t)&&(!OnLeft(q[h],p[t-1])))--t;
28     if(t-h<=1)return 0;
29     p[t]=GetLineIntersection(q[t],q[h]);
30
31     int m=0;
32     for(int i=h;i<=t;++i)poly[m++]=p[i];
33     return m;
34 }

bool OnLeft(Line L,Point p){
return Crs(L.v,p-L.P)>0;
}
Point GetLineIntersection(Line a,Line b){
Vec u=a.P-b.P;
double t=Crs(b.v,u)/Crs(a.v,b.v);
return a.P+a.v*t;
}

Point p[maxn];
Line q[maxn];
int HalfPlaneIntersection(Line *L,int n,Point *poly){
sort(L,L+n);

int h,t;
q[h=t=1]=L[0];
for(int i=1;i<n;++i){
while((h<t)&&(!OnLeft(L[i],p[t-1])))--t;
while((h<t)&&(!OnLeft(L[i],p[h])))++h;
q[++t]=L[i];
if(dcmp(Crs(q[t].v,q[t-1].v))==0){
--t;
if(OnLeft(q[t],L[i].P))q[t]=L[i];
}
if(h<t)p[t-1]=GetLineIntersection(q[t-1],q[t]);
}
while((h<t)&&(!OnLeft(q[h],p[t-1])))--t;
if(t-h<=1)return 0;
p[t]=GetLineIntersection(q[t],q[h]);

int m=0;
for(int i=h;i<=t;++i)poly[m++]=p[i];
return m;
}



posted @ 2018-04-26 17:27  ws_zzy  阅读(218)  评论(0编辑  收藏  举报