diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JFexEleTau.h b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JFexEleTau.h
new file mode 100644
index 0000000000000000000000000000000000000000..4d0669babd4485bad20d591105f956f9e9be937d
--- /dev/null
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JFexEleTau.h
@@ -0,0 +1,32 @@
+#ifndef TRIGT1CALOFEXSIM_JFEXELETAU_H
+#define TRIGT1CALOFEXSIM_JFEXELETAU_H 1
+
+#include "AthAnalysisBaseComps/AthAnalysisAlgorithm.h"
+#include "xAODTrigL1Calo/JGTowerContainer.h"
+#include "xAODTrigger/EmTauRoIContainer.h"
+
+class JFexEleTau: public ::AthAnalysisAlgorithm {
+ public:
+ JFexEleTau( const std::string& name, ISvcLocator* pSvcLocator );
+ virtual ~JFexEleTau();
+
+ ///uncomment and implement methods as required
+ virtual StatusCode initialize(); //once, before any input is loaded
+ virtual StatusCode execute(); //per event
+ virtual StatusCode finalize(); //once, after all events processed
+ virtual StatusCode analyzeSeedRun2(float eta, float phi, const xAOD::JGTowerContainer* jTowers, xAOD::EmTauRoIContainer* clustersForTau); // Use a given seed to produce a candidate
+ virtual StatusCode analyzeSeedRun3(float eta, float phi, const xAOD::JGTowerContainer* jTowers, xAOD::EmTauRoIContainer* clustersForTau); // Use a given seed to produce a candidate
+ virtual StatusCode analyzeSeedEle(float eta, float phi, const xAOD::JGTowerContainer* jTowers, xAOD::EmTauRoIContainer* clustersForEle); // Use a given seed to produce a candidate
+
+ private:
+
+ bool m_regenerateSeeds = false;
+ bool m_applyNoise = false;
+ bool m_checkMax = false;
+ int m_noiseStrategy = 0;
+ bool m_useRun2 = false;
+ bool m_singleTowerSeed = false;
+
+};
+
+#endif //> !TRIGT1CALOFEXSIM_JFEXELETAU_H
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JGTower.h b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JGTower.h
index 05128856571aaf9dcbc8ce6a460403528a9e5d50..d2762d8225d4258e8539673f7fb0fd4fccc8db29 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JGTower.h
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JGTower.h
@@ -59,6 +59,10 @@ class JGTower{
bool inBox(float eta1, float eta2, float deta, float phi1, float phi2, float dphi);
bool withinRadius(float eta1, float eta2, float phi1, float phi2, float dR, bool acceptEqual=false);
float deltaPhi(float phi1, float phi2);
+int GFEX_pFPGA_Int(std::string in);
+std::string GFEX_pFPGA(float eta);
+float GTowerArea(float eta);
+
class TowerHelper{
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JGTowerReader.h b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JGTowerReader.h
index a5c1a56ac4cb9042ba599e27aaf1d54ca1a775a2..de85e43d8e6382d52c1367c47dfb0c664c258c26 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JGTowerReader.h
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JGTowerReader.h
@@ -98,15 +98,18 @@ class JGTowerReader: public ::AthAlgorithm {
bool m_plotSeeds;
bool m_saveSeeds;
- float m_gJet_seed_size;
- float m_gJet_max_r;
float m_gJet_r;
- float m_gJet_seed_tower_noise_multiplier;
- float m_gJet_seed_total_noise_multiplier;
- float m_gJet_seed_min_ET_MeV;
- float m_gJet_jet_tower_noise_multiplier;
- float m_gJet_jet_total_noise_multiplier;
+ float m_gJet_block_tower_noise_multiplier;
+ float m_gJet_block_min_ET_MeV;
+ float m_gJet_tower_noise_multiplier;
float m_gJet_jet_min_ET_MeV;
+ float m_gFEX_pTcone_cut;
+ bool m_gFEX_OnlyPosRho;
+ bool m_gFEX_useNegTowers;
+ bool m_gFEX_Rho_useNegTowers;
+ //job options for gFEX MET algorithms
+ bool m_useRMS;
+ bool m_useMedian;
std::string m_noise_file;
@@ -115,16 +118,6 @@ class JGTowerReader: public ::AthAlgorithm {
float m_jXERHO_rho_up_threshold;
float m_jXERHO_min_noise_cut;
- //job options for gFEX MET algorithms
- bool m_useRMS;
- bool m_useMedian;
- bool m_useNegTowers;
- bool m_developerMET;
- bool m_combine_rhoNoise;
- bool m_combine_skNoise;
- bool m_combine_jwojNoise;
- float m_pTcone_cut;
-
const CaloCell_SuperCell_ID* m_scid;
const JTower_ID* m_jTowerId;
const GTower_ID* m_gTowerId;
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JetAlg.h b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JetAlg.h
index 6ac8f608eb90c53796386de9cc2ff4c78a882851..c15c1d5a37201441954bd07c9a4d41e324083e77 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JetAlg.h
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JetAlg.h
@@ -58,11 +58,11 @@ class JetAlg{
static std::map<TString, std::shared_ptr<Seed>> m_SeedMap;
static StatusCode SeedGrid(const xAOD::JGTowerContainer*towers, TString seedname, bool &m_dumpSeedsEtaPhi);
- static StatusCode SeedFinding(const xAOD::JGTowerContainer*towers, TString seedname, float seed_size, float range, std::vector<float> noise, float seed_tower_noise_multiplier, float seed_total_noise_multiplier, float seed_min_ET_MeV);
+ static StatusCode SeedFinding(const xAOD::JGTowerContainer*towers, TString seedname, float seed_size, float range, std::vector<float> noise, float seed_tower_noise_multiplier, float seed_min_ET_MeV,bool seed_electrons);
- static StatusCode BuildFatJet(const xAOD::JGTowerContainer towers, TString jetname, float jet_r, std::vector<float> noise, float jet_tower_noise_multiplier, float jet_total_noise_multiplier, float jet_min_ET_MeV, float rho);
- static StatusCode BuildJet(const xAOD::JGTowerContainer*towers, TString seedname, TString jetname, float jet_size, std::vector<float> noise, float jet_tower_noise_multiplier, float jet_total_noise_multiplier, float jet_min_ET_MeV, bool m_saveSeeds = false);
- static StatusCode BuildRoundJet(const xAOD::JGTowerContainer*towers, TString seedname, TString jetname, float jet_size, std::vector<float> noise, float jet_tower_noise_multiplier, float jet_total_noise_multiplier, float jet_min_ET_MeV, bool m_saveSeeds = false);
- static StatusCode BuildgBlocksJets(const xAOD::JGTowerContainer* blocks, TString jetname, float rho);
+ static StatusCode BuildFatJet(const xAOD::JGTowerContainer towers, TString jetname, float jet_r, std::vector<float> noise, float jet_tower_noise_multiplier, float seed_threshold, float jet_min_ET_MeV, float rhoA, float rhoB,float rhoC, bool useNegTowers);
+ static StatusCode BuildJet(const xAOD::JGTowerContainer*towers, TString seedname, TString jetname, float jet_size, std::vector<float> noise, float jet_tower_noise_multiplier, float jet_min_ET_MeV, bool m_saveSeeds = false);
+ static StatusCode BuildRoundJet(const xAOD::JGTowerContainer*towers, TString seedname, TString jetname, float jet_size, std::vector<float> noise, float jet_tower_noise_multiplier, float jet_min_ET_MeV, bool m_saveSeeds = false);
+ static StatusCode BuildgBlocksJets(const xAOD::JGTowerContainer* blocks, TString jetname, float rhoA, float rhoB, float rhoC);
};
#endif
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JwoJ.h b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JwoJ.h
index ec8f52ef65dfb734831986cb2d490942ad229b70..480700f2c41dc0b426f70b2d40441bd58741fa8d 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JwoJ.h
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/JwoJ.h
@@ -20,88 +20,59 @@
#include "Objects.h"
#include "GaudiKernel/SystemOfUnits.h"
#include "JGTowerReader.h"
-
+#include "AthenaBaseComps/AthCheckMacros.h"
/**
*@brief Separates Et into hard and soft terms depending on threshold pt. Stores hard and soft MET terms in a vector to be passed.
*@return @c std::vector<float>
*/
-std::vector<float> Run_JwoJ(const xAOD::JGTowerContainer* towers, const std::vector<TowerObject::Block> blocks, float rho, float pTcone_cut, bool useEtaBins, bool useNegTowers){
+std::vector<float> Run_JwoJ(const xAOD::JGTowerContainer* towers, const std::vector<TowerObject::Block> blocks, float pTcone_cut, bool useRho, float RhoA, float RhoB, float RhoC, bool useNegTowers){
/*
By this definition, if we set useNegTowers to true, then we will use all the towers. If we set it to false, then we veto them after pileup subtraction.
*/
+
pTcone_cut*=Gaudi::Units::GeV;
std::vector<float> Et_vals;
- float Ht = 0; float Htx = 0; float Hty = 0;
- float Et = 0; float Etx = 0; float Ety = 0;
+ float Htx = 0; float Hty = 0;
+ float Etx = 0; float Ety = 0;
float Et_tot = 0; float Etx_tot = 0; float Ety_tot = 0;
float scalar_Et = 0;
for(unsigned int b = 0; b < blocks.size(); b++){
- float block_phi = blocks[b].Phi();
- float pt_cone = blocks[b].Pt();
- float seed_Et = (towers->at(blocks[b].seedIndex()))->et();
-
- float block_etx = 0;
- float block_ety = 0;
-
- if(useEtaBins){
- seed_Et -= rho;
- if(!useNegTowers){
- //If we set useNegTowers = False, then we would set all their ET = 0; This is set after pileup subtraction, which is a choice that could be studied, but seems reasonable.
- if(seed_Et<0)seed_Et=0;
- }
- block_etx = seed_Et*TMath::Cos(block_phi);
- block_ety = seed_Et*TMath::Sin(block_phi);
- if(pt_cone > pTcone_cut){
- Ht += seed_Et;
- Htx += block_etx;
- Hty += block_ety;
- }else{
- Et += seed_Et;
- Etx += block_etx;
- Ety += block_ety;
- }
- }else{
- if(TMath::Abs(blocks[b].Eta()) < 2.4){
- seed_Et -= rho;
- if(!useNegTowers){
- //If we set useNegTowers = False, then we would set all their ET = 0;
- if(seed_Et<0)seed_Et=0;
- }
- block_etx = seed_Et*TMath::Cos(block_phi);
- block_ety = seed_Et*TMath::Sin(block_phi);
+ float pt_block = blocks[b].Pt();//pT of the gBlock
+ float seed_Et = (towers->at(blocks[b].seedIndex()))->et();
+ float seed_phi = (towers->at(blocks[b].seedIndex()))->phi();
+ float seed_eta = (towers->at(blocks[b].seedIndex()))->eta();
+ float Area = GTowerArea(seed_eta);
+ std::string FPGA = GFEX_pFPGA(seed_eta);
+ float rho = 0;
+ if(FPGA=="A") rho = RhoA;
+ else if(FPGA=="B") rho = RhoB;
+ else if(FPGA=="C") rho = RhoC;
+ else ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::WARNING,"JwoJ") << "NO FPGA FOUND";
- if(pt_cone > pTcone_cut){
- Ht += seed_Et;
- Htx += block_etx;
- Hty += block_ety;
- }else{
- Et += seed_Et;
- Etx += block_etx;
- Ety += block_ety;
- }
- }else{
- seed_Et -= 4*rho;
- block_etx = seed_Et*TMath::Cos(block_phi);
- block_ety = seed_Et*TMath::Sin(block_phi);
+ if(!useNegTowers && seed_Et < 0)continue;
+
+ //First compute global variables
+ scalar_Et+=seed_Et;
+ Etx_tot += seed_Et*TMath::Cos(seed_phi);
+ Ety_tot += seed_Et*TMath::Sin(seed_phi);
- if(pt_cone > 4*pTcone_cut){
- Ht += seed_Et;
- Htx += block_etx;
- Hty += block_ety;
- }
- else{
- Et += seed_Et;
- Etx += block_etx;
- Ety += block_ety;
- }
+ if(pt_block > pTcone_cut){
+ if(useRho) {
+ //Only apply rho correction for hard term
+ seed_Et -= rho*Area;
+ if(seed_Et < 0 ) seed_Et = 0;
}
+ Htx += seed_Et*TMath::Cos(seed_phi);;
+ Hty += seed_Et*TMath::Sin(seed_phi);;
}
- scalar_Et += seed_Et;
- Etx_tot += block_etx;
- Ety_tot += block_ety;
+ else{
+ Etx += seed_Et*TMath::Cos(seed_phi);;
+ Ety += seed_Et*TMath::Sin(seed_phi);;
+ }
+
}
float MHT = TMath::Sqrt(Htx*Htx + Hty*Hty);
float MET = TMath::Sqrt(Etx*Etx + Ety*Ety);
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/METAlg.h b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/METAlg.h
index 33567280fc2d663d4748e365f8d49e5f6fc68997..69f3822a7fac4da0050abd4a486dc78ba2d44805 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/METAlg.h
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/TrigT1CaloFexSim/METAlg.h
@@ -58,11 +58,11 @@ class METAlg{
/**
*@brief Calculate MET using a fixed 4 sigma noise cut
*/
- static StatusCode Baseline_MET(const xAOD::JGTowerContainer*towers, TString metname, std::vector<float> noise, bool useNegTowers);
+ static StatusCode NoiseCut_MET(const xAOD::JGTowerContainer*towers, TString metname, std::vector<float> noise, bool useNegTowers);
/**
*@brief Calculates MET with pileup subtraction
*/
- static StatusCode SubtractRho_MET(const xAOD::JGTowerContainer* towers, TString metname, bool useEtaBins, bool useRMS, bool useMedian, bool useNegTowers);
+ static StatusCode SubtractRho_MET(const xAOD::JGTowerContainer* towers, TString metname, bool useEtaBins, bool useRMS,bool useNegTowers);
/**
*@brief Calculates MET with pileup subtraction in jFEX
*/
@@ -76,14 +76,14 @@ class METAlg{
/**
*@brief Calculates MET with Jets without Jets
*/
- static StatusCode JwoJ_MET(const xAOD::JGTowerContainer* towers, const std::vector<TowerObject::Block> gBlocks, TString metname, float pTcone_cut, bool useEtaBins, bool useRho, bool useNegTowers);
+ static StatusCode JwoJ_MET(const xAOD::JGTowerContainer* towers, const std::vector<TowerObject::Block> gBlocks, TString metname, float pTcone_cut, bool useRho, float RhoA, float RhoB, float RhoC, bool useNegTowers);
/**
*@brief Calculates MET using PUfit
*/
static StatusCode Pufit_MET(const xAOD::JGTowerContainer* towers, TString metname, bool useNegTowers);
- static float Rho_avg(const xAOD::JGTowerContainer* towers, bool useNegTowers);
- static float Rho_avg_etaRings(const xAOD::JGTowerContainer* towers, int fpga, bool useNegTowers);
+ static float Rho_avg_barrel(const xAOD::JGTowerContainer* towers, bool useNegTowers);
+ static float Rho_avg_etaRings(const xAOD::JGTowerContainer* towers, bool useNegTowers);
};
#endif
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/python/L1SimulationControlFlags.py b/Trigger/TrigT1/TrigT1CaloFexSim/python/L1SimulationControlFlags.py
index 54d0ed094f6d934586939c8775478265e941581b..3aae303adddffea92e50b88d91fb0ef451f33e02 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/python/L1SimulationControlFlags.py
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/python/L1SimulationControlFlags.py
@@ -26,11 +26,14 @@ class SCellType(JobProperty):
"""
statusOn = True
allowedType = ['str']
- StoredValue = "Pulse"
+ StoredValue = "BCID"
+ #We may need this again at some point, so I leave it -- but we are now at the point where we need quality cuts during emulation, B Carlson March 27, 2020
+ '''
def _do_action(self):
if self.get_Value()=="Emulated":
from TrigT1CaloFexSim.L1SimulationControlFlags import L1Phase1SimFlags as simflags
simflags.Calo.ApplySCQual=False
+ '''
_caloflags.append(SCellType)
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/python/L1SimulationSequence.py b/Trigger/TrigT1/TrigT1CaloFexSim/python/L1SimulationSequence.py
index 40653968bf74ba5c1371f3cac8eafd6c2185380f..d8863749a9b38d6f52678f5ec7379124253a0794 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/python/L1SimulationSequence.py
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/python/L1SimulationSequence.py
@@ -164,6 +164,10 @@ def setupRun3L1CaloSimulationSequence(skipCTPEmulation = False, useAlgSequence =
SCBitMask = simflags.Calo.QualBitMask() )
# j/gFEX
l1simAlgSeq += createJGTowerReader(SuperCellType=SCIn) # too much debug output
+
+ #jFEX taus
+ from TrigT1CaloFexSim.TrigT1CaloFexSimConf import JFexEleTau
+ l1simAlgSeq += JFexEleTau(RegenerateSeeds = True, NoiseStrategy=0, ApplyNoise = False, CheckMax = True, UseRun2=False, SingleTowerSeed=True)
#include L1Topo Simulation
if simflags.Topo.RunTopoAlgorithms():
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/python/TrigT1CaloFexSimConfig.py b/Trigger/TrigT1/TrigT1CaloFexSim/python/TrigT1CaloFexSimConfig.py
index b2c43d0e0caf32447b11a47f1be55c751a1a4ec9..0664c6e34f8c5849831414fcaafe84be5f9534a3 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/python/TrigT1CaloFexSimConfig.py
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/python/TrigT1CaloFexSimConfig.py
@@ -30,6 +30,7 @@ def createJGTowerReader( SuperCellType = "SCell", **kwargs ) :
SuperCellType = SuperCellType,
noise_file = "Run3L1CaloSimulation/Noise/noise_r10684_v3.root",
plotSeeds = False,
+ saveSeeds = True,
dumpTowerInfo = False,
makeSquareJets = False,
@@ -75,12 +76,12 @@ def createJGTowerReader( SuperCellType = "SCell", **kwargs ) :
map_jet_total_noise_multiplier = 0.0,
map_jet_min_ET_MeV = 500,
- gJet_seed_size=0.2,#Not for gFEX
- gJet_max_r=1.0,
- gJet_r=1.0,
- gJet_seed_tower_noise_multiplier = 0.0,
- gJet_seed_min_ET_MeV = 500,
+ gJet_r=0.9,
+ gJet_block_tower_noise_multiplier = 0.,
+ gJet_block_min_ET_MeV = 0.,
+ gJet_tower_noise_multiplier = 0.,
gJet_jet_min_ET_MeV = 5000.,
+
)
for arg in kwargs:
towerReader.__setattr__(arg,kwargs[arg])
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/src/JFexEleTau.cxx b/Trigger/TrigT1/TrigT1CaloFexSim/src/JFexEleTau.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..9e6a562ee772bcbaeb2b9d81331a4a070e5707d8
--- /dev/null
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/src/JFexEleTau.cxx
@@ -0,0 +1,494 @@
+// TrigT1CaloFexSim includes
+#include "TrigT1CaloFexSim/JFexEleTau.h"
+#include "xAODTrigL1Calo/JGTowerContainer.h"
+#include "TrigT1CaloFexSim/JetAlg.h"
+//#include "xAODEventInfo/EventInfo.h"
+#include "xAODTrigger/EmTauRoIContainer.h"
+#include "xAODTrigger/EmTauRoIAuxContainer.h"
+#include "xAODTrigger/EmTauRoI.h"
+#include "xAODTrigger/JetRoIContainer.h"
+
+
+
+JFexEleTau::JFexEleTau( const std::string& name, ISvcLocator* pSvcLocator ) : AthAnalysisAlgorithm( name, pSvcLocator ){
+
+
+ declareProperty( "RegenerateSeeds", m_regenerateSeeds = false, "Generate new seeds or use roundJet seeds: note, 'saveSeeds' must be set for the latter" );
+ declareProperty( "ApplyNoise", m_applyNoise = false, "Apply tower noise cut when clustering tau candidate" );
+ declareProperty( "CheckMax", m_checkMax = false, "Require seed to be a local max" );
+ declareProperty( "NoiseStrategy", m_noiseStrategy = 0, "Noise cuts to use (both seeding and clustering, latter only when applyNoise is true): 0 is none, 1 is default, 2 is lowered");
+ declareProperty( "UseRun2", m_useRun2 = false, "Run-II style clustering or Run-III" );
+ declareProperty( "SingleTowerSeed", m_singleTowerSeed = false, "Use 1x1 seed" );
+}
+
+
+JFexEleTau::~JFexEleTau() {}
+
+
+StatusCode JFexEleTau::initialize() {
+ ATH_MSG_INFO ("Initializing " << name() << "...");
+
+ return StatusCode::SUCCESS;
+}
+
+StatusCode JFexEleTau::finalize() {
+ ATH_MSG_INFO ("Finalizing " << name() << "...");
+
+ return StatusCode::SUCCESS;
+}
+
+StatusCode JFexEleTau::execute() {
+ ATH_MSG_DEBUG ("Executing " << name() << "...");
+ setFilterPassed(false); //optional: start with algorithm not passed
+
+ // Prepare output containers (with all xAOD annoying features)
+ xAOD::EmTauRoIContainer* clustersForTau = new xAOD::EmTauRoIContainer();
+ xAOD::EmTauRoIAuxContainer* auxclustersForTau = new xAOD::EmTauRoIAuxContainer();
+ clustersForTau->setStore(auxclustersForTau);
+ clustersForTau->reserve(100);
+ std::string clusterName("jTaus");
+ ATH_CHECK( evtStore()->record(clustersForTau,clusterName) );
+ ATH_CHECK( evtStore()->record(auxclustersForTau,clusterName+"Aux.") );
+
+ // Prepare output containers (with all xAOD annoying features)
+ xAOD::EmTauRoIContainer* clustersForEle = new xAOD::EmTauRoIContainer();
+ xAOD::EmTauRoIAuxContainer* auxclustersForEle = new xAOD::EmTauRoIAuxContainer();
+ clustersForEle->setStore(auxclustersForEle);
+ clustersForEle->reserve(100);
+ std::string clusterNameE("jEles");
+ ATH_CHECK(evtStore()->record(clustersForEle,clusterNameE) );
+ ATH_CHECK ( evtStore()->record(auxclustersForEle,clusterNameE+"Aux.") );
+
+ // Retrieve jTower container
+ const xAOD::JGTowerContainer* jTowers =0;
+ CHECK( evtStore()->retrieve( jTowers,"JTower"));
+
+ // Create noise profile
+ std::vector<float> jT_noise;
+
+ // Control noise cuts
+ for( const auto &jT : *jTowers){
+
+ if(m_noiseStrategy == 0)
+ jT_noise.push_back(0);
+ else if(m_noiseStrategy == 1)
+ {
+ if(jT->sampling()==0) jT_noise.push_back(430);
+ else if(jT->sampling()==1) jT_noise.push_back(2400);
+ else jT_noise.push_back(2000);
+ }
+ else if(m_noiseStrategy == 2)
+ {
+ if(jT->sampling()==0) jT_noise.push_back(125);
+ else if(jT->sampling()==1) jT_noise.push_back(250);
+ else jT_noise.push_back(250);
+ }
+ }
+
+ if(m_regenerateSeeds)
+ {
+ bool dumpPlots = false;
+
+ // We need to regenerate the seed map through JetAlg
+ if( JetAlg::m_SeedMap.find("jTauSeeds") == JetAlg::m_SeedMap.end() )
+ CHECK( JetAlg::SeedGrid(jTowers, "jTauSeeds", dumpPlots) );
+
+ //ATH_MSG_DEBUG("JFexAlg: SeedFinding with jJetSeeds; m_jJet_seed_size = " << m_jJet_seed_size << ", m_jJet_max_r = " << m_jJet_max_r);
+ if(m_singleTowerSeed)
+ CHECK( JetAlg::SeedFinding( jTowers, "jTauSeeds", 0.11, 0.26, jT_noise,
+ 2.0,
+ 0.0,/*electron*/false) );
+ else
+ CHECK( JetAlg::SeedFinding( jTowers, "jTauSeeds", 0.31, 0.26, jT_noise,
+ 2.0,
+ 0.0,/*electron*/false) );
+
+ // Check if we can actually access the 2x2 seeds, since we are technically running after the RoundJets
+ ATH_MSG_DEBUG("Looking for seeds: jFEX Taus");
+
+ std::shared_ptr<JetAlg::Seed> seeds = JetAlg::m_SeedMap["jTauSeeds"];
+
+ ATH_MSG_DEBUG("Looping on seeds: jFEX Taus");
+
+ // Loop over all seeds
+ for(unsigned eta_ind=0; eta_ind<seeds->eta.size(); eta_ind++){
+ for(unsigned phi_ind=0; phi_ind<seeds->phi.at(eta_ind).size(); phi_ind++){
+
+ // Extract seed information
+ float eta = seeds->eta.at(eta_ind);
+ float phi = seeds->phi.at(eta_ind).at(phi_ind);
+
+ // Ignore non-local maxima
+ if(m_checkMax)
+ if(!seeds->local_max.at(eta_ind).at(phi_ind)) continue;
+
+ ATH_MSG_DEBUG("jFEX Taus: got a seed at =" << eta << " " << phi);
+
+ // Use this seed to generate a tau candidate
+ if(m_useRun2)
+ CHECK(analyzeSeedRun2(eta, phi, jTowers, clustersForTau));
+ else
+ CHECK(analyzeSeedRun3(eta, phi, jTowers, clustersForTau));
+
+ }
+ }
+ } else {
+
+ // Retrieve seeds from StoreGate
+ const xAOD::JetRoIContainer* jRoundSeeds = 0;
+ CHECK( evtStore()->retrieve( jRoundSeeds,"jRoundSeeds_localMax"));
+
+ for(unsigned int i=0; i<jRoundSeeds->size(); i++)
+ {
+ const xAOD::JetRoI* thisSeed = jRoundSeeds->at(i);
+
+ // Re-do entire logic
+ float eta = thisSeed->eta();
+ float phi = thisSeed->phi();
+
+ if(m_useRun2)
+ CHECK(analyzeSeedRun2(eta, phi, jTowers, clustersForTau));
+ else
+ CHECK(analyzeSeedRun3(eta, phi, jTowers, clustersForTau));
+
+ }
+
+ // Run electrons
+ // We need to regenerate the seed map through JetAlg
+
+ }
+
+ bool dumpPlots = false;
+ if( JetAlg::m_SeedMap.find("jEleSeeds") == JetAlg::m_SeedMap.end() )
+ CHECK( JetAlg::SeedGrid(jTowers, "jEleSeeds", dumpPlots) );
+
+ // Some notes: set the seed_electrons to true: this will seed on EM only
+ // and will stretch out the 0.21 range window by +0.1 in eta > 2.5
+ CHECK( JetAlg::SeedFinding( jTowers, "jEleSeeds", 0.11, 0.21, jT_noise,
+ 2.0,
+ 0.0, true) );
+
+
+ std::shared_ptr<JetAlg::Seed> seeds = JetAlg::m_SeedMap["jEleSeeds"];
+
+
+ for(unsigned eta_ind=0; eta_ind<seeds->eta.size(); eta_ind++){
+ for(unsigned phi_ind=0; phi_ind<seeds->phi.at(eta_ind).size(); phi_ind++){
+
+ // Extract seed information
+ float eta = seeds->eta.at(eta_ind);
+ float phi = seeds->phi.at(eta_ind).at(phi_ind);
+
+
+ // Ignore non-local maxima
+ if(!seeds->local_max.at(eta_ind).at(phi_ind)) continue;
+
+ ATH_MSG_DEBUG("jFEX Eles: got a seed at =" << eta << " " << phi);
+
+
+ // Use this seed to generate an electron candidate
+ CHECK(analyzeSeedEle(eta, phi, jTowers, clustersForEle));
+
+ }
+ }
+
+ setFilterPassed(true); //if got here, assume that means algorithm passed
+ return StatusCode::SUCCESS;
+}
+
+StatusCode JFexEleTau::analyzeSeedRun2(float eta, float phi, const xAOD::JGTowerContainer* jTowers, xAOD::EmTauRoIContainer* clustersForTau) {
+ ATH_MSG_DEBUG ("analyzing a given seed " << name() << "...");
+
+ // Ignore seeds outside of tracker coverage
+ if(std::abs(eta) > 2.47) return StatusCode::SUCCESS;;
+
+ // Find the 4 EM, 4 Tile towers and 16 isolate EM towers
+ std::vector<const xAOD::JGTower*> EMTowers;
+ std::vector<const xAOD::JGTower*> TileTowers;
+ std::vector<const xAOD::JGTower*> EMIsolation;
+
+ // Loop over constituent towers
+ for(unsigned t=0; t<jTowers->size(); t++){
+
+ const xAOD::JGTower* tower = jTowers->at(t);
+
+ //unused for now
+ //bool isTile = (fabs(tower->eta()) <1.5 && tower->sampling()==1);
+ //bool applyNoise = true;
+ //if(isTile) applyNoise=false;
+
+ float dEta = std::abs(eta - tower->eta());
+ float dPhi = std::abs(TVector2::Phi_mpi_pi(phi - tower->phi()));
+
+ if(dEta < 0.095 && dPhi < 0.095)
+ {
+ // Tower found
+ if(tower->sampling()==1)
+ TileTowers.push_back(tower);
+ else
+ EMTowers.push_back(tower);
+ }
+
+ if(dEta < 0.195 && dPhi < 0.195)
+ {
+ if(tower->sampling()!=1)
+ EMIsolation.push_back(tower);
+ }
+ }
+
+ // Produce tau energy sums
+ float tauEnergy = 0;
+ float tileEnergy = 0;
+ float emEnergy = 0;
+ float fullemEnergy = 0;
+ float isoEnergy = 0;
+
+ // Sum of tile energy: 2x2 sum
+ // Also full 2x2 EM energy: purely for isolation calculation
+ for(unsigned int i=0; i<TileTowers.size(); i++)
+ {
+ tileEnergy += TileTowers[i]->et();
+ fullemEnergy += EMTowers[i]->et();
+ }
+
+ for(unsigned int i=0; i<EMIsolation.size(); i++)
+ isoEnergy += EMIsolation[i]->et();
+
+ // Test a potential problematic case (probably solved by mpi-pi)
+ if(TileTowers.size() == 2 || EMTowers.size() == 2)
+ {
+ ATH_MSG_INFO("Missing towers, eta/phi" << eta << "/" << phi);
+ emEnergy = EMTowers[0]->et() + EMTowers[1]->et();
+ }
+ else if(EMTowers.size() ==4)
+ {
+ // Complicated to sort towers, so try all 2 tower sums in EM
+ // And veto the ones where the towers are not aligned in either eta or phi
+ float biggestSum = 0;
+ for(int i=0; i<4; i++)
+ {
+ for(int j=0; j<4; j++)
+ {
+ if(i==j) continue;
+ bool sameLineEta = std::abs(EMTowers[i]->eta() - EMTowers[j]->eta()) < 0.02;
+ bool sameLinePhi = std::abs(TVector2::Phi_mpi_pi(EMTowers[i]->phi() - EMTowers[j]->phi())) < 0.02;
+
+ float currentSum = EMTowers[i]->et() + EMTowers[j]->et();
+ if(currentSum > biggestSum)
+ {
+ if(sameLineEta || sameLinePhi)
+ biggestSum = currentSum;
+ }
+ }
+ }
+ emEnergy = biggestSum;
+ }
+ else{
+ ATH_MSG_INFO("Wrong number of EMTowers: This is very problematic");
+ }
+
+ // Full tau candidate ET: 1x2 in EM + 2x2 in Tile
+ tauEnergy = emEnergy+tileEnergy;
+ // Tau isolation: 4x4 ring, ignore the 2x2 inner part in EM
+ isoEnergy -= fullemEnergy;
+
+ xAOD::EmTauRoI* clForTau = new xAOD::EmTauRoI();
+ clustersForTau->push_back( clForTau );
+ clForTau->setEta( eta );
+ clForTau->setPhi( phi );
+ clForTau->setTauClus( tauEnergy );
+ clForTau->setEmIsol( isoEnergy );
+
+ return StatusCode::SUCCESS;
+}
+
+StatusCode JFexEleTau::analyzeSeedRun3(float eta, float phi, const xAOD::JGTowerContainer* jTowers, xAOD::EmTauRoIContainer* clustersForTau) {
+ ATH_MSG_DEBUG ("analyzing a given seed " << name() << "...");
+
+ // Ignore seeds outside of tracker coverage
+ if(std::abs(eta) > 2.47) return StatusCode::SUCCESS;;
+
+ // Find the 4 EM, 4 Tile towers and 16 isolate EM towers
+ std::vector<const xAOD::JGTower*> EMTowers;
+ std::vector<const xAOD::JGTower*> TileTowers;
+ std::vector<const xAOD::JGTower*> EMIsolation;
+
+ // Loop over constituent towers
+ for(unsigned t=0; t<jTowers->size(); t++){
+
+ const xAOD::JGTower* tower = jTowers->at(t);
+
+ //unused for now
+ //bool isTile = (fabs(tower->eta()) <1.5 && tower->sampling()==1);
+ //bool applyNoise = true;
+ //if(isTile) applyNoise=false;
+
+ float dEta = std::abs(eta - tower->eta());
+ float dPhi = std::abs(TVector2::Phi_mpi_pi(phi - tower->phi()));
+
+ if(sqrt(dEta*dEta + dPhi*dPhi) < 0.15)
+ {
+ // Tower found
+ if(tower->sampling()==1)
+ TileTowers.push_back(tower);
+ else
+ EMTowers.push_back(tower);
+ }
+
+ if(sqrt(dEta*dEta + dPhi*dPhi) < 0.35)
+ {
+ if(tower->sampling()!=1)
+ EMIsolation.push_back(tower);
+ }
+ }
+
+ // Produce tau energy sums
+ float tauEnergy = 0;
+ float tileEnergy = 0;
+ float fullemEnergy = 0;
+ float isoEnergy = 0;
+
+ // Sum of tile energy: 3x3 sum
+ // Also full 3x3 EM energy
+ for(unsigned int i=0; i<TileTowers.size(); i++)
+ {
+ tileEnergy += TileTowers[i]->et();
+ fullemEnergy += EMTowers[i]->et();
+ }
+
+ for(unsigned int i=0; i<EMIsolation.size(); i++)
+ isoEnergy += EMIsolation[i]->et();
+
+
+ // Full tau candidate ET: 1x2 in EM + 2x2 in Tile
+ tauEnergy = fullemEnergy+tileEnergy;
+
+ // Tau isolation: 4x4 ring, ignore the 2x2 inner part in EM
+ isoEnergy -= fullemEnergy;
+
+ xAOD::EmTauRoI* clForTau = new xAOD::EmTauRoI();
+ clustersForTau->push_back( clForTau );
+ clForTau->setEta( eta );
+ clForTau->setPhi( phi );
+ clForTau->setTauClus( tauEnergy );
+ clForTau->setEmIsol( isoEnergy );
+
+ return StatusCode::SUCCESS;
+}
+
+
+StatusCode JFexEleTau::analyzeSeedEle(float eta, float phi, const xAOD::JGTowerContainer* jTowers, xAOD::EmTauRoIContainer* clustersForEle) {
+ ATH_MSG_DEBUG ("analyzing a given seed " << name() << "...");
+
+ // Ignore seeds outside of tracker coverage
+ if(std::abs(eta) < 2.3) return StatusCode::SUCCESS;;
+
+ // Central tower
+ const xAOD::JGTower* EMTower = 0;
+ const xAOD::JGTower* HADTower = 0;
+
+ // Isolation Ring
+ std::vector<const xAOD::JGTower*> EMIsolation;
+
+ // Loop over constituent towers
+ float nearestdREM= 5.0;
+ float nearestdRHAD = 5.0;
+
+
+ // Look for main tower
+ for(unsigned t=0; t<jTowers->size(); t++){
+
+ const xAOD::JGTower* tower = jTowers->at(t);
+
+ float dEta = std::abs(eta - tower->eta());
+ float dPhi = std::abs(TVector2::Phi_mpi_pi(phi - tower->phi()));
+
+ if(tower->sampling()==1)
+ {
+ if(sqrt(dEta*dEta + dPhi*dPhi) < nearestdRHAD)
+ {
+ HADTower = tower;
+ nearestdRHAD = sqrt(dEta*dEta + dPhi*dPhi);
+ }
+ } else
+ {
+ if(sqrt(dEta*dEta + dPhi*dPhi) < nearestdREM)
+ {
+ EMTower = tower;
+ nearestdREM = sqrt(dEta*dEta + dPhi*dPhi);
+ }
+ }
+
+
+ }
+
+ // Fill isolation
+ for(unsigned t=0; t<jTowers->size(); t++){
+ const xAOD::JGTower* tower = jTowers->at(t);
+
+ float dEta = std::abs(eta - tower->eta());
+ float dPhi = std::abs(TVector2::Phi_mpi_pi(phi - tower->phi()));
+
+ // Ignore central tower
+ if(tower == EMTower)
+ continue;
+
+ if(sqrt(dEta*dEta + dPhi*dPhi) < 0.4)
+ {
+ if(tower->sampling()!=1)
+ EMIsolation.push_back(tower);
+ }
+
+ }
+
+ // Nearest neighbor search
+ float closestdR = 5.0;
+
+ // First, find closest tower in dR
+ for(unsigned int i=0; i< EMIsolation.size(); i++)
+ {
+ const xAOD::JGTower* tower = EMIsolation.at(i);
+
+ float dEta = std::abs(eta - tower->eta());
+ float dPhi = std::abs(TVector2::Phi_mpi_pi(phi - tower->phi()));
+
+ if (sqrt(dEta*dEta + dPhi*dPhi) < closestdR)
+ closestdR = sqrt(dEta*dEta + dPhi*dPhi);
+ }
+
+ // Then, look for largest pT tower that is closer than 1.9 * nearest dR
+ // Which basically implies, approximately, don't go "two towers away" from the center
+ float highestPt = 0;
+ float fullIso = 0;
+ for(unsigned int i=0; i< EMIsolation.size(); i++)
+ {
+ const xAOD::JGTower* tower = EMIsolation.at(i);
+
+ fullIso += tower->et();
+
+ float dEta = std::abs(eta - tower->eta());
+ float dPhi = std::abs(TVector2::Phi_mpi_pi(phi - tower->phi()));
+
+ if (sqrt(dEta*dEta + dPhi*dPhi) < (closestdR*1.9) && tower->et() > highestPt)
+ highestPt = tower->et();
+ }
+
+
+ // Produce electron energy sums
+ float eleEnergy = EMTower->et() + highestPt;
+ float hadEnergy = HADTower->et();
+ float isoEnergy = fullIso - highestPt;
+
+ xAOD::EmTauRoI* clForEle = new xAOD::EmTauRoI();
+ clustersForEle->push_back( clForEle );
+ clForEle->setEta( eta );
+ clForEle->setPhi( phi );
+ clForEle->setEmClus( eleEnergy );
+ clForEle->setTauClus( hadEnergy );
+ clForEle->setEmIsol( isoEnergy );
+
+ return StatusCode::SUCCESS;
+}
+
+
+
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/src/JGTower.cxx b/Trigger/TrigT1/TrigT1CaloFexSim/src/JGTower.cxx
index fcfa0b280e0e05581c091432b456689248e2b757..63c3855c1443ad60a112f6d4fdef5caf96018960 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/src/JGTower.cxx
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/src/JGTower.cxx
@@ -42,11 +42,11 @@ bool withinRadius(float eta1, float eta2, float phi1, float phi2, float dR, bool
//https://root.cern.ch/doc/master/TLorentzVector_8h_source.html#l00456
float dPhi = deltaPhi(phi1, phi2);
float dEta = eta1-eta2;
- float dR2 = dPhi*dPhi + dEta*dEta;
+ float test_dR = TMath::Sqrt(dPhi*dPhi + dEta*dEta);
if(acceptEqual)
- return dR2 <= dR*dR;
+ return test_dR <= dR;
else
- return dR2 < dR*dR;
+ return test_dR < dR;
}
float deltaPhi(float phi1,float phi2){
@@ -54,6 +54,39 @@ float deltaPhi(float phi1,float phi2){
return TVector2::Phi_mpi_pi(phi1-phi2);
}
+int GFEX_pFPGA_Int(std::string in){
+ if(in=="A") return 1;
+ else if(in=="B") return 2;
+ else if(in=="C") return 3;
+
+ return -1;
+}
+
+std::string GFEX_pFPGA(float eta){
+ if(eta <=0 &&eta > -2.5) return "A";
+ else if(eta > 0 &&eta <=2.5) return "B";
+ else if(fabs(eta)>2.5) return "C";
+
+ return "NOTANFPGA";
+}
+
+float GTowerArea(float eta){
+ /*
+ Returns the gTower area in units of 0.2x0.2 towers
+ There are irregular bins from 2.4-2.5 and 3.1-3.2, but the bins from
+ 3.1-3.2 has A = 1
+ */
+ float abs_eta = fabs(eta);
+ if(abs_eta < 2.4) return 1;
+ else if(abs_eta >= 2.4 && abs_eta < 2.5) return 0.5;
+ else if(abs_eta >=2.5 && abs_eta < 3.1) return 1.;
+ else if(abs_eta >= 3.1 && abs_eta < 3.2) return 1.;
+ else if(abs_eta >= 3.2) return 4.;
+
+ return -1;
+
+}
+
//The helpers
TowerHelper::TowerHelper(std::vector<TH2F*> & h_inputs){
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/src/JGTowerReader.cxx b/Trigger/TrigT1/TrigT1CaloFexSim/src/JGTowerReader.cxx
index 0cbc928d3c4d0228d87f56c7071b0661f865425a..0e751066f5502de4c7ca9c222e73fd6554bc2ffe 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/src/JGTowerReader.cxx
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/src/JGTowerReader.cxx
@@ -79,7 +79,8 @@ histSvc("THistSvc",name){
declareProperty("jJetRound_LargeR_max_r", m_jJetRound_LargeR_max_r=0.4);
declareProperty("jJetRound_LargeR_r", m_jJetRound_LargeR_r=0.8);
declareProperty("jJetRound_LargeR_seed_min_ET_MeV", m_jJetRound_LargeR_seed_min_ET_MeV = 2000.0);
- declareProperty("jJetRound_LargeR_jet_min_ET_MeV", m_jJetRound_jet_min_ET_MeV = 2000.0);
+ declareProperty("jJetRound_LargeR_jet_min_ET_MeV", m_jJetRound_LargeR_jet_min_ET_MeV = 2000.0);
+
declareProperty("makeJetsFromMap", m_makeJetsFromMap = false);
declareProperty("towerMap", m_towerMap = "");
@@ -94,21 +95,19 @@ histSvc("THistSvc",name){
declareProperty("saveSeeds", m_saveSeeds = false);
declareProperty("buildgBlockJets", m_buildgBlockJets=true);
- declareProperty("gJet_seed_size", m_gJet_seed_size=0.2);
- declareProperty("gJet_max_r", m_gJet_max_r=1.0); //gFEX constructs large radius jets
+
declareProperty("gJet_r", m_gJet_r=1.0);
- declareProperty("gJet_seed_tower_noise_multiplier", m_gJet_seed_tower_noise_multiplier = 1.0);
- declareProperty("gJet_seed_total_noise_multiplier", m_gJet_seed_total_noise_multiplier = 1.0);
- declareProperty("gJet_seed_min_ET_MeV", m_gJet_seed_min_ET_MeV = 2000.0);
- declareProperty("gJet_jet_tower_noise_multiplier", m_gJet_jet_tower_noise_multiplier = 1.0);
- declareProperty("gJet_jet_total_noise_multiplier", m_gJet_jet_total_noise_multiplier = 0.0);
+ declareProperty("gJet_block_tower_noise_multiplier", m_gJet_block_tower_noise_multiplier = 0.0);
+ declareProperty("gJet_block_min_ET_MeV", m_gJet_block_min_ET_MeV = 0.0);
+ declareProperty("gJet_tower_noise_multiplier", m_gJet_tower_noise_multiplier = 0.0);
declareProperty("gJet_jet_min_ET_MeV", m_gJet_jet_min_ET_MeV = 2000.0);
- declareProperty("developerMET",m_developerMET=false);
declareProperty("useRMS", m_useRMS=true);
declareProperty("useMedian", m_useMedian=false);
- declareProperty("useNegTowers", m_useNegTowers=false);
- declareProperty("pTcone_cut", m_pTcone_cut=25); //cone threshold for Jets without Jets: declared in GeV
+ declareProperty("gFEX_useNegTowers", m_gFEX_useNegTowers=true);
+ declareProperty("gFEX_Rho_useNegTowers",m_gFEX_Rho_useNegTowers=true);
+ declareProperty("gFEX_OnlyPosRho", m_gFEX_OnlyPosRho=false);
+ declareProperty("gFEX_pTcone_cut", m_gFEX_pTcone_cut=25); //cone threshold for Jets without Jets: declared in GeV
declareProperty("jXERHO_correction_file" , m_jXERHO_correction_file="Run3L1CaloSimulation/Noise/jTowerCorrection.20200302.r11364.root"); //correction file for jXERHO
declareProperty("jXERHO_fixed_noise_cut" , m_jXERHO_fixed_noise_cut=0.0);
@@ -388,7 +387,7 @@ StatusCode JGTowerReader::BuildBlocksFromTowers(std::vector<TowerObject::Block>&
std::vector<int> neighbors = grid.neighbors(*seed, blockRows, blockCols);
float seed_Et = seed->et();
- if(!useNegTowers) seed_Et = TMath::Abs(seed_Et);
+ if(!useNegTowers && seed->et() < 0)continue;
double block_area(0.0);
double sum_deta(0.0);
double sum_dphi(0.0);
@@ -401,7 +400,7 @@ StatusCode JGTowerReader::BuildBlocksFromTowers(std::vector<TowerObject::Block>&
sum_deta += neighbor->deta();
sum_dphi += neighbor->dphi();
neighbor_pt = neighbor->et();
- if(!useNegTowers) neighbor_pt = TMath::Abs(neighbor_pt);
+ if(!useNegTowers && neighbor->et() < 0)continue;
block_pt += neighbor_pt;
}
@@ -428,7 +427,7 @@ StatusCode JGTowerReader::JFexAlg(const xAOD::JGTowerContainer* jTs){
ATH_MSG_DEBUG("Found " << jTs->size() << " jTowers");
ATH_MSG_DEBUG("JFexAlg: BuildMET");
- CHECK(METAlg::Baseline_MET(jTs, "jNOISECUT", jT_noise, m_useNegTowers));
+ CHECK(METAlg::NoiseCut_MET(jTs, "jNOISECUT", jT_noise, /*m_useNegTowers*/ false));
if (!buildbins) {
CHECK(METAlg::build_jFEX_bins(jFEX_bins, jFEX_bins_core, jTs ));
buildbins=true;
@@ -448,15 +447,12 @@ StatusCode JGTowerReader::JFexAlg(const xAOD::JGTowerContainer* jTs){
CHECK( JetAlg::SeedGrid(jTs_PUsub, "jRoundSeedsPUsub", m_dumpSeedsEtaPhi) );
ATH_MSG_DEBUG("JFexAlg using JetsPUsub: SeedFinding with jJetSeeds; m_jJet_seed_size = " << m_jJet_seed_size << ", m_jJet_max_r = " << m_jJet_max_r);
- CHECK( JetAlg::SeedFinding( jTs_PUsub, "jRoundSeedsPUsub", m_jJetRound_seed_size, m_jJetRound_max_r, jT_noise,
- m_jJetRound_seed_tower_noise_multiplier, m_jJetRound_seed_total_noise_multiplier,
- m_jJetRound_seed_min_ET_MeV) );
+ CHECK( JetAlg::SeedFinding( jTs_PUsub, "jRoundSeedsPUsub", m_jJetRound_seed_size, m_jJetRound_max_r, jT_noise,
+ /*m_jJetRound_seed_tower_noise_multiplier*/ 0, m_jJetRound_seed_min_ET_MeV, /*electron*/false) );
ATH_MSG_DEBUG("JFexAlg usinf JetsPUsub: BuildRoundJet");
- CHECK( JetAlg::BuildRoundJet(jTs_PUsub, "jRoundSeedsPUsub", "jRoundJetsPUsub", m_jJetRound_r, jT_noise, m_jJetRound_jet_tower_noise_multiplier, m_jJetRound_jet_total_noise_multiplier, m_jJetRound_jet_min_ET_MeV, m_saveSeeds) );
+ CHECK( JetAlg::BuildRoundJet(jTs_PUsub, "jRoundSeedsPUsub", "jRoundJetsPUsub", m_jJetRound_r, jT_noise, /*m_jJetRound_jet_tower_noise_multiplier*/ 0, m_jJetRound_jet_min_ET_MeV, m_saveSeeds) );
}
- delete jTs_PUsub;
- delete jTs_PUsubAux;
}
if(m_makeSquareJets) {
@@ -471,11 +467,11 @@ StatusCode JGTowerReader::JFexAlg(const xAOD::JGTowerContainer* jTs){
ATH_MSG_DEBUG( "JFexAlg: SeedFinding with jSeeds; m_jJet_seed_size = "
<< m_jJet_seed_size << ", m_jJet_max_r = " << m_jJet_max_r);
CHECK( JetAlg::SeedFinding( jTs, "jSeeds", m_jJet_seed_size, m_jJet_max_r, jT_noise,
- m_jJet_seed_tower_noise_multiplier, m_jJet_seed_total_noise_multiplier,
- m_jJet_seed_min_ET_MeV) );
+ m_jJet_seed_tower_noise_multiplier,
+ m_jJet_seed_min_ET_MeV,/*electron*/false) );
ATH_MSG_DEBUG("JFexAlg: BuildJet");
- CHECK( JetAlg::BuildJet(jTs, "jSeeds", "jJets", m_jJet_r, jT_noise, m_jJet_jet_tower_noise_multiplier, m_jJet_jet_total_noise_multiplier, m_jJet_jet_min_ET_MeV, m_saveSeeds) );
+ CHECK( JetAlg::BuildJet(jTs, "jSeeds", "jJets", m_jJet_r, jT_noise, m_jJet_jet_tower_noise_multiplier, m_jJet_jet_min_ET_MeV, m_saveSeeds) );
}
if(m_makeRoundJets) {
if( JetAlg::m_SeedMap.find("jRoundSeeds") == JetAlg::m_SeedMap.end() )
@@ -483,22 +479,17 @@ StatusCode JGTowerReader::JFexAlg(const xAOD::JGTowerContainer* jTs){
ATH_MSG_DEBUG("JFexAlg: SeedFinding with jJetSeeds; m_jJet_seed_size = " << m_jJet_seed_size << ", m_jJet_max_r = " << m_jJet_max_r);
CHECK( JetAlg::SeedFinding( jTs, "jRoundSeeds", m_jJetRound_seed_size, m_jJetRound_max_r, jT_noise,
- m_jJetRound_seed_tower_noise_multiplier, m_jJetRound_seed_total_noise_multiplier,
- m_jJetRound_seed_min_ET_MeV) );
+ m_jJetRound_seed_tower_noise_multiplier,
+ m_jJetRound_seed_min_ET_MeV,/*electron*/false) );
ATH_MSG_DEBUG("JFexAlg: BuildRoundJet");
- CHECK( JetAlg::BuildRoundJet(jTs, "jRoundSeeds", "jRoundJets", m_jJetRound_r, jT_noise, m_jJetRound_jet_tower_noise_multiplier, m_jJetRound_jet_total_noise_multiplier, m_jJetRound_jet_min_ET_MeV, m_saveSeeds) );
+ CHECK( JetAlg::BuildRoundJet(jTs, "jRoundSeeds", "jRoundJets", m_jJetRound_r, jT_noise, m_jJetRound_jet_tower_noise_multiplier, m_jJetRound_jet_min_ET_MeV, m_saveSeeds) );
}
if(m_makeRoundLargeRJets) {
- if( JetAlg::m_SeedMap.find("jRoundLargeRSeeds") == JetAlg::m_SeedMap.end() )
- CHECK( JetAlg::SeedGrid(jTs, "jRoundLargeRSeeds", m_dumpSeedsEtaPhi) );
-
- ATH_MSG_DEBUG("JFexAlg: SeedFinding with jJetLargeRSeeds; m_jJet_LargeR_seed_size = " << m_jJetRound_LargeR_seed_size << ", m_jJetRound_LargeR_max_r = " << m_jJetRound_LargeR_max_r);
- CHECK( JetAlg::SeedFinding(jTs, "jRoundLargeRSeeds", m_jJetRound_LargeR_seed_size, m_jJetRound_LargeR_max_r, jT_noise, m_jJetRound_seed_tower_noise_multiplier, m_jJetRound_seed_total_noise_multiplier, m_jJetRound_LargeR_seed_min_ET_MeV) );
-
ATH_MSG_DEBUG("JFexAlg: BuildRoundLargeRJet");
- CHECK( JetAlg::BuildRoundJet(jTs, "jRoundLargeRSeeds", "jRoundLargeRJets", m_jJetRound_LargeR_r, jT_noise, m_jJetRound_jet_tower_noise_multiplier, m_jJetRound_jet_total_noise_multiplier, m_jJetRound_LargeR_jet_min_ET_MeV, m_saveSeeds) );
+ CHECK( JetAlg::BuildRoundJet(jTs, "jRoundSeeds", "jRoundLargeRJets", m_jJetRound_LargeR_r, jT_noise, m_jJetRound_jet_tower_noise_multiplier, m_jJetRound_LargeR_jet_min_ET_MeV,/*m_saveSeeds*/ false) );
+ CHECK( JetAlg::BuildRoundJet(jTs_PUsub, "jRoundSeedsPUsub", "jRoundLargeRJetsPUsub", m_jJetRound_LargeR_r, jT_noise, /*m_jJetRound_jet_tower_noise_multiplier*/ 0, m_jJetRound_LargeR_jet_min_ET_MeV, false) );
}
if(m_makeJetsFromMap) {
@@ -506,6 +497,9 @@ StatusCode JGTowerReader::JFexAlg(const xAOD::JGTowerContainer* jTs){
CHECK( BuildJetsFromMap(jTs) );
}
+ delete jTs_PUsub;
+ delete jTs_PUsubAux;
+
return StatusCode::SUCCESS;
}
@@ -542,14 +536,19 @@ StatusCode JGTowerReader::GFexAlg(const xAOD::JGTowerContainer* gTs){
if(fabs(eta) < 3.15 && gt_em->sampling() == 0){
const xAOD::JGTower* gt_had = gTs->at(t+544);
+ if(std::abs(eta-gt_had->eta())>0.01)ATH_MSG_WARNING("Unmatched ETA: " << eta << " " << gt_had->eta());
+ if(std::abs(phi-gt_had->phi())>0.01)ATH_MSG_WARNING("Unmatched PHI: " << phi << " "<< gt_had->phi());
totalEt = gt_em->et() + gt_had->et();
}
//if(t > 544) totalEt = 0;
if(fabs(eta)>=3.15){
//For the case where eta > 3.15, the sampling is always 2. Not sure that is quite what we want, but that is how the
//identifiers are configured for now...
+ if(gTs->at(t)->et()>0 && gTs->at(t)->sampling()!=2)ATH_MSG_WARNING("ET: " << gTs->at(t)->et() << " sampling " << gTs->at(t)->sampling());
totalEt=gTs->at(t)->et();
}
+ if(fabs(eta)<3.15 && gt_em->sampling()!=0)continue;
+ //only build 1 layer for the calo towers to avoid double counting;
const float tEt = totalEt;
const std::vector<int> index(2, 0);
@@ -580,18 +579,11 @@ StatusCode JGTowerReader::GFexAlg(const xAOD::JGTowerContainer* gTs){
const xAOD::JGTower* tower = gCaloTowers->at(t);
float eta = tower->eta();
- if(eta > -2.5 && eta < 0){
- temp_a.push_back(tower);
- h_fpga_a->Fill(tower->et());
- }
- if(eta > 0 && eta < 2.5){
- temp_b.push_back(tower);
- h_fpga_b->Fill(tower->et());
- }
- if(fabs(eta) >= 2.5){
- temp_c.push_back(tower);
- h_fpga_c->Fill(tower->et());
- }
+ std::string FPGA = GFEX_pFPGA(eta);
+ if(FPGA=="A")h_fpga_a->Fill(tower->et());
+ if(FPGA=="B")h_fpga_b->Fill(tower->et());
+ if(FPGA=="C")h_fpga_c->Fill(tower->et());
+
}
const xAOD::JGTowerContainer* fpga_a = temp_a.asDataVector();
const xAOD::JGTowerContainer* fpga_b = temp_b.asDataVector();
@@ -606,11 +598,17 @@ StatusCode JGTowerReader::GFexAlg(const xAOD::JGTowerContainer* gTs){
xAOD::JGTowerContainer* pu_sub = new xAOD::JGTowerContainer();
pu_sub->setStore(pu_subAux);
- float rhoA = METAlg::Rho_avg_etaRings(fpga_a, 1, false);
- float rhoB = METAlg::Rho_avg_etaRings(fpga_b, 2, false);
- float rhoC = METAlg::Rho_avg_etaRings(fpga_c, 3, false);
+ float rhoA = METAlg::Rho_avg_etaRings(fpga_a, m_gFEX_Rho_useNegTowers);
+ float rhoB = METAlg::Rho_avg_etaRings(fpga_b, m_gFEX_Rho_useNegTowers);
+ float rhoC = METAlg::Rho_avg_etaRings(fpga_c, m_gFEX_Rho_useNegTowers);
- float rho_barrel = METAlg::Rho_avg(gCaloTowers, false);
+ if(m_gFEX_OnlyPosRho){
+ if(rhoA<0) rhoA=0;
+ if(rhoB<0) rhoB=0;
+ if(rhoC<0) rhoC=0;
+ }
+
+ float rho_barrel = METAlg::Rho_avg_barrel(gCaloTowers, m_gFEX_Rho_useNegTowers);
(*acc_rhoA)(*RhoCont) = rhoA;
(*acc_rhoB)(*RhoCont) = rhoB;
@@ -629,25 +627,17 @@ StatusCode JGTowerReader::GFexAlg(const xAOD::JGTowerContainer* gTs){
const float dphi = tower->dphi();
float Et_sub = tower->et();
- if(eta < 0 && eta >= -2.5){
- if(eta > -2.4) Et_sub -= rhoA;
- else Et_sub -= 0.5*rhoA;
+ float Area = GTowerArea(eta);
+ float rho = -999e6; float thresh = -999e6;
+ std::string FPGA = GFEX_pFPGA(eta);
+ if(FPGA=="A"){rho = rhoA; thresh = thresh_a;}
+ else if(FPGA=="B"){rho = rhoB; thresh = thresh_b;}
+ else if(FPGA=="C"){rho = rhoC; thresh = thresh_c;}
+ else ATH_MSG_WARNING("No valid FPGA found");
- if(Et_sub < thresh_a) Et_sub = 0;
- }
- if(eta > 0 && eta <= 2.5){
- if(eta < 2.4) Et_sub -= rhoB;
- else Et_sub -= 0.5*rhoB;
-
- if(Et_sub < thresh_b) Et_sub = 0;
- }
- if(fabs(eta) > 2.5){
- if(fabs(eta) < 3.2) Et_sub -= rhoC;
- else Et_sub -= 4*rhoC;
+ Et_sub -= Area*rho;
+ if(Et_sub < thresh) Et_sub = 0; //Apply RMS requirement for MET
- if(Et_sub < thresh_c) Et_sub = 0;
- }
-
const std::vector<int> index(2,0);
const int sampling = tower->sampling();
xAOD::JGTower* new_tower = new xAOD::JGTower();
@@ -659,14 +649,14 @@ StatusCode JGTowerReader::GFexAlg(const xAOD::JGTowerContainer* gTs){
new_tower->setSCIndex(index);
new_tower->setTileIndex(index);
new_tower->setSampling(sampling);
- }
+ }//Loop over towers
//build gBlocks
std::vector<TowerObject::Block> temp_gBlocks;
std::vector<TowerObject::Block> temp_puSub_gBlocks;
- CHECK(BuildBlocksFromTowers(temp_gBlocks, *gCaloTowers, 3, 3, m_useNegTowers));
- CHECK(BuildBlocksFromTowers(temp_puSub_gBlocks, *pu_sub, 3, 3, m_useNegTowers));
+ CHECK(BuildBlocksFromTowers(temp_gBlocks, *gCaloTowers, 3, 3, m_gFEX_useNegTowers));
+ CHECK(BuildBlocksFromTowers(temp_puSub_gBlocks, *pu_sub, 3, 3, m_gFEX_useNegTowers));
const std::vector<TowerObject::Block> gBlocks = temp_gBlocks;
const std::vector<TowerObject::Block> puSub_gBlocks = temp_puSub_gBlocks;
@@ -699,46 +689,20 @@ StatusCode JGTowerReader::GFexAlg(const xAOD::JGTowerContainer* gTs){
}
// jet algorithms
- if(JetAlg::m_SeedMap.find("gSeeds") == JetAlg::m_SeedMap.end())
- CHECK(JetAlg::SeedGrid(gTs,"gSeeds",m_dumpSeedsEtaPhi));
- CHECK(JetAlg::SeedFinding(gTs, "gSeeds", m_gJet_seed_size, m_gJet_max_r, gT_noise, m_gJet_seed_tower_noise_multiplier, m_gJet_seed_total_noise_multiplier, m_gJet_seed_min_ET_MeV));
-
- // CHECK(JetAlg::SeedFinding(gTs,gSeeds,m_gJet_seed_size,m_gJet_max_r,gJet_thr)); // the diameter of seed, and its range to be local maximum
- // Careful to ensure the range set to be no tower double counted
- //CHECK(JetAlg::BuildJet(gTs,gSeeds,gL1Jets,m_gJet_r,gJet_thr)); //default gFex jets are cone jets wih radius of 1.0
- CHECK(JetAlg::BuildFatJet(*gCaloTowers, "gL1Jets", m_gJet_r, gT_noise, m_gJet_jet_tower_noise_multiplier, m_gJet_jet_total_noise_multiplier, m_gJet_jet_min_ET_MeV, rho_barrel));
- if(m_buildgBlockJets) CHECK(JetAlg::BuildgBlocksJets(gBs, "gBlockJets",rho_barrel));
+ CHECK(evtStore()->record(RhoCont,"EventVariables"));
+ CHECK(evtStore()->record(RhoContAux,"EventVariablesAux."));
+
+ CHECK(JetAlg::BuildFatJet(*gCaloTowers, "gL1Jets", m_gJet_r, gT_noise, m_gJet_tower_noise_multiplier,m_gJet_block_min_ET_MeV, m_gJet_jet_min_ET_MeV, rhoA, rhoB, rhoC, m_gFEX_useNegTowers));
+ if(m_buildgBlockJets) CHECK(JetAlg::BuildgBlocksJets(gBs, "gBlockJets",rhoA,rhoB, rhoC));
//gFEX MET algorithms
std::vector<float> noNoise;
- if(m_developerMET){
- for(int i=0; i<=1; i++){
- bool NegTowers=false;
- if(i==1)NegTowers=true;
- CHECK(METAlg::Baseline_MET(gCaloTowers,Form("gNoCutsNeg%d",NegTowers),noNoise, NegTowers) );//basic MET, no threshold, pass no noise thresholds, and use negative towers
- CHECK(METAlg::Baseline_MET(gCaloTowers,Form("gBaselineNeg%d",NegTowers),gT_noise, NegTowers)); //basic MET reconstruction with a 4 sigma noise cut applied
- //No RMS
- //CHECK(METAlg::SubtractRho_MET(gTs,Form("RhoSubNeg%dMedian",NegTowers),0,1,NegTowers));//Median
- //CHECK(METAlg::SubtractRho_MET(gTs,Form("RhoSubNeg%dMean",NegTowers),0,0,NegTowers));//Mean
- //Apply 3 sigma RMS
- //CHECK(METAlg::SubtractRho_MET(gTs,Form("RhoSubNeg%dMedianRMS",NegTowers), 1, 1,NegTowers));
- //CHECK(METAlg::SubtractRho_MET(gTs,Form("RhoSubNeg%dMeanRMS",NegTowers), 1, 0,NegTowers));
-
- CHECK(METAlg::Softkiller_MET(pu_sub, Form("SKNeg%d",NegTowers), NegTowers) ); //pileup subtracted SoftKiller (with avg rho)
- CHECK(METAlg::JwoJ_MET(pu_sub, gBlocks,Form("JwoJNeg%d",NegTowers),m_pTcone_cut, false, false, NegTowers) ); //Jets without Jets
- //CHECK(METAlg::Pufit_MET(gTs,Form("PUfitNeg%d",NegTowers), NegTowers) ); //L1 version of PUfit, using gTowers
- }
-
- }//developer met
- else{
- CHECK(METAlg::Baseline_MET(pu_sub, "gXERHO", noNoise, m_useNegTowers));
- CHECK(METAlg::Baseline_MET(gCaloTowers,"gXENOISECUT",gT_noise, m_useNegTowers));
- CHECK(METAlg::JwoJ_MET(pu_sub, puSub_gBlocks,"gXEJWOJRHO",m_pTcone_cut,false, false, m_useNegTowers));
- CHECK(METAlg::JwoJ_MET(gTs, gBlocks, "gXEJWOJ",m_pTcone_cut,false, true, /*m_useNegTowers*/ true));//by default, m_useNegTowers=false, but setting it = true here for consistency
- CHECK(METAlg::Pufit_MET(gCaloTowers,"gXEPUFIT", m_useNegTowers) );
-
-
- }//main definitions for simplicity
+ CHECK(METAlg::NoiseCut_MET(pu_sub, "gXERHO", noNoise, m_gFEX_useNegTowers));
+ CHECK(METAlg::NoiseCut_MET(gCaloTowers,"gXENOISECUT",gT_noise, m_gFEX_useNegTowers));
+ CHECK(METAlg::JwoJ_MET(pu_sub, puSub_gBlocks,"gXEJWOJRHO",m_gFEX_pTcone_cut, /*bool useRho*/ false,rhoA,rhoB,rhoC, /*m_useNegTowers*/ false));//do not use negative towers, so if the subtracted ET is < 0 you remove the tower
+ CHECK(METAlg::JwoJ_MET(gCaloTowers, gBlocks, "gXEJWOJ",m_gFEX_pTcone_cut,/*bool useRho*/ false,rhoA,rhoB,rhoC, /*m_useNegTowers*/ m_gFEX_useNegTowers));
+ CHECK(METAlg::JwoJ_MET(gCaloTowers, gBlocks, "gXEJWOJRHOHT",m_gFEX_pTcone_cut,/*bool useRho*/ true,rhoA,rhoB,rhoC, /*m_useNegTowers*/ m_gFEX_useNegTowers));
+ CHECK(METAlg::Pufit_MET(gCaloTowers,"gXEPUFIT", m_gFEX_useNegTowers) );
//manage conatiners that have been created: save gCaloTowers and pu_sub to SG
CHECK(evtStore()->record(gCaloTowers, "gCaloTowers"));
@@ -748,9 +712,6 @@ StatusCode JGTowerReader::GFexAlg(const xAOD::JGTowerContainer* gTs){
CHECK(evtStore()->record(gBs, "gBlocks"));
CHECK(evtStore()->record(gBAux, "gBlocksAux."));
- CHECK(evtStore()->record(RhoCont,"EventVariables"));
- CHECK(evtStore()->record(RhoContAux,"EventVariablesAux."));
-
delete h_fpga_a;
delete h_fpga_b;
delete h_fpga_c;
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/src/JetAlg.cxx b/Trigger/TrigT1/TrigT1CaloFexSim/src/JetAlg.cxx
index 9b08d34e9803573550c5d38cb473f11e543e6132..a7526077cc361524034430ce606bf22a16faaec9 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/src/JetAlg.cxx
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/src/JetAlg.cxx
@@ -120,7 +120,7 @@ StatusCode JetAlg::SeedGrid(const xAOD::JGTowerContainer*towers, TString seedNam
}
//To find the seeds as local maxima
-StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seedname, float seed_size, float range, std::vector<float> noise, float seed_tower_noise_multiplier, float seed_total_noise_multiplier, float seed_min_ET_MeV){
+StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seedname, float seed_size, float range, std::vector<float> noise, float seed_tower_noise_multiplier, float seed_min_ET_MeV, bool seed_electrons){
// get the energy of each seeds which is defined as 2x2 towers in barrel and endcap, and single tower in fcal
// static MsgStream staticMsg("MyStaticMsgStream"0);
// static MsgStream staticMsg();
@@ -151,6 +151,9 @@ StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seed
float tower_et = tower->et();
float tower_noise = 0;
bool isTile = (fabs(tower_eta) <1.5 && tower->sampling()==1);
+ bool isHad = tower->sampling()==1;
+ // New logic: do not seed on HAD towers
+ if(seed_electrons && isHad) continue;
if(noise.size() < t) {
ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::ERROR,"JetAlg::SeedFinding") << "the noise vector is smaller (at " << noise.size() << " entries) than the tower number " << t << " that you are attempting to use";
@@ -192,7 +195,7 @@ StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seed
nSeeds++;
}
}
-
+ /*
for(unsigned i=0; i<seeds->eta.size(); i++){
for(unsigned ii=0; ii<seeds->phi.at(i).size(); ii++){
if(seeds->et.at(i).at(ii) < seeds->noise.at(i).at(ii)*seed_total_noise_multiplier){
@@ -203,7 +206,7 @@ StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seed
nSeedsAboveThreshold += 1;
}
}
- }
+ }*/
ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::DEBUG,"JetAlg::SeedFinding") << "in total, found " << nSeedsAboveThreshold << " seeds above noise thresholds out of " << nSeeds << " seeds";
@@ -221,15 +224,20 @@ StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seed
// eta_p: et higher than all seeds with larger eta
// eta_0: et higher than the other seeds along the same eta ring
+ // re-define matching range based on eta for electrons
+ float match_range = range;
+ if(seed_electrons && fabs(seeds->eta.at(iseed_eta)) > 2.5)
+ match_range += 0.1;
+
bool eta_n=1, eta_p=1, eta_0=1;
for(unsigned i=iseed_eta+1; ;i++){
if(i>=seeds->eta.size()) break;
- if(fabs(seeds->eta.at(i)-seeds->eta.at(iseed_eta))>range) break;
+ if(fabs(seeds->eta.at(i)-seeds->eta.at(iseed_eta))>match_range) break;
for(unsigned ii=0; ii<seeds->phi.at(i).size(); ii++){
float dphi = fabs(deltaPhi(seeds->phi.at(iseed_eta).at(iseed_phi),seeds->phi.at(i).at(ii)));
- if(dphi>range) continue;
+ if(dphi>match_range) continue;
if(seeds->et.at(iseed_eta).at(iseed_phi) < seeds->et.at(i).at(ii)){
eta_p = false;
break;
@@ -239,10 +247,10 @@ StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seed
for(int i=iseed_eta-1; ;i--){
if(i<0) break;
- if(fabs(seeds->eta.at(iseed_eta)-seeds->eta.at(i))>range) break;
+ if(fabs(seeds->eta.at(iseed_eta)-seeds->eta.at(i))>match_range) break;
for(unsigned ii=0; ii<seeds->phi.at(i).size(); ii++){
float dphi = fabs(deltaPhi(seeds->phi.at(iseed_eta).at(iseed_phi),seeds->phi.at(i).at(ii)));
- if(dphi>range) continue;
+ if(dphi>match_range) continue;
if(seeds->et.at(iseed_eta).at(iseed_phi) < seeds->et.at(i).at(ii)){
eta_n = false;
break;
@@ -253,7 +261,7 @@ StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seed
for(unsigned ii=0; ii<seeds->phi.at(iseed_eta).size(); ii++){
if(ii==iseed_phi) continue;
float dphi = fabs(deltaPhi(seeds->phi.at(iseed_eta).at(iseed_phi),seeds->phi.at(iseed_eta).at(ii)));
- if(dphi>range) continue;
+ if(dphi>match_range) continue;
if(seeds->et.at(iseed_eta).at(iseed_phi) < seeds->et.at(iseed_eta).at(ii)){
eta_0 = false;
break;
@@ -271,6 +279,11 @@ StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seed
for(unsigned iseed_eta=0; iseed_eta<seeds->eta.size(); iseed_eta++){
for(unsigned iseed_phi=0; iseed_phi<seeds->phi.at(iseed_eta).size(); iseed_phi++){
+ // re-define matching range based on eta for electrons
+ float match_range = range;
+ if(seed_electrons && fabs(seeds->eta.at(iseed_eta)) > 2.5)
+ match_range += 0.1;
+
// only look at local max (or equal max) seeds
if(!seeds->local_max.at(iseed_eta).at(iseed_phi))
continue;
@@ -290,12 +303,12 @@ StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seed
// I now have another seed with the same energy. Check it's within range in case of weird coincidence
float deta = seeds->eta.at(iiseed_eta) - seeds->eta.at(iseed_eta);
- if(fabs(deta) > range)
+ if(fabs(deta) > match_range)
continue;
float dphi = deltaPhi(seeds->phi.at(iiseed_eta).at(iiseed_phi), seeds->phi.at(iseed_eta).at(iseed_phi));
//This dphi should not be absolute value, as we need the sign below!
//The correct usage in deltaPhi is the following function from root: TVector2::Phi_mpi_pi(phi1,phi2)
- if(fabs(dphi) > range)
+ if(fabs(dphi) > match_range)
continue;
// I now have another seed that this one has to compete with
@@ -330,7 +343,7 @@ StatusCode JetAlg::SeedFinding(const xAOD::JGTowerContainer*towers, TString seed
return StatusCode::SUCCESS;
}
-StatusCode JetAlg::BuildFatJet(const xAOD::JGTowerContainer towers, TString jetname, float jet_r, std::vector<float> noise, float jet_tower_noise_multiplier, float /*jet_total_noise_multiplier*/, float jet_min_ET_MeV, float rho){
+StatusCode JetAlg::BuildFatJet(const xAOD::JGTowerContainer towers, TString jetname, float jet_r, std::vector<float> noise, float jet_tower_noise_multiplier,float seed_block_threshold, float jet_min_ET_MeV, float rhoA, float rhoB,float rhoC, bool useNegTowers){
std::vector<TowerObject::Block> blocks;
@@ -340,6 +353,8 @@ StatusCode JetAlg::BuildFatJet(const xAOD::JGTowerContainer towers, TString jetn
TowerObject::TowerGrid grid = TowerObject::TowerGrid(towers);
for(const xAOD::JGTower* seed: towers){
+ //if(fabs(seed->eta())>2.4) continue;
+ if(!useNegTowers && seed->et() <0) continue;
int seedIndex = std::find(towers.begin(), towers.end(), seed) - towers.begin();
std::vector<int> neighbors = grid.neighbors(*seed, 3, 3);
@@ -368,35 +383,47 @@ StatusCode JetAlg::BuildFatJet(const xAOD::JGTowerContainer towers, TString jetn
}
std::sort(blocks.rbegin(), blocks.rend());
- float pt_cone_cut = 25*Gaudi::Units::GeV;
-
for(unsigned b = 0; b < blocks.size(); b++){
// const xAOD::JGTower* seed = towers.at(blocks[b].seedIndex());
float block_phi = blocks[b].Phi();
float block_eta = blocks[b].Eta();
- float pt_cone = blocks[b].Pt();
+ float pt_block = blocks[b].Pt();
+ float rho = -999;
+ std::string FPGA = GFEX_pFPGA(block_eta);
+ if(FPGA=="A") rho = rhoA;
+ else if (FPGA=="B") rho = rhoB;
+ else if (FPGA=="C") rho = rhoC;
+ else ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::WARNING,"JetAlg::BuildFatJet") << "NO FPGA FOUND, ETA" << block_eta ;
+ if(pt_block <= seed_block_threshold) continue;
+
+ int NtowersInJet = 0;
float j_Et = 0;
- if(pt_cone > pt_cone_cut){
- for(unsigned int t = 0; t < towers.size(); t++){
- const xAOD::JGTower* tower = towers.at(t);
+
+ for(unsigned int t = 0; t < towers.size(); t++){
+ const xAOD::JGTower* tower = towers.at(t);
+ if(!useNegTowers && tower->et()<0)continue;
+ if(jet_tower_noise_multiplier>0){
+ //only apply noise cuts if the noise multiplier is non-zero
if(fabs(tower->et()) < noise.at(t) * jet_tower_noise_multiplier) continue;
- if(!withinRadius(block_eta, tower->eta(), block_phi, tower->phi(), jet_r, false) ) continue;
- j_Et += tower->et();
- }//looping over all towers
- j_Et -= 64*rho;
- if(j_Et < jet_min_ET_MeV) continue;
- std::shared_ptr<JetAlg::L1Jet> j = std::make_shared<JetAlg::L1Jet>(block_eta, block_phi, j_Et);
- js.push_back(j);
- }//cone above threshold
- }
+ }
+ if(!withinRadius(block_eta, tower->eta(), block_phi, tower->phi(), jet_r, /*acceptEqual*/ false) ) continue;
+ j_Et += tower->et();
+ NtowersInJet++;
+ }//looping over all towers
+ j_Et -= 69*rho;
+ if(j_Et < jet_min_ET_MeV) continue;
+ std::shared_ptr<JetAlg::L1Jet> j = std::make_shared<JetAlg::L1Jet>(block_eta, block_phi, j_Et);
+ js.push_back(j);
+
+ }//blocks
m_JetMap[jetname] = js;
return StatusCode::SUCCESS;
}
-StatusCode JetAlg::BuildJet(const xAOD::JGTowerContainer*towers, TString seedname, TString jetname, float jet_r, std::vector<float> noise, float jet_tower_noise_multiplier, float jet_total_noise_multiplier, float jet_min_ET_MeV, bool m_saveSeeds){
+StatusCode JetAlg::BuildJet(const xAOD::JGTowerContainer*towers, TString seedname, TString jetname, float jet_r, std::vector<float> noise, float jet_tower_noise_multiplier, float jet_min_ET_MeV, bool m_saveSeeds){
std::shared_ptr<JetAlg::Seed> seeds = m_SeedMap[seedname];
std::vector<std::shared_ptr<JetAlg::L1Jet>> js;
@@ -431,18 +458,19 @@ StatusCode JetAlg::BuildJet(const xAOD::JGTowerContainer*towers, TString seednam
}
float j_et = 0;
- float j_totalnoise = 0;
+ TLorentzVector jvector;
for(unsigned t=0; t<towers->size(); t++){
const xAOD::JGTower* tower = towers->at(t);
if(fabs(tower->et()) < noise.at(t)*jet_tower_noise_multiplier) continue;
if(!inBox(eta,tower->eta(),jet_r, phi, tower->phi(),jet_r)) continue;
j_et += tower->et();
- j_totalnoise += noise.at(t);
+ TLorentzVector towerVector;
+ towerVector.SetPtEtaPhiM(tower->et(),tower->eta(),tower->phi(),0);
+ jvector+=towerVector;
ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::DEBUG,"JetAlg::BuildJet") << " adding tower at (eta,phi)=("<<tower->eta()<<","<<tower->phi()<<")";
}
- ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::DEBUG,"JetAlg::BuildJet") << " final jet has et = " << j_et;
+ ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::DEBUG,"JetAlg::BuildJet") << " final jet has et = " << j_et << " and mass " << jvector.M();
if(j_et<jet_min_ET_MeV) continue;
- if(j_et<j_totalnoise*jet_total_noise_multiplier) continue;
std::shared_ptr<JetAlg::L1Jet> j = std::make_shared<JetAlg::L1Jet>(eta, phi, j_et);
js.push_back(j);
}
@@ -460,7 +488,7 @@ StatusCode JetAlg::BuildJet(const xAOD::JGTowerContainer*towers, TString seednam
}
-StatusCode JetAlg::BuildRoundJet(const xAOD::JGTowerContainer*towers, TString seedname, TString jetname, float jet_r, std::vector<float> noise, float jet_tower_noise_multiplier, float jet_total_noise_multiplier, float jet_min_ET_MeV, bool m_saveSeeds){
+StatusCode JetAlg::BuildRoundJet(const xAOD::JGTowerContainer*towers, TString seedname, TString jetname, float jet_r, std::vector<float> noise, float jet_tower_noise_multiplier,float jet_min_ET_MeV, bool m_saveSeeds){
std::vector<std::shared_ptr<JetAlg::L1Jet>> js;
std::shared_ptr<JetAlg::Seed> seeds = m_SeedMap[seedname];
@@ -490,8 +518,7 @@ StatusCode JetAlg::BuildRoundJet(const xAOD::JGTowerContainer*towers, TString se
if(!seeds->local_max.at(eta_ind).at(phi_ind)) continue;
float j_et = 0;
- float j_totalnoise = 0;
- ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::DEBUG,"JetAlg::BuildRoundJet") << "BuildRoundJet: found a local max seed at (eta,phi)=("<<eta<<","<<phi<<")";
+ ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::DEBUG,"JetAlg::BuildRoundJet") << "BuildRoundJet: found a local max seed at (eta,phi)=("<<eta<<","<<phi<<")";
if(m_saveSeeds && et > 0) {
std::shared_ptr<JetAlg::L1Jet> s = std::make_shared<JetAlg::L1Jet>(eta, phi, et);
@@ -506,13 +533,10 @@ StatusCode JetAlg::BuildRoundJet(const xAOD::JGTowerContainer*towers, TString se
if(applyNoise && fabs(tower->et()) < noise.at(t)*jet_tower_noise_multiplier) continue;
if(!withinRadius(eta, tower->eta(), phi, tower->phi(), jet_r,1)) continue;
j_et += tower->et();
- j_totalnoise += noise.at(t);
ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::DEBUG,"JetAlg::BuildRoundJet") << " adding tower at (eta,phi)=("<<tower->eta()<<","<<tower->phi()<<")";
}
ATH_REPORT_MESSAGE_WITH_CONTEXT(MSG::DEBUG,"JetAlg::BuildRoundJet") << " final jet has et = " << j_et;
if(j_et<jet_min_ET_MeV) continue;
- if(j_et<j_totalnoise*jet_total_noise_multiplier) continue;
-
std::shared_ptr<JetAlg::L1Jet> j = std::make_shared<JetAlg::L1Jet>(eta, phi, j_et);
js.push_back(j);
@@ -531,7 +555,7 @@ StatusCode JetAlg::BuildRoundJet(const xAOD::JGTowerContainer*towers, TString se
return StatusCode::SUCCESS;
}
-StatusCode JetAlg::BuildgBlocksJets(const xAOD::JGTowerContainer* gBs, TString jetname, float rho){
+StatusCode JetAlg::BuildgBlocksJets(const xAOD::JGTowerContainer* gBs, TString jetname, float rhoA, float rhoB, float rhoC){
std::vector<std::shared_ptr<JetAlg::L1Jet>> gbJets;
if(m_JetMap.find(jetname)!=m_JetMap.end()) gbJets = m_JetMap[jetname];
@@ -612,17 +636,27 @@ StatusCode JetAlg::BuildgBlocksJets(const xAOD::JGTowerContainer* gBs, TString j
}
}
- rho *= 9; //assuming there are ~9 towers per jet
for(unsigned int i=0; i<max_pt.size(); i++){
int i_max = max_index.at(i);
int i_sec = sec_index.at(i);
float eta_max = gBs->at(i_max)->eta();
float phi_max = gBs->at(i_max)->phi();
+ float rho = -999;
+ std::string FPGA = GFEX_pFPGA(eta_max);
+ if(FPGA=="A") rho = rhoA;
+ else if (FPGA=="B") rho = rhoB;
+ else if (FPGA=="C") rho = rhoC;
+ rho *=9;
float et_max = gBs->at(i_max)->et() - rho;
float eta_sec = gBs->at(i_sec)->eta();
float phi_sec = gBs->at(i_sec)->phi();
+ FPGA = GFEX_pFPGA(eta_sec);
+ if(FPGA=="A") rho = rhoA;
+ else if (FPGA=="B") rho = rhoB;
+ else if (FPGA=="C") rho = rhoC;
+ rho*=9;
float et_sec = gBs->at(i_sec)->et() - rho;
std::shared_ptr<JetAlg::L1Jet> gbJ_max = std::make_shared<JetAlg::L1Jet>(eta_max, phi_max, et_max);
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/src/METAlg.cxx b/Trigger/TrigT1/TrigT1CaloFexSim/src/METAlg.cxx
index 9193e577f87c4a7729c5255aa0342e04c30c5d8d..bf07ab49fc5fc41f46c1100f1c3a1511313df125 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/src/METAlg.cxx
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/src/METAlg.cxx
@@ -11,7 +11,7 @@
#include "TrigT1CaloFexSim/JGTower.h"
//MET algorithms
-#include "TrigT1CaloFexSim/Rho.h"
+
#include "TrigT1CaloFexSim/Softkiller.h"
#include "TrigT1CaloFexSim/JwoJ.h"
#include "TrigT1CaloFexSim/Pufit.h"
@@ -158,7 +158,7 @@ StatusCode METAlg::build_jFEX_bins( std::vector < std::vector < int > > &bins, s
}
//------------------------------------------------------------------------------------------------
-StatusCode METAlg::Baseline_MET(const xAOD::JGTowerContainer*towers, TString metName, std::vector<float> noise, bool useNegTowers){
+StatusCode METAlg::NoiseCut_MET(const xAOD::JGTowerContainer*towers, TString metName, std::vector<float> noise, bool useNegTowers){
float met_x=0;
float met_y=0;
@@ -273,24 +273,19 @@ StatusCode METAlg::jXERHO(const xAOD::JGTowerContainer* towers, TString metName,
return StatusCode::SUCCESS;
}
-StatusCode METAlg::SubtractRho_MET(const xAOD::JGTowerContainer* towers, TString metName, bool useEtaBins, bool useRMS, bool useMedian, bool useNegTowers){
-
+StatusCode METAlg::SubtractRho_MET(const xAOD::JGTowerContainer* towers, TString metName, bool useEtaBins, bool useRMS, bool useNegTowers){
+ /*
+ This function is currently depreciated, as we currently compute a set of towers after pileup
+ subtraction that are passed to the NoiseCut_MET function with a set of noise cuts identical to 0
+ However, because it may be useful in the future for a slightly different algorithm, we leave the function here for now.
+ B. Carlson March 27, 2020
+ */
float EtMiss = 0;
float n = 3.;
float Ex = 0, Ey = 0, Ex_ = 0, Ey_ = 0;
float threshold = 0;
- int fpga = 0;
- //can calculate rho as either the average or median gTower energy in the barrel
- if(metName.Contains("RhoSubA")) fpga = 1;
- if(metName.Contains("RhoSubB")) fpga = 2;
- if(metName.Contains("RhoSubC")) fpga = 3;
-
- //std::cout<<metName<<": "<<fpga<<std::endl;
-
- float rho = Rho_bar(towers, useEtaBins, fpga, 0);
-
- if(useMedian) rho = Rho_med(towers, useNegTowers);
+ float rho = Rho_avg_barrel(towers,useNegTowers);
unsigned int size = towers->size();
@@ -304,11 +299,11 @@ StatusCode METAlg::SubtractRho_MET(const xAOD::JGTowerContainer* towers, TString
float Et = tower->et();
if(!useNegTowers && Et < 0) continue;
+ //if(Et > et_max)continue;
h_Et->Fill(Et);
}
threshold = n*(h_Et->GetRMS());
- //std::cout<<"3 sigma = "<<threshold<<std::endl;
}
for(unsigned int t = 0; t < size; t++){
@@ -318,26 +313,9 @@ StatusCode METAlg::SubtractRho_MET(const xAOD::JGTowerContainer* towers, TString
float phi = tower->phi();
float eta = tower->eta();
- //std::cout<<"Looking at tower ("<<eta<<", "<<phi<<")"<<std::endl;
if(!useNegTowers && Et < 0) continue;
- //float deta = tower->deta();
- //float dphi = tower->dphi();
-
- float area = 0.;//deta*dphi/default_area;
- if(eta >= -2.5 && eta < 0 && fpga == 1){
- if(eta > -2.4) area = 1.;
- else area = 0.5;
- }
- if(eta > 0 && eta <= 2.5 && fpga == 2){
- if(eta < 2.4) area = 1.;
- else area = 0.5;
- }
- if(fabs(eta) > 2.5 && fpga == 3){
- if(fabs(eta) < 3.2) area = 1.;
- else area = 4.;
- }
- //std::cout<<"area = "<<area<<std::endl;
+ float area = GTowerArea(eta);
float Et_sub = 0;
if(useEtaBins) Et_sub = TMath::Abs(Et) - area*rho;
@@ -383,7 +361,12 @@ StatusCode METAlg::SubtractRho_MET(const xAOD::JGTowerContainer* towers, TString
//------------------------------------------------------------------------------------------------
StatusCode METAlg::Softkiller_MET(const xAOD::JGTowerContainer* towers, TString metName, bool useNegTowers){
-
+ /*
+ This function is currently depreciated, but could be updated fairly easily. March 2020
+ This function is designed to compute MET after applying SK
+ First the standard rho correction is applied
+ Then towers that have ET > median are included in the MET calculation
+ */
float median = 0;
unsigned int size = towers->size();
@@ -395,8 +378,7 @@ StatusCode METAlg::Softkiller_MET(const xAOD::JGTowerContainer* towers, TString
//find the true median of the tower energies
median = Et_median_true(&EtMaxPerPatch);
-
- float rho = Rho_bar(towers, false,0,0);
+ float rho = Rho_avg_barrel(towers, useNegTowers);
float Ex = 0, Ey = 0, Ex_ = 0, Ey_ = 0;
for(unsigned int t = 0; t < size; t++){
@@ -436,15 +418,14 @@ StatusCode METAlg::Softkiller_MET(const xAOD::JGTowerContainer* towers, TString
return StatusCode::SUCCESS;
}
-StatusCode METAlg::JwoJ_MET(const xAOD::JGTowerContainer* towers, const std::vector<TowerObject::Block> gBlocks, TString metName, float pTcone_cut, bool useEtaBins, bool useRho, bool useNegTowers){
-
+StatusCode METAlg::JwoJ_MET(const xAOD::JGTowerContainer* towers, const std::vector<TowerObject::Block> gBlocks, TString metName, float pTcone_cut, bool useRho, float RhoA, float RhoB, float RhoC, bool useNegTowers){
- //unsigned int size = towers->size();
+ /*
+ Function retrieve components of jets without jets
+ Then evaluate the calibrated version using the coefficients below
+ */
- float rho = 0;
- if(useRho) rho = Rho_bar(towers, useEtaBins, 0, false);
- else rho = 0;
- std::vector<float> Et_values = Run_JwoJ(towers, gBlocks, rho, pTcone_cut, useEtaBins, useNegTowers);
+ std::vector<float> Et_values = Run_JwoJ(towers, gBlocks, pTcone_cut, useRho, RhoA, RhoB, RhoC, useNegTowers);
//set fit parameters for calculating MET
//Look up table for parameters a,b,c depending on scalar sumEt
@@ -578,7 +559,7 @@ StatusCode METAlg::Pufit_MET(const xAOD::JGTowerContainer*towers, TString metNam
return StatusCode::SUCCESS;
}
-float METAlg::Rho_avg(const xAOD::JGTowerContainer* towers, const bool useNegTowers){
+float METAlg::Rho_avg_barrel(const xAOD::JGTowerContainer* towers, const bool useNegTowers){
float rho = 0;
float et_max = 10*Gaudi::Units::GeV;
@@ -593,6 +574,7 @@ float METAlg::Rho_avg(const xAOD::JGTowerContainer* towers, const bool useNegTow
if(!useNegTowers && Et < 0) continue;
if(TMath::Abs(eta) < 2.4){
+ //estimating rho using the barrel, excluding the one tower that is 0.1x0.2
if(Et < et_max){
rho+=Et;
length++;
@@ -603,44 +585,38 @@ float METAlg::Rho_avg(const xAOD::JGTowerContainer* towers, const bool useNegTow
return rho_bar;
}
-float METAlg::Rho_avg_etaRings(const xAOD::JGTowerContainer* towers, int fpga, bool useNegTowers){
+float METAlg::Rho_avg_etaRings(const xAOD::JGTowerContainer* towers, bool useNegTowers){
+ /*
+ This function accounts for the different area of trigger towers, and
+ returns the average energy per tower, in units of the area for 0.2x0.2 towers.
+ To use this energy for towers that are larger or smaller, multiply by the ratio
+ of towers sizes, e.g., 2, 4, etc.
+ */
+
float rho = 0;
float et_max = 10*Gaudi::Units::GeV;
const unsigned int size = towers->size();
- float area_a = 0;
- float area_b = 0;
- float area_c = 0;
+ float A_covered = 0; //Number of towers that pass selections E < max and >0 (if we allow negative towers)
for(unsigned int i = 0; i < size; i++){
const xAOD::JGTower* tower = towers->at(i);
float eta = tower->eta();
float Et = tower->et();
+ float Area = GTowerArea(eta);
if(!useNegTowers && Et < 0) continue;
+ if(Et > et_max)continue;
- if(eta < 0 && eta >=-2.5){
- if(eta > -2.4) area_a += 1.;
- else area_a += 0.5;
- }
- if(eta <= 2.5 && eta > 0){
- if(eta < 2.4) area_b += 1.;
- else area_b += 0.5;
- }
- if(fabs(eta) > 2.5){
- if(fabs(eta) < 3.2) area_c += 1.;
- else area_c += 4.;
- }
-
- if(Et < et_max) rho+=Et; //sum ET together to compute total
+ A_covered += Area;
+ rho+=Et; //Total ET
}
- //after computing total, divide by total area
- if(fpga == 1) rho=rho/area_a;
- if(fpga == 2) rho=rho/area_b;
- if(fpga == 3) rho=rho/area_c;
+ //Note that this function is designed to compute the area
+ //for a collection of towers of a defined eta range.
+ rho = rho/A_covered;
return rho;
}
diff --git a/Trigger/TrigT1/TrigT1CaloFexSim/src/components/TrigT1CaloFexSim_entries.cxx b/Trigger/TrigT1/TrigT1CaloFexSim/src/components/TrigT1CaloFexSim_entries.cxx
index 5a5d6464a1bd821ff50a86167d4b32ecdfae4ac6..ec10087a78985f9e8bf3b243836d3651222dbed7 100644
--- a/Trigger/TrigT1/TrigT1CaloFexSim/src/components/TrigT1CaloFexSim_entries.cxx
+++ b/Trigger/TrigT1/TrigT1CaloFexSim/src/components/TrigT1CaloFexSim_entries.cxx
@@ -9,8 +9,13 @@ DECLARE_ALGORITHM_FACTORY( JGTowerReader )
#include "TrigT1CaloFexSim/JGTowerMaker.h"
DECLARE_ALGORITHM_FACTORY( JGTowerMaker )
+
+#include "TrigT1CaloFexSim/JFexEleTau.h"
+DECLARE_ALGORITHM_FACTORY( JFexEleTau )
+
DECLARE_FACTORY_ENTRIES( TrigT1CaloFexSim )
{
+ DECLARE_ALGORITHM( JFexEleTau );
DECLARE_ALGORITHM( JGTowerMaker );
DECLARE_ALGORITHM( JGTowerReader );
}