mod2sparse.c

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

/* MOD2SPARSE.C - Procedures for handling sparse mod2 matrices. */

/* 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 mod2sparse.html for documentation on these procedures. */


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

#include "alloc.h"
#include "intio.h"
#include "mod2sparse.h"
#include "rcode.h"


/* ALLOCATE AN ENTRY WITHIN A MATRIX.  This local procedure is used to
   allocate a new entry, representing a non-zero element, within a given
   matrix.  Entries in this matrix that were previously allocated and
   then freed are re-used.  If there are no such entries, a new block
   of entries is allocated. */
//为矩阵分配一个实体，该程序是给一个已给定的矩阵分配一个新的实体表示一个非零元素
//在该矩阵中有预先分配好的空实体，如果没有生成新的实体块
static mod2entry *alloc_entry
( mod2sparse *m
)
{ 
  mod2block *b;
  mod2entry *e;
  int k;
 
  if (m->next_free==0)     //如果该矩阵没有空的实体，则生成新的实体块
  { 
    b = chk_alloc (1, sizeof *b);

    b->next = m->blocks;
    m->blocks = b;

    for (k = 0; k<Mod2sparse_block; k++)
    { b->entry[k].left = m->next_free;
      m->next_free = &b->entry[k];
    }
  }

  e = m->next_free;      //如果矩阵有现成的空实体，把给实体连接上
  m->next_free = e->left;

  e->pr = 0;
  e->lr = 0;
  
  return e;
}


/* ALLOCATE SPACE FOR A SPARSE MOD2 MATRIX.  */
//为一个稀疏矩阵分配内存空间
mod2sparse *mod2sparse_allocate
( int n_rows, 		/* Number of rows in matrix */
  int n_cols		/* Number of columns in matrix */
)
{
  mod2sparse *m;
  mod2entry *e;
  int i, j;

  if (n_rows<=0 || n_cols<=0)   //提示行数或列数出错
  { fprintf(stderr,"mod2sparse_allocate: Invalid number of rows or columns\n");
    exit(1);
  }

  m = chk_alloc (1, sizeof *m);     //分配内存

  m->n_rows = n_rows;           //行数和列数赋值
  m->n_cols = n_cols;    

  m->rows = chk_alloc (n_rows, sizeof *m->rows);   //为行和列链表的头分配空间
  m->cols = chk_alloc (n_cols, sizeof *m->cols);

  m->blocks = 0;     
  m->next_free = 0;

  for (i = 0; i<n_rows; i++)
  { e = &m->rows[i];      //为第i行的行链表的头赋初值
    e->left = e->right = e->up = e->down = e;
    e->row = e->col = -1;
  }

  for (j = 0; j<n_cols; j++)
  { e = &m->cols[j];    //为第j行的列链表的头赋初值
    e->left = e->right = e->up = e->down = e;
    e->row = e->col = -1;
  }

  return m;
}


/* FREE SPACE OCCUPIED BY A SPARSE MOD2 MATRIX.释放一个稀疏矩阵所占的空间 */

void mod2sparse_free
( mod2sparse *m		/* Matrix to free */
)
{ 
  mod2block *b;

  free(m->rows);
  free(m->cols);

  while (m->blocks!=0)
  { b = m->blocks;
    m->blocks = b->next;
    free(b);
  }
}


/* CLEAR A SPARSE MATRIX TO ALL ZEROS.把一个稀疏矩阵全部置0，即释放全部实体 */

void mod2sparse_clear
( mod2sparse *r
)
{
  mod2block *b;
  mod2entry *e;
  int i, j;

  for (i = 0; i<mod2sparse_rows(r); i++)     
  { e = &r->rows[i];                     //把各行链表的头重新赋值
    e->left = e->right = e->up = e->down = e;    
  }

  for (j = 0; j<mod2sparse_cols(r); j++)
  { e = &r->cols[j];               //把各列链表头重新赋值
    e->left = e->right = e->up = e->down = e;
  }

  while (r->blocks!=0)     
  { 
	 b = r->blocks;           //逐个释放实体块
    r->blocks = b->next;
    free(b);
  }
}


/* COPY A SPARSE MATRIX. 复制一个稀疏矩阵*/

