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
|
// Flächeninhalt eines Polygons (nicht selbstschneidend).
// Punkte gegen den Uhrzeigersinn: positiv, sonst negativ.
double area(const vector<pt>& poly) { //poly[0] == poly.back()
double res = 0;
for (int i = 0; i + 1 < sz(poly); i++)
res += cross(poly[i], poly[i + 1]);
return 0.5 * res;
}
// Anzahl drehungen einer Polyline um einen Punkt
// p nicht auf rand und poly[0] == poly.back()
// res != 0 or (res & 1) != 0 um inside zu prüfen bei
// selbstschneidenden Polygonen (definitions Sache)
ll windingNumber(pt p, const vector<pt>& poly) {
ll res = 0;
for (int i = 0; i + 1 < sz(poly); i++) {
pt a = poly[i], b = poly[i + 1];
if (real(a) > real(b)) swap(a, b);
if (real(a) <= real(p) && real(p) < real(b) &&
cross(p, a, b) < 0) {
res += orientation(p, poly[i], poly[i + 1]);
}}
return res;
}
// Testet, ob ein Punkt im Polygon liegt (beliebige Polygone).
// Ändere Zeile 32 falls rand zählt, poly[0] == poly.back()
bool inside(pt p, const vector<pt>& poly) {
bool in = false;
for (int i = 0; i + 1 < sz(poly); i++) {
pt a = poly[i], b = poly[i + 1];
if (pointOnLineSegment(a, b, p)) return false;
if (real(a) > real(b)) swap(a,b);
if (real(a) <= real(p) && real(p) < real(b) &&
cross(p, a, b) < 0) {
in ^= 1;
}}
return in;
}
// convex hull without duplicates, h[0] != h.back()
// apply comments if border counts as inside
bool inside(pt p, const vector<pt>& hull) {
int l = 0, r = sz(hull) - 1;
if (cross(hull[0], hull[r], p) >= 0) return false; // > 0
while (l + 1 < r) {
int m = (l + r) / 2;
if (cross(hull[0], hull[m], p) > 0) l = m; // >= 0
else r = m;
}
return cross(hull[l], hull[r], p) > 0; // >= 0
}
void rotateMin(vector<pt>& hull) {
auto mi = min_element(all(hull), [](const pt& a, const pt& b){
return real(a) == real(b) ? imag(a) < imag(b)
: real(a) < real(b);
});
rotate(hull.begin(), mi, hull.end());
}
// convex hulls without duplicates, h[0] != h.back()
vector<pt> minkowski(vector<pt> ps, vector<pt> qs) {
rotateMin(ps);
rotateMin(qs);
ps.push_back(ps[0]);
qs.push_back(qs[0]);
ps.push_back(ps[1]);
qs.push_back(qs[1]);
vector<pt> res;
for (ll i = 0, j = 0; i + 2 < sz(ps) || j + 2 < sz(qs);) {
res.push_back(ps[i] + qs[j]);
auto c = cross(ps[i + 1] - ps[i], qs[j + 1] - qs[j]);
if(c >= 0) i++;
if(c <= 0) j++;
}
return res;
}
// convex hulls without duplicates, h[0] != h.back()
double dist(const vector<pt>& ps, const vector<pt>& qs) {
for (pt& q : qs) q *= -1;
auto p = minkowski(ps, qs);
p.push_back(p[0]);
double res = 0.0;
//bool intersect = true;
for (ll i = 0; i + 1 < sz(p); i++) {
//intersect &= cross(p[i], p[i+1] - p[i]) <= 0;
res = max(res, cross(p[i], p[i+1]-p[i]) / abs(p[i+1]-p[i]));
}
return res;
}
bool left(pt of, pt p) {return cross(p, of) < 0 ||
(cross(p, of) == 0 && dot(p, of) > 0);}
// convex hulls without duplicates, hull[0] == hull.back() and
// hull[0] must be a convex point (with angle < pi)
// returns index of corner where dot(dir, corner) is maximized
int extremal(const vector<pt>& hull, pt dir) {
dir *= pt(0, 1);
int l = 0, r = sz(hull) - 1;
while (l + 1 < r) {
int m = (l + r) / 2;
pt dm = hull[m+1]-hull[m];
pt dl = hull[l+1]-hull[l];
if (left(dl, dir) != left(dl, dm)) {
if (left(dl, dm)) l = m;
else r = m;
} else {
if (cross(dir, dm) < 0) l = m;
else r = m;
}}
return r;
}
// convex hulls without duplicates, hull[0] == hull.back() and
// hull[0] must be a convex point (with angle < pi)
// {} if no intersection
// {x} if corner is only intersection
// {a, b} segments (a,a+1) and (b,b+1) intersected (if only the
// border is intersected corners a and b are the start and end)
vector<int> intersect(const vector<pt>& hull, pt a, pt b) {
int endA = extremal(hull, (a-b) * pt(0, 1));
int endB = extremal(hull, (b-a) * pt(0, 1));
// cross == 0 => line only intersects border
if (cross(hull[endA], a, b) > 0 ||
cross(hull[endB], a, b) < 0) return {};
int n = sz(hull) - 1;
vector<int> res;
for (auto _ : {0, 1}) {
int l = endA, r = endB;
if (r < l) r += n;
while (l + 1 < r) {
int m = (l + r) / 2;
if (cross(hull[m % n], a, b) <= 0 &&
cross(hull[m % n], a, b) != hull(poly[endB], a, b))
l = m;
else r = m;
}
if (cross(hull[r % n], a, b) == 0) l++;
res.push_back(l % n);
swap(endA, endB);
swap(a, b);
}
if (res[0] == res[1]) res.pop_back();
return res;
}
|