68011dfa395be849ee79dd1c0cdc9b7b334d7697
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / AliAODForwardMult.cxx
1 //
2 // Class that contains the forward multiplicity data per event 
3 //
4 // This class contains a histogram of 
5 // @f[
6 //   \frac{d^2N_{ch}}{d\eta d\phi}\quad,
7 // @f]
8 // as well as a trigger mask for each analysed event.  
9 // 
10 // The eta acceptance of the event is stored in the underflow bins of
11 // the histogram.  So to build the final histogram, one needs to
12 // correct for this acceptance (properly weighted by the events), and
13 // the vertex efficiency.  This simply boils down to defining a 2D
14 // histogram and summing the event histograms in that histogram.  One
15 // should of course also do proper book-keeping of the accepted event.
16 //
17 #include "AliAODForwardMult.h"
18 #include <TBrowser.h>
19 #include <iostream>
20 #include <TMath.h>
21 #include <TObjString.h>
22
23 ClassImp(AliAODForwardMult)
24 #ifdef DOXY_INPUT
25 ; // For Emacs 
26 #endif
27
28 //____________________________________________________________________
29 const Float_t AliAODForwardMult::fgkInvalidIpZ = 1e6;
30
31 //____________________________________________________________________
32 AliAODForwardMult::AliAODForwardMult()
33   : fIsMC(false),
34     fHist(),
35     fTriggers(0),
36     fIpZ(fgkInvalidIpZ), 
37     fCentrality(-1)
38 {
39   // 
40   // Constructor 
41   // 
42 }
43
44 //____________________________________________________________________
45 AliAODForwardMult::AliAODForwardMult(Bool_t isMC) 
46   : fIsMC(isMC),
47     fHist("forwardMult", "d^{2}N_{ch}/d#etad#varphi in the forward regions", 
48           200, -4, 6, 20, 0, 2*TMath::Pi()),
49     fTriggers(0),
50     fIpZ(fgkInvalidIpZ), 
51     fCentrality(-1)
52 {
53   // 
54   // Constructor 
55   // 
56   // Parameters: 
57   //  isMC   If set to true this is for MC data (effects branch name)
58   // 
59   fHist.SetXTitle("#eta");
60   fHist.SetYTitle("#varphi [radians]");
61   fHist.SetZTitle("#frac{d^{2}N_{ch}}{d#etad#varphi}");
62   fHist.SetDirectory(0);
63   fHist.Sumw2();
64 }
65
66 //____________________________________________________________________
67 void
68 AliAODForwardMult::Init(const TAxis& etaAxis)
69 {
70   // Initialize the histogram with an eta axis 
71   // 
72   // Parameters: 
73   //   etaAxis       Eta axis to use 
74   // 
75   fHist.SetBins(etaAxis.GetNbins(), etaAxis.GetXmin(), etaAxis.GetXmax(), 
76                 20, 0, 2*TMath::Pi());
77 }
78
79 //____________________________________________________________________
80 void
81 AliAODForwardMult::Clear(Option_t* option)
82 {
83   // Clear (or reset) internal values 
84   // 
85   // Parameters: 
86   //  option   Passed to TH1::Reset 
87   // 
88   fHist.Reset(option);
89   fTriggers = 0;
90   fIpZ      = fgkInvalidIpZ;
91 }
92 //____________________________________________________________________
93 void
94 AliAODForwardMult::SetSNN(UShort_t snn)
95 {
96   // set the center of mass energy per nucleon pair (GeV). 
97   // This is stored in bin (0,0) of the histogram 
98   // 
99   // Parameters: 
100   //   sNN   Center of mass energy per nuclean 
101   fHist.SetBinContent(0,0,snn);
102 }
103 //____________________________________________________________________
104 void
105 AliAODForwardMult::SetSystem(UShort_t sys)
106 {
107   // set the center of mass energy per nucleon pair (GeV). 
108   // This is stored in bin (N+1,0) of the histogram 
109   // 
110   // Parameters: 
111   //   sys   Collision system number 
112   fHist.SetBinContent(fHist.GetNbinsX()+1,0,sys);
113 }
114
115 //____________________________________________________________________
116 Bool_t
117 AliAODForwardMult::HasIpZ() const
118 {
119   // Check if we have valid z coordinate of the interaction point 
120   // 
121   // Return:
122   //   true if the z coordinate of the interaction point is valid 
123   // 
124   return TMath::Abs(fIpZ - fgkInvalidIpZ) > 1;
125 }
126
127 //____________________________________________________________________
128 UShort_t
129 AliAODForwardMult::GetSNN() const
130 {
131   // set the center of mass energy per nucleon pair (GeV). 
132   // This is stored in bin (0,0) of the histogram 
133   // 
134   // Parameters: 
135   //   sNN   Center of mass energy per nuclean 
136   return UShort_t(fHist.GetBinContent(0,0));
137 }
138
139 //____________________________________________________________________
140 UShort_t
141 AliAODForwardMult::GetSystem() const
142 {
143   // set the center of mass energy per nucleon pair (GeV). 
144   // This is stored in bin (N+1,0) of the histogram 
145   // 
146   // Parameters: 
147   //   sNN   Center of mass energy per nuclean 
148   return UShort_t(fHist.GetBinContent(fHist.GetNbinsX()+1,0));
149 }
150
151 //____________________________________________________________________
152 void
153 AliAODForwardMult::Browse(TBrowser* b)
154 {
155   // Browse this object 
156   // 
157   // Parameters: 
158   //   b   Browser to use 
159   static TObjString ipz;
160   static TObjString trg;
161   static TObjString cnt;
162   ipz = Form("ip_z=%fcm", fIpZ);
163   trg = GetTriggerString(fTriggers);
164   cnt = Form("%+6.1f%%", fCentrality);
165   b->Add(&fHist);
166   b->Add(&ipz);
167   b->Add(&trg);
168   b->Add(&cnt);
169 }
170
171 //____________________________________________________________________
172 const Char_t*
173 AliAODForwardMult::GetTriggerString(UInt_t mask)
174 {
175   // Get a string that describes the triggers 
176   // 
177   // Parameters: 
178   //   mask  Bit pattern of triggers 
179   // Return: 
180   //   Character string representation of mask 
181   static TString trg;
182   trg = "";
183   if ((mask & kInel)        != 0x0) trg.Append("INEL ");
184   if ((mask & kInelGt0)     != 0x0) trg.Append("INEL>0 ");
185   if ((mask & kNSD)         != 0x0) trg.Append("NSD ");
186   if ((mask & kA)           != 0x0) trg.Append("A ");
187   if ((mask & kB)           != 0x0) trg.Append("B ");
188   if ((mask & kC)           != 0x0) trg.Append("C ");
189   if ((mask & kE)           != 0x0) trg.Append("E ");
190   if ((mask & kMCNSD)       != 0x0) trg.Append("MCNSD ");
191   return trg.Data();
192 }
193   
194 //____________________________________________________________________
195 TH1I*
196 AliAODForwardMult::MakeTriggerHistogram(const char* name) 
197 {
198   // 
199   // Make a histogram to record triggers in. 
200   //
201   // The bins defined by the trigger enumeration in this class.  One
202   // can use this enumeration to retrieve the number of triggers for
203   // each class.
204   // 
205   // Parameters:
206   //    name Name of the histogram 
207   // 
208   // Return:
209   //    Newly allocated histogram 
210   //
211   TH1I* ret = new TH1I(name, "Triggers", kAccepted+1, -.5, kAccepted+.5);
212   ret->SetYTitle("Events");
213   ret->SetFillColor(kRed+1);
214   ret->SetFillStyle(3001);
215   ret->GetXaxis()->SetBinLabel(kBinAll,         "All events");
216   ret->GetXaxis()->SetBinLabel(kBinB,           "w/B trigger");
217   ret->GetXaxis()->SetBinLabel(kBinA,           "w/A trigger");
218   ret->GetXaxis()->SetBinLabel(kBinC,           "w/C trigger");
219   ret->GetXaxis()->SetBinLabel(kBinE,           "w/E trigger");
220   ret->GetXaxis()->SetBinLabel(kBinInel,        "Minimum Bias");
221   ret->GetXaxis()->SetBinLabel(kBinInelGt0,     "INEL>0");
222   ret->GetXaxis()->SetBinLabel(kBinNSD,         "NSD");
223   ret->GetXaxis()->SetBinLabel(kBinMCNSD,       "NSD (MC truth)");
224   ret->GetXaxis()->SetBinLabel(kBinPileUp,      "w/Pileup");
225   ret->GetXaxis()->SetBinLabel(kBinOffline,     "w/Offline");
226   ret->GetXaxis()->SetBinLabel(kWithVertex,     "w/Vertex");
227   ret->GetXaxis()->SetBinLabel(kWithTrigger,    "w/Selected trigger");
228   ret->GetXaxis()->SetBinLabel(kAccepted,       "Accepted by cut");
229   ret->GetXaxis()->SetNdivisions(kAccepted, false);
230   ret->SetStats(0);
231
232   return ret;
233 }
234
235 //____________________________________________________________________
236 Bool_t
237 AliAODForwardMult::CheckEvent(Int_t    triggerMask,
238                               Double_t vzMin, Double_t vzMax,
239                               UShort_t cMin,  UShort_t cMax, 
240                               TH1*     hist) const
241 {
242   // 
243   // Check if event meets the passses requirements.   
244   //
245   // It returns true if @e all of the following is true 
246   //
247   // - The trigger is within the bit mask passed.
248   // - The vertex is within the specified limits. 
249   // - The centrality is within the specified limits, or if lower
250   //   limit is equal to or larger than the upper limit.
251   // 
252   // If a histogram is passed in the last parameter, then that
253   // histogram is filled with the trigger bits. 
254   // 
255   // Parameters:
256   //    triggerMask  Trigger mask
257   //    vzMin        Minimum @f$ v_z@f$ (in centimeters)
258   //    vzMax        Maximum @f$ v_z@f$ (in centimeters) 
259   //    cMin         Minimum centrality (in percent)
260   //    cMax         Maximum centrality (in percent)
261   //    hist         Histogram to fill 
262   // 
263   // Return:
264   //    @c true if the event meets the requirements 
265   //
266   if (cMin < cMax && (cMin > fCentrality || cMax <= fCentrality)) return false;
267
268   if (hist) { 
269     hist->AddBinContent(kBinAll);
270     if (IsTriggerBits(kB|triggerMask))  hist->AddBinContent(kBinB);
271     if (IsTriggerBits(kA|triggerMask))  hist->AddBinContent(kBinA);
272     if (IsTriggerBits(kC|triggerMask))  hist->AddBinContent(kBinC);
273     if (IsTriggerBits(kE|triggerMask))  hist->AddBinContent(kBinE);
274     if (IsTriggerBits(kB|kInel))        hist->AddBinContent(kBinInel);
275     if (IsTriggerBits(kB|kInelGt0))     hist->AddBinContent(kBinInelGt0);
276     if (IsTriggerBits(kB|kNSD))         hist->AddBinContent(kBinNSD);
277     if (IsTriggerBits(kPileUp))         hist->AddBinContent(kBinPileUp);
278     if (IsTriggerBits(kMCNSD))          hist->AddBinContent(kBinMCNSD);
279     if (IsTriggerBits(kOffline))        hist->AddBinContent(kBinOffline);
280   }
281   // Check if we have an event of interest. 
282   if (!IsTriggerBits(triggerMask|kB)) return false;
283   
284   // Check for pileup
285   if (IsTriggerBits(kPileUp)) return false;
286   if (hist) hist->AddBinContent(kWithTrigger);
287   
288   // Check that we have a valid vertex
289   if (!HasIpZ()) return false;
290   if (hist) hist->AddBinContent(kWithVertex);
291
292   // Check that vertex is within cuts 
293   if (!InRange(vzMin, vzMax)) return false;
294   if (hist) hist->AddBinContent(kAccepted);
295   
296   return true;
297 }
298
299 //____________________________________________________________________
300 void
301 AliAODForwardMult::Print(Option_t* option) const
302 {
303   // Print this object 
304   // 
305   // Parameters: 
306   //  option   Passed to TH1::Print 
307   fHist.Print(option);
308   UShort_t sys = GetSystem();
309   TString  str = "unknown";
310   switch (sys) { 
311   case 1:  str = "pp"; break;
312   case 2:  str = "PbPb"; break;
313   }
314   std::cout << "Ipz:         " << fIpZ << "cm " << (HasIpZ() ? "" : "in") 
315             << "valid\n"
316             << "Triggers:    " << GetTriggerString(fTriggers) 
317             << "sNN:         " << GetSNN() << "GeV\n" 
318             << "System:      " << str 
319             << "Centrality:  " << fCentrality << "%" 
320             << std::endl;
321 }
322
323 //____________________________________________________________________
324 //
325 // EOF
326 //