Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VEnergyLossProcess Class Referenceabstract

#include <G4VEnergyLossProcess.hh>

Inheritance diagram for G4VEnergyLossProcess:
Collaboration diagram for G4VEnergyLossProcess:

Public Member Functions

 G4VEnergyLossProcess (const G4String &name="EnergyLoss", G4ProcessType type=fElectromagnetic)
 
virtual ~G4VEnergyLossProcess ()
 
virtual G4bool IsApplicable (const G4ParticleDefinition &p) override=0
 
virtual void PrintInfo ()=0
 
virtual void ProcessDescription (std::ostream &outFile) const
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &) override
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &) override
 
G4PhysicsTableBuildDEDXTable (G4EmTableType tType=fRestricted)
 
G4PhysicsTableBuildLambdaTable (G4EmTableType tType=fRestricted)
 
void PrintInfoDefinition (const G4ParticleDefinition &part)
 
virtual void StartTracking (G4Track *) override
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection) override
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &) override
 
G4double SampleSubCutSecondaries (std::vector< G4Track * > &, const G4Step &, G4VEmModel *model, G4int matIdx)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &) override
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false) override
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &directory, G4bool ascii) override
 
G4double GetDEDXDispersion (const G4MaterialCutsCouple *couple, const G4DynamicParticle *dp, G4double length)
 
G4double CrossSectionPerVolume (G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double MeanFreePath (const G4Track &track)
 
G4double ContinuousStepLimit (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety)
 
G4VEmModelSelectModelForMaterial (G4double kinEnergy, size_t &idx) const
 
void AddEmModel (G4int, G4VEmModel *, G4VEmFluctuationModel *fluc=0, const G4Region *region=nullptr)
 
void UpdateEmModel (const G4String &, G4double, G4double)
 
void SetEmModel (G4VEmModel *, G4int index=1)
 
G4VEmModelEmModel (G4int index=1) const
 
G4VEmModelGetModelByIndex (G4int idx=0, G4bool ver=false) const
 
G4int NumberOfModels () const
 
void SetFluctModel (G4VEmFluctuationModel *)
 
G4VEmFluctuationModelFluctModel ()
 
void SetBaseParticle (const G4ParticleDefinition *p)
 
const G4ParticleDefinitionParticle () const
 
const G4ParticleDefinitionBaseParticle () const
 
const G4ParticleDefinitionSecondaryParticle () const
 
void ActivateSubCutoff (G4bool val, const G4Region *region=nullptr)
 
void SetCrossSectionBiasingFactor (G4double f, G4bool flag=true)
 
void ActivateForcedInteraction (G4double length, const G4String &region, G4bool flag=true)
 
void ActivateSecondaryBiasing (const G4String &region, G4double factor, G4double energyLimit)
 
void AddCollaborativeProcess (G4VEnergyLossProcess *)
 
void SetLossFluctuations (G4bool val)
 
void SetIntegral (G4bool val)
 
G4bool IsIntegral () const
 
void SetIonisation (G4bool val)
 
G4bool IsIonisationProcess () const
 
void SetLinearLossLimit (G4double val)
 
void SetStepFunction (G4double v1, G4double v2, G4bool lock=true)
 
void SetLowestEnergyLimit (G4double)
 
G4int NumberOfSubCutoffRegions () const
 
void SetDEDXTable (G4PhysicsTable *p, G4EmTableType tType)
 
void SetCSDARangeTable (G4PhysicsTable *pRange)
 
void SetRangeTableForLoss (G4PhysicsTable *p)
 
void SetSecondaryRangeTable (G4PhysicsTable *p)
 
void SetInverseRangeTable (G4PhysicsTable *p)
 
void SetLambdaTable (G4PhysicsTable *p)
 
void SetSubLambdaTable (G4PhysicsTable *p)
 
void SetDEDXBinning (G4int nbins)
 
void SetMinKinEnergy (G4double e)
 
G4double MinKinEnergy () const
 
void SetMaxKinEnergy (G4double e)
 
G4double MaxKinEnergy () const
 
G4double CrossSectionBiasingFactor () const
 
G4double GetDEDX (G4double &kineticEnergy, const G4MaterialCutsCouple *)
 
G4double GetDEDXForSubsec (G4double &kineticEnergy, const G4MaterialCutsCouple *)
 
G4double GetRange (G4double &kineticEnergy, const G4MaterialCutsCouple *)
 
G4double GetCSDARange (G4double &kineticEnergy, const G4MaterialCutsCouple *)
 
G4double GetRangeForLoss (G4double &kineticEnergy, const G4MaterialCutsCouple *)
 
G4double GetKineticEnergy (G4double &range, const G4MaterialCutsCouple *)
 
G4double GetLambda (G4double &kineticEnergy, const G4MaterialCutsCouple *)
 
G4bool TablesAreBuilt () const
 
G4PhysicsTableDEDXTable () const
 
G4PhysicsTableDEDXTableForSubsec () const
 
G4PhysicsTableDEDXunRestrictedTable () const
 
G4PhysicsTableIonisationTable () const
 
G4PhysicsTableIonisationTableForSubsec () const
 
G4PhysicsTableCSDARangeTable () const
 
G4PhysicsTableSecondaryRangeTable () const
 
G4PhysicsTableRangeTableForLoss () const
 
G4PhysicsTableInverseRangeTable () const
 
G4PhysicsTableLambdaTable () const
 
G4PhysicsTableSubLambdaTable () const
 
const G4ElementGetCurrentElement () const
 
void SetDynamicMassCharge (G4double massratio, G4double charge2ratio)
 
- Public Member Functions inherited from G4VContinuousDiscreteProcess
 G4VContinuousDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VContinuousDiscreteProcess (G4VContinuousDiscreteProcess &)
 
virtual ~G4VContinuousDiscreteProcess ()
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
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)
 
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 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
 
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 &)
 

Protected Member Functions

virtual void InitialiseEnergyLossProcess (const G4ParticleDefinition *, const G4ParticleDefinition *)=0
 
virtual G4double MinPrimaryEnergy (const G4ParticleDefinition *, const G4Material *, G4double cut)
 
virtual G4double GetMeanFreePath (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
 
virtual G4double GetContinuousStepLimit (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety) override
 
G4PhysicsVectorLambdaPhysicsVector (const G4MaterialCutsCouple *, G4double cut)
 
size_t CurrentMaterialCutsCoupleIndex () const
 
void SelectModel (G4double kinEnergy)
 
void SetParticle (const G4ParticleDefinition *p)
 
void SetSecondaryParticle (const G4ParticleDefinition *p)
 
- Protected Member Functions inherited from G4VContinuousDiscreteProcess
void SetGPILSelection (G4GPILSelection selection)
 
G4GPILSelection GetGPILSelection () const
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Protected Attributes

G4ParticleChangeForLoss fParticleChange
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 

Detailed Description

Definition at line 124 of file G4VEnergyLossProcess.hh.

Constructor & Destructor Documentation

G4VEnergyLossProcess::G4VEnergyLossProcess ( const G4String name = "EnergyLoss",
G4ProcessType  type = fElectromagnetic 
)

Definition at line 167 of file G4VEnergyLossProcess.cc.

168  :
169  G4VContinuousDiscreteProcess(name, type),
170  secondaryParticle(nullptr),
171  nSCoffRegions(0),
172  idxSCoffRegions(nullptr),
173  nProcesses(0),
174  theDEDXTable(nullptr),
175  theDEDXSubTable(nullptr),
176  theDEDXunRestrictedTable(nullptr),
177  theIonisationTable(nullptr),
178  theIonisationSubTable(nullptr),
179  theRangeTableForLoss(nullptr),
180  theCSDARangeTable(nullptr),
181  theSecondaryRangeTable(nullptr),
182  theInverseRangeTable(nullptr),
183  theLambdaTable(nullptr),
184  theSubLambdaTable(nullptr),
185  theDensityFactor(nullptr),
186  theDensityIdx(nullptr),
187  baseParticle(nullptr),
188  lossFluctuationFlag(true),
189  rndmStepFlag(false),
190  tablesAreBuilt(false),
191  integral(true),
192  isIon(false),
193  isIonisation(true),
194  useSubCutoff(false),
195  useDeexcitation(false),
196  particle(nullptr),
197  currentCouple(nullptr),
198  mfpKinEnergy(0.0)
199 {
200  theParameters = G4EmParameters::Instance();
201  SetVerboseLevel(1);
202 
203  // low energy limit
204  lowestKinEnergy = theParameters->LowestElectronEnergy();
205  preStepKinEnergy = 0.0;
206  preStepRangeEnergy = 0.0;
207  computedRange = DBL_MAX;
208 
209  // Size of tables assuming spline
210  minKinEnergy = 0.1*keV;
211  maxKinEnergy = 100.0*TeV;
212  nBins = 77;
213  maxKinEnergyCSDA = 1.0*GeV;
214  nBinsCSDA = 35;
215  actMinKinEnergy = actMaxKinEnergy = actBinning = actLinLossLimit
216  = actLossFluc = actIntegral = actStepFunc = false;
217 
218  // default linear loss limit for spline
219  linLossLimit = 0.01;
220  dRoverRange = 0.2;
221  finalRange = CLHEP::mm;
222 
223  // default lambda factor
224  lambdaFactor = 0.8;
225 
226  // cross section biasing
227  biasFactor = 1.0;
228 
229  // particle types
230  theElectron = G4Electron::Electron();
231  thePositron = G4Positron::Positron();
232  theGamma = G4Gamma::Gamma();
233  theGenericIon = nullptr;
234 
235  // run time objects
238  modelManager = new G4EmModelManager();
240  ->GetSafetyHelper();
241  aGPILSelection = CandidateForSelection;
242 
243  // initialise model
244  lManager = G4LossTableManager::Instance();
245  lManager->Register(this);
246  fluctModel = nullptr;
247  currentModel = nullptr;
248  atomDeexcitation = nullptr;
249  subcutProducer = nullptr;
250 
251  biasManager = nullptr;
252  biasFlag = false;
253  weightFlag = false;
254  isMaster = true;
255  lastIdx = 0;
256 
257  idxDEDX = idxDEDXSub = idxDEDXunRestricted = idxIonisation =
258  idxIonisationSub = idxRange = idxCSDA = idxSecRange =
259  idxInverseRange = idxLambda = idxSubLambda = 0;
260 
261  scTracks.reserve(5);
262  secParticles.reserve(5);
263 
264  theCuts = theSubCuts = nullptr;
265  currentMaterial = nullptr;
266  currentCoupleIndex = basedCoupleIndex = 0;
267  massRatio = fFactor = reduceFactor = chargeSqRatio = 1.0;
268  preStepLambda = preStepScaledEnergy = fRange = 0.0;
269 
270  secID = biasID = subsecID = -1;
271 }
G4SafetyHelper * GetSafetyHelper() const
static G4LossTableManager * Instance()
G4double LowestElectronEnergy() const
G4ParticleChangeForLoss fParticleChange
static constexpr double TeV
Definition: G4SIunits.hh:218
void SetSecondaryWeightByProcess(G4bool)
static constexpr double mm
Definition: SystemOfUnits.h:95
void Register(G4VEnergyLossProcess *p)
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
static G4TransportationManager * GetTransportationManager()
static G4Positron * Positron()
Definition: G4Positron.cc:94
static G4EmParameters * Instance()
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
static constexpr double GeV
Definition: G4SIunits.hh:217
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define DBL_MAX
Definition: templates.hh:83
static constexpr double keV
Definition: G4SIunits.hh:216
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:437

Here is the call graph for this function:

G4VEnergyLossProcess::~G4VEnergyLossProcess ( )
virtual

Definition at line 275 of file G4VEnergyLossProcess.cc.

276 {
277  /*
278  G4cout << "** G4VEnergyLossProcess::~G4VEnergyLossProcess() for "
279  << GetProcessName() << " isMaster: " << isMaster
280  << " basePart: " << baseParticle
281  << G4endl;
282  */
283  Clean();
284 
285  // G4cout << " isIonisation " << isIonisation << " "
286  // << theDEDXTable << " " << theIonisationTable << G4endl;
287 
288  if (isMaster && !baseParticle) {
289  if(theDEDXTable) {
290 
291  //G4cout << " theIonisationTable " << theIonisationTable << G4endl;
292  if(theIonisationTable == theDEDXTable) { theIonisationTable = 0; }
293  //G4cout << " delete theDEDXTable " << theDEDXTable << G4endl;
294  theDEDXTable->clearAndDestroy();
295  delete theDEDXTable;
296  theDEDXTable = nullptr;
297  if(theDEDXSubTable) {
298  if(theIonisationSubTable == theDEDXSubTable)
299  { theIonisationSubTable = nullptr; }
300  theDEDXSubTable->clearAndDestroy();
301  delete theDEDXSubTable;
302  theDEDXSubTable = nullptr;
303  }
304  }
305  //G4cout << " theIonisationTable " << theIonisationTable << G4endl;
306  if(theIonisationTable) {
307  //G4cout << " delete theIonisationTable " << theIonisationTable << G4endl;
308  theIonisationTable->clearAndDestroy();
309  delete theIonisationTable;
310  theIonisationTable = nullptr;
311  }
312  if(theIonisationSubTable) {
313  theIonisationSubTable->clearAndDestroy();
314  delete theIonisationSubTable;
315  theIonisationSubTable = nullptr;
316  }
317  if(theDEDXunRestrictedTable && isIonisation) {
318  theDEDXunRestrictedTable->clearAndDestroy();
319  delete theDEDXunRestrictedTable;
320  theDEDXunRestrictedTable = nullptr;
321  }
322  if(theCSDARangeTable && isIonisation) {
323  theCSDARangeTable->clearAndDestroy();
324  delete theCSDARangeTable;
325  theCSDARangeTable = nullptr;
326  }
327  //G4cout << "delete RangeTable: " << theRangeTableForLoss << G4endl;
328  if(theRangeTableForLoss && isIonisation) {
329  theRangeTableForLoss->clearAndDestroy();
330  delete theRangeTableForLoss;
331  theRangeTableForLoss = nullptr;
332  }
333  //G4cout << "delete InvRangeTable: " << theInverseRangeTable << G4endl;
334  if(theInverseRangeTable && isIonisation /*&& !isIon*/) {
335  theInverseRangeTable->clearAndDestroy();
336  delete theInverseRangeTable;
337  theInverseRangeTable = nullptr;
338  }
339  //G4cout << "delete LambdaTable: " << theLambdaTable << G4endl;
340  if(theLambdaTable) {
341  theLambdaTable->clearAndDestroy();
342  delete theLambdaTable;
343  theLambdaTable = nullptr;
344  }
345  if(theSubLambdaTable) {
346  theSubLambdaTable->clearAndDestroy();
347  delete theSubLambdaTable;
348  theSubLambdaTable = nullptr;
349  }
350  }
351 
352  delete modelManager;
353  delete biasManager;
354  lManager->DeRegister(this);
355  //G4cout << "** all removed" << G4endl;
356 }
void DeRegister(G4VEnergyLossProcess *p)
void clearAndDestroy()

Here is the call graph for this function:

Member Function Documentation

void G4VEnergyLossProcess::ActivateForcedInteraction ( G4double  length,
const G4String region,
G4bool  flag = true 
)

Definition at line 2288 of file G4VEnergyLossProcess.cc.

2291 {
2292  if(!biasManager) { biasManager = new G4EmBiasingManager(); }
2293  if(1 < verboseLevel) {
2294  G4cout << "### ActivateForcedInteraction: for "
2295  << " process " << GetProcessName()
2296  << " length(mm)= " << length/mm
2297  << " in G4Region <" << region
2298  << "> weightFlag= " << flag
2299  << G4endl;
2300  }
2301  weightFlag = flag;
2302  biasManager->ActivateForcedInteraction(length, region);
2303 }
void ActivateForcedInteraction(G4double length=0.0, const G4String &r="")
static constexpr double mm
Definition: G4SIunits.hh:115
G4int verboseLevel
Definition: G4VProcess.hh:368
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::ActivateSecondaryBiasing ( const G4String region,
G4double  factor,
G4double  energyLimit 
)

Definition at line 2308 of file G4VEnergyLossProcess.cc.

2311 {
2312  if (0.0 <= factor) {
2313 
2314  // Range cut can be applied only for e-
2315  if(0.0 == factor && secondaryParticle != G4Electron::Electron())
2316  { return; }
2317 
2318  if(!biasManager) { biasManager = new G4EmBiasingManager(); }
2319  biasManager->ActivateSecondaryBiasing(region, factor, energyLimit);
2320  if(1 < verboseLevel) {
2321  G4cout << "### ActivateSecondaryBiasing: for "
2322  << " process " << GetProcessName()
2323  << " factor= " << factor
2324  << " in G4Region <" << region
2325  << "> energyLimit(MeV)= " << energyLimit/MeV
2326  << G4endl;
2327  }
2328  }
2329 }
G4int verboseLevel
Definition: G4VProcess.hh:368
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::ActivateSubCutoff ( G4bool  val,
const G4Region region = nullptr 
)

Definition at line 1053 of file G4VEnergyLossProcess.cc.

1054 {
1055  G4RegionStore* regionStore = G4RegionStore::GetInstance();
1056  const G4Region* reg = r;
1057  if (!reg) {
1058  reg = regionStore->GetRegion("DefaultRegionForTheWorld", false);
1059  }
1060 
1061  // the region is in the list
1062  if (nSCoffRegions > 0) {
1063  for (G4int i=0; i<nSCoffRegions; ++i) {
1064  if (reg == scoffRegions[i]) {
1065  return;
1066  }
1067  }
1068  }
1069  // new region
1070  if(val) {
1071  scoffRegions.push_back(reg);
1072  ++nSCoffRegions;
1073  }
1074 }
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
int G4int
Definition: G4Types.hh:78
static const G4double reg
static G4RegionStore * GetInstance()

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::AddCollaborativeProcess ( G4VEnergyLossProcess p)

Definition at line 2037 of file G4VEnergyLossProcess.cc.

2039 {
2040  G4bool add = true;
2041  if(p->GetProcessName() != "eBrem") { add = false; }
2042  if(add && nProcesses > 0) {
2043  for(G4int i=0; i<nProcesses; ++i) {
2044  if(p == scProcesses[i]) {
2045  add = false;
2046  break;
2047  }
2048  }
2049  }
2050  if(add) {
2051  scProcesses.push_back(p);
2052  ++nProcesses;
2053  if (1 < verboseLevel) {
2054  G4cout << "### The process " << p->GetProcessName()
2055  << " is added to the list of collaborative processes of "
2056  << GetProcessName() << G4endl;
2057  }
2058  }
2059 }
G4int verboseLevel
Definition: G4VProcess.hh:368
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

void G4VEnergyLossProcess::AddEmModel ( G4int  order,
G4VEmModel p,
G4VEmFluctuationModel fluc = 0,
const G4Region region = nullptr 
)

Definition at line 391 of file G4VEnergyLossProcess.cc.

394 {
395  modelManager->AddEmModel(order, p, fluc, region);
396  if(p) { p->SetParticleChange(pParticleChange, fluc); }
397 }
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
Definition: G4VEmModel.cc:418
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283

Here is the call graph for this function:

Here is the caller graph for this function:

G4VParticleChange * G4VEnergyLossProcess::AlongStepDoIt ( const G4Track track,
const G4Step step 
)
overridevirtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 1250 of file G4VEnergyLossProcess.cc.

1252 {
1254  // The process has range table - calculate energy loss
1255  if(!isIonisation || !currentModel->IsActive(preStepScaledEnergy)) {
1256  return &fParticleChange;
1257  }
1258 
1259  // Get the actual (true) Step length
1260  G4double length = step.GetStepLength();
1261  if(length <= 0.0) { return &fParticleChange; }
1262  G4double eloss = 0.0;
1263 
1264  /*
1265  if(-1 < verboseLevel) {
1266  const G4ParticleDefinition* d = track.GetParticleDefinition();
1267  G4cout << "AlongStepDoIt for "
1268  << GetProcessName() << " and particle "
1269  << d->GetParticleName()
1270  << " eScaled(MeV)= " << preStepScaledEnergy/MeV
1271  << " range(mm)= " << fRange/mm
1272  << " s(mm)= " << length/mm
1273  << " rf= " << reduceFactor
1274  << " q^2= " << chargeSqRatio
1275  << " md= " << d->GetPDGMass()
1276  << " status= " << track.GetTrackStatus()
1277  << " " << track.GetMaterial()->GetName()
1278  << G4endl;
1279  }
1280  */
1281 
1282  const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
1283 
1284  // define new weight for primary and secondaries
1286  if(weightFlag) {
1287  weight /= biasFactor;
1289  }
1290 
1291  // stopping
1292  if (length >= fRange || preStepKinEnergy <= lowestKinEnergy) {
1293  eloss = preStepKinEnergy;
1294  if (useDeexcitation) {
1295  atomDeexcitation->AlongStepDeexcitation(scTracks, step,
1296  eloss, currentCoupleIndex);
1297  if(scTracks.size() > 0) { FillSecondariesAlongStep(eloss, weight); }
1298  if(eloss < 0.0) { eloss = 0.0; }
1299  }
1302  return &fParticleChange;
1303  }
1304  //G4cout << theDEDXTable << " idx= " << basedCoupleIndex
1305  // << " " << GetProcessName() << " "<< currentMaterial->GetName()<<G4endl;
1306  //if(particle->GetParticleName() == "e-")G4cout << (*theDEDXTable) <<G4endl;
1307  // Short step
1308  eloss = GetDEDXForScaledEnergy(preStepScaledEnergy)*length;
1309 
1310  //G4cout << "eloss= " << eloss << G4endl;
1311 
1312  // Long step
1313  if(eloss > preStepKinEnergy*linLossLimit) {
1314 
1315  G4double x = (fRange - length)/reduceFactor;
1316  //G4cout << "x= " << x << " " << theInverseRangeTable << G4endl;
1317  eloss = preStepKinEnergy - ScaledKinEnergyForLoss(x)/massRatio;
1318 
1319  /*
1320  if(-1 < verboseLevel)
1321  G4cout << "Long STEP: rPre(mm)= "
1322  << GetScaledRangeForScaledEnergy(preStepScaledEnergy)/mm
1323  << " rPost(mm)= " << x/mm
1324  << " ePre(MeV)= " << preStepScaledEnergy/MeV
1325  << " eloss(MeV)= " << eloss/MeV
1326  << " eloss0(MeV)= "
1327  << GetDEDXForScaledEnergy(preStepScaledEnergy)*length/MeV
1328  << " lim(MeV)= " << preStepKinEnergy*linLossLimit/MeV
1329  << G4endl;
1330  */
1331  }
1332 
1333  /*
1334  G4double eloss0 = eloss;
1335  if(-1 < verboseLevel ) {
1336  G4cout << "Before fluct: eloss(MeV)= " << eloss/MeV
1337  << " e-eloss= " << preStepKinEnergy-eloss
1338  << " step(mm)= " << length/mm
1339  << " range(mm)= " << fRange/mm
1340  << " fluct= " << lossFluctuationFlag
1341  << G4endl;
1342  }
1343  */
1344 
1345  G4double cut = (*theCuts)[currentCoupleIndex];
1346  G4double esec = 0.0;
1347 
1348  //G4cout << "cut= " << cut << " useSubCut= " << useSubCutoff << G4endl;
1349 
1350  // SubCutOff
1351  if(useSubCutoff && !subcutProducer) {
1352  if(idxSCoffRegions[currentCoupleIndex]) {
1353 
1354  G4bool yes = false;
1355  G4StepPoint* prePoint = step.GetPreStepPoint();
1356 
1357  // Check boundary
1358  if(prePoint->GetStepStatus() == fGeomBoundary) { yes = true; }
1359 
1360  // Check PrePoint
1361  else {
1362  G4double preSafety = prePoint->GetSafety();
1363  G4double rcut =
1364  currentCouple->GetProductionCuts()->GetProductionCut(1);
1365 
1366  // recompute presafety
1367  if(preSafety < rcut) {
1368  preSafety = safetyHelper->ComputeSafety(prePoint->GetPosition(),
1369  rcut);
1370  }
1371 
1372  if(preSafety < rcut) { yes = true; }
1373 
1374  // Check PostPoint
1375  else {
1376  G4double postSafety = preSafety - length;
1377  if(postSafety < rcut) {
1378  postSafety = safetyHelper->ComputeSafety(
1379  step.GetPostStepPoint()->GetPosition(), rcut);
1380  if(postSafety < rcut) { yes = true; }
1381  }
1382  }
1383  }
1384 
1385  // Decided to start subcut sampling
1386  if(yes) {
1387 
1388  cut = (*theSubCuts)[currentCoupleIndex];
1389  eloss -= GetSubDEDXForScaledEnergy(preStepScaledEnergy)*length;
1390  esec = SampleSubCutSecondaries(scTracks, step,
1391  currentModel,currentCoupleIndex);
1392  // add bremsstrahlung sampling
1393  /*
1394  if(nProcesses > 0) {
1395  for(G4int i=0; i<nProcesses; ++i) {
1396  (scProcesses[i])->SampleSubCutSecondaries(
1397  scTracks, step, (scProcesses[i])->
1398  SelectModelForMaterial(preStepKinEnergy, currentCoupleIndex),
1399  currentCoupleIndex);
1400  }
1401  }
1402  */
1403  }
1404  }
1405  }
1406 
1407  // Corrections, which cannot be tabulated
1408  if(isIon) {
1409  G4double eadd = 0.0;
1410  G4double eloss_before = eloss;
1411  currentModel->CorrectionsAlongStep(currentCouple, dynParticle,
1412  eloss, eadd, length);
1413  if(eloss < 0.0) { eloss = 0.5*eloss_before; }
1414  }
1415 
1416  // Sample fluctuations
1417  if (lossFluctuationFlag) {
1418  G4VEmFluctuationModel* fluc = currentModel->GetModelOfFluctuations();
1419  if(eloss + esec < preStepKinEnergy) {
1420 
1421  G4double tmax =
1422  std::min(currentModel->MaxSecondaryKinEnergy(dynParticle),cut);
1423  eloss = fluc->SampleFluctuations(currentCouple,dynParticle,
1424  tmax,length,eloss);
1425  /*
1426  if(-1 < verboseLevel)
1427  G4cout << "After fluct: eloss(MeV)= " << eloss/MeV
1428  << " fluc= " << (eloss-eloss0)/MeV
1429  << " ChargeSqRatio= " << chargeSqRatio
1430  << " massRatio= " << massRatio
1431  << " tmax= " << tmax
1432  << G4endl;
1433  */
1434  }
1435  }
1436 
1437  // deexcitation
1438  if (useDeexcitation) {
1439  G4double esecfluo = preStepKinEnergy - esec;
1440  G4double de = esecfluo;
1441  //G4double eloss0 = eloss;
1442  /*
1443  G4cout << "### 1: E(keV)= " << preStepKinEnergy/keV
1444  << " Efluomax(keV)= " << de/keV
1445  << " Eloss(keV)= " << eloss/keV << G4endl;
1446  */
1447  atomDeexcitation->AlongStepDeexcitation(scTracks, step,
1448  de, currentCoupleIndex);
1449 
1450  // sum of de-excitation energies
1451  esecfluo -= de;
1452 
1453  // subtracted from energy loss
1454  if(eloss >= esecfluo) {
1455  esec += esecfluo;
1456  eloss -= esecfluo;
1457  } else {
1458  esec += esecfluo;
1459  eloss = 0.0;
1460  }
1461  /*
1462  if(esecfluo > 0.0) {
1463  G4cout << "### 2: E(keV)= " << preStepKinEnergy/keV
1464  << " Esec(keV)= " << esec/keV
1465  << " Esecf(kV)= " << esecfluo/keV
1466  << " Eloss0(kV)= " << eloss0/keV
1467  << " Eloss(keV)= " << eloss/keV
1468  << G4endl;
1469  }
1470  */
1471  }
1472  if(subcutProducer && idxSCoffRegions[currentCoupleIndex]) {
1473  subcutProducer->SampleSecondaries(step, scTracks, eloss, cut);
1474  }
1475  if(scTracks.size() > 0) { FillSecondariesAlongStep(eloss, weight); }
1476 
1477  // Energy balance
1478  G4double finalT = preStepKinEnergy - eloss - esec;
1479  if (finalT <= lowestKinEnergy) {
1480  eloss += finalT;
1481  finalT = 0.0;
1482  } else if(isIon) {
1484  currentModel->GetParticleCharge(track.GetParticleDefinition(),
1485  currentMaterial,finalT));
1486  }
1487 
1488  eloss = std::max(eloss, 0.0);
1489 
1492  /*
1493  if(-1 < verboseLevel) {
1494  G4double del = finalT + eloss + esec - preStepKinEnergy;
1495  G4cout << "Final value eloss(MeV)= " << eloss/MeV
1496  << " preStepKinEnergy= " << preStepKinEnergy
1497  << " postStepKinEnergy= " << finalT
1498  << " de(keV)= " << del/keV
1499  << " lossFlag= " << lossFluctuationFlag
1500  << " status= " << track.GetTrackStatus()
1501  << G4endl;
1502  }
1503  */
1504  return &fParticleChange;
1505 }
G4double MaxSecondaryKinEnergy(const G4DynamicParticle *dynParticle)
Definition: G4VEmModel.hh:484
virtual G4double SampleFluctuations(const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmax, G4double length, G4double meanLoss)=0
virtual void SampleSecondaries(const G4Step &step, std::vector< G4Track * > &tracks, G4double &eloss, G4double cut) const =0
void AlongStepDeexcitation(std::vector< G4Track * > &tracks, const G4Step &step, G4double &eLoss, G4int coupleIndex)
virtual void CorrectionsAlongStep(const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
Definition: G4VEmModel.cc:370
G4double GetStepLength() const
const G4DynamicParticle * GetDynamicParticle() const
G4double GetProductionCut(G4int index) const
G4StepStatus GetStepStatus() const
G4ParticleChangeForLoss fParticleChange
G4double GetParentWeight() const
tuple x
Definition: test.py:50
G4VEmFluctuationModel * GetModelOfFluctuations()
Definition: G4VEmModel.hh:612
void InitializeForAlongStep(const G4Track &)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
G4StepPoint * GetPreStepPoint() const
void ProposeWeight(G4double finalWeight)
G4double SampleSubCutSecondaries(std::vector< G4Track * > &, const G4Step &, G4VEmModel *model, G4int matIdx)
const G4ThreeVector & GetPosition() const
bool G4bool
Definition: G4Types.hh:79
void SetProposedKineticEnergy(G4double proposedKinEnergy)
const G4ParticleDefinition * GetParticleDefinition() const
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:760
G4double ComputeSafety(const G4ThreeVector &pGlobalPoint, G4double maxRadius=DBL_MAX)
virtual G4double GetParticleCharge(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:362
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4StepPoint * GetPostStepPoint() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double GetSafety() const
void SetProposedCharge(G4double theCharge)
double G4double
Definition: G4Types.hh:76
G4ProductionCuts * GetProductionCuts() const

Here is the call graph for this function:

G4double G4VEnergyLossProcess::AlongStepGetPhysicalInteractionLength ( const G4Track ,
G4double  previousStepSize,
G4double  currentMinimumStep,
G4double currentSafety,
G4GPILSelection selection 
)
overridevirtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 1117 of file G4VEnergyLossProcess.cc.

1120 {
1121  G4double x = DBL_MAX;
1122  *selection = aGPILSelection;
1123  if(isIonisation && currentModel->IsActive(preStepScaledEnergy)) {
1124  fRange = GetScaledRangeForScaledEnergy(preStepScaledEnergy)*reduceFactor;
1125  x = fRange;
1126  G4double finR = finalRange;
1127  if(rndmStepFlag) {
1128  finR = std::min(finR,
1129  currentCouple->GetProductionCuts()->GetProductionCut(1));
1130  }
1131  if(fRange > finR) {
1132  x = fRange*dRoverRange + finR*(1.0 - dRoverRange)*(2.0 - finR/fRange);
1133  }
1134 
1135  // if(particle->GetPDGMass() > 0.9*GeV)
1136  /*
1137  G4cout<<GetProcessName()<<": e= "<<preStepKinEnergy
1138  <<" range= "<<fRange << " idx= " << basedCoupleIndex
1139  << " finR= " << finR
1140  << " limit= " << x <<G4endl;
1141  G4cout << "massRatio= " << massRatio << " Q^2= " << chargeSqRatio
1142  << " finR= " << finR << " dRoverRange= " << dRoverRange
1143  << " finalRange= " << finalRange << G4endl;
1144  */
1145  }
1146  //G4cout<<GetProcessName()<<": e= "<<preStepKinEnergy
1147  //<<" stepLimit= "<<x<<G4endl;
1148  return x;
1149 }
G4double GetProductionCut(G4int index) const
tuple x
Definition: test.py:50
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:760
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double G4double
Definition: G4Types.hh:76
G4ProductionCuts * GetProductionCuts() const
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

Here is the caller graph for this function:

const G4ParticleDefinition * G4VEnergyLossProcess::BaseParticle ( ) const
inline

Definition at line 927 of file G4VEnergyLossProcess.hh.

928 {
929  return baseParticle;
930 }

Here is the caller graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::BuildDEDXTable ( G4EmTableType  tType = fRestricted)

Definition at line 806 of file G4VEnergyLossProcess.cc.

807 {
808  if(1 < verboseLevel ) {
809  G4cout << "G4VEnergyLossProcess::BuildDEDXTable() of type " << tType
810  << " for " << GetProcessName()
811  << " and particle " << particle->GetParticleName()
812  << G4endl;
813  }
814  G4PhysicsTable* table = nullptr;
815  G4double emax = maxKinEnergy;
816  G4int bin = nBins;
817 
818  if(fTotal == tType) {
819  emax = maxKinEnergyCSDA;
820  bin = nBinsCSDA;
821  table = theDEDXunRestrictedTable;
822  } else if(fRestricted == tType) {
823  table = theDEDXTable;
824  } else if(fSubRestricted == tType) {
825  table = theDEDXSubTable;
826  } else {
827  G4cout << "G4VEnergyLossProcess::BuildDEDXTable WARNING: wrong type "
828  << tType << G4endl;
829  }
830 
831  // Access to materials
832  const G4ProductionCutsTable* theCoupleTable=
834  size_t numOfCouples = theCoupleTable->GetTableSize();
835 
836  if(1 < verboseLevel) {
837  G4cout << numOfCouples << " materials"
838  << " minKinEnergy= " << minKinEnergy
839  << " maxKinEnergy= " << emax
840  << " nbin= " << bin
841  << " EmTableType= " << tType
842  << " table= " << table << " " << this
843  << G4endl;
844  }
845  if(!table) { return table; }
846 
847  G4LossTableBuilder* bld = lManager->GetTableBuilder();
848  G4bool splineFlag = theParameters->Spline();
849  G4PhysicsLogVector* aVector = nullptr;
850  G4PhysicsLogVector* bVector = nullptr;
851 
852  for(size_t i=0; i<numOfCouples; ++i) {
853 
854  if(1 < verboseLevel) {
855  G4cout << "G4VEnergyLossProcess::BuildDEDXVector Idx= " << i
856  << " flagTable= " << table->GetFlag(i)
857  << " Flag= " << bld->GetFlag(i) << G4endl;
858  }
859  if(bld->GetFlag(i)) {
860 
861  // create physics vector and fill it
862  const G4MaterialCutsCouple* couple =
863  theCoupleTable->GetMaterialCutsCouple(i);
864  if((*table)[i]) { delete (*table)[i]; }
865  if(bVector) {
866  aVector = new G4PhysicsLogVector(*bVector);
867  } else {
868  bVector = new G4PhysicsLogVector(minKinEnergy, emax, bin);
869  aVector = bVector;
870  }
871  aVector->SetSpline(splineFlag);
872 
873  modelManager->FillDEDXVector(aVector, couple, tType);
874  if(splineFlag) { aVector->FillSecondDerivatives(); }
875 
876  // Insert vector for this material into the table
877  G4PhysicsTableHelper::SetPhysicsVector(table, i, aVector);
878  }
879  }
880 
881  if(1 < verboseLevel) {
882  G4cout << "G4VEnergyLossProcess::BuildDEDXTable(): table is built for "
883  << particle->GetParticleName()
884  << " and process " << GetProcessName()
885  << G4endl;
886  if(2 < verboseLevel) G4cout << (*table) << G4endl;
887  }
888 
889  return table;
890 }
static void SetPhysicsVector(G4PhysicsTable *physTable, size_t idx, G4PhysicsVector *vec)
G4bool Spline() const
G4int verboseLevel
Definition: G4VProcess.hh:368
tuple bin
Definition: plottest35.py:22
G4bool GetFlag(size_t idx) const
int G4int
Definition: G4Types.hh:78
void FillSecondDerivatives()
const G4String & GetParticleName() const
G4LossTableBuilder * GetTableBuilder()
void SetSpline(G4bool)
G4GLOB_DLL std::ostream G4cout
void FillDEDXVector(G4PhysicsVector *, const G4MaterialCutsCouple *, G4EmTableType t=fRestricted)
bool G4bool
Definition: G4Types.hh:79
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
static const G4double emax
static G4ProductionCutsTable * GetProductionCutsTable()
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
#define G4endl
Definition: G4ios.hh:61
G4bool GetFlag(size_t i) const
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::BuildLambdaTable ( G4EmTableType  tType = fRestricted)

Definition at line 894 of file G4VEnergyLossProcess.cc.

895 {
896  G4PhysicsTable* table = nullptr;
897 
898  if(fRestricted == tType) {
899  table = theLambdaTable;
900  } else if(fSubRestricted == tType) {
901  table = theSubLambdaTable;
902  } else {
903  G4cout << "G4VEnergyLossProcess::BuildLambdaTable WARNING: wrong type "
904  << tType << G4endl;
905  }
906 
907  if(1 < verboseLevel) {
908  G4cout << "G4VEnergyLossProcess::BuildLambdaTable() of type "
909  << tType << " for process "
910  << GetProcessName() << " and particle "
911  << particle->GetParticleName()
912  << " EmTableType= " << tType
913  << " table= " << table
914  << G4endl;
915  }
916  if(!table) {return table;}
917 
918  // Access to materials
919  const G4ProductionCutsTable* theCoupleTable=
921  size_t numOfCouples = theCoupleTable->GetTableSize();
922 
923  G4LossTableBuilder* bld = lManager->GetTableBuilder();
924  theDensityFactor = bld->GetDensityFactors();
925  theDensityIdx = bld->GetCoupleIndexes();
926 
927  G4bool splineFlag = theParameters->Spline();
928  G4PhysicsLogVector* aVector = nullptr;
929  G4double scale = G4Log(maxKinEnergy/minKinEnergy);
930 
931  for(size_t i=0; i<numOfCouples; ++i) {
932 
933  if (bld->GetFlag(i)) {
934 
935  // create physics vector and fill it
936  const G4MaterialCutsCouple* couple =
937  theCoupleTable->GetMaterialCutsCouple(i);
938  delete (*table)[i];
939 
940  G4bool startNull = true;
941  G4double emin =
942  MinPrimaryEnergy(particle,couple->GetMaterial(),(*theCuts)[i]);
943  if(minKinEnergy > emin) {
944  emin = minKinEnergy;
945  startNull = false;
946  }
947 
948  G4double emax = maxKinEnergy;
949  if(emax <= emin) { emax = 2*emin; }
950  G4int bin = G4lrint(nBins*G4Log(emax/emin)/scale);
951  bin = std::max(bin, 3);
952  aVector = new G4PhysicsLogVector(emin, emax, bin);
953  aVector->SetSpline(splineFlag);
954 
955  modelManager->FillLambdaVector(aVector, couple, startNull, tType);
956  if(splineFlag) { aVector->FillSecondDerivatives(); }
957 
958  // Insert vector for this material into the table
959  G4PhysicsTableHelper::SetPhysicsVector(table, i, aVector);
960  }
961  }
962 
963  if(1 < verboseLevel) {
964  G4cout << "Lambda table is built for "
965  << particle->GetParticleName()
966  << G4endl;
967  }
968 
969  return table;
970 }
const std::vector< G4double > * GetDensityFactors()
static void SetPhysicsVector(G4PhysicsTable *physTable, size_t idx, G4PhysicsVector *vec)
G4bool Spline() const
G4int verboseLevel
Definition: G4VProcess.hh:368
tuple bin
Definition: plottest35.py:22
G4bool GetFlag(size_t idx) const
void FillLambdaVector(G4PhysicsVector *, const G4MaterialCutsCouple *, G4bool startFromNull=true, G4EmTableType t=fRestricted)
int G4int
Definition: G4Types.hh:78
void FillSecondDerivatives()
const G4String & GetParticleName() const
G4LossTableBuilder * GetTableBuilder()
void SetSpline(G4bool)
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
static const G4double emax
G4double G4Log(G4double x)
Definition: G4Log.hh:230
static G4ProductionCutsTable * GetProductionCutsTable()
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
int G4lrint(double ad)
Definition: templates.hh:163
T max(const T t1, const T t2)
brief Return the largest of the two arguments
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
virtual G4double MinPrimaryEnergy(const G4ParticleDefinition *, const G4Material *, G4double cut)
const G4Material * GetMaterial() const
const std::vector< G4int > * GetCoupleIndexes()

Here is the call graph for this function:

void G4VEnergyLossProcess::BuildPhysicsTable ( const G4ParticleDefinition part)
overridevirtual

Reimplemented from G4VProcess.

Reimplemented in G4ePolarizedIonisation.

Definition at line 697 of file G4VEnergyLossProcess.cc.

698 {
699  if(1 < verboseLevel) {
700 
701  G4cout << "### G4VEnergyLossProcess::BuildPhysicsTable() for "
702  << GetProcessName()
703  << " and particle " << part.GetParticleName()
704  << "; local: " << particle->GetParticleName();
705  if(baseParticle) {
706  G4cout << "; base: " << baseParticle->GetParticleName();
707  }
708  G4cout << " TablesAreBuilt= " << tablesAreBuilt
709  << " isIon= " << isIon << " " << this << G4endl;
710  }
711 
712  if(&part == particle) {
713 
714  G4LossTableBuilder* bld = lManager->GetTableBuilder();
715  if(isMaster) {
716  theDensityFactor = bld->GetDensityFactors();
717  theDensityIdx = bld->GetCoupleIndexes();
718  lManager->BuildPhysicsTable(particle, this);
719 
720  } else {
721 
722  const G4VEnergyLossProcess* masterProcess =
723  static_cast<const G4VEnergyLossProcess*>(GetMasterProcess());
724 
725  // define density factors for worker thread
726  bld->InitialiseBaseMaterials(masterProcess->DEDXTable());
727  theDensityFactor = bld->GetDensityFactors();
728  theDensityIdx = bld->GetCoupleIndexes();
729 
730  // copy table pointers from master thread
731  SetDEDXTable(masterProcess->DEDXTable(),fRestricted);
733  SetDEDXTable(masterProcess->DEDXunRestrictedTable(),fTotal);
734  SetDEDXTable(masterProcess->IonisationTable(),fIsIonisation);
736  SetRangeTableForLoss(masterProcess->RangeTableForLoss());
737  SetCSDARangeTable(masterProcess->CSDARangeTable());
739  SetInverseRangeTable(masterProcess->InverseRangeTable());
740  SetLambdaTable(masterProcess->LambdaTable());
741  SetSubLambdaTable(masterProcess->SubLambdaTable());
742  isIonisation = masterProcess->IsIonisationProcess();
743 
744  tablesAreBuilt = true;
745  // local initialisation of models
746  G4bool printing = true;
747  G4int numberOfModels = modelManager->NumberOfModels();
748  for(G4int i=0; i<numberOfModels; ++i) {
749  G4VEmModel* mod = GetModelByIndex(i, printing);
750  G4VEmModel* mod0= masterProcess->GetModelByIndex(i,printing);
751  mod->InitialiseLocal(particle, mod0);
752  }
753 
754  lManager->LocalPhysicsTables(particle, this);
755  }
756 
757  // needs to be done only once
758  safetyHelper->InitialiseHelper();
759  }
760  // explicitly defined printout by particle name
761  G4String num = part.GetParticleName();
762  if(1 < verboseLevel ||
763  (0 < verboseLevel && (num == "e-" ||
764  num == "e+" || num == "mu+" ||
765  num == "mu-" || num == "proton"||
766  num == "pi+" || num == "pi-" ||
767  num == "kaon+" || num == "kaon-" ||
768  num == "alpha" || num == "anti_proton" ||
769  num == "GenericIon")))
770  {
771  PrintInfoDefinition(part);
772  }
773 
774  // Added tracking cut to avoid tracking artifacts
775  // identify deexcitation flag
776  if(isIonisation) {
777  //VI: seems not needed fParticleChange.SetLowEnergyLimit(lowestKinEnergy);
778  atomDeexcitation = lManager->AtomDeexcitation();
779  if(nSCoffRegions > 0) { subcutProducer = lManager->SubCutProducer(); }
780  if(atomDeexcitation) {
781  if(atomDeexcitation->IsPIXEActive()) { useDeexcitation = true; }
782  }
783  }
784  /*
785  G4cout << "** G4VEnergyLossProcess::BuildPhysicsTable() for "
786  << GetProcessName() << " and " << particle->GetParticleName()
787  << " isMaster: " << isMaster << " isIonisation: "
788  << isIonisation << G4endl;
789  G4cout << " theDEDX: " << theDEDXTable
790  << " theRange: " << theRangeTableForLoss
791  << " theInverse: " << theInverseRangeTable
792  << " theLambda: " << theLambdaTable << G4endl;
793  */
794  //if(1 < verboseLevel || verb) {
795  if(1 < verboseLevel) {
796  G4cout << "### G4VEnergyLossProcess::BuildPhysicsTable() done for "
797  << GetProcessName()
798  << " and particle " << part.GetParticleName();
799  if(isIonisation) { G4cout << " isIonisation flag = 1"; }
800  G4cout << G4endl;
801  }
802 }
G4VSubCutProducer * SubCutProducer()
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
void SetDEDXTable(G4PhysicsTable *p, G4EmTableType tType)
const std::vector< G4double > * GetDensityFactors()
void InitialiseHelper()
G4int verboseLevel
Definition: G4VProcess.hh:368
void PrintInfoDefinition(const G4ParticleDefinition &part)
G4bool IsPIXEActive() const
G4PhysicsTable * SubLambdaTable() const
G4PhysicsTable * RangeTableForLoss() const
G4PhysicsTable * SecondaryRangeTable() const
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel)
Definition: G4VEmModel.cc:218
G4PhysicsTable * CSDARangeTable() const
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const
G4PhysicsTable * IonisationTableForSubsec() const
G4PhysicsTable * IonisationTable() const
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4LossTableBuilder * GetTableBuilder()
G4PhysicsTable * LambdaTable() const
void SetInverseRangeTable(G4PhysicsTable *p)
G4GLOB_DLL std::ostream G4cout
G4PhysicsTable * DEDXTable() const
bool G4bool
Definition: G4Types.hh:79
G4PhysicsTable * DEDXTableForSubsec() const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void SetLambdaTable(G4PhysicsTable *p)
G4int NumberOfModels() const
G4PhysicsTable * InverseRangeTable() const
void BuildPhysicsTable(const G4ParticleDefinition *aParticle)
void InitialiseBaseMaterials(G4PhysicsTable *table)
void SetCSDARangeTable(G4PhysicsTable *pRange)
void LocalPhysicsTables(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
#define G4endl
Definition: G4ios.hh:61
void SetSecondaryRangeTable(G4PhysicsTable *p)
G4VAtomDeexcitation * AtomDeexcitation()
G4PhysicsTable * DEDXunRestrictedTable() const
void SetSubLambdaTable(G4PhysicsTable *p)
void SetRangeTableForLoss(G4PhysicsTable *p)
const std::vector< G4int > * GetCoupleIndexes()
G4bool IsIonisationProcess() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::ContinuousStepLimit ( const G4Track track,
G4double  previousStepSize,
G4double  currentMinimumStep,
G4double currentSafety 
)

Definition at line 1994 of file G4VEnergyLossProcess.cc.

1997 {
1998  G4GPILSelection sel;
1999  return AlongStepGetPhysicalInteractionLength(track, x, y, z, &sel);
2000 }
tuple x
Definition: test.py:50
virtual G4double AlongStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection) override
tuple z
Definition: test.py:28
G4GPILSelection

Here is the call graph for this function:

G4double G4VEnergyLossProcess::CrossSectionBiasingFactor ( ) const
inline

Definition at line 993 of file G4VEnergyLossProcess.hh.

994 {
995  return biasFactor;
996 }
G4double G4VEnergyLossProcess::CrossSectionPerVolume ( G4double  kineticEnergy,
const G4MaterialCutsCouple couple 
)

Definition at line 1962 of file G4VEnergyLossProcess.cc.

1964 {
1965  // Cross section per volume is calculated
1966  DefineMaterial(couple);
1967  G4double cross = 0.0;
1968  if(theLambdaTable) {
1969  cross = GetLambdaForScaledEnergy(kineticEnergy*massRatio);
1970  } else {
1971  SelectModel(kineticEnergy*massRatio);
1972  cross = biasFactor*(*theDensityFactor)[currentCoupleIndex]
1973  *(currentModel->CrossSectionPerVolume(currentMaterial,
1974  particle, kineticEnergy,
1975  (*theCuts)[currentCoupleIndex]));
1976  }
1977  if(cross < 0.0) { cross = 0.0; }
1978  return cross;
1979 }
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:257
void SelectModel(G4double kinEnergy)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::CSDARangeTable ( ) const
inline

Definition at line 1042 of file G4VEnergyLossProcess.hh.

1043 {
1044  return theCSDARangeTable;
1045 }

Here is the caller graph for this function:

size_t G4VEnergyLossProcess::CurrentMaterialCutsCoupleIndex ( ) const
inlineprotected

Definition at line 612 of file G4VEnergyLossProcess.hh.

613 {
614  return currentCoupleIndex;
615 }
G4PhysicsTable * G4VEnergyLossProcess::DEDXTable ( ) const
inline

Definition at line 1007 of file G4VEnergyLossProcess.hh.

1008 {
1009  return theDEDXTable;
1010 }

Here is the caller graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::DEDXTableForSubsec ( ) const
inline

Definition at line 1014 of file G4VEnergyLossProcess.hh.

1015 {
1016  return theDEDXSubTable;
1017 }

Here is the caller graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::DEDXunRestrictedTable ( ) const
inline

Definition at line 1021 of file G4VEnergyLossProcess.hh.

1022 {
1023  return theDEDXunRestrictedTable;
1024 }

Here is the caller graph for this function:

G4VEmModel * G4VEnergyLossProcess::EmModel ( G4int  index = 1) const

Definition at line 418 of file G4VEnergyLossProcess.cc.

419 {
420  G4VEmModel* p = nullptr;
421  if(index >= 0 && index < G4int(emModels.size())) { p = emModels[index]; }
422  return p;
423 }
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78

Here is the caller graph for this function:

G4VEmFluctuationModel * G4VEnergyLossProcess::FluctModel ( )
inline

Definition at line 890 of file G4VEnergyLossProcess.hh.

891 {
892  return fluctModel;
893 }

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::GetContinuousStepLimit ( const G4Track track,
G4double  previousStepSize,
G4double  currentMinimumStep,
G4double currentSafety 
)
overrideprotectedvirtual

Implements G4VContinuousDiscreteProcess.

Definition at line 2016 of file G4VEnergyLossProcess.cc.

2019 {
2020  return DBL_MAX;
2021 }
#define DBL_MAX
Definition: templates.hh:83
G4double G4VEnergyLossProcess::GetCSDARange ( G4double kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 808 of file G4VEnergyLossProcess.hh.

810 {
811  DefineMaterial(couple);
812  G4double x = DBL_MAX;
813  if(theCSDARangeTable) {
814  x=GetLimitScaledRangeForScaledEnergy(kineticEnergy*massRatio)*reduceFactor;
815  }
816  return x;
817 }
tuple x
Definition: test.py:50
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83

Here is the caller graph for this function:

const G4Element * G4VEnergyLossProcess::GetCurrentElement ( ) const

Definition at line 2261 of file G4VEnergyLossProcess.cc.

2262 {
2263  const G4Element* elm = nullptr;
2264  if(currentModel) { elm = currentModel->GetCurrentElement(); }
2265  return elm;
2266 }
const G4Element * GetCurrentElement() const
Definition: G4VEmModel.hh:469

Here is the call graph for this function:

G4double G4VEnergyLossProcess::GetDEDX ( G4double kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 771 of file G4VEnergyLossProcess.hh.

773 {
774  DefineMaterial(couple);
775  return GetDEDXForScaledEnergy(kineticEnergy*massRatio);
776 }

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::GetDEDXDispersion ( const G4MaterialCutsCouple couple,
const G4DynamicParticle dp,
G4double  length 
)

Definition at line 1944 of file G4VEnergyLossProcess.cc.

1948 {
1949  DefineMaterial(couple);
1950  G4double ekin = dp->GetKineticEnergy();
1951  SelectModel(ekin*massRatio);
1952  G4double tmax = currentModel->MaxSecondaryKinEnergy(dp);
1953  tmax = std::min(tmax,(*theCuts)[currentCoupleIndex]);
1954  G4double d = 0.0;
1955  G4VEmFluctuationModel* fm = currentModel->GetModelOfFluctuations();
1956  if(fm) { d = fm->Dispersion(currentMaterial,dp,tmax,length); }
1957  return d;
1958 }
G4double MaxSecondaryKinEnergy(const G4DynamicParticle *dynParticle)
Definition: G4VEmModel.hh:484
G4double GetKineticEnergy() const
void SelectModel(G4double kinEnergy)
G4VEmFluctuationModel * GetModelOfFluctuations()
Definition: G4VEmModel.hh:612
virtual G4double Dispersion(const G4Material *, const G4DynamicParticle *, G4double tmax, G4double length)=0
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::GetDEDXForSubsec ( G4double kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 781 of file G4VEnergyLossProcess.hh.

783 {
784  DefineMaterial(couple);
785  return GetSubDEDXForScaledEnergy(kineticEnergy*massRatio);
786 }

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::GetKineticEnergy ( G4double range,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 837 of file G4VEnergyLossProcess.hh.

839 {
840  DefineMaterial(couple);
841  return ScaledKinEnergyForLoss(range/reduceFactor)/massRatio;
842 }
const G4ParticleDefinition const G4Material *G4double range

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::GetLambda ( G4double kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 847 of file G4VEnergyLossProcess.hh.

849 {
850  DefineMaterial(couple);
851  G4double x = 0.0;
852  if(theLambdaTable) {
853  x = GetLambdaForScaledEnergy(kineticEnergy*massRatio);
854  }
855  return x;
856 }
tuple x
Definition: test.py:50
double G4double
Definition: G4Types.hh:76
G4double G4VEnergyLossProcess::GetMeanFreePath ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
overrideprotectedvirtual

Implements G4VContinuousDiscreteProcess.

Reimplemented in G4ePolarizedIonisation.

Definition at line 2004 of file G4VEnergyLossProcess.cc.

2009 {
2010  *condition = NotForced;
2011  return MeanFreePath(track);
2012 }
G4double condition(const G4ErrorSymMatrix &m)
G4double MeanFreePath(const G4Track &track)

Here is the call graph for this function:

Here is the caller graph for this function:

G4VEmModel * G4VEnergyLossProcess::GetModelByIndex ( G4int  idx = 0,
G4bool  ver = false 
) const

Definition at line 427 of file G4VEnergyLossProcess.cc.

428 {
429  return modelManager->GetModel(idx, ver);
430 }
G4VEmModel * GetModel(G4int idx, G4bool ver=false)

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::GetRange ( G4double kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 791 of file G4VEnergyLossProcess.hh.

793 {
794  G4double x = fRange;
795  DefineMaterial(couple);
796  if(theCSDARangeTable) {
797  x = GetLimitScaledRangeForScaledEnergy(kineticEnergy*massRatio)
798  * reduceFactor;
799  } else if(theRangeTableForLoss) {
800  x = GetScaledRangeForScaledEnergy(kineticEnergy*massRatio)*reduceFactor;
801  }
802  return x;
803 }
tuple x
Definition: test.py:50
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::GetRangeForLoss ( G4double kineticEnergy,
const G4MaterialCutsCouple couple 
)
inline

Definition at line 822 of file G4VEnergyLossProcess.hh.

824 {
825  // G4cout << "GetRangeForLoss: Range from " << GetProcessName() << G4endl;
826  DefineMaterial(couple);
827  G4double x =
828  GetScaledRangeForScaledEnergy(kineticEnergy*massRatio)*reduceFactor;
829  //G4cout << "GetRangeForLoss: Range from " << GetProcessName()
830  // << " e= " << kineticEnergy << " r= " << x << G4endl;
831  return x;
832 }
tuple x
Definition: test.py:50
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

virtual void G4VEnergyLossProcess::InitialiseEnergyLossProcess ( const G4ParticleDefinition ,
const G4ParticleDefinition  
)
protectedpure virtual
G4PhysicsTable * G4VEnergyLossProcess::InverseRangeTable ( ) const
inline

Definition at line 1063 of file G4VEnergyLossProcess.hh.

1064 {
1065  return theInverseRangeTable;
1066 }

Here is the caller graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::IonisationTable ( ) const
inline

Definition at line 1028 of file G4VEnergyLossProcess.hh.

1029 {
1030  return theIonisationTable;
1031 }

Here is the caller graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::IonisationTableForSubsec ( ) const
inline

Definition at line 1035 of file G4VEnergyLossProcess.hh.

1036 {
1037  return theIonisationSubTable;
1038 }

Here is the caller graph for this function:

G4bool G4VEnergyLossProcess::IsIntegral ( ) const
inline

Definition at line 958 of file G4VEnergyLossProcess.hh.

959 {
960  return integral;
961 }
G4bool G4VEnergyLossProcess::IsIonisationProcess ( ) const
inline

Definition at line 965 of file G4VEnergyLossProcess.hh.

966 {
967  return isIonisation;
968 }

Here is the caller graph for this function:

G4PhysicsVector * G4VEnergyLossProcess::LambdaPhysicsVector ( const G4MaterialCutsCouple ,
G4double  cut 
)
protected

Definition at line 2026 of file G4VEnergyLossProcess.cc.

2028 {
2029  G4PhysicsVector* v =
2030  new G4PhysicsLogVector(minKinEnergy, maxKinEnergy, nBins);
2031  v->SetSpline(theParameters->Spline());
2032  return v;
2033 }
G4bool Spline() const
void SetSpline(G4bool)
tuple v
Definition: test.py:18

Here is the call graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::LambdaTable ( ) const
inline

Definition at line 1070 of file G4VEnergyLossProcess.hh.

1071 {
1072  return theLambdaTable;
1073 }

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::MaxKinEnergy ( ) const
inline

Definition at line 986 of file G4VEnergyLossProcess.hh.

987 {
988  return maxKinEnergy;
989 }
G4double G4VEnergyLossProcess::MeanFreePath ( const G4Track track)

Definition at line 1983 of file G4VEnergyLossProcess.cc.

1984 {
1985  DefineMaterial(track.GetMaterialCutsCouple());
1986  preStepLambda = GetLambdaForScaledEnergy(track.GetKineticEnergy()*massRatio);
1987  G4double x = DBL_MAX;
1988  if(0.0 < preStepLambda) { x = 1.0/preStepLambda; }
1989  return x;
1990 }
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
tuple x
Definition: test.py:50
G4double GetKineticEnergy() const
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEnergyLossProcess::MinKinEnergy ( ) const
inline

Definition at line 979 of file G4VEnergyLossProcess.hh.

980 {
981  return minKinEnergy;
982 }
G4double G4VEnergyLossProcess::MinPrimaryEnergy ( const G4ParticleDefinition ,
const G4Material ,
G4double  cut 
)
protectedvirtual

Reimplemented in G4ionIonisation, G4eIonisation, G4MuIonisation, G4hIonisation, G4MuBremsstrahlung, G4ePolarizedIonisation, G4MuPairProduction, G4alphaIonisation, G4hhIonisation, and G4ePairProduction.

Definition at line 382 of file G4VEnergyLossProcess.cc.

385 {
386  return cut;
387 }

Here is the caller graph for this function:

G4int G4VEnergyLossProcess::NumberOfModels ( ) const

Definition at line 434 of file G4VEnergyLossProcess.cc.

435 {
436  return modelManager->NumberOfModels();
437 }
G4int NumberOfModels() const

Here is the call graph for this function:

G4int G4VEnergyLossProcess::NumberOfSubCutoffRegions ( ) const
inline

Definition at line 972 of file G4VEnergyLossProcess.hh.

973 {
974  return nSCoffRegions;
975 }
const G4ParticleDefinition * G4VEnergyLossProcess::Particle ( ) const
inline

Definition at line 920 of file G4VEnergyLossProcess.hh.

921 {
922  return particle;
923 }

Here is the caller graph for this function:

G4VParticleChange * G4VEnergyLossProcess::PostStepDoIt ( const G4Track track,
const G4Step step 
)
overridevirtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 1620 of file G4VEnergyLossProcess.cc.

1622 {
1623  // In all cases clear number of interaction lengths
1625  mfpKinEnergy = currentInteractionLength = DBL_MAX;
1626 
1628  G4double finalT = track.GetKineticEnergy();
1629  if(finalT <= lowestKinEnergy) { return &fParticleChange; }
1630 
1631  G4double postStepScaledEnergy = finalT*massRatio;
1632 
1633  if(!currentModel->IsActive(postStepScaledEnergy)) {
1634  return &fParticleChange;
1635  }
1636  /*
1637  if(-1 < verboseLevel) {
1638  G4cout << GetProcessName()
1639  << "::PostStepDoIt: E(MeV)= " << finalT/MeV
1640  << G4endl;
1641  }
1642  */
1643 
1644  // forced process - should happen only once per track
1645  if(biasFlag) {
1646  if(biasManager->ForcedInteractionRegion(currentCoupleIndex)) {
1647  biasFlag = false;
1648  }
1649  }
1650 
1651  // Integral approach
1652  if (integral) {
1653  G4double lx = GetLambdaForScaledEnergy(postStepScaledEnergy);
1654  /*
1655  if(preStepLambda<lx && 1 < verboseLevel) {
1656  G4cout << "WARNING: for " << particle->GetParticleName()
1657  << " and " << GetProcessName()
1658  << " E(MeV)= " << finalT/MeV
1659  << " preLambda= " << preStepLambda
1660  << " < " << lx << " (postLambda) "
1661  << G4endl;
1662  }
1663  */
1664  if(lx <= 0.0 || preStepLambda*G4UniformRand() > lx) {
1665  return &fParticleChange;
1666  }
1667  }
1668 
1669  SelectModel(postStepScaledEnergy);
1670 
1671  // define new weight for primary and secondaries
1673  if(weightFlag) {
1674  weight /= biasFactor;
1676  }
1677 
1678  const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
1679  G4double tcut = (*theCuts)[currentCoupleIndex];
1680 
1681  // sample secondaries
1682  secParticles.clear();
1683  //G4cout<< "@@@ Eprimary= "<<dynParticle->GetKineticEnergy()/MeV
1684  // << " cut= " << tcut/MeV << G4endl;
1685  currentModel->SampleSecondaries(&secParticles, currentCouple,
1686  dynParticle, tcut);
1687 
1688  G4int num0 = secParticles.size();
1689 
1690  // bremsstrahlung splitting or Russian roulette
1691  if(biasManager) {
1692  if(biasManager->SecondaryBiasingRegion(currentCoupleIndex)) {
1693  G4double eloss = 0.0;
1694  weight *= biasManager->ApplySecondaryBiasing(
1695  secParticles,
1696  track, currentModel,
1697  &fParticleChange, eloss,
1698  currentCoupleIndex, tcut,
1699  step.GetPostStepPoint()->GetSafety());
1700  if(eloss > 0.0) {
1703  }
1704  }
1705  }
1706 
1707  // save secondaries
1708  G4int num = secParticles.size();
1709  if(num > 0) {
1710 
1712  G4double time = track.GetGlobalTime();
1713 
1714  for (G4int i=0; i<num; ++i) {
1715  if(secParticles[i]) {
1716  G4Track* t = new G4Track(secParticles[i], time, track.GetPosition());
1718  t->SetWeight(weight);
1719  if(i < num0) { t->SetCreatorModelIndex(secID); }
1720  else { t->SetCreatorModelIndex(biasID); }
1721 
1722  //G4cout << "Secondary(post step) has weight " << t->GetWeight()
1723  // << ", kenergy " << t->GetKineticEnergy()/MeV << " MeV"
1724  // << " time= " << time/ns << " ns " << G4endl;
1726  }
1727  }
1728  }
1729 
1732  if(particle->GetProcessManager()->GetAtRestProcessVector()->size() > 0)
1735  }
1736 
1737  /*
1738  if(-1 < verboseLevel) {
1739  G4cout << "::PostStepDoIt: Sample secondary; Efin= "
1740  << fParticleChange.GetProposedKineticEnergy()/MeV
1741  << " MeV; model= (" << currentModel->LowEnergyLimit()
1742  << ", " << currentModel->HighEnergyLimit() << ")"
1743  << " preStepLambda= " << preStepLambda
1744  << " dir= " << track.GetMomentumDirection()
1745  << " status= " << track.GetTrackStatus()
1746  << G4endl;
1747  }
1748  */
1749  return &fParticleChange;
1750 }
void InitializeForPostStep(const G4Track &)
G4double ApplySecondaryBiasing(std::vector< G4DynamicParticle * > &, const G4Track &track, G4VEmModel *currentModel, G4ParticleChangeForGamma *pParticleChange, G4double &eloss, G4int coupleIdx, G4double tcut, G4double safety=0.0)
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
const G4DynamicParticle * GetDynamicParticle() const
G4bool SecondaryBiasingRegion(G4int coupleIdx)
const G4ThreeVector & GetPosition() const
G4bool ForcedInteractionRegion(G4int coupleIdx)
G4ParticleChangeForLoss fParticleChange
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
void SelectModel(G4double kinEnergy)
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4double GetParentWeight() const
G4double GetProposedKineticEnergy() const
int G4int
Definition: G4Types.hh:78
void SetWeight(G4double aValue)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void SetCreatorModelIndex(G4int idx)
void ProposeWeight(G4double finalWeight)
G4ProcessVector * GetAtRestProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4double GetKineticEnergy() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4double currentInteractionLength
Definition: G4VProcess.hh:297
G4double GetGlobalTime() const
void AddSecondary(G4Track *aSecondary)
const G4TouchableHandle & GetTouchableHandle() const
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:760
G4int size() const
G4ProcessManager * GetProcessManager() const
void SetNumberOfSecondaries(G4int totSecondaries)
G4double GetLocalEnergyDeposit() const
G4StepPoint * GetPostStepPoint() const
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
G4double GetSafety() const
G4TrackStatus GetTrackStatus() const
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

G4double G4VEnergyLossProcess::PostStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
overridevirtual

Reimplemented from G4VContinuousDiscreteProcess.

Reimplemented in G4ePolarizedIonisation.

Definition at line 1153 of file G4VEnergyLossProcess.cc.

1157 {
1158  // condition is set to "Not Forced"
1159  *condition = NotForced;
1160  G4double x = DBL_MAX;
1161 
1162  // initialisation of material, mass, charge, model
1163  // at the beginning of the step
1164  DefineMaterial(track.GetMaterialCutsCouple());
1165  preStepKinEnergy = track.GetKineticEnergy();
1166  preStepScaledEnergy = preStepKinEnergy*massRatio;
1167  SelectModel(preStepScaledEnergy);
1168 
1169  if(!currentModel->IsActive(preStepScaledEnergy)) {
1172  return x;
1173  }
1174 
1175  // change effective charge of an ion on fly
1176  if(isIon) {
1177  G4double q2 = currentModel->ChargeSquareRatio(track);
1178  if(q2 != chargeSqRatio && q2 > 0.0) {
1179  chargeSqRatio = q2;
1180  fFactor = q2*biasFactor*(*theDensityFactor)[currentCoupleIndex];
1181  reduceFactor = 1.0/(fFactor*massRatio);
1182  }
1183  }
1184  // if(particle->GetPDGMass() > 0.9*GeV)
1185  //G4cout << "q2= "<<chargeSqRatio << " massRatio= " << massRatio << G4endl;
1186 
1187  // forced biasing only for primary particles
1188  if(biasManager) {
1189  if(0 == track.GetParentID()) {
1190  if(biasFlag &&
1191  biasManager->ForcedInteractionRegion(currentCoupleIndex)) {
1192  return biasManager->GetStepLimit(currentCoupleIndex, previousStepSize);
1193  }
1194  }
1195  }
1196 
1197  // compute mean free path
1198  if(preStepScaledEnergy < mfpKinEnergy) {
1199  if (integral) { ComputeLambdaForScaledEnergy(preStepScaledEnergy); }
1200  else { preStepLambda = GetLambdaForScaledEnergy(preStepScaledEnergy); }
1201 
1202  // zero cross section
1203  if(preStepLambda <= 0.0) {
1206  }
1207  }
1208 
1209  // non-zero cross section
1210  if(preStepLambda > 0.0) {
1212 
1213  // beggining of tracking (or just after DoIt of this process)
1216 
1217  } else if(currentInteractionLength < DBL_MAX) {
1218 
1219  // subtract NumberOfInteractionLengthLeft using previous step
1221  previousStepSize/currentInteractionLength;
1222 
1225  }
1226 
1227  // new mean free path and step limit
1228  currentInteractionLength = 1.0/preStepLambda;
1230  }
1231 #ifdef G4VERBOSE
1232  if (verboseLevel>2){
1233  // if(particle->GetPDGMass() > 0.9*GeV){
1234  G4cout << "G4VEnergyLossProcess::PostStepGetPhysicalInteractionLength ";
1235  G4cout << "[ " << GetProcessName() << "]" << G4endl;
1236  G4cout << " for " << track.GetDefinition()->GetParticleName()
1237  << " in Material " << currentMaterial->GetName()
1238  << " Ekin(MeV)= " << preStepKinEnergy/MeV
1239  << " " << track.GetMaterial()->GetName()
1240  <<G4endl;
1241  G4cout << "MeanFreePath = " << currentInteractionLength/cm << "[cm]"
1242  << "InteractionLength= " << x/cm <<"[cm] " <<G4endl;
1243  }
1244 #endif
1245  return x;
1246 }
G4double condition(const G4ErrorSymMatrix &m)
G4ParticleDefinition * GetDefinition() const
G4int GetParentID() const
G4int verboseLevel
Definition: G4VProcess.hh:368
const G4String & GetName() const
Definition: G4Material.hh:178
G4bool ForcedInteractionRegion(G4int coupleIdx)
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
void SelectModel(G4double kinEnergy)
tuple x
Definition: test.py:50
const G4String & GetParticleName() const
G4double GetKineticEnergy() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
static constexpr double cm
Definition: G4SIunits.hh:119
G4double GetStepLimit(G4int coupleIdx, G4double previousStep)
G4double currentInteractionLength
Definition: G4VProcess.hh:297
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
virtual G4double ChargeSquareRatio(const G4Track &)
Definition: G4VEmModel.cc:345
G4Material * GetMaterial() const
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:760
G4double G4Log(G4double x)
Definition: G4Log.hh:230
T max(const T t1, const T t2)
brief Return the largest of the two arguments
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
G4double theInitialNumberOfInteractionLength
Definition: G4VProcess.hh:300
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::PreparePhysicsTable ( const G4ParticleDefinition part)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 442 of file G4VEnergyLossProcess.cc.

443 {
444  if(1 < verboseLevel) {
445  G4cout << "G4VEnergyLossProcess::PreparePhysicsTable for "
446  << GetProcessName() << " for " << part.GetParticleName()
447  << " " << this << G4endl;
448  }
449 
450  const G4VEnergyLossProcess* masterProcess =
451  static_cast<const G4VEnergyLossProcess*>(GetMasterProcess());
452  if(masterProcess && masterProcess != this) { isMaster = false; }
453 
454  currentCouple = nullptr;
455  preStepLambda = 0.0;
456  mfpKinEnergy = DBL_MAX;
457  fRange = DBL_MAX;
458  preStepKinEnergy = 0.0;
459  preStepRangeEnergy = 0.0;
460  chargeSqRatio = 1.0;
461  massRatio = 1.0;
462  reduceFactor = 1.0;
463  fFactor = 1.0;
464  lastIdx = 0;
465 
466  // Are particle defined?
467  if( !particle ) { particle = &part; }
468 
469  if(part.GetParticleType() == "nucleus") {
470 
471  G4String pname = part.GetParticleName();
472  if(pname != "deuteron" && pname != "triton" &&
473  pname != "alpha+" && pname != "helium" &&
474  pname != "hydrogen") {
475 
476  if(!theGenericIon) {
477  theGenericIon =
479  }
480  isIon = true;
481  if(theGenericIon && particle != theGenericIon) {
482  G4ProcessManager* pm = theGenericIon->GetProcessManager();
484  size_t n = v->size();
485  for(size_t j=0; j<n; ++j) {
486  if((*v)[j] == this) {
487  particle = theGenericIon;
488  break;
489  }
490  }
491  }
492  }
493  }
494 
495  if( particle != &part ) {
496  if(!isIon) {
497  lManager->RegisterExtraParticle(&part, this);
498  }
499  if(1 < verboseLevel) {
500  G4cout << "### G4VEnergyLossProcess::PreparePhysicsTable()"
501  << " interrupted for "
502  << part.GetParticleName() << " isIon= " << isIon
503  << " particle " << particle << " GenericIon " << theGenericIon
504  << G4endl;
505  }
506  return;
507  }
508 
509  Clean();
510  lManager->PreparePhysicsTable(&part, this, isMaster);
511  G4LossTableBuilder* bld = lManager->GetTableBuilder();
512 
513  // Base particle and set of models can be defined here
514  InitialiseEnergyLossProcess(particle, baseParticle);
515 
516  const G4ProductionCutsTable* theCoupleTable=
518  size_t n = theCoupleTable->GetTableSize();
519 
520  theDEDXAtMaxEnergy.resize(n, 0.0);
521  theRangeAtMaxEnergy.resize(n, 0.0);
522  theEnergyOfCrossSectionMax.resize(n, 0.0);
523  theCrossSectionMax.resize(n, DBL_MAX);
524 
525  // parameters of the process
526  if(!actIntegral) { integral = theParameters->Integral(); }
527  if(!actLossFluc) { lossFluctuationFlag = theParameters->LossFluctuation(); }
528  rndmStepFlag = theParameters->UseCutAsFinalRange();
529  if(!actMinKinEnergy) { minKinEnergy = theParameters->MinKinEnergy(); }
530  if(!actMaxKinEnergy) { maxKinEnergy = theParameters->MaxKinEnergy(); }
531  if(!actBinning) {
532  nBins = theParameters->NumberOfBinsPerDecade()
533  *G4lrint(std::log10(maxKinEnergy/minKinEnergy));
534  }
535  maxKinEnergyCSDA = theParameters->MaxEnergyForCSDARange();
536  nBinsCSDA = theParameters->NumberOfBinsPerDecade()
537  *G4lrint(std::log10(maxKinEnergyCSDA/minKinEnergy));
538  if(!actLinLossLimit) { linLossLimit = theParameters->LinearLossLimit(); }
539  lambdaFactor = theParameters->LambdaFactor();
540  if(isMaster) { SetVerboseLevel(theParameters->Verbose()); }
541  else { SetVerboseLevel(theParameters->WorkerVerbose()); }
542 
543  G4bool isElec = true;
544  if(particle->GetPDGMass() > CLHEP::MeV) { isElec = false; }
545  theParameters->DefineRegParamForLoss(this, isElec);
546 
547  G4double initialCharge = particle->GetPDGCharge();
548  G4double initialMass = particle->GetPDGMass();
549 
550  if (baseParticle) {
551  massRatio = (baseParticle->GetPDGMass())/initialMass;
552  G4double q = initialCharge/baseParticle->GetPDGCharge();
553  chargeSqRatio = q*q;
554  if(chargeSqRatio > 0.0) { reduceFactor = 1.0/(chargeSqRatio*massRatio); }
555  }
556  if(initialMass < MeV) {
557  lowestKinEnergy = theParameters->LowestElectronEnergy();
558  } else {
559  lowestKinEnergy = theParameters->LowestMuHadEnergy();
560  }
561 
562  // Tables preparation
563  if (isMaster && !baseParticle) {
564 
565  if(theDEDXTable && isIonisation) {
566  if(theIonisationTable && theDEDXTable != theIonisationTable) {
567  theDEDXTable->clearAndDestroy();
568  delete theDEDXTable;
569  theDEDXTable = theIonisationTable;
570  }
571  if(theDEDXSubTable && theIonisationSubTable &&
572  theDEDXSubTable != theIonisationSubTable) {
573  theDEDXSubTable->clearAndDestroy();
574  delete theDEDXSubTable;
575  theDEDXSubTable = theIonisationSubTable;
576  }
577  }
578 
579  theDEDXTable = G4PhysicsTableHelper::PreparePhysicsTable(theDEDXTable);
580  bld->InitialiseBaseMaterials(theDEDXTable);
581 
582  if(theDEDXSubTable) {
583  theDEDXSubTable =
585  }
586 
587  if (theParameters->BuildCSDARange()) {
588  theDEDXunRestrictedTable =
589  G4PhysicsTableHelper::PreparePhysicsTable(theDEDXunRestrictedTable);
590  theCSDARangeTable =
592  }
593 
594  theLambdaTable = G4PhysicsTableHelper::PreparePhysicsTable(theLambdaTable);
595 
596  if(isIonisation) {
597  theRangeTableForLoss =
598  G4PhysicsTableHelper::PreparePhysicsTable(theRangeTableForLoss);
599  theInverseRangeTable =
600  G4PhysicsTableHelper::PreparePhysicsTable(theInverseRangeTable);
601  }
602 
603  if (nSCoffRegions && !lManager->SubCutProducer()) {
604  theDEDXSubTable =
606  theSubLambdaTable =
608  }
609  }
610  /*
611  G4cout << "** G4VEnergyLossProcess::PreparePhysicsTable() for "
612  << GetProcessName() << " and " << particle->GetParticleName()
613  << " isMaster: " << isMaster << " isIonisation: "
614  << isIonisation << G4endl;
615  G4cout << " theDEDX: " << theDEDXTable
616  << " theRange: " << theRangeTableForLoss
617  << " theInverse: " << theInverseRangeTable
618  << " theLambda: " << theLambdaTable << G4endl;
619  */
620  // forced biasing
621  if(biasManager) {
622  biasManager->Initialise(part,GetProcessName(),verboseLevel);
623  biasFlag = false;
624  }
625 
626  // defined ID of secondary particles
627  if(isMaster) {
628  G4String nam1 = GetProcessName();
629  G4String nam4 = nam1 + "_split";
630  G4String nam5 = nam1 + "_subcut";
631  secID = G4PhysicsModelCatalog::Register(nam1);
632  biasID = G4PhysicsModelCatalog::Register(nam4);
633  subsecID= G4PhysicsModelCatalog::Register(nam5);
634  }
635 
636  // initialisation of models
637  G4int nmod = modelManager->NumberOfModels();
638  for(G4int i=0; i<nmod; ++i) {
639  G4VEmModel* mod = modelManager->GetModel(i);
640  mod->SetMasterThread(isMaster);
642  theParameters->UseAngularGeneratorForIonisation());
643  if(mod->HighEnergyLimit() > maxKinEnergy) {
644  mod->SetHighEnergyLimit(maxKinEnergy);
645  }
646  }
647  theCuts = modelManager->Initialise(particle, secondaryParticle,
648  theParameters->MinSubRange(),
649  verboseLevel);
650 
651  // Sub Cutoff
652  if(nSCoffRegions > 0) {
653  if(theParameters->MinSubRange() < 1.0) { useSubCutoff = true; }
654 
655  theSubCuts = modelManager->SubCutoff();
656 
657  idxSCoffRegions = new G4bool[n];
658  for (size_t j=0; j<n; ++j) {
659 
660  const G4MaterialCutsCouple* couple =
661  theCoupleTable->GetMaterialCutsCouple(j);
662  const G4ProductionCuts* pcuts = couple->GetProductionCuts();
663 
664  G4bool reg = false;
665  for(G4int i=0; i<nSCoffRegions; ++i) {
666  if( pcuts == scoffRegions[i]->GetProductionCuts()) {
667  reg = true;
668  break;
669  }
670  }
671  idxSCoffRegions[j] = reg;
672  }
673  }
674 
675  if(1 < verboseLevel) {
676  G4cout << "G4VEnergyLossProcess::PrepearPhysicsTable() is done "
677  << " for local " << particle->GetParticleName()
678  << " isIon= " << isIon;
679  if(baseParticle) {
680  G4cout << "; base: " << baseParticle->GetParticleName();
681  }
682  G4cout << " chargeSqRatio= " << chargeSqRatio
683  << " massRatio= " << massRatio
684  << " reduceFactor= " << reduceFactor << G4endl;
685  if (nSCoffRegions) {
686  G4cout << " SubCutoff Regime is ON for regions: " << G4endl;
687  for (G4int i=0; i<nSCoffRegions; ++i) {
688  const G4Region* r = scoffRegions[i];
689  G4cout << " " << r->GetName() << G4endl;
690  }
691  }
692  }
693 }
G4bool UseCutAsFinalRange() const
G4VSubCutProducer * SubCutProducer()
G4int NumberOfBinsPerDecade() const
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
const G4String & GetName() const
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
G4int WorkerVerbose() const
G4double MaxKinEnergy() const
G4int verboseLevel
Definition: G4VProcess.hh:368
G4VEmModel * GetModel(G4int idx, G4bool ver=false)
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:636
G4double LowestElectronEnergy() const
G4bool Integral() const
void DefineRegParamForLoss(G4VEnergyLossProcess *, G4bool isElectron) const
const G4DataVector * SubCutoff() const
int G4int
Definition: G4Types.hh:78
static G4PhysicsTable * PreparePhysicsTable(G4PhysicsTable *physTable)
G4double MinSubRange() const
const G4String & GetParticleName() const
G4LossTableBuilder * GetTableBuilder()
static const G4double reg
void Initialise(const G4ParticleDefinition &part, const G4String &procName, G4int verbose)
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:732
G4bool BuildCSDARange() const
void SetAngularGeneratorFlag(G4bool)
Definition: G4VEmModel.hh:706
G4double LinearLossLimit() const
G4double LambdaFactor() const
G4GLOB_DLL std::ostream G4cout
G4int Verbose() const
G4double LowestMuHadEnergy() const
bool G4bool
Definition: G4Types.hh:79
static constexpr double MeV
const G4String & GetParticleType() const
G4double MinKinEnergy() const
const G4int n
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
tuple v
Definition: test.py:18
string pname
Definition: eplot.py:33
void RegisterExtraParticle(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
const G4DataVector * Initialise(const G4ParticleDefinition *part, const G4ParticleDefinition *secPart, G4double minSubRange, G4int verb)
void SetMasterThread(G4bool val)
Definition: G4VEmModel.hh:718
G4int size() const
static G4ProductionCutsTable * GetProductionCutsTable()
G4int NumberOfModels() const
G4bool LossFluctuation() const
G4double GetPDGMass() const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
static G4ParticleTable * GetParticleTable()
int G4lrint(double ad)
Definition: templates.hh:163
G4ProcessManager * GetProcessManager() const
void InitialiseBaseMaterials(G4PhysicsTable *table)
G4bool UseAngularGeneratorForIonisation() const
static G4int Register(const G4String &)
void PreparePhysicsTable(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p, G4bool theMaster)
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76
G4ProcessVector * GetAlongStepProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4double MaxEnergyForCSDARange() const
G4double GetPDGCharge() const
G4ProductionCuts * GetProductionCuts() const
#define DBL_MAX
Definition: templates.hh:83
void clearAndDestroy()
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:437
virtual void InitialiseEnergyLossProcess(const G4ParticleDefinition *, const G4ParticleDefinition *)=0
virtual void G4VEnergyLossProcess::PrintInfo ( )
pure virtual
void G4VEnergyLossProcess::PrintInfoDefinition ( const G4ParticleDefinition part)

Definition at line 975 of file G4VEnergyLossProcess.cc.

976 {
977  if(0 < verboseLevel) {
978  G4cout << std::setprecision(6);
979  G4cout << G4endl << GetProcessName() << ": for "
980  << part.GetParticleName()
981  << " SubType= " << GetProcessSubType()
982  << G4endl;
983  G4cout << " dE/dx and range tables from "
984  << G4BestUnit(minKinEnergy,"Energy")
985  << " to " << G4BestUnit(maxKinEnergy,"Energy")
986  << " in " << nBins << " bins" << G4endl
987  << " Lambda tables from threshold to "
988  << G4BestUnit(maxKinEnergy,"Energy")
989  << ", " << theParameters->NumberOfBinsPerDecade()
990  << " bins per decade, spline: "
991  << theParameters->Spline()
992  << G4endl;
993  if(theRangeTableForLoss && isIonisation) {
994  G4cout << " finalRange(mm)= " << finalRange/mm
995  << ", dRoverRange= " << dRoverRange
996  << ", integral: " << integral
997  << ", fluct: " << lossFluctuationFlag
998  << ", linLossLimit= " << linLossLimit
999  << G4endl;
1000  }
1001  PrintInfo();
1002  modelManager->DumpModelList(verboseLevel);
1003  if(theCSDARangeTable && isIonisation) {
1004  G4cout << " CSDA range table up"
1005  << " to " << G4BestUnit(maxKinEnergyCSDA,"Energy")
1006  << " in " << nBinsCSDA << " bins" << G4endl;
1007  }
1008  if(nSCoffRegions>0 && isIonisation) {
1009  G4cout << " Subcutoff sampling in " << nSCoffRegions
1010  << " regions" << G4endl;
1011  }
1012  if(2 < verboseLevel) {
1013  G4cout << " DEDXTable address= " << theDEDXTable << G4endl;
1014  if(theDEDXTable && isIonisation) G4cout << (*theDEDXTable) << G4endl;
1015  G4cout << "non restricted DEDXTable address= "
1016  << theDEDXunRestrictedTable << G4endl;
1017  if(theDEDXunRestrictedTable && isIonisation) {
1018  G4cout << (*theDEDXunRestrictedTable) << G4endl;
1019  }
1020  if(theDEDXSubTable && isIonisation) {
1021  G4cout << (*theDEDXSubTable) << G4endl;
1022  }
1023  G4cout << " CSDARangeTable address= " << theCSDARangeTable
1024  << G4endl;
1025  if(theCSDARangeTable && isIonisation) {
1026  G4cout << (*theCSDARangeTable) << G4endl;
1027  }
1028  G4cout << " RangeTableForLoss address= " << theRangeTableForLoss
1029  << G4endl;
1030  if(theRangeTableForLoss && isIonisation) {
1031  G4cout << (*theRangeTableForLoss) << G4endl;
1032  }
1033  G4cout << " InverseRangeTable address= " << theInverseRangeTable
1034  << G4endl;
1035  if(theInverseRangeTable && isIonisation) {
1036  G4cout << (*theInverseRangeTable) << G4endl;
1037  }
1038  G4cout << " LambdaTable address= " << theLambdaTable << G4endl;
1039  if(theLambdaTable && isIonisation) {
1040  G4cout << (*theLambdaTable) << G4endl;
1041  }
1042  G4cout << " SubLambdaTable address= " << theSubLambdaTable
1043  << G4endl;
1044  if(theSubLambdaTable && isIonisation) {
1045  G4cout << (*theSubLambdaTable) << G4endl;
1046  }
1047  }
1048  }
1049 }
G4int NumberOfBinsPerDecade() const
G4bool Spline() const
static constexpr double mm
Definition: G4SIunits.hh:115
G4int verboseLevel
Definition: G4VProcess.hh:368
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void DumpModelList(G4int verb)
virtual void PrintInfo()=0
#define G4endl
Definition: G4ios.hh:61
G4int GetProcessSubType() const
Definition: G4VProcess.hh:426

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::ProcessDescription ( std::ostream &  outFile) const
virtual

Definition at line 2425 of file G4VEnergyLossProcess.cc.

2426 {
2427  outFile << "Energy loss process <" << GetProcessName() << ">" << G4endl;
2428 }
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::RangeTableForLoss ( ) const
inline

Definition at line 1056 of file G4VEnergyLossProcess.hh.

1057 {
1058  return theRangeTableForLoss;
1059 }

Here is the caller graph for this function:

G4bool G4VEnergyLossProcess::RetrievePhysicsTable ( const G4ParticleDefinition part,
const G4String directory,
G4bool  ascii 
)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 1817 of file G4VEnergyLossProcess.cc.

1820 {
1821  G4bool res = true;
1822  if (!isMaster) return res;
1823  const G4String particleName = part->GetParticleName();
1824 
1825  if(1 < verboseLevel) {
1826  G4cout << "G4VEnergyLossProcess::RetrievePhysicsTable() for "
1827  << particleName << " and process " << GetProcessName()
1828  << "; tables_are_built= " << tablesAreBuilt
1829  << G4endl;
1830  }
1831  if(particle == part) {
1832 
1833  if ( !baseParticle ) {
1834 
1835  G4bool fpi = true;
1836  if(!RetrieveTable(part,theDEDXTable,ascii,directory,"DEDX",fpi))
1837  {fpi = false;}
1838 
1839  // ionisation table keeps individual dEdx and not sum of sub-processes
1840  if(!RetrieveTable(part,theDEDXTable,ascii,directory,"Ionisation",false))
1841  {fpi = false;}
1842 
1843  if(!RetrieveTable(part,theRangeTableForLoss,ascii,directory,"Range",fpi))
1844  {res = false;}
1845 
1846  if(!RetrieveTable(part,theDEDXunRestrictedTable,ascii,directory,
1847  "DEDXnr",false))
1848  {res = false;}
1849 
1850  if(!RetrieveTable(part,theCSDARangeTable,ascii,directory,
1851  "CSDARange",false))
1852  {res = false;}
1853 
1854  if(!RetrieveTable(part,theInverseRangeTable,ascii,directory,
1855  "InverseRange",fpi))
1856  {res = false;}
1857 
1858  if(!RetrieveTable(part,theLambdaTable,ascii,directory,"Lambda",true))
1859  {res = false;}
1860 
1861  G4bool yes = false;
1862  if(nSCoffRegions > 0) {yes = true;}
1863 
1864  if(!RetrieveTable(part,theDEDXSubTable,ascii,directory,"SubDEDX",yes))
1865  {res = false;}
1866 
1867  if(!RetrieveTable(part,theSubLambdaTable,ascii,directory,
1868  "SubLambda",yes))
1869  {res = false;}
1870 
1871  if(!fpi) yes = false;
1872  if(!RetrieveTable(part,theIonisationSubTable,ascii,directory,
1873  "SubIonisation",yes))
1874  {res = false;}
1875  }
1876  }
1877 
1878  return res;
1879 }
G4int verboseLevel
Definition: G4VProcess.hh:368
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

G4double G4VEnergyLossProcess::SampleSubCutSecondaries ( std::vector< G4Track * > &  tracks,
const G4Step step,
G4VEmModel model,
G4int  matIdx 
)

Definition at line 1542 of file G4VEnergyLossProcess.cc.

1546 {
1547  // Fast check weather subcutoff can work
1548  G4double esec = 0.0;
1549  G4double subcut = (*theSubCuts)[idx];
1550  G4double cut = (*theCuts)[idx];
1551  if(cut <= subcut) { return esec; }
1552 
1553  const G4Track* track = step.GetTrack();
1554  const G4DynamicParticle* dp = track->GetDynamicParticle();
1555  G4double e = dp->GetKineticEnergy()*massRatio;
1556  G4double cross = (*theDensityFactor)[idx]*chargeSqRatio
1557  *(((*theSubLambdaTable)[(*theDensityIdx)[idx]])->Value(e, idxSubLambda));
1558  G4double length = step.GetStepLength();
1559 
1560  // negligible probability to get any interaction
1561  if(length*cross < perMillion) { return esec; }
1562  /*
1563  if(-1 < verboseLevel)
1564  G4cout << "<<< Subcutoff for " << GetProcessName()
1565  << " cross(1/mm)= " << cross*mm << ">>>"
1566  << " e(MeV)= " << preStepScaledEnergy
1567  << " matIdx= " << currentCoupleIndex
1568  << G4endl;
1569  */
1570 
1571  // Sample subcutoff secondaries
1572  G4StepPoint* preStepPoint = step.GetPreStepPoint();
1573  G4StepPoint* postStepPoint = step.GetPostStepPoint();
1574  G4ThreeVector prepoint = preStepPoint->GetPosition();
1575  G4ThreeVector dr = postStepPoint->GetPosition() - prepoint;
1576  G4double pretime = preStepPoint->GetGlobalTime();
1577  G4double dt = postStepPoint->GetGlobalTime() - pretime;
1578  G4double fragment = 0.0;
1579 
1580  do {
1581  G4double del = -G4Log(G4UniformRand())/cross;
1582  fragment += del/length;
1583  if (fragment > 1.0) { break; }
1584 
1585  // sample secondaries
1586  secParticles.clear();
1587  model->SampleSecondaries(&secParticles,track->GetMaterialCutsCouple(),
1588  dp,subcut,cut);
1589 
1590  // position of subcutoff particles
1591  G4ThreeVector r = prepoint + fragment*dr;
1592  std::vector<G4DynamicParticle*>::iterator it;
1593  for(it=secParticles.begin(); it!=secParticles.end(); ++it) {
1594 
1595  G4Track* t = new G4Track((*it), pretime + fragment*dt, r);
1597  t->SetCreatorModelIndex(subsecID);
1598  tracks.push_back(t);
1599  esec += t->GetKineticEnergy();
1600  if (t->GetParticleDefinition() == thePositron) {
1601  esec += 2.0*electron_mass_c2;
1602  }
1603 
1604  /*
1605  if(-1 < verboseLevel)
1606  G4cout << "New track "
1607  << t->GetParticleDefinition()->GetParticleName()
1608  << " e(keV)= " << t->GetKineticEnergy()/keV
1609  << " fragment= " << fragment
1610  << G4endl;
1611  */
1612  }
1613  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
1614  } while (fragment <= 1.0);
1615  return esec;
1616 }
static constexpr double perMillion
Definition: G4SIunits.hh:334
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
G4double GetKineticEnergy() const
G4double GetStepLength() const
const G4DynamicParticle * GetDynamicParticle() const
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
void SetCreatorModelIndex(G4int idx)
G4StepPoint * GetPreStepPoint() const
G4double GetKineticEnergy() const
#define G4UniformRand()
Definition: Randomize.hh:97
const G4ThreeVector & GetPosition() const
const G4ParticleDefinition * GetParticleDefinition() const
float electron_mass_c2
Definition: hepunit.py:274
const G4TouchableHandle & GetTouchableHandle() const
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4StepPoint * GetPostStepPoint() const
G4double GetGlobalTime() const
double G4double
Definition: G4Types.hh:76
G4Track * GetTrack() const

Here is the call graph for this function:

Here is the caller graph for this function:

const G4ParticleDefinition * G4VEnergyLossProcess::SecondaryParticle ( ) const
inline

Definition at line 935 of file G4VEnergyLossProcess.hh.

936 {
937  return secondaryParticle;
938 }
G4PhysicsTable * G4VEnergyLossProcess::SecondaryRangeTable ( ) const
inline

Definition at line 1049 of file G4VEnergyLossProcess.hh.

1050 {
1051  return theSecondaryRangeTable;
1052 }

Here is the caller graph for this function:

void G4VEnergyLossProcess::SelectModel ( G4double  kinEnergy)
inlineprotected

Definition at line 619 of file G4VEnergyLossProcess.hh.

620 {
621  currentModel = modelManager->SelectModel(kinEnergy, currentCoupleIndex);
622  currentModel->SetCurrentCouple(currentCouple);
623 }
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:447
G4VEmModel * SelectModel(G4double &energy, size_t &index)

Here is the call graph for this function:

Here is the caller graph for this function:

G4VEmModel * G4VEnergyLossProcess::SelectModelForMaterial ( G4double  kinEnergy,
size_t &  idx 
) const
inline

Definition at line 627 of file G4VEnergyLossProcess.hh.

629 {
630  return modelManager->SelectModel(kinEnergy, idx);
631 }
G4VEmModel * SelectModel(G4double &energy, size_t &index)

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetBaseParticle ( const G4ParticleDefinition p)
inline

Definition at line 913 of file G4VEnergyLossProcess.hh.

914 {
915  baseParticle = p;
916 }
const char * p
Definition: xmltok.h:285

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetCrossSectionBiasingFactor ( G4double  f,
G4bool  flag = true 
)

Definition at line 2270 of file G4VEnergyLossProcess.cc.

2272 {
2273  if(f > 0.0) {
2274  biasFactor = f;
2275  weightFlag = flag;
2276  if(1 < verboseLevel) {
2277  G4cout << "### SetCrossSectionBiasingFactor: for "
2278  << " process " << GetProcessName()
2279  << " biasFactor= " << f << " weightFlag= " << flag
2280  << G4endl;
2281  }
2282  }
2283 }
G4int verboseLevel
Definition: G4VProcess.hh:368
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetCSDARangeTable ( G4PhysicsTable pRange)

Definition at line 2125 of file G4VEnergyLossProcess.cc.

2126 {
2127  theCSDARangeTable = p;
2128 
2129  if(p) {
2130  size_t n = p->length();
2131  G4PhysicsVector* pv;
2132  G4double emax = maxKinEnergyCSDA;
2133 
2134  for (size_t i=0; i<n; ++i) {
2135  pv = (*p)[i];
2136  G4double rmax = 0.0;
2137  if(pv) { rmax = pv->Value(emax, idxCSDA); }
2138  else {
2139  pv = (*p)[(*theDensityIdx)[i]];
2140  if(pv) { rmax = pv->Value(emax, idxCSDA)/(*theDensityFactor)[i]; }
2141  }
2142  theRangeAtMaxEnergy[i] = rmax;
2143  //G4cout << "i= " << i << " Emax(MeV)= " << emax/MeV << " Rmax= "
2144  //<< rmax<< G4endl;
2145  }
2146  }
2147 }
const char * p
Definition: xmltok.h:285
const G4int n
G4double Value(G4double theEnergy, size_t &lastidx) const
static const G4double emax
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetDEDXBinning ( G4int  nbins)

Definition at line 2377 of file G4VEnergyLossProcess.cc.

2378 {
2379  if(2 < n && n < 1000000000) {
2380  nBins = n;
2381  actBinning = true;
2382  } else {
2383  G4double e = (G4double)n;
2384  PrintWarning("SetDEDXBinning", e);
2385  }
2386 }
const G4int n
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetDEDXTable ( G4PhysicsTable p,
G4EmTableType  tType 
)

Definition at line 2064 of file G4VEnergyLossProcess.cc.

2065 {
2066  if(fTotal == tType) {
2067  theDEDXunRestrictedTable = p;
2068  if(p) {
2069  size_t n = p->length();
2070  G4PhysicsVector* pv = (*p)[0];
2071  G4double emax = maxKinEnergyCSDA;
2072 
2073  G4LossTableBuilder* bld = lManager->GetTableBuilder();
2074  theDensityFactor = bld->GetDensityFactors();
2075  theDensityIdx = bld->GetCoupleIndexes();
2076 
2077  for (size_t i=0; i<n; ++i) {
2078  G4double dedx = 0.0;
2079  pv = (*p)[i];
2080  if(pv) {
2081  dedx = pv->Value(emax, idxDEDXunRestricted);
2082  } else {
2083  pv = (*p)[(*theDensityIdx)[i]];
2084  if(pv) {
2085  dedx =
2086  pv->Value(emax, idxDEDXunRestricted)*(*theDensityFactor)[i];
2087  }
2088  }
2089  theDEDXAtMaxEnergy[i] = dedx;
2090  //G4cout << "i= " << i << " emax(MeV)= " << emax/MeV<< " dedx= "
2091  // << dedx << G4endl;
2092  }
2093  }
2094 
2095  } else if(fRestricted == tType) {
2096  /*
2097  G4cout<< "G4VEnergyLossProcess::SetDEDXTable "
2098  << particle->GetParticleName()
2099  << " oldTable " << theDEDXTable << " newTable " << p
2100  << " ion " << theIonisationTable
2101  << " IsMaster " << isMaster
2102  << " " << GetProcessName() << G4endl;
2103  G4cout << (*p) << G4endl;
2104  */
2105  theDEDXTable = p;
2106  } else if(fSubRestricted == tType) {
2107  theDEDXSubTable = p;
2108  } else if(fIsIonisation == tType) {
2109  /*
2110  G4cout<< "G4VEnergyLossProcess::SetIonisationTable "
2111  << particle->GetParticleName()
2112  << " oldTable " << theDEDXTable << " newTable " << p
2113  << " ion " << theIonisationTable
2114  << " IsMaster " << isMaster
2115  << " " << GetProcessName() << G4endl;
2116  */
2117  theIonisationTable = p;
2118  } else if(fIsSubIonisation == tType) {
2119  theIonisationSubTable = p;
2120  }
2121 }
const std::vector< G4double > * GetDensityFactors()
const char * p
Definition: xmltok.h:285
G4LossTableBuilder * GetTableBuilder()
const G4int n
G4double Value(G4double theEnergy, size_t &lastidx) const
static const G4double emax
size_t length() const
double G4double
Definition: G4Types.hh:76
const std::vector< G4int > * GetCoupleIndexes()

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetDynamicMassCharge ( G4double  massratio,
G4double  charge2ratio 
)
inline

Definition at line 652 of file G4VEnergyLossProcess.hh.

654 {
655  massRatio = massratio;
656  fFactor = charge2ratio*biasFactor*(*theDensityFactor)[currentCoupleIndex];
657  chargeSqRatio = charge2ratio;
658  reduceFactor = 1.0/(fFactor*massRatio);
659 }

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetEmModel ( G4VEmModel p,
G4int  index = 1 
)

Definition at line 409 of file G4VEnergyLossProcess.cc.

410 {
411  G4int n = emModels.size();
412  if(index >= n) { for(G4int i=n; i<=index; ++i) {emModels.push_back(0);} }
413  emModels[index] = p;
414 }
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78
const G4int n

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetFluctModel ( G4VEmFluctuationModel p)
inline

Definition at line 883 of file G4VEnergyLossProcess.hh.

884 {
885  fluctModel = p;
886 }
const char * p
Definition: xmltok.h:285

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetIntegral ( G4bool  val)
inline

Definition at line 950 of file G4VEnergyLossProcess.hh.

951 {
952  integral = val;
953  actIntegral = true;
954 }
void G4VEnergyLossProcess::SetInverseRangeTable ( G4PhysicsTable p)

Definition at line 2175 of file G4VEnergyLossProcess.cc.

2176 {
2177  theInverseRangeTable = p;
2178  if(1 < verboseLevel) {
2179  G4cout << "### Set InverseRange table " << p
2180  << " for " << particle->GetParticleName()
2181  << " and process " << GetProcessName() << G4endl;
2182  }
2183 }
G4int verboseLevel
Definition: G4VProcess.hh:368
const char * p
Definition: xmltok.h:285
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetIonisation ( G4bool  val)

Definition at line 2333 of file G4VEnergyLossProcess.cc.

2334 {
2335  isIonisation = val;
2336  if(val) { aGPILSelection = CandidateForSelection; }
2337  else { aGPILSelection = NotCandidateForSelection; }
2338 }

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetLambdaTable ( G4PhysicsTable p)

Definition at line 2187 of file G4VEnergyLossProcess.cc.

2188 {
2189  if(1 < verboseLevel) {
2190  G4cout << "### Set Lambda table " << p
2191  << " for " << particle->GetParticleName()
2192  << " and process " << GetProcessName() << G4endl;
2193  //G4cout << *p << G4endl;
2194  }
2195  theLambdaTable = p;
2196  tablesAreBuilt = true;
2197 
2198  G4LossTableBuilder* bld = lManager->GetTableBuilder();
2199  theDensityFactor = bld->GetDensityFactors();
2200  theDensityIdx = bld->GetCoupleIndexes();
2201 
2202  if(theLambdaTable) {
2203  size_t n = theLambdaTable->length();
2204  G4PhysicsVector* pv = (*theLambdaTable)[0];
2205  G4double e, ss, smax, emax;
2206 
2207  size_t i;
2208 
2209  // first loop on existing vectors
2210  for (i=0; i<n; ++i) {
2211  pv = (*theLambdaTable)[i];
2212  if(pv) {
2213  size_t nb = pv->GetVectorLength();
2214  emax = DBL_MAX;
2215  smax = 0.0;
2216  if(nb > 0) {
2217  for (size_t j=0; j<nb; ++j) {
2218  e = pv->Energy(j);
2219  ss = (*pv)(j);
2220  if(ss > smax) {
2221  smax = ss;
2222  emax = e;
2223  }
2224  }
2225  }
2226  theEnergyOfCrossSectionMax[i] = emax;
2227  theCrossSectionMax[i] = smax;
2228  if(1 < verboseLevel) {
2229  G4cout << "For " << particle->GetParticleName()
2230  << " Max CS at i= " << i << " emax(MeV)= " << emax/MeV
2231  << " lambda= " << smax << G4endl;
2232  }
2233  }
2234  }
2235  // second loop using base materials
2236  for (i=0; i<n; ++i) {
2237  pv = (*theLambdaTable)[i];
2238  if(!pv){
2239  G4int j = (*theDensityIdx)[i];
2240  theEnergyOfCrossSectionMax[i] = theEnergyOfCrossSectionMax[j];
2241  theCrossSectionMax[i] = (*theDensityFactor)[i]*theCrossSectionMax[j];
2242  }
2243  }
2244  }
2245 }
const std::vector< G4double > * GetDensityFactors()
G4int verboseLevel
Definition: G4VProcess.hh:368
const char * p
Definition: xmltok.h:285
size_t GetVectorLength() const
int G4int
Definition: G4Types.hh:78
const G4int smax
const G4String & GetParticleName() const
G4LossTableBuilder * GetTableBuilder()
G4GLOB_DLL std::ostream G4cout
const G4int n
G4double Energy(size_t index) const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
static const G4double emax
size_t length() const
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
const std::vector< G4int > * GetCoupleIndexes()

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetLinearLossLimit ( G4double  val)

Definition at line 2342 of file G4VEnergyLossProcess.cc.

2343 {
2344  if(0.0 < val && val < 1.0) {
2345  linLossLimit = val;
2346  actLinLossLimit = true;
2347  } else { PrintWarning("SetLinearLossLimit", val); }
2348 }

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetLossFluctuations ( G4bool  val)
inline

Definition at line 942 of file G4VEnergyLossProcess.hh.

943 {
944  lossFluctuationFlag = val;
945  actLossFluc = true;
946 }

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetLowestEnergyLimit ( G4double  val)

Definition at line 2369 of file G4VEnergyLossProcess.cc.

2370 {
2371  if(1.e-18 < val && val < 1.e+50) { lowestKinEnergy = val; }
2372  else { PrintWarning("SetLowestEnergyLimit", val); }
2373 }
void G4VEnergyLossProcess::SetMaxKinEnergy ( G4double  e)

Definition at line 2400 of file G4VEnergyLossProcess.cc.

2401 {
2402  if(minKinEnergy < e && e < 1.e+50) {
2403  maxKinEnergy = e;
2404  actMaxKinEnergy = true;
2405  if(e < maxKinEnergyCSDA) { maxKinEnergyCSDA = e; }
2406  } else { PrintWarning("SetMaxKinEnergy", e); }
2407 }

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetMinKinEnergy ( G4double  e)

Definition at line 2390 of file G4VEnergyLossProcess.cc.

2391 {
2392  if(1.e-18 < e && e < maxKinEnergy) {
2393  minKinEnergy = e;
2394  actMinKinEnergy = true;
2395  } else { PrintWarning("SetMinKinEnergy", e); }
2396 }

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetParticle ( const G4ParticleDefinition p)
inlineprotected

Definition at line 897 of file G4VEnergyLossProcess.hh.

898 {
899  particle = p;
900 }
const char * p
Definition: xmltok.h:285
void G4VEnergyLossProcess::SetRangeTableForLoss ( G4PhysicsTable p)

Definition at line 2151 of file G4VEnergyLossProcess.cc.

2152 {
2153  theRangeTableForLoss = p;
2154  if(1 < verboseLevel) {
2155  G4cout << "### Set Range table " << p
2156  << " for " << particle->GetParticleName()
2157  << " and process " << GetProcessName() << G4endl;
2158  }
2159 }
G4int verboseLevel
Definition: G4VProcess.hh:368
const char * p
Definition: xmltok.h:285
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetSecondaryParticle ( const G4ParticleDefinition p)
inlineprotected

Definition at line 905 of file G4VEnergyLossProcess.hh.

906 {
907  secondaryParticle = p;
908 }
const char * p
Definition: xmltok.h:285

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetSecondaryRangeTable ( G4PhysicsTable p)

Definition at line 2163 of file G4VEnergyLossProcess.cc.

2164 {
2165  theSecondaryRangeTable = p;
2166  if(1 < verboseLevel) {
2167  G4cout << "### Set SecondaryRange table " << p
2168  << " for " << particle->GetParticleName()
2169  << " and process " << GetProcessName() << G4endl;
2170  }
2171 }
G4int verboseLevel
Definition: G4VProcess.hh:368
const char * p
Definition: xmltok.h:285
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetStepFunction ( G4double  v1,
G4double  v2,
G4bool  lock = true 
)

Definition at line 2353 of file G4VEnergyLossProcess.cc.

2354 {
2355  if(actStepFunc) { return; }
2356  actStepFunc = lock;
2357  if(0.0 < v1 && 0.0 < v2 && v2 < 1.e+50) {
2358  dRoverRange = std::min(1.0, v1);
2359  finalRange = v2;
2360  } else if(v1 <= 0.0) {
2361  PrintWarning("SetStepFunction", v1);
2362  } else {
2363  PrintWarning("SetStepFunction", v2);
2364  }
2365 }
T min(const T t1, const T t2)
brief Return the smallest of the two arguments

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::SetSubLambdaTable ( G4PhysicsTable p)

Definition at line 2249 of file G4VEnergyLossProcess.cc.

2250 {
2251  theSubLambdaTable = p;
2252  if(1 < verboseLevel) {
2253  G4cout << "### Set SebLambda table " << p
2254  << " for " << particle->GetParticleName()
2255  << " and process " << GetProcessName() << G4endl;
2256  }
2257 }
G4int verboseLevel
Definition: G4VProcess.hh:368
const char * p
Definition: xmltok.h:285
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEnergyLossProcess::StartTracking ( G4Track track)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 1078 of file G4VEnergyLossProcess.cc.

1079 {
1080  /*
1081  G4cout << track->GetDefinition()->GetParticleName()
1082  << " e(MeV)= " << track->GetKineticEnergy()
1083  << " baseParticle " << baseParticle << " proc " << this;
1084  if(particle) G4cout << " " << particle->GetParticleName();
1085  G4cout << " isIon= " << isIon << " dedx " << theDEDXTable <<G4endl;
1086  */
1087  // reset parameters for the new track
1089  currentInteractionLength = mfpKinEnergy = DBL_MAX;
1090  preStepRangeEnergy = 0.0;
1091 
1092  // reset ion
1093  if(isIon) {
1094  chargeSqRatio = 0.5;
1095 
1096  G4double newmass = track->GetDefinition()->GetPDGMass();
1097  if(baseParticle) {
1098  massRatio = baseParticle->GetPDGMass()/newmass;
1099  } else if(theGenericIon) {
1100  massRatio = proton_mass_c2/newmass;
1101  } else {
1102  massRatio = 1.0;
1103  }
1104  }
1105  // forced biasing only for primary particles
1106  if(biasManager) {
1107  if(0 == track->GetParentID()) {
1108  // primary particle
1109  biasFlag = true;
1110  biasManager->ResetForcedInteraction();
1111  }
1112  }
1113 }
G4ParticleDefinition * GetDefinition() const
G4int GetParentID() const
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
G4double currentInteractionLength
Definition: G4VProcess.hh:297
float proton_mass_c2
Definition: hepunit.py:275
G4double GetPDGMass() const
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

G4bool G4VEnergyLossProcess::StorePhysicsTable ( const G4ParticleDefinition part,
const G4String directory,
G4bool  ascii = false 
)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 1754 of file G4VEnergyLossProcess.cc.

1757 {
1758  G4bool res = true;
1759  //G4cout << "G4VEnergyLossProcess::StorePhysicsTable: " << part->GetParticleName()
1760  // << " " << directory << " " << ascii << G4endl;
1761  if (!isMaster || baseParticle || part != particle ) return res;
1762 
1763  if(!StoreTable(part,theDEDXTable,ascii,directory,"DEDX"))
1764  {res = false;}
1765 
1766  if(!StoreTable(part,theDEDXunRestrictedTable,ascii,directory,"DEDXnr"))
1767  {res = false;}
1768 
1769  if(!StoreTable(part,theDEDXSubTable,ascii,directory,"SubDEDX"))
1770  {res = false;}
1771 
1772  if(!StoreTable(part,theIonisationTable,ascii,directory,"Ionisation"))
1773  {res = false;}
1774 
1775  if(!StoreTable(part,theIonisationSubTable,ascii,directory,"SubIonisation"))
1776  {res = false;}
1777 
1778  if(isIonisation &&
1779  !StoreTable(part,theCSDARangeTable,ascii,directory,"CSDARange"))
1780  {res = false;}
1781 
1782  if(isIonisation &&
1783  !StoreTable(part,theRangeTableForLoss,ascii,directory,"Range"))
1784  {res = false;}
1785 
1786  if(isIonisation &&
1787  !StoreTable(part,theInverseRangeTable,ascii,directory,"InverseRange"))
1788  {res = false;}
1789 
1790  if(!StoreTable(part,theLambdaTable,ascii,directory,"Lambda"))
1791  {res = false;}
1792 
1793  if(!StoreTable(part,theSubLambdaTable,ascii,directory,"SubLambda"))
1794  {res = false;}
1795 
1796  if ( !res ) {
1797  if(1 < verboseLevel) {
1798  G4cout << "Physics tables are stored for "
1799  << particle->GetParticleName()
1800  << " and process " << GetProcessName()
1801  << " in the directory <" << directory
1802  << "> " << G4endl;
1803  }
1804  } else {
1805  G4cout << "Fail to store Physics Tables for "
1806  << particle->GetParticleName()
1807  << " and process " << GetProcessName()
1808  << " in the directory <" << directory
1809  << "> " << G4endl;
1810  }
1811  return res;
1812 }
G4int verboseLevel
Definition: G4VProcess.hh:368
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

G4PhysicsTable * G4VEnergyLossProcess::SubLambdaTable ( ) const
inline

Definition at line 1077 of file G4VEnergyLossProcess.hh.

1078 {
1079  return theSubLambdaTable;
1080 }

Here is the caller graph for this function:

G4bool G4VEnergyLossProcess::TablesAreBuilt ( ) const
inline

Definition at line 1000 of file G4VEnergyLossProcess.hh.

1001 {
1002  return tablesAreBuilt;
1003 }
void G4VEnergyLossProcess::UpdateEmModel ( const G4String nam,
G4double  emin,
G4double  emax 
)

Definition at line 401 of file G4VEnergyLossProcess.cc.

403 {
404  modelManager->UpdateEmModel(nam, emin, emax);
405 }
void UpdateEmModel(const G4String &model_name, G4double emin, G4double emax)
static const G4double emax

Here is the call graph for this function:

Member Data Documentation

G4ParticleChangeForLoss G4VEnergyLossProcess::fParticleChange
protected

Definition at line 572 of file G4VEnergyLossProcess.hh.


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