restoration
¶Image restoration module.
skimage.restoration.wiener (image, psf, balance) 
WienerHunt deconvolution 
skimage.restoration.unsupervised_wiener (…) 
Unsupervised WienerHunt deconvolution. 
skimage.restoration.richardson_lucy (image, psf) 
RichardsonLucy deconvolution. 
skimage.restoration.unwrap_phase (image[, …]) 
Recover the original from a wrapped phase image. 
skimage.restoration.denoise_tv_bregman (…) 
Perform totalvariation denoising using splitBregman optimization. 
skimage.restoration.denoise_tv_chambolle (image) 
Perform totalvariation denoising on ndimensional images. 
skimage.restoration.denoise_bilateral (image) 
Denoise image using bilateral filter. 
skimage.restoration.denoise_wavelet (image[, …]) 
Perform wavelet denoising on an image. 
skimage.restoration.denoise_nl_means (image) 
Perform nonlocal means denoising on 2D or 3D grayscale images, and 2D RGB images. 
skimage.restoration.estimate_sigma (image[, …]) 
Robust waveletbased estimator of the (Gaussian) noise standard deviation. 
skimage.restoration.inpaint_biharmonic (…) 
Inpaint masked points in image with biharmonic equations. 
skimage.restoration.cycle_spin (x, func, …) 
Cycle spinning (repeatedly apply func to shifted versions of x). 
skimage.restoration.
wiener
(image, psf, balance, reg=None, is_real=True, clip=True)[source]¶WienerHunt deconvolution
Return the deconvolution with a WienerHunt approach (i.e. with Fourier diagonalisation).
Parameters:  image : (M, N) ndarray
psf : ndarray
balance : float
reg : ndarray, optional
is_real : boolean, optional
clip : boolean, optional


Returns:  im_deconv : (M, N) ndarray

Notes
This function applies the Wiener filter to a noisy and degraded image by an impulse response (or PSF). If the data model is
where \(n\) is noise, \(H\) the PSF and \(x\) the unknown original image, the Wiener filter is
where \(F\) and \(F^\dagger\) are the Fourier and inverse Fourier transfroms respectively, \(\Lambda_H\) the transfer function (or the Fourier transfrom of the PSF, see [Hunt] below) and \(\Lambda_D\) the filter to penalize the restored image frequencies (Laplacian by default, that is penalization of high frequency). The parameter \(\lambda\) tunes the balance between the data (that tends to increase high frequency, even those coming from noise), and the regularization.
These methods are then specific to a prior model. Consequently, the application or the true image nature must corresponds to the prior model. By default, the prior model (Laplacian) introduce image smoothness or pixel correlation. It can also be interpreted as highfrequency penalization to compensate the instability of the solution with respect to the data (sometimes called noise amplification or “explosive” solution).
Finally, the use of Fourier space implies a circulant property of \(H\), see [Hunt].
References
[R406]  François Orieux, JeanFrançois Giovannelli, and Thomas Rodet, “Bayesian estimation of regularization and point spread function parameters for WienerHunt deconvolution”, J. Opt. Soc. Am. A 27, 15931607 (2010) http://www.opticsinfobase.org/josaa/abstract.cfm?URI=josaa2771593 
[R407]  B. R. Hunt “A matrix theory proof of the discrete convolution theorem”, IEEE Trans. on Audio and Electroacoustics, vol. au19, no. 4, pp. 285288, dec. 1971 
Examples
>>> from skimage import color, data, restoration
>>> img = color.rgb2gray(data.astronaut())
>>> from scipy.signal import convolve2d
>>> psf = np.ones((5, 5)) / 25
>>> img = convolve2d(img, psf, 'same')
>>> img += 0.1 * img.std() * np.random.standard_normal(img.shape)
>>> deconvolved_img = restoration.wiener(img, psf, 1100)
skimage.restoration.
unsupervised_wiener
(image, psf, reg=None, user_params=None, is_real=True, clip=True)[source]¶Unsupervised WienerHunt deconvolution.
Return the deconvolution with a WienerHunt approach, where the
hyperparameters are automatically estimated. The algorithm is a
stochastic iterative process (Gibbs sampler) described in the
reference below. See also wiener
function.
Parameters:  image : (M, N) ndarray
psf : ndarray
reg : ndarray, optional
user_params : dict
clip : boolean, optional


Returns:  x_postmean : (M, N) ndarray
chains : dict

Other Parameters:  
The keys of ``user_params`` are: threshold : float
burnin : int
min_iter : int
max_iter : int
callback : callable (None by default)

