new task to find jet trigger candidates (loop over EMCal L1 patches over threshold...
authorloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 31 Jul 2013 07:47:32 +0000 (07:47 +0000)
committerloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 31 Jul 2013 07:47:32 +0000 (07:47 +0000)
PWGJE/EMCALJetTasks/AliJetTriggerSelectionTask.cxx [new file with mode: 0644]
PWGJE/EMCALJetTasks/AliJetTriggerSelectionTask.h [new file with mode: 0644]
PWGJE/EMCALJetTasks/macros/AddTaskJetTriggerSelection.C [new file with mode: 0644]

diff --git a/PWGJE/EMCALJetTasks/AliJetTriggerSelectionTask.cxx b/PWGJE/EMCALJetTasks/AliJetTriggerSelectionTask.cxx
new file mode 100644 (file)
index 0000000..7916447
--- /dev/null
@@ -0,0 +1,185 @@
+// $Id$
+//
+// Jet trigger selection task.
+//
+// Author: S.Aiola
+
+#include <TLorentzVector.h>
+#include <TMath.h>
+#include <TF1.h>
+
+#include "AliVEvent.h"
+#include "AliVCluster.h"
+#include "AliEmcalJet.h"
+#include "AliLog.h"
+#include "AliVVZERO.h"
+#include "AliESDUtils.h"
+
+#include "AliJetTriggerSelectionTask.h"
+
+ClassImp(AliJetTriggerSelectionTask)
+
+//________________________________________________________________________
+AliJetTriggerSelectionTask::AliJetTriggerSelectionTask() : 
+  AliAnalysisTaskEmcalJetDev("AliJetTriggerSelectionTask", kFALSE),
+  fEnergyThreshold(0),
+  fMaxDistance2(0.0225),
+  fTriggerBits(AliVEvent::kEMCEJE),
+  fTaskSettingsOk(kFALSE),
+  fNTriggers(0),
+  fVZERO(0),
+  fV0ATotMult(0),
+  fV0CTotMult(0)
+{
+  // Default constructor.
+  
+  for (Int_t i = 0; i < 999; i++) {
+    fTrigPos[i][0] = -999;
+    fTrigPos[i][1] = -999;
+  }
+}
+
+//________________________________________________________________________
+AliJetTriggerSelectionTask::AliJetTriggerSelectionTask(const char *name) : 
+  AliAnalysisTaskEmcalJetDev(name, kFALSE),
+  fEnergyThreshold(0),
+  fMaxDistance2(0.0225),
+  fTriggerBits(AliVEvent::kEMCEJE),
+  fTaskSettingsOk(kFALSE),
+  fNTriggers(0),
+  fVZERO(0),
+  fV0ATotMult(0),
+  fV0CTotMult(0)
+{
+  // Standard constructor.
+
+  for (Int_t i = 0; i < 999; i++) {
+    fTrigPos[i][0] = -999;
+    fTrigPos[i][1] = -999;
+  }
+}
+
+//________________________________________________________________________
+void AliJetTriggerSelectionTask::ExecOnce()
+{
+  // Initialize the task.
+
+  AliAnalysisTaskEmcalJetDev::ExecOnce();
+
+  fTaskSettingsOk = kTRUE;
+
+  fVZERO = InputEvent()->GetVZEROData();
+  if (!fVZERO) {
+    AliError(Form("%s: AliVVZERO not available, task will not be executed!",GetName()));
+    fTaskSettingsOk = kFALSE;
+  }
+  
+  if (GetClusterArray() == 0) {
+    AliError(Form("%s: No cluster collection provided, task will not be executed!",GetName()));
+    fTaskSettingsOk = kFALSE;
+  }
+
+  if (GetJetArray() == 0) {
+    AliError(Form("%s: No jet collection provided, task will not be executed!",GetName()));
+    fTaskSettingsOk = kFALSE;
+  }
+  if (!fEnergyThreshold) {
+    AliError(Form("%s: No threshold function provided, task will not be executed!",GetName()));
+    fTaskSettingsOk = kFALSE;
+  }
+}
+
+//________________________________________________________________________
+Bool_t AliJetTriggerSelectionTask::RetrieveEventObjects()
+{
+  // Retrieve event objects.
+
+  if (!AliAnalysisTaskEmcalJetDev::RetrieveEventObjects())
+    return kFALSE;
+
+  if (fVZERO) {
+    fV0ATotMult = AliESDUtils::GetCorrV0A(fVZERO->GetMTotV0A(),fVertex[2]);
+    fV0CTotMult = AliESDUtils::GetCorrV0C(fVZERO->GetMTotV0C(),fVertex[2]);
+  }
+
+  return kTRUE;
+}
+
+//________________________________________________________________________
+Bool_t AliJetTriggerSelectionTask::Run()
+{
+  // Run the analysis.
+
+  if (!fTaskSettingsOk) return kFALSE;
+
+  FindTriggers();
+  SelectJets();
+  
+  return kTRUE;
+}
+
+//________________________________________________________________________
+void AliJetTriggerSelectionTask::FindTriggers()
+{
+  fNTriggers = 0;
+
+  Int_t nclusters = GetNClusters();
+
+  Double_t th = fEnergyThreshold->Eval(fV0ATotMult+fV0CTotMult);
+
+  for (Int_t i = 0; i < nclusters; i++) {
+    if (fNTriggers >= 999) {
+      AliError("More than 999 triggers found!");
+      break;
+    }
+
+    AliVCluster *cluster = GetAcceptClusterFromArray(i);
+    if (!cluster)
+      continue;
+
+    //Printf("Cluster energy=%.3f, th=%.3f",cluster->E(), th);
+
+    if (cluster->E() > th) {
+      TLorentzVector vect;
+      cluster->GetMomentum(vect,fVertex);
+      fTrigPos[fNTriggers][0] = vect.Eta();
+      fTrigPos[fNTriggers][1] = vect.Phi();
+      fNTriggers++;
+    }
+  }
+
+  AliDebug(2,Form("%s: %d triggers found among %d candidates (cent=%.1f, mult=%.1f, th=%.2f)!",GetName(),fNTriggers,nclusters,fCent,fV0ATotMult+fV0CTotMult,th));
+}
+
+//________________________________________________________________________
+void AliJetTriggerSelectionTask::SelectJets()
+{
+  for (Int_t c = 0; c < fJetCollArray.GetEntriesFast(); c++) {
+
+    Int_t njets = GetNJets(c);
+
+    for (Int_t i = 0; i < njets; i++) {
+      AliEmcalJet *jet = GetAcceptJetFromArray(i,c);
+      if (IsTriggerJet(jet)) jet->AddTrigger(fTriggerBits);
+    }
+  }
+}
+
+//________________________________________________________________________
+Bool_t AliJetTriggerSelectionTask::IsTriggerJet(AliEmcalJet *jet)
+{
+  if (!jet) return kFALSE;
+
+  for (Int_t i = 0; i < fNTriggers; i++) {
+    Double_t deta = jet->Eta() - fTrigPos[i][0];
+    Double_t dphi = jet->Phi() - fTrigPos[i][1];
+    
+    Double_t d2 = deta * deta + dphi * dphi;
+
+    if (d2 < fMaxDistance2)
+      return kTRUE;
+  }
+  
+  return kFALSE;
+}
diff --git a/PWGJE/EMCALJetTasks/AliJetTriggerSelectionTask.h b/PWGJE/EMCALJetTasks/AliJetTriggerSelectionTask.h
new file mode 100644 (file)
index 0000000..b4f7d7d
--- /dev/null
@@ -0,0 +1,46 @@
+#ifndef ALIJETTRIGGERSELECTIONTASK_H
+#define ALIJETTRIGGERSELECTIONTASK_H
+
+// $Id$
+
+class AliEmcalJet;
+
+#include "AliAnalysisTaskEmcalJetDev.h"
+
+class AliJetTriggerSelectionTask : public AliAnalysisTaskEmcalJetDev {
+ public:
+
+  AliJetTriggerSelectionTask();
+  AliJetTriggerSelectionTask(const char *name);
+  virtual ~AliJetTriggerSelectionTask() {;}
+
+  void                        SetMaxDistance(Double_t d) { fMaxDistance2    = d*d ; }
+  void                        SetEnergyThreshold(TF1 *f) { fEnergyThreshold = f   ; }
+  void                        SetTriggerBits(UInt_t d)   { fTriggerBits     = d   ; }
+
+ protected:
+  Bool_t                      Run();
+  void                        ExecOnce();
+  Bool_t                      RetrieveEventObjects();
+  void                        FindTriggers();
+  void                        SelectJets();
+  Bool_t                      IsTriggerJet(AliEmcalJet *jet);
+
+  TF1                        *fEnergyThreshold;                // energy threshold vs. centrality
+  Double_t                    fMaxDistance2;                   // max distance square between trigger patch and jet
+  UInt_t                      fTriggerBits;                    // trigger bit to be set
+
+  Bool_t                      fTaskSettingsOk;                 //!if false, don't execute task  
+  Int_t                       fNTriggers;                      //!number of triggers in the current event
+  Double_t                    fTrigPos[999][2];                //!(eta,phi) trigger positions in the current event
+  AliVVZERO                  *fVZERO;                          //!Event V0 object
+  Double_t                    fV0ATotMult;                     //!Event V0A total multiplicity
+  Double_t                    fV0CTotMult;                     //!Event V0C total multiplicity
+ private:
+  AliJetTriggerSelectionTask(const AliJetTriggerSelectionTask&);            // not implemented
+  AliJetTriggerSelectionTask &operator=(const AliJetTriggerSelectionTask&); // not implemented
+
+  ClassDef(AliJetTriggerSelectionTask, 1) // jet trigger selection task
+};
+#endif
diff --git a/PWGJE/EMCALJetTasks/macros/AddTaskJetTriggerSelection.C b/PWGJE/EMCALJetTasks/macros/AddTaskJetTriggerSelection.C
new file mode 100644 (file)
index 0000000..71b4b4e
--- /dev/null
@@ -0,0 +1,51 @@
+// $Id$
+
+AliJetTriggerSelectionTask* AddTaskJetTriggerSelection(
+  const char *nclusters          = "CaloClusters",
+  TF1        *eth                = 0,
+  Double_t    maxdistance        = 0.15,
+  UInt_t      triggerbits        = AliVEvent::kEMCEJE,
+  const char *taskname           = "AliJetTriggerSelectionTask"
+)
+{  
+  // Get the pointer to the existing analysis manager via the static access method.
+  //==============================================================================
+  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+  if (!mgr)
+  {
+    ::Error("AddTaskJetTriggerSelection", "No analysis manager to connect to.");
+    return NULL;
+  }  
+  
+  // Check the analysis type using the event handlers connected to the analysis manager.
+  //==============================================================================
+  if (!mgr->GetInputEventHandler())
+  {
+    ::Error("AddTaskJetTriggerSelection", "This task requires an input event handler");
+    return NULL;
+  }
+  
+  //-------------------------------------------------------
+  // Init the task and do settings
+  //-------------------------------------------------------
+
+  TString name(Form("%s_%s",taskname, nclusters));
+
+  AliJetTriggerSelectionTask* jetTask = new AliJetTriggerSelectionTask(name);
+  jetTask->AddClusterContainer(nclusters);
+  jetTask->SetEnergyThreshold(eth);  
+  jetTask->SetMaxDistance(maxdistance);
+  jetTask->SetTriggerBits(triggerbits);
+
+  //-------------------------------------------------------
+  // Final settings, pass to manager and set the containers
+  //-------------------------------------------------------
+  
+  mgr->AddTask(jetTask);
+  
+  // Create containers for input/output
+  AliAnalysisDataContainer *cinput1  = mgr->GetCommonInputContainer()  ;
+  mgr->ConnectInput  (jetTask, 0,  cinput1 );
+  
+  return jetTask;
+}