diff --git a/lib/genqra64.f90 b/lib/genqra64.f90 index babe7d692..fb0c19438 100644 --- a/lib/genqra64.f90 +++ b/lib/genqra64.f90 @@ -1,6 +1,6 @@ subroutine genqra64(msg0,ichk,msgsent,itone,itype) -! Encodes a QRA64 message to yield itone(1:84) +! Encodes a QRA64 message to yield itone(1:84) or a QRA65 msg, itone(1:85) use packjt character*22 msg0 @@ -8,15 +8,15 @@ subroutine genqra64(msg0,ichk,msgsent,itone,itype) character*22 msgsent !Message as it will be received integer itone(85) !QRA64 uses only 84 character*3 cok !' ' or 'OOO' - logical old_qra_sync integer dgen(13) integer sent(63) - integer b11(11) !Barker 11 code + integer isync(22) integer icos7(0:6) data icos7/2,5,6,0,4,1,3/ !Defines a 7x7 Costas array - data b11/1,1,1,0,0,0,1,0,0,1,0/ !Barker 11 definition + data isync/1,9,12,13,15,22,23,26,27,33,35,38,46,50,55,60,62,66,69,74,76,85/ save + print*,'ichk:',ichk if(msg0(1:1).eq.'@') then read(msg0(2:5),*,end=1,err=1) nfreq go to 2 @@ -44,18 +44,17 @@ subroutine genqra64(msg0,ichk,msgsent,itone,itype) if(ichk.eq.1) go to 999 !Return if checking only call qra64_enc(dgen,sent) !Encode using QRA64 - if(ichk.eq.66) then -! Experimental QRA66 (FST66?) mode - j=0 + if(ichk.eq.65) then +! Experimental QRA65 mode + j=1 k=0 do i=1,85 - if(mod(i,4).eq.1) then + if(i.eq.isync(j)) then j=j+1 !Index for next sync symbol - if(j.eq.12) j=1 - itone(i)=b11(j) !Insert a sync symbol + itone(i)=0 !Insert a sync symbol else k=k+1 - itone(i)=sent(k) + 2 + itone(i)=sent(k) + 1 endif enddo else diff --git a/lib/qra/qra66/qra66sim.f90 b/lib/qra/qra66/qra66sim.f90 index 1171542ad..335a37252 100644 --- a/lib/qra/qra66/qra66sim.f90 +++ b/lib/qra/qra66/qra66sim.f90 @@ -65,7 +65,7 @@ program qra66sim mode66=2**(ichar(csubmode) - ichar('A')) print*,csubmode,mode66 - ichk=66 !Flag sent to genqra64 + ichk=65 !Flag sent to genqra64 call genqra64(msg,ichk,msgsent,itone,itype) write(*,1001) itone 1001 format('Channel symbols:'/(20i3)) diff --git a/lib/qra66_decode.f90 b/lib/qra66_decode.f90 index f175fb0ba..f0db7378e 100644 --- a/lib/qra66_decode.f90 +++ b/lib/qra66_decode.f90 @@ -115,6 +115,10 @@ contains endif naptype=maxaptype + call timer('sync66 ',0) + call sync66(iwave,ntrperiod*12000,mode66,nsps,nfqso,ntol,xdt,f0,snr1) + call timer('sync66 ',1) + ! Downsample to give complex data at 6000 S/s fac=2.0/nfft1 c0=fac*iwave(1:nfft1) @@ -122,21 +126,18 @@ contains c0(nfft2/2+1:nfft2)=0. !Zero the top half c0(0)=0.5*c0(0) call four2a(c0,nfft2,1,1,1) !Inverse c2c FFT - call timer('sync66 ',0) - call sync66(iwave,ntrperiod*12000,mode66,nsps,nfqso,ntol,xdt,f0,snr1) - call timer('sync66 ',1) jpk=(xdt+0.5)*6000 - 384 !### Empirical ### if(ntrperiod.ge.60) jpk=(xdt+1.0)*6000 - 384 !### TBD ### if(jpk.lt.0) jpk=0 a=0. - a(1)=-(f0 + 2.0*mode66*baud) !Data tones start 2*mode66 bins higher + a(1)=-(f0 + mode66*baud) !Data tones start mode66 bins higher call twkfreq(c0,c0,ntrperiod*6000,6000.0,a) xdt=jpk/6000.0 - 0.5 LL=64*(mode66+2) NN=63 - call spec66(c0(jpk:),nsps/2,s3,LL,NN) !Compute the synchronized symbol spectra + call spec66(c0(jpk:),nsps/2,s3,LL,NN) !Compute synchronized symbol spectra do j=1,63 !Normalize to symbol baseline call pctile(s3(:,j),LL,40,base) diff --git a/lib/spec66.f90 b/lib/spec66.f90 index 0d049fddb..454f1c18d 100644 --- a/lib/spec66.f90 +++ b/lib/spec66.f90 @@ -6,12 +6,19 @@ subroutine spec66(c0,nsps,s3,LL,NN) complex cs(0:nsps-1) !Complex symbol spectrum real s3(LL,NN) !Synchronized symbol spectra real xbase0(LL),xbase(LL) !Work arrays + integer isync(22) !Indices of sync symbols + data isync/1,9,12,13,15,22,23,26,27,33,35,38,46,50,55,60,62,66,69,74,76,85/ fac=1.0/nsps - ja=-nsps - do j=1,NN - ja=ja+nsps - if(mod(ja/nsps,4).eq.0) ja=ja+nsps + j=0 + n=1 + do k=1,84 + if(k.eq.isync(n)) then + n=n+1 + cycle + endif + j=j+1 + ja=(k-1)*nsps jb=ja+nsps-1 cs=fac*c0(ja:jb) call four2a(cs,nsps,1,-1,1) !c2c FFT to frequency diff --git a/lib/sync66.f90 b/lib/sync66.f90 index bbb72c997..48ed279b4 100644 --- a/lib/sync66.f90 +++ b/lib/sync66.f90 @@ -1,36 +1,32 @@ subroutine sync66(iwave,nmax,mode66,nsps,nfqso,ntol,xdt,f0,snr1) - parameter (NSTEP=4) !Quarter-symbol steps - integer*2 iwave(0:nmax-1) !Raw data - integer b11(11) !Barker 11 code - integer ijpk(2) !Indices i and j at peak of sync_sig - real, allocatable :: s1(:,:) !Symbol spectra - real, allocatable :: x(:) !Work array; demoduated 2FSK sync signal - real sync(4*85) !sync vector - real sync_sig(-64:64,-15:15) + parameter (NSTEP=4) !Quarter-symbol steps + integer*2 iwave(0:nmax-1) !Raw data + integer isync(22) !Indices of sync symbols + integer ijpk(2) !Indices i and j at peak of ccf + real, allocatable :: s1(:,:) !Symbol spectra, quarter-symbol steps + real sync(85) !sync vector + real ccf(-64:64,-15:15) complex, allocatable :: c0(:) !Complex spectrum of symbol - data b11/1,1,1,0,0,0,1,0,0,1,0/ !Barker 11 code + data isync/1,9,12,13,15,22,23,26,27,33,35,38,46,50,55,60,62,66,69,74,76,85/ data sync(1)/99.0/ save sync nfft=2*nsps df=12000.0/nfft istep=nsps/NSTEP - iz=5000.0/df !Uppermost frequency bin, at 5000 Hz + iz=5000.0/df !Uppermost frequency bin, at 5000 Hz txt=85.0*nsps/12000.0 - jz=(txt+1.0)*12000.0/istep !Number of quarter-symbol steps + jz=(txt+1.0)*12000.0/istep !Number of quarter-symbol steps if(nsps.ge.7680) jz=(txt+2.0)*12000.0/istep !For TR 60 s and higher + allocate(s1(iz,jz)) - allocate(x(jz)) - allocate(c0(0:nsps)) + allocate(c0(0:nfft-1)) if(sync(1).eq.99.0) then - sync=0. + sync=-22.0/63.0 !Sync OFF do k=1,22 - kk=k - if(kk.gt.11) kk=k-11 - i=4*NSTEP*(k-1) + 1 - sync(i)=2.0*b11(kk) - 1.0 + sync(isync(k))=1.0 !Sync ON enddo endif @@ -45,7 +41,7 @@ subroutine sync66(iwave,nmax,mode66,nsps,nfqso,ntol,xdt,f0,snr1) k=k+1 c0(k)=fac*cmplx(xx,yy) enddo - c0(k+1:nfft/2)=0. + c0(k+1:)=0. call four2a(c0,nfft,1,-1,0) !r2c FFT do i=1,iz s1(i,j)=real(c0(i))**2 + aimag(c0(i))**2 @@ -58,14 +54,18 @@ subroutine sync66(iwave,nmax,mode66,nsps,nfqso,ntol,xdt,f0,snr1) s1max=20.0 ! Apply AGC - do j=1,jz - x(j)=maxval(s1(i0-64:i0+192,j)) - if(x(j).gt.s1max) s1(i0-64:i0+192,j)=s1(i0-64:i0+192,j)*s1max/x(j) - enddo +! do j=1,jz +! smax=maxval(s1(i0-64:i0+192,j)) +! if(smax.gt.s1max) s1(i0-64:i0+192,j)=s1(i0-64:i0+192,j)*s1max/smax +! enddo - sync_sig=0. +! do i=1,iz +! write(60,3060) i,i*df,sum(s1(i,1:jz)) +!3060 format(i6,f10.3,e12.3) +! enddo + + ccf=0. ia=min(64,nint(ntol/df)) - dt4=nsps/(NSTEP*12000.0) !duration of 1/4 symbol jadd=11 if(nsps.ge.3600) jadd=7 @@ -74,31 +74,35 @@ subroutine sync66(iwave,nmax,mode66,nsps,nfqso,ntol,xdt,f0,snr1) if(nsps.ge.41472) jadd=1 do i=-ia,ia - x=s1(i0+2*mode66+i,:)-s1(i0+i,:) !Do the 2FSK demodulation do lag=-15,15 - do k=1,22 - n=4*NSTEP*(k-1) + 1 + do k=1,85 + n=NSTEP*(k-1) + 1 j=n+lag+jadd - if(j.ge.1 .and. j.le.jz) sync_sig(i,lag)=sync_sig(i,lag) + sync(n)*x(j) + if(j.ge.1 .and. j.le.jz) then + ccf(i,lag)=ccf(i,lag) + sync(k)*s1(i0+i,j) + endif enddo enddo enddo - ijpk=maxloc(sync_sig) - ii=ijpk(1)-65 - jj=ijpk(2)-16 +! do i=-64,64 +! write(61,3061) i,ccf(i,jpk) +!3061 format(i5,e12.3) +! enddo +! do j=-15,15 +! write(62,3061) j,ccf(ipk,j) +! enddo -! Use peakup() here? - f0=nfqso + ii*df - jdt=jj - -! j0=0.5/dt4 -! if(nsps.ge.7680) j0=1.0/dt4 - - tsym=nsps/12000.0 - xdt=jdt*tsym/4.0 - - snr1=maxval(sync_sig)/22.0 + ijpk=maxloc(ccf) + ipk=ijpk(1)-65 + jpk=ijpk(2)-16 + dt4=nsps/(NSTEP*12000.0) !1/4 of symbol duration + if(nsps.ge.7680) j0=1.0/dt4 + f0=nfqso + ipk*df + xdt=jpk*dt4 + snr1=maxval(ccf)/22.0 +! write(*,3100) ipk,jpk,xdt,f0,snr1 +!3100 format(2i5,f7.2,2f10.2) return end subroutine sync66