Garfield++ 3.0
A toolkit for the detailed simulation of particle detectors based on ionisation measurement in gases and semiconductors
Loading...
Searching...
No Matches
Garfield::Medium Class Reference

Abstract base class for media. More...

#include <Medium.hh>

+ Inheritance diagram for Garfield::Medium:

Public Member Functions

 Medium ()
 Constructor.
 
virtual ~Medium ()
 Destructor.
 
int GetId () const
 Return the id number of the class instance.
 
const std::string & GetName () const
 Get the medium name/identifier.
 
virtual bool IsGas () const
 Is this medium a gas?
 
virtual bool IsSemiconductor () const
 Is this medium a semiconductor?
 
virtual bool IsConductor () const
 Is this medium a conductor?
 
void SetTemperature (const double t)
 Set the temperature [K].
 
double GetTemperature () const
 Get the temperature [K].
 
void SetPressure (const double p)
 
double GetPressure () const
 
void SetDielectricConstant (const double eps)
 Set the relative static dielectric constant.
 
double GetDielectricConstant () const
 Get the relative static dielectric constant.
 
unsigned int GetNumberOfComponents () const
 Get number of components of the medium.
 
virtual void GetComponent (const unsigned int i, std::string &label, double &f)
 Get the name and fraction of a given component.
 
virtual void SetAtomicNumber (const double z)
 Set the effective atomic number.
 
virtual double GetAtomicNumber () const
 Get the effective atomic number.
 
virtual void SetAtomicWeight (const double a)
 Set the effective atomic weight.
 
virtual double GetAtomicWeight () const
 Get the effective atomic weight.
 
virtual void SetNumberDensity (const double n)
 Set the number density [cm-3].
 
virtual double GetNumberDensity () const
 Get the number density [cm-3].
 
virtual void SetMassDensity (const double rho)
 Set the mass density [g/cm3].
 
virtual double GetMassDensity () const
 Get the mass density [g/cm3].
 
virtual void EnableDrift (const bool on=true)
 Switch electron/ion/hole on/off.
 
virtual void EnablePrimaryIonisation (const bool on=true)
 Make the medium ionisable or non-ionisable.
 
bool IsDriftable () const
 Is charge carrier transport enabled in this medium?
 
bool IsMicroscopic () const
 Does the medium have electron scattering rates?
 
bool IsIonisable () const
 Is charge deposition by charged particles/photon enabled in this medium?
 
void SetW (const double w)
 Set the W value (average energy to produce an electron/ion or e/h pair).
 
double GetW ()
 Get the W value.
 
void SetFanoFactor (const double f)
 Set the Fano factor.
 
double GetFanoFactor ()
 Get the Fano factor.
 
virtual bool ElectronVelocity (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &vx, double &vy, double &vz)
 Drift velocity [cm / ns].
 
virtual bool ElectronDiffusion (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &dl, double &dt)
 Longitudinal and transverse diffusion coefficients [cm1/2].
 
virtual bool ElectronDiffusion (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double cov[3][3])
 
virtual bool ElectronTownsend (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &alpha)
 Ionisation coefficient [cm-1].
 
virtual bool ElectronAttachment (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &eta)
 Attachment coefficient [cm-1].
 
virtual bool ElectronLorentzAngle (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &lor)
 Lorentz angle.
 
virtual double ElectronMobility ()
 Low-field mobility [cm2 V-1 ns-1].
 
virtual double GetElectronEnergy (const double px, const double py, const double pz, double &vx, double &vy, double &vz, const int band=0)
 Dispersion relation (energy vs. wave vector)
 
virtual void GetElectronMomentum (const double e, double &px, double &py, double &pz, int &band)
 
virtual double GetElectronNullCollisionRate (const int band=0)
 Null-collision rate [ns-1].
 
virtual double GetElectronCollisionRate (const double e, const int band=0)
 Collision rate [ns-1] for given electron energy.
 
virtual bool GetElectronCollision (const double e, int &type, int &level, double &e1, double &dx, double &dy, double &dz, std::vector< std::pair< int, double > > &secondaries, int &ndxc, int &band)
 Sample the collision type. Update energy and direction vector.
 
virtual unsigned int GetNumberOfDeexcitationProducts () const
 
virtual bool GetDeexcitationProduct (const unsigned int i, double &t, double &s, int &type, double &energy) const
 
virtual bool HoleVelocity (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &vx, double &vy, double &vz)
 Drift velocity [cm / ns].
 
virtual bool HoleDiffusion (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &dl, double &dt)
 Longitudinal and transverse diffusion coefficients [cm1/2].
 
virtual bool HoleDiffusion (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double cov[3][3])
 Diffusion tensor.
 
virtual bool HoleTownsend (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &alpha)
 Ionisation coefficient [cm-1].
 
virtual bool HoleAttachment (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &eta)
 Attachment coefficient [cm-1].
 
virtual double HoleMobility ()
 Low-field mobility [cm2 V-1 ns-1].
 
virtual bool IonVelocity (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &vx, double &vy, double &vz)
 Drift velocity [cm / ns].
 
virtual bool IonDiffusion (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &dl, double &dt)
 Longitudinal and transverse diffusion coefficients [cm1/2].
 
virtual bool IonDissociation (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, double &diss)
 Dissociation coefficient.
 
virtual double IonMobility ()
 Low-field mobility [cm2 V-1 ns-1].
 
void SetFieldGrid (double emin, double emax, const size_t ne, bool logE, double bmin=0., double bmax=0., const size_t nb=1, double amin=HalfPi, double amax=HalfPi, const size_t na=1)
 Set the range of fields to be covered by the transport tables.
 
void SetFieldGrid (const std::vector< double > &efields, const std::vector< double > &bfields, const std::vector< double > &angles)
 Set the fields and E-B angles to be used in the transport tables.
 
void GetFieldGrid (std::vector< double > &efields, std::vector< double > &bfields, std::vector< double > &angles)
 Get the fields and E-B angles used in the transport tables.
 
bool SetElectronVelocityE (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double v)
 Set an entry in the table of drift speeds along E.
 
bool GetElectronVelocityE (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &v)
 Get an entry in the table of drift speeds along E.
 
bool SetElectronVelocityExB (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double v)
 Set an entry in the table of drift speeds along ExB.
 
bool GetElectronVelocityExB (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &v)
 Get an entry in the table of drift speeds along ExB.
 
bool SetElectronVelocityB (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double v)
 Set an entry in the table of drift speeds along Btrans.
 
bool GetElectronVelocityB (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &v)
 Get an entry in the table of drift speeds along Btrans.
 
bool SetElectronLongitudinalDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double dl)
 Set an entry in the table of longitudinal diffusion coefficients.
 
bool GetElectronLongitudinalDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &dl)
 Get an entry in the table of longitudinal diffusion coefficients.
 
bool SetElectronTransverseDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double dt)
 Set an entry in the table of transverse diffusion coefficients.
 
bool GetElectronTransverseDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &dt)
 Get an entry in the table of transverse diffusion coefficients.
 
bool SetElectronTownsend (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double alpha)
 Set an entry in the table of Townsend coefficients.
 
bool GetElectronTownsend (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &alpha)
 Get an entry in the table of Townsend coefficients.
 
bool SetElectronAttachment (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double eta)
 Set an entry in the table of attachment coefficients.
 
bool GetElectronAttachment (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &eta)
 Get an entry in the table of attachment coefficients.
 
bool SetElectronLorentzAngle (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double lor)
 Set an entry in the table of Lorentz angles.
 
bool GetElectronLorentzAngle (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &lor)
 Get an entry in the table of Lorentz angles.
 
bool SetHoleVelocityE (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double v)
 Set an entry in the table of drift speeds along E.
 
bool GetHoleVelocityE (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &v)
 Get an entry in the table of drift speeds along E.
 
bool SetHoleVelocityExB (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double v)
 Set an entry in the table of drift speeds along ExB.
 
bool GetHoleVelocityExB (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &v)
 Get an entry in the table of drift speeds along ExB.
 
bool SetHoleVelocityB (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double v)
 Set an entry in the table of drift speeds along Btrans.
 
bool GetHoleVelocityB (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &v)
 Get an entry in the table of drift speeds along Btrans.
 
bool SetHoleLongitudinalDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double dl)
 Set an entry in the table of longitudinal diffusion coefficients.
 
bool GetHoleLongitudinalDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &dl)
 Get an entry in the table of longitudinal diffusion coefficients.
 
bool SetHoleTransverseDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double dt)
 Set an entry in the table of transverse diffusion coefficients.
 
bool GetHoleTransverseDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &dt)
 Get an entry in the table of transverse diffusion coefficients.
 
bool SetHoleTownsend (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double alpha)
 Set an entry in the table of Townsend coefficients.
 
bool GetHoleTownsend (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &alpha)
 Get an entry in the table of Townsend coefficients.
 
bool SetHoleAttachment (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double eta)
 Set an entry in the table of attachment coefficients.
 
bool GetHoleAttachment (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &eta)
 Get an entry in the table of attachment coefficients.
 
bool SetIonMobility (const std::vector< double > &fields, const std::vector< double > &mobilities)
 
bool SetIonMobility (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double mu)
 Set an entry in the table of ion mobilities.
 
bool GetIonMobility (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &mu)
 Get an entry in the table of ion mobilities.
 
bool SetIonLongitudinalDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double dl)
 Set an entry in the table of longitudinal diffusion coefficients.
 
bool GetIonLongitudinalDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &dl)
 Get an entry in the table of longitudinal diffusion coefficients.
 
bool SetIonTransverseDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double dt)
 Set an entry in the table of transverse diffusion coefficients.
 
bool GetIonTransverseDiffusion (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &dt)
 Get an entry in the table of transverse diffusion coefficients.
 
bool SetIonDissociation (const unsigned int ie, const unsigned int ib, const unsigned int ia, const double diss)
 Set an entry in the table of dissociation coefficients.
 
bool GetIonDissociation (const unsigned int ie, const unsigned int ib, const unsigned int ia, double &diss)
 Get an entry in the table of dissociation coefficients.
 
virtual void ResetTables ()
 Reset all tables of transport parameters.
 
void ResetElectronVelocity ()
 
void ResetElectronDiffusion ()
 
void ResetElectronTownsend ()
 
void ResetElectronAttachment ()
 
void ResetElectronLorentzAngle ()
 
void ResetHoleVelocity ()
 
void ResetHoleDiffusion ()
 
void ResetHoleTownsend ()
 
void ResetHoleAttachment ()
 
void ResetIonMobility ()
 
void ResetIonDiffusion ()
 
void ResetIonDissociation ()
 
void SetExtrapolationMethodVelocity (const std::string &extrLow, const std::string &extrHigh)
 
void SetExtrapolationMethodDiffusion (const std::string &extrLow, const std::string &extrHigh)
 
void SetExtrapolationMethodTownsend (const std::string &extrLow, const std::string &extrHigh)
 
void SetExtrapolationMethodAttachment (const std::string &extrLow, const std::string &extrHigh)
 
void SetExtrapolationMethodIonMobility (const std::string &extrLow, const std::string &extrHigh)
 
void SetExtrapolationMethodIonDissociation (const std::string &extrLow, const std::string &extrHigh)
 
void SetInterpolationMethodVelocity (const unsigned int intrp)
 Set the degree of polynomial interpolation (usually 2).
 
void SetInterpolationMethodDiffusion (const unsigned int intrp)
 
void SetInterpolationMethodTownsend (const unsigned int intrp)
 
void SetInterpolationMethodAttachment (const unsigned int intrp)
 
void SetInterpolationMethodIonMobility (const unsigned int intrp)
 
void SetInterpolationMethodIonDissociation (const unsigned int intrp)
 
virtual double ScaleElectricField (const double e) const
 
virtual double UnScaleElectricField (const double e) const
 
virtual double ScaleVelocity (const double v) const
 
virtual double ScaleDiffusion (const double d) const
 
virtual double ScaleDiffusionTensor (const double d) const
 
virtual double ScaleTownsend (const double alpha) const
 
virtual double ScaleAttachment (const double eta) const
 
virtual double ScaleLorentzAngle (const double lor) const
 
virtual double ScaleDissociation (const double diss) const
 
virtual bool GetOpticalDataRange (double &emin, double &emax, const unsigned int i=0)
 Get the energy range [eV] of the available optical data.
 
virtual bool GetDielectricFunction (const double e, double &eps1, double &eps2, const unsigned int i=0)
 Get the complex dielectric function at a given energy.
 
virtual bool GetPhotoAbsorptionCrossSection (const double e, double &sigma, const unsigned int i=0)
 
virtual double GetPhotonCollisionRate (const double e)
 
virtual bool GetPhotonCollision (const double e, int &type, int &level, double &e1, double &ctheta, int &nsec, double &esec)
 
void EnableDebugging ()
 Switch on/off debugging messages.
 
void DisableDebugging ()
 

Protected Member Functions

bool Velocity (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, const std::vector< std::vector< std::vector< double > > > &velE, const std::vector< std::vector< std::vector< double > > > &velB, const std::vector< std::vector< std::vector< double > > > &velX, const double q, double &vx, double &vy, double &vz) const
 
bool Diffusion (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, const std::vector< std::vector< std::vector< double > > > &difL, const std::vector< std::vector< std::vector< double > > > &difT, double &dl, double &dt) const
 
bool Diffusion (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, const std::vector< std::vector< std::vector< std::vector< double > > > > &diff, double cov[3][3]) const
 
bool Alpha (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, const std::vector< std::vector< std::vector< double > > > &tab, unsigned int intp, const unsigned int thr, const std::pair< unsigned int, unsigned int > &extr, double &alpha) const
 
double GetAngle (const double ex, const double ey, const double ez, const double bx, const double by, const double bz, const double e, const double b) const
 
bool Interpolate (const double e, const double b, const double a, const std::vector< std::vector< std::vector< double > > > &table, double &y, const unsigned int intp, const std::pair< unsigned int, unsigned int > &extr) const
 
double Interpolate1D (const double e, const std::vector< double > &table, const std::vector< double > &fields, const unsigned int intpMeth, const std::pair< unsigned int, unsigned int > &extr) const
 
bool SetEntry (const unsigned int i, const unsigned int j, const unsigned int k, const std::string &fcn, std::vector< std::vector< std::vector< double > > > &tab, const double val)
 
bool GetEntry (const unsigned int i, const unsigned int j, const unsigned int k, const std::string &fcn, const std::vector< std::vector< std::vector< double > > > &tab, double &val) const
 
void SetExtrapolationMethod (const std::string &low, const std::string &high, std::pair< unsigned int, unsigned int > &extr, const std::string &fcn)
 
bool GetExtrapolationIndex (std::string str, unsigned int &nb) const
 
unsigned int SetThreshold (const std::vector< std::vector< std::vector< double > > > &tab) const
 
void Clone (std::vector< std::vector< std::vector< double > > > &tab, const std::vector< double > &efields, const std::vector< double > &bfields, const std::vector< double > &angles, const unsigned int intp, const std::pair< unsigned int, unsigned int > &extr, const double init, const std::string &label)
 
void Clone (std::vector< std::vector< std::vector< std::vector< double > > > > &tab, const unsigned int n, const std::vector< double > &efields, const std::vector< double > &bfields, const std::vector< double > &angles, const unsigned int intp, const std::pair< unsigned int, unsigned int > &extr, const double init, const std::string &label)
 
void Init (const size_t nE, const size_t nB, const size_t nA, std::vector< std::vector< std::vector< double > > > &tab, const double val)
 
void Init (const size_t nE, const size_t nB, const size_t nA, const size_t nT, std::vector< std::vector< std::vector< std::vector< double > > > > &tab, const double val)
 

Protected Attributes

std::string m_className = "Medium"
 
int m_id
 
std::string m_name = ""
 
double m_temperature = 293.15
 
double m_pressure = 760.
 
double m_epsilon = 1.
 
unsigned int m_nComponents = 1
 
double m_z = 1.
 
double m_a = 0.
 
double m_density = 0.
 
bool m_driftable = false
 
bool m_microscopic = false
 
bool m_ionisable = false
 
double m_w = 0.
 
double m_fano = 0.
 
bool m_isChanged = true
 
bool m_debug = false
 
std::vector< double > m_eFields
 
std::vector< double > m_bFields
 
std::vector< double > m_bAngles
 
