Merge branch 'hgtd_fix_digi_positions_21p9' into '21.9'

configs, bugfixes and patching of rdo identifiers to have reco running with full digi input

See merge request atlas/athena!46863

HighGranularityTimingDetector/HGTD_DetDescr/HGTD_ReadoutGeometry/src/HGTD_DetectorElement.cxx

@@ -32,14 +32,14 @@ HGTD_DetectorElement::~HGTD_DetectorElement()
     delete m_surface;
 }
 
-void 
+void
 HGTD_DetectorElement::updateCache() const
 {
     SolidStateDetectorElementBase::updateCache();
     m_cacheValid = true;
 }
 
-Identifier 
+Identifier
 HGTD_DetectorElement::identifierFromCellId(const SiCellId & cellId) const
 {
     Identifier id; // Will be initialized in an invalid state.
@@ -65,11 +65,11 @@ HGTD_DetectorElement::cellIdFromIdentifier(const Identifier & identifier) const
             cellId = SiCellId(hgtdIdHelper->phi_index(identifier), hgtdIdHelper->eta_index(identifier));
         }
     }
-  
+
     return cellId;
 }
 
-const std::vector<const Trk::Surface*>& HGTD_DetectorElement::surfaces() const 
+const std::vector<const Trk::Surface*>& HGTD_DetectorElement::surfaces() const
 {
     // This method is needed to satisfy inheritance from TrkDetElementBase
     // so just return the one surface

HighGranularityTimingDetector/HGTD_Digitization/HGTD_Digitization/CMakeLists.txt

@@ -20,3 +20,5 @@ atlas_add_component( HGTD_Digitization
                   InDetSimData InDetSimEvent AthenaKernel AtlasCLHEP_RandomGenerators)
 
 atlas_install_headers( HGTD_Digitization )
+
+atlas_install_python_modules( python/*.py)

HighGranularityTimingDetector/HGTD_Digitization/HGTD_Digitization/HGTD_Digitization/HGTD_DigitizationTool.h

@@ -111,6 +111,7 @@ private:
   const HGTD_DetectorManager* m_det_mgr;
 
   HGTD::HGTD_RDOContainer* m_hgtd_rdo_container;
+  std::string m_output_rdo_cont_name;
 
   InDetSimDataCollection* m_sdo_collection;
   std::string m_output_sdo_coll_name;

HighGranularityTimingDetector/HGTD_Digitization/HGTD_Digitization/python/HGTD_DigitizationConfigDb.py

+# Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
+
+from AthenaCommon.CfgGetter import addTool,addService,addAlgorithm
+addTool("HGTD_Digitization.HGTD_DigitizationConfig.getHGTD_DigitizationTool", "HGTD_DigitizationTool")
+## FIXME HGTD_Amp needs to be added at later point
+# addTool("HGTD_Digitization.HGTD_DigitizationConfig.getHGTD_Amp", "HGTD_Amp" )
+addTool("HGTD_Digitization.HGTD_DigitizationConfig.getHGTD_FrontEndTool", "HGTD_FrontEndTool" )
+addTool("HGTD_Digitization.HGTD_DigitizationConfig.getHGTD_SurfaceChargesGenerator", "HGTD_SurfaceChargesGenerator" )

HighGranularityTimingDetector/HGTD_Digitization/HGTD_Digitization/src/HGTD_DigitizationTool.cxx

@@ -34,6 +34,7 @@ HGTD_DigitizationTool::HGTD_DigitizationTool(const std::string& type,
       m_id_helper(nullptr),
       m_det_mgr(nullptr),
       m_hgtd_rdo_container(nullptr),
+      m_output_rdo_cont_name("HGTD_RDOs"),
       m_sdo_collection(nullptr),
       m_output_sdo_coll_name("HGTD_SDO_Map"),
       m_input_collection_name("HGTD_Hits"),
@@ -48,6 +49,8 @@ HGTD_DigitizationTool::HGTD_DigitizationTool(const std::string& type,
                   "Merge service used in Pixel & HGTD digitization");
   declareProperty("InputCollecionName", m_input_collection_name = "HGTD_Hits",
                   "Name of the collection that is digitized");
+  declareProperty("OutoutRDOContName", m_output_rdo_cont_name = "HGTD_RDOs",
+                  "Name of the retrieved RDO collection");
   declareProperty("ActiveTimeWindow", m_active_time_window = 10000,
                   "Hits within this time window are used for digitization");
   declareProperty("RndmSvc", m_rndm_svc,
@@ -178,7 +181,7 @@ StatusCode HGTD_DigitizationTool::prepareEvent(unsigned int /*index*/) {
       new HGTD::HGTD_RDOContainer(m_id_helper->wafer_hash_max());
 
   // register RDO container with storegate
-  ATH_CHECK(evtStore()->record(m_hgtd_rdo_container, "HGTD_RDOs"));
+  ATH_CHECK(evtStore()->record(m_hgtd_rdo_container, m_output_rdo_cont_name));
 
   // Create a map for the SDO and register it into StoreGate
   m_sdo_collection = new InDetSimDataCollection();
@@ -352,8 +355,14 @@ HGTD_DigitizationTool::createRDOCollection(
     int eta_index = readout_cell.etaIndex();
     int phi_index = readout_cell.phiIndex();
 
+    ATH_MSG_DEBUG("readout_cell ID: " << readout_cell << " eta:" << eta_index
+                                      << " phi:" << phi_index);
+
+    // FIXME the method takes ID, phi, eta
+    // FIXME switching order here to fix upstream issue with orientations
     const Identifier id_readout =
-        m_id_helper->pixel_id(charged_diodes->identify(), eta_index, phi_index);
+        m_id_helper->pixel_id(charged_diodes->identify(), phi_index,
+        eta_index);
 
     SiChargedDiode& diode = (*i_chargedDiode).second;
     const SiTotalCharge& charge = diode.totalCharge();
@@ -449,10 +458,12 @@ void HGTD_DigitizationTool::createAndStoreSDO(
     // particle to which the truth information was kept. Can be HS and PU.
     if (real_particle_hit) {
 
-      auto collection_id = charged_diodes->identify();
-      Identifier id_readout = m_id_helper->pixel_id(
-          collection_id, m_id_helper->eta_index(collection_id),
-          m_id_helper->phi_index(collection_id));
+      InDetDD::SiReadoutCellId readout_cell =
+          (*i_chargedDiode).second.getReadoutCell();
+      int eta_index = readout_cell.etaIndex();
+      int phi_index = readout_cell.phiIndex();
+      const Identifier id_readout = m_id_helper->pixel_id(
+          charged_diodes->identify(), eta_index, phi_index);
 
       m_sdo_collection->insert(std::make_pair(
           id_readout, InDetSimData(deposits, (*i_chargedDiode).second.flag())));

HighGranularityTimingDetector/HGTD_Digitization/HGTD_Digitization/src/HGTD_SurfaceChargesGenerator.cxx

@@ -68,7 +68,7 @@ void HGTD_SurfaceChargesGenerator::createSurfaceChargesFromHit(
     }
   }
 
-  double sensor_thickness = element->design().thickness();
+  float sensor_thickness = element->design().thickness();
   int readout_side = element->design().readoutSide();
 
   float pixel_size_xphi = element->design().phiPitch();
@@ -82,15 +82,15 @@ void HGTD_SurfaceChargesGenerator::createSurfaceChargesFromHit(
                                   << ", xDep=" << start_pos[SiHit::xDep]);
 
   CLHEP::Hep3Vector direction = end_pos - start_pos;
-  double deposit_length = direction.mag();
+  float deposit_length = direction.mag();
   int n_steps = deposit_length / m_small_step_length + 1;
   direction.setMag(deposit_length / static_cast<float>(n_steps));
 
-  double tot_eloss = hit.energyLoss();
+  float tot_eloss = hit.energyLoss();
   // FIXME using the mean ionization energy in Silicon
   const float tot_charge = tot_eloss / (3.62 * CLHEP::eV);
 
-  double charge_per_step = tot_charge / static_cast<float>(n_steps);
+  float charge_per_step = tot_charge / static_cast<float>(n_steps);
 
   // FIXME is this correct? does the eventTime include a "later" truth event and
   // the meanTime is just the TOF?
@@ -106,6 +106,16 @@ void HGTD_SurfaceChargesGenerator::createSurfaceChargesFromHit(
 
   for (int i_step = 0; i_step < n_steps; i_step++) {
     CLHEP::Hep3Vector surface_pos = start_pos + i_step * direction;
+    ATH_MSG_DEBUG("surface_pos x=" << surface_pos.x()
+                                   << ", y=" << surface_pos.y()
+                                   << ", z=" << surface_pos.z());
+    ATH_MSG_DEBUG("surface_pos xEta=" << surface_pos[SiHit::xEta]
+                                      << ", xPhi=" << surface_pos[SiHit::xPhi]
+                                      << ", xDep=" << surface_pos[SiHit::xDep]);
+    // NB! y aka xPhi is the long side of the module!
+    // DEBUG surface_pos x=9.08365,  y=-1.17206, z=-0.025
+    // DEBUG surface_pos xEta=-0.025, xPhi=-1.17206, xDep=9.08365
+    // FIXME: eta, phi andd dep need to be revisited in Rel 22!
 
     // Distance between charge and readout side.  p_design->readoutSide() is
     // +1 if readout side is in +ve depth axis direction and visa-versa.
@@ -117,7 +127,7 @@ void HGTD_SurfaceChargesGenerator::createSurfaceChargesFromHit(
     }
     // diffusion sigma
     // FIXME where is the 0.3 from?
-    double rdif = m_diffusion_constant * std::sqrt(spess / 0.3);
+    float rdif = m_diffusion_constant * std::sqrt(spess / 0.3);
 
     // position at the surface, adding smearing
     // FIXME currently no Lorentz angle considered, can be studied in the future
@@ -163,8 +173,7 @@ void HGTD_SurfaceChargesGenerator::createSurfaceChargesFromHit(
         element->cellIdOfPosition(surface_charge.position());
     ATH_MSG_DEBUG("cell_id x=" << cell_id);
     if (cell_id.isValid()) {
-      // add this charge to the collection (or merge in existing charged
-      // diode)
+      // add this charge to the collection (or merge in existing charged diode)
       diode_coll->add(cell_id, surface_charge.charge());
     }
   } // END LOOP over steps

HighGranularityTimingDetector/HGTD_Reconstruction/HGTD_RecAlgs/HGTD_PRDFormation/CMakeLists.txt

@@ -11,7 +11,7 @@ atlas_depends_on_subdirs( PUBLIC
                           GaudiKernel)
 
 # Component(s) in the package:
-atlas_add_component( TrackTimeExtensionAlg
+atlas_add_component( HGTD_PRDFormation
                      src/*.cxx
                      src/components/*.cxx
                      LINK_LIBRARIES AthenaBaseComps GaudiKernel HGTD_RecToolInterfaces

HighGranularityTimingDetector/HGTD_Reconstruction/HGTD_RecAlgs/HGTD_TrackTimeExtension/CMakeLists.txt

@@ -11,7 +11,7 @@ atlas_depends_on_subdirs( PUBLIC
                           GaudiKernel)
 
 # Component(s) in the package:
-atlas_add_component( TrackTimeExtensionAlgorithm
+atlas_add_component( HGTD_TrackTimeExtension
                      src/*.cxx
                      src/components/*.cxx
                      LINK_LIBRARIES AthenaBaseComps GaudiKernel

HighGranularityTimingDetector/HGTD_Reconstruction/HGTD_RecAlgs/HGTD_TrackTimeExtension/share/postInclude.HGTD_TrkExt.py

+##call clustermakertools before running extension when using RDO input
+
+from HGTD_PadClusterizationTools.HGTD_PadClusterizationToolsConf import HGTD__HGTD_ClusterMakerTool
+
+cluster_maker_tool = HGTD__HGTD_ClusterMakerTool(name="HGTD_ClusterMakerTool")
+cluster_maker_tool.OutputLevel = DEBUG
+ToolSvc += cluster_maker_tool
+
+from HGTD_PadClusterizationTools.HGTD_PadClusterizationToolsConf import HGTD__SinglePadClusterTool
+
+clusterization_tool = HGTD__SinglePadClusterTool(name="SinglePadClusterTool")
+clusterization_tool.ClusterMakerTool = cluster_maker_tool
+clusterization_tool.OutputLevel = DEBUG
+ToolSvc += clusterization_tool
+
+
+##FIXME NAME!! this needs change in the CMAKE
+from HGTD_PRDFormation.PadClusterizationAlgConf import HGTD__PadClusterizationAlg
+
+hgtd_cluster_maker_alg = HGTD__PadClusterizationAlg(name="PadClusterizationAlg")
+hgtd_cluster_maker_alg.OutputLevel = DEBUG
+hgtd_cluster_maker_alg.ClusterizationTool = clusterization_tool
+hgtd_cluster_maker_alg.PRDContainerName = "HGTD_Cluster"
+
+topSequence += hgtd_cluster_maker_alg
+
+##continue with the extension
 
 extrapolator = ToolSvc.AtlasExtrapolator
 
@@ -6,14 +33,14 @@ kalman_updator = ToolSvc.InDetUpdator
 from HGTD_TimeCalibrationTools.HGTD_TimeCalibrationToolsConf import HGTD__StraightLineTOFcorrectionTool
 
 hgtd_tof_corr_tool = HGTD__StraightLineTOFcorrectionTool(name="StraightLineTOFcorrectionTool")
-hgtd_tof_corr_tool.OutputLevel = VERBOSE
+hgtd_tof_corr_tool.OutputLevel = DEBUG
 
 ToolSvc += hgtd_tof_corr_tool
 
 from HGTD_TrackTimeExtensionTools.HGTD_TrackTimeExtensionToolsConf import HGTD__HGTD_IterativeExtensionTool
 
 htgd_ext_tool = HGTD__HGTD_IterativeExtensionTool(name="IterativeExtensionTool")
-htgd_ext_tool.OutputLevel = VERBOSE
+htgd_ext_tool.OutputLevel = DEBUG
 htgd_ext_tool.ExtrapolatorTool = extrapolator
 htgd_ext_tool.UpdatorTool = kalman_updator
 htgd_ext_tool.TOFCorrTool = hgtd_tof_corr_tool
@@ -24,14 +51,17 @@ ToolSvc += htgd_ext_tool
 from HGTD_TruthTools.HGTD_TruthToolsConf import HGTD__ClusterTruthTool
 
 hgtd_cluster_truth_tool = HGTD__ClusterTruthTool(name="ClusterTruthTool")
+hgtd_cluster_truth_tool.OutputLevel = DEBUG
 
 ToolSvc += hgtd_cluster_truth_tool
 
+##FIXME NAME!! this needs change in the CMAKE
 from HGTD_TrackTimeExtension.TrackTimeExtensionAlgConf import HGTD__TrackTimeExtensionAlg
 
 hgtd_alg = HGTD__TrackTimeExtensionAlg(name="TrackTimeExtensionAlg")
-hgtd_alg.OutputLevel = VERBOSE
+hgtd_alg.OutputLevel = DEBUG
 hgtd_alg.TimeExtensionTool = htgd_ext_tool
+hgtd_alg.HGTDClusterContainerName = "HGTD_Cluster"
 hgtd_alg.TruthTool = hgtd_cluster_truth_tool
 
 topSequence += hgtd_alg
@@ -61,7 +91,7 @@ truthList = [ 'xAOD::TruthEventContainer#TruthEvents',
             'PileUpEventInfo#OverlayEvent' ]
 
 # --- create stream
-StreamESD            = AthenaPoolOutputStream ( "StreamESD2", "ESD.HGTD.pool.root",asAlg=True)
+StreamESD            = AthenaPoolOutputStream ( "StreamESD2", "ESD.HGTD.full.pool.root",asAlg=True)
 # --- save MC collections if truth turned on
 StreamESD.ItemList += truthList
 # ---- load list of objects

HighGranularityTimingDetector/HGTD_Reconstruction/HGTD_RecAlgs/HGTD_TrackTimeExtension/src/TrackTimeExtensionAlg.cxx

@@ -196,6 +196,7 @@ StatusCode TrackTimeExtensionAlg::decorateTrackParticle(
        extension) {
 
     if (trk_state) {
+      ATH_MSG_DEBUG("[decorateTrackParticle] extension found");
       has_cluster_vec.emplace_back(true);
 
       chi2_vec.emplace_back(
@@ -236,10 +237,13 @@ StatusCode TrackTimeExtensionAlg::decorateTrackParticle(
         is_merged = truth_info.is_merged;
       }
 
+      ATH_MSG_DEBUG("Truth origin: " << truth_origin);
+
       truth_vec.emplace_back(truth_origin);
       is_shadowed_vec.emplace_back(is_shadowed);
       is_merged_vec.emplace_back(is_merged);
     } else {
+      ATH_MSG_DEBUG("[decorateTrackParticle] NO extension found");
       has_cluster_vec.emplace_back(false);
       if (not skip_deco) {
         chi2_vec.emplace_back(-1.);

HighGranularityTimingDetector/HGTD_Reconstruction/HGTD_RecTools/HGTD_PadClusterizationTools/src/SinglePadClusterTool.cxx

@@ -19,8 +19,10 @@ SinglePadClusterTool::SinglePadClusterTool(const std::string& type,
                                            const std::string& name,
                                            const IInterface* parent)
     : AthAlgTool(type, name, parent),
-      m_cluster_maker(nullptr),
-      m_hgtd_det_mgr(nullptr) {}
+      m_cluster_maker("HGTD::HGTD_ClusterMakerTool/HGTD_ClusterMakerTool"),
+      m_hgtd_det_mgr(nullptr) {
+  declareProperty("ClusterMakerTool", m_cluster_maker);
+}
 
 StatusCode SinglePadClusterTool::initialize() {
 
@@ -71,8 +73,11 @@ SinglePadClusterTool::clusterize(const HGTD_RDOCollection& rdo_coll) const {
     // in col and row direction
     int col_width = 1;
     int row_width = 1;
-    double eta_width = det_design.etaPitch();
-    double phi_width = det_design.phiPitch();
+    //FIXME there is a bug in the detector design! should not be hardcoded here
+    // double eta_width = det_design.etaPitch();
+    // double phi_width = det_design.phiPitch();
+    double eta_width = 1.3;
+    double phi_width = 1.3;
     ATH_MSG_DEBUG("eta_width=" << eta_width << " phi_width=" << phi_width);
 
     InDet::SiWidth si_width(Amg::Vector2D(col_width, row_width),

HighGranularityTimingDetector/HGTD_Reconstruction/HGTD_RecTools/HGTD_PadClusterizationTools/src/components/HGTD_PadClusterizationTools_entries.cxx

 #include "GaudiKernel/DeclareFactoryEntries.h"
 #include "HGTD_PadClusterizationTools/HGTD_ClusterMakerTool.h"
+#include "HGTD_PadClusterizationTools/SinglePadClusterTool.h"
 
 using namespace HGTD;
 
 DECLARE_TOOL_FACTORY(HGTD_ClusterMakerTool)
 
-DECLARE_FACTORY_ENTRIES( HGTD_PadClusterizationTools )
-{
-  DECLARE_TOOL( HGTD_ClusterMakerTool )
+DECLARE_TOOL_FACTORY(SinglePadClusterTool)
+
+DECLARE_FACTORY_ENTRIES(HGTD_PadClusterizationTools) {
+  DECLARE_TOOL(HGTD_ClusterMakerTool)
+  DECLARE_TOOL(SinglePadClusterTool)
 }

HighGranularityTimingDetector/HGTD_Reconstruction/HGTD_RecTools/HGTD_TruthTools/src/ClusterTruthTool.cxx

@@ -30,6 +30,10 @@ HGTD::ClusterTruthInfo HGTD::ClusterTruthTool::classifyCluster(
     const InDetSimDataCollection* sim_data,
     const HepMC::GenEvent* hard_scatter_evnt) const {
 
+  if (not sim_data) {
+    return {ClusterTruthOrigin::UNIDENTIFIED, false, false};
+  }
+
   const std::vector<Identifier>& rdo_id_list = cluster->rdoList();
   // keep record of the cluster origins and if they are shadowed or not
   std::vector<std::pair<HGTD::ClusterTruthOrigin, bool>> shadowed_origins;
@@ -40,10 +44,11 @@ HGTD::ClusterTruthInfo HGTD::ClusterTruthTool::classifyCluster(
     // the InDetSimData contains a std::pair<HepMcParticleLink, float>, where
     // the second entry in the pair holds the time of the SiChargedDiode
     if (pos == sim_data->end()) {
-      ATH_MSG_WARNING(
-          "[ClusterTruthTool::classifyCluster] ID not found in SDO map");
+      ATH_MSG_WARNING("[ClusterTruthTool::classifyCluster] ID not found in SDO "
+                      "map, going to next ID");
       // FIXME I should probably continue here already? otherwise I get an
       // "empty" entry in shadowed_origins
+      continue;
     }
     // collect deposits, sorted with first deposit at start of map bu default
     std::map<float, ClusterTruthOrigin> sorted_deposits;

Simulation/Digitization/python/DigiAlgConfig.py

@@ -173,6 +173,14 @@ def getFastHGTDPileUpTools():
         PileUpToolsList += [ "HGTD_SmearedDigitizationTool" ]
     return PileUpToolsList
 
+def getStandardHGTDPileUpTools():
+    from AthenaCommon.DetFlags import DetFlags
+    from Digitization.DigitizationFlags import digitizationFlags
+    PileUpToolsList = []
+    if DetFlags.digitize.HGTD_on():
+        PileUpToolsList += [ "HGTD_DigitizationTool" ]
+    return PileUpToolsList
+
 def getSplitInDetPileUpTools():
     from AthenaCommon.DetFlags import DetFlags
     PileUpToolsList = []
@@ -315,6 +323,8 @@ def getStandardPileUpToolsList():
     PileUpToolsList += getStandardForwardPileUpTools()
     ## Inner Detector Digitization
     PileUpToolsList += getStandardInDetPileUpTools()
+    ## HGTD Digitization
+    PileUpToolsList += getStandardHGTDPileUpTools()
     ## Calo Digitization
     PileUpToolsList += getStandardCaloPileUpTools()
     ## Muon System Digitization
@@ -331,6 +341,8 @@ def getStandardSignalOnlyTruthPileUpToolsList():
     PileUpToolsList += getStandardForwardPileUpTools()
     ## Inner Detector Digitization
     PileUpToolsList += getStandardInDetPileUpTools()
+    ## HGTD Digitization
+    PileUpToolsList += getStandardHGTDPileUpTools()
     ## Calo Digitization
     PileUpToolsList += getStandardCaloPileUpTools()
     ## Muon System Digitization
@@ -348,7 +360,7 @@ def getStandardInTimeOnlyTruthPileUpToolsList():
     ## Inner Detector Digitization
     PileUpToolsList += getStandardInDetPileUpTools()
     ## HGTD Digitization
-    PileUpToolsList += getFastHGTDPileUpTools()
+    PileUpToolsList += getStandardHGTDPileUpTools()
     ## Calo Digitization
     PileUpToolsList += getStandardCaloPileUpTools()
     ## Muon System Digitization

Simulation/Digitization/python/DigiOutput.py

@@ -133,6 +133,9 @@ def getStreamRDO_ItemList(log):
     if DetFlags.digitize.HGTD_on():
         if 'doHGTD_FastDigi' in digitizationFlags.experimentalDigi():
             StreamRDO_ItemList+=["HGTD::HGTD_ClusterContainer#*"]
+        else:
+            StreamRDO_ItemList+=["HGTD::HGTD_RDOContainer#*"]
+            
     # Calorimeter Output
     if DetFlags.digitize.LAr_on():
         if DetFlags.writeRDOPool.LAr_on():