首页 > ACM题库 > HDU-杭电 > HDU 3982-Harry Potter and J.K.Rowling-计算几何-[解题报告]HOJ
2015
04-14

HDU 3982-Harry Potter and J.K.Rowling-计算几何-[解题报告]HOJ

Harry Potter and J.K.Rowling

问题描述 :

In July 31st, last month, the author of the famous novel series J.K.Rowling celebrated her 46th birthday. Many friends gave their best wishes. They gathered together and shared one large beautiful cake.
Rowling had lots of friends, and she had a knife to cut the cake into many pieces. On the cake there was a cherry. After several cuts, the piece with the cherry was left for Rowling. Before she enjoyed it, she wondered how large this piece was, i.e., she wondered how much percentage of the cake the piece with the only cherry has.

输入:

First line has an integer T, indicating the number of test cases.
T test cases follow. The first line of each test case has one number r (1 <= r <= 10000) and one integer n (0 <= n <= 2000), indicating the radius of the cake and the number of knife cuts. n lines follow. The i-th line has four numbers, x1, y1, x2, y2. (x1, y1) and (x2, y2) are two points on the i-th cut. (-10000 <= x1, y1, x2, y2 <= 10000) The last line of each case has two number x, y, indicating the coordinate(x, y) of the cherry.

Technical Specification

1. R, x1, y2, x2, y2, x, y all have two digits below decimal points.
2. The center of the cake is always on the point (0, 0).
3. The cherry was always on the cake and would not be on the knife cuts.

输出:

First line has an integer T, indicating the number of test cases.
T test cases follow. The first line of each test case has one number r (1 <= r <= 10000) and one integer n (0 <= n <= 2000), indicating the radius of the cake and the number of knife cuts. n lines follow. The i-th line has four numbers, x1, y1, x2, y2. (x1, y1) and (x2, y2) are two points on the i-th cut. (-10000 <= x1, y1, x2, y2 <= 10000) The last line of each case has two number x, y, indicating the coordinate(x, y) of the cherry.

Technical Specification

1. R, x1, y2, x2, y2, x, y all have two digits below decimal points.
2. The center of the cake is always on the point (0, 0).
3. The cherry was always on the cake and would not be on the knife cuts.

样例输入:

1
1.00 2
-1.00 0.00 1.00 0.00
0.00 -1.00 0.00 1.00
0.50 0.50

样例输出:

Case 1: 25.00000%

 昨晚到现在,终于A掉,上午在J2EE课上想清楚了,不过下午还是拍着数据才调出来BUG的 = =。。。这水平。。。到区域赛遇到计算几何神题肿么办吧。

这题一看就有思路啊,很裸的思路,先求半平面交(交得的面积是樱桃所在的那个大块块),然后求得的多边形和圆形求交。


1、半平面交的话,需要把线段方向改下,都变成有效区域为樱桃所在区域。

2、有一个坑就是如果切痕都切不到蛋糕,是需要输出100%的,直接特判下。

3、半平面交我用的N*LOGN的算法,这个就不多说了。但是记得需要初始化区域,我的区域是(-2*r, -2*r) — (2*r, 2*r),这个矩形区域。

关键是求多边形和圆形的交。经过半平面交,交的多边形一定是凸多边形,那么就计算凸多边形和圆形的交。

先求出凸多边形和圆的交点(包括在凸多边形在圆内的顶点),根据樱桃位置极角排序。

然后找到第一个顶点为起始顶点,找相邻顶点和第一个顶点组成的三角形,这个三角形一定是内部的面积。

关键是找弓形面积,如下图,蓝色线表示的是多边形和圆形围城的区域,最终和圆的交点应该是B D E。

当计算到DE这条线段时,做DE的中垂线,交圆与BC(绿色点),判断C或者B是否在多边形内部,如果在的话,说明弓形DE是需要计算的面积。

这里有个计算下面面积和上面面积之分。因为我的计算默认的是逆时针,所以我的计算的有效面积一定是左侧面积,如果需要计算弓形面积,那么计算的一定是DE右侧的面积(这个需要好好想下),然后就需要判断DE是优弧还是劣弧,如果是圆心在DE的右侧,则是优弧,反之为劣弧。这个用叉积判断下就好。

