Garfield 0.3
Toolkit for the detailed simulation of particle detectors based on ionization measurement in gases and semiconductors
Loading...
Searching...
No Matches
Garfield::MediumMagboltz Class Reference

Interface to Magboltz (version 11). More...

#include <MediumMagboltz.hh>

Inheritance diagram for Garfield::MediumMagboltz:
Garfield::MediumGas

Classes

struct  Deexcitation

Public Member Functions

 MediumMagboltz ()
 Default constructor.
 MediumMagboltz (const std::string &gas1, const double f1=1., const std::string &gas2="", const double f2=0., const std::string &gas3="", const double f3=0., const std::string &gas4="", const double f4=0., const std::string &gas5="", const double f5=0., const std::string &gas6="", const double f6=0.)
 Constructor.
virtual ~MediumMagboltz ()
 Destructor.
bool SetMaxElectronEnergy (const double e)
 Set the highest electron energy to be included in the table of scattering rates.
double GetMaxElectronEnergy () const
 Get the highest electron energy in the table of scattering rates.
bool SetMaxPhotonEnergy (const double e)
 Set the highest photon energy to be included in the table of scattering rates.
double GetMaxPhotonEnergy () const
 Get the highest photon energy in the table of scattering rates.
void EnableAnisotropicScattering (const bool on=true)
 Switch on/off anisotropic scattering (enabled by default)
void SetSplittingFunctionOpalBeaty ()
 Sample the secondary electron energy according to the Opal-Beaty model.
void SetSplittingFunctionGreenSawada ()
 Sample the secondary electron energy according to the Green-Sawada model.
void SetSplittingFunctionFlat ()
 Sample the secondary electron energy from a flat distribution.
void EnableDeexcitation ()
 Switch on (microscopic) de-excitation handling.
void DisableDeexcitation ()
 Switch off (microscopic) de-excitation handling.
void EnableRadiationTrapping ()
 Switch on discrete photoabsorption levels.
void DisableRadiationTrapping ()
 Switch off discrete photoabsorption levels.
void SetLineWidth (const double n)
 Set the number of emission line widths within which to apply absorption by discrete lines.
bool EnablePenningTransfer () override
 Switch on simulation of Penning transfers, using pre-implemented parameterisations of the transfer probability (if available).
bool EnablePenningTransfer (const double r, const double lambda) override
 Switch on simulation of Penning transfers by means of transfer probabilities, for all excitation levels in the mixture.
bool EnablePenningTransfer (const double r, const double lambda, std::string gasname) override
 Switch on simulation of Penning transfers by means of transfer probabilities, for all excitations of a given component.
void DisablePenningTransfer () override
 Switch the simulation of Penning transfers off globally.
bool DisablePenningTransfer (std::string gasname) override
 Switch the simulation of Penning transfers off for a given component.
void EnableCrossSectionOutput (const bool on=true)
 Write the gas cross-section table to a file during the initialisation.
void SetExcitationScaling (const double r, std::string gasname)
 Multiply all excitation cross-sections by a uniform scaling factor.
bool Initialise (const bool verbose=false)
 Initialise the table of scattering rates (called internally when a collision rate is requested and the gas mixture or other parameters have changed).
void PrintGas () override
 Print information about the present gas mixture and available data.
double GetElectronNullCollisionRate (const int band) override
double GetElectronCollisionRate (const double e, const int band) override
 Get the (real) collision rate [ns-1] at a given electron energy e [eV].
double GetElectronCollisionRate (const double e, const unsigned int level, const int band)
 Get the collision rate [ns-1] for a specific level.
bool ElectronCollision (const double e, int &type, int &level, double &e1, double &dx, double &dy, double &dz, std::vector< Secondary > &secondaries, int &band) override
 Sample the collision type.
void ComputeDeexcitation (int iLevel, int &fLevel, std::vector< Secondary > &secondaries)
double GetPhotonCollisionRate (const double e) override
bool PhotonCollision (const double e, int &type, int &level, double &e1, double &ctheta, std::vector< Secondary > &secondaries) override
void ResetCollisionCounters ()
 Reset the collision counters.
unsigned int GetNumberOfElectronCollisions () const
 Get the total number of electron collisions.
unsigned int GetNumberOfElectronCollisions (unsigned int &nElastic, unsigned int &nIonising, unsigned int &nAttachment, unsigned int &nInelastic, unsigned int &nExcitation, unsigned int &nSuperelastic) const
 Get the number of collisions broken down by cross-section type.
unsigned int GetNumberOfLevels ()
 Get the number of cross-section terms.
bool GetLevel (const unsigned int i, int &ngas, int &type, std::string &descr, double &e)
 Get detailed information about a given cross-section term i.
bool GetPenningTransfer (const unsigned int i, double &r, double &lambda)
 Get the Penning transfer probability and distance of a specific level.
unsigned int GetNumberOfElectronCollisions (const unsigned int level) const
 Get the number of collisions for a specific cross-section term.
unsigned int GetNumberOfPenningTransfers () const
 Get the number of Penning transfers that occured since the last reset.
unsigned int GetNumberOfPhotonCollisions () const
 Get the total number of photon collisions.
unsigned int GetNumberOfPhotonCollisions (unsigned int &nElastic, unsigned int &nIonising, unsigned int &nInelastic) const
 Get number of photon collisions by collision type.
void EnableThermalMotion (const bool on=true)
 Take the thermal motion of the gas at the selected temperature into account in the calculations done by Magboltz.
void EnableAutoEnergyLimit (const bool on=true)
 Let Magboltz determine the upper energy limit (this is the default) or use the energy limit specified using SetMaxElectronEnergy).
void RunMagboltz (const double e, const double b, const double btheta, const int ncoll, bool verbose, double &vx, double &vy, double &vz, double &wv, double &wr, double &dl, double &dt, double &alpha, double &eta, double &riontof, double &ratttof, double &lor, double &vxerr, double &vyerr, double &vzerr, double &wverr, double &wrerr, double &dlerr, double &dterr, double &alphaerr, double &etaerr, double &riontoferr, double &ratttoferr, double &lorerr, double &alphatof, std::array< double, 6 > &difftens)
 Run Magboltz for a given electric field, magnetic field and angle.
