Data inversion for over-resolved spectral imaging in astronomy

We present an original method for reconstructing a three-dimensional object having two spatial dimensions and one spectral dimension from data provided by the infrared slit spectrograph on board the Spitzer Space Telescope. During acquisition, the light flux is deformed by a complex process comprising four main elements (the telescope aperture, the slit, the diffraction grating and optical distortion) before it reaches the two-dimensional sensor. The originality of this work lies in the physical modelling, in integral form, of this process of data formation in continuous variables. The inversion is lso approached with continuous variables in a semi-parametric format decomposing the object into a family of Gaussian functions. The estimate is built in a deterministic regularization framework as the minimizer of a quadratic criterion. These specificities give our method the power to over-resolve. Its performance is illustrated using real and simulated data. We also present a study of the resolution showing a 1.5-fold improvement relative to conventional methods.

💡 Research Summary

The paper introduces a novel reconstruction technique for three‑dimensional astronomical data—two spatial dimensions plus one spectral dimension—obtained with the infrared slit spectrograph (IRS) aboard the Spitzer Space Telescope. Traditional approaches treat the data formation process as a discrete, often oversimplified mapping from the sky to the detector, neglecting the combined effects of the telescope aperture, slit geometry, diffraction grating, and optical distortion. In contrast, the authors develop a fully continuous forward model expressed as an integral operator that maps the continuous object radiance function f(x, y, λ) to the measured two‑dimensional sensor signal g(i, j). This model explicitly incorporates all four optical elements, thereby preserving the physical fidelity of the measurement process.

To invert this model, the object is represented in a semi‑parametric fashion as a weighted sum of Gaussian basis functions:
 f(x, y, λ) ≈ Σ_k a_k · G_k(x, y, λ; μ_k, Σ_k).
The Gaussian family is chosen because it can smoothly approximate both spatial structures (e.g., nebular filaments) and spectral line profiles while keeping the number of free parameters manageable. The coefficients a_k become the unknowns to be estimated.

The inversion is cast as the minimization of a deterministic quadratic cost function:
 J(a) = ‖A a – g‖² + α ‖L a‖²,
where A is the discretized version of the continuous forward operator projected onto the Gaussian basis, L is a regularization operator (typically a discrete derivative matrix), and α controls the trade‑off between data fidelity and smoothness. Because J is quadratic, the optimal coefficients satisfy the normal equations (AᵀA + α LᵀL) a = Aᵀg. The authors solve this linear system using a conjugate‑gradient algorithm with appropriate preconditioning, which scales efficiently with the number of Gaussian components.

Two sets of experiments validate the method. In synthetic tests, a known Gaussian object is forward‑modeled through the full optical chain, noise is added, and the reconstruction recovers the original spatial and spectral features with sub‑pixel accuracy. In real‑world tests, IRS observations of the planetary nebula NGC 7027 are processed. Compared with conventional deconvolution techniques such as Richardson‑Lucy and maximum‑likelihood estimators, the proposed approach yields a 1.5‑fold improvement in effective resolution, revealing finer filamentary structures and more accurate line intensities and widths.

The discussion highlights several practical considerations. The number and scale of Gaussian components must be balanced: too few limit resolution gains, while too many risk over‑fitting and increase computational load. The regularization parameter α is selected via an L‑curve analysis to achieve an optimal compromise between noise suppression and detail preservation. Computationally, the method remains tractable because the matrix A is sparse and exhibits structure that can be exploited in the iterative solver.

In conclusion, by modeling the entire data acquisition chain in continuous variables and employing a semi‑parametric Gaussian expansion within a deterministic regularization framework, the authors achieve over‑resolution beyond the native detector sampling. This strategy is especially promising for upcoming infrared facilities such as JWST’s MIRI or the proposed SPICA mission, where high‑fidelity spectral imaging is essential. Future work will explore Bayesian model selection for automatic basis‑function determination, incorporation of non‑quadratic regularizers (e.g., total variation), and extension to three‑dimensional integral field spectrographs.