dec.c

用于LDPC编码译码的仿真实现。包括随机生成校验矩阵、由校验矩阵产生生成矩阵、编码、加随机噪声、译码等内容。原作者是老外

/* DEC.C - Decoding procedures. */

/* Copyright (c) 2000, 2001 by Radford M. Neal 
 *
 * Permission is granted for anyone to copy, use, or modify this program 
 * for purposes of research or education, provided this copyright notice 
 * is retained, and note is made of any changes that have been made. 
 *
 * This program is distributed without any warranty, express or implied.
 * As this program was written for research purposes only, it has not been
 * tested to the degree that would be advisable in any important application.
 * All use of this program is entirely at the user's own risk.
 */


/* NOTE:  See decoding.html for general documentation on the decoding methods */


#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#include "alloc.h"
#include "mod2sparse.h"
#include "mod2dense.h"
#include "rand.h"
#include "rcode.h"
#include "check.h"
#include "dec.h"
#include "enc.h"





int max_iter;	/* Maximum number of iteratons of decoding to do */
char *gen_file;	/* Generator file for Enum_block and Enum_bit */



/* DECODE USING PROBABILITY PROPAGATION.  Tries to find the most probable 
   values for the bits of the codeword, given a parity check matrix (H), and
   likelihood ratios (lratio) for each bit.  If max_iter is positive, up to 
   that many iterations of probability propagation are done, stopping before 
   then if the tentative decoding is a valid codeword.  If max_iter is 
   negative, abs(max_iter) iterations are done, regardless of whether a 
   codeword was found earlier.

   Returns the number of iterations done (as an "unsigned" for consistency
   with enum_decode).  Regardless of whether or not a valid codeword was 
   reached, the bit vector from thresholding the bit-by-bit probabilities is 
   stored in dblk, and the resulting parity checks are stored in pchk (all 
   will be zero if the codeword is valid).  The final probabilities for each 
   bit being a 1 are stored in bprb.

   The setup procedure immediately below outputs headers for the detailed trace
   file, if required.
*/



/* INITIALIZE PROBABILITY PROPAGATION.  Stores initial ratios, probabilities,
   and guess at decoding. */

void initprp
( mod2sparse *H,	/* Parity check matrix */
  double *lratio,	/* Likelihood ratios for bits */
  int *dblk		/* Place to store decoding */
)
{ 
  mod2entry *e;
  int N;
  int j;

  N = mod2sparse_cols(H);

  for (j = 0; j<N; j++)
  { for (e = mod2sparse_first_in_col(H,j);
         !mod2sparse_at_end(e);
         e = mod2sparse_next_in_col(e))
    { e->pr = lratio[j];// from bit to check
      e->lr = 0;        //from check to bit
    }
    dblk[j] = lratio[j]>=0;
  }
}


/* DO ONE ITERATION OF PROBABILITY PROPAGATION. */

void iterprp
( mod2sparse *H,	/* Parity check matrix */
  double *lratio,	/* Likelihood ratios for bits */
  int *dblk,		/* Place to store decoding */
  double *LLR
)
{
  double pr, dl, t;
  mod2entry *e;
  int N, M;
  int i, j;

  M = mod2sparse_rows(H);
  N = mod2sparse_cols(H);

   /* Recompute likelihood ratios. */

  for (i = 0; i<M; i++)
  { dl = 1;
    for (e = mod2sparse_first_in_row(H,i);  /* Find the first entry in a row */
         !mod2sparse_at_end(e);             /* See if we've reached the end     */
         e = mod2sparse_next_in_row(e))    /*Move from one entry to the right direction*/
    { e->lr = dl;
      dl *= tanh(-e->pr/2);
    }
    dl = 1;
    for (e = mod2sparse_last_in_row(H,i); /* Find the last entry in a row */
         !mod2sparse_at_end(e);
         e = mod2sparse_prev_in_row(e))   /*Move from one entry to the left direction*/
    { t = e->lr * dl;
		  if (fabs(fabs(t)-1)<1e-12)
		  {
			 e->lr=-t*30;
			// printf("%f,",t);
		  }
		 else
			 e->lr = log((1-t)/(1+t));
      dl *= tanh(-e->pr/2);
    }
  }

 /* Recompute probability ratios.  Also find the next guess based on the
     individually most likely values. */

  for (j = 0; j<N; j++)
  { pr = lratio[j];
    for (e = mod2sparse_first_in_col(H,j); /* Find the first entry in a column */
         !mod2sparse_at_end(e);
         e = mod2sparse_next_in_col(e))   /*Move from one entry to the down direction*/
    { e->pr = pr;
      pr += e->lr;
    }

    LLR[j]=pr;//lxh
	dblk[j] =(pr>=0);
    pr = 0;
    for (e = mod2sparse_last_in_col(H,j); /* Find the last entry in a column */
         !mod2sparse_at_end(e);
         e = mod2sparse_prev_in_col(e))  /* Find the up entry in a column */
    { e->pr += pr;
/*      if (isnan(e->pr)) 
      { e->pr = 1;
      }
*/      pr += e->lr;
    }
  }

}


void WBF_initprp
( mod2sparse *H,	/* Parity check matrix */
  double *lratio,	/* Likelihood ratios for bits */
  int *dblk		/* Place to store decoding */
)
{ 
  mod2entry *e;
  int M,N;
  int i,j;
  double y;

  M = mod2sparse_rows(H);
  N = mod2sparse_cols(H);

  for (j = 0; j<N; j++)
  { for (e = mod2sparse_first_in_col(H,j);
         !mod2sparse_at_end(e);
         e = mod2sparse_next_in_col(e))
    { e->pr = lratio[j];// from bit to check
    }
    dblk[j] = lratio[j]>=0;
  }
  for (i = 0; i<M; i++)
  { e = mod2sparse_first_in_row(H,i);
    y=fabs(e->pr);

    for (e = mod2sparse_first_in_row(H,i);  /* Find the first entry in a row */
         !mod2sparse_at_end(e);             /* See if we've reached the end     */
         e = mod2sparse_next_in_row(e))    /*Move from one entry to the right direction*/
    { 
			 if (fabs(e->pr) < y)
			 y = fabs(e->pr);
    }
    for (e = mod2sparse_first_in_row(H,i);  /* Find the first entry in a row */
         !mod2sparse_at_end(e);             /* See if we've reached the end     */
         e = mod2sparse_next_in_row(e))  
		 { e->lr = y;        //from check to bit
		 }
  }

}




void WBF_iterprp
( mod2sparse *H,	/* Parity check matrix */
  int *dblk,		/* Place to store decoding */
  int *chks
)
{
  double pr, max_pr;
  int max_pos;
  mod2entry *e;
  int N, M;
  int j;

  M = mod2sparse_rows(H);
  N = mod2sparse_cols(H);

   /* Recompute likelihood ratios. */

  
 /* Recompute probability ratios.  Also find the next guess based on the
     individually most likely values. */

    pr = 0;
    for (e = mod2sparse_first_in_col(H,0); /* Find the first entry in a column */
         !mod2sparse_at_end(e);
         e = mod2sparse_next_in_col(e))   /*Move from one entry to the down direction*/
    { pr+= (2*chks[e->row]-1)*(e->lr);
    }

    max_pr=pr;
	max_pos=0;

  
  for (j = 1; j<N; j++)
  { pr = 0;
    for (e = mod2sparse_first_in_col(H,j); /* Find the first entry in a column */
         !mod2sparse_at_end(e);
         e = mod2sparse_next_in_col(e))   /*Move from one entry to the down direction*/
    { pr+= (2*chks[e->row]-1)*(e->lr);
    }

	if (pr>max_pr)
	{
		max_pr=pr;
		max_pos=j;
	}
   
  }
  dblk[max_pos]^=1;

}

