#!/usr/bin/env python3
"""
build_science-explorer.py
Builds the deep-space exploration and science vessel (science-explorer.glb)
from an empty Blender 5.2 scene, bakes CC0 Poly Haven + procedural materials
with Cycles, verifies moving pivots, exports GLB + nodes.json + NOTES.md,
and renders final.png, clay.png, side.png, and paints.png.
"""

import bpy
import bmesh
import math
import os
import sys
import json
import struct
import urllib.request
import numpy as np
from mathutils import Vector, Matrix, Euler

ROOT = os.path.dirname(os.path.abspath(__file__))
TEX_DIR = os.path.join(ROOT, "textures")
os.makedirs(TEX_DIR, exist_ok=True)

# ---------------------------------------------------------------------------
# 1. Poly Haven CC0 Texture Download & Loading
# ---------------------------------------------------------------------------
POLYHAVEN_IDS = [
    "metal_plate_02",
    "metal_plate",
    "painted_metal_shutter",
    "blue_metal_plate",
    "rusty_metal_03",
]

def ensure_polyhaven_textures():
    headers = {"User-Agent": "Mozilla/5.0 (ForgeAssetBuilder/1.0)"}
    for ph_id in POLYHAVEN_IDS:
        for map_type in ("diff", "nor_gl", "arm"):
            fname = f"{ph_id}_{map_type}_1k.jpg"
            fpath = os.path.join(TEX_DIR, fname)
            if os.path.exists(fpath) and os.path.getsize(fpath) > 10000:
                continue
            url = f"https://dl.polyhaven.org/file/ph-assets/Textures/jpg/1k/{ph_id}/{fname}"
            print(f"Downloading {url} -> {fpath}", flush=True)
            req = urllib.request.Request(url, headers=headers)
            with urllib.request.urlopen(req, timeout=30) as resp, open(fpath, "wb") as out:
                out.write(resp.read())

ensure_polyhaven_textures()

# ---------------------------------------------------------------------------
# 2. Scene Reset & Cycles GPU Configuration
# ---------------------------------------------------------------------------
bpy.ops.wm.read_factory_settings(use_empty=True)
sc = bpy.context.scene
sc.unit_settings.system = "METRIC"
sc.unit_settings.scale_length = 1.0

sc.render.engine = "CYCLES"
prefs = bpy.context.preferences.addons["cycles"].preferences
prefs.refresh_devices()
for dev_type in ("OPTIX", "CUDA"):
    try:
        prefs.compute_device_type = dev_type
        for d in prefs.devices:
            d.use = True
        sc.cycles.device = "GPU"
        break
    except Exception:
        pass

sc.cycles.use_denoising = True
sc.view_settings.view_transform = "Standard"
sc.view_settings.look = "None"

def load_image(fname, non_color=False):
    path = os.path.join(TEX_DIR, fname)
    img = bpy.data.images.load(path, check_existing=True)
    img.colorspace_settings.name = "Non-Color" if non_color else "sRGB"
    return img

# ---------------------------------------------------------------------------
# 3. Paint Palette & Procedural + Poly Haven Material Baking (Cycles)
# ---------------------------------------------------------------------------
PAINTS = [
    ("Expedition White", [0.82, 0.82, 0.80]),
    ("Industrial Orange", [0.78, 0.24, 0.04]),
    ("Safety Yellow", [0.85, 0.62, 0.05]),
    ("Deep Space Teal", [0.08, 0.42, 0.48]),
    ("Signal Red", [0.68, 0.07, 0.07]),
    ("Orbital Charcoal", [0.14, 0.16, 0.18]),
]

TEX_FILES = {}

def generate_procedural_baked_maps():
    """
    Combines CC0 Poly Haven 1k maps (diff, nor_gl, arm) with spacecraft panel/trim
    weathering and saves game-weight maps (<=1024px basecolor/normal, <=512px ORM).
    Runs a Cycles GPU bake pass for Paint_basecolor.jpg.
    """
    def read_tex_pixels(fname, target_res):
        img = load_image(fname, non_color=("nor" in fname or "arm" in fname))
        img_copy = img.copy()
        img_copy.scale(target_res, target_res)
        arr = np.array(img_copy.pixels[:], dtype=np.float32).reshape(target_res, target_res, 4)
        bpy.data.images.remove(img_copy)
        return arr[:, :, :3]

    def save_map(mat_name, slot, rgb_or_rgba, is_data=False, use_jpg=True):
        h, w, c = rgb_or_rgba.shape
        ext = "jpg" if (use_jpg and c == 3) else "png"
        fname = f"{mat_name}_{slot}.{ext}"
        out_path = os.path.join(TEX_DIR, fname)
        im = bpy.data.images.new(fname, width=w, height=h, alpha=(c == 4), is_data=is_data)
        im.colorspace_settings.name = "Non-Color" if is_data else "sRGB"
        if c == 3:
            rgba = np.ones((h, w, 4), dtype=np.float32)
            rgba[:, :, :3] = np.clip(rgb_or_rgba, 0.0, 1.0)
        else:
            rgba = np.clip(rgb_or_rgba, 0.0, 1.0)
        im.pixels.foreach_set(rgba.ravel())
        im.filepath_raw = out_path
        if ext == "jpg":
            sc.render.image_settings.file_format = "JPEG"
            sc.render.image_settings.quality = 92
            im.save_render(out_path, scene=sc)
        else:
            im.file_format = "PNG"
            im.save()
        bpy.data.images.remove(im)
        TEX_FILES[(mat_name, slot)] = fname
        return fname

    y1024, x1024 = np.mgrid[0:1024, 0:1024].astype(np.float32) / 1024.0
    y512, x512 = np.mgrid[0:512, 0:512].astype(np.float32) / 512.0

    # --- 1. PAINT (1024 basecolor/normal, 512 orm) ---
    # Light neutral detail ONLY (near white, subtle grey grime streaks, fine panel seams,
    # micro-meteorite pitting; multiplied at runtime by Paint baseColorFactor).
    mp2_diff = read_tex_pixels("metal_plate_02_diff_1k.jpg", 1024)
    mp2_nor = read_tex_pixels("metal_plate_02_nor_gl_1k.jpg", 1024)
    mp2_arm = read_tex_pixels("metal_plate_02_arm_1k.jpg", 512)
    sh_nor = read_tex_pixels("painted_metal_shutter_nor_gl_1k.jpg", 1024)

    lum = (0.299 * mp2_diff[:, :, 0] + 0.587 * mp2_diff[:, :, 1] + 0.114 * mp2_diff[:, :, 2])
    lum_norm = 0.93 + (lum - np.mean(lum)) * 0.14

    # Delicate hair-thin panel lines + subtle plate tonal variation (matches WEATHERED HULL PAINT swatch)
    gx = np.abs(((x1024 * 4.0) % 1.0) - 0.5) * 2.0
    gy = np.abs(((y1024 * 4.0) % 1.0) - 0.5) * 2.0
    seam_main = np.maximum((gx > 0.982).astype(np.float32), (gy > 0.982).astype(np.float32))
    gx2 = np.abs((((x1024 + 0.25 * np.floor(y1024 * 4.0)) * 8.0) % 1.0) - 0.5) * 2.0
    gy2 = np.abs(((y1024 * 8.0) % 1.0) - 0.5) * 2.0
    seam_sub = np.maximum((gx2 > 0.986).astype(np.float32), (gy2 > 0.986).astype(np.float32)) * 0.5

    streak = 0.022 * np.sin(x1024 * 2.0 * math.pi * 7.0) * np.cos(y1024 * 2.0 * math.pi * 2.0)
    plate_var = 0.018 * (((np.floor(x1024 * 4.0) * 7 + np.floor(y1024 * 4.0) * 13) % 5) / 4.0 - 0.5)
    val = np.clip(lum_norm + plate_var - streak - seam_main * 0.065 - seam_sub * 0.035, 0.80, 0.98)
    paint_bc = np.stack([val, val, val], axis=-1)

    dx_seam = np.roll(seam_main + seam_sub, -1, axis=1) - np.roll(seam_main + seam_sub, 1, axis=1)
    dy_seam = np.roll(seam_main + seam_sub, -1, axis=0) - np.roll(seam_main + seam_sub, 1, axis=0)
    paint_nor = mp2_nor * 0.82 + sh_nor * 0.18
    paint_nor[:, :, 0] = np.clip(paint_nor[:, :, 0] + dx_seam * 0.16, 0.0, 1.0)
    paint_nor[:, :, 1] = np.clip(paint_nor[:, :, 1] + dy_seam * 0.16, 0.0, 1.0)

    seam_512 = seam_main[::2, ::2]
    paint_ao = np.clip(mp2_arm[:, :, 0] * (1.0 - seam_512 * 0.15), 0.55, 1.0)
    paint_rough = np.clip(0.38 + (mp2_arm[:, :, 1] - 0.5) * 0.20 + seam_512 * 0.12, 0.24, 0.68)
    paint_metal = np.clip(0.05 + seam_512 * 0.22, 0.0, 0.30)
    paint_orm = np.stack([paint_ao, paint_rough, paint_metal], axis=-1)

    save_map("Paint", "basecolor", paint_bc, is_data=False, use_jpg=True)
    save_map("Paint", "normal", paint_nor, is_data=True, use_jpg=True)
    save_map("Paint", "orm", paint_orm, is_data=True, use_jpg=True)

    # --- 2. METAL (Bright brushed aerospace alloy matching BARE METAL swatch) ---
    mp_diff = read_tex_pixels("metal_plate_diff_1k.jpg", 512)
    mp_nor = read_tex_pixels("metal_plate_nor_gl_1k.jpg", 512)
    mp_arm = read_tex_pixels("metal_plate_arm_1k.jpg", 512)
    # Brushed sheen along U
    brush = 0.04 * np.sin(y512 * 2.0 * math.pi * 64.0)
    metal_bc = np.clip(0.58 + (mp_diff - np.mean(mp_diff)) * 0.35 + brush[:, :, None], 0.42, 0.86)
    metal_ao = np.clip(mp_arm[:, :, 0], 0.55, 1.0)
    metal_rough = np.clip(0.26 + (mp_arm[:, :, 1] - 0.5) * 0.18, 0.16, 0.45)
    metal_met = np.full((512, 512), 0.86, dtype=np.float32)
    save_map("Metal", "basecolor", metal_bc, is_data=False, use_jpg=True)
    save_map("Metal", "normal", mp_nor, is_data=True, use_jpg=True)
    save_map("Metal", "orm", np.stack([metal_ao, metal_rough, metal_met], axis=-1), is_data=True, use_jpg=True)

    # --- 3. SCORCHED (Gunmetal/charcoal-bronze heat-scorched titanium matching HEAT-SCORCHED METAL swatch) ---
    rm_diff = read_tex_pixels("rusty_metal_03_diff_1k.jpg", 512)
    rm_nor = read_tex_pixels("rusty_metal_03_nor_gl_1k.jpg", 512)
    rm_arm = read_tex_pixels("rusty_metal_03_arm_1k.jpg", 512)
    t_grad = y512
    c_gunmetal = np.array([0.22, 0.22, 0.24], dtype=np.float32)
    c_bronze = np.array([0.42, 0.30, 0.20], dtype=np.float32)
    c_violet = np.array([0.25, 0.27, 0.34], dtype=np.float32)
    w_br = np.exp(-((t_grad - 0.48) * 3.5) ** 2)[:, :, None]
    w_vi = np.clip((t_grad - 0.60) * 2.5, 0.0, 1.0)[:, :, None]
    sc_tint = c_gunmetal * (1.0 - w_br) + c_bronze * w_br
    sc_tint = sc_tint * (1.0 - w_vi * 0.55) + c_violet * (w_vi * 0.55)
    scorch_bc = np.clip(sc_tint * (0.72 + 0.45 * rm_diff), 0.12, 0.62)
    scorch_ao = np.clip(rm_arm[:, :, 0], 0.45, 1.0)
    scorch_rough = np.clip(0.32 + (rm_arm[:, :, 1] - 0.5) * 0.22, 0.20, 0.58)
    scorch_met = np.full((512, 512), 0.78, dtype=np.float32)
    save_map("Scorched", "basecolor", scorch_bc, is_data=False, use_jpg=True)
    save_map("Scorched", "normal", rm_nor, is_data=True, use_jpg=True)
    save_map("Scorched", "orm", np.stack([scorch_ao, scorch_rough, scorch_met], axis=-1), is_data=True, use_jpg=True)

    # --- 4. SOLAR (Deep cobalt-blue photovoltaic cells with crisp silver micro-grid) ---
    bm_diff = read_tex_pixels("blue_metal_plate_diff_1k.jpg", 1024)
    bm_nor = read_tex_pixels("blue_metal_plate_nor_gl_1k.jpg", 512)
    bm_arm = read_tex_pixels("blue_metal_plate_arm_1k.jpg", 512)
    cx = np.abs(((x1024 * 16.0) % 1.0) - 0.5) * 2.0
    cy = np.abs(((y1024 * 32.0) % 1.0) - 0.5) * 2.0
    cell_grid = np.maximum((cx > 0.93).astype(np.float32), (cy > 0.93).astype(np.float32))
    wafer_col = np.array([0.06, 0.17, 0.40], dtype=np.float32) + bm_diff * np.array([0.08, 0.16, 0.30], dtype=np.float32)
    silver_col = np.array([0.76, 0.80, 0.86], dtype=np.float32)
    solar_bc = wafer_col * (1.0 - cell_grid[:, :, None] * 0.78) + silver_col * (cell_grid[:, :, None] * 0.78)
    solar_ao = np.clip(bm_arm[:, :, 0], 0.6, 1.0)
    solar_rough = np.clip(0.14 + cell_grid[::2, ::2] * 0.16, 0.10, 0.36)
    solar_met = np.clip(0.55 + cell_grid[::2, ::2] * 0.35, 0.45, 0.90)
    save_map("Solar", "basecolor", solar_bc, is_data=False, use_jpg=True)
    save_map("Solar", "normal", bm_nor, is_data=True, use_jpg=True)
    save_map("Solar", "orm", np.stack([solar_ao, solar_rough, solar_met], axis=-1), is_data=True, use_jpg=True)

    # --- 5. INTERIOR (256 basecolor/normal/orm) ---
    int_bc = np.clip(mp_diff[::2, ::2] * np.array([0.22, 0.25, 0.28], dtype=np.float32), 0.08, 0.35)
    int_nor = mp_nor[::2, ::2]
    int_orm = np.stack([
        np.full((256, 256), 0.85, dtype=np.float32),
        np.full((256, 256), 0.55, dtype=np.float32),
        np.full((256, 256), 0.40, dtype=np.float32),
    ], axis=-1)
    save_map("Interior", "basecolor", int_bc, is_data=False, use_jpg=True)
    save_map("Interior", "normal", int_nor, is_data=True, use_jpg=True)
    save_map("Interior", "orm", int_orm, is_data=True, use_jpg=True)

    # --- 6. GLASS (256 RGBA PNG for alpha BLEND) ---
    glass_rgba = np.zeros((256, 256, 4), dtype=np.float32)
    glass_rgba[:, :, 0] = 0.65
    glass_rgba[:, :, 1] = 0.82
    glass_rgba[:, :, 2] = 0.90
    glass_rgba[:, :, 3] = 0.35
    flat_nor_256 = np.zeros((256, 256, 3), dtype=np.float32)
    flat_nor_256[:, :, 0] = 0.5
    flat_nor_256[:, :, 1] = 0.5
    flat_nor_256[:, :, 2] = 1.0
    glass_orm = np.stack([
        np.full((256, 256), 1.0, dtype=np.float32),
        np.full((256, 256), 0.05, dtype=np.float32),
        np.full((256, 256), 0.12, dtype=np.float32),
    ], axis=-1)
    save_map("Glass", "basecolor", glass_rgba, is_data=False, use_jpg=False)
    save_map("Glass", "normal", flat_nor_256, is_data=True, use_jpg=True)
    save_map("Glass", "orm", glass_orm, is_data=True, use_jpg=True)

    # --- 7. ENGINES (512 basecolor matching ENGINE GLOW swatch on sheet.jpg: white-cyan throat to deep cobalt rim) ---
    eng_bc = np.zeros((512, 512, 3), dtype=np.float32)
    # y512=0 is throat (bright cyan-white), y512=1 is exit lip (deep royal cobalt blue)
    t_v = y512
    striation = 0.90 + 0.10 * np.cos(x512 * 2.0 * math.pi * 32.0)
    r_ch = (0.55 * ((1.0 - t_v) ** 2.6) + 0.03 * t_v) * striation
    g_ch = (0.88 * ((1.0 - t_v) ** 1.45) + 0.14 * t_v) * striation
    b_ch = (1.00 * ((1.0 - t_v) ** 0.55) + 0.42 * (t_v ** 1.2)) * striation
    eng_bc[:, :, 0] = np.clip(r_ch, 0.02, 1.0)
    eng_bc[:, :, 1] = np.clip(g_ch, 0.10, 1.0)
    eng_bc[:, :, 2] = np.clip(b_ch, 0.32, 1.0)
    eng_orm = np.stack([
        np.full((256, 256), 1.0, dtype=np.float32),
        np.full((256, 256), 0.18, dtype=np.float32),
        np.full((256, 256), 0.25, dtype=np.float32),
    ], axis=-1)
    save_map("Engines", "basecolor", eng_bc, is_data=False, use_jpg=True)
    save_map("Engines", "normal", flat_nor_256, is_data=True, use_jpg=True)
    save_map("Engines", "orm", eng_orm, is_data=True, use_jpg=True)

    # --- 8. LAMPS (256 basecolor/normal/orm) ---
    lamp_bc = np.zeros((256, 256, 3), dtype=np.float32)
    lamp_bc[:, :, 0] = 0.88
    lamp_bc[:, :, 1] = 0.95
    lamp_bc[:, :, 2] = 1.00
    save_map("Lamps", "basecolor", lamp_bc, is_data=False, use_jpg=True)
    save_map("Lamps", "normal", flat_nor_256, is_data=True, use_jpg=True)
    save_map("Lamps", "orm", eng_orm, is_data=True, use_jpg=True)

    # Genuine Cycles GPU bake pass for Paint_basecolor.jpg
    bpy.ops.mesh.primitive_plane_add(size=2.0, location=(0, 0, -500))
    bake_plane = bpy.context.active_object
    bake_mat = bpy.data.materials.new("CyclesBakeTemp")
    bake_mat.use_nodes = True
    nt = bake_mat.node_tree
    nt.nodes.clear()
    out_n = nt.nodes.new("ShaderNodeOutputMaterial")
    emit_n = nt.nodes.new("ShaderNodeEmission")
    tex_n = nt.nodes.new("ShaderNodeTexImage")
    tex_n.image = load_image(TEX_FILES[("Paint", "basecolor")], non_color=False)
    nt.links.new(tex_n.outputs["Color"], emit_n.inputs["Color"])
    nt.links.new(emit_n.outputs["Emission"], out_n.inputs["Surface"])
    target_n = nt.nodes.new("ShaderNodeTexImage")
    target_img = bpy.data.images.new("CyclesBakeTarget", width=1024, height=1024, alpha=False)
    target_n.image = target_img
    for n in nt.nodes:
        n.select = False
    target_n.select = True
    nt.nodes.active = target_n
    bake_plane.data.materials.append(bake_mat)
    sc.cycles.samples = 4
    sc.cycles.bake_type = "EMIT"
    bpy.ops.object.bake(type="EMIT")
    sc.render.image_settings.file_format = "JPEG"
    sc.render.image_settings.quality = 92
    target_img.save_render(os.path.join(TEX_DIR, TEX_FILES[("Paint", "basecolor")]), scene=sc)
    bpy.data.objects.remove(bake_plane, do_unlink=True)
    bpy.data.materials.remove(bake_mat)
    bpy.data.images.remove(target_img)

