Generalization of iterative Fourier interpolation algorithms for single frequency estimation

Generalization of iterative Fourier interpolation algorithms for single frequency estimation

Digital Signal Processing 21 (2011) 141–149 Contents lists available at ScienceDirect Digital Signal Processing www.elsevier.com/locate/dsp General...

239KB Sizes 1 Downloads 23 Views

Digital Signal Processing 21 (2011) 141–149

Contents lists available at ScienceDirect

Digital Signal Processing www.elsevier.com/locate/dsp

Generalization of iterative Fourier interpolation algorithms for single frequency estimation Yanhui Liu a,∗ , Zaiping Nie a , Zhiqin Zhao a , Qing Huo Liu b a b

Department of Electronic Engineering, University of Electronic Science and Technology of China, Sichuan 610054, China Department of Electrical and Computer Engineering, Duke University, Durham, NC 27708, USA

a r t i c l e

i n f o

a b s t r a c t

Article history: Available online 1 July 2010

The problem of estimating the frequency of a complex single tone is considered. We generalize two iterative Fourier interpolation algorithms in the literature by introducing an additional parameter to allow for selection of the Fourier interpolation coefficients relative to the true frequency. The generalized algorithms can use more information from FFT results and consequently achieve significant improvement over the original algorithms in either accuracy or efficiency. Simulation results show advantages of the proposed algorithms. © 2010 Elsevier Inc. All rights reserved.

Keywords: Frequency estimation Fourier interpolation Discrete Fourier transform (DFT) Fast Fourier transform (FFT)

1. Introduction The estimation of the frequency of a complex exponential in noise has applications in many areas, such as communications, radar and sonar. For the case of uniformly spaced data, the maximum likelihood (ML) estimation of a single frequency can be obtained by locating the peak of a periodogram [1]. The ML estimator has better threshold property but requires more computation, compared with many time-domain approaches [2]. To improve computational efficiency of the ML estimator, several fast frequency-domain search algorithms have been developed previously, which usually consist of two steps: a coarse search and a fine search [1]. The coarse search returns a coarse estimate of the frequency by searching for the maximum FFT coefficient of the signal. This coarse estimate can be improved by using some fine estimation algorithms, such as nonlinear optimization approaches [3–5] and Fourier interpolation algorithms [6–13]. Compared with the optimization approaches, Fourier interpolation algorithms are usually more efficient since they use prior knowledge about the signal model and consequently require no or fewer calculations of nonstandard DFT coefficients. Quinn, in [7,8] and [9], proposed a number of Fourier interpolators which all interpolate the signal frequency using the maximum FFT coefficient and its two adjacent coefficients. Macleod, in [11] modified one of the three-sample interpolators and also presented a five-sample interpolator for performance improvement. Recently, Aboutanios and Mulgrew, in [12], proposed two iterative Fourier interpolation algorithms (A&M algorithms). One of them uses two complex DFT coefficients located in the middle of neighbor FFT coefficients, and another works on the moduli of the DFT coefficients. The two algorithms both require two iterations to converge with a variance that is asymptotically marginally above the Cramer–Rao bound (CRB) in the entire frequency estimation range. We generalize the A&M iterative interpolation algorithms by introducing an additional parameter to allow for selection of the Fourier interpolation coefficients relative to the true frequency. The two generalized algorithms are both capable of employing more information from the FFT results and consequently improve the performance of the original algorithms in either accuracy or efficiency, depending on the choice of the additional parameter. With a specific choice for this parameter,

*

Corresponding author. E-mail address: [email protected] (Y. Liu).

1051-2004/$ – see front matter doi:10.1016/j.dsp.2010.06.012

©

2010 Elsevier Inc. All rights reserved.

142

Y. Liu et al. / Digital Signal Processing 21 (2011) 141–149

the proposed algorithms save calculations of two DFT coefficients, which would be very useful in real time applications. Numerical results are reported to validate the effectiveness and advantages of the proposed method. It is emphasized that this work deals only with single frequency estimation of uniformly spaced data in additive Gaussian noise. Some early studies in relevant applications include [15,16] in Lomb periodogram estimation for processing nonuniformly spaced data, and [17] in the least-squares solution for the case of colored noise. 2. Methodology The complex single tone under consideration is of the form

z(k) = e

j2π k

f fs



+ w (k),

k = 0, 1 , . . . , N − 1 ,

(1)

