Pentium.Labs

System全家桶:https://zhuanlan.zhihu.com/c_1238468913098731520

导航

[转]计算几何模板

Repost

 

1.

   1 #include<vector>
   2 #include<list>
   3 #include<map>
   4 #include<set>
   5 #include<deque>
   6 #include<queue>
   7 #include<stack>
   8 #include<bitset>
   9 #include<algorithm>
  10 #include<functional>
  11 #include<numeric>
  12 #include<utility>
  13 #include<iostream>
  14 #include<sstream>
  15 #include<iomanip>
  16 #include<cstdio>
  17 #include<cmath>
  18 #include<cstdlib>
  19 #include<cctype>
  20 #include<string>
  21 #include<cstring>
  22 #include<cstdio>
  23 #include<cmath>
  24 #include<cstdlib>
  25 #include<ctime>
  26 #include<climits>
  27 #include<complex>
  28 #define mp make_pair
  29 #define pb push_back
  30 using namespace std;
  31 const double eps=1e-8;
  32 const double pi=acos(-1.0);
  33 const double inf=1e20;
  34 const int maxp=1111;
  35 int dblcmp(double d)
  36 {
  37     if (fabs(d)<eps)return 0;
  38     return d>eps?1:-1;
  39 }
  40 inline double sqr(double x){return x*x;}
  41 struct point
  42 {
  43     double x,y;
  44     point(){}
  45     point(double _x,double _y):
  46     x(_x),y(_y){};
  47     void input()
  48     {
  49         scanf("%lf%lf",&x,&y);
  50     }
  51     void output()
  52     {
  53         printf("%.2f %.2f\n",x,y);
  54     }
  55     bool operator==(point a)const
  56     {
  57         return dblcmp(a.x-x)==0&&dblcmp(a.y-y)==0;
  58     }
  59     bool operator<(point a)const
  60     {
  61         return dblcmp(a.x-x)==0?dblcmp(y-a.y)<0:x<a.x;
  62     }
  63     double len()
  64     {
  65         return hypot(x,y);
  66     }
  67     double len2()
  68     {
  69         return x*x+y*y;
  70     }
  71     double distance(point p)
  72     {
  73         return hypot(x-p.x,y-p.y);
  74     }
  75     point add(point p)
  76     {
  77         return point(x+p.x,y+p.y);
  78     }
  79     point sub(point p)
  80     {
  81         return point(x-p.x,y-p.y);
  82     }
  83     point mul(double b)
  84     {
  85         return point(x*b,y*b);
  86     }
  87     point div(double b)
  88     {
  89         return point(x/b,y/b);
  90     }
  91     double dot(point p)
  92     {
  93         return x*p.x+y*p.y;
  94     }
  95     double det(point p)
  96     {
  97         return x*p.y-y*p.x;
  98     }
  99     double rad(point a,point b)
 100     {
 101         point p=*this;
 102         return fabs(atan2(fabs(a.sub(p).det(b.sub(p))),a.sub(p).dot(b.sub(p))));
 103     }
 104     point trunc(double r)
 105     {
 106         double l=len();
 107         if (!dblcmp(l))return *this;
 108         r/=l;
 109         return point(x*r,y*r);
 110     }
 111     point rotleft()
 112     {
 113         return point(-y,x);
 114     }
 115     point rotright()
 116     {
 117         return point(y,-x);
 118     }
 119     point rotate(point p,double angle)//绕点p逆时针旋转angle角度
 120     {
 121         point v=this->sub(p);
 122         double c=cos(angle),s=sin(angle);
 123         return point(p.x+v.x*c-v.y*s,p.y+v.x*s+v.y*c);
 124     }
 125 };
 126 struct line
 127 {
 128     point a,b;
 129     line(){}
 130     line(point _a,point _b)
 131     {
 132         a=_a;
 133         b=_b;
 134     }
 135     bool operator==(line v)
 136     {
 137         return (a==v.a)&&(b==v.b);
 138     }
 139     //倾斜角angle
 140     line(point p,double angle)
 141     {
 142         a=p;
 143         if (dblcmp(angle-pi/2)==0)
 144         {
 145             b=a.add(point(0,1));
 146         }
 147         else
 148         {
 149             b=a.add(point(1,tan(angle)));
 150         }
 151     }
 152     //ax+by+c=0
 153     line(double _a,double _b,double _c)
 154     {
 155         if (dblcmp(_a)==0)
 156         {
 157             a=point(0,-_c/_b);
 158             b=point(1,-_c/_b);
 159         }
 160         else if (dblcmp(_b)==0)
 161         {
 162             a=point(-_c/_a,0);
 163             b=point(-_c/_a,1);
 164         }
 165         else
 166         {
 167             a=point(0,-_c/_b);
 168             b=point(1,(-_c-_a)/_b);
 169         }
 170     }
 171     void input()
 172     {
 173         a.input();
 174         b.input();
 175     }
 176     void adjust()
 177     {
 178         if (b<a)swap(a,b);
 179     }
 180     double length()
 181     {
 182         return a.distance(b);
 183     }
 184     double angle()//直线倾斜角 0<=angle<180
 185     {
 186         double k=atan2(b.y-a.y,b.x-a.x);
 187         if (dblcmp(k)<0)k+=pi;
 188         if (dblcmp(k-pi)==0)k-=pi;
 189         return k;
 190     }
 191     //点和线段关系
 192     //1 在逆时针
 193     //2 在顺时针
 194     //3 平行
 195     int relation(point p)
 196     {
 197         int c=dblcmp(p.sub(a).det(b.sub(a)));
 198         if (c<0)return 1;
 199         if (c>0)return 2;
 200         return 3;
 201     }
 202     bool pointonseg(point p)
 203     {
 204         return dblcmp(p.sub(a).det(b.sub(a)))==0&&dblcmp(p.sub(a).dot(p.sub(b)))<=0;
 205     }
 206     bool parallel(line v)
 207     {
 208         return dblcmp(b.sub(a).det(v.b.sub(v.a)))==0;
 209     }
 210     //2 规范相交
 211     //1 非规范相交
 212     //0 不相交
 213     int segcrossseg(line v)
 214     {
 215         int d1=dblcmp(b.sub(a).det(v.a.sub(a)));
 216         int d2=dblcmp(b.sub(a).det(v.b.sub(a)));
 217         int d3=dblcmp(v.b.sub(v.a).det(a.sub(v.a)));
 218         int d4=dblcmp(v.b.sub(v.a).det(b.sub(v.a)));
 219         if ((d1^d2)==-2&&(d3^d4)==-2)return 2;
 220         return (d1==0&&dblcmp(v.a.sub(a).dot(v.a.sub(b)))<=0||
 221                 d2==0&&dblcmp(v.b.sub(a).dot(v.b.sub(b)))<=0||
 222                 d3==0&&dblcmp(a.sub(v.a).dot(a.sub(v.b)))<=0||
 223                 d4==0&&dblcmp(b.sub(v.a).dot(b.sub(v.b)))<=0);
 224     }
 225     int linecrossseg(line v)//*this seg v line
 226     {
 227         int d1=dblcmp(b.sub(a).det(v.a.sub(a)));
 228         int d2=dblcmp(b.sub(a).det(v.b.sub(a)));
 229         if ((d1^d2)==-2)return 2;
 230         return (d1==0||d2==0);
 231     }
 232     //0 平行
 233     //1 重合
 234     //2 相交
 235     int linecrossline(line v)
 236     {
 237         if ((*this).parallel(v))
 238         {
 239             return v.relation(a)==3;
 240         }
 241         return 2;
 242     }
 243     point crosspoint(line v)
 244     {
 245         double a1=v.b.sub(v.a).det(a.sub(v.a));
 246         double a2=v.b.sub(v.a).det(b.sub(v.a));
 247         return point((a.x*a2-b.x*a1)/(a2-a1),(a.y*a2-b.y*a1)/(a2-a1));
 248     }
 249     double dispointtoline(point p)
 250     {
 251         return fabs(p.sub(a).det(b.sub(a)))/length();
 252     }
 253     double dispointtoseg(point p)
 254     {
 255         if (dblcmp(p.sub(b).dot(a.sub(b)))<0||dblcmp(p.sub(a).dot(b.sub(a)))<0)
 256         {
 257             return min(p.distance(a),p.distance(b));
 258         }
 259         return dispointtoline(p);
 260     }
 261     point lineprog(point p)
 262     {
 263         return a.add(b.sub(a).mul(b.sub(a).dot(p.sub(a))/b.sub(a).len2()));
 264     }
 265     point symmetrypoint(point p)
 266     {
 267         point q=lineprog(p);
 268         return point(2*q.x-p.x,2*q.y-p.y);
 269     }
 270 };
 271 struct circle
 272 {
 273     point p;
 274     double r;
 275     circle(){}
 276     circle(point _p,double _r):
 277     p(_p),r(_r){};
 278     circle(double x,double y,double _r):
 279     p(point(x,y)),r(_r){};
 280     circle(point a,point b,point c)//三角形的外接圆
 281     {
 282         p=line(a.add(b).div(2),a.add(b).div(2).add(b.sub(a).rotleft())).crosspoint(line(c.add(b).div(2),c.add(b).div(2).add(b.sub(c).rotleft())));
 283         r=p.distance(a);
 284     }
 285     circle(point a,point b,point c,bool t)//三角形的内切圆
 286     {
 287         line u,v;
 288         double m=atan2(b.y-a.y,b.x-a.x),n=atan2(c.y-a.y,c.x-a.x);
 289         u.a=a;
 290         u.b=u.a.add(point(cos((n+m)/2),sin((n+m)/2)));
 291         v.a=b;
 292         m=atan2(a.y-b.y,a.x-b.x),n=atan2(c.y-b.y,c.x-b.x);
 293         v.b=v.a.add(point(cos((n+m)/2),sin((n+m)/2)));
 294         p=u.crosspoint(v);
 295         r=line(a,b).dispointtoseg(p);
 296     }
 297     void input()
 298     {
 299         p.input();
 300         scanf("%lf",&r);
 301     }
 302     void output()
 303     {
 304         printf("%.2lf %.2lf %.2lf\n",p.x,p.y,r);
 305     }
 306     bool operator==(circle v)
 307     {
 308         return ((p==v.p)&&dblcmp(r-v.r)==0);
 309     }
 310     bool operator<(circle v)const
 311     {
 312         return ((p<v.p)||(p==v.p)&&dblcmp(r-v.r)<0);
 313     }
 314     double area()
 315     {
 316         return pi*sqr(r);
 317     }
 318     double circumference()
 319     {
 320         return 2*pi*r;
 321     }
 322     //0 圆外
 323     //1 圆上
 324     //2 圆内
 325     int relation(point b)
 326     {
 327         double dst=b.distance(p);
 328         if (dblcmp(dst-r)<0)return 2;
 329         if (dblcmp(dst-r)==0)return 1;
 330         return 0;
 331     }
 332     int relationseg(line v)
 333     {
 334         double dst=v.dispointtoseg(p);
 335         if (dblcmp(dst-r)<0)return 2;
 336         if (dblcmp(dst-r)==0)return 1;
 337         return 0;
 338     }
 339     int relationline(line v)
 340     {
 341         double dst=v.dispointtoline(p);
 342         if (dblcmp(dst-r)<0)return 2;
 343         if (dblcmp(dst-r)==0)return 1;
 344         return 0;
 345     }
 346     //过a b两点 半径r的两个圆
 347     int getcircle(point a,point b,double r,circle&c1,circle&c2)
 348     {
 349         circle x(a,r),y(b,r);
 350         int t=x.pointcrosscircle(y,c1.p,c2.p);
 351         if (!t)return 0;
 352         c1.r=c2.r=r;
 353         return t;
 354     }
 355     //与直线u相切 过点q 半径r1的圆
 356     int getcircle(line u,point q,double r1,circle &c1,circle &c2)
 357     {
 358         double dis=u.dispointtoline(q);
 359         if (dblcmp(dis-r1*2)>0)return 0;
 360         if (dblcmp(dis)==0)
 361         {
 362             c1.p=q.add(u.b.sub(u.a).rotleft().trunc(r1));
 363             c2.p=q.add(u.b.sub(u.a).rotright().trunc(r1));
 364             c1.r=c2.r=r1;
 365             return 2;
 366         }
 367         line u1=line(u.a.add(u.b.sub(u.a).rotleft().trunc(r1)),u.b.add(u.b.sub(u.a).rotleft().trunc(r1)));
 368         line u2=line(u.a.add(u.b.sub(u.a).rotright().trunc(r1)),u.b.add(u.b.sub(u.a).rotright().trunc(r1)));
 369         circle cc=circle(q,r1);
 370         point p1,p2;
 371         if (!cc.pointcrossline(u1,p1,p2))cc.pointcrossline(u2,p1,p2);
 372         c1=circle(p1,r1);
 373         if (p1==p2)
 374         {
 375             c2=c1;return 1;
 376         }
 377         c2=circle(p2,r1);
 378         return 2;
 379     }
 380     //同时与直线u,v相切 半径r1的圆
 381     int getcircle(line u,line v,double r1,circle &c1,circle &c2,circle &c3,circle &c4)
 382     {
 383         if (u.parallel(v))return 0;
 384         line u1=line(u.a.add(u.b.sub(u.a).rotleft().trunc(r1)),u.b.add(u.b.sub(u.a).rotleft().trunc(r1)));
 385         line u2=line(u.a.add(u.b.sub(u.a).rotright().trunc(r1)),u.b.add(u.b.sub(u.a).rotright().trunc(r1)));
 386         line v1=line(v.a.add(v.b.sub(v.a).rotleft().trunc(r1)),v.b.add(v.b.sub(v.a).rotleft().trunc(r1)));
 387         line v2=line(v.a.add(v.b.sub(v.a).rotright().trunc(r1)),v.b.add(v.b.sub(v.a).rotright().trunc(r1)));
 388         c1.r=c2.r=c3.r=c4.r=r1;
 389         c1.p=u1.crosspoint(v1);
 390         c2.p=u1.crosspoint(v2);
 391         c3.p=u2.crosspoint(v1);
 392         c4.p=u2.crosspoint(v2);
 393         return 4;
 394     }
 395     //同时与不相交圆cx,cy相切 半径为r1的圆
 396     int getcircle(circle cx,circle cy,double r1,circle&c1,circle&c2)
 397     {
 398         circle x(cx.p,r1+cx.r),y(cy.p,r1+cy.r);
 399         int t=x.pointcrosscircle(y,c1.p,c2.p);
 400         if (!t)return 0;
 401         c1.r=c2.r=r1;
 402         return t;
 403     }
 404     int pointcrossline(line v,point &p1,point &p2)//求与线段交要先判断relationseg
 405     {
 406         if (!(*this).relationline(v))return 0;
 407         point a=v.lineprog(p);
 408         double d=v.dispointtoline(p);
 409         d=sqrt(r*r-d*d);
 410         if (dblcmp(d)==0)
 411         {
 412             p1=a;
 413             p2=a;
 414             return 1;
 415         }
 416         p1=a.sub(v.b.sub(v.a).trunc(d));
 417         p2=a.add(v.b.sub(v.a).trunc(d));
 418         return 2;
 419     }
 420     //5 相离
 421     //4 外切
 422     //3 相交
 423     //2 内切
 424     //1 内含
 425     int relationcircle(circle v)
 426     {
 427         double d=p.distance(v.p);
 428         if (dblcmp(d-r-v.r)>0)return 5;
 429         if (dblcmp(d-r-v.r)==0)return 4;
 430         double l=fabs(r-v.r);
 431         if (dblcmp(d-r-v.r)<0&&dblcmp(d-l)>0)return 3;
 432         if (dblcmp(d-l)==0)return 2;
 433         if (dblcmp(d-l)<0)return 1;
 434     }
 435     int pointcrosscircle(circle v,point &p1,point &p2)
 436     {
 437         int rel=relationcircle(v);
 438         if (rel==1||rel==5)return 0;
 439         double d=p.distance(v.p);
 440         double l=(d+(sqr(r)-sqr(v.r))/d)/2;
 441         double h=sqrt(sqr(r)-sqr(l));
 442         p1=p.add(v.p.sub(p).trunc(l).add(v.p.sub(p).rotleft().trunc(h)));
 443         p2=p.add(v.p.sub(p).trunc(l).add(v.p.sub(p).rotright().trunc(h)));
 444         if (rel==2||rel==4)
 445         {
 446             return 1;
 447         }
 448         return 2;
 449     }
 450     //过一点做圆的切线 (先判断点和圆关系)
 451     int tangentline(point q,line &u,line &v)
 452     {
 453         int x=relation(q);
 454         if (x==2)return 0;
 455         if (x==1)
 456         {
 457             u=line(q,q.add(q.sub(p).rotleft()));
 458             v=u;
 459             return 1;
 460         }
 461         double d=p.distance(q);
 462         double l=sqr(r)/d;
 463         double h=sqrt(sqr(r)-sqr(l));
 464         u=line(q,p.add(q.sub(p).trunc(l).add(q.sub(p).rotleft().trunc(h))));
 465         v=line(q,p.add(q.sub(p).trunc(l).add(q.sub(p).rotright().trunc(h))));
 466         return 2;
 467     }
 468     double areacircle(circle v)
 469     {
 470         int rel=relationcircle(v);
 471         if (rel>=4)return 0.0;
 472         if (rel<=2)return min(area(),v.area());
 473         double d=p.distance(v.p);
 474         double hf=(r+v.r+d)/2.0;
 475         double ss=2*sqrt(hf*(hf-r)*(hf-v.r)*(hf-d));
 476         double a1=acos((r*r+d*d-v.r*v.r)/(2.0*r*d));
 477         a1=a1*r*r;
 478         double a2=acos((v.r*v.r+d*d-r*r)/(2.0*v.r*d));
 479         a2=a2*v.r*v.r;
 480         return a1+a2-ss;
 481     }
 482     double areatriangle(point a,point b)
 483     {
 484         if (dblcmp(p.sub(a).det(p.sub(b))==0))return 0.0;
 485         point q[5];
 486         int len=0;
 487         q[len++]=a;
 488         line l(a,b);
 489         point p1,p2;
 490         if (pointcrossline(l,q[1],q[2])==2)
 491         {
 492             if (dblcmp(a.sub(q[1]).dot(b.sub(q[1])))<0)q[len++]=q[1];
 493             if (dblcmp(a.sub(q[2]).dot(b.sub(q[2])))<0)q[len++]=q[2];
 494         }
 495         q[len++]=b;
 496         if (len==4&&(dblcmp(q[0].sub(q[1]).dot(q[2].sub(q[1])))>0))swap(q[1],q[2]);
 497         double res=0;
 498         int i;
 499         for (i=0;i<len-1;i++)
 500         {
 501             if (relation(q[i])==0||relation(q[i+1])==0)
 502             {
 503                 double arg=p.rad(q[i],q[i+1]);
 504                 res+=r*r*arg/2.0;
 505             }
 506             else
 507             {
 508                 res+=fabs(q[i].sub(p).det(q[i+1].sub(p))/2.0);
 509             }
 510         }
 511         return res;
 512     }
 513 };
 514 struct polygon
 515 {
 516     int n;
 517     point p[maxp];
 518     line l[maxp];
 519     void input()
 520     {
 521         n=4;
 522         for (int i=0;i<n;i++)
 523         {
 524             p[i].input();
 525         }
 526     }
 527     void add(point q)
 528     {
 529         p[n++]=q;
 530     }
 531     void getline()
 532     {
 533         for (int i=0;i<n;i++)
 534         {
 535             l[i]=line(p[i],p[(i+1)%n]);
 536         }
 537     }
 538     struct cmp
 539     {
 540         point p;
 541         cmp(const point &p0){p=p0;}
 542         bool operator()(const point &aa,const point &bb)
 543         {
 544             point a=aa,b=bb;
 545             int d=dblcmp(a.sub(p).det(b.sub(p)));
 546             if (d==0)
 547             {
 548                 return dblcmp(a.distance(p)-b.distance(p))<0;
 549             }
 550             return d>0;
 551         }
 552     };
 553     void norm()
 554     {
 555         point mi=p[0];
 556         for (int i=1;i<n;i++)mi=min(mi,p[i]);
 557         sort(p,p+n,cmp(mi));
 558     }
 559     void getconvex(polygon &convex)
 560     {
 561         int i,j,k;
 562         sort(p,p+n);
 563         convex.n=n;
 564         for (i=0;i<min(n,2);i++)
 565         {
 566             convex.p[i]=p[i];
 567         }
 568         if (n<=2)return;
 569         int &top=convex.n;
 570         top=1;
 571         for (i=2;i<n;i++)
 572         {
 573             while (top&&convex.p[top].sub(p[i]).det(convex.p[top-1].sub(p[i]))<=0)
 574                 top--;
 575             convex.p[++top]=p[i];
 576         }
 577         int temp=top;
 578         convex.p[++top]=p[n-2];
 579         for (i=n-3;i>=0;i--)
 580         {
 581             while (top!=temp&&convex.p[top].sub(p[i]).det(convex.p[top-1].sub(p[i]))<=0)
 582                 top--;
 583             convex.p[++top]=p[i];
 584         }
 585     }
 586     bool isconvex()
 587     {
 588         bool s[3];
 589         memset(s,0,sizeof(s));
 590         int i,j,k;
 591         for (i=0;i<n;i++)
 592         {
 593             j=(i+1)%n;
 594             k=(j+1)%n;
 595             s[dblcmp(p[j].sub(p[i]).det(p[k].sub(p[i])))+1]=1;
 596             if (s[0]&&s[2])return 0;
 597         }
 598         return 1;
 599     }
 600     //3 点上
 601     //2 边上
 602     //1 内部
 603     //0 外部
 604     int relationpoint(point q)
 605     {
 606         int i,j;
 607         for (i=0;i<n;i++)
 608         {
 609             if (p[i]==q)return 3;
 610         }
 611         getline();
 612         for (i=0;i<n;i++)
 613         {
 614             if (l[i].pointonseg(q))return 2;
 615         }
 616         int cnt=0;
 617         for (i=0;i<n;i++)
 618         {
 619             j=(i+1)%n;
 620             int k=dblcmp(q.sub(p[j]).det(p[i].sub(p[j])));
 621             int u=dblcmp(p[i].y-q.y);
 622             int v=dblcmp(p[j].y-q.y);
 623             if (k>0&&u<0&&v>=0)cnt++;
 624             if (k<0&&v<0&&u>=0)cnt--;
 625         }
 626         return cnt!=0;
 627     }
 628     //1 在多边形内长度为正
 629     //2 相交或与边平行
 630     //0 无任何交点
 631     int relationline(line u)
 632     {
 633         int i,j,k=0;
 634         getline();
 635         for (i=0;i<n;i++)
 636         {
 637             if (l[i].segcrossseg(u)==2)return 1;
 638             if (l[i].segcrossseg(u)==1)k=1;
 639         }
 640         if (!k)return 0;
 641         vector<point>vp;
 642         for (i=0;i<n;i++)
 643         {
 644             if (l[i].segcrossseg(u))
 645             {
 646                 if (l[i].parallel(u))
 647                 {
 648                     vp.pb(u.a);
 649                     vp.pb(u.b);
 650                     vp.pb(l[i].a);
 651                     vp.pb(l[i].b);
 652                     continue;
 653                 }
 654                 vp.pb(l[i].crosspoint(u));
 655             }
 656         }
 657         sort(vp.begin(),vp.end());
 658         int sz=vp.size();
 659         for (i=0;i<sz-1;i++)
 660         {
 661             point mid=vp[i].add(vp[i+1]).div(2);
 662             if (relationpoint(mid)==1)return 1;
 663         }
 664         return 2;
 665     }
 666     //直线u切割凸多边形左侧
 667     //注意直线方向
 668     void convexcut(line u,polygon &po)
 669     {
 670         int i,j,k;
 671         int &top=po.n;
 672         top=0;
 673         for (i=0;i<n;i++)
 674         {
 675             int d1=dblcmp(p[i].sub(u.a).det(u.b.sub(u.a)));
 676             int d2=dblcmp(p[(i+1)%n].sub(u.a).det(u.b.sub(u.a)));
 677             if (d1>=0)po.p[top++]=p[i];
 678             if (d1*d2<0)po.p[top++]=u.crosspoint(line(p[i],p[(i+1)%n]));
 679         }
 680     }
 681     double getcircumference()
 682     {
 683         double sum=0;
 684         int i;
 685         for (i=0;i<n;i++)
 686         {
 687             sum+=p[i].distance(p[(i+1)%n]);
 688         }
 689         return sum;
 690     }
 691     double getarea()
 692     {
 693         double sum=0;
 694         int i;
 695         for (i=0;i<n;i++)
 696         {
 697             sum+=p[i].det(p[(i+1)%n]);
 698         }
 699         return fabs(sum)/2;
 700     }
 701     bool getdir()//1代表逆时针 0代表顺时针
 702     {
 703         double sum=0;
 704         int i;
 705         for (i=0;i<n;i++)
 706         {
 707             sum+=p[i].det(p[(i+1)%n]);
 708         }
 709         if (dblcmp(sum)>0)return 1;
 710         return 0;
 711     }
 712     point getbarycentre()
 713     {
 714         point ret(0,0);
 715         double area=0;
 716         int i;
 717         for (i=1;i<n-1;i++)
 718         {
 719             double tmp=p[i].sub(p[0]).det(p[i+1].sub(p[0]));
 720             if (dblcmp(tmp)==0)continue;
 721             area+=tmp;
 722             ret.x+=(p[0].x+p[i].x+p[i+1].x)/3*tmp;
 723             ret.y+=(p[0].y+p[i].y+p[i+1].y)/3*tmp;
 724         }
 725         if (dblcmp(area))ret=ret.div(area);
 726         return ret;
 727     }
 728     double areaintersection(polygon po)
 729     {
 730     }
 731     double areaunion(polygon po)
 732     {
 733         return getarea()+po.getarea()-areaintersection(po);
 734     }
 735     double areacircle(circle c)
 736     {
 737         int i,j,k,l,m;
 738         double ans=0;
 739         for (i=0;i<n;i++)
 740         {
 741             int j=(i+1)%n;
 742             if (dblcmp(p[j].sub(c.p).det(p[i].sub(c.p)))>=0)
 743             {
 744                 ans+=c.areatriangle(p[i],p[j]);
 745             }
 746             else
 747             {
 748                 ans-=c.areatriangle(p[i],p[j]);
 749             }
 750         }
 751         return fabs(ans);
 752     }
 753     //多边形和圆关系
 754     //0 一部分在圆外
 755     //1 与圆某条边相切
 756     //2 完全在圆内
 757     int relationcircle(circle c)
 758     {
 759         getline();
 760         int i,x=2;
 761         if (relationpoint(c.p)!=1)return 0;
 762         for (i=0;i<n;i++)
 763         {
 764             if (c.relationseg(l[i])==2)return 0;
 765             if (c.relationseg(l[i])==1)x=1;
 766         }
 767         return x;
 768     }
 769     void find(int st,point tri[],circle &c)
 770     {
 771         if (!st)
 772         {
 773             c=circle(point(0,0),-2);
 774         }
 775         if (st==1)
 776         {
 777             c=circle(tri[0],0);
 778         }
 779         if (st==2)
 780         {
 781             c=circle(tri[0].add(tri[1]).div(2),tri[0].distance(tri[1])/2.0);
 782         }
 783         if (st==3)
 784         {
 785             c=circle(tri[0],tri[1],tri[2]);
 786         }
 787     }
 788     void solve(int cur,int st,point tri[],circle &c)
 789     {
 790         find(st,tri,c);
 791         if (st==3)return;
 792         int i;
 793         for (i=0;i<cur;i++)
 794         {
 795             if (dblcmp(p[i].distance(c.p)-c.r)>0)
 796             {
 797                 tri[st]=p[i];
 798                 solve(i,st+1,tri,c);
 799             }
 800         }
 801     }
 802     circle mincircle()//点集最小圆覆盖
 803     {
 804         random_shuffle(p,p+n);
 805         point tri[4];
 806         circle c;
 807         solve(n,0,tri,c);
 808         return c;
 809     }
 810     int circlecover(double r)//单位圆覆盖
 811     {
 812         int ans=0,i,j;
 813         vector<pair<double,int> >v;
 814         for (i=0;i<n;i++)
 815         {
 816             v.clear();
 817             for (j=0;j<n;j++)if (i!=j)
 818             {
 819                 point q=p[i].sub(p[j]);
 820                 double d=q.len();
 821                 if (dblcmp(d-2*r)<=0)
 822                 {
 823                     double arg=atan2(q.y,q.x);
 824                     if (dblcmp(arg)<0)arg+=2*pi;
 825                     double t=acos(d/(2*r));
 826                     v.push_back(make_pair(arg-t+2*pi,-1));
 827                     v.push_back(make_pair(arg+t+2*pi,1));
 828                 }
 829             }
 830             sort(v.begin(),v.end());
 831             int cur=0;
 832             for (j=0;j<v.size();j++)
 833             {
 834                 if (v[j].second==-1)++cur;
 835                 else --cur;
 836                 ans=max(ans,cur);
 837             }
 838         }
 839         return ans+1;
 840     }
 841     int pointinpolygon(point q)//点在凸多边形内部的判定
 842     {
 843         if (getdir())reverse(p,p+n);
 844         if (dblcmp(q.sub(p[0]).det(p[n-1].sub(p[0])))==0)
 845         {
 846             if (line(p[n-1],p[0]).pointonseg(q))return n-1;
 847             return -1;
 848         }
 849         int low=1,high=n-2,mid;
 850         while (low<=high)
 851         {
 852             mid=(low+high)>>1;
 853             if (dblcmp(q.sub(p[0]).det(p[mid].sub(p[0])))>=0&&dblcmp(q.sub(p[0]).det(p[mid+1].sub(p[0])))<0)
 854             {
 855                 polygon c;
 856                 c.p[0]=p[mid];
 857                 c.p[1]=p[mid+1];
 858                 c.p[2]=p[0];
 859                 c.n=3;
 860                 if (c.relationpoint(q))return mid;
 861                 return -1;
 862             }
 863             if (dblcmp(q.sub(p[0]).det(p[mid].sub(p[0])))>0)
 864             {
 865                 low=mid+1;
 866             }
 867             else
 868             {
 869                 high=mid-1;
 870             }
 871         }
 872         return -1;
 873     }
 874 };
 875 struct polygons
 876 {
 877     vector<polygon>p;
 878     polygons()
 879     {
 880         p.clear();
 881     }
 882     void clear()
 883     {
 884         p.clear();
 885     }
 886     void push(polygon q)
 887     {
 888         if (dblcmp(q.getarea()))p.pb(q);
 889     }
 890     vector<pair<double,int> >e;
 891     void ins(point s,point t,point X,int i)
 892     {
 893         double r=fabs(t.x-s.x)>eps?(X.x-s.x)/(t.x-s.x):(X.y-s.y)/(t.y-s.y);
 894         r=min(r,1.0);r=max(r,0.0);
 895         e.pb(mp(r,i));
 896     }
 897     double polyareaunion()
 898     {
 899         double ans=0.0;
 900         int c0,c1,c2,i,j,k,w;
 901         for (i=0;i<p.size();i++)
 902         {
 903             if (p[i].getdir()==0)reverse(p[i].p,p[i].p+p[i].n);
 904         }
 905         for (i=0;i<p.size();i++)
 906         {
 907             for (k=0;k<p[i].n;k++)
 908             {
 909                 point &s=p[i].p[k],&t=p[i].p[(k+1)%p[i].n];
 910                 if (!dblcmp(s.det(t)))continue;
 911                 e.clear();
 912                 e.pb(mp(0.0,1));
 913                 e.pb(mp(1.0,-1));
 914                 for (j=0;j<p.size();j++)if (i!=j)
 915                 {
 916                     for (w=0;w<p[j].n;w++)
 917                     {
 918                         point a=p[j].p[w],b=p[j].p[(w+1)%p[j].n],c=p[j].p[(w-1+p[j].n)%p[j].n];
 919                         c0=dblcmp(t.sub(s).det(c.sub(s)));
 920                         c1=dblcmp(t.sub(s).det(a.sub(s)));
 921                         c2=dblcmp(t.sub(s).det(b.sub(s)));
 922                         if (c1*c2<0)ins(s,t,line(s,t).crosspoint(line(a,b)),-c2);
 923                         else if (!c1&&c0*c2<0)ins(s,t,a,-c2);
 924                         else if (!c1&&!c2)
 925                         {
 926                             int c3=dblcmp(t.sub(s).det(p[j].p[(w+2)%p[j].n].sub(s)));
 927                             int dp=dblcmp(t.sub(s).dot(b.sub(a)));
 928                             if (dp&&c0)ins(s,t,a,dp>0?c0*((j>i)^(c0<0)):-(c0<0));
 929                             if (dp&&c3)ins(s,t,b,dp>0?-c3*((j>i)^(c3<0)):c3<0);
 930                         }
 931                     }
 932                 }
 933                 sort(e.begin(),e.end());
 934                 int ct=0;
 935                 double tot=0.0,last;
 936                 for (j=0;j<e.size();j++)
 937                 {
 938                     if (ct==2)tot+=e[j].first-last;
 939                     ct+=e[j].second;
 940                     last=e[j].first;
 941                 }
 942                 ans+=s.det(t)*tot;
 943             }
 944         }
 945         return fabs(ans)*0.5;
 946     }
 947 };
 948 const int maxn=500;
 949 struct circles
 950 {
 951     circle c[maxn];
 952     double ans[maxn];//ans[i]表示被覆盖了i次的面积
 953     double pre[maxn];
 954     int n;
 955     circles(){}
 956     void add(circle cc)
 957     {
 958         c[n++]=cc;
 959     }
 960     bool inner(circle x,circle y)
 961     {
 962         if (x.relationcircle(y)!=1)return 0;
 963         return dblcmp(x.r-y.r)<=0?1:0;
 964     }
 965     void init_or()//圆的面积并去掉内含的圆
 966     {
 967         int i,j,k=0;
 968         bool mark[maxn]={0};
 969         for (i=0;i<n;i++)
 970         {
 971             for (j=0;j<n;j++)if (i!=j&&!mark[j])
 972             {
 973                 if ((c[i]==c[j])||inner(c[i],c[j]))break;
 974             }
 975             if (j<n)mark[i]=1;
 976         }
 977         for (i=0;i<n;i++)if (!mark[i])c[k++]=c[i];
 978         n=k;
 979     }
 980     void init_and()//圆的面积交去掉内含的圆
 981     {
 982         int i,j,k=0;
 983         bool mark[maxn]={0};
 984         for (i=0;i<n;i++)
 985         {
 986             for (j=0;j<n;j++)if (i!=j&&!mark[j])
 987             {
 988                 if ((c[i]==c[j])||inner(c[j],c[i]))break;
 989             }
 990             if (j<n)mark[i]=1;
 991         }
 992         for (i=0;i<n;i++)if (!mark[i])c[k++]=c[i];
 993         n=k;
 994     }
 995     double areaarc(double th,double r)
 996     {
 997         return 0.5*sqr(r)*(th-sin(th));
 998     }
 999     void getarea()
1000     {
1001         int i,j,k;
1002         memset(ans,0,sizeof(ans));
1003         vector<pair<double,int> >v;
1004         for (i=0;i<n;i++)
1005         {
1006             v.clear();
1007             v.push_back(make_pair(-pi,1));
1008             v.push_back(make_pair(pi,-1));
1009             for (j=0;j<n;j++)if (i!=j)
1010             {
1011                 point q=c[j].p.sub(c[i].p);
1012                 double ab=q.len(),ac=c[i].r,bc=c[j].r;
1013                 if (dblcmp(ab+ac-bc)<=0)
1014                 {
1015                     v.push_back(make_pair(-pi,1));
1016                     v.push_back(make_pair(pi,-1));
1017                     continue;
1018                 }
1019                 if (dblcmp(ab+bc-ac)<=0)continue;
1020                 if (dblcmp(ab-ac-bc)>0) continue;
1021                 double th=atan2(q.y,q.x),fai=acos((ac*ac+ab*ab-bc*bc)/(2.0*ac*ab));
1022                 double a0=th-fai;
1023                 if (dblcmp(a0+pi)<0)a0+=2*pi;
1024                 double a1=th+fai;
1025                 if (dblcmp(a1-pi)>0)a1-=2*pi;
1026                 if (dblcmp(a0-a1)>0)
1027                 {
1028                     v.push_back(make_pair(a0,1));
1029                     v.push_back(make_pair(pi,-1));
1030                     v.push_back(make_pair(-pi,1));
1031                     v.push_back(make_pair(a1,-1));
1032                 }
1033                 else
1034                 {
1035                     v.push_back(make_pair(a0,1));
1036                     v.push_back(make_pair(a1,-1));
1037                 }
1038             }
1039             sort(v.begin(),v.end());
1040             int cur=0;
1041             for (j=0;j<v.size();j++)
1042             {
1043                 if (cur&&dblcmp(v[j].first-pre[cur]))
1044                 {
1045                     ans[cur]+=areaarc(v[j].first-pre[cur],c[i].r);
1046                     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])).det(point(c[i].p.x+c[i].r*cos(v[j].first),c[i].p.y+c[i].r*sin(v[j].first)));
1047                 }
1048                 cur+=v[j].second;
1049                 pre[cur]=v[j].first;
1050             }
1051         }
1052         for (i=1;i<=n;i++)
1053         {
1054             ans[i]-=ans[i+1];
1055         }
1056     }
1057 };
1058 struct halfplane:public line
1059 {
1060     double angle;
1061     halfplane(){}
1062     //表示向量 a->b逆时针(左侧)的半平面
1063     halfplane(point _a,point _b)
1064     {
1065         a=_a;
1066         b=_b;
1067     }
1068     halfplane(line v)
1069     {
1070         a=v.a;
1071         b=v.b;
1072     }
1073     void calcangle()
1074     {
1075         angle=atan2(b.y-a.y,b.x-a.x);
1076     }
1077     bool operator<(const halfplane &b)const
1078     {
1079         return angle<b.angle;
1080     }
1081 };
1082 struct halfplanes
1083 {
1084     int n;
1085     halfplane hp[maxp];
1086     point p[maxp];
1087     int que[maxp];
1088     int st,ed;
1089     void push(halfplane tmp)
1090     {
1091         hp[n++]=tmp;
1092     }
1093     void unique()
1094     {
1095         int m=1,i;
1096         for (i=1;i<n;i++)
1097         {
1098             if (dblcmp(hp[i].angle-hp[i-1].angle))hp[m++]=hp[i];
1099             else if (dblcmp(hp[m-1].b.sub(hp[m-1].a).det(hp[i].a.sub(hp[m-1].a))>0))hp[m-1]=hp[i];
1100         }
1101         n=m;
1102     }
1103     bool halfplaneinsert()
1104     {
1105         int i;
1106         for (i=0;i<n;i++)hp[i].calcangle();
1107         sort(hp,hp+n);
1108         unique();
1109         que[st=0]=0;
1110         que[ed=1]=1;
1111         p[1]=hp[0].crosspoint(hp[1]);
1112         for (i=2;i<n;i++)
1113         {
1114             while (st<ed&&dblcmp((hp[i].b.sub(hp[i].a).det(p[ed].sub(hp[i].a))))<0)ed--;
1115             while (st<ed&&dblcmp((hp[i].b.sub(hp[i].a).det(p[st+1].sub(hp[i].a))))<0)st++;
1116             que[++ed]=i;
1117             if (hp[i].parallel(hp[que[ed-1]]))return false;
1118             p[ed]=hp[i].crosspoint(hp[que[ed-1]]);
1119         }
1120         while (st<ed&&dblcmp(hp[que[st]].b.sub(hp[que[st]].a).det(p[ed].sub(hp[que[st]].a)))<0)ed--;
1121         while (st<ed&&dblcmp(hp[que[ed]].b.sub(hp[que[ed]].a).det(p[st+1].sub(hp[que[ed]].a)))<0)st++;
1122         if (st+1>=ed)return false;
1123         return true;
1124     }
1125     void getconvex(polygon &con)
1126     {
1127         p[st]=hp[que[st]].crosspoint(hp[que[ed]]);
1128         con.n=ed-st+1;
1129         int j=st,i=0;
1130         for (;j<=ed;i++,j++)
1131         {
1132             con.p[i]=p[j];
1133         }
1134     }
1135 };
1136 struct point3
1137 {
1138     double x,y,z;
1139     point3(){}
1140     point3(double _x,double _y,double _z):
1141     x(_x),y(_y),z(_z){};
1142     void input()
1143     {
1144         scanf("%lf%lf%lf",&x,&y,&z);
1145     }
1146     void output()
1147     {
1148         printf("%.2lf %.2lf %.2lf\n",x,y,z);
1149     }
1150     bool operator==(point3 a)
1151     {
1152         return dblcmp(a.x-x)==0&&dblcmp(a.y-y)==0&&dblcmp(a.z-z)==0;
1153     }
1154     bool operator<(point3 a)const
1155     {
1156         return dblcmp(a.x-x)==0?dblcmp(y-a.y)==0?dblcmp(z-a.z)<0:y<a.y:x<a.x;
1157     }
1158     double len()
1159     {
1160         return sqrt(len2());
1161     }
1162     double len2()
1163     {
1164         return x*x+y*y+z*z;
1165     }
1166     double distance(point3 p)
1167     {
1168         return sqrt((p.x-x)*(p.x-x)+(p.y-y)*(p.y-y)+(p.z-z)*(p.z-z));
1169     }
1170     point3 add(point3 p)
1171     {
1172         return point3(x+p.x,y+p.y,z+p.z);
1173     }
1174     point3 sub(point3 p)
1175     {
1176         return point3(x-p.x,y-p.y,z-p.z);
1177     }
1178     point3 mul(double d)
1179     {
1180         return point3(x*d,y*d,z*d);
1181     }
1182     point3 div(double d)
1183     {
1184         return point3(x/d,y/d,z/d);
1185     }
1186     double dot(point3 p)
1187     {
1188         return x*p.x+y*p.y+z*p.z;
1189     }
1190     point3 det(point3 p)
1191     {
1192         return point3(y*p.z-p.y*z,p.x*z-x*p.z,x*p.y-p.x*y);
1193     }
1194     double rad(point3 a,point3 b)
1195     {
1196         point3 p=(*this);
1197         return acos(a.sub(p).dot(b.sub(p))/(a.distance(p)*b.distance(p)));
1198     }
1199     point3 trunc(double r)
1200     {
1201         r/=len();
1202         return point3(x*r,y*r,z*r);
1203     }
1204     point3 rotate(point3 o,double r)
1205     {
1206     }
1207 };
1208 struct line3
1209 {
1210     point3 a,b;
1211     line3(){}
1212     line3(point3 _a,point3 _b)
1213     {
1214         a=_a;
1215         b=_b;
1216     }
1217     bool operator==(line3 v)
1218     {
1219         return (a==v.a)&&(b==v.b);
1220     }
1221     void input()
1222     {
1223         a.input();
1224         b.input();
1225     }
1226     double length()
1227     {
1228         return a.distance(b);
1229     }
1230     bool pointonseg(point3 p)
1231     {
1232         return dblcmp(p.sub(a).det(p.sub(b)).len())==0&&dblcmp(a.sub(p).dot(b.sub(p)))<=0;
1233     }
1234     double dispointtoline(point3 p)
1235     {
1236         return b.sub(a).det(p.sub(a)).len()/a.distance(b);
1237     }
1238     double dispointtoseg(point3 p)
1239     {
1240         if (dblcmp(p.sub(b).dot(a.sub(b)))<0||dblcmp(p.sub(a).dot(b.sub(a)))<0)
1241         {
1242             return min(p.distance(a),p.distance(b));
1243         }
1244         return dispointtoline(p);
1245     }
1246     point3 lineprog(point3 p)
1247     {
1248         return a.add(b.sub(a).trunc(b.sub(a).dot(p.sub(a))/b.distance(a)));
1249     }
1250     point3 rotate(point3 p,double ang)//p绕此向量逆时针arg角度
1251     {
1252         if (dblcmp((p.sub(a).det(p.sub(b)).len()))==0)return p;
1253         point3 f1=b.sub(a).det(p.sub(a));
1254         point3 f2=b.sub(a).det(f1);
1255         double len=fabs(a.sub(p).det(b.sub(p)).len()/a.distance(b));
1256         f1=f1.trunc(len);f2=f2.trunc(len);
1257         point3 h=p.add(f2);
1258         point3 pp=h.add(f1);
1259         return h.add((p.sub(h)).mul(cos(ang*1.0))).add((pp.sub(h)).mul(sin(ang*1.0)));
1260     }
1261 };
1262 struct plane
1263 {
1264     point3 a,b,c,o;
1265     plane(){}
1266     plane(point3 _a,point3 _b,point3 _c)
1267     {
1268         a=_a;
1269         b=_b;
1270         c=_c;
1271         o=pvec();
1272     }
1273     plane(double _a,double _b,double _c,double _d)
1274     {
1275         //ax+by+cz+d=0
1276         o=point3(_a,_b,_c);
1277         if (dblcmp(_a)!=0)
1278         {
1279             a=point3((-_d-_c-_b)/_a,1,1);
1280         }
1281         else if (dblcmp(_b)!=0)
1282         {
1283             a=point3(1,(-_d-_c-_a)/_b,1);
1284         }
1285         else if (dblcmp(_c)!=0)
1286         {
1287             a=point3(1,1,(-_d-_a-_b)/_c);
1288         }
1289     }
1290     void input()
1291     {
1292         a.input();
1293         b.input();
1294         c.input();
1295         o=pvec();
1296     }
1297     point3 pvec()
1298     {
1299         return b.sub(a).det(c.sub(a));
1300     }
1301     bool pointonplane(point3 p)//点是否在平面上
1302     {
1303         return dblcmp(p.sub(a).dot(o))==0;
1304     }
1305     //0 不在
1306     //1 在边界上
1307     //2 在内部
1308     int pointontriangle(point3 p)//点是否在空间三角形abc上
1309     {
1310         if (!pointonplane(p))return 0;
1311         double s=a.sub(b).det(c.sub(b)).len();
1312         double s1=p.sub(a).det(p.sub(b)).len();
1313         double s2=p.sub(a).det(p.sub(c)).len();
1314         double s3=p.sub(b).det(p.sub(c)).len();
1315         if (dblcmp(s-s1-s2-s3))return 0;
1316         if (dblcmp(s1)&&dblcmp(s2)&&dblcmp(s3))return 2;
1317         return 1;
1318     }
1319     //判断两平面关系
1320     //0 相交
1321     //1 平行但不重合
1322     //2 重合
1323     bool relationplane(plane f)
1324     {
1325         if (dblcmp(o.det(f.o).len()))return 0;
1326         if (pointonplane(f.a))return 2;
1327         return 1;
1328     }
1329     double angleplane(plane f)//两平面夹角
1330     {
1331         return acos(o.dot(f.o)/(o.len()*f.o.len()));
1332     }
1333     double dispoint(point3 p)//点到平面距离
1334     {
1335         return fabs(p.sub(a).dot(o)/o.len());
1336     }
1337     point3 pttoplane(point3 p)//点到平面最近点
1338     {
1339         line3 u=line3(p,p.add(o));
1340         crossline(u,p);
1341         return p;
1342     }
1343     int crossline(line3 u,point3 &p)//平面和直线的交点
1344     {
1345         double x=o.dot(u.b.sub(a));
1346         double y=o.dot(u.a.sub(a));
1347         double d=x-y;
1348         if (dblcmp(fabs(d))==0)return 0;
1349         p=u.a.mul(x).sub(u.b.mul(y)).div(d);
1350         return 1;
1351     }
1352     int crossplane(plane f,line3 &u)//平面和平面的交线
1353     {
1354         point3 oo=o.det(f.o);
1355         point3 v=o.det(oo);
1356         double d=fabs(f.o.dot(v));
1357         if (dblcmp(d)==0)return 0;
1358         point3 q=a.add(v.mul(f.o.dot(f.a.sub(a))/d));
1359         u=line3(q,q.add(oo));
1360         return 1;
1361     }
1362 };
View Code

 

