# The Nonlocal -Laplacian Evolution for Image InterpolationAcademic research paper on "Mathematics" CC BY 0 0
Share paper
Mathematical Problems in Engineering
OECD Field of science
Keywords
{""}

## Academic research paper on topic "The Nonlocal -Laplacian Evolution for Image Interpolation"

﻿Hindawi Publishing Corporation Mathematical Problems in Engineering Volume 2011, Article ID 837426,11 pages doi:10.1155/2011/837426

Research Article

The Nonlocal p-Laplacian Evolution for Image Interpolation

Yi Zhan

College of Mathematics and Statistics, Chongqing Technology and Business University, Chongqing 400067, China

Correspondence should be addressed to Yi Zhan, zhanyimeng@yahoo.com.cn Received 24 June 2011; Accepted 16 August 2011 Academic Editor: P. Liatsis

Copyright © 2011 Yi Zhan. 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.

This paper presents an image interpolation model with nonlocal p-Laplacian regularization. The nonlocal p-Laplacian regularization overcomes the drawback of the partial differential equation (PDE) proposed by Belahmidi and Guichard (2004) that image density diffuses in the directions pointed by local gradient. The grey values of images diffuse along image feature direction not gradient direction under the control of the proposed model, that is, minimal smoothing in the directions across the image features and maximal smoothing in the directions along the image features. The total regularizer combines the advantages of nonlocal p-Laplacian regularization and total variation (TV) regularization (preserving discontinuities and 1D image structures). The derived model efficiently reconstructs the real image, leading to a natural interpolation, with reduced blurring and staircase artifacts. We present experimental results that prove the potential and efficacy of the method.

1. Introduction

Digital image interpolation is an important technology in digital photography, TV, multimedia, advertising, and printing industries, which is applied to obtain higher-resolution image with better perceptual quality. The key task in image interpolation is to remove zigzagging, blurry, and other artifacts producing visually pleasing resulting image. Many literatures are devoted to address these problems [1-10]. These methods, usually known as edge directed, level set based, or isophote oriented, vary considerably.

The methods based on edge direction were proposed to obtain smooth edges of the resulting images [1-3]. However, these methods suffer from degradations on edge-free regions because they rely on local directions estimation, creating false edges in uniform regions. Wang and Ward have developed an interesting technique based on the detection of ridges (straight edges) in images , which allows them to interpolate directionally only

pixels situated on straight edges and avoid the apparition of false edges. In , the zigzagging artifacts are reduced by restricting curvature of the interpolated isophotes (equi-intensity contours). Minimum curvature is required on isophotes of the interpolated images. Cha and Kim proposed a method based on the TV energy in order to remove the so-called checkerboard effect and to form reliable edges . Morse and Schwartzwald  presented a scheme that uses existing interpolation techniques as an initial approximation and then iter-atively reconstructs the isophotes using constrained smoothing. However, they are complex compared to traditional methods and thus computationally expensive. Malgouyres and Guichard  proposed to choose as solution of the interpolation the image that minimizes the TV. This method leads to resulting images without blurring effects, as it allows discontinuities and preserves 1D fine structures. Aly and Dubois  proposed a model-based TV regu-larization image up sampling methods. Image acquisition process is modeled after a lowpass filtering followed by sampling. However, TV minimization is based on the assumption that the desirable image is almost piecewise constant, which yields a result with over smoothed homogeneous regions. Belahmidi and Guichard  have improved the TV-based interpolation by developing a nonlinear anisotropic PDE, hereafter referred to as BG interpolation method. In order to enhance edge preservation, this PDE performs a diffusion with strength and orientation adapted to image structures. This method balances linear zooming on homogeneous regions and anisotropic diffusion near edges, trying to combine the advantages of these two processes. The anisotropic diffusion scheme, including edge-directed or isophote-oriented method, uses the gradient to extract the image feature (edge) direction, that is, the gradient direction is considered to be the direction across the image feature. Nevertheless, the information contained in the gradient is local, not good to determine the edge directions. Nonlocal information should be considered in determining edge directions.

Recently, a nonlocal evolution equation and variations of it have been widely used to model diffusion processes in many areas [11,12]. Let us briefly introduce some references of nonlocal problem considered along this work. A nonlocal evolution equation corresponding to the Laplacian equation is presented as follows

