tcc_dvb-rcs_pc.c

MAP 算法 源程序 第三代 数字通信 turbo码相关

    printf("the e_ber was %1.1e the f_ber was %1.1e \n ", e_ber, f_ber);
    fprintf(fp,"the e_ber was %1.1e the f_ber was %1.1e\n ", e_ber, f_ber);

    fclose (fp);
  } /* end of eb_ovr_n0 loop */
}

float max3(float x, float y, float z)
{
  float maxmum;

  maxmum = x>y?x:y;
  maxmum = maxmum>z?maxmum:z;

  return(maxmum);
}

/*******************/
/* 0/1 Generation  */
/*******************/
void gen01dat( int data_len, int *out_array )
{
  int t;                         /* time */

  rand();                        /* re-seed the random number generator */

  for (t = 0; t < data_len; t++) /* generate random data, write to output array */
    *( out_array + t ) = (int)( rand() / (RAND_MAX / 2) > 0.5 );
}

/***************************/
/* PCCC Turbo code encoder */
/***************************/
void nbrsc_encd(int g[3][K], int *in_array, int *y_array, int *w_array, int s0[K-1], int sN[K-1]);
void deci2bin(int d, int size, int *b);
int bin2deci(int *b, int size);
int start_end_state(int sNdeci);

void turbo_encd( int g[3][K], int M, int *alpha, int *in_array, int *turbo_en_out)
{
  int i,j,temp_data;                 /* loop variables, temporary value */
  int y[2][N];                       /* make matrixs to store the medium output */
  int *pinfo_array;                  /* interleaved input message data */
  int *Y1out_array,*W1out_array;     /* encoder C1 output data */
  int *Y2out_array,*W2out_array;     /* encoder C2 output data */
  int s0deci,sNdeci;                 /* state's number in decimal */
  int s0[K-1],sN[K-1];               /* state's number in binary */

  pinfo_array = malloc( (2*N) * sizeof(int) );
  if (pinfo_array == NULL) {
    printf("\n turbo_encd:  error allocating interleaved info array, aborting!");
    exit(1);
  }

  Y1out_array = malloc( N * sizeof(int) );
  if (Y1out_array == NULL) {
    printf("\n turbo_encd:  error allocating Y array of encoder one, aborting!");
    exit(1);
  }

  W1out_array = malloc( N * sizeof(int) );
  if (W1out_array == NULL) {
    printf("\n turbo_encd:  error allocating W array of encoder one, aborting!");
    exit(1);
  }

  Y2out_array = malloc( N * sizeof(int) );
  if (Y2out_array == NULL) {
    printf("\n turbo_encd:  error allocating Y array of encoder two, aborting!");
    exit(1);
  }

  W2out_array = malloc( N * sizeof(int) );
  if (W2out_array == NULL) {
    printf("\n turbo_encd:  error allocating W array of encoder two, aborting!");
    exit(1);
  }

  s0deci = 0;
  deci2bin(s0deci, K-1, s0);

  nbrsc_encd(g, in_array, Y1out_array, W1out_array, s0, sN);  /* Test the trellis to get the initial state */

  sNdeci = bin2deci(sN, K-1);
  s0deci = start_end_state(sNdeci);
  deci2bin(s0deci, K-1, s0);

  nbrsc_encd(g, in_array, Y1out_array, W1out_array, s0, sN);  /* encoder C1 output */

  /* make a matrix with first and second row corresponding to info sequence */
  /* third and fifth row corresponding to RSC #1 s check/parity bits */
  /* fourth and sixth row corresponding to RSC #2 s check/parity bits */

  for (j=0; j<2*N; j++){         /* get the systematic information array */
    *(pinfo_array + j) = *(in_array + j);
    *(turbo_en_out + j) = *(in_array + j);
  }

  for (j=0; j<N; j++){           /* turbo code permutation Level one */
    if  (j%2 == 0){
      temp_data = *(pinfo_array + 2*j);
      *(pinfo_array + 2*j) = *(pinfo_array + 2*j + 1);
      *(pinfo_array + 2*j + 1) = temp_data;
    }
  }

  for (i=0; i<2; i++){           /* turbo code permutation Level two */
    for (j=0; j<N; j++)
      y[i][j] = *(pinfo_array + 2*j + i);
  }

  for (i=0; i<2; i++){
    for (j=0; j<N; j++)
      *(pinfo_array + 2*j + i) = y[i][*(alpha + j)];
  }

  s0deci = 0;
  deci2bin(s0deci, K-1, s0);

  nbrsc_encd(g, pinfo_array, Y2out_array, W2out_array, s0, sN);  /* Test the trellis to get the initial state */

  sNdeci = bin2deci(sN, K-1);
  s0deci = start_end_state(sNdeci);
  deci2bin(s0deci, K-1, s0);

  nbrsc_encd(g, pinfo_array, Y2out_array, W2out_array, s0, sN);  /* encoder C2 output */

  #if R == 1              /* code rate = 1/3 */
  for (j=0; j<N; j++){
    *(turbo_en_out + 2*j + 2*N) = *(Y1out_array + j);
    *(turbo_en_out + 2*j + 1 + 2*N) = *(Y2out_array + j);
    *(turbo_en_out + 2*j + 4*N) = *(W1out_array + j);
    *(turbo_en_out + 2*j + 1 + 4*N) = *(W2out_array + j);
  }
  #endif

  #if R == 2             /* code rate = 2/5 */
  for (j=0; j<N; j++){
    *(turbo_en_out + 2*j + 2*N) = *(Y1out_array + j);
    *(turbo_en_out + 2*j + 1 + 2*N) = *(Y2out_array + j);
  }
  for (j=0; j<M; j++){
    *(turbo_en_out + 2*j + 4*N) = *(W1out_array + 2*j);
    *(turbo_en_out + 2*j + 1 + 4*N) = *(W2out_array + 2*j);
  } 
  #endif

  #if R == 3             /* code rate = 1/2 */
  for (j=0; j<N; j++){
    *(turbo_en_out + 2*j + 2*N) = *(Y1out_array + j);
    *(turbo_en_out + 2*j + 1 + 2*N) = *(Y2out_array + j);
  }
  #endif

  #if R == 4             /* code rate = 2/3 */
  for (j=0; j<M; j++){
    *(turbo_en_out + 2*j + 2*N) = *(Y1out_array + 2*j);
    *(turbo_en_out + 2*j + 1 + 2*N) = *(Y2out_array + 2*j);
  } 
  #endif

  #if R == 5             /* code rate = 3/4 */
  for (j=0; j<M; j++){
    *(turbo_en_out + 2*j + 2*N) = *(Y1out_array + 3*j);
    *(turbo_en_out + 2*j + 1 + 2*N) = *(Y2out_array + 3*j);
  }
  #endif 

  #if R == 6             /* code rate = 4/5 */
  for (j=0; j<M; j++){
    *(turbo_en_out + 2*j + 2*N) = *(Y1out_array + 4*j);
    *(turbo_en_out + 2*j + 1 + 2*N) = *(Y2out_array + 4*j);
  }
  #endif 

  #if R == 7             /* code rate = 6/7 */
  for (j=0; j<M; j++){
    *(turbo_en_out + 2*j + 2*N) = *(Y1out_array + 6*j);
    *(turbo_en_out + 2*j + 1 + 2*N) = *(Y2out_array + 6*j);
  }
  #endif

  free(pinfo_array);
  free(Y1out_array);
  free(Y2out_array);
  free(W1out_array);
  free(W2out_array);
}