Notes
The estimated image is design as the posterior mean of a probability law (from a Bayesian analysis). The mean is defined as a sum over all the possible images weighted by their respective probability. Given the size of the problem, the exact sum is not tractable. This algorithm use of MCMC to draw image under the posterior law. The practical idea is to only draw highly probable images since they have the biggest contribution to the mean. At the opposite, the less probable images are drawn less often since their contribution is low. Finally the empirical mean of these samples give us an estimation of the mean, and an exact computation with an infinite sample set.
References
[R408]  François Orieux, JeanFrançois Giovannelli, and Thomas Rodet, “Bayesian estimation of regularization and point spread function parameters for WienerHunt deconvolution”, J. Opt. Soc. Am. A 27, 15931607 (2010) http://www.opticsinfobase.org/josaa/abstract.cfm?URI=josaa2771593 
Examples
>>> from skimage import color, data, restoration
>>> img = color.rgb2gray(data.astronaut())
>>> from scipy.signal import convolve2d
>>> psf = np.ones((5, 5)) / 25
>>> img = convolve2d(img, psf, 'same')
>>> img += 0.1 * img.std() * np.random.standard_normal(img.shape)
>>> deconvolved_img = restoration.unsupervised_wiener(img, psf)
skimage.restoration.
richardson_lucy
(image, psf, iterations=50, clip=True)[source]¶RichardsonLucy deconvolution.
Parameters:  image : ndarray
psf : ndarray
iterations : int
clip : boolean, optional


Returns:  im_deconv : ndarray

References
[R409]  http://en.wikipedia.org/wiki/Richardson%E2%80%93Lucy_deconvolution 
Examples
>>> from skimage import color, data, restoration
>>> camera = color.rgb2gray(data.camera())
>>> from scipy.signal import convolve2d
>>> psf = np.ones((5, 5)) / 25
>>> camera = convolve2d(camera, psf, 'same')
>>> camera += 0.1 * camera.std() * np.random.standard_normal(camera.shape)
>>> deconvolved = restoration.richardson_lucy(camera, psf, 5)
skimage.restoration.
unwrap_phase
(image, wrap_around=False, seed=None)[source]¶Recover the original from a wrapped phase image.
From an image wrapped to lie in the interval [pi, pi), recover the original, unwrapped image.
Parameters:  image : 1D, 2D or 3D ndarray of floats, optionally a masked array
wrap_around : bool or sequence of bool, optional
seed : int, optional


Returns:  image_unwrapped : array_like, double

Raises:  ValueError

References
[R410]  Miguel Arevallilo Herraez, David R. Burton, Michael J. Lalor, and Munther A. Gdeisat, “Fast twodimensional phaseunwrapping algorithm based on sorting by reliability following a noncontinuous path”, Journal Applied Optics, Vol. 41, No. 35 (2002) 7437, 
[R411]  AbdulRahman, H., Gdeisat, M., Burton, D., & Lalor, M., “Fast threedimensional phaseunwrapping algorithm based on sorting by reliability following a noncontinuous path. In W. Osten, C. Gorecki, & E. L. Novak (Eds.), Optical Metrology (2005) 32–40, International Society for Optics and Photonics. 
Examples
>>> c0, c1 = np.ogrid[1:1:128j, 1:1:128j]
>>> image = 12 * np.pi * np.exp((c0**2 + c1**2))
>>> image_wrapped = np.angle(np.exp(1j * image))
>>> image_unwrapped = unwrap_phase(image_wrapped)
>>> np.std(image_unwrapped  image) < 1e6 # A constant offset is normal
True
skimage.restoration.
denoise_tv_bregman
(image, weight, max_iter=100, eps=0.001, isotropic=True)[source]¶Perform totalvariation denoising using splitBregman optimization.
Totalvariation denoising (also know as totalvariation regularization) tries to find an image with less totalvariation under the constraint of being similar to the input image, which is controlled by the regularization parameter ([R412], [R413], [R414], [R415]).
Parameters:  image : ndarray
weight : float
eps : float, optional
max_iter : int, optional
isotropic : boolean, optional


Returns:  u : ndarray

