Commit 89a2f0b1 authored by Heinrich Schindler's avatar Heinrich Schindler
Browse files

Change default behaviour of AvalancheMicroscopic such that it uses the B-field...

Change default behaviour of AvalancheMicroscopic such that it uses the B-field stepping algorithm whenever there is a non-zero magnetic field in the sensor.
parent 8b1ac3aa
Pipeline #3160173 passed with stage
in 5 minutes and 26 seconds
......@@ -15,7 +15,7 @@
\vspace{2cm}
\large
Version 2021.4
Version 2021.5
\vspace{2cm}
\large
......@@ -23,7 +23,7 @@
\vfill
September 2021
October 2021
}
\end{titlepage}
......@@ -463,8 +463,8 @@ void GetAvalancheSize(int& ne, int& ni);
Information about the ``history'' of each avalanche electron can be
retrieved by
\begin{lstlisting}
unsigned int GetNumberOfElectronEndpoints() const;
void GetElectronEndpoint(const unsigned int i,
size_t GetNumberOfElectronEndpoints() const;
void GetElectronEndpoint(const size_t i,
double& x0, double& y0, double& z0, double& t0, double& e0,
double& x1, double& y1, double& z1, double& t1, double& e1,
int& status);
......@@ -516,12 +516,14 @@ aval.EnableElectronEnergyHistogramming(&hEnergy);
After each collision,
the histogram is filled with the current electron energy.
The effect of magnetic fields can be included
in the stepping algorithm using the function
If the sensor has a non-zero magnetic field, \texttt{AvalancheMicroscopic}
will by default use a more complicated stepping algorithm which takes
the effect of the $B$ field on the electron trajectory into account.
In order to explicitly switch using magnetic fields on or off
one can use the function
\begin{lstlisting}
void EnableMagneticField();
void EnableMagneticField(const bool on);
\end{lstlisting}
By default, magnetic fields are not taken into account in the calculation.
Using
\begin{lstlisting}
......
......@@ -114,8 +114,11 @@ class AvalancheMicroscopic {
/// Retrieve the currently set size limit.
int GetAvalancheSizeLimit() const { return m_sizeCut; }
/// Enable magnetic field in stepping algorithm (default: off).
void EnableMagneticField(const bool on = true) { m_useBfield = on; }
/// Switch on/off using the magnetic field in the stepping algorithm.
void EnableMagneticField(const bool on = true) {
m_useBfieldAuto = false;
m_useBfield = on;
}
/// Set number of collisions to be skipped for plotting
void SetCollisionSteps(const unsigned int n) { m_nCollSkip = n; }
......@@ -138,7 +141,7 @@ class AvalancheMicroscopic {
/** Return the number of electron trajectories in the last
* simulated avalanche (including captured electrons). */
unsigned int GetNumberOfElectronEndpoints() const {
size_t GetNumberOfElectronEndpoints() const {
return m_endpointsElectrons.size();
}
/** Return the coordinates and time of start and end point of a given
......@@ -149,35 +152,32 @@ class AvalancheMicroscopic {
* \param e0,e1 initial and final energy
* \param status status code (see GarfieldConstants.hh)
*/
void GetElectronEndpoint(const unsigned int i, double& x0, double& y0,
void GetElectronEndpoint(const size_t i, double& x0, double& y0,
double& z0, double& t0, double& e0, double& x1,
double& y1, double& z1, double& t1, double& e1,
int& status) const;
void GetElectronEndpoint(const unsigned int i, double& x0, double& y0,
void GetElectronEndpoint(const size_t i, double& x0, double& y0,
double& z0, double& t0, double& e0, double& x1,
double& y1, double& z1, double& t1, double& e1,
double& dx1, double& dy1, double& dz1,
int& status) const;
unsigned int GetNumberOfElectronDriftLinePoints(
const unsigned int i = 0) const;
unsigned int GetNumberOfHoleDriftLinePoints(const unsigned int i = 0) const;
size_t GetNumberOfElectronDriftLinePoints(const size_t i = 0) const;
size_t GetNumberOfHoleDriftLinePoints(const size_t i = 0) const;
void GetElectronDriftLinePoint(double& x, double& y, double& z, double& t,
const int ip,
const unsigned int iel = 0) const;
void GetHoleDriftLinePoint(double& x, double& y, double& z, double& t,
const int ip, const unsigned int iel = 0) const;
unsigned int GetNumberOfHoleEndpoints() const {
return m_endpointsHoles.size();
}
void GetHoleEndpoint(const unsigned int i, double& x0, double& y0, double& z0,
size_t GetNumberOfHoleEndpoints() const { return m_endpointsHoles.size(); }
void GetHoleEndpoint(const size_t i, double& x0, double& y0, double& z0,
double& t0, double& e0, double& x1, double& y1,
double& z1, double& t1, double& e1, int& status) const;
unsigned int GetNumberOfPhotons() const { return m_photons.size(); }
size_t GetNumberOfPhotons() const { return m_photons.size(); }
// Status codes:
// -2: photon absorbed by gas molecule
void GetPhoton(const unsigned int i, double& e, double& x0, double& y0,
void GetPhoton(const size_t i, double& e, double& x0, double& y0,
double& z0, double& t0, double& x1, double& y1, double& z1,
double& t1, int& status) const;
......@@ -294,6 +294,7 @@ class AvalancheMicroscopic {
bool m_usePhotons = false;
bool m_useBandStructure = true;
bool m_useNullCollisionSteps = false;
bool m_useBfieldAuto = true;
bool m_useBfield = false;
// Transport cuts
......
......@@ -227,12 +227,10 @@ void AvalancheMicroscopic::SetTimeWindow(const double t0, const double t1) {
m_hasTimeWindow = true;
}
void AvalancheMicroscopic::GetElectronEndpoint(const unsigned int i, double& x0,
double& y0, double& z0,
double& t0, double& e0,
double& x1, double& y1,
double& z1, double& t1,
double& e1, int& status) const {
void AvalancheMicroscopic::GetElectronEndpoint(const size_t i,
double& x0, double& y0, double& z0, double& t0, double& e0,
double& x1, double& y1, double& z1, double& t1, double& e1,
int& status) const {
if (i >= m_endpointsElectrons.size()) {
std::cerr << m_className << "::GetElectronEndpoint: Index out of range.\n";
x0 = y0 = z0 = t0 = e0 = 0.;
......@@ -254,9 +252,9 @@ void AvalancheMicroscopic::GetElectronEndpoint(const unsigned int i, double& x0,
status = m_endpointsElectrons[i].status;
}
void AvalancheMicroscopic::GetElectronEndpoint(
const unsigned int i, double& x0, double& y0, double& z0, double& t0,
double& e0, double& x1, double& y1, double& z1, double& t1, double& e1,
void AvalancheMicroscopic::GetElectronEndpoint(const size_t i,
double& x0, double& y0, double& z0, double& t0, double& e0,
double& x1, double& y1, double& z1, double& t1, double& e1,
double& dx1, double& dy1, double& dz1, int& status) const {
if (i >= m_endpointsElectrons.size()) {
std::cerr << m_className << "::GetElectronEndpoint: Index out of range.\n";
......@@ -283,11 +281,10 @@ void AvalancheMicroscopic::GetElectronEndpoint(
status = m_endpointsElectrons[i].status;
}
void AvalancheMicroscopic::GetHoleEndpoint(const unsigned int i, double& x0,
double& y0, double& z0, double& t0,
double& e0, double& x1, double& y1,
double& z1, double& t1, double& e1,
int& status) const {
void AvalancheMicroscopic::GetHoleEndpoint(const size_t i,
double& x0, double& y0, double& z0, double& t0, double& e0,
double& x1, double& y1, double& z1, double& t1, double& e1,
int& status) const {
if (i >= m_endpointsHoles.size()) {
std::cerr << m_className << "::GetHoleEndpoint: Index out of range.\n";
x0 = y0 = z0 = t0 = e0 = 0.;
......@@ -309,11 +306,11 @@ void AvalancheMicroscopic::GetHoleEndpoint(const unsigned int i, double& x0,
status = m_endpointsHoles[i].status;
}
unsigned int AvalancheMicroscopic::GetNumberOfElectronDriftLinePoints(
const unsigned int i) const {
size_t AvalancheMicroscopic::GetNumberOfElectronDriftLinePoints(
const size_t i) const {
if (i >= m_endpointsElectrons.size()) {
std::cerr << m_className << "::GetNumberOfElectronDriftLinePoints:\n";
std::cerr << " Endpoint index (" << i << ") out of range.\n";
std::cerr << m_className << "::GetNumberOfElectronDriftLinePoints: "
<< "Index out of range.\n";
return 0;
}
......@@ -322,11 +319,11 @@ unsigned int AvalancheMicroscopic::GetNumberOfElectronDriftLinePoints(
return m_endpointsElectrons[i].driftLine.size() + 2;
}
unsigned int AvalancheMicroscopic::GetNumberOfHoleDriftLinePoints(
const unsigned int i) const {
size_t AvalancheMicroscopic::GetNumberOfHoleDriftLinePoints(
const size_t i) const {
if (i >= m_endpointsHoles.size()) {
std::cerr << m_className << "::GetNumberOfHoleDriftLinePoints:\n";
std::cerr << " Endpoint index (" << i << ") out of range.\n";
std::cerr << m_className << "::GetNumberOfHoleDriftLinePoints: "
<< "Index out of range.\n";
return 0;
}
......@@ -400,11 +397,9 @@ void AvalancheMicroscopic::GetHoleDriftLinePoint(double& x, double& y,
t = m_endpointsHoles[ih].driftLine[ip - 1].t;
}
void AvalancheMicroscopic::GetPhoton(const unsigned int i, double& e,
double& x0, double& y0, double& z0,
double& t0, double& x1, double& y1,
double& z1, double& t1,
int& status) const {
void AvalancheMicroscopic::GetPhoton(const size_t i, double& e,
double& x0, double& y0, double& z0, double& t0,
double& x1, double& y1, double& z1, double& t1, int& status) const {
if (i >= m_photons.size()) {
std::cerr << m_className << "::GetPhoton: Index out of range.\n";
return;
......@@ -508,6 +503,10 @@ bool AvalancheMicroscopic::TransportElectrons(std::vector<Electron>& stack,
return false;
}
// Do we need to consider the magnetic field?
const bool useBfield = m_useBfieldAuto ? m_sensor->HasMagneticField() :
m_useBfield;
// Loop over the initial set of electrons/holes.
for (auto& p : stack) {
// Make sure that the starting point is inside a medium.
......@@ -641,7 +640,7 @@ bool AvalancheMicroscopic::TransportElectrons(std::vector<Electron>& stack,
// Ratio of transverse electric field component and magnetic field.
double ezovb = 0.;
std::array<std::array<double, 3>, 3> rot;
if (m_useBfield) {
if (useBfield) {
m_sensor->MagneticField(x, y, z, bx, by, bz, status);
const double scale = hole ? Tesla2Internal : -Tesla2Internal;
bx *= scale;
......@@ -720,7 +719,7 @@ bool AvalancheMicroscopic::TransportElectrons(std::vector<Electron>& stack,
double a1 = 0., a2 = 0.;
// Initial velocity.
double vx = 0., vy = 0., vz = 0.;
if (m_useBfield) {
if (useBfield) {
// Calculate the velocity vector in the local frame.
const double vmag = c1 * sqrt(en);
ToLocal(rot, vmag * kx, vmag * ky, vmag * kz, vx, vy, vz);
......@@ -761,7 +760,7 @@ bool AvalancheMicroscopic::TransportElectrons(std::vector<Electron>& stack,
const double r = RndmUniformPos();
dt += -log(r) * fInv;
// Calculate the energy after the proposed step.
if (m_useBfield) {
if (useBfield) {
en1 = en + (a1 + a2 * dt) * dt;
if (omega > Small) {
cphi = cos(omega * dt);
......@@ -814,7 +813,7 @@ bool AvalancheMicroscopic::TransportElectrons(std::vector<Electron>& stack,
// and the proposed new position.
double kx1 = 0., ky1 = 0., kz1 = 0.;
double dx = 0., dy = 0., dz = 0.;
if (m_useBfield) {
if (useBfield) {
// Calculate the new velocity.
double vx1 = vx + 2. * c2 * ex * dt;
double vy1 = vy * cphi + vz * sphi + ezovb;
......@@ -954,7 +953,7 @@ bool AvalancheMicroscopic::TransportElectrons(std::vector<Electron>& stack,
t = t1;
// If switched on, get the magnetic field at the new location.
if (m_useBfield) {
if (useBfield) {
m_sensor->MagneticField(x, y, z, bx, by, bz, status);
const double scale = hole ? Tesla2Internal : -Tesla2Internal;
bx *= scale;
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment