Adding cell FEX code in new style, including noise cuts

Trigger/TrigAlgorithms/TrigEFMissingET/CMakeLists.txt

@@ -80,6 +80,7 @@ atlas_add_component( TrigEFMissingET
   src/FexBase.cxx
   src/MonGroupBuilder.cxx
   src/TrkMHTFex.cxx
+  src/CellFex.cxx
   src/components/TrigEFMissingET_entries.cxx
   src/components/TrigEFMissingET_load.cxx
   INCLUDE_DIRS ${ROOT_INCLUDE_DIRS} ${TDAQ-COMMON_INCLUDE_DIRS} ${FASTJETCONTRIB_INCLUDE_DIRS} ${FASTJET_INCLUDE_DIRS}

Trigger/TrigAlgorithms/TrigEFMissingET/src/CellFex.cxx

+/*
+ * Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+ */
+
+/******************************************************************************
+ * @package Trigger/TrigAlgorithms/TrigEFMissingET
+ * @file CellFex.cxx
+ *
+ * Implementation of cell fex class
+ * @author Jon Burr
+ *****************************************************************************/
+
+#include "CellFex.h"
+#include "TrigEFMissingET/METComponent.h"
+#include "TrigEFMissingET/SignedKinematics.h"
+#include "StoreGate/ReadCondHandle.h"
+#include "CaloGeoHelpers/CaloSampling.h"
+#include <array>
+#include <numeric>
+
+// TODO - we should check a couple of things
+// - What is the 'prefetch' used by the TopoClusterMaker? Do we need to use it?
+// - Should we be using the 2-Gaussian noise estimate or not? (Right now we do)
+// - Should we be skipping 'bad' cells in the loop? (Right now we include them
+// in the sum)
+// - Previously we had the bytestream error flag in the output objects. I can't
+// tell if there's any equivalent in the new access mechanism
+
+namespace {
+  /// The number of samplings that we use in this algorithm. We do not use the
+  /// 'mini FCal' samplings (that appear after FCAL2 in the enum) so we ensure
+  /// that we stop after the FCAL
+  static constexpr std::size_t N_SAMPLINGS{CaloSampling::FCAL2+1};
+} //> end anonymous namespace
+
+namespace HLT { namespace MET {
+  CellFex::CellFex(const std::string& name, ISvcLocator* pSvcLocator) :
+    FexBase(name, pSvcLocator)
+  {}
+
+  StatusCode CellFex::initialize()
+  {
+    CHECK( m_cellsKey.initialize() );
+    CHECK( m_noiseCDOKey.initialize() );
+    CHECK( detStore()->retrieve(m_caloCellID, "CaloCell_ID") );
+    // Build up the list of component names
+    std::vector<std::string> componentNames;
+    componentNames.reserve(N_SAMPLINGS);
+    for (std::size_t ii = 0; ii < N_SAMPLINGS; ++ii)
+      componentNames.push_back(CaloSampling::getSamplingName(ii) );
+    return initializeBase(componentNames);
+  }
+
+  StatusCode CellFex::fillMET(
+      xAOD::TrigMissingET& met,
+      const EventContext& context,
+      MonGroupBuilder&) const
+  {
+    // Retrieve the inputs
+    auto cells = SG::makeHandle(m_cellsKey, context);
+    // NB - there's no makeHandle overload for ReadCondHandle
+    SG::ReadCondHandle noiseCDO(m_noiseCDOKey, context);
+
+    // Prepare the individual components
+    std::array<METComponent, N_SAMPLINGS> sums;
+    // Iterate over the calorimeter cells
+    for (const CaloCell* icell : *cells) {
+      // Get the noise. The two different calls are equivalent for LAr cells,
+      // but do differ for the TileCal. As far as I can see, the 'two Gaussian'
+      // approach is the more recommended one.
+      float noise = m_doTwoGaussianNoise ? 
+        noiseCDO->getEffectiveSigma(icell->ID(), icell->gain(), icell->energy() ) :
+        noiseCDO->getNoise(icell->ID(), icell->gain() );
+      // Noise selections, first |E| < T1*S
+      if (m_absNoiseThreshold > 0 &&
+          std::abs(icell->energy() ) < m_absNoiseThreshold*noise)
+        continue;
+      // Then E > -T2*S
+      if (m_negNoiseThreshold > 0 &&
+          icell->energy() < -m_negNoiseThreshold*noise)
+        continue;
+      // What about bad cells?
+      if (const CaloDetDescrElement* dde = icell->caloDDE() ) {
+        // Get the right component
+        METComponent& sum = sums.at(dde->getSampling() );
+        sum += SignedKinematics::fromEnergyEtaPhi(
+            icell->energy(), dde->eta(), dde->phi() );
+      }
+      else {
+        auto id = icell->ID();
+        METComponent& sum = sums.at(m_caloCellID->sampling(id) );
+        sum += SignedKinematics::fromEnergyEtaPhi(
+            icell->energy(), m_caloCellID->eta(id), m_caloCellID->phi(id) );
+      }
+    }
+    // Save the full sum
+    std::accumulate(sums.begin(), sums.end(), METComponent{}).fillMET(met);
+    // Save each component
+    for (std::size_t ii = 0; ii < N_SAMPLINGS; ++ii)
+      sums.at(ii).fillMETComponent(ii, met);
+
+    return StatusCode::SUCCESS;
+  }
+} } //> end namespace HLT::MET

Trigger/TrigAlgorithms/TrigEFMissingET/src/CellFex.h

+/*
+ * Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+ */
+
+/******************************************************************************
+ * @package Trigger/TrigAlgorithms/TrigEFMissingET
+ * @class   CellFex
+ *
+ * @brief Fex class for the basic cell algorithm
+ * @author Jon Burr
+ *****************************************************************************/
+
+#ifndef TRIGEFMISSINGET_CELLFEX_H
+#define TRIGEFMISSINGET_CELLFEX_H 1
+
+#include "FexBase.h"
+#include "CaloEvent/CaloCellContainer.h"
+#include "CaloConditions/CaloNoise.h"
+#include "StoreGate/ReadCondHandleKey.h"
+#include "CaloIdentifier/CaloCell_ID.h"
+
+
+namespace HLT { namespace MET {
+  /****************************************************************************
+   * @class CellFex
+   *
+   * Class to create output from the cell algorithm
+   *
+   * cell calculates the MET with a sum over all calorimeter cells. A two-sided
+   * noise cut is applied, requiring cells to have |E/S| > T1, E/S > -T2, where
+   * T1 and T2 are thresholds, and S is the 1-sigma noise threshold.
+   ***************************************************************************/
+  class CellFex : public FexBase
+  {
+    public:
+      /// Constructor
+      CellFex(const std::string& name, ISvcLocator* pSvcLocator);
+
+      /// Initialize the fex
+      virtual StatusCode initialize() override;
+
+    private:
+      /************************************************************************
+       * Properties
+       ***********************************************************************/
+      /// Input cells
+      SG::ReadHandleKey<CaloCellContainer> m_cellsKey{
+        this, "CellName", "CaloCells", "Collection containing all input cells"};
+      /// Calorimeter noise CDO (conditions data object)
+      SG::ReadCondHandleKey<CaloNoise> m_noiseCDOKey{
+        this, "CaloNoiseName","totalNoise","SG Key of CaloNoise data object"};
+      /// The threshold on the magnitude of the cell energy
+      Gaudi::Property<float> m_absNoiseThreshold{
+        this, "AbsoluteNoiseThreshold", 2, "Threshold on the magnitude of the "
+          "cell energy (as a multiple of the cell noise level). Selection "
+          "will not be applied if value is negative"};
+      /// The maximum negative cell energy
+      Gaudi::Property<float> m_negNoiseThreshold{
+        this, "NegativeNoiseThreshold", 5, "The maximum negative cell energy. "
+          "Selection will not be applied if value is negative"};
+      /// Use the 'two-gaussian' noise calculation for the TileCal
+      Gaudi::Property<bool> m_doTwoGaussianNoise{
+        this, "TwoGaussianNoise", true,
+          "Whether to use the 'two-Gaussian' noise calculation for the TileCal"};
+
+      /************************************************************************
+       * Internal functions
+       ***********************************************************************/
+      /**
+       * @brief Calculate and fill the output MET value
+       * @param met The object to fill
+       * @param context The event context
+       * @param monitors[out] Any extra monitors to fill
+       */
+      virtual StatusCode fillMET(
+          xAOD::TrigMissingET& met,
+          const EventContext& context,
+          MonGroupBuilder& monitors) const override;
+
+      /************************************************************************
+       * Data members
+       ***********************************************************************/
+      /// Fallback option for calo cells which don't have a detector description
+      /// - we can read the information from this object
+      const CaloCell_ID* m_caloCellID{nullptr};
+  }; //> end class FexBase
+} } //> end namespace HLT::MET
+
+#endif //> !TRIGEFMISSINGET_CELLFEX_H

Trigger/TrigAlgorithms/TrigEFMissingET/src/components/TrigEFMissingET_entries.cxx

@@ -20,6 +20,7 @@
 #include "../EFMissingETFlagsMT.h"
 #include "../EFMissingETComponentCopier.h"
 #include "../TrkMHTFex.h"
+#include "../CellFex.h"
 
 DECLARE_COMPONENT( EFMissingET )
 DECLARE_COMPONENT( EFMissingETBaseTool )
@@ -42,3 +43,4 @@ DECLARE_COMPONENT( EFMissingETFromClustersPufitMT )
 DECLARE_COMPONENT( EFMissingETFromJetsMT )
 DECLARE_COMPONENT( EFMissingETFlagsMT )
 DECLARE_COMPONENT( HLT::MET::TrkMHTFex )
+DECLARE_COMPONENT( HLT::MET::CellFex )

Trigger/TriggerCommon/TriggerMenuMT/python/HLTMenuConfig/MET/METRecoSequences.py

@@ -6,7 +6,8 @@ from TrigEDMConfig.TriggerEDMRun3 import recordable
 from TriggerMenuMT.HLTMenuConfig.Menu.MenuComponents import RecoFragmentsPool
 
 from TrigEFMissingET.TrigEFMissingETConf import (
-        EFMissingETAlgMT, EFMissingETFlagsMT, HLT__MET__TrkMHTFex)
+        EFMissingETAlgMT, EFMissingETFlagsMT, HLT__MET__TrkMHTFex,
+        HLT__MET__CellFex)
 from TrigEFMissingET.TrigEFMissingETMTConfig import getMETMonTool
 
 from TrigT2CaloCommon.CaloDef import clusterFSInputMaker
@@ -26,26 +27,16 @@ def metCellRecoSequence():
 
     from TrigT2CaloCommon.CaloDef import HLTFSCellMakerRecoSequence
     (cellMakerSeq, CellsName) = HLTFSCellMakerRecoSequence()
-
-    #################################################
-    # Add EFMissingETAlg and associated tools
-    #################################################
-    metAlg = EFMissingETAlgMT( name="EFMET_cell" )
-    flagsTool = EFMissingETFlagsMT("theFlagsTool")
-    metAlg.METContainerKey = recordable("HLT_MET_cell")
-    metAlg.MonTool = getMETMonTool()
-
-    #///////////////////////////////////////////
-    # Add EFMissingETFromCells tool
-    #///////////////////////////////////////////
-    from TrigEFMissingET.TrigEFMissingETConf import EFMissingETFromCellsMT
-    cellTool = EFMissingETFromCellsMT( name="METFromCellsTool" )
-    cellTool.CellsCollection = CellsName
-    metAlg.METTools = [cellTool, flagsTool]
-
-    met_recoSequence = seqAND("metCellRecoSequence", [cellMakerSeq, metAlg])
-
-    seqOut = metAlg.METContainerKey
+    alg = HLT__MET__CellFex(
+            name="EFMET_cell",
+            CellName = CellsName,
+            METContainerKey = recordable("HLT_MET_cell"),
+            AbsoluteNoiseThreshold = -1,
+            NegativeNoiseThreshold = -1,
+            MonTool = getMETMonTool() )
+
+    met_recoSequence = seqAND("metCellRecoSequence", [cellMakerSeq, alg])
+    seqOut = alg.METContainerKey
     return (met_recoSequence, seqOut)