#!/usr/bin/env python3
"""
House numbers for address search (Unmapped Search v2, 11 Nov 2026): the points OpenStreetMap tags with
addr:housenumber in South Africa, each put on the named street it belongs to. From OpenStreetMap ONLY, like
scripts/db/address-osm.py: this file is the published method for the house numbers in the derived database offered
under the ODbL at unmapped.co.za/open (ODbL 4.6). Nothing else is joined in (no City of Cape Town address points, no
Municipal Demarcation Board ward, no other source).

    <scratch>/venv/bin/python -I scripts/db/address-houses.py --pbf <south-africa-YYMMDD.osm.pbf> --streets <dir>/osm-streets.ndjson --out <dir>

Writes <dir>/osm-houses.ndjson (one street a line: its id, name and its numbers with their points) and
<dir>/osm-houses-stats.json.

What it keeps: the number and its point. Never a house name, a building name, an owner, an occupant or a contributor.
  - A node, or a building or plot drawn as a way (its middle), tagged addr:housenumber.
  - Numbers as people type them: 1 to 5 digits and at most one letter ("12", "12A"). "12-14" and "12;14" count as both
    numbers at the same point. Anything else ("Unit 3", "Erf 221", a name) is dropped and counted.
  - Its street: addr:street, matched to a street of exactly that name (folded as the street search folds names) no
    more than 400 m away, the nearest piece. A number with no addr:street is put on a street only when one named
    street is within 40 m and no street of another name is within twice that distance; otherwise it is dropped and
    counted (we don't guess which of two streets a corner house faces).
A street here is the street search's street (address-osm.py: same-named ways within about 150 m), so its id is the
piece id without its piece number.
"""
import argparse
import collections
import json
import math
import re
import sys
import unicodedata

import osmium

ABBREV = {'st': 'street', 'str': 'street', 'rd': 'road', 'ave': 'avenue', 'av': 'avenue', 'dr': 'drive', 'cres': 'crescent', 'cr': 'crescent',
          'ln': 'lane', 'pl': 'place', 'blvd': 'boulevard', 'hwy': 'highway', 'cl': 'close', 'tce': 'terrace'}
NUM = re.compile(r'^(\d{1,5})([A-Za-z]?)$')
CELL = 0.01          # about 1 km
WITH_STREET_M = 400  # addr:street given: a street of that name this near
ALONE_M = 40         # no addr:street: one named street this near, and none of another name within twice the distance


def fold(s):
    s = unicodedata.normalize('NFKD', s)
    s = ''.join(c for c in s if not unicodedata.combining(c))
    s = re.sub(r"['’‘`´ʼ!]", '', s.lower())
    words = re.sub(r'[^a-z0-9]+', ' ', s).strip().split(' ')
    return ' '.join(ABBREV.get(w, w) if i > 0 else w for i, w in enumerate(words))


def numbers(v):
    """The house numbers in a tag's value: ['12'], ['12A'], ['12', '14'] for "12-14"; [] when it isn't a plain number."""
    parts = [p.strip() for p in re.split(r'[;,/]|\s+-\s+|(?<=\d)-(?=\d)', v) if p.strip()]
    out = []
    for p in parts[:4]:
        m = NUM.match(p)
        if not m:
            return []
        out.append(f'{int(m.group(1))}{m.group(2).upper()}')
    return [n for n in out if not n.startswith('0')]


def seg_m(px, py, ax, ay, bx, by, c):
    """Metres from p to the segment a-b (degrees in, a flat local approximation: fine under a kilometre)."""
    ax, bx, px = ax * c, bx * c, px * c
    dx, dy = bx - ax, by - ay
    t = 0.0 if dx == 0 and dy == 0 else max(0.0, min(1.0, ((px - ax) * dx + (py - ay) * dy) / (dx * dx + dy * dy)))
    return math.hypot(px - (ax + t * dx), py - (ay + t * dy)) * 110_574


