Bayesian Rain Field Reconstruction using Commercial Microwave Links and Diffusion Model Priors
Summary
This paper presents a Bayesian inverse problem framework for rain field reconstruction using Commercial Microwave Links and Diffusion Model priors, demonstrating improved accuracy over existing baselines.
View Cached Full Text
Cached at: 05/08/26, 07:49 AM
# Bayesian Rain Field Reconstruction using Commercial Microwave Links and Diffusion Model Priors
Source: [https://arxiv.org/html/2605.05520](https://arxiv.org/html/2605.05520)
Albina IlinaHai Victor HabiSalem LahlouYazid JanatiHagit MesserEric Moulines
###### Abstract
Commercial Microwave Links \(CMLs\) offer dense spatial coverage for rainfall sensing but produce path\-integrated measurements that make accurate ground\-level reconstruction challenging\. Existing methods typically oversimplify CMLs as point sensors and neglect line integration relating rainfall to signal attenuation, resulting in degraded performance under heterogeneous precipitation\. In this work, we view rain field reconstruction as a Bayesian inverse problem with Diffusion Models \(DMs\) as high\-fidelity spatial priors\. We show that diffusion models better preserve key rainfall statistics compared to censored Gaussian processes\. Framing rainfall estimation as a Bayesian inverse problem with a DM prior enables training\-free posterior sampling using a broad family of methods, including Plug\-and\-Play, Sequential Monte Carlo, and Replica Exchange methods\. Experiments on synthetic and real\-world datasets demonstrate consistent improvements over established CML\-based reconstruction baselines\.
Machine Learning, ICML
## 1Introduction
Accurate rain fields are central to hydrology, flood warning, and water resource management\. Standard sensing relies on Rain Gauges \(RGs\), weather radars, and satellites\. RGs are accurate but sparse\(Messeret al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib139)\)\. Weather radars provide broad coverage, yet require calibration and are sensitive to clutter and atmospheric effects\(Michelson and Koistinen,[2000](https://arxiv.org/html/2605.05520#bib.bib175); Krajewski and Smith,[2002](https://arxiv.org/html/2605.05520#bib.bib172); Harrisonet al\.,[2000](https://arxiv.org/html/2605.05520#bib.bib167)\)\. Weather satellites extend coverage to data\-scarce regions but often trade resolution for reach\(Funket al\.,[2014](https://arxiv.org/html/2605.05520#bib.bib174)\)\. Hence, no single modality provides dense, accurate, and affordable precipitation monitoring at scale\. Opportunistic remote sensing addresses this gap by reusing existing infrastructure\. In particular,*Commercial Microwave Links*\(CMLs\) in cellular backhaul networks can act as near\-ground rain sensors\(Messeret al\.,[2006](https://arxiv.org/html/2605.05520#bib.bib159); Uijlenhoetet al\.,[2018](https://arxiv.org/html/2605.05520#bib.bib160); Grafet al\.,[2020](https://arxiv.org/html/2605.05520#bib.bib182),[2026](https://arxiv.org/html/2605.05520#bib.bib183)\)\. CMLs form dense networks of transmitters and receivers deployed for telecommunication purposes\. During precipitation, rainfall induces attenuation along each link path\{ℒi\}i=1m\\\{\\mathcal\{L\}\_\{i\}\\\}\_\{i=1\}^\{m\}of the network, which is commonly modeled by the empirical power law\(Leijnseet al\.,[2007](https://arxiv.org/html/2605.05520#bib.bib148); Eshelet al\.,[2021](https://arxiv.org/html/2605.05520#bib.bib138)\)
Yi=a∫𝐬∈ℒiX\(𝐬\)bdℒi\+σiZi,Zi∼𝒩\(0,1\),Y\_\{i\}=a\\int\_\{\\mathbf\{s\}\\in\\mathcal\{L\}\_\{i\}\}X\(\\mathbf\{s\}\)^\{b\}\\,\\mathrm\{d\}\\mathcal\{L\}\_\{i\}\\;\+\\;\\sigma\_\{i\}Z\_\{i\},\\quad Z\_\{i\}\\sim\\mathcal\{N\}\(0,1\),\(1\)whereX\(𝐬\)X\(\\mathbf\{s\}\)is the rain rate \(mm/h\), and\(a,b\)\(a,b\)are constants that depend on the characteristics of the link, whileσiZi\\sigma\_\{i\}Z\_\{i\}is an additive Gaussian noise of standard deviationσi\\sigma\_\{i\}independent of the attenuationYiY\_\{i\}\.
Rain\-field reconstruction from CMLs is an inverse problem: for a region of interestΩ⊂ℝ2\\Omega\\subset\\mathbb\{R\}^\{2\}, and given attenuations\{Yi\}i=1m\\\{Y\_\{i\}\\\}\_\{i=1\}^\{m\}, estimate the gridded field in anH×WH\\times Wdomain\. The problem is ill\-posed because measurements are path\-integrated and nonlinear\. Many reconstruction methods approximate each link as a*Virtual Rain Gauge*\(VRG\) and then interpolate\(Shepard,[1968](https://arxiv.org/html/2605.05520#bib.bib168); Overeemet al\.,[2013](https://arxiv.org/html/2605.05520#bib.bib150); Goldshteinet al\.,[2009a](https://arxiv.org/html/2605.05520#bib.bib147); Eshelet al\.,[2021](https://arxiv.org/html/2605.05520#bib.bib138)\), which can discard valuable spatial information and degrade under heterogeneous rainfall\. More elaborate variants, e\.g\., multiple VRGs per link, still inherit interpolation limitations\(Goldshteinet al\.,[2009a](https://arxiv.org/html/2605.05520#bib.bib147)\)\. The fundamental challenge is that CML measurements provide aggregate information along paths, whereas reconstruction requires resolving fine\-scale spatial structure\. This mismatch is a hallmark of ill\-posed inverse problems: infinitely many rain fields can produce identical observations\. Regularization through prior knowledge is essential, yet classical choices such as isotropic smoothness or Gaussian spatial models fail to capture the complex, multi\-scale, and intermittent structure of real precipitation\. In this work, we instead cast CML\-based rain field reconstruction as Bayesian inference with a learned generative prior\. Recent works show that deep generative models, especially*Diffusion Models*\(DMs\), provide strong priors for inverse problems and are amenable to*training\-free*posterior sampling via inference\-time guidance\(Daraset al\.,[2024](https://arxiv.org/html/2605.05520#bib.bib106); Janatiet al\.,[2025b](https://arxiv.org/html/2605.05520#bib.bib107); Chunget al\.,[2023](https://arxiv.org/html/2605.05520#bib.bib102)\)\. This setting consists in training a DM offline on representative data and then reusing it at inference as a prior to regularize ill\-posed inverse problems\.
Building on this paradigm, we make the following contributions\.1We show that DM priors yield high\-fidelity spatial rain fields and improve key statistics, e\.g\. cumulative precipitation, compared to censored Gaussian Processes\.2We leverage training\-free posterior sampling to enable flexible inference with a fixed pre\-trained DM, covering Plug\-and\-Play, Sequential Monte Carlo, and Replica Exchange approaches\.3We move beyond VRG approximations by explicitly modeling the nonlinear, path\-integrated CML operator governed by the empirical the power\-law\.
We validate the approach on simulated and real datasets, and observe consistent gains over established baselines\. To the best of our knowledge, this is the first work to bring posterior sampling with diffusion priors to opportunistic rain sensing with CMLs\.
## 2Related Work
CMLs were originally proposed as opportunistic rainfall sensors based on the physically grounded relationship between rain\-induced signal attenuation and path\-integrated rainfall intensity\(Messeret al\.,[2006](https://arxiv.org/html/2605.05520#bib.bib159)\)\. A dominant modeling assumption in much of the CML literature is to reduce each path\-integrated CML measurement to a single point observation, typically assigned to the link midpoint, thereby ignoring rainfall variability along the link path\. Under this abstraction, rainfall fields are reconstructed using deterministic point\-based spatial interpolation methods such as Inverse Distance Weighting \(IDW\) and Ordinary Kriging \(OK\)\(Zhanget al\.,[2023](https://arxiv.org/html/2605.05520#bib.bib135); Blettneret al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib137); Messeret al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib139)\)\.IDWestimates values at grid locations as weighted averages of nearby measurements, with weights decaying as a power of the distance, whereasOKexploits second\-order spatial statistics and derives interpolation weights from the spatial correlation structure of the observed data through a variogram model\. While these approaches are simple and computationally efficient, collapsing line\-integral measurements into point estimates inherently discards rainfall variability information contained in the line\-integral nature of CML measurements along the link path\. As a consequence, the reconstructed fields are often overly smooth and fail to reproduce fine\-scale spatial variability and localized extremes\. To mitigate this limitation,Messeret al\.\([2022](https://arxiv.org/html/2605.05520#bib.bib139)\); Eshelet al\.\([2021](https://arxiv.org/html/2605.05520#bib.bib138)\)proposed Goldshtein\-Messer\-Zinevich \(GMZ\) algorithm where multiple VRGs are used to approximate the path\-integrated measurement\. The VRGs are placed along the link, and their intensities are iteratively adjusted to ensure consistency with the original CML observation\. While GMZ better represents sub\-link variability, it ultimately remains within the interpolation framework, as it relies on pseudo\-point measurements rather than directly assimilating the path\-integral observations\.
More advanced stochastic reconstruction methods\(Haeseet al\.,[2017](https://arxiv.org/html/2605.05520#bib.bib136); Hörninget al\.,[2019](https://arxiv.org/html/2605.05520#bib.bib142); Blettneret al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib137)\), such as Random Mixing \(RM\), account for the path\-integrated nature of CML measurements and, unlike interpolation methods, provide a principled treatment of reconstruction uncertainty\. It generates precipitation fields as stochastic linear combinations of unconditional random fields that reproduce the observed spatial dependence structure\. RG measurements are imposed as pointwise constraints, while CML observations enter as nonlinear, path\-integrated constraints\. Repeating this procedure over several steps yields an ensemble of spatial fields, allowing reconstruction uncertainty to be quantified\. However, RM relies on the generation and optimization of large ensembles of unconditional random fields, which is computationally expensive, especially for country\-scale applications\. Besides, it is constrained by the need to specify marginal rainfall distributions and spatial dependence structures\. This makes them rain\-gauge dependent, as the reliable estimation of these statistical properties typically requires a sufficiently dense and gauge network for the marginals\. Consequently, such methods are not well suited for use cases in which only a very limited number of rain gauges are available, or where rainfall information is derived almost exclusively from CMLs\.
In this work, we focus on rain field reconstruction using only CML measurements\. We provide an extended discussion of deterministic and stochastic CML\-based reconstruction methods, including the Random Mixing framework and multi\-sensor fusion approaches, in[AppendixA](https://arxiv.org/html/2605.05520#A1)\.
## 3Rain Field Prior Models
Figure 1:Comparison of rain field statistics between the reference samples and generated samples using the DM prior\. The histograms show densities and are built using50005000samples\.


Figure 2:Examples of generated rain fields using the prior DM compared with reference rain fields\. A total of5,0005\{,\}000samples were generated\. In each subfigure, on the left, samples are randomly drawn from the datasets; on the right, they are drawn from the top90%90\\%wet samples by cumulative rain\.The choice of prior is critical for regularizing the ill\-posed CML inverse problem\. We compare two modeling approaches: \(i\) a censored Gaussian process that explicitly encodes rainfall characteristics such as spatial correlation and non\-negativity, and \(ii\) DM that learns them directly from data\.
### 3\.1Censored Gaussian Processes
Rainfall intensity fields are inherently intermittent: they combine frequent zero\-precipitation events with highly skewed positive intensities, as illustrated in[Figure1](https://arxiv.org/html/2605.05520#S3.F1)\. There exists a rich literature on stochastic precipitation modeling and geostatistical reconstruction\(Wilks and Wilby,[1999](https://arxiv.org/html/2605.05520#bib.bib161); Wheateret al\.,[2000](https://arxiv.org/html/2605.05520#bib.bib162); Yanget al\.,[2005](https://arxiv.org/html/2605.05520#bib.bib163)\)\. A classical and widely used statistical abstraction models this mixed behavior through a latent Gaussian Process that is left censored at zero\(Bardossy and Plate,[1992](https://arxiv.org/html/2605.05520#bib.bib143)\)\. Early frameworks\(Bardossy and Plate,[1992](https://arxiv.org/html/2605.05520#bib.bib143); Ailliotet al\.,[2009](https://arxiv.org/html/2605.05520#bib.bib144)\)represent spatial dependence by a latent Gaussian random field and enforce non\-negativity by censoring, often together with power\-type transformations to better match the heavy right tail of positive rain rates\. Subsequent work refined these ideas by introducing more complex transformations of the censored latent process in order to increase flexibility while retaining analytic and computational convenience\(Baxevani and Lennartsson,[2015](https://arxiv.org/html/2605.05520#bib.bib145); Staufferet al\.,[2017](https://arxiv.org/html/2605.05520#bib.bib146)\)\. Despite their practical appeal, censored Gaussian Process \(GP\) priors remain limited by the structure of Gaussian kernels\. In particular, Gaussian\-based kernels tend to induce comparatively*stiff*priors: they can struggle to reproduce the sharp gradients, multi\-scale variability, and nonstationary spatial organization that commonly arise in rain fields\. In addition, when these models are used as priors within Bayesian inverse problems, inference can become computationally demanding\. The censoring operation introduces inequality constraints and requires manipulating high\-dimensional truncated or censored Gaussian distributions, which in turn involves prohibitive sampling procedures\(Pakman and Paninski,[2014](https://arxiv.org/html/2605.05520#bib.bib164); Botev,[2017](https://arxiv.org/html/2605.05520#bib.bib165)\)\.
We adopt the following censored GP model as a classical baseline prior\. The rain fieldXXis defined on aH×WH\\times Wgrid and is constructed from a latent Gaussian random fieldVV\. Specifically, we assumeV∼𝒢𝒫\(μ,k\)V\\sim\\mathcal\{GP\}\(\\mu,k\)withμ:ℝ2→ℝ\\mu\\mathrel\{\\mathop\{\\ordinarycolon\}\}\\mathbb\{R\}^\{2\}\\rightarrow\\mathbb\{R\}\. To accommodate anisotropy, we use a stationary kernel parameterized by a positive definite matrix𝐐∈ℝ2×2\\mathbf\{Q\}\\in\\mathbb\{R\}^\{2\\times 2\}:
k\(𝐬,𝐬′\)=σ2exp\(−12\(𝐬−𝐬′\)⊤𝐐\(𝐬−𝐬′\)\),k\(\\mathbf\{s\},\\mathbf\{s\}^\{\\prime\}\)=\\sigma^\{2\}\\exp\\left\(\-\\frac\{1\}\{2\}\(\\mathbf\{s\}\-\\mathbf\{s\}^\{\\prime\}\)^\{\\top\}\\mathbf\{Q\}\\ \(\\mathbf\{s\}\-\\mathbf\{s\}^\{\\prime\}\)\\right\),where𝐬,𝐬′∈ℝ2\\mathbf\{s\},\\mathbf\{s\}^\{\\prime\}\\in\\mathbb\{R\}^\{2\}\. The observed rain field is then obtained via a power transformation followed by censoring at zero,
X\(𝐬\)=max\(0,V\(𝐬\)β\),for𝐬∈ℝ2,X\(\\mathbf\{s\}\)=\\max\\big\(0,V\(\\mathbf\{s\}\)^\{\\beta\}\\big\),\\quad\\mathrm\{for\}\\;\\;\\mathbf\{s\}\\in\\mathbb\{R\}^\{2\},which induces a point mass at zero together with a skewed positive component\. The parameterβ\\betacontrols the degree of skewness in the positive intensities, while𝐐\\mathbf\{Q\}andσ\\sigmagovern the spatial correlation structure and marginal scale, respectively\. This model is parameterized by\(μ,𝐐,σ,β\)\(\\mu,\\mathbf\{Q\},\\sigma,\\beta\)\. The Expectation–Maximization algorithm\(Ordoñezet al\.,[2018](https://arxiv.org/html/2605.05520#bib.bib131)\)is used to estimate the parameters of the model\(μ,𝐐,σ,β\)\(\\mu,\\mathbf\{Q\},\\sigma,\\beta\)which treats the latent fieldVVas missing data and accounts for the censoring mechanism during inference; details are provided in[AppendixF](https://arxiv.org/html/2605.05520#A6)\.
### 3\.2Diffusion Models
Letx0∼p0x\_\{0\}\\sim p\_\{0\}denote a sample from the data distribution\. We adopt the*Variance Exploding*\(VE\) formulation of Diffusion Models defined by fixing an increasing noise schedule0=σ0<σ1<⋯<σT0=\\sigma\_\{0\}<\\sigma\_\{1\}<\\cdots<\\sigma\_\{T\}\(Songet al\.,[2021b](https://arxiv.org/html/2605.05520#bib.bib86); Karraset al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib89)\)\. The corruption process is defined as a Markov chain that gradually perturbs the data via Gaussian noise,
qt\+1\|t\(xt\+1∣xt\)=𝒩\(xt,\(σt\+12−σt2\)𝐈\),q\_\{t\+1\|t\}\(x\_\{t\+1\}\\mid x\_\{t\}\)=\\mathcal\{N\}\\big\(x\_\{t\},\(\\sigma\_\{t\+1\}^\{2\}\-\\sigma\_\{t\}^\{2\}\)\\mathbf\{I\}\\big\),\(2\)fort∈⟦0,T−1⟧t\\in\\llbracket 0,T\-1\\rrbracket\. This construction yields the closed\-form marginalqt\|0\(xt∣x0\)=𝒩\(xt;x0,σt2𝐈\)q\_\{t\|0\}\(x\_\{t\}\\mid x\_\{0\}\)=\\mathcal\{N\}\(x\_\{t\};x\_\{0\},\\sigma\_\{t\}^\{2\}\\mathbf\{I\}\), or equivalentlyxt=x0\+σtεx\_\{t\}=x\_\{0\}\+\\sigma\_\{t\}\\varepsilonwithε∼𝒩\(0,𝐈\)\\varepsilon\\sim\\mathcal\{N\}\(0,\\mathbf\{I\}\)\. The resulting joint distribution over the diffusion path factorizes as
p0:T\(x0:T\)=p0\(x0\)∏t=0T−1qt\+1\|t\(xt\+1∣xt\)\.\\textstyle p\_\{0\\mathrel\{\\mathop\{\\ordinarycolon\}\}T\}\(x\_\{0\\mathrel\{\\mathop\{\\ordinarycolon\}\}T\}\)=p\_\{0\}\(x\_\{0\}\)\\prod\_\{t=0\}^\{T\-1\}q\_\{t\+1\|t\}\(x\_\{t\+1\}\\mid x\_\{t\}\)\.The same joint distribution also admits a backward factorization\(Cappéet al\.,[2005](https://arxiv.org/html/2605.05520#bib.bib117)\)
p0:T\(x0:T\)=\[∏t=0T−1pt\|t\+1\(xt∣xt\+1\)\]pT\(xT\),\\textstyle p\_\{0\\mathrel\{\\mathop\{\\ordinarycolon\}\}T\}\(x\_\{0\\mathrel\{\\mathop\{\\ordinarycolon\}\}T\}\)=\\Big\[\\prod\_\{t=0\}^\{T\-1\}p\_\{t\|t\+1\}\(x\_\{t\}\\mid x\_\{t\+1\}\)\\Big\]\\,p\_\{T\}\(x\_\{T\}\),where, for sufficiently largeσT\\sigma\_\{T\}, the terminal distributionpTp\_\{T\}is well approximated by𝒩\(0,σT2𝐈\)\\mathcal\{N\}\(0,\\sigma\_\{T\}^\{2\}\\mathbf\{I\}\)\. The reverse transition kernelspt\|t\+1p\_\{t\|t\+1\}are generally intractable, as they depend on the denoising posteriorp0\|t\(x0∣xt\)∝p0\(x0\)qt\|0\(xt∣x0\)p\_\{0\|t\}\(x\_\{0\}\\mid x\_\{t\}\)\\propto p\_\{0\}\(x\_\{0\}\)\\,q\_\{t\|0\}\(x\_\{t\}\\mid x\_\{0\}\)\. Indeed,pt\|t\+1p\_\{t\|t\+1\}can be expressed using the latter as
pt\|t\+1\(xt\|xt\+1\)=𝔼x0∼p0\|t\(⋅\|xt\+1\)\[qt\|0,t\+1\(xt\|x0,xt\+1\)\],p\_\{t\|t\+1\}\(x\_\{t\}\|x\_\{t\+1\}\)\\\!=\\mathbb\{E\}\_\{x\_\{0\}\\sim p\_\{0\|t\}\(\\cdot\|x\_\{t\+1\}\)\}\\big\[q\_\{t\|0,t\+1\}\(x\_\{t\}\|x\_\{0\},x\_\{t\+1\}\)\\big\],where the conditional transitionqℓ\|0,tq\_\{\\ell\|0,t\}for0<ℓ<t0<\\ell<tadmits a closed\-form expression\(Hoet al\.,[2020](https://arxiv.org/html/2605.05520#bib.bib85); Songet al\.,[2021a](https://arxiv.org/html/2605.05520#bib.bib87)\)\. In the VE setting, it is Gaussian with mean given by a convex combination of\(x0,xt\)\(x\_\{0\},x\_\{t\}\),
qℓ\|0,t\(xℓ∣x0,xt\)=𝒩\(γℓ\|txt\+\(1−γℓ\|t\)x0,σℓ2\(1−γℓ\|t\)𝐈\),q\_\{\\ell\|0,t\}\(x\_\{\\ell\}\\mid x\_\{0\},x\_\{t\}\)\\\!=\\\!\\mathcal\{N\}\\\!\\Big\(\\gamma\_\{\\ell\|t\}x\_\{t\}\+\(1\-\\gamma\_\{\\ell\|t\}\)x\_\{0\},\\;\\sigma\_\{\\ell\}^\{2\}\(1\-\\gamma\_\{\\ell\|t\}\)\\mathbf\{I\}\\\!\\Big\),whereγℓ\|t=σℓ2/σt2\\gamma\_\{\\ell\|t\}=\\sigma\_\{\\ell\}^\{2\}/\\sigma\_\{t\}^\{2\}\. The commonly adopted surrogate for the intractable reverse kernelpt\|t\+1p\_\{t\|t\+1\}replaces the unknown clean samplex0x\_\{0\}in the above transition with the output of a neural network\(t,xt\)↦Dtθ\(xt\)\(t,x\_\{t\}\)\\mapsto D^\{\\theta\}\_\{t\}\(x\_\{t\}\)parameterized byθ\\theta\. This yields the approximate reverse process
pt\|t\+1θ\(xt∣xt\+1\)=qt\|0,t\+1\(xt∣Dt\+1θ\(xt\+1\),xt\+1\)\.p^\{\\theta\}\_\{t\|t\+1\}\(x\_\{t\}\\mid x\_\{t\+1\}\)=q\_\{t\|0,t\+1\}\\big\(x\_\{t\}\\mid D^\{\\theta\}\_\{t\+1\}\(x\_\{t\+1\}\),x\_\{t\+1\}\\big\)\.The neural network is trained to predict the posterior meanDt\(xt\)≔∫x0p0\|t\(x0\|xt\)dx0D\_\{t\}\(x\_\{t\}\)\\coloneqq\\int x\_\{0\}\\,p\_\{0\|t\}\(x\_\{0\}\|x\_\{t\}\)\\mathrm\{d\}x\_\{0\}by minimizing a denoising regression objectives\(Hoet al\.,[2020](https://arxiv.org/html/2605.05520#bib.bib85); Nichol and Dhariwal,[2021](https://arxiv.org/html/2605.05520#bib.bib134); Karraset al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib89)\)\. After training, new samples are generated by initializingxT∼𝒩\(0,σT2𝐈\)x\_\{T\}\\sim\\mathcal\{N\}\(0,\\sigma\_\{T\}^\{2\}\\mathbf\{I\}\)and iteratively drawing samples frompt\|t\+1θ\(⋅\|xt\+1\)p^\{\\theta\}\_\{t\|t\+1\}\(\\cdot\|x\_\{t\+1\}\)down tox0x\_\{0\}\.
### 3\.3Diffusion Model Prior on OpenMRG Data
We train a DM using the EDM framework\(Karraset al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib89)\)on the OpenMRG dataset\(Anderssonet al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib83)\), which provides synchronized CML, radar, and RG measurements over Gothenburg during June–August 2015; refer to[Section5\.2](https://arxiv.org/html/2605.05520#S5.SS2)for details\. Training is performed on21,19621\{,\}196radar\-derived rain\-rate fields, preprocessed withPyNNcml\(Habi,[2026](https://arxiv.org/html/2605.05520#bib.bib170)\)and discretized on a48×3648\\times 36spatial grid\. Additional implementation details are provided in Appendix[B](https://arxiv.org/html/2605.05520#A2)\.
*Evaluation of rain field priors\.*We evaluate the quality of the learned DM prior independently of the inverse problem, focusing on its ability to generate realistic rain fields\.[Figure2](https://arxiv.org/html/2605.05520#S3.F2)presents qualitative samples drawn from the DM alongside reference rain fields, while[Figure1](https://arxiv.org/html/2605.05520#S3.F1)reports histograms of rain rates and cumulative precipitation\. Both statistics show close agreement with the reference data, and the qualitative examples indicate that the trained DM preserves the sparse structure characteristic of rainfall fields\.
To further assess distributional fidelity, we use a classifier two\-sample test\(Lopez\-Paz and Oquab,[2016](https://arxiv.org/html/2605.05520#bib.bib177)\)\. Specifically, it consists of training a classifier to distinguish between a held\-out set of real rain fields from generated samples by the diffusion prior\. Classification accuracy close to50%50\\%indicates that generated samples are difficult to distinguish from real data, and therefore that the learned prior captures the target distribution well\. We use the discriminator architecture ofHabi and Messer \([2021](https://arxiv.org/html/2605.05520#bib.bib180)\)and report classification accuracy on a balanced held\-out test set as the evaluation metric\. The trained discriminator achieves57%57\\%accuracy, where chance level is50%50\\%, indicating only limited distinguishability between generated and real rainfall fields\. We provide more implementation details in[AppendixB](https://arxiv.org/html/2605.05520#A2.SS0.SSS0.Px4)\.
Figure 3:Comparison between the reconstructions of the baselines on the inverse problem with diffusion prior on the setting of GP\. The y\-axis shows intensities while the x\-axis represents a one\-dimensional grid of the\[−5,5\]\[\-5,5\]interval\. The dashed horizontal lines depict values of the integral over the observed intervals\.
## 4Posterior Sampling with Diffusion Models
We now formulate CML\-based rainfall reconstruction as a Bayesian inverse problem and describe training\-free posterior sampling methods that leverage pre\-trained DM prior\.
#### CML Measurement Operator\.
LetX:Ω→ℝ\+X\\mathrel\{\\mathop\{\\ordinarycolon\}\}\\Omega\\to\\mathbb\{R\}\_\{\+\}be a rain field and letx0∈ℝ\+H×Wx\_\{0\}\\in\\mathbb\{R\}\_\{\+\}^\{H\\times W\}denote its discretization on grid points\{𝐬k\}k=1HW\\\{\\mathbf\{s\}\_\{k\}\\\}\_\{k=1\}^\{HW\}, with entries\[x0\]k=X\(𝐬k\)\[x\_\{0\}\]\_\{k\}=X\(\\mathbf\{s\}\_\{k\}\)\. The measurement of the CML at theii\-linkℒi\\mathcal\{L\}\_\{i\}is described by a power\-law parameters\(a,b\)\(a,b\), and a noise levelσi\\sigma\_\{i\}as
Yi=a∫𝐬∈ℒiX\(𝐬\)bdℒi\+σiZi,Zi∼𝒩\(0,1\)\.Y\_\{i\}=a\\int\_\{\\mathbf\{s\}\\in\\mathcal\{L\}\_\{i\}\}X\(\\mathbf\{s\}\)^\{b\}\\,\\mathrm\{d\}\\mathcal\{L\}\_\{i\}\\;\+\\;\\sigma\_\{i\}Z\_\{i\},\\quad Z\_\{i\}\\sim\\mathcal\{N\}\(0,1\)\.\(3\)Under the standard piecewise\-constant assumption on grid cells\(Eshelet al\.,[2021](https://arxiv.org/html/2605.05520#bib.bib138)\), the line integral can be approximated by a weighted sum
Yi=a∑k=1HWΔki\[x0\]kb\+σiZi=⟨aΔi,x0b⟩\+σiZi,Y\_\{i\}=a\\sum\_\{k=1\}^\{HW\}\\Delta^\{i\}\_\{k\}\\,\[x\_\{0\}\]\_\{k\}^\{b\}\+\\sigma\_\{i\}Z\_\{i\}\\ =\\ \\langle a\\Delta^\{i\},x\_\{0\}^\{b\}\\rangle\+\\sigma\_\{i\}Z\_\{i\},\(4\)where*abuse the notation*and usex0bx\_\{0\}^\{b\}to denote the elementwisebb\-th power,Δki\\Delta^\{i\}\_\{k\}is the intersection length between the link pathℒi\\mathcal\{L\}\_\{i\}and cellkk, and⟨⋅,⋅⟩\\langle\\cdot,\\cdot\\rangleis the Frobenius inner product\. Stacking themmlinks gives the observation modely=ℳ\(x0\)\+Σzy=\\mathcal\{M\}\(x\_\{0\}\)\+\\Sigma zwithΣ=diag\(σ1,…,σm\)\\Sigma=\\mathrm\{diag\}\(\\sigma\_\{1\},\\dots,\\sigma\_\{m\}\),z∼𝒩\(0,Im\)z\\sim\\mathcal\{N\}\(0,I\_\{m\}\)andℳ\(x0\)i=⟨aΔi,x0b⟩\\mathcal\{M\}\(x\_\{0\}\)\_\{i\}=\\langle a\\Delta^\{i\},x\_\{0\}^\{b\}\\rangle\. We compute\{Δi\}i=1m\\\{\\Delta^\{i\}\\\}\_\{i=1\}^\{m\}by adapting Siddon’s ray\-tracing algorithm\(Siddon,[1985](https://arxiv.org/html/2605.05520#bib.bib176)\); see Appendix[D](https://arxiv.org/html/2605.05520#A4)for details\. Therefore, the likelihood of the inverse problem is
p\(y∣x0\)∝exp\{−12‖y−ℳ\(x0\)‖Σ−22\},\\displaystyle p\(y\\mid x\_\{0\}\)\\propto\\exp\\\{\-\\tfrac\{1\}\{2\}\\\|y\-\\mathcal\{M\}\(x\_\{0\}\)\\\|\_\{\\Sigma^\{\-2\}\}^\{2\}\\\},\(5\)with∥v∥Σ−22:=v⊤Σ−2v\\\|v\\\|\_\{\\Sigma^\{\-2\}\}^\{2\}\\mathrel\{\\mathop\{\\ordinarycolon\}\}=v^\{\\top\}\\Sigma^\{\-2\}v\. While the operatorℳ\\mathcal\{M\}is nonlinear inx0x\_\{0\}, it is linear in the transformed variableu=x0bu=x\_\{0\}^\{b\}\.
#### Training\-Free Posterior Sampling with DM Priors\.
In training\-free posterior sampling, the goal is, given measurementsyy, sample from the posterior distribution
π0\(x0∣y\)∝p\(y∣x0\)p0\(x0\),\\pi\_\{0\}\(x\_\{0\}\\mid y\)\\propto p\(y\\mid x\_\{0\}\)\\,p\_\{0\}\(x\_\{0\}\),\(6\)wherep0p\_\{0\}is the implicit prior induced by a pre\-trained DM\. Let\{πt\}t=0T\\\{\\pi\_\{t\}\\\}\_\{t=0\}^\{T\}be the marginals obtained by applying the same corruption processqt∣0\(⋅∣x0\)q\_\{t\\mid 0\}\(\\cdot\\mid x\_\{0\}\), as defined in \([2](https://arxiv.org/html/2605.05520#S3.E2)\), toπ0\(⋅∣y\)\\pi\_\{0\}\(\\cdot\\mid y\)\. Then the posterior reverse transitions satisfy
πt∣t\+1\(xt∣xt\+1,y\)∝pt\(y∣xt\)pt∣t\+1\(xt∣xt\+1\),\\displaystyle\\pi\_\{t\\mid t\+1\}\(x\_\{t\}\\mid x\_\{t\+1\},y\)\\propto p\_\{t\}\(y\\mid x\_\{t\}\)\\,p\_\{t\\mid t\+1\}\(x\_\{t\}\\mid x\_\{t\+1\}\),\(7\)wherept\(y∣xt\)=𝔼x0∼p0∣t\(⋅∣xt\)\[p\(y∣x0\)\],\\displaystyle\\text\{where\}\\quad p\_\{t\}\(y\\mid x\_\{t\}\)=\\mathbb\{E\}\_\{x\_\{0\}\\sim p\_\{0\\mid t\}\(\\cdot\\mid x\_\{t\}\)\}\\big\[\\,p\(y\\mid x\_\{0\}\)\\,\\big\],seeJanatiet al\.\([2025b](https://arxiv.org/html/2605.05520#bib.bib107), Eq\. 1\.18\)\. The intermediate likelihoodspt\(y∣xt\)p\_\{t\}\(y\\mid x\_\{t\}\)are typically intractable, as they involve an expectation of the likelihood under the denoising posterior distributionp0∣t\(x0∣xt\)p\_\{0\\mid t\}\(x\_\{0\}\\mid x\_\{t\}\), which is costly to sample\.
A common approximation, introduced and used in the seminal work ofChunget al\.\([2023](https://arxiv.org/html/2605.05520#bib.bib102)\)replaces the full posterior expectation in \([7](https://arxiv.org/html/2605.05520#S4.E7)\) with a point estimate: one plugs the denoiserDt\(xt\)D\_\{t\}\(x\_\{t\}\)into the likelihood;*i\.e\.*
pt\(y∣xt\)≈p0\(y∣Dt\(xt\)\),p\_\{t\}\(y\\mid x\_\{t\}\)\\approx p\_\{0\}\(y\\mid D\_\{t\}\(x\_\{t\}\)\),and uses the resulting gradient to guide sampling, analogous to classifier guidance\(Dhariwal and Nichol,[2021](https://arxiv.org/html/2605.05520#bib.bib158)\)\. This approximation is exact whenp0∣t\(⋅∣xt\)p\_\{0\\mid t\}\(\\cdot\\mid x\_\{t\}\)concentrates at a single point which holds approximately at low noise levels but deteriorates asσt\\sigma\_\{t\}increases\.
Alternatively, several approaches have been proposed to alleviate the problem of intractability using different approximations\.*Score\-based methods*\(Songet al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib133)\)modify the reverse SDE by adding a user\-chosen likelihood score term\.*Projection methods*, such as DDRM\(Kawaret al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib132)\), exploit the special structure in linear Gaussian inverse problems\.*Variational approaches*\(Mardaniet al\.,[2024](https://arxiv.org/html/2605.05520#bib.bib98); Moufadet al\.,[2025](https://arxiv.org/html/2605.05520#bib.bib101); Janatiet al\.,[2025a](https://arxiv.org/html/2605.05520#bib.bib112)\)use different approximation ofpt\(y∣xt\)p\_\{t\}\(y\\mid x\_\{t\}\)together with variational inference to approximate either the posteriorπ0\\pi\_\{0\}in \([6](https://arxiv.org/html/2605.05520#S4.E6)\) \(e\.g\., with a Gaussian\)\(Mardaniet al\.,[2024](https://arxiv.org/html/2605.05520#bib.bib98)\), or the posterior transitions\(Moufadet al\.,[2025](https://arxiv.org/html/2605.05520#bib.bib101); Janatiet al\.,[2025a](https://arxiv.org/html/2605.05520#bib.bib112)\)\.*Sequential Monte Carlo \(SMC\)*approximates the sequence of posterior distributions\{πt\(⋅∣y\)\}t=0T\\\{\\pi\_\{t\}\(\\cdot\\mid y\)\\\}\_\{t=0\}^\{T\}using a population of weighted particles\(Del Moralet al\.,[2006](https://arxiv.org/html/2605.05520#bib.bib105); Chopinet al\.,[2020](https://arxiv.org/html/2605.05520#bib.bib156)\)\. Unlike variational approaches that evolve a single solution, SMC propagates and reweightsNNparticles, thereby maintaining a population\-based representation of the posterior\.*Replica Exchange*constructs a Markov chain targeting the product distribution⨂t=0Tπt\(⋅∣y\)\\bigotimes\_\{t=0\}^\{T\}\\pi\_\{t\}\(\\cdot\\mid y\)\. Each iteration alternates between two types of updates: \(i\)*local moves*, which independently target the marginalsπt\(⋅∣y\)\\pi\_\{t\}\(\\cdot\\mid y\)and can be executed in parallel, and \(ii\)*communication moves*, which propose Metropolis–Hastings swaps between neighboring states\(tm−1,tm\)\(t\_\{m\-1\},t\_\{m\}\)\. InHeet al\.\([2025b](https://arxiv.org/html/2605.05520#bib.bib90)\), the communication step departs from standard swap proposals by instead proposing new states\. Since the resulting acceptance ratio involves the prior marginal density,Heet al\.\([2025b](https://arxiv.org/html/2605.05520#bib.bib90)\)relies on the estimator ofHeet al\.\([2025a](https://arxiv.org/html/2605.05520#bib.bib91)\), which computes density ratios via the forward and the surrogate reverse processes\.
We refer toDaraset al\.\([2024](https://arxiv.org/html/2605.05520#bib.bib106)\); Janatiet al\.\([2025b](https://arxiv.org/html/2605.05520#bib.bib107)\)for detailed survey on the landscape of methods for training\-free posterior sampling with DM priors\.













Figure 4:Comparisons of rain field reconstructions on real CMLs links from OpenMRG\. We depict the network of CMLs in red\.
## 5Experiments
#### Experimental settings\.
We conduct \(i\) simulated and controlled experiments based on GP prior\(Rasmussen and Williams,[2006](https://arxiv.org/html/2605.05520#bib.bib113)\), and \(ii\) real\-world experiments on the OpenMRG dataset\(Anderssonet al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib83); Fenclet al\.,[2024](https://arxiv.org/html/2605.05520#bib.bib169)\)\. In our target setting, CMLs are the primary source of observations\. Hence, we consider as classical baselinesIDW,GMZ, andOK, whereasRMis not applicable because it requires rain\-gauge measurements to estimate the marginal distributions used by the method\.
#### Meteorological baselines\.
We briefly revisit the principles behind the baselines used in rain field reconstruction\.
*Inverse distance weighting \(IDW\)\.*In IDW, the path\-integrated observations of CMLs are converted into VRGs at the middle of the links\. The rain rate at location𝐬\\mathbf\{s\}is estimated as a weighted average of nearby observations
X^\(𝐬\)=∑i=1Nwi\(𝐬\)X\(𝐬i\),wi\(𝐬\)=di−p/∑j=1Ndj−p,\\textstyle\\hat\{X\}\(\\mathbf\{s\}\)=\\sum\_\{i=1\}^\{N\}w\_\{i\}\(\\mathbf\{s\}\)\\,X\(\\mathbf\{s\}\_\{i\}\),\\quad w\_\{i\}\(\\mathbf\{s\}\)=d\_\{i\}^\{\-p\}/\\sum\_\{j=1\}^\{N\}d\_\{j\}^\{\-p\},whereX\(𝐬i\)X\(\\mathbf\{s\}\_\{i\}\)is the rainfall intensity associated with theii\-th VRG positioned at location𝐬i\\mathbf\{s\}\_\{i\},di=‖𝐬−𝐬i‖d\_\{i\}=\\\|\\mathbf\{s\}\-\\mathbf\{s\}\_\{i\}\\\|is the distance to𝐬i\\mathbf\{s\}\_\{i\}, andp\>0p\>0controls the decay of influence with distance\.
*GMZalgorithm\.*Rather than representing each CML by a single VRGs at its midpoint, the GMZ algorithm accounts for the variability of the rainfall rate along the CML link consideringKKVRGs evenly distributed along the CML link\. More precisely, for the CML linkℒi\\mathcal\{L\}\_\{i\}, theKKVRGs are positioned at\{𝐬i,k\}k=1K\\\{\\mathbf\{s\}\_\{i,k\}\\\}\_\{k=1\}^\{K\}, and the rain intensity at the latter points is constrained to match the observationYi=aK∑k=1KXℓ,kbY\_\{i\}=\\frac\{a\}\{K\}\\sum\_\{k=1\}^\{K\}X\_\{\\ell,k\}^\{b\}\. This redistribution is performed iteratively using neighboring links to enforce spatial consistency, and the resulting VRGs are interpolated with IDW\(Goldshteinet al\.,[2009b](https://arxiv.org/html/2605.05520#bib.bib154)\)\.
*Ordinary Kriging \(OK\)\.*Similar to IDW, OK uses the VRG simplification of the path\-integrated measurement\. OK uses GP prior with user\-defined kernel \(variogram\) as an assumption over the rain field to interpolate the observations\. The prior parameters are chosen by minimizing the L1 norm\(Matheron,[1963a](https://arxiv.org/html/2605.05520#bib.bib152); Cressie,[1993a](https://arxiv.org/html/2605.05520#bib.bib153)\)\.
#### Diffusion baselines\.
We select seven training\-free posterior sampling baselines that represent the approach presented in[Section4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2), namely,DPS\(Chunget al\.,[2023](https://arxiv.org/html/2605.05520#bib.bib102)\),DAPS\(Zhanget al\.,[2024](https://arxiv.org/html/2605.05520#bib.bib99)\),RedDiff\(Mardaniet al\.,[2024](https://arxiv.org/html/2605.05520#bib.bib98)\),CREPE\(Heet al\.,[2025b](https://arxiv.org/html/2605.05520#bib.bib90)\),MGPS\(Moufadet al\.,[2025](https://arxiv.org/html/2605.05520#bib.bib101)\),MGDM\(Janatiet al\.,[2025a](https://arxiv.org/html/2605.05520#bib.bib112)\),TDS\(Wuet al\.,[2023](https://arxiv.org/html/2605.05520#bib.bib100)\)\.
We provide details on the implementation of both meteorological and Diffusion baselines as well as hyperparameter choices in[AppendixC](https://arxiv.org/html/2605.05520#A3)\. When reporting results, we highlight the best metric withtxand second besttx\.
WhileCREPEis included in the Gaussian process experiments, we exclude it from the CML experiments\. In this setting, the method exhibits substantially higher computational cost than the other baselines and yields moderate performance\. For these reasons, we omitCREPEfrom both the quantitative results and qualitative comparisons\.
### 5\.1Simulated examples
We consider a one\-dimensional linear inverse problem where the latent field is a Gaussian processX∼𝒢𝒫\(0,k\)X\\sim\\mathcal\{GP\}\(0,k\)on\[−5,5\]\[\-5,5\], with an RBF kernel
k\(s,s′\)=exp\(−\(s−s′\)22ℓ2\),ℓ=0\.6\.k\(s,s^\{\\prime\}\)=\\exp\\\!\\left\(\-\\frac\{\(s\-s^\{\\prime\}\)^\{2\}\}\{2\\ell^\{2\}\}\\right\),\\qquad\\ell=0\.6\.\(8\)The links\{ℒi\}i=1m\\\{\\mathcal\{L\}\_\{i\}\\\}\_\{i=1\}^\{m\}in this setting reduce to intervals\{\[ai,bi\]\}i=1m\\\{\[a\_\{i\},b\_\{i\}\]\\\}\_\{i=1\}^\{m\}111In the GP experiments, we abuse the notation and useaia\_\{i\}andbib\_\{i\}to denote integration intervals instead of the power\-law constants\.,and the observation model
Yi=∫aibiX\(s\)ds\+σZi,Zi∼𝒩\(0,1\),Y\_\{i\}=\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}X\(s\)\\,\\mathrm\{d\}s\+\\sigma Z\_\{i\},\\qquad Z\_\{i\}\\sim\\mathcal\{N\}\(0,1\),\(9\)where\{Zi\}i=1m\\\{Z\_\{i\}\\\}\_\{i=1\}^\{m\}are mutually independent and independent ofXX\. Because the observation operator is linear and the prior is Gaussian, the posterior is a GP with closed\-form mean and covariance\(Rasmussen and Williams,[2006](https://arxiv.org/html/2605.05520#bib.bib113)\):
μX\|y\(s\)\\displaystyle\\mu\_\{X\|y\}\(s\)=𝐤y\(s\)⊤\(𝐊yy\+σ2𝐈\)−1𝐲,\\displaystyle=\\mathbf\{k\}\_\{y\}\(s\)^\{\\top\}\(\\mathbf\{K\}\_\{yy\}\+\\sigma^\{2\}\\mathbf\{I\}\)^\{\-1\}\\mathbf\{y\},kX\|y\(s,s′\)\\displaystyle k\_\{X\|y\}\(s,s^\{\\prime\}\)=k\(s,s′\)−𝐤y\(s\)⊤\(𝐊yy\+σ2𝐈\)−1𝐤y\(s′\),\\displaystyle=k\(s,s^\{\\prime\}\)\-\\mathbf\{k\}\_\{y\}\(s\)^\{\\top\}\(\\mathbf\{K\}\_\{yy\}\+\\sigma^\{2\}\\mathbf\{I\}\)^\{\-1\}\\mathbf\{k\}\_\{y\}\(s^\{\\prime\}\),where𝐤y\(s\)\\mathbf\{k\}\_\{y\}\(s\)is a vector inℝm\\mathbb\{R\}^\{m\}with coordinates\[𝐤y\(s\)\]i=∫aibik\(s,s′\)ds′\[\\mathbf\{k\}\_\{y\}\(s\)\]\_\{i\}=\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}k\(s,s^\{\\prime\}\)\\,\\mathrm\{d\}s^\{\\prime\}, and𝐊yy\\mathbf\{K\}\_\{yy\}is a matrix inℝm×m\\mathbb\{R\}^\{m\\times m\}with elements\[𝐊yy\]ij=∫aibi∫ajbjk\(s,s′\)dsds′\[\\mathbf\{K\}\_\{yy\}\]\_\{ij\}=\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}\\int\_\{a\_\{j\}\}^\{b\_\{j\}\}k\(s,s^\{\\prime\}\)\\,\\mathrm\{d\}s\\,\\mathrm\{d\}s^\{\\prime\}\. In our case, we postpone the derivations of these formulas to[AppendixE](https://arxiv.org/html/2605.05520#A5)\. We discretize the interval\[−5,5\]\[\-5,5\]into5050points\. Evaluating the Gaussian processXXon this grid yields a multivariate Gaussian distribution\. Consequently, for training\-free baselines, we employed a DM prior to targeting this Gaussian distribution\. Such a DM admits a closed\-form construction, as described inMoufadet al\.\([2025](https://arxiv.org/html/2605.05520#bib.bib101), Appendix B\)\.
#### Evaluation\.
This setting provides an oracle posterior, enabling a precise evaluation\. We report \(i\) the sliced Wasserstein distance between the true and approximate posteriors\(Rabinet al\.,[2011](https://arxiv.org/html/2605.05520#bib.bib155); Bonneelet al\.,[2015](https://arxiv.org/html/2605.05520#bib.bib157)\), and \(ii\)ℓ2\\ell\_\{2\}errors between the oracle and the approximate quantities: posterior mean and the5%5\\%\-95%95\\%pointwise quantiles\.
#### Results\.
The results are reported in[Table1](https://arxiv.org/html/2605.05520#S5.T1)\. On the GP benchmark, TDS is the strongest overall method, achieving the best distributional match \(lowest SW\) and the most accurate uncertainty while remaining competitive on the mean\. In contrast,RedDiffachieves the lowest mean error but performs poorly on SW and tails, indicating that mean accuracy alone can mask severe posterior miscalibration\.DAPSandMGPSperform competitively with moderate quantile errors\.MGDMis generally less competitive, with a higher SW and substantially worse tail errors\. We provide qualitative examples in[Figures3](https://arxiv.org/html/2605.05520#S3.F3)and[7](https://arxiv.org/html/2605.05520#A5.F7)\.
Table 1:Quantitative comparison between Diffusion Models baselines on GP experiment\. First row reports Sliced Wasserstein distance\. The last 3 rows reportℓ2\\ell\_\{2\}distance between reference and approximate mentioned quantity\. Lower metrics are better\.MGDMDAPSRedDiffCREPEDPSTDSMGPSSW↓\\downarrow0\.160\.090\.380\.130\.120\.070\.09Mean↓\\downarrow0\.610\.540\.360\.680\.980\.530\.535%5\\%\-Quantile↓\\downarrow2\.211\.125\.701\.821\.170\.850\.9995%95\\%\-Quantile↓\\downarrow2\.201\.295\.731\.511\.301\.021\.35
### 5\.2Real dataset
We evaluate our approach on the OpenMRG dataset\(Anderssonet al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib83)\)and use the Python pacakgePyNNcml\(Habi,[2026](https://arxiv.org/html/2605.05520#bib.bib170); Chwalaet al\.,[2026](https://arxiv.org/html/2605.05520#bib.bib181)\)to download and preprocess it\. The dataset provides synchronized measurements from CMLs, weather radar, and RGs over the Gothenburg metropolitan area \(Sweden\) during June–August 2015\.
#### Commercial microwave links\.
The network contains 364 bidirectional links \(728 sub\-links\) operating at frequencies between 7 and 38 GHz \(mostly above 25 GHz\), with path lengths between 0\.1 km and 15 km \(median≈2\\approx 2km\)\. Transmitted and received signal levels are recorded every 10 seconds for each sub\-link\. Each link is associated with a pair of power\-law parameters\(a,b\)\(a,b\)\(Eq\.[1](https://arxiv.org/html/2605.05520#S1.E1)\), whose dependence on frequency and path length is summarized in[Figure5](https://arxiv.org/html/2605.05520#A2.F5)\.
#### Radar\-derived reference fields\.
Radar data are drawn from the Swedish NORDRAD composite product operated by SMHI\. The composite is generated from multiple Swedish radars using the Rainbow software and provides reflectivity fields on a 2 km grid; coverage over Gothenburg is primarily from the Vara radar \(at≈\\approx78 km\), with partial contributions from nearby radars\. The product is gauge\-corrected using a distance\-dependent gauge\-radar ratio estimation\(Michelson and Koistinen,[2000](https://arxiv.org/html/2605.05520#bib.bib175)\)\. Reflectivity is expressed in dBZ and stored as pseudo\-dBZ values, which are converted back to actual dBZ before processing; values are then mapped to linear units and converted to rainfall intensity via a standardZZ\-RRpower\-law relationship commonly used for stratiform conditions\.
#### Evaluation\.
Since spatially dense ground truth is not directly observed, we use the gauge\-adjusted radar\-derived rain\-rate fields as reference maps for evaluation\. Attenuation observations are generated by applying the CML observation operator to radar rain\-rate fields from the test split, using the OpenMRG network topology, link lengths, and link\-specific power\-law parameters\(ai,bi\)\(a\_\{i\},b\_\{i\}\)\. We consider two noise settings:1*Isotropic noise*where all CML observations are corrupted with the same noise levelσ\\sigmaand2*Heteroscedastic noise*, where the noise level depends on the CML length\. Specifically, for linkiiwith lengthLiL\_\{i\}, we setσi=σ2\(1\+LiLmax\),\\sigma\_\{i\}=\\frac\{\\sigma\}\{2\}\\left\(1\+\\frac\{L\_\{i\}\}\{L\_\{\\max\}\}\\right\),so that longer links have larger noise levels, withσi=σ\\sigma\_\{i\}=\\sigmawhenLi=LmaxL\_\{i\}=L\_\{\\max\}andσ/2≤σi<σ\\sigma/2\\leq\\sigma\_\{i\}<\\sigmaotherwise\. Note that other heteroscedastic noise models could also be considered\. Such settings can be incorporated straightforwardly by specifying the covariance matrixΣ\\Sigmaof the observation operator in[Equation5](https://arxiv.org/html/2605.05520#S4.E5)\.
#### Metrics\.
Existing studies commonly rely on standard pixel\-wise reconstruction metrics\(Grafet al\.,[2021](https://arxiv.org/html/2605.05520#bib.bib140); Zhanget al\.,[2023](https://arxiv.org/html/2605.05520#bib.bib135)\)\. In our assessment, we report the root mean square error \(RMSE\), which quantifies the discrepancy between ground\-truth and reconstructed fields, and the Pearson correlation coefficient \(PCC\), which measures their linear agreement\. We also evaluate the error in cumulative rainfall, defined as the difference between the total precipitation of the reconstructed and ground\-truth fields\. Values closer to zero indicate better agreement, with positive values corresponding to overestimation and negative values to underestimation of total rainfall\. For computational reasons, we restrict the maximum runtime of the diffusion baselines to33minutes\. ForTDS, this requires reducing the number of particles to satisfy the runtime constraint\. Hyperparameter details for all baselines are provided in[Table7](https://arxiv.org/html/2605.05520#A3.T7)\.
#### Results\.
[Table3](https://arxiv.org/html/2605.05520#S5.T3)shows that the diffusion\-model baselines consistently outperform the meteorological baselines across both settings and evaluation metrics\. In particular,MGPSachieves strong performance on all three metrics\. Overall, DM\-based reconstructions exhibit higher agreement with the reference fields, as reflected by their PCC values\. Among the meteorological methods,GMZperforms well in terms of cumulative rainfall, whileOKyields the strongest overall performance within this class of baselines\. The results under heterogeneous noise are slightly better, as expected, since this setting has a lower overall noise level\.
Qualitative comparisons are provided in[Figure4](https://arxiv.org/html/2605.05520#S4.F4), with a detailed comparison including all methods in[Figures8](https://arxiv.org/html/2605.05520#A7.F8)and[9](https://arxiv.org/html/2605.05520#A7.F9)\.
#### Runtime\.
To assess computational cost, we report the runtimes of the meteorological and diffusion\-based baselines in[Table2](https://arxiv.org/html/2605.05520#S5.T2)\. The results highlight a trade\-off: meteorological methods are fast but they \(expectOK\) produce a deterministic reconstruction; the diffusion baselines are more expensive yet they generate ensembles of plausible rain fields\. This enables the assessment of multiple scenarios, including best\- and worst\-case precipitation patterns that may be relevant for flood\-warning applications\.
Table 2:Runtime \(in seconds\) per reconstruction for all baselines\. For diffusion\-based methods, the reported runtime is divided by1010, since a batch of1010samples is generated jointly\.DAPSDPSMGDMMGPSRedDiffTDSIDWGMZOKRuntime6\.065\.276\.875\.225\.9518\.490\.011\.010\.27
In our experiments, generating1010samples with the diffusion baselines takes roughly one minute\. Since CML data is aggregate at intervals of55to1515minutes, these runtimes remain compatible with real\-time meteorological monitoring\.
Table 3:Comparison between DM and meteorological baselines on RMSE, PCC, and difference of Cumulative Rain on real CML links from OpenMRG\. 95%\-confidence intervals are shown\.RMSE↓\\downarrowPCC↑\\uparrowCumulative Rain→0\\rightarrow 0Isotropic noiseDAPS0\.88±\\pm0\.090\.24±\\pm0\.02\-320±\\pm53\.66DPS0\.88±\\pm0\.090\.20±\\pm0\.01\-333±\\pm54\.34MGDM0\.91±\\pm0\.080\.22±\\pm0\.02231±\\pm35\.29MGPS0\.86±\\pm0\.080\.24±\\pm0\.0142±\\pm34\.65RedDiff0\.86±\\pm0\.090\.23±\\pm0\.02\-356±\\pm48\.47TDS0\.84±\\pm0\.090\.21±\\pm0\.02\-216±\\pm35\.89IDW0\.90±\\pm0\.090\.21±\\pm0\.01\-144±\\pm32\.26GMZ0\.92±\\pm0\.090\.21±\\pm0\.01\-65±\\pm29\.26OK0\.87±\\pm0\.090\.23±\\pm0\.02158±\\pm38\.67Heteroscedastic noiseDAPS0\.87±\\pm0\.100\.24±\\pm0\.01\-438\.80±\\pm76\.08DPS0\.85±\\pm0\.090\.21±\\pm0\.01\-349\.25±\\pm53\.49MGDM0\.92±\\pm0\.080\.24±\\pm0\.01300\.00±\\pm36\.24MGPS0\.83±\\pm0\.080\.26±\\pm0\.0170\.94±\\pm36\.40RedDiff0\.84±\\pm0\.090\.24±\\pm0\.02\-370\.92±\\pm48\.78TDS0\.81±\\pm0\.090\.25±\\pm0\.02\-241\.43±\\pm35\.03IDW0\.86±\\pm0\.090\.24±\\pm0\.01\-166\.95±\\pm32\.09GMZ0\.91±\\pm0\.090\.21±\\pm0\.01\-81\.93±\\pm29\.09OK0\.86±\\pm0\.090\.24±\\pm0\.02140\.74±\\pm38\.83
### 5\.3Ablation
#### DM Prior and Nonlinear Observation Model\.
We first disentangle the gains due to the diffusion prior from those due to accurate modeling of the nonlinear path\-integrated observation operator\. To this end, we consider an ablation in which observations are generated under the true nonlinear setting, with power lawbi\>1b\_\{i\}\>1, while inference assumesbi=1b\_\{i\}=1, yielding a linear observation model\. To assess the role of the diffusion prior quality, we repeat this experiment with both a*“Strong”*and a*“Weak”*prior; the latter corresponding to an undertrained diffusion model trained for a small number of epochs\. Results are reported in[Table5](https://arxiv.org/html/2605.05520#S5.T5)\.
*Discussion\.*The results show that the diffusion prior improves the reconstructed rain fields; with accurate modeling of the nonlinear observation operator being essential for strong quantitative performance\.
Table 4:Ablation between DMs and meteorological baselines on RMSE, PCC, and difference of Cumulative Rain \(Cum Rain\)\. Two settings are considered on OpenMRG dataset with simulated CMLs: few\-long and many\-short links, with the power\-lawb=1b=1\. For RMSE lower is better\. For PCC, higher is better\. For Cum Rain, closer to 0 is better\.Few\-long linksMany\-short linksRMSE↓\\downarrowPCC↑\\uparrowCum RainRMSE↓\\downarrowPCC↑\\uparrowCum RainDAPS0\.920\.24\-2431\.160\.29\-348DPS0\.880\.23\-2670\.800\.32\-152MGDM0\.890\.251520\.790\.38153MGPS0\.840\.28110\.750\.4120RedDiff0\.850\.26\-2210\.790\.37\-224TDS0\.850\.23\-1630\.760\.34\-77IDW0\.890\.25\-1640\.750\.423GMZ0\.860\.29\-200\.750\.4315OK0\.830\.28190\.760\.3827
#### CML Network Configuration\.
We next assess the impact of CML network geometry, focusing on the number and spatial configuration of links\. We consider two synthetic configurations: a*few\-long*setting with2525long CMLs and a*many\-short*setting with100100short CMLs\. To isolate the effect of spatial configuration, we restrict this ablation to the case where the power\-law exponent is equal to11, resulting in a linear inverse problem\. Results are reported in[Table4](https://arxiv.org/html/2605.05520#S5.T4)\.
*Discussion\.*In this idealized linear setting, meteorological baselines perform comparably to DM\-based methods, withOKachieving the lowest RMSE\. This suggests that, under linear sensing and sufficiently dense coverage, the rainfall field can be well approximated by a Gaussian\-process model and effectively reconstructed by interpolation\. Meteorological baselines further benefit from increased CML density, particularly when links are well distributed across the domain, as shown in[Table4](https://arxiv.org/html/2605.05520#S5.T4)\.
Real\-world CML observations, however, deviate substantially from this regime\. As shown in[Figure5](https://arxiv.org/html/2605.05520#A2.F5), the power\-law exponentbbvaries over a wide range, inducing nonlinear measurements\. Moreover, operational CML networks are typically spatially heterogeneous, with dense coverage in some regions and sparse or absent coverage in others, as illustrated in[Figure4](https://arxiv.org/html/2605.05520#S4.F4)\. These factors highlight the limitations of interpolation\-based meteorological baselines and motivate the use of expressive priors coupled with physically grounded observation models, as in the proposed DM\-based framework\.
Table 5:Ablation to discern the improvement that comes from the DM prior vs\. the improvement from accounting for the nonlinearity of the measurement operator\.*“Weak\-prior”*is an under\-trained DM, i\.e\. a DM trained with very\-few number epochs\.Linear operatorNonlinear operatorDMSamplerRMSE↓\\downarrowPCC↑\\uparrowCum Rain→0\\rightarrow 0RMSE↓\\downarrowPCC↑\\uparrowCum Rain→0\\rightarrow 0WeakMGPS0\.940\.17\-430\.300\.880\.23132\.21RedDiff1\.010\.16\-661\.380\.900\.19\-456\.52TDS0\.960\.11\-543\.510\.860\.18\-259\.54StrongMGPS0\.960\.16\-436\.110\.860\.2442\.86RedDiff1\.000\.17\-643\.390\.860\.23\-356\.76TDS0\.970\.13\-549\.760\.840\.21\-216\.44
## 6Conclusion
We presented a Bayesian framework for rain field reconstruction from CMLs that leverages DMs as expressive spatial priors\. Formulating the problem as a Bayesian inverse problem enables training\-free posterior sampling using a wide range of modern inference methods\. Crucially, the proposed framework incorporates a physically grounded measurement operator that accounts for the path\-integrated nature of CML observations, moving beyond the VRG assumption\.
#### Limitations and Future Work\.
The DM prior is trained on data from a single geographic region and season\. Its performance under distribution shift to other climates and meteorological regimes has not been assessed\. Moreover, our current formulation treats rainfall fields independently across time, thereby neglecting temporal structure\. Future work will explore spatiotemporal modeling extensions to address these limitations\.
## Impact Statement
The proposed approach may contribute to the monitoring of hydric resources by enabling more effective use of opportunistic sensing data, and may support urban flood and hazard prevention through improved rainfall estimation\. At the same time, because opportunistic sensing relies on telecommunications infrastructure, its deployment may introduce user exposure and privacy\-related risks that warrant careful consideration\. These potential impacts should be taken into account when translating such methods into operational systems\.
## Acknowledgements
The work of Badr Moufad has been supported by Technology Innovation Institute \(TII\), project Fed2Learn\. The work of Eric Moulines has been partly funded by the European Union \(ERC\-2022\-SYG\-OCEAN\-101071601\)\. Views and opinions expressed are however those of the author\(s\) only and do not necessarily reflect those of the European Union or the European Research Council Executive Agency\. Neither the European Union nor the granting authority can be held responsible for them\. This work was granted access to the HPC resources of IDRIS under the allocations 2025\-AD011015980 made by GENCI\.
## References
- P\. Ailliot, C\. Thompson, and P\. Thomson \(2009\)Space–time modelling of precipitation by using a hidden markov model and censored gaussian distributions\.Journal of the Royal Statistical Society Series C: Applied Statistics58\(3\),pp\. 405–426\.Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p1.1)\.
- J\. C\. Andersson, J\. Olsson, R\. Van de Beek, and J\. Hansryd \(2022\)OpenMRG: open data from microwave links, radar, and gauges for rainfall quantification in gothenburg, sweden\.Earth System Science Data14\(12\),pp\. 5411–5426\.Cited by:[§3\.3](https://arxiv.org/html/2605.05520#S3.SS3.p1.2),[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px1.p1.1),[§5\.2](https://arxiv.org/html/2605.05520#S5.SS2.p1.1)\.
- A\. Bardossy and E\. J\. Plate \(1992\)Space\-time model for daily rainfall using atmospheric circulation patterns\.Water resources research28\(5\),pp\. 1247–1259\.Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p1.1)\.
- A\. Baxevani and J\. Lennartsson \(2015\)A spatiotemporal precipitation generator based on a censored latent g aussian field\.Water Resources Research51\(6\),pp\. 4338–4358\.Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p1.1)\.
- N\. Blettner, C\. Chwala, B\. Haese, S\. Hörning, and H\. Kunstmann \(2022\)Combining commercial microwave link and rain gauge observations to estimate countrywide precipitation: a stochastic reconstruction and pattern analysis approach\.Water Resources Research58\(10\),pp\. e2022WR032563\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p5.1),[§2](https://arxiv.org/html/2605.05520#S2.p1.1),[§2](https://arxiv.org/html/2605.05520#S2.p2.1)\.
- N\. Bonneel, J\. Rabin, G\. Peyré, and H\. Pfister \(2015\)Sliced and radon wasserstein barycenters of measures\.Journal of Mathematical Imaging and Vision\.Cited by:[§5\.1](https://arxiv.org/html/2605.05520#S5.SS1.SSS0.Px1.p1.3)\.
- Z\. I\. Botev \(2017\)The normal law under linear restrictions: simulation and estimation via minimax tilting\.Journal of the Royal Statistical Society Series B: Statistical Methodology79\(1\),pp\. 125–148\.Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p1.1)\.
- O\. Cappé, E\. Moulines, and T\. Rydén \(2005\)Inference in hidden markov models\.Springer\.Cited by:[§3\.2](https://arxiv.org/html/2605.05520#S3.SS2.p1.25)\.
- N\. Chopin, O\. Papaspiliopoulos,et al\.\(2020\)An introduction to sequential monte carlo\.Vol\.4,Springer\.Cited by:[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p3.7)\.
- H\. Chung, J\. Kim, M\. T\. Mccann, M\. L\. Klasky, and J\. C\. Ye \(2023\)Diffusion posterior sampling for general noisy inverse problems\.InThe Eleventh International Conference on Learning Representations,External Links:[Link](https://openreview.net/forum?id=OnD9zGAGT0k)Cited by:[§C\.2](https://arxiv.org/html/2605.05520#A3.SS2.SSS0.Px5.p1.1),[§C\.2](https://arxiv.org/html/2605.05520#A3.SS2.SSS0.Px6.p1.3),[§1](https://arxiv.org/html/2605.05520#S1.p2.3),[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p2.1),[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px3.p1.1)\.
- C\. Chwala, A\. Overeem, E\. Øydvin, L\. Petersson Wårdh, J\. Seidel, M\. Graf, B\. Walraven, E\. Covi, H\. V\. Habi, M\. Fencl,et al\.\(2026\)Open\-source tools for processing opportunistic rainfall sensor data: an overview of existing tools and the new opensense software packages poligrain, pypwsqc and mergeplg\.EGUsphere2026,pp\. 1–31\.Cited by:[§5\.2](https://arxiv.org/html/2605.05520#S5.SS2.p1.1)\.
- N\. A\. C\. Cressie \(1993a\)Statistics for spatial data\.Wiley\.Note:Revised editionCited by:[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px2.p4.1)\.
- N\. A\. C\. Cressie \(1993b\)Statistics for spatial data\.Wiley\.Cited by:[§C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px3.p1.1)\.
- G\. Daras, H\. Chung, C\. Lai, Y\. Mitsufuji, J\. C\. Ye, P\. Milanfar, A\. G\. Dimakis, and M\. Delbracio \(2024\)A survey on diffusion models for inverse problems\.arXiv preprint arXiv:2410\.00083\.Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p2.3),[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p4.1)\.
- P\. Del Moral, A\. Doucet, and A\. Jasra \(2006\)Sequential monte carlo samplers\.Journal of the Royal Statistical Society: Series B \(Statistical Methodology\)68\(3\),pp\. 411–436\.External Links:[Document](https://dx.doi.org/https%3A//doi.org/10.1111/j.1467-9868.2006.00553.x),[Link](https://rss.onlinelibrary.wiley.com/doi/abs/10.1111/j.1467-9868.2006.00553.x),https://rss\.onlinelibrary\.wiley\.com/doi/pdf/10\.1111/j\.1467\-9868\.2006\.00553\.xCited by:[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p3.7)\.
- P\. Dhariwal and A\. Nichol \(2021\)Diffusion models beat gans on image synthesis\.Advances in neural information processing systems34,pp\. 8780–8794\.Cited by:[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p2.3)\.
- A\. Eshel, H\. Messer, H\. Kunstmann, P\. Alpert, and C\. Chwala \(2021\)Quantitative analysis of the performance of spatial interpolation methods for rainfall estimation using commercial microwave links\.Journal of Hydrometeorology22\(4\),pp\. 831–843\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p2.1),[§C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px1.p1.5),[§C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px2.p1.3),[§1](https://arxiv.org/html/2605.05520#S1.p1.1),[§1](https://arxiv.org/html/2605.05520#S1.p2.3),[§2](https://arxiv.org/html/2605.05520#S2.p1.1),[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px1.p1.25)\.
- M\. Fencl, R\. Nebuloni, J\. C\. Andersson, V\. Bares, N\. Blettner, G\. Cazzaniga, C\. Chwala, M\. Colli, L\. de Vos, A\. El Hachem,et al\.\(2024\)Data formats and standards for opportunistic rainfall sensors\.Open Research Europe3,pp\. 169\.Cited by:[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px1.p1.1)\.
- C\. C\. Funk, P\. J\. Peterson, M\. F\. Landsfeld, D\. H\. Pedreros, J\. P\. Verdin, J\. D\. Rowland, B\. E\. Romero, G\. J\. Husak, J\. C\. Michaelsen, and A\. P\. Verdin \(2014\)A quasi\-global precipitation time series for drought monitoring\.ReportU\.S\. Geological Survey\(en\)\.External Links:[Document](https://dx.doi.org/10.3133/ds832)Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p1.1)\.
- O\. Goldshtein, H\. Messer, and A\. Zinevich \(2009a\)Rain rate estimation using measurements from commercial telecommunications links\.IEEE Transactions on signal processing57\(4\),pp\. 1616–1625\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p1.1),[§C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px2.p1.3),[§1](https://arxiv.org/html/2605.05520#S1.p2.3)\.
- O\. Goldshtein, H\. Messer, and A\. Zinevich \(2009b\)Rain rate estimation using measurements from commercial telecommunications links\.IEEE Transactions on Signal Processing57\(4\),pp\. 1616–1625\.External Links:[Document](https://dx.doi.org/10.1109/TSP.2009.2012554)Cited by:[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px2.p3.5)\.
- M\. Graf, V\. Bareš, H\. Messer, R\. Nebuloni, M\. Fencl, C\. Chwala, A\. Overeem, R\. Van de Beek, J\. Olsson, J\. Ostrometzky,et al\.\(2026\)The opportunistic precipitation sensing network \(opensense\)\.Bulletin of the American Meteorological Society107\(3\),pp\. E585–E591\.Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p1.1)\.
- M\. Graf, C\. Chwala, J\. Polz, and H\. Kunstmann \(2020\)Rainfall estimation from a german\-wide commercial microwave link network: optimized processing and validation for 1 year of data\.Hydrology and Earth System Sciences24\(6\),pp\. 2931–2950\.Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p1.1)\.
- M\. Graf, A\. El Hachem, M\. Eisele, J\. Seidel, C\. Chwala, H\. Kunstmann, and A\. Bárdossy \(2021\)Rainfall estimates from opportunistic sensors in germany across spatio\-temporal scales\.Journal of Hydrology: Regional Studies37,pp\. 100883\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p3.1),[§5\.2](https://arxiv.org/html/2605.05520#S5.SS2.SSS0.Px4.p1.1)\.
- H\. V\. Habi and H\. Messer \(2021\)RainGAN: a conditional rain fields generator\.In2021 IEEE International Conference on Microwaves, Antennas, Communications and Electronic Systems \(COMCAS\),pp\. 529–532\.Cited by:[Appendix B](https://arxiv.org/html/2605.05520#A2.SS0.SSS0.Px4.p2.5),[§3\.3](https://arxiv.org/html/2605.05520#S3.SS3.p3.3)\.
- H\. V\. Habi \(2026\)PyNNcml: rain estimation and classification from commercial microwave links dataExternal Links:[Link](https://github.com/haihabi/PyNNcml)Cited by:[§C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px2.p1.3),[§3\.3](https://arxiv.org/html/2605.05520#S3.SS3.p1.2),[§5\.2](https://arxiv.org/html/2605.05520#S5.SS2.p1.1)\.
- B\. Haese, S\. Hörning, C\. Chwala, A\. Bárdossy, B\. Schalge, and H\. Kunstmann \(2017\)Stochastic reconstruction and interpolation of precipitation fields using combined information of commercial microwave links and rain gauges\.Water Resources Research53\(12\),pp\. 10740–10756\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p4.1),[§2](https://arxiv.org/html/2605.05520#S2.p2.1)\.
- D\. Harrison, S\. Driscoll, and M\. Kitchen \(2000\)Improving precipitation estimates from weather radar using quality control and correction techniques\.Meteorological Applications7\(2\),pp\. 135–144\.Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p1.1)\.
- J\. He, J\. M\. Hernández\-Lobato, Y\. Du, and F\. Vargas \(2025a\)RNE: a plug\-and\-play framework for diffusion density estimation and inference\-time control\.arXiv preprint arXiv:2506\.05668\.Cited by:[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p3.7)\.
- J\. He, P\. Jeha, P\. Potaptchik, L\. Zhang, J\. M\. Hernández\-Lobato, Y\. Du, S\. Syed, and F\. Vargas \(2025b\)CREPE: controlling diffusion with replica exchange\.arXiv preprint arXiv:2509\.23265\.Cited by:[§C\.2](https://arxiv.org/html/2605.05520#A3.SS2.SSS0.Px4.p1.2),[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p3.7),[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px3.p1.1)\.
- J\. Ho, A\. Jain, and P\. Abbeel \(2020\)Denoising diffusion probabilistic models\.Advances in Neural Information Processing Systems33,pp\. 6840–6851\.Cited by:[§3\.2](https://arxiv.org/html/2605.05520#S3.SS2.p1.15),[§3\.2](https://arxiv.org/html/2605.05520#S3.SS2.p1.24)\.
- S\. Hörning, J\. Sreekanth, and A\. Bárdossy \(2019\)Computational efficient inverse groundwater modeling using random mixing and whittaker–shannon interpolation\.Advances in Water Resources123,pp\. 109–119\.Cited by:[§2](https://arxiv.org/html/2605.05520#S2.p2.1)\.
- International Telecommunication Union \(2005\)Recommendation ITU\-R P\.838\-3: specific attenuation model for rain for use in prediction methods\.Technical reportTechnical ReportP\.838\-3,ITU\-R\.Note:Approved in March 2005\.External Links:[Link](https://www.itu.int/rec/R-REC-P.838-3-200503-I/en)Cited by:[§C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px2.p1.3)\.
- Y\. Janati, B\. Moufad, M\. Abou El Qassime, A\. O\. Durmus, E\. Moulines, and J\. Olsson \(2025a\)A mixture\-based framework for guiding diffusion models\.InForty\-second International Conference on Machine Learning,Cited by:[§C\.2](https://arxiv.org/html/2605.05520#A3.SS2.SSS0.Px1.p1.1),[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p3.7),[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px3.p1.1)\.
- Y\. Janati, E\. Moulines, J\. Olsson, and A\. Oliviero\-Durmus \(2025b\)Bridging diffusion posterior sampling and monte carlo methods: a survey\.Philosophical Transactions A383\(2299\),pp\. 20240331\.Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p2.3),[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p1.7),[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p4.1)\.
- T\. Karras, M\. Aittala, T\. Aila, and S\. Laine \(2022\)Elucidating the design space of diffusion\-based generative models\.Advances in Neural Information Processing Systems35,pp\. 26565–26577\.Cited by:[Appendix B](https://arxiv.org/html/2605.05520#A2.p1.1),[§3\.2](https://arxiv.org/html/2605.05520#S3.SS2.p1.2),[§3\.2](https://arxiv.org/html/2605.05520#S3.SS2.p1.24),[§3\.3](https://arxiv.org/html/2605.05520#S3.SS3.p1.2)\.
- B\. Kawar, M\. Elad, S\. Ermon, and J\. Song \(2022\)Denoising diffusion restoration models\.Advances in Neural Information Processing Systems35,pp\. 23593–23606\.Cited by:[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p3.7)\.
- D\. P\. Kingma and J\. Ba \(2014\)Adam: a method for stochastic optimization\.arXiv preprint arXiv:1412\.6980\.Cited by:[Appendix B](https://arxiv.org/html/2605.05520#A2.SS0.SSS0.Px3.p1.8)\.
- W\.F\. Krajewski and J\.A\. Smith \(2002\)Radar hydrology: rainfall estimation\.Advances in Water Resources25\(8\),pp\. 1387–1394\.External Links:ISSN 0309\-1708,[Document](https://dx.doi.org/https%3A//doi.org/10.1016/S0309-1708%2802%2900062-3),[Link](https://www.sciencedirect.com/science/article/pii/S0309170802000623)Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p1.1)\.
- H\. Leijnse, R\. Uijlenhoet, and J\. Stricker \(2007\)Rainfall measurement using radio links from cellular communication networks\.Water resources research43\(3\)\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p1.1),[§1](https://arxiv.org/html/2605.05520#S1.p1.1)\.
- D\. Lopez\-Paz and M\. Oquab \(2016\)Revisiting classifier two\-sample tests\.arXiv preprint arXiv:1610\.06545\.Cited by:[§3\.3](https://arxiv.org/html/2605.05520#S3.SS3.p3.3)\.
- M\. Mardani, J\. Song, J\. Kautz, and A\. Vahdat \(2024\)A variational perspective on solving inverse problems with diffusion models\.InThe Twelfth International Conference on Learning Representations,External Links:[Link](https://openreview.net/forum?id=1YO4EE3SPB)Cited by:[§C\.2](https://arxiv.org/html/2605.05520#A3.SS2.SSS0.Px3.p1.1),[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p3.7),[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px3.p1.1)\.
- G\. Matheron \(1963a\)Principles of geostatistics\.Economic Geology58\(8\),pp\. 1246–1266\.External Links:[Document](https://dx.doi.org/10.2113/gsecongeo.58.8.1246)Cited by:[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px2.p4.1)\.
- G\. Matheron \(1963b\)Principles of geostatistics\.Economic Geology58\(8\),pp\. 1246–1266\.External Links:[Document](https://dx.doi.org/10.2113/gsecongeo.58.8.1246),[Link](https://doi.org/10.2113/gsecongeo.58.8.1246)Cited by:[§C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px3.p1.1)\.
- H\. Messer, A\. Eshel, H\. Habi, S\. Sagiv, and X\. Zheng \(2022\)Rain field retrieval by ground\-level sensors of various types\.Frontiers in Signal Processing2,pp\. 877336\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p2.1),[§1](https://arxiv.org/html/2605.05520#S1.p1.1),[§2](https://arxiv.org/html/2605.05520#S2.p1.1)\.
- H\. Messer, A\. Zinevich, and P\. Alpert \(2006\)Environmental monitoring by wireless communication networks\.Science312\(5774\),pp\. 713–713\.Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p1.1),[§2](https://arxiv.org/html/2605.05520#S2.p1.1)\.
- D\. Michelson and J\. Koistinen \(2000\)Gauge\-radar network adjustment for the baltic sea experiment\.Physics and Chemistry of the Earth, Part B: Hydrology, Oceans and Atmosphere25\(10\),pp\. 915–920\.Note:First European Conference on Radar MeteorologyExternal Links:ISSN 1464\-1909,[Document](https://dx.doi.org/https%3A//doi.org/10.1016/S1464-1909%2800%2900125-8),[Link](https://www.sciencedirect.com/science/article/pii/S1464190900001258)Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p1.1),[§5\.2](https://arxiv.org/html/2605.05520#S5.SS2.SSS0.Px2.p1.3)\.
- B\. Moufad, Y\. Janati, L\. Bedin, A\. O\. Durmus, randal douc, E\. Moulines, and J\. Olsson \(2025\)Variational diffusion posterior sampling with midpoint guidance\.InThe Thirteenth International Conference on Learning Representations,External Links:[Link](https://openreview.net/forum?id=6EUtjXAvmj)Cited by:[§C\.2](https://arxiv.org/html/2605.05520#A3.SS2.SSS0.Px7.p1.1),[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p3.7),[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px3.p1.1),[§5\.1](https://arxiv.org/html/2605.05520#S5.SS1.p1.15)\.
- B\. Murphy, R\. Yurchak, and S\. Müller \(2025\)GeoStat\-framework/pykrige: v1\.7\.3External Links:[Document](https://dx.doi.org/10.5281/zenodo.17372225),[Link](https://doi.org/10.5281/zenodo.17372225)Cited by:[§C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px3.p1.1)\.
- A\. Q\. Nichol and P\. Dhariwal \(2021\)Improved denoising diffusion probabilistic models\.InProceedings of the 38th International Conference on Machine Learning,M\. Meila and T\. Zhang \(Eds\.\),Proceedings of Machine Learning Research, Vol\.139,pp\. 8162–8171\.External Links:[Link](https://proceedings.mlr.press/v139/nichol21a.html)Cited by:[§3\.2](https://arxiv.org/html/2605.05520#S3.SS2.p1.24)\.
- J\. A\. Ordoñez, D\. Bandyopadhyay, V\. H\. Lachos, and C\. R\. B\. Cabral \(2018\)Geostatistical estimation and prediction for censored responses\.Spatial Statistics23,pp\. 109–123\.External Links:[Document](https://dx.doi.org/10.1016/j.spasta.2017.12.001)Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p2.13)\.
- A\. Overeem, H\. Leijnse, and R\. Uijlenhoet \(2013\)Country\-wide rainfall maps from cellular communication networks\.Proceedings of the National Academy of Sciences110\(8\),pp\. 2741–2745\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p1.1),[§1](https://arxiv.org/html/2605.05520#S1.p2.3)\.
- A\. Overeem, H\. Leijnse, and R\. Uijlenhoet \(2016\)Retrieval algorithm for rainfall mapping from microwave links in a cellular communication network\.Atmospheric Measurement Techniques9\(5\),pp\. 2425–2444\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p1.1)\.
- A\. Pakman and L\. Paninski \(2014\)Exact hamiltonian monte carlo for truncated multivariate gaussians\.Journal of Computational and Graphical Statistics23\(2\),pp\. 518–542\.Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p1.1)\.
- M\. Pasierb, Z\. Bałdysz, J\. Szturc, G\. Nykiel, A\. Jurczyk, K\. Ośródka, M\. Figurski, M\. Wojtczak, and C\. Wojtkowski \(2024\)Application of commercial microwave links \(cmls\) attenuation for quantitative estimation of precipitation\.Meteorological Applications31\(3\),pp\. e2218\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p3.1)\.
- J\. Rabin, G\. Peyré, J\. Delon, and M\. Bernot \(2011\)Wasserstein barycenter and its application to texture mixing\.InScale Space and Variational Methods in Computer Vision \(SSVM\),Cited by:[§5\.1](https://arxiv.org/html/2605.05520#S5.SS1.SSS0.Px1.p1.3)\.
- C\. E\. Rasmussen and C\. K\. Williams \(2006\)Gaussian processes for machine learning\.Cambridge, MA, USA: Massachusetts Institute of Technology Publishing\.Cited by:[Appendix E](https://arxiv.org/html/2605.05520#A5.SS0.SSS0.Px1.p1.6),[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px1.p1.1),[§5\.1](https://arxiv.org/html/2605.05520#S5.SS1.p1.6)\.
- D\. Shepard \(1968\)A two\-dimensional interpolation function for irregularly\-spaced data\.InProceedings of the 1968 23rd ACM national conference,pp\. 517–524\.Cited by:[§C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px1.p1.5),[§1](https://arxiv.org/html/2605.05520#S1.p2.3)\.
- R\. L\. Siddon \(1985\)Fast calculation of the exact radiological path for a three\-dimensional ct array\.Medical physics12\(2\),pp\. 252–255\.Cited by:[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px1.p1.20)\.
- J\. Song, C\. Meng, and S\. Ermon \(2021a\)Denoising diffusion implicit models\.InInternational Conference on Learning Representations,External Links:[Link](https://openreview.net/forum?id=St1giarCHLP)Cited by:[§3\.2](https://arxiv.org/html/2605.05520#S3.SS2.p1.15)\.
- Y\. Song, L\. Shen, L\. Xing, and S\. Ermon \(2022\)Solving inverse problems in medical imaging with score\-based generative models\.InInternational Conference on Learning Representations,External Links:[Link](https://openreview.net/forum?id=vaRCHVj0uGI)Cited by:[§4](https://arxiv.org/html/2605.05520#S4.SS0.SSS0.Px2.p3.7)\.
- Y\. Song, J\. Sohl\-Dickstein, D\. P\. Kingma, A\. Kumar, S\. Ermon, and B\. Poole \(2021b\)Score\-based generative modeling through stochastic differential equations\.InInternational Conference on Learning Representations,Cited by:[Appendix B](https://arxiv.org/html/2605.05520#A2.SS0.SSS0.Px2.p1.2),[§3\.2](https://arxiv.org/html/2605.05520#S3.SS2.p1.2)\.
- R\. Stauffer, G\. J\. Mayr, J\. W\. Messner, N\. Umlauf, and A\. Zeileis \(2017\)Spatio\-temporal precipitation climatology over complex terrain using a censored additive regression model\.International Journal of Climatology37\(7\),pp\. 3264–3275\.Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p1.1)\.
- R\. Uijlenhoet, A\. Overeem, and H\. Leijnse \(2018\)Opportunistic remote sensing of rainfall using microwave links from cellular communication networks\.Wiley Interdisciplinary Reviews: Water5\(4\),pp\. e1289\.Cited by:[§1](https://arxiv.org/html/2605.05520#S1.p1.1)\.
- H\. Wheater, V\. Isham, D\. Cox, R\. Chandler, A\. Kakou, P\. Northrop, L\. Oh, C\. Onof, and I\. Rodriguez\-Iturbe \(2000\)Spatial\-temporal rainfall fields: modelling and statistical aspects\.Hydrology and Earth System Sciences4\(4\),pp\. 581–601\.Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p1.1)\.
- D\. S\. Wilks and R\. L\. Wilby \(1999\)The weather generation game: a review of stochastic weather models\.Progress in physical geography23\(3\),pp\. 329–357\.Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p1.1)\.
- L\. Wu, B\. Trippe, C\. Naesseth, D\. Blei, and J\. P\. Cunningham \(2023\)Practical and asymptotically exact conditional sampling in diffusion models\.Advances in Neural Information Processing Systems36,pp\. 31372–31403\.Cited by:[§C\.2](https://arxiv.org/html/2605.05520#A3.SS2.SSS0.Px6.p1.3),[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px3.p1.1)\.
- C\. Yang, R\. Chandler, V\. Isham, and H\. Wheater \(2005\)Spatial\-temporal rainfall simulation using generalized linear models\.Water Resources Research41\(11\)\.Cited by:[§3\.1](https://arxiv.org/html/2605.05520#S3.SS1.p1.1)\.
- B\. Zhang, W\. Chu, J\. Berner, C\. Meng, A\. Anandkumar, and Y\. Song \(2024\)Improving diffusion inverse problem solving with decoupled noise annealing\.arXiv preprint arXiv:2407\.01521\.Cited by:[§C\.2](https://arxiv.org/html/2605.05520#A3.SS2.SSS0.Px2.p1.1),[§5](https://arxiv.org/html/2605.05520#S5.SS0.SSS0.Px3.p1.1)\.
- P\. Zhang, X\. Liu, and M\. Zou \(2023\)Reconstructing and nowcasting the rainfall field by a cml network\.Earth and Space Science10\(9\),pp\. e2023EA002909\.Cited by:[Appendix A](https://arxiv.org/html/2605.05520#A1.p3.1),[§2](https://arxiv.org/html/2605.05520#S2.p1.1),[§5\.2](https://arxiv.org/html/2605.05520#S5.SS2.SSS0.Px4.p1.1)\.
## Appendix AExtended Related Work
As discussed in Section[2](https://arxiv.org/html/2605.05520#S2), earlier studies on rainfall reconstruction from CMLs\(Leijnseet al\.,[2007](https://arxiv.org/html/2605.05520#bib.bib148); Goldshteinet al\.,[2009a](https://arxiv.org/html/2605.05520#bib.bib147); Overeemet al\.,[2013](https://arxiv.org/html/2605.05520#bib.bib150),[2016](https://arxiv.org/html/2605.05520#bib.bib149)\)typically reduced each path\-integrated link observation to a VRG located at the link midpoint and then produced gridded rain fields via deterministic interpolation \(most commonlyIDWorOK\)\. This abstraction enabled the first large\-scale deployments of CML\-based mapping, but it discards variability along the link and effectively converts line\-integral measurements into point\-scale information, which can yield overly smooth fields and miss localized extremes\.
Recent work has systematically examined the implications of this approximation and the extent to which more faithful representations of CML measurements can improve reconstruction accuracy\.Messeret al\.\([2022](https://arxiv.org/html/2605.05520#bib.bib139)\)quantify the information gain of line\-integral observations relative to point measurements, demonstrating that CMLs contain additional spatial information that is lost when collapsed to a single point\. Building on this insight,Eshelet al\.\([2021](https://arxiv.org/html/2605.05520#bib.bib138)\)evaluate rainfall reconstruction using a dense operational CML network in southern Germany and a semisynthetic framework based on gauge\-adjusted radar data\. Their results confirm that reconstruction accuracy strongly depends on rainfall spatial structure, improving for smoother fields with larger decorrelation distances\. Importantly, they show that representing each CML by multiple VRGs using GMZ yields a consistent, though modest, improvement over the single\-midpoint VRG representation, largely independent of whetherIDWorOKis applied\.
Similar conclusions emerge from other deterministic reconstruction analyses\.Zhanget al\.\([2023](https://arxiv.org/html/2605.05520#bib.bib135)\)investigate two\-dimensional rainfall fields reconstruction by converting attenuation measurements into path\-averaged rain rates via power\-law relationships and then interpolating them withIDWorOK\. Their analysis shows that both methods recover the dominant spatial structures, withOKtypically outperformingIDWdue to its explicit modeling of spatial dependence\. However, reconstruction quality deteriorates in sparsely covered regions and during high\-intensity events, reflecting fundamental limitations imposed by network geometry and path integration\. Moving beyond CML\-only settings,Grafet al\.\([2021](https://arxiv.org/html/2605.05520#bib.bib140)\)demonstrate that jointly assimilating rain gauges, personal weather stations, and CMLs within a block\-kriging framework substantially increases effective observation density and improves representation of spatial rainfall variability\. Crucially, their approach accounts for differences in spatial support and sensor\-specific uncertainties, yielding performance comparable to or exceeding gauge\-adjusted radar at local and regional scales\. Likewise,Pasierbet al\.\([2024](https://arxiv.org/html/2605.05520#bib.bib141)\)emphasize the sensitivity of CML\-based reconstruction to link characteristics \(e\.g\., length, frequency, orientation\) as well as to the strategy used to project path\-averaged measurements onto spatial grids, underscoring that preprocessing and mapping choices play a critical role in reconstruction quality\.
While deterministic interpolation provides a single best estimate, it offer limited uncertainty quantification and may further oversmooth fine\-scale variability when constraints are weak\. To address these limitations, stochastic reconstruction approaches treat rainfall as a random spatial process conditioned on both point and path\-integrated information\. A key advance in this direction is the application of the RM framework byHaeseet al\.\([2017](https://arxiv.org/html/2605.05520#bib.bib136)\), which generates an ensemble of rainfall fields that agree with rain gauges and with the path\-averaged link measurements, so it keeps the fact that CMLs measure rainfall along a path rather than at a point\. Compared toOK, RM can better preserve small\-scale spatial variability and represent higher intensities \(including extreme thresholds\)\. However, because the marginal rainfall distribution is derived primarily from gauge data, RM remains sensitive to gauge representativeness and may underestimate extremes or delay rainfall onset when gauges inadequately sample localized events\.
The RM framework has since been extended to larger spatial domains\.Blettneret al\.\([2022](https://arxiv.org/html/2605.05520#bib.bib137)\)extend stochastic reconstruction to national domains by combining rain gauges and CMLs across an entire country, improving agreement with gauge\-adjusted radar and reducing structural errors in rainfall pattern representation compared to deterministic kriging\. Their uncertainty decomposition \(amplitude, structure, location\) highlights the diagnostic value of ensemble reconstructions\. Nonetheless, increased computational cost, interpretative complexity, and the smoothing implied by path\-integrated constraints continue to limit the recoverable resolution of fine\-scale precipitation features, particularly where observational coverage is sparse\.
## Appendix BPrior Diffusion Model
We base the prior model on the EDM framework\(Karraset al\.,[2022](https://arxiv.org/html/2605.05520#bib.bib89)\)and adapt the publicly released implementation222[https://github\.com/NVlabs/edm](https://github.com/NVlabs/edm)to our setting\.
#### Dataset\.
The OpenMRG dataset contains26,49526\{,\}495rain\-rate fields, which are split chronologically into training and test sets with a0\.8/0\.20\.8/0\.2ratio, yielding21,19621\{,\}196and5,2995\{,\}299samples, respectively\. Unlike the original EDM setup, we operate directly on rain\-rate fields without converting them to PNG images\. Each field has shape\(1,48,36\)\(1,48,36\); the width is cropped from3737to3636for consistency\. To mitigate the long\-tailed distribution of rainfall intensities, we experimented with two approaches:
*Log\-transform\.*We initially experimented with training on a log\-transformed version of the rainfall data, which amounts to training a latent diffusion model where the encoder applies a log\-transform and the decoder applies an exponential transform\. While training in log\-space was significantly faster \(roughly×1\.5\\times 1\.5\), we ultimately opted against using the log\-space model during inference\. The primary reason is that mapping back to rain intensities via an exponential transform introduces an additional nonlinearity into the inverse problem\. As noted in recent work on inverse problems with latent diffusion models, such nonlinearities can introduce numerical instabilities for diffusion samplers; in our case, the exp\-transform led to instabilities for some of the diffusion baselines, namely forDPSandDAPS\.
*Linear\-transform with quantiles\.*All samples are normalized by a fixed constant corresponding to the0\.9990\.999quantile of the training data\. We retained this as our final preprocessing choice: this linear normalization mitigates skewed values while preserving linearity in the input space, which we found more stable for inference\.
#### Architecture\.
We use DDPM\+\+ backbone\(Songet al\.,[2021b](https://arxiv.org/html/2605.05520#bib.bib86)\)\(SongUNet\), a U\-Net encoder–decoder with skip connections and timestep embeddings\. The model usesmodel\_channels=32 anddropout=0\.1\. The denoiser is wrapped with EDM preconditioning \(EDMPrecond\) andσdata=0\.079\\sigma\_\{\\mathrm\{data\}\}=0\.079is used in both the preconditioner and the EDM loss\. A ReLU activation is applied to the denoiser output to enforce the physical constraint that the denoiser expectation𝔼\[X0∣Xt\]\\mathbb\{E\}\[X\_\{0\}\\mid X\_\{t\}\]be non\-negative\.
*Enforcing non\-negative values\.*While the ReLU activation function is a valid heuristic to enforce non\-negative values, it introduces a hard thresholding that can distort gradients near zero\-intensity regions\. We also experimented with a smoother alternative, namely the Softplus activation function333[https://docs\.pytorch\.org/docs/stable/generated/torch\.nn\.Softplus\.html](https://docs.pytorch.org/docs/stable/generated/torch.nn.Softplus.html), which yields smoother gradients around zero\. In practice, however, we found no significant difference in the performance of the diffusion prior between the two variants, indicating that our results are largely insensitive to whether non\-negativity is enforced by a hard or smooth transformation\.
#### Training Setup\.
Training is performed using Adam\(Kingma and Ba,[2014](https://arxiv.org/html/2605.05520#bib.bib32)\)withβ1=0\.9\\beta\_\{1\}=0\.9,β2=0\.999\\beta\_\{2\}=0\.999, and a constant learning rate of10−410^\{\-4\}for approximately5,2845\{,\}284epochs\. Data augmentation is limited to random axis flips, with each sample independently flipped along thexxandyyaxes with probability0\.10\.1\. Exponential moving average updates are disabled\. The training is conducted on two NVIDIA H100 \(80GB\) GPUs, with a total training time of approximately1010hours\.
#### Discriminator Training Details\.
The discriminator is trained as a supervised binary classifier to distinguish real from generated samples, with labels “1” and “0”, respectively\. The two datasets are concatenated, randomly shuffled, and split into training and test sets using a0\.8/0\.20\.8/0\.2ratio\. Optimization is performed using mini\-batches of size256256\.
The discriminator is a convolutional network operating on single\-channel inputs of shape\(1,48,36\)\(1,48,36\)\. Its architecture followsHabi and Messer \([2021](https://arxiv.org/html/2605.05520#bib.bib180)\)and consists of three convolutional blocks\. Each block applies a5×55\\times 5convolution with stride22and padding22, followed by a LeakyReLU activation\. Channel dimensions increase as1→d→2d→4d1\\to d\\to 2d\\to 4d\. The resulting feature maps are flattened and passed through a linear layer to produce a single logit\.
Training minimizes the binary cross\-entropy loss with logits using the Adam optimizer with a learning rate of10−410^\{\-4\}for200200epochs\. At evaluation time, predicted probabilities are thresholded at0\.50\.5to obtain class labels\.


Figure 5:Power\-law parameters of the CML as a function of link’s length and the frequency on the OpenMRG dataset\.
## Appendix CBaselines Implementation
In this section, we provide implementation details for all baseline methods against which we compare our approach\. For the meteorological baselines, hyperparameters are presented in[Table6](https://arxiv.org/html/2605.05520#A3.T6)\. The details about the hyperparameters of DM posterior samplers are summarized in[Table7](https://arxiv.org/html/2605.05520#A3.T7)\. We highlight that all baseline hyperparameters were tuned manually\. For completeness, we also report the runtime required to generate one reconstruction for each baseline in[Table2](https://arxiv.org/html/2605.05520#S5.T2)\. The meteorological baselines are substantially faster than the diffusion baselines, for which the reported runtimes are normalized by batch size, since 10 samples are generated jointly and the total runtime is divided by 10 to obtain a per\-reconstruction comparison\.
### C\.1Meteorological baselines
#### IDW\.
We implement inverse\-distance weighting \(IDW\) followingEshelet al\.\([2021](https://arxiv.org/html/2605.05520#bib.bib138)\); Shepard \([1968](https://arxiv.org/html/2605.05520#bib.bib168)\)\. Each CML link observation is represented as a point measurement located at the link midpoint\. For every grid locationss, we compute a normalized weighted average usingp=2p=2andϵ=10−6\\epsilon=10^\{\-6\}for numerical stability\. A region\-of\-influence radiusroi\\mathrm\{roi\}is set to44in order to discard distant observations\.
#### GMZ\.
We implement the iterative GMZ reconstruction procedure followingGoldshteinet al\.\([2009a](https://arxiv.org/html/2605.05520#bib.bib147)\); Eshelet al\.\([2021](https://arxiv.org/html/2605.05520#bib.bib138)\), that we adapt fromPyNNcml\(Habi,[2026](https://arxiv.org/html/2605.05520#bib.bib170)\)\. We first map each link measurement to the attenuation domain via the standard power lawA=aXbA=aX^\{b\}\(International Telecommunication Union,[2005](https://arxiv.org/html/2605.05520#bib.bib60)\)\. Each link is discretized into a fixed number of uniformly spaced VRGs\. We useK=5K=5virtual points per link and run the algorithm forT=20T=20iterations\. At each iteration, an intermediate field is reconstructed using the same IDW operator as in[C\.1](https://arxiv.org/html/2605.05520#A3.SS1.SSS0.Px1)\. The reconstructed field is bilinearly sampled at the virtual points, and link consistency is enforced by projecting the sampled values such that the mean attenuation along each link matches the original transformed measurement\.
#### OK\.
We apply ordinary kriging algorithm\(Matheron,[1963b](https://arxiv.org/html/2605.05520#bib.bib59); Cressie,[1993b](https://arxiv.org/html/2605.05520#bib.bib129)\)usingOrdinaryKrigingfunction fromPyKrigepackage\(Murphyet al\.,[2025](https://arxiv.org/html/2605.05520#bib.bib171)\)\. We use the Exponential variogram model, with parameters inferred using PyKrige’s default L1\-norm\-based variogram fitting\.
### C\.2Diffusion posterior sampling baselines
#### MGDM\.
We implement theMGDMsampler proposed byJanatiet al\.\([2025a](https://arxiv.org/html/2605.05520#bib.bib112)\)in Algorithm22, by adapting the reference implementation444[https://github\.com/badr\-moufad/mgdm](https://github.com/badr-moufad/mgdm)\. In particular, the code is adapted to the VE setting\. We sample the intermediate steps uniformly and use one Gibbs sampling set\.
#### DAPS\.
We refer toZhanget al\.\([2024](https://arxiv.org/html/2605.05520#bib.bib99), Algo 1\)for the method used in theDAPSsampler implemented based on the released code555[https://github\.com/zhangbingliang2019/DAPS](https://github.com/zhangbingliang2019/DAPS)\. We use the Langevin MCMC sampler for consistency and adapted to the VE setting\.
#### RedDiff\.
We implement theRedDiffsampler as described inMardaniet al\.\([2024](https://arxiv.org/html/2605.05520#bib.bib98), Algo 1\)\. We modify the implementation provided in released repository666[https://github\.com/NVlabs/RED\-diff](https://github.com/NVlabs/RED-diff)to the VE setting\.
#### CREPE\.
We follow the algorithm presented inHeet al\.\([2025b](https://arxiv.org/html/2605.05520#bib.bib90), Algo 1\)combined with the reference code777[https://github\.com/jiajunhe98/CREPE\-Controlling\-diffusion\-with\-REPlica\-Exchange](https://github.com/jiajunhe98/CREPE-Controlling-diffusion-with-REPlica-Exchange)provided by the authors\. We adapted theCREPEsampler to our use case of VE formulation\. Following the authors, we use the prior model as a proposal and use theDPSintermediate likelihoods \(also referred to as intermediate rewards\) to define the annealing path that bridge the GaussianpTp\_\{T\}and the targetπ0\\pi\_\{0\}\.
#### DPS\.
#### TDS\.
We implementWuet al\.\([2023](https://arxiv.org/html/2605.05520#bib.bib100), Algorithm 1\)in the VE setting\. As a proposal distribution, we useDPS\(Chunget al\.,[2023](https://arxiv.org/html/2605.05520#bib.bib102), Algorithm 1\), which differs from the original proposal by normalizing the log\-gradient of the intermediate likelihood\. For the variance of the intermediate likelihoods, we use the surrogateσi2\+\(τσt\)2\\sigma\_\{i\}^\{2\}\+\(\\tau\\sigma\_\{t\}\)^\{2\}, whereτ\\tauis a hyperparameter set toτ=1\\tau=1, as it yields the best empirical performance\.
#### MGPS\.
We implement theMGPS\(Moufadet al\.,[2025](https://arxiv.org/html/2605.05520#bib.bib101), Algo 1\)based on the released code999[https://github\.com/YazidJanati/mgps](https://github.com/YazidJanati/mgps)by adapting it to our VE setting\. We useη=0\.5\\eta=0\.5for selecting the midpoint and we do not use warm start for the algorithm\.
Table 6:The hyperparameters for each meteorological algorithm\. The symbol “–” indicates that the corresponding hyperparameter is not applicable to the given method\.AlgorithmpROI\# Points per link\# IterationsVariogramIDW261––GMZ26520–OK––––exponential
Table 7:The hyperparameters for each posterior sampling algorithm\. We use Karras scheduler withσmin=2×10−3\\sigma\_\{\\min\}=2\\times 10^\{\-3\}andσmax=100\\sigma\_\{\\max\}=100for GP experiments andσmin=2×10−3\\sigma\_\{\\min\}=2\\times 10^\{\-3\}andσmax=80\\sigma\_\{\\max\}=80for rain field experiments\.AlgorithmBase hyperparametersTasksGaussian ProcessReal CML linksDAPSnsteps=100n\_\{\\text\{steps\}\}=100Node=5N\_\{\\text\{ode\}\}=5Min ratio=0\.01\\text\{\{Min ratio\}\}=0\.01MCMC sampler=Langevin\\text\{\{MCMC sampler\}\}=\\text\{\{Langevin\}\}ρ=7\\rho=7MCMC steps=100\\text\{\{MCMC steps\}\}=100η0=5×10−4\\eta\_\{0\}=5\\times 10^\{\-4\}MCMC steps=50\\text\{\{MCMC steps\}\}=50η0=2×10−4\\eta\_\{0\}=2\\times 10^\{\-4\}DPSρ=7\\rho=7nsteps=320n\_\{\\text\{steps\}\}=320γ=4\\gamma=4nsteps=420n\_\{\\text\{steps\}\}=420γ=1\\gamma=1MGDMρ=7\\rho=7lr=3×10−2\\text\{\{lr\}\}=3\\times 10^\{\-2\}n\_gibbs\_repetitions=1\\text\{\{n\\\_gibbs\\\_repetitions\}\}=1n\_gradient\_repetitions=10\\text\{\{n\\\_gradient\\\_repetitions\}\}=10nsteps=32n\_\{\\text\{steps\}\}=32n\_denoising\_steps=3\\text\{\{n\\\_denoising\\\_steps\}\}=3nsteps=50n\_\{\\text\{steps\}\}=50n\_denoising\_steps=5\\text\{\{n\\\_denoising\\\_steps\}\}=5MGPSρ=7\\rho=7lr=3×10−2\\text\{\{lr\}\}=3\\times 10^\{\-2\}n\_gradient\_steps=10\\text\{\{n\\\_gradient\\\_steps\}\}=10nsteps=64n\_\{\\text\{steps\}\}=64nsteps=32n\_\{\\text\{steps\}\}=32RedDiffnsteps=1000n\_\{\\text\{steps\}\}=1000ρ=5\\rho=5obs\_weight=1\\text\{\{obs\\\_weight\}\}=1grad\_term\_weight=1\\text\{\{grad\\\_term\\\_weight\}\}=1lr=10−1\\text\{\{lr\}\}=10^\{\-1\}lr=5×10−3\\text\{\{lr\}\}=5\\times 10^\{\-3\}TDSρ=7\\rho=7nsteps=320n\_\{\\text\{steps\}\}=320γ=4\\gamma=4n\_particles=10\\text\{\{n\\\_particles\}\}=10nsteps=420n\_\{\\text\{steps\}\}=420γ=1\\gamma=1n\_particles=4\\text\{\{n\\\_particles\}\}=4
## Appendix DDetails on Siddon Algorithm

⇒Δi=\[0000000000\.30\.70001\.10\.9000\.31\.20000\]\\Rightarrow\\Delta^\{i\}=\\begin\{bmatrix\}0&0&0&0&0&0\\\\ 0&0&0&0\.3&0\.7&0\\\\ 0&0&1\.1&0\.9&0&0\\\\ 0\.3&1\.2&0&0&0&0\\\\ \\end\{bmatrix\}
Figure 6:Illustration of a line sensor on a6×46\\times 4grid\. Theblue dotsdefine the line sensor whereas thered dotsmark its intersections with the two\-dimensional grid\. The color intensity of each cell is proportional to the intersection length\.We provide implementation details for computing the intersection\-length weightsΔi\\Delta^\{i\}in[Equation4](https://arxiv.org/html/2605.05520#S4.E4), i\.e\., the lengths of intersection between theii\-th CML path and theH×WH\\times Wgrid cells\.
#### Grid parametrization\.
Consider a regularly spaced grid with vertical and horizontal grid lines located at
x=xref\+cΔx,y=yref\+rΔy,c,r∈ℤ,x=x\_\{\\mathrm\{ref\}\}\+c\\,\\Delta x,\\qquad y=y\_\{\\mathrm\{ref\}\}\+r\\,\\Delta y,\\qquad c,r\\in\\mathbb\{Z\},whereΔx,Δy\>0\\Delta x,\\Delta y\>0are the grid spacings and\(xref,yref\)\(x\_\{\\mathrm\{ref\}\},y\_\{\\mathrm\{ref\}\}\)specifies the grid origin\. Under the affine change of variablesx~:=\(x−xref\)/Δx\\tilde\{x\}\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\(x\-x\_\{\\mathrm\{ref\}\}\)/\\Delta xandy~:=\(y−yref\)/Δy\\tilde\{y\}\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\(y\-y\_\{\\mathrm\{ref\}\}\)/\\Delta y, these grid lines map to integers defined by\(x~,y~\)\(\\tilde\{x\},\\tilde\{y\}\)\. Hence, intersection computations can be carried out in the normalized coordinates and mapped back to the original coordinate system\.
*Considered convention\.*In our setting, the radar maps use the origin\(−12,−12\)\(\-\\tfrac\{1\}\{2\},\-\\tfrac\{1\}\{2\}\)and unit spacing in both directions, i\.e\.,Δx=Δy=1\\Delta x=\\Delta y=1\. Indexing cells by the integer coordinates of their centers\(c,r\)∈ℤ2\(c,r\)\\in\\mathbb\{Z\}^\{2\}places the cell boundaries atc±12c\\pm\\tfrac\{1\}\{2\}andr±12r\\pm\\tfrac\{1\}\{2\}; equivalently, cell centers lie at integer coordinates\.
Therefore, the grid lines of the consideredH×WH\\times Wdiscretization are located at
x=c\+12,y=r\+12,c∈\{−1,…,W−1\},andr∈\{−1,…,H−1\}\.x=c\+\\tfrac\{1\}\{2\},\\qquad y=r\+\\tfrac\{1\}\{2\},\\quad c\\in\\\{\-1,\\dots,W\-1\\\},\\ \\text\{and\}\\ r\\in\\\{\-1,\\dots,H\-1\\\}\.
#### Line parameterization\.
For a given link, let𝐬0,𝐬1∈ℝ2\\mathbf\{s\}\_\{0\},\\mathbf\{s\}\_\{1\}\\in\\mathbb\{R\}^\{2\}denote its start and end points, and define the affine interpolation
𝐬t:=𝐬0\+t\(𝐬1−𝐬0\),t∈\[0,1\]\.\\mathbf\{s\}\_\{t\}\\;\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\\;\\mathbf\{s\}\_\{0\}\+t\(\\mathbf\{s\}\_\{1\}\-\\mathbf\{s\}\_\{0\}\),\\qquad t\\in\[0,1\]\.We denote coordinates by\[𝐬t\]x\[\\mathbf\{s\}\_\{t\}\]\_\{x\}and\[𝐬t\]y\[\\mathbf\{s\}\_\{t\}\]\_\{y\}; in particular,\[𝐬0\]x,\[𝐬0\]y\[\\mathbf\{s\}\_\{0\}\]\_\{x\},\[\\mathbf\{s\}\_\{0\}\]\_\{y\}and\[𝐬1\]x,\[𝐬1\]y\[\\mathbf\{s\}\_\{1\}\]\_\{x\},\[\\mathbf\{s\}\_\{1\}\]\_\{y\}for the start and end points respectively\.
#### Algorithm\.
Let𝐝:=𝐬1−𝐬0\\mathbf\{d\}\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\\mathbf\{s\}\_\{1\}\-\\mathbf\{s\}\_\{0\}withdx:=\[𝐝\]xd\_\{x\}\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\[\\mathbf\{d\}\]\_\{x\}anddy:=\[𝐝\]yd\_\{y\}\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\[\\mathbf\{d\}\]\_\{y\}\. We seek parameterst∈\(0,1\)t\\in\(0,1\)such that𝐬t\\mathbf\{s\}\_\{t\}lies on a grid line\. For each vertical grid linex=c\+12x=c\+\\tfrac\{1\}\{2\}withdx≠0d\_\{x\}\\neq 0, the corresponding parameter is
tcx:=\(c\+12\)−\[𝐬0\]xdxt^\{x\}\_\{c\}\\;\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\\;\\frac\{\\big\(c\+\\tfrac\{1\}\{2\}\\big\)\-\[\\mathbf\{s\}\_\{0\}\]\_\{x\}\}\{d\_\{x\}\}and we retain those satisfyingtkx∈\(0,1\)t^\{x\}\_\{k\}\\in\(0,1\)\. Analogously, for each horizontal grid liney=r\+12y=r\+\\tfrac\{1\}\{2\}withdy≠0d\_\{y\}\\neq 0,
try:=\(r\+12\)−\[𝐬0\]ydy,t^\{y\}\_\{r\}\\;\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\\;\\frac\{\\big\(r\+\\tfrac\{1\}\{2\}\\big\)\-\[\\mathbf\{s\}\_\{0\}\]\_\{y\}\}\{d\_\{y\}\},retaining those withtcy∈\(0,1\)t^\{y\}\_\{c\}\\in\(0,1\)\. We then form
T:=\{0,1\}∪\{tcx:tcx∈\(0,1\)\}∪\{try:try∈\(0,1\)\},T\\;\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\\;\\\{0,1\\\}\\;\\cup\\;\\\{t^\{x\}\_\{c\}\\mathrel\{\\mathop\{\\ordinarycolon\}\}\\ t^\{x\}\_\{c\}\\in\(0,1\)\\\}\\;\\cup\\;\\\{t^\{y\}\_\{r\}\\mathrel\{\\mathop\{\\ordinarycolon\}\}\\ t^\{y\}\_\{r\}\\in\(0,1\)\\\},where the set union removes duplicates \(e\.g\., when the line crosses a grid corner\)\. SortingTTyields
0=t0<t1<⋯<tN=1\.0=t\_\{0\}<t\_\{1\}<\\cdots<t\_\{N\}=1\.Consecutive parameters define sub\-segments that lie within single grid cells\. Their Euclidean lengths are
δn:=∥𝐬tn\+1−𝐬tn∥2=∥𝐬1−𝐬0∥2\(tn\+1−tn\),n=0,…,N−1\.\\delta\_\{n\}\\;\\mathrel\{\\mathop\{\\ordinarycolon\}\}=\\;\\\|\\mathbf\{s\}\_\{t\_\{n\+1\}\}\-\\mathbf\{s\}\_\{t\_\{n\}\}\\\|\_\{2\}\\;=\\;\\\|\\mathbf\{s\}\_\{1\}\-\\mathbf\{s\}\_\{0\}\\\|\_\{2\}\\,\(t\_\{n\+1\}\-t\_\{n\}\),\\qquad n=0,\\dots,N\-1\.We store the weights as a sparse matrixΔi∈ℝ≥0H×W\\Delta^\{i\}\\in\\mathbb\{R\}^\{H\\times W\}\_\{\\geq 0\}, where\[Δi\]r,c\[\\Delta^\{i\}\]\_\{r,c\}equals the length of the portion of the link in cell\(c,r\)\(c,r\)\. VectorizingΔi\\Delta^\{i\}yields the coefficients\{Δki\}k=1HW\\\{\\Delta^\{i\}\_\{k\}\\\}\_\{k=1\}^\{HW\}used in[Equation4](https://arxiv.org/html/2605.05520#S4.E4)\. An example of the resulting sparse weight matrix is shown in[Figure6](https://arxiv.org/html/2605.05520#A4.F6)\.
## Appendix EClosed\-form Expression of the Oracle Posterior
#### Setting\.
We recall that the considered setting: a linear inverse problem where the priorXXis a GPX∼𝒢𝒫\(0,k\)X\\sim\\mathcal\{GP\}\(0,k\)on\[−5,5\]\[\-5,5\], with the RBF kernel
k\(s,s′\)=exp\(−\(s−s′\)22ℓ2\),ℓ=0\.6\.k\(s,s^\{\\prime\}\)=\\exp\\\!\\left\(\-\\frac\{\(s\-s^\{\\prime\}\)^\{2\}\}\{2\\ell^\{2\}\}\\right\),\\qquad\\ell=0\.6\.\(10\)The inverse problem is defined by
Yi=∫aibiX\(s\)ds\+σZi,Zi∼𝒩\(0,1\),Y\_\{i\}=\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}X\(s\)\\,\\mathrm\{d\}s\+\\sigma Z\_\{i\},\\qquad Z\_\{i\}\\sim\\mathcal\{N\}\(0,1\),\(11\)where\{\[ai,bi\]\}i=1m\\\{\[a\_\{i\},b\_\{i\}\]\\\}\_\{i=1\}^\{m\}are the considered integration intervals and\{Zi\}i=1m\\\{Z\_\{i\}\\\}\_\{i=1\}^\{m\}are mutually independent and independent ofXX\. Since the inverse problem involves linear transformations of the GP prior, the posterior is also a GP with closed\-form mean and covariance\(Rasmussen and Williams,[2006](https://arxiv.org/html/2605.05520#bib.bib113)\)that writes
μX\|y\(s\)\\displaystyle\\mu\_\{X\|y\}\(s\)=𝐤y\(s\)⊤\(𝐊yy\+σ2𝐈\)−1𝐲,\\displaystyle=\\mathbf\{k\}\_\{y\}\(s\)^\{\\top\}\(\\mathbf\{K\}\_\{yy\}\+\\sigma^\{2\}\\mathbf\{I\}\)^\{\-1\}\\mathbf\{y\},kX\|y\(s,s′\)\\displaystyle k\_\{X\|y\}\(s,s^\{\\prime\}\)=k\(s,s′\)−𝐤y\(s\)⊤\(𝐊yy\+σ2𝐈\)−1𝐤y\(s′\),\\displaystyle=k\(s,s^\{\\prime\}\)\-\\mathbf\{k\}\_\{y\}\(s\)^\{\\top\}\(\\mathbf\{K\}\_\{yy\}\+\\sigma^\{2\}\\mathbf\{I\}\)^\{\-1\}\\mathbf\{k\}\_\{y\}\(s^\{\\prime\}\),where𝐤y\(s\)∈ℝm\\mathbf\{k\}\_\{y\}\(s\)\\in\\mathbb\{R\}^\{m\}is the covariance betweenX\(s\)X\(s\)and theii\-th observationyiy\_\{i\}, namelycov\(X\(s\),yi\)\\mathrm\{cov\}\(X\(s\),y\_\{i\}\); and𝐊yy∈ℝm×m\\mathbf\{K\}\_\{yy\}\\in\\mathbb\{R\}^\{m\\times m\}covariance between the observations\[y1,…,ym\]\[y\_\{1\},\\ldots,y\_\{m\}\], namely\[𝐊yy\]ij\[\\mathbf\{K\}\_\{yy\}\]\_\{ij\}is the covariance between theii\-th andjj\-th observationscov\(yi,yj\)\\mathrm\{cov\}\(y\_\{i\},y\_\{j\}\)\. Getting the closed\-form expression of the posterior GP amounts to deriving the expressions of\[𝐤y\(s\)\]i\[\\mathbf\{k\}\_\{y\}\(s\)\]\_\{i\}and\[𝐊yy\]ij\[\\mathbf\{K\}\_\{yy\}\]\_\{ij\}for alli,ji,jin⟦1,m⟧\\llbracket 1,m\\rrbracket\.
#### Derivation of the covariances\.
The covariance between theii\-th observation andX\(s\)X\(s\)writes as
\[𝐤y\(s\)\]i\\displaystyle\[\\mathbf\{k\}\_\{y\}\(s\)\]\_\{i\}=cov\(X\(s\),∫aibiX\(s′\)ds′\+σZi\)\\displaystyle=\\mathrm\{cov\}\\left\(X\(s\),\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}X\(s^\{\\prime\}\)\\mathrm\{d\}s^\{\\prime\}\+\\sigma Z\_\{i\}\\right\)=∫aibicov\(X\(s\),X\(s′\)\)ds′\\displaystyle=\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}\\mathrm\{cov\}\\big\(X\(s\),X\(s^\{\\prime\}\)\\big\)\\mathrm\{d\}s^\{\\prime\}=∫aibik\(s,s′\)ds′\\displaystyle=\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}k\(s,s^\{\\prime\}\)\\mathrm\{d\}s^\{\\prime\}where the second equality follows from the linearity of the integral and the independence betweenXXandZiZ\_\{i\}\. The value\[𝐤y\(s\)\]i\[\\mathbf\{k\}\_\{y\}\(s\)\]\_\{i\}can be numerically approximated using the Error Functionerf\\mathrm\{erf\}101010[https://docs\.pytorch\.org/docs/stable/generated/torch\.erf\.html](https://docs.pytorch.org/docs/stable/generated/torch.erf.html)
\[𝐤y\(s\)\]i=ℓπ2\[erf\(bi−s2ℓ\)−erf\(ai−s2ℓ\)\]\.\[\\mathbf\{k\}\_\{y\}\(s\)\]\_\{i\}=\\ell\\sqrt\{\\frac\{\\pi\}\{2\}\}\\left\[\\mathrm\{erf\}\\left\(\\frac\{b\_\{i\}\-s\}\{\\sqrt\{2\}\\ell\}\\right\)\-\\mathrm\{erf\}\\left\(\\frac\{a\_\{i\}\-s\}\{\\sqrt\{2\}\\ell\}\\right\)\\right\]\.\(12\)Similarly, the covariance between theii\-th andjj\-th observations writes
\[𝐊yy\]ij\\displaystyle\[\\mathbf\{K\}\_\{yy\}\]\_\{ij\}=cov\(∫aibiX\(s\)ds\+σZi,∫ajbjX\(s′\)ds′\+σZj\)\\displaystyle=\\mathrm\{cov\}\\left\(\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}X\(s\)\\mathrm\{d\}s\+\\sigma Z\_\{i\},\\int\_\{a\_\{j\}\}^\{b\_\{j\}\}X\(s^\{\\prime\}\)\\mathrm\{d\}s^\{\\prime\}\+\\sigma Z\_\{j\}\\right\)=∫aibi∫ajbjcov\(X\(s\),X\(s′\)\)dsds′\+σ2δi=j\\displaystyle=\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}\\int\_\{a\_\{j\}\}^\{b\_\{j\}\}\\mathrm\{cov\}\\big\(X\(s\),X\(s^\{\\prime\}\)\\big\)\\mathrm\{d\}s\\mathrm\{d\}s^\{\\prime\}\+\\sigma^\{2\}\\delta\_\{i=j\}=∫aibi∫ajbjk\(s,s′\)dsds′\+σ2δi=j,\\displaystyle=\\int\_\{a\_\{i\}\}^\{b\_\{i\}\}\\int\_\{a\_\{j\}\}^\{b\_\{j\}\}k\(s,s^\{\\prime\}\)\\mathrm\{d\}s\\mathrm\{d\}s^\{\\prime\}\\ \+\\ \\sigma^\{2\}\\delta\_\{i=j\},whereδi=j\\delta\_\{i=j\}equals11ifi=ji=jand zero otherwise\. The first term on the right\-hand side can be computed by: \(i\) using the result in[Equation12](https://arxiv.org/html/2605.05520#A5.E12), and \(ii\) integrating by leveraging the formula of the primitive ofs↦erf\(s/2ℓ\)s\\mapsto\\mathrm\{erf\}\(\{s\}/\{\\sqrt\{2\\ell\}\}\), namely
H\(z\)=z⋅erf\(z2ℓ\)\+2πℓexp\(−z22ℓ2\)\.H\(z\)=z\\cdot\\text\{erf\}\\left\(\\frac\{z\}\{\\sqrt\{2\}\\ell\}\\right\)\+\\sqrt\{\\frac\{2\}\{\\pi\}\}\\ell\\exp\\left\(\-\\frac\{z^\{2\}\}\{2\\ell^\{2\}\}\\right\)\.Putting it together, the expression of\[𝐊yy\]ij\[\\mathbf\{K\}\_\{yy\}\]\_\{ij\}writes as
\[𝐊yy\]ij=ℓπ2\[H\(bi−aj\)−H\(ai−aj\)−H\(bi−bj\)\+H\(ai−bj\)\]\+σ2δi=j\.\[\\mathbf\{K\}\_\{yy\}\]\_\{ij\}=\\ell\\sqrt\{\\frac\{\\pi\}\{2\}\}\\Big\[H\(b\_\{i\}\-a\_\{j\}\)\-H\(a\_\{i\}\-a\_\{j\}\)\-H\(b\_\{i\}\-b\_\{j\}\)\+H\(a\_\{i\}\-b\_\{j\}\)\\Big\]\+\\ \\sigma^\{2\}\\delta\_\{i=j\}\.
Figure 7:Extended results of[Figure3](https://arxiv.org/html/2605.05520#S3.F3)\. Comparison between the reconstructions of the baselines on the inverse problem with diffusion prior on the setting of GP\. The y\-axis shows intensities while the x\-axis represents a one\-dimensional grid of the\[−5,5\]\[\-5,5\]interval\. The dashed horizontal lines depict values of the integral over the observed intervals\.
## Appendix FCensored Gaussian Process
#### Setting revisted\.
We revisit the setting in[Section3\.1](https://arxiv.org/html/2605.05520#S3.SS1)\. We assume that the rain fieldX:Ω⊂ℝ2→ℝ\+X\\mathrel\{\\mathop\{\\ordinarycolon\}\}\\Omega\\subset\\mathbb\{R\}^\{2\}\\rightarrow\\mathbb\{R\}\_\{\+\}is described by a left\-censored GP with a power transformationβ\>1\\beta\>1
X\(𝐬\)=max\(0,V\(𝐬\)β\),𝐬∈Ω⊂ℝ2,X\(\\mathbf\{s\}\)=\\max\\left\(0,V\(\\mathbf\{s\}\)^\{\\beta\}\\right\),\\qquad\\mathbf\{s\}\\in\\Omega\\subset\\mathbb\{R\}^\{2\},namely,V∼𝒢𝒫\(μ,k\)V\\sim\\mathcal\{GP\}\(\\mu,k\)being a stationary GP with constantμ∈ℝ\\mu\\in\\mathbb\{R\}and kernel
k\(𝐬,𝐬′\)=σ2exp\(−12\(𝐬−𝐬′\)⊤𝐐−2\(𝐬−𝐬′\)\),σ\>0,𝐐=Diag\(ℓ1,ℓ2\)∈ℝ2×2,k\(\\mathbf\{s\},\\mathbf\{s\}^\{\\prime\}\)=\\sigma^\{2\}\\exp\\left\(\-\\frac\{1\}\{2\}\(\\mathbf\{s\}\-\\mathbf\{s\}^\{\\prime\}\)^\{\\top\}\\mathbf\{Q\}^\{\-2\}\\ \(\\mathbf\{s\}\-\\mathbf\{s\}^\{\\prime\}\)\\right\),\\quad\\sigma\>0,\\quad\\mathbf\{Q\}=\\mathrm\{Diag\}\\left\(\{\\ell\}\_\{1\},\{\\ell\}\_\{2\}\\right\)\\in\\mathbb\{R\}^\{2\\times 2\},where the positive definite matrix𝐐\\mathbf\{Q\}is meant to account for anisotropy of the field\. The domain of interestΩ⊂ℝ2\\Omega\\subset\\mathbb\{R\}^\{2\}is discretized in a grid ofH×WH\\times Wat the points\{sk\}k=1HW\\\{s\_\{k\}\\\}\_\{k=1\}^\{HW\}and we assume to have access to the realization ofNNof the discretized fieldxi∈ℝ\+HWx\_\{i\}\\in\\mathbb\{R\}^\{HW\}\_\{\+\}, that is, fori∈⟦1,N⟧i\\in\\llbracket 1,N\\rrbracket, we have
xi=max\(0,viβ\),wherevi∼𝒩\(μ𝟏,Γ\(𝐐\)\)andΓℓ,j=k\(𝐬ℓ,𝐬j\),for0≤ℓ,j≤HWx\_\{i\}=\\max\(0,v\_\{i\}^\{\\beta\}\),\\quad\\text\{ where \}\\quad v\_\{i\}\\sim\\mathcal\{N\}\\left\(\\mu\\mathbf\{1\},\\Gamma\(\\mathbf\{Q\}\)\\right\)\\quad\\text\{ and \}\\quad\\Gamma\_\{\\ell,j\}=k\(\\mathbf\{s\}\_\{\\ell\},\\mathbf\{s\}\_\{j\}\),\\;\\text\{ for \}0\\leq\\ell,j\\leq HWThe goal is to estimate the paramters of the latent model\(μ,𝐐,σ,β\)\(\\mu,\\mathbf\{Q\},\\sigma,\\beta\)\. To do so, we use the EM algorithm and proceed first by deriving the algorithm whenβ=1\\beta=1, and then generalize toβ\>1\\beta\>1\.
#### Parameters estimation of the model\.
One way to estimate the parameters and avoid optimizing jointlyβ\\beta, is to fit the model parameters for different values ofβ\\betaand then select the configuration that gives the best classifier two\-sample test\.
### F\.1EM algorithm forβ=1\\beta=1
The joint model factorizes as
p\(xi,vi\)=δmax\(0,vi\)\(xi\)pV\(vi\),p\(x\_\{i\},v\_\{i\}\)=\\delta\_\{\\max\(0,v\_\{i\}\)\}\(x\_\{i\}\)\\,p\_\{V\}\(v\_\{i\}\),wherepVp\_\{V\}denotes the marginal distribution ofviv\_\{i\}and is given by the Gaussian𝒩\(μ𝟏,Γ\(𝐐\)\)\\mathcal\{N\}\(\\mu\\mathbf\{1\},\\Gamma\(\\mathbf\{Q\}\)\)\. For a given observationxix\_\{i\}, we define the index sets corresponding to positive and zero components as
Ii\+=\{k:\[xi\]k\>0\},Ii0=\{k:\[xi\]k=0\},I\_\{i\}^\{\+\}=\\\{k\\mathrel\{\\mathop\{\\ordinarycolon\}\}\[x\_\{i\}\]\_\{k\}\>0\\\},\\qquad I\_\{i\}^\{0\}=\\\{k\\mathrel\{\\mathop\{\\ordinarycolon\}\}\[x\_\{i\}\]\_\{k\}=0\\\},and denote by\(xi\+,xi0\)\(x\_\{i\}^\{\+\},x\_\{i\}^\{0\}\)the restriction ofxix\_\{i\}toIi\+I\_\{i\}^\{\+\}andIi0I\_\{i\}^\{0\}, respectively\.
*E\-step\.*The E\-step consists in sampling from the posteriorpV\|X\(⋅∣xi\)p\_\{V\|X\}\(\\cdot\\mid x\_\{i\}\), which factorizes as
pV\|X\(vi∣xi\)∝δxi\+\(vi\+\)×𝟏\{vi0≤0\}pV\(vi0∣vi\+=xi\+\)\.p\_\{V\|X\}\(v\_\{i\}\\mid x\_\{i\}\)\\propto\\delta\_\{x\_\{i\}^\{\+\}\}\(v\_\{i\}^\{\+\}\)\\times\\mathbf\{1\}\_\{\\\{v\_\{i\}^\{0\}\\leq 0\\\}\}\\,p\_\{V\}\(v\_\{i\}^\{0\}\\mid v\_\{i\}^\{\+\}=x\_\{i\}^\{\+\}\)\.\(13\)The second factor corresponds to a truncated multivariate Gaussian distribution\. Samples are obtained via Gibbs sampling, leveraging the fact that one\-dimensional truncated Gaussian conditionals can be sampled exactly using inverse transform sampling\.
*M\-step\.*At iterationkk, posterior samplesvi∼pV\|X\(k\)\(⋅∣xi\)v\_\{i\}\\sim p\_\{V\|X\}^\{\(k\)\}\(\\cdot\\mid x\_\{i\}\)are drawn using the current parameter estimates\. The parameters\(μ,Q\)\(\\mu,Q\)are then updated by maximizing the complete\-data log\-likelihood,
argmaxμ,𝐐,σ∑i=1NlogpV\(vi\),\\arg\\max\_\{\\mu,\\mathbf\{Q\},\\sigma\}\\sum\_\{i=1\}^\{N\}\\log p\_\{V\}\(v\_\{i\}\),which corresponds to standard maximum likelihood estimation for a multivariate Gaussian model\.
### F\.2EM algorithm forβ\>1\\beta\>1
Forβ\>1\\beta\>1, the posterior distributionpV\|Xβ\(⋅∣xi\)p\_\{V\|X\}^\{\\beta\}\(\\cdot\\mid x\_\{i\}\)no longer involves a truncated multivariate Gaussian\. Nevertheless, its density admits a closed\-form expression:
pVβ\(vi0∣vi\+=xi\+\)=pV\(\(vi0\)1/β∣vi\+=xi\+\)\|det\(Jf\(vi0\)\)\|,p\_\{V\}^\{\\beta\}\(v\_\{i\}^\{0\}\\mid v\_\{i\}^\{\+\}=x\_\{i\}^\{\+\}\)=p\_\{V\}\\\!\\left\(\(v\_\{i\}^\{0\}\)^\{1/\\beta\}\\mid v\_\{i\}^\{\+\}=x\_\{i\}^\{\+\}\\right\)\\Big\|\\det\\\!\\left\(\\mathrm\{J\}\_\{f\}\(v\_\{i\}^\{0\}\)\\right\)\\Big\|,\(14\)whereJf\(vi0\)\\mathrm\{J\}\_\{f\}\(v\_\{i\}^\{0\}\)denotes the Jacobian of the elementwise power transformationf\(v\)=vβf\(v\)=v^\{\\beta\}\. The Jacobian is diagonal and can be computed exactly and efficiently\.
*E\-step\.*The E\-step is performed using a Metropolis\-within\-Gibbs scheme\. We use the same procedure in E\-step of the caseβ=1\\beta=1as a proposal distribution, and acceptance ratios are computed using the density in \([14](https://arxiv.org/html/2605.05520#A6.E14)\)\.
*M\-step\.*Although the power transformation modifies the posterior density, it does not affect the model parameters\. Consequently, the M\-step mirrors that of theβ=1\\beta=1case, with the exception that posterior samples are transformed back by applying the inverse power\. Specifically, parameter updates are obtained by maximizing
argmaxμ,𝐐,σ∑i=1NlogpV\(vi1/β\)\.\\arg\\max\_\{\\mu,\\mathbf\{Q\},\\sigma\}\\sum\_\{i=1\}^\{N\}\\log p\_\{V\}\\\!\\left\(v\_\{i\}^\{1/\\beta\}\\right\)\.which, similarly, correspond maximum likelihood estimation for a multivariate Gaussian model expect the inputs arevi1/βv\_\{i\}^\{1/\\beta\}\.
### F\.3Practical Computational Challenges
While the procedure described above enables learning from censored processes, it scales poorly in our current high\-dimensional setting, where rain fields have dimension48×36=172848\\times 36=1728and nearly50%50\\%of the entries are censored\. The E\-step requires posterior sampling to impute the censored coordinates of each rain field by sampling from a truncated multivariate Gaussian distribution\. Since the censored coordinates vary across rain fields, this limits the ability to leverage vectorized computation and makes estimation prohibitively dependent on the number of rain fields, which is21,19621\{,\}196in our setting\. In addition, using Gibbs sampling for the truncated multivariate Gaussian requires a long burn\-in period due to the dimensionality of the data\.
## Appendix GExamples of Reconstructions in CML Experiments
We present qualitative examples of rain field reconstructions from CML measurements on the OpenMRG dataset\. The comparison includes all considered baselines, both meteorological and diffusion\-based methods\. Three representative reconstruction examples are shown\.































Figure 8:Qualitative comparisons of rain field reconstructions on real CMLs links from OpenMRG dataset on three reconstruction tasks\. The network of CMLs is depicted in red\.






























Figure 9:Another set of qualitative comparisons of rain field reconstructions on real CMLs links from OpenMRG dataset on three reconstruction tasks\. The network of CMLs is depicted in red\.Similar Articles
Conditional Diffusion Under Linear Constraints: Langevin Mixing and Information-Theoretic Guarantees
This paper analyzes zero-shot conditional sampling with pretrained diffusion models for linear inverse problems, providing information-theoretic guarantees and proposing a projected-Langevin initialization method.
Elucidating the SNR-t Bias of Diffusion Probabilistic Models
This paper identifies a Signal-to-Noise Ratio timestep (SNR-t) bias in diffusion probabilistic models during inference, where SNR-timestep alignment from training is disrupted at inference time. The authors propose a differential correction method that decomposes samples into frequency components and corrects each separately, improving generation quality across models like IDDPM, ADM, DDIM, EDM, and FLUX with minimal computational overhead.
Christoffel-DPS: Optimal sensor placement in diffusion posterior sampling for arbitrary distributions
This paper introduces Christoffel-DPS, a distribution-free framework for optimal sensor placement in diffusion posterior sampling that outperforms classical Gaussian-based methods. It provides theoretical guarantees and practical improvements for reconstructing states from complex, non-Gaussian distributions using generative models.
Multi-Quantile Regression for Extreme Precipitation Downscaling
Introduces Q-srdrn, a multi-quantile super-resolution network using pinball loss to improve extreme precipitation downscaling, achieving dramatic detection rate gains for heavy rainfall events while maintaining overall accuracy.
Modeling Sparse and Bursty Vulnerability Sightings: Forecasting Under Data Constraints
Academic study compares SARIMAX and Poisson regression for forecasting sparse, bursty vulnerability-sighting time-series, finding count-based models more stable.