PaperlayerView exact revision on arXiv

Source-linked AI summary

Recovery of Piecewise Smooth Images from Few Fourier Samples

Greg Ongie, Mathews Jacob

arXiv:1502.00705v1cs.CV

TL;DR

Recovering continuous 2-D piecewise-smooth images from few Fourier samples is difficult for existing methods because image derivatives are supported on curves. This paper extends annihilation to such images and uses structured low-rank matrix completion to jointly recover edges and image amplitudes. It also gives unique edge-recovery conditions for piecewise-constant images and demonstrates MRI phantom recovery from low-resolution Fourier samples.

  • Problem

    Existing Prony-like and off-the-grid methods do not straightforwardly handle 2-D piecewise-smooth images, whose derivatives are supported on separating curves.

    Source: §I · p. 1

  • Method

    The method uses Fourier-domain annihilation for piecewise-smooth signals and a single-step convex structured low-rank matrix completion algorithm to jointly recover the edge set and image amplitudes.

    Source: §I · p. 1, 1) First order partial derivative… · p. 2

  • Results

    The paper derives necessary and sufficient Fourier-sample conditions for unique edge recovery in piecewise-constant images and demonstrates MRI phantom recovery from few low-resolution Fourier samples.

    Source: §IV · p. 3, B. Sufficient conditions · p. 3, §VI · p. 5

  • Takeaways & Limitations

    The annihilation relation supports exact continuous edge-set recovery from few Fourier samples, with the paper demonstrating accurate recovery of MRI phantoms.

    Source: §I · p. 1, §VI · p. 5

    The necessary and sufficient unique-recovery conditions are derived for piecewise-constant signals.

    Source: §IV · p. 3

Abstract

from arXiv · show

We introduce a Prony-like method to recover a continuous domain 2-D piecewise smooth image from few of its Fourier samples. Assuming the discontinuity set of the image is localized to the zero level-set of a trigonometric polynomial, we show the Fourier transform coefficients of partial derivatives of the signal satisfy an annihilation relation. We present necessary and sufficient conditions for unique recovery of piecewise constant images using the above annihilation relation. We pose the recovery of the Fourier coefficients of the signal from the measurements as a convex matrix completion algorithm, which relies on the lifting of the Fourier data to a structured low-rank matrix; this approach jointly estimates the signal and the annihilating filter. Finally, we demonstrate our algorithm on the recovery of MRI phantoms from few low-resolution Fourier samples.

I. INTRODUCTION

Recovering 2-D piecewise-smooth images from few Fourier samples is not a straightforward extension of existing harmonic-retrieval methods, which are not designed for derivatives supported on separating curves. The paper models these images directly and combines edge recovery with image-amplitude recovery in a single convex method.

I. INTRODUCTION

  • Existing Prony-like and off-the-grid methods do not straightforwardly extend to 2-D piecewise-smooth images, whose derivatives are supported on curves separating regions.Those derivatives can be viewed as a continuum of Diracs, which current methods are not designed to handle.Source: §I · p. 1
  • The paper proposes a piecewise-smooth signal model and generalizes annihilation to recover an exact continuous representation of the edge set from few Fourier samples.For piecewise-constant images, it gives necessary and sufficient conditions on the number of samples for perfect edge-set recovery.Source: §I · p. 1
  • A single-step convex structured low-rank matrix completion algorithm jointly estimates the edge set and image amplitudes from Fourier data.The approach lifts ideal Fourier samples to a structured low-rank matrix and is demonstrated on MR phantoms using low-resolution Fourier samples.Source: §I · p. 1

1) Piecewise constant images:

The signal model includes piecewise-constant, piecewise-analytic, and piecewise-harmonic images through choices of differential operators that vanish within each region.

1) Piecewise constant images:

  • For piecewise-constant images, choosing operators whose null space consists of constants reduces the model to the piecewise-constant image model, the paper's primary focus for provable guarantees.The paper identifies this case as its primary focus because of its simplicity and provable guarantees.Source: 1) Piecewise constant images · p. 2

2) Piecewise analytic images:

  • Choosing the complex derivative operator makes the model piecewise analytic, equivalent to the model in, which the paper describes as unrealistic for natural images..Source: 2) Piecewise analytic images · p. 2

3) Piecewise harmonic:

4) Piecewise linear images:

Higher-order differential operators extend the signal model to piecewise-polynomial images, while trigonometric polynomials describe the edge curves and support uniqueness results.

5) Piecewise polynomial:

  • Using all nth-order partial derivatives yields pieces that are polynomials of degree at most n−1.The paper uses derivative annihilation to recover an exact continuous edge-set representation from finitely many Fourier samples.Source: 5) Piecewise polynomial · p. 2

B. Trigonometric polynomials and curves

  • A trigonometric curve is the zero set of a band-limited periodic trigonometric polynomial, and each curve has a unique minimal-degree polynomial up to scaling.The minimal polynomial divides every other trigonometric polynomial defining the same curve.Source: B. Trigonometric polynomials and curves · p. 2
  • If a trigonometric polynomial and its gradient both vanish on a curve, the square of that curve's minimal polynomial divides it.The gradient of the minimal polynomial vanishes at only finitely many points on the curve.Source: B. Trigonometric polynomials and curves · p. 2

III. ANNIHILATION PROPERTY

For a region indicator, first derivatives are supported on its boundary in the distributional sense, enabling a trigonometric polynomial that vanishes there to annihilate their Fourier transforms. The relation extends to piecewise constant signals and piecewise smooth components annihilated by the relevant differential operator.

1) First order partial derivative operators:

  • The derivatives of χΩ can be interpreted as continuous streams of weighted Diracs supported on ∂Ω.The derivatives are defined in the distributional sense because χΩ is non-smooth at the boundary.Source: 1) First order partial derivative… · p. 2
  • A boundary curve C given by the zero set of µ makes the Fourier transform of Df vanish under convolution with µ’s Fourier coefficients for any first-order differential operator D.This applies to f = χΩ and its boundary ∂Ω.Source: 1) First order partial derivative… · p. 2
  • The annihilation relation extends to piecewise constant functions when µ vanishes on the union of their boundaries, and to piecewise smooth components g·χΩ when Dg = 0.It also extends to linear combinations of components satisfying these conditions.Source: 1) First order partial derivative… · p. 3

2) Second order partial derivative operators:

For second-order differential operators, squaring a trigonometric polynomial that vanishes on the boundary gives a sufficient annihilator. When the interior function does not vanish on the boundary, the squared minimal polynomial is also necessary in the stated divisibility sense.

2) Second order partial derivative operators:

  • For f = g·χΩ with Dg = 0, µ² is sufficient to annihilate Df when D is second-order and µ annihilates the partial derivatives of χΩ.The argument uses the product rule and the boundary-zero property of µ.Source: 2) Second order partial derivative… · p. 3
  • When g does not vanish on ∂Ω, every annihilating polynomial ν must be divisible by µ0², where µ0 is the minimal polynomial for ∂Ω.If deg(µ0) = (K, L), then deg(ν) ≥ (2K − 1, 2L − 1).Source: 2) Second order partial derivative… · p. 3

3) Partial derivative operators of arbitrary order:

For an nth-order differential operator D, the annihilator is the nth power of a trigonometric polynomial vanishing on the boundary, provided Dg = 0 within the region. The relation also holds for linear combinations of qualifying components.

3) Partial derivative operators of arbitrary order:

  • If f = g·χΩ, Dg = 0, and ∂Ω lies in the zero set of µ, then µ^n annihilates the Fourier transform of Df by convolution with the Fourier coefficients of µ^n.This holds for any nth-order differential operator D.Source: 3) Partial derivative operators of… · p. 3
  • The nth-order annihilation relation extends by linearity to sums of components gi·χΩi when Dgi = 0 and µ vanishes on every component boundary.Source: 3) Partial derivative operators of… · p. 3

IV. RECOVERY FROM FINITE FOURIER SAMPLES

The section derives a necessary Fourier-sampling condition for recovering the edge-set filter of piecewise constant signals. For square grids and a degree-(K, K) minimal polynomial, the bound corresponds to at least 1.71K × 1.71K samples; experiments suggest, but do not establish, sufficiency.

A. Necessary conditions

  • 2 · (K′ − K + 1) · (L′ − L + 1) ≥ K · L is necessary to recover an edge set with minimal-polynomial degree (K, L) from a rectangular Fourier-sampling grid.The system supplies at most 2 · (K′ − K + 1) · (L′ − L + 1) equations for at most K · L filter coefficients.Source: A. Necessary conditions · p. 3
  • For degree (K, K) and a square sampling region, the necessary bound requires at least 1.71K × 1.71K Fourier samples.Source: A. Necessary conditions · p. 3
  • Simulations indicate the necessary condition might also be sufficient for unique recovery, but the section states this as a hypothesis.Source: A. Necessary conditions · p. 3

B. Sufficient conditions

