1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 /************************************************************************
18 * Implements the abstract base class AliHFEpidBase *
20 * Class contains TPC specific cuts and QA histograms *
21 * Two selection strategies are offered: Selection of certain value *
22 * regions in the TPC dE/dx (by IsSelected), and likelihoods *
26 * Markus Fasel <M.Fasel@gsi.de> *
27 * Markus Heide <mheide@uni-muenster.de> *
30 ************************************************************************/
35 #include "AliESDtrack.h"
36 #include "AliExternalTrackParam.h"
39 #include "AliTPCpidESD.h"
40 #include "AliVParticle.h"
42 #include "AliHFEpidTPC.h"
46 //___________________________________________________________________
47 AliHFEpidTPC::AliHFEpidTPC(const char* name) :
50 , fLineCrossingsEnabled(0)
57 // default constructor
59 memset(fLineCrossingCenter, 0, sizeof(Double_t) * AliPID::kSPECIES);
60 memset(fLineCrossingSigma, 0, sizeof(Double_t) * AliPID::kSPECIES);
62 fPIDtpcESD = new AliTPCpidESD;
65 //___________________________________________________________________
66 AliHFEpidTPC::AliHFEpidTPC(const AliHFEpidTPC &ref) :
68 , fLineCrossingsEnabled(0)
80 //___________________________________________________________________
81 AliHFEpidTPC &AliHFEpidTPC::operator=(const AliHFEpidTPC &ref){
83 // Assignment operator
91 void AliHFEpidTPC::Copy(TObject &o) const{
94 // called in copy constructor and assigment operator
96 AliHFEpidTPC &target = dynamic_cast<AliHFEpidTPC &>(o);
98 target.fLineCrossingsEnabled = fLineCrossingsEnabled;
99 target.fNsigmaTPC = fNsigmaTPC;
100 target.fPID = new AliPID(*fPID);
101 target.fPIDtpcESD = new AliTPCpidESD(*fPIDtpcESD);
102 target.fQAList = dynamic_cast<TList *>(fQAList->Clone());
104 AliHFEpidBase::Copy(target);
107 //___________________________________________________________________
108 AliHFEpidTPC::~AliHFEpidTPC(){
112 if(fPID) delete fPID;
113 if(fPIDtpcESD) delete fPIDtpcESD;
120 //___________________________________________________________________
121 Bool_t AliHFEpidTPC::InitializePID(){
123 // Add TPC dE/dx Line crossings
125 AddTPCdEdxLineCrossing(AliPID::kKaon, 0.3, 0.018);
126 AddTPCdEdxLineCrossing(AliPID::kProton, 0.9, 0.054);
130 //___________________________________________________________________
131 Int_t AliHFEpidTPC::IsSelected(AliVParticle *track)
134 // For the TPC pid we use the 2-sigma band around the bethe bloch curve
136 // exclusion of the crossing points
138 AliESDtrack *esdTrack = 0x0;
139 if(!(esdTrack = dynamic_cast<AliESDtrack *>(track))) return kFALSE;
140 if(IsQAon()) FillTPChistograms(esdTrack);
141 Double_t TPCsignal = esdTrack->GetTPCsignal();
142 // exclude crossing points:
143 // Determine the bethe values for each particle species
144 Double_t p = esdTrack->GetInnerParam()->P();
145 Bool_t isLineCrossing = kFALSE;
146 for(Int_t ispecies = 0; ispecies < AliPID::kSPECIES; ispecies++){
147 if(!(fLineCrossingsEnabled & 1 << ispecies)) continue;
148 if(TMath::Abs(p - fLineCrossingCenter[ispecies]) < fLineCrossingSigma[ispecies]){
149 // Point in a line crossing region, no PID possible
150 isLineCrossing = kTRUE;
154 if(isLineCrossing) return 0;
155 // Check whether distance from the electron line is smaller than n-sigma
156 Double_t beta = p/fPID->ParticleMass(AliPID::kElectron);
157 if(TMath::Abs(TPCsignal - 50*fPIDtpcESD->Bethe(beta)) < GetTPCsigma(p,0)) return 11;
161 //___________________________________________________________________
162 void AliHFEpidTPC::AddTPCdEdxLineCrossing(Int_t species, Double_t p, Double_t sigma_p){
164 // Add exclusion point for the TPC PID where a dEdx line crosses the electron line
165 // Stores line center and line sigma
167 if(species >= AliPID::kSPECIES){
168 AliError("Species doesn't exist");
171 fLineCrossingsEnabled |= 1 << species;
172 fLineCrossingCenter[species] = p;
173 fLineCrossingSigma[species] = sigma_p;
176 //___________________________________________________________________
177 Double_t AliHFEpidTPC::GetTPCsigma(Double_t p, Int_t species){
179 // return the TPC sigma, momentum dependent
181 if(p < 0.1 || p > 20.) return 0.;
182 Double_t beta = p/fPID->ParticleMass(species);
185 return 50*fPIDtpcESD->Bethe(beta) * 0.06;
188 //___________________________________________________________________
189 Double_t AliHFEpidTPC::Likelihood(const AliESDtrack *track, Int_t species, Float_t rsig)
191 //gives probability for track to come from a certain particle species;
192 //no priors considered -> probabilities for equal abundances of all species!
193 //optional restriction to r-sigma bands of different particle species; default: rsig = 2. (see below)
195 //IMPORTANT: Tracks which are judged to be outliers get negative likelihoods -> unusable for combination with further detectors!
197 if(!track) return -1.;
198 Int_t hypo; //marks particle hypotheses for 2-sigma bands
199 Double_t beta;//well o.k., it corresponds to gamma * beta
200 Double_t p = track->GetInnerParam()->P();
201 Double_t TPCsignal = track->GetTPCsignal();
202 Bool_t outlier = kTRUE;
203 // Check whether distance from the respective particle line is smaller than r sigma
204 for(hypo = 0; hypo < 5; hypo++)
206 beta = p/fPID->ParticleMass(hypo);
207 if(TMath::Abs(TPCsignal - (GetTPCsigma(p, hypo))/0.06) > (rsig * GetTPCsigma(p,hypo)))
220 track->GetTPCpid(TPCprob);
222 return TPCprob[species];
224 //___________________________________________________________________
225 Double_t AliHFEpidTPC::Likelihood(const AliESDtrack *track, Int_t species)
228 // for everything else, see above!
230 if(!track) return -1.;
231 Int_t hypo; //marks particle hypotheses for 2-sigma bands
233 Double_t p = track->GetInnerParam()->P();
234 Double_t TPCsignal = track->GetTPCsignal();
235 Bool_t outlier = kTRUE;
236 // Check whether distance from the respective particle line is smaller than 2 sigma
237 for(hypo = 0; hypo < 5; hypo++)
239 beta = p/fPID->ParticleMass(hypo);
241 if(TMath::Abs(TPCsignal - (GetTPCsigma(p, hypo))/0.06) > (2. * GetTPCsigma(p,hypo)))
254 track->GetTPCpid(TPCprob);
256 return TPCprob[species];
258 //___________________________________________________________________
259 Double_t AliHFEpidTPC::Suppression(const AliESDtrack *track, Int_t species)
261 //ratio of likelihoods to be whatever species/to be an electron;
262 //as a cross-check for possible particle type suppression compared to electrons
263 if(!track) return -20;
264 if((Likelihood(track,species) == -2.)||(Likelihood(track,0)== -2.))
266 if(Likelihood(track,species) == 0.)
268 if (Likelihood(track,0) == 0.)
271 return TMath::Log10(Likelihood(track,species)/(Likelihood(track,0)));
275 //___________________________________________________________________
276 void AliHFEpidTPC::FillTPChistograms(const AliESDtrack *track){
281 Double_t tpc_signal = track->GetTPCsignal();
282 Double_t p = track->GetInnerParam() ? track->GetInnerParam()->P() : track->P();
284 switch(TMath::Abs(GetPdgCode(const_cast<AliESDtrack *>(track)))){
285 case 11: (dynamic_cast<TH2I *>(fQAList->At(kHistTPCelectron)))->Fill(p, tpc_signal);
286 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCprobEl)))->Fill(p, Likelihood(track, 0));
287 //histograms with ratio of likelihood to be electron/to be other species (a check for quality of likelihood PID);
288 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCenhanceElPi)))->Fill(p, -Suppression(track, 2));
289 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCenhanceElMu)))->Fill(p, -Suppression(track, 1));
290 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCenhanceElKa)))->Fill(p, -Suppression(track, 3));
291 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCenhanceElPro)))->Fill(p, -Suppression(track, 4));
292 //___________________________________________________________________________________________
293 //Likelihoods for electrons to be other particle species
294 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCElprobPi)))->Fill(p, Likelihood(track, 2));
295 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCElprobMu)))->Fill(p, Likelihood(track, 1));
296 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCElprobKa)))->Fill(p, Likelihood(track, 3));
297 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCElprobPro)))->Fill(p, Likelihood(track, 4));
299 //___________________________________________________________________________________________
300 case 13: (dynamic_cast<TH2I *>(fQAList->At(kHistTPCmuon)))->Fill(p, tpc_signal);
301 //Likelihood of muon to be an electron
302 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCprobMu)))->Fill(p, Likelihood(track, 0));
303 //ratio of likelihood for muon to be a muon/an electron -> indicator for quality of muon suppression
304 //below functions are the same for other species
305 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCsuppressMu)))->Fill(p, Suppression(track, 1));
307 case 211: (dynamic_cast<TH2I *>(fQAList->At(kHistTPCpion)))->Fill(p, tpc_signal);
308 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCprobPi)))->Fill(p, Likelihood(track, 0));
309 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCsuppressPi)))->Fill(p, Suppression(track, 2));
311 case 321: (dynamic_cast<TH2I *>(fQAList->At(kHistTPCkaon)))->Fill(p, tpc_signal);
312 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCprobKa)))->Fill(p, Likelihood(track, 0));
313 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCsuppressKa)))->Fill(p, Suppression(track, 3));
315 case 2212: (dynamic_cast<TH2I *>(fQAList->At(kHistTPCproton)))->Fill(p, tpc_signal);
316 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCprobPro)))->Fill(p, Likelihood(track, 0));
317 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCsuppressPro)))->Fill(p, Suppression(track, 4));
319 default: (dynamic_cast<TH2I *>(fQAList->At(kHistTPCothers)))->Fill(p, tpc_signal);
320 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCprobOth)))->Fill(p, Likelihood(track, 0));
325 //TPC signal and Likelihood to be electron for all tracks (independent of MC information)
326 (dynamic_cast<TH2I *>(fQAList->At(kHistTPCall)))->Fill(p, tpc_signal);
327 (dynamic_cast<TH2F *>(fQAList->At(kHistTPCprobAll)))->Fill(p, Likelihood(track, 0));
333 //___________________________________________________________________
334 void AliHFEpidTPC::AddQAhistograms(TList *qaList){
336 fQAList->SetName("fTPCqaHistos");
338 fQAList->AddAt(new TH2I("fHistTPCelectron","TPC signal for Electrons", 200, 0, 20, 60, 0, 600), kHistTPCelectron);
339 fQAList->AddAt(new TH2I("fHistTPCmuon","TPC signal for Muons", 200, 0, 20, 60, 0, 600), kHistTPCmuon);
340 fQAList->AddAt(new TH2I("fHistTPCpion","TPC signal for Pions", 200, 0, 20, 60, 0, 600), kHistTPCpion);
341 fQAList->AddAt(new TH2I("fHistTPCkaon","TPC signal for Kaons", 200, 0, 20, 60, 0, 600), kHistTPCkaon);
342 fQAList->AddAt(new TH2I("fHistTPCproton","TPC signal for Protons", 200, 0, 20, 60, 0, 600), kHistTPCproton);
343 fQAList->AddAt(new TH2I("fHistTPCothers","TPC signal for other species", 200, 0, 20, 60, 0, 600), kHistTPCothers);
344 fQAList->AddAt(new TH2I("fHistTPCall","TPC signal for all species", 200, 0, 20, 60, 0, 600), kHistTPCall);
346 fQAList->AddAt(new TH2F("fHistTPCprobEl","TPC likelihood for electrons to be an electron vs. p", 200, 0.,20.,200,0.,1.), kHistTPCprobEl);
347 fQAList->AddAt(new TH2F("fHistTPCprobPi","TPC likelihood for pions to be an electron vs. p", 200, 0.,20.,200, 0.,1.), kHistTPCprobPi);
348 fQAList->AddAt(new TH2F("fHistTPCprobMu","TPC likelihood for muons to be an electron vs. p", 200, 0.,20.,200, 0.,1.), kHistTPCprobMu);
349 fQAList->AddAt(new TH2F("fHistTPCprobKa","TPC likelihood for kaons to be an electron vs. p", 200, 0.,20.,200, 0.,1.), kHistTPCprobKa);
350 fQAList->AddAt(new TH2F("fHistTPCprobPro","TPC likelihood for protons to be an electron vs. p", 200, 0.,20.,200, 0.,1.), kHistTPCprobPro);
351 fQAList->AddAt(new TH2F("fHistTPCprobOth","TPC likelihood for other particles to be an electron vs. p", 200, 0.,20.,200, 0.,1.), kHistTPCprobOth);
352 fQAList->AddAt(new TH2F("fHistTPCprobAll","TPC likelihood for all particles to be an electron vs. p", 200, 0.,20.,200, 0.,1.), kHistTPCprobAll);
355 fQAList->AddAt(new TH2F("fHistTPCsuppressPi","log10 of TPC Likelihood(pion)/Likelihood(elec) for pions vs. p", 200, 0.,20.,200,-1.,5.8), kHistTPCsuppressPi);
356 fQAList->AddAt(new TH2F("fHistTPCsuppressMu","log10 of TPC Likelihood(muon)/Likelihood(elec) for muons vs. p", 200, 0.,20.,200,-1.,5.8), kHistTPCsuppressMu);
357 fQAList->AddAt(new TH2F("fHistTPCsuppressKa","log10 of TPC Likelihood(kaon)/Likelihood(elec) for kaons vs. p", 200, 0.,20.,200,-1.,5.8), kHistTPCsuppressKa);
358 fQAList->AddAt(new TH2F("fHistTPCsuppressPro","log10 of TPC Likelihood(proton)/Likelihood(elec)for protons vs. p", 200, 0.,20.,200,-1.,5.8), kHistTPCsuppressPro);
360 fQAList->AddAt(new TH2F("fHistTPCenhanceElPi","log10 of TPC Likelihood(elec)/Likelihood(pion) for electrons vs. p", 200, 0.,20.,200,-1.,5.8), kHistTPCsuppressPi);
361 fQAList->AddAt(new TH2F("fHistTPCenhanceElMu","log10 of TPC Likelihood(elec)/Likelihood(muon) for electrons vs. p", 200, 0.,20.,200,-1.,5.8), kHistTPCsuppressMu);
362 fQAList->AddAt(new TH2F("fHistTPCenhanceElKa","log10 of TPC Likelihood(elec)/Likelihood(kaon) for electrons vs. p", 200, 0.,20.,200,-1.,5.8), kHistTPCsuppressKa);
363 fQAList->AddAt(new TH2F("fHistTPCenhanceElPro","log10 of TPC Likelihood(elec)/Likelihood(proton) for electrons vs. p", 200, 0.,20.,200,-1.,5.8), kHistTPCsuppressPro);
365 fQAList->AddAt(new TH2F("fHistTPCElprobPi","TPC likelihood for electrons to be a pion vs. p", 200, 0.,20.,200,0.,1.), kHistTPCElprobPi);
366 fQAList->AddAt(new TH2F("fHistTPCElprobMu","TPC likelihood for electrons to be a muon vs. p", 200, 0.,20.,200,0.,1.), kHistTPCElprobMu);
367 fQAList->AddAt(new TH2F("fHistTPCElprobKa","TPC likelihood for electrons to be a kaon vs. p", 200, 0.,20.,200,0.,1.), kHistTPCElprobKa);
368 fQAList->AddAt(new TH2F("fHistTPCElprobPro","TPC likelihood for electrons to be a proton vs. p", 200, 0.,20.,200,0.,1.), kHistTPCElprobPro);
371 qaList->AddLast(fQAList);