Update (Andrea)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAODPidHF.cxx
... / ...
CommitLineData
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
17/* $Id$ */
18
19//***********************************************************
20// Class AliAODPidHF
21// class for PID with AliAODRecoDecayHF
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
23//***********************************************************
24#include <TCanvas.h>
25
26#include "AliAODPidHF.h"
27#include "AliAODPid.h"
28#include "AliPID.h"
29#include "AliPIDResponse.h"
30#include "AliAODpidUtil.h"
31#include "AliESDtrack.h"
32
33
34ClassImp(AliAODPidHF)
35
36//------------------------------
37AliAODPidHF::AliAODPidHF():
38 AliAODPid(),
39 fnNSigma(5),
40 fnSigma(),
41 fTOFSigma(160.),
42 fnPriors(5),
43 fPriors(),
44 fnPLimit(2),
45 fPLimit(),
46 fAsym(kFALSE),
47 fTPC(kFALSE),
48 fTOF(kFALSE),
49 fITS(kFALSE),
50 fTRD(kFALSE),
51 fMatch(0),
52 fCompat(kFALSE),
53 fPCompatTOF(1.5),
54 fnNSigmaCompat(2),
55 fnSigmaCompat(),
56 fMC(kFALSE),
57 fOnePad(kFALSE),
58 fMCLowEn2011(kFALSE),
59 fppLowEn2011(kFALSE),
60 fPbPb(kFALSE),
61 fTOFdecide(kFALSE),
62 fOldPid(kTRUE),
63 fPtThresholdTPC(999999.),
64 fPidResponse(0),
65 fPidCombined(new AliPIDCombined())
66{
67 //
68 // Default constructor
69 //
70 fPLimit=new Double_t[fnPLimit];
71 fnSigma=new Double_t[fnNSigma];
72 fPriors=new Double_t[fnPriors];
73 fnSigmaCompat=new Double_t[fnNSigmaCompat];
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 }
84 for(Int_t i=0;i<fnNSigmaCompat;i++){
85 fnSigmaCompat[i]=3.;
86 }
87
88}
89//----------------------
90AliAODPidHF::~AliAODPidHF()
91{
92 // destructor
93 // if(fPLimit) delete fPLimit;
94 // if(fnSigma) delete fnSigma;
95 // if(fPriors) delete fPriors;
96}
97//------------------------
98AliAODPidHF::AliAODPidHF(const AliAODPidHF& pid) :
99 AliAODPid(pid),
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),
107 fAsym(pid.fAsym),
108 fTPC(pid.fTPC),
109 fTOF(pid.fTOF),
110 fITS(pid.fITS),
111 fTRD(pid.fTRD),
112 fMatch(pid.fMatch),
113 fCompat(pid.fCompat),
114 fPCompatTOF(pid.fPCompatTOF),
115 fnNSigmaCompat(pid.fnNSigmaCompat),
116 fnSigmaCompat(pid.fnSigmaCompat),
117 fMC(pid.fMC),
118 fOnePad(pid.fOnePad),
119 fMCLowEn2011(pid.fMCLowEn2011),
120 fppLowEn2011(pid.fppLowEn2011),
121 fPbPb(pid.fPbPb),
122 fTOFdecide(pid.fTOFdecide),
123 fOldPid(pid.fOldPid),
124 fPtThresholdTPC(pid.fPtThresholdTPC),
125 fPidResponse(pid.fPidResponse),
126 fPidCombined(pid.fPidCombined)
127 {
128
129 for(Int_t i=0;i<5;i++){
130 fPriors[i]=pid.fPriors[i];
131 }
132
133 }
134
135//----------------------
136Int_t AliAODPidHF::RawSignalPID(AliAODTrack *track, TString detector) const{
137// raw PID for single detectors, returns the particle type with smaller sigma
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);
142
143 return specie;
144
145}
146//---------------------------
147Bool_t AliAODPidHF::IsKaonRaw(AliAODTrack *track, TString detector) const{
148// checks if the track can be a kaon, raw PID applied for single detectors
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//---------------------------
159Bool_t AliAODPidHF::IsPionRaw (AliAODTrack *track, TString detector) const{
160// checks if the track can be a pion, raw PID applied for single detectors
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//---------------------------
172Bool_t AliAODPidHF::IsProtonRaw (AliAODTrack *track, TString detector) const{
173// checks if the track can be a proton raw PID applied for single detectors
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}
184//--------------------------
185Bool_t AliAODPidHF::IsElectronRaw(AliAODTrack *track, TString detector) const{
186// checks if the track can be an electron raw PID applied for single detectors
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//--------------------------
198Int_t AliAODPidHF::ApplyPidTPCRaw(AliAODTrack *track,Int_t specie) const{
199// n-sigma cut, TPC PID
200
201 if(!CheckStatus(track,"TPC")) return 0;
202 Int_t pid=-1;
203 if(fOldPid){
204 AliAODPid *pidObj = track->GetDetPid();
205
206 Double_t dedx=pidObj->GetTPCsignal();
207 Double_t mom = pidObj->GetTPCmomentum();
208 if(mom>fPtThresholdTPC) return 0;
209 AliTPCPIDResponse tpcResponse;
210 SetBetheBloch(tpcResponse);
211 UShort_t nTPCClus=pidObj->GetTPCsignalN();
212 if(nTPCClus==0) {nTPCClus=track->GetTPCNcls();}
213
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
215 Double_t nsigmaMax=fnSigma[0];
216 for(Int_t ipart=0;ipart<5;ipart++){
217 AliPID::EParticleType type=AliPID::EParticleType(ipart);
218 Double_t nsigma = TMath::Abs(tpcResponse.GetNumberOfSigmas(mom,dedx,nTPCClus,type));
219 if((nsigma<nsigmaMax) && (nsigma<fnSigma[0])) {
220 pid=ipart;
221 nsigmaMax=nsigma;
222 }
223 }
224 }else{ // asks only for one particle specie
225 AliPID::EParticleType type=AliPID::EParticleType(specie);
226 Double_t nsigma = TMath::Abs(tpcResponse.GetNumberOfSigmas(mom,dedx,nTPCClus,type));
227 if (nsigma>fnSigma[0]) {
228 pid=-1;
229 }else{
230 pid=specie;
231 }
232 }
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
255
256 return pid;
257
258}
259//----------------------------
260Int_t AliAODPidHF::ApplyPidITSRaw(AliAODTrack *track,Int_t specie) const{
261// truncated mean, ITS PID
262
263 if(!CheckStatus(track,"ITS")) return 0;
264 Int_t pid=-1;
265
266 if(fOldPid){
267 Double_t mom=track->P();
268 AliAODPid *pidObj = track->GetDetPid();
269
270 Double_t dedx=pidObj->GetITSsignal();
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
280 AliITSPIDResponse itsResponse;
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
282 Double_t nsigmaMax=fnSigma[4];
283 for(Int_t ipart=0;ipart<5;ipart++){
284 AliPID::EParticleType type=AliPID::EParticleType(ipart);
285 Double_t nsigma = TMath::Abs(itsResponse.GetNumberOfSigmas(mom,dedx,type,nPointsForPid,isSA));
286 if((nsigma<nsigmaMax) && (nsigma<fnSigma[4])) {
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));
294 if (nsigma>fnSigma[4]) {
295 pid=-1;
296 }else{
297 pid=specie;
298 }
299 }
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
323 return pid;
324}
325//----------------------------
326Int_t AliAODPidHF::ApplyPidTOFRaw(AliAODTrack *track,Int_t specie) const{
327// n-sigma cut, TOF PID
328
329 if(!CheckStatus(track,"TOF")) return 0;
330
331 Double_t time[AliPID::kSPECIESN];
332 Double_t sigmaTOFPid[AliPID::kSPECIES];
333 AliAODPid *pidObj = track->GetDetPid();
334 pidObj->GetIntegratedTimes(time);
335 Double_t sigTOF=pidObj->GetTOFsignal();
336 pidObj->GetTOFpidResolution(sigmaTOFPid);
337
338 Int_t pid=-1;
339
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];
345 for(Int_t ipart=0;ipart<5;ipart++){
346 Double_t nsigma=TMath::Abs(sigTOF-time[ipart]);
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)) {
350 pid=ipart;
351 nsigmaMax=nsigma;
352 }
353 }
354 }else{ // asks only for one particle specie
355 Double_t nsigma=TMath::Abs(sigTOF-time[specie]);
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 }
364 }
365 return pid;
366
367}
368//------------------------------
369void AliAODPidHF::CombinedProbability(AliAODTrack *track,Bool_t *type) const{
370// combined PID stored inside the AOD track
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//--------------------
386void AliAODPidHF::BayesianProbability(AliAODTrack *track,Double_t *pid) const{
387// bayesian PID for single detectors or combined
388
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);
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 }
406
407 return;
408
409}
410//------------------------------------
411void AliAODPidHF::BayesianProbabilityITS(AliAODTrack *track,Double_t *prob) const{
412
413// bayesian PID for ITS
414 AliAODpidUtil pid;
415 Double_t itspid[AliPID::kSPECIES];
416 pid.MakeITSPID(track,itspid);
417 for(Int_t ind=0;ind<AliPID::kSPECIES;ind++){
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 }
423 }
424 return;
425
426}
427//------------------------------------
428void AliAODPidHF::BayesianProbabilityTPC(AliAODTrack *track,Double_t *prob) const{
429// bayesian PID for TPC
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++){
435 if(fTOF || fITS || fTRD){
436 prob[ind]=tpcpid[ind];
437 }else{
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]);
439 }
440}
441 return;
442
443}
444//------------------------------------
445void AliAODPidHF::BayesianProbabilityTOF(AliAODTrack *track,Double_t *prob) const{
446// bayesian PID for TOF
447
448 AliAODpidUtil pid;
449 Double_t tofpid[AliPID::kSPECIES];
450 pid.MakeTOFPID(track,tofpid);
451 for(Int_t ind=0;ind<AliPID::kSPECIES;ind++){
452 if(fTPC || fITS || fTRD){
453 prob[ind]=tofpid[ind];
454 }else{
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 }
457}
458 return;
459
460}
461//---------------------------------
462void AliAODPidHF::BayesianProbabilityTRD(AliAODTrack *track,Double_t *prob) const{
463// bayesian PID for TRD
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++){
469 if(fTPC || fITS || fTOF){
470 prob[ind]=trdpid[ind];
471 }else{
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]);
473 }
474}
475 return;
476
477 }
478//--------------------------------
479Bool_t AliAODPidHF::CheckStatus(AliAODTrack *track,TString detectors) const{
480
481// Quality cuts on the tracks, detector by detector
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 }
490 if(nPointsForPid<3) return kFALSE;
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;
503 if (!(track->GetStatus()&AliESDtrack::kTOFmismatch)==0) return kFALSE;
504 }
505
506
507 if(detectors.Contains("TRD")){
508 if ((track->GetStatus()&AliESDtrack::kTRDout )==0) return kFALSE;
509 UChar_t ntracklets = track->GetTRDntrackletsPID();
510 if(ntracklets<4) return kFALSE;
511 }
512
513 return kTRUE;
514}
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();
522 if(mom>fPtThresholdTPC) return 0;
523 Double_t nsigma=999.;
524 if(fOldPid){
525 Double_t dedx=pidObj->GetTPCsignal();
526 UShort_t nTPCClus=pidObj->GetTPCsignalN();
527 if(nTPCClus==0) {nTPCClus=track->GetTPCNcls();}
528
529 AliTPCPIDResponse tpcResponse;
530 SetBetheBloch(tpcResponse);
531 AliPID::EParticleType type=AliPID::EParticleType(specie);
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
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.)
550 // convention:
551 // for the single detectors: -1 = kFALSE, 1 = kTRUE, 0 = compatible
552 // the method returns the sum of the response of the 2 detectors
553 if(fTPC && fTOF) {if(!CheckStatus(track,"TPC") && !CheckStatus(track,"TOF")) return 0;}
554
555
556 Int_t tTPCinfo=0;
557 if(fTPC){
558 if(CheckStatus(track,"TPC")) {
559 if(fAsym) {
560 if(TPCRawAsym(track,specie)) {
561 tTPCinfo=1;
562 }else{
563 tTPCinfo=-1;
564 }
565 }else{
566 if(specie==2 && IsPionRaw(track,"TPC")) {
567 tTPCinfo=1;
568 }else{
569 tTPCinfo=-1;
570 }
571 if(specie==3 && IsKaonRaw(track,"TPC")) {
572 tTPCinfo=1;
573 }else{
574 tTPCinfo=-1;
575 }
576 if(specie==4 && IsProtonRaw(track,"TPC")) {
577 tTPCinfo=1;
578 }else{
579 tTPCinfo=-1;
580 }
581
582 }
583
584
585 if(compat && tTPCinfo<0){
586 Double_t sig0tmp=fnSigma[0];
587 SetSigma(0,fnSigmaCompat[0]);
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;
591 SetSigma(0,sig0tmp);
592 }
593
594 }
595 }
596
597 Int_t tTOFinfo=0;
598 if(fTOF){
599 if(!CheckStatus(track,"TOF") && fTPC) return tTPCinfo;
600
601 tTOFinfo=-1;
602
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;
606
607 if(compat && tTOFinfo>0){
608 Double_t ptrack=track->P();
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 }
617 }
618 }
619
620
621 if(tTPCinfo+tTOFinfo==0 && fTOFdecide){
622 if(!CheckStatus(track,"TOF")) return tTPCinfo;
623 return tTOFinfo;
624 }
625
626 if(tTPCinfo+tTOFinfo==0 && fITS){
627 if(!CheckStatus(track,"ITS")) return tTPCinfo+tTOFinfo;
628 Int_t tITSinfo = -1;
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;
632 return tITSinfo;
633 }
634
635 return tTPCinfo+tTOFinfo;
636}
637 if(mode==2){
638 //TPC & TOF (a la' Yifei)
639 // convention: -1 = kFALSE, 1 = kTRUE, 0 = not identified
640 Int_t tTPCinfo=0;
641
642 if(fTPC && CheckStatus(track,"TPC")) {
643 tTPCinfo=1;
644 if(fAsym){
645 if(!TPCRawAsym(track,specie)) tTPCinfo=-1;
646 }else{
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;
650 }
651 }
652
653 Int_t tTOFinfo=1;
654 if(fTOF){
655 if(fTPC && !CheckStatus(track,"TOF")) return tTPCinfo;
656
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;
660 }
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;
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;
670 return tITSinfo;
671 }
672
673 return -1;
674
675}
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
681 if(fTPC && fTOF) if(!CheckStatus(track,"TPC") && !CheckStatus(track,"TOF")) return 0;
682
683 Double_t ptrack=track->P();
684
685
686 Int_t tTPCinfo=-1;
687
688 if(ptrack>=fPLimit[0] && ptrack<fPLimit[1] && fTPC) {
689 if(!CheckStatus(track,"TPC")) return 0;
690 if(fAsym) {
691 if(TPCRawAsym(track,specie)) tTPCinfo=1;
692 }else{
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;
696 }
697 return tTPCinfo;
698 }
699
700 Int_t tTOFinfo=-1;
701 if(ptrack>=fPLimit[1] && fTOF){
702 if(!CheckStatus(track,"TOF")) return 0;
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;
707 }
708
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
718 }
719
720 return -1;
721
722}
723//----------------------------------
724Int_t AliAODPidHF::MakeRawPid(AliAODTrack *track, Int_t specie){
725// general method to compute PID
726 if(fMatch>0){
727 return MatchTPCTOF(track,fMatch,specie,fCompat);
728 }else{
729 if(fTPC && !fTOF && !fITS) {
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 }
740 if(tTPCres==specie){return 1;}else{return tTPCres;};
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) {
746 Int_t tTOFres=ApplyPidTOFRaw(track,specie);
747 if(tTOFres==specie){return 1;}else{return tTOFres;};
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) {
754 Int_t tITSres=ApplyPidITSRaw(track,specie);
755 if(tITSres==specie){return 1;}else{return tITSres;};
756 }else{
757 AliError("You should enable just one detector if you don't want to match");
758 return 0;
759 }
760 }
761
762}
763//--------------------------------------------
764void AliAODPidHF::SetBetheBloch(AliTPCPIDResponse &tpcResp) const{
765
766 Double_t alephParameters[5];
767
768 if(fMC) { // MC
769
770 if(fPbPb) { // PbPb MC
771
772 alephParameters[0] = 1.44405/50.;
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;
777
778 } else { // pp MC
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;
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;
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 }
798 }
799
800 } else { // Real Data
801
802 if(fOnePad) { // pp 1-pad (since LHC10d)
803
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 }
841
842 }
843
844 tpcResp.SetBetheBlochParameters(alephParameters[0],alephParameters[1],alephParameters[2],alephParameters[3],alephParameters[4]);
845
846
847 return;
848
849}
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}
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}
890