2 /**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
17 //-----------------------------------------------------------------
18 // AliAnalysisTaskStrange class
19 // This task is for single strange study from ESD/AOD
20 // Origin: H.Ricaud, Helene.Ricaud@IReS.in2p3.fr
21 //-----------------------------------------------------------------
28 #include "AliAnalysisTaskSE.h"
30 #include "AliESDVertex.h"
31 #include "AliESDEvent.h"
32 #include "AliAODEvent.h"
35 #include "AliESDtrack.h"
37 #include "AliAODTrack.h"
41 #include "AliAnalysisTaskStrange.h"
43 ClassImp(AliAnalysisTaskStrange)
45 //________________________________________________________________________
46 AliAnalysisTaskStrange::AliAnalysisTaskStrange()
47 : AliAnalysisTaskSE(), fAnalysisType("ESD"), fCollidingSystems(0), fOption("no"), fListHist(),
48 fHistPrimaryVertexPosX(0), fHistPrimaryVertexPosY(0), fHistPrimaryVertexPosZ(0),
49 fHistTrackMultiplicity(0), fHistV0Multiplicity(0),
50 fHistDcaPosToPrimVertex(0), fHistDcaNegToPrimVertex(0),
51 fHistDcaPosToPrimVertexZoom(0), fHistDcaNegToPrimVertexZoom(0),
52 fHistRadiusV0(0), fHistDecayLengthV0(0), fHistDcaV0Daughters(0), fHistChi2(0),
53 fHistCosPointAngle(0), fHistCosPointAngleZoom(0),
54 fHistV0MultiplicityOff(0),
55 fHistPtVsYK0sOff(0), fHistPtVsYLambdaOff(0), fHistPtVsYAntiLambdaOff(0),
56 fHistMassK0sOff(0), fHistMassLambdaOff(0), fHistMassAntiLambdaOff(0),
57 fHistMassVsRadiusK0sOff(0), fHistMassVsRadiusLambdaOff(0), fHistMassVsRadiusAntiLambdaOff(0),
58 fHistPtVsMassK0sOff(0), fHistPtVsMassLambdaOff(0), fHistPtVsMassAntiLambdaOff(0),
59 fHistArmenterosPodolanskiOff(0),
60 fHistV0MultiplicityOn(0),
61 fHistPtVsYK0sOn(0), fHistPtVsYLambdaOn(0), fHistPtVsYAntiLambdaOn(0),
62 fHistMassK0sOn(0), fHistMassLambdaOn(0), fHistMassAntiLambdaOn(0),
63 fHistMassVsRadiusK0sOn(0), fHistMassVsRadiusLambdaOn(0), fHistMassVsRadiusAntiLambdaOn(0),
64 fHistPtVsMassK0sOn(0), fHistPtVsMassLambdaOn(0), fHistPtVsMassAntiLambdaOn(0),
65 fHistArmenterosPodolanskiOn(0)
69 //________________________________________________________________________
70 AliAnalysisTaskStrange::AliAnalysisTaskStrange(const char *name, const char *optCuts)
71 : AliAnalysisTaskSE(name), fAnalysisType("ESD"), fCollidingSystems(0), fOption(optCuts), fListHist(),
72 fHistPrimaryVertexPosX(0), fHistPrimaryVertexPosY(0), fHistPrimaryVertexPosZ(0),
73 fHistTrackMultiplicity(0), fHistV0Multiplicity(0),
74 fHistDcaPosToPrimVertex(0), fHistDcaNegToPrimVertex(0),
75 fHistDcaPosToPrimVertexZoom(0), fHistDcaNegToPrimVertexZoom(0),
76 fHistRadiusV0(0), fHistDecayLengthV0(0), fHistDcaV0Daughters(0), fHistChi2(0),
77 fHistCosPointAngle(0), fHistCosPointAngleZoom(0),
78 fHistV0MultiplicityOff(0),
79 fHistPtVsYK0sOff(0), fHistPtVsYLambdaOff(0), fHistPtVsYAntiLambdaOff(0),
80 fHistMassK0sOff(0), fHistMassLambdaOff(0), fHistMassAntiLambdaOff(0),
81 fHistMassVsRadiusK0sOff(0), fHistMassVsRadiusLambdaOff(0), fHistMassVsRadiusAntiLambdaOff(0),
82 fHistPtVsMassK0sOff(0), fHistPtVsMassLambdaOff(0), fHistPtVsMassAntiLambdaOff(0),
83 fHistArmenterosPodolanskiOff(0),
84 fHistV0MultiplicityOn(0),
85 fHistPtVsYK0sOn(0), fHistPtVsYLambdaOn(0), fHistPtVsYAntiLambdaOn(0),
86 fHistMassK0sOn(0), fHistMassLambdaOn(0), fHistMassAntiLambdaOn(0),
87 fHistMassVsRadiusK0sOn(0), fHistMassVsRadiusLambdaOn(0), fHistMassVsRadiusAntiLambdaOn(0),
88 fHistPtVsMassK0sOn(0), fHistPtVsMassLambdaOn(0), fHistPtVsMassAntiLambdaOn(0),
89 fHistArmenterosPodolanskiOn(0)
92 // Define output slots only here
93 // Output slot #1 writes into a TList container
94 DefineOutput(1, TList::Class());
97 //________________________________________________________________________
98 void AliAnalysisTaskStrange::UserCreateOutputObjects()
103 fListHist = new TList();
106 fHistPrimaryVertexPosX = new TH1F("h1PrimaryVertexPosX", "Primary Vertex Position X;Primary Vertex Position X (cm);Events",100,-0.5,0.5);
107 fListHist->Add(fHistPrimaryVertexPosX);
108 fHistPrimaryVertexPosY = new TH1F("h1PrimaryVertexPosY", "Primary Vertex Position Y;Primary Vertex Position Y (cm);Events",100,-0.5,0.5);
109 fListHist->Add(fHistPrimaryVertexPosY);
110 fHistPrimaryVertexPosZ = new TH1F("h1PrimaryVertexPosZ", "Primary Vertex Position Z;Primary Vertex Position Z (cm);Events",200,-2.0,2.0);
111 fListHist->Add(fHistPrimaryVertexPosZ);
114 if (!fHistTrackMultiplicity) {
115 if (fCollidingSystems)
116 fHistTrackMultiplicity = new TH1F("fHistTrackMultiplicity", "Multiplicity distribution;Number of tracks;Events", 200, 0, 40000);
118 fHistTrackMultiplicity = new TH1F("fHistTrackMultiplicity", "Multiplicity distribution;Number of tracks;Events", 250, 0, 250);
119 fListHist->Add(fHistTrackMultiplicity);
121 if (!fHistV0Multiplicity) {
122 if (fCollidingSystems)
123 fHistV0Multiplicity = new TH1F("fHistV0Multiplicity", "Multiplicity distribution;Number of V0s;Events", 200, 0, 40000);
125 fHistV0Multiplicity = new TH1F("fHistV0Multiplicity", "Multiplicity distribution;Number of V0s;Events", 50, 0, 50);
126 fListHist->Add(fHistV0Multiplicity);
130 fHistDcaPosToPrimVertex = new TH2F("h2DcaPosToPrimVertex", "Positive V0 daughter;dca(cm);Status",500,0,5,2,-0.5,1.5);
131 fListHist->Add(fHistDcaPosToPrimVertex);
132 fHistDcaNegToPrimVertex = new TH2F("h2DcaNegToPrimVertex", "Negative V0 daughter;dca(cm);Status",500,0,5,2,-0.5,1.5);
133 fListHist->Add(fHistDcaNegToPrimVertex);
134 fHistDcaPosToPrimVertexZoom = new TH2F("h2DcaPosToPrimVertexZoom", "Positive V0 daughter;dca(cm);Status",100,0,0.1,2,-0.5,1.5);
135 fListHist->Add(fHistDcaPosToPrimVertexZoom);
136 fHistDcaNegToPrimVertexZoom = new TH2F("h2DcaNegToPrimVertexZoom", "Negative V0 daughter;dca(cm);Status",100,0,0.1,2,-0.5,1.5);
137 fListHist->Add(fHistDcaNegToPrimVertexZoom);
138 fHistRadiusV0 = new TH2F("h2RadiusV0", "Radius;Radius(cm);Status",1000,0,100,2,-0.5,1.5);
139 fListHist->Add(fHistRadiusV0);
140 fHistDecayLengthV0 = new TH2F("h2DecayLengthV0", "V0s decay Length;decay length(cm);Status", 200, 0, 100,2,-0.5,1.5);
141 fListHist->Add(fHistDecayLengthV0);
142 fHistDcaV0Daughters = new TH2F("h2DcaV0Daughters", "DCA between daughters;dca(cm);Status", 160, 0, 4,2,-0.5,1.5);
143 fListHist->Add(fHistDcaV0Daughters);
144 fHistChi2 = new TH2F("h2Chi2", "V0s chi2;chi2;Status", 33, 0, 33,2,-0.5,1.5);
145 fListHist->Add(fHistChi2);
146 fHistCosPointAngle = new TH2F("h2CosPointAngle", "Cosine of V0's pointing angle", 100,0,1,2,-0.5,1.5);
147 fListHist->Add(fHistCosPointAngle);
148 fHistCosPointAngleZoom = new TH2F("h2CosPointAngleZoom", "Cosine of V0's pointing angle", 100,0.9,1,2,-0.5,1.5);
149 fListHist->Add(fHistCosPointAngleZoom);
151 // bounds of histograms:
153 const Double_t radius[10] = {0.0,2.5,2.9,3.9,7.6,15.0,23.9,37.8,42.8,100.0};
154 Int_t nBinRadius = 9;
156 // V0 offline distributions
157 if (!fHistV0MultiplicityOff) {
158 if (fCollidingSystems)
159 fHistV0MultiplicityOff = new TH1F("fHistV0MultiplicityOff", "Multiplicity distribution;Number of V0s;Events", 200, 0, 40000);
161 fHistV0MultiplicityOff = new TH1F("fHistV0MultiplicityOff", "Multiplicity distribution;Number of V0s;Events", 50, 0, 50);
162 fListHist->Add(fHistV0MultiplicityOff);
165 fHistPtVsYK0sOff = new TH2F("h2PtVsYK0sOff", "K^{0} Offline candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
166 fListHist->Add(fHistPtVsYK0sOff);
167 fHistPtVsYLambdaOff = new TH2F("h2PtVsYLambdaOff", "#Lambda^{0} Offline candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
168 fListHist->Add(fHistPtVsYLambdaOff);
169 fHistPtVsYAntiLambdaOff = new TH2F("h2PtVsYAntiLambdaOff", "#bar{#Lambda}^{0} Offline candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
170 fListHist->Add(fHistPtVsYAntiLambdaOff);
172 fHistMassK0sOff = new TH1F("h1MassK0sOff", "K^{0} Offline candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});Counts", 100, 0.4, 0.6);
173 fListHist->Add(fHistMassK0sOff);
174 fHistMassLambdaOff = new TH1F("h1MassLambdaOff", "#Lambda^{0} Offline candidates;M(p#pi^{-}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
175 fListHist->Add(fHistMassLambdaOff);
176 fHistMassAntiLambdaOff = new TH1F("h1MassAntiLambdaOff", "#bar{#Lambda}^{0} Offline candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
177 fListHist->Add(fHistMassAntiLambdaOff);
179 fHistMassVsRadiusK0sOff = new TH2F("h2MassVsRadiusK0sOff", "K^{0} Offline candidates;radius (cm);M(#pi^{+}#pi^{-}) (GeV/c^{2})",nBinRadius,radius, 200, 0.4, 0.6);
180 fListHist->Add(fHistMassVsRadiusK0sOff);
181 fHistMassVsRadiusLambdaOff = new TH2F("h2MassVsRadiusLambdaOff", "#Lambda Offline candidates;radius (cm);M(p#pi^{-}) (GeV/c^{2})",nBinRadius,radius, 140, 1.06, 1.2);
182 fListHist->Add(fHistMassVsRadiusLambdaOff);
183 fHistMassVsRadiusAntiLambdaOff = new TH2F("h2MassVsRadiusAntiLambdaOff", "#bar{#Lambda} Offline candidates;radius (cm);M(#bar{p}#pi^{+}) (GeV/c^{2})",nBinRadius,radius, 140, 1.06, 1.2);
184 fListHist->Add(fHistMassVsRadiusAntiLambdaOff);
186 fHistPtVsMassK0sOff = new TH2F("h2PtVsMassK0sOff","K^{0} Offline candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});p_{t} (GeV/c)",200, 0.4, 0.6,100,0,10);
187 fListHist->Add(fHistPtVsMassK0sOff);
188 fHistPtVsMassLambdaOff = new TH2F("h2PtVsMassLambdaOff","#Lambda^{0} Offline candidates;M(p#pi^{-}) (GeV/c^{2});p_{t} (GeV/c)",140, 1.06, 1.2,100,0,10);
189 fListHist->Add(fHistPtVsMassLambdaOff);
190 fHistPtVsMassAntiLambdaOff = new TH2F("h2PtVsMassAntiLambdaOff","#bar{#Lambda}^{0} Offline candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});p_{t} (GeV/c)",140, 1.06, 1.2,100,0,10);
191 fListHist->Add(fHistPtVsMassAntiLambdaOff);
192 //ArmenterosPodolanski:
193 fHistArmenterosPodolanskiOff = new TH2F("h2ArmenterosPodolanskiOff","Armenteros-Podolanski Offline phase space;#alpha;p_{t} arm",100,-1.0,1.0,50,0,0.5);
194 fListHist->Add(fHistArmenterosPodolanskiOff);
196 // V0 on-the-fly distributions
197 if (!fHistV0MultiplicityOn) {
198 if (fCollidingSystems)
199 fHistV0MultiplicityOn = new TH1F("fHistV0MultiplicityOn", "Multiplicity distribution;Number of V0s;Events", 200, 0, 40000);
201 fHistV0MultiplicityOn = new TH1F("fHistV0MultiplicityOn", "Multiplicity distribution;Number of V0s;Events", 50, 0, 50);
202 fListHist->Add(fHistV0MultiplicityOn);
205 fHistPtVsYK0sOn = new TH2F("h2PtVsYK0sOn", "K^{0} Onthefly candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
206 fListHist->Add(fHistPtVsYK0sOn);
207 fHistPtVsYLambdaOn = new TH2F("h2PtVsYLambdaOn", "#Lambda^{0} Onthefly candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
208 fListHist->Add(fHistPtVsYLambdaOn);
209 fHistPtVsYAntiLambdaOn = new TH2F("h2PtVsYAntiLambdaOn", "#bar{#Lambda}^{0} Onthefly candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
210 fListHist->Add(fHistPtVsYAntiLambdaOn);
212 fHistMassK0sOn = new TH1F("h1MassK0sOn", "K^{0} Onthefly candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});Counts", 100, 0.4, 0.6);
213 fListHist->Add(fHistMassK0sOn);
214 fHistMassLambdaOn = new TH1F("h1MassLambdaOn", "#Lambda^{0} Onthefly candidates;M(p#pi^{-}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
215 fListHist->Add(fHistMassLambdaOn);
216 fHistMassAntiLambdaOn = new TH1F("h1MassAntiLambdaOn", "#bar{#Lambda}^{0} Onthefly candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
217 fListHist->Add(fHistMassAntiLambdaOn);
219 fHistMassVsRadiusK0sOn = new TH2F("h2MassVsRadiusK0sOn", "K^{0} Onthefly candidates;radius (cm);M(#pi^{+}#pi^{-}) (GeV/c^{2})",nBinRadius,radius, 200, 0.4, 0.6);
220 fListHist->Add(fHistMassVsRadiusK0sOn);
221 fHistMassVsRadiusLambdaOn = new TH2F("h2MassVsRadiusLambdaOn", "#Lambda Onthefly candidates;radius (cm);M(p#pi^{-}) (GeV/c^{2})",nBinRadius,radius, 140, 1.06, 1.2);
222 fListHist->Add(fHistMassVsRadiusLambdaOn);
223 fHistMassVsRadiusAntiLambdaOn = new TH2F("h2MassVsRadiusAntiLambdaOn", "#bar{#Lambda} Onthefly candidates;radius (cm);M(#bar{p}#pi^{+}) (GeV/c^{2})",nBinRadius,radius, 140, 1.06, 1.2);
224 fListHist->Add(fHistMassVsRadiusAntiLambdaOn);
226 fHistPtVsMassK0sOn = new TH2F("h2PtVsMassK0sOn","K^{0} Onthefly candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});p_{t} (GeV/c)",200, 0.4, 0.6,100,0,10);
227 fListHist->Add(fHistPtVsMassK0sOn);
228 fHistPtVsMassLambdaOn = new TH2F("h2PtVsMassLambdaOn","#Lambda^{0} Onthefly candidates;M(p#pi^{-}) (GeV/c^{2});p_{t} (GeV/c)",140, 1.06, 1.2,100,0,10);
229 fListHist->Add(fHistPtVsMassLambdaOn);
230 fHistPtVsMassAntiLambdaOn = new TH2F("h2PtVsMassAntiLambdaOn","#bar{#Lambda}^{0} Onthefly candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});p_{t} (GeV/c)",140, 1.06, 1.2,100,0,10);
231 fListHist->Add(fHistPtVsMassAntiLambdaOn);
232 //ArmenterosPodolanski:
233 fHistArmenterosPodolanskiOn = new TH2F("h2ArmenterosPodolanskiOn","Armenteros-Podolanski Onthefly phase space;#alpha;p_{t} arm",100,-1.0,1.0,50,0,0.5);
234 fListHist->Add(fHistArmenterosPodolanskiOn);
237 //________________________________________________________________________
238 void AliAnalysisTaskStrange::UserExec(Option_t *)
241 // Called for each event
242 AliVEvent* lEvent = InputEvent();
244 Printf("ERROR: Event not available");
248 if (!(lEvent->GetNumberOfTracks())) {
249 Printf("Strange analysis task: There is no track in this event");
252 fHistTrackMultiplicity->Fill(lEvent->GetNumberOfTracks());
254 Double_t tPrimaryVtxPosition[3];
255 Double_t tPrimaryVtxCov[3];
257 nv0s = lEvent->GetNumberOfV0s();
258 // Printf("Strange analysis task: There are %d v0s in this event",nv0s);
260 Int_t lOnFlyStatus = 0, nv0sOn = 0, nv0sOff = 0;
261 Double_t lChi2V0 = 0;
262 Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0;
263 Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0;
264 Double_t lV0CosineOfPointingAngle = 0;
265 Double_t lV0Radius = 0;
266 Double_t lInvMassK0s = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;
267 Double_t lPt = 0, lRapK0s = 0, lRapLambda = 0;
268 Double_t lAlphaV0 = 0, lPtArmV0 = 0;
270 if(fAnalysisType == "ESD") {
272 const AliESDVertex *primaryVtx = ((AliESDEvent*)lEvent)->GetPrimaryVertex();
273 tPrimaryVtxPosition[0] = primaryVtx->GetXv();
274 tPrimaryVtxPosition[1] = primaryVtx->GetYv();
275 tPrimaryVtxPosition[2] = primaryVtx->GetZv();
277 fHistPrimaryVertexPosX->Fill(tPrimaryVtxPosition[0]);
278 fHistPrimaryVertexPosY->Fill(tPrimaryVtxPosition[1]);
279 fHistPrimaryVertexPosZ->Fill(tPrimaryVtxPosition[2]);
281 Double_t lMagneticField = ((AliESDEvent*)lEvent)->GetMagneticField();
283 primaryVtx->GetCovMatrix(tPrimaryVtxCov);
284 AliAODVertex *primary = new AliAODVertex(tPrimaryVtxPosition, tPrimaryVtxCov, primaryVtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
287 // to get info from ESD files and fill AliAODVertex:
288 Double_t tdcaDaughterToPrimVertex[2]; // ..[0] = Pos and ..[1] = Neg
291 Double_t tV0Position[3];
294 // to fill AliAODtrack:
296 Double_t tTrackPosition[3];
297 Double_t tTrackCovTr[21];
298 Double_t tTrackPid[10];
300 AliAODTrack *myPosAodTrack = new AliAODTrack();
301 AliAODTrack *myNegAodTrack = new AliAODTrack();
302 AliAODVertex *myAODVertex = new AliAODVertex();
303 AliAODv0 *myAODv0 = new AliAODv0();
305 for (Int_t iV0 = 0; iV0 < nv0s; iV0++)
306 {// This is the begining of the V0 loop
307 AliESDv0 *v0 = ((AliESDEvent*)lEvent)->GetV0(iV0);
311 v0->GetXYZ(tV0Position[0], tV0Position[1], tV0Position[2]);
312 v0->GetPosCov(tV0Cov);
313 myAODVertex->SetPosition(tV0Position[0],tV0Position[1],tV0Position[2]);
314 myAODVertex->SetCovMatrix(tV0Cov);
315 myAODVertex->SetChi2perNDF(v0->GetChi2V0());
316 myAODVertex->SetID((Short_t)iV0);
317 myAODVertex->SetParent(primary);
318 myAODVertex->SetType(AliAODVertex::kV0);
320 // AliAODtrack (V0 Daughters)
321 UInt_t lKeyPos = (UInt_t)TMath::Abs(v0->GetPindex());
322 UInt_t lKeyNeg = (UInt_t)TMath::Abs(v0->GetNindex());
324 AliESDtrack *pTrack = ((AliESDEvent*)lEvent)->GetTrack(lKeyPos);
325 AliESDtrack *nTrack = ((AliESDEvent*)lEvent)->GetTrack(lKeyNeg);
326 if (!pTrack || !nTrack) {
327 Printf("ERROR: Could not retreive one of the daughter track");
330 UInt_t lLabelPos = (UInt_t)TMath::Abs(pTrack->GetLabel());
331 UInt_t lLabelNeg = (UInt_t)TMath::Abs(nTrack->GetLabel());
333 myPosAodTrack->SetID((Short_t)(pTrack->GetID()));
334 myPosAodTrack->SetLabel(lLabelPos);
335 pTrack->GetPxPyPz(tTrackP);
336 myPosAodTrack->SetP(tTrackP);
337 pTrack->GetXYZ(tTrackPosition);
338 myPosAodTrack->SetPosition(tTrackPosition,kFALSE);
339 pTrack->GetCovarianceXYZPxPyPz(tTrackCovTr);
340 myPosAodTrack->SetCovMatrix(tTrackCovTr);
341 pTrack->GetESDpid(tTrackPid);
342 myPosAodTrack->SetPID(tTrackPid);
343 myPosAodTrack->SetCharge((Short_t)(pTrack->Charge()));
344 myPosAodTrack->SetITSClusterMap(pTrack->GetITSClusterMap());
345 myPosAodTrack->SetProdVertex(myAODVertex);
346 myPosAodTrack->SetUsedForVtxFit(kTRUE);
347 myPosAodTrack->SetUsedForPrimVtxFit(kFALSE);
348 myPosAodTrack->SetType(AliAODTrack::kSecondary);
349 myPosAodTrack->ConvertAliPIDtoAODPID();
351 myNegAodTrack->SetID((Short_t)(nTrack->GetID()));
352 myNegAodTrack->SetLabel(lLabelNeg);
353 nTrack->GetPxPyPz(tTrackP);
354 myNegAodTrack->SetP(tTrackP);
355 nTrack->GetXYZ(tTrackPosition);
356 myNegAodTrack->SetPosition(tTrackPosition,kFALSE);
357 nTrack->GetCovarianceXYZPxPyPz(tTrackCovTr);
358 myNegAodTrack->SetCovMatrix(tTrackCovTr);
359 nTrack->GetESDpid(tTrackPid);
360 myNegAodTrack->SetPID(tTrackPid);
361 myNegAodTrack->SetCharge((Short_t)(nTrack->Charge()));
362 myNegAodTrack->SetITSClusterMap(nTrack->GetITSClusterMap());
363 myNegAodTrack->SetProdVertex(myAODVertex);
364 myNegAodTrack->SetUsedForVtxFit(kTRUE);
365 myNegAodTrack->SetUsedForPrimVtxFit(kFALSE);
366 myNegAodTrack->SetType(AliAODTrack::kSecondary);
367 myNegAodTrack->ConvertAliPIDtoAODPID();
369 myAODVertex->AddDaughter(myPosAodTrack);
370 myAODVertex->AddDaughter(myNegAodTrack);
373 lDcaV0Daughters = v0->GetDcaV0Daughters();
374 lDcaV0ToPrimVertex = v0->GetD(tPrimaryVtxPosition[0],tPrimaryVtxPosition[1],tPrimaryVtxPosition[2]);
376 if (pTrack) tdcaDaughterToPrimVertex[0] = TMath::Abs(pTrack->GetD(tPrimaryVtxPosition[0],
377 tPrimaryVtxPosition[1],
380 if (nTrack) tdcaDaughterToPrimVertex[1] = TMath::Abs(pTrack->GetD(tPrimaryVtxPosition[0],
381 tPrimaryVtxPosition[1],
384 v0->GetPPxPyPz(tMomPos[0],tMomPos[1],tMomPos[2]);
385 v0->GetNPxPyPz(tMomNeg[0],tMomNeg[1],tMomNeg[2]);
387 myAODv0->Fill(myAODVertex, lDcaV0Daughters, lDcaV0ToPrimVertex, tMomPos, tMomNeg, tdcaDaughterToPrimVertex);
388 myAODv0->SetOnFlyStatus(v0->GetOnFlyStatus());
391 lV0Radius = myAODv0->RadiusV0();
392 lDcaPosToPrimVertex = myAODv0->DcaPosToPrimVertex();
393 lDcaNegToPrimVertex = myAODv0->DcaNegToPrimVertex();
395 lOnFlyStatus = myAODv0->GetOnFlyStatus();
396 lChi2V0 = myAODv0->Chi2V0();
397 lDcaV0Daughters = myAODv0->DcaV0Daughters();
398 lDcaV0ToPrimVertex = myAODv0->DcaV0ToPrimVertex();
399 lV0CosineOfPointingAngle = myAODv0->CosPointingAngle(tPrimaryVtxPosition);
401 lInvMassK0s = myAODv0->MassK0Short();
402 lInvMassLambda = myAODv0->MassLambda();
403 lInvMassAntiLambda = myAODv0->MassAntiLambda();
405 lPt = TMath::Sqrt(myAODv0->Pt2V0());
406 lRapK0s = myAODv0->RapK0Short();
407 lRapLambda = myAODv0->RapLambda();
408 lAlphaV0 = myAODv0->AlphaV0();
409 lPtArmV0 = myAODv0->PtArmV0();
413 if (fOption.Contains("yes")) {
414 if ( (lDcaPosToPrimVertex < 0.036 )||
415 (lDcaPosToPrimVertex < 0.036 )||
416 (lDcaV0Daughters > 0.5 ) ||
417 (lV0CosineOfPointingAngle < 0.99)
421 // Filling histograms
422 fHistDcaPosToPrimVertex->Fill(lDcaPosToPrimVertex,lOnFlyStatus);
423 fHistDcaNegToPrimVertex->Fill(lDcaNegToPrimVertex,lOnFlyStatus);
424 fHistDcaPosToPrimVertexZoom->Fill(lDcaPosToPrimVertex,lOnFlyStatus);
425 fHistDcaNegToPrimVertexZoom->Fill(lDcaNegToPrimVertex,lOnFlyStatus);
426 fHistRadiusV0->Fill(lV0Radius,lOnFlyStatus);
427 fHistDecayLengthV0->Fill(myAODv0->DecayLengthV0(tPrimaryVtxPosition),lOnFlyStatus);
428 fHistDcaV0Daughters->Fill(lDcaV0Daughters,lOnFlyStatus);
429 fHistChi2->Fill(lChi2V0,lOnFlyStatus);
430 fHistCosPointAngle->Fill(lV0CosineOfPointingAngle,lOnFlyStatus);
431 if (lV0CosineOfPointingAngle >= 0.9) fHistCosPointAngleZoom->Fill(lV0CosineOfPointingAngle,lOnFlyStatus);
434 fHistPtVsYK0sOff->Fill(lPt,lRapK0s);
435 fHistPtVsYLambdaOff->Fill(lPt,lRapLambda);
436 fHistPtVsYAntiLambdaOff->Fill(lPt,lRapLambda);
437 fHistArmenterosPodolanskiOff->Fill(lAlphaV0,lPtArmV0);
441 fHistPtVsYK0sOn->Fill(lPt,lRapK0s);
442 fHistPtVsYLambdaOn->Fill(lPt,lRapLambda);
443 fHistPtVsYAntiLambdaOn->Fill(lPt,lRapLambda);
444 fHistArmenterosPodolanskiOn->Fill(lAlphaV0,lPtArmV0);
446 // K0s invariant mass histograms:
447 if (TMath::Abs(lRapK0s) < 1) {
449 fHistMassK0sOff->Fill(lInvMassK0s);
450 fHistMassVsRadiusK0sOff->Fill(lV0Radius,lInvMassK0s);
451 fHistPtVsMassK0sOff->Fill(lInvMassK0s,lPt);
454 fHistMassK0sOn->Fill(lInvMassK0s);
455 fHistMassVsRadiusK0sOn->Fill(lV0Radius,lInvMassK0s);
456 fHistPtVsMassK0sOn->Fill(lInvMassK0s,lPt);
459 // Lambda and AntiLambda invariant mass histograms:
460 if (TMath::Abs(lRapLambda) < 1) {
462 fHistMassLambdaOff->Fill(lInvMassLambda);
463 fHistMassAntiLambdaOff->Fill(lInvMassAntiLambda);
464 fHistMassVsRadiusLambdaOff->Fill(lV0Radius,lInvMassLambda);
465 fHistMassVsRadiusAntiLambdaOff->Fill(lV0Radius,lInvMassAntiLambda);
466 fHistPtVsMassLambdaOff->Fill(lInvMassLambda,lPt);
467 fHistPtVsMassAntiLambdaOff->Fill(lInvMassAntiLambda,lPt);
470 fHistMassLambdaOn->Fill(lInvMassLambda);
471 fHistMassAntiLambdaOn->Fill(lInvMassAntiLambda);
472 fHistMassVsRadiusLambdaOn->Fill(lV0Radius,lInvMassLambda);
473 fHistMassVsRadiusAntiLambdaOn->Fill(lV0Radius,lInvMassAntiLambda);
474 fHistPtVsMassLambdaOn->Fill(lInvMassLambda,lPt);
475 fHistPtVsMassAntiLambdaOn->Fill(lInvMassAntiLambda,lPt);
480 if (primary) delete primary;
481 if (myPosAodTrack) delete myPosAodTrack;
482 if (myNegAodTrack) delete myNegAodTrack;
483 if (myAODVertex) delete myAODVertex;
484 if (myAODv0) delete myAODv0;
486 else if(fAnalysisType == "AOD") {
488 const AliAODVertex *primaryVtx = ((AliAODEvent*)lEvent)->GetPrimaryVertex();
489 tPrimaryVtxPosition[0] = primaryVtx->GetX();
490 tPrimaryVtxPosition[1] = primaryVtx->GetY();
491 tPrimaryVtxPosition[2] = primaryVtx->GetZ();
493 fHistPrimaryVertexPosX->Fill(tPrimaryVtxPosition[0]);
494 fHistPrimaryVertexPosY->Fill(tPrimaryVtxPosition[1]);
495 fHistPrimaryVertexPosZ->Fill(tPrimaryVtxPosition[2]);
497 for (Int_t iV0 = 0; iV0 < nv0s; iV0++)
498 {// This is the begining of the V0 loop
499 AliAODv0 *v0 = ((AliAODEvent*)lEvent)->GetV0(iV0);
503 lV0Radius = v0->RadiusV0();
504 lDcaPosToPrimVertex = v0->DcaPosToPrimVertex();
505 lDcaNegToPrimVertex = v0->DcaNegToPrimVertex();
507 lOnFlyStatus = v0->GetOnFlyStatus();
508 lChi2V0 = v0->Chi2V0();
509 lDcaV0Daughters = v0->DcaV0Daughters();
510 lDcaV0ToPrimVertex = v0->DcaV0ToPrimVertex();
511 lV0CosineOfPointingAngle = v0->CosPointingAngle(tPrimaryVtxPosition);
513 lInvMassK0s = v0->MassK0Short();
514 lInvMassLambda = v0->MassLambda();
515 lInvMassAntiLambda = v0->MassAntiLambda();
517 lPt = TMath::Sqrt(v0->Pt2V0());
518 lRapK0s = v0->RapK0Short();
519 lRapLambda = v0->RapLambda();
520 lAlphaV0 = v0->AlphaV0();
521 lPtArmV0 = v0->PtArmV0();
525 if (fOption.Contains("yes")) {
526 if ( (lDcaPosToPrimVertex < 0.036 )||
527 (lDcaPosToPrimVertex < 0.036 )||
528 (lDcaV0Daughters > 0.5 ) ||
529 (lV0CosineOfPointingAngle < 0.99)
533 // Filling histograms
534 fHistDcaPosToPrimVertex->Fill(lDcaPosToPrimVertex,lOnFlyStatus);
535 fHistDcaNegToPrimVertex->Fill(lDcaNegToPrimVertex,lOnFlyStatus);
536 fHistDcaPosToPrimVertexZoom->Fill(lDcaPosToPrimVertex,lOnFlyStatus);
537 fHistDcaNegToPrimVertexZoom->Fill(lDcaNegToPrimVertex,lOnFlyStatus);
538 fHistRadiusV0->Fill(lV0Radius,lOnFlyStatus);
539 fHistDecayLengthV0->Fill(v0->DecayLengthV0(tPrimaryVtxPosition),lOnFlyStatus);
540 fHistDcaV0Daughters->Fill(lDcaV0Daughters,lOnFlyStatus);
541 fHistChi2->Fill(lChi2V0,lOnFlyStatus);
542 fHistCosPointAngle->Fill(lV0CosineOfPointingAngle,lOnFlyStatus);
543 if (lV0CosineOfPointingAngle >= 0.9) fHistCosPointAngleZoom->Fill(lV0CosineOfPointingAngle,lOnFlyStatus);
546 fHistPtVsYK0sOff->Fill(lPt,lRapK0s);
547 fHistPtVsYLambdaOff->Fill(lPt,lRapLambda);
548 fHistPtVsYAntiLambdaOff->Fill(lPt,lRapLambda);
549 fHistArmenterosPodolanskiOff->Fill(lAlphaV0,lPtArmV0);
553 fHistPtVsYK0sOn->Fill(lPt,lRapK0s);
554 fHistPtVsYLambdaOn->Fill(lPt,lRapLambda);
555 fHistPtVsYAntiLambdaOn->Fill(lPt,lRapLambda);
556 fHistArmenterosPodolanskiOn->Fill(lAlphaV0,lPtArmV0);
558 // K0s invariant mass histograms:
559 if (TMath::Abs(lRapK0s) < 1) {
561 fHistMassK0sOff->Fill(lInvMassK0s);
562 fHistMassVsRadiusK0sOff->Fill(lV0Radius,lInvMassK0s);
563 fHistPtVsMassK0sOff->Fill(lInvMassK0s,lPt);
566 fHistMassK0sOn->Fill(lInvMassK0s);
567 fHistMassVsRadiusK0sOn->Fill(lV0Radius,lInvMassK0s);
568 fHistPtVsMassK0sOn->Fill(lInvMassK0s,lPt);
571 // Lambda and AntiLambda invariant mass histograms:
572 if (TMath::Abs(lRapLambda) < 1) {
574 fHistMassLambdaOff->Fill(lInvMassLambda);
575 fHistMassAntiLambdaOff->Fill(lInvMassAntiLambda);
576 fHistMassVsRadiusLambdaOff->Fill(lV0Radius,lInvMassLambda);
577 fHistMassVsRadiusAntiLambdaOff->Fill(lV0Radius,lInvMassAntiLambda);
578 fHistPtVsMassLambdaOff->Fill(lInvMassLambda,lPt);
579 fHistPtVsMassAntiLambdaOff->Fill(lInvMassAntiLambda,lPt);
582 fHistMassLambdaOn->Fill(lInvMassLambda);
583 fHistMassAntiLambdaOn->Fill(lInvMassAntiLambda);
584 fHistMassVsRadiusLambdaOn->Fill(lV0Radius,lInvMassLambda);
585 fHistMassVsRadiusAntiLambdaOn->Fill(lV0Radius,lInvMassAntiLambda);
586 fHistPtVsMassLambdaOn->Fill(lInvMassLambda,lPt);
587 fHistPtVsMassAntiLambdaOn->Fill(lInvMassAntiLambda,lPt);
593 fHistV0Multiplicity->Fill(nv0s);
594 fHistV0MultiplicityOff->Fill(nv0sOff);
595 fHistV0MultiplicityOn->Fill(nv0sOn);
598 PostData(1, fListHist);
601 //________________________________________________________________________
602 void AliAnalysisTaskStrange::Terminate(Option_t *)
604 // Draw result to the screen
605 // Called once at the end of the query