"""HNX-02 "Copernicus-class" crewed interplanetary transfer vehicle.

Rebuild of the crewed starship with procedural engineering structures
(space trusses, pressure vessels, habitat modules, a spin-gravity
centrifuge, an NTP engine cluster, and a docked Orion-class shuttle as
the scale reference) instead of stacked primitives.

    --step-id <id> --output <glb> --result-output <json> [--render-output <png>]

Design intent (display-scale visual research, NOT closed engineering):
    - Nuclear-thermal propulsion (NERVA-class stainless/epoxy composite
      fuel element heritage) with a shadow shield between reactor and crew.
    - Axial layout keeps the crew section as far from the reactor as the
      truss allows; the water-wall tanks double as radiation shielding.
    - The centrifuge ring provides spin gravity; cars are the pressurized
      crew volumes. No mass/thrust/shielding closure is claimed.

Object ids: hnx02.<part>[.<index>], lowercase dotted, stable per step.
Cumulative object lists must stay in lockstep with Design.shipSteps in
src/SomhairlesDream.Modeling/Definitions.fs.
"""
import argparse
import json
import math
import os
import random
import sys

import bpy  # noqa: F401  (must be imported first: registers bundled bmesh/mathutils)
import bmesh
from mathutils import Vector

# --------------------------------------------------------------- contract

STEP_OBJECTS = {
    "ship-01-keel": [
        "hnx02.hull.core",
        "hnx02.spine.truss",
        "hnx02.node.forward",
    ],
    "ship-02-hab": [
        "hnx02.hull.core",
        "hnx02.spine.truss",
        "hnx02.node.forward",
        "hnx02.hab.a",
        "hnx02.hab.b",
        "hnx02.hab.whipple",
    ],
    "ship-03-radiation": [
        "hnx02.hull.core",
        "hnx02.spine.truss",
        "hnx02.node.forward",
        "hnx02.hab.a",
        "hnx02.hab.b",
        "hnx02.hab.whipple",
        "hnx02.rad.waterwall",
        "hnx02.centrifuge.ring",
    ],
    "ship-04-reactor": [
        "hnx02.hull.core",
        "hnx02.spine.truss",
        "hnx02.node.forward",
        "hnx02.hab.a",
        "hnx02.hab.b",
        "hnx02.hab.whipple",
        "hnx02.rad.waterwall",
        "hnx02.centrifuge.ring",
        "hnx02.prop.tank",
        "hnx02.reactor.core",
        "hnx02.thrust.cone",
        "hnx02.nozzle.cluster",
    ],
    "ship-05-power": [
        "hnx02.hull.core",
        "hnx02.spine.truss",
        "hnx02.node.forward",
        "hnx02.hab.a",
        "hnx02.hab.b",
        "hnx02.hab.whipple",
        "hnx02.rad.waterwall",
        "hnx02.centrifuge.ring",
        "hnx02.prop.tank",
        "hnx02.reactor.core",
        "hnx02.thrust.cone",
        "hnx02.nozzle.cluster",
        "hnx02.power.solar",
        "hnx02.power.radiator",
    ],
    "ship-06-cabin": [
        "hnx02.hull.core",
        "hnx02.spine.truss",
        "hnx02.node.forward",
        "hnx02.hab.a",
        "hnx02.hab.b",
        "hnx02.hab.whipple",
        "hnx02.rad.waterwall",
        "hnx02.centrifuge.ring",
        "hnx02.prop.tank",
        "hnx02.reactor.core",
        "hnx02.thrust.cone",
        "hnx02.nozzle.cluster",
        "hnx02.power.solar",
        "hnx02.power.radiator",
        "hnx02.cabin.airlock",
        "hnx02.shuttle.orion",
    ],
    "ship-07-livery": [
        "hnx02.hull.core",
        "hnx02.spine.truss",
        "hnx02.node.forward",
        "hnx02.hab.a",
        "hnx02.hab.b",
        "hnx02.hab.whipple",
        "hnx02.rad.waterwall",
        "hnx02.centrifuge.ring",
        "hnx02.prop.tank",
        "hnx02.reactor.core",
        "hnx02.thrust.cone",
        "hnx02.nozzle.cluster",
        "hnx02.power.solar",
        "hnx02.power.radiator",
        "hnx02.cabin.airlock",
        "hnx02.shuttle.orion",
        "hnx02.detail.greeble",
        "hnx02.detail.antenna",
        "hnx02.detail.glow",
        "hnx02.detail.glass",
        "hnx02.livery.stripes",
        "hnx02.livery.insignia",
    ],
}

# Material keys per object id.
OBJECT_MATERIAL = {
    "hnx02.hull.core": "hull",
    "hnx02.spine.truss": "hull",
    "hnx02.node.forward": "white",
    "hnx02.hab.a": "white",
    "hnx02.hab.b": "white",
    "hnx02.hab.whipple": "shield",
    "hnx02.rad.waterwall": "hull",
    "hnx02.centrifuge.ring": "hull",
    "hnx02.reactor.core": "dark",
    "hnx02.prop.tank": "gold",
    "hnx02.thrust.cone": "hull",
    "hnx02.nozzle.cluster": "dark",
    "hnx02.power.solar": "pv",
    "hnx02.power.radiator": "radiator",
    "hnx02.cabin.airlock": "white",
    "hnx02.shuttle.orion": "white",
    "hnx02.detail.greeble": "hull",
    "hnx02.detail.antenna": "hull",
    "hnx02.detail.glow": "glow",
    "hnx02.detail.glass": "glass",
    "hnx02.livery.stripes": "gold",
    "hnx02.livery.insignia": "insignia",
}


def utc_now():
    from datetime import datetime, timezone

    return datetime.now(timezone.utc).strftime("%Y-%m-%dT%H:%M:%SZ")


def fail(message, code):
    print(json.dumps({"status": "failed", "message": message}))
    sys.exit(code)


# --------------------------------------------------------------- materials


def make_material(bpy, name, color, *, metallic=0.0, roughness=0.45, emission=None, emission_strength=0.0):
    material = bpy.data.materials.get(name)
    if material is not None:
        return material
    material = bpy.data.materials.new(name)
    material.use_nodes = True
    bsdf = material.node_tree.nodes.get("Principled BSDF")
    bsdf.inputs["Base Color"].default_value = (*color, 1.0)
    bsdf.inputs["Metallic"].default_value = metallic
    bsdf.inputs["Roughness"].default_value = roughness
    if emission is not None:
        bsdf.inputs["Emission Color"].default_value = (*emission, 1.0)
        bsdf.inputs["Emission Strength"].default_value = emission_strength
    return material


def build_materials(bpy):
    return {
        # PLA Space Force livery: "deep-space blue" hull (深空蓝) with gold accents.
        "hull": make_material(bpy, "hnx02.hull", (0.035, 0.075, 0.260), metallic=0.20, roughness=0.46),
        "white": make_material(bpy, "hnx02.white", (0.80, 0.82, 0.85), metallic=0.18, roughness=0.52),
        "shield": make_material(bpy, "hnx02.shield", (0.44, 0.46, 0.50), metallic=0.30, roughness=0.72),
        "dark": make_material(bpy, "hnx02.dark", (0.10, 0.11, 0.13), metallic=0.72, roughness=0.48),
        "radiator": make_material(bpy, "hnx02.radiator", (0.90, 0.91, 0.92), metallic=0.08, roughness=0.30),
        "pv": make_material(bpy, "hnx02.pv", (0.07, 0.10, 0.24), metallic=0.35, roughness=0.28),
        "gold": make_material(bpy, "hnx02.gold", (0.88, 0.64, 0.14), metallic=0.45, roughness=0.32),
        "insignia": make_material(bpy, "hnx02.insignia", (0.72, 0.09, 0.09), metallic=0.15, roughness=0.45),
        "copper": make_material(bpy, "hnx02.copper", (0.70, 0.44, 0.19), metallic=1.0, roughness=0.34),
        "glow": make_material(bpy, "hnx02.glow", (1.0, 0.35, 0.10), metallic=0.0, roughness=0.6,
                              emission=(1.0, 0.42, 0.14), emission_strength=9.0),
        "glass": make_material(bpy, "hnx02.glass", (0.95, 0.84, 0.55), metallic=0.0, roughness=0.25,
                               emission=(0.98, 0.86, 0.52), emission_strength=5.0),
    }


