From 3b3288afa3e1252c3c0c8aa773afdd2b2c5039d4 Mon Sep 17 00:00:00 2001
From: Steven Franke <s.j.franke@icloud.com>
Date: Thu, 25 Feb 2016 03:04:07 +0000
Subject: [PATCH] msksim.f90 works, but needs more work.

git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@6493 ab8295b8-cf94-4d9e-aec4-7959e3be5d79
---
 lib/ldpc/Makefile.msk_osx |  2 +-
 lib/ldpc/dec.c            | 30 +++++++++--------
 lib/ldpc/enc.c            |  9 +++++
 lib/ldpc/msksim.f90       | 69 +++++++++++++++++++++++++++++----------
 lib/ldpc/rcode.c          | 16 +++------
 5 files changed, 83 insertions(+), 43 deletions(-)

diff --git a/lib/ldpc/Makefile.msk_osx b/lib/ldpc/Makefile.msk_osx
index 265871ec6..c5a475bb7 100644
--- a/lib/ldpc/Makefile.msk_osx
+++ b/lib/ldpc/Makefile.msk_osx
@@ -27,7 +27,7 @@ all:	msksim
 
 OBJS = msksim.o alloc.o rcode.o channel.o dec.o enc.o \
 	intio.o blockio.o check.o open.o mod2dense.o \
-	mod2sparse.o mod2convert.o distrib.o rand.o
+	mod2sparse.o mod2convert.o distrib.o rand.o gran.o
 
 msksim:	$(OBJS)
 	$(FC) -o msksim $(OBJS) 
diff --git a/lib/ldpc/dec.c b/lib/ldpc/dec.c
index 2a5e21479..431518d8f 100755
--- a/lib/ldpc/dec.c
+++ b/lib/ldpc/dec.c
@@ -384,19 +384,23 @@ void iterprp
   }
 }
 
-void ldpc_decode_ ( double *lratio, char *dblk, char *pchk, double *bprb )
+void ldpc_decode_ ( double lratio[], char decoded[], int *max_iterations, int *niterations)
 {
-  int i, n;
-  n = prprp_decode ( H, lratio, dblk, pchk, bprb );
-  printf("in ldpc_decode n=%d\n");
-  printf("dblk: ");
-  for( i=0; i<128; i++) {
-    printf("%d",dblk[i]);
+  int i, j, valid;
+  char dblk[N],pchk[M];
+  double bprb[N];
+    
+  max_iter=*max_iterations;
+  *niterations = prprp_decode ( H, lratio, dblk, pchk, bprb );
+  valid = check( H, dblk, pchk )==0;
+  if( !valid ) {
+      *niterations=-1;
+  };
+
+  j=0;
+  for( i=M; i<N; i++ ) {
+    decoded[j]=dblk[cols[i]];
+    j=j+1;
   }
-  printf("\n");
-  printf("pchk: ");
-  for( i=0; i<46; i++) {
-    printf("%d",pchk[i]);
-  }
-  printf("\n");
+
 } 
diff --git a/lib/ldpc/enc.c b/lib/ldpc/enc.c
index dd758017d..74a823abf 100755
--- a/lib/ldpc/enc.c
+++ b/lib/ldpc/enc.c
@@ -159,11 +159,20 @@ char msg[82],
 char cdw[128]
 ){
   int i;
+  char checks[46];
+
   printf("msg: ");
   for (i=0; i<82; i++) {
   printf("%d", msg[i]); }
   printf("\n");
   sparse_encode(msg,cdw);
+  mod2sparse_mulvec (H, cdw, checks);
+  for (i = 0; i<46; i++) {
+    if( checks[i] == 1 ) {
+      printf("Failed to encode.\n");
+      abort();
+    }
+  } 
   printf("cdw: ");
   for (i=0; i<128; i++) {
   printf("%d", cdw[i]); }
diff --git a/lib/ldpc/msksim.f90 b/lib/ldpc/msksim.f90
index 8a1cb3563..1cd2acae0 100644
--- a/lib/ldpc/msksim.f90
+++ b/lib/ldpc/msksim.f90
@@ -5,25 +5,60 @@ use, intrinsic :: iso_c_binding
 parameter (N=128, M=46, K=82) ! M and N are global variables on the C side. 
 integer(1) message(1:K)
 integer(1) codeword(1:N)
-integer(1) decoded(1:N)
-integer(1) pchk(1:M)
-real*8 lratio(N), bitprobs(N)
+integer(1) decoded(1:K)
+real*8 lratio(N)
+character(50) pchk_file,gen_file
 
-write(*,*) "calling init_ldpc"
-call init_ldpc()
+pchk_file="./jtmode_codes/ldpc-128-82-sf11.pchk"
+gen_file="./jtmode_codes/ldpc-128-82-sf11.gen"
 
-message(1:K)=0
-message(10)=1
-write(*,*) "calling ldpc_encode"
+call init_ldpc(trim(pchk_file)//char(0),trim(gen_file)//char(0))
+
+
+message(1:41)=1
+message(42:82)=0
 call ldpc_encode(message,codeword)
-write(*,*) "calling ldpc_decode"
-do i=1,N
-lratio(i)=exp(2.0*(codeword(i)-0.5))
-enddo
-lratio(10)=10.0
-bitprobs(1:N)=0.0 ! shouldn't need this
-call ldpc_decode(lratio, decoded, pchk, bitprobs)
-do i=1,N
-write(*,*) i,bitprobs(i), lratio(i)
+
+max_iterations=10
+ntrials=1000000
+rate=82.0/128.0
+
+write(*,*) "Eb/N0   ngood    nundetected"
+do idb = 0, 11
+  db=idb/2.0-0.5
+  sigma=1/sqrt( 2*rate*(10**(db/10.0)) )
+
+  ngood=0
+  nue=0
+  aviter=0.0
+
+  do itrial=1, ntrials
+
+    do i=1,N
+      rxdata = 2.0*(codeword(i)-0.5) + sigma*gran()
+      lratio(i)=exp(2.0*rxdata/(sigma*sigma))
+    enddo
+
+    call ldpc_decode(lratio, decoded, max_iterations, niterations)
+
+    if( niterations .ge. 0 ) then
+      nueflag=0
+      do i=1,K
+        if( message(i) .ne. decoded(i) ) then
+          nueflag=1
+        endif
+      enddo
+      if( nueflag .eq. 1 ) then
+        nue=nue+1
+      else
+        ngood=ngood+1;
+      endif
+    endif
+
+  enddo
+
+  write(*,"(f4.1,1x,i8,1x,i8)") db,ngood,nue
+
 enddo
+
 end program msksim
diff --git a/lib/ldpc/rcode.c b/lib/ldpc/rcode.c
index 15b091f45..5c7f2e8c8 100755
--- a/lib/ldpc/rcode.c
+++ b/lib/ldpc/rcode.c
@@ -199,19 +199,11 @@ garbled:
   exit(1);
 }
 
-// Fortran interface routine
-void init_ldpc_ ( )
+// Fortran interface routine for WSJT-X
+void init_ldpc_ (char *pfile, char *gfile )
 {
     char *pchk_file,*gen_file;
-    pchk_file="./jtmode_codes/ldpc-128-82-sf11.pchk";
-    gen_file="./jtmode_codes/ldpc-128-82-sf11.gen";
- 
-    printf("pchk_file %s\n",pchk_file);
-    printf("gen_file %s\n",gen_file);
-
-    read_pchk( pchk_file );
-    read_gen( gen_file, 0, 0 );
-
-    printf("N %d M %d\n",N,M);
+    read_pchk( pfile );
+    read_gen( gfile, 0, 0 );
 }