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