void GenerateGasTable (const int numCollisions=10, const bool verbose=true)
 Generate a new gas table (can later be saved to file) by running Magboltz for all electric fields, magnetic fields, and angles in the currently set grid.
void PlotElectronCrossSections (const unsigned int i, TPad *pad)
void PlotElectronCollisionRates (TPad *pad)
void PlotElectronInverseMeanFreePath (TPad *pad)
double CreateGPUTransferObject (MediumGPU *&med_gpu) override
 Create and initialise GPU Transfer class.
Public Member Functions inherited from Garfield::MediumGas
 MediumGas ()
 Constructor.
virtual ~MediumGas ()
 Destructor.
bool SetComposition (const std::string &gas1, const double f1=1., const std::string &gas2="", const double f2=0., const std::string &gas3="", const double f3=0., const std::string &gas4="", const double f4=0., const std::string &gas5="", const double f5=0., const std::string &gas6="", const double f6=0.)
 Set the gas mixture.
void GetComposition (std::string &gas1, double &f1, std::string &gas2, double &f2, std::string &gas3, double &f3, std::string &gas4, double &f4, std::string &gas5, double &f5, std::string &gas6, double &f6) const
 Retrieve the gas mixture.
bool LoadGasFile (const std::string &filename, const bool quiet=false)
 Read table of gas properties (transport parameters) from file.
bool WriteGasFile (const std::string &filename)
 Save the present table of gas properties (transport parameters) to a file.
bool MergeGasFile (const std::string &filename, const bool replaceOld)
 Read table of gas properties from and merge with the existing dataset.
bool GetPenningTransfer (const std::string &gasname, double &r, double &lambda)
 Retrieve the Penning transfer probability and distance for a specific component.
bool LoadIonMobility (const std::string &filename, const bool quiet=false)
 Read a table of (positive) ion mobilities vs. electric field from file.
bool LoadNegativeIonMobility (const std::string &filename, const bool quiet=false)
 Read a table of negative ion mobilities vs. electric field from file.
bool AdjustTownsendCoefficient ()
 Adjust the Townsend coefficient using the excitation and ionisation rates stored in the gas table and the Penning transfer probabilities.
size_t GetNumberOfIonisationLevels () const
 Return the number of ionisation levels in the table.
size_t GetNumberOfExcitationLevels () const
 Return the number of excitation levels in the table.
void GetIonisationLevel (const size_t level, std::string &label, double &energy) const
 Return the identifier and threshold of an ionisation level.
void GetExcitationLevel (const size_t level, std::string &label, double &energy) const
 Return the identifier and energy of an excitation level.
bool GetElectronIonisationRate (const size_t level, const size_t ie, const size_t ib, const size_t ia, double &f) const
 Get an entry in the table of ionisation rates.
bool GetElectronExcitationRate (const size_t level, const size_t ie, const size_t ib, const size_t ia, double &f) const
 Get an entry in the table of excitation rates.
bool IsGas () const override
void GetComponent (const unsigned int i, std::string &label, double &f) override
void SetAtomicNumber (const double z) override
double GetAtomicNumber () const override
void SetAtomicWeight (const double a) override
double GetAtomicWeight () const override
void SetNumberDensity (const double n) override
double GetNumberDensity () const override
void SetMassDensity (const double rho) override
double GetMassDensity () const override
void ResetTables () override
void SetExtrapolationMethodExcitationRates (const std::string &low, const std::string &high)
void SetExtrapolationMethodIonisationRates (const std::string &low, const std::string &high)
void SetInterpolationMethodExcitationRates (const unsigned int intrp)
void SetInterpolationMethodIonisationRates (const unsigned int intrp)
double ScaleElectricField (const double e) const override
double UnScaleElectricField (const double e) const override
double ScaleDiffusion (const double d) const override
double ScaleDiffusionTensor (const double d) const override
double ScaleTownsend (const double alpha) const override
double ScaleAttachment (const double eta) const override
double ScaleLorentzAngle (const double lor) const override
bool GetPhotoAbsorptionCrossSection (const double e, double &sigma, const unsigned int i) override

Static Public Member Functions

static int GetGasNumberMagboltz (const std::string &input)
Static Public Member Functions inherited from Garfield::MediumGas
static void PrintGases ()
 Print a list of all available gases.

Private Member Functions

bool Update (const bool verbose=false)
bool Mixer (const bool verbose=false)
void SetupGreenSawada ()
void GetExcitationIonisationLevels ()
void ComputeDeexcitationTable (const bool verbose)
void AddPenningDeexcitation (Deexcitation &dxc, const double rate, const double pPenning)
double RateConstantWK (const double energy, const double osc, const double pacs, const int igas1, const int igas2) const
double RateConstantHardSphere (const double r1, const double r2, const int igas1, const int igas2) const
double CalcDiscreteLineCf (const Deexcitation &dxc, double e, double cfOth) const
void ComputeDeexcitationInternal (int iLevel, int &fLevel, std::vector< Secondary > &secondaries)
bool ComputePhotonCollisionTable (const bool verbose)

Private Attributes

std::mutex m_mutex
 Mutex.
bool m_useGasMotion = false
 Simulate thermal motion of the gas or not (when running Magboltz).
bool m_autoEnergyLimit = true
 Automatic calculation of the energy limit by Magboltz or not.
double m_eMax
 Max. electron energy in the collision rate tables.
double m_eStep
 Energy spacing in the linear part of the collision rate tables.
double m_eStepInv
 Inverse energy spacing.
double m_eHigh
double m_eHighLog
double m_lnStep
bool m_useCsOutput = false
 Flag enabling/disabling output of cross-section table to file.
unsigned int m_nTerms = 0
 Number of different cross-section types in the current gas mixture.
std::array< double, m_nMaxGasesm_mgas
 Mass.
std::array< double, m_nMaxGasesm_rgas
 Recoil energy parameter.