References
[R412]  (1, 2) http://en.wikipedia.org/wiki/Total_variation_denoising 
[R413]  (1, 2) Tom Goldstein and Stanley Osher, “The Split Bregman Method For L1 Regularized Problems”, ftp://ftp.math.ucla.edu/pub/camreport/cam0829.pdf 
[R414]  (1, 2) Pascal Getreuer, “Rudin–Osher–Fatemi Total Variation Denoising using Split Bregman” in Image Processing On Line on 2012–05–19, http://www.ipol.im/pub/art/2012/gtvd/article_lr.pdf 
[R415]  (1, 2) http://www.math.ucsb.edu/~cgarcia/UGProjects/BregmanAlgorithms_JacquelineBush.pdf 
skimage.restoration.
denoise_tv_chambolle
(image, weight=0.1, eps=0.0002, n_iter_max=200, multichannel=False)[source]¶Perform totalvariation denoising on ndimensional images.
Parameters:  image : ndarray of ints, uints or floats
weight : float, optional
eps : float, optional
n_iter_max : int, optional
multichannel : bool, optional


Returns:  out : ndarray

Notes
Make sure to set the multichannel parameter appropriately for color images.
The principle of total variation denoising is explained in http://en.wikipedia.org/wiki/Total_variation_denoising
The principle of total variation denoising is to minimize the total variation of the image, which can be roughly described as the integral of the norm of the image gradient. Total variation denoising tends to produce “cartoonlike” images, that is, piecewiseconstant images.
This code is an implementation of the algorithm of Rudin, Fatemi and Osher that was proposed by Chambolle in [R416].
References
[R416]  (1, 2) A. Chambolle, An algorithm for total variation minimization and applications, Journal of Mathematical Imaging and Vision, Springer, 2004, 20, 8997. 
Examples
2D example on astronaut image:
>>> from skimage import color, data
>>> img = color.rgb2gray(data.astronaut())[:50, :50]
>>> img += 0.5 * img.std() * np.random.randn(*img.shape)
>>> denoised_img = denoise_tv_chambolle(img, weight=60)
3D example on synthetic data:
>>> x, y, z = np.ogrid[0:20, 0:20, 0:20]
>>> mask = (x  22)**2 + (y  20)**2 + (z  17)**2 < 8**2
>>> mask = mask.astype(np.float)
>>> mask += 0.2*np.random.randn(*mask.shape)
>>> res = denoise_tv_chambolle(mask, weight=100)
skimage.restoration.
denoise_bilateral
(image, win_size=None, sigma_color=None, sigma_spatial=1, bins=10000, mode='constant', cval=0, multichannel=None)[source]¶Denoise image using bilateral filter.
This is an edgepreserving, denoising filter. It averages pixels based on their spatial closeness and radiometric similarity [R417].
Spatial closeness is measured by the Gaussian function of the Euclidean distance between two pixels and a certain standard deviation (sigma_spatial).
Radiometric similarity is measured by the Gaussian function of the Euclidean distance between two color values and a certain standard deviation (sigma_color).
Parameters:  image : ndarray, shape (M, N[, 3])
win_size : int
sigma_color : float
sigma_spatial : float
bins : int
mode : {‘constant’, ‘edge’, ‘symmetric’, ‘reflect’, ‘wrap’}
cval : string
multichannel : bool


Returns:  denoised : ndarray

