#!/usr/bin/env python3
"""Numerical conformance checks for the Paper VII companion exports.

These tests check the distributed binary64 implementation against the paper's
closed formulas.  They are not a substitute for the symbolic proofs in the PDF.
"""

from __future__ import annotations

import json
import math
import tempfile
from pathlib import Path

import numpy as np

import omega_reference_geometry as omega


ROOT = Path(__file__).resolve().parent
REL_TOL = 2.0e-13
ABS_TOL = 2.0e-13


def close(left: float, right: float, *, message: str) -> None:
    if not math.isclose(left, right, rel_tol=REL_TOL, abs_tol=ABS_TOL):
        raise AssertionError(f"{message}: {left!r} != {right!r}")


def distance(left, right) -> float:
    return math.dist(left, right)


def test_constants() -> None:
    close(omega.PHI * omega.PHI, omega.PHI + 1.0, message="golden identity")
    close(omega.SIN_DELTA**2 + omega.COS_DELTA**2, 1.0, message="trigonometric identity")
    close(2.0 * omega.DELTA, omega.EPSILON, message="half-angle identity")


def test_seed_and_extensions() -> None:
    points = omega.canonical_landmarks()
    alpha, d, b = points["alpha"], points["d"], points["b"]
    close(distance(alpha, d), 1.0, message="left seed radius")
    close(distance(alpha, b), 1.0, message="right seed radius")
    # atan2 is stable at one arcsecond; acos(dot) loses digits near cosine 1.
    seed_angle = math.atan2(omega.norm(omega.cross(d, b)), omega.dot(d, b))
    close(seed_angle, omega.EPSILON, message="seed center angle")
    close(distance(points["d"], points["d1"]), 1.0, message="left first extension")
    close(distance(points["b"], points["b1"]), 1.0, message="right first extension")
    close(distance(points["d1"], points["d2"]), omega.PHI, message="left golden extension")
    close(distance(points["b1"], points["b2"]), omega.PHI, message="right golden extension")
    close(distance(points["d"], points["d2"]), omega.PHI**2,
          message="left cumulative golden displacement")
    close(distance(points["b"], points["b2"]), omega.PHI**2,
          message="right cumulative golden displacement")
    close(distance(points["d2"], points["b2"]), 2.0 * (omega.PHI**2 + omega.SIN_DELTA),
          message="total baseline width")


def test_bearings_and_angles() -> None:
    theta0, theta1, theta2 = omega.positive_bearings()
    close(theta0, omega.EPSILON / 2.0, message="theta0")
    close(theta1, math.pi / 4.0 + omega.EPSILON / 4.0, message="theta1")
    gamma = theta2 - theta1
    close(math.tan(gamma), math.tan(math.pi / 4.0 - omega.EPSILON / 4.0) / math.sqrt(5.0),
          message="golden angular increment")
    bearings = list(omega.signed_bearings().values())
    matrix = omega.angle_matrix()
    for row, left in enumerate(bearings):
        for column, right in enumerate(bearings):
            close(matrix[row][column], abs(left - right), message=f"angle[{row},{column}]")


def test_affine_dimensions_gram_and_distances() -> None:
    x0, x1, x2 = omega.X_LEVELS
    z = np.array([-x2, -x1, -x0, x0, x1, x2], dtype=float)
    configuration = np.vstack((z, np.full(6, omega.COS_DELTA), np.zeros(6)))
    gram = configuration.T @ configuration
    if np.linalg.matrix_rank(configuration, tol=1.0e-12) != 2:
        raise AssertionError("configuration rank is not two")
    if np.linalg.matrix_rank(gram, tol=1.0e-12) != 2:
        raise AssertionError("Gram rank is not two")
    baseline = configuration.T
    baseline_differences = baseline[1:] - baseline[0]
    if np.linalg.matrix_rank(baseline_differences, tol=1.0e-12) != 1:
        raise AssertionError("six baseline marks do not have affine dimension one")
    if np.linalg.matrix_rank(configuration, tol=1.0e-12) != 2:
        raise AssertionError("center-augmented configuration does not have affine dimension two")
    for i in range(6):
        for j in range(6):
            squared = gram[i, i] + gram[j, j] - 2.0 * gram[i, j]
            close(squared, (z[i] - z[j]) ** 2, message=f"Gram distance[{i},{j}]")


