//===================================================== file = correl.c ====
//=  Program to compute the correlation coefficient for <x, y> data        =
//==========================================================================
//=  Notes:                                                                =
//=    1) Computes the correlation between a paired x and y series         =
//=    2) Uses malloc() to get memory off of the heap                      =
//=    3) The input file contains <x, y> pairs separated by commas         =
//        and/or whitespace.  No comments are allowed in the input file.   =
//=------------------------------------------------------------------------=
//= Example "in.dat" file:                                                 =
//=                                                                        =
//=    14, 2                                                               =
//=    16, 5                                                               =
//=    27, 7                                                               =
//=    42, 9                                                               =
//=    39, 10                                                              =
//=    50, 13                                                              =
//=    83, 20                                                              =
//=------------------------------------------------------------------------=
//= Example output (for above "in.dat"):                                   =
//=                                                                        =
//=   ------------------------------------------------ correl.c -----      =
//=     Number of <x, y> pairs  = 7                                        =
//=     Correlation Coefficient = 0.985632                                 =
//=   ---------------------------------------------------------------      =
//=------------------------------------------------------------------------=
//=  Build:  bcc32 correl.c                                                =
//=------------------------------------------------------------------------=
//=  Execute: correl < in.dat                                              =
//=------------------------------------------------------------------------=
//=  Contact: Ken Christensen                                              =
//=           University of South Florida                                  =
//=           WWW: http://www.csee.usf.edu/~christen                       =
//=           Email: christen@csee.usf.edu                                 =
//=------------------------------------------------------------------------=
//=  History:  KJC (03/08/07) - Genesis (from lr.c)                        =
//==========================================================================

//----- Include files ------------------------------------------------------
#include <stdio.h>   // Needed for printf() and feof()
#include <stdlib.h>  // Needed for atof() and malloc()
#include <math.h>    // Needed for sqrt() and pow()

//----- Defines ------------------------------------------------------------
#define MAX_SIZE  2000000     // Maximum size of x and y arrays

//==========================================================================
//=  Main program                                                          =
//==========================================================================
void main(void)
{
  double   *x;                // Values for x series
  double   *y;                // Values for y series
  char     instring[80];      // Temporary input string
  int      count;             // Counter for number of <x, y> pairs
  double   mom1_x;            // First moment for x series
  double   mom2_x;            // Second moment for x series
  double   mom1_y;            // First moment for x series
  double   mom2_y;            // Second moment for x series
  double   mu_x;              // Mean of x series
  double   sigma_x;           // Standard deviation of x series
  double   mu_y;              // Mean of x series
  double   sigma_y;           // Standard deviation of x series
  double   cov_xy;            // Covariance for x and y series
  double   correl_xy;         // Correlation coefficient for x and y series
  int      i;                 // Loop counter

  // Output a banner
  printf("------------------------------------------------ correl.c ----- \n");

  // Malloc the memory needed for the x and y arrays
  x = (double *) malloc(sizeof(double)*MAX_SIZE);
  y = (double *) malloc(sizeof(double)*MAX_SIZE);
  if ((x == NULL) || (y == NULL))
  {
    printf("*** ERROR - could not malloc sufficent memory \n");
    return;
  }

  // Main loop to read <x, y> values
  count = 0;
  while (!feof(stdin))
  {
    scanf("%s", instring);
    if (feof(stdin)) break;
    x[count] = atof(instring);
    scanf("%s", instring);
    y[count] = atof(instring);
    count++;
    if (count >= MAX_SIZE)
    {
      printf("*** ERROR - more paired values than MAX_SIZE \n");
      return;
    }
  }

  // Compute first and second moments of x and y series
  mom1_x = mom2_x = mom1_y = mom2_y = 0.0;
  for (i=0; i<count; i++)
  {
    mom1_x = mom1_x + (x[i] / count);
    mom2_x = mom2_x + (pow(x[i], 2.0) / count);

    mom1_y = mom1_y + (y[i] / count);
    mom2_y = mom2_y + (pow(y[i], 2.0) / count);
  }

  // Compute the mean and standard deviation of the x and y series
  mu_x = mom1_x;
  mu_y = mom1_y;
  sigma_x = sqrt(mom2_x - pow(mom1_x, 2.0));
  sigma_y = sqrt(mom2_y - pow(mom1_y, 2.0));

  // Compute the covariance
  cov_xy = 0.0;
  for (i=0; i<count; i++)
    cov_xy = cov_xy + (1.0 / count)*((x[i] - mu_x)*(y[i] - mu_y));

  // Compute the correlation coefficient
  correl_xy = cov_xy / (sigma_x * sigma_y);

  // Output results
  printf("  Number of <x, y> pairs  = %d \n", count);
  printf("  Correlation coefficient = %f \n", correl_xy);
  printf("--------------------------------------------------------------- \n");
}