Revert "Merge PR #9455: raw: develop the default rendering band by band, without the f32 RGB frame"

This reverts commit ddb92be1fa, reversing
changes made to f12be92413.
This commit is contained in:
David Mann
2026-10-25 18:38:22 +01:00
parent 92f592684c
commit 0134d2dc70
3 changed files with 216 additions and 852 deletions
+51 -411
View File
@@ -101,11 +101,6 @@ impl Padded {
Padded { data: vec![0.0; stride * (h + 2 * PAD)], w, h, stride }
}
/// Whether the plane is non-empty and holds all its rows.
fn fits(&self) -> bool {
self.w != 0 && self.h != 0 && self.data.len() >= self.stride * (self.h + 2 * PAD)
}
/// From an unpadded plane (tests and small inputs).
#[cfg(test)]
pub fn from_plane(plane: &[f32], w: usize, h: usize) -> Self {
@@ -195,145 +190,53 @@ fn kernels(phase: [u8; 4], method: Demosaic, stride: usize) -> [[Vec<(isize, f32
pub(crate) fn demosaic(p: &Padded, phase: [u8; 4], method: Demosaic, to_xyz: &Mat3, ahd_jobs: usize) -> Vec<f32> {
let (w, h, s) = (p.w, p.h, p.stride);
let mut out = vec![0.0f32; w * h * 3];
if !p.fits() {
if w == 0 || h == 0 || p.data.len() < s * (h + 2 * PAD) {
return out;
}
if method == Demosaic::Ahd {
let (d, m, labf) = (&p.data, ahd_matrix(to_xyz), LabF::new());
for_bands(
&mut out,
AHD_BAND * w * 3,
ahd_bands(h),
ahd_jobs,
|| AhdScratch::new((AHD_BAND + 2 * M) * s, s),
|band, chunk, sc| {
ahd_band(d, w, s, phase, &m, &labf, band, chunk, sc);
},
);
ahd(p, phase, to_xyz, &mut out, ahd_jobs);
return out;
}
let ks = kernels(phase, method, s);
let rows = par::band_rows(w);
par::chunks_mut(&mut out, rows * w * 3, |band, chunk| linear_band(&p.data, w, s, phase, method, &ks, band * rows, 0, chunk));
let d = &p.data;
par::chunks_mut(&mut out, rows * w * 3, |band, chunk| {
let y0 = band * rows;
for (r, row) in chunk.chunks_exact_mut(w * 3).enumerate() {
let y = y0 + r;
let base = (y + PAD) * s + PAD;
for (x, o) in row.as_chunks_mut::<3>().0.iter_mut().enumerate() {
let pos = (y & 1) * 2 + (x & 1);
let i = base + x;
if method == Demosaic::Bilinear {
// Direct Bayer bilinear (the kernels' result without the tap loops).
let own = usize::from(phase[pos]).min(2);
let (n, so, e, wv) = (d[i - s], d[i + s], d[i + 1], d[i - 1]);
if own == 1 {
let hc = usize::from(phase[(y & 1) * 2 + ((x + 1) & 1)]).min(2);
o[1] = d[i];
o[hc] = (e + wv) * 0.5;
o[2 - hc] = (n + so) * 0.5;
} else {
o[own] = d[i];
o[1] = (n + so + e + wv) * 0.25;
o[2 - own] = (d[i - s - 1] + d[i - s + 1] + d[i + s - 1] + d[i + s + 1]) * 0.25;
}
continue;
}
for (c, k) in ks[pos].iter().enumerate() {
let mut v = 0.0f32;
for &(off, wt) in k {
v += wt * d[(i as isize + off) as usize];
}
o[c] = v.clamp(0.0, 1.0);
}
}
}
});
out
}
/// [`demosaic`] without the full RGB frame and without the padded plane: every band loads the
/// rows it reads (`fill(y, row)` writes the normalized CFA row `y` of the image, `w` samples;
/// the mirrored borders are added here), demosaics into a band buffer, and `emit(rgb, out)`
/// turns it into the band's `out` samples (the same rows, `w * 3` per row). Each pixel depends
/// only on the plane around it, never on how rows are split or ordered, so the result is that
/// of `emit` over [`demosaic`] of the filled plane.
#[allow(clippy::too_many_arguments)]
pub(crate) fn demosaic_fused(
w: usize,
h: usize,
phase: [u8; 4],
method: Demosaic,
to_xyz: &Mat3,
ahd_jobs: usize,
fill: &(impl Fn(usize, &mut [f32]) + Sync),
out: &mut [u16],
emit: &(impl Fn(&[f32], &mut [u16]) + Sync),
) {
if w == 0 || h == 0 || out.len() < w * h * 3 {
return;
}
let s = w + 2 * PAD;
let out = &mut out[..w * h * 3];
// Padded rows of the output rows `y0..y0 + rows` (`rows + 2 * PAD` of them, row `PAD` is `y0`)
// with the mirror reflection `Padded::fill_borders` applies.
let load = |win: &mut [f32], y0: usize, rows: usize| {
for (j, prow) in win.chunks_exact_mut(s).take(rows + 2 * PAD).enumerate() {
fill(reflect((y0 + j) as isize - PAD as isize, h), &mut prow[PAD..PAD + w]);
for px in (0..PAD).chain(PAD + w..s) {
prow[px] = prow[PAD + reflect(px as isize - PAD as isize, w)];
}
}
};
if method == Demosaic::Ahd {
let (m, labf) = (ahd_matrix(to_xyz), LabF::new());
let init = || (AhdScratch::new((AHD_BAND + 2 * M) * s, s), vec![0.0f32; AHD_BAND * w * 3], vec![0.0f32; (AHD_BAND + 2 * PAD) * s]);
for_bands(out, AHD_BAND * w * 3, ahd_bands(h), ahd_jobs, init, |band, chunk, (sc, buf, win)| {
load(win, band * AHD_BAND, chunk.len() / (w * 3));
let buf = &mut buf[..chunk.len()];
// The window starts at a multiple of AHD_BAND (even), so it is band 0 of its own rows.
ahd_band(win, w, s, phase, &m, &labf, 0, buf, sc);
emit(buf, chunk);
});
return;
}
let ks = kernels(phase, method, s);
let init = || (vec![0.0f32; FUSED_ROWS * w * 3], vec![0.0f32; (FUSED_ROWS + 2 * PAD) * s]);
for_bands(out, FUSED_ROWS * w * 3, h.div_ceil(FUSED_ROWS), usize::MAX, init, |band, chunk, (buf, win)| {
load(win, band * FUSED_ROWS, chunk.len() / (w * 3));
let buf = &mut buf[..chunk.len()];
linear_band(win, w, s, phase, method, &ks, band * FUSED_ROWS, band * FUSED_ROWS, buf);
emit(buf, chunk);
});
}
/// Output rows per band of the fused bilinear / MHC path (the window adds `2 * PAD` rows).
pub(crate) const FUSED_ROWS: usize = 32;
/// Rows `y0..` of bilinear or MHC output (`chunk.len() / (w * 3)` of them); `d` is the padded
/// plane from image row `dy` on (`dy` 0: all of it).
#[allow(clippy::too_many_arguments)]
fn linear_band(d: &[f32], w: usize, s: usize, phase: [u8; 4], method: Demosaic, ks: &[[Vec<(isize, f32)>; 3]; 4], y0: usize, dy: usize, chunk: &mut [f32]) {
for (r, row) in chunk.chunks_exact_mut(w * 3).enumerate() {
let y = y0 + r;
let base = (y - dy + PAD) * s + PAD;
for (x, o) in row.as_chunks_mut::<3>().0.iter_mut().enumerate() {
let pos = (y & 1) * 2 + (x & 1);
let i = base + x;
if method == Demosaic::Bilinear {
// Direct Bayer bilinear (the kernels' result without the tap loops).
let own = usize::from(phase[pos]).min(2);
let (n, so, e, wv) = (d[i - s], d[i + s], d[i + 1], d[i - 1]);
if own == 1 {
let hc = usize::from(phase[(y & 1) * 2 + ((x + 1) & 1)]).min(2);
o[1] = d[i];
o[hc] = (e + wv) * 0.5;
o[2 - hc] = (n + so) * 0.5;
} else {
o[own] = d[i];
o[1] = (n + so + e + wv) * 0.25;
o[2 - own] = (d[i - s - 1] + d[i - s + 1] + d[i + s - 1] + d[i + s + 1]) * 0.25;
}
continue;
}
for (c, k) in ks[pos].iter().enumerate() {
let mut v = 0.0f32;
for &(off, wt) in k {
v += wt * d[(i as isize + off) as usize];
}
o[c] = v.clamp(0.0, 1.0);
}
}
}
}
/// Runs `f(band, chunk, scratch)` over the `unit`-sized chunks of `out` (`bands` of them), with
/// at most `jobs` scratch values alive: a few groups of consecutive bands per worker (for
/// balance), each group building its scratch once.
pub(crate) fn for_bands<T: Send, S>(
out: &mut [T],
unit: usize,
bands: usize,
jobs: usize,
init: impl Fn() -> S + Sync,
f: impl Fn(usize, &mut [T], &mut S) + Sync,
) {
let groups = (par::workers().min(jobs.max(1)) * 4).clamp(1, bands.max(1));
let per_group = bands.div_ceil(groups).max(1);
par::chunks_mut_limited(out, per_group * unit, jobs, |group, gchunk| {
let mut sc = init();
for (j, chunk) in gchunk.chunks_mut(unit).enumerate() {
f(group * per_group + j, chunk, &mut sc);
}
});
}
/// A CFA repeat as seen from the plane's origin: colour of plane pixel
/// (x, y) is `colors[(y % h) * w + x % w]`.
#[derive(Debug, Clone)]
@@ -415,36 +318,10 @@ fn phase_taps(pat: &Pattern, px: isize, py: isize) -> Taps {
/// without the edge-directed weights it is no faster and leaves zipper
/// artefacts along every edge of an X-Trans image.
pub(crate) fn demosaic_cfa(p: &Padded, pat: &Pattern) -> Vec<f32> {
let (w, h) = (p.w, p.h);
let mut out = vec![0.0f32; w * h * 3];
let band = par::band_rows(w);
cfa_with(p, pat, |rows| par::chunks_mut(&mut out, band * w * 3, |b, chunk| rows(b * band, chunk)));
out
}
/// [`demosaic_cfa`] without the full RGB frame (the green plane stays): see [`demosaic_fused`].
pub(crate) fn demosaic_cfa_fused(p: &Padded, pat: &Pattern, out: &mut [u16], emit: &(impl Fn(&[f32], &mut [u16]) + Sync)) {
let (w, h) = (p.w, p.h);
if out.len() < w * h * 3 {
return;
}
let band = par::band_rows(w);
let init = || vec![0.0f32; band * w * 3];
cfa_with(p, pat, |rows| {
for_bands(&mut out[..w * h * 3], band * w * 3, h.div_ceil(band), usize::MAX, init, |b, chunk, buf| {
let buf = &mut buf[..chunk.len()];
rows(b * band, buf);
emit(buf, chunk);
});
});
}
/// The generic demosaic's green plane, then `run` with a function that fills the RGB rows
/// `y0..` of a chunk (`chunk.len() / (w * 3)` of them). Not called for an unusable plane or pattern.
fn cfa_with(p: &Padded, pat: &Pattern, run: impl FnOnce(&(dyn Fn(usize, &mut [f32]) + Sync))) {
let (w, h, s) = (p.w, p.h, p.stride);
if !p.fits() || pat.w == 0 || pat.h == 0 || pat.colors.len() < pat.w * pat.h {
return;
let mut out = vec![0.0f32; w * h * 3];
if w == 0 || h == 0 || pat.w == 0 || pat.h == 0 || pat.colors.len() < pat.w * pat.h || p.data.len() < s * (h + 2 * PAD) {
return out;
}
let taps: Vec<Taps> = (0..pat.w * pat.h).map(|i| phase_taps(pat, (i % pat.w) as isize, (i / pat.w) as isize)).collect();
let tap_at = |x: usize, y: usize| &taps[(y % pat.h) * pat.w + x % pat.w];
@@ -512,9 +389,9 @@ fn cfa_with(p: &Padded, pat: &Pattern, run: impl FnOnce(&(dyn Fn(usize, &mut [f3
// 2. Red and blue from colour differences.
let gref = &green;
let gat = |x: usize, y: usize, o: Off| get(gref, w, 0, x, y, o);
run(&|y0: usize, chunk: &mut [f32]| {
par::chunks_mut(&mut out, band * w * 3, |b, chunk| {
for (r, row) in chunk.chunks_exact_mut(w * 3).enumerate() {
let y = y0 + r;
let y = b * band + r;
for (x, o) in row.as_chunks_mut::<3>().0.iter_mut().enumerate() {
let t = tap_at(x, y);
let g = gref.get(y * w + x).copied().unwrap_or(0.0);
@@ -541,6 +418,7 @@ fn cfa_with(p: &Padded, pat: &Pattern, run: impl FnOnce(&(dyn Fn(usize, &mut [f3
}
}
});
out
}
/// Regularizes the inverse-gradient weights (normalized units, about 1/256 of full scale).
@@ -548,10 +426,9 @@ const GRAD_EPS: f32 = 1.0 / 256.0;
/// Regularizes the green-similarity weights of the colour differences.
const DIFF_EPS: f32 = 1.0 / 64.0;
/// CIE f(t) for L*a*b*, tabulated on 0..2 with linear interpolation. Each entry holds the
/// sample and the step to the next one (one load per lookup; the same `a + (b - a) * f`).
/// CIE f(t) for L*a*b*, tabulated on 0..2 with linear interpolation.
struct LabF {
pairs: Vec<[f32; 2]>,
table: Vec<f32>,
}
const LABF_STEPS: usize = 1 << 14;
@@ -559,8 +436,7 @@ const LABF_RANGE: f32 = 2.0;
impl LabF {
fn new() -> Self {
let table: Vec<f32> = (0..=LABF_STEPS + 1).map(|i| Self::exact(i as f32 * LABF_RANGE / LABF_STEPS as f32)).collect();
LabF { pairs: table.windows(2).map(|p| [p[0], p[1] - p[0]]).collect() }
LabF { table: (0..=LABF_STEPS + 1).map(|i| Self::exact(i as f32 * LABF_RANGE / LABF_STEPS as f32)).collect() }
}
fn exact(t: f32) -> f32 {
@@ -574,13 +450,13 @@ impl LabF {
}
let x = t * (LABF_STEPS as f32 / LABF_RANGE);
let i = (x as usize).min(LABF_STEPS);
let [a, step] = self.pairs[i];
a + step * (x - i as f32)
let a = self.table[i];
a + (self.table[i + 1] - a) * (x - i as f32)
}
}
/// Output rows per AHD band.
pub(crate) const AHD_BAND: usize = 32;
const AHD_BAND: usize = 32;
/// Extra band rows above and below the output rows.
const M: usize = 3;
@@ -594,7 +470,7 @@ pub(crate) fn ahd_bands(h: usize) -> usize {
h.div_ceil(AHD_BAND)
}
pub(crate) const AHD_TABLE_BYTES: u64 = (LABF_STEPS as u64 + 1) * 8;
pub(crate) const AHD_TABLE_BYTES: u64 = (LABF_STEPS as u64 + 2) * 4;
/// `v` limited to `lo..=hi` like `f32::clamp`, without its panic on NaN bounds (hostile black
/// levels can put NaN in the plane): same result for every non-NaN bound.
@@ -604,193 +480,7 @@ fn clampf(v: f32, lo: f32, hi: f32) -> f32 {
if v > hi { hi } else { v }
}
/// One AHD band's working planes, reused for every band a worker runs (sized for a full band).
struct AhdScratch {
/// Green interpolated horizontally (0) and vertically (1).
g: [Vec<f32>; 2],
/// Per direction, the RGB candidate (planar).
rgb: [[Vec<f32>; 3]; 2],
/// Per direction, its CIELAB (planar).
lab: [[Vec<f32>; 3]; 2],
hom: [Vec<u8>; 2],
/// Vertical 3-row sums of `hom`, one output row.
vs: [Vec<u32>; 2],
}
impl AhdScratch {
fn new(len: usize, s: usize) -> Self {
let f = || vec![0.0f32; len];
AhdScratch {
g: [f(), f()],
rgb: [[f(), f(), f()], [f(), f(), f()]],
lab: [[f(), f(), f()], [f(), f(), f()]],
hom: [vec![0u8; len], vec![0u8; len]],
vs: [vec![0u32; s], vec![0u32; s]],
}
}
}
/// The camera → XYZ matrix with the rows divided by the white point.
fn ahd_matrix(to_xyz: &Mat3) -> [[f32; 3]; 3] {
let white = [0.9642f32, 1.0, 0.8249];
std::array::from_fn(|i| std::array::from_fn(|j| to_xyz[i][j] as f32 / white[i]))
}
/// One AHD band (see [`ahd_ref`] for the plain statement): the same arithmetic in the same
/// order, per row on slices, planar scratch reused across a worker's bands. Every scratch cell
/// a stage reads was written earlier in the same band, so the planes are never cleared between
/// bands.
#[allow(clippy::too_many_arguments)]
fn ahd_band(d: &[f32], w: usize, s: usize, phase: [u8; 4], m: &[[f32; 3]; 3], labf: &LabF, band: usize, chunk: &mut [f32], sc: &mut AhdScratch) {
let rows = chunk.len() / (w * 3);
let n = rows + 2 * M;
let py = |b: usize| band * AHD_BAND + PAD - M + b;
let drow = |y: usize| &d[y * s..(y + 1) * s];
// 1. Green interpolated horizontally (0) and vertically (1), columns 2..s-2.
{
let [g0, g1] = &mut sc.g;
for b in 0..n {
let y = py(b);
let (rm2, rm1, r0, rp1, rp2) = (drow(y - 2), drow(y - 1), drow(y), drow(y + 1), drow(y + 2));
let gh_row = &mut g0[b * s..(b + 1) * s];
let gv_row = &mut g1[b * s..(b + 1) * s];
let green = [phase[(y & 1) * 2] == 1, phase[(y & 1) * 2 + 1] == 1];
for x in 2..s - 2 {
let c = r0[x];
if green[x & 1] {
gh_row[x] = c;
gv_row[x] = c;
continue;
}
let (l, r) = (r0[x - 1], r0[x + 1]);
let gh = (l + r) * 0.5 + (2.0 * c - r0[x - 2] - r0[x + 2]) * 0.25;
gh_row[x] = clampf(gh, l.min(r), l.max(r));
let (u, dn) = (rm1[x], rp1[x]);
let gv = (u + dn) * 0.5 + (2.0 * c - rm2[x] - rp2[x]) * 0.25;
gv_row[x] = clampf(gv, u.min(dn), u.max(dn));
}
}
}
// 2. RGB candidates and their Lab, band rows 1..n-1, columns 3..s-3.
for dir in 0..2 {
let gd = &sc.g[dir];
let [cr, cg, cb] = &mut sc.rgb[dir];
let [ll, la, lb] = &mut sc.lab[dir];
for b in 1..n - 1 {
let y = py(b);
let (dm1, d0, dp1) = (drow(y - 1), drow(y), drow(y + 1));
let (gm1, g0, gp1) = (&gd[(b - 1) * s..b * s], &gd[b * s..(b + 1) * s], &gd[(b + 1) * s..(b + 2) * s]);
let at = b * s;
let (orr, org, orb) = (&mut cr[at..at + s], &mut cg[at..at + s], &mut cb[at..at + s]);
let (ol, oa, ob) = (&mut ll[at..at + s], &mut la[at..at + s], &mut lb[at..at + s]);
let pr = (y & 1) * 2;
for x in 3..s - 3 {
let own = phase[pr + (x & 1)];
let gv = g0[x];
// `va` goes to channel `ch`, `vb` to `2 - ch` (as `px[ch] = va; px[2 - ch] = vb`,
// built from values: an indexed store into `px` stalls store forwarding).
let (ch, va, vb) = if own == 1 {
(
phase[pr + ((x + 1) & 1)].min(2),
gv + 0.5 * ((d0[x - 1] - g0[x - 1]) + (d0[x + 1] - g0[x + 1])),
gv + 0.5 * ((dm1[x] - gm1[x]) + (dp1[x] - gp1[x])),
)
} else {
(
own.min(2),
d0[x],
gv + 0.25 * ((dm1[x - 1] - gm1[x - 1]) + (dm1[x + 1] - gm1[x + 1]) + (dp1[x - 1] - gp1[x - 1]) + (dp1[x + 1] - gp1[x + 1])),
)
};
let px = match ch {
0 => [va, gv, vb],
2 => [vb, gv, va],
_ => [gv, vb, gv],
};
let px = px.map(|v| v.clamp(0.0, 1.0));
orr[x] = px[0];
org[x] = px[1];
orb[x] = px[2];
let xyz: [f32; 3] = std::array::from_fn(|k| m[k][0] * px[0] + m[k][1] * px[1] + m[k][2] * px[2]);
let (fx, fy, fz) = (labf.f(xyz[0]), labf.f(xyz[1]), labf.f(xyz[2]));
ol[x] = 116.0 * fy - 16.0;
oa[x] = 500.0 * (fx - fy);
ob[x] = 200.0 * (fy - fz);
}
}
}
// 3. Homogeneity maps, band rows 2..n-2, columns 4..s-4.
{
let [l0, a0, b0] = &sc.lab[0];
let [l1, a1, b1] = &sc.lab[1];
let [h0, h1] = &mut sc.hom;
for b in 2..n - 2 {
fn row(pl: &[f32], bb: usize, s: usize) -> &[f32] {
&pl[bb * s..(bb + 1) * s]
}
let (l0u, l0c, l0d) = (row(l0, b - 1, s), row(l0, b, s), row(l0, b + 1, s));
let (a0u, a0c, a0d) = (row(a0, b - 1, s), row(a0, b, s), row(a0, b + 1, s));
let (b0u, b0c, b0d) = (row(b0, b - 1, s), row(b0, b, s), row(b0, b + 1, s));
let (l1u, l1c, l1d) = (row(l1, b - 1, s), row(l1, b, s), row(l1, b + 1, s));
let (a1u, a1c, a1d) = (row(a1, b - 1, s), row(a1, b, s), row(a1, b + 1, s));
let (b1u, b1c, b1d) = (row(b1, b - 1, s), row(b1, b, s), row(b1, b + 1, s));
let (h0r, h1r) = (&mut h0[b * s..(b + 1) * s], &mut h1[b * s..(b + 1) * s]);
for x in 4..s - 4 {
// Neighbours left, right, up, down; (|ΔL|, Δa² + Δb²) from the centre.
let dist = |lc: f32, ac: f32, bc: f32, ln: f32, an: f32, bn: f32| ((lc - ln).abs(), (ac - an) * (ac - an) + (bc - bn) * (bc - bn));
let (lc, ac, bc) = (l0c[x], a0c[x], b0c[x]);
let dh = [
dist(lc, ac, bc, l0c[x - 1], a0c[x - 1], b0c[x - 1]),
dist(lc, ac, bc, l0c[x + 1], a0c[x + 1], b0c[x + 1]),
dist(lc, ac, bc, l0u[x], a0u[x], b0u[x]),
dist(lc, ac, bc, l0d[x], a0d[x], b0d[x]),
];
let (lc, ac, bc) = (l1c[x], a1c[x], b1c[x]);
let dv = [
dist(lc, ac, bc, l1c[x - 1], a1c[x - 1], b1c[x - 1]),
dist(lc, ac, bc, l1c[x + 1], a1c[x + 1], b1c[x + 1]),
dist(lc, ac, bc, l1u[x], a1u[x], b1u[x]),
dist(lc, ac, bc, l1d[x], a1d[x], b1d[x]),
];
let eps_l = dh[0].0.max(dh[1].0).min(dv[2].0.max(dv[3].0));
let eps_c = dh[0].1.max(dh[1].1).min(dv[2].1.max(dv[3].1));
let count = |v: &[(f32, f32); 4]| v.iter().map(|(l, c)| u8::from(*l <= eps_l && *c <= eps_c)).sum::<u8>();
h0r[x] = count(&dh);
h1r[x] = count(&dv);
}
}
}
// 4. Per pixel, the candidate with the larger 3×3 homogeneity sum.
let [vs0, vs1] = &mut sc.vs;
for r in 0..rows {
let b = r + M;
for (hm, vs) in [(&sc.hom[0], &mut *vs0), (&sc.hom[1], &mut *vs1)] {
let (u, c, dn) = (&hm[(b - 1) * s..b * s], &hm[b * s..(b + 1) * s], &hm[(b + 1) * s..(b + 2) * s]);
for x in PAD - 1..PAD + w + 1 {
vs[x] = u32::from(u[x]) + u32::from(c[x]) + u32::from(dn[x]);
}
}
let at = b * s;
let (r0, g0, b0) = (&sc.rgb[0][0][at..at + s], &sc.rgb[0][1][at..at + s], &sc.rgb[0][2][at..at + s]);
let (r1, g1, b1) = (&sc.rgb[1][0][at..at + s], &sc.rgb[1][1][at..at + s], &sc.rgb[1][2][at..at + s]);
let orow = &mut chunk[r * w * 3..(r + 1) * w * 3];
for (x, o) in orow.as_chunks_mut::<3>().0.iter_mut().enumerate() {
let px = x + PAD;
let sh = vs0[px - 1] + vs0[px] + vs0[px + 1];
let sv = vs1[px - 1] + vs1[px] + vs1[px + 1];
let (c0, c1) = ([r0[px], g0[px], b0[px]], [r1[px], g1[px], b1[px]]);
*o = match sh.cmp(&sv) {
std::cmp::Ordering::Greater => c0,
std::cmp::Ordering::Less => c1,
std::cmp::Ordering::Equal => [0.5 * (c0[0] + c1[0]), 0.5 * (c0[1] + c1[1]), 0.5 * (c0[2] + c1[2])],
};
}
}
}
/// The plain AHD statement the optimized [`ahd`] must match bit for bit.
#[cfg(test)]
fn ahd_ref(p: &Padded, phase: [u8; 4], to_xyz: &Mat3, out: &mut [f32], jobs: usize) {
fn ahd(p: &Padded, phase: [u8; 4], to_xyz: &Mat3, out: &mut [f32], jobs: usize) {
let (w, s) = (p.w, p.stride);
let d = &p.data;
// Colour at padded coordinates (PAD is even, so the phase is unchanged).
@@ -924,56 +614,6 @@ mod tests {
assert_eq!(demosaic(&p, [0, 1, 1, 2], Demosaic::Ahd, &IDENTITY, 1).len(), 192);
}
#[test]
fn ahd_matches_reference_bit_for_bit() {
let m = [[0.41, 0.36, 0.18], [0.21, 0.72, 0.07], [0.02, 0.12, 0.95]];
for (w, h) in [(1, 1), (2, 3), (7, 5), (53, 103), (64, 70), (33, 31)] {
let cfa: Vec<f32> = (0..w * h).map(|i| ((i * 7919 + i / w * 104729) % 1013) as f32 / 1000.0 - 0.006).collect();
let p = Padded::from_plane(&cfa, w, h);
for phase in [[0, 1, 1, 2], [2, 1, 1, 0], [1, 0, 2, 1], [1, 2, 0, 1]] {
for jobs in [1, 3, 64] {
let mut want = vec![0.0f32; w * h * 3];
ahd_ref(&p, phase, &m, &mut want, jobs);
let got = demosaic(&p, phase, Demosaic::Ahd, &m, jobs);
assert_eq!(
want.iter().map(|v| v.to_bits()).collect::<Vec<_>>(),
got.iter().map(|v| v.to_bits()).collect::<Vec<_>>(),
"{w}x{h} {phase:?} {jobs}"
);
}
}
}
}
#[test]
fn fused_demosaic_matches_the_full_frame_demosaic() {
let m = [[0.41, 0.36, 0.18], [0.21, 0.72, 0.07], [0.02, 0.12, 0.95]];
// Sizes around the bands (AHD 32 rows, bilinear / MHC 32 rows) and the padding.
for (w, h) in [(1, 1), (2, 3), (7, 5), (53, 103), (64, 64), (33, 31), (5, 70), (31, 65)] {
let cfa: Vec<f32> = (0..w * h).map(|i| ((i * 7919 + i / w * 104729) % 1013) as f32 / 1000.0 - 0.006).collect();
let p = Padded::from_plane(&cfa, w, h);
let fill = |y: usize, row: &mut [f32]| row.copy_from_slice(&cfa[y * w..(y + 1) * w]);
// Any per-sample map stands in for the encode.
let emit = |rgb: &[f32], out: &mut [u16]| {
for (o, v) in out.iter_mut().zip(rgb) {
*o = (v.clamp(0.0, 1.0) * 65535.0 + 0.5) as u16;
}
};
for method in [Demosaic::Bilinear, Demosaic::Mhc, Demosaic::Ahd] {
for phase in [[0, 1, 1, 2], [2, 1, 1, 0], [1, 0, 2, 1], [1, 2, 0, 1]] {
for jobs in [1, 3, 64] {
let full = demosaic(&p, phase, method, &m, jobs);
let mut want = vec![0u16; w * h * 3];
emit(&full, &mut want);
let mut got = vec![0u16; w * h * 3];
demosaic_fused(w, h, phase, method, &m, jobs, &fill, &mut got, &emit);
assert_eq!(want, got, "{w}x{h} {method:?} {phase:?} {jobs}");
}
}
}
}
}
/// Mosaics an RGB image with the given phase.
fn mosaic(rgb: &[[f32; 3]], w: usize, h: usize, phase: [u8; 4]) -> Vec<f32> {
(0..w * h).map(|i| rgb[i][phase[((i / w) & 1) * 2 + ((i % w) & 1)] as usize]).collect()
+159 -423
View File
@@ -5,7 +5,7 @@
//! orientation.
use crate::color::{self, Mat3};
use crate::demosaic::{Demosaic, PAD, Padded, Pattern, demosaic, demosaic_cfa, demosaic_cfa_fused, demosaic_fused, for_bands};
use crate::demosaic::{Demosaic, PAD, Padded, Pattern, demosaic, demosaic_cfa};
use crate::error::{RawError, Result};
use crate::look::{CameraLook, ToneCurve};
use crate::sensor::Sensor;
@@ -194,12 +194,6 @@ pub(crate) fn resolve_neutral(s: &Sensor, opts: &DevelopOptions, sc: &Scale, war
/// Develops sensor data.
pub fn develop_sensor(s: &Sensor, opts: &DevelopOptions) -> Result<Developed> {
develop(s, opts, true)
}
/// [`develop_sensor`]; `allow_fused` false keeps the full-frame path (what the tests compare the
/// fused one against).
fn develop(s: &Sensor, opts: &DevelopOptions, allow_fused: bool) -> Result<Developed> {
let c = s.crop;
let samples = s.width.checked_mul(s.height).and_then(|n| n.checked_mul(s.samples));
if c.is_empty()
@@ -225,12 +219,7 @@ fn develop(s: &Sensor, opts: &DevelopOptions, allow_fused: bool) -> Result<Devel
return Err(RawError::malformed("invalid sensor CFA pattern"));
}
let bayer = s.cfa.as_ref().is_none_or(|cfa| cfa.is_bayer());
// The default rendering (no Dehaze, no local tone, no black point) needs no statistic of
// the whole frame: demosaic and encode band by band, with no full RGB frame (the same
// condition under which `Haze::new` and `ToneMap::new` return `None`).
let fused = allow_fused && opts.tone.is_neutral() && s.black_point <= 0.0;
let orient_copy = opts.orient && (2..=8).contains(&s.orientation);
let ahd_jobs = working_memory(c.width, c.height, s.cfa.is_some(), bayer, Shape { fused, orient_copy }, opts, par::workers())?;
let ahd_jobs = working_memory(c.width, c.height, s.cfa.is_some(), bayer, opts, par::workers())?;
let mut warnings = s.warnings.clone();
let sc = Scale::new(s);
@@ -267,207 +256,161 @@ fn develop(s: &Sensor, opts: &DevelopOptions, allow_fused: bool) -> Result<Devel
let grey_gain = gain as f32;
let (w, h) = (c.width, c.height);
let band = par::band_rows(w);
let mut rgb: Vec<f32>;
let pattern: Option<String>;
match &s.cfa {
Some(cfa) => {
// Normalized, balanced, clipped CFA plane.
let mut plane = Padded::new(w, h);
let stride = plane.stride;
// Fast path: black varies at most per 2×2 position, so black and gain
// repeat every `period` columns (the CFA width, made even).
let bl = &s.black;
let periodic = bl.delta_h.is_empty() && bl.delta_v.is_empty() && matches!(bl.rows, 1 | 2) && matches!(bl.cols, 1 | 2) && cfa.width <= 16;
let period = if cfa.width % 2 == 0 { cfa.width.max(2) } else { cfa.width * 2 };
plane.fill_rows(band, |y0, chunk| {
for (r, prow) in chunk.chunks_exact_mut(stride).enumerate() {
let row = &mut prow[PAD..PAD + w];
let y = c.y + y0 + r;
if periodic {
// Per column of the period: black and the combined scale × white-balance gain.
let k: Vec<(f32, f32)> = (0..period)
.map(|p| {
let x = c.x + p;
let black = bl.at(x.saturating_sub(s.active.x), y.saturating_sub(s.active.y), 0, 1);
let range = s.white[0] - black;
// A non-finite black (a BlackLevel DOUBLE beyond f32) gives 0, as
// `Scale::value` does, instead of (raw - inf) * 0 = NaN.
if !(range >= 1.0 && black.is_finite()) {
return (0.0, 0.0);
}
(black, mult[usize::from(cfa.color(x, y)).min(2)] / range)
})
.collect();
let src = &s.data[y * s.width + c.x..y * s.width + c.x + w];
let ks = k.iter().cycle();
match &sc.lin {
None => {
for ((v, &raw), &(black, gain)) in row.iter_mut().zip(src).zip(ks) {
*v = ((f32::from(raw) - black) * gain).clamp(0.0, 1.0);
}
}
Some(l) => {
for ((v, &raw), &(black, gain)) in row.iter_mut().zip(src).zip(ks) {
let lin = f32::from(l.get(usize::from(raw)).or_else(|| l.last()).copied().unwrap_or(raw));
*v = ((lin - black) * gain).clamp(0.0, 1.0);
}
}
}
} else {
for (i, v) in row.iter_mut().enumerate() {
let x = c.x + i;
let ch = usize::from(cfa.color(x, y)).min(2);
*v = (sc.value(x, y, 0) * mult[ch]).clamp(0.0, 1.0);
}
}
// Lens-shading gain maps, after the white-balance clip so clipped
// highlights stay neutral.
if let Some(g) = gain_row(s, y, c.x, w, 0) {
for (v, g) in row.iter_mut().zip(g) {
*v = (*v * g).clamp(0.0, 1.0);
}
}
}
});
plane.fill_borders();
let names = |p: &[u8]| p.iter().map(|c| ['R', 'G', 'B'][usize::from(*c).min(2)]).collect::<String>();
if bayer {
let p = cfa.phase(c.x, c.y);
pattern = Some(names(&p));
rgb = demosaic(&plane, p, opts.demosaic, &to_xyz, ahd_jobs);
} else {
let pat = Pattern::of(cfa, c.x, c.y);
pattern = Some(pat.colors.chunks(pat.w.max(1)).map(names).collect::<Vec<_>>().join("/"));
rgb = demosaic_cfa(&plane, &pat);
}
drop(plane);
}
None => {
pattern = None;
rgb = vec![0.0f32; w * h * 3];
par::chunks_mut(&mut rgb, band * w * 3, |b, chunk| {
for (r, row) in chunk.chunks_exact_mut(w * 3).enumerate() {
let y = c.y + b * band + r;
for i in 0..w {
for k in 0..3 {
// Monochrome: the one sample in all three channels.
let sample = k.min(s.samples.saturating_sub(1));
row[i * 3 + k] = (sc.value(c.x + i, y, sample) * mult[k]).clamp(0.0, 1.0);
}
}
for k in 0..s.samples.min(3) {
if let Some(g) = gain_row(s, y, c.x, w, k) {
for (i, g) in g.into_iter().enumerate() {
row[i * 3 + k] = (row[i * 3 + k] * g).clamp(0.0, 1.0);
}
}
}
}
});
}
}
// Camera → ProPhoto, exposure, Dehaze, Shadows / Highlights / Clarity / Texture and the black
// point, the look's (or the profile's) tone curve, gamma 1.8, 16 bits.
let y_row = color::rgb_to_xyz(color::ROMM, color::D50_XY)[1].map(|v| v as f32);
let lum = |v: &[f32; 3]| y_row[0] * v[0] + y_row[1] * v[1] + y_row[2] * v[2];
// Monochrome: the one sample, grey; no colour matrix.
let pro = |p: &[f32]| if mono { [grey_gain * p[0]; 3] } else { [0, 1, 2].map(|k| out_m[k][0] * p[0] + out_m[k][1] * p[1] + out_m[k][2] * p[2]) };
let haze = Haze::new(w, h, |i| pro(&rgb[i * 3..i * 3 + 3]), opts.tone.dehaze);
// Linear ProPhoto of pixel `i` after Dehaze.
let linear = |i: usize, p: &[f32]| {
let v = pro(p);
haze.as_ref().map_or(v, |z| z.apply(i % w, i / w, v))
};
let mut map = ToneMap::new(w, h, |i| lum(&linear(i, &rgb[i * 3..i * 3 + 3])), opts.tone, s.black_point);
if let (Some(m), Some(_)) = (&mut map, &haze) {
// The black point belongs to the default rendering: Dehaze does not move it (taken after
// Dehaze it fell to almost 0, and Dehaze darkened less than Photoshop's).
m.black = black_level(w, h, |i| lum(&pro(&rgb[i * 3..i * 3 + 3])), s.black_point);
}
let mut out = vec![0u16; w * h * 3];
let lut = match look {
Some(l) => GammaLut::with_tone(l.tone.as_ref()),
None => GammaLut::with_curve(&s.tone_curve),
};
// The default rendering: matrix then encode, in blocks so the encode's clamp and square root
// run over whole slices.
let encode_plain = |src: &[f32], dst: &mut [u16]| {
let mut buf = [0.0f32; ENCODE_BLOCK];
for (oc, sc) in dst.chunks_mut(ENCODE_BLOCK).zip(src.chunks(ENCODE_BLOCK)) {
let n = sc.len();
for (t, p) in buf[..n].as_chunks_mut::<3>().0.iter_mut().zip(sc.as_chunks::<3>().0) {
*t = pro(p);
// The default rendering (no Dehaze, no local tone): matrix then encode, in blocks so the
// encode's clamp and square root run over whole slices.
let plain = map.is_none() && haze.is_none();
par::chunks_mut(&mut out, band * w * 3, |b, chunk| {
let start = b * band * w * 3;
let src = &rgb[start..start + chunk.len()];
if plain {
let mut buf = [0.0f32; ENCODE_BLOCK];
for (oc, sc) in chunk.chunks_mut(ENCODE_BLOCK).zip(src.chunks(ENCODE_BLOCK)) {
let n = sc.len();
for (t, p) in buf[..n].as_chunks_mut::<3>().0.iter_mut().zip(sc.as_chunks::<3>().0) {
*t = pro(p);
}
lut.encode_into(&buf[..n], oc);
}
lut.encode_into(&buf[..n], oc);
return;
}
};
// Rows `y0..` of a LinearRaw (3-sample or monochrome) frame, balanced and clipped.
let linear_rows = |y0: usize, chunk: &mut [f32]| {
for (r, row) in chunk.chunks_exact_mut(w * 3).enumerate() {
let y = c.y + y0 + r;
for i in 0..w {
for k in 0..3 {
// Monochrome: the one sample in all three channels.
let sample = k.min(s.samples.saturating_sub(1));
row[i * 3 + k] = (sc.value(c.x + i, y, sample) * mult[k]).clamp(0.0, 1.0);
}
for (i, (o, p)) in chunk.as_chunks_mut::<3>().0.iter_mut().zip(src.as_chunks::<3>().0).enumerate() {
let mut v = linear(b * band * w + i, p);
if let Some(m) = &map {
let g = m.gain(i % w, b * band + i / w, lum(&v));
// The black point first, as Camera Raw: lifting before it would flatten the
// deepest shadows' contrast.
v = v.map(|c| m.unblack(c) * g);
}
for k in 0..s.samples.min(3) {
if let Some(g) = gain_row(s, y, c.x, w, k) {
for (i, g) in g.into_iter().enumerate() {
row[i * 3 + k] = (row[i * 3 + k] * g).clamp(0.0, 1.0);
}
}
for k in 0..3 {
o[k] = lut.encode(v[k]);
}
}
};
let mut rgb: Vec<f32> = Vec::new();
// The encoded frame, when the fused path made it.
let mut fused_out: Option<Vec<u16>> = None;
let pattern: Option<String>;
match &s.cfa {
Some(cfa) => {
// One normalized, balanced, clipped CFA row (image row `iy` of the crop).
// Fast path: black varies at most per 2×2 position, so black and gain
// repeat every `period` columns (the CFA width, made even).
let bl = &s.black;
let periodic = bl.delta_h.is_empty() && bl.delta_v.is_empty() && matches!(bl.rows, 1 | 2) && matches!(bl.cols, 1 | 2) && cfa.width <= 16;
let period = if cfa.width % 2 == 0 { cfa.width.max(2) } else { cfa.width * 2 };
let fill_row = |iy: usize, row: &mut [f32]| {
let y = c.y + iy;
if periodic {
// Per column of the period: black and the combined scale × white-balance gain.
let k: Vec<(f32, f32)> = (0..period)
.map(|p| {
let x = c.x + p;
let black = bl.at(x.saturating_sub(s.active.x), y.saturating_sub(s.active.y), 0, 1);
let range = s.white[0] - black;
// A non-finite black (a BlackLevel DOUBLE beyond f32) gives 0, as
// `Scale::value` does, instead of (raw - inf) * 0 = NaN.
if !(range >= 1.0 && black.is_finite()) {
return (0.0, 0.0);
}
(black, mult[usize::from(cfa.color(x, y)).min(2)] / range)
})
.collect();
let src = &s.data[y * s.width + c.x..y * s.width + c.x + w];
let ks = k.iter().cycle();
match &sc.lin {
None => {
for ((v, &raw), &(black, gain)) in row.iter_mut().zip(src).zip(ks) {
*v = ((f32::from(raw) - black) * gain).clamp(0.0, 1.0);
}
}
Some(l) => {
for ((v, &raw), &(black, gain)) in row.iter_mut().zip(src).zip(ks) {
let lin = f32::from(l.get(usize::from(raw)).or_else(|| l.last()).copied().unwrap_or(raw));
*v = ((lin - black) * gain).clamp(0.0, 1.0);
}
}
}
} else {
for (i, v) in row.iter_mut().enumerate() {
let x = c.x + i;
let ch = usize::from(cfa.color(x, y)).min(2);
*v = (sc.value(x, y, 0) * mult[ch]).clamp(0.0, 1.0);
}
}
// Lens-shading gain maps, after the white-balance clip so clipped
// highlights stay neutral.
if let Some(g) = gain_row(s, y, c.x, w, 0) {
for (v, g) in row.iter_mut().zip(g) {
*v = (*v * g).clamp(0.0, 1.0);
}
}
};
let names = |p: &[u8]| p.iter().map(|c| ['R', 'G', 'B'][usize::from(*c).min(2)]).collect::<String>();
if bayer && fused {
// No padded plane either: every band normalizes its own rows.
let p = cfa.phase(c.x, c.y);
pattern = Some(names(&p));
let mut o = vec![0u16; w * h * 3];
demosaic_fused(w, h, p, opts.demosaic, &to_xyz, ahd_jobs, &fill_row, &mut o, &encode_plain);
fused_out = Some(o);
} else {
let mut plane = Padded::new(w, h);
let stride = plane.stride;
plane.fill_rows(band, |y0, chunk| {
for (r, prow) in chunk.chunks_exact_mut(stride).enumerate() {
fill_row(y0 + r, &mut prow[PAD..PAD + w]);
}
});
plane.fill_borders();
if bayer {
let p = cfa.phase(c.x, c.y);
pattern = Some(names(&p));
rgb = demosaic(&plane, p, opts.demosaic, &to_xyz, ahd_jobs);
} else {
let pat = Pattern::of(cfa, c.x, c.y);
pattern = Some(pat.colors.chunks(pat.w.max(1)).map(names).collect::<Vec<_>>().join("/"));
if fused {
let mut o = vec![0u16; w * h * 3];
demosaic_cfa_fused(&plane, &pat, &mut o, &encode_plain);
fused_out = Some(o);
} else {
rgb = demosaic_cfa(&plane, &pat);
}
}
}
}
None => {
pattern = None;
if fused {
let mut o = vec![0u16; w * h * 3];
for_bands(
&mut o,
band * w * 3,
h.div_ceil(band),
usize::MAX,
|| vec![0.0f32; band * w * 3],
|b, chunk, buf| {
let buf = &mut buf[..chunk.len()];
linear_rows(b * band, buf);
encode_plain(buf, chunk);
},
);
fused_out = Some(o);
} else {
rgb = vec![0.0f32; w * h * 3];
par::chunks_mut(&mut rgb, band * w * 3, |b, chunk| linear_rows(b * band, chunk));
}
}
}
let out = match fused_out {
Some(out) => out,
None => {
let haze = Haze::new(w, h, |i| pro(&rgb[i * 3..i * 3 + 3]), opts.tone.dehaze);
// Linear ProPhoto of pixel `i` after Dehaze.
let linear = |i: usize, p: &[f32]| {
let v = pro(p);
haze.as_ref().map_or(v, |z| z.apply(i % w, i / w, v))
};
let mut map = ToneMap::new(w, h, |i| lum(&linear(i, &rgb[i * 3..i * 3 + 3])), opts.tone, s.black_point);
if let (Some(m), Some(_)) = (&mut map, &haze) {
// The black point belongs to the default rendering: Dehaze does not move it (taken after
// Dehaze it fell to almost 0, and Dehaze darkened less than Photoshop's).
m.black = black_level(w, h, |i| lum(&pro(&rgb[i * 3..i * 3 + 3])), s.black_point);
}
let mut out = vec![0u16; w * h * 3];
let plain = map.is_none() && haze.is_none();
par::chunks_mut(&mut out, band * w * 3, |b, chunk| {
let start = b * band * w * 3;
let src = &rgb[start..start + chunk.len()];
if plain {
encode_plain(src, chunk);
return;
}
for (i, (o, p)) in chunk.as_chunks_mut::<3>().0.iter_mut().zip(src.as_chunks::<3>().0).enumerate() {
let mut v = linear(b * band * w + i, p);
if let Some(m) = &map {
let g = m.gain(i % w, b * band + i / w, lum(&v));
// The black point first, as Camera Raw: lifting before it would flatten the
// deepest shadows' contrast.
v = v.map(|c| m.unblack(c) * g);
}
for k in 0..3 {
o[k] = lut.encode(v[k]);
}
}
});
rgb.clear();
rgb.shrink_to_fit();
out
}
};
});
rgb.clear();
rgb.shrink_to_fit();
let (out, ow, oh) = if opts.orient { orient(out, w, h, s.orientation) } else { (out, w, h) };
let nmax = mult.iter().cloned().fold(f32::MAX, f32::min).max(1e-6);
@@ -490,76 +433,29 @@ fn develop(s: &Sensor, opts: &DevelopOptions, allow_fused: bool) -> Result<Devel
})
}
/// How the frame is developed, for [`working_memory`].
#[derive(Clone, Copy)]
struct Shape {
/// Demosaic and encode band by band (no full RGB frame).
fused: bool,
/// The orientation copies the encoded frame.
orient_copy: bool,
}
/// Bytes one worker of the fused path holds for a band of `band_rows` rows: the f32 RGB band
/// and, for a Bayer sensor, the padded CFA window it reads.
fn fused_band_bytes(w: usize, h: usize, band_rows: usize, window: bool) -> u64 {
let rows = band_rows.min(h) as u64;
let rgb = rows.saturating_mul(w as u64).saturating_mul(12);
let window = if window { (rows + 2 * PAD as u64).saturating_mul((w + 2 * PAD) as u64).saturating_mul(4) } else { 0 };
rgb.saturating_add(window)
}
/// Checks the peak of the two development phases. The CFA and non-Bayer green
/// planes are freed before the u16 output is allocated; AHD scratch is bounded
/// by reducing simultaneous bands when needed. The fused path never holds the f32 RGB frame:
/// it has the u16 output (6 bytes a pixel) next to the CFA plane, plus a band buffer and the
/// scratch per worker, and an orientation copy of the output afterwards. Sensor storage belongs
/// to the caller and is checked separately by the decoder.
fn working_memory(w: usize, h: usize, cfa: bool, bayer: bool, shape: Shape, opts: &DevelopOptions, workers: usize) -> Result<usize> {
/// by reducing simultaneous bands when needed. Sensor storage belongs to the
/// caller and is checked separately by the decoder.
fn working_memory(w: usize, h: usize, cfa: bool, bayer: bool, opts: &DevelopOptions, workers: usize) -> Result<usize> {
let limits = &opts.limits;
let Shape { fused, orient_copy } = shape;
limits.check(w as u64, h as u64, if fused { 6 } else { 18 })?;
limits.check(w as u64, h as u64, 18)?;
let pixels = (w as u64).saturating_mul(h as u64);
let table = (GAMMA_STEPS as u64 + 2) * 4;
// Encode phase: the f32 RGB frame and the output, or (fused) the output and its oriented copy.
let output = pixels.saturating_mul(match (fused, orient_copy) {
(false, _) => 18,
(true, false) => 6,
(true, true) => 12,
});
let mut peak = output.saturating_add(table);
let output = pixels.saturating_mul(18).saturating_add((GAMMA_STEPS as u64 + 2) * 4);
let mut peak = output;
let mut jobs = workers.max(1);
// Fused: the u16 output and, per worker running a band, one band of f32 RGB.
let band_of = |rows: usize, window: bool| if fused { fused_band_bytes(w, h, rows, window) } else { 0 };
let band_jobs = |rows: usize| workers.max(1).min(h.div_ceil(rows.max(1))) as u64;
if cfa {
let plane = (w as u64).saturating_add(2 * PAD as u64).saturating_mul((h as u64).saturating_add(2 * PAD as u64)).saturating_mul(4);
// Bytes a pixel next to the plane: the f32 RGB frame (and the non-Bayer green plane), or
// the u16 output (and the green plane, which the fused path keeps).
let frame = match (fused, bayer) {
(false, true) => 12,
(false, false) => 16,
(true, true) => 6,
(true, false) => 10,
};
// The fused Bayer path has no whole plane: each worker's window is counted with its band.
let plane = if fused && bayer { 0 } else { plane };
let mut mosaic = pixels.saturating_mul(frame).saturating_add(plane).saturating_add(table);
let mut mosaic = pixels.saturating_mul(if bayer { 12 } else { 16 }).saturating_add(plane);
if bayer && opts.demosaic == Demosaic::Ahd {
mosaic = mosaic.saturating_add(crate::demosaic::AHD_TABLE_BYTES);
let scratch = crate::demosaic::ahd_scratch(w as u64, h as u64).saturating_add(band_of(crate::demosaic::AHD_BAND, true));
let scratch = crate::demosaic::ahd_scratch(w as u64, h as u64);
let available = limits.max_alloc.saturating_sub(mosaic) / scratch.max(1);
jobs = jobs.min(crate::demosaic::ahd_bands(h)).min(usize::try_from(available).unwrap_or(usize::MAX));
// At least one band must fit, even when no parallelism is available.
mosaic = mosaic.saturating_add(scratch.saturating_mul(jobs.max(1) as u64));
} else {
let rows = if bayer { crate::demosaic::FUSED_ROWS } else { par::band_rows(w) };
mosaic = mosaic.saturating_add(band_of(rows, bayer).saturating_mul(band_jobs(rows)));
}
peak = peak.max(mosaic);
} else if fused {
// LinearRaw: the output and the bands.
let rows = par::band_rows(w);
peak = peak.max(pixels.saturating_mul(6).saturating_add(band_of(rows, false).saturating_mul(band_jobs(rows))));
}
if peak > limits.max_alloc || usize::try_from(peak).is_err() {
return Err(RawError::LimitExceeded(format!("{peak} bytes of RAW working memory exceed {}", limits.max_alloc)));
@@ -738,136 +634,6 @@ pub(crate) fn orient(src: Vec<u16>, w: usize, h: usize, o: u16) -> (Vec<u16>, us
mod tests {
use super::*;
const FULL: Shape = Shape { fused: false, orient_copy: false };
const FUSED: Shape = Shape { fused: true, orient_copy: false };
/// A `w × h` sensor with a pseudo-random gradient, black 64, white 4000.
fn test_sensor(w: usize, h: usize, cfa: Option<crate::Cfa>, samples: usize) -> Sensor {
use crate::*;
let data = (0..w * h * samples).map(|i| 64 + ((i * 7919 + i / w * 104_729 + i / 5 * 31) % 4100) as u16).collect();
Sensor {
format: RawFormat::Dng,
make: None,
model: None,
width: w,
height: h,
samples,
data,
cfa,
linearization: None,
black: BlackLevels::uniform(64.0),
white: [4000.0; 3],
active: Rect::new(0, 0, w, h),
crop: Rect::new(0, 0, w, h),
color: ColorInfo::default(),
camera_wb: None,
orientation: 1,
baseline_exposure: 0.0,
gain_maps: Vec::new(),
tone_curve: vec![[0.18, 0.3], [1.0, 1.0]],
black_point: 0.0,
warnings: Vec::new(),
}
}
/// `f` on a pool of `n` threads (the browser builds have none).
fn with_threads<R: Send>(n: usize, f: impl FnOnce() -> R + Send) -> R {
#[cfg(not(target_arch = "wasm32"))]
return rayon::ThreadPoolBuilder::new().num_threads(n).build().unwrap().install(f);
#[cfg(target_arch = "wasm32")]
{
let _ = n;
f()
}
}
fn assert_fused_matches_full(s: &Sensor, opts: &DevelopOptions, what: &str) {
let full = develop(s, opts, false).unwrap();
for threads in [1, 3, 16] {
let fused = with_threads(threads, || develop(s, opts, true).unwrap());
assert_eq!((full.width, full.height), (fused.width, fused.height), "{what} {threads}");
assert!(full.rgb == fused.rgb, "{what} {threads} threads: fused output differs");
}
}
#[test]
fn fused_development_matches_the_full_frame_path() {
use crate::Cfa;
// Odd sizes, tiny images, band borders (AHD bands are 32 rows) and several bands.
let sizes = [(1, 1), (2, 3), (5, 3), (33, 31), (64, 64), (65, 97), (70, 131), (3, 300)];
for (w, h) in sizes {
for method in [Demosaic::Ahd, Demosaic::Bilinear, Demosaic::Mhc] {
for orientation in [1, 6] {
for pattern in ["RGGB", "BGGR", "GRBG"] {
let mut s = test_sensor(w, h, Cfa::bayer(pattern), 1);
s.orientation = orientation;
let opts = DevelopOptions { demosaic: method, ..Default::default() };
assert_fused_matches_full(&s, &opts, &format!("{w}x{h} {method:?} o{orientation} {pattern}"));
}
}
}
}
}
#[test]
fn fused_development_matches_for_every_orientation_and_look() {
use crate::Cfa;
for orientation in 1..=8 {
let mut s = test_sensor(67, 45, Cfa::bayer("RGGB"), 1);
s.orientation = orientation;
s.baseline_exposure = 0.5;
s.linearization = Some((0..4100u32).map(|v| (v * 15 / 16) as u16).collect());
assert_fused_matches_full(&s, &DevelopOptions::default(), &format!("orientation {orientation}"));
}
}
#[test]
fn fused_development_matches_for_xtrans_and_linear_raw() {
use crate::Cfa;
// 6 × 6 X-Trans repeat.
let rows = ["GGRGGB", "GGBGGR", "BRGRBG", "GGBGGR", "GGRGGB", "RBGBRG"];
let colors: Vec<u8> = rows
.iter()
.flat_map(|r| {
r.bytes().map(|c| match c {
b'R' => 0,
b'G' => 1,
_ => 2,
})
})
.collect();
for (w, h) in [(6, 6), (37, 41), (80, 150)] {
let mut s = test_sensor(w, h, Some(Cfa { width: 6, height: 6, colors: colors.clone(), origin_x: 0, origin_y: 0 }), 1);
s.orientation = 8;
assert_fused_matches_full(&s, &DevelopOptions::default(), &format!("xtrans {w}x{h}"));
}
for samples in [1, 3] {
for (w, h) in [(1, 1), (7, 5), (90, 111)] {
let mut s = test_sensor(w, h, None, samples);
s.orientation = 3;
assert_fused_matches_full(&s, &DevelopOptions::default(), &format!("linear {samples} {w}x{h}"));
}
}
}
#[test]
fn fused_only_where_no_whole_frame_statistic_is_needed() {
use crate::Cfa;
// The condition for the fused path is exactly "Haze and ToneMap are None".
let px = |_: usize| [0.5f32; 3];
for tone in [Tone::NEUTRAL, Tone { shadows: 10.0, ..Tone::NEUTRAL }, Tone { dehaze: 5.0, ..Tone::NEUTRAL }, Tone { texture: 1.0, ..Tone::NEUTRAL }] {
for black_point in [0.0, 0.5] {
let none = Haze::new(8, 8, px, tone.dehaze).is_none() && ToneMap::new(8, 8, |_| 0.5, tone, black_point).is_none();
assert_eq!(none, tone.is_neutral() && black_point <= 0.0, "{tone:?} {black_point}");
}
}
// A tone slider or a black point is still developed (differently from the plain rendering).
let s = test_sensor(40, 40, Cfa::bayer("RGGB"), 1);
let plain = develop_sensor(&s, &DevelopOptions::default()).unwrap();
let toned = develop_sensor(&s, &DevelopOptions { tone: Tone { shadows: 50.0, ..Tone::NEUTRAL }, ..Default::default() }).unwrap();
assert!(plain.rgb != toned.rgb);
}
#[test]
fn block_encode_matches_per_sample_encode() {
let lut = GammaLut::with_curve(&[[0.25, 0.5], [1.0, 1.0]]);
@@ -884,57 +650,27 @@ mod tests {
let mut opts = DevelopOptions::default();
for method in [Demosaic::Bilinear, Demosaic::Mhc, Demosaic::Ahd] {
opts.demosaic = method;
for (fused, orient_copy) in [(false, false), (true, false), (true, true)] {
let shape = Shape { fused, orient_copy };
for workers in [1, 10, 32, 128] {
let jobs = working_memory(11648, 8736, true, true, shape, &opts, workers).unwrap();
assert!(jobs <= workers);
if method == Demosaic::Ahd && workers >= 32 && !fused {
assert!(jobs < workers);
}
for workers in [1, 10, 32, 128] {
let jobs = working_memory(11648, 8736, true, true, &opts, workers).unwrap();
assert!(jobs <= workers);
if method == Demosaic::Ahd && workers >= 32 {
assert!(jobs < workers);
}
}
}
opts.limits.max_alloc = 1 << 30;
assert!(working_memory(11648, 8736, true, true, FULL, &opts, 1).is_err());
// Fused: the u16 output and the CFA plane (+ one band) is about 10 bytes a pixel.
opts.demosaic = Demosaic::Bilinear;
assert!(working_memory(11648, 8736, true, true, FUSED, &opts, 1).is_ok());
assert!(working_memory(11648, 8736, true, true, &opts, 1).is_err());
}
#[test]
fn working_memory_accounts_for_padding_and_at_least_one_ahd_band() {
let mut opts = DevelopOptions { demosaic: Demosaic::Ahd, ..Default::default() };
opts.limits.max_alloc = 100 * 100 * 22;
assert!(working_memory(100, 100, true, true, FULL, &opts, 1).is_err());
assert!(working_memory(100, 100, true, true, &opts, 1).is_err());
opts.demosaic = Demosaic::Mhc;
assert!(working_memory(100, 100, true, true, FULL, &opts, 1).is_ok());
assert!(working_memory(100, 100, true, true, &opts, 1).is_ok());
opts.limits.max_alloc = 100 * 100 * 18;
assert!(working_memory(100, 100, false, true, FULL, &opts, 1).is_err());
}
#[test]
fn fused_working_memory_charges_the_output_the_window_and_the_bands() {
let mut opts = DevelopOptions { demosaic: Demosaic::Mhc, ..Default::default() };
let (w, h) = (1000u64, 1000u64);
let table = (GAMMA_STEPS as u64 + 2) * 4;
let plane = (w + 2 * PAD as u64) * (h + 2 * PAD as u64) * 4;
let pad = PAD as u64;
let exact = |opts: &mut DevelopOptions, cfa: bool, bayer: bool, shape: Shape, need: u64| {
opts.limits.max_alloc = need - 1;
assert!(working_memory(1000, 1000, cfa, bayer, shape, opts, 1).is_err(), "{need} - 1");
opts.limits.max_alloc = need;
assert!(working_memory(1000, 1000, cfa, bayer, shape, opts, 1).is_ok(), "{need}");
};
// Bayer: the u16 output and one band with its padded window, no whole plane.
let band = 32 * w * 12 + (32 + 2 * pad) * (w + 2 * pad) * 4;
exact(&mut opts, true, true, FUSED, w * h * 6 + table + band);
// An oriented frame needs the output twice; a monochrome sensor has no plane.
let mono_band = 65 * w * 12;
exact(&mut opts, false, true, Shape { fused: true, orient_copy: true }, w * h * 12 + table);
exact(&mut opts, false, true, FUSED, w * h * 6 + mono_band);
// X-Trans keeps its plane and its green plane.
exact(&mut opts, true, false, FUSED, w * h * 10 + plane + table + mono_band);
assert!(working_memory(100, 100, false, true, &opts, 1).is_err());
}
#[test]
+6 -18
View File
@@ -8,25 +8,13 @@ fn small_sensor() -> Sensor {
#[test]
fn checks_peak_live_buffers_instead_of_summing_disjoint_phases() {
let (w, h) = (128u64, 96u64);
let sensor = decode(&DngSpec::cfa(w as usize, h as usize, vec![1000; (w * h) as usize]).build(), &Limits::default()).unwrap();
let table = (4096 + 2) * 4;
// Exactly `max_alloc` must work and one byte less must fail.
let check = |opts: &mut DevelopOptions, max_alloc: u64| {
let expected = develop_sensor(&sensor, opts).unwrap();
opts.limits.max_alloc = max_alloc;
assert_eq!(develop_sensor(&sensor, opts).unwrap().rgb, expected.rgb);
opts.limits.max_alloc -= 1;
assert!(matches!(develop_sensor(&sensor, opts), Err(RawError::LimitExceeded(_))));
};
// A tone slider needs the f32 RGB frame next to the output (18 bytes a pixel) ...
let mut opts = DevelopOptions { demosaic: Demosaic::Mhc, tone: Tone { shadows: 20.0, ..Tone::NEUTRAL }, ..Default::default() };
check(&mut opts, w * h * 18 + table);
// ... the default rendering the u16 output and, per worker, one band (32 rows) of f32 RGB with
// the 12 rows more of padded CFA it reads, but no whole plane.
let sensor = decode(&DngSpec::cfa(128, 96, vec![1000; 128 * 96]).build(), &Limits::default()).unwrap();
let mut opts = DevelopOptions { demosaic: Demosaic::Mhc, ..Default::default() };
let band = 32 * w * 12 + (32 + 12) * (w + 12) * 4;
rayon::ThreadPoolBuilder::new().num_threads(1).build().unwrap().install(|| check(&mut opts, w * h * 6 + table + band));
opts.limits.max_alloc = 128 * 96 * 18 + (4096 + 2) * 4;
let expected = develop_sensor(&sensor, &DevelopOptions { demosaic: Demosaic::Mhc, ..Default::default() }).unwrap();
assert_eq!(develop_sensor(&sensor, &opts).unwrap().rgb, expected.rgb);
opts.limits.max_alloc -= 1;
assert!(matches!(develop_sensor(&sensor, &opts), Err(RawError::LimitExceeded(_))));
}
#[test]