Документ взят из кэша поисковой машины. Адрес оригинального документа : http://imaging.cs.msu.ru/pub/2013.PRIA.Mamaev_Lukin_Yurin.HermiteNLM.en.pdf
Дата изменения: Mon Nov 18 12:04:22 2013
Дата индексирования: Thu Feb 27 20:27:57 2014
Кодировка:
HERMITE FUNCTIONS EXPANSION BASED NON-LOCAL MEANS ALGORITHM FOR CT-APPLICATIONS1
N. Mamaev2, A. Lukin3, D. Yurin4, M. Glazkova5, V. Sinitsin6 Laboratory of Mathematical Methods of Image Processing Faculty of Computational Mathematics and Cybernetics Lomonosov Moscow State University Leninskie Gory, Moscow 119991, Russia mamaev.nikolay93@mail.ru, 3lukin@ixbt.com, 4yurin@cs.msu.ru 5,6 Federal Center of Medicine and Rehabilitation Ivan'kovskoye sh., 3, Moscow 125367, Russia 5 mary-ga@yandex.ru, 6vsini@mail.ru
2, 3, 4

2

We propose a new algorithm for noise suppression ­ HeNLM. It is a modification of an existing local jets based method LJNLM-LR. In our method, the neighborhood of a pixel is expanded in a series of Hermite functions, which, unlike Gaussian derivatives, form an orthonormal basis. Testing of our method against Local Jets and NLM has shown better preservation of detail-rich areas, sharp edges in CT scans, without loss of medically-important information.

Introduction Computed tomography (CT) is method of recovery slice images of a patient's body based on Radon transform [5]. Noise is always present in CT images and becomes lager when the power of X-ray radiation is reduced. Therefore, effective noise reduction in obtained CT images allows reducing the patient's radiation dose [6]. Noise in CT-images is close to Gaussian [3]. One class of efficient image denoising algorithms computes weights of pixels depending on similarity of the values characterizing the pixel's patch [2], [4]. First of such class ­ a Non-Local Means (NLM) algorithm ­ is calculating weights depending on the Euclidian distance between patches. Algorithm [4] is a modification of NLM that allows reducing complexity by replacing Euclidian distance with a distance between feature vectors based on Taylor series expansion. The cardinality of feature vectors is much less than a number of pixels in a patch. In this article, a method of non-local filtering based on patch decomposition with Hermite functions (HeNLM) is proposed. HeNLM method allows for better distinguishing of im1

age textures because orthonormality of Hermite functions leads to higher independence of feature vector components. Also HeNLM method better describes high-frequency components of the local neighborhood. Local Jets based Non-Local Means

In NLM algorithm, each output pixel f ( p) , p ( x, y) is produced as a weighed sum of source pixels I ( p) from the neighborhood Q : 1 f ( p) w( p, ) I ( p ) (1) w Q Q Weights w( p, ) depend on similarity of patches (pixel blocks) v(p) around the pixel: v (p, ) v (p ) 2 2 w(p, ) exp (2) 2 2 NLM algorithm provides high quality of filtered image but has a high computational complexity. Furthermore, NLM does not consider possible rotation of pixel patches. For example, weights of pixels lying on the same curved edge, but having different gradient orientation, will be small. These can lead to poor noise reduction along curved or twisting edges. The research was supported by the RFBR grant 13-07-00584


LJNLM-LR method [4] partly overcomes the shortcomings of NLM. In LJNLM-LR the distance between patches in (2) is replaced by a distance between vectors characterizing these patches. Components of the vector are convolutions of a source image I ( x, y ) with derivatives of a Gaussian function:
n m d n m G f nm ( x, y ) I ( x, y ) * 1 n m dx n dy m

(3)

