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