From e4c5811ac365aeafd482c6d869e86578c281bc36 Mon Sep 17 00:00:00 2001 From: Steven Franke Date: Tue, 20 Dec 2016 03:01:43 +0000 Subject: [PATCH] Changes needed to accommodate channel equalization implemented in the new msk144signalquality.f90. Warning - msk144signalquality writes to the executable directory. git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@7412 ab8295b8-cf94-4d9e-aec4-7959e3be5d79 --- lib/analytic.f90 | 47 +++++++++++++++++++++------------------ lib/msk144decodeframe.f90 | 2 +- lib/msk144spd.f90 | 5 +++-- lib/mskrtd.f90 | 36 ++++++++++++++++++++---------- 4 files changed, 53 insertions(+), 37 deletions(-) diff --git a/lib/analytic.f90 b/lib/analytic.f90 index beb7e6f6f..a4f5ef61f 100644 --- a/lib/analytic.f90 +++ b/lib/analytic.f90 @@ -1,29 +1,44 @@ -subroutine analytic(d,npts,nfft,c) +subroutine analytic(d,npts,nfft,c,a,equalize) ! Convert real data to analytic signal parameter (NFFTMAX=1024*1024) real d(npts) - real h(NFFTMAX/2) + complex h(NFFTMAX/2) + real*8 a(5),alast(5),fp + real ac(5),a0(5) complex c(NFFTMAX) + logical*1 equalize data nfft0/0/ - save nfft0,h + data alast/0.0,0.0,0.0,0.0,0.0/ + data ac/1.0,0.05532,0.11438,0.12918,0.09274/ ! amp coeffs for TS2000 + data a0/0.0,0.0,-0.952,0.768,-0.565/ ! baseline phase coeffs for TS2000 + save nfft0,h,alast,a0,ac + +! disable baseline phase correction for commit - should look for a file with coeffs + a0(1:5)=0.0 df=12000.0/nfft nh=nfft/2 - if(nfft.ne.nfft0) then + if( nfft.ne.nfft0 .or. any(alast .ne. a) ) then t=1.0/2000.0 beta=0.1 pi=4.0*atan(1.0) do i=1,nh+1 ff=(i-1)*df f=ff-1500.0 - h(i)=0. - if(abs(f).le.(1-beta)/(2*t)) h(i)=1.0 - if(abs(f).gt.(1-beta)/(2*t) .and. abs(f).le.(1+beta)/(2*t)) then - h(i)=0.5*(1+cos((pi*t/beta )*(abs(f)-(1-beta)/(2*t)))) + fp=f/1000.0 + h(i)=cmplx(1.0,0.0) + if( equalize ) then + phase0=a0(1)+fp*(a0(2)+fp*(a0(3)+fp*(a0(4)+fp*a0(5)))) + phase=a(1)+fp*(a(2)+fp*(a(3)+fp*(a(4)+fp*a(5)))) +! amp=ac(1)+fp*(ac(2)+fp*(ac(3)+fp*(ac(4)+fp*ac(5)))) + amp=1.0 ! no amplitude correction for now + h(i)=amp*cmplx(cos(phase),sin(phase))*cmplx(cos(phase0),sin(phase0)) + endif + if(abs(f).gt.(1-beta)/(2*t) .and. abs(f).le.(1+beta)/(2*t)) then + h(i)=h(i)*0.5*(1+cos((pi*t/beta )*(abs(f)-(1-beta)/(2*t)))) endif -! h(i)=sqrt(h(i)) enddo nfft0=nfft endif @@ -33,22 +48,10 @@ subroutine analytic(d,npts,nfft,c) c(npts+1:nfft)=0. call four2a(c,nfft,1,-1,1) !Forward c2c FFT -! do i=1,nh -! f=(i-1)*df -! s(i)=real(c(i))**2 + aimag(c(i))**2 -! write(12,3001) f,s(i),db(s(i)) -!3001 format(3f12.3) -! enddo - -! ia=700.0/df -! c(1:ia)=0. -! ib=2300.0/df -! c(ib:nfft)=0. - c(1:nh+1)=h(1:nh+1)*c(1:nh+1) c(1)=0.5*c(1) !Half of DC term c(nh+2:nfft)=0. !Zero the negative frequencies call four2a(c,nfft,1,1,1) !Inverse c2c FFT - + alast=a return end subroutine analytic diff --git a/lib/msk144decodeframe.f90 b/lib/msk144decodeframe.f90 index 68abb1e8c..ff8fdac8d 100644 --- a/lib/msk144decodeframe.f90 +++ b/lib/msk144decodeframe.f90 @@ -1,4 +1,4 @@ -subroutine msk144decodeframe(c,msgreceived,nsuccess) +subroutine msk144decodeframe(c,softbits,msgreceived,nsuccess) ! use timer_module, only: timer parameter (NSPM=864) diff --git a/lib/msk144spd.f90 b/lib/msk144spd.f90 index af4021fc7..4506cef71 100644 --- a/lib/msk144spd.f90 +++ b/lib/msk144spd.f90 @@ -1,4 +1,4 @@ -subroutine msk144spd(cbig,n,ntol,nsuccess,msgreceived,fc,fret,tret,navg) +subroutine msk144spd(cbig,n,ntol,nsuccess,msgreceived,fc,fret,tret,navg,ct,softbits) ! MSK144 short-ping-decoder use timer_module, only: timer @@ -23,6 +23,7 @@ subroutine msk144spd(cbig,n,ntol,nsuccess,msgreceived,fc,fret,tret,navg) real rcw(12) real ferrs(MAXCAND) real snrs(MAXCAND) + real softbits(144) real tonespec(NFFT) real tpat(NPATTERNS) real*8 dt, df, fs, pi, twopi @@ -177,7 +178,7 @@ subroutine msk144spd(cbig,n,ntol,nsuccess,msgreceived,fc,fret,tret,navg) if( is.eq.2) ic0=max(1,ic0-1) if( is.eq.3) ic0=min(NSPM,ic0+1) ct=cshift(c,ic0-1) - call msk144decodeframe(ct,msgreceived,ndecodesuccess) + call msk144decodeframe(ct,softbits,msgreceived,ndecodesuccess) if( ndecodesuccess .gt. 0 ) then tret=(nstart(icand)+NSPM/2)/fs diff --git a/lib/mskrtd.f90 b/lib/mskrtd.f90 index a20ae2dbe..05cd40d69 100644 --- a/lib/mskrtd.f90 +++ b/lib/mskrtd.f90 @@ -28,10 +28,14 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & real d(NFFT1) real pow(8) + real softbits(144) real xmc(NPATTERNS) + real*8 pcoeffs(5) logical*1 bshmsg,bcontest logical first + logical*1 equalized + data first/.true./ data iavpatterns/ & 1,1,1,1,0,0,0,0, & @@ -39,18 +43,20 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & 1,1,1,1,1,0,0,0, & 1,1,1,1,1,1,1,0/ data xmc/2.0,4.5,2.5,3.5/ !Used to set time at center of averaging mask - save first,tsec0,nutc00,pnoise,nsnrlast,msglast,cdat + save first,tsec0,nutc00,pnoise,nsnrlast,msglast,cdat,pcoeffs,equalized if(first) then tsec0=tsec nutc00=nutc0 pnoise=-1.0 + pcoeffs(1:5)=0.0 + equalized=.false. first=.false. endif fc=nrxfreq -!!! Dupe checking should probaby be moved to mainwindow.cpp +! Dupe checking setup if(nutc00.ne.nutc0 .or. tsec.lt.tsec0) then ! reset dupe checker msglast=' ' nsnrlast=-99 @@ -68,7 +74,7 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & fac=1.0/rms d(1:NZ)=fac*d(1:NZ) d(NZ+1:NFFT1)=0. - call analytic(d,NZ,NFFT1,cdat) !Convert to analytic signal and filter + call analytic(d,NZ,NFFT1,cdat,pcoeffs,.true.) !Convert to analytic signal and filter ! Calculate average power for each frame and for the entire block. ! If decode is successful, largest power will be taken as signal+noise. @@ -86,10 +92,8 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & ! center a 3-frame analysis window and attempts to decode each of the ! 3 frames along with 2- and 3-frame averages. np=8*NSPM - call msk144spd(cdat,np,ntol,nsuccess,msgreceived,fc,fest,tdec,navg) + call msk144spd(cdat,np,ntol,nsuccess,msgreceived,fc,fest,tdec,navg,ct,softbits) -!############################################################ -!##### hardwired for testing - need to bring in Sh box status if(nsuccess.eq.0 .and. bshmsg) then call msk40spd(cdat,np,ntol,mycall(1:6),hiscall(1:6),nsuccess, & msgreceived,fc,fest,tdec,navg) @@ -98,6 +102,7 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & if( nsuccess .eq. 1 ) then tdec=tsec+tdec decsym=' & ' + if( equalized ) decsym=' ^ ' ipk=0 is=0 goto 900 @@ -120,15 +125,17 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & if( nsyncsuccess .eq. 0 ) cycle do ipk=1,npeaks - do is=1,3 + do is=1,3 ! With equalization, this loop may not be necessary ic0=npkloc(ipk) if(is.eq.2) ic0=max(1,ic0-1) if(is.eq.3) ic0=min(NSPM,ic0+1) ct=cshift(c,ic0-1) - call msk144decodeframe(ct,msgreceived,ndecodesuccess) + call msk144decodeframe(ct,softbits,msgreceived,ndecodesuccess) if(ndecodesuccess .gt. 0) then tdec=tsec+xmc(iavg)*tframe decsym=' & ' + if( equalized ) decsym=' ^ ' + if( equalized .and. is .ne. 1 ) decsym=' ! ' !help decide if is dither is needed goto 900 endif enddo !Slicer dither @@ -149,7 +156,7 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & go to 999 900 continue -! Successful decode - estimate snr !!! noise estimate needs work +! Successful decode - estimate snr if( pnoise .gt. 0.0 ) then snr0=10.0*log10(pmax/pnoise-1.0) else @@ -157,7 +164,10 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & endif nsnr=nint(snr0) -!!!! Temporary - dupe check. Only print if new message, or higher snr. + call msk144signalquality(ct,snr0,fest,tdec,softbits,msgreceived,hiscall, & + ncorrected,eyeopening,equalized,pcoeffs) + +! Dupe check. Only print if new message, or higher snr. if(msgreceived.ne.msglast .or. nsnr.gt.nsnrlast .or. tsec.lt.tsec0) then msglast=msgreceived nsnrlast=nsnr @@ -168,8 +178,8 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & call fix_contest_msg(mycall(1:6),mygrid,hiscall(1:6),msgreceived) endif write(line,1020) nutc0,nsnr,tdec,nint(fest),decsym,msgreceived, & - navg,char(0) -1020 format(i6.6,i4,f5.1,i5,a3,a22,i2,a1) + navg,ncorrected,eyeopening,char(0) +1020 format(i6.6,i4,f5.1,i5,a3,a22,i2,i3,f5.1,a1) endif 999 tsec0=tsec @@ -177,3 +187,5 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & end subroutine mskrtd include 'fix_contest_msg.f90' + +include 'msk144signalquality.f90'