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