# ------------------------------------------------------------ bmesh kit


def basis(axis):
    z = Vector(axis).normalized()
    ref = Vector((0.0, 0.0, 1.0)) if abs(z.z) < 0.9 else Vector((1.0, 0.0, 0.0))
    x = ref.cross(z).normalized()
    y = z.cross(x).normalized()
    return x, y, z


def cyl(bm, p1, p2, r1, r2=None, seg=12, cap1=True, cap2=True):
    p1, p2 = Vector(p1), Vector(p2)
    r2 = r1 if r2 is None else r2
    d = p2 - p1
    if d.length < 1e-7:
        return
    x, y, _ = basis(d)
    ring1, ring2 = [], []
    for i in range(seg):
        a = 2.0 * math.pi * i / seg
        off = x * math.cos(a) + y * math.sin(a)
        ring1.append(bm.verts.new(p1 + off * r1))
        ring2.append(bm.verts.new(p2 + off * r2))
    for i in range(seg):
        j = (i + 1) % seg
        bm.faces.new((ring1[i], ring1[j], ring2[j], ring2[i]))
    if cap1:
        bm.faces.new(tuple(reversed(ring1)))
    if cap2:
        bm.faces.new(tuple(ring2))


def torus(bm, center, R, r, axis=(0.0, 0.0, 1.0), seg=24, rseg=8):
    center = Vector(center)
    x, y, z = basis(axis)
    rings = []
    for i in range(seg):
        a = 2.0 * math.pi * i / seg
        radial = x * math.cos(a) + y * math.sin(a)
        c = center + radial * R
        ring = []
        for j in range(rseg):
            b = 2.0 * math.pi * j / rseg
            ring.append(bm.verts.new(c + radial * (r * math.cos(b)) + z * (r * math.sin(b))))
        rings.append(ring)
    for i in range(seg):
        i2 = (i + 1) % seg
        for j in range(rseg):
            j2 = (j + 1) % rseg
            bm.faces.new((rings[i][j], rings[i][j2], rings[i2][j2], rings[i2][j]))


def dome(bm, base, direction, radius, seg=20, layers=5):
    d = Vector(direction).normalized()
    x, y, _ = basis(d)
    prev = None
    for i in range(layers + 1):
        t = i / layers
        phi = t * math.pi / 2.0
        r = radius * math.cos(phi)
        h = radius * math.sin(phi)
        c = Vector(base) + d * h
        if r < 1e-6:
            tip = bm.verts.new(c)
            if prev:
                for j in range(seg):
                    bm.faces.new((prev[j], prev[(j + 1) % seg], tip))
            break
        ring = []
        for j in range(seg):
            a = 2.0 * math.pi * j / seg
            ring.append(bm.verts.new(c + (x * math.cos(a) + y * math.sin(a)) * r))
        if prev:
            for j in range(seg):
                j2 = (j + 1) % seg
                bm.faces.new((prev[j], prev[j2], ring[j2], ring[j]))
        prev = ring


def bell(bm, throat, direction, length, r_throat, r_exit, seg=32, layers=14, exponent=1.7):
    d = Vector(direction).normalized()
    x, y, _ = basis(d)
    prev = None
    for i in range(layers + 1):
        t = i / layers
        r = r_throat + (r_exit - r_throat) * (t ** exponent)
        c = Vector(throat) + d * (length * t)
        ring = []
        for j in range(seg):
            a = 2.0 * math.pi * j / seg
            ring.append(bm.verts.new(c + (x * math.cos(a) + y * math.sin(a)) * r))
        if prev:
            for j in range(seg):
                j2 = (j + 1) % seg
                bm.faces.new((prev[j], prev[j2], ring[j2], ring[j]))
        prev = ring


def box(bm, center, ex, ey, ez, ax, ay, az):
    c = Vector(center)
    ax, ay, az = Vector(ax), Vector(ay), Vector(az)
    vs = {}
    for sx in (-1, 1):
        for sy in (-1, 1):
            for sz in (-1, 1):
                vs[(sx, sy, sz)] = bm.verts.new(c + ax * (sx * ex) + ay * (sy * ey) + az * (sz * ez))
    quads = [
        (vs[(-1, -1, -1)], vs[(-1, 1, -1)], vs[(-1, 1, 1)], vs[(-1, -1, 1)]),
        (vs[(1, 1, -1)], vs[(1, -1, -1)], vs[(1, -1, 1)], vs[(1, 1, 1)]),
        (vs[(-1, 1, -1)], vs[(1, 1, -1)], vs[(1, 1, 1)], vs[(-1, 1, 1)]),
        (vs[(1, -1, -1)], vs[(-1, -1, -1)], vs[(-1, -1, 1)], vs[(1, -1, 1)]),
        (vs[(-1, -1, 1)], vs[(-1, 1, 1)], vs[(1, 1, 1)], vs[(1, -1, 1)]),
        (vs[(1, -1, -1)], vs[(1, 1, -1)], vs[(-1, 1, -1)], vs[(-1, -1, -1)]),
    ]
    for q in quads:
        bm.faces.new(q)


# ------------------------------------------------------------ structures


def truss(bm, p1, p2, width, bays=6, chord_r=0.09, strut_r=0.05, angles=None):
    p1, p2 = Vector(p1), Vector(p2)
    if angles is None:
        angles = [90.0, 210.0, 330.0]
    x, y, _ = basis(p2 - p1)
    nodes = []
    for deg in angles:
        a = math.radians(deg)
        off = (x * math.cos(a) + y * math.sin(a)) * width
        nodes.append((p1 + off, p2 + off))
    for a, b in nodes:
        cyl(bm, a, b, chord_r, seg=8)
    n = len(nodes)
    for i in range(bays + 1):
        t = i / bays
        pts = [a.lerp(b, t) for a, b in nodes]
        for k in range(n):
            cyl(bm, pts[k], pts[(k + 1) % n], strut_r, seg=6)
        if i < bays:
            t2 = (i + 1) / bays
            pts2 = [a.lerp(b, t2) for a, b in nodes]
            for k in range(n):
                cyl(bm, pts[k], pts2[(k + 1) % n], strut_r, seg=6)


def ring_frame(bm, center, radius, axis, tube=0.07, seg=28, spokes=0, spoke_r=0.05):
    center = Vector(center)
    torus(bm, center, radius, tube, axis=axis, seg=seg, rseg=6)
    if spokes:
        x, y, _ = basis(axis)
        for k in range(spokes):
            a = 2.0 * math.pi * k / spokes
            radial = x * math.cos(a) + y * math.sin(a)
            cyl(bm, center + radial * radius * 0.12, center + radial * radius, spoke_r, seg=6)


def pressure_vessel(bm, p1, p2, radius, rings=5, seg=24, ports=True):
    p1, p2 = Vector(p1), Vector(p2)
    d = (p2 - p1).normalized()
    cyl(bm, p1, p2, radius, seg=seg, cap1=False, cap2=False)
    dome(bm, p1, -d, radius, seg=seg)
    dome(bm, p2, d, radius, seg=seg)
    for i in range(1, rings + 1):
        t = i / (rings + 1)
        torus(bm, p1.lerp(p2, t), radius, 0.06, axis=d, seg=seg, rseg=6)
    if ports:
        x, y, _ = basis(d)
        for t in (0.2, 0.5, 0.8):
            c = p1.lerp(p2, t)
            off = x * radius
            cyl(bm, c + off * 0.88, c + off * 1.45, 0.10, seg=8)
            torus(bm, c + off * 1.45, 0.13, 0.04, axis=x, seg=10, rseg=5)


