Compute analog of "orange sync curve" in MAP65; use it to determind ipol for a decode candidate.

This commit is contained in:
Joe Taylor 2021-04-29 14:35:19 -04:00
parent 88fa43cf77
commit 7f9da0ba90
6 changed files with 128 additions and 42 deletions

View File

@ -1,4 +1,8 @@
set (libm65_FSRCS
# Modules come first:
wideband_sync.f90
# Non-module Fortran routines:
afc65b.f90
astro.f90
astro0.f90

View File

@ -1,8 +1,8 @@
subroutine decode0(dd,ss,savg,nstandalone)
use wideband_sync
use timer_module, only: timer
parameter (NSMAX=60*96000)
parameter (NFFT=32768)
real*4 dd(4,NSMAX),ss(4,322,NFFT),savg(4,NFFT)
real*8 fcenter
@ -17,9 +17,12 @@ subroutine decode0(dd,ss,savg,nstandalone)
data neme0/-99/,mcall3b/1/
save
! write(60) ss,savg
call timer('decode0 ',0)
call timer('wb_sync ',0)
call wb_sync(ss,savg)
call timer('wb_sync ',1)
if(newdat.ne.0) then
nz=52*96000
hist=0

View File

@ -107,7 +107,6 @@ program m65
nch=2
if(nxpol.eq.1) nch=4
! if(ifile.eq.ifile1) call timer('m65 ',0)
do irec=1,9999999
read(10,end=10) i2
do i=1,NREAD,nch
@ -139,9 +138,8 @@ program m65
rejecty,pxdb,pydb,ssz5a,nkhz,ihsym,nzap,slimit,lstrong)
call timer('symspec ',1)
nhsym0=nhsym
! if(ihsym.ge.278) go to 10
endif
enddo
enddo ! irec
10 continue
if(iqadjust.ne.0) write(*,3002) rejectx,rejecty
@ -149,7 +147,7 @@ program m65
nutc=nutc0
nstandalone=1
call decode0(dd,ss,savg,nstandalone)
enddo
enddo ! ifile
call timer('m65 ',1)
call timer('m65 ',101)

View File

