meerkat is hosted by Hepforge, IPPP Durham
Meerkat  v1r3
Multidimensional kernel density estimation package
KernelDensity.hh
Go to the documentation of this file.
1 #ifndef KERNEL_DENSITY
2 #define KERNEL_DENSITY
3 
4 #include "AbsDensity.hh"
5 #include "AbsPhaseSpace.hh"
6 
7 #include "TMath.h"
8 #include "TTree.h"
9 #include "TString.h"
10 #include "TRandom3.h"
11 
12 #include <vector>
13 
14 typedef std::vector<Double_t> TPhspVector;
15 typedef std::vector<TPhspVector> TCell;
16 
18 
19 class KernelDensity : public AbsDensity {
20 
21  public:
22 
23  KernelDensity(const char* pdfName,
24  AbsPhaseSpace* thePhaseSpace,
25  std::vector<Double_t> &width,
26  UInt_t approxSize,
27  AbsDensity* approxDensity = 0);
28 
29  KernelDensity(const char* pdfName,
30  AbsPhaseSpace* thePhaseSpace,
31  UInt_t approxSize,
32  Double_t width1,
33  Double_t width2 = 0,
34  Double_t width3 = 0,
35  Double_t width4 = 0,
36  Double_t width5 = 0);
37 
38  KernelDensity(const char* pdfName,
39  AbsPhaseSpace* thePhaseSpace,
40  UInt_t approxSize,
41  AbsDensity* approxDensity,
42  Double_t width1,
43  Double_t width2 = 0,
44  Double_t width3 = 0,
45  Double_t width4 = 0,
46  Double_t width5 = 0);
47 
48  virtual ~KernelDensity();
49 
50  void setWidth(std::vector<Double_t> &width);
51 
52  Bool_t generateApproximation(UInt_t approxSize);
53 
54  Bool_t readTuple(TTree* tree, std::vector<TString> &vars, UInt_t maxEvents = 0);
55 
56  Bool_t readTuple(TTree* tree, const char* var1, UInt_t maxEvents = 0);
57 
58  Bool_t readTuple(TTree* tree, const char* var1, const char* var2, UInt_t maxEvents = 0);
59 
60  Bool_t readTuple(TTree* tree, const char* var1, const char* var2,
61  const char* var3, UInt_t maxEvents = 0);
62 
63  Bool_t readTuple(TTree* tree, const char* var1, const char* var2,
64  const char* var3, const char* var4, UInt_t maxEvents = 0);
65 
66  Bool_t readTuple(TTree* tree, const char* var1, const char* var2,
67  const char* var3, const char* var4,
68  const char* var5, UInt_t maxEvents = 0);
69 
70  Bool_t readTuple(TTree* tree, const char* var1, const char* var2,
71  const char* var3, const char* var4,
72  const char* var5, const char* var6, UInt_t maxEvents = 0);
73 
74  Double_t density(std::vector<Double_t> &x);
75 
77 
78  private:
79 
80  UInt_t numCells(void);
81 
82  Int_t cellIndex(std::vector<Double_t> &x);
83 
84  Double_t rawDensity(std::vector<Double_t> &x, std::vector<TCell> &vector);
85 
87 
89 
91 
92  std::vector<TCell> m_apprVector;
93 
94  std::vector<TCell> m_dataVector;
95 
96 };
97 
98 #endif
TPhspVector m_width
AbsPhaseSpace * m_phaseSpace
std::vector< TCell > m_dataVector
AbsDensity * m_approxDensity
KernelDensity(const char *pdfName, AbsPhaseSpace *thePhaseSpace, std::vector< Double_t > &width, UInt_t approxSize, AbsDensity *approxDensity=0)
Abstract class which defines probability density interface.
Definition: AbsDensity.hh:16
std::vector< TPhspVector > TCell
virtual ~KernelDensity()
AbsPhaseSpace * phaseSpace()
Return phase space definition for this PDF.
Bool_t generateApproximation(UInt_t approxSize)
Create normalisation vector.
void setWidth(std::vector< Double_t > &width)
Set kernel width.
Double_t density(std::vector< Double_t > &x)
Calculate PDF value at the given point.
Abstract class which defines phase space interface.
Double_t rawDensity(std::vector< Double_t > &x, std::vector< TCell > &vector)
Class that describes the unbinned kernel density.
Bool_t readTuple(TTree *tree, std::vector< TString > &vars, UInt_t maxEvents=0)
UInt_t numCells(void)
std::vector< TCell > m_apprVector
Int_t cellIndex(std::vector< Double_t > &x)
std::vector< Double_t > TPhspVector