再来一道测半平面交模板题 Poj1279 Art Gallery
地址：http://poj.org/problem?id=1279
题目：
Art Gallery
| Time Limit: 1000MS | Memory Limit: 10000K | |
| Total Submissions: 7329 | Accepted: 2938 |
Description
The art galleries of the new and very futuristic building of the Center for Balkan Cooperation have the form of polygons (not necessarily convex). When a big exhibition is organized, watching over all of the pictures is a big security concern. Your task is that for a given gallery to write a program which finds the surface of the area of the floor, from which each point on the walls of the gallery is visible. On the figure 1. a map of a gallery is given in some co-ordinate system. The area wanted is shaded on the figure 2.
Input
The number of tasks T that your program have to solve will be on the first row of the input file. Input data for each task start with an integer N, 5 <= N <= 1500. Each of the next N rows of the input will contain the co-ordinates of a vertex of the polygon ? two integers that fit in 16-bit integer type, separated by a single space. Following the row with the co-ordinates of the last vertex for the task comes the line with the number of vertices for the next test and so on.
Output
For each test you must write on one line the required surface - a number with exactly two digits after the decimal point (the number should be rounded to the second digit after the decimal point).
Sample Input
1
7
0 0
4 4
4 7
9 7
13 -1
8 -6
4 -4
Sample Output
80.00
思路：没什么好说的，和前面几题一样都是用来测模板的题，不过还是wa了两次，因为把题目看成是要四舍五入到第二位小数（盲人acmer）
1 #include <iostream>
2 #include <cstdio>
3 #include <cmath>
4 #include <algorithm>
5
6
7 using namespace std;
8 const double eps = 1e-8;
9 //点
10 class Point
11 {
12 public:
13 double x, y;
14
15 Point(){}
16 Point(double x, double y):x(x),y(y){}
17
18 bool operator < (const Point &_se) const
19 {
20 return x<_se.x || (x==_se.x && y<_se.y);
21 }
22 /*******判断ta与tb的大小关系*******/
23 static int sgn(double ta,double tb)
24 {
25 if(fabs(ta-tb)<eps)return 0;
26 if(ta<tb) return -1;
27 return 1;
28 }
29 static double xmult(const Point &po, const Point &ps, const Point &pe)
30 {
31 return (ps.x - po.x) * (pe.y - po.y) - (pe.x - po.x) * (ps.y - po.y);
32 }
33 friend Point operator + (const Point &_st,const Point &_se)
34 {
35 return Point(_st.x + _se.x, _st.y + _se.y);
36 }
37 friend Point operator - (const Point &_st,const Point &_se)
38 {
39 return Point(_st.x - _se.x, _st.y - _se.y);
40 }
41 //点位置相同(double类型)
42 bool operator == (const Point &_off) const
43 {
44 return Point::sgn(x, _off.x) == 0 && Point::sgn(y, _off.y) == 0;
45 }
46 //点位置不同(double类型)
47 bool operator != (const Point &_Off) const
48 {
49 return ((*this) == _Off) == false;
50 }
51 //两点间距离的平方
52 static double dis2(const Point &_st,const Point &_se)
53 {
54 return (_st.x - _se.x) * (_st.x - _se.x) + (_st.y - _se.y) * (_st.y - _se.y);
55 }
56 //两点间距离
57 static double dis(const Point &_st, const Point &_se)
58 {
59 return sqrt((_st.x - _se.x) * (_st.x - _se.x) + (_st.y - _se.y) * (_st.y - _se.y));
60 }
61 };
62 //两点表示的向量
63 class Line
64 {
65 public:
66 Point s, e;//两点表示，起点[s]，终点[e]
67 double a, b, c;//一般式,ax+by+c=0
68
69 Line(){}
70 Line(const Point &s, const Point &e):s(s),e(e){}
71 Line(double _a,double _b,double _c):a(_a),b(_b),c(_c){}
72
73 //向量与点的叉乘,参数：点[_Off]
74 //[点相对向量位置判断]
75 double operator /(const Point &_Off) const
76 {
77 return (_Off.y - s.y) * (e.x - s.x) - (_Off.x - s.x) * (e.y - s.y);
78 }
79 //向量与向量的叉乘,参数：向量[_Off]
80 friend double operator /(const Line &_st,const Line &_se)
81 {
82 return (_st.e.x - _st.s.x) * (_se.e.y - _se.s.y) - (_st.e.y - _st.s.y) * (_se.e.x - _se.s.x);
83 }
84 friend double operator *(const Line &_st,const Line &_se)
85 {
86 return (_st.e.x - _st.s.x) * (_se.e.x - _se.s.x) - (_st.e.y - _st.s.y) * (_se.e.y - _se.s.y);
87 }
88 //从两点表示转换为一般表示
89 //a=y2-y1,b=x1-x2,c=x2*y1-x1*y2
90 bool pton()
91 {
92 a = e.y - s.y;
93 b = s.x - e.x;
94 c = e.x * s.y - e.y * s.x;
95 return true;
96 }
97
98 //-----------点和直线（向量）-----------
99 //点在向量左边（右边的小于号改成大于号即可,在对应直线上则加上=号）
100 //参数：点[_Off],向量[_Ori]
101 friend bool operator<(const Point &_Off, const Line &_Ori)
102 {
103 return (_Ori.e.y - _Ori.s.y) * (_Off.x - _Ori.s.x)
104 < (_Off.y - _Ori.s.y) * (_Ori.e.x - _Ori.s.x);
105 }
106
107 //点在直线上,参数：点[_Off]
108 bool lhas(const Point &_Off) const
109 {
110 return Point::sgn((*this) / _Off, 0) == 0;
111 }
112 //点在线段上,参数：点[_Off]
113 bool shas(const Point &_Off) const
114 {
115 return lhas(_Off)
116 && Point::sgn(_Off.x - min(s.x, e.x), 0) > 0 && Point::sgn(_Off.x - max(s.x, e.x), 0) < 0
117 && Point::sgn(_Off.y - min(s.y, e.y), 0) > 0 && Point::sgn(_Off.y - max(s.y, e.y), 0) < 0;
118 }
119
120 //点到直线/线段的距离
121 //参数： 点[_Off], 是否是线段[isSegment](默认为直线)
122 double dis(const Point &_Off, bool isSegment = false)
123 {
124 ///化为一般式
125 pton();
126
127 //到直线垂足的距离
128 double td = (a * _Off.x + b * _Off.y + c) / sqrt(a * a + b * b);
129
130 //如果是线段判断垂足
131 if(isSegment)
132 {
133 double xp = (b * b * _Off.x - a * b * _Off.y - a * c) / ( a * a + b * b);
134 double yp = (-a * b * _Off.x + a * a * _Off.y - b * c) / (a * a + b * b);
135 double xb = max(s.x, e.x);
136 double yb = max(s.y, e.y);
137 double xs = s.x + e.x - xb;
138 double ys = s.y + e.y - yb;
139 if(xp > xb + eps || xp < xs - eps || yp > yb + eps || yp < ys - eps)
140 td = min(Point::dis(_Off,s), Point::dis(_Off,e));
141 }
142
143 return fabs(td);
144 }
145
146 //关于直线对称的点
147 Point mirror(const Point &_Off) const
148 {
149 ///注意先转为一般式
150 Point ret;
151 double d = a * a + b * b;
152 ret.x = (b * b * _Off.x - a * a * _Off.x - 2 * a * b * _Off.y - 2 * a * c) / d;
153 ret.y = (a * a * _Off.y - b * b * _Off.y - 2 * a * b * _Off.x - 2 * b * c) / d;
154 return ret;
155 }
156 //计算两点的中垂线
157 static Line ppline(const Point &_a, const Point &_b)
158 {
159 Line ret;
160 ret.s.x = (_a.x + _b.x) / 2;
161 ret.s.y = (_a.y + _b.y) / 2;
162 //一般式
163 ret.a = _b.x - _a.x;
164 ret.b = _b.y - _a.y;
165 ret.c = (_a.y - _b.y) * ret.s.y + (_a.x - _b.x) * ret.s.x;
166 //两点式
167 if(std::fabs(ret.a) > eps)
168 {
169 ret.e.y = 0.0;
170 ret.e.x = - ret.c / ret.a;
171 if(ret.e == ret. s)
172 {
173 ret.e.y = 1e10;
174 ret.e.x = - (ret.c - ret.b * ret.e.y) / ret.a;
175 }
176 }
177 else
178 {
179 ret.e.x = 0.0;
180 ret.e.y = - ret.c / ret.b;
181 if(ret.e == ret. s)
182 {
183 ret.e.x = 1e10;
184 ret.e.y = - (ret.c - ret.a * ret.e.x) / ret.b;
185 }
186 }
187 return ret;
188 }
189
190 //------------直线和直线（向量）-------------
191 //直线重合,参数：直线向量[_st],[_se]
192 static bool equal(const Line &_st, const Line &_se)
193 {
194 return _st.lhas(_se.e) && _se.lhas(_se.s);
195 }
196 //直线平行，参数：直线向量[_st],[_se]
197 static bool parallel(const Line &_st,const Line &_se)
198 {
199 return Point::sgn(_st / _se, 0) == 0;
200 }
201 //两直线（线段）交点，参数：直线向量[_st],[_se],交点
202 //返回-1代表平行，0代表重合，1代表相交
203 static bool crossLPt(const Line &_st,const Line &_se,Point &ret)
204 {
205 if(Line::parallel(_st,_se))
206 {
207 if(Line::equal(_st,_se)) return 0;
208 return -1;
209 }
210 ret = _st.s;
211 double t = (Line(_st.s,_se.s)/_se)/(_st/_se);
212 ret.x += (_st.e.x - _st.s.x) * t;
213 ret.y += (_st.e.y - _st.s.y) * t;
214 return 1;
215 }
216 //------------线段和直线（向量）----------
217 //线段和直线交
218 //参数：直线[_st],线段[_se]
219 friend bool crossSL(const Line &_st,const Line &_se)
220 {
221 return Point::sgn((_st / _se.s) * (_st / _se.e) ,0) <= 0;
222 }
223
224 //------------线段和线段（向量）----------
225 //判断线段是否相交(注意添加eps)，参数：线段[_st],线段[_se]
226 static bool isCrossSS(const Line &_st,const Line &_se)
227 {
228 //1.快速排斥试验判断以两条线段为对角线的两个矩形是否相交
229 //2.跨立试验（等于0时端点重合）
230 return
231 max(_st.s.x, _st.e.x) >= min(_se.s.x, _se.e.x) &&
232 max(_se.s.x, _se.e.x) >= min(_st.s.x, _st.e.x) &&
233 max(_st.s.y, _st.e.y) >= min(_se.s.y, _se.e.y) &&
234 max(_se.s.y, _se.e.y) >= min(_st.s.y, _st.e.y) &&
235 Point::sgn((_st / Line(_st.s, _se.s)) * (_st / Line(_st.s, _se.e)), 0) <= 0 &&
236 Point::sgn((_se / Line(_se.s, _st.s)) * (_se / Line(_se.s, _st.e)), 0) <= 0;
237 }
238 };
239 class Polygon
240 {
241 public:
242 const static int maxpn = 2000;
243 Point pt[maxpn];//点（顺时针或逆时针）
244 int n;//点的个数
245
246 Point& operator[](int _p)
247 {
248 return pt[_p];
249 }
250
251 //求多边形面积，多边形内点必须顺时针或逆时针
252 double area() const
253 {
254 double ans = 0.0;
255 for(int i = 0; i < n; i ++)
256 {
257 int nt = (i + 1) % n;
258 ans += pt[i].x * pt[nt].y - pt[nt].x * pt[i].y;
259 }
260 return fabs(ans / 2.0);
261 }
262 //求多边形重心，多边形内点必须顺时针或逆时针
263 Point gravity() const
264 {
265 Point ans;
266 ans.x = ans.y = 0.0;
267 double area = 0.0;
268 for(int i = 0; i < n; i ++)
269 {
270 int nt = (i + 1) % n;
271 double tp = pt[i].x * pt[nt].y - pt[nt].x * pt[i].y;
272 area += tp;
273 ans.x += tp * (pt[i].x + pt[nt].x);
274 ans.y += tp * (pt[i].y + pt[nt].y);
275 }
276 ans.x /= 3 * area;
277 ans.y /= 3 * area;
278 return ans;
279 }
280 //判断点在凸多边形内，参数：点[_Off]
281 bool chas(const Point &_Off) const
282 {
283 double tp = 0, np;
284 for(int i = 0; i < n; i ++)
285 {
286 np = Line(pt[i], pt[(i + 1) % n]) / _Off;
287 if(tp * np < -eps)
288 return false;
289 tp = (fabs(np) > eps)?np: tp;
290 }
291 return true;
292 }
293 //判断点是否在任意多边形内[射线法]，O(n)
294 bool ahas(const Point &_Off) const
295 {
296 int ret = 0;
297 double infv = 1e-10;//坐标系最大范围
298 Line l = Line(_Off, Point( -infv ,_Off.y));
299 for(int i = 0; i < n; i ++)
300 {
301 Line ln = Line(pt[i], pt[(i + 1) % n]);
302 if(fabs(ln.s.y - ln.e.y) > eps)
303 {
304 Point tp = (ln.s.y > ln.e.y)? ln.s: ln.e;
305 if(fabs(tp.y - _Off.y) < eps && tp.x < _Off.x + eps)
306 ret ++;
307 }
308 else if(Line::isCrossSS(ln,l))
309 ret ++;
310 }
311 return (ret % 2 == 1);
312 }
313 //凸多边形被直线分割,参数：直线[_Off]
314 Polygon split(Line _Off)
315 {
316 //注意确保多边形能被分割
317 Polygon ret;
318 Point spt[2];
319 double tp = 0.0, np;
320 bool flag = true;
321 int i, pn = 0, spn = 0;
322 for(i = 0; i < n; i ++)
323 {
324 if(flag)
325 pt[pn ++] = pt[i];
326 else
327 ret.pt[ret.n ++] = pt[i];
328 np = _Off / pt[(i + 1) % n];
329 if(tp * np < -eps)
330 {
331 flag = !flag;
332 Line::crossLPt(_Off,Line(pt[i], pt[(i + 1) % n]),spt[spn++]);
333 }
334 tp = (fabs(np) > eps)?np: tp;
335 }
336 ret.pt[ret.n ++] = spt[0];
337 ret.pt[ret.n ++] = spt[1];
338 n = pn;
339 return ret;
340 }
341
342
343 /** 卷包裹法求点集凸包，_p为输入点集，_n为点的数量 **/
344 void ConvexClosure(Point _p[],int _n)
345 {
346 sort(_p,_p+_n);
347 n=0;
348 for(int i=0;i<_n;i++)
349 {
350 while(n>1&&Point::sgn(Line(pt[n-2],pt[n-1])/Line(pt[n-2],_p[i]),0)<=0)
351 n--;
352 pt[n++]=_p[i];
353 }
354 int _key=n;
355 for(int i=_n-2;i>=0;i--)
356 {
357 while(n>_key&&Point::sgn(Line(pt[n-2],pt[n-1])/Line(pt[n-2],_p[i]),0)<=0)
358 n--;
359 pt[n++]=_p[i];
360 }
361 if(n>1) n--;//除去重复的点，该点已是凸包凸包起点
362 }
363 // /****** 寻找凸包的graham 扫描法********************/
364 // /****** _p为输入的点集,_n为点的数量****************/
365 // /**使用时需把gmp函数放在类外，并且看情况修改pt[0]**/
366 // bool gcmp(const Point &ta,const Point &tb)/// 选取与最后一条确定边夹角最小的点，即余弦值最大者
367 // {
368 // double tmp=Line(pt[0],ta)/Line(pt[0],tb);
369 // if(Point::sgn(tmp,0)==0)
370 // return Point::dis(pt[0],ta)<Point::dis(pt[0],tb);
371 // else if(tmp>0)
372 // return 1;
373 // return 0;
374 // }
375 // void graham(Point _p[],int _n)
376 // {
377 // int cur=0;
378 // for(int i=1;i<_n;i++)
379 // if(Point::sgn(_p[cur].y,_p[i].y)>0 || (Point::sgn(_p[cur].y,_p[i].y)==0 && Point::sgn(_p[cur].x,_p[i].x)>0))
380 // cur=i;
381 // swap(_p[cur],_p[0]);
382 // n=0,pt[n++]=_p[0];
383 // if(_n==1) return;
384 // sort(_p+1,_p+_n,Polygon::gcmp);
385 // pt[n++]=_p[1],pt[n++]=_p[2];
386 // for(int i=3;i<_n;i++)
387 // {
388 // while(Point::sgn(Line(pt[n-2],pt[n-1])/Line(pt[n-2],_p[i]),0)<0)
389 // n--;
390 // pt[n++]=_p[i];
391 // }
392 // }
393 //凸包旋转卡壳(注意点必须顺时针或逆时针排列)
394 //返回值凸包直径的平方（最远两点距离的平方）
395 double rotating_calipers()
396 {
397 int i = 1;
398 double ret = 0.0;
399 pt[n] = pt[0];
400 for(int j = 0; j < n; j ++)
401 {
402 while(fabs(Point::xmult(pt[i+1],pt[j], pt[j + 1])) > fabs(Point::xmult(pt[i],pt[j], pt[j + 1])) + eps)
403 i = (i + 1) % n;
404 //pt[i]和pt[j],pt[i + 1]和pt[j + 1]可能是对踵点
405 ret = (ret, max(Point::dis(pt[i],pt[j]), Point::dis(pt[i + 1],pt[j + 1])));
406 }
407 return ret;
408 }
409
410 //凸包旋转卡壳(注意点必须逆时针排列)
411 //返回值两凸包的最短距离
412 double rotating_calipers(Polygon &_Off)
413 {
414 int i = 0;
415 double ret = 1e10;//inf
416 pt[n] = pt[0];
417 _Off.pt[_Off.n] = _Off.pt[0];
418 //注意凸包必须逆时针排列且pt[0]是左下角点的位置
419 while(_Off.pt[i + 1].y > _Off.pt[i].y)
420 i = (i + 1) % _Off.n;
421 for(int j = 0; j < n; j ++)
422 {
423 double tp;
424 //逆时针时为 >,顺时针则相反
425 while((tp = Point::xmult(_Off.pt[i + 1],pt[j], pt[j + 1]) - Point::xmult(_Off.pt[i], pt[j], pt[j + 1])) > eps)
426 i = (i + 1) % _Off.n;
427 //(pt[i],pt[i+1])和(_Off.pt[j],_Off.pt[j + 1])可能是最近线段
428 ret = min(ret, Line(pt[j], pt[j + 1]).dis(_Off.pt[i], true));
429 ret = min(ret, Line(_Off.pt[i], _Off.pt[i + 1]).dis(pt[j + 1], true));
430 if(tp > -eps)//如果不考虑TLE问题最好不要加这个判断
431 {
432 ret = min(ret, Line(pt[j], pt[j + 1]).dis(_Off.pt[i + 1], true));
433 ret = min(ret, Line(_Off.pt[i], _Off.pt[i + 1]).dis(pt[j], true));
434 }
435 }
436 return ret;
437 }
438
439 //-----------半平面交-------------
440 //复杂度:O(nlog2(n))
441 //#include <algorithm>
442 //半平面计算极角函数[如果考虑效率可以用成员变量记录]
443 static double hpc_pa(const Line &_Off)
444 {
445 return atan2(_Off.e.y - _Off.s.y, _Off.e.x - _Off.s.x);
446 }
447 //半平面交排序函数[优先顺序: 1.极角 2.前面的直线在后面的左边]
448 static bool hpc_cmp(const Line &l, const Line &r)
449 {
450 double lp = hpc_pa(l), rp = hpc_pa(r);
451 if(fabs(lp - rp) > eps)
452 return lp < rp;
453 return Point::xmult(r.s,l.s, r.e) < -eps;
454 }
455 static int judege(const Line &_lx,const Line &_ly,const Line &_lz)
456 {
457 Point tmp;
458 Line::crossLPt(_lx,_ly,tmp);
459 return Point::sgn(Point::xmult(_lz.s,tmp,_lz.e),0);
460 }
461 //获取半平面交的多边形（多边形的核）
462 //参数：向量集合[l]，向量数量[ln];(半平面方向在向量左边）
463 //函数运行后如果n[即返回多边形的点数量]为0则不存在半平面交的多边形（不存在区域或区域面积无穷大）
464 Polygon& halfPanelCross(Line _Off[], int ln)
465 {
466 Line dequeue[maxpn];//用于计算的双端队列
467 int i, tn;
468 sort(_Off, _Off + ln, hpc_cmp);
469 //平面在向量左边的筛选
470 for(i = tn = 1; i < ln; i ++)
471 if(fabs(hpc_pa(_Off[i]) - hpc_pa(_Off[i - 1])) > eps)
472 _Off[tn ++] = _Off[i];
473 ln = tn,n = 0;
474 int bot = 0, top = 1;
475 dequeue[0] = _Off[0];
476 dequeue[1] = _Off[1];
477 for(i = 2; i < ln; i ++)
478 {
479 while(bot < top && Polygon::judege(dequeue[top],dequeue[top-1],_Off[i]) > 0)
480 top --;
481 while(bot < top && Polygon::judege(dequeue[bot],dequeue[bot+1],_Off[i]) > 0)
482 bot ++;
483 dequeue[++ top] = _Off[i];
484 }
485
486 while(bot < top && Polygon::judege(dequeue[top],dequeue[top-1],dequeue[bot]) > 0)
487 top --;
488 while(bot < top && Polygon::judege(dequeue[bot],dequeue[bot+1],dequeue[top]) > 0)
489 bot ++;
490 //计算交点(注意不同直线形成的交点可能重合)
491 if(top <= bot + 1)
492 return (*this);
493 for(i = bot; i < top; i ++)
494 Line::crossLPt(dequeue[i],dequeue[i + 1],pt[n++]);
495 if(bot < top + 1)
496 Line::crossLPt(dequeue[bot],dequeue[top],pt[n++]);
497 return (*this);
498 }
499 };
500
501 int n,t;
502 Point pt[2000];
503 Line ln[2000];
504 Polygon ans;
505 int main(void)
506 {
507 scanf("%d",&t);
508 while(t--)
509 {
510 scanf("%d",&n);
511 for(int i=0;i<n;i++)
512 scanf("%lf%lf",&pt[i].x,&pt[i].y);
513 pt[n++]=pt[0];
514 for(int i=n-1;i;i--)
515 ln[i-1]=Line(pt[i],pt[i-1]);
516 //for(int i=0;i<n-1;i++)
517 // printf("%.2f %.2f %.2f %.2f\n",ln[i].s.x,ln[i].s.y,ln[i].e.x,ln[i].e.y);
518 ans.halfPanelCross(ln,n-1);
519 double area=0;
520 if(ans.n==0)
521 area=0;
522 else
523 area=ans.area();
524 printf("%.2f\n",area);
525 }
526 return 0;
527 }
作者：weeping
出处：www.cnblogs.com/weeping/
本文版权归作者和博客园共有，欢迎转载，但未经作者同意必须保留此段声明，且在文章页面明显位置给出原文连接，否则保留追究法律责任的权利。
