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

#include <G4EmExtraParameters.hh>

Public Member Functions

 G4EmExtraParameters ()
 
 ~G4EmExtraParameters ()
 
void Initialise ()
 
G4bool GetDirectionalSplitting ()
 
void SetDirectionalSplitting (G4bool v)
 
G4bool QuantumEntanglement ()
 
void SetQuantumEntanglement (G4bool v)
 
void SetDirectionalSplittingRadius (G4double r)
 
G4double GetDirectionalSplittingRadius ()
 
void SetDirectionalSplittingTarget (const G4ThreeVector &v)
 
G4ThreeVector GetDirectionalSplittingTarget () const
 
void SetStepFunction (G4double v1, G4double v2)
 
G4double GetStepFunctionP1 () const
 
G4double GetStepFunctionP2 () const
 
void SetStepFunctionMuHad (G4double v1, G4double v2)
 
G4double GetStepFunctionMuHadP1 () const
 
G4double GetStepFunctionMuHadP2 () const
 
void SetStepFunctionLightIons (G4double v1, G4double v2)
 
G4double GetStepFunctionLightIonsP1 () const
 
G4double GetStepFunctionLightIonsP2 () const
 
void SetStepFunctionIons (G4double v1, G4double v2)
 
G4double GetStepFunctionIonsP1 () const
 
G4double GetStepFunctionIonsP2 () const
 
void FillStepFunction (const G4ParticleDefinition *, G4VEnergyLossProcess *) const
 
void AddPAIModel (const G4String &particle, const G4String &region, const G4String &type)
 
const std::vector< G4String > & ParticlesPAI () const
 
const std::vector< G4String > & RegionsPAI () const
 
const std::vector< G4String > & TypesPAI () const
 
void AddPhysics (const G4String &region, const G4String &type)
 
const std::vector< G4String > & RegionsPhysics () const
 
const std::vector< G4String > & TypesPhysics () const
 
void SetSubCutRegion (const G4String &region)
 
void SetProcessBiasingFactor (const G4String &procname, G4double val, G4bool wflag)
 
void ActivateForcedInteraction (const G4String &procname, const G4String &region, G4double length, G4bool wflag)
 
void ActivateSecondaryBiasing (const G4String &name, const G4String &region, G4double factor, G4double energyLimit)
 
void DefineRegParamForLoss (G4VEnergyLossProcess *) const
 
void DefineRegParamForEM (G4VEmProcess *) const
 
 G4EmExtraParameters (G4EmExtraParameters &)=delete
 
G4EmExtraParametersoperator= (const G4EmExtraParameters &right)=delete
 

Detailed Description

Definition at line 62 of file G4EmExtraParameters.hh.

Constructor & Destructor Documentation

◆ G4EmExtraParameters() [1/2]

G4EmExtraParameters::G4EmExtraParameters ( )
explicit

◆ ~G4EmExtraParameters()

G4EmExtraParameters::~G4EmExtraParameters ( )

Definition at line 63 of file G4EmExtraParameters.cc.

64{
65 delete theMessenger;
66}

◆ G4EmExtraParameters() [2/2]

G4EmExtraParameters::G4EmExtraParameters ( G4EmExtraParameters & )
delete

Member Function Documentation

◆ ActivateForcedInteraction()

void G4EmExtraParameters::ActivateForcedInteraction ( const G4String & procname,
const G4String & region,
G4double length,
G4bool wflag )

Definition at line 323 of file G4EmExtraParameters.cc.

327{
328 const G4String& r = CheckRegion(region);
329 if(length >= 0.0) {
330 std::size_t n = m_procForced.size();
331 for(std::size_t i=0; i<n; ++i) {
332 if(procname == m_procForced[i] && r == m_regnamesForced[i] ) {
333 m_lengthForced[i] = length;
334 m_weightForced[i] = wflag;
335 return;
336 }
337 }
338 m_regnamesForced.push_back(r);
339 m_procForced.push_back(procname);
340 m_lengthForced.push_back(length);
341 m_weightForced.push_back(wflag);
342 } else {
344 ed << "Process: " << procname << " in region " << r
345 << " : forced interacttion length= "
346 << length << " is negative - ignored";
347 PrintWarning(ed);
348 }
349}
std::ostringstream G4ExceptionDescription

