【loj6437】 【PKUSC2018】 PKUSC 計算幾何

題目大意：給你一個m個點的簡單多邊形。對於每一個點i∈[1,n]，做一個以O點爲原點且過點i的圓，求該圓在多邊形內的圓弧長度/圓長。

其中n≤200，m≤500。

 

咱們將n個點分開處理。

 

首先，咱們要判斷需處理的圓，是否被包含在多邊形內，或者圓把多邊形包含了。

咱們顯然能夠從原點出發，向上做一條x=0的射線，判斷該射線與多邊形有多少個交點。

顯然，若交點數量爲奇數個，那麼該點就在多邊形內，不然在多邊形外。

 

若圓與多邊形存在交點，咱們對多邊形的每條邊，求出其與圓，有多少個交點（0個，1個，2個，其實1個點能夠當2個去處理）

這裏就是簡單地推個式子就能夠求出，詳見代碼。

求出這些交點後，進行極角排序。

對於兩個排序後相鄰的點，咱們要判斷，由這兩個點所構成的弧，是否被多邊形包含。

有一種很暴力的思路，就是求出這兩個點在弧上的中點，而後做一條從該點出發，向上做一條平行於x軸的射線，求這條射線與多邊形的相交次數。

顯然，若交點數量爲奇數個，那麼該點就在多邊形內，不然在多邊形外（和上面判斷圓是否被包含的部分相同）。

 

對於每一個圓，最多會出現2m個交點，最後單次判斷兩點所構乘的圓弧是否被多邊形包含，須要O(m)的時間，則總時間複雜度爲O(n*m^2)。

 

下面是一些須要注意的細節：

1，上述判斷奇數偶數個點時，所作射線可能會通過多邊形兩條邊的交點。考慮到輸入的點所有都是整點，故對構成多邊形的點集偏移2eps便可。

2，因爲點要被多邊形嚴格包含，那麼在判斷時，範圍要向兩邊縮減eps。

 

思路簡單，代碼複雜，qwq....

（據說還卡精度）。。。

 

 1 #include<bits/stdc++.h>
 2 #define M 505
 3 #define DB long double 
 4 #define eps 1e-6
 5 #define PI 3.14159265358979323846
 6 using namespace std;  7 
 8 bool zero(DB x){return fabs(x)<eps;}  9 struct pt{  10  DB x,y;  11     pt(){x=y=0;}  12     pt(DB xx,DB yy){x=xx; y=yy;}  13     friend pt operator -(pt a,pt b){return pt(a.x-b.x,a.y-b.y);}  14     friend bool operator <(pt a,pt b){return atan2(a.y,a.x)<atan2(b.y,b.x);}  15     DB mo(){return x*x+y*y;}  16 };  17         
 18 struct line{  19  DB a,b,c;  20     line(){a=b=c=0;}  21     line(DB aa,DB bb,DB cc){a=aa; b=bb; c=cc;}  22  line(pt A,pt B){  23         DB x1=A.x,y1=A.y;  24         DB x2=B.x,y2=B.y;  25         a=y1-y2; b=x2-x1;  26         c=-a*x2-b*y2;  27  }  28 };  29 
 30 pt a[M],b[M],c[M]; int n,m;  31 int type[M]={0};  32 DB sita[M]={0};  33 DB solve(pt hh){  34     memset(c,0,sizeof(c));  35     memset(type,0,sizeof(type));  36     int cnt=0,ok=1;  37     DB R=hh.x*hh.x+hh.y*hh.y,r=sqrt(R);  38     for(int i=1;i<=m;i++){  39         line p=line(b[i],b[i+1]);  40         if(fabs(p.c/sqrt(p.a*p.a+p.b*p.b))-eps<r) ok=0;  41         DB xl=min(b[i].x,b[i+1].x),xr=max(b[i].x,b[i+1].x);  42         DB yl=min(b[i].y,b[i+1].y),yr=max(b[i].y,b[i+1].y);  43         xl-=eps; yl-=eps; xr+=eps; yr+=eps;  44         DB x1,x2,y1,y2; bool is0=0;  45         if(zero(p.a)){  46             y1=y2=-p.c/p.b;  47             DB delta=R-y1*y1;  48             if(delta<-eps) continue;  49             if(zero(delta)) is0=1;  50             x1=-sqrt(delta); x2=sqrt(delta);  51         }else{  52             DB A=p.b*p.b+p.a*p.a;  53             DB B=p.b*p.c*2;  54             DB C=p.c*p.c-R*p.a*p.a;  55             DB delta=B*B-4*A*C;  56             if(delta<-eps) continue;  57             if(zero(delta)) is0=1;  58             y1=(-B+sqrt(delta))/(2*A);  59             y2=(-B-sqrt(delta))/(2*A);  60             x1=(-p.b*y1-p.c)/p.a;  61             x2=(-p.b*y2-p.c)/p.a;  62  }  63         int ok1,ok2;  64         if(xl<=x1&&x1<=xr&&yl<=y1&&y1<=yr){ c[++cnt]=pt(x1,y1); if(is0) continue;}  65         if(xl<=x2&&x2<=xr&&yl<=y2&&y2<=yr) c[++cnt]=pt(x2,y2);  66  }  67     if(ok){  68         //return 1;
 69         c[cnt=1]=pt(0,0);  70     }else{  71         sort(c+1,c+cnt+1);  72         for(int i=1;i<=cnt;i++) sita[i]=atan2(c[i].y,c[i].x);  73  }  74     c[cnt+1]=c[1]; sita[cnt+1]=sita[1]+2*PI;  75     DB ans=0;  76     for(int i=1;i<=cnt;i++){  77         double SITA=(sita[i]+sita[i+1])/2;  78         pt hh=pt(r*cos(SITA)+eps*2,r*sin(SITA));  79         int cnt=0;  80         for(int j=1;j<=m;j++){  81             double l=b[j].x,r=b[j+1].x;  82             if(zero(b[j].x-b[j+1].x)) continue;  83             line p=line(b[j],b[j+1]);  84             double y=-(p.a*hh.x+p.c)/p.b;  85             if(l>r) swap(l,r);  86             l+=eps; r-=eps;  87             if(l<=hh.x&&hh.x<=r&&y+eps>hh.y) cnt++;  88  }  89         if(cnt&1) ans+=sita[i+1]-sita[i];  90  }  91     ans/=(2.*PI);  92     return ans;  93 }  94         
 95 
 96 int main(){  97 // freopen("in.txt","r",stdin);  98 // freopen("out.txt","w",stdout);
 99     scanf("%d%d",&n,&m); 100     for(int i=1;i<=n;i++) cin>>a[i].x>>a[i].y; 101     for(int i=1;i<=m;i++) cin>>b[i].x>>b[i].y; 102     b[m+1]=b[1]; 103     DB ans=0; 104     for(int i=1;i<=n;i++){ 105         ans+=solve(a[i]); 106         //printf("%.5Lf\n",ans);
107  } 108     double res=ans; 109     printf("%.5lf\n",res); 110 }