- added analysis task for comparing offline HLT PHOS properties
authorkkanaki <kkanaki@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 9 Apr 2010 14:04:51 +0000 (14:04 +0000)
committerkkanaki <kkanaki@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 9 Apr 2010 14:04:51 +0000 (14:04 +0000)
HLT/QA/tasks/AliAnalysisTaskHLTPHOS.cxx [new file with mode: 0644]
HLT/QA/tasks/AliAnalysisTaskHLTPHOS.h [new file with mode: 0644]
HLT/QA/tasks/macros/compare-HLT-offline-local.C

diff --git a/HLT/QA/tasks/AliAnalysisTaskHLTPHOS.cxx b/HLT/QA/tasks/AliAnalysisTaskHLTPHOS.cxx
new file mode 100644 (file)
index 0000000..78d91a4
--- /dev/null
@@ -0,0 +1,312 @@
+// $Id$\r
+\r
+//**************************************************************************\r
+//* This file is property of and copyright by the ALICE HLT Project        *\r
+//* ALICE Experiment at CERN, All rights reserved.                         *\r
+//*                                                                        *\r
+//* Primary Authors: Zhongbao Yin <zbyin@mail.ccnu.edu.cn>,                *\r
+//*                  Kalliopi Kanaki <Kalliopi.Kanaki@ift.uib.no>          *\r
+//*                  for The ALICE HLT Project.                            *\r
+//*                                                                        *\r
+//* Permission to use, copy, modify and distribute this software and its   *\r
+//* documentation strictly for non-commercial purposes is hereby granted   *\r
+//* without fee, provided that the above copyright notice appears in all   *\r
+//* copies and that both the copyright notice and this permission notice   *\r
+//* appear in the supporting documentation. The authors make no claims     *\r
+//* about the suitability of this software for any purpose. It is          *\r
+//* provided "as is" without express or implied warranty.                  *\r
+//**************************************************************************\r
+\r
+/** @file   AliAnalysisTaskHLTPHOS.cxx\r
+    @author Zhongbao Yin, Kalliopi Kanaki\r
+    @date\r
+    @brief\r
+*/\r
+\r
+#include <iostream>\r
+\r
+#include "TChain.h"\r
+#include "TTree.h"\r
+#include "TH1F.h"\r
+#include "TH2F.h"\r
+#include "TCanvas.h"\r
+#include "TVector3.h"\r
+#include "TString.h"\r
+#include "TObjArray.h"\r
+#include "TFile.h"\r
+\r
+#include "AliESDEvent.h"\r
+#include "AliESDRun.h"\r
+#include "AliESDInputHandler.h"\r
+\r
+#include "AliAnalysisTask.h"\r
+#include "AliAnalysisManager.h"\r
+#include "AliAnalysisTaskHLTPHOS.h"\r
+\r
+\r
+ClassImp(AliAnalysisTaskHLTPHOS)\r
+\r
+//===========================================================================================\r
+\r
+AliAnalysisTaskHLTPHOS::AliAnalysisTaskHLTPHOS(const char *name)\r
+    : \r
+     AliAnalysisTaskSE(name)\r
+    ,fESDRun(0)\r
+    ,fOutputList(0)\r
+    ,fHistOnlTrk2PHOS(0)\r
+    ,fHistOfflTrk2PHOS(0)\r
+    ,fHistOfflTrk2PHOSTrig(0)\r
+    ,fHistOfflTrk2PHOSNoTrig(0)\r
+    ,fNevt(0)\r
+    ,fTrgClsArray(0)\r
+{\r
+  // Constructor\r
+\r
+  // Define input and output slots here\r
+  // Input slot #0 works with a TChain\r
+  // DefineInput(0, TChain::Class());\r
+  // Output slot #0 writes into a TH1 container\r
+\r
+  DefineOutput(1, TList::Class());\r
+}\r
+\r
+const Float_t AliAnalysisTaskHLTPHOS::fgkEtaMin = -0.12;  \r
+const Float_t AliAnalysisTaskHLTPHOS::fgkEtaMax =  0.12;  \r
+const Float_t AliAnalysisTaskHLTPHOS::fgkPhiMin[5]   = {3.83972, 4.18879, 4.53786, 4.88692, 5.23599};  \r
+const Float_t AliAnalysisTaskHLTPHOS::fgkPhiMax[5]   = {4.18879, 4.53786, 4.88692, 5.23599, 5.58505};  \r
+const Float_t AliAnalysisTaskHLTPHOS::fgkNormX[5]    = {-0.642788, -0.34202, 0, 0.34202, 0.642788};  \r
+const Float_t AliAnalysisTaskHLTPHOS::fgkNormY[5]    = {-0.766044, -0.939693, -1, -0.939693, -0.766044};  \r
+const Float_t AliAnalysisTaskHLTPHOS::fgkInitPosX[5] = {-295.682, -157.329, 0, 157.329, 295.682};  \r
+const Float_t AliAnalysisTaskHLTPHOS::fgkInitPosY[5] = {-352.38, -432.259, -460, -432.259, -352.38};\r
+\r
+//----------------------------------------------------------------------------------------------------\r
+\r
+void AliAnalysisTaskHLTPHOS::UserCreateOutputObjects(){\r
+// Create histograms\r
+\r
+  OpenFile(1);\r
+\r
+  fOutputList = new TList();\r
+  fOutputList->SetName(GetName());\r
+\r
+// --------------- define histograms ---------------------//\r
+\r
+  fHistOfflTrk2PHOS = new TH2F("fHistOfflTrk2PHOS","track intersection point in #eta and #phi (offline)",100, -0.5, 0.5, 100, 240, 340); \r
+  fHistOnlTrk2PHOS  = new TH2F("fHistOnlTrk2PHOS", "track intersection point in #eta and #phi (HLT)",    100, -0.5, 0.5, 100, 240, 340);\r
+  \r
+  fHistOfflTrk2PHOSTrig   = new TH2F("fHistOfflTrk2PHOSTrig",  "track intersection point in #eta and #phi (offline triggered)",    100, -0.5, 0.5, 100, 240, 340);  \r
+  fHistOfflTrk2PHOSNoTrig = new TH2F("fHistOfflTrk2PHOSNoTrig","track intersection point in #eta and #phi (offline not triggered)",100, -0.5, 0.5, 100, 240, 340);   \r
+\r
+  // -------------- add histograms to the output TList -----------------//\r
+  \r
+  fOutputList->Add(fHistOnlTrk2PHOS);\r
+  fOutputList->Add(fHistOfflTrk2PHOS);  \r
+  fOutputList->Add(fHistOfflTrk2PHOSTrig);\r
+  fOutputList->Add(fHistOfflTrk2PHOSNoTrig);\r
+\r
+}\r
+\r
+void AliAnalysisTaskHLTPHOS::NotifyRun(){\r
+// This will not work if the active trigger classes change from run to run.\r
+// Then one has to know all trigger classes before processing the data.\r
+\r
+  AliESDEvent* evESD = dynamic_cast<AliESDEvent*>(InputEvent());\r
+  TString trgClasses = evESD->GetESDRun()->GetActiveTriggerClasses(); \r
\r
+  fTrgClsArray = trgClasses.Tokenize(" ");\r
+  //cout<<fTrgClsArray->GetEntries()<<endl; \r
+    \r
+//   for(Int_t i = 0; i < fTrgClsArray->GetEntries(); i++){ \r
+//     TString str = ((TObjString *)fTrgClsArray->At(i))->GetString(); \r
+//     (fHistTrigger->GetXaxis())->SetBinLabel(i+1, str.Data()); \r
+//     (fHistHLTTrigger->GetXaxis())->SetBinLabel(i+1, str.Data()); \r
+//   } \r
+  \r
+  evESD = NULL;\r
+}\r
+\r
+void AliAnalysisTaskHLTPHOS::UserExec(Option_t *){\r
+\r
+  AliESDEvent* evESD = dynamic_cast<AliESDEvent*>(InputEvent());\r
+  \r
+  if (!evESD) {\r
+    Printf("ERROR: fESD not available");\r
+    return;\r
+  }\r
+  \r
+  AliESDEvent* evHLTESD = 0;\r
+  AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(fInputHandler);\r
+   \r
+  if(esdH) evHLTESD = esdH->GetHLTEvent();\r
+    \r
+  if(!evHLTESD){\r
+      Printf("ERROR: HLTesd not available");\r
+      return;\r
+  }\r
+\r
+  Double_t b = evESD->GetMagneticField();\r
+  \r
+  Double_t pos[] = { 0., 0., 0.};\r
+  AliVertex *vtx = new AliVertex(pos, 0., 0);\r
+    \r
+//   for(Int_t i = 0; i<fTrgClsArray->GetEntries(); i++){\r
+//       if((evESD->GetFiredTriggerClasses()).Contains(((TObjString *)fTrgClsArray->At(i))->GetString()))  fHistTrigger->Fill(i);\r
+//   }\r
+// \r
+//   if(evHLTESD->IsHLTTriggerFired()){\r
+//      //fHistHLTTrigger->Fill(evESD->GetTriggerMask());\r
+//      for(Int_t i = 0; i<fTrgClsArray->GetEntries(); i++){ \r
+//          if((evESD->GetFiredTriggerClasses()).Contains(((TObjString *)fTrgClsArray->At(i))->GetString())) fHistHLTTrigger->Fill(i);\r
+//      } \r
+//   }\r
+\r
+  if(evHLTESD->IsHLTTriggerFired()){\r
+     for(Int_t i = 0; i < evHLTESD->GetNumberOfTracks(); i++){\r
+         AliESDtrack * HLTesdTrk = evHLTESD->GetTrack(i);\r
+      \r
+        TVector3 v; \r
+        if(IsInPHOS(2, HLTesdTrk, b, v)){ \r
+         Float_t phi = v.Phi(); \r
+         if(phi<0) phi += 2.*TMath::Pi(); \r
+         fHistOnlTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg()); \r
+        }else if(IsInPHOS(3, HLTesdTrk, b, v)){ \r
+         Float_t phi = v.Phi();  \r
+         if(phi<0) phi += 2.*TMath::Pi();  \r
+         fHistOnlTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());  \r
+        }else if( IsInPHOS(4, HLTesdTrk, b, v) ){ \r
+         Float_t phi = v.Phi();   \r
+         if(phi<0) phi += 2.*TMath::Pi();   \r
+         fHistOnlTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg()); \r
+        } \r
+\r
+/*\r
+      if(IsInPHOS(2, HLTesdTrk, b) \r
+        || IsInPHOS(3, HLTesdTrk, b)\r
+        || IsInPHOS(4, HLTesdTrk, b) ) cout<<"Good Trigger"<<endl;\r
+*/      \r
+       \r
+     }\r
+  }else{\r
+    for(Int_t i = 0; i < evHLTESD->GetNumberOfTracks(); i++){ \r
+        AliESDtrack * HLTesdTrk = evHLTESD->GetTrack(i); \r
+\r
+        TVector3 v;  \r
+        if(IsInPHOS(2, HLTesdTrk, b, v)){  \r
+          Float_t phi = v.Phi();  \r
+          if(phi<0) phi += 2.*TMath::Pi();  \r
+          fHistOnlTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());  \r
+        }else if(IsInPHOS(3, HLTesdTrk, b, v)){  \r
+          Float_t phi = v.Phi();   \r
+          if(phi<0) phi += 2.*TMath::Pi();   \r
+          fHistOnlTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());   \r
+        }else if( IsInPHOS(4, HLTesdTrk, b, v) ){  \r
+          Float_t phi = v.Phi();    \r
+          if(phi<0) phi += 2.*TMath::Pi();    \r
+          fHistOnlTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());  \r
+        }  \r
+    } \r
+  }\r
+  \r
+  if(evHLTESD->IsHLTTriggerFired()){\r
+\r
+    for(Int_t i = 0; i < evESD->GetNumberOfTracks(); i++){ \r
+      AliESDtrack * esdTrk = evESD->GetTrack(i); \r
+\r
+      TVector3 v;\r
+      if(IsInPHOS(2, esdTrk, b, v)){\r
+       Float_t phi = v.Phi();\r
+       if(phi<0) phi += 2.*TMath::Pi();\r
+       fHistOfflTrk2PHOSTrig->Fill(v.Eta(), phi*TMath::RadToDeg());\r
+       fHistOfflTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());\r
+       //cout<<"Event in PHOS 2: "<<fNevt<<endl;\r
+      }else if(IsInPHOS(3, esdTrk, b, v)){\r
+       Float_t phi = v.Phi(); \r
+        if(phi<0) phi += 2.*TMath::Pi(); \r
+        fHistOfflTrk2PHOSTrig->Fill(v.Eta(), phi*TMath::RadToDeg()); \r
+       fHistOfflTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());\r
+       //cout<<"Event in PHOS 3: "<<fNevt<<endl;\r
+      }else if( IsInPHOS(4, esdTrk, b, v) ){\r
+       Float_t phi = v.Phi();  \r
+        if(phi<0) phi += 2.*TMath::Pi();  \r
+        fHistOfflTrk2PHOSTrig->Fill(v.Eta(), phi*TMath::RadToDeg());\r
+       fHistOfflTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());\r
+       //cout<<"Event in PHOS 4: "<<fNevt<<endl;\r
+      }\r
+    }\r
+  } else {\r
+\r
+    for(Int_t i = 0; i < evESD->GetNumberOfTracks(); i++){ \r
+      AliESDtrack * esdTrk = evESD->GetTrack(i); \r
+\r
+      TVector3 v; \r
+      if(IsInPHOS(2, esdTrk, b, v)){ \r
+        Float_t phi = v.Phi(); \r
+        if(phi<0) phi += 2.*TMath::Pi(); \r
+        fHistOfflTrk2PHOSNoTrig->Fill(v.Eta(), phi*TMath::RadToDeg()); \r
+       fHistOfflTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());\r
+       //cout<<"Event in PHOS 2: "<<fNevt<<endl;\r
+      }else if(IsInPHOS(3, esdTrk, b, v)){ \r
+        Float_t phi = v.Phi();  \r
+        if(phi<0) phi += 2.*TMath::Pi();  \r
+        fHistOfflTrk2PHOSNoTrig->Fill(v.Eta(), phi*TMath::RadToDeg());  \r
+       fHistOfflTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());\r
+       //cout<<"Event in PHOS 3: "<<fNevt<<endl;\r
+      }else if( IsInPHOS(4, esdTrk, b, v) ){ \r
+        Float_t phi = v.Phi();   \r
+        if(phi<0) phi += 2.*TMath::Pi();   \r
+        fHistOfflTrk2PHOSNoTrig->Fill(v.Eta(), phi*TMath::RadToDeg()); \r
+       fHistOfflTrk2PHOS->Fill(v.Eta(), phi*TMath::RadToDeg());\r
+       //cout<<"Event in PHOS 4: "<<fNevt<<endl;\r
+      } \r
+    }    \r
+  }\r
+\r
+  fNevt++;\r
+  delete vtx;\r
+\r
+  // Post output data.\r
+  PostData(1, fOutputList);\r
+\r
+ }\r
+\r
+void AliAnalysisTaskHLTPHOS::Terminate(Option_t *){\r
+  /*\r
+  Printf("Number of tracks thru CE: %d", fNtracksThruZ0);\r
+  Printf("Number of tracks thru CE from triggered events: %d", \r
+        fNtracksThruZ0Trig);\r
+  */\r
+\r
+  // Draw result to the screen\r
+  // Called once at the end of the query\r
+\r
+  //  TCanvas *c1 = new TCanvas("AliAnalysisTaskHLTPHOS","Trigger",10,10,510,510);\r
+  //fHistTrigger->DrawCopy("E");\r
+  \r
+}\r
+\r
+Bool_t AliAnalysisTaskHLTPHOS::IsInPHOS(Int_t iMod, AliESDtrack * trk, Float_t b, TVector3& v){\r
+\r
+  Double_t normVector[3] = {fgkNormX[iMod], fgkNormY[iMod], 0};\r
+\r
+  Double_t point[3] = {fgkInitPosX[iMod], fgkInitPosY[iMod], 0};\r
+  \r
+  if(!trk->Intersect(point, normVector, b)) return kFALSE;\r
+\r
+  TVector3 trackPos(point);\r
+  \r
+  v=trackPos;\r
+\r
+  Double_t phi = 0;\r
+  if(trackPos.Phi() < 0) phi = trackPos.Phi() + 2*TMath::Pi();\r
+  else phi = trackPos.Phi();\r
+\r
+  if(trackPos.Eta() >= fgkEtaMin && \r
+     trackPos.Eta() <= fgkEtaMax &&\r
+     phi >= fgkPhiMin[iMod] &&\r
+     phi <= fgkPhiMax[iMod])\r
+    {\r
+      return kTRUE;\r
+    }\r
+  \r
+  return kFALSE;\r
+}\r
diff --git a/HLT/QA/tasks/AliAnalysisTaskHLTPHOS.h b/HLT/QA/tasks/AliAnalysisTaskHLTPHOS.h
new file mode 100644 (file)
index 0000000..fd43ad4
--- /dev/null
@@ -0,0 +1,71 @@
+// $Id$\r
+\r
+#ifndef ALIANALYSISTASKHLTPHOS_H\r
+#define ALIANALYSISTASKHLTPHOS_H\r
+\r
+//* This file is property of and copyright by the ALICE HLT Project *\r
+//* ALICE Experiment at CERN, All rights reserved.                  *\r
+//* See cxx source for full Copyright notice                        *\r
+\r
+/** @file AliAnalysisTaskHLTTPC.h\r
+    @author Zhongbao Yin, Kalliopi Kanaki\r
+    @date\r
+    @brief An analysis task to compare the offline and HLT esd trees\r
+*/\r
+\r
+\r
+// forward declarations\r
+class TH1F;\r
+class TH2F;\r
+class TList;\r
+class AliESDEvent;\r
+class AliESDtrack;\r
+class AliESDRun;\r
+class TObjArray;\r
+\r
+#include "AliAnalysisTaskSE.h"\r
+\r
+class AliAnalysisTaskHLTPHOS : public AliAnalysisTaskSE {\r
\r
+  public:  \r
+    AliAnalysisTaskHLTPHOS(const char *name);\r
+    virtual ~AliAnalysisTaskHLTPHOS() {}\r
+    virtual void   UserCreateOutputObjects();\r
+    virtual void   UserExec(Option_t *option);\r
+    virtual void   Terminate(Option_t *);\r
+    //virtual Bool_t Notify();\r
+    virtual void NotifyRun();\r
+\r
+ private:\r
+\r
+    AliESDRun *fESDRun;  //!\r
+    TList *fOutputList;\r
+\r
+    TH2F *fHistOnlTrk2PHOS; //! track to PHOS 2,3,4 modules in (eta, phi)\r
+    TH2F *fHistOfflTrk2PHOS; //! \r
+    TH2F *fHistOfflTrk2PHOSTrig; //!\r
+    TH2F *fHistOfflTrk2PHOSNoTrig; //!\r
+\r
+    Int_t fNevt;\r
+    TObjArray *fTrgClsArray;\r
+\r
+    static const Float_t fgkPhiMin[5];\r
+    static const Float_t fgkPhiMax[5];\r
+    static const Float_t fgkEtaMin;\r
+    static const Float_t fgkEtaMax;\r
+    static const Float_t fgkNormX[5];\r
+    static const Float_t fgkNormY[5];\r
+    static const Float_t fgkInitPosX[5];\r
+    static const Float_t fgkInitPosY[5];\r
+\r
+    /** copy constructor */\r
+    AliAnalysisTaskHLTPHOS(const AliAnalysisTaskHLTPHOS&); \r
+    /** assignment operator */\r
+    AliAnalysisTaskHLTPHOS& operator=(const AliAnalysisTaskHLTPHOS&); \r
+\r
+    Bool_t IsInPHOS(Int_t iMod, AliESDtrack * trk, Float_t b, TVector3& v);\r
+\r
+    ClassDef(AliAnalysisTaskHLTPHOS, 0);\r
+};\r
+\r
+#endif\r
index 37d17fd..de5a001 100644 (file)
@@ -84,14 +84,23 @@ void compare_HLT_offline_local(){
  
   //-------------- define the tasks ------------//
   
-  AliAnalysisTaskHLTTPC *task1 = new AliAnalysisTaskHLTTPC("offhlt_comparison");
+  AliAnalysisTaskHLTTPC *task1 = new AliAnalysisTaskHLTTPC("offhlt_comparison_TPC");
   mgr->AddTask(task1);
 
-  AliAnalysisDataContainer *coutput1 =  mgr->CreateContainer("histograms", TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-TPC-comparison.root");  
+  AliAnalysisTaskHLTPHOS *taskPHOS = new AliAnalysisTaskHLTPHOS("offhlt_comparison_PHOS");
+  mgr->AddTask(taskPHOS);
+
+  AliAnalysisDataContainer *coutput1 =  mgr->CreateContainer("tpc_histograms", TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-TPC-comparison.root");  
+  AliAnalysisDataContainer *coutput2 =  mgr->CreateContainer("phos_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-PHOS-comparison.root");  
   
   mgr->ConnectInput(task1,0,mgr->GetCommonInputContainer());
   //mgr->ConnectOutput (task1, 0, mgr->GetCommonOutputContainer());
   mgr->ConnectOutput(task1,1,coutput1);
+  
+  mgr->ConnectInput(taskPHOS,0,mgr->GetCommonInputContainer());
+  mgr->ConnectOutput(taskPHOS,1,coutput2);
+  
+  
   if (!mgr->InitAnalysis()) return;
   mgr->PrintStatus();
   mgr->StartAnalysis("local",chain);