perf(wifi-densepose-sar): incremental phasor rotation in backprojection (~4.4-4.5x, MEASURED)

focus_at_point called Complex64::from_polar (a sin/cos pair) once per
(pose, frequency) term. FrequencySweep::frequencies() produces evenly
spaced frequencies by construction, so the per-term phase is an
arithmetic progression in the frequency index -- the phasor can be
evaluated once per pose and advanced by a fixed complex-multiply step
per frequency instead, turning K trig evaluations into 2.

focus_at_point's signature changes from a raw &[f64] frequency slice
to &FrequencySweep, so the evenly-spaced-frequencies precondition
this optimization depends on is a type-level invariant rather than a
caller-observed one -- an arbitrary non-uniform frequency list is no
longer constructible through this API at all.

MEASURED (criterion regression detection, p < 0.001): ~4.4-4.5x
faster across 512/4096/32768-voxel grids (300us/1.97ms/14.5ms vs the
prior 1.47ms/10.4ms/73.5ms). Proven equivalent, not just faster: a new
test independently reimplements the direct per-frequency computation
as a reference and checks the optimized path against it across four
sweep sizes (incl. the n_steps=1 degenerate case) and both on-target
and off-target points, to <1e-9 relative error.

25 tests (22 unit + 3 integration), 0 failed, clippy-clean.

v2/crates/wifi-densepose-sar/README.md

@@ -42,7 +42,7 @@ cargo bench -p wifi-densepose-sar
 `tests/physics_validation.rs` checks the reconstruction's actual behavior
 against the closed-form formulas in `resolution.rs` (range resolution,
 cross-range/synthetic-aperture resolution, and the antenna-pose coherence
-budget) rather than merely asserting them: 24 tests (21 unit + 3
+budget) rather than merely asserting them: 25 tests (22 unit + 3
 integration), 0 failed, clippy-clean.
 
 ## Performance (MEASURED)
@@ -53,10 +53,22 @@ machine, release profile:
 
 | Voxels | Median time | Throughput |
 |-------:|------------:|-----------:|
-| 512    | 1.47 ms     | ~348,000 voxels/s |
-| 4,096  | 10.4 ms     | ~394,000 voxels/s |
-| 32,768 | 73.5 ms     | ~446,000 voxels/s |
+| 512    | 300 µs      | ~1.71M voxels/s |
+| 4,096  | 1.97 ms     | ~2.08M voxels/s |
+| 32,768 | 14.5 ms     | ~2.26M voxels/s |
 
 Scales as expected: each voxel's cost is independent (`O(poses × freqs)`
 per voxel, embarrassingly parallel), so throughput is roughly constant
 across grid sizes and total time scales linearly with voxel count.
+
+**Optimization (MEASURED, criterion regression detection, p < 0.001): ~4.4-4.5x
+faster** than the first-shipped implementation, across all three grid
+sizes. Frequencies in a [`FrequencySweep`](src/measurement.rs) are evenly
+spaced by construction, so the per-(pose, frequency) phase term is an
+arithmetic progression; `focus_at_point` now evaluates the phasor once per
+pose and advances it by a fixed complex-multiply step per frequency,
+instead of one `sin`/`cos` pair (`Complex64::from_polar`) per frequency —
+K trig evaluations become 2. Proven equivalent (not just faster) to an
+independently-reimplemented direct per-frequency reference in
+`reconstruct::tests::backprojection_incremental_rotation_matches_direct_per_frequency_computation`,
+across several sweep sizes and both on-target and off-target points.

v2/crates/wifi-densepose-sar/src/reconstruct.rs