def test_outer_scaffold() -> None:
    chord = 2.0 * omega.SIN_DELTA
    width = 2.0 * (omega.PHI**2 + omega.SIN_DELTA)
    height = 2.0 * omega.PHI**2
    close(width - height, chord, message="outer width-minus-height")
    center_to_space = omega.PHI**2 - omega.COS_DELTA
    close(center_to_space, omega.PHI + (1.0 - omega.COS_DELTA), message="center-to-Space")


def test_rotated_frame() -> None:
    rotation = omega.rotation_matrix(math.radians(17.0), math.radians(-23.0), math.radians(31.0))
    frame = omega.transformed_frame(rotation)
    close(omega.norm(frame["u"]), 1.0, message="rotated u norm")
    close(omega.norm(frame["v"]), 1.0, message="rotated v norm")
    close(omega.norm(frame["n"]), 1.0, message="rotated n norm")
    close(omega.dot(frame["u"], frame["v"]), 0.0, message="u-v orthogonality")
    computed_normal = omega.cross(frame["u"], frame["v"])
    for index in range(3):
        close(computed_normal[index], frame["n"][index], message=f"normal component {index}")


def test_parameter_validation() -> None:
    for invalid in (0.0, -1.0, float("inf"), float("nan")):
        try:
            omega.validate_parameters(invalid, 4096)
        except ValueError:
            pass
        else:
            raise AssertionError(f"invalid scale accepted: {invalid!r}")
    try:
        omega.validate_parameters(1.0, 63)
    except ValueError:
        pass
    else:
        raise AssertionError("insufficient circle sampling accepted")


def test_distributed_exports() -> None:
    json_path = ROOT / "generated" / "OMEGA_local_kernel_v1_6_0.json"
    obj_path = ROOT / "generated" / "OMEGA_local_kernel_v1_6_0.obj"
    with tempfile.TemporaryDirectory(prefix="omega-v1-6-0-") as directory:
        temp = Path(directory)
        rebuilt_json = temp / json_path.name
        rebuilt_obj = temp / obj_path.name
        omega.export_json(rebuilt_json)
        omega.export_obj(rebuilt_obj)
        if rebuilt_json.read_bytes() != json_path.read_bytes():
            raise AssertionError("distributed JSON differs from a canonical rebuild")
        if rebuilt_obj.read_bytes() != obj_path.read_bytes():
            raise AssertionError("distributed OBJ differs from a canonical rebuild")
    record = json.loads(json_path.read_text(encoding="utf-8"))
    if record.get("schema") != "OMEGA_MARKED_REFERENCE_GEOMETRY/1.2":
        raise AssertionError(f"incorrect export schema: {record.get('schema')!r}")
    if record.get("exact_formulas", {}).get("golden_accumulation") != (
        "a+a*phi=a*(1+phi)=a*phi^2"
    ):
        raise AssertionError("construction-to-OMEGA arithmetic is absent from the export")
    invariants = record["invariants"]
    expected_invariants = {
        "baseline_affine_dimension": 1,
        "center_augmented_affine_dimension": 2,
        "center_relative_matrix_shape": [3, 6],
        "gram_matrix_shape": [6, 6],
        "center_relative_matrix_rank": 2,
        "gram_matrix_rank": 2,
        "ambient_export_dimension": 3,
    }
    for key, expected in expected_invariants.items():
        if invariants.get(key) != expected:
            raise AssertionError(f"incorrect exported invariant {key}: {invariants.get(key)!r}")
    frame = record["frame_approx"]
    computed = omega.cross(frame["u"], frame["v"])
    for index in range(3):
        close(computed[index], frame["n"][index], message=f"JSON normal component {index}")


TESTS = (
    test_constants,
    test_seed_and_extensions,
    test_bearings_and_angles,
    test_affine_dimensions_gram_and_distances,
    test_outer_scaffold,
    test_rotated_frame,
    test_parameter_validation,
    test_distributed_exports,
)


def main() -> None:
    for test in TESTS:
        test()
        print(f"PASS {test.__name__}")
    print(f"PASS all {len(TESTS)} conformance groups")


if __name__ == "__main__":
    main()