Use event specie to identufy laser events
[u/mrichter/AliRoot.git] / TPC / AliTPCtrack.cxx
CommitLineData
73042f01 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
803d1ab0 16/* $Id$ */
0938c7d3 17
73042f01 18//-----------------------------------------------------------------
19// Implementation of the TPC track class
7c27857f 20// This class is used by the AliTPCtracker class
21// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
73042f01 22//-----------------------------------------------------------------
23
19364939 24#include <Riostream.h>
b19a0509 25
73042f01 26#include "AliTPCtrack.h"
7c27857f 27#include "AliCluster.h"
6c94f330 28#include "AliTracker.h"
83d73500 29#include "AliESDtrack.h"
58e536c5 30#include "AliESDVertex.h"
d08b7e45 31#include "TTreeStream.h"
eba5bca1 32#include "AliTPCRecoParam.h"
33#include "AliTPCReconstructor.h"
73042f01 34ClassImp(AliTPCtrack)
b9de75e1 35
73042f01 36//_________________________________________________________________________
e13cc304 37AliTPCtrack::AliTPCtrack():
38 AliKalmanTrack(),
e13cc304 39 fdEdx(0),
e13cc304 40 fSdEdx(1e10),
41 fNFoundable(0),
42 fBConstrain(kFALSE),
43 fLastPoint(-1),
44 fFirstPoint(-1),
45 fRemoval(0),
46 fTrackType(0),
47 fLab2(-1),
48 fNShared(0),
49 fReference()
e3f7105e 50{
2fc0c115 51 //-------------------------------------------------
52 // default constructor
53 //-------------------------------------------------
e13cc304 54 for (Int_t i=0; i<kMaxRow;i++) fIndex[i]=-2;
55 for (Int_t i=0; i<4;i++) fPoints[i]=0.;
56 for (Int_t i=0; i<12;i++) fKinkPoint[i]=0.;
c84a5e9e 57 for (Int_t i=0; i<3;i++) fKinkIndexes[i]=0;
e13cc304 58 for (Int_t i=0; i<3;i++) fV0Indexes[i]=0;
e3f7105e 59}
60
61//_________________________________________________________________________
51ad6848 62
63
64
6c94f330 65AliTPCtrack::AliTPCtrack(Double_t x, Double_t alpha, const Double_t p[5],
66 const Double_t cov[15], Int_t index) :
e13cc304 67 AliKalmanTrack(),
e13cc304 68 fdEdx(0),
e13cc304 69 fSdEdx(1e10),
70 fNFoundable(0),
71 fBConstrain(kFALSE),
72 fLastPoint(0),
73 fFirstPoint(0),
74 fRemoval(0),
75 fTrackType(0),
76 fLab2(0),
77 fNShared(0),
78 fReference()
79{
73042f01 80 //-----------------------------------------------------------------
81 // This is the main track constructor.
82 //-----------------------------------------------------------------
beceeaec 83 Double_t cnv=1./(AliTracker::GetBz()*kB2C);
6c94f330 84
85 Double_t pp[5]={
86 p[0],
87 p[1],
88 x*p[4] - p[2],
89 p[3],
90 p[4]*cnv
91 };
92
93 Double_t c22 = x*x*cov[14] - 2*x*cov[12] + cov[5];
94 Double_t c32 = x*cov[13] - cov[8];
95 Double_t c20 = x*cov[10] - cov[3],
96 c21 = x*cov[11] - cov[4], c42 = x*cov[14] - cov[12];
97
98 Double_t cc[15]={
99 cov[0 ],
100 cov[1 ], cov[2 ],
101 c20, c21, c22,
102 cov[6 ], cov[7 ], c32, cov[9 ],
103 cov[10]*cnv, cov[11]*cnv, c42*cnv, cov[13]*cnv, cov[14]*cnv*cnv
104 };
105
86fc30dc 106 Double_t mostProbablePt=AliExternalTrackParam::GetMostProbablePt();
107 Double_t p0=TMath::Sign(1/mostProbablePt,pp[4]);
beceeaec 108 Double_t w0=cc[14]/(cc[14] + p0*p0), w1=p0*p0/(cc[14] + p0*p0);
109 pp[4] = w0*p0 + w1*pp[4];
110 cc[10]*=w1; cc[11]*=w1; cc[12]*=w1; cc[13]*=w1; cc[14]*=w1;
111
6c94f330 112 Set(x,alpha,pp,cc);
c84a5e9e 113
be9c5115 114 SetNumberOfClusters(1);
e13cc304 115
116 fIndex[0]=index;
9fe5b2ff 117 for (Int_t i=1; i<kMaxRow;i++) fIndex[i]=-2;
118 for (Int_t i=0; i<4;i++) fPoints[i]=0.;
119 for (Int_t i=0; i<12;i++) fKinkPoint[i]=0.;
51ad6848 120 for (Int_t i=0; i<3;i++) fKinkIndexes[i]=0;
81e97e0d 121 for (Int_t i=0; i<3;i++) fV0Indexes[i]=0;
b9de75e1 122}
123
124//_____________________________________________________________________________
d08b7e45 125AliTPCtrack::AliTPCtrack(const AliESDtrack& t, TTreeSRedirector *pcstream) :
e13cc304 126 AliKalmanTrack(),
6c94f330 127 fdEdx(t.GetTPCsignal()),
e13cc304 128 fSdEdx(1e10),
129 fNFoundable(0),
130 fBConstrain(kFALSE),
131 fLastPoint(0),
132 fFirstPoint(0),
133 fRemoval(0),
134 fTrackType(0),
135 fLab2(0),
136 fNShared(0),
137 fReference()
138{
83d73500 139 //-----------------------------------------------------------------
140 // Conversion AliESDtrack -> AliTPCtrack.
141 //-----------------------------------------------------------------
d08b7e45 142 const Double_t kmaxC[4]={10,10,0.1,0.1}; // cuts on the rms /fP0,fP1,fP2,fP3
83d73500 143 SetNumberOfClusters(t.GetTPCclusters(fIndex));
144 SetLabel(t.GetLabel());
145 SetMass(t.GetMass());
9fe5b2ff 146 for (Int_t i=0; i<4;i++) fPoints[i]=0.;
147 for (Int_t i=0; i<12;i++) fKinkPoint[i]=0.;
148 for (Int_t i=0; i<3;i++) fKinkIndexes[i]=0;
149 for (Int_t i=0; i<3;i++) fV0Indexes[i]=0;
d08b7e45 150 //
151 // choose parameters to start
152 //
eba5bca1 153 const AliTPCRecoParam * recoParam = AliTPCReconstructor::GetRecoParam();
d08b7e45 154 Int_t reject=0;
155 AliExternalTrackParam param(t);
e5b3d825 156
d08b7e45 157 const AliExternalTrackParam *tpcout=(t.GetFriendTrack())? ((AliESDfriendTrack*)(t.GetFriendTrack()))->GetTPCOut():0;
158 const AliExternalTrackParam *tpcin = t.GetInnerParam();
159 const AliExternalTrackParam *tpc=(tpcout)?tpcout:tpcin;
160 if (!tpc) tpc=&param;
eba5bca1 161 Bool_t isOK=recoParam->GetUseOuterDetectors();
d08b7e45 162 if (param.GetCovariance()[0]>kmaxC[0]*kmaxC[0]) isOK=kFALSE;
163 if (param.GetCovariance()[2]>kmaxC[1]*kmaxC[1]) isOK=kFALSE;
164 if (param.GetCovariance()[5]>kmaxC[2]*kmaxC[2]) isOK=kFALSE;
165 if (param.GetCovariance()[9]>kmaxC[3]*kmaxC[3]) isOK=kFALSE;
88a42d3e 166 param.Rotate(tpc->GetAlpha());
167 Double_t oldX=param.GetX(), oldY=param.GetY(), oldZ=param.GetZ();
eba5bca1 168 if (!isOK ){
d08b7e45 169 param=*tpc;
170 isOK=kTRUE;
171 reject=1;
172 }
d08b7e45 173 isOK=AliTracker::PropagateTrackToBxByBz(&param,tpc->GetX(),t.GetMass(),2.,kFALSE);
174 if (param.GetCovariance()[0]>kmaxC[0]*kmaxC[0]) isOK=kFALSE;
175 if (param.GetCovariance()[2]>kmaxC[1]*kmaxC[1]) isOK=kFALSE;
176 if (param.GetCovariance()[5]>kmaxC[2]*kmaxC[2]) isOK=kFALSE;
177 if (param.GetCovariance()[9]>kmaxC[3]*kmaxC[3]) isOK=kFALSE;
178 if (!isOK){
179 param=*tpc;
180 isOK=kTRUE;
181 reject=2;
182 }
183 if (reject>0){
184 param.ResetCovariance(4.); // reset covariance if start from backup param
185 }
186 //
187 //
188 if (pcstream){
189 AliExternalTrackParam dummy;
190 AliExternalTrackParam *ptpc=(AliExternalTrackParam *)tpc;
191 if (!ptpc) ptpc=&dummy;
192 AliESDtrack *esd= (AliESDtrack *)&t;
193 (*pcstream)<<"trackP"<<
194 "reject="<<reject<< // flag - rejection of current esd track parameters
195 "esd.="<<esd<< // original esd track
196 "tr.="<<&param<< // starting track parameters
197 "out.="<<ptpc<< // backup tpc parameters
198 "\n";
199 }
200
201 Set(param.GetX(),param.GetAlpha(),param.GetParameter(),param.GetCovariance());
83d73500 202
0798b21e 203 if ((t.GetStatus()&AliESDtrack::kTIME) == 0) return;
83d73500 204 StartTimeIntegral();
205 Double_t times[10]; t.GetIntegratedTimes(times); SetIntegratedTimes(times);
206 SetIntegratedLength(t.GetIntegratedLength());
e5b3d825 207
208 if (GetX()>oldX) {
209 Double_t dX=GetX()-oldX, dY=GetY()-oldY, dZ=GetZ()-oldZ;
210 Double_t d=TMath::Sqrt(dX*dX + dY*dY + dZ*dZ);
211 AddTimeStep(d);
212 }
83d73500 213}
214
215//_____________________________________________________________________________
e13cc304 216AliTPCtrack::AliTPCtrack(const AliTPCtrack& t) :
217 AliKalmanTrack(t),
e13cc304 218 fdEdx(t.fdEdx),
e13cc304 219 fSdEdx(t.fSdEdx),
220 fNFoundable(t.fNFoundable),
221 fBConstrain(t.fBConstrain),
222 fLastPoint(t.fLastPoint),
223 fFirstPoint(t.fFirstPoint),
224 fRemoval(t.fRemoval),
225 fTrackType(t.fTrackType),
226 fLab2(t.fLab2),
227 fNShared(t.fNShared),
228 fReference(t.fReference)
229
230{
73042f01 231 //-----------------------------------------------------------------
232 // This is a track copy constructor.
233 //-----------------------------------------------------------------
6c94f330 234 Set(t.GetX(),t.GetAlpha(),t.GetParameter(),t.GetCovariance());
235
51584047 236 for (Int_t i=0; i<kMaxRow; i++) fIndex[i]=t.fIndex[i];
9fe5b2ff 237 for (Int_t i=0; i<4;i++) fPoints[i]=t.fPoints[i];
238 for (Int_t i=0; i<12;i++) fKinkPoint[i]=t.fKinkPoint[i];
51ad6848 239 for (Int_t i=0; i<3;i++) fKinkIndexes[i]=t.fKinkIndexes[i];
81e97e0d 240 for (Int_t i=0; i<3;i++) fV0Indexes[i]=t.fV0Indexes[i];
be9c5115 241}
242
316c6cd9 243AliTPCtrack& AliTPCtrack::operator=(const AliTPCtrack& o){
244 if(this!=&o){
245 AliKalmanTrack::operator=(o);
246 fdEdx = o.fdEdx;
247 for(Int_t i = 0;i<kMaxRow;++i)fIndex[i] = o.fIndex[i];
248 for(Int_t i = 0;i<4;++i)fPoints[i] = o.fPoints[i];
249 fSdEdx = o.fSdEdx;
250 fNFoundable = o.fNFoundable;
251 fBConstrain = o.fBConstrain;
252 fLastPoint = o.fLastPoint;
253 fFirstPoint = o.fFirstPoint;
254 fTrackType = o.fTrackType;
255 fLab2 = o.fLab2;
256 fNShared = o.fNShared;
257 fReference = o.fReference;
258 for(Int_t i = 0;i<12;++i) fKinkPoint[i] = o.fKinkPoint[i];
259
260 for(Int_t i = 0;i<3;++i){
261 fKinkIndexes[i] = o.fKinkIndexes[i];
262 fV0Indexes[i] = o.fV0Indexes[i];
263 }
264 }
265 return *this;
266
267}
268
269
be9c5115 270//_____________________________________________________________________________
271Int_t AliTPCtrack::Compare(const TObject *o) const {
272 //-----------------------------------------------------------------
273 // This function compares tracks according to the their curvature
274 //-----------------------------------------------------------------
275 AliTPCtrack *t=(AliTPCtrack*)o;
6c23ffed 276 //Double_t co=t->OneOverPt();
277 //Double_t c = OneOverPt();
024a7fe9 278 Double_t co=t->GetSigmaY2()*t->GetSigmaZ2();
279 Double_t c =GetSigmaY2()*GetSigmaZ2();
be9c5115 280 if (c>co) return 1;
281 else if (c<co) return -1;
282 return 0;
283}
284
6c94f330 285Double_t AliTPCtrack::GetPredictedChi2(const AliCluster *c) const {
be9c5115 286 //-----------------------------------------------------------------
287 // This function calculates a predicted chi2 increment.
288 //-----------------------------------------------------------------
6c94f330 289 Double_t p[2]={c->GetY(), c->GetZ()};
290 Double_t cov[3]={c->GetSigmaY2(), 0., c->GetSigmaZ2()};
291 return AliExternalTrackParam::GetPredictedChi2(p,cov);
c630aafd 292}
293
73042f01 294//_____________________________________________________________________________
6c94f330 295Bool_t AliTPCtrack::PropagateTo(Double_t xk,Double_t rho,Double_t x0) {
73042f01 296 //-----------------------------------------------------------------
6c94f330 297 // This function propagates a track to a reference plane x=xk.
298 // rho - density of the crossed matrial (g/cm^3)
299 // x0 - radiation length of the crossed material (g/cm^2)
73042f01 300 //-----------------------------------------------------------------
a326cd74 301 //
302 Double_t bz=GetBz();
303 Double_t zat=0;
304 if (!GetZAt(xk, bz,zat)) return kFALSE;
305 if (TMath::Abs(zat)>250.){
306 // Don't propagate track outside of the fiducial volume - material budget not proper one
307 //
308 //AliWarning("Propagate outside of fiducial volume");
309 return kFALSE;
310 }
311
6c94f330 312 Double_t oldX=GetX(), oldY=GetY(), oldZ=GetZ();
273ccc8a 313 //if (!AliExternalTrackParam::PropagateTo(xk,bz)) return kFALSE;
314 Double_t b[3]; GetBxByBz(b);
315 if (!AliExternalTrackParam::PropagateToBxByBz(xk,b)) return kFALSE;
73042f01 316
6c94f330 317 Double_t d = TMath::Sqrt((GetX()-oldX)*(GetX()-oldX) +
318 (GetY()-oldY)*(GetY()-oldY) +
319 (GetZ()-oldZ)*(GetZ()-oldZ));
320 if (IsStartedTimeIntegral() && GetX()>oldX) AddTimeStep(d);
321
322 if (oldX < xk) d = -d;
0f30c08a 323 if (!AliExternalTrackParam::CorrectForMeanMaterial(d*rho/x0,d*rho,GetMass(),
324 kFALSE,AliExternalTrackParam::BetheBlochGas)) return kFALSE;
e3f7105e 325
6c94f330 326 return kTRUE;
73042f01 327}
328
329//_____________________________________________________________________________
6c94f330 330Bool_t
0f30c08a 331AliTPCtrack::PropagateToVertex(const AliESDVertex *v,Double_t rho,Double_t x0)
73042f01 332{
333 //-----------------------------------------------------------------
0f30c08a 334 // This function propagates tracks to the vertex
335 // rho - density of the crossed matrial (g/cm3)
336 // x0 - radiation length of the crossed material (g/cm2)
73042f01 337 //-----------------------------------------------------------------
0f30c08a 338 Double_t oldX=GetX(), oldY=GetY(), oldZ=GetZ();
339
273ccc8a 340 //Double_t bz=GetBz();
341 //if (!PropagateToDCA(v,bz,kVeryBig)) return kFALSE;
342 Double_t b[3]; GetBxByBz(b);
343 if (!PropagateToDCABxByBz(v,b,kVeryBig)) return kFALSE;
0f30c08a 344
345 Double_t d = TMath::Sqrt((GetX()-oldX)*(GetX()-oldX) +
346 (GetY()-oldY)*(GetY()-oldY) +
347 (GetZ()-oldZ)*(GetZ()-oldZ));
348
349 if (oldX < GetX()) d = -d;
350 if (!AliExternalTrackParam::CorrectForMeanMaterial(d*rho/x0,d*rho,GetMass(),
351 kFALSE,AliExternalTrackParam::BetheBlochGas)) return kFALSE;
352
353 return kTRUE;
73042f01 354}
355
356//_____________________________________________________________________________
6c94f330 357Bool_t AliTPCtrack::Update(const AliCluster *c, Double_t chisq, Int_t index) {
73042f01 358 //-----------------------------------------------------------------
359 // This function associates a cluster with this track.
360 //-----------------------------------------------------------------
6c94f330 361 Double_t p[2]={c->GetY(), c->GetZ()};
362 Double_t cov[3]={c->GetSigmaY2(), 0., c->GetSigmaZ2()};
73042f01 363
6c94f330 364 if (!AliExternalTrackParam::Update(p,cov)) return kFALSE;
73042f01 365
70e20dcf 366 AliTracker::FillResiduals(this,p,cov,c->GetVolumeId());
367
be9c5115 368 Int_t n=GetNumberOfClusters();
369 fIndex[n]=index;
370 SetNumberOfClusters(n+1);
371 SetChi2(GetChi2()+chisq);
372
6c94f330 373 return kTRUE;
73042f01 374}
375
51584047 376////////////////////////////////////////////////////////////////////////
377// MI ADDITION
378
379Float_t AliTPCtrack::Density(Int_t row0, Int_t row1)
380{
381 //
382 // calculate cluster density
383 Int_t good = 0;
384 Int_t found = 0;
385 //if (row0<fFirstPoint) row0 = fFirstPoint;
386 if (row1>fLastPoint) row1 = fLastPoint;
387
388
389 for (Int_t i=row0;i<=row1;i++){
390 // Int_t index = fClusterIndex[i];
391 Int_t index = fIndex[i];
392 if (index!=-1) good++;
393 if (index>0) found++;
394 }
395 Float_t density=0;
396 if (good>0) density = Float_t(found)/Float_t(good);
397 return density;
398}
399
400
401Float_t AliTPCtrack::Density2(Int_t row0, Int_t row1)
402{
403 //
404 // calculate cluster density
405 Int_t good = 0;
406 Int_t found = 0;
407 //
408 for (Int_t i=row0;i<=row1;i++){
409 Int_t index = fIndex[i];
410 if (index!=-1) good++;
411 if (index>0) found++;
412 }
413 Float_t density=0;
414 if (good>0) density = Float_t(found)/Float_t(good);
415 return density;
416}
417
51ad6848 418void AliTPCtrack::UpdatePoints()
419{
420 //--------------------------------------------------
421 //calculates first ,amx dens and last points
422 //--------------------------------------------------
423 Float_t density[160];
424 for (Int_t i=0;i<160;i++) density[i]=-1.;
425 fPoints[0]= 160;
426 fPoints[1] = -1;
427 //
428 Int_t ngood=0;
429 Int_t undeff=0;
430 Int_t nall =0;
431 Int_t range=20;
432 for (Int_t i=0;i<160;i++){
433 Int_t last = i-range;
434 if (nall<range) nall++;
435 if (last>=0){
436 if (fIndex[last]>0&& (fIndex[last]&0x8000)==0) ngood--;
437 if (fIndex[last]==-1) undeff--;
438 }
439 if (fIndex[i]>0&& (fIndex[i]&0x8000)==0) ngood++;
440 if (fIndex[i]==-1) undeff++;
441 if (nall==range &&undeff<range/2) density[i-range/2] = Float_t(ngood)/Float_t(nall-undeff);
442 }
443 Float_t maxdens=0;
444 Int_t indexmax =0;
445 for (Int_t i=0;i<160;i++){
446 if (density[i]<0) continue;
447 if (density[i]>maxdens){
448 maxdens=density[i];
449 indexmax=i;
450 }
451 }
452 //
453 //max dens point
454 fPoints[3] = maxdens;
455 fPoints[1] = indexmax;
456 //
457 // last point
458 for (Int_t i=indexmax;i<160;i++){
459 if (density[i]<0) continue;
460 if (density[i]<maxdens/2.) {
461 break;
462 }
463 fPoints[2]=i;
464 }
465 //
466 // first point
467 for (Int_t i=indexmax;i>0;i--){
468 if (density[i]<0) continue;
469 if (density[i]<maxdens/2.) {
470 break;
471 }
472 fPoints[0]=i;
473 }
474 //
475}
5773defd 476