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::ViewMedium Class Reference

Plot transport coefficients as function of electric and magnetic field. More...

#include <ViewMedium.hh>

+ Inheritance diagram for Garfield::ViewMedium:

Public Types

enum  Property {
  ElectronVelocityE , ElectronTransverseDiffusion , ElectronLongitudinalDiffusion , ElectronTownsend ,
  ElectronAttachment , ElectronLorentzAngle , HoleVelocityE = 10 , HoleTransverseDiffusion ,
  HoleLongitudinalDiffusion , HoleTownsend , HoleAttachment , IonVelocity = 20 ,
  IonTransverseDiffusion , IonLongitudinalDiffusion , ElectronVelocityB , ElectronVelocityExB ,
  HoleVelocityB , HoleVelocityExB
}
 

Public Member Functions

 ViewMedium ()
 Constructor.
 
 ~ViewMedium ()=default
 Destructor.
 
void SetMedium (Medium *m)
 Set the medium from which to retrieve the transport coefficients.
 
void EnableAutoRangeX (const bool on=true)
 Try to choose the x-axis range based on the field grid.
 
void SetRangeE (const double emin, const double emax, const bool logscale)
 Set the limits of the electric field.
 
void SetRangeB (const double bmin, const double bmax, const bool logscale)
 Set the limits of the magnetic field.
 
void SetRangeA (const double amin, const double amax, const bool logscale)
 Set the limits of the angle between electric and magnetic field.
 
void EnableAutoRangeY (const bool on=true)
 Choose the y-axis range based on the function's minima/maxima.
 
void SetRangeY (const double ymin, const double ymax, const bool logscale=false)
 Set the range of the function (velocity etc.) to be plotted.
 
void SetElectricField (const double efield)
 Set the electric field to use when plotting as function of B or angle.
 
void SetMagneticField (const double bfield)
 Set the magnetic field to use when plotting as function of E or angle.
 
void SetAngle (const double angle)
 Set the angle to use when plotting as function of E or B.
 
void PlotElectronVelocity (const char xaxis, const bool same=false)
 
void PlotHoleVelocity (const char xaxis, const bool same=false)
 
void PlotIonVelocity (const char xaxis, const bool same=false)
 
void PlotElectronDiffusion (const char xaxis, const bool same=false)
 
void PlotHoleDiffusion (const char xaxis, const bool same=false)
 
void PlotIonDiffusion (const char xaxis, const bool same=false)
 
void PlotElectronTownsend (const char xaxis, const bool same=false)
 
void PlotHoleTownsend (const char xaxis, const bool same=false)
 
void PlotElectronAttachment (const char xaxis, const bool same=false)
 
void PlotHoleAttachment (const char xaxis, const bool same=false)
 
void PlotElectronLorentzAngle (const char xaxis, const bool same=false)
 
void PlotElectronCrossSections ()
 
double EvaluateFunction (double *pos, double *par)
 
- Public Member Functions inherited from Garfield::ViewBase
 ViewBase ()=delete
 Default constructor.
 
 ViewBase (const std::string &name)
 Constructor.
 
virtual ~ViewBase ()
 Destructor.
 
void SetCanvas (TCanvas *c)
 Set the canvas to be painted on.
 
TCanvas * GetCanvas ()
 Retrieve the canvas.
 
void EnableDebugging (const bool on=true)
 Switch on/off debugging output.
 

Additional Inherited Members

- Protected Member Functions inherited from Garfield::ViewBase
std::string FindUnusedFunctionName (const std::string &s) const
 
std::string FindUnusedHistogramName (const std::string &s) const
 
- Protected Attributes inherited from Garfield::ViewBase
std::string m_className = "ViewBase"
 
bool m_debug = false
 
TCanvas * m_canvas = nullptr
 
bool m_hasExternalCanvas = false
 
double m_proj [3][3]
 

Detailed Description

Plot transport coefficients as function of electric and magnetic field.

Definition at line 19 of file ViewMedium.hh.

Member Enumeration Documentation

◆ Property

Enumerator
ElectronVelocityE 
ElectronTransverseDiffusion 
ElectronLongitudinalDiffusion 
ElectronTownsend 
ElectronAttachment 
ElectronLorentzAngle 
HoleVelocityE 
HoleTransverseDiffusion 
HoleLongitudinalDiffusion 
HoleTownsend 
HoleAttachment 
IonVelocity 
IonTransverseDiffusion 
IonLongitudinalDiffusion 
ElectronVelocityB 
ElectronVelocityExB 
HoleVelocityB 
HoleVelocityExB 