std::array< double, m_nMaxGasesm_s2
std::array< double, Magboltz::nMaxLevelsm_wOpalBeaty
 Opal-Beaty-Peterson splitting parameter [eV].
std::array< std::array< double, 5 >, m_nMaxGasesm_parGreenSawada
 Green-Sawada splitting parameters [eV] (Γs, Γb, Ts, Ta, Tb).
std::array< bool, m_nMaxGasesm_hasGreenSawada
bool m_useOpalBeaty = true
 Sample secondary electron energies using Opal-Beaty parameterisation.
bool m_useGreenSawada = false
 Sample secondary electron energies using Green-Sawada parameterisation.
std::array< double, Magboltz::nMaxLevelsm_energyLoss
 Energy loss.
std::array< int, Magboltz::nMaxLevelsm_csType
 Cross-section type.
std::array< double, Magboltz::nMaxLevelsm_yFluorescence
 Fluorescence yield.
std::array< unsigned int, Magboltz::nMaxLevelsm_nAuger1
 Number of Auger electrons produced in a collision.
std::array< unsigned int, Magboltz::nMaxLevelsm_nAuger2
std::array< double, Magboltz::nMaxLevelsm_eAuger1
 Energy imparted to Auger electrons.
std::array< double, Magboltz::nMaxLevelsm_eAuger2
std::array< unsigned int, Magboltz::nMaxLevelsm_nFluorescence
std::array< double, Magboltz::nMaxLevelsm_eFluorescence
std::vector< std::vector< double > > m_scatPar
std::vector< std::vector< double > > m_scatCut
std::vector< std::vector< double > > m_scatParLog
std::vector< std::vector< double > > m_scatCutLog
std::array< int, Magboltz::nMaxLevelsm_scatModel
std::vector< std::string > m_description
 Level description.
std::vector< double > m_cfTot
std::vector< double > m_cfTotLog
std::vector< std::vector< double > > m_cf
std::vector< std::vector< double > > m_cfLog
bool m_useAnisotropic = true
double m_cfNull = 0.
 Null-collision frequency.
std::array< unsigned int, nCsTypesm_nCollisions
 Collision counters 0: elastic 1: ionisation 2: attachment 3: inelastic 4: excitation 5: super-elastic.
std::vector< unsigned int > m_nCollisionsDetailed
 Number of collisions for each cross-section term.
std::array< double, Magboltz::nMaxLevelsm_rPenning
 Penning transfer probability (by level)
std::array< double, Magboltz::nMaxLevelsm_lambdaPenning
 Mean distance of Penning ionisation (by level)
unsigned int m_nPenning = 0
 Number of Penning ionisations.
bool m_useDeexcitation = false
 Flag enabling/disabling detailed simulation of de-excitation process.
bool m_useRadTrap = true
 Flag enabling/disable radiation trapping (absorption of photons discrete excitation lines)
std::vector< Deexcitationm_deexcitations
std::array< int, Magboltz::nMaxLevelsm_iDeexcitation
double m_nAbsWidths = 1000.
std::array< double, m_nMaxGasesm_ionPot
 Ionisation potentials of each component.
double m_minIonPot = -1.
 Minimum ionisation potential.
std::array< double, m_nMaxGasesm_scaleExc
double m_eFinalGamma
double m_eStepGamma
unsigned int m_nPhotonTerms = 0
std::vector< double > m_cfTotGamma
std::vector< std::vector< double > > m_cfGamma
std::vector< int > csTypeGamma
std::array< unsigned int, nCsTypesGammam_nPhotonCollisions

Static Private Attributes

static constexpr int nEnergyStepsLog = 1000
static constexpr int nEnergyStepsGamma = 5000
static constexpr int nCsTypes = 7
static constexpr int nCsTypesGamma = 4
static const int DxcTypeRad
static const int DxcTypeCollIon
static const int DxcTypeCollNonIon

Additional Inherited Members

Protected Member Functions inherited from Garfield::MediumGas
bool LoadMobility (const std::string &filename, const bool quiet, const bool negative)
bool ReadHeader (std::ifstream &gasfile, int &version, std::bitset< 20 > &gasok, bool &is3d, std::vector< double > &mixture, std::vector< double > &efields, std::vector< double > &bfields, std::vector< double > &angles, std::vector< ExcLevel > &excLevels, std::vector< IonLevel > &ionLevels)
void ReadFooter (std::ifstream &gasfile, std::array< unsigned int, 13 > &extrapH, std::array< unsigned int, 13 > &extrapL, std::array< unsigned int, 13 > &interp, unsigned int &thrAlp, unsigned int &thrAtt, unsigned int &thrDis, double &ionDiffL, double &ionDiffT, double &pgas, double &tgas)
void ReadRecord3D (std::ifstream &gasfile, double &ve, double &vb, double &vx, double &dl, double &dt, double &alpha, double &alpha0, double &eta, double &mu, double &lor, double &dis, std::array< double, 6 > &dif, std::vector< double > &rexc, std::vector< double > &rion, std::bitset< 20 > gasok)
void ReadRecord1D (std::ifstream &gasfile, double &ve, double &vb, double &vx, double &wv, double &wr, double &dl, double &dt, double &alpha, double &alpha0, double &eta, double &riontof, double &ratttof, double &mu, double &lor, double &dis, std::array< double, 6 > &dif, std::vector< double > &rexc, std::vector< double > &rion, std::bitset< 20 > gasok)
void InsertE (const int ie, const int ne, const int nb, const int na)
void InsertB (const int ib, const int ne, const int nb, const int na)
void InsertA (const int ia, const int ne, const int nb, const int na)
void ZeroRowE (const int ie, const int nb, const int na)
void ZeroRowB (const int ib, const int ne, const int na)
void ZeroRowA (const int ia, const int ne, const int nb)
bool GetMixture (const std::vector< double > &mixture, const int version, std::vector< std::string > &gasnames, std::vector< double > &percentages) const
void GetGasBits (std::bitset< 20 > &gasok) const
Static Protected Member Functions inherited from Garfield::MediumGas
static bool GetGasInfo (const std::string &gasname, double &a, double &z, double &w, double &f)
static std::string GetGasName (const int gasnumber, const int version)
static std::string GetGasName (std::string input)
static int GetGasNumberGasFile (const std::string &input)
static const std::vector< std::string > GetAliases (const std::string &gas)
Protected Attributes inherited from Garfield::MediumGas
std::array< std::string, m_nMaxGasesm_gas
std::array< double, m_nMaxGasesm_fraction
std::array< double, m_nMaxGasesm_atWeight
std::array< double, m_nMaxGasesm_atNum
bool m_usePenning = false
double m_rPenningGlobal = 0.
double m_lambdaPenningGlobal = 0.
std::array< double, m_nMaxGasesm_rPenningGas
std::array< double, m_nMaxGasesm_lambdaPenningGas
double m_pressureTable
double m_temperatureTable
std::vector< std::vector< std::vector< double > > > m_eAlp0
std::vector< std::vector< std::vector< std::vector< double > > > > m_excRates
std::vector< std::vector< std::vector< std::vector< double > > > > m_ionRates
std::vector< ExcLevelm_excLevels
std::vector< IonLevelm_ionLevels
std::pair< unsigned int, unsigned int > m_extrExc = {0, 1}
std::pair< unsigned int, unsigned int > m_extrIon = {0, 1}
unsigned int m_intpExc = 2
unsigned int m_intpIon = 2
Static Protected Attributes inherited from Garfield::MediumGas
static constexpr unsigned int m_nMaxGases = 6

