w4kek/WSJT-X
1
0
Fork 0
mirror of https://github.com/saitohirga/WSJT-X.git synced 2026-08-14 15:33:36 -04:00
Code Issues Packages Projects Releases Wiki Activity

Make two sync passes, one for each sync array. Label candidates with their type and decode accordingly.

Browse Source
This commit is contained in:
Steve Franke
2018-06-05 16:37:29 -05:00
parent 772b1bd5fa
commit 48b60eab71
16 changed files with 558 additions and 55 deletions
Download Patch File Download Diff File Expand all files Collapse all files
lib/ft8/.sync8.f90.swp
BIN
View File
Binary file not shown.
lib/ft8/bpdecode174_91.f90
+115
View File
@@ -0,0 +1,115 @@
subroutine bpdecode174_91(llr,apmask,maxiterations,decoded,cw,nharderror,iter)
!
! A log-domain belief propagation decoder for the (174,91) code.
!
integer, parameter:: N=174, K=91, M=N-K
integer*1 codeword(N),cw(N),apmask(N)
integer colorder(N)
integer*1 decoded(K)
integer Nm(7,M)
integer Mn(3,N) ! 3 checks per bit
integer synd(M)
real tov(3,N)
real toc(7,M)
real tanhtoc(7,M)
real zn(N)
real llr(N)
real Tmn
integer nrw(M),ncw
include "ldpc_174_91_c_colorder.f90"
include "ldpc_174_91_c_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
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
! write(*,*) 'number of unsatisfied parity checks ',ncheck
if( ncheck .eq. 0 ) then ! we have a codeword - reorder the columns and return it
codeword=cw(colorder+1)
decoded=codeword(M+1:N)
nerr=0
do i=1,N
if( (2*cw(i)-1)*llr(i) .lt. 0.0 ) nerr=nerr+1
enddo
nharderror=nerr
return
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:7,i)=tanh(-toc(1:7,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 bpdecode174_91
lib/ft8/chkcrc14a.f90
+24
View File
@@ -0,0 +1,24 @@
subroutine chkcrc14a(decoded,nbadcrc)
use crc
integer*1 decoded(91)
integer*1, target:: i1Dec8BitBytes(12)
character*91 cbits
! Write decoded bits into cbits: 77-bit message plus 14-bit CRC
write(cbits,1000) decoded
1000 format(91i1)
read(cbits,1001) i1Dec8BitBytes
1001 format(12b8)
read(cbits,1002) ncrc14 !Received CRC14
1002 format(77x,b14)
i1Dec8BitBytes(10)=iand(i1Dec8BitBytes(10),128+64+32+16+8)
i1Dec8BitBytes(11:12)=0
icrc14=crc14(c_loc(i1Dec8BitBytes),12) !CRC14 computed from 77 msg bits
nbadcrc=1
if(ncrc14.eq.icrc14) nbadcrc=0
return
end subroutine chkcrc14a
lib/ft8/encode174_91.f90
+49
View File
@@ -0,0 +1,49 @@
subroutine encode174_91(message,codeword)
! Encode an 91-bit message and return a 174-bit codeword.
! The generator matrix has dimensions (83,91).
! The code is a (174,91) regular ldpc code with column weight 3.
!
integer, parameter:: N=174, K=91, M=N-K
integer*1 codeword(N)
integer*1 gen(M,K)
integer*1 itmp(N)
integer*1 message(K)
integer*1 pchecks(M)
integer colorder(N)
include "ldpc_174_91_c_generator.f90"
include "ldpc_174_91_c_colorder.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,23
read(g(i)(j:j),"(Z1)") istr
ibmax=4
if(j.eq.23) ibmax=3
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
itmp(1:M)=pchecks
itmp(M+1:N)=message(1:K)
codeword(colorder+1)=itmp(1:N)
return
end subroutine encode174_91
lib/ft8/extractmessage174_91.f90
+40
View File
@@ -0,0 +1,40 @@
subroutine extractmessage174_91(decoded,msgreceived,ncrcflag)
use iso_c_binding, only: c_loc,c_size_t
use crc
use packjt
character*22 msgreceived
character*91 cbits
integer*1 decoded(91)
integer*1, target:: i1Dec8BitBytes(12)
integer*4 i4Dec6BitWords(12)
! Write decoded bits into cbits: 77-bit message plus 14-bit CRC
write(cbits,1000) decoded
1000 format(91i1)
read(cbits,1001) i1Dec8BitBytes
1001 format(12b8)
read(cbits,1002) ncrc14 !Received CRC12
1002 format(77x,b14)
i1Dec8BitBytes(10)=iand(i1Dec8BitBytes(10),128+64+32+16+8)
i1Dec8BitBytes(11:12)=0
icrc14=crc14(c_loc(i1Dec8BitBytes),12) !CRC12 computed from 75 msg bits
if(ncrc14.eq.icrc14 .or. sum(decoded(57:87)).eq.0) then !### Kludge ###
! CRC14 checks out --- unpack 72-bit message
do ibyte=1,12
itmp=0
do ibit=1,6
itmp=ishft(itmp,1)+iand(1,decoded((ibyte-1)*6+ibit))
enddo
i4Dec6BitWords(ibyte)=itmp
enddo
call unpackmsg(i4Dec6BitWords,msgreceived,.false.,' ')
ncrcflag=1
else
msgreceived=' '
ncrcflag=-1
endif
return
end subroutine extractmessage174_91
lib/ft8/ft8b.f90 → lib/ft8/ft8b_1.f90
+5 -5
View File
@@ -1,4 +1,4 @@
subroutine ft8b(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, &
subroutine ft8b_1(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, &
napwid,lsubtract,nagain,iaptype,mycall12,mygrid6,hiscall12,bcontest, &
sync0,f1,xdt,xbase,apsym,nharderrors,dmin,nbadcrc,ipass,iera,msg37,xsnr)
@@ -89,7 +89,7 @@ subroutine ft8b(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, &
i0=nint((xdt+0.5)*fs2) !Initial guess for start of signal
smax=0.0
do idt=i0-8,i0+8 !Search over +/- one quarter symbol
call sync8d(cd0,idt,ctwk,0,sync)
call sync8d(cd0,idt,ctwk,0,1,sync)
if(sync.gt.smax) then
smax=sync
ibest=idt
@@ -108,7 +108,7 @@ subroutine ft8b(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, &
ctwk(i)=cmplx(cos(phi),sin(phi))
phi=mod(phi+dphi,twopi)
enddo
call sync8d(cd0,i0,ctwk,1,sync)
call sync8d(cd0,i0,ctwk,1,1,sync)
if( sync .gt. smax ) then
smax=sync
delfbest=delf
@@ -120,7 +120,7 @@ subroutine ft8b(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, &
xdt=xdt2
f1=f1+delfbest !Improved estimate of DF
call sync8d(cd0,i0,ctwk,2,sync)
call sync8d(cd0,i0,ctwk,2,1,sync)
j=0
do k=1,NN
@@ -437,7 +437,7 @@ subroutine ft8b(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, &
enddo
return
end subroutine ft8b
end subroutine ft8b_1
subroutine normalizebmet(bmet,n)
real bmet(n)
lib/ft8/ft8b_2.f90
+457
View File
@@ -0,0 +1,457 @@
subroutine ft8b_2(dd0,newdat,nQSOProgress,nfqso,nftx,ndepth,lapon,lapcqonly, &
napwid,lsubtract,nagain,iaptype,mycall12,mygrid6,hiscall12,bcontest, &
sync0,f1,xdt,xbase,apsym,nharderrors,dmin,nbadcrc,ipass,iera,msg37,xsnr)
use crc
use timer_module, only: timer
include 'ft8_params.f90'
parameter(NP2=2812)
character*37 msg37
character message*22,msgsent*22
character*12 mycall12,hiscall12
character*6 mycall6,mygrid6,hiscall6,c1,c2
character*87 cbits
logical bcontest
real a(5)
real s1(0:7,ND),s2(0:7,NN),s1sort(8*ND)
real ps(0:7),psl(0:7)
real bmeta(3*ND),bmetb(3*ND),bmetap(3*ND)
real llr(3*ND),llra(3*ND),llr0(3*ND),llr1(3*ND),llrap(3*ND) !Soft symbols
real dd0(15*12000)
integer*1 decoded(91),decoded0(91),apmask(3*ND),cw(3*ND)
integer*1 msgbits(91)
integer apsym(91)
integer mcq(28),mde(28),mrrr(16),m73(16),mrr73(16)
integer itone(NN)
integer indxs1(8*ND)
integer icos7(0:6),ip(1)
integer nappasses(0:5) !Number of decoding passes to use for each QSO state
integer naptypes(0:5,4) ! (nQSOProgress, decoding pass) maximum of 4 passes for now
integer*1, target:: i1hiscall(12)
complex cd0(3200)
complex ctwk(32)
complex csymb(32)
logical first,newdat,lsubtract,lapon,lapcqonly,nagain
equivalence (s1,s1sort)
data icos7/3,1,4,0,6,5,2/ ! Flipped w.r.t. original FT8 sync array
data mcq/1,1,1,1,1,0,1,0,0,0,0,0,1,0,0,0,0,0,1,1,0,0,0,1,1,0,0,1/
data mrrr/0,1,1,1,1,1,1,0,1,1,0,0,1,1,1,1/
data m73/0,1,1,1,1,1,1,0,1,1,0,1,0,0,0,0/
data mde/1,1,1,1,1,1,1,1,0,1,1,0,0,1,0,0,0,0,0,1,1,1,0,1,0,0,0,1/
data mrr73/0,0,0,0,0,0,1,0,0,0,0,1,0,1,0,1/
data first/.true./
save nappasses,naptypes
if(first) then
mcq=2*mcq-1
mde=2*mde-1
mrrr=2*mrrr-1
m73=2*m73-1
mrr73=2*mrr73-1
nappasses(0)=2
nappasses(1)=2
nappasses(2)=2
nappasses(3)=4
nappasses(4)=4
nappasses(5)=3
! iaptype
!------------------------
! 1 CQ ??? ???
! 2 MyCall ??? ???
! 3 MyCall DxCall ???
! 4 MyCall DxCall RRR
! 5 MyCall DxCall 73
! 6 MyCall DxCall RR73
! 7 ??? DxCall ???
naptypes(0,1:4)=(/1,2,0,0/)
naptypes(1,1:4)=(/2,3,0,0/)
naptypes(2,1:4)=(/2,3,0,0/)
naptypes(3,1:4)=(/3,4,5,6/)
naptypes(4,1:4)=(/3,4,5,6/)
naptypes(5,1:4)=(/3,1,2,0/)
first=.false.
endif
max_iterations=30
nharderrors=-1
fs2=12000.0/NDOWN
dt2=1.0/fs2
twopi=8.0*atan(1.0)
delfbest=0.
ibest=0
call timer('ft8_down',0)
call ft8_downsample(dd0,newdat,f1,cd0) !Mix f1 to baseband and downsample
call timer('ft8_down',1)
i0=nint((xdt+0.5)*fs2) !Initial guess for start of signal
smax=0.0
do idt=i0-8,i0+8 !Search over +/- one quarter symbol
call sync8d(cd0,idt,ctwk,0,2,sync)
if(sync.gt.smax) then
smax=sync
ibest=idt
endif
enddo
xdt2=ibest*dt2 !Improved estimate for DT
! Now peak up in frequency
i0=nint(xdt2*fs2)
smax=0.0
do ifr=-5,5 !Search over +/- 2.5 Hz
delf=ifr*0.5
dphi=twopi*delf*dt2
phi=0.0
do i=1,32
ctwk(i)=cmplx(cos(phi),sin(phi))
phi=mod(phi+dphi,twopi)
enddo
call sync8d(cd0,i0,ctwk,1,2,sync)
if( sync .gt. smax ) then
smax=sync
delfbest=delf
endif
enddo
a=0.0
a(1)=-delfbest
call twkfreq1(cd0,NP2,fs2,a,cd0)
xdt=xdt2
f1=f1+delfbest !Improved estimate of DF
call sync8d(cd0,i0,ctwk,2,2,sync)
j=0
do k=1,NN
i1=ibest+(k-1)*32
csymb=cmplx(0.0,0.0)
if( i1.ge.1 .and. i1+31 .le. NP2 ) csymb=cd0(i1:i1+31)
call four2a(csymb,32,1,-1,1)
s2(0:7,k)=abs(csymb(1:8))/1e3
enddo
! sync quality check
is1=0
is2=0
is3=0
do k=1,7
ip=maxloc(s2(:,k))
if(icos7(k-1).eq.(ip(1)-1)) is1=is1+1
ip=maxloc(s2(:,k+36))
if(icos7(k-1).eq.(ip(1)-1)) is2=is2+1
ip=maxloc(s2(:,k+72))
if(icos7(k-1).eq.(ip(1)-1)) is3=is3+1
enddo
! hard sync sum - max is 21
nsync=is1+is2+is3
if(nsync .le. 6) then ! bail out
nbadcrc=1
return
endif
j=0
do k=1,NN
if(k.le.7) cycle
if(k.ge.37 .and. k.le.43) cycle
if(k.gt.72) cycle
j=j+1
s1(0:7,j)=s2(0:7,k)
enddo
call indexx(s1sort,8*ND,indxs1)
xmeds1=s1sort(indxs1(nint(0.5*8*ND)))
s1=s1/xmeds1
do j=1,ND
i4=3*j-2
i2=3*j-1
i1=3*j
! Max amplitude
ps=s1(0:7,j)
r1=max(ps(1),ps(3),ps(5),ps(7))-max(ps(0),ps(2),ps(4),ps(6))
r2=max(ps(2),ps(3),ps(6),ps(7))-max(ps(0),ps(1),ps(4),ps(5))
r4=max(ps(4),ps(5),ps(6),ps(7))-max(ps(0),ps(1),ps(2),ps(3))
bmeta(i4)=r4
bmeta(i2)=r2
bmeta(i1)=r1
bmetap(i4)=r4
bmetap(i2)=r2
bmetap(i1)=r1
! Max log metric
psl=log(ps)
r1=max(psl(1),psl(3),psl(5),psl(7))-max(psl(0),psl(2),psl(4),psl(6))
r2=max(psl(2),psl(3),psl(6),psl(7))-max(psl(0),psl(1),psl(4),psl(5))
r4=max(psl(4),psl(5),psl(6),psl(7))-max(psl(0),psl(1),psl(2),psl(3))
bmetb(i4)=r4
bmetb(i2)=r2
bmetb(i1)=r1
! Metric for Cauchy noise
! r1=log(ps(1)**3+ps(3)**3+ps(5)**3+ps(7)**3)- &
! log(ps(0)**3+ps(2)**3+ps(4)**3+ps(6)**3)
! r2=log(ps(2)**3+ps(3)**3+ps(6)**3+ps(7)**3)- &
! log(ps(0)**3+ps(1)**3+ps(4)**3+ps(5)**3)
! r4=log(ps(4)**3+ps(5)**3+ps(6)**3+ps(7)**3)- &
! log(ps(0)**3+ps(1)**3+ps(2)**3+ps(3)**3)
! Metric for AWGN, no fading
! bscale=2.5
! b0=bessi0(bscale*ps(0))
! b1=bessi0(bscale*ps(1))
! b2=bessi0(bscale*ps(2))
! b3=bessi0(bscale*ps(3))
! b4=bessi0(bscale*ps(4))
! b5=bessi0(bscale*ps(5))
! b6=bessi0(bscale*ps(6))
! b7=bessi0(bscale*ps(7))
! r1=log(b1+b3+b5+b7)-log(b0+b2+b4+b6)
! r2=log(b2+b3+b6+b7)-log(b0+b1+b4+b5)
! r4=log(b4+b5+b6+b7)-log(b0+b1+b2+b3)
if(nQSOProgress .eq. 0 .or. nQSOProgress .eq. 5) then
! When bits 88:115 are set as ap bits, bit 115 lives in symbol 39 along
! with no-ap bits 116 and 117. Take care of metrics for bits 116 and 117.
if(j.eq.39) then ! take care of bits that live in symbol 39
if(apsym(28).lt.0) then
bmetap(i2)=max(ps(2),ps(3))-max(ps(0),ps(1))
bmetap(i1)=max(ps(1),ps(3))-max(ps(0),ps(2))
else
bmetap(i2)=max(ps(6),ps(7))-max(ps(4),ps(5))
bmetap(i1)=max(ps(5),ps(7))-max(ps(4),ps(6))
endif
endif
endif
! When bits 116:143 are set as ap bits, bit 115 lives in symbol 39 along
! with ap bits 116 and 117. Take care of metric for bit 115.
! if(j.eq.39) then ! take care of bit 115
! iii=2*(apsym(29)+1)/2 + (apsym(30)+1)/2 ! known values of bits 116 & 117
! if(iii.eq.0) bmetap(i4)=ps(4)-ps(0)
! if(iii.eq.1) bmetap(i4)=ps(5)-ps(1)
! if(iii.eq.2) bmetap(i4)=ps(6)-ps(2)
! if(iii.eq.3) bmetap(i4)=ps(7)-ps(3)
! endif
! bit 144 lives in symbol 48 and will be 1 if it is set as an ap bit.
! take care of metrics for bits 142 and 143
if(j.eq.48) then ! bit 144 is always 1
bmetap(i4)=max(ps(5),ps(7))-max(ps(1),ps(3))
bmetap(i2)=max(ps(3),ps(7))-max(ps(1),ps(5))
endif
! bit 154 lives in symbol 52 and will be 0 if it is set as an ap bit
! take care of metrics for bits 155 and 156
if(j.eq.52) then ! bit 154 will be 0 if it is set as an ap bit.
bmetap(i2)=max(ps(2),ps(3))-max(ps(0),ps(1))
bmetap(i1)=max(ps(1),ps(3))-max(ps(0),ps(2))
endif
enddo
call normalizebmet(bmeta,3*ND)
call normalizebmet(bmetb,3*ND)
call normalizebmet(bmetap,3*ND)
scalefac=2.83
llr0=scalefac*bmeta
llr1=scalefac*bmetb
llra=scalefac*bmetap ! llr's for use with ap
apmag=scalefac*(maxval(abs(bmetap))*1.01)
! pass #
!------------------------------
! 1 regular decoding
! 2 erase 24
! 3 erase 48
! 4 ap pass 1
! 5 ap pass 2
! 6 ap pass 3
! 7 ap pass 4, etc.
if(lapon) then
if(.not.lapcqonly) then
npasses=4+nappasses(nQSOProgress)
else
npasses=5
endif
else
npasses=4
endif
do ipass=1,npasses
llr=llr0
if(ipass.eq.2) llr=llr1
if(ipass.eq.3) llr(1:24)=0.
if(ipass.eq.4) llr(1:48)=0.
if(ipass.le.4) then
apmask=0
llrap=llr
iaptype=0
endif
if(ipass .gt. 4) then
if(.not.lapcqonly) then
iaptype=naptypes(nQSOProgress,ipass-4)
else
iaptype=1
endif
if(iaptype.ge.3 .and. (abs(f1-nfqso).gt.napwid .and. abs(f1-nftx).gt.napwid) ) cycle
if(iaptype.eq.1 .or. iaptype.eq.2 ) then ! AP,???,???
apmask=0
apmask(88:115)=1 ! first 28 bits are AP
apmask(144)=1 ! not free text
llrap=llr
if(iaptype.eq.1) llrap(88:115)=apmag*mcq
if(iaptype.eq.2) llrap(88:115)=apmag*apsym(1:28)
llrap(116:117)=llra(116:117)
llrap(142:143)=llra(142:143)
llrap(144)=-apmag
endif
if(iaptype.eq.3) then ! mycall, dxcall, ???
apmask=0
apmask(88:115)=1 ! mycall
apmask(116:143)=1 ! hiscall
apmask(144)=1 ! not free text
llrap=llr
llrap(88:143)=apmag*apsym(1:56)
llrap(144)=-apmag
endif
if(iaptype.eq.4 .or. iaptype.eq.5 .or. iaptype.eq.6) then
apmask=0
apmask(88:115)=1 ! mycall
apmask(116:143)=1 ! hiscall
apmask(144:159)=1 ! RRR or 73 or RR73
llrap=llr
llrap(88:143)=apmag*apsym(1:56)
if(iaptype.eq.4) llrap(144:159)=apmag*mrrr
if(iaptype.eq.5) llrap(144:159)=apmag*m73
if(iaptype.eq.6) llrap(144:159)=apmag*mrr73
endif
if(iaptype.eq.7) then ! ???, dxcall, ???
apmask=0
apmask(116:143)=1 ! hiscall
apmask(144)=1 ! not free text
llrap=llr
llrap(115)=llra(115)
llrap(116:143)=apmag*apsym(29:56)
llrap(144)=-apmag
endif
endif
cw=0
call timer('bpd174_91 ',0)
call bpdecode174_91(llrap,apmask,max_iterations,decoded,cw,nharderrors, &
niterations)
call timer('bpd174_91 ',1)
dmin=0.0
if(ndepth.eq.3 .and. nharderrors.lt.0) then
ndeep=3
if(abs(nfqso-f1).le.napwid .or. abs(nftx-f1).le.napwid) then
if((ipass.eq.3 .or. ipass.eq.4) .and. .not.nagain) then
ndeep=3
else
ndeep=4
endif
endif
if(nagain) ndeep=5
call timer('osd174_91 ',0)
call osd174_91(llrap,apmask,ndeep,decoded,cw,nharderrors,dmin)
call timer('osd174_91 ',1)
endif
nbadcrc=1
message=' '
xsnr=-99.0
if(count(cw.eq.0).eq.174) cycle !Reject the all-zero codeword
if(nharderrors.ge.0 .and. nharderrors+dmin.lt.60.0 .and. &
.not.(sync.lt.2.0 .and. nharderrors.gt.35) .and. &
.not.(ipass.gt.2 .and. nharderrors.gt.39) .and. &
.not.(ipass.eq.4 .and. nharderrors.gt.30) &
) then
call chkcrc14a(decoded,nbadcrc)
else
nharderrors=-1
cycle
endif
i3bit=4*decoded(73) + 2*decoded(74) + decoded(75)
iFreeText=decoded(57)
if(nbadcrc.eq.0) then
decoded0=decoded
if(i3bit.eq.1) decoded(57:)=0
call extractmessage174_91(decoded,message,ncrcflag)
decoded=decoded0
! This needs fixing for messages with i3bit=1:
call genft8(message,mygrid6,bcontest,i3bit,msgsent,msgbits,itone)
if(lsubtract) call subtractft8(dd0,itone,f1,xdt2)
xsig=0.0
xnoi=0.0
do i=1,79
xsig=xsig+s2(itone(i),i)**2
ios=mod(itone(i)+4,7)
xnoi=xnoi+s2(ios,i)**2
enddo
xsnr=0.001
if(xnoi.gt.0 .and. xnoi.lt.xsig) xsnr=xsig/xnoi-1.0
xsnr=10.0*log10(xsnr)-27.0
xsnr2=db(xsig/xbase - 1.0) - 32.0
if(.not.nagain) xsnr=xsnr2
if(xsnr .lt. -24.0) xsnr=-24.0
if(i3bit.eq.1) then
do i=1,12
i1hiscall(i)=ichar(hiscall12(i:i))
enddo
icrc10=crc10(c_loc(i1hiscall),12)
write(cbits,1001) decoded
1001 format(87i1)
read(cbits,1002) ncrc10,nrpt
1002 format(56x,b10,b6)
irpt=nrpt-30
i1=index(message,' ')
i2=index(message(i1+1:),' ') + i1
c1=message(1:i1)//' '
c2=message(i1+1:i2)//' '
if(ncrc10.eq.icrc10) msg37=c1//' RR73; '//c2//' <'// &
trim(hiscall12)//'> '
if(ncrc10.ne.icrc10) msg37=c1//' RR73; '//c2//' <...> '
! msg37=c1//' RR73; '//c2//' <...> '
write(msg37(35:37),1010) irpt
1010 format(i3.2)
if(msg37(35:35).ne.'-') msg37(35:35)='+'
iz=len(trim(msg37))
do iter=1,10 !Collapse multiple blanks
ib2=index(msg37(1:iz),' ')
if(ib2.lt.1) exit
msg37=msg37(1:ib2)//msg37(ib2+2:)
iz=iz-1
enddo
else
msg37=message//' '
endif
return
endif
enddo
return
end subroutine ft8b_2
!subroutine normalizebmet(bmet,n)
! real bmet(n)
!
! bmetav=sum(bmet)/real(n)
! bmet2av=sum(bmet*bmet)/real(n)
! var=bmet2av-bmetav*bmetav
! if( var .gt. 0.0 ) then
! bmetsig=sqrt(var)
! else
! bmetsig=sqrt(bmet2av)
! endif
! bmet=bmet/bmetsig
! return
!end subroutine normalizebmet
lib/ft8/ldpc_174_91_c_colorder.f90
+11
View File
@@ -0,0 +1,11 @@
data colorder/ &
0, 1, 2, 3, 28, 4, 5, 6, 7, 8, 9, 10, 11, 34, 12, 32, 13, 14, 15, 16, &
17, 18, 36, 29, 43, 19, 20, 42, 21, 40, 30, 37, 22, 47, 61, 45, 44, 23, 41, 39, &
49, 24, 46, 50, 48, 26, 31, 33, 51, 38, 52, 59, 55, 66, 57, 27, 60, 35, 54, 58, &
25, 56, 62, 64, 67, 69, 63, 68, 70, 72, 65, 73, 75, 74, 71, 77, 78, 76, 79, 80, &
53, 81, 83, 82, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, &
100,101,102,103,104,105,106,107,108,109,110,111,112,113,114,115,116,117,118,119, &
120,121,122,123,124,125,126,127,128,129,130,131,132,133,134,135,136,137,138,139, &
140,141,142,143,144,145,146,147,148,149,150,151,152,153,154,155,156,157,158,159, &
160,161,162,163,164,165,166,167,168,169,170,171,172,173/
lib/ft8/ldpc_174_91_c_generator.f90
+86
View File
@@ -0,0 +1,86 @@
character*23 g(83)
data g/ &
"8329ce11bf31eaf509f27fc", &
"761c264e25c259335493132", &
"dc265902fb277c6410a1bdc", &
"1b3f417858cd2dd33ec7f62", &
"09fda4fee04195fd034783a", &
"077cccc11b8873ed5c3d48a", &
"29b62afe3ca036f4fe1a9da", &
"6054faf5f35d96d3b0c8c3e", &
"e20798e4310eed27884ae90", &
"775c9c08e80e26ddae56318", &
"b0b811028c2bf997213487c", &
"18a0c9231fc60adf5c5ea32", &
"76471e8302a0721e01b12b8", &
"ffbccb80ca8341fafb47b2e", &
"66a72a158f9325a2bf67170", &
"c4243689fe85b1c51363a18", &
"0dff739414d1a1b34b1c270", &
"15b48830636c8b99894972e", &
"29a89c0d3de81d665489b0e", &
"4f126f37fa51cbe61bd6b94", &
"99c47239d0d97d3c84e0940", &
"1919b75119765621bb4f1e8", &
"09db12d731faee0b86df6b8", &
"488fc33df43fbdeea4eafb4", &
"827423ee40b675f756eb5fe", &
"abe197c484cb74757144a9a", &
"2b500e4bc0ec5a6d2bdbdd0", &
"c474aa53d70218761669360", &
"8eba1a13db3390bd6718cec", &
"753844673a27782cc42012e", &
"06ff83a145c37035a5c1268", &
"3b37417858cc2dd33ec3f62", &
"9a4a5a28ee17ca9c324842c", &
"bc29f465309c977e89610a4", &
"2663ae6ddf8b5ce2bb29488", &
"46f231efe457034c1814418", &
"3fb2ce85abe9b0c72e06fbe", &
"de87481f282c153971a0a2e", &
"fcd7ccf23c69fa99bba1412", &
"f0261447e9490ca8e474cec", &
"4410115818196f95cdd7012", &
"088fc31df4bfbde2a4eafb4", &
"b8fef1b6307729fb0a078c0", &
"5afea7acccb77bbc9d99a90", &
"49a7016ac653f65ecdc9076", &
"1944d085be4e7da8d6cc7d0", &
"251f62adc4032f0ee714002", &
"56471f8702a0721e00b12b8", &
"2b8e4923f2dd51e2d537fa0", &
"6b550a40a66f4755de95c26", &
"a18ad28d4e27fe92a4f6c84", &
"10c2e586388cb82a3d80758", &
"ef34a41817ee02133db2eb0", &
"7e9c0c54325a9c15836e000", &
"3693e572d1fde4cdf079e86", &
"bfb2cec5abe1b0c72e07fbe", &
"7ee18230c583cccc57d4b08", &
"a066cb2fedafc9f52664126", &
"bb23725abc47cc5f4cc4cd2", &
"ded9dba3bee40c59b5609b4", &
"d9a7016ac653e6decdc9036", &
"9ad46aed5f707f280ab5fc4", &
"e5921c77822587316d7d3c2", &
"4f14da8242a8b86dca73352", &
"8b8b507ad467d4441df770e", &
"22831c9cf1169467ad04b68", &
"213b838fe2ae54c38ee7180", &
"5d926b6dd71f085181a4e12", &
"66ab79d4b29ee6e69509e56", &
"958148682d748a38dd68baa", &
"b8ce020cf069c32a723ab14", &
"f4331d6d461607e95752746", &
"6da23ba424b9596133cf9c8", &
"a636bcbc7b30c5fbeae67fe", &
"5cb0d86a07df654a9089a20", &
"f11f106848780fc9ecdd80a", &
"1fbb5364fb8d2c9d730d5ba", &
"fcb86bc70a50c9d02a5d034", &
"a534433029eac15f322e34c", &
"c989d9c7c3d3b8c55d75130", &
"7bb38b2f0186d46643ae962", &
"2644ebadeb44b9467d1f42c", &
"608cc857594bfbb55d69600"/
lib/ft8/ldpc_174_91_c_parity.f90
+269
View File
@@ -0,0 +1,269 @@
data Mn/ &
1, 24, 66, &
2, 5, 70, &
3, 31, 65, &
4, 49, 58, &
6, 60, 67, &
7, 32, 75, &
8, 48, 82, &
9, 35, 41, &
10, 39, 62, &
11, 14, 61, &
12, 71, 74, &
13, 23, 78, &
15, 16, 79, &
17, 54, 63, &
18, 50, 57, &
19, 30, 47, &
20, 64, 80, &
21, 28, 69, &
22, 25, 43, &
26, 34, 72, &
27, 36, 37, &
29, 40, 44, &
33, 52, 53, &
38, 55, 83, &
42, 51, 59, &
45, 76, 81, &
46, 68, 77, &
56, 67, 73, &
1, 4, 5, &
2, 47, 51, &
3, 46, 82, &
6, 24, 76, &
7, 9, 16, &
8, 10, 78, &
11, 35, 55, &
12, 38, 64, &
13, 22, 37, &
14, 15, 58, &
17, 50, 53, &
18, 29, 56, &
19, 26, 57, &
20, 39, 77, &
21, 36, 63, &
23, 54, 74, &
25, 69, 73, &
27, 33, 65, &
28, 44, 79, &
30, 43, 52, &
31, 45, 61, &
32, 48, 71, &
34, 60, 62, &
40, 41, 70, &
42, 66, 68, &
49, 75, 83, &
59, 72, 80, &
36, 64, 81, &
1, 65, 74, &
2, 8, 81, &
3, 26, 79, &
4, 22, 72, &
5, 38, 50, &
6, 9, 52, &
7, 18, 77, &
10, 28, 55, &
11, 56, 59, &
12, 67, 68, &
13, 29, 47, &
14, 39, 54, &
15, 60, 70, &
16, 37, 66, &
17, 25, 82, &
19, 21, 71, &
20, 30, 31, &
23, 75, 80, &
24, 69, 75, &
27, 32, 62, &
33, 59, 76, &
34, 53, 61, &
35, 46, 47, &
40, 43, 83, &
41, 42, 63, &
20, 44, 48, &
16, 45, 73, &
42, 49, 57, &
25, 51, 62, &
33, 58, 78, &
1, 44, 45, &
2, 7, 61, &
3, 6, 54, &
4, 35, 48, &
5, 13, 21, &
8, 56, 79, &
9, 64, 69, &
10, 19, 66, &
11, 36, 60, &
12, 37, 58, &
14, 32, 43, &
15, 63, 80, &
17, 28, 77, &
18, 74, 83, &
22, 53, 81, &
23, 30, 34, &
24, 31, 40, &
26, 41, 76, &
27, 57, 70, &
29, 49, 65, &
3, 38, 78, &
5, 39, 82, &
46, 50, 73, &
51, 52, 74, &
55, 71, 72, &
44, 67, 72, &
43, 68, 78, &
1, 32, 59, &
2, 6, 71, &
4, 16, 54, &
7, 65, 67, &
8, 30, 42, &
9, 22, 31, &
10, 18, 76, &
11, 23, 82, &
12, 28, 61, &
13, 52, 79, &
14, 50, 51, &
15, 81, 83, &
17, 29, 60, &
19, 33, 64, &
20, 26, 73, &
21, 34, 40, &
24, 27, 77, &
25, 55, 58, &
35, 53, 66, &
36, 48, 68, &
37, 46, 75, &
38, 45, 47, &
39, 57, 69, &
41, 56, 62, &
20, 49, 53, &
46, 52, 63, &
45, 70, 75, &
27, 35, 80, &
1, 15, 30, &
2, 68, 80, &
3, 36, 51, &
4, 28, 51, &
5, 31, 56, &
6, 20, 37, &
7, 40, 82, &
8, 60, 69, &
9, 10, 49, &
11, 44, 57, &
12, 39, 59, &
13, 24, 55, &
14, 21, 65, &
16, 71, 78, &
17, 30, 76, &
18, 25, 80, &
19, 61, 83, &
22, 38, 77, &
23, 41, 50, &
7, 26, 58, &
29, 32, 81, &
33, 40, 73, &
18, 34, 48, &
13, 42, 64, &
5, 26, 43, &
47, 69, 72, &
54, 55, 70, &
45, 62, 68, &
10, 63, 67, &
14, 66, 72, &
22, 60, 74, &
35, 39, 79, &
1, 46, 64/
data Nm/ &
1, 29, 57, 87, 114, 142, 174, &
2, 30, 58, 88, 115, 143, 0, &
3, 31, 59, 89, 107, 144, 0, &
4, 29, 60, 90, 116, 145, 0, &
2, 29, 61, 91, 108, 146, 166, &
5, 32, 62, 89, 115, 147, 0, &
6, 33, 63, 88, 117, 148, 161, &
7, 34, 58, 92, 118, 149, 0, &
8, 33, 62, 93, 119, 150, 0, &
9, 34, 64, 94, 120, 150, 170, &
10, 35, 65, 95, 121, 151, 0, &
11, 36, 66, 96, 122, 152, 0, &
12, 37, 67, 91, 123, 153, 165, &
10, 38, 68, 97, 124, 154, 171, &
13, 38, 69, 98, 125, 142, 0, &
13, 33, 70, 83, 116, 155, 0, &
14, 39, 71, 99, 126, 156, 0, &
15, 40, 63, 100, 120, 157, 164, &
16, 41, 72, 94, 127, 158, 0, &
17, 42, 73, 82, 128, 138, 147, &
18, 43, 72, 91, 129, 154, 0, &
19, 37, 60, 101, 119, 159, 172, &
12, 44, 74, 102, 121, 160, 0, &
1, 32, 75, 103, 130, 153, 0, &
19, 45, 71, 85, 131, 157, 0, &
20, 41, 59, 104, 128, 161, 166, &
21, 46, 76, 105, 130, 141, 0, &
18, 47, 64, 99, 122, 145, 0, &
22, 40, 67, 106, 126, 162, 0, &
16, 48, 73, 102, 118, 142, 156, &
3, 49, 73, 103, 119, 146, 0, &
6, 50, 76, 97, 114, 162, 0, &
23, 46, 77, 86, 127, 163, 0, &
20, 51, 78, 102, 129, 164, 0, &
8, 35, 79, 90, 132, 141, 173, &
21, 43, 56, 95, 133, 144, 0, &
21, 37, 70, 96, 134, 147, 0, &
24, 36, 61, 107, 135, 159, 0, &
9, 42, 68, 108, 136, 152, 173, &
22, 52, 80, 103, 129, 148, 163, &
8, 52, 81, 104, 137, 160, 0, &
25, 53, 81, 84, 118, 165, 0, &
19, 48, 80, 97, 113, 166, 0, &
22, 47, 82, 87, 112, 151, 0, &
26, 49, 83, 87, 135, 140, 169, &
27, 31, 79, 109, 134, 139, 174, &
16, 30, 67, 79, 135, 167, 0, &
7, 50, 82, 90, 133, 164, 0, &
4, 54, 84, 106, 138, 150, 0, &
15, 39, 61, 109, 124, 160, 0, &
25, 30, 85, 110, 124, 144, 145, &
23, 48, 62, 110, 123, 139, 0, &
23, 39, 78, 101, 132, 138, 0, &
14, 44, 68, 89, 116, 168, 0, &
24, 35, 64, 111, 131, 153, 168, &
28, 40, 65, 92, 137, 146, 0, &
15, 41, 84, 105, 136, 151, 0, &
4, 38, 86, 96, 131, 161, 0, &
25, 55, 65, 77, 114, 152, 0, &
5, 51, 69, 95, 126, 149, 172, &
10, 49, 78, 88, 122, 158, 0, &
9, 51, 76, 85, 137, 169, 0, &
14, 43, 81, 98, 139, 170, 0, &
17, 36, 56, 93, 127, 165, 174, &
3, 46, 57, 106, 117, 154, 0, &
1, 53, 70, 94, 132, 171, 0, &
5, 28, 66, 112, 117, 170, 0, &
27, 53, 66, 113, 133, 143, 169, &
18, 45, 75, 93, 136, 149, 167, &
2, 52, 69, 105, 140, 168, 0, &
11, 50, 72, 111, 115, 155, 0, &
20, 55, 60, 111, 112, 167, 171, &
28, 45, 83, 109, 128, 163, 0, &
11, 44, 57, 100, 110, 172, 0, &
6, 54, 74, 75, 134, 140, 0, &
26, 32, 77, 104, 120, 156, 0, &
27, 42, 63, 99, 130, 159, 0, &
12, 34, 86, 107, 113, 155, 0, &
13, 47, 59, 92, 123, 173, 0, &
17, 55, 74, 98, 141, 143, 157, &
26, 56, 58, 101, 125, 162, 0, &
7, 31, 71, 108, 121, 148, 0, &
24, 54, 80, 100, 125, 158, 0/
data nrw/ &
7,6,6,6,7,6,7,6,6,7,6,6,7,7,6,6,6, &
7,6,7,6,7,6,6,6,7,6,6,6,7,6,6,6,6, &
7,6,6,6,7,7,6,6,6,6,7,7,6,6,6,6,7, &
6,6,6,7,6,6,6,6,7,6,6,6,7,6,6,6,7, &
7,6,6,7,6,6,6,6,6,6,6,7,6,6,6/
ncw=3
lib/ft8/ldpcsim174_91.f90
+202
View File
@@ -0,0 +1,202 @@
program ldpcsim174_91
! End to end test of the (174,91)/crc14 encoder and decoder.
use crc
use packjt
integer, parameter:: N=174, K=91, M=N-K
character*22 msg,msgsent,msgreceived
character*8 arg
character*6 grid
integer*1, allocatable :: codeword(:), decoded(:), message(:)
integer*1, target:: i1Msg8BitBytes(12)
integer*1 msgbits(K)
integer*1 apmask(N), cw(N)
integer*2 checksum
integer*4 i4Msg6BitWords(13)
integer colorder(N)
integer nerrtot(N),nerrdec(N),nmpcbad(K)
logical checksumok
real*8, allocatable :: rxdata(:)
real, allocatable :: llr(:)
include "ldpc_174_91_c_colorder.f90"
nerrtot=0
nerrdec=0
nmpcbad=0 ! Used to collect the number of errors in the message+crc part of the codeword
nargs=iargc()
if(nargs.ne.4) then
print*,'Usage: ldpcsim niter ndepth #trials s '
print*,'eg: ldpcsim 10 2 1000 0.84'
print*,'belief propagation iterations: niter, ordered-statistics depth: ndepth'
print*,'If s is negative, then value is ignored and sigma is calculated from SNR.'
return
endif
call getarg(1,arg)
read(arg,*) max_iterations
call getarg(2,arg)
read(arg,*) ndepth
call getarg(3,arg)
read(arg,*) ntrials
call getarg(4,arg)
read(arg,*) s
! scale Eb/No for a (174,91) code
rate=real(K)/real(N)
write(*,*) "rate: ",rate
write(*,*) "niter= ",max_iterations," s= ",s
allocate ( codeword(N), decoded(K), message(K) )
allocate ( rxdata(N), llr(N) )
msg="K1JT K9AN EN50"
! msg="G4WJS K9AN EN50"
call packmsg(msg,i4Msg6BitWords,itype,.false.) !Pack into 12 6-bit bytes
call unpackmsg(i4Msg6BitWords,msgsent,.false.,grid) !Unpack to get msgsent
write(*,*) "message sent ",msgsent
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
i1Msg8BitBytes(10:12)=0
checksum = crc14 (c_loc (i1Msg8BitBytes), 12)
! For reference, the next 3 lines show how to check the CRC
i1Msg8BitBytes(11)=checksum/256
i1Msg8BitBytes(12)=iand (checksum,255)
checksumok = crc14_check(c_loc (i1Msg8BitBytes), 12)
if( checksumok ) write(*,*) 'Good checksum'
! K=91, For now:
! msgbits(1:72) JT message bits
! msgbits(73:77) 5 free message bits (set to 0)
! msgbits(78:91) CRC14
mbit=0
do i=1, 9
i1=i1Msg8BitBytes(i)
do ibit=1,8
mbit=mbit+1
msgbits(mbit)=iand(1,ishft(i1,ibit-8))
enddo
enddo
msgbits(73:77)=0 ! the five extra message bits go here
i1=i1Msg8BitBytes(11) ! First 6 bits of crc12 are LSB of this byte
do ibit=1,6
msgbits(77+ibit)=iand(1,ishft(i1,ibit-6))
enddo
i1=i1Msg8BitBytes(12) ! Now shift in last 8 bits of the CRC
do ibit=1,8
msgbits(83+ibit)=iand(1,ishft(i1,ibit-8))
enddo
write(*,*) 'message'
write(*,'(12(8i1,1x))') msgbits
call encode174_91(msgbits,codeword)
call init_random_seed()
! call sgran()
write(*,*) 'codeword'
write(*,'(22(8i1,1x))') codeword
write(*,*) "Eb/N0 SNR2500 ngood nundetected nbadcrc sigma"
do idb = 20,-10,-1
!do idb = 0,0,-1
db=idb/2.0-1.0
sigma=1/sqrt( 2*rate*(10**(db/10.0)) )
ngood=0
nue=0
nbadcrc=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
if(nerr.ge.1) nerrtot(nerr)=nerrtot(nerr)+1
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)
nap=0 ! number of AP bits
llr(colorder(174-91+1:174-91+nap)+1)=5*(2.0*msgbits(1:nap)-1.0)
apmask=0
apmask(colorder(174-91+1:174-91+nap)+1)=1
! max_iterations is max number of belief propagation iterations
call bpdecode174_91(llr, apmask, max_iterations, decoded, cw, nharderrors,niterations)
if( ndepth .ge. 0 .and. nharderrors .lt. 0 ) call osd174_91(llr, apmask, ndepth, decoded, cw, nharderrors, dmin)
! If the decoder finds a valid codeword, nharderrors will be .ge. 0.
if( nharderrors .ge. 0 ) then
call extractmessage174_91(decoded,msgreceived,ncrcflag)
if( ncrcflag .ne. 1 ) then
nbadcrc=nbadcrc+1
endif
nueflag=0
nerrmpc=0
do i=1,K ! find number of errors in message+crc part of codeword
if( msgbits(i) .ne. decoded(i) ) then
nueflag=1
nerrmpc=nerrmpc+1
endif
enddo
if(nerrmpc.ge.1) nmpcbad(nerrmpc)=nmpcbad(nerrmpc)+1
if( ncrcflag .eq. 1 ) then
if( nueflag .eq. 0 ) then
ngood=ngood+1
if(nerr.ge.1) nerrdec(nerr)=nerrdec(nerr)+1
else if( nueflag .eq. 1 ) then
nue=nue+1;
endif
endif
endif
enddo
baud=12000/1920
snr2500=db+10.0*log10((baud/2500.0))
pberr=real(nberr)/(real(ntrials*N))
write(*,"(f4.1,4x,f5.1,1x,i8,1x,i8,1x,i8,8x,f5.2,8x,e10.3)") db,SNR2500,ngood,nue,nbadcrc,ss,pberr
enddo
open(unit=23,file='nerrhisto.dat',status='unknown')
do i=1,174
write(23,'(i4,2x,i10,i10,f10.2)') i,nerrdec(i),nerrtot(i),real(nerrdec(i))/real(nerrtot(i)+1e-10)
enddo
close(23)
open(unit=25,file='nmpcbad.dat',status='unknown')
do i=1,87
write(25,'(i4,2x,i10)') i,nmpcbad(i)
enddo
close(25)
end program ldpcsim174_91
lib/ft8/osd174_91.f90
+371
View File
@@ -0,0 +1,371 @@
subroutine osd174_91(llr,apmask,ndeep,decoded,cw,nhardmin,dmin)
!
! An ordered-statistics decoder for the (174,91) code.
!
integer, parameter:: N=174, K=91, M=N-K
integer*1 apmask(N),apmaskr(N)
integer*1 gen(K,N)
integer*1 genmrb(K,N),g2(N,K)
integer*1 temp(K),m0(K),me(K),mi(K),misub(K),e2sub(N-K),e2(N-K),ui(N-K)
integer*1 r2pat(N-K)
integer colorder(N)
integer indices(N),nxor(N)
integer*1 cw(N),ce(N),c0(N),hdec(N)
integer*1 decoded(K)
integer indx(N)
real llr(N),rx(N),absrx(N)
include "ldpc_174_91_c_generator.f90"
include "ldpc_174_91_c_colorder.f90"
logical first,reset
data first/.true./
save first,gen
if( first ) then ! fill the generator matrix
gen=0
do i=1,M
do j=1,23
read(g(i)(j:j),"(Z1)") istr
ibmax=4
if(j.eq.23) ibmax=3
do jj=1, ibmax
irow=(j-1)*4+jj
if( btest(istr,4-jj) ) gen(irow,i)=1
enddo
enddo
enddo
do irow=1,K
gen(irow,M+irow)=1
enddo
first=.false.
endif
! Re-order received vector to place systematic msg bits at the end.
rx=llr(colorder+1)
apmaskr=apmask(colorder+1)
! 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 mrbencode91(m0,c0,g2,N,K)
nxor=ieor(c0,hdec)
nhardmin=sum(nxor)
dmin=sum(nxor*absrx)
cw=c0
ntotal=0
nrejected=0
if(ndeep.eq.0) goto 998 ! norder=0
if(ndeep.gt.5) ndeep=5
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=0
nt=40
ntheta=12
ntau=19
elseif(ndeep.eq.5) then
nord=2
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
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 mrbencode91(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 mrbencode91(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 nextpat91(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 boxit91(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 mrbencode91(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 fetchit91(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 mrbencode91(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 nextpat91(misub,K,nord,iflag)
enddo
endif
998 continue
! Re-order the codeword to place message bits at the end.
cw(indices)=cw
hdec(indices)=hdec
decoded=cw(M+1:N)
cw(colorder+1)=cw ! put the codeword back into received-word order
return
end subroutine osd174_91
subroutine mrbencode91(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 mrbencode91
subroutine nextpat91(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 nextpat91
subroutine boxit91(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 boxit91
subroutine fetchit91(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 fetchit91
lib/ft8/sync8.f90
+23 -18
View File
@@ -10,14 +10,15 @@ subroutine sync8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase)
real x(NFFT1)
real sync2d(NH1,-JZ:JZ)
real red(NH1)
real candidate0(3,200)
real candidate(3,200)
real candidate0(4,200)
real candidate(4,200)
real dd(NMAX)
integer jpeak(NH1)
integer indx(NH1)
integer ii(1)
integer icos7(0:6)
data icos7/2,5,6,0,4,1,3/ !Costas 7x7 tone pattern
integer icos7_1(0:6),icos7_2(0:6),icos7(0:6)
data icos7_1/2,5,6,0,4,1,3/ !Costas 7x7 tone pattern
data icos7_2/3,1,4,0,6,5,2/ !Costas 7x7 tone pattern
equivalence (x,cx)
! Compute symbol spectra, stepping by NSTEP steps.
@@ -49,6 +50,11 @@ subroutine sync8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase)
nfos=NFFT1/NSPS ! # frequency bin oversampling factor
jstrt=0.5/tstep
candidate0=0.
k=0
do itype=1,2
if(itype.eq.1) icos7=icos7_1
if(itype.eq.2) icos7=icos7_2
do i=ia,ib
do j=-JZ,+JZ
ta=0.
@@ -58,16 +64,16 @@ subroutine sync8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase)
t0b=0.
t0c=0.
do n=0,6
k=j+jstrt+nssy*n
if(k.ge.1.and.k.le.NHSYM) then
ta=ta + s(i+nfos*icos7(n),k)
t0a=t0a + sum(s(i:i+nfos*6:nfos,k))
m=j+jstrt+nssy*n
if(m.ge.1.and.m.le.NHSYM) then
ta=ta + s(i+nfos*icos7(n),m)
t0a=t0a + sum(s(i:i+nfos*6:nfos,m))
endif
tb=tb + s(i+nfos*icos7(n),k+nssy*36)
t0b=t0b + sum(s(i:i+nfos*6:nfos,k+nssy*36))
if(k+nssy*72.le.NHSYM) then
tc=tc + s(i+nfos*icos7(n),k+nssy*72)
t0c=t0c + sum(s(i:i+nfos*6:nfos,k+nssy*72))
tb=tb + s(i+nfos*icos7(n),m+nssy*36)
t0b=t0b + sum(s(i:i+nfos*6:nfos,m+nssy*36))
if(m+nssy*72.le.NHSYM) then
tc=tc + s(i+nfos*icos7(n),m+nssy*72)
t0c=t0c + sum(s(i:i+nfos*6:nfos,m+nssy*72))
endif
enddo
t=ta+tb+tc
@@ -98,8 +104,6 @@ subroutine sync8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase)
base=red(ibase)
red=red/base
candidate0=0.
k=0
do i=1,200
n=ia + indx(iz+1-i) - 1
if(red(n).lt.syncmin) exit
@@ -107,7 +111,9 @@ subroutine sync8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase)
candidate0(1,k)=n*df
candidate0(2,k)=(jpeak(n)-1)*tstep
candidate0(3,k)=red(n)
candidate0(4,k)=itype
enddo
enddo
ncand=k
! Put nfqso at top of list, and save only the best of near-dupe freqs.
@@ -121,9 +127,6 @@ subroutine sync8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase)
if(candidate0(3,i).lt.candidate0(3,j)) candidate0(3,i)=0.
endif
enddo
! write(*,3001) i,candidate0(1,i-1),candidate0(1,i),candidate0(3,i-1), &
! candidate0(3,i)
!3001 format(i2,4f8.1)
endif
enddo
@@ -143,6 +146,8 @@ subroutine sync8(dd,nfa,nfb,syncmin,nfqso,s,candidate,ncand,sbase)
candidate(1,k)=abs(candidate0(1,j))
candidate(2,k)=candidate0(2,j)
candidate(3,k)=candidate0(3,j)
candidate(4,k)=candidate0(4,j)
!write(*,*) i,candidate(1:4,k)
k=k+1
endif
enddo
lib/ft8/sync8d.f90
+10 -3
View File
@@ -1,4 +1,4 @@
subroutine sync8d(cd0,i0,ctwk,itwk,sync)
subroutine sync8d(cd0,i0,ctwk,itwk,itype,sync)
! Compute sync power for a complex, downsampled FT8 signal.
@@ -9,8 +9,9 @@ subroutine sync8d(cd0,i0,ctwk,itwk,sync)
complex ctwk(32)
complex z1,z2,z3
logical first
integer icos7(0:6)
data icos7/2,5,6,0,4,1,3/
integer icos7_1(0:6),icos7_2(0:6),icos7(0:6)
data icos7_1/2,5,6,0,4,1,3/
data icos7_2/3,1,4,0,6,5,2/
data first/.true./
save first,twopi,fs2,dt2,taus,baud,csync
@@ -34,6 +35,12 @@ subroutine sync8d(cd0,i0,ctwk,itwk,sync)
first=.false.
endif
if(itype.eq.1) then
icos7=icos7_1
else
icos7=icos7_2
endif
sync=0
do i=0,6 !Sum over 7 Costas frequencies and
i1=i0+i*32 !three Costas arrays