References
[R417]  (1, 2) http://users.soe.ucsc.edu/~manduchi/Papers/ICCV98.pdf 
Examples
>>> from skimage import data, img_as_float
>>> astro = img_as_float(data.astronaut())
>>> astro = astro[220:300, 220:320]
>>> noisy = astro + 0.6 * astro.std() * np.random.random(astro.shape)
>>> noisy = np.clip(noisy, 0, 1)
>>> denoised = denoise_bilateral(noisy, sigma_color=0.05, sigma_spatial=15)
skimage.restoration.
denoise_wavelet
(image, sigma=None, wavelet='db1', mode='soft', wavelet_levels=None, multichannel=False, convert2ycbcr=False, method='BayesShrink')[source]¶Perform wavelet denoising on an image.
Parameters:  image : ndarray ([M[, N[, …P]][, C]) of ints, uints or floats
sigma : float or list, optional
wavelet : string, optional
mode : {‘soft’, ‘hard’}, optional
wavelet_levels : int or None, optional
multichannel : bool, optional
convert2ycbcr : bool, optional
method : {‘BayesShrink’, ‘VisuShrink’}, optional 

Returns:  out : ndarray

Notes
The wavelet domain is a sparse representation of the image, and can be thought of similarly to the frequency domain of the Fourier transform. Sparse representations have most values zero or nearzero and truly random noise is (usually) represented by many small values in the wavelet domain. Setting all values below some threshold to 0 reduces the noise in the image, but larger thresholds also decrease the detail present in the image.
If the input is 3D, this function performs wavelet denoising on each color plane separately. The output image is clipped between either [1, 1] and [0, 1] depending on the input image range.
When YCbCr conversion is done, every color channel is scaled between 0 and 1, and sigma values are applied to these scaled color channels.
Many wavelet coefficient thresholding approaches have been proposed. By
default, denoise_wavelet
applies BayesShrink, which is an adaptive
thresholding method that computes separate thresholds for each wavelet
subband as described in [R418].
If method == "VisuShrink"
, a single “universal threshold” is applied to
all wavelet detail coefficients as described in [R419]. This threshold
is designed to remove all Gaussian noise at a given sigma
with high
probability, but tends to produce images that appear overly smooth.
References
[R418]  (1, 2, 3) Chang, S. Grace, Bin Yu, and Martin Vetterli. “Adaptive wavelet thresholding for image denoising and compression.” Image Processing, IEEE Transactions on 9.9 (2000): 15321546. DOI: 10.1109/83.862633 
[R419]  (1, 2, 3, 4) D. L. Donoho and I. M. Johnstone. “Ideal spatial adaptation by wavelet shrinkage.” Biometrika 81.3 (1994): 425455. DOI: 10.1093/biomet/81.3.425 
Examples
>>> from skimage import color, data
>>> img = img_as_float(data.astronaut())
>>> img = color.rgb2gray(img)
>>> img += 0.1 * np.random.randn(*img.shape)
>>> img = np.clip(img, 0, 1)
>>> denoised_img = denoise_wavelet(img, sigma=0.1)
skimage.restoration.
denoise_nl_means
(image, patch_size=7, patch_distance=11, h=0.1, multichannel=None, fast_mode=True, sigma=0.0)[source]¶Perform nonlocal means denoising on 2D or 3D grayscale images, and 2D RGB images.
Parameters:  image : 2D or 3D ndarray
patch_size : int, optional
patch_distance : int, optional
h : float, optional
multichannel : bool, optional
fast_mode : bool, optional
sigma : float, optional


Returns:  result : ndarray

Notes
The nonlocal means algorithm is well suited for denoising images with specific textures. The principle of the algorithm is to average the value of a given pixel with values of other pixels in a limited neighbourhood, provided that the patches centered on the other pixels are similar enough to the patch centered on the pixel of interest.
In the original version of the algorithm [R420], corresponding to
fast=False
, the computational complexity is:
image.size * patch_size ** image.ndim * patch_distance ** image.ndim
Hence, changing the size of patches or their maximal distance has a strong effect on computing times, especially for 3D images.
However, the default behavior corresponds to fast_mode=True
, for which
another version of nonlocal means [R421] is used, corresponding to a
complexity of:
image.size * patch_distance ** image.ndim
The computing time depends only weakly on the patch size, thanks to
the computation of the integral of patches distances for a given
shift, that reduces the number of operations [R420]. Therefore, this
algorithm executes faster than the classic algorith
(fast_mode=False
), at the expense of using twice as much memory.
This implementation has been proven to be more efficient compared to
other alternatives, see e.g. [R422].
Compared to the classic algorithm, all pixels of a patch contribute to the distance to another patch with the same weight, no matter their distance to the center of the patch. This coarser computation of the distance can result in a slightly poorer denoising performance. Moreover, for small images (images with a linear size that is only a few times the patch size), the classic algorithm can be faster due to boundary effects.
The image is padded using the reflect mode of skimage.util.pad before denoising.
If the noise standard deviation, sigma, is provided a more robust computation of patch weights is used. Subtracting the known noise variance from the computed patch distances improves the estimates of patch similarity, giving a moderate improvement to denoising performance [R423]. It was also mentioned as an option for the fast variant of the algorithm in [R422].
When sigma is provided, a smaller h should typically be used to
avoid oversmoothing. The optimal value for h depends on the image
content and noise level, but a reasonable starting point is
h = 0.8 * sigma
when fast_mode is True, or h = 0.6 * sigma
when
fast_mode is False.
References
[R420]  (1, 2, 3) A. Buades, B. Coll, & JM. Morel. A nonlocal algorithm for image denoising. In CVPR 2005, Vol. 2, pp. 6065, IEEE. DOI: 10.1109/CVPR.2005.38 
[R421]  (1, 2) J. Darbon, A. Cunha, T.F. Chan, S. Osher, and G.J. Jensen, Fast nonlocal filtering applied to electron cryomicroscopy, in 5th IEEE International Symposium on Biomedical Imaging: From Nano to Macro, 2008, pp. 13311334. DOI: 10.1109/ISBI.2008.4541250 
[R422]  (1, 2, 3) Jacques Froment. ParameterFree Fast Pixelwise NonLocal Means Denoising. Image Processing On Line, 2014, vol. 4, pp. 300326. DOI: 10.5201/ipol.2014.120 
[R423]  (1, 2) A. Buades, B. Coll, & JM. Morel. NonLocal Means Denoising. Image Processing On Line, 2011, vol. 1, pp. 208212. DOI: 10.5201/ipol.2011.bcm_nlm 
Examples
>>> a = np.zeros((40, 40))
>>> a[10:10, 10:10] = 1.
>>> a += 0.3 * np.random.randn(*a.shape)
>>> denoised_a = denoise_nl_means(a, 7, 5, 0.1)
skimage.restoration.
estimate_sigma
(image, average_sigmas=False, multichannel=False)[source]¶Robust waveletbased estimator of the (Gaussian) noise standard deviation.
Parameters:  image : ndarray
average_sigmas : bool, optional
multichannel : bool


Returns:  sigma : float or list

Notes
This function assumes the noise follows a Gaussian distribution. The estimation algorithm is based on the median absolute deviation of the wavelet detail coefficients as described in section 4.2 of [R424].
References
[R424]  (1, 2) D. L. Donoho and I. M. Johnstone. “Ideal spatial adaptation by wavelet shrinkage.” Biometrika 81.3 (1994): 425455. DOI:10.1093/biomet/81.3.425 
Examples
>>> import skimage.data
>>> from skimage import img_as_float
>>> img = img_as_float(skimage.data.camera())
>>> sigma = 0.1
>>> img = img + sigma * np.random.standard_normal(img.shape)
>>> sigma_hat = estimate_sigma(img, multichannel=False)
skimage.restoration.
inpaint_biharmonic
(image, mask, multichannel=False)[source]¶Inpaint masked points in image with biharmonic equations.
Parameters:  image : (M[, N[, …, P]][, C]) ndarray
mask : (M[, N[, …, P]]) ndarray
multichannel : boolean, optional


Returns:  out : (M[, N[, …, P]][, C]) ndarray

References
[R425]  N.S.Hoang, S.B.Damelin, “On surface completion and image inpainting by biharmonic functions: numerical aspects”, https://arxiv.org/abs/1707.06567 
[R426]  C. K. Chui and H. N. Mhaskar, MRA ContextualRecovery Extension of Smooth Functions on Manifolds, Appl. and Comp. Harmonic Anal., 28 (2010), 104113, DOI: 10.1016/j.acha.2009.04.004 
Examples
>>> img = np.tile(np.square(np.linspace(0, 1, 5)), (5, 1))
>>> mask = np.zeros_like(img)
>>> mask[2, 2:] = 1
>>> mask[1, 3:] = 1
>>> mask[0, 4:] = 1
>>> out = inpaint_biharmonic(img, mask)
skimage.restoration.
cycle_spin
(x, func, max_shifts, shift_steps=1, num_workers=None, multichannel=False, func_kw={})[source]¶Cycle spinning (repeatedly apply func to shifted versions of x).
Parameters:  x : arraylike
func : function
max_shifts : int or tuple
shift_steps : int or tuple, optional
num_workers : int or None, optional
multichannel : bool, optional
func_kw : dict, optional


Returns:  avg_y : np.ndarray

Notes
Cycle spinning was proposed as a way to approach shiftinvariance via performing several circular shifts of a shiftvariant transform [R427].
For a nlevel discrete wavelet transforms, one may wish to perform all
shifts up to max_shifts = 2**n  1
. In practice, much of the benefit
can often be realized with only a small number of shifts per axis.
For transforms such as the blockwise discrete cosine transform, one may wish to evaluate shifts up to the block size used by the transform.
References
[R427]  (1, 2) R.R. Coifman and D.L. Donoho. “TranslationInvariant DeNoising”. Wavelets and Statistics, Lecture Notes in Statistics, vol.103. Springer, New York, 1995, pp.125150. DOI:10.1007/9781461225447_9 
Examples
>>> import skimage.data
>>> from skimage import img_as_float
>>> from skimage.restoration import denoise_wavelet, cycle_spin
>>> img = img_as_float(skimage.data.camera())
>>> sigma = 0.1
>>> img = img + sigma * np.random.standard_normal(img.shape)
>>> denoised = cycle_spin(img, func=denoise_wavelet, max_shifts=3)