Detailed Description

Interface to Magboltz (version 11).

Definition at line 22 of file MediumMagboltz.hh.

Constructor & Destructor Documentation

◆ MediumMagboltz() [1/2]

Garfield::MediumMagboltz::MediumMagboltz ( )

Default constructor.

◆ MediumMagboltz() [2/2]

Garfield::MediumMagboltz::MediumMagboltz ( const std::string & gas1,
const double f1 = 1.,
const std::string & gas2 = "",
const double f2 = 0.,
const std::string & gas3 = "",
const double f3 = 0.,
const std::string & gas4 = "",
const double f4 = 0.,
const std::string & gas5 = "",
const double f5 = 0.,
const std::string & gas6 = "",
const double f6 = 0. )

Constructor.

◆ ~MediumMagboltz()

virtual Garfield::MediumMagboltz::~MediumMagboltz ( )
inlinevirtual

Destructor.

Definition at line 34 of file MediumMagboltz.hh.

34{}

Member Function Documentation

◆ AddPenningDeexcitation()

void Garfield::MediumMagboltz::AddPenningDeexcitation ( Deexcitation & dxc,
const double rate,
const double pPenning )
inlineprivate

Definition at line 449 of file MediumMagboltz.hh.

450 {
451 dxc.p.push_back(rate * (1. - pPenning));
452 dxc.type.push_back(DxcTypeCollNonIon);
453 dxc.final.push_back(-1);
454 if (pPenning > 0.) {
455 dxc.p.push_back(rate * pPenning);
456 dxc.type.push_back(DxcTypeCollIon);
457 dxc.final.push_back(-1);
458 }
459 }
static const int DxcTypeCollNonIon
static const int DxcTypeCollIon

◆ CalcDiscreteLineCf()

double Garfield::MediumMagboltz::CalcDiscreteLineCf ( const Deexcitation & dxc,
double e,
double cfOth ) const
private

◆ ComputeDeexcitation()

void Garfield::MediumMagboltz::ComputeDeexcitation ( int iLevel,
int & fLevel,
std::vector< Secondary > & secondaries )

◆ ComputeDeexcitationInternal()

void Garfield::MediumMagboltz::ComputeDeexcitationInternal ( int iLevel,
int & fLevel,
std::vector< Secondary > & secondaries )
private

◆ ComputeDeexcitationTable()

void Garfield::MediumMagboltz::ComputeDeexcitationTable ( const bool verbose)
private

◆ ComputePhotonCollisionTable()

bool Garfield::MediumMagboltz::ComputePhotonCollisionTable ( const bool verbose)
private

◆ CreateGPUTransferObject()

double Garfield::MediumMagboltz::CreateGPUTransferObject ( MediumGPU *& med_gpu)
overridevirtual

Create and initialise GPU Transfer class.

Reimplemented from Garfield::MediumGas.

◆ DisableDeexcitation()

void Garfield::MediumMagboltz::DisableDeexcitation ( )
inline

Switch off (microscopic) de-excitation handling.

Definition at line 64 of file MediumMagboltz.hh.

64{ m_useDeexcitation = false; }
bool m_useDeexcitation
Flag enabling/disabling detailed simulation of de-excitation process.

◆ DisablePenningTransfer() [1/2]

void Garfield::MediumMagboltz::DisablePenningTransfer ( )
overridevirtual

Switch the simulation of Penning transfers off globally.

Reimplemented from Garfield::MediumGas.

◆ DisablePenningTransfer() [2/2]

bool Garfield::MediumMagboltz::DisablePenningTransfer ( std::string gasname)
overridevirtual

Switch the simulation of Penning transfers off for a given component.

Reimplemented from Garfield::MediumGas.

◆ DisableRadiationTrapping()

void Garfield::MediumMagboltz::DisableRadiationTrapping ( )
inline

Switch off discrete photoabsorption levels.

Definition at line 68 of file MediumMagboltz.hh.

68{ m_useRadTrap = false; }
bool m_useRadTrap
Flag enabling/disable radiation trapping (absorption of photons discrete excitation lines)

◆ ElectronCollision()

bool Garfield::MediumMagboltz::ElectronCollision ( const double e,
int & type,
int & level,
double & e1,
double & dx,
double & dy,
double & dz,
std::vector< Secondary > & secondaries,
int & band )
override

Sample the collision type.

◆ EnableAnisotropicScattering()

void Garfield::MediumMagboltz::EnableAnisotropicScattering ( const bool on = true)
inline

Switch on/off anisotropic scattering (enabled by default)

