Allow Transient Simulations

This MR adds functionality to perform transient simulations of the sensor response with Allpix Squared by exploiting the Shockley-Ramo theorem @shockley @ramo in different ways as detailed in the following.

Total induced charge

A new transfer module named InducedTransfer is available. It can be used to calculate the total induced charge at any pixel implant after a propagation of both electrons and holes using the GenericPropagation module. With the addition of trapping or recombination (to be added in the future) this module allows to retrieve the partial contribution to the final signal generated by charge carriers which traveled a significant distance through the sensor without finally reaching the implant.

The induced charge on each of the pixel implants is defined as the difference in weighting potential between the end position x_{final} retrieved from the PropagatedCharge and the initial position x_{initial} of the charge carrier obtained from the DepositedCharge object in the history. The total induced charge is calculated by multiplying the potential difference with the charge of the carrier, viz.

Q_n^{ind} = \int_{t_{initial}}^{t_{final}} I_n^{ind} = q \left( \phi (x_{final}) - \phi(x_{initial}) \right)

The resulting induced charge is summed for all propagated charge carriers and returned as a PixelCharge object. The number of neighboring pixels taken into account can be configured using the induction_matrix parameter.

Pulse of induced charge over time

A more intricate simulation is the generation of the full current pulse at the readout electrodes. This can now be performed using the TransientPropagation module which - similar to the InducedTransfer module - calculates the difference in weighting potential for each charge carrier to obtain the induced current. This is repeated for every step of the Runge-Kutta solver, so the individual contributions of each charge carrier during its drift time are added up and result in a current pulse on each pixel implants.

The charge transport is parameterized in time and the time step each simulation step takes can be configured. For each step, the induced charge on the neighboring pixel implants is calculated via the Shockley-Ramo theorem @shockley @ramo by taking the difference in weighting potential between the current position x_1 and the previous position x_0 of the charge carrier

Q_n^{ind} = \int_{t_0}^{t_1} I_n^{ind} = q \left( \phi (x_1) - \phi(x_0) \right)

and multiplying it with the charge. The resulting pulses are stored for every pixel individually.

In order to combine the individual pulses originating from the different PropagatedCharge objects, the PulseTransfer modules should be used. It simply sums the induced current pulses from all PropagatedCharge objects to form one combined pulse per readout electrode.

Changes to the framework core

Some changes to the framework core have been necessary in order to allow these simulations. The main changes are described in the following.

The Detector class now possesses a weighting potential. Since this does not wrap at pixel boundaries like the electric field does, but also needs to be evaluated for neighboring pixels, the DetectorField class has been extended with a new method getRelativeTo() which allows to calculate the field at a position relative to a reference position. For the electric field, this reference position is always the center of the current pixel, while to might be any pixel for the weighting potential.

The objects PropagatedCharge now holds a map of pixel indices and pulses. While this breaks the paradigm of a charge carrier in the sensor not possessing information about the readout channels, this is necessary in order to be able to attribute individual pulse contributions to the charge carriers (keeping the MC history intact) and also to retain separate propagation and transfer steps. Without this addition, the TransientPropagation module would have to directly output PixelCharge objects with the corresponding pulses. While this simplifies things a bit from the code point of view, it also reduces the modularity of the framework and the number of possible module combinations.

The PixelCharge object holds one Pulse, consisting of a time-binned vector with induced charge per time. The default constructor of a PixelCharge object adds a pulse with a single bin and assigns the full charge to this bin. This way, the change is backwards-compatible with algorithms just accessing the total charge of the object, such as the DefaultDigitizer.

Additional modules and tools

The WeightingPotentialReader module allows to read pre-calculated weighting potentials from files, in analogy to what the ElectricFieldReader does for electric fields. In addition, it contains a mode where the weighting field is calculated for the given planar sensor model and position just in time. The weighting potential is given by

\phi_w/V_w = \frac{1}{2\pi}f(x, y, z) - \frac{1}{2\pi}\sum_{n = 1}^N\left[f(x, y, 2nd - z) - g_z(x, y, 2nd + z)\right]

with

f(x, y, u) = \arctan\left( \frac{x_1 y_1}{u\sqrt{x_1^2 + y_1^2 + u^2}} \right) + \arctan\left( \frac{x_2 y_2}{u\sqrt{x_2^2 + y_2^2 + u^2}} \right) - \arctan\left( \frac{x_1 y_2}{u\sqrt{x_1^2 + y_2^2 + u^2}} \right) - \arctan\left( \frac{x_2 y_1}{u\sqrt{x_2^2 + y_1^2 + u^2}} \right)

with x_{1,2} = x \pm \frac{w_x}{2} \qquad y_{1,2} = y \pm \frac{w_y}{2}. The parameters w_{x,y} indicate the size of the collection electrode (i.e. the implant), V_w is the potential of the electrode and d is the thickness of the sensor. Since this uses a numerical method with a Tailor expansion, this is considerably slower than the look-up from pre-calculated field maps.

In order to facilitate the usage of weighting potentials, a generator tool is provided in the repository's tools/ directory. It takes a detector model file as input and generates the corresponding weighting potential using all available cores of the executing machine. The resulting field in INIT or APF format can be read by the WeightingPotentialReader.

This fixes #63 (closed)
This fixes #87 (closed)