DSPRelated.com
Forums

N-pt DFT where n != power of 2

Started by dspdummy August 12, 2009
>dspdummy wrote: > > ... > >> ok. I got it thanks mostly to the cnx reference above and that
wikipedia
>> reference. Thanks a lot for the help! If anyone is interested I can
send
>> you the Matlab M file for the 15-point dft using prime factor
algorithm.
>> Now I can also extend it to 960 or 120 pts. > >My hat is off to you for your enthusiasm and comprehension. >Nevertheless, I doubt that your finished code will be as fast as FFTW. >Why not use the tools that are available for the cost of downloading? > >Jerry >-- >Engineering is the art of making what you want from things you can get. >??????????????????????????????????????????????????????????????????????? >
I just might, but I always like to start with a basic understanding of the principles via use of very simple Matlab scripts before plowing thru some very high performance optimized code. I get about 2 times the no. of complex multiplications needed compared to a standard (power of 2) complex FFT. N*(8+log2(N/15)/2) instead of (N/2)*log2(N). I haven't looked at those FFTW benchmarks yet. Do you have an idea what I could get for a complex input?
On Wed, 12 Aug 2009 17:19:42 -0500, "dspdummy" <chris.glass@orange.fr>
wrote:

>I know that an N-pt DFT can be replaced by an FFT to save operations if N >is a power of 2. What can I do if N is NOT a power of 2 to avoid the DFT? I >need, nevertheless, the exact result. > >Thanks
Hi, Ya' might take a look at: http://www.dsprelated.com/showarticle/63.php Good Luck, [-Rick-]
dspdummy wrote:
>> dspdummy wrote: >> >> ... >> >>> ok. I got it thanks mostly to the cnx reference above and that > wikipedia >>> reference. Thanks a lot for the help! If anyone is interested I can > send >>> you the Matlab M file for the 15-point dft using prime factor > algorithm. >>> Now I can also extend it to 960 or 120 pts. >> My hat is off to you for your enthusiasm and comprehension. >> Nevertheless, I doubt that your finished code will be as fast as FFTW. >> Why not use the tools that are available for the cost of downloading? >> >> Jerry >> -- >> Engineering is the art of making what you want from things you can get. >> ??????????????????????????????????????????????????????????????????????? >> > > I just might, but I always like to start with a basic understanding of the > principles via use of very simple Matlab scripts before plowing thru some > very high performance optimized code.
I applaud that.
> I get about 2 times the no. of complex multiplications needed compared to > a standard (power of 2) complex FFT. N*(8+log2(N/15)/2) instead of > (N/2)*log2(N). I haven't looked at those FFTW benchmarks yet. Do you have > an idea what I could get for a complex input?
Recall that even a real FFT has to deal with complex results after the first set of twiddles. I would expect about N*log2(N). Jerry -- Engineering is the art of making what you want from things you can get. &#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;
On Aug 14, 9:46&#4294967295;am, "dspdummy" <chris.gl...@orange.fr> wrote:
> > Matlab stuff: >
I don&#4294967295;t have MATLAB, so I can&#4294967295;t help much with your code. But you seem to be doing a column of 15 point DFTs, followed by a column of 64 point FFTs (or IDFT and IFFT, if your exponent is different). Although you said that it works, I&#4294967295;d be a bit cautious unless you can compare your results against a &#4294967295;known good&#4294967295; program (e.g.: random number inputs to a double precision &#4294967295;known good&#4294967295; DFT routine with N = 960). Using a 2 stage algorithm is probably a good idea. They&#4294967295;re much easier to visualize, the interstage twiddles are easy to figure out, and for sequential in, bit reversed out, the final bit reversal is a matrix transpose. So for a 960 point, sequential in, bit reverse out 2 stage graph, you&#4294967295;ve got a first column of sixty-four 15 point DFTs, followed by the interstage twiddles, followed by fifteen 64 point FFTs. Since you&#4294967295;re doing things in fixed point, I would think that you could improve accuracy by doing the 64 point FFTs in the first column (use DIT - your round-off error will be less). In general, you want to have the &#4294967295;difficult&#4294967295; twiddles more towards the right side of the graph, as opposed to having them in the earlier stages. Doing the 64 point first, you&#4294967295;d have a first column of fifteen 64 point FFTs, followed by interstage twiddles, followed by a second column of sixty-four 15 point DFTs, followed by a matrix transpose. Using a general purpose scheme to do the 64 point FFTs keeps things simple, but inefficient. You&#4294967295;d be better off with an FFT that was specifically designed for 64 points (e.g.: an 8x8 that uses a straight through 8 point routine and look-up tables for all twiddles). An 8 point FFT only requires internal twiddles of 0, j and the +/-.707... types. The 8x8 requires additional interstage twiddles to stitch the two columns together, so you don&#4294967295;t lose much accuracy with them. And they&#4294967295;ll use fewer multiplies than the general purpose radix 2 approach. You could also break down the 15 point DFTs by using a 5x3 or 3x5 arrangement, but I don&#4294967295;t think (offhand) that it would provide a dramatic improvement (if you profile a 64x15, you could get a better idea). And depending on how you break it down, you could complicate the matrix transpose at the end. You could also improve the efficiency (and numerical accuracy) of the 3 and 5 point DFTs, but I don&#4294967295;t know how important that is for you. Kevin McGee
>On Aug 14, 9:46=A0am, "dspdummy" <chris.gl...@orange.fr> wrote: >> >> Matlab stuff: >> > >I don=92t have MATLAB, so I can=92t help much with your code. But you >seem to be doing a column of 15 point DFTs, followed by a column of 64 >point FFTs (or IDFT and IFFT, if your exponent is different). > >Although you said that it works, I=92d be a bit cautious unless you can >compare your results against a =91known good=92 program (e.g.: random >number inputs to a double precision =91known good=92 DFT routine with N
=3D
>960). > >Using a 2 stage algorithm is probably a good idea. They=92re much >easier to visualize, the interstage twiddles are easy to figure out, >and for sequential in, bit reversed out, the final bit reversal is a >matrix transpose. > >So for a 960 point, sequential in, bit reverse out 2 stage graph, >you=92ve got a first column of sixty-four 15 point DFTs, followed by the >interstage twiddles, followed by fifteen 64 point FFTs. Since you=92re >doing things in fixed point, I would think that you could improve >accuracy by doing the 64 point FFTs in the first column (use DIT - >your round-off error will be less). In general, you want to have the >=91difficult=92 twiddles more towards the right side of the graph, as >opposed to having them in the earlier stages. > >Doing the 64 point first, you=92d have a first column of fifteen 64 >point FFTs, followed by interstage twiddles, followed by a second >column of sixty-four 15 point DFTs, followed by a matrix transpose. > >Using a general purpose scheme to do the 64 point FFTs keeps things >simple, but inefficient. You=92d be better off with an FFT that was >specifically designed for 64 points (e.g.: an 8x8 that uses a straight >through 8 point routine and look-up tables for all twiddles). An 8 >point FFT only requires internal twiddles of 0, j and the +/-.707... >types. The 8x8 requires additional interstage twiddles to stitch the >two columns together, so you don=92t lose much accuracy with them. And >they=92ll use fewer multiplies than the general purpose radix 2 >approach. > >You could also break down the 15 point DFTs by using a 5x3 or 3x5 >arrangement, but I don=92t think (offhand) that it would provide a >dramatic improvement (if you profile a 64x15, you could get a better >idea). And depending on how you break it down, you could complicate >the matrix transpose at the end. > >You could also improve the efficiency (and numerical accuracy) of the >3 and 5 point DFTs, but I don=92t know how important that is for you. > >Kevin McGee
Many thanks to everyone for the help you have furnished for solving my problem! This is a great website and I will recommend it to my colleagues. Bill
On Aug 17, 2:53&#4294967295;pm, "dspdummy" <chris.gl...@orange.fr> wrote:
> Many thanks to everyone for the help you have furnished for solving my > problem! This is a great website and I will recommend it to my colleagues.
(Technically, it's not a website, although many people access it through the Google web front-end. See: http://en.wikipedia.org/wiki/Usenet)