2.

  1 #include<vector>
  2 #include<list>
  3 #include<map>
  4 #include<set>
  5 #include<deque>
  6 #include<queue>
  7 #include<stack>
  8 #include<bitset>
  9 #include<algorithm>
 10 #include<functional>
 11 #include<numeric>
 12 #include<utility>
 13 #include<iostream>
 14 #include<sstream>
 15 #include<iomanip>
 16 #include<cstdio>
 17 #include<cmath>
 18 #include<cstdlib>
 19 #include<cctype>
 20 #include<string>
 21 #include<cstring>
 22 #include<cstdio>
 23 #include<cmath>
 24 #include<cstdlib>
 25 #include<ctime>
 26 #include<climits>
 27 #include<complex>
 28 #define mp make_pair
 29 #define pb push_back
 30 using namespace std;
 31 const double eps=1e-8;//精度
 32 const double pi=acos(-1.0);//π
 33 const double inf=1e20;//无穷大
 34 const int maxp=1111;//最大点数
 35 /*
 36     判断d是否在精度内等于0
 37 */
 38 int dblcmp(double d)
 39 {
 40     if (fabs(d)<eps)return 0;
 41     return d>eps?1:-1;
 42 }
 43 /*
 44     求x的平方
 45 */
 46 inline double sqr(double x){return x*x;}
 47 /*
 48     点/向量
 49 */
 50 struct point
 51 {
 52     double x,y;
 53     point(){}
 54     point(double _x,double _y):x(_x),y(_y){};
 55     //读入一个点
 56     void input()
 57     {
 58         scanf("%lf%lf",&x,&y);
 59     }
 60     //输出一个点
 61     void output()
 62     {
 63         printf("%.2f %.2f\n",x,y);
 64     }
 65     //判断两点是否相等
 66     bool operator==(point a)const
 67     {
 68         return dblcmp(a.x-x)==0&&dblcmp(a.y-y)==0;
 69     }
 70     //判断两点大小
 71     bool operator<(point a)const
 72     {
 73         return dblcmp(a.x-x)==0?dblcmp(y-a.y)<0:x<a.x;
 74     }
 75     //点到源点的距离/向量的长度
 76     double len()
 77     {
 78         return hypot(x,y);
 79     }
 80     //点到源点距离的平方
 81     double len2()
 82     {
 83         return x*x+y*y;
 84     }
 85     //两点间的距离
 86     double distance(point p)
 87     {
 88         return hypot(x-p.x,y-p.y);
 89     }
 90     //向量加
 91     point add(point p)
 92     {
 93         return point(x+p.x,y+p.y);
 94     }
 95     //向量减
 96     point sub(point p)
 97     {
 98         return point(x-p.x,y-p.y);
 99     }
100     //向量乘
101     point mul(double b)
102     {
103         return point(x*b,y*b);
104     }
105     //向量除
106     point div(double b)
107     {
108         return point(x/b,y/b);
109     }
110     //点乘
111     double dot(point p)
112     {
113         return x*p.x+y*p.y;
114     }
115     //叉乘
116     double det(point p)
117     {
118         return x*p.y-y*p.x;
119     }
120     //XXXXXXX
121     double rad(point a,point b)
122     {
123         point p=*this;
124         return fabs(atan2(fabs(a.sub(p).det(b.sub(p))),a.sub(p).dot(b.sub(p))));
125     }
126     //截取长度r
127     point trunc(double r)
128     {
129         double l=len();
130         if (!dblcmp(l))return *this;
131         r/=l;
132         return point(x*r,y*r);
133     }
134     //左转90度
135     point rotleft()
136     {
137         return point(-y,x);
138     }
139     //右转90度
140     point rotright()
141     {
142         return point(y,-x);
143     }
144     //绕点p逆时针旋转angle角度
145     point rotate(point p,double angle)
146     {
147         point v=this->sub(p);
148         double c=cos(angle),s=sin(angle);
149         return point(p.x+v.x*c-v.y*s,p.y+v.x*s+v.y*c);
150     }
151 };
152 /*
153     线段/直线
154 */
155 struct line
156 {
157     point a,b;
158     line(){}
159     line(point _a,point _b)
160     {
161         a=_a;
162         b=_b;
163     }
164     //判断线段相等
165     bool operator==(line v)
166     {
167         return (a==v.a)&&(b==v.b);
168     }
169     //点p做倾斜角为angle的射线
170     line(point p,double angle)
171     {
172         a=p;
173         if (dblcmp(angle-pi/2)==0)
174         {
175             b=a.add(point(0,1));
176         }
177         else
178         {
179             b=a.add(point(1,tan(angle)));
180         }
181     }
182     //直线一般式ax+by+c=0
183     line(double _a,double _b,double _c)
184     {
185         if (dblcmp(_a)==0)
186         {
187             a=point(0,-_c/_b);
188             b=point(1,-_c/_b);
189         }
190         else if (dblcmp(_b)==0)
191         {
192             a=point(-_c/_a,0);
193             b=point(-_c/_a,1);
194         }
195         else
196         {
197             a=point(0,-_c/_b);
198             b=point(1,(-_c-_a)/_b);
199         }
200     }
201     //读入一个线段
202     void input()
203     {
204         a.input();
205         b.input();
206     }
207     //校准线段两点
208     void adjust()
209     {
210         if (b<a)swap(a,b);
211     }
212     //线段长度
213     double length()
214     {
215         return a.distance(b);
216     }
217     //直线倾斜角 0<=angle<180
218     double angle()
219     {
220         double k=atan2(b.y-a.y,b.x-a.x);
221         if (dblcmp(k)<0)k+=pi;
222         if (dblcmp(k-pi)==0)k-=pi;
223         return k;
224     }
225     //点和线段关系
226     //1 在逆时针
227     //2 在顺时针
228     //3 平行
229     int relation(point p)
230     {
231         int c=dblcmp(p.sub(a).det(b.sub(a)));
232         if (c<0)return 1;
233         if (c>0)return 2;
234         return 3;
235     }
236     //点是否在线段上
237     bool pointonseg(point p)
238     {
239         return dblcmp(p.sub(a).det(b.sub(a)))==0&&dblcmp(p.sub(a).dot(p.sub(b)))<=0;
240     }
241     //两线是否平行
242     bool parallel(line v)
243     {
244         return dblcmp(b.sub(a).det(v.b.sub(v.a)))==0;
245     }
246     //线段和线段关系
247     //0 不相交
248     //1 非规范相交
249     //2 规范相交
250     int segcrossseg(line v)
251     {
252         int d1=dblcmp(b.sub(a).det(v.a.sub(a)));
253         int d2=dblcmp(b.sub(a).det(v.b.sub(a)));
254         int d3=dblcmp(v.b.sub(v.a).det(a.sub(v.a)));
255         int d4=dblcmp(v.b.sub(v.a).det(b.sub(v.a)));
256         if ((d1^d2)==-2&&(d3^d4)==-2)return 2;
257         return (d1==0&&dblcmp(v.a.sub(a).dot(v.a.sub(b)))<=0||
258                 d2==0&&dblcmp(v.b.sub(a).dot(v.b.sub(b)))<=0||
259                 d3==0&&dblcmp(a.sub(v.a).dot(a.sub(v.b)))<=0||
260                 d4==0&&dblcmp(b.sub(v.a).dot(b.sub(v.b)))<=0);
261     }
262     //线段和直线v关系
263     int linecrossseg(line v)//*this seg v line
264     {
265         int d1=dblcmp(b.sub(a).det(v.a.sub(a)));
266         int d2=dblcmp(b.sub(a).det(v.b.sub(a)));
267         if ((d1^d2)==-2)return 2;
268         return (d1==0||d2==0);
269     }
270     //直线和直线关系
271     //0 平行
272     //1 重合
273     //2 相交
274     int linecrossline(line v)
275     {
276         if ((*this).parallel(v))
277         {
278             return v.relation(a)==3;
279         }
280         return 2;
281     }
282     //求两线交点
283     point crosspoint(line v)
284     {
285         double a1=v.b.sub(v.a).det(a.sub(v.a));
286         double a2=v.b.sub(v.a).det(b.sub(v.a));
287         return point((a.x*a2-b.x*a1)/(a2-a1),(a.y*a2-b.y*a1)/(a2-a1));
288     }
289     //点p到直线的距离
290     double dispointtoline(point p)
291     {
292         return fabs(p.sub(a).det(b.sub(a)))/length();
293     }
294     //点p到线段的距离
295     double dispointtoseg(point p)
296     {
297         if (dblcmp(p.sub(b).dot(a.sub(b)))<0||dblcmp(p.sub(a).dot(b.sub(a)))<0)
298         {
299             return min(p.distance(a),p.distance(b));
300         }
301         return dispointtoline(p);
302     }
303     //XXXXXXXX
304     point lineprog(point p)
305     {
306         return a.add(b.sub(a).mul(b.sub(a).dot(p.sub(a))/b.sub(a).len2()));
307     }
308     //点p关于直线的对称点
309     point symmetrypoint(point p)
310     {
311         point q=lineprog(p);
312         return point(2*q.x-p.x,2*q.y-p.y);
313     }
314 };
315 
316 /*
317     多边形
318 */
319 struct polygon
320 {
321     int n;//点个数
322     point p[maxp];//顶点
323     line l[maxp];//324     //读入一个多边形
325     void input(int n=4)
326     {
327         for (int i=0;i<n;i++)
328         {
329             p[i].input();
330         }
331     }
332     //添加一个点
333     void add(point q)
334     {
335         p[n++]=q;
336     }
337     //取得边
338     void getline()
339     {
340         for (int i=0;i<n;i++)
341         {
342             l[i]=line(p[i],p[(i+1)%n]);
343         }
344     }
345     struct cmp
346     {
347         point p;
348         cmp(const point &p0){p=p0;}
349         bool operator()(const point &aa,const point &bb)
350         {
351             point a=aa,b=bb;
352             int d=dblcmp(a.sub(p).det(b.sub(p)));
353             if (d==0)
354             {
355                 return dblcmp(a.distance(p)-b.distance(p))<0;
356             }
357             return d>0;
358         }
359     };
360     void norm()
361     {
362         point mi=p[0];
363         for (int i=1;i<n;i++)mi=min(mi,p[i]);
364         sort(p,p+n,cmp(mi));
365     }
366     //求凸包存入多边形convex
367     void getconvex(polygon &convex)
368     {
369         int i,j,k;
370         sort(p,p+n);
371         convex.n=n;
372         for (i=0;i<min(n,2);i++)
373         {
374             convex.p[i]=p[i];
375         }
376         if (n<=2)return;
377         int &top=convex.n;
378         top=1;
379         for (i=2;i<n;i++)
380         {
381             while (top&&convex.p[top].sub(p[i]).det(convex.p[top-1].sub(p[i]))<=0)
382                 top--;
383             convex.p[++top]=p[i];
384         }
385         int temp=top;
386         convex.p[++top]=p[n-2];
387         for (i=n-3;i>=0;i--)
388         {
389             while (top!=temp&&convex.p[top].sub(p[i]).det(convex.p[top-1].sub(p[i]))<=0)
390                 top--;
391             convex.p[++top]=p[i];
392         }
393     }
394     //判断是否凸多边形
395     bool isconvex()
396     {
397         bool s[3];
398         memset(s,0,sizeof(s));
399         int i,j,k;
400         for (i=0;i<n;i++)
401         {
402             j=(i+1)%n;
403             k=(j+1)%n;
404             s[dblcmp(p[j].sub(p[i]).det(p[k].sub(p[i])))+1]=1;
405             if (s[0]&&s[2])return 0;
406         }
407         return 1;
408     }
409     //点与多边形关系
410     //0 外部
411     //1 内部
412     //2 边上
413     //3 点上
414     int relationpoint(point q)
415     {
416         int i,j;
417         for (i=0;i<n;i++)
418         {
419             if (p[i]==q)return 3;
420         }
421         getline();
422         for (i=0;i<n;i++)
423         {
424             if (l[i].pointonseg(q))return 2;
425         }
426         int cnt=0;
427         for (i=0;i<n;i++)
428         {
429             j=(i+1)%n;
430             int k=dblcmp(q.sub(p[j]).det(p[i].sub(p[j])));
431             int u=dblcmp(p[i].y-q.y);
432             int v=dblcmp(p[j].y-q.y);
433             if (k>0&&u<0&&v>=0)cnt++;
434             if (k<0&&v<0&&u>=0)cnt--;
435         }
436         return cnt!=0;
437     }
438     //线段与多边形关系
439     //0 无任何交点
440     //1 在多边形内长度为正
441     //2 相交或与边平行
442     int relationline(line u)
443     {
444         int i,j,k=0;
445         getline();
446         for (i=0;i<n;i++)
447         {
448             if (l[i].segcrossseg(u)==2)return 1;
449             if (l[i].segcrossseg(u)==1)k=1;
450         }
451         if (!k)return 0;
452         vector<point>vp;
453         for (i=0;i<n;i++)
454         {
455             if (l[i].segcrossseg(u))
456             {
457                 if (l[i].parallel(u))
458                 {
459                     vp.pb(u.a);
460                     vp.pb(u.b);
461                     vp.pb(l[i].a);
462                     vp.pb(l[i].b);
463                     continue;
464                 }
465                 vp.pb(l[i].crosspoint(u));
466             }
467         }
468         sort(vp.begin(),vp.end());
469         int sz=vp.size();
470         for (i=0;i<sz-1;i++)
471         {
472             point mid=vp[i].add(vp[i+1]).div(2);
473             if (relationpoint(mid)==1)return 1;
474         }
475         return 2;
476     }
477     //直线u切割凸多边形左侧
478     //注意直线方向
479     void convexcut(line u,polygon &po)
480     {
481         int i,j,k;
482         int &top=po.n;
483         top=0;
484         for (i=0;i<n;i++)
485         {
486             int d1=dblcmp(p[i].sub(u.a).det(u.b.sub(u.a)));
487             int d2=dblcmp(p[(i+1)%n].sub(u.a).det(u.b.sub(u.a)));
488             if (d1>=0)po.p[top++]=p[i];
489             if (d1*d2<0)po.p[top++]=u.crosspoint(line(p[i],p[(i+1)%n]));
490         }
491     }
492     //取得周长
493     double getcircumference()
494     {
495         double sum=0;
496         int i;
497         for (i=0;i<n;i++)
498         {
499             sum+=p[i].distance(p[(i+1)%n]);
500         }
501         return sum;
502     }
503     //取得面积
504     double getarea()
505     {
506         double sum=0;
507         int i;
508         for (i=0;i<n;i++)
509         {
510             sum+=p[i].det(p[(i+1)%n]);
511         }
512         return fabs(sum)/2;
513     }
514     bool getdir()//1代表逆时针 0代表顺时针
515     {
516         double sum=0;
517         int i;
518         for (i=0;i<n;i++)
519         {
520             sum+=p[i].det(p[(i+1)%n]);
521         }
522         if (dblcmp(sum)>0)return 1;
523         return 0;
524     }
525     //取得重心
526     point getbarycentre()
527     {
528         point ret(0,0);
529         double area=0;
530         int i;
531         for (i=1;i<n-1;i++)
532         {
533             double tmp=p[i].sub(p[0]).det(p[i+1].sub(p[0]));
534             if (dblcmp(tmp)==0)continue;
535             area+=tmp;
536             ret.x+=(p[0].x+p[i].x+p[i+1].x)/3*tmp;
537             ret.y+=(p[0].y+p[i].y+p[i+1].y)/3*tmp;
538         }
539         if (dblcmp(area))ret=ret.div(area);
540         return ret;
541     }
542     //点在凸多边形内部的判定
543     int pointinpolygon(point q)
544     {
545         if (getdir())reverse(p,p+n);
546         if (dblcmp(q.sub(p[0]).det(p[n-1].sub(p[0])))==0)
547         {
548             if (line(p[n-1],p[0]).pointonseg(q))return n-1;
549             return -1;
550         }
551         int low=1,high=n-2,mid;
552         while (low<=high)
553         {
554             mid=(low+high)>>1;
555             if (dblcmp(q.sub(p[0]).det(p[mid].sub(p[0])))>=0&&dblcmp(q.sub(p[0]).det(p[mid+1].sub(p[0])))<0)
556             {
557                 polygon c;
558                 c.p[0]=p[mid];
559                 c.p[1]=p[mid+1];
560                 c.p[2]=p[0];
561                 c.n=3;
562                 if (c.relationpoint(q))return mid;
563                 return -1;
564             }
565             if (dblcmp(q.sub(p[0]).det(p[mid].sub(p[0])))>0)
566             {
567                 low=mid+1;
568             }
569             else
570             {
571                 high=mid-1;
572             }
573         }
574         return -1;
575     }
576 };
View Code

 

