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"
43 ClassImp(AliRawReaderChain)
45 TString AliRawReaderChain::fgSearchPath = "/alice/data";
47 AliRawReaderChain::AliRawReaderChain() :
51 // default constructor
54 AliRawReaderChain::AliRawReaderChain(const char* fileName) :
58 // create raw-reader objects which takes as an input a root chain
59 // either from the file list found in 'fileName' (IsCollection = true)
60 // or from entry list found in 'filename' (IsCollection = false)
61 // The entry-list syntax follows root convetion: filename.root/listname
63 fChain = new TChain("RAW");
65 TString fileNameStr = fileName;
66 if (fileNameStr.EndsWith(".xml")) {
68 TGridCollection *collection = NULL;
69 TPluginManager* pluginManager = gROOT->GetPluginManager();
70 TPluginHandler* pluginHandler = pluginManager->FindHandler("TGridCollection", "alice");
72 pluginManager->AddHandler("TGridCollection", "alice",
73 "AliXMLCollection", "ANALYSISalice", "AliXMLCollection(const char*)");
74 pluginHandler = pluginManager->FindHandler("TGridCollection", "alice");
76 gSystem->Load("libANALYSIS");
77 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
78 collection = (TGridCollection*)pluginHandler->ExecPlugin(1,fileNameStr.Data());
85 Bool_t elistsExist = kFALSE;
86 TEntryList *elist = new TEntryList();
87 while (collection->Next()) {
88 fChain->Add(collection->GetTURL(""));
89 TEntryList *list = (TEntryList *)collection->GetEntryList("");
91 list->SetTreeName("RAW");
92 list->SetFileName(collection->GetTURL(""));
98 fChain->SetEntryList(elist,"ne");
101 Info("AliRawReaderChain", "no entry lists found in %s. Using all entries", fileNameStr.Data());
105 else if (fileNameStr.EndsWith(".root")) {
107 TDirectory* dir = gDirectory;
108 TFile *listFile = TFile::Open(fileNameStr.Data());
110 if (!listFile || !listFile->IsOpen()) {
111 Error("AliRawReaderChain", "could not open file %s", fileNameStr.Data());
116 TEntryList *elist = NULL;
118 TIter nextkey(listFile->GetListOfKeys());
119 while ((key=(TKey*)nextkey())){
120 if (strcmp("TEntryList", key->GetClassName())==0){
121 elist = (TEntryList*)key->ReadObj();
125 Error("AliRawReaderChain", "no TEntryList found in %s", fileNameStr.Data());
130 TEntryList *templist = NULL;
131 TList *elists = elist->GetLists();
133 while((templist = (TEntryList*)next())){
134 Info("AliRawReaderChain", "%s added to the chain", templist->GetFileName());
135 fChain->Add(templist->GetFileName());
137 fChain->SetEntryList(elist,"ne");
141 TFileCollection collection("RAW",
142 "Collection with raw-data files",
145 if (!fChain->AddFileInfoList((TCollection*)(collection.GetList()))) {
146 Error("AliRawReaderChain","Bad file list in collection, the chain is empty");
152 fChain->SetBranchStatus("*",1);
153 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
156 AliRawReaderChain::AliRawReaderChain(TFileCollection *collection) :
160 // create raw-reader objects which takes as an input a root chain
161 // from a root file collection
163 fChain = new TChain("RAW");
164 if (!fChain->AddFileInfoList((TCollection*)(collection->GetList()))) {
165 Error("AliRawReaderChain","Bad file list in collection, the chain is empty");
170 fChain->SetBranchStatus("*",1);
171 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
174 AliRawReaderChain::AliRawReaderChain(TChain *chain) :
178 // create raw-reader objects which takes as an input a root chain
179 // from a root file collection
186 fChain->SetBranchStatus("*",1);
187 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
190 AliRawReaderChain::AliRawReaderChain(TEntryList *elist) :
194 // create raw-reader objects which takes as an input a root chain
195 // from a root file collection
202 fChain = new TChain("RAW");
204 TEntryList *templist = NULL;
205 TList *elists = elist->GetLists();
207 while((templist = (TEntryList*)next())){
208 fChain->Add(templist->GetFileName());
210 fChain->SetEntryList(elist,"ne");
212 fChain->SetBranchStatus("*",1);
213 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
216 AliRawReaderChain::AliRawReaderChain(Int_t runNumber) :
220 // create raw-reader objects which takes as an input a root chain
221 // with the raw-data files for a given run
222 // It queries alien FC in order to do that and therefore
223 // it needs alien API to be enabled
225 if (runNumber <= 0) {
226 Error("AliRawReaderChain","Bad run number:%d",runNumber);
230 if (!gGrid) TGrid::Connect("alien://");
236 if (fgSearchPath.IsNull()) fgSearchPath = "/alice/data";
237 TGridResult *res = gGrid->Query(fgSearchPath.Data(),Form("%09d/raw/*%09d*.root",runNumber,runNumber));
238 Int_t nFiles = res->GetEntries();
240 Error("AliRawReaderChain","No raw-data files found for run %d",runNumber);
246 fChain = new TChain("RAW");
247 for (Int_t i = 0; i < nFiles; i++) {
248 TString filename = res->GetKey(i, "turl");
249 if(filename == "") continue;
250 fChain->Add(filename.Data());
254 fChain->SetBranchStatus("*",1);
255 fChain->SetBranchAddress("rawevent",&fEvent,&fBranch);
259 AliRawReaderChain::AliRawReaderChain(const AliRawReaderChain& rawReader) :
260 AliRawReaderRoot(rawReader),
261 fChain(rawReader.fChain)
266 AliRawReaderChain& AliRawReaderChain::operator = (const AliRawReaderChain&
269 // assignment operator
271 this->~AliRawReaderChain();
272 new(this) AliRawReaderChain(rawReader);
276 AliRawReaderChain::~AliRawReaderChain()
278 // delete objects and close root file
286 Bool_t AliRawReaderChain::NextEvent()
288 // go to the next event in the root file
290 if (!fChain || !fChain->GetListOfFiles()->GetEntriesFast()) return kFALSE;
296 Long64_t treeEntry = fChain->LoadTree(fEventIndex+1);
299 if (fBranch->GetEntry(treeEntry) <= 0)
301 fEventHeader = fEvent->GetHeader();
303 } while (!IsEventSelected());
308 Bool_t AliRawReaderChain::RewindEvents()
310 // go back to the beginning of the root file
320 Bool_t AliRawReaderChain::GotoEvent(Int_t event)
322 // go to a particular event
323 // Uses the absolute event index inside the
324 // chain with raw data
326 if (!fChain || !fChain->GetListOfFiles()->GetEntriesFast()) return kFALSE;
331 Long64_t treeEntry = fChain->LoadTree(event);
334 if (fBranch->GetEntry(treeEntry) <= 0)
336 fEventHeader = fEvent->GetHeader();
342 Int_t AliRawReaderChain::GetNumberOfEvents() const
344 // Get the total number of events in the chain
347 if (!fChain) return -1;
349 return fChain->GetEntries();
352 void AliRawReaderChain::SetSearchPath(const char* path)
354 // set alien query search path
355 AliInfoGeneral("SetSearchPath",Form("Setting search path to \"%s\" (was \"%s\")",path,fgSearchPath.Data()));