1 /**************************************************************************
2 * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
17 * Container class for histogram objects. Currenly can handle
22 * Histograms can be stored in groups. For this the parent group is
23 * included inside the histogram name, i.e. /base/inheriting/histogram.
24 * In case just the histogram name is given, it is assumed that the
25 * histogram is stored at the top level.
27 * Author: Markus Fasel
38 #include <THnSparse.h>
39 #include <THashList.h>
40 #include <TObjArray.h>
41 #include <TObjString.h>
44 #include "THistManager.h"
46 ClassImp(THistManager)
48 //______________________________________________________________________________
49 THistManager::THistManager():
55 * Default constructor, only initialising pointers with 0
59 //______________________________________________________________________________
60 THistManager::THistManager(const char *name):
61 TNamed(name, Form("Histogram container %s", name)),
66 * Main constructor, creating also a list for the histograms
68 * @param name: Name of the object (list named accordingly)
70 fHistos = new THashList();
71 fHistos->SetName(Form("histos%s", name));
75 //______________________________________________________________________________
76 THistManager::~THistManager(){
78 * Destructor, deletes the list of histograms if it is the owner
80 if(fHistos && fIsOwner) delete fHistos;
83 //______________________________________________________________________________
84 void THistManager::CreateHistoGroup(const char *groupname, const char *parent) {
86 * Create a new group of histograms within a parent group. Groups are represented as list. The default parent is
87 * always the top list. List name structure accouding to unix paths (i.e. top list /, hirarchies separated by /).
89 * @param groupname: Name of the new group
90 * @param parent (@default "/"): Name of the parent group
91 * @throw HistManagerException
93 THashList *parentgroup = FindGroup(parent);
95 Fatal("THistManager::CreateHistoGroup", "Parent group %s does not exist", parentgroup->GetName());
97 THashList *childgroup = new THashList();
98 childgroup->SetName(groupname);
99 parentgroup->Add(childgroup);
102 //______________________________________________________________________________
103 void THistManager::CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax){
105 * Create a new TH1 within the container. The histogram name also contains the parent group(s) according to the common
108 * @param name: Name of the histogram
109 * @param title: Title of the histogram
110 * @param nbins: number of bins
111 * @param xmin: min. value of the range
112 * @param xmax: max. value of the range
113 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
115 TString dirname(basename(name)), hname(histname(name));
116 THashList *parent(FindGroup(dirname.Data()));
118 Fatal("THistManager::CreateTH1", "Parent %s does not exist", dirname.Data());
119 if(parent->FindObject(hname.Data()))
120 Fatal("THistManager::CreateTH1", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
121 parent->Add(new TH1D(hname.Data(), title, nbins, xmin, xmax));
124 //______________________________________________________________________________
125 void THistManager::CreateTH1(const char *name, const char *title, int nbins, const double *xbins){
127 * Create a new TH1 within the container. The histogram name also contains the parent group(s) according to the common
130 * @param name: Name of the histogram
131 * @param title: Title of the histogram
132 * @param nbins: number of bins
133 * @param xbins: array of bin limits
134 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
136 TString dirname(basename(name)), hname(histname(name));
137 THashList *parent(FindGroup(dirname));
139 Fatal("THistManager::CreateTH1", "Parent %s does not exist", dirname.Data());
140 if(parent->FindObject(hname.Data()))
141 Fatal("THistManager::CreateTH1", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
142 parent->Add(new TH1D(hname.Data(), title, nbins, xbins));
145 //______________________________________________________________________________
146 void THistManager::CreateTH1(const char *name, const char *title, const TArrayD &xbins){
148 * Create a new TH1 within the container. The histogram name also contains the parent group(s) according to the common
151 * @param name: Name of the histogram
152 * @param title: Title of the histogram
153 * @param xbins: array of bin limits (contains also number of bins)
154 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
156 TString dirname(basename(name)), hname(histname(name));
157 THashList *parent(FindGroup(dirname));
159 Fatal("THistManager::CreateTH1", "Parent %s does not exist", dirname.Data());
160 if(parent->FindObject(hname.Data()))
161 Fatal("THistManager::CreateTH1", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
162 parent->Add(new TH1D(hname.Data(), title, xbins.GetSize()-1, xbins.GetArray()));
165 //______________________________________________________________________________
166 void THistManager::CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax){
168 * Create a new TH2 within the container. The histogram name also contains the parent group(s) according to the common
171 * @param name: Name of the histogram
172 * @param title: Title of the histogram
173 * @param nbinsx: number of bins in x-direction
174 * @param xmin: min. value of the range in x-direction
175 * @param xmax: max. value of the range in x-direction
176 * @param nbinsy: number of bins in y-direction
177 * @param ymin: min. value of the range in y-direction
178 * @param ymax: max. value of the range in y-direction
179 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
181 TString dirname(basename(name)), hname(histname(name));
182 THashList *parent(FindGroup(dirname.Data()));
184 Fatal("THistManager::CreateTH2", "Parent %s does not exist", dirname.Data());
185 if(parent->FindObject(hname.Data()))
186 Fatal("THistManager::CreateTH2", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
187 parent->Add(new TH2D(hname.Data(), title, nbinsx, xmin, xmax, nbinsy, ymin, ymax));
190 //______________________________________________________________________________
191 void THistManager::CreateTH2(const char *name, const char *title, int nbinsx, const double *xbins, int nbinsy, const double *ybins){
193 * Create a new TH2 within the container. The histogram name also contains the parent group(s) according to the common
196 * @param name: Name of the histogram
197 * @param title: Title of the histogram
198 * @param nbinsx: number of bins in x-direction
199 * @param xbins: array of bin limits in x-direction
200 * @param nbinsy: number of bins in y-direction
201 * @param ybins: array of bin limits in y-direction
202 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
204 TString dirname(basename(name)), hname(histname(name));
205 THashList *parent(FindGroup(dirname.Data()));
207 Fatal("THistManager::CreateTH2", "Parent %s does not exist", dirname.Data());
208 if(parent->FindObject(hname.Data()))
209 Fatal("THistManager::CreateTH2", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
210 parent->Add(new TH2D(hname.Data(), title, nbinsx, xbins, nbinsy, ybins));
213 //______________________________________________________________________________
214 void THistManager::CreateTH2(const char *name, const char *title, const TArrayD &xbins, const TArrayD &ybins){
216 * Create a new TH2 within the container. The histogram name also contains the parent group(s) according to the common
219 * @param name: Name of the histogram
220 * @param title: Title of the histogram
221 * @param xbins: array of bin limits in x-direction (contains also the number of bins)
222 * @param ybins: array of bin limits in y-direction (contains also the number of bins)
223 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
225 TString dirname(basename(name)), hname(histname(name));
226 THashList *parent(FindGroup(dirname.Data()));
228 Fatal("THistManager::CreateTH2", "Parent %s does not exist", dirname.Data());
229 if(parent->FindObject(hname.Data()))
230 Fatal("THistManager::CreateTH2", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
231 parent->Add(new TH2D(hname.Data(), title, xbins.GetSize() - 1, xbins.GetArray(), ybins.GetSize() - 1, ybins.GetArray()));
234 //______________________________________________________________________________
235 void 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) {
237 * Create a new TH3 within the container. The histogram name also contains the parent group(s) according to the common
240 * @param nbinsx: number of bins in x-direction
241 * @param xmin: min. value of the range in x-direction
242 * @param xmax: max. value of the range in x-direction
243 * @param nbinsy: number of bins in y-direction
244 * @param ymin: min. value of the range in y-direction
245 * @param ymax: max. value of the range in y-direction
246 * @param nbinsz: number of bins in z-direction
247 * @param zmin: min. value of the range in z-direction
248 * @param zmax: max. value of the range in z-direction
249 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
251 TString dirname(basename(name)), hname(histname(name));
252 THashList *parent(FindGroup(dirname.Data()));
254 Fatal("THistManager::CreateTH3", "Parent %s does not exist", dirname.Data());
255 if(parent->FindObject(hname.Data()))
256 Fatal("THistManager::CreateTH3", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
257 parent->Add(new TH3D(hname.Data(), title, nbinsx, xmin, xmax, nbinsy, ymin, ymax, nbinsz, zmin, zmax));
260 //______________________________________________________________________________
261 void THistManager::CreateTH3(const char* name, const char* title, int nbinsx, const double* xbins, int nbinsy, const double* ybins, int nbinsz, const double* zbins) {
263 * Create a new TH3 within the container. The histogram name also contains the parent group(s) according to the common
266 * @param name: Name of the histogram
267 * @param title: Title of the histogram
268 * @param nbinsx: number of bins in x-direction
269 * @param xbins: array of bin limits in x-direction
270 * @param nbinsy: number of bins in y-direction
271 * @param ybins: array of bin limits in y-direction
272 * @param nbinsz: number of bins in z-direction
273 * @param zbins: array of bin limits in z-direction
274 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
276 TString dirname(basename(name)), hname(histname(name));
277 THashList *parent(FindGroup(dirname.Data()));
279 Fatal("THistManager::CreateTH3", "Parent %s does not exist", dirname.Data());
280 if(parent->FindObject(hname.Data()))
281 Fatal("THistManager::CreateTH3", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
282 parent->Add(new TH3D(hname.Data(), title, nbinsx, xbins, nbinsy, ybins, nbinsz, zbins));
285 //______________________________________________________________________________
286 void THistManager::CreateTH3(const char* name, const char* title, const TArrayD& xbins, const TArrayD& ybins, const TArrayD& zbins) {
288 * Create a new TH3 within the container. The histogram name also contains the parent group(s) according to the common
291 * @param name: Name of the histogram
292 * @param title: Title of the histogram
293 * @param xbins: array of bin limits in x-direction (contains also the number of bins)
294 * @param ybins: array of bin limits in y-direction (contains also the number of bins)
295 * @param zbins: array of bin limits in z-direction (contains also the number of bins)
296 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
298 TString dirname(basename(name)), hname(histname(name));
299 THashList *parent(FindGroup(dirname.Data()));
301 Fatal("THistManager::CreateTH3", "Parent %s does not exist", dirname.Data());
302 if(parent->FindObject(hname.Data()))
303 Fatal("THistManager::CreateTH3", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
304 parent->Add(new TH3D(hname.Data(), title, xbins.GetSize()-1, xbins.GetArray(), ybins.GetSize()-1, ybins.GetArray(), zbins.GetSize()-1, zbins.GetArray()));
307 //______________________________________________________________________________
308 void THistManager::CreateTHnSparse(const char *name, const char *title, int ndim, const int *nbins, const double *min, const double *max) {
310 * Create a new THnSparse within the container. The histogram name also contains the parent group(s) according to the common
313 * @param name: Name of the histogram
314 * @param title: Title of the histogram
315 * @param ndim: Number of dimensions
316 * @param nbins: Number of bins per dimension
317 * @param min: min. value of the range for each dimension
318 * @param max: max. value of the range for each dimension
319 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
321 TString dirname(basename(name)), hname(histname(name));
322 THashList *parent(FindGroup(dirname.Data()));
324 Fatal("THistManager::CreateTHnSparse", "Parent %s does not exist", dirname.Data());
325 if(parent->FindObject(hname.Data()))
326 Fatal("THistManager::CreateTHnSparse", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
327 parent->Add(new THnSparseD(hname.Data(), title, ndim, nbins, min, max));
330 //______________________________________________________________________________
331 void THistManager::CreateTHnSparse(const char *name, const char *title, int ndim, const TAxis **axes) {
333 * Create a new THnSparse within the container. The histogram name also contains the parent group(s) according to the common
336 * @param name: Name of the histogram
337 * @param title: Title of the histogram
338 * @param ndim: Number of dimensions
339 * @param axes: Array of pointers to TAxis for containing the axis definition for each dimension
340 * Raises fatals in case the parent group does not exist or the object is attempted to be duplicated within the group
342 TString dirname(basename(name)), hname(histname(name));
343 THashList *parent(FindGroup(dirname.Data()));
345 Fatal("THistManager::CreateTHnSparse", "Parent %s does not exist", dirname.Data());
346 if(parent->FindObject(hname))
347 Fatal("THistManager::CreateTHnSparse", "Object %s already exists in group %s does not exist", hname.Data(), dirname.Data());
348 TArrayD xmin(ndim), xmax(ndim);
350 for(int idim = 0; idim < ndim; ++idim){
351 const TAxis &myaxis = *(axes[idim]);
352 nbins[idim] = myaxis.GetNbins();
353 xmin[idim] = myaxis.GetXmin();
354 xmax[idim] = myaxis.GetXmax();
356 THnSparseD *hsparse = new THnSparseD(hname.Data(), title, ndim, nbins.GetArray(), xmin.GetArray(), xmax.GetArray());
357 for(int id = 0; id < ndim; ++id)
358 *(hsparse->GetAxis(id)) = *(axes[id]);
359 parent->Add(hsparse);
362 //______________________________________________________________________________
363 void THistManager::SetObject(TObject * const o, const char *group) {
365 * Set a new group into the container into the parent group
367 * @param o: the object ot be included
368 * Raises fatals in case the parent group is not found, the object is attempted to be duplicated in the group, or the type of
369 * the object is not a histogram type
371 THashList *parent(FindGroup(group));
373 Fatal("THistManager::SetObject", "Parent %s does not exist", strcmp(group, "/") ? group : "");
374 if(parent->FindObject(o->GetName()))
375 Fatal("THistManager::SetObject", "Parent %s does not exist", strcmp(group, "/") ? group : "");
376 if(!(dynamic_cast<THnBase *>(o) || dynamic_cast<TH1 *>(o)))
377 Fatal("THistManager::SetObject", "Object %s is not of a histogram type",o->GetName());
381 //______________________________________________________________________________
382 void THistManager::FillTH1(const char *name, double x, double weight) {
384 * Fill a 1D histogram within the container. The histogram name also contains the parent group(s) according to the common
387 * @param name: Name of the histogram
388 * @param x: x-coordinate
389 * @param weight (@default 1): optional weight of the entry
390 * Raises fatals in case the parent group is not found or the histogram is not found in the parent group
392 TString dirname(basename(name)), hname(histname(name));
393 THashList *parent(FindGroup(dirname.Data()));
395 Fatal("THistManager::FillTH1", "Parnt group %s does not exist", dirname.Data());
396 TH1 *hist = dynamic_cast<TH1 *>(parent->FindObject(hname.Data()));
398 Fatal("THistManager::FillTH1", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
399 hist->Fill(x, weight);
402 //______________________________________________________________________________
403 void THistManager::FillTH2(const char *name, double x, double y, double weight) {
405 * Fill a 2D histogram within the container. The histogram name also contains the parent group(s) according to the common
408 * @param name: Name of the histogram
409 * @param x: x-coordinate
410 * @param y: y-coordinate
411 * @param weight (@default 1): optional weight of the entry
412 * Raises fatals in case the parent group is not found or the histogram is not found in the parent group
414 TString dirname(basename(name)), hname(histname(name));
415 THashList *parent(FindGroup(dirname.Data()));
417 Fatal("THistManager::FillTH2", "Parnt group %s does not exist", dirname.Data());
418 TH2 *hist = dynamic_cast<TH2 *>(parent->FindObject(hname.Data()));
420 Fatal("THistManager::FillTH2", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
421 hist->Fill(x, y, weight);
424 //______________________________________________________________________________
425 void THistManager::FillTH2(const char *name, double *point, double weight) {
427 * Fill a 2D histogram within the container. The histogram name also contains the parent group(s) according to the common
430 * @param name: Name of the histogram
431 * @param point: coordinates of the data
432 * @param weight (@default 1): optional weight of the entry
433 * Raises fatals in case the parent group is not found or the histogram is not found in the parent group
435 TString dirname(basename(name)), hname(histname(name));
436 THashList *parent(FindGroup(dirname.Data()));
438 Fatal("THistManager::FillTH2", "Parnt group %s does not exist", dirname.Data());
439 TH2 *hist = dynamic_cast<TH2 *>(parent->FindObject(hname.Data()));
441 Fatal("THistManager::FillTH2", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
442 hist->Fill(point[0], point[1], weight);
445 //______________________________________________________________________________
446 void THistManager::FillTH3(const char* name, double x, double y, double z, double weight) {
448 * Fill a 3D histogram within the container. The histogram name also contains the parent group(s) according to the common
451 * @param name: Name of the histogram
452 * @param x: x-coordinate
453 * @param y: y-coordinate
454 * @param z: z-coordinate
455 * @param weight (@default 1): optional weight of the entry
456 * Raises fatals in case the parent group is not found or the histogram is not found in the parent group
458 TString dirname(basename(name)), hname(histname(name));
459 THashList *parent(FindGroup(dirname.Data()));
461 Fatal("THistManager::FillTH3", "Parnt group %s does not exist", dirname.Data());
462 TH3 *hist = dynamic_cast<TH3 *>(parent->FindObject(hname.Data()));
464 Fatal("THistManager::FillTH3", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
465 hist->Fill(x, y, z, weight);
468 //______________________________________________________________________________
469 void THistManager::FillTH3(const char* name, const double* point, double weight) {
471 * Fill a 3D histogram within the container. The histogram name also contains the parent group(s) according to the common
474 * @param name: Name of the histogram
475 * @param point: 3D-coordinate of the point
476 * @param weight (@default 1): optional weight of the entry
477 * Raises fatals in case the parent group is not found or the histogram is not found in the parent group
479 TString dirname(basename(name)), hname(histname(name));
480 THashList *parent(FindGroup(dirname.Data()));
482 Fatal("THistManager::FillTH3", "Parnt group %s does not exist", dirname.Data());
483 TH3 *hist = dynamic_cast<TH3 *>(parent->FindObject(hname.Data()));
485 Fatal("THistManager::FillTH3", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
486 hist->Fill(point[0], point[1], point[2], weight);
490 //______________________________________________________________________________
491 void THistManager::FillTHnSparse(const char *name, const double *x, double weight) {
493 * Fill a nD histogram within the container. The histogram name also contains the parent group(s) according to the common
496 * @param name: Name of the histogram
497 * @param x: coordinates of the data
498 * @param weight (@default 1): optional weight of the entry
499 * Raises fatals in case the parent group is not found or the histogram is not found in the parent group
501 TString dirname(basename(name)), hname(histname(name));
502 THashList *parent(FindGroup(dirname.Data()));
504 Fatal("THistManager::FillTH3", "Parnt group %s does not exist", dirname.Data());
505 THnSparseD *hist = dynamic_cast<THnSparseD *>(parent->FindObject(hname.Data()));
507 Fatal("THistManager::FillTH3", "Histogram %s not found in parent group %s", hname.Data(), dirname.Data());
508 hist->Fill(x, weight);
511 //______________________________________________________________________________
512 TObject *THistManager::FindObject(const char *name) const {
514 * Find an object inside the container. The object can also be within a
515 * histogram group. For this the name has to follow the common notation
517 * @param name: Name of the object to find inside the container
518 * @return: pointer to the object (NULL if not found)
520 TString dirname(basename(name)), hname(histname(name));
521 THashList *parent(FindGroup(dirname.Data()));
522 if(!parent) return NULL;
523 return parent->FindObject(name);
526 //______________________________________________________________________________
527 TObject* THistManager::FindObject(const TObject* obj) const {
529 * Find and object inside the container. The object name is expected to contain the
530 * full path of the histogram object, including parent groups
532 * @param obj: the object to find
533 * @return: pointer to the object (NULL if not found)
535 TString dirname(basename(obj->GetName())), hname(histname(obj->GetName()));
536 THashList *parent(FindGroup(dirname.Data()));
537 if(!parent) return NULL;
538 return parent->FindObject(hname);
541 //______________________________________________________________________________
542 THashList *THistManager::FindGroup(const char *dirname) const {
544 * Find histogram group. Name is using common notation
546 * @param dirname: Path of the group (treat empty path as top node
547 * @return: TList of objects (NULL if group does not exist)
549 if(!strlen(dirname) || !strcmp(dirname, "/")) return fHistos;
550 std::vector<std::string> tokens;
551 TokenizeFilename(dirname, "/", tokens);
552 THashList *currentdir(fHistos);
553 for(std::vector<std::string>::iterator it = tokens.begin(); it != tokens.end(); ++it){
554 currentdir = dynamic_cast<THashList *>(currentdir->FindObject(it->c_str()));
555 if(!currentdir) break;
560 //______________________________________________________________________________
561 void THistManager::TokenizeFilename(const char *name, const char *delim, std::vector<std::string> &listoftokens) const {
563 * Tokenizes a string. Results are stored inside the vector listoftokens
565 * @ param name: string to be tokenised
566 * @ param delim: delimiter string
567 * @ param listoftokens: list of tokens (C++ strings)
570 TObjArray *arr = s.Tokenize(delim);
571 TObjString *ostr(NULL);
573 while((ostr = dynamic_cast<TObjString *>(toks()))){
574 listoftokens.push_back(std::string(ostr->String().Data()));
579 //______________________________________________________________________________
580 const char *THistManager::basename(const char *path) const {
582 * Helper function extracting the basename from a given histogram path.
584 * @param path: histogram path
585 * @return: basename extracted
588 int index = s.Last('/');
589 if(index < 0) return ""; // no directory structure
590 return TString(s(0, index)).Data();
593 //______________________________________________________________________________
594 const char *THistManager::histname(const char *path) const {
596 * Helper function extracting the histogram name from a given histogram path.
598 * @param path: histogram path
599 * @return: basename extracted
602 int index = s.Last('/');
603 if(index < 0) return path; // no directory structure
604 return TString(s(index+1, s.Length() - (index+1))).Data();