The section gives conditions for uniquely recovering edge sets from Fourier-sample equations, including a result for piecewise-constant images with nonintersecting boundaries. It then motivates convex recovery because recovering amplitudes by polynomial factorization is impractical.

B. Sufficient conditions

  • For a simply connected region, the minimal-polynomial coefficients are uniquely recoverable up to scaling as the only non-trivial solution to the stated equations.The proof shows that any other polynomial satisfying the equations vanishes on the boundary and, by degree considerations, is a scalar multiple of the minimal polynomial.Source: B. Sufficient conditions · p. 3
  • For a piecewise-constant image with nonintersecting boundaries, the product of the regions’ minimal polynomials—and thus their combined edge set—is uniquely recoverable up to scaling.The stated result applies when the characteristic functions do not intersect.Source: B. Sufficient conditions · p. 4
  • The sufficient-condition equations require Fourier samples throughout 3Λ, more than the necessary sample count; sharpening the result and extending it to intersecting boundaries remain conjectures.Here, mΛ denotes the set Λ dilated by a factor of m.Source: B. Sufficient conditions · pp. 3–4

V. RECOVERY ALGORITHMS

  • Once the edge set is known, solving for amplitudes through polynomial factorization is impractical, motivating a convex-optimization formulation.The stated obstacle is factoring a high-degree multivariate polynomial into irreducible factors.Source: §V · p. 4

A. Curve-aware recovery

The paper contrasts curve-aware weighted total variation, which assumes a known edge polynomial, with a Fourier-domain structured low-rank completion approach that jointly estimates the image and filter. The latter is demonstrated on an undersampled Shepp–Logan phantom.

A. Curve-aware recovery

  • With the annihilating polynomial known, weighted total variation leaves image gradients along its zero-set curve unpenalized, allowing sharp edges there.The formulation constrains the candidate image to match the measured Fourier samples.Source: A. Curve-aware recovery · p. 4
  • The curve-aware scheme requires discretization, optimization of many parameters, and is sensitive to the estimated polynomial; the alternative jointly recovers the image and polynomial in one stage.The joint method does not rely on an explicit estimate of the annihilating polynomial.Source: A. Curve-aware recovery · p. 4

B. Low-rank recovery

  • The annihilation equations are represented by block Toeplitz matrices whose rows encode distinct patches of the Fourier-domain partial derivatives.The number of rows equals the number of distinct patches, corresponding to the number of annihilation equations.Source: B. Low-rank recovery · p. 4
  • Because the minimal-polynomial degree is unknown, choosing a filter support larger than needed yields a low-rank matrix when a smaller-support solution exists.A polynomial multiple of that smaller-support solution also fits within the larger assumed support.Source: B. Low-rank recovery · pp. 4–5
  • The Fourier-domain formulation casts recovery as structured low-rank matrix completion and can produce a discrete image at a desired resolution by Fourier extrapolation.For noisy measurements or model mismatch, the paper proposes a convex relaxation using the nuclear norm.Source: B. Low-rank recovery · p. 5
  • To avoid the computationally prohibitive repeated SVDs, the method uses an SVD-free algorithm with auxiliary variables and ADMM.The factorization uses the nuclear-norm relation stated in the passage.Source: B. Low-rank recovery · p. 5
  • From 65×49 = 3185 Fourier samples, the method reconstructs a Shepp–Logan phantom on a 256×256 grid at approximately 20-fold undersampling; ringing artifacts are expected and may be reduced by mild post-processing.The reconstruction is obtained by extrapolating the analytical Fourier samples and applying an inverse DFT.Source: B. Low-rank recovery · p. 5

VI. CONCLUSION

The conclusion presents an annihilating-filter extension for piecewise-smooth images and a one-stage Fourier-domain method for piecewise-constant signal recovery. It reports accurate MRI-phantom recovery from few low-resolution Fourier samples.

VI. CONCLUSION

  • The annihilating-filter method extends to 2-D piecewise-smooth functions whose edges lie on the level set of a band-limited function, enabling continuous-domain edge-set recovery from few low-frequency Fourier samples.For piecewise-constant signals, the paper derives necessary and sufficient sample-count conditions for exact edge-set recovery.Source: §VI · p. 5
  • A one-stage algorithm recovers piecewise-constant signals by extrapolating the signal in the Fourier domain.Source: §VI · p. 5
  • The method is demonstrated on a Shepp–Logan phantom using 65×49 = 3185 Fourier samples to recover a 256×256 grid at approximately 20-fold undersampling.The figure shows spatial-domain images in the top row and log-scale Fourier transforms in the bottom row.Source: Figure 2 · p. 5
Loading 1502.00705v1…