diff --git a/CMakeLists.txt b/CMakeLists.txt index 787ac38e9..a8093b100 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -387,6 +387,7 @@ set (wsjt_FSRCS lib/azdist.f90 lib/badmsg.f90 lib/ft8/baseline.f90 + lib/ft4/ft4_baseline.f90 lib/bpdecode40.f90 lib/bpdecode128_90.f90 lib/ft8/bpdecode174_91.f90 @@ -565,7 +566,6 @@ set (wsjt_FSRCS lib/sync64.f90 lib/sync65.f90 lib/ft4/getcandidates4.f90 - lib/ft4/syncft4.f90 lib/ft8/sync8.f90 lib/ft8/sync8d.f90 lib/ft4/sync4d.f90 diff --git a/lib/decoder.f90 b/lib/decoder.f90 index 174b75115..d3f1153f4 100644 --- a/lib/decoder.f90 +++ b/lib/decoder.f90 @@ -152,8 +152,8 @@ subroutine multimode_decoder(ss,id2,params,nfsample) if(params%nmode.eq.5) then call timer('decft4 ',0) call my_ft4%decode(ft4_decoded,id2,params%nQSOProgress,params%nfqso, & - params%nutc,params%nfa,params%nfb,params%ndepth,ncontest, & - mycall,hiscall) + params%nutc,params%nfa,params%nfb,params%ndepth, & + logical(params%lapcqonly),ncontest,mycall,hiscall) call timer('decft4 ',1) go to 800 endif diff --git a/lib/ft4/ft4_baseline.f90 b/lib/ft4/ft4_baseline.f90 new file mode 100644 index 000000000..85a60d9ec --- /dev/null +++ b/lib/ft4/ft4_baseline.f90 @@ -0,0 +1,49 @@ +subroutine ft4_baseline(s,nfa,nfb,sbase) + +! Fit baseline to spectrum +! Input: s(npts) Linear scale in power +! Output: sbase(npts) Baseline + + include 'ft4_params.f90' + implicit real*8 (a-h,o-z) + real*4 s(NH1) + real*4 sbase(NH1) + real*4 base + real*8 x(1000),y(1000),a(5) + data nseg/10/,npct/10/ + + df=12000.0/NFFT1 !5.21 Hz + ia=max(nint(200.0/df),nfa) + ib=min(NH1,nfb) + do i=ia,ib + s(i)=10.0*log10(s(i)) !Convert to dB scale + enddo + + nterms=5 + nlen=(ib-ia+1)/nseg !Length of test segment + i0=(ib-ia+1)/2 !Midpoint + k=0 + do n=1,nseg !Loop over all segments + ja=ia + (n-1)*nlen + jb=ja+nlen-1 + call pctile(s(ja),nlen,npct,base) !Find lowest npct of points + do i=ja,jb + if(s(i).le.base) then + if (k.lt.1000) k=k+1 !Save all "lower envelope" points + x(k)=i-i0 + y(k)=s(i) + endif + enddo + enddo + kz=k + a=0. + call polyfit(x,y,y,kz,nterms,0,a,chisqr) !Fit a low-order polynomial + do i=ia,ib + t=i-i0 + sbase(i)=a(1)+t*(a(2)+t*(a(3)+t*(a(4)+t*(a(5))))) + 0.65 +! write(51,3051) i*df,s(i),sbase(i) +!3051 format(3f12.3) + sbase(i)=10**(sbase(i)/10.0) + enddo + return +end subroutine ft4_baseline diff --git a/lib/ft4/ft4sim.f90 b/lib/ft4/ft4sim.f90 index e35e5ad8a..959cf423c 100644 --- a/lib/ft4/ft4sim.f90 +++ b/lib/ft4/ft4sim.f90 @@ -110,8 +110,10 @@ program ft4sim c0((NN+1)*NSPS:(NN+2)*NSPS-1)=c0((NN+1)*NSPS:(NN+2)*NSPS-1)*(1.0+cos(twopi*(/(i,i=0,NSPS-1)/)/(2.0*NSPS) ))/2.0 c0((NN+2)*NSPS:)=0. - k=nint((xdt+0.5)/dt) + k=nint((xdt+0.5)/dt)-NSPS c0=cshift(c0,-k) + if(k.gt.0) c0(0:k-1)=0.0 + if(k.lt.0) c0(NZZ+k:NZZ-1)=0.0 do ifile=1,nfiles c=c0 diff --git a/lib/ft4/ft4sim_mult.f90 b/lib/ft4/ft4sim_mult.f90 index 582a48f77..d93978e5f 100644 --- a/lib/ft4/ft4sim_mult.f90 +++ b/lib/ft4/ft4sim_mult.f90 @@ -17,8 +17,6 @@ program ft4sim_mult integer itone(NN) integer*1 msgbits(77) integer*2 iwave(NZZ) !Generated full-length waveform - integer icos4(4) - data icos4/0,1,3,2/ ! Get command-line argument(s) nargs=iargc() diff --git a/lib/ft4/getcandidates4.f90 b/lib/ft4/getcandidates4.f90 index 97aeb54bf..dd875336d 100644 --- a/lib/ft4/getcandidates4.f90 +++ b/lib/ft4/getcandidates4.f90 @@ -10,7 +10,6 @@ subroutine getcandidates4(dd,fa,fb,syncmin,nfqso,maxcand,savg,candidate, & complex cx(0:NH1) real candidate(3,maxcand) real dd(NMAX) - integer indx(NH1) integer ipk(1) equivalence (x,cx) logical first @@ -26,7 +25,6 @@ subroutine getcandidates4(dd,fa,fb,syncmin,nfqso,maxcand,savg,candidate, & ! Compute symbol spectra, stepping by NSTEP steps. savg=0. - tstep=NSTEP/12000.0 df=12000.0/NFFT1 fac=1.0/300.0 do j=1,NHSYM @@ -40,27 +38,20 @@ subroutine getcandidates4(dd,fa,fb,syncmin,nfqso,maxcand,savg,candidate, & enddo savg=savg + s(1:NH1,j) !Average spectrum enddo + savg=savg/NHSYM savsm=0. do i=8,NH1-7 savsm(i)=sum(savg(i-7:i+7))/15. enddo + nfa=fa/df - if(nfa.lt.1) nfa=1 + if(nfa.lt.8) nfa=8 nfb=fb/df if(nfb.gt.nint(5000.0/df)) nfb=nint(5000.0/df) - n300=300/df - n2500=2500/df -! np=nfb-nfa+1 - np=n2500-n300+1 - indx=0 - call indexx(savsm(n300:n2500),np,indx) - xn=savsm(n300+indx(nint(0.3*np))) ncand=0 - if(xn.le.1.e-8) return - savsm=savsm/xn -! call ft4_baseline(savg,nfa,nfb,sbase) -! savsm=savsm/sbase - + call ft4_baseline(savg,nfa,nfb,sbase) + if(any(sbase(nfa:nfb).le.0)) return + savsm(nfa:nfb)=savsm(nfa:nfb)/sbase(nfa:nfb) f_offset = -1.5*12000.0/NSPS do i=nfa+1,nfb-1 if(savsm(i).ge.savsm(i-1) .and. savsm(i).ge.savsm(i+1) .and. & diff --git a/lib/ft4/syncft4.f90 b/lib/ft4/syncft4.f90 deleted file mode 100644 index 4841c560e..000000000 --- a/lib/ft4/syncft4.f90 +++ /dev/null @@ -1,145 +0,0 @@ -subroutine syncft4(iwave,nfa,nfb,syncmin,nfqso,maxcand,s,candidate, & - ncand,sbase) - - include 'ft4_params.f90' -! Search over +/- 2.5s relative to 0.5s TX start time. - parameter (JZ=20) - complex cx(0:NH1) - real s(NH1,NHSYM) - real savg(NH1) - real sbase(NH1) - real x(NFFT1) - real sync2d(NH1,-JZ:JZ) - real red(NH1) - real candidate0(3,maxcand) - real candidate(3,maxcand) - real dd(NMAX) - integer jpeak(NH1) - integer indx(NH1) - integer ii(1) - integer*2 iwave(NMAX) - integer icos4(0:3) - data icos4/0,1,3,2/ !Costas 4x4 tone pattern - equivalence (x,cx) - - dd=iwave/1e3 -! Compute symbol spectra, stepping by NSTEP steps. - savg=0. - tstep=NSTEP/12000.0 - df=12000.0/NFFT1 - fac=1.0/300.0 - do j=1,NHSYM - ia=(j-1)*NSTEP + 1 - ib=ia+NSPS-1 - x(1:NSPS)=fac*dd(ia:ib) - x(NSPS+1:)=0. - call four2a(x,NFFT1,1,-1,0) !r2c FFT - do i=1,NH1 - s(i,j)=real(cx(i))**2 + aimag(cx(i))**2 - enddo - savg=savg + s(1:NH1,j) !Average spectrum - enddo - - call baseline(savg,nfa,nfb,sbase) - - ia=max(1,nint(nfa/df)) - ib=nint(nfb/df) - nssy=NSPS/NSTEP ! # steps per symbol - nfos=NFFT1/NSPS ! # frequency bin oversampling factor - jstrt=0.25/tstep - candidate0=0. - k=0 - - do i=ia,ib - do j=-JZ,+JZ - ta=0. - tb=0. - tc=0. - t0a=0. - t0b=0. - t0c=0. - do n=0,3 - m=j+jstrt+nssy*n - if(m.ge.1.and.m.le.NHSYM) then - ta=ta + s(i+nfos*icos4(n),m) - t0a=t0a + sum(s(i:i+nfos*3:nfos,m)) - endif - tb=tb + s(i+nfos*icos4(n),m+nssy*36) - t0b=t0b + sum(s(i:i+nfos*3:nfos,m+nssy*36)) - if(m+nssy*72.le.NHSYM) then - tc=tc + s(i+nfos*icos4(n),m+nssy*72) - t0c=t0c + sum(s(i:i+nfos*3:nfos,m+nssy*72)) - endif - enddo - t=ta+tb+tc - t0=t0a+t0b+t0c - t0=(t0-t)/3.0 - sync_abc=t/t0 - t=tb+tc - t0=t0b+t0c - t0=(t0-t)/3.0 - sync_bc=t/t0 - sync2d(i,j)=max(sync_abc,sync_bc) - enddo - enddo - - red=0. - do i=ia,ib - ii=maxloc(sync2d(i,-JZ:JZ)) - 1 - JZ - j0=ii(1) - jpeak(i)=j0 - red(i)=sync2d(i,j0) - enddo - iz=ib-ia+1 - call indexx(red(ia:ib),iz,indx) - ibase=indx(nint(0.40*iz)) - 1 + ia - if(ibase.lt.1) ibase=1 - if(ibase.gt.nh1) ibase=nh1 - base=red(ibase) - red=red/base - do i=1,min(maxcand,iz) - n=ia + indx(iz+1-i) - 1 - if(red(n).lt.syncmin.or.isnan(red(n)).or.k.eq.maxcand) exit - k=k+1 -! candidate0(1,k)=n*df+37.5*1.5 - candidate0(1,k)=n*df - candidate0(2,k)=(jpeak(n)-1)*tstep - candidate0(3,k)=red(n) - enddo - ncand=k - -! Put nfqso at top of list, and save only the best of near-dupe freqs. - do i=1,ncand - if(abs(candidate0(1,i)-nfqso).lt.10.0) candidate0(1,i)=-candidate0(1,i) - if(i.ge.2) then - do j=1,i-1 - fdiff=abs(candidate0(1,i))-abs(candidate0(1,j)) - if(abs(fdiff).lt.4.0) then - if(candidate0(3,i).ge.candidate0(3,j)) candidate0(3,j)=0. - if(candidate0(3,i).lt.candidate0(3,j)) candidate0(3,i)=0. - endif - enddo - endif - enddo - - fac=20.0/maxval(s) - s=fac*s - -! Sort by sync -! call indexx(candidate0(3,1:ncand),ncand,indx) -! Sort by frequency - call indexx(candidate0(1,1:ncand),ncand,indx) - k=1 -! do i=ncand,1,-1 - do i=1,ncand - j=indx(i) -! if( candidate0(3,j) .ge. syncmin .and. candidate0(2,j).ge.-1.5 ) then - if( candidate0(3,j) .ge. syncmin ) then - candidate(2:3,k)=candidate0(2:3,j) - candidate(1,k)=abs(candidate0(1,j)) - k=k+1 - endif - enddo - ncand=k-1 - return -end subroutine syncft4 diff --git a/lib/ft4_decode.f90 b/lib/ft4_decode.f90 index cf4acaab7..50c9fa47d 100644 --- a/lib/ft4_decode.f90 +++ b/lib/ft4_decode.f90 @@ -24,7 +24,7 @@ module ft4_decode contains subroutine decode(this,callback,iwave,nQSOProgress,nfqso, & - nutc,nfa,nfb,ndepth,ncontest,mycall,hiscall) + nutc,nfa,nfb,ndepth,lapcqonly,ncontest,mycall,hiscall) use timer_module, only: timer use packjt77 include 'ft4/ft4_params.f90' @@ -74,7 +74,8 @@ contains logical nohiscall,unpk77_success logical one(0:255,0:7) ! 256 4-symbol sequences, 8 bits logical first, dobigfft - logical dosubtract + logical dosubtract,doosd + logical, intent(in) :: lapcqonly data icos4a/0,1,3,2/ data icos4b/1,0,2,3/ @@ -210,19 +211,23 @@ contains fb=nfb dd=iwave -! ndepth=3: 2 passes, subtract on each pass -! ndepth=2: 1 pass, no subtraction -! ndepth=1: 1 pass, no subtraction, fewer candidates +! ndepth=3: 3 passes, bp+osd +! ndepth=2: 3 passes, bp only +! ndepth=1: 1 pass, no subtraction max_iterations=40 syncmin=1.2 dosubtract=.true. + doosd=.true. nsp=3 - if(ndepth.lt.3) then + if(ndepth.eq.2) then + doosd=.false. + endif + if(ndepth.eq.1) then nsp=1 dosubtract=.false. + doosd=.false. endif - if(ndepth.eq.1) syncmin=2.0 do isp = 1,nsp if(isp.eq.2) then @@ -255,8 +260,8 @@ contains idfmin=-12 idfmax=12 idfstp=3 - ibmin=-333 - ibmax=1000 + ibmin=-344 + ibmax=1012 ibstp=4 else idfmin=idfbest-4 @@ -403,6 +408,7 @@ contains apmag=maxval(abs(llra))*1.1 npasses=3+nappasses(nQSOProgress) + if(lapcqonly) npasses=4 if(ndepth.eq.1) npasses=3 if(ncontest.ge.5) npasses=3 ! Don't support Fox and Hound do ipass=1,npasses @@ -417,6 +423,7 @@ contains if(ipass .gt. 3) then llrd=llra iaptype=naptypes(nQSOProgress,ipass-3) + if(lapcqonly) iaptype=1 ! ncontest=0 : NONE ! 1 : NA_VHF @@ -482,10 +489,22 @@ contains llr=llrd endif message77=0 + dmin=0.0 call timer('bpdec174',0) call bpdecode174_91(llr,apmask,max_iterations,message77, & cw,nharderror,niterations) call timer('bpdec174',1) + + if(doosd .and. nharderror.lt.0) then + ndeep=3 + if(abs(nfqso-f1).le.napwid) then + ndeep=4 + endif + call timer('osd174_91 ',0) + call osd174_91(llr,apmask,ndeep,message77,cw,nharderror,dmin) + call timer('osd174_91 ',1) + endif + if(sum(message77).eq.0) cycle if( nharderror.ge.0 ) then message77=mod(message77+rvec,2) ! remove rvec scrambling @@ -508,11 +527,11 @@ contains if(snr.gt.0.0) then xsnr=10*log10(snr)-14.8 else - xsnr=-20.0 + xsnr=-21.0 endif - nsnr=nint(max(-20.0,xsnr)) + nsnr=nint(max(-21.0,xsnr)) xdt=ibest/666.67 - 0.5 -!write(21,'(i6.6,i5,2x,f4.1,i6,2x,a37,2x,f4.1,3i3)') nutc,nsnr,xdt,nint(f0),message,sync,iaptype,ipass,isp +!write(21,'(i6.6,i5,2x,f4.1,i6,2x,a37,2x,f4.1,3i3,f5.1)') nutc,nsnr,xdt,nint(f0),message,sync,iaptype,ipass,isp,dmin call this%callback(sync,nsnr,xdt,f0,message,iaptype,qual) exit endif