Definition at line 65 of file ViewMedium.hh.

65 {
72 HoleVelocityE = 10,
77 IonVelocity = 20,
84 };

Constructor & Destructor Documentation

◆ ViewMedium()

Garfield::ViewMedium::ViewMedium ( )

Constructor.

Definition at line 32 of file ViewMedium.cc.

32: ViewBase("ViewMedium") {}
ViewBase()=delete
Default constructor.

◆ ~ViewMedium()

Garfield::ViewMedium::~ViewMedium ( )
default

Destructor.

Member Function Documentation

◆ EnableAutoRangeX()

void Garfield::ViewMedium::EnableAutoRangeX ( const bool  on = true)
inline

Try to choose the x-axis range based on the field grid.

Definition at line 30 of file ViewMedium.hh.

30{ m_autoRangeX = on; }

◆ EnableAutoRangeY()

void Garfield::ViewMedium::EnableAutoRangeY ( const bool  on = true)
inline

Choose the y-axis range based on the function's minima/maxima.

Definition at line 38 of file ViewMedium.hh.

38{ m_autoRangeY = on; }

◆ EvaluateFunction()

double Garfield::ViewMedium::EvaluateFunction ( double *  pos,
double *  par 
)

Definition at line 406 of file ViewMedium.cc.

406 {
407 // to be modified to include B and angle
408
409 if (!m_medium) return 0.;
410 int type = int(par[0]);
411 char xaxis = char(par[1]);
412 const double x = pos[0];
413
414 const double ctheta = xaxis == 'a' ? cos(x) : cos(par[4]);
415 const double stheta = xaxis == 'a' ? sin(x) : sin(par[4]);
416 double ex = xaxis == 'e' ? x : par[2];
417 double ey = 0.;
418 double bx = xaxis == 'b' ? x : par[3];
419 double by = 0.;
420 if (type == ElectronVelocityExB || type == HoleVelocityExB ||
421 type == ElectronVelocityB || type == HoleVelocityB) {
422 ex *= ctheta;
423 ey = xaxis == 'e' ? x * stheta : par[2] * stheta;
424 } else {
425 bx *= ctheta;
426 by = xaxis == 'b' ? x * stheta : par[3] * stheta;
427 }
428 // Return value.
429 double y = 0.;
430 // Auxiliary (dummy) variables.
431 double s = 0., t = 0.;
432 switch (type) {
434 if (!m_medium->ElectronVelocity(ex, 0, 0, bx, by, 0, y, s, t)) return 0.;
435 y = fabs(y);
436 break;
438 if (!m_medium->ElectronDiffusion(ex, 0, 0, bx, by, 0, s, y)) return 0.;
439 break;
441 if (!m_medium->ElectronDiffusion(ex, 0, 0, bx, by, 0, y, s)) return 0.;
442 break;
443 case ElectronTownsend:
444 if (!m_medium->ElectronTownsend(ex, 0, 0, bx, by, 0, y)) return 0.;
445 break;
447 if (!m_medium->ElectronAttachment(ex, 0, 0, bx, by, 0, y)) return 0.;
448 break;
450 if (!m_medium->ElectronLorentzAngle(ex, 0, 0, bx, by, 0, y)) return 0.;
451 break;
452 case HoleVelocityE:
453 if (!m_medium->HoleVelocity(ex, 0, 0, bx, by, 0, y, s, t)) return 0.;
454 break;
456 if (!m_medium->HoleDiffusion(ex, 0, 0, bx, by, 0, s, y)) return 0.;
457 break;
459 if (!m_medium->HoleDiffusion(ex, 0, 0, bx, by, 0, y, s)) return 0.;
460 break;
461 case HoleTownsend:
462 if (!m_medium->HoleTownsend(ex, 0, 0, bx, by, 0, y)) return 0.;
463 break;
464 case HoleAttachment:
465 if (!m_medium->HoleAttachment(ex, 0, 0, bx, by, 0, y)) return 0.;
466 break;
467 case IonVelocity:
468 if (!m_medium->IonVelocity(ex, 0, 0, bx, by, 0, y, s, t)) return 0.;
469 break;
471 if (!m_medium->IonDiffusion(ex, 0, 0, bx, by, 0, s, y)) return 0.;
472 break;
474 if (!m_medium->IonDiffusion(ex, 0, 0, bx, by, 0, y, s)) return 0.;
475 break;
477 if (!m_medium->ElectronVelocity(ex, ey, 0, bx, 0, 0, y, s, t)) return 0.;
478 y = fabs(y);
479 break;
481 if (!m_medium->ElectronVelocity(ex, ey, 0, bx, 0, 0, s, t, y)) return 0.;
482 y = fabs(y);
483 break;
484 case HoleVelocityB:
485 if (!m_medium->HoleVelocity(ex, ey, 0, bx, 0, 0, y, s, t)) return 0.;
486 y = fabs(y);
487 break;
488 case HoleVelocityExB:
489 if (!m_medium->HoleVelocity(ex, ey, 0, bx, 0, 0, s, t, y)) return 0.;
490 y = fabs(y);
491 break;
492 default:
493 std::cerr << m_className << "::EvaluateFunction:\n "
494 << "Unknown type of transport coefficient requested. Bug!\n";
495 return 0.;
496 }
497
498 return y;
499}
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].
Definition: Medium.cc:534
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].
Definition: Medium.cc:513
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.
Definition: Medium.cc:429
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].
Definition: Medium.cc:379
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].
Definition: Medium.cc:565
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].
Definition: Medium.cc:387
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].
Definition: Medium.cc:403
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].
Definition: Medium.cc:416
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].
Definition: Medium.cc:521
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].
Definition: Medium.cc:547
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].
Definition: Medium.cc:600
std::string m_className
Definition: ViewBase.hh:28
DoubleAc cos(const DoubleAc &f)
Definition: DoubleAc.cpp:432
DoubleAc fabs(const DoubleAc &f)
Definition: DoubleAc.h:615
DoubleAc sin(const DoubleAc &f)
Definition: DoubleAc.cpp:384

