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