- About this Journal ·
- Abstracting and Indexing ·
- Advance Access ·
- Aims and Scope ·
- Annual Issues ·
- Article Processing Charges ·
- Articles in Press ·
- Author Guidelines ·
- Bibliographic Information ·
- Citations to this Journal ·
- Contact Information ·
- Editorial Board ·
- Editorial Workflow ·
- Free eTOC Alerts ·
- Publication Ethics ·
- Reviewers Acknowledgment ·
- Submit a Manuscript ·
- Subscription Information ·
- Table of Contents

Computational and Mathematical Methods in Medicine

Volume 2013 (2013), Article ID 902143, 8 pages

http://dx.doi.org/10.1155/2013/902143

## Improving Spatial Adaptivity of Nonlocal Means in Low-Dosed CT Imaging Using Pointwise Fractal Dimension

^{1}College of Computer Science, Sichuan University, No. 29 Jiuyanqiao Wangjiang Road, Chengdu 610064, Sichuan, China^{2}School of Computer Science, Sichuan Normal University, No. 1819 Section 2 of Chenglong Road, Chengdu 610101, Sichuan, China^{3}School of Automation Engineering, University of Electronic Science and Technology of China, No. 2006, Xiyuan Ave, West Hi-Tech Zone, Chengdu 611731, Sichuan, China^{4}School of Information Science and Technology, East China Normal University, No. 500, Dong-Chuan Road, Shanghai 200241, China

Received 25 January 2013; Accepted 6 March 2013

Academic Editor: Shengyong Chen

Copyright © 2013 Xiuqing Zheng et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

NLMs is a state-of-art image denoising method; however, it sometimes oversmoothes anatomical features in low-dose CT (LDCT) imaging. In this paper, we propose a simple way to improve the spatial adaptivity (SA) of NLMs using pointwise fractal dimension (PWFD). Unlike existing fractal image dimensions that are computed on the whole images or blocks of images, the new PWFD, named pointwise box-counting dimension (PWBCD), is computed for each image pixel. PWBCD uses a fixed size local window centered at the considered image pixel to fit the different local structures of images. Then based on PWBCD, a new method that uses PWBCD to improve SA of NLMs directly is proposed. That is, PWBCD is combined with the weight of the difference between local comparison windows for NLMs. Smoothing results for test images and real sinograms show that PWBCD-NLMs with well-chosen parameters can preserve anatomical features better while suppressing the noises efficiently. In addition, PWBCD-NLMs also has better performance both in visual quality and peak signal to noise ratio (PSNR) than NLMs in LDCT imaging.

#### 1. Introduction

Radiation exposure and associated risk of cancer for patients from CT examination have been increasing concerns in recent years. Thus minimizing the radiation exposure to patients has been one of the major efforts in modern clinical X-ray CT radiology [1–8]. However, the presentation of serious noise and many artifacts degrades the quality of low-dose CT images dramatically and decreases the accuracy of diagnosis dose. Although many strategies have been proposed to reduce their noise and artifacts [9–14], filtering noise from clinical scans is still a challenging task, since these scans contain artifacts and consist of many structures with different shape, size, and contrast, which should be preserved for making correct diagnosis.

Recently nonlocal means (NLMs) is proposed for improving the performance of classical adaptive denoising methods [15–17] and shows good performance even in low-dose CT (LDCT) imaging [18–20].

There are two novel ideas for NLMs. One is that the similar points should be found by comparing the difference between their local neighborhoods instead of by comparing their gray levels directly. Since gray levels of LDCT will be polluted seriously by noises and artifacts, finding similar points by local neighborhoods instead of by gray levels directly will help NLMs find correct similar points. The other important idea for NLMs is that the similar points should be searched in large windows to guarantee the reliability of estimation.

Following the previous discussion, the NLMs denoising should be performed in two windows: one is comparison patch and the other is searching window. The sizes of these two windows and the standard deviation of the Gaussian kernel, which is used for computing the distance between two neighborhoods, should be determined according to the standard deviation of noises [15–17], and these three parameters are identical in an image.

Some researchers find that identical sizes of two windows and identical Gaussian kernel in an image are not the best choice for image denoising [21–25]. The straightest motivation is that the parameters should be modified according to the different local structures of images. For example, the parameters near an edge should be different from parameters in a large smooth region.

