perf(pipeline): Highlights/Shadows median about 2x faster on the CPU, bit for bit the same

`median_row` (the running median under process V2's highlights/shadows base) kept a sorted
window of f32 and, per pixel, found the leaving and entering values with two binary searches
on `total_cmp`, then `remove` + `insert`. The searches mispredicted on most steps.

Now the window holds `i32` keys that order exactly like `f32::total_cmp` (the same transform
it applies), both ranks come from one branch-free counting pass the compiler vectorizes, and
one `copy_within` moves the keys between them. Equal keys are equal bits, so the output is
unchanged bit for bit.

24 MP Sony A7 III log luminance, one P-core, interleaved: r = 9 (export) 1977-2027 ms ->
912-958 ms, r = 3 1249-1292 -> 658-724, r = 24 2180-2261 -> 972-1030. All cores, r = 9:
170-188 -> 86-102 ms. The GPU path has its own median and is untouched.

New test: special values (NaN of both signs, signed zeros, infinities), long ties and planes
narrower than the window match a per-pixel sorted median bit for bit.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_0121DZuByhKyGLHw5C9JkUue

crates/pipeline/src/local.rs

@@ -227,22 +227,47 @@ fn median_row(src: &[f32], out: &mut [f32], r: usize) {
     if n == 0 {
         return;
     }
-    let at = |i: isize| src[i.clamp(0, n as isize - 1) as usize];
+    let at = |i: isize| total_key(src[i.clamp(0, n as isize - 1) as usize]);
     let r_ = r as isize;
-    let mut win: Vec<f32> = (-r_..=r_).map(at).collect();
-    win.sort_unstable_by(f32::total_cmp);
+    let mut win: Vec<i32> = (-r_..=r_).map(at).collect();
+    win.sort_unstable();
     for (x, o) in out.iter_mut().enumerate() {
-        *o = win[r];
+        *o = f32::from_bits(total_key_inv(win[r]));
         let (old, new) = (at(x as isize - r_), at(x as isize + r_ + 1));
-        if old.to_bits() != new.to_bits() {
-            let i = win.binary_search_by(|v| v.total_cmp(&old)).unwrap_or_else(|i| i.min(win.len() - 1));
-            win.remove(i);
-            let j = win.binary_search_by(|v| v.total_cmp(&new)).unwrap_or_else(|j| j);
-            win.insert(j, new);
+        if old != new {
+            // Ranks by counting the smaller keys in one branch-free pass the compiler vectorizes.
+            // Binary searches here mispredicted on most steps and, with `remove` + `insert`, cost
+            // 2.4 s per 24 MP at r = 9 on one core.
+            let (mut i, mut j) = (0u32, 0u32);
+            for &v in &win {
+                i += u32::from(v < old);
+                j += u32::from(v < new);
+            }
+            let (i, j) = (i as usize, j as usize);
+            // `old` sits at `i`; shift the keys between the two ranks by one over it
+            if j > i {
+                win.copy_within(i + 1..j, i);
+                win[j - 1] = new;
+            } else {
+                win.copy_within(j..i, j + 1);
+                win[j] = new;
+            }
         }
     }
 }
 
+/// `v`'s bits as an `i32` that orders like [`f32::total_cmp`] (the same transform it applies), so
+/// equal keys are equal bits.
+fn total_key(v: f32) -> i32 {
+    let b = v.to_bits() as i32;
+    b ^ (((b >> 31) as u32) >> 1) as i32
+}
+
+/// The bits [`total_key`] came from (the transform is its own inverse: it keeps the sign bit).
+fn total_key_inv(k: i32) -> u32 {
+    (k ^ (((k >> 31) as u32) >> 1) as i32) as u32
+}
+
 fn transpose(p: &Plane) -> Plane {
     let (w, h) = (p.width, p.height);
     let mut t = Plane::new(h, w);

crates/pipeline/src/tests_tone.rs

@@ -187,6 +187,27 @@ fn median_matches_a_plain_per_pixel_median() {
     }
 }
 
+#[test]
+fn median_orders_special_values_and_ties_like_total_cmp() {
+    // the running window keeps keys in `f32::total_cmp` order: NaNs of either sign, signed zeros,
+    // infinities and long runs of equal values come out bit for bit as a sorted per-pixel window
+    // does, also on planes narrower or shorter than the window
+    let special = [f32::NAN, -f32::NAN, 0.0, -0.0, f32::INFINITY, f32::NEG_INFINITY, f32::MIN_POSITIVE, -1.0, 1.0, 1.0, 1.0];
+    let bits = |p: &Plane| p.data.iter().map(|v| v.to_bits()).collect::<Vec<_>>();
+    for (w, h) in [(1, 1), (2, 3), (5, 1), (19, 7), (61, 9)] {
+        let p = Plane::from_fn(w, h, |x, y| special[(x * 7 + y * 3 + (x * y) % 5) % special.len()]);
+        for r in [1, 2, 9, 24] {
+            let med = |v: &mut Vec<f32>| {
+                v.sort_by(f32::total_cmp);
+                v[r]
+            };
+            let rows = Plane::from_fn(w, h, |x, y| med(&mut (0..=2 * r).map(|k| p.get((x + k).saturating_sub(r).min(w - 1), y)).collect()));
+            let want = Plane::from_fn(w, h, |x, y| med(&mut (0..=2 * r).map(|k| rows.get(x, (y + k).saturating_sub(r).min(h - 1))).collect()));
+            assert_eq!(bits(&local::median_hv(&p, r)), bits(&want), "{w}x{h} r {r}");
+        }
+    }
+}
+
 #[test]
 fn v1_keeps_the_base_it_always_had() {
     // V1's base is the edge-aware filter of log luminance itself, as before process V2: no median