TuneA settings removed from kPyJets
[u/mrichter/AliRoot.git] / STAT / TKDPDF.cxx
CommitLineData
a3408ed3 1#include "TKDPDF.h"
2#include "TKDNodeInfo.h"
3
4#include "TClonesArray.h"
5#include "TTree.h"
6#include "TH2.h"
7#include "TObjArray.h"
8#include "TObjString.h"
9#include "TBox.h"
10#include "TGraph.h"
11#include "TMarker.h"
12
13ClassImp(TKDPDF)
14
15
16
17//_________________________________________________________________
18TKDPDF::TKDPDF() :
50c4eb3a 19 TKDTreeIF()
20 ,TKDInterpolatorBase()
a3408ed3 21{
22// Default constructor. To be used with care since in this case building
23// of data structure is completly left to the user responsability.
24}
25
26//_________________________________________________________________
27TKDPDF::TKDPDF(Int_t npoints, Int_t ndim, UInt_t bsize, Float_t **data) :
50c4eb3a 28 TKDTreeIF(npoints, ndim, bsize, data)
29 ,TKDInterpolatorBase(ndim)
a3408ed3 30{
31// Wrapper constructor for the TKDTree.
32
50c4eb3a 33 Build();
a3408ed3 34}
35
36
37//_________________________________________________________________
38TKDPDF::TKDPDF(TTree *t, const Char_t *var, const Char_t *cut, UInt_t bsize, Long64_t nentries, Long64_t firstentry) :
50c4eb3a 39 TKDTreeIF()
40 ,TKDInterpolatorBase()
a3408ed3 41{
42// Alocate data from a tree. The variables which have to be analysed are
43// defined in the "var" parameter as a colon separated list. The format should
44// be identical to that used by TTree::Draw().
45//
46//
47
50c4eb3a 48 TObjArray *vars = TString(var).Tokenize(":");
49 fNDim = vars->GetEntriesFast(); fNDimm = 2*fNDim;
50 fNSize = fNDim;
eb8908e3 51 if(fNDim > 6/*kDimMax*/) Warning("TKDPDF(TTree*, const Char_t, const Char_t, UInt_t)", "Variable number exceed maximum dimension %d. Results are unpredictable.", 6/*kDimMax*/);
50c4eb3a 52 fBucketSize = bsize;
53
54 Int_t np;
55 Double_t *v;
56 for(int idim=0; idim<fNDim; idim++){
57 if(!(np = t->Draw(((TObjString*)(*vars)[idim])->GetName(), cut, "goff", nentries, firstentry))){
eb8908e3 58 Warning("TKDPDF(TTree*, const Char_t, const Char_t, UInt_t)", "Can not access data for keys %s. Key defined on tree :", ((TObjString*)(*vars)[idim])->GetName());
50c4eb3a 59 TIterator *it = (t->GetListOfLeaves())->MakeIterator();
60 TObject *o;
61 while((o = (*it)())) printf("\t%s\n", o->GetName());
62 continue;
63 }
64 if(!fNPoints){
65 fNPoints = np;
66 //Info("TKDPDF(TTree*, const Char_t, const Char_t, UInt_t)", Form("Allocating %d data points in %d dimensions.", fNpoints, fNDim));
67 fData = new Float_t*[fNDim];
34199116 68 for(int jdim=fNDim; jdim--;) fData[jdim] = new Float_t[fNPoints];
50c4eb3a 69 fDataOwner = kTRUE;
70 }
71 v = t->GetV1();
72 for(int ip=0; ip<fNPoints; ip++) fData[idim][ip] = (Float_t)v[ip];
73 }
50c4eb3a 74 Build();
09d5920f 75 delete vars;
a3408ed3 76}
77
78//_________________________________________________________________
79TKDPDF::~TKDPDF()
80{
81}
82
83//_________________________________________________________________
afb669c1 84Bool_t TKDPDF::Build(Int_t)
a3408ed3 85{
86// Fill interpolator's data array i.e.
87// - estimation points
88// - corresponding PDF values
89
53ee3cad 90 TKDTreeIF::Build();
50c4eb3a 91 if(!fBoundaries) MakeBoundaries();
92 fLambda = 1 + fNDim + (fNDim*(fNDim+1)>>1);
93 //printf("after MakeBoundaries() %d\n", memory());
94
95 // allocate interpolation nodes
afb669c1 96 Int_t fNTNodes = fNPoints/fBucketSize + ((fNPoints%fBucketSize)?1:0);/*TKDTreeIF::GetNTNodes();*/
50c4eb3a 97 TKDInterpolatorBase::Build(fNTNodes);
98
afb669c1 99 TKDNodeInfo *node = NULL;
100 Float_t *bounds = NULL;
50c4eb3a 101 Int_t *indexPoints;
102 for(int inode=0, tnode = fNNodes; inode<fNTNodes-1; inode++, tnode++){
afb669c1 103 node = (TKDNodeInfo*)(*fNodes)[inode];
50c4eb3a 104 node->Val()[0] = Float_t(fBucketSize)/fNPoints;
105 bounds = GetBoundary(tnode);
106 for(int idim=0; idim<fNDim; idim++) node->Val()[0] /= (bounds[2*idim+1] - bounds[2*idim]);
107 node->Val()[1] = node->Val()[0]/TMath::Sqrt(float(fBucketSize));
108
109 indexPoints = GetPointsIndexes(tnode);
110 // loop points in this terminal node
111 for(int idim=0; idim<fNDim; idim++){
112 node->Data()[idim] = 0.;
113 for(int ip = 0; ip<fBucketSize; ip++) node->Data()[idim] += fData[idim][indexPoints[ip]];
114 node->Data()[idim] /= fBucketSize;
115 }
116 memcpy(&(node->Data()[fNDim]), bounds, fNDimm*sizeof(Float_t));
117 }
118
119 // analyze last (incomplete) terminal node
120 Int_t counts = fNPoints%fBucketSize;
121 counts = counts ? counts : fBucketSize;
122 Int_t inode = fNTNodes - 1, tnode = inode + fNNodes;
afb669c1 123 node = (TKDNodeInfo*)(*fNodes)[inode];
50c4eb3a 124 node->Val()[0] = Float_t(counts)/fNPoints;
125 bounds = GetBoundary(tnode);
53ee3cad 126 for(int idim=0; idim<fNDim; idim++){
127 Float_t dx = bounds[2*idim+1]-bounds[2*idim];
128 if(dx < 1.e-30){
eb8908e3 129 Warning("TKDPDF::Build()", "Terminal bucket index[%d] too narrow on the %d dimension.", inode, idim);
53ee3cad 130 continue;
131 }
132 node->Val()[0] /= (bounds[2*idim+1] - bounds[2*idim]);
133 }
50c4eb3a 134 node->Val()[1] = node->Val()[0]/TMath::Sqrt(float(counts));
135
136 // loop points in this terminal node
137 indexPoints = GetPointsIndexes(tnode);
138 for(int idim=0; idim<fNDim; idim++){
139 node->Data()[idim] = 0.;
140 for(int ip = 0; ip<counts; ip++) node->Data()[idim] += fData[idim][indexPoints[ip]];
141 node->Data()[idim] /= counts;
142 }
143 memcpy(&(node->Data()[fNDim]), bounds, fNDimm*sizeof(Float_t));
144
145 delete [] fBoundaries;
afb669c1 146 fBoundaries = NULL;
147
148 return kTRUE;
a3408ed3 149}
150
151
152//_________________________________________________________________
153void TKDPDF::DrawNode(Int_t tnode, UInt_t ax1, UInt_t ax2)
154{
155// Draw node "node" and the data points within.
156//
157// Observation:
158// This function creates some graphical objects
159// but don't delete it. Abusing this function may cause memory leaks !
160
afb669c1 161 if(tnode < 0 || tnode >= GetNTNodes()){
eb8908e3 162 Warning("DrawNode()", "Terminal node %d outside defined range.", tnode);
50c4eb3a 163 return;
164 }
165
166 Int_t inode = tnode;
167 tnode += fNNodes;
168 // select zone of interest in the indexes array
169 Int_t *index = GetPointsIndexes(tnode);
170 Int_t nPoints = (tnode == 2*fNNodes) ? fNPoints%fBucketSize : fBucketSize;
171
172 // draw data points
173 TGraph *g = new TGraph(nPoints);
174 g->SetMarkerStyle(7);
175 for(int ip = 0; ip<nPoints; ip++) g->SetPoint(ip, fData[ax1][index[ip]], fData[ax2][index[ip]]);
176
177 // draw estimation point
afb669c1 178 TKDNodeInfo *node = (TKDNodeInfo*)(*fNodes)[inode];
50c4eb3a 179 TMarker *m=new TMarker(node->Data()[ax1], node->Data()[ax2], 20);
180 m->SetMarkerColor(2);
181 m->SetMarkerSize(1.7);
182
183 // draw node contour
184 Float_t *bounds = GetBoundary(tnode);
185 TBox *n = new TBox(bounds[2*ax1], bounds[2*ax2], bounds[2*ax1+1], bounds[2*ax2+1]);
186 n->SetFillStyle(0);
187
188 g->Draw("ap");
189 m->Draw();
190 n->Draw();
191
192 return;
a3408ed3 193}
194