void nbrsc_encd(int g[3][K], int *in_array, int *y_array, int *w_array, int s0[K-1], int sN[K-1])
{
  int m;                   /* K-1 */
  int i,j;                 /* loop variables */
  int t,tt;                /* bit time, symble time */
  int *a_array,*b_array;   /* pointer to input data array */
  int shift_reg[K];        /* the 1st RSC encoder shift register */
  int sr_head;             /* index to the first element in the shift register */
  int x,y,w;               /* feedback bit and output bits */

  m = K-1;
  a_array = malloc( N * sizeof(int) ); /* allocate space for A input data array */
  if (a_array == NULL) {
    printf("\n nbrsc_encd:  error allocating A array, aborting!");
    exit(1);
  }

  b_array = malloc( N * sizeof(int) ); /* allocate space for B input data array */
  if (b_array == NULL) {
    printf("\n nbrsc_encd:  error allocating B array, aborting!");
    exit(1);
  }

  for (t=0; t<N; t++){     /* read in the data and store it in the array */
    *(a_array+t) = *(in_array+2*t);
    *(b_array+t) = *(in_array+2*t+1);
  }

  shift_reg[0] = 0;        /* Initialize the shift register */
  for (j=0; j<K-1; j++)
    shift_reg[j+1] = s0[j];

  sr_head = 0;             /* Initialize the index of the shift register */
  tt = 0;                  /* Initialize the channel symbol output index */

  for (t=0; t<N; t++){     /* Start the encoding process */
    shift_reg[sr_head] = *(a_array+t);
    x = 0;
    y = 0;
    w = 0;

    for (j=0; j<K; j++){   /* get feedback bit */
      i = (j+sr_head) % K;
      x ^= shift_reg[i] & g[0][j];
    }

    shift_reg[sr_head] = x ^ *(b_array + t);
    for (j=0; j<K; j++){
      i = (j+sr_head) % K;
      y ^= shift_reg[i] & g[1][j];
      w ^= shift_reg[i] & g[2][j];
    }

    for (j=1; j<K-1; j++){
      i = (j+sr_head) % K;
      shift_reg[i] = shift_reg[i] ^ *(b_array+t);
    }

    *(y_array+tt) = y;
    *(w_array+tt) = w;
    tt = tt + 1;

    sr_head -= 1;
    if ( sr_head<0 ) sr_head = m;   /* make sure we adjust pointer modulo K */
  }

  for (j=0; j<(K-1); j++)
    sN[j] = shift_reg[j+1];

  free (a_array);         /* free the dynamically allocated array */
  free (b_array);         /* free the dynamically allocated array */
}