It is called a nonlocal diffusion equation since the diffusion of the density at a point x and time t does not only depend on u(x; t), but on all the values of u in a neighborhood of x through the convolution term J * u. This equation shares many properties with the classical heat equation ut = Au. More precisely, as stated in , if u(t,x) is thought of as the density at the point x at time t, and J (x,y) is thought of as the probability distribution of jumping from location y to location x, then the convolution (J * u)(x,t) = jRN J(y - x)u(y,t)dy is the rate at which individuals are arriving to position x from all other places, and -u(x, t) = - J(y - x)u(y,t)dy is the rate at which they are leaving location x to travel to all other sites. This nonlocal evolution can be thought of as nonlocal isotropic diffusion.

In this paper, we propose a new method for image interpolation based on nonlocal p-Laplacian evolution. The nonlocal p-Laplacian and TV act as a regularizer to restrict edges of resulting image. The evolution is similar to an anisotropic energy dissipation process. The diffusion performs accurately along the direction of edges curves and its orthogonal direction. The magnitude of |u(y, t) - u(x, t) f-2 determines diffusion directions. It suppresses diffusion across the image feature direction and enhances diffusion along the image feature direction.

The rest of the paper is organized as follows. In Section 2, we review the method based on image up sampling with TV regularization in  and BG image interpolation . The proposed nonlocal p-Laplacian evolution for image interpolation is presented in Section 3. In Section 4, we demonstrate the experimental results to verify the effectiveness of our method, and the last section is for conclusion.

2. Background

We can think that a digital lower-resolution image u0 (input image) defined on some lattice is obtained by transforming a high-resolution image u (output image) defined on some better precision lattice, that is,

u0 = Hu, (2.1)

where H is a sparse matrix that combines both filtering and down sampling process. The goal of image interpolation is to solve the inverse problem (2.1), an ill-posed inverse problem. This ill-posed inverse problem is generally approached in a regularization-based framework, which would be formulated as an energy functional ,

E(u) = Jd (u,u0) + XJr (u), (2.2)

where X is a regularization parameter that controls the tradeoff between Jd and Jr. The data fidelity function Jd generally is formulated in the classical least-squares sense as Jd (u,u0) = (1/2)|Hu - u0|2. The TV regularizer Jr is taken as Jr(u) = JQ |Vu(x)|dx. The formula (2.2) is rewritten as follows

E(u) = X| |Vu|dx + 1 |Hu - u0|2. (2.3)

Using Euler's equation, the minimizer is the steady-state solution of the nonlinear PDE given by

ut = X|Vu| divf ^ + HT(Hu - u0). (2.4)

On the other hand, Belahmidi and Guichard solve the ill-posed inverse problem based on the classical heat diffusion model . Let n denote the direction of local gradient, and £ the direction perpendicular to the gradient, namely,