Referenced by G4EmParameters::ActivateForcedInteraction(), and G4EmExtraParametersMessenger::SetNewValue().

◆ ActivateSecondaryBiasing()

void G4EmExtraParameters::ActivateSecondaryBiasing ( const G4String & name,
const G4String & region,
G4double factor,
G4double energyLimit )

Definition at line 352 of file G4EmExtraParameters.cc.

356{
357 const G4String& r = CheckRegion(region);
358 if(factor >= 0.0 && energyLim >= 0.0) {
359 std::size_t n = m_procBiasedSec.size();
360 for(std::size_t i=0; i<n; ++i) {
361 if(procname == m_procBiasedSec[i] && r == m_regnamesBiasedSec[i] ) {
362 m_factBiasedSec[i] = factor;
363 m_elimBiasedSec[i] = energyLim;
364 return;
365 }
366 }
367 m_regnamesBiasedSec.push_back(r);
368 m_procBiasedSec.push_back(procname);
369 m_factBiasedSec.push_back(factor);
370 m_elimBiasedSec.push_back(energyLim);
371 } else {
373 ed << "Process: " << procname << " in region " << r
374 << " : secondary bised factor= "
375 << factor << ", Elim= " << energyLim << " - ignored";
376 PrintWarning(ed);
377 }
378}

Referenced by G4EmParameters::ActivateSecondaryBiasing(), and G4EmExtraParametersMessenger::SetNewValue().

◆ AddPAIModel()

void G4EmExtraParameters::AddPAIModel ( const G4String & particle,
const G4String & region,
const G4String & type )

Definition at line 222 of file G4EmExtraParameters.cc.

225{
226 G4String r = CheckRegion(region);
227 std::size_t nreg = m_regnamesPAI.size();
228
229 // in previously defined region other particles may be already defined
230 // type should be overrided for the same region and particle
231 for(std::size_t i=0; i<nreg; ++i) {
232 if(m_regnamesPAI[i] == r) {
233 if (particle == "all") {
234 m_particlesPAI[i] = particle;
235 m_typesPAI[i] = type;
236 return;
237 } else if(m_particlesPAI[i] == particle || m_particlesPAI[i] == "all") {
238 m_typesPAI[i] = type;
239 return;
240 }
241 }
242 }
243 // new regions and/or particles
244 m_particlesPAI.push_back(particle);
245 m_regnamesPAI.push_back(r);
246 m_typesPAI.push_back(type);
247}

Referenced by G4EmParameters::AddPAIModel(), and G4EmExtraParametersMessenger::SetNewValue().

◆ AddPhysics()

void G4EmExtraParameters::AddPhysics ( const G4String & region,
const G4String & type )

Definition at line 264 of file G4EmExtraParameters.cc.

266{
267 G4String r = CheckRegion(region);
268 std::size_t nreg = m_regnamesPhys.size();
269 for(std::size_t i=0; i<nreg; ++i) {
270 if(r == m_regnamesPhys[i]) { return; }
271 }
272 m_regnamesPhys.push_back(r);
273 m_typesPhys.push_back(type);
274}

Referenced by G4EmParameters::AddPhysics(), and G4EmExtraParametersMessenger::SetNewValue().

◆ DefineRegParamForEM()

void G4EmExtraParameters::DefineRegParamForEM ( G4VEmProcess * ptr) const

Definition at line 416 of file G4EmExtraParameters.cc.

