mirror of
https://github.com/saitohirga/WSJT-X.git
synced 2025-05-30 05:12:26 -04:00
Add routines needed to test (240,101) qso mode.
This commit is contained in:
Steven Franke
parent
c2bcb460e1
commit
73e2aa9573
@ -603,17 +603,17 @@ set (wsjt_FSRCS
|
|||||||
lib/fsk4hf/genwsprcpm.f90
|
lib/fsk4hf/genwsprcpm.f90
|
||||||
lib/fsk4hf/encode204.f90
|
lib/fsk4hf/encode204.f90
|
||||||
lib/fsk4hf/decode174_74.f90
|
lib/fsk4hf/decode174_74.f90
|
||||||
lib/fsk4hf/decode174_101.f90
|
lib/fsk4hf/decode240_101.f90
|
||||||
lib/fsk4hf/ldpcsim174_91.f90
|
lib/fsk4hf/ldpcsim174_91.f90
|
||||||
lib/fsk4hf/ldpcsim174_74.f90
|
lib/fsk4hf/ldpcsim174_74.f90
|
||||||
lib/fsk4hf/ldpcsim174_101.f90
|
lib/fsk4hf/ldpcsim240_101.f90
|
||||||
lib/fsk4hf/get_crc24.f90
|
lib/fsk4hf/get_crc24.f90
|
||||||
lib/fsk4hf/encode174_74.f90
|
lib/fsk4hf/encode174_74.f90
|
||||||
lib/fsk4hf/encode174_101.f90
|
lib/fsk4hf/encode240_101.f90
|
||||||
lib/fsk4hf/bpdecode174_74.f90
|
lib/fsk4hf/bpdecode174_74.f90
|
||||||
lib/fsk4hf/bpdecode174_101.f90
|
lib/fsk4hf/bpdecode240_101.f90
|
||||||
lib/fsk4hf/osd174_74.f90
|
lib/fsk4hf/osd174_74.f90
|
||||||
lib/fsk4hf/osd174_101.f90
|
lib/fsk4hf/osd240_101.f90
|
||||||
lib/fsk4hf/wspr4sim.f90
|
lib/fsk4hf/wspr4sim.f90
|
||||||
lib/fsk4hf/genwspr4.f90
|
lib/fsk4hf/genwspr4.f90
|
||||||
lib/fsk4hf/gen_wspr4wave.f90
|
lib/fsk4hf/gen_wspr4wave.f90
|
||||||
@ -621,7 +621,9 @@ set (wsjt_FSRCS
|
|||||||
lib/fsk4hf/get_wspr4_bitmetrics.f90
|
lib/fsk4hf/get_wspr4_bitmetrics.f90
|
||||||
lib/fsk4hf/ft4slowsim.f90
|
lib/fsk4hf/ft4slowsim.f90
|
||||||
lib/fsk4hf/genft4slow.f90
|
lib/fsk4hf/genft4slow.f90
|
||||||
lib/fsk4hf/decode174_101.f90
|
lib/fsk4hf/decode240_101.f90
|
||||||
|
lib/fsk4hf/ft4sd.f90
|
||||||
|
lib/fsk4hf/get_ft4s_bitmetrics.f90
|
||||||
)
|
)
|
||||||
|
|
||||||
# temporary workaround for a gfortran v7.3 ICE on Fedora 27 64-bit
|
# temporary workaround for a gfortran v7.3 ICE on Fedora 27 64-bit
|
||||||
@ -1368,8 +1370,8 @@ target_link_libraries (ldpcsim174_91 wsjt_fort wsjt_cxx)
|
|||||||
add_executable (ldpcsim174_74 lib/fsk4hf/ldpcsim174_74.f90 wsjtx.rc)
|
add_executable (ldpcsim174_74 lib/fsk4hf/ldpcsim174_74.f90 wsjtx.rc)
|
||||||
target_link_libraries (ldpcsim174_74 wsjt_fort wsjt_cxx)
|
target_link_libraries (ldpcsim174_74 wsjt_fort wsjt_cxx)
|
||||||
|
|
||||||
add_executable (ldpcsim174_101 lib/fsk4hf/ldpcsim174_101.f90 wsjtx.rc)
|
add_executable (ldpcsim240_101 lib/fsk4hf/ldpcsim240_101.f90 wsjtx.rc)
|
||||||
target_link_libraries (ldpcsim174_101 wsjt_fort wsjt_cxx)
|
target_link_libraries (ldpcsim240_101 wsjt_fort wsjt_cxx)
|
||||||
|
|
||||||
add_executable (wspr4sim lib/fsk4hf/wspr4sim.f90 wsjtx.rc)
|
add_executable (wspr4sim lib/fsk4hf/wspr4sim.f90 wsjtx.rc)
|
||||||
target_link_libraries (wspr4sim wsjt_fort wsjt_cxx)
|
target_link_libraries (wspr4sim wsjt_fort wsjt_cxx)
|
||||||
@ -1380,6 +1382,9 @@ target_link_libraries (wspr4d wsjt_fort wsjt_cxx)
|
|||||||
add_executable (ft4slowsim lib/fsk4hf/ft4slowsim.f90 wsjtx.rc)
|
add_executable (ft4slowsim lib/fsk4hf/ft4slowsim.f90 wsjtx.rc)
|
||||||
target_link_libraries (ft4slowsim wsjt_fort wsjt_cxx)
|
target_link_libraries (ft4slowsim wsjt_fort wsjt_cxx)
|
||||||
|
|
||||||
|
add_executable (ft4sd lib/fsk4hf/ft4sd.f90 wsjtx.rc)
|
||||||
|
target_link_libraries (ft4sd wsjt_fort wsjt_cxx)
|
||||||
|
|
||||||
endif(WSJT_BUILD_UTILS)
|
endif(WSJT_BUILD_UTILS)
|
||||||
|
|
||||||
# build the main application
|
# build the main application
|
||||||
|
|||||||
111
lib/fsk4hf/bpdecode240_101.f90
Normal file
111
lib/fsk4hf/bpdecode240_101.f90
Normal file
@ -0,0 +1,111 @@
|
|||||||
|
subroutine bpdecode240_101(llr,apmask,maxiterations,message101,cw,nharderror,iter,ncheck)
|
||||||
|
!
|
||||||
|
! A log-domain belief propagation decoder for the (240,101) code.
|
||||||
|
!
|
||||||
|
integer, parameter:: N=240, K=101, M=N-K
|
||||||
|
integer*1 cw(N),apmask(N)
|
||||||
|
integer*1 decoded(K)
|
||||||
|
integer*1 message101(101)
|
||||||
|
integer nrw(M),ncw
|
||||||
|
integer Nm(6,M)
|
||||||
|
integer Mn(3,N) ! 3 checks per bit
|
||||||
|
integer synd(M)
|
||||||
|
real tov(3,N)
|
||||||
|
real toc(6,M)
|
||||||
|
real tanhtoc(6,M)
|
||||||
|
real zn(N)
|
||||||
|
real llr(N)
|
||||||
|
real Tmn
|
||||||
|
|
||||||
|
include "ldpc_240_101_parity.f90"
|
||||||
|
|
||||||
|
decoded=0
|
||||||
|
toc=0
|
||||||
|
tov=0
|
||||||
|
tanhtoc=0
|
||||||
|
! initialize messages to checks
|
||||||
|
do j=1,M
|
||||||
|
do i=1,nrw(j)
|
||||||
|
toc(i,j)=llr((Nm(i,j)))
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
ncnt=0
|
||||||
|
nclast=0
|
||||||
|
do iter=0,maxiterations
|
||||||
|
! Update bit log likelihood ratios (tov=0 in iteration 0).
|
||||||
|
do i=1,N
|
||||||
|
if( apmask(i) .ne. 1 ) then
|
||||||
|
zn(i)=llr(i)+sum(tov(1:ncw,i))
|
||||||
|
else
|
||||||
|
zn(i)=llr(i)
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! Check to see if we have a codeword (check before we do any iteration).
|
||||||
|
cw=0
|
||||||
|
where( zn .gt. 0. ) cw=1
|
||||||
|
ncheck=0
|
||||||
|
do i=1,M
|
||||||
|
synd(i)=sum(cw(Nm(1:nrw(i),i)))
|
||||||
|
if( mod(synd(i),2) .ne. 0 ) ncheck=ncheck+1
|
||||||
|
! if( mod(synd(i),2) .ne. 0 ) write(*,*) 'check ',i,' unsatisfied'
|
||||||
|
enddo
|
||||||
|
if( ncheck .eq. 0 ) then ! we have a codeword - if crc is good, return it
|
||||||
|
decoded=cw(1:101)
|
||||||
|
call get_crc24(decoded,101,nbadcrc)
|
||||||
|
nharderror=count( (2*cw-1)*llr .lt. 0.0 )
|
||||||
|
if(nbadcrc.eq.0) then
|
||||||
|
message101=decoded(1:101)
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
endif
|
||||||
|
|
||||||
|
if( iter.gt.0 ) then ! this code block implements an early stopping criterion
|
||||||
|
! if( iter.gt.10000 ) then ! this code block implements an early stopping criterion
|
||||||
|
nd=ncheck-nclast
|
||||||
|
if( nd .lt. 0 ) then ! # of unsatisfied parity checks decreased
|
||||||
|
ncnt=0 ! reset counter
|
||||||
|
else
|
||||||
|
ncnt=ncnt+1
|
||||||
|
endif
|
||||||
|
! write(*,*) iter,ncheck,nd,ncnt
|
||||||
|
if( ncnt .ge. 5 .and. iter .ge. 10 .and. ncheck .gt. 15) then
|
||||||
|
nharderror=-1
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
endif
|
||||||
|
nclast=ncheck
|
||||||
|
|
||||||
|
! Send messages from bits to check nodes
|
||||||
|
do j=1,M
|
||||||
|
do i=1,nrw(j)
|
||||||
|
ibj=Nm(i,j)
|
||||||
|
toc(i,j)=zn(ibj)
|
||||||
|
do kk=1,ncw ! subtract off what the bit had received from the check
|
||||||
|
if( Mn(kk,ibj) .eq. j ) then
|
||||||
|
toc(i,j)=toc(i,j)-tov(kk,ibj)
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! send messages from check nodes to variable nodes
|
||||||
|
do i=1,M
|
||||||
|
tanhtoc(1:6,i)=tanh(-toc(1:6,i)/2)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
do j=1,N
|
||||||
|
do i=1,ncw
|
||||||
|
ichk=Mn(i,j) ! Mn(:,j) are the checks that include bit j
|
||||||
|
Tmn=product(tanhtoc(1:nrw(ichk),ichk),mask=Nm(1:nrw(ichk),ichk).ne.j)
|
||||||
|
call platanh(-Tmn,y)
|
||||||
|
! y=atanh(-Tmn)
|
||||||
|
tov(i,j)=2*y
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
enddo
|
||||||
|
nharderror=-1
|
||||||
|
return
|
||||||
|
end subroutine bpdecode240_101
|
||||||
133
lib/fsk4hf/decode240_101.f90
Normal file
133
lib/fsk4hf/decode240_101.f90
Normal file
@ -0,0 +1,133 @@
|
|||||||
|
subroutine decode240_101(llr,Keff,ndeep,apmask,maxsuper,message101,cw,nharderror,iter,ncheck,dmin,isuper)
|
||||||
|
!
|
||||||
|
! A hybrid bp/osd decoder for the (240,101) code.
|
||||||
|
!
|
||||||
|
integer, parameter:: N=240, K=101, M=N-K
|
||||||
|
integer*1 cw(N),apmask(N)
|
||||||
|
integer*1 decoded(K)
|
||||||
|
integer*1 nxor(N),hdec(N)
|
||||||
|
integer*1 message101(101)
|
||||||
|
integer nrw(M),ncw
|
||||||
|
integer Nm(6,M)
|
||||||
|
integer Mn(3,N) ! 3 checks per bit
|
||||||
|
integer synd(M)
|
||||||
|
real tov(3,N)
|
||||||
|
real toc(6,M)
|
||||||
|
real tanhtoc(6,M)
|
||||||
|
real zn(N),zsum(N)
|
||||||
|
real llr(N)
|
||||||
|
real Tmn
|
||||||
|
|
||||||
|
include "ldpc_240_101_parity.f90"
|
||||||
|
|
||||||
|
decoded=0
|
||||||
|
toc=0
|
||||||
|
tov=0
|
||||||
|
tanhtoc=0
|
||||||
|
! initialize messages to checks
|
||||||
|
do j=1,M
|
||||||
|
do i=1,nrw(j)
|
||||||
|
toc(i,j)=llr((Nm(i,j)))
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
ncnt=0
|
||||||
|
nclast=0
|
||||||
|
|
||||||
|
maxiterations=1
|
||||||
|
|
||||||
|
zsum=0.0
|
||||||
|
do isuper=1,maxsuper
|
||||||
|
|
||||||
|
do iter=0,maxiterations
|
||||||
|
! Update bit log likelihood ratios (tov=0 in iteration 0).
|
||||||
|
do i=1,N
|
||||||
|
if( apmask(i) .ne. 1 ) then
|
||||||
|
zn(i)=llr(i)+sum(tov(1:ncw,i))
|
||||||
|
else
|
||||||
|
zn(i)=llr(i)
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
zsum=zsum+zn
|
||||||
|
! Check to see if we have a codeword (check before we do any iteration).
|
||||||
|
cw=0
|
||||||
|
where( zn .gt. 0. ) cw=1
|
||||||
|
ncheck=0
|
||||||
|
do i=1,M
|
||||||
|
synd(i)=sum(cw(Nm(1:nrw(i),i)))
|
||||||
|
if( mod(synd(i),2) .ne. 0 ) ncheck=ncheck+1
|
||||||
|
! if( mod(synd(i),2) .ne. 0 ) write(*,*) 'check ',i,' unsatisfied'
|
||||||
|
enddo
|
||||||
|
! write(*,*) 'number of unsatisfied parity checks ',ncheck
|
||||||
|
if( ncheck .eq. 0 ) then ! we have a codeword - if crc is good, return it
|
||||||
|
decoded=cw(1:K)
|
||||||
|
call get_crc24(decoded,74,nbadcrc)
|
||||||
|
nharderror=count( (2*cw-1)*llr .lt. 0.0 )
|
||||||
|
if(nbadcrc.eq.0) then
|
||||||
|
message101=decoded(1:101)
|
||||||
|
dmin=0.0
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
endif
|
||||||
|
|
||||||
|
! if( iter.gt.0 ) then ! this code block implements an early stopping criterion
|
||||||
|
if( iter.gt.10000 ) then ! this code block implements an early stopping criterion
|
||||||
|
nd=ncheck-nclast
|
||||||
|
if( nd .lt. 0 ) then ! # of unsatisfied parity checks decreased
|
||||||
|
ncnt=0 ! reset counter
|
||||||
|
else
|
||||||
|
ncnt=ncnt+1
|
||||||
|
endif
|
||||||
|
! write(*,*) iter,ncheck,nd,ncnt
|
||||||
|
if( ncnt .ge. 5 .and. iter .ge. 10 .and. ncheck .gt. 15) then
|
||||||
|
nharderror=-1
|
||||||
|
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
endif
|
||||||
|
nclast=ncheck
|
||||||
|
|
||||||
|
! Send messages from bits to check nodes
|
||||||
|
do j=1,M
|
||||||
|
do i=1,nrw(j)
|
||||||
|
ibj=Nm(i,j)
|
||||||
|
toc(i,j)=zn(ibj)
|
||||||
|
do kk=1,ncw ! subtract off what the bit had received from the check
|
||||||
|
if( Mn(kk,ibj) .eq. j ) then
|
||||||
|
toc(i,j)=toc(i,j)-tov(kk,ibj)
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! send messages from check nodes to variable nodes
|
||||||
|
do i=1,M
|
||||||
|
tanhtoc(1:6,i)=tanh(-toc(1:6,i)/2)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
do j=1,N
|
||||||
|
do i=1,ncw
|
||||||
|
ichk=Mn(i,j) ! Mn(:,j) are the checks that include bit j
|
||||||
|
Tmn=product(tanhtoc(1:nrw(ichk),ichk),mask=Nm(1:nrw(ichk),ichk).ne.j)
|
||||||
|
call platanh(-Tmn,y)
|
||||||
|
! y=atanh(-Tmn)
|
||||||
|
tov(i,j)=2*y
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
enddo ! bp iterations
|
||||||
|
|
||||||
|
call osd240_101(zsum,Keff,apmask,ndeep,message101,cw,nharderror,dminosd)
|
||||||
|
if(nharderror.gt.0) then
|
||||||
|
hdec=0
|
||||||
|
where(llr .ge. 0) hdec=1
|
||||||
|
nxor=ieor(hdec,cw)
|
||||||
|
dmin=sum(nxor*abs(llr))
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
enddo ! super iterations
|
||||||
|
|
||||||
|
nharderror=-1
|
||||||
|
|
||||||
|
return
|
||||||
|
end subroutine decode240_101
|
||||||
46
lib/fsk4hf/encode240_101.f90
Normal file
46
lib/fsk4hf/encode240_101.f90
Normal file
@ -0,0 +1,46 @@
|
|||||||
|
subroutine encode240_101(message,codeword)
|
||||||
|
use, intrinsic :: iso_c_binding
|
||||||
|
use iso_c_binding, only: c_loc,c_size_t
|
||||||
|
use crc
|
||||||
|
|
||||||
|
integer, parameter:: N=240, K=101, M=N-K
|
||||||
|
character*24 c24
|
||||||
|
integer*1 codeword(N)
|
||||||
|
integer*1 gen(M,K)
|
||||||
|
integer*1 message(K)
|
||||||
|
integer*1 pchecks(M)
|
||||||
|
integer*4 ncrc24
|
||||||
|
include "ldpc_240_101_generator.f90"
|
||||||
|
logical first
|
||||||
|
data first/.true./
|
||||||
|
save first,gen
|
||||||
|
|
||||||
|
if( first ) then ! fill the generator matrix
|
||||||
|
gen=0
|
||||||
|
do i=1,M
|
||||||
|
do j=1,26
|
||||||
|
read(g(i)(j:j),"(Z1)") istr
|
||||||
|
ibmax=4
|
||||||
|
if(j.eq.26) ibmax=1
|
||||||
|
do jj=1, ibmax
|
||||||
|
icol=(j-1)*4+jj
|
||||||
|
if( btest(istr,4-jj) ) gen(i,icol)=1
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
first=.false.
|
||||||
|
endif
|
||||||
|
|
||||||
|
do i=1,M
|
||||||
|
nsum=0
|
||||||
|
do j=1,K
|
||||||
|
nsum=nsum+message(j)*gen(i,j)
|
||||||
|
enddo
|
||||||
|
pchecks(i)=mod(nsum,2)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
codeword(1:K)=message
|
||||||
|
codeword(K+1:N)=pchecks
|
||||||
|
|
||||||
|
return
|
||||||
|
end subroutine encode240_101
|
||||||
16
lib/fsk4hf/ft4s_params.f90
Normal file
16
lib/fsk4hf/ft4s_params.f90
Normal file
@ -0,0 +1,16 @@
|
|||||||
|
! FT4A
|
||||||
|
! LDPC(240,101)/CRC24 code, four 4x4 Costas arrays for sync, ramp-up and ramp-down symbols
|
||||||
|
|
||||||
|
parameter (KK=77) !Information bits (77 + CRC24)
|
||||||
|
parameter (ND=120) !Data symbols
|
||||||
|
parameter (NS=24) !Sync symbols
|
||||||
|
parameter (NN=NS+ND) !Sync and data symbols (144)
|
||||||
|
parameter (NN2=NS+ND+2) !Total channel symbols (146)
|
||||||
|
parameter (NSPS=9600) !Samples per symbol at 12000 S/s
|
||||||
|
parameter (NZ=NSPS*NN) !Sync and Data samples (1,382,400)
|
||||||
|
parameter (NZ2=NSPS*NN2) !Total samples in shaped waveform (1,397,760)
|
||||||
|
parameter (NMAX=408*3456) !Samples in iwave (1,410,048)
|
||||||
|
parameter (NFFT1=4*NSPS, NH1=NFFT1/2) !Length of FFTs for symbol spectra
|
||||||
|
parameter (NSTEP=NSPS) !Coarse time-sync step size
|
||||||
|
parameter (NHSYM=(NMAX-NFFT1)/NSTEP) !Number of symbol spectra (1/4-sym steps)
|
||||||
|
parameter (NDOWN=32) !Downsample factor
|
||||||
459
lib/fsk4hf/ft4sd.f90
Normal file
459
lib/fsk4hf/ft4sd.f90
Normal file
@ -0,0 +1,459 @@
|
|||||||
|
program ft4sd
|
||||||
|
|
||||||
|
! Decode ft4slow data read from *.c2 or *.wav files.
|
||||||
|
|
||||||
|
use packjt77
|
||||||
|
include 'ft4s_params.f90'
|
||||||
|
parameter (NSPS2=NSPS/32)
|
||||||
|
character arg*8,cbits*50,infile*80,fname*16,datetime*11
|
||||||
|
character ch1*1,ch4*4,cseq*31
|
||||||
|
character*22 decodes(100)
|
||||||
|
character*37 msg
|
||||||
|
character*120 data_dir
|
||||||
|
character*77 c77
|
||||||
|
complex c2(0:NMAX/32-1) !Complex waveform
|
||||||
|
complex cframe(0:144*NSPS2-1) !Complex waveform
|
||||||
|
complex cd(0:144*20-1) !Complex waveform
|
||||||
|
real*8 fMHz
|
||||||
|
real llr(240),llra(240),llrb(240),llrc(240),llrd(240)
|
||||||
|
real candidates(100,2)
|
||||||
|
real bitmetrics(288,4)
|
||||||
|
integer ihdr(11)
|
||||||
|
integer*2 iwave(NMAX) !Generated full-length waveform
|
||||||
|
integer*1 apmask(240),cw(240)
|
||||||
|
integer*1 hbits(288)
|
||||||
|
integer*1 message101(101)
|
||||||
|
logical badsync,unpk77_success
|
||||||
|
|
||||||
|
fs=12000.0/NDOWN !Sample rate
|
||||||
|
dt=1.0/fs !Sample interval (s)
|
||||||
|
tt=NSPS*dt !Duration of "itone" symbols (s)
|
||||||
|
txt=NZ*dt !Transmission length (s)
|
||||||
|
|
||||||
|
nargs=iargc()
|
||||||
|
if(nargs.lt.1) then
|
||||||
|
print*,'Usage: ft4sd [-a <data_dir>] [-f fMHz] file1 [file2 ...]'
|
||||||
|
go to 999
|
||||||
|
endif
|
||||||
|
iarg=1
|
||||||
|
data_dir="."
|
||||||
|
call getarg(iarg,arg)
|
||||||
|
if(arg(1:2).eq.'-a') then
|
||||||
|
call getarg(iarg+1,data_dir)
|
||||||
|
iarg=iarg+2
|
||||||
|
endif
|
||||||
|
call getarg(iarg,arg)
|
||||||
|
if(arg(1:2).eq.'-f') then
|
||||||
|
call getarg(iarg+1,arg)
|
||||||
|
read(arg,*) fMHz
|
||||||
|
iarg=iarg+2
|
||||||
|
endif
|
||||||
|
|
||||||
|
ngood=0
|
||||||
|
do ifile=iarg,nargs
|
||||||
|
call getarg(ifile,infile)
|
||||||
|
open(10,file=infile,status='old',access='stream')
|
||||||
|
j1=index(infile,'.c2')
|
||||||
|
j2=index(infile,'.wav')
|
||||||
|
if(j1.gt.0) then
|
||||||
|
read(10,end=999) fname,ntrmin,fMHz,c2
|
||||||
|
read(fname(8:11),*) nutc
|
||||||
|
write(datetime,'(i11)') nutc
|
||||||
|
else if(j2.gt.0) then
|
||||||
|
read(10,end=999) ihdr,iwave
|
||||||
|
read(infile(j2-4:j2-1),*) nutc
|
||||||
|
datetime=infile(j2-11:j2-1)
|
||||||
|
call ft4s_downsample(iwave,c2)
|
||||||
|
else
|
||||||
|
print*,'Wrong file format?'
|
||||||
|
go to 999
|
||||||
|
endif
|
||||||
|
close(10)
|
||||||
|
|
||||||
|
fa=-100.0
|
||||||
|
fb=100.0
|
||||||
|
fs=12000.0/32.0
|
||||||
|
npts=120*12000.0/32.0
|
||||||
|
|
||||||
|
call getcandidate_ft4s(c2,npts,fs,fa,fb,ncand,candidates) !First approx for freq
|
||||||
|
|
||||||
|
ndecodes=0
|
||||||
|
do icand=1,ncand
|
||||||
|
fc0=candidates(icand,1)
|
||||||
|
xsnr=candidates(icand,2)
|
||||||
|
!write(*,*) 'candidates ',icand,fc0,xsnr
|
||||||
|
do isync=0,1
|
||||||
|
|
||||||
|
del=1.5*fs/300.0
|
||||||
|
if(isync.eq.0) then
|
||||||
|
fc1=fc0-del
|
||||||
|
is0=375
|
||||||
|
ishw=350
|
||||||
|
isst=30
|
||||||
|
ifhw=10
|
||||||
|
df=.1
|
||||||
|
else if(isync.eq.1) then
|
||||||
|
fc1=fc2
|
||||||
|
is0=isbest
|
||||||
|
ishw=100
|
||||||
|
isst=10
|
||||||
|
ifhw=10
|
||||||
|
df=.02
|
||||||
|
endif
|
||||||
|
smax=0.0
|
||||||
|
do if=-ifhw,ifhw
|
||||||
|
fc=fc1+df*if
|
||||||
|
do istart=max(1,is0-ishw),is0+ishw,isst
|
||||||
|
call coherent_sync_ft4s(c2,istart,fc,1,sync)
|
||||||
|
if(sync.gt.smax) then
|
||||||
|
fc2=fc
|
||||||
|
isbest=istart
|
||||||
|
smax=sync
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
! write(*,*) ifile,icand,isync,fc1+del,fc2+del,isbest,smax
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! if(smax .lt. 100.0 ) cycle
|
||||||
|
!isbest=375
|
||||||
|
!fc2=-del
|
||||||
|
do ijitter=0,2
|
||||||
|
if(ijitter.eq.0) ioffset=0
|
||||||
|
if(ijitter.eq.1) ioffset=45
|
||||||
|
if(ijitter.eq.2) ioffset=-45
|
||||||
|
is0=isbest+ioffset
|
||||||
|
if(is0.lt.0) cycle
|
||||||
|
cframe=c2(is0:is0+144*300-1)
|
||||||
|
call downsample_ft4s(cframe,fc2,cd)
|
||||||
|
s2=sum(cd*conjg(cd))/(20*144)
|
||||||
|
cd=cd/sqrt(s2)
|
||||||
|
call get_ft4s_bitmetrics(cd,bitmetrics,badsync)
|
||||||
|
|
||||||
|
hbits=0
|
||||||
|
where(bitmetrics(:,1).ge.0) hbits=1
|
||||||
|
ns1=count(hbits( 1: 8).eq.(/0,0,0,1,1,0,1,1/))
|
||||||
|
ns2=count(hbits( 57: 64).eq.(/0,1,0,0,1,1,1,0/))
|
||||||
|
ns3=count(hbits(113:120).eq.(/1,1,1,0,0,1,0,0/))
|
||||||
|
ns4=count(hbits(169:176).eq.(/1,0,1,1,0,0,0,1/))
|
||||||
|
ns5=count(hbits(225:232).eq.(/0,0,1,1,1,0,0,1/))
|
||||||
|
ns6=count(hbits(281:288).eq.(/0,1,1,1,0,0,1,0/))
|
||||||
|
nsync_qual=ns1+ns2+ns3+ns4+ns5+ns6
|
||||||
|
! if(nsync_qual.lt. 20) cycle
|
||||||
|
|
||||||
|
scalefac=2.83
|
||||||
|
llra( 1: 48)=bitmetrics( 9: 56, 1)
|
||||||
|
llra( 49: 96)=bitmetrics( 65:112, 1)
|
||||||
|
llra( 97:144)=bitmetrics(121:168, 1)
|
||||||
|
llra(145:192)=bitmetrics(177:224, 1)
|
||||||
|
llra(193:240)=bitmetrics(233:280, 1)
|
||||||
|
llra=scalefac*llra
|
||||||
|
llrb( 1: 48)=bitmetrics( 9: 56, 2)
|
||||||
|
llrb( 49: 96)=bitmetrics( 65:112, 2)
|
||||||
|
llrb( 97:144)=bitmetrics(121:168, 2)
|
||||||
|
llrb(145:192)=bitmetrics(177:224, 2)
|
||||||
|
llrb(193:240)=bitmetrics(233:280, 2)
|
||||||
|
llrb=scalefac*llrb
|
||||||
|
llrc( 1: 48)=bitmetrics( 9: 56, 3)
|
||||||
|
llrc( 49: 96)=bitmetrics( 65:112, 3)
|
||||||
|
llrc( 97:144)=bitmetrics(121:168, 3)
|
||||||
|
llrc(145:192)=bitmetrics(177:224, 3)
|
||||||
|
llrc(193:240)=bitmetrics(233:280, 3)
|
||||||
|
llrc=scalefac*llrc
|
||||||
|
llrd( 1: 48)=bitmetrics( 9: 56, 4)
|
||||||
|
llrd( 49: 96)=bitmetrics( 65:112, 4)
|
||||||
|
llrd( 97:144)=bitmetrics(121:168, 4)
|
||||||
|
llrd(145:192)=bitmetrics(177:224, 4)
|
||||||
|
llrd(193:240)=bitmetrics(233:280, 4)
|
||||||
|
llrd=scalefac*llrd
|
||||||
|
apmask=0
|
||||||
|
max_iterations=40
|
||||||
|
|
||||||
|
do itry=4,1,-1
|
||||||
|
if(itry.eq.1) llr=llra
|
||||||
|
if(itry.eq.2) llr=llrb
|
||||||
|
if(itry.eq.3) llr=llrc
|
||||||
|
if(itry.eq.4) llr=llrd
|
||||||
|
nhardbp=0
|
||||||
|
nhardosd=0
|
||||||
|
dmin=0.0
|
||||||
|
call bpdecode240_101(llr,apmask,max_iterations,message101,cw,nhardbp,niterations,nchecks)
|
||||||
|
Keff=91
|
||||||
|
! if(nhardbp.lt.0) call osd240_101(llr,Keff,apmask,5,message101,cw,nhardosd,dmin)
|
||||||
|
maxsuperits=2
|
||||||
|
ndeep=3 ! use ndeep=3 with Keff=91
|
||||||
|
if(nhardbp.lt.0) then
|
||||||
|
! call osd240_101(llr,Keff,apmask,ndeep,message101,cw,nhardosd,dmin)
|
||||||
|
call decode240_101(llr,Keff,ndeep,apmask,maxsuperits,message101,cw,nhardosd,iter,ncheck,dmin,isuper)
|
||||||
|
endif
|
||||||
|
if(nhardbp.ge.0 .or. nhardosd.ge.0) then
|
||||||
|
write(c77,'(77i1)') message101(1:77)
|
||||||
|
call unpack77(c77,0,msg,unpk77_success)
|
||||||
|
if(unpk77_success .and. index(msg,'K9AN').gt.0) then
|
||||||
|
ngood=ngood+1
|
||||||
|
write(*,1100) ifile,icand,xsnr,isbest/375.0-1.0,1500.0+fc2+del,msg(1:14),itry,nhardbp,nhardosd,dmin,ijitter
|
||||||
|
1100 format(i5,2x,i5,2x,f6.1,2x,f6.2,2x,f8.2,2x,a14,i4,i4,i4,f7.2,i6)
|
||||||
|
goto 2002
|
||||||
|
else
|
||||||
|
cycle
|
||||||
|
endif
|
||||||
|
endif
|
||||||
|
enddo ! metrics
|
||||||
|
enddo ! istart jitter
|
||||||
|
enddo !candidate list
|
||||||
|
2002 continue
|
||||||
|
enddo !files
|
||||||
|
nfiles=nargs-iarg+1
|
||||||
|
write(*,*) 'nfiles: ',nfiles,' ngood: ',ngood
|
||||||
|
write(*,1120)
|
||||||
|
1120 format("<DecodeFinished>")
|
||||||
|
|
||||||
|
999 end program ft4sd
|
||||||
|
|
||||||
|
subroutine coherent_sync_ft4s(cd0,i0,f0,itwk,sync)
|
||||||
|
|
||||||
|
! Compute sync power for a complex, downsampled FT4s signal.
|
||||||
|
|
||||||
|
include 'ft4s_params.f90'
|
||||||
|
parameter(NP=NMAX/NDOWN,NSS=NSPS/NDOWN)
|
||||||
|
complex cd0(0:NP-1)
|
||||||
|
complex csynca(4*NSS),csyncb(4*NSS)
|
||||||
|
complex csyncc(4*NSS),csyncd(4*NSS)
|
||||||
|
complex csynce(4*NSS),csyncf(4*NSS)
|
||||||
|
complex csync2(4*NSS)
|
||||||
|
complex ctwk(4*NSS)
|
||||||
|
complex z1,z2,z3,z4,z5,z6
|
||||||
|
logical first
|
||||||
|
integer icos4a(0:3),icos4b(0:3)
|
||||||
|
integer icos4c(0:3),icos4d(0:3)
|
||||||
|
integer icos4e(0:3),icos4f(0:3)
|
||||||
|
data icos4a/0,1,3,2/
|
||||||
|
data icos4b/1,0,2,3/
|
||||||
|
data icos4c/2,3,1,0/
|
||||||
|
data icos4d/3,2,0,1/
|
||||||
|
data icos4e/0,2,3,1/
|
||||||
|
data icos4f/1,2,0,3/
|
||||||
|
data first/.true./
|
||||||
|
save first,twopi,csynca,csyncb,csyncc,csyncd,csynce,csyncf,fac
|
||||||
|
|
||||||
|
p(z1)=(real(z1*fac)**2 + aimag(z1*fac)**2)**0.5 !Statement function for power
|
||||||
|
|
||||||
|
if( first ) then
|
||||||
|
twopi=8.0*atan(1.0)
|
||||||
|
k=1
|
||||||
|
phia=0.0
|
||||||
|
phib=0.0
|
||||||
|
phic=0.0
|
||||||
|
phid=0.0
|
||||||
|
phie=0.0
|
||||||
|
phif=0.0
|
||||||
|
do i=0,3
|
||||||
|
dphia=twopi*icos4a(i)/real(NSS)
|
||||||
|
dphib=twopi*icos4b(i)/real(NSS)
|
||||||
|
dphic=twopi*icos4c(i)/real(NSS)
|
||||||
|
dphid=twopi*icos4d(i)/real(NSS)
|
||||||
|
dphie=twopi*icos4e(i)/real(NSS)
|
||||||
|
dphif=twopi*icos4f(i)/real(NSS)
|
||||||
|
do j=1,NSS
|
||||||
|
csynca(k)=cmplx(cos(phia),sin(phia))
|
||||||
|
csyncb(k)=cmplx(cos(phib),sin(phib))
|
||||||
|
csyncc(k)=cmplx(cos(phic),sin(phic))
|
||||||
|
csyncd(k)=cmplx(cos(phid),sin(phid))
|
||||||
|
csynce(k)=cmplx(cos(phie),sin(phie))
|
||||||
|
csyncf(k)=cmplx(cos(phif),sin(phif))
|
||||||
|
phia=mod(phia+dphia,twopi)
|
||||||
|
phib=mod(phib+dphib,twopi)
|
||||||
|
phic=mod(phic+dphic,twopi)
|
||||||
|
phid=mod(phid+dphid,twopi)
|
||||||
|
phie=mod(phie+dphie,twopi)
|
||||||
|
phif=mod(phif+dphif,twopi)
|
||||||
|
k=k+1
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
first=.false.
|
||||||
|
fac=1.0/(4.0*NSS)
|
||||||
|
endif
|
||||||
|
|
||||||
|
i1=i0 !four Costas arrays
|
||||||
|
i2=i0+28*NSS
|
||||||
|
i3=i0+56*NSS
|
||||||
|
i4=i0+84*NSS
|
||||||
|
i5=i0+112*NSS
|
||||||
|
i6=i0+140*NSS
|
||||||
|
|
||||||
|
z1=0.
|
||||||
|
z2=0.
|
||||||
|
z3=0.
|
||||||
|
z4=0.
|
||||||
|
z5=0.
|
||||||
|
z6=0.
|
||||||
|
|
||||||
|
if(itwk.eq.1) then
|
||||||
|
dt=1/(12000.0/32.0)
|
||||||
|
dphi=twopi*f0*dt
|
||||||
|
phi=0.0
|
||||||
|
do i=1,4*NSS
|
||||||
|
ctwk(i)=cmplx(cos(phi),sin(phi))
|
||||||
|
phi=mod(phi+dphi,twopi)
|
||||||
|
enddo
|
||||||
|
endif
|
||||||
|
|
||||||
|
if(itwk.eq.1) csync2=ctwk*csynca !Tweak the frequency
|
||||||
|
if(i1.ge.0 .and. i1+4*NSS-1.le.NP-1) then
|
||||||
|
z1=sum(cd0(i1:i1+4*NSS-1)*conjg(csync2))
|
||||||
|
elseif( i1.lt.0 ) then
|
||||||
|
npts=(i1+4*NSS-1)/2
|
||||||
|
if(npts.le.40) then
|
||||||
|
z1=0.
|
||||||
|
else
|
||||||
|
z1=sum(cd0(0:i1+4*NSS-1)*conjg(csync2(4*NSS-npts:)))
|
||||||
|
endif
|
||||||
|
endif
|
||||||
|
|
||||||
|
if(itwk.eq.1) csync2=ctwk*csyncb !Tweak the frequency
|
||||||
|
if(i2.ge.0 .and. i2+4*NSS-1.le.NP-1) then
|
||||||
|
z2=sum(cd0(i2:i2+4*NSS-1)*conjg(csync2))
|
||||||
|
endif
|
||||||
|
|
||||||
|
if(itwk.eq.1) csync2=ctwk*csyncc !Tweak the frequency
|
||||||
|
if(i3.ge.0 .and. i3+4*NSS-1.le.NP-1) then
|
||||||
|
z3=sum(cd0(i3:i3+4*NSS-1)*conjg(csync2))
|
||||||
|
endif
|
||||||
|
|
||||||
|
if(itwk.eq.1) csync2=ctwk*csyncd !Tweak the frequency
|
||||||
|
if(i4.ge.0 .and. i4+4*NSS-1.le.NP-1) then
|
||||||
|
z4=sum(cd0(i4:i4+4*NSS-1)*conjg(csync2))
|
||||||
|
endif
|
||||||
|
|
||||||
|
if(itwk.eq.1) csync2=ctwk*csynce !Tweak the frequency
|
||||||
|
if(i5.ge.0 .and. i5+4*NSS-1.le.NP-1) then
|
||||||
|
z5=sum(cd0(i5:i5+4*NSS-1)*conjg(csync2))
|
||||||
|
endif
|
||||||
|
|
||||||
|
if(itwk.eq.1) csync2=ctwk*csyncf !Tweak the frequency
|
||||||
|
if(i6.ge.0 .and. i6+4*NSS-1.le.NP-1) then
|
||||||
|
z6=sum(cd0(i6:i6+4*NSS-1)*conjg(csync2))
|
||||||
|
elseif( i6+4*NSS-1.gt.NP-1 ) then
|
||||||
|
npts=(NP-1-i6+1)
|
||||||
|
if(npts.le.40) then
|
||||||
|
z6=0.
|
||||||
|
else
|
||||||
|
z6=sum(cd0(i6:i6+npts-1)*conjg(csync2(1:npts)))
|
||||||
|
endif
|
||||||
|
endif
|
||||||
|
|
||||||
|
sync = p(z1) + p(z2) + p(z3) + p(z4) + p(z5) + p(z6)
|
||||||
|
|
||||||
|
return
|
||||||
|
end subroutine coherent_sync_ft4s
|
||||||
|
|
||||||
|
subroutine downsample_ft4s(ci,f0,co)
|
||||||
|
parameter(NI=144*300,NH=NI/2,NO=NI/15) ! downsample from 315 samples per symbol to 20
|
||||||
|
complex ci(0:NI-1),ct(0:NI-1)
|
||||||
|
complex co(0:NO-1)
|
||||||
|
fs=12000.0/32.0
|
||||||
|
df=fs/NI
|
||||||
|
ct=ci
|
||||||
|
call four2a(ct,NI,1,-1,1) !c2c FFT to freq domain
|
||||||
|
i0=nint(f0/df)
|
||||||
|
ct=cshift(ct,i0)
|
||||||
|
co=0.0
|
||||||
|
co(0)=ct(0)
|
||||||
|
b=16.0
|
||||||
|
do i=1,NO/2
|
||||||
|
arg=(i*df/b)**2
|
||||||
|
filt=exp(-arg)
|
||||||
|
co(i)=ct(i)*filt
|
||||||
|
co(NO-i)=ct(NI-i)*filt
|
||||||
|
enddo
|
||||||
|
co=co/NO
|
||||||
|
call four2a(co,NO,1,1,1) !c2c FFT back to time domain
|
||||||
|
return
|
||||||
|
end subroutine downsample_ft4s
|
||||||
|
|
||||||
|
subroutine getcandidate_ft4s(c,npts,fs,fa,fb,ncand,candidates)
|
||||||
|
parameter(NFFT1=120*12000/32,NH1=NFFT1/2,NFFT2=120*12000/320,NH2=NFFT2/2)
|
||||||
|
complex c(0:npts-1) !Complex waveform
|
||||||
|
complex cc(0:NFFT1-1)
|
||||||
|
complex csfil(0:NFFT2-1)
|
||||||
|
complex cwork(0:NFFT2-1)
|
||||||
|
real bigspec(0:NFFT2-1)
|
||||||
|
complex c2(0:NFFT1-1) !Short spectra
|
||||||
|
real s(-NH1+1:NH1) !Coarse spectrum
|
||||||
|
real ss(-NH1+1:NH1) !Smoothed coarse spectrum
|
||||||
|
real candidates(100,2)
|
||||||
|
integer indx(NFFT2-1)
|
||||||
|
logical first
|
||||||
|
data first/.true./
|
||||||
|
save first,w,df,csfil
|
||||||
|
|
||||||
|
if(first) then
|
||||||
|
df=10*fs/NFFT1
|
||||||
|
csfil=cmplx(0.0,0.0)
|
||||||
|
do i=0,NFFT2-1
|
||||||
|
! csfil(i)=exp(-((i-NH2)/32.0)**2) ! revisit this
|
||||||
|
csfil(i)=exp(-((i-NH2)/28.0)**2) ! revisit this
|
||||||
|
enddo
|
||||||
|
csfil=cshift(csfil,NH2)
|
||||||
|
call four2a(csfil,NFFT2,1,-1,1)
|
||||||
|
first=.false.
|
||||||
|
endif
|
||||||
|
|
||||||
|
cc=cmplx(0.0,0.0)
|
||||||
|
cc(0:npts-1)=c;
|
||||||
|
call four2a(cc,NFFT1,1,-1,1)
|
||||||
|
cc=abs(cc)**2
|
||||||
|
call four2a(cc,NFFT1,1,-1,1)
|
||||||
|
cwork(0:NH2)=cc(0:NH2)*conjg(csfil(0:NH2))
|
||||||
|
cwork(NH2+1:NFFT2-1)=cc(NFFT1-NH2+1:NFFT1-1)*conjg(csfil(NH2+1:NFFT2-1))
|
||||||
|
|
||||||
|
call four2a(cwork,NFFT2,1,+1,1)
|
||||||
|
bigspec=cshift(real(cwork),-NH2)
|
||||||
|
il=NH2+fa/df
|
||||||
|
ih=NH2+fb/df
|
||||||
|
nnl=ih-il+1
|
||||||
|
call indexx(bigspec(il:il+nnl-1),nnl,indx)
|
||||||
|
xn=bigspec(il-1+indx(nint(0.3*nnl)))
|
||||||
|
bigspec=bigspec/xn
|
||||||
|
ncand=0
|
||||||
|
do i=il,ih
|
||||||
|
if((bigspec(i).gt.bigspec(i-1)).and. &
|
||||||
|
(bigspec(i).gt.bigspec(i+1)).and. &
|
||||||
|
(bigspec(i).gt.1.15).and.ncand.lt.100) then
|
||||||
|
ncand=ncand+1
|
||||||
|
candidates(ncand,1)=df*(i-NH2)
|
||||||
|
candidates(ncand,2)=10*log10(bigspec(i)-1)-28.5
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
return
|
||||||
|
end subroutine getcandidate_ft4s
|
||||||
|
|
||||||
|
subroutine ft4s_downsample(iwave,c)
|
||||||
|
|
||||||
|
! Input: i*2 data in iwave() at sample rate 12000 Hz
|
||||||
|
! Output: Complex data in c(), sampled at 375 Hz
|
||||||
|
|
||||||
|
include 'ft4s_params.f90'
|
||||||
|
parameter (NFFT2=NMAX/32)
|
||||||
|
integer*2 iwave(NMAX)
|
||||||
|
complex c(0:NMAX/32-1)
|
||||||
|
complex c1(0:NFFT2-1)
|
||||||
|
complex cx(0:NMAX/2)
|
||||||
|
real x(NMAX)
|
||||||
|
equivalence (x,cx)
|
||||||
|
|
||||||
|
df=12000.0/NMAX
|
||||||
|
x=iwave
|
||||||
|
call four2a(x,NMAX,1,-1,0) !r2c FFT to freq domain
|
||||||
|
i0=nint(1500.0/df)
|
||||||
|
c1(0)=cx(i0)
|
||||||
|
do i=1,NFFT2/2
|
||||||
|
c1(i)=cx(i0+i)
|
||||||
|
c1(NFFT2-i)=cx(i0-i)
|
||||||
|
enddo
|
||||||
|
c1=c1/NFFT2
|
||||||
|
call four2a(c1,NFFT2,1,1,1) !c2c FFT back to time domain
|
||||||
|
c=c1(0:NMAX/32-1)
|
||||||
|
return
|
||||||
|
end subroutine ft4s_downsample
|
||||||
|
|
||||||
@ -4,7 +4,7 @@ program ft4slowsim
|
|||||||
|
|
||||||
use wavhdr
|
use wavhdr
|
||||||
use packjt77
|
use packjt77
|
||||||
include 'wspr4_params.f90' !Set various constants
|
include 'ft4s_params.f90' !Set various constants
|
||||||
type(hdr) h !Header for .wav file
|
type(hdr) h !Header for .wav file
|
||||||
character arg*12,fname*17
|
character arg*12,fname*17
|
||||||
character msg37*37,msgsent37*37
|
character msg37*37,msgsent37*37
|
||||||
@ -83,8 +83,6 @@ program ft4slowsim
|
|||||||
if(fspread.ne.0.0 .or. delay.ne.0.0) call watterson(c,NMAX,NZ,fs,delay,fspread)
|
if(fspread.ne.0.0 .or. delay.ne.0.0) call watterson(c,NMAX,NZ,fs,delay,fspread)
|
||||||
c=sig*c
|
c=sig*c
|
||||||
wave=real(c)
|
wave=real(c)
|
||||||
peak=maxval(abs(wave))
|
|
||||||
nslots=1
|
|
||||||
if(snrdb.lt.90) then
|
if(snrdb.lt.90) then
|
||||||
do i=1,NMAX !Add gaussian noise at specified SNR
|
do i=1,NMAX !Add gaussian noise at specified SNR
|
||||||
xnoise=gran()
|
xnoise=gran()
|
||||||
|
|||||||
@ -8,13 +8,12 @@ subroutine genft4slow(msg0,ichk,msgsent,msgbits,i4tone)
|
|||||||
! - i4tone array of audio tone values, {0,1,2,3}
|
! - i4tone array of audio tone values, {0,1,2,3}
|
||||||
|
|
||||||
! Frame structure:
|
! Frame structure:
|
||||||
! s16 + 87symbols + 2 ramp up/down = 105 total channel symbols
|
! s4 d24 s4 d24 s4 d24 s4 d24 s4 d24 s4
|
||||||
! r1 + s4 + d29 + s4 + d29 + s4 + d29 + s4 + r1
|
|
||||||
|
|
||||||
! Message duration: TxT = 105*13312/12000 = 116.48 s
|
! Message duration: TxT = 144*9600/12000 = 115.2 s
|
||||||
|
|
||||||
use packjt77
|
use packjt77
|
||||||
include 'wspr4_params.f90'
|
include 'ft4s_params.f90'
|
||||||
character*37 msg0
|
character*37 msg0
|
||||||
character*37 message !Message to be generated
|
character*37 message !Message to be generated
|
||||||
character*37 msgsent !Message as it will be received
|
character*37 msgsent !Message as it will be received
|
||||||
@ -23,13 +22,15 @@ subroutine genft4slow(msg0,ichk,msgsent,msgbits,i4tone)
|
|||||||
integer*4 i4tone(NN),itmp(ND)
|
integer*4 i4tone(NN),itmp(ND)
|
||||||
integer*1 codeword(2*ND)
|
integer*1 codeword(2*ND)
|
||||||
integer*1 msgbits(101),rvec(77)
|
integer*1 msgbits(101),rvec(77)
|
||||||
integer icos4a(4),icos4b(4),icos4c(4),icos4d(4)
|
integer icos4a(4),icos4b(4),icos4c(4),icos4d(4),icos4e(4),icos4f(4)
|
||||||
integer ncrc24
|
integer ncrc24
|
||||||
logical unpk77_success
|
logical unpk77_success
|
||||||
data icos4a/0,1,3,2/
|
data icos4a/0,1,3,2/
|
||||||
data icos4b/1,0,2,3/
|
data icos4b/1,0,2,3/
|
||||||
data icos4c/2,3,1,0/
|
data icos4c/2,3,1,0/
|
||||||
data icos4d/3,2,0,1/
|
data icos4d/3,2,0,1/
|
||||||
|
data icos4e/0,2,3,1/
|
||||||
|
data icos4f/1,2,0,3/
|
||||||
data rvec/0,1,0,0,1,0,1,0,0,1,0,1,1,1,1,0,1,0,0,0,1,0,0,1,1,0,1,1,0, &
|
data rvec/0,1,0,0,1,0,1,0,0,1,0,1,1,1,1,0,1,0,0,0,1,0,0,1,1,0,1,1,0, &
|
||||||
1,0,0,1,0,1,1,0,0,0,0,1,0,0,0,1,0,1,0,0,1,1,1,1,0,0,1,0,1, &
|
1,0,0,1,0,1,1,0,0,0,0,1,0,0,0,1,0,1,0,0,1,1,1,1,0,0,1,0,1, &
|
||||||
0,1,0,1,0,1,1,0,1,1,1,1,1,0,0,0,1,0,1/
|
0,1,0,1,0,1,1,0,1,1,1,1,1,0,0,0,1,0,1/
|
||||||
@ -65,7 +66,7 @@ subroutine genft4slow(msg0,ichk,msgsent,msgbits,i4tone)
|
|||||||
|
|
||||||
entry get_ft4slow_tones_from_101bits(msgbits,i4tone)
|
entry get_ft4slow_tones_from_101bits(msgbits,i4tone)
|
||||||
|
|
||||||
2 call encode174_101(msgbits,codeword)
|
2 call encode240_101(msgbits,codeword)
|
||||||
|
|
||||||
! Grayscale mapping:
|
! Grayscale mapping:
|
||||||
! bits tone
|
! bits tone
|
||||||
@ -82,12 +83,16 @@ entry get_ft4slow_tones_from_101bits(msgbits,i4tone)
|
|||||||
enddo
|
enddo
|
||||||
|
|
||||||
i4tone(1:4)=icos4a
|
i4tone(1:4)=icos4a
|
||||||
i4tone(5:33)=itmp(1:29)
|
i4tone(5:28)=itmp(1:24)
|
||||||
i4tone(34:37)=icos4b
|
i4tone(29:32)=icos4b
|
||||||
i4tone(38:66)=itmp(30:58)
|
i4tone(33:56)=itmp(25:48)
|
||||||
i4tone(67:70)=icos4c
|
i4tone(57:60)=icos4c
|
||||||
i4tone(71:99)=itmp(59:87)
|
i4tone(61:84)=itmp(49:72)
|
||||||
i4tone(100:103)=icos4d
|
i4tone(85:88)=icos4d
|
||||||
|
i4tone(89:112)=itmp(73:96)
|
||||||
|
i4tone(113:116)=icos4e
|
||||||
|
i4tone(117:140)=itmp(97:120)
|
||||||
|
i4tone(141:144)=icos4f
|
||||||
|
|
||||||
999 return
|
999 return
|
||||||
end subroutine genft4slow
|
end subroutine genft4slow
|
||||||
|
|||||||
145
lib/fsk4hf/get_ft4s_bitmetrics.f90
Normal file
145
lib/fsk4hf/get_ft4s_bitmetrics.f90
Normal file
@ -0,0 +1,145 @@
|
|||||||
|
subroutine get_ft4s_bitmetrics(cd,bitmetrics,badsync)
|
||||||
|
|
||||||
|
include 'ft4s_params.f90'
|
||||||
|
parameter (NSS=20)
|
||||||
|
complex cd(0:NN*NSS-1)
|
||||||
|
complex cs(0:3,NN)
|
||||||
|
complex csymb(NSS)
|
||||||
|
integer icos4a(0:3),icos4b(0:3)
|
||||||
|
integer icos4c(0:3),icos4d(0:3)
|
||||||
|
integer icos4e(0:3),icos4f(0:3)
|
||||||
|
integer graymap(0:3)
|
||||||
|
integer ip(1)
|
||||||
|
logical one(0:65535,0:15) ! 65536 8-symbol sequences, 16 bits
|
||||||
|
logical first
|
||||||
|
logical badsync
|
||||||
|
real bitmetrics(2*NN,4)
|
||||||
|
real s2(0:65535)
|
||||||
|
real s4(0:3,NN)
|
||||||
|
|
||||||
|
data icos4a/0,1,3,2/
|
||||||
|
data icos4b/1,0,2,3/
|
||||||
|
data icos4c/2,3,1,0/
|
||||||
|
data icos4d/3,2,0,1/
|
||||||
|
data icos4e/0,2,3,1/
|
||||||
|
data icos4f/1,2,0,3/
|
||||||
|
data graymap/0,1,3,2/
|
||||||
|
data first/.true./
|
||||||
|
save first,one
|
||||||
|
|
||||||
|
if(first) then
|
||||||
|
one=.false.
|
||||||
|
do i=0,65535
|
||||||
|
do j=0,15
|
||||||
|
if(iand(i,2**j).ne.0) one(i,j)=.true.
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
first=.false.
|
||||||
|
endif
|
||||||
|
|
||||||
|
do k=1,NN
|
||||||
|
i1=(k-1)*NSS
|
||||||
|
csymb=cd(i1:i1+NSS-1)
|
||||||
|
call four2a(csymb,NSS,1,-1,1)
|
||||||
|
cs(0:3,k)=csymb(1:4)
|
||||||
|
s4(0:3,k)=abs(csymb(1:4))
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! Sync quality check
|
||||||
|
is1=0
|
||||||
|
is2=0
|
||||||
|
is3=0
|
||||||
|
is4=0
|
||||||
|
is5=0
|
||||||
|
is6=0
|
||||||
|
badsync=.false.
|
||||||
|
ibmax=0
|
||||||
|
|
||||||
|
do k=1,4
|
||||||
|
ip=maxloc(s4(:,k))
|
||||||
|
if(icos4a(k-1).eq.(ip(1)-1)) is1=is1+1
|
||||||
|
ip=maxloc(s4(:,k+28))
|
||||||
|
if(icos4b(k-1).eq.(ip(1)-1)) is2=is2+1
|
||||||
|
ip=maxloc(s4(:,k+56))
|
||||||
|
if(icos4c(k-1).eq.(ip(1)-1)) is3=is3+1
|
||||||
|
ip=maxloc(s4(:,k+84))
|
||||||
|
if(icos4d(k-1).eq.(ip(1)-1)) is4=is4+1
|
||||||
|
ip=maxloc(s4(:,k+112))
|
||||||
|
if(icos4e(k-1).eq.(ip(1)-1)) is5=is5+1
|
||||||
|
ip=maxloc(s4(:,k+140))
|
||||||
|
if(icos4f(k-1).eq.(ip(1)-1)) is6=is6+1
|
||||||
|
enddo
|
||||||
|
nsync=is1+is2+is3+is4+is5+is6 !Number of correct hard sync symbols, 0-24
|
||||||
|
|
||||||
|
badsync=.false.
|
||||||
|
! if(nsync .lt. 8) then
|
||||||
|
! badsync=.true.
|
||||||
|
! return
|
||||||
|
! endif
|
||||||
|
|
||||||
|
do nseq=1,4 !Try coherent sequences of 1, 2, and 4 symbols
|
||||||
|
if(nseq.eq.1) nsym=1
|
||||||
|
if(nseq.eq.2) nsym=2
|
||||||
|
if(nseq.eq.3) nsym=4
|
||||||
|
if(nseq.eq.4) nsym=8
|
||||||
|
nt=2**(2*nsym)
|
||||||
|
do ks=1,NN-nsym+1,nsym !87+16=103 symbols.
|
||||||
|
amax=-1.0
|
||||||
|
do i=0,nt-1
|
||||||
|
! i1=i/64
|
||||||
|
! i2=iand(i,63)/16
|
||||||
|
! i3=iand(i,15)/4
|
||||||
|
! i4=iand(i,3)
|
||||||
|
i1=i/16384
|
||||||
|
i2=iand(i,16383)/4096
|
||||||
|
i3=iand(i,4095)/1024
|
||||||
|
i4=iand(i,1023)/256
|
||||||
|
i5=iand(i,255)/64
|
||||||
|
i6=iand(i,63)/16
|
||||||
|
i7=iand(i,15)/4
|
||||||
|
i8=iand(i,3)
|
||||||
|
if(nsym.eq.1) then
|
||||||
|
s2(i)=abs(cs(graymap(i8),ks))
|
||||||
|
elseif(nsym.eq.2) then
|
||||||
|
s2(i)=abs(cs(graymap(i7),ks)+cs(graymap(i8),ks+1))
|
||||||
|
elseif(nsym.eq.4) then
|
||||||
|
s2(i)=abs(cs(graymap(i5),ks ) + &
|
||||||
|
cs(graymap(i6),ks+1) + &
|
||||||
|
cs(graymap(i7),ks+2) + &
|
||||||
|
cs(graymap(i8),ks+3) &
|
||||||
|
)
|
||||||
|
elseif(nsym.eq.8) then
|
||||||
|
s2(i)=abs(cs(graymap(i1),ks ) + &
|
||||||
|
cs(graymap(i2),ks+1) + &
|
||||||
|
cs(graymap(i3),ks+2) + &
|
||||||
|
cs(graymap(i4),ks+3) + &
|
||||||
|
cs(graymap(i5),ks+4) + &
|
||||||
|
cs(graymap(i6),ks+5) + &
|
||||||
|
cs(graymap(i7),ks+6) + &
|
||||||
|
cs(graymap(i8),ks+7) &
|
||||||
|
)
|
||||||
|
else
|
||||||
|
print*,"Error - nsym must be 1, 2, 4, or 8."
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
ipt=1+(ks-1)*2
|
||||||
|
if(nsym.eq.1) ibmax=1
|
||||||
|
if(nsym.eq.2) ibmax=3
|
||||||
|
if(nsym.eq.4) ibmax=7
|
||||||
|
if(nsym.eq.8) ibmax=15
|
||||||
|
do ib=0,ibmax
|
||||||
|
bm=maxval(s2(0:nt-1),one(0:nt-1,ibmax-ib)) - &
|
||||||
|
maxval(s2(0:nt-1),.not.one(0:nt-1,ibmax-ib))
|
||||||
|
if(ipt+ib.gt.2*NN) cycle
|
||||||
|
bitmetrics(ipt+ib,nseq)=bm
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
call normalizebmet(bitmetrics(:,1),2*NN)
|
||||||
|
call normalizebmet(bitmetrics(:,2),2*NN)
|
||||||
|
call normalizebmet(bitmetrics(:,3),2*NN)
|
||||||
|
call normalizebmet(bitmetrics(:,4),2*NN)
|
||||||
|
return
|
||||||
|
|
||||||
|
end subroutine get_ft4s_bitmetrics
|
||||||
@ -105,6 +105,7 @@ do idb = 20,-10,-1
|
|||||||
! max_iterations is max number of belief propagation iterations
|
! max_iterations is max number of belief propagation iterations
|
||||||
call bpdecode174_91(llr, apmask, max_iterations, message77, cw, nharderrors,niterations,ncheck)
|
call bpdecode174_91(llr, apmask, max_iterations, message77, cw, nharderrors,niterations,ncheck)
|
||||||
if( ndepth .ge. 0 .and. nharderrors .lt. 0 ) call osd174_91(llr, apmask, ndepth, message77, cw, nharderrors, dmin)
|
if( ndepth .ge. 0 .and. nharderrors .lt. 0 ) call osd174_91(llr, apmask, ndepth, message77, cw, nharderrors, dmin)
|
||||||
|
!call decode174_91(llr, apmask, max_iterations, message77, cw, nharderrors,niterations,ncheck,dmin,nsuper)
|
||||||
! If the decoder finds a valid codeword, nharderrors will be .ge. 0.
|
! If the decoder finds a valid codeword, nharderrors will be .ge. 0.
|
||||||
if( nharderrors .ge. 0 ) then
|
if( nharderrors .ge. 0 ) then
|
||||||
call extractmessage77(message77,msgreceived)
|
call extractmessage77(message77,msgreceived)
|
||||||
|
|||||||
144
lib/fsk4hf/ldpcsim240_101.f90
Normal file
144
lib/fsk4hf/ldpcsim240_101.f90
Normal file
@ -0,0 +1,144 @@
|
|||||||
|
program ldpcsim240_101
|
||||||
|
|
||||||
|
! End-to-end test of the (240,101)/crc24 encoder and decoders.
|
||||||
|
|
||||||
|
use packjt77
|
||||||
|
|
||||||
|
parameter(N=240, K=101, M=N-K)
|
||||||
|
character*8 arg
|
||||||
|
character*37 msg0,msg
|
||||||
|
character*77 c77
|
||||||
|
character*24 c24
|
||||||
|
integer*1 msgbits(101)
|
||||||
|
integer*1 apmask(240)
|
||||||
|
integer*1 cw(240)
|
||||||
|
integer*1 codeword(N),message101(101)
|
||||||
|
integer ncrc24
|
||||||
|
real rxdata(N),llr(N)
|
||||||
|
real llrd(240)
|
||||||
|
logical first,unpk77_success
|
||||||
|
data first/.true./
|
||||||
|
|
||||||
|
nargs=iargc()
|
||||||
|
if(nargs.ne.5 .and. nargs.ne.6) then
|
||||||
|
print*,'Usage: ldpcsim niter ndeep #trials s K [msg]'
|
||||||
|
print*,'e.g. ldpcsim240_101 20 5 1000 0.85 91 "K9AN K1JT FN20"'
|
||||||
|
print*,'s : if negative, then value is ignored and sigma is calculated from SNR.'
|
||||||
|
print*,'niter: is the number of BP iterations.'
|
||||||
|
print*,'ndeep: -1 is BP only, ndeep>=0 is OSD order'
|
||||||
|
print*,'K :is the number of message+CRC bits and must be in the range [77,101]'
|
||||||
|
print*,'WSPR-format message is optional'
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
call getarg(1,arg)
|
||||||
|
read(arg,*) max_iterations
|
||||||
|
call getarg(2,arg)
|
||||||
|
read(arg,*) ndeep
|
||||||
|
call getarg(3,arg)
|
||||||
|
read(arg,*) ntrials
|
||||||
|
call getarg(4,arg)
|
||||||
|
read(arg,*) s
|
||||||
|
call getarg(5,arg)
|
||||||
|
read(arg,*) Keff
|
||||||
|
msg0='K9AN K1JT FN20 '
|
||||||
|
if(nargs.eq.6) call getarg(6,msg0)
|
||||||
|
call pack77(msg0,i3,n3,c77)
|
||||||
|
|
||||||
|
rate=real(Keff)/real(N)
|
||||||
|
|
||||||
|
write(*,*) "code rate: ",rate
|
||||||
|
write(*,*) "niter : ",max_iterations
|
||||||
|
write(*,*) "ndeep : ",ndeep
|
||||||
|
write(*,*) "s : ",s
|
||||||
|
write(*,*) "K : ",Keff
|
||||||
|
|
||||||
|
msgbits=0
|
||||||
|
read(c77,'(77i1)') msgbits(1:77)
|
||||||
|
write(*,*) 'message'
|
||||||
|
write(*,'(77i1)') msgbits(1:77)
|
||||||
|
|
||||||
|
call get_crc24(msgbits,101,ncrc24)
|
||||||
|
write(c24,'(b24.24)') ncrc24
|
||||||
|
read(c24,'(24i1)') msgbits(78:101)
|
||||||
|
write(*,'(24i1)') msgbits(78:101)
|
||||||
|
write(*,*) 'message with crc24'
|
||||||
|
write(*,'(101i1)') msgbits(1:101)
|
||||||
|
call encode240_101(msgbits,codeword)
|
||||||
|
call init_random_seed()
|
||||||
|
call sgran()
|
||||||
|
|
||||||
|
write(*,*) 'codeword'
|
||||||
|
write(*,'(77i1,1x,24i1,1x,73i1)') codeword
|
||||||
|
|
||||||
|
write(*,*) "Eb/N0 Es/N0 ngood nundetected sigma symbol error rate"
|
||||||
|
do idb = 8,-3,-1
|
||||||
|
db=idb/2.0-1.0
|
||||||
|
sigma=1/sqrt( 2*rate*(10**(db/10.0)) ) ! to make db represent Eb/No
|
||||||
|
! sigma=1/sqrt( 2*(10**(db/10.0)) ) ! db represents Es/No
|
||||||
|
ngood=0
|
||||||
|
nue=0
|
||||||
|
nberr=0
|
||||||
|
do itrial=1, ntrials
|
||||||
|
! Create a realization of a noisy received word
|
||||||
|
do i=1,N
|
||||||
|
rxdata(i) = 2.0*codeword(i)-1.0 + sigma*gran()
|
||||||
|
enddo
|
||||||
|
nerr=0
|
||||||
|
do i=1,N
|
||||||
|
if( rxdata(i)*(2*codeword(i)-1.0) .lt. 0 ) nerr=nerr+1
|
||||||
|
enddo
|
||||||
|
nberr=nberr+nerr
|
||||||
|
|
||||||
|
rxav=sum(rxdata)/N
|
||||||
|
rx2av=sum(rxdata*rxdata)/N
|
||||||
|
rxsig=sqrt(rx2av-rxav*rxav)
|
||||||
|
rxdata=rxdata/rxsig
|
||||||
|
if( s .lt. 0 ) then
|
||||||
|
ss=sigma
|
||||||
|
else
|
||||||
|
ss=s
|
||||||
|
endif
|
||||||
|
|
||||||
|
llr=2.0*rxdata/(ss*ss)
|
||||||
|
apmask=0
|
||||||
|
! max_iterations is max number of belief propagation iterations
|
||||||
|
call bpdecode240_101(llr,apmask,max_iterations,message101,cw,nharderror,niterations,nchecks)
|
||||||
|
dmin=0.0
|
||||||
|
if( (nharderror .lt. 0) .and. (ndeep .ge. 0) ) then
|
||||||
|
! call osd240_101(llr, Keff, apmask, ndeep, message101, cw, nharderror, dmin)
|
||||||
|
maxsuper=2
|
||||||
|
call decode240_101(llr, Keff, ndeep, apmask, maxsuper, message101, cw, nharderror, iterations, ncheck, dmin, isuper)
|
||||||
|
endif
|
||||||
|
|
||||||
|
if(nharderror.ge.0) then
|
||||||
|
n2err=0
|
||||||
|
do i=1,N
|
||||||
|
if( cw(i)*(2*codeword(i)-1.0) .lt. 0 ) n2err=n2err+1
|
||||||
|
enddo
|
||||||
|
if(n2err.eq.0) then
|
||||||
|
ngood=ngood+1
|
||||||
|
else
|
||||||
|
nue=nue+1
|
||||||
|
endif
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
! snr2500=db+10*log10(200.0/116.0/2500.0)
|
||||||
|
esn0=db+10*log10(rate)
|
||||||
|
pberr=real(nberr)/(real(ntrials*N))
|
||||||
|
write(*,"(f4.1,4x,f5.1,1x,i8,1x,i8,8x,f5.2,8x,e10.3)") db,esn0,ngood,nue,ss,pberr
|
||||||
|
|
||||||
|
if(first) then
|
||||||
|
write(c77,'(77i1)') message101(1:77)
|
||||||
|
write(*,'(101i1)') message101
|
||||||
|
call unpack77(c77,0,msg,unpk77_success)
|
||||||
|
if(unpk77_success) then
|
||||||
|
write(*,1100) msg(1:37)
|
||||||
|
1100 format('Decoded message: ',a37)
|
||||||
|
else
|
||||||
|
print*,'Error unpacking message'
|
||||||
|
endif
|
||||||
|
first=.false.
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
|
||||||
|
end program ldpcsim240_101
|
||||||
403
lib/fsk4hf/osd240_101.f90
Normal file
403
lib/fsk4hf/osd240_101.f90
Normal file
@ -0,0 +1,403 @@
|
|||||||
|
subroutine osd240_101(llr,k,apmask,ndeep,message101,cw,nhardmin,dmin)
|
||||||
|
!
|
||||||
|
! An ordered-statistics decoder for the (240,101) code.
|
||||||
|
! Message payload is 77 bits. Any or all of a 24-bit CRC can be
|
||||||
|
! used for detecting incorrect codewords. The remaining CRC bits are
|
||||||
|
! cascaded with the LDPC code for the purpose of improving the
|
||||||
|
! distance spectrum of the code.
|
||||||
|
!
|
||||||
|
! If p1 (0.le.p1.le.24) is the number of CRC24 bits that are
|
||||||
|
! to be used for bad codeword detection, then the argument k should
|
||||||
|
! be set to 77+p1.
|
||||||
|
!
|
||||||
|
! Valid values for k are in the range [77,101].
|
||||||
|
!
|
||||||
|
character*24 c24
|
||||||
|
integer, parameter:: N=240
|
||||||
|
integer*1 apmask(N),apmaskr(N)
|
||||||
|
integer*1, allocatable, save :: gen(:,:)
|
||||||
|
integer*1, allocatable :: genmrb(:,:),g2(:,:)
|
||||||
|
integer*1, allocatable :: temp(:),m0(:),me(:),mi(:),misub(:),e2sub(:),e2(:),ui(:)
|
||||||
|
integer*1, allocatable :: r2pat(:)
|
||||||
|
integer indices(N),nxor(N)
|
||||||
|
integer*1 cw(N),ce(N),c0(N),hdec(N)
|
||||||
|
integer*1, allocatable :: decoded(:)
|
||||||
|
integer*1 message101(101)
|
||||||
|
integer indx(N)
|
||||||
|
real llr(N),rx(N),absrx(N)
|
||||||
|
|
||||||
|
logical first,reset
|
||||||
|
data first/.true./
|
||||||
|
save first
|
||||||
|
|
||||||
|
allocate( genmrb(k,N), g2(N,k) )
|
||||||
|
allocate( temp(k), m0(k), me(k), mi(k), misub(k), e2sub(N-k), e2(N-k), ui(N-k) )
|
||||||
|
allocate( r2pat(N-k), decoded(k) )
|
||||||
|
|
||||||
|
if( first ) then ! fill the generator matrix
|
||||||
|
!
|
||||||
|
! Create generator matrix for partial CRC cascaded with LDPC code.
|
||||||
|
!
|
||||||
|
! Let p2=101-k and p1+p2=24.
|
||||||
|
!
|
||||||
|
! The last p2 bits of the CRC24 are cascaded with the LDPC code.
|
||||||
|
!
|
||||||
|
! The first p1=k-77 CRC24 bits will be used for error detection.
|
||||||
|
!
|
||||||
|
allocate( gen(k,N) )
|
||||||
|
gen=0
|
||||||
|
do i=1,k
|
||||||
|
message101=0
|
||||||
|
message101(i)=1
|
||||||
|
if(i.le.77) then
|
||||||
|
call get_crc24(message101,101,ncrc24)
|
||||||
|
write(c24,'(b24.24)') ncrc24
|
||||||
|
read(c24,'(24i1)') message101(78:101)
|
||||||
|
message101(78:k)=0
|
||||||
|
endif
|
||||||
|
call encode240_101(message101,cw)
|
||||||
|
gen(i,:)=cw
|
||||||
|
enddo
|
||||||
|
|
||||||
|
first=.false.
|
||||||
|
endif
|
||||||
|
|
||||||
|
rx=llr
|
||||||
|
apmaskr=apmask
|
||||||
|
|
||||||
|
! Hard decisions on the received word.
|
||||||
|
hdec=0
|
||||||
|
where(rx .ge. 0) hdec=1
|
||||||
|
|
||||||
|
! Use magnitude of received symbols as a measure of reliability.
|
||||||
|
absrx=abs(rx)
|
||||||
|
call indexx(absrx,N,indx)
|
||||||
|
|
||||||
|
! Re-order the columns of the generator matrix in order of decreasing reliability.
|
||||||
|
do i=1,N
|
||||||
|
genmrb(1:k,i)=gen(1:k,indx(N+1-i))
|
||||||
|
indices(i)=indx(N+1-i)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! Do gaussian elimination to create a generator matrix with the most reliable
|
||||||
|
! received bits in positions 1:k in order of decreasing reliability (more or less).
|
||||||
|
do id=1,k ! diagonal element indices
|
||||||
|
do icol=id,k+20 ! The 20 is ad hoc - beware
|
||||||
|
iflag=0
|
||||||
|
if( genmrb(id,icol) .eq. 1 ) then
|
||||||
|
iflag=1
|
||||||
|
if( icol .ne. id ) then ! reorder column
|
||||||
|
temp(1:k)=genmrb(1:k,id)
|
||||||
|
genmrb(1:k,id)=genmrb(1:k,icol)
|
||||||
|
genmrb(1:k,icol)=temp(1:k)
|
||||||
|
itmp=indices(id)
|
||||||
|
indices(id)=indices(icol)
|
||||||
|
indices(icol)=itmp
|
||||||
|
endif
|
||||||
|
do ii=1,k
|
||||||
|
if( ii .ne. id .and. genmrb(ii,id) .eq. 1 ) then
|
||||||
|
genmrb(ii,1:N)=ieor(genmrb(ii,1:N),genmrb(id,1:N))
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
exit
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
g2=transpose(genmrb)
|
||||||
|
|
||||||
|
! The hard decisions for the k MRB bits define the order 0 message, m0.
|
||||||
|
! Encode m0 using the modified generator matrix to find the "order 0" codeword.
|
||||||
|
! Flip various combinations of bits in m0 and re-encode to generate a list of
|
||||||
|
! codewords. Return the member of the list that has the smallest Euclidean
|
||||||
|
! distance to the received word.
|
||||||
|
|
||||||
|
hdec=hdec(indices) ! hard decisions from received symbols
|
||||||
|
m0=hdec(1:k) ! zero'th order message
|
||||||
|
absrx=absrx(indices)
|
||||||
|
rx=rx(indices)
|
||||||
|
apmaskr=apmaskr(indices)
|
||||||
|
|
||||||
|
call mrbencode101(m0,c0,g2,N,k)
|
||||||
|
nxor=ieor(c0,hdec)
|
||||||
|
nhardmin=sum(nxor)
|
||||||
|
dmin=sum(nxor*absrx)
|
||||||
|
|
||||||
|
cw=c0
|
||||||
|
ntotal=0
|
||||||
|
nrejected=0
|
||||||
|
npre1=0
|
||||||
|
npre2=0
|
||||||
|
|
||||||
|
if(ndeep.eq.0) goto 998 ! norder=0
|
||||||
|
if(ndeep.gt.6) ndeep=6
|
||||||
|
if( ndeep.eq. 1) then
|
||||||
|
nord=1
|
||||||
|
npre1=0
|
||||||
|
npre2=0
|
||||||
|
nt=40
|
||||||
|
ntheta=12
|
||||||
|
elseif(ndeep.eq.2) then
|
||||||
|
nord=1
|
||||||
|
npre1=1
|
||||||
|
npre2=0
|
||||||
|
nt=40
|
||||||
|
ntheta=12
|
||||||
|
elseif(ndeep.eq.3) then
|
||||||
|
nord=1
|
||||||
|
npre1=1
|
||||||
|
npre2=1
|
||||||
|
nt=40
|
||||||
|
ntheta=12
|
||||||
|
ntau=14
|
||||||
|
elseif(ndeep.eq.4) then
|
||||||
|
nord=2
|
||||||
|
npre1=1
|
||||||
|
npre2=1
|
||||||
|
nt=40
|
||||||
|
ntheta=12
|
||||||
|
ntau=19
|
||||||
|
elseif(ndeep.eq.5) then
|
||||||
|
nord=3
|
||||||
|
npre1=1
|
||||||
|
npre2=1
|
||||||
|
nt=40
|
||||||
|
ntheta=12
|
||||||
|
ntau=19
|
||||||
|
elseif(ndeep.eq.6) then
|
||||||
|
nord=4
|
||||||
|
npre1=1
|
||||||
|
npre2=1
|
||||||
|
nt=40
|
||||||
|
ntheta=12
|
||||||
|
ntau=19
|
||||||
|
endif
|
||||||
|
|
||||||
|
do iorder=1,nord
|
||||||
|
misub(1:k-iorder)=0
|
||||||
|
misub(k-iorder+1:k)=1
|
||||||
|
iflag=k-iorder+1
|
||||||
|
do while(iflag .ge.0)
|
||||||
|
if(iorder.eq.nord .and. npre1.eq.0) then
|
||||||
|
iend=iflag
|
||||||
|
else
|
||||||
|
iend=1
|
||||||
|
endif
|
||||||
|
d1=0.
|
||||||
|
do n1=iflag,iend,-1
|
||||||
|
mi=misub
|
||||||
|
mi(n1)=1
|
||||||
|
if(any(iand(apmaskr(1:k),mi).eq.1)) cycle
|
||||||
|
ntotal=ntotal+1
|
||||||
|
me=ieor(m0,mi)
|
||||||
|
if(n1.eq.iflag) then
|
||||||
|
call mrbencode101(me,ce,g2,N,k)
|
||||||
|
e2sub=ieor(ce(k+1:N),hdec(k+1:N))
|
||||||
|
e2=e2sub
|
||||||
|
nd1kpt=sum(e2sub(1:nt))+1
|
||||||
|
d1=sum(ieor(me(1:k),hdec(1:k))*absrx(1:k))
|
||||||
|
else
|
||||||
|
e2=ieor(e2sub,g2(k+1:N,n1))
|
||||||
|
nd1kpt=sum(e2(1:nt))+2
|
||||||
|
endif
|
||||||
|
if(nd1kpt .le. ntheta) then
|
||||||
|
call mrbencode101(me,ce,g2,N,k)
|
||||||
|
nxor=ieor(ce,hdec)
|
||||||
|
if(n1.eq.iflag) then
|
||||||
|
dd=d1+sum(e2sub*absrx(k+1:N))
|
||||||
|
else
|
||||||
|
dd=d1+ieor(ce(n1),hdec(n1))*absrx(n1)+sum(e2*absrx(k+1:N))
|
||||||
|
endif
|
||||||
|
if( dd .lt. dmin ) then
|
||||||
|
dmin=dd
|
||||||
|
cw=ce
|
||||||
|
nhardmin=sum(nxor)
|
||||||
|
nd1kptbest=nd1kpt
|
||||||
|
endif
|
||||||
|
else
|
||||||
|
nrejected=nrejected+1
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
! Get the next test error pattern, iflag will go negative
|
||||||
|
! when the last pattern with weight iorder has been generated.
|
||||||
|
call nextpat101(misub,k,iorder,iflag)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
if(npre2.eq.1) then
|
||||||
|
reset=.true.
|
||||||
|
ntotal=0
|
||||||
|
do i1=k,1,-1
|
||||||
|
do i2=i1-1,1,-1
|
||||||
|
ntotal=ntotal+1
|
||||||
|
mi(1:ntau)=ieor(g2(k+1:k+ntau,i1),g2(k+1:k+ntau,i2))
|
||||||
|
call boxit101(reset,mi(1:ntau),ntau,ntotal,i1,i2)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
ncount2=0
|
||||||
|
ntotal2=0
|
||||||
|
reset=.true.
|
||||||
|
! Now run through again and do the second pre-processing rule
|
||||||
|
misub(1:k-nord)=0
|
||||||
|
misub(k-nord+1:k)=1
|
||||||
|
iflag=k-nord+1
|
||||||
|
do while(iflag .ge.0)
|
||||||
|
me=ieor(m0,misub)
|
||||||
|
call mrbencode101(me,ce,g2,N,k)
|
||||||
|
e2sub=ieor(ce(k+1:N),hdec(k+1:N))
|
||||||
|
do i2=0,ntau
|
||||||
|
ntotal2=ntotal2+1
|
||||||
|
ui=0
|
||||||
|
if(i2.gt.0) ui(i2)=1
|
||||||
|
r2pat=ieor(e2sub,ui)
|
||||||
|
778 continue
|
||||||
|
call fetchit101(reset,r2pat(1:ntau),ntau,in1,in2)
|
||||||
|
if(in1.gt.0.and.in2.gt.0) then
|
||||||
|
ncount2=ncount2+1
|
||||||
|
mi=misub
|
||||||
|
mi(in1)=1
|
||||||
|
mi(in2)=1
|
||||||
|
if(sum(mi).lt.nord+npre1+npre2.or.any(iand(apmaskr(1:k),mi).eq.1)) cycle
|
||||||
|
me=ieor(m0,mi)
|
||||||
|
call mrbencode101(me,ce,g2,N,k)
|
||||||
|
nxor=ieor(ce,hdec)
|
||||||
|
dd=sum(nxor*absrx)
|
||||||
|
if( dd .lt. dmin ) then
|
||||||
|
dmin=dd
|
||||||
|
cw=ce
|
||||||
|
nhardmin=sum(nxor)
|
||||||
|
endif
|
||||||
|
goto 778
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
call nextpat101(misub,k,nord,iflag)
|
||||||
|
enddo
|
||||||
|
endif
|
||||||
|
|
||||||
|
998 continue
|
||||||
|
! Re-order the codeword to [message bits][parity bits] format.
|
||||||
|
cw(indices)=cw
|
||||||
|
hdec(indices)=hdec
|
||||||
|
message101=cw(1:101)
|
||||||
|
call get_crc24(message101,101,nbadcrc)
|
||||||
|
if(nbadcrc.ne.0) nhardmin=-nhardmin
|
||||||
|
|
||||||
|
return
|
||||||
|
end subroutine osd240_101
|
||||||
|
|
||||||
|
subroutine mrbencode101(me,codeword,g2,N,K)
|
||||||
|
integer*1 me(K),codeword(N),g2(N,K)
|
||||||
|
! fast encoding for low-weight test patterns
|
||||||
|
codeword=0
|
||||||
|
do i=1,K
|
||||||
|
if( me(i) .eq. 1 ) then
|
||||||
|
codeword=ieor(codeword,g2(1:N,i))
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
return
|
||||||
|
end subroutine mrbencode101
|
||||||
|
|
||||||
|
subroutine nextpat101(mi,k,iorder,iflag)
|
||||||
|
integer*1 mi(k),ms(k)
|
||||||
|
! generate the next test error pattern
|
||||||
|
ind=-1
|
||||||
|
do i=1,k-1
|
||||||
|
if( mi(i).eq.0 .and. mi(i+1).eq.1) ind=i
|
||||||
|
enddo
|
||||||
|
if( ind .lt. 0 ) then ! no more patterns of this order
|
||||||
|
iflag=ind
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
ms=0
|
||||||
|
ms(1:ind-1)=mi(1:ind-1)
|
||||||
|
ms(ind)=1
|
||||||
|
ms(ind+1)=0
|
||||||
|
if( ind+1 .lt. k ) then
|
||||||
|
nz=iorder-sum(ms)
|
||||||
|
ms(k-nz+1:k)=1
|
||||||
|
endif
|
||||||
|
mi=ms
|
||||||
|
do i=1,k ! iflag will point to the lowest-index 1 in mi
|
||||||
|
if(mi(i).eq.1) then
|
||||||
|
iflag=i
|
||||||
|
exit
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
return
|
||||||
|
end subroutine nextpat101
|
||||||
|
|
||||||
|
subroutine boxit101(reset,e2,ntau,npindex,i1,i2)
|
||||||
|
integer*1 e2(1:ntau)
|
||||||
|
integer indexes(5000,2),fp(0:525000),np(5000)
|
||||||
|
logical reset
|
||||||
|
common/boxes/indexes,fp,np
|
||||||
|
|
||||||
|
if(reset) then
|
||||||
|
patterns=-1
|
||||||
|
fp=-1
|
||||||
|
np=-1
|
||||||
|
sc=-1
|
||||||
|
indexes=-1
|
||||||
|
reset=.false.
|
||||||
|
endif
|
||||||
|
|
||||||
|
indexes(npindex,1)=i1
|
||||||
|
indexes(npindex,2)=i2
|
||||||
|
ipat=0
|
||||||
|
do i=1,ntau
|
||||||
|
if(e2(i).eq.1) then
|
||||||
|
ipat=ipat+ishft(1,ntau-i)
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
|
||||||
|
ip=fp(ipat) ! see what's currently stored in fp(ipat)
|
||||||
|
if(ip.eq.-1) then
|
||||||
|
fp(ipat)=npindex
|
||||||
|
else
|
||||||
|
do while (np(ip).ne.-1)
|
||||||
|
ip=np(ip)
|
||||||
|
enddo
|
||||||
|
np(ip)=npindex
|
||||||
|
endif
|
||||||
|
return
|
||||||
|
end subroutine boxit101
|
||||||
|
|
||||||
|
subroutine fetchit101(reset,e2,ntau,i1,i2)
|
||||||
|
integer indexes(5000,2),fp(0:525000),np(5000)
|
||||||
|
integer lastpat
|
||||||
|
integer*1 e2(ntau)
|
||||||
|
logical reset
|
||||||
|
common/boxes/indexes,fp,np
|
||||||
|
save lastpat,inext
|
||||||
|
|
||||||
|
if(reset) then
|
||||||
|
lastpat=-1
|
||||||
|
reset=.false.
|
||||||
|
endif
|
||||||
|
|
||||||
|
ipat=0
|
||||||
|
do i=1,ntau
|
||||||
|
if(e2(i).eq.1) then
|
||||||
|
ipat=ipat+ishft(1,ntau-i)
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
index=fp(ipat)
|
||||||
|
|
||||||
|
if(lastpat.ne.ipat .and. index.gt.0) then ! return first set of indices
|
||||||
|
i1=indexes(index,1)
|
||||||
|
i2=indexes(index,2)
|
||||||
|
inext=np(index)
|
||||||
|
elseif(lastpat.eq.ipat .and. inext.gt.0) then
|
||||||
|
i1=indexes(inext,1)
|
||||||
|
i2=indexes(inext,2)
|
||||||
|
inext=np(inext)
|
||||||
|
else
|
||||||
|
i1=-1
|
||||||
|
i2=-1
|
||||||
|
inext=-1
|
||||||
|
endif
|
||||||
|
lastpat=ipat
|
||||||
|
return
|
||||||
|
end subroutine fetchit101
|
||||||
|
|
||||||
Loading…
x
Reference in New Issue
Block a user