# Unsupervised Imaging Inverse Problems with Diffusion Distribution Matching

Giacomo Meanti, Thomas Ryckeboer, Michael Arbel, Julien Mairal

► To cite this version:

Giacomo Meanti, Thomas Ryckeboer, Michael Arbel, Julien Mairal. Unsupervised Imaging Inverse Problems with Diffusion Distribution Matching. International Conference on Computer Vision (ICCV), Oct 2025, Honolulu, United States. hal-05120427

HAL Id: hal-05120427

https://hal.science/hal-05120427v1

Submitted on 19 Jun 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Unsupervised Imaging Inverse Problems
with Diffusion Distribution Matching

Giacomo Meanti Thomas Ryckeboer Michael Arbel Julien Mairal

Univ. Grenoble Alpes, Inria, CNRS, Grenoble INP, LJK, 38000 Grenoble, France

firstname.lastname@inria.fr

Abstract

This work addresses image restoration tasks through the lens of inverse problems using unpaired datasets. In contrast to traditional approaches—which typically assume full knowledge of the forward model or access to paired degraded and ground-truth images—the proposed method operates under minimal assumptions and relies only on small, unpaired datasets. This makes it particularly well-suited for real-world scenarios, where the forward model is often unknown or mis-specified, and collecting paired data is costly or infeasible. The method leverages conditional flow matching to model the distribution of degraded observations, while simultaneously learning the forward model via a distribution-matching loss that arises naturally from the framework. Empirically, it outperforms both single-image blind and unsupervised approaches on deblurring and non-uniform point spread function (PSF) calibration tasks. It also matches state-of-the-art performance on blind super-resolution. We also showcase the effectiveness of our method with a proof of concept for lens calibration: a real-world application traditionally requiring time-consuming experiments and specialized equipment. In contrast, our approach achieves this with minimal data acquisition effort. Code available: https://github.com/inria-thoth/ddm4ip.

<figure><img src="image_2.png" /><figcaption>Figure 1: Learning degradation operators by matching distributions.</figcaption></figure>

# 1 Introduction

Inverse problems arise in many domains where observed data is the result of a noisy and potentially lossy transformation of an underlying signal we wish to recover. More formally, the degradation process is modeled as follows:

$$y = \mathcal{A}(x) + \epsilon, \quad \epsilon \sim \mathcal{N}(0, \sigma^2 I), \qquad (1)$$

where $\mathcal{A}$ is the forward operator and $\epsilon$ denotes additive Gaussian noise. The objective is to recover the original signal $x$ from the measurements $y$. In practice, $\mathcal{A}$ may represent a wide range of linear or nonlinear transformations, such as blur kernels, inpainting masks, the Radon transform, and many others.

To introduce the setting of the current work, we categorize the reconstruction methods based on i) the type and availability of data, ii) the extent of prior knowledge about $\mathcal{A}$. When $\mathcal{A}$ is fully known, two main classes of reconstruction approaches are typically used. The first is plug-and-play
(PnP) methods, which apply iterative algorithms to recover signal $x$ from a single measurement $y$ by maximizing the posterior distribution $p(x | y, \mathcal{A}) \propto p(y | x, \mathcal{A})p(x)$. These methods alternate between a data-fidelity term (dependent on $\mathcal{A}$) and a prior term. The second class leverages supervised learning: if a clean dataset is available, it can be paired with synthetic measurements generated via $\mathcal{A}$ to train a reconstruction model in a supervised fashion.

However, complete knowledge of the forward operator $\mathcal{A}$ is rarely available in practice. For instance, in super-resolution it is common to assume bicubic downsampling to generate low-resolution images from high-resolution ones, yet models trained under this assumption often fail on real-world images where the actual downsampling kernel differs [18, 59].

When $\mathcal{A}$ is unknown or only partially specified, the problem becomes more challenging and is called *blind*. One strategy is to collect (or synthetically generate) a dataset that includes corrupted samples produced from a range of possible $\mathcal{A}$ instances. A supervised model is then trained to generalize across these degradations. An alternative approach is to jointly infer both the clean signal and the unknown operator by maximizing the joint posterior $p(x, \mathcal{A} | y) \propto p(y | x)p(x)p(\mathcal{A})$. This formulation requires priors on both the data distribution $p(x)$ and the corruption process $p(\mathcal{A})$. Both approaches implicitly rely on some prior knowledge or assumptions about $\mathcal{A}$, either through heuristic degradation models or learned probabilistic priors. Consequently, when the true degradation deviates from these assumptions—which is often the case in real-world scenarios—the reconstruction performance degrades significantly. Although these methods typically require very little input data at test time (sometimes even a single image), their robustness to real-world variations is far from guaranteed.

The natural alternative, which we adopt in this work, is to learn information about $\mathcal{A}$ directly from data: specifically, from one dataset of corrupted images and a separate dataset of clean images. Crucially, these datasets do not need to contain corresponding clean-corrupted pairs, which significantly simplifies data collection. During training, the unpaired clean images serve as a reference for what the restored outputs should resemble—effectively providing a model for $p(x)$. This unpaired training setup is often referred to as *unsupervised*, and has been successfully applied to in-the-wild image restoration tasks [31, 45, 51], where the degradation process is unknown and may involve multiple, complex transformations. Most existing algorithms operating in this unpaired regime use an implicit model of the degradation process: a neural network is trained to map corrupted inputs to clean outputs, without explicitly modeling the underlying forward operator.

The method we propose learns the correct degradation by learning an explicit operator $\mathcal{A}$ and needs a relatively small (hundreds to thousands) number of training points. By using diffusion models and an efficient distribution matching algorithm it outperforms both single-image methods and other unsupervised ones. It comprises two distinct steps: the first one is to learn a representation of the noisy distribution, the second is to learn the best degradation operator $\mathcal{A}$ such that the clean samples corrupted by $\mathcal{A}$ match *in distribution* the noisy data prior. Finally the learned corruption can be used in a third step for non-blind restoration. Our contributions are threefold: we begin by devising a principled method for solving imaging inverse problems which only relies on unpaired image data from clean and corrupted distributions, without knowledge of the degradation operator. We prove that under a non-degeneracy assumption on the clean distribution, the true operator can be identified. We then show how this method can provide precise estimates of the degradation operator, without making any distributional assumptions on the operator itself. In particular we focus on the tasks of uniform and non-uniform deblurring. Finally, since the estimates we obtain are considerably more precise than those of single-image methods, we show how our approach can be used as part of a pipeline for camera lens calibration – where accuracy is essential.

The rest of the paper is organized as follows. In Section 2, we provide some necessary background on inverse problems with unknown degradations and on our algorithm. In Section 3, we detail the different steps of our algorithm, whereas the last section is devoted to experiments of increasing complexity on estimating blur kernels.
<table><caption>Table 1: Summary of tasks for solving inverse problems without knowing the forward operator. The headings refer to what is known about the data and the forward operator. Regarding the latter, we note when the method works on data corrupted by a single (unknown) operator, operators from a predefined distribution, or by a collection of unknown operators. The references are not meant to be exhaustive but just to provide examples.</caption><thead><tr><th>Data</th><th>A</th><th>Method</th><th>References</th></tr></thead><tbody><tr><td>prior on <i>A</i></td><td>single known</td><td>PnP</td><td>[42, 47, 57, 62]</td></tr><tr><td>paired</td><td>single unknown</td><td>supervised</td><td>[49]</td></tr><tr><td>prior on <i>A</i></td><td>prior on <i>A</i></td><td>blind PnP</td><td>[12, 23, 44, 52]</td></tr><tr><td>paired</td><td>prior on <i>A</i></td><td>supervised</td><td>[18, 50, 58]</td></tr><tr><td>unpaired</td><td>multiple unknown</td><td>unsupervised</td><td>[34, 45, 51]</td></tr><tr><td>unpaired</td><td>functional form</td><td></td><td><i>this paper</i></td></tr></tbody></table>

# 2 Background

We briefly review the literature on blind and unsupervised inverse problems in imaging, along with the various forms of deblurring addressed in this work.

## 2.1 Unsupervised inverse problems

Solving inverse problems without access to the forward operator is inherently challenging and can be approached from multiple perspectives. In Table 1 we provide a framework to clarify the various settings considered in the literature. Unlike many of the methods discussed below, our approach assumes a fixed degradation operator $A$ which does not vary across measurements $y$. While this restricts generality, it enables higher reconstruction accuracy in specific settings compared to other unpaired methods. We focus on unpaired algorithms, particularly for deblurring and super-resolution tasks, where $A$ typically corresponds to a blur kernel. Data augmentation-based approaches [50, 58] construct synthetic supervised datasets using heuristically designed degradation pipelines which mimic diverse real-world corruptions. To improve adaptability to specific degradations, Zhang et al. [59] propose learning some of the pipeline parameters from a small dataset of noisy reference images. For better alignment with the degradation of each individual image, a line of work tackles the blind MAP inference problem using deep priors over both clean images and degradation operators. These methods are typically coupled with test-time optimization algorithms such as alternating minimization [41], expectation-maximization (EM) [17, 23], or proximal gradient methods [48]. Recent approaches [12, 23, 44] incorporate diffusion models as priors for the blur kernel, jointly estimating $A$ and $x$ via plug-and-play algorithms [47]. Closely related methods include FKP [27], which uses a pretrained normalizing flow as a kernel prior, and DKP [52], which employs MCMC to iteratively refine the blur estimate. GibbsDDRM [38], by contrast, uses a diffusion prior on the clean image and a simpler total variation (TV) prior on the blur kernel.

In contrast to this class of methods, which require only a single degraded image at test time, our approach does not rely on any pretrained priors—neither on clean images nor on the degradation operator. Instead, it learns to adapt directly to the specific corruption process represented in the training data. While it does require a small dataset of degraded images from the same corruption distribution, this targeted adaptation enables significantly improved reconstruction accuracy.

