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

#define MAX_RANDOM        2147483647  // Maximum value of random()
#define NODES             18432   // Maximum of number of code/check nodes
#define J                    11   // Maximum degree of bit nodes
#define K                    8   // Maximum degree of check nodes

int max_size_M;
int max_size_N;

int n;                            // length
int k;                            // dimension
int nk;                           // redundancy
float  rate;                      // code rate

int M,N;                          // Size of parity-check matrix

// ---------------
// NODE STRUCTURES
// ---------------

struct parent_node {
    int size;
    int index[J];                     // indexes of children
    };

struct child_node {
    int size;
    int index[K];                     // indexes of parents
    };

struct parent_node code_node[NODES];
struct child_node check_node[NODES];

double init_snr, final_snr, snr_increment;
double sim, num_sim, ber, fer, ber_rand, amp, tmp;



int error;
int data_int;
int i,j, i2a, j2a, i3a, j3a, i4a, j4a, i5a, j5a, i6a, j6a, i7a, j7a, i8a, j8a, i2b, j2b, i3b, j3b, i4b, j4b, i5b, j5b, i6b, j6b, i7b, j7b, i8b, j8b, l, pa, qa, ta, ua, va, 
wa, pb, qb, tb, ub, vb, wb, i2, i3, j2, j3, pause;
int max_cycle, nb_cycle[100];
char filename[40];
FILE *fp;

// ----------- Functions declaration -------------------

double max(double x,double y);
double min(double x,double y);

// ----------- End of functions declaration -------------------------


