Update master to aliroot
[u/mrichter/AliRoot.git] / STEER / AOD / AliAODpidUtil.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id: AliAODpidUtil.cxx 38329 2010-01-17 19:17:24Z hristov $ */
17
18 //-----------------------------------------------------------------
19 //           Implementation of the combined PID class
20 //           For the AOD Class
21 //           containing information on the particle identification
22 //      Origin: Rosa Romita, GSI, r.romita@gsi.de
23 //-----------------------------------------------------------------
24
25 #include "TRandom.h"
26 #include "AliLog.h"
27 #include "AliPID.h"
28 #include "AliAODpidUtil.h"
29 #include "AliAODEvent.h"
30 #include "AliAODTrack.h"
31 #include "AliAODPid.h"
32 #include "AliTRDPIDResponse.h"
33 #include "AliESDtrack.h"
34 #include "AliAODMCHeader.h"
35 #include "AliAODMCParticle.h"
36 #include "AliTOFPIDParams.h"
37
38 #include <AliDetectorPID.h>
39
40 ClassImp(AliAODpidUtil)
41
42 //_________________________________________________________________________
43 Float_t AliAODpidUtil::GetTPCsignalTunedOnData(const AliVTrack *t) const {
44     AliAODTrack *track = (AliAODTrack *) t;
45     Float_t dedx = track->GetTPCsignalTunedOnData();
46     if(dedx > 0) return dedx;
47
48     dedx = t->GetTPCsignal();
49     track->SetTPCsignalTunedOnData(dedx);
50
51     if(dedx < 20) return dedx;
52
53     
54     AliPID::EParticleType type = AliPID::kPion;
55     
56     AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(track->GetAODEvent()->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
57     if (mcHeader) {
58         
59         TClonesArray *mcArray = (TClonesArray*)track->GetAODEvent()->GetList()->FindObject(AliAODMCParticle::StdBranchName());
60         
61         Bool_t kGood = kTRUE;
62
63         if ( mcArray->At(TMath::Abs(t->GetLabel())) != NULL ) {  // protects against label-0 tracks e.g. the initial proton for Phythia events
64           Int_t iS = TMath::Abs(((AliAODMCParticle*)mcArray->At(TMath::Abs(t->GetLabel())))->GetPdgCode());
65           if(iS==AliPID::ParticleCode(AliPID::kElectron)){
66             type = AliPID::kElectron;
67           }
68           else if(iS==AliPID::ParticleCode(AliPID::kMuon)){
69             type = AliPID::kMuon;
70           }
71           else if(iS==AliPID::ParticleCode(AliPID::kPion)){
72             type = AliPID::kPion;
73           }
74           else if(iS==AliPID::ParticleCode(AliPID::kKaon)){
75             type = AliPID::kKaon;
76           }
77           else if(iS==AliPID::ParticleCode(AliPID::kProton)){
78             type = AliPID::kProton;
79           }
80           else if(iS==AliPID::ParticleCode(AliPID::kDeuteron)){ // d
81             type = AliPID::kDeuteron;
82           }
83           else if(iS==AliPID::ParticleCode(AliPID::kTriton)){ // t
84             type = AliPID::kTriton;
85           }
86           else if(iS==AliPID::ParticleCode(AliPID::kHe3)){ // 3He
87             type = AliPID::kHe3;
88           }
89           else if(iS==AliPID::ParticleCode(AliPID::kAlpha)){ // 4He
90             type = AliPID::kAlpha;
91           }
92         else
93             kGood = kFALSE;
94         } else kGood = kFALSE;
95
96         if(kGood){
97             //TODO maybe introduce different dEdxSources?
98         Double_t bethe = fTPCResponse.GetExpectedSignal(track, type, AliTPCPIDResponse::kdEdxDefault, this->UseTPCEtaCorrection(),
99                                                         this->UseTPCMultiplicityCorrection());
100         Double_t sigma = fTPCResponse.GetExpectedSigma(track, type, AliTPCPIDResponse::kdEdxDefault, this->UseTPCEtaCorrection(),
101                                                        this->UseTPCMultiplicityCorrection());
102         dedx = gRandom->Gaus(bethe,sigma);
103         
104 //          if(iS == AliPID::ParticleCode(AliPID::kHe3) || iS == AliPID::ParticleCode(AliPID::kAlpha)) dedx *= 5;
105         }
106
107     }
108
109     track->SetTPCsignalTunedOnData(dedx);
110     return dedx;
111 }
112 //_________________________________________________________________________
113 Float_t AliAODpidUtil::GetTOFsignalTunedOnData(const AliVTrack *t) const {
114     AliAODTrack *track = (AliAODTrack *) t;
115     Double_t tofSignal = track->GetTOFsignalTunedOnData();
116
117     if(tofSignal <  99999) return (Float_t)tofSignal; // it has been already set
118
119     // read additional mismatch fraction
120     Float_t addmism = GetTOFPIDParams()->GetTOFadditionalMismForMC();
121     if(addmism > 1.){
122       Float_t centr = GetCurrentCentrality();
123       if(centr > 50) addmism *= 0.1667;
124       else if(centr > 20) addmism *= 0.33;
125     }
126
127     AliAODPid *pidObj = track->GetDetPid();
128     tofSignal = pidObj->GetTOFsignal() + fTOFResponse.GetTailRandomValue(t->Pt(),t->Eta(),pidObj->GetTOFsignal(),addmism);
129     track->SetTOFsignalTunedOnData(tofSignal);
130     return (Float_t)tofSignal;
131 }
132
133 //_________________________________________________________________________
134 Float_t AliAODpidUtil::GetSignalDeltaTOFold(const AliVParticle *vtrack, AliPID::EParticleType type, Bool_t ratio/*=kFALSE*/) const
135 {
136   //
137   // Number of sigma implementation for the TOF
138   //
139   
140   AliAODTrack *track=(AliAODTrack*)vtrack;
141   AliAODPid *pidObj = track->GetDetPid();
142   if (!pidObj) return -9999.;
143   Double_t tofTime = 99999;
144   if (fTuneMConData && ((fTuneMConDataMask & kDetTOF) == kDetTOF) ) tofTime = (Double_t)this->GetTOFsignalTunedOnData((AliVTrack*)vtrack);
145   else tofTime=pidObj->GetTOFsignal();
146   const Double_t expTime=fTOFResponse.GetExpectedSignal((AliVTrack*)vtrack,type);
147   Double_t sigmaTOFPid[AliPID::kSPECIES];
148   pidObj->GetTOFpidResolution(sigmaTOFPid);
149   AliAODEvent *event=(AliAODEvent*)track->GetAODEvent();
150   if (event) {  // protection if the user didn't call GetTrack, which sets the internal pointer
151     AliTOFHeader* tofH=(AliTOFHeader*)event->GetTOFHeader();
152     if (tofH && (TMath::Abs(sigmaTOFPid[0]) <= 1.E-16) ) { // new AOD
153         tofTime -= fTOFResponse.GetStartTime(vtrack->P());
154     }
155   } else {
156     AliError("pointer to AliAODEvent not found, please call GetTrack to set it");
157     return -9999.;
158   }
159
160   Double_t delta=-9999.;
161
162   if (!ratio) delta=tofTime-expTime;
163   else if (expTime>1.e-20) delta=tofTime/expTime;
164   
165   return delta;
166 }
167
168 //_________________________________________________________________________
169 Float_t AliAODpidUtil::GetNumberOfSigmasTOFold(const AliVParticle *vtrack, AliPID::EParticleType type) const
170 {
171   //
172   // Number of sigma implementation for the TOF
173   //
174   
175   AliAODTrack *track=(AliAODTrack*)vtrack;
176
177   Bool_t oldAod=kTRUE;
178   Double_t sigTOF=0.;
179   AliAODPid *pidObj = track->GetDetPid();
180   if (!pidObj) return -999.;
181   Double_t tofTime = 99999;
182   if (fTuneMConData && ((fTuneMConDataMask & kDetTOF) == kDetTOF) ) tofTime = (Double_t)this->GetTOFsignalTunedOnData((AliVTrack*)vtrack);
183   else tofTime=pidObj->GetTOFsignal();
184   Double_t expTime=fTOFResponse.GetExpectedSignal((AliVTrack*)vtrack,type);
185   Double_t sigmaTOFPid[AliPID::kSPECIES];
186   pidObj->GetTOFpidResolution(sigmaTOFPid);
187   AliAODEvent *event=(AliAODEvent*)track->GetAODEvent();
188   if (event) {  // protection if the user didn't call GetTrack, which sets the internal pointer
189     AliTOFHeader* tofH=(AliTOFHeader*)event->GetTOFHeader();
190     if (tofH && (TMath::Abs(sigmaTOFPid[0]) <= 1.E-16) ) { // new AOD
191         sigTOF=fTOFResponse.GetExpectedSigma(track->P(),expTime,AliPID::ParticleMassZ(type)); //fTOFResponse is set in InitialiseEvent
192         tofTime -= fTOFResponse.GetStartTime(vtrack->P());
193         oldAod=kFALSE;
194     }
195   } else {
196     AliError("pointer to AliAODEvent not found, please call GetTrack to set it");
197     return -996.;
198   }
199   if (oldAod) { // old AOD
200     if (type <= AliPID::kProton) {
201       sigTOF=sigmaTOFPid[type];
202     } else return -998.;  // light nuclei cannot be supported on old AOD because we don't have timeZero resolution
203   }
204   if (sigTOF>0) return (tofTime - expTime)/sigTOF;
205   else return -997.;
206 }