Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4MicroElecSurface Class Reference

#include <G4MicroElecSurface.hh>

+ Inheritance diagram for G4MicroElecSurface:

Public Member Functions

 G4MicroElecSurface (const G4String &processName="MicroElecSurface", G4ProcessType type=fElectromagnetic)
 
 ~G4MicroElecSurface () override
 
G4bool IsApplicable (const G4ParticleDefinition &aParticleType) override
 
void SetFlagFranchissement ()
 
G4double GetMeanFreePath (const G4Track &, G4double, G4ForceCondition *condition) override
 
G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep) override
 
void BuildPhysicsTable (const G4ParticleDefinition &) override
 
G4MicroElecSurfaceStatus GetStatus () const
 
 G4MicroElecSurface (const G4MicroElecSurface &right)=delete
 
G4MicroElecSurfaceoperator= (const G4MicroElecSurface &right)=delete
 
- Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &aName, G4ProcessType aType=fNotDefined)
 
 G4VDiscreteProcess (G4VDiscreteProcess &)
 
virtual ~G4VDiscreteProcess ()
 
G4VDiscreteProcessoperator= (const G4VDiscreteProcess &)=delete
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4VProcessoperator= (const G4VProcess &)=delete
 
G4bool operator== (const G4VProcess &right) const
 
G4bool operator!= (const G4VProcess &right) const
 
virtual G4VParticleChangePostStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)=0
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)=0
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)=0
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void StartTracking (G4Track *)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
virtual void ProcessDescription (std::ostream &outfile) const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)=0
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double prevStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager = nullptr
 
G4VParticleChangepParticleChange = nullptr
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft = -1.0
 
G4double currentInteractionLength = -1.0
 
G4double theInitialNumberOfInteractionLength = -1.0
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType = fNotDefined
 
G4int theProcessSubType = -1
 
G4double thePILfactor = 1.0
 
G4int verboseLevel = 0
 
G4bool enableAtRestDoIt = true
 
G4bool enableAlongStepDoIt = true
 
G4bool enablePostStepDoIt = true
 

Detailed Description

Definition at line 85 of file G4MicroElecSurface.hh.

Constructor & Destructor Documentation

◆ G4MicroElecSurface() [1/2]

G4MicroElecSurface::G4MicroElecSurface ( const G4String processName = "MicroElecSurface",
G4ProcessType  type = fElectromagnetic 
)

Definition at line 63 of file G4MicroElecSurface.cc.

64 : G4VDiscreteProcess(processName, type),
65 oldMomentum(0.,0.,0.), previousMomentum(0.,0.,0.),
66 theGlobalNormal(0.,0.,0.), theFacetNormal(0.,0.,0.)
67{
68 if ( verboseLevel > 0)
69 {
70 G4cout << GetProcessName() << " is created " << G4endl;
71 }
72
73 isInitialised=false;
75
76 theStatus = UndefinedSurf;
77 material1 = nullptr;
78 material2 = nullptr;
79
81 theParticleMomentum = 0.;
82
83 flag_franchissement_surface = false;
84 flag_normal = false;
85 flag_reflexion = false;
86 teleportToDo = teleportDone = false;
87
88 ekint = thetat = thetaft = energyThreshold = crossingProbability = 0.0;
89}
@ UndefinedSurf
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
G4double GetSurfaceTolerance() const
static G4GeometryTolerance * GetInstance()
G4int verboseLevel
Definition: G4VProcess.hh:356
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:406
const G4String & GetProcessName() const
Definition: G4VProcess.hh:382

◆ ~G4MicroElecSurface()

G4MicroElecSurface::~G4MicroElecSurface ( )
override

Definition at line 93 of file G4MicroElecSurface.cc.

94{}

◆ G4MicroElecSurface() [2/2]

G4MicroElecSurface::G4MicroElecSurface ( const G4MicroElecSurface right)
delete

Member Function Documentation

◆ BuildPhysicsTable()

void G4MicroElecSurface::BuildPhysicsTable ( const G4ParticleDefinition )
overridevirtual

Reimplemented from G4VProcess.

Definition at line 106 of file G4MicroElecSurface.cc.

