多边形可以容纳的半径最大的圆

/*
http://acm.hdu.edu.cn/showproblem.php?pid=3644
枚举，如果一个圆能放在多边形里面，你将这个圆移动，
直到跟多边形“相切",这个相切也不能算是标准的相切定义，
反正有三种情况：一、圆碰到了多边形的两个点；二、圆碰到
了多边形的一边一点；三、圆碰到了多边形的两边。以至于圆
不能在移动了，如此，你可以枚举这些情况，再判断圆是否放
得下。
*/

const double eps = 1e-8;
const double max = 1000000000.0; 
int sign(double d){ return d < -eps ? -1 : d > eps; }
struct point{
	double x, y;
	point(double x=0, double y=0) : x(x), y(y) {}
	void read(){ scanf("%lf%lf", &x, &y);}
	point operator-(point tp){ return point(x-tp.x, y-tp.y);  }
	point operator+(point tp){ return point(x+tp.x, y+tp.y);  }
};
inline double inter_pro(point st1, point ed1, point st2, point ed2){
	return (ed1.x-st1.x)*(ed2.y-st2.y) - (ed1.y-st1.y)*(ed2.x-st2.x);
}
inline double dist(point st, point ed){ return sqrt(pow(st.x-ed.x, 2.0) + pow(st.y-ed.y, 2.0)); }
inline bool isIntersect(point st1, point ed1, point st2, point ed2){  //判断线段是否相交,在端点处相交不算 ??
	if(inter_pro(st1, st2, st1, ed1) * inter_pro(st1, ed2, st1, ed1) >= 0) return false;
	if(inter_pro(st2, st1, st2, ed2) * inter_pro(st2, ed1, st2, ed2) >= 0) return false;
	return true;
}
inline bool isOnSegment(point st, point ed, point tp){ //判断点是否在线段上
	if(inter_pro(st, tp, st, ed) != 0) return false;
	if(fabs(dist(st, tp) + dist(ed, tp) - dist(st, ed)) < eps) return true;
	return false;
}
inline double distPtoL(point st, point ed, point tp){//点到线段的最近距离
	double d1, d2, d3;
	d1 = dist(st, ed); d2 = dist(st, tp); d3 = dist(ed, tp);
	if((d1*d1 + d3*d3 - d2*d2 < 0) || (d1*d1+d2*d2-d3*d3 < 0)) return min(d2, d3);
	return fabs(inter_pro(st, ed, st, tp)) / d1;
}
inline point getMid(point st, point ed){  //求线段的中点
	return point((st.x + ed.x) / 2.0, (st.y + ed.y) / 2.0);
}
inline point goLen(point st, point dir, double len){  //st沿着dir方向走len之后的点
	double tl = sqrt(dir.x * dir.x + dir.y * dir.y);
	if(sign(tl) == 0) return st;
	return point(st.x + len * dir.x / tl, st.y + len * dir.y / tl);
}

