From 60f9a038985b480810bf17b0f22d6d3c08c072a4 Mon Sep 17 00:00:00 2001 From: Jordan Sherer Date: Thu, 30 Apr 2020 17:17:56 -0400 Subject: [PATCH] Working through decoder fixes --- CMakeLists.txt | 1 - lib/ft8/watterson.f90 | 62 ----------------------------- lib/js8/js8_params.f90 | 2 +- lib/js8/js8dec.f90 | 90 +++++++++++++++++++++--------------------- lib/js8/syncjs8.f90 | 18 +++++++-- lib/js8/syncjs8d.f90 | 77 ++++++++++++++++++++---------------- 6 files changed, 103 insertions(+), 147 deletions(-) delete mode 100644 lib/ft8/watterson.f90 diff --git a/CMakeLists.txt b/CMakeLists.txt index acdab63..87fea0c 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -468,7 +468,6 @@ set (wsjt_FSRCS lib/unpackmsg144.f90 lib/update_recent_calls.f90 lib/update_hasharray.f90 - lib/ft8/watterson.f90 lib/wav11.f90 lib/wav12.f90 lib/xcor.f90 diff --git a/lib/ft8/watterson.f90 b/lib/ft8/watterson.f90 deleted file mode 100644 index d047a26..0000000 --- a/lib/ft8/watterson.f90 +++ /dev/null @@ -1,62 +0,0 @@ -subroutine watterson(c,npts,fs,delay,fspread) - - complex c(0:npts-1) - complex c2(0:npts-1) - complex cs1(0:npts-1) - complex cs2(0:npts-1) - - nonzero=0 - df=fs/npts - if(fspread.gt.0.0) then - do i=0,npts-1 - xx=gran() - yy=gran() - cs1(i)=0.707*cmplx(xx,yy) - xx=gran() - yy=gran() - cs2(i)=0.707*cmplx(xx,yy) - enddo - call four2a(cs1,npts,1,-1,1) !To freq domain - call four2a(cs2,npts,1,-1,1) - do i=0,npts-1 - f=i*df - if(i.gt.npts/2) f=(i-npts)*df - x=(f/(0.5*fspread))**2 - a=0. - if(x.le.50.0) then - a=exp(-x) - endif - cs1(i)=a*cs1(i) - cs2(i)=a*cs2(i) - if(abs(f).lt.10.0) then - p1=real(cs1(i))**2 + aimag(cs1(i))**2 - p2=real(cs2(i))**2 + aimag(cs2(i))**2 - if(p1.gt.0.0) nonzero=nonzero+1 -! write(62,3101) f,p1,p2,db(p1+1.e-12)-60,db(p2+1.e-12)-60 -!3101 format(f10.3,2f12.3,2f10.3) - endif - enddo - call four2a(cs1,npts,1,1,1) !Back to time domain - call four2a(cs2,npts,1,1,1) - cs1(0:npts-1)=cs1(0:npts-1)/npts - cs2(0:npts-1)=cs2(0:npts-1)/npts - endif - - nshift=nint(0.001*delay*fs) - c2(0:npts-1)=cshift(c(0:npts-1),nshift) - sq=0. - do i=0,npts-1 - if(nonzero.gt.1) then - c(i)=0.5*(cs1(i)*c(i) + cs2(i)*c2(i)) - else - c(i)=0.5*(c(i) + c2(i)) - endif - sq=sq + real(c(i))**2 + aimag(c(i))**2 -! write(61,3001) i/12000.0,c(i) -!3001 format(3f12.6) - enddo - rms=sqrt(sq/npts) - c=c/rms - - return -end subroutine watterson diff --git a/lib/js8/js8_params.f90 b/lib/js8/js8_params.f90 index 8396358..1ef6896 100644 --- a/lib/js8/js8_params.f90 +++ b/lib/js8/js8_params.f90 @@ -1 +1 @@ -parameter (NWRITELOG=0) !Write log file? +parameter (NWRITELOG=0) !Write log statements to output diff --git a/lib/js8/js8dec.f90 b/lib/js8/js8dec.f90 index 7196ab2..42981b8 100644 --- a/lib/js8/js8dec.f90 +++ b/lib/js8/js8dec.f90 @@ -32,7 +32,7 @@ subroutine js8dec(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, & integer naptypes(0:5,4) ! (nQSOProgress, decoding pass) maximum of 4 passes for now integer*1, target:: i1hiscall(12) complex cd0(0:NP-1) - complex ctwk(4*NSPS/NDOWN) + complex ctwk(7*NSPS/NDOWN) complex csymb(NDOWNSPS) logical first,newdat,lsubtract,lapon,lapcqonly,nagain equivalence (s1,s1sort) @@ -100,8 +100,8 @@ subroutine js8dec(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, & delfbest=0. ibest=0 - if(NWRITELOG.eq.1) then - write(*,*) ' downsampling', fs2, dt2 + if(NWRITELOG.eq.0) then + write(*,*) ' downsampling', f1, fs2, dt2 flush(6) endif @@ -120,7 +120,7 @@ subroutine js8dec(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, & do idt=i0-NQSYMBOL,i0+NQSYMBOL !Search over +/- one quarter symbol call syncjs8d(cd0,idt,ctwk,0,sync) - if(NWRITELOG.eq.1) then + if(NWRITELOG.eq.0) then write(*,*) ' ', 'idt', idt, 'sync', sync flush(6) endif @@ -131,7 +131,7 @@ subroutine js8dec(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, & enddo xdt2=ibest*dt2 !Improved estimate for DT - if(NWRITELOG.eq.1) then + if(NWRITELOG.eq.0) then write(*,*) ' ', 'xdt2', xdt2, 'ibest', ibest flush(6) endif @@ -144,12 +144,12 @@ subroutine js8dec(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, & delf=ifr*0.5 dphi=twopi*delf*dt2 phi=0.0 - do i=1,(4*NSPS/NDOWN) + do i=1,7*NSPS/NDOWN ctwk(i)=cmplx(cos(phi),sin(phi)) phi=mod(phi+dphi,twopi) enddo call syncjs8d(cd0,i0,ctwk,1,sync) - if(NWRITELOG.eq.1) then + if(NWRITELOG.eq.0) then write(*,*) ' ', 'df', delf, 'sync', sync flush(6) endif @@ -164,7 +164,7 @@ subroutine js8dec(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, & xdt=xdt2 f1=f1+delfbest !Improved estimate of DF - if(NWRITELOG.eq.1) then + if(NWRITELOG.eq.0) then write(*,*) ' ', 'twk', xdt, f1, smax flush(6) endif @@ -267,43 +267,43 @@ subroutine js8dec(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, & ! r2=log(b2+b3+b6+b7)-log(b0+b1+b4+b5) ! r4=log(b4+b5+b6+b7)-log(b0+b1+b2+b3) - if(nQSOProgress .eq. 0 .or. nQSOProgress .eq. 5) then -! When bits 88:115 are set as ap bits, bit 115 lives in symbol 39 along -! with no-ap bits 116 and 117. Take care of metrics for bits 116 and 117. - if(j.eq.39) then ! take care of bits that live in symbol 39 - if(apsym(28).lt.0) then - bmetap(i2)=max(ps(2),ps(3))-max(ps(0),ps(1)) - bmetap(i1)=max(ps(1),ps(3))-max(ps(0),ps(2)) - else - bmetap(i2)=max(ps(6),ps(7))-max(ps(4),ps(5)) - bmetap(i1)=max(ps(5),ps(7))-max(ps(4),ps(6)) - endif - endif - endif - -! When bits 116:143 are set as ap bits, bit 115 lives in symbol 39 along -! with ap bits 116 and 117. Take care of metric for bit 115. - if(j.eq.39) then ! take care of bit 115 - iii=2*(apsym(29)+1)/2 + (apsym(30)+1)/2 ! known values of bits 116 & 117 - if(iii.eq.0) bmetap(i4)=ps(4)-ps(0) - if(iii.eq.1) bmetap(i4)=ps(5)-ps(1) - if(iii.eq.2) bmetap(i4)=ps(6)-ps(2) - if(iii.eq.3) bmetap(i4)=ps(7)-ps(3) - endif - -! bit 144 lives in symbol 48 and will be 1 if it is set as an ap bit. -! take care of metrics for bits 142 and 143 - if(j.eq.48) then ! bit 144 is always 1 - bmetap(i4)=max(ps(5),ps(7))-max(ps(1),ps(3)) - bmetap(i2)=max(ps(3),ps(7))-max(ps(1),ps(5)) - endif - -! bit 154 lives in symbol 52 and will be 0 if it is set as an ap bit -! take care of metrics for bits 155 and 156 - if(j.eq.52) then ! bit 154 will be 0 if it is set as an ap bit. - bmetap(i2)=max(ps(2),ps(3))-max(ps(0),ps(1)) - bmetap(i1)=max(ps(1),ps(3))-max(ps(0),ps(2)) - endif +! if(nQSOProgress .eq. 0 .or. nQSOProgress .eq. 5) then +! ! When bits 88:115 are set as ap bits, bit 115 lives in symbol 39 along +! ! with no-ap bits 116 and 117. Take care of metrics for bits 116 and 117. +! if(j.eq.39) then ! take care of bits that live in symbol 39 +! if(apsym(28).lt.0) then +! bmetap(i2)=max(ps(2),ps(3))-max(ps(0),ps(1)) +! bmetap(i1)=max(ps(1),ps(3))-max(ps(0),ps(2)) +! else +! bmetap(i2)=max(ps(6),ps(7))-max(ps(4),ps(5)) +! bmetap(i1)=max(ps(5),ps(7))-max(ps(4),ps(6)) +! endif +! endif +! endif +! +! ! When bits 116:143 are set as ap bits, bit 115 lives in symbol 39 along +! ! with ap bits 116 and 117. Take care of metric for bit 115. +! if(j.eq.39) then ! take care of bit 115 +! iii=2*(apsym(29)+1)/2 + (apsym(30)+1)/2 ! known values of bits 116 & 117 +! if(iii.eq.0) bmetap(i4)=ps(4)-ps(0) +! if(iii.eq.1) bmetap(i4)=ps(5)-ps(1) +! if(iii.eq.2) bmetap(i4)=ps(6)-ps(2) +! if(iii.eq.3) bmetap(i4)=ps(7)-ps(3) +! endif +! +! ! bit 144 lives in symbol 48 and will be 1 if it is set as an ap bit. +! ! take care of metrics for bits 142 and 143 +! if(j.eq.48) then ! bit 144 is always 1 +! bmetap(i4)=max(ps(5),ps(7))-max(ps(1),ps(3)) +! bmetap(i2)=max(ps(3),ps(7))-max(ps(1),ps(5)) +! endif +! +! ! bit 154 lives in symbol 52 and will be 0 if it is set as an ap bit +! ! take care of metrics for bits 155 and 156 +! if(j.eq.52) then ! bit 154 will be 0 if it is set as an ap bit. +! bmetap(i2)=max(ps(2),ps(3))-max(ps(0),ps(1)) +! bmetap(i1)=max(ps(1),ps(3))-max(ps(0),ps(2)) +! endif enddo diff --git a/lib/js8/syncjs8.f90 b/lib/js8/syncjs8.f90 index 222a716..d1c9c96 100644 --- a/lib/js8/syncjs8.f90 +++ b/lib/js8/syncjs8.f90 @@ -114,7 +114,8 @@ subroutine syncjs8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase) sync_bc=t/t0 !sync2d(i,j)=max(max(max(sync_abc, sync_ab), sync_ac), sync_bc) - sync2d(i,j)=max(sync_abc, sync_ab, sync_bc) + !sync2d(i,j)=max(sync_abc, sync_ab, sync_bc) + sync2d(i,j)=max(sync_abc, sync_bc) enddo enddo @@ -128,8 +129,14 @@ subroutine syncjs8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase) iz=ib-ia+1 call indexx(red(ia:ib),iz,indx) + + npctile=nint(0.40*iz) + if(npctile.lt.1) then ! bail + ncand=0 + return; + endif - ibase=indx(nint(0.40*iz)) - 1 + ia + ibase=indx(npctile) - 1 + ia if(ibase.lt.1) ibase=1 if(ibase.gt.nh1) ibase=nh1 base=red(ibase) @@ -157,6 +164,10 @@ subroutine syncjs8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase) do j=1,i-1 fdiff=abs(candidate0(1,i))-abs(candidate0(1,j)) if(abs(fdiff).lt.AZ) then ! note: this dedupe difference is dependent on symbol spacing + if(NWRITELOG.eq.1) then + write(*,*) ' candidate dupe', fdiff, candidate0(1,i), candidate0(1,j) + flush(6) + endif 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 @@ -177,9 +188,8 @@ subroutine syncjs8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase) do i=1,ncand j=indx(i) if( candidate0(3,j) .ge. syncmin ) then + candidate(2:3,k)=candidate0(2:3,j) candidate(1,k)=abs(candidate0(1,j)) - candidate(2,k)=candidate0(2,j) - candidate(3,k)=candidate0(3,j) k=k+1 endif enddo diff --git a/lib/js8/syncjs8d.f90 b/lib/js8/syncjs8d.f90 index fe3c6e6..d255505 100644 --- a/lib/js8/syncjs8d.f90 +++ b/lib/js8/syncjs8d.f90 @@ -3,11 +3,11 @@ subroutine syncjs8d(cd0,i0,ctwk,itwk,sync) !include 'js8_params.f90' - parameter(NP=NMAX/NDOWN,NSS=NSPS/NDOWN) + parameter(NP=NMAX/NDOWN, NP2=NN*NDOWNSPS) complex cd0(0:NP-1) - complex csynca(4*NSS),csyncb(4*NSS),csyncc(4*NSS) - complex csync2(4*NSS) - complex ctwk(4*NSS) + complex csynca(7*NDOWNSPS),csyncb(7*NDOWNSPS),csyncc(7*NDOWNSPS) + complex csync2(7*NDOWNSPS) + complex ctwk(7*NDOWNSPS) complex z1,z2,z3 logical first data first/.true./ @@ -36,24 +36,22 @@ subroutine syncjs8d(cd0,i0,ctwk,itwk,sync) if( first ) then twopi=8.0*atan(1.0) - k=1 + fs2=12000.0/NDOWN !Sample rate after downsampling + dt2=1/fs2 !Corresponding sample interval + taus=NDOWNSPS*dt2 !Symbol duration + baud=1.0/taus !Keying rate + phia=0.0 phib=0.0 phic=0.0 - !fs2=12000.0/NDOWN !Sample rate after downsampling - !dt2=1/fs2 !Corresponding sample interval - !taus=NDOWNSPS*dt2 !Symbol duration - !baud=1.0/taus !Keying rate - do i=0,6 + dphia=twopi*icos7a(i)*baud*dt2 + dphib=twopi*icos7b(i)*baud*dt2 + dphic=twopi*icos7c(i)*baud*dt2 - dphia=2*twopi*icos7a(i)/real(NSS) - dphib=2*twopi*icos7b(i)/real(NSS) - dphic=2*twopi*icos7c(i)/real(NSS) - - do j=1,NSS/2 - + do j=1,NDOWNSPS + k=i*NDOWNSPS+j csynca(k)=cmplx(cos(phia),sin(phia)) !Waveform for Beginning 7x7 Costas array csyncb(k)=cmplx(cos(phib),sin(phib)) !Waveform for Middle 7x7 Costas array csyncc(k)=cmplx(cos(phic),sin(phic)) !Waveform for End 7x7 Costas array @@ -61,35 +59,46 @@ subroutine syncjs8d(cd0,i0,ctwk,itwk,sync) phib=mod(phib+dphib,twopi) phic=mod(phia+dphic,twopi) - k=k+1 - + if(NWRITELOG.eq.1) then + write(*,*) ' computing costas waveforms', k, i, j, phia, phib, phic, dphia, dphib, dphic + flush(6) + endif enddo enddo first=.false. endif - i1=i0 !four Costas arrays - i2=i0+36*NSS - i3=i0+72*NSS + sync=0 - z1=0. - z2=0. - z3=0. + do i=0,6 + i1=i0+i*NDOWNSPS + i2=i1+36*NDOWNSPS + i3=i1+72*NDOWNSPS - csync2=csynca - if(itwk.eq.1) csync2=ctwk*csynca !Tweak the frequency - if(i1.ge.0 .and. i1+8*NSS-1.le.NP-1) z1=sum(cd0(i1:i1+8*NSS-1:2)*conjg(csync2)) + z1=0. + z2=0. + z3=0. - csync2=csyncb - if(itwk.eq.1) csync2=ctwk*csyncb !Tweak the frequency - if(i2.ge.0 .and. i2+8*NSS-1.le.NP-1) z2=sum(cd0(i2:i2+8*NSS-1:2)*conjg(csync2)) + csync2=csynca + if(itwk.eq.1) csync2=ctwk*csynca + if(i1.ge.0 .and. i1+NDOWNSPS-1.le.NP2-1) z1=sum(cd0(i1:i1+7*NDOWNSPS-1)*conjg(csync2)) - csync2=csyncc - if(itwk.eq.1) csync2=ctwk*csyncc !Tweak the frequency - if(i3.ge.0 .and. i3+8*NSS-1.le.NP-1) z3=sum(cd0(i3:i3+8*NSS-1:2)*conjg(csync2)) + csync2=csyncb + if(itwk.eq.1) csync2=ctwk*csyncb + if(i2.ge.0 .and. i2+NDOWNSPS-1.le.NP2-1) z2=sum(cd0(i2:i2+7*NDOWNSPS-1)*conjg(csync2)) - sync = p(z1) + p(z2) + p(z3) + csync2=csyncc + if(itwk.eq.1) csync2=ctwk*csyncc + if(i3.ge.0 .and. i3+NDOWNSPS-1.le.NP2-1) z3=sum(cd0(i3:i3+7*NDOWNSPS-1)*conjg(csync2)) + + sync = sync + p(z1) + p(z2) + p(z3) + + if(NWRITELOG.eq.1) then + write(*,*) ' sync computation', i, sync + flush(6) + endif + enddo return end subroutine syncjs8d