113D_1_EE113D_CR8_FFT - (\/ University of Hertfordshire...

Info iconThis preview shows page 1. Sign up to view the full content.

View Full Document Right Arrow Icon
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview: (\/ University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMA IN LECTURE 4 THE FREQ UENCY DOMAIN OBJECTIVES This lecture should achieve the following: • • • • • • Introduce phasor representation of sinuso ids Introduce Fourier Series and Transforms Introduce Discrete Fourier Transform (DF!) Consider operational complexity ofDFTs Deduce a Radix-2 FFT algorithm Consider some implementation issues of FFTs with DSPs LECTURE 4 4- 1 167 University of Hertfords hire STUDENT'S GUIDE THE FREQUENCY DOMAIN USE OF FREQ UENCY DO MA IN Frequency Respons e ) • SIGNAL FREQUENC IES> IAlt L r r. WILL BE ATTENUATE D A ND DISTORTED SPEECH ANALYSIS • RECOGNITION • IDENTIFICATION . r IAI _ _ IA N' S 'A ' • •• • •• JOHN'S 'A ' • SPECTRUM ANALVZER S CONVERT A TIME DOMAIN SIGNAL INTO FREQUENCY DOMAIN • Frequency Domain In our previous lectures, we have made quite extensive use of frequency domain inform ation about signals. For example, when demonstrating the effect of filtering, we examined the frequency content of signals. Similarly, in deciding the quality of filter performance, we conside red the frequency response of filt ers . In general, such frequency domain information is extremely important and useful in signal processing. Let us consider two more applications where frequency domain information is useful. When we want to transmit signals along a telecommunication channel, it wou ld be usefu l if we knew the frequency resp o nse of the channel. On this information, we can decide what group of frequencies can be sent through the channel without much distortion. . Although our transparency only shows the amplitude information versus frequency, phase information is also necessary to ensure that the sign al will travel through the channel without degradation of signal quality. . Some companies have already implemented security systems that rely on spe aker verification. With speech analysis, sounds are analyzed in orde r to recognize a command. Analysis of the frequency content of each sound plays an important part in its identification. For example, the lower figure in our transparency shows the frequency content of the letter "A" as spo ken by two different people. We can use this information to make a judgment about the letter spoken, or decide by whom it was spoken. In both cases, the frequency content of the sign al provides vital informati on. Television sets that accept voice commands now exist. However, these sets require software which is speaker-independent in order to recognize a word no matter by whom it was spoken. Again , all of this complex analysis must be performed in the frequency domain. LECTURE 4 4-2 168 University of Hertfordshire STUDENrs GUIDE THE FREQUENCY DOMAIN T H E PHASO R MODEL CO MP L EX PL ANE 1m 1m • Imaginary Re . Re. 1 PHASO R - VE C T O R RO TAT ING Speed · co rad per . e c o n d ----~--:.---;-. . Re A m pli tu d e - A 2. PoluFonn 1. Re c tangular Form x(t) " A . I(eDt) where xltl- a + Jb where and A.~ +- rot - tan .t!.a e j(eDt) " COS(eDt) + j s ln (eDt) ro - 2", 1t - 180 degrees Before we can describe the process of converting a time domain representation of a signal into a frequency domain representation, we need some way of modeling the original signal. The method used to describe the signal is called the p hasor model. A phasor is a vector rotating on the complex plane with a magnitude A and a rotational speed of co radians/sec. At any point in time, the signal can be represented by x(t) = a + jb where a is the "real" value of the signal, and b is the complex value of the signal. This is the rectangular representation of the signal. The signal can also be represented as a polar version (an amp litude and a direction) by: x(t) = Aei(rot) where: ei<rot) = cos(cot) + j sin(rot) and: co(t) = ~ = tan-I b/a Using the polar notation x(t) = Aei(lIl l) is more useful than the rectangular version of the equation , as it introduces the variable time I, which allows us to determine the position of the phasor at any point in time with respect to its frequency. Remember that co = 21t fwhere jis the frequency of the signal. Remember that a DSP operates in the time domain, and thus we need some means of converting between the time and frequency domain and vice versa. LECTURE 4 4-3 169 . ,1 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN DISCRETE PHA SOR MOD E L 3. Dls c rele Time Inlllil condi llo n Ihen g ene ral form X(O) x( n) =A e J(no>T.) a xlt) = o ffs el =Ae Ja T•• Sa mp ling Peri od =Ae J(o>l+ a) conllnuou s x (n ) = Ae J(n ",T.+a ) IL..? -' W HA T HAV E WE AC H IEV ED? • INTR ODU CED A PHAS OR MODE L TO ·RE PR ESE NT DISCRETE AND CONTINUOUS SI GNALS • PHASORS . • HAVE AN AMPLITU DE PROPORTIO NA L TO THE SIGNAL AMPLITUDE • ROTATE A T A SPEE D PROPORTION A L TO THE SIGN AL fREQUENCY CAN WE REPRESENT ANY SIGNAL WIT H THIS MODEL? In the last transparency , we showed how the equation could be used to represent the continuous change of the signal with respect to time. While DSPs do not cope with continuously varying signals, they will cope with discrete "snapshots" of a signal. Therefore, we need to express the last equation in a manner that takes this into account First, we introduce two new variables: Ts ' which is our sampling period, and n, which is the number of samples recorded. For example, if we need to represent the signal after IOms, and we have a sampling period of2ms, then T, = 2ms and n = 5. So far, we have assumed that the signal starts with zero phase shifts (i.e., it starts flat along the x-axis) . In reality, not all signals have the same phase, so we have to introduce a phase offset a.. This gives us: Ifwe replace the continuous time variable with discrete time variables, we get: x(n) = Ad(ncoTs+a) We now have a method of representing signals both in a discrete and continuous manner. To reiterate, a p hasor has the following properties: • • An amplitude proportiona l to the signal amplitude A rotation speed proportional to the signal frequency We will now show how any signal can be represented. LECTURE 4 4-4 170 University of Hertfordshire STUDE Nrs GUIDE THE FREQUENCY DOMAIN MODELING SINUSOIDS LET US REWRITE e j lJ>t AS e j(+ ) = cos(. ) + js in(+) AND e· j(+ ) = cos(. ) • jsi n(+ ) LET += (ro t + (X) cos. " OR (nroT. + (X) ---.. THEN e j(+ )+ e • j(+ ) AND 2 DRA WI NG THE PHASORS FOR COS • 1 m IN GENERAL x(l ) =R cos(1Il1+ a ) R j (lJ>1 + a ) (8 +8 OR AS SUM OF TWO PHASORS • j(1Il + a ) 1 ) X(I) = 2 Rf2 We will now demonstrate how a sine or cosine wave can be modeling with the equations we have derived. The main steps are shown in the above diagram, but for clarity, all of these steps are also shown below. From Transparency 2, we know that: ci<CO') = cos(rot) + j sin(rot) and 41 = (rot) or (rot+a) or (ncrT, +a) From this we gain: (I) (2) eK+> = cos(cjI) + j sin( cjI) e-:i<+> = cos(cjI) - j sin( 41) and also by default, Rearranging (l) gives us: (3) cos(cjI) = ci<+> - j sin(41) Rearranging (2) gives us: (4) j sin(cjI) = cos(cjI) - e- < +) j Substituting(4) into (3) gives us: cos(cjI) = ci<+> - (cos(cjI) - e-j<+~ LECTURE 4 4-5 171 u~ 'IJ' University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN Rearranging again gives: 2c o s( cI» = ei(+) + e -j(~) or (5) eM) + e- M) cos( ¢) = - - 2 - Substituting (5) into (4) gives: . . 'f' e +e Jsm("') = - - - - 2 Which fmally gives us: j( ~ - .j(;) (6) sm(¢) = . (e 2j (;») -ei (; ) i So what was the point of defming sin( cI» and cos(cI» ? Suppose we want to model a cosine signal x(t)? This can be defmed as: x(t) = R cos(rot+a) Using (5) we can change this to: x(t) = R/2 (e j t + e -it) Remember that e j t is a phasor, so our cosine signal consists of two phasors. These are complementary phasors of half the total amplitude of the "real" value of the cosine signal, and are rotating in opposite directions at the same velocity. The two imaginary parts of the signal cancel each other out, leaving the real cosine signal. From this we can see that any real sinusoidal signal can be represe nted by a pair of equal but opposite phasors. These are known as a conjugate pair. LECTURE 4 4-6 172 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN ORIGINS OF FOURI ER SERIES As lnCsI) e ,FO U R SIN USOIDS e ADD TH E M TO GIVE A N EW SIGNA L A s ln l b l) As lnlel) Asl n (dl) e CAN THIS B E GENERALIZED? ___ 1/4 Period (e. ))-----'~ A sln(sl) Is Ihe fundsmenlsl frequency s ln C l) , sln(cl) s n d sln(dl) sre multiples of II b 00 s( l) • LA. J( . ... I) k- ·CID fal , .s the Fundamenta. Frequency k Is s multiple (harmonic) of Ihl s frequency We now know that a sinusoidal wav e can be represented by two conjugate phasors, but what about more complex signals? The top graph in the diagram shows four sinusoidal signals of equal amplitude but different frequencies. If we add these signals together, we get the resulting waveform in the lower half of the diagram . We can express this mathematically as: x(t) = Asin(at) + Asin(bt) + Asin(ct) + As in(dt) What is really required is a way of expressing the relationship in a more general form. If we 'look at the top graph , we can see that signal a has the lowest frequency. We will call it the "fundamental" frequency. All of th e other frequencies in the diagram are of a 'higher frequency than a and are 'Called "multiples", or more correctly, harmonics of the fundamental frequency a. To generalize this, we define the fundamental frequen cy as COo and all of the harmonics as kco o> where k is an integer representing the th k harmonic. If we use the above example, a = coo' b = 2coo, C = 4000' x( t) = Ae .Kcvol ) + A e .K 2cvot ) + A e .i< 4cvol ) + Ae j (S4/ot) and d = 8co o If we use phasors in the expression, we get: Instead of just representing four phasors, we need a more general method of expressing thi s ser ies. This can be expressed in the general form of: co x( t) = ~ Ae j( k4/ol ) k =- «> LECTURE 4 4-7 173 u~ 'IJ' University of Hertfordshire STUDENrs GUIDE THE FREQUENCY DOMAIN Our example only had four frequencies , yet this equation shows the summation of an infinite number of frequencies. In reality, our signal also consists of an infmite number of frequencies, but all except four have zero amplitude. So what have we learned? We now have an expression which shows us that a periodic signal can be constructed from the summation of a series of sine waves. Unfortunately , this statement only holds true under certain circumstances, as we will now discover. LECTURE 4 4-8 174 University of Heitfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN FOURIER SERIES , ....".,. .... A Amplitude C . o f a algna' can change with reapect to time. M a them a tic a l Re p re s en ta tio n T p/2 'J\./\} f\ !\ V • • c_....!- J e - /I _", II fX(t) T p dt .T p/2 T p= Re petiti o n Peri o d Tlme_ T p = 21< Illl O The Fourier S e r ie s ll(I)-L C_ell - m, l l . Sum of an Infinite num ber of sinusoid III lind coslnusoldlll terms Elich term has a time dependent coeffi c ie n t r ep rese nte d by C . k- ·ClD The previous equati on ... ao x(t) = LAej(kCl>ol) k =-ao ... has one major flaw. It assum es that all of the composite sinusoidal signals have a constant amplitude. The top half of the transparency shows that we can have a sinusoidal wave with an amplitude that varies over time. Clearly, we need some way to represent this. If we express the amplitude of the phasor with respect to time, we end up with the follo wing equation: Ck =; V2 !x(t)e- j(klllol)dt p - Tp l2 This is a way of expressing that the amplitude of a phasor changes with respect to time. The range (-T f2 to Tfl>"gives us a full range of the fundamental ·time period Tp (also known as the repetition period) as well as the negative values required for mathematical completeness. In other words, the integral gives a constant for a particular freque ncy at a given time. If we substitute this new comp lex constant Ck in our original equation, we end up with: x(t ) = L Ckej(kCl>ol) k =- ao ao This equation is known as the Fourier Series, and represents the sum of an infmite number of sinusoidal and cosinusoidal terms. Each of these has a complex constant (CJ representing its amplitude at anyone point in time. With this equation, we can describe any periodic signal. We will next describe a non-periodic signal. LECTURE 4 4-9 175 , , University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN FOURIER TRANSFO RM 0) Tpl2 C ke x( l) =L J( k CIl, I) k=- WHERE Ck 0) =-'- J' e -II k ". II dt fX(t) Tp -Tpl2 • INCREA SE Tp .. PER IOD INCREA SES T_ 0) p ~ NO REPETITION 1 - T =~ P 211 ~ 211 ID • DISCRETE FREQ UENCY VARIABLE BECOMES CONTINUOUS k ID o _ C(Cil ) ~ • DISCRETE COEFFICIENTS CIO c, BECOME CONTINUOUS _ : d C (CJ» =---!!!. 211 - jx(t) e - j Cll I dt normalize~ ' 0) ~ dCilI 211 • f =:X(ID) .. x(t)e -j01 dt . 0) FT PAIR eo : 'THE INVERSE' - : x(t) =2;" 1 JX( ID ) e j01 d Cll~ r, . On the last diagram, we defmed an equation to repres ent periodic signals as: (1) x( t) = LC GO Tpl2 k e j( k4lol ) and (2) C Ic = ; J x(t)e -j (kCll o ')dt Ic. =- co p -Tp l2 In reality, we need equations to represent non-periodic waveform s such as speech. How can we achieve this ? Think of the example of a single word being spoken. This is as a non-periodi c waveform because it doesn 't repeat itself. However, it could also be viewed as a periodic waveform with an infinite period. If we assume that Tp tends towards infinity, then we can produce equations for non-periodic signals. If Tp tends towards infinity, then IDo tends towards O Because of this, we can replace roo with d ID (which . represents a very small value of to). We know that e =21tfand f = Iff, so we can rearrange these to give us: - 1 Tp =- d eo 2,. Because our fundamental frequency ro o tends towards 0, kroo also tends towards O. Therefore, we mus t remove the harmonic des ignator k from the equations, and C k becomes C(ro), giving us a continuous value. With these changes, equation (2) changes to: dro +<Xl C(ro) =x(t )e - j<Jl1dt 21t J . LECTURE 4 4- 10 176 u$ University of " " Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN If we normalize the equation by dividing both sides by doo /21t , we arrive at: We redefm~ 21t C( oo ) as X(oo ) which gives us: doo +00 (3) X( oo ) = J x(t)e - j<n'dt We now have a function which can tell us the value ofa frequency (00) at any point in tim e. Substituting X( oo ) in equati on (1) and changing the summation to an integral gives: (4) 1 l OJ . x(t)=- X(oo)e )(l)'doo 2n Since the integration is over a frequency range, the integration variable is d(f) = d(oo ) / 21t. We now have one equation (3) to convert from time domain to frequency domain, and another equation (4) to convert from frequency domain to time domain . Together, these two equations are commonly known as the Fourier Pair. LECTURE 4 4- 11 177 . University of Hertfordshire STUDENT' S GUIDE TH E FREQUENCY DOMAIN , DISCRETE FOUR I ER SERIE S ., • FOR DISCRE TE FO U RI E R SERI ES REPLACE t T• • sam pIIng perio d = nr , x(n ). L e ke ' (k "" T' " I • B UT k CO.T.= k ~27t FO R k ~.O . 1.2 • •" • THE SERIES REPEATS ITS EL F AROUN D THE SAMPLING FR EQU ENCY (AN ASSUMPTION WE MA DE IN THE SA M P LI N G LECTURE ) For the Fourier Series to be of use on a DSP, vie must convert it into a discrete version so that it can be applied to a sampled analog input. This conversion is simply performed by replac ing the continuous time variablet with the two variables nTs , where T, is the sampling period and n is the number of samples taken. This gives us the following revised equation: x(n ) = L Ckej(ka>oT,n) k=-co co Ifwe define phase jump (the difference such that the phase is the same) for the kth harmonic component as: kco oT, = 21tDl .. .where m is an integer (such that m = 0, 1,2, ... ) we obtain kcooTs = 0, 21t, 41t, .. . We know that f =..!T to and f = ~ 21r alo _ Therefore, T=21r which means kwoTs = k-21r - 0, 21t, 41t, ... Ws Dividing by 21tgives: kalo _ al s - 0, 1,2 ,... This shows us that·the frequency response of a discrete signal has a period of T s ' In other words, it repeats itself around the sampling frequency. This verifies the assumptions made in the sampling lecture. LECTURE 4 4-12 178 University of Hertfordshire STU DENT'S GUIDE TH E FREQUENCY DOMAIN D IS C R ET E FOURI E R TRAN S FORM FT PAI R FREQUENCY EQUATION <0 X ( CO ) · JIl(I).-J"'1 dl <0 - - --l~ .. X(CO). '" L Il(n). -j" T. n • • • REPLACE I · T. n CONTINUOUS x(l) BECOMES DISCRETE xlnl DO NO T INTEG R A T E BUT SU M AL L DISC RETE SAMPLES FT PA IR TIME EQUA Tl O N <0 X(t)= _1_ 21l JX (CO ) • J m I _~ dm - - - - - i~ .. Il x (n) = _1_ 21l IX (CO ) -n ~ • J (.. T . n) d ( mT.) • • LIM ITS OF IN TEGRA TION D O NOT NE ED TO GO BEYOND REPEA T S IT S EL F OUTSIDE ~ Il fT. (EVERY 2 Il I KEEP INTEG RA TIO N A S X(m liS CONT IN UOU S Il AS TH E SPEC TR U M Once again, we replace the continuous time variable I with the two variables n'T, in the first equation of the Fourier Transform pair, where T, is the sampling period and n is the number of samples taken. Since we no longer have the continuous variable I , x(t) is now referenced as x(n). This gives us the following revised discrete frequency transform: co x( m) = L x(n) e- j aif•D 0= -00 Note that the integration has changed to a summation. This is because we are now using discrete rather than continuous values. We can also change the time-based equation in a similar manner: 1 II x(n) = - fX(co)ej(CI>T.D)d(COT) 2x s - Il In this case, the integration stays because X{co) is a continuous function. How ever, the limits are now reduced to -x to x as we have established that the frequency spectrum repeats itself at intervals of 2x. At this point, we now have two equations which could be implemented on a DSP, which allows us to convert from time domain to frequenc y doma in and back again. LECTURE 4 4-13 179 . . Un iversity of Hertfordshi re STUDENT'S GUIDE THE FREQUENCY DOMA IN SI M P L E R DF T X ( OO ) _ 2: xln). 'J" Tsn n . • INFIN ITE NUMBER OF DISCRETE VA LUES Wil l TAK E A L O N G TIM E TO COMP U TE :~----~--lI"';: Xlk) • L ET US ASSUME WE TAKE ON E PERI OD OF THE LONG SEQUENCE PERIOD - NT, NOW FREQU ENCY "',I N N 6 ~ .. s = 6 =- N lOS . N, ' X(kll ) - I xln). -JkI>Tsn _ . SET k 6 · .. ANDSUMONLYNSAMPLE S 2" 2 lt n-O . X N( k ) -I N·' oj 2ltkn xln). N _ • SINCE I> Ts ~ - N - - ; ; - - N;0 • • • •• •• : ,~ lO s z•••• •••N1 ··· ···· ..: ·, SE T W N -.- ·N ····· ·· ····: :; .kn : J 2" N THEN : : :XN (k ) -L.J x n n- O ~ I ) W kn • N : AND: xln) - - L.. XN(k )W N :: : ... . . ........ . . . . .. ... . Nk_O .: Th e DFT we have previously defmed has two major implementation problems. First, a summation of an infinite number of inputs is not possible in the real world, so we need some way of reduci ng the number of inputs required. Second, we need to reduce the number of frequencies that we are going to compute, as we only have a fmite am ount oftime to perform the calculations. Th e first problem is easily solved by recording a section of the input samples (th is is known as windowing). We kn ow that frequency spectrum is repetitive with respect to oos' so we can deduce how many samples (N) are required to accurately portray our input signal. The number of frequencies (phasors) that we calculate is generally the same as the number of input samples (N). If we say that separation betw een phasors can be represented by the constant 0, then: No = OO s or 0 = ros iN With this equation, we can now digitize the frequency scale of our DFT such that the spectrum can be written in terms of k instead of ro (since ko = e), N- I (I ) X(kO) = Lx( n)e- jk5Tsn n=O From previous diagrams, we know that: and 0= (tJ s N We can then rewnte mT S • 2 =_s -2 1r =- 1r (tJ N (tJs N LECTURE 4 4-14 180 U'IJ' ~ University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN If we substitute this in ( I) we now have: X n (k) = :L x( n )e -J'N n=O - J- N-I .2nkn Th e portion of the equat ion e N is known as the twiddle/actor and is abbreviated to WNo By grouping this part together, we will later be ab le to reduce the amoun t of compu tations required. Incorp orating the tw iddle factors gives us a DFT and inverse (IDF1) defmed as: N-I N- I - .2 1< XN(k)= Lx(n)W:' n=O and x(n)= -..!.- I Xn ( k ) WN N k=O Im LECTURE 4 4- 15 18 1 u4D Universityof Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN P R A C T IC A L CONSID ERATIONS N-I . STA ND ARDDF T XN {k) = Lx.{k )W~.O Sk ~N -I .· 0 • AN EXA MPLE OF AN 8 POINT DF T X N (k) = • I 7 X . { k ) Wr . k = 0 . 1. 2 •...• 7 ••0 WRITING THIS OUT FOR EA CH VALUE OF n X. ( k) = x{O)Wr + X{ I ) W ~ ' + • • • EA CH TERM SUCH AS x(D )Wr + x (7 ) Wr . k = 0. 1•. . .• 7 REQUIRES . MULTIPLICATIONS TOTAL NUMB E R OF MUL TIPL ICA TlONS REQUIRED ••• - 64 DO NOT FOR GET THAT EACH MULTIPLICA TION IS COMPLEX 8 POINT DF T REQUIRES' 2 - 64 M UL TIPLICATIONS 1000 POINT OFT 1 0 0 0 2 - I MILL IO N MULTIPLICATIONS III II AND ALL THESE NEED TO BE SUMMED We will now consider a practical example ofa OFT to give some idea of the computations involved. We will use the example of an 8-point DFT. This consists of eight samp les in the time domain, which will be converted to eight phasors (frequencies) in the frequency domain . K will take values from 0 to 1t to give us eight points. The third equation on the above diagram shows the eight terms of the OFT. Each term consists of a complex multiplication that must be summed with the other terms. This means that we have eight complex multiplications and seven complex additions to be calculated. Each of the terms requires that the eight harmonics be expanded. This is because n takes values from 0 to 7. So we now have 8 x 8 complex multiplications and 8 x 7 complex additions. . This can be generalized such that for any N-point OFT, we have N multiplications and N(N- I) additions. In reality, an 8'-point OFT would not be of any practical use. We would normally deal with around 1024 points. This-means that we would have to perfoim approximately I million multiplications and I million additions. This is too many computations for a practical real-time application. 2 LECTURE 4 4-16 182 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN FAS T FOURIER TRANSFORM SY MME TRY PROPER TY w~ "" = -w: PE RIODICIT Y PROPERTY w~·· = w: i -- SPLI TT IN G THE OFT IN TW O X. (i) = ~ .(2r). W,;ri + ~• • (2r + I). W ~,..". !. . I ~-l OR X .w U ) ,. z., ~ r(2,) .(w:r· + w; ~. %(2, + 1).(W; r .l MANIPU LATING TH E TWIDDLE FAC TO R THE FAST FOUR IER TRANSFORM : o ". ~ !L. 1 XN(k) - , oo . .. ................. . ....... . ............. .. . ........ r.O t !!-. x(2r ) .W : , : + w : r x ( 2r + I) .W:,: ~ .. We have already seen that the DFT requires many computations to arrive at an answe r. It is basically inefficient because it doesn't exploit the symmetry and periodicity properties of the twiddle factor WN. All FFT algorithms exploit these two properties to improve their efficiency. The Symmetry Property states that: W~+N/2 = - W~ This is easy to understand if we remember that the twiddle factor is periodic with respect to N. So an increment of 1/2 N must give us the symmetrical opposite ofN. The Periodicity Property is such that: kN W N + -Wk N Since WN is a periodic function with a limited number of distinct values, if we reduce the amount of times that this is calculated, we can speed up the transform. One possibility is to split the transform .into two series, one consisting of odd sequences and the other consisting of even sequences as shown below: ~-I 2 X N (k) = L x(2r)W~(2 r) + L x(2r + I)W~(2 r+ l ) ·2 ~-I r: O r:O We can increase the amount of identical terms with the following simple manipulation: W~ ( 2 r + J ) = W~ ( 2r ) x W~ L ECTURE 4 4-17 183 u~ University of ...., Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN As W~ does not depend on the term r, we can take it out of the summation. This gives us the following equation : ~ -I 2 X N (k ) = I ~- I 2 x (2r)W~ rk + w~ I x(2r + l)w~rk re O r eO Manipulating the twiddle factor gives us the fmal FIT: ~-I 2 ~- I 2 reO XN( k ) = I X ( 2 r )W:'2 +W~ IX(2 r + I ) W~~2 r eO So what have we achieved? LECTURE 4 4-18 184 University of Hertfordshire STUDENrs GUIDE THE FREQUENCY DOMAIN TI M E SAVINGS ~ ·1 t·, XN(k ) . "" L. x( 2 ,) W N/ 2 ,k • , -0 ~ WN k ~ L. x(2 '+1)W N12 ,k ,"0 "" • "",., ~"oo"' 1 , -0 ~ '"0" ""en."" ," 0 " , N/2 MULTIPLICATION S FOR I POINT FFT => 4 2.4 2.4 - 36 MULTIPLICATIONS SAVING 64·31 - 21 MULTIPLICATIONS FOR 1000 POIN T FFT => 500 2 + 5 00 2 . 500 - 50,500 MUL T IPLIC AT IONS SAVING 1,000,000·50,50 0 - 945, 000 MUL TIPLI CA T IO NS TIME SAVINGS. ASSUME 50n S CY CL E TIME (CLO CK ? I PO IN T FFT SA V E S 1.2 lIS 100 0 POINT FF T SAVES 47 .25mS • • . We now have another way of expressing our OFT, but why is this new method any faster than the previous version? Let us look again at our example 8-point OFT and see if the equivalent 8-point FFT is any faster. The 8-point OFT had to perform 8 x 8 = 64 multiplications. Our FFT performs two 4-point OFTs and adds them tog ether. Unfo rtunately, we must multiply the odd series OFT by W~ before we can add it to the even series OFT. This gives us 4 2 + 4 -t 4 = 36 multiplications to perform, as opposed to the 64 in the original OFT. The savings become even more apparent if we use a toOO-point OFT. With our new equation, we perform 500 2 + 500 2 + 500 = 50,500 multiplications as opposed to I. million multiplications for the original OSP. If we look at time savings using the example of a OSP which has a clock time of SOns, we can save 1.2J.1s on the 8-point OFT and save 47.25ms on a lOoo -point OFT. What is the clock frequency of the above DSP? Such time .savings make it possible to start contemplating performing a Fourier Transform in real time . However, the process can be spe eded up even more as we will see on the next diagram. 2 LECTURE 4 4- 19 I~ TEXAS NSlRUMENrS 185 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN DECIMATION IN TIME • SPLITTING ORIGINA L SERIES INTO lWO IS CALLED DECIMATION IN TIME • L ET US TAKE A SHORT SERIES WHERE N = 8 • DECIMATE ONCE CALLED I RADIX · 2 ] SINC E WE DIVIDE BY lWO n • (0, 1, 2, 3, ... 5, 6, 7) • DECIMA TE AGA IN n =( 0, 2, 4, 6 ) and ( 1, 3, 5, 7 ) n= ( O, 4) ( 2, 6 ) (1,5) an d (3, 7) • GIVES A SAVING OF N2.(NI2)1og,N MULTIPLICATIONS 1024 POINT DFT .',048,576 MULTIPLIC ATIONS 1024 POINT FFT = 5120 MULTIPLICATION • DECIMATION SIMPLI FIES MATHEMATICS BUT THERE ARE MORE lWlDDLE FACTORS TO CALCU LA TE. • A PRACTIC AL FFT INCORPORATES THESE EXTRA FACTORS INTO THE ALGORITH M The process of splitting the OFT into two is known as decimation in time (because we are splitting in the time domain). In the previous example, we only split the series into two once, but there is no reason why we can't continue splitting the sequences until we are left with a series of 2-point OFTs. This particular method of speeding up a OFT is kno wn as a Radix-2 FFT, since all of the final OFTs have two parts. With such a method, the savings become tremendous. We will use the example of a 1024-point OFT on a OSP with a 50ns clock cycle. Such a OFT would require 1,048,576 complex mu ltiplications and 1,047,552 complex additions. Assuming that a complex multiplication and addition can be performed in one clock cycle, the tim e to process the OFT would be 104ms. The computational time for the 1024-point Radix-2 FFT is 768 J.1s. The ratio of OFT multiplications to FFT multiplications is 204.8 , which means that it is possible to pe rform the FFT in real time. It should be noted that the above example is an approximation, and does not take into account the clock cycles required to initially calculate thetw iddle factors . The mathematics of the Radix-2 FFT are considerably more simple than those of the OFT, as we only have to calculate 2-point OFTs. However, the twiddle factors produced will di ffer at each stage of decimation. The goal of a good practical FFT implementation should be to incorporate these extra factors into the mathematics without adding unnecessary computational complexity to the OSP algorithm. LECTURE 4 4-20 186 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMA IN 4 PO INT FFT ; · · · · · · · ·· ·· · · · · · · · · ·· · · · ···· · · · · · · · · · · ·· · · · · · · · ·· · · · · · · 3 · · · ·· · · ~ · · · · · · · · · · · ': : ~. L ET US CONSIDER AN E XA MPL E WHER E : N-4 x.{k ) -I ~~ xlnjW." . .. . . . . .... . . .. . . .. . .... ... . . . .. . . . . . . . . .. . . . .. . . . ... ..Il :' •• • • • • ••• • • • • • • • • ••• • • • • •• • • • • • • • • • • ••• • • ; • • •• •• • • • • • • 'k' .,. ••• •• • • •• • ;~ : ••" : : • DECIMATE IN TIM E INT O 2 SE RIES : n · (0 , 2) and (I, 3) · X.(k) - L x( 2 r j W z ,. 0 k + W. k L xI2r+1)W z ,. 0 k : : : : : : -( ....... .. ........................ . .... xI0) + xI2) W z ) + W. (Xll) +xI3)W z ) ; ...................................... ~ ..··:J¥·k· ·· ····· ···": N : • :• REM EM B ER : WE HAVE TWO TWIDDLE FACTORS. CAN WE RELATE THEM? WN • • - . . . .. . . . . .. .. .. . .. . ... :• NOW OUR FFT BECOMES : k WZ =. . j.l.1L. k =a -j~ . 2k =W.2k 2 4 2k . ) • -(xI0) +xI2)W. ) + W . (x(1)+xI3 ) W . k 2k We will now look in more detail into how FFT works. We will use the examp le of a 4-point FIT to keep things simple, but the same principles apply to all Radix-2 FITs. First we start out with the 4-point DFT, which can be defined as: 3 X 4 (k) = Lx(n)W,f' o ,k = 0,1,2 ,3 We now decimate this into two series ofn = (0, 2) and (1 ,3). This gives us the equation : 1 r=O 1 r=O X 4 (k) = LX(2r)Wf' +WtL x(2r +l)Wf' If we expand the two summations, we have: X4 (k) = [x(O)+ x(2)W;]+ Wt[x(l) + x(3)W;] We have two different twiddle factors; W~ and W; . However, these two can be related to each other as W~ = wl k • From this, our FIT expansion becomes: X 4 (k) = [x(O + x(2)Wl k ] + W; [x(l) + x(3)Wl k ] ) The expansion isn' t too complicated for a 4-point FFT, but as the points increase, so does the complexity of the expansion. To reduce this complexity, the expansion is normally expressed in a different manner, which we will now see. LECTURE 4 4-21 187 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN FLOW DIAGRAM • TWO DFT s ~'.;~)~ ' ; ~ ;o·) ·" ., ~/k~'" ' ~ 4~ :;1" : ;~,' ~ ~;k): ' ~ ~ ~ .·1·. ~"3·1 ') ' :• • WRITING OU T VALUES FOR k - 0 ONL Y o · ·· · ·· ~ · ·· · · · · · ·· ·· ··o .. .... . . . . ... . . . . . . . . .. .. . . . . . . . . .... ... . .. . . . . . . . . . .. -' '': : : X. (O) ..... . . . .. . ... . . .. ..... .... . .. . ... ..... ........ . ... .. REP RESENT WITH A FLOW DIAGRAM - I .(01+ .(2) W 4 ) + W 4 1. (1) + .(3) W 4) . (O) O ~ " o .( 2 ) 70X4 (OI Xu 0 Xu ®-W ~ o .(l )O ~" "(3) • TH IS ONL Y ON E QUAR TER OF TH E F L O W DIAG RAM Instead of expressing the FFT as an equation , we can express it as a signal flow diagram. The flow graph shows the expansion of the equation when k = O. The number in the circles represents the power of W4 required at the stage. If there is no multiplication required, the line has no multiplication factor. The flow graph effectively has time as the input on the left and fr equency as the result on the right. It is eas iest to read the diagram from right to left. From the diagram, we can see that: where Xo2 = x(O) 4- w2 x(2) and xl3 = x(l) + w2 x(3) This gives us: Xi O) = x(O) + W2x(2) + w2 (x(l ) + W x(3» J This , of course, is the original expans ion for k = 0. The above diagram only shows the expansion for the first value (i.e., one quarter) of the equation. However, the rest can be shown in a similar manner. LECTURE 4 4-22 188 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN FULL FLOW DIAGRAM • ~ ( .( 0) + .(2 ) W. (. (0) + . (2) W R ITI NG OU T A LL VA LUES F O R k o 0 ) + W. 0 (.(1) + . (3 ) W. ) NOTICE ~ w. 2 )+ w.' 2 ( . ( 1) + .(3 ) w. 0 2 ) ) - ( . (0) + . (2 ) W. - ( ' (0 ) + ' ( 2) W ; o . ... -e -1 -" ·.. :W ··· ··· ·· · · ·2; · ·· · · · · · · · · · ·: 0: -1- W. ) +W. ( . (1) + . (3) W. 3 ) + W . ( . ( 1)+.(3) W; ) : W..-.. . .. . ... . ... .. . .W. ... . e" -.,- . ... : .. . .. .. . . . . . . . . . . . . . . . . . . . . . . . . 6 -1..11L · 6 2: • (O )O~ .(2) O~o.--~~---.,.OX.(1) SPOT THE BU TT E R FLY? .(l)O~ .( 3 ) ~0'----{ We can now see the full expansion of both the equation for a 4-point FIT and the signal flow graph that represents it The symmetry of the expansion is now visible. This symmetry extends to any-size FIT, which means that implementation is relatively easy. This symmetry is normally referred to as a "butterfly", as the individual "X" parts resemble the wings of a butterfly (assuming that you have a very lucid imagination!). It is important to notice that: W4 4 O W4 and 6 W4 2 W4 These are both examples of symmetry and periodicity that were mentioned earlier. Both of these properties combine to form the key as to why there are so few complex multiplications in an FIT compared to a DFf, since the same twiddle factor can be used in numerous places. LECTURE 4 4-23 189 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN THE BUT T ERFLY TYPI CAL B UTTE RF LY TWIDDLE CONVERS IONS + x • • x, w: x, W1 :. w~. -J W. -1 : W:· J 4 POINT FFT BUTTERFLY 4 POINT FFT EQUATIONS x,. (x. + x J ) X, X, + W: W,' w.' (x,.x . ) X, . (x,. x,) + (x, .x, ) (x ... x.) XJ - (x,+ XI)· W: XI · (x, - x,) - (xt-x l ) X, The butterfl y is more typically implemented as shown above. It consists of two inputs, two outputs, and an optional multiplier. The upper output always consists of the upper input added to the lower input, while the lower output always consists of the lower input subtracted from the upper inp ut. Using this notati on, FFTs can he elegantly represented. A 4-point FFT is shown in the diagram, which consists of two stages each with two butterflies. The first section of the FFT, when k = 0, is derived as follows : Xo= a+ Wfb where and There fore, This appears to be different from the ear lier definition of Xo, which was: (2) X 4(0) = x(O) + w x(2) + wf [x(l) + wf x(3)] f LECTURE 4 4-24 190 u~ '11 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN However, as shown on the transparency, the twiddle factors resolve down to simpler levels. For example: or W4 1 =e 211) - J(- I 4 o = e" 090" = cos(- 90) + j sin(-90) =- j As w = 1 , we can see that both equations (l) and (2) resolve to the following: t Looking at the parts of the FFT equ ation, we can see that they all share common parts. This reduces computation and allows easy implementation. L ECTURE 4 4-25 191 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN A PRACTICAL EXAMPLE TIME DOMA IN FR EQ U EN CY DO MAI N : Am p lit ude : Am plllu de :1 Time In T,) :0 • 0 2 3 Fr eq ue n c y : IkW)kHz: ' : ' • • {l ,O,O, I ) : SA MP L ED AT 10kH z ~ I: ~ ~".•. .J.l .~'.".•. ~~. !~,!~ . . ._ FFT CON VE RS ION x.·•.•• 1 + s,.", • 1+ 0+ 0 +1 • 2 :T . - 10 0 uS : W• .l.IL.. 2n . , S.7k Hz : : Nt 4" 100. 10 " : X l· ..•• ••.J(• •••• ) - 1-0 ••J(O.l ) - 1 • J X.·I' •• '.)· ('.+".1- (1+0). (0+1 ) - 0 X.· I". ·".) - J I".·•• ) - 11. 0) ••JlO. l ) . 1 • J . . . . . ..... .. ... .... . ............ . . . .. .. .. AMPLITUDE OF X.- .[1'+J' · .[2 So far, all of the emphasis has been on how to produce an FFT. At this point, it is probably worth reinforcing the theory with a simple example to demonstrate the conversion from time domain to frequency domain using an FFT. The example we are using is unrealistically simple, but it demonstrates the conversion process quite we ll. In the above example, we are sampling an input signal at 10kHz, which means that we take a sample every 1110kHz = 100J,ls. We are using a 4-point FIT for the conversion, so we take four succ essive samples from the inputs, which have the values: 1, 0 ,0, I Therefore, we have the following input values for the FFT: Xo = 1 Ifwe place these values in the FFT, our four output values are: Xo=2 XI = 1 + j We now have four values in the frequency domain, but how do we relate these values to actua l frequen cies? From our earlier definition of a OFT, we know that: X(k ) =X(w)k Where 21t ro ro = - and f= N 2 1t LECTURE 4 4-26 192 u~ 'IJ' University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN Since we are dealing with a discrete sampled signal, it is very important to remember that (J) is in radians per cycle, and not radians per second. Likewise, f is cycles per sample, not cycles per second . Therefore, if we are interested in frequency, the following is true: X(k) = ( 2~)<k) = (~ )<k)HZ So we know that our output frequencies increase in multiples of 1/4 of the sampling frequency = 205kHz. This means that we now know the amplitude of the signal for frequencies of 0, 2.5, 5, and 7.5kHz, as shown in the above diagram. We already know that the spectrum of an N-point FFT is symmetric around NI2 . Therefore, all values above NI2 are symmetric to earlier values, and can be discarded . We can see that the frequency domain is symmetric around k = 2 and therefore, the value X(3) can be discarded. From this, it can be seen that the N "real" input values are transformed to NI2 complex FFT values of any significance. LECTURE 4 4-27 193 University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN D SP & FFT • F A ST FOU RIER TRAN SFORM IS A GENERA L NAME FOR REDUCING OFT C OMPUTATI ON S W E C ONSID ERED RADIX· 2 HERE. BUT MANY OTHER AL GO RITH MS EXIS T. • THE SIM PLI FIED BUTTE RFLIES CAN BE IMPLE MENTE D WITH A DSP VERY EFFICIENTLY • • • SPECIAL FFT CHIPS I .. PLE .. ENT IT EVE N FA ST ER BU T DSP a ARE PROGRA .... A B L E A ND CA N P ERFORM OTHER OPERAT IO NS O N TH E S IGN A L • F FT REQUIRES ADDRE S S SHUF F LING FOR FASTER DA TA T ABLE ACCESS • MOS T DSP a CAN PER FORM SHUFFLING IN THE BA CKG ROUND • M O DE R N D SPa C A N TAKE FFT OF 1D24 SAMPL ES W EL L UNDER 5mS Fast Fourier Transform is a generic name for any method of reducing the computational complexity of a DFT. There are many d ifferent methods of reducing the complexity (although Radix-2 is the most popular). Many of the other transforms even further-reduce the number of multiplications that must be performed, but at the expense of extra additions. Previously, this was useful because a multiplication took long er to perform than addition. However, most DSPs can perform a multiplication and addition in a single clock cycle, so any increase in additions has an equal penalty on performance. The simple structure of the FFT allows it to be quickly and easily implemented on DSP chips that are specifically designed to perform FFT. Since the algorithms are implemented in hardware, such chips achieve fast FFTs. The maj or disadvantage of such implementations is that they are not programmable and do not have the capability of other signal processing functions. However, a general-purpose DSP can also perform other operations on an input signal apart from the FFT. Ifwe look at the inputs to the butterfly FFT, we can see that the inputs are not in the same order as the output. To perform an FFT quickly, we need a method of shuffling these input data addresses around to the correct. order. This can be done either by reversing the order of the bits that make up the address of the data, or by addition. Many DSPs have special addressing modes that allow them to automatically shuffle the data in the background. Using such meth ods, modem DSPs can easily perform a 1024-point FFT in well under 5ms. LECTURE 4 4-28 194 Univers ity of Hertfordshire STUDENrs GUIDE THE FREQUENCY DOMAIN SUMMARY • • • FR E QUENC Y DOMAIN INF OR M ATI O N FOR A S IGNA L IS IM P O R TAN T FO R P R O C E S SIN G SIN USOIDS CAN BE REPRESENTE D BY PHASORS FOURIER SERIES C A N BE USED TO REPRESENT AN Y PE RI ODIC SIGNA L FOURIER T R A N SFORMS A R E USED TO TRANSFORM S IGNA LS • FRO M TIME TO FREQUENCY DOMAIN • FRO M FRE Q U E N C Y TO TI M E DOM A IN • • • • • OFT ALLOW T R A N S F O R M OPERATIONS ON SAMPLED SIGN A L S OFT COMPUTATIONS CAN BE SP EEDE D UP BY SPLITTING THE OR IG IN AL SERIES INTO TW O OR MORE SE R IES FFT OFFERS CONSIDERABLE SAVINGS IN CO MPUTATION TIME DSPs CAN IMPLEM E N T FFT EFF IC IE NTL Y In this lecture, we discovered ways of working in the frequency domain. The main points in summary are: • Frequency domain information for a signal is important for processing. This information may be used to reveal propertiesof a signal which are not immediately apparent from its time domain picture . Sinusoids can be represented by phaso rs The Fourier Series can be used to represent any periodic signal Fourier Transforms can be used to transform signals from time to frequency and from frequency to time DFTs allow transformations on sampled sign als DFTs can be speeded-up by splitting into sub-series FFTs offer considerable savings in computation time. On a SOns clock cycle DSP, a 1024-point DFT takes l04ms and a l024-point Radix-2 FFT takes less than 5ms. DSPs can implement an FFT efficiently • • • • • • • LECTURE 4 4-29 195 u4D University of Hertfordshire STUDENT'S GUIDE THE FREQUENCY DOMAIN REFERENCES Bateman, A. and Yates, W. [1988J. Digital Signal Processing Design, Pitman Pu blishing, London, UK Burrus, C. ~. and Parks, T. W. [1985J. DFTIFFf and Convolution Algorithms, Wiley Interscience, NY Chan, K. K. [1988J. "Im plementation of Fast Cosine Transforms with Digital Sign al Processors for Image Compression," SPIE Vol 914 Medical Imaging II Chitprasert, B. and Rao, K. R. [ 1990]. "Discrete Cosine Filtering ," Proceedings of the 1990 IEEE International Conference on Acoustic, Speech and Signal Processing Jayant, N . S. and Noll, P. [1984]. Digital Coding of Waveforms, Prentice-Hall, Englewood Cli ffs, NJ Lynn, P. A. and Fuerst, W. [1998]. Introductory Digital Signal Processing, John Wiley and Sons, UK Martin, J. D. [199 IJ. Signals and Processes, Pitman Publishing, London Oppenheim, A. V. and Schafer, R. W. [1975 and 1988]. Digital Signal Processing , Prentice-Hall, Englewood Cliffs, NJ Papamichalis, Panos (ed.) [199 1]. Digital Signal Processing with the TMS320 Family, Volume 3, PrenticeHall, En glewood Cliffs, NJ Pratt, W. K. [1978]. Digital Image Processing, John Wiley and Sons, UK Stafford, R. H. [1980]. Digital Television, Bandwidth Reduction and Commun ication Aspects, John Wiley and Sons, UK LECTURE 4 4-30 ',' 196 E LECTRICA L ENGINEERING 113L D igit a l Signal Processing Lab o ratory @ 1998, Dr. M. Wer t er Le c t u r e 8: Fast Fourier Transfo r m Introd uction The Discrete Fourier Transform (DFT ) is an important tool in many dig ital signal processing systems . It converts informat ion from the ti me domain to a frequency domain r epr esentation. Methods that ar e computationally efficient for computing th e DFT are known as t he Fas t Fourier Transform (F FT) algorithms. In th is experim ent we ar e going t he implement FFTs on t he TMS320C542. Fast I<burier Transform T he DFT of a fin ite length sequence of length N is defined as X(k) =L N-I n=O k x(n)WN k = 0,1 , . . . , N - 1 (1) where W N is the com plex-valued phase factor WN = e- j 2>r/N It is clear from equat ion (1) that a direc t calculat ion of the DFT summation requires N complex mult iplications and (N - 1) complex additions for each of the N output samples. Co nsequently, to com put e all N values of the DFT requires N 2 complex mul t iplicat ions and (N 2 - N) complex additi ons (where each com plex multiplicat ion requires 4 real mult iplicat ions and 2 additions). This direct computation of the DFT is inefficient because it does not exploit the symmetry an d p eriodicity of the phase factors W N. Many methods that efficiently compute the DF T have been developed an d t h ey are collectively called FFT algorithms . In our exp er iments, we will focus on the radix-2 decimat ion-in-tim e (DIT) FFT algor ithm. The decimation-in-time FFT divides the len gth N (where N is a power of 2) input sequence int o two groups, one for the even indexed samples and the oth er for the odd indexed samples, Two length N/2 DFTs are then performed on these subsequences, and their outputs are combined to form the length N DFT. The met hod ology is illus trated by the following equat ions. First the input sequence in equation (1) is divided into the even an d the odd subsequences: X (k ) =L If -I x(2n)W~nk +L .If - I x(2n + I )W~n+l)k n =O =L If -I If - I x(2n) W~n k + wt L x(2n + I)W~n k (2) n =O n=O n= O By the substitutions 2nk _ ( _j2 >r / N ) 2nk _ Wnk WN - e J:!.. 2 and XI (n ) = x(2n) , and x 2(n ) = x (2n + 1) the equa tion becomes 'f - I X(k) = L If -I xI(n)W~k + wt L x2(n)W%-k = XI(k) + Wt X2 (k ) (3) n=O n=O In equation (3), X1(k) and X 2(k) are two Nf2 -point D FTs and t hey can each be further divided into the even and odd sub-sequences to form four N/4 point D FTs. In this way, by recursive decimation , the large DFT is broken into a combination of sm aller DFTs. T his is cont inued until only two-p oint DFTs , which are called radix-2 "butterflies," remain. T he butterfly is the core calculation unit of the FFT. the WN terms appear as coefficients (called twiddle factors ) in the FFT calculation. The entire FFT is performed by computing butterflies in a "but terfly-gro up-stage" fashion in a three-loop algorithm. A complete 8-point DIT FFT is illustrated in Figur e 1. 1 197 xtO) xt4] x{2) xt6) x{1) xtS) xt3] X[O) ~ X[1] X(2) ~ X(3) X[4] ~ X[S) X(6) xfll ~ -1 X(7) Figur e! 1: Eight point orr FFT Each pair of elements in Figure 1, which are combined, represent a butterfly, and the entire FFT is made up of butterflies organized in different groups and stages. The first stage consists of four groups with one butterfly each . The second stage consists of two butterflies each, and the last stage has one group of four butterflies . Each butterfly has one multiplicative coefficient W~ and has two in put points. Every stage contains N J2 = 4 butterflies. For t he more general case where N is a power of 2, the F FT will consist of log2(N ) stages and ~log2 ( N ) butt erHies. The computation of each butterfly requires one complex multiplicat ion and two complex additions . Thus, the entire FFT is computed in ~10g2(N) complex mult iplicat ions and Nlo g 2(N ) complex additions. Compared to the direct comp ut at ion of the DFT, the number of the multiplicat ions and additions is greatly reduced. Not ice that whereas th e output sequence is sequentially ordered, the input sequence is not . This is a consequence of repeatedly div iding the input sequence into subsequences of the even and the odd samples. Therefore, th e inp ut sequence must be scrambled before the butterfly calculations can begin. This is accomplished by a process called bit-reversal. Bit-reversal operates on the binary numb er that represents the index (position) of a sample within the array. The bit-reversed position has the index, which is the tr anspose of the binary num ber representin g the original index ; for example , the transpose of the 3-bit binary number 001 is 100. T he DIT FFT can appear in alt ernate, equivalent forms, one of which is shown in Figure 2, where the input is in normal order b ut the out pu t val ues are in bit-reversed order. xtO) X[O) x[I ) x(2) xt31 x(4) xtS) xt6) X(4) X(2) X(6) ~ X(1) X[S) X(3) ·1 ~ ~ ~ x(7) Figure I 2: 8-point DIT FFT with in put in normal order and ou tput in bit-reversed order . ..: x W=C+j(-S) X(7) I' . <,. jy', x',+jy', x, +jy , Fi gure . ·1 3: Radi x-2 DIT FIT butterfly 2 )98 The fundament al computational structure needed to compute the FF T is t he generalized butterfly graph shown in Figure 3. The var iable x and y represent the real and im aginary parts, respectively, of a sample. Th e twiddle factor is also div ided into its real and imaginary parts as w~ = e - j 21rk / N = cos(2rrkjN) - jsin(2rrk/N) = C - j S The lower input (XI + jyJ) is multiplie d by t he twiddle factor (C - j S ). The res ult of th is mul tiplication is added to t he upper part inp ut (xo + j yo) to prod uce the upp er outp ut (x~ + jy~) and subtracted from the up per input to produce the lower output (x~ + jyD: (4) t l = C X XI + S X YI t2 =C x~ X VI - S X XI x~ = Xo y~ = + tl yo + t2 tl = Xo = yo - V~ t2 T he b utterfly produces two comp lex outputs that bec ome but terfly inputs in the next stage of the FIT . Since each st age has t he same number of but terflies, t he number of butterfly inp uts an d ou tp uts remains const ant from one st age to the next. This makes an "in-place" compu t ation possible, which writes each butterfly output over the corresponding inpu t (x~ over Xo etc.). In t he in-place imp lement ation, the FFT results end up in t he sam e memory locations as the original inpu ts. Implem ent ation of t he r adix-2 FFT The core part of an FFT algorit hm is the computation of the butterfly. As given in equation (4), the computation of each bu t t erfly requires four real multiplications and six real additions. We would like to use four parallel multiply j add inst ruct ions to implement eight of th ese operations . The easiest way to achieve this is to comp ute th e coefficient multiplications from each butt erfly in par allel with the additions from the previous butterfly. Since there are ~log2(N) butt erflies in a length N FFT, each operation removed from the butterfly calculation shaves ~log2(N) cycles of the time to comput e th e entire FFT, which is a significant number. For example, if N = 1024, we can save 5120 instructions . A har d war e featur e of the TMS320C542 that can increase the speed of an FFT , by up to 30%, is the bit-reversed addressing mode . The DIT-FFT requires that the input data must be scrambled in a bit-reversed pattern before the butterfly calculations can start. This process is usually done by swapping data elements between various memory locations and , hence, requires extra instructions. However , the bit-reversed addressing feat ure, which is a special form of the indirect addressing mode, allows the bit-reversal to be done in conjunction with the read-in of the data wit hout any extra cycles. In th e bit-reversed addressing mode, one auxiliary regist er is used to point to the physical location of the data valu e in th e sample array. This auxiliary regis ter is increment ed in a bit-reversed manner each time the r e gist er is access ed. For example, if N = 8, where N is t he FFT length, and base..adr is the ad dr ess of t he first element in the input , then the addr esses will be : base.adr. v base.adr-l-d, base..adr+ 2, base..adr + 6, bas e..adr+ l , base..ad.r+5, base..ad.r+3, base..adr+7. 3 199 ...
View Full Document

This note was uploaded on 11/06/2010 for the course EE 113 taught by Professor Walker during the Spring '08 term at UCLA.

Ask a homework question - tutors are online