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 //-----------------------------------------------------------------
23 #include <Riostream.h>
27 #include <TEventList.h>
28 #include <TEntryList.h>
29 #include <TTreeFormula.h>
33 #include <TGridResult.h>
37 #include "AliRunTag.h"
38 #include "AliEventTag.h"
39 #include "AliTagAnalysis.h"
40 #include "AliEventTagCuts.h"
41 #include "AliDetectorTagCuts.h"
42 #include "AliLHCTagCuts.h"
43 #include "AliRunTagCuts.h"
44 #include "AliXMLCollection.h"
48 ClassImp(AliTagAnalysis)
50 //___________________________________________________________________________
51 AliTagAnalysis::AliTagAnalysis():
58 //Default constructor for a AliTagAnalysis
61 //___________________________________________________________________________
62 AliTagAnalysis::AliTagAnalysis(const char* type):
69 //constructor for a AliTagAnalysis
72 //___________________________________________________________________________
73 AliTagAnalysis::~AliTagAnalysis() {
74 //Default destructor for a AliTagAnalysis
75 if(ftagresult) delete ftagresult;
76 if(fChain) delete fChain;
77 if(fGlobalList) delete fGlobalList;
80 //___________________________________________________________________________
82 AliTagAnalysis::AddTagsFile(const char* alienUrl, Bool_t checkFile)
84 /// Add a single tags file to the chain
86 /// If checkFile=kTRUE (default) the file is opened to check
87 /// it can be and that it contains data.
88 /// It's safer but a lot longer...
90 if (!fChain) fChain = new TChain("T");
94 return ( fChain->AddFile(alienUrl,-1) > 0 );
98 return ( fChain->AddFile(alienUrl) > 0 );
103 //___________________________________________________________________________
104 void AliTagAnalysis::ChainLocalTags(const char *dirname) {
105 //Searches the entries of the provided direcory
106 //Chains the tags that are stored locally
107 fTagDirName = dirname;
108 TString fTagFilename;
110 if (! fChain) fChain = new TChain("T");
111 const char * tagPattern = 0x0;
112 if(fAnalysisType == "ESD") tagPattern = "ESD.tag.root";
113 else if(fAnalysisType == "AOD") tagPattern = "AOD.tag.root";
114 else AliFatal("Only ESD and AOD type is implemented!!!");
116 // Open the working directory
117 void * dirp = gSystem->OpenDirectory(fTagDirName);
118 const char * name = 0x0;
119 // Add all files matching *pattern* to the chain
120 while((name = gSystem->GetDirEntry(dirp))) {
121 if (strstr(name,tagPattern)) {
122 fTagFilename = fTagDirName;
124 fTagFilename += name;
126 fChain->Add(fTagFilename);
127 printf("Tag file %s\n", fTagFilename.Data());
131 //AliInfo(Form("Chained tag files: %d ",fChain->GetEntries()));
132 // AliDebug(Form("Chained tag files: %d ",fChain->GetEntries()));
138 //___________________________________________________________________________
139 TChain * AliTagAnalysis::ChainGridTags(TGridResult *res) {
140 //Loops overs the entries of the TGridResult
141 //Chains the tags that are stored in the GRID
143 Int_t nEntries = ftagresult->GetEntries();
145 if (! fChain) fChain = new TChain("T");
147 TString gridname = "alien://";
150 for(Int_t i = 0; i < nEntries; i++) {
151 alienUrl = ftagresult->GetKey(i,"turl");
152 fChain->Add(alienUrl);
158 //___________________________________________________________________________
159 TChain *AliTagAnalysis::QueryTags(AliRunTagCuts *runTagCuts,
160 AliLHCTagCuts *lhcTagCuts,
161 AliDetectorTagCuts *detTagCuts,
162 AliEventTagCuts *evTagCuts) {
163 //Queries the tag chain using the defined
164 //event tag cuts from the AliEventTagCuts object
165 //and returns a TChain along with the associated TEventList
166 AliInfo(Form("Querying the tags........"));
169 if(fAnalysisType == "ESD") aliceFile = "esdTree";
170 else if(fAnalysisType == "AOD") aliceFile = "aodTree";
171 else AliFatal("Only ESD and AOD type is implemented!!!");
174 TChain *esdChain = new TChain(aliceFile.Data());
176 fGlobalList = new TEntryList();
178 //Defining tag objects
179 AliRunTag *tag = new AliRunTag;
180 // AliEventTag *evTag = 0x0;
181 AliFileTag *flTag = 0x0;
183 fChain->SetBranchAddress("AliTAG",&tag);
189 TEntryList* localList = new TEntryList();
193 for(Int_t iEntry = 0; iEntry < fChain->GetEntries(); iEntry++) {
194 fChain->GetEntry(iEntry);
195 evTagCuts->InitializeTriggerClasses(tag->GetActiveTriggerClasses());
197 if(runTagCuts->IsAccepted(tag)) {
198 if(lhcTagCuts->IsAccepted(tag->GetLHCTag())) {
199 if(detTagCuts->IsAccepted(tag->GetDetectorTags())) {
201 Int_t iEvents = tag->GetNEvents();
203 for (int i = 0; i < iEvents; i++) {
204 // evTag = tag->GetEventTag(i);
205 flTag = tag->GetFileTagForEvent(i);
206 guid = flTag->GetGUID();
207 turl = flTag->GetTURL();
208 path = flTag->GetPath();
209 localList->SetTreeName(aliceFile.Data());
210 if(turl!="") localList->SetFileName(turl.Data());
211 else localList->SetFileName(path.Data());
213 if(evTagCuts->IsAccepted(tag->GetEventTag(i))) localList->Enter(i);
216 // const TClonesArray *tagList = tag->GetEventTags();
217 // for(Int_t i = 0; i < iEvents; i++) {
218 // evTag = (AliEventTag *) tagList->At(i);
219 // guid = evTag->GetGUID();
220 // turl = evTag->GetTURL();
221 // path = evTag->GetPath();
222 // localList->SetTreeName(aliceFile.Data());
223 // if(turl!="") localList->SetFileName(turl.Data());
224 // else localList->SetFileName(path.Data());
226 // if(evTagCuts->IsAccepted(evTag)) localList->Enter(i);
228 iAccepted += localList->GetN();
229 if(turl != "") esdChain->AddFile(turl);
230 else if(path != "") esdChain->AddFile(path);
231 fGlobalList->Add(localList);
237 AliInfo(Form("Accepted events: %d", iAccepted));
239 esdChain->SetEntryList(fGlobalList,"ne");
246 //___________________________________________________________________________
247 TChain *AliTagAnalysis::QueryTags(const char *fRunCut,
249 const char *fDetectorCut,
250 const char *fEventCut) {
251 //Queries the tag chain using the defined
252 //event tag cuts from the AliEventTagCuts object
253 //and returns a TChain along with the associated TEventList
254 AliInfo(Form("Querying the tags........"));
257 if(fAnalysisType == "ESD") aliceFile = "esdTree";
258 else if(fAnalysisType == "AOD") aliceFile = "aodTree";
259 else AliFatal("Only ESD and AOD type is implemented!!!");
263 TChain *esdChain = new TChain(aliceFile.Data());
265 fGlobalList = new TEntryList();
267 //Defining tag objects
268 AliRunTag *tag = new AliRunTag;
269 // AliEventTag *evTag = 0x0;
270 fChain->SetBranchAddress("AliTAG",&tag);
276 TTreeFormula *fRunFormula = new TTreeFormula("fRun",fRunCut,fChain);
277 TTreeFormula *fLHCFormula = new TTreeFormula("fLHC",fLHCCut,fChain);
278 TTreeFormula *fDetectorFormula = new TTreeFormula("fDetector",fDetectorCut,fChain);
279 TTreeFormula *fEventFormula = new TTreeFormula("fEvent",fEventCut,fChain);
281 TEntryList* localList = new TEntryList();
286 for(Int_t iTagFiles = 0; iTagFiles < fChain->GetEntries(); iTagFiles++) {
287 fChain->GetEntry(iTagFiles);
288 if (current != fChain->GetTreeNumber()) {
289 fRunFormula->UpdateFormulaLeaves();
290 fLHCFormula->UpdateFormulaLeaves();
291 fDetectorFormula->UpdateFormulaLeaves();
292 fEventFormula->UpdateFormulaLeaves();
293 current = fChain->GetTreeNumber();
296 if(fRunFormula->EvalInstance(iTagFiles) == 1) {
297 if(fLHCFormula->EvalInstance(iTagFiles) == 1) {
298 if(fDetectorFormula->EvalInstance(iTagFiles) == 1) {
300 // Int_t iEvents = fEventFormula->GetNdata();
303 // const TClonesArray *tagList = tag->GetEventTags();
304 // for(Int_t i = 0; i < iEvents; i++) {
305 // evTag = (AliEventTag *) tagList->At(i);
306 // guid = evTag->GetGUID();
307 // turl = evTag->GetTURL();
308 // path = evTag->GetPath();
309 // localList->SetTreeName(aliceFile.Data());
310 // localList->SetFileName(turl.Data());
311 // if(fEventFormula->EvalInstance(i) == 1) localList->Enter(i);
314 if(path != "") esdChain->AddFile(path);
315 else if(turl != "") esdChain->AddFile(turl);
316 fGlobalList->Add(localList);
317 iAccepted += localList->GetN();
323 AliInfo(Form("Accepted events: %d", iAccepted));
324 esdChain->SetEntryList(fGlobalList,"ne");
331 //___________________________________________________________________________
333 AliTagAnalysis::CreateXMLCollection(const char* name,
334 AliRunTagCuts *runTagCuts,
335 AliLHCTagCuts *lhcTagCuts,
336 AliDetectorTagCuts *detTagCuts,
337 AliEventTagCuts *evTagCuts)
339 /// Queries the tag chain using the defined run, lhc, detector and event tag objects
340 /// and create a XML collection named "name.xml"
341 /// if any of the runTagCuts, lhcTagCuts, detTagCuts or evTagCuts is NULL
342 /// check on that object will be skipped.
344 AliInfo(Form("Creating the collection........"));
348 AliError("fChain is NULL. Cannot make a collection from that !");
352 AliXMLCollection collection;
353 collection.SetCollectionName(name);
354 collection.WriteHeader();
360 TEntryList localList;
363 Int_t iRejectedRun = 0;
364 Int_t iRejectedLHC = 0;
365 Int_t iRejectedDet = 0;
366 Int_t iRejectedEvt = 0;
368 Int_t iTotalEvents = 0;
370 Int_t iAcceptedEvtInFile = 0;
371 Int_t iRejectedEvtInFile = 0;
373 //Defining tag objects
374 AliRunTag* tag = new AliRunTag;
375 fChain->SetBranchAddress("AliTAG",&tag);
379 // AliEventTag *evTag = 0x0;
380 AliFileTag *flTag = 0x0;
382 // for(Int_t iTagFiles = 0; iTagFiles < fChain->GetListOfFiles()->GetEntries(); ++iTagFiles)
383 for(Int_t iRunTags = 0; iRunTags < fChain->GetEntries(); ++iRunTags)
385 fChain->GetEntry(iRunTags);
387 iTotalEvents += tag->GetNEvents();
390 evTagCuts->InitializeTriggerClasses(tag->GetActiveTriggerClasses());
392 if ( !runTagCuts || ( runTagCuts && runTagCuts->IsAccepted(tag) ) )
394 if ( !lhcTagCuts || ( lhcTagCuts && lhcTagCuts->IsAccepted(tag->GetLHCTag())) )
396 if ( !detTagCuts || ( detTagCuts && detTagCuts->IsAccepted(tag->GetDetectorTags())) )
398 for (int iChunk = 0; iChunk < tag->GetNFiles(); iChunk++, iTagFiles++)
400 iRejectedEvtInFile = 0;
401 iAcceptedEvtInFile = 0;
405 flTag = tag->GetFileTag(iChunk);
406 guid = flTag->GetGUID();
407 turl = flTag->GetTURL();
408 lfn = turl(8,turl.Length());
410 for (int i = 0; i<flTag->GetNEvents(); i++)
412 // evTag = flTag->GetEventTag(i);
414 if( !evTagCuts || ( evTagCuts && evTagCuts->IsAccepted(flTag->GetEventTag(i))) )
417 iAcceptedEvtInFile++;
422 ++iRejectedEvtInFile;
428 // TIter next(tag->GetEventTags());
429 // AliEventTag* evTag(0x0);
430 // iRejectedEvtInFile = 0;
431 // iAcceptedEvtInFile = 0;
432 // while ( ( evTag = static_cast<AliEventTag*>(next()) ) )
434 // guid = evTag->GetGUID();
435 // turl = evTag->GetTURL();
436 // lfn = turl(8,turl.Length());
437 // if( !evTagCuts || ( evTagCuts && evTagCuts->IsAccepted(evTag)) )
439 // localList.Enter(i);
440 // iAcceptedEvtInFile++;
445 // ++iRejectedEvtInFile;
449 iAccepted += localList.GetN();
450 collection.WriteBody(iTagFiles+1,guid,lfn,turl,&localList,iAcceptedEvtInFile,iRejectedEvtInFile);
454 iRejectedDet += tag->GetNEvents();
458 iRejectedLHC += tag->GetNEvents();
462 iRejectedRun += tag->GetNEvents();
467 collection.WriteSummary(iTotalEvents, iAccepted, iRejectedRun, iRejectedLHC, iRejectedDet, iRejectedEvt);
473 //___________________________________________________________________________
474 Bool_t AliTagAnalysis::CreateXMLCollection(const char* name,
477 const char *fDetectorCut,
478 const char *fEventCut) {
479 //Queries the tag chain using the defined
480 //event tag cuts from the AliEventTagCuts object
481 //and returns a XML collection
482 AliInfo(Form("Creating the collection........"));
485 AliXMLCollection *collection = new AliXMLCollection();
486 collection->SetCollectionName(name);
487 collection->WriteHeader();
492 TEntryList* localList = new TEntryList();
496 Int_t iRejectedRun = 0;
497 Int_t iRejectedLHC = 0;
498 Int_t iRejectedDet = 0;
499 Int_t iRejectedEvt = 0;
501 Int_t iTotalEvents = 0;
503 Int_t iAcceptedEvtInFile = 0;
504 Int_t iRejectedEvtInFile = 0;
506 //Defining tag objects
507 AliRunTag *tag = new AliRunTag;
508 // AliEventTag *evTag = 0x0;
509 fChain->SetBranchAddress("AliTAG",&tag);
511 TTreeFormula *fRunFormula = new TTreeFormula("fRun",fRunCut,fChain);
512 TTreeFormula *fLHCFormula = new TTreeFormula("fLHC",fLHCCut,fChain);
513 TTreeFormula *fDetectorFormula = new TTreeFormula("fDetector",fDetectorCut,fChain);
514 TTreeFormula *fEventFormula = new TTreeFormula("fEvent",fEventCut,fChain);
518 for(Int_t iTagFiles = 0; iTagFiles < fChain->GetEntries(); iTagFiles++) {
520 fChain->GetEntry(iTagFiles);
521 if (current != fChain->GetTreeNumber()) {
522 fRunFormula->UpdateFormulaLeaves();
523 fLHCFormula->UpdateFormulaLeaves();
524 fDetectorFormula->UpdateFormulaLeaves();
525 fEventFormula->UpdateFormulaLeaves();
526 current = fChain->GetTreeNumber();
530 iTotalEvents += tag->GetNEvents();
532 if(fRunFormula->EvalInstance(iTagFiles) == 1) {
533 if(fLHCFormula->EvalInstance(iTagFiles) == 1) {
534 if(fDetectorFormula->EvalInstance(iTagFiles) == 1) {
535 // Int_t iEvents = fEventFormula->GetNdata();
539 // const TClonesArray *tagList = tag->GetEventTags();
540 // iRejectedEvtInFile = 0;
541 // iAcceptedEvtInFile = 0;
542 // for(Int_t i = 0; i < iEvents; i++) {
543 // evTag = (AliEventTag *) tagList->At(i);
544 // guid = evTag->GetGUID();
545 // turl = evTag->GetTURL();
546 // lfn = turl(8,turl.Length());
547 // if(fEventFormula->EvalInstance(i) == 1) {
548 // localList->Enter(i);
549 // iAcceptedEvtInFile++;
553 // iRejectedEvtInFile++;
557 collection->WriteBody(iTagFiles+1,guid,lfn,turl,localList,iAcceptedEvtInFile, iRejectedEvtInFile);
558 iAccepted += localList->GetN();
561 iRejectedDet += tag->GetNEvents();
565 iRejectedLHC += tag->GetNEvents();
569 iRejectedRun += tag->GetNEvents();
572 collection->WriteSummary(iTotalEvents, iAccepted, iRejectedRun, iRejectedLHC, iRejectedDet, iRejectedEvt);
573 collection->Export();
579 //___________________________________________________________________________
580 TChain *AliTagAnalysis::GetInputChain(const char* system, const char *wn) {
581 //returns the chain+event list - used in batch sessions
582 // this function will be removed once the new root
583 // improvements are committed
584 TString fsystem = system;
587 TChain *fAnalysisChain = 0;
588 if(fAnalysisType == "ESD") fAnalysisChain = new TChain("esdTree");
589 else if(fAnalysisType == "AOD") fAnalysisChain = new TChain("aodTree");
590 else AliFatal("Only ESD and AOD type is implemented!!!");
593 TEventList *fEventList = new TEventList();
594 AliXMLCollection *collection = AliXMLCollection::Open(wn);
597 while (collection->Next()) {
598 AliInfo(Form("Adding: %s",collection->GetTURL("")));
599 fAnalysisChain->Add(collection->GetTURL(""));
600 TEntryList *list = (TEntryList *)collection->GetEventList("");
601 for(Int_t i = 0; i < list->GetN(); i++) fEventList->Enter(iAccepted+list->GetEntry(i));
603 if(fsystem == "pp") iAccepted += 100;
604 else if(fsystem == "PbPb") iAccepted += 1;
607 fAnalysisChain->SetEventList(fEventList);
609 AliInfo(Form("Number of selected events: %d",fEventList->GetN()));
611 return fAnalysisChain;
614 //___________________________________________________________________________
616 AliTagAnalysis::CreateChainFromCollection(const char* collectionname, const char* treename)
618 /// Build a TChain (with its TEntryList object attached) from an XML collection.
619 /// Returned chain must be deleted by the client.
621 TString streename(treename);
622 if ( streename != "esdTree" && streename != "aodTree" )
624 AliErrorClass("Only esdTree and aodTree implemented so far...");
628 TChain* chain = new TChain(streename.Data());
630 // create the event list for the chain. Will be attached to the chain
631 // which thus becomes the owner of it.
632 TEntryList* elist = new TEntryList;
634 AliXMLCollection* collection = AliXMLCollection::Open(collectionname);
636 // Tag selection summary per file
637 TMap* tagCutSummary = new TMap();
638 tagCutSummary->SetName("TagCutSumm");
644 while (collection->Next())
646 AliDebugClass(1,Form("Adding: %s",collection->GetTURL("")));
647 chain->Add(collection->GetTURL(""));
648 TEntryList *list = collection->GetEventList("");
649 list->SetTreeName(streename.Data());
650 list->SetFileName(collection->GetTURL(""));
652 iAccepted += list->GetN();
653 if (collection->GetCutSumm())
655 tagCutSummary->Add(new TObjString(collection->GetTURL("")), new TObjString(collection->GetCutSumm()));
659 chain->SetEntryList(elist,"ne"); // ne => do not expand tree name and/or file names
661 AliDebugClass(1,Form("Number of selected events: %d",iAccepted));
663 TList *aUserInfo = chain->GetUserInfo();
664 aUserInfo->Add(tagCutSummary);
673 collection->GetCollectionSummary(&iTotalEvents, &iAccEv, &iRejRun, &iRejLHC, &iRejDet, &iRejEvt);
677 sprintf(nstr, "TotalEvents=%i", iTotalEvents);
678 TObjString *iTotStr = new TObjString(nstr);
679 aUserInfo->Add(iTotStr);
681 sprintf(nstr, "AcceptedEvents=%i", iAccepted);
682 TObjString *iAccStr = new TObjString(nstr);
683 aUserInfo->Add(iAccStr);
685 sprintf(nstr, "RejectedRun=%i", iRejRun);
686 TObjString *iRejRunStr = new TObjString(nstr);
687 aUserInfo->Add(iRejRunStr);
689 sprintf(nstr, "RejectedLHC=%i", iRejLHC);
690 TObjString *iRejLHCStr = new TObjString(nstr);
691 aUserInfo->Add(iRejLHCStr);
693 sprintf(nstr, "RejectedDet=%i", iRejDet);
694 TObjString *iRejDetStr = new TObjString(nstr);
695 aUserInfo->Add(iRejDetStr);
697 sprintf(nstr, "RejectedEvt=%i", iRejEvt);
698 TObjString *iRejEvtStr = new TObjString(nstr);
699 aUserInfo->Add(iRejEvtStr);