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

@@ -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

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,   &
+       call qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,minsync,mycall,    &
-            hisgrid,sync,nsnr,dtx,nfreq,decoded,nft)
+            hiscall,hisgrid,sync,nsnr,dtx,nfreq,decoded,nft)
        if (associated(this%callback)) then
           ndrift=0
           nflip=1

lib/qra64a.f90

@@ -1,7 +1,9 @@
-subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,mycall_12,hiscall_12, &
+subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,minsync,mycall_12,   &
-     hisgrid_6,sync,nsnr,dtx,nfreq,decoded,nft)
+     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
+           call timer('qra64_de',1)
-!3333       format(i2,i3,2f8.1,i3)
            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

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,22 +59,38 @@ 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
+  jb=6*5000
-  jstep=200
+  jstep=100
   ka=-maxf1
   kb=maxf1
   ipk=0
   kpk=0
-  do iter=1,2
+!  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)
@@ -91,37 +109,46 @@ subroutine sync64(dd,npts,nf1,nf2,nfqso,ntol,mode64,maxf1,dtx,f0,jpk,kpk,   &
         s2(i)=real(c2(i))**2 + aimag(c2(i))**2
         s3(i)=real(c3(i))**2 + aimag(c3(i))**2
      enddo
+     call timer('sync64_1',1)
+
+     call timer('sync64_2',0)
      do k=ka,kb
-           s0b(ia:ib)=s1(ia-k:ib-k) + s2(ia:ib) + s3(ia+k:ib+k)
+        s0(ia:ib)=s1(ia-k:ib-k) + s2(ia:ib) + s3(ia+k:ib+k)
-           s0b(:ia-1)=0.
+        s0(:ia-1)=0.
-           s0b(ib+1:)=0.
+        s0(ib+1:)=0.
-           nadd=(7*mode64)/2
+        if(nadd.ge.3) then
-           if(mod(nadd,2).eq.0) nadd=nadd+1       !Make nadd odd
+           do ii=1,3
-           if(nadd.ge.3) call smo(s0b(ia:ib),iz,s0(ia:ib),nadd)
+              s0b(ia:ib)=s0(ia:ib)
+              call smo(s0b(ia:ib),iz,s0(ia:ib),nadd)
+        
+           enddo
+        endif
         call smo121(s0(ia:ib),iz)
-           nskip=max(14,2*mode64)
+        call averms(s0(ia+id:ib-id),iz-2*id,nskip,ave,rms)
-           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)
+        ip=ipk0(1) + ia - 1
-              f0=(ipk+ia-1)*df3
+        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
-  enddo
 
+  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