main(int argc, char *argv[])
{
  // Command line processing
  if (argc != 5)
    {
      printf("Usage: %s length(n) dimension(k) file_parity-check cycle_max_length\n",
                      argv[0]);
      exit(0);
    }

  sscanf(argv[1],"%d", &n);
  sscanf(argv[2],"%d", &k);
  sscanf(argv[3],"%s", filename);
  sscanf(argv[4],"%d", &max_cycle);

  nk = n-k;
  rate = (float) k / (float) n;
  max_cycle = min(10,2*(int)(max_cycle/2));

  printf("CYCLE ENCOUNTER FOR BINARY SPARSE MATRIX\n");
  printf("Cycles length counter up %d\n", max_cycle);
  printf("\nn=%d, k=%d, n-k=%d, and rate = %lf\n\n",n,k,nk,rate);

  if ((fp = fopen(filename,"r")) != NULL)
    {
      fscanf(fp, "%d %d", &N, &M);
      fscanf(fp, "%d %d", &max_size_N, &max_size_M);
      for (i=0; i<M; i++) fscanf(fp, "%d", &check_node[i].size);
      for (i=0; i<N; i++) {
          fscanf(fp, "%d", &code_node[i].size);
      }

      // Read index sets for check nodes
      for (i=0; i<M; i++)
      {
        for (j=0; j<check_node[i].size; j++)
        {
          fscanf(fp, "%d", &check_node[i].index[j]);
        }
      }

      // Read index sets for code nodes
      for (i=0; i<N; i++)
      {
        for (j=0; j<code_node[i].size; j++)
        {
          fscanf(fp, "%d", &code_node[i].index[j]);
        }
      }
    }
  else
    {
      printf("incorrect input file name ...\n");
      exit(0);
    }

  fclose(fp);


#ifdef PRINT_MATRIX
  printf("%d %d\n", N, M);
  printf("%d %d\n", max_size_N, max_size_M);
  for (i=0; i<M; i++)
    printf("%8d", check_node[i].size);
  printf("\n");
  for (i=0; i<N; i++)
    printf("%8d", code_node[i].size);
  printf("\n");
  for (i=0; i<M; i++)
  {
    for (j=0; j<check_node[i].size; j++)
      printf("%8d", check_node[i].index[j]);
    printf("\n");
  }
  for (i=0; i<N; i++)
  {
    for (j=0; j<code_node[i].size; j++)
      printf("%8d", code_node[i].index[j]);
    printf("\n");
  }
  printf("\n");
#endif

  for (n=0;n<100;n++) nb_cycle[n]=0;

  for (i=1;i<nk;i++){ // Bypass first row
	for (l=0;l<check_node[i].size;l++){
		j = check_node[i].index[l]-1;
		#ifdef DEBUG
		printf("H(%d,%d)\n",i,j);
		#endif
		for (pa=0;pa<check_node[i].size;pa++){
			j2a = check_node[i].index[pa]-1;
			if (j2a<j){
				for (qa=0;qa<code_node[j2a].size;qa++){
					i2a = code_node[j2a].index[qa]-1;
					if (i2a<i){
						if (max_cycle>=4){
						for (pb=0;pb<code_node[j].size;pb++){
							i2b = code_node[j].index[pb]-1;
							if (i2b<i){
								for (qb=0;qb<check_node[i2b].size;qb++){
									j2b = check_node[i2b].index[qb]-1;
									if (j2b<j){
										if (i2a == i2b){ // Test if cycle length 4 occurs
											if (j2a == j2b){
												nb_cycle[4]++;
												#ifdef DEBUG
		printf("cycle of length 4 with H(%d,%d) - H(%d,%d) - H(%d,%d) - H(%d,%d) - H(%d,%d)\n",i2a+1,j2a+1,i+1,j2a+1,i+1,j+1,i2b+1,j+1,i2b+1,j2b+1);
												#endif
											}
										}
										if (max_cycle>=6){
										for (ta=0;ta<check_node[i2a].size;ta++){
											j3a = check_node[i2a].index[ta]-1;
											if (j3a<j && j3a!=j2a){
												for (tb=0;tb<code_node[j2b].size;tb++){
													i3b = code_node[j2b].index[tb]-1;
													if (i3b<i && i3b!=i2b){
														if (i2a == i3b){ // Test if cycle of length 6 occurs
															if (j3a == j2b){
																nb_cycle[6]++;
																#ifdef DEBUG
																printf("cycle of length 6 with H(%d,%d) - H(%d,%d)\n",i2a,j3a,i3b,j2b);
																#endif
															}
														}
														if (max_cycle>=8){
														for (ua=0;ua<code_node[j3a].size;ua++){
															i3a = code_node[j3a].index[ua]-1;
															if (i3a<i && i3a!=i2a){
																for (ub=0;ub<check_node[i3b].size;ub++){
																	j3b = check_node[i3b].index[ub]-1;
																	if (j3b<j && j3b!=j2b){
																		if (i3a == i3b){ // Test if cycle of length 8 occurs
																			if (j3a == j3b){
																				nb_cycle[8]++;
																				#ifdef DEBUG
																				printf("cycle of length 8 with H(%d,%d) - H(%d,%d)\n",i3a,j3a,i3b,j3b);
																				#endif
																			}
																		}
																		if (max_cycle>=10){
																		for (va=0;va<check_node[i3a].size;va++){
																			j4a = check_node[i3a].index[va]-1;
																			if (j4a<j && j4a!=j3a){
																				for (vb=0;vb<code_node[j3b].size;vb++){
																					i4b = code_node[j3b].index[vb]-1;
																					if (i4b<i && i4b!=i3b){
																						if (i3a == i4b){ // Test if cycle of length 10 occurs
																							if (j4a == j3b){
																								nb_cycle[10]++;
																								#ifdef DEBUG
																								printf("cycle of length 10 with H(%d,%d) - H(%d,%d)\n",i3a,j4a,i4b,j3b);
																								#endif
																							}
																						}
																					}
																				}
																			}
																		}
																		}
																	}
																}
															}
														}
														}
													}
												}
											}
										}
										}
									}
								}
							}
						}
						}
					}
				}
			}
		}
	}
  }
  for (n=4;n<=max_cycle;n+=2) printf("%d cycles of length %d\n", nb_cycle[n], n);
}


// -------------- Functions --------------------



double max(double x,double y)
{
	double maximum;
	maximum = y;
	if (x>y) maximum = x;
	return(maximum);
}

double min(double x,double y)
{
        double minimum;    
        minimum = y;
        if (x<y) minimum = x;
        return(minimum);
}