w4kek/WSJT-X
1
0
Fork 0
mirror of https://github.com/saitohirga/WSJT-X.git synced 2024-11-15 08:31:57 -05:00
Code Issues Packages Projects Releases Wiki Activity

Routines for 72ms jtmsk messages. Nothing works yet.

Browse Source 
git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@6704 ab8295b8-cf94-4d9e-aec4-7959e3be5d79
This commit is contained in:
Steven Franke 2016-05-26 21:19:09 +00:00
parent 2bb851d164
commit 9773afcb0a
4 changed files with 662 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

71
lib/JTMSKsim144.f90 Normal file
View File

@ -0,0 +1,71 @@
program JTMSKsim
use wavhdr
parameter (NMAX=15*12000)
real pings(0:NMAX-1)
real waveform(0:864-1)
character arg*8,msg*22,msgsent*22,fname*40
character*3 rpt(0:7)
real wave(0:NMAX-1) !Simulated received waveform
real*8 dt,twopi,freq,phi,dphi0,dphi1,dphi
type(hdr) h !Header for .wav file
integer*2 iwave(0:NMAX-1)
integer itone(144) !Message bits
integer b11(11) !Barker-11 code
data b11/1,1,1,0,0,0,1,0,0,1,0/
data rpt /'26 ','27 ','28 ','R26','R27','R28','RRR','73 '/
nargs=iargc()
if(nargs.ne.5) then
print*,'Usage: JTMSKsim message freq width snr nfiles'
print*,' '
print*,'Examples: JTMSKsim "K1ABC W9XYZ EN37" 1500 0.12 2 1'
print*,' JTMSKsim "<K1ABC W9XYZ> R26" 1500 0.01 1 3'
go to 999
endif
call getarg(1,msg)
call getarg(2,arg)
read(arg,*) freq
call getarg(3,arg)
read(arg,*) width
call getarg(4,arg)
read(arg,*) snrdb
call getarg(5,arg)
read(arg,*) nfiles
sig=sqrt(2.0)*10.0**(0.05*snrdb)
twopi=8.d0*atan(1.d0)
h=default_header(12000,NMAX)
ichk=0
call genmsk(msg,ichk,msgsent,waveform,itype) !Check message type
if(itype.lt.1 .or. itype.gt.7) then
print*,'Illegal message'
go to 999
endif
! call makepings(pings,NMAX,width,sig)
pings=0.0
pings(12345:24000)=sig
do ifile=1,nfiles !Loop over requested number of files
write(fname,1002) ifile !Output filename
1002 format('000000_',i4.4)
open(10,file=fname(1:11)//'.wav',access='stream',status='unknown')
fac=sqrt(6000.0/2500.0)
j=-1
do i=0,NMAX-1
j=mod(j+1,864)
xx=gran()
wave(i)=pings(i)*waveform(j) + fac*xx
write(*,*) pings(i),fac,waveform(j),wave(j)
iwave(i)=30.0*wave(i)
enddo
write(10) h,iwave !Save the .wav file
close(10)
! call jtmsk_short(cwave,NMAX,msg)
enddo
999 end program JTMSKsim

196
lib/genmsk144.f90 Normal file
View File

@ -0,0 +1,196 @@
subroutine genmsk144(msg0,ichk,msgsent,waveform,itype)
!!!!!!!!!!!!!!!!!!!!! Experimental ldpc version
! b11 + 128bit code + 1 zero bit = 140 bits (70ms message duration)
! b13 + 128bit code + 2 slack + 1 zero bit = 144 bits (71ms message duration)
!
! Encode a JTMSK message
! Input:
! - msg0 requested message to be transmitted
! - ichk if ichk=1, return only msgsent
! if ichk.ge.10000, set imsg=ichk-10000 for short msg
! - msgsent message as it will be decoded
! - i4tone array of audio tone values, 0 or 1
! - itype message type
! 1 = standard message "Call_1 Call_2 Grid/Rpt"
! 2 = type 1 prefix
! 3 = type 1 suffix
! 4 = type 2 prefix
! 5 = type 2 suffix
! 6 = free text (up to 13 characters)
! 7 = short message "<Call_1 Call2> Rpt"
use iso_c_binding, only: c_loc,c_size_t
use packjt
use hashing
character*85 pchk_file,gen_file
character*22 msg0
character*22 message !Message to be generated
character*22 msgsent !Message as it will be received
integer*4 i4Msg6BitWords(13) !72-bit message as 6-bit words
integer*4 i4tone(144) !
integer*1, target:: i1Msg8BitBytes(10) !80 bits represented in 10 bytes
integer*1 codeword(128) !Encoded bits before re-ordering
integer*1 msgbits(80) !72-bit message + 8-bit hash
integer*1 reorderedcodeword(128) !Odd bits first, then even
integer*1 bitseq(144) !Tone #s, data and sync (values 0-1)
integer*1 i1hash(4)
integer*1 b7(7)
integer*1 b8(8)
integer*1 b11(11)
integer*1 b13(13)
real*8 pp(12)
real*8 xi(864),xq(864),pi,twopi,phi,dphi
real waveform(864)
data b7/1,1,1,0,0,1,0/
data b8/0,1,1,1,0,0,1,0/
data b11/1,1,1,0,0,0,1,0,0,1,0/ !Barker 11 code
data b13/1,1,1,1,1,0,0,1,1,0,1,0,1/ !Barker 13 code
equivalence (ihash,i1hash)
logical first
data first/.true./
save
if( first ) then
print*,"Initializing ldpc"
first=.false.
nsym=128
pchk_file="peg-128-80-reg3.pchk"
gen_file="peg-128-80-reg3.gen"
call init_ldpc(trim(pchk_file)//char(0),trim(gen_file)//char(0))
pi=4.*atan(1.0)
twopi=8.*atan(1.0)
do i=1,12
pp(i)=sin( (i-1)*pi/12 )
enddo
endif
if(msg0(1:1).eq.'@') then !Generate a fixed tone
read(msg0(2:5),*,end=1,err=1) nfreq !at specified frequency
go to 2
1 nfreq=1000
2 i4tone(1)=nfreq
else
message=msg0
do i=1,22
if(ichar(message(i:i)).eq.0) then
message(i:)=' '
exit
endif
enddo
do i=1,22 !Strip leading blanks
if(message(1:1).ne.' ') exit
message=message(i+1:)
enddo
if(message(1:1).eq.'<') then
call genmsk_short(message,msgsent,ichk,i4tone,itype)
if(itype.lt.0) go to 999
i4tone(36)=-35
go to 999
endif
call packmsg(message,i4Msg6BitWords,itype) !Pack into 12 6-bit bytes
call unpackmsg(i4Msg6BitWords,msgsent) !Unpack to get msgsent
if(ichk.eq.1) go to 999
i4=0
ik=0
im=0
do i=1,12
nn=i4Msg6BitWords(i)
do j=1, 6
ik=ik+1
i4=i4+i4+iand(1,ishft(nn,j-6))
i4=iand(i4,255)
if(ik.eq.8) then
im=im+1
! if(i4.gt.127) i4=i4-256
i1Msg8BitBytes(im)=i4
ik=0
endif
enddo
enddo
ihash=nhash(c_loc(i1Msg8BitBytes),int(9,c_size_t),146)
ihash=2*iand(ihash,32767) !Generate the 8-bit hash
i1Msg8BitBytes(10)=i1hash(1) !CRC to byte 10
mbit=0
do i=1, 10
i1=i1Msg8BitBytes(i)
do ibit=1,8
mbit=mbit+1
msgbits(mbit)=iand(1,ishft(i1,ibit-8))
enddo
enddo
call ldpc_encode(msgbits,codeword)
j=0
do i=1,nsym/2 !Reorder the encoded bits
j=j+1
reorderedcodeword(j)=codeword(2*i-1)
reorderedcodeword(j+nsym/2)=codeword(2*i)
enddo
!Create 144-bit channel vector:
!8-bit sync word + 48 bits + 8-bit sync word + 80 bits
bitseq=0
bitseq(1:8)=b8
bitseq(9:56)=reorderedcodeword(1:48)
bitseq(57:64)=b8
bitseq(65:144)=reorderedcodeword(49:128)
bitseq=2*bitseq-1
xq(1:6)=bitseq(1)*pp(7:12) !first bit is mapped to 1st half-symbol on q
do i=1,71
is=(i-1)*12+7
xq(is:is+11)=bitseq(2*i+1)*pp
enddo
xq(864-5:864)=bitseq(1)*pp(1:6) !first half of last zero bit on q channel
do i=1,72
is=(i-1)*12+1
xi(is:is+11)=bitseq(2*i)*pp
enddo
dphi=2*pi*1500/12000.0;
phi=0.0
do i=1,864
waveform(i)=imag(complex(xi(i),xq(i))*complex(cos(phi),sin(phi)))
phi=mod(phi+dphi,twopi)
enddo
! Map bits to I and Q (in anticipation of coherent demodulation),
! If bit indices run from 1 to 140.
! I bits are 1,3,5,...139
! Q bits are 2,4,6,...140
! Sample the offset qpsk waveform at 140 sampling instants to determine tones.
! Send low tone if I and Q have different sign, send high tone if I and Q have same sign.
! Modulator will take care of ensuring phase continuity.
! i4tone=0
! iq=140
! it=1
! do ii=1,139,2
! do inner=1,2
! if( bitseq(ii) .eq. bitseq(iq) ) then
! i4tone(it)=1
! endif
! print*,it,ii,iq,bitseq(ii),bitseq(iq),bitseq(it),i4tone(it)
! if( mod(it,2) .eq. 1 ) then
! iq=mod(iq+2,140)
! if( iq .eq. 0 ) iq=140
! endif
! it=it+1
! enddo
! enddo
endif
!i4tone(1:11)=b11
! print*,"transmitted tones"
! do i=1,140
! print*,i,i4tone(i)
! enddo
999 return
end subroutine genmsk

125
lib/jtmsk_decode144.f90 Normal file
View File

@ -0,0 +1,125 @@
subroutine jtmsk_decode(id2,narg,line)
! Decoder for JTMSK mode
parameter (NMAX=30*12000)
parameter (NFFTMAX=512*1024)
parameter (NSPM=840) !Samples per JTMSK long message
integer*2 id2(0:NMAX) !Raw i*2 data, up to T/R = 30 s
integer hist(0:32868)
real d(0:NMAX) !Raw r*4 data
real ty(NMAX/512) !Ping times
real yellow(NMAX/512)
complex c(NFFTMAX) !Complex (analytic) data
complex cdat(24000) !Short segments, up to 2 s
complex cdat2(24000)
integer narg(0:14) !Arguments passed from calling pgm
character*22 msg,msg0 !Decoded message
character*80 line(100) !Decodes passed back to caller
equivalence (hist,d)
! Parameters from GUI are in narg():
nutc=narg(0) !UTC
npts=min(narg(1),NMAX) !Number of samples in id2 (12000 Hz)
newdat=narg(3) !1==> new data, compute symbol spectra
minsync=narg(4) !Lower sync limit
npick=narg(5)
t0=0.001*narg(6)
t1=0.001*narg(7)
maxlines=narg(8) !Max # of decodes to return to caller
nmode=narg(9)
nrxfreq=narg(10) !Target Rx audio frequency (Hz)
ntol=narg(11) !Search range, +/- ntol (Hz)
nhashcalls=narg(12)
naggressive=narg(14)
nsnr0=-99
nline=0
line(1:100)(1:1)=char(0)
msg0=' '
msg=msg0
hist=0
do i=0,npts-1
n=abs(id2(i))
hist(n)=hist(n)+1
enddo
ns=0
do n=0,32768
ns=ns+hist(n)
if(ns.gt.npts/2) exit
enddo
fac=1.0/(1.5*n)
d(0:npts-1)=fac*id2(0:npts-1)
! rms=sqrt(dot_product(d(0:npts-1),d(0:npts-1))/npts)
!### Would it be better to set median rms to 1.0 ?
! d(0:npts-1)=d(0:npts-1)/rms !Normalize so that rms=1.0
call mskdt(d,npts,ty,yellow,nyel)
print*,nyel
do i=1,nyel
print*,i,ty(i),yellow(i)
enddo
nyel=min(nyel,5)
n=log(float(npts))/log(2.0) + 1.0
nfft=min(2**n,1024*1024)
call analytic(d,npts,nfft,c) !Convert to analytic signal and filter
nbefore=NSPM
nafter=4*NSPM
! Process ping list (sorted by S/N) from top down.
! do n=1,nyel
do n=1,1
ia=ty(n)*12000.0 - nbefore
if(ia.lt.1) ia=1
ib=ia + nafter
if(ib.gt.NFFTMAX) ib=NFFTMAX
iz=ib-ia+1
cdat2(1:iz)=c(ia:ib) !Select nlen complex samples
ja=ia/NSPM + 1
jb=ib/NSPM
t0=ia/12000.0
! do itry=1,21
do itry=1,1
idf1=(itry/2) * 50
if(mod(itry,2).eq.1) idf1=-idf1
if(abs(idf1).gt.ntol) exit
fpk=idf1 + nrxfreq
call tweak1(cdat2,iz,1500.0-fpk,cdat)
call syncmsk(cdat,iz,jpk,ipk,idf,rmax,snr,metric,msg)
if(metric.eq.-9999) cycle !No output if no significant sync
if(msg(1:1).eq.' ') call jtmsk_short(cdat,iz,narg,tbest,idfpk,msg)
if(msg(1:1).eq.'<' .and. naggressive.eq.0 .and. &
narg(13)/8.ne.narg(12)) msg=' '
if(msg(1:1).ne.' ') then
if(msg.ne.msg0) then
nline=nline+1
nsnr0=-99
endif
freq=fpk+idf
t0=(ia+jpk)/12000.0
y=10.0**(0.1*(yellow(n)-1.5))
nsnr=max(-5,nint(db(y)))
if(nsnr.gt.nsnr0 .and. nline.gt.0) then
call rectify_msk(cdat2(jpk:jpk+NSPM-1),msg,narg(13),freq2)
freq=freq2
if(msg(1:1).eq.'<') freq=freq2+idfpk
!### Check freq values !!!
write(line(nline),1020) nutc,nsnr,t0,nint(freq),msg
1020 format(i6.6,i4,f5.1,i5,' & ',a22)
nsnr0=nsnr
go to 900
endif
msg0=msg
if(nline.ge.maxlines) go to 900
endif
enddo
! print*,'c',nutc,n,nint(yellow(n)-4.0),freq,freq2
enddo
900 continue
! print*,'d',nutc,n,nint(yellow(n)-4.0),freq,freq2
if(line(1)(1:6).eq.' ') line(1)(1:1)=char(0)
return
end subroutine jtmsk_decode

270
lib/syncmsk144.f90 Normal file
View File

@ -0,0 +1,270 @@
subroutine syncmsk(cdat,npts,jpk,ipk,idf,rmax,snr,metric,decoded)
! Attempt synchronization, and if successful decode using Viterbi algorithm.
use iso_c_binding, only: c_loc,c_size_t
use packjt
use hashing
use timer_module, only: timer
parameter (NSPM=840,NSAVE=2000)
character*85 pchk_file,gen_file
complex cdat(npts) !Analytic signal
complex cb(72) !Complex waveform for bit-reversed B11
complex cd(0:11,0:3)
complex c(0:NSPM-1) !Complex data for one message length
complex c2(0:NSPM-1)
complex cb3(1:NSPM,3)
real r(12000)
real rdat(12000)
real ss1(12000)
real symbol(140)
real rdata(140)
real rd2(128)
real rsave(NSAVE)
real xp(29)
real pp(12) !Half-sine pulse shape
real cbi(72),cbq(72)
complex z,z0,z1,z2,z3,cfac
integer*1 e1(128)
integer*1, target :: d8(13)
integer*1 i1hash(4)
integer*1 i1
integer*4 i4Msg6BitWords(12) !72-bit message as 6-bit words
integer mettab(0:255,0:1) !Metric table for BPSK modulation
integer ipksave(NSAVE)
integer jpksave(NSAVE)
integer indx(NSAVE)
integer b11(11) !Barker-11 code
character*22 decoded
character*72 c72
logical first
equivalence (i1,i4)
equivalence (ihash,i1hash)
data first/.true./
data b11/1,1,1,0,0,0,1,0,0,1,0/
save first,cb,cd,pi,twopi,dt,f0,f1,mettab
phi=0.
if(first) then
pchk_file="peg-128-80-reg3.pchk"
gen_file="peg-128-80-reg3.gen"
call init_ldpc(trim(pchk_file)//char(0),trim(gen_file)//char(0))
pi=4.0*atan(1.0)
do i=1,12
angle=(i-1)*pi/12.0
pp(i)=sin(angle)
enddo
cbi(1:12)=(2*b11(11)-1)*pp
cbi(13:24)=(2*b11(9)-1)*pp
cbi(25:36)=(2*b11(7)-1)*pp
cbi(37:48)=(2*b11(5)-1)*pp
cbi(49:60)=(2*b11(3)-1)*pp
cbi(61:72)=(2*b11(1)-1)*pp
cbq(1:6)=0.0
cbq(7:18)=(2*b11(10)-1)*pp
cbq(19:30)=(2*b11(8)-1)*pp
cbq(31:42)=(2*b11(6)-1)*pp
cbq(43:54)=(2*b11(4)-1)*pp
cbq(55:66)=(2*b11(2)-1)*pp
cbq(67:72)=0.0
print*,"cb"
do i=1,72
print*,i,cbi(i),cbq(i)
cb(i)=complex(cbi(i),cbq(i))
enddo
first=.false.
endif
nfft=NSPM
jz=npts-nfft
decoded=" "
ipk=0
jpk=0
metric=-9999
r=0.
open(unit=78,file="blah.txt")
call timer('sync1 ',0)
do j=1,jz !Find the Barker-11 sync vectors
z=0.
ss=0.
do i=1, 72
ss=ss + real(cdat(j+i-1))**2 + aimag(cdat(j+i-1))**2
z=z + cdat(j+i-1)*conjg(cb(i)) !Signal matching Barker 11
enddo
ss=sqrt(ss/72.0)*72.0
r(j)=abs(z)/ss !Goodness-of-fit to Barker 11
ss1(j)=ss
write(78,*) j,r(j),ss1(j)
enddo
z=0
dt=1.0/12000.0
do j=7,jz
z=z+(cdat(j)*conjg(cdat(j-6)))**2
enddo
fff=(atan2(aimag(z),real(z))-pi)/(4*pi*6*dt)
print*,"Frequency error ",fff
call timer('sync1 ',1)
call timer('sync2 ',0)
jz=npts-nfft
rmax=0.
! n1=35, n2=69, n3=94
k=0
do j=1,jz !Find best full-message sync
if(ss1(j).lt.85.0) cycle
r1=r(j) + r(j+282) + r(j+768) ! 6*(12+n1) 6*(24+n1+n2)
r2=r(j) + r(j+486) + r(j+1122) ! 6*(12+n2) 6*(24+n2+n3)
r3=r(j) + r(j+636) + r(j+918) ! 6*(12+n3) 6*(24+n3+n1)
if(r1.gt.rmax) then
rmax=r1
jpk=j
ipk=1
endif
if(r2.gt.rmax) then
rmax=r2
jpk=j
ipk=2
endif
if(r3.gt.rmax) then
rmax=r3
jpk=j
ipk=3
endif
rrmax=max(r1,r2,r3)
if(rrmax.gt.1.9) then
k=min(k+1,NSAVE)
if(r1.eq.rrmax) ipksave(k)=1
if(r2.eq.rrmax) ipksave(k)=2
if(r3.eq.rrmax) ipksave(k)=3
jpksave(k)=j
rsave(k)=rrmax
endif
enddo
call timer('sync2 ',1)
kmax=k
call indexx(rsave,kmax,indx)
call timer('sync3 ',0)
do kk=1,kmax
k=indx(kmax+1-kk)
ipk=ipksave(k)
jpk=jpksave(k)
rmax=rsave(k)
c=conjg(cb3(1:NSPM,ipk))*cdat(jpk:jpk+nfft-1)
smax=0.
dfx=0.
idfbest=0
do itry=1,25
idf=itry/2
if(mod(itry,2).eq.0) idf=-idf
idf=4*idf
twk=idf
call tweak1(c,NSPM,-twk,c2)
z=sum(c2)
if(abs(z).gt.smax) then
dfx=twk
smax=abs(z)
phi=atan2(aimag(z),real(z)) !Carrier phase offset
idfbest=idf
endif
enddo
idf=idfbest
call tweak1(cdat,npts,-dfx,cdat)
cfac=cmplx(cos(phi),-sin(phi))
cdat=cfac*cdat
sig=0.
ref=0.
rdat(1:npts)=cdat
iz=11
do k=1,234 !Compute soft symbols
j=jpk+6*(k-1)
z0=2.0*dot_product(cdat(j:j+iz),cd(0:iz,0))
z1=2.0*dot_product(cdat(j:j+iz),cd(0:iz,1))
z2=2.0*dot_product(cdat(j:j+iz),cd(0:iz,2))
z3=2.0*dot_product(cdat(j:j+iz),cd(0:iz,3))
!### Maybe these should be weighted by yellow() ?
if(j+1404+iz.lt.npts) then
z0=z0 + dot_product(cdat(j+1404:j+1404+iz),cd(0:iz,0))
z1=z1 + dot_product(cdat(j+1404:j+1404+iz),cd(0:iz,1))
z2=z2 + dot_product(cdat(j+1404:j+1404+iz),cd(0:iz,2))
z3=z3 + dot_product(cdat(j+1404:j+1404+iz),cd(0:iz,3))
endif
if(j-1404.ge.1) then
z0=z0 + dot_product(cdat(j-1404:j-1404+iz),cd(0:iz,0))
z1=z1 + dot_product(cdat(j-1404:j-1404+iz),cd(0:iz,1))
z2=z2 + dot_product(cdat(j-1404:j-1404+iz),cd(0:iz,2))
z3=z3 + dot_product(cdat(j-1404:j-1404+iz),cd(0:iz,3))
endif
sym=max(abs(real(z2)),abs(real(z3))) - max(abs(real(z0)),abs(real(z1)))
if(sym.lt.0.0) then
phi=atan2(aimag(z0),real(z0))
sig=sig + real(z0)**2
ref=ref + aimag(z0)**2
else
phi=atan2(aimag(z1),real(z1))
sig=sig + real(z1)**2
ref=ref + aimag(z1)**2
endif
n=k
if(ipk.eq.2) n=k+47
if(ipk.eq.3) n=k+128
if(n.gt.234) n=n-234
ibit=0
if(sym.ge.0) ibit=1
symbol(n)=sym
enddo
snr=db(sig/ref-1.0)
! rdata(1:35)=symbol(12:46)
! rdata(36:104)=symbol(59:127)
! rdata(105:198)=symbol(140:233)
! Re-order the symbols and make them i*1
j=0
do i=1,99
i4=128+rdata(i) !### Should be nint() ??? ###
if(i4.gt.255) i4=255
if(i4.lt.0) i4=0
j=j+1
e1(j)=i1
rd2(j)=rdata(i)
i4=128+rdata(i+99)
if(i4.gt.255) i4=255
if(i4.lt.0) i4=0
j=j+1
e1(j)=i1
rd2(j)=rdata(i+99)
enddo
! Decode the message
nb1=87
call vit213(e1,nb1,mettab,d8,metric)
ihash=nhash(c_loc(d8),int(9,c_size_t),146)
ihash=2*iand(ihash,32767)
decoded=' '
if(d8(10).eq.i1hash(2) .and. d8(11).eq.i1hash(1)) then
write(c72,1012) d8(1:9)
1012 format(9b8.8)
read(c72,1014) i4Msg6BitWords
1014 format(12b6.6)
call unpackmsg(i4Msg6BitWords,decoded) !Unpack to get msgsent
endif
if(decoded.ne.' ') exit
enddo
call timer('sync3 ',1)
return
end subroutine syncmsk