Monday, 30 March 2015

fft - DFT as convolution question


I have tried to make this question as readable and consistent as possible. The short of it, is that I am trying to ascertain how one gets from the math equation shown, (which I understand), to the correct implementation. (Which I also understand). However, how/why one gets from one to the other seems unclear to me.


Here is my setup:


We have an M length 'original' signal x[m]: (Equation 1)


x[m], m=0,1,..M1


We wish to take a 'zoom-DFT' of it, where the DFT length will be, N, where N>M.


We thus now have an N length zero-padded new signal, lets call it xp[n], where: (Equation 2) xp[n], n=0,1,...N1


Our discrete time variable is now n=0,1,...N1, and our discrete frequency variable is k=0,1,...N1.



The N-length DFT of a signal is defined as such: (Equation 3)


X[k]=N1n=0xp[n] ej2πnkN


If the nk in the DFT equation is written as (kn)2+n2+k22, then the DFT can now be written as: (Equation 4)


X[k]=ejπk2N N1n=0xp[n] ejπn2N ejπ(kn)2N


This expresses the full DFT, (i.e, the DFT from 0 to nyquist), as a convolution. To evaluate the 'Zoom-DFT', where we want to interrogate only frequencies between f1 and f2, and so the following minor modification is made: If we let FΔ=(f2f1)fs, and F1=f1fs, then the zoom-DFT, Xz[k] becomes: (Equation 5)


Xz[k]=ejπFΔk2N N1n=0xp[n] ejπn2N ej2π F1n ejπFΔ(kn)2N


We can then simplify, if we let: (Equation 6)


a[u]=xp[u] ejπu2N ej2π F1ub[u]=ejπFΔu2N


, where u is just some general dummy variable, (you will see why I picked it later). Anyway, now finally, we simply have: (Equation 7)


Xz[k]=b[k] N1n=0a[n] b[kn]



This, as we are told here, and here, page 183, is how the DFT can be re-written as a convolution.


Great! So far so good. The problem comes in interpretation vs implementation.


The problem:


My interpretation of (Equation 7), is that a[u] is of length N, and b[u] is also of length N, where we evaluate then for u=0,1,...N1. Then, we do a circular convolution of the series a[u] and b[u], to attain a series also of length N, (which is our desired DFT length), before we finally perform an element-by-element post-multiply with another series, also of length N. I dont so much care about post-multiply, that is not the point of this question. The point here is that a[u] is length N, b[u] is also length N, and u for both of them is evaluated for 0uN1, and then we do a circular convolution, modulo-N.


If you however now take this interpretation and implement it, you will not get the right answer. In fact, you will only ever get the right answer, if the u variable for evaluation of a[u] is taken from u=0,1,...N1, (same as above), BUT, the u variable for evaluation of b[u] is taken from (M1)uN1. So now the series b[u] is not length N anymore, it is of length N+M1. Then, if we do a modulo-(N+M1) circular convolution, and extract the indicies from M to M+N1, we will get the right answer.


...How, from equation-1, does one possibly get this 'correct' interpretation though? This, I cannot seem to explain, however I understand everything else.


FWIW, I have sat down and done the computation 'by hand', so I can 'see' why b[u] MUST be taken from (M1)uN1. (It would otherwise be impossible to compute the DFT for k=0,1,...N1). However I do not understand how one goes from (Equation 7), to this conclusion. Right now, to me, EQ-1 just looks like a garden variety circular convolution, with no hints as to the details for now a[u] and b[u] must be exactly evaluated.


Put another way - AFAIK, when we see equation-1, we simply assume that a[u] and b[u] are evaluated from 0uN1, and then we do a circular convolution, modulo-N - how/why then do we derive the proper evaluations for a and b from this equation??


Thanks!



Answer




You're interpretation that a[u] and b[u] are both N length sequences is incorrect. a[u] is of length N because it is effectively truncated by the length of your input sequence. The period of b[u] is not necessarily N. Depending on the value of N, b[u] may not even be periodic, so what you have in Equation 7 is a linear convolution of two sequences. One sequence (a[u] has length N, the other (b[u]) has length unknown. The limits of summation in the convolution equation are strictly set by the length of a[u] and don't imply circular convolution or anything about the length of sequence b[u].


Next you have to consider what range of values of b[u] are required to compute the convolution. Straight from the form of EQUATION 1, you see that b[u] in general must be defined in the range -[N-1] to N-1. The negative boundary is set by the most negative b[u] value required to calculate X[0] (i.e. you need b[0−(N-1)] for this calculation). The positive boundary is set by the maximum spectral step you want to evaluate in the chirp Z transform (X[K] max). From the way the question has been posed, X[k] is best interpreted as X[N1].


The reason your computations work with a range that only goes back as far as -[M-1] is because you started with an M length sequence and zero padded it up to length N which zeros out any b[u] values going back further than -[M-1].


There is an interesting write up of Chirp-Z along with Goertzel's algorithm that I think puts this in proper perspective: http://www.google.com/url?sa=t&rct=j&q=goertzel%20algorithm&source=web&cd=3&sqi=2&ved=0CEIQFjAC&url=http%3A%2F%2Focw.mit.edu%2Fcourses%2Felectrical-engineering-and-computer-science%2F6-341-discrete-time-signal-processing-fall-2005%2Flecture-notes%2Flec20.pdf&ei=rawjUd6sM-LE0QHJ2IGQAQ&usg=AFQjCNEFEo7R7VSmk5jWYon0f4WqZYpvzg&bvm=bv.42553238,d.dmQ


No comments:

Post a Comment

readings - Appending 内 to a company name is read ない or うち?

For example, if I say マイクロソフト内のパートナーシップは強いです, is the 内 here read as うち or ない? Answer 「内」 in the form: 「Proper Noun + 内」 is always read 「ない...