◆ PlotElectronAttachment()

void Garfield::ViewMedium::PlotElectronAttachment ( const char  xaxis,
const bool  same = false 
)

Definition at line 136 of file ViewMedium.cc.

136 {
137 SetupCanvas();
138 AddFunction(same, ElectronAttachment, xaxis);
139}

◆ PlotElectronCrossSections()

void Garfield::ViewMedium::PlotElectronCrossSections ( )

Definition at line 151 of file ViewMedium.cc.

151 {
152 std::cerr << m_className << "::PlotElectronCrossSections: Not implemented.\n";
153}

◆ PlotElectronDiffusion()

void Garfield::ViewMedium::PlotElectronDiffusion ( const char  xaxis,
const bool  same = false 
)

Definition at line 108 of file ViewMedium.cc.

108 {
109 SetupCanvas();
110 AddFunction(same, ElectronTransverseDiffusion, xaxis);
111 AddFunction(true, ElectronLongitudinalDiffusion, xaxis);
112}

◆ PlotElectronLorentzAngle()

void Garfield::ViewMedium::PlotElectronLorentzAngle ( const char  xaxis,
const bool  same = false 
)

Definition at line 146 of file ViewMedium.cc.

146 {
147 SetupCanvas();
148 AddFunction(same, ElectronLorentzAngle, xaxis);
149}

◆ PlotElectronTownsend()

void Garfield::ViewMedium::PlotElectronTownsend ( const char  xaxis,
const bool  same = false 
)

Definition at line 126 of file ViewMedium.cc.

126 {
127 SetupCanvas();
128 AddFunction(same, ElectronTownsend, xaxis);
129}

◆ PlotElectronVelocity()

void Garfield::ViewMedium::PlotElectronVelocity ( const char  xaxis,
const bool  same = false 
)

Definition at line 89 of file ViewMedium.cc.

89 {
90 SetupCanvas();
91 AddFunction(same, ElectronVelocityE, xaxis);
92 AddFunction(true, ElectronVelocityB, xaxis);
93 AddFunction(true, ElectronVelocityExB, xaxis);
94}

◆ PlotHoleAttachment()

void Garfield::ViewMedium::PlotHoleAttachment ( const char  xaxis,
const bool  same = false 
)

Definition at line 141 of file ViewMedium.cc.

141 {
142 SetupCanvas();
143 AddFunction(same, HoleAttachment, xaxis);
144}

◆ PlotHoleDiffusion()

void Garfield::ViewMedium::PlotHoleDiffusion ( const char  xaxis,
const bool  same = false 
)

Definition at line 114 of file ViewMedium.cc.

114 {
115 SetupCanvas();
116 AddFunction(same, HoleTransverseDiffusion, xaxis);
117 AddFunction(true, HoleLongitudinalDiffusion, xaxis);
118}

◆ PlotHoleTownsend()

void Garfield::ViewMedium::PlotHoleTownsend ( const char  xaxis,
const bool  same = false 
)

