ATLAS Offline Software
FlexErrArray.icc
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 
6 #ifndef NDEBUG
7 #include <iostream>
8 #endif
9 
10 //____________________________________________________________________
11 template <class T>
12 inline FlexErrArray<T>::FlexErrArray( unsigned nbins )
13  : m_indices(reinterpret_cast<char*>(this)+sizeof(*this)),
14  m_fastloop_group2check(0),
15  m_nbins(nbins),
16  m_fastloop_isuper2check(0),
17  m_fastloop_igr2check(0)
18 #ifdef LW_STRICT_ROOT_BEHAVIOUR
19  , m_pretendSumWMode(false)
20 #endif
21 {
22  //Check that the implementation makes sense at all:
23  assert(sizeof(unsigned short)>sizeof(unsigned char));
24  assert(sizeof(T)>sizeof(unsigned short));
25  if (small()) {
26  memset(groups(),0,nGroupsNeeded(nbins)*sizeof(void*));
27  // cppcheck-suppress assertWithSideEffect
28  assert(!groups()[0]);
29  // cppcheck-suppress assertWithSideEffect
30  assert(!groups()[nGroupsNeeded(nbins)-1]);
31  } else {
32  memset(superGroups(),0,nSuperGroups(nbins)*sizeof(void*));
33  // cppcheck-suppress assertWithSideEffect
34  assert(!superGroups()[0]);
35  // cppcheck-suppress assertWithSideEffect
36  assert(!superGroups()[nSuperGroups(nbins)-1]);
37  }
38 }
39 
40 //____________________________________________________________________
41 template <class T>
42 inline FlexErrArray<T>::~FlexErrArray()
43 {
44  if (small()) {
45  const unsigned n(nGroupsNeeded(m_nbins));
46  FlexErrArrayGroup<T> ** grps(groups());
47  for (unsigned igr = 0; igr < n; ++igr) {
48  FlexErrArrayGroup<T> *gr = grps[igr];
49  if (gr)
50  MP_DELETE(gr);
51  }
52  } else {
53  unsigned nsuper(nSuperGroups(m_nbins));
54  FlexErrArrayGroup<T> *** supergroups(superGroups());
55  for (unsigned isuper = 0; isuper<nsuper;++isuper) {
56  FlexErrArrayGroup<T> ** supergroup = supergroups[isuper];
57  if (supergroup) {
58  unsigned upper(isuper+1==nsuper?entriesInSuperGroup(isuper):FLEX1DERRARRAY_NGROUPSPERINDEX);
59  for (unsigned igr = 0; igr<upper;++igr) {
60  FlexErrArrayGroup<T> *gr = supergroup[igr];
61  if (gr)
62  MP_DELETE(gr);
63  }
64  if (upper==FLEX1DERRARRAY_NGROUPSPERINDEX)
65  LWPools::release<FlexErrArrayGroup<T>*,FLEX1DERRARRAY_NGROUPSPERINDEX>(supergroup);
66  else
67  LWPools::release<FlexErrArrayGroup<T>*>(supergroup,upper);
68  }
69  }
70  }
71 }
72 
73 //____________________________________________________________________
74 template <class T>
75 inline unsigned FlexErrArray<T>::getNBins() const
76 {
77  return m_nbins;
78 }
79 
80 //____________________________________________________________________
81 template <class T>
82 inline void FlexErrArray<T>::copyContents(T*__restrict__ cont,double*__restrict__ err) const
83 {
84  //It is faster to null-out the entire arrays at the beginning, and
85  //then not have to null-out the missing pieces later:
86  memset(&(cont[0]),0,m_nbins*sizeof(*cont));
87  if (err)
88  memset(&(err[0]),0,m_nbins*sizeof(*err));
89  //Process group for group:
90  if (small()) {
91  const unsigned n(nGroupsNeeded(m_nbins));
92  const FlexErrArrayGroup<T> * const* grps(groups());
93  for (unsigned igr = 0; igr < n; ++igr) {
94  const FlexErrArrayGroup<T> *gr = grps[igr];
95  if (gr) {
96  unsigned low(FLEXERRARRAYGROUP_MAXBINS*igr);
97  assert(low+gr->getNBins()<=m_nbins);
98  gr->copyContents(&(cont[low]),(err?&(err[low]):0));
99  }
100  }
101  } else {
102  unsigned nsuper(nSuperGroups(m_nbins));
103  const FlexErrArrayGroup<T> * const* const* supergroups(superGroups());
104  for (unsigned isuper = 0; isuper<nsuper;++isuper) {
105  const FlexErrArrayGroup<T> * const* supergroup = supergroups[isuper];
106  if (supergroup) {
107  unsigned upper(isuper+1==nsuper?entriesInSuperGroup(isuper):FLEX1DERRARRAY_NGROUPSPERINDEX);
108  for (unsigned igr = 0; igr<upper;++igr) {
109  const FlexErrArrayGroup<T> *gr = supergroup[igr];
110  if (gr) {
111  unsigned low(FLEXERRARRAYGROUP_MAXBINS*(FLEX1DERRARRAY_NGROUPSPERINDEX*isuper+igr));
112  assert(low+gr->getNBins()<=m_nbins);
113  gr->copyContents(&(cont[low]),(err?&(err[low]):0));
114  }
115  }
116  }
117  }
118  }
119 #ifndef NDEBUG
120  for (unsigned ibin=0;ibin<m_nbins;++ibin)
121  assert(cont[ibin]==getBinContent(ibin));
122  if (err)
123  for (unsigned ibin=0;ibin<m_nbins;++ibin) {
124  const double e1(sqrt(fabs(err[ibin]))),e2(getBinError(ibin));
125  if ((fabs(e1-e2)/(1+std::max(e1,e2)))>1.0e-5&&!bothNaN(e1,e2)) {
126  std::cout<<"e1 = "<<e1<<std::endl;
127  std::cout<<"e2 = "<<e2<<std::endl;
128  std::cout<<"test: "<<(e1-e2)/(1.0+std::max(e1,e2))<<std::endl;
129  assert((e1-e2)/(1.0+std::max(e1,e2))<1.0e-3||bothNaN(e1,e2));//fixme: Investigate if really need to go up to 1.0e-3
130  }
131  }
132 #endif
133 }
134 
135 //____________________________________________________________________
136 template <class T>
137 inline double FlexErrArray<T>::Integral() const
138 {
139  double result(0);
140  if (small()) {
141  const unsigned n(nGroupsNeeded(m_nbins));
142  const FlexErrArrayGroup<T> * const* grps(groups());
143  for (unsigned igr = 0; igr < n; ++igr) {
144  const FlexErrArrayGroup<T> *gr = grps[igr];
145  if (gr)
146  result += gr->Integral();
147  }
148  } else {
149  unsigned nsuper(nSuperGroups(m_nbins));
150  const FlexErrArrayGroup<T> * const* const* supergroups(superGroups());
151  for (unsigned isuper = 0; isuper<nsuper;++isuper) {
152  const FlexErrArrayGroup<T> * const* supergroup = supergroups[isuper];
153  if (supergroup) {
154  unsigned upper(isuper+1==nsuper?entriesInSuperGroup(isuper):FLEX1DERRARRAY_NGROUPSPERINDEX);
155  for (unsigned igr = 0; igr<upper;++igr) {
156  const FlexErrArrayGroup<T> *gr = supergroup[igr];
157  if (gr)
158  result += gr->Integral();
159  }
160  }
161  }
162  }
163  return result;
164 }
165 
166 //____________________________________________________________________
167 template <class T>
168 inline bool FlexErrArray<T>::holdsSeparateSumW2Info() const
169 {
170 #ifdef LW_STRICT_ROOT_BEHAVIOUR
171  return m_pretendSumWMode;
172 #else
173  unsigned nsuper(nSuperGroups(m_nbins));
174  if (small()) {
175  const unsigned n(nGroupsNeeded(m_nbins));
176  const FlexErrArrayGroup<T> ** grps(groups());
177  for (unsigned igr = 0; igr < n; ++igr)
178  if (grps[igr]&&grps[igr]->holdsSeparateSumW2Info())
179  return true;
180  return false;
181  } else {
182  const FlexErrArrayGroup<T> *** supergroups(superGroups());
183  for (unsigned isuper = 0; isuper<nsuper;++isuper) {
184  const FlexErrArrayGroup<T> ** supergroup = supergroups[isuper];
185  if (supergroup) {
186  unsigned upper(isuper+1==nsuper?entriesInSuperGroup(isuper):FLEX1DERRARRAY_NGROUPSPERINDEX);
187  for (unsigned igr = 0; igr<upper;++igr) {
188  const FlexErrArrayGroup<T> *gr = supergroup[igr];
189  if (gr&&gr->holdsSeparateSumW2Info())
190  return true;
191  }
192  }
193  }
194  return false;
195  }
196 #endif
197 }
198 
199 //____________________________________________________________________
200 template <class T>
201 inline T FlexErrArray<T>::getBinContent(unsigned bin) const
202 {
203  assert(bin<m_nbins);
204  const FlexErrArrayGroup<T> * gr = getGroupNoAlloc(bin);
205  return gr ? gr->getBinContent(groupBin(bin)) : 0;
206 }
207 
208 
209 //____________________________________________________________________
210 template <class T>
211 inline double FlexErrArray<T>::getBinError(unsigned bin) const
212 {
213  assert(bin<m_nbins);
214  const FlexErrArrayGroup<T> * gr = getGroupNoAlloc(bin);
215  return gr ? gr->getBinError(groupBin(bin)) : 0;
216 }
217 
218 //____________________________________________________________________
219 template <class T>
220 inline void FlexErrArray<T>::getBinContentAndError(unsigned bin, T& content, double& error ) const
221 {
222  assert(bin<m_nbins);
223  const FlexErrArrayGroup<T> * gr = getGroupNoAlloc(bin);
224  if (gr) {
225  gr->getBinContentAndError(groupBin(bin),content,error);
226  } else {
227  content = 0;
228  error = 0;
229  }
230 }
231 
232 //____________________________________________________________________
233 template <class T>
234 inline void FlexErrArray<T>::setBinContent(unsigned bin, const T& content)
235 {
236  assert(bin<m_nbins);
237  FlexErrArrayGroup<T> * gr;
238  if (content==0) {
239  gr = getGroupNoAlloc(bin);
240  if (!gr)
241  return;
242  } else {
243  gr = getGroup(bin);
244  }
245  assert(gr);
246  gr->setBinContent(groupBin(bin),content STRICT_ROOT_PAR(m_pretendSumWMode));
247  assert(getBinContent(bin)==content||bothNaN(getBinContent(bin),content));
248 }
249 
250 
251 //____________________________________________________________________
252 template <class T>
253 inline void FlexErrArray<T>::setBinError(unsigned bin, const double& error )
254 {
255  assert(bin<m_nbins);
256 #ifdef LW_STRICT_ROOT_BEHAVIOUR
257  m_pretendSumWMode = true;
258 #endif
259  FlexErrArrayGroup<T> * gr;
260  if (error==0) {
261  gr = getGroupNoAlloc(bin);
262  if (!gr)
263  return;
264  } else {
265  gr = getGroup(bin);
266  }
267  assert(gr);
268  gr->setBinError(groupBin(bin),error);
269  assert(fabs(getBinError(bin)-fabs(error))<1.0e-10||bothNaN(error,getBinError(bin)));
270 }
271 
272 //____________________________________________________________________
273 template <class T>
274 inline void FlexErrArray<T>::setBinContentAndError(unsigned bin, const T& content, const double& error )
275 {
276  assert(bin<m_nbins);
277 #ifdef LW_STRICT_ROOT_BEHAVIOUR
278  m_pretendSumWMode = true;
279 #endif
280 
281  FlexErrArrayGroup<T> * gr;
282  if (content==0&&error==0) {
283  //Special case: Maybe we do not need an allocation:
284  gr = getGroupNoAlloc(bin);
285  if (!gr)
286  return;
287  } else {
288  gr = getGroup(bin);
289  }
290  assert(gr);
291  gr->setBinContentAndError( groupBin(bin),
292  content,
293  error );
294 }
295 
296 //____________________________________________________________________
297 template <class T>
298 inline void FlexErrArray<T>::fill(unsigned bin)
299 {
300 #ifndef NDEBUG
301  T expectedval = getBinContent(bin)+1;
302 #endif
303  assert(bin<m_nbins);
304  getGroup(bin)->fill(groupBin(bin) STRICT_ROOT_PAR(m_pretendSumWMode));
305 #ifndef NDEBUG
306  T resval=getBinContent(bin);
307  assert(static_cast<double>(std::abs(resval-expectedval))/(1.0+std::max<T>(std::abs(resval),std::abs(expectedval)))<1.0e-5);
308 #endif
309 }
310 
311 //____________________________________________________________________
312 template <class T>
313 inline void FlexErrArray<T>::fill(unsigned bin, const double& weight)
314 {
315 #ifndef NDEBUG
316  T expectedval = getBinContent(bin)+T(weight);
317 #endif
318  assert(bin<m_nbins);
319 #ifdef LW_STRICT_ROOT_BEHAVIOUR
320  if (weight != 1) m_pretendSumWMode = true;
321 #endif
322  getGroup(bin)->fill(groupBin(bin),weight STRICT_ROOT_PAR(m_pretendSumWMode));
323 #ifndef NDEBUG
324  T resval=getBinContent(bin);
325  assert(static_cast<double>(std::abs(resval-expectedval))/(1.0+std::max<T>(std::abs(resval),std::abs(expectedval)))<1.0e-5);
326 #endif
327 }
328 
329 //____________________________________________________________________
330 template <class T>
331 inline FlexErrArrayGroup<T> * FlexErrArray<T>::getGroup(unsigned bin)
332 {
333  FlexErrArrayGroup<T> ** gr(0);
334  unsigned isuper(0);
335  unsigned igr(0);
336  if (small()) {
337  igr = iGroup(bin);
338  gr = &(groups()[igr]);
339  } else {
340  //1) Get the super-group
341  isuper = iSuperGroup(bin);
342  FlexErrArrayGroup<T> ** supergroup = superGroups()[isuper];
343  if (!supergroup) {
344  unsigned l(entriesInSuperGroup(isuper));
345  if (l==FLEX1DERRARRAY_NGROUPSPERINDEX)
346  supergroup = LWPools::acquire<FlexErrArrayGroup<T>*,FLEX1DERRARRAY_NGROUPSPERINDEX>();
347  else
348  supergroup = LWPools::acquire<FlexErrArrayGroup<T>*>(l);
349  memset(supergroup,0,l*sizeof(supergroup[0]));
350  superGroups()[iSuperGroup(bin)] = supergroup;
351  assert(superGroupBin(bin)<l);
352  }
353  //2) Get the group
354  igr = superGroupBin(bin);
355  gr = &(supergroup[igr]);
356  }
357  assert(gr);
358 
359  if (!(*gr)) {
360  unsigned nbinsGr(FLEXERRARRAYGROUP_MAXBINS);
361  unsigned lastBin(FLEXERRARRAYGROUP_MAXBINS*(isuper*FLEX1DERRARRAY_NGROUPSPERINDEX+igr+1));
362  assert(!(unsigned(igr+1)==nGroupsNeeded(m_nbins))||lastBin>=m_nbins);
363  assert(m_nbins+FLEXERRARRAYGROUP_MAXBINS>lastBin);
364  if (lastBin>m_nbins)
365  nbinsGr -= (lastBin-m_nbins);
366  FlexErrArrayGroup<T> * g = MP_NEW(FlexErrArrayGroup<T>)(nbinsGr);
367  *gr = g;
368  assert((*gr)->getNBins()==nbinsGr);
369  }
370  assert(*gr);
371  return *gr;
372 }
373 
374 
375 
376 //____________________________________________________________________
377 template <class T>
378 inline void FlexErrArray<T>::fastLoop_findAndResetNextGroup( )
379 {
380  //Fixme: We should actually always find the two next, and use
381  //__builtin_prefetch() on the second. We could also call "prefetch
382  //all data" on the group (since we know we will need all chunks/errchunks)
383 
384  assert(!m_fastloop_group2check);
385  if (small()) {
386  const unsigned n(nGroupsNeeded(m_nbins));
387  FlexErrArrayGroup<T> ** grps(groups());
388  for (;m_fastloop_igr2check < n; ++m_fastloop_igr2check) {
389  m_fastloop_group2check = grps[m_fastloop_igr2check];
390  if (m_fastloop_group2check)
391  break;
392  }
393  } else {
394  unsigned nsuper(nSuperGroups(m_nbins));
395  FlexErrArrayGroup<T> *** supergroups(superGroups());
396  for (;m_fastloop_isuper2check<nsuper;++m_fastloop_isuper2check) {
397  FlexErrArrayGroup<T> ** supergroup = supergroups[m_fastloop_isuper2check];
398  if (supergroup) {
399  unsigned upper(m_fastloop_isuper2check+1==nsuper?entriesInSuperGroup(m_fastloop_isuper2check):FLEX1DERRARRAY_NGROUPSPERINDEX);
400  for (;m_fastloop_igr2check < upper; ++m_fastloop_igr2check) {
401  m_fastloop_group2check = supergroup[m_fastloop_igr2check];
402  if (m_fastloop_group2check)
403  break;
404  }
405  if (m_fastloop_group2check)
406  break;
407  }
408  m_fastloop_igr2check = 0;
409  }
410  }
411 
412  if (m_fastloop_group2check)
413  m_fastloop_group2check->resetActiveBinLoop();
414 }
415 
416 
417 //____________________________________________________________________
418 template <class T>
419 inline void FlexErrArray<T>::resetActiveBinLoop()
420 {
421  m_fastloop_isuper2check = 0;
422  m_fastloop_igr2check = 0;
423  m_fastloop_group2check = 0;
424  fastLoop_findAndResetNextGroup();
425 }
426 
427 //____________________________________________________________________
428 template <class T>
429 inline bool FlexErrArray<T>::getNextActiveBin(unsigned& bin, T& content, double& error)
430 {
431  if (m_fastloop_group2check) {
432  //Just check the group once more:
433  if (m_fastloop_group2check->getNextActiveBin(bin, content, error)) {
434  //"bin" is now local. Put into global mode:
435  bin += ( m_fastloop_isuper2check * FLEX1DERRARRAY_NGROUPSPERINDEX + m_fastloop_igr2check) * FLEXERRARRAYGROUP_MAXBINS;//fixme: cache this offset?
436  return true;
437  }
438  }
439  m_fastloop_group2check=0;
440 
441  //Attempt to find the next group (if any):
442  ++m_fastloop_igr2check;
443  if (m_fastloop_igr2check==FLEX1DERRARRAY_NGROUPSPERINDEX&&!small()) {
444  ++m_fastloop_isuper2check;
445  m_fastloop_igr2check = 0;
446  }
447  fastLoop_findAndResetNextGroup();
448 
449  //Repeat:
450  return m_fastloop_group2check ? getNextActiveBin(bin,content,error) : false;
451 }
452 
453 //____________________________________________________________________
454 template <class T>
455 inline void FlexErrArray<T>::scaleContentsAndErrors( const double& fact )
456 {
457  //We simply use the fast-looping functionality here. Just to be
458  //safe, we save the state of the fastlooping temporary variables:
459  unsigned save_fastloop_isuper2check = m_fastloop_isuper2check;
460  unsigned save_fastloop_igr2check = m_fastloop_igr2check;
461 
462  unsigned bin; T cont; double err;
463  resetActiveBinLoop();
464  while (getNextActiveBin(bin, cont, err))
465  setBinContentAndError(bin, static_cast<T>(cont*fact), err*fact );
466 
467  m_fastloop_isuper2check = save_fastloop_isuper2check;
468  m_fastloop_igr2check = save_fastloop_igr2check;
469 }
Generated on Mon Jul 8 2024 21:15:47 for ATLAS Offline Software by   doxygen 1.8.18