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