Messages from Usenet group comp.lang.idl-pvwave, compiled by Paulo Penteado
In article <1193671433.271745.268030@57g2000hsv.googlegroups.com>, knielsen73@gmail.com wrote: > Hi, > > I have a time vs longitude plot and would like to do a 2-D FFT to get > a frequency vs wavenumber plot. I used the IDL FFT but I am not sure > what the output is...I didn't feel the IDL help was very helpful > there. I found a post from 1998 but there were some non-defined > functions and I wasn't exactly sure what they did. > > Does anybody have a good example on how to use 2D FFT in IDL or some > good documentation? > > Cheers, > Kim The main issues with the IDL FFT (or any FFT for that matter) are sorting out the order in which the coefficients are stored and understanding how the sine and cosine components are represented in a complex FFT. My recommendation is that you create an artificial, low-resolution, input data set that contains some simple combination of sines and cosines and that you understand completely. Then look at the FFT. Or you can buy my book (An Introduction to Programming with IDL), which has a chapter that does exactly that for the 1-D case. Ken Bowman
<knielsen73@gmail.com> wrote in message news:1193671433.271745.268030@57g2000hsv.googlegroups.com... > Hi, > > I have a time vs longitude plot and would like to do a 2-D FFT to get > a frequency vs wavenumber plot. I used the IDL FFT but I am not sure > what the output is...I didn't feel the IDL help was very helpful > there. I found a post from 1998 but there were some non-defined > functions and I wasn't exactly sure what they did. > > Does anybody have a good example on how to use 2D FFT in IDL or some > good documentation? a little example:nx = 100 ny = 200 x = (fltarr(ny)+1) ## findgen(nx) y = findgen(ny) ## (fltarr(nx)+1) kx = 20.5/nx ky = 0./ny ;image a =cos(2*!pi*(kx*x + ky*y)) ;spectrum in 2d spe = fft(a) kxvalues = findgen(nx)/nx - 1./2 kyvalues = findgen(ny)/ny - 1./2 !p.multi=[0,1,2] contour,a,/fill,nlevels=14 ; shift to get the normal representation, ; i.e. (neg kx, 0, pos kx) shiftspe = shift(spe,nx/2,ny/2) contour,abs(shiftspe),kxvalues,kyvalues,/fill,nlevels=14 Note: it is symmetric about the origin. Cheers, bob
On Oct 29, 2:41 pm, "R.G. Stockwell" <noem...@please.com> wrote: > <knielse...@gmail.com> wrote in message > > news:1193671433.271745.268030@57g2000hsv.googlegroups.com... > >> Hi, > >> I have a time vs longitude plot and would like to do a 2-D FFT to get >> a frequency vs wavenumber plot. I used the IDL FFT but I am not sure >> what the output is...I didn't feel the IDL help was very helpful >> there. I found a post from 1998 but there were some non-defined >> functions and I wasn't exactly sure what they did. > >> Does anybody have a good example on how to use 2D FFT in IDL or some >> good documentation? > > a little example:nx = 100 > ny = 200 > x = (fltarr(ny)+1) ## findgen(nx) > y = findgen(ny) ## (fltarr(nx)+1) > kx = 20.5/nx > ky = 0./ny > > ;image > a =cos(2*!pi*(kx*x + ky*y)) > > ;spectrum in 2d > spe = fft(a) > > kxvalues = findgen(nx)/nx - 1./2 > kyvalues = findgen(ny)/ny - 1./2 > > !p.multi=[0,1,2] > contour,a,/fill,nlevels=14 > > ; shift to get the normal representation, > ; i.e. (neg kx, 0, pos kx) > shiftspe = shift(spe,nx/2,ny/2) > > contour,abs(shiftspe),kxvalues,kyvalues,/fill,nlevels=14 > > Note: it is symmetric about the origin. > > Cheers, > bob Thanks Bob, I works great with your example. However, when I made a wave with a wavenumber larger than 1, the scale no longer agrees. I can plot it "raw", but again, not sure how to read the values with respect to the real wave number and frequency. Here is the wave I gave as input. ;input a wave number two with an amplitude of 10. ;and eastward phase propagation of 1.7 days amp = 10. wave_number = 2 wavelength = float(jmax)/wave_number period = 6.8333333 for j=0,jmax-1 do begin for n=0,nmax-1 do begin data_in[j,n] = amp*sin(2.*!pi*j/wavelength - 2.*!pi*n/period) +$ amp*.5*cos(2.*!pi*j/wavelength - 2.*!pi*n/period) endfor endfor Cheers, Kim
<knielsen73@gmail.com> wrote in message news:1193692469.369318.51960@50g2000hsm.googlegroups.com... >> kxvalues = findgen(nx)/nx - 1./2 >> kyvalues = findgen(ny)/ny - 1./2 Note, in the above, it shows that the sampling interval is 1. You should divide by your sampling interval there if it is not 1. > amp = 10. > wave_number = 2 > wavelength = float(jmax)/wave_number > period = 6.8333333 > > for j=0,jmax-1 do begin > for n=0,nmax-1 do begin > data_in[j,n] = amp*sin(2.*!pi*j/wavelength - 2.*!pi*n/period) +$ > amp*.5*cos(2.*!pi*j/wavelength - 2.*!pi*n/period) for one thing, you are adding 2 out of phase waves here. Remove the first term (and just use the cos). Also, do just 1D to ensure that you are calculating the correct period and freq. f = 12/len data = amp*cos(2 !pi f t)
"R.G. Stockwell" <noemail@please.com> wrote in message news:47265628$0$497$815e3792@news.qwest.net... > > Also, do just 1D to ensure that you are calculating the correct > period and freq. > > f = 12/len > data = amp*cos(2 !pi f t) oops sorry about that, i've got a 1 yr old climbing on the laptop. Anyways, do the example in 1D first with a simple cos, to ensure that you are calculating the correct freqs. Cheers, bob
On Oct 29, 5:54 pm, "R.G. Stockwell" <noem...@please.com> wrote: > "R.G. Stockwell" <noem...@please.com> wrote in message > > news:47265628$0$497$815e3792@news.qwest.net... > > > >> Also, do just 1D to ensure that you are calculating the correct >> period and freq. > >> f = 12/len >> data = amp*cos(2 !pi f t) > > oops sorry about that, i've got a 1 yr old climbing on the > laptop. > > Anyways, do the example in 1D first with a simple cos, > to ensure that you are calculating the correct freqs. > > Cheers, > bob thanks Bob...got it!