DSPRelated.com
Forums

FFT question

Started by Active8 March 2, 2005
Im' testing a real FFT routine in c++.

I've got a sine wave of 512 samples, norm freg = 1/32, with sample 1
at 0 and the last sample at -1 - a discontinuity.  It's going
through a real FFT routine of length 9.

I'm getting 10.log() spectral components of magntude on the order of
1.27 e37 ?

Something isn't right, eh?

Thanks in advance.
-- 
Best Regards,
Mike
Active8 wrote:
> Im' testing a real FFT routine in c++. > > I've got a sine wave of 512 samples, norm freg = 1/32, with sample 1 > at 0 and the last sample at -1 - a discontinuity. It's going > through a real FFT routine of length 9. > > I'm getting 10.log() spectral components of magntude on the order of > 1.27 e37 ? > > Something isn't right, eh?
Nope... just to clear one point up: What do you mean when you say "FFT routine of length 9"? A 512 pt FFT would require 9 levels of butterflies, is that what you mean? I'd look for wild pointers, or some problem with data formats. Does the FFT routine work with complex numbers internally? Have you allocated eneough enough for the imaginary parts? Do you feed float type numbers to a routine that expects double? Vice versa? These are the types of flaws I'd look for. Rune
On 2 Mar 2005 00:24:35 -0800, Rune Allnor wrote:

> Active8 wrote: >> Im' testing a real FFT routine in c++. >> >> I've got a sine wave of 512 samples, norm freg = 1/32, with sample 1 >> at 0 and the last sample at -1 - a discontinuity. It's going >> through a real FFT routine of length 9. >> >> I'm getting 10.log() spectral components of magntude on the order of >> 1.27 e37 ? >> >> Something isn't right, eh? > > Nope... just to clear one point up: What do you mean when you > say "FFT routine of length 9"? A 512 pt FFT would require 9 > levels of butterflies, is that what you mean?
Yeah. log2(512) Sorry, the code calls it length. It's from a book cd. The author passes log2(512)-1 into the FFT from the RFFT (real) wrapper. That "length", 8 is what determines the iterations of the outer loop. 8. That looks like 8 passes of 9 butterflies, right?
> > I'd look for wild pointers, or some problem with data formats. > Does the FFT routine work with complex numbers internally?
Yup. Trig recombination.
> Have you allocated eneough enough for the imaginary parts?
Yup. eneough is enough ;)
> Do you feed float type numbers to a routine that expects double?
Nope. I'll have to add some more runtime logging. I hate using other people's code ;) TANSTAAFL The Parks/McClellan and FIR code worked well once I got it working :) Thanks. -- Best Regards, Mike
Active8 <reply2group@ndbbm.net> writes:

> I hate using other people's code
You can try my ancient code, http://home.earthlink.net/~yatescr/fftk.c -- Randy Yates Sony Ericsson Mobile Communications Research Triangle Park, NC, USA randy.yates@sonyericsson.com, 919-472-1124
"Randy Yates" <randy.yates@sonyericsson.com> wrote in message 
news:xxp1xayb3o4.fsf@usrts005.corpusers.net...
> Active8 <reply2group@ndbbm.net> writes: > >> I hate using other people's code > > You can try my ancient code, > > http://home.earthlink.net/~yatescr/fftk.c > -- > Randy Yates > Sony Ericsson Mobile Communications > Research Triangle Park, NC, USA > randy.yates@sonyericsson.com, 919-472-1124
Hi Active8, Is you maximum vaue about 5.13? Do you get the right size peak where you expect it and it's only the bins that you expect to be empty that give you the nice low numbers? If so , you may be putting double precission real values , properly sampled, into a length 512 FFT and getting sensible answers out. Best of Luck - Mike
Active8 wrote:

> ... TANSTAAFL
Heinlein, or hearsay? (The Moon is a Harsh Mistress) Jerry -- Anyone who can not cope with mathematics is not fully human. At best he is a tolerable subhuman who has learned to wear shoes, &#4294967295;&#4294967295;&#4294967295; bathe, and not make messes in the house. .. Robert Heinlein &#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;&#4294967295;&#4294967295;&#4294967295;
On 02 Mar 2005 08:03:39 -0500, Randy Yates wrote:

> Active8 <reply2group@ndbbm.net> writes: > >> I hate using other people's code > > You can try my ancient code, > > http://home.earthlink.net/~yatescr/fftk.c
I could probably work with that. It's more readable than the Parks/McClellan code I c->c++ 'd. I think if I switch code, I'll just use Numerical Recipies and built in the other "features" myself. I might look at FFTW eventually. I noticed it said that it "might" be compatible with the VC++ complex class. I forgot there was one, so I looked at it. The += operator for the class only adds the real parts! -- Best Regards, Mike
On Wed, 2 Mar 2005 15:15:24 +0000 (UTC), Mike Yarwood wrote:

> "Randy Yates" <randy.yates@sonyericsson.com> wrote in message > news:xxp1xayb3o4.fsf@usrts005.corpusers.net... >> Active8 <reply2group@ndbbm.net> writes: >> >>> I hate using other people's code >> >> You can try my ancient code, >> >> http://home.earthlink.net/~yatescr/fftk.c >> -- >> Randy Yates >> Sony Ericsson Mobile Communications >> Research Triangle Park, NC, USA >> randy.yates@sonyericsson.com, 919-472-1124 > Hi Active8, > Is you maximum vaue about 5.13?
No. IT's like the one I posted. Here's the first several: 9.328e-041 2.351e-038 -3.087e+010 -1.909e-031 -1.490e+000 8.354e+026 -2.622e-027 1.265e-006 4.799e+010 -2.245e+019 -2.034e+012 4.968e-025 It looks like the author does something to condition the output vector for the screen in the bin axis direction // Get format parameters int center = 128; // 200 //getInput( "Enter center of magnitude", center, 0, ( len / 2 ) - 1 ); int view = 256; // 400 //getInput( "Enter length of magnitude", view, 10, len / 2 ); int begin = center - ( view / 2 ); for( int k = 0; k < view; k++ ) { i = k + begin; if( i < 0 ) i = 0; else if( i >= len / 2 ) i = ( len / 2 ) - 1; vIn[k] = vOut[i].re(); } I think that's why for a norm freq of 1/32, I can get a pip at bin 64 - somewhere around there. The author only runs some wave file for the RFFT routine, he uses a complex fft on the only example with a single sine input - cos, actually. He used 1024 samps of norm freq .05 and ran a 4096 pt fft. He set the display to 400 centered at 200 and got a spectral density centered at 200 for all windowing functions. I can "uncenter" my display and it just looks like it shifts, I get the same hugh numbers. I'm running 512 samps and 512 pt fft. The only windows code in the book is for the plotter and it's a mess of evil MFC code. I don't mess with MFC. I wrote a gui and all that. My scaling and vertical centering of the plot seems to work ok, but I've written the fft values out in text as a sanity check.
> Do you get the right size peak where you > expect it and it's only the bins that you expect to be empty that give you > the nice low numbers?
If I get a correct peak, it's swamped out by the super large ones. The input is 1*sin(2*pi*(1/32)*sampnum) -- Best Regards, Mike
On Wed, 02 Mar 2005 10:28:28 -0500, Jerry Avins wrote:

> Active8 wrote: > >> ... TANSTAAFL > > Heinlein, or hearsay? (The Moon is a Harsh Mistress)
Hearsay. I often find it easier to write my own code. -- Best Regards, Mike
Active8 <reply2group@ndbbm.net> writes:

> On Wed, 02 Mar 2005 10:28:28 -0500, Jerry Avins wrote: > > > Active8 wrote: > > > >> ... TANSTAAFL > > > > Heinlein, or hearsay? (The Moon is a Harsh Mistress) > > Hearsay. I often find it easier to write my own code.
IAOP. -- Randy Yates Sony Ericsson Mobile Communications Research Triangle Park, NC, USA randy.yates@sonyericsson.com, 919-472-1124