Art Gallery
http://poj.org/problem?id=1279
| Time Limit: 1000MS | Memory Limit: 10000K | |
| Total Submissions: 9684 | Accepted: 3747 |
Description
The art galleries of the new and very futuristic building of the Center for Balkan Cooperation have the form of polygons (not necessarily convex). When a big exhibition is organized, watching over all of the pictures is a big security concern. Your task is that for a given gallery to write a program which finds the surface of the area of the floor, from which each point on the walls of the gallery is visible. On the figure 1. a map of a gallery is given in some co-ordinate system. The area wanted is shaded on the figure 2.
Input
The number of tasks T that your program have to solve will be on the first row of the input file. Input data for each task start with an integer N, 5 <= N <= 1500. Each of the next N rows of the input will contain the co-ordinates of a vertex of the polygon ? two integers that fit in 16-bit integer type, separated by a single space. Following the row with the co-ordinates of the last vertex for the task comes the line with the number of vertices for the next test and so on.
Output
For each test you must write on one line the required surface - a number with exactly two digits after the decimal point (the number should be rounded to the second digit after the decimal point).
Sample Input
1 7 0 0 4 4 4 7 9 7 13 -1 8 -6 4 -4
Sample Output
80.00
要特判面积不存在的情况
1 #include<iostream> 2 #include<cstring> 3 #include<cstdio> 4 #include<vector> 5 #include<cmath> 6 #include<algorithm> 7 using namespace std; 8 const double eps=1e-8; 9 const double INF=1e20; 10 const double PI=acos(-1.0); 11 12 13 int sgn(double x){ 14 if(fabs(x)<eps) return 0; 15 if(x<0) return -1; 16 else return 1; 17 } 18 inline double sqr(double x){return x*x;} 19 struct Point{ 20 double x,y; 21 Point(){} 22 Point(double _x,double _y){ 23 x=_x; 24 y=_y; 25 } 26 void input(){ 27 scanf("%lf %lf",&x,&y); 28 } 29 void output(){ 30 printf("%.2f %.2f ",x,y); 31 } 32 bool operator == (const Point &b)const{ 33 return sgn(x-b.x) == 0 && sgn(y-b.y)== 0; 34 } 35 bool operator < (const Point &b)const{ 36 return sgn(x-b.x)==0?sgn(y-b.y)<0:x<b.x; 37 } 38 Point operator - (const Point &b)const{ 39 return Point(x-b.x,y-b.y); 40 } 41 //叉积 42 double operator ^ (const Point &b)const{ 43 return x*b.y-y*b.x; 44 } 45 //点积 46 double operator * (const Point &b)const{ 47 return x*b.x+y*b.y; 48 } 49 //返回长度 50 double len(){ 51 return hypot(x,y); 52 } 53 //返回长度的平方 54 double len2(){ 55 return x*x+y*y; 56 } 57 //返回两点的距离 58 double distance(Point p){ 59 return hypot(x-p.x,y-p.y); 60 } 61 Point operator + (const Point &b)const{ 62 return Point(x+b.x,y+b.y); 63 } 64 Point operator * (const double &k)const{ 65 return Point(x*k,y*k); 66 } 67 Point operator / (const double &k)const{ 68 return Point(x/k,y/k); 69 } 70 71 //计算pa和pb的夹角 72 //就是求这个点看a,b所成的夹角 73 ///LightOJ1202 74 double rad(Point a,Point b){ 75 Point p=*this; 76 return fabs(atan2(fabs((a-p)^(b-p)),(a-p)*(b-p))); 77 } 78 //化为长度为r的向量 79 Point trunc(double r){ 80 double l=len(); 81 if(!sgn(l)) return *this; 82 r/=l; 83 return Point(x*r,y*r); 84 } 85 //逆时针转90度 86 Point rotleft(){ 87 return Point(-y,x); 88 } 89 //顺时针转90度 90 Point rotright(){ 91 return Point(y,-x); 92 } 93 //绕着p点逆时针旋转angle 94 Point rotate(Point p,double angle){ 95 Point v=(*this) -p; 96 double c=cos(angle),s=sin(angle); 97 return Point(p.x+v.x*c-v.y*s,p.y+v.x*s+v.y*c); 98 } 99 }; 100 101 struct Line{ 102 Point s,e; 103 Line(){} 104 Line(Point _s,Point _e){ 105 s=_s; 106 e=_e; 107 } 108 bool operator==(Line v){ 109 return (s==v.s)&&(e==v.e); 110 } 111 //根据一个点和倾斜角angle确定直线，0<=angle<pi 112 Line(Point p,double angle){ 113 s=p; 114 if(sgn(angle-PI/2)==0){ 115 e=(s+Point(0,1)); 116 } 117 else{ 118 e=(s+Point(1,tan(angle))); 119 } 120 } 121 //ax+by+c=0; 122 Line(double a,double b,double c){ 123 if(sgn(a)==0){ 124 s=Point(0,-c/b); 125 e=Point(1,-c/b); 126 } 127 else if(sgn(b)==0){ 128 s=Point(-c/a,0); 129 e=Point(-c/a,1); 130 } 131 else{ 132 s=Point(0,-c/b); 133 e=Point(1,(-c-a)/b); 134 } 135 } 136 void input(){ 137 s.input(); 138 e.input(); 139 } 140 void adjust(){ 141 if(e<s) swap(s,e); 142 } 143 //求线段长度 144 double length(){ 145 return s.distance(e); 146 } 147 //返回直线倾斜角 0<=angle<pi 148 double angle(){ 149 double k=atan2(e.y-s.y,e.x-s.x); 150 if(sgn(k)<0) k+=PI; 151 if(sgn(k-PI)==0) k-=PI; 152 return k; 153 } 154 //点和直线的关系 155 //1 在左侧 156 //2 在右侧 157 //3 在直线上 158 int relation(Point p){ 159 int c=sgn((p-s)^(e-s)); 160 if(c<0) return 1; 161 else if(c>0) return 2; 162 else return 3; 163 } 164 //点在线段上的判断 165 bool pointonseg(Point p){ 166 return sgn((p-s)^(e-s))==0&&sgn((p-s)*(p-e))<=0; 167 } 168 //两向量平行(对应直线平行或重合) 169 bool parallel(Line v){ 170 return sgn((e-s)^(v.e-v.s))==0; 171 } 172 //两线段相交判断 173 //2 规范相交 174 //1 非规范相交 175 //0 不相交 176 int segcrossseg(Line v){ 177 int d1=sgn((e-s)^(v.s-s)); 178 int d2=sgn((e-s)^(v.e-s)); 179 int d3=sgn((v.e-v.s)^(s-v.s)); 180 int d4=sgn((v.e-v.s)^(e-v.s)); 181 if((d1^d2)==-2&&(d3^d4)==-2) return 2; 182 return (d1==0&&sgn((v.s-s)*(v.s-e))<=0|| 183 d2==0&&sgn((v.e-s)*(v.e-e))<=0|| 184 d3==0&&sgn((s-v.s)*(s-v.e))<=0|| 185 d4==0&&sgn((e-v.s)*(e-v.e))<=0); 186 } 187 //直线和线段相交判断 188 //-*this line -v seg 189 //2 规范相交 190 //1 非规范相交 191 //0 不相交 192 int linecrossseg(Line v){ 193 int d1=sgn((e-s)^(v.s-s)); 194 int d2=sgn((e-s)^(v.e-s)); 195 if((d1^d2)==-2) return 2; 196 return (d1==0||d2==0); 197 } 198 //两直线关系 199 //0 平行 200 //1 重合 201 //2 相交 202 int linecrossline(Line v){ 203 if((*this).parallel(v)) 204 return v.relation(s)==3; 205 return 2; 206 } 207 //求两直线的交点 208 //要保证两直线不平行或重合 209 Point crosspoint(Line v){ 210 double a1=(v.e-v.s)^(s-v.s); 211 double a2=(v.e-v.s)^(e-v.s); 212 return Point((s.x*a2-e.x*a1)/(a2-a1),(s.y*a2-e.y*a1)/(a2-a1)); 213 } 214 //点到直线的距离 215 double dispointtoline(Point p){ 216 return fabs((p-s)^(e-s))/length(); 217 } 218 //点到线段的距离 219 double dispointtoseg(Point p){ 220 if(sgn((p-s)*(e-s))<0||sgn((p-e)*(s-e))<0) 221 return min(p.distance(s),p.distance(e)); 222 return dispointtoline(p); 223 } 224 //返回线段到线段的距离 225 //前提是两线段不相交，相交距离就是0了 226 double dissegtoseg(Line v){ 227 return min(min(dispointtoseg(v.s),dispointtoseg(v.e)),min(v.dispointtoseg(s),v.dispointtoseg(e))); 228 } 229 //返回点P在直线上的投影 230 Point lineprog(Point p){ 231 return s+(((e-s)*((e-s)*(p-s)))/((e-s).len2())); 232 } 233 //返回点P关于直线的对称点 234 Point symmetrypoint(Point p){ 235 Point q=lineprog(p); 236 return Point(2*q.x-p.x,2*q.y-p.y); 237 } 238 }; 239 240 struct circle{ 241 Point p; 242 double r; 243 circle(){} 244 circle(Point _p,double _r){ 245 p=_p; 246 r=_r; 247 } 248 249 circle(double x,double y,double _r){ 250 p=Point(x,y); 251 r=_r; 252 } 253 254 circle(Point a,Point b,Point c){///三角形的外接圆 255 Line u=Line((a+b)/2,((a+b)/2)+((b-a).rotleft())); 256 Line v=Line((b+c)/2,((b+c)/2)+((c-b).rotleft())); 257 p=u.crosspoint(v); 258 r=p.distance(a); 259 } 260 261 circle(Point a,Point b,Point c,bool t){///三角形的内切圆 262 Line u,v; 263 double m=atan2(b.y-a.y,b.x-a.x),n=atan2(c.y-a.y,c.x-a.x); 264 u.s=a; 265 u.e=u.s+Point(cos((n+m)/2),sin((n+m)/2)); 266 v.s=b; 267 m=atan2(a.y-b.y,a.x-b.x),n=atan2(c.y-b.y,c.x-b.x); 268 v.e=v.s+Point(cos((n+m)/2),sin((n+m)/2)); 269 p=u.crosspoint(v); 270 r=Line(a,b).dispointtoseg(p); 271 } 272 273 void input(){ 274 p.input(); 275 scanf("%lf",&r); 276 } 277 278 void output(){ 279 printf("%.2f %.2f %.2f ",p.x,p.y,r); 280 } 281 282 bool operator==(circle v){ 283 return (p==v.p)&&sgn(r-v.r)==0; 284 } 285 286 bool operator<(circle v)const{ 287 return ((p<v.p)||((p==v.p)&&sgn(r-v.r)<0)); 288 } 289 290 double area(){ 291 return PI*r*r; 292 } 293 294 double circumference(){ ///周长 295 return 2*PI*r; 296 } 297 298 int relation(Point b){///点和圆的关系 0圆外 1圆上 2圆内 299 double dst=b.distance(p); 300 if(sgn(dst-r)<0) return 2; 301 else if(sgn(dst-r)==0) return 1; 302 return 0; 303 } 304 305 int relationseg(Line v){///线段和圆的关系，比较的是圆心到线段的距离和半径的关系 306 double dst=v.dispointtoseg(p); 307 if(sgn(dst-r)<0) return 2; 308 else if(sgn(dst-r)==0) return 1; 309 return 0; 310 } 311 312 int relationline(Line v){///直线和圆的关系，比较的是圆心到直线的距离和半径的关系 313 double dst=v.dispointtoline(p); 314 if(sgn(dst-r)<0) return 2; 315 else if(sgn(dst-r)==0) return 1; 316 return 0; 317 } 318 319 int relationcircle(circle v){///两圆的关系 5相离 4外切 3相交 2内切 1内含 320 double d=p.distance(v.p); 321 if(sgn(d-r-v.r)>0) return 5; 322 if(sgn(d-r-v.r)==0) return 4; 323 double l=fabs(r-v.r); 324 if(sgn(d-r-v.r)<0&&sgn(d-l)>0) return 3; 325 if(sgn(d-l)==0) return 2; 326 if(sgn(d-l)<0) return 1; 327 } 328 329 int pointcrosscircle(circle v,Point &p1,Point &p2){///求两个圆的交点，0没有交点 1一个交点 2两个交点 330 int rel=relationcircle(v); 331 if(rel == 1 || rel == 5) return 0; 332 double d=p.distance(v.p); 333 double l=(d*d+r*r-v.r*v.r)/2*d; 334 double h=sqrt(r*r-l*l); 335 Point tmp=p+(v.p-p).trunc(l); 336 p1=tmp+((v.p-p).rotleft().trunc(h)); 337 p2=tmp+((v.p-p).rotright().trunc(h)); 338 if(rel == 2 || rel == 4) return 1; 339 return 2; 340 } 341 342 int pointcrossline(Line v,Point &p1,Point &p2){///求直线和圆的交点，返回交点的个数 343 if(!(*this).relationline(v)) return 0; 344 Point a=v.lineprog(p); 345 double d=v.dispointtoline(p); 346 d=sqrt(r*r-d*d); 347 if(sgn(d)==0) { 348 p1=a; 349 p2=a; 350 return 1; 351 } 352 p1=a+(v.e-v.s).trunc(d); 353 p2=a-(v.e-v.s).trunc(d); 354 return 2; 355 } 356 357 int getcircle(Point a,Point b,double r1,circle &c1,circle &c2){///得到过a，b两点，半径为r1的两个圆 358 circle x(a,r1),y(b,r1); 359 int t=x.pointcrosscircle(y,c1.p,c2.p); 360 if(!t) return 0; 361 c1.r=c2.r=r; 362 return t; 363 } 364 365 int getcircle(Line u,Point q,double r1,circle &c1,circle &c2){///得到与直线u相切，过点q，半径为r1的圆 366 double dis = u.dispointtoline(q); 367 if(sgn(dis-r1*2)>0) return 0; 368 if(sgn(dis)==0) { 369 c1.p=q+((u.e-u.s).rotleft().trunc(r1)); 370 c2.p=q+((u.e-u.s).rotright().trunc(r1)); 371 c1.r=c2.r=r1; 372 return 2; 373 } 374 Line u1=Line((u.s+(u.e-u.s).rotleft().trunc(r1)),(u.e+(u.e-u.s).rotleft().trunc(r1))); 375 Line u2=Line((u.s+(u.e-u.s).rotright().trunc(r1)),(u.e+(u.e-u.s).rotright().trunc(r1))); 376 circle cc=circle(q,r1); 377 Point p1,p2; 378 if(!cc.pointcrossline(u1,p1,p2)) cc.pointcrossline(u2,p1,p2); 379 c1=circle(p1,r1); 380 if(p1==p2){ 381 c2=c1; 382 return 1; 383 } 384 c2=circle(p2,r1); 385 return 2; 386 } 387 388 int getcircle(Line u,Line v,double r1,circle &c1,circle &c2,circle &c3,circle &c4){///同时与直线u，v相切，半径为r1的圆 389 if(u.parallel(v)) return 0;///两直线平行 390 Line u1=Line(u.s+(u.e-u.s).rotleft().trunc(r1),u.e+(u.e-u.s).rotleft().trunc(r1)); 391 Line u2=Line(u.s+(u.e-u.s).rotright().trunc(r1),u.e+(u.e-u.s).rotright().trunc(r1)); 392 Line v1=Line(v.s+(v.e-v.s).rotleft().trunc(r1),v.e+(v.e-v.s).rotleft().trunc(r1)); 393 Line v2=Line(v.s+(v.e-v.s).rotright().trunc(r1),v.e+(v.e-v.s).rotright().trunc(r1)); 394 c1.r=c2.r=c3.r=c4.r=r1; 395 c1.p=u1.crosspoint(v1); 396 c2.p=u1.crosspoint(v2); 397 c3.p=u2.crosspoint(v1); 398 c4.p=u2.crosspoint(v2); 399 return 4; 400 } 401 402 int getcircle(circle cx,circle cy,double r1,circle &c1,circle &c2){///同时与不相交圆 cx,cy相切，半径为r1的圆 403 circle x(cx.p,r1+cx.r),y(cy.p,r1+cy.r); 404 int t=x.pointcrosscircle(y,c1.p,c2.p); 405 if(!t) return 0; 406 c1.r=c2.r=r1; 407 return t; 408 } 409 410 int tangentline(Point q,Line &u,Line &v){///过一点作圆的切线（先判断点和圆的关系） 411 int x=relation(q); 412 if(x==2) return 0; 413 if(x==1){ 414 u=Line(q,q+(q-p).rotleft()); 415 v=u; 416 return 1; 417 } 418 double d=p.distance(q); 419 double l=r*r/d; 420 double h=sqrt(r*r-l*l); 421 u=Line(q,p+((q-p).trunc(l)+(q-p).rotleft().trunc(h))); 422 v=Line(q,p+((q-p).trunc(l)+(q-p).rotright().trunc(h))); 423 return 2; 424 } 425 426 double areacircle(circle v){///求两圆相交的面积 427 int rel=relationcircle(v); 428 if(rel >= 4) return 0.0; 429 if(rel <= 2) return min(area(),v.area()); 430 double d=p.distance(v.p); 431 double hf=(r+v.r+d)/2.0; 432 double ss=2*sqrt(hf*(hf-r)*(hf-v.r)*(hf-d)); 433 double a1=acos((r*r+d*d-v.r*v.r)/(2.0*r*d)); 434 a1=a1*r*r; 435 double a2=acos((v.r*v.r+d*d-r*r)/(2.0*v.r*d)); 436 a2=a2*v.r*v.r; 437 return a1+a2-ss; 438 } 439 440 double areatriangle(Point a,Point b){///求圆和三角形pab的相交面积 441 if(sgn((p-a)^(p-b))==0) return 0.0; 442 Point q[5]; 443 int len=0; 444 q[len++]=a; 445 Line l(a,b); 446 Point p1,p2; 447 if(pointcrossline(l,q[1],q[2])==2){ 448 if(sgn((a-q[1])*(b-q[1]))<0) q[len++]=q[1]; 449 if(sgn((a-q[2])*(b-q[2]))<0) q[len++]=q[2]; 450 } 451 q[len++]=b; 452 if(len==4 && sgn((q[0]-q[1])*(q[2]-q[1]))>0) swap(q[1],q[2]); 453 double res=0; 454 for(int i=0;i<len-1;i++){ 455 if(relation(q[i])==0||relation(q[i+1])==0){ 456 double arg=p.rad(q[i],q[i+1]); 457 res+=r*r*arg/2.0; 458 } 459 else{ 460 res+=fabs((q[i]-p)^(q[i+1]-p))/2.0; 461 } 462 } 463 return res; 464 } 465 }; 466 467 struct polygon{ 468 int n; 469 Point p[100010]; 470 Line l[100010]; 471 void input(int _n){ 472 n=_n; 473 for(int i=0;i<n;i++){ 474 p[i].input(); 475 } 476 } 477 478 void add(Point q){ 479 p[n++]=q; 480 } 481 482 void getline(){ 483 for(int i=0;i<n;i++){ 484 l[i]=Line(p[i],p[(i+1)%n]); 485 } 486 } 487 488 struct cmp{ 489 Point p; 490 cmp(const Point &p0){p=p0;} 491 bool operator()(const Point &aa,const Point &bb){ 492 Point a=aa,b=bb; 493 int d=sgn((a-p)^(b-p)); 494 if(d==0){ 495 return sgn(a.distance(p)-b.distance(p))<0; 496 } 497 return d>0; 498 } 499 }; 500 501 void norm(){///进行极角排序 502 Point mi=p[0]; 503 for(int i=1;i<n;i++){ 504 mi=min(mi,p[i]); 505 } 506 sort(p,p+n,cmp(mi)); 507 } 508 509 void getconvex(polygon &convex){///得到第一种凸包的方法，编号为0～n-1，可能要特判所有点共点或共线的特殊情况 510 sort(p,p+n); 511 convex.n=n; 512 for(int i=0;i<min(n,2);i++){ 513 convex.p[i]=p[i]; 514 } 515 if(convex.n==2&&(convex.p[0]==convex.p[1])) convex.n--; 516 if(n<=2) return; 517 int &top=convex.n; 518 top=1; 519 for(int i=2;i<n;i++){ 520 while(top&&sgn((convex.p[top]-p[i])^(convex.p[top-1]-p[i]))<0) top--; 521 convex.p[++top]=p[i]; 522 } 523 int temp=top; 524 convex.p[++top]=p[n-2]; 525 for(int i=n-3;i>=0;i--){ 526 while(top!=temp&&sgn((convex.p[top]-p[i])^(convex.p[top-1]-p[i]))<=0) top--; 527 convex.p[++top]=p[i]; 528 } 529 if(convex.n==2&&(convex.p[0]==convex.p[1])) convex.n--; 530 convex.norm(); 531 } 532 533 void Graham(polygon &convex){///得到凸包的第二种方法 534 norm(); 535 int &top=convex.n; 536 top=0; 537 if(n==1){ 538 top=1; 539 convex.p[0]=p[0]; 540 return; 541 } 542 if(n==2){ 543 top=2; 544 convex.p[0]=p[0]; 545 convex.p[1]=p[1]; 546 if(convex.p[0]==convex.p[1]) top--; 547 return; 548 } 549 convex.p[0]=p[0]; 550 convex.p[1]=p[1]; 551 top=2; 552 for(int i=2;i<n;i++){ 553 while(top>1&&sgn((convex.p[top-1]-convex.p[top-2])^(p[i]-convex.p[top-2]))<=0) top--; 554 convex.p[top++]=p[i]; 555 } 556 if(convex.n==2 && (convex.p[0]==convex.p[1])) convex.n--; 557 } 558 559 bool inconvex(){///判断是不是凸的 560 bool s[3]; 561 memset(s,false,sizeof(s)); 562 for(int i=0;i<n;i++){ 563 int j=(i+1)%n; 564 int k=(j+1)%n; 565 s[sgn((p[j]-p[i])^(p[k]-p[i]))+1]=true; 566 if(s[0]&&s[2]) return false; 567 } 568 return true; 569 } 570 571 int relationpoint(Point q){///判断点和任意多边形的关系 3点上 2边上 1内部 0外部 572 for(int i=0;i<n;i++){ 573 if(p[i]==q) return 3; 574 } 575 getline(); 576 for(int i=0;i<n;i++){ 577 if(l[i].pointonseg(q)) return 2; 578 } 579 int cnt=0; 580 for(int i=0;i<n;i++){ 581 int j=(i+1)%n; 582 int k=sgn((q-p[j])^(p[i]-p[j])); 583 int u=sgn(p[i].y-q.y); 584 int v=sgn(p[j].y-q.y); 585 if(k>0&&u<0&&v>=0) cnt++; 586 if(k<0&&v<0&&u>=0) cnt--; 587 } 588 return cnt!=0; 589 } 590 591 void convexcnt(Line u,polygon &po){///直线u切割凸多边形左侧 注意直线方向 592 int &top=po.n; 593 top=0; 594 for(int i=0;i<n;i++){ 595 int d1=sgn((u.e-u.s)^(p[i]-u.s)); 596 int d2=sgn((u.e-u.s)^(p[(i+1)%n]-u.s)); 597 if(d1>=0) po.p[top++]=p[i]; 598 if(d1*d2<0)po.p[top++]=u.crosspoint(Line(p[i],p[(i+1)%n])); 599 } 600 } 601 602 double getcircumference(){///得到周长 603 double sum=0; 604 for(int i=0;i<n;i++){ 605 sum+=p[i].distance(p[(i+1)%n]); 606 } 607 return sum; 608 } 609 610 double getarea(){///得到面积 611 double sum=0; 612 for(int i=0;i<n;i++){ 613 sum+=(p[i]^p[(i+1)%n]); 614 } 615 return fabs(sum)/2; 616 } 617 618 bool getdir(){///得到方向 1表示逆时针 0表示顺时针 619 double sum=0; 620 for(int i=0;i<n;i++){ 621 sum+=(p[i]^p[(i+1)%n]); 622 } 623 if(sgn(sum)>0) return 1; 624 return 0; 625 } 626 627 Point getbarycentre(){///得到重心 628 Point ret(0,0); 629 double area=0; 630 for(int i=1;i<n-1;i++){ 631 double tmp=(p[i]-p[0])^(p[i+1]-p[0]); 632 if(sgn(tmp)==0) continue; 633 area+=tmp; 634 ret.x+=(p[0].x+p[i].x+p[i+1].x)/3*tmp; 635 ret.y+=(p[0].y+p[i].y+p[i+1].y)/3*tmp; 636 } 637 if(sgn(area)) ret =ret/area; 638 return ret; 639 } 640 641 double areacircle(circle c){///多边形和圆交的面积 642 double ans=0; 643 for(int i=0;i<n;i++){ 644 int j=(i+1)%n; 645 if(sgn((p[j]-c.p)^(p[i]-c.p))>=0) ans+=c.areatriangle(p[i],p[j]); 646 else ans-=c.areatriangle(p[i],p[j]); 647 } 648 return fabs(ans); 649 } 650 651 int relationcircle(circle c){///多边形和圆的关系 2圆完全在多边形内 1圆在多边形里面，碰到了多边形的边界 0其他 652 getline(); 653 int x=2; 654 if(relationpoint(c.p)!=1) return 0; 655 for(int i=0;i<n;i++){ 656 if(c.relationseg(l[i])==2) return 0; 657 if(c.relationseg(l[i])==1) x=1; 658 } 659 return x; 660 } 661 }; 662 663 double cross(Point a,Point b,Point c){///ab x ac 664 return (b-a)^(c-a); 665 } 666 667 double dot(Point a,Point b,Point c){///ab*ac; 668 return (b-a)*(c-a); 669 } 670 671 double minRectangleCover(polygon A){///`最小矩形面积覆盖`,` A 必须是凸包(而且是逆时针顺序)`,` 测试 UVA 10173` 672 //`要特判A.n < 3的情况` 673 if(A.n < 3)return 0.0; 674 A.p[A.n] = A.p[0]; 675 double ans = -1; 676 int r = 1, p = 1, q; 677 for(int i = 0;i < A.n;i++){ 678 //`卡出离边A.p[i] - A.p[i+1]最远的点` 679 while( sgn( cross(A.p[i],A.p[i+1],A.p[r+1]) - cross(A.p[i],A.p[i+1],A.p[r]) ) >= 0 ) 680 r = (r+1)%A.n; 681 //`卡出A.p[i] - A.p[i+1]方向上正向n最远的点` 682 while(sgn( dot(A.p[i],A.p[i+1],A.p[p+1]) - dot(A.p[i],A.p[i+1],A.p[p]) ) >= 0 ) 683 p = (p+1)%A.n; 684 if(i == 0)q = p; 685 //`卡出A.p[i] - A.p[i+1]方向上负向最远的点` 686 while(sgn(dot(A.p[i],A.p[i+1],A.p[q+1]) - dot(A.p[i],A.p[i+1],A.p[q])) <= 0) 687 q = (q+1)%A.n; 688 double d = (A.p[i] - A.p[i+1]).len2(); 689 double tmp = cross(A.p[i],A.p[i+1],A.p[r]) * 690 (dot(A.p[i],A.p[i+1],A.p[p]) - dot(A.p[i],A.p[i+1],A.p[q]))/d; 691 if(ans < 0 || ans > tmp)ans = tmp; 692 } 693 return ans; 694 } 695 696 vector<Point> convexCut(const vector<Point>&ps,Point q1,Point q2){///直线切凸多边形，多边形是逆时针的，在q1q2的左侧 697 vector<Point>qs; 698 int n=ps.size(); 699 for(int i=0;i<n;i++){ 700 Point p1=ps[i],p2=ps[(i+1)%n]; 701 int d1=sgn((q2-q1)^(p1-q1)),d2=sgn((q2-q1)^(p2-q1)); 702 if(d1>=0) qs.push_back(p1); 703 if(d1*d2<0) qs.push_back(Line(p1,p2).crosspoint(Line(q1,q2))); 704 } 705 return qs; 706 } 707 708 struct halfplane:public Line{ 709 double angle; 710 halfplane(){} 711 halfplane(Point _s,Point _e){///表示向量s->e逆时针（左侧）的半平面 712 s=_s; 713 e=_e; 714 } 715 halfplane(Line v){ 716 s=v.s; 717 e=v.e; 718 } 719 void calcangle(){ 720 angle=atan2(e.y-s.y,e.x-s.x); 721 } 722 bool operator<(const halfplane &b)const{ 723 return angle<b.angle; 724 } 725 }; 726 727 struct halfplanes{ 728 int n; 729 halfplane hp[2020]; 730 Point p[2020]; 731 int que[2020]; 732 int st,ed; 733 void push(halfplane tmp){ 734 hp[n++]=tmp; 735 } 736 737 void unique(){///去重 738 int m=1; 739 for(int i=1;i<n;i++){ 740 if(sgn(hp[i].angle-hp[i-1].angle)!=0) hp[m++]=hp[i]; 741 else if(sgn((hp[m-1].e-hp[m-1].s)^(hp[i].s-hp[m-1].s))>0) hp[m-1]=hp[i]; 742 } 743 n=m; 744 } 745 bool halfplaneinsert(){ 746 for(int i=0;i<n;i++) hp[i].calcangle(); 747 sort(hp,hp+n); 748 unique(); 749 que[st=0]=0; 750 que[ed=1]=1; 751 p[1]=hp[0].crosspoint(hp[1]); 752 for(int i=2;i<n;i++){ 753 while(st<ed&&sgn((hp[i].e-hp[i].s)^(p[ed]-hp[i].s))<0) ed--; 754 while(st<ed&&sgn((hp[i].e-hp[i].s)^(p[st+1]-hp[i].s))<0) st++; 755 que[++ed]=i; 756 if(hp[i].parallel(hp[que[ed-1]])) return false; 757 p[ed]=hp[i].crosspoint(hp[que[ed-1]]); 758 } 759 while(st<ed&&sgn((hp[que[st]].e-hp[que[st]].s)^(p[ed]-hp[que[st]].s))<0) ed--; 760 while(st<ed&&sgn((hp[que[ed]].e-hp[que[ed]].s)^(p[st+1]-hp[que[ed]].s))<0) st++; 761 if(st+1>=ed) return false; 762 return true; 763 } 764 765 void getconvex(polygon &con){///得到最后半平面交得到的凸多边形，要先调用halfplaneinsert()且返回true 766 p[st]=hp[que[st]].crosspoint(hp[que[ed]]); 767 con.n=ed-st+1; 768 for(int j=st,i=0;j<=ed;i++,j++){ 769 con.p[i]=p[j]; 770 } 771 } 772 }; 773 774 struct circles{ 775 circle c[1010]; 776 double ans[1010];///ans[i]表示被覆盖了i次的面积 777 double pre[1010]; 778 int n; 779 circles(){} 780 void add(circle cc){ 781 c[n++]=cc; 782 } 783 784 bool inner(circle x,circle y){///x包含在y中 785 if(x.relationcircle(y)!=1) return 0; 786 return sgn(x.r-y.r)<=0?1:0; 787 } 788 789 void init_or(){///圆的面积并去掉内含的圆 790 bool mark[1010]={0}; 791 int i,j,k=0; 792 for(i=0;i<n;i++){ 793 for(j=0;j<n;j++){ 794 if(i!=j&&!mark[j]){ 795 if(c[i]==c[j]||inner(c[i],c[j])) break; 796 } 797 } 798 if(j<n) mark[i]=1; 799 } 800 for(i=0;i<n;i++){ 801 if(!mark[i]) c[k++]=c[i]; 802 } 803 n=k; 804 } 805 806 void init_add(){///圆的面积交去掉内含的圆 807 int i,j,k; 808 bool mark[1010]={0}; 809 for(i=0;i<n;i++){ 810 for(int j=0;j<n;j++){ 811 if(i!=j&&!mark[j]){ 812 if((c[i]==c[j])||inner(c[j],c[i])) break; 813 } 814 } 815 if(j<n) mark[i]=1; 816 } 817 for(i=0;i<n;i++){ 818 if(!mark[i]){ 819 c[k++]=c[i]; 820 } 821 } 822 n=k; 823 } 824 825 double areaarc(double th,double r){///半径为r的圆，弧度为th，对应的弓形的面积 826 return 0.5*r*r*(th-sin(th)); 827 } 828 829 830 void getarea(){ 831 memset(ans,0,sizeof(ans)); 832 vector<pair<double,int> >v; 833 for(int i=0;i<n;i++){ 834 v.clear(); 835 v.push_back(make_pair(-PI,1)); 836 v.push_back(make_pair(PI,-1)); 837 for(int j=0;j<n;j++){ 838 if(i!=j){ 839 Point q=(c[j].p-c[i].p); 840 double ab=q.len(),ac=c[i].r,bc=c[j].r; 841 if(sgn(ab+ac-bc)<=0){ 842 v.push_back(make_pair(-PI,1)); 843 v.push_back(make_pair(PI,-1)); 844 continue; 845 } 846 if(sgn(ab+bc-ac)<=0)continue; 847 if(sgn(ab-ac-bc)>0) continue; 848 double th=atan2(q.y,q.x),fai=acos((ac*ac+ab*ab-bc*bc)/(2.0*ac*ab)); 849 double a0=th-fai; 850 if(sgn(a0+PI)<0) a0+=2*PI; 851 double a1=th+fai; 852 if(sgn(a1-PI)>0) a1-=2*PI; 853 if(sgn(a0-a1)>0){ 854 v.push_back(make_pair(a0,1)); 855 v.push_back(make_pair(PI,-1)); 856 v.push_back(make_pair(-PI,1)); 857 v.push_back(make_pair(a1,-1)); 858 } 859 else{ 860 v.push_back(make_pair(a0,1)); 861 v.push_back(make_pair(a1,-1)); 862 } 863 } 864 sort(v.begin(),v.end()); 865 int cur=0; 866 for(int j=0;j<v.size();j++){ 867 if(cur&&sgn(v[j].first-pre[cur])){ 868 ans[cur]+=areaarc(v[j].first-pre[cur],c[i].r); 869 ans[cur]+=0.5*(Point(c[i].p.x+c[i].r*cos(pre[cur]),c[i].p.y+c[i].r*sin(pre[cur]))^Point(c[i].p.x+c[i].r*cos(v[j].first),c[i].p.y+c[i].r*sin(v[j].first))); 870 } 871 cur+=v[j].second; 872 pre[cur]=v[j].first; 873 } 874 } 875 } 876 for(int i=1;i<n;i++){ 877 ans[i]-=ans[i+1]; 878 } 879 } 880 }; 881 882 883 bool Check(Line a,Line b){ 884 if(sgn((a.s-a.e)^(b.s-a.e))*sgn((a.s-a.e)^(b.e-a.e))>0) return false; 885 if(sgn((b.s-b.e)^(a.s-b.e))*sgn((b.s-b.e)^(a.e-b.e))>0) return false; 886 if(sgn(max(a.s.x,a.e.x)-min(b.s.x,b.e.x))>=0&&sgn(max(b.s.x,b.e.x)-min(a.s.x,a.e.x))>=0 887 &&sgn(max(a.s.y,a.e.y)-min(b.s.y,b.e.y))>=0&&sgn(max(b.s.y,b.e.y)-min(a.s.y,a.e.y))>=0) 888 return true; 889 else return false; 890 } 891 892 halfplanes hps; 893 int n; 894 Point p[105]; 895 polygon poly; 896 897 int main(){ 898 int t; 899 scanf("%d",&t); 900 while(t--){ 901 int n; 902 scanf("%d",&n); 903 for(int i=0;i<n;i++){ 904 p[i].input(); 905 } 906 hps.n=0; 907 halfplane hp; 908 for(int i=0;i<n;i++){ 909 hp.e=p[i]; 910 hp.s=p[(i+1)%n]; 911 hps.push(hp); 912 } 913 if(hps.halfplaneinsert()){ 914 hps.getconvex(poly); 915 printf("%.2f ",poly.getarea()); 916 } 917 else{ 918 printf("0.00 "); 919 } 920 921 } 922 return 0; 923 }