diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/CMakeLists.txt b/Tracker/TrackerDigitization/FaserSCT_Digitization/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..d0965e8dc9ed2520312a1eab0fc8b8cb0ead9e8a
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/CMakeLists.txt
@@ -0,0 +1,55 @@
+################################################################################
+# Package: FaserSCT_Digitization
+################################################################################
+
+# Declare the package name:
+atlas_subdir( FaserSCT_Digitization )
+
+# Declare the package's dependencies:
+atlas_depends_on_subdirs( PUBLIC
+ #Commission/CommissionEvent
+ Control/AthenaBaseComps
+ Control/AthenaKernel
+ Control/PileUpTools
+ DetectorDescription/Identifier
+ Event/xAOD/xAODEventInfo
+ GaudiKernel
+ InnerDetector/InDetConditions/InDetCondTools
+ Tracker/TrackerDigitization/FaserSiDigitization
+ InnerDetector/InDetRawEvent/InDetRawData
+ Tracker/TrackerSimEvent
+ Simulation/HitManagement
+ PRIVATE
+ Control/StoreGate
+ Generators/GeneratorObjects
+ #InnerDetector/InDetConditions/InDetConditionsSummaryService
+ #InnerDetector/InDetConditions/SCT_ConditionsTools
+ InnerDetector/InDetConditions/SiPropertiesTool
+ Tracker/TrackerDetDescr/TrackerIdentifier
+ Tracker/TrackerDetDescr/TrackerReadoutGeometry
+ #InnerDetector/InDetDetDescr/SCT_ModuleDistortions
+ InnerDetector/InDetRawEvent/InDetSimData )
+
+# External dependencies:
+find_package( Boost COMPONENTS filesystem thread system )
+find_package( CLHEP )
+find_package( ROOT COMPONENTS Core Tree MathCore Hist RIO pthread )
+
+# Component(s) in the package:
+atlas_add_component( FaserSCT_Digitization
+ src/*.cxx
+ src/components/*.cxx
+ INCLUDE_DIRS ${ROOT_INCLUDE_DIRS} ${Boost_INCLUDE_DIRS} ${CLHEP_INCLUDE_DIRS}
+ LINK_LIBRARIES ${ROOT_LIBRARIES} ${Boost_LIBRARIES} ${CLHEP_LIBRARIES} AthenaBaseComps AthenaKernel PileUpToolsLib Identifier xAODEventInfo GaudiKernel FaserSiDigitization InDetRawData TrackerSimEvent HitManagement GeneratorObjects SiPropertiesToolLib TrackerIdentifier TrackerReadoutGeometry InDetSimData )
+# LINK_LIBRARIES ${ROOT_LIBRARIES} ${Boost_LIBRARIES} ${CLHEP_LIBRARIES} CommissionEvent AthenaBaseComps AthenaKernel PileUpToolsLib Identifier xAODEventInfo GaudiKernel SiDigitization InDetRawData InDetSimEvent HitManagement GeneratorObjects SiPropertiesToolLib InDetIdentifier InDetReadoutGeometry InDetSimData )
+
+#atlas_add_test( SCT_DigitizationMT_test
+# SCRIPT Digi_tf.py --inputHITSFile /cvmfs/atlas-nightlies.cern.ch/repo/data/data-art/DigitizationTests/HITS.04919495._001041.pool.root.1 --conditionsTag default:OFLCOND-RUN12-SDR-25 --digiSeedOffset1 170 --digiSeedOffset2 170 --geometryVersion ATLAS-R2-2015-03-01-00 --DataRunNumber 222525 --outputRDOFile mc15_2015_ttbar.RDO.pool.root --preInclude HITtoRDO:SimulationJobOptions/preInclude.SCTOnlyConfig.py,Digitization/ForceUseOfAlgorithms.py --postInclude Digitization/FixDataDependenciesForMT.py --skipEvents 0 --maxEvents 100 --athenaopts=--threads=10
+# PROPERTIES TIMEOUT 1200
+# ENVIRONMENT THREADS=10 )
+
+# Install files from the package:
+atlas_install_headers( FaserSCT_Digitization )
+atlas_install_python_modules( python/*.py )
+atlas_install_joboptions( share/*.py )
+
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ATLAS_CHECK_THREAD_SAFETY b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ATLAS_CHECK_THREAD_SAFETY
new file mode 100644
index 0000000000000000000000000000000000000000..c5030c6d4ab80ea56f53ce2ea52e49acec23ed91
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ATLAS_CHECK_THREAD_SAFETY
@@ -0,0 +1 @@
+Tracker/TrackerDigitization/FaserSCT_Digitization
\ No newline at end of file
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_Amp.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_Amp.h
new file mode 100644
index 0000000000000000000000000000000000000000..267840df70550cf4cc63de8a9f5d3de8a9c86fa2
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_Amp.h
@@ -0,0 +1,51 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+/**
+ * ISCT_Amp.h
+ * Header file for abstract base class ISCT_Amp
+ * (c) ATLAS Detector software
+ * Interface for the SiChargedDiode processor classes
+ * 23/08/2007 Kondo.Gnanvo@cern.ch, and others
+ */
+
+#ifndef FASERSIDIGITIZATION_ISCT_AMP_H
+#define FASERSIDIGITIZATION_ISCT_AMP_H
+
+#include "GaudiKernel/IAlgTool.h"
+#include "TrackerSimEvent/SiTotalCharge.h"
+
+#include <vector>
+
+static const InterfaceID IID_ISCT_Amp("ISCT_Amp", 1, 0);
+
+class ISCT_Amp : virtual public IAlgTool {
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+ public:
+ typedef SiTotalCharge::list_t list_t;
+
+ //Retrieve interface ID
+ static const InterfaceID& interfaceID() { return IID_ISCT_Amp; }
+
+ // Destructor:
+ virtual ~ISCT_Amp() {}
+
+ ///////////////////////////////////////////////////////////////////
+ // Pure virtual methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // process the collection of charged diodes
+ /** main purpose: CR-RC^3 response to a list of charges with times */
+ virtual float response(const list_t& Charges, const float timeOverThreshold) const =0;
+ virtual void response(const list_t& Charges, const float timeOverThreshold, std::vector<float>& resp) const =0;
+
+ /** Neighbour strip cross talk response strip to a list of charges with times */
+ virtual float crosstalk(const list_t& Charges, const float timeOverThreshold) const =0;
+ virtual void crosstalk(const list_t& Charges, const float timeOverThreshold, std::vector<float> &resp) const =0;
+};
+
+#endif // FASERSIDIGITIZATION_ISCT_AMP_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_FrontEnd.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_FrontEnd.h
new file mode 100644
index 0000000000000000000000000000000000000000..35a40ddbfa4f78825dae788954202f868c527732
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_FrontEnd.h
@@ -0,0 +1,42 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+/**
+ * ISCT_FrontEnd.h
+ * Header file for interface class for SCT_FrontEnd
+ * (c) ATLAS Detector software
+ */
+
+#ifndef FASERSCT_DIGITIZATION_ISCT_FRONTEND_H
+#define FASERSCT_DIGITIZATION_ISCT_FRONTEND_H
+
+//Inheritance
+#include "FaserSiDigitization/ISiChargedDiodesProcessorTool.h"
+#include "FaserSiDigitization/SiChargedDiode.h"
+
+//methods
+#include "Identifier/Identifier.h"
+class SiChargedDiodeCollection;
+namespace CLHEP {
+ class HepRandomEngine;
+}
+
+static const InterfaceID IID_ISCT_FrontEnd("ISCT_FrontEnd", 1, 0);
+
+class ISCT_FrontEnd : virtual public ISiChargedDiodesProcessorTool {
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+ public:
+
+ //** Retrieve interface ID */
+ static const InterfaceID& interfaceID() { return IID_ISCT_FrontEnd; }
+
+ //** Destructor: */
+ virtual ~ISCT_FrontEnd() {}
+
+};
+
+#endif // FASERSCT_DIGITIZATION_ISCT_FRONTEND_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_RandomDisabledCellGenerator.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_RandomDisabledCellGenerator.h
new file mode 100644
index 0000000000000000000000000000000000000000..d3c69bc251a9402ca1bb6e3a02cbd2c65c45a250
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_RandomDisabledCellGenerator.h
@@ -0,0 +1,42 @@
+// -*- C++ -*-
+
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+/**
+ * ISCT_RandomDisabledCellGenerator.h
+ * Header file for interface class for SCT_RandomDisabledCellGenerator
+ * (c) ATLAS Detector software
+ */
+
+#ifndef FASERSCT_DIGITIZATION_ISCT_RANDOMDISABLEDCELLGENERATOR_H
+#define FASERSCT_DIGITIZATION_ISCT_RANDOMDISABLEDCELLGENERATOR_H
+
+//Inheritance
+#include "FaserSiDigitization/ISiChargedDiodesProcessorTool.h"
+
+//methods
+class SiChargedDiodeCollection;
+namespace CLHEP {
+ class HepRandomEngine;
+}
+
+static const InterfaceID IID_ISCT_RandomDisabledCellGenerator("ISCT_RandomDisabledCellGenerator",1,0);
+
+class ISCT_RandomDisabledCellGenerator : virtual public ISiChargedDiodesProcessorTool {
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+ public:
+
+ //Retrieve interface ID
+ static const InterfaceID& interfaceID() { return IID_ISCT_RandomDisabledCellGenerator; }
+
+ // Destructor:
+ virtual ~ISCT_RandomDisabledCellGenerator() {}
+
+};
+
+#endif // FASERSCT_DIGITIZATION_ISCT_RANDOMDISABLEDCELLGENERATOR_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_SurfaceChargesGenerator.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_SurfaceChargesGenerator.h
new file mode 100644
index 0000000000000000000000000000000000000000..9188844d716ab40eab7ff8fe924a240a0dee6a64
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/FaserSCT_Digitization/ISCT_SurfaceChargesGenerator.h
@@ -0,0 +1,60 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+/**
+ * ISCT_SurfaceChargesGenerator.h
+ * Header file for abstract base class ISCT_SurfaceChargesGenerator
+ * (c) ATLAS Detector software
+ * Interface for the Sct_SurfaceCharge generator classes
+ * Version 23/08/2007 Kondo Gnanvo
+ */
+
+#ifndef FASERSCT_DIGITIZATION_ISCT_SURFACECHARGESGENERATOR_H
+#define FASERSCT_DIGITIZATION_ISCT_SURFACECHARGESGENERATOR_H
+
+// Input/output classes
+#include <list>
+#include "FaserSiDigitization/SiSurfaceCharge.h"
+#include "HitManagement/TimedHitPtr.h"
+
+#include "GaudiKernel/IAlgTool.h"
+
+class FaserSiHit;
+
+namespace TrackerDD {
+ class SiDetectorElement;
+}
+namespace CLHEP {
+ class HepRandomEngine;
+}
+
+class ISiSurfaceChargesInserter
+{
+ public:
+ virtual ~ISiSurfaceChargesInserter() {}
+ virtual void operator() (const SiSurfaceCharge& scharge) const = 0;
+};
+
+static const InterfaceID IID_ISCT_SurfaceChargesGenerator("ISCT_SurfaceChargesGenerator",1,0);
+
+class ISCT_SurfaceChargesGenerator : virtual public IAlgTool {
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+ public:
+
+ //Retrieve interface ID
+ static const InterfaceID& interfaceID() { return IID_ISCT_SurfaceChargesGenerator; }
+
+ // Destructor:
+ virtual ~ISCT_SurfaceChargesGenerator() {}
+
+ virtual void process(const TrackerDD::SiDetectorElement* ele,
+ const TimedHitPtr<FaserSiHit>& phit,
+ const ISiSurfaceChargesInserter& inserter, CLHEP::HepRandomEngine * rndmEngine) const =0;
+ virtual void setFixedTime(float fixedTime) =0;
+};
+
+#endif // FASERSCT_DIGITIZATION_ISCT_SURFACECHARGESGENERATOR_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/python/FaserSCT_DigitizationConfig.py b/Tracker/TrackerDigitization/FaserSCT_Digitization/python/FaserSCT_DigitizationConfig.py
new file mode 100644
index 0000000000000000000000000000000000000000..42b92ef195088b1390cbae4fd21e455f02d1151b
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/python/FaserSCT_DigitizationConfig.py
@@ -0,0 +1,297 @@
+# Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+
+from AthenaCommon import CfgMgr
+# The earliest bunch crossing time for which interactions will be sent
+# to the SCT Digitization code.
+def SCT_FirstXing():
+ return -50
+# The latest bunch crossing time for which interactions will be sent
+# to the SCT Digitization code.
+def SCT_LastXing():
+ return 25
+
+
+
+######################################################################################
+def getSCT_RandomDisabledCellGenerator(name="SCT_RandomDisabledCellGenerator", **kwargs):
+ kwargs.setdefault("TotalBadChannels", 0.01)
+ from FaserSCT_Digitization.FaserSCT_DigitizationConf import SCT_RandomDisabledCellGenerator
+ return SCT_RandomDisabledCellGenerator(name, **kwargs)
+
+
+######################################################################################
+def getSCT_Amp(name="SCT_Amp", **kwargs):
+ kwargs.setdefault("CrossFactor2sides", 0.1)
+ kwargs.setdefault("CrossFactorBack", 0.07)
+ kwargs.setdefault("PeakTime", 21)
+ kwargs.setdefault("deltaT", 1.0)
+ kwargs.setdefault("Tmin", -25.0)
+ kwargs.setdefault("Tmax", 150.0)
+ from FaserSCT_Digitization.FaserSCT_DigitizationConf import SCT_Amp
+ return SCT_Amp(name, **kwargs)
+
+
+######################################################################################
+def getSCT_SurfaceChargesGenerator(name="SCT_SurfaceChargesGenerator", **kwargs):
+ ## Set up services used by SCT_SurfaceChargesGenerator
+ # Set up SCT_DCSConditiosnTool
+ from SCT_ConditionsTools.SCT_DCSConditionsToolSetup import SCT_DCSConditionsToolSetup
+ sct_DCSConditionsToolSetup = SCT_DCSConditionsToolSetup()
+ sct_DCSConditionsToolSetup.setup()
+ # Set up SCT_SiliconConditionsTool
+ from SCT_ConditionsTools.SCT_SiliconConditionsToolSetup import SCT_SiliconConditionsToolSetup
+ sct_SiliconConditionsToolSetup = SCT_SiliconConditionsToolSetup()
+ sct_SiliconConditionsToolSetup.setDcsTool(sct_DCSConditionsToolSetup.getTool())
+ sct_SiliconConditionsToolSetup.setup()
+ # Set up SCT_SiPropertiesTool
+ from SiPropertiesTool.SCT_SiPropertiesToolSetup import SCT_SiPropertiesToolSetup
+ sct_SiPropertiesToolSetup = SCT_SiPropertiesToolSetup()
+ sct_SiPropertiesToolSetup.setSiliconTool(sct_SiliconConditionsToolSetup.getTool())
+ sct_SiPropertiesToolSetup.setup()
+ ## Charge trapping tool - used by SCT_SurfaceChargesGenerator
+ from AthenaCommon.AppMgr import ToolSvc
+ ## SiLorentzAngleTool for SCT_SurfaceChargesGenerator
+ from SiLorentzAngleTool.SCTLorentzAngleToolSetup import SCTLorentzAngleToolSetup
+ sctLorentzAngleToolSetup = SCTLorentzAngleToolSetup()
+
+ kwargs.setdefault("FixedTime", -999)
+ kwargs.setdefault("SubtractTime", -999)
+ kwargs.setdefault("SurfaceDriftTime", 10)
+ kwargs.setdefault("NumberOfCharges", 1)
+ kwargs.setdefault("SmallStepLength", 5)
+ kwargs.setdefault("DepletionVoltage", 70)
+ kwargs.setdefault("BiasVoltage", 150)
+ kwargs.setdefault("SiPropertiesTool", sct_SiPropertiesToolSetup.getTool())
+ kwargs.setdefault("LorentzAngleTool", sctLorentzAngleToolSetup.SCTLorentzAngleTool)
+ from AthenaCommon.GlobalFlags import globalflags
+ kwargs.setdefault("isOverlay", globalflags.isOverlay())
+
+ # kwargs.setdefault("doTrapping", True) # ATL-INDET-INT-2016-019
+
+ from Digitization.DigitizationFlags import digitizationFlags
+ if 'doDetailedSurfChargesGen' in digitizationFlags.experimentalDigi():
+ kwargs.setdefault("ChargeDriftModel", 1)
+ kwargs.setdefault("EFieldModel", 2)
+ kwargs.setdefault("MagneticField", -2.0)
+ kwargs.setdefault("SensorTemperature", 273.15)
+ kwargs.setdefault("TransportTimeStep", 0.25)
+ kwargs.setdefault("TransportTimeMax", 25.0)
+ from FaserSCT_Digitization.FaserSCT_DigitizationConf import SCT_DetailedSurfaceChargesGenerator
+ return SCT_DetailedSurfaceChargesGenerator(name, **kwargs)
+ else:
+ from SCT_ConditionsTools.SCT_ConditionsToolsConf import SCT_RadDamageSummaryTool
+ kwargs.setdefault("RadDamageSummaryTool", SCT_RadDamageSummaryTool(name = "InDetSCT_RadDamageSummaryTool"))
+ from FaserSCT_Digitization.FaserSCT_DigitizationConf import SCT_SurfaceChargesGenerator
+ return SCT_SurfaceChargesGenerator(name, **kwargs)
+
+######################################################################################
+def getSCT_FrontEnd(name="SCT_FrontEnd", **kwargs):
+ from Digitization.DigitizationFlags import digitizationFlags
+ #Setup noise treament in SCT_FrontEnd
+ # To set the mean noise values for the different module types
+ # Default values set at 0 degrees, plus/minus ~5 enc per plus/minus degree
+ kwargs.setdefault("NoiseBarrel", 1500.0)
+ kwargs.setdefault("NoiseBarrel3", 1541.0)
+ kwargs.setdefault("NoiseInners", 1090.0)
+ kwargs.setdefault("NoiseMiddles", 1557.0)
+ kwargs.setdefault("NoiseShortMiddles", 940.0)
+ kwargs.setdefault("NoiseOuters", 1618.0)
+ kwargs.setdefault("NOBarrel", 1.5e-5)
+ kwargs.setdefault("NOBarrel3", 2.1e-5)
+ kwargs.setdefault("NOInners", 5.0e-9)
+ kwargs.setdefault("NOMiddles", 2.7e-5)
+ kwargs.setdefault("NOShortMiddles", 2.0e-9)
+ kwargs.setdefault("NOOuters", 3.5e-5)
+ # If noise is turned off:
+ if not digitizationFlags.doInDetNoise.get_Value():
+ ###kwargs.setdefault("OnlyHitElements", True)
+ print('SCT_Digitization:::: Turned off Noise in SCT_FrontEnd')
+ kwargs.setdefault("NoiseOn", False)
+ kwargs.setdefault("AnalogueNoiseOn", False)
+ else:
+ kwargs.setdefault("NoiseOn", True)
+ kwargs.setdefault("AnalogueNoiseOn", True)
+ # In overlay MC, only analogue noise is on. Noise hits are not added.
+ from AthenaCommon.GlobalFlags import globalflags
+ if globalflags.isOverlay() and globalflags.DataSource == 'geant4':
+ kwargs["NoiseOn"] = False
+ kwargs["AnalogueNoiseOn"] = True
+ # Use Calibration data from Conditions DB, still for testing purposes only
+ kwargs.setdefault("UseCalibData", True)
+
+ # Setup the ReadCalibChip folders and Svc
+ from SCT_ConditionsTools.SCT_ReadCalibChipDataToolSetup import SCT_ReadCalibChipDataToolSetup
+ sct_ReadCalibChipDataToolSetup = SCT_ReadCalibChipDataToolSetup()
+ sct_ReadCalibChipDataToolSetup.setup()
+ kwargs.setdefault("SCT_ReadCalibChipDataTool", sct_ReadCalibChipDataToolSetup.getTool())
+ # DataCompressionMode: 1 is level mode x1x (default), 2 is edge mode 01x, 3 is expanded any hit xxx
+ from AthenaCommon.BeamFlags import jobproperties
+ if digitizationFlags.PileUpPremixing:
+ kwargs.setdefault("DataCompressionMode", 3)
+ elif globalflags.isOverlay() and globalflags.DataSource == 'geant4':
+ kwargs.setdefault("DataCompressionMode", 2)
+ elif (jobproperties.Beam.bunchSpacing() <= 50):
+ kwargs.setdefault("DataCompressionMode", 1)
+ else:
+ kwargs.setdefault("DataCompressionMode", 3)
+ # DataReadOutMode: 0 is condensed mode and 1 is expanded mode
+ if globalflags.isOverlay() and globalflags.DataSource == 'geant4':
+ kwargs.setdefault("DataReadOutMode", 0)
+ else:
+ kwargs.setdefault("DataReadOutMode", 1)
+ from FaserSCT_Digitization.FaserSCT_DigitizationConf import SCT_FrontEnd
+ return SCT_FrontEnd(name, **kwargs)
+
+######################################################################################
+def getPileupSCT_FrontEnd(name="PileupSCT_FrontEnd", **kwargs):
+
+ kwargs.setdefault("NoiseBarrel", 0.0)
+ kwargs.setdefault("NoiseBarrel3", 0.0)
+ kwargs.setdefault("NoiseInners", 0.0)
+ kwargs.setdefault("NoiseMiddles", 0.0)
+ kwargs.setdefault("NoiseShortMiddles", 0.0)
+ kwargs.setdefault("NoiseOuters", 0.0)
+ kwargs.setdefault("NOBarrel", 0.0)
+ kwargs.setdefault("NOBarrel3", 0.0)
+ kwargs.setdefault("NOInners", 0.0)
+ kwargs.setdefault("NOMiddles", 0.0)
+ kwargs.setdefault("NOShortMiddles", 0.0)
+ kwargs.setdefault("NOOuters", 0.0)
+ kwargs.setdefault("NoiseOn", False)
+
+ return getSCT_FrontEnd(name, **kwargs)
+
+######################################################################################
+
+
+def commonSCT_DigitizationConfig(name,**kwargs):
+
+ from Digitization.DigitizationFlags import digitizationFlags
+ # If noise is turned off:
+ if not digitizationFlags.doInDetNoise.get_Value():
+ kwargs.setdefault("OnlyHitElements", True)
+
+ kwargs.setdefault("InputObjectName", "SCT_Hits")
+ kwargs.setdefault("EnableHits", True)
+ kwargs.setdefault("BarrelOnly", False)
+
+ # Use of random disabled cells
+ #kwargs.setdefault("RandomDisabledCells", True)
+
+ # Set FixedTime for cosmics for use in SurfaceChargesGenerator
+ from AthenaCommon.BeamFlags import jobproperties
+ if jobproperties.Beam.beamType == "cosmics" :
+ kwargs.setdefault("CosmicsRun", True)
+ kwargs.setdefault("FixedTime", 10)
+
+ # write out SCT1_RawData
+ #kwargs.setdefault("WriteSCT1_RawData", False)
+
+ if digitizationFlags.doXingByXingPileUp():
+ kwargs.setdefault("FirstXing", SCT_FirstXing())
+ kwargs.setdefault("LastXing", SCT_LastXing() )
+
+ from AthenaCommon import CfgMgr
+ return CfgMgr.SCT_DigitizationTool(name,**kwargs)
+# else:
+# from AthenaCommon import CfgMgr
+# return CfgMgr.SCT_Digitization(name, **kwargs)
+
+######################################################################################
+
+def SCT_DigitizationTool(name="SCT_DigitizationTool", **kwargs):
+ from Digitization.DigitizationFlags import digitizationFlags
+ if digitizationFlags.PileUpPremixing and 'OverlayMT' in digitizationFlags.experimentalDigi():
+ from OverlayCommonAlgs.OverlayFlags import overlayFlags
+ kwargs.setdefault("OutputObjectName", overlayFlags.bkgPrefix() + "SCT_RDOs")
+ kwargs.setdefault("OutputSDOName", overlayFlags.bkgPrefix() + "SCT_SDO_Map")
+ else:
+ kwargs.setdefault("OutputObjectName", "SCT_RDOs")
+ kwargs.setdefault("OutputSDOName", "SCT_SDO_Map")
+
+ kwargs.setdefault("HardScatterSplittingMode", 0)
+ return commonSCT_DigitizationConfig(name,**kwargs)
+
+######################################################################################
+
+def SCT_GeantinoTruthDigitizationTool(name="SCT_GeantinoTruthDigitizationTool", **kwargs):
+ kwargs.setdefault("ParticleBarcodeVeto", 0)
+ return SCT_DigitizationTool(name,**kwargs)
+
+######################################################################################
+
+def SCT_DigitizationToolHS(name="SCT_DigitizationToolHS",**kwargs):
+ kwargs.setdefault("OutputObjectName", "SCT_RDOs")
+ kwargs.setdefault("OutputSDOName", "SCT_SDO_Map")
+ kwargs.setdefault("HardScatterSplittingMode", 1)
+ return commonSCT_DigitizationConfig(name,**kwargs)
+
+######################################################################################
+
+def SCT_DigitizationToolPU(name="SCT_DigitizationToolPU",**kwargs):
+ kwargs.setdefault("OutputObjectName", "SCT_PU_RDOs")
+ kwargs.setdefault("OutputSDOName", "SCT_PU_SDO_Map")
+ kwargs.setdefault("HardScatterSplittingMode", 2)
+ return commonSCT_DigitizationConfig(name,**kwargs)
+
+######################################################################################
+
+def SCT_DigitizationToolSplitNoMergePU(name="SCT_DigitizationToolSplitNoMergePU",**kwargs):
+
+ kwargs.setdefault("InputObjectName", "PileupSCT_Hits")
+ kwargs.setdefault("HardScatterSplittingMode", 0)
+ kwargs.setdefault("OutputObjectName", "SCT_PU_RDOs")
+ kwargs.setdefault("OutputSDOName", "SCT_PU_SDO_Map")
+ kwargs.setdefault("OnlyHitElements", True)
+ kwargs.setdefault("FrontEnd", "PileupSCT_FrontEnd")
+
+ return commonSCT_DigitizationConfig(name,**kwargs)
+
+######################################################################################
+
+def SCT_OverlayDigitizationTool(name="SCT_OverlayDigitizationTool",**kwargs):
+ from OverlayCommonAlgs.OverlayFlags import overlayFlags
+ if overlayFlags.isOverlayMT():
+ kwargs.setdefault("OnlyUseContainerName", False)
+ kwargs.setdefault("OutputObjectName", overlayFlags.sigPrefix() + "SCT_RDOs")
+ kwargs.setdefault("OutputSDOName", overlayFlags.sigPrefix() + "SCT_SDO_Map")
+ else:
+ kwargs.setdefault("OutputObjectName", overlayFlags.evtStore() + "+SCT_RDOs")
+ kwargs.setdefault("OutputSDOName", overlayFlags.evtStore() + "+SCT_SDO_Map")
+ kwargs.setdefault("HardScatterSplittingMode", 0)
+ return commonSCT_DigitizationConfig(name,**kwargs)
+
+######################################################################################
+
+def getSiliconRange(name="SiliconRange" , **kwargs):
+ #this is the time of the xing in ns
+ kwargs.setdefault('FirstXing', SCT_FirstXing() )
+ kwargs.setdefault('LastXing', SCT_LastXing() )
+ kwargs.setdefault('CacheRefreshFrequency', 1.0 ) #default 0 no dataproxy reset
+ kwargs.setdefault('ItemList', ["SiHitCollection#SCT_Hits"] )
+ return CfgMgr.PileUpXingFolder(name, **kwargs)
+
+######################################################################################
+
+def SCT_DigitizationHS(name="SCT_DigitizationHS",**kwargs):
+ kwargs.setdefault("DigitizationTool", "SCT_DigitizationToolHS")
+ from FaserSCT_Digitization.FaserSCT_DigitizationConf import SCT_Digitization
+ return SCT_Digitization(name,**kwargs)
+
+######################################################################################
+
+def SCT_DigitizationPU(name="SCT_DigitizationPU",**kwargs):
+ kwargs.setdefault("DigitizationTool", "SCT_DigitizationToolPU")
+ return CfgMgr.SCT_Digitization(name,**kwargs)
+
+######################################################################################
+
+def SCT_OverlayDigitization(name="SCT_OverlayDigitization",**kwargs):
+ kwargs.setdefault("DigitizationTool", "SCT_OverlayDigitizationTool")
+ # Multi-threading settinggs
+ from AthenaCommon.ConcurrencyFlags import jobproperties as concurrencyProps
+ is_hive = (concurrencyProps.ConcurrencyFlags.NumThreads() > 0)
+ if is_hive:
+ kwargs.setdefault('Cardinality', concurrencyProps.ConcurrencyFlags.NumThreads())
+
+ return CfgMgr.SCT_Digitization(name,**kwargs)
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/python/FaserSCT_DigitizationConfigDb.py b/Tracker/TrackerDigitization/FaserSCT_Digitization/python/FaserSCT_DigitizationConfigDb.py
new file mode 100644
index 0000000000000000000000000000000000000000..11d8480f3c69a0e2c8e0b93a1844e522b24443e1
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/python/FaserSCT_DigitizationConfigDb.py
@@ -0,0 +1,18 @@
+# Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+
+from AthenaCommon.CfgGetter import addTool,addService,addAlgorithm
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.SCT_DigitizationTool" , "SCT_DigitizationTool")
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.SCT_GeantinoTruthDigitizationTool" , "SCT_GeantinoTruthDigitizationTool")
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.SCT_DigitizationToolHS" , "SCT_DigitizationToolHS")
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.SCT_DigitizationToolPU" , "SCT_DigitizationToolPU")
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.SCT_DigitizationToolSplitNoMergePU", "SCT_DigitizationToolSplitNoMergePU")
+addAlgorithm("FaserSCT_Digitization.FaserSCT_DigitizationConfig.SCT_DigitizationHS" , "SCT_DigitizationHS")
+addAlgorithm("FaserSCT_Digitization.FaserSCT_DigitizationConfig.SCT_DigitizationPU" , "SCT_DigitizationPU")
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.getSiliconRange" , "SiliconRange" )
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.getSCT_RandomDisabledCellGenerator", "SCT_RandomDisabledCellGenerator")
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.getSCT_Amp", "SCT_Amp" )
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.getSCT_FrontEnd" , "SCT_FrontEnd" )
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.getPileupSCT_FrontEnd" , "PileupSCT_FrontEnd" )
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.getSCT_SurfaceChargesGenerator", "SCT_SurfaceChargesGenerator" )
+addTool("FaserSCT_Digitization.FaserSCT_DigitizationConfig.SCT_OverlayDigitizationTool", "SCT_OverlayDigitizationTool")
+addAlgorithm("FaserSCT_Digitization.FaserSCT_DigitizationConfig.SCT_OverlayDigitization", "SCT_OverlayDigitization")
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/python/FaserSCT_DigitizationConfigNew.py b/Tracker/TrackerDigitization/FaserSCT_Digitization/python/FaserSCT_DigitizationConfigNew.py
new file mode 100644
index 0000000000000000000000000000000000000000..92b74cad966602cd7e9b9046cba428534805f7fd
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/python/FaserSCT_DigitizationConfigNew.py
@@ -0,0 +1,305 @@
+"""Define methods to construct configured SCT Digitization tools and algorithms
+
+Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+"""
+from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator
+from AthenaConfiguration.ComponentFactory import CompFactory
+from AthenaCommon.Logging import logging
+from FaserSCT_Digitization.FaserSCT_DigitizationConf import (
+ SCT_RandomDisabledCellGenerator,
+ SCT_Amp,
+ SCT_SurfaceChargesGenerator,
+ SCT_FrontEnd,
+ SCT_DigitizationTool,
+ SCT_Digitization,
+)
+PileUpXingFolder=CompFactory.PileUpXingFolder
+from OutputStreamAthenaPool.OutputStreamConfig import OutputStreamCfg
+#SCT_RadDamageSummaryTool=CompFactory.SCT_RadDamageSummaryTool
+from FaserSCT_GeoModel.FaserSCT_GeoModelConfig import FaserSCT_GeometryCfg
+from SCT_ConditionsTools.SCT_DCSConditionsConfig import SCT_DCSConditionsCfg
+from SCT_ConditionsTools.SCT_SiliconConditionsConfig import SCT_SiliconConditionsToolCfg, SCT_SiliconConditionsCfg
+from SCT_ConditionsTools.SCT_ReadCalibChipDataConfig import SCT_ReadCalibChipDataCfg
+from SiPropertiesTool.SCT_SiPropertiesConfig import SCT_SiPropertiesCfg
+from SiLorentzAngleTool.SCT_LorentzAngleConfig import SCT_LorentzAngleCfg
+from Digitization.TruthDigitizationOutputConfig import TruthDigitizationOutputCfg
+from Digitization.PileUpToolsConfig import PileUpToolsCfg
+
+
+# The earliest and last bunch crossing times for which interactions will be sent
+# to the SCT Digitization code
+def SCT_FirstXing():
+ return -50
+
+
+def SCT_LastXing():
+ return 25
+
+
+def SCT_DigitizationCommonCfg(flags, name="SCT_DigitizationToolCommon", **kwargs):
+ """Return ComponentAccumulator with common SCT digitization tool config"""
+ acc = FaserSCT_GeometryCfg(flags)
+ if not flags.Digitization.DoInnerDetectorNoise:
+ kwargs.setdefault("OnlyHitElements", True)
+ kwargs.setdefault("InputObjectName", "SCT_Hits")
+ kwargs.setdefault("EnableHits", True)
+ #kwargs.setdefault("BarrelOnly", False)
+ # Set FixedTime for cosmics for use in SurfaceChargesGenerator
+ if flags.Beam.Type == "cosmics":
+ kwargs.setdefault("CosmicsRun", True)
+ kwargs.setdefault("FixedTime", 10)
+ if flags.Digitization.DoXingByXingPileUp:
+ kwargs.setdefault("FirstXing", SCT_FirstXing())
+ kwargs.setdefault("LastXing", SCT_LastXing() )
+ tool = SCT_DigitizationTool(name, **kwargs)
+ # attach ToolHandles
+ tool.FrontEnd = acc.popToolsAndMerge(SCT_FrontEndCfg(flags))
+ tool.SurfaceChargesGenerator = acc.popToolsAndMerge(SCT_SurfaceChargesGeneratorCfg(flags))
+ tool.RandomDisabledCellGenerator = SCT_RandomDisabledCellGeneratorCfg(flags)
+ acc.setPrivateTools(tool)
+ return acc
+
+
+def SCT_DigitizationToolCfg(flags, name="SCT_DigitizationTool", **kwargs):
+ """Return ComponentAccumulator with configured SCT digitization tool"""
+ if flags.Digitization.PileUpPremixing:
+ kwargs.setdefault("OutputObjectName", flags.Overlay.BkgPrefix + "SCT_RDOs")
+ kwargs.setdefault("OutputSDOName", flags.Overlay.BkgPrefix + "SCT_SDO_Map")
+ else:
+ kwargs.setdefault("OutputObjectName", "SCT_RDOs")
+ kwargs.setdefault("OutputSDOName", "SCT_SDO_Map")
+ #kwargs.setdefault("HardScatterSplittingMode", 0)
+ return SCT_DigitizationCommonCfg(flags, name, **kwargs)
+
+
+def SCT_DigitizationHSToolCfg(flags, name="SCT_DigitizationHSTool", **kwargs):
+ """Return ComponentAccumulator with hard scatter configured SCT digitization tool"""
+ kwargs.setdefault("OutputObjectName", "SCT_RDOs")
+ kwargs.setdefault("OutputSDOName", "SCT_SDO_Map")
+ #kwargs.setdefault("HardScatterSplittingMode", 1)
+ return SCT_DigitizationCommonCfg(flags, name, **kwargs)
+
+
+def SCT_DigitizationPUToolCfg(flags, name="SCT_DigitizationPUTool",**kwargs):
+ """Return ComponentAccumulator with pileup configured SCT digitization tool"""
+ kwargs.setdefault("OutputObjectName", "SCT_PU_RDOs")
+ kwargs.setdefault("OutputSDOName", "SCT_PU_SDO_Map")
+ #kwargs.setdefault("HardScatterSplittingMode", 2)
+ return SCT_DigitizationCommonCfg(flags, name, **kwargs)
+
+
+def SCT_OverlayDigitizationToolCfg(flags, name="SCT_OverlayDigitizationTool",**kwargs):
+ """Return ComponentAccumulator with overlay configured SCT digitization tool"""
+ acc = ComponentAccumulator()
+ kwargs.setdefault("OnlyUseContainerName", False)
+ kwargs.setdefault("OutputObjectName", "StoreGateSvc+" + flags.Overlay.SigPrefix + "SCT_RDOs")
+ kwargs.setdefault("OutputSDOName", "StoreGateSvc+" + flags.Overlay.SigPrefix + "SCT_SDO_Map")
+ #kwargs.setdefault("HardScatterSplittingMode", 0)
+ tool = acc.popToolsAndMerge(SCT_DigitizationCommonCfg(flags, name, **kwargs))
+ acc.setPrivateTools(tool)
+ return acc
+
+
+def SCT_DigitizationToolSplitNoMergePUCfg(flags, name="SCT_DigitizationToolSplitNoMergePU",**kwargs):
+ """Return ComponentAccumulator with merged pileup configured SCT digitization tool"""
+ kwargs.setdefault("InputObjectName", "PileupSCT_Hits")
+ #kwargs.setdefault("HardScatterSplittingMode", 0)
+ kwargs.setdefault("OutputObjectName", "SCT_PU_RDOs")
+ kwargs.setdefault("OutputSDOName", "SCT_PU_SDO_Map")
+ kwargs.setdefault("OnlyHitElements", True)
+ kwargs.setdefault("FrontEnd", "PileupSCT_FrontEnd")
+ return SCT_DigitizationCommonCfg(flags, name, **kwargs)
+
+
+def SCT_DigitizationToolGeantinoTruthCfg(flags, name="SCT_GeantinoTruthDigitizationTool", **kwargs):
+ """Return Geantino truth configured digitization tool"""
+ kwargs.setdefault("ParticleBarcodeVeto", 0)
+ return SCT_DigitizationToolCfg(flags, name, **kwargs)
+
+
+def SCT_RandomDisabledCellGeneratorCfg(flags, name="SCT_RandomDisabledCellGenerator", **kwargs):
+ """Return configured random cell disabling tool"""
+ kwargs.setdefault("TotalBadChannels", 0.01)
+ return SCT_RandomDisabledCellGenerator(name, **kwargs)
+
+
+def SCT_AmpCfg(flags, name="SCT_Amp", **kwargs):
+ """Return configured amplifier and shaper tool"""
+ kwargs.setdefault("CrossFactor2sides", 0.1)
+ kwargs.setdefault("CrossFactorBack", 0.07)
+ kwargs.setdefault("PeakTime", 21)
+ kwargs.setdefault("deltaT", 1.0)
+ kwargs.setdefault("Tmin", -25.0)
+ kwargs.setdefault("Tmax", 150.0)
+ return SCT_Amp(name, **kwargs)
+
+
+def SCT_SurfaceChargesGeneratorCfg(flags, name="SCT_SurfaceChargesGenerator", **kwargs):
+ """Return ComponentAccumulator with configured surface charges tool"""
+ acc = ComponentAccumulator()
+ kwargs.setdefault("FixedTime", -999)
+ kwargs.setdefault("SubtractTime", -999)
+ kwargs.setdefault("SurfaceDriftTime", 10)
+ kwargs.setdefault("NumberOfCharges", 1)
+ kwargs.setdefault("SmallStepLength", 5)
+ kwargs.setdefault("DepletionVoltage", 70)
+ kwargs.setdefault("BiasVoltage", 150)
+ kwargs.setdefault("isOverlay", flags.Detector.Overlay)
+ # kwargs.setdefault("doTrapping", True) # ATL-INDET-INT-2016-019
+ # experimental SCT_DetailedSurfaceChargesGenerator config dropped here
+ tool = SCT_SurfaceChargesGenerator(name, **kwargs)
+ #tool.RadDamageSummaryTool = SCT_RadDamageSummaryTool()
+ DCSCondTool = acc.popToolsAndMerge(SCT_DCSConditionsCfg(flags))
+ SiliCondTool = SCT_SiliconConditionsToolCfg(flags)
+ SiliCondAcc = SCT_SiliconConditionsCfg(flags, DCSConditionsTool=DCSCondTool)
+ SiliPropsAcc = SCT_SiPropertiesCfg(flags, SiConditionsTool=SiliCondTool)
+ acc.merge(SiliCondAcc)
+ #tool.SiConditionsTool = SiliCondTool
+ tool.SiPropertiesTool = acc.popToolsAndMerge(SiliPropsAcc)
+ tool.LorentzAngleTool = acc.popToolsAndMerge(SCT_LorentzAngleCfg(flags))
+ acc.setPrivateTools(tool)
+ return acc
+
+
+def SCT_FrontEndCfg(flags, name="SCT_FrontEnd", **kwargs):
+ """Return ComponentAccumulator with configured front-end electronics tool"""
+ # Setup noise treament in SCT_FrontEnd
+ # To set the mean noise values for the different module types
+ # Default values set at 0 degrees, plus/minus ~5 enc per plus/minus degree
+ kwargs.setdefault("NoiseBarrel", 1500.0)
+ kwargs.setdefault("NoiseBarrel3", 1541.0)
+ kwargs.setdefault("NoiseInners", 1090.0)
+ kwargs.setdefault("NoiseMiddles", 1557.0)
+ kwargs.setdefault("NoiseShortMiddles", 940.0)
+ kwargs.setdefault("NoiseOuters", 1618.0)
+ kwargs.setdefault("NOBarrel", 1.5e-5)
+ kwargs.setdefault("NOBarrel3", 2.1e-5)
+ kwargs.setdefault("NOInners", 5.0e-9)
+ kwargs.setdefault("NOMiddles", 2.7e-5)
+ kwargs.setdefault("NOShortMiddles", 2.0e-9)
+ kwargs.setdefault("NOOuters", 3.5e-5)
+ if not flags.Digitization.DoInnerDetectorNoise:
+ log = logging.getLogger("SCT_FrontEndCfg")
+ log.info("SCT_Digitization:::: Turned off Noise in SCT_FrontEnd")
+ kwargs.setdefault("NoiseOn", False)
+ kwargs.setdefault("AnalogueNoiseOn", False)
+ else:
+ kwargs.setdefault("NoiseOn", True)
+ kwargs.setdefault("AnalogueNoiseOn", True)
+ # In overlay MC, only analogue noise is on. Noise hits are not added.
+ if flags.Detector.Overlay and flags.Input.isMC:
+ kwargs["NoiseOn"] = False
+ kwargs["AnalogueNoiseOn"] = True
+ # Use Calibration data from Conditions DB, still for testing purposes only
+ kwargs.setdefault("UseCalibData", True)
+ # Setup the ReadCalibChip folders and Svc
+ acc = SCT_ReadCalibChipDataCfg(flags)
+ garbage = acc.popPrivateTools()
+ #kwargs.setdefault("SCT_ReadCalibChipDataTool", acc.popPrivateTools())
+ # DataCompressionMode: 1 is level mode x1x (default), 2 is edge mode 01x, 3 is expanded any hit xxx
+ if flags.Digitization.PileUpPremixing:
+ kwargs.setdefault("DataCompressionMode", 3)
+ elif flags.Detector.Overlay and flags.Input.isMC:
+ kwargs.setdefault("DataCompressionMode", 2)
+ elif flags.Beam.BunchSpacing <= 50:
+ kwargs.setdefault("DataCompressionMode", 1)
+ else:
+ kwargs.setdefault("DataCompressionMode", 3)
+ # DataReadOutMode: 0 is condensed mode and 1 is expanded mode
+ if flags.Detector.Overlay and flags.Input.isMC:
+ kwargs.setdefault("DataReadOutMode", 0)
+ else:
+ kwargs.setdefault("DataReadOutMode", 1)
+ acc.setPrivateTools(SCT_FrontEnd(name, **kwargs))
+ return acc
+
+
+def SCT_FrontEndPileupCfg(flags, name="PileupSCT_FrontEnd", **kwargs):
+ """Return ComponentAccumulator with pileup-configured front-end electronics tool"""
+ kwargs.setdefault("NoiseBarrel", 0.0)
+ kwargs.setdefault("NoiseBarrel3", 0.0)
+ kwargs.setdefault("NoiseInners", 0.0)
+ kwargs.setdefault("NoiseMiddles", 0.0)
+ kwargs.setdefault("NoiseShortMiddles", 0.0)
+ kwargs.setdefault("NoiseOuters", 0.0)
+ kwargs.setdefault("NOBarrel", 0.0)
+ kwargs.setdefault("NOBarrel3", 0.0)
+ kwargs.setdefault("NOInners", 0.0)
+ kwargs.setdefault("NOMiddles", 0.0)
+ kwargs.setdefault("NOShortMiddles", 0.0)
+ kwargs.setdefault("NOOuters", 0.0)
+ kwargs.setdefault("NoiseOn", False)
+ return SCT_FrontEndCfg(flags, name, **kwargs)
+
+def SCT_RangeCfg(flags, name="SiliconRange", **kwargs):
+ """Return an SCT configured PileUpXingFolder tool"""
+ kwargs.setdefault("FirstXing", SCT_FirstXing())
+ kwargs.setdefault("LastXing", SCT_LastXing())
+ kwargs.setdefault("CacheRefreshFrequency", 1.0) # default 0 no dataproxy reset
+ kwargs.setdefault("ItemList", ["FaserSiHitCollection#SCT_Hits"] )
+ return PileUpXingFolder(name, **kwargs)
+
+
+def SCT_OutputCfg(flags):
+ """Return ComponentAccumulator with Output for SCT. Not standalone."""
+ acc = ComponentAccumulator()
+ ItemList = ["SCT_RDO_Container#*"]
+ if flags.Digitization.TruthOutput:
+ ItemList += ["InDetSimDataCollection#*"]
+ acc.merge(TruthDigitizationOutputCfg(flags))
+ acc.merge(OutputStreamCfg(flags, "RDO", ItemList))
+ return acc
+
+
+def SCT_DigitizationBasicCfg(flags, **kwargs):
+ """Return ComponentAccumulator for SCT digitization"""
+ acc = ComponentAccumulator()
+ if "PileUpTools" not in kwargs:
+ PileUpTools = acc.popToolsAndMerge(SCT_DigitizationToolCfg(flags))
+ kwargs["PileUpTools"] = PileUpTools
+ acc.merge(PileUpToolsCfg(flags, **kwargs))
+ return acc
+
+
+def SCT_OverlayDigitizationBasicCfg(flags, **kwargs):
+ """Return ComponentAccumulator with SCT Overlay digitization"""
+ acc = ComponentAccumulator()
+ if "DigitizationTool" not in kwargs:
+ tool = acc.popToolsAndMerge(SCT_OverlayDigitizationToolCfg(flags))
+ kwargs["DigitizationTool"] = tool
+ acc.addEventAlgo(SCT_Digitization(**kwargs))
+ return acc
+
+
+# with output defaults
+def SCT_DigitizationCfg(flags, **kwargs):
+ """Return ComponentAccumulator for SCT digitization and Output"""
+ acc = SCT_DigitizationBasicCfg(flags, **kwargs)
+ acc.merge(SCT_OutputCfg(flags))
+ return acc
+
+
+def SCT_OverlayDigitizationCfg(flags, **kwargs):
+ """Return ComponentAccumulator with SCT Overlay digitization and Output"""
+ acc = SCT_OverlayDigitizationBasicCfg(flags, **kwargs)
+ acc.merge(SCT_OutputCfg(flags))
+ return acc
+
+
+# additional specialisations
+def SCT_DigitizationHSCfg(flags, name="SCT_DigitizationHS", **kwargs):
+ """Return ComponentAccumulator for Hard-Scatter-only SCT digitization and Output"""
+ acc = SCT_DigitizationHSToolCfg(flags)
+ kwargs["PileUpTools"] = acc.popPrivateTools()
+ acc = SCT_DigitizationBasicCfg(flags, name=name, **kwargs)
+ acc.merge(SCT_OutputCfg(flags))
+ return acc
+
+
+def SCT_DigitizationPUCfg(flags, name="SCT_DigitizationPU", **kwargs):
+ """Return ComponentAccumulator with Pile-up-only SCT digitization and Output"""
+ acc = SCT_DigitizationPUToolCfg(flags)
+ kwargs["PileUpTools"] = acc.popPrivateTools()
+ acc = SCT_DigitizationBasicCfg(flags, name=name, **kwargs)
+ acc.merge(SCT_OutputCfg(flags))
+ return acc
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Amp.cxx b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Amp.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..22352b26451489fee9f36c56a1a5c756397942ec
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Amp.cxx
@@ -0,0 +1,143 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "SCT_Amp.h"
+
+// CLHEP
+#include "CLHEP/Units/SystemOfUnits.h"
+
+//STD includes
+#include <cmath>
+#include <fstream>
+
+//#define SCT_DIG_DEBUG
+
+// constructor
+SCT_Amp::SCT_Amp(const std::string& type, const std::string& name, const IInterface* parent)
+ : base_class(type, name, parent)
+{
+}
+
+//----------------------------------------------------------------------
+// Initialize
+//----------------------------------------------------------------------
+StatusCode SCT_Amp::initialize() {
+
+ StatusCode sc{AthAlgTool::initialize()};
+ if (sc.isFailure()) {
+ ATH_MSG_FATAL("SCT_Amp::initialize() failed");
+ return sc;
+ }
+ ATH_MSG_DEBUG("SCT_Amp::initialize()");
+
+ /** CHLEP Units */
+ m_PeakTime.setValue(m_PeakTime.value() * CLHEP::ns);
+ m_dt.setValue(m_dt.value() * CLHEP::ns);
+ m_tmin.setValue(m_tmin.value() * CLHEP::ns);
+ m_tmax.setValue(m_tmax.value() * CLHEP::ns);
+
+ m_NormConstCentral = (exp(3.0)/27.0)*(1.0-m_CrossFactor2sides)*(1.0-m_CrossFactorBack);
+ m_NormConstNeigh = exp(3.0-sqrt(3.0))/(6*(2.0*sqrt(3.0)-3.0));
+ m_NormConstNeigh *= (m_CrossFactor2sides/2.0)*(1.0-m_CrossFactorBack);
+
+#ifdef SCT_DIG_DEBUG
+ ATH_MSG_INFO("\tAmp created, PeakTime = " << m_PeakTime);
+ ATH_MSG_INFO("\tResponse will be CR-RC^3 with tp = " << m_PeakTime/3.0);
+ ATH_MSG_INFO("\tCoupling to both neighbours = " << m_CrossFactor2sides);
+ ATH_MSG_INFO("\tCoupling to backplane = " << m_CrossFactorBack);
+ ATH_MSG_INFO("\tNormalization for central " << m_NormConstCentral);
+ ATH_MSG_INFO("\tNormalization for neighbour " << m_NormConstNeigh);
+#endif
+
+ return sc;
+}
+
+//----------------------------------------------------------------------
+// Finalize
+//----------------------------------------------------------------------
+StatusCode SCT_Amp::finalize() {
+ StatusCode sc{AthAlgTool::finalize()};
+ if (sc.isFailure()) {
+ ATH_MSG_FATAL("SCT_Amp::finalize() failed");
+ return sc;
+ }
+ ATH_MSG_DEBUG("SCT_Amp::finalize()");
+ return sc;
+}
+
+//----------------------------------------------------------------------
+// Electronique response is now CR-RC^3 of the charge diode
+//----------------------------------------------------------------------
+float SCT_Amp::response(const list_t& Charges, const float timeOfThreshold) const {
+ float resp{0.0};
+ float tp{static_cast<float>(m_PeakTime/3.0)}; // for CR-RC^3
+ for (const SiCharge& charge: Charges) {
+ float ch{static_cast<float>(charge.charge())};
+ float tC{static_cast<float>(timeOfThreshold - charge.time())};
+ if (tC > 0.0) {
+ tC/=tp; //to avoid doing it four times
+ resp += ch*tC*tC*tC*exp(-tC); //faster than pow
+ }
+ }
+ return resp*m_NormConstCentral;
+}
+
+void SCT_Amp::response(const list_t& Charges, const float timeOfThreshold, std::vector<float>& response) const {
+ short bin_max{static_cast<short>(response.size())};
+ std::fill(response.begin(), response.end(), 0.0);
+ float tp{static_cast<float>(m_PeakTime/3.0)}; // for CR-RC^3
+ for (const SiCharge& charge: Charges) {
+ float ch{static_cast<float>(charge.charge())};
+ float ch_time{static_cast<float>(charge.time())};
+ short bin_end{static_cast<short>(bin_max-1)};
+ for (short bin{-1}; bin<bin_end; ++bin) {
+ float bin_timeOfThreshold{timeOfThreshold + bin*25};//25, fix me
+ float tC{bin_timeOfThreshold - ch_time};
+ if (tC > 0.0) {
+ tC/=tp; //to avoid doing it four times
+ response[bin+1] += ch*tC*tC*tC*exp(-tC); //faster than pow
+ }
+ }
+ }
+ for (short bin{0}; bin<bin_max; ++bin) response[bin] = response[bin]*m_NormConstCentral;
+ return;
+}
+
+//----------------------------------------------------------------------
+// differenciated and scaled pulse on the neighbour strip!
+//----------------------------------------------------------------------
+float SCT_Amp::crosstalk(const list_t& Charges, const float timeOfThreshold) const {
+ float resp{0};
+ float tp{static_cast<float>(m_PeakTime/3.0)}; // for CR-RC^3
+ for (const SiCharge& charge: Charges) {
+ float ch{static_cast<float>(charge.charge())};
+ float tC{static_cast<float>(timeOfThreshold - charge.time())};
+ if (tC > 0.0) {
+ tC/=tp; //to avoid doing it four times
+ resp += ch*tC*tC*exp(-tC)*(3.0-tC); //faster than pow
+ }
+ }
+ return resp*m_NormConstNeigh;
+}
+
+void SCT_Amp::crosstalk(const list_t& Charges, const float timeOfThreshold, std::vector<float>& response) const {
+ short bin_max{static_cast<short>(response.size())};
+ std::fill(response.begin(), response.end(), 0.0);
+ float tp{static_cast<float>(m_PeakTime/3.0)}; // for CR-RC^3
+ for (const SiCharge& charge: Charges) {
+ float ch{static_cast<float>(charge.charge())};
+ float ch_time{static_cast<float>(charge.time())};
+ short bin_end{static_cast<short>(bin_max-1)};
+ for (short bin{-1}; bin<bin_end; ++bin) {
+ float bin_timeOfThreshold{timeOfThreshold + bin*25}; // 25, fix me
+ float tC{bin_timeOfThreshold - ch_time};
+ if (tC > 0.0) {
+ tC/=tp; //to avoid doing it four times
+ response[bin+1] += ch*tC*tC*exp(-tC)*(3.0-tC); //faster than pow
+ }
+ }
+ }
+ for (short bin{0}; bin<bin_max; ++bin) response[bin] = response[bin]*m_NormConstNeigh;
+ return;
+}
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Amp.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Amp.h
new file mode 100644
index 0000000000000000000000000000000000000000..d066e194758804bcce453f97b8b11522831b9aa6
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Amp.h
@@ -0,0 +1,69 @@
+// -*- C++ -*-
+
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+/**
+ * Header file for SCT_Amp Class
+ * @brief A class to model an SCT amplifier and shaper. Gives a CRRC response to a
+ * list of charges with times. Also calculates average input and output for
+ * diagnostic purposes. Questions/comments to Szymon.Gadomski@cern.ch
+ * Name changed from SCTpreamp (misleading) on 09.05.01.
+ * @author Szymon.Gadomski@cern.ch, Awatif.Belymam@cern.ch, Davide.Costanzo@cern.ch,
+ * tgcornel@nikhef.nl, Grant.Gorfine@cern.ch, Paul.Bell@cern.ch, Jorgen.Dalmau@cern.ch,
+ * 23/08/2007 - Kondo.Gnanvo@cern.ch
+ * - Conversion of the SCT_Amp code AlgTool
+ */
+
+#ifndef FASERSCT_DIGITIZATION_SCTAMP_H
+#define FASERSCT_DIGITIZATION_SCTAMP_H
+
+#include "AthenaBaseComps/AthAlgTool.h"
+#include "FaserSCT_Digitization/ISCT_Amp.h"
+
+#include "TrackerSimEvent/SiCharge.h"
+
+class SCT_Amp : public extends<AthAlgTool, ISCT_Amp> {
+ public:
+
+ /** constructor */
+ SCT_Amp(const std::string& type, const std::string& name, const IInterface* parent);
+ /** Destructor */
+ virtual ~SCT_Amp() = default;
+ /** AlgTool initialize */
+ virtual StatusCode initialize() override;
+ /** AlgTool finalize */
+ virtual StatusCode finalize() override;
+
+ /** main purpose: CR-RC^3 response to a list of charges with times */
+ virtual float response(const list_t& Charges, const float timeOverThreshold) const override;
+ virtual void response(const list_t& Charges, const float timeOverThreshold, std::vector<float>& resp) const override;
+
+ /** Neighbour strip cross talk response strip to a list of charges with times */
+ virtual float crosstalk(const list_t& Charges, const float timeOverThreshold) const override;
+ virtual void crosstalk(const list_t& Charges, const float timeOverThreshold, std::vector<float>& resp) const override;
+
+private:
+
+ /** signal peak time */
+ FloatProperty m_PeakTime{this, "PeakTime", 21., "Front End Electronics peaking time"};
+
+ /** Cross factor 2 side */
+ FloatProperty m_CrossFactor2sides{this, "CrossFactor2sides", 0.1, "Loss of charge to neighbour strip constant"};
+
+ /** cross factor */
+ FloatProperty m_CrossFactorBack{this, "CrossFactorBack", 0.07, "Loss of charge to back plane constant"};
+
+ FloatProperty m_tmin{this, "Tmin", -25.0};
+ FloatProperty m_tmax{this, "Tmax", 150.0};
+ FloatProperty m_dt{this, "deltaT", 1.0};
+
+ /** Normalisation factor for the signal response */
+ float m_NormConstCentral{0.};
+
+ /** Normalisation factor for the neighbour strip signal response */
+ float m_NormConstNeigh{0.};
+};
+
+#endif // FASERSCT_DIGITIZATION_SCTAMP_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Digitization.cxx b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Digitization.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..f30bbd62b14fd201ca98a739de49894bafe9c6bc
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Digitization.cxx
@@ -0,0 +1,31 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "SCT_Digitization.h"
+#include "PileUpTools/IPileUpTool.h"
+
+//----------------------------------------------------------------------
+// Constructor with parameters:
+//----------------------------------------------------------------------
+SCT_Digitization::SCT_Digitization(const std::string& name, ISvcLocator* pSvcLocator) :
+ AthAlgorithm(name, pSvcLocator)
+{
+}
+
+//----------------------------------------------------------------------
+// Initialize method:
+//----------------------------------------------------------------------
+StatusCode SCT_Digitization::initialize() {
+ ATH_MSG_DEBUG("SCT_Digitization::initialize()");
+ return StatusCode::SUCCESS ;
+}
+
+//----------------------------------------------------------------------
+// Execute method:
+//----------------------------------------------------------------------
+
+StatusCode SCT_Digitization::execute() {
+ ATH_MSG_DEBUG("execute()");
+ return m_sctDigitizationTool->processAllSubEvents();
+}
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Digitization.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Digitization.h
new file mode 100644
index 0000000000000000000000000000000000000000..3f32cdf0734f555ef6e64b1b476c07218201a534
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_Digitization.h
@@ -0,0 +1,48 @@
+// -*- C++ -*-
+
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+/** @file SCT_Digitization.h Header file for SCT_Digitization class.
+ * @author Szymon.Gadomski@cern.ch, Awatif.Belymam@cern.ch, Davide.Costanzo@cern.ch,
+ * tgcornel@nikhef.nl, Grant.Gorfine@cern.ch, Paul.Bell@cern.ch,
+ * Jorgen.Dalmau@cern.ch, Kondo.Gnanvo@cern.ch, and others
+ * Version 23/08/2007 Kondo.Gnanvo@cern.ch
+ * Conversion of the processors into AlgTool
+ */
+
+// Multiple inclusion protection
+#ifndef FASERSCT_DIGITIZATION_SCT_DIGITIZATION_H
+#define FASERSCT_DIGITIZATION_SCT_DIGITIZATION_H
+
+// Base class
+#include "AthenaBaseComps/AthAlgorithm.h"
+// Gaudi
+#include "GaudiKernel/ToolHandle.h"
+
+class IPileUpTool;
+
+/** Top algorithm class for SCT digitization */
+class SCT_Digitization : public AthAlgorithm {
+
+ public:
+
+ /** Constructor with parameters */
+ SCT_Digitization(const std::string& name, ISvcLocator* pSvcLocator);
+
+ /** Destructor */
+ virtual ~SCT_Digitization() = default;
+
+ /** Basic algorithm methods */
+ virtual StatusCode initialize() override final;
+ virtual StatusCode execute() override final;
+ virtual bool isClonable() const override final { return true; }
+
+ private:
+
+ ToolHandle<IPileUpTool> m_sctDigitizationTool{this, "DigitizationTool", "SCT_DigitizationTool", "SCT_DigitizationTool name"};
+
+};
+
+#endif // FASERSCT_DIGITIZATION_SCT_DIGITIZATION_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_DigitizationTool.cxx b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_DigitizationTool.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..981e3e97d98de8da7c12d9c5f15c91406255d8ee
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_DigitizationTool.cxx
@@ -0,0 +1,795 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "SCT_DigitizationTool.h"
+
+// Mother Package includes
+#include "FaserSiDigitization/SiHelper.h"
+#include "FaserSiDigitization/SiChargedDiodeCollection.h"
+
+// EDM includes
+#include "InDetRawData/SCT1_RawData.h"
+#include "InDetRawData/SCT3_RawData.h"
+
+// Hit class includes
+#include "TrackerSimEvent/FaserSiHit.h"
+#include "Identifier/Identifier.h"
+
+// Det Descr includes
+#include "TrackerReadoutGeometry/SiDetectorElement.h"
+#include "TrackerReadoutGeometry/SCT_ModuleSideDesign.h"
+
+// Data Handle
+#include "StoreGate/ReadCondHandle.h"
+#include "StoreGate/ReadHandle.h"
+
+// Random Number Generation
+#include "AthenaKernel/RNGWrapper.h"
+#include "CLHEP/Random/RandomEngine.h"
+
+// C++ Standard Library
+#include <memory>
+#include <sstream>
+
+// Barcodes at the HepMC level are int
+
+using TrackerDD::SiCellId;
+
+SCT_DigitizationTool::SCT_DigitizationTool(const std::string& type,
+ const std::string& name,
+ const IInterface* parent) :
+ base_class(type, name, parent) {
+ m_WriteSCT1_RawData.declareUpdateHandler(&SCT_DigitizationTool::SetupRdoOutputType, this);
+}
+
+SCT_DigitizationTool::~SCT_DigitizationTool() {
+ delete m_thpcsi;
+ for (FaserSiHitCollection* hit: m_hitCollPtrs) {
+ hit->Clear();
+ delete hit;
+ }
+ m_hitCollPtrs.clear();
+}
+
+// ----------------------------------------------------------------------
+// Initialize method:
+// ----------------------------------------------------------------------
+StatusCode SCT_DigitizationTool::initialize() {
+ ATH_MSG_DEBUG("SCT_DigitizationTool::initialize()");
+
+ // +++ Init the services
+ ATH_CHECK(initServices());
+
+ // +++ Get the Surface Charges Generator tool
+ ATH_CHECK(initSurfaceChargesGeneratorTool());
+
+ // +++ Get the Front End tool
+ ATH_CHECK(initFrontEndTool());
+
+ // +++ Initialise for disabled cells from the random disabled cells tool
+ // +++ Default off, since disabled cells taken form configuration in
+ // reconstruction stage
+ if (m_randomDisabledCells) {
+ ATH_CHECK(initDisabledCells());
+ ATH_MSG_INFO("Use of Random disabled cells");
+ } else {
+ m_sct_RandomDisabledCellGenerator.disable();
+ }
+
+ // check the input object name
+ if (m_hitsContainerKey.key().empty()) {
+ ATH_MSG_FATAL("Property InputObjectName not set !");
+ return StatusCode::FAILURE;
+ }
+ if(m_onlyUseContainerName) m_inputObjectName = m_hitsContainerKey.key();
+ ATH_MSG_DEBUG("Input objects in container : '" << m_inputObjectName << "'");
+
+ // Initialize ReadHandleKey
+ ATH_CHECK(m_hitsContainerKey.initialize(!m_onlyUseContainerName));
+
+ // +++ Initialize WriteHandleKey
+ ATH_CHECK(m_rdoContainerKey.initialize());
+ ATH_CHECK(m_simDataCollMapKey.initialize());
+
+ // Initialize ReadCondHandleKey
+ ATH_CHECK(m_SCTDetEleCollKey.initialize());
+
+ ATH_MSG_DEBUG("SiDigitizationTool::initialize() complete");
+
+ return StatusCode::SUCCESS;
+}
+
+namespace {
+ class SiDigitizationSurfaceChargeInserter
+ : public ISiSurfaceChargesInserter
+ {
+ public:
+ SiDigitizationSurfaceChargeInserter(const TrackerDD::SiDetectorElement* sielement,
+ SiChargedDiodeCollection* chargedDiodes)
+ : m_sielement(sielement),
+ m_chargedDiodes(chargedDiodes) {
+ }
+
+ void operator () (const SiSurfaceCharge& scharge) const;
+ private:
+ const TrackerDD::SiDetectorElement* m_sielement;
+ SiChargedDiodeCollection* m_chargedDiodes;
+ };
+
+
+ void SiDigitizationSurfaceChargeInserter::operator ()
+ (const SiSurfaceCharge& scharge) const {
+ // get the diode in which this charge is
+ SiCellId diode{m_sielement->cellIdOfPosition(scharge.position())};
+
+ if (diode.isValid()) {
+ // add this charge to the collection (or merge in existing charged diode)
+ m_chargedDiodes->add(diode, scharge.charge());
+ }
+ }
+} // anonymous namespace
+
+// ----------------------------------------------------------------------
+// Initialise the surface charge generator Tool
+// ----------------------------------------------------------------------
+StatusCode SCT_DigitizationTool::initSurfaceChargesGeneratorTool() {
+ ATH_CHECK(m_sct_SurfaceChargesGenerator.retrieve());
+
+ if (m_cosmicsRun and m_tfix > -998) {
+ m_sct_SurfaceChargesGenerator->setFixedTime(m_tfix);
+ ATH_MSG_INFO("Use of FixedTime = " << m_tfix << " in cosmics");
+ }
+
+ ATH_MSG_DEBUG("Retrieved and initialised tool " << m_sct_SurfaceChargesGenerator);
+
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+// Initialise the Front End electronics Tool
+// ----------------------------------------------------------------------
+StatusCode SCT_DigitizationTool::initFrontEndTool() {
+ ATH_CHECK(m_sct_FrontEnd.retrieve());
+
+ storeTool(&(*m_sct_FrontEnd));
+
+ ATH_MSG_DEBUG("Retrieved and initialised tool " << m_sct_FrontEnd);
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+// Initialize the different services
+// ----------------------------------------------------------------------
+StatusCode SCT_DigitizationTool::initServices() {
+ // Get SCT ID helper for hash function and Store them using methods from the
+ // SiDigitization.
+ ATH_CHECK(detStore()->retrieve(m_detID, "FaserSCT_ID"));
+
+ ATH_CHECK(m_mergeSvc.retrieve());
+ ATH_CHECK(m_rndmSvc.retrieve());
+
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+// Initialize the disabled cells for cosmics or CTB cases
+// ----------------------------------------------------------------------
+StatusCode SCT_DigitizationTool::initDisabledCells() {
+ // +++ Retrieve the SCT_RandomDisabledCellGenerator
+ ATH_CHECK(m_sct_RandomDisabledCellGenerator.retrieve());
+
+ storeTool(&(*m_sct_RandomDisabledCellGenerator));
+
+ ATH_MSG_INFO("Retrieved the SCT_RandomDisabledCellGenerator tool:" << m_sct_RandomDisabledCellGenerator);
+ return StatusCode::SUCCESS;
+}
+
+StatusCode SCT_DigitizationTool::processAllSubEvents() {
+ if (prepareEvent(0).isFailure()) {
+ return StatusCode::FAILURE;
+ }
+ // Set the RNG to use for this event.
+ ATHRNG::RNGWrapper* rngWrapper = m_rndmSvc->getEngine(this);
+ rngWrapper->setSeed( name(), Gaudi::Hive::currentContext() );
+ CLHEP::HepRandomEngine *rndmEngine = *rngWrapper;
+
+ ATH_MSG_VERBOSE("Begin digitizeAllHits");
+ if (m_enableHits and (not getNextEvent().isFailure())) {
+ digitizeAllHits(&m_rdoContainer, &m_simDataCollMap, &m_processedElements, m_thpcsi, rndmEngine);
+ } else {
+ ATH_MSG_DEBUG("no hits found in event!");
+ }
+ ATH_MSG_DEBUG("Digitized Elements with Hits");
+
+ // loop over elements without hits
+ if (not m_onlyHitElements) {
+ digitizeNonHits(&m_rdoContainer, &m_simDataCollMap, &m_processedElements, rndmEngine);
+ ATH_MSG_DEBUG("Digitized Elements without Hits");
+ }
+
+ delete m_thpcsi;
+ m_thpcsi = nullptr;
+
+ ATH_MSG_VERBOSE("Digitize success!");
+ return StatusCode::SUCCESS;
+}
+
+// ======================================================================
+// prepareEvent
+// ======================================================================
+StatusCode SCT_DigitizationTool::prepareEvent(unsigned int /*index*/) {
+ ATH_MSG_VERBOSE("SCT_DigitizationTool::prepareEvent()");
+ // Create the IdentifiableContainer to contain the digit collections Create
+ // a new RDO container
+ m_rdoContainer = SG::makeHandle(m_rdoContainerKey);
+ ATH_CHECK(m_rdoContainer.record(std::make_unique<SCT_RDO_Container>(m_detID->wafer_hash_max())));
+
+ // Create a map for the SDO and register it into StoreGate
+ m_simDataCollMap = SG::makeHandle(m_simDataCollMapKey);
+ ATH_CHECK(m_simDataCollMap.record(std::make_unique<InDetSimDataCollection>()));
+
+ m_processedElements.clear();
+ m_processedElements.resize(m_detID->wafer_hash_max(), false);
+
+ m_thpcsi = new TimedHitCollection<FaserSiHit>();
+// m_HardScatterSplittingSkipper = false;
+ return StatusCode::SUCCESS;
+}
+
+// =========================================================================
+// mergeEvent
+// =========================================================================
+StatusCode SCT_DigitizationTool::mergeEvent() {
+ ATH_MSG_VERBOSE("SCT_DigitizationTool::mergeEvent()");
+
+ // Set the RNG to use for this event.
+ ATHRNG::RNGWrapper* rngWrapper = m_rndmSvc->getEngine(this);
+ rngWrapper->setSeed( name(), Gaudi::Hive::currentContext() );
+ CLHEP::HepRandomEngine *rndmEngine = *rngWrapper;
+
+ if (m_enableHits) {
+ digitizeAllHits(&m_rdoContainer, &m_simDataCollMap, &m_processedElements, m_thpcsi, rndmEngine);
+ }
+
+ if (not m_onlyHitElements) {
+ digitizeNonHits(&m_rdoContainer, &m_simDataCollMap, &m_processedElements, rndmEngine);
+ }
+
+ for (FaserSiHitCollection* hit: m_hitCollPtrs) {
+ hit->Clear();
+ delete hit;
+ }
+ m_hitCollPtrs.clear();
+
+ delete m_thpcsi;
+ m_thpcsi = nullptr;
+
+ ATH_MSG_DEBUG("Digitize success!");
+ return StatusCode::SUCCESS;
+}
+
+void SCT_DigitizationTool::digitizeAllHits(SG::WriteHandle<SCT_RDO_Container>* rdoContainer, SG::WriteHandle<InDetSimDataCollection>* simDataCollMap, std::vector<bool>* processedElements, TimedHitCollection<FaserSiHit>* thpcsi, CLHEP::HepRandomEngine * rndmEngine) const {
+ /////////////////////////////////////////////////
+ //
+ // In order to process all element rather than just those with hits we
+ // create a vector to keep track of which elements have been processed.
+ // NB. an element is an sct module
+ //
+ /////////////////////////////////////////////////
+ ATH_MSG_DEBUG("Digitizing hits");
+ int hitcount{0}; // First, elements with hits.
+
+ SiChargedDiodeCollection chargedDiodes;
+
+ while (digitizeElement(&chargedDiodes, thpcsi, rndmEngine)) {
+ ATH_MSG_DEBUG("Hit collection ID=" << m_detID->show_to_string(chargedDiodes.identify()));
+
+ hitcount++; // Hitcount will be a number in the hit collection minus
+ // number of hits in missing mods
+
+ ATH_MSG_DEBUG("in digitize elements with hits: station - layer - eta - phi "
+ << m_detID->station(chargedDiodes.identify()) << " - "
+ << m_detID->layer(chargedDiodes.identify()) << " - "
+ << m_detID->eta_module(chargedDiodes.identify()) << " - "
+ << m_detID->phi_module(chargedDiodes.identify()) << " - "
+ << " processing hit number " << hitcount);
+
+ // Have a flag to check if the module is present or not
+ // Generally assume it is:
+
+ IdentifierHash idHash{chargedDiodes.identifyHash()};
+
+ assert(idHash < processedElements->size());
+ (*processedElements)[idHash] = true;
+
+ // create and store RDO and SDO
+
+ if (not chargedDiodes.empty()) {
+ StatusCode sc{createAndStoreRDO(&chargedDiodes, rdoContainer)};
+ if (sc.isSuccess()) { // error msg is given inside
+ // createAndStoreRDO()
+ addSDO(&chargedDiodes, simDataCollMap);
+ }
+ }
+
+ chargedDiodes.clear();
+ }
+ ATH_MSG_DEBUG("hits processed");
+ return;
+}
+
+// digitize elements without hits
+void SCT_DigitizationTool::digitizeNonHits(SG::WriteHandle<SCT_RDO_Container>* rdoContainer, SG::WriteHandle<InDetSimDataCollection>* simDataCollMap, const std::vector<bool>* processedElements, CLHEP::HepRandomEngine * rndmEngine) const {
+ // Get SCT_DetectorElementCollection
+ SG::ReadCondHandle<TrackerDD::SiDetectorElementCollection> sctDetEle(m_SCTDetEleCollKey);
+ const TrackerDD::SiDetectorElementCollection* elements{sctDetEle.retrieve()};
+ if (elements==nullptr) {
+ ATH_MSG_FATAL(m_SCTDetEleCollKey.fullKey() << " could not be retrieved");
+ return;
+ }
+
+ ATH_MSG_DEBUG("processing elements without hits");
+ SiChargedDiodeCollection chargedDiodes;
+
+ for (unsigned int i{0}; i < processedElements->size(); i++) {
+ if (not (*processedElements)[i]) {
+ IdentifierHash idHash{i};
+ if (not idHash.is_valid()) {
+ ATH_MSG_ERROR("SCT Detector element id hash is invalid = " << i);
+ }
+
+ const TrackerDD::SiDetectorElement* element{elements->getDetectorElement(idHash)};
+ if (element) {
+ ATH_MSG_DEBUG("In digitize of untouched elements: layer - phi - eta "
+ << m_detID->layer(element->identify()) << " - "
+ << m_detID->phi_module(element->identify()) << " - "
+ << m_detID->eta_module(element->identify()) << " - "
+ << "size: " << processedElements->size());
+
+ chargedDiodes.setDetectorElement(element);
+ ATH_MSG_DEBUG("calling applyProcessorTools() for NON hits");
+ applyProcessorTools(&chargedDiodes, rndmEngine);
+
+ // Create and store RDO and SDO
+ // Don't create empty ones.
+ if (not chargedDiodes.empty()) {
+ StatusCode sc{createAndStoreRDO(&chargedDiodes, rdoContainer)};
+ if (sc.isSuccess()) {// error msg is given inside
+ // createAndStoreRDO()
+ addSDO(&chargedDiodes, simDataCollMap);
+ }
+ }
+
+ chargedDiodes.clear();
+ }
+ }
+ }
+
+ return;
+}
+
+bool SCT_DigitizationTool::digitizeElement(SiChargedDiodeCollection* chargedDiodes, TimedHitCollection<FaserSiHit>*& thpcsi, CLHEP::HepRandomEngine * rndmEngine) const {
+ if (nullptr == thpcsi) {
+ ATH_MSG_ERROR("thpcsi should not be nullptr!");
+
+ return false;
+ }
+
+ // get the iterator pairs for this DetEl
+
+ TimedHitCollection<FaserSiHit>::const_iterator i, e;
+ if (thpcsi->nextDetectorElement(i, e) == false) { // no more hits
+ return false;
+ }
+
+ // create the identifier for the collection:
+ ATH_MSG_DEBUG("create ID for the hit collection");
+ const TimedHitPtr<FaserSiHit>& firstHit{*i};
+ int barrel{firstHit->getStation()};
+ Identifier id{m_detID->wafer_id(barrel,
+ firstHit->getPlane(),
+ firstHit->getRow(),
+ firstHit->getModule(),
+ firstHit->getSensor())};
+ IdentifierHash waferHash{m_detID->wafer_hash(id)};
+
+ // Get SCT_DetectorElementCollection
+ SG::ReadCondHandle<TrackerDD::SiDetectorElementCollection> sctDetEle(m_SCTDetEleCollKey);
+ const TrackerDD::SiDetectorElementCollection* elements(sctDetEle.retrieve());
+ if (elements==nullptr) {
+ ATH_MSG_FATAL(m_SCTDetEleCollKey.fullKey() << " could not be retrieved");
+ return false;
+ }
+
+ // get the det element from the manager
+ const TrackerDD::SiDetectorElement* sielement{elements->getDetectorElement(waferHash)};
+
+ if (sielement == nullptr) {
+ ATH_MSG_DEBUG("Station=" << barrel << " layer=" << firstHit->getPlane() << " Row=" << firstHit->getRow() << " Module=" << firstHit->getModule() << " Side=" << firstHit->getSensor());
+ ATH_MSG_ERROR("detector manager could not find element with id = " << id);
+ return false;
+ }
+ // create the charged diodes collection
+ chargedDiodes->setDetectorElement(sielement);
+
+ // Loop over the hits and created charged diodes:
+ while (i != e) {
+ TimedHitPtr<FaserSiHit> phit{*i++};
+
+ // skip hits which are more than 10us away
+ if (fabs(phit->meanTime()) < 10000. * CLHEP::ns) {
+ ATH_MSG_DEBUG("HASH = " << m_detID->wafer_hash(m_detID->wafer_id(phit->getStation(),
+ phit->getPlane(),
+ phit->getRow(),
+ phit->getModule(),
+ phit->getSensor())));
+ ATH_MSG_DEBUG("calling process() for all methods");
+ m_sct_SurfaceChargesGenerator->process(sielement, phit, SiDigitizationSurfaceChargeInserter(sielement, chargedDiodes), rndmEngine);
+ ATH_MSG_DEBUG("charges filled!");
+ }
+ }
+ applyProcessorTools(chargedDiodes, rndmEngine); // !< Use of the new AlgTool surface
+ // charges generator class
+ return true;
+}
+
+// -----------------------------------------------------------------------------
+// Applies processors to the current detector element for the current element:
+// -----------------------------------------------------------------------------
+void SCT_DigitizationTool::applyProcessorTools(SiChargedDiodeCollection* chargedDiodes, CLHEP::HepRandomEngine * rndmEngine) const {
+ ATH_MSG_DEBUG("applyProcessorTools()");
+ int processorNumber{0};
+
+ for (ISiChargedDiodesProcessorTool* proc: m_diodeCollectionTools) {
+ proc->process(*chargedDiodes, rndmEngine);
+
+ processorNumber++;
+ ATH_MSG_DEBUG("Applied processor # " << processorNumber);
+ }
+}
+
+StatusCode SCT_DigitizationTool::processBunchXing(int bunchXing,
+ SubEventIterator bSubEvents,
+ SubEventIterator eSubEvents) {
+ ATH_MSG_VERBOSE("SCT_DigitizationTool::processBunchXing() " << bunchXing);
+
+ typedef PileUpMergeSvc::TimedList<FaserSiHitCollection>::type TimedHitCollList;
+ TimedHitCollList hitCollList;
+
+ if ((not (m_mergeSvc->retrieveSubSetEvtData(m_inputObjectName, hitCollList, bunchXing,
+ bSubEvents, eSubEvents).isSuccess())) and
+ hitCollList.size() == 0) {
+ ATH_MSG_ERROR("Could not fill TimedHitCollList");
+ return StatusCode::FAILURE;
+ } else {
+ ATH_MSG_VERBOSE(hitCollList.size() << " SiHitCollections with key " <<
+ m_inputObjectName << " found");
+ }
+
+ TimedHitCollList::iterator endColl{hitCollList.end()};
+ for (TimedHitCollList::iterator iColl{hitCollList.begin()}; iColl != endColl; iColl++) {
+ FaserSiHitCollection *hitCollPtr{new FaserSiHitCollection(*iColl->second)};
+ PileUpTimeEventIndex timeIndex{iColl->first};
+ ATH_MSG_DEBUG("SiHitCollection found with " << hitCollPtr->size() <<
+ " hits");
+ ATH_MSG_VERBOSE("time index info. time: " << timeIndex.time()
+ << " index: " << timeIndex.index()
+ << " type: " << timeIndex.type());
+ m_thpcsi->insert(timeIndex, hitCollPtr);
+ m_hitCollPtrs.push_back(hitCollPtr);
+ }
+
+ return StatusCode::SUCCESS;
+
+}
+
+// =========================================================================
+// property handlers
+// =========================================================================
+void SCT_DigitizationTool::SetupRdoOutputType(Property &) {
+}
+
+// Does nothing, but required by Gaudi
+
+// ----------------------------------------------------------------------
+// Digitisation of non hit elements
+// ----------------------------------------------------------------------
+
+class DigitizeNonHitElementsDebugPrinter
+{
+public:
+ DigitizeNonHitElementsDebugPrinter(const FaserSCT_ID* detID) :
+ m_detID{detID}, m_msgNo{-1} {
+ }
+
+ std::string msg(const TrackerDD::SiDetectorElement* element) {
+ std::ostringstream ost;
+
+ ost << "Digitized unprocessed elements: layer - phi - eta - side "
+ << m_detID->layer(element->identify()) << " - "
+ << m_detID->phi_module(element->identify()) << " - "
+ << m_detID->eta_module(element->identify()) << " - "
+ << m_detID->side(element->identify()) << " - "
+ << " unprocessed hit number: " << ++m_msgNo << '\n';
+
+ return ost.str();
+ }
+
+private:
+ const FaserSCT_ID* m_detID;
+ int m_msgNo;
+};
+
+// ----------------------------------------------------------------------//
+// createAndStoreRDO //
+// ----------------------------------------------------------------------//
+StatusCode SCT_DigitizationTool::createAndStoreRDO(SiChargedDiodeCollection* chDiodeCollection, SG::WriteHandle<SCT_RDO_Container>* rdoContainer) const {
+
+ // Create the RDO collection
+ SCT_RDO_Collection* RDOColl{createRDO(chDiodeCollection)};
+
+ // Add it to storegate
+ if ((*rdoContainer)->addCollection(RDOColl, RDOColl->identifyHash()).isFailure()) {
+ ATH_MSG_FATAL("SCT RDO collection could not be added to container!");
+ delete RDOColl;
+ RDOColl = nullptr;
+ return StatusCode::FAILURE;
+ }
+ return StatusCode::SUCCESS;
+} // SCT_Digitization::createAndStoreRDO()
+
+// ----------------------------------------------------------------------
+// createRDO
+// ----------------------------------------------------------------------
+SCT_RDO_Collection* SCT_DigitizationTool::createRDO(SiChargedDiodeCollection* collection) const {
+
+ // create a new SCT RDO collection
+ SCT_RDO_Collection* p_rdocoll{nullptr};
+
+ // need the DE identifier
+ const Identifier id_de{collection->identify()};
+ IdentifierHash idHash_de{collection->identifyHash()};
+ try {
+ p_rdocoll = new SCT_RDO_Collection(idHash_de);
+ } catch (const std::bad_alloc&) {
+ ATH_MSG_FATAL("Could not create a new SCT_RDORawDataCollection !");
+ }
+ p_rdocoll->setIdentifier(id_de);
+
+ SiChargedDiodeIterator i_chargedDiode{collection->begin()};
+ SiChargedDiodeIterator i_chargedDiode_end{collection->end()};
+ // Choice of producing SCT1_RawData or SCT3_RawData
+ if (m_WriteSCT1_RawData.value()) {
+ for (; i_chargedDiode != i_chargedDiode_end; ++i_chargedDiode) {
+ unsigned int flagmask{static_cast<unsigned int>((*i_chargedDiode).second.flag() & 0xFE)};
+
+ if (!flagmask) { // now check it wasn't masked:
+ // create new SCT RDO, using method 1 for mask:
+ // GroupSize=1: need readout id, make use of
+ // SiTrackerDetDescr
+ TrackerDD::SiReadoutCellId roCell{(*i_chargedDiode).second.getReadoutCell()};
+ int strip{roCell.strip()};
+ if (strip > 0xffff) { // In upgrade layouts strip can be bigger
+ // than 4000
+ ATH_MSG_FATAL("Strip number too big for SCT1 raw data format.");
+ }
+ const Identifier id_readout{m_detID->strip_id(collection->identify(), strip)};
+
+ // build word, masks taken from SiTrackerEvent/SCTRawData.cxx
+ const unsigned int strip_rdo{static_cast<unsigned int>((strip & 0xFFFF) << 16)};
+
+ // user can define what GroupSize is, here 1: TC. Incorrect,
+ // GroupSize >= 1
+ int size{SiHelper::GetStripNum((*i_chargedDiode).second)};
+ unsigned int size_rdo{static_cast<unsigned int>(size & 0xFFFF)};
+
+ // TC. Need to check if there are disabled strips in the cluster
+ int cluscounter{0};
+ if (size > 1) {
+ SiChargedDiodeIterator it2{i_chargedDiode};
+ ++it2;
+ for (; it2 != i_chargedDiode_end; ++it2) {
+ ++cluscounter;
+ if (cluscounter >= size) {
+ break;
+ }
+ if (it2->second.flag() & 0xDE) {
+ int tmp{cluscounter};
+ while ((it2 != i_chargedDiode_end) and (cluscounter < size - 1) and (it2->second.flag() & 0xDE)) {
+ it2++;
+ cluscounter++;
+ }
+ if ((it2 != collection->end()) and !(it2->second.flag() & 0xDE)) {
+ SiHelper::ClusterUsed(it2->second, false);
+ SiHelper::SetStripNum(it2->second, size - cluscounter, &msg());
+ }
+ // groupSize=tmp;
+ size_rdo = tmp & 0xFFFF;
+ break;
+ }
+ }
+ }
+ unsigned int SCT_Word{strip_rdo | size_rdo};
+ SCT1_RawData* p_rdo{new SCT1_RawData(id_readout, SCT_Word)};
+ if (p_rdo) {
+ p_rdocoll->push_back(p_rdo);
+ }
+ }
+ }
+ } else {
+ // Under the current scheme time bin and ERRORS are hard-coded to
+ // default values.
+ int ERRORS{0};
+ static std::vector<int> dummyvector;
+ for (; i_chargedDiode != i_chargedDiode_end; ++i_chargedDiode) {
+ unsigned int flagmask{static_cast<unsigned int>((*i_chargedDiode).second.flag() & 0xFE)};
+
+ if (!flagmask) { // Check it wasn't masked
+ int tbin{SiHelper::GetTimeBin((*i_chargedDiode).second)};
+ // create new SCT RDO
+ TrackerDD::SiReadoutCellId roCell{(*i_chargedDiode).second.getReadoutCell()};
+ int strip{roCell.strip()};
+ Identifier id_readout;
+ id_readout = m_detID->strip_id(collection->identify(), strip);
+
+ // build word (compatible with
+ // SCT_RawDataByteStreamCnv/src/SCT_RodDecoder.cxx)
+ int size{SiHelper::GetStripNum((*i_chargedDiode).second)};
+ int groupSize{size};
+
+ // TC. Need to check if there are disabled strips in the cluster
+ int cluscounter{0};
+ if (size > 1) {
+ SiChargedDiode* diode{i_chargedDiode->second.nextInCluster()};
+ while (diode) {//check if there is a further strip in the cluster
+ ++cluscounter;
+ if (cluscounter >= size) {
+ ATH_MSG_WARNING("Cluster size reached while neighbouring strips still defined.");
+ break;
+ }
+ if (diode->flag() & 0xDE) {//see if it is disabled/below threshold/disconnected/etc (0xDE corresponds to BT_SET | DISABLED_SET | BADTOT_SET | DISCONNECTED_SET | MASKOFF_SET)
+ int tmp{cluscounter};
+ while ((cluscounter < size - 1) and (diode->flag() & 0xDE)) { //check its not the end and still disabled
+ diode = diode->nextInCluster();
+ cluscounter++;
+ }
+ if (diode and !(diode->flag() & 0xDE)) {
+ SiHelper::ClusterUsed(*diode, false);
+ SiHelper::SetStripNum(*diode, size - cluscounter, &msg());
+ }
+ groupSize = tmp;
+ break;
+ }
+ diode = diode->nextInCluster();
+ }
+ }
+
+ int stripIn11bits{strip & 0x7ff};
+ if (stripIn11bits != strip) {
+ ATH_MSG_DEBUG("Strip number " << strip << " doesn't fit into 11 bits - will be truncated");
+ }
+
+ unsigned int SCT_Word{static_cast<unsigned int>(groupSize | (stripIn11bits << 11) | (tbin << 22) | (ERRORS << 25))};
+ SCT3_RawData *p_rdo{new SCT3_RawData(id_readout, SCT_Word, &dummyvector)};
+ if (p_rdo) {
+ p_rdocoll->push_back(p_rdo);
+ }
+ }
+ }
+ }
+ return p_rdocoll;
+} // SCT_Digitization::createRDO()
+
+// ------------------------------------------------------------
+// Get next event and extract collection of hit collections:
+// ------------------------------------------------------------
+StatusCode SCT_DigitizationTool::getNextEvent() {
+ ATH_MSG_DEBUG("SCT_DigitizationTool::getNextEvent()");
+ // get the container(s)
+ typedef PileUpMergeSvc::TimedList<FaserSiHitCollection>::type TimedHitCollList;
+ // this is a list<pair<time_t, DataLink<SiHitCollection> >
+
+ // In case of single hits container just load the collection using read handles
+ if (!m_onlyUseContainerName) {
+ SG::ReadHandle<FaserSiHitCollection> hitCollection(m_hitsContainerKey);
+ if (!hitCollection.isValid()) {
+ ATH_MSG_ERROR("Could not get SCT SiHitCollection container " << hitCollection.name() << " from store " << hitCollection.store());
+ return StatusCode::FAILURE;
+ }
+
+ // create a new hits collection
+ m_thpcsi = new TimedHitCollection<FaserSiHit>{1};
+ m_thpcsi->insert(0, hitCollection.cptr());
+ ATH_MSG_DEBUG("SiHitCollection found with " << hitCollection->size() << " hits");
+
+ return StatusCode::SUCCESS;
+ }
+
+ TimedHitCollList hitCollList;
+ unsigned int numberOfSiHits{0};
+ if (not (m_mergeSvc->retrieveSubEvtsData(m_inputObjectName, hitCollList, numberOfSiHits).isSuccess()) and hitCollList.size() == 0) {
+ ATH_MSG_ERROR("Could not fill TimedHitCollList");
+ return StatusCode::FAILURE;
+ } else {
+ ATH_MSG_DEBUG(hitCollList.size() << " SiHitCollections with key " << m_inputObjectName << " found");
+ }
+ // create a new hits collection
+ m_thpcsi = new TimedHitCollection<FaserSiHit>{numberOfSiHits};
+ // now merge all collections into one
+ TimedHitCollList::iterator endColl{hitCollList.end()};
+ for (TimedHitCollList::iterator iColl{hitCollList.begin()}; iColl != endColl; ++iColl) {
+ const FaserSiHitCollection* p_collection{iColl->second};
+ m_thpcsi->insert(iColl->first, p_collection);
+ ATH_MSG_DEBUG("SiTrackerHitCollection found with " << p_collection->size() << " hits"); // loop on the hit collections
+ }
+ return StatusCode::SUCCESS;
+}
+
+// -----------------------------------------------------------------------------------------------
+// Convert a SiTotalCharge to a InDetSimData, and store it.
+// -----------------------------------------------------------------------------------------------
+void SCT_DigitizationTool::addSDO(SiChargedDiodeCollection* collection, SG::WriteHandle<InDetSimDataCollection>* simDataCollMap) const {
+ typedef SiTotalCharge::list_t list_t;
+ std::vector<InDetSimData::Deposit> deposits;
+ deposits.reserve(5); // no idea what a reasonable number for this would be
+ // with pileup
+ // loop over the charged diodes
+ SiChargedDiodeIterator EndOfDiodeCollection{collection->end()};
+ for (SiChargedDiodeIterator i_chargedDiode{collection->begin()}; i_chargedDiode != EndOfDiodeCollection; ++i_chargedDiode) {
+ deposits.clear();
+ const list_t& charges{(*i_chargedDiode).second.totalCharge().chargeComposition()};
+
+ bool real_particle_hit{false};
+ // loop over the list
+ list_t::const_iterator EndOfChargeList{charges.end()};
+ for (list_t::const_iterator i_ListOfCharges{charges.begin()}; i_ListOfCharges != EndOfChargeList; ++i_ListOfCharges) {
+ const HepMcParticleLink& trkLink{i_ListOfCharges->particleLink()};
+ const int barcode{trkLink.barcode()};
+ if ((barcode == 0) or (barcode == m_vetoThisBarcode)) {
+ continue;
+ }
+ if (not real_particle_hit) {
+ // Types of SiCharges expected from SCT
+ // Noise: barcode==0 and
+ // processType()==SiCharge::noise
+ // Delta Rays: barcode==0 and
+ // processType()==SiCharge::track
+ // Pile Up Tracks With No Truth: barcode!=0 and
+ // processType()==SiCharge::cut_track
+ // Tracks With Truth: barcode!=0 and
+ // processType()==SiCharge::track
+ if (barcode != 0 and i_ListOfCharges->processType() == SiCharge::track) {
+ real_particle_hit = true;
+ }
+ }
+ // check if this track number has been already used.
+ std::vector<InDetSimData::Deposit>::reverse_iterator theDeposit{deposits.rend()}; // dummy value
+ std::vector<InDetSimData::Deposit>::reverse_iterator depositsR_end{deposits.rend()};
+ std::vector<InDetSimData::Deposit>::reverse_iterator i_Deposit{deposits.rbegin()};
+ for (; i_Deposit != depositsR_end; ++i_Deposit) {
+ if ((*i_Deposit).first == trkLink) {
+ theDeposit = i_Deposit;
+ break;
+ }
+ }
+
+ // if the charge has already hit the Diode add it to the deposit
+ if (theDeposit != depositsR_end) {
+ (*theDeposit).second += i_ListOfCharges->charge();
+ } else { // create a new deposit
+ InDetSimData::Deposit deposit(trkLink, i_ListOfCharges->charge());
+ deposits.push_back(deposit);
+ }
+ }
+
+ // add the simdata object to the map:
+ if (real_particle_hit or m_createNoiseSDO) {
+ TrackerDD::SiReadoutCellId roCell{(*i_chargedDiode).second.getReadoutCell()};
+ int strip{roCell.strip()};
+ Identifier id_readout;
+ id_readout = m_detID->strip_id(collection->identify(),strip);
+ (*simDataCollMap)->insert(std::make_pair(id_readout, InDetSimData(deposits, (*i_chargedDiode).second.flag())));
+ }
+ }
+}
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_DigitizationTool.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_DigitizationTool.h
new file mode 100644
index 0000000000000000000000000000000000000000..f34d0619d4bb0538d63e4ded88df576d5b175f44
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_DigitizationTool.h
@@ -0,0 +1,156 @@
+/* -*- C++ -*- */
+
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef FASERSCT_DIGITZATION_SCT_DIGITZATIONTOOL_H
+#define FASERSCT_DIGITZATION_SCT_DIGITZATIONTOOL_H
+/** @file SCT_DigitizationTool.h
+ * @brief Digitize the SCT using an implementation of IPileUpTool
+ * $Id: SCT_DigitizationTool.h,v 1.0 2009-09-22 18:34:42 jchapman Exp $
+ * @author John Chapman - ATLAS Collaboration
+ */
+
+//Base class header
+#include "PileUpTools/PileUpToolBase.h"
+
+// Athena headers
+#include "AthenaKernel/IAthRNGSvc.h"
+#include "HitManagement/TimedHitCollection.h"
+#include "InDetRawData/SCT_RDO_Container.h"
+#include "TrackerReadoutGeometry/SiDetectorElementCollection.h"
+#include "InDetSimData/InDetSimDataCollection.h"
+#include "TrackerSimEvent/FaserSiHitCollection.h"
+#include "PileUpTools/PileUpMergeSvc.h"
+#include "FaserSCT_Digitization/ISCT_FrontEnd.h"
+#include "FaserSCT_Digitization/ISCT_RandomDisabledCellGenerator.h"
+#include "FaserSCT_Digitization/ISCT_SurfaceChargesGenerator.h"
+#include "StoreGate/ReadCondHandleKey.h"
+#include "StoreGate/ReadHandleKey.h"
+#include "StoreGate/WriteHandle.h"
+#include "StoreGate/WriteHandleKey.h"
+
+// Gaudi headers
+#include "GaudiKernel/ServiceHandle.h"
+#include "GaudiKernel/ToolHandle.h"
+
+// STL headers
+#include <limits>
+#include <string>
+
+// Forward declarations
+class ISiChargedDiodesProcessorTool;
+class FaserSCT_ID;
+class SiChargedDiodeCollection;
+
+namespace CLHEP
+{
+ class HepRandomEngine;
+}
+
+class SCT_DigitizationTool : public extends<PileUpToolBase, IPileUpTool>
+{
+public:
+ static const InterfaceID& interfaceID();
+ SCT_DigitizationTool(const std::string& type,
+ const std::string& name,
+ const IInterface* parent);
+ virtual ~SCT_DigitizationTool();
+ /**
+ @brief Called before processing physics events
+ */
+ virtual StatusCode prepareEvent(unsigned int) override final;
+ virtual StatusCode processBunchXing(int bunchXing,
+ SubEventIterator bSubEvents,
+ SubEventIterator eSubEvents) override final;
+ virtual StatusCode mergeEvent() override final;
+
+ virtual StatusCode initialize() override final;
+ virtual StatusCode processAllSubEvents() override final;
+
+protected:
+
+ bool digitizeElement(SiChargedDiodeCollection* chargedDiodes, TimedHitCollection<FaserSiHit>*& thpcsi, CLHEP::HepRandomEngine * rndmEngine) const ; //!
+ void applyProcessorTools(SiChargedDiodeCollection* chargedDiodes, CLHEP::HepRandomEngine * rndmEngine) const; //!
+ void addSDO(SiChargedDiodeCollection* collection, SG::WriteHandle<InDetSimDataCollection>* simDataCollMap) const;
+
+ void storeTool(ISiChargedDiodesProcessorTool* p_processor) {m_diodeCollectionTools.push_back(p_processor);}
+
+private:
+
+ /**
+ @brief initialize the required services
+ */
+ StatusCode initServices();
+ /**
+ @brief Initialize the SCT_FrontEnd AlgTool
+ */
+ StatusCode initFrontEndTool();
+ /**
+ @brief Initialize the SCT_RandomDisabledCellGenerator AlgTool
+ */
+ StatusCode initDisabledCells();
+ /**
+ @brief Initialize the SCT_SurfaceChargesGenerator AlgTool
+ */
+ StatusCode initSurfaceChargesGeneratorTool();
+
+ /** RDO and SDO methods*/
+ /**
+ @brief Create RDOs from the SiChargedDiodeCollection for the current wafer and save to StoreGate
+ @param chDiodeCollection list of the SiChargedDiodes on the current wafer
+ */
+ StatusCode createAndStoreRDO(SiChargedDiodeCollection* chDiodeCollection, SG::WriteHandle<SCT_RDO_Container>* rdoContainer) const;
+ /**
+ @brief Create RDOs from the SiChargedDiodeCollection for the current wafer
+ @param chDiodeCollection list of the SiChargedDiodes on the current wafer
+ */
+ SCT_RDO_Collection* createRDO(SiChargedDiodeCollection* collection) const;
+
+ StatusCode getNextEvent();
+ void digitizeAllHits(SG::WriteHandle<SCT_RDO_Container>* rdoContainer, SG::WriteHandle<InDetSimDataCollection>* simDataCollMap, std::vector<bool>* processedElements, TimedHitCollection<FaserSiHit>* thpcsi, CLHEP::HepRandomEngine * rndmEngine) const; //!< digitize all hits
+ void digitizeNonHits(SG::WriteHandle<SCT_RDO_Container>* rdoContainer, SG::WriteHandle<InDetSimDataCollection>* simDataCollMap, const std::vector<bool>* processedElements, CLHEP::HepRandomEngine * rndmEngine) const; //!< digitize SCT without hits
+
+ /**
+ @brief Called when m_WriteSCT1_RawData is altered. Does nothing, but required by Gaudi.
+ */
+ void SetupRdoOutputType(Property&);
+
+ FloatProperty m_tfix{this, "FixedTime", -999., "Fixed time for Cosmics run selection"};
+ BooleanProperty m_enableHits{this, "EnableHits", true, "Enable hits"};
+ BooleanProperty m_onlyHitElements{this, "OnlyHitElements", false, "Process only elements with hits"};
+ BooleanProperty m_cosmicsRun{this, "CosmicsRun", false, "Cosmics run selection"};
+ BooleanProperty m_randomDisabledCells{this, "RandomDisabledCells", false, "Use Random disabled cells, default no"};
+ BooleanProperty m_createNoiseSDO{this, "CreateNoiseSDO", false, "Create SDOs for strips with only noise hits (huge increase in SDO collection size"};
+ BooleanProperty m_WriteSCT1_RawData{this, "WriteSCT1_RawData", false, "Write out SCT1_RawData rather than SCT3_RawData"};
+
+ BooleanProperty m_onlyUseContainerName{this, "OnlyUseContainerName", true, "Don't use the ReadHandleKey directly. Just extract the container name from it."};
+ SG::ReadHandleKey<FaserSiHitCollection> m_hitsContainerKey{this, "InputObjectName", "SCT_Hits", "Input HITS collection name"};
+ std::string m_inputObjectName{""};
+
+ SG::WriteHandleKey<SCT_RDO_Container> m_rdoContainerKey{this, "OutputObjectName", "SCT_RDOs", "Output Object name"};
+ SG::WriteHandle<SCT_RDO_Container> m_rdoContainer; //!< RDO container handle
+ SG::WriteHandleKey<InDetSimDataCollection> m_simDataCollMapKey{this, "OutputSDOName", "SCT_SDO_Map", "Output SDO container name"};
+ SG::WriteHandle<InDetSimDataCollection> m_simDataCollMap; //!< SDO Map handle
+ SG::ReadCondHandleKey<TrackerDD::SiDetectorElementCollection> m_SCTDetEleCollKey{this, "SCTDetEleCollKey", "SCT_DetectorElementCollection", "Key of SiDetectorElementCollection for SCT"};
+
+ ToolHandle<ISCT_FrontEnd> m_sct_FrontEnd{this, "FrontEnd", "SCT_FrontEnd", "Handle the Front End Electronic tool"};
+ ToolHandle<ISCT_SurfaceChargesGenerator> m_sct_SurfaceChargesGenerator{this, "SurfaceChargesGenerator", "SCT_SurfaceChargesGenerator", "Choice of using a more detailed charge drift model"};
+ ToolHandle<ISCT_RandomDisabledCellGenerator> m_sct_RandomDisabledCellGenerator{this, "RandomDisabledCellGenerator", "SCT_RandomDisabledCellGenerator", ""};
+ ServiceHandle<IAthRNGSvc> m_rndmSvc{this, "RndmSvc", "AthRNGSvc", ""}; //!< Random number service
+ ServiceHandle <PileUpMergeSvc> m_mergeSvc{this, "MergeSvc", "PileUpMergeSvc", "Merge service used in Pixel & SCT digitization"}; //!
+
+ const FaserSCT_ID* m_detID{nullptr}; //!< Handle to the ID helper
+ TimedHitCollection<FaserSiHit>* m_thpcsi{nullptr};
+ std::list<ISiChargedDiodesProcessorTool*> m_diodeCollectionTools;
+ std::vector<bool> m_processedElements; //!< vector of processed elements - set by digitizeHits() */
+ std::vector<FaserSiHitCollection*> m_hitCollPtrs;
+};
+
+static const InterfaceID IID_ISCT_DigitizationTool("SCT_DigitizationTool", 1, 0);
+inline const InterfaceID& SCT_DigitizationTool::interfaceID() {
+ return IID_ISCT_DigitizationTool;
+}
+
+#endif // FASERSCT_DIGITZATION_SCT_DIGITZATIONTOOL_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_FrontEnd.cxx b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_FrontEnd.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..50d06c46c2dceb839ae362a89b6eea8c92320a38
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_FrontEnd.cxx
@@ -0,0 +1,1044 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "SCT_FrontEnd.h"
+
+// Random number
+#include "CLHEP/Random/RandFlat.h"
+#include "CLHEP/Random/RandGaussZiggurat.h" // for RandGaussZiggurat
+#include "CLHEP/Random/RandPoisson.h"
+#include "CLHEP/Random/RandomEngine.h"
+
+#include "FaserSiDigitization/SiHelper.h"
+#include "FaserSCT_Digitization/ISCT_Amp.h"
+
+#include "TrackerIdentifier/FaserSCT_ID.h"
+
+#include "TrackerReadoutGeometry/SCT_DetectorManager.h"
+#include "TrackerReadoutGeometry/SCT_ModuleSideDesign.h"
+
+// C++ Standard Library
+#include <algorithm>
+#include <iostream>
+
+// #define SCT_DIG_DEBUG
+
+using namespace std;
+using namespace TrackerDD;
+
+// constructor
+SCT_FrontEnd::SCT_FrontEnd(const std::string& type, const std::string& name, const IInterface* parent)
+ : base_class(type, name, parent) {
+}
+
+// ----------------------------------------------------------------------
+// Initialize
+// ----------------------------------------------------------------------
+StatusCode SCT_FrontEnd::initialize() {
+ if (m_NoiseOn and (not m_analogueNoiseOn)) {
+ ATH_MSG_FATAL("AnalogueNoiseOn/m_analogueNoiseOn should be true if NoiseOn/m_NoiseOn is true.");
+ return StatusCode::FAILURE;
+ }
+
+ ATH_MSG_DEBUG("SCT_FrontEnd::initialize()");
+ // Get SCT helper
+ ATH_CHECK(detStore()->retrieve(m_sct_id, "FaserSCT_ID"));
+ // Get SCT detector manager
+ ATH_CHECK(detStore()->retrieve(m_SCTdetMgr, "SCT"));
+ // Get the amplifier tool
+ ATH_CHECK(m_sct_amplifier.retrieve());
+ ATH_MSG_DEBUG("SCT Amplifier tool located ");
+
+ // Get the SCT_ReadCaliDataSvc
+ if (m_useCalibData) {
+ // ATH_CHECK(m_ReadCalibChipDataTool.retrieve());
+ // ATH_MSG_DEBUG("CalibChipData Service located ");
+ ATH_MSG_FATAL("Use of calibration data not supported.");
+ } else {
+ // m_ReadCalibChipDataTool.disable();
+ }
+
+ // Get the maximum number of strips of any module
+ m_strip_max = m_SCTdetMgr->numerology().maxNumStrips();
+
+ constexpr float fC = 6242.2;
+ m_Threshold = m_Threshold * fC;
+
+#ifdef SCT_DIG_DEBUG
+ ATH_MSG_INFO("\tNoise factors:");
+ ATH_MSG_INFO("\tBarrel = " << m_NoiseBarrel << " Outer Barrel = " << m_NoiseBarrel3 <<
+ " EC, inners = " << m_NoiseInners << " EC, middles = " << m_NoiseMiddles <<
+ " EC, short middles = " << m_NoiseShortMiddles << " EC, outers = " << m_NoiseOuters);
+ ATH_MSG_INFO("\tThreshold=" << m_Threshold << " fC, time of Threshold=" << m_timeOfThreshold);
+#endif
+
+ ATH_MSG_INFO("m_Threshold (Threshold) = " << m_Threshold);
+ ATH_MSG_INFO("m_timeOfThreshold (TimeOfThreshold) = " << m_timeOfThreshold);
+ ATH_MSG_INFO("m_data_compression_mode (DataCompressionMode) = " << m_data_compression_mode);
+ ATH_MSG_INFO("m_data_readout_mode (DataReadOutMode) = " << m_data_readout_mode);
+
+ return StatusCode::SUCCESS;
+}
+
+StatusCode SCT_FrontEnd::finalize() {
+#ifdef SCT_DIG_DEBUG
+ ATH_MSG_INFO("SCT_FrontEnd::finalize()");
+#endif
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+// Init the class variable vectors
+// ----------------------------------------------------------------------
+StatusCode SCT_FrontEnd::initVectors(int strips, SCT_FrontEndData& data) const {
+ data.m_Offset.assign(strips, 0.0);
+ data.m_GainFactor.assign(strips, 0.0);
+ data.m_NoiseFactor.assign(strips, 0.0);
+
+ if (m_data_compression_mode == 1 and m_data_readout_mode == 0) {
+ data.m_Analogue[0].reserve(1);
+ data.m_Analogue[1].reserve(strips);
+ } else if (m_data_compression_mode == 2 and m_data_readout_mode == 0) {
+ data.m_Analogue[0].reserve(strips);
+ data.m_Analogue[1].reserve(strips);
+ } else if (m_data_compression_mode == 3 or m_data_readout_mode == 1) {
+ data.m_Analogue[0].reserve(strips);
+ data.m_Analogue[1].reserve(strips);
+ data.m_Analogue[2].reserve(strips);
+ } else {
+ return StatusCode::FAILURE;
+ }
+
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+// prepare gain and offset for the strips for a given module
+// ----------------------------------------------------------------------
+StatusCode SCT_FrontEnd::prepareGainAndOffset(SiChargedDiodeCollection& collection, /*const Identifier& moduleId,*/ CLHEP::HepRandomEngine * rndmEngine, SCT_FrontEndData& data) const {
+ // now we need to generate gain and offset channel by channel: some algebra
+ // for generation of partially correlated random numbers
+ float W = m_OGcorr * m_GainRMS * m_Ospread / (m_GainRMS * m_GainRMS - m_Ospread * m_Ospread);
+ float A = 4 * W * W + 1.0;
+ float x1 = (A - sqrt(A)) / (2.0 * A);
+ float sinfi = sqrt(x1);
+ float cosfi = sqrt(1.0 - x1);
+
+ sinfi = sinfi * m_OGcorr / fabs(m_OGcorr);
+ float S = m_GainRMS * m_GainRMS + m_Ospread * m_Ospread;
+ float D = (m_GainRMS * m_GainRMS - m_Ospread * m_Ospread) / (cosfi * cosfi - sinfi * sinfi);
+ float S1 = sqrt((S + D) * 0.5);
+ float S2 = sqrt((S - D) * 0.5);
+ float Noise = 0;
+ int mode = 1;
+
+ // To set noise values for different module types, barrel, EC, inners, middles, short middles, and outers
+ if (m_analogueNoiseOn) {
+ // if (m_sct_id->barrel_ec(moduleId) == 0) { // barrel_ec=0 corresponds to barrel
+ // if (m_sct_id->layer_disk(moduleId) == 3) { // outer barrel layer 10 degrees warmer
+ Noise = m_NoiseBarrel3;
+ // } else {
+ // Noise = m_NoiseBarrel;
+ // }
+ // } else {
+ // int moduleType = m_sct_id->eta_module(moduleId);
+ // switch (moduleType) { // eta = 0, 1, or 2 corresponds to outers, middles and inners?! (at least in the offline world)
+ // case 0: {
+ // Noise = m_NoiseOuters;
+ // break;
+ // }
+ // case 1: {
+ // if (m_sct_id->layer_disk(moduleId) == 7) {
+ // Noise = m_NoiseShortMiddles;
+ // } else {
+ // Noise = m_NoiseMiddles;
+ // }
+ // break;
+ // }
+ // case 2: {
+ // Noise = m_NoiseInners;
+ // break;
+ // }
+ // default: {
+ // Noise = m_NoiseBarrel;
+ // ATH_MSG_ERROR("moduleType(eta): " << moduleType << " unknown, using barrel");
+ // }
+ // }// end of switch structure
+ // }
+ }
+
+ // Loop over collection and setup gain/offset/noise for the hit and neighbouring strips
+ SiChargedDiodeIterator i_chargedDiode = collection.begin();
+ SiChargedDiodeIterator i_chargedDiode_end = collection.end();
+
+ for (; i_chargedDiode != i_chargedDiode_end; ++i_chargedDiode) {
+ SiChargedDiode diode = (*i_chargedDiode).second;
+ // should be const as we aren't trying to change it here - but getReadoutCell() is not a const method...
+ unsigned int flagmask = diode.flag() & 0xFE;
+ // Get the flag for this diode ( if flagmask = 1 If diode is disconnected/disabled skip it)
+ if (!flagmask) { // If the diode is OK (not flagged)
+ const SiReadoutCellId roCell = diode.getReadoutCell();
+ if (roCell.isValid()) {
+ int strip = roCell.strip();
+ int i = std::max(strip - 1, 0);
+ int i_end = std::min(strip + 2, m_strip_max.load());
+
+ // loop over strips
+ for (; i < i_end; i++) {
+ // Need to check if strip is already setup
+ if (data.m_Analogue[1][i] <= 0.0) {
+ float g = CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0.0, S1);
+ float o = CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0.0, S2);
+
+ data.m_GainFactor[i] = 1.0 + (cosfi * g + sinfi * o);
+ data.m_Offset[i] = (cosfi * o - sinfi * g);
+ data.m_NoiseFactor[i] = Noise * mode;
+
+ // Fill the noise and offset values into the Analogue
+ if (m_data_compression_mode == 1 and m_data_readout_mode == 0) { // level mode x1x
+ data.m_Analogue[1][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+ } else if (m_data_compression_mode == 2 and m_data_readout_mode == 0) { // edge mode 01x
+ data.m_Analogue[0][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+ data.m_Analogue[1][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+ } else if (m_data_compression_mode == 3 or m_data_readout_mode == 1) { // any hit mode xxx or expanded read out mode
+ data.m_Analogue[0][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+ data.m_Analogue[1][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+ data.m_Analogue[2][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+ }
+ }
+ }
+ }
+ }
+ }
+
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+// prepare gain and offset for the strips for a given module using
+// Cond Db data to get the chip calibration data
+// ----------------------------------------------------------------------
+// StatusCode SCT_FrontEnd::prepareGainAndOffset(SiChargedDiodeCollection& collection, int side, const Identifier& moduleId, CLHEP::HepRandomEngine * rndmEngine, SCT_FrontEndData& data) const {
+// // Get chip data from calib DB
+// std::vector<float> gainByChipVect = m_ReadCalibChipDataTool->getNPtGainData(moduleId, side, "GainByChip");
+// std::vector<float> gainRMSByChipVect = m_ReadCalibChipDataTool->getNPtGainData(moduleId, side, "GainRMSByChip");
+// std::vector<float> offsetByChipVect = m_ReadCalibChipDataTool->getNPtGainData(moduleId, side, "OffsetByChip");
+// std::vector<float> offsetRMSByChipVect = m_ReadCalibChipDataTool->getNPtGainData(moduleId, side, "OffsetRMSByChip");
+// std::vector<float> noiseByChipVect(6, 0.0);
+
+// if (m_analogueNoiseOn) { // Check if noise should be on or off
+// noiseByChipVect = m_ReadCalibChipDataTool->getNPtGainData(moduleId, side, "NoiseByChip");
+// }
+
+// // Need to check if empty, most should have data, but a few old DEAD modules don't
+// if (gainByChipVect.empty() or noiseByChipVect.empty()) {
+// ATH_MSG_DEBUG("No calibration data in cond DB for module " << moduleId << " using JO values");
+// if (StatusCode::SUCCESS != prepareGainAndOffset(collection, moduleId, rndmEngine, data)) {
+// return StatusCode::FAILURE;
+// } else {
+// return StatusCode::SUCCESS;
+// }
+// }
+
+// // Don't really need to set up values for each chip...
+// float gainMeanValue = meanValue(gainByChipVect);
+// if (gainMeanValue < 0.0) {
+// ATH_MSG_DEBUG("All chip gain values are 0 for module " << moduleId << " using JO values");
+// if (StatusCode::SUCCESS != prepareGainAndOffset(collection, moduleId, rndmEngine, data)) {
+// return StatusCode::FAILURE;
+// } else {
+// return StatusCode::SUCCESS;
+// }
+// }
+
+// std::vector<float> gain(6, 0.0);
+// std::vector<float> offset(6, 0.0);
+// std::vector<float> S1(6, 0.0);
+// std::vector<float> S2(6, 0.0);
+// std::vector<float> sinfi(6, 0.0);
+// std::vector<float> cosfi(6, 0.0);
+// float gainRMS = 0.0;
+// float offsetRMS = 0.0;
+
+// for (int i = 0; i < 6; ++i) {
+// // Some very few chips have 0 values, dead/bypassed/etc, so check and use some fixed values instead
+// if (gainByChipVect[i] > 0.1) {
+// gain[i] = gainByChipVect[i] / gainMeanValue;
+// offset[i] = offsetByChipVect[i] / m_Threshold;
+// gainRMS = gainRMSByChipVect[i] / gainMeanValue;
+// offsetRMS = offsetRMSByChipVect[i] / m_Threshold;
+// } else {
+// gain[i] = 55.0 / gainMeanValue;
+// offset[i] = 42.0 / m_Threshold;
+// gainRMS = 1.3 / gainMeanValue;
+// offsetRMS = 2.0 / m_Threshold;
+// }
+
+// float W = m_OGcorr * gainRMS * offsetRMS / (gainRMS * gainRMS - offsetRMS * offsetRMS);
+// float A = 4 * W * W + 1.0;
+// float x1 = (A - sqrt(A)) / (2.0 * A);
+// sinfi[i] = sqrt(x1);
+// cosfi[i] = sqrt(1.0 - x1);
+// sinfi[i] = sinfi[i] * m_OGcorr / fabs(m_OGcorr);
+// float S = gainRMS * gainRMS + offsetRMS * offsetRMS;
+// float D = (gainRMS * gainRMS - offsetRMS * offsetRMS) / (cosfi[i] * cosfi[i] - sinfi[i] * sinfi[i]);
+// S1[i] = sqrt((S + D) / 2.0);
+// S2[i] = sqrt((S - D) / 2.0);
+// }
+
+// // Loop over collection and setup gain/offset/noise for the hit and neighbouring strips
+// SiChargedDiodeIterator i_chargedDiode = collection.begin();
+// SiChargedDiodeIterator i_chargedDiode_end = collection.end();
+
+// for (; i_chargedDiode != i_chargedDiode_end; ++i_chargedDiode) {
+// SiChargedDiode diode = (*i_chargedDiode).second;
+// // should be const as we aren't trying to change it here - but getReadoutCell() is not a const method...
+// unsigned int flagmask = diode.flag() & 0xFE;
+// // Get the flag for this diode ( if flagmask = 1 If diode is disconnected/disabled skip it)
+// if (!flagmask) { // If the diode is OK (not flagged)
+// const SiReadoutCellId roCell = diode.getReadoutCell();
+
+// if (roCell.isValid()) {
+// int strip = roCell.strip();
+// int i = std::max(strip - 1, 0);
+// int i_end = std::min(strip + 2, m_strip_max.load());
+
+// // loop over strips
+// for (; i < i_end; i++) {
+// // Need to check if strip is already setup
+// if (data.m_Analogue[1][i] <= 0.0) {
+// // Values depends on which chip the strip is on (complex when strip is on chip edge)
+// int chip = i / 128;
+// float g = CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0.0, S1[chip]);
+// float o = CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0.0, S2[chip]);
+
+// data.m_GainFactor[i] = gain[chip] + (cosfi[chip] * g + sinfi[chip] * o);
+// data.m_Offset[i] = offset[chip] + (cosfi[chip] * o - sinfi[chip] * g);
+// data.m_NoiseFactor[i] = noiseByChipVect[chip];
+
+// // Fill the noise and offset values into the Analogue
+// if (m_data_compression_mode == 1 and m_data_readout_mode == 0) { // level mode x1x
+// data.m_Analogue[1][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+// } else if (m_data_compression_mode == 2 and m_data_readout_mode == 0) { // edge mode 01x
+// data.m_Analogue[0][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+// data.m_Analogue[1][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+// } else if (m_data_compression_mode == 3 or m_data_readout_mode == 1) { // any hit mode xxx or expanded read out mode
+// data.m_Analogue[0][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+// data.m_Analogue[1][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+// data.m_Analogue[2][i] = data.m_Offset[i] + data.m_NoiseFactor[i] * CLHEP::RandGaussZiggurat::shoot(rndmEngine);
+// }
+// }
+// }
+// }
+// }
+// }
+
+// return StatusCode::SUCCESS;
+// }
+
+// ----------------------------------------------------------------------
+//
+// ----------------------------------------------------------------------
+StatusCode SCT_FrontEnd::randomNoise(SiChargedDiodeCollection& collection, /*const Identifier& moduleId,*/ CLHEP::HepRandomEngine * rndmEngine, SCT_FrontEndData& data) const {
+ // Add random noise
+
+ double occupancy = 0.0;
+ double NoiseOccupancy = 0.0;
+ float Noise = 0.0;
+ int nNoisyStrips = 0;
+ double mode = 1.;
+
+ const bool noise_expanded_mode = (m_data_compression_mode == 3 and m_data_readout_mode == 1);
+
+ // Will give 3 times as much noise occupancy if running in any hit expanded mode
+ if (noise_expanded_mode) {
+ mode = 3.;
+ }
+
+ // Sets fixed noise occupancy values for different module types, barrel, EC,
+ // inners, middles
+ // short middles, and outers
+// if (m_sct_id->barrel_ec(moduleId) == 0) { // barrel_ec=0 corresponds to barrel
+// if (m_sct_id->layer_disk(moduleId) == 3) { // outer barrel layer 10 degrees warmer
+ NoiseOccupancy = m_NOBarrel3;
+ Noise = m_NoiseBarrel3;
+// } else {
+// NoiseOccupancy = m_NOBarrel;
+// Noise = m_NoiseBarrel;
+// }
+// } else {
+// int moduleType = m_sct_id->eta_module(moduleId);
+// switch (moduleType) { // eta = 0, 1, or 2 corresponds to outers, middles and inners?! (at least in the offline world)
+// case 0: {
+// NoiseOccupancy = m_NOOuters;
+// Noise = m_NoiseOuters;
+// break;
+// }
+// case 1: {
+// if (m_sct_id->layer_disk(moduleId) == 7) {
+// NoiseOccupancy = m_NOShortMiddles;
+// Noise = m_NoiseShortMiddles;
+// } else {
+// NoiseOccupancy = m_NOMiddles;
+// Noise = m_NoiseMiddles;
+// }
+// break;
+// }
+// case 2: {
+// NoiseOccupancy = m_NOInners;
+// Noise = m_NoiseInners;
+// break;
+// }
+// default: {
+// NoiseOccupancy = m_NOBarrel;
+// Noise = m_NoiseBarrel;
+// ATH_MSG_ERROR("moduleType(eta): " << moduleType << " unknown, using barrel");
+// }
+// }// end of switch structure
+// }
+
+ // Calculate the number of "free strips"
+ int nEmptyStrips = 0;
+ std::vector<int> emptyStrips;
+ emptyStrips.reserve(m_strip_max);
+ for (int i = 0; i < m_strip_max; i++) {
+ if (data.m_StripHitsOnWafer[i] == 0) {
+ emptyStrips.push_back(i);
+ ++nEmptyStrips;
+ }
+ }
+
+ if (nEmptyStrips != 0) {
+ // Should randomize the fixed NO values, so we get some differences per
+ // wafer
+ occupancy = CLHEP::RandGaussZiggurat::shoot(rndmEngine, NoiseOccupancy, NoiseOccupancy * 0.1);
+
+ // Modify the occupancy if threshold is not 1.0 fC
+ if (m_Threshold > 6242.3 or m_Threshold < 6242.1) {
+ const float fC = 6242.2;
+ occupancy = occupancy * exp(-(0.5 / (Noise * Noise) * (m_Threshold * m_Threshold - fC * fC)));
+ }
+ nNoisyStrips = CLHEP::RandPoisson::shoot(rndmEngine, m_strip_max * occupancy * mode);
+
+ // Check and adapt the number of noisy strips to the number of free strips
+ if (nEmptyStrips < nNoisyStrips) {
+ nNoisyStrips = nEmptyStrips;
+ }
+
+ // Find random strips to get noise hits
+ for (int i = 0; i < nNoisyStrips; i++) {
+ int index = CLHEP::RandFlat::shootInt(rndmEngine, nEmptyStrips - i); // strip == 10, 12 free strips
+ // have vector [10000100100200211001] 20 strips
+ int strip = emptyStrips.at(index);
+ emptyStrips.erase(emptyStrips.begin()+index); // Erase it not to use it again
+ if (data.m_StripHitsOnWafer[strip]!=0) {
+ ATH_MSG_ERROR(index << "-th empty strip, strip " << strip << " should be empty but is not empty! Something is wrong!");
+ }
+ data.m_StripHitsOnWafer[strip] = 3; // !< Random Noise hit
+ // Add tbin info to noise diode
+ if (noise_expanded_mode) { // !< if any hit mode, any time bin otherwise fixed tbin=2
+ int noise_tbin = CLHEP::RandFlat::shootInt(rndmEngine, 3);
+ // !< random number 0, 1 or 2
+ if (noise_tbin == 0) {
+ noise_tbin = 4; // !< now 1,2 or 4
+ }
+ if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, noise_tbin)) {
+ ATH_MSG_ERROR("Can't add noise hit diode to collection");
+ }
+ } else {
+ if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, 2)) {
+ ATH_MSG_ERROR("Can't add noise hit diode to collection");
+ }
+ }
+ }
+ }
+
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+//
+// ----------------------------------------------------------------------
+// StatusCode SCT_FrontEnd::randomNoise(SiChargedDiodeCollection& collection, const Identifier& moduleId, int side, CLHEP::HepRandomEngine * rndmEngine, SCT_FrontEndData& data) const {
+// const int n_chips = 6;
+// const int chipStripmax = m_strip_max / n_chips;
+// std::vector<float> NOByChipVect(n_chips, 0.0);
+// std::vector<float> ENCByChipVect(n_chips, 0.0);
+// std::vector<int> nNoisyStrips(n_chips, 0);
+// double mode = 1.;
+
+// const bool noise_expanded_mode = (m_data_compression_mode == 3 and m_data_readout_mode == 1);
+
+// // Will give 3 times as much noise occupancy if running in any hit expanded mode
+// if (noise_expanded_mode) {
+// mode = 3.;
+// }
+
+// // Get chip data from calib DB
+// NOByChipVect = m_ReadCalibChipDataTool->getNoiseOccupancyData(moduleId, side, "OccupancyByChip");
+// ENCByChipVect = m_ReadCalibChipDataTool->getNPtGainData(moduleId, side, "NoiseByChip");
+
+// // Need to check if empty, most should have data, but a few old DEAD modules don't, and 9C...
+// if (NOByChipVect.empty()) {
+// ATH_MSG_DEBUG("No calibration data in cond DB for module " << moduleId << " using JO values");
+// if (StatusCode::SUCCESS != randomNoise(collection, moduleId, rndmEngine, data)) {
+// return StatusCode::FAILURE;
+// } else {
+// return StatusCode::SUCCESS;
+// }
+// } else {
+// for (int i = 0; i < n_chips; i++) {
+// // A 0 value can mean two things now, chip out of config for long time and no value was uploaded
+// // or its short middles and inners and the value is for all purposes 0! so ok.
+
+// // Modify the occupancy if threshold is not 1.0 fC
+// if (m_Threshold > 6242.3 or m_Threshold < 6242.1) {
+// constexpr float fC = 6242.2;
+// NOByChipVect[i] = NOByChipVect[i] * exp(-(0.5 / (ENCByChipVect[i]*ENCByChipVect[i]) * (m_Threshold*m_Threshold - fC*fC)));
+// }
+
+// nNoisyStrips[i] = CLHEP::RandPoisson::shoot(rndmEngine, chipStripmax * NOByChipVect[i] * mode);
+// }
+// }
+
+// // Loop over the chips on the wafer
+// for (int chip_index = 0; chip_index < n_chips; ++chip_index) {
+// int chip_strip_offset = chipStripmax * chip_index; // First strip number on chip
+
+// // Calculate the number of "free strips" on this chip
+// int nEmptyStripsOnChip = 0;
+// std::vector<int> emptyStripsOnChip;
+// emptyStripsOnChip.reserve(chipStripmax);
+// for (int i = 0; i < chipStripmax; i++) {
+// if (data.m_StripHitsOnWafer[i + chip_strip_offset] == 0) {
+// emptyStripsOnChip.push_back(i);
+// ++nEmptyStripsOnChip;
+// }
+// }
+
+// // if no empty strips on chip do nothing
+// if (nEmptyStripsOnChip != 0) {
+// // Check and adapt the number of noisy strips to the number of free strips
+// if (nEmptyStripsOnChip < nNoisyStrips[chip_index]) {
+// nNoisyStrips[chip_index] = nEmptyStripsOnChip;
+// }
+
+// // Find random strips to get noise hits
+// for (int i = 0; i < nNoisyStrips[chip_index]; i++) {
+// int index = CLHEP::RandFlat::shootInt(rndmEngine, nEmptyStripsOnChip - i);
+// int strip_on_chip = emptyStripsOnChip.at(index);
+// emptyStripsOnChip.erase(emptyStripsOnChip.begin()+index); // Erase it not to use it again
+// int strip = strip_on_chip + chip_strip_offset;
+// if (data.m_StripHitsOnWafer[strip]!=0) {
+// ATH_MSG_ERROR(index << "-th empty strip, strip " << strip << " should be empty but is not empty! Something is wrong!");
+// }
+// data.m_StripHitsOnWafer[strip] = 3; // !< Random Noise hit
+// // Add tbin info to noise diode
+// if (noise_expanded_mode) { // !< if any hit mode, any time bin
+// // !< otherwise fixed tbin=2
+// int noise_tbin = CLHEP::RandFlat::shootInt(rndmEngine, 3);
+// // !< random number 0, 1 or 2
+// if (noise_tbin == 0) {
+// noise_tbin = 4; // !< now 1, 2 or 4
+// }
+// if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, noise_tbin)) {
+// ATH_MSG_ERROR("Can't add noise hit diode to collection");
+// }
+// } else {
+// if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, 2)) {
+// ATH_MSG_ERROR("Can't add noise hit diode to collection");
+// }
+// }
+// }
+// }
+// }
+
+// return StatusCode::SUCCESS;
+// }
+
+// ----------------------------------------------------------------------
+// process the collection of pre digits this will need to go through
+// all single-strip pre-digits calculate the amplifier response add noise
+// (this could be moved elsewhere later) apply threshold do clustering
+// ----------------------------------------------------------------------
+void SCT_FrontEnd::process(SiChargedDiodeCollection& collection, CLHEP::HepRandomEngine * rndmEngine) const {
+ // get SCT module side design and check it
+ const SCT_ModuleSideDesign *p_design = dynamic_cast<const SCT_ModuleSideDesign*>(&(collection.design()));
+
+ if (p_design==nullptr) {
+ return;
+ }
+
+ SCT_FrontEndData data;
+
+ // Check number of strips in design and from manager(max number of strips on any module)
+ // The design value should always be equal or lower than the manager one
+ // However, no resising is now done in case of a lower value
+ const int strip_max = p_design->cells();
+ // Init vectors
+ if (StatusCode::SUCCESS != initVectors(strip_max, data)) {
+ ATH_MSG_ERROR("Can't resize class variable vectors");
+ return;
+ }
+
+ // Contains strip hit info, reset to 0 for each wafer processed
+ data.m_StripHitsOnWafer.assign(m_strip_max, 0);
+
+ // Containes the charge for each bin on each hit strip
+ if (m_data_compression_mode == 1 and m_data_readout_mode == 0) {
+ for (int i = 0; i < m_strip_max; ++i) {
+ data.m_Analogue[1][i] = 0.0;
+ }
+ } else if (m_data_compression_mode == 2 and m_data_readout_mode == 0) {
+ for (int i = 0; i < m_strip_max; ++i) {
+ data.m_Analogue[0][i] = 0.0;
+ data.m_Analogue[1][i] = 0.0;
+ }
+ } else if (m_data_compression_mode == 3 or m_data_readout_mode == 1) {
+ for (int i = 0; i < m_strip_max; ++i) {
+ data.m_Analogue[0][i] = 0.0;
+ data.m_Analogue[1][i] = 0.0;
+ data.m_Analogue[2][i] = 0.0;
+ }
+ }
+
+ // Get wafer, moduleId and side
+// Identifier waferId = collection.identify();
+// Identifier moduleId = m_sct_id->module_id(waferId);
+// const int side = m_sct_id->side(waferId);
+
+ // Check if collection empty
+ if (not collection.empty()) {
+ // Setup gain/offset/noise to the hit and neighbouring strips
+ if (m_useCalibData) { // Use calib cond DB data
+ ATH_MSG_FATAL("Use of calibration data not supported.");
+ // if (StatusCode::SUCCESS != prepareGainAndOffset(collection, side, moduleId, rndmEngine, data)) {
+ // ATH_MSG_ERROR("\tCan't prepare Gain and Offset");
+ // }
+ } else { // Use JO values
+ if (StatusCode::SUCCESS != prepareGainAndOffset(collection, /*moduleId,*/ rndmEngine, data)) {
+ ATH_MSG_ERROR("\tCan't prepare Gain and Offset");
+ }
+ }
+
+ if (StatusCode::SUCCESS != doSignalChargeForHits(collection, data)) {
+ ATH_MSG_ERROR("\tCan't doSignalChargeForHits");
+ }
+
+ if (StatusCode::SUCCESS != doThresholdCheckForRealHits(collection, data)) {
+ ATH_MSG_ERROR("\tCan't doThresholdCheckForRealHits");
+ }
+
+ if (StatusCode::SUCCESS != doThresholdCheckForCrosstalkHits(collection, data)) {
+ ATH_MSG_ERROR("\tCan't doThresholdCheckForCrosstalkHits");
+ }
+ }
+
+ if (m_NoiseOn) {
+ if (m_useCalibData) { // Check if using DB or not
+ ATH_MSG_FATAL("Use of calibration data not supported.");
+ // if (StatusCode::SUCCESS != randomNoise(collection, moduleId, side, rndmEngine, data)) {
+ // ATH_MSG_ERROR("\tCan't do random noise on wafer?!");
+ // }
+ } else { // Use JO fixed values
+ if (StatusCode::SUCCESS != randomNoise(collection, /*moduleId,*/ rndmEngine, data)) {
+ ATH_MSG_ERROR("\tCan't do random noise on wafer?!");
+ }
+ }
+ }
+
+ // Check for strips above threshold and do clustering
+ if (StatusCode::SUCCESS != doClustering(collection, data)) {
+ ATH_MSG_ERROR("\tCan't cluster the hits?!");
+ }
+}
+
+StatusCode SCT_FrontEnd::doSignalChargeForHits(SiChargedDiodeCollection& collection, SCT_FrontEndData& data) const {
+ typedef SiTotalCharge::list_t list_t;
+
+ // *****************************************************************************
+ // Loop over the diodes (strips ) and for each of them define the total signal
+ // *****************************************************************************
+
+ // set up number of needed bins depending on the compression mode
+ short bin_max = 0;
+ if (m_data_readout_mode == 0) {
+ bin_max = m_data_compression_mode;
+ } else {
+ bin_max = 3;
+ }
+
+ std::vector<float> response(bin_max);
+
+ SiChargedDiodeIterator i_chargedDiode = collection.begin();
+ SiChargedDiodeIterator i_chargedDiode_end = collection.end();
+ for (; i_chargedDiode != i_chargedDiode_end; ++i_chargedDiode) {
+ SiChargedDiode diode = (*i_chargedDiode).second;
+ // should be const as we aren't trying to change it here - but getReadoutCell() is not a const method...
+ unsigned int flagmask = diode.flag() & 0xFE;
+ // Get the flag for this diode ( if flagmask = 1 If diode is disconnected/disabled skip it)
+ if (!flagmask) { // If the diode is OK (not flagged)
+ const SiReadoutCellId roCell = diode.getReadoutCell();
+
+ if (roCell.isValid()) {
+ int strip = roCell.strip();
+
+ const list_t &ChargesOnStrip = diode.totalCharge().chargeComposition();
+
+ if (m_data_compression_mode == 1 and m_data_readout_mode == 0) { // level mode x1x
+ // Amplifier response
+ data.m_Analogue[1][strip] += data.m_GainFactor[strip] * m_sct_amplifier->response(ChargesOnStrip, m_timeOfThreshold);
+
+ // Add Crosstalk signal for neighboring strip
+ response[0] = m_sct_amplifier->crosstalk(ChargesOnStrip, m_timeOfThreshold);
+ if (strip + 1 < m_strip_max) {
+ data.m_Analogue[1][strip + 1] += data.m_GainFactor[strip + 1] * response[0];
+ }
+ if (strip > 0) {
+ data.m_Analogue[1][strip - 1] += data.m_GainFactor[strip - 1] * response[0];
+ }
+ } else if (m_data_compression_mode == 2 and m_data_readout_mode == 0) { // edge mode 01x
+ // Amplifier response
+ m_sct_amplifier->response(ChargesOnStrip, m_timeOfThreshold, response);
+ for (short bin = 0; bin < bin_max; ++bin) {
+ data.m_Analogue[bin][strip] += data.m_GainFactor[strip] * response[bin];
+ }
+ // Add Crosstalk signal for neighboring strip
+ m_sct_amplifier->crosstalk(ChargesOnStrip, m_timeOfThreshold, response);
+ for (short bin = 0; bin < bin_max; ++bin) {
+ if (strip + 1 < m_strip_max) {
+ data.m_Analogue[bin][strip + 1] += data.m_GainFactor[strip + 1] * response[bin];
+ }
+ if (strip > 0) {
+ data.m_Analogue[bin][strip - 1] += data.m_GainFactor[strip - 1] * response[bin];
+ }
+ }
+ } else if (m_data_compression_mode == 3 or m_data_readout_mode == 1) { // any hit mode xxx or expanded read out mode
+ // Amplifier response
+ m_sct_amplifier->response(ChargesOnStrip, m_timeOfThreshold, response);
+ for (short bin = 0; bin < bin_max; ++bin) {
+ data.m_Analogue[bin][strip] += data.m_GainFactor[strip] * response[bin];
+ }
+ // Add Crosstalk signal for neighboring strip
+ m_sct_amplifier->crosstalk(ChargesOnStrip, m_timeOfThreshold, response);
+ for (short bin = 0; bin < bin_max; ++bin) {
+ if (strip + 1 < m_strip_max) {
+ data.m_Analogue[bin][strip + 1] += data.m_GainFactor[strip + 1] * response[bin];
+ }
+ if (strip > 0) {
+ data.m_Analogue[bin][strip - 1] += data.m_GainFactor[strip - 1] * response[bin];
+ }
+ }
+ }
+ } else { // if roCell not valid
+ ATH_MSG_WARNING("\t Cannot get the cell ");
+ }
+ } else {// If diode is disconnected/disabled skip it
+ ATH_MSG_WARNING("\tDisabled or disconnected diode (strip)");
+ }
+ }
+ return StatusCode::SUCCESS;
+}
+
+StatusCode SCT_FrontEnd::doThresholdCheckForRealHits(SiChargedDiodeCollection& collection, SCT_FrontEndData& data) const {
+ // **********************************************************************************
+ // Flag strips below threshold and flag the threshold check into data.m_StripHitsOnWafer
+ // **********************************************************************************
+
+ SiChargedDiodeIterator i_chargedDiode = collection.begin();
+ SiChargedDiodeIterator i_chargedDiode_end = collection.end();
+
+ for (; i_chargedDiode != i_chargedDiode_end; ++i_chargedDiode) {
+ SiChargedDiode& diode = (*i_chargedDiode).second;
+ SiReadoutCellId roCell = diode.getReadoutCell();
+ if (roCell.isValid()) {
+ int strip = roCell.strip();
+ if (strip > -1 and strip < m_strip_max) {
+ if (m_data_compression_mode == 1 and m_data_readout_mode == 0) { // level mode x1x
+ if (data.m_Analogue[1][strip] < m_Threshold) {
+ SiHelper::belowThreshold(diode, true); // Below strip diode signal threshold
+ data.m_StripHitsOnWafer[strip] = -1;
+ } else if (((0x10 & diode.flag()) == 0x10) or ((0x4 & diode.flag()) == 0x4)) {
+ // previously a crazy strip number could have screwed things up here.
+ data.m_StripHitsOnWafer[strip] = -1;
+ } else {
+ data.m_StripHitsOnWafer[strip] = 1;
+ SiHelper::SetTimeBin(diode, 2, &msg()); // set timebin info
+ }
+ } else if (m_data_compression_mode == 2 and m_data_readout_mode == 0) { // edge mode 01x
+ if ((data.m_Analogue[0][strip] >= m_Threshold or data.m_Analogue[1][strip] < m_Threshold)) {
+ SiHelper::belowThreshold(diode, true); // Below strip diode signal threshold
+ data.m_StripHitsOnWafer[strip] = -1;
+ } else if (((0x10 & diode.flag()) == 0x10) or ((0x4 & diode.flag()) == 0x4)) {
+ // previously a crazy strip number could have screwed things up here.
+ data.m_StripHitsOnWafer[strip] = -1;
+ } else {
+ data.m_StripHitsOnWafer[strip] = 1;
+ SiHelper::SetTimeBin(diode, 2, &msg()); // set timebin info
+ }
+ } else if (m_data_compression_mode == 3 or m_data_readout_mode == 1) { // Check hit pattern
+ int have_hit_bin = 0;
+ if (data.m_Analogue[0][strip] >= m_Threshold) {
+ have_hit_bin = 4;
+ }
+ if (data.m_Analogue[1][strip] >= m_Threshold) {
+ have_hit_bin += 2;
+ }
+ if (data.m_Analogue[2][strip] >= m_Threshold) {
+ have_hit_bin += 1;
+ }
+ if (((0x10 & diode.flag()) == 0x10) || ((0x4 & diode.flag()) == 0x4)) {
+ // previously a crazy strip number could have screwed things up here.
+ data.m_StripHitsOnWafer[strip] = -1;
+ } else if (m_data_compression_mode == 1) { // !< level and expanded mode
+ if (have_hit_bin == 2 or have_hit_bin == 3 or have_hit_bin == 6 or have_hit_bin == 7) {
+ data.m_StripHitsOnWafer[strip] = 1;
+ SiHelper::SetTimeBin(diode, have_hit_bin, &msg());
+ } else {
+ SiHelper::belowThreshold(diode, true); // Below strip diode signal threshold
+ data.m_StripHitsOnWafer[strip] = -1;
+ }
+ } else if (m_data_compression_mode == 2) { // !< edge and expanded mode
+ if (have_hit_bin == 2 or have_hit_bin == 3) {
+ data.m_StripHitsOnWafer[strip] = 1;
+ SiHelper::SetTimeBin(diode, have_hit_bin, &msg());
+ } else {
+ SiHelper::belowThreshold(diode, true); // Below strip diode signal threshold
+ data.m_StripHitsOnWafer[strip] = -1;
+ }
+ } else if (m_data_compression_mode == 3) { // !< any hit mode
+ if (have_hit_bin == 0) {
+ SiHelper::belowThreshold(diode, true); // Below strip diode signal threshold
+ data.m_StripHitsOnWafer[strip] = -1;
+ } else {
+ data.m_StripHitsOnWafer[strip] = 1;
+ if (m_data_readout_mode == 1) { // !< check for exp mode or not
+ SiHelper::SetTimeBin(diode, have_hit_bin, &msg());
+ } else {
+ SiHelper::SetTimeBin(diode, 2, &msg());
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+//
+// ----------------------------------------------------------------------
+StatusCode SCT_FrontEnd::doThresholdCheckForCrosstalkHits(SiChargedDiodeCollection& collection, SCT_FrontEndData& data) const {
+ // Check for noise+crosstalk strips above threshold
+ // data.m_StripHitsOnWafer: real hits above threshold == 1 or below/disconnected
+ // == -1
+ // =0 for free strips or strips with charge to be checked (data.m_Analogue[1]!=0)
+ // Set 2 for crosstalk noise hits and -2 for below ones
+
+ for (int strip = 0; strip < m_strip_max; strip++) {
+ // Find strips with data.m_StripHitsOnWafer[strip] == 0
+ if (data.m_StripHitsOnWafer[strip] != 0) { // real hits already checked
+ continue;
+ }
+ if (data.m_Analogue[1][strip] > 0) { // Better way of doing this?!
+ // set data.m_StripHitsOnWafer to x in prepareGainAndOffset
+ if (m_data_compression_mode == 1 and m_data_readout_mode == 0) { // level mode x1x
+ if (data.m_Analogue[1][strip] < m_Threshold) {
+ data.m_StripHitsOnWafer[strip] = -2; // Below threshold
+ } else {
+ data.m_StripHitsOnWafer[strip] = 2; // Crosstalk+Noise hit
+ if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, 2)) {
+ ATH_MSG_ERROR("Can't add noise hit diode to collection");
+ }
+ }
+ } else if (m_data_compression_mode == 2 and m_data_readout_mode == 0) { // edge mode 01x
+ if ((data.m_Analogue[0][strip] >= m_Threshold or data.m_Analogue[1][strip] < m_Threshold)) {
+ data.m_StripHitsOnWafer[strip] = -2; // Below threshold
+ } else {
+ data.m_StripHitsOnWafer[strip] = 2; // Crosstalk+Noise hit
+ if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, 2)) {
+ ATH_MSG_ERROR("Can't add noise hit diode to collection");
+ }
+ }
+ } else if (m_data_compression_mode == 3 or m_data_readout_mode == 1) {
+ int have_hit_bin = 0;
+ if (data.m_Analogue[0][strip] >= m_Threshold) {
+ have_hit_bin = 4;
+ }
+ if (data.m_Analogue[1][strip] >= m_Threshold) {
+ have_hit_bin += 2;
+ }
+ if (data.m_Analogue[2][strip] >= m_Threshold) {
+ have_hit_bin += 1;
+ }
+ if (m_data_compression_mode == 1) { // !< level and expanded mode
+ if (have_hit_bin == 2 or have_hit_bin == 3 or have_hit_bin == 6 or have_hit_bin == 7) {
+ data.m_StripHitsOnWafer[strip] = 2; // Crosstalk+Noise hit
+ if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, have_hit_bin)) {
+ ATH_MSG_ERROR("Can't add noise hit diode to collection");
+ }
+ } else {
+ data.m_StripHitsOnWafer[strip] = -2; // Below threshold
+ }
+ } else if (m_data_compression_mode == 2) { // !< edge and expanded mode
+ if (have_hit_bin == 2 or have_hit_bin == 3) {
+ data.m_StripHitsOnWafer[strip] = 2; // Noise hit
+ if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, have_hit_bin)) {
+ ATH_MSG_ERROR("Can't add noise hit diode to collection");
+ }
+ } else {
+ data.m_StripHitsOnWafer[strip] = -2; // Below threshold
+ }
+ } else if (m_data_compression_mode == 3) { // !< any hit mode
+ if (have_hit_bin == 0) {
+ data.m_StripHitsOnWafer[strip] = -2; // Below threshold
+ } else {
+ data.m_StripHitsOnWafer[strip] = 2; // !< Crosstalk+Noise hit
+ if (m_data_readout_mode == 1) { // !< check for exp mode or not
+ if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, have_hit_bin)) {
+ ATH_MSG_ERROR("Can't add noise hit diode to collection");
+ }
+ } else {
+ if (StatusCode::SUCCESS != addNoiseDiode(collection, strip, 2)) {
+ ATH_MSG_ERROR("Can't add noise hit diode to collection");
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+
+ return StatusCode::SUCCESS;
+}
+
+StatusCode SCT_FrontEnd::doClustering(SiChargedDiodeCollection& collection, SCT_FrontEndData& data) const {
+ // ********************************
+ // now do clustering
+ // ********************************
+ int strip = 0;
+ int clusterSize = 0;
+
+ const SCT_ModuleSideDesign& sctDesign = dynamic_cast<const SCT_ModuleSideDesign&>(collection.design());
+
+ Identifier hitStrip;
+
+ if (m_data_readout_mode == 0) {
+ do {
+ if (data.m_StripHitsOnWafer[strip] > 0) {
+ // ====== First step: Get the cluster size
+ // ===================================================
+ int clusterFirstStrip = strip;
+
+ // Find end of cluster. In multi-row sensors, cluster cannot span rows.
+ int row = sctDesign.row(strip);
+ if (row < 0) {
+ row = 0;
+ }
+
+ int lastStrip1DInRow = 0;
+ for (int i = 0; i < row + 1; ++i) {
+ lastStrip1DInRow += sctDesign.diodesInRow(i);
+ }
+
+ while (strip < lastStrip1DInRow-1 and data.m_StripHitsOnWafer[strip +1] > 0) {
+ ++strip; // !< Find first strip hit and add up the following strips
+ }
+ int clusterLastStrip = strip;
+
+ clusterSize = (clusterLastStrip - clusterFirstStrip) + 1;
+ hitStrip = m_sct_id->strip_id(collection.identify(), clusterFirstStrip);
+ SiChargedDiode& HitDiode = *(collection.find(hitStrip));
+ SiHelper::SetStripNum(HitDiode, clusterSize, &msg());
+
+ SiChargedDiode *PreviousHitDiode = &HitDiode;
+ for (int i = clusterFirstStrip+1; i <= clusterLastStrip; ++i) {
+ hitStrip = m_sct_id->strip_id(collection.identify(), i);
+ SiChargedDiode& HitDiode2 = *(collection.find(hitStrip));
+ SiHelper::ClusterUsed(HitDiode2, true);
+ PreviousHitDiode->setNextInCluster(&HitDiode2);
+ PreviousHitDiode = &HitDiode2;
+ }
+ }
+ ++strip; // !< This is the starting point of the next cluster within this collection
+ } while (strip < m_strip_max);
+ } else {
+ // Expanded read out mode, basically one RDO/strip per cluster
+ // But if consecutively fired strips have the same time bin, those are converted into one cluster.
+ do {
+ clusterSize = 1;
+ if (data.m_StripHitsOnWafer[strip] > 0) {
+ hitStrip = m_sct_id->strip_id(collection.identify(), strip);
+ SiChargedDiode& hitDiode = *(collection.find(hitStrip));
+ int timeBin = SiHelper::GetTimeBin(hitDiode);
+ SiChargedDiode* previousHitDiode = &hitDiode;
+ // Check if consecutively fired strips have the same time bin
+ for (int newStrip=strip+1; newStrip<m_strip_max; newStrip++) {
+ if (not (data.m_StripHitsOnWafer[newStrip]>0)) break;
+ Identifier newHitStrip = m_sct_id->strip_id(collection.identify(), newStrip);
+ SiChargedDiode& newHitDiode = *(collection.find(newHitStrip));
+ if (timeBin!=SiHelper::GetTimeBin(newHitDiode)) break;
+ SiHelper::ClusterUsed(newHitDiode, true);
+ previousHitDiode->setNextInCluster(&newHitDiode);
+ previousHitDiode = &newHitDiode;
+ clusterSize++;
+ }
+ SiHelper::SetStripNum(hitDiode, clusterSize, &msg());
+
+#ifdef SCT_DIG_DEBUG
+ ATH_MSG_DEBUG("RDO Strip = " << strip << ", tbin = " <<
+ SiHelper::GetTimeBin(hitDiode) <<
+ ", HitInfo(1=real, 2=crosstalk, 3=noise): " <<
+ data.m_StripHitsOnWafer[strip]);
+#endif
+ }
+ strip += clusterSize; // If more than one strip fires, those following strips are skipped.
+ } while (strip < m_strip_max);
+ }
+
+ // clusters below threshold, only from pre-digits that existed before no
+ // pure noise clusters below threshold will be made
+ // D. Costanzo. I don't see why we should cluster the strips below
+ // threshold. I'll pass on the info of strips below threshold
+ // to the SDOs. I'll leave the SCT experts the joy to change this if they
+ // don't like what I did or if this requires too much memory/disk
+
+ return StatusCode::SUCCESS;
+}
+
+StatusCode SCT_FrontEnd::addNoiseDiode(SiChargedDiodeCollection& collection, int strip, int tbin) const {
+ const SiCellId ndiode(strip); // !< create a new diode
+ const SiCharge noiseCharge(2 * m_Threshold, 0, SiCharge::noise); // !< add some noise to it
+ collection.add(ndiode, noiseCharge); // !< add it to the collection
+
+ // Get the strip back to check
+ const Identifier idstrip = m_sct_id->strip_id(collection.identify(), strip);
+ SiChargedDiode *NoiseDiode = (collection.find(idstrip));
+ if (NoiseDiode == nullptr) {
+ return StatusCode::FAILURE;
+ }
+ SiHelper::SetTimeBin(*NoiseDiode, tbin, &msg()); // set timebin info
+
+ return StatusCode::SUCCESS;
+}
+
+float SCT_FrontEnd::meanValue(std::vector<float>& calibDataVect) const {
+ float mean_value = 0.0;
+ int nData = 0;
+ const unsigned int vec_size = calibDataVect.size();
+
+ for (unsigned int i = 0; i < vec_size; ++i) {
+ if (calibDataVect[i] > 0.1) {
+ mean_value += calibDataVect[i];
+ ++nData;
+ }
+ }
+
+ if (nData == 0) {
+ return -1;
+ } else {
+ return mean_value / nData;
+ }
+}
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_FrontEnd.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_FrontEnd.h
new file mode 100644
index 0000000000000000000000000000000000000000..03f34b2c9c1681ff32f919bde33c80f6ef7a1b0f
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_FrontEnd.h
@@ -0,0 +1,140 @@
+// -*- C++ -*-
+
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+/**
+ * Header file for SCT_FrontEnd Class
+ * @brief simulation of the SCT front-end electronics
+ * working as a SiPreDigitsProcessor
+ * models response of ABCD chip amplifiers to
+ * collected charges, also does cross-talk, offset
+ * variation and gain variation, in a correlated way
+ * Version 1.0 04.05.2001 Szymon Gadomski
+ * 1.1 07.06.2001 pass amplifier and design as pointers
+ * 1.2, 15/06/2001, compiles with framework and tested
+ * 31.07.2001 changes, constructor now without design
+ * design is only known to process method (from the collection)
+ * 10.08.2001 dedided not to use AlwaysSame flag
+ * Implementation for using it may be provided later
+ * @author Szymon.Gadomski@cern.ch, Awatif.Belymam@cern.ch, Davide.Costanzo@cern.ch,
+ * tgcornel@nikhef.nl, Grant.Gorfine@cern.ch, Paul.Bell@cern.ch, Jorgen.Dalmau@cern.ch,
+ * 23/08/2007 - Kondo.Gnanvo@cern.ch
+ * - Conversion of the SCT_Front code Al gTool
+ */
+
+#ifndef FASERSCT_DIGITIZATION_SCT_FRONTEND_H
+#define FASERSCT_DIGITIZATION_SCT_FRONTEND_H
+
+// Inheritance
+#include "AthenaBaseComps/AthAlgTool.h"
+#include "FaserSCT_Digitization/ISCT_FrontEnd.h"
+
+// Athena
+#include "FaserSiDigitization/SiChargedDiodeCollection.h"
+// #include "SCT_ConditionsTools/ISCT_ReadCalibChipDataTool.h"
+
+// Gaudi
+#include "GaudiKernel/ToolHandle.h"
+
+// STL
+#include <atomic>
+#include <mutex>
+#include <vector>
+
+class ISCT_Amp;
+class FaserSCT_ID;
+
+namespace TrackerDD {
+ class SCT_DetectorManager;
+}
+
+namespace CLHEP {
+ class HepRandomEngine;
+}
+/**
+ * @brief simulation of the SCT front-end electronics
+ * working as a SiPreDigitsProcessor
+ * models response of ABCD chip amplifiers to
+ * collected charges, also does cross-talk, offset
+ * variation and gain variation, in a correlated way
+ */
+
+struct SCT_FrontEndData {
+ std::vector<float> m_Offset; //!< generate offset per channel
+ std::vector<float> m_GainFactor; //!< generate gain per channel (added to the gain per chip from calib data)
+ std::vector<float> m_NoiseFactor; //!< Kondo: 31/08/07 noise per channel (actually noise per chip from calib data)
+ std::vector<double> m_Analogue[3]; //!< To hold the noise and amplifier response
+ std::vector<int> m_StripHitsOnWafer; //!< Info about which strips are above threshold
+};
+
+class SCT_FrontEnd : public extends<AthAlgTool, ISCT_FrontEnd> {
+
+ public:
+
+ /** constructor */
+ SCT_FrontEnd(const std::string& type, const std::string& name, const IInterface* parent);
+ /** Destructor */
+ virtual ~SCT_FrontEnd() = default;
+ //PJ not needed after merging?!
+ /** AlgTool InterfaceID */
+ //static const InterfaceID& interfaceID();
+
+ /** AlgTool initialize */
+ virtual StatusCode initialize() override;
+ /** AlgTool finalize */
+ virtual StatusCode finalize() override;
+
+ /**
+ * process the collection of pre digits: needed to go through all single-strip pre-digits to calculate
+ * the amplifier response add noise (this could be moved elsewhere later) apply threshold do clustering
+ */
+ virtual void process(SiChargedDiodeCollection& collection, CLHEP::HepRandomEngine* rndmEngine) const override;
+ StatusCode doSignalChargeForHits(SiChargedDiodeCollection& collectione, SCT_FrontEndData& data) const;
+ StatusCode doThresholdCheckForRealHits(SiChargedDiodeCollection& collectione, SCT_FrontEndData& data) const;
+ StatusCode doThresholdCheckForCrosstalkHits(SiChargedDiodeCollection& collection, SCT_FrontEndData& data) const;
+ StatusCode doClustering(SiChargedDiodeCollection& collection, SCT_FrontEndData& data) const;
+ StatusCode prepareGainAndOffset(SiChargedDiodeCollection& collection, /*const Identifier& moduleId,*/ CLHEP::HepRandomEngine* rndmEngine, SCT_FrontEndData& data) const;
+// StatusCode prepareGainAndOffset(SiChargedDiodeCollection& collection, int side, const Identifier& moduleId, CLHEP::HepRandomEngine* rndmEngine, SCT_FrontEndData& data) const;
+ StatusCode randomNoise(SiChargedDiodeCollection& collection, /*const Identifier& moduleId,*/ CLHEP::HepRandomEngine* rndmEngine, SCT_FrontEndData& data) const;
+// StatusCode randomNoise(SiChargedDiodeCollection& collection, const Identifier& moduleId, int side, CLHEP::HepRandomEngine* rndmEngine, SCT_FrontEndData& data) const;
+ StatusCode addNoiseDiode(SiChargedDiodeCollection& collection, int strip, int tbin) const;
+ float meanValue(std::vector<float>& calibDataVect) const;
+ StatusCode initVectors(int strips, SCT_FrontEndData& data) const;
+
+ private:
+
+ FloatProperty m_NoiseBarrel{this, "NoiseBarrel", 1500.0, "Noise factor, Barrel (in the case of no use of calibration data)"};
+ FloatProperty m_NoiseBarrel3{this, "NoiseBarrel3", 1541.0, "Noise factor, Barrel3 (in the case of no use of calibration data)"};
+ FloatProperty m_NoiseInners{this, "NoiseInners", 1090.0, "Noise factor, EC Inners (in the case of no use of calibration data)"};
+ FloatProperty m_NoiseMiddles{this, "NoiseMiddles", 1557.0, "Noise factor, EC Middles (in the case of no use of calibration data)"};
+ FloatProperty m_NoiseShortMiddles{this, "NoiseShortMiddles", 940.0, "Noise factor, EC Short Middles (in the case of no use of calibration data)"};
+ FloatProperty m_NoiseOuters{this, "NoiseOuters", 1618.0, "Noise factor, Ec Outers (in the case of no use of calibration data)"};
+ DoubleProperty m_NOBarrel{this, "NOBarrel", 1.5e-5, "Noise factor, Barrel (in the case of no use of calibration data)"};
+ DoubleProperty m_NOBarrel3{this, "NOBarrel3", 2.1e-5, "Noise factor, Barrel3 (in the case of no use of calibration data)"};
+ DoubleProperty m_NOInners{this, "NOInners", 5.0e-9, "Noise Occupancy, EC Inners (in the case of no use of calibration data)"};
+ DoubleProperty m_NOMiddles{this, "NOMiddles", 2.7e-5, "Noise Occupancy, EC Middles (in the case of no use of calibration data)"};
+ DoubleProperty m_NOShortMiddles{this, "NOShortMiddles", 2.0e-9, "Noise Occupancy, EC Short Middles (in the case of no use of calibration data)"};
+ DoubleProperty m_NOOuters{this, "NOOuters", 3.5e-5, "Noise Occupancy, Ec Outers (in the case of no use of calibration data)"};
+ BooleanProperty m_NoiseOn{this, "NoiseOn", true, "To know if noise is on or off when using calibration data"};
+ BooleanProperty m_analogueNoiseOn{this, "AnalogueNoiseOn", true, "To know if analogue noise is on or off"};
+ FloatProperty m_GainRMS{this, "GainRMS", 0.031, "Gain spread parameter within the strips for a given Chip gain"};
+ FloatProperty m_Ospread{this, "Ospread", 0.0001, "offset spread within the strips for a given Chip offset"};
+ FloatProperty m_OGcorr{this, "OffsetGainCorrelation", 0.00001, "Gain/offset correlation for the strips"};
+ FloatProperty m_Threshold{this, "Threshold", 1.0, "Threshold"};
+ FloatProperty m_timeOfThreshold{this, "TimeOfThreshold", 30.0, "Threshold time"};
+ ShortProperty m_data_compression_mode{this, "DataCompressionMode", 1, "Front End Data Compression Mode"};
+ ShortProperty m_data_readout_mode{this, "DataReadOutMode", 0, "Front End Data Read out mode Mode"};
+ BooleanProperty m_useCalibData{this, "UseCalibData", false, "Flag to set the use of calibration data for noise, Gain,offset etc."}; // was true in ATLAS
+
+ ToolHandle<ISCT_Amp> m_sct_amplifier{this, "SCT_Amp", "SCT_Amp", "Handle the Amplifier tool"}; //!< Handle the Amplifier tool
+// ToolHandle<ISCT_ReadCalibChipDataTool> m_ReadCalibChipDataTool{this, "SCT_ReadCalibChipDataTool", "SCT_ReadCalibChipDataTool", "Tool to retrieve chip calibration information"}; //!< Handle to the Calibration ConditionsTool
+
+ const TrackerDD::SCT_DetectorManager* m_SCTdetMgr{nullptr}; //!< Handle to SCT detector manager
+ const FaserSCT_ID* m_sct_id{nullptr}; //!< Handle to SCT ID helper
+
+ mutable std::atomic_int m_strip_max{768}; //!< For SLHC studies
+};
+
+#endif //SCT_FRONTEND_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_RandomDisabledCellGenerator.cxx b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_RandomDisabledCellGenerator.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..b170bd10f5c166e02e8258d5138f0fdd7fc23cac
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_RandomDisabledCellGenerator.cxx
@@ -0,0 +1,40 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "SCT_RandomDisabledCellGenerator.h"
+
+#include "FaserSiDigitization/SiChargedDiodeCollection.h"
+
+#include "FaserSiDigitization/SiHelper.h"
+
+#include "CLHEP/Random/RandomEngine.h"
+#include "CLHEP/Random/RandFlat.h"
+
+// constructor
+SCT_RandomDisabledCellGenerator::SCT_RandomDisabledCellGenerator(const std::string& type, const std::string& name, const IInterface* parent)
+ : base_class(type, name, parent)
+{
+}
+
+StatusCode SCT_RandomDisabledCellGenerator::initialize() {
+ ATH_MSG_DEBUG("SCT_RandomDisabledCellGenerator::initialize()");
+ ATH_MSG_INFO("\tCreating missing bond generator with "<<m_disableProbability<<" probability");
+ return StatusCode::SUCCESS;
+}
+
+StatusCode SCT_RandomDisabledCellGenerator::finalize() {
+ ATH_MSG_INFO("SCT_RandomDisabledCellGenerator::finalize()");
+ return StatusCode::SUCCESS;
+}
+
+// process the collection
+void SCT_RandomDisabledCellGenerator::process(SiChargedDiodeCollection& collection, CLHEP::HepRandomEngine * rndmEngine) const {
+ // disabling is applied to all cells even unconnected or below threshold ones to be able to use these cells as well
+ // loop on all charged diodes
+ for (std::pair<const TrackerDD::SiCellId, SiChargedDiode>& chargedDiode: collection) {
+ if (CLHEP::RandFlat::shoot(rndmEngine)<m_disableProbability) {
+ SiHelper::disconnected(chargedDiode.second, true, false);
+ }
+ }
+}
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_RandomDisabledCellGenerator.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_RandomDisabledCellGenerator.h
new file mode 100644
index 0000000000000000000000000000000000000000..1ea27123d73228fb6fd62fb6960ed8daa24472be
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_RandomDisabledCellGenerator.h
@@ -0,0 +1,67 @@
+// -*- C++ -*-
+
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SCT_RandomDisabledCellGenerator.h
+// Header file for class SCT_RandomDisabledCellGenerator
+///////////////////////////////////////////////////////////////////
+// Class to randomly disable cells for SCT
+///////////////////////////////////////////////////////////////////
+// Version 1.0 25/07/2007 Kondo Gnanvo
+// - Randomly select a cell and modify its flag by using the pointer to
+// a SiHelper function.
+// - Inherits from ISiChargedDiodeProcessor, as SiPreDiodeProcessor
+// has been discontinued
+// This can be used for disabling, disconnecting, and flagging bad_tot
+///////////////////////////////////////////////////////////////////
+
+#ifndef FASERSCT_DIGITIZATION_SCT_RANDOMDISABLEDCELLGENERATOR_H
+#define FASERSCT_DIGITIZATION_SCT_RANDOMDISABLEDCELLGENERATOR_H
+
+//Inheritance
+#include "AthenaBaseComps/AthAlgTool.h"
+#include "FaserSCT_Digitization/ISCT_RandomDisabledCellGenerator.h"
+
+class SiChargedDiode;
+class SiChargedDiodeCollection;
+
+namespace CLHEP {
+ class HepRandomEngine;
+}
+
+class SCT_RandomDisabledCellGenerator : public extends<AthAlgTool, ISCT_RandomDisabledCellGenerator> {
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+
+ public:
+
+ /** constructor */
+ SCT_RandomDisabledCellGenerator(const std::string& type, const std::string& name, const IInterface* parent) ;
+
+ /** Destructor */
+ virtual ~SCT_RandomDisabledCellGenerator() = default;
+
+ /** AlgTool initialize */
+ virtual StatusCode initialize() override;
+ /** AlgTool finalize */
+ virtual StatusCode finalize() override;
+
+ virtual void process(SiChargedDiodeCollection& collection, CLHEP::HepRandomEngine * rndmEngine) const override;
+
+ private:
+
+ ///////////////////////////////////////////////////////////////////
+ // Private data:
+ ///////////////////////////////////////////////////////////////////
+ private:
+
+ FloatProperty m_disableProbability{this, "TotalBadChannels", 0.01, "probability that a cell is disabled"};
+
+};
+
+#endif //FASERSCT_DIGITIZATION_SCT_RANDOMDISABLEDCELLGENERATOR_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_SurfaceChargesGenerator.cxx b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_SurfaceChargesGenerator.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..b05cd7b4468a07ce90afc0346eeb6f07ed77d442
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_SurfaceChargesGenerator.cxx
@@ -0,0 +1,578 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "SCT_SurfaceChargesGenerator.h"
+
+// DD
+#include "TrackerReadoutGeometry/SiDetectorElement.h"
+#include "TrackerReadoutGeometry/SCT_ModuleSideDesign.h"
+
+// Athena
+#include "GeneratorObjects/HepMcParticleLink.h"
+#include "TrackerSimEvent/FaserSiHit.h" // for SiHit, SiHit::::xDep, etc
+#include "HitManagement/TimedHitPtr.h" // for TimedHitPtr
+
+// ROOT
+#include "TH1.h" // for TH1F
+#include "TH2.h" // for TH2F
+#include "TProfile.h" // for TProfile
+
+// CLHEP
+#include "CLHEP/Geometry/Point3D.h"
+#include "CLHEP/Random/RandomEngine.h"
+#include "CLHEP/Random/RandGaussZiggurat.h"
+#include "CLHEP/Units/SystemOfUnits.h"
+
+// STL
+#include <cmath>
+
+using TrackerDD::SiDetectorElement;
+using TrackerDD::SCT_ModuleSideDesign;
+using TrackerDD::SiLocalPosition;
+
+using namespace std;
+
+// constructor
+SCT_SurfaceChargesGenerator::SCT_SurfaceChargesGenerator(const std::string& type,
+ const std::string& name,
+ const IInterface* parent)
+ : base_class(type, name, parent) {
+}
+
+// ----------------------------------------------------------------------
+// Initialize
+// ----------------------------------------------------------------------
+StatusCode SCT_SurfaceChargesGenerator::initialize() {
+ ATH_MSG_DEBUG("SCT_SurfaceChargesGenerator::initialize()");
+
+ // Get ISiPropertiesTool
+ ATH_CHECK(m_siPropertiesTool.retrieve());
+
+ // Get ISiliconConditionsSvc
+// ATH_CHECK(m_siConditionsTool.retrieve());
+
+ if (m_doTrapping) {
+ // -- Get Radiation Damage Tool
+ ATH_MSG_FATAL("Radiation damage tool not supported.");
+// ATH_CHECK(m_radDamageTool.retrieve());
+// } else {
+// m_radDamageTool.disable();
+ }
+
+// if (m_doHistoTrap) {
+// // -- Get Histogram Service
+// ATH_CHECK(m_thistSvc.retrieve());
+
+// m_h_efieldz = new TProfile("efieldz", "", 50, 0., 0.4);
+// ATH_CHECK(m_thistSvc->regHist("/file1/efieldz", m_h_efieldz));
+
+// m_h_efield = new TH1F("efield", "", 100, 200., 800.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/efield", m_h_efield));
+
+// m_h_spess = new TH1F("spess", "", 50, 0., 0.4);
+// ATH_CHECK(m_thistSvc->regHist("/file1/spess", m_h_spess));
+
+// m_h_depD = new TH1F("depD", "", 50, -0.3, 0.3);
+// ATH_CHECK(m_thistSvc->regHist("/file1/depD", m_h_depD));
+
+// m_h_drift_electrode = new TH2F("drift_electrode", "", 50, 0., 20., 50, 0., 20.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/drift_electrode", m_h_drift_electrode));
+
+// m_h_drift_time = new TH1F("drift_time", "", 100, 0., 20.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/drift_time", m_h_drift_time));
+
+// m_h_t_electrode = new TH1F("t_electrode", "", 100, 0., 20.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/t_electrode", m_h_t_electrode));
+
+// m_h_ztrap = new TH1F("ztrap", "", 100, 0., 0.3);
+// ATH_CHECK(m_thistSvc->regHist("/file1/ztrap", m_h_ztrap));
+
+// // More histograms to check what's going on
+// m_h_zhit = new TH1F("zhit", "", 50, -0.2, 0.2);
+// ATH_CHECK(m_thistSvc->regHist("/file1/zhit", m_h_zhit));
+
+// m_h_ztrap_tot = new TH1F("ztrap_tot", "", 100, 0., 0.5);
+// ATH_CHECK(m_thistSvc->regHist("/file1/ztrap_tot", m_h_ztrap_tot));
+
+// m_h_no_ztrap = new TH1F("no_ztrap", "", 100, 0., 0.5);
+// ATH_CHECK(m_thistSvc->regHist("/file1/no_ztrap", m_h_no_ztrap));
+
+// m_h_trap_drift_t = new TH1F("trap_drift_t", "", 100, 0., 20.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/trap_drift_t", m_h_trap_drift_t));
+
+// m_h_notrap_drift_t = new TH1F("notrap_drift_t", "", 100, 0., 20.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/notrap_drift_t", m_h_notrap_drift_t));
+
+// m_h_mob_Char = new TProfile("mob_Char", "", 200, 100., 1000.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/mob_Char", m_h_mob_Char));
+
+// m_h_vel = new TProfile("vel", "", 100, 100., 1000.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/vel", m_h_vel));
+
+// m_h_drift1 = new TProfile("drift1", "", 50, 0., 0.3);
+// ATH_CHECK(m_thistSvc->regHist("/file1/drift1", m_h_drift1));
+
+// m_h_gen = new TProfile("gen", "", 50, 0., 0.3);
+// ATH_CHECK(m_thistSvc->regHist("/file1/gen", m_h_gen));
+
+// m_h_gen1 = new TProfile("gen1", "", 50, 0., 0.3);
+// ATH_CHECK(m_thistSvc->regHist("/file1/gen1", m_h_gen1));
+
+// m_h_gen2 = new TProfile("gen2", "", 50, 0., 0.3);
+// ATH_CHECK(m_thistSvc->regHist("/file1/gen2", m_h_gen2));
+
+// m_h_velocity_trap = new TProfile("velocity_trap", "", 50, 0., 1000.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/velocity_trap", m_h_velocity_trap));
+
+// m_h_mobility_trap = new TProfile("mobility_trap", "", 50, 100., 1000.);
+// ATH_CHECK(m_thistSvc->regHist("/file1/mobility_trap", m_h_mobility_trap));
+
+// m_h_trap_pos = new TH1F("trap_pos", "", 100, 0., 0.3);
+// ATH_CHECK(m_thistSvc->regHist("/file1/trap_pos", m_h_trap_pos));
+// }
+ ///////////////////////////////////////////////////
+
+ m_smallStepLength.setValue(m_smallStepLength.value() * CLHEP::micrometer);
+ m_tSurfaceDrift.setValue(m_tSurfaceDrift.value() * CLHEP::ns);
+
+ // Surface drift time calculation Stuff
+ m_tHalfwayDrift = m_tSurfaceDrift * 0.5;
+ m_distHalfInterStrip = m_distInterStrip * 0.5;
+ if ((m_tSurfaceDrift > m_tHalfwayDrift) and (m_tHalfwayDrift >= 0.0) and
+ (m_distHalfInterStrip > 0.0) and (m_distHalfInterStrip < m_distInterStrip)) {
+ m_SurfaceDriftFlag = true;
+ } else {
+ ATH_MSG_INFO("\tsurface drift still not on, wrong params");
+ }
+
+ // Make sure these two flags are not set simultaneously
+ if (m_tfix>-998. and m_tsubtract>-998.) {
+ ATH_MSG_FATAL("\tCannot set both FixedTime and SubtractTime options!");
+ ATH_MSG_INFO("\tMake sure the two flags are not set simultaneously in jo");
+ return StatusCode::FAILURE;
+ }
+
+ if (m_doDistortions) {
+ ATH_MSG_FATAL("Distortions not supported.");
+// ATH_CHECK(m_distortionsTool.retrieve());
+// } else {
+// m_distortionsTool.disable();
+ }
+
+ ATH_CHECK(m_lorentzAngleTool.retrieve());
+
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+// finalize
+// ----------------------------------------------------------------------
+StatusCode SCT_SurfaceChargesGenerator::finalize() {
+ ATH_MSG_DEBUG("SCT_SurfaceChargesGenerator::finalize()");
+ return StatusCode::SUCCESS;
+}
+
+// ----------------------------------------------------------------------
+// perpandicular Drift time calculation
+// ----------------------------------------------------------------------
+float SCT_SurfaceChargesGenerator::driftTime(float zhit, const SiDetectorElement* element) const {
+ if (element==nullptr) {
+ ATH_MSG_ERROR("SCT_SurfaceChargesGenerator::process element is nullptr");
+ return -2.0;
+ }
+ const SCT_ModuleSideDesign* design{dynamic_cast<const SCT_ModuleSideDesign*>(&(element->design()))};
+ if (design==nullptr) {
+ ATH_MSG_ERROR("SCT_SurfaceChargesGenerator::process can not get " << design);
+ return -2.0;
+ }
+ const double thickness{design->thickness()};
+ const IdentifierHash hashId{element->identifyHash()};
+
+ if ((zhit < 0.0) or (zhit > thickness)) {
+ ATH_MSG_DEBUG("driftTime: hit coordinate zhit=" << zhit / CLHEP::micrometer << " out of range");
+ return -2.0;
+ }
+
+ float depletionVoltage{0.};
+ float biasVoltage{0.};
+ if (m_useSiCondDB) {
+ ATH_MSG_FATAL("Conditions DB voltages not supported.");
+ // depletionVoltage = m_siConditionsTool->depletionVoltage(hashId) * CLHEP::volt;
+ // biasVoltage = m_siConditionsTool->biasVoltage(hashId) * CLHEP::volt;
+ } else {
+ depletionVoltage = m_vdepl * CLHEP::volt;
+ biasVoltage = m_vbias * CLHEP::volt;
+ }
+
+ const float denominator{static_cast<float>(depletionVoltage + biasVoltage - (2.0 * zhit * depletionVoltage / thickness))};
+ if (denominator <= 0.0) {
+ if (biasVoltage >= depletionVoltage) { // Should not happen
+ if(not m_isOverlay) {
+ ATH_MSG_ERROR("driftTime: negative argument X for log(X) " << zhit);
+ }
+ return -1.0;
+ } else {
+ // (m_biasVoltage<m_depletionVoltage) can happen with underdepleted sensors, lose charges in that volume
+ return -10.0;
+ }
+ }
+
+ float t_drift{log((depletionVoltage + biasVoltage) / denominator)};
+ t_drift *= thickness * thickness / (2.0 * m_siPropertiesTool->getSiProperties(hashId).holeDriftMobility() * depletionVoltage);
+ return t_drift;
+}
+
+// ----------------------------------------------------------------------
+// Sigma diffusion calculation
+// ----------------------------------------------------------------------
+float SCT_SurfaceChargesGenerator::diffusionSigma(float zhit, const SiDetectorElement* element) const {
+ if (element==nullptr) {
+ ATH_MSG_ERROR("SCT_SurfaceChargesGenerator::diffusionSigma element is nullptr");
+ return 0.0;
+ }
+ const IdentifierHash hashId{element->identifyHash()};
+ const float t{driftTime(zhit, element)}; // in ns
+
+ if (t > 0.0) {
+ const float sigma{static_cast<float>(sqrt(2. * m_siPropertiesTool->getSiProperties(hashId).holeDiffusionConstant() * t))}; // in mm
+ return sigma;
+ } else {
+ return 0.0;
+ }
+}
+
+// ----------------------------------------------------------------------
+// Maximum drift time
+// ----------------------------------------------------------------------
+float SCT_SurfaceChargesGenerator::maxDriftTime(const SiDetectorElement* element) const {
+ if (element) {
+ const float sensorThickness{static_cast<float>(element->thickness())};
+ return driftTime(sensorThickness, element);
+ } else {
+ ATH_MSG_INFO("Error: SiDetectorElement not set!");
+ return 0.;
+ }
+}
+
+// ----------------------------------------------------------------------
+// Maximum Sigma difusion
+// ----------------------------------------------------------------------
+float SCT_SurfaceChargesGenerator::maxDiffusionSigma(const SiDetectorElement* element) const {
+ if (element) {
+ const float sensorThickness{static_cast<float>(element->thickness())};
+ return diffusionSigma(sensorThickness, element);
+ } else {
+ ATH_MSG_INFO("Error: SiDetectorElement not set!");
+ return 0.;
+ }
+}
+
+// ----------------------------------------------------------------------
+// Calculating the surface drift time but I should confess that
+// I haven't found out yet where the calculation come from
+// ----------------------------------------------------------------------
+float SCT_SurfaceChargesGenerator::surfaceDriftTime(float ysurf) const {
+ if (m_SurfaceDriftFlag) {
+ if ((ysurf >= 0.0) and (ysurf <= m_distInterStrip)) {
+ float t_surfaceDrift{0.};
+ if (ysurf < m_distHalfInterStrip) {
+ const float y{ysurf / m_distHalfInterStrip};
+ t_surfaceDrift= m_tHalfwayDrift * y * y;
+ } else {
+ const float y{(m_distInterStrip - ysurf) / (m_distInterStrip - m_distHalfInterStrip)};
+ t_surfaceDrift = m_tSurfaceDrift + (m_tHalfwayDrift - m_tSurfaceDrift) * y * y;
+ }
+ return t_surfaceDrift;
+ } else {
+ ATH_MSG_INFO(" ysurf out of range " << ysurf);
+ return -1.0;
+ }
+ } else {
+ return 0.0;
+ }
+}
+
+// -------------------------------------------------------------------------------------------
+// create a list of surface charges from a hit - called from SCT_Digitization
+// AthAlgorithm
+// -------------------------------------------------------------------------------------------
+void SCT_SurfaceChargesGenerator::process(const SiDetectorElement* element,
+ const TimedHitPtr<FaserSiHit>& phit,
+ const ISiSurfaceChargesInserter& inserter, CLHEP::HepRandomEngine * rndmEngine) const {
+ ATH_MSG_VERBOSE("SCT_SurfaceChargesGenerator::process starts");
+ processSiHit(element, *phit, inserter, phit.eventTime(), phit.eventId(), rndmEngine);
+ return;
+}
+
+// -------------------------------------------------------------------------------------------
+// create a list of surface charges from a hit - called from both AthAlgorithm
+// and PileUpTool
+// -------------------------------------------------------------------------------------------
+void SCT_SurfaceChargesGenerator::processSiHit(const SiDetectorElement* element,
+ const FaserSiHit& phit,
+ const ISiSurfaceChargesInserter& inserter,
+ float p_eventTime,
+ unsigned short p_eventId, CLHEP::HepRandomEngine * rndmEngine) const {
+ const SCT_ModuleSideDesign* design{dynamic_cast<const SCT_ModuleSideDesign*>(&(element->design()))};
+ if (design==nullptr) {
+ ATH_MSG_ERROR("SCT_SurfaceChargesGenerator::process can not get " << design);
+ return;
+ }
+
+ const double thickness{design->thickness()};
+ const IdentifierHash hashId{element->identifyHash()};
+ const double tanLorentz{m_lorentzAngleTool->getTanLorentzAngle(hashId)};
+
+ // ---**************************************
+ // Time of Flight Calculation - separate method?
+ // ---**************************************
+ // --- Original calculation of Time of Flight of the particle Time needed by the particle to reach the sensor
+ float timeOfFlight{p_eventTime + hitTime(phit)};
+
+ // Kondo 19/09/2007: Use the coordinate of the center of the module to calculate the time of flight
+ timeOfFlight -= (element->center().mag()) / CLHEP::c_light;
+ // !< extract the distance to the origin of the module to Time of flight
+
+ // !< timing set from jo to adjust (subtract) the timing
+ if (m_tsubtract > -998.) {
+ timeOfFlight -= m_tsubtract;
+ }
+ // ---**************************************
+
+ const CLHEP::Hep3Vector pos{phit.localStartPosition()};
+ const float xEta{static_cast<float>(pos[FaserSiHit::xEta])};
+ const float xPhi{static_cast<float>(pos[FaserSiHit::xPhi])};
+ const float xDep{static_cast<float>(pos[FaserSiHit::xDep])};
+
+ const CLHEP::Hep3Vector endPos{phit.localEndPosition()};
+ const float cEta{static_cast<float>(endPos[FaserSiHit::xEta]) - xEta};
+ const float cPhi{static_cast<float>(endPos[FaserSiHit::xPhi]) - xPhi};
+ const float cDep{static_cast<float>(endPos[FaserSiHit::xDep]) - xDep};
+
+ const float LargeStep{sqrt(cEta*cEta + cPhi*cPhi + cDep*cDep)};
+ const int numberOfSteps{static_cast<int>(LargeStep / m_smallStepLength) + 1};
+ const float steps{static_cast<float>(m_numberOfCharges * numberOfSteps)};
+ const float e1{static_cast<float>(phit.energyLoss() / steps)};
+ const float q1{static_cast<float>(e1 * m_siPropertiesTool->getSiProperties(hashId).electronHolePairsPerEnergy())};
+
+ // in the following, to test the code, we will use the original coordinate
+ // system of the SCTtest3SurfaceChargesGenerator x is eta y is phi z is depth
+ float xhit{xEta};
+ float yhit{xPhi};
+ float zhit{xDep};
+
+ if (m_doDistortions) {
+ ATH_MSG_FATAL("Distortions not supported.");
+ // if (element->isBarrel()) {// Only apply disortions to barrel modules
+ // Amg::Vector2D BOW;
+ // BOW[0] = m_distortionsTool->correctSimulation(hashId, xhit, yhit, cEta, cPhi, cDep)[0];
+ // BOW[1] = m_distortionsTool->correctSimulation(hashId, xhit, yhit, cEta, cPhi, cDep)[1];
+ // xhit = BOW.x();
+ // yhit = BOW.y();
+ // }
+ }
+
+ const float StepX{cEta / numberOfSteps};
+ const float StepY{cPhi / numberOfSteps};
+ const float StepZ{cDep / numberOfSteps};
+
+ // check the status of truth information for this SiHit
+ // some Truth information is cut for pile up events
+ const HepMcParticleLink trklink{HepMcParticleLink(phit.trackNumber(), p_eventId)};
+ SiCharge::Process hitproc{SiCharge::track};
+ if (phit.trackNumber() != 0) {
+ if (not trklink.isValid()) {
+ hitproc = SiCharge::cut_track;
+ }
+ }
+
+ float dstep{-0.5};
+ for (int istep{0}; istep < numberOfSteps; ++istep) {
+ dstep += 1.0;
+ float z1{zhit + StepZ * dstep};
+
+ // Distance between charge and readout side.
+ // design->readoutSide() is +1 if readout side is in +ve depth axis direction and visa-versa.
+ float zReadout{static_cast<float>(0.5 * thickness - design->readoutSide() * z1)};
+ // const double spess{zReadout};
+
+ // if (m_doHistoTrap) {
+ // m_h_depD->Fill(z1);
+ // m_h_spess->Fill(spess);
+ // }
+
+ float t_drift{driftTime(zReadout, element)}; // !< t_drift: perpandicular drift time
+ if (t_drift>-2.0000002 and t_drift<-1.9999998) {
+ ATH_MSG_DEBUG("Checking for rounding errors in compression");
+ if ((fabs(z1) - 0.5 * thickness) < 0.000010) {
+ ATH_MSG_DEBUG("Rounding error found attempting to correct it. z1 = " << std::fixed << std::setprecision(8) << z1);
+ if (z1 < 0.0) {
+ z1 = 0.0000005 - 0.5 * thickness;
+ // set new coordinate to be 0.5nm inside wafer volume.
+ } else {
+ z1 = 0.5 * thickness - 0.0000005;
+ // set new coordinate to be 0.5nm inside wafer volume.
+ }
+ zReadout = 0.5 * thickness - design->readoutSide() * z1;
+ t_drift = driftTime(zReadout, element);
+ if (t_drift>-2.0000002 and t_drift<-1.9999998) {
+ ATH_MSG_WARNING("Attempt failed. Making no correction.");
+ } else {
+ ATH_MSG_DEBUG("Correction Successful! z1 = " << std::fixed << std::setprecision(8) << z1 << ", zReadout = " << zReadout << ", t_drift = " << t_drift);
+ }
+ } else {
+ ATH_MSG_DEBUG("No rounding error found. Making no correction.");
+ }
+ }
+ if (t_drift > 0.0) {
+ const float x1{xhit + StepX * dstep};
+ float y1{yhit + StepY * dstep};
+ y1 += tanLorentz * zReadout; // !< Taking into account the magnetic field
+ const float sigma{diffusionSigma(zReadout, element)};
+
+ for (int i{0}; i < m_numberOfCharges; ++i) {
+ const float rx{CLHEP::RandGaussZiggurat::shoot(rndmEngine)};
+ const float xd{x1 + sigma * rx};
+ const float ry{CLHEP::RandGaussZiggurat::shoot(rndmEngine)};
+ const float yd{y1 + sigma * ry};
+
+ // For charge trapping with Ramo potential
+ const double stripPitch{0.080}; // mm
+ double dstrip{y1 / stripPitch}; // mm
+ if (dstrip > 0.) {
+ dstrip -= static_cast<double>(static_cast<int>(dstrip));
+ } else {
+ dstrip -= static_cast<double>(static_cast<int>(dstrip)) + 1;
+ }
+
+ // now y will be x and z will be y ....just to make sure to confuse everebody
+ // double y0{dstrip * stripPitch}; // mm
+ // double z0{thickness - zReadout}; // mm
+
+ // -- Charge Trapping
+ // if (m_doTrapping) {
+ // if (m_doHistoTrap) {
+ // m_h_zhit->Fill(zhit);
+ // }
+ // double trap_pos{-999999.}, drift_time{-999999.}; // FIXME need better default values
+ // if (chargeIsTrapped(spess, element, trap_pos, drift_time)) {
+ // if (not m_doRamo) {
+ // break;
+ // } else { // if we want to take into account also Ramo Potential
+ // double Q_m2{0.}, Q_m1{0.}, Q_00{0.}, Q_p1{0.}, Q_p2{0.}; // Charges
+
+ // dstrip = y1 / stripPitch; // mm
+ // // need the distance from the nearest strips
+ // // edge not centre, xtaka = 1/2*stripPitch
+ // // centre of detector, y1=0, is in the middle of
+ // // an interstrip gap
+ // if (dstrip > 0.) {
+ // dstrip -= static_cast<double>(static_cast<int>(dstrip));
+ // } else {
+ // dstrip -= static_cast<double>(static_cast<int>(dstrip)) + 1;
+ // }
+
+ // // now y will be x and z will be y ....just to make sure to confuse everebody
+ // double yfin{dstrip * stripPitch}; // mm
+ // double zfin{thickness - trap_pos}; // mm
+
+ // m_radDamageTool->holeTransport(y0, z0, yfin, zfin, Q_m2, Q_m1, Q_00, Q_p1, Q_p2);
+ // for (int strip{-2}; strip<=2; strip++) {
+ // const double ystrip{yd + strip * stripPitch}; // mm
+ // const SiLocalPosition position(element->hitLocalToLocal(xd, ystrip));
+ // if (design->inActiveArea(position)) {
+ // double charge{0.};
+ // if (strip == -2) charge = Q_m2;
+ // if (strip == -1) charge = Q_m1;
+ // if (strip == 0) charge = Q_00;
+ // if (strip == 1) charge = Q_p1;
+ // if (strip == 2) charge = Q_p2;
+ // const double time{drift_time};
+ // if (charge != 0.) {
+ // inserter(SiSurfaceCharge(position, SiCharge(q1*charge, time, hitproc, trklink)));
+ // continue;
+ // }
+ // }
+ // }
+ // ATH_MSG_INFO("strip zero charge = " << Q_00); // debug
+ // } // m_doRamo==true
+ // } // chargeIsTrapped()
+ // } // m_doTrapping==true
+
+ if (not m_doRamo) {
+ const SiLocalPosition position{element->hitLocalToLocal(xd, yd)};
+ if (design->inActiveArea(position)) {
+ const float sdist{static_cast<float>(design->scaledDistanceToNearestDiode(position))}; // !< dist on the surface from the hit point to the nearest strip (diode)
+ const float t_surf{surfaceDriftTime(2.0 * sdist)}; // !< Surface drift time
+ const float totaltime{(m_tfix > -998.) ? m_tfix.value() : t_drift + timeOfFlight + t_surf}; // !< Total drift time
+ inserter(SiSurfaceCharge(position, SiCharge(q1, totaltime, hitproc, trklink)));
+ } else {
+ ATH_MSG_VERBOSE(std::fixed << std::setprecision(8) << "Local position (phi, eta, depth): ("
+ << position.xPhi() << ", " << position.xEta() << ", " << position.xDepth()
+ << ") of the element is out of active area, charge = " << q1);
+ }
+ } // end of loop on charges
+ }
+ }
+ }
+ return;
+}
+
+// bool SCT_SurfaceChargesGenerator::chargeIsTrapped(double spess,
+// const SiDetectorElement* element,
+// double& trap_pos,
+// double& drift_time) const {
+// if (element==nullptr) {
+// ATH_MSG_ERROR("SCT_SurfaceChargesGenerator::chargeIsTrapped element is nullptr");
+// return false;
+// }
+// bool isTrapped{false};
+// const IdentifierHash hashId{element->identifyHash()};
+// const SCT_ChargeTrappingCondData condData{m_radDamageTool->getCondData(hashId, spess)};
+// const double electric_field{condData.getElectricField()};
+
+// if (m_doHistoTrap) {
+// const double mobChar{condData.getHoleDriftMobility()};
+// m_h_efieldz->Fill(spess, electric_field);
+// m_h_efield->Fill(electric_field);
+// m_h_mob_Char->Fill(electric_field, mobChar);
+// m_h_vel->Fill(electric_field, electric_field * mobChar);
+// }
+// const double t_electrode{condData.getTimeToElectrode()};
+// drift_time = condData.getTrappingTime();
+// const double z_trap{condData.getTrappingPositionZ()};
+// trap_pos = spess - z_trap;
+// if (m_doHistoTrap) {
+// m_h_drift_time->Fill(drift_time);
+// m_h_t_electrode->Fill(t_electrode);
+// m_h_drift_electrode->Fill(drift_time, t_electrode);
+// m_h_ztrap_tot->Fill(z_trap);
+// }
+// // -- Calculate if the charge is trapped, and at which distance
+// // -- Charge gets trapped before arriving to the electrode
+// if (drift_time < t_electrode) {
+// isTrapped = true;
+// ATH_MSG_INFO("drift_time: " << drift_time << " t_electrode: " << t_electrode << " spess " << spess);
+// ATH_MSG_INFO("z_trap: " << z_trap);
+// if (m_doHistoTrap) {
+// m_h_ztrap->Fill(z_trap);
+// m_h_trap_drift_t->Fill(drift_time);
+// m_h_drift1->Fill(spess, drift_time / t_electrode);
+// m_h_gen->Fill(spess, drift_time);
+// m_h_gen1->Fill(spess, z_trap);
+// m_h_gen2->Fill(spess, z_trap / drift_time * t_electrode);
+// m_h_velocity_trap->Fill(electric_field, z_trap / drift_time);
+// m_h_mobility_trap->Fill(electric_field, z_trap / drift_time / electric_field);
+// m_h_trap_pos->Fill(trap_pos);
+// }
+// } else {
+// isTrapped = false;
+// if (m_doHistoTrap) {
+// const double z_trap{condData.getTrappingPositionZ()};
+// m_h_no_ztrap->Fill(z_trap);
+// m_h_notrap_drift_t->Fill(drift_time);
+// }
+// }
+// return isTrapped;
+// }
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_SurfaceChargesGenerator.h b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_SurfaceChargesGenerator.h
new file mode 100644
index 0000000000000000000000000000000000000000..f9ae5fcda18638ec41925745fc52f710fb0ad804
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/SCT_SurfaceChargesGenerator.h
@@ -0,0 +1,164 @@
+// -*- C++ -*-
+
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+/**
+ * Header file for SCT_SurfaceChargesGenerator Class
+ * @brief Class to drift the charged hits to the sensor surface for SCT
+ * version 1.0a Szymon Gadomski 31.05.2001
+ * @author Szymon.Gadomski@cern.ch, Awatif.Belymam@cern.ch, Davide.Costanzo@cern.ch,
+ * tgcornel@nikhef.nl, Grant.Gorfine@cern.ch, Paul.Bell@cern.ch, Jorgen.Dalmau@cern.ch,
+ * This is based on the SCTtest3SurfaceChargesGenerator
+ * compared to "test3", changes are all related to interfaces
+ * this should work with the new SCT_BarrelModuleSideDesign
+ * Compared to "test2" the following have been added in "test3":
+ * - a possibility to work with long steps and to smear charge along
+ * the step in smaller steps
+ * - Lorentz angle
+ * 19/04/2001 - change of hit class to SiTrackerHit
+ * version 1.0 Szymon Gadomski 09/06/2001, compiles with framework
+ * 26/07/2001 - changed to work with SiDigitization-00-02-00, tested
+ * 20/10/2002 - Awatif Belymam
+ * 29/10/2002 - the thickness of sensors is obtained from the objects
+ * - of SCT_ModuleSideDesign type instead from its jobOptions file.
+ * - Changes are done as well in jobOptions files.
+ * 06/01/2004 - Everything is now in CLHEP units (mm, ns, MeV)
+ * 23/08/2007 - Kondo.Gnanvo@cern.ch
+ * - Conversion of the SCT_SurfaceChargesGenerator code Al gTool
+ */
+
+#ifndef FASERSCT_DIGITIZATION_SCT_SURFACECHARGESGENERATOR_H
+#define FASERSCT_DIGITIZATION_SCT_SURFACECHARGESGENERATOR_H
+
+//Inheritance
+#include "AthenaBaseComps/AthAlgTool.h"
+#include "FaserSCT_Digitization/ISCT_SurfaceChargesGenerator.h"
+
+#include "Identifier/IdentifierHash.h"
+// #include "InDetConditionsSummaryService/ISiliconConditionsTool.h"
+#include "InDetCondTools/ISiLorentzAngleTool.h"
+// #include "SCT_ConditionsTools/ISCT_RadDamageSummaryTool.h"
+// #include "SCT_ModuleDistortions/ISCT_ModuleDistortionsTool.h"
+#include "SiPropertiesTool/ISiPropertiesTool.h"
+
+// Gaudi
+#include "GaudiKernel/ITHistSvc.h" // for ITHistSvc
+#include "GaudiKernel/ServiceHandle.h"
+#include "GaudiKernel/ToolHandle.h"
+
+#include <iostream>
+
+// Charges and hits
+class FaserSiHit;
+
+// -- do charge trapping histos
+class ITHistSvc;
+class TH1F;
+class TH2F;
+class TProfile;
+
+namespace TrackerDD {
+ class SiDetectorElement;
+ class SCT_ModuleSideDesign;
+}
+
+namespace CLHEP {
+ class HepRandomEngine;
+}
+
+template <class HIT> class TimedHitPtr;
+
+class SCT_SurfaceChargesGenerator : public extends<AthAlgTool, ISCT_SurfaceChargesGenerator> {
+ public:
+
+ /** constructor */
+ SCT_SurfaceChargesGenerator(const std::string& type, const std::string& name, const IInterface* parent);
+
+ /** Destructor */
+ virtual ~SCT_SurfaceChargesGenerator() = default;
+
+ /** AlgTool initialize */
+ virtual StatusCode initialize();
+
+ /** AlgTool finalize */
+ virtual StatusCode finalize();
+
+ private:
+
+ virtual void setFixedTime(float fixedTime) {m_tfix = fixedTime;}
+
+ /** create a list of surface charges from a hit */
+ virtual void process(const TrackerDD::SiDetectorElement* element, const TimedHitPtr<FaserSiHit>& phit, const ISiSurfaceChargesInserter& inserter, CLHEP::HepRandomEngine * rndmEngine) const;
+ void processSiHit(const TrackerDD::SiDetectorElement* element, const FaserSiHit& phit, const ISiSurfaceChargesInserter& inserter, float eventTime, unsigned short eventID, CLHEP::HepRandomEngine * rndmEngine) const;
+
+ // some diagnostics methods are needed here too
+ float driftTime(float zhit, const TrackerDD::SiDetectorElement* element) const; //!< calculate drift time perpandicular to the surface for a charge at distance zhit from mid gap
+ float diffusionSigma(float zhit, const TrackerDD::SiDetectorElement* element) const; //!< calculate diffusion sigma from a gaussian dist scattered charge
+ float surfaceDriftTime(float ysurf) const; //!< Calculate of the surface drift time
+ float maxDriftTime(const TrackerDD::SiDetectorElement* element) const; //!< max drift charge equivalent to the detector thickness
+ float maxDiffusionSigma(const TrackerDD::SiDetectorElement* element) const; //!< max sigma diffusion
+
+ // trap_pos and drift_time are updated based on spess.
+// bool chargeIsTrapped(double spess, const TrackerDD::SiDetectorElement* element, double& trap_pos, double& drift_time) const;
+
+ IntegerProperty m_numberOfCharges{this, "NumberOfCharges", 1, "number of charges"};
+ FloatProperty m_smallStepLength{this, "SmallStepLength", 5, "max internal step along the larger G4 step"};
+
+ /** related to the surface drift */
+ FloatProperty m_tSurfaceDrift{this, "SurfaceDriftTime", 10, "max surface drift time"};
+
+ FloatProperty m_tfix{this, "FixedTime", -999., "fixed time"};
+ FloatProperty m_tsubtract{this, "SubtractTime", -999., "subtract drift time from mid gap"};
+
+ BooleanProperty m_doDistortions{this, "doDistortions", false, "Simulation of module distortions"};
+ BooleanProperty m_useSiCondDB{this, "UseSiCondDB", false, "Usage of SiConditions DB values can be disabled to use setable ones"};
+ FloatProperty m_vdepl{this, "DepletionVoltage", 70., "depletion voltage, default 70V"};
+ FloatProperty m_vbias{this, "BiasVoltage", 150., "bias voltage, default 150V"};
+ BooleanProperty m_doTrapping{this, "doTrapping", false, "Flag to set Charge Trapping"};
+// BooleanProperty m_doHistoTrap{this, "doHistoTrap", false, "Histogram the charge trapping effect"};
+ BooleanProperty m_doRamo{this, "doRamo", false, "Ramo Potential for charge trapping effect"};
+ BooleanProperty m_isOverlay{this, "isOverlay", false, "flag for overlay"};
+
+ //ToolHandles
+// ToolHandle<ISCT_ModuleDistortionsTool> m_distortionsTool{this, "SCTDistortionsTool", "SCT_DistortionsTool", "Tool to retrieve SCT distortions"};
+ ToolHandle<ISiPropertiesTool> m_siPropertiesTool{this, "SiPropertiesTool", "SCT_SiPropertiesTool", "Tool to retrieve SCT silicon properties"};
+// ToolHandle<ISCT_RadDamageSummaryTool> m_radDamageTool{this, "RadDamageSummaryTool", "SCT_RadDamageSummaryTool", "Tool to retrieve SCT radiation damages"};
+// ToolHandle<ISiliconConditionsTool> m_siConditionsTool{this, "SiConditionsTool", "SCT_SiliconConditionsTool", "Tool to retrieve SCT silicon information"};
+ ToolHandle<ISiLorentzAngleTool> m_lorentzAngleTool{this, "LorentzAngleTool", "SiLorentzAngleTool/SCTLorentzAngleTool", "Tool to retreive Lorentz angle"};
+
+ ServiceHandle<ITHistSvc> m_thistSvc{this, "THistSvc", "THistSvc"};
+
+ float m_tHalfwayDrift{0.}; //!< Surface drift time
+ float m_distInterStrip{1.0}; //!< Inter strip distance normalized to 1
+ float m_distHalfInterStrip{0.}; //!< Half way distance inter strip
+
+ bool m_SurfaceDriftFlag{false}; //!< surface drift ON/OFF
+
+ // -- Histograms
+ TProfile* m_h_efieldz{nullptr};
+ TH1F* m_h_efield{nullptr};
+ TH1F* m_h_spess{nullptr};
+ TH1F* m_h_depD{nullptr};
+ TH2F* m_h_drift_electrode{nullptr};
+ TH1F* m_h_ztrap{nullptr};
+ TH1F* m_h_drift_time{nullptr};
+ TH1F* m_h_t_electrode{nullptr};
+ TH1F* m_h_zhit{nullptr};
+ TH1F* m_h_ztrap_tot{nullptr};
+ TH1F* m_h_no_ztrap{nullptr};
+ TH1F* m_h_trap_drift_t{nullptr};
+ TH1F* m_h_notrap_drift_t{nullptr};
+ TProfile* m_h_mob_Char{nullptr};
+ TProfile* m_h_vel{nullptr};
+ TProfile* m_h_drift1{nullptr};
+ TProfile* m_h_gen{nullptr};
+ TProfile* m_h_gen1{nullptr};
+ TProfile* m_h_gen2{nullptr};
+ TProfile* m_h_velocity_trap{nullptr};
+ TProfile* m_h_mobility_trap{nullptr};
+ TH1F* m_h_trap_pos{nullptr};
+};
+
+#endif // SCT_SURFACECHARGESGENERATOR_H
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/src/components/FaserSCT_Digitization_entries.cxx b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/components/FaserSCT_Digitization_entries.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..ea833080bce266204197ea69cdccc4c0752bf83d
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/src/components/FaserSCT_Digitization_entries.cxx
@@ -0,0 +1,18 @@
+#include "../SCT_Amp.h"
+#include "../SCT_FrontEnd.h"
+#include "../SCT_Digitization.h"
+#include "../SCT_DigitizationTool.h"
+#include "../SCT_SurfaceChargesGenerator.h"
+#include "../SCT_RandomDisabledCellGenerator.h"
+// Unused by ATLAS
+// #include "../SCT_DetailedSurfaceChargesGenerator.h"
+
+DECLARE_COMPONENT( SCT_Amp )
+DECLARE_COMPONENT( SCT_FrontEnd )
+DECLARE_COMPONENT( SCT_Digitization )
+DECLARE_COMPONENT( SCT_DigitizationTool )
+DECLARE_COMPONENT( SCT_SurfaceChargesGenerator )
+DECLARE_COMPONENT( SCT_RandomDisabledCellGenerator )
+// Unused by ATLAS
+// DECLARE_COMPONENT( SCT_DetailedSurfaceChargesGenerator )
+
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/test/FaserSCT_DigitizationDbg.py b/Tracker/TrackerDigitization/FaserSCT_Digitization/test/FaserSCT_DigitizationDbg.py
new file mode 100644
index 0000000000000000000000000000000000000000..dd51e8f88249c522103f5288ec363256114a02a5
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/test/FaserSCT_DigitizationDbg.py
@@ -0,0 +1,58 @@
+#!/usr/bin/env python
+"""Test various ComponentAccumulator Digitization configuration modules
+
+Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+"""
+import sys
+from AthenaCommon.Logging import log
+from AthenaCommon.Constants import DEBUG
+from AthenaCommon.Configurable import Configurable
+from CalypsoConfiguration.AllConfigFlags import ConfigFlags
+from AthenaConfiguration.AllConfigFlags import ConfigFlags
+from AthenaConfiguration.TestDefaults import defaultTestFiles
+from AthenaConfiguration.MainServicesConfig import MainServicesSerialCfg
+from AthenaPoolCnvSvc.PoolReadConfig import PoolReadCfg
+from AthenaPoolCnvSvc.PoolWriteConfig import PoolWriteCfg
+#from Digitization.DigitizationParametersConfig import writeDigitizationMetadata
+from FaserSCT_Digitization.FaserSCT_DigitizationConfigNew import SCT_DigitizationCfg
+#from MCTruthSimAlgs.RecoTimingConfig import MergeRecoTimingObjCfg
+
+# Set up logging and new style config
+log.setLevel(DEBUG)
+Configurable.configurableRun3Behavior = True
+
+# Configure
+ConfigFlags.Input.Files = ['g4.HITS.root']
+ConfigFlags.Output.RDOFileName = "myRDO.pool.root"
+ConfigFlags.IOVDb.GlobalTag = "OFLCOND-MC16-SDR-14"
+ConfigFlags.GeoModel.Align.Dynamic = False
+ConfigFlags.Concurrency.NumThreads = 1
+ConfigFlags.Beam.NumberOfCollisions = 0.
+
+ConfigFlags.GeoModel.FaserVersion = "FASER-00" # Always needed
+ConfigFlags.GeoModel.AtlasVersion = "ATLAS-R2-2016-01-00-01" # Always needed to fool autoconfig
+
+ConfigFlags.lock()
+
+# Core components
+acc = MainServicesSerialCfg()
+acc.merge(PoolReadCfg(ConfigFlags))
+acc.merge(PoolWriteCfg(ConfigFlags))
+#acc.merge(writeDigitizationMetadata(ConfigFlags))
+
+# Inner Detector
+acc.merge(SCT_DigitizationCfg(ConfigFlags))
+
+# Timing
+#acc.merge(MergeRecoTimingObjCfg(ConfigFlags))
+
+# Dump config
+acc.getService("StoreGateSvc").Dump = True
+acc.getService("ConditionStore").Dump = True
+acc.printConfig(withDetails=True)
+ConfigFlags.dump()
+# Execute and finish
+sc = acc.run(maxEvents=3)
+# Success should be 0
+sys.exit(not sc.isSuccess())
+
diff --git a/Tracker/TrackerDigitization/FaserSCT_Digitization/test/SCT_DigitizationConfigNew_test.py b/Tracker/TrackerDigitization/FaserSCT_Digitization/test/SCT_DigitizationConfigNew_test.py
new file mode 100644
index 0000000000000000000000000000000000000000..4fbd20f575761bc28029c4139df6b5cbd4abd382
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSCT_Digitization/test/SCT_DigitizationConfigNew_test.py
@@ -0,0 +1,37 @@
+#!/usr/bin/env python
+"""Run tests on SCT_DigitizationConfigNew.py
+
+Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+"""
+from AthenaCommon.Logging import log
+from AthenaCommon.Constants import DEBUG
+from AthenaCommon.Configurable import Configurable
+from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator
+from AthenaConfiguration.AllConfigFlags import ConfigFlags
+from AthenaConfiguration.MainServicesConfig import MainServicesSerialCfg
+from AthenaConfiguration.TestDefaults import defaultTestFiles
+from AthenaPoolCnvSvc.PoolReadConfig import PoolReadCfg
+from AtlasGeoModel.InDetGMConfig import InDetGeometryCfg
+from FaserSCT_Digitization.FaserSCT_DigitizationConfigNew import SCT_DigitizationHSCfg
+
+# Set up logging and new style config
+log.setLevel(DEBUG)
+Configurable.configurableRun3Behavior = True
+# Configure
+ConfigFlags.Input.Files = defaultTestFiles.HITS
+ConfigFlags.IOVDb.GlobalTag = "OFLCOND-MC16-SDR-16"
+ConfigFlags.GeoModel.Align.Dynamic = False
+ConfigFlags.Concurrency.NumThreads = 1
+ConfigFlags.lock()
+# Construct our accumulator to run
+acc = MainServicesSerialCfg()
+acc.merge(PoolReadCfg(ConfigFlags))
+acc.merge(SCT_DigitizationHSCfg(ConfigFlags))
+# Dump config
+acc.getService("StoreGateSvc").Dump = True
+acc.getService("ConditionStore").Dump = True
+acc.printConfig(withDetails=True)
+ConfigFlags.dump()
+# Execute and finish
+acc.run(maxEvents=3)
+
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/CMakeLists.txt b/Tracker/TrackerDigitization/FaserSiDigitization/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..afa4ae31542b62d1a6a15b5c9b76b16786d8bb16
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/CMakeLists.txt
@@ -0,0 +1,24 @@
+################################################################################
+# Package: FaserSiDigitization
+################################################################################
+
+# Declare the package name:
+atlas_subdir( FaserSiDigitization )
+
+# Declare the package's dependencies:
+atlas_depends_on_subdirs( PUBLIC
+ Control/AthenaKernel
+ Control/AthAllocators
+ DetectorDescription/Identifier
+ GaudiKernel
+ Tracker/TrackerDetDescr/TrackerReadoutGeometry
+ Tracker/TrackerSimEvent )
+
+# Component(s) in the package:
+atlas_add_library( FaserSiDigitization
+ src/SiChargedDiode.cxx
+ src/SiChargedDiodeCollection.cxx
+ src/SiSurfaceCharge.cxx
+ PUBLIC_HEADERS FaserSiDigitization
+ LINK_LIBRARIES AthenaKernel AthAllocators Identifier GaudiKernel TrackerReadoutGeometry TrackerSimEvent )
+
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/ATLAS_CHECK_THREAD_SAFETY b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/ATLAS_CHECK_THREAD_SAFETY
new file mode 100644
index 0000000000000000000000000000000000000000..da544c9d30a66d7dcb521188c6d363dc77919bd4
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/ATLAS_CHECK_THREAD_SAFETY
@@ -0,0 +1 @@
+Tracker/TrackerDigitization/FaserSiDigitization
\ No newline at end of file
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/ISiChargedDiodesProcessorTool.h b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/ISiChargedDiodesProcessorTool.h
new file mode 100644
index 0000000000000000000000000000000000000000..f45afd9fcd61f695e9d9f5f7b0bcc53d772e6c51
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/ISiChargedDiodesProcessorTool.h
@@ -0,0 +1,47 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+/**
+ * ISiChargedDiodesProcessorTool.h
+ * Header file for abstract base class ISiChargedDiodesProcessorTool
+ * (c) ATLAS Detector software
+ * Interface for the SiChargedDiode processor classes
+ * 23/08/2007 Kondo.Gnanvo@cern.ch, and others
+*/
+
+#ifndef FASERSIDIGITIZATION_ISICHARGEDDIODESPROCESSORTOOL_H
+#define FASERSIDIGITIZATION_ISICHARGEDDIODESPROCESSORTOOL_H
+
+#include "GaudiKernel/IAlgTool.h"
+
+class SiChargedDiodeCollection;
+class NPtGainSummaryData;
+namespace CLHEP {
+ class HepRandomEngine;
+}
+
+static const InterfaceID IID_ISiChargedDiodesProcessorTool("ISiChargedDiodesProcessorTool",1,0);
+
+class ISiChargedDiodesProcessorTool : virtual public IAlgTool {
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+public:
+
+ //Retrieve interface ID
+ static const InterfaceID & interfaceID() {return IID_ISiChargedDiodesProcessorTool; }
+
+ // Destructor:
+ virtual ~ISiChargedDiodesProcessorTool() {}
+
+ ///////////////////////////////////////////////////////////////////
+ // Pure virtual methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // process the collection of charged diodes
+ virtual void process(SiChargedDiodeCollection &collection, CLHEP::HepRandomEngine * rndmEngine) const =0;
+};
+
+#endif // FASERSIDIGITIZATION_ISICHARGEDDIODESPROCESSORTOOL_H
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiChargedDiode.h b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiChargedDiode.h
new file mode 100644
index 0000000000000000000000000000000000000000..fadd13e4555c63976bee38629119e80f179a0755
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiChargedDiode.h
@@ -0,0 +1,148 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SiChargedDiode.h // Header file for class SiChargedDiode
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Class which contains the sum and the composition of all bare
+//charges (SiTotalCharge) deposited in one SiDiode
+///////////////////////////////////////////////////////////////////
+// Version 2.1 09/06/2001 David Calvet
+// Revisited 04-03-03 Davide Costanzo
+// added a int flag as a private data member to store the noise,
+// disconnected, bad_tot information. The relative bunch number is
+// also stored in this word.
+// the word will is meant to be copied as it is in the SDO
+///////////////////////////////////////////////////////////////////
+
+#ifndef FASERSIDIGITIZATION_SICHARGEDDIODE_H
+#define FASERSIDIGITIZATION_SICHARGEDDIODE_H
+
+// Data member classes
+#include "TrackerSimEvent/SiTotalCharge.h"
+#include "TrackerReadoutGeometry/SiCellId.h"
+#include "TrackerReadoutGeometry/SiReadoutCellId.h"
+
+class SiHelper; // used to set the flag word
+//class SiChargedDiode;
+
+class SiChargedDiode {
+
+ friend class SiHelper;
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+public:
+
+ // Constructor with parameters:
+ SiChargedDiode(const SiTotalCharge::alloc_t& alloc,
+ const TrackerDD::SiCellId & diode, const TrackerDD::SiReadoutCellId & roCell, int flagword=0, SiChargedDiode* nextInCluster=NULL);
+
+ ///////////////////////////////////////////////////////////////////
+ // Const methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // Diode which contains this charge:
+ const TrackerDD::SiCellId & diode() const;
+
+ // Readout cell associated to diode
+ // (usually the same id as the diode except for ganged pixels):
+ const TrackerDD::SiReadoutCellId & getReadoutCell();
+
+ // total charge and its composition:
+ const SiTotalCharge & totalCharge() const;
+
+ // total deposited charge:
+ double charge() const;
+
+ // flag, disconnected etc.
+ int flag() const;
+
+ //neighbouring strip for navigation
+ SiChargedDiode * nextInCluster();
+ ///////////////////////////////////////////////////////////////////
+ // Non-const methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // add another charge:
+ void add(const SiCharge &charge);
+ // add a total charge
+ void add(const SiTotalCharge &totcharge);
+ //add a neighbouring strip for navigation
+ void setNextInCluster(SiChargedDiode* nextInCluster);
+
+ ///////////////////////////////////////////////////////////////////
+ // Private methods:
+ ///////////////////////////////////////////////////////////////////
+private:
+
+ SiChargedDiode();
+
+ ///////////////////////////////////////////////////////////////////
+ // Private data:
+ ///////////////////////////////////////////////////////////////////
+private:
+ TrackerDD::SiCellId m_diode; // SiDiode which contains this charge
+ // the pointed SiDiode is owned by the SiChargedDiode
+ SiTotalCharge m_totalCharge; // total charge and its composition
+ TrackerDD::SiReadoutCellId m_readoutCell; //Readout cell associated to this diode
+ int m_word; // a flag for noise etc etc as in InDetSimData
+ SiChargedDiode * m_nextInCluster; //the next strip to navigate to - allows traversing clusters since the SiChargedDiodeCollection is not guaranteed to be contiguous
+};
+
+///////////////////////////////////////////////////////////////////
+// Inline methods:
+///////////////////////////////////////////////////////////////////
+inline const TrackerDD::SiCellId & SiChargedDiode::diode() const
+{
+ return m_diode;
+}
+
+inline int SiChargedDiode::flag() const
+{
+ return m_word;
+}
+
+inline const SiTotalCharge &SiChargedDiode::totalCharge() const
+{
+ return m_totalCharge;
+}
+inline const TrackerDD::SiReadoutCellId & SiChargedDiode::getReadoutCell() {
+ return m_readoutCell;
+}
+
+inline double SiChargedDiode::charge() const
+{
+ return m_totalCharge.charge();
+}
+
+inline SiChargedDiode * SiChargedDiode::nextInCluster()
+{
+ return m_nextInCluster;
+}
+
+inline void SiChargedDiode::add(const SiCharge &charge)
+{
+ m_totalCharge.add(charge);
+}
+
+inline void SiChargedDiode::add(const SiTotalCharge &totcharge)
+{
+ m_totalCharge.add(totcharge);
+}
+
+inline void SiChargedDiode::setNextInCluster(SiChargedDiode* nextInCluster)
+{
+ m_nextInCluster = nextInCluster;
+}
+///////////////////////////////////////////////////////////////////
+// Input/Output stream functions:
+///////////////////////////////////////////////////////////////////
+std::ostream &operator<<(std::ostream &out,const SiChargedDiode &chargedDiode);
+
+#endif // FASERSIDIGITIZATION_SICHARGEDDIODE_H
+
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiChargedDiodeCollection.h b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiChargedDiodeCollection.h
new file mode 100644
index 0000000000000000000000000000000000000000..4f5a9f23734cfc1d48545565e58f5b36e90fe473
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiChargedDiodeCollection.h
@@ -0,0 +1,285 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SiChargedDiodeCollection.h
+// Header file for class SiChargedDiodeCollection
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Class used to store the list of SiChargedDiode objects for one
+// module (Pixel) or half module (SCT) and one event.
+///////////////////////////////////////////////////////////////////
+// Version 3.0 09/07/2001 David Calvet
+// Davide Costanzo. Revisited version. 04-03-03
+// - added a Identifier wafer_id to the private members
+// - added a IdHelper to the private members
+// - constructor modified to initialize wafer_id and IdHelper
+// - replaced <list> with <map> and use the compact id of the
+// SiChargedDiode to map them.
+// - Inherit from Identifiable to enforce the identify() method
+///////////////////////////////////////////////////////////////////
+#ifndef FASERSIDIGITIZATION_SICHARGEDDIODECOLLECTION_H
+#define FASERSIDIGITIZATION_SICHARGEDDIODECOLLECTION_H
+
+// Base class
+#include "Identifier/Identifiable.h"
+
+// Data member classes
+#include <unordered_map>
+#include "FaserSiDigitization/SiChargedDiode.h"
+#include "Identifier/Identifier.h"
+#include "TrackerReadoutGeometry/SiDetectorElement.h"
+
+// Input/output classes
+#include "TrackerSimEvent/FaserSiHit.h"
+
+// STL includes
+#include <atomic>
+#include <set>
+#include <memory>
+
+class FaserDetectorID;
+namespace TrackerDD{
+ class SiDetectorElement;
+ class SiDetectorDesign;
+ class SiCellId;
+}
+
+#include "AthAllocators/ArenaPoolSTLAllocator.h"
+
+struct SiChargedDiodeHash
+{
+ size_t operator() (const TrackerDD::SiCellId & id) const
+ {
+ return id.word();
+ }
+};
+
+typedef std::unordered_map<TrackerDD::SiCellId,
+ SiChargedDiode,
+ SiChargedDiodeHash,
+ std::equal_to<TrackerDD::SiCellId>,
+ SG::ArenaPoolSTLAllocator<
+ std::pair<const TrackerDD::SiCellId, SiChargedDiode> > >
+ SiChargedDiodeMap;
+
+
+// Iterator typedef to make it easier
+typedef SiChargedDiodeMap::iterator SiChargedDiodeIterator;
+
+
+
+//
+// A normal iteration over a SiChargedDiodeCollection will use
+// the unordered_map, so the ordering is not defined. The observed
+// ordering can (and does) change depending on the compiler and library
+// version used. In some cases, though, we are sensitive to the
+// order of iteration, for example in cases where we generate a random
+// number for each element of the collection. To get results which
+// are identical across compilers, we need to instead do the iteration
+// in a well-defined order.
+//
+// This can be done using the methods sortedBegin and sortedEnd instead
+// of begin and end. These work by maintaining a std::set of pointers
+// to the diodes, sorted by the diode number. In order to avoid paying
+// the penalty for maintaining the sorted set when we don't need to, we only
+// start maintaining it the first time that it's requested.
+struct SiChargedDiodeOrderedSetCompare
+{
+ size_t operator() (const SiChargedDiode* a,
+ const SiChargedDiode* b) const
+ {
+ return a->diode().word() < b->diode().word();
+ }
+};
+
+
+
+typedef std::set<SiChargedDiode*,
+ SiChargedDiodeOrderedSetCompare,
+ SG::ArenaPoolSTLAllocator<SiChargedDiode*> >
+ SiChargedDiodeOrderedSet;
+
+
+// Iterator typedef to make it easier
+typedef SiChargedDiodeOrderedSet::iterator SiChargedDiodeOrderedIterator;
+
+class SiChargedDiodeCollection : Identifiable {
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+ public:
+
+
+ // Constructor with parameters:
+ // ref. to the detector element for this collection
+ SiChargedDiodeCollection( );
+
+ SiChargedDiodeCollection(const TrackerDD::SiDetectorElement* );
+
+
+ // Destructor:
+ ~SiChargedDiodeCollection();
+
+ ///////////////////////////////////////////////////////////////////
+ // Const methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // detector element:
+ const TrackerDD::SiDetectorElement * element() const;
+
+ // wafer identifier for this collection
+ virtual Identifier identify() const override final;
+ virtual IdentifierHash identifyHash() const override final;
+
+ // id helper for this collection
+ const FaserDetectorID* id_helper();
+
+ // detector design:
+ const TrackerDD::SiDetectorDesign &design() const;
+
+ // translation from SiReadoutCellId to Identifier
+ Identifier getId(const TrackerDD::SiCellId& id) const
+ {
+ return (element()->identifierFromCellId(id));
+ }
+
+ ///////////////////////////////////////////////////////////////////
+ // Non-const methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // clean up:
+ void clear();
+
+ // read/write access to the collection:
+ SiChargedDiodeMap &chargedDiodes();
+
+ // Set the SiDetectorElement
+ void setDetectorElement(const TrackerDD::SiDetectorElement *SiElement);
+
+ // Add a new SiCharge to the collection
+ // (add or merge in an existing SiChargedDiode):
+ // Diode which contains the new charge
+ // SiCharge to be added.
+ void add(const TrackerDD::SiCellId & diode, const SiCharge &charge);
+ void add(const TrackerDD::SiCellId & diode, const SiTotalCharge &totcharge);
+
+ bool AlreadyHit(const TrackerDD::SiCellId & siId);
+ bool AlreadyHit(const Identifier & id);
+ SiChargedDiodeIterator begin();
+ SiChargedDiodeIterator end();
+ SiChargedDiodeOrderedIterator orderedBegin();
+ SiChargedDiodeOrderedIterator orderedEnd();
+ bool empty() const; // Test if there is anything in the collection.
+
+ // return a Charged diode given its CellId, NULL if doesn't exist
+ SiChargedDiode * find(const TrackerDD::SiCellId & siId);
+ // return a Charged diode given its identifier, NULL if doesn't exist
+ SiChargedDiode * find(Identifier);
+
+ ///////////////////////////////////////////////////////////////////
+ // Private methods:
+ ///////////////////////////////////////////////////////////////////
+ private:
+ SiChargedDiodeCollection(const SiChargedDiodeCollection&);
+ SiChargedDiodeCollection &operator=(const SiChargedDiodeCollection&);
+ void order();
+
+ ///////////////////////////////////////////////////////////////////
+ // Private data:
+ ///////////////////////////////////////////////////////////////////
+ private:
+ //NB m_allocator should always be declared before m_chargedDiodes in
+ //the intialization list. If the allocator is declared after
+ //m_chargedDiodes, when the collection is destroyed, the allocator
+ //will be destroyed (and the memory it manages freed) before the
+ //SiChargedDiodeMap. This will cause a crash unless the
+ //SiChargedDiodeMap is empty.
+ SiTotalCharge::alloc_t m_allocator;
+ SiChargedDiodeMap m_chargedDiodes; // list of SiChargedDiodes
+ SiChargedDiodeOrderedSet m_orderedChargedDiodes; // list of SiChargedDiodes
+ const TrackerDD::SiDetectorElement* m_sielement; // detector element
+};
+
+///////////////////////////////////////////////////////////////////
+// Inline methods:
+///////////////////////////////////////////////////////////////////
+
+// Set the DetectorElement
+inline void SiChargedDiodeCollection::setDetectorElement(const TrackerDD::SiDetectorElement *SiElement)
+{
+ m_sielement=SiElement;
+}
+
+inline SiChargedDiodeMap &SiChargedDiodeCollection::chargedDiodes()
+{
+ return m_chargedDiodes;
+}
+
+// access to the element
+inline const TrackerDD::SiDetectorElement *SiChargedDiodeCollection::element() const
+{
+ return m_sielement;
+}
+
+// access to the design
+inline const TrackerDD::SiDetectorDesign &SiChargedDiodeCollection::design() const
+{
+ return m_sielement->design();
+}
+
+// Get the Identifier of the collection
+inline Identifier SiChargedDiodeCollection::identify() const
+{
+ return m_sielement->identify();
+}
+
+// Get the Identifier of the collection
+inline IdentifierHash SiChargedDiodeCollection::identifyHash() const
+{
+ return m_sielement->identifyHash();
+}
+
+// Get the Id Helper
+inline const FaserDetectorID* SiChargedDiodeCollection::id_helper()
+{
+ return m_sielement->getIdHelper();
+}
+
+
+inline SiChargedDiodeIterator SiChargedDiodeCollection::begin()
+{
+ return m_chargedDiodes.begin();
+}
+
+inline SiChargedDiodeIterator SiChargedDiodeCollection::end()
+{
+ return m_chargedDiodes.end();
+}
+
+inline SiChargedDiodeOrderedIterator SiChargedDiodeCollection::orderedBegin()
+{
+ if (m_orderedChargedDiodes.empty() && !m_chargedDiodes.empty()) {
+ order();
+ }
+ return m_orderedChargedDiodes.begin();
+}
+
+inline SiChargedDiodeOrderedIterator SiChargedDiodeCollection::orderedEnd()
+{
+ if (m_orderedChargedDiodes.empty() && !m_chargedDiodes.empty()) {
+ order();
+ }
+ return m_orderedChargedDiodes.end();
+}
+
+inline bool SiChargedDiodeCollection::empty() const {
+ return m_chargedDiodes.empty();
+}
+
+
+
+#endif // FASERSIDIGITIZATION_SICHARGEDDIODECOLLECTION_H
+
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiHelper.h b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiHelper.h
new file mode 100644
index 0000000000000000000000000000000000000000..1745878defc39fd913365ffba0702e2ff1048af2
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiHelper.h
@@ -0,0 +1,208 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+/*************************************************************************
+ SiHelper for Silicons
+ This is used by the SiDigitization
+ ---------------------------------
+ ATLAS Collaboration
+
+
+ This is a short copy of the PixelSimHelper used by the SDO.
+ It is needed to set the word of the Charged Diodes in the same
+ way as it will be registered in the SDO.
+
+ D. Costanzo 4-3-03
+***************************************************************************/
+
+
+#ifndef FASERSIDIGITIZATION_SIHELPER_H
+#define FASERSIDIGITIZATION_SIHELPER_H
+
+#include "GaudiKernel/MsgStream.h"
+#include "AthenaKernel/getMessageSvc.h"
+#include "FaserSiDigitization/SiChargedDiode.h"
+
+class SiHelper
+{
+ public:
+ // methods to set characteristics of a new object
+
+ static void noise(SiChargedDiode& chDiode, bool flag, bool mask=false);
+ static void belowThreshold(SiChargedDiode& chDiode, bool flag, bool mask=false);
+ static void disabled(SiChargedDiode& chDiode, bool flag, bool mask=false);
+ static void badToT(SiChargedDiode& chDiode, bool flag, bool mask=false);
+ static void disconnected(SiChargedDiode& chDiode, bool flag, bool mask=false);
+ static void maskOut(SiChargedDiode& chDiode, bool flag);
+ static void ClusterUsed(SiChargedDiode& chDiode, bool flag);
+ static void SetBunch(SiChargedDiode& chDiode, int bunch, MsgStream* log=nullptr);
+ static void SetStripNum(SiChargedDiode& chDiode, int nstrip, MsgStream* log=nullptr);
+ static void SetTimeBin(SiChargedDiode& chDiode, int time, MsgStream* log=nullptr);
+
+ static bool isUsable(SiChargedDiode& chDiode);
+ static bool isNoise(SiChargedDiode& chDiode);
+ static bool isBelowThreshold(SiChargedDiode& chDiode);
+ static bool isDisabled(SiChargedDiode& chDiode);
+ static bool isBadToT(SiChargedDiode& chDiode);
+ static bool isDisconnected(SiChargedDiode& chDiode);
+ static bool isMaskOut(SiChargedDiode& chDiode);
+ static bool isClusterUsed(SiChargedDiode &chDiode);
+ static int GetBunch(SiChargedDiode& chDiode);
+ static int GetStripNum(SiChargedDiode& chDiode);
+ static int GetTimeBin(SiChargedDiode& chDiode);
+
+ private:
+ enum {
+ NOISE_SET = 0x1, // main charge is noise
+ BT_SET = 0x2, // below threshold
+ DISABLED_SET = 0x4, // disabled
+ BADTOT_SET = 0x8, // bad TOT
+ DISCONNECTED_SET = 0x10, // disconnected
+ CLUSTERUSED_SET = 0x20, // used in cluster - used by anyone???
+ MASKOUT_SET = 0x40 // charge diode is masked out, not to be used for RDO creation
+ };
+};
+
+inline void SiHelper::maskOut(SiChargedDiode& chDiode, bool flag) {
+ if (flag) {
+ chDiode.m_word |= MASKOUT_SET;
+ } else {
+ chDiode.m_word &= ~MASKOUT_SET;
+ }
+}
+
+inline void SiHelper::noise(SiChargedDiode& chDiode, bool flag, bool mask) {
+ if (flag) {
+ chDiode.m_word |= NOISE_SET;
+ } else {
+ chDiode.m_word &= ~NOISE_SET;
+ }
+ if (mask) SiHelper::maskOut(chDiode,true);
+}
+
+inline void SiHelper::belowThreshold(SiChargedDiode& chDiode, bool flag, bool mask) {
+ if (flag) {
+ chDiode.m_word |= BT_SET;
+ } else {
+ chDiode.m_word &= ~BT_SET;
+ }
+ if (mask) SiHelper::maskOut(chDiode,true);
+}
+
+inline void SiHelper::disabled(SiChargedDiode& chDiode, bool flag, bool mask) {
+ if (flag) {
+ chDiode.m_word |= DISABLED_SET;
+ } else {
+ chDiode.m_word &= ~DISABLED_SET;
+ }
+ if (mask) SiHelper::maskOut(chDiode,true);
+}
+
+inline void SiHelper::badToT(SiChargedDiode& chDiode, bool flag, bool mask) {
+ if (flag) {
+ chDiode.m_word |= BADTOT_SET;
+ } else {
+ chDiode.m_word &= ~BADTOT_SET;
+ }
+ if (mask) SiHelper::maskOut(chDiode,true);
+}
+
+inline void SiHelper::disconnected(SiChargedDiode& chDiode, bool flag, bool mask) {
+ if (flag) {
+ chDiode.m_word |= DISCONNECTED_SET;
+ } else {
+ chDiode.m_word &= ~DISCONNECTED_SET;
+ }
+ if (mask) SiHelper::maskOut(chDiode,true);
+}
+
+
+inline void SiHelper::ClusterUsed(SiChargedDiode& chDiode, bool flag) {
+ if (flag) {
+ chDiode.m_word |= CLUSTERUSED_SET;
+ } else {
+ chDiode.m_word &= ~CLUSTERUSED_SET;
+ }
+}
+
+inline void SiHelper::SetBunch(SiChargedDiode& chDiode, int bunch, MsgStream* log) {
+ //
+ // Code the bunch number in the 8 bits set corresponding to xx in xx00
+ //
+ if (bunch > 0xff) {
+ if (log) (*log) << MSG::ERROR << "Bunch Number not allowed" << endmsg;
+ }
+ chDiode.m_word = chDiode.m_word | ( (bunch&0xff) <<8 ) ;
+}
+
+inline void SiHelper::SetStripNum(SiChargedDiode& chDiode, int nstrip, MsgStream* log) {
+ //
+ // Code the number of strips in the 12 bits set corresponding to xxx in 0xxx0000
+ //
+ if (nstrip > 0xfff) {
+ if (log) (*log) << MSG::ERROR << "Number of strips not allowed" << endmsg;
+ }
+ chDiode.m_word = chDiode.m_word | ((nstrip&0xfff) << 16 ) ;
+}
+
+inline void SiHelper::SetTimeBin(SiChargedDiode& chDiode, int time, MsgStream* log) {
+ //
+ // Code the SCT Timebin number in the 3 bits set corresponding to x in x0000000
+ //
+ if (time > 0xf) {
+ if (log) (*log) << MSG::ERROR << "TimeBin not allowed" << endmsg;
+ }
+ chDiode.m_word = chDiode.m_word | ( (time&0xf) <<28 ) ;
+}
+
+////////////////////////////////////////////////////////////
+
+inline bool SiHelper::isUsable(SiChargedDiode& chDiode) {
+ return !SiHelper::isMaskOut(chDiode);
+ // return (((chDiode.m_word & 0xff) == 0) || // no special status bits sets, either track, xtalk or random noise
+ // (SiHelper::isClusterUsed(chDiode)) ); // cluster used - whatever that means - not set in digitization
+}
+
+inline bool SiHelper::isNoise(SiChargedDiode& chDiode) {
+ return (chDiode.m_word & NOISE_SET);
+}
+
+inline bool SiHelper::isMaskOut(SiChargedDiode& chDiode) {
+ return (chDiode.m_word & MASKOUT_SET);
+}
+
+inline bool SiHelper::isBelowThreshold(SiChargedDiode& chDiode) {
+ return (chDiode.m_word & BT_SET);
+}
+
+inline bool SiHelper::isDisabled(SiChargedDiode& chDiode) {
+ return (chDiode.m_word & DISABLED_SET);
+}
+
+inline bool SiHelper::isBadToT(SiChargedDiode& chDiode) {
+ return (chDiode.m_word & BADTOT_SET);
+}
+
+inline bool SiHelper::isDisconnected(SiChargedDiode& chDiode) {
+ return (chDiode.m_word & DISCONNECTED_SET);
+}
+
+inline bool SiHelper::isClusterUsed(SiChargedDiode& chDiode) {
+ return (chDiode.m_word & CLUSTERUSED_SET);
+}
+
+inline int SiHelper::GetBunch(SiChargedDiode& chDiode) {
+return ( (chDiode.m_word >> 8) & 0xff );
+}
+
+inline int SiHelper::GetStripNum(SiChargedDiode& chDiode) {
+ return ( (chDiode.m_word >> 16) & 0xfff );
+}
+
+inline int SiHelper::GetTimeBin(SiChargedDiode& chDiode) {
+return ( (chDiode.m_word >> 28) & 0xf );
+}
+
+
+#endif
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiSurfaceCharge.h b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiSurfaceCharge.h
new file mode 100644
index 0000000000000000000000000000000000000000..f88065cf62253c5a7a9441649d007307f94406b2
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/FaserSiDigitization/SiSurfaceCharge.h
@@ -0,0 +1,86 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SiSurfaceCharge.h
+// Header file for class SiSurfaceCharge
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Class which contains a charge located at the sensor surface
+///////////////////////////////////////////////////////////////////
+// Version 1.0 09/02/2001 David Calvet
+///////////////////////////////////////////////////////////////////
+
+#ifndef FASERSIDIGITIZATION_SISURFACECHARGE_H
+#define FASERSIDIGITIZATION_SISURFACECHARGE_H
+
+//#include "SiTrackerDetDescr/SiLocalPosition.h"
+#include "TrackerReadoutGeometry/SiLocalPosition.h"
+#include "TrackerSimEvent/SiCharge.h"
+
+
+class SiSurfaceCharge {
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+public:
+
+ // Copy constructor:
+ SiSurfaceCharge(const SiSurfaceCharge &surfaceCharge);
+
+ // Constructor with parameters:
+ // position at the sensor surface
+ // charge located at this position
+ SiSurfaceCharge(const TrackerDD::SiLocalPosition &position,const SiCharge &charge);
+
+ // Destructor:
+ ~SiSurfaceCharge();
+
+ // Assignment operator:
+ SiSurfaceCharge &operator=(const SiSurfaceCharge &surfaceCharge);
+
+ ///////////////////////////////////////////////////////////////////
+ // Const methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // position at the sensor surface:
+ TrackerDD::SiLocalPosition position() const;
+
+ // charge located at this position:
+ SiCharge charge() const;
+
+ ///////////////////////////////////////////////////////////////////
+ // Private methods:
+ ///////////////////////////////////////////////////////////////////
+private:
+
+ SiSurfaceCharge();
+
+ ///////////////////////////////////////////////////////////////////
+ // Private data:
+ ///////////////////////////////////////////////////////////////////
+private:
+ TrackerDD::SiLocalPosition m_position; // position at the sensor surface
+ SiCharge m_charge; // charge located at this position
+};
+
+///////////////////////////////////////////////////////////////////
+// Inline methods:
+///////////////////////////////////////////////////////////////////
+inline SiSurfaceCharge::~SiSurfaceCharge()
+{}
+
+inline TrackerDD::SiLocalPosition SiSurfaceCharge::position() const
+{
+ return m_position;
+}
+
+inline SiCharge SiSurfaceCharge::charge() const
+{
+ return m_charge;
+}
+
+#endif // FASERSIDIGITIZATION_SISURFACECHARGE_H
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/src/SiChargedDiode.cxx b/Tracker/TrackerDigitization/FaserSiDigitization/src/SiChargedDiode.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..622090839b51c6e19b1c04eeb5917b1509810b14
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/src/SiChargedDiode.cxx
@@ -0,0 +1,38 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SiChargedDiode.cxx
+// Implementation file for class SiChargedDiode
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Version 2.1 09/06/2001 David Calvet
+// Davide Costanzo. Revisited version. 04-03-03
+///////////////////////////////////////////////////////////////////
+
+#include "FaserSiDigitization/SiChargedDiode.h"
+
+
+// Constructor with parameters:
+SiChargedDiode::SiChargedDiode(const SiTotalCharge::alloc_t& alloc,
+ const TrackerDD::SiCellId & diode, const TrackerDD::SiReadoutCellId & roCell, int flagword, SiChargedDiode * nextInCluster)
+ : m_diode(diode),
+ m_totalCharge(alloc),
+ m_readoutCell(roCell),
+ m_word(flagword),
+ m_nextInCluster(nextInCluster)
+{}
+
+
+///////////////////////////////////////////////////////////////////
+// Input/Output stream functions:
+///////////////////////////////////////////////////////////////////
+std::ostream &operator<<(std::ostream &out,const SiChargedDiode &chargedDiode)
+{
+ out << "Diode=" << chargedDiode.diode()
+ << " " << chargedDiode.totalCharge();
+ return out;
+}
+
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/src/SiChargedDiodeCollection.cxx b/Tracker/TrackerDigitization/FaserSiDigitization/src/SiChargedDiodeCollection.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..82ae6cb5562854c2225d8e89832128f70dd14b88
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/src/SiChargedDiodeCollection.cxx
@@ -0,0 +1,161 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// sichargeddiodecollection.cxx
+// Implementation file for class SiChargedDiodeCollection
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Version 3.0 09/07/2001 David Calvet
+// 04-03-03 Revisited version. Davide Costanzo
+///////////////////////////////////////////////////////////////////
+// header file
+#include "FaserSiDigitization/SiChargedDiodeCollection.h"
+// member classes
+#include "FaserSiDigitization/SiHelper.h"
+#include "TrackerReadoutGeometry/SiDetectorDesign.h"
+#include "TrackerReadoutGeometry/SiReadoutCellId.h"
+#include "TrackerReadoutGeometry/SiCellId.h"
+#include "GaudiKernel/MsgStream.h"
+#include "AthenaKernel/getMessageSvc.h"
+
+using namespace TrackerDD;
+
+
+
+SiChargedDiodeCollection::SiChargedDiodeCollection( ) :
+ m_chargedDiodes(),
+ m_sielement()
+{
+}
+
+SiChargedDiodeCollection::SiChargedDiodeCollection(const TrackerDD::SiDetectorElement* sielement ) :
+ m_chargedDiodes(),
+ m_sielement(sielement)
+{
+}
+
+
+SiChargedDiodeCollection::~SiChargedDiodeCollection()
+{
+}
+
+// Clean up the collection
+void SiChargedDiodeCollection::clear() {
+ m_sielement = 0;
+ m_chargedDiodes.erase(m_chargedDiodes.begin(), m_chargedDiodes.end() );
+ m_orderedChargedDiodes.clear();
+}
+
+// Add a new SiCharge to the collection
+void SiChargedDiodeCollection::add(const SiCellId & diode,
+ const SiCharge & charge)
+{
+ // check the pointer is correct
+ if (!diode.isValid()) return;
+
+ // find this diode in the charged diode collection
+ // find it by id in the map.
+ //
+
+ SiChargedDiodeIterator the_diode = m_chargedDiodes.find(diode);
+
+ if(the_diode != m_chargedDiodes.end() ) {
+ // Add to existing charge
+ (*the_diode).second.add(charge);
+ } else {
+ // if the new diode has not been found in the collection
+ // get the read out cell from the design.
+ //
+ SiReadoutCellId roCell=design().readoutIdOfCell(diode);
+ if (!roCell.isValid()) { // I don't think this can occur at this stage but cant hurt.
+ MsgStream log(Athena::getMessageSvc(),"SiChargedDiodeCollection");
+ log << MSG::FATAL << "Could not create SiReadoutCellId object !"<< endmsg;
+ }
+ // create a new charged diode
+ SiChargedDiode chargedDiode(m_allocator, diode,roCell);
+ // add the new charge to it
+ chargedDiode.add(charge);
+ if (charge.processType() == SiCharge::extraNoise) SiHelper::noise(chargedDiode,true);
+ // add the new charged diode to the charged diode collection
+ auto p = m_chargedDiodes.emplace(diode,chargedDiode);
+ if (!m_orderedChargedDiodes.empty()) {
+ m_orderedChargedDiodes.insert (&p.first->second);
+ }
+ }
+}
+
+// Add a new SiTotalCharge to the collection
+void SiChargedDiodeCollection::add(const SiCellId & diode,
+ const SiTotalCharge & totcharge)
+{
+ // check the pointer is correct
+ if (!diode.isValid()) return;
+
+ // find this diode in the charged diode collection
+ // find it by id in the map.
+ //
+
+ SiChargedDiodeIterator the_diode = m_chargedDiodes.find(diode);
+
+ if(the_diode != m_chargedDiodes.end() ) {
+ // Add to existing charge
+ (*the_diode).second.add(totcharge);
+ } else {
+ // if the new diode has not been found in the collection
+ // get the read out cell from the design.
+ //
+ SiReadoutCellId roCell=design().readoutIdOfCell(diode);
+ if (!roCell.isValid()) { // I don't think this can occur at this stage but cant hurt.
+ MsgStream log(Athena::getMessageSvc(),"SiChargedDiodeCollection");
+ log << MSG::FATAL << "Could not create SiReadoutCellId object !"<< endmsg;
+ }
+ // create a new charged diode
+ SiChargedDiode chargedDiode(m_allocator, diode,roCell);
+ // add the new charge to it
+ chargedDiode.add(totcharge);
+ // add the new charged diode to the charged diode collection
+ auto p = m_chargedDiodes.emplace(diode,chargedDiode);
+ if (!m_orderedChargedDiodes.empty()) {
+ m_orderedChargedDiodes.insert (&p.first->second);
+ }
+ }
+}
+
+bool SiChargedDiodeCollection::AlreadyHit(const TrackerDD::SiCellId & siId) {
+ if(m_chargedDiodes.find(siId) == m_chargedDiodes.end() ) {
+ return false;
+ } else {
+ return true;
+ }
+}
+
+bool SiChargedDiodeCollection::AlreadyHit(const Identifier & id) {
+
+ const SiCellId cellId = m_sielement->cellIdFromIdentifier(id);
+ return AlreadyHit(cellId);
+}
+
+SiChargedDiode * SiChargedDiodeCollection::find(const TrackerDD::SiCellId & siId) {
+ // get the compact Id to access the map
+ SiChargedDiodeIterator theEl = m_chargedDiodes.find(siId);
+ // if the diode exists return a pointer to it:
+ if (theEl == m_chargedDiodes.end() ) return NULL;
+ else return &( (*theEl).second);
+}
+
+SiChargedDiode * SiChargedDiodeCollection::find(Identifier siId) {
+
+ // Get the key for the map lookup
+ const SiCellId cellId = m_sielement->cellIdFromIdentifier(siId);
+ return find(cellId);
+}
+
+void SiChargedDiodeCollection::order()
+{
+ for (auto& p : m_chargedDiodes) {
+ m_orderedChargedDiodes.insert (&p.second);
+ }
+}
diff --git a/Tracker/TrackerDigitization/FaserSiDigitization/src/SiSurfaceCharge.cxx b/Tracker/TrackerDigitization/FaserSiDigitization/src/SiSurfaceCharge.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..b1d5215f5c9a9d9d1f647eb41a0c47a6613d4fed
--- /dev/null
+++ b/Tracker/TrackerDigitization/FaserSiDigitization/src/SiSurfaceCharge.cxx
@@ -0,0 +1,38 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SiSurfaceCharge.cxx
+// Implementation file for class SiSurfaceCharge
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Version 1.0 01/02/2001 David Calvet
+///////////////////////////////////////////////////////////////////
+
+#include "FaserSiDigitization/SiSurfaceCharge.h"
+using namespace TrackerDD;
+// Copy constructor:
+SiSurfaceCharge::SiSurfaceCharge(const SiSurfaceCharge &surfaceCharge) :
+ m_position(surfaceCharge.m_position),
+ m_charge(surfaceCharge.m_charge)
+{}
+
+// Constructor with parameters:
+SiSurfaceCharge::SiSurfaceCharge(const SiLocalPosition &position,
+ const SiCharge &charge) :
+ m_position(position),
+ m_charge(charge)
+{}
+
+// Assignment operator:
+SiSurfaceCharge &SiSurfaceCharge::operator=(const SiSurfaceCharge &surfaceCharge)
+{
+ if (this!=&surfaceCharge) {
+ m_position=surfaceCharge.m_position;
+ m_charge=surfaceCharge.m_charge;
+ } else {}
+ return *this;
+}
+
diff --git a/Tracker/TrackerSimEvent/TrackerSimEvent/FaserSiHit.h b/Tracker/TrackerSimEvent/TrackerSimEvent/FaserSiHit.h
index 8669b8809f5a98c6c0e579dd4703f72448bcf539..b17acf1db85d63a4fc006c450a34164a08be833d 100644
--- a/Tracker/TrackerSimEvent/TrackerSimEvent/FaserSiHit.h
+++ b/Tracker/TrackerSimEvent/TrackerSimEvent/FaserSiHit.h
@@ -141,8 +141,8 @@ private:
HepMcParticleLink m_partLink;
unsigned int m_ID;
public:
- // enum
- // { xDep = 2, xPhi = 0, xEta = 1};
+ enum
+ { xDep = 2, xPhi = 0, xEta = 1};
};
diff --git a/Tracker/TrackerSimEvent/TrackerSimEvent/SiCharge.h b/Tracker/TrackerSimEvent/TrackerSimEvent/SiCharge.h
new file mode 100644
index 0000000000000000000000000000000000000000..81365d2fc68908e1b002e4f1dbcb0febe080c68e
--- /dev/null
+++ b/Tracker/TrackerSimEvent/TrackerSimEvent/SiCharge.h
@@ -0,0 +1,135 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SiCharge.h
+// Header file for class SiCharge
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Class which contains the bare charge deposited by a single process
+///////////////////////////////////////////////////////////////////
+// Version 1.4 08/06/2001 David Calvet
+///////////////////////////////////////////////////////////////////
+
+#ifndef SITRACKEREVENT_SICHARGE_H
+#define SITRACKEREVENT_SICHARGE_H
+
+#include <iostream>
+
+// Member classes
+#include "GeneratorObjects/HepMcParticleLink.h"
+
+class SiCharge {
+
+public:
+ enum Process {no,track,diodeX_Talk,cellX_Talk,noise,extraNoise,cut_track};
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+public:
+
+ // Copy constructor:
+ SiCharge(const SiCharge &charge);
+
+ // Constructor with parameters:
+ // deposited charge
+ // time of deposition
+ // type of process which produced this charge
+ // Particle Link to the particle generating the Charge
+ SiCharge(const double& charge,const double& time,
+ const Process& processType,const HepMcParticleLink& PL);
+
+ SiCharge(const double& charge,const double& time,
+ const Process& processType);
+
+ // Destructor:
+ ~SiCharge();
+
+ // Assignment operator:
+ SiCharge &operator=(const SiCharge &charge);
+
+ ///////////////////////////////////////////////////////////////////
+ // Const methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // deposited charge:
+ double charge() const;
+
+ // time of deposition:
+ double time() const;
+
+ // type of process which produced this charge:
+ Process processType() const;
+
+ // Barcode of the particle generating the charge:
+ int trackBarcode() const;
+
+ // Particle Link of the particle generating the charge
+ const HepMcParticleLink& particleLink() const;
+
+ ///////////////////////////////////////////////////////////////////
+ // Non-const methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // add another charge, if the process and track are the same:
+ // returns true if the charge was added
+ bool add(const SiCharge &charge);
+
+ ///////////////////////////////////////////////////////////////////
+ // Private methods:
+ ///////////////////////////////////////////////////////////////////
+private:
+
+ SiCharge();
+
+ ///////////////////////////////////////////////////////////////////
+ // Private data:
+ ///////////////////////////////////////////////////////////////////
+private:
+ double m_charge; // deposited charge
+ double m_time; // time of deposition
+ Process m_processType; // type of process which produced this charge
+ // int m_trackNumber; // track number in case of track process
+ HepMcParticleLink m_partLink; //Replace the track number with a PL
+};
+
+///////////////////////////////////////////////////////////////////
+// Inline methods:
+///////////////////////////////////////////////////////////////////
+inline SiCharge::~SiCharge()
+{}
+
+inline double SiCharge::charge() const
+{
+ return m_charge;
+}
+
+inline double SiCharge::time() const
+{
+ return m_time;
+}
+
+inline SiCharge::Process SiCharge::processType() const
+{
+ return m_processType;
+}
+
+inline int SiCharge::trackBarcode() const
+{
+ return m_partLink.barcode();
+}
+
+inline const HepMcParticleLink& SiCharge::particleLink() const
+{
+ return m_partLink;
+}
+
+///////////////////////////////////////////////////////////////////
+// Input/Output stream functions:
+///////////////////////////////////////////////////////////////////
+std::ostream &operator<<(std::ostream &out,const SiCharge &charge);
+
+#endif // SITRACKEREVENT_SICHARGE_H
diff --git a/Tracker/TrackerSimEvent/TrackerSimEvent/SiTotalCharge.h b/Tracker/TrackerSimEvent/TrackerSimEvent/SiTotalCharge.h
new file mode 100644
index 0000000000000000000000000000000000000000..dec3f37333cea86d234560835b884200e1655806
--- /dev/null
+++ b/Tracker/TrackerSimEvent/TrackerSimEvent/SiTotalCharge.h
@@ -0,0 +1,176 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SiTotalCharge.h
+// Header file for class SiTotalCharge
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Class which contains the sum and the composition of all bare charges
+// corresponding to a single element
+///////////////////////////////////////////////////////////////////
+// Version 1.5 08/06/2001 David Calvet
+///////////////////////////////////////////////////////////////////
+
+#ifndef SITRACKEREVENT_SITOTALCHARGE_H
+#define SITRACKEREVENT_SITOTALCHARGE_H
+
+#include <list>
+#include "TrackerSimEvent/SiCharge.h"
+#include "AthAllocators/ArenaSharedHeapSTLAllocator.h"
+
+class SiTotalCharge {
+
+ ///////////////////////////////////////////////////////////////////
+ // Public methods:
+ ///////////////////////////////////////////////////////////////////
+public:
+
+ typedef SG::ArenaSharedHeapSTLAllocator<SiCharge> alloc_t;
+ typedef std::list<SiCharge, alloc_t> list_t;
+
+ // Implicit constructor:
+ SiTotalCharge(const alloc_t& alloc);
+
+ // Copy constructor:
+ SiTotalCharge(const SiTotalCharge &totalCharge);
+
+ // Destructor:
+ ~SiTotalCharge();
+
+ // Assignment operator:
+ SiTotalCharge &operator=(const SiTotalCharge &totalCharge);
+
+ ///////////////////////////////////////////////////////////////////
+ // Const methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // total deposited charge in this element:
+ double charge() const;
+
+ // list of individual charges:
+ const list_t &chargeComposition() const;
+
+ // return true if the charge composition contains more than one charge:
+ bool complexCharge() const;
+
+ // return true if the main charge comes from a track (track/diodeX_Talk):
+ bool fromTrack() const;
+
+ // return true if the main charge is extraNoise:
+ bool extraNoise() const;
+
+ // return the barcode of the main charge:
+ int trackBarcode() const;
+
+ // returns the HepMcParticleLink of the main charge
+ const HepMcParticleLink& particleLink() const;
+
+ // return the time of the main charge
+ double time() const;
+ ///////////////////////////////////////////////////////////////////
+ // Non-const methods:
+ ///////////////////////////////////////////////////////////////////
+
+ // add another SiCharge:
+ void add(const SiCharge &charge);
+
+ // add another SiTotalCharge:
+ void add(const SiTotalCharge &totalCharge);
+
+ // remove time information of the SiCharge objects:
+ void removeTimeInformation();
+
+ // remove small SiCharge objects:
+ // minimum charge allowed to keep a SiCharge in the composition list
+ void removeSmallCharges(const double minimumCharge);
+
+ ///////////////////////////////////////////////////////////////////
+ // Private methods:
+ ///////////////////////////////////////////////////////////////////
+private:
+
+ // add another SiCharge to the charge composition (not the total value):
+ void addSiCharge(const SiCharge &charge);
+
+ // return the SiCharge corresponding to the process which deposited
+ // the biggest amount of charge (returns dummy charge if list is empty)
+ // (used to be a public method, but this caused problems if the list was empty):
+ const SiCharge& mainCharge() const;
+
+ ///////////////////////////////////////////////////////////////////
+ // Private data:
+ ///////////////////////////////////////////////////////////////////
+private:
+ double m_charge; // total deposited charge in this element
+ list_t m_chargeComposition; // list of individual charges
+ HepMcParticleLink m_emptyLink;
+};
+
+///////////////////////////////////////////////////////////////////
+// Inline methods:
+///////////////////////////////////////////////////////////////////
+inline SiTotalCharge::~SiTotalCharge()
+{
+}
+
+inline double SiTotalCharge::charge() const
+{
+ return m_charge;
+}
+
+inline const SiTotalCharge::list_t &SiTotalCharge::chargeComposition() const
+{
+ return m_chargeComposition;
+}
+
+inline bool SiTotalCharge::complexCharge() const
+{
+ return (m_chargeComposition.size()>1);
+}
+
+inline bool SiTotalCharge::extraNoise() const
+{
+ if(m_chargeComposition.empty())
+ {
+ return false;
+ }
+ return (mainCharge().processType()==SiCharge::extraNoise);
+}
+
+inline int SiTotalCharge::trackBarcode() const
+{
+ if(m_chargeComposition.empty())
+ {
+ return 0;
+ }
+ return mainCharge().trackBarcode();
+}
+
+inline double SiTotalCharge::time() const
+{
+ if(m_chargeComposition.empty())
+ {
+ return 0.0;
+ }
+ return mainCharge().time();
+}
+
+inline const HepMcParticleLink& SiTotalCharge::particleLink() const
+{
+ if(m_chargeComposition.empty())
+ {
+ return m_emptyLink;
+ }
+ return mainCharge().particleLink();
+}
+
+
+///////////////////////////////////////////////////////////////////
+// Input/Output stream functions:
+///////////////////////////////////////////////////////////////////
+std::ostream &operator<<(std::ostream &out,const SiTotalCharge &totalCharge);
+
+#endif // SITRACKEREVENT_SITOTALCHARGE_H
diff --git a/Tracker/TrackerSimEvent/src/SiCharge.cxx b/Tracker/TrackerSimEvent/src/SiCharge.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..f24f342adb044ab4f6c79708207efff3b7193f2b
--- /dev/null
+++ b/Tracker/TrackerSimEvent/src/SiCharge.cxx
@@ -0,0 +1,85 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SiCharge.cxx
+// Implementation file for class SiCharge
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Version 1.4 08/06/2001 David Calvet
+///////////////////////////////////////////////////////////////////
+
+#include "TrackerSimEvent/SiCharge.h"
+
+// Copy constructor:
+SiCharge::SiCharge(const SiCharge &charge) :
+ m_charge(charge.m_charge),
+ m_time(charge.m_time),
+ m_processType(charge.m_processType),
+ m_partLink(charge.m_partLink)
+{}
+
+// Constructor with parameters:
+SiCharge::SiCharge(const double& charge,const double& time,
+ const Process& processType,const HepMcParticleLink& PL) :
+ m_charge(charge),
+ m_time(time),
+ m_processType(processType),
+ m_partLink(PL)
+{}
+
+SiCharge::SiCharge(const double& charge,const double& time,
+ const Process& processType) :
+ m_charge(charge),
+ m_time(time),
+ m_processType(processType),
+ m_partLink(HepMcParticleLink((unsigned int) 0, 0))
+{}
+
+
+// Assignment operator:
+SiCharge &SiCharge::operator=(const SiCharge &charge)
+{
+ if (this!=&charge) {
+ m_charge=charge.m_charge;
+ m_time=charge.m_time;
+ m_processType=charge.m_processType;
+ m_partLink=charge.m_partLink;
+ } else {}
+ return *this;
+}
+
+// add another charge, if the process and track are the same:
+bool SiCharge::add(const SiCharge &charge)
+{
+ // check if the two charges are compatible
+ if (charge.m_processType!=m_processType ||
+ charge.m_partLink!=m_partLink ||
+ charge.m_time!=m_time) {
+ return false;
+ } else {
+ m_charge+=charge.m_charge;
+ return true;
+ }
+}
+
+///////////////////////////////////////////////////////////////////
+// Input/Output stream functions:
+///////////////////////////////////////////////////////////////////
+std::ostream &operator<<(std::ostream &out,const SiCharge &charge)
+{
+ out << "Charge=" << charge.charge() << " Time=" << charge.time()
+ << " Process=";
+ if (charge.processType()==SiCharge::no) out << "no";
+ else if (charge.processType()==SiCharge::track) out << "track";
+ else if (charge.processType()==SiCharge::diodeX_Talk) out << "diodeX_Talk";
+ else if (charge.processType()==SiCharge::cellX_Talk) out << "cellX_Talk";
+ else if (charge.processType()==SiCharge::noise) out << "noise";
+ else if (charge.processType()==SiCharge::extraNoise) out << "extraNoise";
+ else if (charge.processType()==SiCharge::cut_track) out << "cut_track";
+ else out << "UNKNOWN !";
+ return (out << " Barcode=" << charge.trackBarcode());
+}
+
diff --git a/Tracker/TrackerSimEvent/src/SiTotalCharge.cxx b/Tracker/TrackerSimEvent/src/SiTotalCharge.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..2d17b4d1144e899f7d098bb3657f2d27ab669b98
--- /dev/null
+++ b/Tracker/TrackerSimEvent/src/SiTotalCharge.cxx
@@ -0,0 +1,152 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// SiTotalCharge.cxx
+// Implementation file for class SiTotalCharge
+///////////////////////////////////////////////////////////////////
+// (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+// Version 1.5 08/06/2001 David Calvet
+///////////////////////////////////////////////////////////////////
+
+#include "TrackerSimEvent/SiTotalCharge.h"
+#include <iterator>
+
+// Implicit constructor:
+SiTotalCharge::SiTotalCharge(const alloc_t& alloc) :
+ m_charge(0),
+ m_chargeComposition(alloc),
+ m_emptyLink((unsigned int) 0, 0)
+{}
+
+// Copy constructor:
+SiTotalCharge::SiTotalCharge(const SiTotalCharge &totalCharge) :
+ m_charge(totalCharge.m_charge),
+ m_chargeComposition(totalCharge.m_chargeComposition),
+ m_emptyLink((unsigned int) 0, 0)
+{}
+
+// Assignment operator:
+SiTotalCharge &SiTotalCharge::operator=(const SiTotalCharge &totalCharge)
+{
+ if (this!=&totalCharge) {
+ m_charge=totalCharge.m_charge;
+ m_chargeComposition=totalCharge.m_chargeComposition;
+ m_emptyLink=totalCharge.m_emptyLink;
+ } else {}
+ return *this;
+}
+
+// give main single process charge:
+const SiCharge& SiTotalCharge::mainCharge() const
+{
+ // start with the first charge
+ list_t::const_iterator p_charge=m_chargeComposition.begin();
+
+ // main charge to be determined, initialized to the first one
+ list_t::const_iterator main_charge = p_charge;
+ ++p_charge;
+ // look for the biggest amount of charge in the remaining charges
+ for( ; p_charge!=m_chargeComposition.end() ; ++p_charge) {
+ if (p_charge->charge()>main_charge->charge()) main_charge=p_charge;
+ }
+ return *main_charge;
+}
+
+bool SiTotalCharge::fromTrack() const
+{
+ if(m_chargeComposition.size()==0)
+ {
+ return false;
+ }
+ SiCharge::Process process=mainCharge().processType();
+ return (process==SiCharge::track ||
+ process==SiCharge::cut_track ||
+ process==SiCharge::diodeX_Talk);
+}
+
+// add another charge:
+void SiTotalCharge::add(const SiCharge &charge)
+{
+ // increase the total deposited charge
+ m_charge+=charge.charge();
+
+ // add the SiCharge to the list of charges
+ addSiCharge(charge);
+}
+
+// add another total charge:
+void SiTotalCharge::add(const SiTotalCharge &totalCharge)
+{
+ // increase the total deposited charge
+ m_charge+=totalCharge.charge();
+
+ // add the new list of charges to the present list of charges
+ for(list_t::const_iterator p_charge=
+ totalCharge.chargeComposition().begin() ;
+ p_charge!=totalCharge.chargeComposition().end() ; ++p_charge) {
+ addSiCharge(*p_charge);
+ }
+}
+
+// remove time information of the SiCharge objects:
+void SiTotalCharge::removeTimeInformation()
+{
+ // save the old charge composition
+ list_t oldComposition;
+ m_chargeComposition.swap(oldComposition);
+
+ // loop on all old charges
+ for(list_t::const_iterator p_charge=oldComposition.begin() ;
+ p_charge!=oldComposition.end() ; ++p_charge) {
+ // add the old charge (without time) to the list
+ addSiCharge(SiCharge(p_charge->charge(),0,
+ p_charge->processType(),p_charge->particleLink()));
+ }
+}
+
+// remove small SiCharge objects:
+void SiTotalCharge::removeSmallCharges(const double minimumCharge)
+{
+ // loop on all charges
+ for(list_t::iterator p_charge=m_chargeComposition.begin() ;
+ p_charge!=m_chargeComposition.end() ; ) {
+ // !!! p_charge is changed in the loop !!!
+
+ // if the charge is too small remove it from list
+ if (p_charge->charge()>-minimumCharge &&
+ p_charge->charge()<minimumCharge) {
+ p_charge=m_chargeComposition.erase(p_charge);
+ } else {
+ ++p_charge;
+ }
+ }
+}
+
+// add another SiCharge to the charge composition (not the total value):
+void SiTotalCharge::addSiCharge(const SiCharge &charge)
+{
+ // try to merge the SiCharge in the list of existing charges
+ for(list_t::iterator p_charge=m_chargeComposition.begin() ;
+ p_charge!=m_chargeComposition.end() ; ++p_charge) {
+ if (p_charge->add(charge)) return;
+ }
+ // add the charge to the list if not merged in existing one
+ m_chargeComposition.push_back(charge);
+}
+
+///////////////////////////////////////////////////////////////////
+// Input/Output stream functions:
+///////////////////////////////////////////////////////////////////
+std::ostream &operator<<(std::ostream &out,const SiTotalCharge &totalCharge)
+{
+ out << "Total charge=" << totalCharge.charge()
+ << " Composition:" << std::endl;
+ copy(totalCharge.chargeComposition().begin(),
+ totalCharge.chargeComposition().end(),
+ std::ostream_iterator<SiCharge>(out,"\n"));
+ return out;
+}
+