void mod2sparse_copy
( mod2sparse *m,	/* Matrix to copy */
  mod2sparse *r		/* Place to store copy of matrix */
)
{
  mod2entry *e, *f;
  int i;

  if (mod2sparse_rows(m)>mod2sparse_rows(r) 
   || mod2sparse_cols(m)>mod2sparse_cols(r))   
  { fprintf(stderr,"mod2sparse_copy: Destination matrix is too small\n");
    exit(1);    //若目标矩阵的行数或列数小与原矩阵，提示出错
  }

  mod2sparse_clear(r);     //先把目标矩阵清零

  for (i = 0; i<mod2sparse_rows(m); i++)
  {
    e = mod2sparse_first_in_row(m,i);  //找到矩阵m第i行的第一个实体

    while (!mod2sparse_at_end(e))
    { f = mod2sparse_insert(r,e->row,e->col);     //在相应的位置插入实体
      f->lr = e->lr;
      f->pr = e->pr;
      e = mod2sparse_next_in_row(e);
    }
  }
}


/* COPY ROWS OF A SPARSE MOD2 MATRIX. */

void mod2sparse_copyrows
( mod2sparse *m,	/* Matrix to copy */
  mod2sparse *r,	/* Place to store copy of matrix */
  int *rows		/* Indexes of rows to copy, from 0 */
)
{ 
  mod2entry *e;
  int i;

  if (mod2sparse_cols(m)>mod2sparse_cols(r))
  { fprintf(stderr,
     "mod2sparse_copyrows: Destination matrix has fewer columns than source\n");
    exit(1);
  }

  mod2sparse_clear(r);

  for (i = 0; i<mod2sparse_rows(r); i++)
  { if (rows[i]<0 || rows[i]>=mod2sparse_rows(m))
    { fprintf(stderr,"mod2sparse_copyrows: Row index out of range\n");
      exit(1);
    }
    e = mod2sparse_first_in_row(m,rows[i]);
    while (!mod2sparse_at_end(e))
    { mod2sparse_insert(r,i,e->col);
      e = mod2sparse_next_in_row(e);
    }
  }
}


/* COPY COLUMNS OF A SPARSE MOD2 MATRIX. */

void mod2sparse_copycols
( mod2sparse *m,	/* Matrix to copy */
  mod2sparse *r,	/* Place to store copy of matrix */
  int *cols		/* Indexes of columns to copy, from 0 */
)
{ 
  mod2entry *e;
  int j;

  if (mod2sparse_rows(m)>mod2sparse_rows(r))
  { fprintf(stderr,
      "mod2sparse_copycols: Destination matrix has fewer rows than source\n");
    exit(1);
  }

  mod2sparse_clear(r);

  for (j = 0; j<mod2sparse_cols(r); j++)
  { if (cols[j]<0 || cols[j]>=mod2sparse_cols(m))
    { fprintf(stderr,"mod2sparse_copycols: Column index out of range\n");
      exit(1);
    }
    e = mod2sparse_first_in_col(m,cols[j]);
    while (!mod2sparse_at_end(e))
    { mod2sparse_insert(r,e->row,j);
      e = mod2sparse_next_in_col(e);
    }
  }
}


/* PRINT A SPARSE MOD2 MATRIX IN HUMAN-READABLE FORM. */

void mod2sparse_print
( FILE *f,
  mod2sparse *m
)
{ 
  int rdigits, cdigits;
  mod2entry *e;
  int i;

  rdigits = mod2sparse_rows(m)<=10 ? 1 
          : mod2sparse_rows(m)<=100 ? 2
          : mod2sparse_rows(m)<=1000 ? 3
          : mod2sparse_rows(m)<=10000 ? 4
          : mod2sparse_rows(m)<=100000 ? 5
          : 6;

  cdigits = mod2sparse_cols(m)<=10 ? 1 
          : mod2sparse_cols(m)<=100 ? 2
          : mod2sparse_cols(m)<=1000 ? 3
          : mod2sparse_cols(m)<=10000 ? 4
          : mod2sparse_cols(m)<=100000 ? 5
          : 6;   //条件运算符右结合

  for (i = 0; i<mod2sparse_rows(m); i++)
  { 
    fprintf(f,"%*d:",rdigits,i);   //如果没有*,则是把rdigits，i两个变量输入到文件指针f所指文件中，加上*以后表示不将rdigits的值输入到文件中，但是她的位置存在

    e = mod2sparse_first_in_row(m,i);  //使e指向i行的第一个元素
    while (!mod2sparse_at_end(e))
    { fprintf(f," %*d",cdigits,mod2sparse_col(e));
      e = mod2sparse_next_in_row(e);
    }

    fprintf(f,"\n");
  }
}


/* WRITE A SPARSE MOD2 MATRIX TO A FILE IN MACHINE-READABLE FORM. 以机器可处理的形式把矩阵写入文件*/

