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 /// This is a class for reading raw data from a root file.
22 /// The root file is expected to contain a tree of name "RAW" with
23 /// a branch of name "rawevent" which contains objects of type
26 /// The file name and the event number are arguments of the constructor
27 /// of AliRawReaderRoot.
29 ///////////////////////////////////////////////////////////////////////////////
33 #include <TTreeIndex.h>
35 #include "AliRawReaderRoot.h"
36 #include "AliRawVEvent.h"
37 #include "AliRawEventHeaderBase.h"
38 #include "AliRawVEquipment.h"
39 #include "AliRawEquipmentHeader.h"
40 #include "AliRawData.h"
43 ClassImp(AliRawReaderRoot)
44 Bool_t AliRawReaderRoot::fgUseOrder = kFALSE;
47 AliRawReaderRoot::AliRawReaderRoot() :
62 // default constructor
66 AliRawReaderRoot::AliRawReaderRoot(const char* fileName, Int_t eventNumber) :
69 fEventIndex(eventNumber),
81 // create an object to read digits from the given input file for the
82 // event with the given number
84 TDirectory* dir = gDirectory;
85 TString flStr = fileName;
86 if (flStr.BeginsWith("alien://") && !gGrid) TGrid::Connect("alien://");
87 fFile = TFile::Open(fileName);
89 if (!fFile || !fFile->IsOpen()) {
90 Error("AliRawReaderRoot", "could not open file %s", fileName);
94 TTree* tree = (TTree*) fFile->Get("RAW");
96 Error("AliRawReaderRoot", "no raw data tree found");
100 fBranch = tree->GetBranch("rawevent");
102 Error("AliRawReaderRoot", "no raw data branch found");
107 fBranch->SetAddress(&fEvent);
108 if (fEventIndex >= 0) {
109 if (fBranch->GetEntry(fEventIndex) <= 0) {
110 Error("AliRawReaderRoot", "no event with number %d found", fEventIndex);
114 fEventHeader = fEvent->GetHeader();
118 AliRawReaderRoot::AliRawReaderRoot(AliRawVEvent* event) :
123 fEventHeader(event->GetHeader()),
133 // create an object to read digits from the given raw event
134 if (!fEvent) fIsValid = kFALSE;
137 AliRawReaderRoot::AliRawReaderRoot(const AliRawReaderRoot& rawReader) :
138 AliRawReader(rawReader),
141 fEventIndex(rawReader.fEventIndex),
144 fSubEventIndex(rawReader.fSubEventIndex),
146 fEquipmentIndex(rawReader.fEquipmentIndex),
155 if (rawReader.fFile) {
156 TDirectory* dir = gDirectory;
157 fFile = TFile::Open(rawReader.fFile->GetName());
159 if (!fFile || !fFile->IsOpen()) {
160 Error("AliRawReaderRoot", "could not open file %s",
161 rawReader.fFile->GetName());
165 TTree* tree = (TTree*) fFile->Get("RAW");
167 Error("AliRawReaderRoot", "no raw data tree found");
171 fBranch = tree->GetBranch("rawevent");
173 Error("AliRawReaderRoot", "no raw data branch found");
178 fBranch->SetAddress(&fEvent);
179 if (fEventIndex >= 0) {
180 if (fBranch->GetEntry(fEventIndex) <= 0) {
181 Error("AliRawReaderRoot", "no event with number %d found",
186 fEventHeader = fEvent->GetHeader();
189 fEvent = rawReader.fEvent;
190 fEventHeader = rawReader.fEventHeader;
193 if (fSubEventIndex > 0) {
194 fSubEvent = fEvent->GetSubEvent(fSubEventIndex-1);
195 fEquipment = fSubEvent->GetEquipment(fEquipmentIndex);
196 fRawData = fEquipment->GetRawData();
198 fHeader = (AliRawDataHeader*) ((UChar_t*) fRawData->GetBuffer() +
199 ((UChar_t*) rawReader.fHeader -
200 (UChar_t*) rawReader.fRawData->GetBuffer()));
201 // Now check the version of the header
203 if (fHeader) version = fHeader->GetVersion();
206 fHeaderV3 = (AliRawDataHeaderV3*) ((UChar_t*) fRawData->GetBuffer() +
207 ((UChar_t*) rawReader.fHeaderV3 -
208 (UChar_t*) rawReader.fRawData->GetBuffer()));
210 fPosition = (UChar_t*) fRawData->GetBuffer() +
211 (rawReader.fPosition - (UChar_t*) rawReader.fRawData->GetBuffer());
212 fEnd = ((UChar_t*) fRawData->GetBuffer()) + fRawData->GetSize();
216 AliRawReaderRoot& AliRawReaderRoot::operator = (const AliRawReaderRoot&
219 // assignment operator
221 this->~AliRawReaderRoot();
222 new(this) AliRawReaderRoot(rawReader);
226 AliRawReaderRoot::~AliRawReaderRoot()
228 // delete objects and close root file
231 if (fEvent) delete fEvent;
237 const AliRawEventHeaderBase* AliRawReaderRoot::GetEventHeader() const
239 // Get the even header
240 // Return NULL in case of failure
244 UInt_t AliRawReaderRoot::GetType() const
246 // get the type from the event header
248 if (!fEventHeader) return 0;
249 return fEventHeader->Get("Type");
252 UInt_t AliRawReaderRoot::GetRunNumber() const
254 // get the run number from the event header
256 if (!fEventHeader) return 0;
257 return fEventHeader->Get("RunNb");
260 const UInt_t* AliRawReaderRoot::GetEventId() const
262 // get the event id from the event header
264 if (!fEventHeader) return NULL;
265 return fEventHeader->GetP("Id");
268 const UInt_t* AliRawReaderRoot::GetTriggerPattern() const
270 // get the trigger pattern from the event header
272 if (!fEventHeader) return NULL;
273 return fEventHeader->GetP("TriggerPattern");
276 const UInt_t* AliRawReaderRoot::GetDetectorPattern() const
278 // get the detector pattern from the event header
280 if (!fEventHeader) return NULL;
281 return fEventHeader->GetP("DetectorPattern");
284 const UInt_t* AliRawReaderRoot::GetAttributes() const
286 // get the type attributes from the event header
288 if (!fEventHeader) return NULL;
289 return fEventHeader->GetP("TypeAttribute");
292 const UInt_t* AliRawReaderRoot::GetSubEventAttributes() const
294 // get the type attributes from the sub event header
296 if (!fSubEvent) return NULL;
297 return fSubEvent->GetHeader()->GetP("TypeAttribute");
300 UInt_t AliRawReaderRoot::GetLDCId() const
302 // get the LDC Id from the event header
304 if (!fEvent || !fSubEvent) return 0;
305 return fSubEvent->GetHeader()->Get("LdcId");
308 UInt_t AliRawReaderRoot::GetGDCId() const
310 // get the GDC Id from the event header
312 if (!fEventHeader) return 0;
313 return fEventHeader->Get("GdcId");
316 UInt_t AliRawReaderRoot::GetTimestamp() const
318 if (!fEventHeader) return 0;
319 return fEventHeader->Get("Timestamp");
322 Int_t AliRawReaderRoot::GetEquipmentSize() const
324 // get the size of the equipment
326 if (!fEvent || !fEquipment || !fEquipment->GetEquipmentHeader()) return 0;
327 return fEquipment->GetEquipmentHeader()->GetEquipmentSize();
330 Int_t AliRawReaderRoot::GetEquipmentType() const
332 // get the type from the equipment header
334 if (!fEvent || !fEquipment || !fEquipment->GetEquipmentHeader()) return -1;
335 return fEquipment->GetEquipmentHeader()->GetEquipmentType();
338 Int_t AliRawReaderRoot::GetEquipmentId() const
340 // get the ID from the equipment header
342 if (!fEvent || !fEquipment || !fEquipment->GetEquipmentHeader()) return -1;
343 return fEquipment->GetEquipmentHeader()->GetId();
346 const UInt_t* AliRawReaderRoot::GetEquipmentAttributes() const
348 // get the attributes from the equipment header
350 if (!fEvent || !fEquipment || !fEquipment->GetEquipmentHeader()) return NULL;
351 return fEquipment->GetEquipmentHeader()->GetTypeAttribute();
354 Int_t AliRawReaderRoot::GetEquipmentElementSize() const
356 // get the basic element size from the equipment header
358 if (!fEvent || !fEquipment || !fEquipment->GetEquipmentHeader()) return 0;
359 return fEquipment->GetEquipmentHeader()->GetBasicSizeType();
362 Int_t AliRawReaderRoot::GetEquipmentHeaderSize() const
364 // get the size of the equipment header (28 bytes by default)
366 if (!fEvent || !fEquipment || !fEquipment->GetEquipmentHeader()) return 0;
367 return fEquipment->GetEquipmentHeader()->HeaderSize();
370 // _________________________________________________________________________
371 void AliRawReaderRoot::SwapData(const void* inbuf, const void* outbuf, UInt_t size) {
372 // The method swaps the contents of the
373 // raw-data event header
374 UInt_t intCount = (size+3)/sizeof(UInt_t);
376 UInt_t* buf = (UInt_t*) inbuf; // temporary integers buffer
377 for (UInt_t i=0; i<intCount; i++, buf++) {
378 UInt_t value = SwapWord(*buf);
380 memcpy((UInt_t*)outbuf+i, &value, size%sizeof(UInt_t));
382 memcpy((UInt_t*)outbuf+i, &value, sizeof(UInt_t));
385 // _________________________________________________________________________
387 Bool_t AliRawReaderRoot::ReadHeader()
389 // read a data header at the current position
390 // returns kFALSE if the data header could not be read
393 if (!fEvent) return kFALSE;
396 // skip payload (if event was not selected)
397 if (fCount > 0) fPosition += fCount;
399 // get the first or the next equipment if at the end of an equipment
400 if (!fEquipment || (fPosition >= fEnd)) {
402 // get the first or the next sub event if at the end of a sub event
403 if (!fSubEvent || (fEquipmentIndex >= fSubEvent->GetNEquipments())) {
405 // check for end of event data
406 if (fSubEventIndex >= fEvent->GetNSubEvents()) return kFALSE;
407 fSubEvent = fEvent->GetSubEvent(fSubEventIndex++);
409 // check the magic word of the sub event
410 if (!fSubEvent->GetHeader()->IsValid()) {
411 Error("ReadHeader", "wrong magic number in sub event!");
412 fSubEvent->GetHeader()->Dump();
413 fErrorCode = kErrMagic;
422 // get the next equipment and raw data
424 if (fEquipmentIndex >= fSubEvent->GetNEquipments()) {
428 fEquipment = fSubEvent->GetEquipment(fEquipmentIndex++);
429 if (!fEquipment) continue;
434 fRawData = fEquipment->GetRawData();
439 fPosition = (UChar_t*) fRawData->GetBuffer();
440 fEnd = ((UChar_t*) fRawData->GetBuffer()) + fRawData->GetSize();
443 // continue with the next equipment if no data left in the payload
444 if (fPosition >= fEnd) continue;
446 if (fRequireHeader) {
447 // check that there are enough bytes left for the data header
448 if (fPosition + sizeof(AliRawDataHeader) > fEnd) {
449 Error("ReadHeader", "could not read data header!");
450 Warning("ReadHeader", "skipping %ld bytes", fEnd - fPosition);
451 fEquipment->GetEquipmentHeader()->Dump();
454 fErrorCode = kErrNoDataHeader;
458 // "read" the data header
459 fHeader = (AliRawDataHeader*) fPosition;
460 // Now check the version of the header
462 if (fHeader) version=fHeader->GetVersion();
465 SwapData((void*) fHeader, (void*) fHeaderSwapped, sizeof(AliRawDataHeader));
466 fHeader=fHeaderSwapped;
468 if ((fPosition + fHeader->fSize) != fEnd) {
469 if (fHeader->fSize != 0xFFFFFFFF)
470 Warning("ReadHeader",
471 "Equipment %d : raw data size found in the header is wrong (%d != %ld)! Using the equipment size instead !",
472 fEquipment->GetEquipmentHeader()->GetId(),fHeader->fSize, fEnd - fPosition);
473 fHeader->fSize = fEnd - fPosition;
475 fPosition += sizeof(AliRawDataHeader);
477 } else if (version==3) {
478 fHeaderV3 = (AliRawDataHeaderV3*) fPosition;
480 SwapData((void*) fHeaderV3, (void*) fHeaderSwapped, sizeof(AliRawDataHeaderV3));
481 fHeaderV3=fHeaderSwappedV3;
483 if ((fPosition + fHeaderV3->fSize) != fEnd) {
484 if (fHeaderV3->fSize != 0xFFFFFFFF)
485 Warning("ReadHeader",
486 "Equipment %d : raw data size found in the header is wrong (%d != %ld)! Using the equipment size instead !",
487 fEquipment->GetEquipmentHeader()->GetId(),fHeader->fSize, fEnd - fPosition);
488 fHeaderV3->fSize = fEnd - fPosition;
490 fPosition += sizeof(AliRawDataHeaderV3);
495 if (fHeader && (fHeader->fSize != 0xFFFFFFFF)) {
496 fCount = fHeader->fSize - sizeof(AliRawDataHeader);
498 // check consistency of data size in the header and in the sub event
499 if (fPosition + fCount > fEnd) {
500 Error("ReadHeader", "size in data header exceeds event size!");
501 Warning("ReadHeader", "skipping %ld bytes", fEnd - fPosition);
502 fEquipment->GetEquipmentHeader()->Dump();
505 fErrorCode = kErrSize;
508 } else if (fHeaderV3 && (fHeaderV3->fSize != 0xFFFFFFFF)) {
509 fCount = fHeaderV3->fSize - sizeof(AliRawDataHeaderV3);
511 // check consistency of data size in the header and in the sub event
512 if (fPosition + fCount > fEnd) {
513 Error("ReadHeader", "size in data header exceeds event size!");
514 Warning("ReadHeader", "skipping %ld bytes", fEnd - fPosition);
515 fEquipment->GetEquipmentHeader()->Dump();
518 fErrorCode = kErrSize;
523 fCount = fEnd - fPosition;
526 } while (!IsSelected());
531 Bool_t AliRawReaderRoot::ReadNextData(UChar_t*& data)
533 // reads the next payload at the current position
534 // returns kFALSE if the data could not be read
537 while (fCount == 0) {
538 if (!ReadHeader()) return kFALSE;
546 Bool_t AliRawReaderRoot::ReadNext(UChar_t* data, Int_t size)
548 // reads the next block of data at the current position
549 // returns kFALSE if the data could not be read
552 if (fPosition + size > fEnd) {
553 Error("ReadNext", "could not read data!");
554 fErrorCode = kErrOutOfBounds;
557 memcpy(data, fPosition, size);
565 Bool_t AliRawReaderRoot::Reset()
567 // reset the current position to the beginning of the event
578 fPosition = fEnd = NULL;
583 Bool_t AliRawReaderRoot::NextEvent()
585 // go to the next event in the root file
587 if (!fBranch) return kFALSE;
589 // check if it uses order or not
590 if (fgUseOrder && !fIndex) MakeIndex();
596 fBranch->SetAddress(&fEvent);
597 Int_t entryToGet = fEventIndex + 1;
600 && entryToGet<fBranch->GetEntries()
601 && entryToGet>-1 ) entryToGet = fIndex[entryToGet];
602 if (fBranch->GetEntry(entryToGet) <= 0)
604 fEventHeader = fEvent->GetHeader();
606 } while (!IsEventSelected());
611 Bool_t AliRawReaderRoot::RewindEvents()
613 // go back to the beginning of the root file
615 if (!fBranch) return kFALSE;
621 fBranch->SetAddress(&fEvent);
626 Bool_t AliRawReaderRoot::GotoEvent(Int_t event)
628 // go to a particular event
629 // Uses the absolute event index inside the
632 if (!fBranch) return kFALSE;
634 // check if it uses order or not
635 if (fgUseOrder && !fIndex) MakeIndex();
640 fBranch->SetAddress(&fEvent);
641 Int_t entryToGet = event;
644 && entryToGet<fBranch->GetEntries()
645 && entryToGet>-1 ) entryToGet = fIndex[entryToGet];
646 if (fBranch->GetEntry(entryToGet) <= 0)
648 fEventHeader = fEvent->GetHeader();
654 Int_t AliRawReaderRoot::GetNumberOfEvents() const
656 // Get the total number of events in
659 if (!fBranch) return -1;
661 return fBranch->GetEntries();
664 Int_t AliRawReaderRoot::CheckData() const
666 // check the consistency of the data
668 if (!fEvent) return 0;
670 AliRawVEvent* subEvent = NULL;
671 Int_t subEventIndex = 0;
672 AliRawVEquipment* equipment = NULL;
673 Int_t equipmentIndex = 0;
674 UChar_t* position = 0;
679 // get the first or the next sub event if at the end of an equipment
680 if (!subEvent || (equipmentIndex >= subEvent->GetNEquipments())) {
682 // check for end of event data
683 if (subEventIndex >= fEvent->GetNSubEvents()) return result;
684 subEvent = fEvent->GetSubEvent(subEventIndex++);
686 // check the magic word of the sub event
687 if (!fSubEvent->GetHeader()->IsValid()) {
695 // get the next equipment and raw data
696 if (equipmentIndex >= subEvent->GetNEquipments()) {
700 equipment = subEvent->GetEquipment(equipmentIndex++);
701 if (!equipment) continue;
702 AliRawData* rawData = equipment->GetRawData();
703 if (!rawData) continue;
704 position = (UChar_t*) rawData->GetBuffer();
705 end = ((UChar_t*) rawData->GetBuffer()) + rawData->GetSize();
707 // continue with the next sub event if no data left in the payload
708 if (position >= end) continue;
710 if (fRequireHeader) {
711 // check that there are enough bytes left for the data header
712 if (position + sizeof(AliRawDataHeader) > end) {
713 result |= kErrNoDataHeader;
717 // check consistency of data size in the header and in the equipment
718 AliRawDataHeader* header = (AliRawDataHeader*) position;
719 if ((position + header->fSize) != end) {
720 if (header->fSize != 0xFFFFFFFF)
721 Warning("ReadHeader",
722 "Equipment %d : raw data size found in the header is wrong (%d != %ld)! Using the equipment size instead !",
723 equipment->GetEquipmentHeader()->GetId(),header->fSize, end - position);
724 header->fSize = end - position;
734 AliRawReader* AliRawReaderRoot::CloneSingleEvent() const
736 // Clones the current event and
737 // creates raw-reader for the cloned event
738 // Can be used in order to make asynchronious
739 // access to the current raw data within
740 // several threads (online event display/reco)
743 // Root formatted raw data
744 AliRawVEvent *gdcRootEvent = (AliRawVEvent*)fEvent->Clone();
745 for (Int_t ldcCounter=0; ldcCounter < gdcRootEvent->GetNSubEvents(); ldcCounter++) {
746 AliRawVEvent *ldcRootEvent = gdcRootEvent->GetSubEvent(ldcCounter);
747 AliRawVEvent *subEvent = fEvent->GetSubEvent(ldcCounter);
748 for (Int_t eqCounter=0; eqCounter < ldcRootEvent->GetNEquipments(); eqCounter++) {
749 AliRawVEquipment *equipment=ldcRootEvent->GetEquipment(eqCounter);
750 AliRawVEquipment *eq = subEvent->GetEquipment(eqCounter);
751 equipment->CloneRawData(eq->GetRawData());
754 // Reset original event and newly
756 gdcRootEvent->GetSubEvent(-1);
757 fEvent->GetSubEvent(-1);
758 return new AliRawReaderRoot(gdcRootEvent);
763 void AliRawReaderRoot::MakeIndex() {
766 TTree * rawTree = fBranch->GetTree();
768 rawTree->BuildIndex("-fEvtHdrs[0].fSize"); // Minus sign to get largest first
769 TTreeIndex * treeInd = (TTreeIndex*)rawTree->GetTreeIndex();
770 if (treeInd) fIndex = treeInd->GetIndex();