Commit 4d3dfe99 authored by Maria Mazza's avatar Maria Mazza Committed by Tadej Novak
Browse files

TrigL0GepPerf: add cone jet algorithm

parent fd1469cd
/*
* Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
*/
#ifndef TRIGL0GEPPERF_CONEJETMAKER_H
#define TRIGL0GEPPERF_CONEJETMAKER_H
#include <string>
#include "xAODTrigger/JetRoI.h"
#include "xAODTrigger/JetRoIAuxContainer.h"
#include "xAODTrigger/JetRoIContainer.h"
#include "TrigL0GepPerf/IJetMaker.h"
#include "TrigL0GepPerf/CustomJet.h"
#include "TrigL0GepPerf/CustomTopoCluster.h"
namespace Gep
{
class ConeJetMaker : public Gep::IJetMaker
{
public:
ConeJetMaker(float jetR, float seedEtThreshold = 5.e3 /*MeV*/, std::string recombScheme = "EScheme") :
m_jetR(jetR),
m_seedEtThreshold(seedEtThreshold),
m_recombScheme(recombScheme)
{}
std::string getName() const { return "ConeJet"; }
std::vector<Gep::CustomJet> makeJet(const std::vector<Gep::CustomTopoCluster> &clusters);
void setSeeds(const xAOD::JetRoIContainer *seeds) { m_seeds = seeds; }
float getJetR() { return m_jetR; }
float getSeedEtThreshold() { return m_seedEtThreshold; }
std::string getRecombScheme() { return m_recombScheme; }
private:
float m_jetR;
float m_seedEtThreshold;
std::string m_recombScheme;
const mutable xAOD::JetRoIContainer *m_seeds;
};
}
#endif //TRIGL0GEPPERF_CONEJETMAKER_H
......@@ -12,8 +12,12 @@ namespace Gep{
{
TLorentzVector vec;
float radius;
std::vector<int> constituentsIndices;
int nConstituents;
float radius {0};
float seedEta {0};
float seedPhi {0};
float seedEt {0};
};
}
......
......@@ -3,19 +3,27 @@ import logging as log
from AthenaCommon.Include import include
#
# The algorithms can be set from the joboptions.
# The current available options are the following.
# The default available options are the following.
#
# topoclustering and default clusters:
# topoclAlgs = ['CaloCal','Calo420','Calo422','CaloWFS']
# Topoclusters:
# 'Calo422' : 42 topoclusters at EM scale
# 'Calo420' : 420 topoclusters at EM scale
# 'CaloCal' : 420 topoclusters LC calibrated
#
# pileup suppression:
# puSupprAlgs = ['', 'Vor', 'SK', 'VorSK']
# Pileup suppression:
# 'Vor' : Voronoi,
# 'SK' : SoftKiller
# 'VorSK' : Voronoi+SoftKiller
#
# jet reconstruction:
# jetAlgs = ['AntiKt4']
# Jet reconstruction:
# 'AntiKt4' : offline AntiKt4
#
# Custom topoclustering and jet reconstruction algorithms can be called
# following the naming scheme in GepClusteringAlg.cxx and GepJetAlg.cxx.
#
def setupL0GepSimulationSequence(
topoclAlgs = ['Calo422'],
puSupprAlgs = [''],
......@@ -49,6 +57,20 @@ def setupL0GepSimulationSequence(
from xAODCaloEventCnv.xAODCaloEventCnvConf import ClusterCreator
topSequence+=ClusterCreator("CaloCluster2xAOD")
# run JFEX algorithms
from TrigT1CaloFexPerf.L1PerfControlFlags import L1Phase1PerfFlags as pflags
pflags.Calo.UseAllCalo=True
pflags.CTP.RunCTPEmulation=False
from RecExConfig.RecFlags import rec
rec.readAOD=True
rec.readESD=True
rec.readRDO=False
rec.doESD=True
rec.doWriteAOD=False
include("TrigT1CaloFexPerf/createL1PerfSequence.py")
for topoMaker in topoclAlgs:
# topoclusters
......
/*
* Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
*/
#include "TrigL0GepPerf/CustomJet.h"
#include "TrigL0GepPerf/ConeJetMaker.h"
double deltaR (double eta_1, double eta_2, double phi_1, double phi_2);
std::vector<Gep::CustomJet> Gep::ConeJetMaker::makeJet( const std::vector<Gep::CustomTopoCluster> &clusters)
{
std::vector<Gep::CustomJet> jets;
for (auto seed: *m_seeds) {
float seedEt = seed->et8x8();
// skip seeds with Et below threshold
if(seedEt < m_seedEtThreshold) continue;
float seedEta = seed->eta();
float seedPhi = seed->phi();
Gep::CustomJet jet;
jet.radius = m_jetR;
jet.seedEta = seedEta;
jet.seedPhi = seedPhi;
jet.seedEt = seedEt;
TLorentzVector jetVec;
float px{0}, py{0}, pz{0};
int clusterIndex {0};
//build jet with clusters within dR from seed
for (const auto &cl: clusters) {
float dR_seed_cl = deltaR(seedEta, cl.vec.Eta(), seedPhi, cl.vec.Phi());
if (dR_seed_cl < m_jetR) {
jetVec += cl.vec;
px += cl.vec.Px();
py += cl.vec.Py();
pz += cl.vec.Pz();
jet.constituentsIndices.push_back(clusterIndex);
}
clusterIndex++;
}
// skip cone jets with 0 constituents
if (jet.constituentsIndices.size()==0) continue;
// recombination scheme
if (m_recombScheme == "EScheme") {
// default option: add four-vectors of constituents
jet.vec = jetVec;
}
else if (m_recombScheme == "SeedScheme") {
// massless jet, correct pt, re-using seed (eta,phi)
float m = 0;
float pt = std::sqrt(px*px + py*py);
jet.vec.SetPtEtaPhiM(pt, seedEta, seedPhi, m);
}
else {
std::cerr << "Error: " << m_recombScheme << " is not a valid recombination scheme.\n";
std::exit(1);
}
jets.push_back(jet);
}
return jets;
}
double deltaR (double eta1, double eta2, double phi1, double phi2) {
double deltaPhi = TVector2::Phi_mpi_pi(phi1 - phi2);
double deltaEta = eta1 - eta2;
return std::sqrt( deltaEta*deltaEta + deltaPhi*deltaPhi );
}
......@@ -9,7 +9,7 @@
// *** Include derived jet algorithm classes ***
#include "TrigL0GepPerf/ModAntikTJetMaker.h"
#include "TrigL0GepPerf/ConeJetMaker.h"
#include "xAODEventInfo/EventInfo.h"
#include "xAODCaloEvent/CaloClusterContainer.h"
......@@ -49,7 +49,7 @@ StatusCode GepJetAlg::execute() {
std::vector<Gep::CustomTopoCluster> customClusters;
for( auto iClus: *clusters){
for( auto iClus: *clusters ){
Gep::CustomTopoCluster clus;
clus.vec.SetPxPyPzE(iClus->p4().Px(), iClus->p4().Py(),
iClus->p4().Pz(), iClus->e());
......@@ -67,8 +67,15 @@ StatusCode GepJetAlg::execute() {
jetMaker = std::move(customJetAlg);
}
else if ( m_jetAlg=="Cone4SeedjRJ" ) {
auto coneJetAlg = std::make_unique<Gep::ConeJetMaker>(0.4);
const xAOD::JetRoIContainer * seeds = 0;
CHECK(evtStore()->retrieve(seeds,"jRoundJetsPerf"));
coneJetAlg->setSeeds(seeds);
jetMaker = std::move(coneJetAlg);
}
if( !jetMaker ){
ATH_MSG_ERROR( "JetMaker is a null pointer." );
return StatusCode::FAILURE;
......@@ -77,6 +84,10 @@ StatusCode GepJetAlg::execute() {
std::vector<Gep::CustomJet> customJets = jetMaker->makeJet( customClusters );
// if no jets were found, skip event
if( customJets.empty() ){
return StatusCode::SUCCESS;
}
// create the new container and its auxiliary store
auto athenaJets = std::make_unique<xAOD::JetContainer>();
......@@ -100,6 +111,15 @@ StatusCode GepJetAlg::execute() {
ijet->setJetP4( P4 );
ijet->setAttribute("RCut", iJet.radius);
ijet->setAttribute("SeedEta", iJet.seedEta); // < custom attributes
ijet->setAttribute("SeedPhi", iJet.seedPhi); //
ijet->setAttribute("SeedEt", iJet.seedEt); //
for (const auto& i: iJet.constituentsIndices) {
ijet->addConstituent(clusters->at(i));
}
}
std::string jetsName = m_jetAlg + m_topoclLabel + "Jets";
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment