//==============================================================================
//  Function for calculation of coefficient of correlation between two vectors:
//  input:
//  nPoints   - number of points for processing
//  pfVect1   - first array with data
//  pfVect2   - second array with data
//  return:
//  coefficient of correlation
//==============================================================================
#include <math.h>
double CalcLinCorrelationCoef(int nPoints, double* pfVectX, double* pfVectY)
{
  int i;
  double fCorr;
  double fDenom;
  double fPoints = (double) nPoints;
  
  double fXSum  = 0., fYSum  = 0.;
  double fX2Sum = 0., fY2Sum = 0.;
  double fXYSum = 0.;
  
  if(nPoints < 2) return -2.; // not enough data !
  
  for(i=0; i<nPoints; i++)
  {
    fXSum += pfVectX[i];
    fYSum += pfVectY[i];
    
    fX2Sum += pfVectX[i] * pfVectX[i];
    fY2Sum += pfVectY[i] * pfVectY[i];
    
    fXYSum += pfVectX[i] * pfVectY[i];
  }

  fDenom = (fPoints*fX2Sum - fXSum*fXSum)*(fPoints*fY2Sum - fYSum*fYSum);
  if(fDenom < 1.E-100) return -3.;
  
  fCorr  = fPoints*fXYSum - fXSum*fYSum;
  fCorr  = fCorr / sqrt(fDenom);
  
  return fCorr;
}

//==============================================================================

#include <stdio.h>

//==============================================================================
//  TEST for function CalcCorrelationCoef :
//==============================================================================
#define  NPOINTS 10
int main(int nArgs, char** ppArgs)
{
  int n = NPOINTS;
  double dat1[NPOINTS] = {0., 1., 2., 3., 4., 5., 6., 7., 8., 9.};
  double dat2[NPOINTS] = {0.1, 1.1, 1.9, 2.95, 4.04, 5.2, 5.9, 7.02, 7.9, 9.0};
  
  double fCorr;
  
  fCorr = CalcLinCorrelationCoef(NPOINTS, dat1, dat2);
  
  printf("correlation=%f\n", fCorr);

  return 0;
}