FOURIER ANALYSIS PART 2: Technicalities, FFT & system analysis

1. FOURIER ANALYSISPART 2: Technicalities, FFT & system analysis Maria Elena Angoletta AB/BDI DISP 2003, 27 February 2003

2. TOPICS 1. DFT windows 2. DFT resolution - improvement 3. Efficient DFT calculation: FFT 4. Hints on system spectral analysis M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 2 / 21

3. (1) Resolution: capability to distinguish different tones. Inversely proportional to main-lobe width. Wish: as high as possible. (2) Peak-sidelobe level: maximum response outside the main lobe. Determines if small signals are hidden by nearby stronger ones. Wish: as low as possible. (1) (3) (2) (3) Sidelobe roll-off: sidelobe decay per decade. Trade-off with (2). Several windows used (application-dependent): Hamming, Hanning, Blackman, Kaiser ... Rectangular window DFT – Window characteristics • Finite discrete sequence spectrum convoluted with rectangular window spectrum. • Leakage amount depends on chosen window & on how signal fits into the window. M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 3 / 21

4. In time it reduces end-points discontinuities. Sampled sequence Non windowed Windowed Some window functions DFT of main windows Windowing reduces leakage by minimising sidelobes magnitude. M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 4 / 21

5. Observed signal Window wish list Far & strong interfering components  high roll-off rate. Near & strong interfering components  small max sidelobe level. Accuracy measure of single tone  wide main-lobe DFT - Window choice Common windows characteristics NB: Strong DC component can shadow nearby small signals. Remove it! M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 5 / 21

6. Solution: sliding (overlapping) DFTs. DFT - Window loss remedial Smooth data-tapering windows cause information loss near edges. • Attenuated inputs get next window’s full gain & leakage reduced. • Usually 50% or 75% overlap (depends on main lobe width). Drawback: increased total processing time. M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 6 / 21

7. After padding bins @ frequencies NS = original samples, L = padded. Zero padding Improves DFT frequency inter-sampling spacing (“resolution”). M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 7 / 21

8. Frequency inter-sampling spacing: increased by zero-padding (DFT “frequency span” unchanged due to same sampling frequency) Capability to distinguish two closely-spaced frequencies: not improved by zero-padding!. • Zero-padding in frequency domain increases sampling rate in time domain. Note: it works only if sampling theorem satisfied! • Additional reason for zero-padding: to reach a power-of-two input samples number (see FFT). Apply zero-padding after windowing (if any)! Otherwise stuffed zeros will partially distort window function. NOTE Zero padding -2 DFT spectral resolution M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 8 / 21

9. We’re lucky here! ~ ~ SL = 20 Log10(|cr+kmax /ckmax|) bin, k kmax SL bin, k May impact on data interpretation (wrong f0!) DFT - scalloping loss (SL) Input frequency f0 btwn. bin centres causes magnitude loss Worst case when f0 falls exactly midway between 2 successive bins (|r|=½) |r|  ½ f0 = (kmax+ r) fS/N Frequency error: f = r fS/N, relative error: R=f /f0 = r/[(kmax+r)] R  1/(1+2 kmax) kmax f0 Note: Non-rectangular windows broaden DFT main lobe  SL less severe Correction depends on window used. M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 9 / 21

10. SL remedial • increasing N (?) • improve windowing, • zero-padding, • interpolation around kmax. DFT - SL Example DC bias correction, Rectang. window, zero padding, FFT DC bias correction, Hanning window, zero padding, FFT M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 10 / 21

11. DFT - parabolic interpolation Rectangular window Hanning window • Parabolic interpolation often enough to find position of peak (i.e. frequency). • Other algorithms available depending on data. M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 11 / 21

12. Incoherent averaging M = No. of DFT to average k = bin number, k=0, 1 .. N-1 Background noise fluctuations reduced, average noise power unchanged. Coherent averaging Background noise power is reduced. DFT averaging M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 12 / 21

13. DFT WNn,k = twiddle factors k = 0,1 .. N-1 Direct DFT calculation redundancy WNkn periodic function calculated many times. Direct DFT calculation requires ~N2 complex multiplications. complexity O(N2) Efficient DFT calculation: FFT VERY BAD ! Algorithms (= Fast Fourier Transform) developed to compute N-points DFT with ~ Nlog2N multiplications (complexity O(Nlog2N) ). M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 13 / 21

14. DFT  N2 FFT  N log2N DSPs & PLDs influenced algorithms design. ‘60s & ‘70s: multiplication counts was “quality factor”. Now: number of additions & memory access (s/w) and communication costs (h/w) also important. FFT advantages NB: Usually you don’t want to write an FFT algorithm, just to “borrow” it !!! Go “shopping” onto the web! M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 14 / 21