An important work to improve the performance of NLMs is quasi-local means (QLMs) proposed by us [21, 22]. We argue that nonlocal searching windows are not necessary for most of image pixels. In fact, for points in smooth regions, which are the majority of image pixels, local searching windows are big enough, while for points near singularities, only the minority of image pixels, nonlocal search windows are necessary. Thus the method is named quasi-local whereit islocal for most of image pixels and nonlocal only for pixels near singularities. The searching windows for quasi-local means (QLMs) are variable for different local structures, and QLMs can get better singularity preservation in image denoising than classical NLMs.

Other important works about improving spatial adaptivity of NLMs are proposed very recently [23–25]. The starting point for these works is that the image pixels are parted into different groups using supervised learning or semisupervised learning and clustering. However, the learning and clustering will waste a lot of computation time and resource, which will hamper them to be applied in medical imaging. Thus we must propose a new method for improving the spatial adaptivity with a simple way.

In this paper we propose a simple and powerful method to improve spatial adaptivity for NLMs in LDCT imaging using pointwise fractal dimension (PWFD) where PWFD is computed pixel by pixel in a fixed-size window centered at the considering pixel. According to the new definition of PWFD, different local structures will be with different local fractal dimensions, for example, pixels near edge regions will be with relatively big PWFDs, while PWFDs of pixels in smooth regions will be zeros. Thus PWFD can provide local structure information for image denoising. After defined PWFD, which can fit different local structures of images well, we design a new weight function by combining the new PWFD difference between two considering pixels with the weight of original NLMs measured by gray level difference between two comparison windows. Thus using this new weight function, the proposed method will not only preserve the gray level adaptivity of NLMs but also improve the SA of NLMs.

The arrangement of this paper is as follows: In Section 2, the backgrounds are introduced, then the new proposed method is presented in Section 3, the experiment results are shown and discussed in Section 4, and the final part is the conclusions and acknowledgment.

#### 2. Backgrounds

In this section, we will introduce related backgrounds of the proposed method.

##### 2.1. Noise Models

Based on repeated phantom experiments, low-mA (or low-dose) CT calibrated projection data after logarithm transform were found to follow approximately a Gaussian distribution with an analytical formula between the sample mean and sample variance; that is, the noise is a signal-dependent Gaussian distribution [11].

The photon noise is due to the limited number of photons collected by the detector. For a given attenuating path in the imaged subject, and denote the incident and the penetrated photon numbers, respectively. Here, denotes the index of detector channel or bin and is the index of projection angle. In the presence of noises, the sinogram should be considered as a random process and the attenuating path is given by where is a constant and is Poisson distribution with mean .

Thus we have

Both its mean value and variance are .

Gaussian distributions of ployenergetic systems were assumed based on limited theorem for high-flux levels and followed many repeated experiments in [11]. We have where is the mean and is the variance of the projection data at detector channel or bin is a scaling parameter, and is a parameter adaptive to different detector bins.

The most common conclusion for the relation between Poisson distribution and Gaussian distribution is that the photon count will obey Gaussian distribution for the case with large incident intensity and Poisson distribution with feeble intensity [11].

##### 2.2. Nonlocal Means (NLMs)

Given a discrete noisy image , the estimated value (), for a pixel , is computed as a weighted nonlocal average: where indicates a neighborhood centered at and size , called searching window, and . The family of weights depend on the similarity between the pixels and and satisfy and .

The similarity between two pixels and , depends on the similarity of the intensity gray level vectors and , where denotes a square window with fixed size and centered at a pixel , named comparison patch: and the weights are computed as where denotes the standard deviation of the noise and is the filtering parameter set depending on the value .

##### 2.3. Box-Counting Dimension

Box-counting dimension, also known as Minkowski dimension or Minkowski-Bouligand dimension, is a way of determining the fractal dimension of a set in a Euclidean space or more generally in a metric space . To calculate this dimension for a fractal , putting this fractal on an evenlyspaced grid and count how many boxes are required to cover the set. The box-counting dimension is calculated by seeing how this number changes as we make the grid finer by applying a box-counting algorithm.

Suppose that is the number of boxes of side length required to cover the set. Then the box-counting dimension is defined as

