Application of Integral Image (II): Non-local mean denoising (Nl-means)

Non-local mean denoising (Nl-means) This paper introduces the basic algorithm of Nl-means, and points out the problem of low efficiency of the algorithm, and uses the integral image technique to accelerate the algorithm.

Assuming that the image is like a vegetarian point, search window size, domain window size, calculate the similarity between the two rectangle neighborhood, for each pixel needs to calculate its similarity with the search window pixels, so nl-means complexity.

After analysis, it can be found that the algorithm can improve the similarity between neighbors, that is, time-consuming operation. In the basic algorithm, it is necessary to traverse two neighbors each time the distance between neighbors is computed, and the difference value is calculated by pixel points.

If we first construct an integral image of the pixel difference:

which

This can be done within a constant time when the distance between the two neighbors is calculated:

In this way, the complexity of the algorithm will be reduced to.

The specific algorithm description can be referenced in [1]:

In order to reduce the space complexity, the above algorithm takes the offset as the outermost loop, that is, the integral image is obtained at one offset direction at a time, and the integral image is processed. Without having to extract all the integral images at once.

Program:

Close all;clear All;clci=double(imread ('lena.tif')); I=i+ *randn (Size (i)); TicO1=nlmeans (i,2,5,); tocticO2= Fastnlmeans (I,2,5,ten); Tocfigure;imshow ([i,o1,o2],[]);

function denoisedimg=Fastnlmeans (i,ds,ds,h)%I: Image with Noise%DS: Neighborhood window Radius%Ds: Search window Radius%h: Gaussian function Smoothing Parameters%denoisedimg: denoising image I=Double(I); [M,n]=size (I); Paddedimg= Padarray (i,[ds+ds+1, ds+ds+1],'symmetric','both'); Paddedv= Padarray (I,[ds,ds],'symmetric','both'); average=zeros (m,n); Sweight=Average;wmax=AVERAGE;H2=h*h;d2=(2*ds+1)^2; fort1=-Ds:ds fort2=-Ds:dsif(t1==0&&t2==0)            Continue; End St=Integralimgsqdiff (PADDEDIMG,DS,T1,T2); V= Paddedv (1+DS+T1:END-DS+T1,1+ds+t2:end-ds+T2); W=zeros (m,n);  forI=1: M forj=1: N i1=i+ds+1; J1=j+ds+1; Dist2=st (I1+ds,j1+ds) +st (i1-ds-1, j1-ds-1)-st (i1+ds,j1-ds-1)-st (i1-ds-1, j1+DS); Dist2=dist2/D2; W (i,j)=exp (-dist2/H2); Sweight (I,J)=sweight (I,J) +W (i,j); Average (i,j)=average (i,j) +w (i,j) *V (i,j); End End Wmax=Max (WMAX,W); Endendaverage=average+wmax.*I;sweight=sweight+wmax;denoisedimg=average./sweight;function Sd=Integralimgsqdiff (PADDEDIMG,DS,T1,T2)%paddedimg: Edge-filled image%Ds: Search window Radius%(T1,T2): Offset%Sd: Integral image [M,n]=size (paddedimg); M1=m-2*ds;n1=n-2*Ds; Sd=zeros (m1,n1);D ist2= (Paddedimg (1+ds:end-ds,1+DS:END-DS)-paddedimg (1+DS+T1:END-DS+T1,1+DS+T2:END-DS+T2)). ^2; forI=1: M1 forj=1: N1ifi==1&& j==1Sd (i,j)=Dist2 (I,J); ElseIf I==1&& j~=1Sd (i,j)=SD (i,j-1)+Dist2 (I,J); ElseIf I~=1&& j==1Sd (i,j)=SD (I-1, j) +Dist2 (I,J); ElseSd (i,j)=dist2 (i,j) +sd (i-1, j) +SD (i,j-1)-SD (I-1, J-1); End EndEnd

Results:

The three images are in turn the original Nl-means algorithm denoising result and the Nl-means algorithm denoising result with integral image acceleration. For the Lena diagram of 256*256, the original algorithm takes 36.251389s, and the algorithm using Integral image acceleration is time consuming 4.647372s.

Of course, for MATLAB, if you take full advantage of its functions and matrix operations, you can further accelerate your programming:

function denoisedimg=fastNLmeans2 (i,ds,ds,h) I=Double(I); [M,n]=size (I); Paddedimg= Padarray (i,[ds+ds+1, ds+ds+1],'symmetric','both'); Paddedv= Padarray (I,[ds,ds],'symmetric','both'); average=zeros (m,n); Wmax=Average;sweight=AVERAGE;H2=h*h;d=(2*ds+1)^2; fort1=-Ds:ds fort2=-Ds:dsif(t1==0&&t2==0)            Continue; End Sd=Integralimgsqdiff (PADDEDIMG,DS,T1,T2); SqDist2=SD (2*ds+2: end-1,2*ds+2: end-1) +SD (1: end-2*ds-2,1: end-2*ds-2)...               -SD (2*ds+2: end-1,1: end-2*ds-2)-SD (1: end-2*ds-2,2*ds+2: end-1); SqDist2=sqdist2/D; W=exp (-sqdist2/H2); V= Paddedv (1+DS+T1:END-DS+T1,1+ds+t2:end-ds+T2); Average=average+w.*v; Wmax=Max (WMAX,W); Sweight=sweight+W; Endendaverage=average+wmax.*I;average=average./(wmax+sweight);D enoisedimg=average;function Sd=Integralimgsqdiff (paddedimg,ds,t1,t2) Dist2= (Paddedimg (1+ds:end-ds,1+DS:END-DS)-paddedimg (1+DS+T1:END-DS+T1,1+DS+T2:END-DS+T2)). ^2; Sd= Cumsum (Dist2,1); Sd= Cumsum (Sd,2);

It takes only 0.416442s to process the Lena diagram using the FastNLmeans2 function described above.

Reference:

[1] FROMENTJ. Parameter-free Fast pixelwise non-local Means denoising[j]. Image Processingon Line, 2014, 4:300-326

Application of Integral Image (II): Non-local mean denoising (Nl-means)