POJ 2079 旋转卡壳

最新推荐文章于 2019-09-30 13:48:14 发布
chm517
最新推荐文章于 2019-09-30 13:48:14 发布
阅读量380
点赞数
分类专栏：【计算几何】
本文链接：https://blog.csdn.net/chm517/article/details/43930065
版权
【计算几何】专栏收录该内容
10 篇文章 0 订阅
订阅专栏
//check if all is d_ ,convenient to double &&long double
#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<vector>
#include<set>
#include<algorithm>
#define sqr(x) (x)*(x)
#define f1 first
#define f2 second
#define pb push_back
#define fr(i,x,y) for(int i=x;i<=y;++i)
#define sqr(x) (x)*(x)
typedef double d_;
const d_ pi=acos(-1.0);
using namespace std;


struct Point{
    d_ x,y;
    Point(d_ x=0,d_ y=0):x(x),y(y){}
    void read(){scanf("%lf%lf",&x,&y);}//lf or LF
    void print(){printf("~~~~~~%lf %lf\n",x,y);}//lf or LF
};
struct Point3{
    d_ x,y,z;
    Point3(d_ x=0,d_ y=0,d_ z=0):x(x),y(y),z(z){}
    void read(){scanf("%lf%lf%lf",&x,&y,&z);}//lf or LF
    void print(){printf("~~~~~~%lf %lf %lf\n",x,y,z);}//lf or LF
};
typedef Point Vector;
typedef Point3 Vector3;
const d_ eps=1e-8;
int dcmp(d_ x){
    if (fabs(x)<eps)  return 0;else return x<0?-1:1;}
inline double ensure_asin(double t){if (fabs(t+1)<eps||fabs(t-1)<eps)t=floor(t+0.5);return t;}
inline double ensure_acos(double t){if (fabs(t+1)<eps||fabs(t-1)<eps)t=floor(t+0.5);return t;}
inline void ensure_atan2(Vector &v){if (fabs(v.y)<eps)v.y=eps;}


inline Vector operator+(Vector a,Vector b) {return Vector(a.x+b.x,a.y+b.y);}
inline Vector operator-(Vector a,Vector b) {return Vector(a.x-b.x,a.y-b.y);}
inline Vector operator*(Vector a,d_ p) {return Vector(a.x*p  ,a.y*p);  }
inline Vector operator*(d_ p,Vector a) {return Vector(a.x*p  ,a.y*p);  }
inline Vector operator/(Vector a,d_ p) {return Vector(a.x/p  ,a.y/p);  }
Vector3 operator+(Vector3 a,Vector3 b) {return Vector3(a.x+b.x,a.y+b.y,a.z+b.z);}
Vector3 operator-(Vector3 a,Vector3 b) {return Vector3(a.x-b.x,a.y-b.y,a.z-b.z);}
Vector3 operator*(Vector3 a,d_ p) {return Vector3(a.x*p  ,a.y*p,a.z*p);  }
Vector3 operator*(d_ p,Vector3 a) {return Vector3(a.x*p  ,a.y*p,a.z*p);  }
Vector3 operator/(Vector3 a,d_ p) {return Vector3(a.x/p  ,a.y/p,a.z/p);  }
d_ cross(Vector a,Vector b){return a.x*b.y-a.y*b.x;}
d_ operator*(const Vector &a,const Vector&b){return a.x*b.y-a.y*b.x;}
d_ cross(Point a,Point b,Point c){return (b-a)*(c-a);}
d_ dot(Vector a,Vector b){return a.x*b.x+a.y*b.y;}
d_ operator%(const Vector &a,const Vector&b){return a.x*b.x+a.y*b.y;}
Vector3 cross(Vector3 a,Vector3 b){
    return Vector3(a.y*b.z-a.z*b.y,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x);
}
Vector3 operator*(const Vector3 &a,const Vector3&b){
    return Vector3(a.y*b.z-a.z*b.y,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x);
}
d_ dot(Vector3 a,Vector3 b){return a.x*b.x+a.y*b.y+a.z*b.z;}
d_ operator%(const Vector3 &a,const Vector3&b){return a.x*b.x+a.y*b.y+a.z*b.z;}
bool operator<(const Point &a,const Point&b)
{
    return dcmp(a.x-b.x)<0||(dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)<0);
}
bool operator==(const Point &a,const Point&b) {return dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)==0;}
bool operator<(const Point3 &a,const Point3 &b)
{
    return dcmp(a.x-b.x)<0||(dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)<0)||(dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)==0&&dcmp(a.z-b.z)<0);
}
bool operator==(const Point3 &a,const Point3&b) {return dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)==0&&dcmp(a.z-b.z)==0;}
//---------------------------------