Given an image whose gray level is G, then the image is part into the grids, which are related to cube grids. If for the th grid, the greatest gray level is in the th box and the smallest is in the th box, then the box number for covering the grid is Therefore the box number for covering the whole image is Selecting different scale , we can get related . Thus we have a group of pairs . The group can be fit with a line using least-squares fitting, the slope of the line is the box-counting dimension.

#### 3. The New Method

In this section, we will present our new proposed algorithm in detail. The motivation for the proposed method is that SA of NLMs should be improved in a simpler way. The new PWFD is introduced firstly to adapt complex image local structures, and then the new weight functions based on PWFD are discussed. At the end of this section, the procedures of the proposed method are shown.

##### 3.1. Pointwise Box-Counting Dimension

In image processing, the fractal dimension usually is used for characterizing roughness and self-similarity of images. However, most of works only focus on how to compute fractal dimensions for images or blocks of images [26–30]. Since fractal dimension can characterize roughness and self-similarity of images, it also can be used for characterizing the local structures of images by generalizing it to PWFD, which is computed pixel by pixel using a fixed-size window centered in the considered pixel. Thus, each pixel in an image has a PWFD and it equals the fractal dimension of the fixed-size window centered in the considered pixel.

Following the previous discussion, the pointwise box-counting dimension (PWBCD) starts from replacing each pixel to a fixed-size window centered at . It is obvious that PWFD can be generalized to all definitions of fractal dimensions. However, in order to make our explanation more clearly, we only extend the new definition to PWBCD.

According to the new PWFD, PWBCD should be computed for each pixel in the image. For each pixel , the PWBCD is computed in a fixed-size window centered at .

The window is parted into the grids, which are related to cube grids. If for the th grid, the greatest gray level is in the th box and the smallest is in the th box, then the box number for covering the grid is Therefore the box number for covering the whole window is Selecting different scale , we can get related . Thus we have a group of pairs . The group can be fit with a line using least-squares fitting; the slope of the line is the box-counting dimension.

Note that each pixel in an image has a PWBCD value. Thus we can test the rationality for PWBCD by showing PWBCD values using an image. In these PWBCD images, high PWBCD values are shown as white points, while low PWBCD values are shown as gray or black points. If PWBCD images are similar to the original images with big PWBCD values near singularities and small PWBCD values in smooth regions, the rationality is testified.

Figure 1 shows PWBCD images for three images: an test image composed by some blocks with different gray levels, a LDCT image, and barbara. The white points signify the pixels with big fractal dimensions, while black points signify the pixels with small fractal dimensions. Here, and . Note that the white parts correspond the texture parts of barbara and soft tissues of the second image in the first row. Moreover, the PWBCD images are very similar to the original images which demonstrate that the PWBCD can be used for characterizing the local structure of images.

##### 3.2. The New Weight Function

After defining the PWBCD, we must find an efficient and powerful way to use the PWBCD in NLMs directly. Just as discussed in the previous subsection, PWBCD can characterize the local structures for images well. Thus PWBCD should be used to weight the points in the searching patch. That is, (6) should be changed as where () is FDBCD value for the considering pixel and is computed according to the method proposed in Section 3.1, denotes the standard deviation of the noise, are the filtering parameters. is the similarity between two pixels and depending on the similarity of the intensity gray level vectors and , where denotes a square window with fixed size and centered at a pixel :

Given a discrete noisy image , the estimated value (), for a pixel is computed as a weighted nonlocal average: where indicates a neighborhood centered at and size , called searching window, and . Note that the family of weights depend on the similarity between the pixels and and satisfy and .

##### 3.3. The Steps of the New Method

The steps of PWBCD-NLMs are as follows.(1)Compute pointwise box-counting dimension for each ofthe pixels. For each of the pixels,given and , compute PWBCD according to Section 3.1, and get a matrix with the same size as the image.(2)Compute weights. determine parameters: , the size of comparison window , and the size of the searching patch . Compute the difference between two comparison windows, , using (13).Compute the weights using (12).(3)Estimate real gray levels: estimate real levels using (14).

#### 4. Experiments and Discussion

The main objective for smoothing LDCT images is to delete the noise while preserving anatomy features for the images.

In order to show the performance of PWBCD-NLMs, a 2-dimensional test phantom is shown in Figure 1(a). The number of bins per view is 888 with 984 views evenly spanned on a circular orbit of . The detector arrays are on an arc concentric to the X-ray source with a distance of 949.075 mm. The distance from the rotation center to the X-ray source is 541 mm. The detector cell spacing is 1.0239 mm.