Domain transfer approaches based on variations of the cycle-consistency loss [61] are conceptually closer to our method. Given a noisy image $y$ and an unpaired clean image $x$, a clean-image generator $G$ and a noisy-image generator $F$ are jointly trained under the constraint that $G(F(x)) \approx x$ and $F(G(y)) \approx y$. For instance, CinCGAN [36, 53] translates images downsampled with an unknown kernel into bicubically downsampled images, which can then be more effectively upscaled using standard super-resolution models. Several related methods [9, 30, 34, 37, 45, 51] employ generative models to synthesize corrupted images from clean ones, thereby enabling the construction of supervised training datasets. This synthesis can be done in a two-stage pipeline or directly in an
end-to-end manner [11, 43]. Notably, Sim et al. [45] address a setting similar to ours, where a model is trained to deblur microscopy images in an unpaired setup. Compared to these approaches, our method leverages a novel loss function derived from diffusion models—a technique not previously applied in this context. Unlike adversarial losses, our formulation is more stable and easier to train. Moreover, by learning an explicit representation of the degradation kernel rather than an implicit one, our method offers significantly improved interpretability. Furthermore, under certain assumptions, we are able to prove the identifiability of the true operator (see Section 3).

We also briefly note that many methods aim to learn the degradation operator in a single-image fashion, but still require paired training data—typically synthetic—for supervision. For instance, IKC [18] iteratively refines kernel estimates in alternation with clean-image predictions. DAN [35] extends this idea by unrolling the refinement process into an end-to-end trainable network. In contrast, KernelGAN [6] learns the degradation operator directly from a single input image, followed by a separate non-blind model to reconstruct the clean image $x$.

## 2.2 Distribution matching

Our method minimizes the distance between two probability distributions over images: the distribution of the observed noisy data $p(y)$ and that of clean data corrupted by a learned degradation operator, i.e., $p(A(x))$. A common strategy for such distribution matching is to use generative adversarial networks (GANs), which optimize an adversarial loss. GANs have been extensively applied to image restoration tasks, as discussed in the previous section [34, 45], and can simultaneously learn both the degradation operator and its inverse within a unified training framework. However, training GANs is notoriously challenging due to instability and sensitivity to hyperparameters [8]. Alternative approaches include normalizing flows, which provide exact likelihoods and invertible mappings. These have been used by DeFlow [51] for matching clean and noisy distributions, and by FKP [27] to generate plausible blur kernels from single images.

More recently, diffusion models – like GANs and normalizing flows – have emerged as powerful tools for modeling empirical distributions and have been used in a range of distribution-matching tasks. For example, DreamFusion [39] learns 3D representations whose 2D projections are consistent with a pretrained diffusion model, while DiffInstruct [32] trains a single-step generator to match the distribution of a multi-step diffusion model. Both approaches rely on a loss function that approximates the KL divergence integrated over all diffusion time steps.

In our work, we use conditional flow matching (CFM) models [28, 29] as a more conceptually straightforward alternative to standard diffusion models. We adapt the integrated KL divergence loss to the CFM framework for learning the degradation operator.

## 2.3 Camera lens calibration

After preliminary experiments on synthetic deblurring, we focus on non-uniform deblurring in the context of a camera calibration pipeline. Camera calibration is a multifaceted task typically broken down into subtasks such as distortion, chromatic aberration, and vignetting corrections among others. In this work, we concentrate solely on compensating for blur induced by lens imperfections, with the goal of obtaining maximally sharp images. These lens aberrations can be characterized by the point spread function (PSF) of the lens-camera system: the blur kernel that transforms an ideal point light source into a spread of colored spots in the image. PSFs are often non-uniform across the image plane and differ across RGB channels, resulting in chromatic aberrations. Correcting such distortions computationally is particularly appealing, as it enables the use of lower-cost lenses without sacrificing image quality.

Traditional non-blind correction methods require precise knowledge of the PSF, which can only be obtained through laborious procedures involving printed calibration patterns or specialized screens [5, 21]. In contrast, blind lens aberration correction remains an understudied problem, although it shares many similarities with non-uniform deblurring. It is more commonly approached under the umbrella of defocus deblurring [2, 24, 40, 54], where blur magnitude varies with scene
depth. In aberration correction, however, the blur is dependent on spatial location in the image plane rather than depth.

Nonetheless, defocus deblurring methods may still be partially effective for lens aberration correction, especially when the induced blur is close to isotropic. As a starting point, we use the PSF dataset from Bauer et al. [5], which provides ground-truth PSFs, and subsequently transition to a realistic calibration scenario using a Panasonic Micro 4/3 camera. In this setting, we acquire a small set of images and aim to learn the lens PSFs without access to ground-truth measurements. Our main comparison will be with the blind method proposed by Eboli et al. [16], which performs deblurring followed by color-fringe removal in a two-stage pipeline [10, 15].

# 3 Method

Given an inverse problem of the form $y = A_{\omega_{*}}(x) + \epsilon$, with known noise level $\sigma$ and an unknown forward operator $A_{\omega_{*}}$ parameterized by a vector $\omega_{*}$, we propose an algorithm to estimate $\omega_{*}$ using unpaired data. On the one hand, we assume access to a clean dataset $\mathcal{X} = \{x_i\}_{i=1}^n$ of images drawn from a distribution $\mathcal{X}$. On the other hand, we have access to a corrupted dataset $\mathcal{Y} = \{y_j\}_{j=1}^m$ of images, generated from unknown clean samples via the forward model $A_{\omega_{*}}$. We denote by $\mathbb{Y}_{\omega_{*}}$ the distribution of such corrupted images, and by $\mathbb{Y}_{\omega}$ the distribution induced by applying an arbitrary forward operator $A_{\omega}$ (with additive noise) to samples from $\mathcal{X}$.

Noting that $\mathbb{Y}_{\omega_{*}}$ is empirically accessible through $\mathcal{Y}$ while $\mathbb{Y}_{\omega}$ can be approximated via $\mathcal{X}$ and a candidate forward operator $A_{\omega}$, we hypothesize that minimizing the distance between $\mathbb{Y}_{\omega_{*}}$ and $\mathbb{Y}_{\omega}$ with respect to $A_{\omega}$ yields a good approximation of the true forward model parameters:

$$ \hat{\omega} = \arg\min_{\omega} \mathcal{D}(\mathbb{Y}_{\omega_{*}}, \mathbb{Y}_{\omega}) \implies A_{\hat{\omega}} \approx A_{\omega_{*}}, \qquad (2) $$

where $\mathcal{D}$ is a distance between distribution that will be specified in detail later. We prove this rigorously in a simplified setting in Proposition 3.1, with the full proof provided in the supplementary material, and we assume that the result generalizes when the distributions are only approximately equal, as is typically the case in practice.

**Proposition 3.1:** For any set of forward model parameters $\omega$, let $p_{\omega}(y) = \int p_{\omega}(y | x)p(x)dx$ where $p_{\omega}(y | x) = \mathcal{N}(y | A_{\omega}x, \sigma^2I)$. Let $\omega_{*}$ be a specific set of parameters which we consider to be the optimal set. Then, assuming the data covariance $\Sigma = \mathbb{E}_{x}[xx^{\top}]$ is invertible, there exists an orthogonal matrix $P$ such that

$$ p_{\omega}(y) = p_{\omega_{*}}(y) \implies A_{\omega} = A_{\omega_{*}}\Sigma^{1/2}P\Sigma^{-1/2} \qquad (3) $$

That is, if the probability distributions $p_{\omega}$ and $p_{\omega_{*}}$ are equal, it is possible to identify $A_{\omega_{*}}$ up to rotations $P$.

Since neither $\mathbb{Y}_{\omega_{*}}$ nor $\mathbb{Y}_{\omega}$ are explicitly available we will use CFM models to act as tractable representations of probability distributions on images.

## 3.1 Distribution matching with diffusion

A CFM model $v_{\theta}(z^{(t)}, t)$ trained on data from a distribution $\mathbb{Y}_{\omega}$ allows one to sample from $\mathbb{Y}_{\omega}$ by solving an ODE between times 0 and 1, defined as

$$ dz^{(t)} = v_{\theta}(z^{(t)}, t)dt, \quad z^{(0)} \sim \mathcal{N}(0, I). $$

such that $z^{(t=1)} \sim \mathbb{Y}_{\omega}$. This conditional flow matching [28, 29] perspective is connected to standard Langevin diffusion since the velocity field $v_{\theta}(x^{(t)}, t)$ is related to the score of a diffusion model: $v_{\theta}(z^{(t)}, t) + z^{(t)} \propto \nabla_z \log p_{\mathbb{Y}_{\omega}^{(t)}}(z_{\omega}^{(t)})$ [60], assuming the neural network model of the velocity field to be exact. Each $z^{(t)}$ follows an intermediate distribution between Gaussian noise and the data $\mathbb{Y}_{\omega}^{(t)}$ for $t \in [0, 1]$. Hence, instead of computing a distance between two distributions, we will compute
it between two sequences of distributions $\mathbb{Y}_{\omega}^{(t)}$ and $\mathbb{Y}_{\omega_{*}}^{(t)}$ which arise from two CFM models. In particular we will use the KL divergence integrated over time, and follow the derivation proposed by DiffInstruct (DI) [32] to compute its gradient with respect to $\omega$. The integrated KL divergence is defined as

$$ \mathrm{IKL}(\mathbb{Y}_{\omega_{*}} \| \mathbb{Y}_{\omega}) = \int_{t=0}^{1} \mathrm{KL}(\mathbb{Y}_{\omega_{*}}^{(t)} \| \mathbb{Y}_{\omega}^{(t)}) dt. \quad (4) $$

Assuming that the score terms for both probability distributions (i.e. $s_{\theta, \omega_{*}}(y, t) \approx \nabla_{y} \log p_{\mathbb{Y}_{\omega_{*}}^{(t)}}(y)$ and $s_{\phi, \omega}(y, t) \approx \nabla_{y} \log p_{\mathbb{Y}_{\omega}^{(t)}}(y)$) can be computed, it is possible to efficiently differentiate the IKL with respect to $\omega$:

$$ \nabla_{\omega} \mathrm{IKL}(\mathbb{Y}_{\omega_{*}} \| \mathbb{Y}_{\omega}) = \int_{t=0}^{1} \mathbb{E} [s_{\theta, \omega_{*}}(y_{\omega}^{(t)}, t) - s_{\phi, \omega}(y_{\omega}^{(t)}, t)]^{\top} \nabla_{\omega} y_{\omega}^{(t)} dt, \quad (5) $$

where the expectation is over $y_{\omega}^{(0)} \sim \mathcal{N}(0, I)$, $x \sim \mathcal{X}$, $y_{\omega}^{(1)} = \mathcal{A}_{\omega}(x) + \epsilon$ and $y_{\omega}^{(t)} = (1-t)y_{\omega}^{(0)} + ty_{\omega}^{(1)}$. Note that eq. (5) has in fact been originally introduced to optimize a 3D scene consistent with a diffusion model's outputs [39] and in [32] was used to learn a distilled model. Here we are showing it can be used effectively on conditional flow matching models as well, and in a completely different setting. The loss gradient shown in eq. (5) needs certain quantities to be computed: the score of $p_{\mathbb{Y}_{\omega_{*}}}$ can be obtained beforehand by training a flow-matching model $v_{\theta, \omega_{*}}$ on the available noisy data $\mathcal{Y}$. The score of $p_{\mathbb{Y}_{\omega}}$ however depends on a distribution which changes with every change in $\omega$. Therefore the strategy we use is to alternately optimize i) an auxiliary flow-matching model $v_{\psi, \omega}$ for a fixed $\omega$ with a standard diffusion loss and ii) forward model parameters $\omega$ with the IKL loss keeping the auxiliary model fixed [32]. By changing the parametrization of $\mathcal{A}_{\omega}$ we can use eq. (5) to learn the degradation of a variety of different inverse problems. We will now go into the details of the algorithm for the task of learning blur operators.

## 3.2 Learning the forward operator

In deblurring, $\mathcal{A}_{\omega}$ corresponds to a convolution with kernel $\omega$. The algorithm we propose proceeds in three distinct steps: i) learn the corrupted data distribution, ii) approximate the forward operator $\mathcal{A}_{\omega_{*}}$ and iii) solve the non-blind inverse problem.

**Step 1: Learning $\mathbb{Y}_{\omega_{*}}$** In more detail, the first step uses the noisy dataset $\{y_j\}_{j=1}^m$ to train a conditional flow matching (CFM) model $v_{\theta, \omega_{*}}$ that transports Gaussian noise samples at time $t=0$, $z^{(0)} \sim \mathcal{N}(0, I)$, to samples from the corrupted data distribution $\mathbb{Y}_{\omega_{*}}$ at time $t=1$, i.e., $z^{(1)} \sim \mathbb{Y}_{\omega_{*}}$. The training objective is the standard conditional flow matching loss [28, 29]:

$$ \mathcal{L}_{\mathrm{CFM}} = \mathbb{E}_{t,z^{(0)},z^{(1)}} [\|v_{\theta,\omega_{*}}(z^{(t)},t) - (z^{(1)} - z^{(0)})\|]. $$

Importantly, for the overall success of our method, it is not necessary for $v_{\theta, \omega_{*}}$ to achieve state-of-the-art generative quality; it only needs to effectively capture the degradation process, which we found to be significantly easier than precisely modeling image content. See Fig. 3 for sample outputs produced by a model at this stage.

**Step 2: distribution matching** The second step uses the clean dataset $\mathcal{X}$ to learn the corruption operator $\hat{\mathcal{A}}_{\omega} \approx \mathcal{A}_{\omega_{*}}$ by minimizing the IKL loss (4). This involves two alternating optimization steps: first, training an auxiliary diffusion model $v_{\phi, \omega}(z^{(t)}, t)$ using the standard flow-matching loss, where $z^{(t)} = (1-t)z^{(0)} + t\mathcal{A}_{\omega}(x)$; second, updating $\hat{\mathcal{A}}_{\omega}$ by following the gradient in eq. (5).

To encourage fast convergence, the auxiliary model $v_{\phi, \omega}$ is initialized with the pretrained weights from $v_{\theta, \omega_{*}}(z^{(t)}, t)$. The parameterization of $\hat{\mathcal{A}}_{\omega}$ is flexible: it may depend on the clean image $x$—which limits certain options for the final step—or be independent of $x$. For instance, we consider modeling non-uniform blur with an operator $\hat{\mathcal{A}}_{\omega}$ that varies with pixel location. The framework is general, but $\hat{\mathcal{A}}_{\omega}$ should not depend on corrupted images, as these are unavailable in the unpaired training setting.
<figure><img src="image_3.png" /><figcaption>Figure 2: Effect of center regularization on reconstruction quality. Without regularization the learned blur filter may be shifted leading to larger reconstruction errors (second row of the plot). A moderate amount of regularization fixes this.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 3: Samples generated by flow matching model with motion-blur kernel. While the model is average at generating faces, it correctly represents the degradation.</figcaption></figure>

The distribution matching step can also be regularized in various ways to introduce prior knowledge about the forward model and improve the quality of results. For example, when $\mathcal{Y}$ and $\mathcal{X}$ consist of patch data, their distribution will be invariant to translation: an image patch translated by some amount in any direction will still follow the same distribution. This invariant can easily lead to learning blur filters which are off-center as shown in Fig. 2. To counter this we add a regularizer which constrains the center of mass of the learned kernels to be in the middle of the filter itself.

**Step 3: solve the inverse problem** The third and final step is to use the learned forward model $\hat{A}_{\omega}$ to solve the inverse problem. At this point the non-blind setting applies, hence multiple strategies can be chosen depending on the specific problem and on the type of data which is available. When a larger clean-image dataset is available, by leveraging $\hat{A}_{\omega}$ it is possible to generate a paired noisy-clean dataset on which to train a supervised image restoration model [49, 54]. The second option is to use a plug and play algorithm [47] which leverages a pretrained prior on clean images (classically this could have been a prior such as total-variation instead) to iteratively convert a noisy image into a clean one. In Section 4 we will use ESRGAN [49] from the first option and DPIR [55], DiffPIR [62] from the second one. Of course classical alternatives such as Wiener deconvolution can be used depending on the specific inverse problem.

# 4 Experiments

## 4.1 Deblurring

We used subsets of the FFHQ dataset [19] to compare with blind deblurring methods. In particular, for each of two degradation operators we train a small CFM model on 1000 images from FFHQ corrupted by the blurring operator $A_{\omega*}$ and subsequently utilize 100 different clean images to learn $\hat{A}_{\omega}$ with distribution matching. We use an isotropic Gaussian blur with standard deviation 1 and a motion blur kernel generated following Borodenko [7]. In both cases, Gaussian noise is added with standard deviation of 0.02. We finally use plug and play algorithm DiffPIR [62] to solve the inverse problem $\hat{A}_{\omega}$. The natural upper baseline for this experiment is to run the same PnP algorithm with the ground-truth kernel. In addition to this, we compare with some lower baselines which learn the correct kernel from a *single image*, and simultaneously perform deblurring: diffusion based methods BlindDPS [12], FastDiffusionEM [23] and KernelDiff [44]. These latter blind algorithms require two priors (in the form of pretrained diffusion models): one on the image dataset and one on the kernels. BlindDPS and FastDiffusionEM use FFHQ as image prior while KernelDiff uses several natural image datasets [14]. All three methods assume the blur kernels to be motion-blur, generated
<table><caption>Table 2: Reconstruction error on FFHQ. KernelDiff was run with no noise. FastEM was run with 16 samples and IGDM.</caption><thead><tr><th></th><th colspan="2">Gaussian (σ = 1)</th><th colspan="2">Motion</th></tr><tr><th></th><th>PSNR ↑</th><th>LPIPS ↓</th><th>PSNR ↑</th><th>LPIPS ↓</th></tr></thead><tbody><tr><td colspan="5" style="text-align: center;">Non-blind (true kernel available)</td></tr><tr><td>DiffPIR</td><td>32.7</td><td>0.031</td><td>29.3</td><td>0.068</td></tr><tr><td colspan="5" style="text-align: center;">Blind, multi-image</td></tr><tr><td>Ours + DiffPIR</td><td>32.7</td><td>0.031</td><td>28.8</td><td>0.069</td></tr><tr><td colspan="5" style="text-align: center;">Blind, single-image</td></tr><tr><td>Blind-DPS [12]</td><td>29.8</td><td>0.077</td><td>25.5</td><td>0.122</td></tr><tr><td>FastEM (n=16) [23]</td><td>—</td><td>—</td><td>24.9</td><td>0.130</td></tr><tr><td>KernelDiff [44]</td><td>—</td><td>—</td><td>21.7</td><td>0.194</td></tr></tbody></table>

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

using the same stochastic procedure [7] as the true one. BlindDPS additionally includes isotropic Gaussian kernels in its kernel prior, which is why in Table 2 we also test it on a Gaussian kernel. It is important to stress that our algorithm works in a different setting from the single-image methods: we require a dataset of clean and noisy images for each degradation, while the single-image algorithms only require a single noisy image. However, note that all single-image algorithms rely strongly on the image and kernel priors on which they were trained: for example we could not successfully run the method from Laroche et al. [23] on a simple Gaussian kernel without first retraining the kernel prior, and KernelDiff which uses a different image prior than FFHQ performs poorly in this setting. The results in Table 2 demonstrate two things: first, even when the problem is purely in-distribution, the gap between blind and non-blind algorithms is large. In second instance, the algorithm we propose significantly reduces this gap by using more data from the same distribution to learn the necessary information about the degradation operator. To increase the robustness of our experiments we compute the standard deviation over 5 different random seeds for the 2nd step of our pipeline (thus keeping the 1st step fixed). For the motion blur kernel, we find a standard deviation of 0.02 for PSNR and 0.0005 for LPIPS which are both negligible.

