/* MAKE-LDPC.C - Make a Low Density Parity Check code's parity check matrix. *//* Copyright (c) 2000, 2001, 2006 by Radford M. Neal and Peter Junteng Liu * * Permission is granted for anyone to copy, use, modify, or distribute this * program and accompanying programs and documents for any purpose, provided  * this copyright notice is retained and prominently displayed, along with * a note saying that the original programs are available from Radford Neal's * web page, and note is made of any changes made to the programs.  The * programs and documents are distributed without any warranty, express or * implied.  As the programs were written for research purposes only, they have * not been tested to the degree that would be advisable in any important * application.  All use of these programs is entirely at the user's own risk. */#include <stdio.h>#include <stdlib.h>#include <string.h>#include <math.h>#include "rand.h"#include "alloc.h"#include "intio.h"#include "open.h"#include "mod2sparse.h"#include "mod2dense.h"#include "mod2convert.h"#include "rcode.h"#include "distrib.h"/* METHODS FOR CONSTRUCTING CODES. */typedef enum { Evencol, 	/* Uniform number of bits per column, with number specified */  Evenboth 	/* Uniform (as possible) over both columns and rows */} make_method; void make_ldpc (int, make_method, distrib *, int);int *column_partition (distrib *, int);void usage (void);/* MAIN PROGRAM. */int main( int argc,  char **argv){  make_method method;  char *file, **meth;  int seed, no4cycle;  distrib *d;  char junk;  FILE *f;    /* Look at initial arguments. */  if (!(file = argv[1])   || !argv[2] || sscanf(argv[2],"%d%c",&M,&junk)!=1 || M<=0   || !argv[3] || sscanf(argv[3],"%d%c",&N,&junk)!=1 || N<=0   || !argv[4] || sscanf(argv[4],"%d%c",&seed,&junk)!=1)  { usage();  }  /* Look at the arguments specifying the method for producing the code. */  meth = argv+5;  if (!meth[0]) usage();  no4cycle = 0;  if (strcmp(meth[0],"evencol")==0 || strcmp(meth[0],"evenboth")==0)  { method = strcmp(meth[0],"evencol")==0 ? Evencol : Evenboth;    if (!meth[1])    { usage();    }    d = distrib_create(meth[1]);    if (d==0)    { usage();    }    if (meth[2])    { if (strcmp(meth[2],"no4cycle")==0)      { no4cycle = 1;        if (meth[3])        { usage();        }      }      else      { usage();      }    }  }  else  { usage();  }  /* Check for some problems. */  if (distrib_max(d)>M)  { fprintf(stderr,      "At least one checks per bit (%d) is greater than total checks (%d)\n",      distrib_max(d), M);    exit(1);  }  if (distrib_max(d)==M && N>1 && no4cycle)  { fprintf(stderr,      "Can't eliminate cycles of length four with this many checks per bit\n");    exit(1);  }   /* Make the parity check matrix. */  make_ldpc(seed,method,d,no4cycle);  /* Write out the parity check matrix. */  f = open_file_std(file,"wb");  if (f==NULL)   { fprintf(stderr,"Can't create parity check file: %s\n",file);    exit(1);  }  intio_write(f,('P'<<8)+0x80);    if (ferror(f) || !mod2sparse_write(f,H) || fclose(f)!=0)  { fprintf(stderr,"Error writing to parity check file %s\n",file);    exit(1);  }  return 0;}/* PRINT USAGE MESSAGE AND EXIT. */void usage(void){ fprintf(stderr,"Usage:  make-ldpc pchk-file n-checks n-bits seed method\n");  fprintf(stderr,"Method: evencol  checks-per-col [ \"no4cycle\" ]\n");  fprintf(stderr,"    or: evencol  checks-distribution [ \"no4cycle\" ]\n");  fprintf(stderr,"    or: evenboth checks-per-col [ \"no4cycle\" ]\n");  fprintf(stderr,"    or: evenboth checks-distribution [ \"no4cycle\" ]\n");  exit(1);}/* CREATE A SPARSE PARITY-CHECK MATRIX.  Of size M by N, stored in H. */void make_ldpc( int seed,		/* Random number seed */  make_method method,	/* How to make it */  distrib *d,		/* Distribution list specified */  int no4cycle		/* Eliminate cycles of length four? */){  mod2entry *e, *f, *g, *h;  int added, uneven, elim4, all_even, n_full, left;  int i, j, k, t, z, cb_N;  int *part, *u;  rand_seed(10*seed+1);  H = mod2sparse_allocate(M,N);  part = column_partition(d,N);  /* Create the initial version of the parity check matrix. */  switch (method)  {     case Evencol:    {       z = 0;      left = part[z];      for (j = 0; j<N; j++)      { while (left==0)        { z += 1;          if (z>distrib_size(d))          { abort();          }          left = part[z];        }        for (k = 0; k<distrib_num(d,z); k++)        { do          { i = rand_int(M);          } while (mod2sparse_find(H,i,j));          mod2sparse_insert(H,i,j);        }        left -= 1;      }      break;    }    case Evenboth:    {      cb_N = 0;      for (z = 0; z<distrib_size(d); z++)      { cb_N += distrib_num(d,z) * part[z];      }            u = chk_alloc (cb_N, sizeof *u);      for (k = cb_N-1; k>=0; k--)      { u[k] = k%M;      }        uneven = 0;      t = 0;      z = 0;      left = part[z];      for (j = 0; j<N; j++)      {         while (left==0)        { z += 1;          if (z>distrib_size(d))          { abort();          }          left = part[z];        }	for (k = 0; k<distrib_num(d,z); k++)        {           for (i = t; i<cb_N && mod2sparse_find(H,u[i],j); i++) ;          if (i==cb_N)          { uneven += 1;            do            { i = rand_int(M);            } while (mod2sparse_find(H,i,j));            mod2sparse_insert(H,i,j);          }          else          { do            { i = t + rand_int(cb_N-t);            } while (mod2sparse_find(H,u[i],j));            mod2sparse_insert(H,u[i],j);            u[i] = u[t];            t += 1;          }        }        left -= 1;      }      if (uneven>0)      { fprintf(stderr,"Had to place %d checks in rows unevenly\n",uneven);      }      break;    }    default: abort();  }  /* Add extra bits to avoid rows with less than two checks. */  added = 0;  for (i = 0; i<M; i++)  { e = mod2sparse_first_in_row(H,i);    if (mod2sparse_at_end(e))    { j = rand_int(N);      e = mod2sparse_insert(H,i,j);      added += 1;    }    e = mod2sparse_first_in_row(H,i);    if (mod2sparse_at_end(mod2sparse_next_in_row(e)) && N>1)    { do       { j = rand_int(N);       } while (j==mod2sparse_col(e));      mod2sparse_insert(H,i,j);      added += 1;    }  }  if (added>0)  { fprintf(stderr,           "Added %d extra bit-checks to make row counts at least two\n",           added);  }  /* Add extra bits to try to avoid problems with even column counts. */  n_full = 0;  all_even = 1;  for (z = 0; z<distrib_size(d); z++)  { if (distrib_num(d,z)==M)     { n_full += part[z];    }    if (distrib_num(d,z)%2==1)    { all_even = 0;    }  }  if (all_even && N-n_full>1 && added<2)  { int a;    for (a = 0; added+a<2; a++)    { do      { i = rand_int(M);        j = rand_int(N);      } while (mod2sparse_find(H,i,j));      mod2sparse_insert(H,i,j);    }    fprintf(stderr, "Added %d extra bit-checks to try to avoid problems from even column counts\n",      a);  }  /* Eliminate cycles of length four, if asked, and if possible. */  if (no4cycle)  {     elim4 = 0;    for (t = 0; t<10; t++)     { k = 0;      for (j = 0; j<N; j++)      { for (e = mod2sparse_first_in_col(H,j);             !mod2sparse_at_end(e);             e = mod2sparse_next_in_col(e))        { for (f = mod2sparse_first_in_row(H,mod2sparse_row(e));               !mod2sparse_at_end(f);               f = mod2sparse_next_in_row(f))          { if (f==e) continue;            for (g = mod2sparse_first_in_col(H,mod2sparse_col(f));                 !mod2sparse_at_end(g);                 g = mod2sparse_next_in_col(g))            { if (g==f) continue;              for (h = mod2sparse_first_in_row(H,mod2sparse_row(g));                   !mod2sparse_at_end(h);                   h = mod2sparse_next_in_row(h))              { if (mod2sparse_col(h)==j)                { do                  { i = rand_int(M);                  } while (mod2sparse_find(H,i,j));                  mod2sparse_delete(H,e);                  mod2sparse_insert(H,i,j);                  elim4 += 1;                  k += 1;                  goto nextj;                }              }            }          }        }      nextj: ;      }      if (k==0) break;    }    if (elim4>0)    { fprintf(stderr,        "Eliminated %d cycles of length four by moving checks within column\n",         elim4);    }    if (t==10)     { fprintf(stderr,        "Couldn't eliminate all cycles of length four in 10 passes\n");    }  }}/* PARTITION THE COLUMNS ACCORDING TO THE SPECIFIED PROPORTIONS.  It   may not be possible to do this exactly.  Returns a pointer to an   array of integers containing the numbers of columns corresponding    to the entries in the distribution passed. */int *column_partition( distrib *d,		/* List of proportions and number of check-bits */  int n			/* Total number of columns to partition */){  double *trunc;  int *part;  int cur, used;  int i, j;  trunc = chk_alloc (distrib_size(d), sizeof(double));  part = chk_alloc (distrib_size(d), sizeof(int));  used = 0;  for (i = 0; i<distrib_size(d); i++)  { cur = floor(distrib_prop(d,i)*n);    part[i] = cur;     trunc[i] = distrib_prop(d,i)*n - cur;     used += cur;   }  if (used>n)   { abort();  }    while (used<n)  { cur = 0;    for (j = 1; j<distrib_size(d); j++)     { if (trunc[j]>trunc[cur])      { cur = j;      }    }    part[cur] += 1;    used += 1;    trunc[cur] = -1;  }  free(trunc);  return part;}