diff --git a/CMakeLists.txt b/CMakeLists.txt index f77a46033..6e295356c 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -299,6 +299,7 @@ set (wsjt_FSRCS lib/badmsg.f90 lib/bpdecode40.f90 lib/bpdecode144.f90 + lib/fsk4hf/bpdecode120.f90 lib/fsk4hf/bpdecode168.f90 lib/baddata.f90 lib/ccf2.f90 @@ -324,6 +325,7 @@ set (wsjt_FSRCS lib/encode4.f90 lib/encode_msk40.f90 lib/encode_msk144.f90 + lib/fsk4hf/encode120.f90 lib/fsk4hf/encode168.f90 lib/entail.f90 lib/ephem.f90 @@ -390,6 +392,7 @@ set (wsjt_FSRCS lib/jt9_decode.f90 lib/jt9fano.f90 lib/ldpcsim144.f90 + lib/fsk4hf/ldpcsim120.f90 lib/fsk4hf/ldpcsim168.f90 lib/ldpcsim40.f90 lib/libration.f90 @@ -1084,10 +1087,13 @@ target_link_libraries (jt65 wsjt_fort wsjt_cxx) add_executable (ldpcsim40 lib/ldpcsim40.f90 wsjtx.rc) target_link_libraries (ldpcsim40 wsjt_fort wsjt_cxx) +add_executable (ldpcsim120 lib/fsk4hf/ldpcsim120.f90 lib/crc.cpp wsjtx.rc) +target_link_libraries (ldpcsim120 wsjt_fort wsjt_cxx) + add_executable (ldpcsim144 lib/ldpcsim144.f90 wsjtx.rc) target_link_libraries (ldpcsim144 wsjt_fort wsjt_cxx) -add_executable (ldpcsim168 lib/fsk4hf/ldpcsim168.f90 wsjtx.rc) +add_executable (ldpcsim168 lib/fsk4hf/ldpcsim168.f90 lib/crc.cpp wsjtx.rc) target_link_libraries (ldpcsim168 wsjt_fort wsjt_cxx) add_executable (msk144sim lib/msk144sim.f90 wsjtx.rc) diff --git a/lib/crc12.cpp b/lib/crc.cpp similarity index 59% rename from lib/crc12.cpp rename to lib/crc.cpp index 7077aaba7..c4dd12b90 100644 --- a/lib/crc12.cpp +++ b/lib/crc.cpp @@ -5,6 +5,8 @@ extern "C" { short crc12 (unsigned char const * data, int length); bool crc12_check (unsigned char const * data, int length); + short crc10 (unsigned char const * data, int length); + bool crc10_check (unsigned char const * data, int length); } // assumes CRC is last 16 bits of the data and is set to zero @@ -18,3 +20,13 @@ bool crc12_check (unsigned char const * data, int length) { return !boost::augmented_crc<12, 0xc06> (data, length); } + +short crc10 (unsigned char const * data, int length) +{ + return boost::augmented_crc<10, 0x08f> (data, length); +} + +bool crc10_check (unsigned char const * data, int length) +{ + return !boost::augmented_crc<10, 0x08f> (data, length); +} diff --git a/lib/crc.f90 b/lib/crc.f90 index 8ac68e051..6e5c262eb 100644 --- a/lib/crc.f90 +++ b/lib/crc.f90 @@ -16,5 +16,22 @@ module crc type (c_ptr), value :: data integer (c_int), value :: length end function crc12_check + + function crc10 (data, length) bind (C, name="crc10") + use, intrinsic :: iso_c_binding, only: c_short, c_ptr, c_int + implicit none + integer (c_short) :: crc10 + type (c_ptr), value :: data + integer (c_int), value :: length + end function crc10 + + function crc10_check (data, length) bind (C, name="crc10_check") + use, intrinsic :: iso_c_binding, only: c_bool, c_ptr, c_int + implicit none + logical (c_bool) :: crc10_check + type (c_ptr), value :: data + integer (c_int), value :: length + end function crc10_check + end interface end module crc diff --git a/lib/fsk4hf/bpdecode120.f90 b/lib/fsk4hf/bpdecode120.f90 new file mode 100644 index 000000000..de6450deb --- /dev/null +++ b/lib/fsk4hf/bpdecode120.f90 @@ -0,0 +1,306 @@ +subroutine bpdecode120(llr,apmask,maxiterations,decoded,niterations) +! +! A log-domain belief propagation decoder for the (120,60) code. +! +integer, parameter:: N=120, K=60, M=N-K +integer*1 codeword(N),cw(N),apmask(N) +integer colorder(N) +integer*1 decoded(K) +integer Nm(7,M) ! 5, 6, or 7 bits per check +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) + +data colorder/ & + 0,1,2,21,3,4,5,6,7,8,20,10,9,11,12,23,13,28,14,31, & + 15,16,22,26,17,30,18,29,25,32,41,34,19,33,27,36,38,43,42,24, & + 37,39,45,40,35,44,47,46,50,51,53,48,52,56,54,57,55,49,58,61, & + 60,59,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76,77,78,79, & + 80,81,82,83,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/ + +data Mn/ & + 1, 18, 48, & + 2, 4, 51, & + 3, 23, 47, & + 5, 36, 42, & + 6, 43, 49, & + 7, 24, 55, & + 8, 35, 60, & + 9, 26, 30, & + 10, 29, 45, & + 11, 13, 46, & + 12, 53, 54, & + 14, 20, 57, & + 15, 16, 58, & + 17, 39, 44, & + 19, 37, 41, & + 21, 28, 34, & + 22, 50, 59, & + 25, 31, 52, & + 27, 32, 38, & + 33, 40, 56, & + 1, 11, 47, & + 2, 10, 16, & + 3, 12, 27, & + 4, 24, 28, & + 5, 23, 60, & + 6, 29, 39, & + 7, 31, 54, & + 8, 50, 56, & + 9, 13, 14, & + 15, 22, 41, & + 17, 26, 40, & + 18, 25, 45, & + 19, 20, 55, & + 21, 30, 36, & + 32, 49, 59, & + 33, 53, 58, & + 34, 38, 46, & + 29, 35, 57, & + 37, 43, 48, & + 42, 51, 52, & + 7, 11, 44, & + 1, 42, 58, & + 2, 13, 49, & + 3, 20, 40, & + 4, 18, 56, & + 5, 45, 55, & + 6, 21, 31, & + 8, 46, 52, & + 9, 12, 48, & + 10, 37, 38, & + 14, 15, 25, & + 16, 17, 60, & + 19, 39, 53, & + 22, 44, 51, & + 23, 28, 41, & + 24, 32, 35, & + 26, 45, 59, & + 27, 33, 36, & + 30, 47, 54, & + 34, 50, 57, & + 33, 43, 55, & + 1, 41, 57, & + 2, 40, 54, & + 3, 6, 24, & + 4, 11, 59, & + 5, 13, 56, & + 7, 16, 34, & + 8, 19, 26, & + 9, 31, 58, & + 10, 21, 53, & + 12, 22, 60, & + 14, 38, 51, & + 15, 43, 46, & + 17, 48, 50, & + 18, 27, 39, & + 20, 28, 44, & + 23, 25, 49, & + 4, 29, 36, & + 30, 32, 52, & + 35, 37, 47, & + 39, 42, 59, & + 1, 21, 40, & + 2, 50, 55, & + 3, 8, 10, & + 5, 31, 37, & + 6, 14, 60, & + 7, 36, 49, & + 9, 34, 39, & + 11, 19, 25, & + 12, 52, 57, & + 13, 22, 29, & + 15, 30, 56, & + 16, 18, 20, & + 17, 24, 46, & + 23, 38, 58, & + 26, 28, 43, & + 2, 27, 41, & + 5, 32, 44, & + 33, 47, 51, & + 35, 48, 53, & + 42, 43, 54, & + 34, 45, 47, & + 1, 8, 49, & + 3, 14, 59, & + 4, 31, 46, & + 6, 20, 50, & + 7, 26, 53, & + 9, 10, 36, & + 11, 58, 60, & + 12, 21, 45, & + 13, 28, 33, & + 15, 17, 35, & + 16, 38, 52, & + 18, 41, 54, & + 19, 23, 32, & + 22, 40, 55, & + 24, 25, 42, & + 26, 27, 56, & + 29, 44, 54, & + 30, 37, 55/ + +data Nm/ & + 1, 21, 42, 62, 82, 103, 0, & + 2, 22, 43, 63, 83, 97, 0, & + 3, 23, 44, 64, 84, 104, 0, & + 2, 24, 45, 65, 78, 105, 0, & + 4, 25, 46, 66, 85, 98, 0, & + 5, 26, 47, 64, 86, 106, 0, & + 6, 27, 41, 67, 87, 107, 0, & + 7, 28, 48, 68, 84, 103, 0, & + 8, 29, 49, 69, 88, 108, 0, & + 9, 22, 50, 70, 84, 108, 0, & + 10, 21, 41, 65, 89, 109, 0, & + 11, 23, 49, 71, 90, 110, 0, & + 10, 29, 43, 66, 91, 111, 0, & + 12, 29, 51, 72, 86, 104, 0, & + 13, 30, 51, 73, 92, 112, 0, & + 13, 22, 52, 67, 93, 113, 0, & + 14, 31, 52, 74, 94, 112, 0, & + 1, 32, 45, 75, 93, 114, 0, & + 15, 33, 53, 68, 89, 115, 0, & + 12, 33, 44, 76, 93, 106, 0, & + 16, 34, 47, 70, 82, 110, 0, & + 17, 30, 54, 71, 91, 116, 0, & + 3, 25, 55, 77, 95, 115, 0, & + 6, 24, 56, 64, 94, 117, 0, & + 18, 32, 51, 77, 89, 117, 0, & + 8, 31, 57, 68, 96, 107, 118, & + 19, 23, 58, 75, 97, 118, 0, & + 16, 24, 55, 76, 96, 111, 0, & + 9, 26, 38, 78, 91, 119, 0, & + 8, 34, 59, 79, 92, 120, 0, & + 18, 27, 47, 69, 85, 105, 0, & + 19, 35, 56, 79, 98, 115, 0, & + 20, 36, 58, 61, 99, 111, 0, & + 16, 37, 60, 67, 88, 102, 0, & + 7, 38, 56, 80, 100, 112, 0, & + 4, 34, 58, 78, 87, 108, 0, & + 15, 39, 50, 80, 85, 120, 0, & + 19, 37, 50, 72, 95, 113, 0, & + 14, 26, 53, 75, 81, 88, 0, & + 20, 31, 44, 63, 82, 116, 0, & + 15, 30, 55, 62, 97, 114, 0, & + 4, 40, 42, 81, 101, 117, 0, & + 5, 39, 61, 73, 96, 101, 0, & + 14, 41, 54, 76, 98, 119, 0, & + 9, 32, 46, 57, 102, 110, 0, & + 10, 37, 48, 73, 94, 105, 0, & + 3, 21, 59, 80, 99, 102, 0, & + 1, 39, 49, 74, 100, 0, 0, & + 5, 35, 43, 77, 87, 103, 0, & + 17, 28, 60, 74, 83, 106, 0, & + 2, 40, 54, 72, 99, 0, 0, & + 18, 40, 48, 79, 90, 113, 0, & + 11, 36, 53, 70, 100, 107, 0, & + 11, 27, 59, 63, 101, 114, 119, & + 6, 33, 46, 61, 83, 116, 120, & + 20, 28, 45, 66, 92, 118, 0, & + 12, 38, 60, 62, 90, 0, 0, & + 13, 36, 42, 69, 95, 109, 0, & + 17, 35, 57, 65, 81, 104, 0, & + 7, 25, 52, 71, 86, 109, 0/ + +data nrw/ & +6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6, & +6,6,6,6,6,7,6,6,6,6,6,6,6,6,6,6,6,6,6,6, & +6,6,6,6,6,6,6,5,6,6,5,6,6,7,7,6,5,6,6,6/ + +ncw=3 + +toc=0 +tov=0 +tanhtoc=0 +!write(*,*) llr +! 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 + niterations=iter + codeword=cw(colorder+1) + decoded=codeword(M+1:N) + return + endif + + if( iter.gt.0 ) 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. 3 .and. iter .ge. 5 .and. ncheck .gt. 10) then + niterations=-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 +niterations=-1 +return +end subroutine bpdecode120 diff --git a/lib/fsk4hf/encode120.f90 b/lib/fsk4hf/encode120.f90 new file mode 100644 index 000000000..577bbe860 --- /dev/null +++ b/lib/fsk4hf/encode120.f90 @@ -0,0 +1,116 @@ +subroutine encode120(message,codeword) +! Encode an 60-bit message and return a 120-bit codeword. +! The generator matrix has dimensions (60,60). +! The code is a (120,60) regular ldpc code with column weight 3. +! The code was generated using the PEG algorithm. +! After creating the codeword, the columns are re-ordered according to +! "colorder" to make the codeword compatible with the parity-check matrix +! +character*15 g(60) +integer*1 codeword(120) +integer colorder(120) +integer*1 gen(60,60) +integer*1 itmp(120) +integer*1 message(60) +integer*1 pchecks(60) +logical first +data first/.true./ +data g/ & + "65541ad98feab6e",& + "27249940a5895a3",& + "c80eac7506bf794",& + "aa50393e3e18d3f",& + "28527e87d47dced",& + "5da0dcaf8db048c",& + "d6509a43ca9b01a",& + "9a7dadd9c94f1d4",& + "bb673d3ba07cf29",& + "65e190f2fbed447",& + "bc2062a4e520969",& + "9e357f3feed059b",& + "aa6b59212036a57",& + "f78a326722d6565",& + "416754bc34c6405",& + "f77000b3f04ff67",& + "d48fbd7d48c5ab9",& + "031ffb5db3a70cb",& + "125964e358c4df5",& + "bd02c32a5a241ea",& + "4c15ecdd8561abd",& + "7f0f1b352c7413e",& + "26edb94dfd0ae79",& + "ca1ba1ee0f8fb24",& + "49878a58cb4544c",& + "3dbcd0ff821b203",& + "c1f4440160d5345",& + "b5ea9dc7a5a70ab",& + "cebcf7d94976be4",& + "0968265f5977c88",& + "c5a36937faa78c3",& + "f0d4fef11e01c10",& + "e35fc0c779bebfe",& + "cf49c3eb41a31d5",& + "3f0b19352c7013e",& + "0e15eccd8521abd",& + "dda8dcaf9d3048c",& + "fee31438fba59ed",& + "ad74a27e939189c",& + "736ac01b439106e",& + "ab5d2729b29bfa1",& + "edf11fb02e5a426",& + "5f38be1c93ecc83",& + "1e4b3b8dc516b3e",& + "84443d8bee614c6",& + "d854d9f355ceac4",& + "a476b5ece51f0ea",& + "831c2b36c4c2f68",& + "f485c97a91615ae",& + "e9376d828ade9ba",& + "cac586f089d3185",& + "b8f8c67613dafe2",& + "1a3142b401b315d",& + "87dbedc43265d2e",& + "bb64ec6e652e7da",& + "e71bfd4c95dfd38",& + "31209af07ad4f75",& + "cff1a8ccc5f4978",& + "742eded1e1dfefd",& + "1cd7154a904dac4"/ + +data colorder/ & + 0,1,2,21,3,4,5,6,7,8,20,10,9,11,12,23,13,28,14,31, & + 15,16,22,26,17,30,18,29,25,32,41,34,19,33,27,36,38,43,42,24, & + 37,39,45,40,35,44,47,46,50,51,53,48,52,56,54,57,55,49,58,61, & + 60,59,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76,77,78,79, & + 80,81,82,83,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/ + +save first,gen + +if( first ) then ! fill the generator matrix + gen=0 + do i=1,60 + do j=1,15 + read(g(i)(j:j),"(Z1)") istr + do jj=1, 4 + icol=(j-1)*4+jj + if( btest(istr,4-jj) ) gen(i,icol)=1 + enddo + enddo + enddo +first=.false. +endif + +do i=1, 60 + nsum=0 + do j=1, 60 + nsum=nsum+message(j)*gen(i,j) + enddo + pchecks(i)=mod(nsum,2) +enddo +itmp(1:60)=pchecks +itmp(61:120)=message(1:60) +codeword(colorder+1)=itmp(1:120) + +return +end subroutine encode120 diff --git a/lib/fsk4hf/ldpcsim120.f90 b/lib/fsk4hf/ldpcsim120.f90 new file mode 100644 index 000000000..ba684f88b --- /dev/null +++ b/lib/fsk4hf/ldpcsim120.f90 @@ -0,0 +1,219 @@ +program ldpcsim120 +! End to end test of the (120,60)/crc10 encoder and decoder. +use crc +use packjt + +parameter(NRECENT=10) +character*12 recent_calls(NRECENT) +character*22 msg,msgsent,msgreceived +character*8 arg +integer*1, allocatable :: codeword(:), decoded(:), message(:) +integer*1, target:: i1Msg8BitBytes(9) +integer*1, target:: i1Dec8BitBytes(9) +integer*1 msgbits(60) +integer*1 apmask(120) +integer*2 checksum +integer colorder(120) +integer nerrtot(120),nerrdec(120),nmpcbad(60) +logical checksumok,fsk,bpsk +real*8, allocatable :: rxdata(:) +real, allocatable :: llr(:) + +data colorder/ & + 0,1,2,21,3,4,5,6,7,8,20,10,9,11,12,23,13,28,14,31, & + 15,16,22,26,17,30,18,29,25,32,41,34,19,33,27,36,38,43,42,24, & + 37,39,45,40,35,44,47,46,50,51,53,48,52,56,54,57,55,49,58,61, & + 60,59,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76,77,78,79, & + 80,81,82,83,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/ + +do i=1,NRECENT + recent_calls(i)=' ' +enddo +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.3) then + print*,'Usage: ldpcsim niter #trials s ' + print*,'eg: ldpcsim 10 1000 0.75' + return +endif +call getarg(1,arg) +read(arg,*) max_iterations +call getarg(2,arg) +read(arg,*) ntrials +call getarg(3,arg) +read(arg,*) s + +fsk=.false. +bpsk=.true. + +! don't count crc bits as data bits +N=120 +K=60 +! scale Eb/No for a (120,50) code +rate=real(50)/real(N) + +write(*,*) "rate: ",rate +write(*,*) "niter= ",max_iterations," s= ",s + +allocate ( codeword(N), decoded(K), message(K) ) +allocate ( rxdata(N), llr(N) ) + +! The message should be packed into the first 7 bytes + i1Msg8BitBytes(1:6)=85 + i1Msg8BitBytes(7)=64 +! The CRC will be put into the last 2 bytes + i1Msg8BitBytes(8:9)=0 + checksum = crc10 (c_loc (i1Msg8BitBytes), 9) +! For reference, the next 3 lines show how to check the CRC + i1Msg8BitBytes(8)=checksum/256 + i1Msg8BitBytes(9)=iand (checksum,255) + checksumok = crc10_check(c_loc (i1Msg8BitBytes), 9) + if( checksumok ) write(*,*) 'Good checksum' +write(*,*) i1Msg8BitBytes(1:9) + + mbit=0 + do i=1, 7 + i1=i1Msg8BitBytes(i) + do ibit=1,8 + mbit=mbit+1 + msgbits(mbit)=iand(1,ishft(i1,ibit-8)) + enddo + enddo + i1=i1Msg8BitBytes(8) ! First 2 bits of crc10 are LSB of this byte + do ibit=1,2 + msgbits(50+ibit)=iand(1,ishft(i1,ibit-2)) + enddo + i1=i1Msg8BitBytes(9) ! Now shift in last 8 bits of the CRC + do ibit=1,8 + msgbits(52+ibit)=iand(1,ishft(i1,ibit-8)) + enddo + + write(*,*) 'message' + write(*,'(9(8i1,1x))') msgbits + + call encode120(msgbits,codeword) + call init_random_seed() + call sgran() + + write(*,*) 'codeword' + write(*,'(15(8i1,1x))') codeword + +write(*,*) "Es/N0 SNR2500 ngood nundetected nbadcrc sigma" +do idb = -10, 24 + 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 + nbadcrc=0 + nberr=0 + do itrial=1, ntrials +! Create a realization of a noisy received word + do i=1,N + if( bpsk ) then + rxdata(i) = 2.0*codeword(i)-1.0 + sigma*gran() + elseif( fsk ) then + if( codeword(i) .eq. 1 ) then + r1=(1.0 + sigma*gran())**2 + (sigma*gran())**2 + r2=(sigma*gran())**2 + (sigma*gran())**2 + elseif( codeword(i) .eq. 0 ) then + r2=(1.0 + sigma*gran())**2 + (sigma*gran())**2 + r1=(sigma*gran())**2 + (sigma*gran())**2 + endif + rxdata(i)=0.35*(sqrt(r1)-sqrt(r2)) +! rxdata(i)=0.35*(exp(r1)-exp(r2)) +! rxdata(i)=0.12*(log(r1)-log(r2)) + endif + enddo + nerr=0 + do i=1,N + if( rxdata(i)*(2*codeword(i)-1.0) .lt. 0 ) nerr=nerr+1 + enddo + nerrtot(nerr)=nerrtot(nerr)+1 + nberr=nberr+nerr + +! Correct signal normalization is important for this decoder. +! rxav=sum(rxdata)/N +! rx2av=sum(rxdata*rxdata)/N +! rxsig=sqrt(rx2av-rxav*rxav) +! rxdata=rxdata/rxsig +! To match the metric to the channel, s should be set to the noise standard deviation. +! For now, set s to the value that optimizes decode probability near threshold. +! The s parameter can be tuned to trade a few tenth's dB of threshold for an order of +! magnitude in UER + 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 bpdecode120(llr, apmask, max_iterations, decoded, niterations) + + +! If the decoder finds a valid codeword, niterations will be .ge. 0. + if( niterations .ge. 0 ) then +! Check the CRC + do ibyte=1,6 + itmp=0 + do ibit=1,8 + itmp=ishft(itmp,1)+iand(1,decoded((ibyte-1)*8+ibit)) + enddo + i1Dec8BitBytes(ibyte)=itmp + enddo + i1Dec8BitBytes(7)=decoded(49)*128+decoded(50)*64 +! Need to pack the received crc into bytes 8 and 9 for crc10_check + i1Dec8BitBytes(8)=decoded(51)*2+decoded(52) + i1Dec8BitBytes(9)=decoded(53)*128+decoded(54)*64+decoded(55)*32+decoded(56)*16 + i1Dec8BitBytes(9)=i1Dec8BitBytes(9)+decoded(57)*8+decoded(58)*4+decoded(59)*2+decoded(60)*1 + ncrcflag=0 + if( crc10_check( c_loc( i1Dec8BitBytes ), 9 ) ) ncrcflag=1 + + 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 + nmpcbad(nerrmpc)=nmpcbad(nerrmpc)+1 ! This histogram should inform our selection of CRC poly + if( ncrcflag .eq. 1 .and. nueflag .eq. 0 ) then + ngood=ngood+1 + nerrdec(nerr)=nerrdec(nerr)+1 + else if( ncrcflag .eq. 1 .and. nueflag .eq. 1 ) then + nue=nue+1; + endif + endif + enddo + snr2500=db+10*log10(0.4166/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,120 + 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,60 + write(25,'(i4,2x,i10)') i,nmpcbad(i) +enddo +close(25) + + + +end program ldpcsim120 diff --git a/lib/fsk4hf/ldpcsim168.f90 b/lib/fsk4hf/ldpcsim168.f90 index 8d0219a41..ca1bee1e4 100644 --- a/lib/fsk4hf/ldpcsim168.f90 +++ b/lib/fsk4hf/ldpcsim168.f90 @@ -14,7 +14,7 @@ integer*1 apmask(168) integer*2 checksum integer*4 i4Msg6BitWords(13) integer colorder(168) -integer nerrtot(168),nerrdec(168) +integer nerrtot(168),nerrdec(168),nmpcbad(84) logical checksumok,fsk,bpsk real*8, allocatable :: rxdata(:) real, allocatable :: llr(:) @@ -33,6 +33,7 @@ do i=1,NRECENT enddo 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.3) then @@ -47,8 +48,8 @@ read(arg,*) ntrials call getarg(3,arg) read(arg,*) s -fsk=.true. -bpsk=.false. +fsk=.false. +bpsk=.true. ! don't count crc bits as data bits N=168 @@ -62,7 +63,8 @@ write(*,*) "niter= ",max_iterations," s= ",s allocate ( codeword(N), decoded(K), message(K) ) allocate ( rxdata(N), llr(N) ) - msg="K1JT K9AN EN50" +! msg="K1JT K9AN EN50" + msg="G4WJS K9AN EN50" call packmsg(msg,i4Msg6BitWords,itype) !Pack into 12 6-bit bytes call unpackmsg(i4Msg6BitWords,msgsent) !Unpack to get msgsent write(*,*) "message sent ",msgsent @@ -180,22 +182,19 @@ do idb = -10, 24 ! If the decoder finds a valid codeword, niterations will be .ge. 0. if( niterations .ge. 0 ) then call extractmessage168(decoded,msgreceived,ncrcflag,recent_calls,nrecent) - if( ncrcflag .eq. 1 ) then - ncrcflag=1 - else - ncrcflag=0 - endif if( ncrcflag .ne. 1 ) then nbadcrc=nbadcrc+1 endif nueflag=0 -! Check the message plus crc against what was sent. - do i=1,K + 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 + nmpcbad(nerrmpc)=nmpcbad(nerrmpc)+1 if( ncrcflag .eq. 1 .and. nueflag .eq. 0 ) then ngood=ngood+1 nerrdec(nerr)=nerrdec(nerr)+1 @@ -211,9 +210,16 @@ do idb = -10, 24 enddo open(unit=23,file='nerrhisto.dat',status='unknown') -do i=1,128 +do i=1,168 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,84 + write(25,'(i4,2x,i10)') i,nmpcbad(i) +enddo +close(25) + + end program ldpcsim168