Updates in PbPb cuts (Andrea Rossi)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAODPidHF.cxx
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
34 ClassImp(AliAODPidHF)
35
36 //------------------------------
37 AliAODPidHF::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   fPbPb(kFALSE),
60   fTOFdecide(kFALSE),
61   fOldPid(kTRUE),
62   fPtThresholdTPC(999999.),
63   fPidResponse(0),
64   fPidCombined(new AliPIDCombined())
65 {
66  //
67  // Default constructor
68  //
69  fPLimit=new Double_t[fnPLimit];
70  fnSigma=new Double_t[fnNSigma];
71  fPriors=new Double_t[fnPriors];
72  fnSigmaCompat=new Double_t[fnNSigmaCompat];
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  }
83  for(Int_t i=0;i<fnNSigmaCompat;i++){
84   fnSigmaCompat[i]=3.;
85  }
86
87 }
88 //----------------------
89 AliAODPidHF::~AliAODPidHF()
90 {
91       // destructor
92  //   if(fPLimit) delete fPLimit;
93  //   if(fnSigma)  delete fnSigma;
94  //   if(fPriors)  delete fPriors;
95 }
96 //------------------------
97 AliAODPidHF::AliAODPidHF(const AliAODPidHF& pid) :
98   AliAODPid(pid),
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),
106   fAsym(pid.fAsym),
107   fTPC(pid.fTPC),
108   fTOF(pid.fTOF),
109   fITS(pid.fITS),
110   fTRD(pid.fTRD),
111   fMatch(pid.fMatch),
112   fCompat(pid.fCompat),
113   fPCompatTOF(pid.fPCompatTOF),
114   fnNSigmaCompat(pid.fnNSigmaCompat),
115   fnSigmaCompat(pid.fnSigmaCompat),
116   fMC(pid.fMC),
117   fOnePad(pid.fOnePad),
118   fMCLowEn2011(pid.fMCLowEn2011),
119   fPbPb(pid.fPbPb),
120   fTOFdecide(pid.fTOFdecide),
121   fOldPid(pid.fOldPid),
122   fPtThresholdTPC(pid.fPtThresholdTPC),
123   fPidResponse(pid.fPidResponse),
124   fPidCombined(pid.fPidCombined)  
125   {
126   
127   for(Int_t i=0;i<5;i++){
128     fPriors[i]=pid.fPriors[i];
129   }
130   
131   }
132
133 //----------------------
134 Int_t AliAODPidHF::RawSignalPID(AliAODTrack *track, TString detector) const{
135 // raw PID for single detectors, returns the particle type with smaller sigma
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);
140
141   return specie;
142
143 }
144 //---------------------------
145 Bool_t AliAODPidHF::IsKaonRaw(AliAODTrack *track, TString detector) const{
146 // checks if the track can be a kaon, raw PID applied for single detectors
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 //---------------------------
157 Bool_t AliAODPidHF::IsPionRaw (AliAODTrack *track, TString detector) const{
158 // checks if the track can be a pion, raw PID applied for single detectors
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 //---------------------------
170 Bool_t AliAODPidHF::IsProtonRaw (AliAODTrack *track, TString detector) const{
171 // checks if the track can be a proton raw PID applied for single detectors
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 }
182 //--------------------------
183 Bool_t AliAODPidHF::IsElectronRaw(AliAODTrack *track, TString detector) const{
184 // checks if the track can be an electron raw PID applied for single detectors
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 //--------------------------
196 Int_t AliAODPidHF::ApplyPidTPCRaw(AliAODTrack *track,Int_t specie) const{
197 // n-sigma cut, TPC PID
198
199   if(!CheckStatus(track,"TPC")) return 0;
200   Int_t pid=-1;
201   if(fOldPid){
202   AliAODPid *pidObj = track->GetDetPid();
203   
204   Double_t dedx=pidObj->GetTPCsignal();
205   Double_t mom = pidObj->GetTPCmomentum();
206   if(mom>fPtThresholdTPC) return 0;
207   AliTPCPIDResponse tpcResponse;
208   SetBetheBloch(tpcResponse); 
209   UShort_t nTPCClus=pidObj->GetTPCsignalN();
210   if(nTPCClus==0) {nTPCClus=track->GetTPCNcls();}
211
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
213    Double_t nsigmaMax=fnSigma[0];
214    for(Int_t ipart=0;ipart<5;ipart++){
215     AliPID::EParticleType type=AliPID::EParticleType(ipart);
216     Double_t nsigma = TMath::Abs(tpcResponse.GetNumberOfSigmas(mom,dedx,nTPCClus,type));
217     if((nsigma<nsigmaMax) && (nsigma<fnSigma[0])) {
218      pid=ipart;
219      nsigmaMax=nsigma;
220     }
221    }
222   }else{ // asks only for one particle specie
223    AliPID::EParticleType type=AliPID::EParticleType(specie);
224     Double_t nsigma = TMath::Abs(tpcResponse.GetNumberOfSigmas(mom,dedx,nTPCClus,type));
225    if (nsigma>fnSigma[0]) {
226     pid=-1; 
227    }else{
228     pid=specie;
229    }
230   }
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
253
254  return pid;
255
256 }
257 //----------------------------
258 Int_t AliAODPidHF::ApplyPidITSRaw(AliAODTrack *track,Int_t specie) const{
259 // truncated mean, ITS PID
260
261   if(!CheckStatus(track,"ITS")) return 0;
262   Int_t pid=-1;
263
264   if(fOldPid){
265   Double_t mom=track->P();
266   AliAODPid *pidObj = track->GetDetPid();
267
268   Double_t dedx=pidObj->GetITSsignal();
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
278   AliITSPIDResponse itsResponse;
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
280    Double_t nsigmaMax=fnSigma[4];
281    for(Int_t ipart=0;ipart<5;ipart++){
282     AliPID::EParticleType type=AliPID::EParticleType(ipart);
283     Double_t nsigma = TMath::Abs(itsResponse.GetNumberOfSigmas(mom,dedx,type,nPointsForPid,isSA));
284     if((nsigma<nsigmaMax) && (nsigma<fnSigma[4])) {
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));
292    if (nsigma>fnSigma[4]) {
293     pid=-1; 
294    }else{
295     pid=specie;
296    }
297   }
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
321  return pid; 
322 }
323 //----------------------------
324 Int_t AliAODPidHF::ApplyPidTOFRaw(AliAODTrack *track,Int_t specie) const{
325 // n-sigma cut, TOF PID
326
327  if(!CheckStatus(track,"TOF")) return 0;
328
329  Double_t time[AliPID::kSPECIESN];
330  Double_t sigmaTOFPid[AliPID::kSPECIES];
331  AliAODPid *pidObj = track->GetDetPid();
332  pidObj->GetIntegratedTimes(time);
333  Double_t sigTOF=pidObj->GetTOFsignal();
334  pidObj->GetTOFpidResolution(sigmaTOFPid);
335
336  Int_t pid=-1;
337
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];
343    for(Int_t ipart=0;ipart<5;ipart++){
344     Double_t nsigma=TMath::Abs(sigTOF-time[ipart]);
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)) {
348      pid=ipart;
349      nsigmaMax=nsigma;
350     }
351    }
352   }else{ // asks only for one particle specie
353     Double_t nsigma=TMath::Abs(sigTOF-time[specie]);
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     }
362   }
363  return pid; 
364
365 }
366 //------------------------------
367 void AliAODPidHF::CombinedProbability(AliAODTrack *track,Bool_t *type) const{
368 // combined PID stored inside the AOD track
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 //--------------------
384 void AliAODPidHF::BayesianProbability(AliAODTrack *track,Double_t *pid) const{
385 // bayesian PID for single detectors or combined
386
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);
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     }
404
405  return;
406
407 }
408 //------------------------------------
409 void AliAODPidHF::BayesianProbabilityITS(AliAODTrack *track,Double_t *prob) const{
410
411 // bayesian PID for ITS
412  AliAODpidUtil pid;
413  Double_t itspid[AliPID::kSPECIES];
414  pid.MakeITSPID(track,itspid);
415  for(Int_t ind=0;ind<AliPID::kSPECIES;ind++){
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   }
421  }
422  return;
423
424 }
425 //------------------------------------
426 void AliAODPidHF::BayesianProbabilityTPC(AliAODTrack *track,Double_t *prob) const{
427 // bayesian PID for TPC
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++){
433   if(fTOF || fITS || fTRD){
434    prob[ind]=tpcpid[ind];
435   }else{
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]);
437  }
438 }
439  return;
440
441 }
442 //------------------------------------
443 void AliAODPidHF::BayesianProbabilityTOF(AliAODTrack *track,Double_t *prob) const{
444 // bayesian PID for TOF
445
446  AliAODpidUtil pid;
447  Double_t tofpid[AliPID::kSPECIES];
448  pid.MakeTOFPID(track,tofpid);
449  for(Int_t ind=0;ind<AliPID::kSPECIES;ind++){
450   if(fTPC || fITS || fTRD){
451    prob[ind]=tofpid[ind];
452   }else{
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  }
455 }
456  return;
457
458 }
459 //---------------------------------
460 void AliAODPidHF::BayesianProbabilityTRD(AliAODTrack *track,Double_t *prob) const{
461 // bayesian PID for TRD
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++){
467   if(fTPC || fITS || fTOF){
468    prob[ind]=trdpid[ind];
469   }else{
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]);
471  }
472 }
473   return;
474
475  }
476 //--------------------------------
477 Bool_t AliAODPidHF::CheckStatus(AliAODTrack *track,TString detectors) const{
478
479 // Quality cuts on the tracks, detector by detector
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   }
488   if(nPointsForPid<3) return kFALSE;
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;
501    if (!(track->GetStatus()&AliESDtrack::kTOFmismatch)==0)    return kFALSE;
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 }
518 //--------------------------------------------
519 Bool_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();
525   if(mom>fPtThresholdTPC) return 0;
526   Double_t nsigma=999.;
527   if(fOldPid){
528   Double_t dedx=pidObj->GetTPCsignal();
529   UShort_t nTPCClus=pidObj->GetTPCsignalN();
530   if(nTPCClus==0) {nTPCClus=track->GetTPCNcls();}
531   
532   AliTPCPIDResponse tpcResponse;
533   SetBetheBloch(tpcResponse); 
534   AliPID::EParticleType type=AliPID::EParticleType(specie);
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
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 //------------------
549 Int_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.)
553  // convention: 
554  // for the single detectors: -1 = kFALSE, 1 = kTRUE, 0 = compatible
555  // the method returns the sum of the response of the 2 detectors
556   if(fTPC && fTOF) {if(!CheckStatus(track,"TPC") && !CheckStatus(track,"TOF")) return 0;}
557
558   
559   Int_t tTPCinfo=0;
560   if(fTPC){
561   if(CheckStatus(track,"TPC")) {
562    if(fAsym) {
563     if(TPCRawAsym(track,specie)) {
564       tTPCinfo=1;
565      }else{
566       tTPCinfo=-1;
567      }
568    }else{
569     if(specie==2 && IsPionRaw(track,"TPC")) {
570      tTPCinfo=1;
571     }else{
572      tTPCinfo=-1;
573     }
574     if(specie==3 && IsKaonRaw(track,"TPC")) {
575      tTPCinfo=1;
576     }else{
577      tTPCinfo=-1;
578     }
579     if(specie==4 && IsProtonRaw(track,"TPC")) {
580      tTPCinfo=1;
581     }else{
582      tTPCinfo=-1;
583     }
584
585    }
586
587
588    if(compat && tTPCinfo<0){
589     Double_t sig0tmp=fnSigma[0];
590     SetSigma(0,fnSigmaCompat[0]);
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;
594     SetSigma(0,sig0tmp);
595    }
596
597   }
598  }
599
600  Int_t tTOFinfo=0;
601  if(fTOF){
602   if(!CheckStatus(track,"TOF") && fTPC) return tTPCinfo;
603
604   tTOFinfo=-1;
605  
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;
609
610   if(compat && tTOFinfo>0){
611    Double_t ptrack=track->P();
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    }
620   }
621  }
622
623  
624  if(tTPCinfo+tTOFinfo==0 && fTOFdecide){
625   if(!CheckStatus(track,"TOF")) return tTPCinfo;
626   return tTOFinfo;
627  }
628
629  if(tTPCinfo+tTOFinfo==0 && fITS){
630   if(!CheckStatus(track,"ITS")) return tTPCinfo+tTOFinfo;
631   Int_t tITSinfo = -1;
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;
635   return tITSinfo;
636  }
637
638  return tTPCinfo+tTOFinfo;
639 }
640  if(mode==2){
641   //TPC & TOF (a la' Yifei)
642  // convention: -1 = kFALSE, 1 = kTRUE, 0 = not identified
643   Int_t tTPCinfo=0; 
644   
645   if(fTPC && CheckStatus(track,"TPC")) {
646    tTPCinfo=1;
647    if(fAsym){
648     if(!TPCRawAsym(track,specie)) tTPCinfo=-1;
649    }else{
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;
653    }
654   }
655
656   Int_t tTOFinfo=1;
657   if(fTOF){
658    if(fTPC && !CheckStatus(track,"TOF")) return tTPCinfo;
659
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;
663   }
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;
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;
673   return tITSinfo;
674  }
675
676    return -1;
677
678 }
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   
684   if(fTPC && fTOF) if(!CheckStatus(track,"TPC") && !CheckStatus(track,"TOF")) return 0;
685
686   Double_t ptrack=track->P();
687   
688
689   Int_t tTPCinfo=-1;
690
691    if(ptrack>=fPLimit[0] && ptrack<fPLimit[1] && fTPC) {  
692     if(!CheckStatus(track,"TPC")) return 0;
693     if(fAsym) {
694      if(TPCRawAsym(track,specie)) tTPCinfo=1;
695     }else{
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;
699     } 
700     return tTPCinfo;
701    }
702
703    Int_t tTOFinfo=-1;
704    if(ptrack>=fPLimit[1] && fTOF){
705     if(!CheckStatus(track,"TOF")) return 0;
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;
710    }
711
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
721  }
722
723  return -1;
724
725 }
726 //----------------------------------   
727 Int_t AliAODPidHF::MakeRawPid(AliAODTrack *track, Int_t specie){
728 // general method to compute PID
729  if(fMatch>0){
730   return MatchTPCTOF(track,fMatch,specie,fCompat); 
731  }else{
732   if(fTPC && !fTOF && !fITS) {
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    }
743    if(tTPCres==specie){return 1;}else{return tTPCres;};
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) {
749    Int_t tTOFres=ApplyPidTOFRaw(track,specie); 
750    if(tTOFres==specie){return 1;}else{return tTOFres;};
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) {
757    Int_t tITSres=ApplyPidITSRaw(track,specie);
758    if(tITSres==specie){return 1;}else{return tITSres;};
759   }else{
760    AliError("You should enable just one detector if you don't want to match");
761    return 0;
762   }
763  } 
764   
765 }
766 //--------------------------------------------
767 void AliAODPidHF::SetBetheBloch(AliTPCPIDResponse &tpcResp) const{
768
769  Double_t alephParameters[5];
770
771  if(fMC) {  // MC
772
773    if(fPbPb) { // PbPb MC
774
775      alephParameters[0] = 1.44405/50.;
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;
780
781    } else {  // pp MC
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;
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;
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      }
801    }
802
803  } else { // Real Data
804
805    if(fOnePad) { // pp 1-pad (since LHC10d)
806
807     alephParameters[0] =1.34490e+00/50.; 
808     alephParameters[1] = 2.69455e+01; 
809     alephParameters[2] = TMath::Exp(-2.97552e+01); 
810     alephParameters[3] = 2.35339e+00; 
811     alephParameters[4] = 5.98079e+00;
812     
813   } else if(fPbPb) { // PbPb 
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     
827   } else {  // pp no 1-pad (LHC10bc)
828
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;
834     
835   }
836   
837  }
838
839  tpcResp.SetBetheBlochParameters(alephParameters[0],alephParameters[1],alephParameters[2],alephParameters[3],alephParameters[4]);
840
841
842  return;
843
844 }
845 //-----------------------
846 Bool_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 }
866 //--------------------------------------------------------------------------
867 void 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 //--------------------------------------------------------------------------
877 void AliAODPidHF::DrawPrior(AliPID::EParticleType type){
878
879         //
880         // Drawing prior distribution for type "type"
881
882         new TCanvas();
883         GetPidCombined()->GetPriorDistribution(type)->Draw();
884 }
885