void modified_WBF_iterprp
( mod2sparse *H,	/* Parity check matrix */
 double *lratio,
  int *dblk,		/* Place to store decoding */
  int *chks
)
{
  double y,pr, max_pr;
  int max_pos;
  mod2entry *e;
  int N, M;
  int i, j;

  M = mod2sparse_rows(H);
  N = mod2sparse_cols(H);

   /* Recompute likelihood ratios. */

  for (i = 0; i<M; i++)
  { e = mod2sparse_first_in_row(H,i);
    y=fabs(e->pr);

    for (e = mod2sparse_first_in_row(H,i);  /* Find the first entry in a row */
         !mod2sparse_at_end(e);             /* See if we've reached the end     */
         e = mod2sparse_next_in_row(e))    /*Move from one entry to the right direction*/
    { 
			 if (fabs(e->pr) < y)
			 y = fabs(e->pr);
    }
    for (e = mod2sparse_first_in_row(H,i);  /* Find the first entry in a row */
         !mod2sparse_at_end(e);             /* See if we've reached the end     */
         e = mod2sparse_next_in_row(e))  
		 { e->lr = y;        //from check to bit
		 }
  }

 /* Recompute probability ratios.  Also find the next guess based on the
     individually most likely values. */

    pr = 0;
    for (e = mod2sparse_first_in_col(H,0); /* Find the first entry in a column */
         !mod2sparse_at_end(e);
         e = mod2sparse_next_in_col(e))   /*Move from one entry to the down direction*/
    { pr+= (2*chks[e->row]-1)*(e->lr);
    }

    max_pr=pr-lratio[0];
	max_pos=0;

  
  for (j = 1; j<N; j++)
  { pr = 0;
    for (e = mod2sparse_first_in_col(H,j); /* Find the first entry in a column */
         !mod2sparse_at_end(e);
         e = mod2sparse_next_in_col(e))   /*Move from one entry to the down direction*/
    { pr+= (2*chks[e->row]-1)*(e->lr);
    }

	if (pr-lratio[j]>max_pr)
	{
		max_pr=pr-lratio[j];
		max_pos=j;
	}
   
  }
  dblk[max_pos]^=1;

}


void WBF_iterprp_lxh
( mod2sparse *H,	/* Parity check matrix */
 double *lratio,
  int *dblk,		/* Place to store decoding */
  int *chks,
  int Num_flip,
  int *encoded_data
)
{
  double y,*pr;
  int *pos;
  mod2entry *e;
  int N, M;
  int i, j;

  M = mod2sparse_rows(H);
  N = mod2sparse_cols(H);

  pr = chk_alloc (N, sizeof *pr);
  pos = chk_alloc (N, sizeof *pos);

  for(i=0;i<N;i++)
	  pos[i]=i;



   /* Recompute likelihood ratios. */

  for (i = 0; i<M; i++)
  { e = mod2sparse_first_in_row(H,i);
    y=fabs(e->pr);

    for (e = mod2sparse_first_in_row(H,i);  /* Find the first entry in a row */
         !mod2sparse_at_end(e);             /* See if we've reached the end     */
         e = mod2sparse_next_in_row(e))    /*Move from one entry to the right direction*/
    { 
			 if (fabs(e->pr) < y)
			 y = fabs(e->pr);
    }
    for (e = mod2sparse_first_in_row(H,i);  /* Find the first entry in a row */
         !mod2sparse_at_end(e);             /* See if we've reached the end     */
         e = mod2sparse_next_in_row(e))  
		 { e->lr = y;        //from check to bit
		 }
  }

 /* Recompute probability ratios.  Also find the next guess based on the
     individually most likely values. */

     
  for (j = 0; j<N; j++)
  { pr[j] = 0;
    for (e = mod2sparse_first_in_col(H,j); /* Find the first entry in a column */
         !mod2sparse_at_end(e);
         e = mod2sparse_next_in_col(e))   /*Move from one entry to the down direction*/
    { pr[j]+= (2*chks[e->row]-1)*(e->lr);
    }

  }

  bubble(pr, N, pos);

  for(i=0;i<Num_flip;i++)
  {
	  dblk[pos[i]]^=1;
	 // printf("%d,%d\n",encoded_data[pos[i]],dblk[pos[i]]);
  }