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