Merge branch 'b-detrapping' into 'v2.4-stable'

\end{minted}

\section{Trapping of Charge Carriers}

\section{Trapping and Detrapping of Charge Carriers}

\label{sec:trapping}

\apsq provides the possibility to simulate the trapping of charge carriers as a consequence of radiation induced lattice defects.

\apsq provides the possibility to simulate the trapping and detrapping of charge carriers as a consequence of radiation induced lattice defects.

Several models exist, that quantify the effective lifetime of electrons and holes, respectively, as a function of the fluence and, partially, the temperature.

The fluence needs to be provided to the corresponding propagation module, and is always interpreted as 1-MeV neutron equivalent fluence~\cite{niel}.

The decision on whether a charge carrier has been trapped during a step during the propagation process is calculated similarly to the recombination precesses, described in \ref{sec:recombination}.

The decision on whether a charge carrier has been trapped during a step during the propagation process is calculated similarly to the recombination processes, described in \ref{sec:recombination}.

It should be noted that the trapping of charge carriers is only one of several effects induced by radiation damage.

In \apsq, these effects are treated independently, i.e. defining the fluence for a propagation module will not affect any other process than trapping.

Until now, no models for de-trapping of charge carriers have been implemented.

In addition, for most modules, the parameters have been extracted under certain annealing conditions.

have been extracted under certain annealing conditions.

A dependency on annealing conditions has not been implemented here.

Please refer to the corresponding reference publications for further details.

The trapping probability is calculated as an exponential decay as a function of the simulation timestep as

\begin{equation*}

p_{e, h} = \left(1 - \exp^{1 \frac{\delta t}{\tau_{e, h}}}\right)

\end{equation*}

where $\delta t$ is the simulation timestep and $\tau{e,h}$ the effective lifetime of electrons and holes, respectively.

At the same time, a total time spent in the trap is calculated if a detrapping model is selected. Here, the time until the

charge carrier is de-trapped is calculated as

\begin{equation*}

\delta t = - \tau_{e.h} \ln{1-p}

\end{equation*}

where $p$ is a probability randomly chosen from a uniform distribution between 0 and 1.

\subsection{Trapping Models}

The following models for trapping of charge carriers can be selected.

\subsection{Ljubljana}

\paragraph{Ljubljana}

\label{sec:trap:ljubljana}

In the Ljubljana (sometimes referred to as \textit{Kramberger}) model~\cite{kramberger}, the trapping time follows the relation

This model can be selected in the configuration file via the parameter \parameter{trapping_model = "ljubljana"}.

\subsection{Dortmund}

\paragraph{Dortmund}

The Dortmund (sometimes referred to as \textit{Krasel}) model~\cite{dortmundTrapping}, describes the effective trapping times as

This model can be selected in the configuration file via the parameter \parameter{trapping_model = "dortmund"}.

\subsection{CMS Tracker}

\paragraph{CMS Tracker}

This effective trapping model has been developed by the CMS Tracker Group.

It follows the results of~\cite{CMSTrackerTrapping}, with measurements at fluences of up to $\Phi_{eq} = \SI{3e15}{n_{eq} \per \cm^2}$, at a temperature of $\SI{-20}{\celsius}$ and an irradiation with protons.

This model can be selected in the configuration file via the parameter \parameter{trapping_model = "cmstracker"}.

\subsection{Mandic}

\paragraph{Mandic}

The Mandi\'{c} model~\cite{Mandic} is an empirical model developed from measurements with high fluences ranging from  $\Phi_{eq} = \SI{5e15}{n_{eq} \per \cm^2}$ to  $\Phi_{eq} = \SI{1e17}{n_{eq} \per \cm^2}$ and describes the lifetime via

This model can be selected in the configuration file via the parameter \parameter{trapping_model = "mandic"}.

\subsection{Custom Trapping Model}

\paragraph{Constant Trapping Model}

For some situations or materials, a constant trapping probability is necessary. This can be achieved with the constant

