Update (Andrea)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAODPidHF.cxx
CommitLineData
7ce8802c 1/**************************************************************************
2 * Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 * *************************************************************************/
16
27de2dfb 17/* $Id$ */
18
7ce8802c 19//***********************************************************
20// Class AliAODPidHF
21// class for PID with AliAODRecoDecayHF
b257e565 22// Authors: D. Caffarri caffarri@pd.infn.it, A.Dainese andrea.dainese@pd.infn.it, S. Dash dash@to.infn.it, F. Prino prino@to.infn.it, R. Romita r.romita@gsi.de, Y. Wang yifei@pi0.physi.uni-heidelberg.de P. Antonioli pietro.antonioli@bo.infn.it
7ce8802c 23//***********************************************************
64e39959 24#include <TCanvas.h>
25
7ce8802c 26#include "AliAODPidHF.h"
27#include "AliAODPid.h"
7ce8802c 28#include "AliPID.h"
b79bfc3e 29#include "AliPIDResponse.h"
7ce8802c 30#include "AliAODpidUtil.h"
f9e3cf8e 31#include "AliESDtrack.h"
7ce8802c 32
33
34ClassImp(AliAODPidHF)
35
36//------------------------------
37AliAODPidHF::AliAODPidHF():
228f3aba 38 AliAODPid(),
39 fnNSigma(5),
40 fnSigma(),
41 fTOFSigma(160.),
42 fnPriors(5),
43 fPriors(),
44 fnPLimit(2),
45 fPLimit(),
cdac45d3 46 fAsym(kFALSE),
47 fTPC(kFALSE),
48 fTOF(kFALSE),
49 fITS(kFALSE),
50 fTRD(kFALSE),
79f3c225 51 fMatch(0),
913d1957 52 fCompat(kFALSE),
ad42e35b 53 fPCompatTOF(1.5),
02baac36 54 fnNSigmaCompat(2),
ad42e35b 55 fnSigmaCompat(),
12141139 56 fMC(kFALSE),
17c4a8d5 57 fOnePad(kFALSE),
5821ea5b 58 fMCLowEn2011(kFALSE),
67a82ef2 59 fppLowEn2011(kFALSE),
439d61b7 60 fPbPb(kFALSE),
b79bfc3e 61 fTOFdecide(kFALSE),
2c9eefd4 62 fOldPid(kTRUE),
63 fPtThresholdTPC(999999.),
11690a06 64 fPidResponse(0),
65 fPidCombined(new AliPIDCombined())
7ce8802c 66{
67 //
68 // Default constructor
69 //
228f3aba 70 fPLimit=new Double_t[fnPLimit];
71 fnSigma=new Double_t[fnNSigma];
72 fPriors=new Double_t[fnPriors];
02baac36 73 fnSigmaCompat=new Double_t[fnNSigmaCompat];
228f3aba 74
75 for(Int_t i=0;i<fnNSigma;i++){
76 fnSigma[i]=0.;
77 }
78 for(Int_t i=0;i<fnPriors;i++){
79 fPriors[i]=0.;
80 }
81 for(Int_t i=0;i<fnPLimit;i++){
82 fPLimit[i]=0.;
83 }
02baac36 84 for(Int_t i=0;i<fnNSigmaCompat;i++){
ad42e35b 85 fnSigmaCompat[i]=3.;
86 }
7ce8802c 87
88}
89//----------------------
90AliAODPidHF::~AliAODPidHF()
91{
92 // destructor
228f3aba 93 // if(fPLimit) delete fPLimit;
94 // if(fnSigma) delete fnSigma;
95 // if(fPriors) delete fPriors;
7ce8802c 96}
97//------------------------
98AliAODPidHF::AliAODPidHF(const AliAODPidHF& pid) :
99 AliAODPid(pid),
228f3aba 100 fnNSigma(pid.fnNSigma),
101 fnSigma(pid.fnSigma),
102 fTOFSigma(pid.fTOFSigma),
103 fnPriors(pid.fnPriors),
104 fPriors(pid.fPriors),
105 fnPLimit(pid.fnPLimit),
106 fPLimit(pid.fPLimit),
cdac45d3 107 fAsym(pid.fAsym),
108 fTPC(pid.fTPC),
109 fTOF(pid.fTOF),
110 fITS(pid.fITS),
111 fTRD(pid.fTRD),
112 fMatch(pid.fMatch),
913d1957 113 fCompat(pid.fCompat),
ad42e35b 114 fPCompatTOF(pid.fPCompatTOF),
02baac36 115 fnNSigmaCompat(pid.fnNSigmaCompat),
ad42e35b 116 fnSigmaCompat(pid.fnSigmaCompat),
12141139 117 fMC(pid.fMC),
17c4a8d5 118 fOnePad(pid.fOnePad),
5821ea5b 119 fMCLowEn2011(pid.fMCLowEn2011),
67a82ef2 120 fppLowEn2011(pid.fppLowEn2011),
439d61b7 121 fPbPb(pid.fPbPb),
b79bfc3e 122 fTOFdecide(pid.fTOFdecide),
123 fOldPid(pid.fOldPid),
2c9eefd4 124 fPtThresholdTPC(pid.fPtThresholdTPC),
11690a06 125 fPidResponse(pid.fPidResponse),
126 fPidCombined(pid.fPidCombined)
7ce8802c 127 {
128
129 for(Int_t i=0;i<5;i++){
130 fPriors[i]=pid.fPriors[i];
131 }
132
133 }
134
135//----------------------
228f3aba 136Int_t AliAODPidHF::RawSignalPID(AliAODTrack *track, TString detector) const{
137// raw PID for single detectors, returns the particle type with smaller sigma
f9e3cf8e 138 Int_t specie=-1;
139 if(detector.Contains("ITS")) return ApplyPidITSRaw(track,specie);
140 if(detector.Contains("TPC")) return ApplyPidTPCRaw(track,specie);
141 if(detector.Contains("TOF")) return ApplyPidTOFRaw(track,specie);
7ce8802c 142
f9e3cf8e 143 return specie;
7ce8802c 144
145}
146//---------------------------
228f3aba 147Bool_t AliAODPidHF::IsKaonRaw(AliAODTrack *track, TString detector) const{
148// checks if the track can be a kaon, raw PID applied for single detectors
7ce8802c 149 Int_t specie=0;
150
151 if(detector.Contains("ITS")) specie=ApplyPidITSRaw(track,3);
152 if(detector.Contains("TPC")) specie=ApplyPidTPCRaw(track,3);
153 if(detector.Contains("TOF")) specie=ApplyPidTOFRaw(track,3);
154
155 if(specie==3) return kTRUE;
156 return kFALSE;
157}
158//---------------------------
228f3aba 159Bool_t AliAODPidHF::IsPionRaw (AliAODTrack *track, TString detector) const{
160// checks if the track can be a pion, raw PID applied for single detectors
7ce8802c 161
162 Int_t specie=0;
163
164 if(detector.Contains("ITS")) specie=ApplyPidITSRaw(track,2);
165 if(detector.Contains("TPC")) specie=ApplyPidTPCRaw(track,2);
166 if(detector.Contains("TOF")) specie=ApplyPidTOFRaw(track,2);
167
168 if(specie==2) return kTRUE;
169 return kFALSE;
170}
171//---------------------------
228f3aba 172Bool_t AliAODPidHF::IsProtonRaw (AliAODTrack *track, TString detector) const{
173// checks if the track can be a proton raw PID applied for single detectors
7ce8802c 174
175 Int_t specie=0;
176 if(detector.Contains("ITS")) specie=ApplyPidITSRaw(track,4);
177 if(detector.Contains("TPC")) specie=ApplyPidTPCRaw(track,4);
178 if(detector.Contains("TOF")) specie=ApplyPidTOFRaw(track,4);
179
180 if(specie==4) return kTRUE;
181
182 return kFALSE;
183}
f9e3cf8e 184//--------------------------
228f3aba 185Bool_t AliAODPidHF::IsElectronRaw(AliAODTrack *track, TString detector) const{
186// checks if the track can be an electron raw PID applied for single detectors
f9e3cf8e 187
188 Int_t specie=-1;
189 if(detector.Contains("ITS")) specie=ApplyPidITSRaw(track,0);
190 if(detector.Contains("TPC")) specie=ApplyPidTPCRaw(track,0);
191 if(detector.Contains("TOF")) specie=ApplyPidTOFRaw(track,0);
192
193 if(specie==0) return kTRUE;
194
195 return kFALSE;
196}
197//--------------------------
228f3aba 198Int_t AliAODPidHF::ApplyPidTPCRaw(AliAODTrack *track,Int_t specie) const{
199// n-sigma cut, TPC PID
7ce8802c 200
cdac45d3 201 if(!CheckStatus(track,"TPC")) return 0;
b79bfc3e 202 Int_t pid=-1;
203 if(fOldPid){
7ce8802c 204 AliAODPid *pidObj = track->GetDetPid();
205
206 Double_t dedx=pidObj->GetTPCsignal();
207 Double_t mom = pidObj->GetTPCmomentum();
2c9eefd4 208 if(mom>fPtThresholdTPC) return 0;
7ce8802c 209 AliTPCPIDResponse tpcResponse;
12141139 210 SetBetheBloch(tpcResponse);
0aba57ed 211 UShort_t nTPCClus=pidObj->GetTPCsignalN();
212 if(nTPCClus==0) {nTPCClus=track->GetTPCNcls();}
12141139 213
f9e3cf8e 214 if(specie<0){ // from RawSignalPID : should return the particle specie to wich the de/dx is closer to the bethe-block curve -> performance to be checked
228f3aba 215 Double_t nsigmaMax=fnSigma[0];
f9e3cf8e 216 for(Int_t ipart=0;ipart<5;ipart++){
7ce8802c 217 AliPID::EParticleType type=AliPID::EParticleType(ipart);
0aba57ed 218 Double_t nsigma = TMath::Abs(tpcResponse.GetNumberOfSigmas(mom,dedx,nTPCClus,type));
228f3aba 219 if((nsigma<nsigmaMax) && (nsigma<fnSigma[0])) {
7ce8802c 220 pid=ipart;
221 nsigmaMax=nsigma;
222 }
223 }
224 }else{ // asks only for one particle specie
225 AliPID::EParticleType type=AliPID::EParticleType(specie);
0aba57ed 226 Double_t nsigma = TMath::Abs(tpcResponse.GetNumberOfSigmas(mom,dedx,nTPCClus,type));
228f3aba 227 if (nsigma>fnSigma[0]) {
f9e3cf8e 228 pid=-1;
7ce8802c 229 }else{
79f3c225 230 pid=specie;
7ce8802c 231 }
232 }
b79bfc3e 233 }else{ //old pid
234 if(specie<0){ // from RawSignalPID : should return the particle specie to wich the de/dx is closer to the bethe-block curve -> performance to be checked
235 Double_t nsigmaMax=fnSigma[0];
236 for(Int_t ipart=0;ipart<5;ipart++){
237 AliPID::EParticleType type=AliPID::EParticleType(ipart);
238 Double_t nsigma = TMath::Abs(fPidResponse->NumberOfSigmasTPC(track,type));
239 if((nsigma<nsigmaMax) && (nsigma<fnSigma[0])) {
240 pid=ipart;
241 nsigmaMax=nsigma;
242 }
243 }
244 }else{ // asks only for one particle specie
245 AliPID::EParticleType type=AliPID::EParticleType(specie);
246 Double_t nsigma = TMath::Abs(fPidResponse->NumberOfSigmasTPC(track,type));
247 if (nsigma>fnSigma[0]) {
248 pid=-1;
249 }else{
250 pid=specie;
251 }
252 }
253
254 } //new pid
7ce8802c 255
256 return pid;
257
258}
259//----------------------------
228f3aba 260Int_t AliAODPidHF::ApplyPidITSRaw(AliAODTrack *track,Int_t specie) const{
261// truncated mean, ITS PID
f9e3cf8e 262
cdac45d3 263 if(!CheckStatus(track,"ITS")) return 0;
b79bfc3e 264 Int_t pid=-1;
f9e3cf8e 265
b79bfc3e 266 if(fOldPid){
7ce8802c 267 Double_t mom=track->P();
268 AliAODPid *pidObj = track->GetDetPid();
269
270 Double_t dedx=pidObj->GetITSsignal();
913d1957 271 UChar_t clumap=track->GetITSClusterMap();
272 Int_t nPointsForPid=0;
273 for(Int_t i=2; i<6; i++){
274 if(clumap&(1<<i)) ++nPointsForPid;
275 }
276
277 Bool_t isSA=kTRUE;
278 if(track->GetStatus() & AliESDtrack::kTPCin) isSA = kFALSE;
279
7ce8802c 280 AliITSPIDResponse itsResponse;
f9e3cf8e 281 if(specie<0){ // from RawSignalPID : should return the particle specie to wich the de/dx is closer to the bethe-block curve -> performance to be checked
228f3aba 282 Double_t nsigmaMax=fnSigma[4];
f9e3cf8e 283 for(Int_t ipart=0;ipart<5;ipart++){
7ce8802c 284 AliPID::EParticleType type=AliPID::EParticleType(ipart);
913d1957 285 Double_t nsigma = TMath::Abs(itsResponse.GetNumberOfSigmas(mom,dedx,type,nPointsForPid,isSA));
228f3aba 286 if((nsigma<nsigmaMax) && (nsigma<fnSigma[4])) {
7ce8802c 287 pid=ipart;
288 nsigmaMax=nsigma;
289 }
290 }
291 }else{ // asks only for one particle specie
292 AliPID::EParticleType type=AliPID::EParticleType(specie);
293 Double_t nsigma = TMath::Abs(itsResponse.GetNumberOfSigmas(mom,dedx,type));
228f3aba 294 if (nsigma>fnSigma[4]) {
f9e3cf8e 295 pid=-1;
7ce8802c 296 }else{
79f3c225 297 pid=specie;
7ce8802c 298 }
299 }
b79bfc3e 300 }else{ // old pid
301
302 if(specie<0){ // from RawSignalPID : should return the particle specie to wich the de/dx is closer to the bethe-block curve -> performance to be checked
303 Double_t nsigmaMax=fnSigma[4];
304 for(Int_t ipart=0;ipart<5;ipart++){
305 AliPID::EParticleType type=AliPID::EParticleType(ipart);
306 Double_t nsigma = TMath::Abs(fPidResponse->NumberOfSigmasITS(track,type));
307 if((nsigma<nsigmaMax) && (nsigma<fnSigma[4])) {
308 pid=ipart;
309 nsigmaMax=nsigma;
310 }
311 }
312 }else{ // asks only for one particle specie
313 AliPID::EParticleType type=AliPID::EParticleType(specie);
314 Double_t nsigma = TMath::Abs(fPidResponse->NumberOfSigmasITS(track,type));
315 if (nsigma>fnSigma[4]) {
316 pid=-1;
317 }else{
318 pid=specie;
319 }
320 }
321 } //new pid
322
7ce8802c 323 return pid;
324}
325//----------------------------
228f3aba 326Int_t AliAODPidHF::ApplyPidTOFRaw(AliAODTrack *track,Int_t specie) const{
327// n-sigma cut, TOF PID
f9e3cf8e 328
cdac45d3 329 if(!CheckStatus(track,"TOF")) return 0;
f9e3cf8e 330
7ce8802c 331 Double_t time[AliPID::kSPECIESN];
b257e565 332 Double_t sigmaTOFPid[AliPID::kSPECIES];
7ce8802c 333 AliAODPid *pidObj = track->GetDetPid();
334 pidObj->GetIntegratedTimes(time);
228f3aba 335 Double_t sigTOF=pidObj->GetTOFsignal();
b257e565 336 pidObj->GetTOFpidResolution(sigmaTOFPid);
337
f9e3cf8e 338 Int_t pid=-1;
7ce8802c 339
b257e565 340 if(specie<0){
341 Double_t sigmaTOFtrack;
342 if (sigmaTOFPid[4]>0) sigmaTOFtrack=sigmaTOFPid[4];
343 else sigmaTOFtrack=fTOFSigma;
344 Double_t nsigmaMax=sigmaTOFtrack*fnSigma[3];
f9e3cf8e 345 for(Int_t ipart=0;ipart<5;ipart++){
228f3aba 346 Double_t nsigma=TMath::Abs(sigTOF-time[ipart]);
b257e565 347 if (sigmaTOFPid[ipart]>0) sigmaTOFtrack=sigmaTOFPid[ipart];
348 else sigmaTOFtrack=fTOFSigma; // backward compatibility for old AODs
349 if((nsigma<nsigmaMax) && (nsigma<fnSigma[3]*sigmaTOFtrack)) {
7ce8802c 350 pid=ipart;
351 nsigmaMax=nsigma;
352 }
353 }
354 }else{ // asks only for one particle specie
228f3aba 355 Double_t nsigma=TMath::Abs(sigTOF-time[specie]);
b257e565 356 Double_t sigmaTOFtrack;
357 if (sigmaTOFPid[specie]>0) sigmaTOFtrack=sigmaTOFPid[specie];
358 else sigmaTOFtrack=fTOFSigma; // backward compatibility for old AODs
359 if (nsigma>fnSigma[3]*sigmaTOFtrack) {
360 pid=-1;
361 }else{
362 pid=specie;
363 }
7ce8802c 364 }
365 return pid;
366
367}
368//------------------------------
228f3aba 369void AliAODPidHF::CombinedProbability(AliAODTrack *track,Bool_t *type) const{
370// combined PID stored inside the AOD track
7ce8802c 371
372 const Double_t *pid=track->PID();
373 Float_t max=0.;
374 Int_t k=-1;
375 for (Int_t i=0; i<10; i++) {
376 if (pid[i]>max) {k=i; max=pid[i];}
377 }
378
379 if(k==2) type[0]=kTRUE;
380 if(k==3) type[1]=kTRUE;
381 if(k==4) type[2]=kTRUE;
382
383 return;
384}
385//--------------------
cdac45d3 386void AliAODPidHF::BayesianProbability(AliAODTrack *track,Double_t *pid) const{
228f3aba 387// bayesian PID for single detectors or combined
7ce8802c 388
cdac45d3 389 if(fITS && !fTPC && !fTOF) {BayesianProbabilityITS(track,pid);return;}
390 if(fTPC && !fITS && !fTOF) {BayesianProbabilityTPC(track,pid);return;}
391 if(fTOF && !fITS && !fTPC) {BayesianProbabilityTOF(track,pid);return;}
392
393 Double_t probITS[5]={1.,1.,1.,1.,1.};
394 Double_t probTPC[5]={1.,1.,1.,1.,1.};
395 Double_t probTOF[5]={1.,1.,1.,1.,1.};
396 if(fITS) BayesianProbabilityITS(track,probITS);
397 if(fTPC) BayesianProbabilityTPC(track,probTPC);
398 if(fTOF) BayesianProbabilityTOF(track,probTOF);
7ce8802c 399 Double_t probTot[5]={0.,0.,0.,0.,0.};
400 for(Int_t i=0;i<5;i++){
401 probTot[i]=probITS[i]*probTPC[i]*probTOF[i];
402 }
403 for(Int_t i2=0;i2<5;i2++){
404 pid[i2]=probTot[i2]*fPriors[i2]/(probTot[0]*fPriors[0]+probTot[1]*fPriors[1]+probTot[2]*fPriors[2]+probTot[3]*fPriors[3]+probTot[4]*fPriors[4]);
405 }
7ce8802c 406
407 return;
408
409}
410//------------------------------------
228f3aba 411void AliAODPidHF::BayesianProbabilityITS(AliAODTrack *track,Double_t *prob) const{
7ce8802c 412
228f3aba 413// bayesian PID for ITS
7ce8802c 414 AliAODpidUtil pid;
415 Double_t itspid[AliPID::kSPECIES];
416 pid.MakeITSPID(track,itspid);
417 for(Int_t ind=0;ind<AliPID::kSPECIES;ind++){
369ad90e 418 if(fTOF || fTPC || fTRD){
419 prob[ind]=itspid[ind];
420 }else{
421 prob[ind]=itspid[ind]*fPriors[ind]/(itspid[0]*fPriors[0]+itspid[1]*fPriors[1]+itspid[2]*fPriors[2]+itspid[3]*fPriors[3]+itspid[4]*fPriors[4]);
422 }
7ce8802c 423 }
424 return;
425
426}
427//------------------------------------
228f3aba 428void AliAODPidHF::BayesianProbabilityTPC(AliAODTrack *track,Double_t *prob) const{
429// bayesian PID for TPC
7ce8802c 430
431 AliAODpidUtil pid;
432 Double_t tpcpid[AliPID::kSPECIES];
433 pid.MakeTPCPID(track,tpcpid);
434 for(Int_t ind=0;ind<AliPID::kSPECIES;ind++){
369ad90e 435 if(fTOF || fITS || fTRD){
436 prob[ind]=tpcpid[ind];
7ce8802c 437 }else{
369ad90e 438 prob[ind]=tpcpid[ind]*fPriors[ind]/(tpcpid[0]*fPriors[0]+tpcpid[1]*fPriors[1]+tpcpid[2]*fPriors[2]+tpcpid[3]*fPriors[3]+tpcpid[4]*fPriors[4]);
7ce8802c 439 }
369ad90e 440}
7ce8802c 441 return;
442
443}
444//------------------------------------
228f3aba 445void AliAODPidHF::BayesianProbabilityTOF(AliAODTrack *track,Double_t *prob) const{
446// bayesian PID for TOF
7ce8802c 447
448 AliAODpidUtil pid;
449 Double_t tofpid[AliPID::kSPECIES];
228f3aba 450 pid.MakeTOFPID(track,tofpid);
7ce8802c 451 for(Int_t ind=0;ind<AliPID::kSPECIES;ind++){
369ad90e 452 if(fTPC || fITS || fTRD){
453 prob[ind]=tofpid[ind];
454 }else{
7ce8802c 455 prob[ind]=tofpid[ind]*fPriors[ind]/(tofpid[0]*fPriors[0]+tofpid[1]*fPriors[1]+tofpid[2]*fPriors[2]+tofpid[3]*fPriors[3]+tofpid[4]*fPriors[4]);
456 }
369ad90e 457}
7ce8802c 458 return;
459
460}
f9e3cf8e 461//---------------------------------
228f3aba 462void AliAODPidHF::BayesianProbabilityTRD(AliAODTrack *track,Double_t *prob) const{
463// bayesian PID for TRD
f9e3cf8e 464
465 AliAODpidUtil pid;
466 Double_t trdpid[AliPID::kSPECIES];
467 pid.MakeTRDPID(track,trdpid);
468 for(Int_t ind=0;ind<AliPID::kSPECIES;ind++){
369ad90e 469 if(fTPC || fITS || fTOF){
470 prob[ind]=trdpid[ind];
f9e3cf8e 471 }else{
369ad90e 472 prob[ind]=trdpid[ind]*fPriors[ind]/(trdpid[0]*fPriors[0]+trdpid[1]*fPriors[1]+trdpid[2]*fPriors[2]+trdpid[3]*fPriors[3]+trdpid[4]*fPriors[4]);
f9e3cf8e 473 }
369ad90e 474}
f9e3cf8e 475 return;
476
477 }
478//--------------------------------
228f3aba 479Bool_t AliAODPidHF::CheckStatus(AliAODTrack *track,TString detectors) const{
f9e3cf8e 480
228f3aba 481// Quality cuts on the tracks, detector by detector
f9e3cf8e 482
483 if(detectors.Contains("ITS")){
484 if ((track->GetStatus()&AliESDtrack::kITSin)==0) return kFALSE;
485 UChar_t clumap=track->GetITSClusterMap();
486 Int_t nPointsForPid=0;
487 for(Int_t i=2; i<6; i++){
488 if(clumap&(1<<i)) ++nPointsForPid;
489 }
913d1957 490 if(nPointsForPid<3) return kFALSE;
f9e3cf8e 491 }
492
493 if(detectors.Contains("TPC")){
494 if ((track->GetStatus()&AliESDtrack::kTPCin )==0) return kFALSE;
495 UShort_t nTPCClus=track->GetTPCClusterMap().CountBits();
496 if (nTPCClus<70) return kFALSE;
497 }
498
499 if(detectors.Contains("TOF")){
500 if ((track->GetStatus()&AliESDtrack::kTOFout )==0) return kFALSE;
501 if ((track->GetStatus()&AliESDtrack::kTIME )==0) return kFALSE;
502 if ((track->GetStatus()&AliESDtrack::kTOFpid )==0) return kFALSE;
b257e565 503 if (!(track->GetStatus()&AliESDtrack::kTOFmismatch)==0) return kFALSE;
f9e3cf8e 504 }
505
506
507 if(detectors.Contains("TRD")){
508 if ((track->GetStatus()&AliESDtrack::kTRDout )==0) return kFALSE;
59a8e853 509 UChar_t ntracklets = track->GetTRDntrackletsPID();
510 if(ntracklets<4) return kFALSE;
f9e3cf8e 511 }
512
513 return kTRUE;
514}
228f3aba 515//--------------------------------------------
516Bool_t AliAODPidHF::TPCRawAsym(AliAODTrack* track,Int_t specie) const{
517// TPC nsigma cut PID, different sigmas in different p bins
518
519 if(!CheckStatus(track,"TPC")) return kFALSE;
520 AliAODPid *pidObj = track->GetDetPid();
521 Double_t mom = pidObj->GetTPCmomentum();
2c9eefd4 522 if(mom>fPtThresholdTPC) return 0;
b79bfc3e 523 Double_t nsigma=999.;
524 if(fOldPid){
228f3aba 525 Double_t dedx=pidObj->GetTPCsignal();
0aba57ed 526 UShort_t nTPCClus=pidObj->GetTPCsignalN();
527 if(nTPCClus==0) {nTPCClus=track->GetTPCNcls();}
228f3aba 528
529 AliTPCPIDResponse tpcResponse;
12141139 530 SetBetheBloch(tpcResponse);
228f3aba 531 AliPID::EParticleType type=AliPID::EParticleType(specie);
b79bfc3e 532 nsigma = TMath::Abs(tpcResponse.GetNumberOfSigmas(mom,dedx,nTPCClus,type));
533
534 }else{ //old pid
535 AliPID::EParticleType type=AliPID::EParticleType(specie);
536 nsigma = TMath::Abs(fPidResponse->NumberOfSigmasTPC(track,type));
537 } //new pid
228f3aba 538
539 if(mom<fPLimit[0] && nsigma<fnSigma[0]) return kTRUE;
540 if(mom<fPLimit[1] && mom>fPLimit[0] && nsigma<fnSigma[1]) return kTRUE;
541 if(mom>fPLimit[1] && nsigma<fnSigma[2]) return kTRUE;
542
543 return kFALSE;
544}
545//------------------
546Int_t AliAODPidHF::MatchTPCTOF(AliAODTrack *track,Int_t mode,Int_t specie,Bool_t compat){
547// combination of the PID info coming from TPC and TOF
548 if(mode==1){
549 //TOF || TPC (a la' Andrea R.)
cdac45d3 550 // convention:
228f3aba 551 // for the single detectors: -1 = kFALSE, 1 = kTRUE, 0 = compatible
552 // the method returns the sum of the response of the 2 detectors
79f3c225 553 if(fTPC && fTOF) {if(!CheckStatus(track,"TPC") && !CheckStatus(track,"TOF")) return 0;}
228f3aba 554
cdac45d3 555
1cb0f658 556 Int_t tTPCinfo=0;
cdac45d3 557 if(fTPC){
228f3aba 558 if(CheckStatus(track,"TPC")) {
559 if(fAsym) {
3779f50a 560 if(TPCRawAsym(track,specie)) {
1cb0f658 561 tTPCinfo=1;
3779f50a 562 }else{
1cb0f658 563 tTPCinfo=-1;
3779f50a 564 }
228f3aba 565 }else{
3779f50a 566 if(specie==2 && IsPionRaw(track,"TPC")) {
1cb0f658 567 tTPCinfo=1;
3779f50a 568 }else{
1cb0f658 569 tTPCinfo=-1;
3779f50a 570 }
571 if(specie==3 && IsKaonRaw(track,"TPC")) {
1cb0f658 572 tTPCinfo=1;
3779f50a 573 }else{
1cb0f658 574 tTPCinfo=-1;
3779f50a 575 }
576 if(specie==4 && IsProtonRaw(track,"TPC")) {
1cb0f658 577 tTPCinfo=1;
3779f50a 578 }else{
1cb0f658 579 tTPCinfo=-1;
3779f50a 580 }
581
228f3aba 582 }
583
584
1cb0f658 585 if(compat && tTPCinfo<0){
cdac45d3 586 Double_t sig0tmp=fnSigma[0];
ad42e35b 587 SetSigma(0,fnSigmaCompat[0]);
1cb0f658 588 if(specie==2 && IsPionRaw(track,"TPC")) tTPCinfo=0;
589 if(specie==3 && IsKaonRaw(track,"TPC")) tTPCinfo=0;
590 if(specie==4 && IsProtonRaw(track,"TPC")) tTPCinfo=0;
cdac45d3 591 SetSigma(0,sig0tmp);
592 }
228f3aba 593
cdac45d3 594 }
228f3aba 595 }
596
1cb0f658 597 Int_t tTOFinfo=0;
cdac45d3 598 if(fTOF){
1cb0f658 599 if(!CheckStatus(track,"TOF") && fTPC) return tTPCinfo;
228f3aba 600
1cb0f658 601 tTOFinfo=-1;
228f3aba 602
1cb0f658 603 if(specie==2 && IsPionRaw(track,"TOF")) tTOFinfo=1;
604 if(specie==3 && IsKaonRaw(track,"TOF")) tTOFinfo=1;
605 if(specie==4 && IsProtonRaw(track,"TOF")) tTOFinfo=1;
228f3aba 606
1cb0f658 607 if(compat && tTOFinfo>0){
cdac45d3 608 Double_t ptrack=track->P();
ad42e35b 609 if(ptrack>fPCompatTOF) {
610 Double_t sig0tmp=fnSigma[3];
611 SetSigma(3,fnSigmaCompat[1]);
612 if(specie==2 && IsPionRaw(track,"TOF")) tTOFinfo=0;
613 if(specie==3 && IsKaonRaw(track,"TOF")) tTOFinfo=0;
614 if(specie==4 && IsProtonRaw(track,"TOF")) tTOFinfo=0;
615 SetSigma(3,sig0tmp);
616 }
cdac45d3 617 }
79f3c225 618 }
913d1957 619
439d61b7 620
621 if(tTPCinfo+tTOFinfo==0 && fTOFdecide){
622 if(!CheckStatus(track,"TOF")) return tTPCinfo;
623 return tTOFinfo;
624 }
625
913d1957 626 if(tTPCinfo+tTOFinfo==0 && fITS){
627 if(!CheckStatus(track,"ITS")) return tTPCinfo+tTOFinfo;
628 Int_t tITSinfo = -1;
2a5813a1 629 if(specie==2 && IsPionRaw(track,"ITS")) tITSinfo=1;
630 if(specie==3 && IsKaonRaw(track,"ITS")) tITSinfo=1;
631 if(specie==4 && IsProtonRaw(track,"ITS")) tITSinfo=1;
913d1957 632 return tITSinfo;
633 }
634
1cb0f658 635 return tTPCinfo+tTOFinfo;
cdac45d3 636}
228f3aba 637 if(mode==2){
638 //TPC & TOF (a la' Yifei)
cdac45d3 639 // convention: -1 = kFALSE, 1 = kTRUE, 0 = not identified
1cb0f658 640 Int_t tTPCinfo=0;
cdac45d3 641
642 if(fTPC && CheckStatus(track,"TPC")) {
1cb0f658 643 tTPCinfo=1;
228f3aba 644 if(fAsym){
1cb0f658 645 if(!TPCRawAsym(track,specie)) tTPCinfo=-1;
228f3aba 646 }else{
1cb0f658 647 if(specie==2 && !IsPionRaw(track,"TPC")) tTPCinfo=-1;
648 if(specie==3 && !IsKaonRaw(track,"TPC")) tTPCinfo=-1;
649 if(specie==4 && !IsProtonRaw(track,"TPC")) tTPCinfo=-1;
228f3aba 650 }
651 }
652
1cb0f658 653 Int_t tTOFinfo=1;
cdac45d3 654 if(fTOF){
1cb0f658 655 if(fTPC && !CheckStatus(track,"TOF")) return tTPCinfo;
228f3aba 656
1cb0f658 657 if(specie==2 && !IsPionRaw(track,"TOF")) tTOFinfo=-1;
658 if(specie==3 && !IsKaonRaw(track,"TOF")) tTOFinfo=-1;
659 if(specie==4 && !IsProtonRaw(track,"TOF")) tTOFinfo=-1;
cdac45d3 660 }
913d1957 661
662 if(tTOFinfo==1 && tTPCinfo==1) return 1;
663
664 if(tTPCinfo+tTOFinfo==0 && fITS){
665 if(!CheckStatus(track,"ITS")) return tTPCinfo+tTOFinfo;
666 Int_t tITSinfo = -1;
2a5813a1 667 if(specie==2 && IsPionRaw(track,"ITS")) tITSinfo=1;
668 if(specie==3 && IsKaonRaw(track,"ITS")) tITSinfo=1;
669 if(specie==4 && IsProtonRaw(track,"ITS")) tITSinfo=1;
913d1957 670 return tITSinfo;
671 }
228f3aba 672
228f3aba 673 return -1;
674
cdac45d3 675}
228f3aba 676
677 if(mode==3){
678 //TPC for p<fPLimit[0], TOF for p>=fPLimit[0] (a la' Andrea A.)
679 // convention (temporary): -1 = kFALSE, 1 = kTRUE, 0 = not identified
680
cdac45d3 681 if(fTPC && fTOF) if(!CheckStatus(track,"TPC") && !CheckStatus(track,"TOF")) return 0;
228f3aba 682
683 Double_t ptrack=track->P();
684
685
1cb0f658 686 Int_t tTPCinfo=-1;
228f3aba 687
913d1957 688 if(ptrack>=fPLimit[0] && ptrack<fPLimit[1] && fTPC) {
228f3aba 689 if(!CheckStatus(track,"TPC")) return 0;
690 if(fAsym) {
1cb0f658 691 if(TPCRawAsym(track,specie)) tTPCinfo=1;
228f3aba 692 }else{
1cb0f658 693 if(specie==2 && IsPionRaw(track,"TPC")) tTPCinfo=1;
694 if(specie==3 && IsKaonRaw(track,"TPC")) tTPCinfo=1;
695 if(specie==4 && IsProtonRaw(track,"TPC")) tTPCinfo=1;
228f3aba 696 }
1cb0f658 697 return tTPCinfo;
228f3aba 698 }
699
1cb0f658 700 Int_t tTOFinfo=-1;
913d1957 701 if(ptrack>=fPLimit[1] && fTOF){
228f3aba 702 if(!CheckStatus(track,"TOF")) return 0;
1cb0f658 703 if(specie==2 && IsPionRaw(track,"TOF")) tTOFinfo=1;
704 if(specie==3 && IsKaonRaw(track,"TOF")) tTOFinfo=1;
705 if(specie==4 && IsProtonRaw(track,"TOF")) tTOFinfo=1;
706 return tTOFinfo;
228f3aba 707 }
708
913d1957 709 Int_t tITSinfo=-1;
710 if(ptrack<fPLimit[0] && fITS){
711 if(!CheckStatus(track,"ITS")) return 0;
712 if(specie==2 && IsPionRaw(track,"ITS")) tITSinfo=1;
713 if(specie==3 && IsKaonRaw(track,"ITS")) tITSinfo=1;
714 if(specie==4 && IsProtonRaw(track,"ITS")) tITSinfo=1;
715 return tITSinfo;
716 }
717
228f3aba 718 }
719
720 return -1;
cdac45d3 721
722}
723//----------------------------------
724Int_t AliAODPidHF::MakeRawPid(AliAODTrack *track, Int_t specie){
1cb0f658 725// general method to compute PID
cdac45d3 726 if(fMatch>0){
1cb0f658 727 return MatchTPCTOF(track,fMatch,specie,fCompat);
cdac45d3 728 }else{
729 if(fTPC && !fTOF && !fITS) {
ad42e35b 730 Int_t tTPCres=0;
731 if(!fAsym){
732 tTPCres=ApplyPidTPCRaw(track,specie);
733 }else{
734 if(TPCRawAsym(track,specie)) {
735 tTPCres=1;
736 }else{
737 tTPCres=-1;
738 }
739 }
1cb0f658 740 if(tTPCres==specie){return 1;}else{return tTPCres;};
cdac45d3 741 }else{
742 AliError("You should enable just one detector if you don't want to match");
743 return 0;
744 }
745 if(fTOF && !fTPC && !fITS) {
1cb0f658 746 Int_t tTOFres=ApplyPidTOFRaw(track,specie);
747 if(tTOFres==specie){return 1;}else{return tTOFres;};
cdac45d3 748 }else{
749 AliError("You should enable just one detector if you don't want to match");
750 return 0;
751 }
752
753 if(fITS && !fTPC && !fTOF) {
1cb0f658 754 Int_t tITSres=ApplyPidITSRaw(track,specie);
755 if(tITSres==specie){return 1;}else{return tITSres;};
cdac45d3 756 }else{
757 AliError("You should enable just one detector if you don't want to match");
758 return 0;
759 }
760 }
761
228f3aba 762}
12141139 763//--------------------------------------------
bad87dce 764void AliAODPidHF::SetBetheBloch(AliTPCPIDResponse &tpcResp) const{
12141139 765
766 Double_t alephParameters[5];
98db0ef2 767
32a978c2 768 if(fMC) { // MC
12141139 769
32a978c2 770 if(fPbPb) { // PbPb MC
771
9ea811e0 772 alephParameters[0] = 1.44405/50.;
32a978c2 773 alephParameters[1] = 2.35409e+01;
774 alephParameters[2] = TMath::Exp(-2.90330e+01);
775 alephParameters[3] = 2.10681e+00;
776 alephParameters[4] = 4.62254e+00;
12141139 777
32a978c2 778 } else { // pp MC
5821ea5b 779 if(fMCLowEn2011){
780 alephParameters[0]=0.0207667;
781 alephParameters[1]=29.9936;
782 alephParameters[2]=3.87866e-11;
783 alephParameters[3]=2.17291;
784 alephParameters[4]=7.1623;
3d4afb8c 785 }else if(fOnePad){
786 alephParameters[0]=0.029021;
787 alephParameters[1]=25.4181;
788 alephParameters[2]=4.66596e-08;
789 alephParameters[3]=1.90008;
790 alephParameters[4]=4.63783;
5821ea5b 791 }else{
792 alephParameters[0] = 2.15898/50.;
793 alephParameters[1] = 1.75295e+01;
794 alephParameters[2] = 3.40030e-09;
795 alephParameters[3] = 1.96178e+00;
796 alephParameters[4] = 3.91720e+00;
797 }
32a978c2 798 }
12141139 799
32a978c2 800 } else { // Real Data
801
802 if(fOnePad) { // pp 1-pad (since LHC10d)
803
67a82ef2 804 alephParameters[0] =1.34490e+00/50.;
805 alephParameters[1] = 2.69455e+01;
806 alephParameters[2] = TMath::Exp(-2.97552e+01);
807 alephParameters[3] = 2.35339e+00;
808 alephParameters[4] = 5.98079e+00;
809
810 } else if(fPbPb) { // PbPb
811
812 // alephParameters[0] = 1.25202/50.;
813 // alephParameters[1] = 2.74992e+01;
814 // alephParameters[2] = TMath::Exp(-3.31517e+01);
815 // alephParameters[3] = 2.46246;
816 // alephParameters[4] = 6.78938;
817
818 alephParameters[0] = 5.10207e+00/50.;
819 alephParameters[1] = 7.94982e+00;
820 alephParameters[2] = TMath::Exp(-9.07942e+00);
821 alephParameters[3] = 2.38808e+00;
822 alephParameters[4] = 1.68165e+00;
823
824 } else if(fppLowEn2011){ // pp low energy
825
826 alephParameters[0]=0.031642;
827 alephParameters[1]=22.353;
828 alephParameters[2]=4.16239e-12;
829 alephParameters[3]=2.61952;
830 alephParameters[4]=5.76086;
831
832 } else { // pp no 1-pad (LHC10bc)
833
834 alephParameters[0] = 0.0283086/0.97;
835 alephParameters[1] = 2.63394e+01;
836 alephParameters[2] = 5.04114e-11;
837 alephParameters[3] = 2.12543e+00;
838 alephParameters[4] = 4.88663e+00;
839
840 }
98db0ef2 841
12141139 842 }
843
844 tpcResp.SetBetheBlochParameters(alephParameters[0],alephParameters[1],alephParameters[2],alephParameters[3],alephParameters[4]);
845
846
847 return;
848
849}
d50b25df 850//-----------------------
851Bool_t AliAODPidHF::IsTOFPiKexcluded(AliAODTrack *track,Double_t nsigmaK){
852
853
854 if(!CheckStatus(track,"TOF")) return 0;
855
856 Double_t time[AliPID::kSPECIESN];
857 Double_t sigmaTOFPid[AliPID::kSPECIES];
858 AliAODPid *pidObj = track->GetDetPid();
859 pidObj->GetIntegratedTimes(time);
860 Double_t sigTOF=pidObj->GetTOFsignal();
861 pidObj->GetTOFpidResolution(sigmaTOFPid);
862 Double_t sigmaTOFtrack;
863 if (sigmaTOFPid[3]>0) sigmaTOFtrack=sigmaTOFPid[3];
864 else sigmaTOFtrack=fTOFSigma; // backward compatibility for old AODs
865
866 if((sigTOF-time[3])>nsigmaK*sigmaTOFtrack)return kTRUE;// K, Pi excluded (->LIKELY A PROTON)
867
868 return kFALSE;
869
870}
64e39959 871//--------------------------------------------------------------------------
872void AliAODPidHF::SetPriorDistribution(AliPID::EParticleType type,TH1F *prior){
873
874 //
875 // method setting the prior distributions to the AliPIDCombined object of the AliAODPidHF data member
876 // all the checks are done directly in the AliPIDCombined object
877 //
878
879 GetPidCombined()->SetPriorDistribution(type,prior);
880}
881//--------------------------------------------------------------------------
882void AliAODPidHF::DrawPrior(AliPID::EParticleType type){
883
884 //
885 // Drawing prior distribution for type "type"
886
887 new TCanvas();
888 GetPidCombined()->GetPriorDistribution(type)->Draw();
889}
d50b25df 890