🔭 astrophysics

Marginal multi-object multi-frame blind deconvolution

This paper introduces a marginal estimator for multi-object multi-frame blind deconvolution that improves upon traditional joint maximum likelihood methods by providing enhanced regularization, robust contrast control with minimal hyperparameter tuning, and simple integration into existing pipelines for high-resolution solar imaging.

Original authors: A. Asensio Ramos (IAC, ULL)

Published 2026-05-13

📖 5 min read🧠 Deep dive

CC BY 4.0

Original authors: A. Asensio Ramos (IAC, ULL)

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

The Problem: Taking a Photo Through a Wobbly Window

Imagine you are trying to take a crystal-clear photo of a distant mountain (the Sun) through a window that is constantly rippling with heat (Earth's atmosphere). Even if you have a super-fast camera that freezes the motion, the air is still wiggling, making the mountain look blurry and distorted.

Astronomers use a technique called Multi-Object Multi-Frame Blind Deconvolution (MOMFBD) to fix this. Think of it as a digital "undo" button. The computer tries to guess two things at the same time:

What the mountain actually looks like (the object).
How the window was wiggling at that exact moment (the atmospheric turbulence).

The Old Way: The "Guess-and-Check" Trap

The traditional method tries to solve for both the mountain and the wiggles simultaneously. It's like trying to solve a puzzle where you don't know the picture and you don't know the shape of the pieces.

The paper explains that this old method has a major flaw: It gets greedy.
When the photo is a bit noisy (grainy), the computer gets confused. It starts thinking that the random graininess is actually part of the mountain's texture. To "fix" the image, it invents fake, high-frequency details (like tiny, sharp spikes) to explain the noise. This is called overfitting.

To stop this, human experts have to manually tune "dials" (hyperparameters) to tell the computer, "Don't get too crazy with the details." If the dial is set too low, the image stays blurry. If it's set too high, the image looks like a cartoon with fake, jagged edges. It's a slow, artisanal process that requires a lot of experience.

The New Way: The "Marginal" Approach

The author, Andrés Asensio Ramos, proposes a smarter way to solve the puzzle. Instead of trying to guess the exact shape of the mountain and the wiggles at the same time, the new method says:

"Let's ignore the specific guess of what the mountain looks like for a moment. Let's just focus on figuring out how the window was wiggling, while acknowledging that we aren't 100% sure what the mountain looks like."

In math terms, this is called marginalization. It treats the exact image of the mountain as a "nuisance parameter"—something we need to account for, but not something we need to pin down perfectly before moving on.

The Analogy:
Imagine you are trying to hear a friend's voice (the object) in a noisy room (the atmosphere).

The Old Way: You try to perfectly reconstruct your friend's face and the noise pattern at the same time. If you hear a cough, you might accidentally think it's part of your friend's voice and start drawing a coughing face.
The New Way: You focus entirely on filtering out the room's noise. You accept that your friend's voice might be slightly unclear, but you use that uncertainty to prevent you from inventing fake words. You let the math "average out" the possibilities of what the voice could be, rather than picking one risky guess.

Why This New Method is Better

The paper highlights three main advantages of this new "Marginal" approach:

It Stops the "Fake Details" (Better Regularization):
Because the new method accounts for the uncertainty of the image, it doesn't get tricked by noise. It won't assign random graininess to high-order atmospheric errors. It naturally smooths out the fake spikes without needing a human to manually turn a dial.
It's "Plug-and-Play" (Robust Contrast):
In the old method, if you turned the "contrast" dial slightly wrong, the image would look terrible. In the new method, the math automatically finds the right balance. The hyperparameters (the dials) can be optimized by the computer itself. You don't need an expert to tweak them; you can just run the code, and it works.
It's Easy to Add:
The author shows that you don't need to rewrite the whole software. You just need to add one small mathematical term (a "log-determinant") to the existing formula. It's like adding a new ingredient to a cake recipe that makes it rise perfectly without changing the oven temperature.

The Results: Clearer Sun Photos

The author tested this on real data from solar telescopes (CRISP and HiFI).

Visuals: The new method produced images that were sharper and had better contrast than the old method, especially when looking at the edge of sunspots (the penumbra).
Automation: The new method automatically figured out how many "modes" (levels of detail) to use. In the old method, using too many modes made the image look messy. In the new method, the math naturally stopped adding detail when it started to just add noise, finding the "sweet spot" automatically.

The Bottom Line

This paper introduces a mathematical upgrade to how we clean up blurry solar images. By changing the strategy from "guessing everything at once" to "accounting for uncertainty," the new method creates clearer images, requires less human tweaking, and is ready to be used by anyone in the solar physics community. The code is already open-source and available for others to use.

Technical Summary: Marginal Multi-Object Multi-Frame Blind Deconvolution

Problem Statement
High-resolution ground-based solar imaging is fundamentally limited by atmospheric turbulence, which introduces refractive index fluctuations that degrade image quality. While Adaptive Optics (AO) systems provide partial correction, residual wavefront errors persist, necessitating post-facto image restoration to approach the theoretical diffraction limit. Traditional restoration methods, specifically Multi-Object Multi-frame Blind Deconvolution (MOMFBD), rely on joint Maximum Likelihood (ML) estimation. In this framework, the solar object and atmospheric aberrations are estimated simultaneously by maximizing the likelihood of observed data.

However, the joint ML approach faces significant challenges, particularly in low signal-to-noise ratio (SNR) regimes. Because it seeks a joint point estimate for both the object and the wavefront, the problem is often ill-posed. In the presence of noise, the algorithm tends to overfit, erroneously assigning noise power to high-order aberrations. This necessitates artisanal Fourier filtering and fine-tuning of regularization parameters. Furthermore, the contrast of the reconstruction is highly sensitive to the hyperparameters dictating the Power Spectral Density (PSD) of the object prior, requiring expert manual tuning that is time-consuming and difficult to standardize.

Methodology
The paper proposes a shift from joint point-estimation to a Marginal Maximum A-Posteriori (MMAP) framework. Instead of solving for the object and aberrations simultaneously, the objects are treated as "nuisance parameters." By defining a suitable Bayesian prior for the objects and marginalizing (integrating) over them, the method isolates the posterior probability of the aberrations alone.

The mathematical formulation proceeds as follows:

Model Setup: The observation model assumes $J$ short-exposure frames for $K$ co-spatial objects with $L$ diversity channels. The observed image is modeled as the convolution of the true object with a Point Spread Function (PSF) representing instantaneous turbulence, plus additive noise.
Bayesian Formulation: The posterior distribution of the unknowns (object $o$ and PSFs $H$ ) is defined. Assuming Gaussian priors for the object and Gaussian noise, the joint posterior factorizes into a product of two Gaussians.
Marginalization: Rather than finding the Maximum A-Posteriori (MAP) estimate of the object and plugging it back into the likelihood (as done in joint estimation), the object is integrated out. This yields a marginalized likelihood function dependent only on the PSFs and hyperparameters.
Implementation: The negative log-marginal likelihood is derived. In the Fourier domain, utilizing the Woodbury matrix identity and the Sherman-Morrison formula, the computation becomes efficient. The resulting merit function is nearly identical to the standard joint estimator but includes an additional log-determinant term.
Optimization: The method utilizes the Adam optimizer with automatic differentiation (PyTorch). A scheduling approach is employed where wavefront expansion modes are gradually increased to avoid local minima. Crucially, the hyperparameters controlling the object prior's PSD and noise variance are optimized directly alongside the wavefront coefficients during the initial iterations, after which they are fixed to prevent divergence.

Key Contributions

Enhanced Regularization: The marginalization process naturally accounts for object uncertainty. This prevents the reconstruction algorithm from assigning noise power to high-order aberrations, a common failure mode in joint estimation.
Robust Contrast Control: The marginal estimator is significantly less sensitive to the hyperparameters governing the object's PSD. This allows these hyperparameters to be optimized directly, enabling a "plug-and-play" deployment that removes the need for expert manual tuning.
Minimal Implementation Overhead: The method requires only the addition of a log-determinant term to the traditional merit function. It has been fully integrated into the open-source torchmfbd package, requiring minimal modifications to existing blind deconvolution pipelines.
Automatic Mode Selection: The behavior of the marginal loss function allows for the automatic selection of the optimal number of wavefront modes. The loss function exhibits a clear minimum, balancing data fit against regularization, whereas the joint loss function decreases monotonically, offering no such guidance.

Results
The proposed method was validated using data from the CRISP instrument (SST) and the HiFI instrument (GREGOR), covering both active regions and quiet Sun observations.

Contrast Stability: In joint estimation, the contrast of the reconstruction is highly sensitive to the PSD hyperparameter ( $\kappa$ ). Small values lead to blurry images, while large values cause artificial contrast enhancement and artifacts. The marginal estimator produces stable, realistic contrast without tuning.
Frame Count Independence: When varying the number of frames ( $M$ ) in the burst, joint estimators require careful adjustment of the PSD parameter to avoid artifacts or blurriness. The marginal estimator automatically produces reconstructions with appropriate contrast and visual quality across different frame counts (e.g., $M=2$ to $M=50$ ).
Mode Count Optimization: Experiments with synthetic data demonstrated that the marginal loss function has a distinct minimum corresponding to the optimal number of KL modes (e.g., 60–80 modes for $r_0=10$ cm). In contrast, the joint loss function decreases monotonically with increasing modes, failing to indicate an optimal stopping point.
Visual Quality: Reconstructions using the marginal estimator displayed slightly more defined structures and higher contrast compared to standard MOMFBD and joint estimators, particularly in penumbral regions, without introducing high-frequency artifacts.

Significance
The paper claims that the marginal estimator addresses the inherent limitations of traditional joint Maximum Likelihood estimation in solar imaging. By providing a stable reconstruction that is less sensitive to high-frequency noise and hyperparameter choices, the method facilitates a more automated and robust workflow. The ability to optimize hyperparameters and select the number of modes directly from the data reduces the reliance on expert intervention, making high-quality solar image restoration more accessible to the broader solar physics community. The integration into torchmfbd ensures immediate availability for routine data processing.