Skip to main content

qualia_core_db/render/
spectral_operator.rs

1//! P7.7 — Unified spectral-operator API surface.
2//!
3//! One entry point for all spectral operations: colour projection, gamut
4//! mapping, metamers, spectral blend, and audio time-frequency surfaces.
5//!
6//! The `SpectralOperator` struct is a zero-allocation facade that dispatches
7//! to the specialised submodules. It is the single API surface that the
8//! renderer, audio engine, and export pipeline call.
9//!
10//! ## Determinism
11//!
12//! All operations are deterministic: identical inputs → bit-identical outputs.
13
14use crate::audio::tf_surface::TfSurface;
15use crate::audio::tf_surface_edit::{
16    apply_gain, copy_patch, crossfade, fade_in, fade_out, pitch_shift, spectral_gate, time_stretch,
17    Region as TfRegion, SurfaceEditError,
18};
19use crate::render::gamut::{gamut_map_clamp, is_in_gamut, linear_srgb_to_xyz};
20use crate::render::metamer::{fibre_spd, is_metameric, metamer_kernel_basis, min_norm_spd_for_xyz};
21use crate::render::spectral_blend::{blend_divergence, spectral_blend_emf, spectral_blend_spd};
22use crate::render::spectral_kernel::{
23    delta_e_76, emf_to_linear_rgb, emf_to_spd, linear_rgb_to_display, spd_to_xyz, xyz_to_lab,
24    xyz_to_linear_srgb, LinearRgb, Spd, Xyz,
25};
26
27// ───────────────────────────────────────────────────────────────────────────
28//  SpectralOperator
29// ───────────────────────────────────────────────────────────────────────────
30
31/// Unified spectral operator — the single entry point for all spectral
32/// operations in the Qualia engine.
33///
34/// This is a zero-allocation facade: it holds no state and dispatches to
35/// the specialised submodules. All methods are `&self` or associated
36/// functions.
37#[derive(Debug, Clone, Copy, Default)]
38pub struct SpectralOperator;
39
40impl SpectralOperator {
41    // ── P7.0: EMF → Colour ───────────────────────────────────────────
42
43    /// Convert an EMF payload `[α, μ, σ]` to a Spectral Power Distribution.
44    #[inline]
45    pub fn emf_to_spd(alpha: f32, mu: f32, sigma: f32) -> Spd {
46        emf_to_spd(alpha, mu, sigma)
47    }
48
49    /// Project an SPD to CIE XYZ tristimulus values.
50    #[inline]
51    pub fn spd_to_xyz(spd: &Spd) -> Xyz {
52        spd_to_xyz(spd)
53    }
54
55    /// Convert CIE XYZ to linear sRGB.
56    #[inline]
57    pub fn xyz_to_linear_srgb(xyz: &Xyz) -> LinearRgb {
58        xyz_to_linear_srgb(xyz)
59    }
60
61    /// Full EMF → linear sRGB pipeline.
62    #[inline]
63    pub fn emf_to_linear_rgb(alpha: f32, mu: f32, sigma: f32) -> LinearRgb {
64        emf_to_linear_rgb(alpha, mu, sigma)
65    }
66
67    /// Convert linear sRGB to 8-bit display sRGB (gamma-encoded).
68    #[inline]
69    pub fn linear_rgb_to_display(rgb: &LinearRgb) -> (u8, u8, u8) {
70        linear_rgb_to_display(rgb)
71    }
72
73    /// Full EMF → 8-bit display sRGB pipeline.
74    #[inline]
75    pub fn emf_to_display_rgb(alpha: f32, mu: f32, sigma: f32) -> (u8, u8, u8) {
76        let rgb = emf_to_linear_rgb(alpha, mu, sigma);
77        linear_rgb_to_display(&rgb)
78    }
79
80    // ── P7.0: Colour difference ──────────────────────────────────────
81
82    /// CIE76 ΔE colour difference between two XYZ values.
83    #[inline]
84    pub fn delta_e(xyz_a: &Xyz, xyz_b: &Xyz) -> f32 {
85        delta_e_76(xyz_a, xyz_b)
86    }
87
88    /// Convert XYZ to CIELAB.
89    #[inline]
90    pub fn xyz_to_lab(xyz: &Xyz) -> (f32, f32, f32) {
91        xyz_to_lab(xyz)
92    }
93
94    // ── P7.1: Metamers ───────────────────────────────────────────────
95
96    /// Compute the kernel basis for metameric-black SPDs.
97    #[inline]
98    pub fn metamer_kernel_basis(
99        out_basis: &mut [f32],
100    ) -> Result<usize, crate::render::metamer::MetamerError> {
101        metamer_kernel_basis(out_basis)
102    }
103
104    /// Compute the minimum-norm SPD for a target XYZ.
105    #[inline]
106    pub fn min_norm_spd_for_xyz(target: &Xyz) -> Spd {
107        min_norm_spd_for_xyz(target)
108    }
109
110    /// Construct a fibre element: particular + Σ c_i · ker_i.
111    #[inline]
112    pub fn fibre_spd(particular: &Spd, basis: &[f32], coeffs: &[f32]) -> Spd {
113        fibre_spd(particular, basis, coeffs)
114    }
115
116    /// Check if an SPD is metameric to a target XYZ.
117    #[inline]
118    pub fn is_metameric(spd: &Spd, target: &Xyz, tolerance: f32) -> bool {
119        is_metameric(spd, target, tolerance)
120    }
121
122    // ── P7.2: Gamut ──────────────────────────────────────────────────
123
124    /// Check if a colour is in the sRGB gamut.
125    #[inline]
126    pub fn is_in_gamut(xyz: &Xyz) -> bool {
127        is_in_gamut(xyz)
128    }
129
130    /// Map an out-of-gamut colour to the closest in-gamut colour.
131    #[inline]
132    pub fn gamut_map(xyz: &Xyz) -> Xyz {
133        gamut_map_clamp(xyz)
134    }
135
136    /// Convert linear sRGB to CIE XYZ.
137    #[inline]
138    pub fn linear_srgb_to_xyz(rgb: &LinearRgb) -> Xyz {
139        linear_srgb_to_xyz(rgb)
140    }
141
142    // ── P7.3: Spectral blend ─────────────────────────────────────────
143
144    /// Blend two SPDs in spectral space.
145    #[inline]
146    pub fn spectral_blend_spd(a: &Spd, b: &Spd, t: f32) -> Spd {
147        spectral_blend_spd(a, b, t)
148    }
149
150    /// Blend two EMF payloads in spectral space and return XYZ.
151    #[inline]
152    pub fn spectral_blend_emf(
153        alpha_a: f32,
154        mu_a: f32,
155        sigma_a: f32,
156        alpha_b: f32,
157        mu_b: f32,
158        sigma_b: f32,
159        t: f32,
160    ) -> Xyz {
161        spectral_blend_emf(alpha_a, mu_a, sigma_a, alpha_b, mu_b, sigma_b, t)
162    }
163
164    /// ΔE divergence between spectral blend and gamma-encoded sRGB lerp.
165    #[inline]
166    pub fn blend_divergence(
167        alpha_a: f32,
168        mu_a: f32,
169        sigma_a: f32,
170        alpha_b: f32,
171        mu_b: f32,
172        sigma_b: f32,
173        t: f32,
174    ) -> f32 {
175        blend_divergence(alpha_a, mu_a, sigma_a, alpha_b, mu_b, sigma_b, t)
176    }
177
178    // ── P7.5: Audio time-frequency surface ───────────────────────────
179
180    /// Create a time-frequency surface view from a raster.
181    #[inline]
182    pub fn tf_surface<'a>(
183        raster: &'a [f32],
184        frame_count: usize,
185        bin_count: usize,
186        sample_rate: u32,
187        hop_size: usize,
188    ) -> TfSurface<'a> {
189        TfSurface::new(raster, frame_count, bin_count, sample_rate, hop_size)
190    }
191
192    // ── P7.6: Audio surface edits ────────────────────────────────────
193
194    /// Apply a gain to a region of the surface.
195    #[inline]
196    pub fn surface_gain(
197        surface: &TfSurface,
198        region: &TfRegion,
199        gain: f32,
200        out: &mut [f32],
201    ) -> Result<usize, SurfaceEditError> {
202        apply_gain(surface, region, gain, out)
203    }
204
205    /// Spectral gate: zero out bins below a threshold.
206    #[inline]
207    pub fn surface_gate(
208        surface: &TfSurface,
209        region: &TfRegion,
210        threshold: f32,
211        out: &mut [f32],
212    ) -> Result<usize, SurfaceEditError> {
213        spectral_gate(surface, region, threshold, out)
214    }
215
216    /// Copy a rectangular patch to a new location.
217    #[inline]
218    pub fn surface_copy_patch(
219        surface: &TfSurface,
220        src_region: &TfRegion,
221        dst_frame: usize,
222        dst_bin: usize,
223        out: &mut [f32],
224    ) -> Result<usize, SurfaceEditError> {
225        copy_patch(surface, src_region, dst_frame, dst_bin, out)
226    }
227
228    /// Time-stretch by resampling along the time axis.
229    #[inline]
230    pub fn surface_time_stretch(
231        surface: &TfSurface,
232        factor: f32,
233        out: &mut [f32],
234    ) -> Result<(usize, usize), SurfaceEditError> {
235        time_stretch(surface, factor, out)
236    }
237
238    /// Pitch-shift by resampling along the frequency axis.
239    #[inline]
240    pub fn surface_pitch_shift(
241        surface: &TfSurface,
242        factor: f32,
243        out: &mut [f32],
244    ) -> Result<usize, SurfaceEditError> {
245        pitch_shift(surface, factor, out)
246    }
247
248    /// Crossfade two surfaces.
249    #[inline]
250    pub fn surface_crossfade(
251        surface_a: &TfSurface,
252        surface_b: &TfSurface,
253        t: f32,
254        out: &mut [f32],
255    ) -> Result<usize, SurfaceEditError> {
256        crossfade(surface_a, surface_b, t, out)
257    }
258
259    /// Fade in over the first `fade_frames` frames.
260    #[inline]
261    pub fn surface_fade_in(
262        surface: &TfSurface,
263        fade_frames: usize,
264        out: &mut [f32],
265    ) -> Result<usize, SurfaceEditError> {
266        fade_in(surface, fade_frames, out)
267    }
268
269    /// Fade out over the last `fade_frames` frames.
270    #[inline]
271    pub fn surface_fade_out(
272        surface: &TfSurface,
273        fade_frames: usize,
274        out: &mut [f32],
275    ) -> Result<usize, SurfaceEditError> {
276        fade_out(surface, fade_frames, out)
277    }
278
279    // ── Batch operations ─────────────────────────────────────────────
280
281    /// Batch EMF → display RGB: process N EMF payloads into N display RGB
282    /// triples. `emf` is `[α, μ, σ]` triples (3*N floats), `out` is N*3 u8.
283    pub fn batch_emf_to_display(emf: &[f32], out: &mut [u8]) -> usize {
284        let n = emf.len() / 3;
285        let n = n.min(out.len() / 3);
286        for i in 0..n {
287            let rgb = emf_to_linear_rgb(emf[i * 3], emf[i * 3 + 1], emf[i * 3 + 2]);
288            let (r, g, b) = linear_rgb_to_display(&rgb);
289            out[i * 3] = r;
290            out[i * 3 + 1] = g;
291            out[i * 3 + 2] = b;
292        }
293        n
294    }
295
296    /// Batch EMF → XYZ: process N EMF payloads into N XYZ triples.
297    /// `emf` is `[α, μ, σ]` triples (3*N floats), `out` is N*3 f32.
298    pub fn batch_emf_to_xyz(emf: &[f32], out: &mut [f32]) -> usize {
299        let n = emf.len() / 3;
300        let n = n.min(out.len() / 3);
301        for i in 0..n {
302            let spd = emf_to_spd(emf[i * 3], emf[i * 3 + 1], emf[i * 3 + 2]);
303            let xyz = spd_to_xyz(&spd);
304            out[i * 3] = xyz.x;
305            out[i * 3 + 1] = xyz.y;
306            out[i * 3 + 2] = xyz.z;
307        }
308        n
309    }
310
311    /// Batch gamut mapping: map N XYZ triples to in-gamut XYZ.
312    /// `out` is N*3 f32.
313    pub fn batch_gamut_map(xyz: &[f32], out: &mut [f32]) -> usize {
314        let n = xyz.len() / 3;
315        let n = n.min(out.len() / 3);
316        for i in 0..n {
317            let mapped = gamut_map_clamp(&Xyz::new(xyz[i * 3], xyz[i * 3 + 1], xyz[i * 3 + 2]));
318            out[i * 3] = mapped.x;
319            out[i * 3 + 1] = mapped.y;
320            out[i * 3 + 2] = mapped.z;
321        }
322        n
323    }
324}
325
326// ───────────────────────────────────────────────────────────────────────────
327//  Tests
328// ───────────────────────────────────────────────────────────────────────────
329
330#[cfg(test)]
331mod tests {
332    use super::*;
333
334    #[test]
335    fn operator_emf_to_display_deterministic() {
336        let (r1, g1, b1) = SpectralOperator::emf_to_display_rgb(1.0, 0.3, 0.5);
337        let (r2, g2, b2) = SpectralOperator::emf_to_display_rgb(1.0, 0.3, 0.5);
338        assert_eq!((r1, g1, b1), (r2, g2, b2));
339    }
340
341    #[test]
342    fn operator_batch_emf_to_display() {
343        let emf = [1.0, 0.0, 0.0, 1.0, 0.0, 0.5, 1.0, 0.0, 1.0];
344        let mut out = [0u8; 9];
345        let n = SpectralOperator::batch_emf_to_display(&emf, &mut out);
346        assert_eq!(n, 3);
347        // Each triple is valid u8, so it's guaranteed to be <= 255.
348    }
349
350    #[test]
351    fn operator_batch_emf_to_xyz() {
352        let emf = [1.0, 0.0, 0.0, 1.0, 0.0, 0.5, 1.0, 0.0, 1.0];
353        let mut out = [0.0f32; 9];
354        let n = SpectralOperator::batch_emf_to_xyz(&emf, &mut out);
355        assert_eq!(n, 3);
356        // All values should be finite.
357        for v in &out {
358            assert!(v.is_finite());
359        }
360    }
361
362    #[test]
363    fn operator_batch_gamut_map() {
364        let xyz = [2.0, 0.5, 0.0, 0.3, 0.4, 0.5];
365        let mut out = [0.0f32; 6];
366        let n = SpectralOperator::batch_gamut_map(&xyz, &mut out);
367        assert_eq!(n, 2);
368        // First was out-of-gamut, should be mapped in.
369        let mapped = Xyz::new(out[0], out[1], out[2]);
370        assert!(SpectralOperator::is_in_gamut(&mapped));
371    }
372
373    #[test]
374    fn operator_metamer_round_trip() {
375        let target = Xyz::new(0.4, 0.5, 0.3);
376        let spd = SpectralOperator::min_norm_spd_for_xyz(&target);
377        let reprojected = SpectralOperator::spd_to_xyz(&spd);
378        assert!((reprojected.x - target.x).abs() < 0.05);
379        assert!((reprojected.y - target.y).abs() < 0.05);
380        assert!((reprojected.z - target.z).abs() < 0.05);
381    }
382
383    #[test]
384    fn operator_blend_pipeline() {
385        let xyz = SpectralOperator::spectral_blend_emf(1.0, 0.1, 0.2, 1.0, 0.1, 0.8, 0.5);
386        assert!(xyz.x.is_finite() && xyz.y.is_finite() && xyz.z.is_finite());
387    }
388
389    #[test]
390    fn operator_full_pipeline_emf_to_display() {
391        // EMF → SPD → XYZ → linear sRGB → display sRGB
392        for i in 0..=10 {
393            let sigma = i as f32 / 10.0;
394            let (_r, _g, _b) = SpectralOperator::emf_to_display_rgb(1.0, 0.2, sigma);
395            // Display RGB returns u8, so it's always <= 255.
396        }
397    }
398
399    #[test]
400    fn operator_delta_e_self_zero() {
401        let xyz = Xyz::new(0.3, 0.5, 0.2);
402        assert!(SpectralOperator::delta_e(&xyz, &xyz) < 1e-6);
403    }
404
405    #[test]
406    fn operator_surface_edits_via_facade() {
407        use crate::audio::audio_spectral_sheet::SPECTRAL_PREVIEW_BINS;
408        let frames = 4;
409        let bins = SPECTRAL_PREVIEW_BINS;
410        let raster = vec![0.5f32; frames * bins];
411        let s = SpectralOperator::tf_surface(&raster, frames, bins, 44100, 512);
412        let region = TfRegion::full(frames, bins);
413        let mut out = vec![0.0f32; frames * bins];
414        SpectralOperator::surface_gain(&s, &region, 2.0, &mut out).unwrap();
415        assert!((out[0] - 1.0).abs() < 1e-6, "gain should double value");
416    }
417}