ROOT logo
#include "TKDPDF.h"
#include "TKDNodeInfo.h"

#include "TClonesArray.h"
#include "TTree.h"
#include "TH2.h"
#include "TObjArray.h"
#include "TObjString.h"
#include "TBox.h"
#include "TGraph.h"
#include "TMarker.h"

ClassImp(TKDPDF)



//_________________________________________________________________
TKDPDF::TKDPDF() :
  TKDTreeIF()
  ,TKDInterpolatorBase()
{
// Default constructor. To be used with care since in this case building
// of data structure is completly left to the user responsability.
}

//_________________________________________________________________
TKDPDF::TKDPDF(Int_t npoints, Int_t ndim, UInt_t bsize, Float_t **data) :
  TKDTreeIF(npoints, ndim, bsize, data)
  ,TKDInterpolatorBase(ndim)
{
// Wrapper constructor for the TKDTree.

  Build();
}


//_________________________________________________________________
TKDPDF::TKDPDF(TTree *t, const Char_t *var, const Char_t *cut, UInt_t bsize, Long64_t nentries, Long64_t firstentry) :
  TKDTreeIF()
  ,TKDInterpolatorBase()
{
// Alocate data from a tree. The variables which have to be analysed are
// defined in the "var" parameter as a colon separated list. The format should
// be identical to that used by TTree::Draw().
//
// 

  TObjArray *vars = TString(var).Tokenize(":");
  fNDim = vars->GetEntriesFast(); fNDimm = 2*fNDim;
  fNSize = fNDim;
  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*/);
  fBucketSize = bsize;

  Int_t np;
  Double_t *v;
  for(int idim=0; idim<fNDim; idim++){
    if(!(np = t->Draw(((TObjString*)(*vars)[idim])->GetName(), cut, "goff", nentries, firstentry))){
      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());
      TIterator *it = (t->GetListOfLeaves())->MakeIterator();
      TObject *o;
      while((o = (*it)())) printf("\t%s\n", o->GetName());
      continue;
    }
    if(!fNPoints){
      fNPoints = np;
      //Info("TKDPDF(TTree*, const Char_t, const Char_t, UInt_t)", Form("Allocating %d data points in %d dimensions.", fNpoints, fNDim));
      fData = new Float_t*[fNDim];
      for(int jdim=fNDim; jdim--;) fData[jdim] = new Float_t[fNPoints];
      fDataOwner = kTRUE;
    }
    v = t->GetV1();
    for(int ip=0; ip<fNPoints; ip++) fData[idim][ip] = (Float_t)v[ip];
  }
  Build();
  delete vars;
}

//_________________________________________________________________
TKDPDF::~TKDPDF()
{
}

//_________________________________________________________________
Bool_t TKDPDF::Build(Int_t)
{
// Fill interpolator's data array i.e.
//  - estimation points 
//  - corresponding PDF values

  TKDTreeIF::Build();
  if(!fBoundaries) MakeBoundaries();
  fLambda = 1 + fNDim + (fNDim*(fNDim+1)>>1);
  //printf("after MakeBoundaries() %d\n", memory());
  
  // allocate interpolation nodes
  Int_t fNTNodes = fNPoints/fBucketSize + ((fNPoints%fBucketSize)?1:0);/*TKDTreeIF::GetNTNodes();*/
  TKDInterpolatorBase::Build(fNTNodes);

  TKDNodeInfo *node = NULL;
  Float_t *bounds = NULL;
  Int_t *indexPoints;
  for(int inode=0, tnode = fNNodes; inode<fNTNodes-1; inode++, tnode++){
    node = (TKDNodeInfo*)(*fNodes)[inode];
    node->Val()[0] =  Float_t(fBucketSize)/fNPoints;
    bounds = GetBoundary(tnode);
    for(int idim=0; idim<fNDim; idim++) node->Val()[0] /= (bounds[2*idim+1] - bounds[2*idim]);
    node->Val()[1] =  node->Val()[0]/TMath::Sqrt(float(fBucketSize));
    
    indexPoints = GetPointsIndexes(tnode);
    // loop points in this terminal node
    for(int idim=0; idim<fNDim; idim++){
      node->Data()[idim] = 0.;
      for(int ip = 0; ip<fBucketSize; ip++) node->Data()[idim] += fData[idim][indexPoints[ip]];
      node->Data()[idim] /= fBucketSize;
    }
    memcpy(&(node->Data()[fNDim]), bounds, fNDimm*sizeof(Float_t));
  }

  // analyze last (incomplete) terminal node
  Int_t counts = fNPoints%fBucketSize;
  counts = counts ? counts : fBucketSize;
  Int_t inode = fNTNodes - 1, tnode = inode + fNNodes;
  node = (TKDNodeInfo*)(*fNodes)[inode];
  node->Val()[0] =  Float_t(counts)/fNPoints;
  bounds = GetBoundary(tnode);
  for(int idim=0; idim<fNDim; idim++){ 
    Float_t dx = bounds[2*idim+1]-bounds[2*idim];
    if(dx < 1.e-30){
      Warning("TKDPDF::Build()", "Terminal bucket index[%d] too narrow on the %d dimension.", inode, idim);
      continue;
    }
    node->Val()[0] /= (bounds[2*idim+1] - bounds[2*idim]);
  }
  node->Val()[1] =  node->Val()[0]/TMath::Sqrt(float(counts));

  // loop points in this terminal node
  indexPoints = GetPointsIndexes(tnode);
  for(int idim=0; idim<fNDim; idim++){
    node->Data()[idim] = 0.;
    for(int ip = 0; ip<counts; ip++) node->Data()[idim] += fData[idim][indexPoints[ip]];
    node->Data()[idim] /= counts;
  }
  memcpy(&(node->Data()[fNDim]), bounds, fNDimm*sizeof(Float_t));

  delete [] fBoundaries;
  fBoundaries = NULL;

  return kTRUE;
}


//_________________________________________________________________
void TKDPDF::DrawNode(Int_t tnode, UInt_t ax1, UInt_t ax2)
{
// Draw node "node" and the data points within.
//
// Observation:
// This function creates some graphical objects
// but don't delete it. Abusing this function may cause memory leaks !

  if(tnode < 0 || tnode >= GetNTNodes()){
    Warning("DrawNode()", "Terminal node %d outside defined range.", tnode);
    return;
  }

  Int_t inode = tnode;
  tnode += fNNodes;
  // select zone of interest in the indexes array
  Int_t *index = GetPointsIndexes(tnode);
  Int_t nPoints = (tnode == 2*fNNodes) ? fNPoints%fBucketSize : fBucketSize;

  // draw data points
  TGraph *g = new TGraph(nPoints);
  g->SetMarkerStyle(7);
  for(int ip = 0; ip<nPoints; ip++) g->SetPoint(ip, fData[ax1][index[ip]], fData[ax2][index[ip]]);

  // draw estimation point
  TKDNodeInfo *node = (TKDNodeInfo*)(*fNodes)[inode];
  TMarker *m=new TMarker(node->Data()[ax1], node->Data()[ax2], 20);
  m->SetMarkerColor(2);
  m->SetMarkerSize(1.7);
  
  // draw node contour
  Float_t *bounds = GetBoundary(tnode);
  TBox *n = new TBox(bounds[2*ax1], bounds[2*ax2], bounds[2*ax1+1], bounds[2*ax2+1]);
  n->SetFillStyle(0);

  g->Draw("ap");
  m->Draw();
  n->Draw();
  
  return;
}

 TKDPDF.cxx:1
 TKDPDF.cxx:2
 TKDPDF.cxx:3
 TKDPDF.cxx:4
 TKDPDF.cxx:5
 TKDPDF.cxx:6
 TKDPDF.cxx:7
 TKDPDF.cxx:8
 TKDPDF.cxx:9
 TKDPDF.cxx:10
 TKDPDF.cxx:11
 TKDPDF.cxx:12
 TKDPDF.cxx:13
 TKDPDF.cxx:14
 TKDPDF.cxx:15
 TKDPDF.cxx:16
 TKDPDF.cxx:17
 TKDPDF.cxx:18
 TKDPDF.cxx:19
 TKDPDF.cxx:20
 TKDPDF.cxx:21
 TKDPDF.cxx:22
 TKDPDF.cxx:23
 TKDPDF.cxx:24
 TKDPDF.cxx:25
 TKDPDF.cxx:26
 TKDPDF.cxx:27
 TKDPDF.cxx:28
 TKDPDF.cxx:29
 TKDPDF.cxx:30
 TKDPDF.cxx:31
 TKDPDF.cxx:32
 TKDPDF.cxx:33
 TKDPDF.cxx:34
 TKDPDF.cxx:35
 TKDPDF.cxx:36
 TKDPDF.cxx:37
 TKDPDF.cxx:38
 TKDPDF.cxx:39
 TKDPDF.cxx:40
 TKDPDF.cxx:41
 TKDPDF.cxx:42
 TKDPDF.cxx:43
 TKDPDF.cxx:44
 TKDPDF.cxx:45
 TKDPDF.cxx:46
 TKDPDF.cxx:47
 TKDPDF.cxx:48
 TKDPDF.cxx:49
 TKDPDF.cxx:50
 TKDPDF.cxx:51
 TKDPDF.cxx:52
 TKDPDF.cxx:53
 TKDPDF.cxx:54
 TKDPDF.cxx:55
 TKDPDF.cxx:56
 TKDPDF.cxx:57
 TKDPDF.cxx:58
 TKDPDF.cxx:59
 TKDPDF.cxx:60
 TKDPDF.cxx:61
 TKDPDF.cxx:62
 TKDPDF.cxx:63
 TKDPDF.cxx:64
 TKDPDF.cxx:65
 TKDPDF.cxx:66
 TKDPDF.cxx:67
 TKDPDF.cxx:68
 TKDPDF.cxx:69
 TKDPDF.cxx:70
 TKDPDF.cxx:71
 TKDPDF.cxx:72
 TKDPDF.cxx:73
 TKDPDF.cxx:74
 TKDPDF.cxx:75
 TKDPDF.cxx:76
 TKDPDF.cxx:77
 TKDPDF.cxx:78
 TKDPDF.cxx:79
 TKDPDF.cxx:80
 TKDPDF.cxx:81
 TKDPDF.cxx:82
 TKDPDF.cxx:83
 TKDPDF.cxx:84
 TKDPDF.cxx:85
 TKDPDF.cxx:86
 TKDPDF.cxx:87
 TKDPDF.cxx:88
 TKDPDF.cxx:89
 TKDPDF.cxx:90
 TKDPDF.cxx:91
 TKDPDF.cxx:92
 TKDPDF.cxx:93
 TKDPDF.cxx:94
 TKDPDF.cxx:95
 TKDPDF.cxx:96
 TKDPDF.cxx:97
 TKDPDF.cxx:98
 TKDPDF.cxx:99
 TKDPDF.cxx:100
 TKDPDF.cxx:101
 TKDPDF.cxx:102
 TKDPDF.cxx:103
 TKDPDF.cxx:104
 TKDPDF.cxx:105
 TKDPDF.cxx:106
 TKDPDF.cxx:107
 TKDPDF.cxx:108
 TKDPDF.cxx:109
 TKDPDF.cxx:110
 TKDPDF.cxx:111
 TKDPDF.cxx:112
 TKDPDF.cxx:113
 TKDPDF.cxx:114
 TKDPDF.cxx:115
 TKDPDF.cxx:116
 TKDPDF.cxx:117
 TKDPDF.cxx:118
 TKDPDF.cxx:119
 TKDPDF.cxx:120
 TKDPDF.cxx:121
 TKDPDF.cxx:122
 TKDPDF.cxx:123
 TKDPDF.cxx:124
 TKDPDF.cxx:125
 TKDPDF.cxx:126
 TKDPDF.cxx:127
 TKDPDF.cxx:128
 TKDPDF.cxx:129
 TKDPDF.cxx:130
 TKDPDF.cxx:131
 TKDPDF.cxx:132
 TKDPDF.cxx:133
 TKDPDF.cxx:134
 TKDPDF.cxx:135
 TKDPDF.cxx:136
 TKDPDF.cxx:137
 TKDPDF.cxx:138
 TKDPDF.cxx:139
 TKDPDF.cxx:140
 TKDPDF.cxx:141
 TKDPDF.cxx:142
 TKDPDF.cxx:143
 TKDPDF.cxx:144
 TKDPDF.cxx:145
 TKDPDF.cxx:146
 TKDPDF.cxx:147
 TKDPDF.cxx:148
 TKDPDF.cxx:149
 TKDPDF.cxx:150
 TKDPDF.cxx:151
 TKDPDF.cxx:152
 TKDPDF.cxx:153
 TKDPDF.cxx:154
 TKDPDF.cxx:155
 TKDPDF.cxx:156
 TKDPDF.cxx:157
 TKDPDF.cxx:158
 TKDPDF.cxx:159
 TKDPDF.cxx:160
 TKDPDF.cxx:161
 TKDPDF.cxx:162
 TKDPDF.cxx:163
 TKDPDF.cxx:164
 TKDPDF.cxx:165
 TKDPDF.cxx:166
 TKDPDF.cxx:167
 TKDPDF.cxx:168
 TKDPDF.cxx:169
 TKDPDF.cxx:170
 TKDPDF.cxx:171
 TKDPDF.cxx:172
 TKDPDF.cxx:173
 TKDPDF.cxx:174
 TKDPDF.cxx:175
 TKDPDF.cxx:176
 TKDPDF.cxx:177
 TKDPDF.cxx:178
 TKDPDF.cxx:179
 TKDPDF.cxx:180
 TKDPDF.cxx:181
 TKDPDF.cxx:182
 TKDPDF.cxx:183
 TKDPDF.cxx:184
 TKDPDF.cxx:185
 TKDPDF.cxx:186
 TKDPDF.cxx:187
 TKDPDF.cxx:188
 TKDPDF.cxx:189
 TKDPDF.cxx:190
 TKDPDF.cxx:191
 TKDPDF.cxx:192
 TKDPDF.cxx:193
 TKDPDF.cxx:194