zoukankan      html  css  js  c++  java
  • 二维计算几何复习

    二维计算几何

     

    声明:

    由于本人较弱,并不能保证以下内容的100%正确

    欢迎大佬来挑错

    基本定义

     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){
    2     return Vec(v.x*cos(rad)-v.y*sin(rad),v.x*sin(rad)+v.y*cos(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 }

    二维凸包

    注意

    输入不能有重复点

    精度高时使用dcmp比较

     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;
    }
    

      

     
  • 相关阅读:
    cannot resolve symbol 'XXX'
    jwt单点登入
    空3
    Hibernate持久化,生命周期
    Hibernate主键生成策略
    Hibernate常用api以及增删改查
    Hibernate配置流程
    Hibernate定义
    Git总结
    spring整合MQ
  • 原文地址:https://www.cnblogs.com/zzyer/p/8953874.html
Copyright © 2011-2022 走看看