107{
108 if (isInitialised) { return; }
109
110 G4ProductionCutsTable* theCoupleTable =
112 G4int numOfCouples = theCoupleTable->GetTableSize();
113 G4cout << "G4MicroElecSurface::Initialise: Ncouples= "
114 << numOfCouples << G4endl;
115
116 for (G4int i = 0; i < numOfCouples; ++i) {
117 const G4Material* material =
118 theCoupleTable->GetMaterialCutsCouple(i)->GetMaterial();
119
120 G4cout << "G4Surface, Material " << i + 1 << " / " << numOfCouples << " : " << material->GetName() << G4endl;
121 if (material->GetName() == "Vacuum") { tableWF[material->GetName()] = 0; continue; }
122 G4String mat = material->GetName();
124 tableWF[mat] = str.GetWorkFunction();
125 }
126 isInitialised = true;
127}
int G4int
Definition: G4Types.hh:85
const G4Material * GetMaterial() const
const G4String & GetName() const
Definition: G4Material.hh:175
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
std::size_t GetTableSize() const
static G4ProductionCutsTable * GetProductionCutsTable()

◆ GetMeanFreePath()

G4double G4MicroElecSurface::GetMeanFreePath ( const G4Track ,
G4double  ,
G4ForceCondition condition 
)
overridevirtual

Implements G4VDiscreteProcess.

Definition at line 393 of file G4MicroElecSurface.cc.

395{
396 *condition = Forced;
397 return DBL_MAX;
398}
G4double condition(const G4ErrorSymMatrix &m)
@ Forced
#define DBL_MAX
Definition: templates.hh:62

◆ GetStatus()

G4MicroElecSurfaceStatus G4MicroElecSurface::GetStatus ( ) const

Definition at line 479 of file G4MicroElecSurface.cc.

480{
481 return theStatus;
482}

◆ IsApplicable()

G4bool G4MicroElecSurface::IsApplicable ( const G4ParticleDefinition aParticleType)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 99 of file G4MicroElecSurface.cc.

100{
101 return ( aParticleType.GetPDGEncoding() == 11 );
102}

◆ operator=()

G4MicroElecSurface & G4MicroElecSurface::operator= ( const G4MicroElecSurface right)
delete

◆ PostStepDoIt()

G4VParticleChange * G4MicroElecSurface::PostStepDoIt ( const G4Track aTrack,
const G4Step aStep 
)
overridevirtual

Reimplemented from G4VDiscreteProcess.

Definition at line 131 of file G4MicroElecSurface.cc.

