I have questions about the spctrm function that calculates the Fourier Power Spectrum.

I have 4096 real data points (actually 640 points with zeroes filling up to 4096). I expected that I would get 2049 (or 2048) data points from the FFT power spectrum. However, without overlap, and k=1, I need to set m=1024 so that there is enough data for the function. Therefore I only get 1024 output points. What am I missing?

Also, the text book says that it returns M frequencies, not M+1 but that to include the P(fc) at the Nyquist freqency would be "straightforward". My math is a little rusty and it's not clear to me how to do this.

Thanks.

Janet

I have read the power spectrum estimation using the fft in the chapter 13. I use the subroutine spctrm to clculate the power spectrum of Sin function, and I get the results having three peaks. It is wrong apparently. I donn't know why. In the calculation, I choose the parameters: m=1024 , k=10 and ovrlap=.true.
I have seen your post and want to know whether you understand the subtoutine clearly. Can you tell me.
Thank you sincerely!

"I have seen your post and want to know whether you understand the subtoutine clearly. Can you tell me."

Sorry, I cannot say that I understand the subroutine. I never got any replies and have not looked at the function since my post. Good luck to you.

Can you find a new method to calculate power Spectrum

Hi,

I'm working with spctrm.c routine also, and as I see, we are having the same problems, that is when setting the parameter m=2048 (so N=2*m=4096) a "one-sided" PSD of 2048 points is expected but I obtain a "two-sided" PSD of 2048 (what means a 1024 "one-sided").

I've been searching in google groups and I found that:

Hi all,


Does anyone knows of bugs in 'spctrm' routine of Numerical Recipes 2nd Ed.?
It says that it uses 2 partitions of 2m REAL data points per segment, then
it reads consecutive data points (real?) into odd positions, and then into even positions.
Then it uses 'four1' routine which is supposed to use COMPLEX data points.
Something doesn't make sense to me here.
Anyone had experience with this routine?


Thanks for any help,


Yair

And the reply:

I have used these routines in a number of programs which have been in use
for years; these routines are trustworthy.


Having said that, the technique of sending a real vector of length 2N to a
routine expecting a complex vector of length N takes some getting used to,
although such use is very common in the world of FFT Fortran code. Read
the book very carefully, and more than once. There are a number of
factors of two in there to keep straight. Also note the interpretation of
the ordering of the returned data; it may not be what you expect.

I don't know what really coulc happend, but my theory is that spctrm returns PSD "two-sided" (althougt in the pdf says it returns p[j] j=1...m

Uri

I forgot to especialy mention the last line (the underlined one):

"Also note the interpretation of
the ordering of the returned data; it may not be what you expect"

It seams like We are not interpreting the result correctly.

Uri