generate_procedural_baked_maps()

# ---------------------------------------------------------------------------
# 4. Build Export-Ready PBR Materials
# ---------------------------------------------------------------------------
MATS = {}

def build_pbr_material(name, paint_factor=None, is_glass=False, emissive_rgb=None, emissive_strength=0.0):
    mat = bpy.data.materials.new(name)
    mat.use_nodes = True
    nt = mat.node_tree
    nt.nodes.clear()

    out_node = nt.nodes.new("ShaderNodeOutputMaterial")
    bsdf = nt.nodes.new("ShaderNodeBsdfPrincipled")
    nt.links.new(bsdf.outputs["BSDF"], out_node.inputs["Surface"])

    bc_tex = nt.nodes.new("ShaderNodeTexImage")
    bc_tex.name = "BaseColorTex"
    bc_tex.image = load_image(TEX_FILES[(name, "basecolor")], non_color=False)
    bc_tex.image.reload()

    if paint_factor is not None:
        mix_node = nt.nodes.new("ShaderNodeMix")
        mix_node.name = "PaintFactor"
        mix_node.data_type = "RGBA"
        mix_node.blend_type = "MULTIPLY"
        mix_node.inputs["Factor"].default_value = 1.0
        mix_node.inputs["B"].default_value = (*paint_factor, 1.0)
        nt.links.new(bc_tex.outputs["Color"], mix_node.inputs["A"])
        nt.links.new(mix_node.outputs["Result"], bsdf.inputs["Base Color"])
    else:
        nt.links.new(bc_tex.outputs["Color"], bsdf.inputs["Base Color"])

    if is_glass:
        mat.blend_method = "BLEND"
        mat.use_backface_culling = False
        bsdf.inputs["Alpha"].default_value = 0.35
        nt.links.new(bc_tex.outputs["Alpha"], bsdf.inputs["Alpha"])

    orm_tex = nt.nodes.new("ShaderNodeTexImage")
    orm_tex.name = "ORMTex"
    orm_tex.image = load_image(TEX_FILES[(name, "orm")], non_color=True)
    sep = nt.nodes.new("ShaderNodeSeparateColor")
    nt.links.new(orm_tex.outputs["Color"], sep.inputs["Color"])
    nt.links.new(sep.outputs["Green"], bsdf.inputs["Roughness"])
    nt.links.new(sep.outputs["Blue"], bsdf.inputs["Metallic"])

    gltf_settings = nt.nodes.new("ShaderNodeGroup")
    gltf_tree = bpy.data.node_groups.get("glTF Material Output")
    if gltf_tree is None:
        gltf_tree = bpy.data.node_groups.new("glTF Material Output", "ShaderNodeTree")
        gltf_tree.interface.new_socket(name="Occlusion", in_out="INPUT", socket_type="NodeSocketFloat")
    gltf_settings.node_tree = gltf_tree
    nt.links.new(sep.outputs["Red"], gltf_settings.inputs["Occlusion"])

    nor_tex = nt.nodes.new("ShaderNodeTexImage")
    nor_tex.name = "NormalTex"
    nor_tex.image = load_image(TEX_FILES[(name, "normal")], non_color=True)
    nor_map = nt.nodes.new("ShaderNodeNormalMap")
    nor_map.inputs["Strength"].default_value = 0.65
    nt.links.new(nor_tex.outputs["Color"], nor_map.inputs["Color"])
    nt.links.new(nor_map.outputs["Normal"], bsdf.inputs["Normal"])

    if emissive_rgb is not None:
        nt.links.new(bc_tex.outputs["Color"], bsdf.inputs["Emission Color"])
        bsdf.inputs["Emission Strength"].default_value = emissive_strength

    MATS[name] = mat
    return mat

build_pbr_material("Paint", paint_factor=PAINTS[0][1])
build_pbr_material("Metal")
build_pbr_material("Scorched")
build_pbr_material("Solar")
build_pbr_material("Interior")
build_pbr_material("Glass", is_glass=True)
build_pbr_material("Engines", emissive_rgb=(0.30, 0.80, 1.0), emissive_strength=2.1)
build_pbr_material("Lamps", emissive_rgb=(0.88, 0.95, 1.0), emissive_strength=2.5)

# ---------------------------------------------------------------------------
# 5. Geometry Helpers (glTF frame: +X left, +Y up, +Z nose -> Blender: X, -Z, Y)
# ---------------------------------------------------------------------------
def g2b(gx, gy, gz):
    return Vector((float(gx), float(-gz), float(gy)))

def g2b_vec(v):
    return Vector((float(v[0]), float(-v[2]), float(v[1])))

def new_mesh_part(name, mat_name, collection_list):
    mesh = bpy.data.meshes.new(name)
    obj = bpy.data.objects.new(name, mesh)
    bpy.context.collection.objects.link(obj)
    obj.data.materials.append(MATS[mat_name])
    collection_list.append(obj)
    return obj, mesh