def hab_module(bm, p1, p2, radius, ribs=8, windows=10, seg=28, hatch=True):
    p1, p2 = Vector(p1), Vector(p2)
    d = (p2 - p1).normalized()
    x, y, _ = basis(d)
    cyl(bm, p1, p2, radius, seg=seg)
    for i in range(1, ribs + 1):
        t = i / (ribs + 1)
        torus(bm, p1.lerp(p2, t), radius, 0.075, axis=d, seg=seg, rseg=6)
    length = (p2 - p1).length
    for row, tt in ((0, 0.34), (1, 0.66)):
        for k in range(windows):
            a = 2.0 * math.pi * k / windows + row * math.pi / windows
            radial = x * math.cos(a) + y * math.sin(a)
            c = p1.lerp(p2, tt) + radial * radius
            cyl(bm, c - radial * 0.08, c + radial * 0.14, 0.17, seg=12, cap1=False)
            torus(bm, c + radial * 0.14, 0.185, 0.04, axis=radial, seg=14, rseg=5)
            box(bm, c + radial * 0.16, 0.10, 0.30, 0.012, radial, d, d.cross(radial))
    if hatch:
        for end, dirv in ((p2, d), (p1, -d)):
            torus(bm, end + dirv * 0.05, radius * 0.52, 0.11, axis=d, seg=seg, rseg=6)
            for k in range(8):
                a = 2.0 * math.pi * k / 8
                radial = x * math.cos(a) + y * math.sin(a)
                cyl(bm, end + radial * radius * 0.52 + dirv * 0.05,
                    end + radial * radius * 0.70 + dirv * 0.05, 0.05, seg=6)


def centrifuge_ring(bm, center, radius, axis=(0.0, 1.0, 0.0), cars=6, seg=44):
    """Rotating gravity ring: twin truss hoops, spokes, suspended cars."""
    center = Vector(center)
    x, y, z = basis(axis)
    for off in (-0.55, 0.55):
        c = center + z * off
        torus(bm, c, radius, 0.085, axis=z, seg=seg, rseg=6)
        torus(bm, c, radius * 0.995, 0.045, axis=z, seg=seg, rseg=5)
    # cross bracing between the two hoops
    for k in range(0, seg, 2):
        a = 2.0 * math.pi * k / seg
        a2 = 2.0 * math.pi * (k + 1) / seg
        radial = x * math.cos(a) + y * math.sin(a)
        radial2 = x * math.cos(a2) + y * math.sin(a2)
        p_a0 = center + radial * radius + z * -0.55
        p_a1 = center + radial * radius + z * 0.55
        p_b0 = center + radial2 * radius + z * -0.55
        p_b1 = center + radial2 * radius + z * 0.55
        cyl(bm, p_a0, p_b1, 0.035, seg=5)
        cyl(bm, p_a1, p_b0, 0.035, seg=5)
        cyl(bm, p_a0, p_a1, 0.035, seg=5)
    # spokes to hub
    for k in range(cars):
        a = 2.0 * math.pi * k / cars
        radial = x * math.cos(a) + y * math.sin(a)
        outer = center + radial * radius
        cyl(bm, center, outer, 0.17, seg=10)
        for tang in (-1.0, 1.0):
            tangent = z * tang
            p0 = center + radial * radius * 0.18 + tangent * 0.7
            p1 = center + radial * radius * 0.55 + tangent * 0.7
            cyl(bm, p0, p1, 0.07, seg=6)
        cyl(bm, center + radial * radius * 0.5 - z * 0.9, center + radial * radius * 0.5 + z * 0.9, 0.09, seg=6)
    # de-spin bearing rings at the hub
    for zz in (-0.9, 0.9):
        torus(bm, center + z * zz, radius * 0.135, 0.11, axis=z, seg=28, rseg=6)
    # suspended cars
    for k in range(cars):
        a = 2.0 * math.pi * k / cars + math.pi / cars
        radial = x * math.cos(a) + y * math.sin(a)
        car = center + radial * (radius + 0.62)
        cyl(bm, car - z * 0.85, car + z * 0.85, 0.62, seg=18)
        torus(bm, car - z * 0.85, 0.64, 0.07, axis=z, seg=18, rseg=5)
        torus(bm, car + z * 0.85, 0.64, 0.07, axis=z, seg=18, rseg=5)
        cyl(bm, center + radial * radius, car, 0.075, seg=7)
        for kk in range(5):
            aa = 2.0 * math.pi * kk / 5
            off = (z * math.cos(aa) + radial * math.sin(aa)) * 0.2
            cyl(bm, car + off, car + off * 2.6, 0.03, seg=5)


def gimbal(bm, pos, direction, radius=0.28, s=1.0):
    pos = Vector(pos)
    d = Vector(direction).normalized()
    x, _, _ = basis(d)
    torus(bm, pos, radius * s, 0.055 * s, axis=d, seg=16, rseg=5)
    torus(bm, pos + d * 0.15 * s, radius * s * 0.9, 0.055 * s, axis=x, seg=16, rseg=5)
    for k in (0, 1):
        side = 1.0 if k == 0 else -1.0
        cyl(bm, pos + x * side * radius * s * 1.35 - d * 0.12 * s,
            pos + x * side * radius * s * 0.8 + d * 0.3 * s, 0.04 * s, seg=6)


def regen_channels(bm, throat, direction, length, r_throat, r_exit, turns=6, tube=0.035, seg_pts=15):
    throat = Vector(throat)
    d = Vector(direction).normalized()
    x, y, _ = basis(d)
    steps = turns * seg_pts
    pts = []
    for i in range(steps + 1):
        t = i / steps
        r = r_throat + (r_exit - r_throat) * (t ** 1.7) + tube * 1.3
        a = 2.0 * math.pi * turns * t
        c = throat + d * (length * t)
        pts.append(c + (x * math.cos(a) + y * math.sin(a)) * r)
    for a, b in zip(pts, pts[1:]):
        cyl(bm, a, b, tube, seg=5)


def ntp_engine(bm, mount, direction, scale=1.0):
    """Nuclear thermal engine: thrust structure, chamber, bell, turbopump."""
    mount = Vector(mount)
    d = Vector(direction).normalized()
    s = scale
    x, y, _ = basis(d)
    cyl(bm, mount, mount + d * 0.5 * s, 0.44 * s, seg=14)
    chamber = mount + d * 0.5 * s
    cyl(bm, chamber, chamber + d * 0.45 * s, 0.36 * s, 0.31 * s, seg=18)
    bell(bm, chamber + d * 0.45 * s, d, 2.3 * s, 0.31 * s, 1.05 * s, seg=32, layers=14)
    for t in (0.3, 0.55, 0.78, 0.95):
        c = chamber + d * (0.45 * s + 2.3 * s * t)
        r = 0.31 * s + (1.05 * s - 0.31 * s) * (t ** 1.7)
        torus(bm, c, r + 0.02 * s, 0.032 * s, axis=d, seg=24, rseg=5)
    regen_channels(bm, chamber + d * 0.45 * s, d, 2.3 * s, 0.31 * s, 1.05 * s)
    # turbopump cluster
    pump = mount + d * 0.75 * s + x * 0.5 * s
    cyl(bm, pump, pump + d * 0.6 * s, 0.21 * s, seg=12)
    cyl(bm, pump - x * 0.22 * s, pump + x * 0.6 * s, 0.10 * s, seg=8)
    cyl(bm, chamber + x * 0.32 * s, pump + d * 0.32 * s + x * 0.27 * s, 0.075 * s, seg=7)
    cyl(bm, mount + d * 0.18 * s + x * 0.48 * s, chamber + x * 0.32 * s, 0.055 * s, seg=6)
    # nozzle exit stiffener ring
    exit_c = chamber + d * (0.45 * s + 2.3 * s)
    torus(bm, exit_c - d * 0.05 * s, 1.05 * s, 0.05 * s, axis=d, seg=32, rseg=6)


def shadow_shield(bm, center, direction, radius, seg=32, ribs=10):
    center = Vector(center)
    d = Vector(direction).normalized()
    x, y, _ = basis(d)
    cyl(bm, center - d * 0.12, center + d * 0.12, radius, seg=seg)
    for k in range(ribs):
        a = 2.0 * math.pi * k / ribs
        radial = x * math.cos(a) + y * math.sin(a)
        cyl(bm, center + radial * radius * 0.12 - d * 0.12,
            center + radial * radius * 0.96 - d * 0.12, 0.05, seg=6)
    torus(bm, center, radius, 0.06, axis=d, seg=seg, rseg=6)
    for k in range(6):
        a = 2.0 * math.pi * k / 6
        radial = x * math.cos(a) + y * math.sin(a)
        cyl(bm, center + radial * radius * 0.55, center + radial * radius * 0.55 - d * 1.2, 0.08, seg=7)


