POJ3608

计算两个凸包之间的最小距离,旋转卡壳法详解在旋转卡壳的用法之计算两个凸 包上的最近距离

#include <iostream>

#include<cstdio>

#include<string.h>

#include<cmath>

using namespace std;

const double eps=1e-10;

const double INF=1e10;

struct point

{

    double x,y;

    point (double a=0,double b=0)

    {

        x=a;

        y=b;

    }

};

double min_val(double a,double b)

{

    return a>b?b:a;

}

int dcmp(double a)

{

    if(fabs(a)<eps)return 0;

    else return a>0?1:-1;

}

point operator -(point a,point b){ return point (a.x-b.x,a.y-b.y); }

bool operator ==(const point &a,const point &b)

{

    return dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)==0;

}

double cross(point a,point b)

{

    return a.x*b.y-a.y*b.x;

}

void readdate(point a[],int n,point b[],int m)

{

    int i;

    for(i=0;i<n;i++)

    scanf("%lf%lf",&a[i].x,&a[i].y);

    for(i=0;i<m;i++)

    scanf("%lf%lf",&b[i].x,&b[i].y);

}

void change(point a[],int b)

{

    int i;

    point temp;

    for(i=0;i<b/2;i++)

    {

        temp=a[i];

        a[i]=a[b-1-i];

        a[b-1-i]=temp;

    }

}

void work(point p[],int n,point q[],int m)

{

    int i;

    for(i=0;i<n-2;i++)

    if(dcmp(cross(p[i+1]-p[i],p[i+2]-p[i])>0)) break;

    else if(dcmp(cross(p[i+1]-p[i],p[i+2]-p[i])<0))

    {

        change(p,n);

        break;

    }

    for(i=0;i<m-2;i++)

    if(dcmp(cross(q[i+1]-q[i],q[i+2]-q[i])>0))break;

    else if(dcmp(cross(q[i+1]-q[i],q[i+2]-q[i])<0))

    {

        change(q,m);

        break;

    }

}

double Dot(point A,point B){return A.x*B.x+A.y*B.y;}

double Length(point A){return sqrt(Dot(A,A));}

double dist(point a,point b,point p)

{

    if(a==b) return Length(p-a);

    point 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);

}

double Mid(point a1,point a2,point b1,point b2)

{

    return min_val(min_val(dist(a1,a2,b1),dist(a1,a2,b2)),

                   min_val(dist(b1,b2,a1),dist(b1,b2,a2)));

}

double solve(point P[],int n,point Q[],int m)

{

    double ans=INF,temp;

    int i,minP=0,maxQ=0;

    for(i=0;i<n;i++)

    if(P[i].y<P[minP].y)

    minP=i;

    for(i=0;i<m;i++)

    if(Q[i].y>Q[maxQ].y)

    maxQ=i;

    Q[m]=Q[0];

    P[n]=P[0];

    for(i=0;i<n;i++)

    {

        while((temp=cross(Q[maxQ+1]-P[minP+1],P[minP]-P[minP+1])

               -cross(Q[maxQ]-P[minP+1],P[minP]-P[minP+1]))>eps)

        maxQ=(maxQ+1)%m;

        if(temp<-eps)

            ans=min_val(ans,dist(P[minP],P[minP+1],Q[maxQ]));

        else ans=min_val(ans,Mid(P[minP],P[minP+1],Q[maxQ],Q[maxQ+1]));

        minP=(minP+1)%n;

    }

    return ans;

}

int main()

{

    point p[10005],q[10005];

    int n,m;

    double a,b;

    while(scanf("%d%d",&n,&m)==2)

    {

        if(m==0&&n==0)break;

        readdate(p,n,q,m);

        work(p,n,q,m);

       a=solve(p,n,q,m);

       b=solve(q,m,p,n);

       printf("%.5lf\n",min_val(a,b));

    }

    return 0;

}