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