def thrust_cone(bm, mount, direction, r_base, r_top, length, struts=10):
    mount = Vector(mount)
    d = Vector(direction).normalized()
    x, y, _ = basis(d)
    top = mount + d * length
    for i, (c, r) in enumerate(((mount, r_base), (top, r_top))):
        torus(bm, c, r, 0.06, axis=d, seg=28, rseg=5)
    for k in range(struts):
        a = 2.0 * math.pi * k / struts
        radial = x * math.cos(a) + y * math.sin(a)
        cyl(bm, mount + radial * r_base, top + radial * r_top, 0.055, seg=6)
    for k in range(struts):
        a = 2.0 * math.pi * k / struts
        b = 2.0 * math.pi * (k + 1) / struts
        ra = x * math.cos(a) + y * math.sin(a)
        rb = x * math.cos(b) + y * math.sin(b)
        mid = mount + d * (length * 0.5)
        cyl(bm, mount + ra * r_base, mid + rb * ((r_base + r_top) * 0.5), 0.04, seg=5)


def mli_bands(bm, p1, p2, radius, bands=10, seg=24, wobble=0.035):
    p1, p2 = Vector(p1), Vector(p2)
    d = (p2 - p1).normalized()
    for i in range(1, bands + 1):
        t = i / (bands + 1)
        c = p1.lerp(p2, t)
        w = 1.0 + wobble * math.sin(i * 2.3)
        torus(bm, c, radius * w, 0.085, axis=d, seg=seg, rseg=7)


def solar_array(bm, root, direction, span, width, cells=72, rows=14, face_normal=(0.0, 0.0, 1.0)):
    """Deployable PV blanket: cell grid, frame rails, torque tube, hinge."""
    root = Vector(root)
    d = Vector(direction).normalized()
    x, _y, _z = basis(d)
    fn = Vector(face_normal).normalized()
    if abs(fn.dot(d)) > 0.9:
        fn = _z
    x = d.cross(fn).normalized()
    cw = span / cells
    ch = width / rows
    for i in range(cells):
        for j in range(rows):
            c = root + d * (cw * (i + 0.5)) + x * ((j - rows * 0.5 + 0.5) * ch)
            box(bm, c, cw * 0.40, ch * 0.38, 0.02, d, x, fn)
    for side in (-1.0, 1.0):
        off = x * side * width * 0.5
        box(bm, root + d * (span * 0.5) + off, span * 0.5, 0.05, 0.05, d, x, fn)
    cyl(bm, root - x * width * 0.5, root + x * width * 0.5, 0.10, seg=10)
    cyl(bm, root, root - d * 0.6, 0.13, seg=10)
    torus(bm, root, 0.24, 0.06, axis=x, seg=16, rseg=5)


def fin_radiator(bm, p1, p2, width, fins=48, fin_h=0.42, coolant=True):
    p1, p2 = Vector(p1), Vector(p2)
    d = (p2 - p1).normalized()
    x, y, _ = basis(d)
    cyl(bm, p1, p2, width * 0.5, seg=4)
    for i in range(fins):
        t = (i + 0.5) / fins
        c = p1.lerp(p2, t)
        for side in (0.5, -0.5):
            off = y * side * 0.05
            box(bm, c + off + y * side * fin_h * 0.5, width * 0.46, fin_h * 0.5, 0.014, x, y, d)
    if coolant:
        for side in (-0.30, 0.0, 0.30):
            off = x * (width * side)
            cyl(bm, p1 + off, p2 + off, 0.035, seg=6)


def pipe_run(bm, points, r=0.06, seg=8):
    for a, b in zip(points, points[1:]):
        cyl(bm, a, b, r, seg=seg)


def rcs_cluster(bm, pos, direction, s=1.0):
    pos = Vector(pos)
    d = Vector(direction).normalized()
    x, y, _ = basis(d)
    cyl(bm, pos, pos + d * 0.22 * s, 0.17 * s, seg=12)
    for k in range(4):
        a = 2.0 * math.pi * k / 4 + math.pi / 4
        off = (x * math.cos(a) + y * math.sin(a)) * 0.14 * s
        cyl(bm, pos + d * 0.22 * s + off, pos + d * 0.40 * s + off * 1.35, 0.055 * s, 0.032 * s, seg=6)


def antenna_dish(bm, base, direction, radius, depth=0.4, ribs=8, s=1.0):
    base = Vector(base)
    d = Vector(direction).normalized()
    x, y, _ = basis(d)
    cyl(bm, base - d * 0.1 * s, base + d * 0.28 * s, 0.16 * s, seg=10)
    gimbal(bm, base + d * 0.28 * s, d, radius=0.2 * s, s=s)
    dish_center = base + d * 0.45 * s
    layers = 5
    prev = None
    for i in range(layers + 1):
        t = i / layers
        r = radius * t
        h = depth * (1.0 - t * t)
        ring = []
        for j in range(24):
            a = 2.0 * math.pi * j / 24
            ring.append(bm.verts.new(dish_center + (x * math.cos(a) + y * math.sin(a)) * r + d * h))
        if prev:
            for j in range(24):
                j2 = (j + 1) % 24
                bm.faces.new((prev[j], prev[j2], ring[j2], ring[j]))
        prev = ring
    for k in range(ribs):
        a = 2.0 * math.pi * k / ribs
        radial = x * math.cos(a) + y * math.sin(a)
        cyl(bm, dish_center + radial * radius * 0.06, dish_center + radial * radius * 0.92 + d * 0.02, 0.03 * s, seg=5)
    feed = dish_center + d * (depth + 1.15 * s * 0.6)
    cyl(bm, dish_center + d * 0.06, feed, 0.035 * s, seg=6)
    cyl(bm, feed - d * 0.06 * s, feed + d * 0.10 * s, 0.09 * s, seg=8)
    for k in range(3):
        a = 2.0 * math.pi * k / 3
        radial = x * math.cos(a) + y * math.sin(a)
        cyl(bm, dish_center + radial * radius * 0.8 + d * depth * 0.9, feed, 0.022 * s, seg=5)


def whipple_panels(bm, p1, p2, radius, rows=6, cols=14, gap=0.04):
    """MMOD shield: arc panels in a regular grid over a cylinder."""
    p1, p2 = Vector(p1), Vector(p2)
    d = (p2 - p1).normalized()
    x, y, _ = basis(d)
    length = (p2 - p1).length
    for i in range(cols):
        t = (i + 0.5) / cols
        for j in range(rows):
            a = 2.0 * math.pi * (j + 0.5) / rows
            radial = x * math.cos(a) + y * math.sin(a)
            tangent = d.cross(radial)
            c = p1 + d * (length * t) + radial * (radius + 0.09)
            box(bm, c, length / cols * 0.44, radius * math.pi / rows * 0.42, 0.035,
                d, tangent, radial)


def equipment_row(bm, start, direction, lateral, count, spacing, s_lo, s_hi, seed, lateral_axis=None):
    """A deliberately arranged row of machinery boxes (not noise)."""
    rng = random.Random(seed)
    start = Vector(start)
    direction = Vector(direction).normalized()
    lateral = Vector(lateral).normalized()
    ax = lateral_axis if lateral_axis is not None else direction.cross(lateral).normalized()
    for i in range(count):
        t = (i - (count - 1) * 0.5) * spacing
        base = start + direction * t
        w = rng.uniform(s_lo, s_hi)
        h = rng.uniform(s_lo, s_hi)
        dep = rng.uniform(s_lo, s_hi * 1.3)
        box(bm, base + lateral * (dep * 0.5), w, h, dep * 0.5, direction, ax, lateral)
        if rng.random() < 0.4:
            cyl(bm, base + lateral * dep * 0.3 - ax * w * 0.6,
                base + lateral * (dep * 0.5) - ax * w * 0.6, w * 0.18, seg=6)