where f is the signal frequency, θ the initial phase, N the size of samples, and f s the sampling frequency. The noise terms w (k) is assumed to be a zero-mean, complex additive white Gaussian process with variance σ 2 . The problem we consider here is to estimate the frequency f that can be written as

f =

m N + δN N

(2)

f s,

where m N indicates the index of the bin that is closest to the true frequency, and δ N is a residual in the interval [−0.5, 0.5]. The subscript N indicates the dependence of the parameters on N. But for simplifying the notation, we will omit this ˆ = m. Thus subscript in the rest of this paper. Assume the coarse search returns the index of the true maximum bin, i.e., m the goal is to estimate the residual δ . The A&M algorithms [12] considered the two DFT coefficients calculated midway between the standard DFT bins

Z ± 0 .5 =

N −1 

z(k)e − j2π k

ˆ ±0.5 m N

(3)

.

k =0

The authors in [12] then derived two Fourier interpolators given below

δˆ = δˆ =





Z 0 .5 + Z − 0 .5 , 2 Z 0 .5 − Z − 0 .5 1 | Z 0 .5 | − | Z − 0 .5 | 1

(4)

Re

2 | Z 0 .5 | + | Z − 0 .5 |

(5)

.

Both interpolators have the same ratio of the asymptotic variance to the CRB, which is given by [12]

R=

π 4 (δ 2 − 0.25)2 (4δ 2 + 1) . 6 cos2 (π δ)

(6)

The dependence of the variance on δ can be eliminated by an iterative procedure as described in Table I of [12]. In this correspondence, we shall show that both interpolators can be generalized by introducing an additional parameter q which allows for selection of DFT coefficients relative to the true frequency. This generalization gives us more freedom to achieve performance improvement. 2.1. Generalization of A&M Fourier interpolators We shall now consider the following DFT coefficients instead

Z q ± 0 .5 =

N −1 

z(k)e − j2π k

ˆ +q±0.5 m N

(7)

,

k =0

where the parameter q is introduced for adjusting the location of the DFT coefficients. Considering (1) and (2), we can rewrite (7) as

Z q ± 0 .5 =

N −1 

z (k)e − j2π k

ˆ ±0.5 m N

(8)

,

k =0

where

z (k) = e j2π k

m+(δ−q) +θ N

q

+ w (k)e − j2π k N .

Clearly, (8) has the same form as (3). In addition, through a shift δ 

(9)

= δ − q, the signal of (9) also has the same form as the original signal described in (1) and (2). Hence, we can follow exactly the same way as in [12] to arrive at

Y. Liu et al. / Digital Signal Processing 21 (2011) 141–149

δˆ = δˆ =

1 2

 Re

Z q + 0 .5 + Z q − 0 .5 Z q + 0 .5 − Z q − 0 .5

1 | Z q + 0 .5 | − | Z q − 0 .5 | 2 | Z q + 0 .5 | + | Z q − 0 .5 |

143

 + q,

+ q.

(10) (11)

It should be noted that the interpolator (11) is valid only under the condition that |δ − q|  0.5. For the case |δ − q|  0.5, such modulus-based interpolators also exist but in other forms that are not given here since in this case the DFT coefficients selected are far away from the signal frequency and contain less signal energy, which degrades the performance for noisy data. Since the noise terms w (k) exp(− j2π kq/ N ) of (9) holds in the same property as w (k), actually we can also arrive at all the performance expressions by simply replacing δ in [12] with (δ − q). Finally, the ratio of the asymptotic variance to the CRB is given by

R=

π 4 [(δ − q)2 − 0.25]2 [4(δ − q)2 + 1] . 6 cos2 [π (δ − q)]

(12)

Note that for the interpolator (10) this ratio asymptotically holds for any q, while for (11) this ratio holds only under the condition |δ − q| < 0.5. 2.2. Choice of q Now consider the problem of how to choose q for (10) and (11). From (12), the closer the parameter q is to δ , the better the performance is. Obviously, the optimum q would be equal to δ (R ≈ 1.1047 in this case). However, δ is of course unknown in practice. A good choice is setting q = 0, so that |δ − q|  0.5 be guaranteed. To achieve a better performance, an elaborate alternative way of choosing q is arrived at by judging the sign of δ using FFT coefficients, as in Quinn’s interpolators [7–9]. Quinn in [7] proposed an estimate of sign(δ ) in order to choose a better one of two interpolators. Macleod [11] studied the performance of Quinn’s test and pointed out that this test is not optimal if the SNR is low and δ is sufficiently close to 0. Instead, he suggested that a better estimate of sign(δ ) is obtained by testing if Re( Z −1 Z 0 ∗ ) − Re( Z 1 Z 0 ∗ ) > 0, as described in (21) of [11]. For clarity, we put the explanation for this test in Appendix A. Finally we suggest the way of determining q as follows