int mod2sparse_write
( FILE *f,
  mod2sparse *m
)
{
  mod2entry *e;
  int i;

  intio_write(f,m->n_rows);
  if (ferror(f)) return 0;

  intio_write(f,m->n_cols);
  if (ferror(f)) return 0;

  for (i = 0; i<mod2sparse_rows(m); i++)
  { 
    e = mod2sparse_first_in_row(m,i);

    if (!mod2sparse_at_end(e))
    {
      intio_write (f, -(i+1));        //行数用负数表示
      if (ferror(f)) return 0;

      while (!mod2sparse_at_end(e))    
      { 
        intio_write (f, mod2sparse_col(e)+1);    //正数输出各列
        if (ferror(f)) return 0;

        e = mod2sparse_next_in_row(e);
      }
    }
  }

  intio_write(f,0);			//最后以0结束
  if (ferror(f)) return 0;

  return 1;
}


/* READ A SPARSE MOD2 MATRIX STORED IN MACHINE-READABLE FORM FROM A FILE. 从文件中读矩阵*/

mod2sparse *mod2sparse_read     
( FILE *f
)
{
  int n_rows, n_cols;
  mod2sparse *m;
  int v, row, col;

  n_rows = intio_read(f);
  if (feof(f) || ferror(f) || n_rows<=0) return 0;

  n_cols = intio_read(f);
  if (feof(f) || ferror(f) || n_cols<=0) return 0;

  m = mod2sparse_allocate(n_rows,n_cols);

  row = -1;

  for (;;)
  {
    v = intio_read(f);
    if (feof(f) || ferror(f)) break;

    if (v==0)      //读到0说明矩阵读完
    { return m;
    }
    else if (v<0)    //v<0说明读到是一行的标记
    { row = -v-1;
      if (row>=n_rows) break;   
    }
    else 
    { col = v-1;
      if (col>=n_cols) break;
      if (row==-1) break;
      mod2sparse_insert(m,row,col);   //在矩阵中插入实体
    }
  }

  /* Error if we get here. */

  mod2sparse_free(m);  //如果没有出错正常地读完整个矩阵会在for循环内结束该函数并返回m
  return 0;   
}


/* LOOK FOR AN ENTRY WITH GIVEN ROW AND COLUMN.在矩阵中找某位置的实体 */

mod2entry *mod2sparse_find
( mod2sparse *m,
  int row,
  int col
)
{ 
  mod2entry *re, *ce;

  if (row<0 || row>=mod2sparse_rows(m) || col<0 || col>=mod2sparse_cols(m))
  { fprintf(stderr,"mod2sparse_find: row or column index out of bounds\n");
    exit(1);     //如果行数和列数超出范围提示出错
  }

  /* Check last entries in row and column. */

  re = mod2sparse_last_in_row(m,row);   //re指向row行的最后一个实体 
  if (mod2sparse_at_end(re) || mod2sparse_col(re)<col) 
  { return 0;               //如果col大于最后一个实体所在的列则没有该实体
  }
  if (mod2sparse_col(re)==col) 
  { return re;             
  }

  ce = mod2sparse_last_in_col(m,col);    //使ce指向col列的最后一个实体
  if (mod2sparse_at_end(ce) || mod2sparse_row(ce)<row)
  { return 0;               //如果row大于最后一个实体所在的行则没有该实体     
  }
  if (mod2sparse_row(ce)==row)
  { return ce;
  }

  /* Search row and column in parallel, from the front. */

  re = mod2sparse_first_in_row(m,row);
  ce = mod2sparse_first_in_col(m,col);

  for (;;)     //并行在row行和col列查找该实体
  { 
    if (mod2sparse_at_end(re) || mod2sparse_col(re)>col)
    { return 0;
    } 
    if (mod2sparse_col(re)==col) 
    { return re;
    }

    if (mod2sparse_at_end(ce) || mod2sparse_row(ce)>row)
    { return 0;
    } 
    if (mod2sparse_row(ce)==row)
    { return ce;
    }

    re = mod2sparse_next_in_row(re);
    ce = mod2sparse_next_in_col(ce);
  }
}


/* INSERT AN ENTRY WITH GIVEN ROW AND COLUMN.在矩阵的某个位置插入实体 */

mod2entry *mod2sparse_insert
( mod2sparse *m,
  int row,
  int col
)
{
  mod2entry *re, *ce, *ne;

  if (row<0 || row>=mod2sparse_rows(m) || col<0 || col>=mod2sparse_cols(m))
  { fprintf(stderr,"mod2sparse_insert: row or column index out of bounds\n");