## 4.2 Space-varying blur

In practice image blur may come from a variety of sources such as camera shake/motion, out-of-focus objects or lens distortions. We focus on the latter which better fits the our framework: it is easy to collect a set of noisy images with equal distribution, but it is not really possible to collect paired datasets (note that paired datasets can be generated synthetically through software camera models [26]). Every imaging system is imperfect, which can be characterized by a spatially varying point spread function (PSF) at every location in the image plane. It can thus be modeled as a per-pixel blur where the blur kernel depends on the image plane location. Since the PSFs are intrinsically connected to the imaging system, multiple pictures taken with the same camera and camera settings will have the same degradation.

For the first experiment we use real PSFs, taken from a subset of those identified in a real camera system [5] but synthetically applied to a clean dataset. In particular we use the green channel PSFs from the Canon EF 24mm f/1.4L II USM lens at 1.4 aperture, and subsample them to a 8x8 grid, shown in Fig. 5 (top-left). This grid is then mapped to images from DIV2K [3] and DPDD [1], and applied as a per-pixel blur by linearly interpolating the kernels to each image location. We use the degraded DIV2K training-set (subdivided into patches 128 pixels wide) to train the diffusion model in the first stage of our algorithm, and the clean DDPD training-set for the second stage – thus ensuring the data is strictly unpaired. For the third and final stage we experiment with two
<table><caption>Table 3: Reconstruction error for non-uniform deblurring on DPDD. The degradation is a spatially varying blur with additive Gaussian noise ($\sigma = 0.01$).</caption><thead><tr><th>Method</th><th>PSNR ↑</th><th>SSIM ↑</th><th>LPIPS ↓</th></tr></thead><tbody><tr><td colspan="4"><em>Non-blind (true kernel available)</em></td></tr><tr><td>ESRGAN [49]</td><td>32.53</td><td>0.906</td><td>0.058</td></tr><tr><td>DPIR [55]</td><td>34.89</td><td>0.933</td><td>0.112</td></tr><tr><td colspan="4"><em>Blind, multi-image</em></td></tr><tr><td>DeFlow+ESRGAN [51]</td><td>27.92</td><td>0.818</td><td>0.206</td></tr><tr><td>Ours+ESRGAN</td><td>32.38</td><td>0.911</td><td>0.061</td></tr><tr><td>Ours+DPIR</td><td>33.84</td><td>0.933</td><td>0.099</td></tr><tr><td colspan="4"><em>Blind, single-image</em></td></tr><tr><td>INIKNet [40]</td><td>29.45</td><td>0.874</td><td>0.160</td></tr><tr><td>Restormer [54]</td><td>28.59</td><td>0.859</td><td>0.203</td></tr></tbody></table>

<figure><img src="image_7.png" /><figcaption>Figure 4: Sample reconstructions on DDPD. INIKNet and Restormer are single image, DPIR is non-blind.</figcaption></figure>

procedures: the plug and play algorithm DPIR [55] which uses a CNN as regularizing prior and can be applied directly to the learned kernels, and supervised method ESRGAN [49] for which we generate a paired clean-noisy dataset using the learned degradation. In order to condition the kernels on image location, we add two positional encoding channels to our images such that the diffusion model learns different distributions based on the patch location. For the degradation prediction we directly learn the 64 kernels without any additional parametrization; both at train and at test time the correct kernels are picked by using the positional encoding and linearly interpolating between the kernels. For this experiment we used the centering regularization and also introduced an isotropic Gaussian regularization which helped stabilize training.

We compare the results obtained on the “target” test-set of DPDD (corrupted with the camera PSFs and not with the original defocus degradations) against two pretrained methods tailored for real-world defocus deblurring: Restormer [54] and INIKNet [40], against unpaired algorithm DeFlow [51] which we retrain for the current task and couple with ESRGAN (DeFlow does not provide explicit degradations hence cannot be coupled with DPIR) and against the upper baselines ESRGAN and DPIR trained with the true degradation. While Restormer and INIKNet were trained on a different degradation domain, defocus blur is also spatially varying and should not be too far from the mostly isotropic kernels of the camera PSF. Nevertheless, their inferior performance compared to our method shows how even small changes to the degradation distribution can have a large impact on reconstruction performance. Table 3 shows the quantitative evaluation on DPDD while the qualitative results are shown in Fig. 4. DeFlow does not manage to learn the correct blur distribution, and simply introduces a small amount of noise in the generated *degraded* images, thus garnering the worst results. Both INIKNet and Restormer perform similarly and succeed at removing some of the blur, but the results are not as sharp as the reconstructions obtained with our method. Both non-blind methods perform very well, with ESRGAN being better under
<figure><img src="image_8.png" /><figcaption>Figure 5: Comparing ground-truth (left) and predicted blur kernels (right), as well as the respective reconstructed images.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Figure 6: Different aberrations in parking lot data. While images taken at f/5.6 are sharp, chromatic aberrations are still present outside of the center portion as evidenced by the middle panel.</figcaption></figure>

perceptual metrics and DPIR under distortion metrics. Importantly, and like we showed for the first round of experiments on FFHQ, our method is very close to its upper limit given by the non-blind counterpart. Note that, as shown in Fig. 5, the predicted kernels are not isotropic despite our regularizer and manage to capture the spread and directional variations of the true kernels. However, there remain a few kernels such as the one at the top-left whose long tails we cannot capture well. This leads to under-compensating the blur as can be seen in the left-most bottom panels of Fig. 5.

## 4.3 Real-world camera lens calibration

As a final experiment we attempt to tackle the same lens-aberrations of the previous task, but this time on real data. We used a Panasonic DC-GX9 camera in aperture priority mode with a Leica DG Summilux 25mm f/1.4 II lens to take 22 pictures in our parking lot at different apertures. Our goal was to learn the PSF of the lens at an extreme aperture (e.g. f/16) in order to correct it using the clean distribution of images taken at a reasonable aperture (e.g. f/5.6). To avoid confounding factors we tried to keep all the image-plane in focus and had ample light to obtain sharp images. Note that while the images at different f-stops were taken using a tripod from the same location no additional care was taken to align them, and they would not be suitable for a supervised learning algorithm. Images were minimally postprocessed, by devignetting and converting to sRGB using lensfunpy [25], in order to preserve the blurry artifacts. Then they were split into patches and fed through our algorithm: an initial diffusion model was trained on the “noisy” images at f/16, and then used as a guide to learn the degradation operator on the center part of f/5.6 images. By using
<figure><img src="image_10.png" /><figcaption>Figure 7: Sample reconstructions on the parking lot dataset. Color differences in DxO are likely due to a different white-balancing algorithm. Eboli [15] correctly removes chromatic aberrations but introduces some noise artifacts. Our method results visually pleasing and significantly sharper than the original image.</figcaption></figure>

the central part of f/5.6 images we should be able to fix the chromatic aberrations which appears around sharp edges but much less in central part of images. By inspecting the data (see Fig. 6) there is a noticeable increase in blurriness between clean and noisy samples, which nevertheless is not very strong. Unlike the previous experiments where we knew the amount of additive Gaussian noise σ used in the forward model, now the noise distribution is completely unknown. For simplicity we again use a Gaussian noise model and treat its standard deviation as a trainable parameter of $\hat{A}_{\omega}$. Having a good noise estimate can be very useful in guiding the final reconstruction with plug-and-play methods. In Fig. 7 we compare our results with the algorithm from Eboli et al. [16] which is a two-step approach to first estimate and remove the blur, and then remove colored fringes using a specialized procedure [15]. For fairness we must note that the algorithm we're comparing against [16] works with single images, and doesn't need retraining for different lenses, partly thanks to strong inductive priors on the blur kernel (Gaussian with 7 parameters) and on the color aberrations. We also compare to commercial solution DxO PhotoLab [22] which exploits information about the specific lens used. We used the web version of the software and applied chromatic aberration filter followed by deblurring with strength 1.23. More comparisons as well as training details are available in the supplementary.

## 4.4 Single image super-resolution

Up to now we have dealt with the task of deblurring, matching clean and corrupted distributions across datasets. Interestingly, there is a very related task whose properties allow us to work on single images, instead of across datasets. Super-resolution is commonly modeled as the composition of blurring with kernel $k$ and subsampling $y = (x \otimes k) \downarrow_s$, hence its close relationship to deblurring. In this setting our algorithm works on single low-resolution images, split into small patches. First, the noisy distribution $Y_{\omega_*}$ is learned with a diffusion model on the patches. Then we learn a kernel $\hat{k}$ such that

$$p((y \otimes \hat{k}) \downarrow_s) \approx Y_{\omega_*} \qquad (6)$$

Note that in eq. (6) we perform a further downscaling of the already low-resolution image $y$ using the kernel $\hat{k}$ which is the target of our learning algorithm. This twice-downscaled image is compared in distribution to the once-downscaled image. Thanks to the scale-invariance of natural images the distributions (once-downscaled and twice-downscaled) will match when the kernels used for downscaling are the same, leading to the recovery of the true kernel. The downscaling step is crucial for this to work: if we omitted it (i.e. pure deblurring), we would always recover the
<table><caption>Table 4: x2 super-resolution on DIV2KRK dataset. Methods are grouped by the respective non-blind solver. NCC is the normalized cross-correlation metric.</caption><thead><tr><th>Method</th><th>Y-PSNR</th><th>SSIM</th><th>LPIPS</th><th>kernel PSNR</th><th>kernel NCC</th></tr></thead><tbody><tr><td colspan="6" style="text-align: center; font-style: italic;">End to end (kernel estimator trained on external data)</td></tr><tr><td>DANv2 [35]</td><td>32.12</td><td>0.8954</td><td>0.148</td><td>56.8</td><td>0.97</td></tr><tr><td>DCLS [33]</td><td>32.31</td><td>0.9006</td><td>0.147</td><td>—</td><td>—</td></tr><tr><td colspan="6" style="text-align: center; font-style: italic;">Non-blind (true kernel available)</td></tr><tr><td>ZSSR</td><td>32.04</td><td>0.8885</td><td>0.196</td><td>∞</td><td>1</td></tr><tr><td>USRNet</td><td>32.06</td><td>0.8923</td><td>0.158</td><td>∞</td><td>1</td></tr><tr><td>DANv2</td><td>31.94</td><td>0.8957</td><td>0.147</td><td>∞</td><td>1</td></tr><tr><td colspan="6" style="text-align: center; font-style: italic;">Two-step (kernel estimate from single image)</td></tr><tr><td>KernelGAN [6]</td><td>29.65</td><td>0.8465</td><td>0.248</td><td>51.5</td><td>0.92</td></tr><tr><td>Ours + ZSSR</td><td>31.61</td><td>0.8837</td><td>0.200</td><td>64.2</td><td>0.96</td></tr><tr><td>DKP [52]</td><td>23.18</td><td>0.6245</td><td>0.226</td><td>53.6</td><td>0.92</td></tr><tr><td>DIP-FKP [27]</td><td>29.52</td><td>0.8530</td><td>0.249</td><td>50.4</td><td>0.95</td></tr><tr><td>Ours + USRNet</td><td>31.11</td><td>0.8838</td><td>0.167</td><td>64.2</td><td>0.96</td></tr><tr><td>Ours + DANv2</td><td>31.86</td><td>0.8938</td><td>0.149</td><td>64.2</td><td>0.96</td></tr></tbody></table>

identity kernel. For a quantitative experiment on single image super-resolution we adopt the setting from Bell-Kligler et al. [6] who also introduced the DIV2KRK dataset consisting of 100 images from DIV2K [3] downscaled with random anisotropic Gaussian kernels. After learning $\hat{k}$ using the proposed algorithm, we tested three non-blind solvers: ZSSR [4] which requires no pretraining and USRNet [56] for which we used the *tiny* pretrained model and DANv2 [35] from which we took just the *Restorer* module plugged with $\hat{k}$. In Table 4 we compare our approach with: i) Two end-to-end solutions [33, 35] which jointly learn super-resolution and kernel-estimation by pre-training on DIV2K + Flickr2K [46] datasets, ii) two non-blind super-resolution algorithms [4, 56] which are natural upper-bounds for our method and iii) three kernel-estimation algorithms [6, 27, 52] which are directly comparable to the proposed method. Note that both DANv2 and DCLS outperform the non-blind methods. We hypothesize this may be due to the strong prior the methods have learned during training on a dataset which is very similar to the one used in testing. To analyze kernel estimation performance in isolation we look at kernel metrics k-PSNR and k-NCC (not available for DCLS which operates in a different kernel space) which show that the proposed algorithm performs the same as DANv2. This can be further corroborated by plugging in the learned kernels with the DANv2 super-resolution module on its own, which significantly reduces the performance gap to the end-to-end DANv2. When compared to the other 2-step methods, the proposed algorithm emerges as the clear winner both on image and kernel metrics.

# 5 Conclusions

Solving blind inverse problems from unpaired data is highly relevant for practical real-world image restoration tasks. In this work, we proposed an algorithm which learns the degradation operator directly from data, making minimal assumptions about the operator itself. Through extensive experiments on uniform and non-uniform deblurring, we have shown that it is possible to reduce the performance gap between blind single-image methods and non-blind methods by leveraging multiple images from the same degradation process. Built-in interpretability and flexibility in its design make the proposed algorithm applicable to solving diverse problems. On the other hand, having multiple images belonging to the same noisy-data distribution may not be realistic for all situations, where single-image methods work best – this single-image setting can be handled with our methods only in specific problems such as super-resolution. Furthermore, while our framework is designed to deal with unpaired data, if the clean and noisy distributions are known through highly unrelated
datasets, the matching procedure will struggle to learn a degradation operator to reconcile the two. The domain shift between clean and noisy distributions hinders training as it becomes larger. Future directions include tackling problems from different domains such as microscopy or medical imaging, as well improving the noise-modeling to be more realistic (e.g. to handle Poisson-Gaussian noise).

## Acknowledgements

This work was supported by ERC grant number 101087696 (APHELEIA project). This work was granted access to the HPC resources of IDRIS under the allocation AD011015445 made by GENCI.

## References

[1] Abdullah Abuolaim and Michael S Brown. Defocus deblurring using dual-pixel data. In *European Conference on Computer Vision*, pages 111–126. Springer, 2020. 8

[2] Abdullah Abuolaim, Mahmoud Afifi, and Michael S Brown. Improving single-image defocus deblurring: How dual-pixel images help through multi-task learning. In *Proceedings of the IEEE/CVF Winter Conference on Applications of Computer Vision*, pages 1231–1239, 2022. 4

[3] Eirikur Agustsson and Radu Timofte. NTIRE 2017 challenge on single image super-resolution: Dataset and study. In *CVPRW*, 2017. 8, 12

[4] Michal Irani Assaf Shocher, Nadav Cohen. "zero-shot" super-resolution using deep internal learning. In *The IEEE Conference on Computer Vision and Pattern Recognition (CVPR)*, 2018. 12

[5] Matthias Bauer, Valentin Volchkov, Michael Hirsch, and Bernhard Schölkopf. Automatic estimation of modulation transfer functions. In *IEEE International Conference on Computational Photography (ICCP)*, pages 1–12, 2018. 4, 5, 8

[6] Sefi Bell-Kligler, Assaf Shocher, and Michal Irani. Blind super-resolution kernel estimation using an internal-gan. In *International Conference on Neural Information Processing Systems*, 2019. 4, 12, 24

[7] Levi Borodenko. Motion blur, 2020. 7, 8, 18

[8] Andrew Brock, Jeff Donahue, and Karen Simonyan. Large scale GAN training for high fidelity natural image synthesis. In *International Conference on Learning Representations*, 2019. 4

[9] Adrian Bulat and Jing Yang and Georgios Tzimiropoulos. To learn image super-resolution, use a GAN to learn how to do image degradation first. In *ECCV*, 2018. 3

[10] Joonyoung Chang, Hee Kang, and Moon Gi Kang. Correction of axial and lateral chromatic aberration with false color filtering. *IEEE Transactions on Image Processing*, 22(3):1186–1198, 2013. 5

[11] Shuaijun Chen, Zhen Han, Enyan Dai, Xu Jia, Ziluan Liu, Xing Liu et al. Unsupervised image super-resolution with an indirect supervised path. In *2020 IEEE/CVF Conference on Computer Vision and Pattern Recognition Workshops (CVPRW)*, pages 1924–1933, 2020. 4

[12] Hyungjin Chung, Jeongsol Kim, Sehui Kim, and Jong Chul Ye. Parallel diffusion models of operator and image for blind inverse problems. In *IEEE/CVF Conference on Computer Vision and Pattern Recognition*, 2023. 3, 7, 8, 18, 19

[13] Dave Coffin. ddraw, 2018. https://dechifro.org/ddraw/ v9.28, accessed: 2025-05-29. 20

[14] Jiangxin Dong, Stefan Roth, and Bernt Schiele. Deep wiener deconvolution: Wiener meets deep learning for image deblurring. *Advances in Neural Information Processing Systems*, 33: 1048–1059, 2020. 7
[15] Thomas Eboli. Fast chromatic aberration correction with 1d filters. Image Process. Line, 13: 198-214, 2023. 5, 11, 22

[16] Thomas Eboli, Jean-Michel Morel, and Gabriele Facciolo. Fast two-step blind optical aberration correction. In European Conference on Computer Vision, pages 693-708, 2022. 5, 11, 20

[17] Angela Gao, Jorge Castellanos, Yisong Yue, Zachary Ross, and Katherine Bouman. Deepgem: Generalized expectation-maximization for blind inversion. In Advances in Neural Information Processing Systems, 2021. 3

[18] Jinjin Gu, Hannan Lu, Wangmeng Zuo, and Chao Dong. Blind super-resolution with iterative kernel correction. In The IEEE Conference on Computer Vision and Pattern Recognition (CVPR), 2019. 2, 3, 4

[19] Tero Karras, Samuli Laine, and Timo Aila. A style-based generator architecture for generative adversarial networks. IEEE Transactions on Pattern Analysis and Machine Intelligence, 43 (12):4217-4228, 2021. 7

[20] Tero Karras, Miika Aittala, Jaakko Lehtinen, Janne Hellsten, Timo Aila, and Samuli Laine. Analyzing and improving the training dynamics of diffusion models. In Proc. CVPR, 2024. 18

[21] Eric Kee, Sylvain Paris, Simon Chen, and Jue Wang. Modeling and removing spatially-varying optical blur. In IEEE international conference on computational photography (ICCP), pages 1-8, 2011. 4

[22] DxO Labs. Dxo photolab 8, 2025. https://www.dxo.com/fr/dxo-photolab/ accessed: 2025-05-29. 11

[23] Charles Laroche, Andrés Almansa, and Eva Coupete. Fast diffusion em: a diffusion model for blind inverse problems with application to deconvolution. In Winter Conference on Applications of Computer Vision, 2024. 3, 7, 8, 18

[24] Junyong Lee, Hyeongseok Son, Jaesung Rim, Sunghyun Cho, and Seungyong Lee. Iterative filter adaptive network for single image defocus deblurring. In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition (CVPR), 2021. 4