def panel_bank(bm, p1, p2, radius, cols=16, rows=8, depth=0.05):
    """Regular rectangular panel grid over a cylinder (machined surface look)."""
    p1, p2 = Vector(p1), Vector(p2)
    d = (p2 - p1).normalized()
    x, y, _ = basis(d)
    length = (p2 - p1).length
    for i in range(cols):
        t = (i + 0.5) / cols
        for j in range(rows):
            a = 2.0 * math.pi * (j + 0.5) / rows
            radial = x * math.cos(a) + y * math.sin(a)
            tangent = d.cross(radial)
            c = p1 + d * (length * t) + radial * radius
            box(bm, c, length / cols * 0.42, radius * math.pi / rows * 0.40, depth,
                d, tangent, radial)


# --------------------------------------------------------------- assembly


def clear_scene(bpy):
    bpy.ops.wm.read_factory_settings(use_empty=True)


def build_object(bpy, name, material, builder):
    bm = bmesh.new()
    builder(bm)
    bmesh.ops.recalc_face_normals(bm, faces=bm.faces[:])
    mesh = bpy.data.meshes.new(name)
    bm.to_mesh(mesh)
    bm.free()
    obj = bpy.data.objects.new(name, mesh)
    bpy.context.collection.objects.link(obj)
    obj.data.materials.append(material)
    return obj


# Geometry constants shared across steps (units: 1 = 2 m).
# Scale philosophy: few large elements, not many small ones. The hull core
# is a continuous 16 m-diameter pressure body; everything else hangs off it.
HULL_P1 = (0.0, -4.0, 0.0)
HULL_P2 = (0.0, 42.0, 0.0)
HULL_R = 4.0                       # 16 m diameter
NODE_C = (0.0, 44.6, 0.0)
CENTRIFUGE_CENTER = (0.0, 6.0, 0.0)
CENTRIFUGE_R = 14.0                # 56 m ring
SPINE_P1 = (0.0, -34.0, 0.0)
SPINE_P2 = (0.0, 20.0, 0.0)
SPINE_R = 1.2
HAB_A = ((5.6, 12.0, 0.0), (5.6, 30.0, 0.0), 2.6)   # side-mounted lab module
HAB_B = ((-5.6, 12.0, 0.0), (-5.6, 30.0, 0.0), 2.6)  # side-mounted hab module
TANK_X = 7.6
TANK_R = 2.5                       # 10 m diameter tanks
REACTOR_C = (0.0, -30.5, 0.0)
NOZZLE_MOUNT_Y = -40.5


# --------------------------------------------------------------- builders


def b_hull_core(bm):
    """Continuous 16 m-diameter pressure hull — the ship's body."""
    p1, p2 = Vector(HULL_P1), Vector(HULL_P2)
    r = HULL_R
    seg = 56
    cyl(bm, p1, p2, r, seg=seg, cap1=False, cap2=False)
    dome(bm, p1, (0.0, -1.0, 0.0), r, seg=seg)
    dome(bm, p2, (0.0, 1.0, 0.0), r, seg=seg)
    frames = 24
    for i in range(1, frames):
        t = i / frames
        torus(bm, p1.lerp(p2, t), r, 0.115, axis=(0, 1, 0), seg=seg, rseg=7)
    for k in range(8):
        a = 2.0 * math.pi * k / 8 + math.pi / 8
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        cyl(bm, p1 + radial * r, p2 + radial * r, 0.10, seg=7)
    for band_t, count in ((0.18, 22), (0.46, 26), (0.74, 22)):
        c0 = p1.lerp(p2, band_t)
        for k in range(count):
            a = 2.0 * math.pi * k / count
            radial = Vector((math.cos(a), 0.0, math.sin(a)))
            wc = c0 + radial * r
            cyl(bm, wc - radial * 0.12, wc + radial * 0.20, 0.28, seg=12, cap1=False)
            torus(bm, wc + radial * 0.20, 0.305, 0.055, axis=radial, seg=14, rseg=5)
    panel_bank(bm, p1 + Vector((0.0, 3.5, 0.0)), p2 - Vector((0.0, 3.5, 0.0)),
               r + 0.035, cols=72, rows=18, depth=0.05)
    for (t0, count, seed) in ((0.10, 26, 601), (0.60, 30, 602), (0.90, 22, 603)):
        base = p1.lerp(p2, t0)
        equipment_row(bm, base + Vector((0.0, 0.0, r + 0.12)), (0, 1, 0), (0, 0, 1), count, 1.35, 0.14, 0.34, seed=seed)
        equipment_row(bm, base + Vector((0.0, 0.0, -r - 0.12)), (0, 1, 0), (0, 0, -1), count, 1.35, 0.14, 0.34, seed=seed + 50)


def b_spine(bm):
    """Outrigger trusses carrying tanks, arrays and radiators."""
    for sx in (-1.0, 1.0):
        x0 = sx * (HULL_R + 1.6)
        truss(bm, (x0, -34.0, 0.0), (x0, 20.0, 0.0), 0.85, bays=36, chord_r=0.085, strut_r=0.05,
              angles=[0.0, 120.0, 240.0] if sx > 0 else [180.0, 60.0, 300.0])
        for st in (-30.0, -22.0, -14.0, -6.0, 2.0, 10.0):
            cyl(bm, Vector((sx * HULL_R, st, 0.0)), Vector((x0, st, 0.0)), 0.10, seg=7)


def b_node(bm):
    c = Vector(NODE_C)
    cyl(bm, c - Vector((0, 1.3, 0)), c + Vector((0, 1.3, 0)), 2.6, seg=24)
    torus(bm, c - Vector((0, 1.3, 0)), 2.63, 0.09, axis=(0, 1, 0), seg=28, rseg=6)
    torus(bm, c + Vector((0, 1.3, 0)), 2.63, 0.09, axis=(0, 1, 0), seg=28, rseg=6)
    for dirv in ((0, 1, 0), (1, 0, 0), (-1, 0, 0), (0, 0, 1), (0, 0, -1)):
        d = Vector(dirv)
        x, y, _ = basis(d)
        base = c + d * 1.25
        cyl(bm, base, base + d * 0.65, 1.05, seg=18)
        torus(bm, base + d * 0.65, 1.12, 0.085, axis=d, seg=20, rseg=6)
        for k in range(10):
            a = 2.0 * math.pi * k / 10
            radial = x * math.cos(a) + y * math.sin(a)
            cyl(bm, base + d * 0.65 + radial * 1.12, base + d * 0.80 + radial * 1.22, 0.055, seg=6)


def b_hab_a(bm):
    p1, p2, r = HAB_A
    hab_module(bm, p1, p2, r, ribs=9, windows=12, seg=32)
    panel_bank(bm, p1, p2, r + 0.02, cols=30, rows=18, depth=0.04)
    equipment_row(bm, (6.6 + r + 0.12, 14.0, 0.0), (0, 1, 0), (1, 0, 0), 12, 1.25, 0.12, 0.28, seed=201)
    equipment_row(bm, (6.6, 14.5, r + 0.12), (0, 1, 0), (0, 0, 1), 10, 1.4, 0.10, 0.26, seed=203)


def b_hab_b(bm):
    p1, p2, r = HAB_B
    hab_module(bm, p1, p2, r, ribs=9, windows=12, seg=32)
    panel_bank(bm, p1, p2, r + 0.02, cols=30, rows=18, depth=0.04)
    equipment_row(bm, (-6.6 - r - 0.12, 15.0, 0.0), (0, 1, 0), (-1, 0, 0), 12, 1.25, 0.12, 0.28, seed=211)
    equipment_row(bm, (-6.6, 15.5, r + 0.12), (0, 1, 0), (0, 0, 1), 10, 1.4, 0.10, 0.26, seed=213)


def b_whipple(bm):
    p1, p2 = Vector(HULL_P1), Vector(HULL_P2)
    whipple_panels(bm, p1 + Vector((0, 3.0, 0)), p2 - Vector((0, 3.0, 0)), HULL_R, rows=14, cols=30)