bool m_tab2d = false
 
std::vector< std::vector< std::vector< double > > > m_eVelE
 
std::vector< std::vector< std::vector< double > > > m_eVelX
 
std::vector< std::vector< std::vector< double > > > m_eVelB
 
std::vector< std::vector< std::vector< double > > > m_eDifL
 
std::vector< std::vector< std::vector< double > > > m_eDifT
 
std::vector< std::vector< std::vector< double > > > m_eAlp
 
std::vector< std::vector< std::vector< double > > > m_eAtt
 
std::vector< std::vector< std::vector< double > > > m_eLor
 
std::vector< std::vector< std::vector< std::vector< double > > > > m_eDifM
 
std::vector< std::vector< std::vector< double > > > m_hVelE
 
std::vector< std::vector< std::vector< double > > > m_hVelX
 
std::vector< std::vector< std::vector< double > > > m_hVelB
 
std::vector< std::vector< std::vector< double > > > m_hDifL
 
std::vector< std::vector< std::vector< double > > > m_hDifT
 
std::vector< std::vector< std::vector< double > > > m_hAlp
 
std::vector< std::vector< std::vector< double > > > m_hAtt
 
std::vector< std::vector< std::vector< std::vector< double > > > > m_hDifM
 
std::vector< std::vector< std::vector< double > > > m_iMob
 
std::vector< std::vector< std::vector< double > > > m_iDifL
 
std::vector< std::vector< std::vector< double > > > m_iDifT
 
std::vector< std::vector< std::vector< double > > > m_iDis
 
unsigned int m_eThrAlp = 0
 
unsigned int m_eThrAtt = 0
 
unsigned int m_hThrAlp = 0
 
unsigned int m_hThrAtt = 0
 
unsigned int m_iThrDis = 0
 
std::pair< unsigned int, unsigned int > m_extrVel = {0, 1}
 
std::pair< unsigned int, unsigned int > m_extrDif = {0, 1}
 
std::pair< unsigned int, unsigned int > m_extrAlp = {0, 1}
 
std::pair< unsigned int, unsigned int > m_extrAtt = {0, 1}
 
std::pair< unsigned int, unsigned int > m_extrLor = {0, 1}
 
std::pair< unsigned int, unsigned int > m_extrMob = {0, 1}
 
std::pair< unsigned int, unsigned int > m_extrDis = {0, 1}
 
unsigned int m_intpVel = 2
 
unsigned int m_intpDif = 2
 
unsigned int m_intpAlp = 2
 
unsigned int m_intpAtt = 2
 
unsigned int m_intpLor = 2
 
unsigned int m_intpMob = 2
 
unsigned int m_intpDis = 2
 

Static Protected Attributes

static int m_idCounter = -1
 

Detailed Description

Abstract base class for media.

Definition at line 13 of file Medium.hh.

Constructor & Destructor Documentation

◆ Medium()

Garfield::Medium::Medium ( )

Constructor.

Definition at line 60 of file Medium.cc.

60 : m_id(++m_idCounter) {
61 // Initialise the transport tables.
62 m_bFields.assign(1, 0.);
63 m_bAngles.assign(1, HalfPi);
64
65 // Set the default grid.
66 SetFieldGrid(100., 100000., 20, true, 0., 0., 1, HalfPi, HalfPi, 1);
67}
std::vector< double > m_bFields
Definition: Medium.hh:547
static int m_idCounter
Definition: Medium.hh:508
void SetFieldGrid(double emin, double emax, const size_t ne, bool logE, double bmin=0., double bmax=0., const size_t nb=1, double amin=HalfPi, double amax=HalfPi, const size_t na=1)
Set the range of fields to be covered by the transport tables.
Definition: Medium.cc:693
std::vector< double > m_bAngles
Definition: Medium.hh:548

◆ ~Medium()

Garfield::Medium::~Medium ( )
virtual

Destructor.

Definition at line 69 of file Medium.cc.

69{}

Member Function Documentation

◆ Alpha()

bool Garfield::Medium::Alpha ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
const std::vector< std::vector< std::vector< double > > > &  tab,
unsigned int  intp,
const unsigned int  thr,
const std::pair< unsigned int, unsigned int > &  extr,
double &  alpha 
) const
protected

Definition at line 348 of file Medium.cc.

353 {
354
355 alpha = 0.;
356 if (tab.empty()) return false;
357
358 // Compute the magnitude of the electric field.
359 const double e = sqrt(ex * ex + ey * ey + ez * ez);
360 const double e0 = ScaleElectricField(e);
361 if (e < Small || e0 < Small) return true;
362
363 // Compute the magnitude of the magnetic field.
364 const double b = m_tab2d ? sqrt(bx * bx + by * by + bz * bz) : 0.;
365 // Compute the angle between B field and E field.
366 const double ebang = m_tab2d ? GetAngle(ex, ey, ez, bx, by, bz, e, b) : 0.;
367
368 // Interpolate.
369 if (e0 < m_eFields[thr]) intp = 1;
370 if (!Interpolate(e0, b, ebang, tab, alpha, intp, extr)) alpha = -30.;
371 if (alpha < -20.) {
372 alpha = 0.;
373 } else {
374 alpha = exp(alpha);
375 }
376 return true;
377}
bool Interpolate(const double e, const double b, const double a, const std::vector< std::vector< std::vector< double > > > &table, double &y, const unsigned int intp, const std::pair< unsigned int, unsigned int > &extr) const
Definition: Medium.cc:1210
double GetAngle(const double ex, const double ey, const double ez, const double bx, const double by, const double bz, const double e, const double b) const
Definition: Medium.cc:1194
virtual double ScaleElectricField(const double e) const
Definition: Medium.hh:476
std::vector< double > m_eFields
Definition: Medium.hh:546
DoubleAc exp(const DoubleAc &f)
Definition: DoubleAc.cpp:377
DoubleAc sqrt(const DoubleAc &f)
Definition: DoubleAc.cpp:314

Referenced by ElectronAttachment(), ElectronTownsend(), HoleAttachment(), HoleTownsend(), and IonDissociation().

◆ Clone() [1/2]

void Garfield::Medium::Clone ( std::vector< std::vector< std::vector< double > > > &  tab,
const std::vector< double > &  efields,
const std::vector< double > &  bfields,
const std::vector< double > &  angles,
const unsigned int  intp,
const std::pair< unsigned int, unsigned int > &  extr,
const double  init,
const std::string &  label 
)
protected

Definition at line 922 of file Medium.cc.

928 {
929 if (m_debug) {
930 std::cout << m_className << "::Clone: Copying " << lbl << " to new grid.\n";
931 }
932
933 if (tab.empty()) {
934 if (m_debug) std::cout << m_className << "::Clone: Table is empty.\n";
935 return;
936 }
937 // Get the dimensions of the new grid.
938 const auto nE = efields.size();
939 const auto nB = bfields.size();
940 const auto nA = angles.size();
941
942 // Create a temporary table to store the values at the new grid points.
943 std::vector<std::vector<std::vector<double> > > tabClone;
944 Init(nE, nB, nA, tabClone, init);
945
946 // Fill the temporary table.
947 for (size_t i = 0; i < nE; ++i) {
948 const double e = efields[i];
949 for (size_t j = 0; j < nB; ++j) {
950 const double b = bfields[j];
951 for (size_t k = 0; k < nA; ++k) {
952 const double a = angles[k];
953 double val = 0.;
954 if (!Interpolate(e, b, a, tab, val, intp, extr)) {
955 std::cerr << m_className << "::Clone:\n"
956 << " Interpolation of " << lbl << " failed.\n"
957 << " Cannot copy value to new grid at E = " << e
958 << ", B = " << b << ", angle: " << a << "\n";
959 continue;
960 }
961 tabClone[k][j][i] = val;
962 }
963 }
964 }
965 // Copy the values to the original table.
966 tab.swap(tabClone);
967 tabClone.clear();
968}
void Init(const size_t nE, const size_t nB, const size_t nA, std::vector< std::vector< std::vector< double > > > &tab, const double val)
Definition: Medium.cc:1304
std::string m_className
Definition: Medium.hh:506

Referenced by SetFieldGrid().

◆ Clone() [2/2]

void Garfield::Medium::Clone ( std::vector< std::vector< std::vector< std::vector< double > > > > &  tab,
const unsigned int  n,
const std::vector< double > &  efields,
const std::vector< double > &  bfields,
const std::vector< double > &  angles,
const unsigned int  intp,
const std::pair< unsigned int, unsigned int > &  extr,
const double  init,
const std::string &  label 
)
protected

Definition at line 970 of file Medium.cc.

975 {
976 // If the table does not exist, do nothing.
977 if (tab.empty()) return;
978
979 // Get the dimensions of the new grid.
980 const unsigned int nE = efields.size();
981 const unsigned int nB = bfields.size();
982 const unsigned int nA = angles.size();
983
984 // Create a temporary table to store the values at the new grid points.
985 std::vector<std::vector<std::vector<std::vector<double> > > > tabClone;
986 Init(nE, nB, nA, n, tabClone, init);
987
988 // Fill the temporary table.
989 for (unsigned int l = 0; l < n; ++l) {
990 for (unsigned int i = 0; i < nE; ++i) {
991 const double e = efields[i];
992 for (unsigned int j = 0; j < nB; ++j) {
993 const double b = bfields[j];
994 for (unsigned int k = 0; k < nA; ++k) {
995 const double a = angles[k];
996 double val = 0.;
997 if (!Interpolate(e, b, a, tab[l], val, intp, extr)) {
998 std::cerr << m_className << "::Clone:\n"
999 << " Interpolation of " << lbl << " failed.\n"
1000 << " Cannot copy value to new grid at index " << l
1001 << ", E = " << e << ", B = " << b << ", angle: " << a
1002 << "\n";
1003 continue;
1004 }
1005 tabClone[l][k][j][i] = val;
1006 }
1007 }
1008 }
1009 }
1010 // Copy the values to the original table.
1011 tab.swap(tabClone);
1012}

◆ Diffusion() [1/2]

bool Garfield::Medium::Diffusion ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
const std::vector< std::vector< std::vector< double > > > &  difL,
const std::vector< std::vector< std::vector< double > > > &  difT,
double &  dl,
double &  dt 
) const
protected

Definition at line 269 of file Medium.cc.

273 {
274
275 dl = dt = 0.;
276 // Compute the magnitude of the electric field.
277 const double e = sqrt(ex * ex + ey * ey + ez * ez);
278 const double e0 = ScaleElectricField(e);
279 if (e < Small || e0 < Small) return true;
280
281 // Compute the magnitude of the magnetic field.
282 const double b = m_tab2d ? sqrt(bx * bx + by * by + bz * bz) : 0.;
283 // Compute the angle between B field and E field.
284 const double ebang = m_tab2d ? GetAngle(ex, ey, ez, bx, by, bz, e, b) : 0.;
285
286 // Interpolate.
287 if (!difL.empty()) {
288 if (!Interpolate(e0, b, ebang, difL, dl, m_intpDif, m_extrDif)) dl = 0.;
289 }
290 if (!difT.empty()) {
291 if (!Interpolate(e0, b, ebang, difT, dt, m_intpDif, m_extrDif)) dt = 0.;
292 }
293
294 // If no data available, calculate
295 // the diffusion coefficients using the Einstein relation
296 if (difL.empty() || difT.empty()) {
297 const double d = sqrt(2. * BoltzmannConstant * m_temperature / e);
298 if (difL.empty()) dl = d;
299 if (difT.empty()) dt = d;
300 }
301 // Verify values and apply scaling.
302 dl = ScaleDiffusion(std::max(dl, 0.));
303 dt = ScaleDiffusion(std::max(dt, 0.));
304 return true;
305}
virtual double ScaleDiffusion(const double d) const
Definition: Medium.hh:479
unsigned int m_intpDif
Definition: Medium.hh:599
std::pair< unsigned int, unsigned int > m_extrDif
Definition: Medium.hh:590
double m_temperature
Definition: Medium.hh:515

Referenced by ElectronDiffusion(), HoleDiffusion(), and IonDiffusion().

◆ Diffusion() [2/2]

bool Garfield::Medium::Diffusion ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
const std::vector< std::vector< std::vector< std::vector< double > > > > &  diff,
double  cov[3][3] 
) const
protected

Definition at line 307 of file Medium.cc.

310 {
311
312 // Initialise the tensor.
313 cov[0][0] = cov[0][1] = cov[0][2] = 0.;
314 cov[1][0] = cov[1][1] = cov[1][2] = 0.;
315 cov[2][0] = cov[2][1] = cov[2][2] = 0.;
316
317 if (diff.empty()) return false;
318
319 // Compute the magnitude of the electric field.
320 const double e = sqrt(ex * ex + ey * ey + ez * ez);
321 const double e0 = ScaleElectricField(e);
322 if (e < Small || e0 < Small) return true;
323
324 // Compute the magnitude of the magnetic field.
325 const double b = m_tab2d ? sqrt(bx * bx + by * by + bz * bz) : 0.;
326 // Compute the angle between B field and E field.
327 const double ebang = m_tab2d ? GetAngle(ex, ey, ez, bx, by, bz, e, b) : 0.;
328
329 for (int j = 0; j < 6; ++j) {
330 // Interpolate.
331 double y = 0.;
332 if (!Interpolate(e0, b, ebang, diff[j], y, m_intpDif, m_extrDif)) y = 0.;
333 // Apply scaling.
335 if (j < 3) {
336 cov[j][j] = y;
337 } else if (j == 3) {
338 cov[0][1] = cov[1][0] = y;
339 } else if (j == 4) {
340 cov[0][2] = cov[2][0] = y;
341 } else if (j == 5) {
342 cov[1][2] = cov[2][1] = y;
343 }
344 }
345 return true;
346}
virtual double ScaleDiffusionTensor(const double d) const
Definition: Medium.hh:480

◆ DisableDebugging()

void Garfield::Medium::DisableDebugging ( )
inline

Definition at line 503 of file Medium.hh.

503{ m_debug = false; }

Referenced by GarfieldPhysics::InitializePhysics().

◆ ElectronAttachment()

bool Garfield::Medium::ElectronAttachment ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  eta 
)
virtual

Attachment coefficient [cm-1].

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumSilicon.

Definition at line 416 of file Medium.cc.

418 {
419
420 if (!Alpha(ex, ey, ez, bx, by, bz, m_eAtt, m_intpAtt, m_eThrAtt, m_extrAtt,
421 eta)) {
422 return false;
423 }
424 // Apply scaling.
425 eta = ScaleAttachment(eta);
426 return true;
427}
virtual double ScaleAttachment(const double eta) const
Definition: Medium.hh:482
std::pair< unsigned int, unsigned int > m_extrAtt
Definition: Medium.hh:592
std::vector< std::vector< std::vector< double > > > m_eAtt
Definition: Medium.hh:559
bool Alpha(const double ex, const double ey, const double ez, const double bx, const double by, const double bz, const std::vector< std::vector< std::vector< double > > > &tab, unsigned int intp, const unsigned int thr, const std::pair< unsigned int, unsigned int > &extr, double &alpha) const
Definition: Medium.cc:348
unsigned int m_intpAtt
Definition: Medium.hh:601
unsigned int m_eThrAtt
Definition: Medium.hh:583

Referenced by Garfield::MediumCdTe::ElectronAttachment(), Garfield::MediumGaAs::ElectronAttachment(), Garfield::MediumSilicon::ElectronAttachment(), and Garfield::ViewMedium::EvaluateFunction().

◆ ElectronDiffusion() [1/2]

bool Garfield::Medium::ElectronDiffusion ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  dl,
double &  dt 
)
virtual

Longitudinal and transverse diffusion coefficients [cm1/2].

Definition at line 387 of file Medium.cc.

390 {
391
392 return Diffusion(ex, ey, ez, bx, by, bz, m_eDifL, m_eDifT, dl, dt);
393}
bool Diffusion(const double ex, const double ey, const double ez, const double bx, const double by, const double bz, const std::vector< std::vector< std::vector< double > > > &difL, const std::vector< std::vector< std::vector< double > > > &difT, double &dl, double &dt) const
Definition: Medium.cc:269
std::vector< std::vector< std::vector< double > > > m_eDifL
Definition: Medium.hh:556
std::vector< std::vector< std::vector< double > > > m_eDifT
Definition: Medium.hh:557