@ -5,10 +5,11 @@ subroutine map65a(dd,ss,savg,newdat,nutc,fcenter,ntol,idphi,nfa,nfb, &
! Processes timf2 data from Linrad to find and decode JT65 signals.
use wideband_sync
use timer_module, only: timer
parameter (MAXMSG=1000) !Size of decoded message list
parameter (NSMAX=60*96000)
parameter (NFFT=32768)
real dd(4,NSMAX)
real*4 ss(4,322,NFFT),savg(4,NFFT)
real tavg(-50:50) !Temp for finding local base level

View File

@ -2,7 +2,9 @@ subroutine q65b(nutc,fcenter,nfcal,nfsample,ikhz,mousedf,ntol,xpol, &
mycall0,hiscall0,hisgrid,mode_q65)
use q65_decode
use wideband_sync
use timer_module, only: timer
parameter (MAXFFT1=5376000) !56*96000
parameter (MAXFFT2=336000) !56*6000 (downsampled by 1/16)
parameter (NMAX=60*12000)
@ -10,6 +12,7 @@ subroutine q65b(nutc,fcenter,nfcal,nfsample,ikhz,mousedf,ntol,xpol, &
complex ca(MAXFFT1),cb(MAXFFT1) !FFTs of raw x,y data
complex cx(0:MAXFFT2-1),cy(0:MAXFFT2-1),cz(0:MAXFFT2)
logical xpol
integer ipk1(1)
real*8 fcenter
character*12 mycall0,hiscall0
character*12 mycall,hiscall
@ -67,43 +70,50 @@ subroutine q65b(nutc,fcenter,nfcal,nfsample,ikhz,mousedf,ntol,xpol, &
! 96000 5376000 0.017857143 336000 6000.000
! 95238 5120000 0.018601172 322560 5999.994
! Get best ipol by referring to the "orange sync curve".
ff=ikhz+0.001*(mousedf+1270.459) !supposed freq of sync tone
ifreq=nint(1000.0*(ff-77.0)*32768.0/96000.0) !Freq index into ss(4,322,32768)
dff=96000.0/32768.0
ia=nint(ifreq-ntol/dff)
ib=nint(ifreq+ntol/dff)
ipk1=maxloc(sync_dat(ia:ib,2))
ipk=ia+ipk1(1)-1
ipol=nint(sync_dat(ipk,4))
nsnr1=-99
npol=1
if(xpol) npol=4
do ipol=1,npol
if(ipol.eq.1) cz(0:MAXFFT2-1)=cx
if(ipol.eq.2) cz(0:MAXFFT2-1)=0.707*(cx+cy)
if(ipol.eq.3) cz(0:MAXFFT2-1)=cy
if(ipol.eq.4) cz(0:MAXFFT2-1)=0.707*(cx-cy)
cz(MAXFFT2)=0.
if(ipol.eq.1) cz(0:MAXFFT2-1)=cx
if(ipol.eq.2) cz(0:MAXFFT2-1)=0.707*(cx+cy)
if(ipol.eq.3) cz(0:MAXFFT2-1)=cy
if(ipol.eq.4) cz(0:MAXFFT2-1)=0.707*(cx-cy)
cz(MAXFFT2)=0.
!Transform to time domain (real), fsample=12000 Hz
call four2a(cz,2*nfft2,1,1,-1)
do i=0,nfft2-1
j=nfft2-1-i
iwave(2*i+1)=nint(real(cz(j)))
iwave(2*i+2)=nint(aimag(cz(j)))
enddo
iwave(2*nfft2+1:)=0
nsubmode=mode_q65-1
nfa=300
nfb=2883
nfqso=1000 + mousedf
newdat=1
nagain=0
nsnr0=-99 !Default snr for no decode
call timer('mmdec ',0)
call map65_mmdec(nutc,iwave,nsubmode,nfa,nfb,nfqso,ntol,newdat,nagain, &
mycall,hiscall,hisgrid)
call timer('mmdec ',1)
if(nsnr0.gt.nsnr1) then
nsnr1=nsnr0
xdt1=xdt0
nfreq1=nfreq0
msg1=msg0
cq1=cq0
ipol1=45*ipol
endif
enddo !ipol
call four2a(cz,2*nfft2,1,1,-1)
do i=0,nfft2-1
j=nfft2-1-i
iwave(2*i+1)=nint(real(cz(j)))
iwave(2*i+2)=nint(aimag(cz(j)))
enddo
iwave(2*nfft2+1:)=0
nsubmode=mode_q65-1
nfa=300
nfb=2883
nfqso=1000 + mousedf
newdat=1
nagain=0
nsnr0=-99 !Default snr for no decode
call timer('mmdec ',0)
call map65_mmdec(nutc,iwave,nsubmode,nfa,nfb,nfqso,ntol,newdat,nagain, &
mycall,hiscall,hisgrid)
call timer('mmdec ',1)
if(nsnr0.gt.nsnr1) then
nsnr1=nsnr0
xdt1=xdt0
nfreq1=nfreq0
msg1=msg0
cq1=cq0
ipol1=45*(ipol-1)
endif
nfreq=nfreq1+mousedf-1000
write(line,1020) ikhz,nfreq,ipol1,nutc,xdt1,nsnr1,msg1(1:27),cq1

View File

@ -0,0 +1,70 @@
module wideband_sync
parameter (NFFT=32768)
integer isync(22)
real sync_dat(NFFT,4) !fkhz, ccfmax, xdt, ipol
contains
subroutine wb_sync(ss,savg)
! Compute "orange sync curve" using the Q65 sync pattern
parameter (LAGMAX=30)
real ss(4,322,NFFT)
real savg(4,NFFT)
logical first
integer isync0(22)
! Q65 sync symbols
data isync0/1,9,12,13,15,22,23,26,27,33,35,38,46,50,55,60,62,66,69,74,76,85/
data first/.true./
save first
tstep=2048.0/11025.0
if(first) then
fac=0.6/tstep
do i=1,22 !Expand the sync stride
isync(i)=nint((isync0(i)-1)*fac) + 1
enddo
first=.false.
endif
df=96000.0/NFFT
ia=nint(7000.0/df) !Flat frequency range for WSE converters
ib=nint(89000.0/df)
lagbest=-1
ipolbest=-1
do i=ia,ib
ccfmax=0.
do ipol=1,4
do lag=0,LAGMAX
ccf=0.
do j=1,22
k=isync(j) + lag
ccf=ccf + ss(ipol,k,i)
enddo
ccf=ccf - savg(ipol,i)*22.0/322.0
if(ccf.gt.ccfmax) then
ipolbest=ipol
lagbest=lag
ccfmax=ccf
endif
enddo ! lag
enddo !ipol
fkhz=0.001*i*df + 125.0 - 48.0
xdt=lagbest*tstep-1.0
sync_dat(i,1)=fkhz
sync_dat(i,2)=ccfmax
sync_dat(i,3)=xdt
sync_dat(i,4)=ipolbest
! write(14,3010) i,sync_dat(i,1:3),nint(sync_dat(i,4))
!3010 format(i6,3f10.3,i5)
enddo
return
end subroutine wb_sync
end module wideband_sync