Phase unwrapping estimates the missing multiples of 2π in measured phase images. For large images, tiling limits the size of local reconstruction problems and enables parallel processing. Adaptive tiling could further reduce the number of local problems and boundaries by retaining large tiles where little refinement is needed. We investigate whether this reduction makes reconstruction faster. We compare complete reconstruction time and accuracy for a regular grid, quadtree, and kd-tree partitions. We also evaluate nine criteria for deciding where quadtree tiles should be subdivided, including residue count, fringe density, and measures of phase variation, at different tile sizes and budgets. In single-threaded experiments on a heterogeneous image dataset, optimized adaptive partitions use fewer tiles but remain slower than the optimized grid, and some reconstructions lose substantial accuracy. Stage measurements explain why: constructing the partition and solving larger retained tiles outweigh the savings at tile boundaries. The criterion comparison also shows that more refinement does not consistently improve accuracy. These results motivate evaluating adaptive partitions by the complete time needed to reach a chosen reconstruction accuracy, including whether limited refinement can provide a faster approximate result.
Figures & tables
Figure 1: Regular and adaptive partitions on three 256×256 phase fields: A, a generated interferometric synthetic aperture radar (InSAR) field; B, a field derived from tree-photograph intensities; C, a generated vortex field. Columns show the known generating phase (ground truth, GT), wrapped phase, signed cell residues (red: +1 , blue: −1 , white: zero), a grid of 8×8 -pixel tiles, and a residue-driven quadtree. Each GT panel has its own color scale in radians; wrapped panels share the cyclic scale [−π,π] . White lines mark tile boundaries. The quadtree has minimum side s=8 pixels and stopping budget B=512 leaves. It retains 268, 514, and four tiles in A–C, versus 1024 grid tiles in each case. A four-way split adds three leaves, allowing the budget to be exceeded by two. These unquantized examples illustrate allocation, not reconstruction accuracy or runtime.
Figure 2: Common pipeline for grid, quadtree, and kd-tree reconstruction. The partition sets the tile boundaries; local Poisson solves use the BRiDCT discrete cosine transform (DCT), with FFTW for unsupported shapes. Seam discrepancies are summarized by their median and median absolute deviation (MAD). A separate reliability tree connects the reconstructed tiles and determines their relative phase offsets. The complete-call timer includes preparation, all reconstruction stages, and temporary-buffer cleanup. Transform plans are reused in the main timing comparison; first-call latency is evaluated separately.
Figure 3: Priority-based quadtree refinement on a 32×32 image, with minimum tile side s=8 pixels and stopping budget B=7 leaves. C is an illustrative nonnegative tile score, and L is the leaf count. (a) The root covers the image. (b) Bisecting both sides creates four 16×16 children; the orange outline selects the highest-scoring eligible child. (c) Splitting that child creates four 8×8 tiles and brings the total to seven leaves. Refinement stops at the budget even though another tile remains eligible with positive score. Zero-score tiles are retained.
Family
Entries
Source and reference
Photographs
45
Image intensities mapped to known phase fields, including USC-SIPI images.
InSAR
69
19 generated fields; 50 InSAR-DLPU entries, including 30 with terrain-model references [ 12 ] .
Synthetic surfaces
18
Analytic or generated fields with known phase.
Published tests
4
Ghiglia–Pritt images and supplied reference arrays [ 5 ] .
Holography
40
20 synthetic and 20 experimental entries; the latter use supplied algorithmic wrap-count references [ 13 ] .
MRI simulations
32
Magnetic resonance imaging (MRI) simulations from the quantitative susceptibility mapping (QSM) challenge, with known phase and brain masks [ 14 ] .
Table 1: Composition of the 208-entry dataset. Counts denote test inputs, not independent acquisitions. Ten entries guide implementation choices; the other 198 provide the main comparison. Reference phases comprise 158 known generating fields, 30 terrain-model references, and 20 algorithmic references. Agreement with the last category does not establish agreement with independent ground truth.
s (px)
Partition
Median t/tG
∑t/∑tG
ΔSπ (pp)
Δ RMSE (rad)
8
Quadtree
1.244
1.371
-0.20
+0.105
8
KD midpoint
1.228
1.381
-0.04
+0.133
8
KD dyadic
1.229
1.294
-0.06
+0.145
16
Quadtree
1.532
1.777
-1.02
+0.184
16
KD midpoint
1.578
1.796
-0.99
+0.213
16
KD dyadic
1.570
1.644
-1.01
+0.225
Table 2: Complete reconstruction time and reference agreement on 198 entries. s is the minimum adaptive tile side in pixels; the grid uses nominal s×s tiles. For each entry, t and tG are the median times over seven repeats for the adaptive method and its grid control. The median of t/tG describes a typical entry, while ∑t/∑tG compares aggregate time across entries. Ratios above one mean slower reconstruction. ΔSπ and Δ RMSE are mean paired changes, adaptive minus grid, in percentage points (pp) and radians. Sπ is the fraction of pixels with absolute error below π after removing a common mean phase offset. Thus positive ΔSπ and negative Δ RMSE indicate better reference agreement. All methods share BRiDCT, FFTW fallback, and reconciliation.
Figure 4: Complete reconstruction time relative to the grid and its stage breakdown. (a) Each adaptive method uses minimum tile side s and is paired with a grid of nominal side s on the same input. Per-entry times are medians over seven repeats. Dots show the median ratio across 198 entries; whiskers show the 10th–90th percentiles across entries, not confidence intervals. The dashed line at one marks equal time. (b) Stage times are summed across entries using their per-entry medians and divided by the summed median complete-call times. Partitioning includes residue computation and tree construction; local solves include data preparation, the Poisson right-hand side, and transforms. Boundary statistics compute seam medians and MAD; joining propagates offsets. The remainder includes cleanup, uninstrumented overhead, and differences caused by summarizing each stage separately.
Figure 5: Effect of the subdivision criterion on a fixed observation. All panels show the same unquantized 256×256 tree-photograph phase field as row B of Fig. 1 , with white quadtree boundaries over the wrapped phase. Panels (a)–(i) use criteria A–I from Sec. 4 ; L is the resulting number of tiles. Minimum side s=8 pixels and stopping budget B=512 leaves are identical throughout. A four-way split can exceed this budget by two leaves. The shared cyclic color scale gives phase in radians. These partitions illustrate how scores allocate a limited tile budget; they are not the unrestricted partitions used in the main runtime comparison.
Full budget B1
Half budget B1/2
Criterion
L
ΔSπ (pp)
Δ RMSE (rad)
L
ΔSπ (pp)
Δ RMSE (rad)
Minimum side s=8 pixels
A. Residue count
302.5
+0.00
+0.000
302.5
+0.00
+0.000
B. Local net charge
419.5
-0.25
+0.266
419.5
+0.50
+0.236
C. Paired cut length
362.5
-0.42
+0.026
362.5
-0.18
+0.070
D. Fringe density
748
+0.07
+0.021
514
+0.37
-0.226
Table 3: Criterion comparison on 208 entries at minimum side s pixels. For original dimensions H×W , the full budget is B1=max(1,⌊H/s⌋⌊W/s⌋) and the half budget is B1/2=max(1,⌊B1/2⌋) . For a 256×256 image at s=8 , these targets are 1024 and 512 tiles; a split may exceed the target by two. L is the median tile count. ΔSπ (percentage points, pp) and Δ RMSE (radians, rad) are mean paired changes relative to residue count at the same s and budget. Positive ΔSπ and negative Δ RMSE improve agreement. The zero control rows do not imply equal accuracy between budgets: absolute control means are given below. Budgets limit tiles, not error, and this table does not compare runtimes.
Family
Entries
Dimensions
Source and reference phase
Photograph-derived
45
256×256 , 512×512 , 1024×1024 ; 250×288
Grayscale test images, including USC-SIPI images, scaled into phase fields; the generating field is the reference. Clean and noisy variants are included.
InSAR
69
256×256
19 locally generated interferograms with known phase; 20 simulated and 30 real TanDEM-X entries from InSAR-DLPU [ 12 ] . The real-data references derive from a digital elevation model.
Synthetic surfaces
18
256×256 , 257×257 ; 458×152
Analytic and generated phase fields, with the generating field retained as the reference.
Published test images
4
257×257 ; 458×157 ; 458×152
Ghiglia–Pritt test images [ 5 ] , using their supplied reference arrays.
Holography
40
256×256
20 synthetic fields with known phase and 20 experimental inputs from Gontarz et al. [ 13 ] . Experimental references are reconstructed from the supplied quality-guided phase-unwrapping (QGPU) wrap counts.
MRI simulations
32
92 – 167 rows; 101 – 155 columns
Cropped slices from four QSM challenge simulation/noise combinations [ 14 ] , two echo times per view. Reference phase is 2πfTE from the supplied frequency f and echo time TE .
Table S1: Dataset families, original image dimensions, and reference-phase construction. Dimensions are height × width in pixels before square padding; ranges cover the selected entries, not the entire source collections. Known generating references include both analytic fields and phases derived from photograph intensities. Terrain-model and algorithmic references are identified separately in the last column.
ID
Input
Size
Origin and role in the comparison
F1
Tree photograph
256×256
Source recorded as USC-SIPI; intensity-derived phase with fine texture.
F2
Portrait
1024×1024
Source recorded as USC-SIPI; larger intensity-derived phase image.
Photograph-derived field retained in the local collection; intermediate subdivision.
Table S2: Ten inputs used for implementation choices and first-call measurements. Dimensions are original height × width in pixels, before padding. F1–F10 identify entries within this report. The last column states their origin and the image property relevant to these tests. Source-collection citations identify provenance; the exact derived benchmark arrays are not publicly available.
Figure S1: Cell-score maps underlying the tree-image partitions in the main report’s Fig. 5 . Panels (a)–(i) correspond to criteria A–I. Color represents ac/amax , where ac is the cell score and amax is the maximum within that panel; a zero map would be displayed as zero. Independent normalization reveals each map’s spatial support but prevents comparison of absolute scores between criteria. Panel labels report the proportion or count of nonzero cells. The uniform control is one everywhere; its inset shows the resulting partition, not a second score map.
Side (px)
SJ baseline
FFTW
Accelerate
Ooura general
Ooura specialized
8
0.048
0.131
—
0.530
0.101
16
0.189
1.189
0.956
1.712
0.489
32
0.882
4.374
3.350
6.451
—
64
6.141
17.668
11.417
25.823
—
128
67.924
94.720
44.160
104.776
—
256
296.262
468.333
186.678
473.060
—
Table S3: Poisson spectral-solve time in μ s per square tile on an Apple M3 Max, using one thread and single precision. Each value is the median of nine randomly ordered timing batches cycling through 16 inputs. The operation comprises a forward DCT, division by nonzero Poisson eigenvalues, and an inverse DCT, including normalization and data rearrangement; setup, allocation, and right-hand-side construction are excluded. SJ baseline is the Shao–Johnson implementation preceding BRiDCT, not a BRiDCT measurement. Bold marks the fastest tested implementation at each size. Dashes indicate sizes unsupported by the tested interface, not by the entire algorithm family. Smaller values are faster; these kernel timings are not complete reconstruction times.
Variant
Time ratio
ΔSπ (pp)
Δ RMSE (rad)
s=8
s=16
s=8
s=16
s=8
s=16
Density 0.001
1.027
1.016
-5.93
-0.78
+0.240
+0.090
Density 0.01
1.061
1.124
-1.52
-1.30
+0.556
+0.016
Low-charge cutoff
0.991
1.023
-1.44
+0.46
+0.873
+0.106
Verified integration
1.039
1.012
-0.04
0.00
+0.002
0.000
Residue-balanced kd-tree
1.525
1.847
-0.70
-3.07
+0.420
-0.190
Table S4: Paired effects of implementation variants on ten inputs at minimum side s pixels. Time ratios compare each variant with the control specified below, using per-input median timings; the table reports the median ratio across inputs. Values above one are slower. ΔSπ and Δ RMSE are mean paired changes, variant minus control, in percentage points (pp) and radians. Sπ is the fraction of pixels with absolute reference error below π after mean-offset removal. Positive ΔSπ and negative Δ RMSE indicate improvement. N is the processed pixel count. Rows from different timing sessions do not establish a ranking between variants.
s (px)
Partition
First/warm
First/grid
8
Grid
2.245
1.000
8
Quadtree
2.151
1.219
8
KD dyadic
2.244
1.244
16
Grid
2.938
1.000
16
Quadtree
2.391
1.378
16
KD dyadic
2.459
1.373
Table S5: First-call overhead and adaptive/grid comparison on the ten exploratory inputs. s is the minimum adaptive side, or nominal grid side, in pixels. Each entry/configuration uses three fresh processes, each with one first call and seven warm calls. First-call time is the median over the three processes; warm time is the median of their seven-call medians. Columns report medians of the resulting per-entry ratios. First/warm compares the same method before and after initialization; first/grid compares its first call with the grid’s first call at the same s . Ratios above one mean longer reconstruction. Process startup, input loading, and encoding are excluded.
Phase unwrapping is a key step in interferometric and coherent imaging, where the physical quantity of interest is carried by a phase that the instrument delivers only modulo 2*pi. The difficulty in two dimensions is to separate the jumps caused by wrapping from those produced by noise, by true discontinuities, by under-sampling or by decorrelation. Spatial-domain and frequency-domain methods have both been studied extensively, each with advantages the other lacks; hybrid schemes combining the two remain scarce. We propose a semi-global tile-based strategy in which every tile is unwrapped in the frequency domain, through the Discrete Cosine Transform (DCT) and the least squares (LS) formalism of Ghiglia et al., the tiles being merged spatially. Unwrapping and noise filtering are performed jointly, which regularizes an otherwise ill-posed inverse problem. An error stays confined to the tile in which it arose; the tiling artifacts are removed by averaging over every shift of the grid and over the symmetries of the square; and that average is weighted by the Poisson residual each pass leaves behind, so that a pass whose tile boundaries fell on a discontinuity does not impose its seam on the result. Experiments on synthetic and real data, against four reference algorithms from four distinct families and under six complementary metrics, show that the proposed method matches or improves on the state of the art. One of these metrics, a corrected cyclic re-wrap residual introduced here, needs no ground truth and therefore remains available on real acquisitions.
Antoine Moevus, Max Mignotte
D´epartement d’informatique et de recherche opérationnelle, Universit´e de Montr´eal, Montr´eal, QC, Canada
Deep unrolling methods for inverse problems augment traditional reconstruction algorithms with learned priors to improve their performance. These methods have been developed under the implicit assumption that the geometry of the measurement operator does not change between training and deployment; when this is no longer the case, the quality of obtained reconstructions can degrade significantly. In this paper, we argue that this is not an inherent limitation of the deep unrolling paradigm and can be mitigated by a careful choice of the underlying iterative algorithm. We propose a deep unrolling framework for the phase retrieval problem that augments the alternating projections method with a learned denoising step. Our method obtains high-quality reconstructions without the need for retraining in settings where conventional model-based alternatives struggle, such as when the dynamic range of the measurement operator changes significantly at test time.
Electrical & Computer Engineering, University of Washington · NSF-Simons National Institute for Theory and Mathematics in Biology · Department of Computer Science, University of Chicago +3
Fourier phase encodes important spatial image structure, but recovering an image without measured spectral magnitude requires additional constraints and leaves absolute intensity ambiguous. We propose a projection-based two-stage framework that combines Fourier-phase and spatial-support constraints with an image-specific neural prior. The first stage alternates constraint enforcement with regularized neural-prior updates, while the second performs phase/support refinement alone with guaranteed convergence. We evaluate two neural-prior implementations on the same 77 microscopy images and compare them with a constraint-only baseline. After 500 final refinement passes, the best-performing variant achieves 31.41 dB pooled PSNR, 35.75 dB mean PSNR, and 0.9531 mean SSIM, improving pooled PSNR by 1.51~dB and reducing pooled MSE by 29.3% relative to the baseline. The results demonstrate the benefit of combining neural guidance with explicit constraint refinement at the evaluated iteration budget, while showing that lower phase residual alone does not guarantee greater reconstruction accuracy.
Ene Meco, Ahmet Enis Cetin
Dept. of Electrical and Computer Engineering, University of Illinois at Chicago, IL, 60607, USA