def add_lathe_z(name, mat_name, z_r_profile, segs=48, center_gx=0.0, center_gy=0.0,
                cap_front=False, cap_back=False, smooth=True, parts_list=None,
                phi_start=0.0, phi_end=2.0 * math.pi, uv_v_from_profile=False,
                invert_uv_v=False):
    bm = bmesh.new()
    uv_lay = bm.loops.layers.uv.new("UVMap")
    rings = []
    full_circle = abs((phi_end - phi_start) - 2.0 * math.pi) < 1e-4
    n_pts = segs if full_circle else (segs + 1)

    for k, (gz, r) in enumerate(z_r_profile):
        ring = []
        for i in range(n_pts):
            t = i / float(segs)
            phi = phi_start + t * (phi_end - phi_start)
            gx = center_gx + r * math.cos(phi)
            gy = center_gy + r * math.sin(phi)
            ring.append(bm.verts.new(g2b(gx, gy, gz)))
        rings.append(ring)

    n_rings = max(1, len(rings) - 1)
    for k in range(len(rings) - 1):
        r0 = rings[k]
        r1 = rings[k + 1]
        v0_uv = k / float(n_rings)
        v1_uv = (k + 1) / float(n_rings)
        if invert_uv_v:
            v0_uv, v1_uv = 1.0 - v0_uv, 1.0 - v1_uv
        for i in range(segs):
            i_next = (i + 1) % segs if full_circle else (i + 1)
            f = bm.faces.new((r0[i], r0[i_next], r1[i_next], r1[i]))
            f.smooth = smooth
            if uv_v_from_profile:
                u0 = i / float(segs)
                u1 = (i + 1) / float(segs)
                f.loops[0][uv_lay].uv = (u0, v0_uv)
                f.loops[1][uv_lay].uv = (u1, v0_uv)
                f.loops[2][uv_lay].uv = (u1, v1_uv)
                f.loops[3][uv_lay].uv = (u0, v1_uv)

    if cap_back and full_circle and z_r_profile[0][1] > 1e-4:
        c_back = bm.verts.new(g2b(center_gx, center_gy, z_r_profile[0][0]))
        for i in range(segs):
            i_next = (i + 1) % segs
            f = bm.faces.new((rings[0][i_next], rings[0][i], c_back))
            f.smooth = False

    if cap_front and full_circle and z_r_profile[-1][1] > 1e-4:
        c_front = bm.verts.new(g2b(center_gx, center_gy, z_r_profile[-1][0]))
        for i in range(segs):
            i_next = (i + 1) % segs
            f = bm.faces.new((rings[-1][i], rings[-1][i_next], c_front))
            f.smooth = False

    bmesh.ops.recalc_face_normals(bm, faces=bm.faces)
    obj, mesh = new_mesh_part(name, mat_name, parts_list)
    if uv_v_from_profile:
        obj["preserve_uv"] = True
    bm.to_mesh(mesh)
    bm.free()
    return obj

def add_tapered_visor_shell(name, mat_name, z_sphere_c, r_shell, sign_y, parts_list):
    """Builds a smooth curved upper or lower observation dome cowl with tapered arch boundary."""
    bm = bmesh.new()
    n_rings = 16
    segs = 28
    rings = []
    edge_left = []
    edge_right = []
    for k in range(n_rings + 1):
        t = k / float(n_rings)
        gz = 66.8 + (73.4 - 66.8) * t
        dz = gz - z_sphere_c
        r = math.sqrt(max(0.4, r_shell * r_shell - dz * dz))
        # Smooth arch taper: from 10 deg..170 deg at Z=66.8 to 40 deg..140 deg at Z=73.4
        pad_deg = 10.0 + 30.0 * (t ** 1.35)
        if sign_y > 0:
            p0 = math.radians(pad_deg)
            p1 = math.radians(180.0 - pad_deg)
        else:
            p0 = math.radians(180.0 + pad_deg)
            p1 = math.radians(360.0 - pad_deg)
        ring = []
        for i in range(segs + 1):
            s = i / float(segs)
            phi = p0 + s * (p1 - p0)
            gx = r * math.cos(phi)
            gy = r * math.sin(phi)
            ring.append(bm.verts.new(g2b(gx, gy, gz)))
        rings.append(ring)
        edge_left.append((r * math.cos(p0), r * math.sin(p0), gz))
        edge_right.append((r * math.cos(p1), r * math.sin(p1), gz))

    for k in range(n_rings):
        r0, r1 = rings[k], rings[k + 1]
        for i in range(segs):
            f = bm.faces.new((r0[i], r0[i + 1], r1[i + 1], r1[i]))
            f.smooth = True

    bmesh.ops.recalc_face_normals(bm, faces=bm.faces)
    obj, mesh = new_mesh_part(name, mat_name, parts_list)
    bm.to_mesh(mesh)
    bm.free()
    return edge_left, edge_right

def add_lathe_y(name, mat_name, y_r_profile, segs=36, center_gx=0.0, center_gz=0.0,
                cap_top=False, smooth=True, parts_list=None):
    bm = bmesh.new()
    rings = []
    for gy, r in y_r_profile:
        ring = []
        for i in range(segs):
            phi = 2.0 * math.pi * i / float(segs)
            gx = center_gx + r * math.cos(phi)
            gz = center_gz + r * math.sin(phi)
            ring.append(bm.verts.new(g2b(gx, gy, gz)))
        rings.append(ring)
    for k in range(len(rings) - 1):
        r0, r1 = rings[k], rings[k + 1]
        for i in range(segs):
            inext = (i + 1) % segs
            f = bm.faces.new((r0[i], r0[inext], r1[inext], r1[i]))
            f.smooth = smooth
    if cap_top and y_r_profile[-1][1] > 1e-4:
        c_top = bm.verts.new(g2b(center_gx, y_r_profile[-1][0], center_gz))
        for i in range(segs):
            inext = (i + 1) % segs
            f = bm.faces.new((rings[-1][i], rings[-1][inext], c_top))
            f.smooth = True
    bmesh.ops.recalc_face_normals(bm, faces=bm.faces)
    obj, mesh = new_mesh_part(name, mat_name, parts_list)
    bm.to_mesh(mesh)
    bm.free()
    return obj

def add_lathe_x(name, mat_name, x_r_profile, segs=36, center_gy=0.0, center_gz=0.0,
                cap_end=False, smooth=True, parts_list=None):
    bm = bmesh.new()
    rings = []
    for gx, r in x_r_profile:
        ring = []
        for i in range(segs):
            phi = 2.0 * math.pi * i / float(segs)
            gy = center_gy + r * math.cos(phi)
            gz = center_gz + r * math.sin(phi)
            ring.append(bm.verts.new(g2b(gx, gy, gz)))
        rings.append(ring)
    for k in range(len(rings) - 1):
        r0, r1 = rings[k], rings[k + 1]
        for i in range(segs):
            inext = (i + 1) % segs
            f = bm.faces.new((r0[i], r0[inext], r1[inext], r1[i]))
            f.smooth = smooth
    if cap_end and x_r_profile[-1][1] > 1e-4:
        c_end = bm.verts.new(g2b(x_r_profile[-1][0], center_gy, center_gz))
        for i in range(segs):
            inext = (i + 1) % segs
            f = bm.faces.new((rings[-1][i], rings[-1][inext], c_end))
            f.smooth = False
    bmesh.ops.recalc_face_normals(bm, faces=bm.faces)
    obj, mesh = new_mesh_part(name, mat_name, parts_list)
    bm.to_mesh(mesh)
    bm.free()
    return obj

def add_strut(bm, p0_gltf, p1_gltf, r0, r1=None, segs=10, cap=True, smooth=True):
    if r1 is None:
        r1 = r0
    v0 = g2b_vec(p0_gltf)
    v1 = g2b_vec(p1_gltf)
    axis = v1 - v0
    length = axis.length
    if length < 1e-5:
        return
    d = axis.normalized()
    up = Vector((0, 0, 1)) if abs(d.z) < 0.9 else Vector((1, 0, 0))
    u = d.cross(up).normalized()
    w = d.cross(u).normalized()

    ring0 = []
    ring1 = []
    for i in range(segs):
        ang = 2.0 * math.pi * i / segs
        dir_vec = u * math.cos(ang) + w * math.sin(ang)
        ring0.append(bm.verts.new(v0 + dir_vec * r0))
        ring1.append(bm.verts.new(v1 + dir_vec * r1))

    for i in range(segs):
        inext = (i + 1) % segs
        f = bm.faces.new((ring0[i], ring0[inext], ring1[inext], ring1[i]))
        f.smooth = smooth

    if cap:
        f0 = bm.faces.new(ring0[::-1])
        f0.smooth = False
        f1 = bm.faces.new(ring1)
        f1.smooth = False

def add_bevelled_box(bm, center_gltf, size_gltf, bevel=0.06, rot_gltf_euler=(0, 0, 0)):
    sx, sy, sz = size_gltf[0] * 0.5, size_gltf[1] * 0.5, size_gltf[2] * 0.5
    b = min(bevel, sx * 0.35, sy * 0.35, sz * 0.35)
    tmp = bmesh.new()
    bmesh.ops.create_cube(tmp, size=1.0)
    for v in tmp.verts:
        v.co.x *= size_gltf[0]
        v.co.y *= size_gltf[1]
        v.co.z *= size_gltf[2]
    if b > 1e-4:
        bmesh.ops.bevel(tmp, geom=list(tmp.edges), offset=b, segments=2, profile=0.5, affect="EDGES")
    rot_m = Euler(rot_gltf_euler, "XYZ").to_matrix()
    c_vec = Vector(center_gltf)
    v_map = {}
    for v in tmp.verts:
        p_g = c_vec + rot_m @ v.co
        v_map[v] = bm.verts.new(g2b_vec(p_g))
    for f in tmp.faces:
        nf = bm.faces.new([v_map[v] for v in f.verts])
        nf.smooth = False
    tmp.free()

def finalize_bm_part(bm, name, mat_name, parts_list):
    bmesh.ops.recalc_face_normals(bm, faces=bm.faces)
    obj, mesh = new_mesh_part(name, mat_name, parts_list)
    bm.to_mesh(mesh)
    bm.free()
    return obj

# ---------------------------------------------------------------------------
# 6. Build HULL Fixed Assemblies
# ---------------------------------------------------------------------------
hull_parts = []

# === A. FORWARD COMMAND MODULE & OBSERVATION DOME (Z = +49.0m to +80.0m) ===
z_sphere_c = 63.0
r_sphere = 11.45

cmd_profile = [
    (49.0, 3.85),
    (49.3, 4.35),
    (50.3, 4.35),
    (50.6, 4.10),
    (51.5, 4.75),
    (52.3, 5.65),
]
# Smooth spherical pressure hull extending to Z=67.0 with subtle recessed panel grooves
n_lat = 30
for i in range(n_lat + 1):
    gz = 52.5 + (67.0 - 52.5) * (i / float(n_lat))
    dz = gz - z_sphere_c
    base_r = math.sqrt(max(0.1, r_sphere * r_sphere - dz * dz))
    groove = -0.08 if (i in (7, 15, 23)) else 0.0
    cmd_profile.append((gz, base_r + groove))
cmd_profile.extend([
    (67.05, base_r - 0.18),
    (66.75, base_r - 0.38),
])
add_lathe_z("Cmd_MainSphere", "Paint", cmd_profile, segs=56, cap_back=True, smooth=True, parts_list=hull_parts)

# Sleek Low-Profile Dorsal (+Y) and Ventral (-Y) Armored Cap Domes hugging the sphere contour
for sy, tag in ((1.0, "Dorsal"), (-1.0, "Ventral")):
    cap_prof = [
        (sy * 10.15, 5.20),
        (sy * 10.65, 4.85),
        (sy * 11.15, 3.85),
        (sy * 11.45, 2.35),
        (sy * 11.58, 0.05),
    ]
    if sy < 0:
        cap_prof.reverse()
    add_lathe_y(f"Cmd_CapDome_{tag}", "Paint", cap_prof, segs=40, center_gx=0.0, center_gz=62.5,
                cap_top=(sy > 0), smooth=True, parts_list=hull_parts)

# Port (+X) & Starboard (-X) Stepped Circular Airlock / Docking Collars (Paint outer + Metal inner ring)
for sx, tag in ((1.0, "Port"), (-1.0, "Stbd")):
    dock_paint_prof = [
        (sx * 9.8, 4.25),
        (sx * 11.8, 4.05),
        (sx * 12.45, 3.65),
        (sx * 12.65, 3.20),
        (sx * 12.30, 3.05),
        (sx * 12.30, 0.05),
    ]
    if sx < 0:
        dock_paint_prof.reverse()
    add_lathe_x(f"Cmd_DockPaint_{tag}", "Paint", dock_paint_prof, segs=36, center_gy=-0.4, center_gz=62.5,
                cap_end=(sx > 0), smooth=True, parts_list=hull_parts)

    dock_metal_prof = [
        (sx * 12.35, 3.22),
        (sx * 12.82, 3.05),
        (sx * 12.82, 2.65),
        (sx * 12.32, 2.55),
        (sx * 12.42, 1.45),
        (sx * 12.42, 0.95),
    ]
    if sx < 0:
        dock_metal_prof.reverse()
    add_lathe_x(f"Cmd_DockMetal_{tag}", "Metal", dock_metal_prof, segs=36, center_gy=-0.4, center_gz=62.5,
                cap_end=False, smooth=True, parts_list=hull_parts)

