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

Starting to build JTMSK short messages into WSJT-X. Beware: this will NOT 

build correctly.


git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@6413 ab8295b8-cf94-4d9e-aec4-7959e3be5d79
This commit is contained in:
Joe Taylor 2016-01-21 17:53:34 +00:00
parent 16cec2e7df
commit eecf44782c
6 changed files with 281 additions and 139 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

175
lib/JTMSKsim.f90
View File

@ -1,134 +1,93 @@
program JTMSKsim
use wavhdr
parameter (NSPM=1404)
parameter (NFFT=256*1024)
parameter (NMAX=15*12000)
type(hdr) h
complex cb11(0:NSPM-1)
complex cmsg(0:NSPM-1)
complex c(0:NFFT-1)
complex c1(0:NSPM-1)
complex c2(0:NSPM-1)
integer i4tone(NSPM)
integer*2 iwave(NMAX)
real w(0:255)
character*8 arg
real pings(0:NMAX-1)
character arg*8,msg*22,msgsent*22,fname*40
character*3 rpt(0:7)
complex cmsg(0:1404-1) !Waveform of message (once)
complex cwave(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(234) !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.2) then
print*,'Usage: JTMSK freq dB'
if(nargs.ne.5) then
print*,'Usage: JTMSKsim message freq width snr nfiles'
print*,' '
print*,'Examples: JTMSKsim "K1ABC W9XYZ EN37" 1500 0.05 0.0 1'
print*,' JTMSKsim "<K1ABC W9XYZ> R26" 1500 0.01 0.0 10'
go to 999
endif
call getarg(1,arg)
read(arg,*) fmid
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
open(10,file='JTMSKcode.out',status='old')
do j=1,234
read(10,*) junk,i4tone(j)
enddo
close(10)
npts=NMAX
h=default_header(12000,npts)
twopi=8.0*atan(1.0)
dt=1.0/12000.0
df=12000.0/NFFT
k=-1
phi=0.
phi0=0.
call getarg(5,arg)
read(arg,*) nfiles
sig=10.0**(0.05*snrdb)
twopi=8.d0*atan(1.d0)
h=default_header(12000,NMAX)
call init_random_seed() ! seed Fortran RANDOM_NUMBER generator
call sgran() ! see C rand generator (used in gran)
ichk=0
call genmsk(msg,ichk,msgsent,itone,itype) !Check message type
if(itype.lt.1 .or. itype.gt.7) then
print*,'Illegal message'
go to 999
endif
dt=1.d0/12000.d0 !Sample interval
dphi0=twopi*(freq-500.d0)*dt !Phase increment, lower tone
dphi1=twopi*(freq+500.d0)*dt !Phase increment, upper tone
! Generate the Tx waveform
cb11=0.
do j=1,234
dphi=twopi*(fmid-500.0)*dt
if(i4tone(j).eq.1) dphi=twopi*(fmid+500.0)*dt
dphi0=twopi*1000.0*dt
if(i4tone(j).eq.1) dphi0=twopi*2000.0*dt
nsym=234
if(itype.eq.7) nsym=35
nspm=6*nsym !Samples per message
k=-1
phi=0.d0
do j=1,nsym
dphi=dphi0
if(itone(j).eq.1) dphi=dphi1
do i=1,6
k=k+1
phi=phi+dphi
phi=phi + dphi
if(phi.gt.twopi) phi=phi-twopi
cmsg(k)=cmplx(cos(phi),sin(phi))
phi0=phi0+dphi0
if(phi0.gt.twopi) phi0=phi0-twopi
! if((k.ge.1 .and. k.le.66) .or. (k.ge.283 .and. k.le.348) .or. &
! (k.ge.769 .and.k.le.834)) cb11(k)=cmplx(cos(phi0),sin(phi0))
if(k.ge.1 .and. k.le.66) cb11(k)=cmplx(cos(phi0),sin(phi0))
! write(11,3001) k*dt,cmsg(k),phi
!3001 format(4f12.6)
xphi=phi
cmsg(k)=cmplx(cos(xphi),sin(xphi))
enddo
enddo
! Generate noise with B=2500 Hz, rms=1.0
c=0.
ia=nint(250.0/df)
ib=nint(2759.0/df)
do i=ia,ib
x=gran()
y=gran()
c(i)=cmplx(x,y)
enddo
call four2a(c,NFFT,1,1,1)
sq=0.
do i=0,npts-1
sq=sq + real(c(i))**2 + aimag(c(i))**2
enddo
rms=sqrt(0.5*sq/npts)
c=c/rms
call makepings(pings,NMAX,width,sig)
i0=3*12000
ncopy=NSPM
! ncopy=0.5*NSPM
c(i0:i0+ncopy-1)=c(i0:i0+ncopy-1) + sig*cmsg(0:ncopy-1)
do i=1,npts
iwave(i)=100.0*real(c(i))
enddo
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')
open(12,file='150901_000000.wav',status='unknown',access='stream')
write(12) h,iwave(1:npts)
smax=0.
fpk=0.
jpk=0
nfft1=256
w=0.
do i=0,65
w(i)=sin(i*twopi/132.0)
enddo
do ia=0,npts-nfft1
c1(0:nfft1-1)=c(ia:ia+nfft1-1)*conjg(cb11(0:nfft1-1))
c2(0:nfft1-1)=w(0:nfft1-1)*c1(0:nfft1-1)
call four2a(c2,nfft1,1,-1,1)
do i=0,nfft1-1
! write(21,1100) i,c1(i)
!1100 format(i6,2f12.3)
fac=sqrt(6000.0/2500.0)
j=-1
do i=0,NMAX-1
j=j+1
if(j.ge.6*nsym) j=j-6*nsym
xx=0.707*gran()
yy=0.707*gran()
cwave(i)=pings(i)*cmsg(j) + fac*cmplx(xx,yy)
iwave(i)=30.0*real(cwave(i))
! write(88,3003) i,i/12000.d0,cwave(i)
!3003 format(i8,f12.6,2f10.3)
enddo
df1=12000.0/nfft1
do i=-20,20
j=i
if(i.lt.0) j=j+nfft1
f=j*df1
if(i.lt.0) f=f-12000.0
s=1.e-3*(real(c2(j))**2 + aimag(c2(j))**2)
if(abs(ia-i0).lt.1404) write(22,1110) f,c2(j),s,ia
1110 format(4f12.3,i6)
if(s.gt.smax) then
smax=s
jpk=ia
fpk=f
endif
enddo
write(10) h,iwave !Save the .wav file
close(10)
call jtmsk_short(cwave,NMAX,msg)
enddo
print*,smax,jpk*dt,fpk
999 end program JTMSKsim

8
lib/Makefile.jt65
View File

@ -36,11 +36,11 @@ OBJS1 = astrosub.o astro0.o astro.o sun.o coord.o tmoonsub.o \
four2a.o grid2deg.o wisdom.o \
symspec.o analytic.o db.o \
encode232.o interleave9.o\
entail.o fano232.o gran.o sync9.o decjt9.o \
fil3.o decoder.o timer.o exp_decode65.o fqso_first.o \
entail.o fano232.o gran.o sync9.o \
fil3.o decoder.o fqso_first.o \
twkfreq.o symspec2.o shell.o sync65.o peakup.o slope.o xcor.o\
fillcom.o chkss2.o zplot9.o flat1.o flat2.o \
jt65a.o symspec65.o flat65.o ccf65.o decode65a.o \
chkss2.o zplot9.o flat1.o flat2.o \
symspec65.o flat65.o ccf65.o decode65a.o \
filbig.o fil6521.o afc65b.o decode65b.o setup65.o \
extract.o fchisq65.o demod64a.o chkhist.o interleave63.o ccf2.o \
move.o indexx.o graycode65.o twkfreq65.o smo.o smo121.o \

12
lib/genmsk.f90
View File

@ -3,16 +3,18 @@ subroutine genmsk(msg0,ichk,msgsent,i4tone,itype)
! Encode a JTMSK message
! Input:
! - msg0 requested message to be transmitted
! - ichk if nonzero, return only msgsent
! - 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 <call1> <call2> <grid/rpt>
! 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
@ -51,15 +53,15 @@ subroutine genmsk(msg0,ichk,msgsent,i4tone,itype)
enddo
if(message(1:1).eq.'<') then
call genmsk_short(message,msgsent,i4tone,itype)
call genmsk_short(message,msgsent,ichk,i4tone,itype)
if(itype.lt.0) go to 999
i4tone(38)=-37
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.ne.0) go to 999
if(ichk.eq.1) go to 999
call entail(i4Msg6BitWords,i1Msg8BitBytes) !Add tail, make 8-bit bytes
ihash=nhash(c_loc(i1Msg8BitBytes),int(9,c_size_t),146)
ihash=2*iand(ihash,32767) !Generate the CRC