Definition at line 131 of file ViewMedium.cc.

131 {
132 SetupCanvas();
133 AddFunction(same, HoleTownsend, xaxis);
134}

◆ PlotHoleVelocity()

void Garfield::ViewMedium::PlotHoleVelocity ( const char  xaxis,
const bool  same = false 
)

Definition at line 96 of file ViewMedium.cc.

96 {
97 SetupCanvas();
98 AddFunction(same, HoleVelocityE, xaxis);
99 AddFunction(true, HoleVelocityB, xaxis);
100 AddFunction(true, HoleVelocityExB, xaxis);
101}

◆ PlotIonDiffusion()

void Garfield::ViewMedium::PlotIonDiffusion ( const char  xaxis,
const bool  same = false 
)

Definition at line 120 of file ViewMedium.cc.

120 {
121 SetupCanvas();
122 AddFunction(same, IonTransverseDiffusion, xaxis);
123 AddFunction(true, IonLongitudinalDiffusion, xaxis);
124}

◆ PlotIonVelocity()

void Garfield::ViewMedium::PlotIonVelocity ( const char  xaxis,
const bool  same = false 
)

Definition at line 103 of file ViewMedium.cc.

103 {
104 SetupCanvas();
105 AddFunction(same, IonVelocity, xaxis);
106}

◆ SetAngle()

void Garfield::ViewMedium::SetAngle ( const double  angle)
inline

Set the angle to use when plotting as function of E or B.

Definition at line 48 of file ViewMedium.hh.

48{ m_angle = angle; }

◆ SetElectricField()

void Garfield::ViewMedium::SetElectricField ( const double  efield)
inline

Set the electric field to use when plotting as function of B or angle.

Definition at line 44 of file ViewMedium.hh.

44{ m_efield = efield; }

◆ SetMagneticField()

void Garfield::ViewMedium::SetMagneticField ( const double  bfield)
inline

Set the magnetic field to use when plotting as function of E or angle.

Definition at line 46 of file ViewMedium.hh.

46{ m_bfield = bfield; }

◆ SetMedium()

void Garfield::ViewMedium::SetMedium ( Medium m)

Set the medium from which to retrieve the transport coefficients.

Definition at line 34 of file ViewMedium.cc.

34 {
35 if (!m) {
36 std::cerr << m_className << "::SetMedium: Null pointer.\n";
37 return;
38 }
39
40 m_medium = m;
41}

◆ SetRangeA()

void Garfield::ViewMedium::SetRangeA ( const double  amin,
const double  amax,
const bool  logscale 
)

Set the limits of the angle between electric and magnetic field.

Definition at line 66 of file ViewMedium.cc.

66 {
67 if (amin >= amax || amin < 0.) {
68 std::cerr << m_className << "::SetRangeA: Incorrect range.\n";
69 return;
70 }
71
72 m_aMin = amin;
73 m_aMax = amax;
74 m_logA = logscale;
75}

◆ SetRangeB()

void Garfield::ViewMedium::SetRangeB ( const double  bmin,
const double  bmax,
const bool  logscale 
)

Set the limits of the magnetic field.

Definition at line 55 of file ViewMedium.cc.

55 {
56 if (bmin >= bmax || bmin < 0.) {
57 std::cerr << m_className << "::SetRangeB: Incorrect range.\n";
58 return;
59 }
60
61 m_bMin = bmin;
62 m_bMax = bmax;
63 m_logB = logscale;
64}

◆ SetRangeE()

void Garfield::ViewMedium::SetRangeE ( const double  emin,
const double  emax,
const bool  logscale 
)

Set the limits of the electric field.

Definition at line 43 of file ViewMedium.cc.

44 {
45 if (emin >= emax || emin < 0.) {
46 std::cerr << m_className << "::SetRangeE: Incorrect range.\n";
47 return;
48 }
49
50 m_eMin = emin;
51 m_eMax = emax;
52 m_logE = logscale;
53}

◆ SetRangeY()

void Garfield::ViewMedium::SetRangeY ( const double  ymin,
const double  ymax,
const bool  logscale = false 
)

Set the range of the function (velocity etc.) to be plotted.

Definition at line 77 of file ViewMedium.cc.

78 {
79 if (ymin >= ymax || ymin < 0.) {
80 std::cerr << m_className << "::SetRangeY: Incorrect range.\n";
81 return;
82 }
83
84 m_yMin = ymin;
85 m_yMax = ymax;
86 m_logY = logscale;
87}

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