00001
00002 #include <math.h>
00003
00004 #include "lregress.h"
00005 #include "autocorr.h"
00006
00007 using namespace std;
00008
00009
00027 void autocorr::acorrSlope( acorrInfo &info )
00028 {
00029
00030 const size_t len = info.points().size();
00031 if (len > 0) {
00032 lregress::lineInfo regressInfo;
00033 lregress lr;
00034 vector<lregress::point> regressPoints;
00035 for (size_t i = 0; i < len; i++) {
00036 double x = log(i+1);
00037 double y = log((info.points())[i]);
00038 regressPoints.push_back( lregress::point(x,y) );
00039 }
00040 lr.lr( regressPoints, regressInfo );
00041 info.slope( regressInfo.slope() );
00042 info.slopeErr( regressInfo.slopeErr() );
00043 }
00044 }
00045
00046
00050 void autocorr::ACF( const double *v,
00051 const size_t N,
00052 acorrInfo &info )
00053 {
00054 if (! info.points().empty() ) {
00055 info.points().clear();
00056 }
00057
00058
00059
00060 double *devMean = new double[ N ];
00061
00062 double sum;
00063 size_t i;
00064
00065
00066 sum = 0;
00067 for (i = 0; i < N; i++) {
00068 sum = sum + v[i];
00069 devMean[i] = v[i];
00070 }
00071 double mean = sum / static_cast<double>(N);
00072
00073
00074
00075 sum = 0;
00076 for (i = 0; i < N; i++) {
00077 devMean[i] = devMean[i] - mean;
00078 sum = sum + (devMean[i]*devMean[i]);
00079 }
00080 double denom = sum / static_cast<double>(N);
00081
00082
00083
00084 double cor = 1.0;
00085 for (size_t shift = 1; shift <= numPoints_ && cor > limit_; shift++) {
00086 info.points().push_back( cor );
00087 size_t n = N - shift;
00088 sum = 0.0;
00089 for (i = 0; i < n; i++) {
00090 sum = sum + (devMean[i] * devMean[i+shift]);
00091 }
00092 double numerator = sum / static_cast<double>(n)
00093 cor = numerator / denom;
00094 }
00095
00096
00097 acorrSlope( info );
00098
00099 delete [] m;
00100 }