AliPhysics  vAN-20150425 (2dcf1b0)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliEMCalHistoContainer.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2014, 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 <cstring>
16 #include <vector>
17 #include <TArrayD.h>
18 #include <TAxis.h>
19 #include <TH1.h>
20 #include <TH2.h>
21 #include <TH3.h>
22 #include <THnSparse.h>
23 #include <THashList.h>
24 #include <TObjArray.h>
25 #include <TObjString.h>
26 #include "AliEMCalHistoContainer.h"
27 #include <TString.h>
28 
29 #include "AliLog.h"
30 
34 
35 namespace EMCalTriggerPtAnalysis {
36 
40 AliEMCalHistoContainer::AliEMCalHistoContainer():
41  TNamed(),
42  fHistos(NULL),
43  fIsOwner(true)
44 {
45 }
46 
53  TNamed(name, Form("Histogram container %s", name)),
54  fHistos(NULL),
55  fIsOwner(true)
56 {
57  fHistos = new THashList();
58  fHistos->SetName(Form("histos%s", name));
59  fHistos->SetOwner();
60 }
61 
66  if(fHistos && fIsOwner) delete fHistos;
67 }
68 
77 void AliEMCalHistoContainer::CreateHistoGroup(const char *groupname, const char *parent) throw(HistoContainerContentException) {
78  THashList *parentgroup = FindGroup(parent);
80  THashList *childgroup = new THashList();
81  childgroup->SetName(groupname);
82  parentgroup->Add(childgroup);
83 }
84 
96 void AliEMCalHistoContainer::CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt) throw(HistoContainerContentException){
97  TString dirname(basename(name)), hname(histname(name));
98  THashList *parent(FindGroup(dirname.Data()));
99  if(!parent)
101  if(parent->FindObject(hname.Data()))
103  TH1 *hist = new TH1D(hname.Data(), title, nbins, xmin, xmax);
104  TString optionstring(opt);
105  optionstring.ToLower();
106  if(optionstring.Contains("s"))
107  hist->Sumw2();
108  parent->Add(hist);
109 }
110 
121 void AliEMCalHistoContainer::CreateTH1(const char *name, const char *title, int nbins, const double *xbins, Option_t *opt) throw(HistoContainerContentException){
122  TString dirname(basename(name)), hname(histname(name));
123  THashList *parent(FindGroup(dirname));
124  if(!parent)
126  if(parent->FindObject(hname.Data()))
128  TH1 *hist = new TH1D(hname.Data(), title, nbins, xbins);
129  TString optionstring(opt);
130  optionstring.ToLower();
131  if(optionstring.Contains("s"))
132  hist->Sumw2();
133  parent->Add(hist);
134 }
135 
145 void AliEMCalHistoContainer::CreateTH1(const char *name, const char *title, const TArrayD &xbins, Option_t *opt) throw(HistoContainerContentException){
146  TString dirname(basename(name)), hname(histname(name));
147  THashList *parent(FindGroup(dirname));
148  if(!parent)
150  if(parent->FindObject(hname.Data()))
152  TH1 *hist = new TH1D(hname.Data(), title, xbins.GetSize()-1, xbins.GetArray());
153  TString optionstring(opt);
154  optionstring.ToLower();
155  if(optionstring.Contains("s"))
156  hist->Sumw2();
157  parent->Add(hist);
158 }
159 
174 void AliEMCalHistoContainer::CreateTH2(const char *name, const char *title,
175  int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt) throw(HistoContainerContentException){
176  TString dirname(basename(name)), hname(histname(name));
177  THashList *parent(FindGroup(dirname.Data()));
178  if(!parent)
180  if(parent->FindObject(hname.Data()))
182  TH2 *hist = new TH2D(hname.Data(), title, nbinsx, xmin, xmax, nbinsy, ymin, ymax);
183  TString optionstring(opt);
184  optionstring.ToLower();
185  if(optionstring.Contains("s"))
186  hist->Sumw2();
187  parent->Add(hist);
188 }
189 
202 void AliEMCalHistoContainer::CreateTH2(const char *name, const char *title,
203  int nbinsx, const double *xbins, int nbinsy, const double *ybins, Option_t *opt) throw(HistoContainerContentException){
204  TString dirname(basename(name)), hname(histname(name));
205  THashList *parent(FindGroup(dirname.Data()));
206  if(!parent)
208  if(parent->FindObject(hname.Data()))
210  TH2 *hist = new TH2D(hname.Data(), title, nbinsx, xbins, nbinsy, ybins);
211  TString optionstring(opt);
212  optionstring.ToLower();
213  if(optionstring.Contains("s"))
214  hist->Sumw2();
215  parent->Add(hist);
216 }
217 
228 void AliEMCalHistoContainer::CreateTH2(const char *name, const char *title, const TArrayD &xbins, const TArrayD &ybins, Option_t *opt) throw(HistoContainerContentException){
229  TString dirname(basename(name)), hname(histname(name));
230  THashList *parent(FindGroup(dirname.Data()));
231  if(!parent)
233  if(parent->FindObject(hname.Data()))
235  TH2 *hist = new TH2D(hname.Data(), title, xbins.GetSize() - 1, xbins.GetArray(), ybins.GetSize() - 1, ybins.GetArray());
236  TString optionstring(opt);
237  optionstring.ToLower();
238  if(optionstring.Contains("s"))
239  hist->Sumw2();
240  parent->Add(hist);
241 }
242 
258 void AliEMCalHistoContainer::CreateTH3(const char* name, const char* title, int nbinsx, double xmin, double xmax,
259  int nbinsy, double ymin, double ymax, int nbinsz, double zmin, double zmax, Option_t *opt) throw (HistoContainerContentException) {
260  TString dirname(basename(name)), hname(histname(name));
261  THashList *parent(FindGroup(dirname.Data()));
262  if(!parent)
264  if(parent->FindObject(hname.Data()))
266  TH3 *hist = new TH3D(hname.Data(), title, nbinsx, xmin, xmax, nbinsy, ymin, ymax, nbinsz, zmin, zmax);
267  TString optionstring(opt);
268  optionstring.ToLower();
269  if(optionstring.Contains("s"))
270  hist->Sumw2();
271  parent->Add(hist);
272 }
273 
288 void AliEMCalHistoContainer::CreateTH3(const char* name, const char* title, int nbinsx, const double* xbins,
289  int nbinsy, const double* ybins, int nbinsz, const double* zbins, Option_t *opt) throw (HistoContainerContentException) {
290  TString dirname(basename(name)), hname(histname(name));
291  THashList *parent(FindGroup(dirname.Data()));
292  if(!parent)
294  if(parent->FindObject(hname.Data()))
296  TH3 *hist = new TH3D(hname.Data(), title, nbinsx, xbins, nbinsy, ybins, nbinsz, zbins);
297  TString optionstring(opt);
298  optionstring.ToLower();
299  if(optionstring.Contains("s"))
300  hist->Sumw2();
301  parent->Add(hist);
302 }
303 
315 void AliEMCalHistoContainer::CreateTH3(const char* name, const char* title, const TArrayD& xbins, const TArrayD& ybins,
316  const TArrayD& zbins, Option_t *opt) throw (HistoContainerContentException) {
317  TString dirname(basename(name)), hname(histname(name));
318  THashList *parent(FindGroup(dirname.Data()));
319  if(!parent)
321  if(parent->FindObject(hname.Data()))
323  TH3 *hist = new TH3D(hname.Data(), title, xbins.GetSize()-1, xbins.GetArray(), ybins.GetSize()-1, ybins.GetArray(), zbins.GetSize()-1, zbins.GetArray());
324  TString optionstring(opt);
325  optionstring.ToLower();
326  if(optionstring.Contains("s"))
327  hist->Sumw2();
328  parent->Add(hist);
329 }
330 
343 void AliEMCalHistoContainer::CreateTHnSparse(const char *name, const char *title,
344  int ndim, const int *nbins, const double *min, const double *max, Option_t *opt) throw(HistoContainerContentException){
345  TString dirname(basename(name)), hname(histname(name));
346  THashList *parent(FindGroup(dirname.Data()));
347  if(!parent)
349  if(parent->FindObject(hname.Data()))
351  THnSparseD *hist = new THnSparseD(hname.Data(), title, ndim, nbins, min, max);
352  TString optionstring(opt);
353  optionstring.ToLower();
354  if(optionstring.Contains("s"))
355  hist->Sumw2();
356  parent->Add(hist);
357 }
358 
369 void AliEMCalHistoContainer::CreateTHnSparse(const char *name, const char *title, int ndim, const TAxis **axes, Option_t *opt) throw(HistoContainerContentException){
370  TString dirname(basename(name)), hname(histname(name));
371  THashList *parent(FindGroup(dirname.Data()));
372  if(!parent)
374  if(parent->FindObject(hname))
376  TArrayD xmin(ndim), xmax(ndim);
377  TArrayI nbins(ndim);
378  for(int idim = 0; idim < ndim; ++idim){
379  const TAxis &myaxis = *(axes[idim]);
380  nbins[idim] = myaxis.GetNbins();
381  xmin[idim] = myaxis.GetXmin();
382  xmax[idim] = myaxis.GetXmax();
383  }
384  THnSparseD *hsparse = new THnSparseD(hname.Data(), title, ndim, nbins.GetArray(), xmin.GetArray(), xmax.GetArray());
385  for(int id = 0; id < ndim; ++id)
386  *(hsparse->GetAxis(id)) = *(axes[id]);
387  TString optionstring(opt);
388  optionstring.ToLower();
389  if(optionstring.Contains("s"))
390  hsparse->Sumw2();
391  parent->Add(hsparse);
392 }
393 
399 void AliEMCalHistoContainer::SetObject(TObject * const o, const char *group) throw(HistoContainerContentException){
400  THashList *parent(FindGroup(group));
401  if(!parent)
402  throw HistoContainerContentException(NULL, strcmp(group, "/") ? group : "", HistoContainerContentException::kGroupException);
403  if(parent->FindObject(o->GetName()))
404  throw HistoContainerContentException(o->GetName(), strcmp(group, "/") ? group : "", HistoContainerContentException::kHistDuplicationException);
405  if(!(dynamic_cast<THnBase *>(o) || dynamic_cast<TH1 *>(o)))
406  throw HistoContainerContentException(o->GetName(), strcmp(group, "/") ? group : "", HistoContainerContentException::kTypeException);
407  fHistos->Add(o);
408 }
409 
419 void AliEMCalHistoContainer::FillTH1(const char *name, double x, double weight) throw(HistoContainerContentException){
420  TString dirname(basename(name)), hname(histname(name));
421  THashList *parent(FindGroup(dirname.Data()));
422  if(!parent)
424  TH1 *hist = dynamic_cast<TH1 *>(parent->FindObject(hname.Data()));
425  if(!hist)
427  hist->Fill(x, weight);
428 }
429 
440 void AliEMCalHistoContainer::FillTH2(const char *name, double x, double y, double weight) throw(HistoContainerContentException){
441  TString dirname(basename(name)), hname(histname(name));
442  THashList *parent(FindGroup(dirname.Data()));
443  if(!parent)
445  TH2 *hist = dynamic_cast<TH2 *>(parent->FindObject(hname.Data()));
446  if(!hist)
448  hist->Fill(x, y, weight);
449 }
450 
460 void AliEMCalHistoContainer::FillTH2(const char *name, double *point, double weight) throw(HistoContainerContentException){
461  TString dirname(basename(name)), hname(histname(name));
462  THashList *parent(FindGroup(dirname.Data()));
463  if(!parent)
465  TH2 *hist = dynamic_cast<TH2 *>(parent->FindObject(hname.Data()));
466  if(!hist)
468  hist->Fill(point[0], point[1], weight);
469 }
470 
482 void AliEMCalHistoContainer::FillTH3(const char* name, double x, double y, double z, double weight) throw (HistoContainerContentException) {
483  TString dirname(basename(name)), hname(histname(name));
484  THashList *parent(FindGroup(dirname.Data()));
485  if(!parent)
487  TH3 *hist = dynamic_cast<TH3 *>(parent->FindObject(hname.Data()));
488  if(!hist)
490  hist->Fill(x, y, z, weight);
491 }
501 void AliEMCalHistoContainer::FillTH3(const char* name, const double* point, double weight) throw (HistoContainerContentException) {
502  TString dirname(basename(name)), hname(histname(name));
503  THashList *parent(FindGroup(dirname.Data()));
504  if(!parent)
506  TH3 *hist = dynamic_cast<TH3 *>(parent->FindObject(hname.Data()));
507  if(!hist)
509  hist->Fill(point[0], point[1], point[2], weight);
510 }
511 
521 void AliEMCalHistoContainer::FillTHnSparse(const char *name, const double *x, double weight) throw(HistoContainerContentException){
522  TString dirname(basename(name)), hname(histname(name));
523  THashList *parent(FindGroup(dirname.Data()));
524  if(!parent)
526  THnSparseD *hist = dynamic_cast<THnSparseD *>(parent->FindObject(hname.Data()));
527  if(!hist)
529  hist->Fill(x, weight);
530 }
531 
539 TObject *AliEMCalHistoContainer::FindObject(const char *name) const {
540  TString dirname(basename(name)), hname(histname(name));
541  THashList *parent(FindGroup(dirname.Data()));
542  if(!parent) return NULL;
543  return parent->FindObject(hname);
544 }
545 
553 TObject* AliEMCalHistoContainer::FindObject(const TObject* obj) const {
554  TString dirname(basename(obj->GetName())), hname(histname(obj->GetName()));
555  THashList *parent(FindGroup(dirname.Data()));
556  if(!parent) return NULL;
557  return parent->FindObject(hname);
558 }
559 
566 THashList *AliEMCalHistoContainer::FindGroup(const char *dirname) const {
567  if(!strlen(dirname) || !strcmp(dirname, "/")) return fHistos;
568  std::vector<std::string> tokens;
569  TokenizeFilename(dirname, "/", tokens);
570  THashList *currentdir(fHistos);
571  for(std::vector<std::string>::iterator it = tokens.begin(); it != tokens.end(); ++it){
572  currentdir = dynamic_cast<THashList *>(currentdir->FindObject(it->c_str()));
573  if(!currentdir) break;
574  }
575  return currentdir;
576 }
577 
585 void AliEMCalHistoContainer::TokenizeFilename(const char *name, const char *delim, std::vector<std::string> &listoftokens) const {
586  TString s(name);
587  TObjArray *arr = s.Tokenize(delim);
588  TObjString *ostr(NULL);
589  TIter toks(arr);
590  while((ostr = dynamic_cast<TObjString *>(toks()))){
591  listoftokens.push_back(std::string(ostr->String().Data()));
592  }
593  delete arr;
594 }
595 
602 const char *AliEMCalHistoContainer::basename(const char *path) const {
603  TString s(path);
604  int index = s.Last('/');
605  if(index < 0) return ""; // no directory structure
606  return TString(s(0, index)).Data();
607 }
608 
615 const char *AliEMCalHistoContainer::histname(const char *path) const {
616  TString s(path);
617  int index = s.Last('/');
618  if(index < 0) return path; // no directory structure
619  return TString(s(index+1, s.Length() - (index+1))).Data();
620 }
621 
622 }
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
void FillTH2(const char *hname, double x, double y, double weight=1.)
void 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="")
void SetObject(TObject *const o, const char *group="/")
void CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Container class for histograms for the high- charged particle analysis.
void FillTHnSparse(const char *name, const double *x, double weight=1.)
void FillTH3(const char *hname, double x, double y, double z, double weight=1.)
Exception thrown by the histogram container in case of problems.
virtual TObject * FindObject(const char *name) const
void CreateTHnSparse(const char *name, const char *title, int ndim, const int *nbins, const double *min, const double *max, Option_t *opt="")
Declarartion of class AliEMCalHistoContainer.
void CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
const char * basename(const char *path) const
void TokenizeFilename(const char *name, const char *delim, std::vector< std::string > &listoftokens) const
const char * histname(const char *path) const
THashList * FindGroup(const char *dirname) const
void FillTH1(const char *hname, double x, double weight=1.)
void CreateHistoGroup(const char *groupname, const char *parent="/")