Definition at line 49 of file MediumMagboltz.hh.

49 {
51 m_isChanged = true;
52 }

◆ EnableAutoEnergyLimit()

void Garfield::MediumMagboltz::EnableAutoEnergyLimit ( const bool on = true)
inline

Let Magboltz determine the upper energy limit (this is the default) or use the energy limit specified using SetMaxElectronEnergy).

Definition at line 164 of file MediumMagboltz.hh.

164{ m_autoEnergyLimit = on; }
bool m_autoEnergyLimit
Automatic calculation of the energy limit by Magboltz or not.

◆ EnableCrossSectionOutput()

void Garfield::MediumMagboltz::EnableCrossSectionOutput ( const bool on = true)
inline

Write the gas cross-section table to a file during the initialisation.

Definition at line 81 of file MediumMagboltz.hh.

81{ m_useCsOutput = on; }
bool m_useCsOutput
Flag enabling/disabling output of cross-section table to file.

◆ EnableDeexcitation()

void Garfield::MediumMagboltz::EnableDeexcitation ( )

Switch on (microscopic) de-excitation handling.

◆ EnablePenningTransfer() [1/3]

bool Garfield::MediumMagboltz::EnablePenningTransfer ( )
overridevirtual

Switch on simulation of Penning transfers, using pre-implemented parameterisations of the transfer probability (if available).

Reimplemented from Garfield::MediumGas.

◆ EnablePenningTransfer() [2/3]

bool Garfield::MediumMagboltz::EnablePenningTransfer ( const double r,
const double lambda )
overridevirtual

Switch on simulation of Penning transfers by means of transfer probabilities, for all excitation levels in the mixture.

Parameters
rtransfer probability [0, 1]
lambdaparameter for sampling the distance of the Penning electron with respect to the excitation.

Reimplemented from Garfield::MediumGas.

◆ EnablePenningTransfer() [3/3]

bool Garfield::MediumMagboltz::EnablePenningTransfer ( const double r,
const double lambda,
std::string gasname )
overridevirtual

Switch on simulation of Penning transfers by means of transfer probabilities, for all excitations of a given component.

Reimplemented from Garfield::MediumGas.

◆ EnableRadiationTrapping()

void Garfield::MediumMagboltz::EnableRadiationTrapping ( )

Switch on discrete photoabsorption levels.

◆ EnableThermalMotion()

void Garfield::MediumMagboltz::EnableThermalMotion ( const bool on = true)
inline

Take the thermal motion of the gas at the selected temperature into account in the calculations done by Magboltz.

By the default, this feature is off (static gas at 0 K).

Definition at line 161 of file MediumMagboltz.hh.

161{ m_useGasMotion = on; }
bool m_useGasMotion
Simulate thermal motion of the gas or not (when running Magboltz).

◆ GenerateGasTable()

void Garfield::MediumMagboltz::GenerateGasTable ( const int numCollisions = 10,
const bool verbose = true )

Generate a new gas table (can later be saved to file) by running Magboltz for all electric fields, magnetic fields, and angles in the currently set grid.

◆ GetElectronCollisionRate() [1/2]

double Garfield::MediumMagboltz::GetElectronCollisionRate ( const double e,
const int band )
override

Get the (real) collision rate [ns-1] at a given electron energy e [eV].

◆ GetElectronCollisionRate() [2/2]

double Garfield::MediumMagboltz::GetElectronCollisionRate ( const double e,
const unsigned int level,
const int band )

Get the collision rate [ns-1] for a specific level.

◆ GetElectronNullCollisionRate()

double Garfield::MediumMagboltz::GetElectronNullCollisionRate ( const int band)
override

◆ GetExcitationIonisationLevels()

void Garfield::MediumMagboltz::GetExcitationIonisationLevels ( )
private

◆ GetGasNumberMagboltz()

int Garfield::MediumMagboltz::GetGasNumberMagboltz ( const std::string & input)
static

◆ GetLevel()

bool Garfield::MediumMagboltz::GetLevel ( const unsigned int i,
int & ngas,
int & type,
std::string & descr,
double & e )

Get detailed information about a given cross-section term i.

◆ GetMaxElectronEnergy()

double Garfield::MediumMagboltz::GetMaxElectronEnergy ( ) const
inline

Get the highest electron energy in the table of scattering rates.

Definition at line 40 of file MediumMagboltz.hh.

40{ return m_eMax; }
double m_eMax
Max. electron energy in the collision rate tables.

◆ GetMaxPhotonEnergy()

double Garfield::MediumMagboltz::GetMaxPhotonEnergy ( ) const
inline

Get the highest photon energy in the table of scattering rates.

Definition at line 46 of file MediumMagboltz.hh.

◆ GetNumberOfElectronCollisions() [1/3]

unsigned int Garfield::MediumMagboltz::GetNumberOfElectronCollisions ( ) const

Get the total number of electron collisions.

◆ GetNumberOfElectronCollisions() [2/3]

unsigned int Garfield::MediumMagboltz::GetNumberOfElectronCollisions ( const unsigned int level) const

Get the number of collisions for a specific cross-section term.

◆ GetNumberOfElectronCollisions() [3/3]

unsigned int Garfield::MediumMagboltz::GetNumberOfElectronCollisions ( unsigned int & nElastic,
unsigned int & nIonising,
unsigned int & nAttachment,
unsigned int & nInelastic,
unsigned int & nExcitation,
unsigned int & nSuperelastic ) const

Get the number of collisions broken down by cross-section type.

◆ GetNumberOfLevels()

unsigned int Garfield::MediumMagboltz::GetNumberOfLevels ( )

Get the number of cross-section terms.

◆ GetNumberOfPenningTransfers()

unsigned int Garfield::MediumMagboltz::GetNumberOfPenningTransfers ( ) const
inline

Get the number of Penning transfers that occured since the last reset.

Definition at line 149 of file MediumMagboltz.hh.

149{ return m_nPenning; }
unsigned int m_nPenning
Number of Penning ionisations.

◆ GetNumberOfPhotonCollisions() [1/2]

