用FFT实现相关函数快速估计
主程序:c----------------------------------------------------------------------
c Main program HCORRE2 : to test subroutine CORRE2
c To compute the auto-correlation function of a simple data by FFT
C Please link subroutine CORRE2,SPLFFT
c----------------------------------------------------------------------
complex x(0:15)
data m/8/,n/16/,icorre/0/
do 10 i=0,m-1
x(i)=i+1.
10 continue
call corre2(x,x,m,n,icorre,ierror)
write(*,*)' ierror=',ierror
if(ierror.ne.0) stop
write(*,*)(x(i),i=0,n-1)
stop
end
子程序:
subroutine corre2(x,y,m,n,icorre,ierror)
c----------------------------------------------------------------------
cRoutinue CORRE2: To Compute the CorrelationFunction of x(i) and
cy(i) by DFT. x(i),y(i),i=0,...,m-1;But the dimension n of x,y must
cbe >=2*m and be the power of 2 ;
cIf: icorre=0:For auto-correlation
c icorre=1:For crros-correlation
c in chapter 11
c----------------------------------------------------------------------
complex x(0:n-1),y(0:n-1)
ierror=1
s=n
1 s=s/2.
if(s-2.)70,5,1
5 do 10 i=m,n-1
10 x(i)=0.
call splfft(x,n,-1)
if(icorre.eq.1)goto 30
do 20 k=0,n-1
x(k)=conjg(x(k))*x(k)/float(m)
20 continue
goto 60
30 do 40 i=m,n-1
40 y(i)=0.
call splfft(y,n,-1)
do 50 k=0,n-1
x(k)=conjg(x(k))*y(k)/float(m)
50 continue
60 call splfft(x,n,1)
ierror=0
70 return
end
页:
[1]