denoise
|
Run benchmarks for module using nose. |
|
Run tests for module using nose. |
denoise.adaptive_soft_matching
|
Adaptive Soft Coefficient Matching |
denoise.gibbs
alias of |
|
partial(func, *args, **keywords) - new function with partial application of the given arguments and keywords. |
|
|
Returns a process pool object |
|
Suppresses Gibbs ringing artefacts of images volumes. |
denoise.localpca
|
Solve a standard or generalized eigenvalue problem for a complex Hermitian or real symmetric matrix. |
|
General function to perform PCA-based denoising of diffusion datasets. |
|
Performs local PCA denoising according to Manjon et al. |
|
Performs PCA-based denoising using the Marcenko-Pastur distribution [1]. |
denoise.nlmeans
|
Non-local means for denoising 3D and 4D images |
Non-local means for denoising 3D images |
denoise.noise_estimate
|
Multidimensional convolution. |
|
Standard deviation estimation from local patches |
|
Probabilistic Identification and Estimation of Noise (PIESNO). |
denoise.non_local_means
Non-Local Means Denoising Using Blockwise Averaging |
|
|
Non-local means for denoising 3D and 4D images, using |
dipy.denoise.
bench
(label='fast', verbose=1, extra_argv=None)Run benchmarks for module using nose.
Identifies the benchmarks to run. This can be a string to pass to the nosetests executable with the ‘-A’ option, or one of several special values. Special values are:
‘fast’ - the default - which corresponds to the nosetests -A
option of ‘not slow’.
‘full’ - fast (as above) and slow benchmarks as in the ‘no -A’ option to nosetests - this is the same as ‘’.
None or ‘’ - run all tests.
attribute_identifier - string passed directly to nosetests as ‘-A’.
Verbosity value for benchmark outputs, in the range 1-10. Default is 1.
List with any extra arguments to pass to nosetests.
Returns True if running the benchmarks works, False if an error occurred.
Notes
Benchmarks are like tests, but have names starting with “bench” instead of “test”, and can be found under the “benchmarks” sub-directory of the module.
Each NumPy module exposes bench in its namespace to run all benchmarks for it.
Examples
>>> success = np.lib.bench()
Running benchmarks for numpy.lib
...
using 562341 items:
unique:
0.11
unique1d:
0.11
ratio: 1.0
nUnique: 56230 == 56230
...
OK
>>> success
True
dipy.denoise.
test
(label='fast', verbose=1, extra_argv=None, doctests=False, coverage=False, raise_warnings=None, timer=False)Run tests for module using nose.
Identifies the tests to run. This can be a string to pass to the nosetests executable with the ‘-A’ option, or one of several special values. Special values are:
‘fast’ - the default - which corresponds to the nosetests -A
option of ‘not slow’.
‘full’ - fast (as above) and slow tests as in the ‘no -A’ option to nosetests - this is the same as ‘’.
None or ‘’ - run all tests.
attribute_identifier - string passed directly to nosetests as ‘-A’.
Verbosity value for test outputs, in the range 1-10. Default is 1.
List with any extra arguments to pass to nosetests.
If True, run doctests in module. Default is False.
If True, report coverage of NumPy code. Default is False. (This requires the coverage module).
This specifies which warnings to configure as ‘raise’ instead of being shown once during the test execution. Valid strings are:
“develop” : equals (Warning,)
“release” : equals ()
, do not raise on any warnings.
Timing of individual tests with nose-timer
(which needs to be
installed). If True, time tests and report on all of them.
If an integer (say N
), report timing results for N
slowest
tests.
Returns the result of running the tests as a
nose.result.TextTestResult
object.
Notes
Each NumPy module exposes test in its namespace to run all tests for it. For example, to run all tests for numpy.lib:
>>> np.lib.test()
Examples
>>> result = np.lib.test()
Running unit tests for numpy.lib
...
Ran 976 tests in 3.933s
OK
>>> result.errors
[]
>>> result.knownfail
[]
dipy.denoise.adaptive_soft_matching.
adaptive_soft_matching
(ima, fimau, fimao, sigma)Adaptive Soft Coefficient Matching
Combines two filtered 3D-images at different resolutions and the original image. Returns the resulting combined image.
filtered image with optimized non-local means using a small block (suggested:3x3), which corresponds to a “high resolution” filter.
filtered image with optimized non-local means using a small block (suggested:5x5), which corresponds to a “low resolution” filter.
that explain the rician noise. Note: In P. Coupe et al. the rician noise was simulated as sqrt((f+x)^2 + (y)^2) where f is the pixel value and x and y are independent realizations of a random variable with Normal distribution, with mean=0 and standard deviation=h
output denoised array which is of the same shape as that of the input
References
Pierrick Coupe, Jose Manjon, Montserrat Robles, Louis Collins. “Multiresolution Non-Local Means Filter for 3D MR Image Denoising” IET Image Processing, Institution of Engineering and Technology, 2011
partial
dipy.denoise.gibbs.
gibbs_removal
(vol, slice_axis=2, n_points=3, inplace=True, num_threads=1)Suppresses Gibbs ringing artefacts of images volumes.
Matrix containing one volume (3D) or multiple (4D) volumes of images.
Data axis corresponding to the number of acquired slices. Default is set to the third axis.
Number of neighbour points to access local TV (see note). Default is set to 3.
If True, the input data is replaced with results. Otherwise, returns a new array. Default is set to True.
Number of threads. Only applies to 3D or 4D data arrays. If None then all available threads will be used. Otherwise, must be a positive integer. Default is set to 1.
Matrix containing one volume (3D) or multiple (4D) volumes of corrected images.
Notes
For 4D matrix last element should always correspond to the number of diffusion gradient directions.
References
Please cite the following articles .. [Rfce522a872fd-1] Neto Henriques, R., 2018. Advanced Methods for Diffusion MRI Data
Analysis and their Application to the Healthy Ageing Brain (Doctoral thesis). https://doi.org/10.17863/CAM.29356
Kellner E, Dhital B, Kiselev VG, Reisert M. Gibbs-ringing artifact removal based on local subvoxel-shifts. Magn Reson Med. 2016 doi: 10.1002/mrm.26054.
dipy.denoise.localpca.
eigh
(a, b=None, lower=True, eigvals_only=False, overwrite_a=False, overwrite_b=False, turbo=True, eigvals=None, type=1, check_finite=True, subset_by_index=None, subset_by_value=None, driver=None)Solve a standard or generalized eigenvalue problem for a complex Hermitian or real symmetric matrix.
Find eigenvalues array w
and optionally eigenvectors array v
of
array a
, where b
is positive definite such that for every
eigenvalue λ (i-th entry of w) and its eigenvector vi
(i-th column of
v
) satisfies:
a @ vi = λ * b @ vi
vi.conj().T @ a @ vi = λ
vi.conj().T @ b @ vi = 1
In the standard problem, b
is assumed to be the identity matrix.
A complex Hermitian or real symmetric matrix whose eigenvalues and eigenvectors will be computed.
A complex Hermitian or real symmetric definite positive matrix in. If omitted, identity matrix is assumed.
Whether the pertinent array data is taken from the lower or upper
triangle of a
and, if applicable, b
. (Default: lower)
Whether to calculate only eigenvalues and no eigenvectors. (Default: both are calculated)
If provided, this two-element iterable defines the start and the end
indices of the desired eigenvalues (ascending order and 0-indexed).
To return only the second smallest to fifth smallest eigenvalues,
[1, 4]
is used. [n-3, n-1]
returns the largest three. Only
available with “evr”, “evx”, and “gvx” drivers. The entries are
directly converted to integers via int()
.
If provided, this two-element iterable defines the half-open interval
(a, b]
that, if any, only the eigenvalues between these values
are returned. Only available with “evr”, “evx”, and “gvx” drivers. Use
np.inf
for the unconstrained ends.
Defines which LAPACK driver should be used. Valid options are “ev”, “evd”, “evr”, “evx” for standard problems and “gv”, “gvd”, “gvx” for generalized (where b is not None) problems. See the Notes section.
For the generalized problems, this keyword specifies the problem type
to be solved for w
and v
(only takes 1, 2, 3 as possible
inputs):
1 => a @ v = w @ b @ v
2 => a @ b @ v = w @ v
3 => b @ a @ v = w @ v
This keyword is ignored for standard problems.
Whether to overwrite data in a
(may improve performance). Default
is False.
Whether to overwrite data in b
(may improve performance). Default
is False.
Whether to check that the input matrices contain only finite numbers. Disabling may give a performance gain, but may result in problems (crashes, non-termination) if the inputs do contain infinities or NaNs.
Deprecated since v1.5.0, use ``driver=gvd`` keyword instead. Use divide and conquer algorithm (faster but expensive in memory, only for generalized eigenvalue problem and if full set of eigenvalues are requested.). Has no significant effect if eigenvectors are not requested.
Deprecated since v1.5.0, use ``subset_by_index`` keyword instead. Indexes of the smallest and largest (in ascending order) eigenvalues and corresponding eigenvectors to be returned: 0 <= lo <= hi <= M-1. If omitted, all eigenvalues and eigenvectors are returned.
The N (1<=N<=M) selected eigenvalues, in ascending order, each repeated according to its multiplicity.
(if eigvals_only == False
)
If eigenvalue computation does not converge, an error occurred, or b matrix is not definite positive. Note that if input matrices are not symmetric or Hermitian, no error will be reported but results will be wrong.
See also
eigvalsh
eigenvalues of symmetric or Hermitian arrays
eig
eigenvalues and right eigenvectors for non-symmetric arrays
eigh_tridiagonal
eigenvalues and right eiegenvectors for symmetric/Hermitian tridiagonal matrices
Notes
This function does not check the input array for being hermitian/symmetric in order to allow for representing arrays with only their upper/lower triangular parts. Also, note that even though not taken into account, finiteness check applies to the whole array and unaffected by “lower” keyword.
This function uses LAPACK drivers for computations in all possible keyword
combinations, prefixed with sy
if arrays are real and he
if
complex, e.g., a float array with “evr” driver is solved via
“syevr”, complex arrays with “gvx” driver problem is solved via “hegvx”
etc.
As a brief summary, the slowest and the most robust driver is the
classical <sy/he>ev
which uses symmetric QR. <sy/he>evr
is seen as
the optimal choice for the most general cases. However, there are certain
occassions that <sy/he>evd
computes faster at the expense of more
memory usage. <sy/he>evx
, while still being faster than <sy/he>ev
,
often performs worse than the rest except when very few eigenvalues are
requested for large arrays though there is still no performance guarantee.
For the generalized problem, normalization with respoect to the given type argument:
type 1 and 3 : v.conj().T @ a @ v = w
type 2 : inv(v).conj().T @ a @ inv(v) = w
type 1 or 2 : v.conj().T @ b @ v = I
type 3 : v.conj().T @ inv(b) @ v = I
Examples
>>> from scipy.linalg import eigh
>>> A = np.array([[6, 3, 1, 5], [3, 0, 5, 1], [1, 5, 6, 2], [5, 1, 2, 2]])
>>> w, v = eigh(A)
>>> np.allclose(A @ v - v @ np.diag(w), np.zeros((4, 4)))
True
Request only the eigenvalues
>>> w = eigh(A, eigvals_only=True)
Request eigenvalues that are less than 10.
>>> A = np.array([[34, -4, -10, -7, 2],
... [-4, 7, 2, 12, 0],
... [-10, 2, 44, 2, -19],
... [-7, 12, 2, 79, -34],
... [2, 0, -19, -34, 29]])
>>> eigh(A, eigvals_only=True, subset_by_value=[-np.inf, 10])
array([6.69199443e-07, 9.11938152e+00])
Request the largest second eigenvalue and its eigenvector
>>> w, v = eigh(A, subset_by_index=[1, 1])
>>> w
array([9.11938152])
>>> v.shape # only a single column is returned
(5, 1)
dipy.denoise.localpca.
genpca
(arr, sigma=None, mask=None, patch_radius=2, pca_method='eig', tau_factor=None, return_sigma=False, out_dtype=None)General function to perform PCA-based denoising of diffusion datasets.
Array of data to be denoised. The dimensions are (X, Y, Z, N), where N are the diffusion gradient directions.
Standard deviation of the noise estimated from the data. If no sigma is given, this will be estimated based on random matrix theory [1],[R4bed7265c934-2]_
A mask with voxels that are true inside the brain and false outside of it. The function denoises within the true part and returns zeros outside of those voxels.
The radius of the local patch to be taken around each voxel (in voxels). Default: 2 (denoise in blocks of 5x5x5 voxels).
Use either eigenvalue decomposition (eig) or singular value decomposition (svd) for principal component analysis. The default method is ‘eig’ which is faster. However, occasionally ‘svd’ might be more accurate.
Thresholding of PCA eigenvalues is done by nulling out eigenvalues that are smaller than:
tau_{factor} can be set to a predefined values (e.g. tau_{factor} = 2.3 [3]), or automatically calculated using random matrix theory (in case that tau_{factor} is set to None). Default: None.
If true, the Standard deviation of the noise will be returned. Default: False.
The dtype for the output array. Default: output has the same dtype as the input.
This is the denoised array of the same size as that of the input data, clipped to non-negative values
References
Veraart J, Novikov DS, Christiaens D, Ades-aron B, Sijbers, Fieremans E, 2016. Denoising of Diffusion MRI using random matrix theory. Neuroimage 142:394-406. doi: 10.1016/j.neuroimage.2016.08.016
Veraart J, Fieremans E, Novikov DS. 2016. Diffusion MRI noise mapping using random matrix theory. Magnetic Resonance in Medicine. doi: 10.1002/mrm.26059.
Manjon JV, Coupe P, Concha L, Buades A, Collins DL (2013) Diffusion Weighted Image Denoising Using Overcomplete Local PCA. PLoS ONE 8(9): e73021. https://doi.org/10.1371/journal.pone.0073021
dipy.denoise.localpca.
localpca
(arr, sigma, mask=None, patch_radius=2, pca_method='eig', tau_factor=2.3, out_dtype=None)Performs local PCA denoising according to Manjon et al. [1].
Array of data to be denoised. The dimensions are (X, Y, Z, N), where N are the diffusion gradient directions.
Standard deviation of the noise estimated from the data.
A mask with voxels that are true inside the brain and false outside of it. The function denoises within the true part and returns zeros outside of those voxels.
The radius of the local patch to be taken around each voxel (in voxels). Default: 2 (denoise in blocks of 5x5x5 voxels).
Use either eigenvalue decomposition (eig) or singular value decomposition (svd) for principal component analysis. The default method is ‘eig’ which is faster. However, occasionally ‘svd’ might be more accurate.
Thresholding of PCA eigenvalues is done by nulling out eigenvalues that are smaller than:
tau_{factor} can be change to adjust the relationship between the noise standard deviation and the threshold tau. If tau_{factor} is set to None, it will be automatically calculated using the Marcenko-Pastur distribution [2]. Default: 2.3 (according to [1])
The dtype for the output array. Default: output has the same dtype as the input.
This is the denoised array of the same size as that of the input data, clipped to non-negative values
References
Manjon JV, Coupe P, Concha L, Buades A, Collins DL (2013) Diffusion Weighted Image Denoising Using Overcomplete Local PCA. PLoS ONE 8(9): e73021. https://doi.org/10.1371/journal.pone.0073021
Veraart J, Novikov DS, Christiaens D, Ades-aron B, Sijbers, Fieremans E, 2016. Denoising of Diffusion MRI using random matrix theory. Neuroimage 142:394-406. doi: 10.1016/j.neuroimage.2016.08.016
dipy.denoise.localpca.
mppca
(arr, mask=None, patch_radius=2, pca_method='eig', return_sigma=False, out_dtype=None)Performs PCA-based denoising using the Marcenko-Pastur distribution [1].
Array of data to be denoised. The dimensions are (X, Y, Z, N), where N are the diffusion gradient directions.
A mask with voxels that are true inside the brain and false outside of it. The function denoises within the true part and returns zeros outside of those voxels.
The radius of the local patch to be taken around each voxel (in voxels). Default: 2 (denoise in blocks of 5x5x5 voxels).
Use either eigenvalue decomposition (eig) or singular value decomposition (svd) for principal component analysis. The default method is ‘eig’ which is faster. However, occasionally ‘svd’ might be more accurate.
If true, a noise standard deviation estimate based on the Marcenko-Pastur distribution is returned [2]. Default: False.
The dtype for the output array. Default: output has the same dtype as the input.
This is the denoised array of the same size as that of the input data, clipped to non-negative values
Estimate of the spatial varying standard deviation of the noise
References
Veraart J, Novikov DS, Christiaens D, Ades-aron B, Sijbers, Fieremans E, 2016. Denoising of Diffusion MRI using random matrix theory. Neuroimage 142:394-406. doi: 10.1016/j.neuroimage.2016.08.016
Veraart J, Fieremans E, Novikov DS. 2016. Diffusion MRI noise mapping using random matrix theory. Magnetic Resonance in Medicine. doi: 10.1002/mrm.26059.
dipy.denoise.nlmeans.
nlmeans
(arr, sigma, mask=None, patch_radius=1, block_radius=5, rician=True, num_threads=None)Non-local means for denoising 3D and 4D images
The array to be denoised
standard deviation of the noise estimated from the data
patch size is 2 x patch_radius + 1
. Default is 1.
block size is 2 x block_radius + 1
. Default is 5.
If True the noise is estimated as Rician, otherwise Gaussian noise is assumed.
Number of threads. If None (default) then all available threads will be used (all CPU cores).
the denoised arr
which has the same shape as arr
.
References
Descoteaux, Maxime and Wiest-Daesslé, Nicolas and Prima, Sylvain and Barillot, Christian and Deriche, Rachid Impact of Rician Adapted Non-Local Means Filtering on HARDI, MICCAI 2008
dipy.denoise.nlmeans.
nlmeans_3d
()Non-local means for denoising 3D images
The array to be denoised
standard deviation of the noise estimated from the data
patch size is 2 x patch_radius + 1
. Default is 1.
block size is 2 x block_radius + 1
. Default is 5.
If True the noise is estimated as Rician, otherwise Gaussian noise is assumed.
Number of threads. If None (default) then all available threads will be used.
the denoised arr
which has the same shape as arr
.
dipy.denoise.noise_estimate.
convolve
(input, weights, output=None, mode='reflect', cval=0.0, origin=0)Multidimensional convolution.
The array is convolved with the given kernel.
The input array.
Array of weights, same number of dimensions as input
The array in which to place the output, or the dtype of the returned array. By default an array of the same dtype as input will be created.
The mode parameter determines how the input array is extended beyond its boundaries. Default is ‘reflect’. Behavior for each valid value is as follows:
The input is extended by reflecting about the edge of the last pixel.
The input is extended by filling all values beyond the edge with the same constant value, defined by the cval parameter.
The input is extended by replicating the last pixel.
The input is extended by reflecting about the center of the last pixel.
The input is extended by wrapping around to the opposite edge.
Value to fill past edges of input if mode is ‘constant’. Default is 0.0
Controls the placement of the filter on the input array’s pixels. A value of 0 (the default) centers the filter over the pixel, with positive values shifting the filter to the left, and negative ones to the right. By passing a sequence of origins with length equal to the number of dimensions of the input array, different shifts can be specified along each axis.
The result of convolution of input with weights.
See also
correlate
Correlate an image with a kernel.
Notes
Each value in result is \(C_i = \sum_j{I_{i+k-j} W_j}\), where W is the weights kernel, j is the N-D spatial index over \(W\), I is the input and k is the coordinate of the center of W, specified by origin in the input parameters.
Examples
Perhaps the simplest case to understand is mode='constant', cval=0.0
,
because in this case borders (i.e., where the weights kernel, centered
on any one value, extends beyond an edge of input) are treated as zeros.
>>> a = np.array([[1, 2, 0, 0],
... [5, 3, 0, 4],
... [0, 0, 0, 7],
... [9, 3, 0, 0]])
>>> k = np.array([[1,1,1],[1,1,0],[1,0,0]])
>>> from scipy import ndimage
>>> ndimage.convolve(a, k, mode='constant', cval=0.0)
array([[11, 10, 7, 4],
[10, 3, 11, 11],
[15, 12, 14, 7],
[12, 3, 7, 0]])
Setting cval=1.0
is equivalent to padding the outer edge of input
with 1.0’s (and then extracting only the original region of the result).
>>> ndimage.convolve(a, k, mode='constant', cval=1.0)
array([[13, 11, 8, 7],
[11, 3, 11, 14],
[16, 12, 14, 10],
[15, 6, 10, 5]])
With mode='reflect'
(the default), outer values are reflected at the
edge of input to fill in missing values.
>>> b = np.array([[2, 0, 0],
... [1, 0, 0],
... [0, 0, 0]])
>>> k = np.array([[0,1,0], [0,1,0], [0,1,0]])
>>> ndimage.convolve(b, k, mode='reflect')
array([[5, 0, 0],
[3, 0, 0],
[1, 0, 0]])
This includes diagonally at the corners.
>>> k = np.array([[1,0,0],[0,1,0],[0,0,1]])
>>> ndimage.convolve(b, k)
array([[4, 2, 0],
[3, 2, 0],
[1, 1, 0]])
With mode='nearest'
, the single nearest value in to an edge in
input is repeated as many times as needed to match the overlapping
weights.
>>> c = np.array([[2, 0, 1],
... [1, 0, 0],
... [0, 0, 0]])
>>> k = np.array([[0, 1, 0],
... [0, 1, 0],
... [0, 1, 0],
... [0, 1, 0],
... [0, 1, 0]])
>>> ndimage.convolve(c, k, mode='nearest')
array([[7, 0, 3],
[5, 0, 2],
[3, 0, 1]])
dipy.denoise.noise_estimate.
estimate_sigma
(arr, disable_background_masking=False, N=0)Standard deviation estimation from local patches
The array to be estimated
If True, uses all voxels for the estimation, otherwise, only non-zeros voxels are used. Useful if the background is masked by the scanner.
Number of coils of the receiver array. Use N = 1 in case of a SENSE reconstruction (Philips scanners) or the number of coils for a GRAPPA reconstruction (Siemens and GE). Use 0 to disable the correction factor, as for example if the noise is Gaussian distributed. See [1] for more information.
standard deviation of the noise, one estimation per volume.
Notes
This function is the same as manually taking the standard deviation of the
background and gives one value for the whole 3D array.
It also includes the coil-dependent correction factor of Koay 2006
(see [1], equation 18) with theta = 0.
Since this function was introduced in [2] for T1 imaging,
it is expected to perform ok on diffusion MRI data, but might oversmooth
some regions and leave others un-denoised for spatially varying noise
profiles. Consider using piesno()
to estimate sigma instead if visual
inaccuracies are apparent in the denoised result.
References
Koay, C. G., & Basser, P. J. (2006). Analytically exact correction
scheme for signal extraction from noisy magnitude MR signals. Journal of Magnetic Resonance), 179(2), 317-22.
Coupe, P., Yger, P., Prima, S., Hellier, P., Kervrann, C., Barillot,
C., 2008. An optimized blockwise nonlocal means denoising filter for 3-D magnetic resonance images, IEEE Trans. Med. Imaging 27, 425-41.
dipy.denoise.noise_estimate.
piesno
(data, N, alpha=0.01, l=100, itermax=100, eps=1e-05, return_mask=False)Probabilistic Identification and Estimation of Noise (PIESNO).
The magnitude signals to analyse. The last dimension must contain the same realisation of the volume, such as dMRI or fMRI data.
The number of phase array coils of the MRI scanner. If your scanner does a SENSE reconstruction, ALWAYS use N=1, as the noise profile is always Rician. If your scanner does a GRAPPA reconstruction, set N as the number of phase array coils.
Probabilistic estimation threshold for the gamma function.
number of initial estimates for sigma to try.
Maximum number of iterations to execute if convergence is not reached.
Tolerance for the convergence criterion. Convergence is reached if two subsequent estimates are smaller than eps.
If True, return a mask identifying all the pure noise voxel that were found.
The estimated standard deviation of the gaussian noise.
A boolean mask indicating the voxels identified as pure noise.
Notes
This function assumes two things : 1. The data has a noisy, non-masked background and 2. The data is a repetition of the same measurements along the last axis, i.e. dMRI or fMRI data, not structural data like T1/T2.
This function processes the data slice by slice, as originally designed in the paper. Use it to get a slice by slice estimation of the noise, as in spinal cord imaging for example.
References
Koay CG, Ozarslan E and Pierpaoli C.
“Probabilistic Identification and Estimation of Noise (PIESNO): A self-consistent approach and its applications in MRI.” Journal of Magnetic Resonance 2009; 199: 94-103.
Koay CG, Ozarslan E and Basser PJ.
“A signal transformational framework for breaking the noise floor and its applications in MRI.” Journal of Magnetic Resonance 2009; 197: 108-119.
dipy.denoise.non_local_means.
nlmeans_block
()Non-Local Means Denoising Using Blockwise Averaging
the input image, corrupted with rician noise
the input mask
similar patches in the non-local means are searched for locally, inside a cube of side 2*v+1 centered at each voxel of interest.
the size of the block to be used (2*f+1)x(2*f+1)x(2*f+1) in the blockwise non-local means implementation (the Coupe’s proposal).
the estimated amount of rician noise in the input image: in P. Coupe et al. the rician noise was simulated as sqrt((f+x)^2 + (y)^2) where f is the pixel value and x and y are independent realizations of a random variable with Normal distribution, with mean=0 and standard deviation=h
If True the noise is estimated as Rician, otherwise Gaussian noise is assumed.
the denoised output which has the same shape as input image.
References
“An Optimized Blockwise Non Local Means Denoising Filter for 3D Magnetic Resonance Images” IEEE Transactions on Medical Imaging, 27(4):425-441, 2008
“Multiresolution Non-Local Means Filter for 3D MR Image Denoising” IET Image Processing, Institution of Engineering and Technology, 2011
dipy.denoise.non_local_means.
non_local_means
(arr, sigma, mask=None, patch_radius=1, block_radius=5, rician=True)blockwise averaging approach
The array to be denoised
standard deviation of the noise estimated from the data
patch size is 2 x patch_radius + 1
. Default is 1.
block size is 2 x block_radius + 1
. Default is 5.
If True the noise is estimated as Rician, otherwise Gaussian noise is assumed.
the denoised arr
which has the same shape as arr
.
References
P. Coupe, P. Yger, S. Prima, P. Hellier, C. Kervrann, C. Barillot, An Optimized Blockwise Non Local Means Denoising Filter for 3D Magnetic Resonance Images, IEEE Transactions on Medical Imaging, 27(4):425-441, 2008
Pierrick Coupe, Jose Manjon, Montserrat Robles, Louis Collins. Adaptive Multiresolution Non-Local Means Filter for 3D MR Image Denoising IET Image Processing, Institution of Engineering and Technology, 2011