Coding rules (Markus)
[u/mrichter/AliRoot.git] / PWG3 / hfe / AliHFEpidTOF.cxx
CommitLineData
809a4336 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**************************************************************************/
50685501 15//
16// Class for TOF PID
17// Implements the abstract base class AliHFEpidBase
18// IsInitialized() does the PID decision
19//
20// Authors:
21// Markus Fasel <M.Fasel@gsi.de>
22// Matus Kalisky <matus.kalisky@cern.ch> (contact)
23//
809a4336 24
25#include <TH2F.h>
26#include <TList.h>
27#include <TMath.h>
28
722347d8 29#include "AliAODTrack.h"
30#include "AliAODMCParticle.h"
809a4336 31#include "AliESDtrack.h"
722347d8 32#include "AliLog.h"
33#include "AliMCParticle.h"
809a4336 34#include "AliPID.h"
35
36#include "AliHFEpidTOF.h"
37#include "AliHFEpidBase.h"
38
39
40ClassImp(AliHFEpidTOF)
41
42//___________________________________________________________________
43AliHFEpidTOF::AliHFEpidTOF(const Char_t *name):
44 AliHFEpidBase(name)
45 , fPID(0x0)
46 , fQAList(0x0)
47{
48 //
49 // Constructor
50 //
51}
52//___________________________________________________________________
53AliHFEpidTOF::AliHFEpidTOF(const AliHFEpidTOF &c):
54 AliHFEpidBase("")
55 , fPID(0x0)
56 , fQAList(0x0)
57{
58 //
59 // Copy operator
60 //
61
62 c.Copy(*this);
63}
64//___________________________________________________________________
65AliHFEpidTOF &AliHFEpidTOF::operator=(const AliHFEpidTOF &ref){
66 //
67 // Assignment operator
68 //
69
70 if(this != &ref){
71 ref.Copy(*this);
72 }
73
74 return *this;
75}
76//___________________________________________________________________
77AliHFEpidTOF::~AliHFEpidTOF(){
78 //
79 // Destructor
80 //
81 if(fPID) delete fPID;
82 if(fQAList){
83 fQAList->Delete();
84 delete fQAList;
85 }
86}
87//___________________________________________________________________
88void AliHFEpidTOF::Copy(TObject &ref) const {
89 //
90 // Performs the copying of the object
91 //
92 AliHFEpidTOF &target = dynamic_cast<AliHFEpidTOF &>(ref);
93
94 target.fPID = fPID;
95 target.fQAList = fQAList;
96
97 AliHFEpidBase::Copy(ref);
98}
99//___________________________________________________________________
100Bool_t AliHFEpidTOF::InitializePID(){
101 //
102 // InitializePID: TOF experts have to implement code here
103 //
104 return kTRUE;
105}
106
107//___________________________________________________________________
722347d8 108Int_t AliHFEpidTOF::IsSelected(AliHFEpidObject *vtrack)
809a4336 109{
110
111 //
112 // as of 22/05/2006 :
113 // returns AliPID based on the ESD TOF PID decision
114 // the ESD PID will be checked and if necessary improved
115 // in the mean time. Best of luck
116 //
117 // returns 10 (== kUnknown)if PID can not be assigned
118 //
722347d8 119 if(vtrack->fAnalysisType == AliHFEpidObject::kESDanalysis){
120 AliESDtrack *esdTrack = dynamic_cast<AliESDtrack *>(vtrack->fRecTrack);
121 if(!esdTrack) return 0;
122 AliMCParticle *mcTrack = dynamic_cast<AliMCParticle *>(vtrack->fMCtrack);
123 return MakePIDesd(esdTrack, mcTrack);
124 } else {
125 AliAODTrack *aodTrack = dynamic_cast<AliAODTrack *>(vtrack->fRecTrack);
126 if(!aodTrack) return 0;
127 AliAODMCParticle *aodmc = dynamic_cast<AliAODMCParticle *>(vtrack->fMCtrack);
128 return MakePIDaod(aodTrack, aodmc);
129 }
130}
809a4336 131
722347d8 132//___________________________________________________________________
133Int_t AliHFEpidTOF::MakePIDesd(AliESDtrack *track, AliMCParticle * /*mcTrack*/){
50685501 134 //
135 // Does particle identification as discribed in IsSelected
136 //
ad75027f 137 Long_t status = 0;
138 status = track->GetStatus();
139
140 if(!(status & AliESDtrack::kTOFout)) return AliPID::kUnknown;
809a4336 141
142 (dynamic_cast<TH1F *>(fQAList->At(kHistTOFpidFlags)))->Fill(0.);
143
ad75027f 144 Double_t tItrackL = track->GetIntegratedLength();
145 Double_t tTOFsignal = track->GetTOFsignal();
146
147 if(tItrackL > 0)
809a4336 148 (dynamic_cast<TH1F *>(fQAList->At(kHistTOFpidFlags)))->Fill(1.);
149
ad75027f 150 if(tTOFsignal > 0)
809a4336 151 (dynamic_cast<TH1F *>(fQAList->At(kHistTOFpidFlags)))->Fill(2.);
152
153
ad75027f 154 if(tItrackL <=0 || tTOFsignal <=0) return AliPID::kUnknown;
809a4336 155
156 (dynamic_cast<TH1F *>(fQAList->At(kHistTOFpidFlags)))->Fill(3.);
ad75027f 157 (dynamic_cast<TH1F *>(fQAList->At(kHistTOFsignal)))->Fill(tTOFsignal/1000.);
158 (dynamic_cast<TH1F *>(fQAList->At(kHistTOFlength)))->Fill(tItrackL);
809a4336 159 // get the TOF pid probabilities
ad75027f 160 Double_t tESDpid[5] = {0., 0., 0., 0., 0.};
75d81601 161 Float_t tTOFpidSum = 0.;
809a4336 162 // find the largest PID probability
ad75027f 163 track->GetTOFpid(tESDpid);
164 Double_t tMAXpid = 0.;
165 Int_t tMAXindex = -1;
809a4336 166 for(Int_t i=0; i<5; ++i){
75d81601 167 tTOFpidSum += tESDpid[i];
ad75027f 168 if(tESDpid[i] > tMAXpid){
169 tMAXpid = tESDpid[i];
170 tMAXindex = i;
809a4336 171 }
172 }
173
ad75027f 174 Double_t p = track->GetOuterParam()->P();
175 Double_t beta = (tItrackL/100.)/(TMath::C()*(tTOFsignal/1e12));
809a4336 176
75d81601 177 if(TMath::Abs(tTOFpidSum - 1) > 0.01) return AliPID::kUnknown;
809a4336 178 else{
179 // should be the same as AliPID flags
180
75d81601 181 (dynamic_cast<TH2F *>(fQAList->At(kHistTOFpid0+tMAXindex)))->Fill(beta, p);
182 (dynamic_cast<TH2F *>(fQAList->At(kHistTOFpidBetavP)))->Fill(beta, p);
ad75027f 183 return tMAXindex;
809a4336 184 }
185}
186
187//___________________________________________________________________
722347d8 188Int_t AliHFEpidTOF::MakePIDaod(AliAODTrack * /*aodTrack*/, AliAODMCParticle * /*mcTrack*/){
189 AliError("AOD PID not yet implemented");
190 return 0;
191}
192
193//___________________________________________________________________
809a4336 194void AliHFEpidTOF::AddQAhistograms(TList *qaList){
195 //
196 // Create QA histograms for TOF PID
197 //
198
199 fQAList = new TList;
200 fQAList->SetName("fTOFqaHistos");
201 fQAList->AddAt(new TH1F("hTOF_flags", "TOF flags;flags (see code for info);counts", 10, -0.25, 4.75), kHistTOFpidFlags);
75d81601 202 fQAList->AddAt(new TH2F("fTOFbeta_v_P_no","beta -v- P; beta;momentum [GeV/c]", 120, 0, 1.2, 200, 0, 20), kHistTOFpidBetavP);
809a4336 203 fQAList->AddAt(new TH1F("hTOF_signal", "TOF signal; TOF signal [ns];counts", 1000, 12, 50), kHistTOFsignal);
204 fQAList->AddAt(new TH1F("hTOF_length", "TOF track length; length [cm];counts", 400, 300, 700), kHistTOFlength);
75d81601 205 fQAList->AddAt(new TH2F("hTOFpid_electron", "TOF reco electron; beta ; momentum [GeV/c]", 120, 0, 1.2, 200, 0, 5), kHistTOFpid0);
206 fQAList->AddAt(new TH2F("hTOFpid_muon", "TOF reco muon; beta ; momentum [GeV/c]", 120, 0, 1.2, 200, 0, 5), kHistTOFpid1);
207 fQAList->AddAt(new TH2F("hTOFpid_pion", "TOF reco pion; beta ; momentum [GeV/c]", 120, 0, 1.2, 200, 0, 5), kHistTOFpid2);
208 fQAList->AddAt(new TH2F("hTOFpid_kaon", "TOF reco kaon; beta ; momentum [GeV/c]", 120, 0, 1.2, 200, 0, 5), kHistTOFpid3);
209 fQAList->AddAt(new TH2F("hTOFpid_proton", "TOF reco proton; beta ; momentum [GeV/c]", 120, 0, 1.2, 200, 0, 5), kHistTOFpid4);
809a4336 210
211 qaList->AddLast(fQAList);
212}