]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/AliTRDclusterizer.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Fix bug due to shadowed variables
[u/mrichter/AliRoot.git] / TRD / AliTRDclusterizer.cxx
CommitLineData
3fe61b77 1
f7336fa3 2/**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
88cb7938 17/* $Id$ */
f7336fa3 18
19///////////////////////////////////////////////////////////////////////////////
20// //
3fe61b77 21// TRD cluster finder //
f7336fa3 22// //
23///////////////////////////////////////////////////////////////////////////////
24
3fe61b77 25#include <TF1.h>
94de3818 26#include <TTree.h>
3fe61b77 27#include <TH1.h>
793ff80c 28#include <TFile.h>
6d50f529 29#include <TObjArray.h>
f7336fa3 30
88cb7938 31#include "AliRunLoader.h"
32#include "AliLoader.h"
3fe61b77 33#include "AliRawReader.h"
6d50f529 34#include "AliLog.h"
3fe61b77 35#include "AliAlignObj.h"
88cb7938 36
f7336fa3 37#include "AliTRDclusterizer.h"
3e1a3ad8 38#include "AliTRDcluster.h"
fc546d21 39#include "AliTRDReconstructor.h"
793ff80c 40#include "AliTRDgeometry.h"
3fe61b77 41#include "AliTRDdataArrayF.h"
42#include "AliTRDdataArrayI.h"
625f5260 43#include "AliTRDdataArrayS.h"
fb44eb8e 44#include "AliTRDdataArrayDigits.h"
3fe61b77 45#include "AliTRDdigitsManager.h"
46#include "AliTRDrawData.h"
3551db50 47#include "AliTRDcalibDB.h"
fc546d21 48#include "AliTRDrecoParam.h"
b43a3e17 49#include "AliTRDCommonParam.h"
3fe61b77 50#include "AliTRDtransform.h"
51#include "AliTRDSignalIndex.h"
dfbb4bb9 52#include "AliTRDrawStreamBase.h"
3fe61b77 53#include "AliTRDfeeParam.h"
54
55#include "Cal/AliTRDCalROC.h"
56#include "Cal/AliTRDCalDet.h"
0e09df31 57#include "Cal/AliTRDCalSingleChamberStatus.h"
f7336fa3 58
59ClassImp(AliTRDclusterizer)
60
61//_____________________________________________________________________________
6d50f529 62AliTRDclusterizer::AliTRDclusterizer()
63 :TNamed()
64 ,fRunLoader(NULL)
65 ,fClusterTree(NULL)
66 ,fRecPoints(NULL)
3fe61b77 67 ,fDigitsManager(NULL)
68 ,fAddLabels(kTRUE)
69 ,fRawVersion(2)
70 ,fIndexesOut(NULL)
71 ,fIndexesMaxima(NULL)
dfbb4bb9 72 ,fTransform(new AliTRDtransform(0))
fc546d21 73 ,fLUTbin(0)
74 ,fLUT(NULL)
f7336fa3 75{
76 //
77 // AliTRDclusterizer default constructor
78 //
79
3fe61b77 80 fRawVersion = AliTRDfeeParam::Instance()->GetRAWversion();
81
f7336fa3 82}
83
84//_____________________________________________________________________________
3fe61b77 85AliTRDclusterizer::AliTRDclusterizer(const Text_t *name, const Text_t *title)
6d50f529 86 :TNamed(name,title)
87 ,fRunLoader(NULL)
88 ,fClusterTree(NULL)
89 ,fRecPoints(NULL)
3fe61b77 90 ,fDigitsManager(new AliTRDdigitsManager())
91 ,fAddLabels(kTRUE)
92 ,fRawVersion(2)
93 ,fIndexesOut(NULL)
94 ,fIndexesMaxima(NULL)
95 ,fTransform(new AliTRDtransform(0))
fc546d21 96 ,fLUTbin(0)
97 ,fLUT(NULL)
f7336fa3 98{
99 //
6d50f529 100 // AliTRDclusterizer constructor
f7336fa3 101 //
102
3fe61b77 103 fDigitsManager->CreateArrays();
104
105 fRawVersion = AliTRDfeeParam::Instance()->GetRAWversion();
106
fc546d21 107 FillLUT();
108
f7336fa3 109}
110
8230f242 111//_____________________________________________________________________________
6d50f529 112AliTRDclusterizer::AliTRDclusterizer(const AliTRDclusterizer &c)
113 :TNamed(c)
114 ,fRunLoader(NULL)
115 ,fClusterTree(NULL)
116 ,fRecPoints(NULL)
3fe61b77 117 ,fDigitsManager(NULL)
118 ,fAddLabels(kTRUE)
119 ,fRawVersion(2)
120 ,fIndexesOut(NULL)
121 ,fIndexesMaxima(NULL)
122 ,fTransform(NULL)
fc546d21 123 ,fLUTbin(0)
124 ,fLUT(0)
8230f242 125{
126 //
127 // AliTRDclusterizer copy constructor
128 //
129
fc546d21 130 FillLUT();
131
8230f242 132}
133
f7336fa3 134//_____________________________________________________________________________
135AliTRDclusterizer::~AliTRDclusterizer()
136{
8230f242 137 //
138 // AliTRDclusterizer destructor
139 //
f7336fa3 140
c10bf383 141 printf("----> AliTRDclusterizer::dtor\n");
142
3fe61b77 143 if (fRecPoints)
144 {
145 fRecPoints->Delete();
146 delete fRecPoints;
147 }
148
149 if (fDigitsManager)
150 {
151 delete fDigitsManager;
152 fDigitsManager = NULL;
153 }
154
155 if (fIndexesOut)
156 {
157 delete fIndexesOut;
158 fIndexesOut = NULL;
159 }
160
161 if (fIndexesMaxima)
162 {
163 delete fIndexesMaxima;
164 fIndexesMaxima = NULL;
165 }
166
167 if (fTransform)
168 {
169 delete fTransform;
170 fTransform = NULL;
171 }
6d50f529 172
fc546d21 173 if (fLUT)
174 {
175 delete [] fLUT;
176 fLUT = NULL;
177 }
178
f7336fa3 179}
180
8230f242 181//_____________________________________________________________________________
dd9a6ee3 182AliTRDclusterizer &AliTRDclusterizer::operator=(const AliTRDclusterizer &c)
183{
184 //
185 // Assignment operator
186 //
187
3fe61b77 188 if (this != &c)
189 {
190 ((AliTRDclusterizer &) c).Copy(*this);
191 }
dd9a6ee3 192 return *this;
193
194}
195
196//_____________________________________________________________________________
e0d47c25 197void AliTRDclusterizer::Copy(TObject &c) const
8230f242 198{
199 //
200 // Copy function
201 //
202
3fe61b77 203 ((AliTRDclusterizer &) c).fClusterTree = NULL;
204 ((AliTRDclusterizer &) c).fRecPoints = NULL;
205 ((AliTRDclusterizer &) c).fDigitsManager = NULL;
206 ((AliTRDclusterizer &) c).fAddLabels = fAddLabels;
207 ((AliTRDclusterizer &) c).fRawVersion = fRawVersion;
208 ((AliTRDclusterizer &) c).fIndexesOut = NULL;
209 ((AliTRDclusterizer &) c).fIndexesMaxima = NULL;
210 ((AliTRDclusterizer &) c).fTransform = NULL;
fc546d21 211 ((AliTRDclusterizer &) c).fLUTbin = 0;
212 ((AliTRDclusterizer &) c).fLUT = NULL;
8230f242 213
214}
215
f7336fa3 216//_____________________________________________________________________________
3e1a3ad8 217Bool_t AliTRDclusterizer::Open(const Char_t *name, Int_t nEvent)
f7336fa3 218{
219 //
3e1a3ad8 220 // Opens the AliROOT file. Output and input are in the same file
f7336fa3 221 //
6d50f529 222
e191bb57 223 TString evfoldname = AliConfig::GetDefaultEventFolderName();
6d50f529 224 fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
225
226 if (!fRunLoader) {
19dd5b2f 227 fRunLoader = AliRunLoader::Open(name);
6d50f529 228 }
229
230 if (!fRunLoader) {
231 AliError(Form("Can not open session for file %s.",name));
232 return kFALSE;
233 }
88cb7938 234
235 OpenInput(nEvent);
236 OpenOutput();
6d50f529 237
3e1a3ad8 238 return kTRUE;
f7336fa3 239
6d50f529 240}
3e1a3ad8 241
242//_____________________________________________________________________________
88cb7938 243Bool_t AliTRDclusterizer::OpenOutput()
3e1a3ad8 244{
245 //
246 // Open the output file
247 //
248
3e1a3ad8 249 TObjArray *ioArray = 0;
88cb7938 250
251 AliLoader* loader = fRunLoader->GetLoader("TRDLoader");
252 loader->MakeTree("R");
6d50f529 253
88cb7938 254 fClusterTree = loader->TreeR();
365d0374 255 fClusterTree->Branch("TRDcluster","TObjArray",&ioArray,32000,0);
3e1a3ad8 256
3e1a3ad8 257 return kTRUE;
258
259}
260
25ca55ce 261//_____________________________________________________________________________
262Bool_t AliTRDclusterizer::OpenOutput(TTree *clusterTree)
263{
264 //
265 // Connect the output tree
266 //
267
268 TObjArray *ioArray = 0;
269
270 fClusterTree = clusterTree;
271 fClusterTree->Branch("TRDcluster","TObjArray",&ioArray,32000,0);
272
273 return kTRUE;
274
275}
276
3e1a3ad8 277//_____________________________________________________________________________
88cb7938 278Bool_t AliTRDclusterizer::OpenInput(Int_t nEvent)
f7336fa3 279{
280 //
281 // Opens a ROOT-file with TRD-hits and reads in the digits-tree
282 //
283
f7336fa3 284 // Import the Trees for the event nEvent in the file
bdbb05bb 285 fRunLoader->GetEvent(nEvent);
88cb7938 286
f7336fa3 287 return kTRUE;
288
289}
290
291//_____________________________________________________________________________
793ff80c 292Bool_t AliTRDclusterizer::WriteClusters(Int_t det)
f7336fa3 293{
294 //
3e1a3ad8 295 // Fills TRDcluster branch in the tree with the clusters
793ff80c 296 // found in detector = det. For det=-1 writes the tree.
a3c76cdc 297 //
793ff80c 298
6d50f529 299 if ((det < -1) ||
300 (det >= AliTRDgeometry::Ndet())) {
301 AliError(Form("Unexpected detector index %d.\n",det));
3e1a3ad8 302 return kFALSE;
793ff80c 303 }
3e1a3ad8 304
3e1a3ad8 305 TBranch *branch = fClusterTree->GetBranch("TRDcluster");
306 if (!branch) {
793ff80c 307 TObjArray *ioArray = 0;
3e1a3ad8 308 branch = fClusterTree->Branch("TRDcluster","TObjArray",&ioArray,32000,0);
793ff80c 309 }
310
6d50f529 311 if ((det >= 0) &&
312 (det < AliTRDgeometry::Ndet())) {
793ff80c 313
bdbb05bb 314 Int_t nRecPoints = RecPoints()->GetEntriesFast();
3e1a3ad8 315 TObjArray *detRecPoints = new TObjArray(400);
316
317 for (Int_t i = 0; i < nRecPoints; i++) {
bdbb05bb 318 AliTRDcluster *c = (AliTRDcluster *) RecPoints()->UncheckedAt(i);
3e1a3ad8 319 if (det == c->GetDetector()) {
320 detRecPoints->AddLast(c);
321 }
322 else {
3fe61b77 323 AliError(Form("Attempt to write a cluster with unexpected detector index: got=%d expected=%d\n"
324 ,c->GetDetector()
325 ,det));
3e1a3ad8 326 }
793ff80c 327 }
328
3e1a3ad8 329 branch->SetAddress(&detRecPoints);
330 fClusterTree->Fill();
331
d9b8978b 332 delete detRecPoints;
ecf39416 333
793ff80c 334 return kTRUE;
3e1a3ad8 335
793ff80c 336 }
337
338 if (det == -1) {
339
6d50f529 340 AliInfo(Form("Writing the cluster tree %s for event %d."
341 ,fClusterTree->GetName(),fRunLoader->GetEventNumber()));
342
a6dd11e9 343 if (fRecPoints) {
344
345 branch->SetAddress(&fRecPoints);
346
347 AliLoader *loader = fRunLoader->GetLoader("TRDLoader");
348 loader->WriteRecPoints("OVERWRITE");
c630aafd 349
a6dd11e9 350 }
351 else {
352
353 AliError("Cluster tree does not exist. Cannot write clusters.\n");
354 return kFALSE;
355
356 }
357
793ff80c 358 return kTRUE;
3e1a3ad8 359
793ff80c 360 }
6d50f529 361
362 AliError(Form("Unexpected detector index %d.\n",det));
3e1a3ad8 363
793ff80c 364 return kFALSE;
88cb7938 365
f7336fa3 366}
793ff80c 367
3551db50 368//_____________________________________________________________________________
3fe61b77 369void AliTRDclusterizer::ResetHelperIndexes(AliTRDSignalIndex *indexesIn)
370{
371 //
372 // Reset the helper indexes
373 //
374
375 if (fIndexesOut)
376 {
377 // carefull here - we assume that only row number may change - most probable
378 if (indexesIn->GetNrow() <= fIndexesOut->GetNrow())
379 fIndexesOut->ResetContent();
380 else
381 fIndexesOut->ResetContentConditional(indexesIn->GetNrow()
382 , indexesIn->GetNcol()
383 , indexesIn->GetNtime());
384 }
385 else
386 {
387 fIndexesOut = new AliTRDSignalIndex(indexesIn->GetNrow()
388 , indexesIn->GetNcol()
389 , indexesIn->GetNtime());
390 }
391
392 if (fIndexesMaxima)
393 {
394 // carefull here - we assume that only row number may change - most probable
395 if (indexesIn->GetNrow() <= fIndexesMaxima->GetNrow())
396 {
397 fIndexesMaxima->ResetContent();
398 }
399 else
400 {
401 fIndexesMaxima->ResetContentConditional(indexesIn->GetNrow()
402 , indexesIn->GetNcol()
403 , indexesIn->GetNtime());
404 }
405 }
406 else
407 {
408 fIndexesMaxima = new AliTRDSignalIndex(indexesIn->GetNrow()
409 , indexesIn->GetNcol()
410 , indexesIn->GetNtime());
411 }
412
413}
414
415//_____________________________________________________________________________
416Bool_t AliTRDclusterizer::ReadDigits()
417{
418 //
419 // Reads the digits arrays from the input aliroot file
420 //
421
422 if (!fRunLoader) {
423 AliError("No run loader available");
424 return kFALSE;
425 }
426
427 AliLoader* loader = fRunLoader->GetLoader("TRDLoader");
428 if (!loader->TreeD()) {
429 loader->LoadDigits();
430 }
431
432 // Read in the digit arrays
433 return (fDigitsManager->ReadDigits(loader->TreeD()));
434
435}
436
437//_____________________________________________________________________________
438Bool_t AliTRDclusterizer::ReadDigits(TTree *digitsTree)
439{
440 //
441 // Reads the digits arrays from the input tree
442 //
443
444 // Read in the digit arrays
445 return (fDigitsManager->ReadDigits(digitsTree));
446
447}
448
449//_____________________________________________________________________________
450Bool_t AliTRDclusterizer::ReadDigits(AliRawReader *rawReader)
451{
452 //
453 // Reads the digits arrays from the ddl file
454 //
455
456 AliTRDrawData raw;
457 fDigitsManager = raw.Raw2Digits(rawReader);
458
459 return kTRUE;
460
461}
462
463//_____________________________________________________________________________
464Bool_t AliTRDclusterizer::MakeClusters()
465{
466 //
467 // Creates clusters from digits
468 //
469
470 // Propagate info from the digits manager
471 if (fAddLabels == kTRUE)
472 {
473 fAddLabels = fDigitsManager->UsesDictionaries();
474 }
475
476 Bool_t fReturn = kTRUE;
477 for (Int_t i = 0; i < AliTRDgeometry::kNdet; i++)
478 {
479
fb44eb8e 480 AliTRDdataArrayDigits *digitsIn = (AliTRDdataArrayDigits*) fDigitsManager->GetDigits(i);
3fe61b77 481 // This is to take care of switched off super modules
625f5260 482 if (!digitsIn->HasData())
3fe61b77 483 {
484 continue;
485 }
486 digitsIn->Expand();
487 AliTRDSignalIndex* indexes = fDigitsManager->GetIndexes(i);
488 if (indexes->IsAllocated() == kFALSE)
489 {
490 fDigitsManager->BuildIndexes(i);
491 }
492
493 Bool_t fR = kFALSE;
494 if (indexes->HasEntry())
495 {
496 if (fAddLabels)
497 {
498 for (Int_t iDict = 0; iDict < AliTRDdigitsManager::kNDict; iDict++)
499 {
500 AliTRDdataArrayI *tracksIn = 0;
625f5260 501 tracksIn = (AliTRDdataArrayI *) fDigitsManager->GetDictionary(i,iDict);
3fe61b77 502 tracksIn->Expand();
503 }
504 }
505 fR = MakeClusters(i);
506 fReturn = fR && fReturn;
507 }
508
509 if (fR == kFALSE)
510 {
511 WriteClusters(i);
512 ResetRecPoints();
513 }
514
515 // No compress just remove
516 fDigitsManager->RemoveDigits(i);
517 fDigitsManager->RemoveDictionaries(i);
518 fDigitsManager->ClearIndexes(i);
519
520 }
521
522 return fReturn;
523
524}
525
526//_____________________________________________________________________________
527Bool_t AliTRDclusterizer::Raw2Clusters(AliRawReader *rawReader)
528{
529 //
530 // Creates clusters from raw data
531 //
532
fc546d21 533 return Raw2ClustersChamber(rawReader);
3fe61b77 534
3fe61b77 535}
536
537//_____________________________________________________________________________
538Bool_t AliTRDclusterizer::Raw2ClustersChamber(AliRawReader *rawReader)
539{
540 //
541 // Creates clusters from raw data
542 //
543
544 // Create the digits manager
545 if (!fDigitsManager)
546 {
547 fDigitsManager = new AliTRDdigitsManager();
548 fDigitsManager->CreateArrays();
549 }
550
551 fDigitsManager->SetUseDictionaries(fAddLabels);
552
dfbb4bb9 553 AliTRDrawStreamBase *pinput = AliTRDrawStreamBase::GetRawStream(rawReader);
554 AliTRDrawStreamBase &input = *pinput;
3fe61b77 555
556 AliInfo(Form("Stream version: %s", input.IsA()->GetName()));
557
558 Int_t det = 0;
559 while ((det = input.NextChamber(fDigitsManager)) >= 0)
560 {
561 Bool_t iclusterBranch = kFALSE;
562 if (fDigitsManager->GetIndexes(det)->HasEntry())
563 {
564 iclusterBranch = MakeClusters(det);
565 }
566 if (iclusterBranch == kFALSE)
567 {
568 WriteClusters(det);
569 ResetRecPoints();
570 }
571 fDigitsManager->RemoveDigits(det);
572 fDigitsManager->RemoveDictionaries(det);
573 fDigitsManager->ClearIndexes(det);
574 }
575
576 delete fDigitsManager;
577 fDigitsManager = NULL;
dfbb4bb9 578
579 delete pinput;
580 pinput = NULL;
3fe61b77 581 return kTRUE;
582
583}
584
fa7427d0 585//_____________________________________________________________________________
586UChar_t AliTRDclusterizer::GetStatus(Short_t &signal)
587{
588 UChar_t status = 0;
589 // check if pad is masked
590 if(signal>0 && TESTBIT(signal, 10)){
591 CLRBIT(signal, 10);
592 for(int ibit=0; ibit<4; ibit++){
593 if(TESTBIT(signal, 11+ibit)){
594 SETBIT(status, ibit);
595 CLRBIT(signal, 11+ibit);
596 }
597 }
598 }
599 return status;
600}
601
3fe61b77 602//_____________________________________________________________________________
603Bool_t AliTRDclusterizer::MakeClusters(Int_t det)
604{
605 //
606 // Generates the cluster.
607 //
608
609 // Get the digits
610 // digits should be expanded beforehand!
611 // digitsIn->Expand();
fb44eb8e 612 AliTRDdataArrayDigits *digitsIn = (AliTRDdataArrayDigits *) fDigitsManager->GetDigits(det);
f5375dcb 613
3fe61b77 614 // This is to take care of switched off super modules
625f5260 615 if (!digitsIn->HasData())
3fe61b77 616 {
617 return kFALSE;
618 }
619
620 AliTRDSignalIndex *indexesIn = fDigitsManager->GetIndexes(det);
621 if (indexesIn->IsAllocated() == kFALSE)
622 {
623 AliError("Indexes do not exist!");
624 return kFALSE;
625 }
626
627 AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
628 if (!calibration)
629 {
630 AliFatal("No AliTRDcalibDB instance available\n");
631 return kFALSE;
632 }
633
3fe61b77 634 // ADC thresholds
635 // There is no ADC threshold anymore, and simParam should not be used in clusterizer. KO
acc49af9 636 Float_t adcThreshold = 0;
3fe61b77 637
c5f589b9 638 if (!AliTRDReconstructor::RecoParam())
639 {
640 AliError("RecoParam does not exist\n");
641 return kFALSE;
642 }
643
3fe61b77 644 // Threshold value for the maximum
fc546d21 645 Float_t maxThresh = AliTRDReconstructor::RecoParam()->GetClusMaxThresh();
3fe61b77 646 // Threshold value for the digit signal
fc546d21 647 Float_t sigThresh = AliTRDReconstructor::RecoParam()->GetClusSigThresh();
3fe61b77 648
649 // Iteration limit for unfolding procedure
650 const Float_t kEpsilon = 0.01;
651 const Int_t kNclus = 3;
652 const Int_t kNsig = 5;
653
654 Int_t iUnfold = 0;
655 Double_t ratioLeft = 1.0;
656 Double_t ratioRight = 1.0;
657
658 Double_t padSignal[kNsig];
659 Double_t clusterSignal[kNclus];
660
661 Int_t icham = indexesIn->GetChamber();
662 Int_t iplan = indexesIn->GetPlane();
663 Int_t isect = indexesIn->GetSM();
664
665 // Start clustering in the chamber
666
667 Int_t idet = AliTRDgeometry::GetDetector(iplan,icham,isect);
668 if (idet != det)
669 {
670 AliError("Strange Detector number Missmatch!");
671 return kFALSE;
672 }
673
674 // TRD space point transformation
675 fTransform->SetDetector(det);
676
677 Int_t ilayer = AliGeomManager::kTRD1 + iplan;
678 Int_t imodule = icham + AliTRDgeometry::Ncham() * isect;
679 UShort_t volid = AliGeomManager::LayerToVolUID(ilayer,imodule);
680
681 Int_t nColMax = digitsIn->GetNcol();
fa7427d0 682 Int_t nRowMax = digitsIn->GetNrow();
3fe61b77 683 Int_t nTimeTotal = digitsIn->GetNtime();
684
685 // Detector wise calibration object for the gain factors
0e09df31 686 const AliTRDCalDet *calGainFactorDet = calibration->GetGainFactorDet();
3fe61b77 687 // Calibration object with pad wise values for the gain factors
0e09df31 688 AliTRDCalROC *calGainFactorROC = calibration->GetGainFactorROC(idet);
3fe61b77 689 // Calibration value for chamber wise gain factor
0e09df31 690 Float_t calGainFactorDetValue = calGainFactorDet->GetValue(idet);
691
3fe61b77 692 Int_t nClusters = 0;
693
fa7427d0 694 AliTRDdataArrayF *digitsOut = new AliTRDdataArrayF(nRowMax, nColMax, nTimeTotal);
695 AliTRDdataArrayS padStatus(nRowMax, nColMax, nTimeTotal);
3fe61b77 696
697 ResetHelperIndexes(indexesIn);
698
699 // Apply the gain and the tail cancelation via digital filter
700 TailCancelation(digitsIn
f5375dcb 701 ,digitsOut
702 ,indexesIn
703 ,fIndexesOut
704 ,nTimeTotal
705 ,adcThreshold
706 ,calGainFactorROC
707 ,calGainFactorDetValue);
3fe61b77 708
709 Int_t row = 0;
710 Int_t col = 0;
711 Int_t time = 0;
712 Int_t iPad = 0;
713
f5375dcb 714 UChar_t status[3]={0, 0, 0}, ipos = 0;
3fe61b77 715 fIndexesOut->ResetCounters();
fa7427d0 716 while (fIndexesOut->NextRCTbinIndex(row, col, time)){
f5375dcb 717 Float_t signalM = TMath::Abs(digitsOut->GetDataUnchecked(row,col,time));
718 status[1] = digitsIn->GetPadStatus(row,col,time);
719 if(status[1]) SETBIT(ipos, AliTRDcluster::kMaskedCenter);
720
721 // Look for the maximum
722 if (signalM >= maxThresh) {
723 if (col + 1 >= nColMax || col-1 < 0) continue;
724
725 Float_t signalL = TMath::Abs(digitsOut->GetDataUnchecked(row,col+1,time));
726 status[0] = digitsIn->GetPadStatus(row,col+1,time);
727 if(status[0]) SETBIT(ipos, AliTRDcluster::kMaskedLeft);
728
729 Float_t signalR = TMath::Abs(digitsOut->GetDataUnchecked(row,col-1,time));
730 status[2] = digitsIn->GetPadStatus(row,col-1,time);
731 if(status[2]) SETBIT(ipos, AliTRDcluster::kMaskedRight);
732
733 // reject candidates with more than 1 problematic pad
734 if(ipos == 3 || ipos > 4) continue;
735
736 if(!status[1]){ // good central pad
737 if(!ipos){ // all pads are OK
738 if ((signalL <= signalM) && (signalR < signalM)) {
739 if ((signalL >= sigThresh) || (signalR >= sigThresh)) {
740 // Maximum found, mark the position by a negative signal
741 digitsOut->SetDataUnchecked(row,col,time,-signalM);
742 fIndexesMaxima->AddIndexTBin(row,col,time);
743 padStatus.SetDataUnchecked(row, col, time, ipos);
744 }
745 }
746 } else { // one of the neighbouring pads are bad
747 if(status[0] && signalR < signalM && signalR >= sigThresh){
748 digitsOut->SetDataUnchecked(row,col,time,-signalM);
749 digitsOut->SetDataUnchecked(row, col, time+1, 0.);
750 fIndexesMaxima->AddIndexTBin(row,col,time);
751 padStatus.SetDataUnchecked(row, col, time, ipos);
752 } else if(status[2] && signalL <= signalM && signalL >= sigThresh){
753 digitsOut->SetDataUnchecked(row,col,time,-signalM);
754 digitsOut->SetDataUnchecked(row, col, time-1, 0.);
755 fIndexesMaxima->AddIndexTBin(row,col,time);
756 padStatus.SetDataUnchecked(row, col, time, ipos);
757 }
758 }
759 } else { // wrong maximum pad
760 if ((signalL >= sigThresh) || (signalR >= sigThresh)) {
761 // Maximum found, mark the position by a negative signal
762 digitsOut->SetDataUnchecked(row,col,time,-maxThresh);
763 fIndexesMaxima->AddIndexTBin(row,col,time);
764 padStatus.SetDataUnchecked(row, col, time, ipos);
765 }
766 }
767 }
768 }
769
770 // The index to the first cluster of a given ROC
771 Int_t firstClusterROC = -1;
772 // The number of cluster in a given ROC
773 Int_t nClusterROC = 0;
774
775 // Now check the maxima and calculate the cluster position
776 fIndexesMaxima->ResetCounters();
777 while (fIndexesMaxima->NextRCTbinIndex(row, col, time)) {
3fe61b77 778
779 // Maximum found ?
f5375dcb 780 if (digitsOut->GetDataUnchecked(row,col,time) < 0.0) {
781 for (iPad = 0; iPad < kNclus; iPad++) {
782 Int_t iPadCol = col - 1 + iPad;
783 clusterSignal[iPad] = TMath::Abs(digitsOut->GetDataUnchecked(row,iPadCol,time));
784 }
3fe61b77 785
f5375dcb 786 // Count the number of pads in the cluster
787 Int_t nPadCount = 0;
788 Int_t ii;
789 // Look to the left
790 ii = 0;
791 while (TMath::Abs(digitsOut->GetDataUnchecked(row,col-ii ,time)) >= sigThresh) {
792 nPadCount++;
793 ii++;
794 if (col-ii < 0) break;
795 }
796 // Look to the right
797 ii = 0;
798 while (TMath::Abs(digitsOut->GetDataUnchecked(row,col+ii+1,time)) >= sigThresh){
799 nPadCount++;
800 ii++;
801 if (col+ii+1 >= nColMax) break;
802 }
803 nClusters++;
804
805 // Look for 5 pad cluster with minimum in the middle
806 Bool_t fivePadCluster = kFALSE;
807 if (col < (nColMax - 3)){
808 if (digitsOut->GetDataUnchecked(row,col+2,time) < 0) {
809 fivePadCluster = kTRUE;
810 }
811 if ((fivePadCluster) && (col < (nColMax - 5))) {
812 if (digitsOut->GetDataUnchecked(row,col+4,time) >= sigThresh){
813 fivePadCluster = kFALSE;
814 }
815 }
816 if ((fivePadCluster) && (col > 1)){
817 if (digitsOut->GetDataUnchecked(row,col-2,time) >= sigThresh){
818 fivePadCluster = kFALSE;
819 }
820 }
821 }
3fe61b77 822
f5375dcb 823 // 5 pad cluster
824 // Modify the signal of the overlapping pad for the left part
825 // of the cluster which remains from a previous unfolding
826 if (iUnfold) {
827 clusterSignal[0] *= ratioLeft;
828 iUnfold = 0;
829 }
3fe61b77 830
f5375dcb 831 // Unfold the 5 pad cluster
832 if (fivePadCluster){
833 for (iPad = 0; iPad < kNsig; iPad++) {
834 padSignal[iPad] = TMath::Abs(digitsOut->GetDataUnchecked(row
835 ,col-1+iPad
836 ,time));
837 }
838 // Unfold the two maxima and set the signal on
839 // the overlapping pad to the ratio
840 ratioRight = Unfold(kEpsilon,iplan,padSignal);
841 ratioLeft = 1.0 - ratioRight;
842 clusterSignal[2] *= ratioRight;
843 iUnfold = 1;
844 }
3fe61b77 845
f5375dcb 846 // The position of the cluster in COL direction relative to the center pad (pad units)
847 Double_t clusterPosCol = 0.0;
848 if (AliTRDReconstructor::RecoParam()->LUTOn()) {
849 // Calculate the position of the cluster by using the
850 // lookup table method
851 clusterPosCol = LUTposition(iplan,clusterSignal[0]
852 ,clusterSignal[1]
853 ,clusterSignal[2]);
854 } else {
855 // Calculate the position of the cluster by using the
856 // center of gravity method
857 for (Int_t i = 0; i < kNsig; i++) {
858 padSignal[i] = 0.0;
859 }
860 padSignal[2] = TMath::Abs(digitsOut->GetDataUnchecked(row,col ,time)); // Central pad
861 padSignal[1] = TMath::Abs(digitsOut->GetDataUnchecked(row,col-1,time)); // Left pad
862 padSignal[3] = TMath::Abs(digitsOut->GetDataUnchecked(row,col+1,time)); // Right pad
863
864 if ((col > 2) &&
865 (TMath::Abs(digitsOut->GetDataUnchecked(row,col-2,time)) < padSignal[1])) {
866 padSignal[0] = TMath::Abs(digitsOut->GetDataUnchecked(row,col-2,time));
867 }
868 if ((col < nColMax - 3) &&
869 (TMath::Abs(digitsOut->GetDataUnchecked(row,col+2,time)) < padSignal[3])){
870 padSignal[4] = TMath::Abs(digitsOut->GetDataUnchecked(row,col+2,time));
871 }
872 clusterPosCol = GetCOG(padSignal);
873 }
3fe61b77 874
f5375dcb 875 // Store the amplitudes of the pads in the cluster for later analysis
876 // and check whether one of these pads is masked in the database
877 Short_t signals[7] = { 0, 0, 0, 0, 0, 0, 0 };
878 for (Int_t jPad = col-3; jPad <= col+3; jPad++) {
879 if ((jPad < 0) ||
880 (jPad >= nColMax-1)) {
881 continue;
882 }
883 signals[jPad-col+3] = TMath::Nint(TMath::Abs(digitsOut->GetDataUnchecked(row,jPad,time)));
884 }
3fe61b77 885
f5375dcb 886 // Transform the local cluster coordinates into calibrated
887 // space point positions defined in the local tracking system.
888 // Here the calibration for T0, Vdrift and ExB is applied as well.
889 Double_t clusterXYZ[6];
890 clusterXYZ[0] = clusterPosCol;
891 clusterXYZ[1] = clusterSignal[0];
892 clusterXYZ[2] = clusterSignal[1];
893 clusterXYZ[3] = clusterSignal[2];
894 clusterXYZ[4] = 0.0;
895 clusterXYZ[5] = 0.0;
896 Int_t clusterRCT[3];
897 clusterRCT[0] = row;
898 clusterRCT[1] = col;
899 clusterRCT[2] = 0;
900
901 Bool_t out = kTRUE;
902 if (fTransform->Transform(clusterXYZ, clusterRCT, ((UInt_t) time), out, 0)) {
903
904 // Add the cluster to the output array
905 // The track indices will be stored later
906 Float_t clusterPos[3];
907 clusterPos[0] = clusterXYZ[0];
908 clusterPos[1] = clusterXYZ[1];
909 clusterPos[2] = clusterXYZ[2];
910 Float_t clusterSig[2];
911 clusterSig[0] = clusterXYZ[4];
912 clusterSig[1] = clusterXYZ[5];
913 Double_t clusterCharge = clusterXYZ[3];
914 Char_t clusterTimeBin = ((Char_t) clusterRCT[2]);
915 AliTRDcluster *cluster = new AliTRDcluster(idet
916 ,clusterCharge
917 ,clusterPos
918 ,clusterSig
919 ,0x0
920 ,((Char_t) nPadCount)
921 ,signals
922 ,((UChar_t) col)
923 ,((UChar_t) row)
924 ,((UChar_t) time)
925 ,clusterTimeBin
926 ,clusterPosCol
927 ,volid);
928 cluster->SetInChamber(!out);
929 UChar_t maskPosition = padStatus.GetDataUnchecked(row, col, time);
930 if(maskPosition){
931 cluster->SetPadMaskedPosition(maskPosition);
932
933 if(maskPosition & AliTRDcluster::kMaskedLeft) cluster->SetPadMaskedStatus(status[0]);
934 else if(maskPosition & AliTRDcluster::kMaskedCenter) cluster->SetPadMaskedStatus(status[1]);
935 else cluster->SetPadMaskedStatus(status[2]);
936 }
3fe61b77 937
f5375dcb 938 // Temporarily store the row, column and time bin of the center pad
939 // Used to later on assign the track indices
940 cluster->SetLabel( row,0);
941 cluster->SetLabel( col,1);
942 cluster->SetLabel(time,2);
943
944 RecPoints()->Add(cluster);
3fe61b77 945
f5375dcb 946 // Store the index of the first cluster in the current ROC
947 if (firstClusterROC < 0) {
948 firstClusterROC = RecPoints()->GetEntriesFast() - 1;
949 }
950
951 // Count the number of cluster in the current ROC
952 nClusterROC++;
953
954 } // if: Transform ok ?
955 } // if: Maximum found ?
956 }
3fe61b77 957
958 delete digitsOut;
959
f5375dcb 960 if (fAddLabels) AddLabels(idet, firstClusterROC, nClusterROC);
3fe61b77 961
962 // Write the cluster and reset the array
963 WriteClusters(idet);
964 ResetRecPoints();
965
966 return kTRUE;
967
968}
969
970//_____________________________________________________________________________
971Bool_t AliTRDclusterizer::AddLabels(Int_t idet, Int_t firstClusterROC, Int_t nClusterROC)
3551db50 972{
973 //
3fe61b77 974 // Add the track indices to the found clusters
3551db50 975 //
976
3fe61b77 977 const Int_t kNclus = 3;
978 const Int_t kNdict = AliTRDdigitsManager::kNDict;
979 const Int_t kNtrack = kNdict * kNclus;
980
981 Int_t iClusterROC = 0;
982
983 Int_t row = 0;
984 Int_t col = 0;
985 Int_t time = 0;
986 Int_t iPad = 0;
987
988 // Temporary array to collect the track indices
989 Int_t *idxTracks = new Int_t[kNtrack*nClusterROC];
990
991 // Loop through the dictionary arrays one-by-one
992 // to keep memory consumption low
993 AliTRDdataArrayI *tracksIn = 0;
994 for (Int_t iDict = 0; iDict < kNdict; iDict++) {
995
996 // tracksIn should be expanded beforehand!
625f5260 997 tracksIn = (AliTRDdataArrayI *) fDigitsManager->GetDictionary(idet,iDict);
3fe61b77 998
999 // Loop though the clusters found in this ROC
1000 for (iClusterROC = 0; iClusterROC < nClusterROC; iClusterROC++) {
1001
1002 AliTRDcluster *cluster = (AliTRDcluster *)
1003 RecPoints()->UncheckedAt(firstClusterROC+iClusterROC);
1004 row = cluster->GetLabel(0);
1005 col = cluster->GetLabel(1);
1006 time = cluster->GetLabel(2);
1007
1008 for (iPad = 0; iPad < kNclus; iPad++) {
1009 Int_t iPadCol = col - 1 + iPad;
1010 Int_t index = tracksIn->GetDataUnchecked(row,iPadCol,time) - 1;
1011 idxTracks[3*iPad+iDict + iClusterROC*kNtrack] = index;
1012 }
1013
1014 }
1015
1016 }
1017
1018 // Copy the track indices into the cluster
1019 // Loop though the clusters found in this ROC
1020 for (iClusterROC = 0; iClusterROC < nClusterROC; iClusterROC++) {
1021
1022 AliTRDcluster *cluster = (AliTRDcluster *)
1023 RecPoints()->UncheckedAt(firstClusterROC+iClusterROC);
1024 cluster->SetLabel(-9999,0);
1025 cluster->SetLabel(-9999,1);
1026 cluster->SetLabel(-9999,2);
1027
1028 cluster->AddTrackIndex(&idxTracks[iClusterROC*kNtrack]);
1029
1030 }
1031
1032 delete [] idxTracks;
1033
1034 return kTRUE;
1035
1036}
1037
1038//_____________________________________________________________________________
acc49af9 1039Double_t AliTRDclusterizer::GetCOG(Double_t signal[5]) const
3fe61b77 1040{
1041 //
1042 // Get COG position
1043 // Used for clusters with more than 3 pads - where LUT not applicable
1044 //
1045
1046 Double_t sum = signal[0]
1047 + signal[1]
1048 + signal[2]
1049 + signal[3]
1050 + signal[4];
1051
1052 Double_t res = (0.0 * (-signal[0] + signal[4])
1053 + (-signal[1] + signal[3])) / sum;
1054
1055 return res;
1056
1057}
1058
1059//_____________________________________________________________________________
1060Double_t AliTRDclusterizer::Unfold(Double_t eps, Int_t plane, Double_t *padSignal)
1061{
1062 //
1063 // Method to unfold neighbouring maxima.
1064 // The charge ratio on the overlapping pad is calculated
1065 // until there is no more change within the range given by eps.
1066 // The resulting ratio is then returned to the calling method.
1067 //
1068
1069 AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
6d50f529 1070 if (!calibration) {
3fe61b77 1071 AliError("No AliTRDcalibDB instance available\n");
1072 return kFALSE;
6d50f529 1073 }
3fe61b77 1074
1075 Int_t irc = 0;
1076 Int_t itStep = 0; // Count iteration steps
1077
1078 Double_t ratio = 0.5; // Start value for ratio
1079 Double_t prevRatio = 0.0; // Store previous ratio
1080
1081 Double_t newLeftSignal[3] = { 0.0, 0.0, 0.0 }; // Array to store left cluster signal
1082 Double_t newRightSignal[3] = { 0.0, 0.0, 0.0 }; // Array to store right cluster signal
1083 Double_t newSignal[3] = { 0.0, 0.0, 0.0 };
1084
1085 // Start the iteration
1086 while ((TMath::Abs(prevRatio - ratio) > eps) && (itStep < 10)) {
1087
1088 itStep++;
1089 prevRatio = ratio;
1090
1091 // Cluster position according to charge ratio
1092 Double_t maxLeft = (ratio*padSignal[2] - padSignal[0])
1093 / (padSignal[0] + padSignal[1] + ratio*padSignal[2]);
1094 Double_t maxRight = (padSignal[4] - (1-ratio)*padSignal[2])
1095 / ((1.0 - ratio)*padSignal[2] + padSignal[3] + padSignal[4]);
1096
1097 // Set cluster charge ratio
1098 irc = calibration->PadResponse(1.0,maxLeft ,plane,newSignal);
1099 Double_t ampLeft = padSignal[1] / newSignal[1];
1100 irc = calibration->PadResponse(1.0,maxRight,plane,newSignal);
1101 Double_t ampRight = padSignal[3] / newSignal[1];
1102
1103 // Apply pad response to parameters
1104 irc = calibration->PadResponse(ampLeft ,maxLeft ,plane,newLeftSignal );
1105 irc = calibration->PadResponse(ampRight,maxRight,plane,newRightSignal);
1106
1107 // Calculate new overlapping ratio
1108 ratio = TMath::Min((Double_t) 1.0
1109 ,newLeftSignal[2] / (newLeftSignal[2] + newRightSignal[0]));
1110
b43a3e17 1111 }
88719a08 1112
3fe61b77 1113 return ratio;
1114
1115}
88719a08 1116
3fe61b77 1117//_____________________________________________________________________________
fb44eb8e 1118void AliTRDclusterizer::TailCancelation(AliTRDdataArrayDigits *digitsIn
3fe61b77 1119 , AliTRDdataArrayF *digitsOut
1120 , AliTRDSignalIndex *indexesIn
1121 , AliTRDSignalIndex *indexesOut
1122 , Int_t nTimeTotal
acc49af9 1123 , Float_t adcThreshold
3fe61b77 1124 , AliTRDCalROC *calGainFactorROC
1125 , Float_t calGainFactorDetValue)
1126{
1127 //
1128 // Applies the tail cancelation and gain factors:
1129 // Transform digitsIn to digitsOut
1130 //
1131
1132 Int_t iRow = 0;
1133 Int_t iCol = 0;
1134 Int_t iTime = 0;
1135
3fe61b77 1136 Double_t *inADC = new Double_t[nTimeTotal]; // ADC data before tail cancellation
1137 Double_t *outADC = new Double_t[nTimeTotal]; // ADC data after tail cancellation
1138 indexesIn->ResetCounters();
1139 while (indexesIn->NextRCIndex(iRow, iCol))
1140 {
1141 Float_t calGainFactorROCValue = calGainFactorROC->GetValue(iCol,iRow);
1142 Double_t gain = calGainFactorDetValue
1143 * calGainFactorROCValue;
1144
f5375dcb 1145 Bool_t corrupted = kFALSE;
3fe61b77 1146 for (iTime = 0; iTime < nTimeTotal; iTime++)
1147 {
1148 // Apply gain gain factor
1149 inADC[iTime] = digitsIn->GetDataUnchecked(iRow,iCol,iTime);
f5375dcb 1150 if(digitsIn->GetPadStatus(iRow, iCol, iTime)) corrupted = kTRUE;
3fe61b77 1151 inADC[iTime] /= gain;
1152 outADC[iTime] = inADC[iTime];
1153 }
f5375dcb 1154 if(!corrupted)
1155 {
1156 // Apply the tail cancelation via the digital filter
1157 // (only for non-coorupted pads)
1158 if (AliTRDReconstructor::RecoParam()->TCOn())
1159 {
1160 DeConvExp(inADC,outADC,nTimeTotal,AliTRDReconstructor::RecoParam()->GetTCnexp());
1161 }
1162 }
3fe61b77 1163
1164 indexesIn->ResetTbinCounter();
1165 while (indexesIn->NextTbinIndex(iTime))
1166 {
1167 // Store the amplitude of the digit if above threshold
acc49af9 1168 if (outADC[iTime] > adcThreshold)
3fe61b77 1169 {
1170 digitsOut->SetDataUnchecked(iRow,iCol,iTime,outADC[iTime]);
1171 indexesOut->AddIndexTBin(iRow,iCol,iTime);
1172 }
1173 } // while itime
1174
1175 } // while irow icol
1176
1177 delete [] inADC;
1178 delete [] outADC;
1179
1180 return;
1181
1182}
1183
1184//_____________________________________________________________________________
1185void AliTRDclusterizer::DeConvExp(Double_t *source, Double_t *target
1186 , Int_t n, Int_t nexp)
1187{
1188 //
1189 // Tail cancellation by deconvolution for PASA v4 TRF
1190 //
1191
1192 Double_t rates[2];
1193 Double_t coefficients[2];
1194
1195 // Initialization (coefficient = alpha, rates = lambda)
1196 Double_t r1 = 1.0;
1197 Double_t r2 = 1.0;
1198 Double_t c1 = 0.5;
1199 Double_t c2 = 0.5;
1200
1201 if (nexp == 1) { // 1 Exponentials
1202 r1 = 1.156;
1203 r2 = 0.130;
1204 c1 = 0.066;
1205 c2 = 0.000;
1206 }
1207 if (nexp == 2) { // 2 Exponentials
1208 r1 = 1.156;
1209 r2 = 0.130;
1210 c1 = 0.114;
1211 c2 = 0.624;
1212 }
1213
1214 coefficients[0] = c1;
1215 coefficients[1] = c2;
1216
1217 Double_t dt = 0.1;
1218
1219 rates[0] = TMath::Exp(-dt/(r1));
1220 rates[1] = TMath::Exp(-dt/(r2));
1221
1222 Int_t i = 0;
1223 Int_t k = 0;
1224
1225 Double_t reminder[2];
1226 Double_t correction = 0.0;
1227 Double_t result = 0.0;
1228
1229 // Attention: computation order is important
1230 for (k = 0; k < nexp; k++) {
1231 reminder[k] = 0.0;
1232 }
1233
1234 for (i = 0; i < n; i++) {
1235
1236 result = (source[i] - correction); // No rescaling
1237 target[i] = result;
1238
1239 for (k = 0; k < nexp; k++) {
1240 reminder[k] = rates[k] * (reminder[k] + coefficients[k] * result);
1241 }
1242
1243 correction = 0.0;
1244 for (k = 0; k < nexp; k++) {
1245 correction += reminder[k];
1246 }
1247
1248 }
6d50f529 1249
1250}
1251
1252//_____________________________________________________________________________
1253void AliTRDclusterizer::ResetRecPoints()
1254{
1255 //
1256 // Resets the list of rec points
1257 //
1258
1259 if (fRecPoints) {
1260 fRecPoints->Delete();
1261 }
1262
1263}
1264
1265//_____________________________________________________________________________
34eaaa7e 1266TObjArray *AliTRDclusterizer::RecPoints()
6d50f529 1267{
1268 //
1269 // Returns the list of rec points
1270 //
1271
1272 if (!fRecPoints) {
1273 fRecPoints = new TObjArray(400);
1274 }
1275
1276 return fRecPoints;
1277
3551db50 1278}
fc546d21 1279
1280//_____________________________________________________________________________
1281void AliTRDclusterizer::FillLUT()
1282{
1283 //
1284 // Create the LUT
1285 //
1286
1287 const Int_t kNlut = 128;
1288
1289 fLUTbin = AliTRDgeometry::kNplan * kNlut;
1290
1291 // The lookup table from Bogdan
1292 Float_t lut[AliTRDgeometry::kNplan][kNlut] = {
1293 {
1294 0.0070, 0.0150, 0.0224, 0.0298, 0.0374, 0.0454, 0.0533, 0.0611,
1295 0.0684, 0.0755, 0.0827, 0.0900, 0.0975, 0.1049, 0.1120, 0.1187,
1296 0.1253, 0.1318, 0.1385, 0.1453, 0.1519, 0.1584, 0.1646, 0.1704,
1297 0.1762, 0.1821, 0.1879, 0.1938, 0.1996, 0.2053, 0.2108, 0.2160,
1298 0.2210, 0.2260, 0.2310, 0.2361, 0.2411, 0.2461, 0.2509, 0.2557,
1299 0.2602, 0.2646, 0.2689, 0.2732, 0.2774, 0.2816, 0.2859, 0.2901,
1300 0.2942, 0.2983, 0.3022, 0.3061, 0.3099, 0.3136, 0.3172, 0.3207,
1301 0.3242, 0.3278, 0.3312, 0.3347, 0.3382, 0.3416, 0.3450, 0.3483,
1302 0.3515, 0.3547, 0.3579, 0.3609, 0.3639, 0.3669, 0.3698, 0.3727,
1303 0.3756, 0.3785, 0.3813, 0.3842, 0.3870, 0.3898, 0.3926, 0.3952,
1304 0.3979, 0.4005, 0.4032, 0.4057, 0.4082, 0.4108, 0.4132, 0.4157,
1305 0.4181, 0.4205, 0.4228, 0.4252, 0.4275, 0.4299, 0.4322, 0.4345,
1306 0.4367, 0.4390, 0.4412, 0.4434, 0.4456, 0.4478, 0.4499, 0.4520,
1307 0.4541, 0.4562, 0.4583, 0.4603, 0.4623, 0.4643, 0.4663, 0.4683,
1308 0.4702, 0.4722, 0.4741, 0.4758, 0.4774, 0.4790, 0.4805, 0.4824,
1309 0.4844, 0.4863, 0.4883, 0.4902, 0.4921, 0.4940, 0.4959, 0.4978
1310 },
1311 {
1312 0.0072, 0.0156, 0.0235, 0.0313, 0.0394, 0.0478, 0.0561, 0.0642,
1313 0.0718, 0.0792, 0.0868, 0.0947, 0.1025, 0.1101, 0.1172, 0.1241,
1314 0.1309, 0.1378, 0.1449, 0.1518, 0.1586, 0.1650, 0.1710, 0.1770,
1315 0.1830, 0.1891, 0.1952, 0.2011, 0.2070, 0.2125, 0.2177, 0.2229,
1316 0.2280, 0.2332, 0.2383, 0.2435, 0.2484, 0.2533, 0.2581, 0.2627,
1317 0.2670, 0.2714, 0.2757, 0.2799, 0.2842, 0.2884, 0.2927, 0.2968,
1318 0.3008, 0.3048, 0.3086, 0.3123, 0.3159, 0.3195, 0.3231, 0.3266,
1319 0.3301, 0.3335, 0.3370, 0.3404, 0.3438, 0.3471, 0.3504, 0.3536,
1320 0.3567, 0.3598, 0.3628, 0.3657, 0.3686, 0.3715, 0.3744, 0.3772,
1321 0.3800, 0.3828, 0.3856, 0.3884, 0.3911, 0.3938, 0.3965, 0.3991,
1322 0.4016, 0.4042, 0.4067, 0.4092, 0.4116, 0.4140, 0.4164, 0.4187,
1323 0.4211, 0.4234, 0.4257, 0.4280, 0.4302, 0.4325, 0.4347, 0.4369,
1324 0.4391, 0.4413, 0.4434, 0.4456, 0.4477, 0.4497, 0.4518, 0.4538,
1325 0.4558, 0.4578, 0.4598, 0.4618, 0.4637, 0.4656, 0.4675, 0.4694,
1326 0.4713, 0.4732, 0.4750, 0.4766, 0.4781, 0.4797, 0.4813, 0.4832,
1327 0.4851, 0.4870, 0.4888, 0.4906, 0.4925, 0.4942, 0.4960, 0.4978
1328 },
1329 {
1330 0.0075, 0.0163, 0.0246, 0.0328, 0.0415, 0.0504, 0.0592, 0.0674,
1331 0.0753, 0.0832, 0.0914, 0.0996, 0.1077, 0.1154, 0.1225, 0.1296,
1332 0.1369, 0.1442, 0.1515, 0.1585, 0.1652, 0.1714, 0.1776, 0.1839,
1333 0.1902, 0.1965, 0.2025, 0.2085, 0.2141, 0.2194, 0.2247, 0.2299,
1334 0.2352, 0.2405, 0.2457, 0.2507, 0.2557, 0.2604, 0.2649, 0.2693,
1335 0.2737, 0.2780, 0.2823, 0.2867, 0.2909, 0.2951, 0.2992, 0.3033,
1336 0.3072, 0.3110, 0.3146, 0.3182, 0.3218, 0.3253, 0.3288, 0.3323,
1337 0.3357, 0.3392, 0.3426, 0.3459, 0.3492, 0.3524, 0.3555, 0.3586,
1338 0.3616, 0.3645, 0.3674, 0.3703, 0.3731, 0.3759, 0.3787, 0.3815,
1339 0.3843, 0.3870, 0.3897, 0.3925, 0.3950, 0.3976, 0.4002, 0.4027,
1340 0.4052, 0.4076, 0.4101, 0.4124, 0.4148, 0.4171, 0.4194, 0.4217,
1341 0.4239, 0.4262, 0.4284, 0.4306, 0.4328, 0.4350, 0.4371, 0.4393,
1342 0.4414, 0.4435, 0.4455, 0.4476, 0.4496, 0.4516, 0.4536, 0.4555,
1343 0.4575, 0.4594, 0.4613, 0.4632, 0.4650, 0.4669, 0.4687, 0.4705,
1344 0.4723, 0.4741, 0.4758, 0.4773, 0.4789, 0.4804, 0.4821, 0.4839,
1345 0.4857, 0.4875, 0.4893, 0.4910, 0.4928, 0.4945, 0.4961, 0.4978
1346 },
1347 {
1348 0.0078, 0.0171, 0.0258, 0.0345, 0.0438, 0.0532, 0.0624, 0.0708,
1349 0.0791, 0.0875, 0.0962, 0.1048, 0.1130, 0.1206, 0.1281, 0.1356,
1350 0.1432, 0.1508, 0.1582, 0.1651, 0.1716, 0.1780, 0.1845, 0.1910,
1351 0.1975, 0.2038, 0.2099, 0.2155, 0.2210, 0.2263, 0.2317, 0.2371,
1352 0.2425, 0.2477, 0.2528, 0.2578, 0.2626, 0.2671, 0.2715, 0.2759,
1353 0.2803, 0.2846, 0.2890, 0.2933, 0.2975, 0.3016, 0.3056, 0.3095,
1354 0.3132, 0.3168, 0.3204, 0.3239, 0.3274, 0.3309, 0.3344, 0.3378,
1355 0.3412, 0.3446, 0.3479, 0.3511, 0.3543, 0.3574, 0.3603, 0.3633,
1356 0.3662, 0.3690, 0.3718, 0.3747, 0.3774, 0.3802, 0.3829, 0.3857,
1357 0.3883, 0.3910, 0.3936, 0.3962, 0.3987, 0.4012, 0.4037, 0.4061,
1358 0.4085, 0.4109, 0.4132, 0.4155, 0.4177, 0.4200, 0.4222, 0.4244,
1359 0.4266, 0.4288, 0.4309, 0.4331, 0.4352, 0.4373, 0.4394, 0.4414,
1360 0.4435, 0.4455, 0.4475, 0.4494, 0.4514, 0.4533, 0.4552, 0.4571,
1361 0.4590, 0.4608, 0.4626, 0.4645, 0.4662, 0.4680, 0.4698, 0.4715,
1362 0.4733, 0.4750, 0.4766, 0.4781, 0.4796, 0.4812, 0.4829, 0.4846,
1363 0.4863, 0.4880, 0.4897, 0.4914, 0.4930, 0.4946, 0.4963, 0.4979
1364 },
1365 {
1366 0.0081, 0.0178, 0.0270, 0.0364, 0.0463, 0.0562, 0.0656, 0.0744,
1367 0.0831, 0.0921, 0.1013, 0.1102, 0.1183, 0.1261, 0.1339, 0.1419,
1368 0.1499, 0.1576, 0.1648, 0.1715, 0.1782, 0.1849, 0.1917, 0.1984,
1369 0.2048, 0.2110, 0.2167, 0.2223, 0.2278, 0.2333, 0.2389, 0.2444,
1370 0.2497, 0.2548, 0.2598, 0.2645, 0.2691, 0.2735, 0.2780, 0.2824,
1371 0.2868, 0.2912, 0.2955, 0.2997, 0.3038, 0.3078, 0.3116, 0.3152,
1372 0.3188, 0.3224, 0.3259, 0.3294, 0.3329, 0.3364, 0.3398, 0.3432,
1373 0.3465, 0.3497, 0.3529, 0.3561, 0.3591, 0.3620, 0.3649, 0.3677,
1374 0.3705, 0.3733, 0.3761, 0.3788, 0.3816, 0.3843, 0.3869, 0.3896,
1375 0.3922, 0.3948, 0.3973, 0.3998, 0.4022, 0.4047, 0.4070, 0.4094,
1376 0.4117, 0.4139, 0.4162, 0.4184, 0.4206, 0.4227, 0.4249, 0.4270,
1377 0.4291, 0.4313, 0.4334, 0.4354, 0.4375, 0.4395, 0.4415, 0.4435,
1378 0.4455, 0.4474, 0.4493, 0.4512, 0.4531, 0.4550, 0.4568, 0.4586,
1379 0.4604, 0.4622, 0.4639, 0.4657, 0.4674, 0.4691, 0.4708, 0.4725,
1380 0.4742, 0.4758, 0.4773, 0.4788, 0.4803, 0.4819, 0.4836, 0.4852,
1381 0.4869, 0.4885, 0.4901, 0.4917, 0.4933, 0.4948, 0.4964, 0.4979
1382 },
1383 {
1384 0.0085, 0.0189, 0.0288, 0.0389, 0.0497, 0.0603, 0.0699, 0.0792,
1385 0.0887, 0.0985, 0.1082, 0.1170, 0.1253, 0.1336, 0.1421, 0.1505,
1386 0.1587, 0.1662, 0.1733, 0.1803, 0.1874, 0.1945, 0.2014, 0.2081,
1387 0.2143, 0.2201, 0.2259, 0.2316, 0.2374, 0.2431, 0.2487, 0.2541,
1388 0.2593, 0.2642, 0.2689, 0.2735, 0.2781, 0.2826, 0.2872, 0.2917,
1389 0.2961, 0.3003, 0.3045, 0.3086, 0.3125, 0.3162, 0.3198, 0.3235,
1390 0.3270, 0.3306, 0.3342, 0.3377, 0.3411, 0.3446, 0.3479, 0.3511,
1391 0.3543, 0.3575, 0.3605, 0.3634, 0.3663, 0.3691, 0.3720, 0.3748,
1392 0.3775, 0.3803, 0.3830, 0.3857, 0.3884, 0.3911, 0.3937, 0.3962,
1393 0.3987, 0.4012, 0.4036, 0.4060, 0.4084, 0.4107, 0.4129, 0.4152,
1394 0.4174, 0.4196, 0.4218, 0.4239, 0.4261, 0.4282, 0.4303, 0.4324,
1395 0.4344, 0.4365, 0.4385, 0.4405, 0.4425, 0.4445, 0.4464, 0.4483,
1396 0.4502, 0.4521, 0.4539, 0.4558, 0.4576, 0.4593, 0.4611, 0.4629,
1397 0.4646, 0.4663, 0.4680, 0.4697, 0.4714, 0.4730, 0.4747, 0.4759,
1398 0.4769, 0.4780, 0.4790, 0.4800, 0.4811, 0.4827, 0.4843, 0.4859,
1399 0.4874, 0.4889, 0.4905, 0.4920, 0.4935, 0.4950, 0.4965, 0.4979
1400 }
1401 };
1402
1403 if (fLUT) {
1404 delete [] fLUT;
1405 }
1406 fLUT = new Double_t[fLUTbin];
1407
1408 for (Int_t iplan = 0; iplan < AliTRDgeometry::kNplan; iplan++) {
1409 for (Int_t ilut = 0; ilut < kNlut; ilut++ ) {
1410 fLUT[iplan*kNlut+ilut] = lut[iplan][ilut];
1411 }
1412 }
1413
1414}
1415
1416//_____________________________________________________________________________
1417Double_t AliTRDclusterizer::LUTposition(Int_t iplane, Double_t ampL
1418 , Double_t ampC, Double_t ampR) const
1419{
1420 //
1421 // Calculates the cluster position using the lookup table.
1422 // Method provided by Bogdan Vulpescu.
1423 //
1424
1425 const Int_t kNlut = 128;
1426
1427 Double_t pos;
1428 Double_t x = 0.0;
1429 Double_t xmin;
1430 Double_t xmax;
1431 Double_t xwid;
1432
1433 Int_t side = 0;
1434 Int_t ix;
1435
1436 Double_t xMin[AliTRDgeometry::kNplan] = { 0.006492, 0.006377, 0.006258
1437 , 0.006144, 0.006030, 0.005980 };
1438 Double_t xMax[AliTRDgeometry::kNplan] = { 0.960351, 0.965870, 0.970445
1439 , 0.974352, 0.977667, 0.996101 };
1440
1441 if (ampL > ampR) {
1442 x = (ampL - ampR) / ampC;
1443 side = -1;
1444 }
1445 else if (ampL < ampR) {
1446 x = (ampR - ampL) / ampC;
1447 side = +1;
1448 }
1449
1450 if (ampL != ampR) {
1451
1452 xmin = xMin[iplane] + 0.000005;
1453 xmax = xMax[iplane] - 0.000005;
1454 xwid = (xmax - xmin) / 127.0;
1455
1456 if (x < xmin) {
1457 pos = 0.0000;
1458 }
1459 else if (x > xmax) {
1460 pos = side * 0.5000;
1461 }
1462 else {
1463 ix = (Int_t) ((x - xmin) / xwid);
1464 pos = side * fLUT[iplane*kNlut+ix];
1465 }
1466
1467 }
1468 else {
1469
1470 pos = 0.0;
1471
1472 }
1473
1474 return pos;
1475
1476}