计算几何板子

        这里放计算几何的板子。关于原理,推荐一篇很好的文章
        长期更新。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
#include <algorithm>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <iostream>
#include <map>
#include <queue>
#include <set>
#include <stack>
#include <vector>
using namespace std;
const double PI = acos(-1.0); //π
const double EPS = 1e-6; //精度控制
const double INF = 1e100; //无穷大
int sgn(double s) //判断浮点数的符号
{
if (fabs(s) < EPS) return 0;
if (s > 0) return 1;
return -1;
}
int cmp(double x, double y) //判断x和y的大小关系,-1表示x<y
{
if (fabs(x - y) < EPS) return 0;
if (x < y) return -1;
return 1;
}
double r2d(double rad) //弧度转角度
{
return rad / PI * 180.0;
}
double d2r(double degree) //角度转弧度
{
return degree / 180.0 * PI;
}
struct Point //点与向量
{
public:
double x, y;
Point(double a = 0, double b = 0)
: x(a), y(b) {}
void input() //输入
{
scanf("%lf%lf", &x, &y);
}
void debug() //输出
{
printf("(%lf,%lf)\n", x, y);
}
Point operator+(Point b) //向量相加
{
return Point(x + b.x, y + b.y);
}
Point operator-(Point b) //向量相减
{
return Point(x - b.x, y - b.y);
}
Point operator*(double b) //向量数乘
{
return Point(x * b, y * b);
}
double operator*(Point b) //向量内积
{
return x * b.x + y * b.y;
}
Point operator/(double b) //除法
{
return Point(x / b, y / b);
}
bool operator<(const Point& b) //大小比较
{
if (cmp(x, b.x) == 0) return y < b.y;
return x < b.x;
}
bool operator==(const Point& b) //判等
{
return cmp(x, b.x) == 0 && cmp(y, b.y) == 0;
}
double cross(Point b) const //向量叉乘
{
return x * b.y - y * b.x;
}
double angle() //极角,弧度(-PI, PI]
{
return atan2(y, x);
}
double length() //模
{
return sqrt((*this) * (*this));
}
double angleTo(Point b) //向量角度,弧度
{
return acos(((*this) * b) / length() / b.length());
}
double angleTo(Point a, Point b) //从该点看其它两点的角度,弧度
{
return Point(a - *this).angleTo(b - *this);
}
double area(Point a, Point b) //求三点组成平行四边形有向面积
{
return (a - *this).cross(b - *this);
}
double dist(Point b) //求两点间欧式距离
{
return (*this - b).length();
}
double manhattanDis(Point b) //求两点间曼哈顿距离
{
return fabs(x - b.x) + fabs(y - b.y);
}
Point rotate(Point p, double ang) //绕点p逆时针旋转ang弧度,不修改本身的值
{
Point v = (*this) - p;
double c = cos(ang), s = sin(ang);
return p + Point(v.x * c - v.y * s, v.x * s + v.y * c);
}
Point rotate(double ang) //直接逆时针旋转ang弧度
{
return rotate(Point(0, 0), ang);
}
Point rotLeft() //左转90
{
return Point(-y, x);
}
Point rotRight() //右转90
{
return Point(y, -x);
}
Point regular() //化为单位向量
{
return *this / length();
}
Point normal() //左转90的单位法向量
{
return rotLeft().regular();
}
Point trunc(double r) //化为长度为r的向量
{
return regular() * r;
}
bool parrelTo(Point b) //向量平行
{
return sgn(cross(b)) == 0;
}
bool verticalTo(Point b) //向量垂直
{
return sgn((*this) * b) == 0;
}
};
typedef Point Vector;
int toLeftTest(Point from, Point to, Point test) //测试test在不在from->to的左边,1在左边,-1在右边,0共线
{
return sgn((to - from).cross(test - to));
}
struct Line
{
Point pos;
Vector to;
Line(Point p, Vector t) //顶点和方向向量
: pos(p), to(t)
{
}
static Line tp(Point a, Point b) //两点构造
{
return Line(a, b - a);
}
Line(double x1, double y1, double x2, double y2) //两点构造
{
pos = Point(x1, y1), to = Point(x2, y2) - Point(x1, y1);
}
Point point(double s) //获取某处的一个点
{
return pos + to * s;
}
double angle() //获得直线倾角,[0,PI)
{
return acos(to.regular().x * sgn(to.regular().y));
}
int onLine(Point s) //判断点与直线的位置关系0在直线外;1在直线上
{
return sgn(to.cross(s - pos)) == 0;
}
bool operator==(Line l) //两条直线是否是同一条(重合)
{
return l.onLine(point(0)) && l.onLine(point(1));
}
bool operator!=(Line l)
{
return !((*this) == l);
}
bool parrelTo(Line l) //直线平行
{
return to.parrelTo(l.to);
}
bool verticalTo(Line l) //直线垂直
{
return to.verticalTo(l.to);
}
Point insWithLine(Line l) //与另一条直线的交点
{
return point(l.to.cross(pos - l.pos) / to.cross(l.to));
}
double disToPoint(Point s) //点到直线的距离
{
return fabs(to.cross(s - pos) / to.length());
}
Point projection(Point s) //点在直线上的投影
{
return point((to * (s - pos)) / (to * to));
}
};
double disToSegment(Point s, Point a, Point b) //点到线段ab最短距离
{
if (a == b) return (s - a).length();
Vector v1 = b - a, v2 = s - a, v3 = s - b;
if (sgn(v1 * v2) < 0) return v2.length();
if (sgn(v1 * v3) > 0) return v3.length();
return Line(a, b - a).disToPoint(s);
}
double area(Point a, Point b, Point c) //三角形面积
{
return fabs((b - a).cross(c - a) / 2);
}
bool onSegment(Point p, Point a, Point b) //点是否在线段上
{
return (sgn((a - p).cross(b - p)) == 0 && sgn((a - p) * (b - p)) < 0) || p == a || p == b;
}
int insWithSegment(Line l, Point a, Point b, Point& ans) //直线与线段相交
{
if (l.parrelTo(Line(a, b - a))) return -1; //不相交
ans = l.insWithLine(Line(a, b - a));
if (!onSegment(ans, a, b)) return -1;
if (ans == a || ans == b) return 0; //与端点相交
return 1; //相交
}
bool segmentIns(Point a1, Point a2, Point b1, Point b2) //判断线段是否相交
{
double c1 = (a2 - a1).cross(b1 - a1), c2 = (a2 - a1).cross(b2 - a1);
double c3 = (b2 - b1).cross(a1 - b1), c4 = (b2 - b1).cross(a2 - b1);
if (!sgn(c1) || !sgn(c2) || !sgn(c3) || !sgn(c4)) { //控制是否可以在顶点处相交
bool f1 = onSegment(b1, a1, a2);
bool f2 = onSegment(b2, a1, a2);
bool f3 = onSegment(a1, b1, b2);
bool f4 = onSegment(a2, b1, b2);
bool f = (f1 | f2 | f3 | f4);
return f;
}
return sgn(c1) * sgn(c2) < 0 && sgn(c3) * sgn(c4) < 0;
}
double polygonArea(Point* p, int n) //计算多边形有向面积,点需要按顺序排列
{
double s = 0;
for (int i = 1; i < n - 1; i++) s += (p[i] - p[0]).cross(p[i + 1] - p[0]) / 2;
return s;
}
double polygonPerimeter(Point* p, int n) //计算多边形周长,点需要按顺序排列
{
double ans = 0;
p[n] = p[0];
for (int i = 0; i < n; i++) ans += p[i].dist(p[i + 1]);
return ans;
}
void sortByPolarAngle(Point at, Point* begin, Point* end) //极角逆时针排序,以at为极点
{
sort(begin, end, [&](Point a, Point b) {
double ans = (a - at).cross(b - at);
return sgn(ans) == 0 ? (a - at).length() <= (b - at).length() : sgn(ans) > 0;
});
}
int gramhamScan(Point* res, int n, Point* ans) //求凸包,O(nlogn),ans存放答案,返回凸包上点的数量
{
if (n < 1) return 0;
Point p = res[0];
int k = 0, top = 1;
for (int i = 1; i < n; i++) {
if (res[i] < p) p = res[i], k = i;
}
swap(res[0], res[k]), sortByPolarAngle(res[0], res + 1, res + n), ans[0] = res[0];
for (int i = 1; i < n; i++) {
while (top > 2 && sgn((ans[top - 1] - ans[top - 2]).cross(res[i] - ans[top - 2])) <= 0) --top;
ans[top++] = res[i];
}
return top;
}
double rotatingCalipers(Point* poly, int n) //旋转卡壳,返回凸包直径
{
double ans = 0;
int j = 2;
poly[n] = poly[0];
for (int i = 0; i < n; i++) {
while (area(poly[i], poly[i + 1], poly[j]) < area(poly[i], poly[i + 1], poly[j + 1])) {
++j;
if (j > n) j = 0;
}
ans = max(ans, max(poly[i].dist(poly[j]), poly[i + 1].dist(poly[j])));
}
return ans;
}
int isPointInPolygon(Point p, Point* poly, int n) //点在多边形内1;在外-1;在边界上0
{
int wn = 0;
for (int i = 0; i < n; ++i) {
if (onSegment(p, poly[i], poly[(i + 1) % n])) return 0;
int k = sgn((poly[(i + 1) % n] - poly[i]).cross(p - poly[i]));
int d1 = sgn(poly[i].y - p.y);
int d2 = sgn(poly[(i + 1) % n].y - p.y);
if (k > 0 && d1 <= 0 && d2 > 0) wn++;
if (k < 0 && d2 <= 0 && d1 > 0) wn--;
}
if (wn != 0) return 1;
return -1;
}

0%