417{
418 std::size_t n = m_procBiasedXS.size();
419 for(std::size_t i=0; i<n; ++i) {
420 if(ptr->GetProcessName() == m_procBiasedXS[i]) {
421 ptr->SetCrossSectionBiasingFactor(m_factBiasedXS[i],
422 m_weightBiasedXS[i]);
423 break;
424 }
425 }
426 n = m_procForced.size();
427 for(std::size_t i=0; i<n; ++i) {
428 if(ptr->GetProcessName() == m_procForced[i]) {
429 ptr->ActivateForcedInteraction(m_lengthForced[i],
430 m_regnamesForced[i],
431 m_weightForced[i]);
432 break;
433 }
434 }
435 n = m_procBiasedSec.size();
436 for(std::size_t i=0; i<n; ++i) {
437 if(ptr->GetProcessName() == m_procBiasedSec[i]) {
438 ptr->ActivateSecondaryBiasing(m_regnamesBiasedSec[i],
439 m_factBiasedSec[i],
440 m_elimBiasedSec[i]);
441 break;
442 }
443 }
444}
void SetCrossSectionBiasingFactor(G4double f, G4bool flag=true)
void ActivateForcedInteraction(G4double length=0.0, const G4String &r="", G4bool flag=true)
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
const G4String & GetProcessName() const

Referenced by G4EmParameters::DefineRegParamForEM().

◆ DefineRegParamForLoss()

void G4EmExtraParameters::DefineRegParamForLoss ( G4VEnergyLossProcess * ptr) const

Definition at line 380 of file G4EmExtraParameters.cc.

381{
382 const G4RegionStore* regionStore = G4RegionStore::GetInstance();
383 std::size_t n = m_regnamesSubCut.size();
384 for(std::size_t i=0; i<n; ++i) {
385 const G4Region* reg = regionStore->GetRegion(m_regnamesSubCut[i], false);
386 if(nullptr != reg) { ptr->ActivateSubCutoff(reg); }
387 }
388 n = m_procBiasedXS.size();
389 for(std::size_t i=0; i<n; ++i) {
390 if(ptr->GetProcessName() == m_procBiasedXS[i]) {
391 ptr->SetCrossSectionBiasingFactor(m_factBiasedXS[i],
392 m_weightBiasedXS[i]);
393 break;
394 }
395 }
396 n = m_procForced.size();
397 for(std::size_t i=0; i<n; ++i) {
398 if(ptr->GetProcessName() == m_procForced[i]) {
399 ptr->ActivateForcedInteraction(m_lengthForced[i],
400 m_regnamesForced[i],
401 m_weightForced[i]);
402 break;
403 }
404 }
405 n = m_procBiasedSec.size();
406 for(std::size_t i=0; i<n; ++i) {
407 if(ptr->GetProcessName() == m_procBiasedSec[i]) {
408 ptr->ActivateSecondaryBiasing(m_regnamesBiasedSec[i],
409 m_factBiasedSec[i],
410 m_elimBiasedSec[i]);
411 break;
412 }
413 }
414}
static G4RegionStore * GetInstance()
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
void SetCrossSectionBiasingFactor(G4double f, G4bool flag=true)
void ActivateForcedInteraction(G4double length, const G4String &region, G4bool flag=true)
void ActivateSubCutoff(const G4Region *region)

Referenced by G4EmParameters::DefineRegParamForLoss().

◆ FillStepFunction()

void G4EmExtraParameters::FillStepFunction ( const G4ParticleDefinition * part,
G4VEnergyLossProcess * proc ) const

Definition at line 202 of file G4EmExtraParameters.cc.

203{
204 // electron and positron
205 if (11 == std::abs(part->GetPDGEncoding())) {
206 proc->SetStepFunction(dRoverRange, finalRange);
207
208 // all heavy ions
209 } else if ("GenericIon" == part->GetParticleName()) {
210 proc->SetStepFunction(dRoverRangeIons, finalRangeIons);
211
212 // light nucleus and anti-nucleus
213 } else if (part->GetParticleType() == "nucleus" || part->GetParticleType() == "anti_nucleus") {
214 proc->SetStepFunction(dRoverRangeLIons, finalRangeLIons);
215
216 // other particles
217 } else {
218 proc->SetStepFunction(dRoverRangeMuHad, finalRangeMuHad);
219 }
220}
const G4String & GetParticleType() const
const G4String & GetParticleName() const
void SetStepFunction(G4double v1, G4double v2)