inline d_ length(Vector a) {return sqrt(dot(a,a));}
d_ angle(Vector a,Vector b){return acos(ensure_acos(dot(a,b)/length(a)/length(b)));}
d_ area2(Point a,Point b,Point c){return cross(b-a,c-a);}
d_ dist(Point &a,Point &b){return length(b-a);}
Vector rotate(Vector a,d_ rad)
{   return Vector(a.x*cos(rad)-a.y*sin(rad),a.x*sin(rad)+a.y*cos(rad));
}
Vector normal(Vector a)//ensure a!=(0,0)
{   d_ l=length(a);
    return Vector(-a.y/l,a.x/l);
}
d_ length(Vector3 a) {return sqrt(dot(a,a));}
d_ dist(Point3 &a,Point3 &b){return length(b-a);}
//--------------------------------------------------------------
struct Line{
    Point p;
    Vector v;
    d_ ang;
    Line(){}
    Line(Point p,Vector v):p(p),v(v){ensure_atan2(v);ang=atan2(v.y,v.x);}
    bool operator<(const Line &l)const{return ang<l.ang;}
    Point getPoint(d_ t)
    {
        return p+v*t;
    }
    d_ getabc(d_ A,d_ B,d_ C)
    {
        A=-v.y;
        B=v.x;
        C=p.x*(v.y+p.y)-p.y*(p.x+v.x);

    }
    //两直线相交求交点
    //第一个值为0表示直线重合，为1表示平行,为2是相交
    //只有第一个值为2时，交点才有意义
    pair<int,Point> operator &(const Line &b)const
    {
        Point res = p;
        if(dcmp(v*b.v) == 0)
        {
            if(dcmp((p-b.p)*b.v) == 0)
                return make_pair(0,res);//重合
            else return make_pair(1,res);//平行
        }
        d_ t = (b.v*(p-b.p))/(v*b.v);
        res=res+t*v;
        return make_pair(2,res);
    }
};

struct Circle{
    Point c;
    d_ r;
    Circle(){}
    Circle(Point c,d_ r):c(c),r(r){}
    Point getPoint(d_ a)
    {
        return Point(c.x+cos(a)*r,c.y+sin(a)*r);
    }
    void read()
    {
        scanf("%lf%lf%lf",&c.x,&c.y,&r);//lf or LF
    }
};
Point get_Line_intersection(Point p,Vector v,Point q,Vector w){
    Vector u=p-q;
    d_ t=cross(w,u)/cross(v,w);
    return p+v*t;
}
d_ distance_to_Line(Point p,Point a,Point b)//p to Line(ab)
{
    Vector v1=b-a,v2=p-a;
    return fabs(cross(v1,v2))/length(v1);
}
inline d_ distance_to_segment(Point p,Point a,Point b)
{
    if (a==b) return length(p-a);
    Vector v1=b-a,v2=p-a,v3=p-b;
    if (dcmp(dot(v1,v2))<0)  return length(v2);
    if (dcmp(dot(v1,v3))>0)  return length(v3);
    return fabs(cross(v1,v2))/length(v1);
}
d_ distance_seg_to_seg(Point A,Point B,Point C,Point D)
{
    return min(min(distance_to_segment(C,A,B),distance_to_segment(D,A,B)),min(distance_to_segment(A,C,D),distance_to_segment(B,C,D)));
}
Point get_Line_projection(Point p,Point a,Point b)
{
    Vector v=b-a;
    return a+v*(dot(v,p-a)/dot(v,v));
}
bool segment_proper_intersection(Point a1,Point a2,Point b1,Point b2)//ensure line(a1,a2)!=line(b1,b2)
{
    d_ c1=cross(a2-a1,b1-a1),c2=cross(a2-a1,b2-a1);
    d_ c3=cross(b2-b1,a1-b1),c4=cross(b2-b1,a2-b1);
    return dcmp(c1)*dcmp(c2)<0&&dcmp(c3)*dcmp(c4)<0;
}

