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