3.

  1 计算几何模板
  2   1 #include<math.h>
  3   2 #define MAXN 1000
  4   3 #define offset 10000
  5   4 #define eps 1e-8
  6   5 #define PI acos(-1.0)//3.14159265358979323846
  7   6 //判断一个数是否为0,是则返回true,否则返回false
  8   7 #define zero(x)(((x)>0?(x):-(x))<eps)
  9   8 //返回一个数的符号,正数返回1,负数返回2,否则返回0
 10   9 #define _sign(x)((x)>eps?1:((x)<-eps?2:0))
 11  10 struct point 
 12  11 {
 13  12     double x,y;
 14  13 };
 15  14 struct line
 16  15 {
 17  16     point a,b;
 18  17 };//直线通过的两个点,而不是一般式的三个系数
 19  18 //求矢量[p0,p1],[p0,p2]的叉积
 20  19 //p0是顶点
 21  20 //若结果等于0,则这三点共线
 22  21 //若结果大于0,则p0p2在p0p1的逆时针方向
 23  22 //若结果小于0,则p0p2在p0p1的顺时针方向
 24  23 double xmult(point p1,point p2,point p0)
 25  24 {
 26  25     return(p1.x-p0.x)*(p2.y-p0.y)-(p2.x-p0.x)*(p1.y-p0.y);
 27  26 }
 28  27 //计算dotproduct(P1-P0).(P2-P0)
 29  28 double dmult(point p1,point p2,point p0)
 30  29 {
 31  30     return(p1.x-p0.x)*(p2.x-p0.x)+(p1.y-p0.y)*(p2.y-p0.y);
 32  31 }
 33  32 //两点距离
 34  33 double distance(point p1,point p2)
 35  34 {
 36  35     return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y));
 37  36 }
 38  37 //判三点共线
 39  38 int dots_inline(point p1,point p2,point p3)
 40  39 {
 41  40     return zero(xmult(p1,p2,p3));
 42  41 }
 43  42 //判点是否在线段上,包括端点
 44  43 int dot_online_in(point p,line l)
 45  44 {
 46  45     return zero(xmult(p,l.a,l.b))&&(l.a.x-p.x)*(l.b.x-p.x)<eps&&(l.a.y-p.y)*(l.b.y-p.y)<eps;
 47  46 }
 48  47 //判点是否在线段上,不包括端点
 49  48 int dot_online_ex(point p,line l)
 50  49 {
 51  50     return dot_online_in(p,l)&&(!zero(p.x-l.a.x)||!zero(p.y-l.a.y))&&(!zero(p.x-l.b.x)||!zero(p.y-l.b.y));
 52  51 }
 53  52 //判两点在线段同侧,点在线段上返回0
 54  53 int same_side(point p1,point p2,line l)
 55  54 {
 56  55     return xmult(l.a,p1,l.b)*xmult(l.a,p2,l.b)>eps;
 57  56 }
 58  57 //判两点在线段异侧,点在线段上返回0
 59  58 int opposite_side(point p1,point p2,line l)
 60  59 {
 61  60     return xmult(l.a,p1,l.b)*xmult(l.a,p2,l.b)<-eps;
 62  61 }
 63  62 //判两直线平行
 64  63 int parallel(line u,line v)
 65  64 {
 66  65     return zero((u.a.x-u.b.x)*(v.a.y-v.b.y)-(v.a.x-v.b.x)*(u.a.y-u.b.y));
 67  66 }
 68  67 //判两直线垂直
 69  68 int perpendicular(line u,line v)
 70  69 {
 71  70     return zero((u.a.x-u.b.x)*(v.a.x-v.b.x)+(u.a.y-u.b.y)*(v.a.y-v.b.y));
 72  71 }
 73  72 //判两线段相交,包括端点和部分重合
 74  73 int intersect_in(line u,line v)
 75  74 {
 76  75     if(!dots_inline(u.a,u.b,v.a)||!dots_inline(u.a,u.b,v.b))
 77  76         return!same_side(u.a,u.b,v)&&!same_side(v.a,v.b,u);
 78  77     return dot_online_in(u.a,v)||dot_online_in(u.b,v)||dot_online_in(v.a,u)||dot_online_in(v.b,u);
 79  78 }
 80  79 //判两线段相交,不包括端点和部分重合
 81  80 int intersect_ex(line u,line v)
 82  81 {
 83  82     return opposite_side(u.a,u.b,v)&&opposite_side(v.a,v.b,u);
 84  83 }
 85  84 //计算两直线交点,注意事先判断直线是否平行!
 86  85 //线段交点请另外判线段相交(同时还是要判断是否平行!)
 87  86 point intersection(line u,line v)
 88  87 {
 89  88     point ret=u.a;
 90  89     double t=((u.a.x-v.a.x)*(v.a.y-v.b.y)-(u.a.y-v.a.y)*(v.a.x-v.b.x))/((u.a.x-u.b.x)*(v.a.y-v.b.y)-(u.a.y-u.b.y)*(v.a.x-v.b.x));
 91  90     ret.x+=(u.b.x-u.a.x)*t;
 92  91     ret.y+=(u.b.y-u.a.y)*t;
 93  92     return ret;
 94  93 }
 95  94 //点到直线上的最近点
 96  95 point ptoline(point p,line l)
 97  96 {
 98  97     point t=p;
 99  98     t.x+=l.a.y-l.b.y,t.y+=l.b.x-l.a.x;
100  99     return intersection(p,t,l.a,l.b);
101 100 }
102 101 //点到直线距离
103 102 double disptoline(point p,line l)
104 103 {
105 104     return fabs(xmult(p,l.a,l.b))/distance(l.a,l.b);
106 105 }
107 106 //点到线段上的最近点
108 107 point ptoseg(point p,line l)
109 108 {
110 109     point t=p;
111 110     t.x+=l.a.y-l.b.y,t.y+=l.b.x-l.a.x;
112 111     if(xmult(l.a,t,p)*xmult(l.b,t,p)>eps)
113 112         return distance(p,l.a)<distance(p,l.b)?l.a:l.b;
114 113     return intersection(p,t,l.a,l.b);
115 114 }
116 115 //点到线段距离
117 116 double disptoseg(point p,line l)
118 117 {
119 118     point t=p;
120 119     t.x+=l.a.y-l.b.y,t.y+=l.b.x-l.a.x;
121 120     if(xmult(l.a,t,p)*xmult(l.b,t,p)>eps)
122 121         return distance(p,l.a)<distance(p,l.b)?distance(p,l.a):distance(p,l.b);
123 122     return fabs(xmult(p,l.a,l.b))/distance(l.a,l.b);
124 123 }
125 124 struct TPoint
126 125 {
127 126     double x,y;
128 127     TPoint operator-(TPoint&a)
129 128     {
130 129         TPoint p1;
131 130         p1.x=x-a.x;
132 131         p1.y=y-a.y;
133 132         return p1;
134 133     }
135 134 };
136 135 
137 136 struct TLine
138 137 {
139 138     double a,b,c;
140 139 };
141 140 
142 141 //求p1关于p2的对称点
143 142 TPoint symmetricalPoint(TPoint p1,TPoint p2)
144 143 {
145 144     TPoint p3;
146 145     p3.x=2*p2.x-p1.x;
147 146     p3.y=2*p2.y-p1.y;
148 147     return p3;
149 148 }
150 149 //p点关于直线L的对称点
151 150 TPoint symmetricalPointofLine(TPoint p,TLine L)
152 151 {
153 152     TPoint p2;
154 153     double d;
155 154     d=L.a*L.a+L.b*L.b;
156 155     p2.x=(L.b*L.b*p.x-L.a*L.a*p.x-2*L.a*L.b*p.y-2*L.a*L.c)/d;
157 156     p2.y=(L.a*L.a*p.y-L.b*L.b*p.y-2*L.a*L.b*p.x-2*L.b*L.c)/d;
158 157     return p2;
159 158 }
160 159 //求线段所在直线,返回直线方程的三个系数
161 160 //两点式化为一般式
162 161 TLine lineFromSegment(TPoint p1,TPoint p2)
163 162 {
164 163     TLine tmp;
165 164     tmp.a=p2.y-p1.y;
166 165     tmp.b=p1.x-p2.x;
167 166     tmp.c=p2.x*p1.y-p1.x*p2.y;
168 167     return tmp;
169 168 }
170 169 //求直线的交点
171 170 //求直线的交点,注意平行的情况无解,避免RE
172 171 TPoint LineInter(TLine l1,TLine l2)
173 172 {
174 173     //求两直线得交点坐标
175 174     TPoint tmp;
176 175     double a1=l1.a;
177 176     double b1=l1.b;
178 177     double c1=l1.c;
179 178     double a2=l2.a;
180 179     double b2=l2.b;
181 180     double c2=l2.c;
182 181     //注意这里b1=0
183 182     if(fabs(b1)<eps){
184 183         tmp.x=-c1/a1;
185 184         tmp.y=(-c2-a2*tmp.x)/b2;
186 185     }
187 186     else{
188 187         tmp.x=(c1*b2-b1*c2)/(b1*a2-b2*a1);
189 188         tmp.y=(-c1-a1*tmp.x)/b1;
190 189     }
191 190     //cout<<"交点坐标"<<endl;
192 191     //cout<<a1*tmp.x+b1*tmp.y+c1<<endl;
193 192     //cout<<a2*tmp.x+b2*tmp.y+c2<<endl;
194 193     return tmp;
195 194 }
196 195 //矢量(点)V以P为顶点逆时针旋转angle(弧度)并放大scale倍
197 196 point rotate(point v,point p,double angle,double scale){
198 197     point ret=p;
199 198     v.x-=p.x,v.y-=p.y;
200 199     p.x=scale*cos(angle);
201 200     p.y=scale*sin(angle);
202 201     ret.x+=v.x*p.x-v.y*p.y;
203 202     ret.y+=v.x*p.y+v.y*p.x;
204 203     return ret;
205 204 }
206 205 //矢量(点)V以P为顶点逆时针旋转angle(弧度)
207 206 point rotate(point v,point p,double angle){
208 207     double cs=cos(angle),sn=sin(angle);
209 208     v.x-=p.x,v.y-=p.y;
210 209     p.x+=v.x*cs-v.y*sn;
211 210     p.y+=v.x*sn+v.y*cs;
212 211     return p;
213 212 }
View Code

 

posted on 2014-12-28 17:41  Pentium.Labs  阅读(241)  评论(0编辑  收藏  举报



Pentium.Lab Since 1998