[25] Lensfun authors. Lensfun, 2022. v0.3.3. 10

[26] Xiu Li, Jinli Suo, Weihang Zhang, Xin Yuan, and Qionghai Dai. Universal and flexible optical aberration correction using deep-prior based deconvolution. In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 2613-2621, 2021. 8

[27] Jingyun Liang, Kai Zhang, Shuhang Gu, Luc Van Gool, and Radu Timofte. Flow-based kernel prior with application to blind super-resolution. In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, pages 10601-10610, 2021. 3, 4, 12

[28] Yaron Lipman, Ricky T. Q. Chen, Heli Ben-Hamu, Maximilian Nickel, and Matt Le. Flow matching for generative modeling. In International Conference on Learning Representations, 2023. 4, 5, 6

[29] Xingchao Liu, Chengyue Gong, and Qiang Liu. Flow straight and fast: Learning to generate and transfer data with rectified flow. In International Conference on Learning Representations, 2022. 4, 5, 6

[30] Andreas Lugmayr, Martin Danelljan, and Radu Timofte. Unsupervised learning for real-world super-resolution. In 2019 IEEE/CVF International Conference on Computer Vision Workshop (ICCVW), pages 3408-3416, 2019. 3
[31] Andreas Lugmayr, Martin Danelljan, Radu Timofte, Namhyuk Ahn, Dongwoon Bai, Jie Cai et al. NTIRE 2020 challenge on real-world image super-resolution: Methods and results. In *The IEEE Conference on Computer Vision and Pattern Recognition (CVPR) Workshops*, 2020. 2

[32] Weijian Luo, Tianyang Hu, Shifeng Zhang, Jiacheng Sun, Zhenguo Li, and Zhihua Zhang. Diff-instruct: A universal approach for transferring knowledge from pre-trained diffusion models. In *Thirty-seventh Conference on Neural Information Processing Systems*, 2023. 4, 6

[33] Ziwei Luo, Haibin Huang, Lei Yu, Youwei Li, Haoqiang Fan, and Shuaicheng Liu. Deep constrained least squares for blind image super-resolution. In *Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition*, pages 17642–17652, 2022. 12, 24

[34] Zhengxiong Luo, Yan Huang, Shang Li, Liang Wang, and Tieniu Tan. Learning the degradation distribution for blind image super-resolution. In *CVPR*, 2022. 3, 4

[35] Zhengxiong Luo, Yan Huang, Shang Li, Liang Wang, and Tieniu Tan. End-to-end alternating optimization for real-world blind super resolution. *Int. J. Comput. Vision*, 131(12):3152–3169, 2023. 4, 12, 24

[36] Shunta Maeda. Unpaired image super-resolution using pseudo-supervision. In *2020 IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR)*, 2020. 3

[37] Subhadip Mukherjee, Marcello Carioni, Ozan Öktem, and Carola-Bibiane Schönlieb. End-to-end reconstruction meets data-driven regularization for inverse problems. In *Advances in Neural Information Processing Systems*, 2021. 3

[38] Naoki Murata, Koichi Saito, Chieh-Hsin Lai, Yuhta Takida, Toshimitsu Uesaka, Yuki Mitsufuji et al. GibbsDDRM: A partially collapsed Gibbs sampler for solving blind inverse problems with denoising diffusion restoration. In *International Conference on Machine Learning*, 2023. 3

[39] Ben Poole, Ajay Jain, Jonathan T. Barron, and Ben Mildenhall. Dreamfusion: text-to-3d using 2d diffusion. In *International Conference on Learning Representations*, 2023. 4, 6

[40] Yuhui Quan, Xin Yao, and Hui Ji. Single image defocus deblurring via implicit neural inverse kernels. In *Proceedings of the IEEE/CVF International Conference on Computer Vision*, pages 12600–12610, 2023. 4, 9

[41] Dongwei Ren, Kai Zhang, Qilong Wang, Qinghua Hu, and Wangmeng Zuo. Neural blind deconvolution using deep priors. In *IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR)*, 2020. 3

[42] Yaniv Romano, Michael Elad, and Peyman Milanfar. The little engine that could: Regularization by denoising (RED). *SIAM Journal on Imaging Sciences*, 10(4):1804–1844, 2017. 3

[43] Andrés Romero, Luc Van Gool, and Radu Timofte. Unpaired real-world super-resolution with pseudo controllable restoration. In *2022 IEEE/CVF Conference on Computer Vision and Pattern Recognition Workshops (CVPRW)*, pages 797–806, 2022. 4

[44] Yash Sanghvi, Yiheng Chi, and Stanley H. Chan. Kernel diffusion: An alternate approach to blind deconvolution. In *European Conference on Computer Vision (ECCV)*, pages 1–20, 2024. 3, 7, 8, 18

[45] Byeongsu Sim, Gyutaek Oh, Jeongsol Kim, Chanyong Jung, and Jong Chul Ye. Optimal transport driven CycleGAN for unsupervised learning in inverse problems. *SIAM Journal on Imaging Sciences*, 13(4):2281–2306, 2020. 2, 3, 4

[46] Radu Timofte, Eirikur Agustsson, Luc Van Gool, Ming-Hsuan Yang, Lei Zhang, Bee Lim et al. NTIRE 2017 challenge on single image super-resolution: Methods and results. In *CVPRW*, 2017. 12
[47] Singanallur V. Venkatakrishnan, Charles A. Bouman, and Brendt Wohlberg. Plug-and-play priors for model based reconstruction. In *Global Conference on Signal and Information Processing*, pages 945–948, 2013. [3, 7]

[48] Ana Fernandez Vidal, Valentin De Bortoli, Marcelo Pereyra, and Alain Durmus. Maximum likelihood estimation of regularization parameters in high-dimensional inverse problems: An empirical bayesian approach part i: Methodology and experiments. *SIAM Journal on Imaging Sciences*, 13(4):1945–1989, 2020. [3]

[49] Xintao Wang, Ke Yu, Shixiang Wu, Jinjin Gu, Yihao Liu, Chao Dong et al. ESRGAN: Enhanced super-resolution generative adversarial networks. In *The European Conference on Computer Vision Workshops (ECCVW)*, 2018. [3, 7, 9]

[50] Xintao Wang, Liangbin Xie, Chao Dong, and Ying Shan. Real-esrgan: Training real-world blind super-resolution with pure synthetic data. In *International Conference on Computer Vision Workshops (ICCVW)*, 2021. [3]

[51] Valentin Wolf, Andreas Lugmayr, Martin Danelljan, Luc Van Gool, and Radu Timofte. DeFlow: Learning complex image degradations from unpaired data with conditional flows. In *IEEE/CVF Conference on Computer Vision and Pattern Recognition*, CVPR, 2021. [2, 3, 4, 9, 20]

[52] Zhixiong Yang, Jingyuan Xia, Shengxi Li, Xinghua Huang, Shuanghui Zhang, Zhen Liu et al. A dynamic kernel prior model for unsupervised blind image super-resolution. In *Proceedings of the IEEE/CVF conference on computer vision and pattern recognition*, pages 26046–26056, 2024. [3, 12]

[53] Yuan Yuan, Siyuan Liu, Jiawei Zhang, Yongbing Zhang, Chao Dong, and Liang Lin. Unsupervised image super-resolution using cycle-in-cycle generative adversarial networks. In *2018 IEEE/CVF Conference on Computer Vision and Pattern Recognition Workshops (CVPRW)*, 2018. [3]

[54] Syed Waqas Zamir, Aditya Arora, Salman Khan, Munawar Hayat, Fahad Shahbaz Khan, and Ming-Hsuan Yang. Restormer: Efficient transformer for high-resolution image restoration. In *CVPR*, 2022. [4, 7, 9]

[55] Kai Zhang, Wangmeng Zuo, Shuhang Gu, and Lei Zhang. Learning deep CNN denoiser prior for image restoration. In *IEEE Conference on Computer Vision and Pattern Recognition*, pages 3929–3938, 2017. [7, 9]

[56] Kai Zhang, Luc Van Gool, and Radu Timofte. Deep unfolding network for image super-resolution. In *IEEE Conference on Computer Vision and Pattern Recognition*, pages 3217–3226, 2020. [12]

[57] Kai Zhang, Yawei Li, Wangmeng Zuo, Lei Zhang, Luc Van Gool, and Radu Timofte. Plug-and-play image restoration with deep denoiser prior. *IEEE Transactions on Pattern Analysis and Machine Intelligence*, 44(10):6360–6376, 2021. [3, 21]

[58] Kai Zhang, Jingyun Liang, Luc Van Gool, and Radu Timofte. Designing a practical degradation model for deep blind image super-resolution. In *IEEE International Conference on Computer Vision*, pages 4791–4800, 2021. [3]

[59] Ruofan Zhang, Jinjin Gu, Haoyu Chen, Chao Dong, Yulun Zhang, and Wenming Yang. Crafting training degradation distribution for the accuracy-generalization trade-off in real-world super-resolution. In *Proceedings of the 40th International Conference on Machine Learning*, 2023. [2, 3]

[60] Yasi Zhang, Peiyu Yu, Yaxuan Zhu, Yingshan Chang, Feng Gao, Ying Nian Wu et al. Flow priors for linear inverse problems via iterative corrupted trajectory matching. In *Neural Information Processing Systems (NeurIPS)*, 2024. [5]
[61] Jun-Yan Zhu, Taesung Park, Phillip Isola, and Alexei A. Efros. Unpaired image-to-image translation using cycle-consistent adversarial networks. In *IEEE International Conference on Computer Vision (ICCV)*, 2017. 3

[62] Yuanzhi Zhu, Kai Zhang, Jingyun Liang, Jiezhang Cao, Bihan Wen, Radu Timofte et al. Denoising diffusion models for plug-and-play image restoration. In *IEEE Conference on Computer Vision and Pattern Recognition Workshops (NTIRE)*, 2023. 3, 7, 18