53
lib/genmsk_short.f90
View File

@ -1,14 +1,21 @@
subroutine genmsk_short(msg,msgsent,itone,itype)
subroutine genmsk_short(msg,msgsent,ichk,itone,itype)
use hashing
character*22 msg0,msg,msgsent
character*22 msg,msgsent
character*3 crpt,rpt(0:7)
integer itone(37)
logical first
integer itone(35)
integer ig24(0:4096-1) !Codewords for Golay (24,12) code
integer b11(11)
integer golay_15_3(0:7)
data b11/1,1,1,0,0,0,1,0,0,1,0/ !Barker 11 code
data rpt /'26 ','27 ','28 ','R26','R27','R28','RRR','73 '/
data golay_15_3/00000,07025,11704,14025,19164,20909,26468,31765/
data rpt /'26 ','27 ','28 ','R26','R27','R28','RRR','73 '/
data first/.true./
save first,ig24
if(first) then
call golay24_table(ig24) !Define the Golay(24,12) codewords
first=.false.
endif
itype=-1
msgsent='*** bad message ***'
@ -22,25 +29,27 @@ subroutine genmsk_short(msg,msgsent,itone,itype)
enddo
go to 900
10 itone(1:11)=b11
irpt=i
ncodeword=golay_15_3(irpt) !Save the 15-bit codeword for report
do i=26,12,-1 !Insert codeword into itone array
itone(i)=iand(ncodeword,1)
ncodeword=ncodeword/2
enddo
call hash(msg(2:i1-1),i1-2,ihash)
ihash=iand(ihash,1023) !10-bit hash for the two callsigns
n=ihash
do i=36,27,-1
itone(i)=iand(n,1)
10 irpt=i !Report index, 0-7
if(ichk.lt.10000) then
call hash(msg(2:i1-1),i1-2,ihash)
ihash=iand(ihash,511) !9-bit hash for the two callsigns
ig=8*ihash + irpt !12-bit message information
else
ig=ichk-10000
endif
ncodeword=ig24(ig)
itone(1:11)=b11 !Insert the Barker-11 code
n=2**24
do i=12,35 !Insert codeword into itone array
n=n/2
enddo
n=count(itone(1:36).eq.0)
itone(37)=1
if(iand(n,1).eq.1) itone(37)=0
itone(i)=0
if(iand(ncodeword,n).ne.0) itone(i)=1
enddo
msgsent=msg
itype=7
n=count(itone(1:35).eq.0)
if(mod(n,2).ne.0) stop 'Parity error in genmsk_short.'
900 return
end subroutine genmsk_short

145
lib/jtmsk_short.f90 Normal file
View File