15. cost(sub-problems) + cost(mapping) < cost(original problem) Different algorithms balance costs differently. Example: Decimation-in-time algorithm time frequency Step 1: Time-domain decomposition. N-points signal  N, 1-point signals (interlace decomposition). Shuffled input data (bit-reversal). log2N stages. (*): only first decomposition shown. (*) FFT philosophy General philosophy (to be applied recursively):divide & conquer. Step 2: 1-point input spectra calculation. (Nothing to do!) Step 3: Frequency-domain synthesis. N spectra synthesised into one. M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 15 / 21

16. N : GCD(*)(N1,N2) = 1 N1, N2 co-prime. Ex: 240 = 16·3·5 N : GCD(N1,N2) <> 1 Ex: N = 2n Cost: SUB-PROBLEMS. Cost: MAPPING. • No twiddle-factors calculations. • Easier mapping (permutations). • Some algorithms: • Good-Thomas, Kolba, • Parks, Winograd. • Twiddle-factors calculations. • Easier sub-problems. • Some algorithms: • Cooley-Tukey, • Decimation-in-time / in-frequency • Radix-2, Radix-4, • Split radix. (*) GCD= Greatest Common Divisor FFT family tree Divide & conquer M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 16 / 21

17. Butterfly: basic FFT calculation element. Dual approach: data to be reordered in time or in frequency! • Decimation-in-timetime data shuffling. • Decimation-in-frequencyfrequency data shuffling. • In-place computation: no auxiliary storage needed, allowed by most algorithms. • DFT pruning: only few bins needed or different from zero  only they get calculated (ex: Goertzel algorithm). • Real-data case: Mirroring effect in DFT coeffs.  only half of them calculated. • N power-of-two: Many common FFT algorithms work with power-of-two number of inputs. When they are not  pad inputs with zeroes. (Some) FFT concepts & notes M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 17 / 21

18. Linear Time Invariant h[t] = impulse response x[n] h[n] X(f) H(f) Y(f) = X(f) · H(f) Transfer function can be estimated by Y(f) / X(f) Systems spectral analysis (hints) System analysis: measure input-output relationship. y[n] predicted from { x[n], h[t] } H(f) : LTI transfer function M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 18 / 21

19. Power Spectral Density of x[t] (FT of autocorrelation). Cross Power Spectrum of x[t] & y[t] (FT of cross-correlation). Transfer Function (ex: beam !) • Coherence function • values in [0,1] • assess degree of linear relationship between x[t] & y[t]. It is a check on H(f) validity! Estimating H(f) (hints) M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 19 / 21

20. References - 1 Papers • Tom, Dick and Mary discover the DFT, J. R. Deller Jr, IEEE Signal Processing Magazine, pg 36 - 50, April 1994. • On the use of windows for harmonic analysis with the Discrete Fourier Transform, F. J. Harris, IEEE Proceedings, Vol. 66, No 1, January 1978. • Some windows with a very good sidelobe behaviour, A. H. Nuttall, IEEE Trans. on acoustics, speech and signal processing, Vol ASSP-29, no. 1, February 1981. • Some novel windows and a concise tutorial comparison of windows families, N. C. Geckinli, D. Yavuz, IEEE Trans. on acoustics, speech and signal processing, Vol ASSP-26, no. 6, December 1978. • Study of the accuracy and computation time requirements of a FFT-based measurement of the frequency, amplitude and phase of betatron oscillations in LEP, H.J. Schmickler, LEP/BI/Note 87-10. • Causes et corrections des erreurs dans la mesure des caracteristiques des oscillations betatroniques obtenues a partir d’une transformation de Fourier, E. Asseo, CERN PS 85-9 (LEA). M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 20 / 21

21. References - 2 • Precise measurements of the betatron tune, R. Bartolini et al.,Particle Accel., 1996, vol. 55, pp 247-256. • How the FFT gained acceptance, J. W. Cooley, IEEE Signal Processing Magazine, January 1992. • A comparative analysis of FFT algorithms, A. Ganapathiraju et al., IEEE Trans.on Signal Processing, December 1997. Books • The Fourier Transform and its applications, R. N. Bracewell, McGraw-Hill, 1986. • A History of scientific computing, edited by S. G. Nash, ACM Press, 1990. • Introduction to Fourier analysis, N. Morrison, John Wiley & Sons, 1994. • The DFT: An owner’s manual for the Discrete Fourier Transform, W. L. Briggs, SIAM, 1995. • The FFT: Fundamentals and concepts, R. W. Ramirez, Prentice Hall, 1985. M. E. Angoletta - DISP2003 - Fourier analysis - Part 2.2: Technicalities 21 / 21