Classes moved to STEERBase.
[u/mrichter/AliRoot.git] / TRD / AliTRDmcmTracklet.cxx
CommitLineData
0ee00e25 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
e3b2b5e5 16///////////////////////////////////////////////////////
17// //
18// //
19// Tracklet object created in the local tracking //
20// //
21// //
22///////////////////////////////////////////////////////
23
0ee00e25 24#include <TGraph.h>
25#include <TMath.h>
26#include <TF1.h>
27
6d50f529 28#include "AliLog.h"
29
0ee00e25 30#include "AliTRDcalibDB.h"
31#include "AliTRDCommonParam.h"
32#include "AliTRDpadPlane.h"
33#include "AliTRDgeometry.h"
0ee00e25 34#include "AliTRDmcmTracklet.h"
35
36ClassImp(AliTRDmcmTracklet)
37
38//_____________________________________________________________________________
39AliTRDmcmTracklet::AliTRDmcmTracklet()
6d50f529 40 :TObject()
f162af62 41 ,fGeo(0)
6d50f529 42 ,fDetector(-1)
43 ,fRow(-1)
44 ,fTrackLabel(-1)
45 ,fNclusters(0)
46 ,fN(0)
47 ,fGPos(0)
48 ,fGAmp(0)
49 ,fTime0(0)
50 ,fRowz(0)
51 ,fSlope(0)
52 ,fOffset(0)
53 ,fPt(0)
54 ,fdQdl(0)
0ee00e25 55{
0ee00e25 56 //
57 // AliTRDmcmTracklet default constructor
58 //
59
0ee00e25 60 for (Int_t time = 0; time < kNtimeBins; time++) {
61 for (Int_t icl = 0; icl < kNclsPads; icl++) {
62 fADC[time][icl] = 0;
63 }
64 for (Int_t it = 0; it < kNdict; it++) {
65 fTrack[time][it] = -1;
66 }
67 fTime[time] = 0;
68 fCol[time] = 0;
69 }
70
f162af62 71 fGeo = new AliTRDgeometry();
72
0ee00e25 73}
74
75//_____________________________________________________________________________
76AliTRDmcmTracklet::AliTRDmcmTracklet(Int_t det, Int_t row, Int_t n)
6d50f529 77 :TObject()
f162af62 78 ,fGeo(0)
6d50f529 79 ,fDetector(det)
80 ,fRow(row)
81 ,fTrackLabel(-1)
82 ,fNclusters(0)
83 ,fN(n)
84 ,fGPos(0)
85 ,fGAmp(0)
86 ,fTime0(0)
87 ,fRowz(0)
88 ,fSlope(0)
89 ,fOffset(0)
90 ,fPt(0)
91 ,fdQdl(0)
0ee00e25 92{
0ee00e25 93 //
94 // AliTRDmcmTracklet default constructor
95 //
96
0ee00e25 97 for (Int_t time = 0; time < kNtimeBins; time++) {
98 for (Int_t icl = 0; icl < kNclsPads; icl++) {
99 fADC[time][icl] = 0;
100 }
101 for (Int_t it = 0; it < kNdict; it++) {
102 fTrack[time][it] = -1;
103 }
104 fTime[time] = 0;
105 fCol[time] = 0;
106 }
107
0ee00e25 108 fGPos = new TGraph(0);
109 fGAmp = new TGraph(0);
110
f162af62 111 fGeo = new AliTRDgeometry();
112
6d50f529 113}
114
115//_____________________________________________________________________________
116AliTRDmcmTracklet::AliTRDmcmTracklet(const AliTRDmcmTracklet &t)
117 :TObject(t)
f162af62 118 ,fGeo(0)
6d50f529 119 ,fDetector(t.fDetector)
120 ,fRow(t.fRow)
121 ,fTrackLabel(t.fTrackLabel)
122 ,fNclusters(t.fNclusters)
123 ,fN(t.fN)
124 ,fGPos(NULL)
125 ,fGAmp(NULL)
126 ,fTime0(t.fTime0)
127 ,fRowz(t.fRowz)
128 ,fSlope(t.fSlope)
129 ,fOffset(t.fOffset)
130 ,fPt(t.fPt)
131 ,fdQdl(t.fdQdl)
132{
133 //
134 // AliTRDmcmTracklet copy constructor
135 //
136
137 for (Int_t time = 0; time < kNtimeBins; time++) {
138 for (Int_t icl = 0; icl < kNclsPads; icl++) {
139 ((AliTRDmcmTracklet &) t).fADC[time][icl] = 0;
140 }
141 for (Int_t it = 0; it < kNdict; it++) {
142 ((AliTRDmcmTracklet &) t).fTrack[time][it] = -1;
143 }
144 ((AliTRDmcmTracklet &) t).fTime[time] = 0;
145 ((AliTRDmcmTracklet &) t).fCol[time] = 0;
146 }
0ee00e25 147
f162af62 148 if (fGeo) {
149 delete fGeo;
150 }
151 fGeo = new AliTRDgeometry();
152
0ee00e25 153}
154
155//_____________________________________________________________________________
156AliTRDmcmTracklet::~AliTRDmcmTracklet()
157{
0ee00e25 158 //
159 // AliTRDmcmTracklet destructor
160 //
3fcb9908 161
162 if (fGPos != 0) delete fGPos;
163 if (fGAmp != 0) delete fGAmp;
164
f162af62 165 if (fGeo) {
166 delete fGeo;
167 }
168
e3b2b5e5 169}
170
171//_____________________________________________________________________________
172AliTRDmcmTracklet &AliTRDmcmTracklet::operator=(const AliTRDmcmTracklet &t)
173{
174 //
6d50f529 175 // Assignment operator
e3b2b5e5 176 //
177
178 if (this != &t) ((AliTRDmcmTracklet &) t).Copy(*this);
179 return *this;
180
181}
182
183//_____________________________________________________________________________
184void AliTRDmcmTracklet::Copy(TObject &t) const
185{
186 //
6d50f529 187 // Copy function
e3b2b5e5 188 //
189
190 ((AliTRDmcmTracklet &) t).fDetector = fDetector;
191 ((AliTRDmcmTracklet &) t).fRow = fRow;
192 ((AliTRDmcmTracklet &) t).fTrackLabel = fTrackLabel;
193 ((AliTRDmcmTracklet &) t).fNclusters = fNclusters;
194 ((AliTRDmcmTracklet &) t).fN = fN;
195 ((AliTRDmcmTracklet &) t).fGPos = NULL;
196 ((AliTRDmcmTracklet &) t).fGAmp = NULL;
197 ((AliTRDmcmTracklet &) t).fTime0 = fTime0;
198 ((AliTRDmcmTracklet &) t).fRowz = fRowz;
199 ((AliTRDmcmTracklet &) t).fSlope = fSlope;
200 ((AliTRDmcmTracklet &) t).fOffset = fOffset;
201 ((AliTRDmcmTracklet &) t).fPt = fPt;
202 ((AliTRDmcmTracklet &) t).fdQdl = fdQdl;
203
204 for (Int_t time = 0; time < kNtimeBins; time++) {
205 for (Int_t icl = 0; icl < kNclsPads; icl++) {
206 ((AliTRDmcmTracklet &) t).fADC[time][icl] = 0;
207 }
208 for (Int_t it = 0; it < kNdict; it++) {
209 ((AliTRDmcmTracklet &) t).fTrack[time][it] = -1;
210 }
211 ((AliTRDmcmTracklet &) t).fTime[time] = 0;
212 ((AliTRDmcmTracklet &) t).fCol[time] = 0;
213 }
214
215}
216
217//_____________________________________________________________________________
3fcb9908 218void AliTRDmcmTracklet::Reset()
e3b2b5e5 219{
220 //
3fcb9908 221 // Reset the tracklet information
222 //
223
224 fDetector = -1;
225 fRow = -1;
226
227 for (Int_t time = 0; time < kNtimeBins; time++) {
228 for (Int_t icl = 0; icl < kNclsPads; icl++) {
229 fADC[time][icl] = 0;
230 }
231 for (Int_t it = 0; it < kNdict; it++) {
232 fTrack[time][it] = -1;
233 }
234 fTime[time] = 0;
235 fCol[time] = 0;
236 }
237
238 fNclusters = 0;
239 fN = 0;
240 fTrackLabel = -1;
241
242 fGPos->Set(0);
243 fGAmp->Set(0);
244
245 fSlope = 0.0;
246 fOffset = 0.0;
247 fTime0 = 0.0;
248 fRowz = 0.0;
249 fPt = 0.0;
250 fdQdl = 0.0;
251
252}
253
254//_____________________________________________________________________________
255void AliTRDmcmTracklet::AddCluster(Int_t icol, Int_t itb, Float_t *adc, Int_t *track)
256{
e3b2b5e5 257 //
3fcb9908 258 // Add a cluster to the tracklet
e3b2b5e5 259 //
260
261 if (fNclusters >= kNtimeBins) return;
262
263 for (Int_t icl = 0; icl < kNclsPads; icl++) {
e3b2b5e5 264 fADC[fNclusters][icl] = adc[icl];
265 }
266
267 fTrack[fNclusters][0] = track[0];
268 fTrack[fNclusters][1] = track[1];
269 fTrack[fNclusters][2] = track[2];
270 fTime[fNclusters] = itb;
271 fCol[fNclusters] = icol;
272
273 fNclusters++;
0ee00e25 274
275}
276
277//_____________________________________________________________________________
278void AliTRDmcmTracklet::MakeTrackletGraph(AliTRDgeometry *geo, Float_t field)
279{
0ee00e25 280 //
281 // Tracklet graph of positions (global coordinates, rotated [cm])
282 //
283
284 if (!geo) {
6d50f529 285 AliError("No geometry.");
0ee00e25 286 return;
287 }
288
f162af62 289 if (!AliTRDCommonParam::Instance()) {
6d50f529 290 AliError("No common parameters.");
0ee00e25 291 return;
292 }
293
294 AliTRDcalibDB* calibration = AliTRDcalibDB::Instance();
6d50f529 295 if (!calibration) {
296 AliError("No instance of AliTRDcalibDB.");
0ee00e25 297 return;
298 }
299
300 Int_t iplan, icham;
301
302 iplan = geo->GetPlane(fDetector);
303 icham = geo->GetChamber(fDetector);
304
f162af62 305 AliTRDpadPlane *padPlane = fGeo->GetPadPlane(iplan,icham);
0ee00e25 306
f162af62 307 Float_t samplFreq = AliTRDCommonParam::Instance()->GetSamplingFrequency();
0ee00e25 308
309 Int_t time, col;
310 Float_t amp[3];
e3b2b5e5 311 Float_t xpos, ypos, xzero, colSize, timeBinSize;
312 Float_t vDrift, omegaTau, lorentzAngle, thetaSlope;
313 Float_t tiltingAngle;
0ee00e25 314 Int_t npg = 0;
315
316 for (Int_t icl = 0; icl < fNclusters; icl++) {
317
318 time = GetClusterTime(icl);
319
320 amp[0] = GetClusterADC(icl)[0];
321 amp[1] = GetClusterADC(icl)[1];
322 amp[2] = GetClusterADC(icl)[2];
323
324 col = GetClusterCol(icl);
325
326 if (amp[0] < 0.0 || amp[1] < 0.0 || amp[2] < 0.0) continue;
327
e3b2b5e5 328 ypos = GetClusY(amp,iplan);
0ee00e25 329
e3b2b5e5 330 colSize = padPlane->GetColSize(col);
0ee00e25 331 vDrift = calibration->GetVdrift(fDetector,col,fRow);
e3b2b5e5 332 timeBinSize = vDrift/samplFreq;
0ee00e25 333
334 // From v4-03-Release to HEAD28Mar06 the sign has changed from "-" to "+"
335 // due to a change in the digitizer
6d50f529 336 omegaTau = TMath::Sign(1.0,(Double_t)field)*GetOmegaTau(vDrift,TMath::Abs(field));
e3b2b5e5 337 lorentzAngle = TMath::ATan(omegaTau)*180.0/TMath::Pi();
0ee00e25 338
e3b2b5e5 339 xpos = (time+0.5) * timeBinSize;
340 xpos = geo->GetTime0(iplan) - xpos;
0ee00e25 341
e3b2b5e5 342 ypos = padPlane->GetColPos(col) - (ypos + 0.5) * colSize;
0ee00e25 343
344 // ExB correction
e3b2b5e5 345 xzero = geo->GetTime0(iplan);
346 ypos = ypos + (xpos-xzero) * omegaTau;
0ee00e25 347
348 // tilted pads correction
e3b2b5e5 349 thetaSlope = - padPlane->GetRowPos(fRow)/geo->GetTime0(iplan);
350 tiltingAngle = padPlane->GetTiltingAngle()/180.0*TMath::Pi();
351 ypos = ypos - (xpos-xzero) * thetaSlope * TMath::Sin(tiltingAngle);
0ee00e25 352
e3b2b5e5 353 fGPos->SetPoint(npg,(Double_t)xpos,(Double_t)ypos);
0ee00e25 354 npg++;
355
356 }
357
358 fGPos->Set(npg);
359
360 fTime0 = geo->GetTime0(iplan) - AliTRDgeometry::CdrHght() - 0.5*AliTRDgeometry::CamHght();
9e79a757 361 fRowz = padPlane->GetRowPos(fRow) - padPlane->GetRowSize(fRow)/2.0;
0ee00e25 362
6d50f529 363 Double_t xMin = 0;
364 Double_t xMax = 0;
365 Double_t x, y;
366 fGPos->GetPoint(0 ,x,y);
367 xMax = x + 0.1;
368 fGPos->GetPoint(npg-1,x,y);
369 xMin = x - 0.1;
0ee00e25 370
371 TF1 *line = new TF1("line","[0]+x*[1]",xMin,xMax);
372 fGPos->Fit(line,"WRQ0");
373
374 fOffset = line->Eval(fTime0);
375 fSlope = TMath::ATan(line->GetParameter(1))*180.0/TMath::Pi();
376
377 line->Delete();
378
e3b2b5e5 379 Float_t fx = fTime0;
380 Float_t fy = fOffset;
0ee00e25 381
e3b2b5e5 382 Float_t infSlope = TMath::ATan(fy/fx)/TMath::Pi()*180.0;
6d50f529 383 Float_t alpha = fSlope - infSlope;
384 Float_t r = TMath::Sqrt(fx*fx + fy*fy)/(2.0*TMath::Sin(alpha/180.0*TMath::Pi()));
0ee00e25 385
e3b2b5e5 386 fPt = 0.3 * field * 0.01 * r;
0ee00e25 387
388 return;
389
390}
391
392//_____________________________________________________________________________
393void AliTRDmcmTracklet::MakeClusAmpGraph()
394{
395 //
396 // Tracklet graph of cluster charges
397 //
398
6d50f529 399 Int_t time;
0ee00e25 400 Float_t amp[3];
6d50f529 401 Int_t npg = 0;
402
0ee00e25 403 fdQdl = 0.0;
404 for (Int_t icl = 0; icl < fNclusters; icl++) {
405
406 time = GetClusterTime(icl);
407
408 amp[0] = GetClusterADC(icl)[0];
409 amp[1] = GetClusterADC(icl)[1];
410 amp[2] = GetClusterADC(icl)[2];
411
412 fGAmp->SetPoint(npg,(Double_t)(time+0.5),(Double_t)(amp[0]+amp[1]+amp[2]));
413 npg++;
414
415 fdQdl += amp[0]+amp[1]+amp[2];
416
417 }
418
419 fGAmp->Set(npg);
420
421 fdQdl /= (Float_t)npg;
422
423 return;
424
425}
426
427//_____________________________________________________________________________
e3b2b5e5 428Float_t AliTRDmcmTracklet::GetClusY(Float_t *adc, Int_t pla) const
0ee00e25 429{
430 //
431 // Cluster position in the phi direction in pad units (relative to the pad border)
432 //
433
e3b2b5e5 434 Float_t ypos = 0.0;
0ee00e25 435
e3b2b5e5 436 Float_t a0 = adc[0];
437 Float_t a1 = adc[1];
438 Float_t a2 = adc[2];
0ee00e25 439
e3b2b5e5 440 Float_t t1, t2, w1 ,w2;
0ee00e25 441
e3b2b5e5 442 Float_t w = 1.0; // pad units
0ee00e25 443
e3b2b5e5 444 Float_t sigma = 0.0;
0ee00e25 445
446 switch(pla) {
447 case 0:
e3b2b5e5 448 sigma = 0.515; break;
0ee00e25 449 case 1:
e3b2b5e5 450 sigma = 0.501; break;
0ee00e25 451 case 2:
e3b2b5e5 452 sigma = 0.491; break;
0ee00e25 453 case 3:
e3b2b5e5 454 sigma = 0.481; break;
0ee00e25 455 case 4:
e3b2b5e5 456 sigma = 0.471; break;
0ee00e25 457 case 5:
e3b2b5e5 458 sigma = 0.463; break;
0ee00e25 459 default:
6d50f529 460 AliError("Wrong plane number.");
0ee00e25 461 return 0.0;
462 }
463
e3b2b5e5 464 sigma *= w;
0ee00e25 465
e3b2b5e5 466 t1 = 0.0;
467 w1 = 0.0;
468 if( a0 > 0 ) {
469 w1 = a0*a0;
470 //w1 = a0;
471 t1 = w1*((sigma*sigma)/w*TMath::Log(a1/a0)-0.5*w);
0ee00e25 472 }
e3b2b5e5 473 t2 = 0.0;
474 w2 = 0.0;
475 if( a2 > 0 ) {
476 w2 = a2*a2;
477 //w2 = a2;
478 t2 = w2*((sigma*sigma)/w*TMath::Log(a2/a1)+0.5*w);
0ee00e25 479 }
480
e3b2b5e5 481 ypos = w*(t1+t2)/(w1+w2); // range: -0.5*w ... +0.5*w
0ee00e25 482
e3b2b5e5 483 return ypos;
0ee00e25 484
485}
486
487//_____________________________________________________________________________
488void AliTRDmcmTracklet::CookLabel(Float_t frac)
489{
490 //
491 // Cook the track label from cluster labels
492 //
493
494 const Int_t kMaxTracks = 10;
495 Int_t trackLabel[kMaxTracks];
496 Int_t trackCount[kMaxTracks];
497 for (Int_t it = 0; it < kMaxTracks; it++) {
498 trackLabel[it] = -1;
499 trackCount[it] = 0;
500 }
501
502 Bool_t counted;
503 Int_t label, nTracks = 0;
504 for (Int_t icl = 0; icl < fNclusters; icl++) {
505
506 for (Int_t id = 0; id < kNdict; id++) {
507
508 if (fTrack[icl][id] == -1) continue;
509
510 label = fTrack[icl][id];
511
512 counted = kFALSE;
513 for (Int_t it = 0; it < nTracks; it++) {
514 if (label == trackLabel[it]) {
515 trackCount[it]++;
516 counted = kTRUE;
517 break;
518 }
519 }
520 if (!counted) {
521 trackLabel[nTracks] = label;
522 trackCount[nTracks]++;
523 nTracks++;
524 if (nTracks == kMaxTracks) {
6d50f529 525 AliWarning("Too many tracks for this tracklet.");
0ee00e25 526 nTracks--;
527 break;
528 }
529 }
530
531 }
532
533 }
534
535 for (Int_t it = 0; it < kMaxTracks; it++) {
536 if (trackCount[it] >= (Int_t)(frac*fNclusters)) {
537 fTrackLabel = trackLabel[it];
538 break;
539 }
540 }
541
542}
543
544//_____________________________________________________________________________
e3b2b5e5 545Float_t AliTRDmcmTracklet::GetOmegaTau(Float_t vdrift, Float_t field) const
0ee00e25 546{
547 //
548 // Returns omega*tau (tan(Lorentz-angle)) for a given drift velocity <vd>
549 // and a B-field <b> for Xe/CO2 (15%).
550 // The values are according to a GARFIELD simulation.
551 //
552
553 //
554 // Copy of the "AliTRDcalibDB" function, taking as argument the magnetic field too
555 //
556
557 const Int_t kNb = 5;
558 Float_t p0[kNb] = { 0.004810, 0.007412, 0.010252, 0.013409, 0.016888 };
559 Float_t p1[kNb] = { 0.054875, 0.081534, 0.107333, 0.131983, 0.155455 };
560 Float_t p2[kNb] = { -0.008682, -0.012896, -0.016987, -0.020880, -0.024623 };
561 Float_t p3[kNb] = { 0.000155, 0.000238, 0.000330, 0.000428, 0.000541 };
562
563 Int_t ib = ((Int_t) (10 * (field - 0.15)));
564 ib = TMath::Max( 0,ib);
565 ib = TMath::Min(kNb,ib);
566
567 Float_t alphaL = p0[ib]
6d50f529 568 + p1[ib] * vdrift
569 + p2[ib] * vdrift*vdrift
570 + p3[ib] * vdrift*vdrift*vdrift;
0ee00e25 571
572 return TMath::Tan(alphaL);
573
574}