注意计算角度不能用asin – -。。。asin的范围是[-π/2, π/2]。。。(哎,我是挫人。。。)

还有就是,如果最后只是两个交点,注意不要计算两次。。。T T。。可以把线段变换下(因为只有两个点的话,需要保证樱桃在线段的右侧)

Harry Potter and J.K.Rowling

#include <set>
#include <map>
#include <queue>
#include <stack>
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <iostream>
#include <limits.h>
#include <string.h>
#include <string>
#include <algorithm>
#define MID(x,y) ( ( x + y ) >> 1 )
#define L(x) ( x << 1 )
#define R(x) ( x << 1 | 1 )
#define FOR(i,s,t) for(int i=s; i<t; i++)
#define BUG puts("here!!!")

using namespace std;

const int MAX = 2010;
struct point { 
	double x,y;
	void P(double xx, double yy){ x = xx; y = yy;}
	void get(){ scanf("%lf%lf",&x, &y);}
};
struct line { point a,b; double ang;};  
point s[MAX*2], ch, p[MAX*2];  
line l[MAX*2],deq[MAX*2];  
const double eps = 1e-6;  
const double pi = acos(-1.0);
bool dy(double x,double y)  {   return x > y + eps;} // x > y   
bool xy(double x,double y)  {   return x < y - eps;} // x < y   
bool dyd(double x,double y) {   return x > y - eps;} // x >= y   
bool xyd(double x,double y) {   return x < y + eps;}     // x <= y   
bool dd(double x,double y)  {   return fabs( x - y ) < eps;}  // x == y  
double crossProduct(point a,point b,point c)//向量 ac 在 ab 的方向   
{  
    return (c.x - a.x)*(b.y - a.y) - (b.x - a.x)*(c.y - a.y);  
}  
double disp2p(point a,point b) //  a b 两点之间的距离 
{
	return sqrt( ( a.x - b.x ) * ( a.x - b.x ) + ( a.y - b.y ) * ( a.y - b.y ) );
}
bool parallel(line u,line v)  
{  
    return dd( (u.a.x - u.b.x)*(v.a.y - v.b.y) - (v.a.x - v.b.x)*(u.a.y - u.b.y) , 0.0 );  
}  
point l2l_inst_p(line l1,line l2)  
{  
    point ans = l1.a;  
    double t = ((l1.a.x - l2.a.x)*(l2.a.y - l2.b.y) - (l1.a.y - l2.a.y)*(l2.a.x - l2.b.x))/  
               ((l1.a.x - l1.b.x)*(l2.a.y - l2.b.y) - (l1.a.y - l1.b.y)*(l2.a.x - l2.b.x));  
    ans.x += (l1.b.x - l1.a.x)*t;  
    ans.y += (l1.b.y - l1.a.y)*t;  
    return ans;  
}  
bool equal_ang(line a,line b)   // 第一次unique的比较函数   
{  
    return dd(a.ang,b.ang);  
}  
bool cmphp(line a,line b)   // 排序的比较函数   
{  
    if( dd(a.ang,b.ang) ) return xy(crossProduct(b.a,b.b,a.a),0.0);  
    return xy(a.ang,b.ang);  
}  
bool equal_p(point a,point b)//第二次unique的比较函数   
{  
    return dd(a.x,b.x) && dd(a.y,b.y);  
}
void makeline_hp(double x1,double y1,double x2,double y2,line &l)
{
	l.a.x = x1; l.a.y = y1; l.b.x = x2; l.b.y = y2;
	l.ang = atan2(y2 - y1,x2 - x1);
}
void inst_hp_nlogn(line *ln,int n,point *s,int &len)  
{  
    len = 0;  
    sort(ln,ln+n,cmphp);  
    n = unique(ln,ln+n,equal_ang) - ln;  
    int bot = 0,top = 1;  
    deq[0] = ln[0]; deq[1] = ln[1];  
    for(int i=2; i<n; i++)  
    {  
        if( parallel(deq[top],deq[top-1]) || parallel(deq[bot],deq[bot+1]) )  
            return ;  
        while( bot < top && dy(crossProduct(ln[i].a,ln[i].b,  
            l2l_inst_p(deq[top],deq[top-1])),0.0) )  
            top--;  
        while( bot < top && dy(crossProduct(ln[i].a,ln[i].b,  
            l2l_inst_p(deq[bot],deq[bot+1])),0.0) )  
            bot++;  
        deq[++top] = ln[i];  
    }  
    while( bot < top && dy(crossProduct(deq[bot].a,deq[bot].b,  
        l2l_inst_p(deq[top],deq[top-1])),0.0) ) top--;  
    while( bot < top && dy(crossProduct(deq[top].a,deq[top].b,  
        l2l_inst_p(deq[bot],deq[bot+1])),0.0) ) bot++;  
    if( top <= bot + 1 ) return ;  
      
    for(int i=bot; i<top; i++)  
        s[len++] = l2l_inst_p(deq[i],deq[i+1]);  
    if( bot < top + 1 ) s[len++] = l2l_inst_p(deq[bot],deq[top]);  
    len = unique(s,s+len,equal_p) - s;  
}
point l2l_inst_p(point u1,point u2,point v1,point v2)
{
	point ans = u1;
	double t = ((u1.x - v1.x)*(v1.y - v2.y) - (u1.y - v1.y)*(v1.x - v2.x))/
				((u1.x - u2.x)*(v1.y - v2.y) - (u1.y - u2.y)*(v1.x - v2.x));
	ans.x += (u2.x - u1.x)*t;
	ans.y += (u2.y - u1.y)*t;
	return ans;
}
void l2c_inst_p(point c,double r,point l1,point l2,point &p1,point &p2)
{
	point p = c;
	double t;
	p.x += l1.y - l2.y;
	p.y += l2.x - l1.x;
	p = l2l_inst_p(p,c,l1,l2);
	t = sqrt(r*r - disp2p(p,c)*disp2p(p,c))/disp2p(l1,l2);
	p1.x = p.x + (l2.x - l1.x)*t;
	p1.y = p.y + (l2.y - l1.y)*t;
	p2.x = p.x - (l2.x - l1.x)*t;
	p2.y = p.y - (l2.y - l1.y)*t;
}
double disp2l(point a,point l1,point l2)
{
	return fabs( crossProduct(a,l1,l2) )/disp2p(l1,l2);
}
bool onSegment(point a, point b, point c)
{
	if( dd(crossProduct(a,b,c),0.0) && dyd(c.x,min(a.x,b.x)) && 
		xyd(c.x,max(a.x,b.x)) && dyd(c.y,min(a.y,b.y)) && xyd(c.y,max(a.y,b.y)) )
		return true;
	return false;
}
bool cmp(point a,point b)
{
	double t1 = atan2(a.y - ch.y, a.x - ch.x);
	double t2 = atan2(b.y - ch.y, b.x - ch.x);
	return xy(t1, t2);
}
double area_triangle(point a,point b,point c)
{
	return fabs( crossProduct(a,b,c) )/2.0;
}
double gong_area(point c,point a,point b,double r)
{
	if( dd(crossProduct(c,a,b), 0.0) )
		return r*r*pi/2;
	double area = area_triangle(a, b, c);
	double ang = acos(1-disp2p(a,b)*disp2p(a,b)/(r*r*2));
	return ang*r*r/2 - area;
}
bool incircle(point a,point c, double r)
{
	return xy(disp2p(a, c), r);
}
bool pin_convexh(point *p,int n,point a)
{
	p[n] = p[0]; p[n+1] = p[1];
	for(int i=0; i<n; i++)
		if( xy(crossProduct(p[i],p[i+1],a)*
			crossProduct(p[i+1],p[i+2],a),0.0) )
			return false;
	return true;
}
point foot_line(point a,point l1,point l2)
{										
	point c;
    c.x = a.x - l2.y + l1.y;
    c.y = a.y + l2.x - l1.x;
	return c;	
}
double solve(point *s, int len, double r)
{
	point c; c.x = c.y = 0;
	point a, b;
	int cnt = 0;
	s[len] = s[0];
	FOR(i, 0, len)	// 求得圆形和多边形的交的多边形的点 
	{	
		if( xy(disp2l(c, s[i], s[i+1]), r) )
		{
			l2c_inst_p(c, r, s[i], s[i+1], a, b);
			if( onSegment(s[i], s[i+1], a) )
				p[cnt++] = a;
			if( onSegment(s[i], s[i+1], b) )
				p[cnt++] = b;
		}
		if( incircle(s[i], c, r) )
			p[cnt++] = s[i];
		if( incircle(s[i+1], c, r) )
			p[cnt++] = s[i+1];
	}
	sort(p, p+cnt, cmp);
	cnt = unique(p, p+cnt, equal_p) - p;
	
	p[cnt] = p[0];

	double area = 0.0;
	if( cnt == 2 ) 	// 如果有两个点的话 是不能计算两次的 
	{
		if( xy(crossProduct(p[0],p[1],ch), 0.0) )
			swap(p[0], p[1]);
		cnt--;
	}
	FOR(i, 0, cnt)
	{
		area += area_triangle(p[0], p[i], p[i+1]);
		if( !dd(disp2p(p[i], c), r) || !dd(disp2p(p[i+1], c), r) )
			continue; 
		point mid;
		mid.P((p[i].x + p[i+1].x)/2, (p[i].y + p[i+1].y)/2);
		point a, b, k;
		
		k = foot_line(c, p[i], p[i+1]);
		l2c_inst_p(c, r, k, c, a, b);

		bool flag = false;
		// 判断两个交点是否在线段右边而且在多边形内 
		if( dyd(crossProduct(p[i], p[i+1], a), 0.0) && pin_convexh(s, len, a) )
			flag = true; 
		if( dyd(crossProduct(p[i], p[i+1], b), 0.0) && pin_convexh(s, len, b) )
			flag = true; 
		if( !flag ) continue;
		
		double area_gong = gong_area(c, p[i], p[i+1], r);
		if( dy(crossProduct(p[i], p[i+1], c), 0.0) )
			area_gong = pi*r*r - area_gong;
		area += area_gong;
	}
	return area/(pi*r*r);
}

