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 // Alice segment manager class //
22 ////////////////////////////////////////////////////////////////////////////
27 #include <TClonesArray.h>
28 #include <TDirectory.h>
34 #include "AliTRDgeometry.h"
35 #include "AliTRDsegmentArray.h"
36 #include "AliTRDsegmentID.h"
37 #include "AliTRDdataArray.h"
38 #include "AliTRDarrayI.h"
40 ClassImp(AliTRDsegmentArray)
42 //_____________________________________________________________________________
43 AliTRDsegmentArray::AliTRDsegmentArray()
53 // Default constructor
58 //_____________________________________________________________________________
59 AliTRDsegmentArray::AliTRDsegmentArray(const char *classname, Int_t n)
69 // Create an array of objects of <classname>. The class must inherit from
70 // AliTRDsegmentID. The second argument sets the number of entries in
74 AliTRDdataArray *dataArray;
78 if (MakeArray(n) == kFALSE) {
79 AliError(Form("Cannot allocate %d segments in memory",n));
83 for (Int_t i = 0; i < n; i++) {
84 dataArray = (AliTRDdataArray *) AddSegment(i);
89 //_____________________________________________________________________________
90 AliTRDsegmentArray::AliTRDsegmentArray(AliTRDsegmentArray &a)
93 ,fTreeIndex(a.fTreeIndex)
94 ,fNSegment(a.fNSegment)
100 // AliTRDsegmentArray copy constructor
107 //_____________________________________________________________________________
108 AliTRDsegmentArray::~AliTRDsegmentArray()
111 // AliTRDsegmentArray destructor
128 //_____________________________________________________________________________
129 AliTRDsegmentArray &AliTRDsegmentArray::operator=(const AliTRDsegmentArray &a)
132 // Assignment operator
135 if (this != &a) ((AliTRDsegmentArray &) a).Copy(*this);
140 //_____________________________________________________________________________
141 void AliTRDsegmentArray::Copy(TObject &a) const
149 fSegment->Copy(*((AliTRDsegmentArray &) a).fSegment);
150 fTreeIndex->Copy(*((AliTRDsegmentArray &) a).fTreeIndex);
151 fClass->Copy(*((AliTRDsegmentArray &) a).fClass);
153 ((AliTRDsegmentArray &) a).fNSegment = fNSegment;
157 //_____________________________________________________________________________
158 Bool_t AliTRDsegmentArray::SetClass(const Char_t *classname)
161 // Sets the classname of the stored object
178 AliFatal("ROOT system not initialized");
182 fClass = gROOT->GetClass(classname);
184 AliError(Form("%s is not a valid class name",classname));
187 if (!fClass->InheritsFrom(AliTRDsegmentID::Class())) {
188 AliError(Form("%s does not inherit from AliTRDsegmentID",classname));
196 //_____________________________________________________________________________
197 AliTRDsegmentID *AliTRDsegmentArray::NewSegment()
200 // Create a new object according to the class information
207 AliTRDsegmentID *segment = (AliTRDsegmentID *) fClass->New();
218 //_____________________________________________________________________________
219 Bool_t AliTRDsegmentArray::AddSegment(AliTRDsegmentID *segment)
222 // Add a segment to the array
235 if (!(segment->IsA()->InheritsFrom(fClass))) {
236 AliError(Form("added class %s is not of proper type"
237 ,segment->IsA()->GetName()));
241 fSegment->AddAt(segment,segment->GetID());
242 fNSegment = fSegment->GetLast() + 1;
248 //_____________________________________________________________________________
249 AliTRDsegmentID *AliTRDsegmentArray::AddSegment(Int_t index)
252 // Add a segment to the array
262 AliTRDsegmentID *segment = NewSegment();
267 fSegment->AddAt(segment,index);
268 segment->SetID(index);
269 fNSegment = fSegment->GetLast() + 1;
275 //_____________________________________________________________________________
276 Bool_t AliTRDsegmentArray::MakeArray(Int_t n)
279 // Create an array of pointers to the segments
286 if (fTreeIndex) delete fTreeIndex;
288 fSegment = new TObjArray(n);
289 fTreeIndex = new AliTRDarrayI();
292 if ((fSegment) && (fTreeIndex)) {
301 //_____________________________________________________________________________
302 void AliTRDsegmentArray::ClearSegment(Int_t index)
305 // Remove a segment from the active memory
308 if (fSegment->At(index)) {
309 fClass->Destructor(fSegment->RemoveAt(index));
314 //_____________________________________________________________________________
315 void AliTRDsegmentArray::MakeTree(char *file)
318 // Create a tree for the segment
321 AliTRDsegmentID *psegment = NewSegment();
327 fTree = new TTree("Segment Tree","Tree with segments");
328 fBranch = fTree->Branch("Segment",psegment->IsA()->GetName(),&psegment,64000);
331 fBranch->SetFile(file);
338 //_____________________________________________________________________________
339 Bool_t AliTRDsegmentArray::ConnectTree(const char *treeName)
342 // Connect a tree from current directory
351 fTree = (TTree *) gDirectory->Get(treeName);
355 fBranch = fTree->GetBranch("Segment");
360 MakeDictionary(TMath::Max(fNSegment,Int_t(fTree->GetEntries())));
366 //_____________________________________________________________________________
367 AliTRDsegmentID *AliTRDsegmentArray::LoadSegment(Int_t index)
370 // Load a segment with index <index> into the memory
373 if (fTreeIndex == 0) {
374 MakeDictionary(3000);
377 // First try to load dictionary
378 if (fTreeIndex == 0) {
384 if (index > fTreeIndex->fN) {
388 AliTRDsegmentID *s = (AliTRDsegmentID *) fSegment->At(index);
395 Int_t treeIndex = (*fTreeIndex)[index];
402 fBranch->SetAddress(&s);
403 fTree->GetEvent(treeIndex);
404 fSegment->AddAt((TObject*) s, index);
414 //_____________________________________________________________________________
415 AliTRDsegmentID *AliTRDsegmentArray::LoadEntry(Int_t index)
418 // Load a segment at position <index> in the tree into the memory
424 if (index > fTree->GetEntries()) {
428 AliTRDsegmentID *s = NewSegment();
430 fBranch->SetAddress(&s);
431 fTree->GetEvent(index);
437 Int_t nindex = s->GetID();
438 ClearSegment(nindex);
439 fSegment->AddAt((TObject *) s, nindex);
445 //_____________________________________________________________________________
446 void AliTRDsegmentArray::StoreSegment(Int_t index)
449 // Make a segment persistent
452 const AliTRDsegmentID *kSegment = (*this)[index];
459 fBranch->SetAddress(&kSegment);
464 //_____________________________________________________________________________
465 Bool_t AliTRDsegmentArray::MakeDictionary(Int_t size)
468 // Create an index table for the tree
478 fTreeIndex = new AliTRDarrayI();
479 fTreeIndex->Set(size);
481 AliTRDsegmentID segment;
482 AliTRDsegmentID *psegment = &segment;
484 fBranch->SetAddress(&psegment);
485 TBranch *brindix = fTree->GetBranch("fSegmentID");
487 Int_t nevent = (Int_t) fTree->GetEntries();
488 for (Int_t i = 0; i < nevent; i++) {
489 brindix->GetEvent(i);
490 Int_t treeIndex = segment.GetID();
491 if (fTreeIndex->fN < treeIndex) {
492 fTreeIndex->Expand(Int_t (Float_t(treeIndex) * 1.5) + 1);
494 (*fTreeIndex)[treeIndex] = i + 1;
501 //_____________________________________________________________________________
502 const AliTRDsegmentID * AliTRDsegmentArray::operator[](Int_t i) const
505 // Returns a segment with the given index <i>
513 return (AliTRDsegmentID *) fSegment->At(i);
517 //_____________________________________________________________________________
518 const AliTRDsegmentID *AliTRDsegmentArray::At(Int_t i) const
521 // Returns a segment with the given index <i>
529 return (AliTRDsegmentID *) fSegment->At(i);
533 //_____________________________________________________________________________
534 void AliTRDsegmentArray::Delete()
537 // Deletes all detector segments from the array
540 for (Int_t iDet = 0; iDet < fNSegment; iDet++) {
546 //_____________________________________________________________________________
547 Bool_t AliTRDsegmentArray::LoadArray(const Char_t *branchname, TTree *tree)
550 // Loads all segments of the array from the branch <branchname> of
551 // the digits tree <tree>
557 AliError("Digits tree is not defined\n");
562 fBranch = fTree->GetBranch(branchname);
564 AliError(Form("Branch %s is not defined\n",branchname));
568 // Loop through all segments and read them from the tree
569 Bool_t status = kTRUE;
570 for (Int_t iSegment = 0; iSegment < fNSegment; iSegment++) {
571 AliTRDdataArray *dataArray = (AliTRDdataArray *) fSegment->At(iSegment);
576 fBranch->SetAddress(&dataArray);
577 fBranch->GetEntry(iSegment);
584 //_____________________________________________________________________________
585 Bool_t AliTRDsegmentArray::StoreArray(const Char_t *branchname, TTree *tree)
588 // Stores all segments of the array in the branch <branchname> of
589 // the digits tree <tree>
595 AliError("Digits tree is not defined\n");
600 fBranch = fTree->GetBranch(branchname);
602 AliError(Form("Branch %s is not defined\n",branchname));
606 // Loop through all segments and fill them into the tree
607 Bool_t status = kTRUE;
608 for (Int_t iSegment = 0; iSegment < fNSegment; iSegment++) {
609 const AliTRDdataArray *kDataArray =
610 (AliTRDdataArray *) AliTRDsegmentArray::At(iSegment);
615 fBranch->SetAddress(&kDataArray);
623 //_____________________________________________________________________________
624 AliTRDdataArray *AliTRDsegmentArray::GetDataArray(Int_t det) const
627 // Returns the data array for a given detector
630 return ((AliTRDdataArray *) AliTRDsegmentArray::At(det));
634 //_____________________________________________________________________________
635 AliTRDdataArray *AliTRDsegmentArray::GetDataArray(Int_t pla
640 // Returns the data array for a given detector
643 Int_t det = AliTRDgeometry::GetDetector(pla,cha,sec);
644 return GetDataArray(det);