Referenced by G4EmParameters::FillStepFunction().

◆ GetDirectionalSplitting()

G4bool G4EmExtraParameters::GetDirectionalSplitting ( )

Definition at line 456 of file G4EmExtraParameters.cc.

456 {
457 return directionalSplitting;
458}

Referenced by G4EmParameters::GetDirectionalSplitting().

◆ GetDirectionalSplittingRadius()

G4double G4EmExtraParameters::GetDirectionalSplittingRadius ( )

Definition at line 481 of file G4EmExtraParameters.cc.

482{
483 return directionalSplittingRadius;
484}

Referenced by G4EmParameters::GetDirectionalSplittingRadius().

◆ GetDirectionalSplittingTarget()

G4ThreeVector G4EmExtraParameters::GetDirectionalSplittingTarget ( ) const

Definition at line 471 of file G4EmExtraParameters.cc.

472{
473 return directionalSplittingTarget;
474}

Referenced by G4EmParameters::GetDirectionalSplittingTarget().

◆ GetStepFunctionIonsP1()

G4double G4EmExtraParameters::GetStepFunctionIonsP1 ( ) const

Definition at line 192 of file G4EmExtraParameters.cc.

193{
194 return dRoverRangeIons;
195}

Referenced by G4EmParameters::StreamInfo().

◆ GetStepFunctionIonsP2()

G4double G4EmExtraParameters::GetStepFunctionIonsP2 ( ) const

Definition at line 197 of file G4EmExtraParameters.cc.

198{
199 return finalRangeIons;
200}

Referenced by G4EmParameters::StreamInfo().

◆ GetStepFunctionLightIonsP1()

G4double G4EmExtraParameters::GetStepFunctionLightIonsP1 ( ) const

Definition at line 169 of file G4EmExtraParameters.cc.

170{
171 return dRoverRangeLIons;
172}

Referenced by G4EmParameters::StreamInfo().

◆ GetStepFunctionLightIonsP2()

G4double G4EmExtraParameters::GetStepFunctionLightIonsP2 ( ) const

Definition at line 174 of file G4EmExtraParameters.cc.

175{
176 return finalRangeLIons;
177}

Referenced by G4EmParameters::StreamInfo().

◆ GetStepFunctionMuHadP1()

G4double G4EmExtraParameters::GetStepFunctionMuHadP1 ( ) const

Definition at line 146 of file G4EmExtraParameters.cc.

147{
148 return dRoverRangeMuHad;
149}

Referenced by G4EmParameters::StreamInfo().

◆ GetStepFunctionMuHadP2()

G4double G4EmExtraParameters::GetStepFunctionMuHadP2 ( ) const

Definition at line 151 of file G4EmExtraParameters.cc.

152{
153 return finalRangeMuHad;
154}

Referenced by G4EmParameters::StreamInfo().

◆ GetStepFunctionP1()

G4double G4EmExtraParameters::GetStepFunctionP1 ( ) const

Definition at line 123 of file G4EmExtraParameters.cc.

124{
125 return dRoverRange;
126}

Referenced by G4EmParameters::StreamInfo().

◆ GetStepFunctionP2()

G4double G4EmExtraParameters::GetStepFunctionP2 ( ) const

Definition at line 128 of file G4EmExtraParameters.cc.

129{
130 return finalRange;
131}

Referenced by G4EmParameters::StreamInfo().

◆ Initialise()

void G4EmExtraParameters::Initialise ( )

Definition at line 70 of file G4EmExtraParameters.cc.

71{
72 quantumEntanglement = false;
73 directionalSplitting = false;
74 directionalSplittingTarget.set(0.,0.,0.);
75 directionalSplittingRadius = 0.;
76
77 dRoverRange = 0.2;
78 finalRange = CLHEP::mm;
79 dRoverRangeMuHad = 0.2;
80 finalRangeMuHad = 0.1*CLHEP::mm;
81 dRoverRangeLIons = 0.2;
82 finalRangeLIons = 0.1*CLHEP::mm;
83 dRoverRangeIons = 0.2;
84 finalRangeIons = 0.1*CLHEP::mm;
85
86 m_regnamesForced.clear();
87 m_procForced.clear();
88 m_lengthForced.clear();
89 m_weightForced.clear();
90 m_regnamesSubCut.clear();
91}
void set(double x, double y, double z)

Referenced by G4EmExtraParameters(), and G4EmParameters::SetDefaults().

◆ operator=()

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

◆ ParticlesPAI()

const std::vector< G4String > & G4EmExtraParameters::ParticlesPAI ( ) const

Definition at line 249 of file G4EmExtraParameters.cc.

250{
251 return m_particlesPAI;
252}

Referenced by G4EmParameters::ParticlesPAI().

◆ QuantumEntanglement()

G4bool G4EmExtraParameters::QuantumEntanglement ( )

Definition at line 446 of file G4EmExtraParameters.cc.

447{
448 return quantumEntanglement;
449}

Referenced by G4EmParameters::QuantumEntanglement(), and G4EmParameters::StreamInfo().

◆ RegionsPAI()

const std::vector< G4String > & G4EmExtraParameters::RegionsPAI ( ) const

Definition at line 254 of file G4EmExtraParameters.cc.

255{
256 return m_regnamesPAI;
257}

Referenced by G4EmParameters::RegionsPAI().

◆ RegionsPhysics()

const std::vector< G4String > & G4EmExtraParameters::RegionsPhysics ( ) const

Definition at line 276 of file G4EmExtraParameters.cc.

277{
278 return m_regnamesPhys;
279}

Referenced by G4EmParameters::RegionsPhysics().

◆ SetDirectionalSplitting()

void G4EmExtraParameters::SetDirectionalSplitting ( G4bool v)

Definition at line 460 of file G4EmExtraParameters.cc.

461{
462 directionalSplitting = v;
463}

Referenced by G4EmParameters::SetDirectionalSplitting(), and G4EmExtraParametersMessenger::SetNewValue().

◆ SetDirectionalSplittingRadius()

void G4EmExtraParameters::SetDirectionalSplittingRadius ( G4double r)

Definition at line 476 of file G4EmExtraParameters.cc.

477{
478 directionalSplittingRadius = r;
479}

Referenced by G4EmParameters::SetDirectionalSplittingRadius(), and G4EmExtraParametersMessenger::SetNewValue().

◆ SetDirectionalSplittingTarget()

void G4EmExtraParameters::SetDirectionalSplittingTarget ( const G4ThreeVector & v)

Definition at line 466 of file G4EmExtraParameters.cc.

467{
468 directionalSplittingTarget = v;
469}

Referenced by G4EmParameters::SetDirectionalSplittingTarget(), and G4EmExtraParametersMessenger::SetNewValue().

◆ SetProcessBiasingFactor()

void G4EmExtraParameters::SetProcessBiasingFactor ( const G4String & procname,
G4double val,
G4bool wflag )

Definition at line 299 of file G4EmExtraParameters.cc.

301{
302 if(val > 0.0) {
303 std::size_t n = m_procBiasedXS.size();
304 for(std::size_t i=0; i<n; ++i) {
305 if(procname == m_procBiasedXS[i]) {
306 m_factBiasedXS[i] = val;
307 m_weightBiasedXS[i]= wflag;
308 return;
309 }
310 }
311 m_procBiasedXS.push_back(procname);
312 m_factBiasedXS.push_back(val);
313 m_weightBiasedXS.push_back(wflag);
314 } else {
316 ed << "Process: " << procname << " XS biasing factor "
317 << val << " is negative - ignored";
318 PrintWarning(ed);
319 }
320}

Referenced by G4EmExtraParametersMessenger::SetNewValue(), and G4EmParameters::SetProcessBiasingFactor().

◆ SetQuantumEntanglement()

void G4EmExtraParameters::SetQuantumEntanglement ( G4bool v)

Definition at line 451 of file G4EmExtraParameters.cc.

452{
453 quantumEntanglement = v;
454}

