Merge branch 'BCMoverlay' into 'master'

BCM overlay

See merge request atlas/athena!33218

Event/EventOverlay/EventOverlayJobTransforms/share/BeamOverlay_jobOptions.py

@@ -9,3 +9,6 @@ if DetFlags.overlay.BCM_on():
 
     from AthenaCommon import CfgGetter
     job += CfgGetter.getAlgorithm("BCM_OverlayDigitization")
+    job += CfgGetter.getAlgorithm("BCMOverlay")
+    if DetFlags.overlay.Truth_on():
+        job += CfgGetter.getAlgorithm("BCMSDOOverlay")
\ No newline at end of file

InnerDetector/InDetDigitization/BCM_Digitization/python/BCM_DigitizationConfig.py

@@ -31,15 +31,26 @@ def BCM_DigitizationTool(name="BCM_DigitizationTool",**kwargs):
  
     if digitizationFlags.doXingByXingPileUp():
         kwargs.setdefault("FirstXing", BCM_FirstXing() )
-        kwargs.setdefault("LastXing",  BCM_LastXing()  ) 
+        kwargs.setdefault("LastXing",  BCM_LastXing()  )
 
-    if digitizationFlags.PileUpPremixing and 'OverlayMT' in digitizationFlags.experimentalDigi():
+    from AthenaCommon.GlobalFlags import globalflags
+    if globalflags.isOverlay():
         from OverlayCommonAlgs.OverlayFlags import overlayFlags
-        kwargs.setdefault("OutputRDOKey", overlayFlags.bkgPrefix() + "BCM_RDOs")
-        kwargs.setdefault("OutputSDOKey", overlayFlags.bkgPrefix() + "BCM_SDO_Map")
+        if overlayFlags.isOverlayMT():
+            kwargs.setdefault("OnlyUseContainerName", False)
+            kwargs.setdefault("OutputRDOKey", overlayFlags.sigPrefix() + "BCM_RDOs")
+            kwargs.setdefault("OutputSDOKey", overlayFlags.sigPrefix() + "BCM_SDO_Map")
+        else:
+            kwargs.setdefault("OutputRDOKey", overlayFlags.evtStore() + "+BCM_RDOs")
+            kwargs.setdefault("OutputSDOKey", overlayFlags.evtStore() + "+BCM_SDO_Map")
     else:
-        kwargs.setdefault("OutputRDOKey", "BCM_RDOs")
-        kwargs.setdefault("OutputSDOKey", "BCM_SDO_Map")
+        if digitizationFlags.PileUpPremixing and 'OverlayMT' in digitizationFlags.experimentalDigi():
+            from OverlayCommonAlgs.OverlayFlags import overlayFlags
+            kwargs.setdefault("OutputRDOKey", overlayFlags.bkgPrefix() + "BCM_RDOs")
+            kwargs.setdefault("OutputSDOKey", overlayFlags.bkgPrefix() + "BCM_SDO_Map")
+        else:
+            kwargs.setdefault("OutputRDOKey", "BCM_RDOs")
+            kwargs.setdefault("OutputSDOKey", "BCM_SDO_Map")
 
     from AthenaCommon import CfgMgr
     return CfgMgr.BCM_DigitizationTool(name,**kwargs)

InnerDetector/InDetRawAlgs/InDetOverlay/CMakeLists.txt

