From 14e13ab8bc0c62367162494ea67ffab77237b1b8 Mon Sep 17 00:00:00 2001 From: mverweij Date: Thu, 27 Nov 2014 14:15:17 +0100 Subject: [PATCH] add class to analyze mock-up trigger avoid problems with incorrect bit settings in ESD/AOD --- PWG/EMCAL/AliEmcalPatchFromCellMaker.cxx | 4 +- PWGJE/CMakelibPWGJEEMCALJetTasks.pkg | 1 + .../AliAnalysisTaskTriggerRejection.cxx | 195 ++++++++++++++++++ .../AliAnalysisTaskTriggerRejection.h | 54 +++++ .../macros/AddTaskTriggerRejection.C | 190 +++++++++++++++++ PWGJE/PWGJEEMCALJetTasksLinkDef.h | 1 + 6 files changed, 444 insertions(+), 1 deletion(-) create mode 100644 PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskTriggerRejection.cxx create mode 100644 PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskTriggerRejection.h create mode 100644 PWGJE/EMCALJetTasks/macros/AddTaskTriggerRejection.C diff --git a/PWG/EMCAL/AliEmcalPatchFromCellMaker.cxx b/PWG/EMCAL/AliEmcalPatchFromCellMaker.cxx index ae998e60102..d26aab17032 100644 --- a/PWG/EMCAL/AliEmcalPatchFromCellMaker.cxx +++ b/PWG/EMCAL/AliEmcalPatchFromCellMaker.cxx @@ -261,7 +261,9 @@ void AliEmcalPatchFromCellMaker::RunSimpleOfflineTrigger() Int_t patchSize = GetDimFastor(); Int_t stepSize = GetSlidingStepSizeFastor(); Int_t maxCol = kPatchCols - patchSize; - Int_t maxRow = kPatchRows - patchSize; + // Int_t maxRow = kPatchRows - patchSize; + Int_t maxRow = fGeom->GetNTotalTRU()*2 - patchSize; //apparently this is the total TRU in phi ?? + // Printf("fGeom->GetNTotalTRU(): %d %d = %d x %d ; %d x %d",fGeom->GetNTRU(),fGeom->GetNTotalTRU(),fGeom->GetNTRUEta(),fGeom->GetNTRUPhi(),fGeom->GetNModulesInTRUEta(),fGeom->GetNModulesInTRUPhi()); for (Int_t i = 0; i <= maxCol; i += stepSize) { for (Int_t j = 0; j <= maxRow; j += stepSize) { diff --git a/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg b/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg index 76781fc5c81..fd7663f12d2 100644 --- a/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg +++ b/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg @@ -94,6 +94,7 @@ set ( SRCS EMCALJetTasks/UserTasks/AliAnalysisTaskSAQA.cxx EMCALJetTasks/UserTasks/AliAnalysisTaskV0sInJetsEmcal.cxx EMCALJetTasks/UserTasks/AliAnalysisTaskSOH.cxx + EMCALJetTasks/UserTasks/AliAnalysisTaskTriggerRejection.cxx EMCALJetTasks/UserTasks/AliEmcalPicoTrackFromJetMaker.cxx EMCALJetTasks/UserTasks/AliCutValueRange.cxx EMCALJetTasks/UserTasks/AliEMCalHistoContainer.cxx diff --git a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskTriggerRejection.cxx b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskTriggerRejection.cxx new file mode 100644 index 00000000000..882d9e60714 --- /dev/null +++ b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskTriggerRejection.cxx @@ -0,0 +1,195 @@ +// +// Jet trigger rejection analysis task. +// +// Author: M.Verweij + +#include +#include +#include +#include +#include +#include +#include +#include + +#include "AliVVZERO.h" +#include "AliEmcalJet.h" +#include "AliRhoParameter.h" +#include "AliLog.h" +#include "AliJetContainer.h" +#include "AliEmcalTriggerPatchInfo.h" +#include "AliPicoTrack.h" + +#include "AliAnalysisTaskTriggerRejection.h" + +ClassImp(AliAnalysisTaskTriggerRejection) + +//________________________________________________________________________ +AliAnalysisTaskTriggerRejection::AliAnalysisTaskTriggerRejection() : + AliAnalysisTaskEmcalJet("AliAnalysisTaskTriggerRejection", kTRUE), + fContainerFull(0), + fContainerCharged(1), + fMaxPatch(0), + fhnTriggerInfo(0) +{ + // Default constructor. + SetMakeGeneralHistograms(kTRUE); +} + +//________________________________________________________________________ +AliAnalysisTaskTriggerRejection::AliAnalysisTaskTriggerRejection(const char *name) : + AliAnalysisTaskEmcalJet(name, kTRUE), + fContainerFull(0), + fContainerCharged(1), + fMaxPatch(0), + fhnTriggerInfo(0) +{ + // Standard constructor. + SetMakeGeneralHistograms(kTRUE); +} + +//________________________________________________________________________ +AliAnalysisTaskTriggerRejection::~AliAnalysisTaskTriggerRejection() +{ + // Destructor. +} + +//________________________________________________________________________ +void AliAnalysisTaskTriggerRejection::UserCreateOutputObjects() +{ + // Create user output. + + AliAnalysisTaskEmcalJet::UserCreateOutputObjects(); + + Bool_t oldStatus = TH1::AddDirectoryStatus(); + TH1::AddDirectory(kFALSE); + + Int_t fgkNCentBins = 21; + Float_t kMinCent = 0.; + Float_t kMaxCent = 105.; + + Int_t fgkNPtBins = 170; + Float_t kMinPt = -50.; + Float_t kMaxPt = 120.; + + Int_t fgkNVZEROBins = 100; + Float_t kMinVZERO = 0.; + Float_t kMaxVZERO = 25000; + + const Int_t fgkNEPatch = 100; + Float_t kMinEPatch = 0.; + Float_t kMaxEPatch = 200.; + + const Int_t fgkNADC = 100; + Float_t kMinADC = 0.; + Float_t kMaxADC = 1500.; + + const Int_t fgkNEta = 10; + const Int_t fgkNPhi = 10; + + const Int_t nDim = 8;//cent;V0mult;ptjet1;ptjet2;Epatch;ADCpatch;EtaPatch;PhiPatch + const Int_t nBins[nDim] = {fgkNCentBins,fgkNVZEROBins,fgkNPtBins,fgkNPtBins,fgkNEPatch,fgkNADC,fgkNEta,fgkNPhi}; + const Double_t xmin0[nDim] = {kMinCent,kMinVZERO,kMinPt,kMinPt,kMinEPatch,kMinADC,-0.7,1.4}; + const Double_t xmax0[nDim] = {kMaxCent,kMaxVZERO,kMaxPt,kMaxPt,kMaxEPatch,kMaxADC, 0.7,3.14}; + fhnTriggerInfo = new THnSparseF("fhnTriggerInfo", + "hnTriggerInfo;cent;V0mult;ptjet1;ptjet2;Epatch;ADCpatch;EtaPatch;PhiPatch",nDim,nBins,xmin0,xmax0); + fOutput->Add(fhnTriggerInfo); + + // =========== Switch on Sumw2 for all histos =========== + for (Int_t i=0; iGetEntries(); ++i) { + TH1 *h1 = dynamic_cast(fOutput->At(i)); + if (h1){ + h1->Sumw2(); + continue; + } + TH2 *h2 = dynamic_cast(fOutput->At(i)); + if (h2){ + h2->Sumw2(); + continue; + } + TH3 *h3 = dynamic_cast(fOutput->At(i)); + if (h3){ + h3->Sumw2(); + continue; + } + THnSparse *hn = dynamic_cast(fOutput->At(i)); + if(hn)hn->Sumw2(); + } + + TH1::AddDirectory(oldStatus); + + PostData(1, fOutput); // Post data for ALL output slots > 0 here. +} + +//________________________________________________________________________ +void AliAnalysisTaskTriggerRejection::ExtractMainPatch() { + + //Find main trigger + if(!fTriggerPatchInfo) + return; + + //number of patches in event + Int_t nPatch = fTriggerPatchInfo->GetEntriesFast(); + + //extract main trigger patch + Double_t emax = -1.; + for (Int_t iPatch = 0; iPatch < nPatch; iPatch++) { + AliEmcalTriggerPatchInfo *patch = (AliEmcalTriggerPatchInfo*)fTriggerPatchInfo->At( iPatch ); + if(patch->GetPatchE()>emax) { + fMaxPatch = patch; + emax = patch->GetPatchE(); + } + } +} +//________________________________________________________________________ +Bool_t AliAnalysisTaskTriggerRejection::FillHistograms() +{ + // Fill histograms. + + if(!GetJetContainer(fContainerFull) || !GetJetContainer(fContainerCharged) || !fMaxPatch) + return kFALSE; + + //Get leading jets + AliEmcalJet *leadJet1 = GetJetContainer(fContainerFull)->GetLeadingJet("rho"); + AliEmcalJet *leadJet2 = GetJetContainer(fContainerCharged)->GetLeadingJet("rho"); + + Double_t ptLeadJet1 = -999; + Double_t ptLeadJet2 = -999; + + if(leadJet1) ptLeadJet1 = leadJet1->Pt() - GetRhoVal(fContainerFull)*leadJet1->Area(); + if(leadJet2) ptLeadJet2 = leadJet2->Pt() - GetRhoVal(fContainerCharged)*leadJet2->Area(); + + Double_t VZEROAmp = (Double_t)(InputEvent()->GetVZEROData()->GetMTotV0A() + InputEvent()->GetVZEROData()->GetMTotV0C()); + + //cent;V0mult;ptjet1;ptjet2;Epatch;ADCpatch;EtaPatch;PhiPatch + Double_t var[8] = { + fCent, + VZEROAmp, + ptLeadJet1, + ptLeadJet2, + fMaxPatch->GetPatchE(), + fMaxPatch->GetADCAmp(), + fMaxPatch->GetEtaGeo(), + fMaxPatch->GetPhiGeo() + }; + fhnTriggerInfo->Fill(var); + + return kTRUE; +} + +//________________________________________________________________________ +Bool_t AliAnalysisTaskTriggerRejection::Run() +{ + // Run analysis code here, if needed. It will be executed before FillHistograms(). + + if(fTriggerPatchInfo) + ExtractMainPatch(); + + return kTRUE; // If return kFALSE FillHistogram() will NOT be executed. +} + +//_______________________________________________________________________ +void AliAnalysisTaskTriggerRejection::Terminate(Option_t *) +{ + // Called once at the end of the analysis. +} diff --git a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskTriggerRejection.h b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskTriggerRejection.h new file mode 100644 index 00000000000..b092ad5064a --- /dev/null +++ b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskTriggerRejection.h @@ -0,0 +1,54 @@ +#ifndef ALIANALYSISTASKTRIGGERREJECTION_H +#define ALIANALYSISTASKTRIGGERREJECTION_H + +class TH1; +class TH2; +class TH3; +class TH3F; +class TProfile; +class THnSparse; +class TClonesArray; +class TArrayI; + +#include +#include +#include + +#include + +#include "AliAnalysisTaskEmcalJet.h" +#include "AliEmcalTriggerPatchInfo.h" + +class AliAnalysisTaskTriggerRejection : public AliAnalysisTaskEmcalJet { + public: + + AliAnalysisTaskTriggerRejection(); + AliAnalysisTaskTriggerRejection(const char *name); + virtual ~AliAnalysisTaskTriggerRejection(); + + void UserCreateOutputObjects(); + void Terminate(Option_t *option); + + //Setters + void SetContainerFull(Int_t c) { fContainerFull = c;} + void SetContainerCharged(Int_t c) { fContainerCharged = c;} + + protected: + Bool_t FillHistograms() ; + Bool_t Run() ; + void ExtractMainPatch(); + + private: + Int_t fContainerFull; // number of container with full jets DET + Int_t fContainerCharged; // number of container with charged jets DET + AliEmcalTriggerPatchInfo *fMaxPatch; // main patch + THnSparse *fhnTriggerInfo; //! correlation between jets, patch energy and event observables + + AliAnalysisTaskTriggerRejection(const AliAnalysisTaskTriggerRejection&); // not implemented + AliAnalysisTaskTriggerRejection &operator=(const AliAnalysisTaskTriggerRejection&); // not implemented + + ClassDef(AliAnalysisTaskTriggerRejection, 1) +}; +#endif + + diff --git a/PWGJE/EMCALJetTasks/macros/AddTaskTriggerRejection.C b/PWGJE/EMCALJetTasks/macros/AddTaskTriggerRejection.C new file mode 100644 index 00000000000..db016feb0db --- /dev/null +++ b/PWGJE/EMCALJetTasks/macros/AddTaskTriggerRejection.C @@ -0,0 +1,190 @@ +enum AlgoType {kKT, kANTIKT}; +enum JetType {kFULLJETS, kCHARGEDJETS, kNEUTRALJETS}; + +AliAnalysisTaskTriggerRejection* AddTaskTriggerRejection(TString kTracksName = "PicoTracks", + TString kClusName = "caloClustersCorr", + Double_t R = 0.4, + Double_t ptminTrack = 0.15, + Double_t etminClus = 0.3, + Int_t rhoType = 1, + const char *CentEst = "V0M", + Int_t pSel = AliVEvent::kAny, + TString kEmcalTriggers = "", + TString kPeriod = "LHC11h", + TString kBeamType = "PbPb", //or pPb or pp + TString tag = "" + ) { + // The following three lines are added for backwards compatibility + kPeriod.ToLower(); + if(kPeriod.EqualTo("lhc10h") || kPeriod.EqualTo("lhc11h")) kBeamType = "PbPb"; + if(kPeriod.EqualTo("lhc13b") || kPeriod.EqualTo("lhc13c") || kPeriod.EqualTo("lhc13d") || kPeriod.EqualTo("lhc13e") || kPeriod.EqualTo("lhc13f")) kBeamType = "pPb"; + + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + if (!mgr) + { + Error("AddTaskTriggerRejection","No analysis manager found."); + return 0; + } + Bool_t ismc=kFALSE; + ismc = (mgr->GetMCtruthEventHandler())?kTRUE:kFALSE; + + // Check the analysis type using the event handlers connected to the analysis manager. + //============================================================================== + if (!mgr->GetInputEventHandler()) + { + ::Error("AddTaskTriggerRejection", "This task requires an input event handler"); + return NULL; + } + + // #### Add necessary jet finder tasks + gROOT->LoadMacro("$ALICE_ROOT/PWGJE/EMCALJetTasks/macros/AddTaskEmcalJet.C"); + + AliEmcalJetTask* jetFinderTask1; + AliEmcalJetTask* jetFinderTask2; + if(kClusName.IsNull()) { //particle level jets + jetFinderTask1 = AddTaskEmcalJet(kTracksName, "", kANTIKT, R, kFULLJETS, ptminTrack, etminClus,0.005,1,"Jet",1.); + jetFinderTask2 = AddTaskEmcalJet(kTracksName, "", kANTIKT, R, kCHARGEDJETS, ptminTrack, etminClus,0.005,1,"Jet",1.); + } + else if(kTracksName.IsNull()) { //neutral/calo jets + jetFinderTask1 = AddTaskEmcalJet(kTracksName, kClusName, kANTIKT, R, kNEUTRALJETS, ptminTrack, etminClus,0.005,1,"Jet",1.); + jetFinderTask2 = AddTaskEmcalJet(kTracksName, kClusName, kANTIKT, R, kFULLJETS, ptminTrack, etminClus,0.005,1,"Jet",1.); + } + else { //full jets + jetFinderTask1 = AddTaskEmcalJet(kTracksName, kClusName, kANTIKT, R, kFULLJETS, ptminTrack, etminClus,0.005,1,"Jet",1.); + jetFinderTask2 = AddTaskEmcalJet(kTracksName, kClusName, kANTIKT, R, kCHARGEDJETS, ptminTrack, etminClus,0.005,1,"Jet",1.); + } + jetFinderTask1->SelectCollisionCandidates(AliVEvent::kAny); + jetFinderTask2->SelectCollisionCandidates(AliVEvent::kAny); + + TString strJets1 = jetFinderTask1->GetName(); + TString strJets2 = jetFinderTask2->GetName(); + + AliAnalysisTaskRhoBase *rhoTask; + if(rhoType==1) { + rhoTask = AttachRhoTask(kBeamType,kTracksName,kClusName,R,ptminTrack,etminClus); + if(rhoTask) { + rhoTask->SetCentralityEstimator(CentEst); + rhoTask->SelectCollisionCandidates(AliVEvent::kAny); + } + else { + Warning("AddTaskTriggerRejection","Asked for rho task but configuration unknown. Continuing configuration without rho task."); + rhoType = 0; + } + } + + TString wagonName = Form("TriggerRejectionQA_%s_%s_%s",strJets1.Data(),strJets2.Data(),tag.Data()); + + //Configure TriggerQA task + AliAnalysisTaskTriggerRejection *task = new AliAnalysisTaskTriggerRejection(wagonName); + AliParticleContainer *trackCont = task->AddParticleContainer(kTracksName.Data()); + AliClusterContainer *clusterCont = task->AddClusterContainer(kClusName.Data()); + + task->SetContainerFull(0); + task->SetContainerCharged(1); + AliJetContainer *jetCont0 = task->AddJetContainer(strJets1.Data(),"EMCAL",R); + if(rhoType==1) task->SetRhoName(rhoTask->GetOutRhoScaledName(),0); + AliJetContainer *jetCont1 = NULL; + if(strJets2.Contains("Charged")) { + jetCont1 = task->AddJetContainer(strJets2.Data(),"TPC",R); + if(rhoType==1) task->SetRhoName(rhoTask->GetOutRhoName(),1); + } + else { + jetCont1 = task->AddJetContainer(strJets2.Data(),"EMCAL",R); + task->SetZLeadingCut(0.98,0.98,1); + } + + task->SetZLeadingCut(0.98,0.98,0); + for(Int_t i=0; i<2; i++) + task->SetPercAreaCut(0.6, i); + + task->SetCaloTriggerPatchInfoName(kEmcalTriggers.Data()); + task->SetCentralityEstimator(CentEst); + task->SelectCollisionCandidates(pSel); + + task->SetUseAliAnaUtils(kFALSE); + task->SetVzRange(-10.,10.); + + mgr->AddTask(task); + + //Connnect input + mgr->ConnectInput (task, 0, mgr->GetCommonInputContainer() ); + + //Connect output + AliAnalysisDataContainer *coutput1 = 0x0; + TString containerName1 = Form("%s",wagonName.Data()); + TString outputfile = Form("%s:%s",AliAnalysisManager::GetCommonFileName(),wagonName.Data()); + coutput1 = mgr->CreateContainer(containerName1, TList::Class(),AliAnalysisManager::kOutputContainer,outputfile); + mgr->ConnectOutput(task,1,coutput1); + + return task; +} + +AliAnalysisTaskRhoBase *AttachRhoTask(TString kBeamType = "pp", + TString kTracksName = "PicoTracks", + TString kClusName = "caloClustersCorr", + Double_t R = 0.4, + Double_t ptminTrack = 0.15, + Double_t etminClus = 0.3 + ) { + + AliAnalysisTaskRhoBase *rhoTaskBase; + + // Add kt jet finder and rho task in case we want background subtraction + Double_t minJetPt = 0.1; + if(kBeamType == "pPb") minJetPt = 0.; + AliEmcalJetTask *jetFinderKt; + AliEmcalJetTask *jetFinderAKt; + jetFinderKt = AddTaskEmcalJet(kTracksName, kClusName, kKT, R, kCHARGEDJETS, ptminTrack, etminClus,0.005,1,"Jet",minJetPt); + jetFinderAKt = AddTaskEmcalJet(kTracksName, kClusName, kANTIKT, R, kCHARGEDJETS, ptminTrack, etminClus,0.005,1,"Jet",1.); + jetFinderKt->SelectCollisionCandidates(AliVEvent::kAny); + jetFinderAKt->SelectCollisionCandidates(AliVEvent::kAny); + + if(kBeamType == "pPb") { + gROOT->LoadMacro("$ALICE_ROOT/PWGJE/EMCALJetTasks/macros/AddTaskRhoSparse.C"); + TF1 *fScale = new TF1("fScale","1.28",0.,100.); //scale factor for pPb + AliAnalysisTaskRhoSparse *rhoTaskSparse = AddTaskRhoSparse( + jetFinderKt->GetName(), + jetFinderAKt->GetName(), + kTracksName, + kClusName, + Form("RhoSparseR%03d",(int)(100*R)), + R, + "TPC", + 0.01, + 0.15, + 0, + fScale, + 0, + kTRUE, + Form("RhoSparseR%03d",(int)(100*R)), + kTRUE + ); + rhoTaskSparse->SetUseAliAnaUtils(kTRUE); + rhoTaskBase = dynamic_castrhoTaskSparse; + } + else if(kBeamType=="PbPb") { + gROOT->LoadMacro("$ALICE_ROOT/PWGJE/EMCALJetTasks/macros/AddTaskRho.C"); + + TF1* sfunc = new TF1("sfunc","[0]*x*x+[1]*x+[2]",-1,100); + sfunc->SetParameter(2,1.76458); + sfunc->SetParameter(1,-0.0111656); + sfunc->SetParameter(0,0.000107296); + AliAnalysisTaskRho *rhoTask = AddTaskRho( + jetFinderKt->GetName(), + kTracksName, + kClusName, + Form("RhoR%03d",(int)(100*R)), + R, + "TPC", + 0.01, + 0, + sfunc, + 2, + kTRUE); + rhoTask->SetHistoBins(100,0,250); + + rhoTaskBase = dynamic_castrhoTask; + } + + return rhoTaskBase; +} diff --git a/PWGJE/PWGJEEMCALJetTasksLinkDef.h b/PWGJE/PWGJEEMCALJetTasksLinkDef.h index 76e679a9afd..53800aba6e1 100644 --- a/PWGJE/PWGJEEMCALJetTasksLinkDef.h +++ b/PWGJE/PWGJEEMCALJetTasksLinkDef.h @@ -75,6 +75,7 @@ #pragma link C++ class AliAnalysisTaskSAJF+; #pragma link C++ class AliAnalysisTaskSAQA+; #pragma link C++ class AliAnalysisTaskSOH+; +#pragma link C++ class AliAnalysisTaskTriggerRejection+; #pragma link C++ class AliAnalysisTaskV0sInJetsEmcal+; #pragma link C++ class AliEmcalPicoTrackFromJetMaker+; #pragma link C++ class AliNtupCumInfo+; -- 2.43.0