diff --git a/CMakeLists.txt b/CMakeLists.txt index 8bb35d268..583eef077 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -435,6 +435,7 @@ set (wsjt_FSRCS lib/softsym9f.f90 lib/softsym9w.f90 lib/shell.f90 + lib/spec64.f90 lib/spec9f.f90 lib/stdmsg.f90 lib/subtract65.f90 diff --git a/lib/qra64a.f90 b/lib/qra64a.f90 index 0892e04ae..9dac72745 100644 --- a/lib/qra64a.f90 +++ b/lib/qra64a.f90 @@ -1,4 +1,4 @@ -subroutine qra64a(dd0,nutc,nf1,nf2,nfqso,ntol,mycall_12,hiscall_12, & +subroutine qra64a(dd,nutc,nf1,nf2,nfqso,ntol,mycall_12,hiscall_12, & hisgrid_6,sync,nsnr,dtx,nfreq,decoded,nft) use packjt @@ -8,106 +8,20 @@ subroutine qra64a(dd0,nutc,nf1,nf2,nfqso,ntol,mycall_12,hiscall_12, & character*6 mycall,hiscall,hisgrid_6 character*4 hisgrid logical ltext - integer*8 count0,count1,clkfreq + complex c00(0:360000) !Complex spectrum of dd() + complex c0(0:360000) !Complex spectrum of dd() +! integer*8 count0,count1,clkfreq integer icos7(0:6) integer dat4(12) - real dd0(NMAX),dd(NMAX) - real s(NZ) - real savg(NZ) - real blue(0:25) - real red(NZ) - real ss(NZ,194) - real ccf(NZ,0:25) + real a(3) + real dd(NMAX) real s3(0:63,1:63) - real s3a(0:63,1:63) - data icos7/2,5,6,0,4,1,3/ !Costas 7x7 pattern + data icos7/2,5,6,0,4,1,3/ !Costas 7x7 pattern data nc1z/-1/,nc2z/-1/,ng2z/-1/ -! common/qra64com/ss(NZ,194),ccf(NZ,0:25),s3(0:63,1:63),s3a(0:63,1:63) save -! write(60) dd0 - - decoded=' ' nft=99 nsnr=-30 - nsps=6912 - istep=nsps/2 - nsteps=52*12000/istep - 2 - df=12000.0/NFFT - call spec64(dd0,-1,s,savg,ss) - fa=max(nf1,nfqso-ntol) - fb=min(nf2,nfqso+ntol) - ia=max(1,nint(fa/df)) - ib=min(NZ,nint(fb/df)) - call pctile(savg(ia),ib-ia+1,45,base) - savg=savg/base - 1.0 - ss=ss/base - - red=-99. - fac=1.0/sqrt(21.0) - sync=0. - do if0=ia,ib - do j=0,25 - t=-3.0 - do n=0,6 - i=if0 + 2*icos7(n) - t=t + ss(i,1+2*n+j) + ss(i,1+2*n+j+78) + ss(i,1+2*n+j+154) - enddo - ccf(if0,j)=fac*t - if(ccf(if0,j).gt.red(if0)) then - red(if0)=ccf(if0,j) - if(red(if0).gt.sync) then - sync=red(if0) - f0=if0*df -! dtx=j*istep/12000.0 - 1.0 - i0=if0 - j0=j - endif - endif - enddo - enddo - - write(17) ia,ib,red(ia:ib) - close(17) - - if0=nint(f0/df) - nfreq=nint(f0) - blue(0:25)=ccf(if0,0:25) - dj=0. - if(j0.ge.1 .and. j0.le.24) call peakup(blue(j0-1),blue(j0),blue(j0+1),dj) - xpk=j0 + dj - dtx=xpk*istep/12000.0 - 1.0 - i=nint(dj*istep) - if(i.lt.0) then - dd(1-i:NMAX)=dd0(1:NMAX+i) - dd(1:-i)=0.0 - else - dd(1:NMAX-i)=dd0(1+i:NMAX) - dd(NMAX-i+1:NMAX)=0.0 - endif - -! Recompute the symbol spectra, this time aligned for best DT estimate. - ss=0. - call spec64(dd,j0,s,savg,ss) - - do i=0,63 !Copy symbol spectra into s3() - k=i0 + 2*i - jj=j0+13 - do j=1,63 - jj=jj+2 - s3(i,j)=ss(k,jj) - if(j.eq.32) jj=jj+14 !Skip over the middle Costas array - enddo - enddo - - if(sync.gt.1.0) snr1=10.0*log10(sync) - 39.0 - nsnr=nint(snr1) -! write(*,5001) dtx,nint(f0),0,snr1 -!5001 format(f6.3,2i6,f7.1) - maxf1=10 -! call sync64(dd0,nf1,nf2,nfqso,ntol,maxf1,dtx,f0,kpk,snr,s3a) -! write(*,5001) dtx,nint(f0),kpk,snr - mycall=mycall_12(1:6) !### May need fixing ### hiscall=hiscall_12(1:6) hisgrid=hisgrid_6(1:4) @@ -124,59 +38,45 @@ subroutine qra64a(dd0,nutc,nf1,nf2,nfqso,ntol,mycall_12,hiscall_12, & ng2z=ng2 endif + maxf1=5 +! maxf1=0 + call sync64(dd,nf1,nf2,nfqso,ntol,maxf1,dtx,f0,jpk,kpk,snr1,c00) + +!### + npts2=216000 snr2=-99. naptype=4 - call system_clock(count0,clkfreq) - call qra64_dec(s3,nc1,nc2,ng2,naptype,0,dat4,snr2,irc) - if(irc.ge.0) then - call unpackmsg(dat4,decoded) !Unpack the user message - call fmtmsg(decoded,iz) - nft=100 + irc - nsnr=nint(snr2) - else - snr2=0. - endif - call system_clock(count1,clkfreq) - tsec=float(count1-count0)/float(clkfreq) - write(78,3900) nutc,sync,snr1,snr2,dtx,nfreq,1,irc,tsec,decoded -3900 format(i4.4,3f6.1,f6.2,i5,i2,i3,f6.3,1x,a22) - flush(78) + do itry0=1,3 + idf0=itry0/2 + if(mod(itry0,2).eq.0) idf0=-idf0 + a(1)=-(f0+0.248*idf0) + nfreq=nint(-a(1)) + a(3)=0. + do itry1=1,1 + idf1=itry1/2 + if(mod(itry1,2).eq.0) idf1=-idf1 + a(2)=0.5*idf1 + call twkfreq(c00,c0,npts2,4000.0,a) + call spec64(c0,npts2,jpk,s3) + call qra64_dec(s3,nc1,nc2,ng2,naptype,0,dat4,snr2,irc) + decoded=' ' + if(irc.ge.0) then + call unpackmsg(dat4,decoded) !Unpack the user message + call fmtmsg(decoded,iz) + nft=100 + irc + nsnr=nint(snr2) + else + snr2=0. + endif + write(78,3900) nutc,snr1,snr2,dtx,nfreq,idf0,idf1,irc,decoded +3900 format(i4.4,2f6.1,f6.2,i5,3i3,1x,a22) + flush(78) +! write(*,3006) idf0,idf1,nfreq,jpk,kpk,irc,decoded +!3006 format(2i4,i6,i7,2i4,2x,a22) + if(irc.ge.0) go to 900 + enddo + enddo +900 continue return end subroutine qra64a - -subroutine spec64(dd,j0,s,savg,ss) - - parameter (NFFT=2*6912,NH=NFFT/2,NZ=5760) - real dd(60*12000) - real s(NZ) - real savg(NZ) - real ss(NZ,194) - real x(NFFT) - complex cx(0:NH) - equivalence (x,cx) - - nsps=6912 - istep=nsps/2 - nsteps=52*12000/istep - 2 - ia=1-istep - savg=0. - df=12000.0/NFFT - mj0=mod(j0,2) - do j=1,nsteps - ia=ia+istep - if(j0.ge.0 .and. mod(j,2).eq.mj0) cycle !Skip half of FFTs in 2nd pass - ib=ia+nsps-1 - x(1:nsps)=1.2e-4*dd(ia:ib) - x(nsps+1:)=0.0 - call four2a(x,nfft,1,-1,0) !r2c FFT - do i=1,NZ - s(i)=real(cx(i))**2 + aimag(cx(i))**2 - enddo - ss(1:NZ,j)=s - savg=savg+s - enddo - savg=savg/nsteps - - return -end subroutine spec64 diff --git a/lib/sync64.f90 b/lib/sync64.f90 index 6c2c5b974..63fe2e029 100644 --- a/lib/sync64.f90 +++ b/lib/sync64.f90 @@ -1,11 +1,9 @@ -subroutine sync64(dd,nf1,nf2,nfqso,ntol,maxf1,dtx,f0,kpk,snrdb,s3a) +subroutine sync64(dd,nf1,nf2,nfqso,ntol,maxf1,dtx,f0,jpk,kpk,snrdb,c0) parameter (NMAX=60*12000) !Max size of raw data at 12000 Hz parameter (NSPS=2304) !Samples per symbol at 4000 Hz parameter (NSPC=7*NSPS) !Samples per Costas array real dd(NMAX) !Raw data - real x(672000) !Up to 56 s at 12000 Hz - real s3a(0:63,1:63) !Synchronized symbol spectra real s1(0:NSPC-1) !Power spectrum of Costas 1 real s2(0:NSPC-1) !Power spectrum of Costas 2 real s3(0:NSPC-1) !Power spectrum of Costas 3 @@ -14,12 +12,11 @@ subroutine sync64(dd,nf1,nf2,nfqso,ntol,maxf1,dtx,f0,kpk,snrdb,s3a) integer icos7(0:6) !Costas 7x7 tones integer ipk0(1) complex cc(0:NSPC-1) !Costas waveform - complex c0(0:336000) !Complex spectrum of dd() + complex c0(0:360000) !Complex spectrum of dd() complex c1(0:NSPC-1) !Complex spectrum of Costas 1 complex c2(0:NSPC-1) !Complex spectrum of Costas 2 complex c3(0:NSPC-1) !Complex spectrum of Costas 3 logical first - equivalence (x,c0) data icos7/2,5,6,0,4,1,3/ !Costas 7x7 tone pattern data first/.true./ save @@ -46,7 +43,10 @@ subroutine sync64(dd,nf1,nf2,nfqso,ntol,maxf1,dtx,f0,kpk,snrdb,s3a) nfft2=nfft1/3 df1=12000.0/nfft1 fac=2.0/nfft1 - x=fac*dd(1:nfft1) +! x=fac*dd(1:nfft1) + do i=0,nfft1/2 + c0(i)=fac*cmplx(dd(1+2*i),dd(2+2*i)) + enddo call four2a(c0,nfft1,1,-1,0) !Forward r2c FFT call four2a(c0,nfft2,1,1,1) !Inverse c2c FFT; c0 is analytic sig npts2=npts0/3 !Downsampled complex data length @@ -94,7 +94,8 @@ subroutine sync64(dd,nf1,nf2,nfqso,ntol,maxf1,dtx,f0,kpk,snrdb,s3a) s=(maxval(s0(ia:ib))-ave)/rms if(s.gt.snr) then jpk=j1 - s0a=s0 +! s0a=(s0-ave)/rms + s0a=s0/rms snr=s dtx=jpk/4000.0 - 1.0 ipk0=maxloc(s0(ia:ib)) @@ -110,11 +111,20 @@ subroutine sync64(dd,nf1,nf2,nfqso,ntol,maxf1,dtx,f0,kpk,snrdb,s3a) ! ka=kpk ! kb=kpk enddo - snrdb=10.0*log10(snr)-39.0 !### -! Now use tweak on the c0() array, then do nfft4=NSPS FFTs to get s3a(), -! the properly synchronized symbol spectra. +! rewind 73 +! do i=ia,ib +! write(73,3201) i*df3,s0a(i),db(s0a(i)) +!3201 format(3f10.3) +! enddo +! flush(73) +!### + + write(17) ia,ib,s0a(ia:ib) !Save red() + close(17) + + snrdb=10.0*log10(snr)-39.0 return end subroutine sync64 diff --git a/plotter.cpp b/plotter.cpp index 0c0add863..1942a4661 100644 --- a/plotter.cpp +++ b/plotter.cpp @@ -132,7 +132,7 @@ void CPlotter::draw(float swide[], bool bScroll, bool bRed) m_line++; if(m_mode=="QRA64" and bRed) { - double df_qra64=12000.0/(2*6912); + double df_qra64=4000.0/(7*2304); int j0,j1; int k=0; float smax,y3max=0; @@ -149,7 +149,7 @@ void CPlotter::draw(float swide[], bool bScroll, bool bRed) } float fac=0.8/qMax(y3max,10.0f); for(int i=1; i0) { + if(FreqfromX(i)>=m_ia*df_qra64 and FreqfromX(i)