Files
DarkRoom/core/dr-film/src/spectrum.rs
T
dtourolleandClaude Opus 5 f14176de29 Give the projector lamp a name instead of a warning
Every bake of a Vision3 stock logged:

    unknown illuminant "K75P", falling back to D55

K75P is a cinema xenon short-arc lamp, and Kodak 2383 and 2393 -- the
projection print films those stocks print onto -- name it as the light their
result is looked at under. It was never implemented, so it fell through to
the unknown branch.

The fallback was the right family: a xenon arc sits near 6000 K, close to
daylight and nothing like the tungsten enlarger above it. So the pixels do not
move. What changes is that D55 is now a documented choice rather than the
consolation prize for an unrecognised string, with the approximation stated --
an arc has line structure a Planckian curve cannot express, and the residue of
that is small here because the viewing step adapts the white point out either
way.

The test is the point of the commit. A profile naming a light nobody
implemented should fail the suite, not whisper into a log that only gets read
when somebody happens to be looking for something else -- which is how this
was found.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
2026-08-26 17:41:34 +02:00

267 lines
10 KiB
Rust
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
//! The fixed colour science: spectra in, tristimulus out.
//!
//! Everything here is per-installation rather than per-stock — one observer,
//! a handful of illuminants, one basis — and all of it is small enough to
//! compile in. See [`crate::tables`] for the numbers themselves.
use crate::tables::{ILLUMINANT_D50, ILLUMINANT_D55, ILLUMINANT_D65, OBSERVER, SPECTRUM};
/// A spectral distribution on the crate's fixed 380–780 nm, 5 nm grid.
pub type Spectrum = [f32; SPECTRUM];
/// Linear sRGB primaries from CIE XYZ, for a D65 white.
///
/// The adaptation to whatever white the picture is actually being viewed under
/// happens before this is applied — see [`Viewing::to_srgb`].
const XYZ_TO_SRGB: [[f32; 3]; 3] = [
[3.2404542, -1.5371385, -0.4985314],
[-0.969266, 1.8760108, 0.041556],
[0.0556434, -0.2040259, 1.0572252],
];
/// sRGB's own white, which is what [`XYZ_TO_SRGB`] expects to be handed.
const D65_WHITE: [f32; 3] = [0.9504559, 1.0, 1.0890578];
/// Bradford cone response, and its inverse.
///
/// A von Kries adaptation done in XYZ — dividing each channel by the white —
/// is the obvious thing and it is wrong: XYZ axes are not cone responses, so
/// the result drifts in hue. Bradford is the transform that makes "the same
/// colour, seen under a different light" mean what a viewer means by it, and
/// it matters here because a print is specified under D50 while sRGB is D65.
const BRADFORD: [[f32; 3]; 3] = [
[0.8951, 0.2664, -0.1614],
[-0.7502, 1.7135, 0.0367],
[0.0389, -0.0685, 1.0296],
];
const BRADFORD_INV: [[f32; 3]; 3] = [
[0.9869929, -0.1470543, 0.1599627],
[0.4323053, 0.5183603, 0.0492912],
[-0.0085287, 0.0400428, 0.9684867],
];
/// Named illuminants a profile may ask for.
///
/// The daylight ones are tabulated because the CIE defines them that way. The
/// tungsten ones are computed from Planck's law, which costs no data at all,
/// and the approximation is unusually safe here: an enlarger's lamp colour is
/// cancelled by the filtration solved in [`crate::bake::print_balance`], the
/// same way a darkroom worker dials it out on the colour head.
pub fn illuminant(name: &str) -> Spectrum {
match name {
"D50" => ILLUMINANT_D50,
"D55" => ILLUMINANT_D55,
"D65" => ILLUMINANT_D65,
// The tungsten-halogen enlarger source, with and without the heat
// filter that a real head carries. Both land on the same blackbody:
// the filter's effect is a colour shift, and a colour shift ahead of
// the balance step is by construction invisible.
"TH-KG3" | "TH-KG3-L" | "T" => blackbody(3400.0),
// A cinema projector's xenon short-arc lamp, which is what a release
// print is *looked at* under — Kodak 2383 and 2393 name it.
//
// Approximated by D55, and named here rather than left to fall through
// to the unknown branch, which logged a warning on every bake of a
// Vision3 stock and implied something was broken. A xenon arc sits near
// 6000 K, close to daylight and nothing like the tungsten above, so
// daylight is the right family to stand in for it. It is not a
// blackbody — the arc has line structure a Planckian curve has no way
// to express — and the residue of that is small here because the
// viewing step adapts the white point out either way.
"K75P" => ILLUMINANT_D55,
other => {
if let Some(kelvin) = other.strip_prefix("BB").and_then(|k| k.parse::<f32>().ok()) {
return blackbody(kelvin);
}
log::warn!("unknown illuminant {other:?}, falling back to D55");
ILLUMINANT_D55
}
}
}
/// Planck's law on the grid, normalised to unit mean.
///
/// Normalised because only the shape matters: absolute level is set by the
/// exposure controls, and leaving it in would make the choice of units a
/// visible parameter.
pub fn blackbody(kelvin: f32) -> Spectrum {
const H: f64 = 6.626_070_15e-34;
const C: f64 = 2.997_924_58e8;
const KB: f64 = 1.380_649e-23;
let mut out = [0.0f32; SPECTRUM];
let mut total = 0.0f64;
for (i, slot) in out.iter_mut().enumerate() {
let lambda =
f64::from(crate::tables::LAMBDA_MIN + crate::tables::LAMBDA_STEP * i as f32) * 1e-9;
let radiance = (2.0 * H * C * C)
/ (lambda.powi(5) * ((H * C / (lambda * KB * f64::from(kelvin))).exp() - 1.0));
*slot = radiance as f32;
total += radiance;
}
let mean = (total / SPECTRUM as f64) as f32;
for slot in &mut out {
*slot /= mean;
}
out
}
/// How a developed image is looked at: an illuminant, and the adaptation it
/// implies.
///
/// Built once per bake rather than per sample, because the white point and the
/// adaptation matrix depend only on the illuminant and computing them inside
/// the LUT loop would be the same work 32 768 times.
pub struct Viewing {
illuminant: Spectrum,
/// The illuminant's own Y, which normalises the integral so that a clear
/// frame comes out at exactly 1.0 rather than at whatever the tabulated
/// units happen to give.
normalisation: f32,
/// XYZ under this illuminant to linear sRGB, adaptation folded in.
matrix: [[f32; 3]; 3],
}
impl Viewing {
pub fn new(illuminant_name: &str) -> Self {
let illuminant = illuminant(illuminant_name);
let mut white = [0.0f32; 3];
let mut normalisation = 0.0f32;
for i in 0..SPECTRUM {
normalisation += illuminant[i] * OBSERVER[i][1];
for c in 0..3 {
white[c] += illuminant[i] * OBSERVER[i][c];
}
}
for c in &mut white {
*c /= normalisation;
}
Self {
illuminant,
normalisation,
matrix: mat3_mul(XYZ_TO_SRGB, bradford_adaptation(white, D65_WHITE)),
}
}
/// Integrate a transmittance (or reflectance) against the illuminant and
/// the observer, and convert to linear sRGB.
pub fn to_srgb(&self, transmittance: &Spectrum) -> [f32; 3] {
let mut xyz = [0.0f32; 3];
for i in 0..SPECTRUM {
let light = transmittance[i] * self.illuminant[i];
for c in 0..3 {
xyz[c] += light * OBSERVER[i][c];
}
}
for c in &mut xyz {
*c /= self.normalisation;
}
mat3_apply(self.matrix, xyz)
}
}
/// The Bradford transform taking `from` white to `to` white.
fn bradford_adaptation(from: [f32; 3], to: [f32; 3]) -> [[f32; 3]; 3] {
let src = mat3_apply(BRADFORD, from);
let dst = mat3_apply(BRADFORD, to);
let scale = [
[dst[0] / src[0], 0.0, 0.0],
[0.0, dst[1] / src[1], 0.0],
[0.0, 0.0, dst[2] / src[2]],
];
mat3_mul(BRADFORD_INV, mat3_mul(scale, BRADFORD))
}
pub(crate) fn mat3_apply(m: [[f32; 3]; 3], v: [f32; 3]) -> [f32; 3] {
[
m[0][0] * v[0] + m[0][1] * v[1] + m[0][2] * v[2],
m[1][0] * v[0] + m[1][1] * v[1] + m[1][2] * v[2],
m[2][0] * v[0] + m[2][1] * v[1] + m[2][2] * v[2],
]
}
pub(crate) fn mat3_mul(a: [[f32; 3]; 3], b: [[f32; 3]; 3]) -> [[f32; 3]; 3] {
let mut out = [[0.0f32; 3]; 3];
for (r, row) in out.iter_mut().enumerate() {
for (c, slot) in row.iter_mut().enumerate() {
*slot = (0..3).map(|k| a[r][k] * b[k][c]).sum();
}
}
out
}
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn a_clear_frame_is_white_under_every_viewing_illuminant() {
// The property the whole viewing step exists to hold: an unexposed
// transparency on a light table is white, whatever the light table is.
// Normalising on luminance alone would leave the illuminant's own
// colour in the result, and a D50 print would render warm.
for name in ["D50", "D55", "D65", "TH-KG3"] {
let rgb = Viewing::new(name).to_srgb(&[1.0; SPECTRUM]);
for c in rgb {
assert!(
(c - 1.0).abs() < 0.01,
"{name}: clear frame gave {rgb:?}, which is not neutral white"
);
}
}
}
#[test]
fn adapting_a_white_to_itself_changes_nothing() {
let m = bradford_adaptation(D65_WHITE, D65_WHITE);
for r in 0..3 {
for c in 0..3 {
let expected = if r == c { 1.0 } else { 0.0 };
assert!(
(m[r][c] - expected).abs() < 1e-5,
"{m:?} is not the identity"
);
}
}
}
#[test]
fn every_illuminant_a_shipped_profile_names_is_known() {
// The warning that found this said "unknown illuminant K75P, falling
// back to D55" on every bake of a Vision3 stock. The fallback was
// reasonable and the silence was not: a profile naming a light nobody
// implemented should be a failing test rather than a line in a log.
for name in ["D50", "D55", "D65", "TH-KG3", "TH-KG3-L", "T", "K75P"] {
let rgb = Viewing::new(name).to_srgb(&[1.0; SPECTRUM]);
assert!(
rgb.iter().all(|c| (c - 1.0).abs() < 0.01),
"{name} is not a working illuminant: clear frame gave {rgb:?}"
);
}
}
#[test]
fn a_hotter_blackbody_is_bluer() {
// Cheap, but it is the one property that catches Planck's law written
// with a sign or a reciprocal wrong, which otherwise produces a
// plausible-looking curve pointing the wrong way.
let cool = blackbody(2800.0);
let hot = blackbody(9000.0);
let blue = 10; // 430nm
let red = 70; // 730nm
assert!(hot[blue] / hot[red] > cool[blue] / cool[red]);
}
#[test]
fn the_observer_integrates_to_a_plausible_white() {
// Guards the generated table against a transposed or misaligned paste:
// D65 through the 1931 observer must land on sRGB's own white.
let v = Viewing::new("D65");
let rgb = v.to_srgb(&[1.0; SPECTRUM]);
assert!(rgb.iter().all(|c| (c - 1.0).abs() < 0.01), "{rgb:?}");
}
}