q = β sign Re ( Z −1 − Z 1 ) Z 0 ∗

 ,

(13)

where β ∈ [0, 0.5]. Asymptotically, the ratio R is given by

R=

π 4 [(|δ| − β)2 − 0.25]2 [4(|δ| − β)2 + 1] . 6 cos2 [π (|δ| − β)]

(14)

Let us consider some special cases of the interpolators (10) and (11).

• Case 1: β = 0. The interpolators (10) and (11) are just the A&M interpolators (4) and (5), respectively. • Case 2: β = 0.25. As indicated by (14), the interpolators (10) and (11) in this case both have the minimum of the maximum variance over δ ∈ [−0.5, 0.5]. This will be validated also by the simulations in Section 3. • Case 3: β = 0.5. The interpolators (10) and (11) become

δˆ = Re



α Zα Zα − Z0

 (15)

and

δˆ =

α| Z α | | Zα | + | Z0|

(16)

respectively, where α = sign{Re[( Z −1 − Z 1 ) Z 0 ∗ ]}. The interpolators (15) and (16) are quite similar to Quinn’s interpolator [7] and the modified Rife and Vincent algorithm [8], respectively. (15) and Quinn’s interpolator have the same asymptotic variance [(14) at β = 0.5]. However, in low SNR, (15) is expected to perform better than Quinn’s interpolator since the polarity test used in (15) is more reliable in this case, as pointed out by Macleod [11]. Also, the same relative advantage can be obtained if we compare (16) and the modified Rife and Vincent algorithm [8]. Such performance improvements will be validated by the simulations in Section 3. Note that compared with Case 1, Case 2 has a lower estimation variance, and Case 3 is much less time-consuming since the DFT coefficients in (15) or (16) are just the FFT coefficients and do not need to be recalculated.

144

Y. Liu et al. / Digital Signal Processing 21 (2011) 141–149

2.3. Iterative Fourier interpolation algorithms A similar iterative interpolation procedure as that used in [12] is introduced here to achieve uniform estimation performance in the whole frequency estimation range. The new iterative Fourier interpolation algorithms are given as follows:

ˆ = arg maxn {| Z (n)|2 }, where n = 0, 1, . . . , N − 1. 1. Let Z = FFT ( z) and m 2. Choose β ∈ [0, 0.5] and let δˆ0 = β sign{Re[( Z −1 − Z 1 ) Z 0 ∗ ]}. 3. Loop: for each i from 1 to Q ( Q  2), do

δˆi =

1

δˆi =

1

2

Z

Re

Z δˆ

i − 1 + 0 .5

or

2

δˆi −1 +0.5

Re

 |Z

+ Z δˆi−1 −0.5  − Z δˆi−1 −0.5

+ δˆi −1 for New Agl.1

δˆi −1 +0.5 | − | Z δˆi −1 −0.5 |

| Z δˆi−1 +0.5 | + | Z δˆi−1 −0.5 |



+ δˆi −1 for New Agl.2.

ˆ +δˆ m 4. Finally ˆf = N Q f s .

