#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>

static void calcSNR(unsigned char *org, unsigned char *new, 
		    int w, int h, double *pv, double *pe);


main(){

  int Max_Frame, i,j,k, p = 0;
  FILE *fd;
  FILE *par;
/******************************************************************************************/
  //darshan
  FILE *YUV;	
  double snrY,snrU,snrV;		
  double sum_snrY,sum_snrU,sum_snrV;		
/******************************************************************************************/
  char org_Y[128], org_U[128], org_V[128];
  char new_Y[128], new_U[128], new_V[128];
  int chrom_hize,chrom_vsize;
  int width, height, chrom_width, chrom_height, size;
  int horizontal, vertical, horizontal1, vertical1;
  unsigned char *org[3];
  unsigned char *new[3];
  int counter =0 ;
  double e, v; 

  if((YUV = fopen("SNR.YUV","w")) == NULL)
    printf("Warnning !!! FILE SNR.YUV could not open\n");
  else{
  printf("FILE SNR.YUV CREATED\n");

   }

  if((par = fopen("SNR.par", "r")) == NULL)
    printf("Warnning !!! FILE SNR.par could not open\n");
  else{
    fscanf(par, "%d %d %d", &width, &height, &Max_Frame);
  }

  fclose(par);

  //  printf(" W: %d \tH: %d\t MAX:%d\n", width, height, Max_Frame);
  Max_Frame = Max_Frame + 1;

  horizontal = width;
  vertical = height;

  chrom_width = width/2;
  chrom_height = height/2;

  horizontal1 = chrom_width;
  vertical1 = chrom_height;
  
  

  //printf(" width = %d height = %d \n", width, height);
  //Dynamic memory allocation
  for(j=0; j< 3; j++){
    size = (j==0) ? width*height : chrom_width*chrom_height; 
    //printf("size = %d\n", size);
    org[j] = (unsigned char *)malloc(size);
    new[j] = (unsigned char *)malloc(size);
  }


  for(k = 0; k< Max_Frame; k++){

    printf("SNR:Frame number %d \n", k);

    /************* Original Y **************/
    sprintf(org_Y, "org%d.Y", k);
    //printf("org:Y = %s\t", org_Y);
    if((fd = fopen(org_Y,"rb")) == NULL){
         printf(" Couldn't open %s\n", org_Y);
    }
    else{
      for(i = 0; i < vertical; i++)
	fread(org[0]+i*width, 1, horizontal, fd);
    }

    fclose(fd);

    /************* Original U **************/
    sprintf(org_U, "org%d.U", k);
    //printf("org:U = %s\t", org_U);
    if((fd = fopen(org_U,"rb")) == NULL){
         printf(" Couldn't open %s\n", org_U);
    }
    else{
      for(i = 0; i < vertical1; i++)
	fread(org[1]+i*chrom_width, 1, horizontal1, fd);
    }

    fclose(fd);

    /************* Original V **************/
    sprintf(org_V, "org%d.V", k);
    //printf("org:V = %s\n", org_V);
    if((fd = fopen(org_V,"rb")) == NULL){
         printf(" Couldn't open %s\n", org_V);
    }
    else{
      for(i = 0; i < vertical1; i++)
	fread(org[2]+i*chrom_width, 1, horizontal1, fd);
    }

    fclose(fd);


    /************* NEW Y **************/
    sprintf(new_Y, "new%d.Y", k);
    //printf("new:Y = %s\t", new_Y);
    if((fd = fopen(new_Y,"rb")) == NULL){
         printf(" Couldn't open %s\n", new_Y);
    }
    else{
      for(i = 0; i < vertical; i++)
	fread(new[0]+i*width, 1, horizontal, fd);
    }

    fclose(fd);
    
    /************* NEW U **************/
    sprintf(new_U, "new%d.U", k);
    //printf("new:U = %s\t", new_U);
    if((fd = fopen(new_U,"rb")) == NULL){
         printf(" Couldn't open %s\n", new_U);
    }
    else{
      for(i = 0; i < vertical1; i++)
	fread(new[1]+i*chrom_width, 1, horizontal1, fd);
    }

    fclose(fd);
    
    /************* NEW V **************/
    sprintf(new_V, "new%d.V", k);
    //printf("new:V = %s\n", new_V);
    if((fd = fopen(new_V,"rb")) == NULL){
         printf(" Couldn't open %s\n", new_V);
    }
    else{
      for(i = 0; i < vertical1; i++)
	fread(new[2]+i*chrom_width, 1, horizontal1, fd);
    }

    fclose(fd);
    
    for(i=0; i<3; i++){
      width = (i==0) ? horizontal : chrom_width; 
      height = (i==0) ? vertical : chrom_height; 
      calcSNR(org[i], new[i], width, height, &v, &e);
      if(i == 0){
//      printf("SNR_Y = %3.3g dB\n", 10.0*log10(255.0*255.0/e));
      snrY = 10.0*log10(255.0*255.0/e);
     // snrY = 10.0*log10(v/e);
//      printf("SNR_Y = %3.3g dB\n", snrY);
      }
      if(i == 1){
//     printf("SNR_U = %3.3g dB\n", 10.0*log10(255.0*255.0/e));
      snrU = 10.0*log10(255.0*255.0/e);
      //snrU = 10.0*log10(v/e);
//    printf("SNR_U = %3.3g dB\n", snrU);
      }
      if(i == 2)
      {
//      printf("SNR_V = %3.3g dB\n\n", 10.0*log10(255.0*255.0/e));
      snrV = 10.0*log10(255.0*255.0/e);
     // snrV = 10.0*log10(v/e);
//      printf("SNR_V = %3.3g dB\n", snrV);
       }



     }

     sum_snrY = sum_snrY + snrY;
     sum_snrU = sum_snrU + snrU;
     sum_snrV = sum_snrV + snrV;
     	
    // printf("%3.3g\t%3.3g\t%3.3g\n",snrY,snrU,snrV);	//darshan
   //  fprintf(YUV,"%3.3g\t%3.3g\t%3.3g\n",snrY,snrU,snrV);	//darshan
    }	//end of for loop for NO.OF FRAMES
//    printf("Max_Frame=%d\n",Max_Frame);		//darshan
    printf("%3.3g\t%3.3g\t%3.3g\n",(sum_snrY/Max_Frame),(sum_snrU/Max_Frame),(sum_snrV/Max_Frame));	//darshan
    fprintf(YUV,"%3.3g\t%3.3g\t%3.3g\n",(sum_snrY/Max_Frame),(sum_snrU/Max_Frame),(sum_snrV/Max_Frame));	//darshan

    for(j=0; j< 3; j++){
      free(org[j]);
      free(new[j]);
    }
 //   printf("\nWidth = %d\tHeight = %d\tTotal Frame=%d\n", horizontal, vertical, Max_Frame);

    fclose(YUV);	//darshan
}


static void calcSNR(org, new, w, h, pv, pe)
unsigned char *org;
unsigned char *new;
int w,h;
double *pv, *pe;
{

  int i , j;
  double v1,s1,s2,e2;

  s1 = s2 = e2 = 0.0;

  for(j=0; j<h; j++){
    for(i=0; i<w;i++){
      v1  = org[i];
      s1 += v1;
      s2 += v1 *v1;
      v1 -= new[i];
      //printf("org[i]=%d\tnew[i]=%d\tv1=%f\n",org[i],new[i],v1);
      e2 += v1 * v1;
    }
  }


  //  printf("width = %d height %d\n", w, h);
  //printf("v1  = %f\ts1 = %f\ts2 = %f\te2 = %f\n",v1,s1,s2,e2 );

  s1 /= w * h;
  s2 /= w * h;
  e2 /= w * h;
  //printf("s1 = %f\ts2 = %f\te2 = %f\n",s1,s2,e2 );
  
  if(e2 == 0.0) e2 = 0.00001;

  *pv = s2 - s1 * s1;
  *pe = e2;

  printf("pv = %f\t pe = %f\n",*pv,*pe );
}