#!/usr/bin/env python3
"""Build base_parks2.json from the paged raw2/*.json Socrata batches.
Preserves MultiPolygon structure (outer + holes); normalises winding:
outer rings CCW (positive shoelace), holes CW (negative shoelace)."""
import json, glob, os

D   = '/sessions/upbeat-quirky-clarke/mnt/outputs/sfmap'
RAW = os.path.join(D, 'raw2')
OUT = os.path.join(D, 'base_parks2.json')
LON0, LON1 = -122.52, -122.35
LAT0, LAT1 = 37.70, 37.84
MIN_ACRES = 0.2

def shoelace(ring):
    s = 0.0
    for i in range(len(ring) - 1):
        x1, y1 = ring[i]; x2, y2 = ring[i+1]
        s += x1*y2 - x2*y1
    return s / 2.0

def clean_ring(ring):
    out = []
    for x, y in ring:
        p = [round(float(x), 6), round(float(y), 6)]
        if not out or out[-1] != p:
            out.append(p)
    if len(out) >= 2 and out[0] != out[-1]:
        out.append(list(out[0]))
    while len(out) >= 3 and out[-1] == out[-2]:
        out.pop()
    return out if len(out) >= 4 else None

rows = []
for f in sorted(glob.glob(os.path.join(RAW, '*.json'))):
    rows.extend(json.load(open(f)))

parks = []
stats = dict(no_geom=0, out_of_bbox=0, small=0, degen_rings=0, degen_polys=0)

for r in rows:
    name, acres, ptype = r['n'], float(r['a']), r.get('t') or ''
    g = r.get('g')
    if not g or not g.get('coordinates'):
        stats['no_geom'] += 1; continue
    if acres < MIN_ACRES:
        stats['small'] += 1; continue
    polys = []
    for poly in g['coordinates']:
        if not poly: continue
        outer = clean_ring(poly[0])
        if outer is None:
            stats['degen_polys'] += 1; continue
        if shoelace(outer) < 0: outer.reverse()          # outer -> CCW
        holes = []
        for hr in poly[1:]:
            h = clean_ring(hr)
            if h is None:
                stats['degen_rings'] += 1; continue
            if shoelace(h) > 0: h.reverse()              # hole -> CW
            holes.append(h)
        polys.append({'outer': outer, 'holes': holes})
    if not polys:
        stats['degen_polys'] += 1; continue
    xs = [p[0] for pg in polys for p in pg['outer']]
    ys = [p[1] for pg in polys for p in pg['outer']]
    if max(xs) < LON0 or min(xs) > LON1 or max(ys) < LAT0 or min(ys) > LAT1:
        stats['out_of_bbox'] += 1
        print('  bbox-dropped:', name, round(min(xs),3), round(min(ys),3)); continue
    parks.append({'n': name, 'acres': acres, 'type': ptype, 'polys': polys})

parks.sort(key=lambda p: -p['acres'])
json.dump({'parks': parks}, open(OUT, 'w'), separators=(',', ':'))

npoly = sum(len(p['polys']) for p in parks)
nhole = sum(len(pg['holes']) for p in parks for pg in p['polys'])
nvert = sum(len(r) for p in parks for pg in p['polys'] for r in [pg['outer']]+pg['holes'])
sz = os.path.getsize(OUT)
print('input rows   :', len(rows))
print('parks        :', len(parks))
print('polygons     :', npoly)
print('holes        :', nhole)
print('vertices     :', nvert)
print('file size    : %d bytes (%.1f KB)' % (sz, sz/1024))
print('stats        :', stats)