# Mako v04: the hull's sections (exec'd by build_mako.py, and by mk_collision.py on its own, after lib2.py, tender_lib.py and
# mk_measure.py). Pure functions, no scene changes.
#
# Each station's section is a PROFILE: the port half, from the keel (lat 0) up to the crest (lat 0), put together from the
# parts the sheet draws (the v04 tables in mk_measure.py), with named FEATURES (points on it) that pin the grid's columns
# (mk_hull.py). Along the hull:
#   * nose (x < 120): the glass bubble's ellipse, whole;
#   * canopy (x 120-286): the jaw below (a flat bottom on its sloped front edge, flanks after FRONT's jaw outline), its dark top
#     rim in to the bubble at the sill, the bubble above (an ellipse: FRONT 1.47 x 1.16 m at the collar), with the lips raised
#     on it (x 205-286) and the hood over its crest (x 197-286);
#   * bezel (x 286-340): the grey sill band, then FRONT's bezel outline (sides, the shoulder over the recess, the peaked top),
#     chamfered 2.5 cm at both faces; the jaw (to x 300) or the lower body below;
#   * body (x 340 on): the lower body up to its widest point (the waist), the dark upper flank (the recess wall at x 340-445)
#     up to the orange shell's edge, the shell (a cape) to the crest.
# The jumps at the bezel's two faces and at the lower body's front face are station pairs 6 mm apart (mk_hull.py).
import math
from mathutils import Vector


def table_f(tab, conv_x=Yi, conv_v=lambda v: v):
    """a measured polyline [(image x, value)] as a function of model y (linear, clamped)"""
    pts = sorted((conv_x(x), conv_v(v)) for x, v in tab)
    def f(y):
        if y <= pts[0][0]:
            return pts[0][1]
        if y >= pts[-1][0]:
            return pts[-1][1]
        for (a, va), (b, vb) in zip(pts, pts[1:]):
            if a <= y <= b:
                return va + (vb - va) * (y - a) / (b - a)
    return f


def _pc(tab, cx=Yi, cv=lambda v: v):
    return pchip(sorted((cx(x), cv(v)) for x, v in tab))


TOPf, BOTf = _pc(TOPY, Yi, Zi), _pc(BOTY, Yi, Zi)                 # the body's crest and keel (x >= 300)
GTf, DWf, DZf = _pc(GTOP, Yi, Zi), _pc(DOMEW, Yi, Lp), _pc(DOMEZC)  # the bubble: crest, half width, centre height
SILLf, JTOPf, JBOTf, JWPf = table_f(SILL, Yi, Zi), table_f(JTOP, Yi, Zi), table_f(JBOT, Yi, Zi), table_f(JW_PLAN)
WJf = pchip(sorted(W_J))                                           # the jaw's half width by height
HWf, HZf = table_f(HOODW), _pc(HOODZ)
LIP0f, LIP1f = table_f([(x, a) for x, a, b in LIPZ]), table_f([(x, b) for x, a, b in LIPZ])
SHWf, SHZf = _pc(SHELL_W), table_f(SHELL_Z)
WLf, WZf = _pc([(x, l) for x, l, z in WAIST]), _pc([(x, z) for x, l, z in WAIST])
GUf = _pc(GLASS_UP)

Y_NOSE, Y_TAIL = Yi(X_TIP), Yi(X_STERN)
Y_JAW, Y_LB = Yi(X_JAW4), Yi(X_LB)
Y_HD0 = Yi(X_HOOD0) if X_HOOD0 is not None else None     # v04 side-first: no hood, no lips
Y_LP0 = Yi(X_LIP0) if X_LIP0 is not None else None
Y_B0, Y_B1 = Yi(X_BZ0), Yi(X_BZ1)
STEP = 0.006                     # station gap of a face (the bezel's front and back faces, the lower body's front face)
CH_B = 0.025                     # the bezel's chamfers: 2.5 cm at 45 deg
LIP_R = 0.015                    # the lips' sloped edges (height on the bubble)
H_SHELL = 0.010                  # the shell's plating over the dark body (mk_hull.py raises it)


def ximg(y):
    return X0 + y / S


