The Sliding DFT

12 downloads 6678 Views 531KB Size Report
This filter computes a single DFT output (the kthbin of an N-point DFT) defined by ..... k = N.fi/fs, when fi is a frequency of interest in Hz. In addition, the SDFT ...
Chapter 18 Sliding Spectrum Analysis Eric Jacobsen

Richard Lyons

Abineau Communications

Besser Associates

The standard method for spectrum analysis in DSP is the discrete Fourier transform (DFT), typically implemented using a fast Fourier transform (FFT) algorithm. However, there are applications that require spectrum analysis only over a subset of the N center frequencies of an N-point DFT. A popular, as well as efficient, technique for computing sparse DFT results is the Goertzel algorithm which computes a single complex DFT spectral bin value for every N input time samples. This chapter describes sliding spectrum analysis techniques whose spectral bin output rates are equal to the input data rate, on a sample-by-sample basis, with the advantage that they require fewer computations than the Goertzel algorithm for real-time spectral analysis. In applications where a new DFT output spectrum is desired every sample, or every few samples, the sliding DFT is computationally simpler than the traditional radix-2 FFT. We'll start our discussion by providing a brief review of the Goertzel algorithm, and use its behavior as a yardstick to evaluate the performance of the sliding DFT technique. Following that, we will examine stability issues regarding the sliding DFT implementation as well as review the process of frequency-domain convolution to accomplish time-domain windowing. Finally, a modified sliding DFT structure is proposed that provides improved computational efficiency.

18.1 Goertzel Algorithm The Goertzel algorithm, used in dual-tone multi-frequency decoding and PSK/FSK modem implementations, is commonly used to compute DFT spectra [1]–[4]. The algorithm is implemented in the form of a 2nd-order IIR filter as shown in Figure 18–1. This filter computes a single DFT output (the kth bin of an N-point DFT) defined by

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

N 1

X (k )   x(n)e j 2πnk / N .

(18–1)

n 0

The filter's y(n) output is equal to the DFT output frequency coefficient, X(k), at the time index n = N. For emphasis, we remind the reader that the filter's y(n) output is not equal to X(k) at any time index when n  N. The frequency-domain index k is an integer in the range 0  k  N–1. The derivation of this filter's structure is readily available in the literature [5]–[7].

y(n)

v(n)

x(n) – 2cos(2k/N)

z–1

–e–j2k/N

y(N) = X(k)

z–1

Figure 18–1 IIR filter implementation of the Goertzel algorithm.

The z-domain transfer function of the Goertzel filter is

HG ( z) 

1  e j 2πk / N z 1 1  2cos(2πk / N ) z 1  z 2

(18–2)

with a single z-domain zero located at z = e–j2k/N and conjugate poles at z = e±j2k/N as shown in Figure 18–2(a). The pole/zero pair at z = e–j2k/N cancel each other. The frequency magnitude response, provided in Figure 18–2(b), shows resonance centered at a normalized frequency of 2k/N, corresponding to a cyclic frequency k.fs/N hertz (where fs is the signal sample rate).

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

Imaginary part

(a)

1

z-plane

0.5 2k/N

0 –0.5 –1 –2

–1

0 Real part

1

2

0

(b)

dB

–10 –20 –30 k–2 k–1

k k+1 k+2 Frequency

Figure 18–2 Goertzel filter: (a) z-domain pole/zero locations; (b) frequency magnitude response.

