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