Data Dependent Noise Models¶
Important: Please read the installation page for details about how to install the toolboxes. $\newcommand{\dotp}[2]{\langle #1, #2 \rangle}$ $\newcommand{\enscond}[2]{\lbrace #1, #2 \rbrace}$ $\newcommand{\pd}[2]{ \frac{ \partial #1}{\partial #2} }$ $\newcommand{\umin}[1]{\underset{#1}{\min}\;}$ $\newcommand{\umax}[1]{\underset{#1}{\max}\;}$ $\newcommand{\umin}[1]{\underset{#1}{\min}\;}$ $\newcommand{\uargmin}[1]{\underset{#1}{argmin}\;}$ $\newcommand{\norm}[1]{\|#1\|}$ $\newcommand{\abs}[1]{\left|#1\right|}$ $\newcommand{\choice}[1]{ \left\{ \begin{array}{l} #1 \end{array} \right. }$ $\newcommand{\pa}[1]{\left(#1\right)}$ $\newcommand{\diag}[1]{{diag}\left( #1 \right)}$ $\newcommand{\qandq}{\quad\text{and}\quad}$ $\newcommand{\qwhereq}{\quad\text{where}\quad}$ $\newcommand{\qifq}{ \quad \text{if} \quad }$ $\newcommand{\qarrq}{ \quad \Longrightarrow \quad }$ $\newcommand{\ZZ}{\mathbb{Z}}$ $\newcommand{\CC}{\mathbb{C}}$ $\newcommand{\RR}{\mathbb{R}}$ $\newcommand{\EE}{\mathbb{E}}$ $\newcommand{\Zz}{\mathcal{Z}}$ $\newcommand{\Ww}{\mathcal{W}}$ $\newcommand{\Vv}{\mathcal{V}}$ $\newcommand{\Nn}{\mathcal{N}}$ $\newcommand{\NN}{\mathcal{N}}$ $\newcommand{\Hh}{\mathcal{H}}$ $\newcommand{\Bb}{\mathcal{B}}$ $\newcommand{\Ee}{\mathcal{E}}$ $\newcommand{\Cc}{\mathcal{C}}$ $\newcommand{\Gg}{\mathcal{G}}$ $\newcommand{\Ss}{\mathcal{S}}$ $\newcommand{\Pp}{\mathcal{P}}$ $\newcommand{\Ff}{\mathcal{F}}$ $\newcommand{\Xx}{\mathcal{X}}$ $\newcommand{\Mm}{\mathcal{M}}$ $\newcommand{\Ii}{\mathcal{I}}$ $\newcommand{\Dd}{\mathcal{D}}$ $\newcommand{\Ll}{\mathcal{L}}$ $\newcommand{\Tt}{\mathcal{T}}$ $\newcommand{\si}{\sigma}$ $\newcommand{\al}{\alpha}$ $\newcommand{\la}{\lambda}$ $\newcommand{\ga}{\gamma}$ $\newcommand{\Ga}{\Gamma}$ $\newcommand{\La}{\Lambda}$ $\newcommand{\si}{\sigma}$ $\newcommand{\Si}{\Sigma}$ $\newcommand{\be}{\beta}$ $\newcommand{\de}{\delta}$ $\newcommand{\De}{\Delta}$ $\newcommand{\phi}{\varphi}$ $\newcommand{\th}{\theta}$ $\newcommand{\om}{\omega}$ $\newcommand{\Om}{\Omega}$
The simplest noise model is Gaussian additive noise, where the variance of the pixel value is independent of the mean (which is the value we look for).
Unfortunately, most real-life data corresponds to noise model that are much more complicated, and in particular the variance of the noise often depends on the parameter of interest (for instance the mean).
warning off
addpath('toolbox_signal')
addpath('toolbox_general')
addpath('solutions/denoisingwav_5_data_dependent')
warning on
Poisson Noise¶
A Poisson model assume that each pixel $x$ of an image $f(x)$ is drawn from a Poisson distribution of parameter $\lambda=f_0(x)$, where $f_0$ is the clean intensity image to recover.
$$ \PP(f(x)=k)=\lambda^k e^{-\lambda}/k! $$Display the Poisson distribution for several value of $\lambda$.
lambda = [4 10 20];
[k,Lambda] = meshgrid(1:50, lambda);
P = Lambda.^k .* exp(-Lambda)./factorial(k);
h = plot(P'); axis('tight');
if using_matlab()
set(h, 'LineWidth', 2);
end
legend('\lambda=2', '\lambda=10', '\lambda=20');
set_label('k', 'P(k)');
This model corresponds to a photon count, where $\lambda$ is proportional to the number of photons that hits the receptor during the exposition time. This is useful to model medical imaging (confocal microscopy), TEP and SPECT tomography, and digital camera noises.
The goal of denoising is to retrieve the value of $f_0(x)=\lambda(x)$, the true image value, which corresponds to a clean image.
Note that $\lambda$ is the mean of the distribution, but is also its variance, so that the noise intensity perturbating the image pixel $f(x)$ is proportional to $f_0(x)$.
We load a clean, unquantized image.
n = 256;
name = 'lena';
f0u = rescale( load_image(name,n) );
Quantize the values to given $\lambda$ range.
lmin = 1;
lmax = 40;
f0 = floor( rescale(f0u,lmin,lmax) );
Generate the noisy image.
f = poissrnd(f0);
Display.
clf;
imageplot(f0, 'Intensity map f0', 1,2,1);
imageplot(f, 'Observed noisy image f', 1,2,2);
Display the difference, which shows that the noise level depends on the intensity. The noise is larger in bright areas.
clf;
imageplot(f0, 'Intensity map f0', 1,2,1);
imageplot(f-f0, 'f-f0', 1,2,2);
Exercise 1
Display noisy image contaminated by Poisson noise of varying range.
exo1()