lipman@fciads.bsd.uchicago.edu (Everett Lipman) (07/23/90)
I am looking for source code for a 2 dimensional Fast Fourier Transform routine. If anyone has one or knows where to get one, please send me e-mail. C is preferred, but any language is fine. Thanks. Everett Lipman (lipman@fciads.bsd.uchicago.edu)
3ksnn64@cadlab1.ecn.purdue.edu (Joe Cychosz) (07/27/90)
In article <1990Jul22.175529.8882@midway.uchicago.edu> lipman@fciads.bsd.uchicago.edu (Everett Lipman) writes: >I am looking for source code for a 2 dimensional Fast Fourier >Transform routine. If anyone has one or knows where to get >one, please send me e-mail. C is preferred, but any language >is fine. Thanks. > > Everett Lipman (lipman@fciads.bsd.uchicago.edu) c Perform two dimensional fast fourier transformation using c routines mfft and fft. c c data = array npix X npix complex numbers c npix = size of image c powof2 = npix must be a power of 2 c transform each column (i.e. along x) CMIC$ PARALLEL SHARED (data,npix,powof2) CMIC$ DO PARALLEL do 4 y = 1,npix 4 call fft (data(1,y),npix,powof2) CMIC$ END DO CMIC$ END PARALLEL c transform each row (i.e. along y) call mfft (data,npix,npix,powof2,cwork) subroutine fft (x,n,ln) c c ------------------------------------------------------------ c c fft - fast fourier transform. c c this subroutines perfroms a fast fourier transform on a c complex vector x. c c on entry: c x complex an array containing the complex vector c to be transformed. c n integer number of complex elements in the c vector. 'n' must be a power of 2. c ln integer log base 2 of n. (i.e. 2**ln = n) c c origin: c c ------------------------------------------------------------ c complex x(n), u, w, t data pi /3.1415926535898/ nv2 = n / 2 nm1 = n - 1 c perform a bit reversal permutation on the x array. j = 1 do 3 i = 1,nm1 if (i .ge. j) goto 1 t = x(j) x(j) = x(i) x(i) = t 1 k = nv2 2 if (k .ge. j) goto 3 j = j - k k = k / 2 goto 2 3 j = j + k do 5 l = 1,ln le = 2**l le1 = le/2 u = (1.0,0.0) w = cmplx (cos (pi/le1),-sin (pi/le1)) do 5 j = 1,le1 do 4 i = j,n,le ip = i + le1 t = x(ip) * u x(ip) = x(i) - t 4 x(i) = x(i) + t 5 u = u * w do 6 i = 1,n 6 x(i) = x(i) / float(n) return end subroutine mfft (x,m,n,ln,t) c c ------------------------------------------------------------ c c mfft - multiple fast fourier transform. c c this subroutines perfroms multiple fast fourier transforms c on a complex array x. c c on entry: c x complex an array of dimensions m x n containing c the complex vectors to be transformed. c m integer number of vectors to transform. c 'm' must not exceed 65535. c n integer number of complex elements in each c vector. 'n' must be a power of 2. c ln integer log base 2 of n. (i.e. 2**ln = n) c t complex a scratch array of size 'm'. c c origin: c c ------------------------------------------------------------ c complex x(m,n), t(m), u, w data pi /3.1415926535898/ nv2 = n / 2 nm1 = n - 1 j = 1 do 3 i = 1,nm1 if (i .ge. j) goto 1 do 11 ll = 1,m t(ll) = x(ll,j) x(ll,j) = x(ll,i) 11 x(ll,i) = t(ll) c t(1;m) = x(1,j;m) c x(1,j;m) = x(1,i;m) c x(1,i;m) = t(1;m) 1 k = nv2 2 if (k .ge. j) goto 3 j = j - k k = k / 2 goto 2 3 j = j + k do 5 l = 1,ln le = 2**l le1 = le/2 u = (1.0,0.0) w = cmplx (cos (pi/le1),-sin (pi/le1)) do 5 j = 1,le1 do 4 i = j,n,le ip = i + le1 do 4 ll = 1,m t(ll) = x(ll,ip) * u x(ll,ip) = x(ll,i) - t(ll) 4 x(ll,i) = x(ll,i) + t(ll) c t(1;m) = x(1,ip;m) * u c x(1,ip;m) = x(1,i;m) - t(1;m) c 4 x(1,i;m) = x(1,i;m) + t(1;m) 5 u = u * w do 6 i = 1,n do 6 l = 1,m 6 x(l,i) = x(l,i) / float(n) return end