updates to TRD trigger analysis

diff --git a/PWGJE/UserTasks/AliAnalysisTaskJetsTriggerTRD.cxx b/PWGJE/UserTasks/AliAnalysisTaskJetsTriggerTRD.cxx
index fe79d78218b4f5033807e5db24cfb015126b8d1a..4a8d791ac964fc9c3a033bbbd8ee53fea941a2ec 100644 (file)
--- a/PWGJE/UserTasks/AliAnalysisTaskJetsTriggerTRD.cxx
+++ b/PWGJE/UserTasks/AliAnalysisTaskJetsTriggerTRD.cxx
@@ -9,6 +9,7 @@
 #include "AliAnalysisManager.h"
 #include "AliInputEventHandler.h"
 #include "AliVEvent.h"
+#include "AliVTrdTrack.h"
 
 // MC stuff
 #include "AliMCEvent.h"
@@ -33,14 +34,13 @@
 
 #include "AliAnalysisTaskJetsTriggerTRD.h"
 
-const char * const histos[] = { "stat", "jetpt", "jetpt_emc", "jetpt_trd" };
-
 AliAnalysisTaskJetsTriggerTRD::AliAnalysisTaskJetsTriggerTRD(const char *name) :
   AliAnalysisTaskSE(name),
+  fTriggerMask(0),
   fOutputList(),
   fHist(),
   fShortTaskId("jets_trg_trd"),
