Most Distant Point from the Sea POJ - 3525 

题意: 求点到多边形的最远距离

思路: 用半平面交,求凸多边形的最大内切圆半径,二分r

#include<cstdio>
#include<vector>
#include<cmath>
#include<string>
#include<string.h>
#include<iostream>
#include<algorithm>
#define PI acos(-1.0)
#define pb push_back
#define F first
#define S second
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
const int N=1e3+3;
const int MOD=1e9+7;
template <class T>
bool sf(T &ret){ //Faster Input
    char c; int sgn; T bit=0.1;
    if(c=getchar(),c==EOF) return 0;
    while(c!='-'&&c!='.'&&(c<'0'||c>'9')) c=getchar();
    sgn=(c=='-')?-1:1;
    ret=(c=='-')?0:(c-'0');
    while(c=getchar(),c>='0'&&c<='9') ret=ret*10+(c-'0');
    if(c==' '||c=='\n'){ ret*=sgn; return 1; }
    while(c=getchar(),c>='0'&&c<='9') ret+=(c-'0')*bit,bit/=10;
    ret*=sgn;
    return 1;
}
int sign(double x) {  //三态函数，减少精度问题
    return abs(x) < 1e-10 ? 0 : x < 0 ? -1 : 1;
}
struct Point { //点的定义
    double x, y;
    Point(double x=0.0, double y=0.0) : x(x), y(y) {}
    Point operator + (const Point &rhs) const {  //向量加法
        return Point(x + rhs.x, y + rhs.y);
    }
    Point operator - (const Point &rhs) const {  //向量减法
        return Point(x - rhs.x, y - rhs.y);
    }
    Point operator * (double p) const {  //向量乘以标量
        return Point(x * p, y * p);
    }
    Point operator / (double p) const {  //向量除以标量
        return Point(x / p, y / p);
    }
    bool operator < (const Point &rhs) const {  //点的坐标排序
        return x < rhs.x || (x == rhs.x && y < rhs.y);
    }
    bool operator == (const Point &rhs) const {  //判断同一个点
        return sign(x - rhs.x) == 0 && sign(y - rhs.y) == 0;
    }

};
typedef Point Vector;
double polar_angle(Vector A) {  //向量极角
    return atan2(A.y, A.x);
}

struct Line {  //有向直线的定义，左边就是对应的半平面
    Point p,s;  //直线上任意一点
    Vector v;  //方向向量
    double ang;  //极角，即从x正方向旋转到向量v所需的角（弧度）
    Line() {}
    Line(Point p, Point s, Vector v) : p(p), s(s) ,v(v) {
        ang = polar_angle(v);
    }
    bool operator < (const Line &rhs) const {
        return ang < rhs.ang;
    }
    Point point(double a) {
        return p + v * a;
    }
};
double cross(Vector A, Vector B) {  //向量叉积，两向量组成的三角形的有向面积的两倍，可判断两向量的方向
    return A.x * B.y - A.y * B.x;
}
bool point_on_left(Point p, Line L) {  //判断点在直线左侧
    return cross(L.v, p-L.p) > 0;
}
Point line_line_inter(Point p, Vector V, Point q, Vector W) {  //两直线交点，参数方程
    Vector U = p - q;
    double t = cross(W, U) / cross(V, W);  //两直线有唯一交点的前提：cross(V, W)非0
    return p + V * t;
}
Point p[N];
Line Q[N];
typedef vector<Point> polygon;
typedef vector<Point> polygon;

bool HPIcmp(Line a,Line b){
	if(fabs(a.ang - b.ang) > 1e-7)return a.ang < b.ang; //斜率排序
	//斜率相同
	return sign(cross((a.p - b.p),(b.s - b.p))) < 0;
}
polygon half_plane_inter(vector<Line> line){
    int n=(int)line.size();
    int tot=n;
    sort(line.begin(),line.end(),HPIcmp);
    polygon ans;
    tot = 1;
	for(int i = 1;i < n;i++)
		if(fabs(line[i].ang - line[i-1].ang) > 1e-7) //去掉斜率重复的
			line[tot++] = line[i];
	int head = 0, tail = 1;
	Q[0] = line[0];
	Q[1] = line[1];
	for(int i = 2; i < tot; i++){ //如果双端队列底端或顶端两个半平面的交点在当前半平面之外，则删除
		if(fabs(cross((Q[tail].s-Q[tail].p),(Q[tail-1].s-Q[tail-1].p))) <1e-7 ||
		fabs(cross((Q[head].s-Q[head].p),(Q[head+1].s-Q[head+1].p))) <1e-7)
			return ans;
		while(head < tail && sign(cross(line_line_inter(Q[tail].p,Q[tail].v,Q[tail-1].p,Q[tail-1].v)-line[i].p,line[i].s-line[i].p))>0)
			tail--;
		while(head < tail && sign(cross(line_line_inter(Q[head].p,Q[head].v,Q[head+1].p,Q[head+1].v)-line[i].p,line[i].s-line[i].p))>0)
			head++;
		Q[++tail] = line[i];
	}
	while(head < tail && sign(cross(line_line_inter(Q[tail].p,Q[tail].v,Q[tail-1].p,Q[tail-1].v)-Q[head].p,Q[head].s-Q[head].p))>0)
		tail--;
	if(tail <= head + 1) return ans;
	//保存点到ans
	for(int i = head; i < tail; i++)
		ans.pb(line_line_inter(Q[i].p,Q[i].v,Q[i+1].p,Q[i+1].v));
	if(head < tail - 1)
		ans.pb(line_line_inter(Q[head].p,Q[head].v,Q[tail].p,Q[tail].v));
    return ans;
}
double dot(Vector a,Vector b){
    return a.x*b.x+a.y*b.y;
}
double length(Vector A) {  //向量长度，点积
    return sqrt(dot(A, A));
}
Vector normal(Vector A) {  //向量的单位法向量，确保A不是零向量
    double len = length(A);
    return Vector(-A.y/len, A.x/len);
}
Vector V1[N];
int main(void){
    int T,n;
    while(cin>>n){
        if(!n)  break;
        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++){
            Line t;
            t={p[i],p[(i+1)%n],p[(i+1)%n]-p[i]};
            V1[i]=normal(p[(i+1)%n]-p[i]);
        }
        double l=0,r=20000;
        for(int i=1;i<=50;++i){
            vector<Line> L;
            double mid=(l+r)/2.0;
            for(int i=0;i<n;i++){
                Line q( (p[i]+V1[i]*mid),(p[(i+1)%n]+V1[i]*mid),p[(i+1)%n]-p[i] );
                L.pb(q);
            }
            vector<Point> res=half_plane_inter(L);

            if((int)res.size()>0)   l=mid;
            else    r=mid;
        }
        printf("%.10f\n",r);

    }

    return 0;
}