added SelectCollisionCandidates J. Otwinowski
[u/mrichter/AliRoot.git] / PWG3 / muon / AliCounterCollection.cxx
CommitLineData
94ef1a28 1/**************************************************************************
2 * Copyright(c) 1998-2007, 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
16//-----------------------------------------------------------------------------
17/// \class AliCounterCollection
18///
19/// generic class to handle a collection of counters
20///
21/// \author Philippe Pillot
22//-----------------------------------------------------------------------------
23
24#include "AliCounterCollection.h"
25
26#include <AliLog.h>
27
28#include <TString.h>
29#include <TObjString.h>
30#include <TObjArray.h>
31#include <THnSparse.h>
32#include <THashList.h>
33#include <TArrayI.h>
34#include <TH1D.h>
35#include <TH2D.h>
36#include <TCollection.h>
37
38ClassImp(AliCounterCollection)
39
40//-----------------------------------------------------------------------
41AliCounterCollection::AliCounterCollection(const char* name) :
42TNamed(name,name),
43fRubrics(new THashList(10)),
44fRubricsSize(new TArrayI(10)),
45fCounters(0x0)
46{
47 /// Constructor
48 fRubrics->SetOwner();
49}
50
51//-----------------------------------------------------------------------
52AliCounterCollection::~AliCounterCollection()
53{
54 /// Destructor
55 delete fRubrics;
56 delete fRubricsSize;
57 delete fCounters;
58}
59
60//-----------------------------------------------------------------------
61void AliCounterCollection::Clear(Option_t*)
62{
63 /// Clear counters
64 fRubrics->Clear();
65 fRubricsSize->Reset();
66 delete fCounters; fCounters = 0x0;
67}
68
69//-----------------------------------------------------------------------
70void AliCounterCollection::AddRubric(TString name, TString listOfKeyWords)
71{
72 /// Add a new rubric with the complete list of related key words separated by "/".
73 /// If the key word "any" is not defined, the overall statistics is
74 /// assumed to be the sum of the statistics under each key word.
75
76 name.ToUpper();
77 listOfKeyWords.ToUpper();
78
79 if (fRubrics->Contains(name.Data())) {
80 AliError(Form("rubric named %s already exist",name.Data()));
81 return;
82 }
83
84 // add the list of autorized key words
85 TObjArray* rubric = listOfKeyWords.Tokenize("/");
86 CleanListOfStrings(rubric);
87 rubric->SetName(name.Data());
88 Int_t nRubrics = fRubrics->GetSize();
89 rubric->SetUniqueID(nRubrics);
90 fRubrics->AddLast(rubric);
91
92 // save the number of autorized key words (expand the array if needed)
93 if (nRubrics+1 > fRubricsSize->GetSize()) fRubricsSize->Set(2*fRubricsSize->GetSize());
94 (*fRubricsSize)[nRubrics] = rubric->GetEntriesFast();
95}
96
97//-----------------------------------------------------------------------
98void AliCounterCollection::AddRubric(TString name, Int_t maxNKeyWords)
99{
100 /// Add a new rubric containing at maximum maxNKeyWords key words.
101 /// Key words will be added as the counters get filled until the maximum is reached.
102 /// If the key word "any" is never defined, the overall statistics is
103 /// assumed to be the sum of the statistics under each key word.
104
105 name.ToUpper();
106
107 if (fRubrics->Contains(name.Data())) {
108 AliError(Form("rubric named %s already exist",name.Data()));
109 return;
110 }
111
112 // create the empty rubric
113 TObjString* rubric = new TObjString(name.Data());
114 Int_t nRubrics = fRubrics->GetSize();
115 rubric->SetUniqueID(nRubrics);
116 fRubrics->AddLast(rubric);
117
118 // save the maximum number of autorized key words
119 if (nRubrics+1 > fRubricsSize->GetSize()) fRubricsSize->Set(2*fRubricsSize->GetSize());
120 (*fRubricsSize)[nRubrics] = maxNKeyWords;
121}
122
123//-----------------------------------------------------------------------
124void AliCounterCollection::Init()
125{
126 /// Initialize the internal counters from the added rubrics.
127
128 // create the counters
129 delete fCounters;
130 fCounters = new THnSparseT<TArrayI>("hCounters", "hCounters", fRubrics->GetSize(), fRubricsSize->GetArray(), 0x0, 0x0);
131
132 // loop over axis
133 TObject* rubric = 0x0;
134 TIter nextRubric(fRubrics);
135 while ((rubric = nextRubric())) {
136 TAxis* axis = fCounters->GetAxis((Int_t)rubric->GetUniqueID());
137
138 // set axis name
139 axis->SetName(rubric->GetName());
140
141 // set labels if already known
142 TObjArray* keyWords = dynamic_cast<TObjArray*>(rubric);
143 if (keyWords) {
144 TObjString* label = 0x0;
145 Int_t bin = 1;
146 TIter nextLabel(keyWords);
147 while ((label = static_cast<TObjString*>(nextLabel()))) axis->SetBinLabel(bin++, label->String().Data());
148 }
149
150 }
151
152}
153
154//-----------------------------------------------------------------------
155const Int_t* AliCounterCollection::FindBins(const TString& externalKey, Bool_t allocate, Int_t& nEmptySlots)
156{
157 /// Return the corresponding bins ordered by rubric or 0x0 if externalKey is not valid.
158 /// The externalKey format must be rubric:keyWord/rubric:keyWord/rubric:keyWord/...
159 /// If allocate = kTRUE, new key words are added to the corresponding rubric if possible.
160 /// If a rubric is not filled in, the coresponding slot contain -1 in the array.
161 /// It is the responsability of the user to delete the returned array.
162
163 // produce an empty array of keys
164 Int_t nRubrics = fRubrics->GetSize();
165 Int_t* bins = new Int_t[nRubrics];
166 memset(bins, -1, sizeof(Int_t) * nRubrics);
167 nEmptySlots = nRubrics;
168 Bool_t isValid = kTRUE;
169
170 // get the list of rubric:keyWord pairs
171 TObjArray* rubricKeyPairs = externalKey.Tokenize("/");
172
173 // loop over each rubric:keyWord pair
174 TObjString* pair = 0x0;
175 TIter next(rubricKeyPairs);
176 while ((pair = static_cast<TObjString*>(next()))) {
177
178 // get both rubric and associated key word
179 TObjArray* rubricKeyPair = pair->String().Tokenize(":");
180
181 // check the format of the pair
182 if (rubricKeyPair->GetEntriesFast() != 2) {
183 AliError("invalid key format");
184 isValid = kFALSE;
185 delete rubricKeyPair;
186 break;
187 }
188
189 // get the axis corresponding to that rubric
190 Int_t dim = FindDim(static_cast<TObjString*>(rubricKeyPair->UncheckedAt(0))->String());
191 if (dim < 0) {
192 isValid = kFALSE;
193 delete rubricKeyPair;
194 break;
195 }
196
197 // find the bin corresponding to that key word
198 Int_t bin = FindBin(dim, static_cast<TObjString*>(rubricKeyPair->UncheckedAt(1))->String(), allocate);
199 if (bin < 0) {
200 isValid = kFALSE;
201 delete rubricKeyPair;
202 break;
203 }
204
205 // check if the array of keys already contains something for that rubric
206 if (bins[dim] >= 0) {
207 AliWarning("key already given for that rubric --> ignored");
208 delete rubricKeyPair;
209 continue;
210 }
211
212 // store the corresponding bin for that slot
213 bins[dim] = bin;
214 nEmptySlots--;
215
216 // clean memory
217 delete rubricKeyPair;
218 }
219
220 // delete the array in case of problem
221 if (!isValid) {
222 delete[] bins;
223 bins = 0x0;
224 nEmptySlots = nRubrics;
225 }
226
227 // clean memory
228 delete rubricKeyPairs;
229
230 return bins;
231}
232
233//-----------------------------------------------------------------------
234Int_t AliCounterCollection::FindDim(const TString& rubricName) const
235{
236 /// Return the dimension corresponding to that rubric (or -1 in case of failure).
237 TObject* rubric = fRubrics->FindObject(rubricName.Data());
238 if (!rubric) {
239 AliError(Form("invalid rubric: %s",rubricName.Data()));
240 return -1;
241 }
242 return (Int_t) rubric->GetUniqueID();
243}
244
245//-----------------------------------------------------------------------
246Int_t AliCounterCollection::FindBin(Int_t dim, const TString& keyWord, Bool_t allocate)
247{
248 /// Return the bin number corresponding to that key word (or -1 in case of failure).
249 /// If allocate = kTRUE, try to add the key word if possible.
250
251 TAxis* axis = fCounters->GetAxis(dim);
252
253 // look for the bin corresponding to keyWord
254 THashList* labels = axis->GetLabels();
255 TObjString* label = (labels) ? static_cast<TObjString*>(labels->FindObject(keyWord.Data())) : 0x0;
256 Int_t bin = (label) ? (Int_t)label->GetUniqueID() : -1;
257
258 // in case the keyWord does not exist, try to add it if required
259 if (bin<0 && allocate) {
260 Int_t nLabels = (labels) ? labels->GetSize() : 0;
261 if (nLabels < axis->GetNbins()) {
262 bin = nLabels+1;
263 axis->SetBinLabel(bin, keyWord.Data());
264 }
265 }
266
267 if (bin<0) AliError(Form("invalid key word: %s:%s",axis->GetName(),keyWord.Data()));
268
269 return bin;
270}
271
272//-----------------------------------------------------------------------
273void AliCounterCollection::CleanListOfStrings(TObjArray* list)
274{
275 /// Make sure all strings appear only once in this list
276
277 // remove multiple-occurrence
278 Int_t nEntries = list->GetEntriesFast();
279 for (Int_t i = 0; i < nEntries; i++) {
280 TObjString* entry1 = static_cast<TObjString*>(list->UncheckedAt(i));
281 if (!entry1) continue;
282 for (Int_t j = i+1; j < nEntries; j++) {
283 TObjString* entry2 = static_cast<TObjString*>(list->UncheckedAt(j));
284 if (entry2 && entry2->IsEqual(entry1)) {
285 AliWarning(Form("multiple-occurence of string \"%s\" --> removed",entry2->String().Data()));
286 list->RemoveAt(j);
287 }
288 }
289 }
290
291 // remove empty slots
292 list->Compress();
293}
294
295//-----------------------------------------------------------------------
296void AliCounterCollection::Count(TString externalKey, Int_t value)
297{
298 /// Add "value" to the counter referenced by "externalKey".
299 /// The externalKey format must be rubric:keyWord/rubric:keyWord/rubric:keyWord/...
300
301 if (value < 1) return;
302 if (!fCounters) {
303 AliError("counters are not initialized");
304 return;
305 }
306
307 externalKey.ToUpper();
308
309 // convert external to internal key
310 Int_t nEmptySlots = 0;
311 const Int_t* bins = FindBins(externalKey, kTRUE, nEmptySlots);
312 if (!bins) return;
313
314 // check for empty slots
315 if (nEmptySlots > 0) {
316 AliError("incomplete key");
317 delete[] bins;
318 return;
319 }
320
321 // increment the corresponding counter
322 fCounters->AddBinContent(bins, (Double_t)value);
323
324 // clean memory
325 delete[] bins;
326}
327
328//-----------------------------------------------------------------------
329void AliCounterCollection::Print(const Option_t* opt) const
330{
331 /// Print every individual counters if opt=="", else call "Print(TString rubrics=opt, TString selections="")".
332
333 if (strcmp(opt,"")) {
334 const_cast<AliCounterCollection*>(this)->Print(opt, "");
335 return;
336 }
337
338 if (!fCounters) {
339 AliError("counters are not initialized");
340 return;
341 }
342
343 if (fCounters->GetNbins() == 0) {
344 printf("\nall counters are empty\n\n");
345 return;
346 }
347
348 Int_t nRubrics = fCounters->GetNdimensions();
349 Int_t* bins = new Int_t[nRubrics];
350
351 // loop over every filled counters
352 for (Long64_t i=0; i<fCounters->GetNbins(); ++i) {
353
354 // get the content of the bin
355 Int_t value = (Int_t) fCounters->GetBinContent(i, bins);
356
357 // build the corresponding counter name
358 TString counter;
359 for (Int_t j=0; j<nRubrics; j++) counter += Form("/%s",fCounters->GetAxis(j)->GetBinLabel(bins[j]));
360 counter += "/";
361
362 // print value
363 printf("\n%s %d", counter.Data(), value);
364 }
365 printf("\n\n");
366
367 // clean memory
368 delete[] bins;
369}
370
371//-----------------------------------------------------------------------
372void AliCounterCollection::PrintKeyWords() const
373{
374 /// Print the full list of key words.
375
376 if (!fCounters) {
377 AliError("counters are not initialized");
378 return;
379 }
380
381 // loop over rubrics
382 Int_t nRubrics = fCounters->GetNdimensions();
383 for (Int_t iDim=0; iDim<nRubrics; iDim++) {
384 TAxis* axis = fCounters->GetAxis(iDim);
385
386 // print rubric's name
387 printf("\n%s:", axis->GetName());
388
389 // loop over key words
390 Bool_t first = kTRUE;
391 TObjString* label = 0x0;
392 TIter nextLabel(axis->GetLabels());
393 while ((label = static_cast<TObjString*>(nextLabel()))) {
394
395 //print key word's name
396 if (first) {
397 printf("%s", label->String().Data());
398 first = kFALSE;
399 } else printf(",%s", label->String().Data());
400
401 }
402 }
403 printf("\n\n");
404}
405
406//-----------------------------------------------------------------------
407void AliCounterCollection::PrintValue(TString selections)
408{
409 /// Print value of selected counter.
410
411 if (!fCounters) {
412 AliError("counters are not initialized");
413 return;
414 }
415
416 selections.ToUpper();
417
418 // convert external to internal key
419 Int_t nEmptySlots = 0;
420 const Int_t* selectedBins = FindBins(selections, kFALSE, nEmptySlots);
421 if (!selectedBins) return;
422
423 // check for empty slots
424 if (nEmptySlots > 0) {
425 AliError("incomplete key");
426 delete[] selectedBins;
427 return;
428 }
429
430 // print value
431 printf("\n%d\n\n", (Int_t) fCounters->GetBinContent(selectedBins));
432
433 // clean memory
434 delete[] selectedBins;
435}
436
437//-----------------------------------------------------------------------
438void AliCounterCollection::Print(TString rubrics, TString selections)
439{
440 /// Print desired rubrics for the given selection:
441 /// - format of "rubrics" is rubric1/rubric2/.. (order matters only for output).
442 /// - format of "selections" is rubric:keyWord/rubric:keyWord/.. (order does not matter).
443 /// If "data" contains 1 rubric, the output will be one counter for each element of that rubric.
444 /// If "data" contains 2 rubrics, the output will be an array of counters, rubric1 vs rubric2.
445 /// If "data" contains 3 rubrics, the output will be an array rubric1 vs rubric2 for each element in rubric3.
446 /// ...
447 /// Results are integrated over rubrics not specified neither in "rubrics" nor in "selections".
448
449 if (!fCounters) {
450 AliError("counters are not initialized");
451 return;
452 }
453
454 rubrics.ToUpper();
455 selections.ToUpper();
456
457 // get the rubrics to print
458 TObjArray* rubricsToPrint = rubrics.Tokenize("/");
459 if (rubricsToPrint->GetEntriesFast() == 0) {
460 delete rubricsToPrint;
461 return;
462 }
463
464 // remove rubrics called twice
465 CleanListOfStrings(rubricsToPrint);
466
467 // project counters in the rubrics to print according to the selections
468 TObject* hist = Projection(*rubricsToPrint, selections);
469 if (!hist) {
470 delete rubricsToPrint;
471 return;
472 }
473
474 // print counters
475 Int_t nRubricsToPrint = rubricsToPrint->GetEntriesFast();
476 if (nRubricsToPrint == 1 && (static_cast<TH1D*>(hist))->Integral() > 0.)
477 PrintList(static_cast<TH1D*>(hist));
478 else if (nRubricsToPrint == 2 && (static_cast<TH2D*>(hist))->Integral() > 0.)
479 PrintArray(static_cast<TH2D*>(hist));
480 else if (nRubricsToPrint > 2 && (static_cast<THnSparse*>(hist))->GetNbins() > 0)
481 PrintListOfArrays(static_cast<THnSparse*>(hist));
482 else
483 printf("\nselected counters are empty\n\n");
484
485 // clean memory
486 delete rubricsToPrint;
487 delete hist;
488}
489
490//-----------------------------------------------------------------------
491void AliCounterCollection::PrintSum(TString rubric, TString selections)
492{
493 /// Print the overall statistics under the given rubric for the given selection:
494 /// - format of "selections" is rubric:keyWord/rubric:keyWord/.. (order does not matter).
495 /// Result is integrated over rubrics not specified neither in "rubric" nor in "selections".
496
497 if (!fCounters) {
498 AliError("counters are not initialized");
499 return;
500 }
501
502 rubric.ToUpper();
503 selections.ToUpper();
504
505 // fill the rubric to sum
506 TObjArray rubricsToSum(1);
507 rubricsToSum.SetOwner();
508 rubricsToSum.AddLast(new TObjString(rubric.Data()));
509
510 // project counters in the rubric to sum according to the selections
511 TH1D* hist = static_cast<TH1D*>(Projection(rubricsToSum, selections));
512 if (!hist) return;
513
514 // check for empty rubric
515 THashList* labels = hist->GetXaxis()->GetLabels();
516 if (!labels) {
517 printf("\n0\n\n");
518 return;
519 }
520
521 // print the sum of counters under that rubric
522 TObjString* any = static_cast<TObjString*>(labels->FindObject("ANY"));
523 if (any) printf("\n%d\n\n", (Int_t) hist->GetBinContent((Int_t)any->GetUniqueID()));
524 else printf("\n%d\n\n", (Int_t) hist->Integral());
525
526 // clean memory
527 delete hist;
528}
529
530//-----------------------------------------------------------------------
531void AliCounterCollection::PrintList(const TH1D* hist) const
532{
533 /// Print the content of 1D histogram as a list.
534
535 // set the format to print labels
536 THashList* labels = hist->GetXaxis()->GetLabels();
537 TString format(Form("\n%%%ds %%9d",GetMaxLabelSize(labels)));
538
539 // print value for each label
540 TObjString* label = 0x0;
541 TIter nextLabel(labels);
542 while ((label = static_cast<TObjString*>(nextLabel()))) {
543 Int_t bin = (Int_t) label->GetUniqueID();
544 printf(format.Data(), label->String().Data(), (Int_t) hist->GetBinContent(bin));
545 }
546 printf("\n\n");
547}
548
549//-----------------------------------------------------------------------
550void AliCounterCollection::PrintArray(const TH2D* hist) const
551{
552 /// Print the content of 2D histogram as an array.
553
554 // set the format to print labels in X direction
555 THashList* labelsX = hist->GetXaxis()->GetLabels();
556 TString formatX(Form("\n%%%ds ",GetMaxLabelSize(labelsX)));
557
558 // set the format to print labels in Y direction and values
559 THashList* labelsY = hist->GetYaxis()->GetLabels();
560 Int_t maxLabelSizeY = TMath::Max(9, GetMaxLabelSize(labelsY));
561 TString formatYs(Form("%%%ds ",maxLabelSizeY));
562 TString formatYd(Form("%%%dd ",maxLabelSizeY));
563
564 // print labels in Y axis
565 printf(formatX.Data()," ");
566 TObjString* labelY = 0x0;
567 TIter nextLabelY(labelsY);
568 while ((labelY = static_cast<TObjString*>(nextLabelY())))
569 printf(formatYs.Data(), labelY->String().Data());
570
571 // fill array for each label in X axis
572 TObjString* labelX = 0x0;
573 TIter nextLabelX(labelsX);
574 while ((labelX = static_cast<TObjString*>(nextLabelX()))) {
575 Int_t binX = (Int_t) labelX->GetUniqueID();
576
577 // print label X
578 printf(formatX.Data(), labelX->String().Data());
579
580 // print value for each label in Y axis
581 nextLabelY.Reset();
582 while ((labelY = static_cast<TObjString*>(nextLabelY()))) {
583 Int_t binY = (Int_t) labelY->GetUniqueID();
584 printf(formatYd.Data(), (Int_t) hist->GetBinContent(binX, binY));
585 }
586 }
587 printf("\n\n");
588}
589
590//-----------------------------------------------------------------------
591void AliCounterCollection::PrintListOfArrays(const THnSparse* hist) const
592{
593 /// Print the content of nD histogram as a list of arrays.
594
595 // set the format to print labels in X direction
596 THashList* labelsX = hist->GetAxis(0)->GetLabels();
597 TString formatX(Form("\n%%%ds ",GetMaxLabelSize(labelsX)));
598
599 // set the format to print labels in Y direction and values
600 THashList* labelsY = hist->GetAxis(1)->GetLabels();
601 Int_t maxLabelSizeY = TMath::Max(9, GetMaxLabelSize(labelsY));
602 TString formatYs(Form("%%%ds ",maxLabelSizeY));
603 TString formatYd(Form("%%%dd ",maxLabelSizeY));
604
605 // create a list containing each combination of labels refering the arrays to be printout
606 TList listOfCombis;
607 listOfCombis.SetOwner();
608
609 // add a first empty combination
610 Int_t nDim = hist->GetNdimensions();
611 listOfCombis.AddLast(new TObjArray(nDim-2));
612
613 // loop over the nDim-2 other rubrics
614 for (Int_t i=2; i<nDim; i++) {
615
616 // save the last label of that rubic
617 THashList* labels = hist->GetAxis(i)->GetLabels();
618 TObjString* lastLabel = (labels) ? static_cast<TObjString*>(labels->Last()) : 0x0;
619 if (!lastLabel) return;
620
621 // prepare iteration over the list of labels
622 TIter nextLabel(labels);
623
624 // loop over existing combinations
625 TObjLink* lnk = listOfCombis.FirstLink();
626 while (lnk) {
627
628 // get the current combination
629 TObjArray* currentCombi = static_cast<TObjArray*>(lnk->GetObject());
630
631 // loop over labels in the current rubric
632 nextLabel.Reset();
633 TObjString* label = 0x0;
634 while ((label = static_cast<TObjString*>(nextLabel()))) {
635
636 // stop at the last one
637 if (label == lastLabel) break;
638
639 // copy the current combination, add the current label to it and add it to the list of combinations
640 TObjArray* combi = new TObjArray(*currentCombi);
641 combi->AddLast(label);
642 listOfCombis.AddBefore(lnk, combi);
643 }
644
645 // add the last label to the current combination
646 currentCombi->AddLast(lastLabel);
647
648 lnk = lnk->Next();
649 }
650
651 }
652
653 // create bin coordinates to access individual counters
654 Int_t* bins = new Int_t[nDim];
655
656 // loop over each combination of labels
657 TObjArray* combi = 0x0;
658 TIter nextCombi(&listOfCombis);
659 while ((combi = static_cast<TObjArray*>(nextCombi()))) {
660
661 // make the name of the combination and fill the corresponding bin coordinates
662 TString combiName = "/";
663 for (Int_t i=2; i<nDim; i++) {
664 TObjString* label = static_cast<TObjString*>(combi->UncheckedAt(i-2));
665 combiName += Form("%s/",label->String().Data());
666 bins[i] = (Int_t)label->GetUniqueID();
667 }
668
669 // skip empty array
670 Bool_t empty = kTRUE;
671 TObjString* labelX = 0x0;
672 TObjString* labelY = 0x0;
673 TIter nextLabelX(labelsX);
674 TIter nextLabelY(labelsY);
675 while ((labelX = static_cast<TObjString*>(nextLabelX()))) {
676 bins[0] = (Int_t) labelX->GetUniqueID();
677 nextLabelY.Reset();
678 while ((labelY = static_cast<TObjString*>(nextLabelY()))) {
679 bins[1] = (Int_t) labelY->GetUniqueID();
680 if (((Int_t) hist->GetBinContent(bins)) > 0) {
681 empty = kFALSE;
682 break;
683 }
684 }
685 if (!empty) break;
686 }
687 if (empty) continue;
688
689 // print the name of the combination of labels refering the incoming array
690 printf("\n%s:\n",combiName.Data());
691
692 // print labels in Y axis
693 printf(formatX.Data()," ");
694 nextLabelY.Reset();
695 while ((labelY = static_cast<TObjString*>(nextLabelY())))
696 printf(formatYs.Data(), labelY->String().Data());
697
698 // fill array for each label in X axis
699 nextLabelX.Reset();
700 while ((labelX = static_cast<TObjString*>(nextLabelX()))) {
701 bins[0] = (Int_t) labelX->GetUniqueID();
702
703 // print label X
704 printf(formatX.Data(), labelX->String().Data());
705
706 // print value for each label in Y axis
707 nextLabelY.Reset();
708 while ((labelY = static_cast<TObjString*>(nextLabelY()))) {
709 bins[1] = (Int_t) labelY->GetUniqueID();
710 printf(formatYd.Data(), (Int_t) hist->GetBinContent(bins));
711 }
712 }
713 printf("\n\n");
714 }
715
716 // clean memory
717 delete[] bins;
718}
719
720//-----------------------------------------------------------------------
721Int_t AliCounterCollection::GetMaxLabelSize(THashList* labels) const
722{
723 /// Return the number of characters of the longest label.
724 Int_t maxLabelSize = 0;
725 TObjString* label = 0x0;
726 TIter nextLabel(labels);
727 while ((label = static_cast<TObjString*>(nextLabel())))
728 maxLabelSize = TMath::Max(maxLabelSize, label->String().Length());
729 return maxLabelSize;
730}
731
732//-----------------------------------------------------------------------
733TH1D* AliCounterCollection::Draw(TString rubric, TString selections)
734{
735 /// Draw counters of the rubric "rubric" for the given "selection".
736 /// Format of "selections" is rubric:keyWord/rubric:keyWord/.. (order does not matter).
737 /// Results are integrated over rubrics not specified neither in "rubric1" nor in "selections".
738 /// It is the responsability of the user to delete the returned histogram.
739
740 if (!fCounters) {
741 AliError("counters are not initialized");
742 return 0x0;
743 }
744
745 rubric.ToUpper();
746 selections.ToUpper();
747
748 // fill the rubrics to print
749 TObjArray rubricsToPrint(1);
750 rubricsToPrint.SetOwner();
751 rubricsToPrint.AddLast(new TObjString(rubric.Data()));
752
753 // project counters in the rubrics to print according to the selections
754 TH1D* hist = static_cast<TH1D*>(Projection(rubricsToPrint, selections));
755
756 // draw counters
757 if (hist) {
758
759 // draw histogram
760 hist->Draw("htext");
761 hist->SetStats(kFALSE);
762
763 // set title
764 TString title = "Selections: ";
765 selections.Remove(TString::kBoth, '/');
766 if (selections.Length() > 0) title += Form("%s/", selections.Data());
767 TObject* rub = 0x0;
768 TIter nextRubric(fRubrics);
769 while ((rub = nextRubric())) {
770 if (selections.Contains(Form("%s:",rub->GetName()))) continue;
771 if (rubricsToPrint.Contains(rub->GetName())) continue;
772 title += Form("%s:ANY/", rub->GetName());
773 }
774 title.ReplaceAll("/", " ");
775 hist->SetTitle(title.Data());
776
777 // draw X axis
778 TAxis* axis = hist->GetXaxis();
779 THashList* labels = axis->GetLabels();
780 Int_t nLabels = (labels) ? labels->GetSize() : 1;
781 axis->SetRange(1,nLabels);
782 axis->SetNdivisions(1,kFALSE);
783 axis->SetTitle(rubric.Data());
784
785 // draw Y axis
786 hist->GetYaxis()->SetTitle("Counts");
787 }
788
789 return hist;
790}
791
792//-----------------------------------------------------------------------
793TH2D* AliCounterCollection::Draw(TString rubric1, TString rubric2, TString selections)
794{
795 /// Draw counters of the "rubric1" vs "rubric2" for the given "selection".
796 /// Format of "selections" is rubric:keyWord/rubric:keyWord/.. (order does not matter).
797 /// Results are integrated over rubrics not specified neither in "rubric1", "rubric2" nor in "selections".
798 /// It is the responsability of the user to delete the returned histogram.
799
800 if (!fCounters) {
801 AliError("counters are not initialized");
802 return 0x0;
803 }
804
805 rubric1.ToUpper();
806 rubric2.ToUpper();
807 selections.ToUpper();
808
809 // fill the rubrics to print
810 TObjArray rubricsToPrint(2);
811 rubricsToPrint.SetOwner();
812 rubricsToPrint.AddLast(new TObjString(rubric2.Data()));
813 rubricsToPrint.AddLast(new TObjString(rubric1.Data()));
814
815 // project counters in the rubrics to print according to the selections
816 TH2D* hist = static_cast<TH2D*>(Projection(rubricsToPrint, selections));
817
818 // draw counters
819 if (hist) {
820
821 // draw histogram
822 hist->Draw("text");
823 hist->SetStats(kFALSE);
824
825 // set title
826 TString title = "Selections: ";
827 selections.Remove(TString::kBoth, '/');
828 if (selections.Length() > 0) title += Form("%s/", selections.Data());
829 TObject* rub = 0x0;
830 TIter nextRubric(fRubrics);
831 while ((rub = nextRubric())) {
832 if (selections.Contains(Form("%s:",rub->GetName()))) continue;
833 if (rubricsToPrint.Contains(rub->GetName())) continue;
834 title += Form("%s:ANY/", rub->GetName());
835 }
836 title.ReplaceAll("/", " ");
837 hist->SetTitle(title.Data());
838
839 // draw X axis
840 TAxis* axisX = hist->GetXaxis();
841 THashList* labelsX = axisX->GetLabels();
842 Int_t nLabelsX = (labelsX) ? labelsX->GetSize() : 1;
843 axisX->SetRange(1,nLabelsX);
844 axisX->SetNdivisions(1,kFALSE);
845 axisX->SetTitle(rubric2.Data());
846
847 // draw Y axis
848 TAxis* axisY = hist->GetYaxis();
849 THashList* labelsY = axisY->GetLabels();
850 Int_t nLabelsY = (labelsY) ? labelsY->GetSize() : 1;
851 axisY->SetRange(1,nLabelsY);
852 axisY->SetNdivisions(1,kFALSE);
853 axisY->SetTitle(rubric1.Data());
854 }
855
856 return hist;
857}
858
859//-----------------------------------------------------------------------
860TObject* AliCounterCollection::Projection(const TObjArray& data, const TString& selections)
861{
862 /// Return desired "data" for the given "selection" stored in a new histogram or 0x0 in case of failure.
863 /// The type of the histogram (TH1D, TH2D or THnSparse) depend on the number of data.
864 /// It is the responsability of the user to delete the returned histogram.
865
866 // get the corresponding dimensions
867 Int_t nTargetDim = data.GetEntriesFast();
868 Int_t* targetDims = new Int_t[nTargetDim];
869 for (Int_t i=0; i<nTargetDim; i++) {
870 targetDims[i] = FindDim(static_cast<TObjString*>(data.UncheckedAt(i))->String());
871 if (targetDims[i] < 0) {
872 delete[] targetDims;
873 return 0x0;
874 }
875 }
876
877 // find bins to select
878 Int_t nEmptySlots = 0;
879 const Int_t* selectedBins = FindBins(selections, kFALSE, nEmptySlots);
880 if (!selectedBins) {
881 delete[] targetDims;
882 return 0x0;
883 }
884
885 // apply selection for each rubric
886 Int_t nRubrics = fCounters->GetNdimensions();
887 for (Int_t iDim=0; iDim<nRubrics; iDim++) {
888 TAxis* axis = fCounters->GetAxis(iDim);
889
890 // select the desired key word
891 if (selectedBins[iDim] >= 0) axis->SetRange(selectedBins[iDim], selectedBins[iDim]);
892
893 // or select all key words
894 else if (data.Contains(axis->GetName())) axis->SetRange();
895
896 // or integrate over all cases
897 else {
898 THashList* labels = axis->GetLabels();
899 TObjString* label = (labels) ? static_cast<TObjString*>(labels->FindObject("ANY")) : 0x0;
900 Int_t binAny = (label) ? (Int_t)label->GetUniqueID() : -1;
901 if (binAny >= 0) axis->SetRange(binAny, binAny);
902 else axis->SetRange();
903 }
904 }
905
906 // do projection
907 TObject* hist = 0x0;
908 if (nTargetDim == 1) {
909
910 // project counters to TH1D
911 hist = fCounters->Projection(targetDims[0]);
912
913 // reset bin labels lost when producing TH1D
914 if (selectedBins[targetDims[0]] >= 0)
915 static_cast<TH1D*>(hist)->GetXaxis()->SetBinLabel(1, fCounters->GetAxis(targetDims[0])->GetBinLabel(selectedBins[targetDims[0]]));
916 else {
917 TObjString* label;
918 TIter nextLabel(fCounters->GetAxis(targetDims[0])->GetLabels());
919 while ((label = static_cast<TObjString*>(nextLabel())))
920 static_cast<TH1D*>(hist)->GetXaxis()->SetBinLabel((Int_t)label->GetUniqueID(), label->String().Data());
921 }
922
923 } else if (nTargetDim == 2) {
924
925 // project counters to TH2D (warning X and Y inverted in THnSparse::Projection(X,Y))
926 hist = fCounters->Projection(targetDims[1], targetDims[0]);
927
928 // reset bin labels in X axis lost when producing TH2D
929 if (selectedBins[targetDims[0]] >= 0)
930 static_cast<TH2D*>(hist)->GetXaxis()->SetBinLabel(1, fCounters->GetAxis(targetDims[0])->GetBinLabel(selectedBins[targetDims[0]]));
931 else {
932 TObjString* label;
933 TIter nextLabel(fCounters->GetAxis(targetDims[0])->GetLabels());
934 while ((label = static_cast<TObjString*>(nextLabel())))
935 static_cast<TH2D*>(hist)->GetXaxis()->SetBinLabel((Int_t)label->GetUniqueID(), label->String().Data());
936 }
937
938 // reset bin labels in Y axis lost when producing TH2D
939 if (selectedBins[targetDims[1]] >= 0)
940 static_cast<TH2D*>(hist)->GetYaxis()->SetBinLabel(1, fCounters->GetAxis(targetDims[1])->GetBinLabel(selectedBins[targetDims[1]]));
941 else {
942 TObjString* label;
943 TIter nextLabel(fCounters->GetAxis(targetDims[1])->GetLabels());
944 while ((label = static_cast<TObjString*>(nextLabel())))
945 static_cast<TH2D*>(hist)->GetYaxis()->SetBinLabel((Int_t)label->GetUniqueID(), label->String().Data());
946 }
947
948 } else {
949
950 // project counters to THnSparse (labels are not lost in that case)
951 hist = fCounters->Projection(nTargetDim, targetDims);
952
953 // reset bin labels in case only one bin has been selected
954 for (Int_t i=0; i<nTargetDim; i++) {
955 if (selectedBins[targetDims[i]] >= 0) {
956 TAxis* axis = static_cast<THnSparse*>(hist)->GetAxis(i);
957 axis->GetLabels()->Clear();
958 axis->SetBinLabel(1, fCounters->GetAxis(targetDims[i])->GetBinLabel(selectedBins[targetDims[i]]));
959 }
960 }
961
962 }
963
964 // clean memory
965 delete[] targetDims;
966 delete[] selectedBins;
967
968 return hist;
969}
970
971//-----------------------------------------------------------------------
972Int_t* AliCounterCollection::CheckConsistency(const AliCounterCollection* c)
973{
974 /// Consistency check of the two counter collections. To be consistent, both counters
975 /// must have the same rubrics with the same list of authorized key words if any.
976 /// Return the correspondence between the local rubric ordering and the one of the other counter,
977 /// or 0x0 in case of problem. It is the responsability of the user to delete the returned array.
978
979 if (!fCounters || !c->fCounters) {
980 AliError("counters are not initialized");
981 return 0x0;
982 }
983
984 // check if the number of rubrics is the same
985 Int_t nRubrics = fRubrics->GetSize();
986 if (c->fRubrics->GetSize() != nRubrics) {
987 AliError("both counters do not contain the same number of rubrics");
988 return 0x0;
989 }
990
991 Int_t* otherDims = new Int_t[nRubrics];
992
993 // loop over local rubrics
994 TObject* rubric1 = 0x0;
995 TIter nextRubric(fRubrics);
996 while ((rubric1 = nextRubric())) {
997
998 // find that rubric in the other counter
999 TObject* rubric2 = c->fRubrics->FindObject(rubric1->GetName());
1000 if (!rubric2) {
1001 AliError(Form("the other counter does not contain the rubric %s", rubric1->GetName()));
1002 delete[] otherDims;
1003 return 0x0;
1004 }
1005
1006 // check the list of authorized key words if any
1007 TObjArray* keyWords1 = dynamic_cast<TObjArray*>(rubric1);
1008 TObjArray* keyWords2 = dynamic_cast<TObjArray*>(rubric2);
1009 if (keyWords1 && keyWords2) {
1010
1011 // check if the number of key words is the same
1012 if (keyWords1->GetEntriesFast() != keyWords2->GetEntriesFast()) {
1013 AliError("that rubric does not contain the same number of authorized key words in both counters");
1014 delete[] otherDims;
1015 return 0x0;
1016 }
1017
1018 // loop over local key words
1019 TObjString* keyWord = 0x0;
1020 TIter nextKeyWord(keyWords1);
1021 while ((keyWord = static_cast<TObjString*>(nextKeyWord()))) {
1022
1023 // find that key word in the corresponding rubric of the other counter
1024 if (!keyWords2->FindObject(keyWord->String().Data())) {
1025 AliError(Form("rubric %s does not contain the key word %s in the other counter", rubric1->GetName(), keyWord->String().Data()));
1026 delete[] otherDims;
1027 return 0x0;
1028 }
1029
1030 }
1031
1032 } else if (keyWords1 || keyWords2) {
1033
1034 // that rubric has not been initialized the same way in both counter
1035 if (keyWords1) {
1036 AliError(Form("rubric %s of the other counter does not contain a list of authorized key words while this does", rubric1->GetName()));
1037 } else {
1038 AliError(Form("rubric %s of this counter does not contain a list of authorized key words while the other does", rubric1->GetName()));
1039 }
1040 delete[] otherDims;
1041 return 0x0;
1042
1043 }
1044
1045 // save the correspondence of rubric IDs in both counters
1046 otherDims[rubric1->GetUniqueID()] = rubric2->GetUniqueID();
1047
1048 }
1049
1050 return otherDims;
1051}
1052
1053//-----------------------------------------------------------------------
1054void AliCounterCollection::Add(const AliCounterCollection* counter)
1055{
1056 /// Add the given AliCounterCollections to this. They must have the
1057 /// same rubrics with the same list of authorized key words if any.
1058
1059 // check the consistency between the other counter and this and get the correspondences between rubric IDs.
1060 Int_t* otherDims = CheckConsistency(counter);
1061 if (!otherDims) return;
1062
1063 Int_t nRubrics = fCounters->GetNdimensions();
1064 Int_t* thisBins = new Int_t[nRubrics];
1065 Int_t* otherBins = new Int_t[nRubrics];
1066
1067 // loop over every filled bins inside the other counter
1068 for (Long64_t i = 0; i < counter->fCounters->GetNbins(); i++) {
1069
1070 // get the content of the bin
1071 Double_t value = counter->fCounters->GetBinContent(i, otherBins);
1072
1073 // convert "other" bin coordinates to "this" bin coordinates
1074 Bool_t ok = kTRUE;
1075 for (Int_t dim = 0; dim < nRubrics; dim++) {
1076 TString label = counter->fCounters->GetAxis(otherDims[dim])->GetBinLabel(otherBins[otherDims[dim]]);
1077 thisBins[dim] = FindBin(dim, label, kTRUE);
1078 if (thisBins[dim] < 0) {
1079 AliError("this counter is full, unable to add that key word");
1080 ok = kFALSE;
1081 break;
1082 }
1083 }
1084 if (!ok) continue;
1085
1086 // increment the corresponding local counter
1087 fCounters->AddBinContent(thisBins, value);
1088 }
1089
1090 // clean memory
1091 delete[] otherDims;
1092 delete[] thisBins;
1093 delete[] otherBins;
1094}
1095
1096//-----------------------------------------------------------------------
1097Long64_t AliCounterCollection::Merge(TCollection* list)
1098{
1099 /// Merge this with a list of AliCounterCollections. All AliCounterCollections provided
1100 /// must have the same rubrics with the same list of authorized key words if any.
1101
1102 if (!list || !fCounters) return 0;
1103 if (list->IsEmpty()) return (Long64_t)fCounters->GetEntries();
1104
1105 TIter next(list);
1106 const TObject* obj = 0x0;
1107 while ((obj = next())) {
1108
1109 // check that "obj" is an object of the class AliCounterCollection
1110 const AliCounterCollection* counter = dynamic_cast<const AliCounterCollection*>(obj);
1111 if (!counter) {
1112 AliError(Form("object named %s is not AliCounterCollection! Skipping it.", counter->GetName()));
1113 continue;
1114 }
1115
1116 // merge counter to this one
1117 Add(counter);
1118
1119 }
1120
1121 return (Long64_t)fCounters->GetEntries();
1122}
1123
1124//-----------------------------------------------------------------------
1125void AliCounterCollection::Sort(Option_t* opt, Bool_t asInt)
1126{
1127 /// Sort rubrics defined without a list of authorized key words or all rubrics if opt=="all".
1128 /// If asInt=kTRUE, key words are ordered as interger instead of alphabetically.
1129
1130 if (!fCounters) {
1131 AliError("counters are not initialized");
1132 return;
1133 }
1134
1135 Bool_t all = (!strcasecmp(opt, "all"));
1136
1137 Bool_t somethingToSort = kFALSE;
1138 Int_t nRubrics = fRubrics->GetSize();
1139 Bool_t* rubricsToSort = new Bool_t[nRubrics];
1140 memset(rubricsToSort, kFALSE, sizeof(Bool_t) * nRubrics);
1141
1142 // choose rubrics to sort
1143 TObject* rubric = 0x0;
1144 TIter nextRubric(fRubrics);
1145 while ((rubric = nextRubric())) {
1146
1147 if (all || dynamic_cast<TObjString*>(rubric)) {
1148
1149 // check if something to sort
1150 THashList* labels = fCounters->GetAxis((Int_t)rubric->GetUniqueID())->GetLabels();
1151 if (!labels || labels->GetSize() < 2) continue;
1152
1153 // select that rubric
1154 rubricsToSort[(Int_t)rubric->GetUniqueID()] = kTRUE;
1155 somethingToSort = kTRUE;
1156
1157 }
1158
1159 }
1160
1161 // sort selected rubrics if any
1162 if (somethingToSort) Sort(rubricsToSort, asInt);
1163
1164 // clean memory
1165 delete[] rubricsToSort;
1166}
1167
1168//-----------------------------------------------------------------------
1169void AliCounterCollection::SortRubric(TString rubric, Bool_t asInt)
1170{
1171 /// Sort only that rubric. If asInt=kTRUE, key words are ordered as interger instead of alphabetically.
1172
1173 if (!fCounters) {
1174 AliError("counters are not initialized");
1175 return;
1176 }
1177
1178 rubric.ToUpper();
1179
1180 // find the rubric to sort
1181 Int_t dim = FindDim(rubric);
1182 if (dim < 0) return;
1183
1184 // check if something to sort
1185 THashList* labels = fCounters->GetAxis(dim)->GetLabels();
1186 if (!labels || labels->GetSize() < 2) return;
1187
1188 // select that rubric
1189 Int_t nRubrics = fRubrics->GetSize();
1190 Bool_t* rubricsToSort = new Bool_t[nRubrics];
1191 memset(rubricsToSort, kFALSE, sizeof(Bool_t) * nRubrics);
1192 rubricsToSort[dim] = kTRUE;
1193
1194 // sort it
1195 Sort(rubricsToSort, asInt);
1196
1197 // clean memory
1198 delete[] rubricsToSort;
1199}
1200
1201//-----------------------------------------------------------------------
1202void AliCounterCollection::Sort(const Bool_t* rubricsToSort, Bool_t asInt)
1203{
1204 /// Sort labels (alphabetically or as integer) in each rubric flagged in "rubricsToSort".
1205
1206 // create a new counter
1207 THnSparse* oldCounters = fCounters;
1208 Int_t nRubrics = fRubrics->GetSize();
1209 fCounters = new THnSparseT<TArrayI>("hCounters", "hCounters", nRubrics, fRubricsSize->GetArray(), 0x0, 0x0);
1210 Int_t** newBins = new Int_t*[nRubrics];
1211
1212 // define the new axes
1213 for (Int_t i=0; i<nRubrics; i++) {
1214 TAxis* oldAxis = oldCounters->GetAxis(i);
1215 TAxis* newAxis = fCounters->GetAxis(i);
1216
1217 // set the name of the new axis
1218 newAxis->SetName(oldAxis->GetName());
1219
1220 // get old labels
1221 THashList* oldLabels = oldAxis->GetLabels();
1222 if (!oldLabels) continue;
1223
1224 // sort them if required
1225 if (rubricsToSort[i]) {
1226 if (asInt) { oldLabels = SortAsInt(oldLabels); }
1227 else { oldLabels->Sort(); }
1228 }
1229
1230 // set labels in the new axis and save the correspondence between new and old bins
1231 newBins[i] = new Int_t[oldLabels->GetSize()+1];
1232 TObjString* label = 0x0;
1233 Int_t bin = 1;
1234 TIter nextLabel(oldLabels);
1235 while ((label = static_cast<TObjString*>(nextLabel()))) {
1236 newAxis->SetBinLabel(bin, label->String().Data());
1237 newBins[i][(Int_t)label->GetUniqueID()] = bin;
1238 bin++;
1239 }
1240
1241 // clean memory
1242 if (rubricsToSort[i] && asInt) delete oldLabels;
1243 }
1244
1245 // fill the new counters
1246 Int_t* oldCoor = new Int_t[nRubrics];
1247 Int_t* newCoor = new Int_t[nRubrics];
1248 for (Long64_t i = 0; i < oldCounters->GetNbins(); i++) {
1249 Double_t value = oldCounters->GetBinContent(i, oldCoor);
1250 for (Int_t dim = 0; dim < nRubrics; dim++) newCoor[dim] = newBins[dim][oldCoor[dim]];
1251 fCounters->AddBinContent(newCoor, value);
1252 }
1253
1254 // clean memory
1255 for (Int_t i=0; i<nRubrics; i++) delete[] newBins[i];
1256 delete[] newBins;
1257 delete[] oldCoor;
1258 delete[] newCoor;
1259 delete oldCounters;
1260}
1261
1262//-----------------------------------------------------------------------
1263THashList* AliCounterCollection::SortAsInt(const THashList* labels)
1264{
1265 /// Return a list (not owner) of labels sorted assuming they are integers.
1266 /// It is the responsability of user to delete the returned list.
1267
1268 THashList* sortedLabels = new THashList(labels->GetSize());
1269 TIter nextSortedLabel(sortedLabels);
1270
1271 // loop over labels
1272 TObjString* label = 0x0;
1273 TIter nextLabel(labels);
1274 while ((label = static_cast<TObjString*>(nextLabel()))) {
1275
1276 // find where to add it
1277 TObjString* sortedLabel = 0x0;
1278 nextSortedLabel.Reset();
1279 while ((sortedLabel = static_cast<TObjString*>(nextSortedLabel())) &&
1280 (sortedLabel->String().Atoi() <= label->String().Atoi())) {}
1281
1282 // add it
1283 if (sortedLabel) sortedLabels->AddBefore(sortedLabel, label);
1284 else sortedLabels->AddLast(label);
1285 }
1286
1287 return sortedLabels;
1288}
1289