code6.c

快速傅立叶变换程序代码,学信号的同学,可要注意了

/*   code6 - software to create low-density parity-check codes

     Copyright   (c) 2002   David J.C. MacKay 

    This program is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation; either version 2 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program; if not, write to the Free Software
    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA

    GNU licenses are here :
    http://www.gnu.org/licenses/licenses.html

    Author contact details are here :
    http://www.inference.phy.cam.ac.uk/mackay/

    If you find this software useful, please feel free to make a donation to
    support David MacKay's research group.
*/
/*
   code6.c
   only change from code5 is that the first matrix has a variable size called
   N1.
   
   code5.c                                   (c) DJCM 96 01 17 from code4r.c

   makes sparse alist matrix with aribtrary REAL "t".

   The matrix is generated in its transpose, from bottom to top.

   You get to choose two average column weights, t1 and t2.
   Thus regular and mildly irregular (highly constrained) codes
   can be made.

   The invertible
   bit of the matrix is at the end (the bottom)

   But actually it's not invertible any more if CHEAP is used

   usage: 

   only args needed are N, M, t_square and t_rest

   code5 -n N -m M -seed seed -t1 5 (-t2 5)
               -o outfile (Default MNC4/N.t1.t2.seed)

   e.g.

   code5 -n 200 -m 100 -t1 2.5

   new option: staircase --- packs a t=2 staircase into the 
   first rows of the matrix

like this

11000000000000
01100000000000
00110000000000
00011000000000
00001100000000
00000110000000

Recommended for rate 1/2 codes.

  Usage

  code5 -stair 1 -n 2000 -m 1000 -t1 2 -t2 3 -seed 1

t1 = 2 is RECOMMENDED


If what is wanted is 
10000000000
10000000000
010000000000
01000000000
00100000000
0010000000

then this can be achieved using standard code5 generator 
eg
code5 -n 60000 -m 40000 -t1 3 -t2 2 -seed 6423

though maybe I should have used
code5 -n 60000 -m 40000 -t1 2 -t2 3 -seed 6423
to make the `2' bit go faster.

No currently by definition the "1" bit is a square matrix.
Except when stair==1 or 2, in which case it's differeeent

*/

#include "./ansi/r.h"
#include "./ansi/rand2.h" 
#include "./ansi/mynr.h"
#include "./ansi/cmatrix.h"

typedef struct { /* code_control */
#include "code5_var_str.h"
} code_creation ;


static void dc_defaults ( code_creation * ) ; 
static int    process_command ( int , char **  , code_creation *  ) ; 
static void   print_usage ( char ** , FILE * ,   code_creation *  );

static int make_sense ( code_creation * ) ; 
static int cheap_make_code ( code_creation * , alist_matrix * ) ;
static int report_code ( code_creation * , alist_matrix * ) ;
static int make_space ( code_creation * , alist_matrix * ) ; 
static int check_alist ( alist_matrix *  , code_creation * ) ;

void   main ( int , char ** ) ;

/*
        MAIN
                     */
void main ( int argc, char *argv[] )
{
  FILE   *fp  ; 

  code_creation dc   ;
  alist_matrix a ; 

  dc_defaults ( &dc ) ; 

  if ( process_command (argc, argv, &dc ) < 0 ) exit (0) ;
  if ( make_sense ( &dc ) < 0 ) exit (0) ;
  if ( make_space ( &dc , &a ) < 0 ) exit (0) ;

  fprintf(stderr,"code6 N=%d, M=%d t1=%f t2=%f\n",
	    dc.N , dc.M , dc.t1 , dc.t2  ) ;
  fflush(stderr);

  if ( 1 || dc.cheap ) { /* no choice now, it's cheap */
    if ( cheap_make_code (  &dc ,  &a ) < 0 ) { 
      printf ( "FAILED\n" ) ; exit ( 0 ) ;
    }
  }
  if ( check_alist ( &a , &dc ) < 0 ) { 
    if ( dc.cheap == 0 ) { report_code ( &dc , &a ) ; }
    printf ( "FAILED\n" ) ; exit ( 0 ) ;
  }

  if ( dc.out ) {
    fp = fopen ( dc.outfile , "w" ) ; 
    if ( !fp ) {
      fprintf ( stderr , "can't open %s\n" , dc.outfile ) , 
      exit (0) ; 
    }
    write_alist_transpose ( fp , &a ) ; 
    fclose ( fp ) ; 
  }
  free_alist ( &a ) ; 
  if ( dc.verbose ) {

    read_allocate_alist ( &a , dc.outfile ) ; 
    if ( dc.out ) {
      fp = fopen ( "tmpout" , "w" ) ; 
      if ( !fp ) {
	fprintf ( stderr , "can't open %s\n" , dc.outfile ) , 
	exit (0) ; 
      }
      write_alist ( fp , &a ) ; 
      fclose ( fp ) ; 
    }
    free_alist ( &a ) ; 

    printf ( "diff tmpout %s\n" , dc.outfile ) ; 
  }
}

