w4kek/WSJT-X
1
0
Fork 0
mirror of https://github.com/saitohirga/WSJT-X.git synced 2025-05-29 04:42:28 -04:00
Code Issues Packages Projects Releases Wiki Activity

Merge branch 'feat-fst280' of bitbucket.org:k1jt/wsjtx into feat-fst280

Browse Source
This commit is contained in:
Nico Palermo/IV3NWV 2020-12-01 03:26:37 +01:00
commit baf427c358
10 changed files with 260 additions and 306 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
CMakeLists.txt
View File

@ -496,10 +496,12 @@ set (wsjt_FSRCS
lib/polyfit.f90 lib/polyfit.f90
lib/prog_args.f90 lib/prog_args.f90
lib/ps4.f90 lib/ps4.f90
lib/q65_sync.f90
lib/qra64a.f90 lib/qra64a.f90
lib/qra_loops.f90 lib/qra_loops.f90
lib/qra/q65/q65_ap.f90 lib/qra/q65/q65_ap.f90
lib/qra/q65/q65_loops.f90 lib/qra/q65/q65_loops.f90
lib/qra/q65/q65_set_list.f90
lib/refspectrum.f90 lib/refspectrum.f90
lib/savec2.f90 lib/savec2.f90
lib/sec0.f90 lib/sec0.f90
@ -529,7 +531,6 @@ set (wsjt_FSRCS
lib/sync4.f90 lib/sync4.f90
lib/sync64.f90 lib/sync64.f90
lib/sync65.f90 lib/sync65.f90
lib/sync_q65.f90
lib/ft4/getcandidates4.f90 lib/ft4/getcandidates4.f90
lib/ft4/get_ft4_bitmetrics.f90 lib/ft4/get_ft4_bitmetrics.f90
lib/ft8/sync8.f90 lib/ft8/sync8.f90

18
lib/decoder.f90
View File

@ -777,8 +777,7 @@ contains
return return
end subroutine fst4_decoded end subroutine fst4_decoded
subroutine q65_decoded (this,nutc,sync,nsnr,dt,freq,decoded,irc, & subroutine q65_decoded (this,nutc,sync,nsnr,dt,freq,decoded,idec,ntrperiod)
qual,ntrperiod,fmid,w50)
use q65_decode use q65_decode
implicit none implicit none
@ -790,22 +789,17 @@ contains
real, intent(in) :: dt real, intent(in) :: dt
real, intent(in) :: freq real, intent(in) :: freq
character(len=37), intent(in) :: decoded character(len=37), intent(in) :: decoded
integer, intent(in) :: irc integer, intent(in) :: idec
real, intent(in) :: qual
integer, intent(in) :: ntrperiod integer, intent(in) :: ntrperiod
real, intent(in) :: fmid
real, intent(in) :: w50
integer navg
navg=irc/100
if(ntrperiod.lt.60) then if(ntrperiod.lt.60) then
write(*,1001) nutc,nsnr,dt,nint(freq),decoded,mod(irc,100),navg write(*,1001) nutc,nsnr,dt,nint(freq),decoded,idec
1001 format(i6.6,i4,f5.1,i5,' + ',1x,a37,1x,i2,i4) 1001 format(i6.6,i4,f5.1,i5,' + ',1x,a37,1x,i2)
write(13,1002) nutc,nint(sync),nsnr,dt,freq,0,decoded write(13,1002) nutc,nint(sync),nsnr,dt,freq,0,decoded
1002 format(i6.6,i4,i5,f6.1,f8.0,i4,3x,a37,' Q65') 1002 format(i6.6,i4,i5,f6.1,f8.0,i4,3x,a37,' Q65')
else else
write(*,1003) nutc,nsnr,dt,nint(freq),decoded,mod(irc,100),navg write(*,1003) nutc,nsnr,dt,nint(freq),decoded,idec
1003 format(i4.4,i4,f5.1,i5,' + ',1x,a37,1x,i2,i4) 1003 format(i4.4,i4,f5.1,i5,' + ',1x,a37,1x,i2)
write(13,1004) nutc,nint(sync),nsnr,dt,freq,0,decoded write(13,1004) nutc,nint(sync),nsnr,dt,freq,0,decoded
1004 format(i4.4,i4,i5,f6.1,f8.0,i4,3x,a37,' Q65') 1004 format(i4.4,i4,i5,f6.1,f8.0,i4,3x,a37,' Q65')

61
lib/q65_decode.f90
View File

@ -8,7 +8,7 @@ module q65_decode
abstract interface abstract interface
subroutine q65_decode_callback (this,nutc,sync,nsnr,dt,freq, & subroutine q65_decode_callback (this,nutc,sync,nsnr,dt,freq, &
decoded,nap,qual,ntrperiod,fmid,w50) decoded,nap,ntrperiod)
import q65_decoder import q65_decoder
implicit none implicit none
class(q65_decoder), intent(inout) :: this class(q65_decoder), intent(inout) :: this
@ -19,10 +19,7 @@ module q65_decode
real, intent(in) :: freq real, intent(in) :: freq
character(len=37), intent(in) :: decoded character(len=37), intent(in) :: decoded
integer, intent(in) :: nap integer, intent(in) :: nap
real, intent(in) :: qual
integer, intent(in) :: ntrperiod integer, intent(in) :: ntrperiod
real, intent(in) :: fmid
real, intent(in) :: w50
end subroutine q65_decode_callback end subroutine q65_decode_callback
end interface end interface
@ -64,6 +61,8 @@ contains
complex, allocatable :: c00(:) !Analytic signal, 6000 Sa/s complex, allocatable :: c00(:) !Analytic signal, 6000 Sa/s
complex, allocatable :: c0(:) !Analytic signal, 6000 Sa/s complex, allocatable :: c0(:) !Analytic signal, 6000 Sa/s
id1=0
id2=0
mode65=2**nsubmode mode65=2**nsubmode
npts=ntrperiod*12000 npts=ntrperiod*12000
nfft1=ntrperiod*12000 nfft1=ntrperiod*12000
@ -90,33 +89,27 @@ contains
this%callback => callback this%callback => callback
if(nutc.eq.-999) print*,lapdx,nfa,nfb,nfqso !Silence warning if(nutc.eq.-999) print*,lapdx,nfa,nfb,nfqso !Silence warning
nFadingModel=1 nFadingModel=1
call q65_set_list(mycall,hiscall,hisgrid,codewords,ncw)
dgen=0
call q65_enc(dgen,codewords) !Initialize Q65
call timer('sync_q65',0) call timer('sync_q65',0)
call sync_q65(iwave,ntrperiod*12000,mode65,nQSOprogress,nsps,nfqso, & call q65_sync(iwave,ntrperiod*12000,mode65,codewords,ncw,nsps, &
ntol,xdt,f0,snr1,width) nfqso,ntol,xdt,f0,snr1,dat4,snr2,id1)
call timer('sync_q65',1) call timer('sync_q65',1)
if(id1.eq.1) then
irc=-9 xdt1=xdt
f1=f0
go to 100
endif
if(snr1.lt.2.8) go to 100 if(snr1.lt.2.8) go to 100
jpk0=(xdt+1.0)*6000 !### Is this OK? jpk0=(xdt+1.0)*6000 !### Is this OK?
if(ntrperiod.le.30) jpk0=(xdt+0.5)*6000 !### if(ntrperiod.le.30) jpk0=(xdt+0.5)*6000 !###
if(jpk0.lt.0) jpk0=0 if(jpk0.lt.0) jpk0=0
fac=1.0/32767.0 fac=1.0/32767.0
dd=fac*iwave(1:npts) dd=fac*iwave(1:npts)
!###
! Optionslly write noise level to LU 56
! sq=dot_product(dd,dd)/npts
! m=nutc
! if(ntrperiod.ge.60) m=100*m
! ihr=m/10000
! imin=mod(m/100,100)
! isec=mod(m,100)
! hours=ihr + imin/60.0 + isec/3600.0
! write(56,3056) m,hours,db(sq)+90.3
!3056 format(i6.6,f10.6,f10.3)
!###
nmode=65 nmode=65
call ana64(dd,npts,c00) call ana64(dd,npts,c00)
call ft8apset(mycall,hiscall,ncontest,apsym0,aph10) call ft8apset(mycall,hiscall,ncontest,apsym0,aph10)
where(apsym0.eq.-1) apsym0=0 where(apsym0.eq.-1) apsym0=0
@ -146,37 +139,29 @@ contains
endif endif
endif endif
call timer('q65loops',0) call timer('q65loops',0)
! call q65_loops(c00,nutc,npts/2,nsps/2,nmode,mode65,nsubmode, & call q65_loops(c00,npts/2,nsps/2,nmode,mode65,nsubmode, &
! nFadingModel,ndepth,jpk0,xdt,f0,width,iaptype,apmask,apsymbols, & nFadingModel,ndepth,jpk0,xdt,f0,iaptype,apmask,apsymbols, &
! snr1,xdt1,f1,snr2,irc,dat4) xdt1,f1,snr2,dat4,id2)
! baud rate required to compute B90TS later
call q65_loops(c00,nutc,npts/2,nsps/2,nmode,mode65,nsubmode, &
nFadingModel,ndepth,jpk0,xdt,f0,width,iaptype,apmask,apsymbols, &
codewords,snr1,xdt1,f1,snr2,irc,dat4)
call timer('q65loops',1) call timer('q65loops',1)
snr2=snr2 + db(6912.0/nsps) snr2=snr2 + db(6912.0/nsps)
if(irc.ge.0) exit if(id2.gt.0) exit
enddo enddo
100 decoded=' ' 100 decoded=' '
! if(irc.lt.0 .and.iaptype.eq.4) print*,'AAA',irc,iaptype if(id1.gt.0 .or. id2.gt.0) then
if(irc.ge.0) then idec=id1+id2
!###
navg=irc/100
! irc=100*navg + ipass
irc=100*navg + iaptype
!###
write(c77,1000) dat4(1:12),dat4(13)/2 write(c77,1000) dat4(1:12),dat4(13)/2
1000 format(12b6.6,b5.5) 1000 format(12b6.6,b5.5)
call unpack77(c77,0,decoded,unpk77_success) !Unpack to get msgsent call unpack77(c77,0,decoded,unpk77_success) !Unpack to get msgsent
nsnr=nint(snr2) nsnr=nint(snr2)
call this%callback(nutc,sync,nsnr,xdt1,f1,decoded, & call this%callback(nutc,sync,nsnr,xdt1,f1,decoded, &
irc,qual,ntrperiod,fmid,w50) idec,ntrperiod)
else else
! Report sync, even if no decode. ! Report sync, even if no decode.
nsnr=db(snr1) - 35.0 nsnr=db(snr1) - 35.0
idec=-1
call this%callback(nutc,sync,nsnr,xdt1,f1,decoded, & call this%callback(nutc,sync,nsnr,xdt1,f1,decoded, &
irc,qual,ntrperiod,fmid,w50) idec,ntrperiod)
endif endif
return return

231
lib/sync_q65.f90 → lib/q65_sync.f90
View File

@ -1,11 +1,11 @@
subroutine sync_q65(iwave,nmax,mode65,nQSOprogress,nsps,nfqso,ntol, & subroutine q65_sync(iwave,nmax,mode_q65,codewords,ncw,nsps,nfqso,ntol, &
xdt,f0,snr1,width) xdt,f0,snr1,dat4,snr2,id1)
! Detect and align with the Q65 sync vector, returning time and frequency ! Detect and align with the Q65 sync vector, returning time and frequency
! offsets and SNR estimate. ! offsets and SNR estimate.
! Input: iwave(0:nmax-1) Raw data ! Input: iwave(0:nmax-1) Raw data
! mode65 Tone spacing 1 2 4 8 16 (A-E) ! mode_q65 Tone spacing 1 2 4 8 16 (A-E)
! nsps Samples per symbol at 12000 Sa/s ! nsps Samples per symbol at 12000 Sa/s
! nfqso Target frequency (Hz) ! nfqso Target frequency (Hz)
! ntol Search range around nfqso (Hz) ! ntol Search range around nfqso (Hz)
@ -15,23 +15,28 @@ subroutine sync_q65(iwave,nmax,mode65,nQSOprogress,nsps,nfqso,ntol, &
parameter (NSTEP=8) !Step size nsps/NSTEP parameter (NSTEP=8) !Step size nsps/NSTEP
parameter (LN=2176*63) !LN=LL*NN; LL=64*(mode_q65+2), NN=63 parameter (LN=2176*63) !LN=LL*NN; LL=64*(mode_q65+2), NN=63
character*37 msg,msgsent
integer*2 iwave(0:nmax-1) !Raw data integer*2 iwave(0:nmax-1) !Raw data
integer isync(22) !Indices of sync symbols integer isync(22) !Indices of sync symbols
integer itone(85) integer itone(85)
real, allocatable :: s1(:,:) !Symbol spectra, quarter-symbol steps integer codewords(63,64)
integer dat4(13)
integer ijpk(2)
real, allocatable :: s1(:,:) !Symbol spectra, 1/8-symbol steps
real, allocatable :: s3(:,:) !Data-symbol energies s3(LL,63)
real, allocatable :: ccf(:,:) !CCF(freq,lag) real, allocatable :: ccf(:,:) !CCF(freq,lag)
real, allocatable :: ccf1(:) !CCF(freq) at best lag real, allocatable :: ccf1(:) !CCF(freq) at best lag
real s3prob(0:63,63) !Symbol-value probabilities
real sync(85) !sync vector real sync(85) !sync vector
real s3(LN) !Symbol spectra
real s3prob(LN) !Symbol-value probabilities
complex, allocatable :: c0(:) !Complex spectrum of symbol complex, allocatable :: c0(:) !Complex spectrum of symbol
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 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/ data sync(1)/99.0/
save sync save sync
nfft=2*nsps id1=0
df=12000.0/nfft !Freq resolution = 0.5*baud dat4=0
LL=64*(2+mode_q65)
nfft=nsps
df=12000.0/nfft !Freq resolution = baud
istep=nsps/NSTEP 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 txt=85.0*nsps/12000.0
@ -40,6 +45,7 @@ subroutine sync_q65(iwave,nmax,mode65,nQSOprogress,nsps,nfqso,ntol, &
ia=ntol/df ia=ntol/df
allocate(s1(iz,jz)) allocate(s1(iz,jz))
allocate(s3(-64:LL-65,63))
allocate(c0(0:nfft-1)) allocate(c0(0:nfft-1))
allocate(ccf(-ia:ia,-53:214)) allocate(ccf(-ia:ia,-53:214))
allocate(ccf1(-ia:ia)) allocate(ccf1(-ia:ia))
@ -68,12 +74,13 @@ subroutine sync_q65(iwave,nmax,mode65,nQSOprogress,nsps,nfqso,ntol, &
s1(i,j)=real(c0(i))**2 + aimag(c0(i))**2 s1(i,j)=real(c0(i))**2 + aimag(c0(i))**2
enddo enddo
! For large Doppler spreads, should we smooth the spectra here? ! For large Doppler spreads, should we smooth the spectra here?
call smo121(s1(1:iz,j),iz) ! call smo121(s1(1:iz,j),iz)
enddo enddo
i0=nint(nfqso/df) !Target QSO frequency i0=nint(nfqso/df) !Target QSO frequency
call pctile(s1(i0-64:i0+192,1:jz),129*jz,40,base) call pctile(s1(i0-64:i0+192,1:jz),129*jz,40,base)
s1=s1/base - 1.0 ! s1=s1/base - 1.0
s1=s1/base
! Apply fast AGC ! Apply fast AGC
s1max=20.0 !Empirical choice s1max=20.0 !Empirical choice
@ -91,8 +98,90 @@ subroutine sync_q65(iwave,nmax,mode65,nQSOprogress,nsps,nfqso,ntol, &
j0=1.0/dtstep !Nominal index for start of signal j0=1.0/dtstep !Nominal index for start of signal
lag2=4.0/dtstep + 0.9999 !Include EME delays lag2=4.0/dtstep + 0.9999 !Include EME delays
endif endif
ccf=0.
if(ncw.lt.1) go to 100
!######################################################################
! Try list decoding via "Deep Likelihood".
ipk=0
jpk=0
ccf_best=0.
do imsg=1,ncw
i=1
k=0
do j=1,85
if(j.eq.isync(i)) then
i=i+1
itone(j)=-1
else
k=k+1
itone(j)=codewords(k,imsg)
endif
enddo
! Compute 2D ccf using all 85 symbols in the list message
ccf=0.
do lag=lag1,lag2
do k=1,85
j=j0 + NSTEP*(k-1) + 1 + lag
if(j.ge.1 .and. j.le.jz) then
do i=-ia,ia
ii=i0+mode_q65*itone(k)+i
ccf(i,lag)=ccf(i,lag) + s1(ii,j)
enddo
endif
enddo
enddo
ccfmax=maxval(ccf)
if(ccfmax.gt.ccf_best) then
ccf_best=ccfmax
ijpk=maxloc(ccf)
ipk=ijpk(1)-ia-1
jpk=ijpk(2)-53-1
f0=nfqso + ipk*df
xdt=jpk*dtstep
endif
enddo ! imsg
ia=i0+ipk-63
ib=ia+LL-1
j=j0+jpk-5
n=0
do k=1,85
j=j+8
if(sync(k).gt.0.0) then
cycle
endif
n=n+1
if(j.ge.1 .and. j.le.jz) s3(-64:LL-65,n)=s1(ia:ib,j)
enddo
nsubmode=0
if(mode_q65.eq.2) nsubmode=1
if(mode_q65.eq.4) nsubmode=2
if(mode_q65.eq.8) nsubmode=3
if(mode_q65.eq.16) nsubmode=4
nFadingModel=1
baud=12000.0/nsps
do ibw=0,10
b90=1.72**ibw
call q65_intrinsics_ff(s3,nsubmode,b90/baud,nFadingModel,s3prob)
call q65_dec_fullaplist(s3,s3prob,codewords,ncw,esnodb,dat4,plog,irc)
if(irc.ge.0) then
snr2=esnodb - db(2500.0/baud)
id1=1
! write(55,3055) nutc,xdt,f0,snr2,plog,irc
!3055 format(i4.4,4f9.2,i5)
go to 900
endif
enddo
!######################################################################
! Establish xdt, f0, and snr1 using sync symbols (and perhaps some AP symbols)
100 ccf=0.
irc=-2
dat4=0
ia=ntol/df
do lag=lag1,lag2 do lag=lag1,lag2
do k=1,85 do k=1,85
n=NSTEP*(k-1) + 1 n=NSTEP*(k-1) + 1
@ -102,23 +191,11 @@ subroutine sync_q65(iwave,nmax,mode65,nQSOprogress,nsps,nfqso,ntol, &
endif endif
enddo enddo
enddo enddo
ijpk=maxloc(ccf)
ic=ntol/df ipk=ijpk(1)-ia-1
ccfmax=0. jpk=ijpk(2)-53-1
ipk=0
jpk=0
do i=-ic,ic
do j=lag1,lag2
if(ccf(i,j).gt.ccfmax) then
ipk=i
jpk=j
ccfmax=ccf(i,j)
endif
enddo
enddo
f0=nfqso + ipk*df f0=nfqso + ipk*df
xdt=jpk*dtstep xdt=jpk*dtstep
sq=0. sq=0.
nsq=0 nsq=0
do j=lag1,lag2 do j=lag1,lag2
@ -130,104 +207,6 @@ subroutine sync_q65(iwave,nmax,mode65,nQSOprogress,nsps,nfqso,ntol, &
rms=sqrt(sq/nsq) rms=sqrt(sq/nsq)
smax=ccf(ipk,jpk) smax=ccf(ipk,jpk)
snr1=smax/rms snr1=smax/rms
! do j=lag1,lag2
! write(55,3055) j,j*dtstep,ccf(ipk,j)/rms
!3055 format(i5,f8.3,f10.3)
! enddo
! do i=-ia,ia
! write(56,3056) i*df,ccf(i,jpk)/rms
!3056 format(2f10.3)
! enddo
! flush(56)
ccf1=ccf(-ia:ia,jpk)
acf0=dot_product(ccf1,ccf1)
do i=1,ia
acf=dot_product(ccf1,cshift(ccf1,i))
if(acf.le.0.5*acf0) exit
enddo
width=i*1.414*df
!### Experimental:
if(nQSOprogress.lt.1) go to 900
! "Deep Likelihood" decode attempt
snr1a_best=0.
do imsg=1,4
ccf=0.
msg='K1ABC W9XYZ RRR'
if(imsg.eq.2) msg='K1ABC W9XYZ RR73'
if(imsg.eq.3) msg='K1ABC W9XYZ 73'
if(imsg.eq.4) msg='CQ K9AN EN50'
call genq65(msg,0,msgsent,itone,i3,n3)
do lag=lag1,lag2
do k=1,85
j=j0 + NSTEP*(k-1) + 1 + lag
if(j.ge.1 .and. j.le.jz) then
do i=-ia,ia
ii=i0+2*itone(k)+i
ccf(i,lag)=ccf(i,lag) + s1(ii,j)
enddo
endif
enddo
enddo
ic=ntol/df
ccfmax=0.
ipk=0
jpk=0
do i=-ic,ic
do j=lag1,lag2
if(ccf(i,j).gt.ccfmax) then
ipk=i
jpk=j
ccfmax=ccf(i,j)
endif
enddo
enddo
f0a=nfqso + ipk*df
xdta=jpk*dtstep
sq=0.
nsq=0
do j=lag1,lag2
if(abs(j-jpk).gt.6) then
sq=sq + ccf(ipk,j)**2
nsq=nsq+1
endif
enddo
rms=sqrt(sq/nsq)
smax=ccf(ipk,jpk)
snr1a=smax/rms
if(snr1a.gt.snr1a_best) then
snr1a_best=snr1a
imsg_best=imsg
xdta_best=xdta
f0a_best=f0a
endif
! write(57,3001) imsg,xdt,xdta,f0,f0a,snr1,snr1a
!3001 format(i1,6f8.2)
! do j=lag1,lag2
! write(55,3055) j,j*dtstep,ccf(ipk,j)/rms
!3055 format(i5,f8.3,f10.3)
! enddo
! do i=-ia,ia
! write(56,3056) i*df,ccf(i,jpk)/rms
!3056 format(2f10.3)
! enddo
enddo
if(snr1a_best.gt.2.0) then
xdt=xdta_best
f0=f0a_best
snr1=1.4*snr1a_best
endif
! write(58,3006) xdta_best,f0a_best,snr1a_best,imsg_best
!3006 format(3f8.2,i3)
900 return 900 return
end subroutine sync_q65 end subroutine q65_sync

3
lib/qra/q65/q65.c
View File

@ -708,10 +708,11 @@ int q65_decode_fullaplist(q65_codec_ds *codec,
maxllh = llh; maxllh = llh;
maxcw = k; maxcw = k;
} }
// printf("BBB %d %f\n",k,llh);
// point to next codeword // point to next codeword
pCw+=nN; pCw+=nN;
} }
q65_llh=maxllh;
if (maxcw<0) // no llh larger than threshold found if (maxcw<0) // no llh larger than threshold found
return Q65_DECODE_FAILED; return Q65_DECODE_FAILED;

2
lib/qra/q65/q65.h
View File

@ -39,7 +39,7 @@
// maximum number of weights for the fast-fading metric evaluation // maximum number of weights for the fast-fading metric evaluation
#define Q65_FASTFADING_MAXWEIGTHS 65 #define Q65_FASTFADING_MAXWEIGTHS 65
float q65_llh;
typedef struct { typedef struct {
const qracode *pQraCode; // qra code to be used by the codec const qracode *pQraCode; // qra code to be used by the codec
float decoderEsNoMetric; // value for which we optimize the decoder metric float decoderEsNoMetric; // value for which we optimize the decoder metric

173
lib/qra/q65/q65_loops.f90
View File

@ -1,33 +1,27 @@
subroutine q65_loops(c00,nutc,npts2,nsps,mode,mode_q65,nsubmode,nFadingModel, & subroutine q65_loops(c00,npts2,nsps,mode,mode_q65,nsubmode,nFadingModel, &
ndepth,jpk0,xdt0,f0,width,iaptype,APmask,APsymbols,codewords,snr1, & ndepth,jpk0,xdt0,f0,iaptype,APmask,APsymbols,xdt1,f1,snr2,dat4,id2)
xdt1,f1,snr2,irc,dat4)
use packjt77 use packjt77
use timer_module, only: timer use timer_module, only: timer
parameter (LN=2176*63) !LN=LL*NN; LL=64*(mode_q65+2), NN=63 parameter (NN=63)
character*37 decoded parameter (LN=1152*63) !LN=LL*NN; LL=64*(mode_q65+2), NN=63
character*77 c77
complex c00(0:npts2-1) !Analytic representation of dd(), 6000 Hz complex c00(0:npts2-1) !Analytic representation of dd(), 6000 Hz
complex ,allocatable :: c0(:) !Ditto, with freq shift complex ,allocatable :: c0(:) !Ditto, with freq shift
real a(3) !twkfreq params f,f1,f2 real a(3) !twkfreq params f,f1,f2
real s3(LN) !Symbol spectra real s3(LN) !Symbol spectra
real s3avg(LN) !Averaged symbol spectra real s3prob(64*NN) !Symbol-value probabilities
real s3prob(LN) !Symbol-value probabilities
logical unpk77_success
integer APmask(13) integer APmask(13)
integer APsymbols(13) integer APsymbols(13)
integer codewords(63,64) integer cw4(63)
! integer cw4(63)
integer dat4(13) !Decoded message (as 13 six-bit integers) integer dat4(13) !Decoded message (as 13 six-bit integers)
integer nap(0:11) !AP return codes integer nap(0:11) !AP return codes
data nap/0,2,3,2,3,4,2,3,6,4,6,6/,nsave/0/ data nap/0,2,3,2,3,4,2,3,6,4,6,6/,nsave/0/
! data cw4/0, 0, 0, 0, 8, 4,60,35,17,48,33,25,34,43,43,43,35,15,46,30, & data cw4/0, 0, 0, 0, 8, 4,60,35,17,48,33,25,34,43,43,43,35,15,46,30, &
! 54,24,26,26,57,57,42, 3,23,11,49,49,16, 2, 6, 6,55,21,39,51, & 54,24,26,26,57,57,42, 3,23,11,49,49,16, 2, 6, 6,55,21,39,51, &
! 51,51,42,42,50,25,31,35,57,30, 1,54,54,10,10,22,44,58,57,40, & 51,51,42,42,50,25,31,35,57,30, 1,54,54,10,10,22,44,58,57,40, &
! 21,21,19/ 21,21,19/
save nsave,s3avg
id2=-1
ircbest=9999 ircbest=9999
allocate(c0(0:npts2-1)) allocate(c0(0:npts2-1))
irc=-99 irc=-99
@ -47,112 +41,69 @@ subroutine q65_loops(c00,nutc,npts2,nsps,mode,mode_q65,nsubmode,nFadingModel, &
ibwmax=5 ibwmax=5
endif endif
LL=64*(mode_q65+2) LL=64*(mode_q65+2)
NN=63
napmin=99 napmin=99
baud=6000.0/nsps baud=6000.0/nsps
xdt1=xdt0 xdt1=xdt0
f1=f0 f1=f0
maxavg=0 do idf=1,idfmax
if(iand(ndepth,16).ne.0) maxavg=1 ndf=idf/2
do iavg=0,maxavg if(mod(idf,2).eq.0) ndf=-ndf
if(iavg.eq.1) then a=0.
idfmax=1 a(1)=-(f0+0.5*baud*ndf)
idtmax=1 call twkfreq(c00,c0,npts2,6000.0,a)
endif do idt=1,idtmax
do idf=1,idfmax ndt=idt/2
ndf=idf/2 if(iaptype.eq.0) then
if(mod(idf,2).eq.0) ndf=-ndf if(mod(idt,2).eq.0) ndt=-ndt
a=0. jpk=jpk0 + nsps*ndt/16 !tsym/16
a(1)=-(f0+0.5*baud*ndf) if(jpk.lt.0) jpk=0
call twkfreq(c00,c0,npts2,6000.0,a) call timer('spec64 ',0)
do idt=1,idtmax call spec64(c0,nsps,mode,mode_q65,jpk,s3,LL,NN)
ndt=idt/2 call timer('spec64 ',1)
if(iaptype.eq.0 .and. iavg.eq.0) then call pctile(s3,LL*NN,40,base)
if(mod(idt,2).eq.0) ndt=-ndt s3=s3/base
jpk=jpk0 + nsps*ndt/16 !tsym/16 where(s3(1:LL*NN)>s3lim) s3(1:LL*NN)=s3lim
if(jpk.lt.0) jpk=0 endif
call timer('spec64 ',0) do ibw=ibwmin,ibwmax
call spec64(c0,nsps,mode,mode_q65,jpk,s3,LL,NN) nbw=ibw
call timer('spec64 ',1) ndist=ndf**2 + ndt**2 + ((nbw-2))**2
call pctile(s3,LL*NN,40,base) if(ndist.gt.maxdist) cycle
s3=s3/base ! b90=1.728**ibw
where(s3(1:LL*NN)>s3lim) s3(1:LL*NN)=s3lim b90=3.0**nbw
endif if(b90.gt.230.0) cycle
if(iavg.eq.1) then call timer('q65_intr',0)
s3(1:LL*NN)=s3avg(1:LL*NN) b90ts = b90/baud
endif call q65_intrinsics_ff(s3,nsubmode,b90ts,nFadingModel,s3prob)
do ibw=ibwmin,ibwmax call timer('q65_intr',1)
nbw=ibw call timer('q65_dec ',0)
ndist=ndf**2 + ndt**2 + ((nbw-2))**2 call q65_dec(s3,s3prob,APmask,APsymbols,esnodb,dat4,irc)
if(ndist.gt.maxdist) cycle call timer('q65_dec ',1)
! b90=1.728**ibw if(irc.ge.0) id2=iaptype+2
b90=3.0**nbw if(irc.ge.0) go to 100
if(b90.gt.230.0) cycle
! if(b90.lt.0.15*width) exit
call timer('q65_intr',0)
b90ts = b90/baud
call q65_intrinsics_ff(s3,nsubmode,b90ts,nFadingModel,s3prob)
call timer('q65_intr',1)
if(iaptype.eq.4) then
! codewords(1:63,4)=cw4
call timer('q65_apli',0)
call q65_dec_fullaplist(s3,s3prob,codewords,3,esnodb,dat4,irc)
call timer('q65_apli',1)
else
call timer('q65_dec ',0)
call q65_dec(s3,s3prob,APmask,APsymbols,esnodb,dat4,irc)
call timer('q65_dec ',1)
endif
if(irc.ge.0) go to 100
! irc > 0 ==> number of iterations required to decode ! irc > 0 ==> number of iterations required to decode
! -1 = invalid params ! -1 = invalid params
! -2 = decode failed ! -2 = decode failed
! -3 = CRC mismatch ! -3 = CRC mismatch
enddo ! ibw (b90 loop) enddo ! ibw (b90 loop)
enddo ! idt (DT loop) enddo ! idt (DT loop)
enddo ! idf (f0 loop) enddo ! idf (f0 loop)
if(iaptype.eq.0 .and. iavg.eq.0) then if(iaptype.eq.0) then
a=0. a=0.
a(1)=-f0 a(1)=-f0
call twkfreq(c00,c0,npts2,6000.0,a) call twkfreq(c00,c0,npts2,6000.0,a)
jpk=3000 !### Are these definitions OK? jpk=3000 !### Are these definitions OK?
if(nsps.ge.3600) jpk=6000 !### TR >= 60 s if(nsps.ge.3600) jpk=6000 !### TR >= 60 s
call spec64(c0,nsps,mode,mode_q65,jpk,s3,LL,NN) call spec64(c0,nsps,mode,mode_q65,jpk,s3,LL,NN)
call pctile(s3,LL*NN,40,base) call pctile(s3,LL*NN,40,base)
s3=s3/base s3=s3/base
where(s3(1:LL*NN)>s3lim) s3(1:LL*NN)=s3lim where(s3(1:LL*NN)>s3lim) s3(1:LL*NN)=s3lim
s3avg(1:LL*NN)=s3avg(1:LL*NN)+s3(1:LL*NN) endif
nsave=nsave+1
endif
if(iavg.eq.0 .and. nsave.lt.2) exit
enddo ! iavg
100 if(irc.ge.0) then 100 if(irc.ge.0) then
navg=nsave
snr2=esnodb - db(2500.0/baud) snr2=esnodb - db(2500.0/baud)
if(iavg.eq.0) navg=0
xdt1=xdt0 + nsps*ndt/(16.0*6000.0) xdt1=xdt0 + nsps*ndt/(16.0*6000.0)
f1=f0 + 0.5*baud*ndf f1=f0 + 0.5*baud*ndf
!### For tests only:
open(53,file='fort.53',status='unknown',position='append')
write(c77,1100) dat4(1:12),dat4(13)/2
1100 format(12b6.6,b5.5)
call unpack77(c77,0,decoded,unpk77_success) !Unpack to get msgsent
m=nutc
if(nsps.ge.3600) m=100*m
ihr=m/10000
imin=mod(m/100,100)
isec=mod(m,100)
hours=ihr + imin/60.0 + isec/3600.0
write(53,3053) m,hours,ndf,ndt,nbw,ndist,irc,iaptype,navg,snr1, &
xdt1,f1,snr2,trim(decoded)
3053 format(i6.6,f8.4,4i3,i4,2i3,f6.1,f6.2,f7.1,f6.1,1x,a)
close(53)
!###
nsave=0
s3avg=0.
irc=irc + 100*navg
endif endif
return return

42
lib/qra/q65/q65_set_list.f90 Normal file
View File

@ -0,0 +1,42 @@
subroutine q65_set_list(mycall,hiscall,hisgrid,codewords,ncw)
character*12 mycall,hiscall
character*6 hisgrid
character*37 msg0,msg,msgsent
integer codewords(63,64)
integer itone(85)
integer isync(22)
data isync/1,9,12,13,15,22,23,26,27,33,35,38,46,50,55,60,62,66,69,74,76,85/
msg0=trim(mycall)//' '//trim(hiscall)
j0=len(trim(msg0))+2
isnr0=-35
do i=1,57
msg=msg0
if(i.eq.2) msg(j0:j0+2)='RRR'
if(i.eq.3) msg(j0:j0+3)='RR73'
if(i.eq.4) msg(j0:j0+1)='73'
if(i.ge.5 .and. i.le.56) then
isnr=isnr0 + (i-5)/2
if(iand(i,1).eq.1) write(msg(j0:j0+2),'(i3.2)') isnr
if(iand(i,1).eq.0) write(msg(j0:j0+3),'("R",i3.2)') isnr
endif
if(i.eq.57) msg='CQ '//trim(hiscall)//' '//hisgrid(1:4)
call genq65(msg,0,msgsent,itone,i3,n3)
i0=1
j=0
do k=1,85
if(k.eq.isync(i0)) then
i0=i0+1
cycle
endif
j=j+1
codewords(j,i)=itone(k) - 1
enddo
ncw=57
! write(*,3001) i,isnr,codewords(1:13,i),trim(msg)
!3001 format(i2,2x,i3.2,2x,13i3,2x,a)
enddo
return
end subroutine q65_set_list

3
lib/qra/q65/q65_subs.c
View File

@ -111,7 +111,7 @@ void q65_dec_(float s3[], float s3prob[], int APmask[], int APsymbols[],
} }
void q65_dec_fullaplist_(float s3[], float s3prob[], int codewords[], void q65_dec_fullaplist_(float s3[], float s3prob[], int codewords[],
int* ncw, float* esnodb0, int xdec[], int* rc0) int* ncw, float* esnodb0, int xdec[], float* plog, int* rc0)
{ {
/* Input: s3[LL,NN] Symbol spectra /* Input: s3[LL,NN] Symbol spectra
* s3prob[LL,NN] Symbol-value intrinsic probabilities * s3prob[LL,NN] Symbol-value intrinsic probabilities
@ -128,6 +128,7 @@ void q65_dec_fullaplist_(float s3[], float s3prob[], int codewords[],
float esnodb; float esnodb;
rc = q65_decode_fullaplist(&codec,ydec,xdec,s3prob,codewords,*ncw); rc = q65_decode_fullaplist(&codec,ydec,xdec,s3prob,codewords,*ncw);
*plog=q65_llh;
*rc0=rc; *rc0=rc;
// rc = -1: Invalid params // rc = -1: Invalid params

30
lib/qra/q65/q65sim.f90
View File

@ -193,21 +193,21 @@ program q65sim
write(10) h,iwave(1:npts) !Save the .wav file write(10) h,iwave(1:npts) !Save the .wav file
close(10) close(10)
if(lsync) then ! if(lsync) then
cd=' ' ! cd=' '
if(ifile.eq.nfiles) cd='d' ! if(ifile.eq.nfiles) cd='d'
nfqso=nint(f0) ! nfqso=nint(f0)
ntol=100 ! ntol=100
call sync_q65(iwave,npts,mode65,nsps,nfqso,ntol,xdt2,f02,snr2) ! call q65_sync(iwave,npts,mode65,nsps,nfqso,ntol,xdt2,f02,snr2)
terr=1.01/(8.0*baud) ! terr=1.01/(8.0*baud)
ferr=1.01*mode65*baud ! ferr=1.01*mode65*baud
if(abs(xdt2-xdt).lt.terr .and. abs(f02-f0).lt.ferr) nsync=nsync+1 ! if(abs(xdt2-xdt).lt.terr .and. abs(f02-f0).lt.ferr) nsync=nsync+1
open(40,file='sync65.out',status='unknown',position='append') ! open(40,file='sync65.out',status='unknown',position='append')
write(40,1030) ifile,65,csubmode,snrdb,fspread,xdt2-xdt,f02-f0, & ! write(40,1030) ifile,65,csubmode,snrdb,fspread,xdt2-xdt,f02-f0, &
snr2,nsync,cd ! snr2,nsync,cd
1030 format(i4,i3,1x,a1,2f7.1,f7.2,2f8.1,i5,1x,a1) !1030 format(i4,i3,1x,a1,2f7.1,f7.2,2f8.1,i5,1x,a1)
close(40) ! close(40)
endif ! endif
enddo enddo
if(lsync) write(*,1040) snrdb,nfiles,nsync if(lsync) write(*,1040) snrdb,nfiles,nsync
1040 format('SNR:',f6.1,' nfiles:',i5,' nsynced:',i5) 1040 format('SNR:',f6.1,' nfiles:',i5,' nsynced:',i5)