AliPhysics  251aa1e (251aa1e)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
THistManager.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2016, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 #include <cfloat>
16 #include <cstring>
17 #include <iostream> // for unit tests
18 #include <sstream>
19 #include <string>
20 #include <exception>
21 #include <vector>
22 #include <TArrayD.h>
23 #include <TAxis.h>
24 #include <TH1.h>
25 #include <TH2.h>
26 #include <TH3.h>
27 #include <THnSparse.h>
28 #include <THashList.h>
29 #include <TMath.h>
30 #include <TObjArray.h>
31 #include <TObjString.h>
32 #include <TProfile.h>
33 #include <TString.h>
34 
35 #include "TBinning.h"
36 #include "THistManager.h"
37 
41 
43  TNamed(),
44  fHistos(NULL),
45  fIsOwner(true)
46 {
47 }
48 
49 THistManager::THistManager(const char *name):
50  TNamed(name, Form("Histogram container %s", name)),
51  fHistos(NULL),
52  fIsOwner(true)
53 {
54  fHistos = new THashList();
55  fHistos->SetName(Form("histos%s", name));
56  fHistos->SetOwner();
57 }
58 
60  if(fHistos && fIsOwner) delete fHistos;
61 }
62 
63 THashList* THistManager::CreateHistoGroup(const char *groupname) {
64  TString parentname = basename(groupname);
65  THashList *parentgroup = FindGroup(parentname);
66  if(!parentgroup) parentgroup = CreateHistoGroup(parentname);
67  THashList *childgroup = new THashList();
68  childgroup->SetName(histname(groupname));
69  childgroup->SetOwner();
70  parentgroup->Add(childgroup);
71  return childgroup;
72 }
73 
74 TH1* THistManager::CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt){
75  TString dirname(basename(name)), hname(histname(name));
76  THashList *parent(FindGroup(dirname));
77  if(!parent) parent = CreateHistoGroup(dirname);
78  if(parent->FindObject(hname)){
79  Fatal("THistManager::CreateTH1", "Object %s already exists in group %s", hname.Data(), dirname.Data());
80  return 0;
81  }
82  TH1* h = new TH1D(hname, title, nbins, xmin, xmax);
83  TString optionstring(opt);
84  optionstring.ToLower();
85  if(optionstring.Contains("s"))
86  h->Sumw2();
87  parent->Add(h);
88  return h;
89 }
90 
91 TH1* THistManager::CreateTH1(const char *name, const char *title, int nbins, const double *xbins, Option_t *opt){
92  TString dirname(basename(name)), hname(histname(name));
93  THashList *parent(FindGroup(dirname));
94  if(!parent) parent = CreateHistoGroup(dirname);
95  if(parent->FindObject(hname)){
96  Fatal("THistManager::CreateTH1", "Object %s already exists in group %s", hname.Data(), dirname.Data());
97  return 0;
98  }
99  TH1* h = new TH1D(hname, title, nbins, xbins);
100  TString optionstring(opt);
101  optionstring.ToLower();
102  if(optionstring.Contains("s"))
103  h->Sumw2();
104  parent->Add(h);
105  return h;
106 }
107 
108 TH1* THistManager::CreateTH1(const char *name, const char *title, const TArrayD &xbins, Option_t *opt){
109  TString dirname(basename(name)), hname(histname(name));
110  THashList *parent(FindGroup(dirname));
111  if(!parent) parent = CreateHistoGroup(dirname);
112  if(parent->FindObject(hname)){
113  Fatal("THistManager::CreateTH1", "Object %s already exists in group %s", hname.Data(), dirname.Data());
114  return 0;
115  }
116  TH1* h = new TH1D(hname, title, xbins.GetSize()-1, xbins.GetArray());
117  TString optionstring(opt);
118  optionstring.ToLower();
119  if(optionstring.Contains("s"))
120  h->Sumw2();
121  parent->Add(h);
122  return h;
123 }
124 
125 TH1* THistManager::CreateTH1(const char *name, const char *title, const TBinning &xbin, Option_t *opt){
126  TArrayD myxbins;
127  try{
128  xbin.CreateBinEdges(myxbins);
129  } catch(std::exception &e){
130  Fatal("THistManager::CreateTH1", "Exception raised: %s", e.what());
131  }
132  return CreateTH1(name, title, myxbins, opt);
133 }
134 
135 TH2* THistManager::CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt){
136  TString dirname(basename(name)), hname(histname(name));
137  THashList *parent(FindGroup(dirname));
138  if(!parent) parent = CreateHistoGroup(dirname);
139  if(parent->FindObject(hname)){
140  Fatal("THistManager::CreateTH2", "Object %s already exists in group %s", hname.Data(), dirname.Data());
141  return 0;
142  }
143  TH2* h = new TH2D(hname, title, nbinsx, xmin, xmax, nbinsy, ymin, ymax);
144  TString optionstring(opt);
145  optionstring.ToLower();
146  if(optionstring.Contains("s"))
147  h->Sumw2();
148  parent->Add(h);
149  return h;
150 }
151 
152 TH2* THistManager::CreateTH2(const char *name, const char *title, int nbinsx, const double *xbins, int nbinsy, const double *ybins, Option_t *opt){
153  TString dirname(basename(name)), hname(histname(name));
154  THashList *parent(FindGroup(dirname));
155  if(!parent) parent = CreateHistoGroup(dirname);
156  if(parent->FindObject(hname)){
157  Fatal("THistManager::CreateTH2", "Object %s already exists in group %s", hname.Data(), dirname.Data());
158  return 0;
159  }
160  TH2* h = new TH2D(hname, title, nbinsx, xbins, nbinsy, ybins);
161  TString optionstring(opt);
162  optionstring.ToLower();
163  if(optionstring.Contains("s"))
164  h->Sumw2();
165  parent->Add(h);
166  return h;
167 }
168 
169 TH2* THistManager::CreateTH2(const char *name, const char *title, const TArrayD &xbins, const TArrayD &ybins, Option_t *opt){
170  TString dirname(basename(name)), hname(histname(name));
171  THashList *parent(FindGroup(dirname));
172  if(!parent) parent = CreateHistoGroup(dirname);
173  if(parent->FindObject(hname)){
174  Fatal("THistManager::CreateTH2", "Object %s already exists in group %s", hname.Data(), dirname.Data());
175  return 0;
176  }
177  TH2* h = new TH2D(hname, title, xbins.GetSize() - 1, xbins.GetArray(), ybins.GetSize() - 1, ybins.GetArray());
178  TString optionstring(opt);
179  optionstring.ToLower();
180  if(optionstring.Contains("s"))
181  h->Sumw2();
182  parent->Add(h);
183  return h;
184 }
185 
186 TH2* THistManager::CreateTH2(const char *name, const char *title, const TBinning &xbins, const TBinning &ybins, Option_t *opt){
187  TArrayD myxbins, myybins;
188  try{
189  xbins.CreateBinEdges(myxbins);
190  } catch (std::exception &e) {
191  Fatal("THistManager::CreateTH2 (x-dir)", "Exception raised: %s", e.what());
192  }
193  try{
194  ybins.CreateBinEdges(myybins);
195  } catch (std::exception &e) {
196  Fatal("THistManager::CreateTH2 (y-dir)", "Exception raised: %s", e.what());
197  }
198 
199  return CreateTH2(name, title, myxbins, myybins, opt);
200 }
201 
202 TH3* THistManager::CreateTH3(const char* name, const char* title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, int nbinsz, double zmin, double zmax, Option_t *opt) {
203  TString dirname(basename(name)), hname(histname(name));
204  THashList *parent(FindGroup(dirname));
205  if(!parent) parent = CreateHistoGroup(dirname);
206  if(parent->FindObject(hname)){
207  Fatal("THistManager::CreateTH3", "Object %s already exists in group %s", hname.Data(), dirname.Data());
208  return 0;
209  }
210  TH3* h = new TH3D(hname, title, nbinsx, xmin, xmax, nbinsy, ymin, ymax, nbinsz, zmin, zmax);
211  TString optionstring(opt);
212  optionstring.ToLower();
213  if(optionstring.Contains("s"))
214  h->Sumw2();
215  parent->Add(h);
216  return h;
217 }
218 
219 TH3* THistManager::CreateTH3(const char* name, const char* title, int nbinsx, const double* xbins, int nbinsy, const double* ybins, int nbinsz, const double* zbins, Option_t *opt) {
220  TString dirname(basename(name)), hname(histname(name));
221  THashList *parent(FindGroup(dirname));
222  if(!parent) parent = CreateHistoGroup(dirname);
223  if(parent->FindObject(hname)){
224  Fatal("THistManager::CreateTH3", "Object %s already exists in group %s", hname.Data(), dirname.Data());
225  return 0;
226  }
227  TH3* h = new TH3D(hname, title, nbinsx, xbins, nbinsy, ybins, nbinsz, zbins);
228  TString optionstring(opt);
229  optionstring.ToLower();
230  if(optionstring.Contains("s"))
231  h->Sumw2();
232  parent->Add(h);
233  return h;
234 }
235 
236 TH3* THistManager::CreateTH3(const char* name, const char* title, const TArrayD& xbins, const TArrayD& ybins, const TArrayD& zbins, Option_t *opt) {
237  TString dirname(basename(name)), hname(histname(name));
238  THashList *parent(FindGroup(dirname));
239  if(!parent) parent = CreateHistoGroup(dirname);
240  if(parent->FindObject(hname)){
241  Fatal("THistManager::CreateTH3", "Object %s already exists in group %s", hname.Data(), dirname.Data());
242  return 0;
243  }
244  TH3* h = new TH3D(hname, title, xbins.GetSize()-1, xbins.GetArray(), ybins.GetSize()-1, ybins.GetArray(), zbins.GetSize()-1, zbins.GetArray());
245  TString optionstring(opt);
246  optionstring.ToLower();
247  if(optionstring.Contains("s"))
248  h->Sumw2();
249  parent->Add(h);
250  return h;
251 }
252 
253 TH3* THistManager::CreateTH3(const char *name, const char *title, const TBinning &xbins, const TBinning &ybins, const TBinning &zbins, Option_t *opt){
254  TArrayD myxbins, myybins, myzbins;
255  try{
256  xbins.CreateBinEdges(myxbins);
257  } catch (std::exception &e) {
258  Fatal("THistManager::CreateTH2 (x-dir)", "Exception raised: %s", e.what());
259  }
260  try{
261  ybins.CreateBinEdges(myybins);
262  } catch (std::exception &e) {
263  Fatal("THistManager::CreateTH2 (y-dir)", "Exception raised: %s", e.what());
264  }
265  try{
266  zbins.CreateBinEdges(myzbins);
267  } catch (std::exception &e) {
268  Fatal("THistManager::CreateTH2 (z-dir)", "Exception raised: %s", e.what());
269  }
270 
271  return CreateTH3(name, title, myxbins, myybins, myzbins);
272 }
273 
274 THnSparse* THistManager::CreateTHnSparse(const char *name, const char *title, int ndim, const int *nbins, const double *min, const double *max, Option_t *opt) {
275  TString dirname(basename(name)), hname(histname(name));
276  THashList *parent(FindGroup(dirname));
277  if(!parent) parent = CreateHistoGroup(dirname);
278  if(parent->FindObject(hname)){
279  Fatal("THistManager::CreateTHnSparse", "Object %s already exists in group %s", hname.Data(), dirname.Data());
280  return 0;
281  }
282  THnSparse* h = new THnSparseD(hname, title, ndim, nbins, min, max);
283  TString optionstring(opt);
284  optionstring.ToLower();
285  if(optionstring.Contains("s"))
286  h->Sumw2();
287  parent->Add(h);
288  return h;
289 }
290 
291 THnSparse* THistManager::CreateTHnSparse(const char *name, const char *title, int ndim, const TAxis **axes, Option_t *opt) {
292  TString dirname(basename(name)), hname(histname(name));
293  THashList *parent(FindGroup(dirname));
294  if(!parent) parent = CreateHistoGroup(dirname);
295  if(parent->FindObject(hname)){
296  Fatal("THistManager::CreateTHnSparse", "Object %s already exists in group %s", hname.Data(), dirname.Data());
297  return 0;
298  }
299  TArrayD xmin(ndim), xmax(ndim);
300  TArrayI nbins(ndim);
301  for(int idim = 0; idim < ndim; ++idim){
302  const TAxis &myaxis = *(axes[idim]);
303  nbins[idim] = myaxis.GetNbins();
304  xmin[idim] = myaxis.GetXmin();
305  xmax[idim] = myaxis.GetXmax();
306  }
307  THnSparseD *hsparse = new THnSparseD(hname, title, ndim, nbins.GetArray(), xmin.GetArray(), xmax.GetArray());
308  for(int id = 0; id < ndim; ++id)
309  *(hsparse->GetAxis(id)) = *(axes[id]);
310  TString optionstring(opt);
311  optionstring.ToLower();
312  if(optionstring.Contains("s"))
313  hsparse->Sumw2();
314  parent->Add(hsparse);
315  return hsparse;
316 }
317 
318 void THistManager::CreateTProfile(const char* name, const char* title, int nbinsX, double xmin, double xmax, Option_t *opt) {
319  TString dirname(basename(name)), hname(histname(name));
320  THashList *parent(FindGroup(dirname));
321  if(!parent) parent = CreateHistoGroup(dirname);
322  if(parent->FindObject(hname))
323  Fatal("THistManager::CreateTProfile", "Object %s already exists in group %s", hname.Data(), dirname.Data());
324  TProfile *hist = new TProfile(hname, title, nbinsX, xmin, xmax, opt);
325  parent->Add(hist);
326 }
327 
328 void THistManager::CreateTProfile(const char* name, const char* title, int nbinsX, const double* xbins, Option_t *opt) {
329  TString dirname(basename(name)), hname(histname(name));
330  THashList *parent(FindGroup(dirname));
331  if(!parent) parent = CreateHistoGroup(dirname);
332  if(parent->FindObject(hname))
333  Fatal("THistManager::CreateTHnSparse", "Object %s already exists in group %s", hname.Data(), dirname.Data());
334  TProfile *hist = new TProfile(hname, title, nbinsX, xbins, opt);
335  parent->Add(hist);
336 }
337 
338 void THistManager::CreateTProfile(const char* name, const char* title, const TArrayD& xbins, Option_t *opt){
339  TString dirname(basename(name)), hname(histname(name));
340  THashList *parent(FindGroup(dirname));
341  if(!parent) parent = CreateHistoGroup(dirname);
342  if(parent->FindObject(hname))
343  Fatal("THistManager::CreateTHnSparse", "Object %s already exists in group %s", hname.Data(), dirname.Data());
344  TProfile *hist = new TProfile(hname.Data(), title, xbins.GetSize()-1, xbins.GetArray(), opt);
345  parent->Add(hist);
346 }
347 
348 void THistManager::CreateTProfile(const char *name, const char *title, const TBinning &xbins, Option_t *opt){
349  TArrayD myxbins;
350  try{
351  xbins.CreateBinEdges(myxbins);
352  } catch (std::exception &e){
353  Fatal("THistManager::CreateProfile", "Exception raised: %s", e.what());
354  }
355  CreateTProfile(name, title, myxbins, opt);
356 }
357 
358 void THistManager::SetObject(TObject * const o, const char *group) {
359  THashList *parent(FindGroup(group));
360  if(!parent) CreateHistoGroup(group);
361  if(parent->FindObject(o->GetName())){
362  Fatal("THistManager::SetObject", "Parent %s does not exist", strcmp(group, "/") ? group : "");
363  return;
364  }
365  if(!(dynamic_cast<THnBase *>(o) || dynamic_cast<TH1 *>(o))){
366  Fatal("THistManager::SetObject", "Object %s is not of a histogram type",o->GetName());
367  return;
368  }
369  fHistos->Add(o);
370 }
371 
372 void THistManager::FillTH1(const char *name, double x, double weight, Option_t *opt) {
373  TString dirname(basename(name)), hname(histname(name));
374  THashList *parent(FindGroup(dirname));
375  if(!parent){
376  Fatal("THistManager::FillTH1", "Parnt group %s does not exist", dirname.Data());
377  return;
378  }
379  TH1 *hist = dynamic_cast<TH1 *>(parent->FindObject(hname));
380  if(!hist){
381  Fatal("THistManager::FillTH1", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
382  return;
383  }
384  TString optionstring(opt);
385  if(optionstring.Contains("w")){
386  // use bin width as weight
387  Int_t bin = hist->GetXaxis()->FindBin(x);
388  // check if not overflow or underflow bin
389  if(bin != 0 && bin != hist->GetXaxis()->GetNbins())
390  weight = 1./hist->GetXaxis()->GetBinWidth(bin);
391  }
392  hist->Fill(x, weight);
393 }
394 
395 void THistManager::FillTH1(const char *name, const char *label, double weight, Option_t *opt) {
396  TString dirname(basename(name)), hname(histname(name));
397  THashList *parent(FindGroup(dirname));
398  if(!parent){
399  Fatal("THistManager::FillTH1", "Parnt group %s does not exist", dirname.Data());
400  return;
401  }
402  TH1 *hist = dynamic_cast<TH1 *>(parent->FindObject(hname));
403  if(!hist){
404  Fatal("THistManager::FillTH1", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
405  return;
406  }
407  TString optionstring(opt);
408  if(optionstring.Contains("w")){
409  // use bin width as weight
410  // get bin for label
411  Int_t bin = hist->GetXaxis()->FindBin(label);
412  // check if not overflow or underflow bin
413  if(bin != 0 && bin != hist->GetXaxis()->GetNbins())
414  weight = 1./hist->GetXaxis()->GetBinWidth(bin);
415  }
416  hist->Fill(label, weight);
417 }
418 
419 void THistManager::FillTH2(const char *name, double x, double y, double weight, Option_t *opt) {
420  TString dirname(basename(name)), hname(histname(name));
421  THashList *parent(FindGroup(dirname));
422  if(!parent){
423  Fatal("THistManager::FillTH2", "Parent group %s does not exist", dirname.Data());
424  return;
425  }
426  TH2 *hist = dynamic_cast<TH2 *>(parent->FindObject(hname));
427  if(!hist){
428  Fatal("THistManager::FillTH2", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
429  return;
430  }
431  TString optstring(opt);
432  Double_t myweight = optstring.Contains("w") ? 1. : weight;
433  if(optstring.Contains("wx")){
434  Int_t binx = hist->GetXaxis()->FindBin(x);
435  if(binx != 0 && binx != hist->GetXaxis()->GetNbins()) myweight *= 1./hist->GetXaxis()->GetBinWidth(binx);
436  }
437  if(optstring.Contains("wy")){
438  Int_t biny = hist->GetYaxis()->FindBin(y);
439  if(biny != 0 && biny != hist->GetYaxis()->GetNbins()) myweight *= 1./hist->GetYaxis()->GetBinWidth(biny);
440  }
441  hist->Fill(x, y, myweight);
442 }
443 
444 void THistManager::FillTH2(const char *name, double *point, double weight, Option_t *opt) {
445  TString dirname(basename(name)), hname(histname(name));
446  THashList *parent(FindGroup(dirname));
447  if(!parent){
448  Fatal("THistManager::FillTH2", "Parent group %s does not exist", dirname.Data());
449  return;
450  }
451  TH2 *hist = dynamic_cast<TH2 *>(parent->FindObject(hname));
452  if(!hist){
453  Fatal("THistManager::FillTH2", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
454  return;
455  }
456  TString optstring(opt);
457  Double_t myweight = optstring.Contains("w") ? 1. : weight;
458  if(optstring.Contains("wx")){
459  Int_t binx = hist->GetXaxis()->FindBin(point[0]);
460  if(binx != 0 && binx != hist->GetXaxis()->GetNbins()) myweight *= 1./hist->GetXaxis()->GetBinWidth(binx);
461  }
462  if(optstring.Contains("wy")){
463  Int_t biny = hist->GetYaxis()->FindBin(point[1]);
464  if(biny != 0 && biny != hist->GetYaxis()->GetNbins()) myweight *= 1./hist->GetYaxis()->GetBinWidth(biny);
465  }
466  hist->Fill(point[0], point[1], weight);
467 }
468 
469 void THistManager::FillTH3(const char* name, double x, double y, double z, double weight, Option_t *opt) {
470  TString dirname(basename(name)), hname(histname(name));
471  THashList *parent(FindGroup(dirname));
472  if(!parent){
473  Fatal("THistManager::FillTH3", "Parent group %s does not exist", dirname.Data());
474  return;
475  }
476  TH3 *hist = dynamic_cast<TH3 *>(parent->FindObject(hname));
477  if(!hist){
478  Fatal("THistManager::FillTH3", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
479  return;
480  }
481  TString optstring(opt);
482  Double_t myweight = optstring.Contains("w") ? 1. : weight;
483  if(optstring.Contains("wx")){
484  Int_t binx = hist->GetXaxis()->FindBin(x);
485  if(binx != 0 && binx != hist->GetXaxis()->GetNbins()) myweight *= 1./hist->GetXaxis()->GetBinWidth(binx);
486  }
487  if(optstring.Contains("wy")){
488  Int_t biny = hist->GetYaxis()->FindBin(y);
489  if(biny != 0 && biny != hist->GetYaxis()->GetNbins()) myweight *= 1./hist->GetYaxis()->GetBinWidth(biny);
490  }
491  if(optstring.Contains("wz")){
492  Int_t binz = hist->GetZaxis()->FindBin(z);
493  if(binz != 0 && binz != hist->GetZaxis()->GetNbins()) myweight *= 1./hist->GetZaxis()->GetBinWidth(binz);
494  }
495  hist->Fill(x, y, z, weight);
496 }
497 
498 void THistManager::FillTH3(const char* name, const double* point, double weight, Option_t *opt) {
499  TString dirname(basename(name)), hname(histname(name));
500  THashList *parent(FindGroup(dirname));
501  if(!parent){
502  Fatal("THistManager::FillTH3", "Parent group %s does not exist", dirname.Data());
503  return;
504  }
505  TH3 *hist = dynamic_cast<TH3 *>(parent->FindObject(hname));
506  if(!hist){
507  Fatal("THistManager::FillTH3", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
508  return;
509  }
510  TString optstring(opt);
511  Double_t myweight = optstring.Contains("w") ? 1. : weight;
512  if(optstring.Contains("wx")){
513  Int_t binx = hist->GetXaxis()->FindBin(point[0]);
514  if(binx != 0 && binx != hist->GetXaxis()->GetNbins()) myweight *= 1./hist->GetXaxis()->GetBinWidth(binx);
515  }
516  if(optstring.Contains("wy")){
517  Int_t biny = hist->GetYaxis()->FindBin(point[1]);
518  if(biny != 0 && biny != hist->GetYaxis()->GetNbins()) myweight *= 1./hist->GetYaxis()->GetBinWidth(biny);
519  }
520  if(optstring.Contains("wz")){
521  Int_t binz = hist->GetZaxis()->FindBin(point[2]);
522  if(binz != 0 && binz != hist->GetZaxis()->GetNbins()) myweight *= 1./hist->GetZaxis()->GetBinWidth(binz);
523  }
524  hist->Fill(point[0], point[1], point[2], weight);
525 }
526 
527 void THistManager::FillTHnSparse(const char *name, const double *x, double weight, Option_t *opt) {
528  TString dirname(basename(name)), hname(histname(name));
529  THashList *parent(FindGroup(dirname));
530  if(!parent){
531  Fatal("THistManager::FillTHnSparse", "Parent group %s does not exist", dirname.Data());
532  return;
533  }
534  THnSparseD *hist = dynamic_cast<THnSparseD *>(parent->FindObject(hname));
535  if(!hist){
536  Fatal("THistManager::FillTHnSparse", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
537  return;
538  }
539  TString optstring(opt);
540  Double_t myweight = optstring.Contains("w") ? 1. : weight;
541  for(Int_t iaxis = 0; iaxis < hist->GetNdimensions(); iaxis++){
542  std::stringstream weighthandler;
543  weighthandler << "w" << iaxis;
544  if(optstring.Contains(weighthandler.str().c_str())){
545  Int_t bin = hist->GetAxis(iaxis)->FindBin(x[iaxis]);
546  if(bin != 0 && bin != hist->GetAxis(iaxis)->GetNbins()) myweight *= hist->GetAxis(iaxis)->GetBinWidth(bin);
547  }
548  }
549 
550  hist->Fill(x, weight);
551 }
552 
553 void THistManager::FillProfile(const char* name, double x, double y, double weight){
554  TString dirname(basename(name)), hname(histname(name));
555  THashList *parent(FindGroup(dirname));
556  if(!parent)
557  Fatal("THistManager::FillTProfile", "Parent group %s does not exist", dirname.Data());
558  TProfile *hist = dynamic_cast<TProfile *>(parent->FindObject(hname));
559  if(!hist)
560  Fatal("THistManager::FillTProfile", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
561  hist->Fill(x, y, weight);
562 }
563 
564 TObject *THistManager::FindObject(const char *name) const {
565  TString dirname(basename(name)), hname(histname(name));
566  THashList *parent(FindGroup(dirname));
567  if(!parent) return NULL;
568  return parent->FindObject(hname);
569 }
570 
572  TString dirname(basename(obj->GetName())), hname(histname(obj->GetName()));
573  THashList *parent(FindGroup(dirname));
574  if(!parent) return NULL;
575  return parent->FindObject(hname);
576 }
577 
578 THashList *THistManager::FindGroup(const char *dirname) const {
579  if(!strlen(dirname) || !strcmp(dirname, "/")) return fHistos;
580  // recursive find - avoids tokenizing filename
581  THashList *parentlist = FindGroup(basename(dirname));
582  if(parentlist) return static_cast<THashList *>(parentlist->FindObject(histname(dirname)));
583  return nullptr;
584 }
585 
587  int index = path.Last('/');
588  if(index < 0) return ""; // no directory structure
589  return TString(path(0, index));
590 }
591 
593  int index = path.Last('/');
594  if(index < 0) return path; // no directory structure
595  return TString(path(index+1, path.Length() - (index+1)));
596 }
597 
603 
605  fkArray(hmgr),
606  fCurrentPos(),
607  fNext(),
608  fDirection(dir)
609 {}
610 
612  fkArray(ref.fkArray),
613  fCurrentPos(ref.fCurrentPos),
614  fNext(ref.fNext),
615  fDirection(ref.fDirection)
616 {}
617 
619  if(this != &ref){
620  fkArray = ref.fkArray;
621  fCurrentPos = ref.fCurrentPos;
622  fNext = ref.fNext;
623  fDirection = ref.fDirection;
624  }
625  return *this;
626 }
627 
629  return fCurrentPos == other.fCurrentPos;
630 }
631 
633  if(fDirection == kTHMIforward)
634  fCurrentPos++;
635  else
636  fCurrentPos--;
637  return *this;
638 }
639 
641  iterator tmp(*this);
642  operator++();
643  return tmp;
644 }
645 
647  if(fDirection == kTHMIforward)
648  fCurrentPos--;
649  else
650  fCurrentPos++;
651  return *this;
652 };
653 
655  iterator tmp(*this);
656  operator--();
657  return tmp;
658 };
659 
661  if(fCurrentPos >=0 && fCurrentPos < fkArray->GetListOfHistograms()->GetEntries())
662  return fkArray->GetListOfHistograms()->At(fCurrentPos);
663  return NULL;
664 }
665 
666 
672 
673 namespace TestTHistManager {
674 
675  int THistManagerTestSuite::TestBuildSimpleHistograms(){
676  THistManager testmgr("testmgr");
677 
678  // Create histogram of each type
679  testmgr.CreateTH1("Test1D", "Test Histogram 1D", 1, 0., 1.);
680  testmgr.CreateTH2("Test2D", "Test Histogram 2D", 2, 0., 2., 10., 0., 10);
681  testmgr.CreateTH3("Test3D", "Test Histogram 3D", 3, 2, 6., 10., 0., 10., 50., 0., 50.);
682  int nbins[3] = {3, 3, 3}; double min[3] = {0., 0., 0}, max[3] = {6, 9, 12}; // dimensions for THnSparseTest
683  testmgr.CreateTHnSparse("TestNSparse", "Test Histogram NSparse", 3, nbins, min, max);
684  testmgr.CreateTProfile("TestProfile", "Test TProfile", 1, 0., 1);
685 
686  // Find histograms in the list (evaluate test)
687  // Tell user why test has failed
688  Bool_t found(true);
689  if(!dynamic_cast<TH1 *>(testmgr.GetListOfHistograms()->FindObject("Test1D"))){
690  std::cout << "Not found: Test1D" << std::endl;
691  found = false;
692  }
693  if(!dynamic_cast<TH2 *>(testmgr.GetListOfHistograms()->FindObject("Test2D"))){
694  std::cout << "Not found: Test2D" << std::endl;
695  found = false;
696  }
697  if(!dynamic_cast<TH3 *>(testmgr.GetListOfHistograms()->FindObject("Test3D"))){
698  std::cout << "Not found: Test3D" << std::endl;
699  found = false;
700  }
701  if(!dynamic_cast<THnSparse *>(testmgr.GetListOfHistograms()->FindObject("TestNSparse"))){
702  std::cout << "Not found: TestNSparse" << std::endl;
703  found = false;
704  }
705  if(!dynamic_cast<TProfile *>(testmgr.GetListOfHistograms()->FindObject("TestProfile"))){
706  std::cout << "Not found: TestProfile" << std::endl;
707  found = false;
708  }
709 
710  return found ? 0 : 1;
711  }
712 
713  int THistManagerTestSuite::TestBuildGroupedHistograms(){
714  THistManager testmgr("testmgr");
715 
716  // Creating 3 groups
717  testmgr.CreateTH1("Group1/Test1", "Test Histogram 1 in group 1", 1, 0., 1.);
718  testmgr.CreateTH1("Group1/Test2", "Test Histogram 2 in group 1", 1, 0., 1.);
719  testmgr.CreateTH2("Group2/Test1", "Test Histogram 1 in group 2", 1, 0., 1., 1, 0., 1.);
720  testmgr.CreateTH2("Group2/Test2", "Test Histogram 2 in group 2", 1, 0., 1., 1, 0., 1.);
721  testmgr.CreateTH3("Group3/Test1", "Test Histogram 1 in group 3", 1, 0., 1., 1, 0., 1., 1, 0., 1.);
722  testmgr.CreateTH3("Group3/Test2", "Test Histogram 2 in group 3", 1, 0., 1., 1, 0., 1., 1, 0., 1.);
723  // Creating histogram in a subgroup
724  testmgr.CreateTProfile("Group4/Subgroup1/Test1", "Test histogram for subgroup handling", 1, 0., 1);
725 
726  // Evalutate test
727  // Tell user why test has failed
728  Bool_t found(true);
729  THashList *currentdir(nullptr);
730  if(!(currentdir = dynamic_cast<THashList *>(testmgr.GetListOfHistograms()->FindObject("Group1")))){
731  std::cout << "Not found: Group1" << std::endl;
732  found = false;
733  } else {
734  if(!dynamic_cast<TH1 *>(currentdir->FindObject("Test1"))){
735  std::cout << "Not found in Group1: Test1" << std::endl;
736  found = false;
737  }
738  if(!dynamic_cast<TH1 *>(currentdir->FindObject("Test2"))){
739  std::cout << "Not found in Group1: Test2" << std::endl;
740  found = false;
741  }
742  }
743  if(!(currentdir = dynamic_cast<THashList *>(testmgr.GetListOfHistograms()->FindObject("Group2")))){
744  std::cout << "Not found: Group2" << std::endl;
745  found = false;
746  } else {
747  if(!dynamic_cast<TH2 *>(currentdir->FindObject("Test1"))){
748  std::cout << "Not found in Group2: Test1" << std::endl;
749  found = false;
750  }
751  if(!dynamic_cast<TH2 *>(currentdir->FindObject("Test2"))){
752  std::cout << "Not found in Group2: Test2" << std::endl;
753  found = false;
754  }
755  }
756  if(!(currentdir = dynamic_cast<THashList *>(testmgr.GetListOfHistograms()->FindObject("Group3")))){
757  std::cout << "Not found: Group3" << std::endl;
758  found = false;
759  } else {
760  if(!static_cast<TH3 *>(currentdir->FindObject("Test1"))){
761  std::cout << "Not found in Group3: Test1" << std::endl;
762  found = false;
763  }
764  if(!static_cast<TH3 *>(currentdir->FindObject("Test2"))){
765  std::cout << "Not found in Group3: Test2" << std::endl;
766  found = false;
767  }
768  }
769  if(!(currentdir = dynamic_cast<THashList *>(testmgr.GetListOfHistograms()->FindObject("Group4")))){
770  std::cout << "Not found: Group4" << std::endl;
771  found = false;
772  } else {
773  if(!(currentdir = dynamic_cast<THashList *>(currentdir->FindObject("Subgroup1")))){
774  std::cout << "Not found in Group4: Subgroup1" << std::endl;
775  found = false;
776  } else {
777  if(!dynamic_cast<TH1 *>(currentdir->FindObject("Test1"))){
778  std::cout << "Not found in Subgroup1: Test1" << std::endl;
779  found = false;
780  }
781  }
782  }
783  return found ? 0 : 1;
784  }
785 
786  int THistManagerTestSuite::TestFillSimpleHistograms(){
787  THistManager testmgr("testmgr");
788 
789  testmgr.CreateTH1("Test1", "Test fill 1D histogram", 1, 0., 1.);
790  testmgr.CreateTH2("Test2", "Test fill 2D histogram", 1, 0., 1., 1, 0., 1.);
791  testmgr.CreateTH3("Test3", "Test fill 3D histogram", 1, 0., 1., 1, 0., 1., 1, 0., 1.);
792  int nbins[4] = {1,1,1,1}; double min[4] = {0.,0.,0.,0.}, max[4] = {1.,1.,1.,1.};
793  testmgr.CreateTHnSparse("TestN", "Test Fill THnSparse", 4, nbins, min, max);
794  testmgr.CreateTProfile("TestProfile", "Test fill Profile histogram", 1, 0., 1.);
795 
796  double point[4] = {0.5, 0.5, 0.5, 0.5};
797  for(int i = 0; i < 100; i++){
798  testmgr.FillTH1("Test1", 0.5);
799  testmgr.FillTH2("Test2", 0.5, 0.5);
800  testmgr.FillTH3("Test3", 0.5, 0.5, 0.5);
801  testmgr.FillProfile("TestProfile", 0.5, 1.);
802  testmgr.FillTHnSparse("TestN", point);
803  }
804 
805  // Evalutate test
806  // tell user why test has failed
807  bool success(true);
808 
809  TH1 *test1 = dynamic_cast<TH1 *>(testmgr.GetListOfHistograms()->FindObject("Test1"));
810  if(test1){
811  if(TMath::Abs(test1->GetBinContent(1) - 100) > DBL_EPSILON){
812  std::cout << "Test1: Mismatch in values, expected 100, found " << test1->GetBinContent(1) << std::endl;
813  success = false;
814  }
815  } else {
816  std::cout << "Not found: Test1" << std::endl;
817  success = false;
818  }
819 
820  TH2 *test2 = dynamic_cast<TH2 *>(testmgr.GetListOfHistograms()->FindObject("Test2"));
821  if(test2){
822  if(TMath::Abs(test2->GetBinContent(1, 1) - 100) > DBL_EPSILON){
823  std::cout << "Test2: Mismatch in values, expected 100, found " << test2->GetBinContent(1,1) << std::endl;
824  success = false;
825  }
826  } else {
827  std::cout << "Not found: Test2" << std::endl;
828  success = false;
829  }
830 
831  TH3 *test3 = dynamic_cast<TH3 *>(testmgr.GetListOfHistograms()->FindObject("Test3"));
832  if(test3){
833  if(TMath::Abs(test3->GetBinContent(1, 1, 1) - 100) > DBL_EPSILON){
834  std::cout << "Test3: Mismatch in values, expected 100, found " << test3->GetBinContent(1,1,1) << std::endl;
835  success = false;
836  }
837  } else {
838  std::cout << "Not found: Test3" << std::endl;
839  success = false;
840  }
841 
842  THnSparse *testN = dynamic_cast<THnSparse *>(testmgr.GetListOfHistograms()->FindObject("TestN"));
843  if(testN){
844  int index[4] = {1,1,1,1};
845  if(TMath::Abs(testN->GetBinContent(index) - 100) > DBL_EPSILON){
846  std::cout << "TestN: Mismatch in values, expected 100, found " << testN->GetBinContent(index) << std::endl;
847  success = false;
848  }
849  } else {
850  std::cout << "Not found: TestN" << std::endl;
851  success = false;
852  }
853 
854  TProfile *testProfile = dynamic_cast<TProfile *>(testmgr.GetListOfHistograms()->FindObject("TestProfile"));
855  if(testProfile){
856  if(TMath::Abs(testProfile->GetBinContent(1) - 1) > DBL_EPSILON){
857  std::cout << "TestProfile: Mismatch in values, expected 1, found " << testProfile->GetBinContent(1) << std::endl;
858  success = false;
859  }
860  } else {
861  std::cout << "Not found: TestProfile" << std::endl;
862  success = false;
863  }
864 
865  return success ? 0 : 1;
866  }
867 
868 
869  int THistManagerTestSuite::TestFillGroupedHistograms(){
870  THistManager testmgr("testmgr");
871 
872  // Creating 3 groups, 1 with 1D and 1 with 2D histograms, and a third with a subgroup and a TProfile
873  testmgr.CreateTH1("Group1/Test1", "Test 1 Group 1D", 1, 0., 1.);
874  testmgr.CreateTH1("Group1/Test2", "Test 2 Group 1D", 1, 0., 1.);
875  testmgr.CreateTH2("Group2/Test1", "Test 1 Group 2D", 1, 0., 1., 1, 0., 1.);
876  testmgr.CreateTH2("Group2/Test2", "Test 2 Group 2D", 1, 0., 1., 1, 0., 1.);
877  testmgr.CreateTProfile("Group3/Subgroup1/Test1", "Test 1 with subgroup", 1, 0., 1.);
878 
879  for(int i = 0; i < 100; i++){
880  testmgr.FillTH1("Group1/Test1", 0.5);
881  testmgr.FillTH1("Group1/Test2", 0.5);
882  testmgr.FillTH2("Group2/Test1", 0.5, 0.5);
883  testmgr.FillTH2("Group2/Test2", 0.5, 0.5);
884  testmgr.FillProfile("Group3/Subgroup1/Test1", 0.5, 1);
885  }
886 
887  // Evaluate test
888  bool success(true);
889 
890  THashList *currentdir = dynamic_cast<THashList *>(testmgr.GetListOfHistograms()->FindObject("Group1"));
891  if(currentdir){
892  TH1 *test1 = dynamic_cast<TH1 *>(currentdir->FindObject("Test1"));
893  if(test1){
894  if(TMath::Abs(test1->GetBinContent(1) - 100) > DBL_EPSILON){
895  std::cout << "Group1/Test1: Value mismatch: expected 100, found " << test1->GetBinContent(1) << std::endl;
896  success = false;
897  }
898  } else {
899  std::cout << "Not found in Group1: Test1" << std::endl;
900  success = false;
901  }
902  test1 = dynamic_cast<TH1 *>(currentdir->FindObject("Test2"));
903  if(test1){
904  if(TMath::Abs(test1->GetBinContent(1) - 100) > DBL_EPSILON){
905  std::cout << "Group1/Test2: Value mismatch: expected 100, found " << test1->GetBinContent(1) << std::endl;
906  success = false;
907  }
908  } else {
909  std::cout << "Not found in Group1: Test2" << std::endl;
910  success = false;
911  }
912  } else {
913  std::cout << "Not found: Group1" << std::endl;
914  success = false;
915  }
916 
917  currentdir = dynamic_cast<THashList *>(testmgr.GetListOfHistograms()->FindObject("Group2"));
918  if(currentdir){
919  TH2 *test2 = dynamic_cast<TH2 *>(currentdir->FindObject("Test1"));
920  if(test2){
921  if(TMath::Abs(test2->GetBinContent(1,1) - 100) > DBL_EPSILON){
922  std::cout << "Group2/Test1: Value mismatch: expected 100, found " << test2->GetBinContent(1,1) << std::endl;
923  success = false;
924  }
925  } else {
926  std::cout << "Not found in Group2: Test1" << std::endl;
927  success = false;
928  }
929  test2 = dynamic_cast<TH2 *>(currentdir->FindObject("Test2"));
930  if(test2){
931  if(TMath::Abs(test2->GetBinContent(1,1) - 100) > DBL_EPSILON){
932  std::cout << "Group2/Test2: Value mismatch: expected 100, found " << test2->GetBinContent(1,1) << std::endl;
933  success = false;
934  }
935  } else {
936  std::cout << "Not found in Group2: Test2" << std::endl;
937  success = false;
938  }
939  } else {
940  std::cout << "Not found: Group2" << std::endl;
941  success = false;
942  }
943 
944  currentdir = dynamic_cast<THashList *>(testmgr.GetListOfHistograms()->FindObject("Group3"));
945  if(currentdir){
946  currentdir = dynamic_cast<THashList *>(currentdir->FindObject("Subgroup1"));
947  if(currentdir){
948  TProfile *testprofile = dynamic_cast<TProfile *>(currentdir->FindObject("Test1"));
949  if(testprofile){
950  if(TMath::Abs(testprofile->GetBinContent(1) - 1) > DBL_EPSILON){
951  std::cout << "Group3/Subgroup1/Test1: Value mismatch: expected 1, found " << testprofile->GetBinContent(1) << std::endl;
952  success = false;
953  }
954  } else {
955  std::cout << "Not found in Group3/Subgroup1: Test1" << std::endl;
956  success = false;
957  }
958  } else {
959  std::cout << "Not found in Group3: Subgroup1" << std::endl;
960  success = false;
961  }
962  } else {
963  std::cout << "Not found: Group3" << std::endl;
964  success = false;
965  }
966  return success ? 0 : 1;
967  }
968 
969  int TestRunAll(){
970  int testresult(0);
971  THistManagerTestSuite testsuite;
972 
973  std::cout << "Running test: Build simple" << std::endl;
974  testresult += testsuite.TestBuildSimpleHistograms();
975  std::cout << "Result after test: " << testresult << std::endl;
976 
977  std::cout << "Running test: Build grouped" << std::endl;
978  testresult += testsuite.TestBuildGroupedHistograms();
979  std::cout << "Result after test: " << testresult << std::endl;
980 
981  std::cout << "Running test: Fill Simple" << std::endl;
982  testresult += testsuite.TestFillSimpleHistograms();
983  std::cout << "Result after test: " << testresult << std::endl;
984 
985  std::cout << "Running test: Fill Grouped" << std::endl;
986  testresult += testsuite.TestFillGroupedHistograms();
987  std::cout << "Result after test: " << testresult << std::endl;
988 
989  return testresult;
990  }
991 
993  THistManagerTestSuite testsuite;
994  return testsuite.TestBuildSimpleHistograms();
995  }
996 
998  THistManagerTestSuite testsuite;
999  return testsuite.TestBuildGroupedHistograms();
1000  }
1001 
1003  THistManagerTestSuite testsuite;
1004  return testsuite.TestFillSimpleHistograms();
1005  }
1006 
1008  THistManagerTestSuite testsuite;
1009  return testsuite.TestFillGroupedHistograms();
1010  }
1011 }
iterator & operator=(const iterator &rhs)
THashList * CreateHistoGroup(const char *groupname)
const Double_t ymax
Definition: AddTaskCFDStar.C:7
double Double_t
Definition: External.C:58
const char * title
Definition: MakeQAPdf.C:26
Definition: External.C:244
bool fIsOwner
Set the ownership.
Definition: THistManager.h:549
Int_t nbinsy
THashList * FindGroup(const char *dirname) const
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
void FillTH3(const char *hname, double x, double y, double z, double weight=1., Option_t *opt="")
virtual void CreateBinEdges(TArrayD &binedges) const =0
const THistManager * fkArray
Definition: THistManager.h:156
void SetObject(TObject *const o, const char *group="/")
Interface for binnings used by the histogram handler.
Definition: TBinning.h:21
THMIDirection_t fDirection
Definition: THistManager.h:159
Bool_t operator!=(const iterator &aIter) const
void FillTHnSparse(const char *name, const double *x, double weight=1., Option_t *opt="")
TH2 * CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
TObject * FindObject(const char *name) const
int Int_t
Definition: External.C:63
void CreateTProfile(const char *name, const char *title, int nbinsX, double xmin, double xmax, Option_t *opt="")
THashList * GetListOfHistograms() const
Definition: THistManager.h:504
Definition: External.C:252
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
const Double_t zmin
Definition: External.C:228
Definition: External.C:212
void FillProfile(const char *name, double x, double y, double weight=1.)
THMIDirection_t
Direction for the iteration.
Definition: THistManager.h:72
stl-iterator for the histogram manager
Definition: THistManager.h:59
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Int_t nbinsx
TString basename(const TString &path) const
Collection of tests for the THistManager.
Definition: THistManager.h:592
TString histname(const TString &path) const
Definition: External.C:220
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
TObject * operator*() const
const Double_t ymin
Definition: AddTaskCFDStar.C:6
Container class for histograms.
Definition: THistManager.h:43
const char Option_t
Definition: External.C:48
const Double_t zmax
const Int_t nbins
bool Bool_t
Definition: External.C:53
THnSparse * CreateTHnSparse(const char *name, const char *title, int ndim, const int *nbins, const double *min, const double *max, Option_t *opt="")
THashList * fHistos
List of histograms.
Definition: THistManager.h:548
Definition: External.C:196
TH3 * CreateTH3(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, int nbinsz, double zmin, double zmax, Option_t *opt="")