AliPhysics  9b6b435 (9b6b435)
AliEmcalList.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2016, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: R. Haake *
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 #include <iostream>
17 #include "TList.h"
18 #include "TH1.h"
19 #include "TH2.h"
20 #include "THnBase.h"
21 #include "AliLog.h"
22 #include "AliEmcalList.h"
23 
24 #ifdef WITH_ROOUNFOLD
25 #include "RooUnfoldResponse.h"
26 #endif
27 
28 ClassImp(AliEmcalList)
29 
30 //________________________________________________________________________
33  fUseScaling(kFALSE),
34  fNameXsec("fHistXsection"),
35  fNameNTrials("fHistTrials")
36 {
37 }
38 
40 //________________________________________________________________________
42 {
43  if(!hlist)
44  return 0;
45 
46  AliInfoStream() << "Scaled merging for list " <<hlist->GetName() << " is " << (fUseScaling ? "" : "not ") << "activated." << std::endl;
47  if(!fUseScaling)
48  return TList::Merge(hlist);
49 
50  // #### Retrieve xsection and ntrials from histograms in this list
51  // NOTE: they must be directly added to the AliEmcalList, not nested in sublists!
52  TH1* xsection = static_cast<TH1*>(FindObject(fNameXsec.Data()));
53  TH1* ntrials = static_cast<TH1*>(FindObject(fNameNTrials.Data()));
54 
55  // #### If xsection or ntrials are not given, go to fatal
56  if(!(xsection && ntrials))
57  AliFatal("Scaled merging is active but AliEmcalList does not contain fHistXsection or fHistTrials. Do not activate scaling for those lists or include these histograms.");
58 
59  // #### Do the scaling only on the last level when we mix several pt hard bins
60  // This is easy to find out checking the std histos in hlist
61  Bool_t isLastLevel = IsLastMergeLevel(hlist);
62 
63  // #### On last level, do the scaling
64  if(isLastLevel)
65  {
66  AliInfoStream() << "===== LAST LEVEL OF MERGING =====" << std::endl;
67 
68  // Scale all histograms in this list
69  Double_t scalingFactor = GetScalingFactor(xsection, ntrials);
70  ScaleAllHistograms(this, scalingFactor);
71 
72  // Scale all histograms in the lists to be merged
73  TIter listIterator(hlist);
74  while (AliEmcalList* tmpList = static_cast<AliEmcalList*>(listIterator()))
75  {
76  xsection = static_cast<TH1*>(tmpList->FindObject(fNameXsec.Data()));
77  ntrials = static_cast<TH1*>(tmpList->FindObject(fNameNTrials.Data()));
78  scalingFactor = GetScalingFactor(xsection, ntrials);
79  ScaleAllHistograms(tmpList, scalingFactor);
80  }
81  }
82 
83  AliInfoStream() << "Merge() done." << std::endl;
84 
85  TList::Merge(hlist);
86  return hlist->GetEntries() + 1;
87 }
88 
90 //________________________________________________________________________
92 {
93  TIter listIterator(hlist);
94  while (TObject* listObject = listIterator())
95  {
96  // Recurse into nested all lists
97  TCollection* sublist = dynamic_cast<TCollection*>(listObject);
98  if(sublist)
99  {
100  ScaleAllHistograms(sublist, scalingFactor);
101  continue;
102  }
103 
104  // Otherwise, scale object if scaling is supported
105  if (!IsScalingSupported(listObject)){
106  AliInfoStream() << "Scaling of objects of type " << listObject->IsA()->GetName() << " unsupported - object " << listObject->GetName() << " will not be scaled" << std::endl;
107  continue;
108  }
109 
110  // Don't scale profiles and histograms used for scaling
111  TString histogram_class (listObject->ClassName());
112  TString histogram_name (listObject->GetName());
113  if (histogram_name.Contains("fHistXsection") || histogram_name.Contains("fHistTrials") || histogram_name.Contains("fHistEvents"))
114  {
115  AliInfoStream() << "Histogram " << listObject->GetName() << " will not be scaled, because a scaling histogram" << std::endl;
116  continue;
117  }
118  if (histogram_class.Contains("TProfile"))
119  {
120  AliInfoStream() << "Histogram " << listObject->GetName() << " will not be scaled, because it is a TProfile" << std::endl;
121  continue;
122  }
123 
124  TH1 *histogram = dynamic_cast<TH1 *>(listObject);
125  if(histogram) {
126  // Handle TH1/TH2/TH3
127  histogram->Sumw2();
128  histogram->Scale(scalingFactor);
129  } else if(listObject->InheritsFrom(THnBase::Class())){
130  // Handle THn
131  THnBase *histogramND = dynamic_cast<THnBase *>(listObject);
132  histogramND->Sumw2();
133  histogramND->Scale(scalingFactor);
134  }
135 #ifdef WITH_ROOUNFOLD
136  else if(listObject->InheritsFrom(RooUnfoldResponse::Class())){
137  RooUnfoldResponse *response = dynamic_cast<RooUnfoldResponse *>(listObject);
138  if(auto measured = response->Hmeasured()) measured->Scale(scalingFactor);
139  if(auto fakes = response->Hfakes()) fakes->Scale(scalingFactor);
140  if(auto truth = response->Htruth()) truth->Scale(scalingFactor);
141  if(auto responseND = response->Hresponse()) responseND->Scale(scalingFactor);
142  }
143 #endif
144  AliInfoStream() << "Histogram " << listObject->GetName() << " (" << histogram_class.Data() << ") was scaled..." << std::endl;
145 
146  }
147 }
148 
150 //________________________________________________________________________
151 Double_t AliEmcalList::GetScalingFactor(const TH1* xsection, const TH1* ntrials) const
152 {
153  Int_t binNumber = GetFilledBinNumber(xsection);
154  if(!binNumber)
155  {
156  AliInfoStream() << "List already scaled or scaling invalid. Scaling factor = 1." << std::endl;
157  return 1.0;
158  }
159 
160  Double_t valNTRIALS = ntrials->GetBinContent(binNumber);
161  Double_t valXSEC = xsection->GetBinContent(binNumber);
162  Double_t scalingFactor = 0;
163  if(valNTRIALS)
164  scalingFactor = valXSEC/valNTRIALS;
165 
166  AliInfoStream() << "## Bin " << binNumber << ": trials=" << valNTRIALS << ", xsec=" << valXSEC << " -> scaling=" << scalingFactor << std::endl;
167  return scalingFactor;
168 }
169 
172 //________________________________________________________________________
174 {
175  // Get the pt hard bin number that is filled in this object
176  TH1* xsection = static_cast<TH1*>(FindObject("fHistXsection"));
177  Int_t binNumberThis = GetFilledBinNumber(xsection);
178 
179  TIter listIterator(collection);
180  while (AliEmcalList* tmpList = static_cast<AliEmcalList*>(listIterator()))
181  {
182  // For every list in the collection, test if they are from different pt hard bins
183  xsection = static_cast<TH1*>(tmpList->FindObject("fHistXsection"));
184  if(binNumberThis != GetFilledBinNumber(xsection))
185  return kTRUE;
186  }
187  return kFALSE;
188 }
189 
199 //________________________________________________________________________
200 bool AliEmcalList::IsScalingSupported(const TObject *scaleobject) const {
201  if(scaleobject->InheritsFrom(TH1::Class()) || scaleobject->InheritsFrom(THnBase::Class())) return true;
202 #ifdef WITH_ROOUNFOLD
203  if(scaleobject->InheritsFrom(RooUnfoldResponse::Class())) return true;
204 #endif
205  return false;
206 }
207 
211 //________________________________________________________________________
213 {
214  AliInfoStream() << hist->GetName() << ": nbinsX=" << hist->GetNbinsX() << std::endl;;
215 
216  Int_t binFound = 0;
217  for(Int_t i=1; i<=hist->GetNbinsX(); i++)
218  {
219  AliInfoStream() << hist->GetName() << ": bin=" << i << ", val=" << hist->GetBinContent(i) << std::endl;
220  if(hist->GetBinContent(i))
221  {
222  if(!binFound)
223  binFound = i;
224  else
225  {
226  return 0; // more than one bin filled (e.g. when merging is partly done)
227  }
228  }
229  }
230 
231  if(!binFound)
232  AliErrorStream() << "No bin filled in scaling histogram." << std::endl;
233 
234  // 0 if no bin found or more than one filled, otherwise returns bin number
235  return binFound;
236 }
double Double_t
Definition: External.C:58
Bool_t fUseScaling
if true, scaling will be done. if false AliEmcalList simplifies to TList
Definition: AliEmcalList.h:44
long long Long64_t
Definition: External.C:43
Double_t GetScalingFactor(const TH1 *xsection, const TH1 *ntrials) const
Helper function scaling factor.
void ScaleAllHistograms(TCollection *hlist, Double_t scalingFactor)
Function that does the scaling of all histograms in hlist recursively.
int Int_t
Definition: External.C:63
TString fNameXsec
Name of the cross section histogram.
Definition: AliEmcalList.h:45
Bool_t IsScalingSupported(const TObject *scaleobject) const
Helper function checking whether type is supported for scaling.
Long64_t Merge(TCollection *hlist)
Overridden Merge function.
TH1 * Merge(const TH1 *cen, const TH1 *fwd, Double_t &xlow, Double_t &xhigh)
TObject * FindObject(int bin, const char *nameH, const TList *lst, Bool_t normPerEvent=kTRUE)
Enhanced TList-derived class that implements correct merging for pt_hard binned production.
Definition: AliEmcalList.h:25
Bool_t IsLastMergeLevel(const TCollection *collection) const
Helper function to determine whether we are in last merge step.
TString fNameNTrials
Name of the histogram with the number of trials.
Definition: AliEmcalList.h:46
bool Bool_t
Definition: External.C:53
Int_t GetFilledBinNumber(const TH1 *hist) const
Helper function that returns the bin in a TH1 that is filled.
Definition: External.C:196