HOA7 Spatial Field Decoder (hoa64 v0.5.0): 7th-order Ambisonics encode/decode, Wigner-D rotation, DOA analysis, vision fuse, diffusion conditioning

README.md

@@ -0,0 +1,126 @@
+---
+language:
+- en
+library_name: other
+pipeline_tag: other
+license: other
+tags:
+- ambisonics
+- spatial-audio
+- hoa
+- spherical-harmonics
+- audio-analysis
+- tool
+- agi-tool
+- wigner-d
+---
+
+# HOA7 Spatial Field Decoder (hoa64)
+
+A **formula-driven spatial field decoder** for **7th-order Ambisonics** — **64 channels** (Ambix **ACN + SN3D**). It encodes/decodes spherical fields, rotates them with **Wigner-D** matrices, analyzes direction-of-arrival / energy, fuses vision detections onto the sphere, and emits **spatial conditioning** for generative pipelines. Pure NumPy; **not an LLM and not a learned model** — deterministic spherical-harmonic geometry, no weights.
+
+| Property | Choice |
+|----------|--------|
+| Version | `0.5.0` |
+| Geometry | Ambix **ACN + SN3D** |
+| Channels | 64 = (7+1)², order 7 |
+| Rotation | **Wigner-D** (~150× vs dense) |
+| Audio | WAV / live mic (ffmpeg·arecord) |
+| Vision | Boxes / YOLO labels / optional torchvision |
+| Agents | HTTP `:8765` + `spatial-report` CLI |
+| Diffusion | Conditioning JSON + optional ComfyUI submit |
+| Language | Python 3, NumPy only (optional torch/torchvision) |
+
+## What it does
+
+- **Encode** point sources / plane waves / scene mixes → HOA-7 coefficients (`encode_points`, `encode_plane_waves`, `encode_scene`).
+- **Decode** coefficients → samples on the sphere, product-grid render, or a single beamform readout (`decode_directions`, `decode_grid`, `beamform`).
+- **Rotate** the field in the listener frame via Wigner-D / zyz (`hoa_rotation_matrix`, `apply_hoa_rotation`, `rotate_yaw_pitch_roll`).
+- **Analyze** DOA (intensity vector + peak search), directional power, field energy, per-STFT-band and per-frame reports (`analysis.py`, `report.py`).
+- **Fuse vision** — object boxes / rays / YOLO detections projected onto the sphere alongside audio (`vision.py`, `detector.py`, `fuse_reports`).
+- **Condition** — spatial reports → plain-text control lines for T2I/T2V prompts, structured JSON for ControlNet-style nodes, or a ComfyUI API payload (`conditioning.py`).
+- **Serve** — HTTP API and iterative agent state (`server.py`, `rnn_stub.py`).
+
+## Coordinates (always)
+
+- **+X** front, **+Y** left, **+Z** up (Ambix listener frame).
+- Azimuth 0° = front, **+90° = left**, −90° = right.
+- Elevation 0° = horizon, +90° = zenith.
+- W (omnidirectional HOA channel) maps to field size / POV: low W → tight / subject-focused / narrow FOV; high W → wide / environmental / immersive FOV.
+
+## Install / run
+
+```bash
+cd spatial-hoa        # this repo root
+export PYTHONPATH="$PWD${PYTHONPATH:+:$PYTHONPATH}"
+python3 -m hoa64 --help
+# or install the CLI on PATH:
+ln -s "$PWD/scripts/spatial-report" ~/.local/bin/spatial-report
+```
+
+**Deps:** NumPy (system `ffmpeg`/`arecord` for live capture; optional torch/torchvision for the real object detector).
+
+## CLI map (`spatial-report`)
+
+| Command | Purpose |
+|---------|---------|
+| `analyze` | Ambix / mono WAV → spatial JSON |
+| `demo-scene` | Synthetic multi-source audio |
+| `vision` | Raw sphere boxes → report |
+| `detect` | Image / YOLO / demo → report |
+| `live` | Mic capture → report |
+| `condition` | Report → diffusion prompt + control vector |
+| `serve` | HTTP API `:8765` |
+
+## HTTP API (`POST /v1/spatial/analyze`)
+
+Modes: `demo_scene` · `ambix_file` · `mono_file` · `vision` · `fuse` · `detect` · `live` · `condition`
+
+```bash
+curl -s -X POST http://127.0.0.1:8765/v1/spatial/analyze \
+  -H 'Content-Type: application/json' \
+  -d '{"mode":"demo_scene","order":3}'
+```
+
+An OpenAI-compatible function schema is provided at `tools/spatial_analyze.openai.json`.
+
+## Agent integration (Qwythos / Pi)
+
+The HOA-7 calculator is **not** part of any LLM's weights. Agents should call it as a tool and treat the returned `one_liner` / `doa_*` / fuse fields as ground-truth geometry:
+
+- Run the CLI: `spatial-report analyze /path/to.wav --ambix -o /tmp/spatial.json && cat /tmp/spatial.json`
+- Or HTTP: `POST http://127.0.0.1:8765/v1/spatial/analyze`
+- See `integrations/qwythos_system_snippet.md` for the exact agent prompt block.
+
+## Quick test pack
+
+```bash
+python3 examples/demo_e2e_testpack.py        # artifacts in /tmp/spatial_hoa_e2e/
+spatial-report analyze /tmp/spatial_hoa_e2e/scene_ambix4.wav --ambix -o /tmp/a.json
+spatial-report detect --demo-image /tmp/frame.png -o /tmp/v.json
+spatial-report condition /tmp/spatial_hoa_e2e/fuse_report.json --prompt 'cinematic interior' -o /tmp/c.json
+```
+
+## Tests
+
+```bash
+python3 tests/test_basis.py
+python3 tests/test_encode_decode.py
+python3 tests/test_rotate_rnn.py
+python3 tests/test_phase1_audio.py
+python3 tests/test_phase2_wigner.py
+python3 tests/test_phase3_vision.py
+python3 tests/test_integration_extras.py
+```
+
+## Limitations
+
+- **Deterministic geometry, not a learned model.** There are no trained weights; accuracy is bounded by the exact spherical-harmonic basis (Farina / Ambix), not by data. The `rnn_stub.py` integrator is an **explicit Euler pose/rotation loop** — learned field dynamics are not implemented yet.
+- **No audio synthesis.** This decodes/analyzes spatial fields and produces conditioning; it does not generate audio content.
+- **Live capture** requires system `ffmpeg`/`arecord`.
+- **Object detection** is optional: the demo path uses synthetic boxes; the real detector downloads torchvision weights on first use.
+- **Coordinate convention is Ambix ACN/SN3D.** Interop with other conventions (FuMa, N3D) requires explicit conversion.
+
+## License
+
+This repository carries no explicit license file. The code is provided as-is by the author; contact `woodfireind` for usage terms. Third-party optional deps (torch/torchvision) keep their own licenses.

examples/demo_e2e_testpack.py

@@ -0,0 +1,85 @@
+#!/usr/bin/env python3
+"""End-to-end pack for manual testing: audio scene, vision detect, fuse, condition.
+
+Writes artifacts under /tmp/spatial_hoa_e2e/ and prints commands for Qwythos/Comfy.
+"""
+
+from __future__ import annotations
+
+import json
+import sys
+from pathlib import Path
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.audio_io import write_wav
+from hoa64.conditioning import build_conditioning, comfy_txt2img_payload, save_conditioning
+from hoa64.detector import detections_to_sphere_boxes, write_demo_image_with_box
+from hoa64.report import report_from_scene
+from hoa64.stream import SourceSpec, encode_scene
+from hoa64.synth import envelope_adsr, tone
+from hoa64.vision import fuse_reports, report_from_boxes
+
+
+def main() -> None:
+    out = Path("/tmp/spatial_hoa_e2e")
+    out.mkdir(parents=True, exist_ok=True)
+
+    sr = 48000
+    dur = 0.5
+    n = int(sr * dur)
+    env = envelope_adsr(n, sr)
+    sources = [
+        SourceSpec(15.0, 0.0, tone(520, dur, sr, amplitude=0.45) * env, "beep"),
+        SourceSpec(-80.0, 5.0, tone(780, dur, sr, amplitude=0.2) * env, "side"),
+    ]
+    audio_rep = report_from_scene(sources, sr, max_order=3)
+    audio_rep.save(out / "audio_report.json")
+    hoa = encode_scene(sources, max_order=3)
+    write_wav(out / "scene_ambix4.wav", hoa[:4], sr)
+
+    img = out / "demo_frame.png"
+    det = write_demo_image_with_box(img)
+    boxes = detections_to_sphere_boxes([det])
+    (out / "boxes.json").write_text(json.dumps(boxes, indent=2) + "\n")
+    vision_rep = report_from_boxes(boxes, max_order=3)
+    vision_rep.save(out / "vision_report.json")
+
+    fused = fuse_reports(
+        {**audio_rep.to_dict(), "one_liner": audio_rep.one_liner()},
+        {**vision_rep.to_dict(), "one_liner": vision_rep.one_liner()},
+    )
+    (out / "fuse_report.json").write_text(json.dumps(fused, indent=2) + "\n")
+
+    cond = build_conditioning(
+        fused,
+        base_prompt="cinematic interior, soft window light, photoreal",
+        style="natural",
+    )
+    save_conditioning(cond, out / "conditioning.json")
+    wf = comfy_txt2img_payload(cond, width=512, height=512, steps=20)
+    (out / "comfy_workflow.json").write_text(json.dumps(wf, indent=2) + "\n")
+
+    print("=== E2E artifacts ===")
+    for p in sorted(out.iterdir()):
+        print(f"  {p}")
+    print("\nAUDIO ", audio_rep.one_liner())
+    print("VISION", vision_rep.one_liner())
+    print("FUSE  ", fused["one_liner"])
+    print("PROMPT", cond["positive_prompt"])
+    print(
+        """
+Test commands:
+  curl -s http://127.0.0.1:8765/health
+  spatial-report analyze /tmp/spatial_hoa_e2e/scene_ambix4.wav --ambix -o /tmp/a.json
+  spatial-report detect --demo-image /tmp/spatial_hoa_e2e/frame2.png -o /tmp/v.json
+  spatial-report condition /tmp/spatial_hoa_e2e/fuse_report.json -o /tmp/c.json --prompt 'moody hall'
+  spatial-report live --duration 1.5 -o /tmp/live.json --write-wav /tmp/live.wav
+  # optional Comfy (if running):
+  spatial-report condition /tmp/spatial_hoa_e2e/fuse_report.json --comfy --write-workflow /tmp/wf.json
+"""
+    )
+
+
+if __name__ == "__main__":
+    main()

examples/demo_phase0.py

@@ -0,0 +1,63 @@
+#!/usr/bin/env python3
+"""Phase 0 demo: encode sources, rotate, iterative walk, print spatial summary."""
+
+from __future__ import annotations
+
+import sys
+from pathlib import Path
+
+import numpy as np
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64 import (
+    N_CHANNELS,
+    channel_names,
+    encode_points,
+    mix,
+    doa_from_intensity,
+    peak_direction,
+    field_energy,
+    beamform,
+)
+from hoa64.rnn_stub import step_rotate, world_from_sources
+
+
+def main() -> None:
+    print(f"hoa64 Phase 0 — channels={N_CHANNELS}")
+    print("names:", ", ".join(channel_names()[:9]), "...")
+
+    # Two sources: front loud, rear quieter
+    field = mix(
+        encode_points([0.0], [0.0], [1.0]),
+        encode_points([180.0], [20.0], [0.4]),
+    )
+    print(f"\nenergy={field_energy(field):.4f}")
+    az, el = doa_from_intensity(field)
+    print(f"intensity DOA: az={az:.1f}° el={el:.1f}°")
+    paz, pel, pv = peak_direction(field)
+    print(f"power peak:    az={paz:.1f}° el={pel:.1f}°  power={pv:.4f}")
+    print(
+        f"beam front={float(beamform(field, 0, 0)):.3f}  "
+        f"rear={float(beamform(field, 180, 20)):.3f}"
+    )
+
+    # Iterative agent motion (order-3 dense for demo speed/quality balance)
+    st = world_from_sources([0.0], [0.0], [1.0])
+    print("\nRNN-stub walk: agent yaws +30° × 3 (order-3 field)")
+    for i in range(3):
+        st = step_rotate(st, d_yaw=30.0, max_order=3, dense=True)
+        snap = st.history[-1]
+        daz, del_ = snap["doa_intensity_az_el"]
+        print(
+            f"  step {i+1}: pose_yaw={st.yaw:.0f}°  "
+            f"head-frame intensity DOA=({daz:.1f},{del_:.1f})  "
+            f"E={snap['energy']:.4f}"
+        )
+
+    print("\nHypothesis check: geometry is formula-driven, state is 64-D, loop integrates pose.")
+    print("Vision tower deferred. Audio encode/analyze/rotate path is live.")
+
+
+if __name__ == "__main__":
+    main()

examples/demo_phase1.py

@@ -0,0 +1,62 @@
+#!/usr/bin/env python3
+"""Phase 1 demo: multi-source scene → Ambix WAV → JSON spatial report for AIs."""
+
+from __future__ import annotations
+
+import json
+import sys
+import tempfile
+from pathlib import Path
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.audio_io import write_wav
+from hoa64.report import report_from_ambix_wav, report_from_scene
+from hoa64.stream import SourceSpec, encode_scene
+from hoa64.synth import envelope_adsr, tone
+
+
+def main() -> None:
+    sr = 48000
+    dur = 0.6
+    n = int(sr * dur)
+    env = envelope_adsr(n, sr, attack=0.02, release=0.08)
+
+    sources = [
+        SourceSpec(0.0, 0.0, tone(440, dur, sr, amplitude=0.45) * env, "front_A4"),
+        SourceSpec(90.0, 10.0, tone(554.37, dur, sr, amplitude=0.3) * env, "left_C#5"),
+        SourceSpec(-120.0, -5.0, tone(659.25, dur, sr, amplitude=0.2) * env, "backright_E5"),
+    ]
+
+    print("=== Synthetic multi-source scene ===")
+    for s in sources:
+        print(f"  {s.label}: az={s.azimuth_deg}° el={s.elevation_deg}°")
+
+    rep = report_from_scene(sources, sr, max_order=3)
+    print("\n" + rep.one_liner())
+    print(f"sources_hint: {json.dumps(rep.sources_hint, indent=2)}")
+    if rep.bands:
+        print("bands (first 3):")
+        for b in rep.bands[:3]:
+            print(
+                f"  {b['band_hz']}: DOA ({b['doa_az_deg']:.1f},{b['doa_el_deg']:.1f}) "
+                f"E={b['energy']:.3g}"
+            )
+
+    out_dir = Path(tempfile.mkdtemp(prefix="hoa64_p1_"))
+    wav_path = out_dir / "scene_ambix4.wav"
+    json_path = out_dir / "spatial_report.json"
+
+    hoa = encode_scene(sources, max_order=3)
+    write_wav(wav_path, hoa[:4], sr)
+    rep.save(json_path)
+
+    # Reload as Ambix and re-analyze
+    rep2 = report_from_ambix_wav(wav_path, max_order=1)
+    print(f"\nFrom written Ambix WAV: {rep2.one_liner()}")
+    print(f"\nArtifacts:\n  {wav_path}\n  {json_path}")
+    print("\nPhase 1: audio path + JSON report for other models — live.")
+
+
+if __name__ == "__main__":
+    main()

examples/demo_phase2.py

@@ -0,0 +1,65 @@
+#!/usr/bin/env python3
+"""Phase 2 demo: Wigner-D vs dense rotation accuracy + timing."""
+
+from __future__ import annotations
+
+import sys
+import time
+from pathlib import Path
+
+import numpy as np
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64 import encode_points, rotate_yaw_pitch_roll, doa_from_intensity
+from hoa64.analysis import angular_error_deg
+from hoa64.rotate import rotate_source_directions
+
+
+def main() -> None:
+    az0, el0 = 25.0, -12.0
+    yaw, pitch, roll = 40.0, 15.0, -10.0
+    a0 = encode_points([az0], [el0], [1.0])
+
+    az1, el1 = rotate_source_directions(az0, el0, yaw=yaw, pitch=pitch, roll=roll)
+    a_gt = encode_points([float(az1)], [float(el1)], [1.0])
+
+    a_w = rotate_yaw_pitch_roll(
+        a0, yaw=yaw, pitch=pitch, roll=roll, method="wigner"
+    )
+    a_d = rotate_yaw_pitch_roll(
+        a0, yaw=yaw, pitch=pitch, roll=roll, method="dense", n_azi=72, n_el=36
+    )
+
+    rel_w = np.linalg.norm(a_w - a_gt) / np.linalg.norm(a_gt)
+    rel_d = np.linalg.norm(a_d - a_gt) / np.linalg.norm(a_gt)
+    print(f"Ground-truth source after rot: az={float(az1):.2f} el={float(el1):.2f}")
+    print(f"Wigner rel error vs re-encode: {rel_w:.3e}")
+    print(f"Dense  rel error vs re-encode: {rel_d:.3e}")
+
+    az_w, el_w = doa_from_intensity(a_w)
+    err = angular_error_deg(float(az1), float(el1), az_w, el_w)
+    print(f"Intensity DOA after Wigner: ({az_w:.2f},{el_w:.2f})  err={err:.3f}°")
+
+    # timing
+    for _ in range(20):
+        rotate_yaw_pitch_roll(a0, yaw=yaw, pitch=pitch, roll=roll, method="wigner")
+    t0 = time.perf_counter()
+    n = 200
+    for _ in range(n):
+        rotate_yaw_pitch_roll(a0, yaw=yaw, pitch=pitch, roll=roll, method="wigner")
+    tw = (time.perf_counter() - t0) / n
+
+    t0 = time.perf_counter()
+    n2 = 5
+    for _ in range(n2):
+        rotate_yaw_pitch_roll(
+            a0, yaw=yaw, pitch=pitch, roll=roll, method="dense", n_azi=48, n_el=24
+        )
+    td = (time.perf_counter() - t0) / n2
+    print(f"Timing: wigner={tw*1e3:.3f} ms  dense={td*1e3:.3f} ms  speedup≈{td/tw:.0f}×")
+    print("Phase 2: Wigner-D rotation is default for pose loops and tools.")
+
+
+if __name__ == "__main__":
+    main()

examples/demo_phase3.py

@@ -0,0 +1,70 @@
+#!/usr/bin/env python3
+"""Phase 3 + Qwythos wiring demo: vision boxes, fuse with audio, HTTP API."""
+
+from __future__ import annotations
+
+import json
+import sys
+import urllib.request
+from pathlib import Path
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.report import report_from_scene
+from hoa64.stream import SourceSpec
+from hoa64.synth import envelope_adsr, tone
+from hoa64.vision import fuse_reports, report_from_boxes
+
+
+def main() -> None:
+    sr = 48000
+    dur = 0.4
+    n = int(sr * dur)
+    env = envelope_adsr(n, sr)
+    audio_rep = report_from_scene(
+        [
+            SourceSpec(20.0, 0.0, tone(500, dur, sr, amplitude=0.5) * env, "beep"),
+        ],
+        sr,
+        max_order=3,
+    )
+    vision_rep = report_from_boxes(
+        [
+            {"az": 25.0, "el": 2.0, "w_deg": 12, "h_deg": 12, "weight": 1.0, "label": "speaker"},
+            {"az": -100.0, "el": 0.0, "w_deg": 8, "h_deg": 8, "weight": 0.3, "label": "clutter"},
+        ],
+        max_order=3,
+    )
+    fused = fuse_reports(audio_rep.to_dict(), {**vision_rep.to_dict(), "one_liner": vision_rep.one_liner()})
+
+    print("AUDIO ", audio_rep.one_liner())
+    print("VISION", vision_rep.one_liner())
+    print("FUSE  ", fused["one_liner"])
+    print(json.dumps({"agreement": fused["agreement"], "sep_deg": fused["angular_separation_deg"]}, indent=2))
+
+    # Optional live API check
+    try:
+        req = urllib.request.Request(
+            "http://127.0.0.1:8765/health",
+            method="GET",
+        )
+        with urllib.request.urlopen(req, timeout=1) as r:
+            health = json.loads(r.read().decode())
+        print("API   ", health)
+        req = urllib.request.Request(
+            "http://127.0.0.1:8765/v1/spatial/analyze",
+            data=json.dumps({"mode": "vision", "order": 3, "boxes": [{"az": 0, "el": 0, "kind": "point"}]}).encode(),
+            headers={"Content-Type": "application/json"},
+            method="POST",
+        )
+        with urllib.request.urlopen(req, timeout=5) as r:
+            api = json.loads(r.read().decode())
+        print("API vision one_liner:", api.get("one_liner"))
+    except Exception as e:
+        print(f"API not running yet ({e}); start: systemctl --user start spatial-hoa")
+
+    print("\nPhase 3 vision + fuse ready. Qwythos: curl :8765 or spatial-report CLI.")
+
+
+if __name__ == "__main__":
+    main()

hoa64/__init__.py

@@ -0,0 +1,92 @@
+"""hoa64 — spatial calculator: 7th-order Ambisonics (Ambix ACN / SN3D).
+
+Fixed spherical-harmonic geometry (Farina / Ambix), not a learned model.
+State dimension is exactly 64 = (7+1)**2.
+
+Phases:
+  0 — basis, encode/decode, rotate, intensity DOA, RNN pose stub
+  1 — WAV I/O, HOA streams, STFT bands, JSON spatial reports for other AIs
+  2 — Fast Wigner-D HOA rotation (default)
+  3 — Vision boxes/rays on sphere + A/V fuse + Qwythos API
+"""
+
+from .basis import (
+    MAX_ORDER,
+    N_CHANNELS,
+    acn_index,
+    acn_nm,
+    channel_names,
+    sh_sn3d,
+    sh_sn3d_batch,
+    unit_vector,
+    az_el_from_unit,
+)
+from .encode import encode_points, encode_plane_waves, mix
+from .decode import decode_directions, decode_grid, beamform
+from .rotate import rotate_yaw_pitch_roll, rotate_matrix_order1
+from .wigner import hoa_rotation_matrix, apply_hoa_rotation, rotation_matrix_zyx
+from .analysis import (
+    field_energy,
+    intensity_vector,
+    doa_from_intensity,
+    directional_power,
+    peak_direction,
+)
+from .report import (
+    SpatialReport,
+    report_from_hoa,
+    report_from_mono_wav,
+    report_from_ambix_wav,
+    report_from_scene,
+)
+from .stream import SourceSpec, encode_mono_plane_wave, encode_scene, analyze_hoa_frames
+from .vision import report_from_boxes, fuse_reports, encode_boxes_to_hoa
+from .detector import detections_to_sphere_boxes, detect_to_sphere
+from .conditioning import build_conditioning, panner_report, spatial_prompt_fragment
+
+__all__ = [
+    "MAX_ORDER",
+    "N_CHANNELS",
+    "acn_index",
+    "acn_nm",
+    "channel_names",
+    "sh_sn3d",
+    "sh_sn3d_batch",
+    "unit_vector",
+    "az_el_from_unit",
+    "encode_points",
+    "encode_plane_waves",
+    "mix",
+    "decode_directions",
+    "decode_grid",
+    "beamform",
+    "rotate_yaw_pitch_roll",
+    "rotate_matrix_order1",
+    "hoa_rotation_matrix",
+    "apply_hoa_rotation",
+    "rotation_matrix_zyx",
+    "field_energy",
+    "intensity_vector",
+    "doa_from_intensity",
+    "directional_power",
+    "peak_direction",
+    "SpatialReport",
+    "report_from_hoa",
+    "report_from_mono_wav",
+    "report_from_ambix_wav",
+    "report_from_scene",
+    "SourceSpec",
+    "encode_mono_plane_wave",
+    "encode_scene",
+    "analyze_hoa_frames",
+    "report_from_boxes",
+    "fuse_reports",
+    "encode_boxes_to_hoa",
+    "detections_to_sphere_boxes",
+    "detect_to_sphere",
+    "build_conditioning",
+    "panner_report",
+    "spatial_prompt_fragment",
+]
+
+__version__ = "0.5.0"

hoa64/__main__.py

@@ -0,0 +1,3 @@
+from .cli import main
+
+raise SystemExit(main())

hoa64/analysis.py

@@ -0,0 +1,118 @@
+"""Spatial field analysis on HOA-7 coefficients (no learning)."""
+
+from __future__ import annotations
+
+import math
+from typing import Tuple
+
+import numpy as np
+
+from .basis import MAX_ORDER, N_CHANNELS, az_el_from_unit, unit_vector
+from .decode import decode_directions, decode_grid
+
+
+def field_energy(hoa: np.ndarray, max_order: int | None = None) -> float | np.ndarray:
+    """Sum of squares of channels (proxy energy; SN3D-weighted optional later).
+
+    hoa: (C,) → scalar; (C,T) → (T,)
+    """
+    a = np.asarray(hoa, dtype=np.float64)
+    if max_order is not None:
+        nch = (max_order + 1) ** 2
+        a = a[..., :nch] if a.ndim == 1 else a[:nch, :]
+    if a.ndim == 1:
+        return float(np.dot(a, a))
+    return np.sum(a * a, axis=0)
+
+
+def intensity_vector(hoa: np.ndarray) -> np.ndarray:
+    """Pseudo-intensity from order-1 SN3D (Ambix).
+
+    For SN3D B-format-like (W,Y,Z,X) = (a0,a1,a2,a3):
+      I ∝ W * (X, Y, Z)   in Cartesian (front, left, up).
+
+    Returns (3,) or (3, T) as [Ix, Iy, Iz].
+    """
+    a = np.asarray(hoa, dtype=np.float64)
+    if a.ndim == 1:
+        W, Y, Z, X = a[0], a[1], a[2], a[3]
+        return np.array([W * X, W * Y, W * Z], dtype=np.float64)
+    W, Y, Z, X = a[0], a[1], a[2], a[3]
+    return np.stack([W * X, W * Y, W * Z], axis=0)
+
+
+def doa_from_intensity(
+    hoa: np.ndarray,
+    *,
+    degrees: bool = True,
+) -> Tuple[float, float] | Tuple[np.ndarray, np.ndarray]:
+    """Direction of arrival from order-1 intensity vector.
+
+    Returns (azimuth, elevation).
+    """
+    I = intensity_vector(hoa)
+    if I.ndim == 1:
+        n = np.linalg.norm(I)
+        if n < 1e-15:
+            return (0.0, 0.0) if degrees else (0.0, 0.0)
+        u = I / n
+        az, el = az_el_from_unit(u, degrees=degrees)
+        return float(az), float(el)
+    # (3, T)
+    n = np.linalg.norm(I, axis=0, keepdims=True)
+    n = np.maximum(n, 1e-15)
+    u = (I / n).T  # (T, 3)
+    az, el = az_el_from_unit(u, degrees=degrees)
+    return az, el
+
+
+def directional_power(
+    hoa: np.ndarray,
+    n_azi: int = 72,
+    n_el: int = 36,
+    *,
+    max_order: int = MAX_ORDER,
+) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
+    """Power map |a · Y(Ω)|² on a sphere grid.
+
+    Returns azi, el, power[azi, el].
+    """
+    azi, el, samp = decode_grid(
+        hoa, n_azi=n_azi, n_el=n_el, degrees=True, max_order=max_order
+    )
+    if samp.ndim == 2:
+        power = samp * samp
+    else:
+        # (A, E, T) → average over time
+        power = np.mean(samp * samp, axis=-1)
+    return azi, el, power
+
+
+def peak_direction(
+    hoa: np.ndarray,
+    n_azi: int = 96,
+    n_el: int = 48,
+    *,
+    max_order: int = MAX_ORDER,
+    degrees: bool = True,
+) -> Tuple[float, float, float]:
+    """Argmax of directional power (azimuth, elevation, peak_value)."""
+    azi, el, power = directional_power(
+        hoa, n_azi=n_azi, n_el=n_el, max_order=max_order
+    )
+    idx = np.unravel_index(int(np.argmax(power)), power.shape)
+    az = float(azi[idx[0]])
+    e = float(el[idx[1]])
+    if not degrees:
+        az, e = math.radians(az), math.radians(e)
+    return az, e, float(power[idx])
+
+
+def angular_error_deg(
+    az0: float, el0: float, az1: float, el1: float
+) -> float:
+    """Great-circle angle between two az/el directions (degrees)."""
+    u = unit_vector(az0, el0, degrees=True)
+    v = unit_vector(az1, el1, degrees=True)
+    c = float(np.clip(np.dot(u, v), -1.0, 1.0))
+    return math.degrees(math.acos(c))

hoa64/audio_io.py

@@ -0,0 +1,120 @@
+"""WAV read/write without scipy/soundfile (stdlib wave + numpy).
+
+Supports mono, multi-channel (e.g. Ambix 4/16/64), and float32/int16.
+"""
+
+from __future__ import annotations
+
+import struct
+import wave
+from pathlib import Path
+from typing import Tuple, Union
+
+import numpy as np
+
+PathLike = Union[str, Path]
+
+
+def read_wav(path: PathLike) -> Tuple[np.ndarray, int]:
+    """Read a WAV file.
+
+    Returns
+    -------
+    audio : np.ndarray
+        Shape (n_channels, n_samples), float64 in roughly [-1, 1].
+    sample_rate : int
+    """
+    path = Path(path)
+    with wave.open(str(path), "rb") as wf:
+        nch = wf.getnchannels()
+        sw = wf.getsampwidth()
+        sr = wf.getframerate()
+        nframes = wf.getnframes()
+        raw = wf.readframes(nframes)
+
+    if sw == 2:
+        mono = np.frombuffer(raw, dtype="<i2").astype(np.float64) / 32768.0
+    elif sw == 4:
+        # Try 32-bit int first; if file is float32 PCM some writers use WAVE_FORMAT_IEEE_FLOAT
+        # which wave may still hand us as bytes — detect by scale.
+        arr_i = np.frombuffer(raw, dtype="<i4")
+        # Heuristic: if max abs > 2^30-ish it's int; else might be float bits misread.
+        # Prefer IEEE float if values look like floats when reinterpreted.
+        arr_f = np.frombuffer(raw, dtype="<f4").astype(np.float64)
+        if np.max(np.abs(arr_f)) <= 8.0 and np.max(np.abs(arr_i)) > 1000:
+            mono = arr_f
+        else:
+            mono = arr_i.astype(np.float64) / 2147483648.0
+    elif sw == 3:
+        # 24-bit packed little-endian
+        a = np.frombuffer(raw, dtype=np.uint8).reshape(-1, 3)
+        vals = (
+            a[:, 0].astype(np.int32)
+            | (a[:, 1].astype(np.int32) << 8)
+            | (a[:, 2].astype(np.int32) << 16)
+        )
+        vals = np.where(vals >= 0x800000, vals - 0x1000000, vals)
+        mono = vals.astype(np.float64) / 8388608.0
+    elif sw == 1:
+        mono = (np.frombuffer(raw, dtype=np.uint8).astype(np.float64) - 128.0) / 128.0
+    else:
+        raise ValueError(f"unsupported sample width {sw}")
+
+    if nch == 1:
+        audio = mono.reshape(1, -1)
+    else:
+        audio = mono.reshape(-1, nch).T.copy()
+    return audio, int(sr)
+
+
+def write_wav(
+    path: PathLike,
+    audio: np.ndarray,
+    sample_rate: int,
+    *,
+    subtype: str = "pcm16",
+) -> None:
+    """Write WAV. audio shape (n_channels, n_samples) or (n_samples,)."""
+    path = Path(path)
+    path.parent.mkdir(parents=True, exist_ok=True)
+    a = np.asarray(audio, dtype=np.float64)
+    if a.ndim == 1:
+        a = a.reshape(1, -1)
+    if a.ndim != 2:
+        raise ValueError("audio must be (C,T) or (T,)")
+    nch, n_samples = a.shape
+    a = np.clip(a, -1.0, 1.0)
+
+    if subtype == "pcm16":
+        pcm = (a.T.reshape(-1) * 32767.0).astype("<i2")
+        sw = 2
+        raw = pcm.tobytes()
+    elif subtype == "float32":
+        pcm = a.T.reshape(-1).astype("<f4")
+        sw = 4
+        raw = pcm.tobytes()
+    else:
+        raise ValueError("subtype must be pcm16 or float32")
+
+    with wave.open(str(path), "wb") as wf:
+        wf.setnchannels(nch)
+        wf.setsampwidth(sw)
+        wf.setframerate(int(sample_rate))
+        # wave module doesn't set IEEE float format tag; pcm16 is portable.
+        if subtype == "float32":
+            # Still write bytes; many tools accept float32 wav with format 3
+            # but stdlib wave always writes PCM. Prefer pcm16 for portability.
+            pass
+        wf.writeframes(raw)
+
+
+def ensure_hoa_channels(audio: np.ndarray, max_order: int = 7) -> np.ndarray:
+    """Pad or truncate multi-channel audio to (max_order+1)**2 Ambix channels."""
+    nch = (max_order + 1) ** 2
+    a = np.asarray(audio, dtype=np.float64)
+    if a.ndim == 1:
+        a = a.reshape(1, -1)
+    out = np.zeros((nch, a.shape[1]), dtype=np.float64)
+    n = min(nch, a.shape[0])
+    out[:n] = a[:n]
+    return out

hoa64/basis.py

@@ -0,0 +1,213 @@
+"""Ambix ACN + SN3D real spherical harmonics up to order 7.
+
+Conventions (Angelo Farina / ISO-style Ambix):
+  - Cartesian: +X front, +Y left, +Z up (ISO 2631-style).
+  - Spherical: azimuth a in XY from +X toward +Y (0=front, +90=left);
+               elevation e from horizontal (+90=zenith, -90=nadir).
+  - Channel order: ACN  n*(n+1)+m  for m = -n .. +n
+  - Normalization: SN3D (Ambix). Order-0 channel W = 1.
+
+Evaluation uses associated Legendre functions (no Condon–Shortley phase)
+with Schmidt/SN3D scaling, matching Farina's explicit Ambix formulas for
+orders 0–5 and 7. Order 6 is generated by the same recurrence (Farina's
+Cartesian transcription for n=6 had inconsistent norms).
+"""
+
+from __future__ import annotations
+
+import math
+from typing import Dict, Tuple
+
+import numpy as np
+
+MAX_ORDER = 7
+N_CHANNELS = (MAX_ORDER + 1) ** 2  # 64
+
+
+def acn_index(n: int, m: int) -> int:
+    """Ambix ACN index for degree n and order m ∈ [-n, n]."""
+    if n < 0 or abs(m) > n:
+        raise ValueError(f"invalid (n,m)=({n},{m})")
+    return n * (n + 1) + m
+
+
+def acn_nm(acn: int) -> Tuple[int, int]:
+    """Inverse of acn_index."""
+    if acn < 0:
+        raise ValueError(acn)
+    n = int(math.floor(math.sqrt(acn)))
+    m = acn - n * (n + 1)
+    return n, m
+
+
+_NAME_N3 = {
+    (0, 0): "W",
+    (1, -1): "Y",
+    (1, 0): "Z",
+    (1, 1): "X",
+    (2, -2): "V",
+    (2, -1): "T",
+    (2, 0): "R",
+    (2, 1): "S",
+    (2, 2): "U",
+    (3, -3): "Q",
+    (3, -2): "O",
+    (3, -1): "M",
+    (3, 0): "K",
+    (3, 1): "L",
+    (3, 2): "N",
+    (3, 3): "P",
+}
+
+
+def channel_names(max_order: int = MAX_ORDER) -> list[str]:
+    names: list[str] = []
+    for n in range(max_order + 1):
+        for m in range(-n, n + 1):
+            names.append(_NAME_N3.get((n, m), f"Y{n}_{m:+d}"))
+    return names
+
+
+def unit_vector(
+    azimuth: float | np.ndarray,
+    elevation: float | np.ndarray,
+    degrees: bool = True,
+) -> np.ndarray:
+    """Direction → unit vector (x, y, z). Broadcasts over arrays. Shape (..., 3)."""
+    a = np.asarray(azimuth, dtype=np.float64)
+    e = np.asarray(elevation, dtype=np.float64)
+    if degrees:
+        a = np.deg2rad(a)
+        e = np.deg2rad(e)
+    ce = np.cos(e)
+    x = np.cos(a) * ce
+    y = np.sin(a) * ce
+    z = np.sin(e)
+    return np.stack([x, y, z], axis=-1)
+
+
+def az_el_from_unit(vec: np.ndarray, degrees: bool = True) -> Tuple[np.ndarray, np.ndarray]:
+    """Unit vector(s) → azimuth, elevation."""
+    v = np.asarray(vec, dtype=np.float64)
+    x, y, z = v[..., 0], v[..., 1], v[..., 2]
+    z = np.clip(z, -1.0, 1.0)
+    el = np.arcsin(z)
+    az = np.arctan2(y, x)
+    if degrees:
+        return np.rad2deg(az), np.rad2deg(el)
+    return az, el
+
+
+def _associated_legendre_no_cs(n_max: int, z: float) -> Dict[Tuple[int, int], float]:
+    """P_n^m(z) without Condon–Shortley phase, 0 ≤ m ≤ n ≤ n_max."""
+    z = float(np.clip(z, -1.0, 1.0))
+    st = math.sqrt(max(0.0, 1.0 - z * z))
+    P: Dict[Tuple[int, int], float] = {(0, 0): 1.0}
+    if n_max >= 1:
+        P[(1, 0)] = z
+        P[(1, 1)] = st
+    for n in range(2, n_max + 1):
+        for m in range(0, n + 1):
+            if m == n:
+                P[(n, n)] = (2 * n - 1) * st * P[(n - 1, n - 1)]
+            elif m == n - 1:
+                P[(n, n - 1)] = (2 * n - 1) * z * P[(n - 1, n - 1)]
+            else:
+                P[(n, m)] = (
+                    (2 * n - 1) * z * P[(n - 1, m)] - (n + m - 1) * P[(n - 2, m)]
+                ) / (n - m)
+    return P
+
+
+def _sn3d_one(x: float, y: float, z: float, max_order: int) -> np.ndarray:
+    """SN3D ACN vector at one unit direction."""
+    r = math.sqrt(x * x + y * y + z * z)
+    if r == 0.0:
+        r = 1.0
+    x, y, z = x / r, y / r, z / r
+    az = math.atan2(y, x)
+    P = _associated_legendre_no_cs(max_order, z)
+    nch = (max_order + 1) ** 2
+    out = np.empty(nch, dtype=np.float64)
+    k = 0
+    for n in range(max_order + 1):
+        for m in range(-n, n + 1):
+            if m == 0:
+                val = P[(n, 0)]
+            else:
+                am = abs(m)
+                norm = math.sqrt(2.0 * math.factorial(n - am) / math.factorial(n + am))
+                if m > 0:
+                    val = norm * P[(n, am)] * math.cos(am * az)
+                else:
+                    val = norm * P[(n, am)] * math.sin(am * az)
+            out[k] = val
+            k += 1
+    return out
+
+
+def _eval_sn3d_cartesian(x: np.ndarray, y: np.ndarray, z: np.ndarray) -> np.ndarray:
+    """Evaluate all 64 SN3D ACN channels at unit direction(s). Output (..., 64)."""
+    x = np.asarray(x, dtype=np.float64)
+    y = np.asarray(y, dtype=np.float64)
+    z = np.asarray(z, dtype=np.float64)
+    shape = np.broadcast(x, y, z).shape
+    x, y, z = np.broadcast_arrays(x, y, z)
+    flat = x.size
+    out = np.empty((flat, N_CHANNELS), dtype=np.float64)
+    xf, yf, zf = x.reshape(-1), y.reshape(-1), z.reshape(-1)
+    for i in range(flat):
+        out[i] = _sn3d_one(float(xf[i]), float(yf[i]), float(zf[i]), MAX_ORDER)
+    return out.reshape(shape + (N_CHANNELS,))
+
+
+def sh_sn3d(
+    azimuth: float | np.ndarray,
+    elevation: float | np.ndarray,
+    degrees: bool = True,
+    max_order: int = MAX_ORDER,
+) -> np.ndarray:
+    """Spherical harmonics Y (Ambix SN3D) at direction(s). Shape (..., n_channels)."""
+    if max_order < 0 or max_order > MAX_ORDER:
+        raise ValueError(f"max_order must be in 0..{MAX_ORDER}")
+    uv = unit_vector(azimuth, elevation, degrees=degrees)
+    y = _eval_sn3d_cartesian(uv[..., 0], uv[..., 1], uv[..., 2])
+    nch = (max_order + 1) ** 2
+    return y[..., :nch]
+
+
+def sh_sn3d_batch(
+    directions_xyz: np.ndarray,
+    max_order: int = MAX_ORDER,
+) -> np.ndarray:
+    """Y at unit vectors. directions_xyz: (..., 3) → (..., n_channels)."""
+    d = np.asarray(directions_xyz, dtype=np.float64)
+    if d.shape[-1] != 3:
+        raise ValueError("last dim must be 3")
+    y = _eval_sn3d_cartesian(d[..., 0], d[..., 1], d[..., 2])
+    nch = (max_order + 1) ** 2
+    return y[..., :nch]
+
+
+def sphere_grid(
+    n_azi: int = 72,
+    n_el: int = 36,
+    degrees: bool = True,
+) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
+    """Product azimuth × elevation grid with integration weights ≈ cos(e) de da."""
+    if degrees:
+        azi = np.linspace(-180.0, 180.0, n_azi, endpoint=False)
+        el = np.linspace(-90.0, 90.0, n_el)
+        el_r = np.deg2rad(el)
+        da = 2.0 * math.pi / n_azi
+        de = math.pi / max(n_el - 1, 1)
+        w_el = np.cos(el_r) * de
+        weights = np.outer(np.full(n_azi, da), np.maximum(w_el, 0.0))
+    else:
+        azi = np.linspace(-math.pi, math.pi, n_azi, endpoint=False)
+        el = np.linspace(-0.5 * math.pi, 0.5 * math.pi, n_el)
+        da = 2.0 * math.pi / n_azi
+        de = math.pi / max(n_el - 1, 1)
+        w_el = np.cos(el) * de
+        weights = np.outer(np.full(n_azi, da), np.maximum(w_el, 0.0))
+    return azi, el, weights

hoa64/cli.py

@@ -0,0 +1,417 @@
+"""CLI: spatial-hoa analyze — produce JSON reports for other models.
+
+Examples
+--------
+  python -m hoa64.cli analyze scene.wav --ambix
+  python -m hoa64.cli analyze mono.wav --az 30 --el 0 -o report.json
+  python -m hoa64.cli demo-scene -o /tmp/scene_report.json
+"""
+
+from __future__ import annotations
+
+import argparse
+import sys
+from pathlib import Path
+
+import numpy as np
+
+from . import __version__
+from .audio_io import write_wav
+from .report import (
+    report_from_ambix_wav,
+    report_from_mono_wav,
+    report_from_scene,
+)
+from .stream import SourceSpec, encode_scene
+from .synth import envelope_adsr, tone
+
+
+def _cmd_analyze(args: argparse.Namespace) -> int:
+    path = Path(args.input)
+    if not path.is_file():
+        print(f"error: file not found: {path}", file=sys.stderr)
+        return 2
+
+    kwargs = dict(
+        max_order=args.order,
+        frame_ms=args.frame_ms,
+        hop_ms=args.hop_ms,
+        include_frames=not args.no_frames,
+        include_bands=not args.no_bands,
+    )
+    if args.ambix:
+        rep = report_from_ambix_wav(path, **kwargs)
+    else:
+        if args.az is None:
+            print(
+                "error: mono plane-wave encode needs --az (or pass --ambix)",
+                file=sys.stderr,
+            )
+            return 2
+        rep = report_from_mono_wav(
+            path, args.az, args.el if args.el is not None else 0.0, **kwargs
+        )
+
+    text = rep.to_json(indent=2 if not args.compact else None)
+    if args.output:
+        Path(args.output).write_text(text + "\n", encoding="utf-8")
+        print(rep.one_liner())
+        print(f"wrote {args.output}")
+    else:
+        print(text)
+    return 0
+
+
+def _cmd_demo_scene(args: argparse.Namespace) -> int:
+    sr = args.sr
+    dur = args.duration
+    n = int(sr * dur)
+    env = envelope_adsr(n, sr)
+    s1 = tone(440.0, dur, sr, amplitude=0.4) * env
+    s2 = tone(660.0, dur, sr, amplitude=0.25) * env
+    sources = [
+        SourceSpec(0.0, 0.0, s1, label="front_A4"),
+        SourceSpec(90.0, 15.0, s2, label="left_E5"),
+    ]
+    rep = report_from_scene(
+        sources,
+        sr,
+        max_order=args.order,
+        frame_ms=args.frame_ms,
+        hop_ms=args.hop_ms,
+    )
+
+    if args.write_wav:
+        hoa = encode_scene(sources, max_order=min(args.order, 3))
+        # Write order-1 B-format-ish (first 4 ch) for portability
+        write_wav(args.write_wav, hoa[:4], sr)
+        print(f"wrote ambix-lite WAV {args.write_wav} (4 ch)")
+
+    text = rep.to_json(indent=2 if not args.compact else None)
+    if args.output:
+        Path(args.output).write_text(text + "\n", encoding="utf-8")
+        print(rep.one_liner())
+        print(f"wrote {args.output}")
+    else:
+        print(text)
+    return 0
+
+
+def _cmd_serve(args: argparse.Namespace) -> int:
+    from .server import serve
+
+    serve(args.host, args.port)
+    return 0
+
+
+def _cmd_detect(args: argparse.Namespace) -> int:
+    import json
+    from .detector import (
+        detect_to_sphere,
+        detections_to_sphere_boxes,
+        load_boxes_json,
+        load_yolo_labels,
+        write_demo_image_with_box,
+        Detection,
+    )
+    from .vision import report_from_boxes
+
+    if args.demo_image:
+        path = Path(args.demo_image)
+        det = write_demo_image_with_box(path)
+        boxes = detections_to_sphere_boxes([det], hfov_deg=args.hfov, vfov_deg=args.vfov)
+        print(f"wrote demo image {path}")
+    elif args.yolo_labels:
+        dets = load_yolo_labels(args.yolo_labels)
+        boxes = detections_to_sphere_boxes(dets, hfov_deg=args.hfov, vfov_deg=args.vfov)
+    elif args.boxes_json:
+        boxes = load_boxes_json(args.boxes_json)
+    elif args.image:
+        img_path = Path(args.image)
+        if not img_path.is_file():
+            print(
+                f"error: image not found: {img_path}\n"
+                "  /path/to/photo.jpg was only an example placeholder.\n"
+                "  Try a real file, or:\n"
+                "    spatial-report detect --demo-image /tmp/demo.png -o /tmp/det.json\n"
+                "    spatial-report detect --image /tmp/spatial_hoa_e2e/demo_frame.png -o /tmp/det.json",
+                file=sys.stderr,
+            )
+            return 2
+        try:
+            boxes = detect_to_sphere(
+                img_path,
+                backend=args.backend,
+                score_thresh=args.score,
+                hfov_deg=args.hfov,
+                vfov_deg=args.vfov,
+            )
+        except Exception as e:
+            print(f"error: detector failed: {e}", file=sys.stderr)
+            print(
+                "  Fallback without neural net:\n"
+                "    spatial-report detect --demo-image /tmp/demo.png -o /tmp/det.json",
+                file=sys.stderr,
+            )
+            return 1
+        if not boxes:
+            print(
+                "warning: no detections above score threshold "
+                "(try --score 0.2, or pass --boxes-json / --demo-image)",
+                file=sys.stderr,
+            )
+    else:
+        print(
+            "error: need --image, --boxes-json, --yolo-labels, or --demo-image\n"
+            "  Example: spatial-report detect --demo-image /tmp/demo.png -o /tmp/det.json",
+            file=sys.stderr,
+        )
+        return 2
+
+    rep = report_from_boxes(boxes, max_order=args.order)
+    text = rep.to_json(indent=2 if not args.compact else None)
+    if args.output:
+        Path(args.output).write_text(text + "\n", encoding="utf-8")
+        print(rep.one_liner())
+        print(f"wrote {args.output}")
+    else:
+        print(text)
+    if args.write_boxes:
+        Path(args.write_boxes).write_text(json.dumps(boxes, indent=2) + "\n")
+        print(f"wrote boxes {args.write_boxes}")
+    return 0
+
+
+def _cmd_live(args: argparse.Namespace) -> int:
+    from .live_audio import live_report, list_pulse_sources
+
+    if args.list_sources:
+        for s in list_pulse_sources():
+            print(s)
+        return 0
+    rep = live_report(
+        duration_sec=args.duration,
+        sample_rate=args.sr,
+        channels=args.channels,
+        source=args.source,
+        az_deg=args.az,
+        el_deg=args.el,
+        max_order=args.order,
+        keep_wav=args.write_wav,
+    )
+    text = rep.to_json(indent=2 if not args.compact else None)
+    if args.output:
+        Path(args.output).write_text(text + "\n", encoding="utf-8")
+        print(rep.one_liner())
+        print(f"wrote {args.output}")
+    else:
+        print(text)
+    return 0
+
+
+def _cmd_condition(args: argparse.Namespace) -> int:
+    import json
+    from .conditioning import (
+        build_conditioning,
+        comfy_txt2img_payload,
+        condition_from_report_file,
+        save_conditioning,
+        submit_comfy_prompt,
+        load_report,
+    )
+
+    rep = load_report(args.report)
+    cond = build_conditioning(
+        rep, base_prompt=args.prompt or "", style=args.style
+    )
+    if args.output:
+        save_conditioning(cond, args.output)
+        print(cond["positive_prompt"])
+        print(f"wrote {args.output}")
+    else:
+        print(json.dumps(cond, indent=2))
+
+    if args.comfy or args.write_workflow:
+        from .conditioning import resolve_comfy_checkpoint, list_comfy_checkpoints
+
+        ckpt = args.checkpoint
+        if not ckpt or args.auto_checkpoint:
+            ckpt = resolve_comfy_checkpoint(ckpt, base_url=args.comfy_url)
+            print(f"comfy checkpoint: {ckpt}")
+        wf = comfy_txt2img_payload(
+            cond,
+            checkpoint=ckpt,
+            width=args.width,
+            height=args.height,
+            steps=args.steps,
+            seed=args.seed,
+            base_url=args.comfy_url,
+            auto_checkpoint=False,
+        )
+        if args.write_workflow:
+            # Write API-ready payload so `curl -d @file` works
+            api_body = {"prompt": wf, "client_id": "spatial-hoa"}
+            Path(args.write_workflow).write_text(
+                json.dumps(api_body, indent=2) + "\n"
+            )
+            print(f"wrote API workflow {args.write_workflow}")
+        if args.comfy:
+            result = submit_comfy_prompt(wf, base_url=args.comfy_url)
+            print(json.dumps(result, indent=2))
+            if result.get("error") and result.get("available_checkpoints"):
+                print(
+                    "available checkpoints:\n  "
+                    + "\n  ".join(result["available_checkpoints"]),
+                    file=sys.stderr,
+                )
+    return 0
+
+
+def _cmd_vision(args: argparse.Namespace) -> int:
+    import json
+    from .vision import report_from_boxes
+
+    if args.boxes_json:
+        boxes = json.loads(Path(args.boxes_json).read_text(encoding="utf-8"))
+    elif args.boxes:
+        boxes = json.loads(args.boxes)
+    else:
+        # demo boxes
+        boxes = [
+            {"az": 0, "el": 0, "w_deg": 10, "h_deg": 10, "weight": 1.0, "label": "front"},
+            {"az": 90, "el": 5, "w_deg": 12, "h_deg": 12, "weight": 0.7, "label": "left"},
+        ]
+    if isinstance(boxes, dict) and "boxes" in boxes:
+        boxes = boxes["boxes"]
+    rep = report_from_boxes(boxes, max_order=args.order)
+    text = rep.to_json(indent=2 if not args.compact else None)
+    if args.output:
+        Path(args.output).write_text(text + "\n", encoding="utf-8")
+        print(rep.one_liner())
+        print(f"wrote {args.output}")
+    else:
+        print(text)
+    return 0
+
+
+def build_parser() -> argparse.ArgumentParser:
+    p = argparse.ArgumentParser(
+        prog="hoa64",
+        description="HOA-7 spatial calculator (audio + vision + JSON for Qwythos)",
+    )
+    p.add_argument("--version", action="version", version=f"hoa64 {__version__}")
+    sub = p.add_subparsers(dest="cmd", required=True)
+
+    a = sub.add_parser("analyze", help="Analyze a WAV → spatial JSON report")
+    a.add_argument("input", help="Path to WAV")
+    a.add_argument(
+        "--ambix",
+        action="store_true",
+        help="Input is multi-channel Ambix ACN (not mono plane-wave)",
+    )
+    a.add_argument("--az", type=float, default=None, help="Plane-wave azimuth (deg)")
+    a.add_argument("--el", type=float, default=0.0, help="Plane-wave elevation (deg)")
+    a.add_argument("--order", type=int, default=7, help="Max HOA order (default 7)")
+    a.add_argument("--frame-ms", type=float, default=40.0)
+    a.add_argument("--hop-ms", type=float, default=20.0)
+    a.add_argument("--no-frames", action="store_true")
+    a.add_argument("--no-bands", action="store_true")
+    a.add_argument("-o", "--output", help="Write JSON to path")
+    a.add_argument("--compact", action="store_true", help="Minified JSON")
+    a.set_defaults(func=_cmd_analyze)
+
+    d = sub.add_parser("demo-scene", help="Synthetic 2-source scene → report")
+    d.add_argument("--sr", type=int, default=48000)
+    d.add_argument("--duration", type=float, default=0.5)
+    d.add_argument("--order", type=int, default=7)
+    d.add_argument("--frame-ms", type=float, default=40.0)
+    d.add_argument("--hop-ms", type=float, default=20.0)
+    d.add_argument("-o", "--output", help="Write JSON to path")
+    d.add_argument("--write-wav", help="Also write 4-ch Ambix-lite WAV")
+    d.add_argument("--compact", action="store_true")
+    d.set_defaults(func=_cmd_demo_scene)
+
+    v = sub.add_parser("vision", help="Vision boxes/rays → spatial JSON report")
+    v.add_argument(
+        "--boxes",
+        help='JSON array of boxes, e.g. \'[{"az":0,"el":0,"weight":1}]\'',
+    )
+    v.add_argument("--boxes-json", help="Path to JSON file with boxes array")
+    v.add_argument("--order", type=int, default=3)
+    v.add_argument("-o", "--output", help="Write JSON to path")
+    v.add_argument("--compact", action="store_true")
+    v.set_defaults(func=_cmd_vision)
+
+    s = sub.add_parser("serve", help="HTTP API for Qwythos/agents (default :8765)")
+    s.add_argument("--host", default="127.0.0.1")
+    s.add_argument("--port", type=int, default=8765)
+    s.set_defaults(func=_cmd_serve)
+
+    det = sub.add_parser("detect", help="Image/YOLO/JSON → vision spatial report")
+    det.add_argument("--image", help="Image path (optional torchvision detector)")
+    det.add_argument("--boxes-json", help="Precomputed boxes or detections JSON")
+    det.add_argument("--yolo-labels", help="YOLO .txt labels for an image")
+    det.add_argument("--demo-image", help="Write synthetic image+box to this path and analyze")
+    det.add_argument("--backend", default="auto", choices=["auto", "torchvision", "none"])
+    det.add_argument("--score", type=float, default=0.5)
+    det.add_argument("--hfov", type=float, default=90.0)
+    det.add_argument("--vfov", type=float, default=60.0)
+    det.add_argument("--order", type=int, default=3)
+    det.add_argument("-o", "--output", help="Spatial report JSON")
+    det.add_argument("--write-boxes", help="Write sphere boxes JSON")
+    det.add_argument("--compact", action="store_true")
+    det.set_defaults(func=_cmd_detect)
+
+    live = sub.add_parser("live", help="Capture mic → spatial report")
+    live.add_argument("--duration", type=float, default=2.0)
+    live.add_argument("--sr", type=int, default=48000)
+    live.add_argument("--channels", type=int, default=1)
+    live.add_argument("--source", help="Pulse source name (pactl list short sources)")
+    live.add_argument("--list-sources", action="store_true")
+    live.add_argument("--az", type=float, default=0.0, help="Plane-wave az if mono")
+    live.add_argument("--el", type=float, default=0.0)
+    live.add_argument("--order", type=int, default=3)
+    live.add_argument("-o", "--output")
+    live.add_argument("--write-wav", help="Keep captured WAV")
+    live.add_argument("--compact", action="store_true")
+    live.set_defaults(func=_cmd_live)
+
+    cond = sub.add_parser("condition", help="Spatial report → diffusion conditioning")
+    cond.add_argument("report", help="Path to spatial/fuse report JSON")
+    cond.add_argument("--prompt", default="cinematic still, photoreal")
+    cond.add_argument("--style", default="natural", choices=["natural", "tags", "technical"])
+    cond.add_argument("-o", "--output", help="Write conditioning JSON")
+    cond.add_argument("--comfy", action="store_true", help="Submit minimal workflow to ComfyUI")
+    cond.add_argument("--comfy-url", default="http://127.0.0.1:8188")
+    cond.add_argument(
+        "--checkpoint",
+        default=None,
+        help="ComfyUI ckpt_name (default: auto-detect from /object_info)",
+    )
+    cond.add_argument(
+        "--auto-checkpoint",
+        action=argparse.BooleanOptionalAction,
+        default=True,
+        help="Resolve checkpoint against ComfyUI's installed list (default: true)",
+    )
+    cond.add_argument("--width", type=int, default=None, help="Latent width (default: 1024 XL / 512 SD)")
+    cond.add_argument("--height", type=int, default=None)
+    cond.add_argument("--steps", type=int, default=20)
+    cond.add_argument("--seed", type=int, default=0)
+    cond.add_argument(
+        "--write-workflow",
+        help="Save ComfyUI API JSON ({prompt: graph}) for curl -d @file",
+    )
+    cond.set_defaults(func=_cmd_condition)
+
+    return p
+
+
+def main(argv: list[str] | None = None) -> int:
+    parser = build_parser()
+    args = parser.parse_args(argv)
+    return int(args.func(args))
+
+
+if __name__ == "__main__":
+    raise SystemExit(main())

hoa64/conditioning.py

@@ -0,0 +1,421 @@
+"""Diffusion / generative conditioning from spatial reports.
+
+Turns HOA calculator output into:
+  * plain-text control lines for T2I / T2V prompts
+  * structured JSON for ControlNet-style / custom nodes
+  * optional ComfyUI API prompt payload (if Comfy is running on :8188)
+"""
+
+from __future__ import annotations
+
+import json
+import urllib.error
+import urllib.request
+from pathlib import Path
+from typing import Any, Mapping, Optional, Union
+
+PathLike = Union[str, Path]
+
+
+def _get(d: Mapping[str, Any], *keys: str, default: Any = None) -> Any:
+    for k in keys:
+        if k in d and d[k] is not None:
+            return d[k]
+    return default
+
+
+def panner_report(
+    az_deg: float = 0.0,
+    el_deg: float = 0.0,
+    w_amplitude: float = 0.5,
+) -> dict:
+    """Synthetic spatial report from a UI spherical panner + W gain.
+
+    W (omnidirectional HOA channel) maps to field size / POV:
+      low W  → tight / subject-focused / narrow FOV
+      high W → wide / environmental / immersive FOV
+    """
+    w = max(0.0, float(w_amplitude))
+    return {
+        "kind": "spatial_panner",
+        "doa_az_deg": float(az_deg),
+        "doa_el_deg": float(el_deg),
+        "w_amplitude": w,
+        "energy": w,
+        "one_liner": (
+            f"panner az={float(az_deg):.0f}° el={float(el_deg):.0f}° W={w:.2f}"
+        ),
+        "meta": {
+            "source": "ui_panner",
+            "field_width_deg": _w_to_field_width_deg(w),
+        },
+    }
+
+
+def _w_to_field_width_deg(w: float) -> float:
+    """Map W amplitude (0..1+) to an angular field width in degrees."""
+    w = max(0.0, min(2.0, float(w)))
+    # 8° pin-point → ~160° ultra-wide at W=1, and beyond at W>1
+    return 8.0 + 152.0 * min(1.0, w) + 40.0 * max(0.0, w - 1.0)
+
+
+def _w_to_field_language(w: float, *, style: str = "natural") -> str:
+    """Natural language for field size / camera POV from W amplitude."""
+    w = max(0.0, float(w))
+    width = _w_to_field_width_deg(w)
+    if style == "tags":
+        if w < 0.25:
+            return f"spatial-fov-tight, spatial-w-{w:.2f}, field-{width:.0f}deg"
+        if w < 0.5:
+            return f"spatial-fov-medium, spatial-w-{w:.2f}, field-{width:.0f}deg"
+        if w < 0.75:
+            return f"spatial-fov-wide, spatial-w-{w:.2f}, field-{width:.0f}deg"
+        return f"spatial-fov-immersive, spatial-w-{w:.2f}, field-{width:.0f}deg"
+    if style == "technical":
+        return f"W={w:.3f} field_width_deg={width:.1f}"
+    if w < 0.2:
+        return (
+            f"tight close-up POV, narrow field of view (~{width:.0f}°), "
+            "subject fills the frame, shallow spatial field"
+        )
+    if w < 0.4:
+        return (
+            f"medium-close framing, moderate field of view (~{width:.0f}°), "
+            "subject-focused with limited environment"
+        )
+    if w < 0.6:
+        return (
+            f"natural mid-shot POV, balanced field of view (~{width:.0f}°), "
+            "subject and surrounding space equally present"
+        )
+    if w < 0.8:
+        return (
+            f"wide environmental framing (~{width:.0f}°), expansive field, "
+            "subject placed in a larger spatial context"
+        )
+    return (
+        f"ultra-wide immersive POV (~{width:.0f}°), large ambient field, "
+        "surrounding space dominates over any single subject"
+    )
+
+
+def spatial_prompt_fragment(
+    report: Mapping[str, Any],
+    *,
+    style: str = "natural",
+) -> str:
+    """Short natural-language spatial control for diffusion prompts.
+
+    style: natural | tags | technical
+
+    When ``w_amplitude`` (or fallback ``energy``) is present, appends field
+    size / POV language so UI panners can drive framing as well as direction.
+    """
+    kind = str(report.get("kind", "spatial_field"))
+    w_raw = _get(report, "w_amplitude", default=None)
+    # Only treat energy as W when the report is from a UI panner (or W is explicit).
+    if w_raw is None and kind == "spatial_panner":
+        w_raw = _get(report, "energy", default=None)
+
+    if kind == "spatial_av_fuse" or "blend_az_deg" in report:
+        a_az = float(_get(report, "audio_doa_az_deg", default=0))
+        a_el = float(_get(report, "audio_doa_el_deg", default=0))
+        v_az = float(_get(report, "vision_doa_az_deg", default=0))
+        v_el = float(_get(report, "vision_doa_el_deg", default=0))
+        sep = float(_get(report, "angular_separation_deg", default=0))
+        agree = bool(_get(report, "agreement", default=False))
+        b_az = float(_get(report, "blend_az_deg", default=a_az))
+        b_el = float(_get(report, "blend_el_deg", default=a_el))
+        if style == "tags":
+            base = (
+                f"spatial-az-{b_az:.0f}, spatial-el-{b_el:.0f}, "
+                f"av-{'aligned' if agree else 'offset'}-{sep:.0f}deg"
+            )
+        elif style == "technical":
+            base = (
+                f"HOA control: blend_az={b_az:.1f} blend_el={b_el:.1f} "
+                f"audio=({a_az:.1f},{a_el:.1f}) vision=({v_az:.1f},{v_el:.1f}) "
+                f"sep={sep:.1f} agree={agree}"
+            )
+        else:
+            side = _az_to_side(b_az)
+            height = _el_to_height(b_el)
+            align = (
+                "sound and subject co-located"
+                if agree
+                else f"sound and subject separated by {sep:.0f} degrees"
+            )
+            base = (
+                f"camera/listener facing forward; primary subject {side}, {height}; "
+                f"{align}; spatial azimuth {b_az:.0f}°, elevation {b_el:.0f}°"
+            )
+        if w_raw is not None:
+            base = f"{base}; {_w_to_field_language(float(w_raw), style=style)}"
+        return base
+
+    az = float(_get(report, "doa_az_deg", "peak_az_deg", default=0))
+    el = float(_get(report, "doa_el_deg", "peak_el_deg", default=0))
+    if style == "tags":
+        base = f"spatial-az-{az:.0f}, spatial-el-{el:.0f}"
+    elif style == "technical":
+        base = f"HOA control: az={az:.1f} el={el:.1f} kind={kind}"
+    else:
+        base = (
+            f"primary direction {_az_to_side(az)}, {_el_to_height(el)}; "
+            f"azimuth {az:.0f}°, elevation {el:.0f}°"
+        )
+    if w_raw is not None:
+        base = f"{base}; {_w_to_field_language(float(w_raw), style=style)}"
+    return base
+
+
+def _az_to_side(az: float) -> str:
+    # Ambix: +az = left
+    if -20 <= az <= 20:
+        return "in front of the camera"
+    if 20 < az <= 70:
+        return "to the front-left"
+    if 70 < az <= 110:
+        return "on the left"
+    if az > 110 or az < -110:
+        return "behind the camera"
+    if -70 <= az < -20:
+        return "to the front-right"
+    return "on the right"
+
+
+def _el_to_height(el: float) -> str:
+    if el > 25:
+        return "above eye level"
+    if el < -25:
+        return "below eye level"
+    return "near eye level"
+
+
+def build_conditioning(
+    report: Mapping[str, Any],
+    *,
+    base_prompt: str = "",
+    negative_prompt: str = "",
+    style: str = "natural",
+) -> dict:
+    """Structured conditioning payload for generative pipelines."""
+    frag = spatial_prompt_fragment(report, style=style)
+    if base_prompt:
+        positive = f"{base_prompt.rstrip(', ').rstrip()}, {frag}"
+    else:
+        positive = frag
+    control = {
+        "schema": "spatial-hoa.conditioning.v1",
+        "spatial_fragment": frag,
+        "positive_prompt": positive,
+        "negative_prompt": negative_prompt,
+        "control_vector": {
+            "az_deg": float(
+                _get(
+                    report,
+                    "blend_az_deg",
+                    "doa_az_deg",
+                    "peak_az_deg",
+                    default=0.0,
+                )
+            ),
+            "el_deg": float(
+                _get(
+                    report,
+                    "blend_el_deg",
+                    "doa_el_deg",
+                    "peak_el_deg",
+                    default=0.0,
+                )
+            ),
+            "energy": float(_get(report, "energy", "audio_energy", default=0.0)),
+            "agreement": _get(report, "agreement", default=None),
+            "angular_separation_deg": _get(
+                report, "angular_separation_deg", default=None
+            ),
+        },
+        "source_report_kind": report.get("kind"),
+        "one_liner": report.get("one_liner") or frag,
+    }
+    return control
+
+
+def list_comfy_checkpoints(base_url: str = "http://127.0.0.1:8188") -> list[str]:
+    """Ask ComfyUI which ckpt_name values are valid."""
+    url = base_url.rstrip("/") + "/object_info/CheckpointLoaderSimple"
+    try:
+        with urllib.request.urlopen(url, timeout=5) as r:
+            data = json.loads(r.read().decode())
+        node = data.get("CheckpointLoaderSimple") or data
+        choices = (
+            node.get("input", {})
+            .get("required", {})
+            .get("ckpt_name", [[]])[0]
+        )
+        # filter non-checkpoint junk (e.g. sam *.pth)
+        return [
+            c
+            for c in choices
+            if isinstance(c, str)
+            and c.endswith((".safetensors", ".ckpt"))
+            and "sam_" not in c.lower()
+        ]
+    except Exception:
+        return []
+
+
+def resolve_comfy_checkpoint(
+    preferred: str | None = None,
+    *,
+    base_url: str = "http://127.0.0.1:8188",
+) -> str:
+    """Pick a checkpoint that exists on this ComfyUI install."""
+    available = list_comfy_checkpoints(base_url)
+    if preferred and preferred in available:
+        return preferred
+    # Prefer SDXL base, then any non-pony XL, then first available
+    for name in available:
+        if name == "sd_xl_base_1.0.safetensors":
+            return name
+    for name in available:
+        if "xl" in name.lower() or "sdxl" in name.lower():
+            return name
+    if available:
+        return available[0]
+    # Offline fallback — may 400 if not installed
+    return preferred or "sd_xl_base_1.0.safetensors"
+
+
+def comfy_txt2img_payload(
+    conditioning: Mapping[str, Any],
+    *,
+    checkpoint: str | None = None,
+    width: int | None = None,
+    height: int | None = None,
+    steps: int = 20,
+    seed: int = 0,
+    cfg: float = 7.0,
+    base_url: str = "http://127.0.0.1:8188",
+    auto_checkpoint: bool = True,
+) -> dict:
+    """Minimal ComfyUI API workflow dict (checkpoint + CLIP + KSampler).
+
+    Load via: POST http://127.0.0.1:8188/prompt  {"prompt": <this>}
+    """
+    if auto_checkpoint or not checkpoint:
+        checkpoint = resolve_comfy_checkpoint(checkpoint, base_url=base_url)
+    # SDXL wants larger latents; SD1.5 512 is fine
+    is_xl = any(t in checkpoint.lower() for t in ("xl", "sdxl", "pony", "zimage"))
+    if width is None:
+        width = 1024 if is_xl else 512
+    if height is None:
+        height = 1024 if is_xl else 512
+
+    positive = conditioning.get("positive_prompt", "")
+    negative = conditioning.get("negative_prompt", "") or (
+        "blurry, low quality, deformed, watermark"
+    )
+    return {
+        "3": {
+            "class_type": "KSampler",
+            "inputs": {
+                "seed": int(seed),
+                "steps": int(steps),
+                "cfg": float(cfg),
+                "sampler_name": "euler",
+                "scheduler": "normal",
+                "denoise": 1.0,
+                "model": ["4", 0],
+                "positive": ["6", 0],
+                "negative": ["7", 0],
+                "latent_image": ["5", 0],
+            },
+        },
+        "4": {
+            "class_type": "CheckpointLoaderSimple",
+            "inputs": {"ckpt_name": checkpoint},
+        },
+        "5": {
+            "class_type": "EmptyLatentImage",
+            "inputs": {"width": int(width), "height": int(height), "batch_size": 1},
+        },
+        "6": {
+            "class_type": "CLIPTextEncode",
+            "inputs": {"text": positive, "clip": ["4", 1]},
+        },
+        "7": {
+            "class_type": "CLIPTextEncode",
+            "inputs": {"text": negative, "clip": ["4", 1]},
+        },
+        "8": {
+            "class_type": "VAEDecode",
+            "inputs": {"samples": ["3", 0], "vae": ["4", 2]},
+        },
+        "9": {
+            "class_type": "SaveImage",
+            "inputs": {"filename_prefix": "spatial_hoa", "images": ["8", 0]},
+        },
+    }
+
+
+def submit_comfy_prompt(
+    workflow: Mapping[str, Any],
+    *,
+    base_url: str = "http://127.0.0.1:8188",
+    client_id: str = "spatial-hoa",
+) -> dict:
+    """POST workflow to ComfyUI. Returns API JSON or error dict with body."""
+    # Accept either raw graph or already-wrapped {"prompt": ...}
+    if "prompt" in workflow and isinstance(workflow.get("prompt"), dict):
+        payload = dict(workflow)
+        payload.setdefault("client_id", client_id)
+    else:
+        payload = {"prompt": dict(workflow), "client_id": client_id}
+
+    url = base_url.rstrip("/") + "/prompt"
+    body = json.dumps(payload).encode("utf-8")
+    req = urllib.request.Request(
+        url, data=body, headers={"Content-Type": "application/json"}, method="POST"
+    )
+    try:
+        with urllib.request.urlopen(req, timeout=30) as r:
+            return json.loads(r.read().decode())
+    except urllib.error.HTTPError as e:
+        err_body = e.read().decode("utf-8", errors="replace")
+        try:
+            detail = json.loads(err_body)
+        except Exception:
+            detail = {"raw": err_body[:2000]}
+        ckpts = list_comfy_checkpoints(base_url)
+        return {
+            "error": f"HTTP {e.code}: {e.reason}",
+            "detail": detail,
+            "available_checkpoints": ckpts,
+            "hint": (
+                "Use --checkpoint <name> from available_checkpoints, "
+                "or omit it to auto-select."
+            ),
+        }
+    except urllib.error.URLError as e:
+        return {"error": str(e), "hint": "Is ComfyUI running on :8188?"}
+    except Exception as e:
+        return {"error": str(e)}
+
+
+def save_conditioning(cond: Mapping[str, Any], path: PathLike) -> None:
+    Path(path).write_text(json.dumps(cond, indent=2) + "\n", encoding="utf-8")
+
+
+def load_report(path: PathLike) -> dict:
+    return json.loads(Path(path).read_text(encoding="utf-8"))
+
+
+def condition_from_report_file(
+    report_path: PathLike,
+    *,
+    base_prompt: str = "cinematic still, photoreal",
+    style: str = "natural",
+) -> dict:
+    rep = load_report(report_path)
+    return build_conditioning(rep, base_prompt=base_prompt, style=style)

hoa64/decode.py

@@ -0,0 +1,70 @@
+"""Decode HOA-7 coefficients to samples on the sphere."""
+
+from __future__ import annotations
+
+import numpy as np
+
+from .basis import MAX_ORDER, sh_sn3d, sh_sn3d_batch, sphere_grid
+
+
+def decode_directions(
+    hoa: np.ndarray,
+    azimuths: np.ndarray | float,
+    elevations: np.ndarray | float,
+    *,
+    degrees: bool = True,
+    max_order: int = MAX_ORDER,
+) -> np.ndarray:
+    """Sample field a · Y(Ω) at given directions.
+
+    hoa: (C,) or (C, T)
+    returns: (...) or (..., T)
+    """
+    a = np.asarray(hoa, dtype=np.float64)
+    nch = (max_order + 1) ** 2
+    # Allow shorter coefficient vectors (truncated order) or longer (ignore tail).
+    if a.ndim == 1:
+        aa = np.zeros(nch, dtype=np.float64)
+        n = min(nch, a.shape[0])
+        aa[:n] = a[:n]
+        Y = sh_sn3d(azimuths, elevations, degrees=degrees, max_order=max_order)
+        return np.einsum("...c,c->...", Y[..., :nch], aa)
+    if a.ndim == 2:
+        aa = np.zeros((nch, a.shape[1]), dtype=np.float64)
+        n = min(nch, a.shape[0])
+        aa[:n] = a[:n, :]
+        Y = sh_sn3d(azimuths, elevations, degrees=degrees, max_order=max_order)
+        return np.einsum("...c,ct->...t", Y[..., :nch], aa)
+    raise ValueError("hoa must be (C,) or (C,T)")
+
+
+def decode_grid(
+    hoa: np.ndarray,
+    n_azi: int = 72,
+    n_el: int = 36,
+    *,
+    degrees: bool = True,
+    max_order: int = MAX_ORDER,
+) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
+    """Decode on a product sphere grid.
+
+    Returns azi, el, samples with samples shape (n_azi, n_el) or (n_azi, n_el, T).
+    """
+    azi, el, _ = sphere_grid(n_azi, n_el, degrees=degrees)
+    AA, EE = np.meshgrid(azi, el, indexing="ij")
+    samp = decode_directions(hoa, AA, EE, degrees=degrees, max_order=max_order)
+    return azi, el, samp
+
+
+def beamform(
+    hoa: np.ndarray,
+    azimuth: float,
+    elevation: float,
+    *,
+    degrees: bool = True,
+    max_order: int = MAX_ORDER,
+) -> float | np.ndarray:
+    """Look / listen toward one direction: scalar (or time series) readout."""
+    return decode_directions(
+        hoa, azimuth, elevation, degrees=degrees, max_order=max_order
+    )

hoa64/detector.py

@@ -0,0 +1,304 @@
+"""Detector → sphere boxes adapter for the vision tower.
+
+Supports:
+  * Precomputed boxes (JSON / dict list) — always available
+  * YOLO-format label files (class x_c y_c w h normalized)
+  * Optional torchvision detection (Faster R-CNN MobileNet) if weights load
+
+No ultralytics/cv2 required. Image I/O via PIL + numpy.
+"""
+
+from __future__ import annotations
+
+import json
+from dataclasses import dataclass, asdict
+from pathlib import Path
+from typing import Any, Mapping, Optional, Sequence, Union
+
+import numpy as np
+
+PathLike = Union[str, Path]
+
+# COCO class names subset for friendlier labels (optional)
+_COCO_NAMES = None
+
+
+def _coco_names() -> list[str]:
+    global _COCO_NAMES
+    if _COCO_NAMES is not None:
+        return _COCO_NAMES
+    # Minimal common set; unknown ids → class_{id}
+    _COCO_NAMES = [
+        "__background__", "person", "bicycle", "car", "motorcycle", "airplane",
+        "bus", "train", "truck", "boat", "traffic light", "fire hydrant",
+        "stop sign", "parking meter", "bench", "bird", "cat", "dog", "horse",
+        "sheep", "cow", "elephant", "bear", "zebra", "giraffe", "backpack",
+        "umbrella", "handbag", "tie", "suitcase", "frisbee", "skis",
+        "snowboard", "sports ball", "kite", "baseball bat", "baseball glove",
+        "skateboard", "surfboard", "tennis racket", "bottle", "wine glass",
+        "cup", "fork", "knife", "spoon", "bowl", "banana", "apple", "sandwich",
+        "orange", "broccoli", "carrot", "hot dog", "pizza", "donut", "cake",
+        "chair", "couch", "potted plant", "bed", "dining table", "toilet",
+        "tv", "laptop", "mouse", "remote", "keyboard", "cell phone", "microwave",
+        "oven", "toaster", "sink", "refrigerator", "book", "clock", "vase",
+        "scissors", "teddy bear", "hair drier", "toothbrush",
+    ]
+    return _COCO_NAMES
+
+
+@dataclass
+class Detection:
+    """Axis-aligned box in normalized image coords (x,y center or corners)."""
+
+    # normalized [0,1] image: x right, y down; origin top-left
+    x1: float
+    y1: float
+    x2: float
+    y2: float
+    score: float = 1.0
+    label: str = ""
+    class_id: int = -1
+
+    @property
+    def cx(self) -> float:
+        return 0.5 * (self.x1 + self.x2)
+
+    @property
+    def cy(self) -> float:
+        return 0.5 * (self.y1 + self.y2)
+
+    @property
+    def w(self) -> float:
+        return max(0.0, self.x2 - self.x1)
+
+    @property
+    def h(self) -> float:
+        return max(0.0, self.y2 - self.y1)
+
+    def to_sphere_box(
+        self,
+        *,
+        hfov_deg: float = 90.0,
+        vfov_deg: float = 60.0,
+        sigma_scale: float = 1.0,
+    ) -> dict:
+        """Map image box → HOA vision box (az/el degrees, Ambix convention)."""
+        # center: image x=0 left → +az, x=1 right → -az (matches vision._as_box)
+        az = (0.5 - self.cx) * hfov_deg
+        el = (0.5 - self.cy) * vfov_deg
+        w_deg = max(2.0, self.w * hfov_deg * sigma_scale)
+        h_deg = max(2.0, self.h * vfov_deg * sigma_scale)
+        return {
+            "az": float(az),
+            "el": float(el),
+            "w_deg": float(w_deg),
+            "h_deg": float(h_deg),
+            "weight": float(self.score),
+            "label": self.label or (f"class_{self.class_id}" if self.class_id >= 0 else "det"),
+            "kind": "box",
+        }
+
+
+def detections_to_sphere_boxes(
+    dets: Sequence[Detection],
+    *,
+    hfov_deg: float = 90.0,
+    vfov_deg: float = 60.0,
+    min_score: float = 0.25,
+) -> list[dict]:
+    out = []
+    for d in dets:
+        if d.score < min_score:
+            continue
+        out.append(d.to_sphere_box(hfov_deg=hfov_deg, vfov_deg=vfov_deg))
+    return out
+
+
+def load_boxes_json(path: PathLike) -> list[dict]:
+    """Load sphere boxes or detections JSON.
+
+    Accepts:
+      [{"az":..., "el":...}, ...]
+      {"boxes": [...]}
+      {"detections": [{"x1","y1","x2","y2",...}, ...]}  # normalized
+    """
+    data = json.loads(Path(path).read_text(encoding="utf-8"))
+    if isinstance(data, list):
+        if data and ("x1" in data[0] or "bbox" in data[0]):
+            dets = []
+            for item in data:
+                if "bbox" in item:
+                    x1, y1, x2, y2 = item["bbox"]
+                else:
+                    x1, y1, x2, y2 = item["x1"], item["y1"], item["x2"], item["y2"]
+                dets.append(
+                    Detection(
+                        float(x1), float(y1), float(x2), float(y2),
+                        score=float(item.get("score", item.get("confidence", 1.0))),
+                        label=str(item.get("label", item.get("class", ""))),
+                        class_id=int(item.get("class_id", -1)),
+                    )
+                )
+            return detections_to_sphere_boxes(dets)
+        return list(data)
+    if isinstance(data, dict):
+        if "boxes" in data:
+            return list(data["boxes"])
+        if "detections" in data:
+            return load_boxes_json_from_obj(data["detections"])
+    raise ValueError(f"unrecognized boxes JSON shape in {path}")
+
+
+def load_boxes_json_from_obj(obj: Any) -> list[dict]:
+    path_like = Path("/tmp/_unused")
+    # reuse logic
+    if isinstance(obj, list) and obj and "x1" in obj[0]:
+        dets = [
+            Detection(
+                float(i["x1"]), float(i["y1"]), float(i["x2"]), float(i["y2"]),
+                score=float(i.get("score", 1.0)),
+                label=str(i.get("label", "")),
+            )
+            for i in obj
+        ]
+        return detections_to_sphere_boxes(dets)
+    if isinstance(obj, list):
+        return list(obj)
+    raise ValueError("bad detections object")
+
+
+def load_yolo_labels(
+    path: PathLike,
+    *,
+    class_names: Optional[Sequence[str]] = None,
+) -> list[Detection]:
+    """YOLO txt: class x_center y_center width height (all normalized)."""
+    dets: list[Detection] = []
+    text = Path(path).read_text(encoding="utf-8").strip()
+    if not text:
+        return dets
+    for line in text.splitlines():
+        parts = line.split()
+        if len(parts) < 5:
+            continue
+        cid = int(float(parts[0]))
+        cx, cy, w, h = map(float, parts[1:5])
+        score = float(parts[5]) if len(parts) > 5 else 1.0
+        x1, y1 = cx - w / 2, cy - h / 2
+        x2, y2 = cx + w / 2, cy + h / 2
+        label = ""
+        if class_names and 0 <= cid < len(class_names):
+            label = class_names[cid]
+        else:
+            label = f"class_{cid}"
+        dets.append(Detection(x1, y1, x2, y2, score=score, label=label, class_id=cid))
+    return dets
+
+
+def detect_torchvision(
+    image_path: PathLike,
+    *,
+    score_thresh: float = 0.5,
+    device: Optional[str] = None,
+    max_dets: int = 32,
+) -> list[Detection]:
+    """Run torchvision Faster R-CNN MobileNet on an image.
+
+    First call may download weights (~50MB). Uses XPU if available.
+    """
+    from PIL import Image
+    import torch
+    import torchvision
+    from torchvision.transforms import functional as F
+
+    img = Image.open(image_path).convert("RGB")
+    w, h = img.size
+    tensor = F.to_tensor(img)
+
+    if device is None:
+        if hasattr(torch, "xpu") and torch.xpu.is_available():
+            device = "xpu"
+        elif torch.cuda.is_available():
+            device = "cuda"
+        else:
+            device = "cpu"
+
+    weights = torchvision.models.detection.FasterRCNN_MobileNet_V3_Large_FPN_Weights.DEFAULT
+    model = torchvision.models.detection.fasterrcnn_mobilenet_v3_large_fpn(weights=weights)
+    model.eval()
+    model.to(device)
+    with torch.no_grad():
+        out = model([tensor.to(device)])[0]
+
+    names = _coco_names()
+    boxes = out["boxes"].detach().cpu().numpy()
+    scores = out["scores"].detach().cpu().numpy()
+    labels = out["labels"].detach().cpu().numpy()
+    dets: list[Detection] = []
+    for box, sc, lab in zip(boxes, scores, labels):
+        if sc < score_thresh:
+            continue
+        x1, y1, x2, y2 = box
+        cid = int(lab)
+        label = names[cid] if cid < len(names) else f"class_{cid}"
+        dets.append(
+            Detection(
+                x1 / w, y1 / h, x2 / w, y2 / h,
+                score=float(sc),
+                label=label,
+                class_id=cid,
+            )
+        )
+        if len(dets) >= max_dets:
+            break
+    return dets
+
+
+def detect_to_sphere(
+    image_path: PathLike,
+    *,
+    backend: str = "auto",
+    score_thresh: float = 0.5,
+    hfov_deg: float = 90.0,
+    vfov_deg: float = 60.0,
+) -> list[dict]:
+    """Image → sphere boxes.
+
+    backend: auto | torchvision | none
+    auto tries torchvision, falls back to empty with a note if unavailable.
+    """
+    if backend in ("auto", "torchvision"):
+        try:
+            dets = detect_torchvision(image_path, score_thresh=score_thresh)
+            return detections_to_sphere_boxes(
+                dets, hfov_deg=hfov_deg, vfov_deg=vfov_deg, min_score=score_thresh
+            )
+        except Exception as e:
+            if backend == "torchvision":
+                raise
+            return []
+    return []
+
+
+def write_demo_image_with_box(
+    path: PathLike,
+    *,
+    size: tuple[int, int] = (640, 480),
+    box_xyxy_norm: tuple[float, float, float, float] = (0.35, 0.35, 0.55, 0.65),
+) -> Detection:
+    """Create a simple synthetic image + known detection (no model)."""
+    from PIL import Image, ImageDraw
+
+    w, h = size
+    img = Image.new("RGB", (w, h), (30, 30, 40))
+    draw = ImageDraw.Draw(img)
+    x1, y1, x2, y2 = box_xyxy_norm
+    px = [x1 * w, y1 * h, x2 * w, y2 * h]
+    draw.rectangle(px, outline=(0, 255, 80), width=4)
+    draw.ellipse(
+        [px[0] + 10, px[1] + 10, px[2] - 10, px[3] - 10],
+        fill=(200, 80, 80),
+    )
+    Path(path).parent.mkdir(parents=True, exist_ok=True)
+    img.save(path)
+    return Detection(x1, y1, x2, y2, score=0.99, label="demo_object", class_id=0)

hoa64/encode.py

@@ -0,0 +1,91 @@
+"""Encode discrete sources into HOA-7 (Ambix SN3D) coefficient vectors."""
+
+from __future__ import annotations
+
+from typing import Sequence
+
+import numpy as np
+
+from .basis import MAX_ORDER, N_CHANNELS, sh_sn3d, sh_sn3d_batch, unit_vector
+
+
+def encode_points(
+    azimuths: Sequence[float] | np.ndarray,
+    elevations: Sequence[float] | np.ndarray,
+    gains: Sequence[float] | np.ndarray | None = None,
+    *,
+    degrees: bool = True,
+    max_order: int = MAX_ORDER,
+) -> np.ndarray:
+    """Encode point / plane-wave sources on the sphere into HOA coeffs.
+
+    Each source at (az, el) with gain g contributes g * Y(az, el).
+    Returns shape (n_channels,) float64. Time-domain: call per sample
+    with complex gains or use encode_plane_waves for signals.
+    """
+    az = np.atleast_1d(np.asarray(azimuths, dtype=np.float64))
+    el = np.atleast_1d(np.asarray(elevations, dtype=np.float64))
+    if az.shape != el.shape:
+        raise ValueError("azimuths and elevations must match")
+    if gains is None:
+        g = np.ones(az.shape, dtype=np.float64)
+    else:
+        g = np.asarray(gains, dtype=np.float64)
+        g = np.broadcast_to(g, az.shape)
+    Y = sh_sn3d(az, el, degrees=degrees, max_order=max_order)  # (N, C)
+    # a = sum_i g_i Y_i
+    return np.einsum("n,nc->c", g.reshape(-1), Y.reshape(-1, Y.shape[-1]))
+
+
+def encode_plane_waves(
+    azimuths: Sequence[float] | np.ndarray,
+    elevations: Sequence[float] | np.ndarray,
+    signals: np.ndarray,
+    *,
+    degrees: bool = True,
+    max_order: int = MAX_ORDER,
+) -> np.ndarray:
+    """Encode multi-channel time signals as plane waves.
+
+    signals: (n_sources, n_samples) or (n_sources,)
+    Returns: (n_channels, n_samples) or (n_channels,)
+    """
+    az = np.atleast_1d(np.asarray(azimuths, dtype=np.float64))
+    el = np.atleast_1d(np.asarray(elevations, dtype=np.float64))
+    sig = np.asarray(signals, dtype=np.float64)
+    if sig.ndim == 1:
+        if sig.shape[0] == az.shape[0]:
+            # (n_sources,) static frame
+            return encode_points(az, el, sig, degrees=degrees, max_order=max_order)
+        raise ValueError("1-D signals length must equal n_sources")
+    if sig.ndim != 2 or sig.shape[0] != az.shape[0]:
+        raise ValueError("signals must be (n_sources, n_samples)")
+    Y = sh_sn3d(az, el, degrees=degrees, max_order=max_order)  # (S, C)
+    # a[c,t] = sum_s Y[s,c] * sig[s,t]
+    return np.einsum("sc,st->ct", Y, sig)
+
+
+def mix(*fields: np.ndarray) -> np.ndarray:
+    """Superpose HOA fields (same channel count)."""
+    if not fields:
+        return np.zeros(N_CHANNELS, dtype=np.float64)
+    out = np.zeros_like(np.asarray(fields[0], dtype=np.float64))
+    for f in fields:
+        out = out + np.asarray(f, dtype=np.float64)
+    return out
+
+
+def encode_xyz(
+    directions_xyz: np.ndarray,
+    gains: np.ndarray | None = None,
+    *,
+    max_order: int = MAX_ORDER,
+) -> np.ndarray:
+    """Encode sources given as unit (or non-unit) Cartesian directions."""
+    d = np.atleast_2d(np.asarray(directions_xyz, dtype=np.float64))
+    Y = sh_sn3d_batch(d, max_order=max_order)
+    if gains is None:
+        g = np.ones(d.shape[0], dtype=np.float64)
+    else:
+        g = np.asarray(gains, dtype=np.float64).reshape(-1)
+    return np.einsum("n,nc->c", g, Y)

hoa64/live_audio.py

@@ -0,0 +1,240 @@
+"""Live microphone / PulseAudio capture → HOA spatial report.
+
+Uses system tools (no sounddevice):
+  * ffmpeg -f pulse -i <source>
+  * arecord (ALSA fallback)
+
+Mono capture is encoded as a plane wave from a configured look direction
+(default: front). Multi-channel (2–4) maps L/R to ±az pseudo-Ambix order-1.
+"""
+
+from __future__ import annotations
+
+import shutil
+import subprocess
+import tempfile
+from pathlib import Path
+from typing import Optional, Sequence, Union
+
+import numpy as np
+
+from .audio_io import read_wav, write_wav
+from .encode import encode_plane_waves
+from .report import SpatialReport, report_from_hoa
+from .stream import encode_mono_plane_wave
+
+PathLike = Union[str, Path]
+
+
+def list_pulse_sources() -> list[str]:
+    """Return Pulse/PipeWire source names via pactl."""
+    if not shutil.which("pactl"):
+        return []
+    try:
+        out = subprocess.check_output(
+            ["pactl", "list", "short", "sources"],
+            text=True,
+            stderr=subprocess.DEVNULL,
+        )
+    except subprocess.CalledProcessError:
+        return []
+    names = []
+    for line in out.splitlines():
+        parts = line.split("\t")
+        if len(parts) >= 2:
+            names.append(parts[1])
+    return names
+
+
+def default_mic_source() -> str:
+    """Prefer a non-monitor mic source."""
+    sources = list_pulse_sources()
+    for s in sources:
+        if "monitor" not in s.lower() and ("mic" in s.lower() or "input" in s.lower() or "Line" in s):
+            return s
+    for s in sources:
+        if "monitor" not in s.lower():
+            return s
+    return "default"
+
+
+def capture_wav(
+    path: PathLike,
+    *,
+    duration_sec: float = 2.0,
+    sample_rate: int = 48000,
+    channels: int = 1,
+    source: Optional[str] = None,
+    backend: str = "auto",
+) -> Path:
+    """Record audio to WAV. backend: auto | ffmpeg | arecord."""
+    path = Path(path)
+    path.parent.mkdir(parents=True, exist_ok=True)
+    source = source or default_mic_source()
+    backend = backend.lower()
+
+    def _ffmpeg() -> None:
+        # pulse device
+        dev = source if source != "default" else "default"
+        cmd = [
+            "ffmpeg", "-y", "-hide_banner", "-loglevel", "error",
+            "-f", "pulse", "-i", dev,
+            "-t", str(duration_sec),
+            "-ar", str(sample_rate),
+            "-ac", str(channels),
+            str(path),
+        ]
+        subprocess.check_call(cmd)
+
+    def _arecord() -> None:
+        cmd = [
+            "arecord",
+            "-d", str(int(max(1, round(duration_sec)))),
+            "-f", "S16_LE",
+            "-r", str(sample_rate),
+            "-c", str(channels),
+            str(path),
+        ]
+        subprocess.check_call(cmd)
+
+    errors = []
+    if backend in ("auto", "ffmpeg") and shutil.which("ffmpeg"):
+        try:
+            _ffmpeg()
+            return path
+        except Exception as e:
+            errors.append(f"ffmpeg: {e}")
+            if backend == "ffmpeg":
+                raise
+    if backend in ("auto", "arecord") and shutil.which("arecord"):
+        try:
+            _arecord()
+            return path
+        except Exception as e:
+            errors.append(f"arecord: {e}")
+            if backend == "arecord":
+                raise
+    raise RuntimeError("capture failed: " + "; ".join(errors) if errors else "no capture backend")
+
+
+def stereo_to_order1_hoa(
+    left: np.ndarray,
+    right: np.ndarray,
+    *,
+    width_az_deg: float = 30.0,
+) -> np.ndarray:
+    """Pseudo order-1 Ambix from stereo: L at +width, R at -width, mid omni."""
+    left = np.asarray(left, dtype=np.float64).reshape(-1)
+    right = np.asarray(right, dtype=np.float64).reshape(-1)
+    n = min(left.shape[0], right.shape[0])
+    left, right = left[:n], right[:n]
+    mid = 0.5 * (left + right)
+    # plane waves
+    L = encode_plane_waves([width_az_deg], [0.0], left[None, :], max_order=1)
+    R = encode_plane_waves([-width_az_deg], [0.0], right[None, :], max_order=1)
+    M = encode_plane_waves([0.0], [0.0], mid[None, :], max_order=1)
+    hoa = L + R + 0.5 * M
+    # pad to 64
+    out = np.zeros((64, n), dtype=np.float64)
+    out[:4] = hoa[:4]
+    return out
+
+
+def audio_to_hoa_stream(
+    audio: np.ndarray,
+    *,
+    az_deg: float = 0.0,
+    el_deg: float = 0.0,
+    max_order: int = 3,
+    stereo_width_az: float = 30.0,
+) -> np.ndarray:
+    """(C,T) HOA from captured audio array (C_in, T) or (T,)."""
+    a = np.asarray(audio, dtype=np.float64)
+    if a.ndim == 1:
+        return encode_mono_plane_wave(a, az_deg, el_deg, max_order=max_order)
+    if a.shape[0] == 1:
+        return encode_mono_plane_wave(a[0], az_deg, el_deg, max_order=max_order)
+    if a.shape[0] >= 4:
+        # assume Ambix ACN already
+        from .audio_io import ensure_hoa_channels
+        return ensure_hoa_channels(a, max_order=max_order)
+    if a.shape[0] == 2:
+        hoa = stereo_to_order1_hoa(a[0], a[1], width_az_deg=stereo_width_az)
+        if max_order < 1:
+            return hoa[:1]
+        return hoa
+    # mixdown
+    mono = np.mean(a, axis=0)
+    return encode_mono_plane_wave(mono, az_deg, el_deg, max_order=max_order)
+
+
+def live_report(
+    *,
+    duration_sec: float = 2.0,
+    sample_rate: int = 48000,
+    channels: int = 1,
+    source: Optional[str] = None,
+    az_deg: float = 0.0,
+    el_deg: float = 0.0,
+    max_order: int = 3,
+    keep_wav: Optional[PathLike] = None,
+) -> SpatialReport:
+    """Capture from mic and return a SpatialReport."""
+    if keep_wav:
+        wav_path = Path(keep_wav)
+    else:
+        tmp = tempfile.NamedTemporaryFile(suffix=".wav", delete=False)
+        wav_path = Path(tmp.name)
+        tmp.close()
+    try:
+        capture_wav(
+            wav_path,
+            duration_sec=duration_sec,
+            sample_rate=sample_rate,
+            channels=channels,
+            source=source,
+        )
+        audio, sr = read_wav(wav_path)
+        hoa = audio_to_hoa_stream(
+            audio, az_deg=az_deg, el_deg=el_deg, max_order=max_order
+        )
+        rep = report_from_hoa(
+            hoa,
+            sr,
+            max_order=max_order,
+            notes=[
+                f"live capture {duration_sec}s ch={audio.shape[0]} src={source or 'default'}",
+                f"encode az={az_deg} el={el_deg}" if audio.shape[0] < 4 else "ambix/multi-ch path",
+            ],
+            meta={
+                "encode": "live_mic",
+                "source": source or default_mic_source(),
+                "wav": str(wav_path) if keep_wav else None,
+            },
+        )
+        return rep
+    finally:
+        if not keep_wav:
+            try:
+                wav_path.unlink(missing_ok=True)
+            except Exception:
+                pass
+
+
+def live_report_from_file(
+    path: PathLike,
+    *,
+    az_deg: float = 0.0,
+    el_deg: float = 0.0,
+    max_order: int = 3,
+) -> SpatialReport:
+    """Same pipeline as live_report but from an existing WAV (offline test)."""
+    audio, sr = read_wav(path)
+    hoa = audio_to_hoa_stream(audio, az_deg=az_deg, el_deg=el_deg, max_order=max_order)
+    return report_from_hoa(
+        hoa,
+        sr,
+        max_order=max_order,
+        notes=[f"live pipeline on file {path}"],
+        meta={"encode": "live_from_file", "path": str(path)},
+    )

hoa64/report.py