/******************************************/
/* Circulation state correspondence table */
/******************************************/
int start_end_state(int sNdeci){
  int s0deci;

  s0deci = 0;

  #if (N%7 == 1)          /* Frame size N modulo 7 = 1 */
  if (sNdeci == 0)
    s0deci = 0;
  else if (sNdeci == 1)
    s0deci = 6;
  else if (sNdeci == 2)
    s0deci = 4;
  else if (sNdeci == 3)
    s0deci = 2;
  else if (sNdeci == 4)
    s0deci = 7;
  else if (sNdeci == 5)
    s0deci = 1;
  else if (sNdeci == 6)
    s0deci = 3;
  else if (sNdeci == 7)
    s0deci = 5;
  #endif

  #if (N%7 == 2)          /* Frame size N modulo 7 = 2 */
  if (sNdeci == 0)
    s0deci = 0;
  else if (sNdeci == 1)
    s0deci = 3;
  else if (sNdeci == 2)
    s0deci = 7;
  else if (sNdeci == 3)
    s0deci = 4;
  else if (sNdeci == 4)
    s0deci = 5;
  else if (sNdeci == 5)
    s0deci = 6;
  else if (sNdeci == 6)
    s0deci = 2;
  else if (sNdeci == 7)
    s0deci = 1;
  #endif

  #if (N%7 == 3)          /* Frame size N modulo 7 = 3 */
  if (sNdeci == 0)
    s0deci = 0;
  else if (sNdeci == 1)
    s0deci = 5;
  else if (sNdeci == 2)
    s0deci = 3;
  else if (sNdeci == 3)
    s0deci = 6;
  else if (sNdeci == 4)
    s0deci = 2;
  else if (sNdeci == 5)
    s0deci = 7;
  else if (sNdeci == 6)
    s0deci = 1;
  else if (sNdeci == 7)
    s0deci = 4;
  #endif

  #if (N%7 == 4)          /* Frame size N modulo 7 = 4 */
  if (sNdeci == 0)
    s0deci = 0;
  else if (sNdeci == 1)
    s0deci = 4;
  else if (sNdeci == 2)
    s0deci = 1;
  else if (sNdeci == 3)
    s0deci = 5;
  else if (sNdeci == 4)
    s0deci = 6;
  else if (sNdeci == 5)
    s0deci = 2;
  else if (sNdeci == 6)
    s0deci = 7;
  else if (sNdeci == 7)
    s0deci = 3;
  #endif

  #if (N%7 == 5)          /* Frame size N modulo 7 = 5 */
  if (sNdeci == 0)
    s0deci = 0;
  else if (sNdeci == 1)
    s0deci = 2;
  else if (sNdeci == 2)
    s0deci = 5;
  else if (sNdeci == 3)
    s0deci = 7;
  else if (sNdeci == 4)
    s0deci = 1;
  else if (sNdeci == 5)
    s0deci = 3;
  else if (sNdeci == 6)
    s0deci = 4;
  else if (sNdeci == 7)
    s0deci = 6;
  #endif

  #if (N%7 == 6)          /* Frame size N modulo 7 = 6 */
  if (sNdeci == 0)
    s0deci = 0;
  else if (sNdeci == 1)
    s0deci = 7;
  else if (sNdeci == 2)
    s0deci = 6;
  else if (sNdeci == 3)
    s0deci = 1;
  else if (sNdeci == 4)
    s0deci = 3;
  else if (sNdeci == 5)
    s0deci = 4;
  else if (sNdeci == 6)
    s0deci = 5;
  else if (sNdeci == 7)
    s0deci = 2;
  #endif

  return(s0deci);
}