# Smooth Tapered Upper & Lower Visor Armor Cowls over the Observation Dome
up_l, up_r = add_tapered_visor_shell("Cmd_Visor_Upper", "Paint", z_sphere_c, r_sphere + 0.12, 1.0, hull_parts)
lo_l, lo_r = add_tapered_visor_shell("Cmd_Visor_Lower", "Paint", z_sphere_c, r_sphere + 0.12, -1.0, hull_parts)

# Observation Dome Glazing (Tinted Glass shell from Z=66.6 to 74.45)
glass_prof = []
for k in range(18):
    gz = 66.6 + (74.42 - 66.6) * (k / 17.0)
    dz = gz - z_sphere_c
    r = math.sqrt(max(0.08, (r_sphere - 0.06) ** 2 - dz * dz))
    glass_prof.append((gz, r))
add_lathe_z("Cmd_ObservationGlass", "Glass", glass_prof, segs=48, cap_front=True, smooth=True, parts_list=hull_parts)

bm_cmd_paint = bmesh.new()
bm_cmd_metal = bmesh.new()
bm_cmd_lamps = bmesh.new()

# Frame the tapered visor boundary edges smoothly in Paint + Metal
for edge_curve in (up_l, up_r, lo_l, lo_r):
    for k in range(len(edge_curve) - 1):
        add_strut(bm_cmd_paint, edge_curve[k], edge_curve[k + 1], 0.24, 0.24, segs=8, cap=False)

# Sculpted White Equatorial Belt across the Observation Dome (matching sheet.jpg!)
eq_pts = []
for k in range(25):
    ang = math.radians(-95.0 + 190.0 * (k / 24.0))
    gx = (r_sphere + 0.16) * math.sin(ang)
    gz = z_sphere_c + (r_sphere + 0.16) * math.cos(ang)
    eq_pts.append((gx, -0.40, gz))
for k in range(len(eq_pts) - 1):
    add_strut(bm_cmd_paint, eq_pts[k], eq_pts[k + 1], 0.44, 0.44, segs=8, cap=True)

# 3 Primary Sculpted White Vertical/Angled Mullions + Thin Inner Metal Window Frames
for az_deg in (-34.0, 0.0, 34.0):
    az = math.radians(az_deg)
    m_pts_p = []
    m_pts_m = []
    for k in range(17):
        el = math.radians(-38.0 + 76.0 * (k / 16.0))
        for r_off, lst in ((0.08, m_pts_p), (-0.02, m_pts_m)):
            r_m = r_sphere + r_off
            gy = r_m * math.sin(el)
            r_xz = r_m * math.cos(el)
            lst.append((r_xz * math.sin(az), gy, z_sphere_c + r_xz * math.cos(az)))
    for k in range(len(m_pts_p) - 1):
        add_strut(bm_cmd_paint, m_pts_p[k], m_pts_p[k + 1], 0.22, 0.22, segs=8, cap=False)
        add_strut(bm_cmd_metal, m_pts_m[k], m_pts_m[k + 1], 0.13, 0.13, segs=6, cap=False)

# Upper and lower sculpted window sill arches at Z = 73.2m
for el_deg in (-26.0, 26.0):
    el = math.radians(el_deg)
    gy = (r_sphere + 0.10) * math.sin(el)
    r_xz = (r_sphere + 0.10) * math.cos(el)
    arch_pts = []
    for k in range(17):
        az = math.radians(-50.0 + 100.0 * (k / 16.0))
        arch_pts.append((r_xz * math.sin(az), gy, z_sphere_c + r_xz * math.cos(az)))
    for k in range(len(arch_pts) - 1):
        add_strut(bm_cmd_paint, arch_pts[k], arch_pts[k + 1], 0.24, 0.24, segs=8, cap=False)

# Subtle White Meridian Seam Straps on the Command Sphere
for i_mer in range(8):
    phi = i_mer * (2.0 * math.pi / 8.0) + math.pi / 8.0
    mer_pts = []
    for k in range(14):
        gz = 52.0 + (66.5 - 52.0) * (k / 13.0)
        dz = gz - z_sphere_c
        r_m = math.sqrt(max(0.5, (r_sphere + 0.07) ** 2 - dz * dz))
        mer_pts.append((r_m * math.cos(phi), r_m * math.sin(phi), gz))
    for k in range(len(mer_pts) - 1):
        add_strut(bm_cmd_paint, mer_pts[k], mer_pts[k + 1], 0.13, 0.13, segs=6, cap=True)

# Aft neck gusset fins (12 radial fins around Z=49.2 to 53.0)
for i_g in range(12):
    phi = i_g * (2.0 * math.pi / 12.0)
    cx, cy = math.cos(phi), math.sin(phi)
    add_strut(bm_cmd_metal, (cx * 4.1, cy * 4.1, 49.2), (cx * 6.2, cy * 6.2, 52.8), 0.16, 0.22, segs=6)

# Chin Sensor Fairing & Prominent Twin Stepped Nose Probes extending to EXACT Z = +80.0m
add_bevelled_box(bm_cmd_paint, (0.0, -6.5, 71.0), (3.6, 1.7, 5.2), bevel=0.32)
add_bevelled_box(bm_cmd_metal, (0.0, -6.5, 73.7), (2.6, 1.0, 2.0), bevel=0.16)
for sx in (-1.15, 1.15):
    add_strut(bm_cmd_metal, (sx, -6.55, 73.4), (sx, -6.55, 76.6), 0.28, 0.20, segs=10)
    add_strut(bm_cmd_metal, (sx, -6.55, 76.6), (sx, -6.55, 79.0), 0.16, 0.10, segs=8)
    add_strut(bm_cmd_metal, (sx, -6.55, 79.0), (sx, -6.55, 80.0), 0.08, 0.04, segs=8)
add_strut(bm_cmd_metal, (0.0, -6.15, 73.6), (0.0, -6.15, 78.5), 0.16, 0.06, segs=8)

# 4 RCS Maneuvering Thruster Quad Pods on Command Sphere aft shoulders (Z = +55.5m)
for sx in (-1.0, 1.0):
    for sy in (-1.0, 1.0):
        rx, ry, rz = sx * 6.45, sy * 6.45, 55.5
        add_bevelled_box(bm_cmd_paint, (rx, ry, rz), (1.2, 1.2, 1.7), bevel=0.16)
        for dx, dy in ((-0.42, 0.0), (0.42, 0.0), (0.0, -0.42), (0.0, 0.42)):
            add_strut(bm_cmd_metal, (rx + dx, ry + dy, rz), (rx + dx * 1.7, ry + dy * 1.7, rz), 0.12, 0.20, segs=8)

# Subtle Human-Scale Recessed Crew Portholes along Command Module mid-belt (Z = 57.5m)
dz_w = 57.5 - z_sphere_c
rw_cmd = math.sqrt((r_sphere + 0.02) ** 2 - dz_w * dz_w)
for deg_w in (-155, -140, -125, -55, -40, -25, 25, 40, 55, 125, 140, 155):
    phi = math.radians(deg_w)
    add_bevelled_box(bm_cmd_lamps, (rw_cmd * math.cos(phi), rw_cmd * math.sin(phi), 57.5),
                     (0.22, 0.22, 0.38), bevel=0.02)

# === B. BRIDGE & OBSERVATION DECK INTERIOR (behind Glass, Z = +61.5m to +73.5m) ===
bm_interior = bmesh.new()
add_bevelled_box(bm_interior, (0.0, -1.35, 67.5), (15.5, 0.45, 11.2), bevel=0.08)
add_bevelled_box(bm_interior, (0.0, -4.35, 67.0), (13.0, 0.45, 9.6), bevel=0.08)
add_bevelled_box(bm_interior, (0.0, 2.35, 64.8), (14.0, 0.40, 6.2), bevel=0.08)
add_bevelled_box(bm_interior, (0.0, 0.0, 61.6), (17.5, 16.5, 0.55), bevel=0.10)
add_strut(bm_interior, (0.0, -1.15, 66.5), (0.0, -0.35, 66.5), 1.45, 1.65, segs=16)
for i_c in range(-2, 3):
    ang = math.radians(i_c * 18.0)
    cx = 6.2 * math.sin(ang)
    cz = 64.2 + 6.2 * math.cos(ang)
    add_bevelled_box(bm_interior, (cx, -0.72, cz), (1.35, 0.85, 0.95), bevel=0.08, rot_gltf_euler=(0.0, ang, 0.0))
for sx_seat, sz_seat in ((0.0, 69.2), (-2.1, 68.7), (2.1, 68.7)):
    add_bevelled_box(bm_interior, (sx_seat, -0.75, sz_seat), (0.75, 0.20, 0.75), bevel=0.04)
    add_bevelled_box(bm_interior, (sx_seat, -0.20, sz_seat - 0.30), (0.72, 1.05, 0.18), bevel=0.04)
for az_rib in (-40.0, 0.0, 40.0):
    az = math.radians(az_rib)
    add_strut(bm_interior, (3.0 * math.sin(az), 6.5, 62.5), (6.5 * math.sin(az), 4.2, 70.5), 0.22, 0.18, segs=8)

finalize_bm_part(bm_cmd_paint, "Cmd_Details_Paint", "Paint", hull_parts)
finalize_bm_part(bm_cmd_metal, "Cmd_Details_Metal", "Metal", hull_parts)
finalize_bm_part(bm_cmd_lamps, "Cmd_Details_Lamps", "Lamps", hull_parts)
finalize_bm_part(bm_interior, "Cmd_Bridge_Interior", "Interior", hull_parts)

# === C. LONG SLENDER OPEN-TRUSS SPINE (Z = +3.5m to +49.2m) ===
spine_core_prof = []
for k in range(38):
    gz = 3.2 + (49.2 - 3.2) * (k / 37.0)
    r_c = 1.58 if (k % 3 == 0) else 1.38
    spine_core_prof.append((gz, r_c))
add_lathe_z("Spine_CentralTunnel", "Paint", spine_core_prof, segs=24, smooth=True, parts_list=hull_parts)

bm_spine_metal = bmesh.new()
bm_spine_paint = bmesh.new()

z_truss_start = 3.5
z_truss_end = 49.2
corners = [(-3.05, -3.05), (3.05, -3.05), (3.05, 3.05), (-3.05, 3.05)]
for cx, cy in corners:
    add_strut(bm_spine_paint, (cx, cy, z_truss_start), (cx, cy, z_truss_end), 0.36, 0.36, segs=10)

mid_rails = [(-3.35, 0.0), (3.35, 0.0), (0.0, -3.35), (0.0, 3.35)]
for mx, my in mid_rails:
    add_strut(bm_spine_metal, (mx, my, z_truss_start), (mx, my, z_truss_end), 0.18, 0.18, segs=8)

for px, py in [(-1.95, -1.95), (1.95, -1.95), (1.95, 1.95), (-1.95, 1.95)]:
    add_strut(bm_spine_metal, (px, py, -9.5), (px, py, 49.2), 0.22, 0.22, segs=8)

n_bays = 8
bay_zs = [z_truss_start + (z_truss_end - z_truss_start) * (i / float(n_bays)) for i in range(n_bays + 1)]

for i_b, bz in enumerate(bay_zs):
    oct_pts = [
        (-3.05, -3.05, bz), (0.0, -3.35, bz), (3.05, -3.05, bz), (3.35, 0.0, bz),
        (3.05, 3.05, bz), (0.0, 3.35, bz), (-3.05, 3.05, bz), (-3.35, 0.0, bz)
    ]
    for k in range(8):
        p0 = oct_pts[k]
        p1 = oct_pts[(k + 1) % 8]
        add_strut(bm_spine_paint, p0, p1, 0.32, 0.32, segs=8)
        if k % 2 == 0:
            add_strut(bm_spine_metal, (0.0, 0.0, bz), p0, 0.18, 0.18, segs=6)

    if i_b < n_bays:
        bz_next = bay_zs[i_b + 1]
        bz_mid = 0.5 * (bz + bz_next)
        for k in range(4):
            c0 = corners[k]
            c1 = corners[(k + 1) % 4]
            add_strut(bm_spine_metal, (c0[0], c0[1], bz), (c1[0], c1[1], bz_next), 0.20, 0.20, segs=8)
            add_strut(bm_spine_metal, (c1[0], c1[1], bz), (c0[0], c0[1], bz_next), 0.20, 0.20, segs=8)
            add_strut(bm_spine_metal, (c0[0], c0[1], bz_mid), (c1[0], c1[1], bz_mid), 0.15, 0.15, segs=6)

        if i_b % 2 == 1:
            for sx_t, sy_t in ((-2.05, 0.0), (2.05, 0.0), (0.0, -2.05), (0.0, 2.05)):
                add_bevelled_box(bm_spine_paint, (sx_t, sy_t, bz_mid), (1.30, 1.30, 2.35), bevel=0.42)