unsigned int Garfield::MediumMagboltz::GetNumberOfPhotonCollisions ( ) const

Get the total number of photon collisions.

◆ GetNumberOfPhotonCollisions() [2/2]

unsigned int Garfield::MediumMagboltz::GetNumberOfPhotonCollisions ( unsigned int & nElastic,
unsigned int & nIonising,
unsigned int & nInelastic ) const

Get number of photon collisions by collision type.

◆ GetPenningTransfer()

bool Garfield::MediumMagboltz::GetPenningTransfer ( const unsigned int i,
double & r,
double & lambda )

Get the Penning transfer probability and distance of a specific level.

◆ GetPhotonCollisionRate()

double Garfield::MediumMagboltz::GetPhotonCollisionRate ( const double e)
override

◆ Initialise()

bool Garfield::MediumMagboltz::Initialise ( const bool verbose = false)

Initialise the table of scattering rates (called internally when a collision rate is requested and the gas mixture or other parameters have changed).

◆ Mixer()

bool Garfield::MediumMagboltz::Mixer ( const bool verbose = false)
private

◆ PhotonCollision()

bool Garfield::MediumMagboltz::PhotonCollision ( const double e,
int & type,
int & level,
double & e1,
double & ctheta,
std::vector< Secondary > & secondaries )
override

◆ PlotElectronCollisionRates()

void Garfield::MediumMagboltz::PlotElectronCollisionRates ( TPad * pad)

◆ PlotElectronCrossSections()

void Garfield::MediumMagboltz::PlotElectronCrossSections ( const unsigned int i,
TPad * pad )

◆ PlotElectronInverseMeanFreePath()

void Garfield::MediumMagboltz::PlotElectronInverseMeanFreePath ( TPad * pad)

◆ PrintGas()

void Garfield::MediumMagboltz::PrintGas ( )
overridevirtual

Print information about the present gas mixture and available data.

Reimplemented from Garfield::MediumGas.

◆ RateConstantHardSphere()

double Garfield::MediumMagboltz::RateConstantHardSphere ( const double r1,
const double r2,
const int igas1,
const int igas2 ) const
private

◆ RateConstantWK()

double Garfield::MediumMagboltz::RateConstantWK ( const double energy,
const double osc,
const double pacs,
const int igas1,
const int igas2 ) const
private

◆ ResetCollisionCounters()

void Garfield::MediumMagboltz::ResetCollisionCounters ( )

Reset the collision counters.

◆ RunMagboltz()

void Garfield::MediumMagboltz::RunMagboltz ( const double e,
const double b,
const double btheta,
const int ncoll,
bool verbose,
double & vx,
double & vy,
double & vz,
double & wv,
double & wr,
double & dl,
double & dt,
double & alpha,
double & eta,
double & riontof,
double & ratttof,
double & lor,
double & vxerr,
double & vyerr,
double & vzerr,
double & wverr,
double & wrerr,
double & dlerr,
double & dterr,
double & alphaerr,
double & etaerr,
double & riontoferr,
double & ratttoferr,
double & lorerr,
double & alphatof,
std::array< double, 6 > & difftens )

Run Magboltz for a given electric field, magnetic field and angle.

Parameters
[in]eelectric field
[in]bmagnetic field
[in]bthetaangle between electric and magnetic field
[in]ncollnumber of collisions (in multiples of 107) to be simulated
[in]verboseverbosity flag
[out]vx,vy,vzdrift velocity vector
[out]wv,wrflux and bulk drift velocity
[out]dl,dtdiffusion cofficients
[out]alphaTownsend cofficient
[out]etaattachment cofficient
[out]riontofTOF ionisation rate
[out]ratttofTOF attachment rate
[out]lorLorentz angle
[out]vxerr,vyerr,vzerrerrors on drift velocity
[out]wverr,wrerrerrors on TOF drifts
[out]dlerr,dterrerrors on diffusion coefficients
[out]alphaerr,etaerrerrors on Townsend/attachment coefficients
[out]riontoferr,ratttoferrerrors on TOF rates
[out]lorerrerror on Lorentz angle
[out]alphatofeffective Townsend coefficient $(\alpha - \eta)$ calculated using time-of-flight method
[out]difftenscomponents of the diffusion tensor (zz, xx, yy, xz, yz, xy)

◆ SetExcitationScaling()

void Garfield::MediumMagboltz::SetExcitationScaling ( const double r,
std::string gasname )

Multiply all excitation cross-sections by a uniform scaling factor.

◆ SetLineWidth()

void Garfield::MediumMagboltz::SetLineWidth ( const double n)
inline

Set the number of emission line widths within which to apply absorption by discrete lines.

Definition at line 71 of file MediumMagboltz.hh.

◆ SetMaxElectronEnergy()

bool Garfield::MediumMagboltz::SetMaxElectronEnergy ( const double e)

Set the highest electron energy to be included in the table of scattering rates.

◆ SetMaxPhotonEnergy()

bool Garfield::MediumMagboltz::SetMaxPhotonEnergy ( const double e)

Set the highest photon energy to be included in the table of scattering rates.

◆ SetSplittingFunctionFlat()

void Garfield::MediumMagboltz::SetSplittingFunctionFlat ( )

Sample the secondary electron energy from a flat distribution.

◆ SetSplittingFunctionGreenSawada()

void Garfield::MediumMagboltz::SetSplittingFunctionGreenSawada ( )

Sample the secondary electron energy according to the Green-Sawada model.

◆ SetSplittingFunctionOpalBeaty()

void Garfield::MediumMagboltz::SetSplittingFunctionOpalBeaty ( )

Sample the secondary electron energy according to the Opal-Beaty model.

◆ SetupGreenSawada()

void Garfield::MediumMagboltz::SetupGreenSawada ( )
private

◆ Update()

bool Garfield::MediumMagboltz::Update ( const bool verbose = false)
private

Member Data Documentation

◆ csTypeGamma

std::vector<int> Garfield::MediumMagboltz::csTypeGamma
private

Definition at line 434 of file MediumMagboltz.hh.

