w4kek
/
WSJT-X
mirror of https://github.com/saitohirga/WSJT-X.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Code cleanup for synctest and wideband2_sync.

This commit is contained in:
Joe Taylor 2021-05-12 14:30:53 -04:00
parent 150c6b179c
commit cd978bae4e
2 changed files with 139 additions and 96 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

77
map65/libm65/synctest.f90
View File

@ -1,13 +1,19 @@
program synctest
! Program to test an algorithm for detecting sync signals for both
! JT65 and Q65-60x signals and rejecting birdies in MAP65 data.
! The important work is done in module wideband_sync.
use timer_module, only: timer
use timer_impl, only: init_timer, fini_timer
use wideband2_sync
parameter (NFFT=32768)
parameter (MAX_CANDIDATES=20)
real ss(4,322,NFFT),savg(4,NFFT)
real candidate(MAX_CANDIDATES,5) !snr1,f0,xdt0,ipol,flip
! real candidate(MAX_CANDIDATES,5) !snr1,f0,xdt0,ipol,flip
character*8 arg
type(candidate) :: cand(MAX_CANDIDATES)
nargs=iargc()
if(nargs.ne.5) then
print*,'Usage: synctest iutc nfa nfb nts_jt65 nts_q65'
@ -36,17 +42,13 @@ program synctest
call timer('synctest',0)
call timer('get_cand',0)
call get_candidates(ss,savg,nfa,nfb,nts_jt65,nts_q65,candidate,ncand)
call get_candidates(ss,savg,nfa,nfb,nts_jt65,nts_q65,cand,ncand)
call timer('get_cand',1)
do k=1,ncand
snr1=candidate(k,1)
f0=candidate(k,2)
xdt0=candidate(k,3)
ipol=nint(candidate(k,4))
iflip=nint(candidate(k,5))
write(*,1010) k,snr1,f0,f0+77,xdt0,ipol,iflip
1010 format(i3,4f10.3,2i3,f8.3)
write(*,1010) k,cand(k)%snr,cand(k)%f,cand(k)%f+77,cand(k)%xdt, &
cand(k)%ipol,cand(k)%iflip
1010 format(i3,4f10.3,2i3)
enddo
998 call timer('synctest',1)
@ -54,58 +56,3 @@ program synctest
call fini_timer()
999 end program synctest
subroutine get_candidates(ss,savg,nfa,nfb,nts_jt65,nts_q65,candidate,ncand)
use wideband2_sync
parameter (MAX0=1000,MAX_CANDIDATES=20)
real ss(4,322,NFFT),savg(4,NFFT)
real candidate(MAX_CANDIDATES,5) !snr1,f0,xdt0,ipol,flip
integer indx(NFFT)
logical skip
call wb2_sync(ss,savg,nfa,nfb)
df3=96000.0/NFFT
ia=nint(1000*nfa/df3)
ib=nint(1000*nfb/df3)
iz=ib-ia+1
do i=ia,ib
f0=0.001*i*df3
write(13,1000) f0,sync_dat(i,2:5)
1000 format(3f10.3,2f5.0)
enddo
call indexx(sync_dat(ia:ib,2),iz,indx)
k=0
do i=1,MAX0
j=indx(iz+1-i) + ia - 1
f0=0.001*j*df3
snr1=sync_dat(j,2)
if(snr1.lt.5.0) exit
flip=sync_dat(j,5)
if(flip.ne.0.0 .and. nts_jt65.eq.0) cycle
if(flip.eq.0.0 .and. nts_q65.eq.0) cycle
if(sync_dat(j,6).lt.3.0) cycle
bw=nts_q65*108.0
if(flip.ne.0) bw=nts_jt65*177.0
skip=.false.
do n=1,k
diffhz=1000.0*(f0-candidate(n,2))
if(diffhz.gt.-10.0 .and. diffhz.lt.bw) skip=.true.
enddo
if(skip) cycle
k=k+1
xdt0=sync_dat(j,3)
ipol=nint(sync_dat(j,4))
candidate(k,1)=snr1 !snr1
candidate(k,2)=f0 !f0
candidate(k,3)=xdt0
candidate(k,4)=ipol
candidate(k,5)=flip
if(k.ge.MAX_CANDIDATES) exit
enddo
ncand=k
return
end subroutine get_candidates