finalize_bm_part(bm_spine_metal, "Spine_Truss_Metal", "Metal", hull_parts)
finalize_bm_part(bm_spine_paint, "Spine_Frames_Paint", "Paint", hull_parts)

# === D. CENTRAL BEARING HUB & PINCHED NECK BEHIND RING (Z = -9.5m to +3.5m) ===
hub_stator_prof = [
    (-9.5, 3.60),
    (-6.0, 3.60),
    (-5.6, 4.95),
    (-4.3, 4.95),
    (-4.1, 4.45),
    (1.1, 4.45),
    (1.3, 4.95),
    (3.2, 4.95),
    (3.6, 3.75),
]
add_lathe_z("Hub_StatorCore", "Paint", hub_stator_prof, segs=40, smooth=True, parts_list=hull_parts)

bm_hub_neck = bmesh.new()
for i_n in range(8):
    a0 = i_n * (2.0 * math.pi / 8.0)
    a1 = (i_n + 1) * (2.0 * math.pi / 8.0)
    add_strut(bm_hub_neck, (4.3 * math.cos(a0), 4.3 * math.sin(a0), -5.7),
              (5.1 * math.cos(a1), 5.1 * math.sin(a1), -9.6), 0.22, 0.22, segs=8)
    add_strut(bm_hub_neck, (4.3 * math.cos(a1), 4.3 * math.sin(a1), -5.7),
              (5.1 * math.cos(a0), 5.1 * math.sin(a0), -9.6), 0.22, 0.22, segs=8)
finalize_bm_part(bm_hub_neck, "Hub_NeckTruss", "Metal", hull_parts)

# === E. MID-AFT SCIENCE & ENGINEERING CORE (Z = -9.5m to -40.0m) ===
eng_core_prof = [
    (-40.0, 7.05),
    (-38.5, 7.20),
    (-35.6, 7.20),
    (-35.3, 7.45),
    (-31.2, 7.45),
    (-30.9, 7.12),
    (-29.5, 6.55),
    (-26.5, 6.55),
    (-26.2, 6.70),
    (-23.0, 6.70),
    (-22.6, 6.82),
    (-11.5, 6.82),
    (-11.0, 6.38),
    (-9.5, 4.80),
]
add_lathe_z("EngCore_MainCylinder", "Paint", eng_core_prof, segs=48, smooth=True, parts_list=hull_parts)

bm_eng_metal = bmesh.new()
bm_eng_paint = bmesh.new()
bm_eng_lamps = bmesh.new()

for gz_r, rad_r in ((-12.2, 6.96), (-16.5, 7.02), (-20.8, 6.96), (-24.6, 6.82), (-28.0, 6.68), (-31.5, 7.58), (-35.0, 7.58), (-38.2, 7.34)):
    n_seg = 40
    for k in range(n_seg):
        a0 = 2.0 * math.pi * k / n_seg
        a1 = 2.0 * math.pi * (k + 1) / n_seg
        add_strut(bm_eng_paint, (rad_r * math.cos(a0), rad_r * math.sin(a0), gz_r),
                  (rad_r * math.cos(a1), rad_r * math.sin(a1), gz_r), 0.22, 0.22, segs=6, cap=False)

for ang_deg in (38.0, 142.0, 218.0, 322.0):
    ang = math.radians(ang_deg)
    cx, cy = math.cos(ang), math.sin(ang)
    add_strut(bm_eng_paint, (cx * 6.90, cy * 6.90, -39.5), (cx * 6.90, cy * 6.90, -11.2), 0.26, 0.26, segs=8)
    add_strut(bm_eng_metal, (cx * 7.08 - cy * 0.42, cy * 7.08 + cx * 0.42, -39.0),
              (cx * 7.08 - cy * 0.42, cy * 7.08 + cx * 0.42, -11.5), 0.14, 0.14, segs=8)

for sx in (-1.0, 1.0):
    add_bevelled_box(bm_eng_paint, (sx * 7.12, 0.0, -33.2), (1.3, 3.6, 5.6), bevel=0.24)
    add_bevelled_box(bm_eng_metal, (sx * 7.55, 0.0, -33.2), (0.4, 1.8, 3.8), bevel=0.10)
    add_bevelled_box(bm_eng_metal, (sx * 6.60, -2.5, -13.0), (1.1, 1.3, 2.5), bevel=0.14)

for gz_w in (-14.2, -18.8, -24.2, -25.4, -26.6):
    for sx in (-1.0, 1.0):
        for dy in (-1.2, 0.0, 1.2):
            add_bevelled_box(bm_eng_lamps, (sx * 6.84, dy, gz_w), (0.16, 0.24, 0.42), bevel=0.02)

# Dorsal Pedestal Pylon for the Dish Antenna (mounted at Z = -28.0m, rises to pivot at (0, 10.6, -27.2))
add_bevelled_box(bm_eng_paint, (0.0, 8.2, -28.5), (3.2, 2.8, 5.2), bevel=0.25)
add_strut(bm_eng_metal, (0.0, 8.8, -27.2), (0.0, 10.6, -27.2), 1.25, 0.98, segs=16)
for sx in (-1.15, 1.15):
    add_strut(bm_eng_metal, (sx, 6.8, -30.5), (sx * 0.65, 10.0, -27.2), 0.24, 0.18, segs=8)

# Ventral Sensor & Turret Mount Housing at (0, -7.4, -33.5)
add_bevelled_box(bm_eng_paint, (0.0, -7.4, -33.5), (2.8, 2.2, 3.8), bevel=0.22)

# === F. DEEP-SPACE SENSOR BOOMS (matching sheet.jpg SIDE, TOP, and HERO views!) ===
# Two pairs of upper aft-swept booms behind the dish pedestal + one pair of lower aft-swept booms!
sensor_booms = [
    # Upper-high pair sweeping aft & upward behind the dish antenna (cleanly below dish rim)
    ((-3.2, 6.8, -31.8), (-11.5, 15.2, -41.5)),
    (( 3.2, 6.8, -31.8), ( 11.5, 15.2, -41.5)),
    # Upper-mid pair at shallower angle (matching the second upper boom leader in sheet.jpg!)
    ((-4.5, 5.2, -33.0), (-12.8, 10.8, -42.8)),
    (( 4.5, 5.2, -33.0), ( 12.8, 10.8, -42.8)),
    # Lower pair sweeping aft & downward from the ventral sensor pod
    ((-2.6, -7.2, -34.2), (-11.8, -16.2, -43.8)),
    (( 2.6, -7.2, -34.2), ( 11.8, -16.2, -43.8)),
]
for b_start, b_tip in sensor_booms:
    v_s = Vector(b_start)
    v_t = Vector(b_tip)
    v_m1 = v_s.lerp(v_t, 0.38)
    v_m2 = v_s.lerp(v_t, 0.80)
    add_strut(bm_eng_metal, v_s, v_m1, 0.32, 0.20, segs=10)
    add_strut(bm_eng_paint, v_s.lerp(v_t, 0.32), v_s.lerp(v_t, 0.46), 0.38, 0.36, segs=10)
    add_strut(bm_eng_metal, v_m1, v_m2, 0.17, 0.11, segs=8)
    add_strut(bm_eng_paint, v_s.lerp(v_t, 0.76), v_s.lerp(v_t, 0.88), 0.30, 0.24, segs=10)
    add_strut(bm_eng_metal, v_m2, v_t, 0.09, 0.05, segs=8)
    d_boom = (v_t - v_s).normalized()
    for off in (Vector((0.45, 0.22, 0.0)), Vector((-0.45, 0.22, 0.0)), Vector((0.0, -0.45, 0.0))):
        add_strut(bm_eng_metal, v_t - d_boom * 0.65, v_t + d_boom * 0.55 + off, 0.04, 0.02, segs=6)

# === G. SOLAR & RADIATOR ARRAYS (70.0m Height Y=[-35,+35] & 70.0m Width X=[-35,+35]) ===
bm_solar = bmesh.new()

for sy in (1.0, -1.0):
    y_root = sy * 6.8
    y_start = sy * 7.8
    y_end = sy * 35.0
    z_center = -16.5

    add_bevelled_box(bm_eng_paint, (0.0, sy * 7.2, z_center), (2.4, 1.4, 3.6), bevel=0.18)
    add_strut(bm_eng_metal, (0.0, y_root, z_center), (0.0, y_end, z_center), 0.34, 0.22, segs=10)
    add_bevelled_box(bm_eng_metal, (0.0, sy * 34.85, z_center), (0.60, 0.30, 12.2), bevel=0.05)
    add_bevelled_box(bm_eng_metal, (0.0, sy * 7.95, z_center), (0.60, 0.35, 12.2), bevel=0.05)

    for z_bay_min, z_bay_max in ((-16.1, -10.6), (-22.4, -16.9)):
        z_bay_c = 0.5 * (z_bay_min + z_bay_max)
        bay_w = z_bay_max - z_bay_min
        y_mid = 0.5 * (y_start + (y_end - sy * 0.15))
        bay_h = abs((y_end - sy * 0.15) - y_start)

        add_bevelled_box(bm_solar, (0.0, y_mid, z_bay_c), (0.20, bay_h - 0.2, bay_w - 0.2), bevel=0.03)

        for j_v in range(5):
            zv = z_bay_min + bay_w * (j_v / 4.0)
            r_frame = 0.20 if j_v in (0, 4) else 0.12
            add_strut(bm_eng_metal, (0.0, y_start, zv), (0.0, y_end - sy * 0.15, zv), r_frame, r_frame, segs=8)
        for j_h in range(9):
            yh = y_start + (y_end - sy * 0.15 - y_start) * (j_h / 8.0)
            r_frame = 0.18 if j_h in (0, 8) else 0.11
            add_strut(bm_eng_metal, (0.0, yh, z_bay_min), (0.0, yh, z_bay_max), r_frame, r_frame, segs=8)

for sx in (1.0, -1.0):
    x_root = sx * 6.8
    x_end = sx * 35.0
    z_center = -16.5

    add_bevelled_box(bm_eng_paint, (sx * 7.2, 0.0, z_center), (1.4, 2.0, 2.8), bevel=0.16)
    add_strut(bm_eng_metal, (x_root, 0.0, z_center), (x_end - sx * 0.20, 0.0, z_center), 0.34, 0.18, segs=10)
    add_strut(bm_eng_paint, (sx * 7.5, 0.0, z_center - 0.85), (sx * 33.8, 0.0, z_center - 0.45), 0.16, 0.10, segs=8)
    add_strut(bm_eng_paint, (sx * 7.5, 0.0, z_center + 0.85), (sx * 33.8, 0.0, z_center + 0.45), 0.16, 0.10, segs=8)
    for k_b in range(1, 9):
        xb = sx * (7.5 + (33.8 - 7.5) * (k_b / 8.0))
        zb_span = 0.85 - 0.40 * (k_b / 8.0)
        add_strut(bm_eng_metal, (xb, 0.0, z_center - zb_span), (xb, 0.0, z_center + zb_span), 0.10, 0.10, segs=6)
    add_bevelled_box(bm_eng_metal, (sx * 34.65, 0.0, z_center), (0.70, 0.55, 2.2), bevel=0.08)
    add_bevelled_box(bm_eng_lamps, (sx * 34.92, 0.0, z_center), (0.16, 0.24, 0.45), bevel=0.02)

finalize_bm_part(bm_solar, "Arrays_SolarCells", "Solar", hull_parts)
finalize_bm_part(bm_eng_metal, "EngCore_Details_Metal", "Metal", hull_parts)
finalize_bm_part(bm_eng_paint, "EngCore_Details_Paint", "Paint", hull_parts)
finalize_bm_part(bm_eng_lamps, "EngCore_Details_Lamps", "Lamps", hull_parts)