The Gaussian function is defined as 2 2 2 1 G ( x, y ) e ( x y ) / 2 (4) 2 2 The multiplier n m in (3) is introduced for producing equal filter response to the same texture at different scales [7]. The denominator 1 n m reflects the number of derivatives of order n m and balances weight contribution from derivatives of different orders. Obtained features are converted to the coordinate system ( g , ) where g is a brightness gradient and is a tangent in the current pixel ( x, y ) , for invariance to rotation. Finally a vector characterizing the pixel patch is defined as follows: ~ v f nm ; n m r , S , (5) ~ where f nm is a feature in the coordinate system ( g , ) , r is the maximal used derivative order, S is set of different scales.

(8) n The multiplier 1 / is introduced for equal filter response at different scales (see a note after formula (4)). Several first Hermite functions and derivatives of a Gaussian function are shown in Fig. 1. Unlike Gaussian derivatives, Hermite functions have more uniform amplitude of oscillations throughout their support interval.

n ( x)

1

x

Fig. 1. Hermite functions (left) and Gaussian derivatives (right).

We define 2D Hermite functions as: nm ( x, y ) n ( x ) m ( y )

(10)





Hermite Functions based Non-Local Means

Hermite Functions Hermite functions are defined as follows:
x2 2 x n
2

1 n ( x) e cn

d n (e dx
x2 2

)

(1) n H n x e cn where cn
n

(6)

2 n!
n
2

d ne x x 2 and H n ( x) (1) e is a dx n Hermite polynomial [1]. They form an orthonormal system L2 (,) [1]:


in



n ( x) m ( x)dx

nm

(7)

Multiscale Hermite functions are defined as follows:

We propose a modification of LJNLM-LR method by using Hermite functions for expansion instead of Gaussian derivatives. Now elements of the feature vector are convolutions of the source image I ( x, y) with Hermite functions: hnm ( x, y ) I ( x, y ) * nm ( x, y ) (11) Similarly to LJNLM-LR method, obtained features are converted in the coordinate system ( g , ) . Consider the rotation of the coordinate system: x cos sin R , R (12) y sin cos The differential operator in a new coordinate system can be expressed as a linear combination of operators in the old coordinate system: d d Rij , where d i j 1, 2 dx j (13) i 1,2; 1 ; 2 ; j 1,2; x1 x; x2 y. A Hermite function is a product of a Gaussian derivative, a radially symmetric multiplier


e( x y ) / 2 , and a coefficient cn . Therefore Hermite function in a new coordinate system can be written as:
cn cm
aj
n ,m

2

2



nm


j 0

a jc jc
j k

n m j




j ,n m j

min{ j , n}

k max{0 , j m}

(

1)

k Cn C

j k m

cos



m j 2 k

sin n

j 2 k



(14)
Fig. 2. Original CT image. White border shows the fragment enlarged in Fig. 3.

n! k!(n k )! We note that the expression for rotation of Gaussian derivatives can be obtained from (14) by replacing n , m with f nm ( x, y ) and setk Cn

ting coefficients ci 1 . By choosing a new basis ( , ) so that matches the gradient of the image T g I ( x , y ) * dG dx , I ( x , y ) * dG dy , we T obtain g g cos , sin in (12) and (14). The use of Hermite functions instead of Gaussian derivatives leads to better characterization of the pixel neighborhood because orthogonality of Hermite functions results in lower de pendence of feature vector components hnm and, accordingly, their greater significance. Also the effective localization region is expanded (see Fig. 1) and, accordingly, the peripheral area of a local neighborhood is better represented in a feature vector.
Results

Fig. 3. Left to right and top to bottom: original CT image and result of filtering with different methods: NLM, LJNLM-LR, and HeNLM. PSNR (in decibels) between noisy and denoised image fragments is shown.