trapping model. Here, the lifetimes are constant and set from the values provided in the configuration file with the

parameters \parameter{trapping_time_electron} and \parameter{trapping_time_hole}:

\begin{minted}[frame=single,framesep=3pt,breaklines=true,tabsize=2,linenos]{ini}

# Constant trapping times for electrons and holes:

trapping_model = "constant"

trapping_time_electron = 5ns

trapping_time_hole = 5ns

\end{minted}

This model can be selected in the configuration file via the parameter \parameter{trapping_model = "constant"}.

\paragraph{Custom Trapping Model}

Similarly to the mobility models described above, \apsq provides the possibility to use fully custom trapping models.

In order to use a custom model, the parameter \parameter{trapping_model = "custom"} needs to be set in the configuration file.

trapping_function_holes = "[0]"

trapping_parameters_holes = 7ns

\end{minted}

\subsection{Detrapping Models}

The detrapping is configured via the \parameter{detrapping_model} parameter. Currently, only \parameter{detrapping_model = "none"} and

\parameter{detrapping_model = "constant"} are supported.

The following models for trapping of charge carriers can be selected:

\paragraph{Constant Detrapping Model}

A constant detrapping probability, with the detrapping time defined separately for electrons and holes, can be implemented via the constant detrapping model.

This model requires the parameters \parameter{detrapping_time_electron} and \parameter{detrapping_time_hole} to be configured.

\begin{minted}[frame=single,framesep=3pt,breaklines=true,tabsize=2,linenos]{ini}

# Constant detrapping times for electrons and holes:

detrapping_model = "constant"

detrapping_time_electron = 10ns

detrapping_time_hole = 10ns

\end{minted}

    config_.setDefault<std::string>("mobility_model", "jacoboni");

    config_.setDefault<std::string>("recombination_model", "none");

    config_.setDefault<std::string>("trapping_model", "none");

    config_.setDefault<std::string>("detrapping_model", "none");

    config_.setDefault<bool>("output_linegraphs", false);

    config_.setDefault<bool>("output_animations", false);

                                  1);

        trapped_histo_ = CreateHistogram<TH1D>(

            "trapping_histo", "Fraction of trapped charge carriers;trapping [N / N_{total}] ;number of events", 100, 0, 1);

            "trapping_histo",

            "Fraction of trapped charge carriers at final state;trapping [N / N_{total}] ;number of events",

            100,

            0,

            1);

        recombination_time_histo_ =

            CreateHistogram<TH1D>("recombination_time_histo",

                                  0,

                                  static_cast<double>(Units::convert(integration_time_, "ns")));

        trapping_time_histo_ = CreateHistogram<TH1D>("trapping_time_histo",

                                                     "Time until trapping of charge carriers;time [ns];charge carriers",

                                                     "Local time of trapping of charge carriers;time [ns];charge carriers",

                                                     static_cast<int>(Units::convert(integration_time_, "ns") * 5),

                                                     0,

                                                     static_cast<double>(Units::convert(integration_time_, "ns")));

        detrapping_time_histo_ =

            CreateHistogram<TH1D>("detrapping_time_histo",

                                  "Time from trapping until detrapping of charge carriers;time [ns];charge carriers",

                                  static_cast<int>(Units::convert(integration_time_, "ns") * 5),

                                  0,

                                  static_cast<double>(Units::convert(integration_time_, "ns")));

    // Prepare trapping model

    trapping_ = Trapping(config_);

    // Prepare trapping model

    detrapping_ = Detrapping(config_);

}

