Document Sample

Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 1 Tampere University of Technology 2 POLYNOMIAL-BASED INTERPOLATION 80558 MULTIRATE SIGNAL PROCESSING FILTERS FOR DSP APPLICATIONS Tapio Saramäki DESIGN, IMPLEMENTATION, AND APPLICATIONS Part III: Polynomial-Based Interpolation for DSP Applications Tapio Saramäki and Jussi Vesma • This pile of lecture notes is mainly based on the research Tampere University of Technology, Finland work done by Jussi Vesma and the lecturer during the last e-mail: saram@vip.fi, or ts@cs.tut.fi three years. http://www.cs.tut.fi/~ts/ • Many thanks to Jussi Vesma for his help in preparing this pile of lecture notes. Contents: • If there is some interest to get the pile of lecture notes of the overall course, please contact Tapio Saramäki using his 1. Interpolation Filters home e-mail adress: saram@vip.fi. 2. Fractional-Delay Filters • The overall course consists of the following parts: 3. Lagrange Interpolation • Part I: Basics of Multirate DSP 4. Analog Model for the Interpolation Filter 5. Polynomial-Based Interpolation Filters • Part II: Design and Implementation of Efficient Decimators and Interpolators 6. Filter Synthesis • Part III: Polynomial-Based Interpolation for DSP Applica- 7. Applications tions 8. Filter Properties • Part IV: Design of FIR Filters Using Multirate DSP and Comment: The following article coming up soon: J. Vesma Complementary Filtering and T. Saramäki, “Polynomial-based Interpolation Fil- • Part V: Multirate Filter Banks Including Discrete-Time − ters−Part I: Filter Synthesis; Part II: Filter Properties and Wavelet Banks Applications. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 3 Tampere University of Technology 4 Interpolation Filters Applications for Interpolation Filters • In many DSP applications there is a need to know the values • Timing adjustment in all-digital receivers (symbol synchro- of the signal also between the exiting discrete-time samples nization) x(n) as shown in Figure 1. • Time delay estimation • Special interpolation filters can be used to compute new • Conversion between arbitrary sampling frequencies sample values y(l) =ya(tl) at arbitrary points t l =(nl +µ l)Tin between the existing samples x(nl ) and x(nl +1). Here, Tin is the • Echo cancellation sampling period. • Phased array antenna systems • Here, ya(t) approximates either the original continuous-time • Speech coding and synthesis signal xa(t) or the signal obtained with the aid of the existing • Derivative approximation of discrete-time signals discrete-time samples x(n) using the sinc interpolation. • The output sample time is determined by nl Tin, the location of • Computer simulation of continuous-time systems the preceding existing sample, and the fractional interval • ML symbol timing recovery in digital receivers µl∈[0,1), the difference between t l and nl Tin as a fraction of Tin. ya(t) = x(n) = y(l) µlTin (nl−3)Tin (nl−2)Tin tl−1 (nl−1)Tin nlTin tl (nl+1)Tin (nl+2)Tin tl+1 (nl+3)Tin Time t + Fig. 1. Interpolation in the time domain. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 5 Tampere University of Technology 6 Discrete-Time Interpolation Interpretation of the Convolution Sum • First, the location of the existing sample preceding or occur- x(n) Interpolation filter y(l) ring at the new sampling instant tl is determined and denoted by nl. h(k, µl) • Based on the location of this sample, the existing samples nl µl located at n =nl −N/2+1, nl −N/2+2,···, nl +N/2 are used. Fig. 2. Simplified block diagram for the interpolation filter. • This means that there are N/2 existing discrete-time samples before and after the desired new time instant tl. • We concentrate on performing the interpolation problem • This gives the best results and explains why N has been se- stated on Page 3 using the system of Figure 2. lected to be an even integer. • The input parameters nl and µ l are used to determine the time • Second, the distance between tl and nlTin is measured as a instant tl for the l th output sample y(l) =ya(tl) as t l = (nl +µ l)Tin. fraction of Tin, giving the fractional interval µ l. • Given tl, these parameters are determined by • In the above convolution, the value of µl determines the im- pulse-response coefficient values h(k, µl). nl = t l / Tin and µ l = tl / Tin − tl / Tin (1) • The purpose is to determine them in such a way that y(l) is • After knowing nl and µ l, the interpolation filter calculates found according to some criterion to be discussed next. y(l) according to the following convolution: N / 2 −1 y (l ) = ∑ x(nl − k )h(k , µl ) (2) k =− N / 2 where N (even) is the filter length and h(k, µ l) is the discrete- time impulse response of the interpolation filter. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 7 Tampere University of Technology 8 Statement of the Interpolation Problem Comments Given N, find the impulse-response coefficients h(k, µ l) for • There exist the following two trivial solutions: k= −N/2+1, −N/2+2,···, N/2 to meet the following two condi- 1. If the ratio is an integer or a ratio of small integers, a con- tions: ventional discrete-time sampling rate alteration is efficient. 1. Optimize them such that y(l) = ya((nl +µl)Tin) for all values of µ l ∈[0,1), where ya(t) approximates according to some time- 2. If the ratio is not a ratio of small integers or an integer, then domain or frequency-domain criterion the signal one can first generate a continuous-time signal with the aid x a (t ) = ∑n = −∞ x(n) sin[π (t − nTin ) / Tin ] [π (t − nTin ) / Tin ] . ∞ of a D/A converter and a reconstruction filter and, then, re- sample this signal using an A/D converter and an anti- 2. The overall system can be implemented digitally using an aliasing filter. efficient structure. • The drawback in the second approach is the fact that if the • In Condition 1, the frequency-domain criteria are usually conversion should be performed with a high accuracy, then preferred for DSP applications. very costly components are needed. • Therefore, it is worth studying whether there exist an effi- cient technique to accomplish the same directly digitally. • This is especially true since all the information is carried by the existing samples. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 9 Tampere University of Technology 10 Various Approaches to Solve the Stated Problem Various Approaches to Solve the Stated Problem • The problem is now to design an interpolation filter that cal- • There exist the following three approaches to solve our prob- culates the interpolated samples y(l) for a given value of µ l lem: based on the above discussion. 1. Fractional delay (FD) filter approach. • Furthermore, there is a need to analyze the performance of 2. Use some classical interpolation method to calculate y(l), the filter, that is, to measure the error between ya(t) and xa(t) e.g., Lagrange or B-spline interpolation (time-domain ap- either in the time domain or frequency domain. proach). • This is not an easy task since the interpolation filter is a 3. Utilize the analog model for the interpolation filter (fre- time-varying system, that is, the coefficients h(k, µl) depend quency-domain approach). on the value of µl. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 11 Tampere University of Technology 12 Fractional Delay (FD) Filter Approach Example FD FIR Filters • There exist several synthesis methods for designing FIR fil- • N =8, passband region =[0, 0.75π], ten filters having 4−µ ters with the transfer function of the form (N is even) with µ = 0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6 ,0.7, 0.8 ,0,9, and (1; N not needed) H ( z , µ ) = ∑ h( k , µ ) z − k k =0 1.1 1 for various values of µ such that in the given passband re- 0.9 gion [0, ωp] 0.8 Amplitude Response 0.7 1. The phase delay approximates Dint−µ, where Dint =N/2 is the 0.6 integer delay and µ∈[0,1) is the fractional delay. 0.5 0.4 2. The amplitude response approximates unity. 0.3 • This approach provides a good solution if µ is fixed. 0.2 0.1 • Otherwise, several FD FIR filters have to be designed for 0 0 0.1π 0.2π 0.3π 0.4π 0.5π 0.6π Angular Frequency ω 0.7π 0.8π 0.9π π numerous values of µ and coefficients have stored in a Fig. 3. Amplitude responses for the example FD filters. lookup table. • This means that if several values of µ are needed to provide 4 a good interpolation, a very large memory is needed. 3.9 3.8 • Furthermore, the analysis of the overall interpolation process Phase Delay Response 3.7 is not so easy to perform. 3.6 3.5 • The next page shows example FD FIR filters. 3.4 3.3 3.2 3.1 3 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Angular Frequency ω Fig. 4. Phase delay responses. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 13 Tampere University of Technology 14 Lagrange Interpolation Example: Lagrange Interpolation: How to generate y(l)? • For the Lagrange interpolation, the approximating continu- ous-time signal ya(t) is formed as follows: ZERO-ORDER HOLD: = x(n) ya(t) = y(l) y a ( kTin ) = x(k ) ≡ x a (kTin ) for − N / 2 + 1 ≤ k ≤ N / 2. (3) • The resulting approximating ya(t) is a polynomial of t and its µlTin degree is M = N − 1. • The next page illustrates how to use the linear (N = 2) and nl (nl+µl)Tin nl+1 the cubic Lagrange interpolation (N = 4) for generating a LINEAR INTERPOLATION: Two samples in use new sample value at an arbitrary point between the existing samples x(nl ) and x(nl +1). ya(t) = x(n) = y(l) • From the curiosity, also the zero-order hold is included. • It should be pointed out that the approximation error µlTin | xa(t) − ya(t)| becomes the smallest in the interval between the existing samples x(nl ) and x(nl +1) if equally many sample nl (nl+µl)Tin nl+1 values are used before and after the new sample instant. • This is why it is preferred to select N to be an even integer. CUBIC LAGRANGE INTERPOLATION: Four samples in use • After introducing a hybrid analog/digital model to be mim- = x(n) ya(t) = y(l) icked digitally, we consider in more details how to realize the Lagrange interpolation using efficient digital implemen- tations. µlTin nl (nl+µl)Tin nl+1 Fig. 5. Lagrange interpolation. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 15 Tampere University of Technology 16 Hybrid Analog/Digital Model to be Mimicked Hybrid Analog/Digital Model to be Mimicked xs(t) ya(t) y(l) • Assuming that ha(t) is zero outside the interval −NTin /2 ≤ x(n) DAC ha(t) t <NTin /2, y(l) obtained by sampling ya(t) at tl is given by N / 2 −1 Resample at the instant y (l ) = y a (t l ) = ∑ x(nl − k )ha ((k + µ l )Tin ) . (5) tl = ( nl + µ l ) Tin k =− N / 2 • By comparing Equations (2) and (5), it can be seen that the Fig. 6. Analog model for the interpolation filter. impulse responses of the analog and discrete-time filters are related as follows: • Interpolation is a reconstruction problem where the ap- proximating signal ya(t) is reconstructed based on the exist- h(k , µ l ) = ha ((k + µ l )Tin ) (6) ing discrete-time samples x(n). for k = −N/2, −N/2+1,···, N/2−1. • Therefore, a useful way of modeling discrete-time interpola- tion filters is to use the analog system shown in Fig. 6. • In the causal case, ha(t) is delayed by NTin/2, i.e., the impulse response is given by ha(t − NTin/2). • In this system, x(n) is first converted with the aid of the ideal D/A converter to the sequence of weighted and shifted • In this case, y(l) obtained NTin/2 time units later is given by ∞ analog impulses x s (t ) = ∑n = −∞ x(n)δ a (t − nTin ) . N −1 y (l ) = ∑ x(n l + N / 2 − k ) ha ((k + µ l − N / 2 )Tin ) . (7) • xs(t) is then filtered using an analog reconstruction filter k =0 with impulse response ha(t) resulting in the following convo- • In the sequel, the non-causal ha(t) [Equation (5)] and the lution: causal ha(t − NTin/2) [Equation (7)] are used for the design ∞ ∞ and interpolation purposes, respectively. y a (t ) = ∫ xs (τ )ha (t − τ )dτ = ∑ x(k )ha (t − kTin ). (4) −∞ k = −∞ Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 17 Tampere University of Technology 18 Why to Use the Analog Model? Why to Use the Analog Model? • Interpolation is generally considered as a time-domain prob- • The use of the analog model converts the interpolation prob- lem of fitting polynomial through the existing samples, lem from the time-domain to the frequency domain in a which is not very practical approach for DSP applications. manner to be considered next. • These include the Lagrange and B-spline interpolations Synthesis problem for interpolation: Determine ha(t) • This is because the time-domain characteristics of the input such that sequence x(n) are not usually known. What is usually 1. It provides the desired filtering performance. known is the frequency-domain performance of the signal. 2. The overall system of Fig. 6 can be implemented digitally • It should be pointed out that recently Atanas Gotchev, using an efficient structure. Karen Egiazarian, and Tapio Saramäki have improved the performance of B-splines in interpolation problems, espe- cially in the case of images, by using modified B-splines consisting of a weighted sum of odd-order B-splines. Con- tact: saram@vip.fi (home e-mail address of Saramäki). • The main idea is to determine the weights in such a manner that the resulting filter effectively preserves the baseband of interest and attenuates the corresponding images. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 19 Tampere University of Technology 20 Frequency-Domain Criteria Role of ha(t) in the Frequency Domain • For the overall system of Fig. 6, the Fourier transform of As shown below, the role of the reconstruction filter with im- ya(t) is related to that of the sequence x(n) or equivalently to pulse response ha(t) is to attenuate the extra images of ∞ ∞ that of the signal x s (t ) = ∑ n = −∞ x( n)δ a (t − nTin ) through xs (t ) = ∑n = −∞ x(n)δ a (t − nTin ) and to preserve the signal components only in the original baseband [0, Fin / 2]. Ya ( j 2πf ) = H a ( j 2πf ) X (e j 2πf / F ) = in ∞ X a( 2 j π f ) (8) = H a ( j 2πf ) Fin ∑X k =− ∞ a ( j 2π ( f − kFin )) αFin/2 where Fin = 1 / Tin is the sampling rate of the input signal f and Ha( j2π f ) is the Fourier transform of the reconstruction filter with impulse response ha(t). Fin/2 • The last form of Equation (8) is for the case where Fig. 7. The spectrum of the original continuous-time signal bandlim- x( n) = xa (nTin ) are samples of a continuous-time signal ited to f ≤ αFin . The sequence is formed as x( n) = x a (nTin ) . xa (t ) with X a ( j 2πf ) being its Fourier transform. 1/Fin Ha( j2π f) X(e j2π f /Fin) Fin (1−α)Fin f Fin/2 Fin 2Fin ∞ Fig. 8. The spectrum of x s (t ) = ∑ n = −∞ x( n)δ a (t − nTin ) , denoted j2π f /Fin by X(e ) and the frequency response of the reconstruction filter with impulse response ha(t), denoted by Ha( j2π f ). Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 21 Tampere University of Technology 22 Criteria for the Uniform Sampling: Interpolation and Practical Criteria Decimation • Like for conventional digital interpolators and decimators, • If y(l) is generated at the time instants t l = lTout , then the criteria can be stated as ∞ 1 − δ p ≤ Fin H a ( j 2πf ) ≤ 1 + δ p for f ∈ [0, f p ] (11a) Y (e j 2πf / F ) = Fout out ∑ Y ( j 2π ( f − kF k =− ∞ a out )) , (9) Fin H a ( j 2πf ) ≤ δ s for f ∈ Ω s , (11b) where Fout = 1 / Tout is the sampling rate of the output signal y(l) and the baseband of interest is [0, Fout/2]. where f p < FC / 2 and • The case β = Fout / Fin > 1 corresponds to the interpolation. • The case β = Fout / Fin < 1 corresponds to the decimation. [F / 2, ∞ ) for Type A • In both cases, the ideal response for Ha( j2π f ) avoiding both C Ω s = [FC − f p , ∞ ) for Type B (11c) imaging and aliasing is given by ∞ r [kFC − f p , kFC + f p ] for Type C. k =1 1 / Fin for 0 ≤ f ≤ FC /2 D( f ) = (10a) 0 for f > FC /2, • For Type A, no aliasing or imaging are allowed. where • For Type C decimation case, aliasing is allowed into the transition band [fp, Fout/2]. For Type B, aliasing into this FC = min ( Fin , Fout ). (10b) band is allowed only from band [ Fout/2, Fout−fp]. • Note that in the interpolation case, it is enough to attenuate • In the interpolation case, Types B and C are useful if most of the images of X(e j2π f /Fin). the energy of the incoming signal is in the range [0, fp]. • In the decimation case, X(e j2π f /Fin) should be band-limited into the range [0, Fout/2], that is, the region [ Fout/2, Fin/2] should be attenuated in order to avoid aliasing. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 23 Tampere University of Technology 24 Polynomial-Based Interpolation Filters Fig. 9. An example of the impulse response of a polynomial-based interpolation filter for N=6 and • In order to arrive at an efficient digital implementation of M=5 (5th order Lagrange). In each interval the system of Fig. 6, it is required that ha(t) is expressible as a polynomial of t in each interval [kTin , ( k+ 1 )Tin ] for + − [kTin , ( k+ 1 ) Tin ] for k= −N/2,···, N/2−1 ha(t) is a k= −N/2,···, N/2−1. polynomial of degree 5. • This is achieved by expressing ha(t) as follows: 1 M ha ((k + µ l ) Tin ) = ∑ c m (k )µ l m ˆ (12) m =0 0.8 Impulse response ha(t) for k= −N/2,···, N/2−1. Here, the c m (k ) ’s are the coeffi- ˆ cients and M is the degree of the polynomials. 0.6 • When µ l varies from 0 to 1, then ha(t) takes in each interval [kTin , ( k +1 ) Tin ] for k = −N/2,···, N/2−1 the following form: 0.4 m M t − kT ∑ c m (k ) T in 0.2 h a (t ) = ˆ , (13) m =0 in 0 as is desired. • An example impulse response ha(t) is shown in Fig. 9 for −3 −2 −1 0 1 2 3 N = 6 and M = 5. Time in T in Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 25 Tampere University of Technology 26 Generation an Efficient Digital Implementation: Interpretation of the Above Equation Original Farrow Structure • According to Equation (15), the l t h output sample y(l) at the • By using the analog model, the digital implementation struc- time instant t = (nl +µl)Tin can be generated based on the N ture for polynomial-based interpolation filters can be derived existing samples x(n) being located at n = nl −N/2+1, by substituting Equation (12) into Equation (7) giving nl −N/2+2,···, nl +N/2 in the following three steps: N −1 M • First, these samples are filtered using M+1 FIR filters having ∑ x(nl − k + N / 2) ∑ c m (k − N / 2)µ l m y (l ) = ˆ . (14) the transfer functions of the form k =0 m=0 N −1 • Alternatively, this equation can be expressed as Cm ( z) = ˆ ∑ c m ( k − N / 2) z − k ˆ for m = 0,1, , M . (16) k =0 M y (l ) = ∑ v m (nl ) µ l m , (15a) • Second, the outputs of these filters, denoted by vm(nl), are multiplied with constants µ l . m m=0 where • Third, the multiplication results are added, leading to the so- called original Farrow structure shown in Figure 10. N −1 vm (nl ) = ∑ x(nl − k + N / 2)cm (k − N / 2) ˆ • In this figure, the input to both structures is denoted by k =0 (15b) x(nl + N/2) to emphasize the fact that this is the last existing N −1 sample value when evaluating y(l). = ∑ cm (k − N / 2)x(nl + N / 2 − k ) . ˆ k =0 • If desired, the transfer functions C m ( z ) can share the delay ˆ elements. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 27 Tampere University of Technology 28 Original Farrow structure for a polynomial-based Characteristics of the Farrow Structure interpolation filter. • The number of FIR filters is M+1. x(nl + N/2) / • The length of these filters is N. • Filter coefficients are directly determined by the polynomial ^ CM(z) ^ C2(z) ^ C1(z) ^ C0(z) coefficients of the impulse response ha(t). • The main advantage of the Farrow structure is that all the fil- vM(nl) v2(nl) v1(nl) v0(nl) ter coefficients are fixed. y(l) µl • The only changeable parameters are the fractional interval µ l as well as nl that depend on the l th output sampling instant. (a) x(nl + N/2) • The control of µ l is easier during the computation than in the ^ − cM(−N/2) ^− c1(−N/2) c − ^0(−N/2) implementation based on the FD filters. • The resolution of µl is limited only by the precision of the −1 −1 −1 Z ^M(−N/2 + 1) c − Z ^1(−N/2 + 1) c − Z ^0(−N/2 + 1) c − arithmetic not by the size of the memory or lookup table. • These characteristics of the Farrow structure make it a very −1 −1 −1 Z Z Z attractive structure to be implemented using a VLSI circuit or a signal processor. ^M(N/2 − 1) c ^1(N/2 − 1) c c0(N/2 − 1) ^ vM(nl) v1(nl) v0(nl) µl y(l) (b) Fig. 10. Farrow structure for a polynomial-based interpolation filter. (a) Basic structure. (b) Details. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 29 Tampere University of Technology 30 What Do We Have Up to Now? Design Strategies 1. The fundamental idea is that the ya(t), the output of our ana- • We consider only polynomial-based interpolation filters be- log model of Figure 6, approximates according to some cri- cause they can be efficiently implemented using the Farrow terion the following ‘ideal’ continuous- time function structure. xa (t ) = ∑ n = −∞ x( n) sin[π (t − nTin ) / Tin ] [π (t − nTin ) / Tin ] . ∞ • The following two design methods for polynomial-based in- 2. By using the analog model, this approximation problem can terpolation filters are considered: be converted to the problem of designing the reconstruction filter ha(t). 1. Conventional time-domain synthesis: Lagrange interpola- tion. 3. If ha(t) is a piecewise polynomial, then the Farrow structure can be used for implementation. 2. Frequency-domain synthesis based on the analog model of Figure 6: Minimax and least-mean-square optimization of the interpolation filters. • The remaining task is to optimize the coefficients c m (k ) for ˆ • These design methods will be compared with the aid of the the polynomial-based impulse response ha(t). See Equations resulting impulse and frequency responses of the corre- (12) and (13). sponding reconstruction filter ha(t). Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 31 Tampere University of Technology 32 Filter responses corresponding to the linear and Lagrange Interpolation cubic Lagrange interpolation • Lagrange interpolation was originally used in mathematics, 1.2 not in signal processing (discovered by Joseph-Louis La- 1 grange 1736-1813). Impulse response h (t) 0.8 ⇒ Does not offer a good filtering characteristics. a 0.6 • Filter design is done in the time domain and the filter coef- 0.4 ficients can be given in a closed form. 0.2 ⇒ No need for optimization. 0 • The synthesis can be accomplished as follows: −0.2 −2 −1 0 1 2 Time in T 1. Choose M, the degree of the interpolation. The length of in the filter is then N = M+1 Fig. 11. Impulse response ha(t) for the cubic (solid line) and linear 2. Based on the time-domain conditions for ya(t) as given by (dashed line) Lagrange interpolation filters. Equation (3), the polynomial coefficients for ha(t) are de- termined from the following equation: 0 M j−x N /2 ∑ c m (k ) x = ∏ − k + j m −10 ˆ (17) Magnitude in dB m =0 j = − N / 2 +1 −20 j≠k −30 for k = −N/2, −N/2+1,···, N/2−1. −40 • As mentioned earlier, it is desired to use an even value for −50 N. This means that M should be an odd integer. −60 0 0.5 1 1.5 2 2.5 3 3.5 4 • The following page show both the time- and frequency- Frequency in Fin domain responses corresponding to the linear (M = 1, Fig. 12. Magnitude responses for the cubic (solid line) and linear N = 2) and cubic (M = 3, N = 4) Langrange interpolation. (dashed line) Lagrange interpolation filters. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 33 Tampere University of Technology 34 Disadvantages of the Lagrange interpolation filters 0 • The frequency response cannot be changed: stopband at- −10 tenuation, passband ripple, and edge frequencies are fixed. −20 • If the degree of the interpolation is M, then the length of the M=1 Magnitude in dB filter is always N = M +1. −30 • The stopband attenuation cannot be improved significantly −40 M=9 by increasing the degree of the interpolation (e.g., 5 or 7), −50 but the number of filter coefficients increases fast. −60 • The characteristics of the resulting reconstruction filter are −70 very poor, at least if the input signal contains frequency components close to half the sampling frequency. −80 • The attenuation of the unwanted images is very good −90 0 0.5 1 1.5 2 2.5 3 3.5 4 around the multiples of Fin. Frequency in F in • Therefore, these filters can be used only if the sampling rate Fig. 13. The magnitude responses |Ha( j 2 π f )| of the Lagrange inter- is increased before using them. polation filters with degree M = 1, 3, 5, 7, and 9 (solid lines) and the • For instance, the third order Lagrange interpolation filter spectrum Xs( j2π f ) of the input signal with the bandwidth of 0.1Fin attenuates the image frequencies at least by 60dB if the (dashed-line). bandwidth of the input signal is 0.1Fin as shown in Fig. 13. • This implies the use of a discrete-time conventional interpo- lation filter that increases the sampling rate by a factor of five. • This makes the overall synthesis more complicated com- pared to the case where a single Farrow structure is used. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 35 Tampere University of Technology 36 General Frequency-Domain Synthesis Technique for Construction of the impulse response ha(t) for the Interpolation Filters modified Farrow structure. • It is desired to design in the analog model of Figure 6 the • In order to arrive at an efficient digital implementation, ha(t) reconstruction filter with impulse response ha(t) only sub- is assumed to satisfy the following conditions: ject to the restriction that the overall system can be imple- 1) ha(t) is nonzero for −NTin/2 ≤ t < NTin/2. mented using the original or modified Farrow structure. 2) The length of the filter N is an even integer. • In this case, M and N can be selected arbitrarily. 3) ha(t) is a piecewise-polynomial of degree M in each inter- • The main goal to is to optimize the coefficients of the Far- val nTin ≤ t < (n +1)Tinfor n = −N/2, −N/2+1,···, N/2−1. row structure and the corresponding impulse response ha(t) in such a manner that the overall system provides the 4) ha(t) is symmetrical, that is ha(−t) = ha(t) except for the desired frequency-domain behavior that depends on the ap- time instants t = nTin for n = −N/2, −N/2+1, ···, N/2. plication. • Conditions 1), 2), and 3) guarantee the corresponding • The purpose is to determine the coefficients in such a man- causal system with impulse response ha(t − NTin /2) can be ner that the amplitude response of the reconstruction filter implemented using an efficient digital implementation. approximates an arbitrary amplitude response either in • The role of Conditions 4) is twofold. the least-mean square sense or in the minimax sense. 1) For the causal system, the phase response is linear. • It is also desired to allow the use of some time-domain 2) In the modified Farrow to be described later, the fixed conditions for the overall impulse response ha(t). FIR filters have either a symmetrical or anti- • Before developing the general synthesis scheme, the modi- symmetrical impulse responses. This enables us to uti- fied Farrow structure is introduced. lize the coefficient symmetry, reducing the number of multipliers in the implementation. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 37 Tampere University of Technology 38 Time-domain conditions How to construct the impulse response ha(t) • There exist the following time-domain constraints of inter- • Instead of using µ l we use 2µ l − 1 as a basis polynomial to est: construct ha(t) as follows: 3) Case A: There are no time-domain conditions. M ha ((n + µ l )Tin ) = ∑ c m (n)(2µ l − 1) m (18) 4) Case B: ha(t) is continuous at t=kTin for k = ±1, ±2,···, m =0 ±N/2−1. f o r n = −N/2, −N/2+1,···, N/2−1(see Fig. 14). Here, the 5) Case C: ha(0) = 1 and ha(kTin) = 0 for k = ±1, ±2,···, cm(n)’s are the unknown polynomial coefficients. ±N/2. 6) Case D: The first derivative of ha(t) is continuous at 1 1 1 1 Amplitude t=kTin for k = 0 and for k = ±1, ±2,···, ±(N/2−1). • Case C guarantees that if the new sampling instant occurs at −1 −1 −1 −1 the instant of the existing sample, then the sample value is 0 1 0 1 0 1 0 1 preserved. Time in T in • Case D is of importance when determining the derivative of a continuous-time signal with the aid of its discrete-time Fig. 14. Polynomials (2µ l − 1)m for m = 0, 1, 2, and 3. samples and a generalized modified Farrow structure. • When µ l varies from 0 to 1, then ha(t) takes in each interval [kTin , ( k +1 ) Tin ] for k = −N/2,···, N/2−1 the following form: 2(t − kTin ) − 1 M m h a (t ) = ∑ c m (k ) , (19) m=0 Tin • The symmetry property ha(−t) = ha(t) is achieved by c m (n) = (−1) m c m (−n − 1) (20) Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 39 Tampere University of Technology 40 for m = 0, 1,···, M and n = 0, 1,···, N/2−1. This condition Example on how to construct ha(t) for N=8 and halves the number of unknowns. M=3. • ha(t) can be now constructed as follows N / 2 −1 M ∑ ∑ c m (n) g (n, m, t ) 0.6 ha (t ) = (21) (a) n =0 m = 0 0 0.5 where cm(n)’s are unknown coefficients and g(n,m,t)’s are (b) 0 the basis functions given by −0.5 2(t − nTin ) m 0.1 (c) − 1 for nTin < t ≤ (n + 1)Tin 0 Tin −0.1 (22) 0.1 m (d) g (n, m, t ) = (−1) 2(t + (n + 1)Tin ) − 1 for − (n + 1)Tin ≤ t < −nTin m 0 −0.1 1 Tin (e) 0 otherwise 0.5 0 1 Amplitude −4 −3 −2 −1 0 1 2 3 4 Time in T 0 in Fig. 16. Construction of the overall impulse response ha(t) for N = 8 and −1 N / 2 −1 ∑c −2 −1 0 1 2 Time in Tin M = 3. The weighted basis functions m (n) g (n, m, t ) for m = 0 (a), n=0 m = 1 (b), m = 2 (c), and m = 3 (d). (e) The resulting impulse response ha(t) Fig. 15. The basis function g(n, m, t) for n = 1 and m = 3. obtained as a sum of these responses. Figure 15 shows an example basis function, whereas Fig. 16 shows how the overall impulse response can be constructed using weighted basis functions. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 41 Tampere University of Technology 42 Modified Farrow structure Modified Farrow Structure • In the above approach, the polynomial coefficients cm(n) x(nl + N/2) / are symmetrical according to Equation (20). • 2µ l −1 is used, instead of µ l . CM(z) C2(z) C1(z) C0(z) • Therefore, interpolation filters based on this approach can be implemented using the modified Farrow structure vM(nl) v2(nl) v1(nl) v0(nl) shown in Fig. 17. y(l) • This modified structure has two differences compared to 2µl-1 the original structure. (a) 1. The output samples vm(nl) of the FIR filters are multiplied x(nl + N/2) with 2µ l − 1 instead of µ l. − cM(−N/2) − c1(−N/2) − c0(−N/2) 2. FIR filters with the transfer functions −1 −1 −1 Z Z Z cM(−N/2 + 1) − c1(−N/2 + 1) − c0(−N/2 + 1) − N −1 Cm ( z) = ∑ c m (k − N / 2) z −k for m = 0, 1, , M (23) −1 −1 −1 k =0 Z Z Z possess the symmetry properties given by Equation (20). • When m is zero or even, c m ( N / 2 − 1 + k ) = c m (− N / 2 − k ) cM(N/2 − 1) c1(N/2 − 1) c0(N/2 − 1) for k = 0, 1,···, N/2−1. vM(nl) v1(nl) v0(nl) • For m odd, cm ( N / 2 − 1 + k ) = −cm (− N / 2 − k ) for k = 0, 2µl − 1 y(l) 1,···, N/2−1. (b) • When exploiting these symmetries, the number of coeffi- Fig. 17. Modified Farrow structure. (a) Basic structure. (b) Details. cients to be implemented can be reduced from (M+1)N to (M+1)N/2. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 43 Tampere University of Technology 44 The frequency response for the reconstruction filter Optimization Problems with impulse response ha(t) • The very attractive property of the above Ha( j2π f ) is that it • After some manipulations, the frequency response can be is linear with with respect its unknowns cm(n). expressed as • These unknowns can be easily found to minimize N / 2 −1 M H a ( j 2πf ) = ∑ ∑ c m (n)G(n, m, f ) , (24) δ ∞ = max W ( f )( H a ( j 2πf ) − D( f ) ) (26) n =0 m =0 f ∈X where G(n, m, f ) is the frequency responses of the basis or function g(n, m, t). δ 2 = ∫X [W ( f )( H a ( j 2π ) − D ( f ) )] df 2 (27) • Since g(n, m, t) is symmetrical around t = 0, G(n, m, f ) is real and is given by subject to the given time-domain conditions of ha(t) stated on 1 sin(πfTin ) Page 37. 2 cos( 2πfTin (n + 2 )) (−1) m!Φ (m, f ) + m/2 πfTin • The first and second criteria correspond to the optimization G ( n , m, f ) = for m even in the minimax and the least-mean-square sense, respec- ( m +1) / 2 2(−1) m!sin( 2πfTin (n + 1 ))Φ ( m, f ) for m odd, 2 tively. • Here X ⊂ [0,∞) is a compact subset and D( f ) is an arbitrary (25a) desired function (continuous) and W( f ) is an arbitrary weighting function (positive). where • Here, the approximation region X consists of a set of pass- ( m −1) / 2 ( −1)k sin(π f Tin ) cos(π f Tin ) band and stopband regions. ∑ 2k − m Φ ( m, f ) = (π f Tin) − . k =0 ( 2k )! π f Tin ( 2k + 1) • The actual optimization can be accomplished in a manner similar to the design of linear-phase FIR filters. (25b) Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 45 Tampere University of Technology 46 • Optimization algorithms have been implemented in Matlab. Parameters for the optimization For minimax problem, linear programming can be used to optimize the filter coefficients. • Design parameters for the optimization programs are the • For both problems, the time-domain conditions are included following: in the problem in such a manner that they become uncon- 1. Edge frequencies for passband(s) and stopband(s). strained problems. 2. Desired amplitude and weight for every band. • This makes the overall optimization algorithms very fast. 3. N, the length of the filter. 4. M, the degree of the interpolation. 5. The number of grid points. 6. Time-domain constraints: 1) Case A: There are no time-domain conditions. 2) Case B: ha(t) is continuous at t=kTin for k= ±1, ±2,···, ±N/2−1. 3) Case C: ha(0) = 1 and ha(kTin) =0 for k =±1, ±2,···, ±N/2. 4) Case D: The first derivative of ha(t) is continuous at t=kTin for k=0 and for k=±1, ±2,···, ±(N/2−1). • Before introducing the applications, two Case A design ex- amples are considered. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 47 Tampere University of Technology 48 Optimized Case A minimax design Optimized Case A least-squared design • M = 7,N = 2 4 , X = [0, 0.4 Fin ] t [0.6 Fin , ∞ ), D( f ) is unity • M = 5, N = 8 , X = [0, 0.35 Fin ] t [0.75 Fin , ∞ ), D ( f ) is unity and zero on the first and second bands, whereas W ( f ) is and zero on the first and second bands, whereas W ( f ) is 0.002 and 1, respectively. 0.02 and 1, respectively. 0 0 −10 −20 −20 Magnitude in dB Magnitude in dB −40 −30 −40 −60 −50 −60 −80 −70 −100 −80 −90 0 0.5 1 1.5 2 2.5 3 3.5 4 0 0.5 1 1.5 2 2.5 3 3.5 4 Frequency in F Frequency in Fin in 1 1 0.8 0.8 Impulse response ha(t) Impulse response h (t) 0.6 a 0.6 0.4 0.4 0.2 0.2 0 0 −0.2 −0.2 −4 −3 −2 −1 0 1 2 3 4 −12 −10 −8 −6 −4 −2 0 2 4 6 8 10 12 Time in T Time in T in in Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 49 Tampere University of Technology 50 Application A: Up-Sampling Between Arbitrary Direct Design Sampling Frequencies • The Farrow structure can be directly used for providing the increse between an arbitrary input sampling rate Fin and an arbitrary output sampling rate Fout. 0 Passband Amplitude Response 1.001 • It is desired that Ha( j 2 πf ) approximates unity for 1.0005 0 ≤ f ≤ 0.45Fin with tolerance of 0.001 and zero for 1 f ≥ 0.5Fin with tolerance of 0.00001 (100-dB attenuation). Amplitude in dB −50 0.9995 • When using the minimax optimization, the given criteria 0.999 are met by N = 92 and M = 6, as shown on Page 50. This 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 implementation requires 7 fixed branch filters of length 92. −100 • The implementation can be simplified using fixed linear- phase FIR interpolators before the Farrow structure, as pro- posed by Saramäki and Ritoniemi. • N = 4 and M = 3 are required if the sampling rate is in- −150 0 2 4 6 8 10 12 14 16 18 20 creased by a factor of six by using a two-stage fixed inter- Frequency f / F in polator with interpolation factors of two and three and FIR filters of order 183 and 11, respectively. See Page 51. • The block diagram for this multistage implementation is shown below. x(n) Modified y(l) 2 2F1(z ) 3 3F2(z ) Farrow Fin 2Fin 6Fin Structure Fout • Note that the same structure can be used for any Fout > Fin. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 51 Tampere University of Technology 52 Design with fixed interpolators before the Farrow Application B: Down-Sampling Between Arbitrary structure: Simultaneous optimization has been used. Sampling Frequencies: First Alternative • There exist two alternatives to perform down-sampling. Solid: 1st interpolator, Dashed: 2nd interpolator, Dot−dashed: Farrow 0 • The first alternative is to increase the sampling rate to a multiple of the output sampling rate and then to decimate to the desired output sampling rate. Amplitude in dB • As an example, consider sampling rate reduction from −50 48 kHz to 44.1 kHz using the a structure shown below −100 x(n) Modified y(l) Farrow F1(z ) 2 F1(z ) 2 Fin Structure 4Fout 2Fout Fout −150 0 2 4 6 8 10 12 14 16 18 20 Frequency as a Fraction of Fin • The passband edge is 20 kHz and aliasing is allowed into 0 the band between 20 kHz and 44.1/2 kHz. Passband Amplitude Response 1.001 • The passband ripple is 0.0001 and the minimum stopband 1.0005 attenuation is 120 dB. 1 Amplitude in dB −50 0.9995 • To meet these criteria L = 4 and N = 6 are required by the 0.999 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 modified Farrow structure, whereas the orders of the first and second decimator are 4 and 126, respectively. • The resulting responses are shown on the next page. −100 −150 0 2 4 6 8 10 12 14 16 18 20 Frequency f / F in Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 53 Tampere University of Technology 54 Three -stage Decimator using the Modified Farrow Application B: Down-Sampling Between Arbitrary structure: Simultaneous optimization has been used. Sampling Frequencies: Second Alternative 0 Passband Amplitude Response • The second alternative is to use the transposed modified 1.0001 Farrow structure together with fixed decimators. −30 • Due to the lack of time, this alternative being more efficient 1 1 −60 is not included here. This is because the transposed struc- Amplitude in dB 1 ture is a very new idea. 0.9999 −90 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 • An up-to-date-version of this pile including the above will −120 be found in two months in http://www.cs.tut.fi/~ts/ called POLYINT. • There is an article: Djordje Babic, Jussi Vesma, Tapio −150 −180 Saramäki, and Markku Renfors, “Implementation of the 0 2 4 6 8 10 12 14 16 18 20 Frequency as a Fraction of Fout transposed Farrrow structure”, submitted to ISCAS 2002. Solid: Farrow, Dashed: 1st decimator, Dot−dashed: 2nd decimator • The main advantage of this structure is that the same struc- ture can be used for any input sampling rate Fin > Fout. 0 −30 Amplitude in dB −60 −90 −120 −150 −180 0 2 4 6 8 10 12 14 16 18 20 Frequency as a Fraction of F out Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 55 Tampere University of Technology 56 Application C: Continuous-Time Signal Processing Generalized Farrow Structure for Determining the Derivative of ya(t) • The Farrow structure can be easily generalized to process digitally the reconstructed signal ya(t). • In the intervals nTi n ≤t <(n+1)Ti n for n = 0, 1, 2,·ya(t) can • These applications include, among others, determining the be expressed as derivative or the integral of ya(t). M • The derivative is widely utilized, for example, in finding y a (t ) t =( n + µ )Tin = p (n, µ ) = ∑ vm (n)[2µ − 1] m (28) m =0 the location of the maximum or minimum of the signal. • The integral can be used to calculate the energy of the sig- where the vm(n)’s are the output samples of the FIR branch nal over the given time interval. filters in the modified Farrow structure. • We concentrate on determining the derivative of ya(t). • The derivative of ya(t) in the intervals is thus given by d y a (t ) d p ( n, µ ) M t = ( n + µ )Tin = = ∑ v m ( n) 2m[2 µ − 1] m −1 . (29) dt dµ m =0 • The derivative of ya(t) at t = (n+µ)Ti n can be determined by multiplying the vm(n)’s by 2m(2µ −1) m−1, instead of (2µ −1) m in the modified Farrow structure. • For estimating the derivative, it is desired that Ha(j2πf) ap- proximates 1/(2πf) in the passband with the weighting equal to 2πf. • In the stopband, it is desired to approximate zero with a constant weight. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 57 Tampere University of Technology 58 Example on the Derivative Approximation Continuous-time processing of an ECG signal • It is desired to estimate the continuous-time derivative of a 1 discrete-time ECG signal shown in Fig. 18(b). Amplitude (a) • Figures 18(b) and 18(c) show the continuous-time interpo- 0.5 lated signal and the derivative signal, respectively. 0 • For ha(t) used for determining the derivative signal, N =8, 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Amplitude M =5, and the passband and stopband edges are located at (b) 0.35Fs and 0.65Fs, respectively. 0.5 • When using the minimax optimization criterion with weig- 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 ting equal to 0.035 in the passband and equal to unity in the 1 stopband, we end up with ha(t) with the amplitude and im- (c) Amplitude 0.5 pulse responses shown in Fig. 19. 0 −0.5 −1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Time t / s Fig.18. Continuous-time processing of an ECG signal. (a) Discrete- time ECG signal. (b) Interpolated continuous-time signal. (c) Con- tinuous-time derivative. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 59 Tampere University of Technology 60 Characteristics of the differentiator Application C: Symbol Synchronization in Digital Impulse response for dh(t)/dt 1.5 Receivers 1 r(t) r(n) x(n) Interpolator y(l) â(l) Matched filter ADC hR(n) h(k,µl) Decision 0.5 Amplitude 0 nl µl −0.5 ∼ Fin=1/ Tin Timing estimation −1 −1.5 −4 −3 −2 −1 0 1 2 3 4 Fig. 20. Digital receiver with non-synchronized sampling. Time / T • The sampling of the received signal is performed by a fixed (a) sampling clock, and thus, sampling is not synchronized to 2.5 the received symbols. 2 ⇒ timing adjustment must be done by digital methods after sampling. Amplitude 1.5 • Can be done by using interpolation filter. 1 • Advantages of nonsynchronized sampling: 0.5 − separates the analog and digital parts − easy to change the sampling rate 0 0 0.5 1 1.5 2 2.5 3 − sampling frequency does not have to be a multiple of Frequency f / F in the symbol frequency (only high enough to avoid (b) aliasing) Fig.19. Optimized differentiator. (a) Impulse response. (b) Ampli- − no need for complex PLL circuit. tude response Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 61 Tampere University of Technology 62 Practical Case • Two polynomial-based interpolators have been designed in the minimax sense to illustrate the use of the above- • When deriving the frequency-domain specifications for the mentioned specifications. anti-imaging filter ha(t), it is assumed that • It is assumed that the received signal has a raised cosine 1) The pulse shape of the transmitted signal has the excess pulse shape with the excess bandwidth of α =0.15 and the bandwidth of α and the ratio between the sampling rate oversampling ratio is R =1.75. Fin and the symbol rate is R. • The passband edge for both filters is fp = βFin = 23/70Fin 2) In order to avoid aliasing, it is required that R ≥ (1+α ). (≈0.33Fin) • Based on these assumptions, the input signal of the interpo- • Furthermore, it is required that the passband distortion is lator x(n) contains the desired component in the frequency less than δp = 0.01 and the minimum stopband attenuation is range [0, βFin ], where β =(1+α )/R/2 and undesired images As = 50dB. in the bands [(k−β)Fin, (k+β)Fin] for k =1, 2,···. • The first filter has a uniform stopband .In order to meet the • Therefore, the desired function for Ha(j2π f ) is specified by specifications, N= 8 and M =3 are required, as shown in Figure 21(a) on the next page 1 / Fin for 0 ≤ f ≤ β Fin D( f ) = (30) • The second filter has a non-uniform stopband as given by 0 for f ∈ Ω s , Eq. (31b) and the spectrum of the raised cosine pulse shape where the stopband region, denoted by Ωs, can be selected as is used as a weighting function. In this case N= 6 and M = 3 meets the requirements giving As =50.0dB and δp = 0.01, as Ω s = [(1 − β ) Fin , ∞ ) (31a) shown in Figure 21(b) on the next page. or ∞ Ω s = t [(k − β ) Fin , (k + β ) Fin ] (31b) k =1 Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 63 Tampere University of Technology 64 Fig.21. The magnitude response of the interpolation filter (solid Properties of Minimax Case A designs line), the spectrum for the raised cosine pulse (dashed line) and for the reconstructed signal ya(t) (dark area). (a)With uniform • Case A: The minimum even value of N can be estimated by stopband. (b) With non-uniform stopbands having the raised cosine weighting. − 20 log10 ( δ pδ s ) − 8.4 N = 2 (32) 7.6( f s − f p ) / Fin 0 where δp and δs are the maximum deviations of the amplitude −10 response from unity for f ∈ [0, fp] and the maximum deviation −20 from zero for f ∈ [ fs, ∞). Magnitude in dB −30 −40 • Here, x stands for the smallest integer that is larger or equal −50 to x. −60 • It has been observed that in most cases the above estimation −70 formula is rather accurate with only a 2 % error. −80 0 0.5 1 1.5 2 2.5 Frequency in Fin 3 3.5 4 • The next problem is to find the minimum value of M to meet the criteria. (a) • To illustrate this the following specifications are considered: 0 Specifications 1: The passband and stopband edges are at −10 fp =0.25Fin and at fs =0.75Fin. −20 Specifications 2: The passband and stopband edges are at Magnitude in dB −30 fp =0.25Fin and at fs =0.5Fin. −40 Specifications 3: The passband and stopband edges are at −50 fp =0.375Fin and at fs =0.675Fin. −60 Specifications 4: The passband and stopband edges are at −70 fp =0.375Fin and at fs =0.5Fin. −80 0 0.5 1 1.5 2 2.5 3 3.5 4 Frequency in Fin (b) Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 65 Tampere University of Technology 66 Properties of Minimax Case A designs Specifications 1 and 2: • In each case, several filters have been designed in the mini- Specifications 1: N= 2, 4,···, 14 Stopband Attenuation in dB Stopband Attenuation in dB 120 max sense with the passband weighting equal to unity and 100 120 100 stopband weightings of Ws =1, Ws =10, Ws =100, and 80 w =1 s 80 ws=10 Ws =1000. 60 40 60 • M, the degree of the polynomial in each subinterval, varies 20 0 40 20 0 2 4 6 8 10 12 0 2 4 6 8 10 12 from 0 to 12. N, the number of intervals varies from 2 to the Degree M Degree M Stopband Attenuation in dB Stopband Attenuation in dB smallest integer for which the stopband ripple for the ampli- 140 tude response is less than or equal to 0.0001 (100 dB) for 120 ws=100 120 w =1000 s 100 Ws =1. 80 100 • For Specifications 1, 2, 3, and 4,, the corresponding small- 60 80 est values of N are 12, 24, 24, and 48, respectively. Recall 40 0 2 4 Degree M 6 8 10 12 60 0 2 4 Degree M 6 8 10 12 that N is an even integer. Specifications 2: N= 2, 4,···, 24 • The following two pages give the results for Case A. Stopband Attenuation in dB Stopband Attenuation in dB 100 100 w =1 ws=10 • For other cases, N is either the same or should be incrased 80 s 80 60 by two. 40 60 • For Case C the passband and stopband edges satisfy 20 40 0 20 0 2 4 6 8 10 12 0 2 4 6 8 10 12 • f p = (1 − ρ ) Fin / 2, f s = (1 + ρ ) Fin / 2 . Degree M Degree M Stopband Attenuation in dB Stopband Attenuation in dB 120 130 120 ws=1000 ws=100 100 110 100 80 90 80 60 70 40 60 0 2 4 6 8 10 12 0 2 4 6 8 10 12 Degree M Degree M Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 67 Tampere University of Technology 68 Specifications 3 and 4: Case A: fs = 0.625Fin , fp = 0.375Fin , δp = 0.01, Specifications 3: N= 2, 4,···, 24 δs = 0.001: N=12; Minimax (solid line): M= 4; Least- squared (dashed line): M= 5 Stopband Attenuation in dB Stopband Attenuation in dB 100 100 80 ws=1 ws=10 80 60 60 Passband amplitude response 40 Amplitude in dB 40 1.01 20 0 0.99 0 20 −20 0.97 0 2 4 6 8 10 12 0 2 4 6 8 10 12 Degree M Degree M −40 0 0.1 0.2 0.3 0.375 −60 Stopband Attenuation in dB Stopband Attenuation in dB 120 −80 w =100 s 120 w =1000 −100 s 100 −120 100 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 80 Frequency as a fraction of Fin 60 80 1 40 60 Impulse response 0 2 4 6 8 10 12 0 2 4 6 8 10 12 0.8 Degree M Degree M 0.6 Specifications 2: N= 2, 4,···, 48 0.4 100 Stopband Attenuation Stopband Attenuation 100 80 ws=1 ws=10 0.2 80 60 0 60 40 −0.2 −6 −5 −4 −3 −2 −1 0 1 2 3 4 5 6 20 40 Time as a fraction of T in 0 20 0 2 4 6 8 10 12 0 2 4 6 8 10 12 Degree M Degree M 120 Stopband Attenuation Stopband Attenuation ws=100 120 w =1000 s 100 100 80 60 80 40 60 0 2 4 6 8 10 12 0 2 4 6 8 10 12 Degree M Degree M Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 69 Tampere University of Technology 70 Case B: fs = 0.625Fin , fp = 0.375Fin , δp = 0.01, Case C: fs = 0.625Fin , fp = 0.375Fin , δp = 0.01, δs = 0.001: N=12; Minimax (solid line): M= 4; Least- δs = 0.001: N=14; Minimax (solid line): M= 5; Least- squared (dashed line): M= 5 squared (dashed line): M= 5 Passband amplitude response Passband amplitude response 1.01 Amplitude in dB Amplitude in dB 1.01 0 0 0.99 1 −20 0.97 −20 −40 −40 0.99 0 0.1 0.2 0.3 0.375 0 0.1 0.2 0.3 0.375 −60 −60 −80 −80 −100 −100 −120 −120 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 Frequency as a fraction of Fin Frequency as a fraction of Fin 1 1 Impulse response 0.8 Impulse response 0.8 0.6 0.6 0.4 0.4 0.2 0.2 0 0 −0.2 −0.2 −6 −5 −4 −3 −2 −1 0 1 2 3 4 5 6 −7 −6 −5 −4 −3 −2 −1 0 1 2 3 4 5 6 7 Time as a fraction of T Time as a fraction of T in in Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 71 Tampere University of Technology 72 Case D: fs = 0.625Fin , fp = 0.375Fin , δp = 0.01, Conclusion δs = 0.001: N=12; Minimax (solid line): M= 5; Least- • An efficient approach has been described for interpolating squared (dashed line): M= 5 new sample values between the existing discrete-time sam- ples. Passband amplitude response • This approach has the following advantages: Amplitude in dB 1.01 0 0.99 • Design directly in the frequency domain is straight- −20 0.97 −40 forward. • The overall system has an efficient implementation 0 0.1 0.2 0.3 0.375 −60 −80 −100 form. −120 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 • The analysis of the system performance is easy. Frequency as a fraction of Fin • There exist several DSP applications. 1 Impulse response 0.8 0.6 0.4 0.2 0 −0.2 −6 −5 −4 −3 −2 −1 0 1 2 3 4 5 6 Time as a fraction of Tin Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 73 Tampere University of Technology 74 References [9] H. Ridha, J. Vesma, T. Saramäki, and M. Renfors, [1] T. I. Laakso, V. Välimäki, M. Karjalainen, and U. K. “Derivative approximations for sampled signals based on Laine, “Splitting the unit delay,” IEEE Signal Processing polynomial interpolation,” in Proc. 13th Int. Conf. on Digi- Magazine, vol. 13, pp. 30-60, Jan. 1996. tal Signal Processing, Santorini, Greece, July 1997, pp. [2] C. W. Farrow, “A continuously variable digital delay 939-942. element,” in Proc. IEEE Int. Symp. Circuits & Syst., [10] H. Ridha, J. Vesma, M. Renfors, and T. Saramäki, Espoo, Finland, June 1988, pp. 2641-2645. “Discrete-time simulation of continuous-time systems us- [3] F. M. Gardner, “Interpolation in digital modems - ing generalized interpolation techniques,” in Proc. 1997 Part I: Fundamentals,” IEEE Trans. Commun., vol. 41, pp. Summer Computer Simulation Conference, Arlington, Vir- 501-507, Mar. 1993. ginia, USA, July 1997, pp. 914-919. [4] L. Erup, F. M. Gardner, and R. A. Harris, “Interpola- [11] V. Tuukkanen, J. Vesma, and M. Renfors, “Com- tion in digital modems - Part II: Implementation and per- bined interpolation and maximum likelihood symbol tim- formance,” IEEE Trans. Commun., vol. 41, pp. 998-1008, ing recovery in digital receivers,” to be presented in 1997 June 1993. IEEE Int. Conference on Universal Personal Communica- tions, San Diego, CA, USA, Oct. 1997. [5] D. Kincaid and W. Cheney, Numerical Analysis. Pa- cific Grove, 1991. [12] T. Saramäki and M. Ritoniemi, "An efficient ap- proach for conversion between arbitrary sampling frequen- [6] J. Vesma and T. Saramäki, “Interpolation filters with cies," in Proc. IEEE Int. Symp. Circuits & Syst., Atlanta, arbitrary frequency response for all-digital receivers,” in GA, May 1996, pp. 285-288. Proc. IEEE Int. Symp. Circuits & Syst., Atlanta, GA, May 1996, pp. 568-571. [13] J. Vesma, R. Hamila, T. Saramäki, and M. Renfors, [7] J. Vesma, M. Renfors, and J. Rinne, “Comparison of “Design of polynomial interpolation filters based on Taylor efficient interpolation techniques for symbol timing recov- series,” in Proc. IX European Signal Processing Conf., ery,” in Proc. IEEE Globecom 96, London, UK, Nov. Rhodes, Greece, Sep. 1998, pp. 283-286. 1996, pp. 953-957. [14] J. Vesma, R. Hamila, M. Renfors, and T. Saramäki, [8] J. Vesma and T. Saramäki, “Optimization and effi- “Continuous-time signal processing based on polynomial cient implementation of FIR filters with adjustable frac- approximation,” in Proc. IEEE Int. Symp. on Circuits and tional delay,” in Proc. IEEE Int. Symp. Circuits & Syst., Systems, Monterey, CA, USA, May 1998, vol. 5, pp. 61- Hong Kong, June 1997, pp. 2256-2259. 65. Jussi Vesma and Tapio Saramäki Jussi Vesma and Tapio Saramäki Tampere University of Technology 75 Tampere University of Technology 76 [15] D. Fu and A. N. Willson, Jr., “Interpolation in timing [22] M. Unser, A. Aldroubi, and M. Eden, “Polynomial recovery using a trigonometric polynomial and its imple- spline signal approximations: Filter design and asymptotic mentation,” in IEEE Globecom 1998 Communications equivalence with Shannon’s sampling theorem,” IEEE Mini Conference Record, Sydney, Australia, Nov. 1998, Trans. Information Theory, vol. 38, pp. 95−103, Jan. 1992. pp. 173−178. [23] J. Vesma, Timing Adjustment in Digital Receivers [16] f. harris, “Performance and design considerations of Using Interpolation. M.Sc. Thesis, Tampere, Finland: Farrow filter used for arbitrary resampling,” in Proc. 13th Tampere University of Tech., Department of Information Int. Conf. on Digital Signal Processing, Santorini, Greece, Technology, Nov. 1995. July 1997, pp. 595−599. [24] V. Välimäki, Discrete-Time Modeling of Acoustic [17] G. Oetken, “A new approach for the design of digital Tubes Using Fractional Delay Filters. Doctoral thesis, interpolation filters,” IEEE Trans. Acoust., Speech, Signal Espoo, Finland: Helsinki University of Technology, Dec. 1995. Process., vol. ASSP−27, pp. 637−643, Dec. 1979. [25] S. R. Dooley and A. K. Nandi, “On explicit time de- [18] T. A. Ramstad, “Digital methods for conversion between lay estimation using the Farrow structure,” Signal Process- arbitrary sampling frequencies,” IEEE Trans. Acoust. ing, vol. 72, pp. 53−57, Jan. 1999. Speech, Signal Processing, vol. ASSP−32, pp. 577−591, [26] J. Vesma, “A frequency-domain approach to poly- June 1984. nomial-based interpolation and the Farrow structure,” to [19] T. A. Ramstad, “Fractional rate decimator and interpola- appear IEEE Trans. on Circuits and Systems II, March tor design,” in Proc. IX European Signal Processing Conf., 2000. Rhodes, Greece, Sep. 1998, pp. 1949−1952. [27] J. Vesma, Optimization and Applications of Polyno- [20] R. W. Schafer and L. R. Rabiner, “A digital signal mial-Based Interpolation Filters. Dr. Tech. Thesis, Tam- processing approach to interpolation,” Proc. IEEE, vol. 61, pere, Finland: Tampere University of Tech., Department of pp. 692−702, June 1973. Information Technology, May 1999 [21] M. Unser, A. Aldroubi, and M. Eden, “Fast B-spline transforms for continuous image representation and inter- polation,” Trans. Pat. Anal., Mach. Int., vol. 13, pp. 277−285, Mar. 1991.

DOCUMENT INFO

Shared By:

Categories:

Tags:
Lagrange Interpolation Example Lagrange Interpolation How to, Lagrange Interpolation, Data Points, Polynomial Interpolation, Linear Interpolation, lagrange indiana, lagrange ky, lagrange illinois, lagrange georgia, lagrange multiplier

Stats:

views: | 745 |

posted: | 6/12/2009 |

language: | Finnish |

pages: | 19 |

OTHER DOCS BY rossmanjerry

How are you planning on using Docstoc?
BUSINESS
PERSONAL

By registering with docstoc.com you agree to our
privacy policy and
terms of service, and to receive content and offer notifications.

Docstoc is the premier online destination to start and grow small businesses. It hosts the best quality and widest selection of professional documents (over 20 million) and resources including expert videos, articles and productivity tools to make every small business better.

Search or Browse for any specific document or resource you need for your business. Or explore our curated resources for Starting a Business, Growing a Business or for Professional Development.

Feel free to Contact Us with any questions you might have.