subroutine fft(s,n,output)
	complex f(8192)
	dimension s(1),output(1)
	call ftrans(s,f,n)
	nd2=n/2
	do 10 j=1,nd2
10	output(j)=cabs(f(j))
	return	
	end
	subroutine ftrans(s,f,n)
	dimension s(1)
	complex f(1)
	call shuff(s,f,n)
	length=2
1	do 2 j=1,n,length
2	call combin(f,j,length)
	length=length+length
	if(length-n) 1,1,3
3	return
	end
	subroutine shuff(s,f,n)
	dimension s(1)
	complex f(1),cmplx
	do 5 ifort=1,n
	i=ifort-1
	j=0
	m2=1
1	m1=m2
	m2=m2+m2
	if(mod(i,m2)-m1) 3,2,2
2	j=j+n/m2
3	if(m2-n) 1,4,4
4	jfort=j+1
5	f(ifort)=cmplx(s(jfort),0.)
	return
	end
	subroutine combin(f,j,n)
	complex f(1),emjt,z,cexp
	emjt=cexp((0.,-1.)*(6.283185/float(n)))
	n2=n/2
	do 1 l=1,n2
	loc1=l+j-1
	loc2=loc1+n2
	z=emjt**(l-1)*f(loc2)
	f(loc2)=f(loc1)-z
1	f(loc1)=f(loc1)+z
	return
	end