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