bool Point_on_line(Point p,Point a,Point b){
    return dcmp(cross(a-p,b-p))==0;
}
bool on_segment(Point p,Point a,Point b)
{
    return dcmp(cross(a-p,b-p))==0&&dcmp(dot(a-p,b-p))<0;
}
bool notstrict_on_segment(Point p,Point a,Point b)
{
    return dcmp(cross(a-p,b-p))==0&&dcmp(dot(a-p,b-p))<=0;
}
//--------------------------------
int get_Line_Circle_intersection(Line l,Circle c,d_ &t1,d_ &t2,vector<Point>&sol)
{
    d_ a=l.v.x,b=l.p.x-c.c.x,cc=l.v.y,d=l.p.y-c.c.y;
    d_ e=a*a+cc*cc,f=2*(a*b+cc*d), g=b*b+d*d-c.r*c.r;
    d_ delta=f*f-4*e*g;
    if (dcmp(delta)<0) return 0;
    if (dcmp(delta)==0) {
        t1=t2=-f/(2*e); sol.pb(l.getPoint(t1));
        return 1;
    }
    t1=(-f-sqrt(delta))/(2*e); sol.pb(l.getPoint(t1));
    t2=(-f+sqrt(delta))/(2*e); sol.pb(l.getPoint(t2));
    return 2;
}
d_ angle(Vector v)
{   return atan2(v.y,v.x);
}

int get_Circle_Cirlce_intersection(Circle c1,Circle c2,vector<Point>&sol)
{
    d_ d=length(c1.c-c2.c);
    if (dcmp(d)==0)
        {
            if (dcmp(c1.r-c2.r)==0) return -1;
            return 0;
        }
    if (dcmp(c1.r+c2.r-d)<0) return 0;
    if (dcmp(fabs(c1.r-c2.r)-d)>0) return 0;

    d_ a=angle(c2.c-c1.c);
    d_ da=acos((c1.r*c1.r+d*d-c2.r*c2.r)/(2*c1.r*d));
    Point p1=c1.getPoint(a-da),p2=c1.getPoint(a+da);
    sol.pb(p1);
    if (p1==p2) return 1;
    sol.pb(p2);
    return 2;
}

int get_tangents1(Point p ,Circle c,vector<Vector>&sol)
{
    Vector u=c.c-p;
    d_ dist=length(u);
         if (dist<c.r) return 0;
    else if (dcmp(dist-c.r)==0)
        {
            sol.pb(rotate(u,pi/2));
            return 1;
        }
    else{   d_ ang=asin(c.r/dist);
            sol.pb(rotate(u,-ang));
            sol.pb(rotate(u,ang));
            return 2;
        }
}

int get_tangents2(Circle ca,Circle cb,Point *a,Point *b)
{
    int cnt=0;
    if (ca.r<cb.r) {swap(ca,cb);
                    swap(a,b);
                   }
    int d2=sqr(ca.c.x-cb.c.x)+sqr(ca.c.y-cb.c.y);
    int rdiff=ca.r-cb.r;
    int rsum=ca.r+cb.r;
    if (d2<rdiff*rdiff) return 0;

    d_ base=atan2(cb.c.y-ca.c.y,cb.c.x-ca.c.x);
    if (d2==0&&ca.r==cb.r) return -1;
    if (d2==rdiff*rdiff)
        {   a[cnt]=ca.getPoint(base); b[cnt]=cb.getPoint(base); cnt++;
            return 1;
        }

    d_ ang=acos((ca.r-cb.r)/sqrt(double(d2)));
    a[cnt]=ca.getPoint(base+ang); b[cnt]=cb.getPoint(base+ang);  cnt++;
    a[cnt]=ca.getPoint(base-ang); b[cnt]=cb.getPoint(base-ang);  cnt++;
    if (d2==rsum*rsum)
        {
            a[cnt]=ca.getPoint(base);
            b[cnt]=cb.getPoint(pi+base);
            cnt++;
        }
    else if (d2>rsum*rsum)
        {
            d_ ang=acos((ca.r+cb.r)/sqrt(double(d2)));
            a[cnt]=ca.getPoint(base+ang); b[cnt]=cb.getPoint(base+pi+ang);  cnt++;
            a[cnt]=ca.getPoint(base-ang); b[cnt]=cb.getPoint(base+pi-ang);  cnt++;
        }
    return cnt;
}


