]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - ANALYSIS/AliAnalysisTask.cxx
AliEveTrackCounter, AliEveTrackCounterEditor
[u/mrichter/AliRoot.git] / ANALYSIS / AliAnalysisTask.cxx
index 6a588faa530ba3f7714ca1b1c4ae83ad7af0114a..c4a123a709c3c561cce5e1ccb251af6947f04122 100644 (file)
 #include <Riostream.h>
 #include <TFile.h>
 #include <TClass.h>
+#include <TTree.h>
+#include <TROOT.h>
 
 #include "AliAnalysisTask.h"
 #include "AliAnalysisDataSlot.h"
@@ -120,7 +122,6 @@ AliAnalysisTask::AliAnalysisTask()
                  fOutputs(NULL)
 {
 // Default constructor.
-   TObject::SetBit(kTaskEvtByEvt, kTRUE);
 }
 
 //______________________________________________________________________________
@@ -136,7 +137,6 @@ AliAnalysisTask::AliAnalysisTask(const char *name, const char *title)
                  fOutputs(NULL)                 
 {
 // Constructor.
-   TObject::SetBit(kTaskEvtByEvt, kTRUE);
    fInputs      = new TObjArray(2);
    fOutputs     = new TObjArray(2);
 }
@@ -233,13 +233,11 @@ void AliAnalysisTask::CheckNotify(Bool_t init)
 // accordingly. This method is called automatically for all tasks connected
 // to a container where the data was published.
    if (init) fInitialized = kFALSE;
