new interface library added to be independent of AliRoot dependecies
[u/mrichter/AliRoot.git] / TRD / AliTRDpidESD.cxx
CommitLineData
b0f03c34 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
b52cb1b6 16////////////////////////////////////////////////////////////////////////////
17// //
18// Implementation of the TRD PID class //
19// //
20// Assigns the electron and pion likelihoods to each ESD track. //
21// The function MakePID(AliESD *event) calculates the probability //
22// of having dedx and a maximum timbin at a given //
23// momentum (mom) and particle type k //
24// from the precalculated distributions. //
25// //
26// Original version: //
27// Prashant Shukla <shukla@pi0.physi.uni-heidelberg.de> //
28// //
29////////////////////////////////////////////////////////////////////////////
b0f03c34 30
b52cb1b6 31#include "AliLog.h"
b0f03c34 32#include "AliESD.h"
33#include "AliESDtrack.h"
b52cb1b6 34
35#include "AliTRDpidESD.h"
36#include "AliTRDgeometry.h"
cc7cef99 37#include "AliTRDcalibDB.h"
7754cd1f 38#include "Cal/AliTRDCalPIDLQ.h"
b0f03c34 39
40ClassImp(AliTRDpidESD)
41
df34b21e 42 Bool_t AliTRDpidESD::fCheckTrackStatus = kTRUE;
43 Bool_t AliTRDpidESD::fCheckKinkStatus = kFALSE;
85815482 44 Int_t AliTRDpidESD::fMinPlane = 0;
df34b21e 45
b52cb1b6 46//_____________________________________________________________________________
47AliTRDpidESD::AliTRDpidESD():TObject()
b0f03c34 48{
49 //
b52cb1b6 50 // Default constructor
b0f03c34 51 //
b52cb1b6 52
b0f03c34 53}
54
b52cb1b6 55//_____________________________________________________________________________
56AliTRDpidESD::AliTRDpidESD(const AliTRDpidESD &p):TObject(p)
eab5961e 57{
b0f03c34 58 //
b52cb1b6 59 // AliTRDpidESD copy constructor
b0f03c34 60 //
eab5961e 61
b52cb1b6 62 ((AliTRDpidESD &) p).Copy(*this);
63
64}
65
66//_____________________________________________________________________________
67AliTRDpidESD &AliTRDpidESD::operator=(const AliTRDpidESD &p)
68{
69 //
70 // Assignment operator
71 //
72
73 if (this != &p) ((AliTRDpidESD &) p).Copy(*this);
74 return *this;
75
76}
77
78//_____________________________________________________________________________
79void AliTRDpidESD::Copy(TObject &p) const
80{
81 //
82 // Copy function
83 //
84
25ca55ce 85 ((AliTRDpidESD &) p).fCheckTrackStatus = fCheckTrackStatus;
86 ((AliTRDpidESD &) p).fCheckKinkStatus = fCheckKinkStatus;
87 ((AliTRDpidESD &) p).fMinPlane = fMinPlane;
eab5961e 88
b0f03c34 89}
90
b52cb1b6 91//_____________________________________________________________________________
b0f03c34 92Int_t AliTRDpidESD::MakePID(AliESD *event)
93{
94 //
b52cb1b6 95 // This function calculates the PID probabilities based on TRD signals
b0f03c34 96 //
b52cb1b6 97 // So far this method produces probabilities based on the total charge
98 // in each layer and the position of the maximum time bin in each layer.
99 // In a final version this should also exploit the charge measurement in
100 // the different slices of a given layer.
101 //
102
103 Double_t p[10];
104 Int_t nSpecies = AliPID::kSPECIES;
105 Int_t nPlanePID = 0;
106 Double_t mom = 0.0;
107 Double_t probTotal = 0.0;
108
109 AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
110 if (!calibration) {
030b4415 111 AliErrorGeneral("AliTRDpidESD::MakePID"
112 ,"No access to calibration data\n");
cc7cef99 113 return -1;
b52cb1b6 114 }
115
116 // Retrieve the CDB container class with the probability distributions
63a700c6 117 const AliTRDCalPIDLQ *pd = calibration->GetPIDLQObject();
b52cb1b6 118 if (!pd) {
030b4415 119 AliErrorGeneral("AliTRDpidESD::MakePID"
120 ,"No access to AliTRDCalPIDLQ\n");
b52cb1b6 121 return -1;
122 }
123
124 // Loop through all ESD tracks
125 Int_t ntrk = event->GetNumberOfTracks();
126 for (Int_t i = 0; i < ntrk; i++) {
127
128 AliESDtrack *t = event->GetTrack(i);
129
130 // Check the ESD track status
131 if (fCheckTrackStatus) {
132 if (((t->GetStatus() & AliESDtrack::kTRDout ) == 0) &&
133 ((t->GetStatus() & AliESDtrack::kTRDrefit) == 0)) {
134 continue;
135 }
136 }
137
138 // Check for ESD kink tracks
139 if (fCheckKinkStatus) {
140 if (t->GetKinkIndex(0) != 0) {
141 continue;
142 }
143 }
144
145 // Skip tracks that have no TRD signal at all
146 if (t->GetTRDsignal() == 0) {
147 continue;
148 }
149
150 mom = t->GetP();
151 probTotal = 0.0;
152 nPlanePID = 0;
153 for (Int_t iSpecies = 0; iSpecies < nSpecies; iSpecies++) {
154 p[iSpecies] = 1.0;
155 }
156
157 // Check the different detector layers
158 for (Int_t iPlan = 0; iPlan < AliTRDgeometry::kNplan; iPlan++) {
159
160 // Use the total charge in a given plane
161 Double_t dedx = t->GetTRDsignals(iPlan,-1);
162 Int_t timebin = t->GetTRDTimBin(iPlan);
163 if ((dedx > 0.0) &&
164 (timebin > -1.0)) {
165
166 nPlanePID++;
167
168 // Get the probabilities for the different particle species
169 for (Int_t iSpecies = 0; iSpecies < nSpecies; iSpecies++) {
170
171 p[iSpecies] *= pd->GetProbability(iSpecies,mom,dedx);
172 p[iSpecies] *= pd->GetProbabilityT(iSpecies,mom,timebin);
173 p[iSpecies] *= 100.0; // ??????????????
174
b52cb1b6 175 }
176
177 }
178
85815482 179 }
180
181 for (Int_t iSpecies = 0; iSpecies < nSpecies; iSpecies++) {
182 probTotal += p[iSpecies];
183 }
b52cb1b6 184
185 for (Int_t iSpecies = 0; iSpecies < nSpecies; iSpecies++) {
186 if ((probTotal > 0.0) &&
187 (nPlanePID > fMinPlane)) {
188 p[iSpecies] /= probTotal;
189 }
190 else {
191 p[iSpecies] = 1.0;
192 }
193 }
194
b0f03c34 195 t->SetTRDpid(p);
b52cb1b6 196
197 }
198
b0f03c34 199 return 0;
b52cb1b6 200
b0f03c34 201}