g4tools  5.4.0
p2
Go to the documentation of this file.
1 // Copyright (C) 2010, Guy Barrand. All rights reserved.
2 // See the file tools.license for terms.
3 
4 #ifndef tools_histo_p2
5 #define tools_histo_p2
6 
7 #include "b2"
8 #include "profile_data"
9 
10 namespace tools {
11 namespace histo {
12 
13 template <class TC,class TO,class TN,class TW,class TH,class TV>
14 class p2 : public b2<TC,TO,TN,TW,TH> {
15  typedef b2<TC,TO,TN,TW,TH> parent;
16 public:
18  typedef typename parent::bn_t bn_t;
19  typedef std::vector<TV> vs_t;
20 protected:
21  virtual TH get_bin_height(TO a_offset) const {
22  return (parent::m_bin_Sw[a_offset] ? (m_bin_Svw[a_offset]/parent::m_bin_Sw[a_offset]):0);
23  }
24 public:
25  bool equals(const p2& a_from,const TW& a_prec,TW(*a_fabs)(TW)) const {
26  if(!parent::equals(a_from,a_prec,a_fabs)) return false;
27  if(m_cut_v!=a_from.m_cut_v) return false;
28  if(!numbers_are_equal(m_min_v,a_from.m_min_v,a_prec,a_fabs)) return false;
29  if(!numbers_are_equal(m_max_v,a_from.m_max_v,a_prec,a_fabs)) return false;
30  if(!vectors_are_equal(m_bin_Svw,a_from.m_bin_Svw,a_prec,a_fabs)) return false;
31  if(!vectors_are_equal(m_bin_Sv2w,a_from.m_bin_Sv2w,a_prec,a_fabs)) return false;
32  return true;
33  }
34  bool equals_TH(const p2& a_from,const TW& a_prec,TW(*a_fabs)(TW)) const {
35  if(!parent::equals_TH(a_from,a_prec,a_fabs,false)) return false;
36  if(m_cut_v!=a_from.m_cut_v) return false;
37  if(!numbers_are_equal(m_min_v,a_from.m_min_v,a_prec,a_fabs)) return false;
38  if(!numbers_are_equal(m_max_v,a_from.m_max_v,a_prec,a_fabs)) return false;
39  if(!vectors_are_equal(m_bin_Svw,a_from.m_bin_Svw,a_prec,a_fabs)) return false;
40  if(!vectors_are_equal(m_bin_Sv2w,a_from.m_bin_Sv2w,a_prec,a_fabs)) return false;
41  return true;
42  }
43 
44  virtual TH bin_error(int aI,int aJ) const { //TH should be the same as TV
45  TO offset;
46  if(!parent::_find_offset(aI,aJ,offset)) return 0;
47 
48  //FIXME Is it correct ?
49  // TProfile::GetBinError with kERRORMEAN mode does :
50  // Stat_t cont = fArray[bin]; //Svw (see TProfile::Fill)
51  // Stat_t sum = parent::m_bin_entries.fArray[bin]; //Sw
52  // Stat_t err2 = fSumw2.fArray[bin]; //Sv2w
53  // if (sum == 0) return 0;
54  // Stat_t eprim;
55  // Stat_t contsum = cont/sum;
56  // Stat_t eprim2 = TMath::Abs(err2/sum - contsum*contsum);
57  // eprim = TMath::Sqrt(eprim2);
58  // ... ???
59  // if (fErrorMode == kERRORMEAN) return eprim/TMath::Sqrt(sum);
60 
61  TW sw = parent::m_bin_Sw[offset]; //ROOT sum
62  if(sw==0) return 0;
63  TV svw = m_bin_Svw[offset]; //ROOT cont
64  TV sv2w = m_bin_Sv2w[offset]; //ROOT err2
65  TV mean = (svw / sw); //ROOT contsum
66  TV rms = ::sqrt(::fabs((sv2w/sw) - mean * mean)); //ROOT eprim
67  // rms = get_bin_rms_value.
68  return rms/::sqrt(sw); //ROOT kERRORMEAN mode returned value
69  }
70 
71 public:
72  bool multiply(TW a_factor){
73  if(!parent::base_multiply(a_factor)) return false;
74  for(bn_t ibin=0;ibin<parent::m_bin_number;ibin++) {
75  m_bin_Svw[ibin] *= a_factor;
76  }
77  return true;
78  }
79  bool scale(TW a_factor) {return multiply(a_factor);}
80 
81  TV bin_Svw(int aI,int aJ) const {
82  TO offset;
83  if(!parent::_find_offset(aI,aJ,offset)) return 0;
84  return m_bin_Svw[offset];
85  }
86  TV bin_Sv2w(int aI,int aJ) const {
87  TO offset;
88  if(!parent::_find_offset(aI,aJ,offset)) return 0;
89  return m_bin_Sv2w[offset];
90  }
91 
92  bool reset() {
94  for(bn_t ibin=0;ibin<parent::m_bin_number;ibin++) {
95  m_bin_Svw[ibin] = 0;
96  m_bin_Sv2w[ibin] = 0;
97  }
98  return true;
99  }
100 
101  void copy_from_data(const pd_t& a_from) {
102  parent::base_from_data(a_from);
103  m_bin_Svw = a_from.m_bin_Svw;
104  m_bin_Sv2w = a_from.m_bin_Sv2w;
105  m_cut_v = a_from.m_cut_v;
106  m_min_v = a_from.m_min_v;
107  m_max_v = a_from.m_max_v;
108  }
110  pd_t hd(parent::dac());
111  hd.m_is_profile = true;
112  hd.m_bin_Svw = m_bin_Svw;
113  hd.m_bin_Sv2w = m_bin_Sv2w;
114  hd.m_cut_v = m_cut_v;
115  hd.m_min_v = m_min_v;
116  hd.m_max_v = m_max_v;
117  return hd;
118  }
119 
120  bool fill(TC aX,TC aY,TV aV,TW aWeight = 1) {
121  if(parent::m_dimension!=2) return false;
122 
123  if(m_cut_v) {
124  if( (aV<m_min_v) || (aV>=m_max_v) ) {
125  return true;
126  }
127  }
128 
129  bn_t ibin,jbin;
130  if(!parent::m_axes[0].coord_to_absolute_index(aX,ibin)) return false;
131  if(!parent::m_axes[1].coord_to_absolute_index(aY,jbin)) return false;
132  bn_t offset = ibin + jbin * parent::m_axes[1].m_offset;
133 
134  parent::m_bin_entries[offset]++;
135  parent::m_bin_Sw[offset] += aWeight;
136  parent::m_bin_Sw2[offset] += aWeight * aWeight;
137 
138  TC xw = aX * aWeight;
139  TC x2w = aX * xw;
140  parent::m_bin_Sxw[offset][0] += xw;
141  parent::m_bin_Sx2w[offset][0] += x2w;
142 
143  TC yw = aY * aWeight;
144  TC y2w = aY * yw;
145  parent::m_bin_Sxw[offset][1] += yw;
146  parent::m_bin_Sx2w[offset][1] += y2w;
147 
148  bool inRange = true;
149  if(ibin==0) inRange = false;
150  else if(ibin==(parent::m_axes[0].m_number_of_bins+1)) inRange = false;
151 
152  if(jbin==0) inRange = false;
153  else if(jbin==(parent::m_axes[1].m_number_of_bins+1)) inRange = false;
154 
156  if(inRange) {
157  parent::m_in_range_plane_Sxyw[0] += aX * aY * aWeight;
158 
159  // fast getters :
161  parent::m_in_range_Sw += aWeight;
162  parent::m_in_range_Sw2 += aWeight*aWeight;
163 
164  parent::m_in_range_Sxw[0] += xw;
165  parent::m_in_range_Sx2w[0] += x2w;
166 
167  parent::m_in_range_Sxw[1] += yw;
168  parent::m_in_range_Sx2w[1] += y2w;
169  }
170 
171  // Profile part :
172  TV vw = aV * aWeight;
173  m_bin_Svw[offset] += vw;
174  m_bin_Sv2w[offset] += aV * vw;
175 
176  return true;
177  }
178 
179  TV bin_rms_value(int aI,int aJ) const {
180  TO offset;
181  if(!parent::_find_offset(aI,aJ,offset)) return 0;
182 
183  TW sw = parent::m_bin_Sw[offset];
184  if(sw==0) return 0;
185  TV svw = m_bin_Svw[offset];
186  TV sv2w = m_bin_Sv2w[offset];
187  TV mean = (svw / sw);
188  return ::sqrt(::fabs((sv2w / sw) - mean * mean));
189  }
190 
191  bool add(const p2& a_histo){
192  parent::base_add(a_histo);
193  for(bn_t ibin=0;ibin<parent::m_bin_number;ibin++) {
194  m_bin_Svw[ibin] += a_histo.m_bin_Svw[ibin];
195  m_bin_Sv2w[ibin] += a_histo.m_bin_Sv2w[ibin];
196  }
197  return true;
198  }
199  bool subtract(const p2& a_histo){
200  parent::base_subtract(a_histo);
201  for(bn_t ibin=0;ibin<parent::m_bin_number;ibin++) {
202  m_bin_Svw[ibin] -= a_histo.m_bin_Svw[ibin];
203  m_bin_Sv2w[ibin] -= a_histo.m_bin_Sv2w[ibin];
204  }
205  return true;
206  }
207 
208  bool cut_v() const {return m_cut_v;}
209  TV min_v() const {return m_min_v;}
210  TV max_v() const {return m_max_v;}
211 public:
212  p2(const std::string& a_title,
213  bn_t aXnumber,TC aXmin,TC aXmax,
214  bn_t aYnumber,TC aYmin,TC aYmax)
215  :parent(a_title,aXnumber,aXmin,aXmax,aYnumber,aYmin,aYmax)
216  ,m_cut_v(false)
217  ,m_min_v(0)
218  ,m_max_v(0)
219  {
222  }
223 
224  p2(const std::string& a_title,
225  bn_t aXnumber,TC aXmin,TC aXmax,
226  bn_t aYnumber,TC aYmin,TC aYmax,
227  TV aVmin,TV aVmax)
228  :parent(a_title,aXnumber,aXmin,aXmax,aYnumber,aYmin,aYmax)
229  ,m_cut_v(true)
230  ,m_min_v(aVmin)
231  ,m_max_v(aVmax)
232  {
235  }
236 
237  p2(const std::string& a_title,
238  const std::vector<TC>& a_edges_x,
239  const std::vector<TC>& a_edges_y)
240  :parent(a_title,a_edges_x,a_edges_y)
241  ,m_cut_v(false)
242  ,m_min_v(0)
243  ,m_max_v(0)
244  {
247  }
248 
249  p2(const std::string& a_title,
250  const std::vector<TC>& a_edges_x,
251  const std::vector<TC>& a_edges_y,
252  TV aVmin,TV aVmax)
253  :parent(a_title,a_edges_x,a_edges_y)
254  ,m_cut_v(true)
255  ,m_min_v(aVmin)
256  ,m_max_v(aVmax)
257  {
260  }
261 
262  virtual ~p2(){}
263 public:
264  p2(const p2& a_from)
265  :parent(a_from)
266  ,m_cut_v(a_from.m_cut_v)
267  ,m_min_v(a_from.m_min_v)
268  ,m_max_v(a_from.m_max_v)
269  ,m_bin_Svw(a_from.m_bin_Svw)
270  ,m_bin_Sv2w(a_from.m_bin_Sv2w)
271  {}
272  p2& operator=(const p2& a_from){
273  parent::operator=(a_from);
274  m_cut_v = a_from.m_cut_v;
275  m_min_v = a_from.m_min_v;
276  m_max_v = a_from.m_max_v;
277  m_bin_Svw = a_from.m_bin_Svw;
278  m_bin_Sv2w = a_from.m_bin_Sv2w;
279  return *this;
280  }
281 public:
282  bool configure(bn_t aXnumber,TC aXmin,TC aXmax,bn_t aYnumber,TC aYmin,TC aYmax){
283  if(!parent::configure(aXnumber,aXmin,aXmax,aYnumber,aYmin,aYmax)) return false;
284  m_bin_Svw.clear();
285  m_bin_Sv2w.clear();
288  m_cut_v = false;
289  m_min_v = 0;
290  m_max_v = 0;
291  return true;
292  }
293  bool configure(const std::vector<TC>& a_edges_x,const std::vector<TC>& a_edges_y) {
294  if(!parent::configure(a_edges_x,a_edges_y)) return false;
295  m_bin_Svw.clear();
296  m_bin_Sv2w.clear();
299  m_cut_v = false;
300  m_min_v = 0;
301  m_max_v = 0;
302  return true;
303  }
304  bool configure(bn_t aXnumber,TC aXmin,TC aXmax,bn_t aYnumber,TC aYmin,TC aYmax,TV aVmin,TV aVmax){
305  if(!parent::configure(aXnumber,aXmin,aXmax,aYnumber,aYmin,aYmax)) return false;
306  m_bin_Svw.clear();
307  m_bin_Sv2w.clear();
310  m_cut_v = true;
311  m_min_v = aVmin;
312  m_max_v = aVmax;
313  return true;
314  }
315  bool configure(const std::vector<TC>& a_edges_x,const std::vector<TC>& a_edges_y,TV aVmin,TV aVmax) {
316  if(!parent::configure(a_edges_x,a_edges_y)) return false;
317  m_bin_Svw.clear();
318  m_bin_Sv2w.clear();
321  m_cut_v = true;
322  m_min_v = aVmin;
323  m_max_v = aVmax;
324  return true;
325  }
326 public:
327  const vs_t& bins_sum_vw() const {return m_bin_Svw;}
328  const vs_t& bins_sum_v2w() const {return m_bin_Sv2w;}
329 
330  TW get_Svw() const {
331  TW sw = 0;
332  for(TO ibin=0;ibin<parent::m_bin_number;ibin++) {
333  if(!histo::is_out(parent::m_axes,ibin)) {
334  sw += m_bin_Svw[ibin];
335  }
336  }
337  return sw;
338  }
339  TW get_Sv2w() const {
340  TW sw = 0;
341  for(TO ibin=0;ibin<parent::m_bin_number;ibin++) {
342  if(!histo::is_out(parent::m_axes,ibin)) {
343  sw += m_bin_Sv2w[ibin];
344  }
345  }
346  return sw;
347  }
348 protected:
349  bool m_cut_v;
354 };
355 
356 }}
357 
358 #endif
359 
tools::histo::histo_data::m_in_range_Sw2
TW m_in_range_Sw2
Definition: histo_data:189
tools::histo::p2::m_min_v
TV m_min_v
Definition: p2:350
tools::histo::profile_data::m_is_profile
bool m_is_profile
Definition: profile_data:73
tools::histo::p2::p2
p2(const std::string &a_title, const std::vector< TC > &a_edges_x, const std::vector< TC > &a_edges_y)
Definition: p2:237
tools::histo::p2::m_bin_Svw
vs_t m_bin_Svw
Definition: p2:352
tools::histo::p2::equals_TH
bool equals_TH(const p2 &a_from, const TW &a_prec, TW(*a_fabs)(TW)) const
Definition: p2:34
tools::histo::p2::bins_sum_vw
const vs_t & bins_sum_vw() const
Definition: p2:327
tools::histo::histo_data::m_bin_Sx2w
std::vector< std::vector< TC > > m_bin_Sx2w
Definition: histo_data:179
tools::histo::p2::p2
p2(const std::string &a_title, bn_t aXnumber, TC aXmin, TC aXmax, bn_t aYnumber, TC aYmin, TC aYmax)
Definition: p2:212
tools::histo::histo_data::m_in_range_Sxw
std::vector< TC > m_in_range_Sxw
Definition: histo_data:190
tools::histo::p2::m_max_v
TV m_max_v
Definition: p2:351
tools::histo::p2::get_Svw
TW get_Svw() const
Definition: p2:330
tools::histo::p2::p2
p2(const std::string &a_title, const std::vector< TC > &a_edges_x, const std::vector< TC > &a_edges_y, TV aVmin, TV aVmax)
Definition: p2:249
tools::histo::base_histo::dac
const hd_t & dac() const
Definition: base_histo:56
tools::histo::p2::configure
bool configure(const std::vector< TC > &a_edges_x, const std::vector< TC > &a_edges_y, TV aVmin, TV aVmax)
Definition: p2:315
tools::histo::p2::m_cut_v
bool m_cut_v
Definition: p2:349
b2
tools::histo::b2::operator=
b2 & operator=(const b2 &a_from)
Definition: b2:317
tools::histo::histo_data::m_bin_entries
std::vector< TN > m_bin_entries
Definition: histo_data:175
tools::histo::p2::m_bin_Sv2w
vs_t m_bin_Sv2w
Definition: p2:353
tools::histo::base_histo::base_add
void base_add(const base_histo &a_histo)
Definition: base_histo:476
tools::histo::p2::bin_Svw
TV bin_Svw(int aI, int aJ) const
Definition: p2:81
tools::histo::p2::bins_sum_v2w
const vs_t & bins_sum_v2w() const
Definition: p2:328
tools::histo::histo_data::m_bin_Sw
std::vector< TW > m_bin_Sw
Definition: histo_data:176
tools::histo::profile_data::m_bin_Sv2w
std::vector< TV > m_bin_Sv2w
Definition: profile_data:75
tools::histo::p2::min_v
TV min_v() const
Definition: p2:209
tools::histo::p2::operator=
p2 & operator=(const p2 &a_from)
Definition: p2:272
tools::histo::p2
Definition: p2:14
tools::histo::p2::p2
p2(const p2 &a_from)
Definition: p2:264
tools::histo::profile_data::m_min_v
TV m_min_v
Definition: profile_data:77
tools::histo::p2::bin_error
virtual TH bin_error(int aI, int aJ) const
Definition: p2:44
tools::histo::p2::configure
bool configure(bn_t aXnumber, TC aXmin, TC aXmax, bn_t aYnumber, TC aYmin, TC aYmax, TV aVmin, TV aVmax)
Definition: p2:304
tools::histo::histo_data::m_in_range_entries
TN m_in_range_entries
Definition: histo_data:187
tools::vectors_are_equal
bool vectors_are_equal(const VEC &a_1, const VEC &a_2, const PREC &a_prec, PREC(*a_fabs)(const PREC &))
Definition: eqT:23
tools::histo::p2::vs_t
std::vector< TV > vs_t
Definition: p2:19
tools::histo::profile_data::m_cut_v
bool m_cut_v
Definition: profile_data:76
tools::histo::base_histo::base_subtract
void base_subtract(const base_histo &a_histo)
Definition: base_histo:493
tools::histo::histo_data::m_in_range_Sx2w
std::vector< TC > m_in_range_Sx2w
Definition: histo_data:191
tools::histo::b2
Definition: b2:15
tools::histo::p2::bin_rms_value
TV bin_rms_value(int aI, int aJ) const
Definition: p2:179
tools::histo::profile_data
Definition: profile_data:13
tools::histo::p2::configure
bool configure(const std::vector< TC > &a_edges_x, const std::vector< TC > &a_edges_y)
Definition: p2:293
tools::histo::profile_data::m_bin_Svw
std::vector< TV > m_bin_Svw
Definition: profile_data:74
tools::histo::is_out
bool is_out(const std::vector< axis< TC, TO > > &a_axes, TO a_offset)
Definition: axes:16
tools::histo::p2::multiply
bool multiply(TW a_factor)
Definition: p2:72
tools::histo::p2::copy_from_data
void copy_from_data(const pd_t &a_from)
Definition: p2:101
tools::histo::base_histo::base_from_data
void base_from_data(const hd_t &a_from)
Definition: base_histo:47
tools::histo::histo_data::m_bin_Sxw
std::vector< std::vector< TC > > m_bin_Sxw
Definition: histo_data:178
tools::histo::b2::configure
bool configure(bn_t aXnumber, TC aXmin, TC aXmax, bn_t aYnumber, TC aYmin, TC aYmax)
Definition: b2:319
profile_data
tools::histo::base_histo< double, unsigned int, unsigned int, double, double >::bn_t
axis_t::bn_t bn_t
Definition: base_histo:37
tools::histo::histo_data::m_axes
std::vector< axis_t > m_axes
Definition: histo_data:181
tools::histo::p2::fill
bool fill(TC aX, TC aY, TV aV, TW aWeight=1)
Definition: p2:120
tools::histo::p2::add
bool add(const p2 &a_histo)
Definition: p2:191
tools::histo::p2::equals
bool equals(const p2 &a_from, const TW &a_prec, TW(*a_fabs)(TW)) const
Definition: p2:25
tools::histo::base_histo::equals
bool equals(const base_histo &a_from, const TW &a_prec, TW(*a_fabs)(TW)) const
Definition: base_histo:83
tools
inlined C code : ///////////////////////////////////
Definition: aida_ntuple:26
tools::histo::base_histo::base_reset
void base_reset()
Definition: base_histo:381
tools::histo::histo_data::m_bin_Sw2
std::vector< TW > m_bin_Sw2
Definition: histo_data:177
tools::histo::p2::reset
bool reset()
Definition: p2:92
tools::histo::histo_data::m_all_entries
TN m_all_entries
Definition: histo_data:186
tools::histo::histo_data::m_in_range_Sw
TW m_in_range_Sw
Definition: histo_data:188
tools::numbers_are_equal
bool numbers_are_equal(const NUMBER &a_left, const NUMBER &a_right, const PREC &a_prec, PREC(*a_fabs)(const NUMBER &))
Definition: eqT:10
tools::histo::p2::max_v
TV max_v() const
Definition: p2:210
tools::histo::p2::p2
p2(const std::string &a_title, bn_t aXnumber, TC aXmin, TC aXmax, bn_t aYnumber, TC aYmin, TC aYmax, TV aVmin, TV aVmax)
Definition: p2:224
tools::histo::p2::get_Sv2w
TW get_Sv2w() const
Definition: p2:339
tools::histo::p2::bn_t
parent::bn_t bn_t
Definition: p2:18
tools::histo::p2::~p2
virtual ~p2()
Definition: p2:262
tools::histo::p2::subtract
bool subtract(const p2 &a_histo)
Definition: p2:199
tools::histo::b2::_find_offset
bool _find_offset(int aI, int aJ, TO &a_offset) const
Definition: b2:340
tools::histo::p2::configure
bool configure(bn_t aXnumber, TC aXmin, TC aXmax, bn_t aYnumber, TC aYmin, TC aYmax)
Definition: p2:282
tools::histo::p2::get_histo_data
pd_t get_histo_data() const
Definition: p2:109
tools::histo::histo_data::m_bin_number
TO m_bin_number
Definition: histo_data:174
tools::histo::base_histo
Definition: base_histo:27
tools::histo::p2::bin_Sv2w
TV bin_Sv2w(int aI, int aJ) const
Definition: p2:86
tools::histo::p2::pd_t
profile_data< TC, TO, TN, TW, TV > pd_t
Definition: p2:17
tools::histo::histo_data::m_in_range_plane_Sxyw
std::vector< TC > m_in_range_plane_Sxyw
Definition: histo_data:183
tools::histo::profile_data::m_max_v
TV m_max_v
Definition: profile_data:78
tools::histo::histo_data::m_dimension
dim_t m_dimension
Definition: histo_data:172
tools::histo::p2::scale
bool scale(TW a_factor)
Definition: p2:79
tools::histo::histo_data::equals_TH
bool equals_TH(const histo_data &a_from, const TW &a_prec, TW(*a_fabs)(TW), bool a_cmp_bin_Sw2) const
Definition: histo_data:144
tools::histo::base_histo::base_multiply
bool base_multiply(const base_histo &a_histo)
Definition: base_histo:513
tools::histo::p2::cut_v
bool cut_v() const
Definition: p2:208
tools::histo::p2::get_bin_height
virtual TH get_bin_height(TO a_offset) const
Definition: p2:21
tools::histo::b2< double, unsigned int, unsigned int, double, double >::bn_t
parent::bn_t bn_t
Definition: b2:22