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