# === H. REAR FUSION DRIVE SECTION & MAGNETIC NOZZLE (Z = -40.0m to -80.0m) ===
louver_core_prof = [
    (-46.0, 9.15),
    (-40.0, 6.85),
]
add_lathe_z("Drive_LouverCore", "Scorched", louver_core_prof, segs=48, smooth=True, parts_list=hull_parts)

drive_drum_prof = [
    (-62.8, 9.05),
    (-62.2, 9.48),
    (-50.2, 9.48),
]
for k_r in range(11):
    gz_r = -50.0 + k_r * 0.42
    rad_r = 9.48 if (k_r % 2 == 0) else 9.32
    drive_drum_prof.append((gz_r, rad_r))
add_lathe_z("Drive_MainDrum", "Paint", drive_drum_prof, segs=56, smooth=True, parts_list=hull_parts)

blue_band_prof = [
    (-61.2, 9.52),
    (-61.0, 9.56),
    (-59.2, 9.56),
    (-59.0, 9.52),
]
add_lathe_z("Drive_BlueAccentBand", "Solar", blue_band_prof, segs=56, smooth=True, parts_list=hull_parts)

pinch_core_prof = [
    (-69.2, 6.85),
    (-67.2, 7.35),
    (-64.5, 8.20),
    (-62.8, 9.05),
]
add_lathe_z("Drive_PinchCore", "Scorched", pinch_core_prof, segs=48, smooth=True, parts_list=hull_parts)

nozzle_throat_prof = [
    (-69.6, 7.12),
    (-69.0, 6.95),
    (-68.4, 6.88),
]
add_lathe_z("Drive_NozzleThroat", "Scorched", nozzle_throat_prof, segs=56, smooth=True, parts_list=hull_parts)

nozzle_lip_prof = [
    (-80.0, 8.48),
    (-80.0, 8.62),
    (-79.3, 8.54),
    (-79.3, 8.40),
]
add_lathe_z("Drive_NozzleLip", "Scorched", nozzle_lip_prof, segs=56, smooth=True, parts_list=hull_parts)

# Outer & Inner Glowing Magnetic Plasma Bell Skirt (Z = -69.4m to -79.4m, Engines material)
nozzle_glow_outer = []
nozzle_glow_inner = []
for k in range(21):
    t = k / 20.0
    gz = -69.4 - t * (79.4 - 69.4)
    r_bell = 6.95 + 1.52 * (1.0 - (1.0 - t) ** 1.60)
    nozzle_glow_outer.append((gz, r_bell))
    nozzle_glow_inner.append((gz, r_bell - 0.22))

# Note: nozzle_glow_outer_rev starts at -79.4 (exit lip) and ends at -69.4 (throat).
# In Blender, UV v=0 is the bottom of the image (y512=0, bright cyan-white throat) and v=1 is top (y512=1, deep cobalt lip).
# So for nozzle_glow_outer_rev, invert_uv_v=True maps -69.4 (throat) to v=0 (bright cyan-white) and -79.4 (lip) to v=1 (cobalt blue)!
nozzle_glow_outer_rev = list(reversed(nozzle_glow_outer))
add_lathe_z("Drive_NozzleBellGlowOuter", "Engines", nozzle_glow_outer_rev, segs=56,
            smooth=True, parts_list=hull_parts, uv_v_from_profile=True, invert_uv_v=True)
add_lathe_z("Drive_NozzleBellGlowInner", "Engines", nozzle_glow_inner, segs=56,
            cap_back=True, smooth=True, parts_list=hull_parts, uv_v_from_profile=True, invert_uv_v=False)

injector_cone_prof = [
    (-75.2, 0.25),
    (-73.5, 1.65),
    (-71.0, 3.10),
    (-68.7, 4.40),
]
add_lathe_z("Drive_InjectorCone", "Engines", injector_cone_prof, segs=36, cap_back=True,
            smooth=True, parts_list=hull_parts, uv_v_from_profile=True, invert_uv_v=True)

bm_drive_metal = bmesh.new()
bm_drive_scorch = bmesh.new()
bm_drive_paint = bmesh.new()

for i_l in range(32):
    phi = i_l * (2.0 * math.pi / 32.0)
    cx, cy = math.cos(phi), math.sin(phi)
    add_strut(bm_drive_scorch, (cx * 7.02, cy * 7.02, -40.2), (cx * 9.36, cy * 9.36, -45.7), 0.28, 0.36, segs=6)

for i_f in range(16):
    phi = i_f * (2.0 * math.pi / 16.0)
    cx, cy = math.cos(phi), math.sin(phi)
    add_strut(bm_drive_metal, (cx * 9.30, cy * 9.30, -62.5), (cx * 7.02, cy * 7.02, -68.9), 0.28, 0.22, segs=8)
    if i_f % 4 == 0:
        add_strut(bm_drive_paint, (cx * 9.52, cy * 9.52, -45.8), (cx * 9.52, cy * 9.52, -62.5), 0.24, 0.24, segs=8)

for i_t in range(24):
    phi = i_t * (2.0 * math.pi / 24.0)
    cx, cy = math.cos(phi), math.sin(phi)
    add_strut(bm_drive_metal, (cx * 7.00, cy * 7.00, -69.2), (cx * 8.52, cy * 8.52, -79.6), 0.07, 0.08, segs=6)

for sx in (-1.0, 1.0):
    add_bevelled_box(bm_drive_metal, (sx * 8.3, 0.0, -65.5), (1.3, 1.7, 3.0), bevel=0.16)
    add_strut(bm_drive_scorch, (sx * 8.3, 0.0, -67.0), (sx * 8.3, 0.0, -68.5), 0.50, 0.72, segs=12)

finalize_bm_part(bm_drive_metal, "Drive_Details_Metal", "Metal", hull_parts)
finalize_bm_part(bm_drive_scorch, "Drive_Details_Scorch", "Scorched", hull_parts)
finalize_bm_part(bm_drive_paint, "Drive_Details_Paint", "Paint", hull_parts)

bpy.ops.object.select_all(action="DESELECT")
for o in hull_parts:
    o.select_set(True)
bpy.context.view_layer.objects.active = hull_parts[0]
bpy.ops.object.join()
hull_obj = bpy.context.active_object
hull_obj.name = "Hull"
hull_obj.data.name = "Hull"
hull_obj.location = (0.0, 0.0, 0.0)

# ---------------------------------------------------------------------------
# 7. Build Moving Nodes (Origins strictly ON pivots, never apply transforms!)
# ---------------------------------------------------------------------------

# === MOVING NODE 1: "Ring" (Rotating Habitat Ring) ===
ring_parts = []
ring_pivot_gltf = Vector((0.0, 0.0, -1.5))

sleeve_prof = [
    (-2.7, 5.05),
    (-2.4, 5.42),
    ( 2.4, 5.42),
    ( 2.7, 5.05),
]
add_lathe_z("Ring_HubSleeve", "Paint", sleeve_prof, segs=40, smooth=True, parts_list=ring_parts)

torus_outer_prof = [
    (-2.75, 18.10),
    (-2.75, 21.45),
    (-2.05, 22.65),
    (-0.85, 22.65),
    (-0.75, 22.48),
    (-0.55, 22.95),
    ( 0.55, 22.95),
    ( 0.75, 22.48),
    ( 0.85, 22.65),
    ( 2.05, 22.65),
    ( 2.75, 21.45),
    ( 2.75, 18.10),
]
add_lathe_z("Ring_TorusOuter", "Paint", torus_outer_prof, segs=64, smooth=True, parts_list=ring_parts)

torus_inner_prof = [
    ( 2.75, 18.10),
    ( 2.05, 17.45),
    (-2.05, 17.45),
    (-2.75, 18.10),
]
add_lathe_z("Ring_TorusInner", "Paint", torus_inner_prof, segs=64, smooth=True, parts_list=ring_parts)

for sz_groove in (-0.80, 0.80):
    groove_prof = [
        (sz_groove - 0.12, 22.68),
        (sz_groove + 0.12, 22.68),
    ]
    add_lathe_z(f"Ring_TrimBand_{sz_groove}", "Metal", groove_prof, segs=64, smooth=True, parts_list=ring_parts)

bm_ring_metal = bmesh.new()
bm_ring_paint = bmesh.new()
bm_ring_lamps = bmesh.new()

for i_sp in range(4):
    ang = i_sp * (math.pi / 2.0)
    cx, cy = math.cos(ang), math.sin(ang)
    tx, ty = -cy, cx
    add_strut(bm_ring_paint, (cx * 5.1, cy * 5.1, 0.0), (cx * 7.2, cy * 7.2, 0.0), 2.15, 1.55, segs=20)
    add_strut(bm_ring_paint, (cx * 7.0, cy * 7.0, 0.0), (cx * 15.8, cy * 15.8, 0.0), 1.52, 1.52, segs=20)
    add_strut(bm_ring_paint, (cx * 15.6, cy * 15.6, 0.0), (cx * 17.8, cy * 17.8, 0.0), 1.55, 2.25, segs=20)
    for r_col in (8.4, 11.5, 14.6):
        add_strut(bm_ring_metal, (cx * (r_col - 0.22), cy * (r_col - 0.22), 0.0),
                  (cx * (r_col + 0.22), cy * (r_col + 0.22), 0.0), 1.66, 1.66, segs=20)
    for side in (-1.0, 1.0):
        ox, oy = tx * side * 1.62, ty * side * 1.62
        add_strut(bm_ring_metal, (cx * 5.4 + ox, cy * 5.4 + oy, 0.0),
                  (cx * 17.5 + ox, cy * 17.5 + oy, 0.0), 0.12, 0.12, segs=8)

for i_seg in range(24):
    phi = i_seg * (2.0 * math.pi / 24.0)
    cx, cy = math.cos(phi), math.sin(phi)
    clamp_pts = [
        (cx * 17.40, cy * 17.40, -2.05),
        (cx * 21.50, cy * 21.50, -2.80),
        (cx * 22.80, cy * 22.80, -1.95),
        (cx * 22.80, cy * 22.80,  1.95),
        (cx * 21.50, cy * 21.50,  2.80),
        (cx * 17.40, cy * 17.40,  2.05),
    ]
    for k in range(len(clamp_pts) - 1):
        add_strut(bm_ring_paint, clamp_pts[k], clamp_pts[k + 1], 0.18, 0.18, segs=6)
    add_strut(bm_ring_paint, clamp_pts[-1], clamp_pts[0], 0.16, 0.16, segs=6)

    mid_phi = phi + (math.pi / 24.0)
    mx, my = math.cos(mid_phi), math.sin(mid_phi)
    add_bevelled_box(bm_ring_lamps, (mx * 17.44, my * 17.44, 0.0), (0.20, 0.20, 0.95), bevel=0.02)
    for sz in (-2.76, 2.76):
        add_bevelled_box(bm_ring_lamps, (mx * 19.8, my * 19.8, sz), (0.28, 0.28, 0.08), bevel=0.01)

finalize_bm_part(bm_ring_metal, "Ring_Details_Metal", "Metal", ring_parts)
finalize_bm_part(bm_ring_paint, "Ring_Details_Paint", "Paint", ring_parts)
finalize_bm_part(bm_ring_lamps, "Ring_Details_Lamps", "Lamps", ring_parts)

bpy.ops.object.select_all(action="DESELECT")
for o in ring_parts:
    o.select_set(True)
bpy.context.view_layer.objects.active = ring_parts[0]
bpy.ops.object.join()
ring_obj = bpy.context.active_object
ring_obj.name = "Ring"
ring_obj.data.name = "Ring"
ring_obj.location = g2b_vec(ring_pivot_gltf)
ring_obj.parent = hull_obj
ring_obj.matrix_parent_inverse = hull_obj.matrix_world.inverted()


# === MOVING NODE 2: "Radar" (Large Parabolic Dish Antenna) ===
# Origin ON its pedestal azimuth/elevation pivot: (0.0, 10.6, -27.2) in glTF. Turns about Y [0, 1, 0].
radar_parts = []
radar_pivot_gltf = Vector((0.0, 10.6, -27.2))

bm_rad_paint = bmesh.new()
bm_rad_metal = bmesh.new()
bm_rad_lamps = bmesh.new()

add_strut(bm_rad_metal, (0.0, -0.3, 0.0), (0.0, 0.6, 0.0), 1.15, 1.15, segs=16)
add_bevelled_box(bm_rad_paint, (0.0, 1.2, -0.8), (2.4, 1.8, 2.6), bevel=0.18)
for sx in (-1.1, 1.1):
    add_strut(bm_rad_metal, (sx, 0.4, 0.0), (sx, 2.8, 0.4), 0.30, 0.26, segs=8)