Referenced by Garfield::ViewMedium::EvaluateFunction().

◆ ElectronDiffusion() [2/2]

bool Garfield::Medium::ElectronDiffusion ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double  cov[3][3] 
)
virtual

Diffusion tensor: diagonal elements are the diffusion coefficients [cm] along e, btrans, e x b, off-diagonal elements are the covariances

Definition at line 395 of file Medium.cc.

398 {
399
400 return Diffusion(ex, ey, ez, bx, by, bz, m_eDifM, cov);
401}
std::vector< std::vector< std::vector< std::vector< double > > > > m_eDifM
Definition: Medium.hh:562

◆ ElectronLorentzAngle()

bool Garfield::Medium::ElectronLorentzAngle ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  lor 
)
virtual

Lorentz angle.

Definition at line 429 of file Medium.cc.

432 {
433 lor = 0.;
434 if (m_eLor.empty()) return false;
435
436 // Compute the magnitude of the electric field.
437 const double e = sqrt(ex * ex + ey * ey + ez * ez);
438 const double e0 = ScaleElectricField(e);
439 if (e < Small || e0 < Small) return true;
440
441 // Compute the magnitude of the magnetic field.
442 const double b = m_tab2d ? sqrt(bx * bx + by * by + bz * bz) : 0.;
443 // Compute the angle between B field and E field.
444 const double ebang = m_tab2d ? GetAngle(ex, ey, ez, bx, by, bz, e, b) : 0.;
445
446 // Interpolate.
447 if (!Interpolate(e0, b, ebang, m_eLor, lor, m_intpLor, m_extrLor)) lor = 0.;
448 // Apply scaling.
449 lor = ScaleLorentzAngle(lor);
450 return true;
451}
virtual double ScaleLorentzAngle(const double lor) const
Definition: Medium.hh:483
std::vector< std::vector< std::vector< double > > > m_eLor
Definition: Medium.hh:560
std::pair< unsigned int, unsigned int > m_extrLor
Definition: Medium.hh:593
unsigned int m_intpLor
Definition: Medium.hh:602

Referenced by Garfield::ViewMedium::EvaluateFunction().

◆ ElectronMobility()

double Garfield::Medium::ElectronMobility ( )
virtual

Low-field mobility [cm2 V-1 ns-1].

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumSilicon.

Definition at line 453 of file Medium.cc.

453 {
454 if (m_eVelE.empty()) return -1.;
455 return m_eVelE[0][0][0] / UnScaleElectricField(m_eFields[0]);
456}
virtual double UnScaleElectricField(const double e) const
Definition: Medium.hh:477
std::vector< std::vector< std::vector< double > > > m_eVelE
Definition: Medium.hh:553

◆ ElectronTownsend()

bool Garfield::Medium::ElectronTownsend ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  alpha 
)
virtual

Ionisation coefficient [cm-1].

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumSilicon.

Definition at line 403 of file Medium.cc.

405 {
406
407 if (!Alpha(ex, ey, ez, bx, by, bz, m_eAlp, m_intpAlp, m_eThrAlp, m_extrAlp,
408 alpha)) {
409 return false;
410 }
411 // Apply scaling.
412 alpha = ScaleTownsend(alpha);
413 return true;
414}
std::vector< std::vector< std::vector< double > > > m_eAlp
Definition: Medium.hh:558
virtual double ScaleTownsend(const double alpha) const
Definition: Medium.hh:481
std::pair< unsigned int, unsigned int > m_extrAlp
Definition: Medium.hh:591
unsigned int m_eThrAlp
Definition: Medium.hh:582
unsigned int m_intpAlp
Definition: Medium.hh:600

Referenced by Garfield::MediumCdTe::ElectronTownsend(), Garfield::MediumGaAs::ElectronTownsend(), Garfield::MediumSilicon::ElectronTownsend(), and Garfield::ViewMedium::EvaluateFunction().

◆ ElectronVelocity()

bool Garfield::Medium::ElectronVelocity ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  vx,
double &  vy,
double &  vz 
)
virtual

Drift velocity [cm / ns].

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumSilicon.

Definition at line 379 of file Medium.cc.

381 {
382
383 return Velocity(ex, ey, ez, bx, by, bz, m_eVelE, m_eVelB, m_eVelX, -1.,
384 vx, vy, vz);
385}
std::vector< std::vector< std::vector< double > > > m_eVelX
Definition: Medium.hh:554
std::vector< std::vector< std::vector< double > > > m_eVelB
Definition: Medium.hh:555
bool Velocity(const double ex, const double ey, const double ez, const double bx, const double by, const double bz, const std::vector< std::vector< std::vector< double > > > &velE, const std::vector< std::vector< std::vector< double > > > &velB, const std::vector< std::vector< std::vector< double > > > &velX, const double q, double &vx, double &vy, double &vz) const
Definition: Medium.cc:160

Referenced by Garfield::MediumCdTe::ElectronVelocity(), Garfield::MediumGaAs::ElectronVelocity(), Garfield::MediumSilicon::ElectronVelocity(), and Garfield::ViewMedium::EvaluateFunction().

◆ EnableDebugging()

void Garfield::Medium::EnableDebugging ( )
inline

Switch on/off debugging messages.

Definition at line 502 of file Medium.hh.

502{ m_debug = true; }

Referenced by GarfieldPhysics::InitializePhysics().

◆ EnableDrift()

virtual void Garfield::Medium::EnableDrift ( const bool  on = true)
inlinevirtual

◆ EnablePrimaryIonisation()

virtual void Garfield::Medium::EnablePrimaryIonisation ( const bool  on = true)
inlinevirtual

◆ GetAngle()

double Garfield::Medium::GetAngle ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
const double  e,
const double  b 
) const
protected

Definition at line 1194 of file Medium.cc.

1196 {
1197 const double eb = e * b;
1198 if (eb <= 0.) return m_bAngles[0];
1199 const double einb = fabs(ex * bx + ey * by + ez * bz);
1200 if (einb > 0.2 * eb) {
1201 const double ebxy = ex * by - ey * bx;
1202 const double ebxz = ex * bz - ez * bx;
1203 const double ebzy = ez * by - ey * bz;
1204 return asin(
1205 std::min(1., sqrt(ebxy * ebxy + ebxz * ebxz + ebzy * ebzy) / eb));
1206 }
1207 return acos(std::min(1., einb / eb));
1208}
DoubleAc acos(const DoubleAc &f)
Definition: DoubleAc.cpp:490
DoubleAc fabs(const DoubleAc &f)
Definition: DoubleAc.h:615
DoubleAc asin(const DoubleAc &f)
Definition: DoubleAc.cpp:470

Referenced by Alpha(), Diffusion(), ElectronLorentzAngle(), IonVelocity(), and Velocity().

◆ GetAtomicNumber()

virtual double Garfield::Medium::GetAtomicNumber ( ) const
inlinevirtual

Get the effective atomic number.

Reimplemented in Garfield::MediumGas.

Definition at line 52 of file Medium.hh.

52{ return m_z; }

◆ GetAtomicWeight()

virtual double Garfield::Medium::GetAtomicWeight ( ) const
inlinevirtual

Get the effective atomic weight.

Reimplemented in Garfield::MediumGas.

Definition at line 56 of file Medium.hh.

56{ return m_a; }

◆ GetComponent()

void Garfield::Medium::GetComponent ( const unsigned int  i,
std::string &  label,
double &  f 
)
virtual

Get the name and fraction of a given component.

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumGas.

Definition at line 105 of file Medium.cc.

105 {
106 if (i >= m_nComponents) {
107 std::cerr << m_className << "::GetComponent: Index out of range.\n";
108 }
109
110 label = m_name;
111 f = 1.;
112}
std::string m_name
Definition: Medium.hh:513
unsigned int m_nComponents
Definition: Medium.hh:521

◆ GetDeexcitationProduct()

bool Garfield::Medium::GetDeexcitationProduct ( const unsigned int  i,
double &  t,
double &  s,
int &  type,
double &  energy 
) const
virtual

Reimplemented in Garfield::MediumMagboltz.

Definition at line 504 of file Medium.cc.

506 {
507 if (m_debug) PrintNotImplemented(m_className, "GetDeexcitationProduct");
508 t = s = energy = 0.;
509 type = 0;
510 return false;
511}

◆ GetDielectricConstant()

double Garfield::Medium::GetDielectricConstant ( ) const
inline

Get the relative static dielectric constant.

Definition at line 42 of file Medium.hh.

42{ return m_epsilon; }
double m_epsilon
Definition: Medium.hh:519

Referenced by Garfield::ComponentNeBem2d::Initialise(), and Garfield::ComponentNeBem3d::Initialise().

◆ GetDielectricFunction()

bool Garfield::Medium::GetDielectricFunction ( const double  e,
double &  eps1,
double &  eps2,
const unsigned int  i = 0 
)
virtual

Get the complex dielectric function at a given energy.

Reimplemented in Garfield::MediumSilicon.

Definition at line 636 of file Medium.cc.

637 {
638 if (i >= m_nComponents) {
639 std::cerr << m_className << "::GetDielectricFunction: Index out of range.\n";
640 return false;
641 }
642
643 if (e < 0.) {
644 std::cerr << m_className << "::GetDielectricFunction: Energy must be > 0.\n";
645 return false;
646 }
647
648 if (m_debug) PrintNotImplemented(m_className, "GetDielectricFunction");
649 eps1 = 1.;
650 eps2 = 0.;
651 return false;
652}

◆ GetElectronAttachment()

bool Garfield::Medium::GetElectronAttachment ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  eta 
)
inline

Get an entry in the table of attachment coefficients.

Definition at line 278 of file Medium.hh.

279 {
280 return GetEntry(ie, ib, ia, "GetElectronAttachment", m_eAtt, eta);
281 }
bool GetEntry(const unsigned int i, const unsigned int j, const unsigned int k, const std::string &fcn, const std::vector< std::vector< std::vector< double > > > &tab, double &val) const
Definition: Medium.cc:888

◆ GetElectronCollision()

bool Garfield::Medium::GetElectronCollision ( const double  e,
int &  type,
int &  level,
double &  e1,
double &  dx,
double &  dy,
double &  dz,
std::vector< std::pair< int, double > > &  secondaries,
int &  ndxc,
int &  band 
)
virtual

Sample the collision type. Update energy and direction vector.

Reimplemented in Garfield::MediumMagboltz, and Garfield::MediumSilicon.

Definition at line 491 of file Medium.cc.

494 {
495 type = level = -1;
496 e1 = e;
497 ndxc = band = 0;
498 RndmDirection(dx, dy, dz);
499
500 if (m_debug) PrintNotImplemented(m_className, "GetElectronCollision");
501 return false;
502}
void RndmDirection(double &dx, double &dy, double &dz, const double length=1.)
Draw a random (isotropic) direction vector.
Definition: Random.hh:107

◆ GetElectronCollisionRate()

double Garfield::Medium::GetElectronCollisionRate ( const double  e,
const int  band = 0 
)
virtual

Collision rate [ns-1] for given electron energy.

Reimplemented in Garfield::MediumMagboltz, and Garfield::MediumSilicon.

Definition at line 485 of file Medium.cc.

486 {
487 if (m_debug) PrintNotImplemented(m_className, "GetElectronCollisionRate");
488 return 0.;
489}

◆ GetElectronEnergy()

double Garfield::Medium::GetElectronEnergy ( const double  px,
const double  py,
const double  pz,
double &  vx,
double &  vy,
double &  vz,
const int  band = 0 
)
virtual

Dispersion relation (energy vs. wave vector)

Reimplemented in Garfield::MediumSilicon.

Definition at line 458 of file Medium.cc.

460 {
461 if (band > 0) {
462 std::cerr << m_className << "::GetElectronEnergy:\n";
463 std::cerr << " Unknown band index.\n";
464 }
465
466 vx = SpeedOfLight * px / ElectronMass;
467 vy = SpeedOfLight * py / ElectronMass;
468 vz = SpeedOfLight * pz / ElectronMass;
469
470 return 0.5 * (px * px + py * py + pz * pz) / ElectronMass;
471}

◆ GetElectronLongitudinalDiffusion()

bool Garfield::Medium::GetElectronLongitudinalDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  dl 
)
inline

Get an entry in the table of longitudinal diffusion coefficients.

Definition at line 241 of file Medium.hh.

243 {
244 return GetEntry(ie, ib, ia, "GetElectronLongitudinalDiffusion",
245 m_eDifL, dl);
246 }

◆ GetElectronLorentzAngle()

bool Garfield::Medium::GetElectronLorentzAngle ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  lor 
)
inline

Get an entry in the table of Lorentz angles.

Definition at line 289 of file Medium.hh.

290 {
291 return GetEntry(ie, ib, ia, "GetElectronLorentzAngle", m_eLor, lor);
292 }

◆ GetElectronMomentum()

void Garfield::Medium::GetElectronMomentum ( const double  e,
double &  px,
double &  py,
double &  pz,
int &  band 
)
virtual

Sample the momentum vector for a given energy (only meaningful in semiconductors).

Reimplemented in Garfield::MediumSilicon.

Definition at line 473 of file Medium.cc.

474 {
475 const double p = sqrt(2. * ElectronMass * e) / SpeedOfLight;
476 RndmDirection(px, py, pz, p);
477 band = -1;
478}

◆ GetElectronNullCollisionRate()

double Garfield::Medium::GetElectronNullCollisionRate ( const int  band = 0)
virtual

Null-collision rate [ns-1].

Reimplemented in Garfield::MediumMagboltz, and Garfield::MediumSilicon.

Definition at line 480 of file Medium.cc.

480 {
481 if (m_debug) PrintNotImplemented(m_className, "GetElectronNullCollisionRate");
482 return 0.;
483}

◆ GetElectronTownsend()

bool Garfield::Medium::GetElectronTownsend ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  alpha 
)
inline

Get an entry in the table of Townsend coefficients.

Definition at line 268 of file Medium.hh.

269 {
270 return GetEntry(ie, ib, ia, "GetElectronTownsend", m_eAlp, alpha);
271 }

◆ GetElectronTransverseDiffusion()

bool Garfield::Medium::GetElectronTransverseDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  dt 
)
inline

Get an entry in the table of transverse diffusion coefficients.

Definition at line 256 of file Medium.hh.

258 {
259 return GetEntry(ie, ib, ia, "GetElectronTransverseDiffusion",
260 m_eDifT, dt);
261 }

◆ GetElectronVelocityB()

bool Garfield::Medium::GetElectronVelocityB ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  v 
)
inline

Get an entry in the table of drift speeds along Btrans.

Definition at line 228 of file Medium.hh.

229 {
230 return GetEntry(ie, ib, ia, "GetElectronVelocityB", m_eVelB, v);
231 }

◆ GetElectronVelocityE()

bool Garfield::Medium::GetElectronVelocityE ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  v 
)
inline

Get an entry in the table of drift speeds along E.

Definition at line 208 of file Medium.hh.

209 {
210 return GetEntry(ie, ib, ia, "GetElectronVelocityE", m_eVelE, v);
211 }

◆ GetElectronVelocityExB()

bool Garfield::Medium::GetElectronVelocityExB ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  v 
)
inline

Get an entry in the table of drift speeds along ExB.

Definition at line 218 of file Medium.hh.

219 {
220 return GetEntry(ie, ib, ia, "GetElectronVelocityExB", m_eVelX, v);
221 }

◆ GetEntry()

bool Garfield::Medium::GetEntry ( const unsigned int  i,
const unsigned int  j,
const unsigned int  k,
const std::string &  fcn,
const std::vector< std::vector< std::vector< double > > > &  tab,
double &  val 
) const
protected

Definition at line 888 of file Medium.cc.

891 {
892 val = 0.;
893 if (i >= m_eFields.size() || j >= m_bFields.size() || k >= m_bAngles.size()) {
894 PrintOutOfRange(m_className, fcn, i, j, k);
895 return false;
896 }
897 if (tab.empty()) {
898 if (m_debug) PrintDataNotAvailable(m_className, fcn);
899 return false;
900 }
901 val = tab[k][j][i];
902 return true;
903}

Referenced by GetElectronAttachment(), GetElectronLongitudinalDiffusion(), GetElectronLorentzAngle(), GetElectronTownsend(), GetElectronTransverseDiffusion(), GetElectronVelocityB(), GetElectronVelocityE(), GetElectronVelocityExB(), GetHoleAttachment(), GetHoleLongitudinalDiffusion(), GetHoleTownsend(), GetHoleTransverseDiffusion(), GetHoleVelocityB(), GetHoleVelocityE(), GetHoleVelocityExB(), GetIonDissociation(), GetIonLongitudinalDiffusion(), GetIonMobility(), and GetIonTransverseDiffusion().

◆ GetExtrapolationIndex()

bool Garfield::Medium::GetExtrapolationIndex ( std::string  str,
unsigned int &  nb 
) const
protected

Definition at line 1127 of file Medium.cc.

1127 {
1128 // Convert to upper-case
1129 for (unsigned int i = 0; i < str.length(); ++i) {
1130 str[i] = toupper(str[i]);
1131 }
1132
1133 if (str == "CONST" || str == "CONSTANT") {
1134 nb = 0;
1135 } else if (str == "LIN" || str == "LINEAR") {
1136 nb = 1;
1137 } else if (str == "EXP" || str == "EXPONENTIAL") {
1138 nb = 2;
1139 } else {
1140 return false;
1141 }
1142
1143 return true;
1144}

Referenced by SetExtrapolationMethod().

◆ GetFanoFactor()

double Garfield::Medium::GetFanoFactor ( )
inline

Get the Fano factor.

Definition at line 87 of file Medium.hh.

87{ return m_fano; }
double m_fano
Definition: Medium.hh:537

◆ GetFieldGrid()

void Garfield::Medium::GetFieldGrid ( std::vector< double > &  efields,
std::vector< double > &  bfields,
std::vector< double > &  angles 
)

Get the fields and E-B angles used in the transport tables.

Definition at line 864 of file Medium.cc.

866 {
867 efields = m_eFields;
868 bfields = m_bFields;
869 angles = m_bAngles;
870}

◆ GetHoleAttachment()

bool Garfield::Medium::GetHoleAttachment ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  eta 
)
inline

Get an entry in the table of attachment coefficients.

Definition at line 364 of file Medium.hh.

365 {
366 return GetEntry(ie, ib, ia, "GetHoleAttachment", m_hAtt, eta);
367 }
std::vector< std::vector< std::vector< double > > > m_hAtt
Definition: Medium.hh:571

◆ GetHoleLongitudinalDiffusion()

bool Garfield::Medium::GetHoleLongitudinalDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  dl 
)
inline

Get an entry in the table of longitudinal diffusion coefficients.

Definition at line 332 of file Medium.hh.

334 {
335 return GetEntry(ie, ib, ia, "GetHoleLongitudinalDiffusion", m_hDifL, dl);
336 }
std::vector< std::vector< std::vector< double > > > m_hDifL
Definition: Medium.hh:568

◆ GetHoleTownsend()

bool Garfield::Medium::GetHoleTownsend ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  alpha 
)
inline

Get an entry in the table of Townsend coefficients.

Definition at line 354 of file Medium.hh.

355 {
356 return GetEntry(ie, ib, ia, "GetHoleTownsend", m_hAlp, alpha);
357 }
std::vector< std::vector< std::vector< double > > > m_hAlp
Definition: Medium.hh:570

◆ GetHoleTransverseDiffusion()

bool Garfield::Medium::GetHoleTransverseDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  dt 
)
inline

Get an entry in the table of transverse diffusion coefficients.

Definition at line 344 of file Medium.hh.

345 {
346 return GetEntry(ie, ib, ia, "GetHoleTransverseDiffusion", m_hDifT, dt);
347 }
std::vector< std::vector< std::vector< double > > > m_hDifT
Definition: Medium.hh:569

◆ GetHoleVelocityB()

bool Garfield::Medium::GetHoleVelocityB ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  v 
)
inline

Get an entry in the table of drift speeds along Btrans.

Definition at line 320 of file Medium.hh.

321 {
322 return GetEntry(ie, ib, ia, "GetHoleVelocityB", m_hVelB, v);
323 }
std::vector< std::vector< std::vector< double > > > m_hVelB
Definition: Medium.hh:567

◆ GetHoleVelocityE()

bool Garfield::Medium::GetHoleVelocityE ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  v 
)
inline

Get an entry in the table of drift speeds along E.

Definition at line 300 of file Medium.hh.

301 {
302 return GetEntry(ie, ib, ia, "GetHoleVelocityE", m_hVelE, v);
303 }
std::vector< std::vector< std::vector< double > > > m_hVelE
Definition: Medium.hh:565

◆ GetHoleVelocityExB()

bool Garfield::Medium::GetHoleVelocityExB ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  v 
)
inline

Get an entry in the table of drift speeds along ExB.

Definition at line 310 of file Medium.hh.

311 {
312 return GetEntry(ie, ib, ia, "GetHoleVelocityExB", m_hVelX, v);
313 }
std::vector< std::vector< std::vector< double > > > m_hVelX
Definition: Medium.hh:566

◆ GetId()

int Garfield::Medium::GetId ( ) const
inline

Return the id number of the class instance.

Definition at line 21 of file Medium.hh.

21{ return m_id; }

Referenced by Garfield::GeometrySimple::AddSolid(), and Garfield::ViewGeometry::Plot().

◆ GetIonDissociation()

bool Garfield::Medium::GetIonDissociation ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  diss 
)
inline

Get an entry in the table of dissociation coefficients.

Definition at line 410 of file Medium.hh.

411 {
412 return GetEntry(ie, ib, ia, "GetIonDissociation", m_iDis, diss);
413 }
std::vector< std::vector< std::vector< double > > > m_iDis
Definition: Medium.hh:579

◆ GetIonLongitudinalDiffusion()

bool Garfield::Medium::GetIonLongitudinalDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  dl 
)
inline

Get an entry in the table of longitudinal diffusion coefficients.

Definition at line 390 of file Medium.hh.

391 {
392 return GetEntry(ie, ib, ia, "GetIonLongitudinalDiffusion", m_iDifL, dl);
393 }
std::vector< std::vector< std::vector< double > > > m_iDifL
Definition: Medium.hh:577

◆ GetIonMobility()

bool Garfield::Medium::GetIonMobility ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  mu 
)
inline

Get an entry in the table of ion mobilities.

Definition at line 379 of file Medium.hh.

380 {
381 return GetEntry(ie, ib, ia, "GetIonMobility", m_iMob, mu);
382 }
std::vector< std::vector< std::vector< double > > > m_iMob
Definition: Medium.hh:576

◆ GetIonTransverseDiffusion()

bool Garfield::Medium::GetIonTransverseDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
double &  dt 
)
inline

Get an entry in the table of transverse diffusion coefficients.

Definition at line 400 of file Medium.hh.

401 {
402 return GetEntry(ie, ib, ia, "GetIonTransverseDiffusion", m_iDifT, dt);
403 }
std::vector< std::vector< std::vector< double > > > m_iDifT
Definition: Medium.hh:578

◆ GetMassDensity()

double Garfield::Medium::GetMassDensity ( ) const
virtual

Get the mass density [g/cm3].

Reimplemented in Garfield::MediumGas.

Definition at line 101 of file Medium.cc.

101 {
102 return m_density * AtomicMassUnit * m_a;
103}
double m_density
Definition: Medium.hh:527

Referenced by Garfield::TrackHeed::NewTrack(), Garfield::GeometryRoot::SetMedium(), Garfield::TrackHeed::TransportDeltaElectron(), and Garfield::TrackHeed::TransportPhoton().

◆ GetName()

◆ GetNumberDensity()

virtual double Garfield::Medium::GetNumberDensity ( ) const
inlinevirtual

Get the number density [cm-3].

Reimplemented in Garfield::MediumGas.

Definition at line 60 of file Medium.hh.

60{ return m_density; }

Referenced by Garfield::TrackElectron::GetCluster(), Garfield::TrackPAI::GetCluster(), and Garfield::TrackPAI::NewTrack().

◆ GetNumberOfComponents()

unsigned int Garfield::Medium::GetNumberOfComponents ( ) const
inline

Get number of components of the medium.

Definition at line 45 of file Medium.hh.

45{ return m_nComponents; }

◆ GetNumberOfDeexcitationProducts()

virtual unsigned int Garfield::Medium::GetNumberOfDeexcitationProducts ( ) const
inlinevirtual

Reimplemented in Garfield::MediumMagboltz.

Definition at line 144 of file Medium.hh.

144{ return 0; }

◆ GetOpticalDataRange()

bool Garfield::Medium::GetOpticalDataRange ( double &  emin,
double &  emax,
const unsigned int  i = 0 
)
virtual

Get the energy range [eV] of the available optical data.

Reimplemented in Garfield::MediumSilicon.

Definition at line 624 of file Medium.cc.

625 {
626 if (i >= m_nComponents) {
627 std::cerr << m_className << "::GetOpticalDataRange: Index out of range.\n";
628 return false;
629 }
630
631 if (m_debug) PrintNotImplemented(m_className, "GetOpticalDataRange");
632 emin = emax = 0.;
633 return false;
634}

◆ GetPhotoAbsorptionCrossSection()

bool Garfield::Medium::GetPhotoAbsorptionCrossSection ( const double  e,
double &  sigma,
const unsigned int  i = 0 
)
virtual

Reimplemented in Garfield::MediumGas.

Definition at line 654 of file Medium.cc.

655 {
656 if (i >= m_nComponents) {
657 std::cerr << m_className << "::GetPhotoAbsorptionCrossSection:\n";
658 std::cerr << " Component " << i << " does not exist.\n";
659 return false;
660 }
661
662 if (e < 0.) {
663 std::cerr << m_className << "::GetPhotoAbsorptionCrossSection:\n";
664 std::cerr << " Energy must be > 0.\n";
665 return false;
666 }
667
668 if (m_debug) {
669 PrintNotImplemented(m_className, "GetPhotoAbsorptionCrossSection");
670 }
671 sigma = 0.;
672 return false;
673}

Referenced by GetPhotonCollisionRate().

◆ GetPhotonCollision()

bool Garfield::Medium::GetPhotonCollision ( const double  e,
int &  type,
int &  level,
double &  e1,
double &  ctheta,
int &  nsec,
double &  esec 
)
virtual

Reimplemented in Garfield::MediumMagboltz.

Definition at line 682 of file Medium.cc.

684 {
685 type = level = -1;
686 e1 = e;
687 ctheta = 1.;
688 nsec = 0;
689 esec = 0.;
690 return false;
691}

◆ GetPhotonCollisionRate()

double Garfield::Medium::GetPhotonCollisionRate ( const double  e)
virtual

Reimplemented in Garfield::MediumMagboltz.

Definition at line 675 of file Medium.cc.

675 {
676 double sigma = 0.;
677 if (!GetPhotoAbsorptionCrossSection(e, sigma)) return 0.;
678
679 return sigma * m_density * SpeedOfLight;
680}
virtual bool GetPhotoAbsorptionCrossSection(const double e, double &sigma, const unsigned int i=0)
Definition: Medium.cc:654

◆ GetPressure()

double Garfield::Medium::GetPressure ( ) const
inline

Definition at line 38 of file Medium.hh.

38{ return m_pressure; }
double m_pressure
Definition: Medium.hh:517

◆ GetTemperature()

double Garfield::Medium::GetTemperature ( ) const
inline

Get the temperature [K].

Definition at line 34 of file Medium.hh.

34{ return m_temperature; }

◆ GetW()

double Garfield::Medium::GetW ( )
inline

Get the W value.

Definition at line 83 of file Medium.hh.

83{ return m_w; }

◆ HoleAttachment()

bool Garfield::Medium::HoleAttachment ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  eta 
)
virtual

Attachment coefficient [cm-1].

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumSilicon.

Definition at line 547 of file Medium.cc.

549 {
550
551 if (!Alpha(ex, ey, ez, bx, by, bz, m_hAtt, m_intpAtt, m_hThrAtt, m_extrAtt,
552 eta)) {
553 return false;
554 }
555 // Apply scaling.
556 eta = ScaleAttachment(eta);
557 return true;
558}
unsigned int m_hThrAtt
Definition: Medium.hh:585

Referenced by Garfield::ViewMedium::EvaluateFunction(), Garfield::MediumCdTe::HoleAttachment(), Garfield::MediumGaAs::HoleAttachment(), and Garfield::MediumSilicon::HoleAttachment().

◆ HoleDiffusion() [1/2]

bool Garfield::Medium::HoleDiffusion ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  dl,
double &  dt 
)
virtual

Longitudinal and transverse diffusion coefficients [cm1/2].

Definition at line 521 of file Medium.cc.

523 {
524 return Diffusion(ex, ey, ez, bx, by, bz, m_hDifL, m_hDifT, dl, dt);
525}

Referenced by Garfield::ViewMedium::EvaluateFunction().

◆ HoleDiffusion() [2/2]

bool Garfield::Medium::HoleDiffusion ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double  cov[3][3] 
)
virtual

Diffusion tensor.

Definition at line 527 of file Medium.cc.

529 {
530
531 return Diffusion(ex, ey, ez, bx, by, bz, m_hDifM, cov);
532}
std::vector< std::vector< std::vector< std::vector< double > > > > m_hDifM
Definition: Medium.hh:573

◆ HoleMobility()

double Garfield::Medium::HoleMobility ( )
virtual

Low-field mobility [cm2 V-1 ns-1].

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumSilicon.

Definition at line 560 of file Medium.cc.

560 {
561 if (m_hVelE.empty()) return -1.;
562 return m_hVelE[0][0][0] / UnScaleElectricField(m_eFields[0]);
563}

◆ HoleTownsend()

bool Garfield::Medium::HoleTownsend ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  alpha 
)
virtual

Ionisation coefficient [cm-1].

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumSilicon.

Definition at line 534 of file Medium.cc.

536 {
537
538 if (!Alpha(ex, ey, ez, bx, by, bz, m_hAlp, m_intpAlp, m_hThrAlp, m_extrAlp,
539 alpha)) {
540 return false;
541 }
542 // Apply scaling.
543 alpha = ScaleTownsend(alpha);
544 return true;
545}
unsigned int m_hThrAlp
Definition: Medium.hh:584

Referenced by Garfield::ViewMedium::EvaluateFunction(), Garfield::MediumCdTe::HoleTownsend(), Garfield::MediumGaAs::HoleTownsend(), and Garfield::MediumSilicon::HoleTownsend().

◆ HoleVelocity()

bool Garfield::Medium::HoleVelocity ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  vx,
double &  vy,
double &  vz 
)
virtual

Drift velocity [cm / ns].

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumSilicon.

Definition at line 513 of file Medium.cc.

515 {
516
517 return Velocity(ex, ey, ez, bx, by, bz, m_hVelE, m_hVelB, m_hVelX, +1.,
518 vx, vy, vz);
519}

Referenced by Garfield::ViewMedium::EvaluateFunction(), Garfield::MediumCdTe::HoleVelocity(), Garfield::MediumGaAs::HoleVelocity(), and Garfield::MediumSilicon::HoleVelocity().

◆ Init() [1/2]

void Garfield::Medium::Init ( const size_t  nE,
const size_t  nB,
const size_t  nA,
const size_t  nT,
std::vector< std::vector< std::vector< std::vector< double > > > > &  tab,
const double  val 
)
protected

Definition at line 1315 of file Medium.cc.

1318 {
1319 if (nE == 0 || nB == 0 || nA == 0 || nT == 0) {
1320 std::cerr << m_className << "::Init: Invalid grid.\n";
1321 return;
1322 }
1323
1324 tab.assign(nT, std::vector<std::vector<std::vector<double> > >(
1325 nA, std::vector<std::vector<double> >(
1326 nB, std::vector<double>(nE, val))));
1327}

◆ Init() [2/2]

void Garfield::Medium::Init ( const size_t  nE,
const size_t  nB,
const size_t  nA,
std::vector< std::vector< std::vector< double > > > &  tab,
const double  val 
)
protected

Definition at line 1304 of file Medium.cc.

1306 {
1307 if (nE == 0 || nB == 0 || nA == 0) {
1308 std::cerr << m_className << "::Init: Invalid grid.\n";
1309 return;
1310 }
1311 tab.assign(
1312 nA, std::vector<std::vector<double> >(nB, std::vector<double>(nE, val)));
1313}

