diff --git a/Tracking/TrkConfig/python/AtlasTrackingGeometrySvcConfig.py b/Tracking/TrkConfig/python/AtlasTrackingGeometrySvcConfig.py
index 39d77143c6cfacbd8a10fe501be656d9a56527b1..94012d79769cf23ed1a285f79e215ec956e6ab8a 100644
--- a/Tracking/TrkConfig/python/AtlasTrackingGeometrySvcConfig.py
+++ b/Tracking/TrkConfig/python/AtlasTrackingGeometrySvcConfig.py
@@ -24,11 +24,9 @@ def _setupCondDB(flags, CoolDataBaseFolder, quiet=True):
if flags.Detector.GeometryITk:
AtlasMaterialTag = flags.ITk.trackingGeometry.materialTag+str(flags.ITk.trackingGeometry.version)+'_'
cfolder = CoolDataBaseFolder +'<tag>TagInfoMajor/'+AtlasMaterialTag+'/GeoAtlas</tag>'
-
-
+
if flags.Detector.GeometryITk and flags.ITk.trackingGeometry.loadLocalDbForMaterialMaps:
DataBaseName=flags.ITk.trackingGeometry.localDatabaseName
-
from IOVDbSvc.IOVDbSvcConfig import addFolders
result.merge(addFolders(flags,"/GLOBAL/TrackingGeo/LayerMaterialITK",detDb=DataBaseName, tag=AtlasMaterialTag))
cfolder = CoolDataBaseFolder +'<tag>TagInfoMajor/'+AtlasMaterialTag+'</tag>'
diff --git a/Tracking/TrkDetDescr/TrkDetDescrAlgs/CMakeLists.txt b/Tracking/TrkDetDescr/TrkDetDescrAlgs/CMakeLists.txt
index 05f41201d75f58181ba70bd78ae710b3f04276dc..9dffe2e8ce16673dc68f813682e3ad803d542386 100644
--- a/Tracking/TrkDetDescr/TrkDetDescrAlgs/CMakeLists.txt
+++ b/Tracking/TrkDetDescr/TrkDetDescrAlgs/CMakeLists.txt
@@ -15,4 +15,4 @@ atlas_add_component( TrkDetDescrAlgs
# Install files from the package:
atlas_install_headers( TrkDetDescrAlgs )
-
+atlas_install_python_modules( python/*.py POST_BUILD_CMD ${ATLAS_FLAKE8} )
diff --git a/Tracking/TrkDetDescr/TrkDetDescrAlgs/TrkDetDescrAlgs/MaterialMapping.h b/Tracking/TrkDetDescr/TrkDetDescrAlgs/TrkDetDescrAlgs/MaterialMapping.h
index 8695db946d8efb9b2c40fd85b330a6b87cb05d36..ab803dd50ca77e2421df88f4b9050a383fbd5b02 100644
--- a/Tracking/TrkDetDescr/TrkDetDescrAlgs/TrkDetDescrAlgs/MaterialMapping.h
+++ b/Tracking/TrkDetDescr/TrkDetDescrAlgs/TrkDetDescrAlgs/MaterialMapping.h
@@ -1,13 +1,14 @@
-/*
- Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
-*/
-
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
//////////////////////////////////////////////////////////////////
// MaterialMapping.h, (c) ATLAS Detector software
///////////////////////////////////////////////////////////////////
#ifndef TRKDETDESCRALGS_MATERIALMAPPING_H
#define TRKDETDESCRALGS_MATERIALMAPPING_H
+#define LEGACY_TRKGEOM
// Athena & Gaudi includes
#include "AthenaBaseComps/AthAlgorithm.h"
@@ -19,8 +20,19 @@
#include <string>
#include <map>
#include "GeoPrimitives/GeoPrimitives.h"
-#include "StoreGate/ReadHandleKey.h"
-#include "TrkGeometry/MaterialStepCollection.h"
+#include "StoreGate/ReadHandleKey.h"
+#include "TrkGeometry/MaterialStepCollection.h"
+//TrkDetDescr Algs, Interfaces, Utils
+#include "TrkDetDescrInterfaces/IMaterialMapper.h"
+#include "TrkDetDescrInterfaces/ILayerMaterialCreator.h"
+#include "TrkDetDescrInterfaces/ILayerMaterialAnalyser.h"
+// TrkExtrapolation
+#include "TrkExInterfaces/IExtrapolationEngine.h"
+// TrkGeometry
+#include "TrkGeometry/TrackingGeometry.h"
+#ifdef LEGACY_TRKGEOM
+#include "TrkDetDescrInterfaces/ITrackingGeometrySvc.h"
+#endif
#ifdef TRKDETDESCR_MEMUSAGE
#include "TrkDetDescrUtils/MemoryLogger.h"
@@ -43,7 +55,6 @@ namespace Trk {
class Layer;
class TrackingVolume;
- class TrackingGeometry;
class SurfaceMaterialRecord;
class LayerMaterialRecord;
class LayerMaterialMap;
@@ -52,11 +63,6 @@ namespace Trk {
class BinnedLayerMaterial;
class CompressedLayerMaterial;
class ElementTable;
- class IExtrapolationEngine;
- class IMaterialMapper;
- class ILayerMaterialAnalyser;
- class ILayerMaterialCreator;
- class ITrackingGeometrySvc;
/** @class MaterialMapping
@@ -85,7 +91,7 @@ namespace Trk {
/** standard Athena-Algorithm method */
StatusCode finalize();
-
+
private:
/** Associate the Step to the Layer */
@@ -104,15 +110,13 @@ namespace Trk {
//!< create the LayerMaterialRecord */
void insertLayerMaterialRecord(const Trk::Layer& lay);
- /** Retrieve the TrackingGeometry */
- StatusCode retrieveTrackingGeometry();
-
- /** Tracking Geometry */
- ServiceHandle<Trk::ITrackingGeometrySvc> m_trackingGeometrySvc; //!< Name of the TrackingGeometrySvc
- mutable const TrackingGeometry* m_trackingGeometry; //!< The underlying TrackingGeometry
+ /** Retrieve the TrackingGeometry and its informations */
+ StatusCode handleTrackingGeometry();
+
+ const TrackingGeometry& trackingGeometry() const;
bool m_checkForEmptyHits; //!< use extrapoaltion engine to check for empty hits
- ToolHandle<IExtrapolationEngine> m_extrapolationEngine; //!< cross-check for empty hit scaling
+ ToolHandle<IExtrapolationEngine> m_extrapolationEngine {this, "ExtrapolationEngine", "", "Extrapolation Engine"};
std::string m_mappingVolumeName;
const Trk::TrackingVolume* m_mappingVolume;
@@ -127,13 +131,13 @@ namespace Trk {
bool m_useLayerThickness; //!< use the actual layer thickness
int m_associationType;
- ToolHandle<ILayerMaterialAnalyser> m_layerMaterialRecordAnalyser; //!< the layer material analyser for the layer material record
- ToolHandleArray<ILayerMaterialAnalyser> m_layerMaterialAnalysers; //!< the layer material analysers per creator (if wanted)
- ToolHandleArray<ILayerMaterialCreator> m_layerMaterialCreators; //!< the layer material creators
+ ToolHandle<ILayerMaterialAnalyser> m_layerMaterialRecordAnalyser {this, "LayerMaterialRecordAnalyser", "", "Layer material analyser for the layer material record"};
+ ToolHandleArray<ILayerMaterialAnalyser> m_layerMaterialAnalysers {this, "LayerMaterialAnalysers", {}, "Layer material analysers per creator (if wanted)"};
+ ToolHandleArray<ILayerMaterialCreator> m_layerMaterialCreators {this, "LayerMaterialCreators", {}, "Layer material creators"};
/** Mapper and Inspector */
bool m_mapMaterial;
- ToolHandle<IMaterialMapper> m_materialMapper; //!< Pointer to an IMaterialMapper algTool
+ ToolHandle<IMaterialMapper> m_materialMapper {this, "MaterialMapper", "" , "IMaterialMapper algTool"};
bool m_mapComposition; //!< map the composition of the material
double m_minCompositionFraction; //!< minimal fraction to be accounted for the composition recording
@@ -154,12 +158,44 @@ namespace Trk {
int m_layerMaterialScreenOutput;
+
+#ifdef LEGACY_TRKGEOM
+ ServiceHandle<ITrackingGeometrySvc> m_trackingGeometrySvc {this, "TrackingGeometrySvc", "",""};
+#endif
+ void throwFailedToGetTrackingGeometry() const;
+ const TrackingGeometry* retrieveTrackingGeometry(const EventContext& ctx) const {
+#ifdef LEGACY_TRKGEOM
+ if (m_trackingGeometryReadKey.key().empty()) {
+ return m_trackingGeometrySvc->trackingGeometry();
+ }
+#endif
+ SG::ReadCondHandle<TrackingGeometry> handle(m_trackingGeometryReadKey,ctx);
+ if (!handle.isValid()) {
+ ATH_MSG_FATAL("Could not load TrackingGeometry with name '" << m_trackingGeometryReadKey.key() << "'. Aborting." );
+ throwFailedToGetTrackingGeometry();
+ }
+ return handle.cptr();
+ }
+
+ SG::ReadCondHandleKey<TrackingGeometry> m_trackingGeometryReadKey
+ {this, "TrackingGeometryReadKey", "", "Key of the TrackingGeometry conditions data."};
+
+
#ifdef TRKDETDESCR_MEMUSAGE
MemoryLogger m_memoryLogger; //!< in case the memory is logged
#endif
};
+
+ inline const Trk::TrackingGeometry& Trk::MaterialMapping::trackingGeometry() const {
+ const Trk::TrackingGeometry *tracking_geometry = retrieveTrackingGeometry(Gaudi::Hive::currentContext());
+ if (!tracking_geometry){
+ ATH_MSG_FATAL("Did not get valid TrackingGeometry. Aborting." );
+ throw GaudiException("MaterialMapping", "Problem with TrackingGeometry loading.", StatusCode::FAILURE);
+ }
+ return *tracking_geometry;
+ }
}
#endif
diff --git a/Tracking/TrkDetDescr/TrkDetDescrAlgs/TrkDetDescrAlgs/MaterialValidation.h b/Tracking/TrkDetDescr/TrkDetDescrAlgs/TrkDetDescrAlgs/MaterialValidation.h
index 35ce8535d7d858827488ed0909be587733942d92..0da48192283e0d2feadaf4098b7f4d584d54b931 100755
--- a/Tracking/TrkDetDescr/TrkDetDescrAlgs/TrkDetDescrAlgs/MaterialValidation.h
+++ b/Tracking/TrkDetDescr/TrkDetDescrAlgs/TrkDetDescrAlgs/MaterialValidation.h
@@ -8,6 +8,7 @@
#ifndef TRKDETDESCRALGS_MATERIALVALIDATION_H
#define TRKDETDESCRALGS_MATERIALVALIDATION_H
+#define LEGACY_TRKGEOM
// Gaudi includes
#include "AthenaBaseComps/AthAlgorithm.h"
@@ -16,12 +17,15 @@
#include "GaudiKernel/ToolHandle.h"
//Eigen
#include "GeoPrimitives/GeoPrimitives.h"
+// TrkGeometry
+#include "TrkGeometry/TrackingGeometry.h"
+#ifdef LEGACY_TRKGEOM
+#include "TrkDetDescrInterfaces/ITrackingGeometrySvc.h"
+#endif
namespace Trk {
- class ITrackingGeometrySvc;
class IMaterialMapper;
- class TrackingGeometry;
class TrackingVolume;
class Surface;
@@ -72,16 +76,32 @@ namespace Trk {
const Amg::Vector3D& position,
const Amg::Vector3D& direction);
- /** Retrieve the TrackingGeometry */
- StatusCode retrieveTrackingGeometry();
-
- /** Tracking Geometry */
- ServiceHandle<Trk::ITrackingGeometrySvc> m_trackingGeometrySvc; //!< Name of the TrackingGeometrySvc
- mutable const TrackingGeometry* m_trackingGeometry; //!< The underlying TrackingGeometry
+ const TrackingGeometry& trackingGeometry() const;
+
+ #ifdef LEGACY_TRKGEOM
+ ServiceHandle<ITrackingGeometrySvc> m_trackingGeometrySvc {this, "TrackingGeometrySvc", "",""};
+#endif
+ void throwFailedToGetTrackingGeometry() const;
+ const TrackingGeometry* retrieveTrackingGeometry(const EventContext& ctx) const {
+#ifdef LEGACY_TRKGEOM
+ if (m_trackingGeometryReadKey.key().empty()) {
+ return m_trackingGeometrySvc->trackingGeometry();
+ }
+#endif
+ SG::ReadCondHandle<TrackingGeometry> handle(m_trackingGeometryReadKey,ctx);
+ if (!handle.isValid()) {
+ ATH_MSG_FATAL("Could not load TrackingGeometry with name '" << m_trackingGeometryReadKey.key() << "'. Aborting." );
+ throwFailedToGetTrackingGeometry();
+ }
+ return handle.cptr();
+ }
+
+ SG::ReadCondHandleKey<TrackingGeometry> m_trackingGeometryReadKey
+ {this, "TrackingGeometryReadKey", "", "Key of the TrackingGeometry conditions data."};
/** Mapper and Inspector */
ToolHandle<IMaterialMapper> m_materialMapper; //!< Pointer to an IMaterialMapper algTool
- int m_maxMaterialMappingEvents; //!< limit the number of validation records to avoid 2G files
+ int m_maxMaterialValidationEvents; //!< limit the number of validation records to avoid 2G files
Rndm::Numbers* m_flatDist; //!< Random generator for flat distribution
@@ -91,9 +111,16 @@ namespace Trk {
double m_accTinX0; //!< accumulated t in X0
-
-
};
+
+ inline const Trk::TrackingGeometry& Trk::MaterialValidation::trackingGeometry() const {
+ const Trk::TrackingGeometry *tracking_geometry = retrieveTrackingGeometry(Gaudi::Hive::currentContext());
+ if (!tracking_geometry){
+ ATH_MSG_FATAL("Did not get valid TrackingGeometry. Aborting." );
+ throw GaudiException("MaterialValidation", "Problem with TrackingGeometry loading.", StatusCode::FAILURE);
+ }
+ return *tracking_geometry;
+ }
}
#endif
diff --git a/Tracking/TrkDetDescr/TrkDetDescrAlgs/python/TrkDetDescrAlgsConfig.py b/Tracking/TrkDetDescr/TrkDetDescrAlgs/python/TrkDetDescrAlgsConfig.py
new file mode 100644
index 0000000000000000000000000000000000000000..dfb7f27d9fc334d8c859717d8f68a5af1887de48
--- /dev/null
+++ b/Tracking/TrkDetDescr/TrkDetDescrAlgs/python/TrkDetDescrAlgsConfig.py
@@ -0,0 +1,142 @@
+# Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
+
+"""Define methods to configure TrkDetDescrAlgs"""
+
+from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator
+from AthenaConfiguration.ComponentFactory import CompFactory
+
+
+def ITkMaterialMappingCfg(flags, name="ITkMaterialMapping", **kwargs):
+ """Return configured ComponentAccumulator and tool for MaterialMapping"""
+
+ result=ComponentAccumulator()
+ LayerMaterialName = 'LayerMaterialITK'
+ LayerMaterialDirectory = '/GLOBAL/TrackingGeo/'
+
+ from PixelGeoModelXml.ITkPixelGeoModelConfig import ITkPixelGeometryCfg
+ itkPixel = ITkPixelGeometryCfg(flags)
+ result.merge(itkPixel)
+
+ from StripGeoModelXml.ITkStripGeoModelConfig import ITkStripGeometryCfg
+ itkStrip = ITkStripGeometryCfg(flags)
+ result.merge(itkStrip)
+
+ from BeamPipeGeoModel.BeamPipeGMConfig import BeamPipeGeometryCfg
+ result.merge(BeamPipeGeometryCfg(flags))
+
+ # get the correct TrackingGeometry setup
+ geom_svc=None
+ geom_cond_key=''
+ from InDetRecExample.TrackingCommon import use_tracking_geometry_cond_alg
+ if not use_tracking_geometry_cond_alg :
+ from TrkConfig.AtlasTrackingGeometrySvcConfig import TrackingGeometrySvcCfg
+ acc = TrackingGeometrySvcCfg(flags)
+ geom_svc = acc.getPrimary()
+ result.merge(acc)
+ kwargs.setdefault("TrackingGeometrySvc", geom_svc)
+ else :
+ from TrackingGeometryCondAlg.AtlasTrackingGeometryCondAlgConfig import TrackingGeometryCondAlgCfg
+ result.merge( TrackingGeometryCondAlgCfg(flags) )
+ geom_cond_key = 'AtlasTrackingGeometry'
+ kwargs.setdefault("TrackingGeometryReadKey", geom_cond_key)
+
+ if 'MappingVolumeName' not in kwargs :
+ kwargs.setdefault("MappingVolumeName", 'InDet::Containers::InnerDetector')
+
+ if 'ExtrapolationEngine' not in kwargs :
+ from TrkConfig.AtlasExtrapolationEngineConfig import AtlasExtrapolationEngineCfg
+ kwargs.setdefault("ExtrapolationEngine", result.getPrimaryAndMerge(AtlasExtrapolationEngineCfg(flags)))
+
+ if 'MaterialMapper' not in kwargs :
+ MaterialMapper = CompFactory.Trk.MaterialMapper("MaterialMapper")
+ kwargs.setdefault("MaterialMapper", MaterialMapper)
+
+ if 'LayerMaterialRecordAnalyser' not in kwargs :
+ LayerMaterialAnalyser = CompFactory.Trk.LayerMaterialAnalyser(name="LayerMaterialRecordAnalyser")
+ kwargs.setdefault("LayerMaterialRecordAnalyser", LayerMaterialAnalyser)
+
+ if 'LayerMaterialCreators' not in kwargs :
+ BinnedLayerMaterialCreator = CompFactory.Trk.BinnedLayerMaterialCreator(name="BinnedLayerMaterialCreator")
+ BinnedLayerMaterialCreator.LayerMaterialName = LayerMaterialName
+ BinnedLayerMaterialCreator.LayerMaterialDirectory = LayerMaterialDirectory
+ LayerMaterialCreators = [BinnedLayerMaterialCreator]
+ kwargs.setdefault("LayerMaterialCreators", LayerMaterialCreators)
+
+ if 'LayerMaterialAnalysers' not in kwargs :
+ BinnedLayerMaterialAnalyser = CompFactory.Trk.LayerMaterialAnalyser(name="BinnedLayerMaterialAnalyser")
+ BinnedLayerMaterialAnalyser.LayerMaterialName = LayerMaterialName
+ BinnedLayerMaterialAnalyser.ValidationTreeName = 'BinnedLayerMaterialAnalyser'
+ BinnedLayerMaterialAnalyser.ValidationTreeDescription = 'Output of the BinnedLayerMaterialAnalyser'
+ BinnedLayerMaterialAnalyser.ValidationTreeFolder = '/val/BinnedLayerMaterialAnalyser'
+ LayerMaterialAnalysers = [ BinnedLayerMaterialAnalyser ]
+ kwargs.setdefault("LayerMaterialAnalysers", LayerMaterialAnalysers)
+
+ histSvc = CompFactory.THistSvc(Output = ["val DATAFILE='AtlasGeant4Geometry.root' TYPE='ROOT' OPT='RECREATE'"])
+ result.addService( histSvc )
+
+ algo = CompFactory.Trk.MaterialMapping(name=name, **kwargs)
+ result.addEventAlgo(algo, primary = True)
+
+ from RegistrationServices.OutputConditionsAlgConfig import OutputConditionsAlgCfg
+ result.merge(OutputConditionsAlgCfg(flags, name = "CondAlg_Material",
+ outputFile="AtlasLayerMaterial.pool.root",
+ ObjectList=['Trk::LayerMaterialMap#'+LayerMaterialDirectory+LayerMaterialName],
+ WriteIOV=True,IOVTagList=[flags.ITk.trackingGeometry.materialTag] ))
+
+ result.addService(CompFactory.IOVRegistrationSvc(RecreateFolders = True))
+
+ return result
+
+
+def ITkMaterialValidationCfg(flags, name="MaterialValidation", **kwargs):
+ """Return configured ComponentAccumulator and tool for MaterialMapping"""
+
+ result=ComponentAccumulator()
+
+ from PixelGeoModelXml.ITkPixelGeoModelConfig import ITkPixelGeometryCfg
+ itkPixel = ITkPixelGeometryCfg(flags)
+ result.merge(itkPixel)
+
+ from StripGeoModelXml.ITkStripGeoModelConfig import ITkStripGeometryCfg
+ itkStrip = ITkStripGeometryCfg(flags)
+ result.merge(itkStrip)
+
+ from BeamPipeGeoModel.BeamPipeGMConfig import BeamPipeGeometryCfg
+ result.merge(BeamPipeGeometryCfg(flags))
+
+ # get the correct TrackingGeometry setup
+ geom_svc=None
+ geom_cond_key=''
+ from InDetRecExample.TrackingCommon import use_tracking_geometry_cond_alg
+ if not use_tracking_geometry_cond_alg :
+ from TrkConfig.AtlasTrackingGeometrySvcConfig import TrackingGeometrySvcCfg
+ acc = TrackingGeometrySvcCfg(flags)
+ geom_svc = acc.getPrimary()
+ result.merge(acc)
+ kwargs.setdefault("TrackingGeometrySvc", geom_svc)
+ else :
+ from TrackingGeometryCondAlg.AtlasTrackingGeometryCondAlgConfig import TrackingGeometryCondAlgCfg
+ result.merge( TrackingGeometryCondAlgCfg(flags) )
+ geom_cond_key = 'AtlasTrackingGeometry'
+ kwargs.setdefault("TrackingGeometryReadKey", geom_cond_key)
+
+
+ if 'MaterialMapper' not in kwargs :
+ MaterialMapper = CompFactory.Trk.MaterialMapper("MaterialMapper")
+ kwargs.setdefault("MaterialMapper", MaterialMapper)
+
+ if 'MinEta' not in kwargs :
+ minEta = -6.
+ kwargs.setdefault("MinEta", minEta)
+
+ if 'MaxEta' not in kwargs :
+ maxEta = 6.
+ kwargs.setdefault("MaxEta", maxEta)
+
+ histSvc = CompFactory.THistSvc(Output = ["val DATAFILE='AtlasTrackingGeometry.root' TYPE='ROOT' OPT='RECREATE'"])
+ result.addService( histSvc )
+
+ algo = CompFactory.Trk.MaterialValidation(name=name, **kwargs)
+ result.addEventAlgo(algo, primary = True)
+
+ return result
diff --git a/Tracking/TrkDetDescr/TrkDetDescrAlgs/scripts/RunMaterialMappingITk.py b/Tracking/TrkDetDescr/TrkDetDescrAlgs/scripts/RunMaterialMappingITk.py
new file mode 100644
index 0000000000000000000000000000000000000000..3af26e28265b4cbb6c4fae51f81a94e79c303e9a
--- /dev/null
+++ b/Tracking/TrkDetDescr/TrkDetDescrAlgs/scripts/RunMaterialMappingITk.py
@@ -0,0 +1,116 @@
+#!/usr/bin/env python
+"""
+
+Run material mapping for tracking geometry.
+Uses as MaterialStepsCollections as input.
+
+"""
+from AthenaCommon.Configurable import Configurable
+from AthenaCommon.Logging import log
+from argparse import ArgumentParser
+from AthenaConfiguration.AllConfigFlags import ConfigFlags
+
+Configurable.configurableRun3Behavior = True
+
+# Argument parsing
+parser = ArgumentParser("RunMaterialMappingITk.py")
+parser.add_argument("detectors", metavar="detectors", type=str, nargs="*",
+ help="Specify the list of detectors")
+parser.add_argument("--localgeo", default=False, action="store_true",
+ help="Use local geometry Xml files")
+parser.add_argument("-V", "--verboseAccumulators", default=False,
+ action="store_true",
+ help="Print full details of the AlgSequence")
+parser.add_argument("-S", "--verboseStoreGate", default=False,
+ action="store_true",
+ help="Dump the StoreGate(s) each event iteration")
+parser.add_argument("--maxEvents",default=10, type=int,
+ help="The number of events to run. 0 skips execution")
+parser.add_argument("--geometrytag",default="ATLAS-P2-ITK-24-00-00", type=str,
+ help="The geometry tag to use")
+parser.add_argument("--inputfile",
+ default="MaterialStepCollection.root",
+ help="The input material step file to use")
+args = parser.parse_args()
+
+# Some info about the job
+print("----MaterialMapping for ITk geometry----")
+print()
+print("Using Geometry Tag: "+args.geometrytag)
+if args.localgeo:
+ print("...overridden by local Geometry Xml files")
+print("Input File:"+args.inputfile)
+if not args.detectors:
+ print("Running complete detector")
+else:
+ print("Running with: {}".format(", ".join(args.detectors)))
+print()
+
+LocalDataBaseName = ConfigFlags.ITk.trackingGeometry.localDatabaseName
+ConfigFlags.IOVDb.DBConnection='sqlite://;schema='+LocalDataBaseName+';dbname=OFLP200'
+
+# necessity to create a new PoolFileCatalog
+import os
+if os.path.exists('./PoolFileCatalog.xml') :
+ print('[!] PoolFileCatalog exists in the run directory (may use old PFN!)')
+ print('[>] Deleting it now !')
+ os.remove('./PoolFileCatalog.xml')
+
+ConfigFlags.Input.isMC = True
+
+if args.localgeo:
+ configFlags.GeoModel.useLocalGeometry = True
+
+from AthenaConfiguration.DetectorConfigFlags import setupDetectorsFromList
+detectors = args.detectors if 'detectors' in args and args.detectors else ['ITkPixel', 'ITkStrip']
+detectors.append('Bpipe') # always run with beam pipe
+setupDetectorsFromList(ConfigFlags, detectors, toggle_geometry=True)
+
+ConfigFlags.GeoModel.AtlasVersion = args.geometrytag
+ConfigFlags.IOVDb.GlobalTag = "OFLCOND-SIM-00-00-00"
+ConfigFlags.GeoModel.Align.Dynamic = False
+ConfigFlags.TrackingGeometry.MaterialSource = "None"
+ConfigFlags.Beam.Type =''
+
+ConfigFlags.Detector.GeometryCalo = False
+ConfigFlags.Detector.GeometryMuon = False
+
+# This should run serially for the moment.
+ConfigFlags.Concurrency.NumThreads = 1
+ConfigFlags.Concurrency.NumConcurrentEvents = 1
+
+import glob
+FileList = glob.glob(args.inputfile)
+ConfigFlags.Input.Files = FileList
+
+log.debug('Lock config flags now.')
+ConfigFlags.lock()
+
+from AthenaConfiguration.MainServicesConfig import MainServicesCfg
+cfg=MainServicesCfg(ConfigFlags)
+
+### setup dumping of additional information
+if args.verboseAccumulators:
+ cfg.printConfig(withDetails=True)
+if args.verboseStoreGate:
+ cfg.getService("StoreGateSvc").Dump = True
+
+log.debug('Dumping of ConfigFlags now.')
+ConfigFlags.dump()
+
+from AthenaPoolCnvSvc.PoolReadConfig import PoolReadCfg
+cfg.merge(PoolReadCfg(ConfigFlags))
+
+from TrkDetDescrAlgs.TrkDetDescrAlgsConfig import ITkMaterialMappingCfg
+cfg.merge(ITkMaterialMappingCfg(ConfigFlags))
+
+cfg.printConfig(withDetails = True, summariseProps = True)
+
+events = args.maxEvents
+if events<=0:
+ events = 10000000000
+cfg.run(maxEvents=events)
+f=open("MaterialMappingITk.pkl","wb")
+cfg.store(f)
+f.close()
+
diff --git a/Tracking/TrkDetDescr/TrkDetDescrAlgs/scripts/RunMaterialValidationITk.py b/Tracking/TrkDetDescr/TrkDetDescrAlgs/scripts/RunMaterialValidationITk.py
new file mode 100644
index 0000000000000000000000000000000000000000..38c1610bc4b45233e9991c10259144962ad2fe1c
--- /dev/null
+++ b/Tracking/TrkDetDescr/TrkDetDescrAlgs/scripts/RunMaterialValidationITk.py
@@ -0,0 +1,93 @@
+#!/usr/bin/env python
+"""
+
+Run material validation to check material maps for tracking geometry.
+
+"""
+
+from AthenaCommon.Configurable import Configurable
+from AthenaCommon.Logging import log
+from argparse import ArgumentParser
+from AthenaConfiguration.AllConfigFlags import ConfigFlags
+
+Configurable.configurableRun3Behavior = True
+
+# Argument parsing
+parser = ArgumentParser("RunMaterialValidationITk.py")
+parser.add_argument("detectors", metavar="detectors", type=str, nargs="*",
+ help="Specify the list of detectors")
+parser.add_argument("--localgeo", default=False, action="store_true",
+ help="Use local geometry Xml files")
+parser.add_argument("-V", "--verboseAccumulators", default=False,
+ action="store_true",
+ help="Print full details of the AlgSequence")
+parser.add_argument("-S", "--verboseStoreGate", default=False,
+ action="store_true",
+ help="Dump the StoreGate(s) each event iteration")
+parser.add_argument("--maxEvents",default=10, type=int,
+ help="The number of events to run. 0 skips execution")
+parser.add_argument("--geometrytag",default="ATLAS-P2-ITK-24-00-00", type=str,
+ help="The geometry tag to use")
+args = parser.parse_args()
+
+# Some info about the job
+print("----MaterialValidation for ITk geometry----")
+print()
+print("Using Geometry Tag: "+args.geometrytag)
+if args.localgeo:
+ print("...overridden by local Geometry Xml files")
+print()
+
+ConfigFlags.Input.isMC = True
+
+if args.localgeo:
+ configFlags.GeoModel.useLocalGeometry = True
+
+from AthenaConfiguration.DetectorConfigFlags import setupDetectorsFromList
+detectors = args.detectors if 'detectors' in args and args.detectors else ['ITkPixel', 'ITkStrip']
+detectors.append('Bpipe') # always run with beam pipe
+setupDetectorsFromList(ConfigFlags, detectors, toggle_geometry=True)
+
+ConfigFlags.GeoModel.AtlasVersion = args.geometrytag
+ConfigFlags.IOVDb.GlobalTag = "OFLCOND-SIM-00-00-00"
+ConfigFlags.GeoModel.Align.Dynamic = False
+ConfigFlags.TrackingGeometry.MaterialSource = "COOL"
+ConfigFlags.Beam.Type =''
+
+ConfigFlags.Detector.GeometryCalo = False
+ConfigFlags.Detector.GeometryMuon = False
+
+# This should run serially for the moment.
+ConfigFlags.Concurrency.NumThreads = 1
+ConfigFlags.Concurrency.NumConcurrentEvents = 1
+
+ConfigFlags.ITk.trackingGeometry.loadLocalDbForMaterialMaps=True
+LocalDataBaseName = ConfigFlags.ITk.trackingGeometry.localDatabaseName
+ConfigFlags.IOVDb.DBConnection='sqlite://;schema='+LocalDataBaseName+';dbname=OFLP200'
+
+log.debug('Lock config flags now.')
+ConfigFlags.lock()
+
+from AthenaConfiguration.MainServicesConfig import MainServicesCfg
+cfg=MainServicesCfg(ConfigFlags)
+
+### setup dumping of additional information
+if args.verboseAccumulators:
+ cfg.printConfig(withDetails=True)
+if args.verboseStoreGate:
+ cfg.getService("StoreGateSvc").Dump = True
+
+log.debug('Dumping of ConfigFlags now.')
+ConfigFlags.dump()
+
+from TrkDetDescrAlgs.TrkDetDescrAlgsConfig import ITkMaterialValidationCfg
+cfg.merge(ITkMaterialValidationCfg(ConfigFlags))
+
+cfg.printConfig(withDetails = True, summariseProps = True)
+
+## as no input file, maxEvents should be a valid number
+cfg.run(maxEvents=args.maxEvents)
+f=open("MaterialValidationITk.pkl","wb")
+cfg.store(f)
+f.close()
+
diff --git a/Tracking/TrkDetDescr/TrkDetDescrAlgs/src/MaterialMapping.cxx b/Tracking/TrkDetDescr/TrkDetDescrAlgs/src/MaterialMapping.cxx
index a8d949b23efc8866294e03fce3665dcfea17d556..0ecfab62c7a8d92ddd7a68fef64968f204cf2216 100755
--- a/Tracking/TrkDetDescr/TrkDetDescrAlgs/src/MaterialMapping.cxx
+++ b/Tracking/TrkDetDescr/TrkDetDescrAlgs/src/MaterialMapping.cxx
@@ -10,19 +10,11 @@
#include "GaudiKernel/SystemOfUnits.h"
//TrkDetDescr Algs, Interfaces, Utils
#include "TrkDetDescrAlgs/MaterialMapping.h"
-#include "TrkDetDescrInterfaces/ITrackingGeometrySvc.h"
-#include "TrkDetDescrInterfaces/IMaterialMapper.h"
-#include "TrkDetDescrInterfaces/ILayerMaterialCreator.h"
-#include "TrkDetDescrInterfaces/ILayerMaterialAnalyser.h"
#include "TrkDetDescrUtils/GeometryStatics.h"
#include "TrkDetDescrUtils/LayerIndex.h"
#include "TrkDetDescrUtils/BinUtility.h"
-// TrkExtrapolation
-#include "TrkExInterfaces/IExtrapolationEngine.h"
// TrkGeometry
#include "TrkGeometry/LayerMaterialRecord.h"
-//#include "TrkGeometry/ElementTable.h"
-#include "TrkGeometry/TrackingGeometry.h"
#include "TrkGeometry/TrackingVolume.h"
#include "TrkGeometry/MaterialStep.h"
#include "TrkGeometry/MaterialProperties.h"
@@ -43,8 +35,6 @@
Trk::MaterialMapping::MaterialMapping(const std::string& name, ISvcLocator* pSvcLocator)
: AthAlgorithm(name,pSvcLocator),
- m_trackingGeometrySvc("AtlasTrackingGeometrySvc",name),
- m_trackingGeometry(nullptr),
m_checkForEmptyHits(true),
m_mappingVolumeName("Atlas"),
m_mappingVolume(nullptr),
@@ -53,10 +43,8 @@ Trk::MaterialMapping::MaterialMapping(const std::string& name, ISvcLocator* pSvc
m_etaSide(0),
m_useLayerThickness(false),
m_associationType(1),
- m_layerMaterialRecordAnalyser(""),
m_mapMaterial(true),
- m_materialMapper(""),
- m_mapComposition(true),
+ m_mapComposition(false),
m_minCompositionFraction(0.005),
m_elementTable(nullptr),
m_inputEventElementTable("ElementTable"),
@@ -70,26 +58,18 @@ Trk::MaterialMapping::MaterialMapping(const std::string& name, ISvcLocator* pSvc
,m_memoryLogger()
#endif
{
- // the name of the TrackingGeometry to be retrieved
- declareProperty("TrackingGeometrySvc" , m_trackingGeometrySvc);
+ // the name of the volume to map
declareProperty("MappingVolumeName" , m_mappingVolumeName);
// the extrapolation engine
declareProperty("CheckForEmptyHits" , m_checkForEmptyHits);
- declareProperty("ExtrapolationEngine" , m_extrapolationEngine);
// general steering
declareProperty("EtaCutOff" , m_etaCutOff);
declareProperty("EtaSide" , m_etaSide);
- // for the analysis of the material
- declareProperty("LayerMaterialRecordAnalyser" , m_layerMaterialRecordAnalyser);
- // for the creation of the material
- declareProperty("LayerMaterialCreators" , m_layerMaterialCreators);
- declareProperty("LayerMaterialAnalysers" , m_layerMaterialAnalysers);
// the toolhandle of the MaterialMapper to be used
declareProperty("MapMaterial" , m_mapMaterial);
- declareProperty("MaterialMapper" , m_materialMapper);
- // Composition related parameters
- declareProperty("MapComposition" , m_mapComposition);
- declareProperty("MinCompositionFraction" , m_minCompositionFraction);
+ // Composition related parameters
+ declareProperty("MapComposition" , m_mapComposition);
+ declareProperty("MinCompositionFraction" , m_minCompositionFraction);
// Steer the layer thickness
declareProperty("UseActualLayerThicknesss" , m_useLayerThickness);
// some job setup
@@ -110,23 +90,26 @@ StatusCode Trk::MaterialMapping::initialize()
ATH_MSG_INFO("initialize()");
- if ( (m_trackingGeometrySvc.retrieve()).isFailure() )
- ATH_MSG_WARNING("Could not retrieve TrackingGeometrySvc");
+#ifdef LEGACY_TRKGEOM
+ if (!m_trackingGeometrySvc.empty()) {
+ ATH_CHECK( m_trackingGeometrySvc.retrieve());
+ }
+#endif
+ ATH_CHECK( m_trackingGeometryReadKey.initialize(!m_trackingGeometryReadKey.key().empty()) );
- if (m_extrapolationEngine.retrieve().isFailure())
- ATH_MSG_WARNING("Could not retrieve ExtrapolationEngine - needed for emopty hit scaling");
+ ATH_CHECK(m_extrapolationEngine.retrieve());
- if (!m_materialMapper.empty() && (m_materialMapper.retrieve()).isFailure() )
- ATH_MSG_WARNING("Could not retrieve MaterialMapper");
+ if ( !m_materialMapper.empty() )
+ ATH_CHECK( m_materialMapper.retrieve() );
- if ( !m_layerMaterialRecordAnalyser.empty() && m_layerMaterialRecordAnalyser.retrieve().isFailure() )
- ATH_MSG_WARNING("Could not retrieve LayerMaterialAnalyser");
+ if ( !m_layerMaterialRecordAnalyser.empty() )
+ ATH_CHECK( m_layerMaterialRecordAnalyser.retrieve() );
- if ( !m_layerMaterialCreators.empty() && m_layerMaterialCreators.retrieve().isFailure() )
- ATH_MSG_WARNING("Could not retrieve any LayerMaterialCreators");
+ if ( !m_layerMaterialCreators.empty() )
+ ATH_CHECK( m_layerMaterialCreators.retrieve() );
- if ( !m_layerMaterialAnalysers.empty() && m_layerMaterialAnalysers.retrieve().isFailure() )
- ATH_MSG_WARNING("Could not retrieve any LayerMaterialAnalysers");
+ if ( !m_layerMaterialAnalysers.empty() )
+ ATH_CHECK( m_layerMaterialAnalysers.retrieve() );
ATH_CHECK( m_inputMaterialStepCollection.initialize() );
ATH_CHECK( m_inputEventElementTable.initialize() );
@@ -140,189 +123,190 @@ StatusCode Trk::MaterialMapping::execute()
ATH_MSG_VERBOSE("MaterialMapping execute() start");
// ------------------------------- get the trackingGeometry at first place
- if (!m_trackingGeometry) {
- StatusCode retrieveCode = retrieveTrackingGeometry();
+ if (!m_mappingVolume) {
+ StatusCode retrieveCode = handleTrackingGeometry();
if (retrieveCode.isFailure()){
- ATH_MSG_INFO("Could not retrieve TrackingGeometry. Exiting.");
+ ATH_MSG_INFO("Could not retrieve mapping volume from tracking geometry. Exiting.");
return retrieveCode;
}
}
- if (m_trackingGeometry)
- ATH_MSG_VERBOSE("TrackingGeometry sucessfully retrieved");
+ if (m_mappingVolume)
+ ATH_MSG_VERBOSE("Mapping volume correctly retrieved from tracking geometry");
+
SG::ReadHandle<MaterialStepCollection> materialStepCollection(m_inputMaterialStepCollection);
-
+
// --------- prepare the element table ---------------------------------------------------
- SG::ReadHandle<Trk::ElementTable> eTableEvent(m_inputEventElementTable);
+ if (m_mapComposition) {
+ SG::ReadHandle<Trk::ElementTable> eTableEvent(m_inputEventElementTable);
+ (*m_elementTable) += (*eTableEvent); // accummulate the table
+ }
- (*m_elementTable) += (*eTableEvent); // accummulate the table
- m_mapComposition = eTableEvent.isValid();
-
-
- // event parameters - associated asteps, and layers hit per event
- int associatedSteps = 0;
- m_accumulatedMaterialXX0 = 0.;
- m_accumulatedRhoS = 0.;
- m_layersRecordedPerEvent.clear();
- // clearing the recorded layers per event
- if (materialStepCollection.isValid() && !materialStepCollection->empty()){
-
- // get the number of material steps
- size_t materialSteps = materialStepCollection->size();
- ATH_MSG_DEBUG("[+] Successfully read "<< materialSteps << " geantino steps");
-
- // create a direction out of the last material step
- double dirx = (*materialStepCollection)[materialSteps-1]->hitX();
- double diry = (*materialStepCollection)[materialSteps-1]->hitY();
- double dirz = (*materialStepCollection)[materialSteps-1]->hitZ();
- Amg::Vector3D direction = Amg::Vector3D(dirx,diry,dirz).unit();
-
- double eta = direction.eta();
- // skip the event if the eta cut is not met
- if ( fabs(eta) > m_etaCutOff || (m_etaSide && m_etaSide*eta < 0.) ) {
- ATH_MSG_VERBOSE("[-] Event is outside eta acceptance of " << m_etaCutOff << ". Skipping it.");
- return StatusCode::SUCCESS;
- }
-
- // now propagate through the full detector and collect the layers
- Trk::NeutralCurvilinearParameters ncP(Amg::Vector3D(0.,0.,0.), direction, 0.);
- // create a neutral extrapolation cell
- Trk::ExtrapolationCell<Trk::NeutralParameters> ecc(ncP);
- ecc.navigationCurvilinear = false;
- ecc.addConfigurationMode(Trk::ExtrapolationMode::StopAtBoundary);
- ecc.addConfigurationMode(Trk::ExtrapolationMode::CollectPassive);
- ecc.addConfigurationMode(Trk::ExtrapolationMode::CollectBoundary);
-
- // let's extrapolate through the detector and remember which layers (with material) should have been hit
- std::vector< std::pair<const Trk::Layer*, Amg::Vector3D> > layersAndHits;
- // call the extrapolation engine
- Trk::ExtrapolationCode eCode = m_extrapolationEngine->extrapolate(ecc);
- // end the parameters if there
- if (eCode.isSuccess()){
- // name of passive surfaces found
- size_t nLayersHit = ecc.extrapolationSteps.size();
- ATH_MSG_VERBOSE("[+] Extrapolation to layers did succeed and found " << nLayersHit << " layers.");
- // reserve the size of the vectors
- layersAndHits.reserve(nLayersHit);
- // for screen output
- size_t ilayer = 0;
- // find all the intersected material - remember the last parameters
- const Trk::NeutralParameters* parameters = nullptr;
- // loop over the collected information
- for (auto& es : ecc.extrapolationSteps){
- // continue if we have parameters
- parameters = es.parameters;
- if (parameters){
- const Trk::Surface& pSurface = parameters->associatedSurface();
- // get the surface with associated layer (that has material)
- ATH_MSG_VERBOSE("[L] Testing layer with associatedLayer() " << pSurface.associatedLayer() << " and materialLayer() " << pSurface.materialLayer() );
- //
- if ( ( pSurface.associatedLayer() && pSurface.associatedLayer()->layerMaterialProperties() ) || pSurface.materialLayer() ){
- // material layer
- const Trk::Layer* mLayer = pSurface.materialLayer() ? pSurface.materialLayer() : pSurface.associatedLayer();
- // record that one
- std::pair<const Trk::Layer*, Amg::Vector3D> layerHitPair(mLayer, parameters->position());
- ATH_MSG_VERBOSE("[L] Layer " << ++ilayer << " with index " << mLayer->layerIndex().value() << " hit at " << Amg::toString(parameters->position()));
- layersAndHits.push_back(layerHitPair);
- }
- delete parameters;
+
+ // event parameters - associated asteps, and layers hit per event
+ int associatedSteps = 0;
+ m_accumulatedMaterialXX0 = 0.;
+ m_accumulatedRhoS = 0.;
+ m_layersRecordedPerEvent.clear();
+ // clearing the recorded layers per event
+ if (materialStepCollection.isValid() && !materialStepCollection->empty()){
+
+ // get the number of material steps
+ size_t materialSteps = materialStepCollection->size();
+ ATH_MSG_DEBUG("[+] Successfully read "<< materialSteps << " geantino steps");
+
+ // create a direction out of the last material step
+ double dirx = (*materialStepCollection)[materialSteps-1]->hitX();
+ double diry = (*materialStepCollection)[materialSteps-1]->hitY();
+ double dirz = (*materialStepCollection)[materialSteps-1]->hitZ();
+ Amg::Vector3D direction = Amg::Vector3D(dirx,diry,dirz).unit();
+
+ double eta = direction.eta();
+ // skip the event if the eta cut is not met
+ if ( fabs(eta) > m_etaCutOff || (m_etaSide && m_etaSide*eta < 0.) ) {
+ ATH_MSG_VERBOSE("[-] Event is outside eta acceptance of " << m_etaCutOff << ". Skipping it.");
+ return StatusCode::SUCCESS;
+ }
+
+ // now propagate through the full detector and collect the layers
+ Trk::NeutralCurvilinearParameters ncP(Amg::Vector3D(0.,0.,0.), direction, 0.);
+ // create a neutral extrapolation cell
+ Trk::ExtrapolationCell<Trk::NeutralParameters> ecc(ncP);
+ ecc.navigationCurvilinear = false;
+ ecc.addConfigurationMode(Trk::ExtrapolationMode::StopAtBoundary);
+ ecc.addConfigurationMode(Trk::ExtrapolationMode::CollectPassive);
+ ecc.addConfigurationMode(Trk::ExtrapolationMode::CollectBoundary);
+
+ // let's extrapolate through the detector and remember which layers (with material) should have been hit
+ std::vector< std::pair<const Trk::Layer*, Amg::Vector3D> > layersAndHits;
+ // call the extrapolation engine
+ Trk::ExtrapolationCode eCode = m_extrapolationEngine->extrapolate(ecc);
+ // end the parameters if there
+ if (eCode.isSuccess()){
+ // name of passive surfaces found
+ size_t nLayersHit = ecc.extrapolationSteps.size();
+ ATH_MSG_VERBOSE("[+] Extrapolation to layers did succeed and found " << nLayersHit << " layers.");
+ // reserve the size of the vectors
+ layersAndHits.reserve(nLayersHit);
+ // for screen output
+ size_t ilayer = 0;
+ // find all the intersected material - remember the last parameters
+ const Trk::NeutralParameters* parameters = nullptr;
+ // loop over the collected information
+ for (auto& es : ecc.extrapolationSteps){
+ // continue if we have parameters
+ parameters = es.parameters;
+ if (parameters){
+ const Trk::Surface& pSurface = parameters->associatedSurface();
+ // get the surface with associated layer (that has material)
+ ATH_MSG_VERBOSE("[L] Testing layer with associatedLayer() " << pSurface.associatedLayer() << " and materialLayer() " << pSurface.materialLayer() );
+ //
+ if ( ( pSurface.associatedLayer() && pSurface.associatedLayer()->layerMaterialProperties() ) || pSurface.materialLayer() ){
+ // material layer
+ const Trk::Layer* mLayer = pSurface.materialLayer() ? pSurface.materialLayer() : pSurface.associatedLayer();
+ // record that one
+ std::pair<const Trk::Layer*, Amg::Vector3D> layerHitPair(mLayer, parameters->position());
+ ATH_MSG_VERBOSE("[L] Layer " << ++ilayer << " with index " << mLayer->layerIndex().value() << " hit at " << Amg::toString(parameters->position()));
+ layersAndHits.push_back(layerHitPair);
}
- }
- // cleanup of the final hits
- if (ecc.endParameters != parameters) delete ecc.endParameters;
-
- // we have no layers and Hits
- if (layersAndHits.empty()){
- ATH_MSG_VERBOSE("[!] No Layer was intersected - skipping.");
- return StatusCode::SUCCESS;
- }
-
- // layers are ordered, hence you can move the starting point along
- size_t currentLayer = 0;
- // loop through hits and find the closest layer, the start point moves outwards as we go
- for ( const Trk::MaterialStep* step : *materialStepCollection ) {
- // verbose output
- ATH_MSG_VERBOSE("[L] starting from layer " << currentLayer << " from layer collection for this step.");
- // step length and position
- double t = step->steplength();
- Amg::Vector3D pos(step->hitX(), step->hitY(), step->hitZ());
- // skip if :
- // -- 0) no mapping volume exists
- // -- 1) outside the mapping volume
- // -- 2) outside the eta acceptance
- if (!m_mappingVolume || !(m_mappingVolume->inside(pos)) || fabs(pos.eta()) > m_etaCutOff ){
- ++m_skippedOutside;
- continue;
- }
- // now find the closest layer
- // (a) if the currentLayer is the last layer and the hit is still inside -> assign
- if (currentLayer < nLayersHit-1) {
- // search through the layers - this is the reference distance for projection
- double currentDistance = (pos-layersAndHits[currentLayer].second).mag();
- ATH_MSG_VERBOSE("- current distance is " << currentDistance << " from " << Amg::toString(pos) << " and " << Amg::toString(layersAndHits[currentLayer].second) );
- for (size_t testLayer = (currentLayer+1); testLayer < nLayersHit; ++testLayer){
- // calculate teh distance to the testLayer
- double testDistance = (pos-layersAndHits[testLayer].second).mag();
- ATH_MSG_VERBOSE("[L] Testing layer " << testLayer << " from layer collection for this step.");
- ATH_MSG_VERBOSE("- test distance is " << testDistance << " from " << Amg::toString(pos) << " and " << Amg::toString(layersAndHits[testLayer].second) );
- if ( testDistance < currentDistance ){
- // screen output
- ATH_MSG_VERBOSE("[L] Skipping over to current layer " << testLayer << " because " << testDistance << " < " << currentDistance);
- // the test distance did shrink - update currentLayer
- currentLayer = testLayer;
- currentDistance = testDistance;
- } else {
- // stick to the layer you have
- break;
- }
- }
- }
- // the currentLayer *should* be correct now
- const Trk::Layer* assignedLayer = layersAndHits[currentLayer].first;
- Amg::Vector3D assignedPosition = layersAndHits[currentLayer].second;
- // associate the hit
- // (1) count it
- ++associatedSteps;
- // (2) associate it
- associateHit(*assignedLayer, pos, assignedPosition, t, step->fullMaterial());
- } // loop over material Steps
+ delete parameters;
+ }
+ }
+ // cleanup of the final hits
+ if (ecc.endParameters != parameters) delete ecc.endParameters;
- // check for the empty hits - they need to be taken into account
- ATH_MSG_VERBOSE("Found " << layersAndHits.size() << " intersected layers - while having " << m_layersRecordedPerEvent.size() << " recorded ones.");
+ // we have no layers and Hits
+ if (layersAndHits.empty()){
+ ATH_MSG_VERBOSE("[!] No Layer was intersected - skipping.");
+ return StatusCode::SUCCESS;
+ }
- // now - cross-chek if you have additional layers
- for ( auto& lhp : layersAndHits){
- // check if you find the layer int he already done record-map : not found - we need to do an empty hit scaling
- if (m_layersRecordedPerEvent.find(lhp.first) == m_layersRecordedPerEvent.end()){
- // try to find the layer material record
- auto clIter = m_layerRecords.find(lhp.first);
- if (clIter != m_layerRecords.end() ){
- (*clIter).second.associateEmptyHit(lhp.second);
- ATH_MSG_VERBOSE("- to layer with index "<< lhp.first->layerIndex().value() << " with empty hit detected.");
- } else
- ATH_MSG_WARNING("- no Layer found in the associated map! Should not happen.");
- }
+ // layers are ordered, hence you can move the starting point along
+ size_t currentLayer = 0;
+ // loop through hits and find the closest layer, the start point moves outwards as we go
+ for ( const Trk::MaterialStep* step : *materialStepCollection ) {
+ // verbose output
+ ATH_MSG_VERBOSE("[L] starting from layer " << currentLayer << " from layer collection for this step.");
+ // step length and position
+ double t = step->steplength();
+ Amg::Vector3D pos(step->hitX(), step->hitY(), step->hitZ());
+ // skip if :
+ // -- 0) no mapping volume exists
+ // -- 1) outside the mapping volume
+ // -- 2) outside the eta acceptance
+ if (!m_mappingVolume || !(m_mappingVolume->inside(pos)) || fabs(pos.eta()) > m_etaCutOff ){
+ ++m_skippedOutside;
+ continue;
}
-
- // check whether the event was good for at least one hit
- if (associatedSteps) {
- ATH_MSG_VERBOSE("There are associated steps, need to call finalizeEvent() & record to the MaterialMapper.");
- // finalize the event --------------------- Layers ---------------------------------------------
- for (auto& lRecord : m_layerRecords ) {
- // associated material
- Trk::AssociatedMaterial* assMatHit = lRecord.second.finalizeEvent((*lRecord.first));
- // record the full layer hit
- if (assMatHit && !m_materialMapper.empty()) m_materialMapper->recordLayerHit(*assMatHit, true);
- delete assMatHit;
- // call the material mapper finalize method
- ATH_MSG_VERBOSE("Calling finalizeEvent on the MaterialMapper ...");
- }
- } // the event had at least one associated hit
-
- } // end of eCode.success : needed for new mapping schema
-
- } // material steps existed
+ // now find the closest layer
+ // (a) if the currentLayer is the last layer and the hit is still inside -> assign
+ if (currentLayer < nLayersHit-1) {
+ // search through the layers - this is the reference distance for projection
+ double currentDistance = (pos-layersAndHits[currentLayer].second).mag();
+ ATH_MSG_VERBOSE("- current distance is " << currentDistance << " from " << Amg::toString(pos) << " and " << Amg::toString(layersAndHits[currentLayer].second) );
+ for (size_t testLayer = (currentLayer+1); testLayer < nLayersHit; ++testLayer){
+ // calculate teh distance to the testLayer
+ double testDistance = (pos-layersAndHits[testLayer].second).mag();
+ ATH_MSG_VERBOSE("[L] Testing layer " << testLayer << " from layer collection for this step.");
+ ATH_MSG_VERBOSE("- test distance is " << testDistance << " from " << Amg::toString(pos) << " and " << Amg::toString(layersAndHits[testLayer].second) );
+ if ( testDistance < currentDistance ){
+ // screen output
+ ATH_MSG_VERBOSE("[L] Skipping over to current layer " << testLayer << " because " << testDistance << " < " << currentDistance);
+ // the test distance did shrink - update currentLayer
+ currentLayer = testLayer;
+ currentDistance = testDistance;
+ } else {
+ // stick to the layer you have
+ break;
+ }
+ }
+ }
+ // the currentLayer *should* be correct now
+ const Trk::Layer* assignedLayer = layersAndHits[currentLayer].first;
+ Amg::Vector3D assignedPosition = layersAndHits[currentLayer].second;
+ // associate the hit
+ // (1) count it
+ ++associatedSteps;
+ // (2) associate it
+ associateHit(*assignedLayer, pos, assignedPosition, t, step->fullMaterial());
+ } // loop over material Steps
+
+ // check for the empty hits - they need to be taken into account
+ ATH_MSG_VERBOSE("Found " << layersAndHits.size() << " intersected layers - while having " << m_layersRecordedPerEvent.size() << " recorded ones.");
+
+ // now - cross-chek if you have additional layers
+ for ( auto& lhp : layersAndHits){
+ // check if you find the layer int he already done record-map : not found - we need to do an empty hit scaling
+ if (m_layersRecordedPerEvent.find(lhp.first) == m_layersRecordedPerEvent.end()){
+ // try to find the layer material record
+ auto clIter = m_layerRecords.find(lhp.first);
+ if (clIter != m_layerRecords.end() ){
+ (*clIter).second.associateEmptyHit(lhp.second);
+ ATH_MSG_VERBOSE("- to layer with index "<< lhp.first->layerIndex().value() << " with empty hit detected.");
+ } else
+ ATH_MSG_WARNING("- no Layer found in the associated map! Should not happen.");
+ }
+ }
+
+ // check whether the event was good for at least one hit
+ if (associatedSteps) {
+ ATH_MSG_VERBOSE("There are associated steps, need to call finalizeEvent() & record to the MaterialMapper.");
+ // finalize the event --------------------- Layers ---------------------------------------------
+ for (auto& lRecord : m_layerRecords ) {
+ // associated material
+ Trk::AssociatedMaterial* assMatHit = lRecord.second.finalizeEvent((*lRecord.first));
+ // record the full layer hit
+ if (assMatHit && !m_materialMapper.empty()) m_materialMapper->recordLayerHit(*assMatHit, true);
+ delete assMatHit;
+ // call the material mapper finalize method
+ ATH_MSG_VERBOSE("Calling finalizeEvent on the MaterialMapper ...");
+ }
+ } // the event had at least one associated hit
+
+ } // end of eCode.success : needed for new mapping schema
+
+ } // material steps existed
return StatusCode::SUCCESS;
@@ -341,7 +325,7 @@ bool Trk::MaterialMapping::associateHit( const Trk::Layer& associatedLayer,
const Trk::Layer* layer = &associatedLayer;
// get the associated volume
- const Trk::TrackingVolume* associatedVolume = m_trackingGeometry->lowestTrackingVolume(pos);
+ const Trk::TrackingVolume* associatedVolume = trackingGeometry().lowestTrackingVolume(pos);
// try to find the layer material record
auto clIter = m_layerRecords.find(layer);
@@ -525,30 +509,23 @@ StatusCode Trk::MaterialMapping::finalize()
}
-StatusCode Trk::MaterialMapping::retrieveTrackingGeometry()
+StatusCode Trk::MaterialMapping::handleTrackingGeometry()
{
-
- // Retrieve the TrackingGeometry from the DetectorStore
- if ((detStore()->retrieve(m_trackingGeometry, m_trackingGeometrySvc->trackingGeometryName())).isFailure()) {
- ATH_MSG_FATAL("Could not retrieve TrackingGeometry from DetectorStore!");
- return StatusCode::FAILURE;
- }
-
// either get a string volume or the highest one
- const Trk::TrackingVolume* trackingVolume = m_trackingGeometry->highestTrackingVolume();
+ const Trk::TrackingVolume* trackingVolume = trackingGeometry().highestTrackingVolume();
// prepare the mapping volume
- m_mappingVolume = m_trackingGeometry->trackingVolume(m_mappingVolumeName);
+ m_mappingVolume = trackingGeometry().trackingVolume(m_mappingVolumeName);
// register the confined layers from the TrackingVolume
registerVolume(*trackingVolume, 0);
ATH_MSG_INFO("Add "<< m_layerRecords.size() << " confined volume layers to mapping setup.");
- ATH_MSG_INFO("Add "<< m_trackingGeometry->boundaryLayers().size() << " boundary layers to mapping setup.");
+ ATH_MSG_INFO("Add "<< trackingGeometry().boundaryLayers().size() << " boundary layers to mapping setup.");
// register the layers from boundary surfaces
- auto bLayerIter = m_trackingGeometry->boundaryLayers().begin();
- for (; bLayerIter != m_trackingGeometry->boundaryLayers().end(); ++bLayerIter)
+ auto bLayerIter = trackingGeometry().boundaryLayers().begin();
+ for (; bLayerIter != trackingGeometry().boundaryLayers().end(); ++bLayerIter)
insertLayerMaterialRecord(*bLayerIter->first);
ATH_MSG_INFO("Map for "<< m_layerRecords.size() << " layers booked & prepared for mapping procedure");
@@ -632,3 +609,10 @@ void Trk::MaterialMapping::insertLayerMaterialRecord(const Trk::Layer& lay){
}
}
+void Trk::MaterialMapping::throwFailedToGetTrackingGeometry() const {
+ std::stringstream msg;
+ msg << "Failed to get conditions data " << m_trackingGeometryReadKey.key() << ".";
+ throw std::runtime_error(msg.str());
+}
+
+
diff --git a/Tracking/TrkDetDescr/TrkDetDescrAlgs/src/MaterialValidation.cxx b/Tracking/TrkDetDescr/TrkDetDescrAlgs/src/MaterialValidation.cxx
index 93f10c158c4d17db03c8668233ffbc3659023519..cc027352e509f357bc9b85e2e7e596da6879b079 100755
--- a/Tracking/TrkDetDescr/TrkDetDescrAlgs/src/MaterialValidation.cxx
+++ b/Tracking/TrkDetDescr/TrkDetDescrAlgs/src/MaterialValidation.cxx
@@ -11,12 +11,10 @@
// Trk
#include "TrkDetDescrAlgs/MaterialValidation.h"
#include "TrkGeometry/TrackingVolume.h"
-#include "TrkGeometry/TrackingGeometry.h"
#include "TrkGeometry/Layer.h"
#include "TrkGeometry/LayerMaterialProperties.h"
#include "TrkGeometry/MaterialProperties.h"
#include "TrkGeometry/AssociatedMaterial.h"
-#include "TrkDetDescrInterfaces/ITrackingGeometrySvc.h"
#include "TrkDetDescrInterfaces/IMaterialMapper.h"
#include "TrkNeutralParameters/NeutralParameters.h"
@@ -29,10 +27,8 @@
Trk::MaterialValidation::MaterialValidation(const std::string& name, ISvcLocator* pSvcLocator)
: AthAlgorithm(name,pSvcLocator) ,
- m_trackingGeometrySvc("AtlasTrackingGeometrySvc",name),
- m_trackingGeometry(nullptr),
m_materialMapper("Trk::MaterialMapper/MappingMaterialMapper"),
- m_maxMaterialMappingEvents(25000),
+ m_maxMaterialValidationEvents(25000),
m_flatDist(nullptr),
m_etaMin(-3.),
m_etaMax(3.),
@@ -40,12 +36,10 @@ Trk::MaterialValidation::MaterialValidation(const std::string& name, ISvcLocator
m_accTinX0(0)
{
- // ---------------------- The TrackingGeometrySvc ------------------------ //
- declareProperty("TrackingGeometrySvc" , m_trackingGeometrySvc);
// ---------------------- The Material Mapping -------------------------- //
// the toolhandle of the MaterialMapper to be used
declareProperty("MaterialMapper" , m_materialMapper);
- declareProperty("MaximumMappingEvents" , m_maxMaterialMappingEvents);
+ declareProperty("MaximumMappingEvents" , m_maxMaterialValidationEvents);
// ---------------------- Range setup ----------------------------------- //
declareProperty("MinEta" , m_etaMin);
declareProperty("MaxEta" , m_etaMax);
@@ -63,11 +57,12 @@ StatusCode Trk::MaterialValidation::initialize()
{
// Get the TrackingGeometry from StoreGate
- // initialize the TrackingGeometrySvc
- if (m_trackingGeometrySvc.retrieve().isFailure()) {
- ATH_MSG_FATAL( "Cannot retrieve TrackingGeometrySvc. Abort job. " );
- return StatusCode::FAILURE;
+#ifdef LEGACY_TRKGEOM
+ if (!m_trackingGeometrySvc.empty()) {
+ ATH_CHECK( m_trackingGeometrySvc.retrieve());
}
+#endif
+ ATH_CHECK( m_trackingGeometryReadKey.initialize(!m_trackingGeometryReadKey.key().empty()) );
if ( (m_materialMapper.retrieve()).isFailure() )
ATH_MSG_WARNING("Could not retrieve MaterialMapper");
@@ -82,15 +77,6 @@ StatusCode Trk::MaterialValidation::execute()
{
ATH_MSG_VERBOSE( "MaterialValidation execute() start ================================================" );
- // ------------------------------- get the trackingGeometry at first place
- if (!m_trackingGeometry) {
- StatusCode retrieveCode = retrieveTrackingGeometry();
- if (retrieveCode.isFailure()){
- ATH_MSG_INFO( "Could not retrieve TrackingGeometry. Exiting." );
- return retrieveCode;
- }
- }
-
// create the random direction - flat in eta
double eta = m_etaMin + (m_etaMax-m_etaMin)*m_flatDist->shoot();
double theta = 2.*atan(exp(-eta));
@@ -104,7 +90,7 @@ StatusCode Trk::MaterialValidation::execute()
ATH_MSG_DEBUG("[>] Start mapping event with phi | eta = " << phi << " | " << direction.eta());
// find the start TrackingVolume
- const Trk::TrackingVolume* sVolume = m_trackingGeometry->lowestTrackingVolume(position);
+ const Trk::TrackingVolume* sVolume = trackingGeometry().lowestTrackingVolume(position);
const Trk::TrackingVolume* nVolume = sVolume;
while (nVolume ) {
Trk::PositionAtBoundary paB = collectMaterialAndExit(*nVolume, position, direction);
@@ -312,13 +298,10 @@ StatusCode Trk::MaterialValidation::finalize()
return StatusCode::SUCCESS;
}
-StatusCode Trk::MaterialValidation::retrieveTrackingGeometry()
-{
-
- // Retrieve the TrackingGeometry from the DetectorStore
- if ((detStore()->retrieve(m_trackingGeometry, m_trackingGeometrySvc->trackingGeometryName())).isFailure()) {
- ATH_MSG_FATAL( "Could not retrieve TrackingGeometry from DetectorStore!" );
- return StatusCode::FAILURE;
- }
- return StatusCode::SUCCESS;
+void Trk::MaterialValidation::throwFailedToGetTrackingGeometry() const {
+ std::stringstream msg;
+ msg << "Failed to get conditions data " << m_trackingGeometryReadKey.key() << ".";
+ throw std::runtime_error(msg.str());
}
+
+
diff --git a/Tracking/TrkDetDescr/TrkDetDescrTestTools/TrkDetDescrTestTools/LayerMaterialAnalyser.h b/Tracking/TrkDetDescr/TrkDetDescrTestTools/TrkDetDescrTestTools/LayerMaterialAnalyser.h
index 66d41538fa2470baf06e0341af4e10e27cae67e9..c0a87deeb88b5c3f685aae83f4bbbe0ec502715d 100644
--- a/Tracking/TrkDetDescr/TrkDetDescrTestTools/TrkDetDescrTestTools/LayerMaterialAnalyser.h
+++ b/Tracking/TrkDetDescr/TrkDetDescrTestTools/TrkDetDescrTestTools/LayerMaterialAnalyser.h
@@ -74,25 +74,25 @@ namespace Trk {
std::string m_validationTreeDescription; //!< validation tree description - second argument in TTree
std::string m_validationTreeFolder; //!< stream/folder to for the TTree to be written out
- mutable int m_layerIndex; //!< the layer index given by the TrackingGeometry
- mutable int m_layerType; //!< the type of the layer 1 - cylinder, 2 - disk
- mutable float m_layerTranslation[3]{}; //!< center of the transform
- mutable float m_layerRotation[9]{}; //!< orientation of the layer
- mutable float m_layerDimension0; //!< dimension 0 : cylinder r, disk r_min
- mutable float m_layerDimension1; //!< dimension 1 : cylinder z, disk r_max
- mutable int m_layerBins; //!< total number of bins - loc0 * loc 1
- mutable int m_layerBins0; //!< total number of bins - loc 0
- mutable int m_layerBins1; //!< total number of bins - loc 0
- mutable std::vector<int> m_bin0 {LAYERMAXBINS}; //!< bin 0
- mutable std::vector<int> m_bin1 {LAYERMAXBINS}; //!< bin 1
- mutable std::vector<float> m_thickness {LAYERMAXBINS}; //!< gathered thickness from material mapping/material properties
- mutable std::vector<float> m_X0 {LAYERMAXBINS}; //!< gathered X0 from material mapping/material properties
- mutable std::vector<float> m_L0 {LAYERMAXBINS}; //!< gathered L0 from material mapping/material properties
- mutable std::vector<float> m_A {LAYERMAXBINS}; //!< gathered A from material mapping/material properties
- mutable std::vector<float> m_Z {LAYERMAXBINS}; //!< gathered Z from material mapping/material properties
- mutable std::vector<float> m_Rho {LAYERMAXBINS}; //!< gathered rho from material mapping/material properties
- mutable std::vector<int> m_elements {LAYERMAXBINS}; //!< gathered number of elements from material mapping/material properties
- mutable std::vector<int> m_binCounter {LAYERMAXBINS}; //!< how often was this bin hit / used
+ mutable int m_layerIndex; //!< the layer index given by the TrackingGeometry
+ mutable int m_layerType; //!< the type of the layer 1 - cylinder, 2 - disk
+ std::vector<float>* m_layerTranslation; //!< center of the transform
+ std::vector<float>* m_layerRotation; //!< orientation of the layer
+ mutable float m_layerDimension0; //!< dimension 0 : cylinder r, disk r_min
+ mutable float m_layerDimension1; //!< dimension 1 : cylinder z, disk r_max
+ mutable int m_layerBins; //!< total number of bins - loc0 * loc 1
+ mutable int m_layerBins0; //!< total number of bins - loc 0
+ mutable int m_layerBins1; //!< total number of bins - loc 0
+ std::vector<int>* m_bin0; //!< bin 0
+ std::vector<int>* m_bin1; //!< bin 1
+ std::vector<float>* m_thickness; //!< gathered thickness from material mapping/material properties
+ std::vector<float>* m_X0; //!< gathered X0 from material mapping/material properties
+ std::vector<float>* m_L0; //!< gathered L0 from material mapping/material properties
+ std::vector<float>* m_A; //!< gathered A from material mapping/material properties
+ std::vector<float>* m_Z; //!< gathered Z from material mapping/material properties
+ std::vector<float>* m_Rho; //!< gathered rho from material mapping/material properties
+ std::vector<int>* m_elements; //!< gathered number of elements from material mapping/material properties
+ std::vector<int>* m_binCounter; //!< how often was this bin hit / used
};
diff --git a/Tracking/TrkDetDescr/TrkDetDescrTestTools/src/LayerMaterialAnalyser.cxx b/Tracking/TrkDetDescr/TrkDetDescrTestTools/src/LayerMaterialAnalyser.cxx
index 4fccf295713b8ef034bfe16844dad1b2b006fa15..971ea4ddb37182f3bf148a1a2699296702065c08 100644
--- a/Tracking/TrkDetDescr/TrkDetDescrTestTools/src/LayerMaterialAnalyser.cxx
+++ b/Tracking/TrkDetDescr/TrkDetDescrTestTools/src/LayerMaterialAnalyser.cxx
@@ -37,15 +37,25 @@ Trk::LayerMaterialAnalyser::LayerMaterialAnalyser(const std::string& t, const st
m_validationTreeName("LayerMaterialAnalyser"),
m_validationTreeDescription("LayerMaterialAnalyser information"),
m_validationTreeFolder("/val/LayerMaterialAnalyser"),
- m_layerIndex{},
- m_layerType{},
- m_layerDimension0{},
- m_layerDimension1{},
- m_layerBins{},
- m_layerBins0{},
- m_layerBins1{}
-
-
+ m_layerIndex(0),
+ m_layerType(0),
+ m_layerTranslation(nullptr),
+ m_layerRotation(nullptr),
+ m_layerDimension0(0.),
+ m_layerDimension1(0.),
+ m_layerBins(0),
+ m_layerBins0(0),
+ m_layerBins1(0),
+ m_bin0(nullptr),
+ m_bin1(nullptr),
+ m_thickness(nullptr),
+ m_X0(nullptr),
+ m_L0(nullptr),
+ m_A(nullptr),
+ m_Z(nullptr),
+ m_Rho(nullptr),
+ m_elements(nullptr),
+ m_binCounter(nullptr)
{
declareInterface<Trk::ILayerMaterialAnalyser>(this);
// give the map a name
@@ -63,7 +73,20 @@ Trk::LayerMaterialAnalyser::~LayerMaterialAnalyser()
// initialize
StatusCode Trk::LayerMaterialAnalyser::initialize()
{
-
+
+ m_layerTranslation = new std::vector<float>(3, 0.);
+ m_layerRotation = new std::vector<float>(9, 0.);
+ m_bin0 = new std::vector<int>(LAYERMAXBINS, 0);
+ m_bin1 = new std::vector<int>(LAYERMAXBINS, 0);
+ m_thickness = new std::vector<float>(LAYERMAXBINS, 0.);
+ m_X0 = new std::vector<float>(LAYERMAXBINS, 0.);
+ m_L0 = new std::vector<float>(LAYERMAXBINS, 0.);
+ m_A = new std::vector<float>(LAYERMAXBINS, 0.);
+ m_Z = new std::vector<float>(LAYERMAXBINS, 0.);
+ m_Rho = new std::vector<float>(LAYERMAXBINS, 0.);
+ m_elements = new std::vector<int>(LAYERMAXBINS, 0);
+ m_binCounter = new std::vector<int>(LAYERMAXBINS, 0);
+
// now register the Tree
ITHistSvc* tHistSvc = nullptr;
@@ -71,25 +94,25 @@ StatusCode Trk::LayerMaterialAnalyser::initialize()
m_validationTree = new TTree(m_validationTreeName.c_str(), m_validationTreeDescription.c_str());
// position coordinates of the update
- m_validationTree->Branch("LayerIndex", &m_layerIndex, "lIndex/I");
- m_validationTree->Branch("LayerType", &m_layerType, "lType/I");
- m_validationTree->Branch("LayerTranslation", m_layerTranslation, "lTranslation[3]/F");
- m_validationTree->Branch("LayerRotation", m_layerRotation, "lRotation[9]/F");
- m_validationTree->Branch("LayerDimension0", &m_layerDimension0, "lDimension0/F");
- m_validationTree->Branch("LayerDimension1", &m_layerDimension1, "lDimension1/F");
- m_validationTree->Branch("LayerBins", &m_layerBins, "lBins/I");
- m_validationTree->Branch("LayerBins0", &m_layerBins0, "lBins0/I");
- m_validationTree->Branch("LayerBins1", &m_layerBins1, "lBins1/I");
- m_validationTree->Branch("LayerBin0", m_bin0.data(), "lBin0[lBins]/I");
- m_validationTree->Branch("LayerBin1", m_bin1.data(), "lBin1[lBins]/I");
- m_validationTree->Branch("LayerBinCounter", m_binCounter.data(), "lBinC[lBins]/I");
- m_validationTree->Branch("LayerThickness", m_thickness.data(), "lt[lBins]/F");
- m_validationTree->Branch("LayerX0", m_X0.data(), "lX0[lBins]/F");
- m_validationTree->Branch("LayerL0", m_L0.data(), "lL0[lBins]/F");
- m_validationTree->Branch("LayerA", m_A.data(), "lA[lBins]/F");
- m_validationTree->Branch("LayerZ", m_Z.data(), "lZ[lBins]/F");
- m_validationTree->Branch("LayerRo", m_Rho.data(), "lRho[lBins]/F");
- m_validationTree->Branch("LayerElements", m_elements.data(), "lElements[lBins]/I");
+ m_validationTree->Branch("LayerIndex", &m_layerIndex );
+ m_validationTree->Branch("LayerType", &m_layerType );
+ m_validationTree->Branch("LayerTranslation", &m_layerTranslation );
+ m_validationTree->Branch("LayerRotation", &m_layerRotation );
+ m_validationTree->Branch("LayerDimension0", &m_layerDimension0 );
+ m_validationTree->Branch("LayerDimension1", &m_layerDimension1 );
+ m_validationTree->Branch("LayerBins", &m_layerBins );
+ m_validationTree->Branch("LayerBins0", &m_layerBins0 );
+ m_validationTree->Branch("LayerBins1", &m_layerBins1 );
+ m_validationTree->Branch("LayerBin0", &m_bin0 );
+ m_validationTree->Branch("LayerBin1", &m_bin1 );
+ m_validationTree->Branch("LayerBinCounter", &m_binCounter );
+ m_validationTree->Branch("LayerThickness", &m_thickness );
+ m_validationTree->Branch("LayerX0", &m_X0 );
+ m_validationTree->Branch("LayerL0", &m_L0 );
+ m_validationTree->Branch("LayerA", &m_A );
+ m_validationTree->Branch("LayerZ", &m_Z );
+ m_validationTree->Branch("LayerRo", &m_Rho );
+ m_validationTree->Branch("LayerElements", &m_elements );
// now register the Tree
if (service("THistSvc",tHistSvc).isFailure()) {
@@ -109,6 +132,18 @@ StatusCode Trk::LayerMaterialAnalyser::initialize()
// finalize
StatusCode Trk::LayerMaterialAnalyser::finalize()
{
+ delete m_layerTranslation ;
+ delete m_layerRotation ;
+ delete m_bin0 ;
+ delete m_bin1 ;
+ delete m_thickness ;
+ delete m_X0 ;
+ delete m_L0 ;
+ delete m_A ;
+ delete m_Z ;
+ delete m_Rho ;
+ delete m_elements ;
+ delete m_binCounter ;
return StatusCode::SUCCESS;
}
@@ -167,24 +202,24 @@ StatusCode Trk::LayerMaterialAnalyser::analyse(const Trk::Layer& layer,
const Trk::MaterialPropertiesMatrix& mpMatrix,
const std::vector< std::vector< unsigned int > >* bCounter ) const
{
-
+
// general layer information
m_layerIndex = layer.layerIndex().value();
const Trk::Surface& lSurface = layer.surfaceRepresentation();
- m_layerTranslation[0] = lSurface.center().x();
- m_layerTranslation[1] = lSurface.center().y();
- m_layerTranslation[2] = lSurface.center().z();
+ m_layerTranslation->at(0) = lSurface.center().x();
+ m_layerTranslation->at(1) = lSurface.center().y();
+ m_layerTranslation->at(2) = lSurface.center().z();
AmgMatrix(3,3) rMatrix = lSurface.transform().rotation();
- m_layerRotation[0] = rMatrix(0,0);
- m_layerRotation[1] = rMatrix(1,0);
- m_layerRotation[2] = rMatrix(2,0);
- m_layerRotation[3] = rMatrix(0,1);
- m_layerRotation[4] = rMatrix(1,1);
- m_layerRotation[5] = rMatrix(2,1);
- m_layerRotation[6] = rMatrix(0,2);
- m_layerRotation[7] = rMatrix(1,2);
- m_layerRotation[8] = rMatrix(2,2);
+ m_layerRotation->at(0) = rMatrix(0,0);
+ m_layerRotation->at(1) = rMatrix(1,0);
+ m_layerRotation->at(2) = rMatrix(2,0);
+ m_layerRotation->at(3) = rMatrix(0,1);
+ m_layerRotation->at(4) = rMatrix(1,1);
+ m_layerRotation->at(5) = rMatrix(2,1);
+ m_layerRotation->at(6) = rMatrix(0,2);
+ m_layerRotation->at(7) = rMatrix(1,2);
+ m_layerRotation->at(8) = rMatrix(2,2);
// cylinder bounds
if ( lSurface.type() == Trk::SurfaceType::Cylinder ){
@@ -213,29 +248,29 @@ StatusCode Trk::LayerMaterialAnalyser::analyse(const Trk::Layer& layer,
for (const auto & outerIter : mpMatrix){
int bin0 = 0;
for (const auto & innerIter : outerIter ){
- m_bin0[m_layerBins] = bin0;
- m_bin1[m_layerBins] = bin1;
+ m_bin0->at(m_layerBins) = bin0;
+ m_bin1->at(m_layerBins) = bin1;
// get the material
const Trk::MaterialProperties* mProperties = innerIter;
if (mProperties){
- m_thickness[m_layerBins] = mProperties->thickness();
- m_X0[m_layerBins] = mProperties->x0();
- m_L0[m_layerBins] = mProperties->l0();
- m_A[m_layerBins] = mProperties->averageA();
- m_Z[m_layerBins] = mProperties->averageZ();
- m_Rho[m_layerBins] = mProperties->averageRho();
- m_elements[m_layerBins] = mProperties->material().composition ? mProperties->material().composition->size() : 0;
+ m_thickness->at(m_layerBins) = mProperties->thickness();
+ m_X0->at(m_layerBins) = mProperties->x0();
+ m_L0->at(m_layerBins) = mProperties->l0();
+ m_A->at(m_layerBins) = mProperties->averageA();
+ m_Z->at(m_layerBins) = mProperties->averageZ();
+ m_Rho->at(m_layerBins) = mProperties->averageRho();
+ m_elements->at(m_layerBins) = mProperties->material().composition ? mProperties->material().composition->size() : 0;
} else {
- m_thickness[m_layerBins] = 0.;
- m_X0[m_layerBins] = 0.;
- m_L0[m_layerBins] = 0.;
- m_A[m_layerBins] = 0.;
- m_Z[m_layerBins] = 0.;
- m_Rho[m_layerBins] = 0.;
- m_elements[m_layerBins] = 0.;
+ m_thickness->at(m_layerBins) = 0.;
+ m_X0->at(m_layerBins) = 0.;
+ m_L0->at(m_layerBins) = 0.;
+ m_A->at(m_layerBins) = 0.;
+ m_Z->at(m_layerBins) = 0.;
+ m_Rho->at(m_layerBins) = 0.;
+ m_elements->at(m_layerBins) = 0.;
}
// set the bin Counter
- m_binCounter[m_layerBins] = bCounter ? (*bCounter)[bin1][bin0] : 1;
+ m_binCounter->at(m_layerBins) = bCounter ? (*bCounter)[bin1][bin0] : 1;
//
++bin0;
if (!bin1) ++m_layerBins0;
diff --git a/Tracking/TrkG4Components/TrkG4UserActions/python/TrkG4UserActionsConfigNew.py b/Tracking/TrkG4Components/TrkG4UserActions/python/TrkG4UserActionsConfigNew.py
new file mode 100644
index 0000000000000000000000000000000000000000..4c3962d5cd51dce5bb9deb1f71abe3aea7ec664d
--- /dev/null
+++ b/Tracking/TrkG4Components/TrkG4UserActions/python/TrkG4UserActionsConfigNew.py
@@ -0,0 +1,37 @@
+# Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
+
+"""Define methods to configure TrkG4UserActions"""
+
+def MaterialStepRecorderCfg(configFlags, name="G4UA::UserActionSvc.MaterialStepRecorderTool", **kwargs):
+ from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator
+ from AthenaConfiguration.ComponentFactory import CompFactory
+ result = ComponentAccumulator()
+ result.setPrivateTools(CompFactory.G4UA.MaterialStepRecorderTool(name, **kwargs))
+ return result
+
+
+def MaterialStepRecorderUserActionSvcCfg(configFlags, name="G4UA::MaterialStepRecorderUserActionSvc", **kwargs):
+ from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator
+ from AthenaConfiguration.ComponentFactory import CompFactory
+
+ result = ComponentAccumulator()
+
+ #Setting up the CA for the MaterialStepRecorder
+ actionAcc = ComponentAccumulator()
+ actions = []
+ actions += [actionAcc.popToolsAndMerge(MaterialStepRecorderCfg(configFlags))]
+ actionAcc.setPrivateTools(actions)
+ MaterialStepRecorderAction = result.popToolsAndMerge(actionAcc)
+
+ #Retrieving the default action list
+ from G4AtlasServices.G4AtlasUserActionConfigNew import getDefaultActions
+ defaultActions = result.popToolsAndMerge(getDefaultActions(configFlags))
+
+ #Adding LengthIntegrator to defaults
+ actionList = (defaultActions + MaterialStepRecorderAction)
+
+ #Setting up UserActionsService
+ kwargs.setdefault("UserActionTools",actionList)
+ result.addService(CompFactory.G4UA.UserActionSvc(name, **kwargs))
+
+ return result
diff --git a/Tracking/TrkG4Components/TrkG4UserActions/python/__init__.py b/Tracking/TrkG4Components/TrkG4UserActions/python/__init__.py
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/Tracking/TrkG4Components/TrkG4UserActions/scripts/RunGeantinoStepRecordingITk.py b/Tracking/TrkG4Components/TrkG4UserActions/scripts/RunGeantinoStepRecordingITk.py
new file mode 100644
index 0000000000000000000000000000000000000000..a5be31adac988d0ea4952660d37ea30270b96836
--- /dev/null
+++ b/Tracking/TrkG4Components/TrkG4UserActions/scripts/RunGeantinoStepRecordingITk.py
@@ -0,0 +1,141 @@
+#!/usr/bin/env python
+"""
+
+Run geantino processing for material step creation
+
+"""
+
+from argparse import ArgumentParser
+from AthenaCommon.Configurable import Configurable
+from AthenaCommon.Logging import log
+from AthenaConfiguration.AllConfigFlags import ConfigFlags
+from AthenaConfiguration.MainServicesConfig import MainServicesCfg
+from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator
+from AthenaConfiguration.ComponentFactory import CompFactory
+
+# Set up logging and new style config
+Configurable.configurableRun3Behavior = True
+
+# Argument parsing
+parser = ArgumentParser("RunGeantinoStepRecordingITk.py")
+parser.add_argument("detectors", metavar="detectors", type=str, nargs="*",
+ help="Specify the list of detectors")
+parser.add_argument("--simulate", default=True, action="store_true",
+ help="Run Simulation")
+parser.add_argument("--localgeo", default=False, action="store_true",
+ help="Use local geometry Xml files")
+parser.add_argument("-V", "--verboseAccumulators", default=False,
+ action="store_true",
+ help="Print full details of the AlgSequence")
+parser.add_argument("-S", "--verboseStoreGate", default=False,
+ action="store_true",
+ help="Dump the StoreGate(s) each event iteration")
+parser.add_argument("--maxEvents",default=10, type=int,
+ help="The number of events to run. 0 skips execution")
+parser.add_argument("--skipEvents",default=0, type=int,
+ help="The number of events to skip")
+parser.add_argument("--geometrytag",default="ATLAS-P2-ITK-24-00-00", type=str,
+ help="The geometry tag to use")
+parser.add_argument("--inputevntfile",
+ default="/cvmfs/atlas-nightlies.cern.ch/repo/data/data-art/InDetSLHC_Example/inputs/pgun_2M_10GeV_geantinos_Eta6_v2_EVNT.root",
+ help="The input EVNT file to use")
+parser.add_argument("--outputhitsfile",default="myHITS.pool.root", type=str,
+ help="The output HITS filename")
+args = parser.parse_args()
+
+
+# Some info about the job
+print("----GeantinoStepRecording for ITk geometry----")
+print()
+print("Using Geometry Tag: "+args.geometrytag)
+if args.localgeo:
+ print("...overridden by local Geometry Xml files")
+print("Input EVNT File:"+args.inputevntfile)
+if not args.detectors:
+ print("Running complete detector")
+else:
+ print("Running with: {}".format(", ".join(args.detectors)))
+print()
+
+# Configure
+if args.localgeo:
+ ConfigFlags.GeoModel.useLocalGeometry = True
+
+ConfigFlags.Input.Files = [args.inputevntfile]
+ConfigFlags.Output.HITSFileName = args.outputhitsfile
+
+ConfigFlags.GeoModel.AtlasVersion = args.geometrytag
+ConfigFlags.IOVDb.GlobalTag = "OFLCOND-SIM-00-00-00"
+ConfigFlags.GeoModel.Align.Dynamic = False
+
+ConfigFlags.Exec.SkipEvents = args.skipEvents
+
+from AthenaConfiguration.DetectorConfigFlags import setupDetectorsFromList
+detectors = args.detectors if 'detectors' in args and args.detectors else ['ITkPixel', 'ITkStrip']
+detectors.append('Bpipe') # always run with beam pipe
+setupDetectorsFromList(ConfigFlags, detectors, toggle_geometry=True)
+
+log.debug('Lock config flags now.')
+ConfigFlags.lock()
+
+# Construct our accumulator to run
+acc = MainServicesCfg(ConfigFlags)
+
+### setup dumping of additional information
+if args.verboseAccumulators:
+ acc.printConfig(withDetails=True)
+if args.verboseStoreGate:
+ acc.getService("StoreGateSvc").Dump = True
+
+log.debug('Dumping of ConfigFlags now.')
+ConfigFlags.dump()
+
+from AthenaPoolCnvSvc.PoolReadConfig import PoolReadCfg
+from AthenaPoolCnvSvc.PoolWriteConfig import PoolWriteCfg
+acc.merge(PoolReadCfg(ConfigFlags))
+acc.merge(PoolWriteCfg(ConfigFlags))
+
+# add BeamEffectsAlg
+from BeamEffects.BeamEffectsAlgConfig import BeamEffectsAlgCfg
+acc.merge(BeamEffectsAlgCfg(ConfigFlags))
+
+kwargs = {}
+
+svcName = "G4UA::MaterialStepRecorderUserActionSvc"
+from TrkG4UserActions.TrkG4UserActionsConfigNew import MaterialStepRecorderUserActionSvcCfg
+acc.merge(MaterialStepRecorderUserActionSvcCfg(ConfigFlags,svcName,**kwargs))
+kwargs.update(UserActionSvc=svcName)
+
+if args.simulate:
+ from G4AtlasAlg.G4AtlasAlgConfigNew import G4AtlasAlgCfg
+ acc.merge(G4AtlasAlgCfg(ConfigFlags, "ITkG4AtlasAlg", **kwargs))
+
+AthenaOutputStream=CompFactory.AthenaOutputStream
+AthenaOutputStreamTool=CompFactory.AthenaOutputStreamTool
+writingTool = AthenaOutputStreamTool( "MaterialStepCollectionStreamTool" )
+
+outputStream = AthenaOutputStream(name = "MaterialStepCollectionStream",
+ WritingTool = writingTool,
+ ItemList=['EventInfo#*', 'Trk::MaterialStepCollection#*'],
+ MetadataItemList = [ "EventStreamInfo#MaterialStepCollectionStream", "IOVMetaDataContainer#*" ],
+ OutputFile = "MaterialStepCollection.root")
+
+StoreGateSvc=CompFactory.StoreGateSvc
+acc.addService(StoreGateSvc("MetaDataStore"))
+outputStream.MetadataStore = acc.getService("MetaDataStore")
+
+MakeEventStreamInfo=CompFactory.MakeEventStreamInfo
+streamInfoTool = MakeEventStreamInfo( "MaterialStepCollectionStream_MakeEventStreamInfo" )
+streamInfoTool.Key = "MaterialStepCollectionStream"
+streamInfoTool.EventInfoKey = "EventInfo"
+outputStream.HelperTools.append(streamInfoTool)
+
+acc.addEventAlgo(outputStream)
+
+acc.printConfig(withDetails = True, summariseProps = True)
+
+acc.run(maxEvents=args.maxEvents)
+
+f=open("GeantinoStepRecordingITk.pkl","wb")
+acc.store(f)
+f.close()