96cad837 |
1 | /************************************************************************** |
2 | * Author: Boris Hippolyte. * |
3 | * Contributors are mentioned in the code where appropriate. * |
4 | * * |
5 | * Permission to use, copy, modify and distribute this software and its * |
6 | * documentation strictly for non-commercial purposes is hereby granted * |
7 | * without fee, provided that the above copyright notice appears in all * |
8 | * copies and that both the copyright notice and this permission notice * |
9 | * appear in the supporting documentation. The authors make no claims * |
10 | * about the suitability of this software for any purpose. It is * |
11 | * provided "as is" without express or implied warranty. * |
12 | **************************************************************************/ |
13 | |
14 | //----------------------------------------------------------------- |
c8eada4f |
15 | // AliAnalysisTaskCheckV0 class |
16 | // This task is for QAing the V0s from ESD/AOD |
96cad837 |
17 | // Origin: B.H. Nov2007, hippolyt@in2p3.fr |
18 | //----------------------------------------------------------------- |
19 | #include "TChain.h" |
20 | #include "TTree.h" |
21 | #include "TList.h" |
22 | #include "TH1F.h" |
23 | #include "TCanvas.h" |
c8eada4f |
24 | #include "TLegend.h" |
96cad837 |
25 | |
c8eada4f |
26 | #include "AliAnalysisTaskSE.h" |
96cad837 |
27 | #include "AliAnalysisManager.h" |
28 | |
29 | #include "AliESDEvent.h" |
30 | #include "AliESDInputHandler.h" |
c8eada4f |
31 | #include "AliAODEvent.h" |
32 | #include "AliAODInputHandler.h" |
96cad837 |
33 | |
34 | #include "AliESDv0.h" |
35 | |
c8eada4f |
36 | #include "AliAnalysisTaskCheckV0.h" |
96cad837 |
37 | |
c8eada4f |
38 | ClassImp(AliAnalysisTaskCheckV0) |
96cad837 |
39 | |
40 | //________________________________________________________________________ |
c8eada4f |
41 | AliAnalysisTaskCheckV0::AliAnalysisTaskCheckV0(const char *name) |
42 | : AliAnalysisTaskSE(name), fESD(0), fAOD(0), |
43 | fAnalysisType("ESD"), fCollidingSystems(0), fListHist(), |
96cad837 |
44 | fHistTrackMultiplicity(0), fHistV0Multiplicity(0), fHistV0OnFlyStatus(0), |
45 | fHistV0MultiplicityOff(0), fHistV0Chi2Off(0), |
46 | fHistDcaV0DaughtersOff(0), fHistV0CosineOfPointingAngleOff(0), |
47 | fHistV0RadiusOff(0),fHistDcaV0ToPrimVertexOff(0), |
48 | fHistDcaPosToPrimVertexOff(0),fHistDcaNegToPrimVertexOff(0), |
49 | fHistMassK0Off(0),fHistMassLambdaOff(0),fHistMassAntiLambdaOff(0), |
50 | fHistV0MultiplicityOn(0), fHistV0Chi2On(0), |
51 | fHistDcaV0DaughtersOn(0), fHistV0CosineOfPointingAngleOn(0), |
52 | fHistV0RadiusOn(0),fHistDcaV0ToPrimVertexOn(0), |
53 | fHistDcaPosToPrimVertexOn(0),fHistDcaNegToPrimVertexOn(0), |
54 | fHistMassK0On(0),fHistMassLambdaOn(0),fHistMassAntiLambdaOn(0) |
55 | { |
56 | // Constructor |
57 | |
58 | // Define input and output slots here |
59 | // Input slot #0 works with a TChain |
60 | DefineInput(0, TChain::Class()); |
61 | // Output slot #0 writes into a TH1 container |
62 | DefineOutput(0, TH1F::Class()); |
63 | // Output slot #1 writes into a TList container |
64 | DefineOutput(1, TList::Class()); |
65 | } |
66 | |
67 | //________________________________________________________________________ |
c8eada4f |
68 | void AliAnalysisTaskCheckV0::ConnectInputData(Option_t *) |
96cad837 |
69 | { |
70 | // Connect ESD or AOD here |
71 | // Called once |
72 | |
73 | TTree* tree = dynamic_cast<TTree*> (GetInputData(0)); |
74 | if (!tree) { |
75 | Printf("ERROR: Could not read chain from input slot 0"); |
76 | } else { |
c8eada4f |
77 | // tree->Print(); |
78 | if(fAnalysisType == "ESD") { |
79 | tree->SetBranchStatus("*", kFALSE); |
80 | tree->SetBranchStatus("Tracks.*", kTRUE); |
81 | tree->SetBranchStatus("V0s.*", kTRUE); |
82 | // tree->SetBranchStatus("fTracks.*", kTRUE); |
83 | // tree->SetBranchStatus("fV0s.*", kTRUE); |
84 | |
85 | AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()); |
86 | |
87 | if (!esdH) { |
88 | Printf("ERROR: Could not get ESDInputHandler"); |
89 | } else |
90 | fESD = esdH->GetEvent(); |
91 | } |
92 | else if(fAnalysisType == "AOD") { |
93 | AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()); |
94 | |
95 | if (!aodH) { |
96 | Printf("ERROR: Could not get AODInputHandler"); |
97 | } else |
98 | fAOD = aodH->GetEvent(); |
99 | } |
100 | else |
101 | Printf("Wrong analysis type: Only ESD and AOD types are allowed! (actually only ESD types is supported for the moment!)"); |
96cad837 |
102 | } |
c8eada4f |
103 | |
96cad837 |
104 | } |
105 | |
106 | //________________________________________________________________________ |
c8eada4f |
107 | void AliAnalysisTaskCheckV0::CreateOutputObjects() |
96cad837 |
108 | { |
109 | // Create histograms |
110 | // Called once |
111 | |
112 | // Distinguish Track and V0 Multiplicity ! |
113 | |
114 | fListHist = new TList(); |
115 | |
116 | if (!fHistTrackMultiplicity) { |
c8eada4f |
117 | if (fCollidingSystems) |
118 | fHistTrackMultiplicity = new TH1F("fHistTrackMultiplicity", "Multiplicity distribution;Number of tracks;Events", 200, 0, 40000); |
119 | else |
120 | fHistTrackMultiplicity = new TH1F("fHistTrackMultiplicity", "Multiplicity distribution;Number of tracks;Events", 250, 0, 250); |
96cad837 |
121 | fListHist->Add(fHistTrackMultiplicity); |
122 | } |
123 | if (!fHistV0Multiplicity) { |
c8eada4f |
124 | if (fCollidingSystems) |
125 | fHistV0Multiplicity = new TH1F("fHistV0Multiplicity", "Multiplicity distribution;Number of V0s;Events", 200, 0, 40000); |
126 | else |
127 | fHistV0Multiplicity = new TH1F("fHistV0Multiplicity", "Multiplicity distribution;Number of V0s;Events", 50, 0, 50); |
96cad837 |
128 | fListHist->Add(fHistV0Multiplicity); |
129 | } |
130 | if (!fHistV0OnFlyStatus) { |
131 | fHistV0OnFlyStatus = new TH1F("fHistV0OnFlyStatus", "V0 On fly status;status;Number of V0s", 3, 0, 3); |
132 | fListHist->Add(fHistV0OnFlyStatus); |
133 | } |
134 | |
135 | // V0 offline distributions |
136 | if (!fHistV0MultiplicityOff) { |
c8eada4f |
137 | if (fCollidingSystems) |
138 | fHistV0MultiplicityOff = new TH1F("fHistV0MultiplicityOff", "Multiplicity distribution;Number of V0s;Events", 200, 0, 40000); |
139 | else |
140 | fHistV0MultiplicityOff = new TH1F("fHistV0MultiplicityOff", "Multiplicity distribution;Number of V0s;Events", 50, 0, 50); |
96cad837 |
141 | fListHist->Add(fHistV0MultiplicityOff); |
142 | } |
143 | if (!fHistV0Chi2Off) { |
144 | fHistV0Chi2Off = new TH1F("fHistV0Chi2Off", "V0 chi2;chi2;Number of V0s", 33, 0, 33); |
145 | fListHist->Add(fHistV0Chi2Off); |
146 | } |
147 | if (!fHistDcaV0DaughtersOff) { |
148 | fHistDcaV0DaughtersOff = new TH1F("fHistDcaV0DaughtersOff", "DCA between V0 daughters;DCA (cm);Number of V0s", 300, 0, 3); |
149 | fListHist->Add(fHistDcaV0DaughtersOff); |
150 | } |
151 | if (!fHistV0CosineOfPointingAngleOff) { |
152 | fHistV0CosineOfPointingAngleOff = new TH1F("fHistV0CosineOfPointingAngleOff", "V0 Cosine of Pointing Angle;Number of V0s", 200, 0, 1); |
153 | fListHist->Add(fHistV0CosineOfPointingAngleOff); |
154 | } |
155 | if (!fHistV0RadiusOff) { |
156 | fHistV0RadiusOff = new TH1F("fHistV0RadiusOff", "V0 decay radius;Radius (cm);Number of V0s", 33, 0, 33); |
157 | fListHist->Add(fHistV0RadiusOff); |
158 | } |
159 | if (!fHistDcaV0ToPrimVertexOff) { |
160 | fHistDcaV0ToPrimVertexOff = new TH1F("fHistDcaV0ToPrimVertexOff", "DCA of V0 to Prim. Vertex;DCA (cm);Number of V0s", 300, 0, 3); |
161 | fListHist->Add(fHistDcaV0ToPrimVertexOff); |
162 | } |
163 | if (!fHistDcaPosToPrimVertexOff) { |
164 | fHistDcaPosToPrimVertexOff = new TH1F("fHistDcaPosToPrimVertexOff", "DCA of V0 neg daughter to Prim. Vertex;DCA (cm);Number of V0s", 300, 0, 3); |
165 | fListHist->Add(fHistDcaPosToPrimVertexOff); |
166 | } |
167 | if (!fHistDcaNegToPrimVertexOff) { |
168 | fHistDcaNegToPrimVertexOff = new TH1F("fHistDcaNegToPrimVertexOff", "DCA of V0 pos daughter to Prim. Vertex;DCA (cm);Number of V0s", 300, 0, 3); |
169 | fListHist->Add(fHistDcaNegToPrimVertexOff); |
170 | } |
171 | |
172 | if (!fHistMassK0Off) { |
c8eada4f |
173 | fHistMassK0Off = new TH1F("fHistMassK0Off","K^{0} candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});Counts",100,0.4,0.6); |
96cad837 |
174 | fListHist->Add(fHistMassK0Off); |
175 | } |
176 | if (!fHistMassLambdaOff) { |
c8eada4f |
177 | fHistMassLambdaOff = new TH1F("fHistMassLambdaOff","#Lambda^{0} candidates;M(p#pi^{-}) (GeV/c^{2});Counts",75,1.05,1.2); |
96cad837 |
178 | fListHist->Add(fHistMassLambdaOff); |
179 | } |
180 | if (!fHistMassAntiLambdaOff) { |
c8eada4f |
181 | fHistMassAntiLambdaOff = new TH1F("fHistMassAntiLambdaOff","#bar{#Lambda}^{0} candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts",75,1.05,1.2); |
96cad837 |
182 | fListHist->Add(fHistMassAntiLambdaOff); |
183 | } |
184 | |
185 | // V0 on-the-fly distributions |
186 | if (!fHistV0MultiplicityOn) { |
c8eada4f |
187 | if (fCollidingSystems) |
188 | fHistV0MultiplicityOn = new TH1F("fHistV0MultiplicityOn", "Multiplicity distribution;Number of V0s;Events", 200, 0, 40000); |
189 | else |
190 | fHistV0MultiplicityOn = new TH1F("fHistV0MultiplicityOn", "Multiplicity distribution;Number of V0s;Events", 50, 0, 50); |
96cad837 |
191 | fListHist->Add(fHistV0MultiplicityOn); |
192 | } |
193 | if (!fHistV0Chi2On) { |
194 | fHistV0Chi2On = new TH1F("fHistV0Chi2On", "V0 chi2;chi2;Number of V0s", 33, 0, 33); |
195 | fListHist->Add(fHistV0Chi2On); |
196 | } |
197 | if (!fHistDcaV0DaughtersOn) { |
198 | fHistDcaV0DaughtersOn = new TH1F("fHistDcaV0DaughtersOn", "DCA between V0 daughters;DCA (cm);Number of V0s", 300, 0, 3); |
199 | fListHist->Add(fHistDcaV0DaughtersOn); |
200 | } |
201 | if (!fHistV0CosineOfPointingAngleOn) { |
202 | fHistV0CosineOfPointingAngleOn = new TH1F("fHistV0CosineOfPointingAngleOn", "V0 Cosine of Pointing Angle;Number of V0s", 200, 0, 1); |
203 | fListHist->Add(fHistV0CosineOfPointingAngleOn); |
204 | } |
205 | if (!fHistV0RadiusOn) { |
206 | fHistV0RadiusOn = new TH1F("fHistV0RadiusOn", "V0 decay radius;Radius (cm);Number of V0s", 33, 0, 33); |
207 | fListHist->Add(fHistV0RadiusOn); |
208 | } |
209 | if (!fHistDcaV0ToPrimVertexOn) { |
210 | fHistDcaV0ToPrimVertexOn = new TH1F("fHistDcaV0ToPrimVertexOn", "DCA of V0 to Prim. Vertex;DCA (cm);Number of V0s", 300, 0, 3); |
211 | fListHist->Add(fHistDcaV0ToPrimVertexOn); |
212 | } |
213 | if (!fHistDcaPosToPrimVertexOn) { |
214 | fHistDcaPosToPrimVertexOn = new TH1F("fHistDcaPosToPrimVertexOn", "DCA of V0 neg daughter to Prim. Vertex;DCA (cm);Number of V0s", 300, 0, 3); |
215 | fListHist->Add(fHistDcaPosToPrimVertexOn); |
216 | } |
217 | if (!fHistDcaNegToPrimVertexOn) { |
218 | fHistDcaNegToPrimVertexOn = new TH1F("fHistDcaNegToPrimVertexOn", "DCA of V0 pos daughter to Prim. Vertex;DCA (cm);Number of V0s", 300, 0, 3); |
219 | fListHist->Add(fHistDcaNegToPrimVertexOn); |
220 | } |
221 | |
222 | if (!fHistMassK0On) { |
c8eada4f |
223 | fHistMassK0On = new TH1F("fHistMassK0On","K^{0} candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});Counts",100,0.4,0.6); |
96cad837 |
224 | fListHist->Add(fHistMassK0On); |
225 | } |
226 | if (!fHistMassLambdaOn) { |
c8eada4f |
227 | fHistMassLambdaOn = new TH1F("fHistMassLambdaOn","#Lambda^{0} candidates;M(p#pi^{-}) (GeV/c^{2});Counts",75,1.05,1.2); |
96cad837 |
228 | fListHist->Add(fHistMassLambdaOn); |
229 | } |
230 | if (!fHistMassAntiLambdaOn) { |
c8eada4f |
231 | fHistMassAntiLambdaOn = new TH1F("fHistMassAntiLambdaOn","#bar{#Lambda}^{0} candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts",75,1.05,1.2); |
96cad837 |
232 | fListHist->Add(fHistMassAntiLambdaOn); |
233 | } |
234 | |
235 | } |
236 | |
237 | //________________________________________________________________________ |
c8eada4f |
238 | void AliAnalysisTaskCheckV0::Exec(Option_t *) |
96cad837 |
239 | { |
240 | // Main loop |
241 | // Called for each event |
242 | |
c8eada4f |
243 | if(fAnalysisType == "ESD") { |
244 | |
245 | if (!fESD) { |
246 | Printf("ERROR: fESD not available"); |
247 | return; |
248 | } |
249 | // Printf("There are %d tracks in this event", fESD->GetNumberOfTracks()); |
250 | fHistTrackMultiplicity->Fill(fESD->GetNumberOfTracks()); |
251 | |
252 | Int_t nv0s = 0; |
253 | nv0s = fESD->GetNumberOfV0s(); |
254 | |
255 | Int_t lOnFlyStatus = 0, nv0sOn = 0, nv0sOff = 0; |
256 | Double_t lChi2V0 = 0; |
257 | Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0; |
258 | Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0; |
259 | Double_t lV0CosineOfPointingAngle = 0; |
260 | Double_t lV0Radius = 0; |
261 | Double_t lInvMassK0 = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0; |
262 | |
263 | for (Int_t iV0 = 0; iV0 < nv0s; iV0++) |
264 | {// This is the begining of the V0 loop |
265 | AliESDv0 *v0 = fESD->GetV0(iV0); |
266 | if (!v0) continue; |
267 | |
268 | Double_t tDecayVertexV0[3]; v0->GetXYZ(tDecayVertexV0[0],tDecayVertexV0[1],tDecayVertexV0[2]); |
269 | |
270 | lV0Radius = TMath::Sqrt(tDecayVertexV0[0]*tDecayVertexV0[0]+tDecayVertexV0[1]*tDecayVertexV0[1]); |
271 | |
272 | UInt_t lKeyPos = (UInt_t)TMath::Abs(v0->GetPindex()); |
273 | UInt_t lKeyNeg = (UInt_t)TMath::Abs(v0->GetNindex()); |
274 | |
275 | Double_t lMomPos[3]; v0->GetPPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]); |
276 | Double_t lMomNeg[3]; v0->GetNPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]); |
277 | |
278 | AliESDtrack *pTrack=fESD->GetTrack(lKeyPos); |
279 | AliESDtrack *nTrack=fESD->GetTrack(lKeyNeg); |
280 | if (!pTrack || !nTrack) { |
281 | Printf("ERROR: Could not retreive one of the daughter track"); |
282 | continue; |
283 | } |
284 | |
285 | Float_t tDcaPosToPrimVertex[2]; |
286 | if(pTrack) pTrack->GetImpactParameters(tDcaPosToPrimVertex[0],tDcaPosToPrimVertex[1]); |
287 | else { tDcaPosToPrimVertex[0]=999.; tDcaPosToPrimVertex[1]=999.;} |
288 | lDcaPosToPrimVertex = TMath::Sqrt(tDcaPosToPrimVertex[0]*tDcaPosToPrimVertex[0]+tDcaPosToPrimVertex[1]*tDcaPosToPrimVertex[1]); |
289 | |
290 | Float_t tDcaNegToPrimVertex[2]; |
291 | if(nTrack) nTrack->GetImpactParameters(tDcaNegToPrimVertex[0],tDcaNegToPrimVertex[1]); |
292 | else { tDcaNegToPrimVertex[0]=999.; tDcaNegToPrimVertex[1]=999.;} |
293 | lDcaNegToPrimVertex = TMath::Sqrt(tDcaNegToPrimVertex[0]*tDcaNegToPrimVertex[0]+tDcaNegToPrimVertex[1]*tDcaNegToPrimVertex[1]); |
294 | |
295 | |
296 | lOnFlyStatus = v0->GetOnFlyStatus(); |
297 | lChi2V0 = v0->GetChi2V0(); |
298 | lDcaV0Daughters = v0->GetDcaV0Daughters(); |
299 | lDcaV0ToPrimVertex = v0->GetD(); |
300 | lV0CosineOfPointingAngle = v0->GetV0CosineOfPointingAngle(); |
301 | |
302 | // Getting invariant mass infos directly from ESD |
303 | v0->ChangeMassHypothesis(310); |
304 | lInvMassK0 = v0->GetEffMass(); |
305 | v0->ChangeMassHypothesis(3122); |
306 | lInvMassLambda = v0->GetEffMass(); |
307 | v0->ChangeMassHypothesis(-3122); |
308 | lInvMassAntiLambda = v0->GetEffMass(); |
309 | |
310 | fHistV0OnFlyStatus->Fill(lOnFlyStatus); |
311 | if(!lOnFlyStatus){ |
312 | nv0sOff++; |
313 | fHistV0Chi2Off->Fill(lChi2V0); |
314 | fHistDcaV0ToPrimVertexOff->Fill(lDcaV0ToPrimVertex); |
315 | fHistDcaV0DaughtersOff->Fill(lDcaV0Daughters); |
316 | fHistV0CosineOfPointingAngleOff->Fill(lV0CosineOfPointingAngle); |
317 | |
318 | fHistV0RadiusOff->Fill(lV0Radius); |
319 | fHistDcaPosToPrimVertexOff->Fill(lDcaPosToPrimVertex); |
320 | fHistDcaNegToPrimVertexOff->Fill(lDcaNegToPrimVertex); |
321 | |
322 | // Filling invariant mass histos for all candidates |
323 | fHistMassK0Off->Fill(lInvMassK0); |
324 | fHistMassLambdaOff->Fill(lInvMassLambda); |
325 | fHistMassAntiLambdaOff->Fill(lInvMassAntiLambda); |
326 | } |
327 | else { |
328 | nv0sOn++; |
329 | fHistV0Chi2On->Fill(lChi2V0); |
330 | fHistDcaV0ToPrimVertexOn->Fill(lDcaV0ToPrimVertex); |
331 | fHistDcaV0DaughtersOn->Fill(lDcaV0Daughters); |
332 | fHistV0CosineOfPointingAngleOn->Fill(lV0CosineOfPointingAngle); |
333 | |
334 | fHistV0RadiusOn->Fill(lV0Radius); |
335 | fHistDcaPosToPrimVertexOn->Fill(lDcaPosToPrimVertex); |
336 | fHistDcaNegToPrimVertexOn->Fill(lDcaNegToPrimVertex); |
337 | |
338 | // Filling invariant mass histos for all candidates |
339 | fHistMassK0On->Fill(lInvMassK0); |
340 | fHistMassLambdaOn->Fill(lInvMassLambda); |
341 | fHistMassAntiLambdaOn->Fill(lInvMassAntiLambda); |
342 | } |
343 | }// This is the end of the V0 loop |
344 | |
345 | fHistV0Multiplicity->Fill(nv0s); |
346 | fHistV0MultiplicityOff->Fill(nv0sOff); |
347 | fHistV0MultiplicityOn->Fill(nv0sOn); |
348 | } // end of "ESD" analysis |
349 | |
350 | else if(fAnalysisType == "AOD") { |
351 | |
352 | if (!fAOD) { |
353 | Printf("ERROR: fAOD not available"); |
354 | return; |
355 | } |
356 | Printf("Well: not fully implemented yet"); |
357 | |
358 | // Printf("There are %d tracks in this event", fAOD->GetNumberOfTracks()); |
359 | fHistTrackMultiplicity->Fill(fAOD->GetNumberOfTracks()); |
360 | |
361 | Int_t nv0s = 0; |
362 | nv0s = fAOD->GetNV0s(); |
363 | |
364 | Int_t lOnFlyStatus = 0, nv0sOn = 0, nv0sOff = 0; |
365 | Double_t lChi2V0 = 0; |
366 | Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0; |
367 | Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0; |
368 | Double_t lV0CosineOfPointingAngle = 0; |
369 | Double_t lV0Radius = 0; |
370 | Double_t lInvMassK0 = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0; |
371 | |
372 | for (Int_t iV0 = 0; iV0 < nv0s; iV0++) |
373 | {// This is the begining of the V0 loop |
374 | AliAODv0 *v0 = fAOD->GetV0(iV0); |
375 | if (!v0) continue; |
376 | |
377 | lV0Radius = v0->RadiusV0(); |
378 | lDcaPosToPrimVertex = v0->DcaPosToPrimVertex(); |
379 | lDcaNegToPrimVertex = v0->DcaNegToPrimVertex(); |
380 | |
381 | |
382 | lOnFlyStatus = 0; // v0->GetOnFlyStatus(); TBD |
383 | lChi2V0 = v0->Chi2V0(); |
384 | lDcaV0Daughters = v0->DcaV0Daughters(); |
385 | lDcaV0ToPrimVertex = v0->DcaV0ToPrimVertex(); |
386 | lV0CosineOfPointingAngle = 0; // v0->GetV0CosineOfPointingAngle(); TBD |
387 | |
388 | lInvMassK0 = v0->MassK0Short(); |
389 | lInvMassLambda = v0->MassLambda(); |
390 | lInvMassAntiLambda = v0->MassAntiLambda(); |
391 | |
392 | fHistV0OnFlyStatus->Fill(lOnFlyStatus); |
393 | if(!lOnFlyStatus){ |
394 | nv0sOff++; |
395 | fHistV0Chi2Off->Fill(lChi2V0); |
396 | fHistDcaV0ToPrimVertexOff->Fill(lDcaV0ToPrimVertex); |
397 | fHistDcaV0DaughtersOff->Fill(lDcaV0Daughters); |
398 | fHistV0CosineOfPointingAngleOff->Fill(lV0CosineOfPointingAngle); |
399 | |
400 | fHistV0RadiusOff->Fill(lV0Radius); |
401 | fHistDcaPosToPrimVertexOff->Fill(lDcaPosToPrimVertex); |
402 | fHistDcaNegToPrimVertexOff->Fill(lDcaNegToPrimVertex); |
403 | |
404 | // Filling invariant mass histos for all candidates |
405 | fHistMassK0Off->Fill(lInvMassK0); |
406 | fHistMassLambdaOff->Fill(lInvMassLambda); |
407 | fHistMassAntiLambdaOff->Fill(lInvMassAntiLambda); |
408 | } |
409 | else { |
410 | nv0sOn++; |
411 | fHistV0Chi2On->Fill(lChi2V0); |
412 | fHistDcaV0ToPrimVertexOn->Fill(lDcaV0ToPrimVertex); |
413 | fHistDcaV0DaughtersOn->Fill(lDcaV0Daughters); |
414 | fHistV0CosineOfPointingAngleOn->Fill(lV0CosineOfPointingAngle); |
415 | |
416 | fHistV0RadiusOn->Fill(lV0Radius); |
417 | fHistDcaPosToPrimVertexOn->Fill(lDcaPosToPrimVertex); |
418 | fHistDcaNegToPrimVertexOn->Fill(lDcaNegToPrimVertex); |
419 | |
420 | // Filling invariant mass histos for all candidates |
421 | fHistMassK0On->Fill(lInvMassK0); |
422 | fHistMassLambdaOn->Fill(lInvMassLambda); |
423 | fHistMassAntiLambdaOn->Fill(lInvMassAntiLambda); |
424 | } |
425 | }// This is the end of the V0 loop |
426 | |
427 | fHistV0Multiplicity->Fill(nv0s); |
428 | fHistV0MultiplicityOff->Fill(nv0sOff); |
429 | fHistV0MultiplicityOn->Fill(nv0sOn); |
430 | } // end of "AOD" analysis |
96cad837 |
431 | |
432 | // Post output data. |
433 | PostData(0, fHistTrackMultiplicity); |
434 | PostData(1, fListHist); |
435 | } |
436 | |
437 | //________________________________________________________________________ |
c8eada4f |
438 | void AliAnalysisTaskCheckV0::Terminate(Option_t *) |
96cad837 |
439 | { |
440 | // Draw result to the screen |
441 | // Called once at the end of the query |
442 | |
443 | fHistTrackMultiplicity = dynamic_cast<TH1F*> (GetOutputData(0)); |
444 | if (!fHistTrackMultiplicity) { |
445 | Printf("ERROR: fHistTrackMultiplicity not available"); |
446 | return; |
447 | } |
448 | fHistV0Multiplicity = dynamic_cast<TH1F*> (((TList*)GetOutputData(1))->FindObject("fHistV0Multiplicity")); |
449 | if (!fHistV0Multiplicity) { |
450 | Printf("ERROR: fHistV0Multiplicity not available"); |
451 | return; |
452 | } |
453 | fHistV0MultiplicityOff = dynamic_cast<TH1F*> (((TList*)GetOutputData(1))->FindObject("fHistV0MultiplicityOff")); |
454 | if (!fHistV0MultiplicityOff) { |
455 | Printf("ERROR: fHistV0MultiplicityOff not available"); |
456 | return; |
457 | } |
458 | fHistV0MultiplicityOn = dynamic_cast<TH1F*> (((TList*)GetOutputData(1))->FindObject("fHistV0MultiplicityOn")); |
459 | if (!fHistV0MultiplicityOn) { |
460 | Printf("ERROR: fHistV0MultiplicityOn not available"); |
461 | return; |
462 | } |
463 | |
c8eada4f |
464 | TCanvas *canCheckV0 = new TCanvas("AliAnalysisTaskCheckV0","Multiplicity",10,10,510,510); |
465 | canCheckV0->Divide(2,2); |
466 | canCheckV0->cd(1)->SetLogy(); |
96cad837 |
467 | fHistTrackMultiplicity->SetMarkerStyle(22); |
468 | fHistTrackMultiplicity->DrawCopy("E"); |
469 | fHistV0Multiplicity->SetMarkerStyle(26); |
470 | fHistV0Multiplicity->DrawCopy("ESAME"); |
471 | fHistV0MultiplicityOff->SetMarkerStyle(24); |
472 | fHistV0MultiplicityOff->DrawCopy("ESAME"); |
473 | fHistV0MultiplicityOn->SetMarkerStyle(20); |
474 | fHistV0MultiplicityOn->DrawCopy("ESAME"); |
c8eada4f |
475 | |
476 | TLegend *legendMultiplicity = new TLegend(0.4,0.4,0.65,0.6); |
477 | legendMultiplicity->AddEntry(fHistTrackMultiplicity,"tracks"); |
478 | legendMultiplicity->AddEntry(fHistV0Multiplicity,"all V^{0}"); |
479 | legendMultiplicity->AddEntry(fHistV0MultiplicityOff,"offline V^{0}"); |
480 | legendMultiplicity->AddEntry(fHistV0MultiplicityOn,"onthefly V^{0}"); |
481 | legendMultiplicity->Draw(); |
482 | |
483 | fHistMassK0Off = dynamic_cast<TH1F*> (((TList*)GetOutputData(1))->FindObject("fHistMassK0Off")); |
484 | if (!fHistMassK0Off) { |
485 | Printf("ERROR: fHistMassK0Off not available"); |
486 | return; |
487 | } |
488 | fHistMassK0On = dynamic_cast<TH1F*> (((TList*)GetOutputData(1))->FindObject("fHistMassK0On")); |
489 | if (!fHistMassK0On) { |
490 | Printf("ERROR: fHistMassK0On not available"); |
491 | return; |
492 | } |
493 | fHistMassLambdaOff = dynamic_cast<TH1F*> (((TList*)GetOutputData(1))->FindObject("fHistMassLambdaOff")); |
494 | if (!fHistMassLambdaOff) { |
495 | Printf("ERROR: fHistMassLambdaOff not available"); |
496 | return; |
497 | } |
498 | fHistMassLambdaOn = dynamic_cast<TH1F*> (((TList*)GetOutputData(1))->FindObject("fHistMassLambdaOn")); |
499 | if (!fHistMassLambdaOn) { |
500 | Printf("ERROR: fHistMassLambdaOn not available"); |
501 | return; |
502 | } |
503 | fHistMassAntiLambdaOff = dynamic_cast<TH1F*> (((TList*)GetOutputData(1))->FindObject("fHistMassAntiLambdaOff")); |
504 | if (!fHistMassAntiLambdaOff) { |
505 | Printf("ERROR: fHistMassAntiLambdaOff not available"); |
506 | return; |
507 | } |
508 | fHistMassAntiLambdaOn = dynamic_cast<TH1F*> (((TList*)GetOutputData(1))->FindObject("fHistMassAntiLambdaOn")); |
509 | if (!fHistMassAntiLambdaOn) { |
510 | Printf("ERROR: fHistMassAntiLambdaOn not available"); |
511 | return; |
512 | } |
513 | |
514 | canCheckV0->cd(2); |
515 | fHistMassK0Off->SetMarkerStyle(24); |
516 | fHistMassK0Off->DrawCopy("E"); |
517 | fHistMassK0On->SetMarkerStyle(20); |
518 | fHistMassK0On->DrawCopy("ESAME"); |
519 | |
520 | canCheckV0->cd(3); |
521 | fHistMassLambdaOff->SetMarkerStyle(24); |
522 | fHistMassLambdaOff->DrawCopy("E"); |
523 | fHistMassLambdaOn->SetMarkerStyle(20); |
524 | fHistMassLambdaOn->DrawCopy("ESAME"); |
525 | |
526 | canCheckV0->cd(4); |
527 | fHistMassAntiLambdaOff->SetMarkerStyle(24); |
528 | fHistMassAntiLambdaOff->DrawCopy("E"); |
529 | fHistMassAntiLambdaOn->SetMarkerStyle(20); |
530 | fHistMassAntiLambdaOn->DrawCopy("ESAME"); |
96cad837 |
531 | } |