Provable and Robust Wavefront Sensing via Self-Reference Interferometry¶
Conference: ECCV 2026
Paper: ECCV Official Page
Area: Medical Imaging
Keywords: Wavefront Sensing, Self-Reference Interferometry, Phase Retrieval, Co-Prime Shifts, Graph Propagation Algorithm
TL;DR¶
This paper introduces an analytical, non-iterative wavefront sensing framework based on spatially shifted self-reference interferometry, proving that consecutive co-prime shifts create a connected measurement graph with minimal hop diameter to bound worst-case error accumulation, and combining hop-constrained path averaging with closed-form Fourier least-squares refinement to achieve fast, provably robust phase recovery from as few as 8 to 16 measurements.
Background & Motivation¶
Wavefront sensing aims to recover the full complex-valued optical field (both intensity and phase), serving as a foundational enabler across adaptive optics, quantitative phase microscopy, imaging through scattering media, and digital holography. However, conventional photoelectric sensors respond strictly to time-averaged optical intensity, entirely discarding phase information. Consequently, recovering phase from intensity-only measurements represents a notoriously ill-posed, non-convex inverse problem. Prevailing computational phase retrieval approaches rely on alternating projection heuristics (e.g., Gerchberg-Saxton variants) or non-convex gradient descent optimization; these methods remain fundamentally plagued by high sensitivity to initial guesses, protracted convergence times, and susceptibility to local minima. While deep-learning priors and plug-and-play denoisers have improved visual fidelity, they operate as black boxes without rigorous error bounds or physical provability.
Phase-shifting interferometry (PSI) can yield analytical, closed-form reconstructions by superimposing the unknown wavefront with a calibrated reference beam. Nevertheless, separate-path reference beams are exceedingly fragile, suffering severe fringe degradation and phase-step errors in the presence of mechanical vibrations, thermal drift, or in vivo biological fluctuations. Common-path architectures such as point diffraction interferometry (PDI) spatially filter a bright spot of the incident wave through a sub-micron pinhole to form a self-reference; however, the extreme spatial filtering severely throttles reference beam irradiance, drastically diminishing the measurement signal-to-noise ratio (SNR). Lateral shearing interferometry (LSI) avoids pinhole attenuation by interfering the beam with its laterally shifted replica, but standard numerical difference integration is hypersensitive to noise, suffering from linear error accumulation over long pixel baselines while lacking theoretical principles for optimal shift selection.
This work attacks the problem from a fresh angle: by modeling self-interference between an incident wave and its spatially shifted replicas as a discrete phasor difference integration on a graph, and co-designing the measurement shift vectors using number-theoretic and graph-topological principles. Core Idea: employ consecutive co-prime shifts to construct a fully connected pixel graph that provably achieves the theoretical lower bound on maximum hop distance to cap worst-case error accumulation, coupled with non-iterative hop-constrained phasor propagation and closed-form Fourier least-squares refinement for high-throughput, provably robust wavefront reconstruction.
Method¶
Overall Architecture¶
The input consists of the unknown complex wavefront \(x = |x| \odot e^{j\phi}\) interfered with its spatially shifted replicas under four-step quadrature phase stepping \(\varphi_q\) on a spatial light modulator (SLM), yielding intensity measurements \(y_{k,q}\). The output is the globally reconstructed, high-fidelity 2D phase profile \(\hat{\phi}_{LS}\). Bypassing iterative forward-adjoint operator evaluations entirely, the reconstruction pipeline executes across four streamlined phases: self-referenced phasor difference extraction, co-prime shift graph construction, hop-constrained phasor propagation, and Fourier-domain global least-squares refinement.
%%{init: {'flowchart': {'rankSpacing': 24, 'nodeSpacing': 28, 'padding': 6, 'wrappingWidth': 400}}}%%
flowchart TD
In["Input Incident Wavefront & Interferograms<br/>x, y_{k,q}"] --> D1["Self-Referenced Phasor Difference Extraction<br/>Quadrature stepping to isolate relative phasor p_k"]
D1 --> D2["Co-Prime Shift Graph Construction<br/>Consecutive co-prime shifts s, t minimizing diameter h*"]
D2 --> D3["Hop-Constrained Phasor Propagation<br/>BFS traversal with equal-hop phasor averaging"]
D3 --> D4["Fourier Global Least-Squares Refinement<br/>Integer cycle unwrapping & closed-form FFT solver"]
D4 --> Out["Output High-Fidelity Wavefront Phase<br/>phi_{LS}"]
Key Designs¶
1. Self-Referenced Phasor Difference Extraction: Bypassing External Beams and Decoupling Intensity Variations Standard common-path point diffraction setups suffer from catastrophic SNR drops due to pinhole filtering. Here, a full-aperture self-referenced architecture is formulated where the incident wavefront \(x\) interferes directly with its shifted replica under shift operator \(\mathcal{S}_{\Delta_k}\). For each 2D shift vector \(\Delta_k\), four quadrature phase steps \(\varphi_q = q\frac{\pi}{2}\) (\(q \in \{0,1,2,3\}\)) are introduced, yielding intensity captures: $\(y_{k,q} = |x + e^{j\varphi_q}\mathcal{S}_{\Delta_k}(x)|^2 + \eta_{k,q} = |x|^2 + |\mathcal{S}_{\Delta_k}(x)|^2 + 2\text{Re}\left(|x| \odot |\mathcal{S}_{\Delta_k}(x)| e^{j(\phi - \mathcal{S}_{\Delta_k}(\phi) - \varphi_q)}\right) + \eta_{k,q}\)$ Taking the quadrature-weighted linear combination \(\sum_q y_{k,q} e^{-j\varphi_q}\) perfectly cancels out the static zero-frequency background and auto-intensity components, isolating the complex mutual interference term \(4 |x| \odot |\mathcal{S}_{\Delta_k}(x)| e^{j(\mathcal{S}_{\Delta_k}(\phi) - \phi)}\). Normalizing this quantity to unit magnitude strips away spatially non-uniform amplitude modulation, analytically isolating the unit phase-difference phasor map \(p_k(i) \approx e^{j(\phi(i+\Delta_k) - \phi(i))}\). This extraction requires zero external optical paths, rendering the system impervious to external vibrations while preserving maximum optical throughput.
2. Co-Prime Shift Graph Construction: Guaranteeing Full Connectivity and Minimum Hop Bounds In conventional shearing interferometry or single-pixel neighbor difference schemes (\(\Delta=1\)), the pixel interaction topology degenerates into a linear chain graph. Pixels remote from the reference node require \(O(N)\) hops, causing independent edge noise variances to accumulate and diverge across long baselines. This work formalizes shift-based measurements over a 1D grid \(V\) (extended to 2D via tensor product) with edge set \(E = \{(i, i+s_k)\}\). The authors prove that whenever two shifts \(s, t \in \mathbb{N}\) are co-prime (\(\gcd(s,t)=1\)) and satisfy \(s+t \le N\), full graph reachability from a central reference pixel is guaranteed via modulo \((s+t)\) residue coverage and sliding window overlap. Crucially, to prevent noise explosion along traversal paths, the paper establishes a universal lower bound on the maximum hop distance \(h^\star\): $\(h^\star \ge \left\lceil \frac{-1 + \sqrt{2N - 1}}{2} \right\rceil\)$ The authors rigorously prove that selecting consecutive co-prime integers \(s = \lfloor \sqrt{N/2} \rfloor\) and \(t = \lfloor \sqrt{N/2} \rfloor + 1\) achieves this theoretical lower bound. For an image dimension of \(N=512\), optimal shifts \(s=16, t=17\) collapse the maximum hop distance from 256 hops down to merely 16 hops, provably capping worst-case error accumulation.
3. Hop-Constrained Phasor Propagation: Suppressing Variance via Non-Iterative Wavefront Expansion Integrating differential phases in the complex phasor domain avoids trigonometric branch cuts and wraps. A breadth-first search (BFS) queue is initiated at the reference pixel with known phase \(\hat{p}(i_0) = 1\). As the expansion reaches candidate pixels through multiple distinct paths, naive averaging across all arriving paths would corrupt the estimate, because paths with longer hop counts carry substantially higher accumulated noise variance. To counter this, a hop-constrained path averaging rule is enforced: only paths that reach a node with identical minimal hop length \(h(i') = h(i) + 1\) are averaged vectorially; any longer path arriving subsequently is discarded outright. This maintains an unbiased phase estimate while driving down estimation variance, and executes with high efficiency via vectorized GPU tensor operations.
4. Fourier Global Least-Squares Refinement: Closed-Form Cycle Correction across Redundant Loops While tree-based BFS propagation is exceptionally fast, it leaves cross-loop redundant constraints underutilized. Representing the forward spatial difference operator for shift \(\Delta_k\) as \(D_{\Delta_k} = \mathcal{S}_{\Delta_k} - I\), the preliminary propagation phase estimate \(\hat{\phi}\) is first leveraged to compute exact integer wrapping jump offsets: $\(m_k = \text{round}\left( \frac{D_{\Delta_k}\hat{\phi} - \nabla_k\phi}{2\pi} \right)\)$ This unwraps raw periodic differences \(\nabla_k\phi = \angle p_k \in (-\pi, \pi]\) into continuous offsets \(\nabla_k\phi^\text{unwrap} = \nabla_k\phi + 2\pi m_k\). A global least-squares objective \(\min_\phi \sum_k \|D_{\Delta_k}\phi - \nabla_k\phi^\text{unwrap}\|_2^2\) is formulated over all measurements. Because \(D_{\Delta_k}\) corresponds to a spatial convolution kernel \(d_{\Delta_k}\), the entire least-squares system is diagonalized in the 2D spatial frequency domain, producing a closed-form analytical solution: $\(\hat{\phi}_{LS} = \mathcal{F}^{-1}\left( \frac{\sum_k \overline{\mathcal{F}(d_{\Delta_k})} \odot \mathcal{F}(\nabla_k\phi^\text{unwrap})}{\sum_k |\mathcal{F}(d_{\Delta_k})|^2 + \lambda} \right)\)$ where \(\mathcal{F}\) denotes the 2D discrete Fourier transform and \(\lambda\) is a regularizer. Solvable in a fraction of a second via 2D FFTs, this refinement eliminates residual traversal streaks and pushes reconstruction accuracy to the theoretical limit.
Key Experimental Results¶
Main Results¶
Simulations were evaluated on 100 intensity patterns from the DIV2K validation dataset, combined with three distinct phase profiles: quadratic, random, and smooth (MATLAB peaks) profiles under combined Poisson and Gaussian noise at 22 dB SNR. Baselines include gradient descent with random initialization (GD-Rand), spectral initialization (GD-Spec), plug-and-play FISTA with DRUNet prior (PnP-FISTA), wavefront imaging sensor with high resolution (WISH with phase ramps), and Deep Image Prior (DIP).
| Method | Quadratic (16 Meas.) | Quadratic (32 Meas.) | Random (16 Meas.) | Random (32 Meas.) | Smooth (16 Meas.) | Smooth (32 Meas.) |
|---|---|---|---|---|---|---|
| GD-Rand | 0.695 ± 0.365 | 0.470 ± 0.434 | 0.138 ± 0.148 | 0.079 ± 0.152 | 0.791 ± 0.361 | 0.570 ± 0.404 |
| GD-Spec | 0.737 ± 0.400 | 0.400 ± 0.432 | 0.691 ± 0.367 | 0.449 ± 0.440 | 0.686 ± 0.397 | 0.371 ± 0.374 |
| WISH | 0.852 ± 0.299 | 0.423 ± 0.366 | 0.854 ± 0.320 | 0.453 ± 0.377 | 0.898 ± 0.330 | 0.425 ± 0.367 |
| PnP-FISTA | 0.669 ± 0.341 | 0.456 ± 0.304 | 0.735 ± 0.329 | 0.425 ± 0.264 | 0.663 ± 0.341 | 0.479 ± 0.338 |
| DIP | 1.486 ± 0.031 | 1.491 ± 0.051 | 1.495 ± 0.037 | 1.447 ± 0.065 | 1.442 ± 0.063 | 1.423 ± 0.081 |
| Ours | 0.158 ± 0.104 | 0.099 ± 0.083 | 0.156 ± 0.099 | 0.094 ± 0.075 | 0.155 ± 0.096 | 0.097 ± 0.071 |
| Ours with LS | 0.183 ± 0.082 | 0.126 ± 0.065 | 0.158 ± 0.075 | 0.121 ± 0.055 | 0.205 ± 0.063 | 0.159 ± 0.050 |
Ablation Study¶
The impact of shift selection on graph topology, hop diameter, and phase error accumulation was investigated on an \(N=512\) grid, alongside execution runtime and speedup comparisons on an NVIDIA RTX 6000 Ada GPU.
| Shift Config / Method | Max Hops \(h\) | Graph Connectivity | Error Accumulation Behavior | Runtime (s) | Relative Speedup |
|---|---|---|---|---|---|
| Single Shift (\(\Delta_1=1\)) | 256 | Linear chain topology | Linear growth with coordinate distance | - | - |
| Non-Co-Prime (e.g., 16, 18) | \(\infty\) | Disconnected (2 isolated components) | Disconnected blind zones; failure | - | - |
| Suboptimal Co-Prime (e.g., 2, 3) | 86 | Sparse long-reach connectivity | Large loop diameters, noise amplification | - | - |
| Optimal Co-Prime (16, 17) | 16 | Fully connected (achieves lower bound \(h^\star\)) | Tightly concentrated hops, bounded error | - | - |
| PnP-FISTA Iterative | - | - | Heavy deep denoiser iteration | 10.7 | 1.0× (Baseline) |
| DIP Neural Fitting | - | - | Test-time gradient backprop | 8.7 | 1.2× |
| WISH Alternating Proj. | - | - | Multi-round spatial-Fourier projection | 3.6 | 3.0× |
| Gradient Descent (GD) | - | - | Non-convex loss optimization | 1.5 | 7.1× |
| Ours / Ours with LS | - | - | Non-iterative propagation + closed FFT | 1.2 | 8.9× |
Key Findings¶
- Under the tight budget of 16 measurements (2 optimal shifts \(\times\) 2 axes \(\times\) 4 phase steps), all prior iterative and learning baselines experience severe phase degradation (mean errors \(>0.60\) rad on smooth and quadratic wavefronts), whereas the proposed method reconstructs fine phase details with low error (\(\approx 0.155\) rad).
- Across SNR levels from 10 dB to 30 dB, the least-squares refinement step (Ours + LS) delivers dramatic robustness gains at lower SNRs, successfully suppressing random phase deviations.
- Hardware bench experiments confirmed physical fidelity: recovering the continuous concentric profile of a biconvex lens (\(f=75\text{ mm}\)), enabling purely numerical refocusing via Fresnel propagation to resolve a sharp star target at \(z=67\text{ mm}\), and peering through a diffuser (Thorlabs DG10-120-A) to restore previously obscured target objects.
Highlights & Insights¶
- Mapping physical interferometric sampling onto an optimal graph diameter problem bridges optical physics, number theory, and discrete algorithms, offering a provable paradigm that departs from black-box heuristics.
- Selecting consecutive co-prime shifts \((s, s+1)\) near \(\sqrt{N/2}\) achieves an optimal balance between step reach and local density, shrinking graph diameter from \(O(N)\) to \(O(\sqrt{N})\) and locking error variance into a narrow bound.
- The hop-constrained averaging principle provides a valuable general insight: when aggregating noisy relational measurements, rejecting longer paths entirely rather than averaging them prevents error leakage from degraded estimates.
Limitations & Future Work¶
- The current algorithm applies uniform weights to all paths sharing equal hop counts; incorporating quality-guided metrics (e.g., local fringe visibility or intensity-weighted confidence) could further enhance robustness under localized shadowing or high-density speckle.
- On uncorrelated random phase maps, where phase gradients exhibit sharp discontinuities without spatial smoothness, the advantages of graph propagation over unconstrained gradient descent become less pronounced.
- The hardware prototype currently implements phase stepping and shifts sequentially via an SLM; adapting the framework to single-shot snapshot architectures (such as polarization-division multiplexing) represents an important direction for dynamic adaptive optics and in vivo imaging.
Related Work & Insights¶
- vs Lateral Shearing Interferometry (LSI): Classical LSI relies on small fixed shears and numerical gradient integration, which suffers from boundary divergence and cumulative noise; this method proves co-prime shear optimality, shrinking maximum hops to \(O(\sqrt{N})\) with analytical recovery.
- vs WISH & Modulated Phase Retrieval: WISH requires random phase masks and dozens of alternating projection iterations; the proposed framework replaces iterative projections with analytical graph traversal, cutting reconstruction latency from 3.6 s to 1.2 s.
- vs Deep Learning & Plug-and-Play (PnP / DIP): PnP-FISTA and DIP demand heavy computation (8–10 s) and are vulnerable to hallucination or smoothing artifacts under scarce measurements; the proposed approach is fully deterministic, computationally lightweight, and requires no training data.
Rating¶
- Novelty: ⭐⭐⭐⭐⭐ Elegant synthesis of number theory, graph shortest paths, and self-reference interferometry with formal proofs.
- Experimental Thoroughness: ⭐⭐⭐⭐⭐ Comprehensive coverage from mathematical bounds and synthetic benchmarks to a physical 4f optical prototype demonstrating lens sensing, auto-refocusing, and scattering medium descattering.
- Writing Quality: ⭐⭐⭐⭐⭐ Rigorous lemma and theorem statements, clean narrative structure, and compelling experimental visualizations.
- Value: ⭐⭐⭐⭐⭐ Delivers an ultra-fast, reference-beam-free wavefront sensing solution with immediate utility across biological microscopy and adaptive optics.