def b_waterwall(bm):
    """Water tanks ringed around the forward hull as radiation shielding."""
    for k in range(14):
        a = 2.0 * math.pi * k / 14
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        p1 = Vector((0.0, 22.0, 0.0)) + radial * (HULL_R + 1.15)
        p2 = Vector((0.0, 38.0, 0.0)) + radial * (HULL_R + 1.15)
        pressure_vessel(bm, p1, p2, 1.15, rings=3, seg=16, ports=False)
        mli_bands(bm, p1, p2, 1.18, bands=7, seg=18, wobble=0.02)
        for t in (0.25, 0.75):
            cco = p1.lerp(p2, t)
            cyl(bm, cco - radial * 1.1, cco - radial * 0.4, 0.08, seg=6)


def b_centrifuge(bm):
    centrifuge_ring(bm, CENTRIFUGE_CENTER, CENTRIFUGE_R, axis=(0, 1, 0), cars=8, seg=72)
    hub = Vector(CENTRIFUGE_CENTER)
    cyl(bm, hub - Vector((0, 2.2, 0)), hub + Vector((0, 2.2, 0)), 1.9, seg=24)
    torus(bm, hub + Vector((0, 2.2, 0)), 1.93, 0.09, axis=(0, 1, 0), seg=28, rseg=6)
    torus(bm, hub - Vector((0, 2.2, 0)), 1.93, 0.09, axis=(0, 1, 0), seg=28, rseg=6)
    for k in range(8):
        a = 2.0 * math.pi * k / 8
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        cyl(bm, hub + radial * 1.9, hub + radial * 3.0, 0.11, seg=7)
    for k in range(32):
        a = 2.0 * math.pi * k / 32
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        c = Vector(CENTRIFUGE_CENTER) + radial * (CENTRIFUGE_R - 1.9)
        cyl(bm, c - Vector((0, 1.1, 0)), c + Vector((0, 1.1, 0)), 0.62, seg=16)
        torus(bm, c - Vector((0, 1.1, 0)), 0.64, 0.055, axis=(0, 1, 0), seg=16, rseg=5)
        torus(bm, c + Vector((0, 1.1, 0)), 0.64, 0.055, axis=(0, 1, 0), seg=16, rseg=5)
    for k in range(28):
        a = 2.0 * math.pi * k / 28 + 0.22
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        c = Vector(CENTRIFUGE_CENTER) + radial * (CENTRIFUGE_R + 0.4)
        box(bm, c, 0.38, 0.62, 0.38, radial, Vector((0, 1, 0)), Vector((-math.sin(a), 0.0, math.cos(a))))


def b_prop_tank(bm):
    for sx in (-1.0, 1.0):
        for yc in (-14.0, -23.0):
            x0 = sx * TANK_X
            p1 = (x0, yc - 4.5, 0.0)
            p2 = (x0, yc + 4.5, 0.0)
            pressure_vessel(bm, p1, p2, TANK_R, rings=5, seg=30)
            mli_bands(bm, p1, p2, TANK_R + 0.05, bands=13, seg=30, wobble=0.03)
            panel_bank(bm, Vector(p1) + Vector((0, 1.4, 0)), Vector(p2) - Vector((0, 1.4, 0)),
                       TANK_R + 0.05, cols=22, rows=16, depth=0.035)
            for st in (yc - 3.0, yc, yc + 3.0):
                cyl(bm, Vector((sx * (HULL_R + 1.6), st, 0.0)), Vector((x0, st, 0.0)), 0.16, seg=8)
                cyl(bm, Vector((sx * (HULL_R + 1.6), st - 0.5, 0.0)), Vector((x0, st, 0.6)), 0.07, seg=6)
                torus(bm, Vector((x0, st, 0.0)), TANK_R + 0.07, 0.06, axis=(0, 1, 0), seg=24, rseg=5)


def b_reactor(bm):
    base = Vector(REACTOR_C)
    cyl(bm, base - Vector((0, 1.8, 0)), base + Vector((0, 1.8, 0)), 2.1, seg=26)
    for yy in (-1.8, 1.8):
        torus(bm, base + Vector((0, yy, 0)), 2.13, 0.09, axis=(0, 1, 0), seg=28, rseg=6)
    for k in range(14):
        a = 2.0 * math.pi * k / 14
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        cyl(bm, base + radial * 2.0 - Vector((0, 1.5, 0)), base + radial * 2.0 + Vector((0, 1.5, 0)), 0.11, seg=6)
    for k in range(8):
        a = 2.0 * math.pi * k / 8
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        cyl(bm, base + radial * 2.1, base + radial * 3.1 - Vector((0, 3.4, 0)), 0.10, seg=7)


def b_thrust_cone(bm):
    # Aft transition structure: heavy beams + ring frames that visibly tie the
    # hull to the propulsion section instead of leaving a bare truss gap.
    for k in range(4):
        a = math.pi / 4.0 + k * math.pi / 2.0
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        cyl(bm, Vector((0.0, -3.0, 0.0)) + radial * (HULL_R * 0.92),
            Vector((0.0, -30.0, 0.0)) + radial * 2.6, 0.26, seg=10)
    for t in (0.12, 0.32, 0.52, 0.72, 0.9):
        yy = -3.0 - 27.0 * t
        rr = HULL_R * 0.92 + (2.6 - HULL_R * 0.92) * t
        torus(bm, Vector((0.0, yy, 0.0)), rr, 0.11, axis=(0, 1, 0), seg=24, rseg=6)
    # propellant feed lines: tanks -> reactor -> engines
    for sx in (-1.0, 1.0):
        for dz in (0.6, -0.6):
            pipe_run(bm, [
                (sx * TANK_X, -25.0, dz * TANK_R),
                (sx * 4.6, -31.0, dz * 2.3),
                (sx * 2.2, -36.0, dz * 1.5),
                (sx * 1.3, -40.0, dz * 0.9),
            ], 0.15)
    thrust_cone(bm, (0.0, -33.5, 0.0), (0, -1, 0), 3.6, 2.9, 1.8, struts=14)
    shadow_shield(bm, (0.0, -35.6, 0.0), (0, -1, 0), 4.4, seg=36, ribs=12)


def b_nozzle_cluster(bm):
    for ex in (-2.4, 2.4):
        mount = Vector((ex, NOZZLE_MOUNT_Y, 0.0))
        gimbal(bm, mount, (0, -1, 0), radius=0.55, s=2.0)
        ntp_engine(bm, mount, (0, -1, 0), scale=2.0)


def b_solar(bm):
    for sx in (-1.0, 1.0):
        for yy in (24.0, -18.0):
            root = Vector((sx * (HULL_R + 2.4), yy, 0.0))
            solar_array(bm, root, (sx, 0, 0), 30.0, 13.0, cells=84, rows=18)
            d = Vector((sx, 0, 0))
            x, _, _ = basis(d)
            for k in range(16):
                t = (k + 0.5) / 16
                c = root + d * (30.0 * t)
                cyl(bm, c - x * 6.5, c + x * 6.5, 0.05, seg=6)
            for side in (-1.0, 1.0):
                off = x * 6.5 * side
                cyl(bm, root + off, root + d * 30.0 + off, 0.07, seg=6)
            cyl(bm, root, root + Vector((sx * 2.4, 0, 0)), 0.20, seg=10)
            for k in range(4):
                a = 2.0 * math.pi * k / 4
                cyl(bm, Vector((sx * (HULL_R + 0.4), yy, 0.0)),
                    root + Vector((sx * 1.2, 0, 0)) + x * 0.8 * math.cos(a), 0.06, seg=6)


def b_radiator(bm):
    for sx in (-1.0, 1.0):
        for yy in (-32.0, -12.0, 32.0):
            p1 = Vector((sx * (HULL_R + 3.4), yy - 6.0, 0.0))
            p2 = Vector((sx * (HULL_R + 3.4), yy + 6.0, 0.0))
            fin_radiator(bm, p1, p2, 3.4, fins=64, fin_h=0.55)
            mid = p1.lerp(p2, 0.5)
            cyl(bm, Vector((sx * (HULL_R + 0.6), mid.y, 0.0)), mid, 0.11, seg=8)
            torus(bm, mid, 1.7, 0.06, axis=(1, 0, 0), seg=20, rseg=5)


