No Spectral Leakage Time Domain Frequency Domain Signal
- Slides: 55
No Spectral Leakage Time Domain Frequency Domain • Signal: sin(2π * 3 fs) • Three periods per sampling time period – When f=3, DFT correlation squares each value; result >0 – Result is N/2 times the time domain amplitude • The four periods per time period – Positives and negatives cancel each other – FFT spectral value for f=4 totals zero
Spectral Leakage Time Domain Frequency Domain • Signal: sin(2π * 3. 5 fs) • Three periods per sampling time period – When f=3, DFT correlation attenuates the total – Some of the amplitude spills across the spectrum
DFT Window-sync Distribution • Window sync function: sin(x)/x – Period is the same as sin(x) – Amplitude degrades as x increases – Sin(π(k-m))/{π(k-m)}; k is signal frequency, m is DFT basis function • Integer k points compute to zero because sin(πk) = 0 • Non integer k points ≠ 0 because sin(π(k-m)) is not a multiple of π
Spectral Leakage Examples
Leakage Symmetry 3. 4 cycles period 64 point DFT Circular view of the frequency spectrum
Leakage Symmetry (cont. ) Linear view of the frequency spectrum
More Examples
Scalloping Loss • • • Leakage means power loss in non-center frequencies Loss amount greatest for frequencies at the intersections Loss can be significant (0. 637 in this example) Result is the scalloped pattern of frequency amplitudes Solutions: – Zero DFT creates finer frequency resolution, and therefor causes the scalloping low points to be closer together – Applying a windowing function alters the scalloping curve shape
Windowing • Issue: In practice, we apply DFT to finite signal pieces (frames) • Impact: – The start of the truncated signal does not typically line up with the end (rough edge) – The frequency domain experiences distortion • Using a window – First Create the window: Create an array of frame length with values calculated with a windowing formula – Apply the window: multiply these values by each time domain frame amplitude
Speech Application Speech signal characteristics change as we speak • Goal: Extract short-term signal features from speech signal frames • Approach – – – Block the speech into short overlapping segments Duration of each block varies from 10 to 30 ms Assumption: Speech within a single frame is unchanging Apply a window function to smooth the signal Apply DFT after applying the window Perform frequency analysis to extract speech characteristics
Window Function Application of windows to sinusoid signal 1. rectangular (no) window 2. Trangular window 3. Hanning window
Popular Window Types • Rectangular: wk = 1 where k = 0 … M – Advantage: Easy to calculate, array elements unchanged – Disadvantage: Messes up the frequency domain • Ideal (window-sync): sin( 2 π f i) / (πi) – Disadvantage: Must be infinitely long – Truncation causes ripples and overshoots • Hamming: wk = 0. 54 – 0. 46 cos(2 kπ/M) – Advantage: Fast roll-off in frequency domain – Disadvantage: worse attenuation • Blackman: wk = 0. 42 – 0. 5 cos(2 kπ/M) + 0. 08 cos(4 kπ/M) – Advantage: better attenuation – Disadvantage: slower roll-off Multiply the window, point by point, to the audio signal
Other Windowing Formulae • Hanning: w[n] = 0. 5 -0. 5 cos(2πn/(N-1)) – Advantage: Adding 50% overlapped frames produces no gain – Advantage: Excellent for pitch manipulation algorithms • Bartlett: w[n] = 2/(N-1)/2 - |(n–(N-1)/2|) • Triangle: w[n] = 2/N. (N/2 - |(n–(N-1)/2|) • Blackman: w[n] = a 0 – a 1 cos(2πn/(N-1)) – a 2 cos(4πn/(N-1)) where: a 0 = (1 -α)/2 ; a 1 = ½ ; a 2 = α /2; α =0. 16 • Blackman-Harris: w[n] = a 0–a 1 cos(2πn/(N-1))–a 2 cos(4πn/(N-1)) - a 3 cos(6πn/(N-1)) where: a 0=0. 35875; a 1=0. 48829; a 2=0. 14128; a 3=0. 01168 Note: Small changes can cause large frequency domain impact
Code: Hamming Window Create Window double[] window = new double[window. Size]; double c = 2*Math. PI / (window. Size - 1); for (int h=0; h<window. Size; h++) window[h] = 0. 54 - 0. 46*Math. cos(c*h); Apply window for(int i=0; i<window. length; i++) frame[i]=frame[i]*window[i];
Rectangular Window Frequency Response Time Domain Filter Top Right: Amplitude linear scale, Bottom Right: Power logarithmic scale
Understanding Digital Signal Processing, Third Edition, Richard Lyons (0 -13 -261480 -4) © Pearson Education, 2011.
Windowing Frequency Response • Main Lobe: Narrow implies better frequency resolution. As the window length grows, the main lobe width narrows (less initial spectral leakage). • Side lobe: Higher for abrupt time domain window transitions from 1 to 0. • Roll-off rate: Slower for abrupt window discontinuities. Side Lobe Roll-off Rate Main Lobe
Rectangular Window • Narrow main lobe width: 4π/M, Side lobe width: 2 π /M • Spectral leakage directed to places that hurt analysis • Poor roll off and high initial lobe
Convoluting Rectangular Window with a particular frequency Note: Windows with more points narrow the central lobe
Triangular (Bartlett Window)
Compare: Hamming to Hanning Hamming Hanning: Faster roll of; Hamming: lower first lobe
Blackman & Hamming Frequency Response
Comparing Windows
Spectral Leakage using Hanning Window Understanding Digital Signal Processing, Third Edition, Richard Lyons (0 -13 -261480 -4) © Pearson Education, 2011.
Spectral Leakage using Hamming Window Understanding Digital Signal Processing, Third Edition, Richard Lyons (0 -13 -261480 -4) © Pearson Education, 2011.
Efficiency Issues • Calculation of Sin(x), Cos(x), ex – Sin(x) = x – x 3/3! + x 5/5! – x 7/7! … – Cos(x) = x – x 2/2! + x 4/4! – x 6/6! … – ex = x + x 2/2! + x 3/3! + x 3/4! … • Suggestions – Use table look up – x 3 = x * x – Consider the use of the cache when using loops
Fourier Analysis Real-time DSP Algorithms must be as efficient as possible • The Naïve DFT requires a double loop – Outer loop is T long (length of the time domain frame) – Inner loop is N long (number of FFT frequency bins) – Total calculations: O(N*T), or O(N 2) when N = T • Fast Fourier Transform (FFT) – N must be an even power of 2, if not, pad with zeros – Total calculations O(N lg. N) calculations • Compare efficiencies N N 2 lg N N lg N Improvement 128 16, 384 7 896 18. 2857 1, 024 1, 048, 576 10 10, 240 102. 4000 8, 192 67, 108, 864 13 106, 496 630. 1538
Understanding Digital Signal Processing, Third Edition, Richard Lyons (0 -13 -261480 -4) © Pearson Education, 2011.
Optimization of DFT 1. Replace the double loop by a recursive approach 2. Eliminate the recursion for further efficiency 3. Perform obscure optimizations to extract even more speed How do we know if the optimized version works? 1. Observe the frequency domain results a. b. c. d. Value at index N/2 -1 = value at index N/2+1 Value at index N/2 -k = value at index N/2+k The 0 index contains a constant multiplier that when removed causes the signal’s mean to be zero Value at N/2 is unrelated to the rest of the spectrum 2. Perform an inverse FFT to restore the original signal 3. Compare the results to the slower DFT algorithm 4. Apply algorithm to simple sinusoid signals and check results
Why is FFT fast • Evaluate polynomial at n roots of unity • Why is this important? – Dealing with exponents is easy • eaeb=ea+b, (ea)b = eab – The ith root of unity ωin = e 2πi/n. – For even numbers, wnn/2=w 21=-1 – For even numbers, squares of the n complex roots of unity are the same as the n/2 complex roots of unity – Summing all n roots of unity equals zero – The calculation of e 4πi/n=e 2 π i/n because of the periodicity • We can avoid many calculations, and use a divide and conquer recursive solution – Solve smaller parts of the linear equations and easily combine results to solve the larger equation
The Fast Fourier Transform • • The FFT is an efficient divide and conquer algorithm For polynomials, If n is even, we can divide into two smaller polynomials Peven(x) = a 0 + a 2 x 2 + a 4 x 4 + … + an-2 xn-2 Podd(x) = a 1 + a 3 x 2 + … + an-1 xn-1 and we can write the equation as follows P(x) = peven(x) + x * podd(x) • This suggests recursion of lg(n) steps. The above equation is the recursive relationship between levels
Theory for Optimization Base Case: x[0] Recursive Relationship (2 t)/N = t/(N/2) Fk = ∑t=0 N-1 x[t] e-j 2πkt/N = ∑t=0 N/2 -1 x[2 t] e-j 2πk(2 t)/N + ∑ -j 2πk(2 t+1)/N t=0 N/2 -1 x[2 t+1] e = ∑t=0 N/2 -1 x[2 t] e-j 2πkt/(N/2) + e-j 2πk/N∑ -j 2πkt/(N/2) x[2 t+1] e t=0 N/2 -1 = Fkeven + e-j 2πk/N * Fkodd
The FFT Algorithm The recursive step has O(N) computations The recursion requires O(lg(N)) levels Total calculations are O(N lg(N))
FFT Recursive Efficiency Figure to the right – FFT of 16 points – O(N) work at each level – lg N + 1 = 5 levels – Total calculations: O(N lg(N)
The FFT Algorithm Overview • Initial Recursive Call (FFT(time) time contains signal amplitudes, and is a power of 2 length (N) • Base Case: n=1, then return the time array • Divide Step separate time into even (timeeven) and odd (timeodd) sub-arrays • Recursive Step: yeven = FFT(timeeven), yodd= FFT(timeodd) • Combine Step For loop of complex number multiplications
Simple Recursive Solution public static Complex[] FFT(Complex[] time) { int N = time. length; Complex[] fourier = new Complex[N]; if (N==1) { fourier[0] = time[0]; return fourier; } Complex[] even = new Complex[N/2], odd = new Complex[N/2]; for (int m=0; m<N/2; m++) { even[m] = time[2*m]; odd[m] = time[2*m+1]; } } Complex[] q = FFT(even), r = FFT(odd); for (int k=0; k<N/2; k++) { double exp = -2*k* Math. PI /N; Complex wk = new Complex(Math. cos(exp), Math. sin(exp)); fourier[k] = q[k]. plus(wk. times(r[k])); fourier[k+N/2] = q[k]. minus(wk. times(r[k])); } return fourier; -2 kπ/N+N/2 Note: e = -e
Inefficiencies Computations still are an order of magnitude slower than necessary • Repetitive calculations of sines and cosines are extremely slow • Recursive activation record overhead is huge • Declaring and copying arrays at every step slows things down at least by half • N>>1 is ten times faster than N/2 • The Complex class methods cause machine code jumps that put pressure on the hardware cache
Eliminating the Recursion Butterfly algorithm 000 000 001 010 100 011 100 101 110 111 101 111 100 110 001 010 110 001 101 011 111 • The numbers in the rectangles are the array indices • Notice how the indices change during the recursion • Can you see a pattern ?
Bit flipping implementation • Determine if a number is an even power of two: n & (n-1) • Determine the upper most bit number, B int B=1, x = n-1; while (x&(x-1)) B++; Result = 0 FOR each bit b from 0 to B IF bit set, shift B – b and OR into result
Butterfly Code Flip bits from left to right int j = N>>1, k; for (int i=1; i<N-1; i++) { if (i < j) { swap (x[i], x[j]); } k = N>>1; while (k>=2 & j>=k) { j -= k; k >>= 1; } j += k; } • Most Significant Bit Swap. Bit ( x, x + lg. N) • Second most significant bit Swap. Bit(x, x + lg(N/2) • Third most significant bit Swap. Bit(x, x + lg(N/4) • kth most significant bit Swap. Bit(x, x + lg(N/2 k))
Sin and Cosine Table Look Up Compute the values ahead of time and save repetitive calculations • ei 2πk/N = cos(2πk/N) + i sin(2πk/N) • We can store in an array (sin. X[]) sin(2π0/N), sin(2π1/N), sin(2π2/N) sin(2π3/N), … sin(2π(N-1)/N) • cos(2πk/N) = cos. X[(k+(N>>2))%N]
Wk. N= e 2πk/N This relationship allows table lookups for precomputed sines and cosines Understanding Digital Signal Processing, Third Edition, Richard Lyons (0 -13 -261480 -4) © Pearson Education, 2011.
Optimized FFT // Perform the fft calculations. for (int stage=1; stage<=M; stage++) // M = lg N { // Remember that complex numbers require pairs of doubles fft. Sub. Group. Gap = 2<<stage; // 4, 8, 16, . . . – subgroup distance gap = fft. Sub. Group. Gap>>1; // 2, 4, 8, . . . – odd/even distance k. Inc = N>>1; // Number of 2 PIki/N steps for odd/even entries. // Outer loop: each sub-fft group; inner loop: combine group elements for (int even=0; even<complex. length; even+=fft. Sub. Group. Gap) { k = 0; // Index into the trigonometric lookup table. for (int element=even; element<(even+gap); element+=2) { // ***** See Next Slide ***** k += k. Inc; } } k. Inc >>= 1; } // position for next look up.
Multiplication Portion // Look up e^2 PIki/N avoiding trig calculations here. real. W = sines[(k+(N>>2))%N]; // cos(2 PIk/N); imag. W = -sines[k%N]; // -sin(2 PIk/N); // Complex multiplication of the odd entry of the subgroup // with (e^2 PIi/N)^k = (cos(2 PI/N) - i * sin(2*PI/N)^k j = (element + gap); temp. Real = real. W * complex[j] - imag. W * complex[j+1]; temp. Imag = real. W * complex[j+1] + imag. W * complex[j]; // Adjust the odd entry (subtract: the fft is periodic). complex[j] = complex[element] - temp. Real; complex[j+1] = complex[element+1] -temp. Imag; //Adjust the even entry. complex[element] += temp. Real; complex[element+1] += temp. Imag;
Final Notes • The above analysis – Assumes that N is a power of 2 – If not, pad the array with zeroes – For speech recognition, we always use power of 2 frames • What if not a power of 2 – The same analysis is possible with other radices – When hit a prime, fall back to the slower O(n 2) algorithm
Final Notes • Standard Fast Fourier Transform – requires N to be a power of 2 for recursion to work – Can pad the array with zeroes to extend frequency domain • Can it work if N is not a power of 2? – Yes, but special slower processing is needed • How do we know if it works? – Point N/2 -1 = Point N/2+1, Point N/2 -2 = Point N/2+2, Point N/2 -k = Point N/2 + k, etc. – Note: Points 0 and N/2 don't match, so don’t check these – The FFT Inverse should restore the time domain signal – Compare to the slower correlation DFT calculation – Try some simple impulses and check the results
- Spectral leakage
- Spectral leakage
- Frquency domain
- "gulf coast college" -casino -events -jobs
- Z domain to frequency domain
- Z transform time reversal
- Z transform integrator
- Impulse erp
- Laser synchronization
- Time frequency
- Baseband signal and bandpass signal
- Baseband signal and bandpass signal
- Digital signal as a composite analog signal
- The product of two odd signals is
- Row conditional relative frequency
- Probability with relative frequency
- Constant phase difference
- Vmax shm
- Frequency vs relative frequency
- Marginal frequency distribution
- Joint relative frequency vs conditional relative frequency
- Physiological effect of electricity
- Identify leakages on claims
- Subthreshold leakage current equation
- Foreign trade multiplier formula
- Short note on elcb
- Iec 60601 leakage current limits
- Microshock and macroshock
- Eli p llc clean air
- Megger dcm305e earth leakage clamp meter
- Cmos leakage current
- Leakage current and temperature relationship
- The air leakage rate for a combination vehicle
- Gate-induced drain leakage
- Touch current vs leakage current
- Gidl dibl
- Craneware.com
- Symptoms of preterm labor
- Pericatheter leak meaning
- Leakage point
- Eec micro project pdf
- Front wheel brakes are good under all conditions
- Nfpa 99 leakage current limits
- Spectral regrowth
- Spectral regrowth
- Spectral classification
- Profil spectral rigel
- Spectral normalization
- Spectral hashing
- Spectral graph theory and its applications
- Spectral efficiency
- Séquence principale
- Spectral bands
- Eric xing
- Spectral clustering
- Foster and freeman vsc 80 price