1 / 39

ELEN E4810: Digital Signal Processing Topic 10: The Fast Fourier Transform

ELEN E4810: Digital Signal Processing Topic 10: The Fast Fourier Transform. Calculation of the DFT The Fast Fourier Transform algorithm Short-Time Fourier Transform. 1. Calculation of the DFT. Filter design so far has been oriented to time-domain processing - cheaper!

kita
Download Presentation

ELEN E4810: Digital Signal Processing Topic 10: The Fast Fourier Transform

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. ELEN E4810: Digital Signal ProcessingTopic 10: The Fast Fourier Transform • Calculation of the DFT • The Fast Fourier Transform algorithm • Short-Time Fourier Transform Dan Ellis

  2. 1. Calculation of the DFT • Filter design so far has been oriented to time-domain processing - cheaper! • But: frequency-domain processing makes some problems very simple: • Need an efficient way to calculate DFT! X[k] Fourier domain processing Y[k] DFT IDFT x[n] y[n] Dan Ellis

  3. ( ) The DFT • Recall the DFT: • discrete transform of discrete sequence • Matrix form: WN@2p/N WNrhas only N distinct values Lots of structure  opportunities for efficient algorithms Dan Ellis

  4. Computational Complexity • N cpx multiplies + N-1 cpx adds per ptN points (k = 0..N-1) • cpx mult: (a+jb)(c+jd) = ac - bd + j(ad+bc)= 4 real mults + 2 real adds • cpx add = 2 real adds • Total: 4N2 real mults, 4N2-2N real adds Dan Ellis

  5. Goertzel’s Algorithm • Now: • i.e.where looks like a convolution x[n] 0 ≤ n < N xe[n] = { 0 n = N WN-knn ≥ 0 hk[n] = { 0 n < 0 xe[n]xe[N] = 0 yk[n]yk[-1] = 0 yk[N] = X[k] + z-1 WN-k Dan Ellis

  6. Goertzel’s Algorithm • Separate ‘filters’ for each X[k] • can calculate for just a few samples of k • No large buffer, no coefficient table • Same complexity for full X[k](4N2mults, 4N2 - 2N adds) • but: can halve multiplies by making the denominator real: evaluate only for last step 2 real mults per step Dan Ellis

  7. 2. Fast Fourier Transform FFT • Reduce complexity of DFT from O(N2)to O(N·logN) • grows more slowly with larger N • Works by decomposing large DFT into several stages of smaller DFTs • Often provided as a highly optimized library Dan Ellis

  8. Decimation in Time (DIT) FFT • Can rearrange DFT formula in 2 halves: Arrange terms in pairs Group terms from each pair X0[<k>N/2] X1[<k>N/2] N/2 pt DFT of x for evenn N/2 pt DFT of x for oddn Dan Ellis

  9. Decimation in Time (DIT) FFT • Finite sequence x[n], 0 ≤ n < N, N = 2M • i.e. length is a power of 2 • Divide z-transform into parts coming from even and odd values of n: ZT of N/2 point sequence formed from even pts of x[n] ZT of N/2 point sequence formed from odd pts of x[n] Dan Ellis

  10. Decimation in Time (DIT) FFT • Now, but • Hence: N/2point DFT defined only for k = 0..N/2-1 k = 0..N-1 Dan Ellis

  11. Decimation in Time (DIT) FFT • We can evaluate an N-pt DFT as two N/2-pt DFTs (plus a few mults/adds) • But if DFTN{•} ~ O(N2)thenDFTN/2{•} ~ O((N/2)2) = 1/4 O(N2)  Total computation ~ 2  1/4 O(N2)= 1/2 the computation (+e) of direct DFT Dan Ellis

  12. One-Stage DIT Flowgraph “twiddle factors” always apply to odd-terms ouptut NOT mirror-image Even points from x[n] Same as X[0..3] except for factors on X1[•]terms Odd points from x[n] Classic FFT structure Dan Ellis

  13. Multiple DIT Stages • If decomposing one DFTN into two smaller DFTN/2’s speeds things up...Why not further divide into DFTN/4’s ?i.e.so have • Similarly, 0 ≤ k < N 0 ≤ k < N/2 N/4-pt DFT of odd points from even subset N/4-pt DFT of even points in even subset of x[n] Dan Ellis

  14. Two-Stage DIT Flowgraph Dan Ellis

  15. Multi-stage DIT FFT • Can keep doing this until we get down to 2-pt DFTs: N = 2M-pt DFT reduces to M stages of twiddle factors & summation (O(N2) part vanishes)  real mults < M·4N , real adds < 2M·2N  complexity ~O(N·M) = O(N·log2N) “butterfly” element X[0] = x[0] + x[1]  DFT2 1 = W20 X[1] = x[0] - x[1] -1 = W21 Dan Ellis

  16. FFT Implementation Details • Basic butterfly (at any stage): • Can simplify: XX0[r] XX[r] • • 2 cpx mults WNr • • • • XX1[r] XX[r+N/2] WNr+N/2 XX0[r] XX[r] just one cpx mult! XX1[r] XX[r+N/2] WNr -1 i.e. SUB rather than ADD Dan Ellis

  17. 8-pt DIT FFT Flowgraph • -1’s absorbed into summation nodes • WN0 disappears • ‘in-place’ algorithm: sequential stages bit-reversed indexing Dan Ellis

  18. FFT for Other Values of N • Having N = 2M meant we could divide each stage into 2 halves = “radix-2 FFT” • Same approach works for: • N = 3M radix-3 • N = 4M radix-4 - more optimized radix-2 • etc... • Composite N = a·b·c·d mixed radix(different N/rpoint FFTs at each stage) • .. or just zero-pad to make N = 2M Dan Ellis

  19. 1/N Re{X[k]} Re{x[n]} Re Re DFT Im{X[k]} Im{x[n]} Im Im -1 -1/N only differencesfrom forward DFT Inverse FFT • Recall IDFT: • Thus: • Hence, use FFT to calculate IFFT: Forward DFT of x'[n] = X*[k]|k=n i.e. time sequence made from spectrum Dan Ellis

  20. DFT of Real Sequences • If x[n] is pure-real, DFT wastes mult’s • Realx[n]Conj. symm.X[k] = X*[-k] • Given two real sequences, x[n] and w[n] call y[n] = j·w[n] , v[n] = x[n]+y[n] • N-pt DFTV[k] = X[k]+Y[k]but: V[k]+V*[-k] = X[k]+X*[-k]+Y[k]+Y*[-k]  X[k]=1/2(V[k]+V*[-k]) , W[k]=-j/2(V[k]-V*[-k]) • i.e. compute DFTs of twoN-pt real sequences with a singleN-pt DFT -Y[k] X[k] M Dan Ellis

  21. 3. Short-Time Fourier Transform (STFT) • Fourier Transform (e.g. DTFT) gives spectrum of an entire sequence: • How to see a time-varying spectrum? • e.g. slow AM of a sinusoid carrier: Dan Ellis

  22. Fourier Transform of AM Sine • Spectrum of whole sequence indicates modulation indirectly... • ... as cancellationbetween closely-tuned sines 2cAcB = cA+B +cA-B Nsin2pkn N -Nsin2p(k-1)n 2 N -Nsin2p(k+1)n 2 N Dan Ellis

  23. A[n] w w0 Fourier Transform of AM Sine • Sometimes we’d rather separate modulation and carrier:x[n] = A[n]cosw0n • A[n]varies on a different (slower) timescale • One approach: • chop x[n] into short sub-sequences • .. each with ~ constant modulation scale • DFT spectrum of pieces  show variation Dan Ellis

  24. FT of Short Segments • Break up x[n] into successive, shorter chunks of length NFT, then DFT each: Shows amplitude modulation of w0 energy Dan Ellis

  25. The Spectrogram • Plot successive DFTs in time-frequency: • This image is called the Spectrogram time hopsize (between successive frames)= 128 points Dan Ellis

  26. Short-Time Fourier Transform • Spectrogram = STFT magnitudeplotted on time-frequency plane • STFT is (DFT form): • intensity as a function of time & frequency window frequency index DFTkernel time index NFT points of x starting at n Dan Ellis

  27. w[n] n STFT Window Shape • w[n] provides ‘time localization’ of STFT • e.g. rectangularselects x[n], n0 ≤ n < n0+NW • But: resulting spectrum has same problems as windowing for FIR design: DTFT form of STFT spectrum of short-time window is convolved with (twisted) parent spectrum Dan Ellis

  28. STFT Window Shape • e.g. if x[n] is a pure sinusoid, • Hence, use tapered window for w[n] blurring (mainlobe)+ ghosting (sidelobes) e.g. Hamming sidelobes < -40 dB Dan Ellis

  29. STFT Window Length • Length of w[n] sets temporal resolution • Window length  1/(Mainlobe width) • moretime detail lessfrequency detail short window measuresonly local properties longer window averagesspectral character shorter window  more blurred spectrum Dan Ellis

  30. STFT Window Length • Can illustrate time-frequency tradeoffon the time-frequency plane: • Alternate tilingsof time-freq: • disks show ‘blurring’ • due to window length • area of disks is constant • Uncertainty principle: df·dt ≥ k half-length window  half as many DFT samples Dan Ellis

  31. Spectrograms of Real Sounds time-domain successive short DFTs individual t-f cells merge into continuous image Dan Ellis

  32. Narrowband vs. Wideband • Effect of varying window length: Dan Ellis

  33. Spectrogram in Matlab >> [d,sr]=wavread('sx419.wav'); >> Nw=256; >> specgram(d,Nw,sr) >> colorbar (hann) window length actual sampling rate(to label time axis) dB Dan Ellis

  34. STFT as a Filterbank • Consider one ‘row’ of STFT: where • Each STFT row is output of a filter (subsampled by the STFT hop size) convolutionwith complex IR just one freq. Dan Ellis

  35. STFT as a Filterbank • Ifthen • Each STFT row is the same bandpass response defined byW(ejw), frequency-shifted to a given DFT bin: shift-in-w A bank of identical, frequency-shifted bandpass filters: “filterbank” Dan Ellis

  36. STFT Analysis-Synthesis • IDFT of STFT frames can reconstruct (part of) original waveform • e.g. ifthen • Can shift by n0, combine, to get x[n]: • Could divide by w[n-n0] to recover x[n]... ^ Dan Ellis

  37. STFT Analysis-Synthesis • Dividing by small values of w[n] is bad • Prefer to overlap windows: i.e. sample X[k,n0]at n0 = r·H where H = N/2 (for example) • Then hopsize window length if Dan Ellis

  38. STFT Analysis-Synthesis • Hann or Hamming windowswith 50% overlap sum to constant • Can also modify individual frames of and reconstruct • complex, time-varying modifications • tapered overlap makes things OK Dan Ellis

  39. STFT Analysis-Synthesis • e.g. Noise reduction: STFT of original speech Speech corrupted by white noise Energy threshold mask M Dan Ellis

More Related