]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/AliTRDclusterizer.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
setting the name for non-std file
[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
66f6bfd9 24#include <TClonesArray.h>
6d50f529 25#include <TObjArray.h>
f7336fa3 26
88cb7938 27#include "AliRunLoader.h"
28#include "AliLoader.h"
3fe61b77 29#include "AliAlignObj.h"
88cb7938 30
f7336fa3 31#include "AliTRDclusterizer.h"
3e1a3ad8 32#include "AliTRDcluster.h"
fc546d21 33#include "AliTRDReconstructor.h"
793ff80c 34#include "AliTRDgeometry.h"
b65e5048 35#include "AliTRDarrayDictionary.h"
36#include "AliTRDarrayADC.h"
3fe61b77 37#include "AliTRDdigitsManager.h"
e570c2a5 38#include "AliTRDdigitsParam.h"
3fe61b77 39#include "AliTRDrawData.h"
3551db50 40#include "AliTRDcalibDB.h"
3fe61b77 41#include "AliTRDtransform.h"
42#include "AliTRDSignalIndex.h"
dfbb4bb9 43#include "AliTRDrawStreamBase.h"
3fe61b77 44#include "AliTRDfeeParam.h"
a5b99acd 45#include "AliTRDtrackletWord.h"
3fe61b77 46
29f95561 47#include "TTreeStream.h"
48
3fe61b77 49#include "Cal/AliTRDCalROC.h"
50#include "Cal/AliTRDCalDet.h"
0e09df31 51#include "Cal/AliTRDCalSingleChamberStatus.h"
f7336fa3 52
53ClassImp(AliTRDclusterizer)
54
55//_____________________________________________________________________________
486c2339 56AliTRDclusterizer::AliTRDclusterizer(const AliTRDReconstructor *const rec)
6d50f529 57 :TNamed()
3a039a31 58 ,fReconstructor(rec)
6d50f529 59 ,fRunLoader(NULL)
60 ,fClusterTree(NULL)
61 ,fRecPoints(NULL)
a5b99acd 62 ,fTracklets(NULL)
9c7c9ec1 63 ,fTrackletTree(NULL)
3d0c7d6d 64 ,fDigitsManager(new AliTRDdigitsManager())
9c7c9ec1 65 ,fTrackletContainer(NULL)
3fe61b77 66 ,fRawVersion(2)
dfbb4bb9 67 ,fTransform(new AliTRDtransform(0))
1b3a719f 68 ,fDigits(NULL)
064d808d 69 ,fIndexes(NULL)
064d808d 70 ,fMaxThresh(0)
71 ,fSigThresh(0)
72 ,fMinMaxCutSigma(0)
73 ,fMinLeftRightCutSigma(0)
74 ,fLayer(0)
75 ,fDet(0)
76 ,fVolid(0)
77 ,fColMax(0)
78 ,fTimeTotal(0)
79 ,fCalGainFactorROC(NULL)
80 ,fCalGainFactorDetValue(0)
81 ,fCalNoiseROC(NULL)
82 ,fCalNoiseDetValue(0)
1b3a719f 83 ,fCalPadStatusROC(NULL)
064d808d 84 ,fClusterROC(0)
85 ,firstClusterROC(0)
3d0c7d6d 86 ,fNoOfClusters(0)
615f0826 87 ,fBaseline(0)
8dbc94c7 88 ,fRawStream(NULL)
f7336fa3 89{
90 //
91 // AliTRDclusterizer default constructor
92 //
93
b72f4eaf 94 SetBit(kLabels, kTRUE);
8cbe253a 95 SetBit(knewDM, kFALSE);
3d0c7d6d 96
a7ac01d2 97 AliTRDcalibDB *trd = 0x0;
98 if (!(trd = AliTRDcalibDB::Instance())) {
99 AliFatal("Could not get calibration object");
100 }
101
3fe61b77 102 fRawVersion = AliTRDfeeParam::Instance()->GetRAWversion();
103
3a039a31 104 // Initialize debug stream
105 if(fReconstructor){
a2fbb6ec 106 if(fReconstructor->GetRecoParam()->GetStreamLevel(AliTRDrecoParam::kClusterizer) > 1){
3a039a31 107 TDirectory *savedir = gDirectory;
29f95561 108 //fgGetDebugStream = new TTreeSRedirector("TRD.ClusterizerDebug.root");
3a039a31 109 savedir->cd();
110 }
111 }
eb52b657 112
f7336fa3 113}
114
115//_____________________________________________________________________________
486c2339 116AliTRDclusterizer::AliTRDclusterizer(const Text_t *name, const Text_t *title, const AliTRDReconstructor *const rec)
6d50f529 117 :TNamed(name,title)
3a039a31 118 ,fReconstructor(rec)
6d50f529 119 ,fRunLoader(NULL)
120 ,fClusterTree(NULL)
121 ,fRecPoints(NULL)
a5b99acd 122 ,fTracklets(NULL)
9c7c9ec1 123 ,fTrackletTree(NULL)
3d0c7d6d 124 ,fDigitsManager(new AliTRDdigitsManager())
9c7c9ec1 125 ,fTrackletContainer(NULL)
3fe61b77 126 ,fRawVersion(2)
3fe61b77 127 ,fTransform(new AliTRDtransform(0))
1b3a719f 128 ,fDigits(NULL)
064d808d 129 ,fIndexes(NULL)
064d808d 130 ,fMaxThresh(0)
131 ,fSigThresh(0)
132 ,fMinMaxCutSigma(0)
133 ,fMinLeftRightCutSigma(0)
134 ,fLayer(0)
135 ,fDet(0)
136 ,fVolid(0)
137 ,fColMax(0)
138 ,fTimeTotal(0)
139 ,fCalGainFactorROC(NULL)
140 ,fCalGainFactorDetValue(0)
141 ,fCalNoiseROC(NULL)
142 ,fCalNoiseDetValue(0)
1b3a719f 143 ,fCalPadStatusROC(NULL)
064d808d 144 ,fClusterROC(0)
145 ,firstClusterROC(0)
3d0c7d6d 146 ,fNoOfClusters(0)
615f0826 147 ,fBaseline(0)
8dbc94c7 148 ,fRawStream(NULL)
f7336fa3 149{
150 //
6d50f529 151 // AliTRDclusterizer constructor
f7336fa3 152 //
153
b72f4eaf 154 SetBit(kLabels, kTRUE);
8cbe253a 155 SetBit(knewDM, kFALSE);
3d0c7d6d 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
b72f4eaf 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)
a5b99acd 177 ,fTracklets(NULL)
9c7c9ec1 178 ,fTrackletTree(NULL)
3fe61b77 179 ,fDigitsManager(NULL)
66f6bfd9 180 ,fTrackletContainer(NULL)
3fe61b77 181 ,fRawVersion(2)
3fe61b77 182 ,fTransform(NULL)
1b3a719f 183 ,fDigits(NULL)
064d808d 184 ,fIndexes(NULL)
064d808d 185 ,fMaxThresh(0)
186 ,fSigThresh(0)
187 ,fMinMaxCutSigma(0)
188 ,fMinLeftRightCutSigma(0)
189 ,fLayer(0)
190 ,fDet(0)
191 ,fVolid(0)
192 ,fColMax(0)
193 ,fTimeTotal(0)
194 ,fCalGainFactorROC(NULL)
195 ,fCalGainFactorDetValue(0)
196 ,fCalNoiseROC(NULL)
197 ,fCalNoiseDetValue(0)
1b3a719f 198 ,fCalPadStatusROC(NULL)
064d808d 199 ,fClusterROC(0)
200 ,firstClusterROC(0)
3d0c7d6d 201 ,fNoOfClusters(0)
615f0826 202 ,fBaseline(0)
8dbc94c7 203 ,fRawStream(NULL)
8230f242 204{
205 //
206 // AliTRDclusterizer copy constructor
207 //
208
b72f4eaf 209 SetBit(kLabels, kTRUE);
8cbe253a 210 SetBit(knewDM, kFALSE);
3d0c7d6d 211
b72f4eaf 212 //FillLUT();
fc546d21 213
8230f242 214}
215
f7336fa3 216//_____________________________________________________________________________
217AliTRDclusterizer::~AliTRDclusterizer()
218{
8230f242 219 //
220 // AliTRDclusterizer destructor
221 //
f7336fa3 222
48f8adf3 223 if (fRecPoints/* && IsClustersOwner()*/){
66f6bfd9 224 fRecPoints->Delete();
225 delete fRecPoints;
226 }
3fe61b77 227
a5b99acd 228 if (fTracklets){
229 fTracklets->Delete();
230 delete fTracklets;
231 }
232
66f6bfd9 233 if (fDigitsManager) {
234 delete fDigitsManager;
235 fDigitsManager = NULL;
236 }
9c7c9ec1 237
66f6bfd9 238 if (fTrackletContainer){
8cbe253a 239 delete [] fTrackletContainer[0];
240 delete [] fTrackletContainer[1];
241 delete [] fTrackletContainer;
66f6bfd9 242 fTrackletContainer = NULL;
243 }
3fe61b77 244
66f6bfd9 245 if (fTransform){
246 delete fTransform;
8dbc94c7 247 fTransform = NULL;
248 }
249
250 if (fRawStream){
251 delete fRawStream;
252 fRawStream = NULL;
66f6bfd9 253 }
fc546d21 254
f7336fa3 255}
256
8230f242 257//_____________________________________________________________________________
dd9a6ee3 258AliTRDclusterizer &AliTRDclusterizer::operator=(const AliTRDclusterizer &c)
259{
260 //
261 // Assignment operator
262 //
263
3fe61b77 264 if (this != &c)
265 {
266 ((AliTRDclusterizer &) c).Copy(*this);
267 }
053767a4 268
dd9a6ee3 269 return *this;
270
271}
272
273//_____________________________________________________________________________
0ab5ff49 274void AliTRDclusterizer::Copy(TObject &c) const
8230f242 275{
276 //
277 // Copy function
278 //
279
3fe61b77 280 ((AliTRDclusterizer &) c).fClusterTree = NULL;
281 ((AliTRDclusterizer &) c).fRecPoints = NULL;
9c7c9ec1 282 ((AliTRDclusterizer &) c).fTrackletTree = NULL;
3fe61b77 283 ((AliTRDclusterizer &) c).fDigitsManager = NULL;
9c7c9ec1 284 ((AliTRDclusterizer &) c).fTrackletContainer = NULL;
3fe61b77 285 ((AliTRDclusterizer &) c).fRawVersion = fRawVersion;
3fe61b77 286 ((AliTRDclusterizer &) c).fTransform = NULL;
1b3a719f 287 ((AliTRDclusterizer &) c).fDigits = NULL;
064d808d 288 ((AliTRDclusterizer &) c).fIndexes = NULL;
064d808d 289 ((AliTRDclusterizer &) c).fMaxThresh = 0;
290 ((AliTRDclusterizer &) c).fSigThresh = 0;
291 ((AliTRDclusterizer &) c).fMinMaxCutSigma= 0;
292 ((AliTRDclusterizer &) c).fMinLeftRightCutSigma = 0;
293 ((AliTRDclusterizer &) c).fLayer = 0;
294 ((AliTRDclusterizer &) c).fDet = 0;
295 ((AliTRDclusterizer &) c).fVolid = 0;
296 ((AliTRDclusterizer &) c).fColMax = 0;
297 ((AliTRDclusterizer &) c).fTimeTotal = 0;
298 ((AliTRDclusterizer &) c).fCalGainFactorROC = NULL;
299 ((AliTRDclusterizer &) c).fCalGainFactorDetValue = 0;
300 ((AliTRDclusterizer &) c).fCalNoiseROC = NULL;
301 ((AliTRDclusterizer &) c).fCalNoiseDetValue = 0;
1b3a719f 302 ((AliTRDclusterizer &) c).fCalPadStatusROC = NULL;
064d808d 303 ((AliTRDclusterizer &) c).fClusterROC = 0;
304 ((AliTRDclusterizer &) c).firstClusterROC= 0;
3d0c7d6d 305 ((AliTRDclusterizer &) c).fNoOfClusters = 0;
615f0826 306 ((AliTRDclusterizer &) c).fBaseline = 0;
8dbc94c7 307 ((AliTRDclusterizer &) c).fRawStream = NULL;
615f0826 308
8230f242 309}
310
f7336fa3 311//_____________________________________________________________________________
3e1a3ad8 312Bool_t AliTRDclusterizer::Open(const Char_t *name, Int_t nEvent)
f7336fa3 313{
314 //
3e1a3ad8 315 // Opens the AliROOT file. Output and input are in the same file
f7336fa3 316 //
6d50f529 317
e191bb57 318 TString evfoldname = AliConfig::GetDefaultEventFolderName();
6d50f529 319 fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
320
321 if (!fRunLoader) {
19dd5b2f 322 fRunLoader = AliRunLoader::Open(name);
6d50f529 323 }
324
325 if (!fRunLoader) {
326 AliError(Form("Can not open session for file %s.",name));
327 return kFALSE;
328 }
88cb7938 329
330 OpenInput(nEvent);
331 OpenOutput();
6d50f529 332
3e1a3ad8 333 return kTRUE;
f7336fa3 334
6d50f529 335}
3e1a3ad8 336
337//_____________________________________________________________________________
88cb7938 338Bool_t AliTRDclusterizer::OpenOutput()
3e1a3ad8 339{
340 //
341 // Open the output file
342 //
343
66f6bfd9 344 if (!fReconstructor->IsWritingClusters()) return kTRUE;
345
346 TObjArray *ioArray = 0x0;
88cb7938 347
348 AliLoader* loader = fRunLoader->GetLoader("TRDLoader");
349 loader->MakeTree("R");
6d50f529 350
88cb7938 351 fClusterTree = loader->TreeR();
66f6bfd9 352 fClusterTree->Branch("TRDcluster", "TObjArray", &ioArray, 32000, 0);
3e1a3ad8 353
3e1a3ad8 354 return kTRUE;
355
356}
357
25ca55ce 358//_____________________________________________________________________________
e7295a3a 359Bool_t AliTRDclusterizer::OpenOutput(TTree *const clusterTree)
25ca55ce 360{
361 //
362 // Connect the output tree
363 //
364
66f6bfd9 365 // clusters writing
366 if (fReconstructor->IsWritingClusters()){
367 TObjArray *ioArray = 0x0;
368 fClusterTree = clusterTree;
369 fClusterTree->Branch("TRDcluster", "TObjArray", &ioArray, 32000, 0);
370 }
11d0be11 371 return kTRUE;
372}
373
374//_____________________________________________________________________________
375Bool_t AliTRDclusterizer::OpenTrackletOutput()
376{
377 //
378 // Tracklet writing
379 //
9c7c9ec1 380
3a039a31 381 if (fReconstructor->IsWritingTracklets()){
9c7c9ec1 382 TString evfoldname = AliConfig::GetDefaultEventFolderName();
383 fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
384
385 if (!fRunLoader) {
386 fRunLoader = AliRunLoader::Open("galice.root");
387 }
388 if (!fRunLoader) {
389 AliError(Form("Can not open session for file galice.root."));
390 return kFALSE;
391 }
392
393 UInt_t **leaves = new UInt_t *[2];
394 AliDataLoader *dl = fRunLoader->GetLoader("TRDLoader")->GetDataLoader("tracklets");
395 if (!dl) {
396 AliError("Could not get the tracklets data loader!");
90cf7133 397 dl = new AliDataLoader("TRD.Tracklets.root","tracklets", "tracklets");
9c7c9ec1 398 fRunLoader->GetLoader("TRDLoader")->AddDataLoader(dl);
399 }
11d0be11 400 fTrackletTree = dl->Tree();
401 if (!fTrackletTree)
402 {
403 dl->MakeTree();
404 fTrackletTree = dl->Tree();
405 }
406 TBranch *trkbranch = fTrackletTree->GetBranch("trkbranch");
407 if (!trkbranch)
408 fTrackletTree->Branch("trkbranch",leaves[0],"det/i:side/i:tracklets[256]/i");
9c7c9ec1 409 }
410
25ca55ce 411 return kTRUE;
25ca55ce 412}
413
3e1a3ad8 414//_____________________________________________________________________________
88cb7938 415Bool_t AliTRDclusterizer::OpenInput(Int_t nEvent)
f7336fa3 416{
417 //
418 // Opens a ROOT-file with TRD-hits and reads in the digits-tree
419 //
420
f7336fa3 421 // Import the Trees for the event nEvent in the file
bdbb05bb 422 fRunLoader->GetEvent(nEvent);
88cb7938 423
f7336fa3 424 return kTRUE;
425
426}
427
428//_____________________________________________________________________________
793ff80c 429Bool_t AliTRDclusterizer::WriteClusters(Int_t det)
f7336fa3 430{
431 //
3e1a3ad8 432 // Fills TRDcluster branch in the tree with the clusters
793ff80c 433 // found in detector = det. For det=-1 writes the tree.
a3c76cdc 434 //
793ff80c 435
6d50f529 436 if ((det < -1) ||
437 (det >= AliTRDgeometry::Ndet())) {
438 AliError(Form("Unexpected detector index %d.\n",det));
3e1a3ad8 439 return kFALSE;
793ff80c 440 }
317b5644 441
66f6bfd9 442 TObjArray *ioArray = new TObjArray(400);
3e1a3ad8 443 TBranch *branch = fClusterTree->GetBranch("TRDcluster");
444 if (!branch) {
3e1a3ad8 445 branch = fClusterTree->Branch("TRDcluster","TObjArray",&ioArray,32000,0);
66f6bfd9 446 } else branch->SetAddress(&ioArray);
66f6bfd9 447
448 Int_t nRecPoints = RecPoints()->GetEntriesFast();
449 if(det >= 0){
3e1a3ad8 450 for (Int_t i = 0; i < nRecPoints; i++) {
bdbb05bb 451 AliTRDcluster *c = (AliTRDcluster *) RecPoints()->UncheckedAt(i);
66f6bfd9 452 if(det != c->GetDetector()) continue;
453 ioArray->AddLast(c);
793ff80c 454 }
3e1a3ad8 455 fClusterTree->Fill();
839f6b9b 456 ioArray->Clear();
66f6bfd9 457 } else {
839f6b9b 458 Int_t detOld = -1, nw(0);
66f6bfd9 459 for (Int_t i = 0; i < nRecPoints; i++) {
460 AliTRDcluster *c = (AliTRDcluster *) RecPoints()->UncheckedAt(i);
461 if(c->GetDetector() != detOld){
839f6b9b 462 nw += ioArray->GetEntriesFast();
66f6bfd9 463 fClusterTree->Fill();
464 ioArray->Clear();
465 detOld = c->GetDetector();
466 }
467 ioArray->AddLast(c);
a6dd11e9 468 }
839f6b9b 469 if(ioArray->GetEntriesFast()){
470 nw += ioArray->GetEntriesFast();
471 fClusterTree->Fill();
472 ioArray->Clear();
473 }
474 AliDebug(2, Form("Clusters FOUND[%d] WRITTEN[%d] STATUS[%s]", nRecPoints, nw, nw==nRecPoints?"OK":"FAILED"));
793ff80c 475 }
66f6bfd9 476 delete ioArray;
6d50f529 477
66f6bfd9 478 return kTRUE;
f7336fa3 479}
793ff80c 480
9c7c9ec1 481//_____________________________________________________________________________
482Bool_t AliTRDclusterizer::WriteTracklets(Int_t det)
483{
eb52b657 484 //
485 // Write the raw data tracklets into seperate file
486 //
9c7c9ec1 487
488 UInt_t **leaves = new UInt_t *[2];
489 for (Int_t i=0; i<2 ;i++){
317b5644 490 leaves[i] = new UInt_t[258];
491 leaves[i][0] = det; // det
492 leaves[i][1] = i; // side
493 memcpy(leaves[i]+2, fTrackletContainer[i], sizeof(UInt_t) * 256);
9c7c9ec1 494 }
495
9c7c9ec1 496 if (!fTrackletTree){
497 AliDataLoader *dl = fRunLoader->GetLoader("TRDLoader")->GetDataLoader("tracklets");
498 dl->MakeTree();
499 fTrackletTree = dl->Tree();
500 }
501
502 TBranch *trkbranch = fTrackletTree->GetBranch("trkbranch");
503 if (!trkbranch) {
504 trkbranch = fTrackletTree->Branch("trkbranch",leaves[0],"det/i:side/i:tracklets[256]/i");
505 }
506
507 for (Int_t i=0; i<2; i++){
317b5644 508 if (leaves[i][2]>0) {
509 trkbranch->SetAddress(leaves[i]);
510 fTrackletTree->Fill();
511 }
9c7c9ec1 512 }
513
514 AliDataLoader *dl = fRunLoader->GetLoader("TRDLoader")->GetDataLoader("tracklets");
515 dl->WriteData("OVERWRITE");
516 //dl->Unload();
517 delete [] leaves;
518
eb52b657 519 return kTRUE;
520
9c7c9ec1 521}
522
3fe61b77 523//_____________________________________________________________________________
524Bool_t AliTRDclusterizer::ReadDigits()
525{
526 //
527 // Reads the digits arrays from the input aliroot file
528 //
529
530 if (!fRunLoader) {
531 AliError("No run loader available");
532 return kFALSE;
533 }
534
535 AliLoader* loader = fRunLoader->GetLoader("TRDLoader");
536 if (!loader->TreeD()) {
537 loader->LoadDigits();
538 }
539
540 // Read in the digit arrays
541 return (fDigitsManager->ReadDigits(loader->TreeD()));
542
543}
544
545//_____________________________________________________________________________
546Bool_t AliTRDclusterizer::ReadDigits(TTree *digitsTree)
547{
548 //
549 // Reads the digits arrays from the input tree
550 //
551
552 // Read in the digit arrays
553 return (fDigitsManager->ReadDigits(digitsTree));
554
555}
556
557//_____________________________________________________________________________
558Bool_t AliTRDclusterizer::ReadDigits(AliRawReader *rawReader)
559{
560 //
561 // Reads the digits arrays from the ddl file
562 //
563
564 AliTRDrawData raw;
565 fDigitsManager = raw.Raw2Digits(rawReader);
566
567 return kTRUE;
568
569}
570
571//_____________________________________________________________________________
572Bool_t AliTRDclusterizer::MakeClusters()
573{
574 //
575 // Creates clusters from digits
576 //
577
578 // Propagate info from the digits manager
b72f4eaf 579 if (TestBit(kLabels)){
580 SetBit(kLabels, fDigitsManager->UsesDictionaries());
66f6bfd9 581 }
582
3fe61b77 583 Bool_t fReturn = kTRUE;
66f6bfd9 584 for (Int_t i = 0; i < AliTRDgeometry::kNdet; i++){
585
b65e5048 586 AliTRDarrayADC *digitsIn = (AliTRDarrayADC*) fDigitsManager->GetDigits(i); //mod
66f6bfd9 587 // This is to take care of switched off super modules
588 if (!digitsIn->HasData()) continue;
589 digitsIn->Expand();
0d64b05f 590 digitsIn->DeleteNegatives(); // Restore digits array to >=0 values
66f6bfd9 591 AliTRDSignalIndex* indexes = fDigitsManager->GetIndexes(i);
592 if (indexes->IsAllocated() == kFALSE){
593 fDigitsManager->BuildIndexes(i);
594 }
595
596 Bool_t fR = kFALSE;
597 if (indexes->HasEntry()){
b72f4eaf 598 if (TestBit(kLabels)){
66f6bfd9 599 for (Int_t iDict = 0; iDict < AliTRDdigitsManager::kNDict; iDict++){
b65e5048 600 AliTRDarrayDictionary *tracksIn = 0; //mod
601 tracksIn = (AliTRDarrayDictionary *) fDigitsManager->GetDictionary(i,iDict); //mod
66f6bfd9 602 tracksIn->Expand();
3fe61b77 603 }
66f6bfd9 604 }
605 fR = MakeClusters(i);
606 fReturn = fR && fReturn;
3fe61b77 607 }
66f6bfd9 608
609 //if (fR == kFALSE){
610 // if(IsWritingClusters()) WriteClusters(i);
611 // ResetRecPoints();
612 //}
613
614 // No compress just remove
615 fDigitsManager->RemoveDigits(i);
616 fDigitsManager->RemoveDictionaries(i);
617 fDigitsManager->ClearIndexes(i);
618 }
619
620 if(fReconstructor->IsWritingClusters()) WriteClusters(-1);
3fe61b77 621
66f6bfd9 622 AliInfo(Form("Number of found clusters : %d", RecPoints()->GetEntriesFast()));
3fe61b77 623
66f6bfd9 624 return fReturn;
eb52b657 625
3fe61b77 626}
627
628//_____________________________________________________________________________
629Bool_t AliTRDclusterizer::Raw2Clusters(AliRawReader *rawReader)
630{
631 //
632 // Creates clusters from raw data
633 //
634
fc546d21 635 return Raw2ClustersChamber(rawReader);
3fe61b77 636
3fe61b77 637}
638
639//_____________________________________________________________________________
640Bool_t AliTRDclusterizer::Raw2ClustersChamber(AliRawReader *rawReader)
641{
642 //
643 // Creates clusters from raw data
644 //
645
646 // Create the digits manager
66f6bfd9 647 if (!fDigitsManager){
8cbe253a 648 SetBit(knewDM, kTRUE);
3d0c7d6d 649 fDigitsManager = new AliTRDdigitsManager(kTRUE);
66f6bfd9 650 fDigitsManager->CreateArrays();
651 }
3fe61b77 652
b72f4eaf 653 fDigitsManager->SetUseDictionaries(TestBit(kLabels));
3fe61b77 654
317b5644 655 // tracklet container for raw tracklet writing
a5b99acd 656 if (!fTrackletContainer && ( fReconstructor->IsWritingTracklets() || fReconstructor->IsProcessingTracklets() )) {
317b5644 657 // maximum tracklets for one HC
658 const Int_t kTrackletChmb=256;
659 fTrackletContainer = new UInt_t *[2];
660 fTrackletContainer[0] = new UInt_t[kTrackletChmb];
661 fTrackletContainer[1] = new UInt_t[kTrackletChmb];
662 }
663
8dbc94c7 664 if(!fRawStream)
665 fRawStream = AliTRDrawStreamBase::GetRawStream(rawReader);
666 else
667 fRawStream->SetReader(rawReader);
668
17e0e535 669 if(fReconstructor->IsHLT()){
8dbc94c7 670 fRawStream->SetSharedPadReadout(kFALSE);
17e0e535 671 fRawStream->SetNoErrorWarning();
672 }
3fe61b77 673
17e0e535 674 AliDebug(1,Form("Stream version: %s", fRawStream->IsA()->GetName()));
3fe61b77 675
676 Int_t det = 0;
8dbc94c7 677 while ((det = fRawStream->NextChamber(fDigitsManager,fTrackletContainer)) >= 0){
66f6bfd9 678 Bool_t iclusterBranch = kFALSE;
679 if (fDigitsManager->GetIndexes(det)->HasEntry()){
486c2339 680 iclusterBranch = MakeClusters(det);
3fe61b77 681 }
66f6bfd9 682
66f68968 683 fDigitsManager->ClearArrays(det);
684
317b5644 685 if (!fReconstructor->IsWritingTracklets()) continue;
686 if (*(fTrackletContainer[0]) > 0 || *(fTrackletContainer[1]) > 0) WriteTracklets(det);
9c7c9ec1 687 }
486c2339 688
a5b99acd 689 if (fTrackletContainer){
317b5644 690 delete [] fTrackletContainer[0];
691 delete [] fTrackletContainer[1];
692 delete [] fTrackletContainer;
693 fTrackletContainer = NULL;
694 }
695
66f6bfd9 696 if(fReconstructor->IsWritingClusters()) WriteClusters(-1);
3fe61b77 697
8cbe253a 698 if(!TestBit(knewDM)){
699 delete fDigitsManager;
700 fDigitsManager = NULL;
701 }
dfbb4bb9 702
3d0c7d6d 703 AliInfo(Form("Number of found clusters : %d", fNoOfClusters));
66f6bfd9 704 return kTRUE;
eb52b657 705
3fe61b77 706}
707
fa7427d0 708//_____________________________________________________________________________
709UChar_t AliTRDclusterizer::GetStatus(Short_t &signal)
710{
023b669c 711 //
712 // Check if a pad is masked
713 //
714
317b5644 715 UChar_t status = 0;
716
717 if(signal>0 && TESTBIT(signal, 10)){
718 CLRBIT(signal, 10);
719 for(int ibit=0; ibit<4; ibit++){
720 if(TESTBIT(signal, 11+ibit)){
721 SETBIT(status, ibit);
722 CLRBIT(signal, 11+ibit);
723 }
724 }
725 }
726 return status;
fa7427d0 727}
728
6d9e79cc 729//_____________________________________________________________________________
e7295a3a 730void AliTRDclusterizer::SetPadStatus(const UChar_t status, UChar_t &out) const {
6d9e79cc 731 //
732 // Set the pad status into out
733 // First three bits are needed for the position encoding
734 //
064d808d 735 out |= status << 3;
6d9e79cc 736}
737
738//_____________________________________________________________________________
064d808d 739UChar_t AliTRDclusterizer::GetPadStatus(UChar_t encoding) const {
6d9e79cc 740 //
741 // return the staus encoding of the corrupted pad
742 //
743 return static_cast<UChar_t>(encoding >> 3);
744}
745
746//_____________________________________________________________________________
064d808d 747Int_t AliTRDclusterizer::GetCorruption(UChar_t encoding) const {
6d9e79cc 748 //
749 // Return the position of the corruption
750 //
751 return encoding & 7;
752}
753
3fe61b77 754//_____________________________________________________________________________
755Bool_t AliTRDclusterizer::MakeClusters(Int_t det)
756{
757 //
758 // Generates the cluster.
759 //
760
761 // Get the digits
615f0826 762 fDigits = (AliTRDarrayADC *) fDigitsManager->GetDigits(det); //mod
763 fBaseline = fDigitsManager->GetDigitsParam()->GetADCbaseline();
f5375dcb 764
3fe61b77 765 // This is to take care of switched off super modules
b72f4eaf 766 if (!fDigits->HasData()) return kFALSE;
3fe61b77 767
064d808d 768 fIndexes = fDigitsManager->GetIndexes(det);
b72f4eaf 769 if (fIndexes->IsAllocated() == kFALSE) {
770 AliError("Indexes do not exist!");
771 return kFALSE;
772 }
064d808d 773
3fe61b77 774 AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
b72f4eaf 775 if (!calibration) {
776 AliFatal("No AliTRDcalibDB instance available\n");
777 return kFALSE;
778 }
3fe61b77 779
3a039a31 780 if (!fReconstructor){
781 AliError("Reconstructor not set\n");
782 return kFALSE;
783 }
c5f589b9 784
839f6b9b 785
064d808d 786 fMaxThresh = fReconstructor->GetRecoParam()->GetClusMaxThresh();
787 fSigThresh = fReconstructor->GetRecoParam()->GetClusSigThresh();
788 fMinMaxCutSigma = fReconstructor->GetRecoParam()->GetMinMaxCutSigma();
789 fMinLeftRightCutSigma = fReconstructor->GetRecoParam()->GetMinLeftRightCutSigma();
3fe61b77 790
064d808d 791 Int_t istack = fIndexes->GetStack();
792 fLayer = fIndexes->GetLayer();
793 Int_t isector = fIndexes->GetSM();
3fe61b77 794
795 // Start clustering in the chamber
796
064d808d 797 fDet = AliTRDgeometry::GetDetector(fLayer,istack,isector);
798 if (fDet != det) {
b65e5048 799 AliError("Strange Detector number Missmatch!");
800 return kFALSE;
801 }
3fe61b77 802
839f6b9b 803 AliDebug(2, Form("Det[%d] @ Sec[%d] Stk[%d] Ly[%d]", fDet, isector, istack, fLayer));
804
3fe61b77 805 // TRD space point transformation
806 fTransform->SetDetector(det);
807
064d808d 808 Int_t iGeoLayer = AliGeomManager::kTRD1 + fLayer;
053767a4 809 Int_t iGeoModule = istack + AliTRDgeometry::Nstack() * isector;
064d808d 810 fVolid = AliGeomManager::LayerToVolUID(iGeoLayer,iGeoModule);
3fe61b77 811
a5b99acd 812 if(fReconstructor->IsProcessingTracklets() && fTrackletContainer)
813 AddTrackletsToArray();
814
1b3a719f 815 fColMax = fDigits->GetNcol();
816 //Int_t nRowMax = fDigits->GetNrow();
817 fTimeTotal = fDigits->GetNtime();
3fe61b77 818
819 // Detector wise calibration object for the gain factors
064d808d 820 const AliTRDCalDet *calGainFactorDet = calibration->GetGainFactorDet();
3fe61b77 821 // Calibration object with pad wise values for the gain factors
064d808d 822 fCalGainFactorROC = calibration->GetGainFactorROC(fDet);
3fe61b77 823 // Calibration value for chamber wise gain factor
064d808d 824 fCalGainFactorDetValue = calGainFactorDet->GetValue(fDet);
0e09df31 825
df83a620 826 // Detector wise calibration object for the noise
064d808d 827 const AliTRDCalDet *calNoiseDet = calibration->GetNoiseDet();
df83a620 828 // Calibration object with pad wise values for the noise
064d808d 829 fCalNoiseROC = calibration->GetNoiseROC(fDet);
df83a620 830 // Calibration value for chamber wise noise
064d808d 831 fCalNoiseDetValue = calNoiseDet->GetValue(fDet);
1b3a719f 832
833 // Calibration object with the pad status
834 fCalPadStatusROC = calibration->GetPadStatusROC(fDet);
6b0d4ad5 835
a2fbb6ec 836 SetBit(kLUT, fReconstructor->GetRecoParam()->UseLUT());
837 SetBit(kGAUS, fReconstructor->GetRecoParam()->UseGAUS());
b72f4eaf 838 SetBit(kHLT, fReconstructor->IsHLT());
6b0d4ad5 839
064d808d 840 firstClusterROC = -1;
841 fClusterROC = 0;
3fe61b77 842
b72f4eaf 843 // Apply the gain and the tail cancelation via digital filter
a2fbb6ec 844 if(fReconstructor->GetRecoParam()->UseTailCancelation()) TailCancelation();
023b669c 845
486c2339 846 MaxStruct curr, last;
064d808d 847 Int_t nMaximas = 0, nCorrupted = 0;
f5375dcb 848
064d808d 849 // Here the clusterfining is happening
850
615f0826 851 for(curr.time = 0; curr.time < fTimeTotal; curr.time++){
852 while(fIndexes->NextRCIndex(curr.row, curr.col)){
853 if(IsMaximum(curr, curr.padStatus, &curr.signals[0])){
854 if(last.row>-1){
855 if(curr.time==last.time && curr.row==last.row && curr.col==last.col+2) FivePadCluster(last, curr);
b72f4eaf 856 CreateCluster(last);
857 }
615f0826 858 last=curr; curr.fivePad=kFALSE;
b72f4eaf 859 }
b72f4eaf 860 }
861 }
615f0826 862 if(last.row>-1) CreateCluster(last);
df83a620 863
a2fbb6ec 864 if(fReconstructor->GetRecoParam()->GetStreamLevel(AliTRDrecoParam::kClusterizer) > 2 && fReconstructor->IsDebugStreaming()){
865 TTreeSRedirector* fDebugStream = fReconstructor->GetDebugStream(AliTRDrecoParam::kClusterizer);
064d808d 866 (*fDebugStream) << "MakeClusters"
b72f4eaf 867 << "Detector=" << det
868 << "NMaxima=" << nMaximas
869 << "NClusters=" << fClusterROC
870 << "NCorrupted=" << nCorrupted
871 << "\n";
df83a620 872 }
b72f4eaf 873 if (TestBit(kLabels)) AddLabels();
f5375dcb 874
064d808d 875 return kTRUE;
f5375dcb 876
064d808d 877}
3fe61b77 878
064d808d 879//_____________________________________________________________________________
1b3a719f 880Bool_t AliTRDclusterizer::IsMaximum(const MaxStruct &Max, UChar_t &padStatus, Short_t *const Signals)
064d808d 881{
882 //
883 // Returns true if this row,col,time combination is a maximum.
884 // Gives back the padStatus and the signals of the center pad and the two neighbouring pads.
885 //
c96d03ba 886
615f0826 887 Float_t gain = fCalGainFactorDetValue * fCalGainFactorROC->GetValue(Max.col,Max.row);
888 Signals[1] = (Short_t)((fDigits->GetData(Max.row, Max.col, Max.time) - fBaseline) / gain + 0.5f);
c96d03ba 889 if(Signals[1] < fMaxThresh) return kFALSE;
3fe61b77 890
615f0826 891 Float_t noiseMiddleThresh = fMinMaxCutSigma*fCalNoiseDetValue*fCalNoiseROC->GetValue(Max.col, Max.row);
c96d03ba 892 if (Signals[1] < noiseMiddleThresh) return kFALSE;
f5375dcb 893
615f0826 894 if (Max.col + 1 >= fColMax || Max.col < 1) return kFALSE;
3fe61b77 895
b72f4eaf 896 UChar_t status[3]={
615f0826 897 fCalPadStatusROC->GetStatus(Max.col-1, Max.row)
898 ,fCalPadStatusROC->GetStatus(Max.col, Max.row)
899 ,fCalPadStatusROC->GetStatus(Max.col+1, Max.row)
b72f4eaf 900 };
486c2339 901
615f0826 902 gain = fCalGainFactorDetValue * fCalGainFactorROC->GetValue(Max.col-1,Max.row);
19dd8eb9 903 Signals[0] = (Short_t)((fDigits->GetData(Max.row, Max.col-1, Max.time) - fBaseline) / gain + 0.5f);
904 gain = fCalGainFactorDetValue * fCalGainFactorROC->GetValue(Max.col+1,Max.row);
905 Signals[2] = (Short_t)((fDigits->GetData(Max.row, Max.col+1, Max.time) - fBaseline) / gain + 0.5f);
c96d03ba 906
486c2339 907 if(!(status[0] | status[1] | status[2])) {//all pads are good
c96d03ba 908 if ((Signals[2] <= Signals[1]) && (Signals[0] < Signals[1])) {
909 if ((Signals[2] >= fSigThresh) || (Signals[0] >= fSigThresh)) {
5a8db426 910 if(Signals[0]<0)Signals[0]=0;
911 if(Signals[2]<0)Signals[2]=0;
b72f4eaf 912 Float_t noiseSumThresh = fMinLeftRightCutSigma
913 * fCalNoiseDetValue
615f0826 914 * fCalNoiseROC->GetValue(Max.col, Max.row);
c96d03ba 915 if ((Signals[2]+Signals[0]+Signals[1]) < noiseSumThresh) return kFALSE;
b72f4eaf 916 padStatus = 0;
917 return kTRUE;
eb52b657 918 }
064d808d 919 }
b72f4eaf 920 } else { // at least one of the pads is bad, and reject candidates with more than 1 problematic pad
5a8db426 921 if(Signals[0]<0)Signals[0]=0;
922 if(Signals[2]<0)Signals[2]=0;
c96d03ba 923 if (status[2] && (!(status[0] || status[1])) && Signals[1] > Signals[0] && Signals[0] >= fSigThresh) {
924 Signals[2]=0;
486c2339 925 SetPadStatus(status[2], padStatus);
926 return kTRUE;
927 }
c96d03ba 928 else if (status[0] && (!(status[1] || status[2])) && Signals[1] >= Signals[2] && Signals[2] >= fSigThresh) {
929 Signals[0]=0;
486c2339 930 SetPadStatus(status[0], padStatus);
931 return kTRUE;
932 }
c96d03ba 933 else if (status[1] && (!(status[0] || status[2])) && ((Signals[2] >= fSigThresh) || (Signals[0] >= fSigThresh))) {
615f0826 934 Signals[1] = (Short_t)(fMaxThresh + 0.5f);
486c2339 935 SetPadStatus(status[1], padStatus);
064d808d 936 return kTRUE;
937 }
938 }
939 return kFALSE;
940}
3fe61b77 941
064d808d 942//_____________________________________________________________________________
1b3a719f 943Bool_t AliTRDclusterizer::FivePadCluster(MaxStruct &ThisMax, MaxStruct &NeighbourMax)
064d808d 944{
945 //
946 // Look for 5 pad cluster with minimum in the middle
947 // Gives back the ratio
948 //
615f0826 949
950 if (ThisMax.col >= fColMax - 3) return kFALSE;
951 Float_t gain;
952 if (ThisMax.col < fColMax - 5){
953 gain = fCalGainFactorDetValue * fCalGainFactorROC->GetValue(ThisMax.col+4,ThisMax.row);
954 if (fDigits->GetData(ThisMax.row, ThisMax.col+4, ThisMax.time) - fBaseline >= fSigThresh * gain)
486c2339 955 return kFALSE;
064d808d 956 }
615f0826 957 if (ThisMax.col > 1) {
958 gain = fCalGainFactorDetValue * fCalGainFactorROC->GetValue(ThisMax.col-2,ThisMax.row);
959 if (fDigits->GetData(ThisMax.row, ThisMax.col-2, ThisMax.time) - fBaseline >= fSigThresh * gain)
486c2339 960 return kFALSE;
961 }
962
486c2339 963 const Float_t kEpsilon = 0.01;
615f0826 964 Double_t padSignal[5] = {ThisMax.signals[0], ThisMax.signals[1], ThisMax.signals[2],
965 NeighbourMax.signals[1], NeighbourMax.signals[2]};
486c2339 966
967 // Unfold the two maxima and set the signal on
968 // the overlapping pad to the ratio
1b3a719f 969 Float_t ratio = Unfold(kEpsilon,fLayer,padSignal);
615f0826 970 ThisMax.signals[2] = (Short_t)(ThisMax.signals[2]*ratio + 0.5f);
971 NeighbourMax.signals[0] = (Short_t)(NeighbourMax.signals[0]*(1-ratio) + 0.5f);
972 ThisMax.fivePad=kTRUE;
973 NeighbourMax.fivePad=kTRUE;
486c2339 974 return kTRUE;
c96d03ba 975
064d808d 976}
eb52b657 977
064d808d 978//_____________________________________________________________________________
486c2339 979void AliTRDclusterizer::CreateCluster(const MaxStruct &Max)
064d808d 980{
981 //
3d0c7d6d 982 // Creates a cluster at the given position and saves it in fRecPoints
064d808d 983 //
eb52b657 984
c96d03ba 985 Int_t nPadCount = 1;
615f0826 986 Short_t signals[7] = { 0, 0, Max.signals[0], Max.signals[1], Max.signals[2], 0, 0 };
c96d03ba 987 if(!TestBit(kHLT)) CalcAdditionalInfo(Max, signals, nPadCount);
988
615f0826 989 AliTRDcluster cluster(fDet, ((UChar_t) Max.col), ((UChar_t) Max.row), ((UChar_t) Max.time), signals, fVolid);
c96d03ba 990 cluster.SetNPads(nPadCount);
b72f4eaf 991 if(TestBit(kLUT)) cluster.SetRPhiMethod(AliTRDcluster::kLUT);
992 else if(TestBit(kGAUS)) cluster.SetRPhiMethod(AliTRDcluster::kGAUS);
993 else cluster.SetRPhiMethod(AliTRDcluster::kCOG);
3fe61b77 994
615f0826 995 cluster.SetFivePad(Max.fivePad);
b72f4eaf 996 // set pads status for the cluster
997 UChar_t maskPosition = GetCorruption(Max.padStatus);
998 if (maskPosition) {
999 cluster.SetPadMaskedPosition(maskPosition);
1000 cluster.SetPadMaskedStatus(GetPadStatus(Max.padStatus));
1001 }
064d808d 1002
1003 // Transform the local cluster coordinates into calibrated
1004 // space point positions defined in the local tracking system.
1005 // Here the calibration for T0, Vdrift and ExB is applied as well.
b72f4eaf 1006 if(!fTransform->Transform(&cluster)) return;
1007 // Temporarily store the Max.Row, column and time bin of the center pad
1008 // Used to later on assign the track indices
615f0826 1009 cluster.SetLabel(Max.row, 0);
1010 cluster.SetLabel(Max.col, 1);
1011 cluster.SetLabel(Max.time,2);
c96d03ba 1012
b72f4eaf 1013 //needed for HLT reconstruction
1014 AddClusterToArray(&cluster);
1015
1016 // Store the index of the first cluster in the current ROC
1017 if (firstClusterROC < 0) firstClusterROC = fNoOfClusters;
1018
1019 fNoOfClusters++;
1020 fClusterROC++;
3fe61b77 1021}
1022
3d0c7d6d 1023//_____________________________________________________________________________
1024void AliTRDclusterizer::CalcAdditionalInfo(const MaxStruct &Max, Short_t *const signals, Int_t &nPadCount)
1025{
1026 // Look to the right
1027 Int_t ii = 1;
615f0826 1028 while (fDigits->GetData(Max.row, Max.col-ii, Max.time) >= fSigThresh) {
3d0c7d6d 1029 nPadCount++;
1030 ii++;
615f0826 1031 if (Max.col < ii) break;
3d0c7d6d 1032 }
1033 // Look to the left
1034 ii = 1;
615f0826 1035 while (fDigits->GetData(Max.row, Max.col+ii, Max.time) >= fSigThresh) {
3d0c7d6d 1036 nPadCount++;
1037 ii++;
615f0826 1038 if (Max.col+ii >= fColMax) break;
3d0c7d6d 1039 }
1040
1041 // Store the amplitudes of the pads in the cluster for later analysis
1042 // and check whether one of these pads is masked in the database
615f0826 1043 signals[2]=Max.signals[0];
1044 signals[3]=Max.signals[1];
1045 signals[4]=Max.signals[2];
1046 Float_t gain;
3d0c7d6d 1047 for(Int_t i = 0; i<2; i++)
1048 {
615f0826 1049 if(Max.col+i >= 3){
1050 gain = fCalGainFactorDetValue * fCalGainFactorROC->GetValue(Max.col-3+i,Max.row);
1051 signals[i] = (Short_t)((fDigits->GetData(Max.row, Max.col-3+i, Max.time) - fBaseline) / gain + 0.5f);
1052 }
1053 if(Max.col+3-i < fColMax){
1054 gain = fCalGainFactorDetValue * fCalGainFactorROC->GetValue(Max.col+3-i,Max.row);
1055 signals[6-i] = (Short_t)((fDigits->GetData(Max.row, Max.col+3-i, Max.time) - fBaseline) / gain + 0.5f);
1056 }
3d0c7d6d 1057 }
1058 /*for (Int_t jPad = Max.Col-3; jPad <= Max.Col+3; jPad++) {
1059 if ((jPad >= 0) && (jPad < fColMax))
1060 signals[jPad-Max.Col+3] = TMath::Nint(fDigits->GetData(Max.Row,jPad,Max.Time));
1061 }*/
1062}
1063
1064//_____________________________________________________________________________
e7295a3a 1065void AliTRDclusterizer::AddClusterToArray(AliTRDcluster* cluster)
3d0c7d6d 1066{
1067 //
1068 // Add a cluster to the array
1069 //
1070
1071 Int_t n = RecPoints()->GetEntriesFast();
1072 if(n!=fNoOfClusters)AliError(Form("fNoOfClusters != RecPoints()->GetEntriesFast %i != %i \n", fNoOfClusters, n));
1073 new((*RecPoints())[n]) AliTRDcluster(*cluster);
1074}
1075
a5b99acd 1076//_____________________________________________________________________________
1077void AliTRDclusterizer::AddTrackletsToArray()
1078{
1079 //
1080 // Add the online tracklets of this chamber to the array
1081 //
1082
1083 UInt_t* trackletword;
1084 for(Int_t side=0; side<2; side++)
1085 {
1086 Int_t trkl=0;
1087 trackletword=fTrackletContainer[side];
97be9934 1088 while(trackletword[trkl]>0){
a5b99acd 1089 Int_t n = TrackletsArray()->GetEntriesFast();
35943b1e 1090 AliTRDtrackletWord tmp(trackletword[trkl]);
1091 new((*TrackletsArray())[n]) AliTRDcluster(&tmp,fDet,fVolid);
a5b99acd 1092 trkl++;
97be9934 1093 }
a5b99acd 1094 }
1095}
1096
3fe61b77 1097//_____________________________________________________________________________
6cd9a21f 1098Bool_t AliTRDclusterizer::AddLabels()
3551db50 1099{
1100 //
3fe61b77 1101 // Add the track indices to the found clusters
3551db50 1102 //
1103
3fe61b77 1104 const Int_t kNclus = 3;
1105 const Int_t kNdict = AliTRDdigitsManager::kNDict;
1106 const Int_t kNtrack = kNdict * kNclus;
1107
1108 Int_t iClusterROC = 0;
1109
1110 Int_t row = 0;
1111 Int_t col = 0;
1112 Int_t time = 0;
1113 Int_t iPad = 0;
1114
1115 // Temporary array to collect the track indices
6cd9a21f 1116 Int_t *idxTracks = new Int_t[kNtrack*fClusterROC];
3fe61b77 1117
1118 // Loop through the dictionary arrays one-by-one
1119 // to keep memory consumption low
b65e5048 1120 AliTRDarrayDictionary *tracksIn = 0; //mod
3fe61b77 1121 for (Int_t iDict = 0; iDict < kNdict; iDict++) {
1122
1123 // tracksIn should be expanded beforehand!
6cd9a21f 1124 tracksIn = (AliTRDarrayDictionary *) fDigitsManager->GetDictionary(fDet,iDict);
3fe61b77 1125
1126 // Loop though the clusters found in this ROC
6cd9a21f 1127 for (iClusterROC = 0; iClusterROC < fClusterROC; iClusterROC++) {
317b5644 1128
3fe61b77 1129 AliTRDcluster *cluster = (AliTRDcluster *)
b65e5048 1130 RecPoints()->UncheckedAt(firstClusterROC+iClusterROC);
3fe61b77 1131 row = cluster->GetLabel(0);
1132 col = cluster->GetLabel(1);
1133 time = cluster->GetLabel(2);
1134
1135 for (iPad = 0; iPad < kNclus; iPad++) {
b65e5048 1136 Int_t iPadCol = col - 1 + iPad;
1137 Int_t index = tracksIn->GetData(row,iPadCol,time); //Modification of -1 in Track
1138 idxTracks[3*iPad+iDict + iClusterROC*kNtrack] = index;
3fe61b77 1139 }
1140
1141 }
1142
1143 }
1144
1145 // Copy the track indices into the cluster
1146 // Loop though the clusters found in this ROC
6cd9a21f 1147 for (iClusterROC = 0; iClusterROC < fClusterROC; iClusterROC++) {
317b5644 1148
3fe61b77 1149 AliTRDcluster *cluster = (AliTRDcluster *)
1150 RecPoints()->UncheckedAt(firstClusterROC+iClusterROC);
1151 cluster->SetLabel(-9999,0);
1152 cluster->SetLabel(-9999,1);
1153 cluster->SetLabel(-9999,2);
1154
1155 cluster->AddTrackIndex(&idxTracks[iClusterROC*kNtrack]);
1156
1157 }
1158
1159 delete [] idxTracks;
1160
1161 return kTRUE;
1162
1163}
1164
3fe61b77 1165//_____________________________________________________________________________
e7295a3a 1166Float_t AliTRDclusterizer::Unfold(Double_t eps, Int_t layer, const Double_t *const padSignal) const
3fe61b77 1167{
1168 //
1169 // Method to unfold neighbouring maxima.
1170 // The charge ratio on the overlapping pad is calculated
1171 // until there is no more change within the range given by eps.
1172 // The resulting ratio is then returned to the calling method.
1173 //
1174
1175 AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
6d50f529 1176 if (!calibration) {
3fe61b77 1177 AliError("No AliTRDcalibDB instance available\n");
1178 return kFALSE;
6d50f529 1179 }
3fe61b77 1180
1181 Int_t irc = 0;
1182 Int_t itStep = 0; // Count iteration steps
1183
1184 Double_t ratio = 0.5; // Start value for ratio
1185 Double_t prevRatio = 0.0; // Store previous ratio
1186
1187 Double_t newLeftSignal[3] = { 0.0, 0.0, 0.0 }; // Array to store left cluster signal
1188 Double_t newRightSignal[3] = { 0.0, 0.0, 0.0 }; // Array to store right cluster signal
1189 Double_t newSignal[3] = { 0.0, 0.0, 0.0 };
1190
1191 // Start the iteration
1192 while ((TMath::Abs(prevRatio - ratio) > eps) && (itStep < 10)) {
1193
1194 itStep++;
1195 prevRatio = ratio;
1196
1197 // Cluster position according to charge ratio
1198 Double_t maxLeft = (ratio*padSignal[2] - padSignal[0])
064d808d 1199 / (padSignal[0] + padSignal[1] + ratio * padSignal[2]);
3fe61b77 1200 Double_t maxRight = (padSignal[4] - (1-ratio)*padSignal[2])
1201 / ((1.0 - ratio)*padSignal[2] + padSignal[3] + padSignal[4]);
1202
1203 // Set cluster charge ratio
064d808d 1204 irc = calibration->PadResponse(1.0, maxLeft, layer, newSignal);
3fe61b77 1205 Double_t ampLeft = padSignal[1] / newSignal[1];
064d808d 1206 irc = calibration->PadResponse(1.0, maxRight, layer, newSignal);
3fe61b77 1207 Double_t ampRight = padSignal[3] / newSignal[1];
1208
1209 // Apply pad response to parameters
053767a4 1210 irc = calibration->PadResponse(ampLeft ,maxLeft ,layer,newLeftSignal );
1211 irc = calibration->PadResponse(ampRight,maxRight,layer,newRightSignal);
3fe61b77 1212
1213 // Calculate new overlapping ratio
1214 ratio = TMath::Min((Double_t) 1.0
1215 ,newLeftSignal[2] / (newLeftSignal[2] + newRightSignal[0]));
1216
b43a3e17 1217 }
88719a08 1218
3fe61b77 1219 return ratio;
1220
1221}
88719a08 1222
3fe61b77 1223//_____________________________________________________________________________
064d808d 1224void AliTRDclusterizer::TailCancelation()
3fe61b77 1225{
1226 //
1227 // Applies the tail cancelation and gain factors:
1b3a719f 1228 // Transform fDigits to fDigits
3fe61b77 1229 //
1230
1231 Int_t iRow = 0;
1232 Int_t iCol = 0;
1233 Int_t iTime = 0;
1234
1b3a719f 1235 Double_t *inADC = new Double_t[fTimeTotal]; // ADC data before tail cancellation
064d808d 1236 Double_t *outADC = new Double_t[fTimeTotal]; // ADC data after tail cancellation
1b3a719f 1237
a2fbb6ec 1238 TTreeSRedirector *fDebugStream = fReconstructor->GetDebugStream(AliTRDrecoParam::kClusterizer);
cb5522da 1239 Bool_t debugStreaming = fReconstructor->GetRecoParam()->GetStreamLevel(AliTRDrecoParam::kClusterizer) > 7 && fReconstructor->IsDebugStreaming();
064d808d 1240 while(fIndexes->NextRCIndex(iRow, iCol))
3fe61b77 1241 {
f5375dcb 1242 Bool_t corrupted = kFALSE;
5a8db426 1243 if (fCalPadStatusROC->GetStatus(iCol, iRow)) corrupted = kTRUE;
1244
cb5522da 1245 // Save data into the temporary processing array and substract the baseline,
1246 // since DeConvExp does not expect a baseline
064d808d 1247 for (iTime = 0; iTime < fTimeTotal; iTime++)
cb5522da 1248 inADC[iTime] = fDigits->GetData(iRow,iCol,iTime)-fBaseline;
5a8db426 1249
cb5522da 1250 if(debugStreaming){
1251 for (iTime = 0; iTime < fTimeTotal; iTime++)
5a8db426 1252 (*fDebugStream) << "TailCancellation"
1253 << "col=" << iCol
1254 << "row=" << iRow
1255 << "time=" << iTime
1256 << "inADC=" << inADC[iTime]
1257 << "outADC=" << outADC[iTime]
1258 << "corrupted=" << corrupted
1259 << "\n";
cb5522da 1260 }
5a8db426 1261
b65e5048 1262 if (!corrupted)
1263 {
1264 // Apply the tail cancelation via the digital filter
1265 // (only for non-coorupted pads)
1b3a719f 1266 DeConvExp(&inADC[0],&outADC[0],fTimeTotal,fReconstructor->GetRecoParam() ->GetTCnexp());
b65e5048 1267 }
5a8db426 1268 else memcpy(&outADC[0],&inADC[0],fTimeTotal*sizeof(inADC[0]));
3fe61b77 1269
cb5522da 1270 // Save tailcancalled data and add the baseline
5a8db426 1271 for(iTime = 0; iTime < fTimeTotal; iTime++)
cb5522da 1272 fDigits->SetData(iRow,iCol,iTime,(Short_t)(outADC[iTime] + fBaseline + 0.5));
5a8db426 1273
3fe61b77 1274 } // while irow icol
1b3a719f 1275
3fe61b77 1276 delete [] inADC;
1277 delete [] outADC;
1278
1279 return;
1280
1281}
1282
1283//_____________________________________________________________________________
064d808d 1284void AliTRDclusterizer::DeConvExp(const Double_t *const source, Double_t *const target
1285 ,const Int_t n, const Int_t nexp)
3fe61b77 1286{
1287 //
1288 // Tail cancellation by deconvolution for PASA v4 TRF
1289 //
1290
1291 Double_t rates[2];
1292 Double_t coefficients[2];
1293
1294 // Initialization (coefficient = alpha, rates = lambda)
1295 Double_t r1 = 1.0;
1296 Double_t r2 = 1.0;
1297 Double_t c1 = 0.5;
1298 Double_t c2 = 0.5;
1299
1300 if (nexp == 1) { // 1 Exponentials
1301 r1 = 1.156;
1302 r2 = 0.130;
1303 c1 = 0.066;
1304 c2 = 0.000;
1305 }
1306 if (nexp == 2) { // 2 Exponentials
181c7f7e 1307 Double_t par[4];
a2fbb6ec 1308 fReconstructor->GetRecoParam()->GetTCParams(par);
181c7f7e 1309 r1 = par[0];//1.156;
1310 r2 = par[1];//0.130;
1311 c1 = par[2];//0.114;
1312 c2 = par[3];//0.624;
3fe61b77 1313 }
1314
1315 coefficients[0] = c1;
1316 coefficients[1] = c2;
1317
1318 Double_t dt = 0.1;
1319
1320 rates[0] = TMath::Exp(-dt/(r1));
1321 rates[1] = TMath::Exp(-dt/(r2));
1322
1323 Int_t i = 0;
1324 Int_t k = 0;
1325
1326 Double_t reminder[2];
1327 Double_t correction = 0.0;
1328 Double_t result = 0.0;
1329
1330 // Attention: computation order is important
1331 for (k = 0; k < nexp; k++) {
1332 reminder[k] = 0.0;
1333 }
1334
1335 for (i = 0; i < n; i++) {
1336
1337 result = (source[i] - correction); // No rescaling
1338 target[i] = result;
1339
1340 for (k = 0; k < nexp; k++) {
1341 reminder[k] = rates[k] * (reminder[k] + coefficients[k] * result);
1342 }
1343
1344 correction = 0.0;
1345 for (k = 0; k < nexp; k++) {
1346 correction += reminder[k];
1347 }
1348
1349 }
6d50f529 1350
1351}
1352
1353//_____________________________________________________________________________
1354void AliTRDclusterizer::ResetRecPoints()
1355{
1356 //
1357 // Resets the list of rec points
1358 //
1359
1360 if (fRecPoints) {
1361 fRecPoints->Delete();
48f8adf3 1362 delete fRecPoints;
6d50f529 1363 }
6d50f529 1364}
1365
1366//_____________________________________________________________________________
66f6bfd9 1367TClonesArray *AliTRDclusterizer::RecPoints()
6d50f529 1368{
1369 //
1370 // Returns the list of rec points
1371 //
1372
1373 if (!fRecPoints) {
48f8adf3 1374 if(!(fRecPoints = AliTRDReconstructor::GetClusters())){
8ae98148 1375 // determine number of clusters which has to be allocated
1376 Float_t nclusters = fReconstructor->GetRecoParam()->GetNClusters();
8ae98148 1377
1378 fRecPoints = new TClonesArray("AliTRDcluster", Int_t(nclusters));
48f8adf3 1379 }
1380 //SetClustersOwner(kTRUE);
1381 AliTRDReconstructor::SetClusters(0x0);
6d50f529 1382 }
6d50f529 1383 return fRecPoints;
eb52b657 1384
3551db50 1385}
fc546d21 1386
a5b99acd 1387//_____________________________________________________________________________
1388TClonesArray *AliTRDclusterizer::TrackletsArray()
1389{
1390 //
1391 // Returns the list of rec points
1392 //
1393
1394 if (!fTracklets && fReconstructor->IsProcessingTracklets()) {
94ac01ef 1395 fTracklets = new TClonesArray("AliTRDcluster", 2*MAXTRACKLETSPERHC);
a5b99acd 1396 //SetClustersOwner(kTRUE);
1397 //AliTRDReconstructor::SetTracklets(0x0);
1398 }
1399 return fTracklets;
1400
1401}
1402