def dome(y):
    a = max(1e-4, DWf(y))
    zc = DZf(y)
    return a, zc, max(1e-4, GTf(y) - zc)


def _th_at(zc, b, z):
    return math.asin(max(-1.0, min(1.0, (z - zc) / b)))


class Prof:
    """a profile under construction: points (lat, z) from the keel to the crest and named features (point indices)"""
    def __init__(self):
        self.p, self.f = [], {}

    def add(self, pts):
        for q in pts:
            if self.p and abs(q[0] - self.p[-1][0]) < 1e-7 and abs(q[1] - self.p[-1][1]) < 1e-7:
                continue
            self.p.append((q[0], q[1]))

    def mark(self, name):
        self.f[name] = len(self.p) - 1

    def line(self, a, b, step=0.006):
        n = max(1, int(math.ceil(math.hypot(b[0] - a[0], b[1] - a[1]) / step)))
        self.add([(a[0] + (b[0] - a[0]) * k / n, a[1] + (b[1] - a[1]) * k / n) for k in range(n + 1)])


def _inset_poly(poly, d):
    """a 2D polyline (port half, keel side first) moved d inward (to the left of its direction)"""
    out = []
    n = len(poly)
    for i in range(n):
        a = Vector(poly[max(0, i - 1)]); b = Vector(poly[min(n - 1, i + 1)])
        t = (b - a)
        t = t.normalized() if t.length > 1e-9 else Vector((0, 1))
        nin = Vector((-t.y, t.x))
        q = Vector(poly[i]) + nin * d
        if i == n - 1:
            q.x = 0.0
        out.append((max(0.0, q.x), q.y))
    return out


def _bezel_outline(y):
    """FRONT's bezel outline at y, its top scaled to the hull's crest there; densified, the shoulder's tip rounded"""
    top = TOPf(y)
    pts = []
    for (l, z) in BEZEL:
        if z > 0.75:
            z = 0.75 + (z - 0.75) * (top - 0.75) / (1.10 - 0.75)
        pts.append((l, z))
    return pts


def _dome_upper(P, y, th_from):
    """the bubble from the sill (angle th_from) to its crest, with the lips raised on it and the hood over its crest"""
    a, zc, b = dome(y)
    x = ximg(y)
    brk, lip = [], None
    if X_LIP0 is not None and X_LIP0 - 1e-6 <= x <= X_BZ0 + 1e-6:
        h = LIP_H * smoothstep(X_LIP0, X_BZ0, x)
        z0, z1 = LIP0f(y), LIP1f(y)
        lip = (h, [_th_at(zc, b, z0 - LIP_R), _th_at(zc, b, z0), _th_at(zc, b, z1), _th_at(zc, b, z1 + LIP_R)])
    th_end = math.pi / 2
    hood = None
    if X_HOOD0 is not None and X_HOOD0 - 1e-6 <= x <= X_BZ0 + 1e-6:
        wh = min(HWf(y), 0.97 * a)
        th_end = math.acos(wh / a)
        hood = wh

    def raise_at(th):
        if lip is None:
            return 0.0
        h, (t0, t1, t2, t3) = lip
        if th <= t0 or th >= t3:
            return 0.0
        if th < t1:
            return h * (th - t0) / max(1e-9, t1 - t0)
        if th > t2:
            return h * (t3 - th) / max(1e-9, t3 - t2)
        return h

    def pt(th):
        c, s_ = math.cos(th), math.sin(th)
        nx, nz = c / a, s_ / b
        L = math.hypot(nx, nz)
        r = raise_at(th)
        return (a * c + r * nx / L, zc + b * s_ + r * nz / L)

    marks = []
    zgu = GUf(y)                                               # v04 round 2: the line between the pale upper glass and the clear
    if zc - b < zgu < zc + b:
        t_gu = _th_at(zc, b, zgu)
        if th_from < t_gu < th_end:
            marks.append((t_gu, "GU"))
    if lip is not None:
        for nm, t in zip(("L0a", "L0b", "L1a", "L1b"), lip[1]):
            if th_from < t < th_end:
                marks.append((t, nm))
    marks.sort(key=lambda m: m[0])
    marks.append((th_end, "H" if hood is not None else None))
    t_prev = th_from
    P.add([pt(th_from)])
    for t_m, nm in marks:
        arc = abs(t_m - t_prev) * max(a, b)
        n = max(1, int(math.ceil(arc / 0.006)))
        P.add([pt(t_prev + (t_m - t_prev) * k / n) for k in range(1, n + 1)])
        if nm:
            P.mark(nm)
        t_prev = t_m
    if hood is not None:                       # the hood's roof: steep sides off the bubble, a flat top
        le, ze = P.p[-1]
        zh = max(HZf(y), ze + 0.002)
        n = 40
        P.add([(le * (1 - u), ze + (zh - ze) * (1 - (1 - u) ** 6) ** (1 / 6)) for u in [k / n for k in range(1, n + 1)]])


