|
ATLAS Offline Software
|
Go to the documentation of this file.
15 #ifndef SCT_DIGITIZATION_SCT_DETAILEDSURFACECHARGESGENERATOR_H
16 #define SCT_DIGITIZATION_SCT_DETAILEDSURFACECHARGESGENERATOR_H
28 #include "GaudiKernel/ITHistSvc.h"
29 #include "GaudiKernel/ServiceHandle.h"
30 #include "GaudiKernel/ToolHandle.h"
33 #include "CLHEP/Units/PhysicalConstants.h"
46 class SiDetectorElement;
50 class HepRandomEngine;
108 void EField(
double x,
double y,
double& Ex,
double& Ey)
const;
109 double induced (
int istrip,
double x,
double y)
const;
110 bool electron(
double x_e,
double y_e,
double& vx_e,
double& vy_e,
double& D_e)
const;
111 bool hole(
double x_h,
double y_h,
double& vx_h,
double& vy_h,
double& D_h);
112 double mud_e(
double E)
const;
113 double mud_h(
double E)
const;
117 void holeTransport(
double& x0,
double& y0,
double* Q_m2,
double* Q_m1,
double* Q_00,
double* Q_p1,
double* Q_p2, CLHEP::HepRandomEngine * rndmEngine);
118 void electronTransport(
double& x0,
double& y0,
double* Q_m2,
double* Q_m1,
double* Q_00,
double* Q_p1,
double* Q_p2, CLHEP::HepRandomEngine * rndmEngine)
const;
125 FloatProperty
m_tSurfaceDrift{
this,
"SurfaceDriftTime", 10.,
"max surface drift time"};
127 FloatProperty
m_tfix{
this,
"FixedTime", -999.,
"ixed time"};
128 FloatProperty
m_tsubtract{
this,
"SubtractTime", -999.,
"subtract drift time from mid gap"};
131 BooleanProperty
m_doHistoTrap{
this,
"doHistoTrap",
false,
"Allow filling of histos for charge trapping effect"};
132 BooleanProperty
m_doTrapping{
this,
"doTrapping",
false,
"Simulation of charge trapping effect"};
133 DoubleProperty
m_Fluence{
this,
"Fluence", 0.,
"Fluence for charge trapping effect"};
147 BooleanProperty
m_isOverlay{
this,
"isOverlay",
false,
"flag for overlay"};
150 ToolHandle<ISiPropertiesTool>
m_siPropertiesTool{
this,
"SiPropertiesTool",
"SCT_SiPropertiesTool",
"Tool to retrieve SCT silicon properties"};
151 ToolHandle<ISiliconConditionsTool>
m_siConditionsTool{
this,
"SiConditionsTool",
"SCT_SiliconConditionsTool",
"Tool to retrieve SCT silicon information"};
152 ToolHandle<ISiLorentzAngleTool>
m_lorentzAngleTool{
this,
"LorentzAngleTool",
"SiLorentzAngleTool/SCTLorentzAngleTool",
"Tool to retreive Lorentz angle"};
205 #endif // SCT_DETAILEDSURFACECHARGESGENERATOR_H
IntegerProperty m_numberOfCharges
number of charges
float m_distInterStrip
Inter strip distance normalized to 1.
void initTransportModel()
BooleanProperty m_isOverlay
double mud_h(double E) const
FloatProperty m_smallStepLength
max internal step along the larger G4 step
void init_mud_h(double T)
DoubleProperty m_depletionVoltage
BooleanProperty m_doHistoTrap
bool hole(double x_h, double y_h, double &vx_h, double &vy_h, double &D_h)
void initPotentialValue()
bool electron(double x_e, double y_e, double &vx_e, double &vy_e, double &D_e) const
bool m_SurfaceDriftFlag
surface drift ON/OFF
double mud_e(double E) const
IntegerProperty m_chargeDriftModel
float m_distHalfInterStrip
Half way distance inter strip.
DoubleProperty m_transportTimeMax
FloatProperty m_tsubtract
float MaxDiffusionSigma(const InDetDD::SiDetectorElement *element, const EventContext &ctx) const
max sigma diffusion
void processSiHit(const InDetDD::SiDetectorElement *element, const SiHit &phit, ISiSurfaceChargesInserter &inserter, const float eventTime, const unsigned short eventID, CLHEP::HepRandomEngine *rndmEngine, const EventContext &ctx)
void init_mud_e(double T)
double inducedCharge(int &istrip, double &x, double &y, double &t) const
static double GetPotentialValue(int ix, int iy)
SCT_DetailedSurfaceChargesGenerator(const std::string &type, const std::string &name, const IInterface *parent)
constructor
void EField(double x, double y, double &Ex, double &Ey) const
static double EyValue150(int ix, int iy)
ToolHandle< ISiliconConditionsTool > m_siConditionsTool
float m_tHalfwayDrift
Surface drift time.
double m_PotentialValue[81][115]
def TProfile(*args, **kwargs)
virtual StatusCode finalize() override
AlgTool finalize.
::StatusCode StatusCode
StatusCode definition for legacy code.
Class to take calculate Charge Transport in the SCT with a detailed charge transport model.
float SurfaceDriftTime(float ysurf) const
Calculate of the surface drift time.
double m_ExValue150[17][115]
double induced(int istrip, double x, double y) const
ServiceHandle< ITHistSvc > m_thistSvc
static double ExValue150(int ix, int iy)
ToolHandle< ISiPropertiesTool > m_siPropertiesTool
float DriftTime(float zhit, const InDetDD::SiDetectorElement *element, const EventContext &ctx) const
calculate drift time perpandicular to the surface for a charge at distance zhit from mid gap
TProfile2D * m_h_yzEfield
virtual void process(const InDetDD::SiDetectorElement *element, const TimedHitPtr< SiHit > &phit, ISiSurfaceChargesInserter &inserter, CLHEP::HepRandomEngine *rndmEngine, const EventContext &ctx) override
create a list of surface charges from a hit
virtual StatusCode initialize() override
AlgTool initialize.
void electronTransport(double &x0, double &y0, double *Q_m2, double *Q_m1, double *Q_00, double *Q_p1, double *Q_p2, CLHEP::HepRandomEngine *rndmEngine) const
DoubleProperty m_transportTimeStep
DoubleProperty m_biasVoltage
FloatProperty m_tSurfaceDrift
related to the surface drift
float DiffusionSigma(float zhit, const InDetDD::SiDetectorElement *element, const EventContext &ctx) const
calculate diffusion sigma from a gaussian dist scattered charge
double m_EyValue150[17][115]
DoubleProperty m_sensorTemperature
virtual void setFixedTime(float fixedTime) override
virtual ~SCT_DetailedSurfaceChargesGenerator()=default
Destructor.
void holeTransport(double &x0, double &y0, double *Q_m2, double *Q_m1, double *Q_00, double *Q_p1, double *Q_p2, CLHEP::HepRandomEngine *rndmEngine)
IntegerProperty m_eFieldModel
DoubleProperty m_magneticField
BooleanProperty m_doTrapping
ToolHandle< ISiLorentzAngleTool > m_lorentzAngleTool
float MaxDriftTime(const InDetDD::SiDetectorElement *element, const EventContext &ctx) const
max drift charge equivalent to the detector thickness