ce706320d9e39eefc4c8925f1cd8595fd97fba53
[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 //   AliAnalysysManager - 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 <Riostream.h>
29
30 #include <TClass.h>
31 #include <TFile.h>
32 #include <TMethodCall.h>
33 #include <TChain.h>
34 #include <TSystem.h>
35 #include <TROOT.h>
36
37 #include "AliAnalysisTask.h"
38 #include "AliAnalysisDataContainer.h"
39 #include "AliAnalysisDataSlot.h"
40 #include "AliVirtualEventHandler.h"
41 #include "AliAnalysisManager.h"
42
43 ClassImp(AliAnalysisManager)
44
45 AliAnalysisManager *AliAnalysisManager::fgAnalysisManager = NULL;
46
47 //______________________________________________________________________________
48 AliAnalysisManager::AliAnalysisManager() 
49                    :TNamed(),
50                     fTree(NULL),
51                     fOutputEventHandler(NULL),
52                     fMCtruthEventHandler(NULL),
53                     fCurrentEntry(-1),
54                     fMode(kLocalAnalysis),
55                     fInitOK(kFALSE),
56                     fDebug(0),
57                     fTasks(NULL),
58                     fTopTasks(NULL),
59                     fZombies(NULL),
60                     fContainers(NULL),
61                     fInputs(NULL),
62                     fOutputs(NULL)
63 {
64 // Dummy constructor.
65    fgAnalysisManager = this;
66    SetEventLoop(kTRUE);
67 }
68
69 //______________________________________________________________________________
70 AliAnalysisManager::AliAnalysisManager(const char *name, const char *title)
71                    :TNamed(name,title),
72                     fTree(NULL),
73                     fOutputEventHandler(NULL),
74                     fMCtruthEventHandler(NULL),
75                     fCurrentEntry(-1),
76                     fMode(kLocalAnalysis),
77                     fInitOK(kFALSE),
78                     fDebug(0),
79                     fTasks(NULL),
80                     fTopTasks(NULL),
81                     fZombies(NULL),
82                     fContainers(NULL),
83                     fInputs(NULL),
84                     fOutputs(NULL)
85 {
86 // Default constructor.
87    fgAnalysisManager = this;
88    fTasks      = new TObjArray();
89    fTopTasks   = new TObjArray();
90    fZombies    = new TObjArray();
91    fContainers = new TObjArray();
92    fInputs     = new TObjArray();
93    fOutputs    = new TObjArray();
94    SetEventLoop(kTRUE);
95 }
96
97 //______________________________________________________________________________
98 AliAnalysisManager::AliAnalysisManager(const AliAnalysisManager& other)
99                    :TNamed(other),
100                     fTree(NULL),
101                     fOutputEventHandler(NULL),
102                     fMCtruthEventHandler(NULL),
103                     fCurrentEntry(-1),
104                     fMode(other.fMode),
105                     fInitOK(other.fInitOK),
106                     fDebug(other.fDebug),
107                     fTasks(NULL),
108                     fTopTasks(NULL),
109                     fZombies(NULL),
110                     fContainers(NULL),
111                     fInputs(NULL),
112                     fOutputs(NULL)
113 {
114 // Copy constructor.
115    fTasks      = new TObjArray(*other.fTasks);
116    fTopTasks   = new TObjArray(*other.fTopTasks);
117    fZombies    = new TObjArray(*other.fZombies);
118    fContainers = new TObjArray(*other.fContainers);
119    fInputs     = new TObjArray(*other.fInputs);
120    fOutputs    = new TObjArray(*other.fOutputs);
121    fgAnalysisManager = this;
122 }
123    
124 //______________________________________________________________________________
125 AliAnalysisManager& AliAnalysisManager::operator=(const AliAnalysisManager& other)
126 {
127 // Assignment
128    if (&other != this) {
129       TNamed::operator=(other);
130       fOutputEventHandler  = other.fOutputEventHandler;
131       fMCtruthEventHandler = other.fMCtruthEventHandler;
132       fTree       = NULL;
133       fCurrentEntry = -1;
134       fMode       = other.fMode;
135       fInitOK     = other.fInitOK;
136       fDebug      = other.fDebug;
137       fTasks      = new TObjArray(*other.fTasks);
138       fTopTasks   = new TObjArray(*other.fTopTasks);
139       fZombies    = new TObjArray(*other.fZombies);
140       fContainers = new TObjArray(*other.fContainers);
141       fInputs     = new TObjArray(*other.fInputs);
142       fOutputs    = new TObjArray(*other.fOutputs);
143       fgAnalysisManager = this;
144    }
145    return *this;
146 }
147
148 //______________________________________________________________________________
149 AliAnalysisManager::~AliAnalysisManager()
150 {
151 // Destructor.
152    if (fTasks) {fTasks->Delete(); delete fTasks;}
153    if (fTopTasks) delete fTopTasks;
154    if (fZombies) delete fZombies;
155    if (fContainers) {fContainers->Delete(); delete fContainers;}
156    if (fInputs) delete fInputs;
157    if (fOutputs) delete fOutputs;
158    if (fgAnalysisManager==this) fgAnalysisManager = NULL;
159 }
160
161 //______________________________________________________________________________
162 Int_t AliAnalysisManager::GetEntry(Long64_t entry, Int_t getall)
163 {
164 // Read one entry of the tree or a whole branch.
165    if (fDebug > 1) {
166       cout << "== AliAnalysisManager::GetEntry()" << endl;
167    }   
168    fCurrentEntry = entry;
169    return fTree ? fTree->GetTree()->GetEntry(entry, getall) : 0;
170 }
171    
172 //______________________________________________________________________________
173 void AliAnalysisManager::Init(TTree *tree)
174 {
175   // The Init() function is called when the selector needs to initialize
176   // a new tree or chain. Typically here the branch addresses of the tree
177   // will be set. It is normaly not necessary to make changes to the
178   // generated code, but the routine can be extended by the user if needed.
179   // Init() will be called many times when running with PROOF.
180    if (!tree) return;
181    if (fDebug > 1) {
182       printf("->AliAnalysisManager::Init(%s)\n", tree->GetName());
183    }
184    if (!fInitOK) InitAnalysis();
185    if (!fInitOK) return;
186    fTree = tree;
187    AliAnalysisDataContainer *top = (AliAnalysisDataContainer*)fInputs->At(0);
188    if (!top) {
189       cout<<"Error: No top input container !" <<endl;
190       return;
191    }
192    top->SetData(tree);
193    if (fDebug > 1) {
194       printf("<-AliAnalysisManager::Init(%s)\n", tree->GetName());
195    }
196 }
197
198 //______________________________________________________________________________
199 void AliAnalysisManager::Begin(TTree *tree)
200 {
201   // The Begin() function is called at the start of the query.
202   // When running with PROOF Begin() is only called on the client.
203   // The tree argument is deprecated (on PROOF 0 is passed).
204    if (fDebug > 1) {
205       cout << "AliAnalysisManager::Begin()" << endl;
206    }   
207    Init(tree);
208 }
209
210 //______________________________________________________________________________
211 void AliAnalysisManager::SlaveBegin(TTree *tree)
212 {
213   // The SlaveBegin() function is called after the Begin() function.
214   // When running with PROOF SlaveBegin() is called on each slave server.
215   // The tree argument is deprecated (on PROOF 0 is passed).
216    if (fDebug > 1) {
217       cout << "->AliAnalysisManager::SlaveBegin()" << endl;
218    }
219    // Call InitIO of EventHandler
220    if (fOutputEventHandler) {
221        if (fMode == kProofAnalysis) {
222            fOutputEventHandler->InitIO("proof");
223        } else {
224            fOutputEventHandler->InitIO("local");
225        }
226    }
227    //
228    TIter next(fTasks);
229    AliAnalysisTask *task;
230    // Call CreateOutputObjects for all tasks
231    while ((task=(AliAnalysisTask*)next())) {
232       TDirectory *curdir = gDirectory;
233       task->CreateOutputObjects();
234       if (curdir) curdir->cd();
235    }   
236    if (fMode == kLocalAnalysis) Init(tree);   
237    if (fDebug > 1) {
238       cout << "<-AliAnalysisManager::SlaveBegin()" << endl;
239    }
240 }
241
242 //______________________________________________________________________________
243 Bool_t AliAnalysisManager::Notify()
244 {
245    // The Notify() function is called when a new file is opened. This
246    // can be either for a new TTree in a TChain or when when a new TTree
247    // is started when using PROOF. It is normaly not necessary to make changes
248    // to the generated code, but the routine can be extended by the
249    // user if needed. The return value is currently not used.
250    if (fTree) {
251       TFile *curfile = fTree->GetCurrentFile();
252       if (curfile && fDebug>1) printf("AliAnalysisManager::Notify() file: %s\n", curfile->GetName());
253       TIter next(fTasks);
254       AliAnalysisTask *task;
255       // Call Notify for all tasks
256       while ((task=(AliAnalysisTask*)next())) 
257          task->Notify();
258       // Call Notify of the MC truth handler
259       if (fMCtruthEventHandler) fMCtruthEventHandler->Notify();
260    }
261    return kTRUE;
262 }    
263
264 //______________________________________________________________________________
265 Bool_t AliAnalysisManager::Process(Long64_t entry)
266 {
267   // The Process() function is called for each entry in the tree (or possibly
268   // keyed object in the case of PROOF) to be processed. The entry argument
269   // specifies which entry in the currently loaded tree is to be processed.
270   // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
271   // to read either all or the required parts of the data. When processing
272   // keyed objects with PROOF, the object is already loaded and is available
273   // via the fObject pointer.
274   //
275   // This function should contain the "body" of the analysis. It can contain
276   // simple or elaborate selection criteria, run algorithms on the data
277   // of the event and typically fill histograms.
278
279   // WARNING when a selector is used with a TChain, you must use
280   //  the pointer to the current TTree to call GetEntry(entry).
281   //  The entry is always the local entry number in the current tree.
282   //  Assuming that fChain is the pointer to the TChain being processed,
283   //  use fChain->GetTree()->GetEntry(entry).
284    if (fDebug > 1) {
285       cout << "->AliAnalysisManager::Process()" << endl;
286    }
287
288    if (fOutputEventHandler)  fOutputEventHandler ->BeginEvent();
289    if (fMCtruthEventHandler) fMCtruthEventHandler->BeginEvent();
290    
291    GetEntry(entry);
292    ExecAnalysis();
293    if (fDebug > 1) {
294       cout << "<-AliAnalysisManager::Process()" << endl;
295    }
296    return kTRUE;
297 }
298
299 //______________________________________________________________________________
300 void AliAnalysisManager::PackOutput(TList *target)
301 {
302   // Pack all output data containers in the output list. Called at SlaveTerminate
303   // stage in PROOF case for each slave.
304    if (fDebug > 1) {
305       cout << "->AliAnalysisManager::PackOutput()" << endl;
306    }   
307    if (!target) {
308       Error("PackOutput", "No target. Aborting.");
309       return;
310    }
311
312    if (fOutputEventHandler)  fOutputEventHandler ->Terminate();
313    if (fMCtruthEventHandler) fMCtruthEventHandler->Terminate();
314    
315    if (fMode == kProofAnalysis) {
316       TIter next(fOutputs);
317       AliAnalysisDataContainer *output;
318       while ((output=(AliAnalysisDataContainer*)next())) {
319          if (fDebug > 1) printf("   Packing container %s...\n", output->GetName());
320          if (output->GetData()) target->Add(output->ExportData());
321       }
322    } 
323    if (fDebug > 1) {
324       printf("<-AliAnalysisManager::PackOutput: output list contains %d containers\n", target->GetSize());
325    }
326 }
327
328 //______________________________________________________________________________
329 void AliAnalysisManager::ImportWrappers(TList *source)
330 {
331 // Import data in output containers from wrappers coming in source.
332    if (fDebug > 1) {
333       cout << "->AliAnalysisManager::ImportWrappers()" << endl;
334    }   
335    TIter next(fOutputs);
336    AliAnalysisDataContainer *cont;
337    AliAnalysisDataWrapper   *wrap;
338    Int_t icont = 0;
339    while ((cont=(AliAnalysisDataContainer*)next())) {
340       wrap = (AliAnalysisDataWrapper*)source->FindObject(cont->GetName());
341       if (!wrap && fDebug>1) {
342          printf("(WW) ImportWrappers: container %s not found in analysis output !\n", cont->GetName());
343          continue;
344       }
345       icont++;
346       if (fDebug > 1) printf("   Importing data for container %s\n", wrap->GetName());
347       if (cont->GetFileName()) printf("    -> %s\n", cont->GetFileName());
348       cont->ImportData(wrap);
349    }         
350    if (fDebug > 1) {
351       cout << "<-AliAnalysisManager::ImportWrappers(): "<< icont << " containers imported" << endl;
352    }   
353 }
354
355 //______________________________________________________________________________
356 void AliAnalysisManager::UnpackOutput(TList *source)
357 {
358   // Called by AliAnalysisSelector::Terminate. Output containers should
359   // be in source in the same order as in fOutputs.
360    if (fDebug > 1) {
361       cout << "->AliAnalysisManager::UnpackOutput()" << endl;
362    }   
363    if (!source) {
364       Error("UnpackOutput", "No target. Aborting.");
365       return;
366    }
367    if (fDebug > 1) {
368       printf("   Source list contains %d containers\n", source->GetSize());
369    }   
370
371    if (fMode == kProofAnalysis) ImportWrappers(source);
372
373    TIter next(fOutputs);
374    AliAnalysisDataContainer *output;
375    while ((output=(AliAnalysisDataContainer*)next())) {
376       if (!output->GetData()) continue;
377       // Check if there are client tasks that run post event loop
378       if (output->HasConsumers()) {
379          // Disable event loop semaphore
380          output->SetPostEventLoop(kTRUE);
381          TObjArray *list = output->GetConsumers();
382          Int_t ncons = list->GetEntriesFast();
383          for (Int_t i=0; i<ncons; i++) {
384             AliAnalysisTask *task = (AliAnalysisTask*)list->At(i);
385             task->CheckNotify(kTRUE);
386             // If task is active, execute it
387             if (task->IsPostEventLoop() && task->IsActive()) {
388                if (fDebug > 1) {
389                   cout << "== Executing post event loop task " << task->GetName() << endl;
390                }                  
391                task->ExecuteTask();
392             }   
393          }
394       }   
395       // Check if the output need to be written to a file.
396       const char *filename = output->GetFileName();
397       if (!(strcmp(filename, "default"))) {
398           if (fOutputEventHandler) filename = fOutputEventHandler->GetOutputFileName();
399       }
400       
401       if (!filename || !strlen(filename)) continue;
402       TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
403       if (file) file->cd();
404       else      file = new TFile(filename, "RECREATE");
405       if (file->IsZombie()) continue;
406       // Reparent data to this file
407       TMethodCall callEnv;
408       if (output->GetData()->IsA())
409          callEnv.InitWithPrototype(output->GetData()->IsA(), "SetDirectory", "TDirectory*");
410       if (callEnv.IsValid()) {
411          callEnv.SetParam((Long_t) file);
412          callEnv.Execute(output->GetData());
413       }
414       output->GetData()->Write();
415    }
416    if (fDebug > 1) {
417       cout << "<-AliAnalysisManager::UnpackOutput()" << endl;
418    }   
419 }
420
421 //______________________________________________________________________________
422 void AliAnalysisManager::Terminate()
423 {
424   // The Terminate() function is the last function to be called during
425   // a query. It always runs on the client, it can be used to present
426   // the results graphically.
427    if (fDebug > 1) {
428       cout << "->AliAnalysisManager::Terminate()" << endl;
429    }   
430    AliAnalysisTask *task;
431    TIter next(fTasks);
432    // Call Terminate() for tasks
433    while ((task=(AliAnalysisTask*)next())) task->Terminate();
434    if (fDebug > 1) {
435       cout << "<-AliAnalysisManager::Terminate()" << endl;
436    }   
437    //
438    if (fOutputEventHandler) fOutputEventHandler->TerminateIO();
439 }
440
441 //______________________________________________________________________________
442 void AliAnalysisManager::AddTask(AliAnalysisTask *task)
443 {
444 // Adds a user task to the global list of tasks.
445    task->SetActive(kFALSE);
446    fTasks->Add(task);
447 }  
448
449 //______________________________________________________________________________
450 AliAnalysisTask *AliAnalysisManager::GetTask(const char *name) const
451 {
452 // Retreive task by name.
453    if (!fTasks) return NULL;
454    return (AliAnalysisTask*)fTasks->FindObject(name);
455 }
456
457 //______________________________________________________________________________
458 AliAnalysisDataContainer *AliAnalysisManager::CreateContainer(const char *name, 
459                                 TClass *datatype, EAliAnalysisContType type, const char *filename)
460 {
461 // Create a data container of a certain type. Types can be:
462 //   kExchangeContainer  = 0, used to exchange date between tasks
463 //   kInputContainer   = 1, used to store input data
464 //   kOutputContainer  = 2, used for posting results
465    if (fContainers->FindObject(name)) {
466       Error("CreateContainer","A container named %s already defined !\n",name);
467       return NULL;
468    }   
469    AliAnalysisDataContainer *cont = new AliAnalysisDataContainer(name, datatype);
470    fContainers->Add(cont);
471    switch (type) {
472       case kInputContainer:
473          fInputs->Add(cont);
474          break;
475       case kOutputContainer:
476          fOutputs->Add(cont);
477          if (filename && strlen(filename)) cont->SetFileName(filename);
478          break;
479       case kExchangeContainer:
480          break;   
481    }
482    return cont;
483 }
484          
485 //______________________________________________________________________________
486 Bool_t AliAnalysisManager::ConnectInput(AliAnalysisTask *task, Int_t islot,
487                                         AliAnalysisDataContainer *cont)
488 {
489 // Connect input of an existing task to a data container.
490    if (!fTasks->FindObject(task)) {
491       AddTask(task);
492       Warning("ConnectInput", "Task %s not registered. Now owned by analysis manager", task->GetName());
493    } 
494    Bool_t connected = task->ConnectInput(islot, cont);
495    return connected;
496 }   
497
498 //______________________________________________________________________________
499 Bool_t AliAnalysisManager::ConnectOutput(AliAnalysisTask *task, Int_t islot,
500                                         AliAnalysisDataContainer *cont)
501 {
502 // Connect output of an existing task to a data container.
503    if (!fTasks->FindObject(task)) {
504       AddTask(task);
505       Warning("ConnectOutput", "Task %s not registered. Now owned by analysis manager", task->GetName());
506    } 
507    Bool_t connected = task->ConnectOutput(islot, cont);
508    return connected;
509 }   
510                                
511 //______________________________________________________________________________
512 void AliAnalysisManager::CleanContainers()
513 {
514 // Clean data from all containers that have already finished all client tasks.
515    TIter next(fContainers);
516    AliAnalysisDataContainer *cont;
517    while ((cont=(AliAnalysisDataContainer *)next())) {
518       if (cont->IsOwnedData() && 
519           cont->IsDataReady() && 
520           cont->ClientsExecuted()) cont->DeleteData();
521    }
522 }
523
524 //______________________________________________________________________________
525 Bool_t AliAnalysisManager::InitAnalysis()
526 {
527 // Initialization of analysis chain of tasks. Should be called after all tasks
528 // and data containers are properly connected
529    // Check for input/output containers
530    fInitOK = kFALSE;
531    // Check for top tasks (depending only on input data containers)
532    if (!fTasks->First()) {
533       Error("InitAnalysis", "Analysis has no tasks !");
534       return kFALSE;
535    }   
536    TIter next(fTasks);
537    AliAnalysisTask *task;
538    AliAnalysisDataContainer *cont;
539    Int_t ntop = 0;
540    Int_t nzombies = 0;
541    Bool_t iszombie = kFALSE;
542    Bool_t istop = kTRUE;
543    Int_t i;
544    while ((task=(AliAnalysisTask*)next())) {
545       istop = kTRUE;
546       iszombie = kFALSE;
547       Int_t ninputs = task->GetNinputs();
548       for (i=0; i<ninputs; i++) {
549          cont = task->GetInputSlot(i)->GetContainer();
550          if (!cont) {
551             if (!iszombie) {
552                task->SetZombie();
553                fZombies->Add(task);
554                nzombies++;
555                iszombie = kTRUE;
556             }   
557             Error("InitAnalysis", "Input slot %d of task %s has no container connected ! Declared zombie...", 
558                   i, task->GetName()); 
559          }
560          if (iszombie) continue;
561          // Check if cont is an input container
562          if (istop && !fInputs->FindObject(cont)) istop=kFALSE;
563          // Connect to parent task
564       }
565       if (istop) {
566          ntop++;
567          fTopTasks->Add(task);
568       }
569    }
570    if (!ntop) {
571       Error("InitAnalysis", "No top task defined. At least one task should be connected only to input containers");
572       return kFALSE;
573    }                        
574    // Check now if there are orphan tasks
575    for (i=0; i<ntop; i++) {
576       task = (AliAnalysisTask*)fTopTasks->At(i);
577       task->SetUsed();
578    }
579    Int_t norphans = 0;
580    next.Reset();
581    while ((task=(AliAnalysisTask*)next())) {
582       if (!task->IsUsed()) {
583          norphans++;
584          Warning("InitAnalysis", "Task %s is orphan", task->GetName());
585       }   
586    }          
587    // Check the task hierarchy (no parent task should depend on data provided
588    // by a daughter task)
589    for (i=0; i<ntop; i++) {
590       task = (AliAnalysisTask*)fTopTasks->At(i);
591       if (task->CheckCircularDeps()) {
592          Error("InitAnalysis", "Found illegal circular dependencies between following tasks:");
593          PrintStatus("dep");
594          return kFALSE;
595       }   
596    }
597    // Check that all containers feeding post-event loop tasks are in the outputs list
598    TIter nextcont(fContainers); // loop over all containers
599    while ((cont=(AliAnalysisDataContainer*)nextcont())) {
600       if (!cont->IsPostEventLoop() && !fOutputs->FindObject(cont)) {
601          if (cont->HasConsumers()) {
602          // Check if one of the consumers is post event loop
603             TIter nextconsumer(cont->GetConsumers());
604             while ((task=(AliAnalysisTask*)nextconsumer())) {
605                if (task->IsPostEventLoop()) {
606                   fOutputs->Add(cont);
607                   break;
608                }
609             }
610          }
611       }
612    }   
613    fInitOK = kTRUE;
614    return kTRUE;
615 }   
616
617 //______________________________________________________________________________
618 void AliAnalysisManager::PrintStatus(Option_t *option) const
619 {
620 // Print task hierarchy.
621    TIter next(fTopTasks);
622    AliAnalysisTask *task;
623    while ((task=(AliAnalysisTask*)next()))
624       task->PrintTask(option);
625 }
626
627 //______________________________________________________________________________
628 void AliAnalysisManager::ResetAnalysis()
629 {
630 // Reset all execution flags and clean containers.
631    CleanContainers();
632 }
633
634 //______________________________________________________________________________
635 void AliAnalysisManager::StartAnalysis(const char *type, TTree *tree)
636 {
637 // Start analysis for this manager. Analysis task can be: LOCAL, PROOF or GRID.
638    if (!fInitOK) {
639       Error("StartAnalysis","Analysis manager was not initialized !");
640       return;
641    }
642    if (fDebug>1) {
643       cout << "StartAnalysis: " << GetName() << endl;   
644    }   
645    TString anaType = type;
646    anaType.ToLower();
647    fMode = kLocalAnalysis;
648    if (tree) {
649       if (anaType.Contains("proof"))     fMode = kProofAnalysis;
650       else if (anaType.Contains("grid")) fMode = kGridAnalysis;
651    }
652    if (fMode == kGridAnalysis) {
653       Warning("StartAnalysis", "GRID analysis mode not implemented. Running local.");
654       fMode = kLocalAnalysis;
655    }
656    char line[128];
657    SetEventLoop(kFALSE);
658    // Disable all branches if requested and set event loop mode
659    if (tree) {
660       if (TestBit(kDisableBranches)) {
661          printf("Disabling all branches...\n");
662 //         tree->SetBranchStatus("*",0); // not yet working
663       }   
664       SetEventLoop(kTRUE);
665    }   
666
667    TChain *chain = dynamic_cast<TChain*>(tree);
668
669    // Initialize locally all tasks
670    TIter next(fTasks);
671    AliAnalysisTask *task;
672    while ((task=(AliAnalysisTask*)next())) {
673       task->LocalInit();
674    }
675    
676    switch (fMode) {
677       case kLocalAnalysis:
678          if (!tree) {
679             TIter next(fTasks);
680             AliAnalysisTask *task;
681             // Call CreateOutputObjects for all tasks
682             while ((task=(AliAnalysisTask*)next())) {
683                TDirectory *curdir = gDirectory;
684                task->CreateOutputObjects();
685                if (curdir) curdir->cd();
686             }   
687             ExecAnalysis();
688             Terminate();
689             return;
690          } 
691          // Run tree-based analysis via AliAnalysisSelector  
692 //         gROOT->ProcessLine(".L $ALICE_ROOT/ANALYSIS/AliAnalysisSelector.cxx+");
693          cout << "===== RUNNING LOCAL ANALYSIS " << GetName() << " ON TREE " << tree->GetName() << endl;
694          sprintf(line, "AliAnalysisSelector *selector = new AliAnalysisSelector((AliAnalysisManager*)0x%lx);",(ULong_t)this);
695          gROOT->ProcessLine(line);
696          sprintf(line, "((TTree*)0x%lx)->Process(selector);",(ULong_t)tree);
697          gROOT->ProcessLine(line);
698          break;
699       case kProofAnalysis:
700          if (!gROOT->GetListOfProofs() || !gROOT->GetListOfProofs()->GetEntries()) {
701             printf("StartAnalysis: no PROOF!!!\n");
702             return;
703          }   
704          sprintf(line, "gProof->AddInput((TObject*)0x%lx);", (ULong_t)this);
705          gROOT->ProcessLine(line);
706          if (chain) {
707             chain->SetProof();
708             cout << "===== RUNNING PROOF ANALYSIS " << GetName() << " ON CHAIN " << chain->GetName() << endl;
709             chain->Process("AliAnalysisSelector");
710          } else {
711             printf("StartAnalysis: no chain\n");
712             return;
713          }      
714          break;
715       case kGridAnalysis:
716          Warning("StartAnalysis", "GRID analysis mode not implemented. Running local.");
717    }   
718 }   
719
720 //______________________________________________________________________________
721 void AliAnalysisManager::ExecAnalysis(Option_t *option)
722 {
723 // Execute analysis.
724    if (!fInitOK) {
725      Error("ExecAnalysis", "Analysis manager was not initialized !");
726       return;
727    }   
728    AliAnalysisTask *task;
729    // Check if the top tree is active.
730    if (fTree) {
731       if (fDebug>1) {
732          printf("AliAnalysisManager::ExecAnalysis\n");
733       }   
734       TIter next(fTasks);
735    // De-activate all tasks
736       while ((task=(AliAnalysisTask*)next())) task->SetActive(kFALSE);
737       AliAnalysisDataContainer *cont = (AliAnalysisDataContainer*)fInputs->At(0);
738       if (!cont) {
739               Error("ExecAnalysis","Cannot execute analysis in TSelector mode without at least one top container");
740          return;
741       }   
742       cont->SetData(fTree); // This will notify all consumers
743 //
744 //    Call BeginEvent() for optional output and MC services 
745       if (fOutputEventHandler)  fOutputEventHandler ->BeginEvent();
746       if (fMCtruthEventHandler) fMCtruthEventHandler->BeginEvent();
747 //
748 //    Execute the tasks
749       TIter next1(cont->GetConsumers());
750       while ((task=(AliAnalysisTask*)next1())) {
751          if (fDebug >1) {
752             cout << "    Executing task " << task->GetName() << endl;
753          }   
754          
755          task->ExecuteTask(option);
756       }
757 //
758 //    Call FinishEvent() for optional output and MC services 
759       if (fOutputEventHandler)  fOutputEventHandler ->FinishEvent();
760       if (fMCtruthEventHandler) fMCtruthEventHandler->FinishEvent();
761 //
762       return;
763    }   
764    // The event loop is not controlled by TSelector   
765 //
766 //  Call BeginEvent() for optional output and MC services 
767    if (fOutputEventHandler)  fOutputEventHandler ->BeginEvent();
768    if (fMCtruthEventHandler) fMCtruthEventHandler->BeginEvent();
769    TIter next2(fTopTasks);
770    while ((task=(AliAnalysisTask*)next2())) {
771       task->SetActive(kTRUE);
772       if (fDebug > 1) {
773          cout << "    Executing task " << task->GetName() << endl;
774       }   
775       task->ExecuteTask(option);
776    }   
777 //
778 // Call FinishEvent() for optional output and MC services 
779    if (fOutputEventHandler)  fOutputEventHandler->FinishEvent();
780    if (fMCtruthEventHandler) fMCtruthEventHandler->FinishEvent();
781 }
782
783 //______________________________________________________________________________
784 void AliAnalysisManager::FinishAnalysis()
785 {
786 // Finish analysis.
787 }