Added AliL3Stopwatch.
[u/mrichter/AliRoot.git] / HLT / comp / AliL3ModelTrack.cxx
CommitLineData
735e167e 1//$Id$
2
3// Author: Anders Vestbo <mailto:vestbo$fi.uib.no>
4//*-- Copyright &copy ASV
5
6#include <stream.h>
7#include <string.h>
8487f697 8#include <stdlib.h>
95a00d93 9#include <math.h>
735e167e 10
11#include "AliL3ModelTrack.h"
95a00d93 12#include "AliL3Transform.h"
735e167e 13
4a838220 14//_____________________________________________________________
15// AliL3ModelTrack
16//
17//
18
735e167e 19ClassImp(AliL3ModelTrack)
20
21AliL3ModelTrack::AliL3ModelTrack()
22{
23 fNClusters = 0;
24 fClusters = 0;
029912b7 25 fOverlap = 0;
735e167e 26 fPad=0;
27 fTime=0;
28 fClusterCharge=0;
95a00d93 29 fTrackModel=0;
8487f697 30 fLabel=0;
735e167e 31}
32
33
34AliL3ModelTrack::~AliL3ModelTrack()
35{
36 if(fClusters)
37 delete [] fClusters;
38 if(fPad)
39 delete [] fPad;
40 if(fTime)
41 delete [] fTime;
029912b7 42 if(fOverlap)
43 delete [] fOverlap;
95a00d93 44 if(fTrackModel)
45 delete fTrackModel;
735e167e 46}
47
48void AliL3ModelTrack::Init(Int_t slice,Int_t patch)
49{
50 fNClusters = 0;
029912b7 51 fSlice=slice;
52 fPatch=patch;
4a838220 53 Int_t nrows = AliL3Transform::GetNRows(fPatch);
95a00d93 54 fClusters = new AliL3ClusterModel[nrows];
029912b7 55 fPad = new Float_t[nrows];
56 fTime = new Float_t[nrows];
95a00d93 57 fTrackModel = new AliL3TrackModel;
029912b7 58 fOverlap = new Int_t[nrows];
8487f697 59
029912b7 60 memset(fClusters,0,nrows*sizeof(AliL3ClusterModel));
61 memset(fPad,0,nrows*sizeof(Float_t));
62 memset(fTime,0,nrows*sizeof(Float_t));
63 memset(fTrackModel,0,sizeof(AliL3TrackModel));
64 for(Int_t i=0; i<nrows; i++)
65 fOverlap[i]=-1;
8487f697 66#ifdef do_mc
67 for(Int_t i=0; i<nrows; i++)
68 fClusters[i].fTrackID[0]=fClusters[i].fTrackID[1]=fClusters[i].fTrackID[2]=-2;
69#endif
4a838220 70 fClusterCharge = 100;
029912b7 71
4a838220 72 // 100 micrometers:
73 fXYResidualQ = 0.01/AliL3Transform::GetPadPitchWidth(patch);
74 fZResidualQ = 0.01/AliL3Transform::GetPadPitchWidth(patch);
92a876e2 75
4a838220 76 fXYWidthQ = 0.005/AliL3Transform::GetPadPitchWidth(patch);
77 fZWidthQ = 0.005/AliL3Transform::GetPadPitchWidth(patch);
735e167e 78}
79
80
4a838220 81void AliL3ModelTrack::SetCluster(Int_t row,Float_t fpad,Float_t ftime,Float_t charge,Float_t sigmaY2,Float_t sigmaZ2,Int_t npads)
735e167e 82{
4a838220 83 Int_t index = row - AliL3Transform::GetFirstRow(fPatch);
029912b7 84 if(index != fNClusters)
85 cout<<"AliL3ModelTrack::SetCluster() : Mismatch ; index: "<<index<<" nclusters "<<fNClusters<<endl;
86
4a838220 87 if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
029912b7 88 {
89 cerr<<"AliL3ModelTrack::SetCluster() : Wrong index: "<<index<<" row "<<row<<endl;
90 return;
91 }
92 AliL3ClusterModel *cl = GetClusterModel(row);
735e167e 93 if(!charge)
95a00d93 94 cl->fPresent = kFALSE;
735e167e 95 else
96 {
95a00d93 97 cl->fPresent = kTRUE;
029912b7 98 cl->fDTime = (ftime - GetTimeHit(row))/fXYResidualQ;
99 cl->fDPad = (fpad - GetPadHit(row))/fZResidualQ;
100 cl->fDCharge = charge - fClusterCharge;
4a838220 101 cl->fDSigmaY2 = (sigmaY2 - GetParSigmaY2(row))/fXYWidthQ;
102 cl->fDSigmaZ2 = (sigmaZ2 - GetParSigmaZ2(row))/fZWidthQ;
103 cl->fNPads = npads;
735e167e 104 }
029912b7 105
735e167e 106 fNClusters++;
107}
108
4a838220 109Int_t AliL3ModelTrack::CheckClustersQuality(UInt_t npads=3)
110{
111
112 //Check the quality of clusters,- remove clusters with less than
113 //npads.
114 //Returns the number of good clusters left.
115
116 Int_t count=0;
117
118 for(Int_t i=AliL3Transform::GetFirstRow(fPatch); i<=AliL3Transform::GetLastRow(fPatch); i++)
119 {
120 AliL3ClusterModel *cl = GetClusterModel(i);
121 if(cl->fNPads < npads)
122 cl->fPresent = kFALSE;
123 if(cl->fPresent)
124 count++;
125 }
126
127 return count;
128}
95a00d93 129
95a00d93 130void AliL3ModelTrack::FillModel()
131{
029912b7 132 //Fill the track structure
133
134 if(!fTrackModel)
135 {
136 cerr<<"AliL3ModelTrack::FillModel() : No trackmodel "<<endl;
137 return;
138 }
95a00d93 139 fTrackModel->fKappa = GetKappa();
140 fTrackModel->fFirstPointX = GetFirstPointX();
141 fTrackModel->fFirstPointY = GetFirstPointY();
142 fTrackModel->fFirstPointZ = GetFirstPointZ();
143 fTrackModel->fTgl = GetTgl();
144 fTrackModel->fPsi = GetPsi();
029912b7 145 fTrackModel->fLength = (Short_t)GetLength();
95a00d93 146 fTrackModel->fClusterCharge = fClusterCharge;
147 fTrackModel->fNClusters = fNClusters;
148
149}
150
029912b7 151void AliL3ModelTrack::FillTrack()
95a00d93 152{
029912b7 153 //Fill the track parameters from the structure.
154
155 if(!fTrackModel)
156 {
157 cerr<<"AliL3ModelTrack::FillTrack() : No data!!"<<endl;
158 return;
159 }
160 SetKappa(fTrackModel->fKappa);
161 SetCharge((-1*(Int_t)copysign(1.,GetKappa())));
162 SetFirstPoint(fTrackModel->fFirstPointX,fTrackModel->fFirstPointY,fTrackModel->fFirstPointZ);
163 SetTgl(fTrackModel->fTgl);
164 SetPsi(fTrackModel->fPsi);
165 SetLength(fTrackModel->fLength);
166 fClusterCharge=fTrackModel->fClusterCharge;
167 fNClusters = fTrackModel->fNClusters;
4a838220 168 SetPt((BFACT*AliL3Transform::GetBField())/fabs(GetKappa()));
029912b7 169
170 CalculateHelix();
171
029912b7 172 Float_t hit[3];
173 Int_t sector,row;
4a838220 174 for(Int_t i=AliL3Transform::GetFirstRow(fPatch); i<=AliL3Transform::GetLastRow(fPatch); i++)
029912b7 175 {
176 AliL3ClusterModel *cl = GetClusterModel(i);
177 if(!cl) continue;
178 GetCrossingPoint(i,hit);
4a838220 179 AliL3Transform::Slice2Sector(fSlice,i,sector,row);
180 AliL3Transform::Local2Raw(hit,sector,row);
029912b7 181 SetPadHit(i,hit[1]);
182 SetTimeHit(i,hit[2]);
183 }
184}
95a00d93 185
8487f697 186
187void AliL3ModelTrack::SetTrackID(Int_t row,Int_t *trackID)
188{
189#ifdef do_mc
190 AliL3ClusterModel *cluster = GetClusterModel(row);
191 cluster->fTrackID[0] = trackID[0];
192 cluster->fTrackID[1] = trackID[1];
193 cluster->fTrackID[2] = trackID[2];
194 return;
195#endif
196 cerr<<"AliL3ModelTrack::SetTrackID : Compile with do_mc flag"<<endl;
197}
198
199
029912b7 200void AliL3ModelTrack::SetPadHit(Int_t row,Float_t pad)
201{
4a838220 202 Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
203 if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
95a00d93 204 {
029912b7 205 cerr<<"AliL3ModelTrack::SetPadHit() : Wrong index: "<<index<<endl;
206 return;
95a00d93 207 }
029912b7 208 fPad[index]=pad;
209
210}
211
212void AliL3ModelTrack::SetTimeHit(Int_t row,Float_t time)
213{
4a838220 214 Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
215 if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
029912b7 216 {
217 cerr<<"AliL3ModelTrack::SetTimeHit() : Wrong index: "<<index<<endl;
218 return;
219 }
220 fTime[index]=time;
221}
222
223void AliL3ModelTrack::SetOverlap(Int_t row,Int_t id)
224{
4a838220 225 Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
226 if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
029912b7 227 {
228 cerr<<"AliL3ModelTrack::SetOverlap() : Wrong index: "<<index<<endl;
229 return;
230 }
231 fOverlap[index]=id;
232}
233
8487f697 234
235Int_t AliL3ModelTrack::GetTrackID(Int_t row,Int_t index)
236{
237
238#ifdef do_mc
239 AliL3ClusterModel *cl = GetClusterModel(row);
240 return cl->fTrackID[index];
241#endif
242 cerr<<"AliL3ModelTrack::GetTrackID : Compile with do_mc flag"<<endl;
243}
244
245
4a838220 246Int_t AliL3ModelTrack::GetNPads(Int_t row)
247{
248 AliL3ClusterModel *cl = GetClusterModel(row);
249 return cl->fNPads;
250}
251
029912b7 252Bool_t AliL3ModelTrack::GetPad(Int_t row,Float_t &pad)
253{
254 //(ftime - GetTimeHit(fNClusters))/fXYResidualQ;
255 //(fpad - GetPadHit(fNClusters))/fZResidualQ;
256
257 AliL3ClusterModel *cl = GetClusterModel(row);
258 pad = cl->fDPad*fXYResidualQ + GetPadHit(row);
259
260 return (Bool_t)cl->fPresent;
261}
262
263Bool_t AliL3ModelTrack::GetTime(Int_t row,Float_t &time)
264{
265 AliL3ClusterModel *cl = GetClusterModel(row);
266 time = cl->fDTime*fZResidualQ + GetTimeHit(row);
267
268 return (Bool_t)cl->fPresent;
269}
270
271Bool_t AliL3ModelTrack::GetClusterCharge(Int_t row,Int_t &charge)
272{
273 AliL3ClusterModel *cl = GetClusterModel(row);
274 charge = (Int_t)cl->fDCharge + fClusterCharge;
275
276 return (Bool_t)cl->fPresent;
277}
278
279Bool_t AliL3ModelTrack::GetXYWidth(Int_t row,Float_t &width)
280{
281 AliL3ClusterModel *cl = GetClusterModel(row);
4a838220 282 width = cl->fDSigmaY2*fXYWidthQ + GetParSigmaY2(row);
029912b7 283
284 return (Bool_t)cl->fPresent;
285}
286
287Bool_t AliL3ModelTrack::GetZWidth(Int_t row,Float_t &width)
288{
289 AliL3ClusterModel *cl = GetClusterModel(row);
4a838220 290 width = cl->fDSigmaZ2*fZWidthQ + GetParSigmaZ2(row);
029912b7 291
292 return (Bool_t)cl->fPresent;
95a00d93 293}
294
295Bool_t AliL3ModelTrack::GetPadResidual(Int_t row,Float_t &res)
296{
4a838220 297 AliL3ClusterModel *cl = GetClusterModel(row);
298 res = cl->fDPad;
299 return cl->fPresent;
95a00d93 300}
301
302Bool_t AliL3ModelTrack::GetTimeResidual(Int_t row,Float_t &res)
303{
4a838220 304 AliL3ClusterModel *cl = GetClusterModel(row);
305 res = cl->fDTime;
306 return cl->fPresent;
307}
308
309Bool_t AliL3ModelTrack::GetXYWidthResidual(Int_t row,Float_t &res)
310{
311 AliL3ClusterModel *cl = GetClusterModel(row);
312 res = cl->fDSigmaY2;
313 return cl->fPresent;
314}
315
316Bool_t AliL3ModelTrack::GetZWidthResidual(Int_t row,Float_t &res)
317{
318 AliL3ClusterModel *cl = GetClusterModel(row);
319 res = cl->fDSigmaZ2;
320 return cl->fPresent;
95a00d93 321}
322
029912b7 323Float_t AliL3ModelTrack::GetPadHit(Int_t row)
324{
4a838220 325 Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
326 if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
029912b7 327 {
328 cerr<<"AliL3ModelTrack::GetPadHit() : Wrong index: "<<index<<" row "<<row<<endl;
329 return 0;
330 }
331 return fPad[index];
332}
333
334Float_t AliL3ModelTrack::GetTimeHit(Int_t row)
335{
4a838220 336 Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
337 if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
029912b7 338 {
339 cerr<<"AliL3ModelTrack::GetTimeHit() : Wrong index: "<<index<<" row "<<row<<endl;
340 return 0;
341 }
342 return fTime[index];
343}
344
345Int_t AliL3ModelTrack::GetOverlap(Int_t row)
346{
4a838220 347 Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
348 if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
029912b7 349 {
350 cerr<<"AliL3ModelTrack::GetOverlap() : Wrong index: "<<index<<endl;
351 return 0;
352 }
353 return fOverlap[index];
354}
355
356AliL3ClusterModel *AliL3ModelTrack::GetClusterModel(Int_t row)
357{
358 if(!fClusters) return 0;
4a838220 359 Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
360 if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
029912b7 361 {
362 cerr<<"AliL3ModelTrack::GetClusterModel() : Wrong index: "<<index<<endl;
363 return 0;
364 }
365 return &fClusters[index];
366}
367
368void AliL3ModelTrack::Print()
369{
370 //Print info
371
4a838220 372 cout<<"----Slice "<<fSlice<<" Patch "<<fPatch<<"----"<<endl;
029912b7 373 cout<<"First point "<<GetFirstPointX()<<" "<<GetFirstPointY()<<" "<<GetFirstPointZ()<<endl;
374 cout<<"Last point "<<GetLastPointX()<<" "<<GetLastPointY()<<" "<<GetLastPointZ()<<endl;
375 cout<<"Pt "<<GetPt()<<" kappa "<<GetKappa()<<" tgl "<<GetTgl()<<" psi "<<GetPsi()<<" charge "<<GetCharge()<<endl;
376 cout<<"Center "<<GetCenterX()<<" "<<GetCenterY()<<endl<<endl;
377 cout<<"NHits "<<GetNClusters()<<endl;
378 cout<<"Clusters:"<<endl;
379
4a838220 380 for(Int_t i=AliL3Transform::GetFirstRow(fPatch); i<=AliL3Transform::GetLastRow(fPatch); i++)
029912b7 381 {
382 AliL3ClusterModel *cl = GetClusterModel(i);
383
384 if(!cl->fPresent)
385 cout<<i<<" Empty"<<" Padcrossing "<<GetPadHit(i)<<" Timecrossing "<<GetTimeHit(i)<<" ";
386 else
387 {
388 cout<<i<<" Dpad "<<cl->fDPad<<" Dtime "<<cl->fDTime<<" Dcharge "<<cl->fDCharge;
4a838220 389 cout<<" DsigmaY2 "<<cl->fDSigmaY2<<" DsigmaZ2 "<<cl->fDSigmaZ2;
029912b7 390 cout<<" Padcrossing "<<GetPadHit(i)<<" Timecrossing "<<GetTimeHit(i)<<" ";
4a838220 391 cout<<"Number of pads "<<GetNPads(i)<<endl;
029912b7 392 }
393 cout<<"Overlapping index "<<GetOverlap(i)<<endl;
394 }
395}
396
4a838220 397Double_t AliL3ModelTrack::GetParSigmaY2(Int_t row)
95a00d93 398{
4a838220 399 //Calculate the expected cluster width, based on the track parameters and drift distance.
400
92a876e2 401 Float_t pad,time;
402 if(!GetTime(row,time) || !GetPad(row,pad))
403 return -1;
404
405 Float_t xyz[3];
406 Int_t sector,padrow;
4a838220 407 AliL3Transform::Slice2Sector(fSlice,row,sector,padrow);
408 AliL3Transform::Raw2Local(xyz,sector,padrow,pad,time);
92a876e2 409
4a838220 410 //Calculate the drift length:
411 Double_t drift;
412 if(xyz[2] > 0)
413 drift = AliL3Transform::GetZLength() - xyz[2];
414 else
415 drift = AliL3Transform::GetZLength() + xyz[2];
416
417 Double_t prf = AliL3Transform::GetPRFSigma(fPatch);
418 Double_t diffT = AliL3Transform::GetDiffT();
419 Double_t padlength = AliL3Transform::GetPadLength(fPatch);
420 Double_t anode = AliL3Transform::GetAnodeWireSpacing();
421 Double_t beta = GetCrossingAngle(row);
95a00d93 422
4a838220 423 Double_t sigmaY2 = prf*prf + diffT*diffT*drift + padlength*padlength*tan(beta)*tan(beta)/12 + anode*anode*pow( tan(beta)-0.15, 2)/12;
95a00d93 424
4a838220 425 //Convert back to raw coordinates.
426 sigmaY2 = sigmaY2/pow(AliL3Transform::GetPadPitchWidth(fPatch),2);
427 return sigmaY2;
95a00d93 428}
429
4a838220 430Double_t AliL3ModelTrack::GetParSigmaZ2(Int_t row)
95a00d93 431{
4a838220 432 //Calculate the expected cluster width, based on the track parameters and drift distance.
95a00d93 433
92a876e2 434 Float_t pad,time;
435 if(!GetTime(row,time) || !GetPad(row,pad))
436 return -1;
437
438 Float_t xyz[3];
439 Int_t sector,padrow;
4a838220 440 AliL3Transform::Slice2Sector(fSlice,row,sector,padrow);
441 AliL3Transform::Raw2Local(xyz,sector,padrow,pad,time);
92a876e2 442
4a838220 443 //Calculate the drift length:
444 Double_t drift;
445 if(xyz[2] > 0)
446 drift = AliL3Transform::GetZLength() - xyz[2];
447 else
448 drift = AliL3Transform::GetZLength() + xyz[2];
449
450 Double_t sigma0 = AliL3Transform::GetTimeSigma();
451 Double_t diffL = AliL3Transform::GetDiffL();
452 Double_t padlength = AliL3Transform::GetPadLength(fPatch);
453 Double_t tanl = GetTgl();
454
455 Double_t sigmaZ2 = sigma0*sigma0 + diffL*diffL*drift + padlength*padlength * tanl*tanl/12;
456
457 //Convert back to raw coodinates:
458 sigmaZ2 = sigmaZ2/pow(AliL3Transform::GetZWidth(),2);
459 return sigmaZ2;
95a00d93 460
95a00d93 461}
8487f697 462
463void AliL3ModelTrack::AssignTrackID(Float_t wrong=0.10)
464{
465 //Assign a track ID to the track, corresponding to the MC TParticle ID.
466 //Can only be done if you compiled with do_mc flag, of course.
467 //The function loops over the assigned clusters, and finds the label (ID)
468 //of each clusters, and assigns the ID with the most hits to the track.
469 //If there are more than wrong% clusters of a different ID, the track is
470 //considered to be fake, and label will be assigned as negative.
471
472#ifdef do_mc
473 Int_t *lb = new Int_t[GetNClusters()];
474 Int_t *mx = new Int_t[GetNClusters()];
475
476 Int_t i,j;
477 for(Int_t i=0; i<GetNClusters(); i++)
478 lb[i]=mx[i]=0;
479
480 Int_t lab=123456789;
481
482 for(i=0; i<GetNClusters(); i++)
483 {
484 lab = abs(GetTrackID(i,0));
485 for (j=0; j<GetNClusters(); j++)
486 if (lb[j]==lab || mx[j]==0) break;
487 lb[j]=lab;
488 (mx[j])++;
489 }
490
491 Int_t max=0;
492 for (i=0; i<GetNClusters(); i++)
493 {
494 if(mx[i] > max)
495 {
496 max=mx[i];
497 lab=lb[i];
498 }
499 }
500
501 for (i=0; i<GetNClusters(); i++)
502 {
503 if(abs(GetTrackID(i,1)) == lab ||
504 abs(GetTrackID(i,2)) == lab)
505 max++;
506 }
507
508 if ((1.- Float_t(max)/GetNClusters()) > wrong) lab=-lab;
509
510 SetLabel(lab);
511
512 delete[] lb;
513 delete[] mx;
514 return;
515#endif
516 cerr<<"AliL3ModelTrack::AssignTrackID : Compile with do_mc flag"<<endl;
517}