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