158
map65/libm65/wideband2_sync.f90
View File

@ -1,25 +1,122 @@
module wideband2_sync
parameter (NFFT=32768)
integer isync(22)
integer nkhz_center
real sync_dat(NFFT,6) !fkhz, ccfmax, xdt, ipol, flip
real savg_med(4)
type candidate
real :: snr !Relative S/N of sync detection
real :: f !Freq of sync tone, 0 to 96000 Hz
real :: xdt !DT of matching sync pattern, -1.0 to +4.0 s
integer :: ipol !Polarization angle, 1 to 4 ==> 0, 45, 90, 135 deg
integer :: iflip !Sync type: JT65 = +/- 1, Q65 = 0
end type candidate
type sync_dat
real :: ccfmax
real :: xdt
integer :: ipol
integer :: iflip
logical :: birdie
end type sync_dat
contains
subroutine wb2_sync(ss,savg,nfa,nfb)
subroutine get_candidates(ss,savg,nfa,nfb,nts_jt65,nts_q65,cand,ncand)
! Search symbol spectra ss() over frequency range nfa to nfb (in kHz) for
! JT65 and Q65 sync patterns. The nts_* variables are the submode tone
! spacings: 1 2 4 8 16 for A B C D E. Birdies are detected and
! excised. Candidates are returned in the structure array cand().
parameter (NFFT=32768)
parameter (MAX_PEAKS=300,MAX_CANDIDATES=20)
real ss(4,322,NFFT),savg(4,NFFT)
real pavg(-20:20)
integer indx(NFFT)
logical skip
type(candidate) :: cand(MAX_CANDIDATES)
type(sync_dat) :: sync(NFFT)
do j=322,1,-1 !Find end of data in ss()
if(sum(ss(1,j,1:NFFT)).gt.0.0) exit
enddo
jz=j
call wb2_sync(ss,savg,jz,nfa,nfb,sync)
tstep=2048.0/11025.0 !0.185760 s: 0.5*tsym_jt65, 0.3096*tsym_q65
df3=96000.0/NFFT
ia=nint(1000*nfa/df3) + 1
ib=nint(1000*nfb/df3) + 1
iz=ib-ia+1
call indexx(sync(ia:ib)%ccfmax,iz,indx) !Sort by relative snr
k=0
do i=1,MAX_PEAKS
n=indx(iz+1-i) + ia - 1
f0=0.001*(n-1)*df3
snr1=sync(n)%ccfmax
if(snr1.lt.5.0) exit
flip=sync(n)%iflip
if(flip.ne.0.0 .and. nts_jt65.eq.0) cycle
if(flip.eq.0.0 .and. nts_q65.eq.0) cycle
if(sync(n)%birdie) cycle
! Test for signal outside of TxT range and set bw for this signal type
j1=(sync(n)%xdt + 1.0)/tstep - 1.0
j2=(sync(n)%xdt + 52.0)/tstep + 1.0
bw=nts_q65*108.0
if(flip.ne.0) then
bw=nts_jt65*177.0
j2=(sync(n)%xdt + 47.811)/tstep + 1.0
endif
ipol=sync(n)%ipol
pavg=0.
do j=1,j1
pavg=pavg + ss(ipol,j,n-20:n+20)
enddo
do j=j2,jz
pavg=pavg + ss(ipol,j,n-20:n+20)
enddo
jsum=j1 + (jz-j2+1)
pmax=maxval(pavg(-2:2)) !### Why not just pavg(0) ?
base=(sum(pavg)-pmax)/jsum
pmax=pmax/base
if(pmax.gt.5.0) cycle
skip=.false.
do m=1,k !Skip false syncs within signal bw
diffhz=1000.0*(f0-cand(m)%f)
if(diffhz.gt.-10.0 .and. diffhz.lt.bw) skip=.true.
enddo
if(skip) cycle
k=k+1
cand(k)%snr=snr1
cand(k)%f=f0
cand(k)%xdt=sync(n)%xdt
cand(k)%ipol=sync(n)%ipol
cand(k)%iflip=nint(flip)
if(k.ge.MAX_CANDIDATES) exit
enddo
ncand=k
return
end subroutine get_candidates
subroutine wb2_sync(ss,savg,jz,nfa,nfb,sync)
! Compute "orange sync curve" using the Q65 sync pattern
parameter (NFFT=32768)
parameter (LAGMAX=30)
real ss(4,322,NFFT)
real savg(4,NFFT)
real savg_med(4)
logical first
integer isync0(22)
integer isync(22)
integer jsync0(63),jsync1(63)
type(sync_dat) :: sync(NFFT)
! 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 isync/1,9,12,13,15,22,23,26,27,33,35,38,46,50,55,60,62,66,69,74,76,85/
data jsync0/ &
1, 4, 5, 9, 10, 11, 12, 13, 14, 16, 18, 22, 24, 25, 28, 32, &
33, 34, 37, 38, 39, 40, 42, 43, 45, 46, 47, 48, 52, 53, 55, 57, &
@ -31,13 +128,13 @@ subroutine wb2_sync(ss,savg,nfa,nfb)
72, 74, 75, 76, 77, 78, 79, 82, 83, 84, 85, 86, 87, 88, 91, 93, &
94, 96, 99,101,103,105,106,109,110,112,113,115,116,117,118/
data first/.true./
save first
save first,isync,jsync0,jsync1
tstep=2048.0/11025.0 !0.185760 s: 0.5*tsym_jt65, 0.3096*tsym_q65
if(first) then
fac=0.6/tstep
do i=1,22 !Expand the Q65 sync stride
isync(i)=nint((isync0(i)-1)*fac) + 1
isync(i)=nint((isync(i)-1)*fac) + 1
enddo
do i=1,63
jsync0(i)=2*(jsync0(i)-1) + 1
@ -46,21 +143,15 @@ subroutine wb2_sync(ss,savg,nfa,nfb)
first=.false.
endif
do j=322,1,-1
if(sum(ss(1,j,1:NFFT)).gt.0.0) exit
enddo
jz=j
df3=96000.0/NFFT
ia=nint(1000*nfa/df3) !Flat frequency range for WSE converters
ib=nint(1000*nfb/df3)
ia=nint(1000*nfa/df3) + 1 !Flat frequency range for WSE converters
ib=nint(1000*nfb/df3) + 1
do i=1,4
call pctile(savg(i,ia:ib),ib-ia+1,50,savg_med(i))
enddo
do i=ia,ib
write(14,3014) 0.001*i*df3,savg(1:4,i)
write(14,3014) 0.001*(i-1)*df3,savg(1:4,i)
3014 format(5f10.3)
enddo
@ -115,18 +206,23 @@ subroutine wb2_sync(ss,savg,nfa,nfb)
enddo ! lag
enddo !ipol
fkhz=0.001*i*df3 + nkhz_center - 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
sync_dat(i,5)=flip
sync_dat(i,6)=ccfmax/(savg(ipolbest,i)/savg_med(ipolbest))
enddo
sync(i)%ccfmax=ccfmax
sync(i)%xdt=lagbest*tstep-1.0
sync(i)%ipol=ipolbest
sync(i)%iflip=flip
sync(i)%birdie=.false.
if(ccfmax/(savg(ipolbest,i)/savg_med(ipolbest)).lt.3.0) sync(i)%birdie=.true.
enddo ! i (frequency bin)
call pctile(sync_dat(ia:ib,2),ib-ia+1,50,base)
sync_dat(ia:ib,2)=sync_dat(ia:ib,2)/base
! do i=ia,ib
! write(15,3015) 0.001*(i-1)*df3,sync(i)%ccfmax,sync(i)%xdt,sync(i)%ipol, &
! sync(i)%iflip,sync(i)%birdie
!3015 format(3f10.3,2i6,L5)
! enddo
call pctile(sync(ia:ib)%ccfmax,ib-ia+1,50,base)
sync(ia:ib)%ccfmax=sync(ia:ib)%ccfmax/base
! print*,base
return
end subroutine wb2_sync