Referenced by Clone(), Garfield::MediumMagboltz::GenerateGasTable(), Garfield::MediumGas::LoadGasFile(), Garfield::MediumGas::MergeGasFile(), SetEntry(), and SetIonMobility().

◆ Interpolate()

bool Garfield::Medium::Interpolate ( const double  e,
const double  b,
const double  a,
const std::vector< std::vector< std::vector< double > > > &  table,
double &  y,
const unsigned int  intp,
const std::pair< unsigned int, unsigned int > &  extr 
) const
protected

Definition at line 1210 of file Medium.cc.

1214 {
1215 if (table.empty()) {
1216 y = 0.;
1217 return false; // TODO: true!
1218 }
1219
1220 if (m_tab2d) {
1222 m_bAngles.size(), m_bFields.size(), m_eFields.size(),
1223 a, b, e, y, intp);
1224 } else {
1225 y = Interpolate1D(e, table[0][0], m_eFields, intp, extr);
1226 }
1227 return true;
1228}
double Interpolate1D(const double e, const std::vector< double > &table, const std::vector< double > &fields, const unsigned int intpMeth, const std::pair< unsigned int, unsigned int > &extr) const
Definition: Medium.cc:1230
bool Boxin3(const std::vector< std::vector< std::vector< double > > > &value, const std::vector< double > &xAxis, const std::vector< double > &yAxis, const std::vector< double > &zAxis, const int nx, const int ny, const int nz, const double xx, const double yy, const double zz, double &f, const int iOrder)
Definition: Numerics.cc:1495

Referenced by Alpha(), Clone(), Diffusion(), ElectronLorentzAngle(), IonVelocity(), and Velocity().

◆ Interpolate1D()

double Garfield::Medium::Interpolate1D ( const double  e,
const std::vector< double > &  table,
const std::vector< double > &  fields,
const unsigned int  intpMeth,
const std::pair< unsigned int, unsigned int > &  extr 
) const
protected

Definition at line 1230 of file Medium.cc.

1233 {
1234 // This function is a generalized version of the Fortran functions
1235 // GASVEL, GASVT1, GASVT2, GASLOR, GASMOB, GASDFT, and GASDFL
1236 // for the case of a 1D table. All variables are generic.
1237
1238 const int nSizeTable = fields.size();
1239
1240 if (e < 0. || nSizeTable < 1) return 0.;
1241
1242 double result = 0.;
1243
1244 if (nSizeTable == 1) {
1245 // Only one point
1246 result = table[0];
1247 } else if (e < fields[0]) {
1248 // Extrapolation towards small fields
1249 if (fields[0] >= fields[1]) {
1250 if (m_debug) {
1251 std::cerr << m_className << "::Interpolate1D:\n";
1252 std::cerr << " First two field values coincide.\n";
1253 std::cerr << " No extrapolation to lower fields.\n";
1254 }
1255 result = table[0];
1256 } else if (extr.first == 1) {
1257 // Linear extrapolation
1258 const double extr4 = (table[1] - table[0]) / (fields[1] - fields[0]);
1259 const double extr3 = table[0] - extr4 * fields[0];
1260 result = extr3 + extr4 * e;
1261 } else if (extr.first == 2) {
1262 // Logarithmic extrapolation
1263 const double extr4 = log(table[1] / table[0]) / (fields[1] - fields[0]);
1264 const double extr3 = log(table[0] - extr4 * fields[0]);
1265 result = std::exp(std::min(50., extr3 + extr4 * e));
1266 } else {
1267 result = table[0];
1268 }
1269 } else if (e > fields[nSizeTable - 1]) {
1270 // Extrapolation towards large fields
1271 if (fields[nSizeTable - 1] <= fields[nSizeTable - 2]) {
1272 if (m_debug) {
1273 std::cerr << m_className << "::Interpolate1D:\n";
1274 std::cerr << " Last two field values coincide.\n";
1275 std::cerr << " No extrapolation to higher fields.\n";
1276 }
1277 result = table[nSizeTable - 1];
1278 } else if (extr.second == 1) {
1279 // Linear extrapolation
1280 const double extr2 = (table[nSizeTable - 1] - table[nSizeTable - 2]) /
1281 (fields[nSizeTable - 1] - fields[nSizeTable - 2]);
1282 const double extr1 =
1283 table[nSizeTable - 1] - extr2 * fields[nSizeTable - 1];
1284 result = extr1 + extr2 * e;
1285 } else if (extr.second == 2) {
1286 // Logarithmic extrapolation
1287 const double extr2 = log(table[nSizeTable - 1] / table[nSizeTable - 2]) /
1288 (fields[nSizeTable - 1] - fields[nSizeTable - 2]);
1289 const double extr1 =
1290 log(table[nSizeTable - 1]) - extr2 * fields[nSizeTable - 1];
1291 result = exp(std::min(50., extr1 + extr2 * e));
1292 } else {
1293 result = table[nSizeTable - 1];
1294 }
1295 } else {
1296 // Intermediate points, spline interpolation (not implemented).
1297 // Intermediate points, Newtonian interpolation
1298 result = Numerics::Divdif(table, fields, nSizeTable, e, intpMeth);
1299 }
1300
1301 return result;
1302}
double Divdif(const std::vector< double > &f, const std::vector< double > &a, int nn, double x, int mm)
Definition: Numerics.cc:1206

Referenced by Interpolate(), and SetIonMobility().

◆ IonDiffusion()

bool Garfield::Medium::IonDiffusion ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  dl,
double &  dt 
)
virtual

Longitudinal and transverse diffusion coefficients [cm1/2].

Definition at line 600 of file Medium.cc.

602 {
603
604 return Diffusion(ex, ey, ez, bx, by, bz, m_iDifL, m_iDifT, dl, dt);
605}

Referenced by Garfield::ViewMedium::EvaluateFunction().

◆ IonDissociation()

bool Garfield::Medium::IonDissociation ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  diss 
)
virtual

Dissociation coefficient.

Definition at line 607 of file Medium.cc.

609 {
610
611 if (!Alpha(ex, ey, ez, bx, by, bz, m_iDis, m_intpDis, m_iThrDis, m_extrDis,
612 diss)) {
613 return false;
614 }
615 // Apply scaling.
616 diss = ScaleDissociation(diss);
617 return true;
618}
unsigned int m_intpDis
Definition: Medium.hh:604
virtual double ScaleDissociation(const double diss) const
Definition: Medium.hh:484
std::pair< unsigned int, unsigned int > m_extrDis
Definition: Medium.hh:595
unsigned int m_iThrDis
Definition: Medium.hh:586

◆ IonMobility()

double Garfield::Medium::IonMobility ( )
virtual

Low-field mobility [cm2 V-1 ns-1].

Definition at line 620 of file Medium.cc.

620 {
621 return m_iMob.empty() ? -1. : m_iMob[0][0][0];
622}

◆ IonVelocity()

bool Garfield::Medium::IonVelocity ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
double &  vx,
double &  vy,
double &  vz 
)
virtual

Drift velocity [cm / ns].

Definition at line 565 of file Medium.cc.

567 {
568 vx = vy = vz = 0.;
569 if (m_iMob.empty()) return false;
570 // Compute the magnitude of the electric field.
571 const double e = sqrt(ex * ex + ey * ey + ez * ez);
572 const double e0 = ScaleElectricField(e);
573 if (e < Small || e0 < Small) return true;
574 // Compute the magnitude of the electric field.
575 const double b = sqrt(bx * bx + by * by + bz * bz);
576
577 // Compute the angle between B field and E field.
578 const double ebang = m_tab2d ? GetAngle(ex, ey, ez, bx, by, bz, e, b) : 0.;
579 double mu = 0.;
580 if (!Interpolate(e0, b, ebang, m_iMob, mu, m_intpMob, m_extrMob)) mu = 0.;
581
582 constexpr double q = 1.;
583 mu *= q;
584 if (b < Small) {
585 vx = mu * ex;
586 vy = mu * ey;
587 vz = mu * ez;
588 } else {
589 const double eb = bx * ex + by * ey + bz * ez;
590 const double mu2 = mu * mu;
591 const double f = mu / (1. + mu2 * b * b);
592 vx = f * (ex + mu * (ey * bz - ez * by) + mu2 * bx * eb);
593 vy = f * (ey + mu * (ez * bx - ex * bz) + mu2 * by * eb);
594 vz = f * (ez + mu * (ex * by - ey * bx) + mu2 * bz * eb);
595 }
596
597 return true;
598}
unsigned int m_intpMob
Definition: Medium.hh:603
std::pair< unsigned int, unsigned int > m_extrMob
Definition: Medium.hh:594

Referenced by Garfield::ViewMedium::EvaluateFunction().

◆ IsConductor()

virtual bool Garfield::Medium::IsConductor ( ) const
inlinevirtual

Is this medium a conductor?

Reimplemented in Garfield::MediumConductor.

Definition at line 29 of file Medium.hh.

29{ return false; }

Referenced by Garfield::ComponentNeBem2d::ElectricField().

◆ IsDriftable()

◆ IsGas()

virtual bool Garfield::Medium::IsGas ( ) const
inlinevirtual

Is this medium a gas?

Reimplemented in Garfield::MediumGas.

Definition at line 25 of file Medium.hh.

25{ return false; }

Referenced by Garfield::TrackElectron::NewTrack(), and Garfield::ViewGeometry::Plot().

◆ IsIonisable()

◆ IsMicroscopic()

bool Garfield::Medium::IsMicroscopic ( ) const
inline

Does the medium have electron scattering rates?

Definition at line 76 of file Medium.hh.

76{ return m_microscopic; }
bool m_microscopic
Definition: Medium.hh:531

◆ IsSemiconductor()

virtual bool Garfield::Medium::IsSemiconductor ( ) const
inlinevirtual

Is this medium a semiconductor?

Reimplemented in Garfield::MediumCdTe, Garfield::MediumGaAs, and Garfield::MediumSilicon.

Definition at line 27 of file Medium.hh.

27{ return false; }

Referenced by Garfield::ViewGeometry::Plot().

◆ ResetElectronAttachment()

void Garfield::Medium::ResetElectronAttachment ( )
inline

Definition at line 429 of file Medium.hh.

429{ m_eAtt.clear(); }

Referenced by ResetTables().

◆ ResetElectronDiffusion()

void Garfield::Medium::ResetElectronDiffusion ( )
inline

Definition at line 423 of file Medium.hh.

423 {
424 m_eDifL.clear();
425 m_eDifT.clear();
426 m_eDifM.clear();
427 }

Referenced by ResetTables().

◆ ResetElectronLorentzAngle()

void Garfield::Medium::ResetElectronLorentzAngle ( )
inline

Definition at line 430 of file Medium.hh.

430{ m_eLor.clear(); }

Referenced by ResetTables().

◆ ResetElectronTownsend()

void Garfield::Medium::ResetElectronTownsend ( )
inline

Definition at line 428 of file Medium.hh.

428{ m_eAlp.clear(); }

Referenced by ResetTables().

◆ ResetElectronVelocity()

void Garfield::Medium::ResetElectronVelocity ( )
inline

Definition at line 418 of file Medium.hh.

418 {
419 m_eVelE.clear();
420 m_eVelB.clear();
421 m_eVelX.clear();
422 }

Referenced by ResetTables().

◆ ResetHoleAttachment()

void Garfield::Medium::ResetHoleAttachment ( )
inline

Definition at line 443 of file Medium.hh.

443{ m_hAtt.clear(); }

Referenced by ResetTables().

◆ ResetHoleDiffusion()

void Garfield::Medium::ResetHoleDiffusion ( )
inline

Definition at line 437 of file Medium.hh.

437 {
438 m_hDifL.clear();
439 m_hDifT.clear();
440 m_hDifM.clear();
441 }

Referenced by ResetTables().

◆ ResetHoleTownsend()

void Garfield::Medium::ResetHoleTownsend ( )
inline

Definition at line 442 of file Medium.hh.

442{ m_hAlp.clear(); }

Referenced by ResetTables().

◆ ResetHoleVelocity()

void Garfield::Medium::ResetHoleVelocity ( )
inline

Definition at line 432 of file Medium.hh.

432 {
433 m_hVelE.clear();
434 m_hVelB.clear();
435 m_hVelX.clear();
436 }

Referenced by ResetTables().

◆ ResetIonDiffusion()

void Garfield::Medium::ResetIonDiffusion ( )
inline

Definition at line 446 of file Medium.hh.

446 {
447 m_iDifL.clear();
448 m_iDifT.clear();
449 }

Referenced by ResetTables().

◆ ResetIonDissociation()

void Garfield::Medium::ResetIonDissociation ( )
inline

Definition at line 450 of file Medium.hh.

450{ m_iDis.clear(); }

Referenced by ResetTables().

◆ ResetIonMobility()

void Garfield::Medium::ResetIonMobility ( )
inline

Definition at line 445 of file Medium.hh.

445{ m_iMob.clear(); }

Referenced by ResetTables(), and SetIonMobility().

◆ ResetTables()

void Garfield::Medium::ResetTables ( )
virtual

Reset all tables of transport parameters.

Reimplemented in Garfield::MediumGas.

Definition at line 905 of file Medium.cc.

905 {
911
916
920}
void ResetHoleAttachment()
Definition: Medium.hh:443
void ResetIonMobility()
Definition: Medium.hh:445
void ResetElectronAttachment()
Definition: Medium.hh:429
void ResetIonDiffusion()
Definition: Medium.hh:446
void ResetElectronLorentzAngle()
Definition: Medium.hh:430
void ResetHoleVelocity()
Definition: Medium.hh:432
void ResetIonDissociation()
Definition: Medium.hh:450
void ResetElectronDiffusion()
Definition: Medium.hh:423
void ResetHoleDiffusion()
Definition: Medium.hh:437
void ResetElectronTownsend()
Definition: Medium.hh:428
void ResetHoleTownsend()
Definition: Medium.hh:442
void ResetElectronVelocity()
Definition: Medium.hh:418

Referenced by Garfield::MediumGas::ResetTables().

◆ ScaleAttachment()

virtual double Garfield::Medium::ScaleAttachment ( const double  eta) const
inlinevirtual

Reimplemented in Garfield::MediumGas.

Definition at line 482 of file Medium.hh.

482{ return eta; }

Referenced by ElectronAttachment(), and HoleAttachment().

◆ ScaleDiffusion()

virtual double Garfield::Medium::ScaleDiffusion ( const double  d) const
inlinevirtual

Reimplemented in Garfield::MediumGas.

Definition at line 479 of file Medium.hh.

479{ return d; }

Referenced by Diffusion().

◆ ScaleDiffusionTensor()

virtual double Garfield::Medium::ScaleDiffusionTensor ( const double  d) const
inlinevirtual

Reimplemented in Garfield::MediumGas.

Definition at line 480 of file Medium.hh.

480{ return d; }

Referenced by Diffusion().

◆ ScaleDissociation()

virtual double Garfield::Medium::ScaleDissociation ( const double  diss) const
inlinevirtual

Definition at line 484 of file Medium.hh.

484{ return diss; }

Referenced by IonDissociation().

◆ ScaleElectricField()

virtual double Garfield::Medium::ScaleElectricField ( const double  e) const
inlinevirtual

Reimplemented in Garfield::MediumGas.

Definition at line 476 of file Medium.hh.

476{ return e; }

Referenced by Alpha(), Diffusion(), ElectronLorentzAngle(), IonVelocity(), and Velocity().

◆ ScaleLorentzAngle()

virtual double Garfield::Medium::ScaleLorentzAngle ( const double  lor) const
inlinevirtual

Reimplemented in Garfield::MediumGas.

Definition at line 483 of file Medium.hh.

483{ return lor; }

Referenced by ElectronLorentzAngle().

◆ ScaleTownsend()

virtual double Garfield::Medium::ScaleTownsend ( const double  alpha) const
inlinevirtual

Reimplemented in Garfield::MediumGas.

Definition at line 481 of file Medium.hh.

481{ return alpha; }

Referenced by ElectronTownsend(), and HoleTownsend().

◆ ScaleVelocity()

virtual double Garfield::Medium::ScaleVelocity ( const double  v) const
inlinevirtual

Definition at line 478 of file Medium.hh.

478{ return v; }

◆ SetAtomicNumber()

void Garfield::Medium::SetAtomicNumber ( const double  z)
virtual

Set the effective atomic number.

