meerkat is hosted by Hepforge, IPPP Durham
Meerkat  v1r3
Multidimensional kernel density estimation package
AdaptiveKernelDensity.hh
Go to the documentation of this file.
1 #ifndef ADAPTIVE_KERNEL_DENSITY
2 #define ADAPTIVE_KERNEL_DENSITY
3 
4 #include "AbsDensity.hh"
5 
6 #include "TMath.h"
7 
8 #include <vector>
9 
15 
16  public:
17 
19 
32  AdaptiveKernelDensity(const char* pdfName,
33  AbsPhaseSpace* thePhaseSpace,
34  TTree* tree,
35  const char *vars1,
36  UInt_t bins1,
37  Double_t width1,
38  AbsDensity* widthScale,
39  AbsDensity* approx = 0,
40  UInt_t toyEvents = 0,
41  UInt_t maxEvents = 0,
42  UInt_t skipEvents = 0
43  );
44 
46 
60  AdaptiveKernelDensity(const char* pdfName,
61  AbsPhaseSpace* thePhaseSpace,
62  TTree* tree,
63  const char *vars1,
64  const char* weight,
65  UInt_t bins1,
66  Double_t width1,
67  AbsDensity* widthScale,
68  AbsDensity* approx = 0,
69  UInt_t toyEvents = 0,
70  UInt_t maxEvents = 0,
71  UInt_t skipEvents = 0
72  );
73 
75 
91  AdaptiveKernelDensity(const char* pdfName,
92  AbsPhaseSpace* thePhaseSpace,
93  TTree* tree,
94  const char *vars1, const char *vars2,
95  UInt_t bins1, UInt_t bins2,
96  Double_t width1, Double_t width2,
97  AbsDensity* widthScale,
98  AbsDensity* approx = 0,
99  UInt_t toyEvents = 0,
100  UInt_t maxEvents = 0,
101  UInt_t skipEvents = 0
102  );
103 
105 
122  AdaptiveKernelDensity(const char* pdfName,
123  AbsPhaseSpace* thePhaseSpace,
124  TTree* tree,
125  const char *vars1, const char *vars2,
126  const char *weight,
127  UInt_t bins1, UInt_t bins2,
128  Double_t width1, Double_t width2,
129  AbsDensity* widthScale,
130  AbsDensity* approx = 0,
131  UInt_t toyEvents = 0,
132  UInt_t maxEvents = 0,
133  UInt_t skipEvents = 0
134  );
135 
137 
156  AdaptiveKernelDensity(const char* pdfName,
157  AbsPhaseSpace* thePhaseSpace,
158  TTree* tree,
159  const char *vars1, const char *vars2, const char *vars3,
160  UInt_t bins1, UInt_t bins2, UInt_t bins3,
161  Double_t width1, Double_t width2, Double_t width3,
162  AbsDensity* widthScale,
163  AbsDensity* approx = 0,
164  UInt_t toyEvents = 0,
165  UInt_t maxEvents = 0,
166  UInt_t skipEvents = 0
167  );
168 
170 
190  AdaptiveKernelDensity(const char* pdfName,
191  AbsPhaseSpace* thePhaseSpace,
192  TTree* tree,
193  const char *vars1, const char *vars2, const char *vars3,
194  const char *weight,
195  UInt_t bins1, UInt_t bins2, UInt_t bins3,
196  Double_t width1, Double_t width2, Double_t width3,
197  AbsDensity* widthScale,
198  AbsDensity* approx = 0,
199  UInt_t toyEvents = 0,
200  UInt_t maxEvents = 0,
201  UInt_t skipEvents = 0
202  );
203 
205 
227  AdaptiveKernelDensity(const char* pdfName,
228  AbsPhaseSpace* thePhaseSpace,
229  TTree* tree,
230  const char *vars1, const char *vars2, const char *vars3, const char *vars4,
231  UInt_t bins1, UInt_t bins2, UInt_t bins3, UInt_t bins4,
232  Double_t width1, Double_t width2, Double_t width3, Double_t width4,
233  AbsDensity* widthScale,
234  AbsDensity* approx = 0,
235  UInt_t toyEvents = 0,
236  UInt_t maxEvents = 0,
237  UInt_t skipEvents = 0
238  );
239 
241 
264  AdaptiveKernelDensity(const char* pdfName,
265  AbsPhaseSpace* thePhaseSpace,
266  TTree* tree,
267  const char *vars1, const char *vars2, const char *vars3, const char *vars4,
268  const char *weight,
269  UInt_t bins1, UInt_t bins2, UInt_t bins3, UInt_t bins4,
270  Double_t width1, Double_t width2, Double_t width3, Double_t width4,
271  AbsDensity* widthScale,
272  AbsDensity* approx = 0,
273  UInt_t toyEvents = 0,
274  UInt_t maxEvents = 0,
275  UInt_t skipEvents = 0
276  );
277 
279 
304  AdaptiveKernelDensity(const char* pdfName,
305  AbsPhaseSpace* thePhaseSpace,
306  TTree* tree,
307  const char *vars1, const char *vars2, const char *vars3, const char *vars4, const char *vars5,
308  UInt_t bins1, UInt_t bins2, UInt_t bins3, UInt_t bins4, UInt_t bins5,
309  Double_t width1, Double_t width2, Double_t width3, Double_t width4, Double_t width5,
310  AbsDensity* widthScale,
311  AbsDensity* approx = 0,
312  UInt_t toyEvents = 0,
313  UInt_t maxEvents = 0,
314  UInt_t skipEvents = 0
315  );
316 
318 
344  AdaptiveKernelDensity(const char* pdfName,
345  AbsPhaseSpace* thePhaseSpace,
346  TTree* tree,
347  const char *vars1, const char *vars2, const char *vars3, const char *vars4, const char *vars5,
348  const char *weight,
349  UInt_t bins1, UInt_t bins2, UInt_t bins3, UInt_t bins4, UInt_t bins5,
350  Double_t width1, Double_t width2, Double_t width3, Double_t width4, Double_t width5,
351  AbsDensity* widthScale,
352  AbsDensity* approx = 0,
353  UInt_t toyEvents = 0,
354  UInt_t maxEvents = 0,
355  UInt_t skipEvents = 0
356  );
357 
359 
372  AdaptiveKernelDensity(const char* pdfName,
373  AbsPhaseSpace* thePhaseSpace,
374  TTree* tree,
375  std::vector<TString> &vars,
376  std::vector<UInt_t> &binning,
377  std::vector<Double_t> &width,
378  AbsDensity* widthScale,
379  AbsDensity* approx = 0,
380  UInt_t toyEvents = 0,
381  UInt_t maxEvents = 0,
382  UInt_t skipEvents = 0
383  );
384 
386  virtual ~AdaptiveKernelDensity();
387 
390 
393  void writeToFile(const char* fileName);
394 
396 
399  void writeToTextFile(const char* fileName);
400 
402 
405  void writeToRootFile(const char* fileName);
406 
408 
412  Double_t density(std::vector<Double_t> &x);
413 
415 
419 
421 
424  void setFractionalMode(Bool_t mode = true) { m_fractionalMode = mode; }
425 
427  void normalise(void);
428 
429  private:
430 
432 
444  void init(AbsPhaseSpace* thePhaseSpace,
445  TTree* tree,
446  std::vector<TString> &vars,
447  std::vector<UInt_t> &binning,
448  std::vector<Double_t> &width,
449  AbsDensity* widthScale,
450  AbsDensity* approx = 0,
451  UInt_t toyEvents = 0,
452  UInt_t maxEvents = 0,
453  UInt_t skipEvents = 0
454  );
455 
457 
461  UInt_t iterToIndex( std::vector<UInt_t> &iter );
462 
464 
470  void addToMap(std::vector<Double_t> &map, std::vector<Double_t> &point,
471  Double_t widthScale = 1., Double_t weight = 1.);
472 
473  void fillMapFromTree( TTree* tree, std::vector<TString> &vars,
474  UInt_t maxEvents = 0, UInt_t skipEvents = 0);
475 
476  void fillMapFromDensity(AbsDensity* theDensity, UInt_t toyEvents = 0);
477 
479 
483  Double_t mapDensity(std::vector<Double_t> &map, std::vector<Double_t> &x);
484 
486  std::vector<Double_t> m_map;
487 
489  std::vector<Double_t> m_approxMap;
490 
492  std::vector<UInt_t> m_binning;
493 
496 
499 
502 
504  std::vector<Double_t> m_width;
505 
507  UInt_t m_dim;
508 
511 
513  Double_t m_minValue;
514 
516  Double_t m_maxValue;
517 
519  Double_t m_minScale;
520 
522  Double_t m_maxScale;
523 
524 };
525 
526 #endif
void normalise(void)
Normalise the PDF such that the average PDF value over the allowed phase space equals to 1...
Abstract class which defines probability density interface.
Definition: AbsDensity.hh:16
void addToMap(std::vector< Double_t > &map, std::vector< Double_t > &point, Double_t widthScale=1., Double_t weight=1.)
Add a kernel density of a single data point to the binned map.
AbsDensity * m_widthScaleDensity
Reference to PDF used for kernel width scaling.
AbsPhaseSpace * phaseSpace()
Return the phase space.
std::vector< Double_t > m_width
Vector of kernel widths.
void fillMapFromTree(TTree *tree, std::vector< TString > &vars, UInt_t maxEvents=0, UInt_t skipEvents=0)
std::vector< Double_t > m_approxMap
Bin map of approximation PDF convolved with the kernel.
virtual ~AdaptiveKernelDensity()
Destructor.
Double_t density(std::vector< Double_t > &x)
Calculate PDF density at the point. Multilinear binned interpolation is used which is calculated in t...
Bool_t m_fractionalMode
Fractional mode flag.
Double_t mapDensity(std::vector< Double_t > &map, std::vector< Double_t > &x)
Calculate the raw density using the binned map (estimated or approximation) at a given point...
Abstract class which defines phase space interface.
std::vector< Double_t > m_map
Bin map of estimated PDF.
AbsPhaseSpace * m_phaseSpace
Reference to phase space.
Double_t m_maxScale
Maximum width scale factor.
UInt_t iterToIndex(std::vector< UInt_t > &iter)
Convert an N-dimensional iterator vector into a linear bin index in the bin map.
AbsDensity * m_approxDensity
Reference to approximation PDF.
void writeToTextFile(const char *fileName)
Write the binned PDF to a text file.
void init(AbsPhaseSpace *thePhaseSpace, TTree *tree, std::vector< TString > &vars, std::vector< UInt_t > &binning, std::vector< Double_t > &width, AbsDensity *widthScale, AbsDensity *approx=0, UInt_t toyEvents=0, UInt_t maxEvents=0, UInt_t skipEvents=0)
Common initialise method used by all constructors.
Double_t m_maxValue
Maximum value of the width scale PDF to be used for scaling.
void writeToRootFile(const char *fileName)
Write the binned density into a ROOT file.
Double_t m_minScale
Minimum width scale factor.
Double_t m_minValue
Minimum value of the width scale PDF to be used for scaling.
void fillMapFromDensity(AbsDensity *theDensity, UInt_t toyEvents=0)
void writeToFile(const char *fileName)
Write density map to file depending on extension.
std::vector< UInt_t > m_binning
Vector of bin numbers in each variable.
UInt_t m_dim
Cached dimensionality of the phase space.
AdaptiveKernelDensity(const char *pdfName, AbsPhaseSpace *thePhaseSpace, TTree *tree, const char *vars1, UInt_t bins1, Double_t width1, AbsDensity *widthScale, AbsDensity *approx=0, UInt_t toyEvents=0, UInt_t maxEvents=0, UInt_t skipEvents=0)
Constructor for 1-dimensional adaptive kernel PDF from the sample of points in an NTuple...
void setFractionalMode(Bool_t mode=true)
Set fractional mode of operation.