Referenced by G4EmExtraParametersMessenger::SetNewValue(), and G4EmParameters::SetQuantumEntanglement().

◆ SetStepFunction()

void G4EmExtraParameters::SetStepFunction ( G4double v1,
G4double v2 )

Definition at line 110 of file G4EmExtraParameters.cc.

111{
112 if(v1 > 0.0 && v1 <= 1.0 && v2 > 0.0) {
113 dRoverRange = v1;
114 finalRange = v2;
115 } else {
117 ed << "Values of step function are out of range: "
118 << v1 << ", " << v2/CLHEP::mm << " mm - are ignored";
119 PrintWarning(ed);
120 }
121}

Referenced by G4EmExtraParametersMessenger::SetNewValue(), and G4EmParameters::SetStepFunction().

◆ SetStepFunctionIons()

void G4EmExtraParameters::SetStepFunctionIons ( G4double v1,
G4double v2 )

Definition at line 179 of file G4EmExtraParameters.cc.

180{
181 if(v1 > 0.0 && v1 <= 1.0 && v2 > 0.0) {
182 dRoverRangeIons = v1;
183 finalRangeIons = v2;
184 } else {
186 ed << "Values of step function are out of range: "
187 << v1 << ", " << v2/CLHEP::mm << " mm - are ignored";
188 PrintWarning(ed);
189 }
190}

Referenced by G4EmExtraParametersMessenger::SetNewValue(), and G4EmParameters::SetStepFunctionIons().

◆ SetStepFunctionLightIons()

void G4EmExtraParameters::SetStepFunctionLightIons ( G4double v1,
G4double v2 )

Definition at line 156 of file G4EmExtraParameters.cc.

157{
158 if(v1 > 0.0 && v1 <= 1.0 && v2 > 0.0) {
159 dRoverRangeLIons = v1;
160 finalRangeLIons = v2;
161 } else {
163 ed << "Values of step function are out of range: "
164 << v1 << ", " << v2/CLHEP::mm << " mm - are ignored";
165 PrintWarning(ed);
166 }
167}

Referenced by G4EmExtraParametersMessenger::SetNewValue(), and G4EmParameters::SetStepFunctionLightIons().

◆ SetStepFunctionMuHad()

void G4EmExtraParameters::SetStepFunctionMuHad ( G4double v1,
G4double v2 )

Definition at line 133 of file G4EmExtraParameters.cc.

134{
135 if(v1 > 0.0 && v1 <= 1.0 && v2 > 0.0) {
136 dRoverRangeMuHad = v1;
137 finalRangeMuHad = v2;
138 } else {
140 ed << "Values of step function are out of range: "
141 << v1 << ", " << v2/CLHEP::mm << " mm - are ignored";
142 PrintWarning(ed);
143 }
144}

Referenced by G4EmExtraParametersMessenger::SetNewValue(), and G4EmParameters::SetStepFunctionMuHad().

◆ SetSubCutRegion()

void G4EmExtraParameters::SetSubCutRegion ( const G4String & region)

Definition at line 286 of file G4EmExtraParameters.cc.

287{
288 const G4String& r = CheckRegion(region);
289 std::size_t nreg = m_regnamesSubCut.size();
290 for(std::size_t i=0; i<nreg; ++i) {
291 if(r == m_regnamesSubCut[i]) {
292 return;
293 }
294 }
295 m_regnamesSubCut.push_back(r);
296}

Referenced by G4EmExtraParametersMessenger::SetNewValue(), and G4EmParameters::SetSubCutRegion().

◆ TypesPAI()

const std::vector< G4String > & G4EmExtraParameters::TypesPAI ( ) const

Definition at line 259 of file G4EmExtraParameters.cc.

260{
261 return m_typesPAI;
262}

Referenced by G4EmParameters::TypesPAI().

◆ TypesPhysics()

const std::vector< G4String > & G4EmExtraParameters::TypesPhysics ( ) const

Definition at line 281 of file G4EmExtraParameters.cc.

282{
283 return m_typesPhys;
284}

Referenced by G4EmParameters::TypesPhysics().


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