w4kek/WSJT-X
1
0
Fork 0
mirror of https://github.com/saitohirga/WSJT-X.git synced 2024-11-12 07:06:17 -05:00
Code Issues Packages Projects Releases Wiki Activity

Test program

Browse Source 
git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@6607 ab8295b8-cf94-4d9e-aec4-7959e3be5d79
This commit is contained in:
Joe Taylor 2016-04-11 15:24:55 +00:00
parent 03d74a6208
commit 3541d12a67
1 changed files with 30 additions and 29 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

59
lib/t3.f90
View File

@ -19,30 +19,32 @@ program t3
call filter(x0(ia:ib),x1(ia:ib)) call filter(x0(ia:ib),x1(ia:ib))
enddo enddo
do i=1,NZ x1(1:NZ-NBLK)=x1(NBLK+1:NZ)
do i=1,NZ-NBLK
write(13,1001) i,x0(i),x1(i),x1(i)-x0(i) write(13,1001) i,x0(i),x1(i),x1(i)-x0(i)
1001 format(i6,3e12.3) 1001 format(i6,3f13.9)
enddo enddo
end program t3 end program t3
subroutine filter(x0,x1) subroutine filter(x0,x1)
! Process time-domain data sequentially, optionally using 'refspec.dat' ! Process time-domain data sequentially, optionally using a frequency-domain
! to flatten the spectrum. ! filter to alter the spectrum.
! NB: sin^2 window with 50% overlap; sin^2 + cos^2 = 1.0. ! NB: uses a sin^2 window with 50% overlap.
parameter (NFFT=6912,NH=NFFT/2) parameter (NFFT=6912,NH=NFFT/2)
real x0(0:NH-1) !Real input data real x0(0:NH-1) !Input samples
real x1(0:NH-1) !Real output data real x1(0:NH-1) !Output samples (delayed by one block)
real xov(0:NH-1) real x0s(0:NH-1) !Saved upper half of input samples
real x1s(0:NH-1) !Saved upper half of output samples
real x(0:NFFT-1) real x(0:NFFT-1) !Work array
complex cx(0:NH) real*4 w(0:NFFT-1) !Window function
real*4 w(0:NFFT-1) real f(0:NH) !Filter to be applied
real*4 s(0:NH) real*4 s(0:NH) !Average spectrum
logical first logical first
complex cx(0:NH) !Complex frequency-domain work array
equivalence (x,cx) equivalence (x,cx)
data first/.true./ data first/.true./
save save
@ -50,26 +52,25 @@ subroutine filter(x0,x1)
if(first) then if(first) then
pi=4.0*atan(1.0) pi=4.0*atan(1.0)
do i=0,NFFT-1 do i=0,NFFT-1
w(i)=(sin(i*pi/NFFT))**2 ww=sin(i*pi/NFFT)
w(i)=ww*ww/NFFT
enddo enddo
s=0. s=0.0
fac=1.0/NFFT f=1.0
x0s=0.
x1s=0.
first=.false. first=.false.
xov=0.
endif endif
x(:NH-1)=xov !Previous 2nd half to new 1st half x(0:NH-1)=x0s !Previous 2nd half to new 1st half
x(NH:)=x0 !New 2nd half x(NH:NFFT-1)=x0 !New 2nd half
x=x*w !Apply window x0s=x0 !Save the new 2nd half
call four2a(x,NFFT,1,-1,0) !r2c FFT: to freq domain x=w*x !Apply window
call four2a(x,NFFT,1,-1,0) !r2c FFT (to frequency domain)
! Apply filter to cx() cx=f*cx
call four2a(cx,NFFT,1,1,-1) !c2r FFT (back to time domain)
call four2a(cx,NFFT,1,1,-1) !c2r FFT: back to time domain x1=x1s + x(0:NH-1) !Add previous segment's 2nd half
x1s=x(NH:NFFT-1) !Save the new 2nd half
x(0:NH-1)=x(0:NH-1)+xov(0:NH-1) !Add previous segment's 2nd half
x1=x
return return
end subroutine filter end subroutine filter