Logo ROOT   6.08/07
Reference Guide
BinarySearchTree.cxx
Go to the documentation of this file.
1 // @(#)root/tmva $Id$
2 // Author: Andreas Hoecker, Joerg Stelzer, Helge Voss, Kai Voss
3 
4 /**********************************************************************************
5  * Project: TMVA - a Root-integrated toolkit for multivariate data analysis *
6  * Package: TMVA *
7  * Class : BinarySearchTree *
8  * Web : http://tmva.sourceforge.net *
9  * *
10  * Description: *
11  * Implementation (see header file for description) *
12  * *
13  * Authors (alphabetical): *
14  * Andreas Hoecker <Andreas.Hocker@cern.ch> - CERN, Switzerland *
15  * Joerg Stelzer <stelzer@cern.ch> - DESY, Germany *
16  * Helge Voss <Helge.Voss@cern.ch> - MPI-K Heidelberg, Germany *
17  * Kai Voss <Kai.Voss@cern.ch> - U. of Victoria, Canada *
18  * *
19  * Copyright (c) 2005: *
20  * CERN, Switzerland *
21  * U. of Victoria, Canada *
22  * MPI-K Heidelberg, Germany *
23  * LAPP, Annecy, France *
24  * *
25  * Redistribution and use in source and binary forms, with or without *
26  * modification, are permitted according to the terms listed in LICENSE *
27  * (http://tmva.sourceforge.net/LICENSE) *
28  * *
29  **********************************************************************************/
30 
31 //////////////////////////////////////////////////////////////////////////
32 // //
33 // BinarySearchTree //
34 // //
35 // A simple Binary search tree including a volume search method //
36 // //
37 //////////////////////////////////////////////////////////////////////////
38 
39 #include <stdexcept>
40 #include <cstdlib>
41 #include <queue>
42 #include <algorithm>
43 
44 // #if ROOT_VERSION_CODE >= 364802
45 // #ifndef ROOT_TMathBase
46 // #include "TMathBase.h"
47 // #endif
48 // #else
49 // #ifndef ROOT_TMath
50 #include "TMath.h"
51 // #endif
52 // #endif
53 
54 #include "TMatrixDBase.h"
55 #include "TObjString.h"
56 #include "TTree.h"
57 
58 #ifndef ROOT_TMVA_MsgLogger
59 #include "TMVA/MsgLogger.h"
60 #endif
61 #ifndef ROOT_TMVA_MethodBase
62 #include "TMVA/MethodBase.h"
63 #endif
64 #ifndef ROOT_TMVA_Tools
65 #include "TMVA/Tools.h"
66 #endif
67 #ifndef ROOT_TMVA_Event
68 #include "TMVA/Event.h"
69 #endif
70 #ifndef ROOT_TMVA_BinarySearchTree
71 #include "TMVA/BinarySearchTree.h"
72 #endif
73 
74 #include "TMVA/BinaryTree.h"
75 #include "TMVA/Types.h"
76 #include "TMVA/Node.h"
77 
79 
80 ////////////////////////////////////////////////////////////////////////////////
81 /// default constructor
82 
84 BinaryTree(),
85  fPeriod ( 1 ),
86  fCurrentDepth( 0 ),
87  fStatisticsIsValid( kFALSE ),
88  fSumOfWeights( 0 ),
89  fCanNormalize( kFALSE )
90 {
91  fNEventsW[0]=fNEventsW[1]=0.;
92 }
93 
94 ////////////////////////////////////////////////////////////////////////////////
95 /// copy constructor that creates a true copy, i.e. a completely independent tree
96 
98  : BinaryTree(),
99  fPeriod ( b.fPeriod ),
100  fCurrentDepth( 0 ),
101  fStatisticsIsValid( kFALSE ),
102  fSumOfWeights( b.fSumOfWeights ),
103  fCanNormalize( kFALSE )
104 {
105  fNEventsW[0]=fNEventsW[1]=0.;
106  Log() << kFATAL << " Copy constructor not implemented yet " << Endl;
107 }
108 
109 ////////////////////////////////////////////////////////////////////////////////
110 /// destructor
111 
113 {
114  for(std::vector< std::pair<Double_t, const TMVA::Event*> >::iterator pIt = fNormalizeTreeTable.begin();
115  pIt != fNormalizeTreeTable.end(); pIt++) {
116  delete pIt->second;
117  }
118 }
119 
120 ////////////////////////////////////////////////////////////////////////////////
121 /// re-create a new tree (decision tree or search tree) from XML
122 
124  std::string type("");
125  gTools().ReadAttr(node,"type", type);
127  bt->ReadXML( node, tmva_Version_Code );
128  return bt;
129 }
130 
131 ////////////////////////////////////////////////////////////////////////////////
132 /// insert a new "event" in the binary tree
133 
135 {
136  fCurrentDepth=0;
138 
139  if (this->GetRoot() == NULL) { // If the list is empty...
140  this->SetRoot( new BinarySearchTreeNode(event)); //Make the new node the root.
141  // have to use "s" for start as "r" for "root" would be the same as "r" for "right"
142  this->GetRoot()->SetPos('s');
143  this->GetRoot()->SetDepth(0);
144  fNNodes = 1;
145  fSumOfWeights = event->GetWeight();
146  ((BinarySearchTreeNode*)this->GetRoot())->SetSelector((UInt_t)0);
147  this->SetPeriode(event->GetNVariables());
148  }
149  else {
150  // sanity check:
151  if (event->GetNVariables() != (UInt_t)this->GetPeriode()) {
152  Log() << kFATAL << "<Insert> event vector length != Periode specified in Binary Tree" << Endl
153  << "--- event size: " << event->GetNVariables() << " Periode: " << this->GetPeriode() << Endl
154  << "--- and all this when trying filling the "<<fNNodes+1<<"th Node" << Endl;
155  }
156  // insert a new node at the propper position
157  this->Insert(event, this->GetRoot());
158  }
159 
160  // normalise the tree to speed up searches
161  if (fCanNormalize) fNormalizeTreeTable.push_back( std::make_pair(0.0,new const Event(*event)) );
162 }
163 
164 ////////////////////////////////////////////////////////////////////////////////
165 /// private internal function to insert a event (node) at the proper position
166 
168  Node *node )
169 {
170  fCurrentDepth++;
172 
173  if (node->GoesLeft(*event)){ // If the adding item is less than the current node's data...
174  if (node->GetLeft() != NULL){ // If there is a left node...
175  // Add the new event to the left node
176  this->Insert(event, node->GetLeft());
177  }
178  else { // If there is not a left node...
179  // Make the new node for the new event
180  BinarySearchTreeNode* current = new BinarySearchTreeNode(event);
181  fNNodes++;
182  fSumOfWeights += event->GetWeight();
183  current->SetSelector(fCurrentDepth%((Int_t)event->GetNVariables()));
184  current->SetParent(node); // Set the new node's previous node.
185  current->SetPos('l');
186  current->SetDepth( node->GetDepth() + 1 );
187  node->SetLeft(current); // Make it the left node of the current one.
188  }
189  }
190  else if (node->GoesRight(*event)) { // If the adding item is less than or equal to the current node's data...
191  if (node->GetRight() != NULL) { // If there is a right node...
192  // Add the new node to it.
193  this->Insert(event, node->GetRight());
194  }
195  else { // If there is not a right node...
196  // Make the new node.
197  BinarySearchTreeNode* current = new BinarySearchTreeNode(event);
198  fNNodes++;
199  fSumOfWeights += event->GetWeight();
200  current->SetSelector(fCurrentDepth%((Int_t)event->GetNVariables()));
201  current->SetParent(node); // Set the new node's previous node.
202  current->SetPos('r');
203  current->SetDepth( node->GetDepth() + 1 );
204  node->SetRight(current); // Make it the left node of the current one.
205  }
206  }
207  else Log() << kFATAL << "<Insert> neither left nor right :)" << Endl;
208 }
209 
210 ////////////////////////////////////////////////////////////////////////////////
211 ///search the tree to find the node matching "event"
212 
214 {
215  return this->Search( event, this->GetRoot() );
216 }
217 
218 ////////////////////////////////////////////////////////////////////////////////
219 /// Private, recursive, function for searching.
220 
222 {
223  if (node != NULL) { // If the node is not NULL...
224  // If we have found the node...
225  if (((BinarySearchTreeNode*)(node))->EqualsMe(*event))
226  return (BinarySearchTreeNode*)node; // Return it
227  if (node->GoesLeft(*event)) // If the node's data is greater than the search item...
228  return this->Search(event, node->GetLeft()); //Search the left node.
229  else //If the node's data is less than the search item...
230  return this->Search(event, node->GetRight()); //Search the right node.
231  }
232  else return NULL; //If the node is NULL, return NULL.
233 }
234 
235 ////////////////////////////////////////////////////////////////////////////////
236 ///return the sum of event (node) weights
237 
239 {
240  if (fSumOfWeights <= 0) {
241  Log() << kWARNING << "you asked for the SumOfWeights, which is not filled yet"
242  << " I call CalcStatistics which hopefully fixes things"
243  << Endl;
244  }
245  if (fSumOfWeights <= 0) Log() << kFATAL << " Zero events in your Search Tree" <<Endl;
246 
247  return fSumOfWeights;
248 }
249 
250 ////////////////////////////////////////////////////////////////////////////////
251 ///return the sum of event (node) weights
252 
254 {
255  if (fSumOfWeights <= 0) {
256  Log() << kWARNING << "you asked for the SumOfWeights, which is not filled yet"
257  << " I call CalcStatistics which hopefully fixes things"
258  << Endl;
259  }
260  if (fSumOfWeights <= 0) Log() << kFATAL << " Zero events in your Search Tree" <<Endl;
261 
262  return fNEventsW[ ( theType == Types::kSignal) ? 0 : 1 ];
263 }
264 
265 ////////////////////////////////////////////////////////////////////////////////
266 /// create the search tree from the event collection
267 /// using ONLY the variables specified in "theVars"
268 
269 Double_t TMVA::BinarySearchTree::Fill( const std::vector<Event*>& events, const std::vector<Int_t>& theVars,
270  Int_t theType )
271 {
272  fPeriod = theVars.size();
273  return Fill(events, theType);
274 }
275 
276 ////////////////////////////////////////////////////////////////////////////////
277 /// create the search tree from the events in a TTree
278 /// using ALL the variables specified included in the Event
279 
280 Double_t TMVA::BinarySearchTree::Fill( const std::vector<Event*>& events, Int_t theType )
281 {
282  UInt_t n=events.size();
283 
284  UInt_t nevents = 0;
285  if (fSumOfWeights != 0) {
286  Log() << kWARNING
287  << "You are filling a search three that is not empty.. "
288  << " do you know what you are doing?"
289  << Endl;
290  }
291  for (UInt_t ievt=0; ievt<n; ievt++) {
292  // insert event into binary tree
293  if (theType == -1 || (Int_t(events[ievt]->GetClass()) == theType) ) {
294  this->Insert( events[ievt] );
295  nevents++;
296  fSumOfWeights += events[ievt]->GetWeight();
297  }
298  } // end of event loop
299  CalcStatistics(0);
300 
301  return fSumOfWeights;
302 }
303 
304 ////////////////////////////////////////////////////////////////////////////////
305 
306 void TMVA::BinarySearchTree::NormalizeTree ( std::vector< std::pair<Double_t, const TMVA::Event*> >::iterator leftBound,
307  std::vector< std::pair<Double_t, const TMVA::Event*> >::iterator rightBound,
308  UInt_t actDim )
309 {
310  // normalises the binary-search tree to reduce the branch length and hence speed up the
311  // search procedure (on average)
312  if (leftBound == rightBound) return;
313 
314  if (actDim == fPeriod) actDim = 0;
315  for (std::vector< std::pair<Double_t, const TMVA::Event*> >::iterator i=leftBound; i!=rightBound; i++) {
316  i->first = i->second->GetValue( actDim );
317  }
318 
319  std::sort( leftBound, rightBound );
320 
321  std::vector< std::pair<Double_t, const TMVA::Event*> >::iterator leftTemp = leftBound;
322  std::vector< std::pair<Double_t, const TMVA::Event*> >::iterator rightTemp = rightBound;
323 
324  // meet in the middle
325  while (true) {
326  rightTemp--;
327  if (rightTemp == leftTemp ) {
328  break;
329  }
330  leftTemp++;
331  if (leftTemp == rightTemp) {
332  break;
333  }
334  }
335 
336  std::vector< std::pair<Double_t, const TMVA::Event*> >::iterator mid = leftTemp;
337  std::vector< std::pair<Double_t, const TMVA::Event*> >::iterator midTemp = mid;
338 
339  if (mid!=leftBound) midTemp--;
340 
341  while (mid != leftBound && mid->second->GetValue( actDim ) == midTemp->second->GetValue( actDim )) {
342  mid--;
343  midTemp--;
344  }
345 
346  Insert( mid->second );
347 
348  // Print(std::cout);
349  // std::cout << std::endl << std::endl;
350 
351  NormalizeTree( leftBound, mid, actDim+1 );
352  mid++;
353  // Print(std::cout);
354  // std::cout << std::endl << std::endl;
355  NormalizeTree( mid, rightBound, actDim+1 );
356 
357 
358  return;
359 }
360 
361 ////////////////////////////////////////////////////////////////////////////////
362 /// Normalisation of tree
363 
365 {
366  SetNormalize( kFALSE );
367  Clear( NULL );
368  this->SetRoot(NULL);
370 }
371 
372 ////////////////////////////////////////////////////////////////////////////////
373 /// clear nodes
374 
376 {
377  BinarySearchTreeNode* currentNode = (BinarySearchTreeNode*)(n == NULL ? this->GetRoot() : n);
378 
379  if (currentNode->GetLeft() != 0) Clear( currentNode->GetLeft() );
380  if (currentNode->GetRight() != 0) Clear( currentNode->GetRight() );
381 
382  if (n != NULL) delete n;
383 
384  return;
385 }
386 
387 ////////////////////////////////////////////////////////////////////////////////
388 /// search the whole tree and add up all weigths of events that
389 /// lie within the given voluem
390 
392  std::vector<const BinarySearchTreeNode*>* events )
393 {
394  return SearchVolume( this->GetRoot(), volume, 0, events );
395 }
396 
397 ////////////////////////////////////////////////////////////////////////////////
398 /// recursively walk through the daughter nodes and add up all weigths of events that
399 /// lie within the given volume
400 
402  std::vector<const BinarySearchTreeNode*>* events )
403 {
404  if (t==NULL) return 0; // Are we at an outer leave?
405 
407 
408  Double_t count = 0.0;
409  if (InVolume( st->GetEventV(), volume )) {
410  count += st->GetWeight();
411  if (NULL != events) events->push_back( st );
412  }
413  if (st->GetLeft()==NULL && st->GetRight()==NULL) {
414 
415  return count; // Are we at an outer leave?
416  }
417 
418  Bool_t tl, tr;
419  Int_t d = depth%this->GetPeriode();
420  if (d != st->GetSelector()) {
421  Log() << kFATAL << "<SearchVolume> selector in Searchvolume "
422  << d << " != " << "node "<< st->GetSelector() << Endl;
423  }
424  tl = (*(volume->fLower))[d] < st->GetEventV()[d]; // Should we descend left?
425  tr = (*(volume->fUpper))[d] >= st->GetEventV()[d]; // Should we descend right?
426 
427  if (tl) count += SearchVolume( st->GetLeft(), volume, (depth+1), events );
428  if (tr) count += SearchVolume( st->GetRight(), volume, (depth+1), events );
429 
430  return count;
431 }
432 
433 Bool_t TMVA::BinarySearchTree::InVolume(const std::vector<Float_t>& event, Volume* volume ) const
434 {
435  // test if the data points are in the given volume
436 
437  Bool_t result = false;
438  for (UInt_t ivar=0; ivar< fPeriod; ivar++) {
439  result = ( (*(volume->fLower))[ivar] < event[ivar] &&
440  (*(volume->fUpper))[ivar] >= event[ivar] );
441  if (!result) break;
442  }
443  return result;
444 }
445 
446 ////////////////////////////////////////////////////////////////////////////////
447 /// calculate basic statistics (mean, rms for each variable)
448 
450 {
451  if (fStatisticsIsValid) return;
452 
453  BinarySearchTreeNode * currentNode = (BinarySearchTreeNode*)n;
454 
455  // default, start at the tree top, then descend recursively
456  if (n == NULL) {
457  fSumOfWeights = 0;
458  for (Int_t sb=0; sb<2; sb++) {
459  fNEventsW[sb] = 0;
460  fMeans[sb] = std::vector<Float_t>(fPeriod);
461  fRMS[sb] = std::vector<Float_t>(fPeriod);
462  fMin[sb] = std::vector<Float_t>(fPeriod);
463  fMax[sb] = std::vector<Float_t>(fPeriod);
464  fSum[sb] = std::vector<Double_t>(fPeriod);
465  fSumSq[sb] = std::vector<Double_t>(fPeriod);
466  for (UInt_t j=0; j<fPeriod; j++) {
467  fMeans[sb][j] = fRMS[sb][j] = fSum[sb][j] = fSumSq[sb][j] = 0;
468  fMin[sb][j] = FLT_MAX;
469  fMax[sb][j] = -FLT_MAX;
470  }
471  }
472  currentNode = (BinarySearchTreeNode*) this->GetRoot();
473  if (currentNode == NULL) return; // no root-node
474  }
475 
476  const std::vector<Float_t> & evtVec = currentNode->GetEventV();
477  Double_t weight = currentNode->GetWeight();
478  // Int_t type = currentNode->IsSignal();
479  // Int_t type = currentNode->IsSignal() ? 0 : 1;
480  Int_t type = Int_t(currentNode->GetClass())== Types::kSignal ? 0 : 1;
481 
482  fNEventsW[type] += weight;
483  fSumOfWeights += weight;
484 
485  for (UInt_t j=0; j<fPeriod; j++) {
486  Float_t val = evtVec[j];
487  fSum[type][j] += val*weight;
488  fSumSq[type][j] += val*val*weight;
489  if (val < fMin[type][j]) fMin[type][j] = val;
490  if (val > fMax[type][j]) fMax[type][j] = val;
491  }
492 
493  if ( (currentNode->GetLeft() != NULL) ) CalcStatistics( currentNode->GetLeft() );
494  if ( (currentNode->GetRight() != NULL) ) CalcStatistics( currentNode->GetRight() );
495 
496  if (n == NULL) { // i.e. the root node
497  for (Int_t sb=0; sb<2; sb++) {
498  for (UInt_t j=0; j<fPeriod; j++) {
499  if (fNEventsW[sb] == 0) { fMeans[sb][j] = fRMS[sb][j] = 0; continue; }
500  fMeans[sb][j] = fSum[sb][j]/fNEventsW[sb];
501  fRMS[sb][j] = TMath::Sqrt(fSumSq[sb][j]/fNEventsW[sb] - fMeans[sb][j]*fMeans[sb][j]);
502  }
503  }
505  }
506 
507  return;
508 }
509 
510 Int_t TMVA::BinarySearchTree::SearchVolumeWithMaxLimit( Volume *volume, std::vector<const BinarySearchTreeNode*>* events,
511  Int_t max_points )
512 {
513  // recursively walk through the daughter nodes and add up all weigths of events that
514  // lie within the given volume a maximum number of events can be given
515  if (this->GetRoot() == NULL) return 0; // Are we at an outer leave?
516 
517  std::queue< std::pair< const BinarySearchTreeNode*, Int_t > > queue;
518  std::pair< const BinarySearchTreeNode*, Int_t > st = std::make_pair( (const BinarySearchTreeNode*)this->GetRoot(), 0 );
519  queue.push( st );
520 
521  Int_t count = 0;
522 
523  while ( !queue.empty() ) {
524  st = queue.front(); queue.pop();
525 
526  if (count == max_points)
527  return count;
528 
529  if (InVolume( st.first->GetEventV(), volume )) {
530  count++;
531  if (NULL != events) events->push_back( st.first );
532  }
533 
534  Bool_t tl, tr;
535  Int_t d = st.second;
536  if ( d == Int_t(this->GetPeriode()) ) d = 0;
537 
538  if (d != st.first->GetSelector()) {
539  Log() << kFATAL << "<SearchVolume> selector in Searchvolume "
540  << d << " != " << "node "<< st.first->GetSelector() << Endl;
541  }
542 
543  tl = (*(volume->fLower))[d] < st.first->GetEventV()[d] && st.first->GetLeft() != NULL; // Should we descend left?
544  tr = (*(volume->fUpper))[d] >= st.first->GetEventV()[d] && st.first->GetRight() != NULL; // Should we descend right?
545 
546  if (tl) queue.push( std::make_pair( (const BinarySearchTreeNode*)st.first->GetLeft(), d+1 ) );
547  if (tr) queue.push( std::make_pair( (const BinarySearchTreeNode*)st.first->GetRight(), d+1 ) );
548  }
549 
550  return count;
551 }
Int_t SearchVolumeWithMaxLimit(TMVA::Volume *, std::vector< const TMVA::BinarySearchTreeNode *> *events=0, Int_t=-1)
std::vector< Double_t > * fLower
Definition: Volume.h:75
std::vector< Double_t > fSum[2]
BinarySearchTreeNode * Search(Event *event) const
search the tree to find the node matching "event"
std::vector< Float_t > fMax[2]
std::vector< Double_t > fSumSq[2]
MsgLogger & Endl(MsgLogger &ml)
Definition: MsgLogger.h:162
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t mid
Definition: TRolke.cxx:630
virtual ~BinarySearchTree(void)
destructor
Bool_t InVolume(const std::vector< Float_t > &, Volume *) const
float Float_t
Definition: RtypesCore.h:53
Double_t Fill(const std::vector< TMVA::Event *> &events, const std::vector< Int_t > &theVars, Int_t theType=-1)
create the search tree from the event collection using ONLY the variables specified in "theVars" ...
std::vector< Double_t > * fUpper
Definition: Volume.h:76
UInt_t GetDepth() const
Definition: Node.h:118
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
virtual void SetRight(Node *r)
Definition: Node.h:97
const Bool_t kFALSE
Definition: Rtypes.h:92
Double_t GetSumOfWeights(void) const
return the sum of event (node) weights
UInt_t GetPeriode(void) const
void NormalizeTree()
Normalisation of tree.
void SetDepth(UInt_t d)
Definition: Node.h:115
std::vector< std::pair< Double_t, const TMVA::Event * > > fNormalizeTreeTable
TClass * GetClass(T *)
Definition: TClass.h:555
Tools & gTools()
Definition: Tools.cxx:79
MsgLogger & Log() const
Definition: BinaryTree.cxx:236
virtual Bool_t GoesLeft(const Event &) const =0
virtual void SetLeft(Node *l)
Definition: Node.h:96
virtual Node * GetRight() const
Definition: Node.h:92
virtual Node * GetLeft() const
Definition: Node.h:91
void Clear(TMVA::Node *n=0)
clear nodes
void CalcStatistics(TMVA::Node *n=0)
calculate basic statistics (mean, rms for each variable)
BinarySearchTree(void)
default constructor
static BinarySearchTree * CreateFromXML(void *node, UInt_t tmva_Version_Code=TMVA_VERSION_CODE)
re-create a new tree (decision tree or search tree) from XML
unsigned int UInt_t
Definition: RtypesCore.h:42
void SetNormalize(Bool_t norm)
void ReadAttr(void *node, const char *, T &value)
Definition: Tools.h:296
void Insert(const Event *)
insert a new "event" in the binary tree
virtual void ReadXML(void *node, UInt_t tmva_Version_Code=TMVA_VERSION_CODE)
read attributes from XML
Definition: BinaryTree.cxx:145
virtual void SetParent(Node *p)
Definition: Node.h:98
UInt_t GetNVariables() const
accessor to the number of variables
Definition: Event.cxx:305
#define ClassImp(name)
Definition: Rtypes.h:279
double Double_t
Definition: RtypesCore.h:55
int type
Definition: TGX11.cxx:120
std::vector< Float_t > fMeans[2]
const std::vector< Float_t > & GetEventV() const
Abstract ClassifierFactory template that handles arbitrary types.
virtual Bool_t GoesRight(const Event &) const =0
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
#define NULL
Definition: Rtypes.h:82
virtual Node * GetRoot() const
Definition: BinaryTree.h:89
double result[121]
void SetRoot(Node *r)
Definition: BinaryTree.h:86
Double_t Sqrt(Double_t x)
Definition: TMath.h:464
std::vector< Float_t > fMin[2]
Double_t SearchVolume(Volume *, std::vector< const TMVA::BinarySearchTreeNode *> *events=0)
search the whole tree and add up all weigths of events that lie within the given voluem ...
std::vector< Float_t > fRMS[2]
const Bool_t kTRUE
Definition: Rtypes.h:91
const Int_t n
Definition: legend1.C:16
void SetPos(char s)
Definition: Node.h:121