/*	correlation.c
**	this program generate correlation
**      coefficient for two input images!
**  COPY RIGHT: BRANDT TSO
*/

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


void main(void)
{
  FILE 		*fin1, *fin2;
  unsigned char store1[2],store2[2],in_file[50];
  unsigned int pat_length;
  unsigned int count1, img_row,img_col;
  double total1,total2,ex1 = 2,ex2 = 0.5;
  double mean1,mean2, var1,var2, diff1,diff2,div1,div2,r;




      printf("Please input the number of image rows:");
      scanf("%d",&img_row);
	  printf("Please input the number of image columns:");
      scanf("%d",&img_col);

      pat_length = img_row*img_col;

	  printf("Please input image file 1:");
	  scanf("%s",in_file);
      if((fin1 = fopen(in_file,"rb")) == NULL)
	{
	  fprintf(stderr,"Failed to open %s\n",in_file);
	  exit(EXIT_FAILURE);
	}

	  printf("Please input image file 2:");
	  scanf("%s",in_file);
      if((fin2 = fopen(in_file,"rb")) == NULL)
       {
         fprintf(stderr,"Failed to open %s\n",in_file);
         exit(EXIT_FAILURE);
       }




      mean1=0; mean2=0;
      total1=0; total2=0;


      for(count1 = 0;count1 < pat_length;count1++)
	{
          fread(store1,1,1,fin1);
          fread(store2,1,1,fin2);
          total1 = total1 + store1[0];
          total2 = total2 + store2[0];
	}

      mean1 = (double) total1/pat_length;
      mean2 = (double) total2/pat_length;
      fseek(fin1,0,0);
      fseek(fin2,0,0);
      total1=0; total2=0;

      for(count1 = 0;count1 < pat_length;count1++)
	{
          fread(store1,1,1,fin1);
          fread(store2,1,1,fin2);
          diff1 = store1[0] - mean1;
          diff2 = store2[0] - mean2;
          if(diff1 < 0)
             diff1 = -diff1;
          if(diff2 < 0)
             diff2 = -diff2;
          total1 = total1 + pow(diff1,ex1);
          total2 = total2 + pow(diff2,ex1);
	}

       var1 = (double) total1/pat_length;
       var2 = (double) total2/pat_length;
       div1 = pow(var1,ex2);
       div2 = pow(var2,ex2);

       fseek(fin1,0,0);
       fseek(fin2,0,0);
       total1=0; total2=0;
       r =0;

       for(count1 = 0;count1 < pat_length;count1++)
	{
          fread(store1,1,1,fin1);
          fread(store2,1,1,fin2);
          diff1 = store1[0] - mean1;
          diff2 = store2[0] - mean2;
          r = r + (double) (diff1*diff2)/(div1*div2);
	}

        r = (double) r/pat_length;
        printf("\n**************Done****************\n");
        printf("mean for image 1:%f image 2:%f\n",mean1,mean2);
        printf("deviation1:%f deviation2:%f\n",div1,div2);
        printf("Correlation coefficient:%f\n",r);


        fclose(fin1);
        fclose(fin2);



}