def _jaw_lower(P, y):
    """the jaw from the keel to its top outer corner (Q): a flat bottom on the front edge, rounded corner, FRONT's flanks"""
    a, zc, b = dome(y)
    zt = JTOPf(y)
    zb = min(JBOTf(y), zt)
    wt = max(0.01, min(WJf(zt), JWPf(y), 0.93 * a))
    h = zt - zb
    if h < 0.02:
        P.line((0.0, zb), (wt, zt))
    else:
        r = min(0.045, h / 3)
        wb = min(WJf(zb + r), wt)
        P.line((0.0, zb), (max(0.0, wb - r), zb))
        P.add([(wb - r + r * math.sin(math.pi / 2 * k / 8), zb + r - r * math.cos(math.pi / 2 * k / 8)) for k in range(1, 9)])
        n = max(2, int(math.ceil((zt - zb - r) / 0.008)))
        P.add([(min(WJf(z), wt), z) for z in [zb + r + (zt - zb - r) * k / n for k in range(1, n + 1)]])
    P.mark("Q")


def nlb(x):
    return 3.0 - 0.7 * smoothstep(450, 884, x)


def _body_lower(P, y):
    """the lower body from the keel to its widest point (Q): a superellipse quarter, boxy under the cockpit"""
    wq, zq, zb = WLf(y), WZf(y), BOTf(y)
    n = nlb(ximg(y))
    N = 40
    for k in range(N + 1):
        t = math.pi / 2 * (1 - k / N)
        P.add([(wq * math.cos(t) ** (2 / n), zq - (zq - zb) * math.sin(t) ** (2 / n))])
    P.mark("Q")


def ncape(x):
    return 2.2 - 0.4 * smoothstep(380, 620, x)


