Introduced tree caching and async reading for data (ESD and AOD) and MC. An read...
[u/mrichter/AliRoot.git] / ANALYSIS / AliAnalysisManager.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
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  **************************************************************************/
15
16 /* $Id$ */
17 // Author: Andrei Gheata, 31/05/2006
18
19 //==============================================================================
20 //   AliAnalysisManager - Manager analysis class. Allows creation of several
21 // analysis tasks and data containers storing their input/output. Allows
22 // connecting/chaining tasks via shared data containers. Serializes the current
23 // event for all tasks depending only on initial input data.
24 //==============================================================================
25 //
26 //==============================================================================
27
28 #include "AliAnalysisManager.h"
29
30 #include <cerrno>
31 #include <Riostream.h>
32 #include <TError.h>
33 #include <TMap.h>
34 #include <TClass.h>
35 #include <TFile.h>
36 #include <TTreeCache.h>
37 #include <TEnv.h>
38 #include <TMath.h>
39 #include <TH1.h>
40 #include <TMethodCall.h>
41 #include <TChain.h>
42 #include <TSystem.h>
43 #include <TROOT.h>
44 #include <TCanvas.h>
45 #include <TStopwatch.h>
46
47 #include "AliLog.h"
48 #include "AliAnalysisSelector.h"
49 #include "AliAnalysisGrid.h"
50 #include "AliAnalysisTask.h"
51 #include "AliAnalysisDataContainer.h"
52 #include "AliAnalysisDataSlot.h"
53 #include "AliVEventHandler.h"
54 #include "AliVEventPool.h"
55 #include "AliSysInfo.h"
56 #include "AliAnalysisStatistics.h"
57
58 using std::ofstream;
59 using std::ios;
60 using std::cout;
61 using std::endl;
62 ClassImp(AliAnalysisManager)
63
64 AliAnalysisManager *AliAnalysisManager::fgAnalysisManager = NULL;
65 TString AliAnalysisManager::fgCommonFileName = "";
66 Int_t AliAnalysisManager::fPBUpdateFreq = 1;
67
68 //______________________________________________________________________________
69 AliAnalysisManager::AliAnalysisManager(const char *name, const char *title)
70                    :TNamed(name,title),
71                     fTree(0),
72                     fInputEventHandler(0),
73                     fOutputEventHandler(0),
74                     fMCtruthEventHandler(0),
75                     fEventPool(0),
76                     fCurrentEntry(-1),
77                     fNSysInfo(0),
78                     fMode(kLocalAnalysis),
79                     fInitOK(kFALSE),
80                     fMustClean(kFALSE),
81                     fIsRemote(kFALSE),
82                     fLocked(kFALSE),
83                     fDebug(0),
84                     fSpecialOutputLocation(""), 
85                     fTasks(0),
86                     fTopTasks(0),
87                     fZombies(0),
88                     fContainers(0),
89                     fInputs(0),
90                     fOutputs(0),
91                     fParamCont(0),
92                     fDebugOptions(0),
93                     fFileDescriptors(new TObjArray()),
94                     fCurrentDescriptor(0),
95                     fCommonInput(0),
96                     fCommonOutput(0),
97                     fSelector(0),
98                     fGridHandler(0),
99                     fExtraFiles(),
100                     fFileInfoLog(),
101                     fAutoBranchHandling(kTRUE), 
102                     fAsyncReading(kTRUE), // default prefetching on
103                     fTable(),
104                     fRunFromPath(0),
105                     fNcalls(0),
106                     fMaxEntries(0),
107                     fCacheSize(100000000), // default 100 MB
108                     fStatisticsMsg(),
109                     fRequestedBranches(),
110                     fStatistics(0),
111                     fGlobals(0),
112                     fIOTimer(0),
113                     fCPUTimer(0),
114                     fInitTimer(0),
115                     fIOTime(0),
116                     fCPUTime(0),
117                     fInitTime(0)
118 {
119 // Default constructor.
120    fgAnalysisManager = this;
121    fgCommonFileName  = "AnalysisResults.root";
122    if (TClass::IsCallingNew() != TClass::kDummyNew) {
123      fTasks      = new TObjArray();
124      fTopTasks   = new TObjArray();
125      fZombies    = new TObjArray();
126      fContainers = new TObjArray();
127      fInputs     = new TObjArray();
128      fOutputs    = new TObjArray();
129      fParamCont  = new TObjArray();
130      fGlobals    = new TMap();
131    }
132    fIOTimer = new TStopwatch();
133    fCPUTimer = new TStopwatch();
134    fInitTimer = new TStopwatch();
135    SetEventLoop(kTRUE);
136 }
137
138 //______________________________________________________________________________
139 AliAnalysisManager::AliAnalysisManager(const AliAnalysisManager& other)
140                    :TNamed(other),
141                     fTree(NULL),
142                     fInputEventHandler(NULL),
143                     fOutputEventHandler(NULL),
144                     fMCtruthEventHandler(NULL),
145                     fEventPool(NULL),
146                     fCurrentEntry(-1),
147                     fNSysInfo(0),
148                     fMode(other.fMode),
149                     fInitOK(other.fInitOK),
150                     fMustClean(other.fMustClean),
151                     fIsRemote(other.fIsRemote),
152                     fLocked(other.fLocked),
153                     fDebug(other.fDebug),
154                     fSpecialOutputLocation(""), 
155                     fTasks(NULL),
156                     fTopTasks(NULL),
157                     fZombies(NULL),
158                     fContainers(NULL),
159                     fInputs(NULL),
160                     fOutputs(NULL),
161                     fParamCont(NULL),
162                     fDebugOptions(NULL),
163                     fFileDescriptors(new TObjArray()),
164                     fCurrentDescriptor(0),
165                     fCommonInput(NULL),
166                     fCommonOutput(NULL),
167                     fSelector(NULL),
168                     fGridHandler(NULL),
169                     fExtraFiles(other.fExtraFiles),
170                     fFileInfoLog(other.fFileInfoLog),
171                     fAutoBranchHandling(other.fAutoBranchHandling), 
172                     fAsyncReading(other.fAsyncReading),
173                     fTable(),
174                     fRunFromPath(0),
175                     fNcalls(other.fNcalls),
176                     fMaxEntries(other.fMaxEntries),
177                     fCacheSize(other.fCacheSize),
178                     fStatisticsMsg(other.fStatisticsMsg),
179                     fRequestedBranches(other.fRequestedBranches),
180                     fStatistics(other.fStatistics),
181                     fGlobals(other.fGlobals),
182                     fIOTimer(new TStopwatch()),
183                     fCPUTimer(new TStopwatch()),
184                     fInitTimer(new TStopwatch()),
185                     fIOTime(0),
186                     fCPUTime(0),
187                     fInitTime(0)
188 {
189 // Copy constructor.
190    fTasks      = new TObjArray(*other.fTasks);
191    fTopTasks   = new TObjArray(*other.fTopTasks);
192    fZombies    = new TObjArray(*other.fZombies);
193    fContainers = new TObjArray(*other.fContainers);
194    fInputs     = new TObjArray(*other.fInputs);
195    fOutputs    = new TObjArray(*other.fOutputs);
196    fParamCont  = new TObjArray(*other.fParamCont);
197    fgCommonFileName  = "AnalysisResults.root";
198    fgAnalysisManager = this;
199 }
200    
201 //______________________________________________________________________________
202 AliAnalysisManager& AliAnalysisManager::operator=(const AliAnalysisManager& other)
203 {
204 // Assignment
205    if (&other != this) {
206       TNamed::operator=(other);
207       fInputEventHandler   = other.fInputEventHandler;
208       fOutputEventHandler  = other.fOutputEventHandler;
209       fMCtruthEventHandler = other.fMCtruthEventHandler;
210       fEventPool           = other.fEventPool;
211       fTree       = NULL;
212       fCurrentEntry = -1;
213       fNSysInfo   = other.fNSysInfo;
214       fMode       = other.fMode;
215       fInitOK     = other.fInitOK;
216       fIsRemote   = other.fIsRemote;
217       fLocked     = other.fLocked;
218       fDebug      = other.fDebug;
219       fTasks      = new TObjArray(*other.fTasks);
220       fTopTasks   = new TObjArray(*other.fTopTasks);
221       fZombies    = new TObjArray(*other.fZombies);
222       fContainers = new TObjArray(*other.fContainers);
223       fInputs     = new TObjArray(*other.fInputs);
224       fOutputs    = new TObjArray(*other.fOutputs);
225       fParamCont  = new TObjArray(*other.fParamCont);
226       fDebugOptions = NULL;
227       fFileDescriptors = new TObjArray();
228       fCurrentDescriptor = 0;
229       fCommonInput = NULL;
230       fCommonOutput = NULL;
231       fSelector   = NULL;
232       fGridHandler = NULL;
233       fExtraFiles = other.fExtraFiles;
234       fFileInfoLog = other.fFileInfoLog;
235       fgCommonFileName = "AnalysisResults.root";
236       fgAnalysisManager = this;
237       fAutoBranchHandling = other.fAutoBranchHandling;
238       fAsyncReading = other.fAsyncReading;
239       fTable.Clear("nodelete");
240       fRunFromPath = other.fRunFromPath;
241       fNcalls     = other. fNcalls;
242       fMaxEntries = other.fMaxEntries;
243       fCacheSize = other.fCacheSize;
244       fStatisticsMsg = other.fStatisticsMsg;
245       fRequestedBranches = other.fRequestedBranches;
246       fStatistics = other.fStatistics;
247       fGlobals = new TMap();
248       fIOTimer = new TStopwatch();
249       fCPUTimer = new TStopwatch();
250       fInitTimer = new TStopwatch();
251       fIOTime = 0.;
252       fCPUTime = 0.;
253       fInitTime = 0.;
254    }
255    return *this;
256 }
257
258 //______________________________________________________________________________
259 AliAnalysisManager::~AliAnalysisManager()
260 {
261 // Destructor.
262    if (fTasks) {fTasks->Delete(); delete fTasks;}
263    if (fTopTasks) delete fTopTasks;
264    if (fZombies) delete fZombies;
265    if (fContainers) {fContainers->Delete(); delete fContainers;}
266    if (fInputs) delete fInputs;
267    if (fOutputs) delete fOutputs;
268    if (fParamCont) delete fParamCont;
269    if (fDebugOptions) delete fDebugOptions;
270    if (fGridHandler) delete fGridHandler;
271    if (fInputEventHandler) delete fInputEventHandler;
272    if (fOutputEventHandler) delete fOutputEventHandler;
273    if (fMCtruthEventHandler) delete fMCtruthEventHandler;
274    if (fEventPool) delete fEventPool;
275    if (fgAnalysisManager==this) fgAnalysisManager = NULL;
276    if (fGlobals) {fGlobals->DeleteAll(); delete fGlobals;}
277    if (fFileDescriptors) {fFileDescriptors->Delete(); delete fFileDescriptors;}
278    delete fIOTimer;
279    delete fCPUTimer;
280    delete fInitTimer;
281 }
282
283 //______________________________________________________________________________
284 void AliAnalysisManager::CreateReadCache()
285 {
286 // Create cache for reading according fCacheSize and fAsyncReading.
287    if (!fTree || !fTree->GetCurrentFile()) {
288       Error("CreateReadCache","Current tree or tree file not yet defined");
289       return;
290    }   
291    if (!fCacheSize) {
292       if (fDebug) Info("CreateReadCache","=== Read caching disabled ===");
293       return;
294    }
295 //   gEnv->SetValue("TFile.AsyncPrefetching",(Int_t)fAsyncReading);
296 //   if (fAsyncReading) gEnv->SetValue("Cache.Directory",Form("file://%s/cache", gSystem->WorkingDirectory()));
297    if (fAsyncReading) gEnv->SetValue("TFile.AsyncReading",1);
298    fTree->SetCacheSize(fCacheSize);
299    TTreeCache::SetLearnEntries(1);  //<<< we can take the decision after 1 entry
300    fTree->AddBranchToCache("*",kTRUE);    //<<< add all branches to the cache
301    if (fDebug) {
302       Info("CreateReadCache","Read cache enabled %lld bytes with async reading=%d",fCacheSize, (Int_t)fAsyncReading);
303    }
304    return;
305 }   
306       
307 //______________________________________________________________________________
308 Int_t AliAnalysisManager::GetEntry(Long64_t entry, Int_t getall)
309 {
310 // Read one entry of the tree or a whole branch.
311    fCurrentEntry = entry;
312    if (!fAutoBranchHandling)
313      return 123456789;
314    if (!fTree) return -1; 
315    fIOTimer->Start(kTRUE); 
316    Long64_t readbytes = fTree->GetTree()->GetEntry(entry, getall);
317    fIOTimer->Stop();
318    fIOTime += fIOTimer->RealTime();
319    return (Int_t)readbytes;
320 }
321
322 //______________________________________________________________________________
323 Int_t AliAnalysisManager::GetRunFromAlienPath(const char *path)
324 {
325 // Attempt to extract run number from input data path. Works only for paths to
326 // alice data in alien.
327 //    sim:  /alice/sim/<production>/run_no/...
328 //    data: /alice/data/year/period/000run_no/... (ESD or AOD)
329    TString type = "unknown";
330    TString s(path);
331    if (s.Contains("/alice/data")) type = "real";
332    else if (s.Contains("/alice/sim")) type = "simulated";
333    TString srun;
334    Int_t ind1, ind2;
335    ind1 = s.Index("/00");
336    if (ind1>0) {
337       ind2 = s.Index("/",ind1+1);
338       if (ind2-ind1>8) srun = s(ind1+1, ind2-ind1-1);
339    }   
340    if (srun.IsNull()) {
341       ind1 = s.Index("/LHC");
342       if (ind1>0) {
343          ind1 = s.Index("/",ind1+1);
344          if (ind1>0) {
345             ind2 = s.Index("/",ind1+1);
346             if (ind2>0) srun = s(ind1+1, ind2-ind1-1);
347          }
348       }
349    }         
350    Int_t run = srun.Atoi();
351    if (run>0) printf("=== GetRunFromAlienPath: run %d of %s data ===\n", run, type.Data());
352    return run;
353 }   
354
355 //______________________________________________________________________________
356 Bool_t AliAnalysisManager::Init(TTree *tree)
357 {
358   // The Init() function is called when the selector needs to initialize
359   // a new tree or chain. Typically here the branch addresses of the tree
360   // will be set. It is normaly not necessary to make changes to the
361   // generated code, but the routine can be extended by the user if needed.
362   // Init() will be called many times when running with PROOF.
363    Bool_t init = kFALSE;
364    if (!tree) return kFALSE; // Should not happen - protected in selector caller
365    if (fDebug > 1) {
366       printf("->AliAnalysisManager::Init(%s)\n", tree->GetName());
367    }
368    // Call InitTree of EventHandler
369    if (fOutputEventHandler) {
370       if (fMode == kProofAnalysis) {
371          init = fOutputEventHandler->Init(0x0, "proof");
372       } else {
373          init = fOutputEventHandler->Init(0x0, "local");
374       }
375       if (!init) {
376          Error("Init", "Output event handler failed to initialize");
377          return kFALSE;
378       }         
379    }
380    
381    if (fInputEventHandler) {
382       if (fMode == kProofAnalysis) {
383          init = fInputEventHandler->Init(tree, "proof");
384       } else {
385          init = fInputEventHandler->Init(tree, "local");
386       }
387       if (!init) {
388          Error("Init", "Input event handler failed to initialize tree"); 
389          return kFALSE;
390       }         
391    } else {
392       // If no input event handler we need to get the tree once
393       // for the chain
394       if(!tree->GetTree()) {
395          Long64_t readEntry = tree->LoadTree(0);
396          if (readEntry == -2) {
397             Error("Init", "Input tree has no entry. Exiting");
398             return kFALSE;
399          }
400       }   
401    }
402
403    if (fMCtruthEventHandler) {
404       if (fMode == kProofAnalysis) {
405          init = fMCtruthEventHandler->Init(0x0, "proof");
406       } else {
407          init = fMCtruthEventHandler->Init(0x0, "local");
408       }
409       if (!init) {
410          Error("Init", "MC event handler failed to initialize"); 
411          return kFALSE;
412       }         
413    }
414
415    if (!fInitOK) InitAnalysis();
416    if (!fInitOK) return kFALSE;
417    fTree = tree;
418    CreateReadCache();
419    fTable.Rehash(100);
420    AliAnalysisDataContainer *top = fCommonInput;
421    if (!top) top = (AliAnalysisDataContainer*)fInputs->At(0);
422    if (!top) {
423       Error("Init","No top input container !");
424       return kFALSE;
425    }
426    top->SetData(tree);
427    CheckBranches(kFALSE);
428    if (fDebug > 1) {
429       printf("<-AliAnalysisManager::Init(%s)\n", tree->GetName());
430    }
431    return kTRUE;
432 }
433
434 //______________________________________________________________________________
435 void AliAnalysisManager::SlaveBegin(TTree *tree)
436 {
437   // The SlaveBegin() function is called after the Begin() function.
438   // When running with PROOF SlaveBegin() is called on each slave server.
439   // The tree argument is deprecated (on PROOF 0 is passed).
440    if (fDebug > 1) printf("->AliAnalysisManager::SlaveBegin()\n");
441    // Init timer should be already started
442    // Apply debug options
443    ApplyDebugOptions();
444    if (fCacheSize && fMCtruthEventHandler) fMCtruthEventHandler->SetCacheSize(fCacheSize);
445    if (!CheckTasks()) Fatal("SlaveBegin", "Not all needed libraries were loaded");
446    static Bool_t isCalled = kFALSE;
447    Bool_t init = kFALSE;
448    Bool_t initOK = kTRUE;
449    TString msg;
450    TDirectory *curdir = gDirectory;
451    // Call SlaveBegin only once in case of mixing
452    if (isCalled && fMode==kMixingAnalysis) return;
453    gROOT->cd();
454    // Call Init of EventHandler
455    if (fOutputEventHandler) {
456       if (fMode == kProofAnalysis) {
457          // Merging AOD's in PROOF via TProofOutputFile
458          if (fDebug > 1) printf("   Initializing AOD output file %s...\n", fOutputEventHandler->GetOutputFileName());
459          init = fOutputEventHandler->Init("proof");
460          if (!init) msg = "Failed to initialize output handler on worker";
461       } else {
462          init = fOutputEventHandler->Init("local");
463          if (!init) msg = "Failed to initialize output handler";
464       }
465       initOK &= init;
466       if (!fSelector) Error("SlaveBegin", "Selector not set");
467       else if (!init) {fSelector->Abort(msg); fSelector->SetStatus(-1);}
468    }
469    gROOT->cd();
470    if (fInputEventHandler) {
471       fInputEventHandler->SetInputTree(tree);
472       if (fMode == kProofAnalysis) {
473          init = fInputEventHandler->Init("proof");
474          if (!init) msg = "Failed to initialize input handler on worker";
475       } else {
476          init = fInputEventHandler->Init("local");
477          if (!init) msg = "Failed to initialize input handler";
478       }
479       initOK &= init;
480       if (!fSelector) Error("SlaveBegin", "Selector not set");      
481       else if (!init) {fSelector->Abort(msg); fSelector->SetStatus(-1);}
482    }
483    gROOT->cd();
484    if (fMCtruthEventHandler) {
485       if (fMode == kProofAnalysis) {
486          init = fMCtruthEventHandler->Init("proof");
487          if (!init) msg = "Failed to initialize MC handler on worker";
488       } else {
489          init = fMCtruthEventHandler->Init("local");
490          if (!init) msg = "Failed to initialize MC handler";
491       }
492       initOK &= init;
493       if (!fSelector) Error("SlaveBegin", "Selector not set");      
494       else if (!init) {fSelector->Abort(msg); fSelector->SetStatus(-1);}
495    }
496    if (curdir) curdir->cd();
497    isCalled = kTRUE;
498    if (!initOK) return;   
499    TIter next(fTasks);
500    AliAnalysisTask *task;
501    // Call CreateOutputObjects for all tasks
502    Bool_t getsysInfo = ((fNSysInfo>0) && (fMode==kLocalAnalysis))?kTRUE:kFALSE;
503    Bool_t dirStatus = TH1::AddDirectoryStatus();
504    Int_t itask = 0;
505    while ((task=(AliAnalysisTask*)next())) {
506       gROOT->cd();
507       // Start with memory as current dir and make sure by default histograms do not get attached to files.
508       TH1::AddDirectory(kFALSE);
509       task->CreateOutputObjects();
510       if (!task->CheckPostData()) {
511          Error("SlaveBegin","####### IMPORTANT! ####### \n\n\n\
512                 Task %s (%s) did not call PostData() for all its outputs in (User)CreateOutputObjects()\n\n\
513                 ####### FIX YOUR CODE, THIS WILL PRODUCE A FATAL ERROR IN FUTURE! ##########", task->GetName(), task->ClassName());
514       }
515       if (getsysInfo) AliSysInfo::AddStamp(Form("%s_CREATEOUTOBJ",task->ClassName()), 0, itask, 0);
516       itask++;
517    }
518    TH1::AddDirectory(dirStatus);
519    if (curdir) curdir->cd();
520    fInitTimer->Stop();
521    fInitTime += fInitTimer->RealTime();
522    fInitTimer->Continue();
523    printf("Initialization time: %g [sec]\n", fInitTime);
524    if (fDebug > 1) printf("<-AliAnalysisManager::SlaveBegin()\n");
525 }
526
527 //______________________________________________________________________________
528 Bool_t AliAnalysisManager::Notify()
529 {
530    // The Notify() function is called when a new file is opened. This
531    // can be either for a new TTree in a TChain or when when a new TTree
532    // is started when using PROOF. It is normaly not necessary to make changes
533    // to the generated code, but the routine can be extended by the
534    // user if needed. The return value is currently not used.
535    fIOTimer->Start(kTRUE); 
536    if (!fTree) return kFALSE;
537    if (!TObject::TestBit(AliAnalysisManager::kTrueNotify)) return kFALSE;
538
539    fTable.Clear("nodelete"); // clearing the hash table may not be needed -> C.L.
540    if (fMode == kProofAnalysis) fIsRemote = kTRUE;
541
542    TFile *curfile = fTree->GetCurrentFile();
543    if (!curfile) {
544       Error("Notify","No current file");
545       return kFALSE;
546    }  
547    if (IsCollectThroughput()) {
548       if (fCurrentDescriptor) fCurrentDescriptor->Done();
549       fCurrentDescriptor = new AliAnalysisFileDescriptor(curfile);
550       fFileDescriptors->Add(fCurrentDescriptor);
551    }   
552    
553    if (fDebug > 1) printf("->AliAnalysisManager::Notify() file: %s\n", curfile->GetName());
554    Int_t run = AliAnalysisManager::GetRunFromAlienPath(curfile->GetName());
555    if (run && (run != fRunFromPath)) {
556       fRunFromPath = run;
557       if (fDebug > 1) printf("   ### run found from path: %d\n", run);
558    }
559    TIter next(fTasks);
560    AliAnalysisTask *task;
561         
562    // Call Notify of the event handlers
563    if (fInputEventHandler) {
564        fInputEventHandler->Notify(curfile->GetName());
565    }
566
567    if (fOutputEventHandler) {
568        fOutputEventHandler->Notify(curfile->GetName());
569    }
570
571    if (fMCtruthEventHandler) {
572        fMCtruthEventHandler->Notify(curfile->GetName());
573    }
574
575    // Call Notify for all tasks
576    while ((task=(AliAnalysisTask*)next())) 
577       task->Notify();
578
579    if (fDebug > 1) printf("<-AliAnalysisManager::Notify()\n");
580    fIOTimer->Stop();
581    fIOTime += fIOTimer->RealTime();
582    return kTRUE;
583 }    
584
585 //______________________________________________________________________________
586 Bool_t AliAnalysisManager::Process(Long64_t)
587 {
588   // The Process() function is called for each entry in the tree (or possibly
589   // keyed object in the case of PROOF) to be processed. The entry argument
590   // specifies which entry in the currently loaded tree is to be processed.
591   // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
592   // to read either all or the required parts of the data. When processing
593   // keyed objects with PROOF, the object is already loaded and is available
594   // via the fObject pointer.
595   //
596   // This function should contain the "body" of the analysis. It can contain
597   // simple or elaborate selection criteria, run algorithms on the data
598   // of the event and typically fill histograms.
599
600   // WARNING when a selector is used with a TChain, you must use
601   //  the pointer to the current TTree to call GetEntry(entry).
602   //  The entry is always the local entry number in the current tree.
603   //  Assuming that fChain is the pointer to the TChain being processed,
604   //  use fChain->GetTree()->GetEntry(entry).
605
606    // This method is obsolete. ExecAnalysis is called instead.
607    return kTRUE;
608 }
609
610 //______________________________________________________________________________
611 void AliAnalysisManager::PackOutput(TList *target)
612 {
613   // Pack all output data containers in the output list. Called at SlaveTerminate
614   // stage in PROOF case for each slave.
615    if (fDebug > 1) printf("->AliAnalysisManager::PackOutput()\n");
616    fIOTimer->Start(kTRUE);
617    std::ofstream out;
618    if (IsCollectThroughput()) {
619       if (fCurrentDescriptor) fCurrentDescriptor->Done();
620       fFileDescriptors->Print();
621       if (fFileInfoLog.IsNull()) fFileInfoLog = "fileinfo.log";
622       out.open(fFileInfoLog, std::ios::app);
623       if (out.bad()) Error("SavePrimitive", "Bad file name: %s", fFileInfoLog.Data());
624       else {
625          TIter nextflog(fFileDescriptors);
626          TObject *log;
627          while ((log=nextflog())) log->SavePrimitive(out,"");
628       }
629    }   
630    if (!target) {
631       Error("PackOutput", "No target. Exiting.");
632       return;
633    }
634    TDirectory *cdir = gDirectory;
635    gROOT->cd();
636    if (fInputEventHandler)   fInputEventHandler  ->Terminate();
637    if (fOutputEventHandler)  fOutputEventHandler ->Terminate();
638    if (fMCtruthEventHandler) fMCtruthEventHandler->Terminate();
639    gROOT->cd();
640
641    // Call FinishTaskOutput() for each event loop task (not called for 
642    // post-event loop tasks - use Terminate() fo those)
643    TIter nexttask(fTasks);
644    AliAnalysisTask *task;
645    while ((task=(AliAnalysisTask*)nexttask())) {
646       if (!task->IsPostEventLoop()) {
647          if (fDebug > 1) printf("->FinishTaskOutput: task %s\n", task->GetName());
648          task->FinishTaskOutput();
649          gROOT->cd();
650          if (fDebug > 1) printf("<-FinishTaskOutput: task %s\n", task->GetName());
651       }
652    }
653    // Write statistics message on the workers.
654    if (fStatistics) WriteStatisticsMsg(fNcalls);
655    
656    if (fMode == kProofAnalysis) {
657       TIter next(fOutputs);
658       AliAnalysisDataContainer *output;
659       Bool_t isManagedByHandler = kFALSE;
660       TList filestmp;
661       filestmp.SetOwner();
662       while ((output=(AliAnalysisDataContainer*)next())) {
663          // Do not consider outputs of post event loop tasks
664          isManagedByHandler = kFALSE;
665          if (output->GetProducer() && output->GetProducer()->IsPostEventLoop()) continue;
666          const char *filename = output->GetFileName();
667          if (!(strcmp(filename, "default")) && fOutputEventHandler) {
668             isManagedByHandler = kTRUE;
669             printf("#### Handler output. Extra: %s\n", fExtraFiles.Data());
670             filename = fOutputEventHandler->GetOutputFileName();
671          }
672          // Check if data was posted to this container. If not, issue an error.
673          if (!output->GetData() && !isManagedByHandler) {
674             Error("PackOutput", "No data for output container %s. Forgot to PostData ?", output->GetName());
675             continue;
676          }   
677          if (!output->IsSpecialOutput()) {
678             // Normal outputs
679             if (strlen(filename) && !isManagedByHandler) {
680                // Backup current folder
681                TDirectory *opwd = gDirectory;
682                // File resident outputs. 
683                // Check first if the file exists.
684                TString openoption = "RECREATE";
685                Bool_t firsttime = kTRUE;
686                if (filestmp.FindObject(output->GetFileName())) {
687                   firsttime = kFALSE;
688                } else {   
689                   filestmp.Add(new TNamed(output->GetFileName(),""));
690                }   
691                if (!gSystem->AccessPathName(output->GetFileName()) && !firsttime) openoption = "UPDATE";
692 //               TFile *file = AliAnalysisManager::OpenFile(output, openoption, kTRUE);
693                // Save data to file, then close.
694                if (output->GetData()->InheritsFrom(TCollection::Class())) {
695                   // If data is a collection, we set the name of the collection 
696                   // as the one of the container and we save as a single key.
697                   TCollection *coll = (TCollection*)output->GetData();
698                   coll->SetName(output->GetName());
699 //                  coll->Write(output->GetName(), TObject::kSingleKey);
700                } else {
701                   if (output->GetData()->InheritsFrom(TTree::Class())) {
702                      TFile *file = AliAnalysisManager::OpenFile(output, openoption, kTRUE);
703                      // Save data to file, then close.
704                      TTree *tree = (TTree*)output->GetData();
705                      // Check if tree is in memory
706                      if (tree->GetDirectory()==gROOT) tree->SetDirectory(gDirectory);
707                      tree->AutoSave();
708                      file->Close();
709                   } else {
710 //                     output->GetData()->Write();
711                   }   
712                }      
713                if (fDebug > 1) printf("PackOutput %s: memory merge, file resident output\n", output->GetName());
714 //               if (fDebug > 2) {
715 //                  printf("   file %s listing content:\n", filename);
716 //                  file->ls();
717 //               }   
718                // Clear file list to release object ownership to user.
719 //               file->Clear();
720 //               file->Close();
721                output->SetFile(NULL);
722                // Restore current directory
723                if (opwd) opwd->cd();
724             } else {
725                // Memory-resident outputs   
726                if (fDebug > 1) printf("PackOutput %s: memory merge memory resident output\n", filename);
727             }   
728             AliAnalysisDataWrapper *wrap = 0;
729             if (isManagedByHandler) {
730                wrap = new AliAnalysisDataWrapper(fOutputEventHandler->GetTree());
731                wrap->SetName(output->GetName());
732             }   
733             else                    wrap =output->ExportData();
734             // Output wrappers must NOT delete data after merging - the user owns them
735             wrap->SetDeleteData(kFALSE);
736             target->Add(wrap);
737          } else {
738          // Special outputs. The file must be opened and connected to the container.
739             TDirectory *opwd = gDirectory;
740             TFile *file = output->GetFile();
741             if (!file) {
742                AliAnalysisTask *producer = output->GetProducer();
743                Fatal("PackOutput", 
744                      "File %s for special container %s was NOT opened in %s::CreateOutputObjects !!!",
745                      output->GetFileName(), output->GetName(), producer->ClassName());
746                continue;
747             }   
748             TString outFilename = file->GetName();
749             if (fDebug > 1) printf("PackOutput %s: special output\n", output->GetName());
750             if (isManagedByHandler) {
751                // Terminate IO for files managed by the output handler
752                // file->Write() moved to AOD handler (A.G. 11.01.10)
753 //               if (file) file->Write();
754                if (file && fDebug > 2) {
755                   printf("   handled file %s listing content:\n", file->GetName());
756                   file->ls();
757                }   
758                fOutputEventHandler->TerminateIO();
759             } else {               
760                file->cd();
761                // Release object ownership to users after writing data to file
762                if (output->GetData()->InheritsFrom(TCollection::Class())) {
763                   // If data is a collection, we set the name of the collection 
764                   // as the one of the container and we save as a single key.
765                   TCollection *coll = (TCollection*)output->GetData();
766                   coll->SetName(output->GetName());
767                   coll->Write(output->GetName(), TObject::kSingleKey);
768                } else {
769                   if (output->GetData()->InheritsFrom(TTree::Class())) {
770                      TTree *tree = (TTree*)output->GetData();
771                      tree->SetDirectory(file);
772                      tree->AutoSave();
773                   } else {
774                      output->GetData()->Write();
775                   }   
776                }      
777                if (fDebug > 2) {
778                   printf("   file %s listing content:\n", output->GetFileName());
779                   file->ls();
780                }
781                // Clear file list to release object ownership to user.
782 //               file->Clear();
783                file->Close();
784                output->SetFile(NULL);
785             }
786             // Restore current directory
787             if (opwd) opwd->cd();
788             // Check if a special output location was provided or the output files have to be merged
789             if (strlen(fSpecialOutputLocation.Data())) {
790                TString remote = fSpecialOutputLocation;
791                remote += "/";
792                Int_t gid = gROOT->ProcessLine("gProofServ->GetGroupId();");
793                if (remote.BeginsWith("alien:")) {
794                   gROOT->ProcessLine("TGrid::Connect(\"alien:\", gProofServ->GetUser());");
795                   remote += outFilename;
796                   remote.ReplaceAll(".root", Form("_%d.root", gid));
797                } else {   
798                   remote += Form("%s_%d_", gSystem->HostName(), gid);
799                   remote += outFilename;
800                }   
801                if (fDebug > 1) 
802                   Info("PackOutput", "Output file for container %s to be copied \n   at: %s. No merging.",
803                        output->GetName(), remote.Data());
804                TFile::Cp ( outFilename.Data(), remote.Data() );
805                // Copy extra outputs
806                if (fExtraFiles.Length() && isManagedByHandler) {
807                   TObjArray *arr = fExtraFiles.Tokenize(" ");
808                   TObjString *os;
809                   TIter nextfilename(arr);
810                   while ((os=(TObjString*)nextfilename())) {
811                      outFilename = os->GetString();
812                      remote = fSpecialOutputLocation;
813                      remote += "/";
814                      if (remote.BeginsWith("alien://")) {
815                         remote += outFilename;
816                         remote.ReplaceAll(".root", Form("_%d.root", gid));
817                      } else {   
818                         remote += Form("%s_%d_", gSystem->HostName(), gid);
819                         remote += outFilename;
820                      }   
821                      if (fDebug > 1) 
822                         Info("PackOutput", "Extra AOD file %s to be copied \n   at: %s. No merging.",
823                              outFilename.Data(), remote.Data());
824                      TFile::Cp ( outFilename.Data(), remote.Data() );
825                   }   
826                   delete arr;
827                }   
828             } else {
829             // No special location specified-> use TProofOutputFile as merging utility
830             // The file at this output slot must be opened in CreateOutputObjects
831                if (fDebug > 1) printf("   File for container %s to be merged via file merger...\n", output->GetName());
832             }
833          }      
834       }
835    } 
836    fIOTime += fIOTimer->RealTime();
837    if ((fDebug || IsCollectThroughput())) {
838       fInitTimer->Stop();
839       fInitTime = fInitTimer->RealTime()-fIOTime-fCPUTime;
840       printf("=Analysis %s= init time:       %g[sec]\
841             \n              I/O & data mng.: %g [sec]\
842             \n              task execution: %g [sec]\
843             \n              total time:     CPU=%g [sec]  REAL=%g[sec]\n",
844             GetName(), fInitTime, fIOTime, fCPUTime, fInitTimer->CpuTime(), fInitTimer->RealTime());
845       if (IsCollectThroughput()) {
846          out << "#summary#########################################################" << endl;
847          out << "train_name   " << GetName() << endl;
848          out << "root_time    " << fInitTimer->RealTime() << endl;
849          out << "root_cpu     " << fInitTimer->CpuTime() << endl;
850          out << "init_time    " << fInitTime << endl;
851          out << "io_mng_time  " << fIOTime << endl;
852          out << "exec_time    " << fCPUTime << endl;
853          TString aliensite = gSystem->Getenv("ALIEN_SITE");
854          out << "alien_site   " << aliensite << endl;
855          out << "host_name    ";
856          TString hostname = gSystem->Getenv("ALIEN_HOSTNAME");
857          if (hostname.IsNull()) {
858             out.close();
859             gSystem->Exec(Form("hostname -f >> %s", fFileInfoLog.Data()));
860          } else {
861             out << hostname << endl;
862          }   
863       }
864    }
865               
866    if (cdir) cdir->cd();
867    if (fDebug > 1) printf("<-AliAnalysisManager::PackOutput: output list contains %d containers\n", target->GetSize());
868 }
869
870 //______________________________________________________________________________
871 void AliAnalysisManager::ImportWrappers(TList *source)
872 {
873 // Import data in output containers from wrappers coming in source.
874    if (fDebug > 1) printf("->AliAnalysisManager::ImportWrappers()\n");
875    fIOTimer->Start(kTRUE);
876    TIter next(fOutputs);
877    AliAnalysisDataContainer *cont;
878    AliAnalysisDataWrapper   *wrap;
879    Int_t icont = 0;
880    Bool_t inGrid = (fMode == kGridAnalysis)?kTRUE:kFALSE;
881    TDirectory *cdir = gDirectory;
882    while ((cont=(AliAnalysisDataContainer*)next())) {
883       wrap = 0;
884       if (cont->GetProducer() && cont->GetProducer()->IsPostEventLoop() && !inGrid) continue;
885       if (cont->IsRegisterDataset()) continue;
886       const char *filename = cont->GetFileName();
887       Bool_t isManagedByHandler = kFALSE;
888       if (!(strcmp(filename, "default")) && fOutputEventHandler) {
889          isManagedByHandler = kTRUE;
890          filename = fOutputEventHandler->GetOutputFileName();
891       }
892       if (cont->IsSpecialOutput() || inGrid) {
893          if (strlen(fSpecialOutputLocation.Data())) continue;
894          // Copy merged file from PROOF scratch space. 
895          // In case of grid the files are already in the current directory.
896          if (!inGrid) {
897             if (isManagedByHandler && fExtraFiles.Length()) {
898                // Copy extra registered dAOD files.
899                TObjArray *arr = fExtraFiles.Tokenize(" ");
900                TObjString *os;
901                TIter nextfilename(arr);
902                while ((os=(TObjString*)nextfilename())) GetFileFromWrapper(os->GetString(), source);
903                delete arr;
904             }
905             if (!GetFileFromWrapper(filename, source)) continue;
906          }   
907          // Normally we should connect data from the copied file to the
908          // corresponding output container, but it is not obvious how to do this
909          // automatically if several objects in file...
910          TFile *f = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
911          if (!f) f = TFile::Open(filename, "READ");
912          if (!f) {
913             Error("ImportWrappers", "Cannot open file %s in read-only mode", filename);
914             continue;
915          }   
916          f->cd();
917          TObject *obj = 0;
918          // Cd to the directory pointed by the container
919          TString folder = cont->GetFolderName();
920          if (!folder.IsNull()) f->cd(folder);
921          // Try to fetch first an object having the container name.
922          obj = gDirectory->Get(cont->GetName());
923          if (!obj) {
924             Warning("ImportWrappers", "Could not import object of type:%s for container %s in file %s:%s.\n Object will not be available in Terminate(). Try if possible to name the output object as the container (%s) or to embed it in a TList", 
925                     cont->GetType()->GetName(), cont->GetName(), filename, cont->GetFolderName(), cont->GetName());
926             continue;
927          }  
928          wrap = new AliAnalysisDataWrapper(obj);
929          wrap->SetDeleteData(kFALSE);
930       }   
931       if (!wrap) wrap = (AliAnalysisDataWrapper*)source->FindObject(cont->GetName());
932       if (!wrap) {
933          Error("ImportWrappers","Container %s not found in analysis output !", cont->GetName());
934          continue;
935       }
936       icont++;
937       if (fDebug > 1) {
938          printf("   Importing data for container %s\n", cont->GetName());
939          if (strlen(filename)) printf("    -> file %s\n", filename);
940          else printf("\n");
941       }   
942       cont->ImportData(wrap);
943    }
944    if (cdir) cdir->cd();
945    fIOTimer->Stop();
946    fIOTime += fIOTimer->RealTime();
947    if (fDebug > 1) printf("<-AliAnalysisManager::ImportWrappers(): %d containers imported\n", icont);
948 }
949
950 //______________________________________________________________________________
951 void AliAnalysisManager::UnpackOutput(TList *source)
952 {
953   // Called by AliAnalysisSelector::Terminate only on the client.
954    fIOTimer->Start(kTRUE);
955    if (fDebug > 1) printf("->AliAnalysisManager::UnpackOutput()\n");
956    if (!source) {
957       Error("UnpackOutput", "No target. Exiting.");
958       return;
959    }
960    if (fDebug > 1) printf("   Source list contains %d containers\n", source->GetSize());
961
962    if (fMode == kProofAnalysis) ImportWrappers(source);
963
964    TIter next(fOutputs);
965    AliAnalysisDataContainer *output;
966    while ((output=(AliAnalysisDataContainer*)next())) {
967       if (!output->GetData()) continue;
968       // Check if there are client tasks that run post event loop
969       if (output->HasConsumers()) {
970          // Disable event loop semaphore
971          output->SetPostEventLoop(kTRUE);
972          TObjArray *list = output->GetConsumers();
973          Int_t ncons = list->GetEntriesFast();
974          for (Int_t i=0; i<ncons; i++) {
975             AliAnalysisTask *task = (AliAnalysisTask*)list->At(i);
976             task->CheckNotify(kTRUE);
977             // If task is active, execute it
978             if (task->IsPostEventLoop() && task->IsActive()) {
979                if (fDebug > 1) printf("== Executing post event loop task %s\n", task->GetName());
980                task->ExecuteTask();
981             }   
982          }
983       }   
984    }
985    fIOTimer->Stop();
986    fIOTime += fIOTimer->RealTime();
987    if (fDebug > 1) printf("<-AliAnalysisManager::UnpackOutput()\n");
988 }
989
990 //______________________________________________________________________________
991 void AliAnalysisManager::Terminate()
992 {
993   // The Terminate() function is the last function to be called during
994   // a query. It always runs on the client, it can be used to present
995   // the results graphically.
996    if (fDebug > 1) printf("->AliAnalysisManager::Terminate()\n");
997    fInitTimer->Start(kTRUE);
998    TDirectory *cdir = gDirectory;
999    gROOT->cd();
1000    AliAnalysisTask *task;
1001    AliAnalysisDataContainer *output;
1002    TIter next(fTasks);
1003    TStopwatch timer;
1004    Bool_t getsysInfo = ((fNSysInfo>0) && (fMode==kLocalAnalysis))?kTRUE:kFALSE;
1005    // Call Terminate() for tasks
1006    Int_t itask = 0;
1007    while (!IsSkipTerminate() && (task=(AliAnalysisTask*)next())) {
1008       // Save all the canvases produced by the Terminate
1009       TString pictname = Form("%s_%s", task->GetName(), task->ClassName());
1010       task->Terminate();
1011       gROOT->cd();
1012       if (getsysInfo) 
1013          AliSysInfo::AddStamp(Form("%s_TERMINATE",task->ClassName()),0, itask, 2);
1014       itask++;   
1015       if (TObject::TestBit(kSaveCanvases)) {
1016          if (!gROOT->IsBatch()) {
1017             if (fDebug>1) printf("Waiting 5 sec for %s::Terminate() to finish drawing ...\n", task->ClassName());
1018             timer.Start(kTRUE);
1019             while (timer.RealTime()<5) {
1020                timer.Continue();
1021                gSystem->ProcessEvents();
1022             }
1023          }
1024          Int_t iend = gROOT->GetListOfCanvases()->GetEntries();
1025          if (iend==0) continue;
1026          TCanvas *canvas;
1027          for (Int_t ipict=0; ipict<iend; ipict++) {
1028             canvas = (TCanvas*)gROOT->GetListOfCanvases()->At(ipict);
1029             if (!canvas) continue;         
1030             canvas->SaveAs(Form("%s_%02d.gif", pictname.Data(),ipict));
1031          } 
1032          gROOT->GetListOfCanvases()->Delete(); 
1033       }
1034    }   
1035    //
1036    if (fInputEventHandler)   fInputEventHandler  ->TerminateIO();
1037    if (fOutputEventHandler)  fOutputEventHandler ->TerminateIO();
1038    if (fMCtruthEventHandler) fMCtruthEventHandler->TerminateIO();
1039    gROOT->cd();
1040    TObjArray *allOutputs = new TObjArray();
1041    Int_t icont;
1042    for (icont=0; icont<fOutputs->GetEntriesFast(); icont++) allOutputs->Add(fOutputs->At(icont));
1043    if (!IsSkipTerminate())
1044       for (icont=0; icont<fParamCont->GetEntriesFast(); icont++) allOutputs->Add(fParamCont->At(icont));
1045    TIter next1(allOutputs);
1046    TString handlerFile = "";
1047    TString extraOutputs = "";
1048    if (fOutputEventHandler) {
1049       handlerFile = fOutputEventHandler->GetOutputFileName();
1050       extraOutputs = fOutputEventHandler->GetExtraOutputs();
1051    }
1052    icont = 0;
1053    TList filestmp;
1054    while ((output=(AliAnalysisDataContainer*)next1())) {
1055       // Special outputs or grid files have the files already closed and written.
1056       icont++;
1057       if (fMode == kGridAnalysis && icont<=fOutputs->GetEntriesFast()) continue;
1058       if (fMode == kProofAnalysis) {
1059         if (output->IsSpecialOutput() || output->IsRegisterDataset()) continue;
1060       }  
1061       const char *filename = output->GetFileName();
1062       TString openoption = "RECREATE";
1063       if (!(strcmp(filename, "default"))) continue;
1064       if (!strlen(filename)) continue;
1065       if (!output->GetData()) continue;
1066       TDirectory *opwd = gDirectory;
1067       TFile *file = output->GetFile();
1068       if (!file) file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
1069       if (!file) {
1070               //if (handlerFile == filename && !gSystem->AccessPathName(filename)) openoption = "UPDATE";
1071          Bool_t firsttime = kTRUE;
1072          if (filestmp.FindObject(filename) || extraOutputs.Contains(filename)) {
1073             firsttime = kFALSE;
1074          } else {   
1075             filestmp.Add(new TNamed(filename,""));
1076          }   
1077          if (!gSystem->AccessPathName(filename) && !firsttime) openoption = "UPDATE";
1078               if (fDebug>1) printf("Opening file: %s  option=%s\n",filename, openoption.Data());
1079          file = new TFile(filename, openoption);
1080       } else {
1081          if (fDebug>1) printf("File <%s> already opened with option: <%s> \n", filename, file->GetOption());
1082          openoption = file->GetOption();
1083          if (openoption == "READ") {
1084             if (fDebug>1) printf("...reopening in UPDATE mode\n");
1085             file->ReOpen("UPDATE");            
1086          }
1087       }   
1088       if (file->IsZombie()) {
1089          Error("Terminate", "Cannot open output file %s", filename);
1090          continue;
1091       }   
1092       output->SetFile(file);
1093       file->cd();
1094       // Check for a folder request
1095       TString dir = output->GetFolderName();
1096       if (!dir.IsNull()) {
1097          if (!file->GetDirectory(dir)) file->mkdir(dir);
1098          file->cd(dir);
1099       }  
1100       if (fDebug > 1) printf("...writing container %s to file %s:%s\n", output->GetName(), file->GetName(), output->GetFolderName());
1101       if (output->GetData()->InheritsFrom(TCollection::Class())) {
1102       // If data is a collection, we set the name of the collection 
1103       // as the one of the container and we save as a single key.
1104          TCollection *coll = (TCollection*)output->GetData();
1105          coll->SetName(output->GetName());
1106          coll->Write(output->GetName(), TObject::kSingleKey);
1107       } else {
1108          if (output->GetData()->InheritsFrom(TTree::Class())) {
1109             TTree *tree = (TTree*)output->GetData();
1110             tree->SetDirectory(gDirectory);
1111             tree->AutoSave();
1112          } else {
1113             output->GetData()->Write();
1114          }   
1115       }      
1116       if (opwd) opwd->cd();
1117    }
1118    gROOT->cd();
1119    next1.Reset();
1120    TString copiedFiles;
1121    while ((output=(AliAnalysisDataContainer*)next1())) {
1122       // Close all files at output
1123       TDirectory *opwd = gDirectory;
1124       if (output->GetFile()) {
1125          // Clear file list to release object ownership to user.
1126 //         output->GetFile()->Clear();
1127          output->GetFile()->Close();
1128          // Copy merged outputs in alien if requested
1129          if (fSpecialOutputLocation.BeginsWith("alien://")) {
1130             if (copiedFiles.Contains(output->GetFile()->GetName())) {
1131                if (opwd) opwd->cd();
1132                output->SetFile(NULL);
1133                continue;
1134             } 
1135             Info("Terminate", "Copy file %s to %s", output->GetFile()->GetName(),fSpecialOutputLocation.Data()); 
1136             gROOT->ProcessLine("if (!gGrid) TGrid::Connect(\"alien:\");");
1137             TFile::Cp(output->GetFile()->GetName(), 
1138                       Form("%s/%s", fSpecialOutputLocation.Data(), output->GetFile()->GetName()));
1139             copiedFiles += output->GetFile()->GetName();
1140          }             
1141          output->SetFile(NULL);
1142       }   
1143       if (opwd) opwd->cd();
1144    }   
1145    delete allOutputs;
1146    //Write statistics information on the client
1147    if (fStatistics) WriteStatisticsMsg(fNcalls);
1148    if (getsysInfo) {
1149       TDirectory *crtdir = gDirectory;
1150       TFile f("syswatch.root", "RECREATE");
1151       TH1 *hist;
1152       TString cut;
1153       if (!f.IsZombie()) {
1154          TTree *tree = AliSysInfo::MakeTree("syswatch.log");
1155          tree->SetName("syswatch");
1156          tree->SetMarkerStyle(kCircle);
1157          tree->SetMarkerColor(kBlue);
1158          tree->SetMarkerSize(0.5);
1159          if (!gROOT->IsBatch()) {
1160             tree->SetAlias("event", "id0");
1161             tree->SetAlias("task",  "id1");
1162             tree->SetAlias("stage", "id2");
1163             // Already defined aliases
1164             // tree->SetAlias("deltaT","stampSec-stampOldSec");
1165             // tree->SetAlias("T","stampSec-first");
1166             // tree->SetAlias("deltaVM","(pI.fMemVirtual-pIOld.fMemVirtual)");
1167             // tree->SetAlias("VM","pI.fMemVirtual");
1168             TCanvas *canvas = new TCanvas("SysInfo","SysInfo",10,10,1200,800);
1169             Int_t npads = 1 /*COO plot for all tasks*/ +
1170                           fTopTasks->GetEntries() /*Exec plot per task*/ +
1171                           1 /*Terminate plot for all tasks*/ +
1172                           1; /*vm plot*/
1173                           
1174             Int_t iopt = (Int_t)TMath::Sqrt((Double_t)npads);
1175             if (npads<iopt*(iopt+1))
1176                canvas->Divide(iopt, iopt+1, 0.01, 0.01);
1177             else
1178                canvas->Divide(iopt+1, iopt+1, 0.01, 0.01);
1179             Int_t ipad = 1;
1180             // draw the plot of deltaVM for Exec for each task
1181             for (itask=0; itask<fTopTasks->GetEntriesFast(); itask++) {
1182                task = (AliAnalysisTask*)fTopTasks->At(itask);
1183                canvas->cd(ipad++);
1184                cut = Form("task==%d && stage==1", itask);
1185                tree->Draw("deltaVM:event",cut,"", 1234567890, 0);
1186                hist = (TH1*)gPad->GetListOfPrimitives()->FindObject("htemp");            
1187                if (hist) {
1188                   hist->SetTitle(Form("%s: Exec dVM[MB]/event", task->GetName()));
1189                   hist->GetYaxis()->SetTitle("deltaVM [MB]");
1190                }   
1191             }
1192             // Draw the plot of deltaVM for CreateOutputObjects for all tasks
1193             canvas->cd(ipad++);
1194             tree->SetMarkerStyle(kFullTriangleUp);
1195             tree->SetMarkerColor(kRed);
1196             tree->SetMarkerSize(0.8);
1197             cut = "task>=0 && task<1000 && stage==0";
1198             tree->Draw("deltaVM:sname",cut,"", 1234567890, 0);
1199             hist = (TH1*)gPad->GetListOfPrimitives()->FindObject("htemp");            
1200             if (hist) {
1201                hist->SetTitle("Memory in CreateOutputObjects()");
1202                hist->GetYaxis()->SetTitle("deltaVM [MB]");
1203                hist->GetXaxis()->SetTitle("task");
1204             }   
1205             // draw the plot of deltaVM for Terminate for all tasks
1206             canvas->cd(ipad++);
1207             tree->SetMarkerStyle(kOpenSquare);
1208             tree->SetMarkerColor(kMagenta);
1209             cut = "task>=0 && task<1000 && stage==2";
1210             tree->Draw("deltaVM:sname",cut,"", 1234567890, 0);
1211             hist = (TH1*)gPad->GetListOfPrimitives()->FindObject("htemp");
1212             if (hist) {
1213                hist->SetTitle("Memory in Terminate()");
1214                hist->GetYaxis()->SetTitle("deltaVM [MB]");
1215                hist->GetXaxis()->SetTitle("task");
1216             }   
1217             // Full VM profile
1218             canvas->cd(ipad++);
1219             tree->SetMarkerStyle(kFullCircle);
1220             tree->SetMarkerColor(kGreen);
1221             cut = Form("task==%d && stage==1",fTopTasks->GetEntriesFast()-1);            
1222             tree->Draw("VM:event",cut,"", 1234567890, 0);
1223             hist = (TH1*)gPad->GetListOfPrimitives()->FindObject("htemp");
1224             if (hist) {
1225                hist->SetTitle("Virtual memory");
1226                hist->GetYaxis()->SetTitle("VM [MB]");
1227             }
1228             canvas->Modified();   
1229          }   
1230          tree->SetMarkerStyle(kCircle);
1231          tree->SetMarkerColor(kBlue);
1232          tree->SetMarkerSize(0.5);
1233          tree->Write();
1234          f.Close();
1235          delete tree;
1236       }
1237       if (crtdir) crtdir->cd();
1238    }
1239    // Validate the output files
1240    if (ValidateOutputFiles() && fIsRemote && fMode!=kProofAnalysis) {
1241       ofstream out;
1242       out.open("outputs_valid", ios::out);
1243       out.close();
1244    }
1245    if (cdir) cdir->cd();  
1246    fInitTimer->Stop();
1247    if (fDebug || IsCollectThroughput()) {
1248       printf("=Analysis %s= Terminate time:  %g[sec]\n", GetName(), fInitTimer->RealTime());
1249    }
1250    if (fDebug > 1) printf("<-AliAnalysisManager::Terminate()\n");
1251 }
1252 //______________________________________________________________________________
1253 void AliAnalysisManager::ProfileTask(Int_t itop, const char *option) const
1254 {
1255 // Profiles the task having the itop index in the list of top (first level) tasks.
1256    AliAnalysisTask *task = (AliAnalysisTask*)fTopTasks->At(itop);
1257    if (!task) {
1258       Error("ProfileTask", "There are only %d top tasks in the manager", fTopTasks->GetEntries());
1259       return;
1260    }
1261    ProfileTask(task->GetName(), option);
1262 }      
1263
1264 //______________________________________________________________________________
1265 void AliAnalysisManager::ProfileTask(const char *name, const char */*option*/) const
1266 {
1267 // Profile a managed task after the execution of the analysis in case NSysInfo
1268 // was used.
1269    if (gSystem->AccessPathName("syswatch.root")) {
1270       Error("ProfileTask", "No file syswatch.root found in the current directory");
1271       return;
1272    }
1273    if (gROOT->IsBatch()) return;
1274    AliAnalysisTask *task = (AliAnalysisTask*)fTopTasks->FindObject(name);
1275    if (!task) {
1276       Error("ProfileTask", "No top task named %s known by the manager.", name);
1277       return;
1278    }
1279    Int_t itop = fTopTasks->IndexOf(task);
1280    Int_t itask = fTasks->IndexOf(task);
1281    // Create canvas with 2 pads: first draw COO + Terminate, second Exec
1282    TDirectory *cdir = gDirectory;
1283    TFile f("syswatch.root");
1284    TTree *tree = (TTree*)f.Get("syswatch");
1285    if (!tree) {
1286       Error("ProfileTask", "No tree named <syswatch> found in file syswatch.root");
1287       return;
1288    }   
1289    if (fDebug > 1) printf("=== Profiling task %s (class %s)\n", name, task->ClassName());
1290    TCanvas *canvas = new TCanvas(Form("profile_%d",itop),Form("Profile of task %s (class %s)",name,task->ClassName()),10,10,800,600);
1291    canvas->Divide(2, 2, 0.01, 0.01);
1292    Int_t ipad = 1;
1293    TString cut;
1294    TH1 *hist;
1295    // VM profile for COO and Terminate methods
1296    canvas->cd(ipad++);
1297    cut = Form("task==%d && (stage==0 || stage==2)",itask);
1298    tree->Draw("deltaVM:sname",cut,"", 1234567890, 0);
1299    hist = (TH1*)gPad->GetListOfPrimitives()->FindObject("htemp");
1300    if (hist) {
1301       hist->SetTitle("Alocated VM[MB] for COO and Terminate");
1302       hist->GetYaxis()->SetTitle("deltaVM [MB]");
1303       hist->GetXaxis()->SetTitle("method");
1304    }   
1305    // CPU profile per event
1306    canvas->cd(ipad++);
1307    cut = Form("task==%d && stage==1",itop);
1308    tree->Draw("deltaT:event",cut,"", 1234567890, 0);
1309    hist = (TH1*)gPad->GetListOfPrimitives()->FindObject("htemp");
1310    if (hist) {
1311       hist->SetTitle("Execution time per event");
1312       hist->GetYaxis()->SetTitle("CPU/event [s]");
1313    }   
1314    // VM profile for Exec
1315    canvas->cd(ipad++);
1316    cut = Form("task==%d && stage==1",itop);
1317    tree->Draw("deltaVM:event",cut,"", 1234567890, 0);
1318    hist = (TH1*)gPad->GetListOfPrimitives()->FindObject("htemp");
1319    if (hist) {
1320       hist->SetTitle("Alocated VM[MB] per event");
1321       hist->GetYaxis()->SetTitle("deltaVM [MB]");
1322    }   
1323    canvas->Modified();
1324    delete tree;
1325    f.Close();
1326    if (cdir) cdir->cd();
1327 }     
1328
1329 //______________________________________________________________________________
1330 void AliAnalysisManager::AddTask(AliAnalysisTask *task)
1331 {
1332 // Adds a user task to the global list of tasks.
1333    if (fInitOK) {
1334       Error("AddTask", "Cannot add task %s since InitAnalysis was already called", task->GetName());
1335       return;
1336    }   
1337       
1338    if (fTasks->FindObject(task)) {
1339       Warning("AddTask", "Task %s: the same object already added to the analysis manager. Not adding.", task->GetName());
1340       return;
1341    }   
1342    task->SetActive(kFALSE);
1343    fTasks->Add(task);
1344 }  
1345
1346 //______________________________________________________________________________
1347 AliAnalysisTask *AliAnalysisManager::GetTask(const char *name) const
1348 {
1349 // Retreive task by name.
1350    if (!fTasks) return NULL;
1351    return (AliAnalysisTask*)fTasks->FindObject(name);
1352 }
1353
1354 //______________________________________________________________________________
1355 AliAnalysisDataContainer *AliAnalysisManager::CreateContainer(const char *name, 
1356                                 TClass *datatype, EAliAnalysisContType type, const char *filename)
1357 {
1358 // Create a data container of a certain type. Types can be:
1359 //   kExchangeContainer  = 0, used to exchange data between tasks
1360 //   kInputContainer   = 1, used to store input data
1361 //   kOutputContainer  = 2, used for writing result to a file
1362 // filename: composed by file#folder (e.g. results.root#INCLUSIVE) - will write
1363 // the output object to a folder inside the output file
1364    if (fContainers->FindObject(name)) {
1365       Error("CreateContainer","A container named %s already defined !",name);
1366       return NULL;
1367    }   
1368    AliAnalysisDataContainer *cont = new AliAnalysisDataContainer(name, datatype);
1369    fContainers->Add(cont);
1370    switch (type) {
1371       case kInputContainer:
1372          fInputs->Add(cont);
1373          break;
1374       case kOutputContainer:
1375          fOutputs->Add(cont);
1376          if (filename && strlen(filename)) {
1377             cont->SetFileName(filename);
1378             cont->SetDataOwned(kFALSE);  // data owned by the file
1379          }   
1380          break;
1381       case kParamContainer:
1382          fParamCont->Add(cont);
1383          if (filename && strlen(filename)) {
1384             cont->SetFileName(filename);
1385             cont->SetDataOwned(kFALSE);  // data owned by the file
1386          }   
1387          break;
1388       case kExchangeContainer:
1389          break;   
1390    }
1391    return cont;
1392 }
1393          
1394 //______________________________________________________________________________
1395 Bool_t AliAnalysisManager::ConnectInput(AliAnalysisTask *task, Int_t islot,
1396                                         AliAnalysisDataContainer *cont)
1397 {
1398 // Connect input of an existing task to a data container.
1399    if (!task) {
1400       Error("ConnectInput", "Task pointer is NULL");
1401       return kFALSE;
1402    }   
1403    if (!fTasks->FindObject(task)) {
1404       AddTask(task);
1405       Info("ConnectInput", "Task %s was not registered. Now owned by analysis manager", task->GetName());
1406    } 
1407    Bool_t connected = task->ConnectInput(islot, cont);
1408    return connected;
1409 }   
1410
1411 //______________________________________________________________________________
1412 Bool_t AliAnalysisManager::ConnectOutput(AliAnalysisTask *task, Int_t islot,
1413                                         AliAnalysisDataContainer *cont)
1414 {
1415 // Connect output of an existing task to a data container.
1416    if (!task) {
1417       Error("ConnectOutput", "Task pointer is NULL");
1418       return kFALSE;
1419    }   
1420    if (!fTasks->FindObject(task)) {
1421       AddTask(task);
1422       Warning("ConnectOutput", "Task %s not registered. Now owned by analysis manager", task->GetName());
1423    } 
1424    Bool_t connected = task->ConnectOutput(islot, cont);
1425    return connected;
1426 }   
1427                                
1428 //______________________________________________________________________________
1429 void AliAnalysisManager::CleanContainers()
1430 {
1431 // Clean data from all containers that have already finished all client tasks.
1432    TIter next(fContainers);
1433    AliAnalysisDataContainer *cont;
1434    while ((cont=(AliAnalysisDataContainer *)next())) {
1435       if (cont->IsOwnedData() && 
1436           cont->IsDataReady() && 
1437           cont->ClientsExecuted()) cont->DeleteData();
1438    }
1439 }
1440
1441 //______________________________________________________________________________
1442 Bool_t AliAnalysisManager::InitAnalysis()
1443 {
1444 // Initialization of analysis chain of tasks. Should be called after all tasks
1445 // and data containers are properly connected
1446    // Reset flag and remove valid_outputs file if exists
1447    if (fInitOK) return kTRUE;
1448    if (!gSystem->AccessPathName("outputs_valid"))
1449       gSystem->Unlink("outputs_valid");
1450    // Check for top tasks (depending only on input data containers)
1451    if (!fTasks->First()) {
1452       Error("InitAnalysis", "Analysis has no tasks !");
1453       return kFALSE;
1454    }   
1455    TIter next(fTasks);
1456    AliAnalysisTask *task;
1457    AliAnalysisDataContainer *cont;
1458    Int_t ntop = 0;
1459    Int_t nzombies = 0;
1460    Bool_t iszombie = kFALSE;
1461    Bool_t istop = kTRUE;
1462    Int_t i;
1463    while ((task=(AliAnalysisTask*)next())) {
1464       istop = kTRUE;
1465       iszombie = kFALSE;
1466       Int_t ninputs = task->GetNinputs();
1467       for (i=0; i<ninputs; i++) {
1468          cont = task->GetInputSlot(i)->GetContainer();
1469          if (!cont) {
1470             if (!iszombie) {
1471                task->SetZombie();
1472                fZombies->Add(task);
1473                nzombies++;
1474                iszombie = kTRUE;
1475             }   
1476             Error("InitAnalysis", "Input slot %d of task %s has no container connected ! Declared zombie...", 
1477                   i, task->GetName()); 
1478          }
1479          if (iszombie) continue;
1480          // Check if cont is an input container
1481          if (istop && !fInputs->FindObject(cont)) istop=kFALSE;
1482          // Connect to parent task
1483       }
1484       if (istop) {
1485          ntop++;
1486          fTopTasks->Add(task);
1487       }
1488    }
1489    if (!ntop) {
1490       Error("InitAnalysis", "No top task defined. At least one task should be connected only to input containers");
1491       return kFALSE;
1492    }                        
1493    // Check now if there are orphan tasks
1494    for (i=0; i<ntop; i++) {
1495       task = (AliAnalysisTask*)fTopTasks->At(i);
1496       task->SetUsed();
1497    }
1498    Int_t norphans = 0;
1499    next.Reset();
1500    while ((task=(AliAnalysisTask*)next())) {
1501       if (!task->IsUsed()) {
1502          norphans++;
1503          Warning("InitAnalysis", "Task %s is orphan", task->GetName());
1504       }   
1505    }          
1506    // Check the task hierarchy (no parent task should depend on data provided
1507    // by a daughter task)
1508    for (i=0; i<ntop; i++) {
1509       task = (AliAnalysisTask*)fTopTasks->At(i);
1510       if (task->CheckCircularDeps()) {
1511          Error("InitAnalysis", "Found illegal circular dependencies between following tasks:");
1512          PrintStatus("dep");
1513          return kFALSE;
1514       }   
1515    }
1516    // Check that all containers feeding post-event loop tasks are in the outputs list
1517    TIter nextcont(fContainers); // loop over all containers
1518    while ((cont=(AliAnalysisDataContainer*)nextcont())) {
1519       if (!cont->IsPostEventLoop() && !fOutputs->FindObject(cont)) {
1520          if (cont->HasConsumers()) {
1521          // Check if one of the consumers is post event loop
1522             TIter nextconsumer(cont->GetConsumers());
1523             while ((task=(AliAnalysisTask*)nextconsumer())) {
1524                if (task->IsPostEventLoop()) {
1525                   fOutputs->Add(cont);
1526                   break;
1527                }
1528             }
1529          }
1530       }
1531    }   
1532    // Check if all special output containers have a file name provided
1533    TIter nextout(fOutputs);
1534    while ((cont=(AliAnalysisDataContainer*)nextout())) {
1535       if (cont->IsSpecialOutput() && !strlen(cont->GetFileName())) {
1536          Error("InitAnalysis", "Wrong container %s : a file name MUST be provided for special outputs", cont->GetName());
1537          return kFALSE;
1538       }
1539    }
1540    // Initialize requested branch list if needed
1541    if (!fAutoBranchHandling) {
1542       next.Reset();
1543       while ((task=(AliAnalysisTask*)next())) {
1544          if (!task->HasBranches()) {
1545             Error("InitAnalysis", "Manual branch loading requested but task %s of type %s does not define branches.\nUse: fBranchNames = \"ESD:br1,br2,...,brN AOD:bra1,bra2,...,braM\"",
1546                   task->GetName(), task->ClassName());
1547             return kFALSE;
1548          }
1549          if (!fInputEventHandler || !strlen(fInputEventHandler->GetDataType())) {
1550             Error("InitAnalysis", "Manual branch loading requested but no input handler defined or handler does not define data type.");
1551             return kFALSE;
1552          }
1553          TString taskbranches;
1554          task->GetBranches(fInputEventHandler->GetDataType(), taskbranches);
1555          if (taskbranches.IsNull()) {
1556             Error("InitAnalysis", "Manual branch loading requested but task %s of type %s does not define branches of type %s:",
1557                   task->GetName(), task->ClassName(), fInputEventHandler->GetDataType());
1558             return kFALSE;      
1559          }
1560          AddBranches(taskbranches);
1561       }         
1562    }
1563    fInitOK = kTRUE;
1564    return kTRUE;
1565 }   
1566
1567 //______________________________________________________________________________
1568 void AliAnalysisManager::AddBranches(const char *branches)
1569 {
1570 // Add branches to the existing fRequestedBranches.
1571    TString br(branches);
1572    TObjArray *arr = br.Tokenize(",");
1573    TIter next(arr);
1574    TObject *obj;
1575    while ((obj=next())) {
1576       if (!fRequestedBranches.Contains(obj->GetName())) {
1577          if (!fRequestedBranches.IsNull()) fRequestedBranches += ",";
1578          fRequestedBranches += obj->GetName();
1579       }
1580    }
1581   delete arr;
1582 }   
1583
1584 //______________________________________________________________________________
1585 void AliAnalysisManager::CheckBranches(Bool_t load)
1586 {
1587 // The method checks the input branches to be loaded during the analysis.
1588    if (fAutoBranchHandling || fRequestedBranches.IsNull() || !fTree) return;   
1589    TObjArray *arr = fRequestedBranches.Tokenize(",");
1590    TIter next(arr);
1591    TObject *obj;
1592    while ((obj=next())) {
1593       TBranch *br = dynamic_cast<TBranch*>(fTable.FindObject(obj->GetName()));
1594       if (!br) {
1595          br = fTree->GetBranch(obj->GetName());
1596          if (!br) {
1597             Error("CheckBranches", "Could not find branch %s",obj->GetName());
1598             continue;
1599          }
1600       }   
1601       fTable.Add(br);
1602       if (load && br->GetReadEntry()!=GetCurrentEntry()) {
1603          br->GetEntry(GetCurrentEntry());
1604       }      
1605    }
1606   delete arr;
1607 }
1608
1609 //______________________________________________________________________________
1610 Bool_t AliAnalysisManager::CheckTasks() const
1611 {
1612 // Check consistency of tasks.
1613    Int_t ntasks = fTasks->GetEntries();
1614    if (!ntasks) {
1615       Error("CheckTasks", "No tasks connected to the manager. This may be due to forgetting to compile the task or to load their library.");
1616       return kFALSE;
1617    }
1618    // Get the pointer to AliAnalysisTaskSE::Class()
1619    TClass *badptr = (TClass*)gROOT->ProcessLine("AliAnalysisTaskSE::Class()");
1620    // Loop all tasks to check if their corresponding library was loaded
1621    TIter next(fTasks);
1622    TObject *obj;
1623    while ((obj=next())) {
1624       if (obj->IsA() == badptr) {
1625          Error("CheckTasks", "##################\n \
1626          Class for task %s NOT loaded. You probably forgot to load the library for this task (or compile it dynamically).\n###########################\n",obj->GetName());
1627          return kFALSE;
1628       }
1629    }
1630    return kTRUE;      
1631 }   
1632
1633 //______________________________________________________________________________
1634 void AliAnalysisManager::PrintStatus(Option_t *option) const
1635 {
1636 // Print task hierarchy.
1637    if (!fInitOK) {
1638       Info("PrintStatus", "Analysis manager %s not initialized : call InitAnalysis() first", GetName());
1639       return;
1640    }   
1641    Bool_t getsysInfo = ((fNSysInfo>0) && (fMode==kLocalAnalysis))?kTRUE:kFALSE;
1642    if (getsysInfo)
1643       Info("PrintStatus", "System information will be collected each %lld events", fNSysInfo);
1644    TIter next(fTopTasks);
1645    AliAnalysisTask *task;
1646    while ((task=(AliAnalysisTask*)next()))
1647       task->PrintTask(option);
1648   
1649    if (!fAutoBranchHandling && !fRequestedBranches.IsNull()) 
1650       printf("Requested input branches:\n%s\n", fRequestedBranches.Data());
1651   
1652   TString sopt(option);
1653   sopt.ToUpper();
1654   
1655   if (sopt.Contains("ALL"))
1656   {
1657     if ( fOutputEventHandler )
1658     {
1659       cout << TString('_',78) << endl;
1660       cout << "OutputEventHandler:" << endl;
1661       fOutputEventHandler->Print("   ");
1662     }
1663   }
1664 }
1665
1666 //______________________________________________________________________________
1667 void AliAnalysisManager::ResetAnalysis()
1668 {
1669 // Reset all execution flags and clean containers.
1670    CleanContainers();
1671 }
1672
1673 //______________________________________________________________________________
1674 void AliAnalysisManager::RunLocalInit()
1675 {
1676 // Run LocalInit method for all tasks.
1677    TDirectory *cdir = gDirectory;
1678    if (IsTrainInitialized()) return;
1679    TIter nextTask(fTasks);
1680    AliAnalysisTask *task;
1681    while ((task=(AliAnalysisTask*)nextTask())) {
1682       gROOT->cd();
1683       task->LocalInit();
1684    }
1685    if (cdir) cdir->cd();
1686    TObject::SetBit(kTasksInitialized, kTRUE);
1687 }   
1688
1689 //______________________________________________________________________________
1690 Long64_t AliAnalysisManager::StartAnalysis(const char *type, Long64_t nentries, Long64_t firstentry)
1691 {
1692 // Start analysis having a grid handler.
1693    if (!fGridHandler) {
1694       Error("StartAnalysis", "Cannot start analysis providing just the analysis type without a grid handler.");
1695       Info("===", "Add an AliAnalysisAlien object as plugin for this manager and configure it.");
1696       return -1;
1697    }
1698    TTree *tree = NULL;
1699    return StartAnalysis(type, tree, nentries, firstentry);
1700 }
1701
1702 //______________________________________________________________________________
1703 Long64_t AliAnalysisManager::StartAnalysis(const char *type, TTree * const tree, Long64_t nentries, Long64_t firstentry)
1704 {
1705 // Start analysis for this manager. Analysis task can be: LOCAL, PROOF, GRID or
1706 // MIX. Process nentries starting from firstentry
1707    Long64_t retv = 0;
1708    // Backup current directory and make sure gDirectory points to gROOT
1709    TDirectory *cdir = gDirectory;
1710    gROOT->cd();
1711    if (!fInitOK) {
1712       Error("StartAnalysis","Analysis manager was not initialized !");
1713       if (cdir) cdir->cd();
1714       return -1;
1715    }
1716    if (!CheckTasks()) Fatal("StartAnalysis", "Not all needed libraries were loaded");
1717    if (fDebug > 1) {
1718       printf("StartAnalysis %s\n",GetName());
1719       AliLog::SetGlobalLogLevel(AliLog::kInfo);
1720    }
1721    fMaxEntries = nentries;
1722    fIsRemote = kFALSE;
1723    TString anaType = type;
1724    anaType.ToLower();
1725    fMode = kLocalAnalysis;
1726    if (anaType.Contains("file"))      fIsRemote = kTRUE;
1727    if (anaType.Contains("proof"))     fMode = kProofAnalysis;
1728    else if (anaType.Contains("grid")) fMode = kGridAnalysis;
1729    else if (anaType.Contains("mix"))  fMode = kMixingAnalysis;
1730
1731    if (fMode == kGridAnalysis) {
1732       fIsRemote = kTRUE;
1733       if (!anaType.Contains("terminate")) {
1734          if (!fGridHandler) {
1735             Error("StartAnalysis", "Cannot start grid analysis without a grid handler.");
1736             Info("===", "Add an AliAnalysisAlien object as plugin for this manager and configure it.");
1737             if (cdir) cdir->cd();
1738             return -1;
1739          }
1740          // Write analysis manager in the analysis file
1741          cout << "===== RUNNING GRID ANALYSIS: " << GetName() << endl;
1742          // run local task configuration
1743          RunLocalInit();
1744          if (!fGridHandler->StartAnalysis(nentries, firstentry)) {
1745             Info("StartAnalysis", "Grid analysis was stopped and cannot be terminated");
1746             if (cdir) cdir->cd();
1747             return -1;
1748          }   
1749
1750          // Terminate grid analysis
1751          if (fSelector && fSelector->GetStatus() == -1) {if (cdir) cdir->cd(); return -1;}
1752          if (fGridHandler->GetRunMode() == AliAnalysisGrid::kOffline) {if (cdir) cdir->cd(); return 0;}
1753          cout << "===== MERGING OUTPUTS REGISTERED BY YOUR ANALYSIS JOB: " << GetName() << endl;
1754          if (!fGridHandler->MergeOutputs()) {
1755             // Return if outputs could not be merged or if it alien handler
1756             // was configured for offline mode or local testing.
1757             if (cdir) cdir->cd();
1758             return 0;
1759          }
1760       }   
1761       cout << "===== TERMINATING GRID ANALYSIS JOB: " << GetName() << endl;
1762       if (cdir) cdir->cd();
1763       ImportWrappers(NULL);
1764       Terminate();
1765       if (cdir) cdir->cd();
1766       return 0;
1767    }
1768    TString line;
1769    SetEventLoop(kFALSE);
1770    // Enable event loop mode if a tree was provided
1771    if (tree || fGridHandler || fMode==kMixingAnalysis) SetEventLoop(kTRUE);
1772
1773    TChain *chain = 0;
1774    TString ttype = "TTree";
1775    if (tree && tree->IsA() == TChain::Class()) {
1776       chain = (TChain*)tree;
1777       if (!chain || !chain->GetListOfFiles()->First()) {
1778          Error("StartAnalysis", "Cannot process null or empty chain...");
1779          if (cdir) cdir->cd();
1780          return -1;
1781       }   
1782       ttype = "TChain";
1783    }   
1784
1785    Bool_t getsysInfo = ((fNSysInfo>0) && (fMode==kLocalAnalysis))?kTRUE:kFALSE;
1786    if (getsysInfo) AliSysInfo::AddStamp("Start", 0);
1787    // Initialize locally all tasks (happens for all modes)
1788    TIter next(fTasks);
1789    AliAnalysisTask *task;
1790    RunLocalInit();
1791    
1792    switch (fMode) {
1793       case kLocalAnalysis:
1794          if (!tree && !fGridHandler) {
1795             TIter nextT(fTasks);
1796             // Call CreateOutputObjects for all tasks
1797             Int_t itask = 0;
1798             Bool_t dirStatus = TH1::AddDirectoryStatus();
1799             while ((task=(AliAnalysisTask*)nextT())) {
1800                TH1::AddDirectory(kFALSE);
1801                task->CreateOutputObjects();
1802                if (!task->CheckPostData()) {
1803                   Error("SlaveBegin","####### IMPORTANT! ####### \n\n\n\
1804                         Task %s (%s) did not call PostData() for all its outputs in (User)CreateOutputObjects()\n\n\
1805                         ########### FIX YOUR CODE, THIS WILL PRODUCE A FATAL ERROR IN FUTURE! ###########", task->GetName(), task->ClassName());
1806                }
1807                if (getsysInfo) AliSysInfo::AddStamp(Form("%s_CREATEOUTOBJ",task->ClassName()), 0, itask, 0);
1808                gROOT->cd();
1809                itask++;
1810             }   
1811             TH1::AddDirectory(dirStatus);
1812             if (IsExternalLoop()) {
1813                Info("StartAnalysis", "Initialization done. Event loop is controlled externally.\
1814                      \nSetData for top container, call ExecAnalysis in a loop and then Terminate manually");
1815                return 0;
1816             }
1817             ExecAnalysis();
1818             Terminate();
1819             return 0;
1820          } 
1821          fSelector = new AliAnalysisSelector(this);
1822          // Check if a plugin handler is used
1823          if (fGridHandler) {
1824             // Get the chain from the plugin
1825             TString dataType = "esdTree";
1826             if (fInputEventHandler) {
1827                dataType = fInputEventHandler->GetDataType();
1828                dataType.ToLower();
1829                dataType += "Tree";
1830             }   
1831             chain = fGridHandler->GetChainForTestMode(dataType);
1832             if (!chain) {
1833                Error("StartAnalysis", "No chain for test mode. Aborting.");
1834                return -1;
1835             }
1836             cout << "===== RUNNING LOCAL ANALYSIS" << GetName() << " ON CHAIN " << chain->GetName() << endl;
1837             retv = chain->Process(fSelector, "", nentries, firstentry);
1838             break;
1839          }
1840          // Run tree-based analysis via AliAnalysisSelector  
1841          cout << "===== RUNNING LOCAL ANALYSIS " << GetName() << " ON TREE " << tree->GetName() << endl;
1842          retv = tree->Process(fSelector, "", nentries, firstentry);
1843          break;
1844       case kProofAnalysis:
1845          fIsRemote = kTRUE;
1846          // Check if the plugin is used
1847          if (fGridHandler) {
1848             return StartAnalysis(type, fGridHandler->GetProofDataSet(), nentries, firstentry);
1849          }
1850          if (!gROOT->GetListOfProofs() || !gROOT->GetListOfProofs()->GetEntries()) {
1851             Error("StartAnalysis", "No PROOF!!! Exiting.");
1852             if (cdir) cdir->cd();
1853             return -1;
1854          }   
1855          line = Form("gProof->AddInput((TObject*)%p);", this);
1856          gROOT->ProcessLine(line);
1857          if (chain) {
1858             chain->SetProof();
1859             cout << "===== RUNNING PROOF ANALYSIS " << GetName() << " ON CHAIN " << chain->GetName() << endl;
1860             retv = chain->Process("AliAnalysisSelector", "", nentries, firstentry);
1861          } else {
1862             Error("StartAnalysis", "No chain!!! Exiting.");
1863             if (cdir) cdir->cd();
1864             return -1;
1865          }      
1866          break;
1867       case kGridAnalysis:
1868          fIsRemote = kTRUE;
1869          if (!anaType.Contains("terminate")) {
1870             if (!fGridHandler) {
1871                Error("StartAnalysis", "Cannot start grid analysis without a grid handler.");
1872                Info("===", "Add an AliAnalysisAlien object as plugin for this manager and configure it.");
1873                if (cdir) cdir->cd();
1874                return -1;
1875             }
1876             // Write analysis manager in the analysis file
1877             cout << "===== RUNNING GRID ANALYSIS: " << GetName() << endl;
1878             // Start the analysis via the handler
1879             if (!fGridHandler->StartAnalysis(nentries, firstentry)) {
1880                Info("StartAnalysis", "Grid analysis was stopped and cannot be terminated");
1881                if (cdir) cdir->cd();
1882                return -1;
1883             }   
1884
1885             // Terminate grid analysis
1886             if (fSelector && fSelector->GetStatus() == -1) {if (cdir) cdir->cd(); return -1;}
1887             if (fGridHandler->GetRunMode() == AliAnalysisGrid::kOffline) {if (cdir) cdir->cd(); return 0;}
1888             cout << "===== MERGING OUTPUTS REGISTERED BY YOUR ANALYSIS JOB: " << GetName() << endl;
1889             if (!fGridHandler->MergeOutputs()) {
1890                // Return if outputs could not be merged or if it alien handler
1891                // was configured for offline mode or local testing.
1892                if (cdir) cdir->cd();
1893                return 0;
1894             }
1895          }   
1896          cout << "===== TERMINATING GRID ANALYSIS JOB: " << GetName() << endl;
1897          ImportWrappers(NULL);
1898          Terminate();
1899          if (cdir) cdir->cd();
1900          return 0;
1901       case kMixingAnalysis:   
1902          // Run event mixing analysis
1903          if (!fEventPool) {
1904             Error("StartAnalysis", "Cannot run event mixing without event pool");
1905             if (cdir) cdir->cd();
1906             return -1;
1907          }
1908          cout << "===== RUNNING EVENT MIXING ANALYSIS " << GetName() << endl;
1909          fSelector = new AliAnalysisSelector(this);
1910          while ((chain=fEventPool->GetNextChain())) {
1911             next.Reset();
1912             // Call NotifyBinChange for all tasks
1913             while ((task=(AliAnalysisTask*)next()))
1914                if (!task->IsPostEventLoop()) task->NotifyBinChange();
1915             retv = chain->Process(fSelector);
1916             if (retv < 0) {
1917                Error("StartAnalysis", "Mixing analysis failed");
1918                if (cdir) cdir->cd();
1919                return retv;
1920             }   
1921          }
1922          PackOutput(fSelector->GetOutputList());
1923          Terminate();
1924    }
1925    if (cdir) cdir->cd();
1926    return retv;
1927 }   
1928
1929 //______________________________________________________________________________
1930 Long64_t AliAnalysisManager::StartAnalysis(const char *type, const char *dataset, Long64_t nentries, Long64_t firstentry)
1931 {
1932 // Start analysis for this manager on a given dataset. Analysis task can be: 
1933 // LOCAL, PROOF or GRID. Process nentries starting from firstentry.
1934    if (!fInitOK) {
1935       Error("StartAnalysis","Analysis manager was not initialized !");
1936       return -1;
1937    }
1938    fIsRemote = kTRUE;
1939    if (fDebug > 1) printf("StartAnalysis %s\n",GetName());
1940    TString anaType = type;
1941    anaType.ToLower();
1942    if (!anaType.Contains("proof")) {
1943       Error("StartAnalysis", "Cannot process datasets in %s mode. Try PROOF.", type);
1944       return -1;
1945    }   
1946    fMode = kProofAnalysis;
1947    TString line;
1948    SetEventLoop(kTRUE);
1949    // Set the dataset flag
1950    TObject::SetBit(kUseDataSet);
1951    fTree = 0;
1952    if (fGridHandler) {
1953       // Start proof analysis using the grid handler
1954       if (!fGridHandler->StartAnalysis(nentries, firstentry)) {
1955          Error("StartAnalysis", "The grid plugin could not start PROOF analysis");
1956          return -1;
1957       }
1958       // Check if the plugin is in test mode
1959       if (fGridHandler->GetRunMode() == AliAnalysisGrid::kTest) {
1960          dataset = "test_collection";
1961       } else {
1962          dataset = fGridHandler->GetProofDataSet();
1963       }
1964    }   
1965
1966    if (!gROOT->GetListOfProofs() || !gROOT->GetListOfProofs()->GetEntries()) {
1967       Error("StartAnalysis", "No PROOF!!! Exiting.");
1968       return -1;
1969    }   
1970
1971    // Initialize locally all tasks
1972    RunLocalInit();
1973       
1974    line = Form("gProof->AddInput((TObject*)%p);", this);
1975    gROOT->ProcessLine(line);
1976    Long_t retv;
1977    line = Form("gProof->Process(\"%s\", \"AliAnalysisSelector\", \"\", %lld, %lld);",
1978                dataset, nentries, firstentry);
1979    cout << "===== RUNNING PROOF ANALYSIS " << GetName() << " ON DATASET " << dataset << endl;
1980    retv = (Long_t)gROOT->ProcessLine(line);
1981    return retv;
1982 }   
1983
1984 //______________________________________________________________________________
1985 TFile *AliAnalysisManager::OpenFile(AliAnalysisDataContainer *cont, const char *option, Bool_t ignoreProof)
1986 {
1987 // Opens according the option the file specified by cont->GetFileName() and changes
1988 // current directory to cont->GetFolderName(). If the file was already opened, it
1989 // checks if the option UPDATE was preserved. File open via TProofOutputFile can
1990 // be optionally ignored.
1991   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1992   TString filename = cont->GetFileName();
1993   TFile *f = NULL;
1994   if (filename.IsNull()) {
1995     ::Error("AliAnalysisManager::OpenFile", "No file name specified for container %s", cont->GetName());
1996     return NULL;
1997   }
1998   if (mgr->GetAnalysisType()==AliAnalysisManager::kProofAnalysis && cont->IsSpecialOutput()
1999       && !ignoreProof)
2000     f = mgr->OpenProofFile(cont,option);
2001   else {
2002     // Check first if the file is already opened
2003     f = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
2004     if (f) {
2005       // Check if option "UPDATE" was preserved 
2006       TString opt(option);
2007       opt.ToUpper();
2008       if ((opt=="UPDATE") && (opt!=f->GetOption())) 
2009         ::Info("AliAnalysisManager::OpenFile", "File %s already opened in %s mode!", cont->GetFileName(), f->GetOption());
2010     } else {
2011       f = TFile::Open(filename, option);
2012     }    
2013   }   
2014   if (f && !f->IsZombie() && !f->TestBit(TFile::kRecovered)) {
2015     cont->SetFile(f);
2016     // Cd to file
2017     f->cd();
2018     // Check for a folder request
2019     TString dir = cont->GetFolderName(); 
2020     if (!dir.IsNull()) {
2021       if (!f->GetDirectory(dir)) f->mkdir(dir);
2022       f->cd(dir);
2023     }
2024     return f;
2025   }
2026   ::Fatal("AliAnalysisManager::OpenFile", "File %s could not be opened", filename.Data());
2027   cont->SetFile(NULL);
2028   return NULL;
2029 }    
2030  
2031 //______________________________________________________________________________
2032 TFile *AliAnalysisManager::OpenProofFile(AliAnalysisDataContainer *cont, const char *option, const char *extaod)
2033 {
2034 // Opens a special output file used in PROOF.
2035   TString line;
2036   TString filename = cont->GetFileName();
2037   if (cont == fCommonOutput) {
2038      if (fOutputEventHandler) {
2039         if (strlen(extaod)) filename = extaod;
2040         filename = fOutputEventHandler->GetOutputFileName();
2041      }   
2042      else Fatal("OpenProofFile","No output container. Exiting.");
2043   }   
2044   TFile *f = NULL;
2045   if (fMode!=kProofAnalysis || !fSelector) {
2046     Fatal("OpenProofFile","Cannot open PROOF file %s: no PROOF or selector",filename.Data());
2047     return NULL;
2048   } 
2049   if (fSpecialOutputLocation.Length()) {
2050     f = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
2051     if (f) {
2052       // Check if option "UPDATE" was preserved 
2053       TString opt(option);
2054       opt.ToUpper();
2055       if ((opt=="UPDATE") && (opt!=f->GetOption()))
2056         ::Info("OpenProofFile", "File %s already opened in %s mode!", cont->GetFileName(), f->GetOption());
2057     } else {
2058       f = new TFile(filename, option);
2059     }
2060     if (f && !f->IsZombie() && !f->TestBit(TFile::kRecovered)) {
2061       cont->SetFile(f);
2062       // Cd to file
2063       f->cd();
2064       // Check for a folder request
2065       TString dir = cont->GetFolderName(); 
2066       if (dir.Length()) {
2067         if (!f->GetDirectory(dir)) f->mkdir(dir);
2068         f->cd(dir);
2069       }      
2070       return f;
2071     }
2072     Fatal("OpenProofFile", "File %s could not be opened", cont->GetFileName());
2073     cont->SetFile(NULL);
2074     return NULL;       
2075   }
2076   // Check if there is already a proof output file in the output list
2077   TObject *pof = fSelector->GetOutputList()->FindObject(filename);
2078   if (pof) {
2079     // Get the actual file
2080     line = Form("((TProofOutputFile*)%p)->GetFileName();", pof);
2081     filename = (const char*)gROOT->ProcessLine(line);
2082     if (fDebug>1) {
2083       printf("File: %s already booked via TProofOutputFile\n", filename.Data());
2084     }  
2085     f = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
2086     if (!f) {
2087        Fatal("OpenProofFile", "Proof output file found but no file opened for %s", filename.Data());
2088        return NULL;
2089     }   
2090     // Check if option "UPDATE" was preserved 
2091     TString opt(option);
2092     opt.ToUpper();
2093     if ((opt=="UPDATE") && (opt!=f->GetOption())) 
2094       Fatal("OpenProofFile", "File %s already opened, but not in UPDATE mode!", cont->GetFileName());
2095   } else {
2096     if (cont->IsRegisterDataset()) {
2097       TString dsetName = filename;
2098       dsetName.ReplaceAll(".root", cont->GetTitle());
2099       dsetName.ReplaceAll(":","_");
2100       if (fDebug>1) printf("Booking dataset: %s\n", dsetName.Data());
2101       line = Form("TProofOutputFile *pf = new TProofOutputFile(\"%s\", \"DROV\", \"%s\");", filename.Data(), dsetName.Data());
2102     } else {
2103       if (fDebug>1) printf("Booking TProofOutputFile: %s to be merged\n", filename.Data());
2104       line = Form("TProofOutputFile *pf = new TProofOutputFile(\"%s\");", filename.Data());
2105     }
2106     if (fDebug > 1) printf("=== %s\n", line.Data());
2107     gROOT->ProcessLine(line);
2108     line = Form("pf->OpenFile(\"%s\");", option);
2109     gROOT->ProcessLine(line);
2110     f = gFile;
2111     if (fDebug > 1) {
2112       gROOT->ProcessLine("pf->Print()");
2113       printf(" == proof file name: %s", f->GetName());
2114     }   
2115     // Add to proof output list
2116     line = Form("((TList*)%p)->Add(pf);",fSelector->GetOutputList());
2117     if (fDebug > 1) printf("=== %s\n", line.Data());
2118     gROOT->ProcessLine(line);
2119   }
2120   if (f && !f->IsZombie() && !f->TestBit(TFile::kRecovered)) {
2121     cont->SetFile(f);
2122     // Cd to file
2123     f->cd();
2124     // Check for a folder request
2125     TString dir = cont->GetFolderName(); 
2126     if (!dir.IsNull()) {
2127       if (!f->GetDirectory(dir)) f->mkdir(dir);
2128       f->cd(dir);
2129     }
2130     return f;
2131   }
2132   Fatal("OpenProofFile", "File %s could not be opened", cont->GetFileName());
2133   cont->SetFile(NULL);  
2134   return NULL;
2135 }   
2136
2137 //______________________________________________________________________________
2138 void AliAnalysisManager::ExecAnalysis(Option_t *option)
2139 {
2140 // Execute analysis.
2141    static Long64_t nentries = 0;
2142    static TTree *lastTree = 0;
2143    static TStopwatch *timer = new TStopwatch();
2144    // Only the first call to Process will trigger a true Notify. Other Notify
2145    // coming before is ignored.
2146    if (!TObject::TestBit(AliAnalysisManager::kTrueNotify)) {
2147       TObject::SetBit(AliAnalysisManager::kTrueNotify);
2148       Notify();
2149    }   
2150    if (fDebug > 0) printf("MGR: Processing event #%d\n", fNcalls);
2151    else {
2152       if (fTree && (fTree != lastTree)) {
2153          nentries += fTree->GetEntries();
2154          lastTree = fTree;
2155       }   
2156       if (!fNcalls) timer->Start();
2157       if (!fIsRemote && TObject::TestBit(kUseProgressBar)) ProgressBar("Processing event", fNcalls, TMath::Min(fMaxEntries,nentries), timer, kFALSE);
2158    }
2159    fIOTimer->Start(kTRUE);
2160    gROOT->cd();
2161    TDirectory *cdir = gDirectory;
2162    Bool_t getsysInfo = ((fNSysInfo>0) && (fMode==kLocalAnalysis))?kTRUE:kFALSE;
2163    if (getsysInfo && ((fNcalls%fNSysInfo)==0)) AliSysInfo::AddStamp("Exec_start", (Int_t)fNcalls);
2164    if (!fInitOK) {
2165       Error("ExecAnalysis", "Analysis manager was not initialized !");
2166       if (cdir) cdir->cd();
2167       return;
2168    }
2169    fNcalls++;
2170    AliAnalysisTask *task;
2171    // Check if the top tree is active.
2172    if (fTree) {
2173       if (getsysInfo && ((fNcalls%fNSysInfo)==0)) 
2174          AliSysInfo::AddStamp("Handlers_BeginEventGroup",fNcalls, 1002, 0);
2175       TIter next(fTasks);
2176    // De-activate all tasks
2177       while ((task=(AliAnalysisTask*)next())) task->SetActive(kFALSE);
2178       AliAnalysisDataContainer *cont = fCommonInput;
2179       if (!cont) cont = (AliAnalysisDataContainer*)fInputs->At(0);
2180       if (!cont) {
2181               Error("ExecAnalysis","Cannot execute analysis in TSelector mode without at least one top container");
2182          if (cdir) cdir->cd();
2183          return;
2184       }   
2185       cont->SetData(fTree); // This will notify all consumers
2186       Long64_t entry = fTree->GetTree()->GetReadEntry();      
2187 //
2188 //    Call BeginEvent() for optional input/output and MC services 
2189       if (fInputEventHandler)   fInputEventHandler  ->BeginEvent(entry);
2190       if (fOutputEventHandler)  fOutputEventHandler ->BeginEvent(entry);
2191       if (fMCtruthEventHandler) fMCtruthEventHandler->BeginEvent(entry);
2192       gROOT->cd();
2193       if (getsysInfo && ((fNcalls%fNSysInfo)==0)) 
2194          AliSysInfo::AddStamp("Handlers_BeginEvent",fNcalls, 1000, 0);
2195 //
2196 //    Execute the tasks
2197 //      TIter next1(cont->GetConsumers());
2198       fIOTimer->Stop();
2199       fIOTime += fIOTimer->RealTime();
2200       fCPUTimer->Start(kTRUE);
2201       TIter next1(fTopTasks);
2202       Int_t itask = 0;
2203       while ((task=(AliAnalysisTask*)next1())) {
2204          if (fDebug >1) {
2205             cout << "    Executing task " << task->GetName() << endl;
2206          }       
2207          task->ExecuteTask(option);
2208          gROOT->cd();
2209          if (getsysInfo && ((fNcalls%fNSysInfo)==0)) 
2210             AliSysInfo::AddStamp(task->ClassName(), fNcalls, itask, 1);
2211          itask++;   
2212       }
2213       fCPUTimer->Stop();
2214       fCPUTime += fCPUTimer->RealTime();
2215       fIOTimer->Start(kTRUE);
2216 //
2217 //    Call FinishEvent() for optional output and MC services 
2218       if (fInputEventHandler)   fInputEventHandler  ->FinishEvent();
2219       if (fOutputEventHandler)  fOutputEventHandler ->FinishEvent();
2220       if (fMCtruthEventHandler) fMCtruthEventHandler->FinishEvent();
2221       // Gather system information if requested
2222       if (getsysInfo && ((fNcalls%fNSysInfo)==0)) 
2223          AliSysInfo::AddStamp("Handlers_FinishEvent",fNcalls, 1001, 1);
2224       if (cdir) cdir->cd();   
2225       fIOTimer->Stop();
2226       fIOTime += fIOTimer->RealTime();
2227       return;
2228    }   
2229    // The event loop is not controlled by TSelector   
2230 //
2231 //  Call BeginEvent() for optional input/output and MC services 
2232    fIOTimer->Start(kTRUE);
2233    if (fInputEventHandler)   fInputEventHandler  ->BeginEvent(-1);
2234    if (fOutputEventHandler)  fOutputEventHandler ->BeginEvent(-1);
2235    if (fMCtruthEventHandler) fMCtruthEventHandler->BeginEvent(-1);
2236    fIOTimer->Stop();
2237    fIOTime += fIOTimer->RealTime();
2238    gROOT->cd();
2239    if (getsysInfo && ((fNcalls%fNSysInfo)==0)) 
2240       AliSysInfo::AddStamp("Handlers_BeginEvent",fNcalls, 1000, 0);
2241    fCPUTimer->Start(kTRUE);
2242    TIter next2(fTopTasks);
2243    while ((task=(AliAnalysisTask*)next2())) {
2244       task->SetActive(kTRUE);
2245       if (fDebug > 1) {
2246          cout << "    Executing task " << task->GetName() << endl;
2247       }   
2248       task->ExecuteTask(option);
2249       gROOT->cd();
2250    }   
2251    fCPUTimer->Stop();
2252    fCPUTime += fCPUTimer->RealTime();
2253 //
2254 // Call FinishEvent() for optional output and MC services 
2255    fIOTimer->Start(kTRUE);
2256    if (fInputEventHandler)   fInputEventHandler  ->FinishEvent();
2257    if (fOutputEventHandler)  fOutputEventHandler ->FinishEvent();
2258    if (fMCtruthEventHandler) fMCtruthEventHandler->FinishEvent();
2259    if (getsysInfo && ((fNcalls%fNSysInfo)==0)) 
2260       AliSysInfo::AddStamp("Handlers_FinishEvent",fNcalls, 1000, 1);
2261    if (cdir) cdir->cd();   
2262    fIOTimer->Stop();
2263    fIOTime += fIOTimer->RealTime();
2264 }
2265
2266 //______________________________________________________________________________
2267 Bool_t AliAnalysisManager::IsPipe(std::ostream &out)
2268 {
2269 // Check if the stdout is connected to a pipe (C.Holm)
2270   Bool_t ispipe = kFALSE;
2271   out.seekp(0, std::ios_base::cur);
2272   if (out.fail()) {
2273     out.clear();
2274     if (errno == ESPIPE) ispipe = kTRUE;
2275   }
2276   return ispipe;
2277 }
2278    
2279 //______________________________________________________________________________
2280 void AliAnalysisManager::SetInputEventHandler(AliVEventHandler* const handler)
2281 {
2282 // Set the input event handler and create a container for it.
2283    Changed();
2284    fInputEventHandler   = handler;
2285    if (!fCommonInput) fCommonInput = CreateContainer("cAUTO_INPUT", TChain::Class(), AliAnalysisManager::kInputContainer);
2286 }
2287
2288 //______________________________________________________________________________
2289 void AliAnalysisManager::SetOutputEventHandler(AliVEventHandler* const handler)
2290 {
2291 // Set the input event handler and create a container for it.
2292    Changed();
2293    fOutputEventHandler   = handler;
2294    if (!fCommonOutput) fCommonOutput = CreateContainer("cAUTO_OUTPUT", TTree::Class(), AliAnalysisManager::kOutputContainer, "default");
2295    fCommonOutput->SetSpecialOutput();
2296 }
2297
2298 //______________________________________________________________________________
2299 void AliAnalysisManager::SetDebugLevel(UInt_t level)
2300 {
2301 // Set verbosity of the analysis manager. If the progress bar is used, the call is ignored
2302    if (TObject::TestBit(kUseProgressBar)) {
2303       Info("SetDebugLevel","Ignored. Disable the progress bar first.");
2304       return;
2305    }
2306    fDebug = level;
2307 }
2308    
2309 //______________________________________________________________________________
2310 void AliAnalysisManager::SetUseProgressBar(Bool_t flag, Int_t freq)
2311 {
2312 // Enable a text mode progress bar. Resets debug level to 0.
2313    Info("SetUseProgressBar", "Progress bar enabled, updated every %d events.\n  ### NOTE: Debug level reset to 0 ###", freq);
2314    TObject::SetBit(kUseProgressBar,flag);
2315    fPBUpdateFreq = freq;
2316    fDebug = 0;
2317 }   
2318
2319 //______________________________________________________________________________
2320 void AliAnalysisManager::RegisterExtraFile(const char *fname)
2321 {
2322 // This method is used externally to register output files which are not
2323 // connected to any output container, so that the manager can properly register,
2324 // retrieve or merge them when running in distributed mode. The file names are
2325 // separated by blancs. The method has to be called in MyAnalysisTask::LocalInit().
2326    if (fExtraFiles.Contains(fname)) return;
2327    if (fExtraFiles.Length()) fExtraFiles += " ";
2328    fExtraFiles += fname;
2329 }
2330
2331 //______________________________________________________________________________
2332 Bool_t AliAnalysisManager::GetFileFromWrapper(const char *filename, const TList *source)
2333 {
2334 // Copy a file from the location specified ina the wrapper with the same name from the source list.
2335    char fullPath[512];
2336    char chUrl[512];
2337    char tmp[1024];
2338    TObject *pof =  source->FindObject(filename);
2339    if (!pof || !pof->InheritsFrom("TProofOutputFile")) {
2340       Error("GetFileFromWrapper", "TProofOutputFile object not found in output list for file %s", filename);
2341       return kFALSE;
2342    }
2343    gROOT->ProcessLine(Form("sprintf((char*)%p, \"%%s\", ((TProofOutputFile*)%p)->GetOutputFileName());", fullPath, pof));
2344    gROOT->ProcessLine(Form("sprintf((char*)%p, \"%%s\", gProof->GetUrl());",chUrl));
2345    TString clientUrl(chUrl);
2346    TString fullPath_str(fullPath);
2347    if (clientUrl.Contains("localhost")){
2348       TObjArray* array = fullPath_str.Tokenize ( "//" );
2349       TObjString *strobj = ( TObjString *)array->At(1);
2350       TObjArray* arrayPort = strobj->GetString().Tokenize ( ":" );
2351       TObjString *strobjPort = ( TObjString *) arrayPort->At(1);
2352       fullPath_str.ReplaceAll(strobj->GetString().Data(),"localhost:PORT");
2353       fullPath_str.ReplaceAll(":PORT",Form(":%s",strobjPort->GetString().Data()));
2354       if (fDebug > 1) Info("GetFileFromWrapper","Using tunnel from %s to %s",fullPath_str.Data(),filename);
2355       delete arrayPort;
2356       delete array;
2357    }
2358    else if (clientUrl.Contains("__lite__")) { 
2359      // Special case for ProofLite environement - get file info and copy. 
2360      gROOT->ProcessLine(Form("sprintf((char*)%p,\"%%s\",((TProofOutputFile*)%p)->GetDir());", tmp, pof));
2361      fullPath_str = Form("%s/%s", tmp, fullPath);
2362    }
2363    if (fDebug > 1) 
2364      Info("GetFileFromWrapper","Copying file %s from PROOF scratch space to %s", fullPath_str.Data(),filename);
2365    Bool_t gotit = TFile::Cp(fullPath_str.Data(), filename); 
2366    if (!gotit)
2367       Error("GetFileFromWrapper", "Could not get file %s from proof scratch space", filename);
2368    return gotit;
2369 }
2370
2371 //______________________________________________________________________________
2372 void AliAnalysisManager::GetAnalysisTypeString(TString &type) const
2373 {
2374 // Fill analysis type in the provided string.
2375    switch (fMode) {
2376       case kLocalAnalysis:
2377          type = "local";
2378          return;
2379       case kProofAnalysis:
2380          type = "proof";
2381          return;
2382       case kGridAnalysis:
2383          type = "grid";
2384          return;
2385       case kMixingAnalysis:
2386          type = "mix";
2387    }
2388 }
2389
2390 //______________________________________________________________________________
2391 Bool_t AliAnalysisManager::ValidateOutputFiles() const
2392 {
2393 // Validate all output files.
2394    TIter next(fOutputs);
2395    AliAnalysisDataContainer *output;
2396    TDirectory *cdir = gDirectory;
2397    TString openedFiles;
2398    while ((output=(AliAnalysisDataContainer*)next())) {
2399       if (output->IsRegisterDataset()) continue;
2400       TString filename = output->GetFileName();
2401       if (filename == "default") {
2402          if (!fOutputEventHandler) continue;
2403          filename = fOutputEventHandler->GetOutputFileName();
2404          // Main AOD may not be there
2405          if (gSystem->AccessPathName(filename)) continue;
2406       }
2407       // Check if the file is closed
2408       if (openedFiles.Contains(filename)) continue;;
2409       TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
2410       if (file) {
2411          Warning("ValidateOutputs", "File %s was not closed. Closing.", filename.Data());
2412          // Clear file list to release object ownership to user.
2413 //         file->Clear();
2414          file->Close();
2415       }
2416       file = TFile::Open(filename);
2417       if (!file || file->IsZombie() || file->TestBit(TFile::kRecovered)) {
2418          Error("ValidateOutputs", "Output file <%s> was not created or invalid", filename.Data());
2419          if (cdir) cdir->cd();
2420          return kFALSE;
2421       }
2422       file->Close();
2423       openedFiles += filename;
2424       openedFiles += " ";
2425    }
2426    if (cdir) cdir->cd();
2427    return kTRUE;
2428 }   
2429
2430 //______________________________________________________________________________
2431 void AliAnalysisManager::ProgressBar(const char *opname, Long64_t current, Long64_t size, TStopwatch * const watch, Bool_t last, Bool_t refresh)
2432 {
2433 // Implements a nice text mode progress bar.
2434    static Long64_t icount = 0;
2435    static TString oname;
2436    static TString nname;
2437    static Long64_t ocurrent = 0;
2438    static Long64_t osize = 0;
2439    static Int_t oseconds = 0;
2440    static TStopwatch *owatch = 0;
2441    static Bool_t oneoftwo = kFALSE;
2442    static Int_t nrefresh = 0;
2443    static Int_t nchecks = 0;
2444    static char lastChar = 0;
2445    const char symbol[4] = {'-','\\','|','/'}; 
2446    
2447    if (!lastChar) lastChar = (IsPipe(std::cerr))?'\r':'\n';
2448    if (!refresh) {
2449       nrefresh = 0;
2450       if (!size) return;
2451       owatch = watch;
2452       oname = opname;
2453       ocurrent = TMath::Abs(current);
2454       osize = TMath::Abs(size);
2455       if (ocurrent > osize) ocurrent=osize;
2456    } else {
2457       nrefresh++;
2458       if (!osize) return;
2459    }     
2460    if ((current % fPBUpdateFreq) != 0) return;
2461    icount++;
2462    char progress[11] = "          ";
2463    Int_t ichar = icount%4;
2464    Double_t time = 0.;
2465    Int_t hours = 0;
2466    Int_t minutes = 0;
2467    Int_t seconds = 0;
2468    if (owatch && !last) {
2469       owatch->Stop();
2470       time = owatch->RealTime();
2471       seconds   = int(time) % 60;
2472       minutes   = (int(time) / 60) % 60;
2473       hours     = (int(time) / 60 / 60);
2474       if (refresh)  {
2475          if (oseconds==seconds) {
2476             owatch->Continue();
2477             return;
2478          }
2479          oneoftwo = !oneoftwo;   
2480       }
2481       oseconds = seconds;   
2482    }
2483    if (refresh && oneoftwo) {
2484       nname = oname;
2485       if (nchecks <= 0) nchecks = nrefresh+1;
2486       Int_t pctdone = (Int_t)(100.*nrefresh/nchecks);
2487       oname = Form("     == %d%% ==", pctdone);
2488    }         
2489    Double_t percent = 100.0*ocurrent/osize;
2490    Int_t nchar = Int_t(percent/10);
2491    if (nchar>10) nchar=10;
2492    Int_t i;
2493    for (i=0; i<nchar; i++)  progress[i] = '=';
2494    progress[nchar] = symbol[ichar];
2495    for (i=nchar+1; i<10; i++) progress[i] = ' ';
2496    progress[10] = '\0';
2497    oname += "                    ";
2498    oname.Remove(20);
2499    if(size<10000) fprintf(stderr, "%s [%10s] %4lld ", oname.Data(), progress, ocurrent);
2500    else if(size<100000) fprintf(stderr, "%s [%10s] %5lld ",oname.Data(), progress, ocurrent);
2501    else fprintf(stderr, "%s [%10s] %7lld ",oname.Data(), progress, ocurrent);
2502    if (time>0.) {
2503      Int_t full   = Int_t(ocurrent > 0 ? 
2504                           time * (float(osize)/ocurrent) + .5 : 
2505                           99*3600+59*60+59); 
2506      Int_t remain = Int_t(full - time);
2507      Int_t rsec   = remain % 60;
2508      Int_t rmin   = (remain / 60) % 60;
2509      Int_t rhour  = (remain / 60 / 60);
2510      fprintf(stderr, "[%6.2f %%]   TIME %.2d:%.2d:%.2d  ETA %.2d:%.2d:%.2d%c",
2511              percent, hours, minutes, seconds, rhour, rmin, rsec, lastChar);
2512    }
2513    else fprintf(stderr, "[%6.2f %%]%c", percent, lastChar);
2514    if (refresh && oneoftwo) oname = nname;
2515    if (owatch) owatch->Continue();
2516    if (last) {
2517       icount = 0;
2518       owatch = 0;
2519       ocurrent = 0;
2520       osize = 0;
2521       oseconds = 0;
2522       oneoftwo = kFALSE;
2523       nrefresh = 0;
2524       fprintf(stderr, "\n");
2525    }   
2526 }
2527
2528 //______________________________________________________________________________
2529 void AliAnalysisManager::DoLoadBranch(const char *name) 
2530 {
2531   // Get tree and load branch if needed.
2532   static Long64_t crtEntry = -100;
2533
2534   if (fAutoBranchHandling || !fTree)
2535     return;
2536
2537   TBranch *br = dynamic_cast<TBranch*>(fTable.FindObject(name));
2538   if (!br) {
2539     br = fTree->GetBranch(name);
2540     if (!br) {
2541       Error("DoLoadBranch", "Could not find branch %s",name);
2542       return;
2543     }
2544     fTable.Add(br);
2545   }
2546   if (br->GetReadEntry()==fCurrentEntry) return;
2547   Long64_t readbytes = br->GetEntry(GetCurrentEntry());
2548   if (readbytes<0) {
2549     Error("DoLoadBranch", "Could not load entry %lld from branch %s",GetCurrentEntry(), name);
2550     if (crtEntry != fCurrentEntry) {
2551       CountEvent(1,0,1,0);
2552       crtEntry = fCurrentEntry;
2553     }  
2554   } else {
2555     if (crtEntry != fCurrentEntry) {
2556       CountEvent(1,1,0,0);
2557       crtEntry = fCurrentEntry;
2558     }
2559   }
2560 }
2561
2562 //______________________________________________________________________________
2563 void AliAnalysisManager::AddStatisticsTask(UInt_t offlineMask)
2564 {
2565 // Add the statistics task to the manager.
2566   if (fStatistics) {
2567      Info("AddStatisticsTask", "Already added");
2568      return;
2569   }
2570   TString line = Form("AliAnalysisTaskStat::AddToManager(%u);", offlineMask);
2571   gROOT->ProcessLine(line);
2572 }  
2573
2574 //______________________________________________________________________________
2575 void AliAnalysisManager::CountEvent(Int_t ninput, Int_t nprocessed, Int_t nfailed, Int_t naccepted)
2576 {
2577 // Bookkeep current event;
2578    if (!fStatistics) return;
2579    fStatistics->AddInput(ninput);
2580    fStatistics->AddProcessed(nprocessed);
2581    fStatistics->AddFailed(nfailed);
2582    fStatistics->AddAccepted(naccepted);
2583 }   
2584
2585 //______________________________________________________________________________
2586 void AliAnalysisManager::AddStatisticsMsg(const char *line)
2587 {
2588 // Add a line in the statistics message. If available, the statistics message is written
2589 // at the end of the SlaveTerminate phase on workers AND at the end of Terminate
2590 // on the client.
2591    if (!strlen(line)) return;
2592    if (!fStatisticsMsg.IsNull()) fStatisticsMsg += "\n";
2593    fStatisticsMsg += line;
2594 }
2595
2596 //______________________________________________________________________________
2597 void AliAnalysisManager::WriteStatisticsMsg(Int_t)
2598 {
2599 // If fStatistics is present, write the file in the format ninput_nprocessed_nfailed_naccepted.stat
2600    static Bool_t done = kFALSE;
2601    if (done) return;
2602    done = kTRUE;
2603    if (!fStatistics) return;
2604    ofstream out;
2605    AddStatisticsMsg(Form("Number of input events:        %lld",fStatistics->GetNinput()));
2606    AddStatisticsMsg(Form("Number of processed events:    %lld",fStatistics->GetNprocessed()));      
2607    AddStatisticsMsg(Form("Number of failed events (I/O): %lld",fStatistics->GetNfailed()));
2608    AddStatisticsMsg(Form("Number of accepted events for mask %s: %lld", AliAnalysisStatistics::GetMaskAsString(fStatistics->GetOfflineMask()), fStatistics->GetNaccepted()));
2609    out.open(Form("%lld_%lld_%lld_%lld.stat",fStatistics->GetNinput(),
2610                  fStatistics->GetNprocessed(),fStatistics->GetNfailed(),
2611                  fStatistics->GetNaccepted()), ios::out);      
2612    out << fStatisticsMsg << endl;
2613    out.close();
2614 }
2615
2616 //______________________________________________________________________________
2617 const char* AliAnalysisManager::GetOADBPath()
2618 {
2619 // returns the path of the OADB
2620 // this static function just depends on environment variables
2621
2622    static TString oadbPath;
2623
2624    if (gSystem->Getenv("OADB_PATH"))
2625       oadbPath = gSystem->Getenv("OADB_PATH");
2626    else if (gSystem->Getenv("ALICE_ROOT"))
2627       oadbPath.Form("%s/OADB", gSystem->Getenv("ALICE_ROOT"));
2628    else
2629       ::Fatal("AliAnalysisManager::GetOADBPath", "Cannot figure out AODB path. Define ALICE_ROOT or OADB_PATH!");
2630       
2631    return oadbPath;
2632 }
2633
2634 //______________________________________________________________________________
2635 void AliAnalysisManager::SetGlobalStr(const char *key, const char *value)
2636 {
2637 // Define a custom string variable mapped to a global unique name. The variable
2638 // can be then retrieved by a given analysis macro via GetGlobalStr(key).
2639    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2640    if (!mgr) {
2641       ::Error("AliAnalysisManager::SetGlobalStr", "No analysis manager defined");
2642       return;
2643    }   
2644    Bool_t valid = kFALSE;
2645    TString existing = AliAnalysisManager::GetGlobalStr(key, valid);
2646    if (valid) {
2647       ::Error("AliAnalysisManager::SetGlobalStr", "Global %s = %s already defined.", key, existing.Data());
2648       return;
2649    }
2650    mgr->GetGlobals()->Add(new TObjString(key), new TObjString(value));
2651 }
2652
2653 //______________________________________________________________________________
2654 const char *AliAnalysisManager::GetGlobalStr(const char *key, Bool_t &valid)
2655 {
2656 // Static method to retrieve a global variable defined via SetGlobalStr.
2657    valid = kFALSE;
2658    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2659    if (!mgr) return 0;
2660    TObject *value = mgr->GetGlobals()->GetValue(key);
2661    if (!value) return 0;
2662    valid = kTRUE;
2663    return value->GetName();
2664 }
2665
2666 //______________________________________________________________________________
2667 void AliAnalysisManager::SetGlobalInt(const char *key, Int_t value)
2668 {
2669 // Define a custom integer variable mapped to a global unique name. The variable
2670 // can be then retrieved by a given analysis macro via GetGlobalInt(key).
2671    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2672    if (!mgr) {
2673       ::Error("AliAnalysisManager::SetGlobalStr", "No analysis manager defined");
2674       return;
2675    }   
2676    Bool_t valid = kFALSE;
2677    Int_t existing = AliAnalysisManager::GetGlobalInt(key, valid);
2678    if (valid) {
2679       ::Error("AliAnalysisManager::SetGlobalInt", "Global %s = %i already defined.", key, existing);
2680       return;
2681    }
2682    mgr->GetGlobals()->Add(new TObjString(key), new TObjString(TString::Format("%i",value)));
2683 }
2684
2685 //______________________________________________________________________________
2686 Int_t AliAnalysisManager::GetGlobalInt(const char *key, Bool_t &valid)
2687 {
2688 // Static method to retrieve a global variable defined via SetGlobalInt.
2689    valid = kFALSE;
2690    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2691    if (!mgr) return 0;
2692    TObject *value = mgr->GetGlobals()->GetValue(key);
2693    if (!value) return 0;
2694    valid = kTRUE;
2695    TString s = value->GetName();
2696    return s.Atoi();
2697 }
2698
2699 //______________________________________________________________________________
2700 void AliAnalysisManager::SetGlobalDbl(const char *key, Double_t value)
2701 {
2702 // Define a custom double precision variable mapped to a global unique name. The variable
2703 // can be then retrieved by a given analysis macro via GetGlobalInt(key).
2704    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2705    if (!mgr) {
2706       ::Error("AliAnalysisManager::SetGlobalStr", "No analysis manager defined");
2707       return;
2708    }   
2709    Bool_t valid = kFALSE;
2710    Double_t existing = AliAnalysisManager::GetGlobalDbl(key, valid);
2711    if (valid) {
2712       ::Error("AliAnalysisManager::SetGlobalInt", "Global %s = %g already defined.", key, existing);
2713       return;
2714    }
2715    mgr->GetGlobals()->Add(new TObjString(key), new TObjString(TString::Format("%f.16",value)));
2716 }
2717
2718 //______________________________________________________________________________
2719 Double_t AliAnalysisManager::GetGlobalDbl(const char *key, Bool_t &valid)
2720 {
2721 // Static method to retrieve a global variable defined via SetGlobalDbl.
2722    valid = kFALSE;
2723    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2724    if (!mgr) return 0;
2725    TObject *value = mgr->GetGlobals()->GetValue(key);
2726    if (!value) return 0;
2727    valid = kTRUE;
2728    TString s = value->GetName();
2729    return s.Atof();
2730 }
2731
2732 //______________________________________________________________________________
2733 void AliAnalysisManager::AddClassDebug(const char *className, Int_t debugLevel)
2734 {
2735 // Sets Class debug level
2736
2737    if (!fDebugOptions) {
2738       fDebugOptions = new TObjArray();
2739       fDebugOptions->SetOwner(kTRUE);
2740    }
2741
2742    // substracting DebugOffset, beacuse of AliLog::SetClassDebugLevel()
2743    debugLevel -= AliLog::kDebug-1;
2744
2745    TNamed *debugOpt = (TNamed*)fDebugOptions->FindObject(className);
2746    if (!debugOpt) {
2747      AliInfo(TString::Format("Adding debug level %d for class %s",debugLevel+AliLog::kDebug-1,className).Data());
2748      fDebugOptions->Add(new TNamed(className,TString::Format("%d",debugLevel).Data()));
2749    } else {
2750       TString oldDebugStr = debugOpt->GetTitle();
2751       Int_t oldDebug = oldDebugStr.Atoi();
2752       if (debugLevel > oldDebug) {
2753          AliWarning(TString::Format("Overwriting debug level to %d class %s, because it is higher then previously set (%d).",debugLevel+AliLog::kDebug-1,className,oldDebug+AliLog::kDebug-1).Data());
2754          debugOpt->SetTitle(TString::Format("%d",debugLevel).Data());
2755       } else {
2756          AliWarning(TString::Format("Ignoring debug level to %d class %s, because it is smaller then previously set (%d).",debugLevel+AliLog::kDebug-1,className,oldDebug+AliLog::kDebug-1).Data());
2757       }
2758    }
2759 }
2760
2761 //______________________________________________________________________________
2762 void AliAnalysisManager::ApplyDebugOptions()
2763 {
2764 // Apply debug options
2765
2766    if (!fDebugOptions) return;
2767    
2768    TIter next(fDebugOptions);
2769    TNamed *debug;
2770    TString debugLevel;
2771    while ((debug=dynamic_cast<TNamed*>(next()))) {
2772       debugLevel = debug->GetTitle();
2773       AliInfo(TString::Format("ApplyDebugOptions : Class=%s debulLevel=%d",debug->GetName(),debugLevel.Atoi()+AliLog::kDebug-1).Data());
2774       AliLog::SetClassDebugLevel(debug->GetName(), debugLevel.Atoi());
2775    }
2776 }
2777
2778 //______________________________________________________________________________
2779 void AliAnalysisManager::Lock()
2780 {
2781 // Security lock. This is to detect NORMAL user errors and not really to
2782 // protect against intentional hacks.
2783    if (fLocked) return;
2784    fLocked = kTRUE;
2785    if (fInputEventHandler)  fInputEventHandler->Lock();
2786    if (fOutputEventHandler) fOutputEventHandler->Lock();
2787    if (fMCtruthEventHandler) fMCtruthEventHandler->Lock();
2788    Info("Lock","====== ANALYSIS MANAGER LOCKED ======");
2789 }
2790
2791 //______________________________________________________________________________
2792 void AliAnalysisManager::UnLock()
2793 {
2794 // Verbose unlocking. Hackers will be punished ;-) ... 
2795    if (!fLocked) return;
2796    fLocked = kFALSE;
2797    if (fInputEventHandler)  fInputEventHandler->UnLock();
2798    if (fOutputEventHandler) fOutputEventHandler->UnLock();
2799    if (fMCtruthEventHandler) fMCtruthEventHandler->UnLock();
2800    Info("UnLock", "====== ANALYSIS MANAGER UNLOCKED ======");
2801 }
2802
2803 //______________________________________________________________________________
2804 void AliAnalysisManager::Changed()
2805 {
2806 // All critical setters pass through the Changed method that throws an exception 
2807 // in case the lock was set.
2808    if (fLocked) Fatal("Changed","Critical setter called in locked mode");
2809 }