WSJT-X/lib/ft8/foxgen.f90

172 lines
4.0 KiB
Fortran

subroutine foxgen()
use crc
parameter (NN=79,ND=58,KK=87,NSPS=4*1920)
parameter (NWAVE=NN*NSPS,NFFT=614400,NH=NFFT/2)
character*32 cmsg
character*22 msg,msgsent
character*6 mygrid
character*87 cbits
character*88 cb88
logical bcontest,checksumok
integer itone(NN)
integer icos7(0:6)
integer*1 msgbits(KK),codeword(3*ND)
integer*1, target:: i1Msg8BitBytes(11)
integer*1, target:: mycall
real x(NFFT),y(NFFT)
real*8 dt,twopi,f0,fstep,dfreq,phi,dphi
complex cx(0:NH),cy(0:NH)
common/foxcom/wave(NWAVE),nslots,i3bit(5),cmsg(5),mycall(6)
equivalence (x,cx),(y,cy)
data icos7/2,5,6,0,4,1,3/ !Costas 7x7 tone pattern
bcontest=.false.
fstep=60.d0
dfreq=6.25d0
dt=1.d0/48000.d0
twopi=8.d0*atan(1.d0)
wave=0.
mygrid=' '
nrpt=0
do n=1,nslots
!###
! if(n.eq.1) then
! cmsg(n)='W0AAA W3DDD'
! i3bit(n)=0
! endif
! if(n.eq.2) then
! cmsg(n)='W0AAA RR73; W3DDD <K1JT> -12'
! i3bit(n)=1
! endif
!###
i3b=i3bit(n)
if(i3b.eq.0) then
msg=cmsg(n)(1:22)
else
i1=index(cmsg(n),' ')
i2=index(cmsg(n),';')
i3=index(cmsg(n),'<')
i4=index(cmsg(n),'>')
msg=cmsg(n)(1:i1)//cmsg(n)(i2+1:i3-2)//' '
read(cmsg(n)(i4+2:i4+4),*) nrpt
endif
call genft8(msg,mygrid,bcontest,0,msgsent,msgbits,itone)
if(i3b.eq.1) then
icrc10=crc10(c_loc(mycall),6)
ng16=64*icrc10 + nrpt+30
write(cbits,1001) msgbits(1:56),ng16,i3b,0
1001 format(56b1.1,b16.16,b3.3,b12.12)
read(cbits,1002) msgbits
1002 format(87i1)
cb88=cbits//'0'
read(cb88,1003) i1Msg8BitBytes(1:11)
1003 format(11b8)
icrc12=crc12(c_loc(i1Msg8BitBytes),11)
write(cbits,1001) msgbits(1:56),ng16,i3b,icrc12
read(cbits,1002) msgbits
call encode174(msgbits,codeword) !Encode the test message
! Message structure: S7 D29 S7 D29 S7
itone(1:7)=icos7
itone(36+1:36+7)=icos7
itone(NN-6:NN)=icos7
k=7
do j=1,ND
i=3*j -2
k=k+1
if(j.eq.30) k=k+7
itone(k)=codeword(i)*4 + codeword(i+1)*2 + codeword(i+2)
enddo
endif
!###
! call chkcrc12a(msgbits,nbadcrc)
! i3bb=4*msgbits(73) + 2*msgbits(74) + msgbits(75)
! iFreeText=msgbits(57)
! write(*,3001) i3b,i3bb,icrc10,icrc12,nbadcrc,msgsent
!3001 format(5i6,2x,a22)
!###
f0=1500.d0 + fstep*(n-1)
phi=0.d0
k=0
do j=1,NN
f=f0 + dfreq*itone(j)
dphi=twopi*f*dt
do ii=1,NSPS
k=k+1
phi=phi+dphi
xphi=phi
wave(k)=wave(k)+sin(xphi)
enddo
enddo
enddo
sqx=0.
do i=1,NWAVE
sqx=sqx + wave(i)*wave(i)
enddo
sigmax=sqrt(sqx/NWAVE)
wave=wave/sigmax !Force rms=1.0
do i=1,NWAVE
wave(i)=h1(wave(i)) !Compress the waveform
enddo
fac=1.0/maxval(abs(wave)) !Set maxval = 1.0
wave=fac*wave
if(NWAVE.ne.-99) go to 100 !### Omit filtering, for now ###
x(1:k)=wave
x(k+1:)=0.
call four2a(x,nfft,1,-1,0)
nadd=64
k=0
df=48000.0/NFFT
rewind(29)
do i=1,NH/nadd - 1
sx=0.
! sy=0.
do j=1,nadd
k=k+1
sx=sx + real(cx(k))**2 + aimag(cx(k))**2
! sy=sy + real(cy(k))**2 + aimag(cy(k))**2
enddo
freq=df*(k-nadd/2+0.5)
write(29,1022) freq,sx,sy,db(sx)-90.0,db(sy)-90.0
1022 format(f10.3,2e12.3,2f10.3)
if(freq.gt.3000.0) exit
enddo
flush(29)
100 continue
return
end subroutine foxgen
real function h1(x)
! sigma=1.0/sqrt(2.0)
sigma=1.0
xlim=sigma/sqrt(6.0)
ax=abs(x)
sgnx=1.0
if(x.lt.0) sgnx=-1.0
if(ax.le.xlim) then
h1=x
else
z=exp(1.0/6.0 - (ax/sigma)**2)
h1=sgnx*sqrt(6.0)*sigma*(2.0/3.0 - 0.5*z)
endif
return
end function h1