SolidStateDetectors.SSD_examples
— ConstantSSD_examples::Dict{Symbol,String}
Dictionary with the paths to the example detector configuration files provided by the package.
Find the possible keys of the dictionary with keys(SSD_examples)
.
The example detector configuration files can be loaded via
path_to_config_file = SSD_examples[:InvertedCoax]
sim = Simulation(path_to_config_file)
SolidStateDetectors.ADLChargeDriftModel
— TypeADLChargeDriftModel{T <: SSDFloat, M <: AbstractDriftMaterial, N, TM <: AbstractTemperatureModel{T}} <: AbstractChargeDriftModel{T}
Charge drift model for electrons and holes based on the AGATA Detector Library. Find a detailed description of the calculations in ADL Charge Drift Model.
Fields
electrons::CarrierParameters{T}
: Parameters to describe the electron drift along the <100> and <111> axes.holes::CarrierParameters{T}
: Parameters to describe the hole drift along the <100> and <111> axes.crystal_orientation::SMatrix{3,3,T,9}
: Rotation matrix that transforms the global coordinate system to the crystal coordinate system given by the <100>, <010> and <001> axes of the crystal.γ::SVector{N,SMatrix{3,3,T,9}}
: Reciprocal mass tensors to theN
valleys of the conduction band.parameters::ADLParameters{T}
: Parameters needed for the calculation of the electron drift velocity.temperaturemodel::TM
: Models to scale the resulting drift velocities with respect to temperature
See also CarrierParameters
.
SolidStateDetectors.AbstractChargeDensity
— Typeabstract type AbstractChargeDensity{T <: SSDFloat} end
Struct defining the charge density inside a Passive
.
For each charge density, there should be a method get_charge_density
which returns the charge density in SI units (C/m³) at a given point pt
.
Examples
SolidStateDetectors.AbstractImpurityDensity
— Typeabstract type AbstractImpurityDensity{T <: SSDFloat} end
Struct defining an impurity density inside a Semiconductor
.
For each impurity density, there should be a method get_impurity_density
which returns the impurity density in SI units (1/m³) at a given point.
Examples
The sign of the impurity density is important. It is taken into account in the conversion to a charge density and, thus, defines where the semiconductor is n-type or p-type.
SolidStateDetectors.CarrierParameters
— Typestruct CarrierParameters{T <: SSDFloat}
Parameters needed to describe the electron or hole drift along the <100> and <111> axes.
Fields
axis100::VelocityParameters{T}
: Parameters to describe the charge drift along the <100> axis.axis111::VelocityParameters{T}
: Parameters to describe the charge drift along the <111> axis.
See also VelocityParameters
.
SolidStateDetectors.ConstantChargeDensity
— Typestruct ConstantChargeDensity{T <: SSDFloat} <: AbstractChargeDensity{T}
Charge density model that assumes a constant charge density everywhere.
Fields
ρ::T
: the constant value of the charge density.
Definition in Configuration File
A ConstantChargeDensity
is defined in the configuration file through the field charge_density
(of a passive
or surrounding
) with name: constant
and a value
for ρ
.
An example definition of a constant charge density looks like this:
charge_density:
name: constant
value: 1.0e-10 # C/m³
SolidStateDetectors.ConstantImpurityDensity
— Typestruct ConstantImpurityDensity{T <: SSDFloat} <: AbstractImpurityDensity{T}
Impurity density model that assumes a constant impurity density everywhere.
Fields
ρ::T
: the constant value of the impurity density.
Definition in Configuration File
A ConstantImpurityDensity
is defined in the configuration file through the field impurity_density
(of a semiconductor
) with name: constant
and a value
for ρ
.
An example definition of a constant impurity density looks like this:
impurity_density:
name: constant
value: 1.0e10 # 1/m³
SolidStateDetectors.Contact
— Typemutable struct Contact{T, G, MT} <: AbstractContact{T}
Contact of a SolidStateDetector
.
For the simulation of the ElectricPotential
, all contacts are fixed to a constant potential value.
Parametric types:
T
: Precision type.G
: Type ofgeometry
.MT
: Type ofmaterial
.
Fields
potential::T
: Potential (in V) to which the contact will be fixed during the calculation of theElectricPotential
.material::MT
: Material of the contact.id::Int
: Unique id that will unambiguously identify the contact.name::String
: Custom name for the contact, relevant for plotting.geometry::G
: Geometry of the contact, see Constructive Solid Geometry (CSG).
Definition in Configuration File
A Contact
is defined in the configuration file through an entry in the contacts
array of a detector. It needs id
, potential
and geometry
and can optionally be given a name
and material
.
An example definition of contacts looks like this:
contacts:
- name: "n+ contact"
id: 1
potential: 5000V
material: HPGe # optional
geometry: # ....
- name: "p+ contact"
id: 2
potential: 0
material: HPGe #optional
geometry: # ....
SolidStateDetectors.CylindricalChargeDensity
— Typestruct CylindricalChargeDensity{T <: SSDFloat} <: AbstractChargeDensity{T}
Charge density model which assumes a linear gradient in charge density in each spatial dimension of a cylindrical coordinate system.
Fields
offsets::NTuple{3,T}
: charge density values at the origin of each dimension.gradients::NTuple{3,T}
: linear slopes inr
andz
direction.
Definition in Configuration File
A CylindricalChargeDensity
is defined in the configuration file through the field charge_density
(of a passive
or surrounding
) with name: cylindrical
and optional fields r
and z
that can each contain init
for initial values at 0 and gradient
for gradients in that dimension.
An example definition of a cylindrical charge density looks like this:
charge_density:
name: cylindrical
r: # impurity profile with linear gradient in r
init: 1.0e-10 # C/m³
gradient: 1.0e-11 # C/m⁴
SolidStateDetectors.CylindricalImpurityDensity
— Typestruct CylindricalImpurityDensity{T <: SSDFloat} <: AbstractImpurityDensity{T}
Impurity density model which assumes a linear gradient in impurity density in each spatial dimension of a cylindrical coordinate system.
Fields
offsets::NTuple{3,T}
: impurity density values at the origin of each dimension.gradients::NTuple{3,T}
: linear slopes inr
andz
direction.
Definition in Configuration File
A CylindricalImpurityDensity
is defined in the configuration file through the field impurity_density
(of a passive
or surrounding
) with name: cylindrical
and optional fields r
and z
that can each contain init
for initial values at 0 and gradient
for gradients in that dimension.
An example definition of a cylindrical impurity density looks like this:
impurity_density:
name: cylindrical
r: # impurity profile with linear gradient in r
init: 1.0e10 # 1/m³
gradient: 1.0e11 # 1/m⁴
SolidStateDetectors.DeadVolume
— Typestruct DeadVolume{T, G} <: AbstractVirtualVolume{T}
Volume inside which the charge drift is set to zero.
Parametric types
T
: Precision type.G
: Type ofgeometry
.
Fields
name::String
: Name of the dead volume, relevant for plotting.geometry::G
: Geometry of the dead volume, see Constructive Solid Geometry (CSG).
Definition in Configuration File
A DeadVolume
is defined through an entry in the virtual_drift_volumes
array of a detector with model: dead
. It needs a geometry
and can optionally be given a name
.
An example definition of dead volumes looks like this:
virtual_drift_volume:
- name: Volume 1
model: dead
geometry: # ...
- name: Volume 2
model: dead
geometry: # ...
SolidStateDetectors.DielectricDistribution
— Typestruct DielectricDistribution{T, N, S, AT} <: AbstractArray{T, N}
Dielectric distribution, or distribution of the relative permittivity, $\epsilon_r$, needed to calculate the ElectricPotential
. The dielectric distribution is unitless.
Parametric types
T
: Element type ofdata
.N
: Dimension of thegrid
anddata
array.S
: Coordinate system (Cartesian
orCylindrical
).AT
: Axes type.
Fields
data::Array{T, N}
: Array containing the values of the dielectric distribution at the discrete points of thegrid
.grid::Grid{T, N, S, AT}
:Grid
defining the discrete points at which the electric potential is determined.
The data
array contains the values of the dielectric distribution at the discrete points between the points defined by the axes ticks of the extended grid of grid
.
Thus, size(data) == size(grid) .+ 1
!
SolidStateDetectors.DiscreteAxis
— Typestruct DiscreteAxis{T, BL, BR, I} <: AbstractAxis{T, BL, BR, I}
Axis with discrete ticks which is used to define a dimension of a Grid
.
Parametric types
T
: Type of ticksBL
: Boundary condition at the left endpoint.BR
: Boundary condition at the right endpoint.I
: IntervalSets.Interval (closed or open boundaries)
The boundary conditions of a DiscreteAxis
can be BL, BR ∈ {:periodic, :reflecting, :infinite, :r0, :fixed}
.
Fields
interval::I
: Interval that defines the range of the axis.ticks::Vector{T}
: Array of values that correspond to the discrete ticks of the axis.
See also Grid
.
SolidStateDetectors.DiscreteAxis
— MethodDiscreteAxis(left_endpoint::T, right_endpoint::T, BL::Symbol, BR::Symbol, L::Symbol, R::Symbol, ticks::AbstractVector{T}) where {T}
Constructor of a DiscreteAxis
.
Arguments
left_endpoint::T
: Left endpoint of the interval of theDiscreteAxis
.right_endpoint::T
: Right endpoint of the interval of theDiscreteAxis
.BL::Symbol
: Boundary condition at the left endpoint.BR::Symbol
: Boundary condition at the right endpoint.L::Symbol
: Boundary type of the left endpoint.R::Symbol
: Boundary type of the right endpoint.ticks::AbstractVector{T}
: Array of values that correspond to the discrete ticks of the axis.
The boundary conditions of a DiscreteAxis
can be BL, BR ∈ {:periodic, :reflecting, :infinite, :r0, :fixed}
.
The boundary types of a DiscreteAxis
can be L, R ∈ {:closed, :open}
.
Examples
DiscreteAxis(-2.0, 2.0, :infinite, :infinite, :closed, :closed, collect(-2:0.1:2))
SolidStateDetectors.EffectiveChargeDensity
— Typestruct EffectiveChargeDensity{T, N, S, AT} <: AbstractArray{T, N}
Effective charge density needed to calculate the ElectricPotential
. The effective charge density is the charge density (in C/m³) times the volume of the voxel of the respective grid point (in m³). Thus, the unit of the effective charge density is Coulomb (C).
Parametric types
T
: Element type ofdata
.N
: Dimension of thegrid
anddata
array.S
: Coordinate system (Cartesian
orCylindrical
).AT
: Axes type.
Fields
data::Array{T, N}
: Array containing the values of the effective charge density at the discrete points of thegrid
.grid::Grid{T, N, S, AT}
:Grid
defining the discrete points at which the electric potential is determined.
SolidStateDetectors.ElectricField
— Typestruct ElectricField{T, N, S, AT} <: AbstractArray{T, N}
Electric field of the simulation in units of volt per meter (V/m).
Parametric types
T
: Element type ofgrid.axes
.N
: Dimension of thegrid
anddata
array.S
: Coordinate system (Cartesian
orCylindrical
).AT
: Axes type.
Fields
data::Array{<:StaticArray{Tuple{N}, T}, N}
: Array containing the field vectors of the electric field at the discrete points of thegrid
.grid::Grid{T, N, S, AT}
:Grid
defining the discrete points for which the electric field is determined.
SolidStateDetectors.ElectricFieldChargeDriftModel
— Typestruct ElectricFieldChargeDriftModel{T <: SSDFloat} <: AbstractChargeDriftModel{T}
Charge drift model in which the electrons and holes drift along the electric field with a mobility of ± 1m²/Vs.
This model is the default when no charge drift model is defined in the configuration file.
SolidStateDetectors.ElectricPotential
— Typestruct ElectricPotential{T, N, S, AT} <: AbstractArray{T, N}
Electric potential of the simulation in units of volt (V).
Parametric types
T
: Element type ofdata
.N
: Dimension of thegrid
anddata
array.S
: Coordinate system (Cartesian
orCylindrical
).AT
: Axes type.
Fields
data::Array{T, N}
: Array containing the values of the electric potential at the discrete points of thegrid
.grid::Grid{T, N, S, AT}
:Grid
defining the discrete points for which the electric potential is determined.
SolidStateDetectors.Event
— Typemutable struct Event{T <: SSDFloat}
Collection struct for individual events. This (mutable) struct is meant to be used to look at individual events, not to process a huge amount of events.
Fields
locations::Union{Vector{<:AbstractCoordinatePoint{T}}, Missing}
: Vector of the positions of all hits of the event.energies::Union{Vector{T}, Missing}
: Vector of energies corresponding to the hits of the event.drift_paths::Union{Vector{EHDriftPath{T}}, Missing}
: Calculated drift paths of each hit position.waveforms::Union{Vector{<:Any}, Missing}
: Generated signals (waveforms) of the event.
SolidStateDetectors.Grid
— Typestruct Grid{T, N, S <: AbstractCoordinateSystem, AT} <: AbstractGrid{T, N}
Collection of N
axes that define the dimensions of the grid needed to calculate ElectricPotential
, ElectricField
or WeightingPotential
.
Parametric types
T
: Tick type (element type) of the axes.N
: Dimension of the grid.S
: Coordinate system (Cartesian
orCylindrical
).AT
: Axes type.
Fields
axes::AT
: Tuple of lengthN
containingDiscreteAxis
for each dimension of the grid.
See also DiscreteAxis
.
SolidStateDetectors.Grid
— MethodGrid(sim::Simulation{T, Cartesian}; kwargs...)
Grid(sim::Simulation{T, Cylindrical}; kwargs...)
Initializes a Grid
based on the objects defined in a Simulation
.
The important points of all objects are sampled and added to the ticks of the grid. The grid initialization can be tuned using a set of keyword arguments listed below.
Arguments
sim::Simulation{T, S}
:Simulation
for which the grid will be defined.
Keywords
max_tick_distance = missing
: Maximum distance between neighbouring ticks of the grid. Additional grid ticks will be added if two neighbouring ticks are too far apart.max_tick_distance
can either be aQuantity
, e.g.1u"mm"
, or a Tuple ofQuantity
, e.g.(1u"mm", 15u"°", 3u"mm")
, to set it for each axis of theGrid
separately. Note that aCartesianGrid3D
requires aTuple{LengthQuantity, LengthQuantity, LengthQuantity}
while aCylindricalGrid
requires aTuple{LengthQuantity, AngleQuantity, LengthQuantity}
. Ifmax_tick_distance
ismissing
, one fourth of the axis length is used.max_distance_ratio::Real = 5
: If the ratio between a tick and its left and right neighbour is greater thanmax_distance_ratio
, additional ticks are added between the ticks that are further apart. This prevents the ticks from being too unevenly spaced.add_points_between_important_point::Bool = true
: If set totrue
, additional points will be added in between the important points obtained from sampling the objects of the simulation. If some objects are too close together, this will ensure a noticeable gap between them in the calculation of potentials and fields.for_weighting_potential::Bool = false
: Grid will be optimized for the calculation of anElectricPotential
if set totrue
, and of aWeightingPotential
if set tofalse
.
Additional Keyword for a CylindricalGrid
full_2π::Bool = false
: Grid will be extended to2π
if set totrue
and be left as is if set tofalse
.
SolidStateDetectors.LinearChargeDensity
— Typestruct LinearChargeDensity{T <: SSDFloat} <: AbstractChargeDensity{T}
Charge density model which assumes a linear gradient in charge density in each dimension of a Cartesian coordinate system.
Fields
offsets::NTuple{3,T}
: charge density values at the origin of each dimension.gradients::NTuple{3,T}
: linear slopes inx
,y
andz
direction.
Definition in Configuration File
A LinearChargeDensity
is defined in the configuration file through the field charge_density
(of a passive
or surrounding
) with name: linear
and optional fields x
, y
and z
that can each contain init
for initial values at 0 and gradient
for gradients in that dimension.
An example definition of a linear charge density looks like this:
charge_density:
name: cylindrical
x: # impurity profile with linear gradient in x
init: 1.0e-10 # C/m³
gradient: 1.0e-11 # C/m⁴
SolidStateDetectors.LinearImpurityDensity
— Typestruct LinearImpurityDensity{T <: SSDFloat} <: AbstractImpurityDensity{T}
Impurity density model which assumes a linear gradient in impurity density in each dimension of a Cartesian coordinate system.
Fields
offsets::NTuple{3,T}
: impurity density values at the origin of each dimension.gradients::NTuple{3,T}
: linear slopes inx
,y
andz
direction.
Definition in Configuration File
A LinearImpurityDensity
is defined in the configuration file through the field impurity_density
(of a passive
or surrounding
) with name: linear
and optional fields x
, y
and z
that can each contain init
for initial values at 0 and gradient
for gradients in that dimension.
An example definition of a linear impurity density looks like this:
impurity_density:
name: cylindrical
x: # impurity profile with linear gradient in x
init: 1.0e10 # 1/m³
gradient: 1.0e11 # 1/m⁴
SolidStateDetectors.Passive
— Typemutable struct Passive{T,G,MT,CDM} <: AbstractPassive{T}
Passive object, assigned to a SolidStateDetector
.
For the calculation of the ElectricPotential
and WeightingPotential
, passives can be fixed to a constant potential. They can additionally have a charge density profile that has an influence on the ElectricPotential
.
Parametric types
T
: Precision type.G
: Type ofgeometry
.MT
: Type ofmaterial
.CDM
: Type ofcharge_density_model
.
Fields
name::String
: Custom name for the passive, relevant for plotting.id::Int
: Unique id that will unambiguously identify the passive.potential::T
: Potential (in V) to which the passive will be fixed during the calculation of the electric potential. For floating passives, thepotential
value isNaN
.temperature::T
: Temperature (in K) of the passive.material::MT
: Material of the passive.charge_density_model::CDM
: Charge density model for the points inside the passive.geometry::G
: Geometry of the passive, see Constructive Solid Geometry (CSG).
Definition in Configuration File
A Passive
is defined through an entry in the passives
array of a detector or an entry in the surroundings
array in the configuration file. It needs material
and geometry
and can optionally be given a name
, id
, potential
, temperature
and charge_density
.
An example definition of passives looks like this:
passives:
- name: Passivated Surface
material: HPGe
charge_density: # ...
geometry: # ...
- name: Cryostat
id: 3
potential: 0
temperature: 293K
material: Al
geometry: # ...
SolidStateDetectors.PointType
— Typeconst PointType = UInt8
Stores certain information about a grid point via bit-flags.
Right now, there are:
const update_bit = 0x01
const undepleted_bit = 0x02
const pn_junction_bit = 0x04
Examples
How to get information out of a PointType
variable point_type
:
point_type & update_bit == 0
-> do not update this point (for fixed points)point_type & update_bit > 0
-> do update this pointpoint_type & undepleted_bit > 0
-> this point is undepletedpoint_type & pn_junction_bit > 0
-> this point belongs to the solid state detector. So it is in the volume of the n-type or p-type material.
SolidStateDetectors.PointTypes
— Typestruct PointTypes{T, N, S, AT} <: AbstractArray{T, N}
Information about the grid points used to calculate the ElectricPotential
stored via bit-flags. Data is stored as PointType
which is an UInt8
.
Parametric types
T
: Element type ofgrid.axes
.N
: Dimension of thegrid
anddata
array.S
: Coordinate system (Cartesian
orCylindrical
).AT
: Axes type.
Fields
data::Array{PointType, N}
: Array containing the point type values at the discrete points of thegrid
.grid::Grid{T, N, S, AT}
:Grid
defining the discrete points for which the point types are determined.
See also PointType
.
SolidStateDetectors.SSDInterval
— Typestruct SSDInterval{T <: SSDFloat, L, R, BL, BR} <: IntervalSets.TypedEndpointsInterval{L,R,T}
Interval containing boundary conditions of left and right boundary as parametric type (BL
and BR
).
Parametric types
T
: Precision type.L
: Boundary type of the left endpoint.R
: Boundary type of the right endpoint.BL
: Boundary condition at the left endpoint.BR
: Boundary condition at the right endpoint.
The boundary types of an SSDInterval
can be L, R ∈ {:closed, :open}
.
The boundary conditions of an SSDInterval
can be BL, BR ∈ {:periodic, :reflecting, :infinite, :r0, :fixed}
.
Fields
left::T
: Value of the left endpoint.right::T
: Value of the right endpoint.
SolidStateDetectors.Semiconductor
— Typestruct Semiconductor{T,G,MT,CDM,IDM} <: AbstractSemiconductor{T}
Semiconductor bulk of a SolidStateDetector
.
This is the volume in which electrons and holes will drift during the signal development.
Parametric types
T
: Precision type.G
: Type ofgeometry
.MT
: Type ofmaterial
.CDM
: Type ofcharge_drift_model
.IDM
: Type ofimpurity_density_model
.
Fields
temperature::T
: Temperature (in K) of the semiconductor.material::MT
: Material of the semiconductor.impurity_density_model::IDM
: Impurity density model for the points inside the semiconductor.charge_drift_model::CDM
: Model that describes the drift of electrons and holes inside the semiconductor.geometry::G
: Geometry of the semiconductor, see Constructive Solid Geometry (CSG).
Definition in Configuration File
A Semiconductor
is defined through the semiconductor
field of a detector. It needs material
and geometry
, and can optionally be given a temperature
, impurity_density
and charge_drift_model
.
An example definition of a semiconductor looks like this:
semiconductor:
material: HPGe
temperature: 78
impurity_density: # ...
charge_drift_model: # ...
geometry: # ...
SolidStateDetectors.Simulation
— Typemutable struct Simulation{T <: SSDFloat, CS <: AbstractCoordinateSystem} <: AbstractSimulation{T}
Collection of all parts of a simulation of a SolidStateDetector
.
Parametric types
T
: Precision type.CS
: Coordinate system (Cartesian
orCylindrical
).
Fields
config_dict::Dict
: Dictionary (parsed configuration file) which initialized the simulation.input_units::NamedTuple
: Units with which theconfig_dict
should be parsed.medium::NamedTuple
: Medium of the world.detector::Union{SolidStateDetector{T}, Missing}
: TheSolidStateDetector
of the simulation.world::World{T, 3, CS}
: TheWorld
of the simulation.q_eff_imp::Union{EffectiveChargeDensity{T}, Missing}
: Effective charge resulting from the impurites in theSemiconductor
of thedetector
.q_eff_fix::Union{EffectiveChargeDensity{T}, Missing}
: Fixed charge resulting from fixed space charges inPassive
of thedetector
.ϵ_r::Union{DielectricDistribution{T}, Missing}
: TheDielectricDistribution
of the simulation.point_types::Union{PointTypes{T}, Missing}
: ThePointTypes
of the simulation.electric_potential::Union{ElectricPotential{T}, Missing}
: TheElectricPotential
of the simulation.weighting_potentials::Vector{Any}
: TheWeightingPotential
for eachContact
of thedetector
in the simulation.electric_field::Union{ElectricField{T}, Missing}
: TheElectricField
of the simulation.electron_drift_field::Union{ElectricField{T}, Missing}
: The electron drift field of the simulation.hole_drift_field::Union{ElectricField{T}, Missing}
: The hole drift field of the simulation.
SolidStateDetectors.SolidStateDetector
— Typestruct SolidStateDetector{T,SC,CT,PT,VDM} <: AbstractConfig{T}
Struct to describe all parts of a solid state detector, i.e. the Semiconductor
, a set of Contact
and (optionally) Passive
and virtual drift volumes.
The properties of the parts (charge densities, fixed potentials, relative permittivity of the materials) will be used as input to the calculation of ElectricPotential
and WeightingPotential
in the Simulation
.
Parametric types
T
: Precision type.SC
: Type of thesemiconductor
.CT
: Type of thecontacts
.PT
: Type of thepassives
.VDM
: Type of thevirtual_drift_volumes
.
Fields
name::String
: Name of the detector.semiconductor::SC
:Semiconductor
of the detector.contacts::CT
: Vector ofContact
of the detector.passives::PT
: Vector ofPassive
objects, e.g. holding structures around the detector.virtual_drift_volumes::VDM
: Vector of virtual drift volumes in which the drift can be modulated by user-defined methods formodulate_driftvector
, e.g.DeadVolume
.
See also Semiconductor
, Contact
, Passive
and DeadVolume
.
SolidStateDetectors.VelocityParameters
— Typestruct VelocityParameters{T <: SSDFloat}
Values needed to parametrize the longitudinal drift velocity of electrons or hole along a crystal axis as a function of the electric field strength.
Background information
The parameterization for the longitudinal drift velocity, $v_l$, as a function of the electric field strength, $E$, was proposed by D.M. Caughey and R.E. Thomas and later expanded by L. Mihailescu et al.:
\[v_l = \frac{\mu_0 E}{(1 + (E/E_0 )^{\beta})^{1/ \beta}} - \mu_{n} E.\]
with the four parameters, $\mu_0$, $\beta$, $E_0$ and $\mu_n$, which are different for electrons and holes and for the different crystal axes.
The parameter $\mu_n$ accounts for the Gunn effects for electrons and should be 0 for holes.
Fields
mu0::T
: Parameter $\mu_0$ in the parameterization shown above.beta::T
: Parameter $\beta$ in the parameterization shown above.E0::T
: Parameter $E_0$ in the parameterization shown above.mun::T
: Parameter $\mu_n$ in the parameterization shown above.
SolidStateDetectors.WeightingPotential
— Typestruct WeightingPotential{T, N, S, AT} <: AbstractArray{T, N}
Weighting potential for a given Contact
which is a unitless potential.
Parametric types
T
: Element type ofdata
.N
: Dimension of thegrid
anddata
array.S
: Coordinate system (Cartesian
orCylindrical
).AT
: Axes type.
Fields
data::Array{T, N}
: Array containing the values of the weighting potential at the discrete points of thegrid
.grid::Grid{T, N, S, AT}
:Grid
defining the discrete points for which the weighting potential is determined.
SolidStateDetectors.World
— Typestruct World{T <: SSDFloat, N, S} <: AbstractWorld{T, N}
Definition of the finite volume on which a Simulation
is performed.
Parametric types
T
: Precision type.N
: Dimensions of the world.S
: Coordinate system (Cartesian
orCylindrical
).
Fields
intervals::NTuple{N, SSDInterval{T}}
: A set ofSSDInterval
defining the dimensions of the world.
SolidStateDetectors.add_baseline_and_extend_tail
— Methodadd_baseline_and_extend_tail(wv::RadiationDetectorSignals.RDWaveform{T,U,TV,UV}, n_baseline_samples::Int, total_waveform_length::Int) where {T,U,TV,UV}
Adds a zero-valued baseline in front of the waveform wv
and extends (or cuts off) the waveform at the end with the last value of wv
. A waveform of length total_waveform_length
is returned.
Arguments
wv::RadiationDetectorSignals.RDWaveform{T,U,TV,UV}
: A waveform (signal over time).n_baseline_samples::Int
: Number of samples added in front of the waveform with values 0.total_waveform_length::Int
: Number of samples of the extended waveform which is returned.
Examples
add_baseline_and_extend_tail(wv, 1000, 5000)
This functions assumes that the time steps between the samples of the input waveform wv
are the same. Thus, that the input waveform is sampled with a fixed frequency.
SolidStateDetectors.add_fano_noise
— Functionadd_fano_noise(E_dep::RealQuantity, E_ionisation::RealQuantity, f_fano::Real)::RealQuantity
Adds Fano noise to an energy deposition E_dep
, assuming a detector material ionisation energy E_ionisation
and a Fano factor f_fano
.
Arguments
E_dep::RealQuantity
: Energy deposited in a semiconductor material.E_ionisation
: Energy needed to create one electron-hole-pair in the semiconductor material.f_fano
: Fano factor of the material.
Example
add_fano_noise(100u"keV", 2.95u"eV", 0.129)
Some material properties are stored in SolidStateDetectors.material_properties
and can be used here:
material = SolidStateDetectors.material_properties[:HPGe]
add_fano_noise(100u"keV", material.E_ionisation, material.f_fano)
Using values with units for E_dep
or E_ionisation
requires the package Unitful.jl.
SolidStateDetectors.apply_initial_state!
— Methodapply_initial_state!(sim::Simulation{T}, ::Type{ElectricPotential}, grid::Grid{T} = Grid(sim);
not_only_paint_contacts::Bool = true, paint_contacts::Bool = true)::Nothing where {T <: SSDFloat}
Applies the initial state for the calculation of the ElectricPotential
. It overwrites sim.electric_potential
, sim.q_eff_imp
, sim.q_eff_fix
, sim.ϵ
and sim.point_types
with the material properties and fixed potentials defined in sim.detector
.
Arguments
sim::Simulation{T}
:Simulation
for which the initial state should be applied.grid::Grid{T}
:Grid
to apply the initial state on. If nogrid
is given, a defaultGrid
is determined fromsim
.
Keywords
not_only_paint_contacts::Bool = true
: Whether to only use the painting algorithm of the surfaces ofContact
without checking if points are actually inside them. Setting it tofalse
should improve the performance but the points inside ofContact
are not fixed anymore.paint_contacts::Bool = true
: Enable or disable the painting of the surfaces of theContact
onto thegrid
.
Examples
apply_initial_state!(sim, ElectricPotential, paint_contacts = false)
SolidStateDetectors.apply_initial_state!
— Methodapply_initial_state!(sim::Simulation{T}, ::Type{WeightingPotential}, contact_id::Int, grid::Grid{T} = Grid(sim))::Nothing
Applies the initial state for the calculation of the WeightingPotential
for the [Contact
}(@ref) with the id contact_id
. It overwrites sim.weighting_potentials[contact_id]
with the fixed values on the [Contact
}(@ref).
Arguments
sim::Simulation{T}
:Simulation
for which the initial state should be applied.contact_id::Int
: Theid
of theContact
for which theWeightingPotential
is to be calculated.grid::Grid{T}
:Grid
to apply the initial state on. If nogrid
is given, a defaultGrid
is determined fromsim
.
Keywords
not_only_paint_contacts::Bool = true
: Whether to only use the painting algorithm of the surfaces ofContact
without checking if points are actually inside them. Setting it tofalse
should improve the performance but the points inside ofContact
are not fixed anymore.paint_contacts::Bool = true
: Enable or disable the painting of the surfaces of theContact
onto thegrid
.
Examples
apply_initial_state!(sim, WeightingPotential, 1) # => applies initial state for weighting potential of contact with id 1
SolidStateDetectors.calculate_capacitance
— Methodcalculate_capacitance(sim::Simulation{T}) where {T <: SSDFloat}
Calculates and returns the capacitance of a SolidStateDetector
in a given Simulation
in units of pF.
Arguments
sim::Simulation{T}
:Simulation
withsim.detector
for which the capacitance is calculated.
This method only works if sim.electric_field
has already been calculated and is not missing
.
SolidStateDetectors.calculate_drift_fields!
— Methodcalculate_drift_fields!(sim::Simulation{T}; use_nthreads::Int = Base.Threads.nthreads())
Calculates the drift fields for electrons and holes from the ElectricField
, sim.electric_field
and the charge drift model sim.detector.semiconductor.charge_drift_model
and stores them in sim.electron_drift_field
and sim.hole_drift_field
.
Arguments
sim::Simulation{T}
:Simulation
for whichsim.electric_field
has already been calculated.
Keywords
use_nthreads::Int = Base.Threads.nthreads()
: Number of threads that should be used when calculating the drift fields.
Examples
calculate_drift_fields!(sim, use_nthreads = 4)
This method only works if sim.electric_field
has already been calculated and is not missing
.
SolidStateDetectors.calculate_electric_field!
— Methodcalculate_electric_field!(sim::Simulation{T}; n_points_in_φ::Union{Missing, Int} = missing)::Nothing
Calculates the ElectricField
from the ElectricPotential
stored in sim.electric_potential
and stores it in sim.electric_field
.
Arguments
sim::Simulation{T}
:Simulation
for whichsim.electric_potential
has already been calculated.
Keywords
n_points_in_φ::Union{Missing, Int}
: For a 2DElectricPotential
(cylindrical coordinates and symmetric inφ
),sim.electric_potential
is extended ton_points_in_φ
"layers" inφ
in order to calculate a 3D [ElectricField
]. Ifn_points_in_φ
ismissing
, the default value is36
.
Examples
calculate_electric_field!(sim, n_points_in_φ = 32)
This method only works if sim.electric_potential
has already been calculated and is not missing
.
SolidStateDetectors.calculate_electric_potential!
— Methodcalculate_electric_potential!(sim::Simulation{T}; kwargs...)::Nothing
Calculates the ElectricPotential
for a given Simulation
sim
on an adaptive grid through successive over relaxation and stores it in sim.electric_potential
.
There are several keyword arguments which can be used to tune the calculation.
Arguments
sim::Simulation{T}
:Simulation
for which theElectricPotential
is calculated.
Keywords
convergence_limit::Real
:convergence_limit
times the bias voltage sets the convergence limit of the relaxation. The convergence value is the absolute maximum difference of the potential between two iterations of all grid points. Default ofconvergence_limit
is1e-7
(times bias voltage).refinement_limits
: Defines the maximum relative (to applied bias voltage) allowed differences of the potential value of neighbored grid points in each dimension for each refinement.rl::Real
-> One refinement withrl
equal in all 3 dimensions.rl::Tuple{<:Real,<:Real,<:Real}
-> One refinement withrl
set individual for each dimension.rl::Vector{<:Real}
->length(l)
refinements withrl[i]
being the limit for the i-th refinement.rl::Vector{<:Real,<:Real,<:Real}}
->length(rl)
refinements withrl[i]
being the limits for thei
-th refinement.
min_tick_distance::Tuple{<:Quantity, <:Quantity, <:Quantity}
: Tuple of the minimum allowed distance between two grid ticks for each dimension. It prevents the refinement to make the grid too fine. Default is1e-5
for linear axes and1e-5 / (0.25 * r_max)
for the polar axis in case of a cylindricalgrid
.max_tick_distance::Tuple{<:Quantity, <:Quantity, <:Quantity}
: Tuple of the maximum allowed distance between two grid ticks for each dimension used in the initialization of the grid. Default is 1/4 of size of the world of the respective dimension.max_distance_ratio::Real
: Maximum allowed ratio between the two distances in any dimension to the two neighbouring grid points. If the ratio is too large, additional ticks are generated such that the new ratios are smaller thanmax_distance_ratio
. Default is5
.grid::Grid
: Initial grid used to start the simulation. Default isGrid(sim)
.depletion_handling::Bool
: Enables the handling of undepleted regions. Default isfalse
.use_nthreads::Union{Int, Vector{Int}}
: If<:Int
,use_nthreads
defines the maximum number of threads to be used in the computation. Fewer threads might be used depending on the current grid size due to threading overhead. Default isBase.Threads.nthreads()
. If<:Vector{Int}
,use_nthreads[i]
defines the number of threads used for each grid (refinement) stage of the field simulation. The environment variableJULIA_NUM_THREADS
must be set appropriately before the Julia session was started (e.g.export JULIA_NUM_THREADS=8
in case of bash).sor_consts::Union{<:Real, NTuple{2, <:Real}}
: Two element tuple in case of cylindrical coordinates. First element contains the SOR constant forr
= 0. Second contains the constant at the outer most grid point inr
. A linear scaling is applied in between. First element should be smaller than the second one and both should be∈ [1.0, 2.0]
. Default is[1.4, 1.85]
. In case of Cartesian coordinates, only one value is taken.max_n_iterations::Int
: Set the maximum number of iterations which are performed after each grid refinement. Default is10000
. If set to-1
there will be no limit.not_only_paint_contacts::Bool = true
: Whether to only use the painting algorithm of the surfaces ofContact
without checking if points are actually inside them. Setting it tofalse
should improve the performance but the points inside ofContact
are not fixed anymore.paint_contacts::Bool = true
: Enable or disable the painting of the surfaces of theContact
onto thegrid
.verbose::Bool=true
: Boolean whether info output is produced or not.
Example
calculate_electric_potential!(sim, refinement_limits = [0.3, 0.1, 0.05], max_distance_ratio = 4, max_n_iterations = 20000)
SolidStateDetectors.calculate_stored_energy
— Methodcalculate_stored_energy(sim::Simulation{T}) where {T <: SSDFloat}
Calculates and returns the energy stored in the ElectricField
of a SolidStateDetector
in a given Simulation
in units of J.
Arguments
sim::Simulation{T}
:Simulation
withsim.detector
for which the stored energy is calculated.
This method only works if sim.electric_field
has already been calculated and is not missing
.
SolidStateDetectors.calculate_weighting_potential!
— Methodcalculate_weighting_potential!(sim::Simulation{T}, contact_id::Int; kwargs...)::Nothing
Calculates the WeightingPotential
for a Contact
with contact_id
given Simulation
sim
on an adaptive grid through successive over relaxation and stores it in sim.weighting_potentials[contact_id]
.
There are several keyword arguments which can be used to tune the calculation.
Keywords
convergence_limit::Real
:convergence_limit
sets the convergence limit of the relaxation. The convergence value is the absolute maximum difference of the potential between two iterations of all grid points. Default ofconvergence_limit
is1e-7
.refinement_limits
: Defines the maximum relative allowed differences of the potential value of neighbored grid points in each dimension for each refinement.rl::Real
-> One refinement withrl
equal in all 3 dimensions.rl::Tuple{<:Real,<:Real,<:Real}
-> One refinement withrl
set individual for each dimension.rl::Vector{<:Real}
->length(l)
refinements withrl[i]
being the limit for the i-th refinement.rl::Vector{<:Real,<:Real,<:Real}}
->length(rl)
refinements withrl[i]
being the limits for thei
-th refinement.
min_tick_distance::Tuple{<:Quantity, <:Quantity, <:Quantity}
: Tuple of the minimum allowed distance between two grid ticks for each dimension. It prevents the refinement to make the grid too fine. Default is1e-5
for linear axes and1e-5 / (0.25 * r_max)
for the polar axis in case of a cylindricalgrid
.max_tick_distance::Tuple{<:Quantity, <:Quantity, <:Quantity}
: Tuple of the maximum allowed distance between two grid ticks for each dimension used in the initialization of the grid. Default is 1/4 of size of the world of the respective dimension.max_distance_ratio::Real
: Maximum allowed ratio between the two distances in any dimension to the two neighbouring grid points. If the ratio is too large, additional ticks are generated such that the new ratios are smaller thanmax_distance_ratio
. Default is5
.grid::Grid
: Initial grid used to start the simulation. Default isGrid(sim)
.use_nthreads::Union{Int, Vector{Int}}
: If<:Int
,use_nthreads
defines the maximum number of threads to be used in the computation. Fewer threads might be used depending on the current grid size due to threading overhead. Default isBase.Threads.nthreads()
. If<:Vector{Int}
,use_nthreads[i]
defines the number of threads used for each grid (refinement) stage of the field simulation. The environment variableJULIA_NUM_THREADS
must be set appropriately before the Julia session was started (e.g.export JULIA_NUM_THREADS=8
in case of bash).sor_consts::Union{<:Real, NTuple{2, <:Real}}
: Two element tuple in case of cylindrical coordinates. First element contains the SOR constant forr
= 0. Second contains the constant at the outer most grid point inr
. A linear scaling is applied in between. First element should be smaller than the second one and both should be∈ [1.0, 2.0]
. Default is[1.4, 1.85]
. In case of Cartesian coordinates, only one value is taken.max_n_iterations::Int
: Set the maximum number of iterations which are performed after each grid refinement. Default is10000
. If set to-1
there will be no limit.not_only_paint_contacts::Bool = true
: Whether to only use the painting algorithm of the surfaces ofContact
without checking if points are actually inside them. Setting it tofalse
should improve the performance but the points inside ofContact
are not fixed anymore.paint_contacts::Bool = true
: Enable or disable the painting of the surfaces of theContact
onto thegrid
.verbose::Bool=true
: Boolean whether info output is produced or not.
Example
calculate_weighting_potential!(sim, 1, refinement_limits = [0.3, 0.1, 0.05], max_distance_ratio = 4, max_n_iterations = 20000)
SolidStateDetectors.drift_charges!
— Methoddrift_charges!(evt::Event{T}, sim::Simulation{T}; max_nsteps::Int = 1000, Δt::RealQuantity = 5u"ns", verbose::Bool = true)::Nothing where {T <: SSDFloat}
Calculates the electron and hole drift paths for the given Event
and Simulation
and stores them in evt.drift_paths
.
Arguments
evt::Event{T}
:Event
for which the charges should be drifted.sim::Simulation{T}
:Simulation
which defines the setup in which the charges inevt
should drift.
Keywords
max_nsteps::Int = 1000
: Maximum number of steps in the drift of each hit.Δt::RealQuantity = 5u"ns"
: Time step used for the drift.verbose = true
: Activate or deactivate additional info output.
Example
drift_charges!(evt, sim, Δt = 1u"ns", verbose = false)
Using values with units for Δt
requires the package Unitful.jl.
SolidStateDetectors.get_active_volume
— Methodget_active_volume(point_types::PointTypes{T}) where {T}
Returns an approximation of the active volume of the detector by summing up the cell volumes of all cells marked as depleted.
Arguments
point_types::PointTypes{T}
: Point types of aSimulation
.
Examples
get_active_volume(sim.point_types)
Only φ
-symmetries are taken into account.
SolidStateDetectors.get_charge_density
— Functionget_charge_density(cd::AbstractChargeDensity, pt::AbstractCoordinatePoint)
Returns the charge density at a given point, pt
, based on the charge density model cd
.
Arguments
cd::AbstractChargeDensity
: TheAbstractChargeDensity
defining the charge density inside aPassive
.pt::AbstractCoordinatePoint
: The point at whichcd
is to be evaluated.
The value returned by get_charge_density
is in units of C/m³ (SI units).
SolidStateDetectors.get_electron_drift_field
— Methodget_electron_drift_field(ef::Array{SVector{3, T},3}, chargedriftmodel::AbstractChargeDriftModel)::Array{SVector{3,T},3} where {T <: SSDFloat}
Applies the electron charge drift model onto the electric field vectors and returns the electron drift field.
The field vectors in ef
have to be in Cartesian coordinates.
Arguments
ef::Array{SVector{3, T},3}
: Three-dimensional array withElectricField
vectors, e.g.sim.electric_field.data
.chargedriftmodel::AbstractChargeDriftModel
: Model that describes the electron drift in theSemiconductor
.
Keywords
use_nthreads::Int = Base.Threads.nthreads()
: Number of threads that should be used when calculating the drift fields.
Example
get_electron_drift_field(sim.electric_field.data, sim.detector.semiconductor.charge_drift_model)
This method only works if sim.electric_field
has already been calculated and is not missing
.
SolidStateDetectors.get_hole_drift_field
— Methodget_hole_drift_field(ef::Array{SVector{3, T},3}, chargedriftmodel::AbstractChargeDriftModel)::Array{SVector{3,T},3} where {T <: SSDFloat}
Applies the hole charge drift model onto the electric field vectors and returns the hole drift field.
The field vectors in ef
have to be in Cartesian coordinates.
Arguments
ef::Array{SVector{3, T},3}
: Three-dimensional array withElectricField
vectors, e.g.sim.electric_field.data
.chargedriftmodel::AbstractChargeDriftModel
: Model that describes the hole drift in theSemiconductor
.
Keywords
use_nthreads::Int = Base.Threads.nthreads()
: Number of threads that should be used when calculating the drift fields.
Example
get_hole_drift_field(sim.electric_field.data, sim.detector.semiconductor.charge_drift_model)
This method only works if sim.electric_field
has already been calculated and is not missing
.
SolidStateDetectors.get_impurity_density
— Functionget_impurity_density(id::AbstractImpurityDensity, pt::AbstractCoordinatePoint)
Returns the impurity density at a given point, pt
, based on the impurity density model id
.
Arguments
id::AbstractImpurityDensity
: TheAbstractImpurityDensity
defining the impurity density inside aSemiconductor
.pt::AbstractCoordinatePoint
: The point at whichid
is to be evaluated.
The value returned by get_impurity_density
is in units of 1/m³ (SI Units).
SolidStateDetectors.get_path_to_example_config_files
— Methodget_path_to_example_config_files()::String
Returns the path to the example detector configuration files provided by the package.
See also SSD_examples
.
SolidStateDetectors.get_signals!
— Methodget_signals!(evt::Event{T}, sim::Simulation{T}; Δt::RealQuantity = 5u"ns")::Nothing where {T <: SSDFloat}
Generates the signals/waveforms from the drift paths of an Event
for each Contact
, for which a WeightingPotential
is specified in sim.weighting_potentials
.
The output is stored in evt.waveforms
.
Arguments
evt::Event{T}
:Event
for which the waveforms should be generated.sim::Simulation{T}
:Simulation
which defines the setup in which the waveforms are generated.
Keywords
Δt::RealQuantity = 5u"ns"
: Time steps with which the drift paths were calculated.
Example
get_signals!(evt, sim, Δt = 1u"ns") # if evt.drift_paths were calculated in time steps of 1ns
This method only works if evt.drift_paths
has already been calculated and is not missing
.
SolidStateDetectors.is_depleted
— Methodis_depleted(point_types::PointTypes)::Bool
Returns true
if all PointType
values of the PointTypes
of a Simulation
are marked as depleted and false
if any point in the PointTypes
is marked as undepleted.
It can be used to determine whether the SolidStateDetector
is depleted at the provided bias voltage.
Arguments
point_types::PointTypes
:PointTypes
of aSimulation
.
Examples
is_depleted(sim.point_types)
SolidStateDetectors.parse_config_file
— Methodparse_config_file(filename::AbstractString)::Dict where {T <: SSDFloat}
Reads in a configuration file and returns a parsed dictionary which holds all the information specified in the configuration file.
Find detailed information on configuration files in Configuration Files.
Arguments
filename::AbstractString
: File name of the configuration file. If the file is not in the same directory, a path to the file is required.
SolidStateDetectors.readsiggen
— Methodreadsiggen(file_path::String; T::Type)
Reads a SigGen configuration file (ending in .config
) in 'file_path' and returns a dictionary of all parameters. Non-existing parameteres are set to 0
.
Arguments
file_path::String
: File path leading to the SigGen configuration file.
Keywords
T::Type
: Type of the parameters in the output dictionary. Default isFloat64
.
SolidStateDetectors.refine!
— Methodrefine!(sim::Simulation{T}, ::Type{ElectricPotential}, max_diffs::Tuple, minimum_distances::Tuple, kwargs...)
Takes the current state of sim.electric_potential
and refines it with respect to the input arguments max_diffs
and minimum_distances
by
- extending the
grid
ofsim.electric_potential
to be a closed grid in all dimensions, - refining the axis of the grid based on
max_diffs
andminimum_distances
: Insert new ticks between two existing ticks such that the potential difference between each tick becomes smaller thanmax_diff[i]
(i
-> dimension) but that the distances between the ticks stays larger thanminimum_distances[i]
, and - creating the new data array for the refined grid and fill it by interpolation of the the initial
grid
.
Arguments
sim::Simulation{T}
:Simulation
for whichsim.electric_potential
will be refined.max_diffs::Tuple{<:Real,<:Real,<:Real}
: Maximum potential difference between two discrete ticks ofsim.electric_potential.grid
after refinement.minimum_distances::Tuple{<:Real,<:Real,<:Real}
: Minimum distance (in SI Units) between two discrete ticks ofsim.electric_potential.grid
after refinement.
Examples
SolidStateDetectors.refine!(sim, ElectricPotential, max_diffs = (100, 100, 100), minimum_distances = (0.01, 0.02, 0.01))
SolidStateDetectors.refine!
— Methodrefine!(sim::Simulation{T}, ::Type{WeightingPotential}, max_diffs::Tuple{<:Real,<:Real,<:Real}, minimum_distances::Tuple{<:Real,<:Real,<:Real})
Takes the current state of sim.weighting_potentials[contact_id]
and refines it with respect to the input arguments max_diffs
and minimum_distances
by
- extending the
grid
ofsim.weighting_potentials[contact_id]
to be a closed grid in all dimensions, - refining the axis of the grid based on
max_diffs
andminimum_distances
: Insert new ticks between two existing ticks such that the potential difference between each tick becomes smaller thanmax_diff[i]
(i
-> dimension) but that the distances between the ticks stays larger thanminimum_distances[i]
, and - creating the new data array for the refined grid and fill it by interpolation of the the initial
grid
.
Arguments
sim::Simulation{T}
:Simulation
for whichsim.weighting_potentials[contact_id]
will be refined.contact_id::Int
: Theid
of theContact
for which theWeightingPotential
is refined.max_diffs::Tuple{<:Real,<:Real,<:Real}
: Maximum potential difference between two discrete ticks ofsim.weighting_potentials[contact_id].grid
after refinement.minimum_distances::Tuple{<:Real,<:Real,<:Real}
: Minimum distance (in SI Units) between two discrete ticks ofsim.weighting_potentials[contact_id].grid
after refinement.
Examples
SolidStateDetectors.refine!(sim, WeightingPotential, 1, max_diffs = (0.01, 0.01, 0.01), minimum_distances = (0.01, 0.02, 0.01))
SolidStateDetectors.siggentodict
— Methodsiggentodict(config::Dict; units::Dict)
Converts the dictionary containing the parameters from a SigGen configuration file to a dictionary that can be understood by SolidStateDetectors.jl.
Arguments
config::Dict
: Dictionary containing SigGen parameters (output ofreadsiggen()
`).
Keywords
units::Dict
: Units used in SigGen configuration file (set to"mm"
,"deg"
,"V"
and"K"
). The dictionary needs the fields"length"
,"angle"
,"potential"
and"temperature"
.
SolidStateDetectors.simulate!
— Methodsimulate!( sim::Simulation{T}; kwargs...) where {T, S}
Performs a full chain simulation for a given Simulation
by
- calculating the
ElectricPotential
, - calculating the
ElectricField
, - calculating the electron and hole drift fields and
- calculating the
WeightingPotential
for eachContact
.
The output is stored in sim.electric_potential
, sim.electric_field
, sim.electron_drift_field
, sim.hole_drift_field
and sim.weighting_potentials
, respectively.
There are several keyword arguments which can be used to tune the simulation.
Arguments
sim::Simulation{T}
:Simulation
for which the full chain simulation should be performed.
Keywords
convergence_limit::Real
:convergence_limit
times the bias voltage sets the convergence limit of the relaxation. The convergence value is the absolute maximum difference of the potential between two iterations of all grid points. Default ofconvergence_limit
is1e-7
(times bias voltage).refinement_limits
: Defines the maximum relative (to applied bias voltage) allowed differences of the potential value of neighboured grid points in each dimension for each refinement.rl::Real
-> One refinement withrl
equal in all 3 dimensions.rl::Tuple{<:Real,<:Real,<:Real}
-> One refinement withrl
set individual for each dimension.rl::Vector{<:Real}
->length(l)
refinements withrl[i]
being the limit for the i-th refinement.rl::Vector{<:Real,<:Real,<:Real}}
->length(rl)
refinements withrl[i]
being the limits for thei
-th refinement.
min_tick_distance::Tuple{<:Quantity, <:Quantity, <:Quantity}
: Tuple of the minimum allowed distance between two grid ticks for each dimension. It prevents the refinement to make the grid too fine. Default is1e-5
for linear axes and1e-5 / (0.25 * r_max)
for the polar axis in case of a cylindricalgrid
.max_tick_distance::Tuple{<:Quantity, <:Quantity, <:Quantity}
: Tuple of the maximum allowed distance between two grid ticks for each dimension used in the initialization of the grid. Default is 1/4 of size of the world of the respective dimension.max_distance_ratio::Real
: Maximum allowed ratio between the two distances in any dimension to the two neighbouring grid points. If the ratio is too large, additional ticks are generated such that the new ratios are smaller thanmax_distance_ratio
. Default is5
.grid::Grid
: Initial grid used to start the simulation. Default isGrid(sim)
.depletion_handling::Bool
: Enables the handling of undepleted regions. Default isfalse
.use_nthreads::Int
: Number of threads to use in the computation. Default isBase.Threads.nthreads()
. The environment variableJULIA_NUM_THREADS
must be set appropriately before the Julia session was started (e.g.export JULIA_NUM_THREADS=8
in case of bash).sor_consts::Union{<:Real, NTuple{2, <:Real}}
: Two element tuple in case of cylindrical coordinates. First element contains the SOR constant forr
= 0. Second contains the constant at the outer most grid point inr
. A linear scaling is applied in between. First element should be smaller than the second one and both should be∈ [1.0, 2.0]
. Default is[1.4, 1.85]
. In case of Cartesian coordinates, only one value is taken.max_n_iterations::Int
: Set the maximum number of iterations which are performed after each grid refinement. Default is10000
. If set to-1
there will be no limit.not_only_paint_contacts::Bool = true
: Whether to only use the painting algorithm of the surfaces ofContact
without checking if points are actually inside them. Setting it tofalse
should improve the performance but the points inside ofContact
are not fixed anymore.paint_contacts::Bool = true
: Enable or disable the painting of the surfaces of theContact
onto thegrid
.verbose::Bool=true
: Boolean whether info output is produced or not.
See also calculate_electric_potential!
, calculate_electric_field!
, calculate_drift_fields!
and calculate_weighting_potential!
.
Example
simulate!(sim, refinement_limits = [0.3, 0.1, 0.05], max_distance_ratio = 4, max_n_iterations = 20000)
SolidStateDetectors.simulate!
— Methodsimulate!(evt::Event{T}, sim::Simulation{T}; kwargs...)::Nothing where {T <: SSDFloat}
Simulates the waveforms for the Event
for a given Simulation
by
- calculating the drift paths of all energy hits, at
evt.locations
and - generating the waveforms for each
Contact
, for which aWeightingPotential
is specified insim.weighting_potentials
.
The output is stored in evt.drift_paths
and evt.waveforms
.
Arguments
evt::Event{T}
:Event
for which the charges should be drifted.sim::Simulation{T}
:Simulation
which defines the setup in which the charges inevt
should drift.
Keywords
max_nsteps::Int = 1000
: Maximum number of steps in the drift of each hit.Δt::RealQuantity = 5u"ns"
: Time step used for the drift.verbose = true
: Activate or deactivate additional info output.
Example
simulate!(evt, sim, Δt = 1u"ns", verbose = false)
See also drift_charges!
and get_signals!
.
SolidStateDetectors.simulate_waveforms
— Methodsimulate_waveforms( mcevents::TypedTables.Table, sim::Simulation{T}; kwargs...)
simulate_waveforms( mcevents::TypedTables.Table, sim::Simulation{T}, output_dir::AbstractString, output_base_name::AbstractString; kwargs...)
Simulates the waveforms for all events defined in mcevents
for a given Simulation
by
- calculating the drift paths of all energy hits defined in
mcevents
based on the drift fields for electrons and holes stored insim.electron_drift_field
andsim.hole_drift_field
, - determining the signal (waveforms) for each
Contact
, for which aWeightingPotential
is specified insim.weighting_potentials
.
Arguments
mcevents::TypedTables.Table
: Table with information about events in the simulated setup.sim::Simulation{T}
:Simulation
which defines the setup in which the charges inmcevents
should drift.
If HDF5.jl is loaded, this function has additional arguments.
Additional Arguments (HDF5)
output_dir::AbstractString
: Directory where the HDF5 output file is saved.output_base_name::AbstractString
: Basename of the HDF5 output file, default is"generated_waveforms"
.
Keywords
max_nsteps::Int = 1000
: Maximum number of steps in the drift of each hit.Δt::RealQuantity = 4u"ns"
: Time step used for the drift.verbose = false
: Activate or deactivate additional info output.chunk_n_physics_events::Int = 1000
(HDF5 only): Number of events that should be saved in a single HDF5 output file.
Examples
simulate_waveforms(mcevents, sim, Δt = 1u"ns", verbose = false)
# => returns the input table `mcevents` with an additional column `waveform` in which the generated waveforms are stored
import HDF5
simulate_waveforms(mcevents, sim, "output_dir", "my_basename", Δt = 1u"ns", verbose = false)
# => simulates the charge drift and saves the output to "output_dir/my_basename_evts_xx.h5"
The drift paths are just calculated temporarily and not returned.
Using values with units for Δt
requires the package Unitful.jl.
SolidStateDetectors.ssd_read
— Functionssd_read(filename::AbstractString, ::Type{Simulation})
Reads a Simulation
from a HDF5 file with a given filename
using LegendHDF5IO.jl.
Arguments
filename::AbstractString
: Filename of the HDF5 file.
Example
using HDF5
using LegendHDF5IO
using SolidStateDetectors
sim = ssd_read("example_sim.h5", Simulation)
In order to use this method, the packages HDF5.jl and LegendHDF5IO.jl have to be loaded before loading SolidStateDetectors.jl.
See also ssd_write
.
SolidStateDetectors.ssd_write
— Functionssd_write(filename::AbstractString, sim::Simulation)
Converts a Simulation
to a NamedTuple
and writes it to a HDF5 file with a given filename
using LegendHDF5IO.jl.
Arguments
filename::AbstractString
: Filename of the HDF5 file.sim::Simulation
:Simulation
that should be written to the HDF5 file.
Example
using HDF5
using LegendHDF5IO
using SolidStateDetectors
sim = Simulation(SSD_examples[:InvertedCoax])
simulate!(sim)
ssd_write("example_sim.h5", sim)
If a file with filename
already exists, it will be overwritten by this method.
In order to use this method, the packages HDF5.jl and LegendHDF5IO.jl have to be loaded before loading SolidStateDetectors.jl.
See also ssd_read
.
SolidStateDetectors.update_till_convergence!
— Methodupdate_till_convergence!( sim::Simulation{T} ::Type{ElectricPotential}, convergence_limit::Real; kwargs...)::T
Takes the current state of sim.electric_potential
and updates it until it has converged.
There are several keyword arguments which can be used to tune the simulation.
Arguments
sim::Simulation{T}
:Simulation
for whichsim.electric_potential
will be updated.convergence_limit::Real
:convergence_limit
times the bias voltage sets the convergence limit of the relaxation. The convergence value is the absolute maximum difference of the potential between two iterations of all grid points. Default is1e-7
.
Keywords
n_iterations_between_checks::Int
: Number of iterations between checks. Default is set to500
.max_n_iterations::Int
: Set the maximum number of iterations which are performed after each grid refinement. Default is-1
. If set to-1
there will be no limit.depletion_handling::Bool
: Enables the handling of undepleted regions. Default isfalse
.use_nthreads::Int
: Number of threads to use in the computation. Default isBase.Threads.nthreads()
. The environment variableJULIA_NUM_THREADS
must be set appropriately before the Julia session was started (e.g.export JULIA_NUM_THREADS=8
in case of bash).not_only_paint_contacts::Bool = true
: Whether to only use the painting algorithm of the surfaces ofContact
without checking if points are actually inside them. Setting it tofalse
should improve the performance but the points inside ofContact
are not fixed anymore.paint_contacts::Bool = true
: Enable or disable the painting of the surfaces of theContact
onto thegrid
.sor_consts::Union{<:Real, NTuple{2, <:Real}}
: Two element tuple in case of cylindrical coordinates. First element contains the SOR constant forr
= 0. Second contains the constant at the outer most grid point inr
. A linear scaling is applied in between. First element should be smaller than the second one and both should be∈ [1.0, 2.0]
. Default is[1.4, 1.85]
. In case of Cartesian coordinates, only one value is taken.verbose::Bool=true
: Boolean whether info output is produced or not.
Example
SolidStateDetectors.update_till_convergence!(sim, ElectricPotential, 1e-6, depletion_handling = true)
SolidStateDetectors.update_till_convergence!
— Methodupdate_till_convergence!( sim::Simulation{T} ::Type{WeightingPotential}, contact_id::Int, convergence_limit::Real; kwargs...)::T
Takes the current state of sim.weighting_potentials[contact_id]
and updates it until it has converged.
There are several keyword arguments which can be used to tune the simulation.
Arguments
sim::Simulation{T}
:Simulation
for whichsim.weighting_potentials[contact_id]
will be updated.contact_id::Int
: Theid
of theContact
for which theWeightingPotential
is to be calculated.convergence_limit::Real
:convergence_limit
times the bias voltage sets the convergence limit of the relaxation. The convergence value is the absolute maximum difference of the potential between two iterations of all grid points. Default is1e-7
.
Keywords
n_iterations_between_checks::Int
: Number of iterations between checks. Default is set to500
.max_n_iterations::Int
: Set the maximum number of iterations which are performed after each grid refinement. Default is-1
. If set to-1
there will be no limit.depletion_handling::Bool
: Enables the handling of undepleted regions. Default isfalse
.use_nthreads::Int
: Number of threads to use in the computation. Default isBase.Threads.nthreads()
. The environment variableJULIA_NUM_THREADS
must be set appropriately before the Julia session was started (e.g.export JULIA_NUM_THREADS=8
in case of bash).not_only_paint_contacts::Bool = true
: Whether to only use the painting algorithm of the surfaces ofContact
without checking if points are actually inside them. Setting it tofalse
should improve the performance but the points inside ofContact
are not fixed anymore.paint_contacts::Bool = true
: Enable or disable the painting of the surfaces of theContact
onto thegrid
.sor_consts::Union{<:Real, NTuple{2, <:Real}}
: Two element tuple in case of cylindrical coordinates. First element contains the SOR constant forr
= 0. Second contains the constant at the outer most grid point inr
. A linear scaling is applied in between. First element should be smaller than the second one and both should be∈ [1.0, 2.0]
. Default is[1.4, 1.85]
. In case of Cartesian coordinates, only one value is taken.verbose::Bool=true
: Boolean whether info output is produced or not.
Example
SolidStateDetectors.update_till_convergence!(sim, WeightingPotential, 1, 1e-6, use_nthreads = 4)