function spectra_2d,data,yfac=yfac,nofold=nofold,taper_=taper_

; Purpose: compute 2-dimensional fft (usually: x=longitude, y=time)
; Category: numerical statistics
; Output: anonymous structure {spec,w,f}
;	where:
;	spec=array (nx,ny/2+1) where:
;  	  y(0)=zero-freq; y(1)=lowest-freq; y(ny/2)=highest-freq
;	  x(0)=fastest-w-prop; x(nx/2-1)=slowest-w-prop; x(nx/2)=stationary; x(nx/2+1)=slowest-e
;	w=array(nx) of wavenumbers from -nx/2+1 to nx/2
;	f=array(ny/2+1) of frequencies from 0 to 0.5
; See-also: spectra_2d_plot
; Keywords:
;	yfac=number - multiply f by this; eg if y is months and you want
;	              cycles/year, yfac=12.

@comm_error

; Check the data dimensions are even (not a real requirement but simplifies)
s=size(data)
if (s(0) ne 2) then message,'Require 2-d array'
if ((s(1) mod 2 ne 0) or (s(2) mod 2 ne 0)) then message,'Require both dimensions to be even'

; Keywords
if (n_elements(yfac) eq 0) then yfac=1

; Define some useful numbers
m=s(1)
m1=m/2
n=s(2)
n1=n/2

; Taper as required
data0=data
if (keyword_set(taper_)) then $
  for i=0,m-1 do data0(i,*)=taper(data(i,*),per=taper_)

; Compute fft
data1=abs(fft(data0))

; Make so that x-direction is highest-w -> lowest-w -> zero -> lowest-e -> highest-e
data2=data1; data2(*,*)=-999
data2(m1-1+indgen(m1+1),*)=data1(indgen(m1+1),*)
data2(indgen(m1-1),*)=data1(m1+1+indgen(m1-1),*)
i=where(data2 eq -999,count) & if (count gt 0) then message,'Internal error'

; Fold over
if (not keyword_set(nofold)) then begin
  data3=data2(*,indgen(n1+1))
  for i=1,n1-1 do data3(*,i)=data2(*,i)+shift(reverse(data2(*,n-i)),-1)
  f=indgen(n1+1)/float(n)*yfac
endif else begin
  data3=data2
  f=indgen(n)/float(n)*yfac
endelse

; Return
return,{spec:data3,w:indgen(m)-m1+1,f:f}
end