diff --git a/lib/superfox/qpc_decode2.f90 b/lib/superfox/qpc_decode2.f90 index 2a260e021..7900cb7ee 100644 --- a/lib/superfox/qpc_decode2.f90 +++ b/lib/superfox/qpc_decode2.f90 @@ -149,16 +149,17 @@ subroutine smo121a(x,nz,a,b) return end subroutine smo121a -subroutine remove_tone(c0,npts) +subroutine remove_tone(c0,fsync) + parameter (NMAX=15*12000) parameter (NFILT=8000) - complex c0(npts) + complex c0(NMAX) complex cwindow(15*12000) - complex cref(npts) - complex cfilt(npts) + complex cref(NMAX) + complex cfilt(NMAX) real window(-NFILT/2:NFILT/2) ! real endcorrection(NFILT/2+1) - real s(npts/4) + real s(NMAX/4) integer ipk(1) logical first data first/.true./ @@ -166,7 +167,7 @@ subroutine remove_tone(c0,npts) if(first) then pi=4.0*atan(1.0) - fac=1.0/float(npts) + fac=1.0/float(NMAX) sumw=0.0 do j=-NFILT/2,NFILT/2 window(j)=cos(pi*j/NFILT)**2 @@ -175,37 +176,60 @@ subroutine remove_tone(c0,npts) cwindow=0. cwindow(1:NFILT+1)=window/sumw cwindow=cshift(cwindow,NFILT/2+1) - call four2a(cwindow,npts,1,-1,1) + call four2a(cwindow,NMAX,1,-1,1) cwindow=cwindow*fac ! frequency domain smoothing filter first=.false. endif fsample=12000.0 - df=fsample/npts - fac=1.0/npts + df=fsample/NMAX + fac=1.0/NMAX cfilt=fac*c0 - call four2a(cfilt,npts,1,-1,1) ! fourier transform of input data - iz=npts/4 + call four2a(cfilt,NMAX,1,-1,1) ! fourier transform of input data + iz=NMAX/4 do i=1,iz s(i)=real(cfilt(i))**2 + aimag(cfilt(i))**2 enddo - ia=nint(700.0/df) - ib=nint(800.0/df) + ia=nint((fsync-100.0)/df) + ib=nint((fsync+100.0)/df) ipk=maxloc(s(ia:ib)) i0=ipk(1) + ia - 1 - f2=df*i0 - write(*,*) 'remove_tone - frequency: ',f2 + + i10=nint(10.0/df) + ia=i0-i10 + ib=i0+i10 + s0=0.0 + s1=0.0 + s2=0.0 + do i=ia,ib + s0=s0+s(i) + s1=s1+(i-i0)*s(i) + enddo + delta=s1/s0 + i0=nint(i0+delta) + f2=i0*df + + ia=i0-i10 + ib=i0+i10 + do i=ia,ib + s2=s2 + s(i)*(i-i0)**2 + enddo + sigma=sqrt(s2/s0)*df + +! write(61,*) 'frequency, spectral width ',f2,sigma + if(sigma .gt. 2.5) go to 999 +! write(61,*) 'remove_tone - frequency: ',f2 dt=1.0/fsample - do i=1, npts + do i=1, NMAX arg=2*pi*f2*i*dt cref(i)=cmplx(cos(arg),sin(arg)) enddo cfilt=c0*conjg(cref) ! baseband to be filtered - call four2a(cfilt,npts,1,-1,1) + call four2a(cfilt,NMAX,1,-1,1) cfilt=cfilt*cwindow - call four2a(cfilt,npts,1,1,1) + call four2a(cfilt,NMAX,1,1,1) nframe=50*3456 do i=1,nframe @@ -213,6 +237,6 @@ subroutine remove_tone(c0,npts) c0(i)=c0(i)-cref(i) enddo - return +999 return end subroutine remove_tone diff --git a/lib/superfox/sfrx_sub.f90 b/lib/superfox/sfrx_sub.f90 index b0e1dcb92..0779b23ef 100644 --- a/lib/superfox/sfrx_sub.f90 +++ b/lib/superfox/sfrx_sub.f90 @@ -36,6 +36,8 @@ subroutine sfrx_sub(nyymmdd,nutc,nfqso,ntol,iwave) call sfox_ana(dd,npts,c0,npts) +! call remove_tone(c0,fsync) ! Needs testing + ndepth=3 dth=0.5 damp=1.0