Fig. 2 shows the original CT image. Fig. 3 shows its enlarged fragment and the results of filtering with different methods: NLM, LJNLM-LR and HeNLM. The size of neighborhood for all methods has been set to 21в21 pixels; the size of patch for NLM has been set to 7в7 pixels. The maximal order of Hermite functions and Gaussian derivatives has been set to 4; the scale for Gaussian functions has been set to 1.0 and 4.0. Thus, the cardinality of a feature vector in LJNLM-LR and HeNLM is 30. The parameter has been manually selected for each method so that PSNR between noisy and filtered images (i.e. filtering strength) is the same for all compared methods.

Fig. 3 shows that LJNLM-LR and HeNLM methods perform filtering stronger than NLM in fragments with small details. Also HeNLM preserves edges better than LJNLM-LR. Another comparison of the methods was performed with phantom images described in [6] (results are available at our laboratory website: http://imaging.cs.msu.ru/en/research/ct_fil tering/pria2013add). Fig. 4 shows one of these images. A clean image was contaminated by a Gaussian noise with 3000 and then filtered with different algorithms. The parameter was manually chosen so that PSNR between original and filtered images is maximal for each algorithm separately. Table 1 shows some results for phantom images.


Fig. 4. Example of phantom image used for comparison of the methods. White borders show fragments in which PSNR was computed.

Fig. 5. Measurement of the contrast level (HU), image noise, and signal-to-noise ratio (SNR) at the root of ascending aorta. The region of interest (ROI) was defined as 100 mm2.

Table Fragment NLM LJNLM HeNLM

1. Results on a phantom image Over1 2 3 4 all 31.98 30.00 42.25 49.03 40.86 34.12 31.52 43.15 46.73 41.90 34.17 31.80 43.09 47.41 42.10

Conclusion

The proposed HeNLM method provides an improved image quality without the loss of medically important diagnostic information.
Acknowledgements

The table shows that LJNLM-LR and HeNLM perform filtering better than NLM on the whole image (Table 1, column "Overall"). HeNLM performs filtering slightly better than LJNLM-LR in fragments where there are many small details. However NLM performs filtering better in flat areas. The image quality was assessed by two independent radiologists. They measured contrast level, image noise and signal-to-noise ratio (SNR) of ascending aorta, myocardium, left ventricular cavity and pulmonary trunk (Fig. 5). The quality of images was also evaluated subjectively (visualization of small coronary arteries and coronary stents). Contrast level, image noise, SNR inside ascending aorta, myocardium left ventricular cavity and pulmonary trunk were significantly better after filtering, in comparison with the original CT data (Friedman test p-value p=0.001). Visualization of small arteries and coronary stents was not statistically different (p>0.05) between the original CT data and the result of filtering with different methods.

The research was supported by the RFBR grant 13-07-00584.
References
1. M. Abramowitz, I.A. Stegun. "Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables" // U.S. National Bureau of Standards (1964); Dover, New York 1965. 2. A. Buades. "A non-local algorithm for image denoising" // IEEE Computer Society Conference on Computer Vision and Pattern Recognition, 2005, Volume 2, Pages 60­65. 3. Kijewski M.F., Judy P.F. "The noise power spectrum of CT images" // Phys. Med. Biol., 1987, Vol. 32, No 5, pp. 565­575. Printed in the UK. 4. A. Manzanera. "Local Jet based similarity for NLMeans filtering" // Pattern Recognition (ICPR), 20th International Conference on Computer Vision and Pattern Recognition, 2010, Pages 2668­2671. 5. F. Natterer. "The Mathematics of Computerized Tomography" // SIAM: Society for Industrial and Applied Mathematics, Philadelphia, PA, USA, July 2001, 222 p. 6. M. V. Storozhilova, A. S. Lukin, D. V. Yurin, V. E. Sinitsyn. "Two approaches for noise filtering in 3D medical CT-images" // Proc. of 22nd Int. Conf. on Computer Graphics GraphiCon'2012, Moscow, Russia, 2012, pp. 68-72. 7. T. Lindeberg. "Scale­Space Theory in Computer Vision" // Kluwer Academic Publishers, Dordrecht, 1994.