◆ DxcTypeCollIon

const int Garfield::MediumMagboltz::DxcTypeCollIon
staticprivate

Definition at line 225 of file MediumMagboltz.hh.

◆ DxcTypeCollNonIon

const int Garfield::MediumMagboltz::DxcTypeCollNonIon
staticprivate

Definition at line 226 of file MediumMagboltz.hh.

◆ DxcTypeRad

const int Garfield::MediumMagboltz::DxcTypeRad
staticprivate

Definition at line 224 of file MediumMagboltz.hh.

◆ m_autoEnergyLimit

bool Garfield::MediumMagboltz::m_autoEnergyLimit = true
private

Automatic calculation of the energy limit by Magboltz or not.

Definition at line 234 of file MediumMagboltz.hh.

◆ m_cf

std::vector<std::vector<double> > Garfield::MediumMagboltz::m_cf
private

Definition at line 332 of file MediumMagboltz.hh.

◆ m_cfGamma

std::vector<std::vector<double> > Garfield::MediumMagboltz::m_cfGamma
private

Definition at line 433 of file MediumMagboltz.hh.

◆ m_cfLog

std::vector<std::vector<double> > Garfield::MediumMagboltz::m_cfLog
private

Definition at line 333 of file MediumMagboltz.hh.

◆ m_cfNull

double Garfield::MediumMagboltz::m_cfNull = 0.
private

Null-collision frequency.

Definition at line 338 of file MediumMagboltz.hh.

◆ m_cfTot

std::vector<double> Garfield::MediumMagboltz::m_cfTot
private

Definition at line 329 of file MediumMagboltz.hh.

◆ m_cfTotGamma

std::vector<double> Garfield::MediumMagboltz::m_cfTotGamma
private

Definition at line 431 of file MediumMagboltz.hh.

◆ m_cfTotLog

std::vector<double> Garfield::MediumMagboltz::m_cfTotLog
private

Definition at line 330 of file MediumMagboltz.hh.

◆ m_csType

std::array<int, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_csType
private

Cross-section type.

Definition at line 305 of file MediumMagboltz.hh.

◆ m_deexcitations

std::vector<Deexcitation> Garfield::MediumMagboltz::m_deexcitations
private

Definition at line 410 of file MediumMagboltz.hh.

◆ m_description

std::vector<std::string> Garfield::MediumMagboltz::m_description
private

Level description.

Definition at line 326 of file MediumMagboltz.hh.

◆ m_eAuger1

std::array<double, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_eAuger1
private

Energy imparted to Auger electrons.

Definition at line 313 of file MediumMagboltz.hh.

◆ m_eAuger2

std::array<double, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_eAuger2
private

Definition at line 314 of file MediumMagboltz.hh.

◆ m_eFinalGamma

double Garfield::MediumMagboltz::m_eFinalGamma
private

Definition at line 427 of file MediumMagboltz.hh.

◆ m_eFluorescence

std::array<double, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_eFluorescence
private

Definition at line 316 of file MediumMagboltz.hh.

◆ m_eHigh

double Garfield::MediumMagboltz::m_eHigh
private

Definition at line 243 of file MediumMagboltz.hh.

◆ m_eHighLog

double Garfield::MediumMagboltz::m_eHighLog
private

Definition at line 243 of file MediumMagboltz.hh.

◆ m_eMax

double Garfield::MediumMagboltz::m_eMax
private

Max. electron energy in the collision rate tables.

Definition at line 238 of file MediumMagboltz.hh.

◆ m_energyLoss

std::array<double, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_energyLoss
private

Energy loss.

Definition at line 303 of file MediumMagboltz.hh.

◆ m_eStep

double Garfield::MediumMagboltz::m_eStep
private

Energy spacing in the linear part of the collision rate tables.

Definition at line 240 of file MediumMagboltz.hh.

◆ m_eStepGamma

double Garfield::MediumMagboltz::m_eStepGamma
private

Definition at line 427 of file MediumMagboltz.hh.

◆ m_eStepInv

double Garfield::MediumMagboltz::m_eStepInv
private

Inverse energy spacing.

Definition at line 242 of file MediumMagboltz.hh.

◆ m_hasGreenSawada

std::array<bool, m_nMaxGases> Garfield::MediumMagboltz::m_hasGreenSawada
private

Definition at line 261 of file MediumMagboltz.hh.

◆ m_iDeexcitation

std::array<int, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_iDeexcitation
private

Definition at line 412 of file MediumMagboltz.hh.

◆ m_ionPot

std::array<double, m_nMaxGases> Garfield::MediumMagboltz::m_ionPot
private

Ionisation potentials of each component.

Definition at line 419 of file MediumMagboltz.hh.

◆ m_lambdaPenning

std::array<double, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_lambdaPenning
private

Mean distance of Penning ionisation (by level)

Definition at line 371 of file MediumMagboltz.hh.

◆ m_lnStep

double Garfield::MediumMagboltz::m_lnStep
private

Definition at line 244 of file MediumMagboltz.hh.

◆ m_mgas

std::array<double, m_nMaxGases> Garfield::MediumMagboltz::m_mgas
private

Mass.

Definition at line 252 of file MediumMagboltz.hh.

◆ m_minIonPot

double Garfield::MediumMagboltz::m_minIonPot = -1.
private

Minimum ionisation potential.

Definition at line 421 of file MediumMagboltz.hh.

◆ m_mutex

std::mutex Garfield::MediumMagboltz::m_mutex
private

Mutex.

Definition at line 229 of file MediumMagboltz.hh.

◆ m_nAbsWidths

double Garfield::MediumMagboltz::m_nAbsWidths = 1000.
private

Definition at line 416 of file MediumMagboltz.hh.

◆ m_nAuger1

std::array<unsigned int, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_nAuger1
private

Number of Auger electrons produced in a collision.

Definition at line 310 of file MediumMagboltz.hh.

◆ m_nAuger2

std::array<unsigned int, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_nAuger2
private

Definition at line 311 of file MediumMagboltz.hh.

◆ m_nCollisions

