POJ 3384 Feng Shui （半平面交）

题目链接：http://poj.org/problem?id=3384

题目：给出一个凸多边形的房间，根据风水要求，把两个圆形地毯铺在房间里，不能折叠，不能切割，可以重叠。问最多能覆盖多大空间，输出两个地毯的圆心坐标。多组解输出其中一个

思路：

题目保证至少可以放入一个圆，上一题中判断过在一个多边形内是否能放入一个半径为r的圆。

同样将多边形的边内移R之后，半平面交区域便是可以放入圆的可行区域。

题目要求覆盖的面积最大，也就是两个圆的半径相同，圆心越远，面积就越大。

在半平面交区域内找到最远点对便是题目要求的解。

#include <stdio.h>
#include <string.h>
#include <iostream>
#include <algorithm>
#include <vector>
#include <queue>
#include <set>
#include <map>
#include <string>
#include <math.h>
#include <stdlib.h>
#include <time.h>
using namespace std;
const double eps = 1e-8;
const double PI = acos(-1.0);
int sgn(double x)
{
    if(fabs(x) < eps) return 0;
    if(x < 0) return -1;
    else return 1;
}
struct Point
{
    double x,y;
    Point(){}
    Point(double _x,double _y)
    {
        x = _x; y = _y;
    }
    Point operator -(const Point &b)const
    {
        return Point(x - b.x, y - b.y);
    }
    double operator ^(const Point &b)const
    {
        return x*b.y - y*b.x;
    }
    double operator *(const Point &b)const
    {
        return x*b.x + y*b.y;
    }
};
struct Line
{
    Point s,e;
    double k;
    Line(){}
    Line(Point _s,Point _e)
    {
        s = _s; e = _e;
        k = atan2(e.y - s.y,e.x - s.x);
    }
    Point operator &(const Line &b)const
    {
        Point res = s;
        double t = ((s - b.s)^(b.s - b.e))/((s - e)^(b.s - b.e));
        res.x += (e.x - s.x)*t;
        res.y += (e.y - s.y)*t;
        return res;
    }
};
//半平面交，直线的左边代表有效区域
bool HPIcmp(Line a,Line b)
{
    if(fabs(a.k - b.k) > eps)return a.k < b.k;
    return ((a.s - b.s)^(b.e - b.s)) < 0;
}
Line Q[1010];
void HPI(Line line[], int n, Point res[], int &resn)
{
    int tot = n;
    sort(line,line+n,HPIcmp);
    tot = 1;
    for(int i = 1;i < n;i++)
        if(fabs(line[i].k - line[i-1].k) > eps)
            line[tot++] = line[i];
    int head = 0, tail = 1;
    Q[0] = line[0];
    Q[1] = line[1];
    resn = 0;
    for(int i = 2; i < tot; i++)
    {
        if(fabs((Q[tail].e-Q[tail].s)^(Q[tail-1].e-Q[tail-1].s)) < eps || fabs((Q[head].e-Q[head].s)^(Q[head+1].e-Q[head+1].s)) < eps)
            return;
        while(head < tail && (((Q[tail]&Q[tail-1]) - line[i].s)^(line[i].e-line[i].s)) > eps)
            tail--;
        while(head < tail && (((Q[head]&Q[head+1]) - line[i].s)^(line[i].e-line[i].s)) > eps)
            head++;
        Q[++tail] = line[i];
    }
    while(head < tail && (((Q[tail]&Q[tail-1]) - Q[head].s)^(Q[head].e-Q[head].s)) > eps)
        tail--;
    while(head < tail && (((Q[head]&Q[head-1]) - Q[tail].s)^(Q[tail].e-Q[tail].e)) > eps)
        head++;
    if(tail <= head + 1)return;
    for(int i = head; i < tail; i++)
        res[resn++] = Q[i]&Q[i+1];
    if(head < tail - 1)
        res[resn++] = Q[head]&Q[tail];
}
Point p[1010];
Line line[1010];
//*两点间距离
double dist(Point a,Point b)
{
    return sqrt((a-b)*(a-b));
}
void change(Point a,Point b,Point &c,Point &d,double p)//将线段ab往左移动距离p
{
    double len = dist(a,b);
    double dx = (a.y - b.y)*p/len;
    double dy = (b.x - a.x)*p/len;
    c.x = a.x + dx; c.y = a.y + dy;
    d.x = b.x + dx; d.y = b.y + dy;
}
int main()
{
    //freopen("in.txt","r",stdin);
    //freopen("out.txt","w",stdout);
    int n;
    double r;
    while(scanf("%d%lf",&n,&r) == 2)
    {
        for(int i = 0;i < n;i++)
            scanf("%lf%lf",&p[i].x,&p[i].y);
        reverse(p,p+n);
        for(int i = 0;i < n;i++)
        {
            Point t1,t2;
            change(p[i],p[(i+1)%n],t1,t2,r);
            line[i] = Line(t1,t2);
        }
        int resn;
        HPI(line,n,p,resn);
        //resn为半平面的交点个数
        int res1 = 0, res2 = 0;
        for(int i = 0;i < resn;i++)
            for(int j = i;j < resn;j++)
                if(dist(p[i],p[j]) > dist(p[res1],p[res2]))
                    res1 = i, res2 = j;
        printf("%.5f %.5f %.5f %.5f\n",p[res1].x,p[res1].y,p[res2].x,p[res2].y);
    }
    return 0;
}