Reimplemented in Garfield::MediumGas.

Definition at line 114 of file Medium.cc.

114 {
115 if (z < 1.) {
116 std::cerr << m_className << "::SetAtomicNumber:\n"
117 << " Atomic number must be >= 1.\n";
118 return;
119 }
120 m_z = z;
121 m_isChanged = true;
122}
bool m_isChanged
Definition: Medium.hh:540

Referenced by Garfield::MediumCdTe::MediumCdTe(), Garfield::MediumGaAs::MediumGaAs(), and Garfield::MediumSilicon::MediumSilicon().

◆ SetAtomicWeight()

void Garfield::Medium::SetAtomicWeight ( const double  a)
virtual

Set the effective atomic weight.

Reimplemented in Garfield::MediumGas.

Definition at line 124 of file Medium.cc.

124 {
125 if (a <= 0.) {
126 std::cerr << m_className << "::SetAtomicWeight:\n"
127 << " Atomic weight must be greater than zero.\n";
128 return;
129 }
130 m_a = a;
131 m_isChanged = true;
132}

Referenced by Garfield::MediumCdTe::MediumCdTe(), Garfield::MediumGaAs::MediumGaAs(), and Garfield::MediumSilicon::MediumSilicon().

◆ SetDielectricConstant()

void Garfield::Medium::SetDielectricConstant ( const double  eps)

Set the relative static dielectric constant.

Definition at line 91 of file Medium.cc.

91 {
92 if (eps < 1.) {
93 std::cerr << m_className << "::SetDielectricConstant:\n"
94 << " Dielectric constant must be >= 1.\n";
95 return;
96 }
97 m_epsilon = eps;
98 m_isChanged = true;
99}

Referenced by Garfield::MediumCdTe::MediumCdTe(), Garfield::MediumGaAs::MediumGaAs(), and Garfield::MediumSilicon::MediumSilicon().

◆ SetElectronAttachment()

bool Garfield::Medium::SetElectronAttachment ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  eta 
)
inline

Set an entry in the table of attachment coefficients.

Definition at line 273 of file Medium.hh.

274 {
275 return SetEntry(ie, ib, ia, "SetElectronAttachment", m_eAtt, eta);
276 }
bool SetEntry(const unsigned int i, const unsigned int j, const unsigned int k, const std::string &fcn, std::vector< std::vector< std::vector< double > > > &tab, const double val)
Definition: Medium.cc:872

◆ SetElectronLongitudinalDiffusion()

bool Garfield::Medium::SetElectronLongitudinalDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  dl 
)
inline

Set an entry in the table of longitudinal diffusion coefficients.

Definition at line 233 of file Medium.hh.

236 {
237 return SetEntry(ie, ib, ia, "SetElectronLongitudinalDiffusion",
238 m_eDifL, dl);
239 }

◆ SetElectronLorentzAngle()

bool Garfield::Medium::SetElectronLorentzAngle ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  lor 
)
inline

Set an entry in the table of Lorentz angles.

Definition at line 284 of file Medium.hh.

285 {
286 return SetEntry(ie, ib, ia, "SetElectronLorentzAngle", m_eLor, lor);
287 }

◆ SetElectronTownsend()

bool Garfield::Medium::SetElectronTownsend ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  alpha 
)
inline

Set an entry in the table of Townsend coefficients.

Definition at line 263 of file Medium.hh.

264 {
265 return SetEntry(ie, ib, ia, "SetElectronTownsend", m_eAlp, alpha);
266 }

◆ SetElectronTransverseDiffusion()

bool Garfield::Medium::SetElectronTransverseDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  dt 
)
inline

Set an entry in the table of transverse diffusion coefficients.

Definition at line 248 of file Medium.hh.

251 {
252 return SetEntry(ie, ib, ia, "SetElectronTransverseDiffusion",
253 m_eDifT, dt);
254 }

◆ SetElectronVelocityB()

bool Garfield::Medium::SetElectronVelocityB ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  v 
)
inline

Set an entry in the table of drift speeds along Btrans.

Definition at line 223 of file Medium.hh.

224 {
225 return SetEntry(ie, ib, ia, "SetElectronVelocityB", m_eVelB, v);
226 }

◆ SetElectronVelocityE()

bool Garfield::Medium::SetElectronVelocityE ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  v 
)
inline

Set an entry in the table of drift speeds along E.

Definition at line 203 of file Medium.hh.

204 {
205 return SetEntry(ie, ib, ia, "SetElectronVelocityE", m_eVelE, v);
206 }

◆ SetElectronVelocityExB()

bool Garfield::Medium::SetElectronVelocityExB ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  v 
)
inline

Set an entry in the table of drift speeds along ExB.

Definition at line 213 of file Medium.hh.

214 {
215 return SetEntry(ie, ib, ia, "SetElectronVelocityExB", m_eVelX, v);
216 }

◆ SetEntry()

bool Garfield::Medium::SetEntry ( const unsigned int  i,
const unsigned int  j,
const unsigned int  k,
const std::string &  fcn,
std::vector< std::vector< std::vector< double > > > &  tab,
const double  val 
)
protected

◆ SetExtrapolationMethod()

void Garfield::Medium::SetExtrapolationMethod ( const std::string &  low,
const std::string &  high,
std::pair< unsigned int, unsigned int > &  extr,
const std::string &  fcn 
)
protected

Definition at line 1107 of file Medium.cc.

1110 {
1111 unsigned int i = 0;
1112 if (GetExtrapolationIndex(low, i)) {
1113 extr.first = i;
1114 } else {
1115 std::cerr << m_className << "::SetExtrapolationMethod" << fcn << ":\n"
1116 << " Unknown extrapolation method (" << low << ")\n";
1117 }
1118 unsigned int j = 0;
1119 if (GetExtrapolationIndex(high, j)) {
1120 extr.second = j;
1121 } else {
1122 std::cerr << m_className << "::SetExtrapolationMethod" << fcn << ":\n"
1123 << " Unknown extrapolation method (" << high << ")\n";
1124 }
1125}
bool GetExtrapolationIndex(std::string str, unsigned int &nb) const
Definition: Medium.cc:1127

Referenced by SetExtrapolationMethodAttachment(), SetExtrapolationMethodDiffusion(), Garfield::MediumGas::SetExtrapolationMethodExcitationRates(), SetExtrapolationMethodIonDissociation(), Garfield::MediumGas::SetExtrapolationMethodIonisationRates(), SetExtrapolationMethodIonMobility(), SetExtrapolationMethodTownsend(), and SetExtrapolationMethodVelocity().

◆ SetExtrapolationMethodAttachment()

void Garfield::Medium::SetExtrapolationMethodAttachment ( const std::string &  extrLow,
const std::string &  extrHigh 
)

Definition at line 1092 of file Medium.cc.

1093 {
1094 SetExtrapolationMethod(low, high, m_extrAtt, "Attachment");
1095}
void SetExtrapolationMethod(const std::string &low, const std::string &high, std::pair< unsigned int, unsigned int > &extr, const std::string &fcn)
Definition: Medium.cc:1107

◆ SetExtrapolationMethodDiffusion()

void Garfield::Medium::SetExtrapolationMethodDiffusion ( const std::string &  extrLow,
const std::string &  extrHigh 
)

Definition at line 1082 of file Medium.cc.

1083 {
1084 SetExtrapolationMethod(low, high, m_extrDif, "Diffusion");
1085}

◆ SetExtrapolationMethodIonDissociation()

void Garfield::Medium::SetExtrapolationMethodIonDissociation ( const std::string &  extrLow,
const std::string &  extrHigh 
)

Definition at line 1102 of file Medium.cc.

1103 {
1104 SetExtrapolationMethod(low, high, m_extrDis, "IonDissociation");
1105}

◆ SetExtrapolationMethodIonMobility()

void Garfield::Medium::SetExtrapolationMethodIonMobility ( const std::string &  extrLow,
const std::string &  extrHigh 
)

Definition at line 1097 of file Medium.cc.

1098 {
1099 SetExtrapolationMethod(low, high, m_extrMob, "IonMobility");
1100}

◆ SetExtrapolationMethodTownsend()

void Garfield::Medium::SetExtrapolationMethodTownsend ( const std::string &  extrLow,
const std::string &  extrHigh 
)

Definition at line 1087 of file Medium.cc.

1088 {
1089 SetExtrapolationMethod(low, high, m_extrAlp, "Townsend");
1090}

◆ SetExtrapolationMethodVelocity()

void Garfield::Medium::SetExtrapolationMethodVelocity ( const std::string &  extrLow,
const std::string &  extrHigh 
)

Select the extrapolation method for fields below/above the table range. Possible options are "constant", "linear", and "exponential".

Definition at line 1077 of file Medium.cc.

1078 {
1079 SetExtrapolationMethod(low, high, m_extrVel, "Velocity");
1080}
std::pair< unsigned int, unsigned int > m_extrVel
Definition: Medium.hh:589

◆ SetFanoFactor()

void Garfield::Medium::SetFanoFactor ( const double  f)
inline

Set the Fano factor.

Definition at line 85 of file Medium.hh.

85{ m_fano = f; }

◆ SetFieldGrid() [1/2]

void Garfield::Medium::SetFieldGrid ( const std::vector< double > &  efields,
const std::vector< double > &  bfields,
const std::vector< double > &  angles 
)

Set the fields and E-B angles to be used in the transport tables.

Definition at line 788 of file Medium.cc.

790 {
791 const std::string hdr = m_className + "::SetFieldGrid";
792 if (!CheckFields(efields, hdr, "E-fields")) return;
793 if (!CheckFields(bfields, hdr, "B-fields")) return;
794 if (!CheckFields(angles, hdr, "angles")) return;
795
796 if (m_debug) {
797 std::cout << m_className << "::SetFieldGrid:\n E-fields:\n";
798 for (const auto efield : efields) std::cout << " " << efield << "\n";
799 std::cout << " B-fields:\n";
800 for (const auto bfield : bfields) std::cout << " " << bfield << "\n";
801 std::cout << " Angles:\n";
802 for (const auto angle : angles) std::cout << " " << angle << "\n";
803 }
804
805 // Clone the existing tables.
806 // Electrons
807 Clone(m_eVelE, efields, bfields, angles, m_intpVel, m_extrVel, 0.,
808 "electron velocity along E");
809 Clone(m_eVelB, efields, bfields, angles, m_intpVel, m_extrVel, 0.,
810 "electron velocity along Bt");
811 Clone(m_eVelX, efields, bfields, angles, m_intpVel, m_extrVel, 0.,
812 "electron velocity along ExB");
813 Clone(m_eDifL, efields, bfields, angles, m_intpDif, m_extrDif, 0.,
814 "electron longitudinal diffusion");
815 Clone(m_eDifT, efields, bfields, angles, m_intpDif, m_extrDif, 0.,
816 "electron transverse diffusion");
817 Clone(m_eAlp, efields, bfields, angles, m_intpAlp, m_extrAlp, -30.,
818 "electron Townsend coefficient");
819 Clone(m_eAtt, efields, bfields, angles, m_intpAtt, m_extrAtt, -30.,
820 "electron attachment coefficient");
821 Clone(m_eLor, efields, bfields, angles, m_intpLor, m_extrLor, 0.,
822 "electron Lorentz angle");
823 if (!m_eDifM.empty()) {
824 Clone(m_eDifM, 6, efields, bfields, angles, m_intpDif, m_extrDif, 0.,
825 "electron diffusion tensor");
826 }
827
828 // Holes
829 Clone(m_hVelE, efields, bfields, angles, m_intpVel, m_extrVel, 0.,
830 "hole velocity along E");
831 Clone(m_hVelB, efields, bfields, angles, m_intpVel, m_extrVel, 0.,
832 "hole velocity along Bt");
833 Clone(m_hVelX, efields, bfields, angles, m_intpVel, m_extrVel, 0.,
834 "hole velocity along ExB");
835 Clone(m_hDifL, efields, bfields, angles, m_intpDif, m_extrDif, 0.,
836 "hole longitudinal diffusion");
837 Clone(m_hDifT, efields, bfields, angles, m_intpDif, m_extrDif, 0.,
838 "hole transverse diffusion");
839 Clone(m_hAlp, efields, bfields, angles, m_intpAlp, m_extrAlp, -30.,
840 "hole Townsend coefficient");
841 Clone(m_hAtt, efields, bfields, angles, m_intpAtt, m_extrAtt, -30.,
842 "hole attachment coefficient");
843 if (!m_hDifM.empty()) {
844 Clone(m_hDifM, 6, efields, bfields, angles, m_intpDif, m_extrDif, 0.,
845 "hole diffusion tensor");
846 }
847
848 // Ions
849 Clone(m_iMob, efields, bfields, angles, m_intpMob, m_extrMob, 0.,
850 "ion mobility");
851 Clone(m_iDifL, efields, bfields, angles, m_intpDif, m_extrDif, 0.,
852 "ion longitudinal diffusion");
853 Clone(m_iDifT, efields, bfields, angles, m_intpDif, m_extrDif, 0.,
854 "ion transverse diffusion");
855 Clone(m_iDis, efields, bfields, angles, m_intpDis, m_extrDis, -30.,
856 "ion dissociation");
857
858 if (bfields.size() > 1 || angles.size() > 1) m_tab2d = true;
859 m_eFields = efields;
860 m_bFields = bfields;
861 m_bAngles = angles;
862}
unsigned int m_intpVel
Definition: Medium.hh:598
void Clone(std::vector< std::vector< std::vector< double > > > &tab, const std::vector< double > &efields, const std::vector< double > &bfields, const std::vector< double > &angles, const unsigned int intp, const std::pair< unsigned int, unsigned int > &extr, const double init, const std::string &label)
Definition: Medium.cc:922

◆ SetFieldGrid() [2/2]

void Garfield::Medium::SetFieldGrid ( double  emin,
double  emax,
const size_t  ne,
bool  logE,
double  bmin = 0.,
double  bmax = 0.,
const size_t  nb = 1,
double  amin = HalfPi,
double  amax = HalfPi,
const size_t  na = 1 
)

Set the range of fields to be covered by the transport tables.

Definition at line 693 of file Medium.cc.

695 {
696 // Check if the requested E-field range makes sense.
697 if (ne <= 0) {
698 std::cerr << m_className << "::SetFieldGrid:\n"
699 << " Number of E-fields must be > 0.\n";
700 return;
701 }
702
703 if (emin < 0. || emax < 0.) {
704 std::cerr << m_className << "::SetFieldGrid:\n"
705 << " Electric fields must be positive.\n";
706 return;
707 }
708
709 if (emax < emin) {
710 std::cerr << m_className << "::SetFieldGrid: Swapping min./max. E-field.\n";
711 std::swap(emin, emax);
712 }
713
714 double estep = 0.;
715 if (logE) {
716 // Logarithmic scale
717 if (emin < Small) {
718 std::cerr << m_className << "::SetFieldGrid:\n"
719 << " Min. E-field must be non-zero for log. scale.\n";
720 return;
721 }
722 if (ne == 1) {
723 std::cerr << m_className << "::SetFieldGrid:\n"
724 << " Number of E-fields must be > 1 for log. scale.\n";
725 return;
726 }
727 estep = pow(emax / emin, 1. / (ne - 1.));
728 } else {
729 // Linear scale
730 if (ne > 1) estep = (emax - emin) / (ne - 1.);
731 }
732
733 // Check if the requested B-field range makes sense.
734 if (nb <= 0) {
735 std::cerr << m_className << "::SetFieldGrid:\n"
736 << " Number of B-fields must be > 0.\n";
737 return;
738 }
739 if (bmax < 0. || bmin < 0.) {
740 std::cerr << m_className << "::SetFieldGrid:\n"
741 << " Magnetic fields must be positive.\n";
742 return;
743 }
744 if (bmax < bmin) {
745 std::cerr << m_className << "::SetFieldGrid: Swapping min./max. B-field.\n";
746 std::swap(bmin, bmax);
747 }
748
749 const double bstep = nb > 1 ? (bmax - bmin) / (nb - 1.) : 0.;
750
751 // Check if the requested angular range makes sense.
752 if (na <= 0) {
753 std::cerr << m_className << "::SetFieldGrid:\n"
754 << " Number of angles must be > 0.\n";
755 return;
756 }
757 if (amax < 0. || amin < 0.) {
758 std::cerr << m_className << "::SetFieldGrid:"
759 << " Angles must be positive.\n";
760 return;
761 }
762 if (amax < amin) {
763 std::cerr << m_className << "::SetFieldGrid: Swapping min./max. angle.\n";
764 std::swap(amin, amax);
765 }
766 const double astep = na > 1 ? (amax - amin) / (na - 1.) : 0;
767
768 // Setup the field grids.
769 std::vector<double> eFieldsNew(ne);
770 std::vector<double> bFieldsNew(nb);
771 std::vector<double> bAnglesNew(na);
772 for (size_t i = 0; i < ne; ++i) {
773 eFieldsNew[i] = logE ? emin * pow(estep, i) : emin + i * estep;
774 }
775 for (size_t i = 0; i < nb; ++i) {
776 bFieldsNew[i] = bmin + i * bstep;
777 }
778 if (na == 1 && nb == 1 && fabs(bmin) < 1.e-4) {
779 bAnglesNew[0] = HalfPi;
780 } else {
781 for (size_t i = 0; i < na; ++i) {
782 bAnglesNew[i] = amin + i * astep;
783 }
784 }
785 SetFieldGrid(eFieldsNew, bFieldsNew, bAnglesNew);
786}
DoubleAc pow(const DoubleAc &f, double p)
Definition: DoubleAc.cpp:337