@@ -0,0 +1,299 @@
+"""JSON spatial report — compact payload for other models / tools."""
+
+from __future__ import annotations
+
+import json
+import time
+from dataclasses import asdict, dataclass, field
+from pathlib import Path
+from typing import Any, Optional, Sequence, Union
+
+import numpy as np
+
+from .analysis import angular_error_deg, doa_from_intensity, field_energy, peak_direction
+from .audio_io import ensure_hoa_channels, read_wav
+from .basis import MAX_ORDER, N_CHANNELS, channel_names
+from .stream import (
+    FrameAnalysis,
+    SourceSpec,
+    analyze_hoa_frames,
+    analyze_hoa_stft_bands,
+    encode_mono_plane_wave,
+    encode_scene,
+    hoa_rms,
+)
+
+PathLike = Union[str, Path]
+
+REPORT_SCHEMA_VERSION = "spatial-hoa.report.v1"
+
+
+@dataclass
+class SpatialReport:
+    """Machine-readable spatial summary for LLM / diffusion conditioning."""
+
+    schema: str = REPORT_SCHEMA_VERSION
+    kind: str = "spatial_field"
+    sample_rate: int = 0
+    n_samples: int = 0
+    duration_sec: float = 0.0
+    max_order: int = MAX_ORDER
+    n_channels: int = N_CHANNELS
+    # Global broadband
+    energy: float = 0.0
+    doa_az_deg: float = 0.0
+    doa_el_deg: float = 0.0
+    peak_az_deg: float = 0.0
+    peak_el_deg: float = 0.0
+    peak_power: float = 0.0
+    # Optional detail
+    bands: list[dict] = field(default_factory=list)
+    frames: list[dict] = field(default_factory=list)
+    sources_hint: list[dict] = field(default_factory=list)
+    channel_rms: list[float] = field(default_factory=list)
+    notes: list[str] = field(default_factory=list)
+    meta: dict = field(default_factory=dict)
+
+    def to_dict(self) -> dict[str, Any]:
+        return asdict(self)
+
+    def to_json(self, *, indent: int | None = 2) -> str:
+        return json.dumps(self.to_dict(), indent=indent)
+
+    def save(self, path: PathLike) -> None:
+        Path(path).write_text(self.to_json() + "\n", encoding="utf-8")
+
+    def one_liner(self) -> str:
+        """Short string suitable for tool_result / system hints."""
+        return (
+            f"spatial: DOA az={self.doa_az_deg:.1f}° el={self.doa_el_deg:.1f}° "
+            f"peak=({self.peak_az_deg:.1f},{self.peak_el_deg:.1f}) "
+            f"E={self.energy:.4g} T={self.duration_sec:.2f}s"
+        )
+
+
+def _frames_to_dicts(frames: Sequence[FrameAnalysis]) -> list[dict]:
+    return [
+        {
+            "t_sec": f.t_center_sec,
+            "energy": f.energy,
+            "doa_az_deg": f.doa_az_deg,
+            "doa_el_deg": f.doa_el_deg,
+            "peak_az_deg": f.peak_az_deg,
+            "peak_el_deg": f.peak_el_deg,
+            "order1_energy": f.order1_energy,
+        }
+        for f in frames
+    ]
+
+
+def report_from_hoa(
+    hoa: np.ndarray,
+    sample_rate: int,
+    *,
+    max_order: int = MAX_ORDER,
+    frame_ms: float = 40.0,
+    hop_ms: float = 20.0,
+    include_frames: bool = True,
+    include_bands: bool = True,
+    include_peak_map: bool = True,
+    max_frames: int = 64,
+    meta: Optional[dict] = None,
+    notes: Optional[list[str]] = None,
+) -> SpatialReport:
+    """Build a SpatialReport from HOA coefficients (C,T) or (C,)."""
+    a = np.asarray(hoa, dtype=np.float64)
+    if a.ndim == 1:
+        a = a.reshape(-1, 1)
+    a = ensure_hoa_channels(a, max_order=max_order)
+    n_samples = int(a.shape[1])
+    duration = n_samples / float(sample_rate) if sample_rate else 0.0
+
+    # Broadband DOA must be AC-safe (time-mean of coeffs → 0 for audio).
+    # 1) global intensity from products W*X, W*Y, W*Z
+    # 2) energy-weighted blend of per-frame DOAs
+    from .basis import az_el_from_unit, unit_vector
+
+    W, Yc, Zc, Xc = a[0], a[1], a[2], a[3]
+    I_global = np.array(
+        [
+            float(np.mean(W * Xc)),
+            float(np.mean(W * Yc)),
+            float(np.mean(W * Zc)),
+        ]
+    )
+    nI = float(np.linalg.norm(I_global))
+    if nI > 1e-18:
+        az, el = az_el_from_unit(I_global / nI, degrees=True)
+        az, el = float(az), float(el)
+    else:
+        az, el = 0.0, 0.0
+
+    frames = analyze_hoa_frames(
+        a,
+        sample_rate,
+        frame_ms=frame_ms,
+        hop_ms=hop_ms,
+        max_order=max_order,
+        peak_grid=False,
+    )
+    if frames:
+        vecs = []
+        for f in frames:
+            if f.energy <= 1e-18:
+                continue
+            u = unit_vector(f.doa_az_deg, f.doa_el_deg, degrees=True)
+            vecs.append(u * f.energy)
+        if vecs:
+            v = np.sum(vecs, axis=0)
+            n = float(np.linalg.norm(v))
+            if n > 1e-15:
+                az_f, el_f = az_el_from_unit(v / n, degrees=True)
+                # Prefer frame blend when global intensity is weak
+                if nI < 1e-12:
+                    az, el = float(az_f), float(el_f)
+                else:
+                    # average unit vectors of global + frames
+                    ug = unit_vector(az, el, degrees=True)
+                    uf = unit_vector(float(az_f), float(el_f), degrees=True)
+                    um = ug + uf
+                    nm = float(np.linalg.norm(um))
+                    if nm > 1e-15:
+                        az, el = az_el_from_unit(um / nm, degrees=True)
+                        az, el = float(az), float(el)
+
+    energy = float(np.mean(np.sum(a * a, axis=0)))
+
+    if include_peak_map:
+        # Pseudo-static field from RMS * sign(corr with W) for map peak
+        rms = np.sqrt(np.mean(a * a, axis=1) + 1e-30)
+        sign = np.sign(np.mean(a * a[0:1, :], axis=1) + 1e-30)
+        pseudo = rms * sign
+        paz, pel, pv = peak_direction(
+            pseudo, n_azi=72, n_el=36, max_order=min(max_order, 3)
+        )
+    else:
+        paz, pel, pv = az, el, energy
+
+    bands: list[dict] = []
+    if include_bands and n_samples >= 256:
+        bands = analyze_hoa_stft_bands(a, sample_rate)
+
+    frame_dicts: list[dict] = []
+    if include_frames and frames:
+        step = max(1, len(frames) // max_frames)
+        frame_dicts = _frames_to_dicts(frames[::step][:max_frames])
+
+    return SpatialReport(
+        sample_rate=int(sample_rate),
+        n_samples=n_samples,
+        duration_sec=duration,
+        max_order=max_order,
+        n_channels=(max_order + 1) ** 2,
+        energy=energy,
+        doa_az_deg=float(az),
+        doa_el_deg=float(el),
+        peak_az_deg=float(paz),
+        peak_el_deg=float(pel),
+        peak_power=float(pv),
+        bands=bands,
+        frames=frame_dicts,
+        channel_rms=[float(x) for x in hoa_rms(a)[:8]],  # first 8 for brevity
+        notes=list(notes or []),
+        meta=dict(meta or {}),
+    )
+
+
+def report_from_mono_wav(
+    path: PathLike,
+    azimuth_deg: float,
+    elevation_deg: float = 0.0,
+    *,
+    max_order: int = MAX_ORDER,
+    **kwargs,
+) -> SpatialReport:
+    """Load mono (or mixdown) WAV, encode as plane wave from known direction, report."""
+    audio, sr = read_wav(path)
+    if audio.shape[0] == 1:
+        mono = audio[0]
+    else:
+        mono = np.mean(audio, axis=0)
+        notes_extra = ["mixed multi-channel input down to mono before plane-wave encode"]
+    hoa = encode_mono_plane_wave(mono, azimuth_deg, elevation_deg, max_order=max_order)
+    notes = list(kwargs.pop("notes", []) or [])
+    if audio.shape[0] > 1:
+        notes.extend(notes_extra)
+    notes.append(
+        f"encoded mono as plane wave az={azimuth_deg}° el={elevation_deg}°"
+    )
+    meta = dict(kwargs.pop("meta", None) or {})
+    meta.update(
+        {
+            "source_path": str(path),
+            "encode": "mono_plane_wave",
+            "encode_az_deg": azimuth_deg,
+            "encode_el_deg": elevation_deg,
+            "created_unix": time.time(),
+        }
+    )
+    return report_from_hoa(
+        hoa, sr, max_order=max_order, notes=notes, meta=meta, **kwargs
+    )
+
+
+def report_from_ambix_wav(
+    path: PathLike,
+    *,
+    max_order: int = MAX_ORDER,
+    **kwargs,
+) -> SpatialReport:
+    """Load multi-channel Ambix WAV (ACN order) and analyze."""
+    audio, sr = read_wav(path)
+    hoa = ensure_hoa_channels(audio, max_order=max_order)
+    notes = list(kwargs.pop("notes", []) or [])
+    notes.append(f"loaded Ambix-style multi-channel WAV with {audio.shape[0]} ch")
+    meta = dict(kwargs.pop("meta", None) or {})
+    meta.update(
+        {
+            "source_path": str(path),
+            "encode": "ambix_wav",
+            "input_channels": int(audio.shape[0]),
+            "created_unix": time.time(),
+        }
+    )
+    return report_from_hoa(
+        hoa, sr, max_order=max_order, notes=notes, meta=meta, **kwargs
+    )
+
+
+def report_from_scene(
+    sources: Sequence[SourceSpec],
+    sample_rate: int,
+    *,
+    max_order: int = MAX_ORDER,
+    **kwargs,
+) -> SpatialReport:
+    """Analyze a synthetic multi-source scene."""
+    hoa = encode_scene(sources, max_order=max_order)
+    hints = [
+        {
+            "label": s.label or f"src{i}",
+            "azimuth_deg": s.azimuth_deg,
+            "elevation_deg": s.elevation_deg,
+            "rms": float(np.sqrt(np.mean(np.asarray(s.signal) ** 2) + 1e-30)),
+        }
+        for i, s in enumerate(sources)
+    ]
+    notes = list(kwargs.pop("notes", []) or [])
+    notes.append(f"synthetic scene with {len(sources)} plane-wave source(s)")
+    meta = dict(kwargs.pop("meta", None) or {})
+    meta["encode"] = "synthetic_scene"
+    rep = report_from_hoa(
+        hoa, sample_rate, max_order=max_order, notes=notes, meta=meta, **kwargs
+    )
+    rep.sources_hint = hints
+    return rep
+
+
+def load_report(path: PathLike) -> dict:
+    return json.loads(Path(path).read_text(encoding="utf-8"))

hoa64/rnn_stub.py

@@ -0,0 +1,101 @@
+"""RNN integrator stub — Phase 0 interface for iterative motion in a field.
+
+No learned weights yet: explicit Euler integration of pose + field rotation.
+Validates the "calculator + loop" hypothesis before training dynamics.
+"""
+
+from __future__ import annotations
+
+from dataclasses import dataclass, field
+
+import numpy as np
+
+from .basis import N_CHANNELS
+from .encode import encode_points
+from .rotate import rotate_matrix_order1, rotate_yaw_pitch_roll
+from .analysis import doa_from_intensity, peak_direction, field_energy
+
+
+def _pad64(a: np.ndarray) -> np.ndarray:
+    a = np.asarray(a, dtype=np.float64).reshape(-1)
+    out = np.zeros(N_CHANNELS, dtype=np.float64)
+    out[: min(N_CHANNELS, a.shape[0])] = a[:N_CHANNELS]
+    return out
+
+
+@dataclass
+class SpatialState:
+    """Agent + field state for iterative spatial calculation."""
+
+    hoa: np.ndarray  # (64,) world field in listener frame after last step
+    yaw: float = 0.0  # degrees, agent heading
+    pitch: float = 0.0
+    roll: float = 0.0
+    history: list = field(default_factory=list)
+
+    def __post_init__(self) -> None:
+        self.hoa = _pad64(self.hoa)
+
+    def snapshot(self) -> dict:
+        az, el = doa_from_intensity(self.hoa)
+        paz, pel, pval = peak_direction(self.hoa, n_azi=72, n_el=36)
+        return {
+            "yaw": self.yaw,
+            "pitch": self.pitch,
+            "roll": self.roll,
+            "energy": field_energy(self.hoa),
+            "doa_intensity_az_el": (az, el),
+            "doa_peak_az_el": (paz, pel),
+            "peak_power": pval,
+        }
+
+
+def step_rotate(
+    state: SpatialState,
+    d_yaw: float = 0.0,
+    d_pitch: float = 0.0,
+    d_roll: float = 0.0,
+    *,
+    max_order: int = 7,
+    dense: bool = True,
+) -> SpatialState:
+    """Integrate a pose increment: rotate the field opposite agent turn.
+
+    If the agent yaws +θ (turns left), the world field in head frame yaws −θ.
+    """
+    # Agent pose update
+    yaw = state.yaw + d_yaw
+    pitch = state.pitch + d_pitch
+    roll = state.roll + d_roll
+    # Field in listener frame: rotate by -d_*
+    if dense and max_order > 1:
+        hoa = rotate_yaw_pitch_roll(
+            state.hoa,
+            yaw=-d_yaw,
+            pitch=-d_pitch,
+            roll=-d_roll,
+            degrees=True,
+            max_order=max_order,
+        )
+    else:
+        hoa = rotate_matrix_order1(
+            state.hoa, yaw=-d_yaw, pitch=-d_pitch, roll=-d_roll, degrees=True
+        )
+    hoa = _pad64(hoa)
+    new = SpatialState(hoa=hoa, yaw=yaw, pitch=pitch, roll=roll, history=list(state.history))
+    new.history.append(new.snapshot())
+    return new
+
+
+def world_from_sources(
+    azimuths,
+    elevations,
+    gains=None,
+    *,
+    max_order: int = 7,
+) -> SpatialState:
+    """Build initial state from world-frame sources (listener at origin, identity pose)."""
+    a = encode_points(azimuths, elevations, gains, degrees=True, max_order=max_order)
+    st = SpatialState(hoa=a)
+    st.history.append(st.snapshot())
+    return st

hoa64/rotate.py

@@ -0,0 +1,186 @@
+"""Rotate HOA fields.
+
+Phase 2 default: Wigner-D block rotation (fast, all orders).
+Phase 0 reference: dense spherical re-projection (``method="dense"``).
+Order-1 shortcut: ``rotate_matrix_order1`` (exact Cartesian).
+"""
+
+from __future__ import annotations
+
+import math
+
+import numpy as np
+
+from .basis import MAX_ORDER, N_CHANNELS, sh_sn3d_batch, sphere_grid, unit_vector, acn_index
+from .wigner import (
+    apply_hoa_rotation,
+    hoa_rotation_matrix,
+    rotation_matrix_zyx,
+)
+
+
+def _rotation_matrix_zyx(yaw: float, pitch: float, roll: float, degrees: bool) -> np.ndarray:
+    return rotation_matrix_zyx(yaw, pitch, roll, degrees=degrees)
+
+
+def rotate_matrix_order1(
+    hoa: np.ndarray,
+    yaw: float = 0.0,
+    pitch: float = 0.0,
+    roll: float = 0.0,
+    *,
+    degrees: bool = True,
+) -> np.ndarray:
+    """Exact rotation for order-1 channels (Y,Z,X); W and n≥2 unchanged."""
+    a = np.array(hoa, dtype=np.float64, copy=True)
+    R = _rotation_matrix_zyx(yaw, pitch, roll, degrees)
+    if a.ndim == 1:
+        cart = np.array([a[3], a[1], a[2]])
+        Xp, Yp, Zp = R @ cart
+        a[3], a[1], a[2] = Xp, Yp, Zp
+        return a
+    cart = np.stack([a[3], a[1], a[2]], axis=0)
+    rot = R @ cart
+    a[3], a[1], a[2] = rot[0], rot[1], rot[2]
+    return a
+
+
+def _sn3d_channel_scale(max_order: int) -> np.ndarray:
+    nch = (max_order + 1) ** 2
+    scale = np.empty(nch, dtype=np.float64)
+    for n in range(max_order + 1):
+        s = float(2 * n + 1)
+        for m in range(-n, n + 1):
+            scale[acn_index(n, m)] = s
+    return scale
+
+
+def rotate_yaw_pitch_roll(
+    hoa: np.ndarray,
+    yaw: float = 0.0,
+    pitch: float = 0.0,
+    roll: float = 0.0,
+    *,
+    degrees: bool = True,
+    max_order: int = MAX_ORDER,
+    method: str = "wigner",
+    n_azi: int = 96,
+    n_el: int = 48,
+) -> np.ndarray:
+    """Rotate HOA field (active): values move with R.
+
+    Parameters
+    ----------
+    method :
+        ``"wigner"`` (default, Phase 2) — exact band-limited SH rotation.
+        ``"dense"`` — Phase 0 spherical sample / re-encode (slow reference).
+    """
+    if method == "wigner":
+        R = _rotation_matrix_zyx(yaw, pitch, roll, degrees)
+        M = hoa_rotation_matrix(R, max_order=max_order)
+        return apply_hoa_rotation(hoa, M)
+    if method == "dense":
+        return _rotate_dense(
+            hoa,
+            yaw,
+            pitch,
+            roll,
+            degrees=degrees,
+            max_order=max_order,
+            n_azi=n_azi,
+            n_el=n_el,
+        )
+    raise ValueError("method must be 'wigner' or 'dense'")
+
+
+def _rotate_dense(
+    hoa: np.ndarray,
+    yaw: float,
+    pitch: float,
+    roll: float,
+    *,
+    degrees: bool,
+    max_order: int,
+    n_azi: int,
+    n_el: int,
+) -> np.ndarray:
+    """Phase 0 dense re-encode (reference)."""
+    a = np.asarray(hoa, dtype=np.float64)
+    if a.ndim == 2:
+        out = np.zeros(((max_order + 1) ** 2, a.shape[1]), dtype=np.float64)
+        for t in range(a.shape[1]):
+            out[:, t] = _rotate_dense(
+                a[:, t],
+                yaw,
+                pitch,
+                roll,
+                degrees=degrees,
+                max_order=max_order,
+                n_azi=n_azi,
+                n_el=n_el,
+            )
+        full = np.zeros((N_CHANNELS, a.shape[1]), dtype=np.float64)
+        full[: out.shape[0]] = out
+        return full
+
+    from .decode import decode_directions
+
+    nch = (max_order + 1) ** 2
+    a = a[:nch]
+    R = _rotation_matrix_zyx(yaw, pitch, roll, degrees)
+
+    azi, el, weights = sphere_grid(n_azi, n_el, degrees=True)
+    AA, EE = np.meshgrid(azi, el, indexing="ij")
+    dirs = unit_vector(AA, EE, degrees=True).reshape(-1, 3)
+    w = weights.reshape(-1)
+    f = decode_directions(
+        a, AA.reshape(-1), EE.reshape(-1), degrees=True, max_order=max_order
+    )
+    dirs_rot = (R @ dirs.T).T
+    Y_rot = sh_sn3d_batch(dirs_rot, max_order=max_order)
+    a_new = np.einsum("i,i,ic->c", w, f, Y_rot)
+    a_new *= _sn3d_channel_scale(max_order) / (4.0 * math.pi)
+    full = np.zeros(N_CHANNELS, dtype=np.float64)
+    full[: a_new.shape[0]] = a_new
+    return full
+
+
+def rotate_source_directions(
+    azimuths: np.ndarray,
+    elevations: np.ndarray,
+    yaw: float = 0.0,
+    pitch: float = 0.0,
+    roll: float = 0.0,
+    *,
+    degrees: bool = True,
+) -> tuple[np.ndarray, np.ndarray]:
+    """Rotate source az/el by the same R (exact plane-wave ground truth)."""
+    dirs = unit_vector(azimuths, elevations, degrees=degrees)
+    flat = np.atleast_2d(dirs.reshape(-1, 3))
+    R = _rotation_matrix_zyx(yaw, pitch, roll, degrees)
+    rot = (R @ flat.T).T
+    from .basis import az_el_from_unit
+
+    az, el = az_el_from_unit(rot, degrees=degrees)
+    return az.reshape(np.shape(azimuths)), el.reshape(np.shape(elevations))
+
+
+def rotate_plane_wave_field(
+    hoa: np.ndarray,
+    yaw: float = 0.0,
+    pitch: float = 0.0,
+    roll: float = 0.0,
+    *,
+    degrees: bool = True,
+    max_order: int = MAX_ORDER,
+    method: str = "wigner",
+) -> np.ndarray:
+    return rotate_yaw_pitch_roll(
+        hoa,
+        yaw,
+        pitch,
+        roll,
+        degrees=degrees,
+        max_order=max_order,
+        method=method,
+    )

hoa64/server.py

@@ -0,0 +1,369 @@
+"""Minimal HTTP API so Qwythos / agents can request spatial reports.
+
+Stdlib only (no Flask). Default bind: 127.0.0.1:8765
+
+Endpoints
+---------
+GET  /health
+GET  /v1/spatial/schema          — OpenAI-style tool descriptor
+POST /v1/spatial/analyze        — JSON body → SpatialReport
+POST /v1/spatial/analyze_file   — {path, mode, az, el, ...}
+POST /v1/spatial/demo_scene     — synthetic multi-source report
+POST /v1/spatial/vision         — Phase 3: boxes/rays → spatial report
+POST /v1/spatial/fuse           — merge audio + vision report dicts
+
+OpenAI tools (for agents with bash): also install ``spatial-report`` on PATH.
+"""
+
+from __future__ import annotations
+
+import json
+import traceback
+from http.server import BaseHTTPRequestHandler, ThreadingHTTPServer
+from pathlib import Path
+from typing import Any
+from urllib.parse import urlparse
+
+from . import __version__
+from .report import (
+    REPORT_SCHEMA_VERSION,
+    report_from_ambix_wav,
+    report_from_hoa,
+    report_from_mono_wav,
+    report_from_scene,
+)
+from .stream import SourceSpec
+from .synth import envelope_adsr, tone
+
+DEFAULT_HOST = "127.0.0.1"
+DEFAULT_PORT = 8765
+
+TOOL_SCHEMA = {
+    "type": "function",
+    "function": {
+        "name": "spatial_analyze",
+        "description": (
+            "Analyze spatial audio (Ambix HOA or mono plane-wave) or visual "
+            "bounding boxes on the sphere. Returns a compact spatial report "
+            "(DOA, bands, frames) from the HOA-7 calculator — not an LLM."
+        ),
+        "parameters": {
+            "type": "object",
+            "properties": {
+                "mode": {
+                    "type": "string",
+                    "enum": ["ambix_file", "mono_file", "demo_scene", "vision", "fuse"],
+                    "description": "Analysis mode",
+                },
+                "path": {
+                    "type": "string",
+                    "description": "WAV path for ambix_file or mono_file",
+                },
+                "az": {
+                    "type": "number",
+                    "description": "Azimuth degrees for mono plane-wave encode (0=front, +90=left)",
+                },
+                "el": {
+                    "type": "number",
+                    "description": "Elevation degrees (0=horizon, +90=zenith)",
+                },
+                "order": {
+                    "type": "integer",
+                    "description": "Max HOA order 0..7 (default 3 for speed, 7 full)",
+                    "default": 3,
+                },
+                "boxes": {
+                    "type": "array",
+                    "description": "Vision boxes: [{az, el, w_deg?, h_deg?, weight?, label?}]",
+                    "items": {"type": "object"},
+                },
+                "audio_report": {
+                    "type": "object",
+                    "description": "Existing audio SpatialReport dict for fuse mode",
+                },
+                "vision_report": {
+                    "type": "object",
+                    "description": "Existing vision SpatialReport dict for fuse mode",
+                },
+            },
+            "required": ["mode"],
+        },
+    },
+}
+
+
+def _json_response(handler: BaseHTTPRequestHandler, code: int, obj: Any) -> None:
+    body = json.dumps(obj, indent=None).encode("utf-8")
+    handler.send_response(code)
+    handler.send_header("Content-Type", "application/json")
+    handler.send_header("Content-Length", str(len(body)))
+    handler.send_header("Access-Control-Allow-Origin", "*")
+    handler.end_headers()
+    handler.wfile.write(body)
+
+
+def _read_json(handler: BaseHTTPRequestHandler) -> dict:
+    n = int(handler.headers.get("Content-Length", "0"))
+    raw = handler.rfile.read(n) if n else b"{}"
+    if not raw:
+        return {}
+    return json.loads(raw.decode("utf-8"))
+
+
+def handle_analyze(body: dict) -> dict:
+    mode = body.get("mode", "demo_scene")
+    order = int(body.get("order", 3))
+    order = max(0, min(7, order))
+
+    if mode == "demo_scene":
+        sr = int(body.get("sample_rate", 48000))
+        dur = float(body.get("duration", 0.4))
+        n = int(sr * dur)
+        env = envelope_adsr(n, sr)
+        sources = [
+            SourceSpec(0.0, 0.0, tone(440, dur, sr, amplitude=0.4) * env, "front"),
+            SourceSpec(90.0, 10.0, tone(660, dur, sr, amplitude=0.25) * env, "left"),
+        ]
+        rep = report_from_scene(sources, sr, max_order=order)
+        d = rep.to_dict()
+        d["one_liner"] = rep.one_liner()
+        return d
+
+    if mode == "ambix_file":
+        path = body.get("path") or body.get("file")
+        if not path:
+            raise ValueError("path required for ambix_file")
+        rep = report_from_ambix_wav(path, max_order=order)
+        d = rep.to_dict()
+        d["one_liner"] = rep.one_liner()
+        return d
+
+    if mode == "mono_file":
+        path = body.get("path") or body.get("file")
+        if not path:
+            raise ValueError("path required for mono_file")
+        if body.get("az") is None:
+            raise ValueError("az required for mono_file plane-wave encode")
+        rep = report_from_mono_wav(
+            path,
+            float(body["az"]),
+            float(body.get("el", 0.0)),
+            max_order=order,
+        )
+        d = rep.to_dict()
+        d["one_liner"] = rep.one_liner()
+        return d
+
+    if mode == "vision":
+        from .vision import report_from_boxes
+
+        boxes = body.get("boxes") or []
+        rep = report_from_boxes(boxes, max_order=order)
+        d = rep.to_dict()
+        d["one_liner"] = rep.one_liner()
+        return d
+
+    if mode == "fuse":
+        from .vision import fuse_reports
+
+        ar = body.get("audio_report") or {}
+        vr = body.get("vision_report") or {}
+        return fuse_reports(ar, vr)
+
+    if mode == "detect":
+        from .detector import detect_to_sphere, load_boxes_json
+        from .vision import report_from_boxes
+
+        if body.get("boxes_path"):
+            boxes = load_boxes_json(body["boxes_path"])
+        elif body.get("image"):
+            boxes = detect_to_sphere(
+                body["image"],
+                backend=body.get("backend", "auto"),
+                score_thresh=float(body.get("score", 0.5)),
+                hfov_deg=float(body.get("hfov", 90)),
+                vfov_deg=float(body.get("vfov", 60)),
+            )
+        else:
+            boxes = body.get("boxes") or []
+        rep = report_from_boxes(boxes, max_order=order)
+        d = rep.to_dict()
+        d["one_liner"] = rep.one_liner()
+        d["boxes"] = boxes
+        return d
+
+    if mode == "live":
+        from .live_audio import live_report
+
+        rep = live_report(
+            duration_sec=float(body.get("duration", 2.0)),
+            sample_rate=int(body.get("sample_rate", 48000)),
+            channels=int(body.get("channels", 1)),
+            source=body.get("source"),
+            az_deg=float(body.get("az", 0.0)),
+            el_deg=float(body.get("el", 0.0)),
+            max_order=order,
+            keep_wav=body.get("write_wav"),
+        )
+        d = rep.to_dict()
+        d["one_liner"] = rep.one_liner()
+        return d
+
+    if mode == "panner":
+        # UI spherical panner: az/el + W amplitude → spatial report (+ optional condition)
+        from .conditioning import build_conditioning, panner_report
+
+        rep = panner_report(
+            float(body.get("az", body.get("az_deg", 0.0))),
+            float(body.get("el", body.get("el_deg", 0.0))),
+            float(body.get("w", body.get("w_amplitude", body.get("energy", 0.5)))),
+        )
+        if body.get("condition") or body.get("prompt"):
+            return build_conditioning(
+                rep,
+                base_prompt=str(body.get("prompt", "")),
+                style=str(body.get("style", "natural")),
+            )
+        return rep
+
+    if mode == "condition":
+        from .conditioning import build_conditioning, panner_report
+
+        rep = body.get("report")
+        if rep is None and (
+            "az" in body or "az_deg" in body or "w" in body or "w_amplitude" in body
+        ):
+            # Convenience: condition directly from panner knobs
+            rep = panner_report(
+                float(body.get("az", body.get("az_deg", 0.0))),
+                float(body.get("el", body.get("el_deg", 0.0))),
+                float(body.get("w", body.get("w_amplitude", body.get("energy", 0.5)))),
+            )
+        else:
+            rep = rep or body
+        return build_conditioning(
+            rep,
+            base_prompt=str(body.get("prompt", "")),
+            style=str(body.get("style", "natural")),
+        )
+
+    if mode == "hoa_vector":
+        # raw coefficients list
+        import numpy as np
+
+        coeffs = body.get("hoa") or body.get("coefficients")
+        if coeffs is None:
+            raise ValueError("hoa coefficients required")
+        sr = int(body.get("sample_rate", 48000))
+        a = np.asarray(coeffs, dtype=np.float64)
+        if a.ndim == 1:
+            a = a.reshape(-1, 1)
+        rep = report_from_hoa(a, sr, max_order=order)
+        d = rep.to_dict()
+        d["one_liner"] = rep.one_liner()
+        return d
+
+    raise ValueError(f"unknown mode: {mode}")
+
+
+class SpatialHandler(BaseHTTPRequestHandler):
+    server_version = f"spatial-hoa/{__version__}"
+
+    def log_message(self, fmt: str, *args) -> None:
+        # quieter default
+        sys_stderr = __import__("sys").stderr
+        sys_stderr.write("%s - %s\n" % (self.address_string(), fmt % args))
+
+    def do_OPTIONS(self) -> None:
+        self.send_response(204)
+        self.send_header("Access-Control-Allow-Origin", "*")
+        self.send_header("Access-Control-Allow-Methods", "GET, POST, OPTIONS")
+        self.send_header("Access-Control-Allow-Headers", "Content-Type")
+        self.end_headers()
+
+    def do_GET(self) -> None:
+        path = urlparse(self.path).path
+        if path in ("/health", "/v1/health"):
+            _json_response(
+                self,
+                200,
+                {
+                    "status": "ok",
+                    "service": "spatial-hoa",
+                    "version": __version__,
+                    "schema": REPORT_SCHEMA_VERSION,
+                },
+            )
+            return
+        if path in ("/v1/spatial/schema", "/v1/tools"):
+            _json_response(
+                self,
+                200,
+                {
+                    "tools": [TOOL_SCHEMA],
+                    "report_schema": REPORT_SCHEMA_VERSION,
+                },
+            )
+            return
+        _json_response(self, 404, {"error": "not found", "path": path})
+
+    def do_POST(self) -> None:
+        path = urlparse(self.path).path
+        try:
+            body = _read_json(self)
+        except Exception as e:
+            _json_response(self, 400, {"error": f"invalid json: {e}"})
+            return
+        try:
+            if path in (
+                "/v1/spatial/analyze",
+                "/v1/spatial/analyze_file",
+                "/v1/spatial/demo_scene",
+                "/v1/spatial/vision",
+                "/v1/spatial/fuse",
+            ):
+                # map path to mode if not set
+                if path.endswith("demo_scene") and "mode" not in body:
+                    body["mode"] = "demo_scene"
+                elif path.endswith("vision") and "mode" not in body:
+                    body["mode"] = "vision"
+                elif path.endswith("fuse") and "mode" not in body:
+                    body["mode"] = "fuse"
+                elif path.endswith("analyze_file") and "mode" not in body:
+                    body["mode"] = "ambix_file" if body.get("ambix") else "mono_file"
+                elif "mode" not in body:
+                    body["mode"] = "demo_scene"
+                result = handle_analyze(body)
+                _json_response(self, 200, result)
+                return
+            _json_response(self, 404, {"error": "not found", "path": path})
+        except Exception as e:
+            _json_response(
+                self,
+                400,
+                {"error": str(e), "trace": traceback.format_exc()[-800:]},
+            )
+
+
+def serve(host: str = DEFAULT_HOST, port: int = DEFAULT_PORT) -> None:
+    httpd = ThreadingHTTPServer((host, port), SpatialHandler)
+    print(f"spatial-hoa API http://{host}:{port}  (hoa64 {__version__})")
+    print("  GET  /health")
+    print("  GET  /v1/spatial/schema")
+    print("  POST /v1/spatial/analyze")
+    httpd.serve_forever()
+
+
+def main(argv: list[str] | None = None) -> int:
+    import argparse
+
+    p = argparse.ArgumentParser(description="spatial-hoa HTTP API for Qwythos/agents")
+    p.add_argument("--host", default=DEFAULT_HOST)
+    p.add_argument("--port", type=int, default=DEFAULT_PORT)
+    args = p.parse_args(argv)
+    serve(args.host, args.port)
+    return 0
+
+
+if __name__ == "__main__":
+    raise SystemExit(main())

hoa64/stft.py

@@ -0,0 +1,96 @@
+"""Minimal STFT / ISTFT in pure NumPy (no scipy)."""
+
+from __future__ import annotations
+
+import numpy as np
+
+
+def hann_window(n: int) -> np.ndarray:
+    if n <= 1:
+        return np.ones(n, dtype=np.float64)
+    return 0.5 - 0.5 * np.cos(2.0 * np.pi * np.arange(n, dtype=np.float64) / n)
+
+
+def stft(
+    x: np.ndarray,
+    *,
+    n_fft: int = 1024,
+    hop: int = 256,
+    window: np.ndarray | None = None,
+) -> tuple[np.ndarray, np.ndarray]:
+    """Short-time FFT.
+
+    Parameters
+    ----------
+    x : (n_samples,) real
+    Returns
+    -------
+    freqs_bins : (n_fft//2+1,)  (normalized later by caller with sr)
+    S : complex64 (n_bins, n_frames)
+    """
+    x = np.asarray(x, dtype=np.float64).reshape(-1)
+    if window is None:
+        window = hann_window(n_fft)
+    window = np.asarray(window, dtype=np.float64)
+    if window.shape[0] != n_fft:
+        raise ValueError("window length must equal n_fft")
+
+    if x.shape[0] < n_fft:
+        x = np.pad(x, (0, n_fft - x.shape[0]))
+
+    n_frames = 1 + (x.shape[0] - n_fft) // hop
+    n_bins = n_fft // 2 + 1
+    S = np.empty((n_bins, n_frames), dtype=np.complex128)
+    for i in range(n_frames):
+        start = i * hop
+        frame = x[start : start + n_fft] * window
+        spec = np.fft.rfft(frame, n=n_fft)
+        S[:, i] = spec
+    return S
+
+
+def stft_freqs(n_fft: int, sample_rate: int) -> np.ndarray:
+    return np.fft.rfftfreq(n_fft, d=1.0 / sample_rate)
+
+
+def frame_signal(
+    x: np.ndarray,
+    *,
+    frame_len: int,
+    hop: int,
+    window: np.ndarray | None = None,
+) -> np.ndarray:
+    """Slice a 1-D signal into overlapping frames (n_frames, frame_len)."""
+    x = np.asarray(x, dtype=np.float64).reshape(-1)
+    if window is None:
+        window = hann_window(frame_len)
+    if x.shape[0] < frame_len:
+        x = np.pad(x, (0, frame_len - x.shape[0]))
+    n_frames = 1 + (x.shape[0] - frame_len) // hop
+    out = np.empty((n_frames, frame_len), dtype=np.float64)
+    for i in range(n_frames):
+        start = i * hop
+        out[i] = x[start : start + frame_len] * window
+    return out
+
+
+def frame_multichannel(
+    audio: np.ndarray,
+    *,
+    frame_len: int,
+    hop: int,
+) -> np.ndarray:
+    """Frame multi-channel audio (C, T) → (n_frames, C, frame_len) rectangular (no window)."""
+    a = np.asarray(audio, dtype=np.float64)
+    if a.ndim != 2:
+        raise ValueError("audio must be (C,T)")
+    C, T = a.shape
+    if T < frame_len:
+        a = np.pad(a, ((0, 0), (0, frame_len - T)))
+        T = a.shape[1]
+    n_frames = 1 + (T - frame_len) // hop
+    out = np.empty((n_frames, C, frame_len), dtype=np.float64)
+    for i in range(n_frames):
+        start = i * hop
+        out[i] = a[:, start : start + frame_len]
+    return out

hoa64/stream.py

@@ -0,0 +1,235 @@
+"""Time-domain HOA streams: encode mono sources, frame-wise analysis."""
+
+from __future__ import annotations
+
+from dataclasses import dataclass
+from typing import Sequence
+
+import numpy as np
+
+from .analysis import (
+    angular_error_deg,
+    doa_from_intensity,
+    field_energy,
+    peak_direction,
+)
+from .audio_io import ensure_hoa_channels
+from .basis import MAX_ORDER, N_CHANNELS, sh_sn3d
+from .encode import encode_plane_waves, mix
+from .stft import frame_multichannel, hann_window, stft, stft_freqs
+
+
+@dataclass
+class SourceSpec:
+    """One plane-wave source in a synthetic or annotated scene."""
+
+    azimuth_deg: float
+    elevation_deg: float
+    signal: np.ndarray  # (n_samples,)
+    label: str = ""
+
+
+def encode_mono_plane_wave(
+    signal: np.ndarray,
+    azimuth_deg: float,
+    elevation_deg: float = 0.0,
+    *,
+    max_order: int = MAX_ORDER,
+) -> np.ndarray:
+    """Mono signal from a known direction → Ambix HOA stream (C, T)."""
+    sig = np.asarray(signal, dtype=np.float64).reshape(1, -1)
+    return encode_plane_waves(
+        [azimuth_deg],
+        [elevation_deg],
+        sig,
+        degrees=True,
+        max_order=max_order,
+    )
+
+
+def encode_scene(
+    sources: Sequence[SourceSpec],
+    *,
+    max_order: int = MAX_ORDER,
+) -> np.ndarray:
+    """Superpose multiple plane-wave sources into one HOA stream (C, T)."""
+    if not sources:
+        raise ValueError("sources must be non-empty")
+    lengths = [int(np.asarray(s.signal).reshape(-1).shape[0]) for s in sources]
+    T = max(lengths)
+    fields = []
+    for s in sources:
+        sig = np.asarray(s.signal, dtype=np.float64).reshape(-1)
+        if sig.shape[0] < T:
+            sig = np.pad(sig, (0, T - sig.shape[0]))
+        elif sig.shape[0] > T:
+            sig = sig[:T]
+        fields.append(
+            encode_mono_plane_wave(
+                sig, s.azimuth_deg, s.elevation_deg, max_order=max_order
+            )
+        )
+    out = fields[0]
+    for f in fields[1:]:
+        out = out + f
+    return out
+
+
+def hoa_rms(hoa: np.ndarray) -> np.ndarray:
+    """Per-channel RMS. hoa (C,T) → (C,)."""
+    a = np.asarray(hoa, dtype=np.float64)
+    return np.sqrt(np.mean(a * a, axis=-1) + 1e-30)
+
+
+@dataclass
+class FrameAnalysis:
+    t_center_sec: float
+    energy: float
+    doa_az_deg: float
+    doa_el_deg: float
+    peak_az_deg: float
+    peak_el_deg: float
+    peak_power: float
+    order1_energy: float
+
+
+def analyze_hoa_frames(
+    hoa: np.ndarray,
+    sample_rate: int,
+    *,
+    frame_ms: float = 40.0,
+    hop_ms: float = 20.0,
+    max_order: int = MAX_ORDER,
+    peak_grid: bool = False,
+) -> list[FrameAnalysis]:
+    """Short-time spatial analysis of an HOA stream.
+
+    Uses rectangular frames; DOA from order-1 intensity on frame-averaged
+    (or energy-weighted) coefficients. Optional dense peak per frame is slower.
+    """
+    a = ensure_hoa_channels(hoa, max_order=max_order)
+    frame_len = max(1, int(round(sample_rate * frame_ms / 1000.0)))
+    hop = max(1, int(round(sample_rate * hop_ms / 1000.0)))
+    frames = frame_multichannel(a, frame_len=frame_len, hop=hop)
+    # Energy-weighted mean coefficient per frame: sum_t a[c,t]*|a_w| style —
+    # use simple mean of coeffs (works for quasi-stationary plane waves).
+    win = hann_window(frame_len)
+    win = win / (np.sum(win) + 1e-30)
+
+    out: list[FrameAnalysis] = []
+    nch_o1 = 4
+    for i in range(frames.shape[0]):
+        block = frames[i]  # (C, L)
+        # AC-safe: do NOT average coeffs (→0 for audio). Use intensity products.
+        w = win  # (L,)
+        # weighted instantaneous intensity ~ W*X etc.
+        W = block[0] * w
+        Yc = block[1] * w
+        Zc = block[2] * w
+        Xc = block[3] * w
+        I = np.array(
+            [
+                float(np.sum(W * Xc)),
+                float(np.sum(W * Yc)),
+                float(np.sum(W * Zc)),
+            ],
+            dtype=np.float64,
+        )
+        nrm = float(np.linalg.norm(I))
+        if nrm < 1e-18:
+            az, el = 0.0, 0.0
+        else:
+            from .basis import az_el_from_unit
+
+            az, el = az_el_from_unit(I / nrm, degrees=True)
+            az, el = float(az), float(el)
+
+        # RMS energy of the frame
+        energy = float(np.mean(np.sum(block * block, axis=0)))
+        o1e = float(np.mean(np.sum(block[:nch_o1] ** 2, axis=0)))
+
+        if peak_grid:
+            # Build a pseudo-static vector: sign-stable energy-weighted mean
+            # via sqrt of mean squares * sign of correlation with W
+            rms = np.sqrt(np.mean(block * block, axis=1) + 1e-30)
+            sign = np.sign(np.mean(block * block[0:1, :], axis=1) + 1e-30)
+            pseudo = rms * sign
+            paz, pel, pv = peak_direction(
+                pseudo, n_azi=48, n_el=24, max_order=min(max_order, 3)
+            )
+        else:
+            paz, pel, pv = float(az), float(el), energy
+        t_center = (i * hop + 0.5 * frame_len) / float(sample_rate)
+        out.append(
+            FrameAnalysis(
+                t_center_sec=t_center,
+                energy=energy,
+                doa_az_deg=float(az),
+                doa_el_deg=float(el),
+                peak_az_deg=float(paz),
+                peak_el_deg=float(pel),
+                peak_power=float(pv),
+                order1_energy=o1e,
+            )
+        )
+    return out
+
+
+def analyze_hoa_stft_bands(
+    hoa: np.ndarray,
+    sample_rate: int,
+    *,
+    n_fft: int = 1024,
+    hop: int = 256,
+    band_edges_hz: Sequence[float] | None = None,
+) -> list[dict]:
+    """Per-frequency-band intensity DOA using order-1 HOA channels only.
+
+    Returns list of {band_hz: [lo,hi], doa_az, doa_el, energy}.
+    """
+    a = ensure_hoa_channels(hoa, max_order=1)
+    if band_edges_hz is None:
+        band_edges_hz = [0, 250, 500, 1000, 2000, 4000, 8000, sample_rate / 2]
+
+    # STFT of W,Y,Z,X
+    specs = []
+    for c in range(4):
+        S = stft(a[c], n_fft=n_fft, hop=hop)
+        specs.append(S)
+    freqs = stft_freqs(n_fft, sample_rate)
+    # Time-average power-weighted intensity per bin then fold into bands
+    W, Y, Z, X = specs
+    # Use complex conjugate product for active intensity-like measure
+    # I_x ~ Re(W * conj(X)), etc., averaged over frames
+    Ix = np.mean(np.real(W * np.conj(X)), axis=1)
+    Iy = np.mean(np.real(W * np.conj(Y)), axis=1)
+    Iz = np.mean(np.real(W * np.conj(Z)), axis=1)
+    Ew = np.mean(np.abs(W) ** 2, axis=1)
+
+    edges = list(band_edges_hz)
+    reports = []
+    for lo, hi in zip(edges[:-1], edges[1:]):
+        mask = (freqs >= lo) & (freqs < hi)
+        if not np.any(mask):
+            continue
+        I = np.array(
+            [np.sum(Ix[mask]), np.sum(Iy[mask]), np.sum(Iz[mask])],
+            dtype=np.float64,
+        )
+        n = np.linalg.norm(I)
+        if n < 1e-15:
+            az, el = 0.0, 0.0
+        else:
+            from .basis import az_el_from_unit
+
+            az, el = az_el_from_unit(I / n, degrees=True)
+            az, el = float(az), float(el)
+        reports.append(
+            {
+                "band_hz": [float(lo), float(hi)],
+                "doa_az_deg": az,
+                "doa_el_deg": el,
+                "energy": float(np.sum(Ew[mask])),
+            }
+        )
+    return reports

hoa64/synth.py

@@ -0,0 +1,46 @@
+"""Synthetic test signals for Phase 1 audio pipeline."""
+
+from __future__ import annotations
+
+import numpy as np
+
+
+def tone(
+    freq_hz: float,
+    duration_sec: float,
+    sample_rate: int,
+    *,
+    amplitude: float = 0.5,
+    phase: float = 0.0,
+) -> np.ndarray:
+    t = np.arange(int(round(duration_sec * sample_rate)), dtype=np.float64) / sample_rate
+    return amplitude * np.sin(2.0 * np.pi * freq_hz * t + phase)
+
+
+def noise(
+    duration_sec: float,
+    sample_rate: int,
+    *,
+    amplitude: float = 0.2,
+    seed: int = 0,
+) -> np.ndarray:
+    rng = np.random.default_rng(seed)
+    n = int(round(duration_sec * sample_rate))
+    return amplitude * rng.standard_normal(n)
+
+
+def envelope_adsr(
+    n: int,
+    sample_rate: int,
+    *,
+    attack: float = 0.01,
+    release: float = 0.05,
+) -> np.ndarray:
+    env = np.ones(n, dtype=np.float64)
+    a = min(n, int(attack * sample_rate))
+    r = min(n, int(release * sample_rate))
+    if a > 0:
+        env[:a] = np.linspace(0.0, 1.0, a)
+    if r > 0:
+        env[-r:] = np.linspace(1.0, 0.0, r)
+    return env

hoa64/vision.py

@@ -0,0 +1,237 @@
+"""Phase 3 — vision tower: project boxes/rays onto the HOA-7 sphere.
+
+No neural net. Detections become soft spherical Gaussians in Ambix SN3D
+space (same 64-D geometry as audio). Suitable as conditioning / tool output
+for language models and later diffusion control.
+"""
+
+from __future__ import annotations
+
+import math
+from typing import Any, Iterable, Mapping, Optional, Sequence
+
+import numpy as np
+
+from .basis import MAX_ORDER, N_CHANNELS, sh_sn3d, unit_vector
+from .encode import mix
+from .report import SpatialReport, report_from_hoa
+
+
+def _as_box(b: Mapping[str, Any]) -> dict:
+    """Normalize a box/ray dict."""
+    if "az" in b or "azimuth" in b or "azimuth_deg" in b:
+        az = float(b.get("az", b.get("azimuth", b.get("azimuth_deg", 0.0))))
+        el = float(b.get("el", b.get("elevation", b.get("elevation_deg", 0.0))))
+    elif "x" in b and "y" in b:
+        # normalized image coords: x in [0,1] left→right, y in [0,1] top→bottom
+        # map to az in [-180,180] around center, el in [-90,90]
+        x = float(b["x"])
+        y = float(b["y"])
+        # simple equirectangular pinhole-ish: center (0.5,0.5) = front
+        # az: left of image = +az (left), right = -az if we want Ambix left=+Y
+        # Convention: image x=0 left → az=+hfov/2, x=1 right → az=-hfov/2
+        hfov = float(b.get("hfov_deg", 90.0))
+        vfov = float(b.get("vfov_deg", 60.0))
+        az = (0.5 - x) * hfov
+        el = (0.5 - y) * vfov
+    else:
+        raise ValueError(f"box needs az/el or x/y: {b}")
+
+    return {
+        "az": az,
+        "el": el,
+        "w_deg": float(b.get("w_deg", b.get("width_deg", b.get("sigma_deg", 8.0)))),
+        "h_deg": float(b.get("h_deg", b.get("height_deg", b.get("sigma_deg", 8.0)))),
+        "weight": float(b.get("weight", b.get("score", b.get("confidence", 1.0)))),
+        "label": str(b.get("label", b.get("class", b.get("name", "")))),
+        "kind": str(b.get("kind", "box")),  # box | ray | point
+    }
+
+
+def angular_gaussian_weights(
+    center_az: float,
+    center_el: float,
+    sample_az: np.ndarray,
+    sample_el: np.ndarray,
+    sigma_az_deg: float,
+    sigma_el_deg: float,
+) -> np.ndarray:
+    """Soft lobe on the sphere (product of wrapped azimuth + elevation Gaussians)."""
+    # great-circle-ish separation using unit vectors is better
+    c = unit_vector(center_az, center_el, degrees=True)
+    # sample_az, sample_el may be mesh
+    u = unit_vector(sample_az, sample_el, degrees=True)
+    # cos gamma = c · u
+    cosg = np.clip(np.sum(u * c, axis=-1), -1.0, 1.0)
+    gamma = np.rad2deg(np.arccos(cosg))
+    # isotropic sigma from geometric mean of width axes
+    sig = max(0.5, math.sqrt(max(sigma_az_deg, 0.5) * max(sigma_el_deg, 0.5)))
+    return np.exp(-0.5 * (gamma / sig) ** 2)
+
+
+def encode_boxes_to_hoa(
+    boxes: Sequence[Mapping[str, Any]],
+    *,
+    max_order: int = MAX_ORDER,
+    n_azi: int = 72,
+    n_el: int = 36,
+) -> np.ndarray:
+    """Project one or more vision boxes/rays into a static HOA-7 field (C,).
+
+    Each box is a soft spherical Gaussian; the field is the weighted sum of
+    SN3D plane-wave encodings of the discrete sphere samples of that lobe.
+    """
+    if not boxes:
+        return np.zeros((max_order + 1) ** 2, dtype=np.float64)
+
+    azi = np.linspace(-180.0, 180.0, n_azi, endpoint=False)
+    el = np.linspace(-90.0, 90.0, n_el)
+    AA, EE = np.meshgrid(azi, el, indexing="ij")
+    # solid-angle-ish weights
+    w_el = np.cos(np.deg2rad(el))
+    w_el = np.maximum(w_el, 0.0)
+    dA = (2.0 * math.pi / n_azi) * (math.pi / max(n_el - 1, 1))
+    sa = dA * w_el[None, :]  # (1,E) broadcast to (A,E)
+
+    field = np.zeros((max_order + 1) ** 2, dtype=np.float64)
+    for raw in boxes:
+        b = _as_box(raw)
+        if b["kind"] == "ray" or b["kind"] == "point":
+            # delta: single direction plane wave
+            Y = sh_sn3d(b["az"], b["el"], degrees=True, max_order=max_order)
+            field = field + b["weight"] * Y
+            continue
+        lobe = angular_gaussian_weights(
+            b["az"], b["el"], AA, EE, b["w_deg"], b["h_deg"]
+        )
+        # discrete HOA analysis of the lobe function f(Ω)
+        # a = sum f(Ω) Y(Ω) sa(Ω) * (2n+1)/4π   — use same scale as rotate dense
+        Y = sh_sn3d(AA, EE, degrees=True, max_order=max_order)  # (A,E,C)
+        f = lobe * b["weight"]
+        # a_c = sum_{ae} f_ae * sa_ae * Y_ae,c
+        weighted = (f * sa)[..., None] * Y
+        a = np.sum(weighted, axis=(0, 1))
+        # SN3D analysis scale per order
+        nch = a.shape[0]
+        scale = np.empty(nch, dtype=np.float64)
+        for n in range(max_order + 1):
+            s = float(2 * n + 1)
+            for m in range(-n, n + 1):
+                scale[n * (n + 1) + m] = s
+        a = a * scale / (4.0 * math.pi)
+        field = field + a
+
+    # pad to 64
+    out = np.zeros(N_CHANNELS, dtype=np.float64)
+    out[: field.shape[0]] = field
+    return out
+
+
+def report_from_boxes(
+    boxes: Sequence[Mapping[str, Any]],
+    *,
+    max_order: int = MAX_ORDER,
+    sample_rate: int = 48000,
+    n_samples: int = 1,
+    meta: Optional[dict] = None,
+) -> SpatialReport:
+    """Vision-only spatial report from detection boxes/rays."""
+    hoa = encode_boxes_to_hoa(boxes, max_order=max_order)
+    # report_from_hoa expects (C,T) — use a constant "frame"
+    stream = np.tile(hoa.reshape(-1, 1), (1, max(1, n_samples)))
+    notes = [f"vision tower: {len(boxes)} box/ray event(s) on HOA sphere"]
+    hints = []
+    for i, raw in enumerate(boxes):
+        try:
+            b = _as_box(raw)
+            hints.append(
+                {
+                    "label": b["label"] or f"box{i}",
+                    "azimuth_deg": b["az"],
+                    "elevation_deg": b["el"],
+                    "weight": b["weight"],
+                    "kind": b["kind"],
+                }
+            )
+        except Exception:
+            continue
+    m = dict(meta or {})
+    m["encode"] = "vision_boxes"
+    m["n_boxes"] = len(boxes)
+    rep = report_from_hoa(
+        stream,
+        sample_rate,
+        max_order=max_order,
+        include_frames=False,
+        include_bands=False,
+        include_peak_map=True,
+        notes=notes,
+        meta=m,
+    )
+    rep.kind = "spatial_vision"
+    rep.sources_hint = hints
+    return rep
+
+
+def fuse_reports(audio: Mapping[str, Any], vision: Mapping[str, Any]) -> dict:
+    """Merge audio + vision spatial reports into one agent-facing payload.
+
+    Strategy:
+      - Keep both DOAs
+      - angular_separation_deg between them
+      - agreement flag if within 15°
+      - combined one_liner
+    """
+    from .analysis import angular_error_deg
+
+    a_az = float(audio.get("doa_az_deg", 0.0))
+    a_el = float(audio.get("doa_el_deg", 0.0))
+    v_az = float(vision.get("doa_az_deg", vision.get("peak_az_deg", 0.0)))
+    v_el = float(vision.get("doa_el_deg", vision.get("peak_el_deg", 0.0)))
+    sep = angular_error_deg(a_az, a_el, v_az, v_el)
+    agree = sep <= 15.0
+
+    # energy-weighted blend of unit vectors if both have energy
+    from .basis import az_el_from_unit, unit_vector
+
+    ea = float(audio.get("energy", 0.0)) + 1e-12
+    ev = float(vision.get("energy", 0.0)) + 1e-12
+    ua = unit_vector(a_az, a_el, degrees=True)
+    uv = unit_vector(v_az, v_el, degrees=True)
+    um = (ea * ua + ev * uv) / (ea + ev)
+    nm = float(np.linalg.norm(um))
+    if nm > 1e-15:
+        baz, bel = az_el_from_unit(um / nm, degrees=True)
+        baz, bel = float(baz), float(bel)
+    else:
+        baz, bel = a_az, a_el
+
+    one = (
+        f"spatial-fuse: audio=({a_az:.1f},{a_el:.1f}) vision=({v_az:.1f},{v_el:.1f}) "
+        f"sep={sep:.1f}° agree={agree} blend=({baz:.1f},{bel:.1f})"
+    )
+    return {
+        "schema": "spatial-hoa.fuse.v1",
+        "kind": "spatial_av_fuse",
+        "audio_doa_az_deg": a_az,
+        "audio_doa_el_deg": a_el,
+        "vision_doa_az_deg": v_az,
+        "vision_doa_el_deg": v_el,
+        "angular_separation_deg": sep,
+        "agreement": agree,
+        "blend_az_deg": baz,
+        "blend_el_deg": bel,
+        "audio_energy": float(audio.get("energy", 0.0)),
+        "vision_energy": float(vision.get("energy", 0.0)),
+        "audio_one_liner": audio.get("one_liner")
+        or f"audio DOA ({a_az:.1f},{a_el:.1f})",
+        "vision_one_liner": vision.get("one_liner")
+        or f"vision DOA ({v_az:.1f},{v_el:.1f})",
+        "vision_sources": vision.get("sources_hint") or [],
+        "audio_bands": audio.get("bands") or [],
+        "one_liner": one,
+        "notes": [
+            "Fused audio intensity DOA with vision spherical presence peak.",
+            "agreement true if separation ≤ 15°.",
+        ],
+    }

hoa64/wigner.py

@@ -0,0 +1,266 @@
+"""Fast HOA rotation via Wigner-D matrices (real Ambix ACN / SN3D).
+
+Within each order n the (2n+1) coefficients transform by a real matrix built
+from complex Wigner D functions. Normalization (SN3D vs N3D) cancels inside
+an order, so the same matrices apply to Ambix SN3D.
+
+Public:
+  rotation_matrix_zyx     — 3×3 active rotation (yaw/pitch/roll)
+  hoa_rotation_matrix     — full ((N+1)²)×((N+1)²) block-diagonal matrix
+  apply_hoa_rotation      — a' = M @ a  (also (C,T) streams)
+"""
+
+from __future__ import annotations
+
+import math
+from functools import lru_cache
+from typing import Tuple
+
+import numpy as np
+
+from .basis import MAX_ORDER, N_CHANNELS, acn_index
+
+
+def rotation_matrix_zyx(
+    yaw: float, pitch: float, roll: float, *, degrees: bool = True
+) -> np.ndarray:
+    """R = Rz(yaw) @ Ry(pitch) @ Rx(roll); active rotation of column vectors."""
+    if degrees:
+        yaw, pitch, roll = map(math.radians, (yaw, pitch, roll))
+    cy, sy = math.cos(yaw), math.sin(yaw)
+    cp, sp = math.cos(pitch), math.sin(pitch)
+    cr, sr = math.cos(roll), math.sin(roll)
+    Rz = np.array([[cy, -sy, 0.0], [sy, cy, 0.0], [0.0, 0.0, 1.0]])
+    Ry = np.array([[cp, 0.0, sp], [0.0, 1.0, 0.0], [-sp, 0.0, cp]])
+    Rx = np.array([[1.0, 0.0, 0.0], [0.0, cr, -sr], [0.0, sr, cr]])
+    return Rz @ Ry @ Rx
+
+
+def _fact(n: int) -> float:
+    if n < 0:
+        return 0.0
+    return float(math.factorial(n))
+
+
+@lru_cache(maxsize=8192)
+def wigner_d(j: int, mp: int, m: int, beta: float) -> float:
+    """Wigner small-d matrix element d^j_{mp,m}(beta), real.
+
+    Explicit finite sum (safe for j ≤ 7 used here).
+    """
+    if abs(mp) > j or abs(m) > j:
+        return 0.0
+    # Numerical stability: beta in [0, pi] preferred but general ok
+    cb = math.cos(beta * 0.5)
+    sb = math.sin(beta * 0.5)
+    # Avoid 0**negative
+    s_min = max(0, m - mp)
+    s_max = min(j + m, j - mp)
+    total = 0.0
+    for s in range(s_min, s_max + 1):
+        den = (
+            _fact(j + m - s)
+            * _fact(j - mp - s)
+            * _fact(s)
+            * _fact(s + mp - m)
+        )
+        if den == 0.0:
+            continue
+        num = _fact(j + m) * _fact(j - m) * _fact(j + mp) * _fact(j - mp)
+        pref = ((-1.0) ** (mp - m + s)) * math.sqrt(num) / den
+        cpow = 2 * j + m - mp - 2 * s
+        spow = mp - m + 2 * s
+        # handle base cases
+        cterm = 1.0 if cpow == 0 else (cb ** cpow if abs(cb) > 1e-15 or cpow > 0 else 0.0)
+        sterm = 1.0 if spow == 0 else (sb ** spow if abs(sb) > 1e-15 or spow > 0 else 0.0)
+        if cpow < 0 and abs(cb) < 1e-15:
+            cterm = 0.0
+        if spow < 0 and abs(sb) < 1e-15:
+            sterm = 0.0
+        total += pref * cterm * sterm
+    return total
+
+
+def wigner_D_complex(
+    j: int, alpha: float, beta: float, gamma: float
+) -> np.ndarray:
+    """Complex Wigner-D matrix for order j, indices m',m ∈ [-j..j].
+
+    Ordering: row/col index = m + j  (m from -j to +j).
+    D_{mp,m} = e^{-i mp α} d_{mp,m}(β) e^{-i m γ}
+    """
+    dim = 2 * j + 1
+    D = np.zeros((dim, dim), dtype=np.complex128)
+    for imp, mp in enumerate(range(-j, j + 1)):
+        for im, m in enumerate(range(-j, j + 1)):
+            d = wigner_d(j, mp, m, beta)
+            D[imp, im] = (
+                math.cos(mp * alpha)
+                - 1j * math.sin(mp * alpha)
+            ) * d * (
+                math.cos(m * gamma) - 1j * math.sin(m * gamma)
+            )
+            # e^{-i θ} = cosθ - i sinθ
+    return D
+
+
+def _real_to_complex_matrix(j: int) -> np.ndarray:
+    """Unitary map U: real ACN coeffs (m=-j..j) → complex m=-j..j.
+
+    Convention (common in Ambisonics / real SH):
+      c_0 = r_0
+      c_{+m} = (-1)^m / √2 * (r_{+m} - i r_{-m})
+      c_{-m} =         1 / √2 * (r_{+m} + i r_{-m})
+    so r = U^H c  and c = U r  with U unitary.
+    """
+    dim = 2 * j + 1
+    U = np.zeros((dim, dim), dtype=np.complex128)
+    # index helper: m -> i = m + j
+    def idx(m: int) -> int:
+        return m + j
+
+    U[idx(0), idx(0)] = 1.0 + 0.0j
+    s2 = 1.0 / math.sqrt(2.0)
+    for m in range(1, j + 1):
+        sign = (-1.0) ** m
+        # c_{+m} from r_{+m}, r_{-m}
+        U[idx(m), idx(m)] = sign * s2
+        U[idx(m), idx(-m)] = -1j * sign * s2
+        # c_{-m} from r_{+m}, r_{-m}
+        U[idx(-m), idx(m)] = s2
+        U[idx(-m), idx(-m)] = 1j * s2
+    return U
+
+
+def real_sh_rotation_block(
+    j: int, alpha: float, beta: float, gamma: float
+) -> np.ndarray:
+    """Real (2j+1)×(2j+1) rotation matrix for ACN order-j block (m=-j..j)."""
+    if j == 0:
+        return np.array([[1.0]], dtype=np.float64)
+    U = _real_to_complex_matrix(j)
+    D = wigner_D_complex(j, alpha, beta, gamma)
+    # real coeffs: r' = U^H D U r
+    R_c = U.conj().T @ D @ U
+    # Should be real symmetric orthogonal (numerically tiny imag)
+    return np.real(R_c)
+
+
+def rotation_matrix_to_zyz(R: np.ndarray) -> Tuple[float, float, float]:
+    """Extract ZYZ Euler angles (α, β, γ) from a right-handed rotation matrix.
+
+    R = Rz(α) @ Ry(β) @ Rz(γ)  (active, column vectors).
+    """
+    R = np.asarray(R, dtype=np.float64)
+    # β ∈ [0, π]
+    # R[2,2] = cos β
+    cbeta = float(np.clip(R[2, 2], -1.0, 1.0))
+    beta = math.acos(cbeta)
+    sb = math.sin(beta)
+    if abs(sb) > 1e-10:
+        # α = atan2(R[1,2]/sinβ, R[0,2]/sinβ)
+        alpha = math.atan2(R[1, 2] / sb, R[0, 2] / sb)
+        # γ = atan2(R[2,1]/sinβ, -R[2,0]/sinβ)
+        gamma = math.atan2(R[2, 1] / sb, -R[2, 0] / sb)
+    else:
+        # Gimbal: β ≈ 0 or π — only α+γ or α-γ determined
+        alpha = math.atan2(-R[0, 1], R[0, 0])
+        gamma = 0.0
+        if cbeta < 0:
+            # β = π
+            alpha = math.atan2(R[0, 1], -R[0, 0])
+    return alpha, beta, gamma
+
+
+def hoa_rotation_matrix(
+    R3: np.ndarray,
+    *,
+    max_order: int = MAX_ORDER,
+) -> np.ndarray:
+    """Full Ambix ACN rotation matrix for orders 0..max_order.
+
+    Applies the *same* geometric rotation as R3 does to Cartesian vectors.
+    """
+    max_order = int(max_order)
+    nch = (max_order + 1) ** 2
+    M = np.zeros((nch, nch), dtype=np.float64)
+    alpha, beta, gamma = rotation_matrix_to_zyz(R3)
+
+    # Order 0
+    M[0, 0] = 1.0
+
+    # Order 1: direct Cartesian for numerical exactness / convention lock.
+    # ACN: [Y, Z, X] = indices 1,2,3 ; cart = [X,Y,Z]
+    # v' = R3 @ v  ⇒  [X',Y',Z'] = R3 @ [X,Y,Z]
+    if max_order >= 1:
+        # Build 3×3 block mapping [Y,Z,X] -> [Y',Z',X']
+        # [X']   [R00 R01 R02] [X]
+        # [Y'] = [R10 R11 R12] [Y]
+        # [Z']   [R20 R21 R22] [Z]
+        # Y' = R10 X + R11 Y + R12 Z
+        # Z' = R20 X + R21 Y + R22 Z
+        # X' = R00 X + R01 Y + R02 Z
+        # Coefficients of (Y,Z,X):
+        # Y' = R11 Y + R12 Z + R10 X
+        # Z' = R21 Y + R22 Z + R20 X
+        # X' = R01 Y + R02 Z + R00 X
+        B = np.array(
+            [
+                [R3[1, 1], R3[1, 2], R3[1, 0]],
+                [R3[2, 1], R3[2, 2], R3[2, 0]],
+                [R3[0, 1], R3[0, 2], R3[0, 0]],
+            ],
+            dtype=np.float64,
+        )
+        M[1:4, 1:4] = B
+
+    for n in range(2, max_order + 1):
+        block = real_sh_rotation_block(n, alpha, beta, gamma)
+        # Verify det ~ 1; if complex conversion convention is flipped we may
+        # need transpose — tests catch this against plane-wave re-encode.
+        i0 = n * n  # ACN start for order n is n^2 (m=-n → n*(n+1)+(-n)=n^2)
+        # Wait: n*(n+1)+(-n) = n^2 + n - n = n^2. Yes.
+        dim = 2 * n + 1
+        M[i0 : i0 + dim, i0 : i0 + dim] = block
+
+    return M
+
+
+def apply_hoa_rotation(
+    hoa: np.ndarray,
+    M: np.ndarray,
+) -> np.ndarray:
+    """Apply precomputed rotation matrix to (C,) or (C,T) coefficients."""
+    a = np.asarray(hoa, dtype=np.float64)
+    nch = M.shape[0]
+    if a.ndim == 1:
+        aa = np.zeros(nch, dtype=np.float64)
+        n = min(nch, a.shape[0])
+        aa[:n] = a[:n]
+        out = M @ aa
+        if a.shape[0] > nch:
+            full = np.zeros_like(a)
+            full[:nch] = out
+            return full
+        if a.shape[0] < N_CHANNELS:
+            full = np.zeros(N_CHANNELS, dtype=np.float64)
+            full[:nch] = out
+            return full
+        return out
+    if a.ndim == 2:
+        aa = np.zeros((nch, a.shape[1]), dtype=np.float64)
+        n = min(nch, a.shape[0])
+        aa[:n] = a[:n]
+        out = M @ aa
+        if a.shape[0] >= N_CHANNELS:
+            full = np.zeros((max(a.shape[0], N_CHANNELS), a.shape[1]), dtype=np.float64)
+            full[:nch] = out
+            return full[: a.shape[0]]
+        full = np.zeros((N_CHANNELS, a.shape[1]), dtype=np.float64)
+        full[:nch] = out
+        return full
+    raise ValueError("hoa must be (C,) or (C,T)")
+
+
+def clear_wigner_cache() -> None:
+    wigner_d.cache_clear()

integrations/qwythos_system_snippet.md

@@ -0,0 +1,50 @@
+# Qwythos / Pi bridge — spatial tool wiring
+
+The HOA-7 spatial calculator is **not** part of the LLM weights. Qwythos should
+call it as a tool and treat the JSON / one-liner as ground-truth geometry.
+
+## Runtime
+
+| Piece | Address |
+|-------|---------|
+| Qwythos (llama-server) | `http://127.0.0.1:8000/v1` |
+| Olympus Pi bridge | `http://127.0.0.1:8642` |
+| **spatial-hoa API** | `http://127.0.0.1:8765` |
+| CLI on PATH | `spatial-report` |
+
+Enable API at login:
+
+```bash
+systemctl --user enable --now spatial-hoa.service
+curl -s http://127.0.0.1:8765/health
+```
+
+## For Pi bridge (bash tool)
+
+Pi already has `bash`. Instruct the agent:
+
+```text
+When the user asks about spatial audio/location/direction of sound or objects,
+run one of:
+  spatial-report demo-scene -o /tmp/spatial.json && cat /tmp/spatial.json
+  spatial-report analyze /path/to.wav --ambix -o /tmp/spatial.json && cat /tmp/spatial.json
+  spatial-report analyze /path/to.wav --az DEG --el DEG -o /tmp/spatial.json
+  spatial-report vision --boxes '[{"az":0,"el":0,"label":"obj"}]' -o /tmp/v.json
+Or HTTP:
+  curl -s -X POST http://127.0.0.1:8765/v1/spatial/analyze \
+    -H 'Content-Type: application/json' \
+    -d '{"mode":"demo_scene","order":3}'
+Use the returned one_liner and doa_* fields as factual spatial state.
+```
+
+Optional: append that block to the Pi system prompt / AGENTS notes.
+
+## OpenAI tool schema
+
+See `tools/spatial_analyze.openai.json` for a function-calling descriptor
+compatible with agents that support tools against a custom executor.
+
+## Coordinates (always)
+
+- Azimuth 0° = front, +90° = left, −90° = right  
+- Elevation 0° = horizon, +90° = zenith  

scripts/spatial-report

@@ -0,0 +1,11 @@
+#!/usr/bin/env bash
+# spatial-report — PATH entry for Qwythos / Pi bash tool / humans.
+# Wraps hoa64 CLI. Examples:
+#   spatial-report demo-scene -o /tmp/r.json
+#   spatial-report analyze /path/file.wav --ambix
+#   spatial-report vision --boxes '[{"az":30,"el":0,"label":"obj"}]'
+#   spatial-report serve
+set -euo pipefail
+ROOT="$(cd "$(dirname "$0")/.." && pwd)"
+export PYTHONPATH="${ROOT}${PYTHONPATH:+:$PYTHONPATH}"
+exec python3 -m hoa64 "$@"

tests/test_basis.py

@@ -0,0 +1,107 @@
+"""Basis / Farina formula sanity checks."""
+
+from __future__ import annotations
+
+import math
+import sys
+from pathlib import Path
+
+import numpy as np
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.basis import (
+    MAX_ORDER,
+    N_CHANNELS,
+    acn_index,
+    acn_nm,
+    channel_names,
+    sh_sn3d,
+    unit_vector,
+)
+
+
+def test_layout():
+    assert N_CHANNELS == 64
+    assert MAX_ORDER == 7
+    assert acn_index(0, 0) == 0
+    assert acn_index(1, -1) == 1
+    assert acn_index(1, 0) == 2
+    assert acn_index(1, 1) == 3
+    assert acn_index(7, 7) == 63
+    for i in range(64):
+        n, m = acn_nm(i)
+        assert acn_index(n, m) == i
+    names = channel_names()
+    assert names[0] == "W"
+    assert names[1] == "Y"
+    assert names[2] == "Z"
+    assert names[3] == "X"
+    assert len(names) == 64
+
+
+def test_unit_vector_cardinals():
+    # front, left, up
+    f = unit_vector(0.0, 0.0)
+    l = unit_vector(90.0, 0.0)
+    u = unit_vector(0.0, 90.0)
+    np.testing.assert_allclose(f, [1, 0, 0], atol=1e-12)
+    np.testing.assert_allclose(l, [0, 1, 0], atol=1e-12)
+    np.testing.assert_allclose(u, [0, 0, 1], atol=1e-12)
+
+
+def test_order1_cartesian_identity():
+    """Order-1 SN3D: Y=y, Z=z, X=x at unit directions."""
+    for az, el, xyz in [
+        (0, 0, (1, 0, 0)),
+        (90, 0, (0, 1, 0)),
+        (0, 90, (0, 0, 1)),
+        (-90, 0, (0, -1, 0)),
+    ]:
+        y = sh_sn3d(az, el)
+        assert y.shape == (64,)
+        assert abs(y[0] - 1.0) < 1e-12  # W
+        x, yy, z = xyz
+        assert abs(y[1] - yy) < 1e-12
+        assert abs(y[2] - z) < 1e-12
+        assert abs(y[3] - x) < 1e-12
+
+
+def test_front_source_signs():
+    y = sh_sn3d(0.0, 0.0)
+    # front: X dominant positive among order-1 dipoles
+    assert y[3] > 0.9
+    assert abs(y[1]) < 1e-12
+    assert abs(y[2]) < 1e-12
+
+
+def test_finite_all_channels():
+    y = sh_sn3d(33.0, -12.0)
+    assert y.shape == (64,)
+    assert np.all(np.isfinite(y))
+    # batch
+    yb = sh_sn3d([0, 90, 180], [0, 0, 0])
+    assert yb.shape == (3, 64)
+
+
+def test_zenith_only_z_orders():
+    y = sh_sn3d(0.0, 90.0)
+    # at north pole, only m=0 channels should be nonzero (real SH)
+    for n in range(0, 8):
+        for m in range(-n, n + 1):
+            idx = acn_index(n, m)
+            if m != 0:
+                assert abs(y[idx]) < 1e-8, f"n={n} m={m} val={y[idx]}"
+
+
+if __name__ == "__main__":
+    for fn in [
+        test_layout,
+        test_unit_vector_cardinals,
+        test_order1_cartesian_identity,
+        test_front_source_signs,
+        test_finite_all_channels,
+        test_zenith_only_z_orders,
+    ]:
+        fn()
+        print("OK", fn.__name__)

tests/test_encode_decode.py

@@ -0,0 +1,95 @@
+"""Encode / decode / DOA hypothesis tests."""
+
+from __future__ import annotations
+
+import sys
+from pathlib import Path
+
+import numpy as np
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.analysis import (
+    angular_error_deg,
+    doa_from_intensity,
+    field_energy,
+    peak_direction,
+)
+from hoa64.decode import beamform, decode_directions
+from hoa64.encode import encode_plane_waves, encode_points, mix
+from hoa64.basis import sh_sn3d
+
+
+def test_encode_matches_basis():
+    a = encode_points([45.0], [10.0], [1.0])
+    y = sh_sn3d(45.0, 10.0)
+    np.testing.assert_allclose(a, y, atol=1e-12)
+
+
+def test_beamform_peaks_at_source():
+    a = encode_points([30.0], [-15.0], [1.0])
+    on = float(beamform(a, 30.0, -15.0))
+    off = float(beamform(a, -150.0, 40.0))
+    assert on > off
+    assert on > 0.5  # SN3D self-inner-product roughly order-dependent
+
+
+def test_intensity_doa_cardinal():
+    for az, el in [(0, 0), (90, 0), (-90, 0), (0, 45)]:
+        a = encode_points([az], [el], [1.0])
+        az_h, el_h = doa_from_intensity(a)
+        err = angular_error_deg(az, el, az_h, el_h)
+        assert err < 2.0, f"src=({az},{el}) hat=({az_h},{el_h}) err={err}"
+
+
+def test_peak_direction_near_source():
+    a = encode_points([120.0], [25.0], [1.0])
+    paz, pel, pval = peak_direction(a, n_azi=120, n_el=60)
+    err = angular_error_deg(120.0, 25.0, paz, pel)
+    assert err < 5.0, f"peak=({paz},{pel}) err={err}"
+    assert pval > 0
+
+
+def test_superposition_linearity():
+    a1 = encode_points([0], [0], [1.0])
+    a2 = encode_points([90], [0], [0.5])
+    m = mix(a1, a2)
+    np.testing.assert_allclose(m, a1 + a2)
+    assert field_energy(m) > field_energy(a1)
+
+
+def test_plane_wave_time_series():
+    t = np.linspace(0, 1, 100, endpoint=False)
+    sig = np.sin(2 * np.pi * 5 * t)[None, :]  # (1, T)
+    a = encode_plane_waves([0.0], [0.0], sig)
+    assert a.shape == (64, 100)
+    # W channel tracks the signal
+    np.testing.assert_allclose(a[0], sig[0], atol=1e-12)
+    # X (front) also tracks for front source
+    np.testing.assert_allclose(a[3], sig[0], atol=1e-12)
+
+
+def test_two_source_separation_energy():
+    a = mix(
+        encode_points([0], [0], [1.0]),
+        encode_points([180], [0], [1.0]),
+    )
+    front = float(beamform(a, 0, 0))
+    back = float(beamform(a, 180, 0))
+    side = float(beamform(a, 90, 0))
+    assert front > side
+    assert back > side
+
+
+if __name__ == "__main__":
+    for fn in [
+        test_encode_matches_basis,
+        test_beamform_peaks_at_source,
+        test_intensity_doa_cardinal,
+        test_peak_direction_near_source,
+        test_superposition_linearity,
+        test_plane_wave_time_series,
+        test_two_source_separation_energy,
+    ]:
+        fn()
+        print("OK", fn.__name__)

tests/test_integration_extras.py

@@ -0,0 +1,108 @@
+"""Detector, live pipeline (file), conditioning tests."""
+
+from __future__ import annotations
+
+import json
+import sys
+import tempfile
+from pathlib import Path
+
+import numpy as np
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.analysis import angular_error_deg
+from hoa64.audio_io import write_wav
+from hoa64.conditioning import build_conditioning, comfy_txt2img_payload, spatial_prompt_fragment
+from hoa64.detector import (
+    Detection,
+    detections_to_sphere_boxes,
+    load_yolo_labels,
+    write_demo_image_with_box,
+)
+from hoa64.live_audio import audio_to_hoa_stream, live_report_from_file, stereo_to_order1_hoa
+from hoa64.synth import tone
+from hoa64.vision import fuse_reports, report_from_boxes
+
+
+def test_demo_image_detect_sphere():
+    with tempfile.TemporaryDirectory() as td:
+        img = Path(td) / "demo.png"
+        det = write_demo_image_with_box(img)
+        boxes = detections_to_sphere_boxes([det])
+        assert boxes[0]["label"] == "demo_object"
+        # center of box ~ (0.45, 0.5) → slight left az positive? cx=0.45 → az=(0.5-0.45)*90=+4.5
+        assert abs(boxes[0]["az"] - 4.5) < 1.0
+        rep = report_from_boxes(boxes, max_order=3)
+        assert rep.energy > 0
+
+
+def test_yolo_labels():
+    with tempfile.TemporaryDirectory() as td:
+        p = Path(td) / "labels.txt"
+        # class 0 center front-ish
+        p.write_text("0 0.5 0.5 0.2 0.2 0.9\n")
+        dets = load_yolo_labels(p)
+        assert len(dets) == 1
+        boxes = detections_to_sphere_boxes(dets)
+        assert abs(boxes[0]["az"]) < 1.0  # center → az≈0
+
+
+def test_live_pipeline_from_synth_wav():
+    with tempfile.TemporaryDirectory() as td:
+        wav = Path(td) / "m.wav"
+        sr = 16000
+        sig = tone(800, 0.3, sr, amplitude=0.4)
+        write_wav(wav, sig, sr)
+        rep = live_report_from_file(wav, az_deg=30.0, el_deg=0.0, max_order=1)
+        err = angular_error_deg(30.0, 0.0, rep.doa_az_deg, rep.doa_el_deg)
+        assert err < 8.0, f"err={err}"
+
+
+def test_stereo_pseudo_hoa():
+    sr = 8000
+    L = tone(400, 0.2, sr, amplitude=0.5)
+    R = tone(400, 0.2, sr, amplitude=0.1)
+    hoa = stereo_to_order1_hoa(L, R, width_az_deg=40)
+    assert hoa.shape[0] == 64
+    # left-dominant → positive Y energy
+    assert np.mean(hoa[1] ** 2) > 0
+
+
+def test_conditioning_from_fuse():
+    audio = {"doa_az_deg": 15.0, "doa_el_deg": 0.0, "energy": 1.0, "kind": "spatial_field"}
+    vision = {
+        "doa_az_deg": 18.0,
+        "doa_el_deg": 2.0,
+        "energy": 0.5,
+        "kind": "spatial_vision",
+        "peak_az_deg": 18.0,
+        "peak_el_deg": 2.0,
+    }
+    fused = fuse_reports(audio, vision)
+    cond = build_conditioning(fused, base_prompt="cinematic room", style="natural")
+    assert "cinematic room" in cond["positive_prompt"]
+    assert cond["schema"] == "spatial-hoa.conditioning.v1"
+    assert "control_vector" in cond
+    frag = spatial_prompt_fragment(fused, style="tags")
+    assert "spatial-az" in frag
+    # offline payload with explicit checkpoint (no Comfy required)
+    wf = comfy_txt2img_payload(
+        cond, checkpoint="sd_xl_base_1.0.safetensors", auto_checkpoint=False
+    )
+    assert "3" in wf and wf["6"]["inputs"]["text"]
+    assert wf["4"]["inputs"]["ckpt_name"] == "sd_xl_base_1.0.safetensors"
+    assert wf["5"]["inputs"]["width"] == 1024  # XL default size
+
+
+if __name__ == "__main__":
+    test_demo_image_detect_sphere()
+    print("OK test_demo_image_detect_sphere")
+    test_yolo_labels()
+    print("OK test_yolo_labels")
+    test_live_pipeline_from_synth_wav()
+    print("OK test_live_pipeline_from_synth_wav")
+    test_stereo_pseudo_hoa()
+    print("OK test_stereo_pseudo_hoa")
+    test_conditioning_from_fuse()
+    print("OK test_conditioning_from_fuse")

tests/test_phase1_audio.py

@@ -0,0 +1,116 @@
+"""Phase 1: WAV I/O, streams, JSON spatial reports."""
+
+from __future__ import annotations
+
+import json
+import sys
+import tempfile
+from pathlib import Path
+
+import numpy as np
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.analysis import angular_error_deg
+from hoa64.audio_io import read_wav, write_wav
+from hoa64.report import (
+    REPORT_SCHEMA_VERSION,
+    report_from_ambix_wav,
+    report_from_mono_wav,
+    report_from_scene,
+)
+from hoa64.stream import SourceSpec, encode_mono_plane_wave, encode_scene
+from hoa64.synth import envelope_adsr, tone
+
+
+def test_wav_roundtrip(tmp_path: Path | None = None):
+    base = Path(tmp_path) if tmp_path else Path(tempfile.mkdtemp())
+    sr = 16000
+    x = tone(440, 0.2, sr, amplitude=0.3)
+    path = base / "mono.wav"
+    write_wav(path, x, sr)
+    audio, sr2 = read_wav(path)
+    assert sr2 == sr
+    assert audio.shape[0] == 1
+    assert audio.shape[1] == x.shape[0]
+    assert np.corrcoef(audio[0], x)[0, 1] > 0.99
+
+
+def test_mono_encode_report_doa():
+    sr = 24000
+    sig = tone(1000, 0.4, sr, amplitude=0.5) * envelope_adsr(
+        int(0.4 * sr), sr
+    )
+    hoa = encode_mono_plane_wave(sig, 45.0, 0.0, max_order=3)
+    rep = report_from_hoa_local(hoa, sr, az_true=45.0)
+    err = angular_error_deg(45.0, 0.0, rep.doa_az_deg, rep.doa_el_deg)
+    assert err < 5.0, f"DOA err={err} got ({rep.doa_az_deg},{rep.doa_el_deg})"
+    assert rep.schema == REPORT_SCHEMA_VERSION
+    assert rep.duration_sec > 0.3
+    d = rep.to_dict()
+    assert "bands" in d and "frames" in d
+    # JSON serializable
+    json.dumps(d)
+
+
+def report_from_hoa_local(hoa, sr, az_true=None):
+    from hoa64.report import report_from_hoa
+
+    return report_from_hoa(hoa, sr, max_order=3, include_peak_map=True)
+
+
+def test_scene_two_sources_json(tmp_path: Path | None = None):
+    base = Path(tmp_path) if tmp_path else Path(tempfile.mkdtemp())
+    sr = 16000
+    n = int(0.35 * sr)
+    env = envelope_adsr(n, sr)
+    sources = [
+        SourceSpec(-30.0, 0.0, tone(500, 0.35, sr, amplitude=0.5) * env, "A"),
+        SourceSpec(120.0, 10.0, tone(900, 0.35, sr, amplitude=0.35) * env, "B"),
+    ]
+    rep = report_from_scene(sources, sr, max_order=3)
+    assert len(rep.sources_hint) == 2
+    assert rep.energy > 0
+    out = base / "report.json"
+    rep.save(out)
+    loaded = json.loads(out.read_text())
+    assert loaded["schema"] == REPORT_SCHEMA_VERSION
+    assert loaded["n_channels"] == 16  # order 3
+    assert "one_liner" not in loaded
+    assert rep.one_liner().startswith("spatial:")
+
+
+def test_ambix_wav_analyze(tmp_path: Path | None = None):
+    base = Path(tmp_path) if tmp_path else Path(tempfile.mkdtemp())
+    sr = 16000
+    sig = tone(700, 0.25, sr, amplitude=0.4)
+    hoa = encode_mono_plane_wave(sig, -60.0, 5.0, max_order=1)
+    path = base / "bformat.wav"
+    write_wav(path, hoa[:4], sr)
+    rep = report_from_ambix_wav(path, max_order=1, include_bands=True)
+    err = angular_error_deg(-60.0, 5.0, rep.doa_az_deg, rep.doa_el_deg)
+    assert err < 8.0, f"ambix DOA err={err}"
+
+
+def test_mono_wav_cli_path(tmp_path: Path | None = None):
+    base = Path(tmp_path) if tmp_path else Path(tempfile.mkdtemp())
+    sr = 16000
+    sig = tone(800, 0.2, sr, amplitude=0.4)
+    wav = base / "m.wav"
+    write_wav(wav, sig, sr)
+    rep = report_from_mono_wav(wav, 0.0, 0.0, max_order=1)
+    err = angular_error_deg(0.0, 0.0, rep.doa_az_deg, rep.doa_el_deg)
+    assert err < 5.0
+
+
+if __name__ == "__main__":
+    test_wav_roundtrip()
+    print("OK test_wav_roundtrip")
+    test_mono_encode_report_doa()
+    print("OK test_mono_encode_report_doa")
+    test_scene_two_sources_json()
+    print("OK test_scene_two_sources_json")
+    test_ambix_wav_analyze()
+    print("OK test_ambix_wav_analyze")
+    test_mono_wav_cli_path()
+    print("OK test_mono_wav_cli_path")

tests/test_phase2_wigner.py

@@ -0,0 +1,117 @@
+"""Phase 2: Wigner-D HOA rotation accuracy and speed."""
+
+from __future__ import annotations
+
+import sys
+import time
+from pathlib import Path
+
+import numpy as np
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.analysis import angular_error_deg, doa_from_intensity
+from hoa64.encode import encode_points
+from hoa64.rotate import rotate_matrix_order1, rotate_source_directions, rotate_yaw_pitch_roll
+from hoa64.wigner import hoa_rotation_matrix, rotation_matrix_zyx, apply_hoa_rotation
+
+
+def test_order1_wigner_matches_cartesian():
+    a = encode_points([25.0], [-12.0], [1.0])
+    for yaw, pitch, roll in [(90, 0, 0), (0, 45, 0), (30, -20, 15), (180, 0, 0)]:
+        a_w = rotate_yaw_pitch_roll(
+            a, yaw=yaw, pitch=pitch, roll=roll, max_order=1, method="wigner"
+        )
+        a_c = rotate_matrix_order1(a, yaw=yaw, pitch=pitch, roll=roll)
+        np.testing.assert_allclose(a_w[:4], a_c[:4], atol=1e-9, rtol=1e-9)
+
+
+def test_wigner_matches_plane_wave_reencode():
+    """Ground truth: rotate source direction, re-encode."""
+    cases = [
+        (0.0, 0.0, 90.0, 0.0, 0.0),
+        (40.0, -15.0, 35.0, 0.0, 0.0),
+        (10.0, 20.0, 0.0, 40.0, 0.0),
+        (-70.0, 5.0, 20.0, -25.0, 30.0),
+    ]
+    for az0, el0, yaw, pitch, roll in cases:
+        a0 = encode_points([az0], [el0], [1.0], max_order=7)
+        az1, el1 = rotate_source_directions(az0, el0, yaw=yaw, pitch=pitch, roll=roll)
+        a_gt = encode_points([float(az1)], [float(el1)], [1.0], max_order=7)
+        a_w = rotate_yaw_pitch_roll(
+            a0, yaw=yaw, pitch=pitch, roll=roll, max_order=7, method="wigner"
+        )
+        # Relative error on full 64-vector
+        denom = np.linalg.norm(a_gt) + 1e-15
+        rel = np.linalg.norm(a_w[:64] - a_gt[:64]) / denom
+        assert rel < 1e-6, (
+            f"rel={rel:.3e} for src=({az0},{el0}) rot=({yaw},{pitch},{roll})"
+        )
+
+
+def test_wigner_orthogonal_blocks():
+    R = rotation_matrix_zyx(33.0, -17.0, 8.0, degrees=True)
+    M = hoa_rotation_matrix(R, max_order=7)
+    # Each order block should be orthogonal (rotation)
+    for n in range(0, 8):
+        i0 = n * n
+        dim = 2 * n + 1
+        B = M[i0 : i0 + dim, i0 : i0 + dim]
+        I = B.T @ B
+        np.testing.assert_allclose(I, np.eye(dim), atol=1e-8)
+
+
+def test_wigner_faster_than_dense():
+    a = encode_points([15.0], [10.0], [1.0], max_order=7)
+    # warm-up
+    rotate_yaw_pitch_roll(a, yaw=20.0, pitch=10.0, roll=5.0, method="wigner")
+    rotate_yaw_pitch_roll(
+        a, yaw=20.0, pitch=10.0, roll=5.0, method="dense", n_azi=48, n_el=24
+    )
+
+    t0 = time.perf_counter()
+    for _ in range(50):
+        rotate_yaw_pitch_roll(a, yaw=20.0, pitch=10.0, roll=5.0, method="wigner")
+    t_w = time.perf_counter() - t0
+
+    t0 = time.perf_counter()
+    for _ in range(5):
+        rotate_yaw_pitch_roll(
+            a, yaw=20.0, pitch=10.0, roll=5.0, method="dense", n_azi=48, n_el=24
+        )
+    t_d = time.perf_counter() - t0
+    # per-call times
+    tw = t_w / 50
+    td = t_d / 5
+    # Wigner should be substantially faster (typically 50–1000×)
+    assert tw < td, f"wigner {tw:.4f}s not faster than dense {td:.4f}s"
+    print(f"  timing: wigner={tw*1e3:.3f} ms/call  dense={td*1e3:.3f} ms/call  speedup={td/tw:.0f}x")
+
+
+def test_stream_rotation_CT():
+    # (C,T) path
+    from hoa64.encode import encode_plane_waves
+
+    t = np.linspace(0, 1, 32, endpoint=False)
+    sig = np.sin(2 * np.pi * 3 * t)[None, :]
+    hoa = encode_plane_waves([0.0], [0.0], sig, max_order=3)
+    out = rotate_yaw_pitch_roll(hoa, yaw=90.0, max_order=3, method="wigner")
+    assert out.shape[0] >= 16
+    # DOA of first frame energy via products
+    W, Y, Z, X = out[0], out[1], out[2], out[3]
+    I = np.array([np.mean(W * X), np.mean(W * Y), np.mean(W * Z)])
+    n = np.linalg.norm(I)
+    assert n > 1e-9
+
+
+if __name__ == "__main__":
+    test_order1_wigner_matches_cartesian()
+    print("OK test_order1_wigner_matches_cartesian")
+    test_wigner_matches_plane_wave_reencode()
+    print("OK test_wigner_matches_plane_wave_reencode")
+    test_wigner_orthogonal_blocks()
+    print("OK test_wigner_orthogonal_blocks")
+    test_wigner_faster_than_dense()
+    print("OK test_wigner_faster_than_dense")
+    test_stream_rotation_CT()
+    print("OK test_stream_rotation_CT")

tests/test_phase3_vision.py

@@ -0,0 +1,96 @@
+"""Phase 3 vision + Qwythos API smoke tests."""
+
+from __future__ import annotations
+
+import json
+import sys
+import threading
+import time
+from pathlib import Path
+
+import numpy as np
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.analysis import angular_error_deg
+from hoa64.server import TOOL_SCHEMA, handle_analyze
+from hoa64.vision import encode_boxes_to_hoa, fuse_reports, report_from_boxes
+
+
+def test_vision_point_doa():
+    rep = report_from_boxes(
+        [{"az": 45.0, "el": -10.0, "kind": "point", "weight": 1.0, "label": "p"}],
+        max_order=3,
+    )
+    err = angular_error_deg(45.0, -10.0, rep.doa_az_deg, rep.doa_el_deg)
+    # intensity on static field works; peak should be tight
+    perr = angular_error_deg(45.0, -10.0, rep.peak_az_deg, rep.peak_el_deg)
+    assert perr < 8.0, f"peak err={perr} doa=({rep.doa_az_deg},{rep.doa_el_deg})"
+    assert rep.kind == "spatial_vision"
+    assert rep.sources_hint[0]["label"] == "p"
+
+
+def test_vision_box_lobe_near_center():
+    hoa = encode_boxes_to_hoa(
+        [{"az": 0.0, "el": 0.0, "w_deg": 15, "h_deg": 15, "weight": 1.0}],
+        max_order=3,
+        n_azi=48,
+        n_el=24,
+    )
+    assert hoa.shape[0] == 64
+    assert np.linalg.norm(hoa) > 0
+    # W channel should dominate for broad front lobe
+    assert abs(hoa[0]) > 0
+
+
+def test_fuse_agreement():
+    audio = {"doa_az_deg": 10.0, "doa_el_deg": 0.0, "energy": 1.0}
+    vision = {"doa_az_deg": 12.0, "doa_el_deg": 1.0, "energy": 0.8}
+    f = fuse_reports(audio, vision)
+    assert f["agreement"] is True
+    assert f["angular_separation_deg"] < 15
+    assert "one_liner" in f
+
+
+def test_fuse_disagreement():
+    audio = {"doa_az_deg": 0.0, "doa_el_deg": 0.0, "energy": 1.0}
+    vision = {"doa_az_deg": 90.0, "doa_el_deg": 0.0, "energy": 1.0}
+    f = fuse_reports(audio, vision)
+    assert f["agreement"] is False
+    assert f["angular_separation_deg"] > 80
+
+
+def test_handle_analyze_vision_and_demo():
+    d = handle_analyze({"mode": "demo_scene", "order": 2})
+    assert d["schema"].startswith("spatial-hoa")
+    assert "one_liner" in d
+    v = handle_analyze(
+        {
+            "mode": "vision",
+            "order": 3,
+            "boxes": [{"az": -30, "el": 5, "kind": "ray", "label": "x"}],
+        }
+    )
+    assert v["kind"] == "spatial_vision"
+    err = angular_error_deg(-30, 5, v["peak_az_deg"], v["peak_el_deg"])
+    assert err < 10.0
+
+
+def test_tool_schema_shape():
+    assert TOOL_SCHEMA["type"] == "function"
+    assert TOOL_SCHEMA["function"]["name"] == "spatial_analyze"
+
+
+if __name__ == "__main__":
+    test_vision_point_doa()
+    print("OK test_vision_point_doa")
+    test_vision_box_lobe_near_center()
+    print("OK test_vision_box_lobe_near_center")
+    test_fuse_agreement()
+    print("OK test_fuse_agreement")
+    test_fuse_disagreement()
+    print("OK test_fuse_disagreement")
+    test_handle_analyze_vision_and_demo()
+    print("OK test_handle_analyze_vision_and_demo")
+    test_tool_schema_shape()
+    print("OK test_tool_schema_shape")

tests/test_rotate_rnn.py

@@ -0,0 +1,81 @@
+"""Rotation + iterative RNN-stub motion tests."""
+
+from __future__ import annotations
+
+import sys
+from pathlib import Path
+
+import numpy as np
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from hoa64.analysis import angular_error_deg, doa_from_intensity, peak_direction
+from hoa64.encode import encode_points
+from hoa64.rotate import (
+    rotate_matrix_order1,
+    rotate_source_directions,
+    rotate_yaw_pitch_roll,
+)
+from hoa64.rnn_stub import step_rotate, world_from_sources
+
+
+def test_order1_yaw_matches_reencode():
+    a = encode_points([0.0], [0.0], [1.0])
+    # rotate field +90° yaw: front → left
+    a_fast = rotate_matrix_order1(a, yaw=90.0)
+    az, el = doa_from_intensity(a_fast)
+    err = angular_error_deg(90.0, 0.0, az, el)
+    assert err < 1.0, f"intensity DOA after yaw90: ({az},{el}) err={err}"
+
+
+def test_source_dir_rotation_exact_reencode():
+    """Ground truth: rotate source coordinates, re-encode."""
+    az0, el0 = 20.0, -10.0
+    a0 = encode_points([az0], [el0], [1.0])
+    az1, el1 = rotate_source_directions(az0, el0, yaw=45.0)
+    a1 = encode_points([float(az1)], [float(el1)], [1.0])
+    # dense rotate of a0 by same yaw should ≈ a1 for low orders
+    a_rot = rotate_yaw_pitch_roll(a0, yaw=45.0, max_order=3, n_azi=120, n_el=60)
+    # compare order ≤ 3 channels
+    nch = 16
+    rel = np.linalg.norm(a_rot[:nch] - a1[:nch]) / (np.linalg.norm(a1[:nch]) + 1e-12)
+    assert rel < 0.15, f"relative L2 order≤3 after rotate: {rel}"
+
+
+def test_rnn_loop_keeps_tracking_after_turns():
+    """Agent turns through 360° in steps; peak should stay world-stable in intensity
+    when we interpret DOA in head frame correctly.
+
+    After total yaw +90 (agent turns left), a world-front source is at head-right
+    (−90) in listener frame if field is counter-rotated.
+    """
+    st = world_from_sources([0.0], [0.0], [1.0])
+    # 6 steps of +15° agent yaw
+    for _ in range(6):
+        st = step_rotate(st, d_yaw=15.0, max_order=3, dense=True)
+    assert abs(st.yaw - 90.0) < 1e-9
+    az, el = doa_from_intensity(st.hoa)
+    # front source, agent yawed +90 → source appears at az=-90 (right)
+    err = angular_error_deg(-90.0, 0.0, az, el)
+    assert err < 8.0, f"after +90 agent yaw, head-frame DOA=({az},{el}) err={err}"
+
+
+def test_iterative_path_energy_stable_order1():
+    st = world_from_sources([45.0], [0.0], [1.0])
+    e0 = float(np.dot(st.hoa[:4], st.hoa[:4]))
+    for _ in range(24):
+        st = step_rotate(st, d_yaw=15.0, max_order=1, dense=False)
+    e1 = float(np.dot(st.hoa[:4], st.hoa[:4]))
+    # full 360° of order-1 rotations should preserve energy exactly
+    assert abs(e1 - e0) / e0 < 1e-9
+
+
+if __name__ == "__main__":
+    for fn in [
+        test_order1_yaw_matches_reencode,
+        test_source_dir_rotation_exact_reencode,
+        test_rnn_loop_keeps_tracking_after_turns,
+        test_iterative_path_energy_stable_order1,
+    ]:
+        fn()
+        print("OK", fn.__name__)

tools/spatial_analyze.openai.json

@@ -0,0 +1,65 @@
+{
+  "tools": [
+    {
+      "type": "function",
+      "function": {
+        "name": "spatial_analyze",
+        "description": "HOA-7 spatial calculator (not an LLM). Analyze Ambix audio, mono plane-wave audio, synthetic scenes, or vision boxes on the sphere. Returns DOA, energy, bands, frames. Call via HTTP POST http://127.0.0.1:8765/v1/spatial/analyze or bash: spatial-report ...",
+        "parameters": {
+          "type": "object",
+          "properties": {
+            "mode": {
+              "type": "string",
+              "enum": [
+                "ambix_file",
+                "mono_file",
+                "demo_scene",
+                "vision",
+                "fuse",
+                "detect",
+                "live",
+                "condition"
+              ]
+            },
+            "image": { "type": "string", "description": "Image path for detect mode" },
+            "duration": { "type": "number", "description": "Seconds for live mic capture" },
+            "prompt": { "type": "string", "description": "Base diffusion prompt for condition mode" },
+            "report": { "type": "object", "description": "Spatial report object for condition mode" },
+            "path": { "type": "string" },
+            "az": { "type": "number" },
+            "el": { "type": "number" },
+            "order": { "type": "integer", "default": 3 },
+            "boxes": {
+              "type": "array",
+              "items": {
+                "type": "object",
+                "properties": {
+                  "az": { "type": "number" },
+                  "el": { "type": "number" },
+                  "w_deg": { "type": "number" },
+                  "h_deg": { "type": "number" },
+                  "weight": { "type": "number" },
+                  "label": { "type": "string" },
+                  "kind": { "type": "string", "enum": ["box", "ray", "point"] }
+                }
+              }
+            }
+          },
+          "required": ["mode"]
+        }
+      }
+    }
+  ],
+  "http": {
+    "base": "http://127.0.0.1:8765",
+    "analyze": "POST /v1/spatial/analyze",
+    "health": "GET /health",
+    "schema": "GET /v1/spatial/schema"
+  },
+  "bash": {
+    "demo": "spatial-report demo-scene -o /tmp/spatial.json",
+    "ambix": "spatial-report analyze FILE.wav --ambix -o /tmp/spatial.json",
+    "mono": "spatial-report analyze FILE.wav --az 30 --el 0 -o /tmp/spatial.json",
+    "vision": "spatial-report vision --boxes '[{\"az\":0,\"el\":0,\"label\":\"obj\"}]' -o /tmp/vision.json"
+  }
+}