without malloc
[u/mrichter/AliRoot.git] / STAT / TKDInterpolator.h
CommitLineData
f2040a8f 1#ifndef ROOT_TKDInterpolator
2#define ROOT_TKDInterpolator
3
4#ifndef ROOT_TKDTree
5#include "TKDTree.h"
6#endif
316a7f5a 7
1adc956a 8#ifndef ROOT_TVectorD
9#include "TVectorD.h"
10#endif
11#ifndef ROOT_TMatrixD
12#include "TMatrixD.h"
13#endif
14
117c62ab 15// Non parametric interpolation class based on local polinomial regression.
16// The class can be used to approximate PDF together with TKDTree or for
17// general regression when the data points are given.
f2040a8f 18
1adc956a 19// template <typename Value> class TVectorT;
20// typedef struct TVectorT<Double_t> TVectorD;
21// template <typename Value> class TMatrixT;
22// typedef class TMatrixT<Double_t> TMatrixD;
23// template <typename Index, typename Value> class TKDTree;
24// typedef class TKDTree<Int_t, Float_t> TKDTreeIF;
f2040a8f 25class TTree;
26class TLinearFitter;
27class TKDInterpolator : public TKDTreeIF
28{
29public:
30 TKDInterpolator();
316a7f5a 31 TKDInterpolator(TTree *t, const Char_t *var, const Char_t *cut = 0, UInt_t bsize = 100, Long64_t nentries = 1000000000, Long64_t firstentry = 0);
f2040a8f 32 TKDInterpolator(Int_t npoints, Int_t ndim, UInt_t bsize, Float_t **data);
33 ~TKDInterpolator();
e0b38166 34
117c62ab 35 Double_t Eval(const Double_t *point, Double_t &result, Double_t &error, Bool_t force = kFALSE);
36 Float_t GetAlpha() const {return fAlpha;}
37 inline Bool_t GetCOGPoint(Int_t node, Float_t *coord, Float_t &val, Float_t &error) const;
e0b38166 38 inline Bool_t GetDataPoint(Int_t n, Float_t *p) const;
316a7f5a 39 TKDTreeIF* GetHelper() {return fKDhelper;}
117c62ab 40 Bool_t GetInterpolationMethod() const {return fStatus&1;}
41 Int_t GetNTNodes() const {return fNTNodes;}
42 Bool_t GetStore() const {return fStatus&2;}
43 Bool_t GetWeights() const {return fStatus&4;}
316a7f5a 44
45 void DrawNodes(UInt_t ax1 = 0, UInt_t ax2 = 1, Int_t depth = -1);
e0b38166 46 void DrawNode(Int_t tnode, UInt_t ax1 = 0, UInt_t ax2 = 1);
316a7f5a 47 void GetStatus();
117c62ab 48 void SetAlpha(const Float_t a){if(a>0.) fAlpha = a;}
49 void SetInterpolationMethod(const Bool_t on = kTRUE);
50 void SetStore(const Bool_t on = kTRUE);
51 void SetWeights(const Bool_t on = kTRUE);
316a7f5a 52
f2040a8f 53private:
a9c20b1f 54 TKDInterpolator(const TKDInterpolator &);
55 TKDInterpolator& operator=(const TKDInterpolator &);
e0b38166 56 void Build();
316a7f5a 57
117c62ab 58public:
59 class TKDNodeInfo
60 {
61 public:
62 TKDNodeInfo(const Int_t ndim = 0);
63 virtual ~TKDNodeInfo();
64 void Build(const Int_t ndim);
65 Double_t CookPDF(const Double_t *point, Double_t &result, Double_t &error);
66 void Store(const TVectorD &par, const TMatrixD &cov);
67
68 Int_t fNDim; // data dimension
69 Float_t *fRefPoint; //[fNDim] node's COG
70 Float_t fRefValue; // measured value for node
71 TMatrixD *fCov; // interpolator covariance matrix
72 TVectorD *fPar; // interpolator parameters
73
74 private:
75 TKDNodeInfo(const TKDNodeInfo &);
76 TKDNodeInfo& operator=(const TKDNodeInfo &);
77
78 ClassDef(TKDNodeInfo, 1) // node info for interpolator
79 };
80
f2040a8f 81protected:
82 Int_t fNTNodes; //Number of evaluation data points
5f38a39d 83 TKDNodeInfo *fTNodes; //[fNTNodes] interpolation node
f2040a8f 84
85private:
316a7f5a 86 UChar_t fStatus; // status of the interpolator
117c62ab 87 UChar_t fLambda; // number of parameters in polynom
88 Short_t fDepth; //! depth of the KD Tree structure used
89 Float_t fAlpha; // alpha parameter
5f38a39d 90 Float_t **fRefPoints; //! temporary storage of COG data
316a7f5a 91 Double_t *fBuffer; //! working space [2*fLambda]
92 TKDTreeIF *fKDhelper; //! kNN finder
93 TLinearFitter *fFitter; //! linear fitter
f2040a8f 94
95 ClassDef(TKDInterpolator, 1) // data interpolator based on KD tree
96};
97
98//__________________________________________________________________
5f38a39d 99Bool_t TKDInterpolator::GetCOGPoint(Int_t node, Float_t *coord, Float_t &val, Float_t &err) const
f2040a8f 100{
e0b38166 101 if(node < 0 || node > fNTNodes) return kFALSE;
f2040a8f 102
5f38a39d 103 for(int idim=0; idim<fNDim; idim++) coord[idim] = fTNodes[node].fRefPoint[idim];
104 val = fTNodes[node].fRefValue;
105 err = fTNodes[node].fRefValue/TMath::Sqrt(fBucketSize);
f2040a8f 106 return kTRUE;
107}
108
e0b38166 109//__________________________________________________________________
110Bool_t TKDInterpolator::GetDataPoint(Int_t n, Float_t *p) const
111{
112 if(n < 0 || n >= fNpoints) return kFALSE;
117c62ab 113 if(!fData) return kFALSE;
114
e0b38166 115 for(int i=0; i<fNDim; i++) p[i] = fData[i][n];
116 return kTRUE;
117}
118
119
f2040a8f 120#endif
121