VoxelForge/Standalone/Research/WorldAssemblerM3/Linux37/RunMouthCreaseAuditL37.py

#!/usr/bin/env python3
"""Freeze a 20-degree mouth-crease gate and project its sample spacing at 1440p."""
import ast
import contextlib
import importlib.util
import io
import math
import re
import argparse
from pathlib import Path

ROOT = Path(__file__).resolve().parents[4]
BASE_PATH = ROOT / "Standalone/Research/WorldAssemblerM3/Linux37/RunMouthFragmentAuditL37.py"
SPEC = importlib.util.spec_from_file_location("mouth_fragment_l37", BASE_PATH)
BASE = importlib.util.module_from_spec(SPEC)
SPEC.loader.exec_module(BASE)

# Frozen before the L37 mouth edit: 20 degrees over the audit's 0.125 m samples is
# the visible-crease limit. At the 6.016 m standard camera and 80-degree horizontal
# FOV, each sample gap spans multiple pixels at 2560x1440.
BASE.NORMAL_LIMIT_DEG = 20.0
CAMERA = (1.691820, -3.424982, -54.224198)
TARGET = (6.298233, 1.137233, -55.468264)
FOV_DEG = 80.0
WIDTH_1440P = 2560
HEIGHT_1440P = 1440
MOUTH = (5.600, 0.445, -55.156)
OVERVIEW_CAMERA = (MOUTH[0] - 9.75, MOUTH[1] - 11.4, MOUTH[2] + 8.0)


def normalize(vector):
    magnitude = math.sqrt(sum(value * value for value in vector))
    return tuple(value / magnitude for value in vector)


def view_basis(camera, target):
    forward = normalize(tuple(target[i] - camera[i] for i in range(3)))
    right = normalize((forward[1], -forward[0], 0.0))
    up = (right[1] * forward[2] - right[2] * forward[1],
          right[2] * forward[0] - right[0] * forward[2],
          right[0] * forward[1] - right[1] * forward[0])
    return forward, right, up


def visible_crease(path, camera, target):
    cx, cy = BASE.CENTER
    bounds = (cx - BASE.OUTER, cx + BASE.OUTER, cy - BASE.OUTER, cy + BASE.OUTER)
    voxel_um, triangles, _ = BASE.RING.read_scene(
        str(path), 0, BASE.HEIGHT_MIN, BASE.HEIGHT_MAX, bounds)
    cell = max(2.0, 16.0 * voxel_um / 1_000_000.0)
    buckets = BASE.RING.build_buckets(triangles, True, bounds, cell)
    forward, right, up = view_basis(camera, target)
    tan_h = math.tan(math.radians(FOV_DEG / 2.0))
    tan_v = math.tan(math.atan(tan_h * HEIGHT_1440P / WIDTH_1440P))
    focal = WIDTH_1440P / (2.0 * tan_h)
    lines = []
    count = round((BASE.OUTER - BASE.INNER) / BASE.STEP)
    for degrees in (135.0, 180.0, 225.0, 270.0, 315.0):
        angle = math.radians(degrees)
        lines.append([(cx + radius * math.cos(angle), cy + radius * math.sin(angle))
                      for radius in (BASE.INNER + index * BASE.STEP
                                     for index in range(count + 1))])
    tangent_count = round(8.0 / BASE.STEP)
    for degrees in (225.0, 315.0):
        angle = math.radians(degrees)
        normal = (math.cos(angle), math.sin(angle))
        tangent = (-normal[1], normal[0])
        lines.append([(cx + 6.5 * normal[0] + (-4.0 + index * BASE.STEP) * tangent[0],
                       cy + 6.5 * normal[1] + (-4.0 + index * BASE.STEP) * tangent[1])
                      for index in range(tangent_count + 1)])

    maximum, largest_pixels, worst_depth, visible_pairs = 0.0, 0.0, 0.0, 0
    worst_pair_pixels = 0.0
    worst_pair = None

    def project(point):
        offset = tuple(point[i] - camera[i] for i in range(3))
        depth = sum(offset[i] * forward[i] for i in range(3))
        x = sum(offset[i] * right[i] for i in range(3))
        y = sum(offset[i] * up[i] for i in range(3))
        return depth, x, y

    for points in lines:
        samples = [BASE.RING.sample_surface(buckets, x, y, True, cell, BASE.PREFERRED)
                   for x, y in points]
        for index in range(1, len(samples)):
            before, current = samples[index - 1], samples[index]
            if before is None or current is None:
                continue
            p0 = (points[index - 1][0], points[index - 1][1], before[0])
            p1 = (points[index][0], points[index][1], current[0])
            midpoint = tuple((p0[axis] + p1[axis]) * 0.5 for axis in range(3))
            depth, horizontal, vertical = project(midpoint)
            if (depth <= 0.0 or abs(horizontal) > depth * tan_h
                    or abs(vertical) > depth * tan_v):
                continue
            visible_pairs += 1
            angle = math.degrees(math.acos(max(-1.0, min(1.0,
                sum(a * b for a, b in zip(before[1], current[1]))))))
            d0, x0, y0 = project(p0)
            d1, x1, y1 = project(p1)
            pixels = math.hypot(focal * x0 / d0 - focal * x1 / d1,
                                focal * y0 / d0 - focal * y1 / d1)
            largest_pixels = max(largest_pixels, pixels)
            if angle > maximum:
                maximum, worst_pair, worst_depth = angle, (p0, p1), depth
                worst_pair_pixels = pixels
    return maximum, worst_depth, largest_pixels, visible_pairs, worst_pair, worst_pair_pixels