The LDCT projection data (sinogram) is simulated by adding Gaussian-dependent noise (GDN) whose analytic form between its mean and variance has been shown in (3) with and . The projection data is reconstructed by standard Filtered Back Projection (FBP). Since both the original projection data and sinogram have been provided, the evaluation is based on peak signal to noise ration (PSNR) between the ideal reconstructed image and reconstructed image.

The PWBCDs for images are computed according to Section 3.1, and the parameters are and . The new proposed method is compared with NLMs, and their common parameters includes the standard deviation of noise ; the size of comparison window is , while the size of searching patch is . The other parameter for NLMswhickis the Gaussian kernel for weights defined on (13) is and the parameters for the new method are the sizes of Gaussian kernel for two weights defined on (12): for the weights of difference between comparison window and for the weights between two PWBCDs. All parameters are chosen by hand with many experiments, which has the best performance.

Table 1 summarized PSNR between the ideal reconstructed image and filtered reconstructed image. The PWBCD-NLMs has better performance in different noise levels in the term of PSNR than NLMs.

Figure 2 shows noisy test images and their reconstructed images using NLMs and the proposed method. Although the reconstructed images are very similar to each other, the reconstructed images using the new method also show better performance in edge preservation especially in weak and curve edge preserving than the NLMs. Since PWBCD-NLMs provides a more flexible way for handling different local image structures, it has much good performance in denoising while preserving structures.

One abdominal CT images of a 62-year-old woman were scanned from a 16 multidetector row CT unit (Somatom Sensation 16; Siemens Medical Solutions) using 120 kVp and 5 mm slice thickness. Other remaining scanning parameters are gantry rotation time, 0.5 second; detector configuration (number of detector rows section thickness), mm; table feed per gantry rotation, 24 mm; pitch, 1 : 1; and reconstruction method, Filtered Back Projection (FBP) algorithm with the soft-tissue convolution kernel “B30f”. Different CT doses were controlled by using two different fixed tube currents 60 mAs for LDCT and 150 mAs (60 mA or 300 mAs) for SDCT, resp.). The CT dose index volumes (CTDIvol) for LDCT images and SDCT images are in positive linear correlation to the tube current and are calculated to be approximately ranged between 15.32 mGy and 3.16 mGy [18].

On sinogram space, the PWBCDs for images are computed according to Section 3.1 and the parameters are and . The new proposed method is compared with NLMs and their common parameters includes the standard deviation of noise ; the size of comparison window is , while the size of searching patch is . The other parameter for NLMswhichis the Gaussian kernel for weights defined on (13) is and the parameters for the new method are the sizes of Gaussian kernel for two weights defined on (12): for the weights of difference between comparison window and for the weights between two PWBCDs.

Comparing the original SDCT images and LDCT images in Figure 3, we found that the LDCT images were severely degraded by nonstationary noise and streak artifacts. In Figure 3(d), for the proposed approach, experiments obtain more smooth images. Both in Figures 3(c) and 3(d), we can observe better noise/artifacts suppression and edge preservation than the LDCT image. Especially, compared to their corresponding original SDCT images, the fine features representing the hepatic cyst were well restored by using the proposed method. We can observe that the noise grains and artifacts were significantly reduced for the NLMs and PWBCD-NLMs processed LDCT images with suitable parameters both in Figures 3(c) and 3(d). The fine anatomical/pathological features can be well preserved compared to the original SDCT images (Figure 3(a)) under standard dose conditions.

#### 5. Conclusions

In this paper, we propose a new PWBCD-NLMs method for LDCT imaging based on pointwise boxing-counting dimension and its new weight function. Since PWBCD can characterize the local structures of image well and also can be combined with NLMs easily, it provides a more flexible way to balance the noise reduction and anatomical details preservation. Smoothing results for phantoms and real sinograms show that PWBCD-NLMs with suitable parameters has good performance in visual quality and PSNR.

#### Acknowledgments

This paper is supported by the National Natural Science Foundation of China (no 60873102), Major State Basic Research Development Program (no. 2010CB732501), and Open Foundation of Visual Computing and Virtual Reality Key Laboratory Of Sichuan Province (no. J2010N03). Ming Li also acknowledges the supports by the NSFC under the Project Grant nos. 61272402, 61070214 and 60873264 and the 973 plan under the Project Grant no. 2011CB302800.