-  fNoJetPtBins(40),
+  fNoJetPtBins(80),
   fJetPtBinMax(400)
 {
   // default ctor
@@ -64,15 +64,34 @@ void AliAnalysisTaskJetsTriggerTRD::UserCreateOutputObjects()
   fOutputList->SetOwner();
 
   // setup histograms
-  AddHistogram(HISTID(kHistStat), "event statistics",
-              10, -.5, 9.5);
-  AddHistogram(HISTID(kHistJetPt), "jet spectrum (|#eta| < 0.5);p_{T} (GeV/c);counts",
+  TH1 *histStat = AddHistogram(ID(kHistStat), "event statistics;;counts",
+                               kStatLast-1, .5, kStatLast-.5);
+  histStat->GetXaxis()->SetBinLabel(ID(kStatSeen));
+  histStat->GetXaxis()->SetBinLabel(ID(kStatUsed));
+  histStat->GetXaxis()->SetBinLabel(ID(kStatMB));
+
+  AddHistogram(ID(kHistNoJets), "number of jets;N^{jet};counts",
+               100, -.5, 99.5);
+
+  AddHistogram(ID(kHistTrackGTU), "GTU track p_{T};p_{T} (GeV/c);counts",
+               100, 0., 25.);
+
+  AddHistogram(ID(kHistJetPt), "jet spectrum (|#eta| < 0.5);p_{T} (GeV/c);counts",
+              fNoJetPtBins, 0., fJetPtBinMax);
+  AddHistogram(ID(kHistJetPtITS), "jet spectrum (|#eta| < 0.5);p_{T} (GeV/c);counts",
               fNoJetPtBins, 0., fJetPtBinMax);
-  AddHistogram(HISTID(kHistJetPtEMC), "jet spectrum (|#eta| < 0.5);p_{T} (GeV/c);counts",
+  AddHistogram(ID(kHistJetPt3x3), "jet spectrum (|#eta| < 0.5);p_{T} (GeV/c);counts",
               fNoJetPtBins, 0., fJetPtBinMax);
-  AddHistogram(HISTID(kHistJetPtTRD), "jet spectrum (|#eta| < 0.5);p_{T} (GeV/c);counts",
+
+  AddHistogram(ID(kHistJetPtEMC), "jet spectrum (|#eta| < 0.5);p_{T} (GeV/c);counts",
+              fNoJetPtBins, 0., fJetPtBinMax);
+  AddHistogram(ID(kHistJetPtHJT), "jet spectrum (|#eta| < 0.5);p_{T} (GeV/c);counts",
               fNoJetPtBins, 0., fJetPtBinMax);
 
+  AddHistogram(ID(kHistJetPtNoTracks3), "number of tracks above 3 GeV;p_{T}^{jet};no. of tracks",
+               fNoJetPtBins, 0., fJetPtBinMax,
+               20, -.5, 19.5);
+
   PostData(1, fOutputList);
 }
 
@@ -104,7 +123,7 @@ void AliAnalysisTaskJetsTriggerTRD::UserExec(Option_t * /* option */)
 
   // Int_t nTracksMC  = 0; // no. of MC tracks
   // Int_t nTracksESD = 0; // no. of global ESD tracks
-  // Int_t nTracksGTU = 0; // no. of GTU tracks
+  Int_t nTracksGTU = 0; // no. of GTU tracks
 
   // Int_t nTracks[6][90]; // tracks above lower pt threshold, counted stack-wise
   // memset(nTracks, 0, sizeof(Int_t)*6*90);
@@ -117,27 +136,21 @@ void AliAnalysisTaskJetsTriggerTRD::UserExec(Option_t * /* option */)
   // Float_t leadingJetPtMC = 0.; // leading jet energy from MC information
   Float_t leadingJetPtRec = 0.; // leading jet energy from AOD information
 
-  // number of sampled events
-  FillH1(kHistStat, 1);
+  // record number of sampled events and detect trigger contributions
+  FillH1(kHistStat, kStatSeen);
+  if (!DetectTriggers()) {
+    AliError("Failed to detect the triggers");
+    return;
+  }
 
-  // event selection
-  AliVEvent::EOfflineTriggerTypes physSel = (AliVEvent::EOfflineTriggerTypes) ((AliInputEventHandler*) (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
-  if (fDebug > 2)
-    printf("trg: %8s %8s %8s %8s %8s %8s (%s)\n",
-           (physSel & AliVEvent::kAnyINT)  ? "kAnyINT"  : " ",
-           (physSel & AliVEvent::kCINT5)   ? "kCINT5" : " ",
-           (physSel & AliVEvent::kINT7)    ? "kINT7" : " ",
-           (physSel & AliVEvent::kINT8)    ? "kINT8" : " ",
-           (physSel & AliVEvent::kMB)      ? "kMB"   : " ",
-           (physSel & AliVEvent::kEMC7)    ? "kEMC7" : " ",
-           esdEvent ? esdEvent->GetFiredTriggerClasses().Data() : "<no ESD event>"
-           );
-  // Int_t triggerMask = 0;
-  if (!(physSel & (AliVEvent::kAnyINT)))
+  if (!IsTrigger(kTrgInt))
     return;
+  FillH1(kHistStat, kStatUsed);
 
-  // store no. of sampled events
-  FillH1(kHistStat, 2);
+  for (Int_t iTrack = 0; iTrack < nTracksGTU; ++iTrack) {
+    AliVTrdTrack *trk = 0x0;
+    FillH1(kHistTrackGTU, trk->Pt());
+  }
 
   AliAODJet *leadJet = 0x0;
   AliAODJet *subleadJet = 0x0;
@@ -146,29 +159,33 @@ void AliAnalysisTaskJetsTriggerTRD::UserExec(Option_t * /* option */)
     TClonesArray *jetArray = dynamic_cast<TClonesArray*> (aodEvent->FindListObject(fJetBranchName));
     if (jetArray) {
       Int_t nJets = jetArray->GetEntriesFast();
+      FillH1(kHistNoJets, nJets);
       for (Int_t iJet = 0; iJet < nJets; ++iJet) {
         AliAODJet *jet = (AliAODJet*) (*jetArray)[iJet];
 
         // check contributing tracks
         Int_t nJetTracks = jet->GetRefTracks()->GetEntriesFast();
+        Int_t nJetTracks3 = 0;
         Int_t iLeadingTrack = -1;
         Float_t ptLeadingTrack = 0.;
         for (Int_t iTrack=0; iTrack < nJetTracks; ++iTrack) {
           AliAODTrack *track = (AliAODTrack*) jet->GetRefTracks()->At(iTrack);
+          // count tracks above 3 GeV/c
+          if (track->Pt() > 3.)
+            ++nJetTracks3;
+          // find the leading track
           if (track->Pt() > ptLeadingTrack) {
             ptLeadingTrack = track->Pt();
             iLeadingTrack = iTrack;
           }
         }
 
-        // check tracking flags for leading track
-        // discard the jet if the leading track has no ITS contribution with a hit in any SPD layer
+        // retrieve the leading track
         AliAODTrack *leadingTrack = (AliAODTrack*) jet->GetRefTracks()->At(iLeadingTrack);
-        if (!(leadingTrack->GetFlags() & AliVTrack::kITSrefit) ||
-            !(leadingTrack->HasPointOnITSLayer(0) || leadingTrack->HasPointOnITSLayer(1)))
-          continue;
 
         if (TMath::Abs(jet->Eta()) < 0.5) {
+
+          // find leading jet
           if (TMath::Abs(jet->Pt()) > leadingJetPtRec) {
             leadingJetPtRec = TMath::Abs(jet->Pt());
             subleadJet = leadJet;
@@ -176,11 +193,22 @@ void AliAnalysisTaskJetsTriggerTRD::UserExec(Option_t * /* option */)
           }
 
          FillH1(kHistJetPt, jet->Pt());
-          // if (esdEvent && esdEvent->IsTriggerClassFired("CEMC7WU-B-NOPF-ALL"))
-         //   FillH1("jetpt_emc", jet->Pt());
-          // // if (esdEvent && esdEvent->IsTriggerClassFired("CINT7WU-I-NOPF-ALL"))
-          // if (triggerMask & (1 << kMBHJT))
-         //   FillH1("jetpt_trd", jet->Pt());
+
+          // discard the jet if the leading track has no ITS contribution with a hit in any SPD layer
+          if (leadingTrack &&
+              (leadingTrack->GetFlags() & AliVTrack::kITSrefit) &&
+              (leadingTrack->HasPointOnITSLayer(0) || leadingTrack->HasPointOnITSLayer(1)))
+            FillH1(kHistJetPtITS, jet->Pt());
+
+          // discard the jet if it does not have 3 tracks above 3 GeV/c
+          FillH2(kHistJetPtNoTracks3, jet->Pt(), nJetTracks3);
+          if (nJetTracks3 > 2)
+            FillH1(kHistJetPt3x3, jet->Pt());
+
+          if (IsTrigger(kTrgEMC))
+            FillH1(kHistJetPtEMC, jet->Pt());
+          if (IsTrigger(kTrgHJT))
+            FillH1(kHistJetPtHJT, jet->Pt());
         }
       }
     }
@@ -198,6 +226,38 @@ void AliAnalysisTaskJetsTriggerTRD::Terminate(const Option_t * /* option */)
 
 }
 
+Bool_t AliAnalysisTaskJetsTriggerTRD::DetectTriggers()
+{
+  fTriggerMask = 0;
+
+  AliVEvent::EOfflineTriggerTypes physSel = (AliVEvent::EOfflineTriggerTypes) ((AliInputEventHandler*) (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
+  TString trgClasses = InputEvent()->GetFiredTriggerClasses();
+
+  if (fDebug > 2)
+    printf("trg: %8s %8s %8s %8s %8s %8s %8s (%s)\n",
+           (physSel & AliVEvent::kAnyINT)  ? "kAnyINT"  : " ",
+           (physSel & AliVEvent::kCINT5)   ? "kCINT5" : " ",
+           (physSel & AliVEvent::kINT7)    ? "kINT7" : " ",
+           (physSel & AliVEvent::kINT8)    ? "kINT8" : " ",
+           (physSel & AliVEvent::kMB)      ? "kMB"   : " ",
+           (physSel & AliVEvent::kEMC7)    ? "kEMC7" : " ",
+           (physSel & AliVEvent::kTRD)     ? "kTRD" : " ",
+           trgClasses.Data()
+           );
+
+  // physics selection
+  if (physSel & (AliVEvent::kAnyINT))
+    MarkTrigger(kTrgInt);
+
+  // trigger classes
+  if (trgClasses.Contains("CINT7WUHJT") || trgClasses.Contains("CINT8WUHJT"))
+    MarkTrigger(kTrgHJT);
+  if (trgClasses.Contains("CEMC7WU-") || trgClasses.Contains("CEMC8WU-"))
+    MarkTrigger(kTrgEMC);
+
+  return kTRUE;
+}
+
 // ----- histogram management -----
 TH1* AliAnalysisTaskJetsTriggerTRD::AddHistogram(Hist_t hist, const char *hid, TString title,
                                                 Int_t xbins, Float_t xmin, Float_t xmax,

diff --git a/PWGJE/UserTasks/AliAnalysisTaskJetsTriggerTRD.h b/PWGJE/UserTasks/AliAnalysisTaskJetsTriggerTRD.h
index 49ce90157ad8c3c91ad2f0434ed14ad9ee155abe..7e1549fb363a4a77ce011595693653556f2dadd0 100644 (file)
--- a/PWGJE/UserTasks/AliAnalysisTaskJetsTriggerTRD.h
+++ b/PWGJE/UserTasks/AliAnalysisTaskJetsTriggerTRD.h
@@ -6,7 +6,7 @@
 
 #include "AliAnalysisTaskSE.h"
 
-#define HISTID(x) x, #x
+#define ID(x) x, #x
 
 class TList;
 
@@ -30,14 +30,46 @@ public:
   void    SetJetPtBinMax(Float_t ptmax) { fJetPtBinMax = ptmax; }
   Float_t GetJetPtBinMax() const { return fJetPtBinMax; }
 
-  void SetJetBranchName(const char* const branchName) { strncpy(fJetBranchName, branchName, fgkStringLength); }
+  void SetJetBranchName(const char* const branchName) { strncpy(fJetBranchName, branchName, fgkStringLength-1); }
   const char* GetJetBranchName() const { return fJetBranchName; }
 
-  // // histograms
-  enum Hist_t { kHistStat = 0, kHistJetPt, kHistJetPtEMC, kHistJetPtTRD, kHistLast };
+  // histograms
+  enum Hist_t {
+      kHistStat = 0,
+      kHistNoJets,
+      kHistTrackGTU,
+      kHistJetPt, kHistJetPtITS, kHistJetPt3x3,
+      kHistJetPtEMC, kHistJetPtHJT,
+      kHistJetPtNoTracks3,
+      kHistLast
+  };
+
+  // statistics
+  enum Stat_t {
+      kStatSeen = 1,
+      kStatUsed,
+      kStatMB,
+      kStatLast
+  };
+
+  // trigger conditions
+  enum Trigger_t { 
+      kTrgMB = 0,
+      kTrgInt,
+      kTrgInt78,
+      kTrgHJT,
+      kTrgEMC,
+      kTrgLast
+  };      
 
 protected:
-  // // output objects
+  UInt_t fTriggerMask;         // internal representation of trigger conditions
+
+  Bool_t DetectTriggers();
+  void   MarkTrigger(Trigger_t trg) { fTriggerMask |= (1 << trg); }
+  Bool_t IsTrigger(Trigger_t trg) const { return (fTriggerMask & (1 << trg)); }
+
+  // output objects
   TList *fOutputList;          // list of output objects
 
   // histogram management

diff --git a/PWGJE/macros/AddTaskJetsTriggerTRD.C b/PWGJE/macros/AddTaskJetsTriggerTRD.C
index 26e271176d06dcd6cd089b0a21097fad23e0473c..3b5e7c12c148bf218d6543575c44e8d153b3222e 100644 (file)
--- a/PWGJE/macros/AddTaskJetsTriggerTRD.C
+++ b/PWGJE/macros/AddTaskJetsTriggerTRD.C
@@ -34,8 +34,8 @@ AliAnalysisTaskJetsTriggerTRD* AddTaskJetsTriggerTRD(const char *name = "jets_tr
   mgr->AddTask(task);
 
   AliAnalysisDataContainer *coutput =
-    mgr->CreateContainer("histograms", TList::Class(), AliAnalysisManager::kOutputContainer,
-                         Form("%s:PWGJE_trg_trd", AliAnalysisManager::GetCommonFileName()));
+    mgr->CreateContainer(Form("hist_%s", name), TList::Class(), AliAnalysisManager::kOutputContainer,
+                         Form("%s:PWGJE_jets_trg_trd", AliAnalysisManager::GetCommonFileName()));
 
   mgr->ConnectInput (task, 0, mgr->GetCommonInputContainer());
   // mgr->ConnectOutput(task, 0, mgr->GetCommonOutputContainer());