diff --git a/lib/badmsg.f90 b/lib/badmsg.f90 index 56366f9d2..578270dce 100644 --- a/lib/badmsg.f90 +++ b/lib/badmsg.f90 @@ -18,8 +18,29 @@ subroutine badmsg(irc,dat,nc1,nc2,ng2) iand(ishft(dat(10),-4),3) ig=ishft(iand(dat(10),15),12) + ishft(dat(11),6) + dat(12) + +! Test for blank, -01 to -30, R-01 to R-30, RO, RRR, 73 + if(ig.ge.32401 .and. ig.le.32461) return - if(ic1.eq.nc1 .and. ic2.eq.nc2 .and. ng.ne.32401 .and. ig.ne.ng) irc=-1 + if(ig.ge.14220 .and. ig.le.14229) return !-41 to -50 + if(ig.ge.14040 .and. ig.le.14049) return !-31 to -40 + + if(ig.ge.13320 .and. ig.le.13329) return !+00 to +09 + if(ig.ge.13140 .and. ig.le.13149) return !+10 to +19 + if(ig.ge.12960 .and. ig.le.12969) return !+20 to +29 + if(ig.ge.12780 .and. ig.le.12789) return !+30 to +39 + if(ig.ge.12600 .and. ig.le.12609) return !+40 to +49 + + if(ig.ge.12420 .and. ig.le.12429) return !R-41 to R-50 + if(ig.ge.12240 .and. ig.le.12249) return !R-31 to R-40 + + if(ig.ge.11520 .and. ig.le.11529) return !R+00 to R+09 + if(ig.ge.11340 .and. ig.le.11349) return !R+10 to R+19 + if(ig.ge.11160 .and. ig.le.11169) return !R+20 to R+29 + if(ig.ge.10980 .and. ig.le.10989) return !R+30 to R+39 + if(ig.ge.10800 .and. ig.le.10809) return !R+40 to R+49 + + if(ic1.eq.nc1 .and. ic2.eq.nc2 .and. ng2.ne.32401 .and. ig.ne.ng2) irc=-1 return end subroutine badmsg diff --git a/lib/pctile.f90 b/lib/pctile.f90 index 867070c52..92253cddf 100644 --- a/lib/pctile.f90 +++ b/lib/pctile.f90 @@ -1,6 +1,6 @@ subroutine pctile(x,npts,npct,xpct) - parameter (NMAX=32768) + parameter (NMAX=100000) real*4 x(npts) real*4 tmp(NMAX) diff --git a/lib/qra64a.f90 b/lib/qra64a.f90 index d37f17ba4..c30554784 100644 --- a/lib/qra64a.f90 +++ b/lib/qra64a.f90 @@ -62,7 +62,12 @@ subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,minsync,ndepth, & 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 + +! sync=5 +! dtx=0. +! f0=1000. + +! if(sync.lt.minsync+3.5) go to 900 npts2=npts/2 itz=11 @@ -86,9 +91,21 @@ subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,minsync,ndepth, & a(2)=-0.67*(idf1 + 0.67*kpk) call twkfreq(c00,c0,npts2,6000.0,a) call spec64(c0,npts2,mode64,jpk,s3a,LL,NN) + +!### Parameters 3000.0 and 10.0 might be better optimized? + call pctile(s3a,LL*NN,45,base) + s3max=3000.0 + if(sync.gt.5.0) s3max=15000.0/sync + s3max=max(10.0,s3max) + do i=1,LL*NN + if(s3a(i).gt.s3max*base) s3a(i)=s3max*base + enddo + s3pk=maxval(s3a(1:LL*NN)) + napmin=99 do iter=itz,0,-2 b90=1.728**iter + if(b90.gt.230.0) cycle s3(1:LL*NN)=s3a(1:LL*NN) call timer('qra64_de',0) call qra64_dec(s3,nc1,nc2,ng2,naptype,0,nSubmode,b90, & @@ -114,6 +131,7 @@ subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,minsync,ndepth, & irc=irckeep endif 10 decoded=' ' + if(irc.ge.0) then call unpackmsg(dat4,decoded) !Unpack the user message call fmtmsg(decoded,iz) @@ -131,7 +149,7 @@ subroutine qra64a(dd,npts,nutc,nf1,nf2,nfqso,ntol,mode64,minsync,ndepth, & irc=-1 endif if(irc.lt.0) then - sy=max(1.0,sync+1.0) + sy=max(1.0,sync-2.5) 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 diff --git a/lib/sync64.f90 b/lib/sync64.f90 index 2fee84022..44e079dc4 100644 --- a/lib/sync64.f90 +++ b/lib/sync64.f90 @@ -13,6 +13,7 @@ subroutine sync64(dd,npts,nf1,nf2,nfqso,ntol,mode64,maxf1,dtx,f0,jpk,kpk, & real s0(0:NSPC-1) !Sum of s1+s2+s3 real s0a(0:NSPC-1) !Best synchromized spectrum (saved) real s0b(0:NSPC-1) !tmp + real s0c(0:NSPC-1) !tmp logical old_qra_sync integer icos7(0:6) !Costas 7x7 tones integer ipk0(1) @@ -78,17 +79,9 @@ subroutine sync64(dd,npts,nf1,nf2,nfqso,ntol,mode64,maxf1,dtx,f0,jpk,kpk, & ja=0 jb=6*5000 jstep=100 - ka=-maxf1 - kb=maxf1 ipk=0 - kpk=0 -! 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) - + nskip=max(14,10*mode64) + do j1=ja,jb,jstep call timer('sync64_1',0) j2=j1 + 39*NSPS @@ -112,39 +105,40 @@ subroutine sync64(dd,npts,nf1,nf2,nfqso,ntol,mode64,maxf1,dtx,f0,jpk,kpk, & 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) + s0(ia:ib)=s1(ia:ib) + s2(ia:ib) + s3(ia:ib) + s0(:ia-1)=0. + s0(ib+1:)=0. + smax=0. + do na=0,5 + nadd=7*(2**na) + if(nadd.gt.161) nadd=161 + if(mod(nadd,2).eq.0) nadd=nadd+1 + call smo(s0(ia:ib)/nadd,iz,s0b(ia:ib),nadd) + call smo(s0b(ia:ib)/nadd,iz,s0(ia:ib),nadd) 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 + s=0. + if(rms.gt.0.0) s=(maxval(s0(ia:ib))-ave)/rms + if(s.gt.smax) then + smax=s + nabest=na + avebest=ave + rmsbest=rms + s0c(ia:ib)=s0(ia:ib) endif enddo + s0=s0c + ipk0=maxloc(s0(ia:ib)) + ip=ipk0(1) + ia - 1 + if(smax.gt.sync .and. ip.ge.iaa .and. ip.le.ibb) then + jpk=j1 + s0a=(s0-avebest)/rmsbest + sync=smax + dtx=jpk/6000.0 - 1.0 + ipk=ip + f0=ip*df3 + endif 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