void GenericPropagationModule::run(Event* event) {

            }

            // Propagate a single charge deposit

            auto [final_position, time, alive, trapped] = propagate(

                initial_position, deposit.getType(), deposit.getLocalTime(), event->getRandomEngine(), output_plot_points);

            auto [final_position, time, alive, trapped] = propagate(initial_position,

                                                                    deposit.getType(),

                                                                    deposit.getLocalTime(),

                                                                    event->getRandomEngine(),

                                                                    output_plot_points,

                                                                    charge_per_step);

            if(!alive) {

                LOG(DEBUG) << " Recombined " << charge_per_step << " at " << Units::display(final_position, {"mm", "um"})

                LOG(DEBUG) << " Trapped " << charge_per_step << " at " << Units::display(final_position, {"mm", "um"})

                           << " in " << Units::display(time, "ns") << " time, removing";

                trapped_charges_count += charge_per_step;

                if(output_plots_) {

                    trapping_time_histo_->Fill(static_cast<double>(Units::convert(time, "ns")), charge_per_step);

                }

                continue;

            }

                                    const CarrierType& type,

                                    const double initial_time,

                                    RandomNumberGenerator& random_generator,

                                    OutputPlotPoints& output_plot_points) const {

                                    OutputPlotPoints& output_plot_points,

                                    const unsigned int charge) const {

    // Create a runge kutta solver using the electric field as step function

    Eigen::Vector3d position(pos.x(), pos.y(), pos.z());

        return diffusion;

    };

    // Survival probability of this charge carrier package, evaluated at every step

    allpix::uniform_real_distribution<double> survival(0, 1);

    // Survival or detrap probability of this charge carrier package, evaluated at every step

    allpix::uniform_real_distribution<double> uniform_distribution(0, 1);

    // Define lambda functions to compute the charge carrier velocity with or without magnetic field

    std::function<Eigen::Vector3d(double, const Eigen::Vector3d&)> carrier_velocity_noB =

        // Check if charge carrier is still alive:

        is_alive = !recombination_(type,

                                   detector_->getDopingConcentration(static_cast<ROOT::Math::XYZPoint>(position)),

                                   survival(random_generator),

                                   uniform_distribution(random_generator),

                                   timestep);

        // Check if the charge carrier has been trapped:

        auto [trapped, traptime] = trapping_(type, survival(random_generator), timestep, std::sqrt(efield.Mag2()));

        auto trapped = trapping_(type, uniform_distribution(random_generator), timestep, std::sqrt(efield.Mag2()));

        if(trapped) {

            if((initial_time + runge_kutta.getTime() + traptime) < integration_time_) {

            if(output_plots_) {

                trapping_time_histo_->Fill(static_cast<double>(Units::convert(runge_kutta.getTime(), "ns")), charge);

            }

            auto detrap_time = detrapping_(type, uniform_distribution(random_generator), std::sqrt(efield.Mag2()));

            if((initial_time + runge_kutta.getTime() + detrap_time) < integration_time_) {

                LOG(DEBUG) << "De-trapping charge carrier after " << Units::display(detrap_time, {"ns", "us"});

                // De-trap and advance in time if still below integration time

                LOG(DEBUG) << "De-trapping charge carrier after " << Units::display(traptime, {"ns", "us"});

                runge_kutta.advanceTime(traptime);

                runge_kutta.advanceTime(detrap_time);

                if(output_plots_) {

                    detrapping_time_histo_->Fill(static_cast<double>(Units::convert(detrap_time, "ns")), charge);

                }

            } else {

                // Mark as trapped otherwise

                is_trapped = true;

        trapped_histo_->Write();

        recombination_time_histo_->Write();

        trapping_time_histo_->Write();

        detrapping_time_histo_->Write();

    }

    long double average_time = static_cast<long double>(total_time_picoseconds_) / 1e3 /

#include "objects/DepositedCharge.hpp"

#include "objects/PropagatedCharge.hpp"

#include "physics/Detrapping.hpp"

#include "physics/Mobility.hpp"

#include "physics/Recombination.hpp"

#include "physics/Trapping.hpp"

         * @param initial_time Initial time passed before propagation starts in local time coordinates

         * @param random_generator Reference to the random number engine to be used

         * @param output_plot_points Reference to vector to hold points for line graph output plots

         * @param charge Total charge of the observed charge carrier set

         * @return Tuple with the point where the deposit ended after propagation, the time the propagation took and a flag

         * whether it has recombined or was trapped

         */

                                                                       const CarrierType& type,

                                                                       const double initial_time,

                                                                       RandomNumberGenerator& random_generator,

                                                                       OutputPlotPoints& output_plot_points) const;

                                                                       OutputPlotPoints& output_plot_points,

                                                                       const unsigned int charge) const;

        // Local copies of configuration parameters to avoid costly lookup:

        double temperature_{}, timestep_min_{}, timestep_max_{}, timestep_start_{}, integration_time_{},

        Mobility mobility_;

        Recombination recombination_;

        Trapping trapping_;

        Detrapping detrapping_;

        // Precalculated value for Boltzmann constant:

        double boltzmann_kT_;

        Histogram<TH1D> trapped_histo_;

        Histogram<TH1D> recombination_time_histo_;

        Histogram<TH1D> trapping_time_histo_;

        Histogram<TH1D> detrapping_time_histo_;

        std::mutex stats_mutex_;

    };