132{
133 theStatus = UndefinedSurf;
134
135 //Definition of the parameters for the particle
138
139 G4StepPoint* pPreStepPoint = aStep.GetPreStepPoint();
140 G4StepPoint* pPostStepPoint = aStep.GetPostStepPoint();
141
142 material1 = pPreStepPoint -> GetMaterial();
143 material2 = pPostStepPoint -> GetMaterial();
144
145 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
146
147 theParticleMomentum = aParticle->GetTotalMomentum();
148 previousMomentum = oldMomentum;
149 oldMomentum = aParticle->GetMomentumDirection();
150
151 //First case: not a boundary
152 if (pPostStepPoint->GetStepStatus() != fGeomBoundary ||
153 pPostStepPoint->GetPhysicalVolume() == pPreStepPoint->GetPhysicalVolume())
154 {
155 theStatus = NotAtBoundarySurf;
156 flag_franchissement_surface = false;
157 flag_reflexion = false;
158 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
159
160 }
161 theStatus = UndefinedSurf;
162
163 //Third case: same material
164 if (material1 == material2)
165 {
166 theStatus = SameMaterialSurf;
167 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
168
169 }
170 if (verboseLevel > 0)
171 {
172 G4cout << G4endl << " Electron at Boundary! " << G4endl;
173 G4VPhysicalVolume* thePrePV = pPreStepPoint->GetPhysicalVolume();
174 G4VPhysicalVolume* thePostPV = pPostStepPoint->GetPhysicalVolume();
175 if (thePrePV) G4cout << " thePrePV: " << thePrePV->GetName() << G4endl;
176 if (thePostPV) G4cout << " thePostPV: " << thePostPV->GetName() << G4endl;
177 G4cout << " Old Momentum Direction: " << oldMomentum << G4endl;
178 }
179
180 //Definition of the parameters for the surface
181 G4ThreeVector theGlobalPoint = pPostStepPoint->GetPosition();
182
183 G4Navigator* theNavigator =
185 GetNavigatorForTracking();
186
187 G4bool valid;
188
189 theGlobalNormal = theNavigator->GetGlobalExitNormal(theGlobalPoint, &valid);
190 // G4cout << "Global exit normal = " << theGlobalNormal << " valid = " << valid << G4endl;
191
192 if (valid)
193 {
194 theGlobalNormal = -theGlobalNormal;
195 }
196 else
197 {
199 ed << " G4MicroElecSurface/PostStepDoIt(): "
200 << " The Navigator reports that it returned an invalid normal.\n"
201 << "PV: " << pPreStepPoint->GetPhysicalVolume()->GetName()
202 << " TrackID= " << aTrack.GetTrackID()
203 << " Ekin(MeV)= " << aTrack.GetKineticEnergy()
204 << " position: " << theGlobalPoint
205 << " direction: " << oldMomentum
206 << G4endl;
207 G4Exception("G4MuElecSurf::PostStepDoIt", "OpBoun01",
208 FatalException, ed,
209 "Invalid Surface Normal - Geometry must return valid surface normal");
210 return 0;
211 }
212
213 //Exception: the particle is not in the right direction
214 if (oldMomentum * theGlobalNormal > 0.0)
215 {
216 theGlobalNormal = -theGlobalNormal;
217 }
218
219 //Second case: step too small
220 //Corrections bug rotation + réflexion
221 if (aTrack.GetStepLength()<=kCarTolerance)
222 {
223 theStatus = StepTooSmallSurf;
224 if (pPostStepPoint) {
225
226 WorkFunctionTable::iterator postStepWF;
227 postStepWF = tableWF.find(pPostStepPoint->GetMaterial()->GetName());
228 WorkFunctionTable::iterator preStepWF;
229 preStepWF = tableWF.find(pPreStepPoint->GetMaterial()->GetName());
230
231 if (postStepWF == tableWF.end()) {
232 G4String str = "Material ";
233 str += pPostStepPoint->GetMaterial()->GetName() + " not found!";
234 G4Exception("G4Surface::G4Surface", "em0002", FatalException, str);
235 return 0;
236 }
237
238 else if (preStepWF == tableWF.end()) {
239 G4String str = "Material ";
240 str += pPreStepPoint->GetMaterial()->GetName() + " not found!";
241 G4Exception("G4Surface::G4Surface", "em0002", FatalException, str);
242 return 0;
243 }
244
245 if (pPreStepPoint->GetMaterial() != pPostStepPoint->GetMaterial()) {
246
247 flag_franchissement_surface = false;
248
249 if (flag_reflexion == true && flag_normal == true) {
251 flag_reflexion = false;
252 flag_normal = false;
253 }
254 }
255 }
256 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
257 }
258
259 flag_normal = (theGlobalNormal.x() == 0.0 && theGlobalNormal.y() == 0.0);
260
261 G4LogicalSurface* Surface = nullptr;
262
264 (pPreStepPoint ->GetPhysicalVolume(),
265 pPostStepPoint->GetPhysicalVolume());
266
267 if (Surface == nullptr)
268 {
269 G4bool enteredDaughter=(pPostStepPoint->GetPhysicalVolume()
270 ->GetMotherLogical() ==
271 pPreStepPoint->GetPhysicalVolume()
272 ->GetLogicalVolume());
273 if(enteredDaughter)
274 {
276 (pPostStepPoint->GetPhysicalVolume()->
277 GetLogicalVolume());
278
279 if(Surface == nullptr)
281 (pPreStepPoint->GetPhysicalVolume()->
282 GetLogicalVolume());
283 }
284 else
285 {
287 (pPreStepPoint->GetPhysicalVolume()->
288 GetLogicalVolume());
289
290 if(Surface == nullptr)
292 (pPostStepPoint->GetPhysicalVolume()->
293 GetLogicalVolume());
294 }
295 }
296
297
298 G4VPhysicalVolume* thePrePV = pPreStepPoint->GetPhysicalVolume();
299 G4VPhysicalVolume* thePostPV = pPostStepPoint->GetPhysicalVolume();
300
301 if (thePostPV)
302 {
303 WorkFunctionTable::iterator postStepWF;
304 postStepWF = tableWF.find(thePostPV->GetLogicalVolume()->GetMaterial()->GetName());
305 WorkFunctionTable::iterator preStepWF;
306 preStepWF = tableWF.find(thePrePV->GetLogicalVolume()->GetMaterial()->GetName());
307
308 if (postStepWF == tableWF.end()) {
309 G4String str = "Material ";
310 str += thePostPV->GetLogicalVolume()->GetMaterial()->GetName() + " not found!";
311 G4Exception("G4Surface::G4Surface", "em0002", FatalException, str);
312 return 0;
313 }
314
315 else if (preStepWF == tableWF.end()) {
316 G4String str = "Material ";
317 str += thePrePV->GetLogicalVolume()->GetMaterial()->GetName() + " not found!";
318 G4Exception("G4Surface::G4Surface", "em0002", FatalException, str);
319 return 0;
320 }
321
322 else
323 {
324 G4double thresholdNew = postStepWF->second;
325 G4double thresholdOld = preStepWF->second;
326 energyThreshold = thresholdNew - thresholdOld;
327 }
328 }
329
330 ekint = pPreStepPoint->GetKineticEnergy();
331 thetat= GetIncidentAngle(); //angle d'incidence
332 G4double ekinNormalt=ekint*std::cos(thetat)*std::cos(thetat);
333
334 G4double atet = std::sqrt(ekint/(ekint+energyThreshold))*std::sin(thetat);
335 thetaft = (atet > 1.0) ? pi*0.5 : std::asin(atet);//Angle de réfraction
336
337 G4double aleat=G4UniformRand();
338
339 const G4double waveVectort=std::sqrt(2*9.1093826E-31*1.602176487E-19)/(6.6260755E-34/(2.0*pi));
340
341 //Parameter for an exponential barrier of potential (Thèse P68)
342 const G4double at=0.5E-10;
343
344 //G4double modif already declared in .hh
345 crossingProbability=0;
346
347 G4double kft=waveVectort*std::sqrt(ekint+energyThreshold)*std::cos(thetaft);
348 G4double kit=waveVectort*std::sqrt(ekinNormalt);
349
350 G4double yy = std::sinh(pi*at*(kit-kft))/std::sinh(pi*at*(kit+kft));
351 crossingProbability = 1 - yy*yy;
352
353 //First case: the electron crosses the surface
354 if((aleat<=crossingProbability)&&(ekint>std::abs(energyThreshold)))
355 {
356 if (pPreStepPoint->GetMaterial() != pPostStepPoint->GetMaterial()) {
357 flag_franchissement_surface = true;
358 }
359
360 thetaft=std::abs(thetaft-thetat);
361
363 G4ThreeVector xVerst = zVerst.orthogonal();
364 G4ThreeVector yVerst = zVerst.cross(xVerst);
365
366 G4double cost = std::cos(thetaft);
367 G4double xDirt = std::sqrt(1. - cost*cost);
368 G4double yDirt = xDirt;
369
370 G4ThreeVector zPrimeVerst = xDirt*xVerst + yDirt*yVerst + cost*zVerst;
371
373 }
374 else if ((aleat > crossingProbability) && (ekint>std::abs(energyThreshold)))
375 {
376 flag_reflexion = true;
377 if (flag_normal) { aParticleChange.ProposeMomentumDirection(-oldMomentum.unit()); }
378 else { aParticleChange.ProposeMomentumDirection(Reflexion(aStep.GetPostStepPoint())); }
379
380 }
381 else {
382
383 if (flag_normal) { aParticleChange.ProposeMomentumDirection(-oldMomentum.unit()); }
384 else { aParticleChange.ProposeMomentumDirection(Reflexion(aStep.GetPostStepPoint())); }
385 flag_reflexion = true;
386 }
387
388 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
389}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
@ StepTooSmallSurf
@ SameMaterialSurf
@ NotAtBoundarySurf
@ fGeomBoundary
Definition: G4StepStatus.hh:43
double G4double
Definition: G4Types.hh:83
bool G4bool
Definition: G4Types.hh:86
#define G4UniformRand()
Definition: Randomize.hh:52
Hep3Vector unit() const
Hep3Vector orthogonal() const
double x() const
double y() const
Hep3Vector cross(const Hep3Vector &) const
const G4ThreeVector & GetMomentumDirection() const
G4double GetTotalMomentum() const
static G4LogicalBorderSurface * GetSurface(const G4VPhysicalVolume *vol1, const G4VPhysicalVolume *vol2)
static G4LogicalSkinSurface * GetSurface(const G4LogicalVolume *vol)
G4Material * GetMaterial() const
virtual G4ThreeVector GetGlobalExitNormal(const G4ThreeVector &point, G4bool *valid)
void ProposeVelocity(G4double finalVelocity)
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
virtual void Initialize(const G4Track &)
G4StepStatus GetStepStatus() const
G4Material * GetMaterial() const
const G4ThreeVector & GetPosition() const
const G4ThreeVector & GetMomentumDirection() const
G4VPhysicalVolume * GetPhysicalVolume() const
G4double GetKineticEnergy() const
G4StepPoint * GetPreStepPoint() const
G4StepPoint * GetPostStepPoint() const
G4double GetVelocity() const
G4int GetTrackID() const
const G4DynamicParticle * GetDynamicParticle() const
G4double GetKineticEnergy() const
G4double GetStepLength() const
static G4TransportationManager * GetTransportationManager()
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
G4LogicalVolume * GetMotherLogical() const
G4LogicalVolume * GetLogicalVolume() const
const G4String & GetName() const
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:327
const G4double pi

◆ SetFlagFranchissement()

void G4MicroElecSurface::SetFlagFranchissement ( )

Definition at line 487 of file G4MicroElecSurface.cc.

488{
489 flag_franchissement_surface = false;
490}

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