* provided "as is" without express or implied warranty. *
**************************************************************************/
+//-------------------------------------------------------------------------
+// Implementation of the AliSplineFit class
+// The class performs a spline fit on an incoming TGraph. The graph is
+// divided into several parts (identified by knots between each part).
+// Spline fits are performed on each part. According to user parameters,
+// the function, first and second derivative are requested to be continuous
+// at each knot.
+// Origin: Marian Ivanov, CERN, Marian.Ivanov@cern.ch
+// Adjustments by Haavard Helstrup, Haavard.Helstrup@cern.ch
+//-------------------------------------------------------------------------
+
#include "AliSplineFit.h"
-
-ClassImp(AliSplineFit);
-TLinearFitter AliSplineFit::fitterStatic = TLinearFitter(4,"pol3","");
+ClassImp(AliSplineFit)
+
+TLinearFitter* AliSplineFit::fitterStatic()
+{
+ static TLinearFitter* fit = new TLinearFitter(4,"pol3","");
+ return fit;
+}
AliSplineFit::AliSplineFit() :
fBDump(kFALSE),
fGraph (0),
fNmin (0),
+ fMinPoints(0),
fSigma (0),
fMaxDelta (0),
fN0 (0),
fBDump (source.fBDump),
fGraph (source.fGraph),
fNmin (source.fNmin),
+ fMinPoints (source.fMinPoints),
fSigma (source.fSigma),
fMaxDelta (source.fMaxDelta),
fN0 (source.fN0),
+ fParams (0),
+ fCovars (0),
+ fIndex (0),
fN (source.fN),
- fChi2 (source.fChi2)
+ fChi2 (0.0),
+ fX (0),
+ fY0 (0),
+ fY1 (0),
+ fChi2I (source.fChi2I)
{
//
// Copy constructor
fParams = (TClonesArray*)source.fParams->Clone();
fCovars = (TClonesArray*)source.fCovars->Clone();
for (Int_t i=0; i<fN0; i++) fIndex[i] = source.fIndex[i];
-
+
fX = new Double_t[fN];
fY0 = new Double_t[fN];
fY1 = new Double_t[fN];
fX[i] = source.fX[i];
fY0[i] = source.fY0[i];
fY1[i] = source.fY1[i];
+ fChi2I[i] = source.fChi2I[i];
}
}
AliSplineFit& AliSplineFit::operator=(const AliSplineFit& source){
//
// assignment operator
-//
+//
if (&source == this) return *this;
//
// reassign memory as previous fit could have a different size
//
- if ( fN0 != source.fN0) {
+ if ( fN0 != source.fN0) {
delete fParams;
delete fCovars;
fParams = new TClonesArray("TVectorD",fN0);
fCovars = new TClonesArray("TMatrixD",fN0);
}
- if ( fN != source.fN) {
+ if ( fN != source.fN) {
delete []fX;
delete []fY0;
delete []fY1;
delete []fChi2I;
- fN = source.fN;
+ fN = source.fN;
fX = new Double_t[fN];
fY0 = new Double_t[fN];
fY1 = new Double_t[fN];
}
// use copy constructor (without reassigning memory) to copy values
-
+
new (this) AliSplineFit(source);
-
+
return *this;
}
-
+
AliSplineFit::~AliSplineFit(){
//
// destructor. Don't delete fGraph, as this normally comes as input parameter
// = 3: 3rd derivative
//
// a2 = -(3*a0 -3*b0 + 2*a1*dx +b1*dx)/(dx*dx)
- // a3 = -(-2*a0+2*b0 - a1*dx - b1*dx)/(dx*dx*dx)
+ // a3 = -(-2*a0+2*b0 - a1*dx - b1*dx)/(dx*dx*dx)
Int_t index = TMath::BinarySearch(fN,fX,x);
if (index<0) index =0;
}
TGraph * graph = new TGraph(npoints,time,value);
-
+
delete [] value;
delete [] time;
return graph;
delete [] time;
return graph;
}
-
+
TGraph * AliSplineFit::MakeGraph(Double_t xmin, Double_t xmax, Int_t npoints, Int_t deriv) const {
//
TGraph * AliSplineFit::MakeDiff(TGraph * graph0) const {
//
- // Make graph of difference to reference graph
+ // Make graph of difference to reference graph
//
Int_t npoints=graph0->GetN();
TH1F * AliSplineFit::MakeDiffHisto(TGraph * graph0) const {
//
- // Make histogram of difference to reference graph
+ // Make histogram of difference to reference graph
//
Int_t npoints=graph0->GetN();
}else{
if (y<min) min=y;
if (y>max) max=y;
- }
+ }
}
TH1F *his = new TH1F("hdiff","hdiff", 100, min, max);
void AliSplineFit::InitKnots(TGraph * graph, Int_t min, Int_t iter, Double_t maxDelta){
//
- // initialize knots + estimate sigma of noise + make initial parameters
+ // initialize knots + estimate sigma of noise + make initial parameters
//
//
fNmin = min;
fMaxDelta = maxDelta;
Int_t npoints = fGraph->GetN();
+
+ // Make simple spline if too few points in graph
+
+ if (npoints < fMinPoints ) {
+ CopyGraph();
+ return;
+ }
+
fN0 = (npoints/fNmin)+1;
Float_t delta = Double_t(npoints)/Double_t(fN0-1);
for (Int_t iKnot=0; iKnot<fN0; iKnot++){
TMatrixD & cov = *((TMatrixD*)fCovars->At(iKnot));
cov*=fSigma*fSigma;
- }
+ }
OptimizeKnots(iter);
fN = 0;
fChi2I = new Double_t[fN];
Int_t iKnot=0;
for (Int_t i=0; i<fN0; i++){
- if (fIndex[i]<0) continue;
+ if (fIndex[i]<0) continue;
if (iKnot>=fN) {
printf("AliSplineFit::InitKnots: Knot number > Max knot number\n");
break;
TVectorD * param = (TVectorD*) fParams->At(i);
fX[iKnot] = fGraph->GetX()[fIndex[i]];
fY0[iKnot] = (*param)(0);
- fY1[iKnot] = (*param)(1);
+ fY1[iKnot] = (*param)(1);
fChi2I[iKnot] = 0;
iKnot++;
}
//
//
//
- const Double_t kEpsilon = 1.e-7;
Int_t iPrevious=(iKnot>0) ?iKnot-1: 0;
Int_t iNext =(iKnot<fN0)?iKnot+1: fN0-1;
if (iNext>=fN0) iNext=fN0-1;
Double_t startX = fGraph->GetX()[fIndex[iKnot]];
- AliSplineFit::fitterStatic.ClearPoints();
+ AliSplineFit::fitterStatic()->ClearPoints();
Int_t indPrev = fIndex[iPrevious];
Int_t indNext = fIndex[iNext];
Double_t *graphX = fGraph->GetX();
// ePoint[indVec] = fSigma+TMath::Abs(y)*kEpsilon;
// AliSplineFit::fitterStatic.AddPoint(&dxl,y,fSigma+TMath::Abs(y)*kEpsilon);
}
- AliSplineFit::fitterStatic.AssignData(nPoints,1,xPoint,yPoint,ePoint);
- AliSplineFit::fitterStatic.Eval();
+ AliSplineFit::fitterStatic()->AssignData(nPoints,1,xPoint,yPoint,ePoint);
+ AliSplineFit::fitterStatic()->Eval();
// delete temporary arrays
TMatrixD * covar = (TMatrixD*)fCovars->At(iKnot);
TVectorD * param = (TVectorD*)fParams->At(iKnot);
- AliSplineFit::fitterStatic.GetParameters(*param);
- AliSplineFit::fitterStatic.GetCovarianceMatrix(*covar);
+ AliSplineFit::fitterStatic()->GetParameters(*param);
+ AliSplineFit::fitterStatic()->GetCovarianceMatrix(*covar);
return 0;
}
TVectorD sPrevious = tPrevious*pPrevious+tNext*pNext;
TVectorD sKnot = tpKnot*sPrevious;
- TVectorD sNext = tpNext*sPrevious;
+ TVectorD sNext = tpNext*sPrevious;
TMatrixD csPrevious00(tPrevious, TMatrixD::kMult,cPrevious);
csPrevious00 *= tPrevious.T();
}
-
+
void AliSplineFit::SplineFit(Int_t nder){
//
// Cubic spline fit of graph
- //
+ //
// nder
// nder<0 - no continuity requirement
// =0 - continous 0 derivative
// =1 - continous 1 derivative
- // >1 - continous 2 derivative
+ // >1 - continous 2 derivative
//
if (!fGraph) return;
TGraph * graph = fGraph;
if (nder>1) nder=2;
Int_t nknots = fN;
+ if (nknots < 2 ) return;
Int_t npoints = graph->GetN();
//
//
// spline fit
// each knot 4 parameters
- //
+ //
TMatrixD *pmatrix = 0;
TVectorD *pvalues = 0;
if (nder>1){
delete h2;
delete hS;
}
- delete pcstream;
+ delete pcstream;
+}
+
+void AliSplineFit::Cleanup(){
+ //
+ // deletes extra information to reduce amount of information stored on the data
+ // base
+
+ delete fGraph; fGraph=0;
+ delete fParams; fParams=0;
+ delete fCovars; fCovars=0;
+ delete [] fIndex; fIndex=0;
+ delete [] fChi2I; fChi2I=0;
+}
+
+
+void AliSplineFit::CopyGraph() {
+ //
+ // enter graph points directly to fit parameters
+ // (to bee used when too few points are available)
+ //
+ fN = fGraph->GetN();
+ fX = new Double_t[fN];
+ fY0 = new Double_t[fN];
+ fY1 = new Double_t[fN];
+ for (Int_t i=0; i<fN; i++ ) {
+ fX[i] = fGraph->GetX()[i];
+ fY0[i] = fGraph->GetY()[i];
+ fY1[i] = 0;
+ }
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff