code wiki / _hdl_build / nx_geo_buffer.nx

nx_geo_buffer.nx source

↩ module page · 37 lines · 1843 B

1// nx_geo_buffer.nx -- LIB: GEO-017 BUFFER / OFFSET MEASURE (Minkowski sum of a convex polygon with a 2// disk of radius r). For a CONVEX polygon the buffered area has the exact closed form 3// area(P (+) disk_r) = area(P) + perimeter(P)*r + pi*r^2 4// (the straight edges sweep rectangles of total area perimeter*r; the corners sweep wedges that sum to 5// one full disk pi*r^2). We compute this MEASURE in integer arithmetic, doubled to match GEO-010: 6// buffer_area2 = geo_area2(P) + 2*perimeter*r + 2*pi*r^2, 2*pi*r^2 ~ 710*r^2/113 (pi ~ 355/113). 7// 8// HONEST SCOPE: this returns the buffer COVERAGE MEASURE (doubled microdeg^2), not the rounded buffer 9// BOUNDARY polygon. It is DETERMINISTIC integer arithmetic; the only inexactness is pi via the 355/113 10// rational (error ~2.7e-7) and nx_isqrt's floor on non-perfect-square edge lengths -- no float drift. 11// This is the "expand a geofence by R meters" coverage query. Convention: poly = flat [lat0,lon0,...], 12// r in the same microdeg units. license_tier: ORIGINAL 13import "nx_geo_area.nx" 14import "nx_isqrt.nx" 15import "nx_syscalls.nx" 16 17// polygon perimeter = sum of edge lengths (microdeg), via nx_isqrt of each edge's squared length. 18func geo_perimeter(poly: *i64, npts: i64) -> i64 { 19 var per: i64 = 0 20 var i: i64 = 0 21 var j: i64 = npts - 1 22 while i < npts { 23 let dx: i64 = poly[i * 2 + 1] - poly[j * 2 + 1] 24 let dy: i64 = poly[i * 2] - poly[j * 2] 25 per = per + nx_isqrt(dx * dx + dy * dy) 26 j = i 27 i = i + 1 28 } 29 return per 30} 31 32// doubled area (microdeg^2) of the convex polygon buffered outward by radius r. r==0 -> geo_area2(P). 33func geo_buffer_area2(poly: *i64, npts: i64, r: i64) -> i64 { 34 let a2: i64 = geo_area2(poly, npts) 35 let per: i64 = geo_perimeter(poly, npts) 36 return a2 + 2 * per * r + 710 * r * r / 113 37}