code wiki / _hdl_build / nx_geo_area.nx
nx_geo_area.nx source
↩ module page · 52 lines · 2486 B
1// nx_geo_area.nx -- LIB: GEO-010 polygon AREA + ORIENTATION, INTEGER-EXACT (shoelace / Gauss).
2//
3// THE EXCEED ANGLE (measured, not asserted): the area is the shoelace sum of INTEGER cross products
4// over microdegree coordinates -- ZERO floating point. So both the doubled-area magnitude AND the
5// orientation (its sign) are EXACT and DETERMINISTIC, and collinear / degenerate vertices contribute
6// EXACTLY 0. Those are precisely the cases where Turf.js / Shapely / PostGIS (all float) accumulate
7// epsilon error or flip orientation. Resolution is 1 microdegree (~0.11 m); a 1-microdegree apex
8// yields its exact area, not noise.
9//
10// Convention (matches nx_geo): poly = flat [lat0,lon0, lat1,lon1, ...], npts vertices, implicitly
11// closed (last->first); lat = Y, lon = X. The result is the DOUBLED (un-halved) area in microdegree^2
12// (deg^2 * 1e12) so it stays a pure integer -- halving is the caller's choice and exact when the
13// doubled value is even. Overflow-safe for Earth coords in i64 (|lon*lat| <= 1.8e8 * 9e7 = 1.6e16;
14// the sum of hundreds of such terms fits in i64's 9.2e18).
15//
16// Foundation rung: the area-weighted CENTROID and polygon clipping/union compose THIS. license_tier: ORIGINAL
17import "nx_syscalls.nx"
18
19// Signed DOUBLED area via the shoelace formula: sum over edges (prev=j -> curr=i) of
20// (x_j*y_i - x_i*y_j), with x=lon, y=lat. Sign encodes orientation: > 0 counter-clockwise,
21// < 0 clockwise, == 0 degenerate (zero-area / all-collinear). Integer-exact.
22func geo_signed_area2(poly: *i64, npts: i64) -> i64 {
23 var sum: i64 = 0
24 var i: i64 = 0
25 var j: i64 = npts - 1
26 while i < npts {
27 let yi: i64 = poly[i * 2]
28 let xi: i64 = poly[i * 2 + 1]
29 let yj: i64 = poly[j * 2]
30 let xj: i64 = poly[j * 2 + 1]
31 sum = sum + (xj * yi - xi * yj)
32 j = i
33 i = i + 1
34 }
35 return sum
36}
37
38// Absolute DOUBLED area (always >= 0), in microdegree^2. The unsigned magnitude for size queries.
39func geo_area2(poly: *i64, npts: i64) -> i64 {
40 var s: i64 = geo_signed_area2(poly, npts)
41 if s < 0 { s = 0 - s }
42 return s
43}
44
45// Polygon ORIENTATION as an exact integer verdict: +1 counter-clockwise, -1 clockwise, 0 degenerate.
46// This is the integer-exact win -- no epsilon, no platform-dependent flip on near-collinear vertices.
47func geo_orientation(poly: *i64, npts: i64) -> i64 {
48 let s: i64 = geo_signed_area2(poly, npts)
49 if s > 0 { return 1 }
50 if s < 0 { return 0 - 1 }
51 return 0
52}