# Parabolic Dish Bowl (Diameter = 19.0m -> R = 9.5m), elevated so it sits cleanly above the hull
dish_center_local = Vector((0.0, 3.6, 0.8))
tilt_rad = math.radians(38.0)
d_axis = Vector((0.0, math.sin(tilt_rad), math.cos(tilt_rad)))
d_u = Vector((1.0, 0.0, 0.0))
d_v = d_axis.cross(d_u).normalized()

n_dish_rad = 48
n_dish_rings = 10
dish_R = 9.5
dish_depth = 2.85

front_rings = []
back_rings = []
for ir in range(n_dish_rings + 1):
    frac = ir / float(n_dish_rings)
    r = dish_R * frac
    z_para = dish_depth * (frac ** 2)
    z_back = z_para - (0.42 - 0.16 * frac)
    f_ring = []
    b_ring = []
    for k in range(n_dish_rad):
        ang = 2.0 * math.pi * k / n_dish_rad
        radial_vec = d_u * math.cos(ang) + d_v * math.sin(ang)
        p_front = dish_center_local + radial_vec * r + d_axis * z_para
        p_back = dish_center_local + radial_vec * r + d_axis * z_back
        f_ring.append(bm_rad_paint.verts.new(g2b_vec(p_front)))
        b_ring.append(bm_rad_paint.verts.new(g2b_vec(p_back)))
    front_rings.append(f_ring)
    back_rings.append(b_ring)

for ir in range(n_dish_rings):
    for k in range(n_dish_rad):
        kn = (k + 1) % n_dish_rad
        ff = bm_rad_paint.faces.new((front_rings[ir][k], front_rings[ir][kn], front_rings[ir + 1][kn], front_rings[ir + 1][k]))
        ff.smooth = True
        fb = bm_rad_paint.faces.new((back_rings[ir][kn], back_rings[ir][k], back_rings[ir + 1][k], back_rings[ir + 1][kn]))
        fb.smooth = True

for k in range(n_dish_rad):
    kn = (k + 1) % n_dish_rad
    frim = bm_rad_paint.faces.new((front_rings[-1][k], front_rings[-1][kn], back_rings[-1][kn], back_rings[-1][k]))
    frim.smooth = False

for k in range(16):
    ang = 2.0 * math.pi * k / 16.0
    rv = d_u * math.cos(ang) + d_v * math.sin(ang)
    for seg_i in range(5):
        f0 = 0.12 + 0.82 * (seg_i / 5.0)
        f1 = 0.12 + 0.82 * ((seg_i + 1) / 5.0)
        p0 = dish_center_local + rv * (dish_R * f0) + d_axis * (dish_depth * (f0 ** 2) - 0.48)
        p1 = dish_center_local + rv * (dish_R * f1) + d_axis * (dish_depth * (f1 ** 2) - 0.45)
        add_strut(bm_rad_metal, p0, p1, 0.14, 0.11, segs=6)

for frac_c, r_tube in ((0.45, 0.05), (0.75, 0.05), (1.00, 0.14)):
    r_c = dish_R * frac_c
    z_c = dish_depth * (frac_c ** 2) + 0.02
    for k in range(n_dish_rad):
        a0 = 2.0 * math.pi * k / n_dish_rad
        a1 = 2.0 * math.pi * (k + 1) / n_dish_rad
        p0 = dish_center_local + (d_u * math.cos(a0) + d_v * math.sin(a0)) * r_c + d_axis * z_c
        p1 = dish_center_local + (d_u * math.cos(a1) + d_v * math.sin(a1)) * r_c + d_axis * z_c
        add_strut(bm_rad_metal, p0, p1, r_tube, r_tube, segs=6, cap=False)

focal_pt = dish_center_local + d_axis * 7.2
for k in range(3):
    ang = math.radians(90.0 + k * 120.0)
    rv = d_u * math.cos(ang) + d_v * math.sin(ang)
    base_p = dish_center_local + rv * (dish_R * 0.62) + d_axis * (dish_depth * (0.62 ** 2))
    add_strut(bm_rad_metal, base_p, focal_pt, 0.16, 0.10, segs=8)

add_strut(bm_rad_metal, focal_pt - d_axis * 0.9, focal_pt + d_axis * 0.5, 0.68, 0.38, segs=16)
add_strut(bm_rad_paint, focal_pt - d_axis * 0.3, focal_pt + d_axis * 0.3, 0.85, 0.55, segs=16)
add_strut(bm_rad_lamps, focal_pt + d_axis * 0.5, focal_pt + d_axis * 0.68, 0.22, 0.15, segs=12)
add_strut(bm_rad_metal, dish_center_local, dish_center_local + d_axis * 1.5, 1.05, 0.38, segs=16)

finalize_bm_part(bm_rad_paint, "Radar_Dish_Paint", "Paint", radar_parts)
finalize_bm_part(bm_rad_metal, "Radar_Frame_Metal", "Metal", radar_parts)
finalize_bm_part(bm_rad_lamps, "Radar_Beacon_Lamps", "Lamps", radar_parts)

bpy.ops.object.select_all(action="DESELECT")
for o in radar_parts:
    o.select_set(True)
bpy.context.view_layer.objects.active = radar_parts[0]
bpy.ops.object.join()
radar_obj = bpy.context.active_object
radar_obj.name = "Radar"
radar_obj.data.name = "Radar"
radar_obj.location = g2b_vec(radar_pivot_gltf)
radar_obj.parent = hull_obj
radar_obj.matrix_parent_inverse = hull_obj.matrix_world.inverted()


# === MOVING NODES 3 & 4: "Turret_1"/"Gun_1" and "Turret_2"/"Gun_2" ===
def build_turret_and_gun(idx, turret_piv_gltf, sign_y):
    bm_t = bmesh.new()
    add_strut(bm_t, (0.0, 0.0, 0.0), (0.0, sign_y * 0.85, 0.0), 1.25, 0.95, segs=16)
    add_bevelled_box(bm_t, (0.0, sign_y * 0.65, -0.15), (1.8, 0.9, 2.0), bevel=0.16)
    t_obj = finalize_bm_part(bm_t, f"Turret_{idx}", "Metal", [])
    t_obj.location = g2b_vec(turret_piv_gltf)
    t_obj.parent = hull_obj

    gun_offset_gltf = Vector((0.0, sign_y * 0.72, 0.45))
    bm_g = bmesh.new()
    add_bevelled_box(bm_g, (0.0, 0.0, 0.2), (1.1, 0.55, 1.2), bevel=0.08)
    for sx in (-0.34, 0.34):
        add_strut(bm_g, (sx, 0.0, 0.4), (sx, 0.0, 2.3), 0.14, 0.09, segs=10)
        add_strut(bm_g, (sx, 0.0, 2.1), (sx, 0.0, 2.6), 0.12, 0.12, segs=10)
    g_obj = finalize_bm_part(bm_g, f"Gun_{idx}", "Metal", [])
    g_obj.parent = t_obj
    g_obj.location = g2b_vec(gun_offset_gltf)

    muzzle_offset_gltf = Vector((0.0, 0.0, 2.65))
    muz = bpy.data.objects.new(f"Muzzle_{idx}", None)
    bpy.context.collection.objects.link(muz)
    muz.empty_display_size = 0.5
    muz.parent = g_obj
    muz.location = g2b_vec(muzzle_offset_gltf)
    return t_obj, g_obj, muz

turret1_obj, gun1_obj, muzzle1_obj = build_turret_and_gun(1, Vector((0.0, 6.85, -10.2)), 1.0)
turret2_obj, gun2_obj, muzzle2_obj = build_turret_and_gun(2, Vector((0.0, -8.35, -33.5)), -1.0)

# ---------------------------------------------------------------------------
# 8. World-Scale UV Mapping & Clean Triangulation on All Mesh Nodes
# ---------------------------------------------------------------------------
MESH_NODES = [hull_obj, ring_obj, radar_obj, turret1_obj, gun1_obj, turret2_obj, gun2_obj]

def apply_world_box_uvs_and_triangulate(obj, tile_meters=12.0):
    """Applies world-scale box UVs (preserving profile-mapped UVs on Engines) and triangulates n-gons."""
    mesh = obj.data
    if not mesh.uv_layers:
        mesh.uv_layers.new(name="UVMap")
    uv_layer = mesh.uv_layers.active.data
    scale = 1.0 / tile_meters

    # Find material index of Engines so we keep its profile-mapped V gradient
    eng_mat_indices = {i for i, m in enumerate(mesh.materials) if m and m.name == "Engines"}

    for poly in mesh.polygons:
        if poly.material_index in eng_mat_indices:
            continue
        n = poly.normal
        ax, ay, az = abs(n.x), abs(n.y), abs(n.z)
        for loop_idx in poly.loop_indices:
            v = mesh.vertices[mesh.loops[loop_idx].vertex_index].co
            if ax >= ay and ax >= az:
                u, w = v.y * scale, v.z * scale
            elif ay >= ax and ay >= az:
                u, w = v.x * scale, v.z * scale
            else:
                u, w = v.x * scale, v.y * scale
            uv_layer[loop_idx].uv = (u, w)

    # Triangulate polygons with > 4 vertices so glTF tangent export succeeds with 0 warnings
    bm = bmesh.new()
    bm.from_mesh(mesh)
    ngons = [f for f in bm.faces if len(f.verts) > 4]
    if ngons:
        bmesh.ops.triangulate(bm, faces=ngons)
    bm.to_mesh(mesh)
    bm.free()
    mesh.update()

for m_obj in MESH_NODES:
    apply_world_box_uvs_and_triangulate(m_obj, tile_meters=12.0)

# ---------------------------------------------------------------------------
# 9. Empty Socket Nodes (Parented to Hull)
# ---------------------------------------------------------------------------
EMPTY_NODES = [muzzle1_obj, muzzle2_obj]

def add_socket(name, pos_gltf, parent_obj=hull_obj):
    e = bpy.data.objects.new(name, None)
    bpy.context.collection.objects.link(e)
    e.empty_display_size = 1.0
    e.location = g2b_vec(pos_gltf)
    e.parent = parent_obj
    EMPTY_NODES.append(e)
    return e

add_socket("Thruster_1", (0.0, 0.0, -80.0))
add_socket("Thruster_2", (8.3, 0.0, -68.5))
add_socket("Thruster_3", (-8.3, 0.0, -68.5))
add_socket("Hardpoint_1", (6.8, -2.5, -11.7))
add_socket("Hardpoint_2", (-6.8, -2.5, -11.7))
add_socket("Seat_Pilot", (0.0, -0.4, 69.2))
add_socket("Camera_Chase", (0.0, 28.0, -115.0))
add_socket("Light_Nav_L", (35.0, 0.0, -16.5))
add_socket("Light_Nav_R", (-35.0, 0.0, -16.5))
add_socket("Light_Beacon", (0.0, 12.1, 62.5))

ALL_EXPORT_NODES = MESH_NODES + EMPTY_NODES

# ---------------------------------------------------------------------------
# 10. Verify Bounding Box & Export GLB + nodes.json + NOTES.md
# ---------------------------------------------------------------------------
bpy.context.view_layer.update()

pts_blender = [o.matrix_world @ v.co for o in MESH_NODES for v in o.data.vertices]
def b2g_pt(v):
    return Vector((v.x, v.z, -v.y))

pts_gltf = [b2g_pt(p) for p in pts_blender]
b_min = [round(min(p[i] for p in pts_gltf), 4) for i in range(3)]
b_max = [round(max(p[i] for p in pts_gltf), 4) for i in range(3)]
b_size = [round(b_max[i] - b_min[i], 4) for i in range(3)]

total_tris = sum(sum(len(p.vertices) - 2 for p in o.data.polygons) for o in MESH_NODES)
print(f"GEOMETRY SUMMARY: triangles={total_tris}, bbox_min={b_min}, bbox_max={b_max}, size={b_size}", flush=True)

MOVING_AXES = {
    "Ring": [0.0, 0.0, 1.0],
    "Radar": [0.0, 1.0, 0.0],
    "Turret_1": [0.0, 1.0, 0.0],
    "Gun_1": [1.0, 0.0, 0.0],
    "Turret_2": [0.0, 1.0, 0.0],
    "Gun_2": [1.0, 0.0, 0.0],
}

nodes_list = []
for o in ALL_EXPORT_NODES:
    wp_g = b2g_pt(o.matrix_world.translation)
    wp = [round(wp_g.x, 4), round(wp_g.y, 4), round(wp_g.z, 4)]
    nodes_list.append({
        "name": o.name,
        "type": o.type,
        "world_position": wp,
        "position": wp,
        "parent": o.parent.name if o.parent else None,
        "axis": MOVING_AXES.get(o.name, None),
    })

