From 46a1693e835b6f96c43100d3f25b7128592ae72a Mon Sep 17 00:00:00 2001 From: Joe Taylor Date: Tue, 15 Nov 2016 21:39:55 +0000 Subject: [PATCH] Working on QRA64 sync, snr estimates, etc., for QRA64 decoder. git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@7321 ab8295b8-cf94-4d9e-aec4-7959e3be5d79 --- lib/decoder.f90 | 2 +- lib/jt65_decode.f90 | 4 +- lib/qra64a.f90 | 34 +++++++++--- lib/sync64.f90 | 129 ++++++++++++++++++++++++++------------------ 4 files changed, 109 insertions(+), 60 deletions(-) diff --git a/lib/decoder.f90 b/lib/decoder.f90 index 8cd5067c6..6049bb49e 100644 --- a/lib/decoder.f90 +++ b/lib/decoder.f90 @@ -249,7 +249,7 @@ contains if(ft.ge.80) then nft=ft-100 csync=': ' - if(snr.gt.-30 .or. nft.ge.0) csync=':*' + if(sync.ge.float(minsync) .or. nft.ge.0) csync=':*' if(nft.lt.0) then write(*,1009) params%nutc,snr,dt,freq,csync,decoded else diff --git a/lib/jt65_decode.f90 b/lib/jt65_decode.f90 index 71d0b42a2..2ac0c6d65 100644 --- a/lib/jt65_decode.f90 +++ b/lib/jt65_decode.f90 @@ -99,8 +99,8 @@ contains if(nsubmode.ge.100) then ! This is QRA64 mode mode64=2**(nsubmode-100) - call qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,mycall,hiscall, & - hisgrid,sync,nsnr,dtx,nfreq,decoded,nft) + call qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,minsync,mycall, & + hiscall,hisgrid,sync,nsnr,dtx,nfreq,decoded,nft) if (associated(this%callback)) then ndrift=0 nflip=1 diff --git a/lib/qra64a.f90 b/lib/qra64a.f90 index a0bb1cdab..fb8210595 100644 --- a/lib/qra64a.f90 +++ b/lib/qra64a.f90 @@ -1,7 +1,9 @@ -subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,mycall_12,hiscall_12, & - hisgrid_6,sync,nsnr,dtx,nfreq,decoded,nft) +subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,minsync,mycall_12, & + hiscall_12,hisgrid_6,sync,nsnr,dtx,nfreq,decoded,nft) use packjt + use timer_module, only: timer + parameter (NMAX=60*12000,LN=1152*63) character decoded*22 character*12 mycall_12,hiscall_12 @@ -20,6 +22,8 @@ subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,mycall_12,hiscall_12, & data nc1z/-1/,nc2z/-1/,ng2z/-1/ save + call timer('qra64a ',0) + decoded=' ' if(nfqso.lt.nf1 .or. nfqso.gt.nf2) go to 900 nft=99 nsnr=-30 @@ -47,7 +51,13 @@ subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,mycall_12,hiscall_12, & endif maxf1=0 - call sync64(dd,npts,nf1,nf2,nfqso,ntol,mode64,maxf1,dtx,f0,jpk,kpk,snr1,c00) +! write(60) npts,dd(1:npts),nf1,nf2,nfqso,ntol,mode64,maxf1 + call timer('sync64 ',0) + call sync64(dd,npts,nf1,nf2,nfqso,ntol,mode64,maxf1,dtx,f0,jpk,kpk,sync,c00) + call timer('sync64 ',1) + irc=-99 + if(sync.lt.float(minsync)) go to 900 + npts2=npts/2 itz=10 if(mode64.eq.4) itz=9 @@ -75,10 +85,10 @@ subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,mycall_12,hiscall_12, & do iter=itz,0,-1 b90=1.728**iter s3(1:LL*NN)=s3a(1:LL*NN) + call timer('qra64_de',0) call qra64_dec(s3,nc1,nc2,ng2,naptype,0,nSubmode,b90, & nFadingModel,dat4,snr2,irc) -! if(irc.ge.0) write(*,3333) iter,idf0,-a(1),b90,irc -!3333 format(i2,i3,2f8.1,i3) + call timer('qra64_de',1) if(abs(snr2).gt.30.) snr2=-30.0 if(irc.eq.0) go to 10 if(irc.gt.0 .and. irc.le.ircmin) then @@ -104,7 +114,6 @@ subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,mycall_12,hiscall_12, & else snr2=0. endif -! if(irc.ge.0) go to 900 enddo enddo 900 continue @@ -113,6 +122,19 @@ subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,mycall_12,hiscall_12, & decoded=' ' irc=-1 endif + if(irc.lt.0) then + sy=max(1.0,sync+1.0) + if(nSubmode.eq.0) nsnr=nint(10.0*log10(sy)-38.0) !A + if(nSubmode.eq.1) nsnr=nint(10.0*log10(sy)-36.0) !B + if(nSubmode.eq.2) nsnr=nint(10.0*log10(sy)-34.0) !C + if(nSubmode.eq.3) nsnr=nint(10.0*log10(sy)-29.0) !D + if(nSubmode.eq.4) nsnr=nint(10.0*log10(sy)-24.0) !E + endif + +! write(70,3303) sync,snr2,nsnr,irc +!3303 format(2f8.1,2i5) + + call timer('qra64a ',1) return end subroutine qra64a diff --git a/lib/sync64.f90 b/lib/sync64.f90 index 079934838..2fee84022 100644 --- a/lib/sync64.f90 +++ b/lib/sync64.f90 @@ -1,5 +1,7 @@ subroutine sync64(dd,npts,nf1,nf2,nfqso,ntol,mode64,maxf1,dtx,f0,jpk,kpk, & - snrdb,c0) + sync,c0) + + use timer_module, only: timer parameter (NMAX=60*12000) !Max size of raw data at 12000 Hz parameter (NSPS=3456) !Samples per symbol at 6000 Hz @@ -57,71 +59,96 @@ subroutine sync64(dd,npts,nf1,nf2,nfqso,ntol,mode64,maxf1,dtx,f0,jpk,kpk, & nfft3=NSPC nh3=nfft3/2 df3=6000.0/nfft3 + fa=max(nf1,nfqso-ntol) fb=min(nf2,nfqso+ntol) + iaa=max(maxf1,nint(fa/df3)) + ibb=min(NSPC-1-maxf1,nint(fb/df3)) + + maxtol=max(ntol,500) + fa=max(nf1,nfqso-maxtol) + fb=min(nf2,nfqso+maxtol) ia=max(maxf1,nint(fa/df3)) ib=min(NSPC-1-maxf1,nint(fb/df3)) + id=0.1*(ib-ia) + iz=ib-ia+1 - snr=0. + sync=0. jpk=0 ja=0 - jb=6*6000 - jstep=200 + jb=6*5000 + jstep=100 ka=-maxf1 kb=maxf1 ipk=0 kpk=0 - do iter=1,2 - do j1=ja,jb,jstep - j2=j1 + 39*NSPS - j3=j1 + 77*NSPS - c1=1.e-4*c0(j1:j1+NSPC-1) * conjg(cc) - call four2a(c1,nfft3,1,-1,1) - c2=1.e-4*c0(j2:j2+NSPC-1) * conjg(cc) - call four2a(c2,nfft3,1,-1,1) - c3=1.e-4*c0(j3:j3+NSPC-1) * conjg(cc) - call four2a(c3,nfft3,1,-1,1) - s0=0. - s1=0. - s2=0. - s3=0. - do i=ia,ib - freq=i*df3 - s1(i)=real(c1(i))**2 + aimag(c1(i))**2 - s2(i)=real(c2(i))**2 + aimag(c2(i))**2 - s3(i)=real(c3(i))**2 + aimag(c3(i))**2 - enddo - do k=ka,kb - s0b(ia:ib)=s1(ia-k:ib-k) + s2(ia:ib) + s3(ia+k:ib+k) - s0b(:ia-1)=0. - s0b(ib+1:)=0. - nadd=(7*mode64)/2 - if(mod(nadd,2).eq.0) nadd=nadd+1 !Make nadd odd - if(nadd.ge.3) call smo(s0b(ia:ib),iz,s0(ia:ib),nadd) - call smo121(s0(ia:ib),iz) - nskip=max(14,2*mode64) - call averms(s0(ia:ib),iz,nskip,ave,rms) - s=(maxval(s0(ia:ib))-ave)/rms - if(s.gt.snr) then - jpk=j1 - s0a=s0/rms - snr=s - dtx=jpk/6000.0 - 1.0 - ipk0=maxloc(s0(ia:ib)) - ipk=ipk0(1) - f0=(ipk+ia-1)*df3 - kpk=k - endif - enddo +! nadd=(7*mode64)/2 +! nadd=7*mode64 + nadd=10*mode64 + if(mod(nadd,2).eq.0) nadd=nadd+1 !Make nadd odd +! nskip=max(14,2*mode64) + nskip=max(14,nadd) + + do j1=ja,jb,jstep + call timer('sync64_1',0) + j2=j1 + 39*NSPS + j3=j1 + 77*NSPS + c1=1.e-4*c0(j1:j1+NSPC-1) * conjg(cc) + call four2a(c1,nfft3,1,-1,1) + c2=1.e-4*c0(j2:j2+NSPC-1) * conjg(cc) + call four2a(c2,nfft3,1,-1,1) + c3=1.e-4*c0(j3:j3+NSPC-1) * conjg(cc) + call four2a(c3,nfft3,1,-1,1) + s0=0. + s1=0. + s2=0. + s3=0. + do i=ia,ib + freq=i*df3 + s1(i)=real(c1(i))**2 + aimag(c1(i))**2 + s2(i)=real(c2(i))**2 + aimag(c2(i))**2 + s3(i)=real(c3(i))**2 + aimag(c3(i))**2 enddo - ja=max(0,jpk-2*jstep) - jb=min(336000-NSPC,jpk+2*jstep) - jstep=10 - enddo + call timer('sync64_1',1) + call timer('sync64_2',0) + do k=ka,kb + s0(ia:ib)=s1(ia-k:ib-k) + s2(ia:ib) + s3(ia+k:ib+k) + s0(:ia-1)=0. + s0(ib+1:)=0. + if(nadd.ge.3) then + do ii=1,3 + s0b(ia:ib)=s0(ia:ib) + call smo(s0b(ia:ib),iz,s0(ia:ib),nadd) + + enddo + endif + call smo121(s0(ia:ib),iz) + call averms(s0(ia+id:ib-id),iz-2*id,nskip,ave,rms) + s=(maxval(s0(ia:ib))-ave)/rms + ipk0=maxloc(s0(ia:ib)) + ip=ipk0(1) + ia - 1 + if(s.gt.sync .and. ip.ge.iaa .and. ip.le.ibb) then + jpk=j1 + s0a=(s0-ave)/rms + sync=s + dtx=jpk/6000.0 - 1.0 + ipk=ip + f0=ip*df3 + kpk=k + endif + enddo + call timer('sync64_2',1) + enddo + sync=sync-3.5 + + ja=max(0,jpk-2*jstep) + jb=min(336000-NSPC,jpk+2*jstep) + jstep=10 + + s0a=s0a+2.0 write(17) ia,ib,s0a(ia:ib) !Save data for red curve close(17) - snrdb=10.0*log10(snr)-39.0 return end subroutine sync64