bool segment_intersection(Point a,Point b,Point c,Point d)//have and only one intersection
{
    Line l1(a,b-a),l2(c,d-c);
    if ((l1&l2).first==0)
    {
        if (a==c&&!on_segment(b,c,d)&&!on_segment(d,a,b))return 1;
        if (a==d&&!on_segment(b,c,d)&&!on_segment(c,a,b))return 1;
        if (b==c&&!on_segment(a,c,d)&&!on_segment(d,a,b))return 1;
        if (b==d&&!on_segment(a,c,d)&&!on_segment(c,a,b))return 1;
        return 0;
    }
   if (segment_proper_intersection(a,b,c,d)||on_segment(a,c,d)||on_segment(b,c,d)||on_segment(c,a,b)||on_segment(d,a,b)||(a==c)||(a==d)||(b==c)||(b==d))
   return 1;else return 0;
}
bool segment_notproper_intersection(Point a,Point b,Point c,Point d)//have at least one intersection
{
    Line l1(a,b-a),l2(c,d-c);
    if ((l1&l2).first==0)
    {
        if (a==c||a==d||b==c||b==d) return 1;
        if (on_segment(a,c,d)||on_segment(b,c,d)||on_segment(c,a,b)||on_segment(d,a,b)==1)return 1;
        return 0;
    }
   if (segment_proper_intersection(a,b,c,d)||on_segment(a,c,d)||on_segment(b,c,d)||on_segment(c,a,b)||on_segment(d,a,b)||(a==c)||(a==d)||(b==c)||(b==d))
   return 1;else return 0;
}
int n;

d_ polygon_area(Point *p,int n)
{
    d_ re=0;
    fr(i,1,n-2)
        re+=(p[i]-p[0])*(p[i+1]-p[0]);
    return re;
}

Point bcenter1(Point *p,int n)
{
    d_ area=0,tt;
    Point t,s,ans(0,0);
    t=p[0];
    fr(i,1,n)
        {
            s=p[(i==n)?0:i];
            tt=t*s;area+=tt;
            ans=ans+(t+s)*tt;
            t=s;
        }
    ans=ans/(area*3);
    return ans;
}
Point bcenter2(Point *p,int n)
{
    d_ area=0,tt;
    Point t,s,ans(0,0);
    fr(i,1,n-2)
        {
            tt=(p[i]-p[0])*(p[i+1]-p[0]);
            ans=ans+(p[i]+p[i+1]+p[0])*tt;
            area+=tt;
        }
    ans=ans/(area*3);
    return ans;
}
int is_point_in(Point p,Point poly[],int n)
{
    int w=0;
    fr (i,0,n-1)
        {
            if (notstrict_on_segment(p,poly[i],poly[(i+1)%n])) return -1;
            int k=dcmp((poly[(i+1)%n]-poly[i])*(p-poly[i]));
            int d1=dcmp(poly[i].y-p.y);
            int d2=dcmp(poly[(i+1)%n].y-p.y);
            if (k>0&&d1<=0&&d2>0)w++;
            if (k<0&&d2<=0&&d1>0)w--;
        }
    return (w!=0);
}
//判断直线ab是否与线段cd相交
//0 不相交
//1 规范相交
//2 不规范相交
int seg(const Point& a,const Point& b,const Point& c,const Point& d)
{
    int d1,d2;
    d1 = dcmp((b-a)*(c-a));
    d2 = dcmp((b-a)*(d-a));
    if(d1*d2<0)
        return 1;
    if(d1==0||d2==0)
        return 2;
    return 0;
}


