2.161: Lecture 21 Notes

Massachusetts Institute of Technology Department of Mechanical Engineering 2.161 Signal Processing - Continuous and Discrete Lecture 21 Reading: • Class Handout: Interpolation (Up-sampling) and Decimation (Down-sampling) • Proakis and Manolakis: Secs. 11.1 – 11.5, 12.1 • Oppenheim, Schafer, and Buck: Sec. 4.6, Appendix A • Stearns and Hush, Ch. 9 1 1.1 Interpolation and Decimation Up-Sampling (Interpolation) by an Integer Factor Consider a data set {fn } of length N , where fn = f (nΔT ), n = 0 . . . N − 1 and ΔT is the sampling interval. The task is to resample the data at a higher rate so as to create a new data set {fˆn }, of length KN , representing samples of the same continuous waveform f (t), sampled at intervals ΔT /K. The following figure shows a cosinusoidal data set with N = 8 samples, and resampled with the same data set interpolated by a factor K = 4. 1 1 (a) original data set (8 samples) 0.5 0.5 0 0 −0.5 −0.5 −1 2 4 6 8 21–1 −1 (b) interpolated data set (32 samples) 10 20 30 1.1.1 A Frequency Domain Method This method is useful for a finite-sized data record. Consider the DFTs {Fm } and {F̂m } of a pair of sample sets {fn } and {fˆn }, both recorded from f (t) from 0 ≤ t < T , but with sampling intervals ΔT and ΔT /K respectively. Let N and KN be the corresponding sample sizes. It is assumed that ΔT has been chosen to satisfy the Nyquist criterion: • Let F (jΩ) = F {f (t)} be the Fourier transform of f (t), and let f (t) be sampled at intervals ΔT to produce f ∗ (t). Then    ∞ 1  2πn F (jΩ) = F j Ω− ΔT n=0 ΔT ∗ (1) is periodic with period 2π/ΔT , and consists of scaled and shifted replicas of F (jΩ). Let the total sampling interval be T to produce N = T /ΔT samples. • If the same waveform f (t) is sampled at intervals ΔT /K to produce fˆ∗ (t) the period of its Fourier transform F̂ ∗ (jΩ) is 2πK/ΔT and    ∞ K  2πKn F̂ (jΩ) = F j Ω− ΔT n=0 ΔT ∗ (2) which differs only by a scale factor, and an increase in the period. Let the total sampling period be T as above, to generate KN samples. • We consider the DFTs to be sampled representations of a single period of F ∗ (jΩ) and F̂ ∗ (jΩ). The equivalent line spacing in the DFT depends only on the total duration of the sample set T , and is ΔΩ = 2π/T in each case:    2πm ∗ Fm = F j , m = 0, 1, . . . N − 1 T    2πm ∗ , m = 0, 1, . . . KN − 1. F̂m = F̂ j T   From Eqs. (1) and (2) the two DFTs {Fm } and F̂m are related: ⎧ m = 0, 1, . . . , N/2 − 1 ⎨ KFm m = N/2, . . . , N K − N/2 − 1 0 F̂m = ⎩ KFm−(K−1)N m = N K − N/2, . . . , KN − 1 • The effect of increasing N (or decreasing ΔT ) in the sample set, while maintaining T = N ΔT constant, is to increase the length of the DFT by raising the effective Nyquist frequency ΩN . π Nπ = ΩN = ΔT T 21–2 F m N u m b e r o f s a m p le s : N S a m p lin g in te r v a l: D T N y q u is t fr e q u e n c y 0 -N /2 N /2 m 2 N 3 N /2 N (a ) F m N u m b e r o f s a m p le s : 2 N S a m p lin g in te r v a l: D T /2 N y q u is t fr e q u e n c y 0 m 2 N N (b ) The above figure demonstrates these effects by schematically, by comparing the DFT of (a) a data set of length N derived by sampling at intervals ΔT , and (b) a data set of length 2N resulting from sampling at intervals ΔT /2. The low frequency region of both spectra are similar, except for a scale factor, and the primary difference lies in the “high frequency” region, centered around the Nyquist frequency, in which all data points are zero. The above leads to an algorithm for the interpolation of additional points into a data set, by a constant factor K: 1. Take the DFT of the original data set to create {Fm } of length N . 2. Insert (K − 1)N zeros into the center of the DFT to create a length KN array. 3. Take the IDFT of the expanded array, and scale the sequence by a factor K. 1.1.2 A Time-Domain Method We now examine an interpolation scheme that may implemented on a sample-by-sample basis in real-time using time domain processing alone. As before, assume that the process f (t) is sampled at intervals ΔT , generating a sequence {fn } = {f (nΔT )}. Now assume that K − 1 zeros are inserted between the samples to form a sequence {fˆk } at intervals ΔT /K. f n 3 2 1 0 f n n 0 1 2 3 4 (a) Initial data set, N=8 5 6 7 3 2 1 0 n 0 5 10 15 20 (b) Data set with two samples interpolated between samples, N=24 21–3 This is illustrated above, where a data record with N = 8 samples has been expanded by a factor K = 3 to form a new data record of length N = 24 formed by inserting two zero samples between each of the original data points. We now examine the effect of inserting K − 1 samples with amplitude 0 after each sample. The DFT of the original data set is Fm = N −1  fn e−j 2πmn N m = 0...N − 1 , n=0 and for the extended data set {fˆn }, n = 0 . . . KN − 1 F̂m = KN −1  2πmn fˆn e−j KN , m = 0 . . . KN − 1 n=0 However, only the original N samples contribute to the sum, so that we can write F̂m = N −1  2πmk fˆKk e−j N k=0 m = 0 . . . KN − 1 = Fm , since fˆKk = fk . We note that {Fm } is periodic with period N , and {F̂m } is periodic with period KN , so that {F̂m } will contain K repetitions of {Fm }. F ^ m lo w - p a s s filte r p a s s b a n d r e p lic a tio n s o f F N y q u is t fre q u e n c y m 0 w F ^ K N /2 = K p /, T K N w = 2 K p /, T lo w - p a s s d ig ita l filte r m m 0 w K N /2 = K p /, T w K N = 2 K p /, T The magnitude of the DFTs of the two waveforms is shown above. The effect of inserting the K −1 zeros between the original samples has been to generate a waveform with an equivalent sampling interval of ΔT /K s, and a Nyquist frequency of Kπ/ΔT rad/s. The line resolution is unchanged, and the original DFT {Fm } is replicated K times within the frequency span of 2Kπ/ΔT rad/s. 21–4 F 10 m 8 6 4 2 0 m 0 1 2 3 4 (a) DFT of Initial data set, N=8 5 6 7 ^ F 10 m 8 6 4 2 0 m 0 5 10 15 20 (b) DFT of data set with two samples interpolated between samples, N=24 The fully interpolated waveform may be reconstructed by elimination of the replications of the original spectral components. While this might be done in the frequency domain, the most common method is to low-pass filter the padded data sequence to retain only the base-band portion of the spectrum as shown below. , T /K , T f(t) F ( jW ) A n ti- a lia s in g lo w - p a s s filte r H ( jW ) = 0 , |W |> p /, T ~ f(t) F ( jW ) H ( jW ) c o n tin u o u s d o m a in ^ K S a m p le r , T { fn } {fn } In s e rt K -1 z e ro s b e tw e e n s a m p le s L o w - p a s s d ig ita l filte r H (z ) d is c r e te d o m a in in te r p o la te d w a v e fo rm , T /K 1.2 Down-Sampling (Decimation) by an Integer Factor Decimation by an integer factor K is the reverse of interpolation, that is increasing the sampling interval T by an an integer factor (or decreasing the sampling frequency). At first glance the process seems to be simple: simply retain every Kth sample from the data sequence so that f˜n = fnK where {f˜n } is the down-sampled sequence, as shown in below. L D T D T t D o w n - s a m p le r t L Caution must be taken however to prevent aliasing in the decimated sequence. It is not valid to directly down-sample a sequence directly unless it is known a-priori that the spectrum of the data set is identically zero at frequencies at and above the Nyquist frequency defined by the lower sampling frequency. 21–5 In discussing sampling of continuous waveforms, we described the use of a pre-aliasing filter to eliminate (or at least significantly reduce) spectral components that would introduce aliasing into the sampled data set. When down-sampling, the digital equivalent is required: a digital low-pass filter is used to eliminate all spectral components that would cause aliasing in the resampled data set. The complete down-sampling scheme is: D T D T L D T t t D o w n - s a m p le r D ig ita l lo w - p a s s a n ti- a lia s in g filte r 1.3 L Resampling with a non-integer factor Assume that the goal is to re-sample a sequence by a non-integer factor p that can be expressed as a rational fraction, that is P = N M where N and M are positive integers. This can be achieved by (1) interpolation by a factor N , followed by (2) decimation by a factor M , as shown below. in te r p o la tio n filte r u p - s a m p le r { fn } H N u p a n ti- a lia s in g filte r (z ) H d n d o w n - s a m p le r (z ) ^ { fm } M d e c im a tio n in te r p o la tio n However, since the two low-pass filters are cascaded, they may be replaced with a single filter with a cut-off frequency that is the lower of the two filters, as is shown below { fn } 2 u p - s a m p le r N lo w - p a s s filte r H lp (z ) d o w n - s a m p le r M ^ { fm } Introduction to Random Signals In dealing with physical phenomena and systems we are frequently confronted with non­ deterministic, (stochastic, or random) signals, where the temporal function can not be de­ scribed explicitly, nor predicted. Some simple examples are wind loading on a structure, additive noise in a communication system, and speech waveforms. In this brief examina­ tion of random phenomena we concentrate on common statistical descriptors of stochastic waveforms and input-output relationships of linear systems excited by random waveforms. L in e a r s y s te m n o n - d e r te r m in is tic in p u t n o n - d e r te r m in is tic o u tp u t 21–6 Since we cannot describe f (t), we must use statistical descriptors that capture the essence of the waveform. There are two basic methods of doing this: (a) Describe the waveform based on temporal measurements, for example define the mean μ of the waveform as 1 T /2 μ = lim f (t) dt T →∞ T −T /2 (b) Conjecture an ensemble of random processes fi (t), i = 1, . . . N , with identical statistics and define the descriptors my measurements made across the ensemble at a given time, for example ∞  fi (t). μ(t) = lim N →∞ N =1 • A stationary process is defined as one whose ensemble statistics are indepen­ dent of time. • An ergodic process is one in which the temporal statistics are identical to the ensemble statistics. Clearly, ergodicity impies stationarity. f1 (t) tim e - b a s e d s ta tis tic s d e r iv e d fr o m a s in g le r a n d o m p ro c e s s f2 (t) f3 (t) e n s e m b le - b a s e d s ta tis tic s d e r iv e d fr o m a h y p o th e tic a l c o lle c tio n o f r a n d o m p r o c e s s e s w ith id e n tic a l d e s c r ip to r s In practice statistical descriptors are usually derived experimentally from measurements. For example, the mean of of a waveform might be estimated from a set of 1000 samples of a waveform and computed as 1000 1  fi . μ̂1 = 1000 i=1 But μ̂ is not the mean, it is simply an estimator of the true mean. If we repeated the experiment, we would come up with a different value μ̂2 . In statistical descriptions we use the terms expected value, or expectation, designated E {x}, and say E {μ̂i } = μ to indicate that our experimental estimates μ̂i will be clustered around the true mean μ. 21–7 2.1 Ensemble Based Statistics 2.1.1 The Probability Density Function (pdf ) The pdf is strictly an ensemble statistic that describes the distribution of samples x across the amplitude scale. Its definition is 1 Prob {xa ≤ x ≤ xa + Δx} Δx→0 Δx p(xa ) = lim so that the probability that a single sample lies in the range a ≤ x ≤ b is Prob {a ≤ x ≤ b} = and we note Prob {−∞ ≤ x ≤ ∞} = b a p(x) dx, ∞ −∞ p(x) dx = 1. For an ergodic process, the pdf can also be described from a single time series x x (t) a re a = P ro b {a £ x £ b } b b a t a p (x ) and may be interpreted as the fraction of time that the waveform “dwells” in the range b ≤ f (t) ≤ b. Two common pdfs are (a) The Uniform distribution a random sample taken from a uniformly distributed ran­ dom process is equally likely to be found anywhere between a minimum and maximum value. 1 a≤x≤b p(X) = b−a 0 elsewhwere. p (x ) 1 b - a a b 21–8 x (b) The normal (or gaussian) distribution] The normal distribution defines the well known “bell-shaped-curve” of elementary statistics p(x) = √ (x−μ)2 1 e− 2σ2 2πσ where μ is the mean of the distribution, and σ 2 is the variance. p (x ) m - s m x m + s Note: The central limit theorem of statistics states that any random process that is the sum of a large number of underlying independent random processes, regardless of their distributions will be described by a gaussian distribution. Many ensemble based statistical descriptors may be described in terms of the pdf, for example The mean E {x} = μ = The variance 2 E (x − μ) 2 =σ = ∞ −∞ xp(x)dx ∞ −∞ (x − μ)2 p(x)dx and by expanding the integral 2 σ = ∞ −∞ 2 x p(x)dx − 2μ ∞ −∞ xp(x)dx + μ 2 ∞ −∞ p(x)dx = E x2 − μ2 Example 1 Find the mean and variance a random variable that is uniformly distributed between x0 and x0 + Δ. The pdf is 1 x0 ≤ x ≤ x0 + Δ p(X) = Δ 0 elsewhere. 21–9 The mean is ∞ μ= −∞ x0+Δ xp(x)dx = x0 x Δ dx = x0 + Δ 2 which is the mid-point of the range. The variance is 1 σ 2 = E x2 − μ2 = Δ2 12 2.2 Time-based Statistics Two stochastic waveforms may have identical pdfs (and hence equal means and variances), but be very different qualitatively, for example f1 (t) t f2 (t) t t o to + t These two waveforms obviously differ in • spectral content, or • self-similarity between themselves at a time t0 and some time τ later. A random waveform cannot be predicted exactly at any time, but clearly in the above figure the upper waveform (with greater high frequency spectral content) has less self-similarity, or correlation with itself, after a delay of τ . The correlation functions are a measure of the degree to which the value of a function depends upon its past. For infinite duration waveforms the auto-correlation function φf f (τ ) is defined as 1 T /2 f (t)f (t + τ ) dt φf f (τ ) = lim T →∞ T −T /2 21–10 and is a measure of the self-similarity of the function f (t) at time t and at a time τ later. The cross-correlation function φf g (τ ) measures the similarity between two different functions f (t) and g(t) at two times τ apart. 1 φf g (τ ) = lim T →∞ T � T /2 f (t)g(t + τ ) dt −T /2 Note that these definitions must be modified for finite duration waveforms, and if f (t) exists in the interval T1 ≤ t ≤ T2 , we define � T2 φf f (τ ) = f (t)f (t + τ ) dt T1 21–11