#### References

- D. J. Brenner and E. J. Hall, “Computed tomography-an increasing source of radiation exposure,”
*New England Journal of Medicine*, vol. 357, no. 22, pp. 2277–2284, 2007. View at Publisher · View at Google Scholar · View at Scopus - J. Hansen and A. G. Jurik, “Survival and radiation risk in patients obtaining more than six CT examinations during one year,”
*Acta Oncologica*, vol. 48, no. 2, pp. 302–307, 2009. View at Publisher · View at Google Scholar · View at Scopus - H. J. Brisse, J. Brenot, N. Pierrat et al., “The relevance of image quality indices for dose optimization in abdominal multi-detector row CT in children: experimental assessment with pediatric phantoms,”
*Physics in Medicine and Biology*, vol. 54, no. 7, pp. 1871–1892, 2009. View at Publisher · View at Google Scholar · View at Scopus - L. Yu, “Radiation dose reduction in computed tomography: techniques and future perspective,”
*Imaging in Medicine*, vol. 1, no. 1, pp. 65–84, 2009. View at Publisher · View at Google Scholar - J. Weidemann, G. Stamm, M. Galanski, and M. Keberle, “Comparison of the image quality of various fixed and dose modulated protocols for soft tissue neck CT on a GE Lightspeed scanner,”
*European Journal of Radiology*, vol. 69, no. 3, pp. 473–477, 2009. View at Publisher · View at Google Scholar · View at Scopus - W. Qi, J. Li, and X. Du, “Method for automatic tube current selection for obtaining a consistent image quality and dose optimization in a cardiac multidetector CT,”
*Korean Journal of Radiology*, vol. 10, no. 6, pp. 568–574, 2009. View at Publisher · View at Google Scholar · View at Scopus - A. Kuettner, B. Gehann, J. Spolnik et al., “Strategies for dose-optimized imaging in pediatric cardiac dual source CT,”
*RoFo*, vol. 181, no. 4, pp. 339–348, 2009. View at Publisher · View at Google Scholar · View at Scopus - P. Kropil, R. S. Lanzman, C. Walther et al., “Dose reduction and image quality in MDCT of the upper abdomen: potential of an adaptive post-processing filter,”
*RoFo*, vol. 182, no. 3, pp. 248–253, 2009. - H. B. Lu, X. Li, L. Li et al., “Adaptive noise reduction toward low-dose computed tomography,” in
*Proceedings of the Medical Imaging 2003: Physics of Medical Imaging, parts 1 and 2*, vol. 5030, pp. 759–766, February 2003. View at Publisher · View at Google Scholar · View at Scopus - J. C. Giraldo, Z. S. Kelm, L. S. Guimaraes et al., “Comparative study of two image space noise reduction methods for computed tomography: bilateral filter and nonlocal means,” in
*Proceedings of the Annual International Conference of the IEEE Engineering in Medicine and Biology Society*, vol. 1, pp. 3529–3532, 2009. View at Scopus - H. B. Lu, I. T. Hsiao, X. Li, and Z. Liang, “Noise properties of low-dose CT projections and noise treatment by scale transformations,” in
*Proceedings of the IEEE Nuclear Science Symposium Conference Record*, vol. 1–4, pp. 1662–1666, November 2002. View at Scopus - P. J. La Rivière, “Penalized-likelihood sinogram smoothing for low-dose CT,”
*Medical Physics*, vol. 32, no. 6, pp. 1676–1683, 2005. View at Publisher · View at Google Scholar · View at Scopus - S. Hu, Z. Liao, and W. Chen, “Reducing noises and artifacts simultaneously of low-dosed X-ray computed tomography using bilateral filter weighted by Gaussian filtered sinogram,”
*Mathematical Problems in Engineering*, vol. 2012, Article ID 138581, 14 pages, 2012. View at Publisher · View at Google Scholar - S. Hu, Z. Liao, and W. Chen, “Sinogram restoration for low-dosed X-ray computed tomography using fractional-order Perona-Malik diffusion,”
*Mathematical Problems in Engineering*, vol. 2012, Article ID 391050, 13 pages, 2012. View at Publisher · View at Google Scholar - A. Buades, B. Coll, and J. M. Morel, “A review of image denoising algorithms, with a new one,”
*Multiscale Modeling and Simulation*, vol. 4, no. 2, pp. 490–530, 2005. View at Publisher · View at Google Scholar · View at Scopus - A. Buades, B. Coll, and J. M. Morel, “A non-local algorithm for image denoising,” in
*Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition (CVPR '05)*, vol. 2, pp. 60–65, June 2005. View at Scopus - A. Buades, B. Coll, and J. M. Morel, “Nonlocal image and movie denoising,”
*International Journal of Computer Vision*, vol. 76, no. 2, pp. 123–139, 2008. View at Publisher · View at Google Scholar · View at Scopus - C. Yang, C. Wufan, Y. Xindao et al., “Improving low-dose abdominal CT images by weighted intensity averaging over large-scale neighborhoods,”
*European Journal of Radiology*, vol. 80, no. 2, pp. e42–e49, 2011. View at Publisher · View at Google Scholar - Y. Chen, Z. Yang, W. Chen, et al., “Thoracic low-dose CT image processing using an artifact suppressed largescale nonlocal means,”
*Physics in Medicine and Biology*, vol. 57, no. 9, pp. 2667–2688, 2012. View at Publisher · View at Google Scholar - Y. Chen, D. Gao, C. Nie et al., “Bayesian statistical reconstruction for low-dose X-ray computed tomography using an adaptive-weighting nonlocal prior,”
*Computerized Medical Imaging and Graphics*, vol. 33, no. 7, pp. 495–500, 2009. View at Scopus - Z. Liao, S. Hu, and W. Chen, “Determining neighborhoods of image pixels automatically for adaptive image denoising using nonlinear time series analysis,”
*Mathematical Problems in Engineering*, vol. 2010, Article ID 914564, 2010. View at Publisher · View at Google Scholar · View at Scopus - Z. Liao, S. Hu, M. Li, and W. Chen, “Noise estimation for single- slice sinogram of low-dose X-ray computed tomography using homogenous patch,”
*Mathematical Problems in Engineering*, vol. 2012, Article ID 696212, 16 pages, 2012. View at Publisher · View at Google Scholar - T. Thaipanich, B. T. Oh, P.-H. Wu, and C.-J. Kuo, “Adaptive nonlocal means algorithm for image denoising,” in
*Proceedings of the IEEE International Conference on Consumer Electronics (ICCE '10)*, 2010. - T. Thaipanich and C.-C. J. Kuo, “An adaptive nonlocal means scheme for medical image denoising,” in
*Proceedings of the SPIE Medical Imaging 2010: Image Processing*, vol. 7623, March 2010. View at Publisher · View at Google Scholar · View at Scopus - R. Yan, L. Shao, S. D. Cvetkovic, and J. Klijn, “Improved nonlocal means based on pre-classification and invariant block matching,”
*Journal of Display Technology*, vol. 8, no. 4, pp. 212–218, 2012. View at Publisher · View at Google Scholar - A. K. Bisoi and J. Mishra, “On calculation of fractal dimension of images,”
*Pattern Recognition Letters*, vol. 22, no. 6-7, pp. 631–637, 2001. View at Publisher · View at Google Scholar · View at Scopus - R. Creutzberg and E. Ivanov, “Computing fractal dimension of image segments,” in
*Proceedings of the 3rd International Conference of Computer Analysis of Images and Patterns (CAIP '89)*, 1989. - M. Ghazel, G. H. Freeman, and E. R. Vrscay, “Fractal image denoising,”
*IEEE Transactions on Image Processing*, vol. 12, no. 12, pp. 1560–1578, 2003. View at Publisher · View at Google Scholar · View at Scopus - M. Ghazel, G. H. Freeman, and E. R. Vrscay, “Fractal-wavelet image denoising revisited,”
*IEEE Transactions on Image Processing*, vol. 15, no. 9, pp. 2669–2675, 2006. View at Publisher · View at Google Scholar · View at Scopus - B. Pesquet-Popescu and J. L. Vehel, “Stochastic fractal models for image processing,”
*IEEE Signal Processing Magazine*, vol. 19, no. 5, pp. 48–62, 2002. View at Publisher · View at Google Scholar