Logo ROOT   6.08/07
Reference Guide
PDEFoam.cxx
Go to the documentation of this file.
1 // @(#)root/tmva $Id$
2 // Author: S.Jadach, Tancredi Carli, Dominik Dannheim, Alexander Voigt
3 
4 /**********************************************************************************
5  * Project: TMVA - a Root-integrated toolkit for multivariate data analysis *
6  * Package: TMVA *
7  * Classes: PDEFoam *
8  * Web : http://tmva.sourceforge.net *
9  * *
10  * Description: *
11  * Implementations *
12  * *
13  * Authors (alphabetical): *
14  * Tancredi Carli - CERN, Switzerland *
15  * Dominik Dannheim - CERN, Switzerland *
16  * S. Jadach - Institute of Nuclear Physics, Cracow, Poland *
17  * Alexander Voigt - TU Dresden, Germany *
18  * Peter Speckmayer - CERN, Switzerland *
19  * *
20  * Copyright (c) 2008, 2010: *
21  * CERN, Switzerland *
22  * MPI-K Heidelberg, Germany *
23  * *
24  * Redistribution and use in source and binary forms, with or without *
25  * modification, are permitted according to the terms listed in LICENSE *
26  * (http://tmva.sourceforge.net/LICENSE) *
27  **********************************************************************************/
28 
29 //_____________________________________________________________________
30 //
31 // Implementation of PDEFoam
32 //
33 // The PDEFoam method is an extension of the PDERS method, which uses
34 // self-adapting binning to divide the multi-dimensional phase space
35 // in a finite number of hyper-rectangles (boxes).
36 //
37 // For a given number of boxes, the binning algorithm adjusts the size
38 // and position of the boxes inside the multidimensional phase space,
39 // minimizing the variance of the signal and background densities inside
40 // the boxes. The binned density information is stored in binary trees,
41 // allowing for a very fast and memory-efficient classification of
42 // events.
43 //
44 // The implementation of the PDEFoam is based on the monte-carlo
45 // integration package TFoam included in the analysis package ROOT.
46 //
47 // The class TMVA::PDEFoam defines the default interface for the
48 // PDEFoam variants:
49 //
50 // - PDEFoamEvent
51 // - PDEFoamDiscriminant
52 // - PDEFoamTarget
53 // - PDEFoamMultiTarget
54 // - PDEFoamDecisionTree
55 //
56 // Per default PDEFoam stores in the cells the number of events (event
57 // weights) and therefore acts as an event density estimator.
58 // However, the above listed derived classes override this behaviour
59 // to implement certain PDEFoam variations.
60 //
61 // In order to use PDEFoam the user has to set the density estimator
62 // of the type TMVA::PDEFoamDensityBase, which is used to during the foam
63 // build-up. The default PDEFoam should be used with
64 // PDEFoamEventDensity.
65 // _____________________________________________________________________
66 
67 #include "TMVA/PDEFoam.h"
68 
69 #include "TMVA/Event.h"
70 #include "TMVA/MsgLogger.h"
71 #include "TMVA/PDEFoamKernelBase.h"
72 #include "TMVA/Timer.h"
73 #include "TMVA/Tools.h"
74 #include "TMVA/Types.h"
75 
76 #include "TStyle.h"
77 #include "TObject.h"
78 #include "TH1D.h"
79 #include "TMath.h"
80 #include "TVectorT.h"
81 #include "TRandom3.h"
82 #include "TColor.h"
83 #include "TDirectory.h"
84 #include "TObjArray.h"
85 
86 #include <cassert>
87 #include <fstream>
88 #include <iostream>
89 #include <iomanip>
90 #include <limits>
91 #include <sstream>
92 
94 
95 static const Float_t kHigh= FLT_MAX;
96 static const Float_t kVlow=-FLT_MAX;
97 
98 using namespace std;
99 
100 ////////////////////////////////////////////////////////////////////////////////
101 /// Default constructor for streamer, user should not use it.
102 
104  fName("PDEFoam"),
105  fDim(0),
106  fNCells(0),
107  fNBin(5),
108  fNSampl(2000),
109  fEvPerBin(0),
110  fMaskDiv(0),
111  fInhiDiv(0),
112  fNoAct(1),
113  fLastCe(-1),
114  fCells(0),
115  fHistEdg(0),
116  fRvec(0),
117  fPseRan(new TRandom3(4356)),
118  fAlpha(0),
119  fFoamType(kSeparate),
120  fXmin(0),
121  fXmax(0),
122  fNElements(0),
123  fNmin(100),
124  fMaxDepth(0),
125  fVolFrac(1.0/15.0),
126  fFillFoamWithOrigWeights(kFALSE),
127  fDTSeparation(kFoam),
128  fPeekMax(kTRUE),
129  fDistr(NULL),
130  fTimer(new Timer(0, "PDEFoam", kTRUE)),
131  fVariableNames(new TObjArray()),
132  fLogger(new MsgLogger("PDEFoam"))
133 {
134  // fVariableNames may delete it's heap-based content
135  if (fVariableNames)
136  fVariableNames->SetOwner(kTRUE);
137 }
138 
139 ////////////////////////////////////////////////////////////////////////////////
140 /// User constructor, to be employed by the user
141 
142 TMVA::PDEFoam::PDEFoam(const TString& name) :
143  fName(name),
144  fDim(0),
145  fNCells(1000),
146  fNBin(5),
147  fNSampl(2000),
148  fEvPerBin(0),
149  fMaskDiv(0),
150  fInhiDiv(0),
151  fNoAct(1),
152  fLastCe(-1),
153  fCells(0),
154  fHistEdg(0),
155  fRvec(0),
156  fPseRan(new TRandom3(4356)),
157  fAlpha(0),
158  fFoamType(kSeparate),
159  fXmin(0),
160  fXmax(0),
161  fNElements(0),
162  fNmin(100),
163  fMaxDepth(0),
164  fVolFrac(1.0/15.0),
165  fFillFoamWithOrigWeights(kFALSE),
166  fDTSeparation(kFoam),
167  fPeekMax(kTRUE),
168  fDistr(NULL),
169  fTimer(new Timer(1, "PDEFoam", kTRUE)),
170  fVariableNames(new TObjArray()),
171  fLogger(new MsgLogger("PDEFoam"))
172 {
173  if(strlen(name) > 128)
174  Log() << kFATAL << "Name too long " << name.Data() << Endl;
175 
176  // fVariableNames may delete it's heap-based content
177  if (fVariableNames)
179 }
180 
181 ////////////////////////////////////////////////////////////////////////////////
182 /// Default destructor
183 
185 {
186  delete fVariableNames;
187  delete fTimer;
188  if (fDistr) delete fDistr;
189  if (fPseRan) delete fPseRan;
190  if (fXmin) { delete [] fXmin; fXmin=0; }
191  if (fXmax) { delete [] fXmax; fXmax=0; }
192 
194  if(fCells!= 0) {
195  for(Int_t i=0; i<fNCells; i++) delete fCells[i]; // PDEFoamCell*[]
196  delete [] fCells;
197  }
198  delete [] fRvec; //double[]
199  delete [] fAlpha; //double[]
200  delete [] fMaskDiv; //int[]
201  delete [] fInhiDiv; //int[]
202 
203  delete fLogger;
204 }
205 
206 ////////////////////////////////////////////////////////////////////////////////
207 /// Copy Constructor NOT IMPLEMENTED (NEVER USED)
208 
210  TObject(from)
211  , fDim(0)
212  , fNCells(0)
213  , fNBin(0)
214  , fNSampl(0)
215  , fEvPerBin(0)
216  , fMaskDiv(0)
217  , fInhiDiv(0)
218  , fNoAct(0)
219  , fLastCe(0)
220  , fCells(0)
221  , fHistEdg(0)
222  , fRvec(0)
223  , fPseRan(0)
224  , fAlpha(0)
226  , fXmin(0)
227  , fXmax(0)
228  , fNElements(0)
229  , fNmin(0)
230  , fMaxDepth(0)
231  , fVolFrac(1.0/15.0)
234  , fPeekMax(kTRUE)
235  , fDistr(0)
236  , fTimer(0)
237  , fVariableNames(0)
238  , fLogger(new MsgLogger(*from.fLogger))
239 {
240  Log() << kFATAL << "COPY CONSTRUCTOR NOT IMPLEMENTED" << Endl;
241 
242  // fVariableNames may delete it's heap-based content
243  if (fVariableNames)
245 }
246 
247 ////////////////////////////////////////////////////////////////////////////////
248 /// Sets dimension of cubical space
249 
251 {
252  if (kDim < 1)
253  Log() << kFATAL << "<SetDim>: Dimension is zero or negative!" << Endl;
254 
255  fDim = kDim;
256  if (fXmin) delete [] fXmin;
257  if (fXmax) delete [] fXmax;
258  fXmin = new Double_t[GetTotDim()];
259  fXmax = new Double_t[GetTotDim()];
260 }
261 
262 ////////////////////////////////////////////////////////////////////////////////
263 /// set lower foam bound in dimension idim
264 
266 {
267  if (idim<0 || idim>=GetTotDim())
268  Log() << kFATAL << "<SetXmin>: Dimension out of bounds!" << Endl;
269 
270  fXmin[idim]=wmin;
271 }
272 
273 ////////////////////////////////////////////////////////////////////////////////
274 /// set upper foam bound in dimension idim
275 
277 {
278  if (idim<0 || idim>=GetTotDim())
279  Log() << kFATAL << "<SetXmax>: Dimension out of bounds!" << Endl;
280 
281  fXmax[idim]=wmax;
282 }
283 
284 ////////////////////////////////////////////////////////////////////////////////
285 /// Basic initialization of FOAM invoked by the user.
286 /// IMPORTANT: Random number generator and the distribution object has to be
287 /// provided using SetPseRan and SetRho prior to invoking this initializator!
288 ///
289 /// After the foam is grown, space for 2 variables is reserved in
290 /// every cell. They are used for filling the foam cells.
291 
293 {
294  Bool_t addStatus = TH1::AddDirectoryStatus();
296 
297  if(fPseRan==0) Log() << kFATAL << "Random number generator not set" << Endl;
298  if(fDistr==0) Log() << kFATAL << "Distribution function not set" << Endl;
299  if(fDim==0) Log() << kFATAL << "Zero dimension not allowed" << Endl;
300 
301  /////////////////////////////////////////////////////////////////////////
302  // ALLOCATE SMALL LISTS //
303  // it is done globally, not for each cell, to save on allocation time //
304  /////////////////////////////////////////////////////////////////////////
305  fRvec = new Double_t[fDim]; // Vector of random numbers
306  if(fRvec==0) Log() << kFATAL << "Cannot initialize buffer fRvec" << Endl;
307 
308  if(fDim>0){
309  fAlpha = new Double_t[fDim]; // sum<1 for internal parametrization of the simplex
310  if(fAlpha==0) Log() << kFATAL << "Cannot initialize buffer fAlpha" << Endl;
311  }
312 
313  //====== List of directions inhibited for division
314  if(fInhiDiv == 0){
315  fInhiDiv = new Int_t[fDim];
316  for(Int_t i=0; i<fDim; i++) fInhiDiv[i]=0;
317  }
318  //====== Dynamic mask used in Explore for edge determination
319  if(fMaskDiv == 0){
320  fMaskDiv = new Int_t[fDim];
321  for(Int_t i=0; i<fDim; i++) fMaskDiv[i]=1;
322  }
323  //====== Initialize list of histograms
324  fHistEdg = new TObjArray(fDim); // Initialize list of histograms
325  for(Int_t i=0; i<fDim; i++){
326  TString hname, htitle;
327  hname = fName+TString("_HistEdge_");
328  hname += i;
329  htitle = TString("Edge Histogram No. ");
330  htitle += i;
331  (*fHistEdg)[i] = new TH1D(hname.Data(),htitle.Data(),fNBin,0.0, 1.0); // Initialize histogram for each edge
332  ((TH1D*)(*fHistEdg)[i])->Sumw2();
333  }
334 
335  // ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| //
336  // BUILD-UP of the FOAM //
337  // ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| //
338 
339  // prepare PDEFoam for growing
340  ResetCellElements(); // reset all cell elements
341 
342  // Define and explore root cell(s)
343  InitCells();
344  Grow();
345 
346  TH1::AddDirectory(addStatus);
347 
348  // prepare PDEFoam for the filling with events
349  ResetCellElements(); // reset all cell elements
350 } // Create
351 
352 ////////////////////////////////////////////////////////////////////////////////
353 /// Internal subprogram used by Create.
354 /// It initializes "root part" of the FOAM of the tree of cells.
355 
357 {
358  fLastCe =-1; // Index of the last cell
359  if(fCells!= 0) {
360  for(Int_t i=0; i<fNCells; i++) delete fCells[i];
361  delete [] fCells;
362  }
363 
364  fCells = new(nothrow) PDEFoamCell*[fNCells];
365  if (!fCells) {
366  Log() << kFATAL << "not enough memory to create " << fNCells
367  << " cells" << Endl;
368  }
369  for(Int_t i=0; i<fNCells; i++){
370  fCells[i]= new PDEFoamCell(fDim); // Allocate BIG list of cells
371  fCells[i]->SetSerial(i);
372  }
373 
374  /////////////////////////////////////////////////////////////////////////////
375  // Single Root Hypercube //
376  /////////////////////////////////////////////////////////////////////////////
377  CellFill(1, 0); // 0-th cell ACTIVE
378 
379  // Exploration of the root cell(s)
380  for(Long_t iCell=0; iCell<=fLastCe; iCell++){
381  Explore( fCells[iCell] ); // Exploration of root cell(s)
382  }
383 }//InitCells
384 
385 ////////////////////////////////////////////////////////////////////////////////
386 /// Internal subprogram used by Create.
387 /// It initializes content of the newly allocated active cell.
388 
390 {
391  PDEFoamCell *cell;
392  if (fLastCe==fNCells){
393  Log() << kFATAL << "Too many cells" << Endl;
394  }
395  fLastCe++; // 0-th cell is the first
396 
397  cell = fCells[fLastCe];
398 
399  cell->Fill(status, parent, 0, 0);
400 
401  cell->SetBest( -1); // pointer for planning division of the cell
402  cell->SetXdiv(0.5); // factor for division
403  Double_t xInt2,xDri2;
404  if(parent!=0){
405  xInt2 = 0.5*parent->GetIntg();
406  xDri2 = 0.5*parent->GetDriv();
407  cell->SetIntg(xInt2);
408  cell->SetDriv(xDri2);
409  }else{
410  cell->SetIntg(0.0);
411  cell->SetDriv(0.0);
412  }
413  return fLastCe;
414 }
415 
416 ////////////////////////////////////////////////////////////////////////////////
417 /// Internal subprogram used by Create.
418 /// It explores newly defined cell with help of special short MC sampling.
419 /// As a result, estimates of kTRUE and drive volume is defined/determined
420 /// Average and dispersion of the weight distribution will is found along
421 /// each edge and the best edge (minimum dispersion, best maximum weight)
422 /// is memorized for future use.
423 /// The optimal division point for eventual future cell division is
424 /// determined/recorded. Recorded are also minimum and maximum weight etc.
425 /// The volume estimate in all (inactive) parent cells is updated.
426 /// Note that links to parents and initial volume = 1/2 parent has to be
427 /// already defined prior to calling this routine.
428 ///
429 /// If fNmin > 0 then the total number of (training) events found in
430 /// the cell during the exploration is stored in the cell. This
431 /// information is used withing PeekMax() to avoid splitting cells
432 /// which contain less than fNmin events.
433 
435 {
436  Double_t wt, dx, xBest=0, yBest;
437  Double_t intOld, driOld;
438 
439  Long_t iev;
440  Double_t nevMC;
441  Int_t i, j, k;
442  Int_t nProj, kBest;
443  Double_t ceSum[5], xproj;
444 
445  Double_t event_density = 0;
446  Double_t totevents = 0;
447  Double_t toteventsOld = 0;
448 
449  PDEFoamVect cellSize(fDim);
450  PDEFoamVect cellPosi(fDim);
451 
452  cell->GetHcub(cellPosi,cellSize);
453 
454  PDEFoamCell *parent;
455 
456  Double_t *xRand = new Double_t[fDim];
457 
458  Double_t *volPart=0;
459 
460  // calculate volume scale
461  Double_t vol_scale = 1.0;
462  for (Int_t idim = 0; idim < fDim; ++idim)
463  vol_scale *= fXmax[idim] - fXmin[idim];
464 
465  cell->CalcVolume();
466  dx = cell->GetVolume() * vol_scale;
467  intOld = cell->GetIntg(); //memorize old values,
468  driOld = cell->GetDriv(); //will be needed for correcting parent cells
469  toteventsOld = GetCellElement(cell, 0);
470 
471  /////////////////////////////////////////////////////
472  // Special Short MC sampling to probe cell //
473  /////////////////////////////////////////////////////
474  ceSum[0]=0;
475  ceSum[1]=0;
476  ceSum[2]=0;
477  ceSum[3]=kHigh; //wtmin
478  ceSum[4]=kVlow; //wtmax
479 
480  for (i=0;i<fDim;i++) ((TH1D *)(*fHistEdg)[i])->Reset(); // Reset histograms
481 
482  Double_t nevEff=0.;
483  // ||||||||||||||||||||||||||BEGIN MC LOOP|||||||||||||||||||||||||||||
484  for (iev=0;iev<fNSampl;iev++){
485  MakeAlpha(); // generate uniformly vector inside hypercube
486 
487  if (fDim>0) for (j=0; j<fDim; j++) xRand[j]= cellPosi[j] +fAlpha[j]*(cellSize[j]);
488 
489  wt = dx*Eval(xRand, event_density);
490  totevents += event_density;
491 
492  nProj = 0;
493  if (fDim>0) {
494  for (k=0; k<fDim; k++) {
495  xproj =fAlpha[k];
496  ((TH1D *)(*fHistEdg)[nProj])->Fill(xproj,wt);
497  nProj++;
498  }
499  }
500 
501  ceSum[0] += wt; // sum of weights
502  ceSum[1] += wt*wt; // sum of weights squared
503  ceSum[2]++; // sum of 1
504  if (ceSum[3]>wt) ceSum[3]=wt; // minimum weight;
505  if (ceSum[4]<wt) ceSum[4]=wt; // maximum weight
506  // test MC loop exit condition
507  if (ceSum[1]>0) nevEff = ceSum[0]*ceSum[0]/ceSum[1];
508  else nevEff = 0;
509  if ( nevEff >= fNBin*fEvPerBin) break;
510  } // ||||||||||||||||||||||||||END MC LOOP|||||||||||||||||||||||||||||
511  totevents *= dx;
512 
513  if (fNSampl>0) totevents /= fNSampl;
514 
515  // make shure that, if root cell is explored, more than zero
516  // events were found.
517  if (cell==fCells[0] && ceSum[0]<=0.0){
518  if (ceSum[0]==0.0)
519  Log() << kFATAL << "No events were found during exploration of "
520  << "root cell. Please check PDEFoam parameters nSampl "
521  << "and VolFrac." << Endl;
522  else
523  Log() << kWARNING << "Negative number of events found during "
524  << "exploration of root cell" << Endl;
525  }
526 
527  //------------------------------------------------------------------
528  //--- predefine logics of searching for the best division edge ---
529  for (k=0; k<fDim;k++){
530  fMaskDiv[k] =1; // default is all
531  if ( fInhiDiv[k]==1) fMaskDiv[k] =0; // inhibit some...
532  }
533  kBest=-1;
534 
535  /////////////////////////////////////////////////////////////////////////////
536 
537  nevMC = ceSum[2];
538  Double_t intTrue = ceSum[0]/(nevMC+0.000001);
539  Double_t intDriv=0.;
540 
541  if (kBest == -1) Varedu(ceSum,kBest,xBest,yBest); // determine the best edge,
542  intDriv =sqrt(ceSum[1]/nevMC) -intTrue; // Foam build-up, sqrt(<w**2>) -<w>
543 
544  //=================================================================================
545  cell->SetBest(kBest);
546  cell->SetXdiv(xBest);
547  cell->SetIntg(intTrue);
548  cell->SetDriv(intDriv);
549  SetCellElement(cell, 0, totevents);
550 
551  // correct/update integrals in all parent cells to the top of the tree
552  Double_t parIntg, parDriv;
553  for (parent = cell->GetPare(); parent!=0; parent = parent->GetPare()){
554  parIntg = parent->GetIntg();
555  parDriv = parent->GetDriv();
556  parent->SetIntg( parIntg +intTrue -intOld );
557  parent->SetDriv( parDriv +intDriv -driOld );
558  SetCellElement( parent, 0, GetCellElement(parent, 0) + totevents - toteventsOld);
559  }
560  delete [] volPart;
561  delete [] xRand;
562 }
563 
564 ////////////////////////////////////////////////////////////////////////////////
565 /// Internal subrogram used by Create.
566 /// In determines the best edge candidate and the position of the cell division plane
567 /// in case of the variance reduction for future cell division,
568 /// using results of the MC exploration run stored in fHistEdg
569 
570 void TMVA::PDEFoam::Varedu(Double_t ceSum[5], Int_t &kBest, Double_t &xBest, Double_t &yBest)
571 {
572  Double_t nent = ceSum[2];
573  // Double_t swAll = ceSum[0];
574  Double_t sswAll = ceSum[1];
575  Double_t ssw = sqrt(sswAll)/sqrt(nent);
576  //
577  Double_t sswIn,sswOut,xLo,xUp;
578  kBest =-1;
579  xBest =0.5;
580  yBest =1.0;
581  Double_t maxGain=0.0;
582  // Now go over all projections kProj
583  for(Int_t kProj=0; kProj<fDim; kProj++) {
584  if( fMaskDiv[kProj]) {
585  // initialize search over bins
586  // Double_t sigmIn =0.0; Double_t sigmOut =0.0;
587  Double_t sswtBest = kHigh;
588  Double_t gain =0.0;
589  Double_t xMin=0.0; Double_t xMax=0.0;
590  // Double loop over all pairs jLo<jUp
591  for(Int_t jLo=1; jLo<=fNBin; jLo++) {
592  Double_t aswIn=0; Double_t asswIn=0;
593  for(Int_t jUp=jLo; jUp<=fNBin;jUp++) {
594  aswIn += ((TH1D *)(*fHistEdg)[kProj])->GetBinContent(jUp);
595  asswIn += Sqr(((TH1D *)(*fHistEdg)[kProj])->GetBinError( jUp));
596  xLo=(jLo-1.0)/fNBin;
597  xUp=(jUp*1.0)/fNBin;
598  // swIn = aswIn/nent;
599  // swOut = (swAll-aswIn)/nent;
600  if ( (xUp-xLo) < std::numeric_limits<double>::epsilon()) sswIn=0.;
601  else sswIn = sqrt(asswIn) /sqrt(nent*(xUp-xLo)) *(xUp-xLo);
602  if ( (1.0-xUp+xLo) < std::numeric_limits<double>::epsilon()) sswOut=0.;
603  else if ( sswAll-asswIn < std::numeric_limits<double>::epsilon()) sswOut=0.;
604  else sswOut= sqrt(sswAll-asswIn)/sqrt(nent*(1.0-xUp+xLo)) *(1.0-xUp+xLo);
605  if( (sswIn+sswOut) < sswtBest) {
606  sswtBest = sswIn+sswOut;
607  gain = ssw-sswtBest;
608  // sigmIn = sswIn -swIn; // Debug
609  // sigmOut = sswOut-swOut; // Debug
610  xMin = xLo;
611  xMax = xUp;
612  }
613  }//jUp
614  }//jLo
615  Int_t iLo = (Int_t) (fNBin*xMin);
616  Int_t iUp = (Int_t) (fNBin*xMax);
617 
618  if(gain>=maxGain) {
619  maxGain=gain;
620  kBest=kProj; // <--- !!!!! The best edge
621  xBest=xMin;
622  yBest=xMax;
623  if(iLo == 0) xBest=yBest; // The best division point
624  if(iUp == fNBin) yBest=xBest; // this is not really used
625  }
626  }
627  } //kProj
628 
629  if( (kBest >= fDim) || (kBest<0) )
630  Log() << kFATAL << "Something wrong with kBest" << Endl;
631 } //PDEFoam::Varedu
632 
633 ////////////////////////////////////////////////////////////////////////////////
634 /// Internal subrogram used by Create.
635 /// Provides random vector Alpha 0< Alpha(i) < 1
636 
638 {
639  // simply generate and load kDim uniform random numbers
640  fPseRan->RndmArray(fDim,fRvec); // kDim random numbers needed
641  for(Int_t k=0; k<fDim; k++) fAlpha[k] = fRvec[k];
642 } //MakeAlpha
643 
644 ////////////////////////////////////////////////////////////////////////////////
645 /// Internal subprogram used by Create. It finds cell with maximal
646 /// driver integral for the purpose of the division. This function
647 /// is overridden by the PDEFoam Class to apply cuts on the number
648 /// of events in the cell (fNmin) and the cell tree depth
649 /// (GetMaxDepth() > 0) during cell buildup.
650 
652 {
653  Long_t iCell = -1;
654 
655  Long_t i;
656  Double_t drivMax, driv, xDiv;
657  Bool_t bCutNmin = kTRUE;
658  Bool_t bCutMaxDepth = kTRUE;
659  // drivMax = kVlow;
660  drivMax = 0; // only split cells if gain>0 (this also avoids splitting at cell boundary)
661  for(i=0; i<=fLastCe; i++) {//without root
662  if( fCells[i]->GetStat() == 1 ) {
663  // if driver integral < numeric limit, skip cell
664  driv = fCells[i]->GetDriv();
666  continue;
667 
668  // do not split cell at the edges
669  xDiv = TMath::Abs(fCells[i]->GetXdiv());
672  continue;
673 
674  // apply cut on depth
675  if (GetMaxDepth() > 0)
676  bCutMaxDepth = fCells[i]->GetDepth() < GetMaxDepth();
677 
678  // apply Nmin-cut
679  if (GetNmin() > 0)
680  bCutNmin = GetCellElement(fCells[i], 0) > GetNmin();
681 
682  // choose cell
683  if(driv > drivMax && bCutNmin && bCutMaxDepth) {
684  drivMax = driv;
685  iCell = i;
686  }
687  }
688  }
689 
690  if (iCell == -1){
691  if (!bCutNmin)
692  Log() << kVERBOSE << "Warning: No cell with more than "
693  << GetNmin() << " events found!" << Endl;
694  else if (!bCutMaxDepth)
695  Log() << kVERBOSE << "Warning: Maximum depth reached: "
696  << GetMaxDepth() << Endl;
697  else
698  Log() << kWARNING << "<PDEFoam::PeekMax>: no more candidate cells (drivMax>0) found for further splitting." << Endl;
699  }
700 
701  return(iCell);
702 }
703 
704 ////////////////////////////////////////////////////////////////////////////////
705 /// Internal subrogram used by Create.
706 /// It divides cell iCell into two daughter cells.
707 /// The iCell is retained and tagged as inactive, daughter cells are appended
708 /// at the end of the buffer.
709 /// New vertex is added to list of vertices.
710 /// List of active cells is updated, iCell removed, two daughters added
711 /// and their properties set with help of MC sampling (PDEFoam_Explore)
712 /// Returns Code RC=-1 of buffer limit is reached, fLastCe=fnBuf.
713 
715 {
716  // Double_t xdiv;
717  Int_t kBest;
718 
719  if(fLastCe+1 >= fNCells) Log() << kFATAL << "Buffer limit is reached, fLastCe=fnBuf" << Endl;
720 
721  cell->SetStat(0); // reset cell as inactive
722  fNoAct++;
723 
724  // xdiv = cell->GetXdiv();
725  kBest = cell->GetBest();
726  if( kBest<0 || kBest>=fDim ) Log() << kFATAL << "Wrong kBest" << Endl;
727 
728  //////////////////////////////////////////////////////////////////
729  // define two daughter cells (active) //
730  //////////////////////////////////////////////////////////////////
731 
732  Int_t d1 = CellFill(1, cell);
733  Int_t d2 = CellFill(1, cell);
734  cell->SetDau0((fCells[d1]));
735  cell->SetDau1((fCells[d2]));
736 
737  Explore( (fCells[d1]) );
738  Explore( (fCells[d2]) );
739 
740  return 1;
741 } // PDEFoam_Divide
742 
743 ////////////////////////////////////////////////////////////////////////////////
744 /// Internal subprogram.
745 /// Evaluates (training) distribution.
746 
748 {
749  // Transform variable xRand, since Foam boundaries are [0,1] and
750  // fDistr is filled with events which range in [fXmin,fXmax]
751  //
752  // Transformation: [0, 1] --> [xmin, xmax]
753  std::vector<Double_t> xvec;
754  xvec.reserve(GetTotDim());
755  for (Int_t idim = 0; idim < GetTotDim(); ++idim)
756  xvec.push_back( VarTransformInvers(idim, xRand[idim]) );
757 
758  return GetDistr()->Density(xvec, event_density);
759 }
760 
761 ////////////////////////////////////////////////////////////////////////////////
762 /// Internal subrogram used by Create.
763 /// It grow new cells by the binary division process.
764 /// This function is overridden by the PDEFoam class to stop the foam buildup process
765 /// if one of the cut conditions stop the cell split.
766 
768 {
769  fTimer->Init(fNCells);
770 
771  Long_t iCell;
772  PDEFoamCell* newCell;
773 
774  while ( (fLastCe+2) < fNCells ) { // this condition also checked inside Divide
775  iCell = PeekMax(); // peek up cell with maximum driver integral
776 
777  if ( (iCell<0) || (iCell>fLastCe) ) {
778  Log() << kVERBOSE << "Break: "<< fLastCe+1 << " cells created" << Endl;
779  // remove remaining empty cells
780  for (Long_t jCell=fLastCe+1; jCell<fNCells; jCell++)
781  delete fCells[jCell];
782  fNCells = fLastCe+1;
783  break;
784  }
785  newCell = fCells[iCell];
786 
787  OutputGrow();
788 
789  if ( Divide( newCell )==0) break; // and divide it into two
790  }
791  OutputGrow( kTRUE );
792  CheckAll(1); // set arg=1 for more info
793 
794  Log() << kVERBOSE << GetNActiveCells() << " active cells created" << Endl;
795 }// Grow
796 
797 ////////////////////////////////////////////////////////////////////////////////
798 /// This can be called before Create, after setting kDim
799 /// It defines which variables are excluded in the process of the cell division.
800 /// For example 'FoamX->SetInhiDiv(1, 1);' inhibits division of y-variable.
801 
803 {
804  if(fDim==0) Log() << kFATAL << "SetInhiDiv: fDim=0" << Endl;
805  if(fInhiDiv == 0) {
806  fInhiDiv = new Int_t[ fDim ];
807  for(Int_t i=0; i<fDim; i++) fInhiDiv[i]=0;
808  }
809  //
810  if( ( 0<=iDim) && (iDim<fDim)) {
811  fInhiDiv[iDim] = inhiDiv;
812  } else
813  Log() << kFATAL << "Wrong iDim" << Endl;
814 }//SetInhiDiv
815 
816 ////////////////////////////////////////////////////////////////////////////////
817 /// User utility, miscellaneous and debug.
818 /// Checks all pointers in the tree of cells. This is useful autodiagnostic.
819 /// level=0, no printout, failures causes STOP
820 /// level=1, printout, failures lead to WARNINGS only
821 
823 {
824  Int_t errors, warnings;
825  PDEFoamCell *cell;
826  Long_t iCell;
827 
828  errors = 0; warnings = 0;
829  if (level==1) Log() << kVERBOSE << "Performing consistency checks for created foam" << Endl;
830  for(iCell=1; iCell<=fLastCe; iCell++) {
831  cell = fCells[iCell];
832  // checking general rules
833  if( ((cell->GetDau0()==0) && (cell->GetDau1()!=0) ) ||
834  ((cell->GetDau1()==0) && (cell->GetDau0()!=0) ) ) {
835  errors++;
836  if (level==1) Log() << kFATAL << "ERROR: Cell's no %d has only one daughter " << iCell << Endl;
837  }
838  if( (cell->GetDau0()==0) && (cell->GetDau1()==0) && (cell->GetStat()==0) ) {
839  errors++;
840  if (level==1) Log() << kFATAL << "ERROR: Cell's no %d has no daughter and is inactive " << iCell << Endl;
841  }
842  if( (cell->GetDau0()!=0) && (cell->GetDau1()!=0) && (cell->GetStat()==1) ) {
843  errors++;
844  if (level==1) Log() << kFATAL << "ERROR: Cell's no %d has two daughters and is active " << iCell << Endl;
845  }
846 
847  // checking parents
848  if( (cell->GetPare())!=fCells[0] ) { // not child of the root
849  if ( (cell != cell->GetPare()->GetDau0()) && (cell != cell->GetPare()->GetDau1()) ) {
850  errors++;
851  if (level==1) Log() << kFATAL << "ERROR: Cell's no %d parent not pointing to this cell " << iCell << Endl;
852  }
853  }
854 
855  // checking daughters
856  if(cell->GetDau0()!=0) {
857  if(cell != (cell->GetDau0())->GetPare()) {
858  errors++;
859  if (level==1) Log() << kFATAL << "ERROR: Cell's no %d daughter 0 not pointing to this cell " << iCell << Endl;
860  }
861  }
862  if(cell->GetDau1()!=0) {
863  if(cell != (cell->GetDau1())->GetPare()) {
864  errors++;
865  if (level==1) Log() << kFATAL << "ERROR: Cell's no %d daughter 1 not pointing to this cell " << iCell << Endl;
866  }
867  }
868  if(cell->GetVolume()<1E-50) {
869  errors++;
870  if(level==1) Log() << kFATAL << "ERROR: Cell no. " << iCell << " has Volume of <1E-50" << Endl;
871  }
872  }// loop after cells;
873 
874  // Check for cells with Volume=0
875  for(iCell=0; iCell<=fLastCe; iCell++) {
876  cell = fCells[iCell];
877  if( (cell->GetStat()==1) && (cell->GetVolume()<1E-11) ) {
878  errors++;
879  if(level==1) Log() << kFATAL << "ERROR: Cell no. " << iCell << " is active but Volume is 0 " << Endl;
880  }
881  }
882  // summary
883  if(level==1){
884  Log() << kVERBOSE << "Check has found " << errors << " errors and " << warnings << " warnings." << Endl;
885  }
886  if(errors>0){
887  Info("CheckAll","Check - found total %d errors \n",errors);
888  }
889 } // Check
890 
891 ////////////////////////////////////////////////////////////////////////////////
892 /// Prints geometry of and elements of 'iCell', as well as relations
893 /// to parent and daughter cells.
894 
896 {
897  if (iCell < 0 || iCell > fLastCe) {
898  Log() << kWARNING << "<PrintCell(iCell=" << iCell
899  << ")>: cell number " << iCell << " out of bounds!"
900  << Endl;
901  return;
902  }
903 
904  PDEFoamVect cellPosi(fDim), cellSize(fDim);
905  fCells[iCell]->GetHcub(cellPosi,cellSize);
906  Int_t kBest = fCells[iCell]->GetBest();
907  Double_t xBest = fCells[iCell]->GetXdiv();
908 
909  Log() << "Cell[" << iCell << "]={ ";
910  Log() << " " << fCells[iCell] << " " << Endl; // extra DEBUG
911  Log() << " Xdiv[abs. coord.]="
912  << VarTransformInvers(kBest,cellPosi[kBest] + xBest*cellSize[kBest])
913  << Endl;
914  Log() << " Abs. coord. = (";
915  for (Int_t idim=0; idim<fDim; idim++) {
916  Log() << "dim[" << idim << "]={"
917  << VarTransformInvers(idim,cellPosi[idim]) << ","
918  << VarTransformInvers(idim,cellPosi[idim] + cellSize[idim])
919  << "}";
920  if (idim < fDim-1)
921  Log() << ", ";
922  }
923  Log() << ")" << Endl;
924  fCells[iCell]->Print("1");
925  // print the cell elements
926  Log() << "Elements: [";
927  TVectorD *vec = (TVectorD*)fCells[iCell]->GetElement();
928  if (vec != NULL){
929  for (Int_t i=0; i<vec->GetNrows(); i++){
930  if (i>0) Log() << ", ";
931  Log() << GetCellElement(fCells[iCell], i);
932  }
933  } else
934  Log() << "not set";
935  Log() << "]" << Endl;
936  Log()<<"}"<<Endl;
937 }
938 
939 ////////////////////////////////////////////////////////////////////////////////
940 /// Prints geometry of ALL cells of the FOAM
941 
943 {
944  for(Long_t iCell=0; iCell<=fLastCe; iCell++)
945  PrintCell(iCell);
946 }
947 
948 ////////////////////////////////////////////////////////////////////////////////
949 /// This function fills a weight 'wt' into the PDEFoam cell, which
950 /// corresponds to the given event 'ev'. Per default cell element 0
951 /// is filled with the weight 'wt', and cell element 1 is filled
952 /// with the squared weight. This function can be overridden by a
953 /// subclass in order to change the values stored in the foam cells.
954 
956 {
957  // find corresponding foam cell
958  std::vector<Float_t> values = ev->GetValues();
959  std::vector<Float_t> tvalues = VarTransform(values);
960  PDEFoamCell *cell = FindCell(tvalues);
961 
962  // 0. Element: Sum of weights 'wt'
963  // 1. Element: Sum of weights 'wt' squared
964  SetCellElement(cell, 0, GetCellElement(cell, 0) + wt);
965  SetCellElement(cell, 1, GetCellElement(cell, 1) + wt*wt);
966 }
967 
968 ////////////////////////////////////////////////////////////////////////////////
969 /// Remove the cell elements from all cells.
970 
972 {
973  if (!fCells) return;
974 
975  Log() << kVERBOSE << "Delete cell elements" << Endl;
976  for (Long_t iCell = 0; iCell < fNCells; ++iCell) {
977  TObject* elements = fCells[iCell]->GetElement();
978  if (elements) {
979  delete elements;
980  fCells[iCell]->SetElement(NULL);
981  }
982  }
983 }
984 
985 ////////////////////////////////////////////////////////////////////////////////
986 /// Returns true, if the value of the given cell is undefined.
987 /// Default value: kFALSE. This function can be overridden by
988 /// sub-classes.
989 
991 {
992  return kFALSE;
993 }
994 
995 ////////////////////////////////////////////////////////////////////////////////
996 /// This function finds the cell, which corresponds to the given
997 /// untransformed event vector 'xvec' and return its value, which is
998 /// given by the parameter 'cv'. If kernel != NULL, then
999 /// PDEFoamKernelBase::Estimate() is called on the transformed event
1000 /// variables.
1001 ///
1002 /// Parameters:
1003 ///
1004 /// - xvec - event vector (untransformed, [fXmin,fXmax])
1005 ///
1006 /// - cv - the cell value to return
1007 ///
1008 /// - kernel - PDEFoam kernel estimator. If NULL is given, than the
1009 /// pure cell value is returned
1010 ///
1011 /// Return:
1012 ///
1013 /// The cell value, corresponding to 'xvec', estimated by the given
1014 /// kernel.
1015 
1016 Float_t TMVA::PDEFoam::GetCellValue(const std::vector<Float_t> &xvec, ECellValue cv, PDEFoamKernelBase *kernel)
1017 {
1018  std::vector<Float_t> txvec(VarTransform(xvec));
1019  if (kernel == NULL)
1020  return GetCellValue(FindCell(txvec), cv);
1021  else
1022  return kernel->Estimate(this, txvec, cv);
1023 }
1024 
1025 ////////////////////////////////////////////////////////////////////////////////
1026 /// This function finds all cells, which corresponds to the given
1027 /// (incomplete) untransformed event vector 'xvec' and returns the
1028 /// cell values, according to the parameter 'cv'.
1029 ///
1030 /// Parameters:
1031 ///
1032 /// - xvec - map for the untransformed vector. The key (Int_t) is
1033 /// the dimension, and the value (Float_t) is the event
1034 /// coordinate. Note that not all coordinates have to be
1035 /// specified.
1036 ///
1037 /// - cv - cell values to return
1038 ///
1039 /// Return:
1040 ///
1041 /// cell values from all cells that were found
1042 
1043 std::vector<Float_t> TMVA::PDEFoam::GetCellValue( const std::map<Int_t,Float_t>& xvec, ECellValue cv )
1044 {
1045  // transformed event
1046  std::map<Int_t,Float_t> txvec;
1047  for (std::map<Int_t,Float_t>::const_iterator it=xvec.begin(); it!=xvec.end(); ++it)
1048  txvec.insert(std::pair<Int_t, Float_t>(it->first, VarTransform(it->first, it->second)));
1049 
1050  // find all cells, which correspond to the transformed event
1051  std::vector<PDEFoamCell*> cells = FindCells(txvec);
1052 
1053  // get the cell values
1054  std::vector<Float_t> cell_values;
1055  cell_values.reserve(cells.size());
1056  for (std::vector<PDEFoamCell*>::const_iterator cell_it=cells.begin();
1057  cell_it != cells.end(); ++cell_it)
1058  cell_values.push_back(GetCellValue(*cell_it, cv));
1059 
1060  return cell_values;
1061 }
1062 
1063 ////////////////////////////////////////////////////////////////////////////////
1064 /// Find cell that contains 'xvec' (in foam coordinates [0,1]).
1065 ///
1066 /// Loop to find cell that contains 'xvec' starting at root cell,
1067 /// and traversing binary tree to find the cell quickly. Note, that
1068 /// if 'xvec' lies outside the foam, the cell which is nearest to
1069 /// 'xvec' is returned. (The returned pointer should never be
1070 /// NULL.)
1071 ///
1072 /// Parameters:
1073 ///
1074 /// - xvec - event vector (in foam coordinates [0,1])
1075 ///
1076 /// Return:
1077 ///
1078 /// PDEFoam cell corresponding to 'xvec'
1079 
1080 TMVA::PDEFoamCell* TMVA::PDEFoam::FindCell( const std::vector<Float_t> &xvec ) const
1081 {
1082  PDEFoamVect cellPosi0(GetTotDim()), cellSize0(GetTotDim());
1083  PDEFoamCell *cell, *cell0;
1084 
1085  cell=fCells[0]; // start with root cell
1086  Int_t idim=0;
1087  while (cell->GetStat()!=1) { //go down binary tree until cell is found
1088  idim=cell->GetBest(); // dimension that changed
1089  cell0=cell->GetDau0();
1090  cell0->GetHcub(cellPosi0,cellSize0);
1091 
1092  if (xvec.at(idim)<=cellPosi0[idim]+cellSize0[idim])
1093  cell=cell0;
1094  else
1095  cell=(cell->GetDau1());
1096  }
1097  return cell;
1098 }
1099 
1100 ////////////////////////////////////////////////////////////////////////////////
1101 /// This is a helper function for std::vector<PDEFoamCell*>
1102 /// FindCells(...) and a generalisation of PDEFoamCell* FindCell().
1103 /// It saves in 'cells' all cells, which contain the coordinates
1104 /// specifies in 'txvec'. Note, that not all coordinates have to be
1105 /// specified in 'txvec'.
1106 ///
1107 /// Parameters:
1108 ///
1109 /// - txvec - event vector in foam coordinates [0,1]. The key is
1110 /// the dimension and the value is the event coordinate. Note,
1111 /// that not all coordinates have to be specified.
1112 ///
1113 /// - cell - cell to start searching with (usually root cell
1114 /// fCells[0])
1115 ///
1116 /// - cells - list of cells that were found
1117 
1118 void TMVA::PDEFoam::FindCells(const std::map<Int_t, Float_t> &txvec, PDEFoamCell* cell, std::vector<PDEFoamCell*> &cells) const
1119 {
1120  PDEFoamVect cellPosi0(GetTotDim()), cellSize0(GetTotDim());
1121  PDEFoamCell *cell0;
1122  Int_t idim=0;
1123 
1124  while (cell->GetStat()!=1) { //go down binary tree until cell is found
1125  idim=cell->GetBest(); // dimension that changed
1126 
1127  // check if dimension 'idim' is specified in 'txvec'
1128  map<Int_t, Float_t>::const_iterator it = txvec.find(idim);
1129 
1130  if (it != txvec.end()){
1131  // case 1: cell is splitten in a dimension which is specified
1132  // in txvec
1133  cell0=cell->GetDau0();
1134  cell0->GetHcub(cellPosi0,cellSize0);
1135  // check, whether left daughter cell contains txvec
1136  if (it->second <= cellPosi0[idim] + cellSize0[idim])
1137  cell=cell0;
1138  else
1139  cell=cell->GetDau1();
1140  } else {
1141  // case 2: cell is splitten in target dimension
1142  FindCells(txvec, cell->GetDau0(), cells);
1143  FindCells(txvec, cell->GetDau1(), cells);
1144  return;
1145  }
1146  }
1147  cells.push_back(cell);
1148 }
1149 
1150 ////////////////////////////////////////////////////////////////////////////////
1151 /// Find all cells, that contain txvec. This function can be used,
1152 /// when the dimension of the foam is greater than the dimension of
1153 /// txvec. E.g. this is the case for multi-target regression.
1154 ///
1155 /// Parameters:
1156 ///
1157 /// - txvec - event vector of variables, transformed into foam
1158 /// coordinates [0,1]. The size of txvec can be smaller than the
1159 /// dimension of the foam.
1160 ///
1161 /// Return value:
1162 ///
1163 /// - vector of cells, that fit txvec
1164 
1165 std::vector<TMVA::PDEFoamCell*> TMVA::PDEFoam::FindCells(const std::vector<Float_t> &txvec) const
1166 {
1167  // copy the coordinates from 'txvec' into a map
1168  std::map<Int_t, Float_t> txvec_map;
1169  for (UInt_t i=0; i<txvec.size(); ++i)
1170  txvec_map.insert(std::pair<Int_t, Float_t>(i, txvec.at(i)));
1171 
1172  // the cells found
1173  std::vector<PDEFoamCell*> cells(0);
1174 
1175  // loop over all target dimensions
1176  FindCells(txvec_map, fCells[0], cells);
1177 
1178  return cells;
1179 }
1180 
1181 ////////////////////////////////////////////////////////////////////////////////
1182 /// Find all cells, that contain the coordinates specified in txvec.
1183 /// The key in 'txvec' is the dimension, and the corresponding value
1184 /// is the coordinate. Note, that not all coordinates have to be
1185 /// specified in txvec.
1186 ///
1187 /// Parameters:
1188 ///
1189 /// - txvec - map of coordinates (transformed into foam coordinates
1190 /// [0,1])
1191 ///
1192 /// Return value:
1193 ///
1194 /// - vector of cells, that fit txvec
1195 
1196 std::vector<TMVA::PDEFoamCell*> TMVA::PDEFoam::FindCells(const std::map<Int_t, Float_t> &txvec) const
1197 {
1198  // the cells found
1199  std::vector<PDEFoamCell*> cells(0);
1200 
1201  // loop over all target dimensions
1202  FindCells(txvec, fCells[0], cells);
1203 
1204  return cells;
1205 }
1206 
1207 ////////////////////////////////////////////////////////////////////////////////
1208 /// Draws 1-dimensional foam (= histogram)
1209 ///
1210 /// Parameters:
1211 ///
1212 /// - cell_value - the cell value to draw
1213 ///
1214 /// - nbin - number of bins of result histogram
1215 ///
1216 /// - kernel - a PDEFoam kernel.
1217 
1219 {
1220  // avoid plotting of wrong dimensions
1221  if ( GetTotDim()!=1 )
1222  Log() << kFATAL << "<Draw1Dim>: function can only be used for 1-dimensional foams!"
1223  << Endl;
1224 
1225  TString hname("h_1dim");
1226  TH1D* h1=(TH1D*)gDirectory->Get(hname);
1227  if (h1) delete h1;
1228  h1= new TH1D(hname, "1-dimensional Foam", nbin, fXmin[0], fXmax[0]);
1229 
1230  if (!h1) Log() << kFATAL << "ERROR: Can not create histo" << hname << Endl;
1231 
1232  // loop over all bins
1233  for (Int_t ibinx=1; ibinx<=h1->GetNbinsX(); ++ibinx) {
1234  // get event vector corresponding to bin
1235  std::vector<Float_t> txvec;
1236  txvec.push_back( VarTransform(0, h1->GetBinCenter(ibinx)) );
1237  Float_t val = 0;
1238  if (kernel != NULL) {
1239  // get cell value using the kernel
1240  val = kernel->Estimate(this, txvec, cell_value);
1241  } else {
1242  val = GetCellValue(FindCell(txvec), cell_value);
1243  }
1244  // fill value to histogram
1245  h1->SetBinContent(ibinx, val + h1->GetBinContent(ibinx));
1246  }
1247 
1248  return h1;
1249 }
1250 
1251 ////////////////////////////////////////////////////////////////////////////////
1252 /// Project foam variable idim1 and variable idim2 to histogram.
1253 ///
1254 /// Parameters:
1255 ///
1256 /// - idim1, idim2 - dimensions to project to
1257 ///
1258 /// - cell_value - the cell value to draw
1259 ///
1260 /// - kernel - a PDEFoam kernel (optional). If NULL is given, the
1261 /// kernel is ignored and the pure cell values are
1262 /// plotted.
1263 ///
1264 /// - nbin - number of bins in x and y direction of result histogram
1265 /// (optional, default is 50).
1266 ///
1267 /// Returns:
1268 /// a 2-dimensional histogram
1269 
1270 TH2D* TMVA::PDEFoam::Project2( Int_t idim1, Int_t idim2, ECellValue cell_value, PDEFoamKernelBase *kernel, UInt_t nbin )
1271 {
1272  // avoid plotting of wrong dimensions
1273  if ((idim1>=GetTotDim()) || (idim1<0) ||
1274  (idim2>=GetTotDim()) || (idim2<0) ||
1275  (idim1==idim2) )
1276  Log() << kFATAL << "<Project2>: wrong dimensions given: "
1277  << idim1 << ", " << idim2 << Endl;
1278 
1279  // root can not handle too many bins in one histogram --> catch this
1280  // Furthermore, to have more than 1000 bins in the histogram doesn't make
1281  // sense.
1282  if (nbin>1000){
1283  Log() << kWARNING << "Warning: number of bins too big: " << nbin
1284  << " Using 1000 bins for each dimension instead." << Endl;
1285  nbin = 1000;
1286  } else if (nbin<1) {
1287  Log() << kWARNING << "Wrong bin number: " << nbin
1288  << "; set nbin=50" << Endl;
1289  nbin = 50;
1290  }
1291 
1292  // create result histogram
1293  TString hname(Form("h_%d_vs_%d",idim1,idim2));
1294 
1295  // if histogram with this name already exists, delete it
1296  TH2D* h1=(TH2D*)gDirectory->Get(hname.Data());
1297  if (h1) delete h1;
1298  h1= new TH2D(hname.Data(), Form("var%d vs var%d",idim1,idim2), nbin, fXmin[idim1], fXmax[idim1], nbin, fXmin[idim2], fXmax[idim2]);
1299 
1300  if (!h1) Log() << kFATAL << "ERROR: Can not create histo" << hname << Endl;
1301 
1302  // ============== start projection algorithm ================
1303  // loop over all histogram bins (2-dim)
1304  for (Int_t xbin = 1; xbin <= h1->GetNbinsX(); ++xbin) {
1305  for (Int_t ybin = 1; ybin <= h1->GetNbinsY(); ++ybin) {
1306  // calculate the phase space point, which corresponds to this
1307  // bin combination
1308  std::map<Int_t, Float_t> txvec;
1309  txvec[idim1] = VarTransform(idim1, h1->GetXaxis()->GetBinCenter(xbin));
1310  txvec[idim2] = VarTransform(idim2, h1->GetYaxis()->GetBinCenter(ybin));
1311 
1312  // find the cells, which corresponds to this phase space
1313  // point
1314  std::vector<TMVA::PDEFoamCell*> cells = FindCells(txvec);
1315 
1316  // loop over cells and fill the histogram with the cell
1317  // values
1318  Float_t sum_cv = 0; // sum of the cell values
1319  for (std::vector<TMVA::PDEFoamCell*>::const_iterator it = cells.begin();
1320  it != cells.end(); ++it) {
1321  // get cell position and size
1322  PDEFoamVect cellPosi(GetTotDim()), cellSize(GetTotDim());
1323  (*it)->GetHcub(cellPosi,cellSize);
1324  // Create complete event vector from txvec. The missing
1325  // coordinates of txvec are set to the cell center.
1326  std::vector<Float_t> tvec;
1327  for (Int_t i=0; i<GetTotDim(); ++i) {
1328  if ( i != idim1 && i != idim2 )
1329  tvec.push_back(cellPosi[i] + 0.5*cellSize[i]);
1330  else
1331  tvec.push_back(txvec[i]);
1332  }
1333  if (kernel != NULL) {
1334  // get the cell value using the kernel
1335  sum_cv += kernel->Estimate(this, tvec, cell_value);
1336  } else {
1337  sum_cv += GetCellValue(FindCell(tvec), cell_value);
1338  }
1339  }
1340 
1341  // fill the bin content
1342  h1->SetBinContent(xbin, ybin, sum_cv + h1->GetBinContent(xbin, ybin));
1343  }
1344  }
1345 
1346  return h1;
1347 }
1348 
1349 ////////////////////////////////////////////////////////////////////////////////
1350 /// Returns the cell value of 'cell' corresponding to the given
1351 /// option 'cv'. This function should be overridden by the subclass
1352 /// in order to specify which cell elements to return for a given
1353 /// cell value 'cv'. By default kValue returns cell element 0, and
1354 /// kValueError returns cell element 1.
1355 
1357 {
1358  // calculate cell value (depending on the given option 'cv')
1359  switch (cv) {
1360 
1361  case kValue:
1362  return GetCellElement(cell, 0);
1363 
1364  case kValueError:
1365  return GetCellElement(cell, 1);
1366 
1367  case kValueDensity: {
1368 
1369  Double_t volume = cell->GetVolume();
1370  if (volume > numeric_limits<double>::epsilon()) {
1371  return GetCellValue(cell, kValue)/volume;
1372  } else {
1373  if (volume<=0){
1374  cell->Print("1"); // debug output
1375  Log() << kWARNING << "<GetCellDensity(cell)>: ERROR: cell volume"
1376  << " negative or zero!"
1377  << " ==> return cell density 0!"
1378  << " cell volume=" << volume
1379  << " cell entries=" << GetCellValue(cell, kValue) << Endl;
1380  } else {
1381  Log() << kWARNING << "<GetCellDensity(cell)>: WARNING: cell volume"
1382  << " close to zero!"
1383  << " cell volume: " << volume << Endl;
1384  }
1385  }
1386  }
1387  return 0;
1388 
1389  case kMeanValue:
1390  return cell->GetIntg();
1391 
1392  case kRms:
1393  return cell->GetDriv();
1394 
1395  case kRmsOvMean:
1396  if (cell->GetIntg() != 0)
1397  return cell->GetDriv()/cell->GetIntg();
1398  else
1399  return 0;
1400 
1401  case kCellVolume:
1402  return cell->GetVolume();
1403 
1404  default:
1405  Log() << kFATAL << "<GetCellValue>: unknown cell value" << Endl;
1406  return 0;
1407  }
1408 
1409  return 0;
1410 }
1411 
1412 ////////////////////////////////////////////////////////////////////////////////
1413 /// Returns cell element i of cell 'cell'. If the cell has no
1414 /// elements or the index 'i' is out of range, than 0 is returned.
1415 
1417 {
1418  // dynamic_cast doesn't seem to work here ?!
1419  TVectorD *vec = (TVectorD*)cell->GetElement();
1420 
1421  // if vec is not set or index out of range, return 0
1422  if (!vec || i >= (UInt_t) vec->GetNrows())
1423  return 0;
1424 
1425  return (*vec)(i);
1426 }
1427 
1428 ////////////////////////////////////////////////////////////////////////////////
1429 /// Set cell element i of cell to value. If the cell element i does
1430 /// not exist, it is created.
1431 
1433 {
1434  TVectorD *vec = NULL;
1435 
1436  // if no cell elements are set, create TVectorD with i+1 entries,
1437  // ranging from [0,i]
1438  if (cell->GetElement() == NULL) {
1439  vec = new TVectorD(i+1);
1440  vec->Zero(); // set all values to zero
1441  (*vec)(i) = value; // set element i to value
1442  cell->SetElement(vec);
1443  } else {
1444  // dynamic_cast doesn't seem to work here ?!
1445  vec = (TVectorD*)cell->GetElement();
1446  if (!vec)
1447  Log() << kFATAL << "<SetCellElement> ERROR: cell element is not a TVectorD*" << Endl;
1448  // check vector size and resize if necessary
1449  if (i >= (UInt_t) vec->GetNrows())
1450  vec->ResizeTo(0,i);
1451  // set element i to value
1452  (*vec)(i) = value;
1453  }
1454 }
1455 
1456 ////////////////////////////////////////////////////////////////////////////////
1457 /// Overridden function of PDEFoam to avoid native foam output.
1458 /// Draw TMVA-process bar instead.
1459 
1461 {
1462  if (finished) {
1463  Log() << kINFO << "Elapsed time: " + fTimer->GetElapsedTime()
1464  << " " << Endl;
1465  return;
1466  }
1467 
1468  Int_t modulo = 1;
1469 
1470  if (fNCells >= 100) modulo = Int_t(fNCells/100);
1471  if (fLastCe%modulo == 0) fTimer->DrawProgressBar( fLastCe );
1472 }
1473 
1474 ////////////////////////////////////////////////////////////////////////////////
1475 /// Debugging tool which plots the cells of a 2-dimensional PDEFoam
1476 /// as rectangles in C++ format readable for ROOT.
1477 ///
1478 /// Parameters:
1479 /// - filename - filename of ouput root macro
1480 ///
1481 /// - opt - cell_value, rms, rms_ov_mean
1482 /// If cell_value is set, the following values will be filled into
1483 /// the result histogram:
1484 /// - number of events - in case of classification with 2 separate
1485 /// foams or multi-target regression
1486 /// - discriminator - in case of classification with one
1487 /// unified foam
1488 /// - target - in case of mono-target regression
1489 /// If none of {cell_value, rms, rms_ov_mean} is given, the cells
1490 /// will not be filled.
1491 /// If 'opt' contains the string 'cellnumber', the index of
1492 /// each cell is draw in addition.
1493 ///
1494 /// - createCanvas - whether to create a new canvas or not
1495 ///
1496 /// - colors - whether to fill cells with colors or shades of grey
1497 ///
1498 /// Example:
1499 ///
1500 /// The following commands load a mono-target regression foam from
1501 /// file 'foam.root' and create a ROOT macro 'output.C', which
1502 /// draws all PDEFoam cells with little boxes. The latter are
1503 /// filled with colors according to the target value stored in the
1504 /// cell. Also the cell number is drawn.
1505 ///
1506 /// TFile file("foam.root");
1507 /// TMVA::PDEFoam *foam = (TMVA::PDEFoam*) gDirectory->Get("MonoTargetRegressionFoam");
1508 /// foam->RootPlot2dim("output.C","cell_value,cellnumber");
1509 /// gROOT->Macro("output.C");
1510 
1511 void TMVA::PDEFoam::RootPlot2dim( const TString& filename, TString opt,
1512  Bool_t createCanvas, Bool_t colors )
1513 {
1514  if (GetTotDim() != 2)
1515  Log() << kFATAL << "RootPlot2dim() can only be used with "
1516  << "two-dimensional foams!" << Endl;
1517 
1518  // select value to plot
1519  ECellValue cell_value = kValue;
1520  Bool_t plotcellnumber = kFALSE;
1521  Bool_t fillcells = kTRUE;
1522  if (opt.Contains("cell_value")){
1523  cell_value = kValue;
1524  } else if (opt.Contains("rms_ov_mean")){
1525  cell_value = kRmsOvMean;
1526  } else if (opt.Contains("rms")){
1527  cell_value = kRms;
1528  } else {
1529  fillcells = kFALSE;
1530  }
1531  if (opt.Contains("cellnumber"))
1532  plotcellnumber = kTRUE;
1533 
1534  // open file (root macro)
1535  std::ofstream outfile(filename, std::ios::out);
1536 
1537  outfile<<"{" << std::endl;
1538 
1539  // declare boxes and set the fill styles
1540  if (!colors) { // define grayscale colors from light to dark,
1541  // starting from color index 1000
1542  outfile << "TColor *graycolors[100];" << std::endl;
1543  outfile << "for (Int_t i=0.; i<100; i++)" << std::endl;
1544  outfile << " graycolors[i]=new TColor(1000+i, 1-(Float_t)i/100.,1-(Float_t)i/100.,1-(Float_t)i/100.);"<< std::endl;
1545  }
1546  if (createCanvas)
1547  outfile << "cMap = new TCanvas(\"" << fName << "\",\"Cell Map for "
1548  << fName << "\",600,600);" << std::endl;
1549 
1550  outfile<<"TBox*a=new TBox();"<<std::endl;
1551  outfile<<"a->SetFillStyle(0);"<<std::endl; // big frame
1552  outfile<<"a->SetLineWidth(4);"<<std::endl;
1553  outfile<<"TBox *b1=new TBox();"<<std::endl; // single cell
1554  outfile<<"TText*t=new TText();"<<std::endl; // text for numbering
1555  if (fillcells) {
1556  outfile << (colors ? "gStyle->SetPalette(1, 0);" : "gStyle->SetPalette(0);")
1557  << std::endl;
1558  outfile <<"b1->SetFillStyle(1001);"<<std::endl;
1559  outfile<<"TBox *b2=new TBox();"<<std::endl; // single cell
1560  outfile <<"b2->SetFillStyle(0);"<<std::endl;
1561  }
1562  else {
1563  outfile <<"b1->SetFillStyle(0);"<<std::endl;
1564  }
1565 
1566  if (fillcells)
1567  (colors ? gStyle->SetPalette(1, 0) : gStyle->SetPalette(0) );
1568 
1569  Float_t zmin = 1E8; // minimal value (for color calculation)
1570  Float_t zmax = -1E8; // maximal value (for color calculation)
1571 
1572  // if cells shall be filled, calculate minimal and maximal plot
1573  // value --> store in zmin and zmax
1574  if (fillcells) {
1575  for (Long_t iCell=1; iCell<=fLastCe; iCell++) {
1576  if ( fCells[iCell]->GetStat() == 1) {
1577  Float_t value = GetCellValue(fCells[iCell], cell_value);
1578  if (value<zmin)
1579  zmin=value;
1580  if (value>zmax)
1581  zmax=value;
1582  }
1583  }
1584  outfile << "// observed minimum and maximum of distribution: " << std::endl;
1585  outfile << "// Float_t zmin = "<< zmin << ";" << std::endl;
1586  outfile << "// Float_t zmax = "<< zmax << ";" << std::endl;
1587  }
1588 
1589  outfile << "// used minimum and maximum of distribution (taking into account log scale if applicable): " << std::endl;
1590  outfile << "Float_t zmin = "<< zmin << ";" << std::endl;
1591  outfile << "Float_t zmax = "<< zmax << ";" << std::endl;
1592 
1593  Float_t x1,y1,x2,y2,x,y; // box and text coordintates
1594  Float_t offs = 0.01;
1595  Float_t lpag = 1-2*offs;
1596  Int_t ncolors = colors ? gStyle->GetNumberOfColors() : 100;
1597  Float_t scale = (ncolors-1)/(zmax - zmin);
1598  PDEFoamVect cellPosi(GetTotDim()), cellSize(GetTotDim());
1599 
1600  // loop over cells and draw a box for every cell (and maybe the
1601  // cell number as well)
1602  outfile << "// =========== Rectangular cells ==========="<< std::endl;
1603  for (Long_t iCell=1; iCell<=fLastCe; iCell++) {
1604  if ( fCells[iCell]->GetStat() == 1) {
1605  fCells[iCell]->GetHcub(cellPosi,cellSize);
1606  x1 = offs+lpag*(cellPosi[0]);
1607  y1 = offs+lpag*(cellPosi[1]);
1608  x2 = offs+lpag*(cellPosi[0]+cellSize[0]);
1609  y2 = offs+lpag*(cellPosi[1]+cellSize[1]);
1610 
1611  if (fillcells) {
1612  // get cell value
1613  Float_t value = GetCellValue(fCells[iCell], cell_value);
1614 
1615  // calculate fill color
1616  Int_t color;
1617  if (colors)
1618  color = gStyle->GetColorPalette(Int_t((value-zmin)*scale));
1619  else
1620  color = 1000+(Int_t((value-zmin)*scale));
1621 
1622  // set fill color of box b1
1623  outfile << "b1->SetFillColor(" << color << ");" << std::endl;
1624  }
1625 
1626  // cell rectangle
1627  outfile<<"b1->DrawBox("<<x1<<","<<y1<<","<<x2<<","<<y2<<");"<<std::endl;
1628  if (fillcells)
1629  outfile<<"b2->DrawBox("<<x1<<","<<y1<<","<<x2<<","<<y2<<");"<<std::endl;
1630 
1631  // cell number
1632  if (plotcellnumber) {
1633  outfile<<"t->SetTextColor(4);"<<std::endl;
1634  if(fLastCe<51)
1635  outfile<<"t->SetTextSize(0.025);"<<std::endl; // text for numbering
1636  else if(fLastCe<251)
1637  outfile<<"t->SetTextSize(0.015);"<<std::endl;
1638  else
1639  outfile<<"t->SetTextSize(0.008);"<<std::endl;
1640  x = offs+lpag*(cellPosi[0]+0.5*cellSize[0]);
1641  y = offs+lpag*(cellPosi[1]+0.5*cellSize[1]);
1642  outfile<<"t->DrawText("<<x<<","<<y<<","<<"\""<<iCell<<"\""<<");"<<std::endl;
1643  }
1644  }
1645  }
1646  outfile<<"// ============== End Rectangles ==========="<< std::endl;
1647 
1648  outfile << "}" << std::endl;
1649  outfile.flush();
1650  outfile.close();
1651 }
1652 
1653 ////////////////////////////////////////////////////////////////////////////////
1654 /// Insert event to internal foam's density estimator
1655 /// PDEFoamDensityBase.
1656 
1658 {
1660 }
1661 
1662 ////////////////////////////////////////////////////////////////////////////////
1663 /// Delete the foam's density estimator, which contains the binary
1664 /// search tree.
1665 
1667 {
1668  if(fDistr) delete fDistr;
1669  fDistr = NULL;
1670 }
Double_t * fRvec
Definition: PDEFoam.h:121
Bool_t fFillFoamWithOrigWeights
Definition: PDEFoam.h:134
PDEFoamCell * GetDau1() const
Definition: PDEFoamCell.h:101
void Varedu(Double_t [], Int_t &, Double_t &, Double_t &)
Internal subrogram used by Create.
Definition: PDEFoam.cxx:570
virtual void FillFoamCells(const Event *ev, Float_t wt)
This function fills a weight &#39;wt&#39; into the PDEFoam cell, which corresponds to the given event &#39;ev&#39;...
Definition: PDEFoam.cxx:955
An array of TObjects.
Definition: TObjArray.h:39
TVectorT< Element > & ResizeTo(Int_t lwb, Int_t upb)
Resize the vector to [lwb:upb] .
Definition: TVectorT.cxx:292
MsgLogger * fLogger
Definition: PDEFoam.h:140
UInt_t fNmin
Definition: PDEFoam.h:131
Random number generator class based on M.
Definition: TRandom3.h:29
virtual Double_t GetBinCenter(Int_t bin) const
Return bin center for 1D histogram.
Definition: TH1.cxx:8251
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition: TObject.cxx:899
MsgLogger & Endl(MsgLogger &ml)
Definition: MsgLogger.h:162
void SetElement(TObject *fobj)
Definition: PDEFoamCell.h:112
void Print(Option_t *option) const
Printout of the cell geometry parameters for the debug purpose.
Long_t PeekMax()
Internal subprogram used by Create.
Definition: PDEFoam.cxx:651
TObjArray * fVariableNames
timer for graphical output
Definition: PDEFoam.h:139
void Init(Int_t ncounts)
Definition: Timer.cxx:104
Int_t * fInhiDiv
[fDim] Dynamic Mask for cell division
Definition: PDEFoam.h:114
Float_t fVolFrac
Definition: PDEFoam.h:133
Int_t fDim
Definition: PDEFoam.h:106
float Float_t
Definition: RtypesCore.h:53
EFoamType fFoamType
Definition: PDEFoam.h:127
void OutputGrow(Bool_t finished=false)
message logger
Definition: PDEFoam.cxx:1460
Float_t VarTransform(Int_t idim, Float_t x) const
Definition: PDEFoam.h:303
R__EXTERN TStyle * gStyle
Definition: TStyle.h:418
Double_t GetVolume() const
Definition: PDEFoamCell.h:91
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
void Fill(Int_t, PDEFoamCell *, PDEFoamCell *, PDEFoamCell *)
Fills in certain data into newly allocated cell.
std::vector< TMVA::PDEFoamCell * > FindCells(const std::vector< Float_t > &) const
Find all cells, that contain txvec.
Definition: PDEFoam.cxx:1165
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition: TH1.cxx:4638
void Grow()
Internal subrogram used by Create.
Definition: PDEFoam.cxx:767
void DrawProgressBar(Int_t, const TString &comment="")
draws progress bar in color or B&W caution:
Definition: Timer.cxx:186
Double_t * fXmax
Definition: PDEFoam.h:129
static Bool_t AddDirectoryStatus()
Static function: cannot be inlined on Windows/NT.
Definition: TH1.cxx:700
void SetXmin(Int_t idim, Double_t wmin)
set lower foam bound in dimension idim
Definition: PDEFoam.cxx:265
Int_t GetNrows() const
Definition: TVectorT.h:81
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
const Bool_t kFALSE
Definition: Rtypes.h:92
Double_t GetCellElement(const PDEFoamCell *cell, UInt_t i) const
Returns cell element i of cell &#39;cell&#39;.
Definition: PDEFoam.cxx:1416
void SetCellElement(PDEFoamCell *cell, UInt_t i, Double_t value)
Set cell element i of cell to value.
Definition: PDEFoam.cxx:1432
Int_t Divide(PDEFoamCell *)
Internal subrogram used by Create.
Definition: PDEFoam.cxx:714
STL namespace.
void SetBest(Int_t Best)
Definition: PDEFoamCell.h:85
TString GetElapsedTime(Bool_t Scientific=kTRUE)
Definition: Timer.cxx:129
TH1D * Draw1Dim(ECellValue cell_value, Int_t nbin, PDEFoamKernelBase *kernel=NULL)
Draws 1-dimensional foam (= histogram)
Definition: PDEFoam.cxx:1218
PDEFoamCell * GetPare() const
Definition: PDEFoamCell.h:99
void CheckAll(Int_t)
User utility, miscellaneous and debug.
Definition: PDEFoam.cxx:822
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
static void AddDirectory(Bool_t add=kTRUE)
Sets the flag controlling the automatic add of histograms in memory.
Definition: TH1.cxx:1218
Int_t GetTotDim() const
Definition: PDEFoam.h:219
Int_t CellFill(Int_t, PDEFoamCell *)
Internal subprogram used by Create.
Definition: PDEFoam.cxx:389
void MakeAlpha()
Internal subrogram used by Create.
Definition: PDEFoam.cxx:637
Timer * fTimer
distribution of training events
Definition: PDEFoam.h:138
double sqrt(double)
virtual Bool_t CellValueIsUndefined(PDEFoamCell *)
Returns true, if the value of the given cell is undefined.
Definition: PDEFoam.cxx:990
static const double x2[5]
virtual Float_t Estimate(PDEFoam *, std::vector< Float_t > &, ECellValue)=0
Double_t x[n]
Definition: legend1.C:17
virtual Double_t Density(std::vector< Double_t > &Xarg, Double_t &event_density)=0
void RootPlot2dim(const TString &filename, TString opt, Bool_t createCanvas=kTRUE, Bool_t colors=kTRUE)
Debugging tool which plots the cells of a 2-dimensional PDEFoam as rectangles in C++ format readable ...
Definition: PDEFoam.cxx:1511
UInt_t GetNmin()
Definition: PDEFoam.h:228
PDEFoamCell * FindCell(const std::vector< Float_t > &) const
Find cell that contains &#39;xvec&#39; (in foam coordinates [0,1]).
Definition: PDEFoam.cxx:1080
PDEFoamDensityBase * fDistr
Definition: PDEFoam.h:137
virtual Double_t GetBinCenter(Int_t bin) const
Return center of bin.
Definition: TAxis.cxx:464
void SetDau1(PDEFoamCell *Daug)
Definition: PDEFoamCell.h:103
TH1F * h1
Definition: legend1.C:5
Int_t GetColorPalette(Int_t i) const
Return color number i in current palette.
Definition: TStyle.cxx:735
TRandom3 * fPseRan
Definition: PDEFoam.h:123
void ResetCellElements()
Remove the cell elements from all cells.
Definition: PDEFoam.cxx:971
virtual void Explore(PDEFoamCell *Cell)
Internal subprogram used by Create.
Definition: PDEFoam.cxx:434
static const Float_t kHigh
Definition: PDEFoam.cxx:95
TVectorT< Double_t > TVectorD
Definition: TVectorDfwd.h:24
void FillBinarySearchTree(const Event *ev)
This method inserts the given event &#39;ev&#39; it into the binary search tree.
Int_t fNBin
Definition: PDEFoam.h:109
void InitCells()
Internal subprogram used by Create.
Definition: PDEFoam.cxx:356
Int_t GetBest() const
Definition: PDEFoamCell.h:84
void PrintCells()
Prints geometry of ALL cells of the FOAM.
Definition: PDEFoam.cxx:942
Float_t VarTransformInvers(Int_t idim, Float_t x) const
Definition: PDEFoam.h:320
void SetXdiv(Double_t Xdiv)
Definition: PDEFoamCell.h:86
void SetXmax(Int_t idim, Double_t wmax)
set upper foam bound in dimension idim
Definition: PDEFoam.cxx:276
UInt_t fNElements
Definition: PDEFoam.h:130
void CalcVolume()
Calculates volume of the cell using size params which are calculated.
Int_t fEvPerBin
Definition: PDEFoam.h:111
PDEFoamDensityBase * GetDistr() const
Definition: PDEFoam.h:180
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content see convention for numbering bins in TH1::GetBin In case the bin number is greater th...
Definition: TH1.cxx:8323
void SetStat(Int_t Stat)
Definition: PDEFoamCell.h:98
TVectorT< Element > & Zero()
Set vector elements to zero.
Definition: TVectorT.cxx:451
unsigned int UInt_t
Definition: RtypesCore.h:42
Bool_t fPeekMax
Definition: PDEFoam.h:136
char * Form(const char *fmt,...)
Double_t E()
Definition: TMath.h:54
UInt_t GetDepth()
Get depth of cell in binary tree, where the root cell has depth 1.
void SetDriv(Double_t Driv)
Definition: PDEFoamCell.h:95
TObject * GetElement() const
Definition: PDEFoamCell.h:113
TAxis * GetYaxis()
Definition: TH1.h:325
void SetDim(Int_t kDim)
Sets dimension of cubical space.
Definition: PDEFoam.cxx:250
tomato 1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:618
TObjArray * fHistEdg
Definition: PDEFoam.h:120
void PrintCell(Long_t iCell=0)
Prints geometry of and elements of &#39;iCell&#39;, as well as relations to parent and daughter cells...
Definition: PDEFoam.cxx:895
MsgLogger & Log() const
Definition: PDEFoam.h:262
Int_t GetNumberOfColors() const
Return number of colors in the color palette.
Definition: TStyle.cxx:801
EDTSeparation fDTSeparation
Definition: PDEFoam.h:135
Double_t * fXmin
Definition: PDEFoam.h:128
REAL epsilon
Definition: triangle.c:617
Double_t GetIntg() const
Definition: PDEFoamCell.h:92
virtual void RndmArray(Int_t n, Float_t *array)
Return an array of n random numbers uniformly distributed in ]0,1].
Definition: TRandom3.cxx:138
Int_t fNSampl
Definition: PDEFoam.h:110
long Long_t
Definition: RtypesCore.h:50
Color * colors
Definition: X3DBuffer.c:19
Double_t Eval(Double_t *xRand, Double_t &event_density)
Internal subprogram.
Definition: PDEFoam.cxx:747
static const double x1[5]
#define ClassImp(name)
Definition: Rtypes.h:279
double Double_t
Definition: RtypesCore.h:55
Double_t GetXdiv() const
Definition: PDEFoamCell.h:83
ECellValue
Definition: PDEFoam.h:85
Int_t fNoAct
[fDim] Flags for inhibiting cell division
Definition: PDEFoam.h:116
void Create()
Basic initialization of FOAM invoked by the user.
Definition: PDEFoam.cxx:292
PDEFoamCell * GetDau0() const
Definition: PDEFoamCell.h:100
Double_t y[n]
Definition: legend1.C:17
Double_t * fAlpha
Definition: PDEFoam.h:125
PDEFoamCell ** fCells
Definition: PDEFoam.h:118
virtual TH2D * Project2(Int_t idim1, Int_t idim2, ECellValue cell_value=kValue, PDEFoamKernelBase *kernel=NULL, UInt_t nbin=50)
Project foam variable idim1 and variable idim2 to histogram.
Definition: PDEFoam.cxx:1270
Int_t * fMaskDiv
Definition: PDEFoam.h:113
THist< 2, double, THistStatContent, THistStatUncertainty > TH2D
Definition: THist.hxx:307
virtual Float_t GetCellValue(const std::vector< Float_t > &xvec, ECellValue cv, PDEFoamKernelBase *)
This function finds the cell, which corresponds to the given untransformed event vector &#39;xvec&#39; and re...
Definition: PDEFoam.cxx:1016
void SetSerial(Int_t Serial)
Definition: PDEFoamCell.h:105
UInt_t fMaxDepth
Definition: PDEFoam.h:132
Mother of all ROOT objects.
Definition: TObject.h:37
void GetHcub(PDEFoamVect &, PDEFoamVect &) const
Provides size and position of the cell These parameter are calculated by analyzing information in all...
void SetIntg(Double_t Intg)
Definition: PDEFoamCell.h:94
Abstract ClassifierFactory template that handles arbitrary types.
void DeleteBinarySearchTree()
Delete the foam&#39;s density estimator, which contains the binary search tree.
Definition: PDEFoam.cxx:1666
std::vector< Float_t > & GetValues()
Definition: Event.h:96
void FillBinarySearchTree(const Event *ev)
Insert event to internal foam&#39;s density estimator PDEFoamDensityBase.
Definition: PDEFoam.cxx:1657
T Sqr(T x) const
Definition: PDEFoam.h:183
Int_t GetStat() const
Definition: PDEFoamCell.h:97
Double_t GetDriv() const
Definition: PDEFoamCell.h:93
virtual ~PDEFoam()
Default destructor.
Definition: PDEFoam.cxx:184
#define NULL
Definition: Rtypes.h:82
Int_t fLastCe
Definition: PDEFoam.h:117
THist< 1, double, THistStatContent, THistStatUncertainty > TH1D
Definition: THist.hxx:301
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition: TH2.h:90
Int_t fNCells
Definition: PDEFoam.h:107
#define gDirectory
Definition: TDirectory.h:221
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content.
Definition: TH2.cxx:2475
void SetDau0(PDEFoamCell *Daug)
Definition: PDEFoamCell.h:102
virtual Int_t GetNbinsX() const
Definition: TH1.h:301
static const Float_t kVlow
Definition: PDEFoam.cxx:96
const Bool_t kTRUE
Definition: Rtypes.h:91
TString fName
Definition: PDEFoam.h:105
void SetPalette(Int_t ncolors=kBird, Int_t *colors=0, Float_t alpha=1.)
See TColor::SetPalette.
Definition: TStyle.cxx:1450
UInt_t GetNActiveCells() const
Definition: PDEFoam.h:221
char name[80]
Definition: TGX11.cxx:109
TAxis * GetXaxis()
Definition: TH1.h:324
void SetInhiDiv(Int_t, Int_t)
This can be called before Create, after setting kDim It defines which variables are excluded in the p...
Definition: PDEFoam.cxx:802
virtual Int_t GetNbinsY() const
Definition: TH1.h:302
PDEFoam()
Default constructor for streamer, user should not use it.
Definition: PDEFoam.cxx:103
tomato 2-D histogram with a double per channel (see TH1 documentation)}
Definition: TH2.h:296
UInt_t GetMaxDepth() const
Definition: PDEFoam.h:230