Blind Deconvolution: Estimating the Blur Width

We are given a noiseless observation x0 = hσ ∗ x, blurred by an isotropic Gaussian kernel of unknown width σ. For a candidate σ, the sharp image is recovered exactly via the closed-form normal equation (Tikhonov-regularized least squares) x̂σ = argminx ‖ hσ ∗ x − x0 ‖22 + λ ‖ x ‖22 — no gradient descent needed for that inner solve. "Run" then searches over σ itself: starting from a random guess (gradient descent) or a random population (differential evolution), it repeatedly (a) solves x̂σ, (b) scores the candidate, and (c) updates σ, until it converges. Naively scoring by reblur error alone would fail here — it is minimized by less blur regardless of the truth, since blurring can only shrink a signal's energy (the classic "no-blur bias" in blind deconvolution), and with the observation noiseless, Morozov's discrepancy principle (matching that error to a known noise level) needs noise > 0 to work. Instead the score is a gradient- sparsity ratio: (mean |∇x̂σ|) / √(mean ∇x̂σ2), small when x̂'s edges are few and strong (a sharp, natural-looking image) and large when its gradient energy is smeared over many small ones (still blurry, or ringing from over-correction). Being a ratio, it needs no hand-tuned trade-off weight, unlike scoring by sharpness alone (which nothing stops from growing forever). Gradient descent moves along the exact analytic ∂loss/∂σ (hand-derived through the kernel, the normal-equation solve, and the smoothed sparsity ratio — every step differentiable, not a finite-difference approximation), with its step size set by a backtracking line search so it lands precisely rather than overshooting; differential evolution runs its usual mutate/select generations with an annealed mutation weight, then polishes its best candidate with that same gradient descent for a precise final estimate.

 | 

Target

The sharp image x we want to recover.
Resize: 100%

Observation

x0 = hσ ∗ x (noiseless)
True σ (hidden from estimator): 0.5px

Deconvolved (at current σ estimate)

x̂σ = argminx ‖ hσ ∗ x − x0 ‖2 + λ ‖ x ‖2, solved exactly
Regularization λ: 0.020
(click "Run")

Pixel-wise Difference

|x̂ − x| (Target vs. Deconvolved)
mean |xhat - x|: —

σ Estimation

Search range for σ (green = ground truth, blue = trajectory from the random initial guess to the current estimate, orange = current).
(click "Run")