|
From: glowkeeper on 24 Apr 2008 12:32 The planning function I'm using is fftwf_plan_r2r_1d, and I use it this way: fftwf_plan_r2r_1d( 2048, m_pFourierIn, m_pFourierOut, FFTW_R2HC, FFTW_BACKWARD ); But nonetheless, I'm still getting the power algorithm wrong it seems.... >>I suppose I could also do something like taking the square root of the >>square, but that seems expensive (although granted, I don't know what >fabs >>is doing under the hood). > >First of all, your Fourier coefficients should be complex, even if the >input data are real. From FFTW's website, the planning function is: > >fftw_plan fftw_plan_dft_r2c_1d(int n, double *in, fftw_complex *out, > unsigned flags); > >where *in is a pointer to your real input array, and *out is a pointer to >your complex output array. > >Next, fabs() performs an absolute value of a real data point. For power, >you actually want the absolute value squared. I.e., given a complex >coefficient of (a+jb), the power is (a+jb)(a-jb), or abs(a+jb)^2. Since >you're using FFTW, the output array is defined as type fftw_complex (or >fftwf_complex for single-precision), unless you've included the complex.h >header from the C99 math library. If so, then it will be of type >double/float complex at the output. I've never done it this way but maybe >Steven Johnson can chime in on how it works better (I prefer their >convention since the complex math doesn't work as well on the processors >I've used). Either way, the appropriate function at that point would be >cabs(x)*cabs(x) (well, y = cabs(x), then y2 = y*y, if you want to out-think >the compiler), which is the complex magnitude squared. > >Mark >
From: glowkeeper on 24 Apr 2008 12:38 Actually, I should say that I'm using fftwf_plan_r2r_1d( 2048, m_pFourierIn, m_pFourierOut, FFTW_R2HC, FFTW_BACKWARD ); and then when I come to do the power computations I work on the first 1024 entries in the returned transform (ignoring the imaginary parts). Is this incorrect? >The planning function I'm using is fftwf_plan_r2r_1d, and I use it this >way: > >fftwf_plan_r2r_1d( 2048, m_pFourierIn, m_pFourierOut, FFTW_R2HC, >FFTW_BACKWARD ); > >But nonetheless, I'm still getting the power algorithm wrong it >seems.... > >>>I suppose I could also do something like taking the square root of the >>>square, but that seems expensive (although granted, I don't know what >>fabs >>>is doing under the hood). >> >>First of all, your Fourier coefficients should be complex, even if the >>input data are real. From FFTW's website, the planning function is: >> >>fftw_plan fftw_plan_dft_r2c_1d(int n, double *in, fftw_complex *out, >> unsigned flags); >> >>where *in is a pointer to your real input array, and *out is a pointer >to >>your complex output array. >> >>Next, fabs() performs an absolute value of a real data point. For >power, >>you actually want the absolute value squared. I.e., given a complex >>coefficient of (a+jb), the power is (a+jb)(a-jb), or abs(a+jb)^2. Since >>you're using FFTW, the output array is defined as type fftw_complex (or >>fftwf_complex for single-precision), unless you've included the >complex.h >>header from the C99 math library. If so, then it will be of type >>double/float complex at the output. I've never done it this way but >maybe >>Steven Johnson can chime in on how it works better (I prefer their >>convention since the complex math doesn't work as well on the processors >>I've used). Either way, the appropriate function at that point would be >>cabs(x)*cabs(x) (well, y = cabs(x), then y2 = y*y, if you want to >out-think >>the compiler), which is the complex magnitude squared. >> >>Mark >> >
From: markt on 24 Apr 2008 12:59 >Actually, I should say that I'm using > >fftwf_plan_r2r_1d( 2048, m_pFourierIn, m_pFourierOut, FFTW_R2HC, >FFTW_BACKWARD ); Ah, gotcha. OK, the way I read the documentation on this is that the complex data still exists, but in the second half of the transform output, which you would need to incorporate into your absolute value function, i.e. you still get a complex output, but the data are ordered all of the inphase first, then quadrature last. I'll defer to SJ for a better explanation since I've never used the R2HC option. Mark
From: Ron N. on 24 Apr 2008 14:15 On Apr 24, 8:38 am, "glowkeeper" <glowkee...(a)yahoo.co.uk> wrote: > Actually, I should say that I'm using > > fftwf_plan_r2r_1d( 2048, m_pFourierIn, m_pFourierOut, FFTW_R2HC, > FFTW_BACKWARD ); > > and then when I come to do the power computations I work on the first 1024 > entries in the returned transform (ignoring the imaginary parts). Is this > incorrect? Yes. If you ignore the "imaginary" parts, then a phase shift of 90 degrees could cause the fundamental or any harmonic to disappear from your composite "spectrum". If your fft frames are not phase locked to the signal generator, then this could easily cause the random type of results you may be seeing.
From: Ron N. on 24 Apr 2008 14:24 On Apr 23, 11:33 pm, "glowkeeper" <glowkee...(a)yahoo.co.uk> wrote: > >If you downsample your fft, shouldn't that change the bin > >width? Or are you throwing away spectral energy? That > >loss of information could certainly effect the results. > > Ah, now this was something I wasn't sure of to be honest. Maybe? I output > my product to an array the same size as the fft, and therefore I thought > the bin width was the same. > > But I think you're suggestingt that my 2x bin width is ~5Hz/2, 3x ~5Hz/3, > 4x ~5Hz/4, ... Actually, the opposite. When you downsample a spectrum, each bin should end up representing a wider frequency band, e.g. 5 Hz * 2, etc., or you will be throwing away energy. IMHO. YMMV.
First
|
Prev
|
Pages: 1 2 Prev: Zeroforcing precoding/beamforming Next: Retry: Java Audio question - PLEASE HELP |