/********************************************/
/* QPSK simulator for non-binary turbo code */
/********************************************/
void gngauss(float mean, float sigma, int encd_block_len, float *gauss1, float *gauss2);

void addnoise(float eb_ovr_n0, int encd_block_len, float rate, int *in_array, float *out_array)
{
  int t;
  float mean, eb, sn_ratio, sigma, signal_1, signal_2;
  float *gauss1, *gauss2;

  mean = 0;
  eb = 1;
  sn_ratio = (float) pow(10, ( eb_ovr_n0 / 10) );
  sigma =  (float) sqrt (eb / (2 * sn_ratio * rate) );

  gauss1 = malloc( encd_block_len * sizeof(float) );
  if (gauss1 == NULL) {
    printf("\n addnoise:  error allocating cos I noise array, aborting!");
    exit(1);
  }

  gauss2 = malloc( encd_block_len * sizeof(float) );
  if (gauss2 == NULL) {
    printf("\n addnoise:  error allocating sin Q noise array, aborting!");
    exit(1);
  }

  gngauss(mean, sigma, encd_block_len, gauss1, gauss2);

  for (t = 0; t < encd_block_len; t++) {       /* transform the data from 0/1 to +1/-1 and add noise */
    signal_1 = 1 - 2 * *( in_array + 2*t );
    *( out_array + 2*t ) = signal_1 + *(gauss1+t);
    signal_2 = 1 - 2 * *( in_array + 2*t + 1 );
    *( out_array + 2*t + 1) = signal_2 + *(gauss2+t);
  }

  free(gauss1);
  free(gauss2);
}

void gngauss(float mean, float sigma, int encd_block_len, float *gauss1, float *gauss2) 
{ 
  int i;
  double u, r;                         /* uniform and Rayleigh random variables */

  for (i=0; i<encd_block_len; i++){
    u = (double)rand() / RAND_MAX;     /* uniformly distributed random number u between 0 and 1 - 1E-6*/
    if (u == 1.0) u = 0.999999999;

    r = sigma * sqrt( 2.0 * log( 1.0 / (1.0 - u) ) ); /* generate a Rayleigh-distributed random number r using u */

    u = (double)rand() / RAND_MAX;             /* generate another uniformly-distributed random number u as before*/
    if (u == 1.0) u = 0.999999999;

    *(gauss1+i) = mean + r * cos(2 * PI * u);  /* generate two Gaussian-distributed random number using r and u */
    *(gauss2+i) = mean + r * sin(2 * PI * u);
  } 
}

/**************************************************/
/* Serial to parallel demultiplex at the receiver */
/**************************************************/
void de_multiplex(int msg_length, int M, int *alpha, float *in_array, float *sys_array,
		  float *pinfo_array, float *parity1_array, float *parity2_array)
{
  int i,j;                    /* loop variables */
  float permu[2][N]={0.0};    /* make a matrix store systematic info array */
  float yy[2][N]={0.0};       /* make matrix to store first parity info of encoder1&2 temporarily */
  float ww[2][N]={0.0};       /* make matrix to store second parity info of encoder1&2 temporarily */
  float temp_data;            /* temporary value */

  for (j=0; j<msg_length; j++){       /* get systematic information array */
    *(sys_array + j) = *(in_array + j);
    *(pinfo_array + j) = *(in_array + j);
  }

  for (j=0; j<N; j++){                /* turbo code permutation Level one */
    if  (j%2 == 0){
      temp_data = *(pinfo_array + 2*j);
      *(pinfo_array + 2*j) = *(pinfo_array + 2*j + 1);
      *(pinfo_array + 2*j + 1) = temp_data;