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 **************************************************************************/
18 #include "AliESDtrack.h"
20 #include "AliTRDtrackV1.h"
21 #include "AliTRDcluster.h"
22 #include "AliTRDcalibDB.h"
23 #include "AliTRDrecoParam.h"
25 ClassImp(AliTRDtrackV1)
27 ///////////////////////////////////////////////////////////////////////////////
29 // Represents a reconstructed TRD track //
30 // Local TRD Kalman track //
33 // Alex Bercuci <A.Bercuci@gsi.de> //
34 // Markus Fasel <M.Fasel@gsi.de> //
36 ///////////////////////////////////////////////////////////////////////////////
38 //_______________________________________________________________
39 AliTRDtrackV1::AliTRDtrackV1()
43 // Default constructor
46 for(int ip=0; ip<6; ip++){
47 fTrackletIndex[ip] = -1;
48 fTracklet[ip].Reset();
52 //_______________________________________________________________
53 AliTRDtrackV1::AliTRDtrackV1(const AliESDtrack &t)
57 // Constructor from AliESDtrack
60 //AliInfo(Form("alpha %f", GetAlpha()));
61 t.GetTRDtracklets(&fTrackletIndex[0]);
62 for(int ip=0; ip<6; ip++) fTracklet[ip].Reset();
65 //_______________________________________________________________
66 AliTRDtrackV1::AliTRDtrackV1(const AliTRDtrackV1 &ref)
73 for(int ip=0; ip<6; ip++){
74 fTrackletIndex[ip] = ref.fTrackletIndex[ip];
75 fTracklet[ip] = ref.fTracklet[ip];
79 //_______________________________________________________________
80 // AliTRDtrackV1::~AliTRDtrackV1()
85 //_______________________________________________________________
86 AliTRDtrackV1::AliTRDtrackV1(AliTRDseedV1 *trklts, const Double_t p[5], const Double_t cov[15]
87 , Double_t x, Double_t alpha)
91 // The stand alone tracking constructor
92 // TEMPORARY !!!!!!!!!!!
94 // 1. covariance matrix
95 // 2. dQdl calculation
98 Double_t cnv = 1.0 / (GetBz() * kB2C);
100 Double_t pp[5] = { p[0]
106 Double_t c22 = x*x*cov[14] - 2*x*cov[12] + cov[ 5];
107 Double_t c32 = x*cov[13] - cov[ 8];
108 Double_t c20 = x*cov[10] - cov[ 3];
109 Double_t c21 = x*cov[11] - cov[ 4];
110 Double_t c42 = x*cov[14] - cov[12];
112 Double_t cc[15] = { cov[ 0]
115 , cov[ 6], cov[ 7], c32, cov[ 9]
116 , cov[10]*cnv, cov[11]*cnv, c42*cnv, cov[13]*cnv, cov[14]*cnv*cnv };
119 Double_t s = GetSnp();
120 Double_t t = GetTgl();
122 Int_t ncl = 0, nclPlane; AliTRDcluster *c = 0x0;
123 for(int iplane=0; iplane<6; iplane++){
124 fTrackletIndex[iplane] = -1;
125 fTracklet[iplane] = trklts[iplane];
127 for(int ic = 0; ic<AliTRDseed::knTimebins; ic++){
128 if(!fTracklet[iplane].IsUsable(ic)) continue;
129 if(!(c = fTracklet[iplane].GetClusters(ic))) continue;
131 fIndex[ncl] = fTracklet[iplane].GetIndexes(ic);
132 Double_t q = TMath::Abs(c->GetQ());
135 // This is not correctly. Has to be updated in the AliTRDtrackerV1::FollowBackProlonagation()
136 fdQdl[ncl] = q * (s*s < 1.) ? TMath::Sqrt((1-s*s)/(1+t*t)) : 1.;
140 //printf("%d N clusters plane %d [%d %d].\n", iplane, nclPlane, fTracklet[iplane].GetN2(), trklts[iplane].GetN());
142 //printf("N clusters in AliTRDtrackV1 %d.\n", ncl);
143 SetNumberOfClusters(/*ncl*/);
146 //_______________________________________________________________
147 Bool_t AliTRDtrackV1::CookLabel(Float_t wrong)
149 // set MC label for this tracklet
151 Int_t s[kMAXCLUSTERSPERTRACK][2];
152 for (Int_t i = 0; i < kMAXCLUSTERSPERTRACK; i++) {
159 AliTRDcluster *c = 0x0;
160 for (Int_t ip = 0; ip < AliESDtrack::kNPlane; ip++) {
161 if(fTrackletIndex[ip] < 0) continue;
162 for (Int_t ic = 0; ic < AliTRDseed::knTimebins; ic++) {
163 if(!(c = fTracklet[ip].GetClusters(ic))) continue;
164 for (Int_t k = 0; k < 3; k++) {
165 label = c->GetLabel(k);
169 while ((!labelAdded) && (j < kMAXCLUSTERSPERTRACK)) {
170 if ((s[j][0] == label) ||
185 for (Int_t i = 0; i < kMAXCLUSTERSPERTRACK; i++) {
186 if (s[i][1] <= max) continue;
191 if ((1. - Float_t(max)/GetNumberOfClusters()) > wrong) label = -label;
198 //_______________________________________________________________
199 Bool_t AliTRDtrackV1::CookPID()
202 // Cook the PID information
205 // CookdEdx(); // truncated mean ... do we still need it ?
207 // CookdEdxTimBin(seed->GetID());
209 // Sets the a priori probabilities
210 for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) {
211 fPID[ispec] = 1.0 / AliPID::kSPECIES;
214 // steer PID calculation @ tracklet level
215 Double_t *prob = 0x0;
217 for(int itrklt=0; itrklt<AliESDtrack::kNPlane; itrklt++){
218 //for (Int_t iSlice = 0; iSlice < AliESDtrack::kNSlice; iSlice++) fdEdxPlane[itrklt][iSlice] = -1.;
220 if(fTrackletIndex[itrklt]<0) continue;
221 if(!fTracklet[itrklt].IsOK()) continue;
222 if(!(prob = fTracklet[itrklt].GetProbability())) return kFALSE;
224 Int_t nspec = 0; // quality check of tracklet dEdx
225 for(int ispec=0; ispec<AliPID::kSPECIES; ispec++){
226 if(prob[ispec] < 0.) continue;
227 fPID[ispec] *= prob[ispec];
235 // no tracklet found for PID calculations
236 if(!fPIDquality) return kTRUE;
238 // slot for PID calculation @ track level
241 // normalize probabilities
242 Double_t probTotal = 0.0;
243 for (Int_t is = 0; is < AliPID::kSPECIES; is++) probTotal += fPID[is];
246 if (probTotal <= 0.0) {
247 AliWarning("The total probability over all species <= 0. This may be caused by some error in the reference data.");
251 for (Int_t iSpecies = 0; iSpecies < AliPID::kSPECIES; iSpecies++) fPID[iSpecies] /= probTotal;
256 //_______________________________________________________________
257 Float_t AliTRDtrackV1::GetMomentum(Int_t plane) const
260 // Get the momentum at a given plane
263 return plane >=0 && plane < 6 && fTrackletIndex[plane] >= 0 ? fTracklet[plane].GetMomentum() : -1.;
266 //_______________________________________________________________
267 Double_t AliTRDtrackV1::GetPredictedChi2(const AliTRDseedV1 *trklt) const
270 // Get the predicted chi2
273 Double_t x = trklt->GetX0();
274 Double_t p[2] = { trklt->GetYat(x)
275 , trklt->GetZat(x) };
277 trklt->GetCovAt(x, cov);
279 return AliExternalTrackParam::GetPredictedChi2(p, cov);
283 //_______________________________________________________________
284 Bool_t AliTRDtrackV1::IsOwner() const
287 // Check whether track owns the tracklets
290 for (Int_t ip = 0; ip < AliESDtrack::kNPlane; ip++) {
291 if(fTrackletIndex[ip] < 0) continue;
292 if(!fTracklet[ip].IsOwner()) return kFALSE;
297 //___________________________________________________________
298 void AliTRDtrackV1::SetNumberOfClusters()
300 // Calculate the number of clusters attached to this track
303 for(int ip=0; ip<6; ip++){
304 if(fTrackletIndex[ip] >= 0) ncls += fTracklet[ip].GetN();
306 AliKalmanTrack::SetNumberOfClusters(ncls);
310 //_______________________________________________________________
311 void AliTRDtrackV1::SetOwner(Bool_t own)
314 // Toggle ownership of tracklets
317 for (Int_t ip = 0; ip < AliESDtrack::kNPlane; ip++) {
318 if(fTrackletIndex[ip] < 0) continue;
319 //AliInfo(Form("p[%d] index[%d]", ip, fTrackletIndex[ip]));
320 fTracklet[ip].SetOwner(own);
324 //_______________________________________________________________
325 void AliTRDtrackV1::SetTracklet(AliTRDseedV1 *trklt, Int_t plane, Int_t index)
331 if(plane < 0 || plane >= AliESDtrack::kNPlane) return;
332 fTracklet[plane] = (*trklt);
333 fTrackletIndex[plane] = index;
336 //_______________________________________________________________
337 Bool_t AliTRDtrackV1::Update(AliTRDseedV1 *trklt, Double_t chisq)
340 // Update track parameters
343 Double_t x = trklt->GetX0();
344 Double_t p[2] = { trklt->GetYat(x)
345 , trklt->GetZat(x) };
347 trklt->GetCovAt(x, cov);
350 // AliInfo(Form("cov[%f %f %f]", cov[0], cov[1], cov[2]));
352 if(!AliExternalTrackParam::Update(p, cov)) return kFALSE;
355 // Register info to track
356 // Int_t n = GetNumberOfClusters();
357 // fIndex[n] = index;
359 SetNumberOfClusters(/*GetNumberOfClusters()+trklt->GetN()*/);
360 SetChi2(GetChi2() + chisq);
363 trklt->SetMomentum(GetP());
364 trklt->SetSnp(GetSnp());
365 trklt->SetTgl(GetTgl());
369 //_______________________________________________________________
370 void AliTRDtrackV1::UpdateESDtrack(AliESDtrack *track)
373 // Update the ESD track
378 for (Int_t ip = 0; ip < AliESDtrack::kNPlane; ip++) {
379 if(fTrackletIndex[ip] < 0) continue;
380 fTracklet[ip].CookdEdx(AliESDtrack::kNSlice);
381 dedx = fTracklet[ip].GetdEdx();
382 for (Int_t js = 0; js < AliESDtrack::kNSlice; js++) track->SetTRDsignals(dedx[js], ip, js);
383 //track->SetTRDTimBin(fTimBinPlane[i], i);
387 track->SetTRDpid(fPID);
388 track->SetTRDpidQuality(fPIDquality);