int convex(Point* st,int n,Point *ps)
{
    int k=0;
    sort(st,st+n);
    for (int i=0;i<n;i++)
        {   while (k>1&&(ps[k-1]-ps[k-2])*(st[i]-ps[k-2])<=0) k--;
            //若想含边上的点，可以将<=改成<
            ps[k++]=st[i];
            //printf("~~i=%d k=%d\n",i,k);
        }
    for (int i=n-2,t=k;i>=0;i--)
        {   while (k>t&&(ps[k-1]-ps[k-2])*(st[i]-ps[k-2])<=0) k--;
            ps[k++]=st[i];
            //printf("!!i=%d k=%d\n",i,k);
        }
    if (n>1)k--;
    return k;
}
bool onleft(Line l,Point p)
{
    return dcmp(l.v*(p-l.p))>0;
}
int half_plane_intersection(Line *l,int n,Point*poly)
{
    sort(l,l+n);
    int fir,last;
    Point *p=new Point[n];
    Line *q=new Line[n];
    q[fir=last=0]=l[0];
    fr(i,1,n-1)
        {
            while (fir<last&&!onleft(l[i],p[last-1])) last--;
            while (fir<last&&!onleft(l[i],p[fir])) fir++;
            q[++last]=l[i];
            if (fabs(q[last].v*q[last-1].v)<eps)
                {     if (dcmp(q[last].v%q[last-1].v)<=0){last--;continue;}
                      last--;
                      if(onleft(q[last],l[i].p) )q[last]=l[i];
                }
            if (fir<last)
                p[last-1]=get_Line_intersection(q[last-1].p,q[last-1].v,q[last].p,q[last].v);
        }
    while (fir<last&&!onleft(q[fir],p[last-1]))last--;
    if (last-fir<=1){delete p;delete q;return 0;}
    p[last]=get_Line_intersection(q[last].p,q[last].v,q[fir].p,q[fir].v);
    if (!onleft(q[last-1],p[last])){delete p;delete q;return -1;}
    int m=0;
    fr(i,fir,last) poly[m++]=p[i];
    delete p;
    delete q;
    return m;
}
struct Poly {
    vector<Point> pt;
    Poly() { pt.clear();}
    Poly(vector<Point> &pt) : pt(pt) {}
    Point operator [] (int x) const { return pt[x];}
    int size() { return pt.size();}
    double area()
        {   int n=size();
            d_ re=0;
            fr(i,1,n-2)
                re+=(pt[i]-pt[0])*(pt[i+1]-pt[0]);
            return fabs(re/2.0);
        }
 };
int ptInPoly(Point p, Poly &poly) {
 int wn = 0, sz = poly.size();
 for (int i = 0; i < sz; i++) {
 if (on_segment(p, poly[i], poly[(i + 1) % sz])) return -1;
 int k = dcmp(cross(poly[(i + 1) % sz] - poly[i], p - poly[i]));
 int d1 = dcmp(poly[i].y - p.y);
 int d2 = dcmp(poly[(i + 1) % sz].y - p.y);
 if (k > 0 && d1 <= 0 && d2 > 0) wn++;
 if (k < 0 && d2 <= 0 && d1 > 0) wn--;
 }
 if (wn != 0) return 1;
 return 0;
 }
bool isIntersect(Point a, Point b, Poly &poly)
{   int sz=poly.size();
    fr (i,0,sz-1)
    if (segment_intersection(a, b, poly[i],poly[(i+1)%n])) return 1;
    return 0;
 }