We remind the reader that while the typical Goertzel algorithm description in the literature specifies the frequency resonance variable k in (18–2) and Figure 18–1 to be an integer (making the Goertzel filter's output equivalent to an N-point DFT bin output). Variable k can in fact be any value between 0 and N–1 giving us full flexibility in specifying a Goertzel filter's resonance frequency. While the Goertzel algorithm is derived from the standard DFT equation, it's important to realize that the filter's frequency magnitude response is not the sin(x)/(x)like response of a single-bin DFT. The Goertzel filter is a complex resonator having an infinite-length unit impulse response, h(n) = ej2nk/N, and that's why its magnitude response is so narrow. The time-domain difference equations for the Goertzel filter are

v(n) = 2cos(2k/N)v(n–1) –v(n–2) + x(n).

(18–3a)

y(n) = v(n) –e–j2k/Nv(n–1).

(18–3b)

An advantage of the Goertzel filter in calculating an N-point X(k) DFT bin is that (18–3a) is implemented N times while (18–3b), the feedforward path in Figure 18– 1, need only be computed once after the arrival of the Nth input sample. Thus for real From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

x(n) the filter requires N + 2 real multiplies and 2N + 1 real adds to compute an N-point X(k). However, when modeling the Goertzel filter if the time index begins at n = 0, the filter must process N+1 time samples with x(N) = 0 to compute X(k). Now let's look at the sliding DFT process.

18.2 Sliding Discrete Fourier Transform (SDFT) The sliding DFT (SDFT) algorithm performs an N-point DFT on time samples within a sliding-window as shown in Figure 18–3. In this example the SDFT initially computes the DFT of the N = 16 time samples in Figure 18–3(a). The time window is then advanced one sample, as in Figure 18–3(b), and a new N-point DFT is calculated. The value of this process is that each new DFT is efficiently computed directly from the results of the previous DFT. The incremental advance of the time window for each output computation is what leads to the name sliding DFT or sliding-window DFT.

Window 1 (a)

0 –1 Time

x(–15)

x(0)

Shifted window 1 (b)

0 –1 x(–16)

x(–15)

Time

x(0)

Figure 18–3 Signal windowing for two 16-point DFTs: (a) data samples in the first computation; (b) second computation samples.

The principle used for the SDFT is known as the DFT shifting theorem, or the circular shift property [8]. It states that if the DFT of a windowed (finite-length) timedomain sequence is X(k), then the DFT of that sequence, circularly shifted by one sample is X(k)ej2k/N. Thus the spectral components of a shifted time sequence are the original (unshifted) spectral components multiplied by ej2k/N, where k is the DFT bin of From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

interest. We use this shift principal to express our sliding DFT process as

Sk(n) = ej2k/N[Sk(n–1) + x(n) –x(n–N)]

(18–4)

where Sk(n) is the new spectral component and Sk(n–1) is the previous spectral component. The subscript k reminds us that the spectra are those associated with the kth DFT bin. Equation (18–4), whose derivation is provided in the Appendix, reveals the value of this process in computing real-time spectra. We calculate Sk(n) by phase shifting the sum of the previous Sk(n–1) with the difference between the current x(n) sample and the x(n–N) sample. The difference between the x(n) and x(n–N) samples can be computed once for each n and used for each Sk(n) computation. So the SDFT requires only one complex multiply and two real adds per output sample. The computational complexity of each successive N-point output is then O(N) for the sliding DFT compared to O(N2) for the DFT and O[Nlog2(N)] for the FFT. Unlike the DFT or FFT, however, due to its recursive nature the sliding DFT output must be computed for each new input sample. If a new N-point DFT output is required only every N inputs, the sliding DFT requires O(N2) computations and is equivalent to the DFT. When output computations are required every M input samples, and M is less than log2(N), the sliding DFT can be computationally superior to traditional FFT implementations even when all N DFT outputs are required. Equation (18–4) leads to the single-bin SDFT filter structure shown in Figure 18–4.

e j2k/N

x(n)

Sk(n)

– z–N

z–1 Sk(n–1)

Figure 18–4 Single-bin sliding DFT filter structure. The single-bin SDFT algorithm is implemented as an IIR filter with a comb filter followed by a complex resonator [9]. (If you want to compute all N DFT spectral From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

components, N resonators with k = 0 to N–1 will be needed, all driven by a single comb filter.) The comb filter delay of N samples forces the filter's transient response to be N– 1 samples in length, so the output will not reach steady state until the Sk(N) sample. In practical applications the algorithm can be initialized with zero-input and zero-output. The output will not be valid, or equivalent to (18–1)'s X(k), until N input samples have been processed. The z-domain transfer function for the kth bin of the sliding DFT filter is

HSDFT ( z ) 

e j 2πk / N (1  z  N ) . 1  e j 2πk / N z 1

(18–5)

This complex filter has N zeros equally-spaced around the z-domain's unit circle, due to the N-delay comb filter, as well as a single pole canceling the zero at z = ej2k/N as shown in Figure 18–5(a). The SDFT filter's complex h(n) unit impulse response is shown in Figure 18–5(b) for the example where k = 2 and N = 20.

Imaginary part

1

(a)

0.5 2k/N

0 –0.5 k=2 –1 N = 20 –1

0 Real part

1

1 0 –1 0

2

4

6

8 10 12 Real[h(n)]

14

16

18

20

2

4

6

8 10 12 Imag[h(n)]

14

16

18

20

(b) 1 0 –1 0

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

Figure 18–5 Sliding DFT characteristics for k = 2 and N = 20: (a) impulse response; (b) pole/zero locations.

Because of the comb subfilter, the SDFT filter's complex sinusoidal unit impulse response is finite in length—truncated in time to N samples—and that property makes the frequency magnitude response of the SDFT filter identical to the sin(Nx)/sin(x) response of a single DFT bin centered at a normalized frequency of 2k/N. We've encountered a useful property of the SDFT that's not widely known, but is important. If we change the SDFT's comb filter feedforward coefficient (in Figure 18– 4) from –1 to +1, the comb's zeros will be rotated counterclockwise around the unit circle by an angle of /N radians. This situation, for N = 8, is shown on the right side of Figure 18–6(a). The zeros are located at angles of 2(k + 1/2)/N radians. The k = 0 zeros are shown as solid dots. Figure 18–6(b) shows the zeros locations for an N = 9 SDFT under the two conditions of the comb filter's feedforward coefficient being –1 and +1.

2/N

0 –1

Imaginary part

(b)

1 N=8 3/N

0 –1

–1

N odd

Imaginary part

(a)

1 N=8

Comb coefficient = +1

0 Real part

–1

1

1 N=9

Imaginary part

N even

Imaginary part

Comb coefficient = –1

4/N 0 –1

0 Real part

1

1 N=9 5/N 0 –1

–1

0 Real part

1

–1

0 Real part

1

Figure 18–6 Four possible orientations of comb filter zeros on the unit circle.

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

This alternate situation is useful, we can now expand our set of spectrum analysis center frequencies to more than just N angular frequency points around the unit circle. The analysis frequencies can be either 2k/N or 2(k+1/2)/N , where integer k is in the range 0  k  N–1. Thus we can build an SDFT analyzer that resonates at any one of 2N frequencies between 0 and fs Hz. Of course, if the comb filter's feedforward coefficient is set to +1, the resonator's feedforward coefficient must be ej2(k+1/2)/N to achieve pole/zero cancellation. One of the attributes of the SDFT is that once an Sk(n–1) is obtained, the number of computations to calculate Sk(n) is fixed and independent of N. A computational workload comparison between the Goertzel and SDFT filters is provided later in this discussion. Unlike the radix-2 FFT, the SDFT's N can be any positive integer giving us greater flexibility to tune the SDFT's center frequency by defining integer k such that k = N.fi/fs, when fi is a frequency of interest in Hz. In addition, the SDFT requires no bitreversal processing as does the FFT. Like Goertzel, the SDFT is especially efficient for narrowband spectrum analysis. For completeness, we mention that a radix-2 sliding FFT technique exists for computing all N bins of X(k) in (18–1) [10], [11]. This method is computationally attractive because it requires only N complex multiplies to update the N-point FFT for all N bins, however it requires 3N memory locations (2N for data and N for twiddle coefficients). Unlike the SDFT, the radix-2 sliding FFT scheme requires address bitreversal processing and restricts N to be an integer power of two.

18.3 SDFT Stability The SDFT filter is only marginally stable because its pole resides on the z-domain's unit circle. If filter coefficient numerical rounding error is not severe, the SDFT is boundedinput-bounded-output stable. Filter instability can be a problem, however, if numerical coefficient rounding causes the filter's pole to move outside the unit circle. We can use a damping factor r to force the pole to be at a radius of r inside the unit circle and guarantee stability using a transfer function of

HSDFT,gs ( z ) 

re j 2πk / N (1  r N z  N ) 1  re j 2πk / N z 1

(18–6)

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

with the subscript "gs" meaning guaranteed-stable. The stabilized feedforward and feedback coefficients become –rN and rej2k/N respectively. The difference equation for the stable SDFT filter becomes

Sk(n) = rej2k/N [Sk(n–1) + x(n) –rNx(n–N)]

(18–7)

with the stabilized-filter structure shown in Figure 18–7.

re j2k/N

x(n)

Sk(n)

– z–N

rN

z–1 Sk(n–1)

Figure 18–7 Guaranteed-stable sliding DFT filter structure.

Using a damping factor as in Figure 18–7 guarantees stability, but the Sk(n) output, defined by

N 1

X r 1 (k )   x(n)r ( N n ) e j 2πnk / N n 0

(18–8)

is no longer exactly equal to the kth bin of an N-point DFT in (18–1). While the error is reduced by making r very close to (but less than) unity, a scheme does exist for eliminating that error completely once every N output samples at the expense of additional conditional logic operations [12]. Determining if the damping factor r is necessary for a particular SDFT application requires careful empirical investigation. Another stabilization method worth consideration for fixed-point applications is decrementing the largest component (either real or imaginary) of the filter's ej2k/N feedback coefficient by one least significant bit. This technique can be applied selectively to problematic output bins and is effective in combating instability due to rounding errors which result in finite-precision ej2k/N coefficients having magnitudes From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

greater than unity. Like the DFT, the SDFT's output is proportional to N, so in fixed-point binary implementations the designer must allocate sufficiently-wide registers to hold the computed results.

18.4 Time-Domain Windowing In the Frequency Domain The spectral leakage of the SDFT can be reduced by the standard concept of windowing the x(n) input time samples. However, windowing by time-domain multiplication would compromise the computational simplicity of the SDFT. Alternatively, we can implement a time-domain window by means of frequency-domain convolution. Spectral leakage reduction performed in the frequency domain is accomplished by convolving adjacent Sk(n) values with the DFT of a window function. For example, the DFT of a Hanning window comprises only three non-zero values, –0.25, 0.5, and –0.25. As such we can compute a Hanning-windowed Sk(n), the kth DFT bin, with a threepoint convolution using

Hanning-windowed Sk(n) = –0.25.Sk–1(n) +0.5.Sk(n) –0.25.Sk+1(n).

(18–9)

Figure 18–8 shows this process where the comb filter stage need only be implemented once. Thus a Hanning window can be implemented by binary right shifts and two complex adds for each SDFT bin, making the Hanning window attractive in ASIC and FPGA implementations where single-cycle hardware multiplies are costly. If a gain of four is acceptable, then only two left shifts (one for the real part and one for the imaginary parts of Sk(n)) and two complex adds are required using

Hanning-windowed Sk(n) = –Sk–1(n) +2.Sk(n) –Sk+1(n).

(18–10)

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

Resonator k–1 x(n) – z–N

Sk–1(n)

Resonator Sk(n) k Resonator Sk+1(n) k+1

–0.25

0.5

Windowed Sk(n) output

–0.25

Figure 18–8 Three-resonator structure to compute three SDFT bin results, and a three-point convolution. The Hanning window is a member of a category called cos(x) window functions [13], [14]. These functions are also known as generalized cosine windows because their N-point time-domain samples are defined as

w(n) 

 1

(1)m am cos(2πmn / N )  m 0

(18–11)

where n = 0, 1, 2, ..., N–1, and the integer  specifies the number of terms in the window's time function. These window functions are attractive for frequency domain convolution because their DFTs contain only a few non-zero samples. The frequency domain vectors of various cos(x) window functions follow the form (1/2)·(a2, – a1, 2a0, –a1, a2), with a few examples presented in Table 18–1. Additional cos(x) window functions are described in the literature [14]. Table 18–1. cos(x) windows, frequency domain coefficients

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

Window function:

a0

a1

a2

Rectangular

1.0





Hanning, ( = 2)

0.5

0.5



Blackman, ( = 3)

0.42

0.5

0.08

Exact Blackman,

0.4265907 0.4965606 0.0768487

( = 3)

Hamming,

0.54



0.46

18.5 Sliding Goertzel DFT We can reduce the number of multiplications required in the SDFT by creating a new pole/zero pair in its HSDFT(z) system function [15]. This is done by multiplying the numerator and denominator of HSDFT(z) in (18–5) by the factor (1 –e–j2k/Nz–1) yielding

HSG ( z ) 

(1  e j 2πk / N z 1 )(1  z  N ) (1  e j 2πk / N z 1 )(1  e j 2πk / N z 1 )



(1  e j 2πk / N z 1 )(1  z  N ) 1  2cos(2πk / N ) z 1  z 2

(18–12)

where the subscript "SG" means sliding Goertzel. The filter block diagram for HSG(z) is shown in Figure 18–9 where this new filter is recognized as the standard Goertzel filter preceded by a comb filter. The sliding Goertzel DFT filter, unlike the standard Goertzel filter, has a finite-duration impulse response identical to that shown in Figure 18–5(b), for k = 2 and N = 20.

y(n)

v(n)

x(n) –



–N

z

2cos(2k/N)

z–1

–e–j2k/N

z–1

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

Figure 18–9 Structure of the sliding Goertzel DFT filter.

Of course, unlike the traditional Goertzel filter in Figure 18–1, the sliding Goertzel DFT filter's complex feedforward computations must be performed for each input time sample. The sliding Goertzel filter's sin(Nx)/sin(x) frequency magnitude response, for k = 2 and N = 20, is provided in Figure 18–10(a). The asymmetrical frequency response is defined by the filter's N zeros equally-spaced around the z-domain's unit circle in Figure 18–10(b) due to the N-delay comb filter, as well as an additional (uncanceled) zero located at z = e–j2k/N on account of the (1 –e–j2k/Nz–1) factor in the HSG(z) transfer function's numerator. In addition, the filter has conjugate poles canceling zeros at z = e±j2k/N. 0 k=2 N = 20

dB

–5

(a)

–10 –15 –20 –25 –30 0

fs/4

fs/2 Frequency

3fs/4

fs

(b)

Imaginary part

1 0.5 2k/N

0 –0.5 k=2 –1 N = 20 –1

two zeros 0 Real part

1

Figure 18–10 Sliding Goertzel filter for N = 20 and k = 2: (a) frequency magnitude response; (b) z-domain pole/zero locations.

The sliding Goertzel DFT filter is of interest because its computational workload is less than that of the SDFT. This is because the v(n) samples in Figure 18–9 are real-only due to the real-only feedback coefficients. A single-bin DFT computational comparison, for real-only inputs, is provided in Table 18–2. For real-time processing requiring From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

spectral updates on a sample by sample basis the sliding Goertzel method requires fewer multiplies than either the SDFT or the traditional Goertzel algorithm.

Table 18–2. Single-bin DFT comparison

Method

DFT Goertzel Sliding DFT Sliding Goertzel

Single Sk(n)

Next Sk(n+1)

computation:

computation:

Real

Real

Real

Real

multiplies:

adds:

multiplies:

adds:

2N

2N

2N

2N

N+2

2N + 1

N+2

2N + 1

4N

4N

4

4

N+2

3N + 1

3

4

18.6 Conclusions The sliding DFT process for spectrum analysis was presented, and shown to be more efficient than the popular Goertzel algorithm for sample-by-sample DFT bin computations. The sliding DFT provides computational advantages over the traditional DFT or FFT for many applications requiring successive output calculations, especially when only a subset of the DFT output bins are required. Methods for output stabilization as well as time-domain data windowing by means of frequency-domain convolution were also discussed. A modified sliding DFT algorithm, called the sliding Goertzel DFT, was proposed to further reduce computational workload.

18.7 References [1] M. Felder, J. Mason, and B. Evans, "Efficient dual-tone multi-frequency detection using the non-uniform discrete Fourier transform," IEEE Signal Processing Lett., Vol. 5, pp. 160-163, July 1998. [2] Analog Devices Inc., ADSP-2100 Family User’s Manual, (3rd Ed), Chapter 14; Vol. 1, 1995, Norwood, MA, http://www.analog.com/Analog_Root/static/library/dspManuals/Using_ADSP2100_Vol1_books.html

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

[3] S. Gay, J. Hartung, and G. Smith, “Algorithms for multi-channel DTMF detection for the WERDSP32 family,” in Proceedings on the International Conference on ASSP, pp. 1134-1137, 1989. [4] K. Banks, “The Goertzel algorithm,” Embedded Systems Programming Magazine, pp. 34-42, Sept. 2002. [5] G. Goertzel, “An algorithm for the evaluation of finite trigonometric series,” American Math. Monthly, Vol. 65, pp. 34-35, 1958. [6] J. Proakis and D. Manolakis, Digital Signal Processing-Principles, Algorithms, and Applications, Third Edition, Prentice Hall, Upper Saddle River, New Jersey, 1996, pp. 480-481. [7] A. Oppenheim, R. Schafer, and J. Buck, "Discrete-Time Signal Processing", 2nd ed. Prentice Hall, Upper Saddle River, New Jersey, 1996, pp. 633-634. [8] T. Springer, "Sliding FFT computes frequency spectra in real time", EDN Magazine, pp. 161-170, Sept. 29, 1988. [9] L. Rabiner and B. Gold, Theory and Application of Digital Signal Processing, Prentice Hall, Upper Saddle River, New Jersey, 1975, pp. 382-383. [10] B. Farhang-Boroujeny and Y. Lim, "A comment on the computational complexity of sliding FFT," IEEE Trans. Circuits and Syst. II, Vol. 39, No. 12, pp. 875-876, Dec. 1992. [11] B. Farhang-Boroujeny and S. Gazor, "Generalized sliding FFT and its application to implementation of block LMS adaptive filters," IEEE Trans. Sig. Proc., Vol. 42, No. 3, pp. 532-538, Mar. 1994. [12] S. Douglas and J. Soh, "A numerically-stable sliding-window estimator and its application to adaptive filters," Proc. 31st Annual Asilomar Conf. on Signals, Systems, and Computers, Pacific Grove, CA, vol. 1, pp. 111-115, November 1997. [13] F. Harris, "On the use of windows for harmonic analysis with the discrete Fourier transform," Proc. IEEE, Vol. 66, pp. 51-84, Jan. 1978. [14] A. Nuttall, "Some windows with very good sidelobe behavior," IEEE Trans., Vol. 29, No. 1, pp. 84-91, Feb. 1981. [15] K. Larson, Texas Instruments Inc., Private communication, Nov. 2002.

18.8 Appendix The derivation of the general SDFT can be understood starting with the definition of the DFT of a contiguous subset of a longer sequence. For this derivation we consider a DFT From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

of length N computed on a subset of an input sequence with length of at least N+q+1 where q is the index of the start of the DFT window in the input sequence. The additional time unit increment is merely to accommodate the slide to the next window. We start the derivation with the definition of the kth DFT bin for a window starting at the qth element of the input sequence: N 1

X (k , q)   x(n  q)e j 2πnk / N . n 0

(18–A1)

The transform of the (q+1)th window, the next DFT in the sliding sequence, is then: N 1

X (k , q  1)   x(n  q  1)e j 2πnk / N . n 0

(18–A2)

Substituting p = n+1, so the range of p is 1 to N and we have:

N

X (k , q  1)   x( p  q)e j 2π( p 1) k / N . p 1

(18–A3)

Next we change the summation by formally expressing the Nth component separately and adding the p = 0 case to the summation, and then subtracting it formally:

X (k , q  1) 

N 1

 x ( p  q )e p 0

 j 2π( p 1) k / N

+ x(q+N)e–j2(N–1)k/N – x(q)ej2k/N .

(18–A4)

The exponential terms can be factored as follows:  N 1 X (k , q  1)  e j 2πk / N   x( p  q)e  j 2πpk / N  p 0  x(q  N )e j 2πNk / N  x(q)  .

(18–A5)

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012

Because e–j2Nk/N = 1, (18–A5) can be rewritten as:

X(k,q+1) =  N 1  e j 2πk / N   x( p  q)e j 2πpk / N  x(q  N )  x(q).  p 0 

(18–A6)

Notice that the summation is the DFT of the qth time window, but (18–A6) is the DFT of the (q+1)th window. The sliding DFT can therefore be expressed as:

X(k,q+1) = ej2k/N [X(k,q) + x(q+N) – x(q)].

(18–A7)

The individual, kth, frequency bins for the (q+1)th SDFT window can therefore be computed from the qth window bins and the time-domain inputs by:

Sk(n) = ej2k/N[Sk(n–1) + x(n) –x(n–N)].

(18–A8)

The above derivation provides exact equivalence to the traditional DFT computation with a recursive algorithm that reduces the computational for successive sliding windows. ___________________________________________________________ ___________________________________________________________

Editor Comment: A variation of this sliding DFT process is discussed in Chapter 21.

From: "Streamlining Digital Signal Processing", © Copyright IEEE Press-John Wiley & Sons, 2012