Referenced by Medium(), and SetFieldGrid().

◆ SetHoleAttachment()

bool Garfield::Medium::SetHoleAttachment ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  eta 
)
inline

Set an entry in the table of attachment coefficients.

Definition at line 359 of file Medium.hh.

360 {
361 return SetEntry(ie, ib, ia, "SetHoleAttachment", m_hAtt, eta);
362 }

◆ SetHoleLongitudinalDiffusion()

bool Garfield::Medium::SetHoleLongitudinalDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  dl 
)
inline

Set an entry in the table of longitudinal diffusion coefficients.

Definition at line 326 of file Medium.hh.

328 {
329 return SetEntry(ie, ib, ia, "SetHoleLongitudinalDiffusion", m_hDifL, dl);
330 }

◆ SetHoleTownsend()

bool Garfield::Medium::SetHoleTownsend ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  alpha 
)
inline

Set an entry in the table of Townsend coefficients.

Definition at line 349 of file Medium.hh.

350 {
351 return SetEntry(ie, ib, ia, "SetHoleTownsend", m_hAlp, alpha);
352 }

◆ SetHoleTransverseDiffusion()

bool Garfield::Medium::SetHoleTransverseDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  dt 
)
inline

Set an entry in the table of transverse diffusion coefficients.

Definition at line 338 of file Medium.hh.

340 {
341 return SetEntry(ie, ib, ia, "SetHoleTransverseDiffusion", m_hDifT, dt);
342 }

◆ SetHoleVelocityB()

bool Garfield::Medium::SetHoleVelocityB ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  v 
)
inline

Set an entry in the table of drift speeds along Btrans.

Definition at line 315 of file Medium.hh.

316 {
317 return SetEntry(ie, ib, ia, "SetHoleVelocityB", m_hVelB, v);
318 }

◆ SetHoleVelocityE()

bool Garfield::Medium::SetHoleVelocityE ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  v 
)
inline

Set an entry in the table of drift speeds along E.

Definition at line 295 of file Medium.hh.

296 {
297 return SetEntry(ie, ib, ia, "SetHoleVelocityE", m_hVelE, v);
298 }

◆ SetHoleVelocityExB()

bool Garfield::Medium::SetHoleVelocityExB ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  v 
)
inline

Set an entry in the table of drift speeds along ExB.

Definition at line 305 of file Medium.hh.

306 {
307 return SetEntry(ie, ib, ia, "SetHoleVelocityExB", m_hVelX, v);
308 }

◆ SetInterpolationMethodAttachment()

void Garfield::Medium::SetInterpolationMethodAttachment ( const unsigned int  intrp)

Definition at line 1182 of file Medium.cc.

1182 {
1183 if (intrp > 0) m_intpAtt = intrp;
1184}

◆ SetInterpolationMethodDiffusion()

void Garfield::Medium::SetInterpolationMethodDiffusion ( const unsigned int  intrp)

Definition at line 1174 of file Medium.cc.

1174 {
1175 if (intrp > 0) m_intpDif = intrp;
1176}

◆ SetInterpolationMethodIonDissociation()

void Garfield::Medium::SetInterpolationMethodIonDissociation ( const unsigned int  intrp)

Definition at line 1190 of file Medium.cc.

1190 {
1191 if (intrp > 0) m_intpDis = intrp;
1192}

◆ SetInterpolationMethodIonMobility()

void Garfield::Medium::SetInterpolationMethodIonMobility ( const unsigned int  intrp)

Definition at line 1186 of file Medium.cc.

1186 {
1187 if (intrp > 0) m_intpMob = intrp;
1188}

◆ SetInterpolationMethodTownsend()

void Garfield::Medium::SetInterpolationMethodTownsend ( const unsigned int  intrp)

Definition at line 1178 of file Medium.cc.

1178 {
1179 if (intrp > 0) m_intpAlp = intrp;
1180}

◆ SetInterpolationMethodVelocity()

void Garfield::Medium::SetInterpolationMethodVelocity ( const unsigned int  intrp)

Set the degree of polynomial interpolation (usually 2).

Definition at line 1170 of file Medium.cc.

1170 {
1171 if (intrp > 0) m_intpVel = intrp;
1172}

◆ SetIonDissociation()

bool Garfield::Medium::SetIonDissociation ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  diss 
)
inline

Set an entry in the table of dissociation coefficients.

Definition at line 405 of file Medium.hh.

406 {
407 return SetEntry(ie, ib, ia, "SetIonDissociation", m_iDis, diss);
408 }

◆ SetIonLongitudinalDiffusion()

bool Garfield::Medium::SetIonLongitudinalDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  dl 
)
inline

Set an entry in the table of longitudinal diffusion coefficients.

Definition at line 385 of file Medium.hh.

386 {
387 return SetEntry(ie, ib, ia, "SetIonLongitudinalDiffusion", m_iDifL, dl);
388 }

◆ SetIonMobility() [1/2]

bool Garfield::Medium::SetIonMobility ( const std::vector< double > &  fields,
const std::vector< double > &  mobilities 
)

Initialise the table of ion mobilities from a list of electric fields and corresponding mobilities. The mobilities will be interpolated at the electric fields of the currently set grid.

Definition at line 1044 of file Medium.cc.

1045 {
1046 const int ne = efields.size();
1047 const int nm = mobs.size();
1048 if (ne != nm) {
1049 std::cerr << m_className << "::SetIonMobility:\n"
1050 << " E-field and mobility arrays have different sizes.\n";
1051 return false;
1052 }
1053
1055 const unsigned int nEfields = m_eFields.size();
1056 const unsigned int nBfields = m_bFields.size();
1057 const unsigned int nAngles = m_bAngles.size();
1058 Init(nEfields, nBfields, nAngles, m_iMob, 0.);
1059 for (unsigned int i = 0; i < nEfields; ++i) {
1060 const double e = m_eFields[i];
1061 const double mu = Interpolate1D(e, mobs, efields, m_intpMob, m_extrMob);
1062 m_iMob[0][0][i] = mu;
1063 }
1064
1065 if (m_tab2d) {
1066 for (unsigned int i = 0; i < nAngles; ++i) {
1067 for (unsigned int j = 0; j < nBfields; ++j) {
1068 for (unsigned int k = 0; k < nEfields; ++k) {
1069 m_iMob[i][j][k] = m_iMob[0][0][k];
1070 }
1071 }
1072 }
1073 }
1074 return true;
1075}

Referenced by Garfield::MediumGas::LoadIonMobility().

◆ SetIonMobility() [2/2]

bool Garfield::Medium::SetIonMobility ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  mu 
)

Set an entry in the table of ion mobilities.

Definition at line 1014 of file Medium.cc.

1015 {
1016 // Check the index.
1017 if (ie >= m_eFields.size() || ib >= m_bFields.size() ||
1018 ia >= m_bAngles.size()) {
1019 PrintOutOfRange(m_className, "SetIonMobility", ie, ib, ia);
1020 return false;
1021 }
1022
1023 if (m_iMob.empty()) {
1024 std::cerr << m_className << "::SetIonMobility:\n"
1025 << " Ion mobility table not initialised.\n";
1026 return false;
1027 }
1028
1029 if (mu == 0.) {
1030 std::cerr << m_className << "::SetIonMobility: Zero value not allowed.\n";
1031 return false;
1032 }
1033
1034 m_iMob[ia][ib][ie] = mu;
1035 if (m_debug) {
1036 std::cout << m_className << "::SetIonMobility:\n Ion mobility at E = "
1037 << m_eFields[ie] << " V/cm, B = "
1038 << m_bFields[ib] << " T, angle "
1039 << m_bAngles[ia] << " set to " << mu << " cm2/(V ns).\n";
1040 }
1041 return true;
1042}

◆ SetIonTransverseDiffusion()

bool Garfield::Medium::SetIonTransverseDiffusion ( const unsigned int  ie,
const unsigned int  ib,
const unsigned int  ia,
const double  dt 
)
inline

Set an entry in the table of transverse diffusion coefficients.

Definition at line 395 of file Medium.hh.

396 {
397 return SetEntry(ie, ib, ia, "SetIonTransverseDiffusion", m_iDifT, dt);
398 }

◆ SetMassDensity()

void Garfield::Medium::SetMassDensity ( const double  rho)
virtual

Set the mass density [g/cm3].

Reimplemented in Garfield::MediumGas.

Definition at line 144 of file Medium.cc.

144 {
145 if (rho <= 0.) {
146 std::cerr << m_className << "::SetMassDensity:\n"
147 << " Density [g/cm3] must be greater than zero.\n";
148 return;
149 }
150
151 if (m_a <= 0.) {
152 std::cerr << m_className << "::SetMassDensity:\n"
153 << " Atomic weight is not defined.\n";
154 return;
155 }
156 m_density = rho / (AtomicMassUnit * m_a);
157 m_isChanged = true;
158}

Referenced by Garfield::MediumCdTe::MediumCdTe(), Garfield::MediumGaAs::MediumGaAs(), and Garfield::MediumSilicon::MediumSilicon().

◆ SetNumberDensity()

void Garfield::Medium::SetNumberDensity ( const double  n)
virtual

Set the number density [cm-3].

Reimplemented in Garfield::MediumGas.

Definition at line 134 of file Medium.cc.

134 {
135 if (n <= 0.) {
136 std::cerr << m_className << "::SetNumberDensity:\n"
137 << " Density [cm-3] must be greater than zero.\n";
138 return;
139 }
140 m_density = n;
141 m_isChanged = true;
142}

◆ SetPressure()

void Garfield::Medium::SetPressure ( const double  p)

Definition at line 81 of file Medium.cc.

81 {
82 if (p <= 0.) {
83 std::cerr << m_className << "::SetPressure:\n"
84 << " Pressure [Torr] must be greater than zero.\n";
85 return;
86 }
87 m_pressure = p;
88 m_isChanged = true;
89}

Referenced by GarfieldPhysics::InitializePhysics(), and main().

◆ SetTemperature()

void Garfield::Medium::SetTemperature ( const double  t)

Set the temperature [K].

Definition at line 71 of file Medium.cc.

71 {
72 if (t <= 0.) {
73 std::cerr << m_className << "::SetTemperature:\n"
74 << " Temperature [K] must be greater than zero.\n";
75 return;
76 }
77 m_temperature = t;
78 m_isChanged = true;
79}

Referenced by GarfieldPhysics::InitializePhysics(), main(), Garfield::MediumCdTe::MediumCdTe(), Garfield::MediumGaAs::MediumGaAs(), and Garfield::MediumSilicon::MediumSilicon().

◆ SetThreshold()

unsigned int Garfield::Medium::SetThreshold ( const std::vector< std::vector< std::vector< double > > > &  tab) const
protected

Definition at line 1146 of file Medium.cc.

1147 {
1148
1149 if (tab.empty()) return 0;
1150 const unsigned int nE = m_eFields.size();
1151 const unsigned int nB = m_bFields.size();
1152 const unsigned int nA = m_bAngles.size();
1153 for (unsigned int i = 0; i < nE; ++i) {
1154 bool below = false;
1155 for (unsigned int k = 0; k < nA; ++k) {
1156 for (unsigned int j = 0; j < nB; ++j) {
1157 if (tab[k][j][i] < -20.) {
1158 below = true;
1159 break;
1160 }
1161 }
1162 if (below) break;
1163 }
1164 if (below) continue;
1165 return i;
1166 }
1167 return nE - 1;
1168}

Referenced by Garfield::MediumGas::AdjustTownsendCoefficient(), Garfield::MediumMagboltz::GenerateGasTable(), and Garfield::MediumGas::MergeGasFile().

◆ SetW()

void Garfield::Medium::SetW ( const double  w)
inline

Set the W value (average energy to produce an electron/ion or e/h pair).

Definition at line 81 of file Medium.hh.

81{ m_w = w; }

◆ UnScaleElectricField()

virtual double Garfield::Medium::UnScaleElectricField ( const double  e) const
inlinevirtual

Reimplemented in Garfield::MediumGas.

Definition at line 477 of file Medium.hh.

477{ return e; }

Referenced by ElectronMobility(), and HoleMobility().

◆ Velocity()

bool Garfield::Medium::Velocity ( const double  ex,
const double  ey,
const double  ez,
const double  bx,
const double  by,
const double  bz,
const std::vector< std::vector< std::vector< double > > > &  velE,
const std::vector< std::vector< std::vector< double > > > &  velB,
const std::vector< std::vector< std::vector< double > > > &  velX,
const double  q,
double &  vx,
double &  vy,
double &  vz 
) const
protected

Definition at line 160 of file Medium.cc.

165 {
166
167 vx = vy = vz = 0.;
168 // Make sure there is at least a table of velocities along E.
169 if (velE.empty()) return false;
170
171 // Compute the magnitude of the electric field.
172 const double e = sqrt(ex * ex + ey * ey + ez * ez);
173 const double e0 = ScaleElectricField(e);
174 if (e < Small || e0 < Small) return false;
175
176 // Compute the magnitude of the magnetic field.
177 const double b = sqrt(bx * bx + by * by + bz * bz);
178 // Compute the angle between B field and E field.
179 const double ebang = GetAngle(ex, ey, ez, bx, by, bz, e, b);
180
181 // Calculate the velocity along E.
182 double ve = 0.;
183 if (!Interpolate(e0, b, ebang, velE, ve, m_intpVel, m_extrVel)) {
184 std::cerr << m_className << "::Velocity: Interpolation along E failed.\n";
185 return false;
186 }
187 if (b < Small) {
188 // No magnetic field.
189 const double mu = q * ve / e;
190 vx = mu * ex;
191 vy = mu * ey;
192 vz = mu * ez;
193 return true;
194 } else if (velX.empty() || velB.empty()) {
195 // Magnetic field, velocities along ExB, Bt not available.
196 const double mu = q * ve / e;
197 const double mu2 = mu * mu;
198 const double eb = bx * ex + by * ey + bz * ez;
199 const double f = mu / (1. + mu2 * b * b);
200 vx = f * (ex + mu * (ey * bz - ez * by) + mu2 * bx * eb);
201 vy = f * (ey + mu * (ez * bx - ex * bz) + mu2 * by * eb);
202 vz = f * (ez + mu * (ex * by - ey * bx) + mu2 * bz * eb);
203 return true;
204 }
205
206 // Magnetic field, velocities along ExB and Bt available.
207 // Compute unit vectors along E, E x B and Bt.
208 double ue[3] = {ex / e, ey / e, ez / e};
209 double uexb[3] = {ey * bz - ez * by, ez * bx - ex * bz, ex * by - ey * bx};
210 const double exb =
211 sqrt(uexb[0] * uexb[0] + uexb[1] * uexb[1] + uexb[2] * uexb[2]);
212 if (exb > 0.) {
213 uexb[0] /= exb;
214 uexb[1] /= exb;
215 uexb[2] /= exb;
216 } else {
217 uexb[0] = ue[0];
218 uexb[1] = ue[1];
219 uexb[2] = ue[2];
220 }
221
222 double ubt[3] = {uexb[1] * ez - uexb[2] * ey, uexb[2] * ex - uexb[0] * ez,
223 uexb[0] * ey - uexb[1] * ex};
224 const double bt = sqrt(ubt[0] * ubt[0] + ubt[1] * ubt[1] + ubt[2] * ubt[2]);
225 if (bt > 0.) {
226 ubt[0] /= bt;
227 ubt[1] /= bt;
228 ubt[2] /= bt;
229 } else {
230 ubt[0] = ue[0];
231 ubt[1] = ue[1];
232 ubt[2] = ue[2];
233 }
234
235 if (m_debug) {
236 std::cout << std::setprecision(5);
237 std::cout << m_className << "::Velocity:\n"
238 << " unit vector along E: (" << ue[0] << ", " << ue[1]
239 << ", " << ue[2] << ")\n";
240 std::cout << " unit vector along E x B: (" << uexb[0] << ", "
241 << uexb[1] << ", " << uexb[2] << ")\n";
242 std::cout << " unit vector along Bt: (" << ubt[0] << ", " << ubt[1]
243 << ", " << ubt[2] << ")\n";
244 }
245
246 // Calculate the velocities in all directions.
247 double vexb = 0.;
248 if (!Interpolate(e0, b, ebang, velX, vexb, m_intpVel, m_extrVel)) {
249 std::cerr << m_className << "::Velocity: Interpolation along ExB failed.\n";
250 return false;
251 }
252 double vbt = 0.;
253 if (!Interpolate(e0, b, ebang, velB, vbt, m_intpVel, m_extrVel)) {
254 std::cerr << m_className << "::Velocity: Interpolation along Bt failed.\n";
255 return false;
256 }
257 if (ex * bx + ey * by + ez * bz > 0.) {
258 vbt = fabs(vbt);
259 } else {
260 vbt = -fabs(vbt);
261 }
262 vx = q * (ve * ue[0] + q * q * vbt * ubt[0] + q * vexb * uexb[0]);
263 vy = q * (ve * ue[1] + q * q * vbt * ubt[1] + q * vexb * uexb[1]);
264 vz = q * (ve * ue[2] + q * q * vbt * ubt[2] + q * vexb * uexb[2]);
265
266 return true;
267}