@ -0,0 +1,145 @@
subroutine jtmsk_short(cdat,npts,msg)
parameter (NMAX=15*12000,NSAVE=100)
character*22 msg,decoded,msgsent
character*3 rpt(0:7)
complex cdat(0:npts-1)
complex cw(0:209,0:4096) !Waveforms of possible messages
complex cb11(0:65) !Complex waveform of Barker 11
complex z1,z2a,z2b
real*8 dt,twopi,freq,phi,dphi0,dphi1,dphi
real r1(0:NMAX-1)
real r2(0:4096)
real r1save(NSAVE)
integer*8 count0,count1,clkfreq
integer itone(234) !Message bits
integer jgood(NSAVE)
integer indx(NSAVE)
logical first
data rpt /'26 ','27 ','28 ','R26','R27','R28','RRR','73 '/
data first/.true./
save first,cw
! msg='<W2DEF K3GHI> 73 '
if(first) then
dt=1.d0/12000.d0
twopi=8.d0*atan(1.d0)
freq=1500.d0
dphi0=twopi*(freq-500.d0)*dt !Phase increment, lower tone
dphi1=twopi*(freq+500.d0)*dt !Phase increment, upper tone
nsym=35 !Number of symbols
nspm=6*nsym !Samples per message
do imsg=0,4096 !Generate all possible message waveforms
ichk=0
if(imsg.lt.4096) ichk=10000+imsg
call genmsk(msg,ichk,msgsent,itone,itype) !Encode the message
k=-1
phi=0.d0
do j=1,nsym
dphi=dphi0
if(itone(j).eq.1) dphi=dphi1
do i=1,6
k=k+1
phi=phi + dphi
if(phi.gt.twopi) phi=phi-twopi
xphi=phi
cw(k,imsg)=cmplx(cos(xphi),sin(xphi))
enddo
enddo
enddo
first=.false.
endif
r1thresh=0.40
r2thresh=0.50
rmax=0.9
ngood=0
nbad=0
maxdecodes=999
cb11=cw(0:65,0)
call system_clock(count0,clkfreq)
r1max=0.
do j=0,npts-210 !Find the B11 sync vectors
z1=0.
ss=0.
do i=0,65
ss=ss + real(cdat(j+i))**2 + aimag(cdat(j+i))**2
z1=z1 + cdat(j+i)*conjg(cb11(i)) !Signal matching B11
enddo
ss=sqrt(ss/66.0)*66.0
r1(j)=abs(z1)/(0.908*ss) !Goodness-of-fit to B11
if(r1(j).gt.r1max) then
r1max=r1(j)
jpk=j
endif
enddo
k=0
do j=1,npts-211
if(r1(j).gt.r1thresh .and. r1(j).ge.r1(j-1) .and. r1(j).ge.r1(j+1) ) then
k=k+1
jgood(k)=j
r1save(k)=r1(j)
if(k.ge.NSAVE) exit
endif
enddo
kmax=k
call indexx(r1save,kmax,indx)
r2bad=0.
do kk=1,kmax
k=indx(kmax+1-kk)
j=jgood(k)
u1=0.
u2=0.
r2max=0.
do imsg=0,4096
ssa=0.
ssb=0.
do i=0,209
ssa=ssa + real(cdat(j+i))**2 + aimag(cdat(j+i))**2
ssb=ssb + real(cdat(j+i-144))**2 + aimag(cdat(j+i-144))**2
enddo
z2a=dot_product(cw(0:209,imsg),cdat(j:j+209))
z2b=dot_product(cw(0:65,imsg),cdat(j:j+65)) + &
dot_product(cw(66:209,imsg),cdat(j-144:j-1))
ssa=sqrt(ssa/210.0)*210.0
ssb=sqrt(ssb/210.0)*210.0
r2(imsg)=max(abs(z2a)/(0.908*ssa),abs(z2b)/(0.908*ssb))
if(r2(imsg).gt.r2max) then
r2max=r2(imsg)
ibest=imsg
u2=u1
u1=r2max
endif
enddo
r1or2=r1(j)/r2max
if(r2max.ge.r2thresh .and. u2/u1.lt.rmax .and. r1or2.ge.0.91) then
t=j/12000.0
n=0
irpt=iand(ibest,7)
decoded="<...> "//rpt(irpt)
if(r2max.eq.r2(4096)) then
n=1
decoded=msg(1:14)//rpt(irpt)
endif
if(n.eq.0) nbad=nbad+1
if(n.eq.1) ngood=ngood+1
if(n.eq.0 .and. r2max.gt.r2bad) r2bad=r2max
write(52,3020) k,t,ibest,r1(j),r2max,u2/u1,r1or2,n,decoded
3020 format(i3,f9.4,i5,4f7.2,i2,1x,a22)
if(ngood+nbad.ge.maxdecodes) exit
endif
enddo
call system_clock(count1,clkfreq)
t=float(count1-count0)/float(clkfreq)
print "('Decode: ',f6.3)",t
print "('Worst false decode:',f6.3)",r2bad
print "('Good:',i3,' Bad:',i3)",ngood,nbad
return
end subroutine jtmsk_short

27
lib/makepings.f90 Normal file
View File

@ -0,0 +1,27 @@
subroutine makepings(pings,npts,width,sig)
real pings(npts)
real*8 t
real t0(14)
iping0=-999
dt=1.0/12000.0
do i=1,14
t0(i)=i !Make pings at t=1, 2, ... 14 s.
enddo
w=width
amp=sig
do i=1,npts
iping=min(max(1,i/12000),14)
t=(i*dt-t0(iping))/w
if(t.lt.0.d0 .and. t.lt.10.0) then
fac=0.
else
fac=2.718*t*dexp(-t)
endif
pings(i)=fac*amp
enddo
return
end subroutine makepings