@@ -120,8 +120,7 @@ pub fn backproject(
     grid: &VoxelGrid,
 ) -> ReflectivityImage {
     assert_eq!(poses.len(), measurement.n_poses, "pose count must match measurement");
-    let freqs = sweep.frequencies();
-    assert_eq!(freqs.len(), measurement.n_freqs, "frequency count must match measurement");
+    assert_eq!(sweep.n_steps, measurement.n_freqs, "frequency count must match measurement");
 
     let n = grid.len();
 
@@ -130,7 +129,7 @@ pub fn backproject(
         .map(|linear| {
             let (i, j, k) = grid.unflatten(linear);
             let voxel = grid.voxel_center(i, j, k);
-            focus_at_point(measurement, poses, &freqs, &voxel)
+            focus_at_point(measurement, poses, sweep, &voxel)
         })
         .collect();
 
@@ -143,11 +142,30 @@ pub fn backproject(
 /// (and tests) can measure focus quality exactly at a location of
 /// interest -- e.g. a known target position -- rather than only at
 /// whatever grid points happen to be sampled.
-pub fn focus_at_point(measurement: &Measurement, poses: &[AntennaPose], freqs: &[f64], point: &Point3) -> f64 {
+///
+/// Takes `sweep` rather than a raw frequency slice specifically so the
+/// evenly-spaced-frequencies guarantee ([`FrequencySweep::frequencies`])
+/// is a type-level invariant, not a caller-observed precondition: the
+/// implementation below relies on it (see the comment inside the pose
+/// loop). Passing an arbitrary non-uniform frequency list is not possible
+/// through this signature.
+pub fn focus_at_point(measurement: &Measurement, poses: &[AntennaPose], sweep: &FrequencySweep, point: &Point3) -> f64 {
     assert_eq!(poses.len(), measurement.n_poses, "pose count must match measurement");
-    assert_eq!(freqs.len(), measurement.n_freqs, "frequency count must match measurement");
+    assert_eq!(sweep.n_steps, measurement.n_freqs, "frequency count must match measurement");
 
     let n_terms = (measurement.n_poses * measurement.n_freqs) as f64;
+    let k = sweep.n_steps;
+    // Frequencies are evenly spaced by construction: f_kf = start_hz + kf *
+    // delta_f. That makes the per-term phase phase_kf = 4*pi*f_kf*r/c an
+    // arithmetic progression in kf, so instead of K trig evaluations
+    // (Complex64::from_polar per frequency step) the phasor is evaluated
+    // once and advanced by a fixed per-step rotation -- one complex
+    // multiply per step instead of a sin/cos pair. Proven equivalent to
+    // the direct per-frequency computation (independently reimplemented,
+    // not reusing this code) in
+    // `backprojection_incremental_rotation_matches_direct_per_frequency_computation`.
+    let delta_f = if k > 1 { (sweep.stop_hz - sweep.start_hz) / (k - 1) as f64 } else { 0.0 };
+
     let mut acc = Complex64::new(0.0, 0.0);
     for (m, pose) in poses.iter().enumerate() {
         let r = pose.position.distance(point);
@@ -155,9 +173,15 @@ pub fn focus_at_point(measurement: &Measurement, poses: &[AntennaPose], freqs: &
             continue;
         }
         let gain_compensation = r * r;
-        for (kf, &f) in freqs.iter().enumerate() {
-            let phase = 4.0 * PI * f * r / SPEED_OF_LIGHT_M_PER_S;
-            acc += measurement.get(m, kf) * gain_compensation * Complex64::from_polar(1.0, phase);
+        let base_phase = 4.0 * PI * sweep.start_hz * r / SPEED_OF_LIGHT_M_PER_S;
+        let step_phase = 4.0 * PI * delta_f * r / SPEED_OF_LIGHT_M_PER_S;
+        let step = Complex64::from_polar(1.0, step_phase);
+        let mut rot = Complex64::from_polar(1.0, base_phase);
+        for kf in 0..k {
+            acc += measurement.get(m, kf) * gain_compensation * rot;
+            if kf + 1 < k {
+                rot *= step;
+            }
         }
     }
     acc.norm() / n_terms
@@ -211,4 +235,65 @@ mod tests {
 
         assert!(peak_mag > mean_mag * 5.0, "coherent focus at the target must dominate the incoherent background: peak={peak_mag}, mean={mean_mag}");
     }
+
+    /// Independent reference: the direct per-frequency computation
+    /// `focus_at_point` used before the incremental-phasor-rotation
+    /// optimization (one `Complex64::from_polar` per (pose, frequency)
+    /// term, no recurrence). Deliberately reimplemented here rather than
+    /// calling any shared helper, so this test cannot pass by construction.
+    fn focus_at_point_direct_reference(
+        measurement: &Measurement,
+        poses: &[AntennaPose],
+        sweep: &FrequencySweep,
+        point: &Point3,
+    ) -> f64 {
+        let freqs = sweep.frequencies();
+        let n_terms = (measurement.n_poses * measurement.n_freqs) as f64;
+        let mut acc = Complex64::new(0.0, 0.0);
+        for (m, pose) in poses.iter().enumerate() {
+            let r = pose.position.distance(point);
+            if r < 1e-6 {
+                continue;
+            }
+            let gain_compensation = r * r;
+            for (kf, &f) in freqs.iter().enumerate() {
+                let phase = 4.0 * PI * f * r / SPEED_OF_LIGHT_M_PER_S;
+                acc += measurement.get(m, kf) * gain_compensation * Complex64::from_polar(1.0, phase);
+            }
+        }
+        acc.norm() / n_terms
+    }
+
+    /// PERF PROOF: the incremental-phasor-rotation `focus_at_point` (2
+    /// trig evaluations/pose instead of K) matches the direct
+    /// per-frequency reference to within f64 rounding, across several
+    /// sweep sizes, ranges, and off-axis points (not just the on-target
+    /// case, where errors could cancel).
+    #[test]
+    fn backprojection_incremental_rotation_matches_direct_per_frequency_computation() {
+        let poses = linear_aperture(Point3::new(-0.7, 0.0, 0.0), Point3::new(0.6, 0.1, 0.0), 17);
+        let targets = vec![
+            ScatteringTarget::new(Point3::new(0.1, 2.3, -0.2), 1.0),
+            ScatteringTarget::new(Point3::new(-0.4, 1.9, 0.3), 0.6),
+        ];
+        let test_points = [
+            Point3::new(0.1, 2.3, -0.2),  // on a target
+            Point3::new(-0.4, 1.9, 0.3),  // on the other target
+            Point3::new(0.0, 2.0, 0.0),   // off-target
+            Point3::new(-0.55, 2.6, 0.4), // off-target, far corner
+        ];
+        for &(n_steps, start_hz, stop_hz) in &[(1usize, 3.0e9, 3.0e9), (2, 2.0e9, 6.0e9), (8, 1.0e9, 9.0e9), (64, 2.4e9, 2.5e9)] {
+            let sweep = FrequencySweep::new(start_hz, stop_hz, n_steps);
+            let measurement = simulate_measurement(&poses, &sweep, &targets, 0.0, 42);
+            for point in test_points {
+                let fast = focus_at_point(&measurement, &poses, &sweep, &point);
+                let reference = focus_at_point_direct_reference(&measurement, &poses, &sweep, &point);
+                let scale = reference.max(1e-12);
+                assert!(
+                    (fast - reference).abs() / scale < 1e-9,
+                    "incremental rotation diverged from the direct reference at n_steps={n_steps}, point={point:?}: fast={fast}, reference={reference}"
+                );
+            }
+        }
+    }
 }

v2/crates/wifi-densepose-sar/tests/physics_validation.rs

@@ -179,7 +179,6 @@ fn phase_error_from_pose_jitter_degrades_focus_beyond_pose_budget() {
             .collect()
     };
 
-    let freqs = sweep.frequencies();
     let focus_at_epsilon = |epsilon: f64| -> f64 {
         let true_poses = jittered_poses_at(epsilon);
         // The measurement is recorded at the (jittered) TRUE antenna
@@ -191,7 +190,7 @@ fn phase_error_from_pose_jitter_degrades_focus_beyond_pose_budget() {
         // happens to focus slightly better and mask the coherence loss
         // this test is measuring).
         let measurement = simulate_measurement(&true_poses, &sweep, &[target], 0.0, 1);
-        focus_at_point(&measurement, &nominal_poses, &freqs, &target.position)
+        focus_at_point(&measurement, &nominal_poses, &sweep, &target.position)
     };
 
     let levels = [0.0, 0.5 * budget, budget, 2.0 * budget, 4.0 * budget];