The trapping probability is calculated at each step of the propagation by drawing a random number from an uniform distribution with $`0 \leq r \leq 1`$ and comparing it to the expression $`1 - e^{-dt/\tau_{eff}}`$, where $`dt`$ is the time step of the last charge carrier movement and $`\tau_{eff}`$ the effective trapping time constant.

A list of available models can be found in the user manual.

Detrapping of charge carriers can be enabled by setting a detrapping model via the parameter `detrapping_model`.

The default value is `none`, corresponding to no charge carrier detrapping being simulated.

A list of available models can be found in the user manual.

The propagation module also produces a variety of output plots. These include a 3D line plot of the path of all separately propagated charge carrier sets from their point of deposition to the end of their drift, with nearby paths having different colors. In this coloring scheme, electrons are marked in blue colors, while holes are presented in different shades of orange.

In addition, a 3D GIF animation for the drift of all individual sets of charges (with the size of the point proportional to the number of charges in the set) can be produced. Finally, the module produces 2D contour animations in all the planes normal to the X, Y and Z axis, showing the concentration flow in the sensor.

It should be noted that generating the animations is very time-consuming and should be switched off even when investigating drift behavior.

* `recombination_model`: Charge carrier lifetime model to be used for the propagation. Defaults to `none`, a list of available models can be found in the documentation. This feature requires a doping concentration to be present for the detector.

* `trapping_model`: Model for simulating charge carrier trapping from radiation-induced damage. Defaults to `none`, a list of available models can be found in the documentation. All models require explicitly setting a fluence parameter.

* `fluence`: 1MeV-neutron equivalent fluence the sensor has been exposed to.

* `detrapping_model`: Model for simulating charge carrier detrapping from radiation-induced damage. Defaults to `none`, a list of available models can be found in the documentation.

* `charge_per_step` : Maximum number of charge carriers to propagate together. Divides the total number of deposited charge carriers at a specific point into sets of this number of charge carriers and a set with the remaining charge carriers. A value of 10 charges per step is used by default if this value is not specified.

* `max_charge_groups`: Maximum number of charge groups to propagate from a single deposit point. Temporarily increases the value of `charge_per_step` to reduce the number of propagated groups if the deposit is larger than the value `max_charge_groups`*`charge_per_step`, thus reducing the negative performance impact of unexpectedly large deposits. If it is set to 0, there is no upper limit on the number of charge groups propagated. The default value is 0 charge groups.

* `spatial_precision` : Spatial precision to aim for. The timestep of the Runge-Kutta propagation is adjusted to reach this spatial precision after calculating the uncertainty from the fifth-order error method. Defaults to 0.25nm.

trapping_function_holes = "[0]"

trapping_parameters_holes = 10ns

#PASS Trapped 1385 charges during transport

#PASS Trapped 1365 charges during transport