static int make_sense ( code_creation *dc ) 
{ /* routine to correct silly control parameters */
  int status = 0 ; 
  char junk[100] ; 

  /* M is number of rows of the real matrix A and N is the number of cols.
     In this program the transpose is generated, row by row (from 1 to N)
     Confused? So am I.

     New convention for filename is N.M
*/

  dc->NS = dc->N - dc->M ; /*  signal block length, K -- not that we care */
  if ( dc->NS < 0 ) { fprintf ( stderr , "ERROR Need N > M \n" ) ; exit(0) ;}
  else if ( dc->NS == 0 )  { fprintf ( stderr , "warning N = M \n" ) ; }

  if ( dc->t2 == 0.0 ) dc->t2 = dc->t1 ;
  sprintf ( junk , "MNC%d%s/%d.%d.%3.1f.%ld" , 
	   4 , ( dc->cheap ? "RC" : "R" ) , dc->N , dc->M , dc->t1 ,
	   dc->seed 
	   ) ; 
  if ( dc->out==1 )  
    sprintf ( dc->outfile , "%s" ,  junk ) ; 
  printf ( "code will be written to %s\n" , dc->outfile ) ; 
  return status ; 
}

static int make_space ( code_creation *dc  , alist_matrix *a) {
  int status = 0 ; 
  int tr , tc , M , N ; 

  if ( dc->verbose ) printf ( "Making space\n" ) ; 

    N = dc->N ; M = dc->M ;  
    tc = (int)( MAX ( dc->t1 , dc->t2 ) ) + 5 ; 
    tr = (int) ( tc * ( (double) N / (double) M ) ) ;
  if ( dc->verbose ) 
    printf ( "for up to %d per col and %d per row\n" , tc , tr ) ; 

  initialize_alist ( a , M , N , tr , tc ) ;
/* this allocates the memory */
/* NB, this is transposed relative to the normal use of alist */

  return status ; 
}

/* calls uneven_sparse_rectangular_alist which is found in cmatrix.c */
static int cheap_make_code ( code_creation *dc , alist_matrix *a ) {
  int status = 0 ; 

  if ( dc->verbose ) printf ( "Making code\n" ) ; 
  ran_seed ( dc->seed ) ; 

  if ( !( dc->stair ) ) { /* the standard program. Modified in code6 from code5 version */
    status += uneven_sparse_rectangular_alist ( a , 0 , dc->N1 , dc->t1 , 
					       dc->M , dc->N , dc->hd , dc->maxattempts ) ;
    if ( dc->verbose )  { printf ( "done 1st matrix\n" ) ; fflush(stdout) ; }

    status += uneven_sparse_rectangular_alist ( a , dc->N1 , dc->N , dc->t2 , 
					       dc->M , dc->N , dc->hd , dc->maxattempts ) ;
    /* goes back to the last column for HD purposes, but only 
       generates from dc-> back */
    if ( dc->verbose )  { printf ( "done 2nd matrix\n" ) ; fflush(stdout) ; }

  }  else if (dc->stair==1) { /* make a staircase of height NS (K) */
/* NB, if t1=2 and NS > M then this is pointless for making good codes */
    status += staircase_alist ( a , 0 , dc->NS , (int)(dc->t1) , 
			       dc->M , dc->N ) ;

    status += uneven_sparse_rectangular_alist ( a , dc->NS , dc->N , dc->t2 , 
					       dc->M , dc->N , dc->hd , dc->maxattempts ) ;
  } else  { /* dc->stair == 2 */
    status += slope_alist ( a , 0 , dc->M / 2 , (int)(dc->t1) , 
			       dc->M , dc->N ) ;

    status += uneven_sparse_rectangular_alist ( a , dc->M/2 , dc->N , dc->t2 , 
					       dc->M , dc->N , dc->hd , dc->maxattempts ) ;
  }

  status += finish_off_alist ( a ) ;

  return status ; 
}