# Appendices

## A Proof of Proposition 1.

We restate the statement of the proposition for completeness.

**Proposition A.1:** For any set of forward model parameters $\omega$, let $p_\omega(y) = \int p_\omega(y | x)p(x)dx$ where $p_\omega(y | x) = N(y | A_\omega x, \sigma^2 I)$. Let $\omega_*$ be a specific set of parameters which we consider to be the optimal set. Then, assuming the data covariance $\Sigma = \mathbb{E}_x[xx^\top]$ is invertible, there exists an orthogonal matrix $P$ such that

$$p_{\omega}(y) = p_{\omega_{*}}(y) \implies \mathcal{A}_{\omega} = \Sigma^{-1/2} P \Sigma^{1/2} \mathcal{A}_{\omega_{*}} \quad (7)$$

That is, if the probability distributions $p_\omega$ and $p_{\omega_*}$ are equal, it is possible to identify $\omega_*$ up to rotations $P$.

*Proof.* Let $f$ be a quadratic function of the form $f(x) = x^\top C x$. Consider the form of the conditional, and perform a change of variables from $y$ to $Ax + \epsilon$ with $\epsilon \sim N(0, \sigma^2 I)$ independent of $x$ to write:

$$\begin{aligned}
0 &= \mathbb{E}_{y \sim p_{\omega}}[f(y)] - \mathbb{E}_{y \sim p_{\omega_*}}[f(y)] \\
&= \int f(y) p_{\omega}(y | x) p(x) dx dy - \int f(y) p_{\omega_*}(y | x) p(x) dx dy \\
&= \int f(\mathcal{A}_{\omega} x + \epsilon) p(\epsilon) p(x) dx d\epsilon - \int f(\mathcal{A}_{\omega_*} x + \epsilon) p(\epsilon) p(x) dx d\epsilon \\
&= \mathbb{E}_{x \sim p(x), \epsilon}[f(\mathcal{A}_{\omega} x + \epsilon) - f(\mathcal{A}_{\omega_*} x + \epsilon)].
\end{aligned}$$

Now

$$\begin{aligned}
\mathbb{E}_{x \sim p(x), \epsilon}[f(\mathcal{A}_{\omega} x + \epsilon) - f(\mathcal{A}_{\omega_*} x + \epsilon)] &= \mathbb{E}_{x \sim p(x), \epsilon}[(\mathcal{A}_{\omega} x + \epsilon)^{\top} C (\mathcal{A}_{\omega} x + \epsilon) - (\mathcal{A}_{\omega_*} x + \epsilon)^{\top} C (\mathcal{A}_{\omega_*} x + \epsilon)] \\
&= \mathbb{E}_{x, \epsilon}[x^{\top} \mathcal{A}_{\omega}^{\top} C \mathcal{A}_{\omega} x - x^{\top} \mathcal{A}_{\omega_*}^{\top} C \mathcal{A}_{\omega_*} x + 2x^{\top} \mathcal{A}_{\omega}^{\top} C \epsilon - 2x^{\top} \mathcal{A}_{\omega_*}^{\top} C \epsilon] \\
&= \mathbb{E}_{x}[x^{\top} \mathcal{A}_{\omega}^{\top} C \mathcal{A}_{\omega} x - x^{\top} \mathcal{A}_{\omega_*}^{\top} C \mathcal{A}_{\omega_*} x] \\
&= \mathbb{E}_{x}[\text{tr}(x^{\top} \mathcal{A}_{\omega}^{\top} C \mathcal{A}_{\omega} x) - \text{tr}(x^{\top} \mathcal{A}_{\omega_*}^{\top} C \mathcal{A}_{\omega_*} x)] \\
&= \mathbb{E}_{x}[\text{tr}(\mathcal{A}_{\omega} x x^{\top} \mathcal{A}_{\omega}^{\top} C) - \text{tr}(\mathcal{A}_{\omega_*} x x^{\top} \mathcal{A}_{\omega_*}^{\top} C)] \\
&= \langle (\mathcal{A}_{\omega} \Sigma \mathcal{A}_{\omega}^{\top} - \mathcal{A}_{\omega_*} \Sigma \mathcal{A}_{\omega_*}^{\top}), C \rangle_F = 0 \\
&\implies \mathcal{A}_{\omega} \Sigma \mathcal{A}_{\omega}^{\top} - \mathcal{A}_{\omega_*} \Sigma \mathcal{A}_{\omega_*}^{\top} = 0 \quad (8)
\end{aligned}$$

where $\Sigma = \mathbb{E}_x[xx^\top]$ is the data covariance. Eq. (8) holds since $C$ can be any matrix. Now since $\Sigma$ is invertible, we can take $\tilde{\mathcal{A}}_\omega = \mathcal{A}_\omega \Sigma^{1/2}$ and $\tilde{\mathcal{A}}_{\omega_*} = \mathcal{A}_{\omega_*} \Sigma^{1/2}$ and we have that $\tilde{\mathcal{A}}_\omega \tilde{\mathcal{A}}_\omega^\top = \tilde{\mathcal{A}}_{\omega_*} \tilde{\mathcal{A}}_{\omega_*}^\top$. Now looking at the SVD of $\tilde{\mathcal{A}}_\omega = U_\omega S_\omega V_\omega^\top$ and $\tilde{\mathcal{A}}_{\omega_*} = U_{\omega_*} S_{\omega_*} V_{\omega_*}^\top$ it holds that

$$U_{\omega} S_{\omega}^{2} U_{\omega}^{\top} = U_{\omega_{*}} S_{\omega_{*}}^{2} U_{\omega_{*}}^{\top} \quad (9)$$
<figure><img src="image_11.png" /><figcaption>Figure 8: Sample reconstructions on FFHQ (motion-blur kernel). All methods single-image apart from ours.</figcaption></figure>

which implies that both $U_\omega = U_{\omega_*}$ and $S_\omega = S_{\omega_*}$. Take $P$ the orthonormal change of basis from $V_\omega$ to $V_{\omega_*}$, then

$$ \tilde{A}_\omega = \sum_i \sigma_\omega^{(i)} u_\omega^{(i)} v_\omega^{(i)*} P = \sum_i \sigma_{\omega_*}^{(i)} u_{\omega_*}^{(i)} v_\omega^{(i)*} P = \sum_i \sigma_{\omega_*}^{(i)} u_{\omega_*}^{(i)} v_{\omega_*}^{(i)*} = \tilde{A}_{\omega_*} \quad (10) $$

and going back to the original operators, $A_\omega = A_{\omega_*} \Sigma^{1/2} P \Sigma^{-1/2}$ with $P$ an orthogonal matrix. $\square$

# B Experiment Details

For all experiments we started from the diffusion model implementation EDM2 [20], and adapted it to train flow-matching models. The EMA for model weights based on a power-function instead of an exponential and is parameterized in an unconventional way through the relative standard deviation (see Karras et al. [20] for a precise description).

## B.1 FFHQ experiment