Poly cutPoly(Poly &poly, Point a, Point b)
{
    int n=poly.size();
    Poly npoly=Poly();
    fr(i,0,n-1)
        {
            Point c=poly[i],d=poly[(i+1)%n];
            if (dcmp((b-a)*(c-a))>=0) npoly.pt.pb(c);
            if (dcmp((b-a)*(c-d))!=0)
                {
                    Point ip=get_Line_intersection(a,b-a,c,d-c);
                    if (on_segment(ip,c,d))npoly.pt.pb(ip);
                }
        }
    return npoly;
}


 vector<Poly> cutPolies(Point s, Point t, vector<Poly> polies) {
    vector<Poly> ret;
    ret.clear();
    for (int i = 0, sz = polies.size(); i < sz; i++)
    {
    Poly tmp;
    tmp = cutPoly(polies[i], s, t);
    if (tmp.size() >= 3 && tmp.area() > eps) ret.push_back(tmp);
    tmp = cutPoly(polies[i], t, s);
    if (tmp.size() >= 3 && tmp.area() > eps) ret.push_back(tmp);
    }
    return ret;
 }





double rotating_calipers(Point *p1,Point *p2,int top1,int top2)
{
    int s1=0,s2=0;
    for(int i=1;i<top1;i++)
        if(dcmp(p1[i].y-p1[s1].y)<0||(dcmp(p1[i].y-p1[s1].y)==0)&&dcmp(p1[i].x-p1[s1].x)<0) s1=i;
    for(int i=1;i<top2;i++)
        if(dcmp(p2[i].y-p2[s2].y)>0||(dcmp(p2[i].y-p2[s2].y)==0)&&dcmp(p2[i].x-p2[s2].x)>0) s2=i;
    int t1=s1,t2=s2;
    double ans=length(p1[s1]-p2[s2]);
    int nn=0;
    do
    {
        double af=(p1[(s1+1)%top1]-p1[s1])*(p2[(s2+1)%top2]-p2[s2]);
        if(dcmp(af)==0)//卡壳到两条边
        {
            ans=min(ans,distance_to_segment(p1[s1],p2[s2],p2[(s2+1)%top2]));
            ans=min(ans,distance_to_segment(p1[(s1+1)%top1],p2[s2],p2[(s2+1)%top2]));
            ans=min(ans,distance_to_segment(p2[s2],p1[s1],p1[(s1+1)%top1]));
            ans=min(ans,distance_to_segment(p2[(s2+1)%top2],p1[s1],p1[(s1+1)%top1]));
            s1=(s1+1)%top1; s2=(s2+1)%top2;
        }
        else if(dcmp(af)<0)//卡壳到第一个的边，第二个的点
        {
            ans=min(ans,distance_to_segment(p2[s2],p1[s1],p1[(s1+1)%top1]));
            s1=(s1+1)%top1;
        }
        else//卡壳到第二个的边，第一个的点
        {
            ans=min(ans,distance_to_segment(p1[s1],p2[s2],p2[(s2+1)%top2]));
            s2=(s2+1)%top2;
        }
    }while(t1!=s1||t2!=s2);
    return ans;
}


double ro(int n,Point *p)
{
    int i,j,k;
    double ar=0;
    for (i=0;i<n;i++)
        {
            j=(i+1)%n;k=(j+1)%n;
            double c1,c2;
            while (k!=i)
            {
                double c1=cross(p[i],p[j],p[k]),c2=cross(p[i],p[j],p[(k+1)%n]);
                if (dcmp(c1-c2)<0) k=(k+1)%n;else break;
            }
            if (k==i)continue;
            int kk =(k+1)%n;
            while (j!=kk&&k!=i)
                {
                   double c1=cross(p[i],p[j],p[k]);
                    ar=max(ar,c1);
                    while (k!=i)
                    {
                        double c1=cross(p[i],p[j],p[k]),c2=cross(p[i],p[j],p[(k+1)%n]);
                        if (dcmp(c1-c2)<0) k=(k+1)%n;else break;
                    }
                    j=(j+1)%n;
                }
        }
        return ar/2;
}
#define N 50005
Point aa[N],a[N];
void doit()
{
    fr(i,0,n-1)aa[i].read();
    n=convex(aa,n,a);
    printf("%.2lf\n",ro(n,a));
}
int main()
{
    int cas;
    while (scanf("%d",&n),n!=-1) doit();
}