int main()
{
	int ncases, ind = 1, n;
	double r;
	scanf("%d", &ncases);
	
	while( ncases-- )
	{
		scanf("%lf%d",&r, &n);
		FOR(i, 0, n)
		{
			l[i].a.get();
			l[i].b.get();
		}
		ch.get();

		bool f = false;
		point c; c.P(0,0);
		FOR(i, 0, n)
		{
			if( dy(crossProduct(l[i].a, l[i].b, ch), 0.0) )
				swap(l[i].a, l[i].b);
			l[i].ang = atan2(l[i].b.y - l[i].a.y, l[i].b.x - l[i].a.x);
			if( xy(disp2l(c,l[i].a, l[i].b), r) )
				f = true;		 
		}
		if( !f )		// 如果输入的切痕都切不到蛋糕,要输出100%
		{
			printf("Case %d: %.5lf%%\n",ind++, 100.0);
			continue;
		}
		
		makeline_hp(-2*r, -2*r, 2*r, -2*r, l[n++]); 
		makeline_hp(2*r, -2*r, 2*r, 2*r, l[n++]);
		makeline_hp(2*r, 2*r, -2*r, 2*r, l[n++]);
		makeline_hp(-2*r, 2*r, -2*r, -2*r, l[n++]);

		int len;
		inst_hp_nlogn(l, n, s, len); 
		
		double ans = solve(s, len, r);
		
		printf("Case %d: %.5lf%%\n",ind++, ans*100);
	}

return 0;
}

版权声明:本文为博主原创文章,未经博主允许不得转载。

参考:http://blog.csdn.net/zxy_snow/article/details/6739561


  1. 5.1处,反了;“上一个操作符的优先级比操作符ch的优先级大,或栈是空的就入栈。”如代码所述,应为“上一个操作符的优先级比操作符ch的优先级小,或栈是空的就入栈。”