/* MAKE-LDPC.C - Make a Low Density Parity Check code's parity check matrix. */

/* 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.
 */

//#include <stdafx.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <time.h>

#include "rand.h"
#include "alloc.h"
#include "mod2sparse.h"
#include "rcode.h"
#include "make-ldpc.h"
#include "make-gen.h"
#include "encode.h"
#include "transmit.h"
#include "dec.h"
#include "check.h"
#include "extract.h"




seed = 15963;	//initiation number for random
block_num=1000;		//number of blocks to tansmitting要傳送的分組數
max_iter=200;   //最大迭代次數  
M=252;N=512;   //校驗矩陣的維數


/* METHODS FOR CONSTRUCTING CODES. */

/* MAIN PROGRAM. */
void main()
{
	int i;

	mod2sparse_strategy strategy;  //枚舉型變量

	/*parameter for decode解碼參數 */
	int iter,c;   //iter迭代次數，c用于標實解得的碼字是否滿足校驗矩陣
	double tot_iter=0;
	int nb;     //傳送的碼塊數
	double *lratio,*LLR;   /*demodulator output解調輸出*/

	FILE *stream;	/* point to "simu_report.txt" which is used to record the results */
	int err_bit,totbit;
 
	clock_t start, end;
	
	/* = = = = parameter  for channel = = = = = */
	float Eb_N0_dB[7] = {0.5,1,1.5,2,2.5,3};  //信道的幾種信噪比的情況
	float EbN0dB;
	int ien;

	double sgma;// = 1/sqrt(pow( 10.0, Eb_N0_dB/10.0 )*2*(N-M)/N);	//Eb / N0 = 1 / 2Rs2
	
	int *source_data, *encoded_data, *decoded_data, *decoded_source,*pchk;
	double *received_data;
	int err_num,err_blck;  //錯碼數和錯分組數
	double blck_err_prop[7],bit_err_prop[7];  //對應于每種信道參數的誤碼率和誤分組率

	/* = = = =生成 LDPC 碼的輸入參數 = = = = = = */

	char *pchk_file="pchk_file",*gen_file="gen_file";
	double  *lamita;	//distribution of check bits in each col
	double  *lao;		//lao[i]/i is distribution of check bits in each row

    int Dlmax = 10,Drmax = 10;			/* max number of check bits in rows */
  			/* max number of check bits in cols */


   lamita = chk_alloc (Dlmax+1, sizeof *lamita);
    for (i=0; i<Dlmax+1; i++)
       lamita[i] = 0;

    lao = chk_alloc (Drmax+1, sizeof *lao);
    for (i=0; i<Drmax+1; i++)
       lao[i] = 0;

	lamita[3] = 1;

	lao[6] = 1;		//intiation the distribution of check bit in each col
   

   /*====*/


	
	received_data = chk_alloc (N, sizeof *received_data);
 	
	source_data = chk_alloc (N-M, sizeof *source_data);
	encoded_data = chk_alloc (N, sizeof *encoded_data);
	decoded_data = chk_alloc (N, sizeof *decoded_data);
	pchk = chk_alloc (M, sizeof *pchk);
	decoded_source = chk_alloc (N-M, sizeof *decoded_source);
	lratio = chk_alloc (N, sizeof *lratio);
	LLR = chk_alloc (N, sizeof *LLR);


	if ((stream=fopen( "test_ldpc_report.txt", "a+" )) == NULL)
    {
	  printf("\nError! Can not open file simu_report.txt\n");
	  exit(1);
    }


	strategy = Mod2sparse_minprod;
	
	/*print the ldpc code structure to screen*/

	printf("/*=================================*/\n");
	printf("LDPC code structure\n");

	printf("M=%d, N=%d, blck_num=%d\n", M, N, block_num);
	

	/*print the ldpc code structure to simu_report.txt*/

	fprintf(stream,"/*=================================*/\n");
	fprintf(stream,"LDPC code structure KTS\n");

	fprintf(stream,"M=%d, N=%d, blck_num=%d\n", M, N, block_num);

    /*cycle for different signal to noise ratio*/

	for (ien=0;ien<6;ien++)
	{
		EbN0dB = (float)Eb_N0_dB[ien];
		sgma = 1/sqrt(pow( 10.0, EbN0dB/10.0 )*2*(N-M)/N);
		err_num=0,err_blck=0,totbit=0;

		/*====   output the simulation parameters	=====*/
		/*	to screen	*/
		printf("/*=================================*/\n");
		printf("LDPC code simulation:\n");
		printf("/*=================================*/\n");
		printf("====signal to noise ratio====\n");
		printf("Eb_N0_dB=%f, sgma=%f\n",EbN0dB, sgma);

		/*encode and decode for each block*/
		for (nb=0; nb<block_num;nb++)
		{

		    start=clock();
		    err_bit=0;

            make_ldpc_l(pchk_file,lamita, Dlmax,lao, Drmax, seed);
	
	        make_gen(pchk_file,gen_file, strategy);
		

         rand01(source_data, N-M,seed);

		 encode(pchk_file,gen_file, source_data,encoded_data);
		

		 /*modulation mode in communication system*/
		 bpsk_transmit(encoded_data, received_data, N ,seed, sgma);
		 
		 for (i=0;i<N;i++)
		 {
		    lratio[i] = (2*received_data[i]/(sgma*sgma)) ;//log(P(x=1)/P(x=0))
		 }
	
		 /* Read parity check matrix. */

         read_pchk(pchk_file);

		 initprp(H,lratio,decoded_data);   //初始化概率解碼
		   
		 for (iter = 0; ; iter++)
		 {
			   c = check(H,decoded_data,pchk);

				   
               if (iter==max_iter || (max_iter>0 && c==0))
			   { break; 
			   }

		    iterprp(H, lratio, decoded_data, LLR);
		}

   		   mod2sparse_free(H);

		   extract(gen_file, decoded_data, decoded_source);
		   err_bit=0;
		   /* = = = statistic error= = = */
		   for (i = 0; i<(N-M); i++)
		   {
			   err_bit += decoded_source[i] ^ source_data[i];

		   }
		   if (err_bit!=0)
		   {
			   err_num +=err_bit;
		       err_blck ++;
		       err_bit = 0;
		   }
		   printf("%d %s", nb+1, "block transmitted\n");
		   end=clock();
		   printf( "\n%5f seconds\n", (double)(end - start) / CLOCKS_PER_SEC);
		   totbit+=N-M;
		}

		blck_err_prop[ien] = (double) err_blck /(nb+1);
	    bit_err_prop[ien] = (double) err_num / totbit;
        fprintf(stream,"perr_num=%d,  err_blck=%d,  bit_err_prop=%f,  blck_err_pro=%f\n",
 		  err_num, err_blck, bit_err_prop[ien], blck_err_prop[ien]);

	}

	fclose(stream);

    free(received_data);
    free(source_data);
	free(encoded_data);
	free(decoded_data);
	free(decoded_source);
	free(pchk);
	free(lratio);
 
}