//! TRACES: FR-DEV-3e //! Camera input profiles — deciding what a sensor's numbers are supposed to //! mean. //! //! A RAW file is a count of electrons per photosite. Turning that into a //! colour needs two things the file cannot supply on its own: a **matrix** //! saying how this sensor's three responses relate to the CIE observer, and a //! **rendering** saying what to do with the resulting scene-referred values so //! that a photograph looks like a photograph. This module supplies the first //! and looks up the second ([`crate::base_curve`]). //! //! # What is extracted, and from where //! //! | Quantity | DNG tag | Where rawler keeps it | //! |---|---|---| //! | `ColorMatrix1/2` | 50721 / 50722 | `RawImage::color_matrix`, keyed by illuminant | //! | `CalibrationIlluminant1/2` | 50778 / 50779 | *is* that map's key | //! | `ForwardMatrix1/2` | 50964 / 50965 | nowhere — read here, from the IFD | //! | `AsShotNeutral` | 50728 | reciprocated into `RawImage::wb_coeffs` | //! //! The interesting row is the third. rawler 0.7.2 parses `ColorMatrix1/2` and //! their illuminants for DNG files and drops `ForwardMatrix1/2` on the floor, //! so that pair is read straight out of the root IFD here — the decoder hands //! out a parsed `IFD` through [`rawler::decoders::Decoder::ifd`], which means //! no TIFF parser of our own and no second interpretation of the same bytes. //! //! The second row is the one worth not misreading. rawler's map is keyed by //! `Illuminant`, so `CalibrationIlluminant1/2` are not *lost* by being absent //! from the table above — they are the keys. That is also true for bodies that //! are not DNG at all: rawler's camera database stores the same per-illuminant //! matrices for native formats, so a Canon CR2 arrives with an `A` matrix and //! a `D65` matrix exactly as an Adobe-converted DNG of the same frame would. //! Everything below therefore applies to every supported body, not only to //! DNGs. //! //! # Two rawler fields that are traps, and stay traps //! //! Both measured on a Canon 6D CR2 (2026-08-09), and neither is used here: //! //! 1. `RawImage::xyz_to_cam` is **all zeros**. It carries an upstream //! deprecation note and 0.7.2 no longer fills it. Reading it silently //! yields no colour transform at all — a black frame, from code that looks //! exactly right. //! 2. `cam_to_xyz_normalized()` divides each of four rows by its own sum, and //! the fourth row (emerald/white, unused on any Bayer body) sums to zero. //! Every element came back `NaN`. Inverting the 3×3 ourselves avoids the //! fourth channel entirely, which is what [`crate::invert3`] is for. //! //! # The rendering this replaces //! //! Before this module, the decoder took whichever matrix was keyed `D65`, //! ignored the other, and inverted it. That is the dcraw default and it is the //! rendering FR-DEV-3e exists to get away from: correct in the abstract, flat //! and poor on skin in practice, because a single daylight matrix is asked to //! describe a sensor under every illuminant it will ever see. //! //! It also rendered some bodies with no matrix at all. The lookup asked for //! `D65`, then `A`, and gave up — so a Phase One IQ3, calibrated at `D55` and //! `D75`, went through the pipeline uncalibrated despite carrying two perfectly //! good matrices. Placing every illuminant on a temperature axis fixes that as //! a side effect of making the interpolation possible at all. //! //! # What is deliberately not here yet //! //! FR-DEV-3e defers full `.dcp` support — `HueSatDeltas` and //! `ProfileLookTable` — and requires that they arrive as *additions* rather //! than as a pipeline reordering. They would: both are lookups applied to a //! colour after this matrix and before, or alongside, the base curve, so they //! extend [`CameraProfile`] with more calibration data and extend the shader's //! camera-profile stage with more work. Nothing above would move. use crate::{cam_to_srgb_from, invert3}; use rawler::imgop::xyz::Illuminant; /// One calibration: a matrix, and the light it was measured under. /// /// `forward` is the DNG `ForwardMatrix`, present only where the file carries /// one. It is a different and better-conditioned statement of the same /// relationship — see [`CameraProfile::cam_to_srgb`] — and its absence is /// ordinary rather than exceptional, because no native RAW format has the tag /// and rawler's camera database does not carry one either. #[derive(Debug, Clone, Copy, PartialEq)] pub struct Calibration { /// Correlated colour temperature of the calibration illuminant, in kelvin. pub temperature: f32, /// CIE XYZ → camera RGB, as the tag defines it. /// /// The DNG specification states `ColorMatrix` against D50, its connection /// space. [`crate::cam_to_srgb_from`] instead normalises the rows by the /// camera's response to *sRGB* white and then inverts, which is the dcraw /// treatment: it performs the adaptation in camera space rather than in /// XYZ. The two differ by a chromatic adaptation that the row /// normalisation absorbs, and the property that has to hold — camera /// neutral lands on sRGB neutral — holds either way. Doing it the other /// way round is a colour-management concern (ARCH §5.2) and would want an /// ICC pipeline behind it rather than a different constant here. pub xyz_to_cam: [[f32; 3]; 3], /// White-balanced camera RGB → CIE XYZ (D50), where the file says so. pub forward: Option<[[f32; 3]; 3]>, } /// TRACES: FR-DEV-3e /// Everything known about how this body sees colour. /// /// Held as *calibrations plus a neutral* rather than as a finished matrix /// because the finished matrix depends on the scene: a two-illuminant profile /// describes the sensor under tungsten and under daylight, and the frame in /// hand was shot under neither. [`Self::cam_to_srgb`] is where that gets /// resolved. #[derive(Debug, Clone, PartialEq)] pub struct CameraProfile { /// Ordered coolest last; two in practice, since that is all rawler /// surfaces. Never empty — [`Self::extract`] returns `None` rather than /// build an empty profile, so that "no profile" and "a profile that says /// nothing" cannot be confused downstream. calibrations: Vec, /// `AsShotNeutral`: the camera-space colour the camera considered white. /// /// Optional because a body that wrote no white balance leaves nothing to /// estimate the scene's illuminant from. The interpolation then falls back /// to the coolest calibration, which is the daylight one on every body /// this has been seen on. neutral: Option<[f32; 3]>, } /// The DNG tags rawler parses but does not surface, read from the root IFD. /// /// Kept as a separate step from [`CameraProfile::extract`] because it needs /// the *decoder*, which exists only for as long as the file is open, while the /// profile is built from the decoded `RawImage` and outlives it. #[derive(Debug, Clone, Copy, Default, PartialEq)] pub struct DngMatrices { /// `(CalibrationIlluminant1, ForwardMatrix1)`, and the same for 2. /// /// The illuminant is carried alongside because the pairing is *positional* /// in the file — `ForwardMatrix2` belongs to `CalibrationIlluminant2` — /// and rawler's illuminant-keyed map has already thrown that ordering /// away. Reading the illuminant here is what lets a forward matrix be /// matched back to the colour matrix it was measured with. pub forward: [Option<(Illuminant, [[f32; 3]; 3])>; 2], /// `AsShotNeutral`, in camera RGB. /// /// rawler already reciprocates this into `wb_coeffs`, so this is a /// cross-check rather than the only source — but it is the source that is /// exactly what the file said, with no green normalisation applied and no /// `NaN` fourth channel to strip. pub as_shot_neutral: Option<[f32; 3]>, } impl CameraProfile { /// TRACES: FR-DEV-3e /// Build a profile from a decoded image and whatever the file's IFD added. /// /// `None` where the body has no usable matrix at all — an unknown camera, /// or rawler's all-zero placeholder. That is a real answer and the caller /// must keep treating it as one: rendering uncalibrated is recoverable, /// rendering through a matrix that is secretly zeros is a black frame. pub fn extract(image: &rawler::RawImage, dng: &DngMatrices) -> Option { let mut calibrations: Vec = Vec::new(); for (illuminant, flat) in &image.color_matrix { if flat.len() < 9 { // A four-row matrix (RGBE sensors) would be 12; anything // shorter than nine is not a 3×3 and there is nothing // sensible to do with it. continue; } let xyz_to_cam = [ [flat[0], flat[1], flat[2]], [flat[3], flat[4], flat[5]], [flat[6], flat[7], flat[8]], ]; if !usable(&xyz_to_cam) { continue; } let Some(temperature) = illuminant_temperature(*illuminant) else { // A calibration under a light with no defined temperature — // `Flash`, `Unknown` — cannot be placed on the axis the // interpolation runs along. Skipping it is better than // inventing a temperature, which would drag the interpolation // toward a point that means nothing. continue; }; // The forward matrix measured under *this* illuminant, if the file // carried one. Matched by illuminant rather than by position, // because position is what the map key already discarded. let forward = dng .forward .iter() .flatten() .find(|(i, _)| i == illuminant) .map(|(_, m)| *m) .filter(usable); calibrations.push(Calibration { temperature, xyz_to_cam, forward, }); } if calibrations.is_empty() { return None; } // Coolest last, so "the daylight one" is `.last()` and the // interpolation below can talk about a low and a high end without // re-sorting. `HashMap` iteration order is not defined, so without // this the same file would produce different profiles between runs. calibrations.sort_by(|a, b| a.temperature.total_cmp(&b.temperature)); let neutral = dng .as_shot_neutral .or_else(|| neutral_from_coefficients(image.wb_coeffs)) .filter(|n| n.iter().all(|v| v.is_finite() && *v > 0.0)); Some(Self { calibrations, neutral, }) } /// Build a profile directly from calibrations, for tests and for callers /// that already hold the numbers. pub fn new(calibrations: Vec, neutral: Option<[f32; 3]>) -> Option { if calibrations.is_empty() { return None; } let mut calibrations = calibrations; calibrations.sort_by(|a, b| a.temperature.total_cmp(&b.temperature)); Some(Self { calibrations, neutral, }) } /// TRACES: FR-DEV-3e /// The scene's estimated colour temperature, in kelvin. /// /// # Why this is a fixed point rather than a calculation /// /// The temperature is read off the as-shot neutral — the camera-space /// colour the camera called white — by carrying it into XYZ and asking /// which black body it sits nearest. But carrying it into XYZ needs a /// matrix, and *which* matrix is exactly what the temperature is being /// computed to decide. The two definitions are circular. /// /// Adobe's DNG SDK resolves it by iterating, and so does this: guess /// daylight, interpolate, re-read the temperature, repeat. It converges in /// two or three rounds because the interpolation weight is a smooth, /// shallow function of the temperature — a hundred kelvin of error in the /// guess moves the next estimate by far less than a hundred kelvin. /// /// Three rounds, fixed, rather than a convergence test. The remaining /// disagreement after three is smaller than the difference between two /// bodies' idea of the same illuminant, and a loop that ran until it /// settled would spend its time proving that. pub fn scene_temperature(&self) -> f32 { // A single calibration has nothing to interpolate between, so its own // illuminant is the answer and the iteration would be busy-work. if self.calibrations.len() < 2 { return self.calibrations[0].temperature; } let Some(neutral) = self.neutral else { // No white balance to read: assume the light the coolest // calibration was measured under, which is D65 on every // dual-illuminant body rawler knows. Not the midpoint — a midpoint // is a scene lit by nothing anyone has photographed, whereas // daylight is the single most likely answer, and it is also what // this decoder assumed before it could interpolate at all. return self.calibrations.last().expect("checked above").temperature; }; let mut temperature = self.calibrations.last().expect("checked above").temperature; for _ in 0..3 { let m = self.interpolated_xyz_to_cam_at(temperature); let Some((x, y)) = neutral_to_xy(&neutral, &m) else { // A singular interpolant. Keep the previous estimate rather // than propagating a NaN into the matrix that renders the // frame. break; }; temperature = cct_from_xy(x, y); } temperature } /// TRACES: FR-DEV-3e /// The XYZ→camera matrix for the light this frame was actually shot under. /// /// Exposed because two callers need the *same* answer: the colour /// conversion, and the white-balance fallback that runs when a file /// carries no coefficients. Two readers disagreeing about which matrix a /// body uses would balance to one white and convert from another, which /// looks like a decode fault rather than like a disagreement. pub fn xyz_to_cam(&self) -> [[f32; 3]; 3] { self.interpolated_xyz_to_cam_at(self.scene_temperature()) } /// The interpolation itself, at a stated temperature. /// /// **In mireds, not kelvin.** Reciprocal temperature is the axis on which /// equal steps are equally visible — 2000 K to 2500 K is an enormous /// change in the light and 9000 K to 9500 K is barely one — so a linear /// blend in kelvin would put nearly all of its travel in the tungsten end /// and leave daylight to a sliver at the top. Adobe's DNG specification /// interpolates in mireds for the same reason, so this also agrees with /// what the profile's author was picturing. fn interpolated_xyz_to_cam_at(&self, temperature: f32) -> [[f32; 3]; 3] { if let [only] = self.calibrations.as_slice() { return only.xyz_to_cam; } let (low, high) = self.bracketing(temperature); lerp3( &low.xyz_to_cam, &high.xyz_to_cam, blend(low, high, temperature), ) } /// The two calibrations this temperature sits between. /// /// Two is the whole story for every file rawler can produce — the DNG /// decoder reads `ColorMatrix1/2` and stops, and no camera in its database /// carries more than two usable illuminants. Bracketing rather than /// hard-coding the pair anyway, because `ColorMatrix3` exists in DNG 1.6 /// and the day rawler starts reading it, a profile with a tungsten, a /// daylight and a shade calibration must interpolate between the two the /// scene actually falls between rather than across the whole range. /// /// Outside the calibrated range the outermost pair is returned and /// [`blend`] clamps, which is the same answer either way. fn bracketing(&self, temperature: f32) -> (&Calibration, &Calibration) { let cals = &self.calibrations; debug_assert!(cals.len() >= 2, "callers handle the single-matrix case"); let last_pair = cals.len() - 2; let i = cals .iter() .rposition(|c| c.temperature <= temperature) .unwrap_or(0) .min(last_pair); (&cals[i], &cals[i + 1]) } /// TRACES: FR-DEV-3e /// Camera RGB → linear sRGB (D65), row-major, ready for the shader. /// /// # Two routes, and why the forward matrix wins where it exists /// /// The colour-matrix route inverts a measurement of the *sensor* to get a /// transform of the *scene*: it says what camera values a given XYZ would /// produce, and rendering needs the other direction. Inversion is exact /// arithmetic on an inexact measurement, and it amplifies whatever error /// the measurement had — most visibly in saturated reds and in the /// near-neutrals that skin is made of. /// /// `ForwardMatrix` is the manufacturer's answer to that: the same /// relationship measured in the direction rendering actually wants, from /// white-balanced camera RGB straight to XYZ, and constrained so that a /// neutral maps exactly onto the D50 white point. Where a file carries /// one, using it is strictly better than inverting the colour matrix, and /// it costs one matrix multiply less. /// /// Both routes end by normalising the rows to sum to one. That is not /// cosmetic: the shader applies green-normalised as-shot multipliers /// before this matrix, so camera neutral arrives as (1, 1, 1), and a /// matrix whose rows did not sum to one would render it as something other /// than white. It is also what makes the two routes interchangeable — they /// agree about neutral by construction, so a body gaining a forward matrix /// changes its colour rendering without changing its exposure. pub fn cam_to_srgb(&self) -> Option<[f32; 9]> { let temperature = self.scene_temperature(); // The forward route needs a forward matrix at *both* ends, or at the // only end there is. Interpolating a forward matrix against an // inverted colour matrix would mix two different statements of the // same thing and land between them, which is worse than either. if let Some(forward) = self.interpolated_forward_at(temperature) { return forward_to_srgb(&forward); } cam_to_srgb_from(&self.interpolated_xyz_to_cam_at(temperature)) } /// The forward matrix for this temperature, where every calibration has /// one. fn interpolated_forward_at(&self, temperature: f32) -> Option<[[f32; 3]; 3]> { if let [only] = self.calibrations.as_slice() { return only.forward; } let (low, high) = self.bracketing(temperature); let (a, b) = (low.forward?, high.forward?); Some(lerp3(&a, &b, blend(low, high, temperature))) } /// The calibrations this profile was built from, coolest first. pub fn calibrations(&self) -> &[Calibration] { &self.calibrations } /// The as-shot neutral, in camera RGB, where the file carried one. pub fn neutral(&self) -> Option<[f32; 3]> { self.neutral } } /// Where `temperature` sits between two calibrations, in mireds, clamped. /// /// Returns the weight of `high` — 0 at `low`'s illuminant, 1 at `high`'s. /// Clamped because extrapolating a two-point calibration past either end /// produces matrices that are not measurements of anything: candlelight is /// well below illuminant A and open shade well above D65, and both happen /// often enough to matter. fn blend(low: &Calibration, high: &Calibration, temperature: f32) -> f32 { let mired = |k: f32| 1.0e6 / k.max(1.0); let (lo, hi) = (mired(low.temperature), mired(high.temperature)); if (lo - hi).abs() < 1e-6 { // Two calibrations under the same light. Nothing to interpolate; take // the second, which is the one a DNG writer would have meant as the // primary. return 1.0; } ((mired(temperature) - lo) / (hi - lo)).clamp(0.0, 1.0) } /// Linear blend of two matrices. `t` is the weight of `b`. fn lerp3(a: &[[f32; 3]; 3], b: &[[f32; 3]; 3], t: f32) -> [[f32; 3]; 3] { let mut out = [[0.0f32; 3]; 3]; for i in 0..3 { for j in 0..3 { out[i][j] = a[i][j] + (b[i][j] - a[i][j]) * t; } } out } /// The chromaticity of a camera-space neutral, under a given matrix. /// /// `xyz_to_cam` says what the camera would report for a colour; the neutral is /// what it *did* report, so the colour is the inverse applied to it. fn neutral_to_xy(neutral: &[f32; 3], xyz_to_cam: &[[f32; 3]; 3]) -> Option<(f32, f32)> { let cam_to_xyz = invert3(xyz_to_cam)?; let mut xyz = [0.0f32; 3]; for (i, row) in cam_to_xyz.iter().enumerate() { for k in 0..3 { xyz[i] += row[k] * neutral[k]; } } let sum: f32 = xyz.iter().sum(); if !sum.is_finite() || sum.abs() < 1e-9 { return None; } let (x, y) = (xyz[0] / sum, xyz[1] / sum); if !x.is_finite() || !y.is_finite() { return None; } Some((x, y)) } /// TRACES: FR-DEV-3e /// Correlated colour temperature of a chromaticity, by McCamy's approximation. /// /// The exact answer is the nearest point on the Planckian locus measured in a /// uniform chromaticity space, which means either a minimisation or Robertson's /// thirty-one-row isotemperature table. McCamy's cubic reproduces it to within /// a couple of kelvin between 2856 K and 6504 K — which is not a coincidence, /// it is the range every dual-illuminant profile is calibrated across — and it /// is four multiplies. /// /// Accuracy here buys very little in any case. The value feeds an /// interpolation weight that moves smoothly and slowly: on a Canon 6D, being /// 200 K wrong about the scene shifts the resulting matrix by well under a /// thousandth in every element. What would matter is being wrong by /// *thousands*, which is what the clamp below prevents. fn cct_from_xy(x: f32, y: f32) -> f32 { // The epicentre McCamy fits the isotemperature lines around. A neutral // landing exactly on it would divide by zero, which is a chromaticity no // camera reports but not one to leave undefended. let denominator = 0.1858 - y; if denominator.abs() < 1e-6 { return 6504.0; } let n = (x - 0.3320) / denominator; let cct = 449.0 * n * n * n + 3525.0 * n * n + 6823.3 * n + 5520.33; if !cct.is_finite() { return 6504.0; } // The DNG specification's own working range. Outside it the fit is not a // temperature so much as a number, and a profile has nothing calibrated // there to interpolate toward anyway. cct.clamp(1667.0, 25000.0) } /// TRACES: FR-DEV-3e /// The nominal colour temperature of a DNG calibration illuminant. /// /// `None` for the ones that do not have one. `Flash` is a real illuminant with /// a real spectrum and no agreed temperature; `Unknown` is the tag's way of /// saying nothing. Neither can be placed on the interpolation axis, and giving /// them a plausible-looking number would put a calibration somewhere its author /// never claimed it belonged. /// /// The fluorescent entries are the CIE F-series' nominal correlated /// temperatures. They are correlated temperatures of decidedly non-Planckian /// light, so interpolating toward one is approximate in a way daylight is not — /// but a profile calibrated under fluorescent light is describing a sensor /// under fluorescent light, and placing it at roughly the right colour is much /// better than discarding it. fn illuminant_temperature(illuminant: Illuminant) -> Option { Some(match illuminant { // CIE standard illuminant A: a tungsten filament at 2856 K. The low // end of essentially every dual-illuminant profile ever written. Illuminant::A | Illuminant::Tungsten => 2856.0, Illuminant::IsoStudioTungsten => 3200.0, Illuminant::WhiteFluorescent => 3500.0, Illuminant::CoolWhiteFluorescent => 4150.0, Illuminant::B => 4874.0, Illuminant::D50 => 5003.0, Illuminant::DaylightWhiteFluorescent => 5000.0, Illuminant::D55 | Illuminant::Daylight | Illuminant::FineWeather => 5503.0, Illuminant::Fluorescent => 4230.0, Illuminant::DaylightFluorescent => 6430.0, // The sRGB white point, and the high end of essentially every // dual-illuminant profile. Illuminant::D65 => 6504.0, Illuminant::C => 6774.0, Illuminant::CloudyWeather => 6504.0, Illuminant::D75 | Illuminant::Shade => 7504.0, Illuminant::Flash | Illuminant::Unknown => return None, }) } /// Compose a forward matrix into camera RGB → linear sRGB. /// /// `forward` takes white-balanced camera RGB to XYZ under D50, which is the /// DNG connection space and *not* the space sRGB is defined in. So: adapt D50 /// to D65 with Bradford, convert to sRGB primaries, and normalise the rows so /// that camera neutral lands on sRGB white. /// /// The normalisation should be very close to a no-op — a forward matrix's rows /// sum to the D50 white point by definition, which survives the adaptation as /// D65 white and the conversion as (1, 1, 1). It is applied anyway, because /// "should be" is doing load-bearing work in that sentence for a tag written /// by a hundred different converters, and because it also absorbs the green /// normalisation the pipeline's white balance applies and this matrix would /// otherwise have to know about. fn forward_to_srgb(forward: &[[f32; 3]; 3]) -> Option<[f32; 9]> { // Bradford-adapted D50 → D65. The chromatic adaptation everything else in // colour management uses, so that a DarkRoom render and an ICC-aware // application disagree about nothing here. #[allow(clippy::excessive_precision)] const D50_TO_D65: [[f32; 3]; 3] = [ [0.9555766, -0.0230393, 0.0631636], [-0.0282895, 1.0099416, 0.0210077], [0.0122982, -0.0204830, 1.3299098], ]; // XYZ (D65) → linear sRGB. Written out rather than shared with // `cam_to_srgb_from`'s copy for the reason `daylight_wb` states: two // matrices that must agree are better checked against the specification // than against each other. #[allow(clippy::excessive_precision)] const XYZ_TO_SRGB: [[f32; 3]; 3] = [ [3.2404542, -1.5371385, -0.4985314], [-0.9692660, 1.8760108, 0.0415560], [0.0556434, -0.2040259, 1.0572252], ]; if !usable(forward) { return None; } let m = mul3(&XYZ_TO_SRGB, &mul3(&D50_TO_D65, forward)); let mut out = [0.0f32; 9]; for (i, row) in m.iter().enumerate() { let sum: f32 = row.iter().sum(); // A row summing to zero would send neutral to black in that channel. // Something is very wrong with the tag; say so rather than render it. if !sum.is_finite() || sum.abs() < 1e-6 { return None; } for j in 0..3 { out[i * 3 + j] = row[j] / sum; } } if out.iter().any(|v| !v.is_finite()) { return None; } Some(out) } /// Ordinary 3×3 multiply, `a · b`. fn mul3(a: &[[f32; 3]; 3], b: &[[f32; 3]; 3]) -> [[f32; 3]; 3] { let mut out = [[0.0f32; 3]; 3]; for i in 0..3 { for j in 0..3 { for k in 0..3 { out[i][j] += a[i][k] * b[k][j]; } } } out } /// Whether a matrix is a measurement rather than a placeholder. /// /// All-zero is how rawler represents "this body has no matrix", and inverting /// it renders black. A non-finite element is a corrupt tag. Both are absence, /// and absence has to reach the caller as absence. fn usable(m: &[[f32; 3]; 3]) -> bool { m.iter().flatten().any(|v| v.abs() > f32::EPSILON) && m.iter().flatten().all(|v| v.is_finite()) } /// Recover `AsShotNeutral` from rawler's white balance coefficients. /// /// They are its reciprocal — `wb_coeffs[i] == 1.0 / neutral[i]` — for DNG /// files by construction, and for native formats because a maker's white /// balance tag means the same thing. The fourth channel is `NaN` on every /// three-colour sensor and is not consulted. fn neutral_from_coefficients(coefficients: [f32; 4]) -> Option<[f32; 3]> { let usable = |v: f32| v.is_finite() && v > 0.0; if !(usable(coefficients[0]) && usable(coefficients[1]) && usable(coefficients[2])) { return None; } Some([ 1.0 / coefficients[0], 1.0 / coefficients[1], 1.0 / coefficients[2], ]) } /// TRACES: FR-DEV-3e /// Read the DNG tags rawler parses but does not expose. /// /// Returns an empty set for every non-DNG file, and that is not a failure: /// [`rawler::decoders::Decoder::ifd`] is implemented by the DNG decoder alone /// and defaults to `Ok(None)` everywhere else, so a CR2 or an ARW simply has /// no forward matrix to find. Those bodies still get dual-illuminant /// interpolation from rawler's own camera database — the forward matrix is the /// only thing this adds that they cannot have. pub fn read_dng_matrices(decoder: &dyn rawler::decoders::Decoder) -> DngMatrices { use rawler::decoders::WellKnownIFD; use rawler::tags::DngTag; let Ok(Some(ifd)) = decoder.ifd(WellKnownIFD::Root) else { return DngMatrices::default(); }; // Recursive, because a DNG converter is free to put the colour tags in a // sub-IFD beside the raw image rather than in IFD0, and both layouts are // in the wild. let matrix_at = |tag: DngTag| -> Option<[[f32; 3]; 3]> { let entry = ifd.get_entry_recursive(tag)?; if entry.count() < 9 { return None; } let at = |i: usize| entry.force_f32(i); Some([ [at(0), at(1), at(2)], [at(3), at(4), at(5)], [at(6), at(7), at(8)], ]) }; let illuminant_at = |tag: DngTag| -> Illuminant { // 21 is D65, which is the DNG specification's default for a missing // `CalibrationIlluminant` and is what rawler assumes for the colour // matrices. Assuming anything else here would pair a forward matrix // with the wrong colour matrix. ifd.get_entry_recursive(tag) .map(|e| e.force_u16(0)) .unwrap_or(21) .try_into() .unwrap_or(Illuminant::D65) }; let forward = [ matrix_at(DngTag::ForwardMatrix1) .map(|m| (illuminant_at(DngTag::CalibrationIlluminant1), m)), matrix_at(DngTag::ForwardMatrix2) .map(|m| (illuminant_at(DngTag::CalibrationIlluminant2), m)), ]; let as_shot_neutral = ifd .get_entry_recursive(DngTag::AsShotNeutral) .filter(|e| e.count() >= 3) .map(|e| [e.force_f32(0), e.force_f32(1), e.force_f32(2)]) .filter(|n| n.iter().all(|v| v.is_finite() && *v > 0.0)); DngMatrices { forward, as_shot_neutral, } } #[cfg(test)] mod tests { use super::*; /// The Canon EOS 6D's two matrices, as rawler's camera database holds /// them. A real dual-illuminant body, and the one every other measurement /// in this crate was taken on. const SIX_D_A: [[f32; 3]; 3] = [ [0.7546, -0.1435, -0.0929], [-0.3846, 1.1488, 0.2692], [-0.0332, 0.1209, 0.6370], ]; const SIX_D_D65: [[f32; 3]; 3] = [ [0.7034, -0.0804, -0.1014], [-0.4420, 1.2564, 0.2058], [-0.0851, 0.1994, 0.5758], ]; fn six_d(neutral: Option<[f32; 3]>) -> CameraProfile { CameraProfile::new( vec![ Calibration { temperature: 2856.0, xyz_to_cam: SIX_D_A, forward: None, }, Calibration { temperature: 6504.0, xyz_to_cam: SIX_D_D65, forward: None, }, ], neutral, ) .expect("two calibrations") } /// The camera-space neutral of a scene at roughly `kelvin`, derived from /// the body's own matrices so the test is not asserting against a number /// somebody typed. fn neutral_at(profile: &CameraProfile, kelvin: f32) -> [f32; 3] { // The daylight locus is close enough to Planckian over this range for // a test fixture; what matters is that warm scenes and cool scenes // produce distinguishable neutrals, not that they are exact. let m = profile.interpolated_xyz_to_cam_at(kelvin); let xy = planckian_xy(kelvin); let xyz = [xy.0 / xy.1, 1.0, (1.0 - xy.0 - xy.1) / xy.1]; let mut n = [0.0f32; 3]; for (i, row) in m.iter().enumerate() { for k in 0..3 { n[i] += row[k] * xyz[k]; } } n } /// Kim et al.'s cubic fit to the daylight/Planckian locus, good from /// 1667 K to 25000 K. Test scaffolding only — nothing outside this module /// needs to go from a temperature back to a chromaticity. // // Coefficients at their published precision rather than trimmed to what // f32 can hold, for the reason `cam_to_srgb_from` states about its own: // a fit checked against its paper is worth more than one checked against // a linter, and the rounding happens identically either way. #[allow(clippy::excessive_precision)] fn planckian_xy(kelvin: f32) -> (f32, f32) { let t = kelvin.clamp(1667.0, 25000.0); let (t1, t2, t3) = (1.0e3 / t, 1.0e6 / (t * t), 1.0e9 / (t * t * t)); let x = if t <= 4000.0 { -0.2661239 * t3 - 0.2343589 * t2 + 0.8776956 * t1 + 0.179910 } else { -3.0258469 * t3 + 2.1070379 * t2 + 0.2226347 * t1 + 0.240390 }; let y = if t <= 2222.0 { -1.1063814 * x * x * x - 1.34811020 * x * x + 2.18555832 * x - 0.20219683 } else if t <= 4000.0 { -0.9549476 * x * x * x - 1.37418593 * x * x + 2.09137015 * x - 0.16748867 } else { 3.0817580 * x * x * x - 5.87338670 * x * x + 3.75112997 * x - 0.37001483 }; (x, y) } #[test] fn a_single_matrix_body_is_unchanged_by_the_interpolation() { // The regression that matters most: most bodies in rawler's database // carry one matrix, and dual-illuminant support must not perturb what // they rendered as yesterday. let profile = CameraProfile::new( vec![Calibration { temperature: 6504.0, xyz_to_cam: SIX_D_D65, forward: None, }], Some([0.5, 1.0, 0.7]), ) .expect("one calibration"); assert_eq!(profile.xyz_to_cam(), SIX_D_D65); assert_eq!(profile.cam_to_srgb(), cam_to_srgb_from(&SIX_D_D65)); } #[test] fn a_tungsten_scene_lands_nearer_the_tungsten_matrix() { // The whole point of two matrices. If the weighting ran the wrong way // this test is the only thing between that and a warm scene rendered // through a daylight calibration, which is the flat, orange-skinned // look the requirement exists to avoid. let profile = six_d(None); let warm = neutral_at(&profile, 2900.0); let profile = six_d(Some(warm)); let t = profile.scene_temperature(); assert!( (2500.0..3600.0).contains(&t), "a tungsten neutral was read as {t} K" ); let m = profile.xyz_to_cam(); let to_a: f32 = distance(&m, &SIX_D_A); let to_d65: f32 = distance(&m, &SIX_D_D65); assert!( to_a < to_d65, "interpolated matrix sits {to_a} from A and {to_d65} from D65" ); } #[test] fn a_daylight_scene_lands_nearer_the_daylight_matrix() { let profile = six_d(None); let cool = neutral_at(&profile, 6500.0); let profile = six_d(Some(cool)); let t = profile.scene_temperature(); assert!( (5200.0..8000.0).contains(&t), "a daylight neutral was read as {t} K" ); let m = profile.xyz_to_cam(); assert!(distance(&m, &SIX_D_D65) < distance(&m, &SIX_D_A)); } #[test] fn the_interpolated_matrix_still_keeps_neutral_neutral() { // A blend of two valid matrices is not automatically a valid matrix. // This is the property the render depends on: whatever the scene, the // composed transform must send camera neutral to sRGB neutral, or the // frame carries a cast that looks like a bug in the decoder. for kelvin in [2500.0, 3000.0, 4500.0, 6500.0, 9000.0] { let base = six_d(None); let profile = six_d(Some(neutral_at(&base, kelvin))); let m = profile.cam_to_srgb().expect("a real body inverts"); for (i, row) in m.chunks(3).enumerate() { let sum: f32 = row.iter().sum(); assert!( (sum - 1.0).abs() < 1e-4, "at {kelvin} K row {i} sums to {sum}" ); } } } #[test] fn extrapolation_is_clamped_to_what_was_measured() { // Candlelight is far below illuminant A and open shade far above D65. // A linear blend run past either end produces matrices that are not // measurements of anything, and the failure is silent — a plausible // matrix with implausible primaries. let profile = six_d(None); let candle = six_d(Some(neutral_at(&profile, 1700.0))); let shade = six_d(Some(neutral_at(&profile, 20000.0))); assert_eq!(candle.xyz_to_cam(), SIX_D_A, "clamped to the warm end"); assert_eq!(shade.xyz_to_cam(), SIX_D_D65, "clamped to the cool end"); } #[test] fn with_no_white_balance_the_daylight_calibration_is_assumed() { // Not the midpoint: a midpoint is a scene lit by nothing anyone has // photographed. Daylight is the single most likely answer and is the // behaviour every previous version of this decoder had. let profile = six_d(None); assert_eq!(profile.xyz_to_cam(), SIX_D_D65); } #[test] fn a_forward_matrix_is_preferred_and_agrees_about_white() { // The property that lets the two routes coexist: a body gaining a // forward matrix must change its colour rendering without changing its // exposure or its neutral. Rows summing to one is what says so. // // The matrix is the D50-adapted forward matrix of a body whose // primaries are close to sRGB's, which is enough to exercise the path // — the arithmetic under test is the adaptation and the normalisation, // not anyone's measurement. let forward = [ [0.7500, 0.1800, 0.0347], [0.2900, 0.6800, 0.0300], [0.0200, 0.0900, 0.7147], ]; let profile = CameraProfile::new( vec![Calibration { temperature: 6504.0, xyz_to_cam: SIX_D_D65, forward: Some(forward), }], None, ) .expect("one calibration"); let m = profile.cam_to_srgb().expect("a forward matrix composes"); for (i, row) in m.chunks(3).enumerate() { let sum: f32 = row.iter().sum(); assert!((sum - 1.0).abs() < 1e-4, "row {i} sums to {sum}"); } assert_ne!( Some(m), cam_to_srgb_from(&SIX_D_D65), "the forward matrix must actually be used, not merely accepted" ); } #[test] fn a_forward_matrix_at_only_one_end_is_not_mixed_with_an_inverted_one() { // Interpolating a forward matrix against an inverted colour matrix // would blend two different statements of the same relationship and // land between them — worse than either, and invisible. let profile = CameraProfile::new( vec![ Calibration { temperature: 2856.0, xyz_to_cam: SIX_D_A, forward: Some([[0.75, 0.18, 0.03], [0.29, 0.68, 0.03], [0.02, 0.09, 0.71]]), }, Calibration { temperature: 6504.0, xyz_to_cam: SIX_D_D65, forward: None, }, ], None, ) .expect("two calibrations"); assert_eq!( profile.cam_to_srgb(), cam_to_srgb_from(&SIX_D_D65), "a half-populated forward pair must fall back to the colour matrices" ); } #[test] fn a_third_calibration_brackets_rather_than_widening_the_blend() { // Not reachable through rawler today, which reads `ColorMatrix1/2` and // stops. Asserted anyway because the failure it prevents is silent: a // profile with tungsten, daylight and shade calibrations interpolated // across the *extremes* would ignore the daylight matrix entirely and // render every ordinary photograph through a blend of candlelight and // open shade. let shade: [[f32; 3]; 3] = [ [0.68, -0.07, -0.10], [-0.46, 1.28, 0.20], [-0.09, 0.21, 0.56], ]; let profile = CameraProfile::new( vec![ Calibration { temperature: 2856.0, xyz_to_cam: SIX_D_A, forward: None, }, Calibration { temperature: 6504.0, xyz_to_cam: SIX_D_D65, forward: None, }, Calibration { temperature: 7504.0, xyz_to_cam: shade, forward: None, }, ], None, ) .expect("three calibrations"); // 5000 K falls between A and D65, so the shade matrix contributes // nothing at all. let m = profile.interpolated_xyz_to_cam_at(5000.0); let between = |a: f32, b: f32, v: f32| (v - a) * (v - b) <= 1e-6; assert!( between(SIX_D_A[0][0], SIX_D_D65[0][0], m[0][0]), "5000 K landed at {}, outside the A..D65 span", m[0][0] ); // And 7000 K falls in the upper pair, where A contributes nothing. let m = profile.interpolated_xyz_to_cam_at(7000.0); assert!(between(SIX_D_D65[0][0], shade[0][0], m[0][0])); } #[test] fn an_illuminant_with_no_temperature_is_dropped_rather_than_guessed() { // `Flash` has a real spectrum and no agreed temperature, `Unknown` is // the tag saying nothing. Placing either on the interpolation axis // would drag the result toward a point its author never claimed. assert_eq!(illuminant_temperature(Illuminant::Flash), None); assert_eq!(illuminant_temperature(Illuminant::Unknown), None); assert_eq!(illuminant_temperature(Illuminant::A), Some(2856.0)); assert_eq!(illuminant_temperature(Illuminant::D65), Some(6504.0)); } #[test] fn the_temperature_estimate_is_stable_under_iteration() { // The fixed point is only worth three rounds if three rounds is where // it has stopped moving. If this failed, the render would depend on an // iteration count nobody chose deliberately. let base = six_d(None); for kelvin in [3000.0, 4500.0, 6500.0] { let profile = six_d(Some(neutral_at(&base, kelvin))); let once = profile.scene_temperature(); // Re-seeding from the answer must not move it. let again = { let mut t = once; for _ in 0..3 { let m = profile.interpolated_xyz_to_cam_at(t); let (x, y) = neutral_to_xy(&profile.neutral.unwrap(), &m).unwrap(); t = cct_from_xy(x, y); } t }; assert!( (once - again).abs() < 25.0, "at {kelvin} K the estimate moved from {once} to {again}" ); } } #[test] fn mccamy_recovers_the_temperatures_the_profiles_are_calibrated_at() { // The only range the approximation has to be right over: the two ends // of a dual-illuminant profile, and the daylight between them. for kelvin in [2856.0f32, 3200.0, 4000.0, 5003.0, 5503.0, 6504.0] { let (x, y) = planckian_xy(kelvin); let got = cct_from_xy(x, y); assert!( (got - kelvin).abs() / kelvin < 0.03, "{kelvin} K read back as {got} K" ); } } #[test] fn an_all_zero_matrix_never_becomes_a_calibration() { // rawler's placeholder for "this body has no matrix". Accepting it // would render black, which is much harder to diagnose than // uncalibrated. assert!(!usable(&[[0.0; 3]; 3])); assert!(!usable(&[ [f32::NAN, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0] ])); assert!(usable(&SIX_D_D65)); } #[test] fn as_shot_neutral_is_the_reciprocal_of_the_coefficients() { // Confusing the two inverts the estimate: a tungsten scene would be // read as daylight and rendered through the wrong end of the profile. let n = neutral_from_coefficients([2.0, 1.0, 4.0, f32::NAN]).expect("three usable"); assert!((n[0] - 0.5).abs() < 1e-6); assert!((n[1] - 1.0).abs() < 1e-6); assert!((n[2] - 0.25).abs() < 1e-6); } #[test] fn absent_coefficients_yield_no_neutral_rather_than_a_neutral_one() { // (1, 1, 1) is a claim about the scene — that it was lit by something // this sensor happens to see as grey — and it is one no file made. assert_eq!(neutral_from_coefficients([0.0, 0.0, 0.0, 0.0]), None); assert_eq!(neutral_from_coefficients([f32::NAN, 1.0, 1.0, 1.0]), None); } /// Sum of absolute differences between two matrices. fn distance(a: &[[f32; 3]; 3], b: &[[f32; 3]; 3]) -> f32 { a.iter() .flatten() .zip(b.iter().flatten()) .map(|(x, y)| (x - y).abs()) .sum() } }