Referenced by ElectronVelocity(), and HoleVelocity().

Member Data Documentation

◆ m_a

double Garfield::Medium::m_a = 0.
protected

Definition at line 525 of file Medium.hh.

Referenced by GetAtomicWeight(), GetMassDensity(), SetAtomicWeight(), and SetMassDensity().

◆ m_bAngles

◆ m_bFields

◆ m_className

std::string Garfield::Medium::m_className = "Medium"
protected

Definition at line 506 of file Medium.hh.

Referenced by Garfield::MediumGas::AdjustTownsendCoefficient(), Clone(), Garfield::MediumMagboltz::ComputeDeexcitation(), Garfield::MediumGas::DisablePenningTransfer(), Garfield::MediumMagboltz::DisablePenningTransfer(), Garfield::MediumSilicon::ElectronAttachment(), Garfield::MediumSilicon::ElectronTownsend(), Garfield::MediumSilicon::ElectronVelocity(), Garfield::MediumMagboltz::EnableDeexcitation(), Garfield::MediumGas::EnablePenningTransfer(), Garfield::MediumMagboltz::EnablePenningTransfer(), Garfield::MediumMagboltz::EnableRadiationTrapping(), Garfield::MediumMagboltz::GenerateGasTable(), GetComponent(), Garfield::MediumCdTe::GetComponent(), Garfield::MediumGas::GetComponent(), Garfield::MediumSilicon::GetConductionBandDensityOfStates(), GetDeexcitationProduct(), GetDielectricFunction(), Garfield::MediumSilicon::GetDielectricFunction(), Garfield::MediumSilicon::GetElectronBandPopulation(), GetElectronCollision(), Garfield::MediumMagboltz::GetElectronCollision(), Garfield::MediumSilicon::GetElectronCollision(), Garfield::MediumMagboltz::GetElectronCollisionRate(), Garfield::MediumSilicon::GetElectronCollisionRate(), GetElectronCollisionRate(), GetElectronEnergy(), Garfield::MediumSilicon::GetElectronEnergy(), Garfield::MediumSilicon::GetElectronMomentum(), Garfield::MediumMagboltz::GetElectronNullCollisionRate(), Garfield::MediumSilicon::GetElectronNullCollisionRate(), GetElectronNullCollisionRate(), GetEntry(), Garfield::MediumGas::GetGasName(), Garfield::MediumGas::GetGasNumberGasFile(), Garfield::MediumMagboltz::GetLevel(), Garfield::MediumGas::GetMixture(), Garfield::MediumMagboltz::GetNumberOfElectronCollisions(), Garfield::MediumSilicon::GetNumberOfElectronCollisions(), Garfield::MediumMagboltz::GetNumberOfLevels(), GetOpticalDataRange(), Garfield::MediumSilicon::GetOpticalDataRange(), Garfield::MediumGas::GetPhotoAbsorptionCrossSection(), GetPhotoAbsorptionCrossSection(), Garfield::MediumMagboltz::GetPhotonCollision(), Garfield::MediumMagboltz::GetPhotonCollisionRate(), Garfield::MediumSilicon::GetValenceBandDensityOfStates(), Garfield::MediumSilicon::HoleAttachment(), Garfield::MediumSilicon::HoleTownsend(), Garfield::MediumSilicon::HoleVelocity(), Init(), Garfield::MediumSilicon::Initialise(), Garfield::MediumMagboltz::Initialise(), Interpolate1D(), Garfield::MediumGas::LoadGasFile(), Garfield::MediumGas::LoadIonMobility(), Garfield::MediumCdTe::MediumCdTe(), Garfield::MediumConductor::MediumConductor(), Garfield::MediumGaAs::MediumGaAs(), Garfield::MediumGas::MediumGas(), Garfield::MediumMagboltz::MediumMagboltz(), Garfield::MediumPlastic::MediumPlastic(), Garfield::MediumSilicon::MediumSilicon(), Garfield::MediumGas::MergeGasFile(), Garfield::MediumGas::PrintGas(), Garfield::MediumGas::ReadHeader(), Garfield::MediumMagboltz::RunMagboltz(), SetAtomicNumber(), Garfield::MediumGas::SetAtomicNumber(), SetAtomicWeight(), Garfield::MediumGas::SetAtomicWeight(), Garfield::MediumGas::SetComposition(), SetDielectricConstant(), Garfield::MediumSilicon::SetDoping(), SetEntry(), Garfield::MediumMagboltz::SetExcitationScaling(), SetExtrapolationMethod(), SetFieldGrid(), SetIonMobility(), Garfield::MediumCdTe::SetLowFieldMobility(), Garfield::MediumGaAs::SetLowFieldMobility(), Garfield::MediumSilicon::SetLowFieldMobility(), SetMassDensity(), Garfield::MediumGas::SetMassDensity(), Garfield::MediumMagboltz::SetMaxElectronEnergy(), Garfield::MediumSilicon::SetMaxElectronEnergy(), Garfield::MediumMagboltz::SetMaxPhotonEnergy(), SetNumberDensity(), Garfield::MediumGas::SetNumberDensity(), SetPressure(), Garfield::MediumSilicon::SetSaturationVelocity(), Garfield::MediumMagboltz::SetSplittingFunctionGreenSawada(), SetTemperature(), Garfield::MediumSilicon::SetTrapCrossSection(), Garfield::MediumSilicon::SetTrapDensity(), Garfield::MediumSilicon::SetTrappingTime(), Velocity(), and Garfield::MediumGas::WriteGasFile().

◆ m_debug

◆ m_density

double Garfield::Medium::m_density = 0.
protected

◆ m_driftable

bool Garfield::Medium::m_driftable = false
protected

Definition at line 530 of file Medium.hh.

Referenced by EnableDrift(), and IsDriftable().

◆ m_eAlp

◆ m_eAtt

◆ m_eDifL

◆ m_eDifM

◆ m_eDifT

◆ m_eFields

◆ m_eLor

◆ m_epsilon

double Garfield::Medium::m_epsilon = 1.
protected

Definition at line 519 of file Medium.hh.

Referenced by GetDielectricConstant(), and SetDielectricConstant().

◆ m_eThrAlp

unsigned int Garfield::Medium::m_eThrAlp = 0
protected

◆ m_eThrAtt

unsigned int Garfield::Medium::m_eThrAtt = 0
protected

◆ m_eVelB

◆ m_eVelE

◆ m_eVelX

◆ m_extrAlp

std::pair<unsigned int, unsigned int> Garfield::Medium::m_extrAlp = {0, 1}
protected

◆ m_extrAtt

◆ m_extrDif

std::pair<unsigned int, unsigned int> Garfield::Medium::m_extrDif = {0, 1}
protected

◆ m_extrDis

std::pair<unsigned int, unsigned int> Garfield::Medium::m_extrDis = {0, 1}
protected

◆ m_extrLor

std::pair<unsigned int, unsigned int> Garfield::Medium::m_extrLor = {0, 1}
protected

◆ m_extrMob

std::pair<unsigned int, unsigned int> Garfield::Medium::m_extrMob = {0, 1}
protected

◆ m_extrVel

std::pair<unsigned int, unsigned int> Garfield::Medium::m_extrVel = {0, 1}
protected

◆ m_fano

double Garfield::Medium::m_fano = 0.
protected

◆ m_hAlp

std::vector<std::vector<std::vector<double> > > Garfield::Medium::m_hAlp
protected

◆ m_hAtt

std::vector<std::vector<std::vector<double> > > Garfield::Medium::m_hAtt
protected

◆ m_hDifL

std::vector<std::vector<std::vector<double> > > Garfield::Medium::m_hDifL
protected

◆ m_hDifM

std::vector<std::vector<std::vector<std::vector<double> > > > Garfield::Medium::m_hDifM
protected

Definition at line 573 of file Medium.hh.

Referenced by HoleDiffusion(), ResetHoleDiffusion(), and SetFieldGrid().

◆ m_hDifT

std::vector<std::vector<std::vector<double> > > Garfield::Medium::m_hDifT
protected

◆ m_hThrAlp

unsigned int Garfield::Medium::m_hThrAlp = 0
protected

Definition at line 584 of file Medium.hh.

Referenced by HoleTownsend().

◆ m_hThrAtt

unsigned int Garfield::Medium::m_hThrAtt = 0
protected

Definition at line 585 of file Medium.hh.

Referenced by HoleAttachment().

◆ m_hVelB

std::vector<std::vector<std::vector<double> > > Garfield::Medium::m_hVelB
protected

◆ m_hVelE

std::vector<std::vector<std::vector<double> > > Garfield::Medium::m_hVelE
protected

◆ m_hVelX

std::vector<std::vector<std::vector<double> > > Garfield::Medium::m_hVelX
protected

◆ m_id

int Garfield::Medium::m_id
protected

Definition at line 511 of file Medium.hh.

Referenced by GetId().

◆ m_idCounter

int Garfield::Medium::m_idCounter = -1
staticprotected

Definition at line 508 of file Medium.hh.

◆ m_iDifL

◆ m_iDifT

◆ m_iDis

◆ m_iMob

◆ m_intpAlp

◆ m_intpAtt

◆ m_intpDif

◆ m_intpDis

◆ m_intpLor

◆ m_intpMob

◆ m_intpVel

◆ m_ionisable

bool Garfield::Medium::m_ionisable = false
protected

Definition at line 532 of file Medium.hh.

Referenced by EnablePrimaryIonisation(), and IsIonisable().

◆ m_isChanged

bool Garfield::Medium::m_isChanged = true
protected

Definition at line 540 of file Medium.hh.

Referenced by Garfield::MediumMagboltz::ComputeDeexcitation(), Garfield::MediumSilicon::ElectronAttachment(), Garfield::MediumGaAs::ElectronTownsend(), Garfield::MediumSilicon::ElectronTownsend(), Garfield::MediumCdTe::ElectronVelocity(), Garfield::MediumGaAs::ElectronVelocity(), Garfield::MediumSilicon::ElectronVelocity(), Garfield::MediumMagboltz::EnableAnisotropicScattering(), Garfield::MediumMagboltz::EnableDeexcitation(), Garfield::MediumMagboltz::EnablePenningTransfer(), Garfield::MediumMagboltz::EnableRadiationTrapping(), Garfield::MediumMagboltz::GetElectronCollision(), Garfield::MediumSilicon::GetElectronCollision(), Garfield::MediumMagboltz::GetElectronCollisionRate(), Garfield::MediumSilicon::GetElectronCollisionRate(), Garfield::MediumMagboltz::GetElectronNullCollisionRate(), Garfield::MediumSilicon::GetElectronNullCollisionRate(), Garfield::MediumMagboltz::GetLevel(), Garfield::MediumMagboltz::GetNumberOfLevels(), Garfield::MediumMagboltz::GetPhotonCollision(), Garfield::MediumMagboltz::GetPhotonCollisionRate(), Garfield::MediumSilicon::HoleAttachment(), Garfield::MediumGaAs::HoleTownsend(), Garfield::MediumSilicon::HoleTownsend(), Garfield::MediumCdTe::HoleVelocity(), Garfield::MediumGaAs::HoleVelocity(), Garfield::MediumSilicon::HoleVelocity(), Garfield::MediumSilicon::Initialise(), Garfield::MediumMagboltz::Initialise(), Garfield::MediumGas::LoadGasFile(), Garfield::MediumGas::MediumGas(), Garfield::MediumMagboltz::MediumMagboltz(), Garfield::MediumMagboltz::PrintGas(), SetAtomicNumber(), SetAtomicWeight(), Garfield::MediumGas::SetComposition(), SetDielectricConstant(), Garfield::MediumSilicon::SetDoping(), Garfield::MediumSilicon::SetDopingMobilityModelMasetti(), Garfield::MediumSilicon::SetDopingMobilityModelMinimos(), Garfield::MediumMagboltz::SetExcitationScaling(), Garfield::MediumSilicon::SetHighFieldMobilityModelCanali(), Garfield::MediumSilicon::SetHighFieldMobilityModelMinimos(), Garfield::MediumSilicon::SetHighFieldMobilityModelReggiani(), Garfield::MediumSilicon::SetImpactIonisationModelGrant(), Garfield::MediumSilicon::SetImpactIonisationModelMassey(), Garfield::MediumSilicon::SetImpactIonisationModelVanOverstraetenDeMan(), Garfield::MediumSilicon::SetLatticeMobilityModelMinimos(), Garfield::MediumSilicon::SetLatticeMobilityModelReggiani(), Garfield::MediumSilicon::SetLatticeMobilityModelSentaurus(), Garfield::MediumCdTe::SetLowFieldMobility(), Garfield::MediumGaAs::SetLowFieldMobility(), Garfield::MediumSilicon::SetLowFieldMobility(), SetMassDensity(), Garfield::MediumMagboltz::SetMaxElectronEnergy(), Garfield::MediumSilicon::SetMaxElectronEnergy(), Garfield::MediumMagboltz::SetMaxPhotonEnergy(), SetNumberDensity(), SetPressure(), Garfield::MediumSilicon::SetSaturationVelocity(), Garfield::MediumSilicon::SetSaturationVelocityModelCanali(), Garfield::MediumSilicon::SetSaturationVelocityModelMinimos(), Garfield::MediumSilicon::SetSaturationVelocityModelReggiani(), Garfield::MediumMagboltz::SetSplittingFunctionGreenSawada(), SetTemperature(), Garfield::MediumSilicon::SetTrapCrossSection(), Garfield::MediumSilicon::SetTrapDensity(), Garfield::MediumSilicon::SetTrappingTime(), Garfield::MediumCdTe::UnsetLowFieldMobility(), and Garfield::MediumGaAs::UnsetLowFieldMobility().

◆ m_iThrDis

unsigned int Garfield::Medium::m_iThrDis = 0
protected

◆ m_microscopic

◆ m_name

◆ m_nComponents

◆ m_pressure

◆ m_tab2d

◆ m_temperature

◆ m_w

double Garfield::Medium::m_w = 0.
protected

◆ m_z

double Garfield::Medium::m_z = 1.
protected

Definition at line 523 of file Medium.hh.

Referenced by GetAtomicNumber(), and SetAtomicNumber().


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