!------------------------------------------------- !interface module of fftlib !------------------------------------------------- module fft_interface interface !--------------------------------------------------------- !周期実数値データのフーリエ変換 subroutine rffti(n,wsave) !初期化をおこなう integer, intent(in) :: n !処理するデータの長さ real, intent(out), dimension(*) :: wsave !作業用配列.長さは少なくとも2n+15以上 end subroutine subroutine rfftf(n,r,wsave) !フーリエ順変換をおこなう integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列.長さは少なくとも2n+15以上 real, intent(inout), dimension(*) :: r !処理する実数型配列 end subroutine subroutine rfftb(n,r,wsave) !フーリエ逆変換をおこなう. integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列.長さは少なくとも2n+15以上 real, intent(inout), dimension(*) :: r !処理する実数型配列 end subroutine !--------------------------------------------------------- !rffti, rfftf, rfftbの簡易型サブルーチン subroutine ezffti(n,wsave) !初期化をおこなう; integer, intent(in) :: n !処理するデータの長さ real, intent(out), dimension(*) :: wsave !作業用配列 end subroutine subroutine ezfftf(n,r,a0,a,b,wsave) !フーリエ順変換をおこなう; integer, intent(in) :: n !処理するデータの長さ real, intent(in), dimension(*) :: wsave !作業用配列 real, intent(in), dimension(*) :: r !処理する実数型配列 real, intent(out) :: a0 !上記定義における a_0 および a_0 real, intent(out), dimension(*) :: a !nが偶数のとき n/2, nが奇数のとき(n-1)/2の長さの実数型配列 real, intent(out), dimension(*) :: b ! end subroutine subroutine ezfftb(n,r,a0,a,b,wsave) !ezfftbはフーリエ逆変換をおこなう. integer, intent(in) :: n !処理するデータの長さ real, intent(in), dimension(*) :: wsave !作業用配列 real, intent(out), dimension(*) :: r !処理する実数型配列 real, intent(in) :: a0 !上記定義における a_0 および a_0 real, intent(in), dimension(*) :: a !nが偶数のとき n/2, nが奇数のとき(n-1)/2の長さの実数型配列 real, intent(in), dimension(*) :: b end subroutine !--------------------------------------------------------- !奇の周期データのsine変換をおこなう. subroutine sinti(n,wsave) !初期化をおこなう; integer, intent(in) :: n !処理するデータの長さ real, intent(out), dimension(*) :: wsave !作業用配列 end subroutine subroutine sint(n,r,wsave) !sine変換をおこなう. integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列 real, intent(inout), dimension(*) :: r !処理する実数型配列 end subroutine !--------------------------------------------------------- !偶の周期データのcosine変換をおこなう subroutine costi(n,wsave) !初期化をおこなう integer, intent(in) :: n !処理するデータの長さ real, intent(out), dimension(*) :: wsave !作業用配列 end subroutine subroutine cost(n,r,wsave) !cosine変換をおこなう. integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列 real, intent(inout), dimension(*) :: r !処理する実数型配列 end subroutine !--------------------------------------------------------- !奇数波数成分のみのsin変換をおこなう. subroutine sinqi(n,wsave) !初期化をおこなう integer, intent(in) :: n !処理するデータの長さ real, intent(out), dimension(*) :: wsave !作業用配列 end subroutine subroutine sinqf(n,r,wsave) !sine順変換をおこなう integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列 real, intent(inout), dimension(*) :: r !処理する実数型配列 end subroutine subroutine sinqb(n,r,wsave) !逆変換をおこなう. integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列 real, intent(inout), dimension(*) :: r !処理する実数型配列 end subroutine !--------------------------------------------------------- !偶数波数成分のみのcossine変換をおこなう. subroutine cosqi(n,wsave) !初期化をおこなう integer, intent(in) :: n !処理するデータの長さ real, intent(out), dimension(*) :: wsave !作業用配列 end subroutine subroutine cosqf(n,r,wsave) !cosine順変換をおこなう integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列 real, intent(inout), dimension(*) :: r !処理する実数型配列 end subroutine subroutine cosqb(n,r,wsave) !cosine逆変換をおこなう. integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列 real, intent(inout), dimension(*) :: r !処理する実数型配列 end subroutine !--------------------------------------------------------- !周期複素数データのフーリエ変換をおこなう subroutine cffti(n,wsave) !初期化をおこなう integer, intent(in) :: n !処理するデータの長さ real, intent(out), dimension(*) :: wsave !作業用配列 end subroutine subroutine cfftf(n,c,wsave) !フーリエ順変換をおこなう integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列 complex, intent(inout), dimension(*) :: c !処理する複素数型配列 end subroutine subroutine cfftb(n,c,wsave) !フーリエ逆変換をおこなう. integer, intent(in) :: n !処理するデータの長さ real, intent(inout), dimension(*) :: wsave !作業用配列 complex, intent(inout), dimension(*) :: c !処理する複素数型配列 end subroutine !--------------------------------------------------------- !下位ルーチン subroutine cfftb1(n,c,ch,wa,ifac) integer, intent(in) :: n complex, intent(inout), dimension(*) :: c complex, intent(inout), dimension(*) :: ch real, intent(in), dimension(*) :: wa integer, intent(in), dimension(*) :: ifac end subroutine subroutine cfftf1(n,c,ch,wa,ifac) integer, intent(in) :: n complex, intent(inout), dimension(*) :: c complex, intent(inout), dimension(*) :: ch real, intent(in), dimension(*) :: wa integer, intent(in), dimension(*) :: ifac end subroutine subroutine cffti1(n,wa,ifac) integer, intent(in) :: n real, intent(out), dimension(*) :: wa integer, intent(out), dimension(*) :: ifac end subroutine subroutine cosqb1(n,x,w,xh) integer, intent(in) :: n real, intent(inout), dimension(*) :: x real, intent(in), dimension(*) :: w real, intent(in), dimension(*) :: xh end subroutine subroutine cosqf1(n,x,w,xh) integer, intent(in) :: n real, intent(inout), dimension(*) :: x real, intent(in), dimension(*) :: w real, intent(inout), dimension(*) :: xh end subroutine subroutine ezfft1(n,wa,ifac) integer, intent(in) :: n real, intent(out), dimension(*) :: wa integer, intent(out), dimension(*) :: ifac end subroutine subroutine passb(nac,ido,ip,l1,idl1,cc,c1,c2,ch,ch2,wa) integer, intent(in) :: nac integer, intent(in) :: ido integer, intent(in) :: ip integer, intent(in) :: l1 integer, intent(in) :: idl1 complex, intent(in), dimension(ido,ip,l1) :: cc complex, intent(inout), dimension(ido,l1,ip) :: c1 complex, intent(inout), dimension(idl1,ip) :: c2 complex, intent(inout), dimension(ido,l1,ip) :: ch complex, intent(inout), dimension(idl1,ip) :: ch2 real, intent(in), dimension(*) :: wa end subroutine subroutine passb2(ido,l1,cc,ch,wa1) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,2,l1) :: cc complex, intent(out), dimension(ido,l1,2) :: ch real, intent(in), dimension(*) :: wa1 end subroutine subroutine passb3(ido,l1,cc,ch,wa1,wa2) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,3,l1) :: cc complex, intent(out), dimension(ido,l1,3) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 end subroutine subroutine passb4(ido,l1,cc,ch,wa1,wa2,wa3) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,4,l1) :: cc complex, intent(out), dimension(ido,l1,4) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 real, intent(in), dimension(*) :: wa3 end subroutine subroutine passb5(ido,l1,cc,ch,wa1,wa2,wa3,wa4) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,5,l1) :: cc complex, intent(out), dimension(ido,l1,5) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 real, intent(in), dimension(*) :: wa3 real, intent(in), dimension(*) :: wa4 end subroutine subroutine passf(nac,ido,ip,l1,idl1,cc,c1,c2,ch,ch2,wa) integer, intent(in) :: nac integer, intent(in) :: ido integer, intent(in) :: ip integer, intent(in) :: l1 integer, intent(in) :: idl1 complex, intent(in), dimension(ido,ip,l1) :: cc complex, intent(inout), dimension(ido,l1,ip) :: c1 complex, intent(inout), dimension(idl1,ip) :: c2 complex, intent(inout), dimension(ido,l1,ip) :: ch complex, intent(inout), dimension(idl1,ip) :: ch2 real, intent(in), dimension(*) :: wa end subroutine subroutine passf2(ido,l1,cc,ch,wa1) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,2,l1) :: cc complex, intent(out), dimension(ido,l1,2) :: ch real, intent(in), dimension(*) :: wa1 end subroutine subroutine passf3(ido,l1,cc,ch,wa1,wa2) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,3,l1) :: cc complex, intent(out), dimension(ido,l1,3) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 end subroutine subroutine passf4(ido,l1,cc,ch,wa1,wa2,wa3) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,4,l1) :: cc complex, intent(out), dimension(ido,l1,4) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 real, intent(in), dimension(*) :: wa3 end subroutine subroutine passf5(ido,l1,cc,ch,wa1,wa2,wa3,wa4) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,5,l1) :: cc complex, intent(out), dimension(ido,l1,5) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 real, intent(in), dimension(*) :: wa3 real, intent(in), dimension(*) :: wa4 end subroutine subroutine radb2 (ido,l1,cc,ch,wa1) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,2,l1) :: cc complex, intent(out), dimension(ido,l1,2) :: ch real, intent(in), dimension(*) :: wa1 end subroutine subroutine radb3(ido,l1,cc,ch,wa1,wa2) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,3,l1) :: cc complex, intent(out), dimension(ido,l1,3) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 end subroutine subroutine radb4(ido,l1,cc,ch,wa1,wa2,wa3) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,4,l1) :: cc complex, intent(out), dimension(ido,l1,4) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 real, intent(in), dimension(*) :: wa3 end subroutine subroutine radb5(ido,l1,cc,ch,wa1,wa2,wa3,wa4) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,5,l1) :: cc complex, intent(out), dimension(ido,l1,5) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 real, intent(in), dimension(*) :: wa3 real, intent(in), dimension(*) :: wa4 end subroutine subroutine radbg (ido,ip,l1,idl1,cc,c1,c2,ch,ch2,wa) integer, intent(in) :: ido integer, intent(in) :: ip integer, intent(in) :: l1 integer, intent(in) :: idl1 complex, intent(in), dimension(ido,ip,l1) :: cc complex, intent(inout), dimension(ido,l1,ip) :: c1 complex, intent(inout), dimension(idl1,ip) :: c2 complex, intent(inout), dimension(ido,l1,ip) :: ch complex, intent(inout), dimension(idl1,ip) :: ch2 real, intent(in), dimension(*) :: wa end subroutine subroutine radf2 (ido,l1,cc,ch,wa1) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,2,l1) :: cc complex, intent(out), dimension(ido,l1,2) :: ch real, intent(in), dimension(*) :: wa1 end subroutine subroutine radf3(ido,l1,cc,ch,wa1,wa2) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,3,l1) :: cc complex, intent(out), dimension(ido,l1,3) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 end subroutine subroutine radf4(ido,l1,cc,ch,wa1,wa2,wa3) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,4,l1) :: cc complex, intent(out), dimension(ido,l1,4) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 real, intent(in), dimension(*) :: wa3 end subroutine subroutine radf5(ido,l1,cc,ch,wa1,wa2,wa3,wa4) integer, intent(in) :: ido integer, intent(in) :: l1 complex, intent(in), dimension(ido,5,l1) :: cc complex, intent(out), dimension(ido,l1,5) :: ch real, intent(in), dimension(*) :: wa1 real, intent(in), dimension(*) :: wa2 real, intent(in), dimension(*) :: wa3 real, intent(in), dimension(*) :: wa4 end subroutine subroutine radfg(ido,ip,l1,idl1,cc,c1,c2,ch,ch2,wa) integer, intent(in) :: ido integer, intent(in) :: ip integer, intent(in) :: l1 integer, intent(in) :: idl1 complex, intent(in), dimension(ido,ip,l1) :: cc complex, intent(inout), dimension(ido,l1,ip) :: c1 complex, intent(inout), dimension(idl1,ip) :: c2 complex, intent(inout), dimension(ido,l1,ip) :: ch complex, intent(inout), dimension(idl1,ip) :: ch2 real, intent(in), dimension(*) :: wa end subroutine subroutine rfftb1(n,c,ch,wa,ifac) integer, intent(in) :: n complex, intent(inout), dimension(*) :: c complex, intent(inout), dimension(*) :: ch real, intent(in), dimension(*) :: wa integer, intent(in), dimension(*) :: ifac end subroutine subroutine rfftf1(n,c,ch,wa,ifac) integer, intent(in) :: n complex, intent(inout), dimension(*) :: c complex, intent(inout), dimension(*) :: ch real, intent(in), dimension(*) :: wa integer, intent(in), dimension(*) :: ifac end subroutine subroutine rffti1(n,wa,ifac) integer, intent(in) :: n real, intent(out), dimension(*) :: wa integer, intent(out), dimension(*) :: ifac end subroutine subroutine sint1(n,war,was,xh,x,ifac) integer, intent(in) :: n real, intent(inout), dimension(*) :: war real, intent(in), dimension(*) :: was real, intent(out), dimension(*) :: xh real, intent(inout), dimension(*) :: x integer, intent(in), dimension(*) :: ifac end subroutine end interface end module !fftlib library end ----