]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/AntiprotonToProton/AliProtonAnalysisBase.cxx
changes for Vertex and Tracks classes
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / AntiprotonToProton / AliProtonAnalysisBase.cxx
CommitLineData
0ab648ea 1/**************************************************************************
2 * Author: Panos Christakoglou. *
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/* $Id: AliProtonAnalysisBase.cxx 31056 2009-02-16 14:31:41Z pchrist $ */
15
16//-----------------------------------------------------------------
17// AliProtonAnalysisBase class
18// This is the class to deal with the proton analysis
19// Origin: Panos Christakoglou | Panos.Christakoglou@cern.ch
20//-----------------------------------------------------------------
21#include <Riostream.h>
0ab648ea 22#include <TCanvas.h>
23#include <TLatex.h>
73aba974 24#include <TF1.h>
735cc63d 25#include <TList.h>
26#include <TH1F.h>
0ab648ea 27
28#include <AliExternalTrackParam.h>
29#include <AliESDEvent.h>
0ab648ea 30#include <AliPID.h>
0ab648ea 31#include <AliVertexerTracks.h>
10d100d4 32#include <AliESDpid.h>
f81f1b19 33#include <AliTPCPIDResponse.h>
73aba974 34class AliLog;
35class AliESDVertex;
36
37#include "AliProtonAnalysisBase.h"
38
0ab648ea 39ClassImp(AliProtonAnalysisBase)
40
41//____________________________________________________________________//
42AliProtonAnalysisBase::AliProtonAnalysisBase() :
790140ac 43 TObject(), fProtonAnalysisLevel("ESD"), fAnalysisMC(kFALSE),
e56f08ed 44 fTriggerMode(kMB2), kUseOnlineTrigger(kFALSE), kUseOfflineTrigger(kFALSE),
45 fPhysicsSelection(0),
76161fe0 46 fProtonAnalysisMode(kTPC), fProtonPIDMode(kBayesian),
0ab648ea 47 fAnalysisEtaMode(kFALSE),
2f6dae44 48 fRunQAAnalysis(kFALSE),
e56f08ed 49 fVxMax(100.), fVyMax(100.), fVzMax(100.), fMinNumOfContributors(0),
1b4f5688 50 fNBinsX(0), fMinX(0.), fMaxX(0.),
51 fNBinsY(0), fMinY(0.), fMaxY(0.),
0ab648ea 52 fMinTPCClusters(0), fMinITSClusters(0),
1b4f5688 53 fMaxChi2PerTPCCluster(0.), fMaxChi2PerITSCluster(0.),
54 fMaxCov11(0.), fMaxCov22(0.), fMaxCov33(0.), fMaxCov44(0.), fMaxCov55(0.),
55 fMaxSigmaToVertex(0.), fMaxSigmaToVertexTPC(0.),
56 fMaxDCAXY(0.), fMaxDCAXYTPC(0.),
57 fMaxDCAZ(0.), fMaxDCAZTPC(0.),
58 fMaxDCA3D(0.), fMaxDCA3DTPC(0.),
8e317d0d 59 fMaxConstrainChi2(0.), fMinTPCdEdxPoints(0),
0ab648ea 60 fMinTPCClustersFlag(kFALSE), fMinITSClustersFlag(kFALSE),
61 fMaxChi2PerTPCClusterFlag(kFALSE), fMaxChi2PerITSClusterFlag(kFALSE),
62 fMaxCov11Flag(kFALSE), fMaxCov22Flag(kFALSE),
63 fMaxCov33Flag(kFALSE), fMaxCov44Flag(kFALSE), fMaxCov55Flag(kFALSE),
64 fMaxSigmaToVertexFlag(kFALSE), fMaxSigmaToVertexTPCFlag(kFALSE),
65 fMaxDCAXYFlag(kFALSE), fMaxDCAXYTPCFlag(kFALSE),
66 fMaxDCAZFlag(kFALSE), fMaxDCAZTPCFlag(kFALSE),
67 fMaxDCA3DFlag(kFALSE), fMaxDCA3DTPCFlag(kFALSE),
68 fMaxConstrainChi2Flag(kFALSE),
69 fITSRefitFlag(kFALSE), fTPCRefitFlag(kFALSE),
f62e9410 70 fESDpidFlag(kFALSE), fTPCpidFlag(kFALSE), fTOFpidFlag(kFALSE),
4787e0ca 71 fPointOnSPDLayersFlag(0), fPointOnSDDLayersFlag(0), fPointOnSSDLayersFlag(0),
0ab648ea 72 fPointOnITSLayer1Flag(0), fPointOnITSLayer2Flag(0),
73 fPointOnITSLayer3Flag(0), fPointOnITSLayer4Flag(0),
74 fPointOnITSLayer5Flag(0), fPointOnITSLayer6Flag(0),
87a55728 75 fMinTPCdEdxPointsFlag(kFALSE),
8e317d0d 76 fPtDependentDcaXY(0), fPtDependentDcaXYFlag(kFALSE), fNSigmaDCAXY(0),
0ab648ea 77 fFunctionProbabilityFlag(kFALSE),
8e317d0d 78 fNSigma(0), fNRatio(0.),
0ab648ea 79 fElectronFunction(0), fMuonFunction(0),
80 fPionFunction(0), fKaonFunction(0), fProtonFunction(0),
bc554c73 81 fDebugMode(kFALSE), fListVertexQA(0) {
0ab648ea 82 //Default constructor
83 for(Int_t i = 0; i < 5; i++) fPartFrac[i] = 0.0;
c6909683 84 /*for(Int_t i = 0; i < 24; i++) {
87a55728 85 fdEdxMean[i] = 0.0;
86 fdEdxSigma[i] = 0.0;
c6909683 87 }*/
bc554c73 88
89 fListVertexQA = new TList();
735cc63d 90 fListVertexQA->SetName("fListVertexQA");
91 TH1F *gHistVx = new TH1F("gHistVx",
92 "Vx distribution;V_{x} [cm];Entries",
93 100,-5.,5.);
94 gHistVx->SetFillColor(kRed-2);
95 fListVertexQA->Add(gHistVx);
96 TH1F *gHistVxAccepted = new TH1F("gHistVxaccepted",
97 "Vx distribution;V_{x} [cm];Entries",
98 100,-5.,5.);
99 fListVertexQA->Add(gHistVxAccepted);
100 TH1F *gHistVy = new TH1F("gHistVy",
101 "Vy distribution;V_{y} [cm];Entries",
102 100,-5.,5.);
103 gHistVy->SetFillColor(kRed-2);
104 fListVertexQA->Add(gHistVy);
105 TH1F *gHistVyAccepted = new TH1F("gHistVyaccepted",
106 "Vy distribution;V_{y} [cm];Entries",
107 100,-5.,5.);
108 fListVertexQA->Add(gHistVyAccepted);
109 TH1F *gHistVz = new TH1F("gHistVz",
110 "Vz distribution;V_{z} [cm];Entries",
111 100,-25.,25.);
112 gHistVz->SetFillColor(kRed-2);
113 fListVertexQA->Add(gHistVz);
114 TH1F *gHistVzAccepted = new TH1F("gHistVzaccepted",
115 "Vz distribution;V_{z} [cm];Entries",
116 100,-25.,25.);
117 fListVertexQA->Add(gHistVzAccepted);
118
e56f08ed 119 TH1F *gHistNumberOfContributors = new TH1F("gHistNumberOfContributors",
120 "Number of contributors;N_{contr.};Entries",
121 100,0.,100.);
122 fListVertexQA->Add(gHistNumberOfContributors);
123
735cc63d 124
0ab648ea 125}
126
127//____________________________________________________________________//
128AliProtonAnalysisBase::~AliProtonAnalysisBase() {
129 //Default destructor
130 if(fElectronFunction) delete fElectronFunction;
131 if(fMuonFunction) delete fMuonFunction;
132 if(fPionFunction) delete fPionFunction;
133 if(fKaonFunction) delete fKaonFunction;
134 if(fProtonFunction) delete fProtonFunction;
735cc63d 135 if(fListVertexQA) delete fListVertexQA;
bc554c73 136 if(fPtDependentDcaXY) delete fPtDependentDcaXY;
137 if(fPhysicsSelection) delete fPhysicsSelection;
0ab648ea 138}
139
140//____________________________________________________________________//
141Double_t AliProtonAnalysisBase::GetParticleFraction(Int_t i, Double_t p) {
142 //Return the a priori probs
143 Double_t partFrac=0;
144 if(fFunctionProbabilityFlag) {
145 if(i == 0) partFrac = fElectronFunction->Eval(p);
146 if(i == 1) partFrac = fMuonFunction->Eval(p);
147 if(i == 2) partFrac = fPionFunction->Eval(p);
148 if(i == 3) partFrac = fKaonFunction->Eval(p);
149 if(i == 4) partFrac = fProtonFunction->Eval(p);
150 }
151 else partFrac = fPartFrac[i];
152
153 return partFrac;
154}
155
156//____________________________________________________________________//
157Bool_t AliProtonAnalysisBase::IsInPhaseSpace(AliESDtrack* const track) {
158 // Checks if the track is outside the analyzed y-Pt phase space
df201289 159 Double_t gP = 0.0, gPt = 0.0, gPx = 0.0, gPy = 0.0, gPz = 0.0;
0ab648ea 160 Double_t eta = 0.0;
161
afc2ac17 162 if((fProtonAnalysisMode == kTPC) || (fProtonAnalysisMode == kHybrid) || (fProtonAnalysisMode == kFullHybrid)) {
0ab648ea 163 AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
164 if(!tpcTrack) {
df201289 165 gP = 0.0; gPt = 0.0; gPx = 0.0; gPy = 0.0; gPz = 0.0; eta = -10.0;
0ab648ea 166 }
167 else {
df201289 168 gP = tpcTrack->P();
0ab648ea 169 gPt = tpcTrack->Pt();
170 gPx = tpcTrack->Px();
171 gPy = tpcTrack->Py();
172 gPz = tpcTrack->Pz();
173 eta = tpcTrack->Eta();
174 }
a84984fa 175 }//standalone TPC or Hybrid TPC approaches
0ab648ea 176 else {
df201289 177 gP = track->P();
0ab648ea 178 gPt = track->Pt();
179 gPx = track->Px();
180 gPy = track->Py();
181 gPz = track->Pz();
182 eta = track->Eta();
183 }
184
e56f08ed 185 if((gPt < fMinY) || (gPt > fMaxY)) {
0ab648ea 186 if(fDebugMode)
187 Printf("IsInPhaseSpace: Track rejected because it has a Pt value of %lf (accepted interval: %lf - %lf)",gPt,fMinY,fMaxY);
188 return kFALSE;
189 }
e56f08ed 190 if((gP < fMinY) || (gP > fMaxY)) {
191 if(fDebugMode)
192 Printf("IsInPhaseSpace: Track rejected because it has a P value of %lf (accepted interval: %lf - %lf)",gP,fMinY,fMaxY);
193 return kFALSE;
194 }
0ab648ea 195 if(fAnalysisEtaMode) {
196 if((eta < fMinX) || (eta > fMaxX)) {
197 if(fDebugMode)
198 Printf("IsInPhaseSpace: Track rejected because it has an eta value of %lf (accepted interval: %lf - %lf)",eta,fMinX,fMaxX);
199 return kFALSE;
200 }
201 }
202 else {
203 if((Rapidity(gPx,gPy,gPz) < fMinX) || (Rapidity(gPx,gPy,gPz) > fMaxX)) {
204 if(fDebugMode)
205 Printf("IsInPhaseSpace: Track rejected because it has a y value of %lf (accepted interval: %lf - %lf)",Rapidity(gPx,gPy,gPz),fMinX,fMaxX);
206 return kFALSE;
207 }
208 }
209
210 return kTRUE;
211}
212
213//____________________________________________________________________//
30f87a5b 214Bool_t AliProtonAnalysisBase::IsAccepted(AliESDtrack* track) {
0ab648ea 215 // Checks if the track is excluded from the cuts
e0919e87 216 //Int_t fIdxInt[200];
217 //Int_t nClustersITS = track->GetITSclusters(fIdxInt);
218 //Int_t nClustersTPC = track->GetTPCclusters(fIdxInt);
219 Int_t nClustersITS = track->GetITSclusters(0x0);
220 Int_t nClustersTPC = track->GetTPCclusters(0x0);
0ab648ea 221
222 Float_t chi2PerClusterITS = -1;
223 if (nClustersITS!=0)
224 chi2PerClusterITS = track->GetITSchi2()/Float_t(nClustersITS);
225 Float_t chi2PerClusterTPC = -1;
226 if (nClustersTPC!=0)
227 chi2PerClusterTPC = track->GetTPCchi2()/Float_t(nClustersTPC);
228
229 Double_t extCov[15];
230 track->GetExternalCovariance(extCov);
231
4787e0ca 232 if(fPointOnSPDLayersFlag) {
233 if((!track->HasPointOnITSLayer(0))&&(!track->HasPointOnITSLayer(1))) {
234 if(fDebugMode)
235 Printf("IsAccepted: Track rejected because it doesn't have a point on either SPD layers");
236 return kFALSE;
237 }
238 }
239 if(fPointOnSDDLayersFlag) {
240 if((!track->HasPointOnITSLayer(2))&&(!track->HasPointOnITSLayer(3))) {
241 if(fDebugMode)
242 Printf("IsAccepted: Track rejected because it doesn't have a point on either SDD layers");
243 return kFALSE;
244 }
245 }
246 if(fPointOnSSDLayersFlag) {
247 if((!track->HasPointOnITSLayer(4))&&(!track->HasPointOnITSLayer(5))) {
248 if(fDebugMode)
249 Printf("IsAccepted: Track rejected because it doesn't have a point on either SSD layers");
250 return kFALSE;
251 }
252 }
0ab648ea 253 if(fPointOnITSLayer1Flag) {
254 if(!track->HasPointOnITSLayer(0)) {
255 if(fDebugMode)
256 Printf("IsAccepted: Track rejected because it doesn't have a point on the 1st ITS layer");
257 return kFALSE;
258 }
259 }
260 if(fPointOnITSLayer2Flag) {
261 if(!track->HasPointOnITSLayer(1)) {
262 if(fDebugMode)
263 Printf("IsAccepted: Track rejected because it doesn't have a point on the 2nd ITS layer");
264 return kFALSE;
265 }
266 }
267 if(fPointOnITSLayer3Flag) {
268 if(!track->HasPointOnITSLayer(2)) {
269 if(fDebugMode)
270 Printf("IsAccepted: Track rejected because it doesn't have a point on the 3rd ITS layer");
271 return kFALSE;
272 }
273 }
274 if(fPointOnITSLayer4Flag) {
275 if(!track->HasPointOnITSLayer(3)) {
276 if(fDebugMode)
277 Printf("IsAccepted: Track rejected because it doesn't have a point on the 4th ITS layer");
278 return kFALSE;
279 }
280 }
281 if(fPointOnITSLayer5Flag) {
282 if(!track->HasPointOnITSLayer(4)) {
283 if(fDebugMode)
284 Printf("IsAccepted: Track rejected because it doesn't have a point on the 5th ITS layer");
285 return kFALSE;
286 }
287 }
288 if(fPointOnITSLayer6Flag) {
289 if(!track->HasPointOnITSLayer(5)) {
290 if(fDebugMode)
291 Printf("IsAccepted: Track rejected because it doesn't have a point on the 6th ITS layer");
292 return kFALSE;
293 }
294 }
295 if(fMinITSClustersFlag) {
296 if(nClustersITS < fMinITSClusters) {
297 if(fDebugMode)
298 Printf("IsAccepted: Track rejected because it has %d ITS points (min. requested: %d)",nClustersITS,fMinITSClusters);
299 return kFALSE;
300 }
301 }
302 if(fMaxChi2PerITSClusterFlag) {
303 if(chi2PerClusterITS > fMaxChi2PerITSCluster) {
304 if(fDebugMode)
305 Printf("IsAccepted: Track rejected because it has a chi2 per ITS cluster %lf (max. requested: %lf)",chi2PerClusterITS,fMaxChi2PerITSCluster);
306 return kFALSE;
307 }
308 }
309 if(fMinTPCClustersFlag) {
310 if(nClustersTPC < fMinTPCClusters) {
311 if(fDebugMode)
312 Printf("IsAccepted: Track rejected because it has %d TPC clusters (min. requested: %d)",nClustersTPC,fMinTPCClusters);
313 return kFALSE;
314 }
315 }
316 if(fMaxChi2PerTPCClusterFlag) {
317 if(chi2PerClusterTPC > fMaxChi2PerTPCCluster) {
318 if(fDebugMode)
319 Printf("IsAccepted: Track rejected because it has a chi2 per TPC cluster %lf (max. requested: %lf)",chi2PerClusterTPC,fMaxChi2PerTPCCluster);
320 return kFALSE;
321 }
322 }
323 if(fMaxCov11Flag) {
324 if(extCov[0] > fMaxCov11) {
325 if(fDebugMode)
326 Printf("IsAccepted: Track rejected because it has a cov11 value of %lf (max. requested: %lf)",extCov[0],fMaxCov11);
327 return kFALSE;
328 }
329 }
330 if(fMaxCov22Flag) {
331 if(extCov[2] > fMaxCov22) {
332 if(fDebugMode)
333 Printf("IsAccepted: Track rejected because it has a cov22 value of %lf (max. requested: %lf)",extCov[2],fMaxCov22);
334 return kFALSE;
335 }
336 }
337 if(fMaxCov33Flag) {
338 if(extCov[5] > fMaxCov33) {
339 if(fDebugMode)
340 Printf("IsAccepted: Track rejected because it has a cov33 value of %lf (max. requested: %lf)",extCov[5],fMaxCov33);
341 return kFALSE;
342 }
343 }
344 if(fMaxCov44Flag) {
345 if(extCov[9] > fMaxCov44) {
346 if(fDebugMode)
347 Printf("IsAccepted: Track rejected because it has a cov44 value of %lf (max. requested: %lf)",extCov[9],fMaxCov44);
348 return kFALSE;
349 }
350 }
351 if(fMaxCov55Flag) {
352 if(extCov[14] > fMaxCov55) {
353 if(fDebugMode)
354 Printf("IsAccepted: Track rejected because it has a cov55 value of %lf (max. requested: %lf)",extCov[14],fMaxCov55);
355 return kFALSE;
356 }
357 }
87a55728 358 if(fMinTPCdEdxPointsFlag) {
359 if(track->GetTPCsignalN() < fMinTPCdEdxPoints) {
360 if(fDebugMode)
361 Printf("IsAccepted: Track rejected because it has %d TPC points for the calculation of the energy loss (min. requested: %d)",track->GetTPCsignalN(),fMinTPCdEdxPoints);
362 return kFALSE;
363 }
364 }
0ab648ea 365 if(fITSRefitFlag) {
366 if ((track->GetStatus() & AliESDtrack::kITSrefit) == 0) {
367 if(fDebugMode)
368 Printf("IsAccepted: Track rejected because it has no ITS refit flag");
369 return kFALSE;
370 }
371 }
372 if(fTPCRefitFlag) {
373 if ((track->GetStatus() & AliESDtrack::kTPCrefit) == 0) {
374 if(fDebugMode)
375 Printf("IsAccepted: Track rejected because it has no TPC refit flag");
376 return kFALSE;
377 }
378 }
379 if(fESDpidFlag) {
380 if ((track->GetStatus() & AliESDtrack::kESDpid) == 0) {
381 if(fDebugMode)
382 Printf("IsAccepted: Track rejected because it has no ESD pid flag");
383 return kFALSE;
384 }
385 }
386 if(fTPCpidFlag) {
387 if ((track->GetStatus() & AliESDtrack::kTPCpid) == 0) {
388 if(fDebugMode)
389 Printf("IsAccepted: Track rejected because it has no TPC pid flag");
390 return kFALSE;
391 }
392 }
f62e9410 393 if(fTOFpidFlag) {
394 if ((track->GetStatus() & AliESDtrack::kTOFpid) == 0) {
395 if(fDebugMode)
396 Printf("IsAccepted: Track rejected because it has no TOF pid flag");
397 return kFALSE;
398 }
399 }
0ab648ea 400
401 return kTRUE;
402}
403
ff1c9f70 404//____________________________________________________________________//
405void AliProtonAnalysisBase::SetPtDependentDCAxy(Int_t nSigma,
406 Double_t p0,
407 Double_t p1,
408 Double_t p2) {
409 //Pt dependent dca cut in xy
410 fNSigmaDCAXY = nSigma;
411 fPtDependentDcaXY = new TF1("fPtDependentDcaXY","[0]+[1]/x^[2]",0.1,10.1);
412 fPtDependentDcaXY->SetParameter(0,p0);
413 fPtDependentDcaXY->SetParameter(1,p1);
414 fPtDependentDcaXY->SetParameter(2,p2);
415 fPtDependentDcaXYFlag = kTRUE;
416}
417
e56f08ed 418//____________________________________________________________________//
419Bool_t AliProtonAnalysisBase::IsPrimary(AliESDEvent *esd,
420 const AliESDVertex *vertex,
421 AliESDtrack* track) {
422 // Checks if the track is a primary-like candidate
842c574e 423 const Double_t kMicrometer2Centimeter = 0.0001;
e56f08ed 424 Double_t gPt = 0.0, gPx = 0.0, gPy = 0.0, gPz = 0.0;
425 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0}; //The impact parameters and their covariance.
426 Double_t dca3D = 0.0;
afc2ac17 427 Float_t dcaXY = 0.0, dcaZ = 0.0;
428
e56f08ed 429 if((fProtonAnalysisMode == kTPC)||(fProtonAnalysisMode == kHybrid)) {
430 AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
431 if(!tpcTrack) {
432 gPt = 0.0; gPx = 0.0; gPy = 0.0; gPz = 0.0;
433 dca[0] = -100.; dca[1] = -100.; dca3D = -100.;
434 cov[0] = -100.; cov[1] = -100.; cov[2] = -100.;
435 }
436 else {
437 gPt = tpcTrack->Pt();
438 gPx = tpcTrack->Px();
439 gPy = tpcTrack->Py();
440 gPz = tpcTrack->Pz();
441 tpcTrack->PropagateToDCA(vertex,
442 esd->GetMagneticField(),
443 100.,dca,cov);
444 }
445 }//standalone TPC or hybrid TPC approaches
afc2ac17 446 if(fProtonAnalysisMode == kFullHybrid) {
447 AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
448 AliExternalTrackParam cParam;
449 if(!tpcTrack) {
450 gPt = 0.0; gPx = 0.0; gPy = 0.0; gPz = 0.0;
451 dca[0] = -100.; dca[1] = -100.; dca3D = -100.;
452 cov[0] = -100.; cov[1] = -100.; cov[2] = -100.;
453 }
454 else {
455 gPt = tpcTrack->Pt();
456 gPx = tpcTrack->Px();
457 gPy = tpcTrack->Py();
458 gPz = tpcTrack->Pz();
459 track->RelateToVertex(vertex,
460 esd->GetMagneticField(),
461 100.,&cParam);
462 track->GetImpactParameters(dcaXY,dcaZ);
463 dca[0] = dcaXY; dca[1] = dcaZ;
464 }
465 }//standalone TPC or hybrid TPC approaches
e56f08ed 466 else {
467 gPt = track->Pt();
468 gPx = track->Px();
469 gPy = track->Py();
470 gPz = track->Pz();
471 track->PropagateToDCA(vertex,
472 esd->GetMagneticField(),
473 100.,dca,cov);
474 }
475 dca3D = TMath::Sqrt(TMath::Power(dca[0],2) +
476 TMath::Power(dca[1],2));
477
e56f08ed 478 if(fMaxSigmaToVertexFlag) {
479 if(GetSigmaToVertex(track) > fMaxSigmaToVertex) {
480 if(fDebugMode)
481 Printf("IsPrimary: Track rejected because it has a %lf sigmas to vertex (max. requested: %lf)",GetSigmaToVertex(track),fMaxSigmaToVertex);
482 return kFALSE;
483 }
484 }
485 if(fMaxSigmaToVertexTPCFlag) {
486 if(GetSigmaToVertex(track) > fMaxSigmaToVertexTPC) {
487 if(fDebugMode)
488 Printf("IsPrimary: Track rejected because it has a %lf sigmas to vertex TPC (max. requested: %lf)",GetSigmaToVertex(track),fMaxSigmaToVertexTPC);
489 return kFALSE;
490 }
491 }
492 if(fMaxDCAXYFlag) {
493 if(TMath::Abs(dca[0]) > fMaxDCAXY) {
494 if(fDebugMode)
495 Printf("IsPrimary: Track rejected because it has a value of dca(xy) of %lf (max. requested: %lf)",TMath::Abs(dca[0]),fMaxDCAXY);
496 return kFALSE;
497 }
498 }
499 if(fMaxDCAXYTPCFlag) {
500 if(TMath::Abs(dca[0]) > fMaxDCAXYTPC) {
501 if(fDebugMode)
502 Printf("IsPrimary: Track rejected because it has a value of dca(xy) (TPC) of %lf (max. requested: %lf)",TMath::Abs(dca[0]),fMaxDCAXYTPC);
503 return kFALSE;
504 }
505 }
506 if(fMaxDCAZFlag) {
507 if(TMath::Abs(dca[1]) > fMaxDCAZ) {
508 if(fDebugMode)
509 Printf("IsPrimary: Track rejected because it has a value of dca(z) of %lf (max. requested: %lf)",TMath::Abs(dca[1]),fMaxDCAZ);
510 return kFALSE;
511 }
512 }
513 if(fMaxDCAZTPCFlag) {
514 if(TMath::Abs(dca[1]) > fMaxDCAZTPC) {
515 if(fDebugMode)
516 Printf("IsPrimary: Track rejected because it has a value of dca(z) (TPC) of %lf (max. requested: %lf)",TMath::Abs(dca[1]),fMaxDCAZTPC);
517 return kFALSE;
518 }
519 }
520 if(fMaxDCA3DFlag) {
521 if(TMath::Abs(dca3D) > fMaxDCA3D) {
522 if(fDebugMode)
523 Printf("IsPrimary: Track rejected because it has a value of dca(3D) of %lf (max. requested: %lf)",TMath::Abs(dca3D),fMaxDCA3D);
524 return kFALSE;
525 }
526 }
527 if(fMaxDCA3DTPCFlag) {
528 if(TMath::Abs(dca3D) > fMaxDCA3DTPC) {
529 if(fDebugMode)
530 Printf("IsPrimary: Track rejected because it has a value of dca(3D) (TPC) of %lf (max. requested: %lf)",TMath::Abs(dca3D),fMaxDCA3DTPC);
531 return kFALSE;
532 }
533 }
534 if(fMaxConstrainChi2Flag) {
535 if(track->GetConstrainedChi2() > 0)
536 if(TMath::Log(track->GetConstrainedChi2()) > fMaxConstrainChi2) {
537 if(fDebugMode)
538 Printf("IsPrimary: Track rejected because it has a value of the constrained chi2 to the vertex of %lf (max. requested: %lf)",TMath::Log(track->GetConstrainedChi2()),fMaxConstrainChi2);
539 return kFALSE;
540 }
541 }
ff1c9f70 542 if(fPtDependentDcaXYFlag) {
842c574e 543 if(TMath::Abs(dca[0]) > kMicrometer2Centimeter*fNSigmaDCAXY*fPtDependentDcaXY->Eval(gPt)) {
ff1c9f70 544 if(fDebugMode)
9e5b51c2 545 Printf("IsPrimary: Track rejected because it has a value of the dca(xy) higher than the %d sigma pt dependent cut: %lf (max. requested: %lf)",fNSigmaDCAXY,TMath::Abs(dca[0]),fNSigmaDCAXY*fPtDependentDcaXY->Eval(gPt));
ff1c9f70 546 return kFALSE;
547 }
548 }
e56f08ed 549
550 return kTRUE;
551}
552
0ab648ea 553//____________________________________________________________________//
554Float_t AliProtonAnalysisBase::GetSigmaToVertex(AliESDtrack* esdTrack) const {
555 // Calculates the number of sigma to the vertex.
0ab648ea 556 Float_t b[2];
557 Float_t bRes[2];
558 Float_t bCov[3];
afc2ac17 559 if((fProtonAnalysisMode == kTPC)&&(fProtonAnalysisMode != kHybrid)&&(fProtonAnalysisMode != kFullHybrid))
0ab648ea 560 esdTrack->GetImpactParametersTPC(b,bCov);
561 else
562 esdTrack->GetImpactParameters(b,bCov);
563
564 if (bCov[0]<=0 || bCov[2]<=0) {
565 //AliDebug(1, "Estimated b resolution lower or equal zero!");
566 bCov[0]=0; bCov[2]=0;
567 }
568 bRes[0] = TMath::Sqrt(bCov[0]);
569 bRes[1] = TMath::Sqrt(bCov[2]);
570
571 if (bRes[0] == 0 || bRes[1] ==0) return -1;
572
573 Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
574
575 if (TMath::Exp(-d * d / 2) < 1e-10) return 1000;
576
577 d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
578
579 return d;
580}
581
582//____________________________________________________________________//
583Double_t AliProtonAnalysisBase::Rapidity(Double_t gPx,
584 Double_t gPy,
585 Double_t gPz) const {
586 //returns the rapidity of the proton - to be removed
587 Double_t fMass = 9.38270000000000048e-01;
588
589 Double_t gP = TMath::Sqrt(TMath::Power(gPx,2) +
590 TMath::Power(gPy,2) +
591 TMath::Power(gPz,2));
592 Double_t energy = TMath::Sqrt(gP*gP + fMass*fMass);
593 Double_t y = -999;
594 if(energy != gPz)
595 y = 0.5*TMath::Log((energy + gPz)/(energy - gPz));
596
597 return y;
598}
599
600//________________________________________________________________________
601const AliESDVertex* AliProtonAnalysisBase::GetVertex(AliESDEvent* esd,
602 AnalysisMode mode,
603 Double_t gVxMax,
604 Double_t gVyMax,
605 Double_t gVzMax) {
606 // Get the vertex from the ESD and returns it if the vertex is valid
607 // Second argument decides which vertex is used (this selects
608 // also the quality criteria that are applied)
609 const AliESDVertex* vertex = 0;
afc2ac17 610 if((mode == kHybrid)||(mode == kFullHybrid))
0ab648ea 611 vertex = esd->GetPrimaryVertexSPD();
612 else if(mode == kTPC){
613 Double_t kBz = esd->GetMagneticField();
614 AliVertexerTracks vertexer(kBz);
615 vertexer.SetTPCMode();
616 AliESDVertex *vTPC = vertexer.FindPrimaryVertex(esd);
617 esd->SetPrimaryVertexTPC(vTPC);
618 for (Int_t i=0; i<esd->GetNumberOfTracks(); i++) {
619 AliESDtrack *t = esd->GetTrack(i);
620 t->RelateToVertexTPC(vTPC, kBz, kVeryBig);
621 }
622 delete vTPC;
623 vertex = esd->GetPrimaryVertexTPC();
624 }
625 else if(mode == kGlobal)
626 vertex = esd->GetPrimaryVertex();
627 else
628 Printf("GetVertex: ERROR: Invalid second argument %d", mode);
629
630 if(!vertex) {
631 if(fDebugMode)
632 Printf("GetVertex: Event rejected because there is no valid vertex object");
633 return 0;
634 }
635
636 // check Ncontributors
637 if(vertex->GetNContributors() <= 0) {
638 if(fDebugMode)
639 Printf("GetVertex: Event rejected because the number of contributors for the vertex determination is <= 0");
640 return 0;
641 }
642
643 // check resolution
644 Double_t zRes = vertex->GetZRes();
645 if(zRes == 0) {
646 if(fDebugMode)
647 Printf("GetVertex: Event rejected because the value of the vertex resolution in z is 0");
648 return 0;
649 }
e690d4d0 650 ((TH1F *)(fListVertexQA->At(0)))->Fill(vertex->GetX());
651 ((TH1F *)(fListVertexQA->At(2)))->Fill(vertex->GetY());
652 ((TH1F *)(fListVertexQA->At(4)))->Fill(vertex->GetZ());
e56f08ed 653
0ab648ea 654 //check position
e690d4d0 655 if(TMath::Abs(vertex->GetX()) > gVxMax) {
0ab648ea 656 if(fDebugMode)
e690d4d0 657 Printf("GetVertex: Event rejected because it has a Vx value of %lf cm (accepted interval: -%lf - %lf)",TMath::Abs(vertex->GetX()),gVxMax,gVxMax);
0ab648ea 658 return 0;
659 }
e690d4d0 660 if(TMath::Abs(vertex->GetY()) > gVyMax) {
0ab648ea 661 if(fDebugMode)
e690d4d0 662 Printf("GetVertex: Event rejected because it has a Vy value of %lf cm (accepted interval: -%lf - %lf)",TMath::Abs(vertex->GetY()),gVyMax,gVyMax);
0ab648ea 663 return 0;
664 }
e690d4d0 665 if(TMath::Abs(vertex->GetZ()) > gVzMax) {
0ab648ea 666 if(fDebugMode)
e690d4d0 667 Printf("GetVertex: Event rejected because it has a Vz value of %lf cm (accepted interval: -%lf - %lf)",TMath::Abs(vertex->GetZ()),gVzMax,gVzMax);
0ab648ea 668 return 0;
669 }
e690d4d0 670 ((TH1F *)(fListVertexQA->At(1)))->Fill(vertex->GetX());
671 ((TH1F *)(fListVertexQA->At(3)))->Fill(vertex->GetY());
672 ((TH1F *)(fListVertexQA->At(5)))->Fill(vertex->GetZ());
e56f08ed 673 ((TH1F *)(fListVertexQA->At(6)))->Fill(vertex->GetNContributors());
674
675 //check number of contributors
676 if(fMinNumOfContributors > 0) {
677 if(fMinNumOfContributors > vertex->GetNContributors()) {
678 if(fDebugMode)
679 Printf("GetVertex: Event rejected because it has %d number of contributors (requested minimum: %d)",vertex->GetNContributors(),fMinNumOfContributors);
680
681 return 0;
682 }
683 }
0ab648ea 684
685 return vertex;
686}
687
688//________________________________________________________________________
689Bool_t AliProtonAnalysisBase::IsEventTriggered(const AliESDEvent *esd,
690 TriggerMode trigger) {
691 // check if the event was triggered
692 ULong64_t triggerMask = esd->GetTriggerMask();
42270c4c 693 TString firedTriggerClass = esd->GetFiredTriggerClasses();
0ab648ea 694
42270c4c 695 if(fAnalysisMC) {
696 // definitions from p-p.cfg
697 ULong64_t spdFO = (1 << 14);
698 ULong64_t v0left = (1 << 11);
699 ULong64_t v0right = (1 << 12);
700
701 switch (trigger) {
702 case kMB1: {
703 if (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
704 return kTRUE;
705 break;
706 }
707 case kMB2: {
708 if (triggerMask & spdFO && ((triggerMask & v0left) || (triggerMask & v0right)))
709 return kTRUE;
710 break;
711 }
712 case kSPDFASTOR: {
713 if (triggerMask & spdFO)
714 return kTRUE;
715 break;
716 }
717 }//switch
0ab648ea 718 }
42270c4c 719 else {
4de4661f 720 if(kUseOnlineTrigger) {
721 if(firedTriggerClass.Contains("CINT1B-ABCE-NOPF-ALL"))
722 return kTRUE;
723 }
724 else if(!kUseOnlineTrigger)
0ab648ea 725 return kTRUE;
0ab648ea 726 }
0ab648ea 727
728 return kFALSE;
729}
730
731//________________________________________________________________________
732TCanvas *AliProtonAnalysisBase::GetListOfCuts() {
733 // return the list of cuts and their cut values for reference
734 TLatex l;
735 l.SetTextAlign(12);
736 l.SetTextSize(0.04);
737
738 TString listOfCuts;
739
740 TCanvas *c = new TCanvas("cListOfCuts","List of cuts",0,0,900,600);
741 c->SetFillColor(10); c->SetHighLightColor(41);
742 c->Divide(3,2);
743
744 c->cd(1)->SetFillColor(10);
745 l.DrawLatex(0.3,0.9,"Analysis details: List of cuts\n\n");
746
747 listOfCuts = "Analysis level: "; listOfCuts += fProtonAnalysisLevel;
748 l.DrawLatex(0.1,0.82,listOfCuts.Data());
749 listOfCuts = "Analysis mode: ";
750 if(fProtonAnalysisMode == kTPC) listOfCuts += "TPC standalone";
751 if(fProtonAnalysisMode == kHybrid) listOfCuts += "Hybrid TPC";
afc2ac17 752 if(fProtonAnalysisMode == kFullHybrid) listOfCuts += "Full Hybrid TPC";
0ab648ea 753 if(fProtonAnalysisMode == kGlobal) listOfCuts += "Global tracking";
754 l.DrawLatex(0.1,0.74,listOfCuts.Data());
755 listOfCuts = "Trigger mode: ";
756 if(fTriggerMode == kMB1) listOfCuts += "Minimum bias 1";
757 if(fTriggerMode == kMB2) listOfCuts += "Minimum bias 2";
758 if(fTriggerMode == kSPDFASTOR) listOfCuts += "FastOR (SPD)";
759 l.DrawLatex(0.1,0.66,listOfCuts.Data());
760 listOfCuts = "PID mode: ";
761 if(fProtonPIDMode == kBayesian) listOfCuts += "Bayesian PID";
c6909683 762 if(fProtonPIDMode == kRatio) {
763 listOfCuts += "Z = ln((dE/dx)_{exp.}/(dE/dx)_{theor.}) > ";
764 listOfCuts += fNRatio;
87a55728 765 }
c6909683 766 if(fProtonPIDMode == kSigma) {
767 listOfCuts += "N_{#sigma} area: "; listOfCuts += fNSigma;
87a55728 768 listOfCuts += " #sigma";
769 }
c6909683 770 //if(fProtonPIDMode == kSigma2) {
771 //listOfCuts += "N_{#sigma}(2) area: "; listOfCuts += fNSigma;
772 //listOfCuts += " #sigma";
773 //}
0ab648ea 774 l.DrawLatex(0.1,0.58,listOfCuts.Data());
775 listOfCuts = "Accepted vertex diamond: ";
c6909683 776 l.DrawLatex(0.1,0.52,listOfCuts.Data());
0ab648ea 777 listOfCuts = "|V_{x}| < "; listOfCuts += fVxMax; listOfCuts += " [cm]";
c6909683 778 l.DrawLatex(0.6,0.52,listOfCuts.Data());
0ab648ea 779 listOfCuts = "|V_{y}| < "; listOfCuts += fVyMax; listOfCuts += " [cm]";
c6909683 780 l.DrawLatex(0.6,0.45,listOfCuts.Data());
0ab648ea 781 listOfCuts = "|V_{z}| < "; listOfCuts += fVzMax; listOfCuts += " [cm]";
c6909683 782 l.DrawLatex(0.6,0.38,listOfCuts.Data());
783 listOfCuts = "N_{contributors} > "; listOfCuts += fMinNumOfContributors;
784 l.DrawLatex(0.6,0.31,listOfCuts.Data());
0ab648ea 785 listOfCuts = "Phase space: ";
786 l.DrawLatex(0.1,0.2,listOfCuts.Data());
787 if(fAnalysisEtaMode) listOfCuts = "|#eta| < ";
788 else listOfCuts = "|y| < ";
789 listOfCuts += TMath::Abs(fMinX);
790 l.DrawLatex(0.35,0.2,listOfCuts.Data());
791 listOfCuts = "N_{bins} = ";
792 listOfCuts += fNBinsX; listOfCuts += " (binning: ";
793 listOfCuts += (fMaxX-fMinX)/fNBinsX; listOfCuts += ")";
794 l.DrawLatex(0.35,0.15,listOfCuts.Data());
795 listOfCuts = ""; listOfCuts += fMinY; listOfCuts += " < P_{T} < ";
796 listOfCuts += fMaxY; listOfCuts += "GeV/c";
797 l.DrawLatex(0.35,0.1,listOfCuts.Data());
798 listOfCuts = "N_{bins} = ";
799 listOfCuts += fNBinsY; listOfCuts += " (binning: ";
800 listOfCuts += (fMaxY-fMinY)/fNBinsY; listOfCuts += ")";
801 l.DrawLatex(0.35,0.05,listOfCuts.Data());
802
803 c->cd(2)->SetFillColor(2);
4787e0ca 804 l.DrawLatex(0.3,0.95,"ITS related cuts");
805 listOfCuts = "Request a cluster on either of the SPD layers: ";
806 listOfCuts += fPointOnSPDLayersFlag;
807 l.DrawLatex(0.1,0.9,listOfCuts.Data());
808 listOfCuts = "Request a cluster on either of the SDD layers: ";
809 listOfCuts += fPointOnSDDLayersFlag;
810 l.DrawLatex(0.1,0.83,listOfCuts.Data());
811 listOfCuts = "Request a cluster on either of the SSD layers: ";
812 listOfCuts += fPointOnSSDLayersFlag;
813 l.DrawLatex(0.1,0.76,listOfCuts.Data());
814
0ab648ea 815 listOfCuts = "Request a cluster on SPD1: ";
816 listOfCuts += fPointOnITSLayer1Flag;
4787e0ca 817 l.DrawLatex(0.1,0.69,listOfCuts.Data());
0ab648ea 818 listOfCuts = "Request a cluster on SPD2: ";
819 listOfCuts += fPointOnITSLayer2Flag;
4787e0ca 820 l.DrawLatex(0.1,0.62,listOfCuts.Data());
0ab648ea 821 listOfCuts = "Request a cluster on SDD1: ";
822 listOfCuts += fPointOnITSLayer3Flag;
4787e0ca 823 l.DrawLatex(0.1,0.55,listOfCuts.Data());
0ab648ea 824 listOfCuts = "Request a cluster on SDD2: ";
825 listOfCuts += fPointOnITSLayer4Flag;
4787e0ca 826 l.DrawLatex(0.1,0.48,listOfCuts.Data());
0ab648ea 827 listOfCuts = "Request a cluster on SSD1: ";
828 listOfCuts += fPointOnITSLayer5Flag;
4787e0ca 829 l.DrawLatex(0.1,0.41,listOfCuts.Data());
0ab648ea 830 listOfCuts = "Request a cluster on SSD2: ";
831 listOfCuts += fPointOnITSLayer6Flag;
4787e0ca 832 l.DrawLatex(0.1,0.34,listOfCuts.Data());
0ab648ea 833 listOfCuts = "Minimum number of ITS clusters: ";
834 if(fMinITSClustersFlag) listOfCuts += fMinITSClusters;
835 else listOfCuts += "Not used";
4787e0ca 836 l.DrawLatex(0.1,0.27,listOfCuts.Data());
0ab648ea 837 listOfCuts = "Maximum #chi^{2} per ITS cluster: ";
838 if(fMaxChi2PerITSClusterFlag) listOfCuts += fMaxChi2PerITSCluster;
839 else listOfCuts += "Not used";
4787e0ca 840 l.DrawLatex(0.1,0.2,listOfCuts.Data());
0ab648ea 841
842 c->cd(3)->SetFillColor(3);
843 l.DrawLatex(0.3,0.9,"TPC related cuts");
844 listOfCuts = "Minimum number of TPC clusters: ";
845 if(fMinTPCClustersFlag) listOfCuts += fMinTPCClusters;
846 else listOfCuts += "Not used";
847 l.DrawLatex(0.1,0.7,listOfCuts.Data());
848 listOfCuts = "Maximum #chi^{2} per TPC cluster: ";
849 if(fMaxChi2PerTPCClusterFlag) listOfCuts += fMaxChi2PerTPCCluster;
850 else listOfCuts += "Not used";
851 l.DrawLatex(0.1,0.5,listOfCuts.Data());
57e749bb 852 listOfCuts = "Minimum number of TPC points for the dE/dx: ";
853 if(fMinTPCdEdxPointsFlag) listOfCuts += fMinTPCdEdxPoints;
854 else listOfCuts += "Not used";
855 l.DrawLatex(0.1,0.3,listOfCuts.Data());
0ab648ea 856
857 c->cd(4)->SetFillColor(4);
858 l.DrawLatex(0.3,0.9,"Tracking related cuts");
859 listOfCuts = "Maximum cov11: ";
860 if(fMaxCov11Flag) listOfCuts += fMaxCov11;
861 else listOfCuts += "Not used";
862 l.DrawLatex(0.1,0.75,listOfCuts.Data());
863 listOfCuts = "Maximum cov22: ";
864 if(fMaxCov22Flag) listOfCuts += fMaxCov22;
865 else listOfCuts += "Not used";
866 l.DrawLatex(0.1,0.6,listOfCuts.Data());
867 listOfCuts = "Maximum cov33: ";
868 if(fMaxCov33Flag) listOfCuts += fMaxCov33;
869 else listOfCuts += "Not used";
870 l.DrawLatex(0.1,0.45,listOfCuts.Data());
871 listOfCuts = "Maximum cov44: ";
872 if(fMaxCov44Flag) listOfCuts += fMaxCov44;
873 else listOfCuts += "Not used";
874 l.DrawLatex(0.1,0.3,listOfCuts.Data());
875 listOfCuts = "Maximum cov55: ";
876 if(fMaxCov55Flag) listOfCuts += fMaxCov55;
877 else listOfCuts += "Not used";
878 l.DrawLatex(0.1,0.15,listOfCuts.Data());
879
880 c->cd(5)->SetFillColor(5);
881 l.DrawLatex(0.3,0.9,"DCA related cuts");
882 listOfCuts = "Maximum sigma to vertex: ";
883 if(fMaxSigmaToVertexFlag) listOfCuts += fMaxSigmaToVertex;
884 else listOfCuts += "Not used";
885 l.DrawLatex(0.1,0.8,listOfCuts.Data());
886 listOfCuts = "Maximum sigma to vertex (TPC): ";
887 if(fMaxSigmaToVertexTPCFlag) listOfCuts += fMaxSigmaToVertexTPC;
888 else listOfCuts += "Not used";
889 l.DrawLatex(0.1,0.72,listOfCuts.Data());
890 listOfCuts = "Maximum DCA in xy: ";
891 if(fMaxDCAXYFlag) listOfCuts += fMaxDCAXY;
892 else listOfCuts += "Not used";
893 l.DrawLatex(0.1,0.64,listOfCuts.Data());
894 listOfCuts = "Maximum DCA in z: ";
895 if(fMaxDCAZFlag) listOfCuts += fMaxDCAZ;
896 else listOfCuts += "Not used";
897 l.DrawLatex(0.1,0.56,listOfCuts.Data());
898 listOfCuts = "Maximum DCA in 3D: ";
899 if(fMaxDCA3DFlag) listOfCuts += fMaxDCA3D;
900 else listOfCuts += "Not used";
901 l.DrawLatex(0.1,0.48,listOfCuts.Data());
902 listOfCuts = "Maximum DCA in xy (TPC): ";
903 if(fMaxDCAXYFlag) listOfCuts += fMaxDCAXYTPC;
904 else listOfCuts += "Not used";
905 l.DrawLatex(0.1,0.4,listOfCuts.Data());
906 listOfCuts = "Maximum DCA in z (TPC): ";
907 if(fMaxDCAZFlag) listOfCuts += fMaxDCAZTPC;
908 else listOfCuts += "Not used";
909 l.DrawLatex(0.1,0.32,listOfCuts.Data());
910 listOfCuts = "Maximum DCA in 3D (TPC): ";
911 if(fMaxDCA3DFlag) listOfCuts += fMaxDCA3DTPC;
912 else listOfCuts += "Not used";
913 l.DrawLatex(0.1,0.24,listOfCuts.Data());
914 listOfCuts = "Maximum constrained #chi^{2} (vertex): ";
915 if(fMaxConstrainChi2Flag) listOfCuts += fMaxConstrainChi2;
916 else listOfCuts += "Not used";
917 l.DrawLatex(0.1,0.16,listOfCuts.Data());
918
919 c->cd(6)->SetFillColor(6);
920 l.DrawLatex(0.3,0.9,"Tracking bits related cuts");
921 listOfCuts = "Request the ITS refit bit: "; listOfCuts += fITSRefitFlag;
922 l.DrawLatex(0.1,0.7,listOfCuts.Data());
923 listOfCuts = "Request the TPC refit bit: "; listOfCuts += fTPCRefitFlag;
924 l.DrawLatex(0.1,0.5,listOfCuts.Data());
925 listOfCuts = "Request the TPC pid bit: "; listOfCuts += fTPCpidFlag;
926 l.DrawLatex(0.1,0.3,listOfCuts.Data());
927 listOfCuts = "Request the ESD pid bit: "; listOfCuts += fESDpidFlag;
928 l.DrawLatex(0.1,0.1,listOfCuts.Data());
929
930 return c;
931}
932
933//________________________________________________________________________
934Bool_t AliProtonAnalysisBase::IsProton(AliESDtrack *track) {
935 //Function that checks if a track is a proton
99cdfc57 936 Double_t probability[5] = {0.,0.,0.,0.,0.};
87a55728 937 Double_t gPt = 0.0, gP = 0.0, gEta = 0.0;
0ab648ea 938 Long64_t fParticleType = 0;
939
940 //Bayesian approach for the PID
941 if(fProtonPIDMode == kBayesian) {
afc2ac17 942 if((fProtonAnalysisMode == kTPC)||(fProtonAnalysisMode == kHybrid)||(fProtonAnalysisMode == kFullHybrid)) {
0ab648ea 943 AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
944 if(tpcTrack) {
945 gPt = tpcTrack->Pt();
946 gP = tpcTrack->P();
947 track->GetTPCpid(probability);
948 }
99cdfc57 949 else return 0;
0ab648ea 950 }//TPC standalone or Hybrid TPC
951 else if(fProtonAnalysisMode == kGlobal) {
952 gPt = track->Pt();
953 gP = track->P();
954 track->GetESDpid(probability);
955 }//Global tracking
956
957 Double_t rcc = 0.0;
958 for(Int_t i = 0; i < AliPID::kSPECIES; i++)
959 rcc += probability[i]*GetParticleFraction(i,gP);
960 if(rcc != 0.0) {
961 Double_t w[5];
962 for(Int_t i = 0; i < AliPID::kSPECIES; i++)
963 w[i] = probability[i]*GetParticleFraction(i,gP)/rcc;
964 fParticleType = TMath::LocMax(AliPID::kSPECIES,w);
965 }
966 if(fParticleType == 4)
967 return kTRUE;
99cdfc57 968 }//Bayesian pid
0ab648ea 969 //Ratio of the measured over the theoretical dE/dx a la STAR
970 else if(fProtonPIDMode == kRatio) {
c6909683 971 AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
972 if(tpcTrack) {
973 gPt = tpcTrack->Pt();
b8cca027 974 gP = track->GetInnerParam()->P();
c6909683 975 gEta = tpcTrack->Eta();
976 }
99cdfc57 977 Double_t fAlephParameters[5] = {0.,0.,0.,0.,0.};
c6909683 978 if(fAnalysisMC) {
979 fAlephParameters[0] = 2.15898e+00/50.;
980 fAlephParameters[1] = 1.75295e+01;
981 fAlephParameters[2] = 3.40030e-09;
982 fAlephParameters[3] = 1.96178e+00;
983 fAlephParameters[4] = 3.91720e+00;
984 }
985 else {
986 fAlephParameters[0] = 0.0283086;
987 fAlephParameters[1] = 2.63394e+01;
988 fAlephParameters[2] = 5.04114e-11;
989 fAlephParameters[3] = 2.12543e+00;
990 fAlephParameters[4] = 4.88663e+00;
991 }
992
bc554c73 993 AliTPCPIDResponse *tpcResponse = new AliTPCPIDResponse();
994 tpcResponse->SetBetheBlochParameters(fAlephParameters[0],fAlephParameters[1],fAlephParameters[2],fAlephParameters[3],fAlephParameters[4]);
9a3204e2 995
996 Double_t normalizeddEdx = -10.;
bc554c73 997 if((track->GetTPCsignal() > 0.0) && (tpcResponse->GetExpectedSignal(gP,AliPID::kProton) > 0.0))
998 normalizeddEdx = TMath::Log(track->GetTPCsignal()/tpcResponse->GetExpectedSignal(gP,AliPID::kProton));
999
1000 delete tpcResponse;
1001
c6909683 1002 if(normalizeddEdx >= fNRatio)
1003 return kTRUE;
1004 }//kRatio PID mode
0ab648ea 1005 //Definition of an N-sigma area around the dE/dx vs P band
c6909683 1006 else if(fProtonPIDMode == kSigma) {
2c080079 1007 Double_t fAlephParameters[5];
0bb8f6b4 1008 if(fAnalysisMC) {
e56f08ed 1009 fAlephParameters[0] = 2.15898e+00/50.;
1010 fAlephParameters[1] = 1.75295e+01;
1011 fAlephParameters[2] = 3.40030e-09;
1012 fAlephParameters[3] = 1.96178e+00;
1013 fAlephParameters[4] = 3.91720e+00;
0bb8f6b4 1014 }
1015 else {
e56f08ed 1016 fAlephParameters[0] = 0.0283086;
1017 fAlephParameters[1] = 2.63394e+01;
1018 fAlephParameters[2] = 5.04114e-11;
1019 fAlephParameters[3] = 2.12543e+00;
1020 fAlephParameters[4] = 4.88663e+00;
0bb8f6b4 1021 }
2c080079 1022
30f87a5b 1023 Double_t nsigma = 100.0;
bc554c73 1024 AliTPCPIDResponse *tpcResponse = new AliTPCPIDResponse();
1025 tpcResponse->SetBetheBlochParameters(fAlephParameters[0],fAlephParameters[1],fAlephParameters[2],fAlephParameters[3],fAlephParameters[4]);
2c080079 1026
bc554c73 1027 Double_t mom = track->GetP();
1028 const AliExternalTrackParam *in = track->GetInnerParam();
1029 if (in)
1030 mom = in->GetP();
1031
1032 nsigma = TMath::Abs(tpcResponse->GetNumberOfSigmas(mom,track->GetTPCsignal(),track->GetTPCsignalN(),AliPID::kProton));
30f87a5b 1033
bc554c73 1034 delete tpcResponse;
2c080079 1035 if(nsigma <= fNSigma)
87a55728 1036 return kTRUE;
c6909683 1037 }//kSigma PID method
0ab648ea 1038
1039 return kFALSE;
1040}
1041
1042
1043