def main():
    ap = argparse.ArgumentParser()
    ap.add_argument('--pbf', required=True)
    ap.add_argument('--streets', required=True)
    ap.add_argument('--out', required=True)
    a = ap.parse_args()
    stats = collections.Counter()

    # The street search's pieces, in a grid: each piece as the line end, middle, end.
    grid = collections.defaultdict(list)
    pieces = []
    for line in open(a.streets, encoding='utf8'):
        r = json.loads(line)
        ends = r.get('ends') or []
        pts = ([ends[0]] if ends else []) + [r['at']] + ([ends[1]] if len(ends) > 1 else [])
        i = len(pieces)
        pieces.append((r['id'], fold(r['name']), pts, r['name']))
        lats = [p[0] for p in pts]
        lons = [p[1] for p in pts]
        pad = 0.005
        for gy in range(int(math.floor((min(lats) - pad) / CELL)), int(math.floor((max(lats) + pad) / CELL)) + 1):
            for gx in range(int(math.floor((min(lons) - pad) / CELL)), int(math.floor((max(lons) + pad) / CELL)) + 1):
                grid[(gy, gx)].append(i)
    print(f'{len(pieces)} street pieces', file=sys.stderr)
    base = lambda pid: re.sub(r'-\d+$', '', pid)
    stats['streets'] = len({base(p[0]) for p in pieces})

    def near(lat, lon):
        """Each piece near a point with its distance in metres, nearest first."""
        c = math.cos(math.radians(lat))
        out = []
        for i in grid.get((int(math.floor(lat / CELL)), int(math.floor(lon / CELL))), ()):
            pts = pieces[i][2]
            d = min((seg_m(lon, lat, p[1], p[0], q[1], q[0], c) for p, q in zip(pts, pts[1:])), default=None)
            if d is None:
                d = seg_m(lon, lat, pts[0][1], pts[0][0], pts[0][1], pts[0][0], c)
            out.append((d, i))
        out.sort()
        return out

    houses = collections.defaultdict(dict)  # street id -> number -> (lat, lon, piece id)
    fp = osmium.FileProcessor(a.pbf).with_locations().with_filter(osmium.filter.KeyFilter('addr:housenumber'))
    for o in fp:
        if o.is_node():
            if not o.location.valid():
                continue
            lat, lon = o.location.lat, o.location.lon
        elif o.is_way():
            ll = [(n.location.lat, n.location.lon) for n in o.nodes if n.location.valid()]
            if not ll:
                continue
            lat, lon = sum(p[0] for p in ll) / len(ll), sum(p[1] for p in ll) / len(ll)
        else:
            stats['dropped: a relation'] += 1
            continue
        stats['tagged'] += 1
        nums = numbers(o.tags.get('addr:housenumber') or '')
        if not nums:
            stats['dropped: not a plain number'] += 1
            continue
        street = o.tags.get('addr:street')
        cands = near(lat, lon)
        pick = None
        if street:
            f = fold(street)
            same = [(d, i) for d, i in cands if pieces[i][1] == f and d <= WITH_STREET_M]
            if same:
                pick = same[0][1]
                stats['placed: by addr:street'] += 1
            else:
                stats['dropped: no street of that name within 400 m'] += 1
        else:
            if cands and cands[0][0] <= ALONE_M:
                d0, i0 = cands[0]
                other = next((d for d, i in cands[1:] if pieces[i][1] != pieces[i0][1]), None)
                if other is None or other >= 2 * d0 + 15:
                    pick = i0
                    stats['placed: no addr:street, the one street within 40 m'] += 1
                else:
                    stats['dropped: no addr:street, two streets near'] += 1
            else:
                stats['dropped: no addr:street, no street within 40 m'] += 1
        if pick is None:
            continue
        pid = pieces[pick][0]
        h = houses[base(pid)]
        for n in nums:
            # The same number again on the same piece of a street is a duplicate (a building and its entrance); on
            # another piece, about a kilometre on, it is another house (a long road's numbers start again by suburb).
            if (n, pid) in h:
                stats['same number twice on one piece of a street (first kept)'] += 1
            else:
                h[(n, pid)] = (round(lat, 6), round(lon, 6))

    names = {}
    for pid, f, pts, name in pieces:
        names.setdefault(base(pid), name)
    key = lambda k: (int(re.match(r'\d+', k[0]).group()), k[0], k[1])
    with open(f'{a.out}/osm-houses.ndjson', 'w', encoding='utf8') as out:
        for sid in sorted(houses):
            h = houses[sid]
            out.write(json.dumps({'id': sid, 'name': names[sid], 'nums': [[k[0], h[k][0], h[k][1], k[1][len(sid):].lstrip('-')] for k in sorted(h, key=key)]}, ensure_ascii=False) + '\n')
    stats['streets with a number'] = len(houses)
    stats['numbers kept'] = sum(len(h) for h in houses.values())
    stats['streets with two numbers or more'] = sum(1 for h in houses.values() if len({k[0] for k in h}) >= 2)
    json.dump({'counts': dict(stats)}, open(f'{a.out}/osm-houses-stats.json', 'w'), indent=1)
    print(json.dumps(dict(stats), indent=1), file=sys.stderr)


if __name__ == '__main__':
    main()
