1/*
2 * 기하 알고리즘 (Computational Geometry)
3 * CCW, 볼록 껍질, 선분 교차, 다각형
4 *
5 * 2D 평면에서의 기하학적 연산입니다.
6 */
7
8#include <stdio.h>
9#include <stdlib.h>
10#include <math.h>
11#include <stdbool.h>
12
13#define EPS 1e-9
14#define PI 3.14159265358979323846
15
16/* =============================================================================
17 * 1. 점과 벡터
18 * ============================================================================= */
19
20typedef struct {
21 double x, y;
22} Point;
23
24Point point_create(double x, double y) {
25 Point p = {x, y};
26 return p;
27}
28
29Point point_add(Point a, Point b) {
30 return point_create(a.x + b.x, a.y + b.y);
31}
32
33Point point_sub(Point a, Point b) {
34 return point_create(a.x - b.x, a.y - b.y);
35}
36
37Point point_scale(Point p, double k) {
38 return point_create(p.x * k, p.y * k);
39}
40
41double point_dot(Point a, Point b) {
42 return a.x * b.x + a.y * b.y;
43}
44
45double point_cross(Point a, Point b) {
46 return a.x * b.y - a.y * b.x;
47}
48
49double point_norm(Point p) {
50 return sqrt(p.x * p.x + p.y * p.y);
51}
52
53double point_dist(Point a, Point b) {
54 return point_norm(point_sub(a, b));
55}
56
57/* =============================================================================
58 * 2. CCW (Counter-Clockwise)
59 * ============================================================================= */
60
61/* 반시계 방향: 1, 시계 방향: -1, 일직선: 0 */
62int ccw(Point a, Point b, Point c) {
63 double cross = point_cross(point_sub(b, a), point_sub(c, a));
64 if (cross > EPS) return 1;
65 if (cross < -EPS) return -1;
66 return 0;
67}
68
69/* =============================================================================
70 * 3. 선분 교차 판정
71 * ============================================================================= */
72
73bool on_segment(Point p, Point a, Point b) {
74 double min_x = (a.x < b.x) ? a.x : b.x;
75 double max_x = (a.x > b.x) ? a.x : b.x;
76 double min_y = (a.y < b.y) ? a.y : b.y;
77 double max_y = (a.y > b.y) ? a.y : b.y;
78
79 return p.x >= min_x - EPS && p.x <= max_x + EPS &&
80 p.y >= min_y - EPS && p.y <= max_y + EPS;
81}
82
83bool segments_intersect(Point a, Point b, Point c, Point d) {
84 int d1 = ccw(a, b, c);
85 int d2 = ccw(a, b, d);
86 int d3 = ccw(c, d, a);
87 int d4 = ccw(c, d, b);
88
89 if (((d1 > 0 && d2 < 0) || (d1 < 0 && d2 > 0)) &&
90 ((d3 > 0 && d4 < 0) || (d3 < 0 && d4 > 0))) {
91 return true;
92 }
93
94 /* 일직선 상의 교차 */
95 if (d1 == 0 && on_segment(c, a, b)) return true;
96 if (d2 == 0 && on_segment(d, a, b)) return true;
97 if (d3 == 0 && on_segment(a, c, d)) return true;
98 if (d4 == 0 && on_segment(b, c, d)) return true;
99
100 return false;
101}
102
103/* 선분 교차점 계산 */
104bool line_intersection(Point a, Point b, Point c, Point d, Point* result) {
105 double a1 = b.y - a.y;
106 double b1 = a.x - b.x;
107 double c1 = a1 * a.x + b1 * a.y;
108
109 double a2 = d.y - c.y;
110 double b2 = c.x - d.x;
111 double c2 = a2 * c.x + b2 * c.y;
112
113 double det = a1 * b2 - a2 * b1;
114 if (fabs(det) < EPS) return false; /* 평행 */
115
116 result->x = (b2 * c1 - b1 * c2) / det;
117 result->y = (a1 * c2 - a2 * c1) / det;
118 return true;
119}
120
121/* =============================================================================
122 * 4. 볼록 껍질 (Convex Hull)
123 * ============================================================================= */
124
125int compare_points(const void* a, const void* b) {
126 Point* p1 = (Point*)a;
127 Point* p2 = (Point*)b;
128 if (fabs(p1->x - p2->x) > EPS)
129 return (p1->x < p2->x) ? -1 : 1;
130 return (p1->y < p2->y) ? -1 : 1;
131}
132
133/* Graham Scan */
134int convex_hull_graham(Point points[], int n, Point hull[]) {
135 if (n < 3) return 0;
136
137 /* 가장 아래, 왼쪽 점 찾기 */
138 int min_idx = 0;
139 for (int i = 1; i < n; i++) {
140 if (points[i].y < points[min_idx].y ||
141 (fabs(points[i].y - points[min_idx].y) < EPS &&
142 points[i].x < points[min_idx].x)) {
143 min_idx = i;
144 }
145 }
146
147 Point pivot = points[min_idx];
148 points[min_idx] = points[0];
149 points[0] = pivot;
150
151 /* 각도 기준 정렬 */
152 for (int i = 1; i < n - 1; i++) {
153 for (int j = i + 1; j < n; j++) {
154 int c = ccw(pivot, points[i], points[j]);
155 if (c < 0 || (c == 0 &&
156 point_dist(pivot, points[i]) > point_dist(pivot, points[j]))) {
157 Point temp = points[i];
158 points[i] = points[j];
159 points[j] = temp;
160 }
161 }
162 }
163
164 /* 스택 사용 */
165 int hull_size = 0;
166 for (int i = 0; i < n; i++) {
167 while (hull_size >= 2 &&
168 ccw(hull[hull_size - 2], hull[hull_size - 1], points[i]) <= 0) {
169 hull_size--;
170 }
171 hull[hull_size++] = points[i];
172 }
173
174 return hull_size;
175}
176
177/* Monotone Chain */
178int convex_hull_monotone(Point points[], int n, Point hull[]) {
179 if (n < 3) return 0;
180
181 qsort(points, n, sizeof(Point), compare_points);
182
183 int hull_size = 0;
184
185 /* 아래 껍질 */
186 for (int i = 0; i < n; i++) {
187 while (hull_size >= 2 &&
188 ccw(hull[hull_size - 2], hull[hull_size - 1], points[i]) <= 0) {
189 hull_size--;
190 }
191 hull[hull_size++] = points[i];
192 }
193
194 /* 위 껍질 */
195 int lower_size = hull_size;
196 for (int i = n - 2; i >= 0; i--) {
197 while (hull_size > lower_size &&
198 ccw(hull[hull_size - 2], hull[hull_size - 1], points[i]) <= 0) {
199 hull_size--;
200 }
201 hull[hull_size++] = points[i];
202 }
203
204 return hull_size - 1; /* 마지막 점은 첫 점과 동일 */
205}
206
207/* =============================================================================
208 * 5. 다각형 연산
209 * ============================================================================= */
210
211/* 다각형 넓이 (신발끈 공식) */
212double polygon_area(Point poly[], int n) {
213 double area = 0;
214 for (int i = 0; i < n; i++) {
215 int j = (i + 1) % n;
216 area += point_cross(poly[i], poly[j]);
217 }
218 return fabs(area) / 2.0;
219}
220
221/* 점이 다각형 내부에 있는지 (Ray Casting) */
222bool point_in_polygon(Point p, Point poly[], int n) {
223 int crossings = 0;
224 for (int i = 0; i < n; i++) {
225 int j = (i + 1) % n;
226 if ((poly[i].y <= p.y && p.y < poly[j].y) ||
227 (poly[j].y <= p.y && p.y < poly[i].y)) {
228 double x = (poly[j].x - poly[i].x) * (p.y - poly[i].y) /
229 (poly[j].y - poly[i].y) + poly[i].x;
230 if (p.x < x) crossings++;
231 }
232 }
233 return (crossings % 2) == 1;
234}
235
236/* 다각형 무게중심 */
237Point polygon_centroid(Point poly[], int n) {
238 double cx = 0, cy = 0, area = 0;
239 for (int i = 0; i < n; i++) {
240 int j = (i + 1) % n;
241 double cross = point_cross(poly[i], poly[j]);
242 cx += (poly[i].x + poly[j].x) * cross;
243 cy += (poly[i].y + poly[j].y) * cross;
244 area += cross;
245 }
246 area /= 2.0;
247 return point_create(cx / (6.0 * area), cy / (6.0 * area));
248}
249
250/* =============================================================================
251 * 6. 가장 가까운 점 쌍 (Closest Pair)
252 * ============================================================================= */
253
254double min_double(double a, double b) {
255 return (a < b) ? a : b;
256}
257
258int compare_by_y(const void* a, const void* b) {
259 Point* p1 = (Point*)a;
260 Point* p2 = (Point*)b;
261 return (p1->y < p2->y) ? -1 : 1;
262}
263
264double closest_pair_recursive(Point px[], Point py[], int n) {
265 if (n <= 3) {
266 double min_dist = 1e18;
267 for (int i = 0; i < n; i++) {
268 for (int j = i + 1; j < n; j++) {
269 double d = point_dist(px[i], px[j]);
270 if (d < min_dist) min_dist = d;
271 }
272 }
273 return min_dist;
274 }
275
276 int mid = n / 2;
277 Point mid_point = px[mid];
278
279 Point* pyl = malloc(mid * sizeof(Point));
280 Point* pyr = malloc((n - mid) * sizeof(Point));
281 int li = 0, ri = 0;
282
283 for (int i = 0; i < n; i++) {
284 if (py[i].x <= mid_point.x && li < mid)
285 pyl[li++] = py[i];
286 else
287 pyr[ri++] = py[i];
288 }
289
290 double dl = closest_pair_recursive(px, pyl, mid);
291 double dr = closest_pair_recursive(px + mid, pyr, n - mid);
292 double d = min_double(dl, dr);
293
294 free(pyl);
295 free(pyr);
296
297 /* 스트립 내 점들 */
298 Point* strip = malloc(n * sizeof(Point));
299 int strip_n = 0;
300 for (int i = 0; i < n; i++) {
301 if (fabs(py[i].x - mid_point.x) < d) {
302 strip[strip_n++] = py[i];
303 }
304 }
305
306 for (int i = 0; i < strip_n; i++) {
307 for (int j = i + 1; j < strip_n && (strip[j].y - strip[i].y) < d; j++) {
308 double dist = point_dist(strip[i], strip[j]);
309 if (dist < d) d = dist;
310 }
311 }
312
313 free(strip);
314 return d;
315}
316
317double closest_pair(Point points[], int n) {
318 Point* px = malloc(n * sizeof(Point));
319 Point* py = malloc(n * sizeof(Point));
320 for (int i = 0; i < n; i++) {
321 px[i] = py[i] = points[i];
322 }
323
324 qsort(px, n, sizeof(Point), compare_points);
325 qsort(py, n, sizeof(Point), compare_by_y);
326
327 double result = closest_pair_recursive(px, py, n);
328
329 free(px);
330 free(py);
331 return result;
332}
333
334/* =============================================================================
335 * 7. 회전 캘리퍼스 (Rotating Calipers)
336 * ============================================================================= */
337
338/* 볼록 다각형의 지름 (가장 먼 두 점 사이 거리) */
339double rotating_calipers(Point hull[], int n) {
340 if (n < 2) return 0;
341 if (n == 2) return point_dist(hull[0], hull[1]);
342
343 int j = 1;
344 double max_dist = 0;
345
346 for (int i = 0; i < n; i++) {
347 Point edge = point_sub(hull[(i + 1) % n], hull[i]);
348
349 while (1) {
350 Point next_edge = point_sub(hull[(j + 1) % n], hull[j]);
351 double cross = point_cross(edge, point_sub(hull[(j + 1) % n], hull[j]));
352 if (cross <= 0) break;
353 j = (j + 1) % n;
354 }
355
356 double d1 = point_dist(hull[i], hull[j]);
357 double d2 = point_dist(hull[(i + 1) % n], hull[j]);
358 if (d1 > max_dist) max_dist = d1;
359 if (d2 > max_dist) max_dist = d2;
360 }
361
362 return max_dist;
363}
364
365/* =============================================================================
366 * 테스트
367 * ============================================================================= */
368
369int main(void) {
370 printf("============================================================\n");
371 printf("기하 알고리즘 예제\n");
372 printf("============================================================\n");
373
374 /* 1. CCW */
375 printf("\n[1] CCW (반시계 방향 판정)\n");
376 Point a = {0, 0}, b = {4, 0}, c = {2, 2};
377 printf(" A(0,0), B(4,0), C(2,2)\n");
378 int result = ccw(a, b, c);
379 printf(" CCW: %d (%s)\n", result,
380 result > 0 ? "반시계" : (result < 0 ? "시계" : "일직선"));
381
382 /* 2. 선분 교차 */
383 printf("\n[2] 선분 교차 판정\n");
384 Point p1 = {0, 0}, p2 = {4, 4}, p3 = {0, 4}, p4 = {4, 0};
385 printf(" 선분1: (0,0)-(4,4), 선분2: (0,4)-(4,0)\n");
386 printf(" 교차: %s\n", segments_intersect(p1, p2, p3, p4) ? "예" : "아니오");
387
388 Point intersection;
389 if (line_intersection(p1, p2, p3, p4, &intersection)) {
390 printf(" 교차점: (%.1f, %.1f)\n", intersection.x, intersection.y);
391 }
392
393 /* 3. 볼록 껍질 */
394 printf("\n[3] 볼록 껍질\n");
395 Point points[] = {{0, 0}, {1, 1}, {2, 2}, {4, 4}, {0, 4}, {4, 0}, {2, 1}, {1, 2}};
396 int n = 8;
397 Point hull[10];
398
399 printf(" 점들: ");
400 for (int i = 0; i < n; i++) {
401 printf("(%.0f,%.0f) ", points[i].x, points[i].y);
402 }
403 printf("\n");
404
405 int hull_size = convex_hull_monotone(points, n, hull);
406 printf(" 볼록 껍질 (%d개): ", hull_size);
407 for (int i = 0; i < hull_size; i++) {
408 printf("(%.0f,%.0f) ", hull[i].x, hull[i].y);
409 }
410 printf("\n");
411
412 /* 4. 다각형 넓이 */
413 printf("\n[4] 다각형 넓이\n");
414 Point triangle[] = {{0, 0}, {4, 0}, {2, 3}};
415 printf(" 삼각형: (0,0), (4,0), (2,3)\n");
416 printf(" 넓이: %.1f\n", polygon_area(triangle, 3));
417
418 Point square[] = {{0, 0}, {4, 0}, {4, 4}, {0, 4}};
419 printf(" 정사각형: (0,0), (4,0), (4,4), (0,4)\n");
420 printf(" 넓이: %.1f\n", polygon_area(square, 4));
421
422 /* 5. 점의 다각형 내부 판정 */
423 printf("\n[5] 점의 다각형 내부 판정\n");
424 Point test1 = {2, 2}, test2 = {5, 5};
425 printf(" 정사각형: (0,0), (4,0), (4,4), (0,4)\n");
426 printf(" 점 (2,2): %s\n", point_in_polygon(test1, square, 4) ? "내부" : "외부");
427 printf(" 점 (5,5): %s\n", point_in_polygon(test2, square, 4) ? "내부" : "외부");
428
429 /* 6. 가장 가까운 점 쌍 */
430 printf("\n[6] 가장 가까운 점 쌍\n");
431 Point closest_points[] = {{2, 3}, {12, 30}, {40, 50}, {5, 1}, {12, 10}, {3, 4}};
432 printf(" 점들: (2,3), (12,30), (40,50), (5,1), (12,10), (3,4)\n");
433 printf(" 최소 거리: %.4f\n", closest_pair(closest_points, 6));
434
435 /* 7. 복잡도 */
436 printf("\n[7] 복잡도\n");
437 printf(" | 알고리즘 | 시간복잡도 |\n");
438 printf(" |-----------------|---------------|\n");
439 printf(" | CCW | O(1) |\n");
440 printf(" | 선분 교차 | O(1) |\n");
441 printf(" | 볼록 껍질 | O(n log n) |\n");
442 printf(" | 다각형 넓이 | O(n) |\n");
443 printf(" | 점 내부 판정 | O(n) |\n");
444 printf(" | 가장 가까운 쌍 | O(n log n) |\n");
445 printf(" | 회전 캘리퍼스 | O(n) |\n");
446
447 printf("\n============================================================\n");
448
449 return 0;
450}