The above iterative algorithms are quite similar to the A&M algorithms [12] but with more generalization which allows for selection of the DFT coefficients at the first iteration. The dominant time cost comes from the N-point FFT and the calculations of four interpolated DFT coefficients (when Q = 2). Therefore the proposed algorithms both require about 0.5N log2 N + 4N complex multiplications and N log2 N + 4N complex additions (suppose that the radix-2 algorithm [14] is used for the FFT implementation). This computational load is obviously the same as that of the A&M algorithms. However, for a specific choice of β = 0.5, the proposed algorithms do not need to calculate the two DFT coefficients at the first iteration, which saves 2N complex multiplications and 2N − 2 complex additions. This saving is very useful, especially in fast data-rate applications. 3. Simulation results Two simulation experiments, with N = 1024 and N = 16, were conducted. Both experiments used f s = 1 without loss of generality in the results. For each experiment, we scanned the frequency offset δ from −0.5 to 0.5 with increment of 0.025 and also varied the SNR with increment of 1 dB. The proposed Fourier interpolation algorithms with different choices of β were used to estimate frequency. In addition, several other interpolation algorithms in [12,7] and [8], and the ML estimator as well, were also used for comparison. All these interpolation algorithms can be classified as the following two types: (1) The first-type algorithms (based on complex DFT coefficients) • A&M: the first algorithm of Aboutanios and Mulgrew [12]; • β = 0.25: the proposed first algorithm with β = 0.25; • β = 0.5: the proposed first algorithm with β = 0.5; • Quinn, [7]: Quinn’s interpolator presented in [7]. (2) The second-type algorithms (based on moduli of DFT coefficients) • A&M: the second algorithm of Aboutanios and Mulgrew [12]; • β = 0.25: the proposed second algorithm with β = 0.25; • β = 0.5: the proposed second algorithm with β = 0.5; • Quinn, [8]: the modified Rife and Vincent algorithm presented by Quinn in [8]. Note that both A&M algorithms are special cases of the proposed algorithms when β = 0. In the rest of this paper, we will show the simulation results of the first-type algorithms in Figs. 1–4(a), and the second-type algorithms in Figs. 1–4(b), respectively. The results for the ML estimator are also included for comparison. Fig. 1 shows the mean squared frequency errors (MSE) averaged over δ for N = 1024 at SNR from −14 dB to 0 dB. First, consider the first iteration ( Q = 1). We can see, in Fig. 1(a), that compared with the first A&M algorithm [12], the proposed first algorithm with β = 0.5 has almost the same estimation accuracy (but requires no calculations of the DFT coefficients), and this algorithm with β = 0.25 gives much better accuracy (its performance curve is already very close to the CRB). Also in Fig. 1(b), the same relative performance is observed between the proposed second algorithm and the second A&M algorithm. We can also compare the proposed algorithms with Quinn’s algorithms. In Fig. 1(a), Quinn’s interpolator in [7] obtains the same accuracy as the proposed first algorithm with β = 0.5 at high SNR but deteriorates at low SNR. And in Fig. 1(b), Quinn’s second algorithm (the modified R&V algorithm in [8]) is also less accurate than the proposed second algorithm, especially at low SNR. Then let us check the performance at the second iteration ( Q = 2). As can be seen, there is no any visible difference in the final accuracy performance between the proposed two algorithms (for both β = 0.25 and β = 0.5) and the A&M algorithms. The MSEs for all these algorithms at the second iteration are approximately equal to the

Y. Liu et al. / Digital Signal Processing 21 (2011) 141–149

145

(a)

(b) Fig. 1. The estimator MSE averaged over δ for N = 1024. (a) and (b) show the results of the first-type algorithms and the second type algorithms, respectively. Q denotes the number of iterations. There were an average of 2000 trials.

CRB. It should be noted that the threshold effect of the proposed algorithms come from the coarse search stage, which is characteristic of the ML estimator [1]. Remember all the results in Fig. 1 are shown in the average performance over δ . We also check the performance at different δ in detail for the proposed algorithms and other algorithms. Fig. 2 shows the ratio of the estimation MSE to the CRB over δ ∈ [−0.5, 0.5] at SNR = 0 dB. The results at the first iteration shed some light on how the choice of β affects the accuracy of the proposed interpolators, which also confirm the conclusion in the relative average performance we made from Fig. 1. The results at the second iteration confirm that our generalization of the A&M algorithms still converges with a final performance uniformly distributed over the whole range of δ . In addition, Fig. 3 shows the situation at SNR = −13 dB that is (close to) the threshold SNR of the ML estimator. We can see that at such SNR, the proposed algorithms still maintain almost the same performance relative to the CRB, except for a slight degradation in a range of δ close to 0. At such SNR, this degradation should be acceptable. However, Quinn’s two algorithms unfortunately get worse at such SNR condition, especially for δ > 0. This is because the polarity test used in these two algorithms is not reliable enough at low SNR, as mentioned previously. Another simulation experiment is for checking the performance of the proposed algorithm for small N. Note that in this simulation we artificially removed the threshold of the coarse search by using the true maximum bin index for all the tested algorithms, in order to check the performance of their fine interpolators at an SNR below the threshold of the ML estimator. Definitely, this procedure does not affect the performance at high SNR. Fig. 4 shows the results on the MSE averaged over δ for N = 16 at SNR varying from −2 dB to 14 dB. As can be seen, the proposed algorithms for this small N maintain