static int report_code ( code_creation *dc , alist_matrix *a ) {
  int status = 0 ; 

  report_alist ( a , dc->verbose + 1 ) ; 

  return status ; 
}

static int check_alist ( alist_matrix *alist ,  code_creation *dc )
{
  int status = 0 ; 
  int tmpi = 0 , n , m ;

  if ( dc->verbose ) printf ( "Checking alist\n" ) ; 

    /* check: every col should have more than 1 1 in it */
    for ( n = 1 ; n <= alist->N ; n ++ ) {
      if ( alist->num_nlist[n] <= 1 ) {
	tmpi ++ ; 
	fprintf ( stderr , "%d:" , n ) ; 
      }
    }
    if ( tmpi > 0 ) {
      fprintf ( stderr , "Warning unchecked bits n : %d\n", tmpi ) ;
      status -- ; 
    }

    /* check: every row should have more than 1 1 in it */

    tmpi = 0 ; 
    for ( m = 1 ; m <= alist->M ; m ++ )  {
      if ( alist->num_mlist[m] <= 1 ) {
	tmpi ++ ; 
	fprintf ( stderr , "%d:" , m ) ; 
      }
    }
    if ( tmpi > 0 ) {
      fprintf ( stderr , "Warning unchecked bits n : %d\n", tmpi ) ;
      status -= 10 ; 
    }

  if ( status < 0 && dc->verbose ) {
    fprintf ( stderr , "Failing ... here is the offending alist\n" ) ;
    report_alist ( alist , 2 ) ; 
  }
  else if ( dc->verbose >= 2 )   {
    report_alist ( alist , dc->verbose ) ; 
  }
  return status ; 
}

static void dc_defaults ( code_creation *dc ) 
{
#include "code5_var_def.c"
}

static int process_command ( int argc , char **argv , code_creation *dc ) {

  int p_usage = 0 ;
  int status = 0 ;
  int cs , i , fake ;

  if ( argc < 1 )     {
    p_usage = 1 ; 
    status -- ;
  }

#define ERROR1 fprintf ( stderr , "arg to `%s' missing\n" , argv[i] ) ; \
               status --
#define ERROR2 fprintf ( stderr , "args to `%s' missing\n" , argv[i] ) ; \
               status --
  for (i = 1 ; i < argc; i++)    {
    cs = 1 ;
    if ( strcmp (argv[i], "-fake") == 0 ) {
      if ( i + 1 == argc ) { ERROR1;      }
      else cs *= sscanf(argv[++i], "%d", &(fake)); 
    }  
#include "code5_var_clr.c"
    else {
      fprintf ( stderr , "arg `%s' not recognised\n" , argv[i] ) ; 
      p_usage = 1 ;
      status -- ;
    }
    if ( cs == 0 ) {
      fprintf ( stderr , "arg at or before `%s' has incorrect format\n" , 
	       argv[i] ) ;
      p_usage = 1 ;
      status -- ;
    }
  }
  if ( p_usage ) print_usage ( argv , stderr , dc ) ;
  return ( status ) ;
}
#undef ERROR1
#undef ERROR2

#define DNT fprintf( fp, "\n        ")
#define NLNE  fprintf( fp, "\n")

static void print_usage ( char **argv , FILE * fp ,
			  code_creation *dc
)
{
  fprintf( fp, "Usage: %s ",argv[0]);
  fprintf( fp, " [optional arguments]");

#include "code5_var_usg.c"
  return ;
}
#undef DNT
#undef NLNE


/*
<!-- hhmts start -->
Last modified: Sun Mar 24 21:38:05 1996
<!-- hhmts end -->
*/