def profile(y, case=None, lower=None, inset=None):
    """the section at y: Prof with points from the keel to the crest and its features. case ('canopy', 'bezel', 'body') and
    lower ('jaw', 'lowerbody') force which parts it is made of, for the stations the collar's arcs warp along the hull
    (mk_hull.py); by default y decides"""
    x = ximg(y)
    P = Prof()
    if case == "nose" or (case is None and y < Y_JAW):         # the nose: the bubble's ellipse, whole
        a, zc, b = dome(y)
        N = 90
        t_gu = _th_at(zc, b, GUf(y))
        ts = sorted(set([-math.pi / 2 + math.pi * k / N for k in range(N + 1)] + [t_gu]))
        for t in ts:
            P.add([(a * math.cos(t), zc + b * math.sin(t))])
            if abs(t - t_gu) < 1e-12 and -math.pi / 2 < t < math.pi / 2:
                P.mark("GU")
        P.f["crest"] = len(P.p) - 1
        return P
    if case is None:
        case = "canopy" if y < Y_B0 + STEP / 2 else ("bezel" if y < Y_B1 - STEP / 2 else "body")
    if lower is None:
        lower = "jaw" if y < Y_LB + STEP / 2 else "lowerbody"
    in_bezel = case == "bezel"
    if lower == "jaw":
        _jaw_lower(P, y)
    else:
        _body_lower(P, y)
    q = P.p[-1]
    if case == "canopy":                                        # canopy: the jaw's top rim in to the bubble, the bubble
        a, zc, b = dome(y)
        zs = max(SILLf(y), zc - b + 0.002)
        th_s = _th_at(zc, b, zs)
        P.line(q, (a * math.cos(th_s), zc + b * math.sin(th_s)))
        P.mark("S")
        _dome_upper(P, y, th_s)
    elif in_bezel:                                              # the bezel: grey sill band, then its outline (chamfered)
        if inset is None:                                       # 2.5 cm chamfers at both faces
            d = (y - Y_B0) if y < 0.5 * (Y_B0 + Y_B1) else (Y_B1 - y)
            inset = CH_B * max(0.0, min(1.0, 1.0 - (d - STEP) / (CH_B - STEP)))
        out = _bezel_outline(y)
        if inset > 1e-6:
            out = _inset_poly(out, inset)
        P.line(q, out[0])
        P.mark("S")
        for i in range(1, len(out)):
            P.line(out[i - 1], out[i], 0.008)
            if i == 5:
                P.mark("SU")
            if i == 8:
                P.mark("SH")
    else:                                                       # the body: dark upper flank, the shell's edge, the cape
        wq, zq = q
        we, ze = min(SHWf(y), wq - 0.004), SHZf(y)
        N = 36
        P.add([(we + (wq - we) * math.cos(t), zq + (ze - zq) * math.sin(t)) for t in [math.pi / 2 * k / N for k in range(1, N + 1)]])
        P.mark("E")
        zt = TOPf(y)
        n = ncape(x)
        N = 48
        P.add([(we * math.cos(t) ** (2 / n), ze + (zt - ze) * math.sin(t) ** (2 / n)) for t in [math.pi / 2 * k / N for k in range(1, N + 1)]])
    P.f["crest"] = len(P.p) - 1
    if P.p[-1][0] != 0.0:
        P.add([(0.0, P.p[-1][1])])
        P.f["crest"] = len(P.p) - 1
    return P


def arc_of(P):
    arc = [0.0]
    for i in range(1, len(P.p)):
        arc.append(arc[-1] + math.hypot(P.p[i][0] - P.p[i - 1][0], P.p[i][1] - P.p[i - 1][1]))
    return arc


def at_arc(P, arc, s):
    """the profile's point at arc length s (linear between its points)"""
    if s <= 0:
        return P.p[0]
    if s >= arc[-1]:
        return P.p[-1]
    lo, hi = 0, len(arc) - 1
    while hi - lo > 1:
        m = (lo + hi) // 2
        if arc[m] <= s:
            lo = m
        else:
            hi = m
    t = (s - arc[lo]) / max(1e-12, arc[hi] - arc[lo])
    (a0, b0), (a1, b1) = P.p[lo], P.p[hi]
    return (a0 + (a1 - a0) * t, b0 + (b1 - b0) * t)


def section_pts(y, n=16):
    """n points round the whole section at y (port and starboard) with their outward 2D normals, for collision hulls"""
    P = profile(y)
    arc = arc_of(P)
    L = arc[-1]
    out = []
    half = n // 2
    for k in range(half + 1):
        s = L * k / half
        p = at_arc(P, arc, s)
        q = at_arc(P, arc, min(L, s + 0.01)); o = at_arc(P, arc, max(0.0, s - 0.01))
        t = (q[0] - o[0], q[1] - o[1])
        nn = Vector((t[1], -t[0]))
        nn = nn.normalized() if nn.length > 1e-9 else Vector((1, 0))
        out.append((Vector((p[0], y, p[1])), Vector((nn.x, 0, nn.y))))
    full = out + [(Vector((-p.x, p.y, p.z)), Vector((-m.x, 0, m.z))) for p, m in reversed(out[1:-1])]
    return full


def side_poly(img):
    return [(Yi(x), Zi(v)) for x, v in img]


def pip(pt, poly):
    x, y = pt; inside = False; n = len(poly)
    for i in range(n):
        x1, y1 = poly[i]; x2, y2 = poly[(i + 1) % n]
        if (y1 > y) != (y2 > y) and x < (x2 - x1) * (y - y1) / (y2 - y1) + x1:
            inside = not inside
    return inside


P_BAY = side_poly(BAY)