146

Y. Liu et al. / Digital Signal Processing 21 (2011) 141–149

(a)

(b) Fig. 2. The ratio of the estimator MSE to the CRB as a function of δ at SNR = 0 dB. (a) and (b) show the results for the first-type algorithms and for the second type algorithms, respectively. The legend is the same as that used in Fig. 1. There were an average of 5000 trials.

almost the same performance relative to the CRB as them for N = 1024. We also see that compared with Quinn’s two algorithms, the proposed interpolaters show much lower thresholds. It should be noted that we did many other simulations with N varying from 16 to 4096, and found that for larger N, the threshold of the proposed interpolators themselves can be much lower than that of the coarse search stage (the characteristic of ML estimator [1]). This property is useful for some applications where prior information about the signal frequency is available and the search for the maximum FFT coefficients can be therefore restricted to a small frequency range [11]. 4. Conclusion and discussion In this paper, we have proposed a generalization for the two state-of-the-art iterative Fourier interpolation algorithms (A&M algorithms [12]) by introducing an additional parameter. By choosing different values of β , the proposed algorithms allow for selecting different DFT coefficients to calculate the signal frequency. With β = 0, the proposed algorithms reduce to the original A&M algorithms. With β = 0.25, the proposed algorithms have much better accuracy performance than that of the original A&M algorithms at the first iteration. This accuracy is even already acceptable. With β = 0.5, the proposed algorithms maintain almost the same performance as the A&M algorithms in both accuracy and threshold but save the calculations of two DFT coefficients at the first iteration. This leads to saving of 2N complex multiplications and 2N − 2

Y. Liu et al. / Digital Signal Processing 21 (2011) 141–149

147

(a)

(b) Fig. 3. The ratio of the estimator MSE to the CRB as a function of δ at SNR = −13 dB. (a) and (b) show the results for the first-type algorithms and for the second type algorithms, respectively. The legend is the same as that used in Fig. 1. There were an average of 5000 trials.

complex additions. In addition, the proposed algorithms (with any β ∈ [−0.5, 0.5]) are much more accurate and more robust than the Quinn’s two algorithms at low SNR or for small N. Acknowledgments This work was supported in part by U.S. National Science Foundation under Grant CCF-0621862, in part by the 111 project of China under Contract No. B07046, in part by the Natural Science Foundation of China under Grant No. 60431010, and in part by the Joint Ph.D. Fellowship Program of the China Scholarship Council. Appendix A. Explanation of polarity test We give the explanation for the polarity test of δ in (13). Substituting (1) and (2) into (7) and carrying out some necessary manipulations (the approximation that e x ≈ 1 + x for |x|  1 is used), we obtain

Z q ± 0 .5 = b q ± 0 .5 where bq is given by

δ + W q ± 0 .5 , (δ − q) ∓ 0.5

(17)

148

Y. Liu et al. / Digital Signal Processing 21 (2011) 141–149

(a)

(b) Fig. 4. The estimator MSE averaged over δ for N = 16. (a) and (b) show the results of the first-type algorithms and the second type algorithms, respectively. Q denotes the number of iterations. There were an average of 2000 trials.

bq±0.5 = − Ne j θ

1 + e j2π (δ−q) j2π δ

and W q±0.5 are the noise DFT coefficients. By setting p = ±0.5, the coefficients Z l (l = −1, 0, 1) is given by

Z l ≈ b0

δ . δ −l

(18)

Hence, we have





Re ( Z −1 − Z 1 ) Z 0 ∗ ≈ 2|b0 |2

δ 1 − δ2

.

Since 1 − δ 2 > 0 for |δ|  0.5, the sign of (19) is asymptotically equal to sign(δ ). References [1] D.C. Rife, R.R. Boorstyn, Single-tone parameter estimation from discrete-time observations, IEEE Trans. Inform. Theory 20 (1974) 591–598. [2] M. Morelli, U. Mengali, Feedforwad frequency estimation for PSK: a tutorial review, Eur. Trans. Commun. Related Technol. 9 (1998) 103–116.

(19)

Y. Liu et al. / Digital Signal Processing 21 (2011) 141–149

149

