]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/AliTRDclusterizerV1.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
SetSolenoidField method to set the L3 field strength. 2 kG is default.
[u/mrichter/AliRoot.git] / TRD / AliTRDclusterizerV1.cxx
CommitLineData
f7336fa3 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
16/*
17$Log$
26edf6a4 18Revision 1.11 2001/05/21 16:45:47 hristov
19Last minute changes (C.Blume)
20
db30bf0f 21Revision 1.10 2001/05/07 08:06:44 cblume
22Speedup of the code. Create only AliTRDcluster
23
3e1a3ad8 24Revision 1.9 2000/11/01 14:53:20 cblume
25Merge with TRD-develop
26
793ff80c 27Revision 1.1.4.5 2000/10/15 23:40:01 cblume
28Remove AliTRDconst
29
30Revision 1.1.4.4 2000/10/06 16:49:46 cblume
31Made Getters const
32
33Revision 1.1.4.3 2000/10/04 16:34:58 cblume
34Replace include files by forward declarations
35
36Revision 1.1.4.2 2000/09/22 14:49:49 cblume
37Adapted to tracking code
38
39Revision 1.8 2000/10/02 21:28:19 fca
40Removal of useless dependecies via forward declarations
41
94de3818 42Revision 1.7 2000/06/27 13:08:50 cblume
43Changed to Copy(TObject &A) to appease the HP-compiler
44
43da34c0 45Revision 1.6 2000/06/09 11:10:07 cblume
46Compiler warnings and coding conventions, next round
47
dd9a6ee3 48Revision 1.5 2000/06/08 18:32:58 cblume
49Make code compliant to coding conventions
50
8230f242 51Revision 1.4 2000/06/07 16:27:01 cblume
52Try to remove compiler warnings on Sun and HP
53
9d0b222b 54Revision 1.3 2000/05/08 16:17:27 cblume
55Merge TRD-develop
56
6f1e466d 57Revision 1.1.4.1 2000/05/08 15:09:01 cblume
58Introduce AliTRDdigitsManager
59
c0dd96c3 60Revision 1.1 2000/02/28 18:58:54 cblume
61Add new TRD classes
62
f7336fa3 63*/
64
65///////////////////////////////////////////////////////////////////////////////
66// //
67// TRD cluster finder for the slow simulator.
68// //
69///////////////////////////////////////////////////////////////////////////////
70
71#include <TF1.h>
94de3818 72#include <TTree.h>
793ff80c 73#include <TH1.h>
f7336fa3 74
793ff80c 75#include "AliRun.h"
76
77#include "AliTRD.h"
f7336fa3 78#include "AliTRDclusterizerV1.h"
79#include "AliTRDmatrix.h"
80#include "AliTRDgeometry.h"
81#include "AliTRDdigitizer.h"
6f1e466d 82#include "AliTRDdataArrayF.h"
793ff80c 83#include "AliTRDdataArrayI.h"
84#include "AliTRDdigitsManager.h"
f7336fa3 85
86ClassImp(AliTRDclusterizerV1)
87
88//_____________________________________________________________________________
89AliTRDclusterizerV1::AliTRDclusterizerV1():AliTRDclusterizer()
90{
91 //
92 // AliTRDclusterizerV1 default constructor
93 //
94
6f1e466d 95 fDigitsManager = NULL;
f7336fa3 96
3e1a3ad8 97 fClusMaxThresh = 0;
98 fClusSigThresh = 0;
99
db30bf0f 100 fUseLUT = kFALSE;
101
f7336fa3 102}
103
104//_____________________________________________________________________________
105AliTRDclusterizerV1::AliTRDclusterizerV1(const Text_t* name, const Text_t* title)
106 :AliTRDclusterizer(name,title)
107{
108 //
109 // AliTRDclusterizerV1 default constructor
110 //
111
6f1e466d 112 fDigitsManager = new AliTRDdigitsManager();
f7336fa3 113
114 Init();
115
116}
117
8230f242 118//_____________________________________________________________________________
dd9a6ee3 119AliTRDclusterizerV1::AliTRDclusterizerV1(const AliTRDclusterizerV1 &c)
8230f242 120{
121 //
122 // AliTRDclusterizerV1 copy constructor
123 //
124
dd9a6ee3 125 ((AliTRDclusterizerV1 &) c).Copy(*this);
8230f242 126
127}
128
f7336fa3 129//_____________________________________________________________________________
130AliTRDclusterizerV1::~AliTRDclusterizerV1()
131{
8230f242 132 //
133 // AliTRDclusterizerV1 destructor
134 //
f7336fa3 135
6f1e466d 136 if (fDigitsManager) {
137 delete fDigitsManager;
f7336fa3 138 }
139
140}
141
dd9a6ee3 142//_____________________________________________________________________________
143AliTRDclusterizerV1 &AliTRDclusterizerV1::operator=(const AliTRDclusterizerV1 &c)
144{
145 //
146 // Assignment operator
147 //
148
149 if (this != &c) ((AliTRDclusterizerV1 &) c).Copy(*this);
150 return *this;
151
152}
153
8230f242 154//_____________________________________________________________________________
43da34c0 155void AliTRDclusterizerV1::Copy(TObject &c)
8230f242 156{
157 //
158 // Copy function
159 //
160
db30bf0f 161 ((AliTRDclusterizerV1 &) c).fUseLUT = fUseLUT;
43da34c0 162 ((AliTRDclusterizerV1 &) c).fClusMaxThresh = fClusMaxThresh;
163 ((AliTRDclusterizerV1 &) c).fClusSigThresh = fClusSigThresh;
43da34c0 164 ((AliTRDclusterizerV1 &) c).fDigitsManager = NULL;
db30bf0f 165 for (Int_t ilut = 0; ilut < kNlut; ilut++) {
166 ((AliTRDclusterizerV1 &) c).fLUT[ilut] = fLUT[ilut];
167 }
8230f242 168
169 AliTRDclusterizer::Copy(c);
170
171}
172
f7336fa3 173//_____________________________________________________________________________
174void AliTRDclusterizerV1::Init()
175{
176 //
177 // Initializes the cluster finder
178 //
179
180 // The default parameter for the clustering
3e1a3ad8 181 fClusMaxThresh = 3;
182 fClusSigThresh = 1;
f7336fa3 183
db30bf0f 184 // Use the lookup table for the position determination
185 fUseLUT = kTRUE;
186
187 // The lookup table from Bogdan
188 Float_t lut[128] = {
189 0.0068, 0.0198, 0.0318, 0.0432, 0.0538, 0.0642, 0.0742, 0.0838,
190 0.0932, 0.1023, 0.1107, 0.1187, 0.1268, 0.1347, 0.1423, 0.1493,
191 0.1562, 0.1632, 0.1698, 0.1762, 0.1828, 0.1887, 0.1947, 0.2002,
192 0.2062, 0.2118, 0.2173, 0.2222, 0.2278, 0.2327, 0.2377, 0.2428,
193 0.2473, 0.2522, 0.2567, 0.2612, 0.2657, 0.2697, 0.2743, 0.2783,
194 0.2822, 0.2862, 0.2903, 0.2943, 0.2982, 0.3018, 0.3058, 0.3092,
195 0.3128, 0.3167, 0.3203, 0.3237, 0.3268, 0.3302, 0.3338, 0.3368,
196 0.3402, 0.3433, 0.3462, 0.3492, 0.3528, 0.3557, 0.3587, 0.3613,
197 0.3643, 0.3672, 0.3702, 0.3728, 0.3758, 0.3783, 0.3812, 0.3837,
198 0.3862, 0.3887, 0.3918, 0.3943, 0.3968, 0.3993, 0.4017, 0.4042,
199 0.4067, 0.4087, 0.4112, 0.4137, 0.4157, 0.4182, 0.4207, 0.4227,
200 0.4252, 0.4272, 0.4293, 0.4317, 0.4338, 0.4358, 0.4383, 0.4403,
201 0.4423, 0.4442, 0.4462, 0.4482, 0.4502, 0.4523, 0.4543, 0.4563,
202 0.4582, 0.4602, 0.4622, 0.4638, 0.4658, 0.4678, 0.4697, 0.4712,
203 0.4733, 0.4753, 0.4767, 0.4787, 0.4803, 0.4823, 0.4837, 0.4857,
204 0.4873, 0.4888, 0.4908, 0.4922, 0.4942, 0.4958, 0.4972, 0.4988
205 };
206 for (Int_t ilut = 0; ilut < kNlut; ilut++) {
207 fLUT[ilut] = lut[ilut];
208 }
209
f7336fa3 210}
211
212//_____________________________________________________________________________
213Bool_t AliTRDclusterizerV1::ReadDigits()
214{
215 //
216 // Reads the digits arrays from the input aliroot file
217 //
218
219 if (!fInputFile) {
220 printf("AliTRDclusterizerV1::ReadDigits -- ");
221 printf("No input file open\n");
222 return kFALSE;
223 }
224
f7336fa3 225 // Read in the digit arrays
6f1e466d 226 return (fDigitsManager->ReadDigits());
f7336fa3 227
228}
229
230//_____________________________________________________________________________
793ff80c 231Bool_t AliTRDclusterizerV1::MakeClusters()
f7336fa3 232{
233 //
234 // Generates the cluster.
235 //
236
237 Int_t row, col, time;
238
3e1a3ad8 239 if (fTRD->IsVersion() != 1) {
f7336fa3 240 printf("AliTRDclusterizerV1::MakeCluster -- ");
241 printf("TRD must be version 1 (slow simulator).\n");
242 return kFALSE;
243 }
244
245 // Get the geometry
3e1a3ad8 246 AliTRDgeometry *geo = fTRD->GetGeometry();
f7336fa3 247
248 printf("AliTRDclusterizerV1::MakeCluster -- ");
249 printf("Start creating clusters.\n");
250
8230f242 251 AliTRDdataArrayI *digits;
793ff80c 252 AliTRDdataArrayI *track0;
253 AliTRDdataArrayI *track1;
254 AliTRDdataArrayI *track2;
f7336fa3 255
3e1a3ad8 256 // Threshold value for the maximum
257 Int_t maxThresh = fClusMaxThresh;
258 // Threshold value for the digit signal
259 Int_t sigThresh = fClusSigThresh;
f7336fa3 260
261 // Iteration limit for unfolding procedure
8230f242 262 const Float_t kEpsilon = 0.01;
f7336fa3 263
8230f242 264 const Int_t kNclus = 3;
265 const Int_t kNsig = 5;
3e1a3ad8 266 const Int_t kNtrack = 3 * kNclus;
267
db30bf0f 268 // For the LUT
269 const Float_t kLUTmin = 0.106113;
270 const Float_t kLUTmax = 0.995415;
271
272 Int_t iType = 0;
273 Int_t iUnfold = 0;
274
275 Float_t ratioLeft = 1.0;
276 Float_t ratioRight = 1.0;
277
3e1a3ad8 278 Float_t padSignal[kNsig];
279 Float_t clusterSignal[kNclus];
280 Float_t clusterPads[kNclus];
281 Int_t clusterDigit[kNclus];
282 Int_t clusterTracks[kNtrack];
f7336fa3 283
284 Int_t chamBeg = 0;
793ff80c 285 Int_t chamEnd = AliTRDgeometry::Ncham();
3e1a3ad8 286 if (fTRD->GetSensChamber() >= 0) {
287 chamBeg = fTRD->GetSensChamber();
6f1e466d 288 chamEnd = chamBeg + 1;
f7336fa3 289 }
290 Int_t planBeg = 0;
793ff80c 291 Int_t planEnd = AliTRDgeometry::Nplan();
3e1a3ad8 292 if (fTRD->GetSensPlane() >= 0) {
293 planBeg = fTRD->GetSensPlane();
f7336fa3 294 planEnd = planBeg + 1;
295 }
296 Int_t sectBeg = 0;
793ff80c 297 Int_t sectEnd = AliTRDgeometry::Nsect();
f7336fa3 298
3e1a3ad8 299 // Start clustering in every chamber
f7336fa3 300 for (Int_t icham = chamBeg; icham < chamEnd; icham++) {
301 for (Int_t iplan = planBeg; iplan < planEnd; iplan++) {
302 for (Int_t isect = sectBeg; isect < sectEnd; isect++) {
303
3e1a3ad8 304 if (fTRD->GetSensSector() >= 0) {
305 Int_t sens1 = fTRD->GetSensSector();
306 Int_t sens2 = sens1 + fTRD->GetSensSectorRange();
793ff80c 307 sens2 -= ((Int_t) (sens2 / AliTRDgeometry::Nsect()))
308 * AliTRDgeometry::Nsect();
dd9a6ee3 309 if (sens1 < sens2) {
9d0b222b 310 if ((isect < sens1) || (isect >= sens2)) continue;
dd9a6ee3 311 }
312 else {
9d0b222b 313 if ((isect < sens1) && (isect >= sens2)) continue;
dd9a6ee3 314 }
9d0b222b 315 }
316
8230f242 317 Int_t idet = geo->GetDetector(iplan,icham,isect);
f7336fa3 318
db30bf0f 319 Int_t nClusters = 0;
320 Int_t nClusters2pad = 0;
321 Int_t nClusters3pad = 0;
322 Int_t nClusters4pad = 0;
323 Int_t nClusters5pad = 0;
324 Int_t nClustersLarge = 0;
3e1a3ad8 325
f7336fa3 326 printf("AliTRDclusterizerV1::MakeCluster -- ");
327 printf("Analyzing chamber %d, plane %d, sector %d.\n"
3e1a3ad8 328 ,icham,iplan,isect);
f7336fa3 329
3e1a3ad8 330 Int_t nRowMax = geo->GetRowMax(iplan,icham,isect);
331 Int_t nColMax = geo->GetColMax(iplan);
332 Int_t nTimeBefore = geo->GetTimeBefore();
333 Int_t nTimeTotal = geo->GetTimeTotal();
f7336fa3 334
3e1a3ad8 335 // Get the digits
8230f242 336 digits = fDigitsManager->GetDigits(idet);
3e1a3ad8 337 digits->Expand();
793ff80c 338 track0 = fDigitsManager->GetDictionary(idet,0);
3e1a3ad8 339 track0->Expand();
793ff80c 340 track1 = fDigitsManager->GetDictionary(idet,1);
3e1a3ad8 341 track1->Expand();
793ff80c 342 track2 = fDigitsManager->GetDictionary(idet,2);
3e1a3ad8 343 track2->Expand();
344
345 // Loop through the chamber and find the maxima
346 for ( row = 0; row < nRowMax; row++) {
347 for ( col = 2; col < nColMax; col++) {
348 for (time = 0; time < nTimeTotal; time++) {
349
350 Int_t signalL = digits->GetDataUnchecked(row,col ,time);
351 Int_t signalM = digits->GetDataUnchecked(row,col-1,time);
352 Int_t signalR = digits->GetDataUnchecked(row,col-2,time);
353
354 // Look for the maximum
db30bf0f 355 if (signalM >= maxThresh) {
356 if (((signalL >= sigThresh) &&
357 (signalL < signalM)) ||
358 ((signalR >= sigThresh) &&
359 (signalR < signalM))) {
3e1a3ad8 360 // Maximum found, mark the position by a negative signal
361 digits->SetDataUnchecked(row,col-1,time,-signalM);
362 }
363 }
364
365 }
366 }
367 }
368
369 // Now check the maxima and calculate the cluster position
370 for ( row = 0; row < nRowMax ; row++) {
db30bf0f 371 for (time = 0; time < nTimeTotal; time++) {
372 for ( col = 1; col < nColMax-1; col++) {
3e1a3ad8 373
374 // Maximum found ?
375 if (digits->GetDataUnchecked(row,col,time) < 0) {
f7336fa3 376
9d0b222b 377 Int_t iPad;
8230f242 378 for (iPad = 0; iPad < kNclus; iPad++) {
3e1a3ad8 379 Int_t iPadCol = col - 1 + iPad;
380 clusterSignal[iPad] = TMath::Abs(digits->GetDataUnchecked(row
381 ,iPadCol
382 ,time));
383 clusterDigit[iPad] = digits->GetIndexUnchecked(row,iPadCol,time);
384 clusterTracks[3*iPad ] = track0->GetDataUnchecked(row,iPadCol,time) - 1;
385 clusterTracks[3*iPad+1] = track1->GetDataUnchecked(row,iPadCol,time) - 1;
386 clusterTracks[3*iPad+2] = track2->GetDataUnchecked(row,iPadCol,time) - 1;
f7336fa3 387 }
388
db30bf0f 389 // Count the number of pads in the cluster
390 Int_t nPadCount = 0;
391 Int_t ii = 0;
392 while (TMath::Abs(digits->GetDataUnchecked(row,col-ii ,time))
393 >= sigThresh) {
394 nPadCount++;
395 ii++;
396 if (col-ii < 0) break;
397 }
398 ii = 0;
399 while (TMath::Abs(digits->GetDataUnchecked(row,col+ii+1,time))
400 >= sigThresh) {
401 nPadCount++;
402 ii++;
403 if (col+ii+1 >= nColMax) break;
404 }
405
406 nClusters++;
407 switch (nPadCount) {
408 case 2:
409 iType = 0;
410 nClusters2pad++;
411 break;
412 case 3:
413 iType = 1;
414 nClusters3pad++;
415 break;
416 case 4:
417 iType = 2;
418 nClusters4pad++;
419 break;
420 case 5:
421 iType = 3;
422 nClusters5pad++;
423 break;
424 default:
425 iType = 4;
426 nClustersLarge++;
427 break;
428 };
429
430 // Don't analyze large clusters
431 //if (iType == 4) continue;
432
433 // Look for 5 pad cluster with minimum in the middle
434 Bool_t fivePadCluster = kFALSE;
3e1a3ad8 435 if (col < nColMax-3) {
436 if (digits->GetDataUnchecked(row,col+2,time) < 0) {
db30bf0f 437 fivePadCluster = kTRUE;
438 }
439 if ((fivePadCluster) && (col < nColMax-5)) {
440 if (digits->GetDataUnchecked(row,col+4,time) >= sigThresh) {
441 fivePadCluster = kFALSE;
442 }
443 }
444 if ((fivePadCluster) && (col > 1)) {
445 if (digits->GetDataUnchecked(row,col-2,time) >= sigThresh) {
446 fivePadCluster = kFALSE;
447 }
448 }
449 }
450
451 // 5 pad cluster
452 // Modify the signal of the overlapping pad for the left part
453 // of the cluster which remains from a previous unfolding
454 if (iUnfold) {
455 clusterSignal[0] *= ratioLeft;
456 iType = 3;
457 iUnfold = 0;
458 }
459
460 // Unfold the 5 pad cluster
461 if (fivePadCluster) {
462 for (iPad = 0; iPad < kNsig; iPad++) {
463 padSignal[iPad] = TMath::Abs(digits->GetDataUnchecked(row
464 ,col-1+iPad
465 ,time));
f7336fa3 466 }
db30bf0f 467 // Unfold the two maxima and set the signal on
468 // the overlapping pad to the ratio
469 ratioRight = Unfold(kEpsilon,padSignal);
470 ratioLeft = 1.0 - ratioRight;
471 clusterSignal[2] *= ratioRight;
472 iType = 3;
473 iUnfold = 1;
f7336fa3 474 }
f7336fa3 475
3e1a3ad8 476 Float_t clusterCharge = clusterSignal[0]
477 + clusterSignal[1]
478 + clusterSignal[2];
479
db30bf0f 480 // The position of the cluster
3e1a3ad8 481 clusterPads[0] = row + 0.5;
3e1a3ad8 482 // Take the shift of the additional time bins into account
483 clusterPads[2] = time - nTimeBefore + 0.5;
484
db30bf0f 485 if (fUseLUT) {
486
487 // Calculate the position of the cluster by using the
488 // lookup table method
489 Float_t ratioLUT;
490 Float_t signLUT;
491 Float_t lut = 0.0;
492 if (clusterSignal[0] > clusterSignal[2]) {
493 ratioLUT = clusterSignal[0] / clusterSignal[1];
494 signLUT = -1.0;
495 }
496 else {
497 ratioLUT = clusterSignal[2] / clusterSignal[1];
498 signLUT = 1.0;
499 }
500 if (ratioLUT < kLUTmin) {
501 lut = 0.0;
502 }
503 else if (ratioLUT > kLUTmax) {
504 lut = 0.5;
505 }
506 else {
507 Int_t indexLUT = TMath::Nint ((kNlut-1) * (ratioLUT - kLUTmin)
508 / (kLUTmax - kLUTmin));
509 lut = fLUT[indexLUT];
510 }
511 clusterPads[1] = col + 0.5 + signLUT * lut;
512
513 }
514 else {
515
516 // Calculate the position of the cluster by using the
517 // center of gravity method
518 clusterPads[1] = col + 0.5
519 + (clusterSignal[2] - clusterSignal[0])
520 / clusterCharge;
521
522 }
523
3e1a3ad8 524 if (fVerbose) {
525 printf("-----------------------------------------------------------\n");
526 printf("Create cluster no. %d\n",nClusters);
527 printf("Position: row = %f, col = %f, time = %f\n",clusterPads[0]
528 ,clusterPads[1]
529 ,clusterPads[2]);
530 printf("Indices: %d, %d, %d\n",clusterDigit[0]
531 ,clusterDigit[1]
532 ,clusterDigit[2]);
533 printf("Total charge = %f\n",clusterCharge);
534 printf("Tracks: pad0 %d, %d, %d\n",clusterTracks[0]
535 ,clusterTracks[1]
536 ,clusterTracks[2]);
537 printf(" pad1 %d, %d, %d\n",clusterTracks[3]
538 ,clusterTracks[4]
539 ,clusterTracks[5]);
540 printf(" pad2 %d, %d, %d\n",clusterTracks[6]
541 ,clusterTracks[7]
542 ,clusterTracks[8]);
db30bf0f 543 printf("Type = %d, Number of pads = %d\n",iType,nPadCount);
f7336fa3 544 }
545
f7336fa3 546 // Add the cluster to the output array
3e1a3ad8 547 fTRD->AddCluster(clusterPads
548 ,clusterDigit
549 ,idet
550 ,clusterCharge
551 ,clusterTracks
552 ,iType);
f7336fa3 553
554 }
3e1a3ad8 555 }
556 }
557 }
f7336fa3 558
3e1a3ad8 559 // Compress the arrays
560 digits->Compress(1,0);
561 track0->Compress(1,0);
562 track1->Compress(1,0);
563 track2->Compress(1,0);
f7336fa3 564
3e1a3ad8 565 // Write the cluster and reset the array
793ff80c 566 WriteClusters(idet);
3e1a3ad8 567 fTRD->ResetRecPoints();
793ff80c 568
3e1a3ad8 569 printf("AliTRDclusterizerV1::MakeCluster -- ");
db30bf0f 570 printf("Found %d clusters in total.\n"
571 ,nClusters);
572 printf(" 2pad: %d\n",nClusters2pad);
573 printf(" 3pad: %d\n",nClusters3pad);
574 printf(" 4pad: %d\n",nClusters4pad);
575 printf(" 5pad: %d\n",nClusters5pad);
576 printf(" Large: %d\n",nClustersLarge);
f7336fa3 577
3e1a3ad8 578 }
579 }
580 }
f7336fa3 581
f7336fa3 582 printf("AliTRDclusterizerV1::MakeCluster -- ");
583 printf("Done.\n");
584
585 return kTRUE;
586
587}
588
589//_____________________________________________________________________________
590Float_t AliTRDclusterizerV1::Unfold(Float_t eps, Float_t* padSignal)
591{
592 //
593 // Method to unfold neighbouring maxima.
594 // The charge ratio on the overlapping pad is calculated
595 // until there is no more change within the range given by eps.
596 // The resulting ratio is then returned to the calling method.
597 //
598
3e1a3ad8 599 Int_t itStep = 0; // Count iteration steps
f7336fa3 600
3e1a3ad8 601 Float_t ratio = 0.5; // Start value for ratio
602 Float_t prevRatio = 0; // Store previous ratio
f7336fa3 603
3e1a3ad8 604 Float_t newLeftSignal[3] = {0}; // Array to store left cluster signal
605 Float_t newRightSignal[3] = {0}; // Array to store right cluster signal
f7336fa3 606
3e1a3ad8 607 // Start the iteration
f7336fa3 608 while ((TMath::Abs(prevRatio - ratio) > eps) && (itStep < 10)) {
609
610 itStep++;
611 prevRatio = ratio;
612
3e1a3ad8 613 // Cluster position according to charge ratio
614 Float_t maxLeft = (ratio*padSignal[2] - padSignal[0])
615 / (padSignal[0] + padSignal[1] + ratio*padSignal[2]);
616 Float_t maxRight = (padSignal[4] - (1-ratio)*padSignal[2])
617 / ((1-ratio)*padSignal[2] + padSignal[3] + padSignal[4]);
f7336fa3 618
3e1a3ad8 619 // Set cluster charge ratio
f7336fa3 620 Float_t ampLeft = padSignal[1];
621 Float_t ampRight = padSignal[3];
622
3e1a3ad8 623 // Apply pad response to parameters
624 newLeftSignal[0] = ampLeft * PadResponse(-1 - maxLeft);
625 newLeftSignal[1] = ampLeft * PadResponse( 0 - maxLeft);
626 newLeftSignal[2] = ampLeft * PadResponse( 1 - maxLeft);
f7336fa3 627
3e1a3ad8 628 newRightSignal[0] = ampRight * PadResponse(-1 - maxRight);
629 newRightSignal[1] = ampRight * PadResponse( 0 - maxRight);
630 newRightSignal[2] = ampRight * PadResponse( 1 - maxRight);
f7336fa3 631
3e1a3ad8 632 // Calculate new overlapping ratio
26edf6a4 633 ratio = TMath::Min((Float_t)1.0,newLeftSignal[2] /
db30bf0f 634 (newLeftSignal[2] + newRightSignal[0]));
f7336fa3 635
636 }
637
638 return ratio;
639
640}
641
642//_____________________________________________________________________________
643Float_t AliTRDclusterizerV1::PadResponse(Float_t x)
644{
645 //
646 // The pad response for the chevron pads.
647 // We use a simple Gaussian approximation which should be good
648 // enough for our purpose.
3e1a3ad8 649 // Updated for new PRF 1/5/01.
f7336fa3 650 //
651
652 // The parameters for the response function
3e1a3ad8 653 const Float_t kA = 0.8303;
654 const Float_t kB = -0.00392;
655 const Float_t kC = 0.472 * 0.472;
656 const Float_t kD = 2.19;
f7336fa3 657
3e1a3ad8 658 Float_t pr = kA * (kB + TMath::Exp(-TMath::Power(x*x,kD) / (2.*kC)));
f7336fa3 659
660 return (pr);
661
662}