Jpp  19.1.0
the software that should make you happy
JRootTestkit.hh
Go to the documentation of this file.
1 #ifndef __JROOT__JROOTTESTKIT__
2 #define __JROOT__JROOTTESTKIT__
3 
4 #include <cstddef>
5 #include <vector>
6 #include <algorithm>
7 
8 #pragma GCC diagnostic push
9 #pragma GCC diagnostic ignored "-Wall"
10 #include "TH1.h"
11 #include "TH2.h"
12 #include "TH3.h"
13 #include "TRandom3.h"
14 #pragma GCC diagnostic pop
15 
16 
17 /**
18  * \author mdejong
19  */
20 
21 namespace JROOT {}
22 namespace JPP { using namespace JROOT; }
23 
24 namespace JROOT {
25 
26  /**
27  * Fill 1D histogram according Poisson statistics with expectation values from given 1D function.
28  *
29  * \param h1 histogram
30  * \param f1 function
31  */
32  template<class T>
33  inline void FillRandom(TH1* h1, const T& f1)
34  {
35  for (Int_t ix = 1; ix <= h1->GetXaxis()->GetNbins(); ++ix){
36 
37  const double x = h1->GetXaxis()->GetBinCenter(ix);
38 
39  h1->SetBinContent(ix, gRandom->Poisson(f1(x)));
40 
41  if (h1->GetSumw2N() != 0) {
42  h1->SetBinError(ix, sqrt(h1->GetBinContent(ix)));
43  }
44  }
45  }
46 
47 
48  /**
49  * Fill 2D histogram according Poisson statistics with expectation values from given 2D function.
50  *
51  * \param h2 histogram
52  * \param f2 function
53  */
54  template<class T>
55  inline void FillRandom(TH2* h2, const T& f2)
56  {
57  for (Int_t ix = 1; ix <= h2->GetXaxis()->GetNbins(); ++ix){
58  for (Int_t iy = 1; iy <= h2->GetYaxis()->GetNbins(); ++iy){
59 
60  const double x = h2->GetXaxis()->GetBinCenter(ix);
61  const double y = h2->GetYaxis()->GetBinCenter(iy);
62 
63  h2->SetBinContent(ix, iy, gRandom->Poisson(f2(x,y)));
64 
65  if (h2->GetSumw2N() != 0) {
66  h2->SetBinError(ix, iy, sqrt(h2->GetBinContent(ix,iy)));
67  }
68  }
69  }
70  }
71 
72 
73  /**
74  * Fill 3D histogram according Poisson statistics with expectation values from given 3D function.
75  *
76  * \param h3 histogram
77  * \param f3 function
78  */
79  template<class T>
80  inline void FillRandom(TH3* h3, const T& f3)
81  {
82  for (Int_t ix = 1; ix <= h3->GetXaxis()->GetNbins(); ++ix){
83  for (Int_t iy = 1; iy <= h3->GetYaxis()->GetNbins(); ++iy){
84  for (Int_t iz = 1; iz <= h3->GetZaxis()->GetNbins(); ++iz){
85 
86  const double x = h3->GetXaxis()->GetBinCenter(ix);
87  const double y = h3->GetYaxis()->GetBinCenter(iy);
88  const double z = h3->GetZaxis()->GetBinCenter(iz);
89 
90  h3->SetBinContent(ix, iy, iz, gRandom->Poisson(f3(x,y,z)));
91 
92  if (h3->GetSumw2N() != 0) {
93  h3->SetBinError(ix, iy, iy, sqrt(h3->GetBinContent(ix,iy,iz)));
94  }
95  }
96  }
97  }
98  }
99 
100 
101  /**
102  * Fill 1D histogram according PDF as given 1D function.
103  *
104  * \param h1 histogram
105  * \param f1 function
106  * \param ns number of entries
107  */
108  template<class T>
109  inline void FillRandom(TH1* h1, const T& f1, const size_t ns)
110  {
111  using namespace std;
112 
113  struct tuple {
114 
115  Double_t x;
116  double W;
117 
118  inline bool operator<(const double W) { return this->W < W; }
119  };
120 
121  vector<tuple> buffer;
122 
123  double W = 0.0;
124 
125  for (Int_t ix = 1; ix <= h1->GetXaxis()->GetNbins(); ++ix){
126 
127  const Double_t x = h1->GetXaxis()->GetBinCenter(ix);
128 
129  W += f1(x);
130 
131  buffer.push_back({x, W});
132  }
133 
134  for (size_t i = 0; i != ns; ++i) {
135 
136  typename vector<tuple>::const_iterator p = lower_bound(buffer.begin(), buffer.end(), W * gRandom->Rndm());
137 
138  h1->Fill(p->x);
139  }
140  }
141 
142 
143  /**
144  * Fill 2D histogram according PDF as given 2D function.
145  *
146  * \param h2 histogram
147  * \param f2 function
148  * \param ns number of entries
149  */
150  template<class T>
151  inline void FillRandom(TH2* h2, const T& f2, const size_t ns)
152  {
153  using namespace std;
154 
155  struct tuple {
156 
157  Double_t x;
158  Double_t y;
159  double W;
160 
161  inline bool operator<(const double W) { return this->W < W; }
162  };
163 
164  vector<tuple> buffer;
165 
166  double W = 0.0;
167 
168  for (Int_t ix = 1; ix <= h2->GetXaxis()->GetNbins(); ++ix){
169  for (Int_t iy = 1; iy <= h2->GetYaxis()->GetNbins(); ++iy){
170 
171  const Double_t x = h2->GetXaxis()->GetBinCenter(ix);
172  const Double_t y = h2->GetYaxis()->GetBinCenter(iy);
173 
174  W += f2(x,y);
175 
176  buffer.push_back({x, y, W});
177  }
178  }
179 
180  for (size_t i = 0; i != ns; ++i) {
181 
182  typename vector<tuple>::const_iterator p = lower_bound(buffer.begin(), buffer.end(), W * gRandom->Rndm());
183 
184  h2->Fill(p->x, p->y);
185  }
186  }
187 
188 
189  /**
190  * Fill 3D histogram according PDF as given 3D function.
191  *
192  * \param h3 histogram
193  * \param f3 function
194  * \param ns number of entries
195  */
196  template<class T>
197  inline void FillRandom(TH3* h3, const T& f3, const size_t ns)
198  {
199  using namespace std;
200 
201  struct tuple {
202 
203  Double_t x;
204  Double_t y;
205  Double_t z;
206  double W;
207 
208  inline bool operator<(const double W) { return this->W < W; }
209  };
210 
211  vector<tuple> buffer;
212 
213  double W = 0.0;
214 
215  for (Int_t ix = 1; ix <= h3->GetXaxis()->GetNbins(); ++ix){
216  for (Int_t iy = 1; iy <= h3->GetYaxis()->GetNbins(); ++iy){
217  for (Int_t iz = 1; iz <= h3->GetYaxis()->GetNbins(); ++iz){
218 
219  const Double_t x = h3->GetXaxis()->GetBinCenter(ix);
220  const Double_t y = h3->GetYaxis()->GetBinCenter(iy);
221  const Double_t z = h3->GetZaxis()->GetBinCenter(iz);
222 
223  W += f3(x,y,z);
224 
225  buffer.push_back({x, y, z, W});
226  }
227  }
228  }
229 
230  for (size_t i = 0; i != ns; ++i) {
231 
232  typename vector<tuple>::const_iterator p = lower_bound(buffer.begin(), buffer.end(), W * gRandom->Rndm());
233 
234  h3->Fill(p->x, p->y, p->z);
235  }
236  }
237 }
238 
239 #endif
240 
double f3(const double x, const double y, const double z)
3D function.
Definition: JPolynome3D.cc:23
const JPolynome f1(1.0, 2.0, 3.0)
Function.
bool operator<(const Head &first, const Head &second)
Less than operator.
Definition: JHead.hh:1817
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Auxiliary classes and methods for ROOT I/O.
void FillRandom(TH3 *h3, const T &f3, const size_t ns)
Fill 3D histogram according PDF as given 3D function.
std::vector< size_t > ns
Definition: JSTDTypes.hh:14