1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
17 // Author: Andrei Gheata, 31/05/2006
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 //==============================================================================
26 //==============================================================================
33 #include "AliAnalysisManager.h"
34 #include "AliAnalysisTask.h"
35 #include "AliAnalysisDataContainer.h"
36 #include "AliAnalysisDataSlot.h"
38 ClassImp(AliAnalysisManager)
40 //______________________________________________________________________________
41 AliAnalysisManager::AliAnalysisManager() : TSelector(),
51 // Default constructor.
52 if (TClass::IsCallingNew() != TClass::kDummyNew) {
53 fContainers = new TObjArray();
54 fInputs = new TObjArray();
55 fOutputs = new TObjArray();
56 fTasks = new TObjArray();
57 fTopTasks = new TObjArray();
58 fZombies = new TObjArray();
59 // fStatus = new AliAnalysisInfo(this);
64 //______________________________________________________________________________
65 AliAnalysisManager::AliAnalysisManager(const AliAnalysisManager& other)
77 fInitOK = other.fInitOK;
78 fContainers = new TObjArray(*other.fContainers);
79 fInputs = new TObjArray(*other.fInputs);
80 fOutputs = new TObjArray(*other.fOutputs);
81 fTasks = new TObjArray(*other.fTasks);
82 fTopTasks = new TObjArray(*other.fTopTasks);
83 fZombies = new TObjArray(*other.fZombies);
84 // fStatus = new AliAnalysisInfo(this);
87 //______________________________________________________________________________
88 AliAnalysisManager& AliAnalysisManager::operator=(const AliAnalysisManager& other)
92 TSelector::operator=(other);
94 fInitOK = other.fInitOK;
95 fContainers = new TObjArray(*other.fContainers);
96 fInputs = new TObjArray(*other.fInputs);
97 fOutputs = new TObjArray(*other.fOutputs);
98 fTasks = new TObjArray(*other.fTasks);
99 fTopTasks = new TObjArray(*other.fTopTasks);
100 fZombies = new TObjArray(*other.fZombies);
101 // fStatus = new AliAnalysisInfo(this);
107 //______________________________________________________________________________
108 AliAnalysisManager::~AliAnalysisManager()
111 if (fContainers) {fContainers->Delete(); delete fContainers;}
112 if (fInputs) delete fInputs;
113 if (fOutputs) delete fOutputs;
114 if (fTasks) {fTasks->Delete(); delete fTasks;}
115 if (fTopTasks) delete fTopTasks;
116 if (fZombies) delete fZombies;
118 //______________________________________________________________________________
119 Int_t AliAnalysisManager::GetEntry(Long64_t entry, Int_t getall)
121 // Read one entry of the tree or a whole branch.
122 return fTree ? fTree->GetTree()->GetEntry(entry, getall) : 0;
125 //______________________________________________________________________________
126 void AliAnalysisManager::Init(TTree *tree)
128 // The Init() function is called when the selector needs to initialize
129 // a new tree or chain. Typically here the branch addresses of the tree
130 // will be set. It is normaly not necessary to make changes to the
131 // generated code, but the routine can be extended by the user if needed.
132 // Init() will be called many times when running with PROOF.
133 printf("AliAnalysisManager::Init(%s)\n", tree->GetName());
135 AliError("You have to call InitAnalysis first");
140 AliAnalysisDataContainer *top = (AliAnalysisDataContainer*)fInputs->At(0);
144 //______________________________________________________________________________
145 void AliAnalysisManager::Begin(TTree *tree)
147 // The Begin() function is called at the start of the query.
148 // When running with PROOF Begin() is only called on the client.
149 // The tree argument is deprecated (on PROOF 0 is passed).
150 printf("AliAnalysisManager::Begin(%s)\n", tree->GetName());
154 //______________________________________________________________________________
155 void AliAnalysisManager::SlaveBegin(TTree *tree)
157 // The SlaveBegin() function is called after the Begin() function.
158 // When running with PROOF SlaveBegin() is called on each slave server.
159 // The tree argument is deprecated (on PROOF 0 is passed).
160 printf("AliAnalysisManager::SlaveBegin(%s)\n", tree->GetName());
164 //______________________________________________________________________________
165 Bool_t AliAnalysisManager::Notify()
167 // The Notify() function is called when a new file is opened. This
168 // can be either for a new TTree in a TChain or when when a new TTree
169 // is started when using PROOF. It is normaly not necessary to make changes
170 // to the generated code, but the routine can be extended by the
171 // user if needed. The return value is currently not used.
173 TFile *curfile = fTree->GetCurrentFile();
174 if (curfile) printf("AliAnalysisManager::Notify() file: %s\n", curfile->GetName());
179 //______________________________________________________________________________
180 Bool_t AliAnalysisManager::Process(Long64_t entry)
182 // The Process() function is called for each entry in the tree (or possibly
183 // keyed object in the case of PROOF) to be processed. The entry argument
184 // specifies which entry in the currently loaded tree is to be processed.
185 // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
186 // to read either all or the required parts of the data. When processing
187 // keyed objects with PROOF, the object is already loaded and is available
188 // via the fObject pointer.
190 // This function should contain the "body" of the analysis. It can contain
191 // simple or elaborate selection criteria, run algorithms on the data
192 // of the event and typically fill histograms.
194 // WARNING when a selector is used with a TChain, you must use
195 // the pointer to the current TTree to call GetEntry(entry).
196 // The entry is always the local entry number in the current tree.
197 // Assuming that fChain is the pointer to the TChain being processed,
198 // use fChain->GetTree()->GetEntry(entry).
200 // printf("AliAnalysisManager::Process(%lld)\n", entry);
206 //______________________________________________________________________________
207 void AliAnalysisManager::SlaveTerminate()
209 // The SlaveTerminate() function is called after all entries or objects
210 // have been processed. When running with PROOF SlaveTerminate() is called
211 // on each slave server.
213 printf("AliAnalysisManager::SlaveTerminate()\n");
216 AliError("ERROR: Output list not initialized.");
219 TIter next(fOutputs);
220 AliAnalysisDataContainer *output;
221 while ((output=(AliAnalysisDataContainer *)next())) {
222 output->SetDataOwned(kFALSE);
223 fOutput->Add(output->GetData());
227 //______________________________________________________________________________
228 void AliAnalysisManager::Terminate()
230 // The Terminate() function is the last function to be called during
231 // a query. It always runs on the client, it can be used to present
232 // the results graphically or save the results to file.
233 printf("AliAnalysisManager::Terminate()\n");
234 AliAnalysisDataContainer *output;
235 AliAnalysisTask *task;
236 TIter next(fOutputs);
237 while ((output=(AliAnalysisDataContainer *)next())) output->WriteData();
239 // Call Terminate() for tasks
240 while ((task=(AliAnalysisTask*)next1())) task->Terminate();
243 //______________________________________________________________________________
244 void AliAnalysisManager::AddTask(AliAnalysisTask *task)
246 // Adds a user task to the global list of tasks.
247 task->SetActive(kFALSE);
251 //______________________________________________________________________________
252 AliAnalysisTask *AliAnalysisManager::GetTask(const char *name) const
254 // Retreive task by name.
255 if (!fTasks) return NULL;
256 return (AliAnalysisTask*)fTasks->FindObject(name);
259 //______________________________________________________________________________
260 AliAnalysisDataContainer *AliAnalysisManager::CreateContainer(const char *name,
261 TClass *datatype, EAliAnalysisContType type)
263 // Create a data container of a certain type. Types can be:
264 // kNormalContainer = 0, used to exchange date between tasks
265 // kInputContainer = 1, used to store input data
266 // kOutputContainer = 2, used for posting results
267 AliAnalysisDataContainer *cont = new AliAnalysisDataContainer(name, datatype);
268 fContainers->Add(cont);
270 case kInputContainer:
273 case kOutputContainer:
276 case kNormalContainer:
282 //______________________________________________________________________________
283 Bool_t AliAnalysisManager::ConnectInput(AliAnalysisTask *task, Int_t islot,
284 AliAnalysisDataContainer *cont)
286 // Connect input of an existing task to a data container.
287 if (!fTasks->FindObject(task)) {
289 AliInfo(Form("Task %s not registered. Now owned by analysis manager", task->GetName()));
291 Bool_t connected = task->ConnectInput(islot, cont);
295 //______________________________________________________________________________
296 Bool_t AliAnalysisManager::ConnectOutput(AliAnalysisTask *task, Int_t islot,
297 AliAnalysisDataContainer *cont)
299 // Connect output of an existing task to a data container.
300 if (!fTasks->FindObject(task)) {
302 AliInfo(Form("Task %s not registered. Now owned by analysis manager", task->GetName()));
304 Bool_t connected = task->ConnectOutput(islot, cont);
308 //______________________________________________________________________________
309 void AliAnalysisManager::CleanContainers()
311 // Clean data from all containers that have already finished all client tasks.
312 TIter next(fContainers);
313 AliAnalysisDataContainer *cont;
314 while ((cont=(AliAnalysisDataContainer *)next())) {
315 if (cont->IsOwnedData() &&
316 cont->IsDataReady() &&
317 cont->ClientsExecuted()) cont->DeleteData();
321 //______________________________________________________________________________
322 Bool_t AliAnalysisManager::InitAnalysis()
324 // Initialization of analysis chain of tasks. Should be called after all tasks
325 // and data containers are properly connected
326 // Check for input/output containers
328 // if (!fInputs->GetEntriesFast()) {
329 // AliError("No input container defined. At least one container should store input data");
332 // if (!fOutputs->GetEntriesFast()) {
333 // AliError("No output container defined. At least one container should store output data");
336 // Check for top tasks (depending only on input data containers)
337 if (!fTasks->First()) {
338 AliError("Analysis have no tasks !");
342 AliAnalysisTask *task;
343 AliAnalysisDataContainer *cont;
346 Bool_t iszombie = kFALSE;
347 Bool_t istop = kTRUE;
349 while ((task=(AliAnalysisTask*)next())) {
352 Int_t ninputs = task->GetNinputs();
354 // task->SetZombie();
355 // fZombies->Add(task);
357 // AliWarning(Form("Task %s has no input slots defined ! Declared zombie...",task->GetName()));
360 for (i=0; i<ninputs; i++) {
361 cont = task->GetInputSlot(i)->GetContainer();
369 AliWarning(Form("Input slot %i of task %s has no container connected ! Declared zombie...",
372 if (iszombie) continue;
373 // Check if cont is an input container
374 if (istop && !fInputs->FindObject(cont)) istop=kFALSE;
375 // Connect to parent task
379 fTopTasks->Add(task);
383 AliError("No top task defined. At least one task should be connected only to input containers");
386 // Check now if there are orphan tasks
387 for (i=0; i<ntop; i++) {
388 task = (AliAnalysisTask*)fTopTasks->At(i);
393 while ((task=(AliAnalysisTask*)next())) {
394 if (!task->IsUsed()) {
396 AliWarning(Form("Task %s is orphan",task->GetName()));
399 // Check the task hierarchy (no parent task should depend on data provided
400 // by a daughter task)
401 for (i=0; i<ntop; i++) {
402 task = (AliAnalysisTask*)fTopTasks->At(i);
403 if (task->CheckCircularDeps()) {
404 AliError("Found illegal circular dependencies between following tasks:");
413 //______________________________________________________________________________
414 void AliAnalysisManager::PrintStatus(Option_t *option) const
416 // Print task hierarchy.
417 TIter next(fTopTasks);
418 AliAnalysisTask *task;
419 while ((task=(AliAnalysisTask*)next()))
420 task->PrintTask(option);
423 //______________________________________________________________________________
424 void AliAnalysisManager::ResetAnalysis()
426 // Reset all execution flags and clean containers.
430 //______________________________________________________________________________
431 void AliAnalysisManager::ExecAnalysis(Option_t *option)
435 AliError("Analysis manager was not initialized !");
438 AliAnalysisTask *task;
439 // Check if the top tree is active.
442 // De-activate all tasks
443 while ((task=(AliAnalysisTask*)next())) task->SetActive(kFALSE);
444 AliAnalysisDataContainer *cont = (AliAnalysisDataContainer*)fInputs->At(0);
446 AliError("Cannot execute analysis in TSelector mode without at least one top container");
449 cont->SetData(fTree); // This will notify all consumers
450 TIter next1(cont->GetConsumers());
451 while ((task=(AliAnalysisTask*)next1())) {
452 // task->SetActive(kTRUE);
453 task->ExecuteTask(option);
457 // The event loop is not controlled by TSelector
458 TIter next2(fTopTasks);
459 while ((task=(AliAnalysisTask*)next2())) {
460 task->SetActive(kTRUE);
461 printf("executing %s\n", task->GetName());
462 task->ExecuteTask(option);
466 //______________________________________________________________________________
467 void AliAnalysisManager::FinishAnalysis()