1 /**************************************************************************
2 * Author: Panos Christakoglou. *
3 * Contributors are mentioned in the code where appropriate. *
5 * Permission to use, copy, modify and distribute this software and its *
6 * documentation strictly for non-commercial purposes is hereby granted *
7 * without fee, provided that the above copyright notice appears in all *
8 * copies and that both the copyright notice and this permission notice *
9 * appear in the supporting documentation. The authors make no claims *
10 * about the suitability of this software for any purpose. It is *
11 * provided "as is" without express or implied warranty. *
12 **************************************************************************/
16 //-----------------------------------------------------------------
17 // AliTagAnalysis class
18 // This is the class to deal with the tag analysis
19 // Origin: Panos Christakoglou, UOA-CERN, Panos.Christakoglou@cern.ch
20 //-----------------------------------------------------------------
26 #include <TEventList.h>
27 #include <TEntryList.h>
28 #include <TTreeFormula.h>
31 #include <TGridResult.h>
35 #include "AliRunTag.h"
36 #include "AliEventTag.h"
37 #include "AliTagAnalysis.h"
38 #include "AliEventTagCuts.h"
39 #include "AliDetectorTagCuts.h"
40 #include "AliLHCTagCuts.h"
41 #include "AliRunTagCuts.h"
42 #include "AliXMLCollection.h"
46 ClassImp(AliTagAnalysis)
48 //___________________________________________________________________________
49 AliTagAnalysis::AliTagAnalysis():
56 //Default constructor for a AliTagAnalysis
59 //___________________________________________________________________________
60 AliTagAnalysis::AliTagAnalysis(const char* type):
67 //constructor for a AliTagAnalysis
70 //___________________________________________________________________________
71 AliTagAnalysis::~AliTagAnalysis() {
72 //Default destructor for a AliTagAnalysis
75 //___________________________________________________________________________
76 Bool_t AliTagAnalysis::AddTagsFile(const char *alienUrl) {
78 // Add a single tags file to the chain
82 if (! fChain) fChain = new TChain("T");
84 TFile *f = TFile::Open(alienUrl,"READ");
85 fChain->Add(alienUrl);
86 AliInfo(Form("Chained tag files: %d ",fChain->GetEntries()));
89 if (fChain->GetEntries() == 0 )
95 //___________________________________________________________________________
96 void AliTagAnalysis::ChainLocalTags(const char *dirname) {
97 //Searches the entries of the provided direcory
98 //Chains the tags that are stored locally
99 fTagDirName = dirname;
100 TString fTagFilename;
102 if (! fChain) fChain = new TChain("T");
103 const char * tagPattern = 0x0;
104 if(fAnalysisType == "ESD") tagPattern = "ESD.tag.root";
105 else if(fAnalysisType == "AOD") tagPattern = "AOD.tag.root";
106 else AliFatal("Only ESD and AOD type is implemented!!!");
108 // Open the working directory
109 void * dirp = gSystem->OpenDirectory(fTagDirName);
110 const char * name = 0x0;
111 // Add all files matching *pattern* to the chain
112 while((name = gSystem->GetDirEntry(dirp))) {
113 if (strstr(name,tagPattern)) {
114 fTagFilename = fTagDirName;
116 fTagFilename += name;
118 fChain->Add(fTagFilename);
121 AliInfo(Form("Chained tag files: %d ",fChain->GetEntries()));
125 //___________________________________________________________________________
126 void AliTagAnalysis::ChainGridTags(TGridResult *res) {
127 //Loops overs the entries of the TGridResult
128 //Chains the tags that are stored in the GRID
130 Int_t nEntries = ftagresult->GetEntries();
132 if (! fChain) fChain = new TChain("T");
134 TString gridname = "alien://";
137 for(Int_t i = 0; i < nEntries; i++) {
138 alienUrl = ftagresult->GetKey(i,"turl");
139 fChain->Add(alienUrl);
144 //___________________________________________________________________________
145 TChain *AliTagAnalysis::QueryTags(AliRunTagCuts *runTagCuts, AliLHCTagCuts *lhcTagCuts, AliDetectorTagCuts *detTagCuts, AliEventTagCuts *evTagCuts) {
146 //Queries the tag chain using the defined
147 //event tag cuts from the AliEventTagCuts object
148 //and returns a TChain along with the associated TEventList
149 AliInfo(Form("Querying the tags........"));
152 if(fAnalysisType == "ESD") fAliceFile = "esdTree";
153 else if(fAnalysisType == "AOD") fAliceFile = "aodTree";
154 else AliFatal("Only ESD and AOD type is implemented!!!");
157 TChain *fESDchain = new TChain(fAliceFile.Data());
159 TEntryList *fGlobalList = new TEntryList();
161 //Defining tag objects
162 AliRunTag *tag = new AliRunTag;
163 AliEventTag *evTag = new AliEventTag;
164 fChain->SetBranchAddress("AliTAG",&tag);
171 for(Int_t iTagFiles = 0; iTagFiles < fChain->GetEntries(); iTagFiles++) {
172 fChain->GetEntry(iTagFiles);
173 if(runTagCuts->IsAccepted(tag)) {
174 if(lhcTagCuts->IsAccepted(tag->GetLHCTag())) {
175 if(detTagCuts->IsAccepted(tag->GetDetectorTags())) {
176 TEntryList *fLocalList = new TEntryList();
177 Int_t iEvents = tag->GetNEvents();
178 const TClonesArray *tagList = tag->GetEventTags();
179 for(Int_t i = 0; i < iEvents; i++) {
180 evTag = (AliEventTag *) tagList->At(i);
181 guid = evTag->GetGUID();
182 turl = evTag->GetTURL();
183 path = evTag->GetPath();
184 fLocalList->SetTreeName(fAliceFile.Data());
185 if(turl!="") fLocalList->SetFileName(turl.Data());
186 else fLocalList->SetFileName(path.Data());
187 if(evTagCuts->IsAccepted(evTag)) fLocalList->Enter(i);
190 if(path != "") fESDchain->AddFile(path);
191 else if(turl != "") fESDchain->AddFile(turl);
192 fGlobalList->Add(fLocalList);
193 iAccepted += fLocalList->GetN();
198 AliInfo(Form("Accepted events: %d",iAccepted));
199 fESDchain->SetEntryList(fGlobalList,"ne");
204 //___________________________________________________________________________
205 TChain *AliTagAnalysis::QueryTags(const char *fRunCut, const char *fLHCCut, const char *fDetectorCut, const char *fEventCut) {
206 //Queries the tag chain using the defined
207 //event tag cuts from the AliEventTagCuts object
208 //and returns a TChain along with the associated TEventList
209 AliInfo(Form("Querying the tags........"));
212 if(fAnalysisType == "ESD") fAliceFile = "esdTree";
213 else if(fAnalysisType == "AOD") fAliceFile = "aodTree";
214 else AliFatal("Only ESD and AOD type is implemented!!!");
217 TChain *fESDchain = new TChain(fAliceFile.Data());
219 TEntryList *fGlobalList = new TEntryList();
221 //Defining tag objects
222 AliRunTag *tag = new AliRunTag;
223 AliEventTag *evTag = new AliEventTag;
224 fChain->SetBranchAddress("AliTAG",&tag);
230 TTreeFormula *fRunFormula = new TTreeFormula("fRun",fRunCut,fChain);
231 TTreeFormula *fLHCFormula = new TTreeFormula("fLHC",fLHCCut,fChain);
232 TTreeFormula *fDetectorFormula = new TTreeFormula("fDetector",fDetectorCut,fChain);
233 TTreeFormula *fEventFormula = new TTreeFormula("fEvent",fEventCut,fChain);
237 for(Int_t iTagFiles = 0; iTagFiles < fChain->GetEntries(); iTagFiles++) {
238 fChain->GetEntry(iTagFiles);
239 if (current != fChain->GetTreeNumber()) {
240 fRunFormula->UpdateFormulaLeaves();
241 fLHCFormula->UpdateFormulaLeaves();
242 fDetectorFormula->UpdateFormulaLeaves();
243 fEventFormula->UpdateFormulaLeaves();
244 current = fChain->GetTreeNumber();
246 if(fRunFormula->EvalInstance(iTagFiles) == 1) {
247 if(fLHCFormula->EvalInstance(iTagFiles) == 1) {
248 if(fDetectorFormula->EvalInstance(iTagFiles) == 1) {
249 TEntryList *fLocalList = new TEntryList();
250 Int_t iEvents = fEventFormula->GetNdata();
251 const TClonesArray *tagList = tag->GetEventTags();
252 for(Int_t i = 0; i < iEvents; i++) {
253 evTag = (AliEventTag *) tagList->At(i);
254 guid = evTag->GetGUID();
255 turl = evTag->GetTURL();
256 path = evTag->GetPath();
257 fLocalList->SetTreeName(fAliceFile.Data());
258 fLocalList->SetFileName(turl.Data());
259 if(fEventFormula->EvalInstance(i) == 1) fLocalList->Enter(i);
261 iAccepted += fLocalList->GetN();
263 if(path != "") fESDchain->AddFile(path);
264 else if(turl != "") fESDchain->AddFile(turl);
265 fGlobalList->Add(fLocalList);
266 iAccepted += fLocalList->GetN();
271 AliInfo(Form("Accepted events: %d",iAccepted));
272 fESDchain->SetEntryList(fGlobalList,"ne");
277 //___________________________________________________________________________
278 Bool_t AliTagAnalysis::CreateXMLCollection(const char* name, AliRunTagCuts *runTagCuts, AliLHCTagCuts *lhcTagCuts, AliDetectorTagCuts *detTagCuts, AliEventTagCuts *evTagCuts) {
279 //Queries the tag chain using the defined
280 //event tag cuts from the AliEventTagCuts object
281 //and returns a XML collection
282 AliInfo(Form("Creating the collection........"));
284 AliXMLCollection *collection = new AliXMLCollection();
285 collection->SetCollectionName(name);
286 collection->WriteHeader();
292 //Defining tag objects
293 AliRunTag *tag = new AliRunTag;
294 AliEventTag *evTag = new AliEventTag;
295 fChain->SetBranchAddress("AliTAG",&tag);
297 for(Int_t iTagFiles = 0; iTagFiles < fChain->GetEntries(); iTagFiles++) {
299 TEntryList *fList = new TEntryList();
300 fChain->GetEntry(iTagFiles);
301 if(runTagCuts->IsAccepted(tag)) {
302 if(lhcTagCuts->IsAccepted(tag->GetLHCTag())) {
303 if(detTagCuts->IsAccepted(tag->GetDetectorTags())) {
304 Int_t iEvents = tag->GetNEvents();
305 const TClonesArray *tagList = tag->GetEventTags();
306 for(Int_t i = 0; i < iEvents; i++) {
307 evTag = (AliEventTag *) tagList->At(i);
308 guid = evTag->GetGUID();
309 turl = evTag->GetTURL();
310 lfn = turl(8,turl.Length());
311 if(evTagCuts->IsAccepted(evTag)) fList->Enter(i);
313 collection->WriteBody(iTagFiles+1,guid,lfn,turl,fList);
318 collection->Export();
323 //___________________________________________________________________________
324 Bool_t AliTagAnalysis::CreateXMLCollection(const char* name, const char *fRunCut, const char *fLHCCut, const char *fDetectorCut, const char *fEventCut) {
325 //Queries the tag chain using the defined
326 //event tag cuts from the AliEventTagCuts object
327 //and returns a XML collection
328 AliInfo(Form("Creating the collection........"));
330 AliXMLCollection *collection = new AliXMLCollection();
331 collection->SetCollectionName(name);
332 collection->WriteHeader();
338 //Defining tag objects
339 AliRunTag *tag = new AliRunTag;
340 AliEventTag *evTag = new AliEventTag;
341 fChain->SetBranchAddress("AliTAG",&tag);
343 TTreeFormula *fRunFormula = new TTreeFormula("fRun",fRunCut,fChain);
344 TTreeFormula *fLHCFormula = new TTreeFormula("fLHC",fLHCCut,fChain);
345 TTreeFormula *fDetectorFormula = new TTreeFormula("fDetector",fDetectorCut,fChain);
346 TTreeFormula *fEventFormula = new TTreeFormula("fEvent",fEventCut,fChain);
349 for(Int_t iTagFiles = 0; iTagFiles < fChain->GetEntries(); iTagFiles++) {
351 TEntryList *fList = new TEntryList();
352 fChain->GetEntry(iTagFiles);
353 if (current != fChain->GetTreeNumber()) {
354 fRunFormula->UpdateFormulaLeaves();
355 fLHCFormula->UpdateFormulaLeaves();
356 fDetectorFormula->UpdateFormulaLeaves();
357 fEventFormula->UpdateFormulaLeaves();
358 current = fChain->GetTreeNumber();
360 if(fRunFormula->EvalInstance(iTagFiles) == 1) {
361 if(fLHCFormula->EvalInstance(iTagFiles) == 1) {
362 if(fDetectorFormula->EvalInstance(iTagFiles) == 1) {
363 Int_t iEvents = fEventFormula->GetNdata();
364 const TClonesArray *tagList = tag->GetEventTags();
365 for(Int_t i = 0; i < iEvents; i++) {
366 evTag = (AliEventTag *) tagList->At(i);
367 guid = evTag->GetGUID();
368 turl = evTag->GetTURL();
369 lfn = turl(8,turl.Length());
370 if(fEventFormula->EvalInstance(i) == 1) fList->Enter(i);
372 collection->WriteBody(iTagFiles+1,guid,lfn,turl,fList);
377 collection->Export();
382 //___________________________________________________________________________
383 TChain *AliTagAnalysis::GetInputChain(const char* system, const char *wn) {
384 //returns the chain+event list - used in batch sessions
385 // this function will be removed once the new root
386 // improvements are committed
387 TString fsystem = system;
390 TChain *fAnalysisChain = 0;
391 if(fAnalysisType == "ESD") fAnalysisChain = new TChain("esdTree");
392 else if(fAnalysisType == "AOD") fAnalysisChain = new TChain("aodTree");
393 else AliFatal("Only ESD and AOD type is implemented!!!");
396 TEventList *fEventList = new TEventList();
397 AliXMLCollection *collection = AliXMLCollection::Open(wn);
400 while (collection->Next()) {
401 AliInfo(Form("Adding: %s",collection->GetTURL("")));
402 fAnalysisChain->Add(collection->GetTURL(""));
403 TEntryList *list = (TEntryList *)collection->GetEventList("");
404 for(Int_t i = 0; i < list->GetN(); i++) fEventList->Enter(iAccepted+list->GetEntry(i));
406 if(fsystem == "pp") iAccepted += 100;
407 else if(fsystem == "PbPb") iAccepted += 1;
410 fAnalysisChain->SetEventList(fEventList);
412 AliInfo(Form("Number of selected events: %d",fEventList->GetN()));
414 return fAnalysisChain;
417 //___________________________________________________________________________
418 TChain *AliTagAnalysis::GetChainFromCollection(const char* collectionname, const char* treename) {
419 //returns the TChain+TEntryList object- used in batch sessions
420 TString fAliceFile = treename;
422 TChain *fAnalysisChain = 0;
423 if(fAliceFile == "esdTree") fAnalysisChain = new TChain("esdTree");
424 else if(fAliceFile == "aodTree") fAnalysisChain = new TChain("aodTree");
425 else AliFatal("Inconsistent tree name - use esdTree or aodTree!");
428 TEntryList *fGlobalList = new TEntryList();
429 AliXMLCollection *collection = AliXMLCollection::Open(collectionname);
432 while (collection->Next()) {
433 AliInfo(Form("Adding: %s",collection->GetTURL("")));
434 fAnalysisChain->Add(collection->GetTURL(""));
435 TEntryList *list = (TEntryList *)collection->GetEventList("");
436 list->SetTreeName(fAliceFile.Data());
437 list->SetFileName(collection->GetTURL(""));
438 fGlobalList->Add(list);
439 iAccepted += list->GetN();
442 fAnalysisChain->SetEntryList(fGlobalList,"ne");
444 AliInfo(Form("Number of selected events: %d",iAccepted));
446 return fAnalysisChain;