- Email: [email protected]
Sparse deconvolution method for improving the time-resolution of ultrasonic NDE signals
ARTICLE IN PRESS NDT&E International 42 (2009) 430–434
Contents lists available at ScienceDirect
NDT&E International journal homepage: www.elsevier.com/locate/ndteint
Sparse deconvolution method for improving the time-resolution of ultrasonic NDE signals Liang Wei a,, Zuo-ying Huang b, Que Pei-wenb a b
School of Electromechanical Engineering, University of Electronic Science and Technology of China, Chengdu 610054, PR China Departments of Information Measurement Technology and Instruments, Shanghai Jiaotong University, Shanghai 200240, PR China
a r t i c l e in f o
a b s t r a c t
Article history: Received 7 October 2008 Received in revised form 3 December 2008 Accepted 26 January 2009 Available online 6 February 2009
The estimation of the time-of-arrival (TOA) and/or time-of-flight (TOF) of the ultrasonic echoes is essential in ultrasonic non-destructive evaluation. In this paper, a sparse deconvolution method is proposed for deconvolving ultrasonic signals. To obtain the transducer pulse-echo wavelet, matching pursuit (MP) method has been used to analyze noisy ultrasonic pulse-echo wavelet and decompose the noisy pulse-echo wavelet into unit-norm vectors, then, the approximation pulse-echo wavelet obtained by the reconstruction result with several large coefficient decomposed unit-norm vectors. A weighting matrix is the central to the efficiency of the sparse deconvolution method. It increases the sparsity of the deconvolution result and decreases the influence of the added noise. The deconvolution gives highresolution ultrasonic reflectivity function. We can obtain the accurate TOA or TOF from the ultrasonic reflectivity function. Results using real ultrasonic data show the sparse deconvolution method yields good estimates of the spikes associated with the reverberation echoes and the thickness of a steel sample can be calculated by the ultrasonic reflectivity function with a reasonable accuracy. & 2009 Elsevier Ltd. All rights reserved.
Keywords: MP Over-complete dictionary TOA Sparse deconvolution
1. Introduction Many ultrasonic testing applications are based on the estimation of the time-of-arrival (TOA) and/or time-of-flight (TOF) of the ultrasonic echoes. It produces a peak at the TOA of the received signal. Unfortunately, the received signals are also contaminated by noise originating from both the measurement system and specimen. The propagation of the ultrasonic wave is influenced by reflections at all boundary surfaces in the media. These reflections contained in the ultrasonic wave leaving the material are suitable to obtain the necessary information about the properties of the test material, i.e., thickness, defects, and layers of different materials. Sometimes, however, a small defect is masked by the echo from a larger nearby reflector, and it will be difficult to be detected. To discover the internal structure of the test sample as precisely as possible, various groups have explored different techniques to restore the signal leaving the material surface containing only the reflection information from the detected ultrasonic signal [1]. In ultrasonic inspection, the received ultrasonic signal can be modeled as a convolution between the transducer’s impulse response and a reflection sequence [2]. Deconvolution is one method that can be utilized to restore the reflection sequence or rather to find a good estimate of it.
Corresponding author.
E-mail address: [email protected] (L. Wei). 0963-8695/$ - see front matter & 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.ndteint.2009.01.009
Deconvolution, which is used mainly for improving the resolution of the signals, can be particularly useful for revealing a small defect that has been masked by the echo from a larger nearby reflector [3]. In conventional deconvolution, the system and the system’s output must be known, and the problem is to estimate the input. In blind deconvolution both the system and the system’s input are not known, yet the system’s input is desired [4]. Wiggins [5] proposed a method called minimum entropy deconvolution (MED). The approach was to use the knowledge that the source signal had a sparse distribution, and try to find the deconvolution filter whose output distribution was as sparse as possible. Based on this idea, a number of algorithms using higher-order statistics are developed [1,6,7]. These algorithms seek the smallest number of large spikes consistent with the data. Thus, the order in the data is maximized, or equivalently, the entropy or disorder in the data is minimized. However, the reflection sequence restored by the above-mentioned deconvolution methods is not sparse enough. In other words, the resolution of the reflection sequence is not high enough for computing the exact time between ultrasound impulse signals. Furthermore, the high computational cost and, in some cases, compromised convergence, is another drawback of these methods. In this study, we use matching pursuit (MP) method to analyze noisy ultrasonic pulse-echo and decompose the noisy pulse-echo into unit-norm vectors. The error ratio of the mean square is used to monitor the acceptability of the unit-norm vector obtained by
ARTICLE IN PRESS L. Wei et al. / NDT&E International 42 (2009) 430–434
the matching pursuit method. The approximation of the ultrasonic pulse-echo can be obtained by the reconstruction result with several large coefficient unit-norm vectors. Using the reconstruction transducer pulse-echo, a sparse deconvolution method is proposed for deconvolving ultrasonic signals. A weighting matrix is central to the efficiency of the sparse deconvolution method. It increases the sparsity of the deconvolution result and decreases the influence of the added noise. The deconvolution gives high-resolution ultrasonic reflectivity function. Results using real ultrasonic data show the sparse deconvolution method yields good estimates of the spikes associated with the reverberation echoes and the thickness of a steel sample can be calculated by the ultrasonic reflectivity function with a reasonable accuracy. The remainder of this paper is organized as follows. Section 2 is devoted to the study of convolutional model of ultrasound signal. In Section 3, we present the modified blind deconvolution algorithm. In Section 4, the definition and basic properties of the matching pursuit are recalled. Meanwhile, the error ratio of the mean square is developed in this section. In Section 5, results obtained from the analysis of actual ultrasonic signals are presented. Conclusions are given in the last.
2. Convolutional model of ultrasound signal In ultrasonic non-destructive testing, the measured ultrasound signal y(n) can be modeled as 1 X
yðnÞ ¼
hðn kÞxðkÞ
(1)
k¼1
where h(n) represents the transmitted wavelet which is influenced by the propagation as well as the transducer impulse response (distortion function). x(k) is the reflectivity function, it is known a priori that x(k) is sparse which consists of a very low number of noncontiguous nonzero entries. Due to the sparsity assumption, x(k) is more efficiently represented by two vectors: t giving the time of the spikes which is related to the location of the reflector as the distance from the transducer over the speed of ultrasound in the propagation path, and b giving their amplitudes which is primarily governed by the impedance, size, or orientation of the reflector. The relationship between x(k),t, and b is xðkÞ ¼
M X
bi dðn ti Þ;
n ¼ 1; 2; . . . N
(2)
i¼1
where M is the number of spikes, and d( ) is the Dirac delta function, which is one when its argument is zero and zero otherwise.
3. Sparse deconvolution method In ultrasonic inspection, in order to improve the timeresolution of the measured ultrasonic signal, it is effective to restore the reflectivity function x(n) or rather to find a good estimate of it, i.e., find a good estimate of the reflectivity vector c ¼ [b,t]. The Eq. (1) can be written in matrix notation y ¼ Hx
(3)
where H depends on t and is given by Hnk ¼ h(nk) for k ¼ 1,2, y, N, and N is the length of the ultrasound signal y(n). Thus, each column of H contains a copy of the transducer pulse-echo wavelet h(n) that is shifted to the corresponding spike position. The minimum norm solution is the most widely used estimate for Eq. (3) and is found by assuming the minimum Euclidian or
431
P l2-norm criterion minJxJ ¼ ( nN¼ 1x2(n))1/2 on the solution. This solution is unique and is computed as x ¼ Hþ y
(4)
where H+ ¼ H0 (HH0 )1 denotes the Moore–Penrose inverse [8,9]. The solution has a number of computational advantages, but it does not provide sparse solutions. In the focal underdetermined system solver (FOCUSS), the sparse solution to Eq. (3) is achieved as follows: x ¼ WðHWÞþ y
(5)
where W is the weighting matrix defined as W ¼ DiagðWð1Þ; Wð2Þ; . . . ; WðNÞÞ
(6) +
P
Then, the cost objective becomes JW xJ ¼ ( W(n))2)1/2. So the iteration resolution solution in the ith step is xi ¼ W i ðHW i Þþ y
N n ¼ 1(x(n)/
(7)
The iteration operation process can be terminated once the convergence pattern becomes clear. In order to obtain the sparse solutions and decrease the influence of the added noise, in this paper, the ith weighting matrix Wi can be selected as follows [10]: W i ðnÞ ¼ "
pffiffiffiffi N
1 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi, N P x2i ðkÞ k¼1;kan
N P
!# jxi ðkÞj
k¼1;kan
1 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi, " !# N N pffiffiffiffi P P 2 xi ðkÞ jxi ðkÞj N k¼1
(8)
k¼1
According to Eq. (8), if the value in xi was large, the pffiffiffiffiqffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN 2 corresponding element of the term ½ N k¼1;kan xi ðkÞ PN = k¼1;kan jxi ðkÞj is small. In other words, the relatively large pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN PN 2 values of the term ½ N ð k¼1 xi ðkÞ= k¼1 jxi ðkÞjÞ reduce the contribution of the corresponding elements of xi to the cost Eq. (8), and vice versa. Thus, larger entries in xi(k) result in larger corresponding entries in xi(n)/Wi(n), the sparsity of xi(n)/Wi(n) is larger than that of xi(n). This indicates that this method can obtain the sparse solutions, at the same time, weakens the influence of the added noise. A simulated reflection sequence x(n), which is shown in Fig. 1(a), consisted of three spikes, inserted in zero mean white Gaussian noise sequence with scalar variance of 0.002. Fig. 1(b) and (c) shows the results processed by the non-linear functions pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN PN 2 PN PN 2 ½ Nð k¼1;kan xi ðkÞ= k¼1;kan jxi ðkÞjÞ and ½ N ð k¼1 xi ðkÞ= k¼1 jxi ðkÞjÞ, respectively. It is shown in Fig. 1(b) and (c) that two processed results are almost the same if x(n) was the added noise, but result processed by the non-linear functions pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN PN 2 ½ Nð k¼1 xi ðkÞ= k¼1 jxi ðkÞjÞ is larger than that processed by pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN PN 2 the non-linear functions ½ Nð k¼1;kan xi ðkÞ= k¼1;kan jxi ðkÞjÞ if x(n) was the large spikes. Fig. 1(d) shows the weighting matrix W(n). This indicates that this method can obtain the sparse solutions, at the same time, weakens the influence of the added noise and it almost holds only the already prominent entries. However, in order to obtain the sparse solution to Eq. (3), an accurate estimate of the reflectivity vector is on the premise that the high-resolution pulse-echo wavelet (transducer impulse response or the reference echo) estimation is obtained firstly.
ARTICLE IN PRESS L. Wei et al. / NDT&E International 42 (2009) 430–434
Amplitude
432
4 2 0
Amplitude
30 20 10 0 -10
Amplitude
30 20 10 0 -10
Amplitude
-2 0
100
200
300
400
500
600
700
800
900
1000
0
100
200
300
400
500
600
700
800
900
1000
0
100
200
300
400
500
600
700
800
900
1000
0
100
200
300
400
500 600 Time index
700
800
900
1000
4 2 0
Fig. 1. Enhance the sparsity of the reflection sequence by the weighting matrix W(n): (a) a simulated reflection sequence; (b) result processed by the non-linear functions pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN 2 ffi PN PN PN 2 N ; (c) result processed by the non-linear functions N ; (d) the weighting matrix W(n). k¼1 xi ðkÞ k¼1 jxi ðkÞj k¼1;kan xi ðkÞ k¼1;kan jxi ðkÞj
4. Matching pursuit for utrasonic pulse-echo estimation
compute all the inner products with the residual rm
Let H be the Hilbert space H ¼ RN. The inner product of vectors P f,gAH is defined as /f,gS nf(n)g(n) and the norm of a vector f is 1/2 JfJ ¼ /f,fS . Let D ¼ fg i ; i 2 Zg, where gi is unit-norm vector (JgiJ ¼ 1), span(D) ¼ H and card(D)bN be an over-complete basis for H. D can be called over-complete dictionary. If a measured noisy ultrasonic pulse-echo can be expressed as follows:
g iðmÞ ¼
xðtÞ ¼ hðtÞ þ eðtÞ
(9)
where h(t) is the ultrasonic pulse-echo, e(t) can be characterized as white Gaussian noise (WGN). We look for decompositions of the signal x(t) on the overcomplete dictionary D of indexed unit-norm vector gi, in the form x¼
m X
ai g i þ rm
(10)
i¼1
where ai is the weight associated to the unit-norm vector gi. Our goal is to look for the appropriate number m to make the Pm term i ¼ 1ai gi be the approximation of the ultrasonic pulseecho h(t) with a higher quality [13]. Matching pursuit [11] is an iterative algorithm that offers a suboptimal solution for decomposing a signal f in terms of unitnorm vector gi chosen from an over-complete dictionary D. At the first iteration, the unit-norm vector gi which has the largest inner product with analyzed signal f is chosen. Then the contribution of this vector is then subtracted from the signal, and the process is repeated on the residual. At the mth iteration, the residue is ( f; m¼0 m (11) r ¼ rmþ1 þ aiðmÞ g iðmÞ ; ma0 where ai(m) is the weight associated with the optimum unit-norm vector gi(m) at the mth iteration. The weight ai(m) associated with each unit-norm vector giAD at the mth iteration is introduced to
hr m ; g i i hr m ; g i i ¼ ¼ hrm ; g i i hg i ; g i i kg i k2
(12)
The optimum unit-norm vector gi(m) (and its weight ai(m)) at the mth iteration can be obtained g iðmÞ ¼ arg min kr mþ1 k2 ¼ arg max jaiðmÞ j2 g i 2D
g i 2D
(13)
The correlation updating procedure is performed as follows: hr m ; g i i ¼ hrm ; g i i aiðmÞ hg iðmÞ ; g i i
(14)
The correlations /gi(m),giS can be precalculated and stored, once over-complete dictionary D has been determined [12]. The criterion chosen to stop the iteration is: for two given residuals rm and rm+1, if the error ratio of the mean square e ¼ (rm+1mrm)/mrm is less than a specified tolerance, then the matching pursuit stops. Where m ¼ (E((rm)2)/E((rm+1)2))1/2 and E( ) is the expected value. m makes the two given residuals rm and rm+1 have the same expected value, and it also can guarantee the convergence of the ratio e ¼ (rm+1mrm)/mrm in the mean square. The error criterion is used to monitor the acceptability of the unit-norm vector obtained by the matching pursuit performance. And the matching pursuit operation will continue until the error ratio of the mean square (rm+1mrm)/mrm is less than a predetermined threshold. The error criterion is an important factor for the pulse-echo approximation while noise level is unknown; furthermore, it is generally different for different ultrasonic signals. Compared to that of the residual rm, when the mean square of the term (rm+1mrm) can be neglected, the optimum unit-norm vector gi(m) at the mth iteration will make little contribution to the pulse-echo, it can be regarded as a partial signal come from the additive noise. Thus, the matching pursuit operation process can be terminated.
ARTICLE IN PRESS L. Wei et al. / NDT&E International 42 (2009) 430–434
According to the MP algorithm, our method has to define the unit-norm vectors in the dictionary D and the way to approximate pulse-echo. Since, a real time implementation is claimed for usefulness, the complexity of this method must be reduced as much as possible. Because MP is an iterative algorithm, the goal can be achieved if both the amount of operations at each iteration and the number of iterations to estimate pulse-echo are reduced. In ultrasonic non-destructive testing, the magnitude spectrum of the transducer pulse-echo has band-pass characteristics. Hence, the over-complete dictionary D can compose of unit-norm vectors gi,k that mimic the ultrasonic pulse-echo to be estimated, and in turn depend on the impulsive response of the transducer used in the experiments. The unit-norm vector in the dictionary D can be defined as the time-shifted Gaussian echo wavelet, and can be expressed as 1 g i;k ½n ¼ pffiffiffi exp½gk ðn ti Þ2 cosð2pf c ðn ti Þ jk Þ
ti ¼ 1; . . . ; L; L ¼ M þ N; n ¼ ti
M M þ 1; . . . ; ti þ 2 2
(15)
where gk, fc and jk are the Gaussian echo wavelet parameters. The pffiffiffi normalizes the energy of Gaussian echo wavelet, and termp1= ffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ p=2gk . They are related to the frequency and the bandwidth of the transducer response. L is the number of time shifts to take into account all possible flaw echo positions, N is the ultrasonic data length and M is the length of the unit-norm vectors gi,k.
5. Experimental results Many ultrasonic non-destructive testing applications require the accurate estimation of TOA or/and TOF of the ultrasonic echoes. In this paper, ultrasonic detection signals are analyzed in order to validate the proposed method. In the first experimental scheme, a steel sample block with approximately 1.0 cm thickness is tested using the unfocussed longitudinal transducer with 2.5 MHz central frequency. The ultrasound signals were recorded in the format of A-scan with a sample frequency of 100 MHz. Fig. 2(a) shows the measured
433
echoes reflected from the front and back surfaces of the steel sample block. Three echoes can be found in the signal. The first and second echoes are visible at the delay of approximately 1.13 and 4.56 ms, respectively while the third echo can be assumed to arrive at 8.04 ms. The TOFs of two neighboring ultrasonic pulseecho are 3.43, 3.48 ms, respectively. The average error is about 0.05 ms. It is shown from Fig. 2(a) that the ultrasound signal has a low time-resolution. Fig. 2(b) plots the time domain representation of the measured ultrasonic pulse-echo in solid line reflected from the inner surface of the offshore pipeline sample, and it also plots the reconstruction result in dotted line with the former fourteen large coefficient unit-norm vectors. It is shown in Fig. 2(b) that the reconstruction result achieves a reasonable quality. Fig. 2(c) shows the result of the sparse deconvolution method proposed in this paper. The delay between the first and the second echo is 3.44 ms, and the delay between the second and third echo is also 3.45 ms. The average error is only about 0.01 ms. If the sound velocity for longitudinal waves in ordinary steel is 5.9 mm/ms, the thickness of the steel sample block can be computed as 10.148 and 10.178 mm. The average error is only about 0.03 mm. This result makes evidence that the good performance of the sparse deconvolution method which is in much improved accuracy in the estimation of echo arrival time. In the second experimental scheme, an offshore pipeline specimen was applied in which a man-made 10 1 5.6 mm3 flaw was fabricated close to the outer wall. The experiment was done inside the pipe. The ultrasonic probe with 10 MHz central frequency and 6 mm diameter was used and water was selected as couplant. The flaw echoes were recorded in the format of A-scan with a sample frequency of 100 MHz. Fig. 3(a) shows the measured ultrasonic signal. Fig. 3(b) plots the time domain representation of the measured ultrasonic pulse-echo reflected from the inner surface of the offshore pipeline sample (in solid line) and the reconstruction result in dotted line with the former fourteen large coefficient unit-norm vectors. Fig. 3(c) shows the deconvolution result of the sparse deconvolution method. It is shown from Fig. 3(c) that the major peaks are visible in the deconvolution result. Three TOFs, 0.81, 2.95, and 5.01 ms can be identified clearly in Fig. 3(c). The depth of the flaw can be calculated: 6.077 mm.
Fig. 2. The sparse deconvolution method for the estimation of echo arrival time: (a) the measured ultrasonic pulse-echo reflected from the front surface of the steel sample block; (b) comparison between the approximation result (dotted line) and the measured ultrasonic pulse-echo (solid line); (c) deconvolution result using the sparse deconvolution method.
ARTICLE IN PRESS 434
L. Wei et al. / NDT&E International 42 (2009) 430–434
Fig. 3. The sparse deconvolution method for the estimation of echo arrival time: (a) the measured ultrasonic signal received from a 10 1 5.6 mm3 flaw close to the outer wall; (b) comparison between the approximation result (dotted line) and the measured ultrasonic pulse-echo (solid line); (c) deconvolution result using the sparse deconvolution method.
6. Conclusion Many ultrasonic non-destructive testing applications require the accurate estimation of TOA or/and TOF of the ultrasonic echoes. However, the received signals might have been contaminated by noise originating from both the measurement system and the specimen. In these cases, it is impossible to obtain the accurate TOA or TOF from the received signals. Based the characteristics of the transducer pulse-echo, the over-complete dictionary D has been defined. We use the matching pursuit method to analyze noisy ultrasonic pulse-echo and decompose the noisy ultrasonic pulse-echo into elementary components. Because, the over-complete dictionary D is composed of the unitnorm vectors that mimic the transducer impulse response, the complexity of the matching pursuit method is reduced. The method is very efficient in eliminating the noise and the reconstruction result achieves a reasonable quality. Using the reconstruction transducer pulse-echo wavelet, a sparse deconvolution method is proposed for deconvolving ultrasonic signals. A weighting matrix is central to the efficiency of the sparse deconvolution method. It increases the sparsity of the deconvolution result and decreases the influence of the added noise. The deconvolution gives high-resolution ultrasonic reflectivity function. We can obtain the accurate TOA or TOF from the ultrasonic reflectivity function. Results using real ultrasonic data show the sparse deconvolution method that yields good estimates of the spikes associated with the reverberation echoes and the thickness of a steel sample can be calculated by the ultrasonic reflectivity function with a reasonable accuracy.
Acknowledgement This work was supported by the China National High-Tech Research and Development Program, the National Science
Foundation for Young Scientists of China and the JIN CHUAN Advance Research Program of Science and Technology. References [1] Nandi AK, Mampel D, Roscher B. Blind deconvolution of ultrasonic signals in non-destructive testing applications. IEEE Transactions on Signal Processing 1997;45:1382–9. [2] Cassereau D, Guyomar D, Fink M. Time deconvolution of diffraction effects—application to calibration and prediction of transducer waveforms. Journal of the Acoustical Society of America 1988;84:1073–85. [3] Olofsson T. Deconvolution and model-based restoration of clipped ultrasonic signals. IEEE transactions on instrumentation and measurement 2005; 54:1235–40. [4] Lee JY, Nandi AK. Extraction of impacting signals using blind deconvolution. Journal of Sound and Vibration 2000;232:945–62. [5] Wiggins RA. Minimum entropy deconvolution. Geophysical Exploration 1978;16:21–35. [6] Nandi Asoke K. Blind identification of FIR systems using third order cumulates. Signal Processing 1994;39:131–47. [7] Nandi AK, Mehlan R. Parameter estimation and phase reconstruction of moving average processes using third-order cumulants. Mechanical Systems and Signal Processing 1994;8:421–36. [8] Gorodnitsky Irina F, Rao Bhaskar D. Sparse signal reconstruction from limited data using FOCUSS: a re-weighted minimum norm algorithm. IEEE Transactions on Signal Processing 1997;45:600–16. [9] Liu HH, Gao XR, Schimpf Paul H. A recursive algorithm for the threedimensional imaging of brain electric activity. Shrinking LORETA–FOCUSS 2004;51:1794–802. [10] Liang W, Lei HM, Que PW. Sparsity enhancement for blind deconvolution of ultrasonic signals in nondestructive testing application. Review of Scientific Instruments 2008;79:014901. [11] Mallat S, Zhang Z. Matching pursuits with time-frequency dictionaries. IEEE Transactions on Signal Processing 1993;41:3397. [12] Ruiz-Reyes N, Vera-Candeas P, Curpia´n-Alonso J. New matching pursuit-based algorithm for SNR improvement in ultrasonic NDT. NDT & E International 2005;38:453. [13] Liang W, Que PW, Lei HM, et al. Matching pursuit for decomposition and approximation of ultrasonic pulse-echo wavelet and its application in ultrasonic nondestructive evaluation. Review of Scientific Instruments 2008;79(7):075105.
Contents lists available at ScienceDirect
NDT&E International journal homepage: www.elsevier.com/locate/ndteint
Sparse deconvolution method for improving the time-resolution of ultrasonic NDE signals Liang Wei a,, Zuo-ying Huang b, Que Pei-wenb a b
School of Electromechanical Engineering, University of Electronic Science and Technology of China, Chengdu 610054, PR China Departments of Information Measurement Technology and Instruments, Shanghai Jiaotong University, Shanghai 200240, PR China
a r t i c l e in f o
a b s t r a c t
Article history: Received 7 October 2008 Received in revised form 3 December 2008 Accepted 26 January 2009 Available online 6 February 2009
The estimation of the time-of-arrival (TOA) and/or time-of-flight (TOF) of the ultrasonic echoes is essential in ultrasonic non-destructive evaluation. In this paper, a sparse deconvolution method is proposed for deconvolving ultrasonic signals. To obtain the transducer pulse-echo wavelet, matching pursuit (MP) method has been used to analyze noisy ultrasonic pulse-echo wavelet and decompose the noisy pulse-echo wavelet into unit-norm vectors, then, the approximation pulse-echo wavelet obtained by the reconstruction result with several large coefficient decomposed unit-norm vectors. A weighting matrix is the central to the efficiency of the sparse deconvolution method. It increases the sparsity of the deconvolution result and decreases the influence of the added noise. The deconvolution gives highresolution ultrasonic reflectivity function. We can obtain the accurate TOA or TOF from the ultrasonic reflectivity function. Results using real ultrasonic data show the sparse deconvolution method yields good estimates of the spikes associated with the reverberation echoes and the thickness of a steel sample can be calculated by the ultrasonic reflectivity function with a reasonable accuracy. & 2009 Elsevier Ltd. All rights reserved.
Keywords: MP Over-complete dictionary TOA Sparse deconvolution
1. Introduction Many ultrasonic testing applications are based on the estimation of the time-of-arrival (TOA) and/or time-of-flight (TOF) of the ultrasonic echoes. It produces a peak at the TOA of the received signal. Unfortunately, the received signals are also contaminated by noise originating from both the measurement system and specimen. The propagation of the ultrasonic wave is influenced by reflections at all boundary surfaces in the media. These reflections contained in the ultrasonic wave leaving the material are suitable to obtain the necessary information about the properties of the test material, i.e., thickness, defects, and layers of different materials. Sometimes, however, a small defect is masked by the echo from a larger nearby reflector, and it will be difficult to be detected. To discover the internal structure of the test sample as precisely as possible, various groups have explored different techniques to restore the signal leaving the material surface containing only the reflection information from the detected ultrasonic signal [1]. In ultrasonic inspection, the received ultrasonic signal can be modeled as a convolution between the transducer’s impulse response and a reflection sequence [2]. Deconvolution is one method that can be utilized to restore the reflection sequence or rather to find a good estimate of it.
Corresponding author.
E-mail address: [email protected] (L. Wei). 0963-8695/$ - see front matter & 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.ndteint.2009.01.009
Deconvolution, which is used mainly for improving the resolution of the signals, can be particularly useful for revealing a small defect that has been masked by the echo from a larger nearby reflector [3]. In conventional deconvolution, the system and the system’s output must be known, and the problem is to estimate the input. In blind deconvolution both the system and the system’s input are not known, yet the system’s input is desired [4]. Wiggins [5] proposed a method called minimum entropy deconvolution (MED). The approach was to use the knowledge that the source signal had a sparse distribution, and try to find the deconvolution filter whose output distribution was as sparse as possible. Based on this idea, a number of algorithms using higher-order statistics are developed [1,6,7]. These algorithms seek the smallest number of large spikes consistent with the data. Thus, the order in the data is maximized, or equivalently, the entropy or disorder in the data is minimized. However, the reflection sequence restored by the above-mentioned deconvolution methods is not sparse enough. In other words, the resolution of the reflection sequence is not high enough for computing the exact time between ultrasound impulse signals. Furthermore, the high computational cost and, in some cases, compromised convergence, is another drawback of these methods. In this study, we use matching pursuit (MP) method to analyze noisy ultrasonic pulse-echo and decompose the noisy pulse-echo into unit-norm vectors. The error ratio of the mean square is used to monitor the acceptability of the unit-norm vector obtained by
ARTICLE IN PRESS L. Wei et al. / NDT&E International 42 (2009) 430–434
the matching pursuit method. The approximation of the ultrasonic pulse-echo can be obtained by the reconstruction result with several large coefficient unit-norm vectors. Using the reconstruction transducer pulse-echo, a sparse deconvolution method is proposed for deconvolving ultrasonic signals. A weighting matrix is central to the efficiency of the sparse deconvolution method. It increases the sparsity of the deconvolution result and decreases the influence of the added noise. The deconvolution gives high-resolution ultrasonic reflectivity function. Results using real ultrasonic data show the sparse deconvolution method yields good estimates of the spikes associated with the reverberation echoes and the thickness of a steel sample can be calculated by the ultrasonic reflectivity function with a reasonable accuracy. The remainder of this paper is organized as follows. Section 2 is devoted to the study of convolutional model of ultrasound signal. In Section 3, we present the modified blind deconvolution algorithm. In Section 4, the definition and basic properties of the matching pursuit are recalled. Meanwhile, the error ratio of the mean square is developed in this section. In Section 5, results obtained from the analysis of actual ultrasonic signals are presented. Conclusions are given in the last.
2. Convolutional model of ultrasound signal In ultrasonic non-destructive testing, the measured ultrasound signal y(n) can be modeled as 1 X
yðnÞ ¼
hðn kÞxðkÞ
(1)
k¼1
where h(n) represents the transmitted wavelet which is influenced by the propagation as well as the transducer impulse response (distortion function). x(k) is the reflectivity function, it is known a priori that x(k) is sparse which consists of a very low number of noncontiguous nonzero entries. Due to the sparsity assumption, x(k) is more efficiently represented by two vectors: t giving the time of the spikes which is related to the location of the reflector as the distance from the transducer over the speed of ultrasound in the propagation path, and b giving their amplitudes which is primarily governed by the impedance, size, or orientation of the reflector. The relationship between x(k),t, and b is xðkÞ ¼
M X
bi dðn ti Þ;
n ¼ 1; 2; . . . N
(2)
i¼1
where M is the number of spikes, and d( ) is the Dirac delta function, which is one when its argument is zero and zero otherwise.
3. Sparse deconvolution method In ultrasonic inspection, in order to improve the timeresolution of the measured ultrasonic signal, it is effective to restore the reflectivity function x(n) or rather to find a good estimate of it, i.e., find a good estimate of the reflectivity vector c ¼ [b,t]. The Eq. (1) can be written in matrix notation y ¼ Hx
(3)
where H depends on t and is given by Hnk ¼ h(nk) for k ¼ 1,2, y, N, and N is the length of the ultrasound signal y(n). Thus, each column of H contains a copy of the transducer pulse-echo wavelet h(n) that is shifted to the corresponding spike position. The minimum norm solution is the most widely used estimate for Eq. (3) and is found by assuming the minimum Euclidian or
431
P l2-norm criterion minJxJ ¼ ( nN¼ 1x2(n))1/2 on the solution. This solution is unique and is computed as x ¼ Hþ y
(4)
where H+ ¼ H0 (HH0 )1 denotes the Moore–Penrose inverse [8,9]. The solution has a number of computational advantages, but it does not provide sparse solutions. In the focal underdetermined system solver (FOCUSS), the sparse solution to Eq. (3) is achieved as follows: x ¼ WðHWÞþ y
(5)
where W is the weighting matrix defined as W ¼ DiagðWð1Þ; Wð2Þ; . . . ; WðNÞÞ
(6) +
P
Then, the cost objective becomes JW xJ ¼ ( W(n))2)1/2. So the iteration resolution solution in the ith step is xi ¼ W i ðHW i Þþ y
N n ¼ 1(x(n)/
(7)
The iteration operation process can be terminated once the convergence pattern becomes clear. In order to obtain the sparse solutions and decrease the influence of the added noise, in this paper, the ith weighting matrix Wi can be selected as follows [10]: W i ðnÞ ¼ "
pffiffiffiffi N
1 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi, N P x2i ðkÞ k¼1;kan
N P
!# jxi ðkÞj
k¼1;kan
1 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi, " !# N N pffiffiffiffi P P 2 xi ðkÞ jxi ðkÞj N k¼1
(8)
k¼1
According to Eq. (8), if the value in xi was large, the pffiffiffiffiqffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN 2 corresponding element of the term ½ N k¼1;kan xi ðkÞ PN = k¼1;kan jxi ðkÞj is small. In other words, the relatively large pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN PN 2 values of the term ½ N ð k¼1 xi ðkÞ= k¼1 jxi ðkÞjÞ reduce the contribution of the corresponding elements of xi to the cost Eq. (8), and vice versa. Thus, larger entries in xi(k) result in larger corresponding entries in xi(n)/Wi(n), the sparsity of xi(n)/Wi(n) is larger than that of xi(n). This indicates that this method can obtain the sparse solutions, at the same time, weakens the influence of the added noise. A simulated reflection sequence x(n), which is shown in Fig. 1(a), consisted of three spikes, inserted in zero mean white Gaussian noise sequence with scalar variance of 0.002. Fig. 1(b) and (c) shows the results processed by the non-linear functions pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN PN 2 PN PN 2 ½ Nð k¼1;kan xi ðkÞ= k¼1;kan jxi ðkÞjÞ and ½ N ð k¼1 xi ðkÞ= k¼1 jxi ðkÞjÞ, respectively. It is shown in Fig. 1(b) and (c) that two processed results are almost the same if x(n) was the added noise, but result processed by the non-linear functions pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN PN 2 ½ Nð k¼1 xi ðkÞ= k¼1 jxi ðkÞjÞ is larger than that processed by pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN PN 2 the non-linear functions ½ Nð k¼1;kan xi ðkÞ= k¼1;kan jxi ðkÞjÞ if x(n) was the large spikes. Fig. 1(d) shows the weighting matrix W(n). This indicates that this method can obtain the sparse solutions, at the same time, weakens the influence of the added noise and it almost holds only the already prominent entries. However, in order to obtain the sparse solution to Eq. (3), an accurate estimate of the reflectivity vector is on the premise that the high-resolution pulse-echo wavelet (transducer impulse response or the reference echo) estimation is obtained firstly.
ARTICLE IN PRESS L. Wei et al. / NDT&E International 42 (2009) 430–434
Amplitude
432
4 2 0
Amplitude
30 20 10 0 -10
Amplitude
30 20 10 0 -10
Amplitude
-2 0
100
200
300
400
500
600
700
800
900
1000
0
100
200
300
400
500
600
700
800
900
1000
0
100
200
300
400
500
600
700
800
900
1000
0
100
200
300
400
500 600 Time index
700
800
900
1000
4 2 0
Fig. 1. Enhance the sparsity of the reflection sequence by the weighting matrix W(n): (a) a simulated reflection sequence; (b) result processed by the non-linear functions pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN 2 ffi PN PN PN 2 N ; (c) result processed by the non-linear functions N ; (d) the weighting matrix W(n). k¼1 xi ðkÞ k¼1 jxi ðkÞj k¼1;kan xi ðkÞ k¼1;kan jxi ðkÞj
4. Matching pursuit for utrasonic pulse-echo estimation
compute all the inner products with the residual rm
Let H be the Hilbert space H ¼ RN. The inner product of vectors P f,gAH is defined as /f,gS nf(n)g(n) and the norm of a vector f is 1/2 JfJ ¼ /f,fS . Let D ¼ fg i ; i 2 Zg, where gi is unit-norm vector (JgiJ ¼ 1), span(D) ¼ H and card(D)bN be an over-complete basis for H. D can be called over-complete dictionary. If a measured noisy ultrasonic pulse-echo can be expressed as follows:
g iðmÞ ¼
xðtÞ ¼ hðtÞ þ eðtÞ
(9)
where h(t) is the ultrasonic pulse-echo, e(t) can be characterized as white Gaussian noise (WGN). We look for decompositions of the signal x(t) on the overcomplete dictionary D of indexed unit-norm vector gi, in the form x¼
m X
ai g i þ rm
(10)
i¼1
where ai is the weight associated to the unit-norm vector gi. Our goal is to look for the appropriate number m to make the Pm term i ¼ 1ai gi be the approximation of the ultrasonic pulseecho h(t) with a higher quality [13]. Matching pursuit [11] is an iterative algorithm that offers a suboptimal solution for decomposing a signal f in terms of unitnorm vector gi chosen from an over-complete dictionary D. At the first iteration, the unit-norm vector gi which has the largest inner product with analyzed signal f is chosen. Then the contribution of this vector is then subtracted from the signal, and the process is repeated on the residual. At the mth iteration, the residue is ( f; m¼0 m (11) r ¼ rmþ1 þ aiðmÞ g iðmÞ ; ma0 where ai(m) is the weight associated with the optimum unit-norm vector gi(m) at the mth iteration. The weight ai(m) associated with each unit-norm vector giAD at the mth iteration is introduced to
hr m ; g i i hr m ; g i i ¼ ¼ hrm ; g i i hg i ; g i i kg i k2
(12)
The optimum unit-norm vector gi(m) (and its weight ai(m)) at the mth iteration can be obtained g iðmÞ ¼ arg min kr mþ1 k2 ¼ arg max jaiðmÞ j2 g i 2D
g i 2D
(13)
The correlation updating procedure is performed as follows: hr m ; g i i ¼ hrm ; g i i aiðmÞ hg iðmÞ ; g i i
(14)
The correlations /gi(m),giS can be precalculated and stored, once over-complete dictionary D has been determined [12]. The criterion chosen to stop the iteration is: for two given residuals rm and rm+1, if the error ratio of the mean square e ¼ (rm+1mrm)/mrm is less than a specified tolerance, then the matching pursuit stops. Where m ¼ (E((rm)2)/E((rm+1)2))1/2 and E( ) is the expected value. m makes the two given residuals rm and rm+1 have the same expected value, and it also can guarantee the convergence of the ratio e ¼ (rm+1mrm)/mrm in the mean square. The error criterion is used to monitor the acceptability of the unit-norm vector obtained by the matching pursuit performance. And the matching pursuit operation will continue until the error ratio of the mean square (rm+1mrm)/mrm is less than a predetermined threshold. The error criterion is an important factor for the pulse-echo approximation while noise level is unknown; furthermore, it is generally different for different ultrasonic signals. Compared to that of the residual rm, when the mean square of the term (rm+1mrm) can be neglected, the optimum unit-norm vector gi(m) at the mth iteration will make little contribution to the pulse-echo, it can be regarded as a partial signal come from the additive noise. Thus, the matching pursuit operation process can be terminated.
ARTICLE IN PRESS L. Wei et al. / NDT&E International 42 (2009) 430–434
According to the MP algorithm, our method has to define the unit-norm vectors in the dictionary D and the way to approximate pulse-echo. Since, a real time implementation is claimed for usefulness, the complexity of this method must be reduced as much as possible. Because MP is an iterative algorithm, the goal can be achieved if both the amount of operations at each iteration and the number of iterations to estimate pulse-echo are reduced. In ultrasonic non-destructive testing, the magnitude spectrum of the transducer pulse-echo has band-pass characteristics. Hence, the over-complete dictionary D can compose of unit-norm vectors gi,k that mimic the ultrasonic pulse-echo to be estimated, and in turn depend on the impulsive response of the transducer used in the experiments. The unit-norm vector in the dictionary D can be defined as the time-shifted Gaussian echo wavelet, and can be expressed as 1 g i;k ½n ¼ pffiffiffi exp½gk ðn ti Þ2 cosð2pf c ðn ti Þ jk Þ
ti ¼ 1; . . . ; L; L ¼ M þ N; n ¼ ti
M M þ 1; . . . ; ti þ 2 2
(15)
where gk, fc and jk are the Gaussian echo wavelet parameters. The pffiffiffi normalizes the energy of Gaussian echo wavelet, and termp1= ffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ p=2gk . They are related to the frequency and the bandwidth of the transducer response. L is the number of time shifts to take into account all possible flaw echo positions, N is the ultrasonic data length and M is the length of the unit-norm vectors gi,k.
5. Experimental results Many ultrasonic non-destructive testing applications require the accurate estimation of TOA or/and TOF of the ultrasonic echoes. In this paper, ultrasonic detection signals are analyzed in order to validate the proposed method. In the first experimental scheme, a steel sample block with approximately 1.0 cm thickness is tested using the unfocussed longitudinal transducer with 2.5 MHz central frequency. The ultrasound signals were recorded in the format of A-scan with a sample frequency of 100 MHz. Fig. 2(a) shows the measured
433
echoes reflected from the front and back surfaces of the steel sample block. Three echoes can be found in the signal. The first and second echoes are visible at the delay of approximately 1.13 and 4.56 ms, respectively while the third echo can be assumed to arrive at 8.04 ms. The TOFs of two neighboring ultrasonic pulseecho are 3.43, 3.48 ms, respectively. The average error is about 0.05 ms. It is shown from Fig. 2(a) that the ultrasound signal has a low time-resolution. Fig. 2(b) plots the time domain representation of the measured ultrasonic pulse-echo in solid line reflected from the inner surface of the offshore pipeline sample, and it also plots the reconstruction result in dotted line with the former fourteen large coefficient unit-norm vectors. It is shown in Fig. 2(b) that the reconstruction result achieves a reasonable quality. Fig. 2(c) shows the result of the sparse deconvolution method proposed in this paper. The delay between the first and the second echo is 3.44 ms, and the delay between the second and third echo is also 3.45 ms. The average error is only about 0.01 ms. If the sound velocity for longitudinal waves in ordinary steel is 5.9 mm/ms, the thickness of the steel sample block can be computed as 10.148 and 10.178 mm. The average error is only about 0.03 mm. This result makes evidence that the good performance of the sparse deconvolution method which is in much improved accuracy in the estimation of echo arrival time. In the second experimental scheme, an offshore pipeline specimen was applied in which a man-made 10 1 5.6 mm3 flaw was fabricated close to the outer wall. The experiment was done inside the pipe. The ultrasonic probe with 10 MHz central frequency and 6 mm diameter was used and water was selected as couplant. The flaw echoes were recorded in the format of A-scan with a sample frequency of 100 MHz. Fig. 3(a) shows the measured ultrasonic signal. Fig. 3(b) plots the time domain representation of the measured ultrasonic pulse-echo reflected from the inner surface of the offshore pipeline sample (in solid line) and the reconstruction result in dotted line with the former fourteen large coefficient unit-norm vectors. Fig. 3(c) shows the deconvolution result of the sparse deconvolution method. It is shown from Fig. 3(c) that the major peaks are visible in the deconvolution result. Three TOFs, 0.81, 2.95, and 5.01 ms can be identified clearly in Fig. 3(c). The depth of the flaw can be calculated: 6.077 mm.
Fig. 2. The sparse deconvolution method for the estimation of echo arrival time: (a) the measured ultrasonic pulse-echo reflected from the front surface of the steel sample block; (b) comparison between the approximation result (dotted line) and the measured ultrasonic pulse-echo (solid line); (c) deconvolution result using the sparse deconvolution method.
ARTICLE IN PRESS 434
L. Wei et al. / NDT&E International 42 (2009) 430–434
Fig. 3. The sparse deconvolution method for the estimation of echo arrival time: (a) the measured ultrasonic signal received from a 10 1 5.6 mm3 flaw close to the outer wall; (b) comparison between the approximation result (dotted line) and the measured ultrasonic pulse-echo (solid line); (c) deconvolution result using the sparse deconvolution method.
6. Conclusion Many ultrasonic non-destructive testing applications require the accurate estimation of TOA or/and TOF of the ultrasonic echoes. However, the received signals might have been contaminated by noise originating from both the measurement system and the specimen. In these cases, it is impossible to obtain the accurate TOA or TOF from the received signals. Based the characteristics of the transducer pulse-echo, the over-complete dictionary D has been defined. We use the matching pursuit method to analyze noisy ultrasonic pulse-echo and decompose the noisy ultrasonic pulse-echo into elementary components. Because, the over-complete dictionary D is composed of the unitnorm vectors that mimic the transducer impulse response, the complexity of the matching pursuit method is reduced. The method is very efficient in eliminating the noise and the reconstruction result achieves a reasonable quality. Using the reconstruction transducer pulse-echo wavelet, a sparse deconvolution method is proposed for deconvolving ultrasonic signals. A weighting matrix is central to the efficiency of the sparse deconvolution method. It increases the sparsity of the deconvolution result and decreases the influence of the added noise. The deconvolution gives high-resolution ultrasonic reflectivity function. We can obtain the accurate TOA or TOF from the ultrasonic reflectivity function. Results using real ultrasonic data show the sparse deconvolution method that yields good estimates of the spikes associated with the reverberation echoes and the thickness of a steel sample can be calculated by the ultrasonic reflectivity function with a reasonable accuracy.
Acknowledgement This work was supported by the China National High-Tech Research and Development Program, the National Science
Foundation for Young Scientists of China and the JIN CHUAN Advance Research Program of Science and Technology. References [1] Nandi AK, Mampel D, Roscher B. Blind deconvolution of ultrasonic signals in non-destructive testing applications. IEEE Transactions on Signal Processing 1997;45:1382–9. [2] Cassereau D, Guyomar D, Fink M. Time deconvolution of diffraction effects—application to calibration and prediction of transducer waveforms. Journal of the Acoustical Society of America 1988;84:1073–85. [3] Olofsson T. Deconvolution and model-based restoration of clipped ultrasonic signals. IEEE transactions on instrumentation and measurement 2005; 54:1235–40. [4] Lee JY, Nandi AK. Extraction of impacting signals using blind deconvolution. Journal of Sound and Vibration 2000;232:945–62. [5] Wiggins RA. Minimum entropy deconvolution. Geophysical Exploration 1978;16:21–35. [6] Nandi Asoke K. Blind identification of FIR systems using third order cumulates. Signal Processing 1994;39:131–47. [7] Nandi AK, Mehlan R. Parameter estimation and phase reconstruction of moving average processes using third-order cumulants. Mechanical Systems and Signal Processing 1994;8:421–36. [8] Gorodnitsky Irina F, Rao Bhaskar D. Sparse signal reconstruction from limited data using FOCUSS: a re-weighted minimum norm algorithm. IEEE Transactions on Signal Processing 1997;45:600–16. [9] Liu HH, Gao XR, Schimpf Paul H. A recursive algorithm for the threedimensional imaging of brain electric activity. Shrinking LORETA–FOCUSS 2004;51:1794–802. [10] Liang W, Lei HM, Que PW. Sparsity enhancement for blind deconvolution of ultrasonic signals in nondestructive testing application. Review of Scientific Instruments 2008;79:014901. [11] Mallat S, Zhang Z. Matching pursuits with time-frequency dictionaries. IEEE Transactions on Signal Processing 1993;41:3397. [12] Ruiz-Reyes N, Vera-Candeas P, Curpia´n-Alonso J. New matching pursuit-based algorithm for SNR improvement in ultrasonic NDT. NDT & E International 2005;38:453. [13] Liang W, Que PW, Lei HM, et al. Matching pursuit for decomposition and approximation of ultrasonic pulse-echo wavelet and its application in ultrasonic nondestructive evaluation. Review of Scientific Instruments 2008;79(7):075105.