1 /**************************************************************************
2 * Author: Boris Hippolyte. *
3 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
14 //-----------------------------------------------------------------
15 // AliAnalysisTaskESDCheckV0 class
16 // This task is for QAing the dE/dx from the ESD
17 // Origin: B.H. Nov2007, hippolyt@in2p3.fr
18 //-----------------------------------------------------------------
20 //-----------------------------------------------------------------
22 // 1) extent the use of fMidPseudoRapidityFlag to primaries
23 // 2) discuss the GetSigmaToVertex(AliESDtrack*)
24 // 3) decide the kind of refit to be used
25 //-----------------------------------------------------------------
35 #include "AliAnalysisTask.h"
36 #include "AliAnalysisManager.h"
38 #include "AliESDEvent.h"
39 #include "AliESDInputHandler.h"
41 #include "AliESDVertex.h"
44 #include "AliAnalysisTaskESDDedx.h"
46 ClassImp(AliAnalysisTaskESDDedx)
48 //________________________________________________________________________
49 AliAnalysisTaskESDDedx::AliAnalysisTaskESDDedx(const char *rName,
50 const Bool_t rAllConstrainedFlag,
51 const Bool_t rMidPseudoRapidityFlag,
52 const Bool_t rSelTrackRemoveKink,
53 const Bool_t rSelTrackWithOnTheFlyV0,
54 const Int_t rSelTrackMinClustersTPC,
55 const Int_t rSelTrackMinClustersITS,
56 const Float_t rSelTrackMaxChi2PerClusterTPC,
57 const Float_t rSelTrackMaxChi2PerClusterITS,
58 const Double_t rSelTrackMaxCov11,
59 const Double_t rSelTrackMaxCov22,
60 const Double_t rSelTrackMaxCov33,
61 const Double_t rSelTrackMaxCov44,
62 const Double_t rSelTrackMaxCov55,
63 const Double_t rSelV0MaxDcaDaughters,
64 const Double_t rSelV0MinDecayLength)
66 : AliAnalysisTask(rName, ""), fESD(0), fListHist(), fHistPtot(0),
68 fHistMultiplicity(0), fHistTPCDedxVsMomentum(0), fHistITSDedxVsMomentum(0),
69 fHistMassK0(0), fHistMassLambda(0), fHistMassAntiLambda(0),
70 fHistTPCDedxVsMomPosK0(0), fHistTPCDedxVsMomNegK0(0),
71 fHistTPCDedxVsMomPosLambda(0), fHistTPCDedxVsMomNegLambda(0),
72 fHistTPCDedxVsMomPosAntiLambda(0), fHistTPCDedxVsMomNegAntiLambda(0),
73 fHistDiffInOutMomentum(0), fHistDiffPrimOutMomentum(0),
74 fHistDiffPrimMeanMomentum(0), fHistPercPrimMeanMomentum(0), fHistPrimEta(0),
75 fHistPercPrimMeanMomentumVsEta(0), fHistPercPrimMeanMomentumVsPrim(0),
77 fHistMultiplicityCuts(0), fHistTPCDedxVsMomentumCuts(0), fHistITSDedxVsMomentumCuts(0),
78 fHistMassK0Cuts(0), fHistMassLambdaCuts(0), fHistMassAntiLambdaCuts(0),
79 fHistTPCDedxVsMomPosK0Cuts(0), fHistTPCDedxVsMomNegK0Cuts(0),
80 fHistTPCDedxVsMomPosLambdaCuts(0), fHistTPCDedxVsMomNegLambdaCuts(0),
81 fHistTPCDedxVsMomPosAntiLambdaCuts(0), fHistTPCDedxVsMomNegAntiLambdaCuts(0),
82 fHistDiffInOutMomentumCuts(0), fHistDiffPrimOutMomentumCuts(0),
83 fHistDiffPrimMeanMomentumCuts(0), fHistPercPrimMeanMomentumCuts(0), fHistPrimEtaCuts(0),
84 fHistPercPrimMeanMomentumVsEtaCuts(0), fHistPercPrimMeanMomentumVsPrimCuts(0),
86 fAllConstrainedFlag(rAllConstrainedFlag), fMidPseudoRapidityFlag(rMidPseudoRapidityFlag),
88 fSelTrackRemoveKink(rSelTrackRemoveKink), fSelTrackWithOnTheFlyV0(rSelTrackWithOnTheFlyV0),
89 fSelTrackMinClustersTPC(rSelTrackMinClustersTPC), fSelTrackMinClustersITS(rSelTrackMinClustersITS),
90 fSelTrackMaxChi2PerClusterTPC(rSelTrackMaxChi2PerClusterTPC), fSelTrackMaxChi2PerClusterITS(rSelTrackMaxChi2PerClusterITS),
91 fSelTrackMaxCov11(rSelTrackMaxCov11), fSelTrackMaxCov22(rSelTrackMaxCov22), fSelTrackMaxCov33(rSelTrackMaxCov33),
92 fSelTrackMaxCov44(rSelTrackMaxCov44), fSelTrackMaxCov55(rSelTrackMaxCov55),
93 fSelV0MaxDcaDaughters(rSelV0MaxDcaDaughters), fSelV0MinDecayLength(rSelV0MinDecayLength)
97 // Define input and output slots here
98 // Input slot #0 works with a TChain
99 DefineInput(0, TChain::Class());
100 // Output slot #0 writes into a TH1 container
101 DefineOutput(0, TH1F::Class());
102 // Output slot #1 writes into a TList container
103 DefineOutput(1, TList::Class());
106 //________________________________________________________________________
107 void AliAnalysisTaskESDDedx::ConnectInputData(Option_t *)
109 // Connect ESD or AOD here
112 TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
114 Printf("ERROR: Could not read chain from input slot 0");
116 // Disable all branches and enable only the needed ones
117 // The next two lines are different when data produced as AliESDEvent is read
118 tree->SetBranchStatus("*", kFALSE);
119 tree->SetBranchStatus("fTracks.*", kTRUE);
120 tree->SetBranchStatus("fV0s.*", kTRUE);
122 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
125 Printf("ERROR: Could not get ESDInputHandler");
127 fESD = esdH->GetEvent();
131 //________________________________________________________________________
132 void AliAnalysisTaskESDDedx::CreateOutputObjects()
136 fListHist = new TList();
138 fHistPtot = new TH1F("fHistPtot", "P_{tot} distribution;P_{tot} (GeV/c);dN/dP_{tot} (c/GeV)", 15, 0.1, 3.1);
139 fHistPtot->SetMarkerStyle(kFullCircle);
140 fListHist->Add(fHistPtot);
142 if (!fHistMultiplicity) {
143 fHistMultiplicity = new TH1F("fHistMultiplicity", "Multiplicity distribution;Number of tracks;Events", 250, 0, 250);
144 fHistMultiplicity->SetMarkerStyle(kFullCircle);
145 fListHist->Add(fHistMultiplicity);
147 if (!fHistTPCDedxVsMomentum) {
148 fHistTPCDedxVsMomentum = new TH2F("h2TPCDedxVsMomentum","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
149 fListHist->Add(fHistTPCDedxVsMomentum);
151 if (!fHistITSDedxVsMomentum) {
152 fHistITSDedxVsMomentum = new TH2F("h2ITSDedxVsMomentum","Bethe-Bloch Distribution for ITS;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
153 fListHist->Add(fHistITSDedxVsMomentum);
156 fHistMassK0 = new TH1F("h1MassK0","K^{0} candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});Counts",100,0.4,0.6);
157 fListHist->Add(fHistMassK0);
159 if (!fHistMassLambda) {
160 fHistMassLambda = new TH1F("h1MassLambda","#Lambda^{0} candidates;M(p#pi^{-}) (GeV/c^{2});Counts",75,1.05,1.2);
161 fListHist->Add(fHistMassLambda);
163 if (!fHistMassAntiLambda) {
164 fHistMassAntiLambda = new TH1F("h1MassAntiLambda","#bar{#Lambda}^{0} candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts",75,1.05,1.2);
165 fListHist->Add(fHistMassAntiLambda);
167 if (!fHistTPCDedxVsMomPosK0) {
168 fHistTPCDedxVsMomPosK0 = new TH2F("h2TPCDedxVsMomPosK0","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
169 fListHist->Add(fHistTPCDedxVsMomPosK0);
171 if (!fHistTPCDedxVsMomNegK0) {
172 fHistTPCDedxVsMomNegK0 = new TH2F("h2TPCDedxVsMomNegK0","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
173 fListHist->Add(fHistTPCDedxVsMomNegK0);
175 if (!fHistTPCDedxVsMomPosLambda) {
176 fHistTPCDedxVsMomPosLambda = new TH2F("h2TPCDedxVsMomPosLambda","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
177 fListHist->Add(fHistTPCDedxVsMomPosLambda);
179 if (!fHistTPCDedxVsMomNegLambda) {
180 fHistTPCDedxVsMomNegLambda = new TH2F("h2TPCDedxVsMomNegLambda","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
181 fListHist->Add(fHistTPCDedxVsMomNegLambda);
183 if (!fHistTPCDedxVsMomPosAntiLambda) {
184 fHistTPCDedxVsMomPosAntiLambda = new TH2F("h2TPCDedxVsMomPosAntiLambda","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
185 fListHist->Add(fHistTPCDedxVsMomPosAntiLambda);
187 if (!fHistTPCDedxVsMomNegAntiLambda) {
188 fHistTPCDedxVsMomNegAntiLambda = new TH2F("h2TPCDedxVsMomNegAntiLambda","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
189 fListHist->Add(fHistTPCDedxVsMomNegAntiLambda);
191 if (!fHistDiffInOutMomentum) {
192 fHistDiffInOutMomentum = new TH1F("h1DiffInOutMomentum","Momentum Difference In-Out;Momentum (GeV/c);Counts",200,-0.1,0.1);
193 fListHist->Add(fHistDiffInOutMomentum);
195 if (!fHistDiffPrimOutMomentum) {
196 fHistDiffPrimOutMomentum = new TH1F("h1DiffPrimOutMomentum","Momentum Difference Prim-Out;Momentum (GeV/c);Counts",200,-0.1,0.1);
197 fListHist->Add(fHistDiffPrimOutMomentum);
199 if (!fHistDiffPrimMeanMomentum) {
200 fHistDiffPrimMeanMomentum = new TH1F("h1DiffPrimMeanMomentum","Momentum Difference Prim-Mean;Momentum (GeV/c);Counts",200,-0.1,0.1);
201 fListHist->Add(fHistDiffPrimMeanMomentum);
203 if (!fHistPercPrimMeanMomentum) {
204 fHistPercPrimMeanMomentum = new TH1F("h1PercPrimMeanMomentum","Momentum Percentage (Prim-Mean)/Prim;(Primary-Mean)/Primary (%);Counts",200,-10,10);
205 fListHist->Add(fHistPercPrimMeanMomentum);
208 fHistPrimEta = new TH1F("h1PrimEta","Pseudorapidity Distribution (Primaries);#eta;Counts",300,-1.5,1.5);
209 fListHist->Add(fHistPrimEta);
211 if (!fHistPercPrimMeanMomentumVsEta) {
212 fHistPercPrimMeanMomentumVsEta = new TH2F("h2PercPrimMeanMomentumVsEta","Momentum Percentage (Prim-Mean)/Prim vs #eta;#eta;(Primary-Mean)/Primary (%)",100,-1.5,1.5,100,-10,10);
213 fListHist->Add(fHistPercPrimMeanMomentumVsEta);
215 if (!fHistPercPrimMeanMomentumVsPrim) {
216 fHistPercPrimMeanMomentumVsPrim = new TH2F("h2PercPrimMeanMomentumVsPrim","Momentum Percentage (Prim-Mean)/Prim vs Prim;Momentum (GeV/c);(Primary-Mean)/Primary (%)",1500,0,15,100,-10,10);
217 fListHist->Add(fHistPercPrimMeanMomentumVsPrim);
219 // Histograms after selections
220 if (!fHistMultiplicityCuts) {
221 fHistMultiplicityCuts = new TH1F("h1HistMultiplicityCuts","Number of Good tracks;Number of tracks;Event",201,-0.5,200.5);
222 fListHist->Add(fHistMultiplicityCuts);
224 if (!fHistTPCDedxVsMomentumCuts) {
225 fHistTPCDedxVsMomentumCuts = new TH2F("h2TPCDedxVsMomentumCuts","Bethe-Bloch Distribution for TPC w/cuts;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
226 fListHist->Add(fHistTPCDedxVsMomentumCuts);
228 if (!fHistITSDedxVsMomentumCuts) {
229 fHistITSDedxVsMomentumCuts = new TH2F("h2ITSDedxVsMomentumCuts","Bethe-Bloch Distribution for ITS w/cuts;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
230 fListHist->Add(fHistITSDedxVsMomentumCuts);
232 if (!fHistMassK0Cuts) {
233 fHistMassK0Cuts = new TH1F("h1MassK0Cuts","K^{0} candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});Counts",100,0.4,0.6);
234 fListHist->Add(fHistMassK0Cuts);
236 if (!fHistMassLambdaCuts) {
237 fHistMassLambdaCuts = new TH1F("h1MassLambdaCuts","#Lambda^{0} candidates;M(p#pi^{-}) (GeV/c^{2});Counts",75,1.05,1.2);
238 fListHist->Add(fHistMassLambdaCuts);
240 if (!fHistMassAntiLambdaCuts) {
241 fHistMassAntiLambdaCuts = new TH1F("h1MassAntiLambdaCuts","#bar{#Lambda}^{0} candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts",75,1.05,1.2);
242 fListHist->Add(fHistMassAntiLambdaCuts);
244 if (!fHistTPCDedxVsMomPosK0Cuts) {
245 fHistTPCDedxVsMomPosK0Cuts = new TH2F("h2TPCDedxVsMomPosK0Cuts","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
246 fListHist->Add(fHistTPCDedxVsMomPosK0Cuts);
248 if (!fHistTPCDedxVsMomNegK0Cuts) {
249 fHistTPCDedxVsMomNegK0Cuts = new TH2F("h2TPCDedxVsMomNegK0Cuts","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
250 fListHist->Add(fHistTPCDedxVsMomNegK0Cuts);
252 if (!fHistTPCDedxVsMomPosLambdaCuts) {
253 fHistTPCDedxVsMomPosLambdaCuts = new TH2F("h2TPCDedxVsMomPosLambdaCuts","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
254 fListHist->Add(fHistTPCDedxVsMomPosLambdaCuts);
256 if (!fHistTPCDedxVsMomNegLambdaCuts) {
257 fHistTPCDedxVsMomNegLambdaCuts = new TH2F("h2TPCDedxVsMomNegLambdaCuts","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
258 fListHist->Add(fHistTPCDedxVsMomNegLambdaCuts);
260 if (!fHistTPCDedxVsMomPosAntiLambdaCuts) {
261 fHistTPCDedxVsMomPosAntiLambdaCuts = new TH2F("h2TPCDedxVsMomPosAntiLambdaCuts","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
262 fListHist->Add(fHistTPCDedxVsMomPosAntiLambdaCuts);
264 if (!fHistTPCDedxVsMomNegAntiLambdaCuts) {
265 fHistTPCDedxVsMomNegAntiLambdaCuts = new TH2F("h2TPCDedxVsMomNegAntiLambdaCuts","Bethe-Bloch Distribution for TPC;Momentum (GeV/c);dE/dx (A.U.)",1500,0,15,100,0,100);
266 fListHist->Add(fHistTPCDedxVsMomNegAntiLambdaCuts);
268 if (!fHistDiffInOutMomentumCuts) {
269 fHistDiffInOutMomentumCuts = new TH1F("h1DiffInOutMomentumCuts","Momentum Difference In-Out;Momentum (GeV/c);Counts",200,-0.1,0.1);
270 fListHist->Add(fHistDiffInOutMomentumCuts);
272 if (!fHistDiffPrimOutMomentumCuts) {
273 fHistDiffPrimOutMomentumCuts = new TH1F("h1DiffPrimOutMomentumCuts","Momentum Difference Prim-Out;Momentum (GeV/c);Counts",200,-0.1,0.1);
274 fListHist->Add(fHistDiffPrimOutMomentumCuts);
276 if (!fHistDiffPrimMeanMomentumCuts) {
277 fHistDiffPrimMeanMomentumCuts = new TH1F("h1DiffPrimMeanMomentumCuts","Momentum Difference Prim-Mean;Momentum (GeV/c);Counts",200,-0.1,0.1);
278 fListHist->Add(fHistDiffPrimMeanMomentumCuts);
280 if (!fHistPercPrimMeanMomentumCuts) {
281 fHistPercPrimMeanMomentumCuts = new TH1F("h1PercPrimMeanMomentumCuts","Momentum Percentage (Prim-Mean)/Prim;(Primary-Mean)/Primary (%);Counts",200,-10,10);
282 fListHist->Add(fHistPercPrimMeanMomentumCuts);
284 if (!fHistPrimEtaCuts) {
285 fHistPrimEtaCuts = new TH1F("h1PrimEtaCuts","Pseudorapidity Distribution (Primaries);#eta;Counts",300,-1.5,1.5);
286 fListHist->Add(fHistPrimEtaCuts);
288 if (!fHistPercPrimMeanMomentumVsEtaCuts) {
289 fHistPercPrimMeanMomentumVsEtaCuts = new TH2F("h2PercPrimMeanMomentumVsEtaCuts","Momentum Percentage (Prim-Mean)/Prim vs #eta;#eta;(Primary-Mean)/Primary (%)",100,-1.5,1.5,100,-10,10);
290 fListHist->Add(fHistPercPrimMeanMomentumVsEtaCuts);
292 if (!fHistPercPrimMeanMomentumVsPrimCuts) {
293 fHistPercPrimMeanMomentumVsPrimCuts = new TH2F("h2PercPrimMeanMomentumVsPrimCuts","Momentum Percentage (Prim-Mean)/Prim vs Prim;Momentum (GeV/c);(Primary-Mean)/Primary (%)",1500,0,15,100,-10,10);
294 fListHist->Add(fHistPercPrimMeanMomentumVsPrimCuts);
298 //________________________________________________________________________
299 void AliAnalysisTaskESDDedx::Exec(Option_t *)
302 // Called for each event
305 Printf("ERROR: fESD not available");
309 Int_t nGetTracks = fESD->GetNumberOfTracks();
310 // Printf("There are %d tracks in this event", nGetTracks);
312 Int_t nAllTracks = 0, nGoodTracks = 0;
313 Double_t lTPCDedx = 0, lITSDedx = 0;
316 for (Int_t iTracks = 0; iTracks < nGetTracks; iTracks++) {
317 AliESDtrack* track = fESD->GetTrack(iTracks);
319 Printf("ERROR: Could not retrieve track %d", iTracks);
322 fHistPtot->Fill(track->P()); // Fill a ptot spectrum before any selection
324 Double_t lPrimvtxMomemtum[3];
325 Bool_t lIsTrackConstrained = track->GetConstrainedPxPyPz(lPrimvtxMomemtum);
326 if (fAllConstrainedFlag && (!lIsTrackConstrained)) continue; // Constrained or not to Prim Vertex.
328 const AliExternalTrackParam *lInnerETP = track->GetInnerParam();
329 const AliExternalTrackParam *lOuterETP = track->GetOuterParam();
330 if ((!lInnerETP)||(!lOuterETP)) continue; // No inner nor outer params for this track.
332 Double_t lInnerMomentum = 0, lOuterMomentum = 0, lMeanMomentum = 0;
333 if (lInnerETP) lInnerMomentum = lInnerETP->GetP(); // The "inner" momentum.
334 if (lOuterETP) lOuterMomentum = lOuterETP->GetP(); // The "outer" momentum.
335 if (lInnerETP&&lOuterETP) lMeanMomentum = (lInnerMomentum+lOuterMomentum)/2.; // The "mean" momentum
337 lTPCDedx = track->GetTPCsignal();
338 lITSDedx = track->GetITSsignal();
339 // This needs to be discussed and implemented
340 // nSigmaToVertex = GetSigmaToVertex(track);
342 Float_t lPrimPseudoRap = 0, lPrimMomentum = 0;
343 if (lIsTrackConstrained){
344 TVector3 pPrimvtx(lPrimvtxMomemtum);
345 lPrimPseudoRap = pPrimvtx.Eta();
346 lPrimMomentum = pPrimvtx.Mag();
350 ((!fAllConstrainedFlag)||(TMath::Abs(lPrimPseudoRap) < 0.5)) &&
353 // fAllConstrained = 0 -> !fAllConstrained = 1 so ok anyway
354 // fAllConstrained = 1 -> !fAllConstrained = 0 so it depends on eta condition
357 fHistTPCDedxVsMomentum->Fill(lInnerMomentum,lTPCDedx);
358 fHistDiffInOutMomentum->Fill(lInnerMomentum-lOuterMomentum);
359 if (lIsTrackConstrained){
360 fHistITSDedxVsMomentum->Fill(lPrimMomentum,lITSDedx);
361 fHistDiffPrimOutMomentum->Fill(lPrimMomentum-lOuterMomentum);
362 fHistDiffPrimMeanMomentum->Fill(lPrimMomentum-lMeanMomentum);
364 fHistPercPrimMeanMomentum->Fill(100.*(lPrimMomentum-lMeanMomentum)/lPrimMomentum);
365 fHistPercPrimMeanMomentumVsEta->Fill(lPrimPseudoRap,100.*(lPrimMomentum-lMeanMomentum)/lPrimMomentum);
366 fHistPercPrimMeanMomentumVsPrim->Fill(lPrimMomentum,100.*(lPrimMomentum-lMeanMomentum)/lPrimMomentum);
368 fHistPrimEta->Fill(lPrimPseudoRap);
371 // (!fAllConstrainedFlag || (nSigmaToVertex < 3)) &&
374 fHistTPCDedxVsMomentumCuts->Fill(lInnerMomentum,lTPCDedx);
375 fHistDiffInOutMomentumCuts->Fill(lInnerMomentum-lOuterMomentum);
376 if (lIsTrackConstrained){
377 fHistITSDedxVsMomentumCuts->Fill(lPrimMomentum,lITSDedx);
378 fHistDiffPrimOutMomentumCuts->Fill(lPrimMomentum-lOuterMomentum);
379 fHistDiffPrimMeanMomentumCuts->Fill(lPrimMomentum-lMeanMomentum);
381 fHistPercPrimMeanMomentumCuts->Fill(100.*(lPrimMomentum-lMeanMomentum)/lPrimMomentum);
382 fHistPercPrimMeanMomentumVsEtaCuts->Fill(lPrimPseudoRap,100.*(lPrimMomentum-lMeanMomentum)/lPrimMomentum);
383 fHistPercPrimMeanMomentumVsPrimCuts->Fill(lPrimMomentum,100.*(lPrimMomentum-lMeanMomentum)/lPrimMomentum);
385 fHistPrimEtaCuts->Fill(lPrimPseudoRap);
391 fHistMultiplicity->Fill(nGetTracks);
392 fHistMultiplicityCuts->Fill(nGoodTracks);
395 nv0s = fESD->GetNumberOfV0s();
397 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
398 Double_t lTPCDedxPos = 0, lTPCDedxNeg = 0;
399 Float_t lInvMassK0 = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;
400 for (Int_t iV0 = 0; iV0 < nv0s; iV0++)
401 {// This is the V0 loop
402 AliESDv0 *v0 = fESD->GetV0(iV0);
404 if (v0->GetOnFlyStatus() != fSelTrackWithOnTheFlyV0) continue;
406 const AliESDVertex *lPrimaryVertex = fESD->GetVertex();
407 Double_t tPositionPrimaryVertex[3]; lPrimaryVertex->GetXYZ(tPositionPrimaryVertex);
408 Double_t v0PositionX = 0, v0PositionY = 0, v0PositionZ = 0; v0->GetXYZ(v0PositionX,v0PositionY,v0PositionZ);
409 Double_t v0DecayLength = TMath::Sqrt(
410 (v0PositionX-tPositionPrimaryVertex[0])*(v0PositionX-tPositionPrimaryVertex[0])+
411 (v0PositionY-tPositionPrimaryVertex[1])*(v0PositionY-tPositionPrimaryVertex[1])+
412 (v0PositionZ-tPositionPrimaryVertex[2])*(v0PositionZ-tPositionPrimaryVertex[2])
414 if (v0->GetDcaV0Daughters() > fSelV0MaxDcaDaughters) continue;
415 if (v0DecayLength < fSelV0MinDecayLength) continue;
417 // Getting invariant mass infos directly from ESD
418 v0->ChangeMassHypothesis(310);
419 lInvMassK0 = v0->GetEffMass();
420 v0->ChangeMassHypothesis(3122);
421 lInvMassLambda = v0->GetEffMass();
422 v0->ChangeMassHypothesis(-3122);
423 lInvMassAntiLambda = v0->GetEffMass();
425 // Filling invariant mass histos for all candidates
426 fHistMassK0->Fill(lInvMassK0);
427 fHistMassLambda->Fill(lInvMassLambda);
428 fHistMassAntiLambda->Fill(lInvMassAntiLambda);
430 // Accessing the daughter track infos
431 lIndexTrackPos = TMath::Abs(v0->GetPindex());
432 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
433 AliESDtrack* trackPos = fESD->GetTrack(lIndexTrackPos);
434 AliESDtrack* trackNeg = fESD->GetTrack(lIndexTrackNeg);
436 Double_t lPrimvtxMomPos[3];
437 Bool_t lIsPosConstrained = trackPos->GetConstrainedPxPyPz(lPrimvtxMomPos);
438 if (fAllConstrainedFlag && (!lIsPosConstrained)) continue; // Constrained or not to Prim Vertex.
440 Double_t lPrimvtxMomNeg[3];
441 Bool_t lIsNegConstrained = trackNeg->GetConstrainedPxPyPz(lPrimvtxMomNeg);
442 if (fAllConstrainedFlag && (!lIsNegConstrained)) continue; // Constrained or not to Prim Vertex.
444 const AliExternalTrackParam *lInnerPosETP = trackPos->GetInnerParam();
445 const AliExternalTrackParam *lInnerNegETP = trackNeg->GetInnerParam();
446 if ((!lInnerPosETP)||(!lInnerNegETP)) continue; // No inner params for at least one track.
448 Double_t lInnerMomPos = 0, lInnerMomNeg = 0;
449 if (lInnerPosETP) lInnerMomPos = lInnerPosETP->GetP();
450 if (lInnerNegETP) lInnerMomNeg = lInnerNegETP->GetP();
452 lTPCDedxPos = trackPos->GetTPCsignal();
453 lTPCDedxNeg = trackNeg->GetTPCsignal();
455 Float_t lPrimPseudoRapPos = 0, lPrimPseudoRapNeg = 0;
456 if (lIsPosConstrained && lIsNegConstrained){
457 TVector3 pPrimPos(lPrimvtxMomPos);
458 TVector3 pPrimNeg(lPrimvtxMomNeg);
459 lPrimPseudoRapPos = pPrimPos.Eta();
460 lPrimPseudoRapNeg = pPrimNeg.Eta();
463 if (TMath::Abs(lInvMassK0-0.497)<0.01) {
464 fHistMassK0Cuts->Fill(lInvMassK0);
465 fHistTPCDedxVsMomPosK0->Fill(lInnerMomPos,lTPCDedxPos);
466 fHistTPCDedxVsMomNegK0->Fill(lInnerMomNeg,lTPCDedxNeg);
468 (!fMidPseudoRapidityFlag) || ( TMath::Abs(lPrimPseudoRapPos) < 0.5) &&
469 (IsAccepted(trackPos))
471 fHistTPCDedxVsMomPosK0Cuts->Fill(lInnerMomPos,lTPCDedxPos);
473 (!fMidPseudoRapidityFlag) || ( TMath::Abs(lPrimPseudoRapNeg) < 0.5) &&
474 (IsAccepted(trackNeg))
476 fHistTPCDedxVsMomNegK0Cuts->Fill(lInnerMomNeg,lTPCDedxNeg);
478 if (TMath::Abs(lInvMassLambda-1.115)<0.01) {
479 fHistMassLambdaCuts->Fill(lInvMassLambda);
480 fHistTPCDedxVsMomPosLambda->Fill(lInnerMomPos,lTPCDedxPos);
481 fHistTPCDedxVsMomNegLambda->Fill(lInnerMomNeg,lTPCDedxNeg);
483 (!fMidPseudoRapidityFlag) || ( TMath::Abs(lPrimPseudoRapPos) < 0.5) &&
484 (IsAccepted(trackPos))
486 fHistTPCDedxVsMomPosLambdaCuts->Fill(lInnerMomPos,lTPCDedxPos);
488 (!fMidPseudoRapidityFlag) || ( TMath::Abs(lPrimPseudoRapNeg) < 0.5) &&
489 (IsAccepted(trackNeg))
491 fHistTPCDedxVsMomNegLambdaCuts->Fill(lInnerMomNeg,lTPCDedxNeg);
493 if (TMath::Abs(lInvMassAntiLambda-1.115)<0.01) {
494 fHistMassAntiLambdaCuts->Fill(lInvMassAntiLambda);
495 fHistTPCDedxVsMomPosAntiLambda->Fill(lInnerMomPos,lTPCDedxPos);
496 fHistTPCDedxVsMomNegAntiLambda->Fill(lInnerMomNeg,lTPCDedxNeg);
498 (!fMidPseudoRapidityFlag) || ( TMath::Abs(lPrimPseudoRapPos) < 0.5) &&
499 (IsAccepted(trackPos))
501 fHistTPCDedxVsMomPosAntiLambdaCuts->Fill(lInnerMomPos,lTPCDedxPos);
503 (!fMidPseudoRapidityFlag) || ( TMath::Abs(lPrimPseudoRapNeg) < 0.5) &&
504 (IsAccepted(trackNeg))
506 fHistTPCDedxVsMomNegAntiLambdaCuts->Fill(lInnerMomNeg,lTPCDedxNeg);
511 PostData(0, fHistPtot);
512 PostData(1, fListHist);
515 //________________________________________________________________________
516 void AliAnalysisTaskESDDedx::Terminate(Option_t *)
518 // Draw result to the screen
519 // Called once at the end of the query
521 Printf("Reminder for options and selections: fAllConstrainedFlag = %d fMidPseudoRapidityFlag = %d fSelTrackRemoveKink = %d fSelTrackWithOnTheFlyV0 = %d fSelTrackMinClustersTPC = %d fSelTrackMinClustersITS = %d fSelTrackMaxChi2PerClusterTPC = %.2f fSelTrackMaxChi2PerClusterITS = %.2f fSelTrackMaxCov11 = %.2f fSelTrackMaxCov22 = %.2f fSelTrackMaxCov33 = %.2f fSelTrackMaxCov44 = %.2f fSelTrackMaxCov55 = %.2f fSelV0MaxDcaDaughters = %.2f fSelV0MinDecayLength= %.2f",
523 fMidPseudoRapidityFlag,fSelTrackRemoveKink,fSelTrackWithOnTheFlyV0,
524 fSelTrackMinClustersTPC,fSelTrackMinClustersITS,
525 fSelTrackMaxChi2PerClusterTPC,fSelTrackMaxChi2PerClusterITS,
526 fSelTrackMaxCov11,fSelTrackMaxCov22,fSelTrackMaxCov33,fSelTrackMaxCov44,fSelTrackMaxCov55,
527 fSelV0MaxDcaDaughters,fSelV0MinDecayLength);
529 fHistPtot = dynamic_cast<TH1F*> (GetOutputData(0));
531 Printf("ERROR: fHistPtot not available");
535 TCanvas *cTaskDedx = new TCanvas("AliAnalysisTaskESDDedx","Ptot",10,10,510,510);
536 cTaskDedx->cd(1)->SetLogy();
537 fHistPtot->DrawCopy("E");
539 //____________________________________________________________________//
540 Bool_t AliAnalysisTaskESDDedx::IsAccepted(AliESDtrack* track) {
541 // Checks if the track is excluded from the cuts
544 Int_t nClustersITS = track->GetITSclusters(fIdxInt);
545 Int_t nClustersTPC = track->GetTPCclusters(fIdxInt);
547 Float_t chi2PerClusterITS = -1;
548 Float_t chi2PerClusterTPC = -1;
550 chi2PerClusterTPC = track->GetTPCchi2()/Float_t(nClustersTPC);
553 track->GetExternalCovariance(extCov);
555 if (fSelTrackRemoveKink)
556 if (track->GetKinkIndex(0)>0) return kFALSE;
557 if (fSelTrackMinClustersTPC)
558 if (nClustersTPC < fSelTrackMinClustersTPC) return kFALSE;
559 if (fSelTrackMinClustersITS)
560 if (nClustersITS < fSelTrackMinClustersITS) return kFALSE;
561 if (fSelTrackMaxChi2PerClusterTPC)
562 if (chi2PerClusterTPC > fSelTrackMaxChi2PerClusterTPC) return kFALSE;
563 if (fSelTrackMaxChi2PerClusterITS)
564 if (chi2PerClusterITS > fSelTrackMaxChi2PerClusterITS) return kFALSE;
566 if (fSelTrackMaxCov11)
567 if (extCov[0] > fSelTrackMaxCov11) return kFALSE;
568 if (fSelTrackMaxCov22)
569 if (extCov[2] > fSelTrackMaxCov22) return kFALSE;
570 if (fSelTrackMaxCov33)
571 if (extCov[5] > fSelTrackMaxCov33) return kFALSE;
572 if (fSelTrackMaxCov44)
573 if (extCov[9] > fSelTrackMaxCov44) return kFALSE;
574 if(fSelTrackMaxCov55)
575 if (extCov[14] > fSelTrackMaxCov55) return kFALSE;
577 if(fMaxSigmaToVertexFlag)
578 if(GetSigmaToVertex(track) > fMaxSigmaToVertex) return kFALSE;
580 if ((track->GetStatus() & AliESDtrack::kITSrefit) == 0) return kFALSE;
582 if ((track->GetStatus() & AliESDtrack::kTPCrefit) == 0) return kFALSE;
583 if((Pt < fMinPt) || (Pt > fMaxPt)) return kFALSE;