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}