We used a very small CFM model with 4M parameters to train on the first 1000 samples from the FFHQ dataset, downscaled to 256x256, degraded using a fixed motion-blur kernel generated with the code from Borodenko [7] with intensity 0.5 and Gaussian noise with standard deviation 0.02. We trained the model from scratch for 4000 epochs, at the end of which it could generate samples like those in Fig. 3. We then proceeded to the second step of our algorithm by learning a 28x28 blur kernel from 100 different clean samples from FFHQ. We enforced the kernel to sum to 1, and used a sparsity promoting regularizer (simply average of the absolute values of the kernel) with weight 1. We trained the second step for 5000 epochs with a low learning rate ($10^{-5}$ for the blur kernel and $4 \times 10^{-5}$ for the auxiliary CFM model). Finally we used DiffPIR [62] in the implementation of DeepInv (https://github.com/deepinv/deepinv) for reconstructing the full FFHQ test set (1000 images). The baselines were all run from their respective repositories. The only difficulty encountered was with KernelDiff [44] which cannot handle noisy measurements. We thus ran KernelDiff without noise. Figures 8 and 9 compare the image and kernel reconstructions of the different methods for both motion-blur and isotropic Gaussian kernels. Note that the same regularization and other hyperparameters were used to train our method in both experiments.

### B.1.1 FFHQ with varying noise levels

The experiments described above and in the main text all used the same random Gaussian noise with standard deviation of 0.02. To determine whether our method is robust to different levels of noise we ran the same experiment with noise-levels 0, 0.04 and 0.08. To compare the models fairly, we then ran DiffPIR using the learned kernels (learned with different noise levels) on data which has been corrupted with the original noise level of 0.02. The training hyper-parameters for all three stages were kept fixed from the previous section. The results of this experiment show a small decrease in accuracy as the noise-level increases, in Section B.1.1 support a small deterioration in both kernel estimates and consequently in the reconstruction quality, while remaining fairly robust.
<figure><img src="image_12.png" /><figcaption>Figure 9: Sample reconstructions on FFHQ (gaussian-blur kernel). All methods single-image apart from ours.</figcaption></figure>

<figure><img src="image_13.png" /></figure>

<table><caption>Figure 10: Qualitative and quantitative results for motion-deblurring with increasing amounts of noise. The figure on the left shows the learned kernel and a sample reconstruction with the different noise levels. On the right, the table provides reconstruction metrics and kernel-accuracy metrics for the FFHQ test-set (1000 random images). Note that kernel NCC is the normalized cross-correlation of the learned kernel with the ground-truth kernel.</caption><thead><tr><th>σ<sub>train</sub></th><th>PSNR</th><th>LPIPS</th><th>kernel NCC</th></tr></thead><tbody><tr><td>0</td><td>28.93</td><td>0.068</td><td>0.98</td></tr><tr><td>0.02</td><td>28.76</td><td>0.069</td><td>0.94</td></tr><tr><td>0.04</td><td>28.84</td><td>0.070</td><td>0.96</td></tr><tr><td>0.08</td><td>28.01</td><td>0.079</td><td>0.91</td></tr></tbody></table>

Note that the noise is applied to 0-1 normalized images such that a noise level of 0.08 corresponds to 8%, which is very clearly noticeable in the corrupted images.

## B.2 Space-varying blur on DPDD

Here we used a larger CFM model with 52M parameters, pretrained on the DIV2K clean training set. This then helps speed up the first step of our algorithm which consists in fine-tuning the model on the same dataset, degraded using the per-pixel blur kernels, as linearly interpolated from the 8x8 PSFs shown in Fig. 5 and Gaussian noise with standard deviation 0.01. The model was trained on patches of size 128x128 from the degraded DIV2K dataset with learning rate $10^{-4}$ and exponential moving average over the weights. For the second step, we once again directly optimized the parameters of an 8x8 grid of kernels (the same as the ground-truth). A learning rate scheduler is used for the auxiliary flow matching model (square root scheduling with warmup). A batch size of 20 was used, and we found the batch size to significantly affect convergence speed (i.e. convergence depended on how many times each of the two models was updated, so higher batch sizes did not speed up much). For this step, we use a regularization composed of three different terms, the first one is a center term which avoid the kernels to have centering offset, the second one is a sparsity one which is basically a 11 regularization which limit the number of non zero pixels in the kernels and the last one is a Gaussian prior which help with the overall shape of kernels. Each term has the same weight of 0.003. For the third step we used the DPDD dataset and perform the forward model inversion in two ways:
<figure><img src="image_14.png" /><figcaption>Figure 11: Full restored DPDD image for each method analyzed</figcaption></figure>

* ESRGAN, used by generating a paired dataset on the “train_c” split of DPDD. We used the implementation at https://github.com/XPixelGroup/BasicSR, changed the upscaling factor to 1 since we were not interested in super-resolution and reduced the number of iterations to 60000 since convergence was fast.

* DPIR, again using the implementation from DeepInv, we only need to process the DPDD test-set.

To train the Deflow [51] we use the official repository which is available on github at https://github.com/volflow/DeFlow. Deflow requires a bit of preprocessing to create the training datasets. Indeed, we use both DIV2K and DPDD to train this model. So we use DPDD as a clean dataset with the low quality dataset being the DPDD training set downsampled by 4 using BICUBIC interpolation. And we use DIV2K for the noisy dataset applying a per pixel kernel interpolated from the 8x8 grid on which we add a gaussian noise of standard deviation 0.01. We also downsample the noisy DIV2K to have noisy low quality dataset. Concerning the parameters of the training we only change the batch size in order to make it fit our GPUs and we change the scale to 4 because we only found a pretrained model for this downsampling scale.

In Fig. 11, we show the full reconstructed images used for Fig. 4. The first 2 images are methods which use the real degradation which are our baselines. The 3 following images come from the other methods we compare ours to and the last image is our method.

## B.3 Parking-lot experiment

**Dataset** We used a Panasonic DC-GX9 camera in aperture priority mode with a Leica DG Summilux 25mm f/1.4 II lens to take 22 distinct photographs in a parking lot at apertures f/16 and f/5.6. We used a tripod for stabilisation but the two images are still not perfectly aligned (for this reason we do not provide any quantitative metric on experiments using this dataset). We process the raw images using ddraw [13] to fix white balance, image highlights, demosaicking, and a small amount of noise reduction. We then converted images to png format. A low-resolution sample of the images in the dataset is shown in Fig. 12. The high-resolution counterparts have 5200x3904 pixels.

**Training** We started from a scratch CFM model as the size of the one for DDPD experiment, we train it on 64-sized patches from the 18 parking lot images at f-stop 16. For the distribution-matching step we used the central 768x768 pixels of f/5.6 images as the clean distribution and learn a 8x8 grid of 13x13 RGB kernels. We compare predicted kernels with Eboli et al. [16] in Figs. 13 and 14.
<figure><img src="image_15.png" /><figcaption>Figure 12: Samples from our parking-lot dataset. All 6 images are at f-stop 5.6</figcaption></figure>

Our PSFs are consistent across images by design, while the competing method adapts between different images. The PSFs predicted by our method more accurately recover the true blur, as seen in Figs. 7 and 15.

**Hyperparameters** The first step consists in learning a diffusion model on a single, low-resolution image. We used a flow-matching model with about 4M parameters, trained on patches of size 64x64. We used a batch-size of 1024 patches, without any added noise, learning rate of 0.01 and 0.05 as an EMA value. We trained for 16M patches. To learn a different blur filter for different locations in the image plane, we condition the diffusion model on the original patch's location in the full image. We do this by concatenating the input image with two positioning-channels (for the x and y axes). Since the dataset is small, to avoid completely overfitting the diffusion model to the location-specific conditioning, we replace the conditioning for 20% of the patches with a randomly chosen location. For the second step we used a smaller batch size of 128, learning rate of 0.00001 for the auxiliary diffusion model and 0.00004 for the kernel network (described in the next paragraph). We added sparsity regularization (11) with strength 0.01 and gaussian prior regularization with strength 0.1. The second step was trained for 6M patches. Here we wished to learn a location-conditioned degradation which, when applied to any part of the central portion of a clean image, it would replicate the degradations of a specific location of the corrupted images. To achieve this, we conditioned the step-2 model (whose inputs are patches from the central portion of clean images) on random patch locations (spanning the full image size). Like for DPDD we parameterized the degradations by a 8x8 grid of kernels, which unlike DPDD are different for the three RGB channels (to represent chromatic aberrations accurately). Final filters for reconstruction (which is patch-wise) are linearly interpolated from this 8x8 grid. The learned PSF grid is shown in Fig. 13. We additionally added a trainable parameter to learn the noise standard deviation of the images, which was assumed to be fixed throughout the image. The noise estimate was low (around 0.003) because similar amounts of noise were present in the clean and corrupted images. For the third step, we used plug-and-play method DPIR [57] which uses a deep natural image prior. We found the results were very similar to classical Wiener deconvolution, apart for some extra noise removal with DPIR. Additional qualitative reconstructions are shown in Fig. 15.

## B.4 Super Resolution

**Hyperparameters** The first step consists in learning a diffusion model on a single, low-resolution image. We used a small flow-matching model with about 4M parameters, trained on patches of size 32x32. We used a batch-size of 512 patches, without any added noise, learning rate of 0.01 and a
<figure><img src="image_16.png" /><figcaption>Figure 13: Lens PSF predicted by our model brightened x2. Note how the chromatic aberration (different position of red and blue channels) varies smoothly between the top, bottom and left, right parts of the lens.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Figure 14: Lens PSF by Eboli [15] for two different images. Note that this method corrects for color aberration in a separate successive step. The PSFs here are not brightened, and hence much smaller than those in Fig. 13.</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>Figure 15: Additional reconstructions on parking lot data (best viewed zoomed-in).</figcaption></figure>
small EMA value. We trained for 1M images (or approximately 2000 steps). For the second step we used a smaller batch size of 64, learning rate of 0.00001 for the auxiliary diffusion model and 0.00008 for the kernel network (described in the next paragraph). We added center-regularization with strength 0.01, sparsity regularization (l1) with strength 0.1 and sum-to-one regularization with strength 0.1. Note that the kernel network was not forced to output normalized kernels (unlike all previous experiments), so the sum-to-one regularization was important. The second step was trained for 512k images (or approximately 8000 steps).

**Kernel parameterization** Similarly to KernelGAN [6] and DCLS [33] we use a linear convolutional network to model the blur kernel. This network can be used directly as forward model $A_\omega$, and can be collapsed into an explicit kernel by convolving a dirac function (i.e. an image with a single non-zero pixel in the center). We could confirm the findings of Bell-Kligler et al. [6] that this parameterization is easier to train than an equivalent CNN with a single layer. We used 4 layers with kernel-sizes 7, 5, 5 and 1, (for an overall receptive field of 15) no bias and 64 channels.

**Additional results** In Fig. 16 we show zoom-ins on four different images (35, 36, 37 and 39) from the DIV2KRK dataset as reconstructed with different methods, as well as the blur kernels predicted by the respective methods. We can see that while DANv2 [35] has better reconstruction, the predicted kernels do not show significant differences to the ones predicted with our method. The better reconstruction is due to the relatively lower accuracy of ZSSR compared to the integrated super-resolver in DANv2. FKP has decent reconstruction because it always predicts blur kernels which are close to being isotropic Gaussians with a small standard deviation, while DKP predicts kernels with a higher variance which are not very close to the true ones.

In Fig. 17 we plot the distribution of the difference in PSNR between blind SR methods and their non-blind counterpart. Note that this can be greater than zero when due to random stochasticity of the non-blind method, the estimated kernel results in a better reconstruction than the true one.

**Metrics** We specify here the parameters used to calculate metrics in Table 4. Firstly, we found USRNet to produce very strong border artifacts which – if included – would have completely skewed the results. For this reason we crop 60 pixels off each edge for every image before computing metrics. For the PSNR, in order to maintain consistency with results reported in the literature, we used the Y-PSNR (i.e. computed on just the luminance channel, after converting RGB images to YCbCr). For the LPIPS metric we used the "alex-net" network. The PSNR on kernels was computed after padding the kernels with zeros up to size 25x25 and after shifting them so they were appropriately centered. The NCC metric is the 2d normalized cross-correlation, which can better account for small centering errors in the kernels.
<figure><img src="image_19.png" /><figcaption>Figure 16: Super-resolution: qualitative results. Note that DCLS kernels are not comparable as they live in a different space.</figcaption></figure>
<figure><img src="image_20.png" /><figcaption>Figure 17: Super-resolution: Performance gap of two-step methods compared to their natural upper-bound. The metric considered is the Y-PSNR (PSNR on the luminance only).</figcaption></figure>