2 // See implementation or Doxygen comments for more information
4 #ifndef ALIAODFORWARDMULT_H
5 #define ALIAODFORWARDMULT_H
10 * Class that contains the forward multiplicity data per event
12 * This class contains a histogram of
14 * \frac{d^2N_{ch}}{d\eta d\phi}\quad,
16 * as well as a trigger mask for each analysed event.
18 * The eta acceptance of the event is stored in the underflow bins of
19 * the histogram. So to build the final histogram, one needs to
20 * correct for this acceptance (properly weighted by the events), and
21 * the vertex efficiency. This simply boils down to defining a 2D
22 * histogram and summing the event histograms in that histogram. One
23 * should of course also do proper book-keeping of the accepted event.
28 * TFile* file = TFile::Open("AliAODs.root","READ");
29 * TTree* tree = static_cast<TTree*>(file->Get("aodTree"));
35 * TH2D* sum = 0; // Summed hist
36 * TTree* tree = GetAODTree(); // AOD tree
37 * AliAODForwardMult* mult = 0; // AOD object
38 * Int_t nTriggered = 0; // # of triggered ev.
39 * Int_t nWithVertex= 0; // # of ev. w/vertex
40 * Int_t nAccepted = 0; // # of ev. used
41 * Int_t nAvailable = tree->GetEntries(); // How many entries
42 * Float_t vzLow = -10; // Lower ip cut
43 * Float_t vzHigh = 10; // Upper ip cut
44 * Int_t mask = AliAODForwardMult::kInel;// Trigger mask
45 * tree->SetBranchAddress("forward", &forward); // Set the address
47 * for (int i = 0; i < nAvailable; i++) {
48 * // Create sum histogram on first event - to match binning to input
49 * if (!sum) sum = static_cast<TH2D*>(mult->Clone("d2ndetadphi"));
53 * // Other trigger/event requirements could be defined
54 * if (!mult->IsTriggerBits(mask)) continue;
57 * // Check if we have vertex
58 * if (!mult->HasIpZ()) continue;
61 * // Select vertex range (in centimeters)
62 * if (!mult->InRange(vzLow, vzHigh) continue;
65 * // Add contribution from this event
66 * sum->Add(&(mult->GetHistogram()));
69 * // Get acceptance normalisation from underflow bins
70 * TH1D* norm = sum->Projection("norm", 0, 1, "");
71 * // Project onto eta axis - _ignoring_underflow_bins_!
72 * TH1D* dndeta = sum->Projection("dndeta", 1, -1, "e");
73 * // Normalize to the acceptance
74 * dndeta->Divide(norm);
75 * // Scale by the vertex efficiency
76 * dndeta->Scale(Double_t(nWithVertex)/nTriggered, "width");
77 * // And draw the result
82 * The above code will draw the final @f$ dN_{ch}/d\eta@f$ for the
83 * selected event class and vertex range
85 * The histogram can be used as input for other kinds of analysis too,
86 * like flow, event-plane, centrality, and so on.
88 * @ingroup pwg2_forward
90 class AliAODForwardMult : public TObject
94 * Bits of the trigger pattern
97 /** In-elastic collision */
99 /** In-elastic collision with at least one SPD tracklet */
101 /** Non-single diffractive collision */
103 /** Empty bunch crossing */
105 /** A-side trigger */
107 /** B(arrel) trigger */
109 /** C-side trigger */
113 /** pileup from SPD */
115 /** true NSD from MC */
120 * Default constructor
122 * Used by ROOT I/O sub-system - do not use
128 * @param isMC Whether this was from MC or not
130 AliAODForwardMult(Bool_t isMC);
134 virtual ~AliAODForwardMult() {} // Destructor
138 * @param etaAxis Pseudo-rapidity axis
140 void Init(const TAxis& etaAxis);
142 * Get the @f$ d^2N_{ch}/d\eta d\phi@f$ histogram,
144 * @return @f$ d^2N_{ch}/d\eta d\phi@f$ histogram,
146 const TH2D& GetHistogram() const { return fHist; } // Get histogram
148 * Get the @f$ d^2N_{ch}/d\eta d\phi@f$ histogram,
150 * @return @f$ d^2N_{ch}/d\eta d\phi@f$ histogram,
152 TH2D& GetHistogram() { return fHist; } // Get histogram
154 * Get the trigger mask
156 * @return Trigger mask
158 UInt_t GetTriggerMask() const { return fTriggers; } // Get triggers
160 * Set the trigger mask
162 * @param trg Trigger mask
164 void SetTriggerMask(UInt_t trg) { fTriggers = trg; } // Set triggers
166 * Set bit(s) in trigger mask
168 * @param bits bit(s) to set
170 void SetTriggerBits(UInt_t bits) { fTriggers |= bits; } // Set trigger bits
172 * Check if bit(s) are set in the trigger mask
174 * @param bits Bits to test for
178 Bool_t IsTriggerBits(UInt_t bits) const;
180 * Whether we have any trigger bits
182 Bool_t HasTrigger() const { return fTriggers != 0; } // Check for triggers
186 * @param option Passed on to TH2::Reset verbatim
188 void Clear(Option_t* option="");
194 void Browse(TBrowser* b);
198 * @return Always true
200 Bool_t IsFolder() const { return kTRUE; } // Always true
204 * @param option Passed verbatim to TH2::Print
206 void Print(Option_t* option="") const;
208 * Set the z coordinate of the interaction point
210 * @param ipZ Interaction point z coordinate
212 void SetIpZ(Float_t ipZ) { fIpZ = ipZ; } // Set Ip's Z coordinate
214 * Set the center of mass energy per nucleon-pair. This is stored
215 * in the (0,0) of the histogram
217 * @param sNN Center of mass energy per nucleon pair (GeV)
219 void SetSNN(UShort_t sNN);
221 * Get the collision system number
226 * @param sys Collision system number
228 void SetSystem(UShort_t sys);
230 * Set the event centrality
232 * @param c Centrality
234 void SetCentrality(Float_t c) { fCentrality = c; }
236 * Set the z coordinate of the interaction point
238 * @return Interaction point z coordinate
240 Float_t GetIpZ() const { return fIpZ; } // Get Ip's Z coordinate
242 * Check if we have a valid z coordinate of the interaction point
244 * @return True if we have a valid interaction point z coordinate
246 Bool_t HasIpZ() const;
248 * Get the center of mass energy per nucleon pair (GeV)
250 * @return Center of mass energy per nucleon pair (GeV)
252 UShort_t GetSNN() const;
254 * Get the collision system number
259 * @return Collision system number
261 UShort_t GetSystem() const;
263 * Check if the z coordinate of the interaction point is within the
264 * given limits. Note that the convention used corresponds to the
265 * convention used in ROOTs TAxis.
267 * @param low Lower cut (inclusive)
268 * @param high Upper cut (exclusive)
270 * @return true if @f$ low \ge ipz < high@f$
272 Bool_t InRange(Float_t low, Float_t high) const;
274 * Get the event centrality
279 Float_t GetCentrality() const { return fCentrality; }
281 * Check if we have a valid centrality
286 Bool_t HasCentrality() const { return !(fCentrality < 0); }
289 * Get the name of the object
291 * @return Name of object
293 const Char_t* GetName() const { return (fIsMC ? "ForwardMC" : "Forward"); }
295 * Get a string correspondig to the trigger mask
297 * @param mask Trigger mask
299 * @return Static string (copy before use)
301 static const Char_t* GetTriggerString(UInt_t mask);
303 Bool_t fIsMC; // Whether this is from MC
304 TH2D fHist; // Histogram of d^2N_{ch}/(deta dphi) for this event
305 UInt_t fTriggers; // Trigger bit mask
306 Float_t fIpZ; // Z coordinate of the interaction point
307 Float_t fCentrality; // Event centrality
309 static const Float_t fgkInvalidIpZ; // Invalid IpZ value
310 ClassDef(AliAODForwardMult,2); // AOD forward multiplicity
313 //____________________________________________________________________
315 AliAODForwardMult::InRange(Float_t low, Float_t high) const
317 return HasIpZ() && fIpZ >= low && fIpZ < high;
320 //____________________________________________________________________
322 AliAODForwardMult::IsTriggerBits(UInt_t bits) const
324 return HasTrigger() && ((fTriggers & bits) != 0);