1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // TRD cluster finder //
22 ///////////////////////////////////////////////////////////////////////////////
28 #include <TClonesArray.h>
29 #include <TObjArray.h>
31 #include "AliRunLoader.h"
32 #include "AliLoader.h"
33 #include "AliRawReader.h"
35 #include "AliAlignObj.h"
37 #include "AliTRDclusterizer.h"
38 #include "AliTRDcluster.h"
39 #include "AliTRDReconstructor.h"
40 #include "AliTRDgeometry.h"
41 #include "AliTRDarrayDictionary.h"
42 #include "AliTRDarrayADC.h"
43 #include "AliTRDdigitsManager.h"
44 #include "AliTRDrawData.h"
45 #include "AliTRDcalibDB.h"
46 #include "AliTRDrecoParam.h"
47 #include "AliTRDCommonParam.h"
48 #include "AliTRDtransform.h"
49 #include "AliTRDSignalIndex.h"
50 #include "AliTRDrawStreamBase.h"
51 #include "AliTRDfeeParam.h"
53 #include "TTreeStream.h"
55 #include "Cal/AliTRDCalROC.h"
56 #include "Cal/AliTRDCalDet.h"
57 #include "Cal/AliTRDCalSingleChamberStatus.h"
59 ClassImp(AliTRDclusterizer)
61 //_____________________________________________________________________________
62 AliTRDclusterizer::AliTRDclusterizer(const AliTRDReconstructor *const rec)
69 ,fDigitsManager(new AliTRDdigitsManager())
70 ,fTrackletContainer(NULL)
72 ,fTransform(new AliTRDtransform(0))
81 ,fMinLeftRightCutSigma(0)
87 ,fCalGainFactorROC(NULL)
88 ,fCalGainFactorDetValue(0)
91 ,fCalPadStatusROC(NULL)
97 // AliTRDclusterizer default constructor
100 SetBit(kAddLabels, kTRUE);
102 AliTRDcalibDB *trd = 0x0;
103 if (!(trd = AliTRDcalibDB::Instance())) {
104 AliFatal("Could not get calibration object");
107 fRawVersion = AliTRDfeeParam::Instance()->GetRAWversion();
109 // Initialize debug stream
111 if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kClusterizer) > 1){
112 TDirectory *savedir = gDirectory;
113 //fgGetDebugStream = new TTreeSRedirector("TRD.ClusterizerDebug.root");
120 //_____________________________________________________________________________
121 AliTRDclusterizer::AliTRDclusterizer(const Text_t *name, const Text_t *title, const AliTRDReconstructor *const rec)
128 ,fDigitsManager(new AliTRDdigitsManager())
129 ,fTrackletContainer(NULL)
131 ,fTransform(new AliTRDtransform(0))
140 ,fMinLeftRightCutSigma(0)
146 ,fCalGainFactorROC(NULL)
147 ,fCalGainFactorDetValue(0)
149 ,fCalNoiseDetValue(0)
150 ,fCalPadStatusROC(NULL)
156 // AliTRDclusterizer constructor
159 SetBit(kAddLabels, kTRUE);
161 AliTRDcalibDB *trd = 0x0;
162 if (!(trd = AliTRDcalibDB::Instance())) {
163 AliFatal("Could not get calibration object");
166 fDigitsManager->CreateArrays();
168 fRawVersion = AliTRDfeeParam::Instance()->GetRAWversion();
174 //_____________________________________________________________________________
175 AliTRDclusterizer::AliTRDclusterizer(const AliTRDclusterizer &c)
177 ,fReconstructor(c.fReconstructor)
182 ,fDigitsManager(NULL)
183 ,fTrackletContainer(NULL)
194 ,fMinLeftRightCutSigma(0)
200 ,fCalGainFactorROC(NULL)
201 ,fCalGainFactorDetValue(0)
203 ,fCalNoiseDetValue(0)
204 ,fCalPadStatusROC(NULL)
210 // AliTRDclusterizer copy constructor
213 SetBit(kAddLabels, kTRUE);
219 //_____________________________________________________________________________
220 AliTRDclusterizer::~AliTRDclusterizer()
223 // AliTRDclusterizer destructor
226 if (fRecPoints/* && IsClustersOwner()*/){
227 fRecPoints->Delete();
231 if (fDigitsManager) {
232 delete fDigitsManager;
233 fDigitsManager = NULL;
236 if (fTrackletContainer){
237 delete fTrackletContainer;
238 fTrackletContainer = NULL;
253 //_____________________________________________________________________________
254 AliTRDclusterizer &AliTRDclusterizer::operator=(const AliTRDclusterizer &c)
257 // Assignment operator
262 ((AliTRDclusterizer &) c).Copy(*this);
269 //_____________________________________________________________________________
270 void AliTRDclusterizer::Copy(TObject &c) const
276 ((AliTRDclusterizer &) c).fClusterTree = NULL;
277 ((AliTRDclusterizer &) c).fRecPoints = NULL;
278 ((AliTRDclusterizer &) c).fTrackletTree = NULL;
279 ((AliTRDclusterizer &) c).fDigitsManager = NULL;
280 ((AliTRDclusterizer &) c).fTrackletContainer = NULL;
281 ((AliTRDclusterizer &) c).fRawVersion = fRawVersion;
282 ((AliTRDclusterizer &) c).fTransform = NULL;
283 ((AliTRDclusterizer &) c).fLUTbin = 0;
284 ((AliTRDclusterizer &) c).fLUT = NULL;
285 ((AliTRDclusterizer &) c).fDigits = NULL;
286 ((AliTRDclusterizer &) c).fIndexes = NULL;
287 ((AliTRDclusterizer &) c).fADCthresh = 0;
288 ((AliTRDclusterizer &) c).fMaxThresh = 0;
289 ((AliTRDclusterizer &) c).fSigThresh = 0;
290 ((AliTRDclusterizer &) c).fMinMaxCutSigma= 0;
291 ((AliTRDclusterizer &) c).fMinLeftRightCutSigma = 0;
292 ((AliTRDclusterizer &) c).fLayer = 0;
293 ((AliTRDclusterizer &) c).fDet = 0;
294 ((AliTRDclusterizer &) c).fVolid = 0;
295 ((AliTRDclusterizer &) c).fColMax = 0;
296 ((AliTRDclusterizer &) c).fTimeTotal = 0;
297 ((AliTRDclusterizer &) c).fCalGainFactorROC = NULL;
298 ((AliTRDclusterizer &) c).fCalGainFactorDetValue = 0;
299 ((AliTRDclusterizer &) c).fCalNoiseROC = NULL;
300 ((AliTRDclusterizer &) c).fCalNoiseDetValue = 0;
301 ((AliTRDclusterizer &) c).fCalPadStatusROC = NULL;
302 ((AliTRDclusterizer &) c).fClusterROC = 0;
303 ((AliTRDclusterizer &) c).firstClusterROC= 0;
304 ((AliTRDclusterizer &) c).fNoOfClusters = 0;
307 //_____________________________________________________________________________
308 Bool_t AliTRDclusterizer::Open(const Char_t *name, Int_t nEvent)
311 // Opens the AliROOT file. Output and input are in the same file
314 TString evfoldname = AliConfig::GetDefaultEventFolderName();
315 fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
318 fRunLoader = AliRunLoader::Open(name);
322 AliError(Form("Can not open session for file %s.",name));
333 //_____________________________________________________________________________
334 Bool_t AliTRDclusterizer::OpenOutput()
337 // Open the output file
340 if (!fReconstructor->IsWritingClusters()) return kTRUE;
342 TObjArray *ioArray = 0x0;
344 AliLoader* loader = fRunLoader->GetLoader("TRDLoader");
345 loader->MakeTree("R");
347 fClusterTree = loader->TreeR();
348 fClusterTree->Branch("TRDcluster", "TObjArray", &ioArray, 32000, 0);
354 //_____________________________________________________________________________
355 Bool_t AliTRDclusterizer::OpenOutput(TTree *clusterTree)
358 // Connect the output tree
362 if (fReconstructor->IsWritingClusters()){
363 TObjArray *ioArray = 0x0;
364 fClusterTree = clusterTree;
365 fClusterTree->Branch("TRDcluster", "TObjArray", &ioArray, 32000, 0);
369 if (fReconstructor->IsWritingTracklets()){
370 TString evfoldname = AliConfig::GetDefaultEventFolderName();
371 fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
374 fRunLoader = AliRunLoader::Open("galice.root");
377 AliError(Form("Can not open session for file galice.root."));
381 UInt_t **leaves = new UInt_t *[2];
382 AliDataLoader *dl = fRunLoader->GetLoader("TRDLoader")->GetDataLoader("tracklets");
384 AliError("Could not get the tracklets data loader!");
385 dl = new AliDataLoader("TRD.Tracklets.root","tracklets", "tracklets");
386 fRunLoader->GetLoader("TRDLoader")->AddDataLoader(dl);
389 fTrackletTree = dl->Tree();
393 fTrackletTree = dl->Tree();
395 TBranch *trkbranch = fTrackletTree->GetBranch("trkbranch");
397 fTrackletTree->Branch("trkbranch",leaves[0],"det/i:side/i:tracklets[256]/i");
405 //_____________________________________________________________________________
406 Bool_t AliTRDclusterizer::OpenInput(Int_t nEvent)
409 // Opens a ROOT-file with TRD-hits and reads in the digits-tree
412 // Import the Trees for the event nEvent in the file
413 fRunLoader->GetEvent(nEvent);
419 //_____________________________________________________________________________
420 Bool_t AliTRDclusterizer::WriteClusters(Int_t det)
423 // Fills TRDcluster branch in the tree with the clusters
424 // found in detector = det. For det=-1 writes the tree.
428 (det >= AliTRDgeometry::Ndet())) {
429 AliError(Form("Unexpected detector index %d.\n",det));
433 TObjArray *ioArray = new TObjArray(400);
434 TBranch *branch = fClusterTree->GetBranch("TRDcluster");
436 branch = fClusterTree->Branch("TRDcluster","TObjArray",&ioArray,32000,0);
437 } else branch->SetAddress(&ioArray);
439 Int_t nRecPoints = RecPoints()->GetEntriesFast();
441 for (Int_t i = 0; i < nRecPoints; i++) {
442 AliTRDcluster *c = (AliTRDcluster *) RecPoints()->UncheckedAt(i);
443 if(det != c->GetDetector()) continue;
446 fClusterTree->Fill();
450 for (Int_t i = 0; i < nRecPoints; i++) {
451 AliTRDcluster *c = (AliTRDcluster *) RecPoints()->UncheckedAt(i);
452 if(c->GetDetector() != detOld){
453 fClusterTree->Fill();
455 detOld = c->GetDetector();
466 //_____________________________________________________________________________
467 Bool_t AliTRDclusterizer::WriteTracklets(Int_t det)
470 // Write the raw data tracklets into seperate file
473 UInt_t **leaves = new UInt_t *[2];
474 for (Int_t i=0; i<2 ;i++){
475 leaves[i] = new UInt_t[258];
476 leaves[i][0] = det; // det
477 leaves[i][1] = i; // side
478 memcpy(leaves[i]+2, fTrackletContainer[i], sizeof(UInt_t) * 256);
482 AliDataLoader *dl = fRunLoader->GetLoader("TRDLoader")->GetDataLoader("tracklets");
484 fTrackletTree = dl->Tree();
487 TBranch *trkbranch = fTrackletTree->GetBranch("trkbranch");
489 trkbranch = fTrackletTree->Branch("trkbranch",leaves[0],"det/i:side/i:tracklets[256]/i");
492 for (Int_t i=0; i<2; i++){
493 if (leaves[i][2]>0) {
494 trkbranch->SetAddress(leaves[i]);
495 fTrackletTree->Fill();
499 AliDataLoader *dl = fRunLoader->GetLoader("TRDLoader")->GetDataLoader("tracklets");
500 dl->WriteData("OVERWRITE");
508 //_____________________________________________________________________________
509 Bool_t AliTRDclusterizer::ReadDigits()
512 // Reads the digits arrays from the input aliroot file
516 AliError("No run loader available");
520 AliLoader* loader = fRunLoader->GetLoader("TRDLoader");
521 if (!loader->TreeD()) {
522 loader->LoadDigits();
525 // Read in the digit arrays
526 return (fDigitsManager->ReadDigits(loader->TreeD()));
530 //_____________________________________________________________________________
531 Bool_t AliTRDclusterizer::ReadDigits(TTree *digitsTree)
534 // Reads the digits arrays from the input tree
537 // Read in the digit arrays
538 return (fDigitsManager->ReadDigits(digitsTree));
542 //_____________________________________________________________________________
543 Bool_t AliTRDclusterizer::ReadDigits(AliRawReader *rawReader)
546 // Reads the digits arrays from the ddl file
550 fDigitsManager = raw.Raw2Digits(rawReader);
556 //_____________________________________________________________________________
557 Bool_t AliTRDclusterizer::MakeClusters()
560 // Creates clusters from digits
563 // Propagate info from the digits manager
564 if (TestBit(kAddLabels)){
565 SetBit(kAddLabels, fDigitsManager->UsesDictionaries());
568 Bool_t fReturn = kTRUE;
569 for (Int_t i = 0; i < AliTRDgeometry::kNdet; i++){
571 AliTRDarrayADC *digitsIn = (AliTRDarrayADC*) fDigitsManager->GetDigits(i); //mod
572 // This is to take care of switched off super modules
573 if (!digitsIn->HasData()) continue;
575 digitsIn->DeleteNegatives(); // Restore digits array to >=0 values
576 AliTRDSignalIndex* indexes = fDigitsManager->GetIndexes(i);
577 if (indexes->IsAllocated() == kFALSE){
578 fDigitsManager->BuildIndexes(i);
582 if (indexes->HasEntry()){
583 if (TestBit(kAddLabels)){
584 for (Int_t iDict = 0; iDict < AliTRDdigitsManager::kNDict; iDict++){
585 AliTRDarrayDictionary *tracksIn = 0; //mod
586 tracksIn = (AliTRDarrayDictionary *) fDigitsManager->GetDictionary(i,iDict); //mod
590 fR = MakeClusters(i);
591 fReturn = fR && fReturn;
595 // if(IsWritingClusters()) WriteClusters(i);
599 // No compress just remove
600 fDigitsManager->RemoveDigits(i);
601 fDigitsManager->RemoveDictionaries(i);
602 fDigitsManager->ClearIndexes(i);
605 if(fReconstructor->IsWritingClusters()) WriteClusters(-1);
607 AliInfo(Form("Number of found clusters : %d", RecPoints()->GetEntriesFast()));
613 //_____________________________________________________________________________
614 Bool_t AliTRDclusterizer::Raw2Clusters(AliRawReader *rawReader)
617 // Creates clusters from raw data
620 return Raw2ClustersChamber(rawReader);
624 //_____________________________________________________________________________
625 Bool_t AliTRDclusterizer::Raw2ClustersChamber(AliRawReader *rawReader)
628 // Creates clusters from raw data
631 // Create the digits manager
632 if (!fDigitsManager){
633 fDigitsManager = new AliTRDdigitsManager(kTRUE);
634 fDigitsManager->CreateArrays();
637 fDigitsManager->SetUseDictionaries(TestBit(kAddLabels));
639 // tracklet container for raw tracklet writing
640 if (!fTrackletContainer && fReconstructor->IsWritingTracklets()) {
641 // maximum tracklets for one HC
642 const Int_t kTrackletChmb=256;
643 fTrackletContainer = new UInt_t *[2];
644 fTrackletContainer[0] = new UInt_t[kTrackletChmb];
645 fTrackletContainer[1] = new UInt_t[kTrackletChmb];
648 AliTRDrawStreamBase *input = AliTRDrawStreamBase::GetRawStream(rawReader);
650 AliInfo(Form("Stream version: %s", input->IsA()->GetName()));
653 while ((det = input->NextChamber(fDigitsManager,fTrackletContainer)) >= 0){
654 Bool_t iclusterBranch = kFALSE;
655 if (fDigitsManager->GetIndexes(det)->HasEntry()){
656 iclusterBranch = MakeClusters(det);
659 fDigitsManager->ResetArrays(det);
661 if (!fReconstructor->IsWritingTracklets()) continue;
662 if (*(fTrackletContainer[0]) > 0 || *(fTrackletContainer[1]) > 0) WriteTracklets(det);
665 if (fReconstructor->IsWritingTracklets()){
666 delete [] fTrackletContainer[0];
667 delete [] fTrackletContainer[1];
668 delete [] fTrackletContainer;
669 fTrackletContainer = NULL;
672 if(fReconstructor->IsWritingClusters()) WriteClusters(-1);
674 delete fDigitsManager;
675 fDigitsManager = NULL;
680 AliInfo(Form("Number of found clusters : %d", fNoOfClusters));
685 //_____________________________________________________________________________
686 UChar_t AliTRDclusterizer::GetStatus(Short_t &signal)
689 // Check if a pad is masked
694 if(signal>0 && TESTBIT(signal, 10)){
696 for(int ibit=0; ibit<4; ibit++){
697 if(TESTBIT(signal, 11+ibit)){
698 SETBIT(status, ibit);
699 CLRBIT(signal, 11+ibit);
706 //_____________________________________________________________________________
707 void AliTRDclusterizer::SetPadStatus(const UChar_t status, UChar_t &out){
709 // Set the pad status into out
710 // First three bits are needed for the position encoding
715 //_____________________________________________________________________________
716 UChar_t AliTRDclusterizer::GetPadStatus(UChar_t encoding) const {
718 // return the staus encoding of the corrupted pad
720 return static_cast<UChar_t>(encoding >> 3);
723 //_____________________________________________________________________________
724 Int_t AliTRDclusterizer::GetCorruption(UChar_t encoding) const {
726 // Return the position of the corruption
731 //_____________________________________________________________________________
732 Bool_t AliTRDclusterizer::MakeClusters(Int_t det)
735 // Generates the cluster.
739 // digits should be expanded beforehand!
740 // digitsIn->Expand();
741 fDigits = (AliTRDarrayADC *) fDigitsManager->GetDigits(det); //mod
743 // This is to take care of switched off super modules
744 if (!fDigits->HasData())
749 fIndexes = fDigitsManager->GetIndexes(det);
750 if (fIndexes->IsAllocated() == kFALSE)
752 AliError("Indexes do not exist!");
756 AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
759 AliFatal("No AliTRDcalibDB instance available\n");
765 if (!fReconstructor){
766 AliError("Reconstructor not set\n");
770 TTreeSRedirector *fDebugStream = fReconstructor->GetDebugStream(AliTRDReconstructor::kClusterizer);
772 fMaxThresh = fReconstructor->GetRecoParam()->GetClusMaxThresh();
773 fSigThresh = fReconstructor->GetRecoParam()->GetClusSigThresh();
774 fMinMaxCutSigma = fReconstructor->GetRecoParam()->GetMinMaxCutSigma();
775 fMinLeftRightCutSigma = fReconstructor->GetRecoParam()->GetMinLeftRightCutSigma();
777 Int_t istack = fIndexes->GetStack();
778 fLayer = fIndexes->GetLayer();
779 Int_t isector = fIndexes->GetSM();
781 // Start clustering in the chamber
783 fDet = AliTRDgeometry::GetDetector(fLayer,istack,isector);
785 AliError("Strange Detector number Missmatch!");
789 // TRD space point transformation
790 fTransform->SetDetector(det);
792 Int_t iGeoLayer = AliGeomManager::kTRD1 + fLayer;
793 Int_t iGeoModule = istack + AliTRDgeometry::Nstack() * isector;
794 fVolid = AliGeomManager::LayerToVolUID(iGeoLayer,iGeoModule);
796 fColMax = fDigits->GetNcol();
797 //Int_t nRowMax = fDigits->GetNrow();
798 fTimeTotal = fDigits->GetNtime();
800 // Detector wise calibration object for the gain factors
801 const AliTRDCalDet *calGainFactorDet = calibration->GetGainFactorDet();
802 // Calibration object with pad wise values for the gain factors
803 fCalGainFactorROC = calibration->GetGainFactorROC(fDet);
804 // Calibration value for chamber wise gain factor
805 fCalGainFactorDetValue = calGainFactorDet->GetValue(fDet);
807 // Detector wise calibration object for the noise
808 const AliTRDCalDet *calNoiseDet = calibration->GetNoiseDet();
809 // Calibration object with pad wise values for the noise
810 fCalNoiseROC = calibration->GetNoiseROC(fDet);
811 // Calibration value for chamber wise noise
812 fCalNoiseDetValue = calNoiseDet->GetValue(fDet);
814 // Calibration object with the pad status
815 fCalPadStatusROC = calibration->GetPadStatusROC(fDet);
817 SetBit(kIsLUT, fReconstructor->GetRecoParam()->IsLUT());
818 SetBit(kIsHLT, fReconstructor->IsHLT());
820 firstClusterROC = -1;
823 if(fReconstructor->GetRecoParam()->IsTailCancelation()){
824 // Apply the gain and the tail cancelation via digital filter
828 MaxStruct curr, last;
829 Int_t nMaximas = 0, nCorrupted = 0;
831 // Here the clusterfining is happening
833 for(curr.Time = 0; curr.Time < fTimeTotal; curr.Time++)
834 while(fIndexes->NextRCIndex(curr.Row, curr.Col))
835 if(IsMaximum(curr, curr.padStatus, &curr.Signals[0]))
839 if(curr.Time==last.Time && curr.Row==last.Row && curr.Col==last.Col+2)
840 FivePadCluster(last, curr);
843 last=curr; curr.FivePad=kFALSE;
848 if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kClusterizer) > 2){
849 (*fDebugStream) << "MakeClusters"
850 << "Detector=" << det
851 << "NMaxima=" << nMaximas
852 << "NClusters=" << fClusterROC
853 << "NCorrupted=" << nCorrupted
857 if (TestBit(kAddLabels)) {
858 AddLabels(fDet,firstClusterROC,fClusterROC);
865 //_____________________________________________________________________________
866 Bool_t AliTRDclusterizer::IsMaximum(const MaxStruct &Max, UChar_t &padStatus, Short_t *const Signals)
869 // Returns true if this row,col,time combination is a maximum.
870 // Gives back the padStatus and the signals of the center pad and the two neighbouring pads.
873 Signals[1] = fDigits->GetData(Max.Row, Max.Col, Max.Time);
874 if(Signals[1] < fMaxThresh) return kFALSE;
876 Float_t noiseMiddleThresh = fMinMaxCutSigma*fCalNoiseDetValue*fCalNoiseROC->GetValue(Max.Col, Max.Row);
877 if (Signals[1] < noiseMiddleThresh) return kFALSE;
879 if (Max.Col + 1 >= fColMax || Max.Col < 1) return kFALSE;
881 UChar_t status[3]={fCalPadStatusROC->GetStatus(Max.Col-1, Max.Row),
882 fCalPadStatusROC->GetStatus(Max.Col, Max.Row),
883 fCalPadStatusROC->GetStatus(Max.Col+1, Max.Row)};
885 Signals[0] = fDigits->GetData(Max.Row, Max.Col-1, Max.Time);
886 Signals[2] = fDigits->GetData(Max.Row, Max.Col+1, Max.Time);
888 if(!(status[0] | status[1] | status[2])) {//all pads are good
889 if ((Signals[2] <= Signals[1]) && (Signals[0] < Signals[1])) {
890 if ((Signals[2] >= fSigThresh) || (Signals[0] >= fSigThresh)) {
891 Float_t noiseSumThresh = fMinLeftRightCutSigma
893 * fCalNoiseROC->GetValue(Max.Col, Max.Row);
894 if ((Signals[2]+Signals[0]+Signals[1]) < noiseSumThresh) return kFALSE;
900 else { // at least one of the pads is bad, and reject candidates with more than 1 problematic pad
901 if (status[2] && (!(status[0] || status[1])) && Signals[1] > Signals[0] && Signals[0] >= fSigThresh) {
903 SetPadStatus(status[2], padStatus);
906 else if (status[0] && (!(status[1] || status[2])) && Signals[1] >= Signals[2] && Signals[2] >= fSigThresh) {
908 SetPadStatus(status[0], padStatus);
911 else if (status[1] && (!(status[0] || status[2])) && ((Signals[2] >= fSigThresh) || (Signals[0] >= fSigThresh))) {
912 Signals[1]=fMaxThresh;
913 SetPadStatus(status[1], padStatus);
920 //_____________________________________________________________________________
921 Bool_t AliTRDclusterizer::FivePadCluster(MaxStruct &ThisMax, MaxStruct &NeighbourMax)
924 // Look for 5 pad cluster with minimum in the middle
925 // Gives back the ratio
928 if (ThisMax.Col >= fColMax - 3) return kFALSE;
929 if (ThisMax.Col < fColMax - 5){
930 if (fDigits->GetData(ThisMax.Row, ThisMax.Col+4, ThisMax.Time) >= fSigThresh)
933 if (ThisMax.Col > 1) {
934 if (fDigits->GetData(ThisMax.Row, ThisMax.Col-2, ThisMax.Time) >= fSigThresh)
938 //if (fSignalsThisMax[1] >= 0){ //TR: mod
940 const Float_t kEpsilon = 0.01;
941 Double_t padSignal[5] = {ThisMax.Signals[0], ThisMax.Signals[1], ThisMax.Signals[2],
942 NeighbourMax.Signals[1], NeighbourMax.Signals[2]};
944 // Unfold the two maxima and set the signal on
945 // the overlapping pad to the ratio
946 Float_t ratio = Unfold(kEpsilon,fLayer,padSignal);
947 ThisMax.Signals[2] = TMath::Nint(ThisMax.Signals[2]*ratio);
948 NeighbourMax.Signals[0] = TMath::Nint(NeighbourMax.Signals[0]*(1-ratio));
949 ThisMax.FivePad=kTRUE;
950 NeighbourMax.FivePad=kTRUE;
954 //_____________________________________________________________________________
955 void AliTRDclusterizer::CreateCluster(const MaxStruct &Max)
958 // Creates a cluster at the given position and saves it in fRecPoints
961 // The position of the cluster in COL direction relative to the center pad (pad units)
962 Double_t clusterPosCol = 0.0;
963 if (TestBit(kIsLUT)) {
964 // Calculate the position of the cluster by using the
965 // lookup table method
966 clusterPosCol = LUTposition(fLayer,Max.Signals[0]
971 // Calculate the position of the cluster by using the
972 // center of gravity method
973 const Int_t kNsig = 5;
974 Double_t padSignal[kNsig];
975 padSignal[1] = Max.Signals[0];
976 padSignal[2] = Max.Signals[1];
977 padSignal[3] = Max.Signals[2];
979 padSignal[0] = fDigits->GetData(Max.Row, Max.Col-2, Max.Time);
980 if(padSignal[0]>= padSignal[1])
983 if(Max.Col < fColMax - 3){
984 padSignal[4] = fDigits->GetData(Max.Row, Max.Col+2, Max.Time);
985 if(padSignal[4]>= padSignal[3])
988 clusterPosCol = GetCOG(padSignal);
991 // Count the number of pads in the cluster
993 Short_t signals[7] = { 0, 0, 0, 0, 0, 0, 0 };
995 if(!TestBit(kIsHLT))CalcAdditionalInfo(Max, signals, nPadCount);
997 // Transform the local cluster coordinates into calibrated
998 // space point positions defined in the local tracking system.
999 // Here the calibration for T0, Vdrift and ExB is applied as well.
1000 Double_t clusterXYZ[6];
1001 clusterXYZ[0] = clusterPosCol;
1002 clusterXYZ[1] = Max.Signals[2];
1003 clusterXYZ[2] = Max.Signals[1];
1004 clusterXYZ[3] = Max.Signals[0];
1005 clusterXYZ[4] = 0.0;
1006 clusterXYZ[5] = 0.0;
1007 Int_t clusterRCT[3];
1008 clusterRCT[0] = Max.Row;
1009 clusterRCT[1] = Max.Col;
1013 if (fTransform->Transform(clusterXYZ,clusterRCT,((UInt_t) Max.Time),out,0)) {
1015 Char_t clusterTimeBin = ((Char_t) clusterRCT[2]);
1016 Float_t clusterPos[3];
1017 clusterPos[0] = clusterXYZ[0];
1018 clusterPos[1] = clusterXYZ[1];
1019 clusterPos[2] = clusterXYZ[2];
1020 Float_t clusterSig[2];
1021 clusterSig[0] = clusterXYZ[4];
1022 clusterSig[1] = clusterXYZ[5];
1023 Float_t clusterCharge = clusterXYZ[3];
1025 AliTRDcluster cluster(
1027 clusterCharge, clusterPos, clusterSig,
1029 ((Char_t) nPadCount),
1031 ((UChar_t) Max.Col), ((UChar_t) Max.Row), ((UChar_t) Max.Time),
1032 clusterTimeBin, clusterPosCol,
1035 cluster.SetInChamber(!out);
1036 cluster.SetFivePad(Max.FivePad);
1038 UChar_t maskPosition = GetCorruption(Max.padStatus);
1040 cluster.SetPadMaskedPosition(maskPosition);
1041 cluster.SetPadMaskedStatus(GetPadStatus(Max.padStatus));
1044 // Temporarily store the Max.Row, column and time bin of the center pad
1045 // Used to later on assign the track indices
1046 cluster.SetLabel(Max.Row, 0);
1047 cluster.SetLabel(Max.Col, 1);
1048 cluster.SetLabel(Max.Time,2);
1050 AddClusterToArray(&cluster); //needs to be like that because HLT does things differently
1052 //AddCluster(Max,clusterXYZ,clusterTimeBin,signals,nPadCount,out,clusterPosCol);
1053 // Store the index of the first cluster in the current ROC
1054 if (firstClusterROC < 0) {
1055 firstClusterROC = fNoOfClusters;
1064 //_____________________________________________________________________________
1065 void AliTRDclusterizer::CalcAdditionalInfo(const MaxStruct &Max, Short_t *const signals, Int_t &nPadCount)
1067 // Look to the right
1069 while (fDigits->GetData(Max.Row, Max.Col-ii, Max.Time) >= fSigThresh) {
1072 if (Max.Col < ii) break;
1076 while (fDigits->GetData(Max.Row, Max.Col+ii, Max.Time) >= fSigThresh) {
1079 if (Max.Col+ii >= fColMax) break;
1082 // Store the amplitudes of the pads in the cluster for later analysis
1083 // and check whether one of these pads is masked in the database
1084 signals[2]=Max.Signals[0];
1085 signals[3]=Max.Signals[1];
1086 signals[4]=Max.Signals[2];
1087 for(Int_t i = 0; i<2; i++)
1090 signals[i] = fDigits->GetData(Max.Row, Max.Col-3+i, Max.Time);
1091 if(Max.Col+3-i < fColMax)
1092 signals[6-i] = fDigits->GetData(Max.Row, Max.Col+3-i, Max.Time);
1094 /*for (Int_t jPad = Max.Col-3; jPad <= Max.Col+3; jPad++) {
1095 if ((jPad >= 0) && (jPad < fColMax))
1096 signals[jPad-Max.Col+3] = TMath::Nint(fDigits->GetData(Max.Row,jPad,Max.Time));
1100 //_____________________________________________________________________________
1101 void AliTRDclusterizer::AddClusterToArray(AliTRDcluster *cluster)
1104 // Add a cluster to the array
1107 Int_t n = RecPoints()->GetEntriesFast();
1108 if(n!=fNoOfClusters)AliError(Form("fNoOfClusters != RecPoints()->GetEntriesFast %i != %i \n", fNoOfClusters, n));
1109 new((*RecPoints())[n]) AliTRDcluster(*cluster);
1112 //_____________________________________________________________________________
1113 Bool_t AliTRDclusterizer::AddLabels(const Int_t idet, const Int_t firstClusterROC, const Int_t nClusterROC)
1116 // Add the track indices to the found clusters
1119 const Int_t kNclus = 3;
1120 const Int_t kNdict = AliTRDdigitsManager::kNDict;
1121 const Int_t kNtrack = kNdict * kNclus;
1123 Int_t iClusterROC = 0;
1130 // Temporary array to collect the track indices
1131 Int_t *idxTracks = new Int_t[kNtrack*nClusterROC];
1133 // Loop through the dictionary arrays one-by-one
1134 // to keep memory consumption low
1135 AliTRDarrayDictionary *tracksIn = 0; //mod
1136 for (Int_t iDict = 0; iDict < kNdict; iDict++) {
1138 // tracksIn should be expanded beforehand!
1139 tracksIn = (AliTRDarrayDictionary *) fDigitsManager->GetDictionary(idet,iDict);
1141 // Loop though the clusters found in this ROC
1142 for (iClusterROC = 0; iClusterROC < nClusterROC; iClusterROC++) {
1144 AliTRDcluster *cluster = (AliTRDcluster *)
1145 RecPoints()->UncheckedAt(firstClusterROC+iClusterROC);
1146 row = cluster->GetLabel(0);
1147 col = cluster->GetLabel(1);
1148 time = cluster->GetLabel(2);
1150 for (iPad = 0; iPad < kNclus; iPad++) {
1151 Int_t iPadCol = col - 1 + iPad;
1152 Int_t index = tracksIn->GetData(row,iPadCol,time); //Modification of -1 in Track
1153 idxTracks[3*iPad+iDict + iClusterROC*kNtrack] = index;
1160 // Copy the track indices into the cluster
1161 // Loop though the clusters found in this ROC
1162 for (iClusterROC = 0; iClusterROC < nClusterROC; iClusterROC++) {
1164 AliTRDcluster *cluster = (AliTRDcluster *)
1165 RecPoints()->UncheckedAt(firstClusterROC+iClusterROC);
1166 cluster->SetLabel(-9999,0);
1167 cluster->SetLabel(-9999,1);
1168 cluster->SetLabel(-9999,2);
1170 cluster->AddTrackIndex(&idxTracks[iClusterROC*kNtrack]);
1174 delete [] idxTracks;
1180 //_____________________________________________________________________________
1181 Double_t AliTRDclusterizer::GetCOG(Double_t signal[5]) const
1185 // Used for clusters with more than 3 pads - where LUT not applicable
1188 Double_t sum = signal[0]
1195 // Go to 3 pad COG ????
1197 Double_t res = (0.0 * (-signal[0] + signal[4])
1198 + (-signal[1] + signal[3])) / sum;
1204 //_____________________________________________________________________________
1205 Float_t AliTRDclusterizer::Unfold(Double_t eps, Int_t layer, Double_t *padSignal) const
1208 // Method to unfold neighbouring maxima.
1209 // The charge ratio on the overlapping pad is calculated
1210 // until there is no more change within the range given by eps.
1211 // The resulting ratio is then returned to the calling method.
1214 AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
1216 AliError("No AliTRDcalibDB instance available\n");
1221 Int_t itStep = 0; // Count iteration steps
1223 Double_t ratio = 0.5; // Start value for ratio
1224 Double_t prevRatio = 0.0; // Store previous ratio
1226 Double_t newLeftSignal[3] = { 0.0, 0.0, 0.0 }; // Array to store left cluster signal
1227 Double_t newRightSignal[3] = { 0.0, 0.0, 0.0 }; // Array to store right cluster signal
1228 Double_t newSignal[3] = { 0.0, 0.0, 0.0 };
1230 // Start the iteration
1231 while ((TMath::Abs(prevRatio - ratio) > eps) && (itStep < 10)) {
1236 // Cluster position according to charge ratio
1237 Double_t maxLeft = (ratio*padSignal[2] - padSignal[0])
1238 / (padSignal[0] + padSignal[1] + ratio * padSignal[2]);
1239 Double_t maxRight = (padSignal[4] - (1-ratio)*padSignal[2])
1240 / ((1.0 - ratio)*padSignal[2] + padSignal[3] + padSignal[4]);
1242 // Set cluster charge ratio
1243 irc = calibration->PadResponse(1.0, maxLeft, layer, newSignal);
1244 Double_t ampLeft = padSignal[1] / newSignal[1];
1245 irc = calibration->PadResponse(1.0, maxRight, layer, newSignal);
1246 Double_t ampRight = padSignal[3] / newSignal[1];
1248 // Apply pad response to parameters
1249 irc = calibration->PadResponse(ampLeft ,maxLeft ,layer,newLeftSignal );
1250 irc = calibration->PadResponse(ampRight,maxRight,layer,newRightSignal);
1252 // Calculate new overlapping ratio
1253 ratio = TMath::Min((Double_t) 1.0
1254 ,newLeftSignal[2] / (newLeftSignal[2] + newRightSignal[0]));
1262 //_____________________________________________________________________________
1263 void AliTRDclusterizer::TailCancelation()
1266 // Applies the tail cancelation and gain factors:
1267 // Transform fDigits to fDigits
1274 Double_t *inADC = new Double_t[fTimeTotal]; // ADC data before tail cancellation
1275 Double_t *outADC = new Double_t[fTimeTotal]; // ADC data after tail cancellation
1277 fIndexes->ResetCounters();
1278 TTreeSRedirector *fDebugStream = fReconstructor->GetDebugStream(AliTRDReconstructor::kClusterizer);
1279 while(fIndexes->NextRCIndex(iRow, iCol))
1281 Float_t fCalGainFactorROCValue = fCalGainFactorROC->GetValue(iCol,iRow);
1282 Double_t gain = fCalGainFactorDetValue
1283 * fCalGainFactorROCValue;
1285 Bool_t corrupted = kFALSE;
1286 for (iTime = 0; iTime < fTimeTotal; iTime++)
1288 // Apply gain gain factor
1289 inADC[iTime] = fDigits->GetData(iRow,iCol,iTime);
1290 if (fCalPadStatusROC->GetStatus(iCol, iRow)) corrupted = kTRUE;
1291 inADC[iTime] /= gain;
1292 outADC[iTime] = inADC[iTime];
1293 if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kClusterizer) > 7){
1294 (*fDebugStream) << "TailCancellation"
1298 << "inADC=" << inADC[iTime]
1300 << "outADC=" << outADC[iTime]
1301 << "corrupted=" << corrupted
1307 // Apply the tail cancelation via the digital filter
1308 // (only for non-coorupted pads)
1309 DeConvExp(&inADC[0],&outADC[0],fTimeTotal,fReconstructor->GetRecoParam() ->GetTCnexp());
1312 for(iTime = 0; iTime < fTimeTotal; iTime++)//while (fIndexes->NextTbinIndex(iTime))
1314 // Store the amplitude of the digit if above threshold
1315 if (outADC[iTime] > fADCthresh)
1316 fDigits->SetData(iRow,iCol,iTime,TMath::Nint(outADC[iTime]));
1318 fDigits->SetData(iRow,iCol,iTime,0);
1321 } // while irow icol
1330 //_____________________________________________________________________________
1331 void AliTRDclusterizer::DeConvExp(const Double_t *const source, Double_t *const target
1332 ,const Int_t n, const Int_t nexp)
1335 // Tail cancellation by deconvolution for PASA v4 TRF
1339 Double_t coefficients[2];
1341 // Initialization (coefficient = alpha, rates = lambda)
1347 if (nexp == 1) { // 1 Exponentials
1353 if (nexp == 2) { // 2 Exponentials
1355 fReconstructor->GetTCParams(par);
1356 r1 = par[0];//1.156;
1357 r2 = par[1];//0.130;
1358 c1 = par[2];//0.114;
1359 c2 = par[3];//0.624;
1362 coefficients[0] = c1;
1363 coefficients[1] = c2;
1367 rates[0] = TMath::Exp(-dt/(r1));
1368 rates[1] = TMath::Exp(-dt/(r2));
1373 Double_t reminder[2];
1374 Double_t correction = 0.0;
1375 Double_t result = 0.0;
1377 // Attention: computation order is important
1378 for (k = 0; k < nexp; k++) {
1382 for (i = 0; i < n; i++) {
1384 result = (source[i] - correction); // No rescaling
1387 for (k = 0; k < nexp; k++) {
1388 reminder[k] = rates[k] * (reminder[k] + coefficients[k] * result);
1392 for (k = 0; k < nexp; k++) {
1393 correction += reminder[k];
1400 //_____________________________________________________________________________
1401 void AliTRDclusterizer::ResetRecPoints()
1404 // Resets the list of rec points
1408 fRecPoints->Delete();
1413 //_____________________________________________________________________________
1414 TClonesArray *AliTRDclusterizer::RecPoints()
1417 // Returns the list of rec points
1421 if(!(fRecPoints = AliTRDReconstructor::GetClusters())){
1422 // determine number of clusters which has to be allocated
1423 Float_t nclusters = fReconstructor->GetRecoParam()->GetNClusters();
1425 fRecPoints = new TClonesArray("AliTRDcluster", Int_t(nclusters));
1427 //SetClustersOwner(kTRUE);
1428 AliTRDReconstructor::SetClusters(0x0);
1434 //_____________________________________________________________________________
1435 void AliTRDclusterizer::FillLUT()
1441 const Int_t kNlut = 128;
1443 fLUTbin = AliTRDgeometry::kNlayer * kNlut;
1445 // The lookup table from Bogdan
1446 Float_t lut[AliTRDgeometry::kNlayer][kNlut] = {
1448 0.0070, 0.0150, 0.0224, 0.0298, 0.0374, 0.0454, 0.0533, 0.0611,
1449 0.0684, 0.0755, 0.0827, 0.0900, 0.0975, 0.1049, 0.1120, 0.1187,
1450 0.1253, 0.1318, 0.1385, 0.1453, 0.1519, 0.1584, 0.1646, 0.1704,
1451 0.1762, 0.1821, 0.1879, 0.1938, 0.1996, 0.2053, 0.2108, 0.2160,
1452 0.2210, 0.2260, 0.2310, 0.2361, 0.2411, 0.2461, 0.2509, 0.2557,
1453 0.2602, 0.2646, 0.2689, 0.2732, 0.2774, 0.2816, 0.2859, 0.2901,
1454 0.2942, 0.2983, 0.3022, 0.3061, 0.3099, 0.3136, 0.3172, 0.3207,
1455 0.3242, 0.3278, 0.3312, 0.3347, 0.3382, 0.3416, 0.3450, 0.3483,
1456 0.3515, 0.3547, 0.3579, 0.3609, 0.3639, 0.3669, 0.3698, 0.3727,
1457 0.3756, 0.3785, 0.3813, 0.3842, 0.3870, 0.3898, 0.3926, 0.3952,
1458 0.3979, 0.4005, 0.4032, 0.4057, 0.4082, 0.4108, 0.4132, 0.4157,
1459 0.4181, 0.4205, 0.4228, 0.4252, 0.4275, 0.4299, 0.4322, 0.4345,
1460 0.4367, 0.4390, 0.4412, 0.4434, 0.4456, 0.4478, 0.4499, 0.4520,
1461 0.4541, 0.4562, 0.4583, 0.4603, 0.4623, 0.4643, 0.4663, 0.4683,
1462 0.4702, 0.4722, 0.4741, 0.4758, 0.4774, 0.4790, 0.4805, 0.4824,
1463 0.4844, 0.4863, 0.4883, 0.4902, 0.4921, 0.4940, 0.4959, 0.4978
1466 0.0072, 0.0156, 0.0235, 0.0313, 0.0394, 0.0478, 0.0561, 0.0642,
1467 0.0718, 0.0792, 0.0868, 0.0947, 0.1025, 0.1101, 0.1172, 0.1241,
1468 0.1309, 0.1378, 0.1449, 0.1518, 0.1586, 0.1650, 0.1710, 0.1770,
1469 0.1830, 0.1891, 0.1952, 0.2011, 0.2070, 0.2125, 0.2177, 0.2229,
1470 0.2280, 0.2332, 0.2383, 0.2435, 0.2484, 0.2533, 0.2581, 0.2627,
1471 0.2670, 0.2714, 0.2757, 0.2799, 0.2842, 0.2884, 0.2927, 0.2968,
1472 0.3008, 0.3048, 0.3086, 0.3123, 0.3159, 0.3195, 0.3231, 0.3266,
1473 0.3301, 0.3335, 0.3370, 0.3404, 0.3438, 0.3471, 0.3504, 0.3536,
1474 0.3567, 0.3598, 0.3628, 0.3657, 0.3686, 0.3715, 0.3744, 0.3772,
1475 0.3800, 0.3828, 0.3856, 0.3884, 0.3911, 0.3938, 0.3965, 0.3991,
1476 0.4016, 0.4042, 0.4067, 0.4092, 0.4116, 0.4140, 0.4164, 0.4187,
1477 0.4211, 0.4234, 0.4257, 0.4280, 0.4302, 0.4325, 0.4347, 0.4369,
1478 0.4391, 0.4413, 0.4434, 0.4456, 0.4477, 0.4497, 0.4518, 0.4538,
1479 0.4558, 0.4578, 0.4598, 0.4618, 0.4637, 0.4656, 0.4675, 0.4694,
1480 0.4713, 0.4732, 0.4750, 0.4766, 0.4781, 0.4797, 0.4813, 0.4832,
1481 0.4851, 0.4870, 0.4888, 0.4906, 0.4925, 0.4942, 0.4960, 0.4978
1484 0.0075, 0.0163, 0.0246, 0.0328, 0.0415, 0.0504, 0.0592, 0.0674,
1485 0.0753, 0.0832, 0.0914, 0.0996, 0.1077, 0.1154, 0.1225, 0.1296,
1486 0.1369, 0.1442, 0.1515, 0.1585, 0.1652, 0.1714, 0.1776, 0.1839,
1487 0.1902, 0.1965, 0.2025, 0.2085, 0.2141, 0.2194, 0.2247, 0.2299,
1488 0.2352, 0.2405, 0.2457, 0.2507, 0.2557, 0.2604, 0.2649, 0.2693,
1489 0.2737, 0.2780, 0.2823, 0.2867, 0.2909, 0.2951, 0.2992, 0.3033,
1490 0.3072, 0.3110, 0.3146, 0.3182, 0.3218, 0.3253, 0.3288, 0.3323,
1491 0.3357, 0.3392, 0.3426, 0.3459, 0.3492, 0.3524, 0.3555, 0.3586,
1492 0.3616, 0.3645, 0.3674, 0.3703, 0.3731, 0.3759, 0.3787, 0.3815,
1493 0.3843, 0.3870, 0.3897, 0.3925, 0.3950, 0.3976, 0.4002, 0.4027,
1494 0.4052, 0.4076, 0.4101, 0.4124, 0.4148, 0.4171, 0.4194, 0.4217,
1495 0.4239, 0.4262, 0.4284, 0.4306, 0.4328, 0.4350, 0.4371, 0.4393,
1496 0.4414, 0.4435, 0.4455, 0.4476, 0.4496, 0.4516, 0.4536, 0.4555,
1497 0.4575, 0.4594, 0.4613, 0.4632, 0.4650, 0.4669, 0.4687, 0.4705,
1498 0.4723, 0.4741, 0.4758, 0.4773, 0.4789, 0.4804, 0.4821, 0.4839,
1499 0.4857, 0.4875, 0.4893, 0.4910, 0.4928, 0.4945, 0.4961, 0.4978
1502 0.0078, 0.0171, 0.0258, 0.0345, 0.0438, 0.0532, 0.0624, 0.0708,
1503 0.0791, 0.0875, 0.0962, 0.1048, 0.1130, 0.1206, 0.1281, 0.1356,
1504 0.1432, 0.1508, 0.1582, 0.1651, 0.1716, 0.1780, 0.1845, 0.1910,
1505 0.1975, 0.2038, 0.2099, 0.2155, 0.2210, 0.2263, 0.2317, 0.2371,
1506 0.2425, 0.2477, 0.2528, 0.2578, 0.2626, 0.2671, 0.2715, 0.2759,
1507 0.2803, 0.2846, 0.2890, 0.2933, 0.2975, 0.3016, 0.3056, 0.3095,
1508 0.3132, 0.3168, 0.3204, 0.3239, 0.3274, 0.3309, 0.3344, 0.3378,
1509 0.3412, 0.3446, 0.3479, 0.3511, 0.3543, 0.3574, 0.3603, 0.3633,
1510 0.3662, 0.3690, 0.3718, 0.3747, 0.3774, 0.3802, 0.3829, 0.3857,
1511 0.3883, 0.3910, 0.3936, 0.3962, 0.3987, 0.4012, 0.4037, 0.4061,
1512 0.4085, 0.4109, 0.4132, 0.4155, 0.4177, 0.4200, 0.4222, 0.4244,
1513 0.4266, 0.4288, 0.4309, 0.4331, 0.4352, 0.4373, 0.4394, 0.4414,
1514 0.4435, 0.4455, 0.4475, 0.4494, 0.4514, 0.4533, 0.4552, 0.4571,
1515 0.4590, 0.4608, 0.4626, 0.4645, 0.4662, 0.4680, 0.4698, 0.4715,
1516 0.4733, 0.4750, 0.4766, 0.4781, 0.4796, 0.4812, 0.4829, 0.4846,
1517 0.4863, 0.4880, 0.4897, 0.4914, 0.4930, 0.4946, 0.4963, 0.4979
1520 0.0081, 0.0178, 0.0270, 0.0364, 0.0463, 0.0562, 0.0656, 0.0744,
1521 0.0831, 0.0921, 0.1013, 0.1102, 0.1183, 0.1261, 0.1339, 0.1419,
1522 0.1499, 0.1576, 0.1648, 0.1715, 0.1782, 0.1849, 0.1917, 0.1984,
1523 0.2048, 0.2110, 0.2167, 0.2223, 0.2278, 0.2333, 0.2389, 0.2444,
1524 0.2497, 0.2548, 0.2598, 0.2645, 0.2691, 0.2735, 0.2780, 0.2824,
1525 0.2868, 0.2912, 0.2955, 0.2997, 0.3038, 0.3078, 0.3116, 0.3152,
1526 0.3188, 0.3224, 0.3259, 0.3294, 0.3329, 0.3364, 0.3398, 0.3432,
1527 0.3465, 0.3497, 0.3529, 0.3561, 0.3591, 0.3620, 0.3649, 0.3677,
1528 0.3705, 0.3733, 0.3761, 0.3788, 0.3816, 0.3843, 0.3869, 0.3896,
1529 0.3922, 0.3948, 0.3973, 0.3998, 0.4022, 0.4047, 0.4070, 0.4094,
1530 0.4117, 0.4139, 0.4162, 0.4184, 0.4206, 0.4227, 0.4249, 0.4270,
1531 0.4291, 0.4313, 0.4334, 0.4354, 0.4375, 0.4395, 0.4415, 0.4435,
1532 0.4455, 0.4474, 0.4493, 0.4512, 0.4531, 0.4550, 0.4568, 0.4586,
1533 0.4604, 0.4622, 0.4639, 0.4657, 0.4674, 0.4691, 0.4708, 0.4725,
1534 0.4742, 0.4758, 0.4773, 0.4788, 0.4803, 0.4819, 0.4836, 0.4852,
1535 0.4869, 0.4885, 0.4901, 0.4917, 0.4933, 0.4948, 0.4964, 0.4979
1538 0.0085, 0.0189, 0.0288, 0.0389, 0.0497, 0.0603, 0.0699, 0.0792,
1539 0.0887, 0.0985, 0.1082, 0.1170, 0.1253, 0.1336, 0.1421, 0.1505,
1540 0.1587, 0.1662, 0.1733, 0.1803, 0.1874, 0.1945, 0.2014, 0.2081,
1541 0.2143, 0.2201, 0.2259, 0.2316, 0.2374, 0.2431, 0.2487, 0.2541,
1542 0.2593, 0.2642, 0.2689, 0.2735, 0.2781, 0.2826, 0.2872, 0.2917,
1543 0.2961, 0.3003, 0.3045, 0.3086, 0.3125, 0.3162, 0.3198, 0.3235,
1544 0.3270, 0.3306, 0.3342, 0.3377, 0.3411, 0.3446, 0.3479, 0.3511,
1545 0.3543, 0.3575, 0.3605, 0.3634, 0.3663, 0.3691, 0.3720, 0.3748,
1546 0.3775, 0.3803, 0.3830, 0.3857, 0.3884, 0.3911, 0.3937, 0.3962,
1547 0.3987, 0.4012, 0.4036, 0.4060, 0.4084, 0.4107, 0.4129, 0.4152,
1548 0.4174, 0.4196, 0.4218, 0.4239, 0.4261, 0.4282, 0.4303, 0.4324,
1549 0.4344, 0.4365, 0.4385, 0.4405, 0.4425, 0.4445, 0.4464, 0.4483,
1550 0.4502, 0.4521, 0.4539, 0.4558, 0.4576, 0.4593, 0.4611, 0.4629,
1551 0.4646, 0.4663, 0.4680, 0.4697, 0.4714, 0.4730, 0.4747, 0.4759,
1552 0.4769, 0.4780, 0.4790, 0.4800, 0.4811, 0.4827, 0.4843, 0.4859,
1553 0.4874, 0.4889, 0.4905, 0.4920, 0.4935, 0.4950, 0.4965, 0.4979
1560 fLUT = new Double_t[fLUTbin];
1562 for (Int_t ilayer = 0; ilayer < AliTRDgeometry::kNlayer; ilayer++) {
1563 for (Int_t ilut = 0; ilut < kNlut; ilut++ ) {
1564 fLUT[ilayer*kNlut+ilut] = lut[ilayer][ilut];
1570 //_____________________________________________________________________________
1571 Double_t AliTRDclusterizer::LUTposition(Int_t ilayer
1574 , Double_t ampR) const
1577 // Calculates the cluster position using the lookup table.
1578 // Method provided by Bogdan Vulpescu.
1581 const Int_t kNlut = 128;
1592 Double_t xMin[AliTRDgeometry::kNlayer] = { 0.006492, 0.006377, 0.006258
1593 , 0.006144, 0.006030, 0.005980 };
1594 Double_t xMax[AliTRDgeometry::kNlayer] = { 0.960351, 0.965870, 0.970445
1595 , 0.974352, 0.977667, 0.996101 };
1598 x = (ampL - ampR) / ampC;
1601 else if (ampL < ampR) {
1602 x = (ampR - ampL) / ampC;
1608 xmin = xMin[ilayer] + 0.000005;
1609 xmax = xMax[ilayer] - 0.000005;
1610 xwid = (xmax - xmin) / 127.0;
1615 else if (x > xmax) {
1616 pos = side * 0.5000;
1619 ix = (Int_t) ((x - xmin) / xwid);
1620 pos = side * fLUT[ilayer*kNlut+ix];