-//   Bool_t execperevent = IsExecPerEvent();
+   Bool_t single_shot = IsPostEventLoop();
    AliAnalysisDataContainer *cinput;
    for (Int_t islot=0; islot<fNinputs; islot++) {
       cinput = GetInputSlot(islot)->GetContainer();
-      if (!cinput) return;
-      if (!cinput->GetData()) {
-//      if (!cinput->GetData() || execperevent!=cinput->IsEventByEvent()) {
+      if (!cinput->GetData() || (single_shot && !cinput->IsPostEventLoop())) {
          SetActive(kFALSE);
          return;
       }   
@@ -290,6 +288,8 @@ Bool_t AliAnalysisTask::ConnectOutput(Int_t islot, AliAnalysisDataContainer *con
    }            
    // Connect the slot to the container as output         
    if (!output->ConnectContainer(cont)) return kFALSE;
+   // Set event loop type the same as for the task
+   cont->SetPostEventLoop(IsPostEventLoop());
    // Declare this as the data producer
    cont->SetProducer(this, islot);
    AreSlotsConnected();
@@ -385,6 +385,34 @@ Bool_t AliAnalysisTask::SetBranchAddress(Int_t islot, const char *branch, void *
    return GetInputSlot(islot)->SetBranchAddress(branch, address);
 }   
 
+//______________________________________________________________________________
+void AliAnalysisTask::EnableBranch(Int_t islot, const char *bname) const
+{
+// Call this in ConnectInputData() to enable only the branches needed by this 
+// task. "*" will enable everything.
+   AliAnalysisDataSlot *input = GetInputSlot(islot);
+   if (!input || !input->GetType()->InheritsFrom(TTree::Class())) {
+      Error("EnableBranch", "Wrong slot type #%d for task %s: not TTree-derived type", islot, GetName());
+      return;
+   }   
+   TTree *tree = (TTree*)input->GetData();
+   if (!strcmp(bname, "*")) {
+      tree->SetBranchStatus("*",1);
+      return;
+   }
+   AliAnalysisDataSlot::EnableBranch(bname, tree);
+}
+
+//______________________________________________________________________________
+void AliAnalysisTask::FinishTaskOutput()
+{
+// Optional method that is called in SlaveTerminate phase. 
+// Used for calling aditional methods just after the last event was processed ON
+// THE WORKING NODE. The call is made also in local case.
+// Do NOT delete output objects here since they will have to be sent for 
+// merging in PROOF mode - use class destructor for cleanup.
+}
+      
 //______________________________________________________________________________
 void AliAnalysisTask::ConnectInputData(Option_t *)
 {
@@ -409,7 +437,7 @@ void AliAnalysisTask::CreateOutputObjects()
 }
 
 //______________________________________________________________________________
-void AliAnalysisTask::OpenFile(Int_t iout, Option_t *option) const
+TFile *AliAnalysisTask::OpenFile(Int_t iout, Option_t *option) const
 {
 // This method has to be called INSIDE the user redefined CreateOutputObjects
 // method, before creating each object corresponding to the output containers
@@ -431,11 +459,31 @@ void AliAnalysisTask::OpenFile(Int_t iout, Option_t *option) const
 //    fHist2 = new TH2F("my quality check hist2",...);
 // }
    
-   if (iout<0 || iout>=fNoutputs) return;
+   if (iout<0 || iout>=fNoutputs) {
+      Error("OpenFile", "No output slot for task %s with index %d", GetName(), iout);
+      return NULL;
+   }   
+   // We allow file opening also on the slaves (AG)
    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
-   if (!mgr || mgr->GetAnalysisType()==AliAnalysisManager::kProofAnalysis) return;
    AliAnalysisDataContainer *cont = GetOutputSlot(iout)->GetContainer();
-   if (strlen(cont->GetFileName())) new TFile(cont->GetFileName(), option);
+   TFile *f = NULL;
+   if (!strlen(cont->GetFileName())) {
+      Error("OpenFile", "No file name specified for container %s", cont->GetName());
+      return f;
+   }   
+   if (mgr->GetAnalysisType()==AliAnalysisManager::kProofAnalysis && cont->IsSpecialOutput())
+      f = mgr->OpenProofFile(cont->GetFileName(),option);
+   else {
+      // Check first if the file is already opened
+      f = (TFile*)gROOT->GetListOfFiles()->FindObject(cont->GetFileName());
+      if (!f) f = new TFile(cont->GetFileName(), option);
+   }   
+   if (f && !f->IsZombie()) {
+      cont->SetFile(f);
+      return f;
+   }
+   cont->SetFile(NULL);
+   return NULL;
 }
 
 //______________________________________________________________________________
@@ -445,6 +493,13 @@ Bool_t AliAnalysisTask::Notify()
    return kTRUE;
 }
 
+//______________________________________________________________________________
+Bool_t AliAnalysisTask::NotifyBinChange()
+{
+// Overload this IF you need to treat bin change in event mixing.
+   return kTRUE;
+}
+
 //______________________________________________________________________________
 void AliAnalysisTask::Terminate(Option_t *)
 {
@@ -522,7 +577,7 @@ void AliAnalysisTask::PrintTask(Option_t *option, Int_t indent) const
    AliAnalysisDataContainer *cont;
    for (Int_t i=0; i<indent; i++) ind += " ";
    if (!dep || (dep && IsChecked())) {
-      printf("%s\n", Form("%stask: %s  ACTIVE=%i", ind.Data(), GetName(),IsActive()));
+      printf("%s\n", Form("%stask: %s  ACTIVE=%i POST_LOOP=%i", ind.Data(), GetName(),IsActive(),IsPostEventLoop()));
       if (dep) thistask->SetChecked(kFALSE);
       else {
          for (islot=0; islot<fNinputs; islot++) {
@@ -557,16 +612,16 @@ void AliAnalysisTask::PrintContainers(Option_t *option, Int_t indent) const
 }
 
 //______________________________________________________________________________
-void AliAnalysisTask::SetExecPerEvent(Bool_t flag)
+void AliAnalysisTask::SetPostEventLoop(Bool_t flag)
 {
-// Set the task execution mode - run in a event loop or single shot. All output
+// Set the task execution mode - run after event loop or not. All output
 // containers of this task will get the same type.
-   TObject::SetBit(kTaskEvtByEvt,flag);
+   TObject::SetBit(kTaskPostEventLoop,flag);
    AliAnalysisDataContainer *cont;
    Int_t islot;
    for (islot=0; islot<fNoutputs; islot++) {
       cont = GetOutputSlot(islot)->GetContainer();
-      if (cont) cont->SetEventByEvent(flag);
+      if (cont) cont->SetPostEventLoop(flag);
    }   
 }