As in image enhancement the goal of restoration is to improve an image for further processing. In contrast to image enhancement that was subjective and largely based on heuristics, restoration attempts to reconstruct or recover an image that has been distorted by a known degradation phenomenon. Restoration techniques thus try to model the degradation process and apply the inverse process in order to reconstruct the original image.
The degradation process is generally modelled as a degradation function and the additive noise term
together they yield
.
Given , some knowledge about
, and
it is the objective of restoration to estimate
. Of course this estimate should be as close as possible to
.
![]() Fig 6.1 Image degradation model |
If the degradation function is a linear, shift-invariant process, then the degraded image is given in the spatial domain by
![]() | (6.1) |
where is the spatial representation of the degradation function, and
indicates convolution.
As we know from the Convolution Theorem this can be rewritten in the frequency domain
![]() | (6.2) |
where the capital letters are the Fourier transforms of the respective function in Eq 6.1.
The principle sources of noise in digital images arise during image acquisition and/or transmission.
The performance of imaging sensors are affected by a variety of factors during acquisition, such as
Images can also be corrupted during transmission due to interference in the channel for example
Depending on the specific noise source, a different model must be selected that accurately reproduces the spatial characteristics of the noise.
Because of its mathematical tractability in both the spatial and frequency domain, Gaussian noise models (aka normal distribution) are used frequently in practice. The PDF (Probability density function) of a Gaussian random variable
where |
![]() Fig 6.2: Gaussian probability density function ![]() |
The application of the Gaussian model is so convenient that it is often used in situations in which they are marginally applicable at best.
The PDF of Rayleigh noise is defined by
the mean and variance are given by
|
![]() Fig 6.3: Rayleigh probability density function ![]() |
As the shape of the Rayleigh density function is skewed it is useful for approximating skewed histograms.
The PDF of the Erlang noise is given by
where
|
![]() Fig 6.4: Erlang probability density function ![]() |
The PDF of the exponential noise is given by
where
|
![]() Fig 6.5: Probability density function ![]() |
The exponential noise model is a special case of the Erlang noise model with .
The PDF of the uniform noise is given by
The mean and variance of the density function are given by
|
![]() Fig 6.6: Probability density function ![]() |
The PDF of bipolar impulse noise model is given by
|
![]() Fig 6.7: Probability density function ![]() |
If , grey-level
appears as a light dot (salt) in the image. Conversely,
will appear as dark dot (pepper). If either
is zero, the PDF is called unipolar.
Because impulse corruption is generally large compared to the signal strength, the assumption is usually that and
are digitised as saturated values thus black (pepper) and white (salt).
The presented noise models provide a useful tool for approximating a broad range of noise corruption situations found in practice. For example:
Although the corresponding histograms show close resemblance to the PDF of the respective noise model, it is virtually impossible to select the type of noise causing the degradation from the image alone.
|
![]() |
![]() |
|
![]() |
![]() |
|
![]() |
![]() |
|
![]() |
![]() |
|
![]() |
![]() |
(Salt & Pepper) |
![]() |
![]() |
Periodic noise typically arises from electrical or electromechanical interference during image acquisition and is spatial dependent. |
|
![]() Fig 6.8: Periodic noise on a satellite image of Pompeii |
![]() Fig 6.9: Spectrum of the Pompeii image with periodic noise (peaks of the periodic noise visually enhanced) |
In rare cases the parameters of noise PDFs may be known from sensor specifications, but it is often necessary to estimate them for a particular imaging arrangement. If the imaging system is available, one simple way to study noise characteristics would be to capture multiple images of homogeneous environments, such as solid grey cards. |
![]() Fig 6.10: Sample histogram of a small patch of Erlang noise |
In many cases the noise parameters have to be estimated from already captured images. The simplest way is to estimate the mean and variance of the grey levels from small patches of reasonably constant grey level.
The parameters of periodic noise are typically estimated by inspection of the Fourier spectrum. The frequency spikes can often be detected by visual analysis. Automatic analysis is possible if the noise spikes are either exceptionally pronounced or a priori knowledge about their location is known. |
![]() Fig. 6.11: |
![]() Fig. 6.12: |
When noise is the only image degradation present in an image, thus , Eqs 6.1 & 6.2 become
![]() | (6.13) |
and
![]() | (6.14) |
In case of periodic noise, it is usually possible to estimate from the spectrum of
and subtract it to obtain an estimate of the original image
.
But estimating the noise terms is unreasonable, so subtracting them from
is impossible. Spatial filtering is the method of choice in situations when only additive noise is present.
Bandreject filters remove or attenuate a band of frequencies around the origin in the Fourier domain.
|
![]() Fig 6.13: Ideal bandreject filter |
||
|
![]() Fig 6.14: Butterworth bandreject filter |
||
|
![]() Fig 6.15: Gaussian bandreject filter |
![]() Fig 6.16: Aerial image of Pompeii with periodic noise |
![]() Fig 6.17: Ideal bandreject filter ![]() |
![]() Fig 6.18: Butterworth bandreject filter ![]() |
![]() Fig 6.19: Gaussian bandreject filter ![]() |
The Bandpass filter performs the opposite of a bandreject filter and thus lets pass the frequencies in a narrow band of width around
. The transfer function
can be obtained from a corresponding bandreject filter by using the equation
![]() | (6.18) |
The perspective plots of the corresponding filters are depicted in Fig 6.20.
![]() (a) |
![]() (b) |
![]() (c) |
Fig 6.20 (a) Ideal bandpass filter, (b) Butterworth bandpass filter, and (c) Gaussian bandpass filter |
![]() Fig 6.21: Aerial image of Pompeii with periodic noise |
![]() Fig 6.22: Ideal bandpass filter ![]() |
![]() Fig 6.23: Butterworth bandpass filter ![]() |
![]() Fig 6.24: Gaussian bandpass filter ![]() |
A Notch filter rejects or passes frequencies in predefined neighbourhoods about a centre frequency. Due to the symmetry of the Fourier transform, notch filters must appear in symmetric pairs about the origin (except the one at the origin).
The transfer functions of the ideal, Butterworth, and Gaussian notch reject filter of radius (width) with centres at
and, by symmetry, at
, are
|
![]() Fig 6.25: Ideal notch filter |
||
|
![]() Fig 6.26: Butterworth notch filter |
||
|
![]() Fig 6.27: Gaussian notch filter |
where
![]() | (6.22) |
Similar to the bandpass/bandreject filters, the Notch reject filters can be turned into Notch pass filters with the relation
![]() | (6.23) |
where is the transfer function of the notch pass filter, and
the corresponding Notch reject filter.
If![]()
- the Notch reject filter becomes a highpass filter and
- the Notch pass filter becomes a lowpass filter.
![]() Fig 6.28: Pompeii aerial image with periodic noise |
|
![]() Fig 6.29: Notch reject filtered image |
![]() Fig 6.30: Notch pass filtered image (contrast enhanced) |
![]() Fig 6.31: Pompeii aerial image with periodic noise |
|
![]() Fig 6.32: Notch reject filtered image |
![]() Fig 6.33: Notch pass filtered image (contrast enhanced) |
![]() Fig 6.34: Pompeii aerial image with periodic noise |
|
![]() Fig 6.35: Notch reject filtered image |
![]() Fig 6.36: Notch pass filtered image (contrast enhanced) |
![]() Fig 6.37: Original Mariner 6 martian image |
![]() Fig 6.38: Log Fourier spectra of the image |
![]() Fig 6.39: Notch filtered log spectra |
![]() Fig 6.40: Notch filtered image |
Images derived from electro-optical scanners, such as those used in space and aerial imaging, are sometimes corrupted by coupling and amplification of low-level signals in the scanners' circuitry. The resulting images tend to contain substantial periodic noise. The interference pattern are, however, in practice not always as clearly defined as in the samples shown thus far.
The first step in Optimum notch filtering is to find the principal frequency components and placing notch pass filters at the location of each spike, yielding . The Fourier transform of the interference pattern is thus given by
![]() | (6.24) |
where is the Fourier transform of the corrupted image. The corresponding interference pattern in the spatial domain is obtained with the inverse Fourier transform
![]() | (6.25) |
As the corrupted image is assumed to be formed by the addition of the uncorrupted image
and the interference noise
![]() | (6.26) |
If the noise term were known completely, subtracting it from
would yield
. The problem is of course that the estimated interference pattern
is only an approximation. The effect of components not present in the estimate
can be minimised be subtracting from
a weighted portion of
to obtain an estimate for
![]() | (6.27) |
where is the estimate of
, and
is to be determined and optimises
in some meaningful way.
One common approach is to select so that the variance of
is minimised over a specified neighbourhood of size
about every point
.
The additional assumption that remains constant over the neighbourhood yields the approximation
![]() | (6.28) |
and we can solve the minimisation problem with
![]() | (6.29) |
where is the average value of
in the neighbourhood; that is
![]() | (6.30) |
Substituting Eq 6.27 into Eq 6.29 yields
![]() | (6.31) |
Assuming that remains constant over the neighbourhood gives the approximation
![]() | (6.32) |
for and
. This assumption results in the expression
![]() | (6.33) |
in the neighbourhood. With these approximations Eq 6.30 becomes
![]() | (6.34) |
To minimise , we solve
![]() | (6.35) |
for . The result is
![]() | (6.36) |
To obtain an estimated restored image ,
must be computed from Eq 6.36 and substituted into Eq 6.27 to determine
.
![]() Fig 6.41 Principle of the optimum notch filter |
![]() Fig 6.42 Aerial image of Pompeii with heavy periodic noise distortions |
![]() Fig 6.43 Noise estimation with Gaussian notch pass filter |
![]() Fig 6.44 Notch reject filtered image |
![]() Fig 6.45 Result filtered with an optimum notch filter of size ![]() |
![]() ![]()
![]() |
Fig 6.46 Aerial image of Pompeii corrected with the optimal notch filter of size with the window shifted pixelwise (mouseover) shows the result when shifting the window by 50
![]() |
Fig 6.47 of the example image in Fig 6.46
So far we concentrated mainly at the additive noise term . In the remainder of this presentation we concentrate on ways to remove the degradation function
.
![]() Fig 6.48 Image degradation model |
There are three principle methods to estimate the degradation function
The process of restoring a corrupted image using the estimated degradation function is sometimes called blind deconvolution, as the true degradation is seldom known completely.
In recent years a forth method based on maximum-likelihood estimation (MLE) often called blind deconvolution has been proposed.
Beware: All multiplications/divisions in the equations in this part of the lecture are componentwise.
A degraded image without any information about the degradation function .
Equipment similar to the equipment used to acquire the degraded image
Degradation modelling has been used for many years. It, however, requires an indepth understanding of the involved physical phenomenon. In some cases it can even incorporate environmental conditions that cause degradations. Hufnagel and Stanly proposed in 1964 a model to characterise atmospheric turbulence.
![]() | (6.37) |
where is a constant that depends on the nature of the turbulence.
Except of the 5/6 power the above equation has the same form as a Gaussian lowpass filter.
![]() Fig 6.49: Aerial image of Pompeii |
![]() Fig 6.50: Pompeii with mild turbulence ![]() |
![]() Fig 6.51: Pompeii with medium turbulence ![]() |
![]() Fig 6.52 Pompeii with heavy turbulence ![]() |
The simplest approach to restore a degraded image is to form an estimate of the form
![]() | (6.38) |
and then obtain the corresponding image by inverse Fourier transform. This method is called direct inverse filtering. With the degradation model Equation 6.2 we get
![]() | (6.39) |
This simple equation tells us that, even if we knew exactly we could not recover the undegraded image
completely because:
Direct Inverse filtering is seldom a suitable approach in practical applications
![]() | (6.40) |
where is a low-pass filter that should eliminiate the very low (or even zero) values often experienced in the high frequencies.
![]() | (6.41) |
gives values, where is lower than a threshold
a special treatment.
These appraoches are also known as Pseudo-Inverse Filtering.
![]() Fig 6.53: Checkerboard sample |
![]() Fig 6.54: Motion blur of ![]() ![]() |
![]() Fig 6.55: Additive noise with ![]() |
![]() Fig 6.56: Motion and noise degraded checkerboard |
Although very tempting from its simplicity, the direct inverse filtering approach does not work for practical applications. Even if the degradation is known exactly.
![]() Fig 6.57: Direct inverse filtering result |
![]() Fig 6.58: Inverse filtering of Figure 6.52 |
![]() Fig 6.59: Pseudo Inverse filtering for freq. radius ![]() |
![]() Fig 6.60: Pseudo Inverse filtering for frequencies radius ![]() |
![]() Fig 6.61: Pseudo Inverse filtering for frequencies radius ![]() |
The best known improvement to inverse filtering is the Wiener Filter. The Wiener filter seeks an estimate that minimises the statistical error function
![]() | (6.42) |
where is the expected value and
the undegraded image. The solution to this equation in the frequency domain is
![]() | (6.43) |
where
If the noise spectrum is zero, the noise-to-signal power ratio
vanishes and the Wiener Filter reduces to the inverse filter.
![]() | (6.44) |
This is no problem, as the inverse filter works fine if no noise is present.
The main problem with the Wiener filter is that the power spectrum of the undegraded image is seldom known.
A frequently used approach when these quantities cannot be estimated is to approximate Eq 6.43 by the expression
![]() | (6.45) |
where is a user selected constant.
This simplification can be partly justified when dealing with spectrally white noise, where the noise spectrum is constant. However, the problem still remains that the power spectrum
of the undegraded image is unknown and must be estimated.
Even if the actual ratiois not known, it becomes a simple matter to experiment interactively varying the constant and viewing the results.
This example shows the performance of the Wiener filter on the example Fig 6.56 using (1) the approximation of Eq 6.45 with in Fig 6.62 and (2) the Wiener filter using full knowledge of the noise and undegraded image's power spectra in Fig 6.63.
![]() Fig 6.62: Parametric Wiener filter using a constant ratio |
![]() Fig 6.63: Wiener filter knowing the noise and undegraded signal power spectra |
![]() Fig 6.64: Pompeii blurred with Gaussian kernel ![]() ![]() |
![]() Fig 6.65: Result of direct inverse filtering |
![]() Fig 6.66: Result with Parametric Wiener filtering ![]() |
![]() Fig 6.67: Result with Wiener filtering |
![]() Fig 6.68: Pompeii blurred with Gaussian kernel ![]() ![]() |
![]() Fig 6.69: Result of direct inverse filtering |
![]() Fig 6.70: Result with Parametric Wiener filtering ![]() |
![]() Fig 6.71: Result with Wiener filtering |
![]() Fig 6.72: Pompeii blurred with Gaussian kernel ![]() ![]() |
![]() Fig 6.73: Result of direct inverse filtering |
![]() Fig 6.74: Result with Parametric Wiener filtering ![]() |
![]() Fig 6.75: Result with Wiener filtering |
![]() ![]()
|
![]() (a) |
![]() (b) |
Fig 6.76 Aerial image of Pompeii with Gaussian blur () and noise (
) (a) restored with the Parametric Wiener (
), (b) (mouseover) same image restored with Wiener filter.
The problem of having to know something about the degradation function is common to all methods discussed in this part.
The Wiener filter has the additional difficulty, that
must be known.
Although an approximation by a constant is possible as in the Parametric Wiener filter, see Fig 6.62, this approach is not always suitable.
The Geometric Mean Filter is a generalisation of inverse filtering and the Wiener filtering. Through its parameters it allows access to an entire family of filters.
![]() | (6.46) |
Depending on the parameters the filter characteristics can be adjusted
The Constrained Least Squares Filtering approach only requires knowledge of the mean and variance of the noise. As shown in previous lectures these parameters can be usually estimated from the degraded image.
The definition of the 2D discrete convolution is
![]() | (6.47) |
Using this equation we can express the linear degradation model in vector notation as
![]() | (6.48) |
It seems obvious that the restoration problem is now reduced to simple matrix manipulations. Unfortunately this is not the case. The problem with the matrix formulation is that the matrix is very large
and that its inverse is firstly very sensitive to noise and secondly does not necessarily exist.
One way to deal with these issues is to base optimality of restoration on a measure of smoothness, such as the second derivative of the image, e.g. the Laplacian.
![]() | (6.49) |
subject to the constraint
![]() | (6.50) |
where is the Euclidian norm, and
an estimate of the undegraded image.
The Frequency domain solution to this optimisation problem is given by
![]() | (6.51) |
where is a parameter that must be adjusted so that the constraint in Eq 6.50 is satisfied, and
is the Fourier transform of the Laplace operator
![]() | (6.52) |
Note: must be properly padded with zeros prior to computing the Fourier transform.
It is possible to adjust the parameterinteractively until acceptable results are obtained.
If we are interested in optimality, the parameter must be adjusted so that the constraint
in Eq 6.50 is satisfied.
It turns out that
![]() | (6.53) |
is the optimal selection given the above optimisation equation where are the mean and variance of the noise.
It is important to understand, that optimum restoration in the sense of constrained least squares does not necessarily imply best in the visual sense. In general, automatic determined restoration filters yield inferior results to manual adjustment of filter parameters.
This example shows the performance of the Constrained Least Squares Filter on the example Fig 6.56 using (1) the theoretical optimal value for in Fig 6.77 and (2) the constrained least squares filter with a manually selected
of
6.78.
![]() Fig 6.77: Result with optimal ![]() |
![]() Fig 6.78: Results with manually selected ![]() |
![]() Fig 6.79: Pompeii blurred with Gaussian kernel ![]() |
![]() Fig 6.80: Result with low Gaussian noise of ![]() |
![]() Fig 6.81: Result with medium Gaussian noise of ![]() |
![]() Fig 6.82: Result with high Gaussian noise of ![]() |