// (st->ed)的左侧为多边形的内侧
point change(point st, point ed, point next, double l){  //l为线段移动的长度
        //next为和返回的点相对应的点
	point dd;
	dd.x = -(ed - st).y;
	dd.y = (ed - st).x;
	double len = sqrt(dd.x * dd.x + dd.y * dd.y);
	dd.x /= len, dd.y /= len;
	dd.x *= l, dd.y *= l;
	dd = dd + next;
	return dd;
}
point intersectPoint(point st1, point ed1, point st2, point ed2){  //求相交直线的交点
	double t = inter_pro(st2, st1, st2, ed2) / ((ed1.y-st1.y) * (ed2.x-st2.x) - (ed1.x-st1.x) * (ed2.y-st2.y));
	return point(st1.x+(ed1.x-st1.x)*t, st1.y+(ed1.y-st1.y)*t);
}
point get_point(point st, point ed, point p){  //求p在直线（st，ed）上的垂点
	point ans;
	if(fabs(st.x-ed.x) < eps){
		ans.x = st.x;
		ans.y = p.y;
	}else if(fabs(st.y-ed.y) < eps){
		ans.x = p.x;
		ans.y = st.y;
	}else{
		double k1, d1, k2, d2;
		k1 = (ed.y - st.y) / (ed.x - st.x);
		d1 = st.y - st.x*(ed.y - st.y) / (ed.x - st.x);
		k2 = -1 / k1;
		d2 = p.y + p.x / k1;
		ans.x = (d2 - d1) * k1 / (k1*k1 + 1);
		ans.y = k1 * ans.x + d1;
	}
	return ans;
}
struct poly{
	static const int N = 10005; //点数目的最大值
	point ps[N];
	int pn;
	poly(){ pn = 0; }
	void push(point tp){ ps[pn++] = tp; }
	point& operator[](int k){ return ps[k];  }
	bool inPoly(point tp){  //判断点是否在多边形内
		int i;
		for(i = 0; i < 3; i++) ps[pn+i] = ps[i];
		for(i = 0; i < pn; i++){
			if(isOnSegment(ps[i], ps[i+1], tp)) return true;
		}
		int c = 0;  //c表示交点的个数
		point d1 = tp, d2;
		d2.y = d1.y;
		d2.x = maxx;   //构造一条足够长的线段,d2.x要根据数据的范围而定
		for(i = 0; i < pn; i++){
			if(isIntersect(d1, d2, ps[i], ps[i+1])
				||(isOnSegment(d1, d2, ps[i+1]) && (inter_pro(d1, d2, d1, ps[i]) *  inter_pro(d1, d2, d1, ps[i+2]) < 0))
				||(isOnSegment(d1, d2, ps[i+1]) && isOnSegment(d1, d2, ps[i+2]) && (inter_pro(d1, d2, d1, ps[i]) * inter_pro(d1, d2, d1, ps[i+3]) < 0)))
				c++;
		}
		return c % 2;
	}
	double getArea(){  //求面积
		ps[pn] = ps[0];
		int i;
		double ans;
		for(ans = i = 0; i < pn; i++) ans += ps[i].x * ps[i+1].y - ps[i].y * ps[i+1].x;
		return ans;
	}
	void reverse(){  //翻转
		int l = 0, r = pn - 1;
		point tp;
		while(l < r){ tp = ps[l]; ps[l] = ps[r]; ps[r] = tp;  l++, r--; }
	}
	bool isInnerCir(double R, point cen){  //判断半径为R，圆心为cen的圆是否在多边形内
		if(!inPoly(cen)) return false;
		ps[pn] = ps[0];
		int i;
		double dd;
		for(i = 0; i < pn; i++){ if(sign(distPtoL(ps[i], ps[i+1], cen) - R) < 0) return false; }
		return true; //圆在多边形内返回true
	}
	bool canHold(double R){  //判断该多边形是否可以容纳半径为R的圆
		if(sign(getArea()) < 0) reverse();
		ps[pn] = ps[0];
		int i, j;
		double d;
		point p1, p2, p3, p4, cen;
		for(i = 0; i < pn; i++){ //每举两个点
			for(j = i + 1; j < pn; j++){
				if(sign((d = dist(ps[i], ps[j])) - 2*R) <= 0){
					p1 = getMid(ps[i], ps[j]);
					p2 = ps[j] - ps[i];
					swap(p2.x, p2.y);  p2.x *= -1;   //向量逆时针旋转90度
					d = sqrt(R * R - d * d / 4.0);
					cen = goLen(p1, p2,d);
					if(isInnerCir(R, cen)) return true;
					p2.x *= -1; p2.y *= -1;  //反向
					cen = goLen(p1, p2,d);
					if(isInnerCir(R, cen)) return true;
				}
				//枚举线段和线段卡住圆的情况，把两条线段内移R的距离，将两条新的线段的交点作为圆心
				p1 = change(ps[i], ps[i+1], ps[i], R);
				p2 = change(ps[i], ps[i+1], ps[i+1], R);
				p3 = change(ps[j], ps[j+1], ps[j], R);
				p4 = change(ps[j], ps[j+1], ps[j+1], R);
				if(isIntersect(p1, p2, p3, p4)){
					cen = intersectPoint(p1, p2, p3, p4);
					if(isInnerCir(R, cen)) return true;
				}
				//枚举点和线段卡住圆的情况
				p1 = get_point(ps[i], ps[i+1], ps[j]);
				p2 = ps[i+1] - ps[i];
				if(sign(2 * R - dist(p1, ps[j])) < 0) continue;
				d = sqrt(R * R - (ps[j].x - p1.x) * (ps[j].y - p1.y));
				cen = goLen(p1, p2, d);
				if(isInnerCir(R, cen)) return true;
				p2.x *= -1; p2.y *= -1;
				cen = goLen(p1, p2, d);
				if(isInnerCir(R, cen)) return true;
			}
			int k;
		}
		return false;
	}
	double getMax(){ //计算多边形可以容纳的最大半径的圆
		double low, high, mid, step;
		low = 0;
		high = 10000;
		step = 0.001;
		while(low <= high){
			mid = (low + high) * 0.5;
			if(canHold(mid)) low = mid + step;
			else high = mid - step;
		}
		return low - step;
	}
};