n = |vU| {Ux>u)> 5 = jvU|(~%'Ux)'

The second-order directional derivatives of the image u along the directions of n and & are easily computed as follows:

df2 " f

Uxx Uxy

d2U (Uxx Uxy

UxxUn 2UxyUxUy + UyyUx

*-yyUx

UXXUx + 2uxyuxuy + UWUi/

■lyyUy

The interpolation scheme based on heat diffusion model is formulated as follows:

d2u d2u

u = jVuj^ + g(jVuj)- HT(Hu - uo).

In this equation, the function g (s) is typically defined as

g(s) =

1 + (s/K)2

with K > 0 is a constant to be tuned for a particular application. The role of the diffusion coefficient g (jVuj) is to control the smoothing adaptively.

When g = 0, (2.7) reduces to (2.4) with X = 1. All these models can be viewed as interpolation schemes based on nonlinear diffusion model. The two regularizers jVuj div(Vu/jVuj) and jVuj(52u/d^2)+g(jVuj)(52u/5n2), respectively, in (2.4) and (2.7) result in different interpolation effects. In fact, jVuj div(Vu/jVuj) = jVuj(d2u/d&2) is the second-order directional derivative in the direction that is orthogonal to the gradient jVuj, and d1u/dn1 is the second-order directional derivative in the direction of the gradient jVuj.

From the viewpoint of geometry, the evolution processes in the artificial time t given by these models are seen as energy dissipation processes in two orthogonal directions n and & . The diffusion process of u(x,t) along I will preserve the location and the intensity transitions of the contours, while smoothing along them maintaining their crispness. This diffusion term is used to maintain edges with smooth isophotes in [9,10]. The diffusion of the grey values along n walks across both sides of the local image contour. This process blurs crisp contours as in the case of linear interpolators. Two divided means are used to deal with it. The diffusion process along n is cast aside in , while controlled by edge-stopping function g(jVuj) to balance the two diffusion terms in .

TV regularization in (2.3) does an excellent job at preserving edges while reconstructing images . This phenomenon can also be explained physically, since the resulting diffusion is strictly orthogonal to the gradient of the image. But TV-based interpolation favors solutions that are piecewise constant. This sometimes causes a staircasing effect in homogeneous regions, which are long observed in the literature in denoising, for example, [15,16]. Not only having blocky solutions, but they can also develop false edges in resulting image.

The function g in (2.7) is to be chosen with values between 0 and 1. The energy dissipation process (2.7) is adaptively controlled in the direction n and &, that is, minimal smoothing in the directions n (across the image features) preserving sharp edges, and maximal smoothing in the directions & (along the image features) obtaining smooth contours.

Uxy Uyy

The anisotropic diffusion scheme that uses the gradient to extract the image feature direction can mistakenly give maximal smoothing to the across feature direction and severely damage the image features, especially the image lines and textures . And blurry and/or oscillatory edges are introduced in interpolated image. The drawback of this model is that the gradient used to extract the image feature direction is too local. The information contained in the gradient is limited to a point and its immediate neighbors, while the edge curve that determines the edge directions is not a local event. The interpolation direction extraction should base on a larger neighborhood.

3. Nonlocal p-Laplacian Image Interpolation

In this section, we adopt nonlocal p-Laplacian evolution to overcome the local limit. Our proposed energy functional for regularized image interpolation is given by

E(u) = a( \Vu\dx + J(x - y)|u(y) - u(x)\dydx + 1 \Hu - u0\2. (3.1)

Jq 2p J JQ 2

The first part (the sum of the first term and the second) of right-hand side is regularizer, and the other is data fidelity function. The gradient flow associated to the functional E(u) is

ut(x,t) = a div(\Vu\-1Vu) + P(u) - HT(Hu - uo),

u(x, 0) = HTu0,

PJ(u) = J{x,y) \u(y,t) - u(x,t)\p 2(u{y,t) - u(x,t))dy. (3.3)

The kernel J : Q ^ R is assumed to be nonnegative, bounded continuous radial function, with sup(J) c B(0,d) and J"Q J(z)dz = 1.

The nonlocal energy dissipation is implemented mostly by PpJ(u) in our model. It is necessary to investigate the relation between the heat diffusion equation related to (2.7) and the p-Laplacian equation. The p-Laplacian evolution equation ut = is well

studied in image processing, which can be represented as 

ut = |V«r2+ (p - l)|V«r20• (3-4)

When p = 1, it is the TV flow keeping the edges but suffering from the staircase effect. When p = 2, ut = Au, this is isotropic diffusion because of the same diffusion coefficients. This model can smooth image, while bluring sharp edges. When 1 < p < 2, the grey values of u(x,t) in (3.4) diffuse along the directions n and I, respectively, as in the following equation:

ut = | Vu | ^ + g( | Vu | ) ^, (3-5)

which has different control factor compared with (3.4). That is to say, (3.4) and (3.5) all preserve advantage of adaptive smoothing.

A nonlocal improvement to p-Laplacian equation was studied in ,

ut(x,t) = J(x,y)\u(y,t) - u(x,t)\p 2(u{y,t) - u(x,t))dy. (3.6)

Moreover, nonlocal problems of type (3.6) have been used recently in the study of deblurring and denoising of images . With Neumann boundary conditions, the solutions to problem (3.6) converge to the solution of the classical p-Laplacian if p > 1 . The nonlocal p-La-placian evolution (3.6) improves the limit of the diffusion direction extraction depending on gradient (local information). The diffusion of the density at a point x and time t depends on all the values of u in a larger neighborhood of x. More precisely, if u(t,x) is thought of as the density at the point x at time t, and J(x,y) is thought of as the probability distribution of jumping from location y to location x, then PJ (u) is the rate at which individuals are arriving to position x from all other places. In image interpolation, PJ (u) is also the rate at which individuals are devoting to interpolated pixel x from all other pixels. The evolution process in the artificial time t given by (3.6) is seen as an anisotropic energy dissipation process. The direction of anisotropic diffusion is indicated by the magnitude of |u(y, t) - u(x, t)|p-2 in a larger neighborhood. It approximates to the direction of edge curve more accurate than the direction indicated by gradient. When 1 < p < 2, the energy dissipation process is adaptively controlled by |u(y, t) - u(x, t)|p-2 along the direction of edge curves and the orthogonal direction to edge curves. The diffusion process along the direction of edge curves is suppressed for small |u(y, t) - u(x, t)|p-2, and the diffusion along the orthogonal direction is enhanced for larger |u(y, t) - u(x, t)|p-2. This results in minimal smoothing in the directions across the image features preserving sharp edges and maximal smoothing in the directions along the image features reducing zigzagging artifacts and oscillatory.

4. Numerical Algorithm and Experimental Results

In this section, we develop a fully discrete numerical method to approximate problem (3.2). We recall first the notations in the finite differences scheme used in our paper. Let h and At be the space and time steps, respectively, and let (x^; x2j) = (ih; jh) be the grid points. Let un(i; j) be an approximation of the function u(nAt; x1i; x2j), with n > 0. Equation (3.2) can be discretized as follows:

un+1( i,j) - un(i,j) ( /Vu

At =(adiv( rV^)- HT(Hu - i

+ £ £ [J( (k,l) - j\un(k,l) - un{i,j)\p-2(un(k,l) - un(i,j)j\.

(fc,/)eQ

Figure 1: Barbara reduced and expanded by 2 x 2 (portion). (a) NEDI interpolation; (b) EGI interpolation; (c) BG method; (d) proposed method.

In all numerical experiments, we choose the following kernel function:

J (x) :=

Cexpi 2 0

|x|2 - d2 0 if x>d.

if x < d,

The constant C > 0 is selected such that JQ J (x)dx = 1.

We tested the proposed interpolation method on a variety of images. Some of the results are shown in Figures 1-3. Images are expanded by a factor of 2 x 2 in Figures 1 and 2 and by a factor of 10 x 10 in Figure 3. For comparison, we also show images interpolated using the BG image interpolation proposed in , the edge-guided image interpolation (EGI) in , and the edge-directed interpolation (NEDI) proposed in . The choice of the parameters is based on subjective quality of the results assessed informally by our personal preference as human viewers in terms of edge sharpness, contour crispness, no ringing in smooth regions, and no ringing near edges. We use the following parameters: a = 0.5, fi = 0.0001, and p = 1.8 for the proposed interpolation method, k = 0.0001 for the heat diffusion model, and time step At = 0.15 for these experiments. There is a non visible improvement on subjective or objective quality of the results when the parameters are not badly

changed. The iteration is terminated, when iterations.

In the first experiments, Barbara image with a size of 512 x512 was lowpass filtered and subsampled by a factor of 2 x 2, then the subsampled image was interpolated to the original image size. The interpolation was performed by four different methods, and a portion of the results is shown in Figure 1. Figure 2 shows a portion of a result interpolating a given flower image with a size of 320 x 240 by a factor of 2 x 2 without subsample. From the two examples, the NEDI interpolation method tends to introduce the zigzagging artifacts (Figures 1(a) and 2(a)). The BG results produce slight blurry edges (stripes in Figure 1(c)) and ringing in smooth regions (as shown in Figure 2(c)), because the direction decided by gradient misses their real directions as stated in Section 2. Our proposed method and the EGI method produce sharp edges and smooth contours, but the EGI method and the NEDI interpolation are applied only by a factor of 2 x 2.

Figure 2: Flower expanded by 2x2 (portion). (a) NEDI interpolation; (b) EGI interpolation; (c) BG method; (d) proposed method.

(d) (e) (f)

Figure 3: A portion of Barbara, Mandrill, and house images expanded by 10 x 10. Top row: BG method; bottom row: proposed method.

Table 1: Comparison of different interpolation algorithms using PSNR values.

Image The proposed BD EGI NEDI

Lena 34.6483 33.5438 30.5031 35.1024

Cameraman 26.7846 26.2239 24.3987 27.5413

Barche 28.0214 27.4137 25.7585 29.6238

Peppers 27.9167 27.4998 25.9046 29.8716

Mandrill 23.7041 23.2248 22.4871 24.8313

Resolution test 21.6378 21.0777 18.7030 22.7104

Barbara 25.2955 24.9608 24.3268 26.9344

Table 2: Comparison of different interpolation algorithms using MSSIM values.

Image The proposed BD EGI NEDI

Lena 0.9908 0.9836 0.9640 0.9609

Cameraman 0.8728 0.8575 0.8202 0.8105

Barche 0.8447 0.8237 0.7790 0.7671

Peppers 0.9144 0.9036 0.8704 0.8616

Mandrill 0.9528 0.9169 0.8851 0.8718

Resolution test 0.9147 0.9042 0.8547 0.8388

Barbara 0.9438 0.9204 0.9011 0.8952

The second experiment directly interpolates a portion of images (Barbara, mandrill, and house) by a factor of 10 x 10 shown in Figure 3. This experiment is performed by the BG method and our proposed method since the other methods only resize image by a factor of 2 x 2. It is clear from the figures that results obtained with the proposed approach are better than the results by the BG method in . The proposed method generates sharper and crisper stripes in Figure 3(d) compared to the result in Figure 3(a). The BG method produces blurry edges (the stripes in Figure 3(a), the beard in Figure 3(b)), zigzagging artifacts, and oscillatory (the line in Figure 3(c)), and images tend to be less natural. In the BG interpolation result, the boundaries of the text suffer from artifacts that make visualization difficult. It can be seen that our method results in an interpolated image with the fewer spurious patterns.

We use two measures, the classic PSNR and the mean structural similarity (MSSIM) index , to characterize the difference between the reference image and the output of a method. The MSSIM seems to approximate the perceived visual quality of an image better than PSNR or various other measures . MSSIM index takes values in [0,1] and increases as the quality increases. It is calculated by the code available at http://www .cns.nyu.edu/lcv/ssim/, using the default parameters. We use several test images of size 512x512 including Lena, mandrill, and Barbara of size 256x256 including cameraman, barche, peppers, and resolution test. To show the true power of the interpolation algorithms, we first downsampled each image by a factor of 2 x 2 and then interpolated the result back to its original size. The PSNR is shown on Table 1 and the MSSIM on Table 2. From the two tables, the proposed method yields improved PSNR (except for the NEDI), and MSSIM results in all the experiments. This improvement may be attributed to the fact that the nonlocal p-Lapla-cian evolution works better than other methods.

5. Conclusion

In this paper, a new image interpolation model based on TV and nonlocal p-Laplacian reg-ularization is proposed. It combines the advantages of TV regularizer and nonlocal p-Lapla-cian regularizer, that is, allowing discontinuities and preserving 1D image structures and the diffusion of the grey values of images along image feature direction. It overcomes the drawbacks of the anisotropic diffusion proposed by Belahmidi and Guichard. The direction of anisotropic diffusion is indicated by the information of image feature in a larger neighborhood. This results in minimal smoothing in the directions across the image features preserving sharp edges and maximal smoothing in the directions along the image features reducing zigzagging artifacts and oscillatory. We have shown improvement over nonlocal p-Laplacian on a subjective scale, and in many cases with an improvement in PSNR and MSSIM. We expect to prove convergence of the evolution equation in future work.

Acknowledgments

This research was partially supported by Natural Science Foundation Project of CQ CSTC (Grant no. cstcjjA40012) and the National Natural Science Foundation of China (Grant no. 10871217).

References

[1 ] S. Carrato, G. Ramponi, and S. Marsi, "Simple edge-sensitive image interpolation filter," in Proceedings of IEEE International Conference on Image Processing (ICIP '96), pp. 711-714, September 1996.

 H. Jiang and C. Moloney, "A new direction adaptive scheme for image interpolation," in Proceedings of IEEE International Conference on Image Processing (ICIP '02), pp. 369-372, September 2002.

 X. Li and M. T. Orchard, "New edge-directed interpolation," IEEE Transactions on Image Processing, vol. 10, no. 10, pp. 1521-1527, 2001.

 Q. Wang and R. Ward, "A new edge-directed image expansion scheme," in Proceedings of IEEE International Conference on Image Processing (ICIP '01), pp. 899-902, October 2001.

 B. S. Morse and D. Schwartzwald, "Image magnification using level-set reconstruction," in Proceedings ofIEEE Computer Society Conference on Computer Vision and Pattern Recognition, pp. 333-340, December 2001.

 Y. Cha and S. Kim, "Edge-forming methods for image zooming," Journal of Mathematical Imaging and Vision, vol. 25, no. 3, pp. 353-364, 2006.

 B. S. Morse and D. Schwartzwald, "Isophote-based interpolation," in Proceedings of IEEE International Conference on Image Processing (ICIP '98), vol. 3, pp. 227-231, Chicago, 111, USA, October 1998.

 F. Malgouyres and F. Guichard, "Edge direction preserving image zooming: a mathematical and numerical analysis," SIAM Journal on Numerical Analysis, vol. 39, no. 1, pp. 1-37, 2002.

 H. A. Aly and E. Dubois, "Image up-sampling using total-variation regularization with a new observation model," IEEE Transactions on Image Processing, vol. 14, no. 10, pp. 1647-1659, 2005.

 A. Belahmidi and F. Guichard, "A partial differential equation approach to image zoom," in Proceedings of IEEE International Conference on Image Processing (ICIP '04), pp. 649-652, October 2004.

 E. Chasseigne, M. Chaves, and J. D. Rossi, "Asymptotic behavior for nonlocal diffusion equations," Journal des Mathematiques Pures et Appliquees, vol. 86, no. 3, pp. 271-291, 2006.

 C. Cortazar, M. Elgueta, J. D. Rossi, and N. Wolanski, "Boundary fluxes for nonlocal diffusion," Journal of Differential Equations, vol. 234, no. 2, pp. 360-390, 2007.

 P. Fife, "Some nonclassical trends in parabolic and parabolic-like evolutions," in Trends in Nonlinear Analysis, pp. 153-191, Springer, Berlin, Germany, 2003.

 Y. Chen, S. Levine, and M. Rao, "Variable exponent, linear growth functionals in image restoration," SIAM Journal on Applied Mathematics, vol. 66, no. 4, pp. 1383-1406, 2006.

 W. Ring, "Structural properties of solutions to total variation regularization problems," Mathematical Modelling and Numerical Analysis, vol. 34, no. 4, pp. 799-810, 2000.

 Y. L. You, W. Xu, A. Tannenbaum, and M. Kaveh, "Behavioral analysis of anisotropic diffusion in image processing," IEEE Transactions on Image Processing, vol. 5, no. 11, pp. 1539-1553,1996.

 R. A. Carmona and S. Zhong, "Adaptive smoothing respecting feature directions," IEEE Transactions on Image Processing, vol. 7, no. 3, pp. 353-358,1998.

 H. Y. Zhang, Q. C. Peng, and Y. D. Wu, "Wavelet inpainting based on p-Laplace operator," Acta Automatica Sinica, vol. 33, no. 5, pp. 546-549, 2007.

 F. Andreu, J. M. Mazon, J. D. Rossi, and J. Toledo, "A nonlocal p-Laplacian evolution equation with Neumann boundary conditions," Journal des Mathematiques Pures et Appliquees, vol. 90, no. 2, pp. 201227, 2008.

 S. Kindermann, S. Osher, and P. W. Jones, "Deblurring and denoising of images by nonlocal function-als," Multiscale Modeling and Simulation, vol. 4, no. 4, pp. 1091-1115, 2005.

 L. Zhang and X. Wu, "An edge-guided image interpolation algorithm via directional filtering and data fusion," IEEE Transactions on Image Processing, vol. 15, no. 8, pp. 2226-2238, 2006.

 Z. Wang, A. C. Bovik, H. R. Sheikh, and E. P. Simoncelli, "Image quality assessment: from error visibility to structural similarity," IEEE Transactions on Image Processing, vol. 13, no. 4, pp. 600-612, 2004.

 A. Roussos and P. Maragos, "Vector-valued image interpolation by an anisotropic diffusion-projection PDE," in Proceedings of the 1st International Conference on Scale Space and Variational Methods in Computer Vision (SSVM '07), vol. 4485 of Lecture Notes in Computer Science, pp. 104-115, May 2007.

Copyright of Mathematical Problems in Engineering is the property of Hindawi Publishing Corporation and its content may not be copied or emailed to multiple sites or posted to a listserv without the copyright holder's express written permission. However, users may print, download, or email articles for individual use.