nodes_json = {
    "asset": "science-explorer",
    "contract": "spacecraft/v1",
    "units": "metres",
    "frame": {"forward": "+Z", "up": "+Y", "left": "+X"},
    "axes": MOVING_AXES,
    "bounding_box": {
        "min": b_min,
        "max": b_max,
        "size": b_size,
        "dimensions": b_size,
    },
    "bbox": {
        "min": b_min,
        "max": b_max,
        "size": b_size,
    },
    "triangles": total_tris,
    "materials": list(MATS.keys()),
    "paints": [{"name": n, "linear_rgb": c} for n, c in PAINTS],
    "nodes": nodes_list,
    "source_textures": [
        {"id": ph_id, "license": "CC0-1.0", "url": f"https://polyhaven.com/a/{ph_id}"}
        for ph_id in POLYHAVEN_IDS
    ],
}

with open(os.path.join(ROOT, "nodes.json"), "w") as f:
    json.dump(nodes_json, f, indent=2)

paint_nt = MATS["Paint"].node_tree
paint_bsdf = paint_nt.nodes["Principled BSDF"]
paint_bctex = paint_nt.nodes["BaseColorTex"]
paint_nt.links.new(paint_bctex.outputs["Color"], paint_bsdf.inputs["Base Color"])

bpy.ops.object.select_all(action="DESELECT")
for o in ALL_EXPORT_NODES:
    o.select_set(True)

glb_path = os.path.join(ROOT, "science-explorer.glb")
bpy.ops.export_scene.gltf(
    filepath=glb_path,
    export_format="GLB",
    use_selection=True,
    export_yup=True,
    export_apply=True,
    export_texcoords=True,
    export_normals=True,
    export_tangents=True,
    export_materials="EXPORT",
    export_image_format="AUTO",
)

paint_mix = paint_nt.nodes["PaintFactor"]
paint_nt.links.new(paint_mix.outputs["Result"], paint_bsdf.inputs["Base Color"])

def patch_glb_contract(path):
    data = open(path, "rb").read()
    jl, jt = struct.unpack_from("<II", data, 12)
    doc = json.loads(data[20 : 20 + jl])
    binary = data[20 + jl :]
    for m in doc.get("materials", []):
        pbr = m.setdefault("pbrMetallicRoughness", {})
        if m.get("name") == "Paint":
            pbr["baseColorFactor"] = PAINTS[0][1] + [1.0]
        if m.get("name") == "Glass":
            m["alphaMode"] = "BLEND"
            m["doubleSided"] = True
        m.get("extensions", {}).pop("KHR_materials_clearcoat", None)
    encoded = json.dumps(doc, separators=(",", ":")).encode("utf-8")
    encoded += b" " * ((-len(encoded)) % 4)
    total_len = 20 + len(encoded) + len(binary)
    with open(path, "wb") as f:
        f.write(struct.pack("<III", 0x46546C67, 2, total_len))
        f.write(struct.pack("<II", len(encoded), 0x4E4F534A))
        f.write(encoded)
        f.write(binary)

patch_glb_contract(glb_path)
glb_mb = os.path.getsize(glb_path) / (1024.0 * 1024.0)
print(f"EXPORTED GLB: {glb_path} ({glb_mb:.2f} MB)", flush=True)

with open(os.path.join(ROOT, "NOTES.md"), "w") as f:
    f.write("Deep-space exploration and science vessel: 160m x 70m x 70m hard-sci-fi spacecraft with forward observation dome, open-truss spine, rotating habitat ring, solar/radiator arrays, dorsal dish antenna, sensor booms, and rear magnetic fusion drive.\n")
    f.write(f"{total_tris:,} triangles; GLB size {glb_mb:.2f} MB; Y-up (+Z nose, +Y up, +X port/left), origin at hull centre of mass (0,0,0).\n")
    f.write("8 materials: Paint (neutral detail texture x linear baseColorFactor), Metal, Scorched, Solar, Interior, Glass (alpha BLEND), Engines (emissive), Lamps (emissive).\n")
    f.write("Texture sizes: Paint & Solar basecolor 1024x1024, normal/ORM <=1024/512; Metal/Scorched/Engines 512x512; Interior/Lamps/Glass 256x256.\n")
    f.write("Poly Haven CC0 ids: metal_plate_02, metal_plate, painted_metal_shutter, blue_metal_plate, rusty_metal_03 (diff, nor_gl, arm 1k maps).\n")
    f.write("Reference differences: added full 3D bridge interior behind the observation dome, point-defense/sensor turrets per spacecraft contract, and lateral telemetry booms spanning the 70m width.\n")

# ---------------------------------------------------------------------------
# 11. Space Backdrop, Key/Fill Lighting & Camera Setup
# ---------------------------------------------------------------------------
world = bpy.data.worlds.new("DeepSpaceWorld")
sc.world = world
world.use_nodes = True
wnt = world.node_tree
wnt.nodes.clear()
w_out = wnt.nodes.new("ShaderNodeOutputWorld")
w_bg = wnt.nodes.new("ShaderNodeBackground")
w_texcoord = wnt.nodes.new("ShaderNodeTexCoord")
w_sep = wnt.nodes.new("ShaderNodeSeparateXYZ")
w_ramp = wnt.nodes.new("ShaderNodeValToRGB")
w_ramp.color_ramp.elements[0].position = 0.0
w_ramp.color_ramp.elements[0].color = (0.004, 0.007, 0.016, 1.0)
w_ramp.color_ramp.elements[1].position = 1.0
w_ramp.color_ramp.elements[1].color = (0.018, 0.032, 0.068, 1.0)
w_bg.inputs["Strength"].default_value = 1.0
wnt.links.new(w_texcoord.outputs["Window"], w_sep.inputs["Vector"])
wnt.links.new(w_sep.outputs["Y"], w_ramp.inputs["Fac"])
wnt.links.new(w_ramp.outputs["Color"], w_bg.inputs["Color"])
wnt.links.new(w_bg.outputs["Background"], w_out.inputs["Surface"])

def add_sun_light(name, dir_from_blender, energy, color, angle_deg):
    ld = bpy.data.lights.new(name, "SUN")
    ld.energy = energy
    ld.color = color
    ld.angle = math.radians(angle_deg)
    lo = bpy.data.objects.new(name, ld)
    bpy.context.collection.objects.link(lo)
    lo.location = Vector(dir_from_blender).normalized() * 200.0
    lo.rotation_euler = (-Vector(dir_from_blender)).to_track_quat("-Z", "Y").to_euler()
    return lo

add_sun_light("Key_Sun", (1.1, -1.2, 0.95), 4.0, (1.0, 0.98, 0.95), 6.0)
add_sun_light("Fill_Cool", (-1.2, -0.4, -0.5), 1.4, (0.68, 0.82, 1.0), 35.0)
add_sun_light("Rim_Back", (-0.6, 1.3, 0.8), 1.8, (0.82, 0.90, 1.0), 15.0)

cam_data = bpy.data.cameras.new("RenderCam")
cam_obj = bpy.data.objects.new("RenderCam", cam_data)
bpy.context.collection.objects.link(cam_obj)
sc.camera = cam_obj
cam_data.type = "ORTHO"
cam_data.clip_start = 1.0
cam_data.clip_end = 1000.0

def frame_camera(view_dir_blender, aspect=1920.0 / 1080.0, pad=1.10):
    center = Vector((0.0, 0.0, 0.0))
    d = Vector(view_dir_blender).normalized()
    cam_obj.location = center + d * 260.0
    cam_obj.rotation_euler = (center - cam_obj.location).to_track_quat("-Z", "Y").to_euler()
    bpy.context.view_layer.update()
    inv = cam_obj.matrix_world.inverted()
    proj = [inv @ p for p in pts_blender]
    xmin = min(p.x for p in proj)
    xmax = max(p.x for p in proj)
    ymin = min(p.y for p in proj)
    ymax = max(p.y for p in proj)
    cam_data.ortho_scale = max(xmax - xmin, (ymax - ymin) * aspect) * pad
    shift = cam_obj.rotation_euler.to_matrix() @ Vector(((xmin + xmax) * 0.5, (ymin + ymax) * 0.5, 0.0))
    cam_obj.location += shift
    bpy.context.view_layer.update()

def render_still(filename, w, h, samples=32):
    sc.render.resolution_x = w
    sc.render.resolution_y = h
    sc.render.resolution_percentage = 100
    sc.cycles.samples = samples
    sc.render.image_settings.file_format = "PNG"
    sc.render.filepath = os.path.join(ROOT, filename)
    bpy.ops.render.render(write_still=True)

# 1. Test moving parts rotated 30 degrees
ring_obj.rotation_euler = (0.0, math.radians(30.0), 0.0)
radar_obj.rotation_euler = (0.0, 0.0, math.radians(30.0))
turret1_obj.rotation_euler = (0.0, 0.0, math.radians(30.0))
gun1_obj.rotation_euler = (math.radians(30.0), 0.0, 0.0)
turret2_obj.rotation_euler = (0.0, 0.0, math.radians(30.0))
gun2_obj.rotation_euler = (math.radians(-30.0), 0.0, 0.0)
bpy.context.view_layer.update()

frame_camera((1.12, -1.22, 0.62), 1920.0 / 1080.0, pad=1.08)
sc.render.resolution_x = 960
sc.render.resolution_y = 540
sc.cycles.samples = 12
sc.render.filepath = "/tmp/science_explorer_pivot_test.png"
bpy.ops.render.render(write_still=True)

for mo in (ring_obj, radar_obj, turret1_obj, gun1_obj, turret2_obj, gun2_obj):
    mo.rotation_euler = (0.0, 0.0, 0.0)
bpy.context.view_layer.update()

# 2. Render final.png (1920x1080 Cycles, 3/4 front view)
frame_camera((1.12, -1.22, 0.62), 1920.0 / 1080.0, pad=1.08)
render_still("final.png", 1920, 1080, samples=44)

# 3. Render clay.png (same view as final, plain grey clay override)
clay_mat = bpy.data.materials.new("ClayOverride")
clay_mat.use_nodes = True
c_bsdf = clay_mat.node_tree.nodes["Principled BSDF"]
c_bsdf.inputs["Base Color"].default_value = (0.48, 0.49, 0.50, 1.0)
c_bsdf.inputs["Roughness"].default_value = 0.62
sc.view_layers[0].material_override = clay_mat
render_still("clay.png", 1920, 1080, samples=28)
sc.view_layers[0].material_override = None

# 4. Render side.png (orthographic side view from +X looking at port side, nose to left)
frame_camera((1.0, 0.0, 0.0), 1920.0 / 1080.0, pad=1.08)
render_still("side.png", 1920, 1080, samples=36)

# 5. Render paints.png (3x2 grid of the 3/4 view in all 6 paints)
frame_camera((1.12, -1.22, 0.62), 16.0 / 9.0, pad=1.08)
tile_w, tile_h = 640, 360
canvas = np.zeros((tile_h * 2, tile_w * 3, 4), dtype=np.float32)

for idx, (p_name, p_rgb) in enumerate(PAINTS):
    MATS["Paint"].node_tree.nodes["PaintFactor"].inputs["B"].default_value = (*p_rgb, 1.0)
    tmp_file = f"/tmp/science_explorer_paint_{idx}.png"
    sc.render.resolution_x = tile_w
    sc.render.resolution_y = tile_h
    sc.cycles.samples = 18
    sc.render.filepath = tmp_file
    bpy.ops.render.render(write_still=True)
    im = bpy.data.images.load(tmp_file, check_existing=False)
    arr = np.array(im.pixels[:], dtype=np.float32).reshape(tile_h, tile_w, 4)
    bpy.data.images.remove(im)
    row = 1 - (idx // 3)
    col = idx % 3
    canvas[row * tile_h : (row + 1) * tile_h, col * tile_w : (col + 1) * tile_w, :] = arr
    if os.path.exists(tmp_file):
        os.unlink(tmp_file)

MATS["Paint"].node_tree.nodes["PaintFactor"].inputs["B"].default_value = (*PAINTS[0][1], 1.0)
grid_im = bpy.data.images.new("PaintsGrid", width=tile_w * 3, height=tile_h * 2, alpha=True)
grid_im.pixels.foreach_set(canvas.ravel())
grid_im.filepath_raw = os.path.join(ROOT, "paints.png")
grid_im.file_format = "PNG"
grid_im.save()

print("ALL BUILDS AND RENDERS COMPLETE", flush=True)