@@ -15,6 +15,7 @@ atlas_depends_on_subdirs( PUBLIC
                           Control/StoreGate
                           DetectorDescription/IdDictParser
                           Generators/GeneratorObjects
+                          InnerDetector/InDetRawEvent/InDetBCM_RawData
                           InnerDetector/InDetDetDescr/InDetIdentifier
                           InnerDetector/InDetRawEvent/InDetSimData
                           InnerDetector/InDetRecTools/TRT_ElectronPidTools
@@ -48,7 +49,7 @@ atlas_add_component( InDetOverlay
                      src/*.cxx
                      src/components/*.cxx
                      INCLUDE_DIRS ${CLHEP_INCLUDE_DIRS}
-                     LINK_LIBRARIES ${CLHEP_LIBRARIES} AthenaBaseComps IDC_OverlayBase GaudiKernel InDetRawData StoreGateLib SGtests GeneratorObjects InDetIdentifier InDetSimData TrkTrack TRT_ConditionsServicesLib)
+                     LINK_LIBRARIES ${CLHEP_LIBRARIES} AthenaBaseComps IDC_OverlayBase GaudiKernel InDetRawData StoreGateLib SGtests GeneratorObjects InDetBCM_RawData InDetIdentifier InDetSimData TrkTrack TRT_ConditionsServicesLib)
 
 # Install files from the package:
 atlas_install_headers( InDetOverlay )
@@ -67,3 +68,8 @@ atlas_add_test( SCTOverlayConfig_test
 atlas_add_test( TRTOverlayConfig_test
                 SCRIPT test/TRTOverlayConfig_test.py
                 PROPERTIES TIMEOUT 300 )
+
+atlas_add_test( BCMOverlay_test
+                SOURCES test/BCMOverlay_test.cxx src/BCMOverlay.cxx
+                INCLUDE_DIRS ${CLHEP_INCLUDE_DIRS} ${GTEST_INCLUDE_DIRS}
+                LINK_LIBRARIES ${CLHEP_LIBRARIES} AthenaBaseComps GaudiKernel StoreGateLib SGtests GeneratorObjects InDetBCM_RawData InDetSimData ${GTEST_LIBRARIES} )

InnerDetector/InDetRawAlgs/InDetOverlay/InDetOverlay/BCMOverlay.h

+/*
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
+*/
+
+/** @class BCMOverlay
+
+  BCM overlay algorithm
+      @author  Jakob Novak <jakob.novak@cern.ch>
+
+*/
+
+#ifndef INDETOVERLAY_BCMOVERLAY_H
+#define INDETOVERLAY_BCMOVERLAY_H
+
+#include "AthenaBaseComps/AthReentrantAlgorithm.h"
+#include "InDetBCM_RawData/BCM_RDO_Container.h"
+
+struct BCM_Pulse {
+  BCM_Pulse(unsigned int p_, unsigned int w_) {p = p_; w = w_;};
+  unsigned int p;
+  unsigned int w;
+};
+
+class BCMOverlay : public AthReentrantAlgorithm
+{
+public:
+
+  BCMOverlay(const std::string &name, ISvcLocator *pSvcLocator);
+
+  virtual StatusCode initialize() override final;
+  virtual StatusCode execute(const EventContext& ctx) const override final;
+
+  StatusCode overlayContainer(const BCM_RDO_Container *bkgContainer, 
+                              const BCM_RDO_Container *signalContainer, 
+                              BCM_RDO_Container *outputContainer) const;
+
+private:
+
+  BCM_RawData *mergeChannel(const BCM_RawData *bkgRDO, 
+                            const BCM_RawData *signalRDO) const;
+
+  void overlayPulses(std::vector<std::unique_ptr<BCM_Pulse>>& merged_pulses) const;
+
+  std::pair<BCM_Pulse*, BCM_Pulse*> timeOrder(BCM_Pulse* pulse1, BCM_Pulse* pulse2) const;
+
+  static bool compare(const std::unique_ptr<BCM_Pulse>& a, const std::unique_ptr<BCM_Pulse>& b);
+
+  SG::ReadHandleKey<BCM_RDO_Container> m_bkgInputKey{ this, "BkgInputKey", "Bkg_BCM_RDOs", "ReadHandleKey for Background Input BCM_RDO_Container" };
+  SG::ReadHandleKey<BCM_RDO_Container> m_signalInputKey{ this, "SignalInputKey", "Sig_BCM_RDOs", "ReadHandleKey for Signal Input BCM_RDO_Container" };
+  SG::WriteHandleKey<BCM_RDO_Container> m_outputKey{ this, "OutputKey", "BCM_RDOs", "WriteHandleKey for Output BCM_RDO_Container" };
+
+};
+
+#endif // INDETOVERLAY_BCMOVERLAY_H

InnerDetector/InDetRawAlgs/InDetOverlay/python/InDetOverlayConfig.py

@@ -110,3 +110,35 @@ def getTRTSDOOverlay(name="TRTSDOOverlay", **kwargs):
         kwargs.setdefault("OutputKey", overlayFlags.outputStore() + "+TRT_SDO_Map");
 
     return CfgMgr.InDetSDOOverlay(name, **kwargs)
+
+
+def getBCMOverlay(name="BCMOverlay", **kwargs):
+    from OverlayCommonAlgs.OverlayFlags import overlayFlags
+
+    if overlayFlags.isOverlayMT():
+        kwargs.setdefault("BkgInputKey", overlayFlags.bkgPrefix() + "BCM_RDOs");
+        kwargs.setdefault("SignalInputKey", overlayFlags.sigPrefix() + "BCM_RDOs");
+        kwargs.setdefault("OutputKey", "BCM_RDOs");
+    else:
+        kwargs.setdefault("BkgInputKey", overlayFlags.dataStore() + "+BCM_RDOs");
+        kwargs.setdefault("SignalInputKey", overlayFlags.evtStore() + "+BCM_RDOs");
+        kwargs.setdefault("OutputKey", overlayFlags.outputStore() + "+BCM_RDOs");
+
+    return CfgMgr.BCMOverlay(name, **kwargs)
+
+
+def getBCMSDOOverlay(name="BCMSDOOverlay", **kwargs):
+    from OverlayCommonAlgs.OverlayFlags import overlayFlags
+
+    # We do not need background pixel SDOs
+    kwargs.setdefault("BkgInputKey", "");
+
+    if overlayFlags.isOverlayMT():
+        kwargs.setdefault("SignalInputKey", overlayFlags.sigPrefix() + "BCM_SDO_Map");
+        kwargs.setdefault("OutputKey", "BCM_SDO_Map");
+    else:
+        kwargs.setdefault("SignalInputKey", overlayFlags.evtStore() + "+BCM_SDO_Map");
+        kwargs.setdefault("OutputKey", overlayFlags.outputStore() + "+BCM_SDO_Map");
+
+    return CfgMgr.InDetSDOOverlay(name, **kwargs)
+

InnerDetector/InDetRawAlgs/InDetOverlay/python/InDetOverlayConfigDb.py

@@ -5,6 +5,9 @@ from AthenaCommon.CfgGetter import addAlgorithm
 addAlgorithm("InDetOverlay.InDetOverlayConfig.getPixelOverlay", "PixelOverlay")
 addAlgorithm("InDetOverlay.InDetOverlayConfig.getPixelSDOOverlay", "PixelSDOOverlay")
 
+addAlgorithm("InDetOverlay.InDetOverlayConfig.getBCMOverlay", "BCMOverlay")
+addAlgorithm("InDetOverlay.InDetOverlayConfig.getBCMSDOOverlay", "BCMSDOOverlay")
+
 addAlgorithm("InDetOverlay.InDetOverlayConfig.getSCTOverlay", "SCTOverlay")
 addAlgorithm("InDetOverlay.InDetOverlayConfig.getSCTSDOOverlay", "SCTSDOOverlay")
 

InnerDetector/InDetRawAlgs/InDetOverlay/src/BCMOverlay.cxx

+/*
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "InDetOverlay/BCMOverlay.h"
+
+#include "StoreGate/ReadHandle.h"
+#include "StoreGate/WriteHandle.h"
+
+
+BCMOverlay::BCMOverlay(const std::string &name, ISvcLocator *pSvcLocator)
+  : AthReentrantAlgorithm(name, pSvcLocator)
+{
+}
+
+StatusCode BCMOverlay::initialize()
+{
+  ATH_MSG_DEBUG("Initializing...");
+
+  // Check and initialize keys
+  ATH_CHECK( m_bkgInputKey.initialize(!m_bkgInputKey.empty()) );
+  ATH_MSG_VERBOSE("Initialized ReadHandleKey: " << m_bkgInputKey);
+  ATH_CHECK( m_signalInputKey.initialize() );
+  ATH_MSG_VERBOSE("Initialized ReadHandleKey: " << m_signalInputKey);
+  ATH_CHECK( m_outputKey.initialize() );
+  ATH_MSG_VERBOSE("Initialized WriteHandleKey: " << m_outputKey);
+
+  return StatusCode::SUCCESS;
+}
+
+StatusCode BCMOverlay::execute(const EventContext& ctx) const
+{
+  ATH_MSG_DEBUG("execute() begin");
+
+  // Reading the input RDOs
+  ATH_MSG_VERBOSE("Retrieving input RDO containers");
+
+  const BCM_RDO_Container *bkgContainerPtr = nullptr;
+  if (!m_bkgInputKey.empty()) {
+    SG::ReadHandle<BCM_RDO_Container> bkgContainer(m_bkgInputKey, ctx);
+    if (!bkgContainer.isValid()) {
+      ATH_MSG_ERROR("Could not get background BCM RDO container " << bkgContainer.name() << " from store " << bkgContainer.store());
+      return StatusCode::FAILURE;
+    }
+    bkgContainerPtr = bkgContainer.cptr();
+
+    ATH_MSG_DEBUG("Found background BCM RDO container " << bkgContainer.name() << " in store " << bkgContainer.store());
+  }
+
+  SG::ReadHandle<BCM_RDO_Container> signalContainer(m_signalInputKey, ctx);
+  if (!signalContainer.isValid()) {
+    ATH_MSG_ERROR("Could not get signal BCM RDO container " << signalContainer.name() << " from store " << signalContainer.store());
+    return StatusCode::FAILURE;
+  }
+  ATH_MSG_DEBUG("Found signal BCM RDO container " << signalContainer.name() << " in store " << signalContainer.store());
+
+  // Creating output RDO container
+  SG::WriteHandle<BCM_RDO_Container> outputContainer(m_outputKey, ctx);
+  ATH_CHECK(outputContainer.record(std::make_unique<BCM_RDO_Container>()));
+  if (!outputContainer.isValid()) {
+    ATH_MSG_ERROR("Could not record output BCM RDO container " << outputContainer.name() << " to store " << outputContainer.store());
+    return StatusCode::FAILURE;
+  }
+  ATH_MSG_DEBUG("Recorded output BCM RDO container " << outputContainer.name() << " in store " << outputContainer.store());
+
+  ATH_CHECK(overlayContainer(bkgContainerPtr, signalContainer.cptr(), outputContainer.ptr()));
+
+  ATH_MSG_DEBUG("execute() end");
+  return StatusCode::SUCCESS;
+}
+
+StatusCode BCMOverlay::overlayContainer(const BCM_RDO_Container *bkgContainer, 
+                                        const BCM_RDO_Container *signalContainer, 
+                                        BCM_RDO_Container *outputContainer) const
+{
+
+  if (!bkgContainer) {
+    for (const BCM_RDO_Collection *coll : *signalContainer) {
+
+      std::unique_ptr<BCM_RDO_Collection> outputColl = std::make_unique<BCM_RDO_Collection>();
+      
+      for (const BCM_RawData *rdo : *coll) {
+        outputColl->push_back(new BCM_RawData(rdo->getWord1(), rdo->getWord2()));
+      }
+      outputContainer->push_back(outputColl.release());
+    }
+
+    return StatusCode::SUCCESS;
+  }
+
+  size_t containerSize = signalContainer->size();
+
+  for (size_t i = 0; i < containerSize; i++) {
+
+    const BCM_RDO_Collection *sigColl = signalContainer->at(i);
+    const BCM_RDO_Collection *bkgColl = bkgContainer->at(i);
+
+    std::unique_ptr<BCM_RDO_Collection> outputColl = std::make_unique<BCM_RDO_Collection>();
+    size_t collectionSize = sigColl->size();
+    if (collectionSize != bkgColl->size()) {
+      ATH_MSG_ERROR ("BCM signal and background collection size mismatch");
+      return StatusCode::FAILURE;
+    }
+
+    if (bkgColl->getChannel() == sigColl->getChannel()) {
+      outputColl->setChannel(sigColl->getChannel());
+    } else {
+      ATH_MSG_ERROR ("BCM signal and background channel mismatch");
+      return StatusCode::FAILURE;
+    }
+
+    for (size_t j = 0; j < collectionSize; j++) {
+
+      if (bkgColl->at(j)->getChannel() == sigColl->at(j)->getChannel()) {
+        BCM_RawData *mergedRDO = mergeChannel(bkgColl->at(j),sigColl->at(j));
+        if (mergedRDO) outputColl->push_back(mergedRDO);
+        else {
+          ATH_MSG_ERROR ("BCM channel merging failed");
+          return StatusCode::FAILURE;
+        }
+      }
+      else {
+        ATH_MSG_ERROR ("BCM signal and background channel mismatch");
+        return StatusCode::FAILURE;
+      }
+
+    }
+    outputContainer->push_back(outputColl.release());
+  }
+
+  return StatusCode::SUCCESS;
+}
+
+BCM_RawData *BCMOverlay::mergeChannel(const BCM_RawData *bkgRDO, 
+                                      const BCM_RawData *signalRDO) const
+{
+
+  if (bkgRDO->getPulse1Width()==0) {
+    return new BCM_RawData(signalRDO->getWord1(), signalRDO->getWord2());
+  } else if (signalRDO->getPulse1Width()==0) {
+    return new BCM_RawData(bkgRDO->getWord1(), bkgRDO->getWord2());
+  }
+
+  unsigned int bkg_p1 = bkgRDO->getPulse1Position();
+  unsigned int bkg_w1 = bkgRDO->getPulse1Width();
+  unsigned int bkg_p2 = bkgRDO->getPulse2Position();
+  unsigned int bkg_w2 = bkgRDO->getPulse2Width();
+  unsigned int sig_p1 = signalRDO->getPulse1Position();
+  unsigned int sig_w1 = signalRDO->getPulse1Width();
+  unsigned int sig_p2 = signalRDO->getPulse2Position();
+  unsigned int sig_w2 = signalRDO->getPulse2Width();
+
+  std::vector<std::unique_ptr<BCM_Pulse>> merged_pulses;
+
+  if (bkg_w1 > 0) merged_pulses.push_back(std::make_unique<BCM_Pulse>(bkg_p1,bkg_w1));
+  if (bkg_w2 > 0) merged_pulses.push_back(std::make_unique<BCM_Pulse>(bkg_p2,bkg_w2));
+  if (sig_w1 > 0) merged_pulses.push_back(std::make_unique<BCM_Pulse>(sig_p1,sig_w1));
+  if (sig_w2 > 0) merged_pulses.push_back(std::make_unique<BCM_Pulse>(sig_p2,sig_w2));
+
+  overlayPulses(merged_pulses);
+  std::sort(merged_pulses.begin(), merged_pulses.end(), compare);
+
+  // Check whether some of the pulses merged
+  for (size_t i = 0; i < merged_pulses.size()-1; i++) {
+
+    BCM_Pulse *early = merged_pulses.at(i).get();
+    BCM_Pulse *later = merged_pulses.at(i+1).get();
+
+    if ( (early->p + early->w - 1) >= later->p ) {
+      early->w = later->p - early->p + later->w;      // Enlarge early pulse
+      merged_pulses.erase(merged_pulses.begin()+i+1, 
+                          merged_pulses.begin()+i+2); // Omit later pulse
+      i--;
+    }
+  }
+
+  unsigned int merged_p1;
+  unsigned int merged_w1;
+  unsigned int merged_p2;
+  unsigned int merged_w2;
+
+  if (merged_pulses.size() > 0) {
+    merged_p1 = merged_pulses.at(0)->p;
+    merged_w1 = merged_pulses.at(0)->w;
+  } else {
+    merged_p1 = 0;
+    merged_w1 = 0;
+  }
+
+  if (merged_pulses.size() > 1) {
+    merged_p2 = merged_pulses.at(1)->p;
+    merged_w2 = merged_pulses.at(1)->w;
+  } else {
+    merged_p2 = 0;
+    merged_w2 = 0;
+  }
+
+  // Record two earliest pulses into the output RDO
+  return new BCM_RawData(signalRDO->getChannel(),
+                         merged_p1, merged_w1,
+                         merged_p2, merged_w2,
+                         signalRDO->getLVL1A(),
+                         signalRDO->getBCID(),
+                         signalRDO->getLVL1ID());
+}
+
+void BCMOverlay::overlayPulses(std::vector<std::unique_ptr<BCM_Pulse>>& merged_pulses) const
+{
+
+  constexpr double fullPulseWidth{15.};      // Analogue pulse width
+                                             // before application of
+                                             // time over threshold
+  constexpr double slopeUpFraction{1./3.};   // Pulse rise time,
+                                             // expressed in the fraction
+                                             // of full pulse width
+  constexpr double slopeDownFraction{2./3.}; // Pulse fall time,
+                                             // expressed in the fraction
+                                             // of full pulse width
+
+  for (size_t i = 0; i < merged_pulses.size(); i++) {
+    BCM_Pulse* pulse_1 = merged_pulses.at(i).get();
+
+    for (size_t j = 0; j < i; j++) {
+      BCM_Pulse* pulse_2 = merged_pulses.at(j).get();
+      auto[early,later] = timeOrder(pulse_1, pulse_2);
+
+      double slope_up = 1./slopeUpFraction/(fullPulseWidth - later->w);
+      double slope_down = 1./slopeDownFraction/(fullPulseWidth - early->w);
+      int bin_min = early->p + early->w;
+      int bin_max = later->p;
+
+      // Widen the pulses, if they lie close enough to each other
+      for (int bin_iter=bin_min; bin_iter < bin_max; bin_iter++) {
+        if (slope_up*(bin_iter - bin_max) - slope_down*(bin_iter - bin_min) > -1) {
+          early->w++;
+        }
+        else break;
+      }
+      for (int bin_iter=bin_max-1; bin_iter >= bin_min; bin_iter--) {
+        if (slope_up*(bin_iter - bin_max) - slope_down*(bin_iter - bin_min) > -1) {
+          later->p = bin_iter;
+          later->w++;
+        }
+        else break;
+      }
+
+    }
+  }
+
+}
+
+std::pair<BCM_Pulse*, BCM_Pulse*> BCMOverlay::timeOrder(BCM_Pulse* pulse1, 
+                                                        BCM_Pulse* pulse2) const
+{
+
+  if (pulse2->p > pulse1->p) return std::pair(pulse1,pulse2);
+  else return std::pair(pulse2,pulse1);
+
+}
+
+bool BCMOverlay::compare(const std::unique_ptr<BCM_Pulse>& a, 
+                         const std::unique_ptr<BCM_Pulse>& b)
+{
+
+  return a->p < b->p;
+  
+}
\ No newline at end of file

InnerDetector/InDetRawAlgs/InDetOverlay/src/components/InDetOverlay_entries.cxx

+#include "InDetOverlay/BCMOverlay.h"
 #include "InDetOverlay/PixelOverlay.h"
 #include "InDetOverlay/SCTOverlay.h"
 #include "InDetOverlay/TRTOverlay.h"
@@ -6,6 +7,7 @@
 
 DECLARE_COMPONENT( InDetSDOOverlay )
   
+DECLARE_COMPONENT( BCMOverlay )
 DECLARE_COMPONENT( PixelOverlay )
 DECLARE_COMPONENT( SCTOverlay )
 DECLARE_COMPONENT( TRTOverlay )

Tools/PROCTools/python/RunTier0TestsTools.py

@@ -28,7 +28,7 @@ ciRefFileMap = {
                 's3505-22.0'           : 'v3',
                 # OverlayTier0Test_required-test
                 'overlay-d1498-21.0'   : 'v2',
-                'overlay-d1498-22.0'   : 'v30',
+                'overlay-d1498-22.0'   : 'v31',
                 'overlay-bkg-21.0'     : 'v1',
                 'overlay-bkg-22.0'     : 'v4',
                }