ROOT logo
KDetSim » KDetector

class KDetector: public KGeometry, public KMaterial


KDetector

The base class for all detectors. It incorporates the calculation of
electric and weithing field as well as simulation of drift.

Calculation of the drift:

#frac{d^{2}f(x)}{dx^{2}}=#frac{f(x+h1)-2*f(x+0.5(h1-h2))+f(x-h2)}{(0.5 h1+ 0.5 h2)^{2}}

Function Members (Methods)

public:
KDetector()
KDetector(const KDetector&)
~KDetector()
voidCalField(Int_t)
voidCalM(KStruct* seg, Double_t* data, Int_t = -1)
voidCalPhyField()
voidCalRamoField()
static TClass*Class()
voidDeclaration(Int_t)
static Double_tKMaterial::dEdx(Double_t)
static Float_tKMaterial::dEX(Double_t, Double_t*, Double_t*, Double_t)
TH2F*Draw(Char_t*, Float_t = 1)
TH1F*Draw1D(Char_t*, Float_t, Int_t, Float_t)
voidDrift(Double_t, Double_t, Double_t, Float_t, KStruct*, Double_t = 0)
voidKGeometry::ElCylinder(Float_t* Pos, Float_t R, Float_t L, Int_t O, Int_t Wei, Int_t Mat)
voidKGeometry::ElLine(Float_t* r0, Float_t* r1, Float_t* W, Int_t Wei, Int_t Mat)
voidKGeometry::ElRectangle(Float_t* Pos, Float_t* Size, Int_t Wei, Int_t Mat)
voidKGeometry::ElRectangle(Float_t* Pos, Float_t* Size, TH3F*, Float_t)
voidGaussBeam(Int_t, Float_t, Float_t, Float_t, Float_t, Float_t)
TH3F*KGeometry::GetGeom()
voidKGeometry::GetGrid(TH3I*, Short_t = 0)
Float_tKGeometry::GetLowEdge(Int_t)
Double_tGetPrecision()
Double_tKGeometry::GetStepSize(Int_t, Int_t)
Double_tKGeometry::GetStepSize(Int_t, Float_t)
Float_tKGeometry::GetUpEdge(Int_t)
virtual TClass*IsA() const
Double_tkappa(int, int, int, int)
TH3F*KGeometry::MapToGeometry(Double_t*, Double_t = 1)
voidMipIR(Int_t = 20, Float_t = 0)
static Int_tKMaterial::MobMod()
KDetector&operator=(const KDetector&)
static Float_tKMaterial::Perm(Int_t = 1)
TFile*Read(Char_t*, Char_t*)
voidKGeometry::Reset(Int_t = 0, Int_t = 0)
voidResetRnd(Int_t seed)
static Float_tKMaterial::Rho()
voidSave(Char_t*, Char_t*)
Int_tKGeometry::SetBoundaryConditions()
voidSetCalculationParameters(Double_t x, Int_t y)
voidSetDebug(Short_t x)
voidSetDriftHisto(Float_t x, Int_t = 200)
intKGeometry::SetElecVolt(int i)
voidSetEntryPoint(Float_t x, Float_t y, Float_t z)
voidSetExitPoint(Float_t x, Float_t y, Float_t z)
voidSetNeff(TF3* neff, Int_t calnow = 1)
voidSetNeff(TH3F* neff, Int_t calnow = 1)
voidSetPrecision(Double_t x)
voidSetVoltage(Float_t x, Int_t calnow = 1)
voidShowGaussBeam(Int_t, Float_t, Float_t, Float_t, Int_t, Int_t)
virtual voidShowMembers(TMemberInspector&)
voidShowMipIR(Int_t, Int_t = 14, Int_t = 1)
virtual voidStreamer(TBuffer&)
voidStreamerNVirtual(TBuffer& ClassDef_StreamerNVirtual_b)
Double_tV(int, int)

Data Members

public:
Float_tB[3]magnetic field
Float_tBDTreshhole multiplication - break down treshold
Int_tBreakDownif break down occurs it goes to 1 otherwise is 0
TH3I*KGeometry::DMdetector material
TH3I*KGeometry::EGelectrode geometry
static Int_tKMaterial::ImpactIonizationimpact ionization model
Float_tMTreshtreshold for taking multiplication into account
static Int_tKMaterial::MatMaterial index
static Int_tKMaterial::Mobilitymobility model for each material
TF3*NeffFeffective dopping concentration function
TH3F*NeffHeffective dopping concentration histogram
KField*Ramoramo field
KField*Realelectric field
Float_tSStepSimulation step size;
static Float_tKMaterial::TemperatureTemperature
Float_tVoltageVoltage
Float_tVoltage2Voltage2
TArrayFVoltagesArray of voltages
Int_taverageAverage (over how many events)
Int_tdiffDiffusion simulation (yes=1, no=0)
Float_tenp[3]entry point for the charge drift
Float_texp[3]exit point for the cahrge drift
TH1F*negcontribution of the electrons to the total drift current
Int_tKGeometry::nxx-divisions
Int_tKGeometry::nyy-divisions
Int_tKGeometry::nzz-divisions
TH1F*poscontribution of the holes to the total drift current
TH1F*sumtotal drift current
Float_ttaueeffective trapping time constants - used if Multiplication is ON
Float_ttauheffective trapping time constants - used if Multiplication is ON
private:
Double_tCalErrError of the solver
Short_tDebugPrint information of drift calculation etc.
Double_tDeps
Int_tMaxIterMaximum number of iterations in eq solver
TRandom*ranrandom number generator

Class Charts

Inheritance Chart:
KGeometry
KMaterial
KDetector
K3D
KPad
KPixel
KStrip

Function documentation

double V(int , int )
Double_t kappa(int , int , int , int )
Sets the effective space charge values for given point in the mesh!
void Declaration(Int_t )
 New declaration for a general detector class
void CalField(Int_t )
void Drift(Double_t , Double_t , Double_t , Float_t , KStruct* , Double_t = 0)
Drift simulation for a point charge (Float_t charg;)
starting from ( sx,sy, sz)
KStruct *seg is the structure where the  the drift paths, drift times and induced cahrges are stored
TH1F * Draw1D(Char_t* , Float_t , Int_t , Float_t )
 Draws the 1D projection of the Field
 Char_t option;  see ::Draw()
 Char_t proj;    see ::Draw()
 Int_t axis;     0=x, 1=y, 2=z;
 Float_t pos;    position along the missing coordinate if 2D z=0.5;
TH2F * Draw(Char_t* , Float_t = 1)
The function draws weighting and electric field and also the event display
 Float_t proj; position along the axis of projection
 Char_t *option:
           Which field do you want to plot?
		  W -weighting
 		  E -electric
              G -geometry
              M -material
              N -Neff = space charge
           What do do want to plot?
		  P - potential
 		  F - |field|
              X - x component of the field
              Y - y component of the field
              Z - z component of the field
          In case of 3D which plane?
              yz - the cross section is along the x value of proj
              xz - the cross section is along the y value of proj
              xy - the cross section is along the z value of proj
void MipIR(Int_t = 20, Float_t = 0)
 The simulation of the drift for the laser induced carriers - attenuation of the light.
 A track is devided into Int_ div buckets. Each bucket is drifted in the field. The
 induced currents for each carrier is calculated as the sum  all buckets.
	Int_t MobMod; mobility model
	Float_t B; magnetic field
void ShowMipIR(Int_t , Int_t = 14, Int_t = 1)
 The simulation of the drift for the minimum ionizing particles.
 A track is devided into Int_ div buckets. Each bucket is drifted in the field. The
 induced currents for each carrier is calculated as the sum  all buckets.
KDetector()
Author: Gregor Kramberger
Default constructor for KDetector class
Default = no space charge -> default space charge is function

~KDetector()
Destructor of the detector class
void CalM(KStruct* seg, Double_t* data, Int_t = -1)
void SetDriftHisto(Float_t x, Int_t = 200)
void Save(Char_t* , Char_t* )
TFile * Read(Char_t* , Char_t* )
void GaussBeam(Int_t , Float_t , Float_t , Float_t , Float_t , Float_t )
The simulation of the drift for Gaussian Beam.
A track is divided into Int_ div buckets. Each bucket is drifted in the field.
The induced currents for each carrier is calculated as the sum all buckets.

	- Int_t MobMod; mobility model                                                  
	- Float_t B; magnetic field                                                     

- Lambda [um]
- w0 [um]


void ShowGaussBeam(Int_t , Float_t , Float_t , Float_t , Int_t , Int_t )
Gaussian Beam for strip detector

A track is divided into Int_ div buckets. Each bucket is drifted in the field.
The induced currents for each carrier is calculated as the sum all buckets.

- Lambda [um]
- w0 [um]

KDetector()
 Constructors and destructor
void ResetRnd(Int_t seed)
Configuration functions
{delete ran; ran=new TRandom(seed);}
void SetCalculationParameters(Double_t x, Int_t y)
{CalErr=x; MaxIter=y;}
void CalPhyField()
{CalField(0);}
void CalRamoField()
{CalField(1);}
void SetVoltage(Float_t x, Int_t calnow = 1)
 Calculation in case of any changes
{Voltage=x; if(calnow) CalPhyField(); }
void SetNeff(TF3* neff, Int_t calnow = 1)
{NeffF=neff; if(calnow) CalPhyField(); }
void SetNeff(TH3F* neff, Int_t calnow = 1)
{NeffH=neff; if(calnow) CalPhyField(); }
void SetEntryPoint(Float_t x, Float_t y, Float_t z)
 Simulation of drift
{enp[0]=x; enp[1]=y; enp[2]=z;}
void SetExitPoint(Float_t x, Float_t y, Float_t z)
{exp[0]=x; exp[1]=y; exp[2]=z;}
void SetDebug(Short_t x)
 precision of drift
{Debug=x;}
void SetPrecision(Double_t x)
{Deps=x;}
Double_t GetPrecision()
The main detector class = the basis for simulation
{return Deps;}