[3] T.J. Abatzoglou, A fast maximum likelihood algorithm for frequency estimation of a sinusoid based on Newton’s method, IEEE Trans. Acoust. Speech Signal Process. 33 (1985) 77–78. [4] Y.V. Zakharov, T.C. Tozer, Frequency estimator with dichotomous search of periodogram peak, Electron. Lett. 35 (1999) 1608–1609. [5] E. Aboutanios, A modified dichotomous search frequency estimator, IEEE Signal Process. Lett. 11 (2004) 186–188. [6] D.C. Rife, G.A. Vincent, Use of the discrete Fourier transform in the measurement of frequencies and levels of tones, Bell Syst. Tech. J. 49 (1970) 197–228. [7] B.G. Quinn, Estimating frequency by interpolation using Fourier coefficients, IEEE Trans. Signal Process. 42 (1994) 1264–1268. [8] B.G. Quinn, Estimation of frequency, amplitude, and phase from the DFT of a time series, IEEE Trans. Signal Process. 45 (1997) 814–817. [9] B.G. Quinn, E.J. Hannan, The Estimation and Tracking of Frequency, Cambridge University Press, New York, 2001. [10] B.G. Quinn, Recent advances in rapid frequency estimation, Digital Signal Process. 19 (6) (2008) 942–948, doi:10.1016/j.dsp.2008.04.004. [11] M.D. Macleod, Fast nearly ML estimation of the parameters of real or complex single tones or resolved multiple tones, IEEE Trans. Signal Process. 46 (1998) 141–148. [12] E. Aboutanios, B. Mulgrew, Iterative frequency estimation by interpolation on Fourier coefficients, IEEE Trans. Signal Process. 53 (2005) 1237–1241. [13] I. Perisa, J. Lindner, Employing simple FFT-interpolation for improved complex tone detection and fine estimation, in: 3rd ISWCS, 2006, pp. 744–748. [14] A.V. Oppenheim, R.W. Schafer, Digital Signal Processing, Prentice Hall, Englewood Cliffs, 1975. [15] N.R. Lomb, Least-square frequency analysis of unequally spaced data, Astrophys. Space Sci. 39 (1976) 447–462. [16] J.D. Scargle, Studies in astronomical time series analysis. II. Statistical aspects of spectral analysis of unevenly spaced data, Astrophys. J. 263 (1982) 835–853. [17] P. Stoica, A. Jakobsson, Jian Li, Cisoid parameter estimation in the colored noise case: Asymptotic Cramer–Rao bound, maximum likelihood, and nonlinear least-squares, IEEE Trans. Signal Process. 45 (1997) 2048–2059.

Yanhui Liu received the B.S. degree in electronic engineering from the University of Electronic Science and Technology of China (UESTC) in 2004, where he is currently working toward the Ph.D. degree. From Sep. 2007 to Jun. 2009, he was a Visiting Researcher in the Department of Electrical Engineering, Duke University. His research interests include spectral analysis and frequency estimation, radar imaging, electromagnetic system design. Zaiping Nie received the B.S. degree in radio engineering and the M.S. degree from the Chengdu Institute of Radio Engineering (now UESTC) in 1968 and 1981, respectively. From 1987 to 1989, he was a Visiting Scholar with the Electromagnetics Laboratory, University of Illinois. Currently, he is a Professor with the Department of Electromagnetic Engineering, UESTC. His research interests include electromagnetic propagation and scattering, antenna theory and technique. Zhiqin Zhao received his Ph.D. degrees in electronic engineering from the University of Electronic Science and Technology of China (UESTC) in 1996, and his second Ph.D. degree in electrical engineering from Oklahoma State University in 2002. From Jan. 2003 to Mar. 2006, he worked as a Post-doctoral Research Associate, later a Research Scientist, with Duke University. Since Apr. 2006, he has joined UESTC as a professor. His current research interests include signal processing, remote sensing and biomedical imaging. Qing Huo Liu received the Ph.D. degree in electrical engineering from the University of Illinois at Urbana-Champaign in 1989. He was a Research Scientist and Program Leader with Schlumberger-Doll Research from 1990 to 1995, and then was an Associate Professor with New Mexico State University from 1996 to May 1999. Since June 1999, he has been with Duke University where he is now a Professor of Electrical and Computer Engineering. His research interests include computational electromagnetics and acoustics, inverse problems, biomedical imaging and electronic packaging.