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 **************************************************************************/
16 ///////////////////////////////////////////////////////////////////////////////
18 /// This is a class for reading raw data from a root chain.
19 /// There are two constructors available - one from a text file containing the
20 /// list of root raw-data files to be processed and one directly from
23 /// cvetan.cheshkov@cern.ch 29/07/2008
25 ///////////////////////////////////////////////////////////////////////////////
28 #include <TFileCollection.h>
29 #include <TEntryList.h>
30 #include "TGridCollection.h"
31 #include <TPluginManager.h>
37 #include <TGridResult.h>
39 #include "AliRawReaderChain.h"
40 #include "AliRawVEvent.h"
42 ClassImp(AliRawReaderChain)
44 AliRawReaderChain::AliRawReaderChain() :
48 // default constructor
51 AliRawReaderChain::AliRawReaderChain(const char* fileName) :
55 // create raw-reader objects which takes as an input a root chain
56 // either from the file list found in 'fileName' (IsCollection = true)
57 // or from entry list found in 'filename' (IsCollection = false)
58 // The entry-list syntax follows root convetion: filename.root/listname
60 fChain = new TChain("RAW");
62 TString fileNameStr = fileName;
63 if (fileNameStr.EndsWith(".xml")) {
65 TGridCollection *collection = NULL;
66 TPluginManager* pluginManager = gROOT->GetPluginManager();
67 TPluginHandler* pluginHandler = pluginManager->FindHandler("TGridCollection", "alice");
69 pluginManager->AddHandler("TGridCollection", "alice",
70 "AliXMLCollection", "ANALYSISalice", "AliXMLCollection(const char*)");
71 pluginHandler = pluginManager->FindHandler("TGridCollection", "alice");
73 gSystem->Load("libANALYSIS");
74 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
75 collection = (TGridCollection*)pluginHandler->ExecPlugin(1,fileNameStr.Data());
82 Bool_t elistsExist = kFALSE;
83 TEntryList *elist = new TEntryList();
84 while (collection->Next()) {
85 fChain->Add(collection->GetTURL(""));
86 TEntryList *list = (TEntryList *)collection->GetEntryList("");
88 list->SetTreeName("RAW");
89 list->SetFileName(collection->GetTURL(""));
95 fChain->SetEntryList(elist,"ne");
98 Info("AliRawReaderChain", "no entry lists found in %s. Using all entries", fileNameStr.Data());
102 else if (fileNameStr.EndsWith(".root")) {
104 TDirectory* dir = gDirectory;
105 TFile *listFile = TFile::Open(fileNameStr.Data());
107 if (!listFile || !listFile->IsOpen()) {
108 Error("AliRawReaderChain", "could not open file %s", fileNameStr.Data());
113 TEntryList *elist = NULL;
115 TIter nextkey(listFile->GetListOfKeys());
116 while ((key=(TKey*)nextkey())){
117 if (strcmp("TEntryList", key->GetClassName())==0){
118 elist = (TEntryList*)key->ReadObj();
122 Error("AliRawReaderChain", "no TEntryList found in %s", fileNameStr.Data());
127 TEntryList *templist = NULL;
128 TList *elists = elist->GetLists();
130 while((templist = (TEntryList*)next())){
131 Info("AliRawReaderChain", "%s added to the chain", templist->GetFileName());
132 fChain->Add(templist->GetFileName());
134 fChain->SetEntryList(elist,"ne");
138 TFileCollection collection("RAW",
139 "Collection with raw-data files",
142 if (!fChain->AddFileInfoList((TCollection*)(collection.GetList()))) {
143 Error("AliRawReaderChain","Bad file list in collection, the chain is empty");
149 fChain->SetBranchStatus("*",1);
150 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
153 AliRawReaderChain::AliRawReaderChain(TFileCollection *collection) :
157 // create raw-reader objects which takes as an input a root chain
158 // from a root file collection
160 fChain = new TChain("RAW");
161 if (!fChain->AddFileInfoList((TCollection*)(collection->GetList()))) {
162 Error("AliRawReaderChain","Bad file list in collection, the chain is empty");
167 fChain->SetBranchStatus("*",1);
168 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
171 AliRawReaderChain::AliRawReaderChain(TChain *chain) :
175 // create raw-reader objects which takes as an input a root chain
176 // from a root file collection
178 if (!fChain) fIsValid = kFALSE;
180 fChain->SetBranchStatus("*",1);
181 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
184 AliRawReaderChain::AliRawReaderChain(TEntryList *elist) :
188 // create raw-reader objects which takes as an input a root chain
189 // from a root file collection
191 if (!elist) fIsValid = kFALSE;
193 fChain = new TChain("RAW");
195 TEntryList *templist = NULL;
196 TList *elists = elist->GetLists();
198 while((templist = (TEntryList*)next())){
199 fChain->Add(templist->GetFileName());
201 fChain->SetEntryList(elist,"ne");
203 fChain->SetBranchStatus("*",1);
204 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
207 AliRawReaderChain::AliRawReaderChain(Int_t runNumber) :
211 // create raw-reader objects which takes as an input a root chain
212 // with the raw-data files for a given run
213 // It queries alien FC in order to do that and therefore
214 // it needs alien API to be enabled
216 if (runNumber <= 0) {
217 Error("AliRawReaderChain","Bad run number:%d",runNumber);
221 if (!gGrid) TGrid::Connect("alien://");
227 TGridResult *res = gGrid->Query("/alice/data",Form("%09d/raw/*%09d*0.root",runNumber,runNumber));
228 Int_t nFiles = res->GetEntries();
230 Error("AliRawReaderChain","No raw-data files found for run %d",runNumber);
236 fChain = new TChain("RAW");
237 for (Int_t i = 0; i < nFiles; i++) {
238 TString filename = res->GetKey(i, "turl");
239 if(filename == "") continue;
240 fChain->Add(filename.Data());
244 fChain->SetBranchStatus("*",1);
245 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
249 AliRawReaderChain::AliRawReaderChain(const AliRawReaderChain& rawReader) :
250 AliRawReaderRoot(rawReader),
251 fChain(rawReader.fChain)
256 AliRawReaderChain& AliRawReaderChain::operator = (const AliRawReaderChain&
259 // assignment operator
261 this->~AliRawReaderChain();
262 new(this) AliRawReaderChain(rawReader);
266 AliRawReaderChain::~AliRawReaderChain()
268 // delete objects and close root file
276 Bool_t AliRawReaderChain::NextEvent()
278 // go to the next event in the root file
280 if (!fChain || !fChain->GetListOfFiles()->GetEntriesFast()) return kFALSE;
286 Long64_t treeEntry = fChain->LoadTree(fEventIndex+1);
289 if (fBranch->GetEntry(treeEntry) <= 0)
291 fEventHeader = fEvent->GetHeader();
293 } while (!IsEventSelected());
298 Bool_t AliRawReaderChain::RewindEvents()
300 // go back to the beginning of the root file
310 Bool_t AliRawReaderChain::GotoEvent(Int_t event)
312 // go to a particular event
313 // Uses the absolute event index inside the
314 // chain with raw data
316 if (!fChain || !fChain->GetListOfFiles()->GetEntriesFast()) return kFALSE;
321 Long64_t treeEntry = fChain->LoadTree(event);
324 if (fBranch->GetEntry(treeEntry) <= 0)
326 fEventHeader = fEvent->GetHeader();
332 Int_t AliRawReaderChain::GetNumberOfEvents() const
334 // Get the total number of events in the chain
337 if (!fChain) return -1;
339 return fChain->GetEntries();