def b_airlock(bm):
    p1 = Vector((0.0, 16.0, HULL_R - 0.5))
    p2 = Vector((0.0, 16.0, HULL_R + 4.6))
    cyl(bm, p1, p2, 1.6, seg=22)
    for zz in (HULL_R - 0.5, HULL_R + 4.6):
        torus(bm, Vector((0.0, 16.0, zz)), 1.63, 0.08, axis=(0, 0, 1), seg=24, rseg=6)
    hatch = Vector((0.0, 16.0, HULL_R + 4.6 + 0.4))
    cyl(bm, hatch - Vector((0, 0, 0.4)), hatch + Vector((0, 0, 0.25)), 1.0, seg=18)
    torus(bm, hatch + Vector((0, 0, 0.25)), 1.06, 0.075, axis=(0, 0, 1), seg=20, rseg=6)
    for k in range(8):
        a = 2.0 * math.pi * k / 8
        off = Vector((math.cos(a), math.sin(a), 0.0)) * 1.06
        cyl(bm, hatch + off + Vector((0, 0, 0.25)), hatch + off + Vector((0, 0, 0.42)), 0.055, seg=6)
    equipment_row(bm, (1.9, 16.0, HULL_R + 2.4), (0, 0, 1), (1, 0, 0), 8, 0.5, 0.10, 0.22, seed=401)


def b_orion(bm):
    """Orion-class shuttle docked forward — the human-scale reference."""
    base = Vector((0.0, 47.6, 0.0))
    cyl(bm, base, base + Vector((0, 5.4, 0)), 1.5, seg=20)
    for t in (0.3, 0.62):
        torus(bm, base + Vector((0, 5.4 * t, 0)), 1.53, 0.065, axis=(0, 1, 0), seg=22, rseg=5)
    cyl(bm, base + Vector((0, 5.4, 0)), base + Vector((0, 7.8, 0)), 1.5, 0.55, seg=20)
    torus(bm, base + Vector((0, 5.4, 0)), 1.53, 0.08, axis=(0, 1, 0), seg=22, rseg=6)
    cyl(bm, base + Vector((0, 7.8, 0)), base + Vector((0, 8.5, 0)), 0.55, seg=14)
    torus(bm, base + Vector((0, 8.5, 0)), 0.65, 0.07, axis=(0, 1, 0), seg=16, rseg=6)
    for k in range(8):
        a = 2.0 * math.pi * k / 8
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        c = base + Vector((0, 6.9, 0)) + radial * 1.06
        cyl(bm, c - radial * 0.05, c + radial * 0.11, 0.13, seg=10, cap1=False)
    for k in range(4):
        a = 2.0 * math.pi * k / 4 + math.pi / 4
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        root = base + Vector((0, 5.4, 0)) + radial * 1.5
        solar_array(bm, root, radial, 4.6, 2.0, cells=22, rows=8)
    bell(bm, base - Vector((0, 0.6, 0)), (0, -1, 0), 1.5, 0.22, 0.68, seg=22, layers=9)


def five_star(bm, center, normal, radius, depth=0.09):
    """Filled five-pointed star stamped on the hull (PLA insignia)."""
    n = Vector(normal).normalized()
    x, y, _ = basis(n)
    c = Vector(center)
    front, back = [], []
    for i in range(10):
        a = math.pi / 2.0 + i * math.pi / 5.0
        rr = radius if i % 2 == 0 else radius * 0.382
        p = c + (x * math.cos(a) + y * math.sin(a)) * rr
        front.append(bm.verts.new(p + n * depth * 0.5))
        back.append(bm.verts.new(p - n * depth * 0.5))
    cf = bm.verts.new(c + n * depth * 0.5)
    cb = bm.verts.new(c - n * depth * 0.5)
    for i in range(10):
        j = (i + 1) % 10
        bm.faces.new((cf, front[i], front[j]))
        bm.faces.new((cb, back[j], back[i]))
        bm.faces.new((front[i], back[i], back[j], front[j]))


def b_stripes(bm):
    """Twin gold stripes down both flanks — PLA Space Force colour band."""
    p1, p2 = Vector(HULL_P1), Vector(HULL_P2)
    length = (p2 - p1).length
    for base_deg in (0.0, 180.0):
        for offset in (-20.0, 20.0):
            a = math.radians(base_deg + offset)
            radial = Vector((math.cos(a), 0.0, math.sin(a)))
            tangent = radial.cross(Vector((0.0, 1.0, 0.0))).normalized()
            c = (p1 + p2) * 0.5 + radial * (HULL_R + 0.14)
            box(bm, c, 0.035, length * 0.455, 0.26, radial, Vector((0.0, 1.0, 0.0)), tangent)


def b_insignia(bm):
    """Red star insignia placed where nothing occludes them."""
    for sx in (-1.0, 1.0):
        a = 0.0 if sx > 0.0 else math.pi
        radial = Vector((math.cos(a), 0.0, math.sin(a)))
        # forward hull, in the clear band between the water wall and the node
        five_star(bm, Vector((0.0, 40.2, 0.0)) + radial * (HULL_R + 0.06), radial, 1.9, 0.13)
        # outboard face of the side-mounted modules
        five_star(bm, Vector((sx * (6.6 + 2.6 + 0.05), 21.0, 0.0)), (sx, 0.0, 0.0), 1.5, 0.11)


def b_glow(bm):
    """Heat glow rings inside the engine bells."""
    for ex in (-2.4, 2.4):
        mount = Vector((ex, NOZZLE_MOUNT_Y, 0.0))
        chamber = mount + Vector((0.0, -1.0, 0.0))
        throat = chamber + Vector((0.0, -0.9, 0.0))
        for t in (0.32, 0.58, 0.84):
            c = throat + Vector((0.0, -4.6 * t, 0.0))
            r = 0.62 + (2.10 - 0.62) * (t ** 1.7)
            torus(bm, c, r * 0.90, 0.07, axis=(0, 1, 0), seg=26, rseg=5)


def b_glass(bm):
    """Lit viewports along the hull bands."""
    p1, p2 = Vector(HULL_P1), Vector(HULL_P2)
    for band_t, count in ((0.18, 22), (0.46, 26), (0.74, 22)):
        c0 = p1.lerp(p2, band_t)
        for k in range(count):
            a = 2.0 * math.pi * k / count
            radial = Vector((math.cos(a), 0.0, math.sin(a)))
            wc = c0 + radial * (HULL_R + 0.18)
            cyl(bm, wc - radial * 0.02, wc + radial * 0.06, 0.23, seg=10)


def b_greeble(bm):
    """Organised machinery, cable looms, nav lights, hatch plates."""
    p1, p2 = Vector(HULL_P1), Vector(HULL_P2)
    r = HULL_R
    for rx, (yy, ang) in ((0, (28.0, 0.9)), (0, (36.0, 2.5)), (0, (10.0, 4.2))):
        radial = Vector((math.cos(ang), 0.0, math.sin(ang)))
        base = Vector((0.0, yy, 0.0)) + radial * r
        cyl(bm, base - radial * 0.05, base + radial * 0.13, 0.55, seg=14)
        torus(bm, base + radial * 0.13, 0.60, 0.055, axis=radial, seg=16, rseg=5)
    for sx in (-1.0, 1.0):
        x0 = sx * (HULL_R + 1.6)
        for dz in (0.95, -0.95):
            pipe_run(bm, [(x0, -34.0, dz), (x0, 20.0, dz)], 0.055)
    for (pos, dirv, rad) in (((-4.0, 40.0, 0.0), (-1, 0, 0), 0.20),
                             ((4.0, 40.0, 0.0), (1, 0, 0), 0.20),
                             ((0.0, 40.0, 4.0), (0, 0, 1), 0.20),
                             ((-4.0, -2.0, 0.0), (-1, 0, 0), 0.16),
                             ((4.0, -2.0, 0.0), (1, 0, 0), 0.16)):
        d = Vector(dirv)
        cyl(bm, Vector(pos), Vector(pos) + d * 0.15, rad, seg=10)
        cyl(bm, Vector(pos) + d * 0.15, Vector(pos) + d * 0.34, rad * 0.7, seg=10)
    equipment_row(bm, (0.0, 34.0, r + 0.15), (0, 1, 0), (0, 0, 1), 12, 1.0, 0.10, 0.24, seed=501)
    # EVA handrail runs (crew translation paths along the hull)
    for sx in (-1.0, 1.0):
        for (yy0, count) in ((24.0, 10), (4.0, 12)):
            a = math.radians(-14.0 if sx > 0.0 else 194.0)
            radial = Vector((math.cos(a), 0.0, math.sin(a)))
            tangent = Vector((0.0, 1.0, 0.0))
            for k in range(count):
                yy = yy0 + k * 1.35
                base = Vector((0.0, yy, 0.0)) + radial * HULL_R
                tip = base + radial * 0.42
                cyl(bm, base, tip, 0.045, seg=6)
                if k < count - 1:
                    cyl(bm, tip, tip + tangent * 1.35, 0.035, seg=5)
    equipment_row(bm, (0.0, 6.0, -r - 0.15), (0, 1, 0), (0, 0, -1), 14, 1.0, 0.10, 0.24, seed=502)


