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