def run(path, label):
    stream = io.StringIO()
    with contextlib.redirect_stdout(stream):
        passed, maximum = BASE.audit(path)
    output = stream.getvalue()
    print(output, end="")
    match = re.search(r"worst_pair=(\([^\n]*\)) status=", output)
    worst = ast.literal_eval(match.group(1)) if match else None
    path = Path(path)
    near = visible_crease(path, CAMERA, TARGET)
    far = visible_crease(path, OVERVIEW_CAMERA, TARGET)
    print(f"[MouthCreasePixels] build={label} sample_step_m={BASE.STEP:.3f} "
          f"horizontal_fov_deg={FOV_DEG:.1f} resolution=2560x1440 "
          f"near_camera_to_mouth_m=6.016 near_visible_pairs={near[3]} "
          f"near_max_visible_deg={near[0]:.6f} near_depth_m={near[1]:.3f} "
          f"near_sample_spacing_px={near[2]:.2f} "
          f"near_worst_pair_spacing_px={near[5]:.2f} "
          f"overview_camera_to_mouth_m={math.dist(OVERVIEW_CAMERA, MOUTH):.3f} "
          f"overview_visible_pairs={far[3]} overview_max_visible_deg={far[0]:.6f} "
          f"overview_depth_m={far[1]:.3f} overview_sample_spacing_px={far[2]:.2f} "
          f"overview_worst_pair_spacing_px={far[5]:.2f} "
          f"overall_max_deg={maximum:.6f}")
    return passed, maximum


def main():
    parser = argparse.ArgumentParser(description=__doc__)
    parser.add_argument("--candidate", type=Path,
                        default=ROOT / "Saved/Linux/l37/owner-fi-ground-mouth.vfmesh")
    args = parser.parse_args()
    paths = {
        "L35": ROOT / "Saved/Linux/l35/owner-fi-ground-mouth.vfmesh",
        "L36": ROOT / "Saved/Linux/l36/owner-fi-ground-mouth.vfmesh",
        "L37": args.candidate if args.candidate.is_absolute()
            else ROOT / args.candidate,
    }
    results = {label: run(path, label) for label, path in paths.items()}
    controls = (not results["L35"][0] and not results["L36"][0]
                and results["L35"][1] > 20.0 and results["L36"][1] > 20.0)
    final = results["L37"][0] and results["L37"][1] <= 20.0
    print(f"[MouthCreaseGate] normal_limit_deg=20 frozen_before_fix=true "
          f"L35_42deg_fail_control={'PASS' if controls else 'FAIL'} "
          f"L36_27deg_fail_control={'PASS' if controls else 'FAIL'} "
          f"L37_candidate={'PASS' if final else 'FAIL'} "
          f"status={'PASS' if controls and final else 'FAIL'}")
    raise SystemExit(0 if controls and final else 1)


if __name__ == "__main__":
    main()