std::array<unsigned int, nCsTypes> Garfield::MediumMagboltz::m_nCollisions
private

Collision counters 0: elastic 1: ionisation 2: attachment 3: inelastic 4: excitation 5: super-elastic.

Definition at line 363 of file MediumMagboltz.hh.

◆ m_nCollisionsDetailed

std::vector<unsigned int> Garfield::MediumMagboltz::m_nCollisionsDetailed
private

Number of collisions for each cross-section term.

Definition at line 365 of file MediumMagboltz.hh.

◆ m_nFluorescence

std::array<unsigned int, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_nFluorescence
private

Definition at line 315 of file MediumMagboltz.hh.

◆ m_nPenning

unsigned int Garfield::MediumMagboltz::m_nPenning = 0
private

Number of Penning ionisations.

Definition at line 373 of file MediumMagboltz.hh.

◆ m_nPhotonCollisions

std::array<unsigned int, nCsTypesGamma> Garfield::MediumMagboltz::m_nPhotonCollisions
private

Definition at line 440 of file MediumMagboltz.hh.

◆ m_nPhotonTerms

unsigned int Garfield::MediumMagboltz::m_nPhotonTerms = 0
private

Definition at line 429 of file MediumMagboltz.hh.

◆ m_nTerms

unsigned int Garfield::MediumMagboltz::m_nTerms = 0
private

Number of different cross-section types in the current gas mixture.

Definition at line 249 of file MediumMagboltz.hh.

◆ m_parGreenSawada

std::array<std::array<double, 5>, m_nMaxGases> Garfield::MediumMagboltz::m_parGreenSawada
private

Green-Sawada splitting parameters [eV] (Γs, Γb, Ts, Ta, Tb).

Definition at line 260 of file MediumMagboltz.hh.

◆ m_rgas

std::array<double, m_nMaxGases> Garfield::MediumMagboltz::m_rgas
private

Recoil energy parameter.

Definition at line 254 of file MediumMagboltz.hh.

◆ m_rPenning

std::array<double, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_rPenning
private

Penning transfer probability (by level)

Definition at line 369 of file MediumMagboltz.hh.

◆ m_s2

std::array<double, m_nMaxGases> Garfield::MediumMagboltz::m_s2
private

Definition at line 255 of file MediumMagboltz.hh.

◆ m_scaleExc

std::array<double, m_nMaxGases> Garfield::MediumMagboltz::m_scaleExc
private

Definition at line 424 of file MediumMagboltz.hh.

◆ m_scatCut

std::vector<std::vector<double> > Garfield::MediumMagboltz::m_scatCut
private

Definition at line 320 of file MediumMagboltz.hh.

◆ m_scatCutLog

std::vector<std::vector<double> > Garfield::MediumMagboltz::m_scatCutLog
private

Definition at line 322 of file MediumMagboltz.hh.

◆ m_scatModel

std::array<int, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_scatModel
private

Definition at line 323 of file MediumMagboltz.hh.

◆ m_scatPar

std::vector<std::vector<double> > Garfield::MediumMagboltz::m_scatPar
private

Definition at line 319 of file MediumMagboltz.hh.

◆ m_scatParLog

std::vector<std::vector<double> > Garfield::MediumMagboltz::m_scatParLog
private

Definition at line 321 of file MediumMagboltz.hh.

◆ m_useAnisotropic

bool Garfield::MediumMagboltz::m_useAnisotropic = true
private

Definition at line 336 of file MediumMagboltz.hh.

◆ m_useCsOutput

bool Garfield::MediumMagboltz::m_useCsOutput = false
private

Flag enabling/disabling output of cross-section table to file.

Definition at line 247 of file MediumMagboltz.hh.

◆ m_useDeexcitation

bool Garfield::MediumMagboltz::m_useDeexcitation = false
private

Flag enabling/disabling detailed simulation of de-excitation process.

Definition at line 377 of file MediumMagboltz.hh.

◆ m_useGasMotion

bool Garfield::MediumMagboltz::m_useGasMotion = false
private

Simulate thermal motion of the gas or not (when running Magboltz).

Definition at line 232 of file MediumMagboltz.hh.

◆ m_useGreenSawada

bool Garfield::MediumMagboltz::m_useGreenSawada = false
private

Sample secondary electron energies using Green-Sawada parameterisation.

Definition at line 266 of file MediumMagboltz.hh.

◆ m_useOpalBeaty

bool Garfield::MediumMagboltz::m_useOpalBeaty = true
private

Sample secondary electron energies using Opal-Beaty parameterisation.

Definition at line 264 of file MediumMagboltz.hh.

◆ m_useRadTrap

bool Garfield::MediumMagboltz::m_useRadTrap = true
private

Flag enabling/disable radiation trapping (absorption of photons discrete excitation lines)

Definition at line 380 of file MediumMagboltz.hh.

◆ m_wOpalBeaty

std::array<double, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_wOpalBeaty
private

Opal-Beaty-Peterson splitting parameter [eV].

Definition at line 257 of file MediumMagboltz.hh.

◆ m_yFluorescence

std::array<double, Magboltz::nMaxLevels> Garfield::MediumMagboltz::m_yFluorescence
private

Fluorescence yield.

Definition at line 308 of file MediumMagboltz.hh.

◆ nCsTypes

int Garfield::MediumMagboltz::nCsTypes = 7
staticconstexprprivate

Definition at line 219 of file MediumMagboltz.hh.

◆ nCsTypesGamma

int Garfield::MediumMagboltz::nCsTypesGamma = 4
staticconstexprprivate

Definition at line 222 of file MediumMagboltz.hh.

◆ nEnergyStepsGamma

int Garfield::MediumMagboltz::nEnergyStepsGamma = 5000
staticconstexprprivate

Definition at line 218 of file MediumMagboltz.hh.

◆ nEnergyStepsLog

int Garfield::MediumMagboltz::nEnergyStepsLog = 1000
staticconstexprprivate

Definition at line 217 of file MediumMagboltz.hh.


The documentation for this class was generated from the following file: