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