Signal Processing ] (]]]]) ]]]–]]]
1
Contents lists available at ScienceDirect
3
Signal Processing
5
journal homepage: www.elsevier.com/locate/sigpro
7 9
Fast communication
11
Modified BM3D algorithm for image denoising using nonlocal centralization prior
13 15 Q1
Hua Zhong, Ke Ma, Yang Zhou
17 Q2
Key Laboratory of Intelligent Perception and Image Understanding of Ministry of Education, International Research Center for Intelligent Perception and Computation, Xidian University, Mailbox 224, No. 2 South TaiBai Road, Xi'an 710071, Shaanxi Province, China
19 21 23 25
a r t i c l e i n f o
abstract
Article history: Received 9 April 2014 Received in revised form 23 July 2014 Accepted 7 August 2014
Block matching and 3D collaborative filtering (BM3D) has shown great power for image denoising. For grouped image blocks, this letter proposes to remove the 1D transform inter-blocks and introduce the nonlocal centralization prior to better utilize both local sparsity of wavelet coefficients and nonlocal similarity of grouped blocks. Three nonlocal shrinkage functions are developed under different norm restrictions of wavelet coefficients intra- and inter-blocks. Such shrinkage functions are verified to be competitive or better than the BM3D algorithm and other state-of-the-art denoising methods. & 2014 Elsevier B.V. All rights reserved.
27 Keywords: Denoising Wavelet Shrinkage Block matching Nonlocal Sparsity
29 31 33 35 37 39 41 43 45 47 49 51 53 55
1. Introduction Q3
In recent years, nonlocal similarity between patches in one image has attracted a lot of attention. Many excellent denoising methods have been developed since the nonlocal means (NLM) proposed by Buads et al. [1]. Among them, to name a few, methods such as BM3D [2], centralized sparse representation (CSR) [3] and learned simultaneous sparse coding (LSSC) [4], etc., achieve state of the art denoising performance by combining both local and nonlocal similarity of grouped blocks. As the first one that integrates the advantage of wavelet shrinkage [5] and nonlocal similarity, the BM3D algorithm succeeds in combining two critical operations: the block matching and 3D collaborative filtering. Based on these two operations, its several variants [6] are developed with improved performance. Other attempts replace wavelet transform with the methods of dictionary learning [3,4,7],
57 E-mail address:
[email protected] (H. Zhong).
59
which can provide non-fixed, adaptive atoms (bases) for the representation of the images. The combination of dictionary learning and nonlocal similarity forms another prospective research interest but with increased complexity. The focus of this letter is on the combination of wavelet shrinkage and nonlocal similarity. The BM3D algorithm provides an effective and efficient framework to achieve this purpose. Through effective block matching, the 3D transform can achieve highly sparse representation of grouped blocks. This sparsity makes the so called 3D collaborative filtering very effective in preserving the details shared by similar blocks and also the difference between them. The essence of the BM3D algorithm is to simultaneously utilize both local sparsity of wavelet coefficients within one block, and the nonlocal similarity between group blocks. However, the 3D collaborative filtering employs only simple thresholding function, which is originally derived for the cases of 1D and 2D wavelet shrinkage under the assumption of local sparsity of wavelet coefficients. The nonlocal similarity might not be fully
http://dx.doi.org/10.1016/j.sigpro.2014.08.014 0165-1684/& 2014 Elsevier B.V. All rights reserved.
61 Please cite this article as: H. Zhong, et al., Modified BM3D algorithm for image denoising using nonlocal centralization prior, Signal Processing (2014), http://dx.doi.org/10.1016/j.sigpro.2014.08.014i
H. Zhong et al. / Signal Processing ] (]]]]) ]]]–]]]
2
1 3 5 7 9 11 13 15 17 19 21
exploited. In fact, image details shared by similar blocks are possible to be weakened or blurred by the BM3D algorithm, particularly under the circumstance of strong noise. In this letter, nonlocal centralization prior [3] is introduced into the BM3D framework. Only 2D wavelet transform is applied to blocks and different nonlocal shrinkage functions are derived for the filtering task instead of the simple thresholding function. This modification can further exploit the intra- and inter-block correlation, providing a promising new way of integrating wavelet shrinkage with nonlocal similarity. The derived nonlocal shrinkage functions use different norm restrictions on wavelet coefficients intra- and inter-blocks. Experimental results are shown to demonstrate the potential of the proposed method when compared with the BM3D algorithm and other state-of-theart denoising methods. This letter is organized as follows. Section 2 briefly reviews the BM3D algorithm and the wavelet shrinkage model. Section 3 introduces the proposed method and the detailed implementation. Results and conclusions are presented in Sections 4 and 5, respectively.
23 25
2. BM3D and wavelet shrinkage
27
45
As one kind of nonlocal filter, the BM3D algorithm consists of two main operations: block matching and 3D collaborative filtering. The purpose of block matching is to group image blocks similar to a given reference block, and then stack them together into a 3D array. Let Z be the reference block, the Euclidian distance between Z and all its candidate blocks in the fixed neighborhood will be calculated. The m nearest blocks will be selected to form the group, denoted as S ¼ fZ j ; j ¼ 1; :::; mg. The 3D collaborative filtering will be performed on the group S, which can be seen as a 3D array of these blocks. In BM3D algorithm, 2D wavelet transform will be applied to each block at first for the representation of local sparsity. the 3rd transform such as 1D wavelet and PCA, which utilizes nonlocal similarity within the group S, will be performed across the coefficients of different blocks. Let 1 T 3D be the above inverse 3D transformation and α be the noisy coefficient, the following shrinkage model is introduced for denoising task [5]:
47
1 α^ ¼ arg min 12 jjZ T 3D αjj þ λjjαjj
29 31 33 35 37 39 41 43
49 51 53 55 57 59 61
α
ð1Þ
where λ is the regularization multiplier and jj Ujjp is the p norm, also written as Lp norm. Finally, such block-wise estimate α^ will be inversely transformed and aggregated to obtain the denoised image. More implementation details can be obtained in [2]. Notice that though (1) is in fact a kind of nonlocal shrinkage, the employed hard or soft shrinkage function is still a local way derived under the assumption of local sparsity [5] and might not be optimal in exploiting the nonlocal similarity within the group S. As a result, it is preferred to derive a shrinkage function in a nonlocal fashion, which provides another way of integrating local sparsity and nonlocal similarity.
3. The proposed method In this letter, we remove the 1D transform across grouped blocks and introduce into (1) the nonlocal centralization prior [3]. The modified shrinkage model is expressed as 1 α^ ¼ arg min 12 jjZ T 2D αjj22 þ λ1 jjαjjp þ λ2 jjα βjjq α
ð2Þ
1 and α are the inverse 2D wavelet transform and where T 2D the coefficient, respectively. The term jjα βjjq originates from the CSR model [3] but generalizing the L1 norm to Lq norm. β denotes the nonlocal mean of the coefficients at the same position across similar blocks. λ1 and λ2 are the Lagrangian multipliers to balance between local sparsity restriction and nonlocal restriction. Though (2) is inspired by the CSR model [3], the 2D wavelet transform T 2D in (2) will bring some other advantages such as efficiency and also effectiveness when compared with the iteratively implemented PCA transform in the CSR model. More importantly, there exist numerous wavelet shrinkage methods, such as bivariate shrinkage [8], Bayesian shrinkage [9], Gaussian scale mixtures (GSM) [10], etc., which might be integrated into this nonlocal shrinkage model. In the below, additive Gaussian noise environment is supposed and three nonlocal shrinkage functions are derived based on different selections for the norms Lp and Lq. Note that we adopt the commonly used L1 norm and L2 norm for simplicity, while nonlocal shrinkage functions based on other norms and shrinkage functions can also be obtained similarly.
3.1. L1 and L2 shrinkage function We first set p¼ 1 and q¼ 2, which leads to the following expression: 1 _ α ¼ arg min 12 jjZ T 2D αjj22 þ λ1 jjαjj1 þλ2 jjα βjj2 : α
ð3Þ
In (3), L1 norm is chosen for the term jjαjjp , which represents local sparsity of wavelet coefficients within the same band and same block. It implies that α is characterized by i.i.d. Laplacian distribution and the classical soft shrinkage function will be adopted. The term jjα βjjq models centralization prior of the coefficients across different blocks in the same group. Let A ¼ fαj;k ; j ¼ 1; :::; mg be the coefficient set, where j denotes the block index and k denotes position index within the block, the value of q will depend on the coefficient distribution in Ai and consequently on the effect of block matching. Usually, the underlying clean image blocks in the group are not the same, and the matching error in the transform domain can be defined as E ¼ fαj;k β; j ¼ 1; :::; mg. Here we set q¼ 2 in (3) with the assumption that the matching error α β follows a Gaussian distribution. Similar with the derivation of the wavelet soft shrinkage function [5], one needs to translate (3) into its scalar version gðtÞ ¼ 12 ðt t 0 Þ2 þτ1 jtj þ τ2 ðt bÞ2 ;
ð4Þ
Please cite this article as: H. Zhong, et al., Modified BM3D algorithm for image denoising using nonlocal centralization prior, Signal Processing (2014), http://dx.doi.org/10.1016/j.sigpro.2014.08.014i
H. Zhong et al. / Signal Processing ] (]]]]) ]]]–]]]
1 3 5 7 9 11 13 15 17
where τ1 ¼ λ1 =c and τ2 ¼ λ2 =c are scaled relaxation parameters by the constant c, t0 and b are the scalar components of T 2D Z and β, respectively. Then, the minimization of (4) leads to the nonlocal shrinkage function, namely L1 þL2 shrinkage function, as 8 ατ þτ β 1 2 α 4 τ1 τ2 β > > < 1 þ τ2 0 else ð5Þ α^ ¼ Sτ1 ;τ2 ;β ðαÞ ¼ > > : α þ τ1 þ τ2 β α o τ1 τ2 β 1 þ τ2 where pffiffiffi Sτ1 ;τ2 ;β ð:Þ denotes the shrinkage function, τ1 ¼ c1 2σ 2n =σ, τ2 ¼ c2 σ 2n =δ2 , c1 and c2 are constant, σ 2n is the variance of additive noise, σ and δ are the standard deviations of the sets A and E, respectively. We can find that (5) (and the followed double L1 shrinkage function) is more general and wavelet soft shrinkage function is one case of (5) by setting both τ2 and β as zero.
19 3.2. Double L1 shrinkage function 21 23 25 27 29
The CSR model suggests to use double L1 norms to solve the same problem in sparse coding, where the distributions of the code of a block and also across similar blocks are assumed Laplacian, respectively. Similarly, using wavelet transform instead of the codebook, double l1-norms can also be applied and leads to the following optimal problem: 1 _ α ¼ arg min 12 jjZ T 2D αjj22 þ λ1 jjαjj1 þλ2 jjα βjj1 : α
31 33 35 37 39 41 43 45 47 49
ð6Þ
where both p and q are set as one. In fact, when block matching is fairly good, blocks in the same group will have little difference so that the empirical distribution of the matching error α β can be well modeled by the Laplacian distribution. The double L1 shrinkage function can be derived similarly, which has been given in the CSR model [3] ( Sτ1 ;τ2 ;β ðαÞ β Z0 α^ ¼ ð7Þ Sτ1 ;τ2 ; β ðαÞ β o0 pffiffiffi pffiffiffi where τ1 ¼ c1 2σ 2n =σ, τ2 ¼ c2 2σ 2n =δ, and the shrinkage function Sτ1 ;τ2 ;β ð:Þ has the form [3] 8 αþ τ1 þ τ2 α o τ1 τ2 > > > > > τ1 τ2 r α rτ1 τ2 > <0 ð8Þ Sτ1 ;τ2 ;β ðαÞ ¼ α τ1 þ τ2 τ1 τ2 r αr τ1 τ2 þ β : > > > τ þ β r α rτ þτ þ β β τ > 1 2 1 2 > > : α τ τ α 4τ þ τ þ β 1
2
1
2
51 53
3
norm: 1 α^ ¼ arg min 12 jjZ T 2D αjj22 þ λ2 jjα βjj1 : α
ð9Þ
For (9), the following shrinkage function can be derived: α^ ¼ Sτ2 ðα βÞ þ β; ð10Þ pffiffiffi 2 where τ2 ¼ c2 2σ n =δ, c2 is constant, and Sτ2 ð:Þis the soft shrinkage function but works in the way of intra-block ( 0 jαj r τ2 Sτ2 ðαÞ ¼ ð11Þ α sgnðαÞτ2 jαj 4 τ2
3.4. Implementation details The proposed algorithm differs from the BM3D algorithm only at its first stage. The reference blocks will be sampled from the image per 3 pixels horizontally and vertically, and bior1.5 wavelet is adopted to transform each block. Then, for each group we replace the original shrinkage function with the above nonlocal shrinkage functions. The parameter β is defined as the nonlocal means [1] with its kth component calculated as , m
βk ¼ ∑ ωj αj;k j¼1
m
∑ ωj ;
ð12Þ
j¼1
where αj;k is the kth wavelet coefficient of block Z j in the group, the weight ωj ¼ expð jjZ Z j jj22 =N2 hÞ, the patch size N ¼8 and the smoothing factor h ¼ 12σ n . The size of search neighborhood for similar blocks is 39 39, and the group size m is set as 16 when σ n r30 and 20 when σ n 430. Two factors c1 and c2 are set as 0.1 and 0.9, respectively, following the way in the CSR model [3]. After the shrinkage and inverse transform, block-wised estimates are obtained and then aggregated to get the estimate of the whole image, as done in [2]. Notice that the first stage needs to be iterated once (σ n r50) or twice (σ n 450) again to get the refined estimate using the technique of residual feedback [3]. That is, we first compute the image y ¼ y þδðy y'Þ, where y and y0 are the noisy observation and the estimate of last iteration, respectively, and the factor δ¼0.2. Then, use y as input and repeat the first stage again with the newly pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi estimated noise standard deviation σ^ n ¼ λ σ n 2 δn 2 , where λ is set as 0.23 and δ2n is the variance of the difference image y y0 . The final output of the first stage is considered as the basic image estimate. The second stage is the same as in BM3D [2], where the basic image estimate and the noisy image are as the input, and wiener filtering will be performed to get the final estimate.
3.3. Nonlocal L1 shrinkage function 55 57 59 61
In practice, jj! α jj1 is used to guarantee the sparsity of the wavelet coefficients in each block. However, an interesting fact is that this local sparsity restriction can sometimes be neglected due to the great power of nonlocal term jjα βjjq . Consequently, we remove the regularization term jj! α jj1 from (2) and obtain the following nonlocal L1
4. Results 4.1. Comparison of three nonlocal shrinkage functions The purpose of this part is to compare the three nonlocal shrinkage functions whose local and nonlocal norms
Please cite this article as: H. Zhong, et al., Modified BM3D algorithm for image denoising using nonlocal centralization prior, Signal Processing (2014), http://dx.doi.org/10.1016/j.sigpro.2014.08.014i
H. Zhong et al. / Signal Processing ] (]]]]) ]]]–]]]
4
1 3 5 7 9 11 13 15 17
are different. Table 1 gives the PSNR values of denoised images with the modified BM3D algorithm. We can observe that the differences between them are very small when σ n r50, and become a little larger when σ n 4 50. For images with large smooth regions such as House and Cameraman, L1 þL2 function works slightly worse, at least 0.1–0.2 dB lower than the others. For image House, this difference can reach 0.7 dB when σ n ¼ 75. However, for images with large amount of texture such as Baboon, L1þ L2 function obtains the best results of PSNR values, which demonstrate the role of L2 norm in some aspect Q4 when representing nonlocal similarity. In comparison, double L1 function holds the best denoising power as a whole, which can be concluded from the PSNR values in Table 1. Nonlocal L1 function, which has the simplest form, achieves competitive performance with double L1 function.
19 21 23 25 27 29
Table 1 PSNR comparison of nonlocal shrinkage functions under different noise standard deviation. Function
House
L1þ L2 32.64 Double L1 32.85 Nonlocal L1 32.82
31.21 31.38 31.35
29.41 29.66 29.64
26.94 27.59 27.59
25.29 26.00 26.00
Cameraman
L1þ L2 29.30 Double L1 29.41 Nonlocal L1 29.41
27.78 27.91 27.90
26.19 26.29 26.29
24.26 24.39 24.43
22.91 23.13 23.16
Lena
L1þ L2 31.94 Double L1 32.10 Nonlocal L1 32.07
30.50 30.62 30.61
28.85 29.04 29.03
26.73 27.27 27.28
25.35 25.97 25.96
Man
L1þ L2 29.64 Double L1 29.67 Nonlocal L1 29.67
28.25 28.26 28.27
26.76 26.84 26.83
25.01 25.30 25.31
23.87 24.25 24.25
Barbara
L1þ L2 30.53 Double L1 30.78 Nonlocal L1 30.75
28.97 29.15 29.13
27.19 27.37 27.36
24.95 25.13 25.15
23.42 23.58 23.61
Baboon
L1þ L2 25.52 Double L1 25.50 Nonlocal L1 25.50
23.98 23.92 23.93
22.59 22.48 22.51
21.14 21.09 20.78
20.29 20.29 20.08
31 33 35 37 39 41 43
σ n ¼25 σ n ¼35 σ n ¼50 σ n ¼ 75 σ n ¼ 100
Image
Fig. 1 gives the visual comparison of these functions when σ n ¼ 75. It is clear that L1þL2 function has the worst smoothing power, producing many ring-like artifacts in smooth regions. Double L1 function can smooth the smooth regions very well with less ring-like artifacts, but the detail is blurred very slightly. The result of nonlocal L1 function is shown in Fig. 1(c), where the detail is similar with that of double L1 function, but appearing many snow-like artifacts. The reason may be resorted to the lack of local sparsity restriction. These observations can be used to explain why three functions behave differently in Table 1. 4.2. Comparison with other state of the art methods This part compares the denoising performance of our method with three recently developed denoising methods, including BM3D [2], CSR [3], and LSSC [4]. All the results of these benchmark methods are obtained using the source codes released by the original authors, while our method works with the double L1 function. Table 2 presents the PSNR comparison and the best PSNR value in each cell is bolded. From Table 2, we can conclude that our method is very competitive. As a whole, our method always outperforms the other methods in averaged PSNR values. When σ n Z50, the averaged PSNR gains of our methods become obvious, about 0.3 dB over the BM3D algorithm, and 0.1 dB over the CSR and LSSC algorithms. Visual inspection is also given in Figs. 2 and 3, where σ n takes a large value of 100 and 50, respectively. We can see in Fig. 2 that our method introduces fewer artifacts in the uniform regions of House and preserves the details better than other methods. Fig. 3 demonstrates that our method has competitive power in preserving the texture of Barbara. We can find that the texture in the circled area of Fig. 3(d) is slightly better than the others. 5. Conclusion In this letter, a modified BM3D algorithm using nonlocal centralization prior is presented and verified to be effective when compared with the original BM3D algorithm and other state of the art methods. More importantly, the presented
45 47 49 51 53 55 57 59 61
Fig. 1. Denoised images with different nonlocal shrinkage functions when σ n ¼ 75: (a) L1þ L2 (PSNR ¼ 24.26 dB); (b) double L1 (PSNR ¼24.39 dB) and (c) nonlocal L1 (PSNR¼ 24.43 dB).
Please cite this article as: H. Zhong, et al., Modified BM3D algorithm for image denoising using nonlocal centralization prior, Signal Processing (2014), http://dx.doi.org/10.1016/j.sigpro.2014.08.014i
H. Zhong et al. / Signal Processing ] (]]]]) ]]]–]]]
1 3
Table 2 PSNR comparison of different methods under different noise standard deviation. In each cell, top left: BM3D; top right: CSR; bottom left: LSSC; bottom right: our method. Image
σ n ¼ 25
House
32.86 33.00
32.89 32.85
31.34 31.54
31.41 31.38
29.37 29.88
29.55 29.66
27.20 27.60
27.29 27.59
25.50 25.68
25.68 26.00
Cameraman
29.45 29.47
29.43 29.41
27.84 27.91
27.75 27.91
25.85 26.31
26.18 26.29
24.06 24.33
24.21 24.39
22.82 23.09
22.88 23.13
Lena
32.08 31.85
31.97 32.10
30.50 30.50
30.54 30.62
28.86 28.84
28.94 29.04
27.02 27.07
27.16 27.27
25.57 25.73
25.88 25.97
Man
29.62 29.61
29.58 29.67
28.19 28.14
28.12 28.26
26.59 26.62
26.74 26.84
25.10 24.97
25.24 25.30
23.97 23.81
24.19 24.25
Barbara
30.72 30.35
30.69 30.78
28.81 28.76
29.08 29.15
27.17 26.90
27.18 27.37
25.10 24.90
25.06 25.13
23.49 23.45
23.62 23.58
Baboon
25.35 25.54
25.41 25.50
23.55 23.96
23.81 23.92
22.06 22.45
22.35 22.48
20.68 21.00
20.99 21.09
19.92 20.13
20.19 20.29
Average
30.01 29.97
30.00 30.05
28.37 28.47
28.45 28.54
26.65 26.83
26.82 26.95
24.86 24.98
24.99 25.13
23.55 23.65
23.74 23.87
5 7 9 11 13 15 17 19
5
σ n ¼ 35
σ n ¼ 50
σ n ¼ 75
σ n ¼100
21 23 25 27 29 31 33 35
Fig. 2. Denoised image for House whenσ n ¼ 100: (a) BM3D (PSNR¼ 25.50 dB); (b) CSR (PSNR ¼25.68 dB); (c) LSSC (PSNR ¼ 25.68 dB) and (d) our method (PSNR ¼26.00 dB).
37 39 41 43 45 47 49 51 53
Fig. 3. Denoised image for Barbara when σ n ¼ 50: (a) BM3D (PSNR¼ 27.17 dB); (b) CSR (PSNR ¼ 27.18 dB); (c) LSSC (PSNR ¼26.90 dB) and (d) our method (PSNR ¼27.37 dB).
55 57 59 61
work might stimulate more meaningful nonlocal shrinkage functions based on numerous studies in wavelet literature. In future work, proper wavelet shrinkage functions will be considered. Another attempt might be to exploit the relationship between local and nonlocal restrictions for different image types.
Acknowledgments This work was supported by the National Basic Q5 Research Program (973 Program) of China (No. 2013 CB329402), the National Natural Science Foundation of China (Nos. 61072106, 61271302 and 61271298), the Fund
Please cite this article as: H. Zhong, et al., Modified BM3D algorithm for image denoising using nonlocal centralization prior, Signal Processing (2014), http://dx.doi.org/10.1016/j.sigpro.2014.08.014i
H. Zhong et al. / Signal Processing ] (]]]]) ]]]–]]]
6
1 3
for Foreign Scholars in University Research and Teaching Programs (the 111 Project) (Nos. B07048), and the Fundamental Research Funds for the Central Universities (No. K5051302004).
5 References 7 9 11 13 15
[1] A. Buads, B. Coll, and J.M. Morel, A non-local algorithm for image denoising, in: Proceedings of IEEE CVPR, 2005, pp. 60–65. [2] K. Dabov, A. Foi, V. Katkovnik, K. Egiazarian, Image denoising by sparse 3D transform-domain collaborative filtering, IEEE Trans. Image Process. 16 (8) (2007) 2080–2095. [3] W. Dong, L. Zhang, G. Shi, Centralized sparse representation for image restoration, in: Proceedings of IEEE ICCV, 2011, pp. 1259–1266. [4] J. Mairal, F. Bach, J. Ponce, G. Sapiro, A. Zisserman, Non-local sparse models for image restoration, in: Proceedings of IEEE ICCV, 2009, pp. 2272–2279.
[5] D.L. Donoho, I.M. Johnstone, Idea spatial adaptation via wavelet shrinkage, Biometrika 81 (1994) 425–455. [6] V. Katkovnik, A. Foi, K. Egiazarian, J. Astola, From local kernel to nonlocal multiple-model image denoising, Int. J. Comput. Vis. 86 (1) (2010) 1–32. [7] W. Dong, L. Zhang, G. Shi, Nonlocally centralized sparse representation for image restoration, IEEE Trans. Image Process. 22 (4) (2013) 1618–1628. [8] L. Sendur, I.W. Selesnick, Bivariate shrinkage functions for waveletbased denoising exploiting inter-scale dependency, IEEE Trans. Signal Process. 50 (11) (2002) 2744–2756. [9] H.A. Chipman, E.D. Kolaczyk, R.E. McCulloch, Adaptive Bayesian wavelet shrinkage, J. Am. Stat. Assoc. 92 (440) (1997) 1413–1421. [10] J. Portilla, V. Strela, M.J. Wainwright, E.P. Simoncelli, Image denoising using scale mixture of gaussians in the wavelet domain, IEEE Trans. Image Process. 12 (11) (2003) 1338–1351.
17
Please cite this article as: H. Zhong, et al., Modified BM3D algorithm for image denoising using nonlocal centralization prior, Signal Processing (2014), http://dx.doi.org/10.1016/j.sigpro.2014.08.014i