def b_antenna(bm):
    antenna_dish(bm, (0.0, 30.0, HULL_R + 1.4), (0.0, 0.3, 1.0), 2.6, depth=0.7, ribs=12)
    antenna_dish(bm, (HULL_R + 1.2, 38.0, 0.0), (1.0, 0.25, 0.1), 1.4, depth=0.4, ribs=9)
    for (bx, by, bz, ll) in ((0.0, 42.5, 3.0, 2.2), (0.0, 8.0, -3.0, 2.0), (3.0, 41.0, 0.0, 1.6)):
        cyl(bm, Vector((bx, by, bz)), Vector((bx, by + 0.2, bz + ll)), 0.055, seg=6)
        torus(bm, Vector((bx, by + 0.2, bz + ll)), 0.4, 0.045, axis=(0, 0, 1), seg=16, rseg=5)


BUILDERS = {
    "hnx02.hull.core": b_hull_core,
    "hnx02.spine.truss": b_spine,
    "hnx02.node.forward": b_node,
    "hnx02.hab.a": b_hab_a,
    "hnx02.hab.b": b_hab_b,
    "hnx02.hab.whipple": b_whipple,
    "hnx02.rad.waterwall": b_waterwall,
    "hnx02.centrifuge.ring": b_centrifuge,
    "hnx02.prop.tank": b_prop_tank,
    "hnx02.reactor.core": b_reactor,
    "hnx02.thrust.cone": b_thrust_cone,
    "hnx02.nozzle.cluster": b_nozzle_cluster,
    "hnx02.power.solar": b_solar,
    "hnx02.power.radiator": b_radiator,
    "hnx02.cabin.airlock": b_airlock,
    "hnx02.shuttle.orion": b_orion,
    "hnx02.detail.greeble": b_greeble,
    "hnx02.detail.antenna": b_antenna,
    "hnx02.detail.glow": b_glow,
    "hnx02.detail.glass": b_glass,
    "hnx02.livery.stripes": b_stripes,
    "hnx02.livery.insignia": b_insignia,
}


def build_geometry(bpy, materials, step_id):
    objects = []
    for object_id in STEP_OBJECTS[step_id]:
        material = materials[OBJECT_MATERIAL[object_id]]
        objects.append(build_object(bpy, object_id, material, BUILDERS[object_id]))
    return objects


# ------------------------------------------------------------ render / io


def add_render_setup(bpy):
    scene = bpy.context.scene
    scene.render.resolution_x = 1280
    scene.render.resolution_y = 860
    scene.render.resolution_percentage = 100
    scene.render.image_settings.file_format = "PNG"
    if scene.world is None:
        scene.world = bpy.data.worlds.new("hnx02.world")
    scene.world.color = (0.02, 0.03, 0.04)

    camera_data = bpy.data.cameras.new("hnx02.camera")
    camera_data.type = "ORTHO"
    camera_data.ortho_scale = 118.0
    camera = bpy.data.objects.new("hnx02.camera", camera_data)
    bpy.context.collection.objects.link(camera)
    camera.location = (52.0, -66.0, 40.0)
    camera.rotation_euler = (Vector((0.0, 2.0, 0.0)) - camera.location).to_track_quat("-Z", "Y").to_euler()
    scene.camera = camera

    def add_light(name, location, energy, size, color=(1.0, 1.0, 1.0)):
        light_data = bpy.data.lights.new(name, type="AREA")
        light_data.energy = energy
        light_data.shape = "DISK"
        light_data.size = size
        light_data.color = color
        light = bpy.data.objects.new(name, light_data)
        bpy.context.collection.objects.link(light)
        light.location = location
        light.rotation_euler = (Vector((0.0, 2.0, 0.0)) - Vector(location)).to_track_quat("-Z", "Y").to_euler()
        return light

    add_light("hnx02.key", (34.0, -30.0, 40.0), 16000, 26.0, (1.0, 0.97, 0.92))
    add_light("hnx02.fill", (-30.0, -42.0, 10.0), 5200, 30.0, (0.68, 0.79, 1.0))
    add_light("hnx02.rim", (-36.0, 26.0, -16.0), 9000, 16.0, (1.0, 0.85, 0.58))


def export_glb(bpy, objects, output):
    for obj in bpy.context.selected_objects:
        obj.select_set(False)
    for obj in objects:
        obj.select_set(True)
    bpy.context.view_layer.objects.active = objects[-1]
    result = bpy.ops.export_scene.gltf(
        filepath=output,
        export_format="GLB",
        use_selection=True,
        export_apply=True,
        export_draco_mesh_compression_enable=True,
        export_draco_mesh_compression_level=6,
        export_draco_position_quantization=14,
        export_draco_normal_quantization=10,
        export_draco_texcoord_quantization=12,
    )
    if "FINISHED" not in result or not os.path.isfile(output) or os.path.getsize(output) <= 0:
        fail("GLB export failed", 4)


def write_result(path, response):
    path = os.path.abspath(path)
    os.makedirs(os.path.dirname(path), exist_ok=True)
    temporary_path = path + ".tmp"
    with open(temporary_path, "w", encoding="utf-8") as stream:
        json.dump(response, stream, separators=(",", ":"))
        stream.flush()
        os.fsync(stream.fileno())
    os.replace(temporary_path, path)


def verify_reimport(bpy, output, expected):
    clear_scene(bpy)
    result = bpy.ops.import_scene.gltf(filepath=output)
    if "FINISHED" not in result:
        fail("GLB reimport failed", 4)
    names = [obj.name for obj in bpy.context.scene.objects if obj.type == "MESH"]
    if len(names) != len(expected) or set(names) != set(expected):
        fail("GLB reimport object IDs did not match", 4)
    return [object_id for object_id in expected if object_id in names]


def main():
    parser = argparse.ArgumentParser(add_help=False)
    parser.add_argument("--step-id", required=True)
    parser.add_argument("--output", required=True)
    parser.add_argument("--render-output")
    parser.add_argument("--result-output", required=True)
    args = parser.parse_args()

    step_id = args.step_id
    if step_id not in STEP_OBJECTS:
        fail("unknown step id: " + step_id, 3)

    clear_scene(bpy)
    materials = build_materials(bpy)
    objects = build_geometry(bpy, materials, step_id)

    render_path = None
    if args.render_output:
        add_render_setup(bpy)
        scene = bpy.context.scene
        scene.render.filepath = args.render_output
        bpy.ops.render.render(write_still=True)
        if os.path.isfile(args.render_output):
            render_path = os.path.abspath(args.render_output)

    exported_at = utc_now()
    export_glb(bpy, objects, args.output)
    reimported = verify_reimport(bpy, args.output, STEP_OBJECTS[step_id])

    response = {
        "ok": True,
        "blenderVersion": ".".join(str(part) for part in bpy.app.version),
        "exportedAt": exported_at,
        "objectIds": STEP_OBJECTS[step_id],
        "verifiedObjectIds": reimported,
        "glbBytes": os.path.getsize(args.output),
    }
    if render_path is not None:
        response["renderPath"] = render_path
        response["renderedAt"] = utc_now()
    write_result(args.result_output, response)
    print("hnx02 result written: " + args.result_output)


if __name__ == "__main__":
    main()