HDU3644 2010 Asia Regional Hangzhou Site Online Contest D A Chocolate Manufacturer's Problem
1、判断多边形有向边顺、逆时针：
取最右点p[i]，p[i-1]->p[i]与p[i]->p[i+1]成右手关系则为逆时针。
2、判断点与简单多边形位置关系：
参考文献：王学军,沈连婠,朱绍源等.基于左边的点在简单多边形内的判别算法[J].机械工程师,2006,(2):53-54.
给定一个简单多边形，判别点u在多边形G内外的判断算法步骤：
Step1：过u作一水平射线；
Step1：过u作一水平射线；
Step2：求出射线与简单多边形G的交点；
Step3：若交点为0，则u在G外，结束；
Step4：若交点数大于0，求出与u点距离最近交点v；
Step5：找到交点v所在边的两顶点，记p1为p[i]，p2为p[i+1]；
Step6：判断交点v是否是顶点，若v与p[i]重合，则p1为p[i-1]，若v与p[i+1]重合，则p2为p[i+2]；
Step7：判断 vu 到 vp1 沿逆时针方向的角度是否小于 vp2到vp1沿逆时针方向的角度，若是，则u在G内，否则在外，结束。
3、模拟退火，在边上选采样点可避免直接采样到多边形外，其他方面在方法上没有特别之处。参数调整比较纠结，改天编译器变了或者HDU数据变了这代码就不一定能AC了。
1 #include<stdio.h> 2 #include<string.h> 3 #include<stdlib.h> 4 #include<math.h> 5 #include<time.h> 6 const int maxn = 55; 7 const double eps = 1e-8; 8 const double pi = acos(-1.0); 9 int dcmp(double x) 10 { 11 if(x > eps) return 1; 12 return x < -eps ? -1 : 0; 13 } 14 inline double min(double a, double b) 15 {return a < b ? a : b;} 16 inline double max(double a, double b) 17 {return a > b ? a : b;} 18 inline double Sqr(double x) 19 {return x * x;} 20 /**************************************************/ 21 //点及关于点、线的操作 22 struct Point 23 { 24 double x, y; 25 Point(){x = y = 0;} 26 Point(double a, double b) 27 {x = a, y = b;} 28 inline Point operator-(const Point &b)const 29 {return Point(x - b.x, y - b.y);} 30 inline Point operator+(const Point &b)const 31 {return Point(x + b.x, y + b.y);} 32 inline Point operator-() 33 {return Point(-x, -y);} 34 inline Point operator*(const double &b)const 35 {return Point(x * b, y * b);} 36 inline double dot(const Point &b)const 37 {return x * b.x + y * b.y;} 38 inline double cross(const Point &b, const Point &c)const 39 {return (b.x - x) * (c.y - y) - (c.x - x) * (b.y - y);} 40 inline double Dis(const Point &b)const 41 {return sqrt((*this - b).dot(*this - b));} 42 inline bool operator<(const Point &b)const 43 { 44 if(!dcmp(x - b.x)) return y < b.y; 45 return x < b.x; 46 } 47 inline bool operator>(const Point &b)const 48 {return b < *this;} 49 inline bool operator==(const Point &b)const 50 {return !dcmp(x - b.x) && !dcmp(y - b.y);} 51 inline bool InLine(const Point &b, const Point &c)const 52 {return !dcmp(cross(b, c));} 53 inline bool OnSeg(const Point &b, const Point &c)const//包括端点 54 {return InLine(b, c) && (*this - c).dot(*this - b) < eps;} 55 inline bool InSeg(const Point &b, const Point &c)const//不包括端点 56 {return InLine(b, c) && (*this - c).dot(*this - b) < -eps;} 57 double ToSeg(const Point&, const Point&)const; 58 }; 59 /**************************************************/ 60 bool Parallel(const Point &a, const Point &b, const Point &c, const Point &d) 61 {return !dcmp(a.cross(b, a + d - c));} 62 Point LineCross(const Point &a, const Point &b, const Point &c, const Point &d) 63 { 64 double u = a.cross(b, c), v = b.cross(a, d); 65 return Point((c.x * v + d.x * u) / (u + v), (c.y * v + d.y * u) / (u + v)); 66 } 67 double Point::ToSeg(const Point &b, const Point &c)const 68 { 69 Point t(x + b.y - c.y, y + c.x - b.x); 70 if(cross(t, b) * cross(t, c) > eps) 71 return min(Dis(b), Dis(c)); 72 return Dis(LineCross(*this, t, b, c)); 73 } 74 bool SegCross(const Point &a, const Point &b, 75 const Point &c, const Point &d, Point &p)//包括端点 76 { 77 if(Parallel(a, b, c, d)) return false; 78 if(a.cross(b, c) * a.cross(b, d) <= 0 && c.cross(d, a) * c.cross(d, b) <= 0) 79 { 80 p = LineCross(a, b, c, d); 81 return true; 82 } 83 return false; 84 } 85 /**************************************************/ 86 Point p[maxn]; 87 int n; 88 double r; 89 /**************************************************/ 90 double AngCounterClock(double start, double end) 91 {return end - start + (end > start - eps ? 2 * pi : 0);} 92 bool InSimplePolygon(Point u, Point p[], int n/*double neg_inf*/) 93 //判断点在简单多边形内，不包括边界，多边形点为逆时针 94 { 95 double neg_inf = -1e20, angvu, angvp1, angvp2; 96 Point v(0, u.y), p1, p2, tmp; 97 int i, id;//距离u最近交点对应线段起始点id 98 bool flag = false; 99 p[n] = p[0]; 100 //设u->v为水平负向射线 101 /* for(i = 0, neg_inf = 0; i < n; ++ i) neg_inf = min(neg_inf, p[i].x); 102 neg_inf -= 100.0; */ 103 v.x = neg_inf; 104 for(i = 0; i < n; ++ i) 105 { 106 if(u.OnSeg(p[i], p[i + 1])) return false; 107 if(!SegCross(u, v, p[i], p[i + 1], tmp)) continue; 108 flag = true; 109 if(tmp.x - v.x > eps) v.x = tmp.x, id = i; 110 } 111 if(!flag) return false; 112 p1 = v == p[id] ? p[(id + n - 1) % n] : p[id]; 113 p2 = v == p[id + 1] ? p[(id + 2) % n] : p[id + 1]; 114 angvu = atan2(u.y - v.y, u.x - v.x); 115 angvp1 = atan2(p1.y - v.y, p1.x - v.x); 116 angvp2 = atan2(p2.y - v.y, p2.x - v.x); 117 return AngCounterClock(angvu, angvp1) < 118 AngCounterClock(angvp2, angvp1) - eps; 119 } 120 /**************************************************/ 121 double PtoPolygon(Point u, Point p[], int n)//点到多边形边、点的最近距离 122 { 123 int i; 124 double res = 1e120; 125 p[n] = p[0]; 126 for(i = 0; i < n; ++ i) 127 res = min(res, u.ToSeg(p[i], p[i + 1])); 128 return res; 129 } 130 /**************************************************/ 131 //模拟退火 132 const int maxsam = 20; 133 const int pacelen = 5; 134 const double depace = 0.57; 135 const double endeps = 1e-3; 136 Point sam[maxsam]; 137 double mindis[maxsam]; 138 bool SA(double r) 139 { 140 int i, j, samnum; 141 Point tmp; 142 double tmpdis; 143 double maxx = -1e30, maxy = -1e30, minx = 1e30, miny = 1e30; 144 for(i = 0; i < n; ++ i) 145 { 146 maxx = max(p[i].x, maxx); 147 maxy = max(p[i].y, maxy); 148 minx = min(p[i].x, minx); 149 miny = min(p[i].y, miny); 150 } 151 samnum = maxsam < n ? maxsam : n; 152 double pace = sqrt(Sqr(maxx - minx) + Sqr(maxy - miny)); 153 for(i = 0; i < samnum; ++ i)//随机选边，随机比例选择边上的点作为采样点 154 { 155 j = rand() % n; 156 sam[i] = p[j] + (p[j + 1] - p[j]) * (rand() / 32767.0); 157 mindis[i] = 0; 158 } 159 for(; pace > endeps; pace *= depace) 160 for(i = 0; i < samnum; ++ i) 161 for(j = 0; j < pacelen; ++ j) 162 { 163 tmp.x = sam[i].x + cos((double)rand()) * pace; 164 tmp.y = sam[i].y + cos((double)rand()) * pace; 165 if(!InSimplePolygon(tmp, p, n)) continue; 166 tmpdis = PtoPolygon(tmp, p, n); 167 if(tmpdis > r - endeps) return true; 168 if(tmpdis > mindis[i]) 169 { 170 sam[i] = tmp; 171 mindis[i] = tmpdis; 172 } 173 } 174 return false; 175 } 176 /**************************************************/ 177 void MakeCounterClock(Point p[], int n)//顺时针则反转多边形的有向边 178 { 179 int i, id = 0; 180 p[n] = p[0]; 181 for(i = 0; i < n; ++ i) if(p[i].x > p[id].x) id = i; 182 if(p[id].cross(p[id + 1], p[(id + n - 1) % n]) > eps) return; 183 Point tmp; 184 for(i = n - 1 >> 1; i >= 0; -- i) 185 tmp = p[i], p[i] = p[n - i - 1], p[n - i - 1] = tmp; 186 } 187 int main() 188 { 189 srand(419); 190 while(scanf("%d", &n) && n) 191 { 192 for(int i = 0; i < n; ++ i) 193 scanf("%lf%lf", &p[i].x, &p[i].y); 194 MakeCounterClock(p, n); 195 scanf("%lf", &r); 196 printf(SA(r) ? "Yes\n" : "No\n"); 197 } 198 return 0; 199 }
