Geant4  10.02.p02
G4VEnergyLossProcess.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4VEnergyLossProcess.cc 93264 2015-10-14 09:30:04Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4VEnergyLossProcess
34 //
35 // Author: Vladimir Ivanchenko
36 //
37 // Creation date: 03.01.2002
38 //
39 // Modifications:
40 //
41 // 13-11-02 Minor fix - use normalised direction (V.Ivanchenko)
42 // 04-12-02 Minor change in PostStepDoIt (V.Ivanchenko)
43 // 23-12-02 Change interface in order to move to cut per region (V.Ivanchenko)
44 // 26-12-02 Secondary production moved to derived classes (V.Ivanchenko)
45 // 04-01-03 Fix problem of very small steps for ions (V.Ivanchenko)
46 // 20-01-03 Migrade to cut per region (V.Ivanchenko)
47 // 24-01-03 Temporarily close a control on usage of couples (V.Ivanchenko)
48 // 24-01-03 Make models region aware (V.Ivanchenko)
49 // 05-02-03 Fix compilation warnings (V.Ivanchenko)
50 // 06-02-03 Add control on tmax in PostStepDoIt (V.Ivanchenko)
51 // 13-02-03 SubCutoffProcessors defined for regions (V.Ivanchenko)
52 // 15-02-03 Lambda table can be scaled (V.Ivanchenko)
53 // 17-02-03 Fix problem of store/restore tables (V.Ivanchenko)
54 // 18-02-03 Add control on CutCouple usage (V.Ivanchenko)
55 // 26-02-03 Simplify control on GenericIons (V.Ivanchenko)
56 // 06-03-03 Control on GenericIons using SubType+ update verbose (V.Ivanchenko)
57 // 10-03-03 Add Ion registration (V.Ivanchenko)
58 // 22-03-03 Add Initialisation of cash (V.Ivanchenko)
59 // 26-03-03 Remove finalRange modification (V.Ivanchenko)
60 // 09-04-03 Fix problem of negative range limit for non integral (V.Ivanchenko)
61 // 26-04-03 Fix retrieve tables (V.Ivanchenko)
62 // 06-05-03 Set defalt finalRange = 1 mm (V.Ivanchenko)
63 // 12-05-03 Update range calculations + lowKinEnergy (V.Ivanchenko)
64 // 13-05-03 Add calculation of precise range (V.Ivanchenko)
65 // 23-05-03 Remove tracking cuts (V.Ivanchenko)
66 // 03-06-03 Fix initialisation problem for STD ionisation (V.Ivanchenko)
67 // 21-07-03 Add UpdateEmModel method (V.Ivanchenko)
68 // 03-11-03 Fix initialisation problem in RetrievePhysicsTable (V.Ivanchenko)
69 // 04-11-03 Add checks in RetrievePhysicsTable (V.Ivanchenko)
70 // 12-11-03 G4EnergyLossSTD -> G4EnergyLossProcess (V.Ivanchenko)
71 // 21-01-04 Migrade to G4ParticleChangeForLoss (V.Ivanchenko)
72 // 27-02-04 Fix problem of loss in low presure gases, cleanup precise range
73 // calculation, use functions ForLoss in AlongStepDoIt (V.Ivanchenko)
74 // 10-03-04 Fix a problem of Precise Range table (V.Ivanchenko)
75 // 19-03-04 Fix a problem energy below lowestKinEnergy (V.Ivanchenko)
76 // 31-03-04 Fix a problem of retrieve tables (V.Ivanchenko)
77 // 21-07-04 Check weather AtRest are active or not (V.Ivanchenko)
78 // 03-08-04 Add pointer of DEDX table to all processes (V.Ivanchenko)
79 // 06-08-04 Clear up names of member functions (V.Ivanchenko)
80 // 06-08-04 Clear up names of member functions (V.Ivanchenko)
81 // 27-08-04 Add NeedBuildTables method (V.Ivanchneko)
82 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivantchenko)
83 // 11-03-05 Shift verbose level by 1 (V.Ivantchenko)
84 // 08-04-05 Major optimisation of internal interfaces (V.Ivantchenko)
85 // 11-04-05 Use MaxSecondaryEnergy from a model (V.Ivanchenko)
86 // 25-07-05 Add extra protection PostStep for non-integral mode (V.Ivanchenko)
87 // 12-08-05 Integral=false; SetStepFunction(0.2, 0.1*mm) (mma)
88 // 18-08-05 Return back both AlongStep and PostStep from 7.0 (V.Ivanchenko)
89 // 02-09-05 Default StepFunction 0.2 1 mm + integral (V.Ivanchenko)
90 // 04-09-05 default lambdaFactor 0.8 (V.Ivanchenko)
91 // 05-10-05 protection against 0 energy loss added (L.Urban)
92 // 17-10-05 protection above has been removed (L.Urban)
93 // 06-01-06 reset currentCouple when StepFunction is changed (V.Ivanchenko)
94 // 10-01-06 PreciseRange -> CSDARange (V.Ivantchenko)
95 // 18-01-06 Clean up subcutoff including recalculation of presafety (VI)
96 // 20-01-06 Introduce G4EmTableType and reducing number of methods (VI)
97 // 22-03-06 Add control on warning printout AlongStep (VI)
98 // 23-03-06 Use isIonisation flag (V.Ivanchenko)
99 // 07-06-06 Do not reflect AlongStep in subcutoff regime (V.Ivanchenko)
100 // 14-01-07 add SetEmModel(index) and SetFluctModel() (mma)
101 // 16-01-07 add IonisationTable and IonisationSubTable (V.Ivanchenko)
102 // 16-02-07 set linLossLimit=1.e-6 (V.Ivanchenko)
103 // 13-03-07 use SafetyHelper instead of navigator (V.Ivanchenko)
104 // 10-04-07 use unique SafetyHelper (V.Ivanchenko)
105 // 12-04-07 Add verbosity at destruction (V.Ivanchenko)
106 // 25-04-07 move initialisation of safety helper to BuildPhysicsTable (VI)
107 // 27-10-07 Virtual functions moved to source (V.Ivanchenko)
108 // 24-06-09 Removed hidden bin in G4PhysicsVector (V.Ivanchenko)
109 // 01-25-09 (Xin Dong) Phase II change for Geant4 multi-threading:
110 // New methods SlavePreparePhysicsTable, SlaveBuildPhysicsTable
111 // Worker threads share physics tables with the master thread for
112 // this kind of process. This member function is used by worker
113 // threads to achieve the partial effect of the master thread when
114 // it builds physcis tables.
115 // 15-10-10 Fixed 4-momentum balance if deexcitation is active (L.Pandola)
116 // 30-05-12 Call new ApplySecondaryBiasing so 2ries may be unique (D. Sawkey)
117 // 30-05-12 Fix bug in forced biasing: now called on first step (D. Sawkey)
118 // 04-06-13 Adoptation to MT mode, adding internal cache to GetRangeForLoss,
119 // more accurate initialisation for ions (V.Ivanchenko)
120 //
121 // Class Description:
122 //
123 // It is the unified energy loss process it calculates the continuous
124 // energy loss for charged particles using a set of Energy Loss
125 // models valid for different energy regions. There are a possibility
126 // to create and access to dE/dx and range tables, or to calculate
127 // that information on fly.
128 // -------------------------------------------------------------------
129 //
130 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
131 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
132 
133 #include "G4VEnergyLossProcess.hh"
134 #include "G4PhysicalConstants.hh"
135 #include "G4SystemOfUnits.hh"
136 #include "G4ProcessManager.hh"
137 #include "G4LossTableManager.hh"
138 #include "G4LossTableBuilder.hh"
139 #include "G4Step.hh"
140 #include "G4ParticleDefinition.hh"
141 #include "G4ParticleTable.hh"
142 #include "G4VEmModel.hh"
143 #include "G4VEmFluctuationModel.hh"
144 #include "G4DataVector.hh"
145 #include "G4PhysicsLogVector.hh"
146 #include "G4VParticleChange.hh"
147 #include "G4Gamma.hh"
148 #include "G4Electron.hh"
149 #include "G4Positron.hh"
150 #include "G4ProcessManager.hh"
151 #include "G4UnitsTable.hh"
152 #include "G4ProductionCutsTable.hh"
153 #include "G4Region.hh"
154 #include "G4RegionStore.hh"
155 #include "G4PhysicsTableHelper.hh"
156 #include "G4SafetyHelper.hh"
158 #include "G4EmConfigurator.hh"
159 #include "G4VAtomDeexcitation.hh"
160 #include "G4VSubCutProducer.hh"
161 #include "G4EmBiasingManager.hh"
162 #include "G4Log.hh"
163 #include <iostream>
164 
165 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
166 
168  G4ProcessType type):
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  nWarnings(0),
199  mfpKinEnergy(0.0)
200 {
202  SetVerboseLevel(1);
203 
204  // low energy limit
206  preStepKinEnergy = 0.0;
207  preStepRangeEnergy = 0.0;
209 
210  // Size of tables assuming spline
211  minKinEnergy = 0.1*keV;
212  maxKinEnergy = 10.0*TeV;
213  nBins = 77;
214  maxKinEnergyCSDA = 1.0*GeV;
215  nBinsCSDA = 35;
217  = actLossFluc = false;
218 
219  // default linear loss limit for spline
220  linLossLimit = 0.01;
221 
222  // default dRoverRange and finalRange
223  SetStepFunction(0.2, 1.0*mm);
224 
225  // default lambda factor
226  lambdaFactor = 0.8;
227 
228  // cross section biasing
229  biasFactor = 1.0;
230 
231  // particle types
235  theGenericIon = nullptr;
236 
237  // run time objects
242  ->GetSafetyHelper();
244 
245  // initialise model
247  lManager->Register(this);
248  fluctModel = nullptr;
249  currentModel = nullptr;
250  atomDeexcitation = nullptr;
251  subcutProducer = nullptr;
252 
253  biasManager = nullptr;
254  biasFlag = false;
255  weightFlag = false;
256  isMaster = true;
257  lastIdx = 0;
258 
262 
263  scTracks.reserve(5);
264  secParticles.reserve(5);
265 
266  theCuts = theSubCuts = nullptr;
267  currentMaterial = nullptr;
271 
272  secID = biasID = subsecID = -1;
273 }
274 
275 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
276 
278 {
279  /*
280  G4cout << "** G4VEnergyLossProcess::~G4VEnergyLossProcess() for "
281  << GetProcessName() << " isMaster: " << isMaster
282  << " basePart: " << baseParticle
283  << G4endl;
284  */
285  Clean();
286 
287  // G4cout << " isIonisation " << isIonisation << " "
288  // << theDEDXTable << " " << theIonisationTable << G4endl;
289 
290  if (isMaster && !baseParticle) {
291  if(theDEDXTable) {
292 
293  //G4cout << " theIonisationTable " << theIonisationTable << G4endl;
295  //G4cout << " delete theDEDXTable " << theDEDXTable << G4endl;
297  delete theDEDXTable;
298  theDEDXTable = nullptr;
299  if(theDEDXSubTable) {
301  { theIonisationSubTable = nullptr; }
303  delete theDEDXSubTable;
304  theDEDXSubTable = nullptr;
305  }
306  }
307  //G4cout << " theIonisationTable " << theIonisationTable << G4endl;
308  if(theIonisationTable) {
309  //G4cout << " delete theIonisationTable " << theIonisationTable << G4endl;
311  delete theIonisationTable;
312  theIonisationTable = nullptr;
313  }
316  delete theIonisationSubTable;
317  theIonisationSubTable = nullptr;
318  }
322  theDEDXunRestrictedTable = nullptr;
323  }
326  delete theCSDARangeTable;
327  theCSDARangeTable = nullptr;
328  }
329  //G4cout << "delete RangeTable: " << theRangeTableForLoss << G4endl;
332  delete theRangeTableForLoss;
333  theRangeTableForLoss = nullptr;
334  }
335  //G4cout << "delete InvRangeTable: " << theInverseRangeTable << G4endl;
336  if(theInverseRangeTable && isIonisation /*&& !isIon*/) {
338  delete theInverseRangeTable;
339  theInverseRangeTable = nullptr;
340  }
341  //G4cout << "delete LambdaTable: " << theLambdaTable << G4endl;
342  if(theLambdaTable) {
344  delete theLambdaTable;
345  theLambdaTable = nullptr;
346  }
347  if(theSubLambdaTable) {
349  delete theSubLambdaTable;
350  theSubLambdaTable = nullptr;
351  }
352  }
353 
354  delete modelManager;
355  delete biasManager;
356  lManager->DeRegister(this);
357  //G4cout << "** all removed" << G4endl;
358 }
359 
360 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
361 
363 {
364  /*
365  if(1 < verboseLevel) {
366  G4cout << "G4VEnergyLossProcess::Clear() for " << GetProcessName()
367  << G4endl;
368  }
369  */
370  delete [] idxSCoffRegions;
371 
372  tablesAreBuilt = false;
373 
374  scProcesses.clear();
375  nProcesses = 0;
376 
380 }
381 
382 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
383 
385  const G4Material*,
386  G4double cut)
387 {
388  return cut;
389 }
390 
391 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
392 
394  G4VEmFluctuationModel* fluc,
395  const G4Region* region)
396 {
397  modelManager->AddEmModel(order, p, fluc, region);
398  if(p) { p->SetParticleChange(pParticleChange, fluc); }
399 }
400 
401 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
402 
404  G4double emin, G4double emax)
405 {
406  modelManager->UpdateEmModel(nam, emin, emax);
407 }
408 
409 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
410 
412 {
413  G4int n = emModels.size();
414  if(index >= n) { for(G4int i=n; i<=index; ++i) {emModels.push_back(0);} }
415  emModels[index] = p;
416 }
417 
418 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
419 
421 {
422  G4VEmModel* p = nullptr;
423  if(index >= 0 && index < G4int(emModels.size())) { p = emModels[index]; }
424  return p;
425 }
426 
427 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
428 
430 {
431  return modelManager->GetModel(idx, ver);
432 }
433 
434 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
435 
437 {
438  return modelManager->NumberOfModels();
439 }
440 
441 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
442 
443 void
445 {
446  if(1 < verboseLevel) {
447  G4cout << "G4VEnergyLossProcess::PreparePhysicsTable for "
448  << GetProcessName() << " for " << part.GetParticleName()
449  << " " << this << G4endl;
450  }
451 
452  const G4VEnergyLossProcess* masterProcess =
453  static_cast<const G4VEnergyLossProcess*>(GetMasterProcess());
454  if(masterProcess && masterProcess != this) { isMaster = false; }
455 
456  currentCouple = nullptr;
457  preStepLambda = 0.0;
459  fRange = DBL_MAX;
460  preStepKinEnergy = 0.0;
461  preStepRangeEnergy = 0.0;
462  chargeSqRatio = 1.0;
463  massRatio = 1.0;
464  reduceFactor = 1.0;
465  fFactor = 1.0;
466  lastIdx = 0;
467 
468  // Are particle defined?
469  if( !particle ) { particle = &part; }
470 
471  if(part.GetParticleType() == "nucleus") {
472 
473  G4String pname = part.GetParticleName();
474  if(pname != "deuteron" && pname != "triton" &&
475  pname != "alpha+" && pname != "helium" &&
476  pname != "hydrogen") {
477 
478  if(!theGenericIon) {
479  theGenericIon =
481  }
482  isIon = true;
486  size_t n = v->size();
487  for(size_t j=0; j<n; ++j) {
488  if((*v)[j] == this) {
490  break;
491  }
492  }
493  }
494  }
495  }
496 
497  if( particle != &part ) {
498  if(!isIon) {
499  lManager->RegisterExtraParticle(&part, this);
500  }
501  if(1 < verboseLevel) {
502  G4cout << "### G4VEnergyLossProcess::PreparePhysicsTable()"
503  << " interrupted for "
504  << part.GetParticleName() << " isIon= " << isIon
505  << " particle " << particle << " GenericIon " << theGenericIon
506  << G4endl;
507  }
508  return;
509  }
510 
511  Clean();
512  lManager->PreparePhysicsTable(&part, this, isMaster);
514 
515  // Base particle and set of models can be defined here
517 
518  const G4ProductionCutsTable* theCoupleTable=
520  size_t n = theCoupleTable->GetTableSize();
521 
522  theDEDXAtMaxEnergy.resize(n, 0.0);
523  theRangeAtMaxEnergy.resize(n, 0.0);
524  theEnergyOfCrossSectionMax.resize(n, 0.0);
525  theCrossSectionMax.resize(n, DBL_MAX);
526 
527  // parameters of the process
535  *G4lrint(std::log10(maxKinEnergyCSDA/minKinEnergy));
540 
541  G4double initialCharge = particle->GetPDGCharge();
542  G4double initialMass = particle->GetPDGMass();
543 
544  if (baseParticle) {
545  massRatio = (baseParticle->GetPDGMass())/initialMass;
546  G4double q = initialCharge/baseParticle->GetPDGCharge();
547  chargeSqRatio = q*q;
548  if(chargeSqRatio > 0.0) { reduceFactor = 1.0/(chargeSqRatio*massRatio); }
549  }
550  if(initialMass < MeV) {
552  } else {
554  }
555 
556  // Tables preparation
557  if (isMaster && !baseParticle) {
558 
559  if(theDEDXTable && isIonisation) {
562  delete theDEDXTable;
564  }
568  delete theDEDXSubTable;
570  }
571  }
572 
575 
576  if(theDEDXSubTable) {
577  theDEDXSubTable =
579  }
580 
581  if (lManager->BuildCSDARange()) {
586  }
587 
589 
590  if(isIonisation) {
595  }
596 
598  theDEDXSubTable =
602  }
603  }
604  /*
605  G4cout << "** G4VEnergyLossProcess::PreparePhysicsTable() for "
606  << GetProcessName() << " and " << particle->GetParticleName()
607  << " isMaster: " << isMaster << " isIonisation: "
608  << isIonisation << G4endl;
609  G4cout << " theDEDX: " << theDEDXTable
610  << " theRange: " << theRangeTableForLoss
611  << " theInverse: " << theInverseRangeTable
612  << " theLambda: " << theLambdaTable << G4endl;
613  */
614  // forced biasing
615  if(biasManager) {
617  biasFlag = false;
618  }
619 
620  // defined ID of secondary particles
621  if(isMaster) {
622  G4String nam1 = GetProcessName();
623  G4String nam4 = nam1 + "_split";
624  G4String nam5 = nam1 + "_subcut";
628  }
629 
630  // initialisation of models
632  for(G4int i=0; i<nmod; ++i) {
633  G4VEmModel* mod = modelManager->GetModel(i);
637  if(mod->HighEnergyLimit() > maxKinEnergy) {
639  }
640  }
643  verboseLevel);
644 
645  // Sub Cutoff
646  if(nSCoffRegions > 0) {
647  if(theParameters->MinSubRange() < 1.0) { useSubCutoff = true; }
648 
650 
651  idxSCoffRegions = new G4bool[n];
652  for (size_t j=0; j<n; ++j) {
653 
654  const G4MaterialCutsCouple* couple =
655  theCoupleTable->GetMaterialCutsCouple(j);
656  const G4ProductionCuts* pcuts = couple->GetProductionCuts();
657 
658  G4bool reg = false;
659  for(G4int i=0; i<nSCoffRegions; ++i) {
660  if( pcuts == scoffRegions[i]->GetProductionCuts()) {
661  reg = true;
662  break;
663  }
664  }
665  idxSCoffRegions[j] = reg;
666  }
667  }
668 
669  if(1 < verboseLevel) {
670  G4cout << "G4VEnergyLossProcess::PrepearPhysicsTable() is done "
671  << " for local " << particle->GetParticleName()
672  << " isIon= " << isIon;
673  if(baseParticle) {
674  G4cout << "; base: " << baseParticle->GetParticleName();
675  }
676  G4cout << " chargeSqRatio= " << chargeSqRatio
677  << " massRatio= " << massRatio
678  << " reduceFactor= " << reduceFactor << G4endl;
679  if (nSCoffRegions) {
680  G4cout << " SubCutoff Regime is ON for regions: " << G4endl;
681  for (G4int i=0; i<nSCoffRegions; ++i) {
682  const G4Region* r = scoffRegions[i];
683  G4cout << " " << r->GetName() << G4endl;
684  }
685  }
686  }
687 }
688 
689 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
690 
692 {
693  if(1 < verboseLevel) {
694 
695  G4cout << "### G4VEnergyLossProcess::BuildPhysicsTable() for "
696  << GetProcessName()
697  << " and particle " << part.GetParticleName()
698  << "; local: " << particle->GetParticleName();
699  if(baseParticle) {
700  G4cout << "; base: " << baseParticle->GetParticleName();
701  }
702  G4cout << " TablesAreBuilt= " << tablesAreBuilt
703  << " isIon= " << isIon << " " << this << G4endl;
704  }
705 
706  if(&part == particle) {
707 
709  if(isMaster) {
713 
714  } else {
715 
716  const G4VEnergyLossProcess* masterProcess =
717  static_cast<const G4VEnergyLossProcess*>(GetMasterProcess());
718 
719  // define density factors for worker thread
720  bld->InitialiseBaseMaterials(masterProcess->DEDXTable());
723 
724  // copy table pointers from master thread
725  SetDEDXTable(masterProcess->DEDXTable(),fRestricted);
727  SetDEDXTable(masterProcess->DEDXunRestrictedTable(),fTotal);
728  SetDEDXTable(masterProcess->IonisationTable(),fIsIonisation);
730  SetRangeTableForLoss(masterProcess->RangeTableForLoss());
731  SetCSDARangeTable(masterProcess->CSDARangeTable());
733  SetInverseRangeTable(masterProcess->InverseRangeTable());
734  SetLambdaTable(masterProcess->LambdaTable());
735  SetSubLambdaTable(masterProcess->SubLambdaTable());
736  isIonisation = masterProcess->IsIonisationProcess();
737 
738  tablesAreBuilt = true;
739  // local initialisation of models
740  G4bool printing = true;
741  G4int numberOfModels = modelManager->NumberOfModels();
742  for(G4int i=0; i<numberOfModels; ++i) {
743  G4VEmModel* mod = GetModelByIndex(i, printing);
744  G4VEmModel* mod0= masterProcess->GetModelByIndex(i,printing);
745  mod->InitialiseLocal(particle, mod0);
746  }
747 
749  }
750 
751  // needs to be done only once
753  }
754  // explicitly defined printout by particle name
755  G4String num = part.GetParticleName();
756  if(1 < verboseLevel ||
757  (0 < verboseLevel && (num == "e-" ||
758  num == "e+" || num == "mu+" ||
759  num == "mu-" || num == "proton"||
760  num == "pi+" || num == "pi-" ||
761  num == "kaon+" || num == "kaon-" ||
762  num == "alpha" || num == "anti_proton" ||
763  num == "GenericIon")))
764  {
765  PrintInfoDefinition(part);
766  }
767 
768  // Added tracking cut to avoid tracking artifacts
769  // identify deexcitation flag
770  if(isIonisation) {
771  //VI: seems not needed fParticleChange.SetLowEnergyLimit(lowestKinEnergy);
774  if(atomDeexcitation) {
775  if(atomDeexcitation->IsPIXEActive()) { useDeexcitation = true; }
776  }
777  }
778  /*
779  G4cout << "** G4VEnergyLossProcess::BuildPhysicsTable() for "
780  << GetProcessName() << " and " << particle->GetParticleName()
781  << " isMaster: " << isMaster << " isIonisation: "
782  << isIonisation << G4endl;
783  G4cout << " theDEDX: " << theDEDXTable
784  << " theRange: " << theRangeTableForLoss
785  << " theInverse: " << theInverseRangeTable
786  << " theLambda: " << theLambdaTable << G4endl;
787  */
788  //if(1 < verboseLevel || verb) {
789  if(1 < verboseLevel) {
790  G4cout << "### G4VEnergyLossProcess::BuildPhysicsTable() done for "
791  << GetProcessName()
792  << " and particle " << part.GetParticleName();
793  if(isIonisation) { G4cout << " isIonisation flag = 1"; }
794  G4cout << G4endl;
795  }
796 }
797 
798 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
799 
801 {
802  if(1 < verboseLevel ) {
803  G4cout << "G4VEnergyLossProcess::BuildDEDXTable() of type " << tType
804  << " for " << GetProcessName()
805  << " and particle " << particle->GetParticleName()
806  << G4endl;
807  }
808  G4PhysicsTable* table = nullptr;
810  G4int bin = nBins;
811 
812  if(fTotal == tType) {
813  emax = maxKinEnergyCSDA;
814  bin = nBinsCSDA;
815  table = theDEDXunRestrictedTable;
816  } else if(fRestricted == tType) {
817  table = theDEDXTable;
818  } else if(fSubRestricted == tType) {
819  table = theDEDXSubTable;
820  } else {
821  G4cout << "G4VEnergyLossProcess::BuildDEDXTable WARNING: wrong type "
822  << tType << G4endl;
823  }
824 
825  // Access to materials
826  const G4ProductionCutsTable* theCoupleTable=
828  size_t numOfCouples = theCoupleTable->GetTableSize();
829 
830  if(1 < verboseLevel) {
831  G4cout << numOfCouples << " materials"
832  << " minKinEnergy= " << minKinEnergy
833  << " maxKinEnergy= " << emax
834  << " nbin= " << bin
835  << " EmTableType= " << tType
836  << " table= " << table << " " << this
837  << G4endl;
838  }
839  if(!table) { return table; }
840 
842  G4bool splineFlag = theParameters->Spline();
843  G4PhysicsLogVector* aVector = nullptr;
844  G4PhysicsLogVector* bVector = nullptr;
845 
846  for(size_t i=0; i<numOfCouples; ++i) {
847 
848  if(1 < verboseLevel) {
849  G4cout << "G4VEnergyLossProcess::BuildDEDXVector Idx= " << i
850  << " flagTable= " << table->GetFlag(i)
851  << " Flag= " << bld->GetFlag(i) << G4endl;
852  }
853  if(bld->GetFlag(i)) {
854 
855  // create physics vector and fill it
856  const G4MaterialCutsCouple* couple =
857  theCoupleTable->GetMaterialCutsCouple(i);
858  if((*table)[i]) { delete (*table)[i]; }
859  if(bVector) {
860  aVector = new G4PhysicsLogVector(*bVector);
861  } else {
862  bVector = new G4PhysicsLogVector(minKinEnergy, emax, bin);
863  aVector = bVector;
864  }
865  aVector->SetSpline(splineFlag);
866 
867  modelManager->FillDEDXVector(aVector, couple, tType);
868  if(splineFlag) { aVector->FillSecondDerivatives(); }
869 
870  // Insert vector for this material into the table
871  G4PhysicsTableHelper::SetPhysicsVector(table, i, aVector);
872  }
873  }
874 
875  if(1 < verboseLevel) {
876  G4cout << "G4VEnergyLossProcess::BuildDEDXTable(): table is built for "
878  << " and process " << GetProcessName()
879  << G4endl;
880  if(2 < verboseLevel) G4cout << (*table) << G4endl;
881  }
882 
883  return table;
884 }
885 
886 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
887 
889 {
890  G4PhysicsTable* table = nullptr;
891 
892  if(fRestricted == tType) {
893  table = theLambdaTable;
894  } else if(fSubRestricted == tType) {
895  table = theSubLambdaTable;
896  } else {
897  G4cout << "G4VEnergyLossProcess::BuildLambdaTable WARNING: wrong type "
898  << tType << G4endl;
899  }
900 
901  if(1 < verboseLevel) {
902  G4cout << "G4VEnergyLossProcess::BuildLambdaTable() of type "
903  << tType << " for process "
904  << GetProcessName() << " and particle "
906  << " EmTableType= " << tType
907  << " table= " << table
908  << G4endl;
909  }
910  if(!table) {return table;}
911 
912  // Access to materials
913  const G4ProductionCutsTable* theCoupleTable=
915  size_t numOfCouples = theCoupleTable->GetTableSize();
916 
920 
921  G4bool splineFlag = theParameters->Spline();
922  G4PhysicsLogVector* aVector = nullptr;
924 
925  for(size_t i=0; i<numOfCouples; ++i) {
926 
927  if (bld->GetFlag(i)) {
928 
929  // create physics vector and fill it
930  const G4MaterialCutsCouple* couple =
931  theCoupleTable->GetMaterialCutsCouple(i);
932  delete (*table)[i];
933 
934  G4bool startNull = true;
935  G4double emin =
936  MinPrimaryEnergy(particle,couple->GetMaterial(),(*theCuts)[i]);
937  if(minKinEnergy > emin) {
938  emin = minKinEnergy;
939  startNull = false;
940  }
941 
943  if(emax <= emin) { emax = 2*emin; }
944  G4int bin = G4lrint(nBins*G4Log(emax/emin)/scale);
945  bin = std::max(bin, 3);
946  aVector = new G4PhysicsLogVector(emin, emax, bin);
947  aVector->SetSpline(splineFlag);
948 
949  modelManager->FillLambdaVector(aVector, couple, startNull, tType);
950  if(splineFlag) { aVector->FillSecondDerivatives(); }
951 
952  // Insert vector for this material into the table
953  G4PhysicsTableHelper::SetPhysicsVector(table, i, aVector);
954  }
955  }
956 
957  if(1 < verboseLevel) {
958  G4cout << "Lambda table is built for "
960  << G4endl;
961  }
962 
963  return table;
964 }
965 
966 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
967 
968 void
970 {
971  if(0 < verboseLevel) {
972  G4cout << std::setprecision(6);
973  G4cout << G4endl << GetProcessName() << ": for "
974  << part.GetParticleName()
975  << " SubType= " << GetProcessSubType()
976  << G4endl;
977  G4cout << " dE/dx and range tables from "
978  << G4BestUnit(minKinEnergy,"Energy")
979  << " to " << G4BestUnit(maxKinEnergy,"Energy")
980  << " in " << nBins << " bins" << G4endl
981  << " Lambda tables from threshold to "
982  << G4BestUnit(maxKinEnergy,"Energy")
983  << ", " << theParameters->NumberOfBinsPerDecade()
984  << " bins per decade, spline: "
985  << theParameters->Spline()
986  << G4endl;
988  G4cout << " finalRange(mm)= " << finalRange/mm
989  << ", dRoverRange= " << dRoverRange
990  << ", integral: " << integral
991  << ", fluct: " << lossFluctuationFlag
992  << ", linLossLimit= " << linLossLimit
993  << G4endl;
994  }
995  PrintInfo();
998  G4cout << " CSDA range table up"
999  << " to " << G4BestUnit(maxKinEnergyCSDA,"Energy")
1000  << " in " << nBinsCSDA << " bins" << G4endl;
1001  }
1002  if(nSCoffRegions>0 && isIonisation) {
1003  G4cout << " Subcutoff sampling in " << nSCoffRegions
1004  << " regions" << G4endl;
1005  }
1006  if(2 < verboseLevel) {
1007  G4cout << " DEDXTable address= " << theDEDXTable << G4endl;
1008  if(theDEDXTable && isIonisation) G4cout << (*theDEDXTable) << G4endl;
1009  G4cout << "non restricted DEDXTable address= "
1012  G4cout << (*theDEDXunRestrictedTable) << G4endl;
1013  }
1014  if(theDEDXSubTable && isIonisation) {
1015  G4cout << (*theDEDXSubTable) << G4endl;
1016  }
1017  G4cout << " CSDARangeTable address= " << theCSDARangeTable
1018  << G4endl;
1020  G4cout << (*theCSDARangeTable) << G4endl;
1021  }
1022  G4cout << " RangeTableForLoss address= " << theRangeTableForLoss
1023  << G4endl;
1025  G4cout << (*theRangeTableForLoss) << G4endl;
1026  }
1027  G4cout << " InverseRangeTable address= " << theInverseRangeTable
1028  << G4endl;
1030  G4cout << (*theInverseRangeTable) << G4endl;
1031  }
1032  G4cout << " LambdaTable address= " << theLambdaTable << G4endl;
1033  if(theLambdaTable && isIonisation) {
1034  G4cout << (*theLambdaTable) << G4endl;
1035  }
1036  G4cout << " SubLambdaTable address= " << theSubLambdaTable
1037  << G4endl;
1039  G4cout << (*theSubLambdaTable) << G4endl;
1040  }
1041  }
1042  }
1043 }
1044 
1045 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1046 
1048 {
1049  G4RegionStore* regionStore = G4RegionStore::GetInstance();
1050  const G4Region* reg = r;
1051  if (!reg) {
1052  reg = regionStore->GetRegion("DefaultRegionForTheWorld", false);
1053  }
1054 
1055  // the region is in the list
1056  if (nSCoffRegions > 0) {
1057  for (G4int i=0; i<nSCoffRegions; ++i) {
1058  if (reg == scoffRegions[i]) {
1059  return;
1060  }
1061  }
1062  }
1063 
1064  // new region
1065  if(val) {
1066  scoffRegions.push_back(reg);
1067  ++nSCoffRegions;
1068  }
1069 }
1070 
1071 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1072 
1074 {
1075  /*
1076  G4cout << track->GetDefinition()->GetParticleName()
1077  << " e(MeV)= " << track->GetKineticEnergy()
1078  << " baseParticle " << baseParticle << " proc " << this;
1079  if(particle) G4cout << " " << particle->GetParticleName();
1080  G4cout << " isIon= " << isIon << " dedx " << theDEDXTable <<G4endl;
1081  */
1082  // reset parameters for the new track
1085  preStepRangeEnergy = 0.0;
1086 
1087  // reset ion
1088  if(isIon) {
1089  chargeSqRatio = 0.5;
1090 
1091  G4double newmass = track->GetDefinition()->GetPDGMass();
1092  if(baseParticle) {
1093  massRatio = baseParticle->GetPDGMass()/newmass;
1094  } else if(theGenericIon) {
1095  massRatio = proton_mass_c2/newmass;
1096  } else {
1097  massRatio = 1.0;
1098  }
1099  }
1100  // forced biasing only for primary particles
1101  if(biasManager) {
1102  if(0 == track->GetParentID()) {
1103  // primary particle
1104  biasFlag = true;
1106  }
1107  }
1108 }
1109 
1110 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1111 
1114  G4GPILSelection* selection)
1115 {
1116  G4double x = DBL_MAX;
1117  *selection = aGPILSelection;
1120  x = fRange;
1121  G4double finR = finalRange;
1122  if(rndmStepFlag) {
1123  finR = std::min(finR,
1125  }
1126  if(fRange > finR) {
1127  x = fRange*dRoverRange + finR*(1.0 - dRoverRange)*(2.0 - finR/fRange);
1128  }
1129 
1130  // if(particle->GetPDGMass() > 0.9*GeV)
1131  /*
1132  G4cout<<GetProcessName()<<": e= "<<preStepKinEnergy
1133  <<" range= "<<fRange << " idx= " << basedCoupleIndex
1134  << " finR= " << finR
1135  << " limit= " << x <<G4endl;
1136  G4cout << "massRatio= " << massRatio << " Q^2= " << chargeSqRatio
1137  << " finR= " << finR << " dRoverRange= " << dRoverRange
1138  << " finalRange= " << finalRange << G4endl;
1139  */
1140  }
1141  //G4cout<<GetProcessName()<<": e= "<<preStepKinEnergy
1142  //<<" stepLimit= "<<x<<G4endl;
1143  return x;
1144 }
1145 
1146 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1147 
1149  const G4Track& track,
1150  G4double previousStepSize,
1152 {
1153  // condition is set to "Not Forced"
1154  *condition = NotForced;
1155  G4double x = DBL_MAX;
1156 
1157  // initialisation of material, mass, charge, model
1158  // at the beginning of the step
1161  preStepScaledEnergy = preStepKinEnergy*massRatio;
1163 
1167  return x;
1168  }
1169 
1170  // change effective charge of an ion on fly
1171  if(isIon) {
1173  if(q2 != chargeSqRatio && q2 > 0.0) {
1174  chargeSqRatio = q2;
1175  fFactor = q2*biasFactor*(*theDensityFactor)[currentCoupleIndex];
1176  reduceFactor = 1.0/(fFactor*massRatio);
1177  }
1178  }
1179  // if(particle->GetPDGMass() > 0.9*GeV)
1180  //G4cout << "q2= "<<chargeSqRatio << " massRatio= " << massRatio << G4endl;
1181 
1182  // forced biasing only for primary particles
1183  if(biasManager) {
1184  if(0 == track.GetParentID()) {
1185  if(biasFlag &&
1187  return biasManager->GetStepLimit(currentCoupleIndex, previousStepSize);
1188  }
1189  }
1190  }
1191 
1192  // compute mean free path
1196 
1197  // zero cross section
1198  if(preStepLambda <= 0.0) {
1201  }
1202  }
1203 
1204  // non-zero cross section
1205  if(preStepLambda > 0.0) {
1207 
1208  // beggining of tracking (or just after DoIt of this process)
1211 
1212  } else if(currentInteractionLength < DBL_MAX) {
1213 
1214  // subtract NumberOfInteractionLengthLeft using previous step
1216  previousStepSize/currentInteractionLength;
1217 
1220  }
1221 
1222  // new mean free path and step limit
1225  }
1226 #ifdef G4VERBOSE
1227  if (verboseLevel>2){
1228  // if(particle->GetPDGMass() > 0.9*GeV){
1229  G4cout << "G4VEnergyLossProcess::PostStepGetPhysicalInteractionLength ";
1230  G4cout << "[ " << GetProcessName() << "]" << G4endl;
1231  G4cout << " for " << track.GetDefinition()->GetParticleName()
1232  << " in Material " << currentMaterial->GetName()
1233  << " Ekin(MeV)= " << preStepKinEnergy/MeV
1234  << " " << track.GetMaterial()->GetName()
1235  <<G4endl;
1236  G4cout << "MeanFreePath = " << currentInteractionLength/cm << "[cm]"
1237  << "InteractionLength= " << x/cm <<"[cm] " <<G4endl;
1238  }
1239 #endif
1240  return x;
1241 }
1242 
1243 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1244 
1246  const G4Step& step)
1247 {
1249  // The process has range table - calculate energy loss
1251  return &fParticleChange;
1252  }
1253 
1254  // Get the actual (true) Step length
1255  G4double length = step.GetStepLength();
1256  if(length <= 0.0) { return &fParticleChange; }
1257  G4double eloss = 0.0;
1258 
1259  /*
1260  if(-1 < verboseLevel) {
1261  const G4ParticleDefinition* d = track.GetParticleDefinition();
1262  G4cout << "AlongStepDoIt for "
1263  << GetProcessName() << " and particle "
1264  << d->GetParticleName()
1265  << " eScaled(MeV)= " << preStepScaledEnergy/MeV
1266  << " range(mm)= " << fRange/mm
1267  << " s(mm)= " << length/mm
1268  << " rf= " << reduceFactor
1269  << " q^2= " << chargeSqRatio
1270  << " md= " << d->GetPDGMass()
1271  << " status= " << track.GetTrackStatus()
1272  << " " << track.GetMaterial()->GetName()
1273  << G4endl;
1274  }
1275  */
1276 
1277  const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
1278 
1279  // define new weight for primary and secondaries
1281  if(weightFlag) {
1282  weight /= biasFactor;
1284  }
1285 
1286  // stopping
1287  if (length >= fRange || preStepKinEnergy <= lowestKinEnergy) {
1288  eloss = preStepKinEnergy;
1289  if (useDeexcitation) {
1291  eloss, currentCoupleIndex);
1292  if(scTracks.size() > 0) { FillSecondariesAlongStep(eloss, weight); }
1293  if(eloss < 0.0) { eloss = 0.0; }
1294  }
1297  return &fParticleChange;
1298  }
1299  //G4cout << theDEDXTable << " idx= " << basedCoupleIndex
1300  // << " " << GetProcessName() << " "<< currentMaterial->GetName()<<G4endl;
1301  //if(particle->GetParticleName() == "e-")G4cout << (*theDEDXTable) <<G4endl;
1302  // Short step
1304 
1305  //G4cout << "eloss= " << eloss << G4endl;
1306 
1307  // Long step
1308  if(eloss > preStepKinEnergy*linLossLimit) {
1309 
1310  G4double x = (fRange - length)/reduceFactor;
1311  //G4cout << "x= " << x << " " << theInverseRangeTable << G4endl;
1313 
1314  /*
1315  if(-1 < verboseLevel)
1316  G4cout << "Long STEP: rPre(mm)= "
1317  << GetScaledRangeForScaledEnergy(preStepScaledEnergy)/mm
1318  << " rPost(mm)= " << x/mm
1319  << " ePre(MeV)= " << preStepScaledEnergy/MeV
1320  << " eloss(MeV)= " << eloss/MeV
1321  << " eloss0(MeV)= "
1322  << GetDEDXForScaledEnergy(preStepScaledEnergy)*length/MeV
1323  << " lim(MeV)= " << preStepKinEnergy*linLossLimit/MeV
1324  << G4endl;
1325  */
1326  }
1327 
1328  /*
1329  G4double eloss0 = eloss;
1330  if(-1 < verboseLevel ) {
1331  G4cout << "Before fluct: eloss(MeV)= " << eloss/MeV
1332  << " e-eloss= " << preStepKinEnergy-eloss
1333  << " step(mm)= " << length/mm
1334  << " range(mm)= " << fRange/mm
1335  << " fluct= " << lossFluctuationFlag
1336  << G4endl;
1337  }
1338  */
1339 
1340  G4double cut = (*theCuts)[currentCoupleIndex];
1341  G4double esec = 0.0;
1342 
1343  //G4cout << "cut= " << cut << " useSubCut= " << useSubCutoff << G4endl;
1344 
1345  // SubCutOff
1346  if(useSubCutoff && !subcutProducer) {
1348 
1349  G4bool yes = false;
1350  G4StepPoint* prePoint = step.GetPreStepPoint();
1351 
1352  // Check boundary
1353  if(prePoint->GetStepStatus() == fGeomBoundary) { yes = true; }
1354 
1355  // Check PrePoint
1356  else {
1357  G4double preSafety = prePoint->GetSafety();
1358  G4double rcut =
1360 
1361  // recompute presafety
1362  if(preSafety < rcut) {
1363  preSafety = safetyHelper->ComputeSafety(prePoint->GetPosition(),
1364  rcut);
1365  }
1366 
1367  if(preSafety < rcut) { yes = true; }
1368 
1369  // Check PostPoint
1370  else {
1371  G4double postSafety = preSafety - length;
1372  if(postSafety < rcut) {
1373  postSafety = safetyHelper->ComputeSafety(
1374  step.GetPostStepPoint()->GetPosition(), rcut);
1375  if(postSafety < rcut) { yes = true; }
1376  }
1377  }
1378  }
1379 
1380  // Decided to start subcut sampling
1381  if(yes) {
1382 
1383  cut = (*theSubCuts)[currentCoupleIndex];
1385  esec = SampleSubCutSecondaries(scTracks, step,
1386  currentModel,currentCoupleIndex);
1387  // add bremsstrahlung sampling
1388  /*
1389  if(nProcesses > 0) {
1390  for(G4int i=0; i<nProcesses; ++i) {
1391  (scProcesses[i])->SampleSubCutSecondaries(
1392  scTracks, step, (scProcesses[i])->
1393  SelectModelForMaterial(preStepKinEnergy, currentCoupleIndex),
1394  currentCoupleIndex);
1395  }
1396  }
1397  */
1398  }
1399  }
1400  }
1401 
1402  // Corrections, which cannot be tabulated
1403  if(isIon) {
1404  G4double eadd = 0.0;
1405  G4double eloss_before = eloss;
1407  eloss, eadd, length);
1408  if(eloss < 0.0) { eloss = 0.5*eloss_before; }
1409  }
1410 
1411  // Sample fluctuations
1412  if (lossFluctuationFlag) {
1414  if(eloss + esec < preStepKinEnergy) {
1415 
1416  G4double tmax =
1417  std::min(currentModel->MaxSecondaryKinEnergy(dynParticle),cut);
1418  eloss = fluc->SampleFluctuations(currentCouple,dynParticle,
1419  tmax,length,eloss);
1420  /*
1421  if(-1 < verboseLevel)
1422  G4cout << "After fluct: eloss(MeV)= " << eloss/MeV
1423  << " fluc= " << (eloss-eloss0)/MeV
1424  << " ChargeSqRatio= " << chargeSqRatio
1425  << " massRatio= " << massRatio
1426  << " tmax= " << tmax
1427  << G4endl;
1428  */
1429  }
1430  }
1431 
1432  // deexcitation
1433  if (useDeexcitation) {
1434  G4double esecfluo = preStepKinEnergy - esec;
1435  G4double de = esecfluo;
1436  //G4double eloss0 = eloss;
1437  /*
1438  G4cout << "### 1: E(keV)= " << preStepKinEnergy/keV
1439  << " Efluomax(keV)= " << de/keV
1440  << " Eloss(keV)= " << eloss/keV << G4endl;
1441  */
1443  de, currentCoupleIndex);
1444 
1445  // sum of de-excitation energies
1446  esecfluo -= de;
1447 
1448  // subtracted from energy loss
1449  if(eloss >= esecfluo) {
1450  esec += esecfluo;
1451  eloss -= esecfluo;
1452  } else {
1453  esec += esecfluo;
1454  eloss = 0.0;
1455  }
1456  /*
1457  if(esecfluo > 0.0) {
1458  G4cout << "### 2: E(keV)= " << preStepKinEnergy/keV
1459  << " Esec(keV)= " << esec/keV
1460  << " Esecf(kV)= " << esecfluo/keV
1461  << " Eloss0(kV)= " << eloss0/keV
1462  << " Eloss(keV)= " << eloss/keV
1463  << G4endl;
1464  }
1465  */
1466  }
1468  subcutProducer->SampleSecondaries(step, scTracks, eloss, cut);
1469  }
1470  if(scTracks.size() > 0) { FillSecondariesAlongStep(eloss, weight); }
1471 
1472  // Energy balance
1473  G4double finalT = preStepKinEnergy - eloss - esec;
1474  if (finalT <= lowestKinEnergy) {
1475  eloss += finalT;
1476  finalT = 0.0;
1477  } else if(isIon) {
1480  currentMaterial,finalT));
1481  }
1482 
1483  eloss = std::max(eloss, 0.0);
1484 
1487  /*
1488  if(-1 < verboseLevel) {
1489  G4double del = finalT + eloss + esec - preStepKinEnergy;
1490  G4cout << "Final value eloss(MeV)= " << eloss/MeV
1491  << " preStepKinEnergy= " << preStepKinEnergy
1492  << " postStepKinEnergy= " << finalT
1493  << " de(keV)= " << del/keV
1494  << " lossFlag= " << lossFluctuationFlag
1495  << " status= " << track.GetTrackStatus()
1496  << G4endl;
1497  }
1498  */
1499  return &fParticleChange;
1500 }
1501 
1502 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1503 
1504 void
1506 {
1507  G4int n0 = scTracks.size();
1508 
1509  // weight may be changed by biasing manager
1510  if(biasManager) {
1512  weight *=
1514  }
1515  }
1516 
1517  // fill secondaries
1518  G4int n = scTracks.size();
1520 
1521  for(G4int i=0; i<n; ++i) {
1522  G4Track* t = scTracks[i];
1523  if(t) {
1524  t->SetWeight(weight);
1526  if(i >= n0) { t->SetCreatorModelIndex(biasID); }
1527  //G4cout << "Secondary(along step) has weight " << t->GetWeight()
1528  //<< ", kenergy " << t->GetKineticEnergy()/MeV << " MeV" <<G4endl;
1529  }
1530  }
1531  scTracks.clear();
1532 }
1533 
1534 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1535 
1536 G4double
1537 G4VEnergyLossProcess::SampleSubCutSecondaries(std::vector<G4Track*>& tracks,
1538  const G4Step& step,
1539  G4VEmModel* model,
1540  G4int idx)
1541 {
1542  // Fast check weather subcutoff can work
1543  G4double esec = 0.0;
1544  G4double subcut = (*theSubCuts)[idx];
1545  G4double cut = (*theCuts)[idx];
1546  if(cut <= subcut) { return esec; }
1547 
1548  const G4Track* track = step.GetTrack();
1549  const G4DynamicParticle* dp = track->GetDynamicParticle();
1551  G4double cross = (*theDensityFactor)[idx]*chargeSqRatio
1552  *(((*theSubLambdaTable)[(*theDensityIdx)[idx]])->Value(e, idxSubLambda));
1553  G4double length = step.GetStepLength();
1554 
1555  // negligible probability to get any interaction
1556  if(length*cross < perMillion) { return esec; }
1557  /*
1558  if(-1 < verboseLevel)
1559  G4cout << "<<< Subcutoff for " << GetProcessName()
1560  << " cross(1/mm)= " << cross*mm << ">>>"
1561  << " e(MeV)= " << preStepScaledEnergy
1562  << " matIdx= " << currentCoupleIndex
1563  << G4endl;
1564  */
1565 
1566  // Sample subcutoff secondaries
1567  G4StepPoint* preStepPoint = step.GetPreStepPoint();
1568  G4StepPoint* postStepPoint = step.GetPostStepPoint();
1569  G4ThreeVector prepoint = preStepPoint->GetPosition();
1570  G4ThreeVector dr = postStepPoint->GetPosition() - prepoint;
1571  G4double pretime = preStepPoint->GetGlobalTime();
1572  G4double dt = postStepPoint->GetGlobalTime() - pretime;
1573  //G4double dt = length/preStepPoint->GetVelocity();
1574  G4double fragment = 0.0;
1575 
1576  do {
1577  G4double del = -G4Log(G4UniformRand())/cross;
1578  fragment += del/length;
1579  if (fragment > 1.0) { break; }
1580 
1581  // sample secondaries
1582  secParticles.clear();
1584  dp,subcut,cut);
1585 
1586  // position of subcutoff particles
1587  G4ThreeVector r = prepoint + fragment*dr;
1588  std::vector<G4DynamicParticle*>::iterator it;
1589  for(it=secParticles.begin(); it!=secParticles.end(); ++it) {
1590 
1591  G4Track* t = new G4Track((*it), pretime + fragment*dt, r);
1594  tracks.push_back(t);
1595  esec += t->GetKineticEnergy();
1596  if (t->GetParticleDefinition() == thePositron) {
1597  esec += 2.0*electron_mass_c2;
1598  }
1599 
1600  /*
1601  if(-1 < verboseLevel)
1602  G4cout << "New track "
1603  << t->GetParticleDefinition()->GetParticleName()
1604  << " e(keV)= " << t->GetKineticEnergy()/keV
1605  << " fragment= " << fragment
1606  << G4endl;
1607  */
1608  }
1609  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
1610  } while (fragment <= 1.0);
1611  return esec;
1612 }
1613 
1614 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1615 
1617  const G4Step& step)
1618 {
1619  // In all cases clear number of interaction lengths
1622 
1624  G4double finalT = track.GetKineticEnergy();
1625  if(finalT <= lowestKinEnergy) { return &fParticleChange; }
1626 
1627  G4double postStepScaledEnergy = finalT*massRatio;
1628 
1629  if(!currentModel->IsActive(postStepScaledEnergy)) {
1630  return &fParticleChange;
1631  }
1632  /*
1633  if(-1 < verboseLevel) {
1634  G4cout << GetProcessName()
1635  << "::PostStepDoIt: E(MeV)= " << finalT/MeV
1636  << G4endl;
1637  }
1638  */
1639 
1640  // forced process - should happen only once per track
1641  if(biasFlag) {
1643  biasFlag = false;
1644  }
1645  }
1646 
1647  // Integral approach
1648  if (integral) {
1649  G4double lx = GetLambdaForScaledEnergy(postStepScaledEnergy);
1650  /*
1651  if(preStepLambda<lx && 1 < verboseLevel && nWarnings<200) {
1652  G4cout << "WARNING: for " << particle->GetParticleName()
1653  << " and " << GetProcessName()
1654  << " E(MeV)= " << finalT/MeV
1655  << " preLambda= " << preStepLambda
1656  << " < " << lx << " (postLambda) "
1657  << G4endl;
1658  ++nWarnings;
1659  }
1660  */
1661  if(lx <= 0.0 || preStepLambda*G4UniformRand() > lx) {
1662  return &fParticleChange;
1663  }
1664  }
1665 
1666  SelectModel(postStepScaledEnergy);
1667 
1668  // define new weight for primary and secondaries
1670  if(weightFlag) {
1671  weight /= biasFactor;
1673  }
1674 
1675  const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
1676  G4double tcut = (*theCuts)[currentCoupleIndex];
1677 
1678  // sample secondaries
1679  secParticles.clear();
1680  //G4cout<< "@@@ Eprimary= "<<dynParticle->GetKineticEnergy()/MeV
1681  // << " cut= " << tcut/MeV << G4endl;
1683  dynParticle, tcut);
1684 
1685  G4int num0 = secParticles.size();
1686 
1687  // bremsstrahlung splitting or Russian roulette
1688  if(biasManager) {
1690  G4double eloss = 0.0;
1692  secParticles,
1693  track, currentModel,
1694  &fParticleChange, eloss,
1695  currentCoupleIndex, tcut,
1696  step.GetPostStepPoint()->GetSafety());
1697  if(eloss > 0.0) {
1700  }
1701  }
1702  }
1703 
1704  // save secondaries
1705  G4int num = secParticles.size();
1706  if(num > 0) {
1707 
1709  G4double time = track.GetGlobalTime();
1710 
1711  for (G4int i=0; i<num; ++i) {
1712  if(secParticles[i]) {
1713  G4Track* t = new G4Track(secParticles[i], time, track.GetPosition());
1715  t->SetWeight(weight);
1716  if(i < num0) { t->SetCreatorModelIndex(secID); }
1717  else { t->SetCreatorModelIndex(biasID); }
1718 
1719  //G4cout << "Secondary(post step) has weight " << t->GetWeight()
1720  // << ", kenergy " << t->GetKineticEnergy()/MeV << " MeV"
1721  // << " time= " << time/ns << " ns " << G4endl;
1723  }
1724  }
1725  }
1726 
1732  }
1733 
1734  /*
1735  if(-1 < verboseLevel) {
1736  G4cout << "::PostStepDoIt: Sample secondary; Efin= "
1737  << fParticleChange.GetProposedKineticEnergy()/MeV
1738  << " MeV; model= (" << currentModel->LowEnergyLimit()
1739  << ", " << currentModel->HighEnergyLimit() << ")"
1740  << " preStepLambda= " << preStepLambda
1741  << " dir= " << track.GetMomentumDirection()
1742  << " status= " << track.GetTrackStatus()
1743  << G4endl;
1744  }
1745  */
1746  return &fParticleChange;
1747 }
1748 
1749 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1750 
1752  const G4ParticleDefinition* part, const G4String& directory,
1753  G4bool ascii)
1754 {
1755  G4bool res = true;
1756  //G4cout << "G4VEnergyLossProcess::StorePhysicsTable: " << part->GetParticleName()
1757  // << " " << directory << " " << ascii << G4endl;
1758  if (!isMaster || baseParticle || part != particle ) return res;
1759 
1760  if(!StoreTable(part,theDEDXTable,ascii,directory,"DEDX"))
1761  {res = false;}
1762 
1763  if(!StoreTable(part,theDEDXunRestrictedTable,ascii,directory,"DEDXnr"))
1764  {res = false;}
1765 
1766  if(!StoreTable(part,theDEDXSubTable,ascii,directory,"SubDEDX"))
1767  {res = false;}
1768 
1769  if(!StoreTable(part,theIonisationTable,ascii,directory,"Ionisation"))
1770  {res = false;}
1771 
1772  if(!StoreTable(part,theIonisationSubTable,ascii,directory,"SubIonisation"))
1773  {res = false;}
1774 
1775  if(isIonisation &&
1776  !StoreTable(part,theCSDARangeTable,ascii,directory,"CSDARange"))
1777  {res = false;}
1778 
1779  if(isIonisation &&
1780  !StoreTable(part,theRangeTableForLoss,ascii,directory,"Range"))
1781  {res = false;}
1782 
1783  if(isIonisation &&
1784  !StoreTable(part,theInverseRangeTable,ascii,directory,"InverseRange"))
1785  {res = false;}
1786 
1787  if(!StoreTable(part,theLambdaTable,ascii,directory,"Lambda"))
1788  {res = false;}
1789 
1790  if(!StoreTable(part,theSubLambdaTable,ascii,directory,"SubLambda"))
1791  {res = false;}
1792 
1793  if ( !res ) {
1794  if(1 < verboseLevel) {
1795  G4cout << "Physics tables are stored for "
1797  << " and process " << GetProcessName()
1798  << " in the directory <" << directory
1799  << "> " << G4endl;
1800  }
1801  } else {
1802  G4cout << "Fail to store Physics Tables for "
1804  << " and process " << GetProcessName()
1805  << " in the directory <" << directory
1806  << "> " << G4endl;
1807  }
1808  return res;
1809 }
1810 
1811 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
1812 
1813 G4bool
1815  const G4String& directory,
1816  G4bool ascii)
1817 {
1818  G4bool res = true;
1819  if (!isMaster) return res;
1820  const G4String particleName = part->GetParticleName();
1821 
1822  if(1 < verboseLevel) {
1823  G4cout << "G4VEnergyLossProcess::RetrievePhysicsTable() for "
1824  << particleName << " and process " << GetProcessName()
1825  << "; tables_are_built= " << tablesAreBuilt
1826  << G4endl;
1827  }
1828  if(particle == part) {
1829 
1830  if ( !baseParticle ) {
1831 
1832  G4bool fpi = true;
1833  if(!RetrieveTable(part,theDEDXTable,ascii,directory,"DEDX",fpi))
1834  {fpi = false;}
1835 
1836  // ionisation table keeps individual dEdx and not sum of sub-processes
1837  if(!RetrieveTable(part,theDEDXTable,ascii,directory,"Ionisation",false))
1838  {fpi = false;}
1839 
1840  if(!RetrieveTable(part,theRangeTableForLoss,ascii,directory,"Range",fpi))
1841  {res = false;}
1842 
1843  if(!RetrieveTable(part,theDEDXunRestrictedTable,ascii,directory,
1844  "DEDXnr",false))
1845  {res = false;}
1846 
1847  if(!RetrieveTable(part,theCSDARangeTable,ascii,directory,
1848  "CSDARange",false))
1849  {res = false;}
1850 
1851  if(!RetrieveTable(part,theInverseRangeTable,ascii,directory,
1852  "InverseRange",fpi))
1853  {res = false;}
1854 
1855  if(!RetrieveTable(part,theLambdaTable,ascii,directory,"Lambda",true))
1856  {res = false;}
1857 
1858  G4bool yes = false;
1859  if(nSCoffRegions > 0) {yes = true;}
1860 
1861  if(!RetrieveTable(part,theDEDXSubTable,ascii,directory,"SubDEDX",yes))
1862  {res = false;}
1863 
1864  if(!RetrieveTable(part,theSubLambdaTable,ascii,directory,
1865  "SubLambda",yes))
1866  {res = false;}
1867 
1868  if(!fpi) yes = false;
1869  if(!RetrieveTable(part,theIonisationSubTable,ascii,directory,
1870  "SubIonisation",yes))
1871  {res = false;}
1872  }
1873  }
1874 
1875  return res;
1876 }
1877 
1878 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
1879 
1881  G4PhysicsTable* aTable, G4bool ascii,
1882  const G4String& directory,
1883  const G4String& tname)
1884 {
1885  //G4cout << "G4VEnergyLossProcess::StoreTable: " << aTable
1886  // << " " << directory << " " << tname << G4endl;
1887  G4bool res = true;
1888  if ( aTable ) {
1889  const G4String name = GetPhysicsTableFileName(part,directory,tname,ascii);
1890  G4cout << name << G4endl;
1891  //G4cout << *aTable << G4endl;
1892  if( !aTable->StorePhysicsTable(name,ascii)) res = false;
1893  }
1894  return res;
1895 }
1896 
1897 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
1898 
1899 G4bool
1901  G4PhysicsTable* aTable,
1902  G4bool ascii,
1903  const G4String& directory,
1904  const G4String& tname,
1905  G4bool mandatory)
1906 {
1907  G4bool isRetrieved = false;
1908  G4String filename = GetPhysicsTableFileName(part,directory,tname,ascii);
1909  if(aTable) {
1910  if(aTable->ExistPhysicsTable(filename)) {
1911  if(G4PhysicsTableHelper::RetrievePhysicsTable(aTable,filename,ascii)) {
1912  isRetrieved = true;
1913  if(theParameters->Spline()) {
1914  size_t n = aTable->length();
1915  for(size_t i=0; i<n; ++i) {
1916  if((*aTable)[i]) { (*aTable)[i]->SetSpline(true); }
1917  }
1918  }
1919  if (0 < verboseLevel) {
1920  G4cout << tname << " table for " << part->GetParticleName()
1921  << " is Retrieved from <" << filename << ">"
1922  << G4endl;
1923  }
1924  }
1925  }
1926  }
1927  if(mandatory && !isRetrieved) {
1928  if(0 < verboseLevel) {
1929  G4cout << tname << " table for " << part->GetParticleName()
1930  << " from file <"
1931  << filename << "> is not Retrieved"
1932  << G4endl;
1933  }
1934  return false;
1935  }
1936  return true;
1937 }
1938 
1939 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1940 
1942  const G4MaterialCutsCouple *couple,
1943  const G4DynamicParticle* dp,
1944  G4double length)
1945 {
1946  DefineMaterial(couple);
1947  G4double ekin = dp->GetKineticEnergy();
1948  SelectModel(ekin*massRatio);
1950  tmax = std::min(tmax,(*theCuts)[currentCoupleIndex]);
1951  G4double d = 0.0;
1953  if(fm) { d = fm->Dispersion(currentMaterial,dp,tmax,length); }
1954  return d;
1955 }
1956 
1957 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1958 
1960  G4double kineticEnergy, const G4MaterialCutsCouple* couple)
1961 {
1962  // Cross section per volume is calculated
1963  DefineMaterial(couple);
1964  G4double cross = 0.0;
1965  if(theLambdaTable) {
1966  cross = GetLambdaForScaledEnergy(kineticEnergy*massRatio);
1967  } else {
1968  SelectModel(kineticEnergy*massRatio);
1969  cross = biasFactor*(*theDensityFactor)[currentCoupleIndex]
1971  particle, kineticEnergy,
1973  }
1974  if(cross < 0.0) { cross = 0.0; }
1975  return cross;
1976 }
1977 
1978 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1979 
1981 {
1984  G4double x = DBL_MAX;
1985  if(0.0 < preStepLambda) { x = 1.0/preStepLambda; }
1986  return x;
1987 }
1988 
1989 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1990 
1992  G4double x, G4double y,
1993  G4double& z)
1994 {
1995  G4GPILSelection sel;
1996  return AlongStepGetPhysicalInteractionLength(track, x, y, z, &sel);
1997 }
1998 
1999 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2000 
2002  const G4Track& track,
2003  G4double,
2005 
2006 {
2007  *condition = NotForced;
2008  return MeanFreePath(track);
2009 }
2010 
2011 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2012 
2014  const G4Track&,
2016 {
2017  return DBL_MAX;
2018 }
2019 
2020 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2021 
2024  G4double)
2025 {
2026  G4PhysicsVector* v =
2029  return v;
2030 }
2031 
2032 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2033 
2036 {
2037  G4bool add = true;
2038  if(p->GetProcessName() != "eBrem") { add = false; }
2039  if(add && nProcesses > 0) {
2040  for(G4int i=0; i<nProcesses; ++i) {
2041  if(p == scProcesses[i]) {
2042  add = false;
2043  break;
2044  }
2045  }
2046  }
2047  if(add) {
2048  scProcesses.push_back(p);
2049  ++nProcesses;
2050  if (1 < verboseLevel) {
2051  G4cout << "### The process " << p->GetProcessName()
2052  << " is added to the list of collaborative processes of "
2053  << GetProcessName() << G4endl;
2054  }
2055  }
2056 }
2057 
2058 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2059 
2060 void
2062 {
2063  if(fTotal == tType) {
2065  if(p) {
2066  size_t n = p->length();
2067  G4PhysicsVector* pv = (*p)[0];
2069 
2073 
2074  for (size_t i=0; i<n; ++i) {
2075  G4double dedx = 0.0;
2076  pv = (*p)[i];
2077  if(pv) {
2078  dedx = pv->Value(emax, idxDEDXunRestricted);
2079  } else {
2080  pv = (*p)[(*theDensityIdx)[i]];
2081  if(pv) {
2082  dedx =
2083  pv->Value(emax, idxDEDXunRestricted)*(*theDensityFactor)[i];
2084  }
2085  }
2086  theDEDXAtMaxEnergy[i] = dedx;
2087  //G4cout << "i= " << i << " emax(MeV)= " << emax/MeV<< " dedx= "
2088  // << dedx << G4endl;
2089  }
2090  }
2091 
2092  } else if(fRestricted == tType) {
2093  /*
2094  G4cout<< "G4VEnergyLossProcess::SetDEDXTable "
2095  << particle->GetParticleName()
2096  << " oldTable " << theDEDXTable << " newTable " << p
2097  << " ion " << theIonisationTable
2098  << " IsMaster " << isMaster
2099  << " " << GetProcessName() << G4endl;
2100  G4cout << (*p) << G4endl;
2101  */
2102  theDEDXTable = p;
2103  } else if(fSubRestricted == tType) {
2104  theDEDXSubTable = p;
2105  } else if(fIsIonisation == tType) {
2106  /*
2107  G4cout<< "G4VEnergyLossProcess::SetIonisationTable "
2108  << particle->GetParticleName()
2109  << " oldTable " << theDEDXTable << " newTable " << p
2110  << " ion " << theIonisationTable
2111  << " IsMaster " << isMaster
2112  << " " << GetProcessName() << G4endl;
2113  */
2114  theIonisationTable = p;
2115  } else if(fIsSubIonisation == tType) {
2117  }
2118 }
2119 
2120 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2121 
2123 {
2124  theCSDARangeTable = p;
2125 
2126  if(p) {
2127  size_t n = p->length();
2128  G4PhysicsVector* pv;
2130 
2131  for (size_t i=0; i<n; ++i) {
2132  pv = (*p)[i];
2133  G4double rmax = 0.0;
2134  if(pv) { rmax = pv->Value(emax, idxCSDA); }
2135  else {
2136  pv = (*p)[(*theDensityIdx)[i]];
2137  if(pv) { rmax = pv->Value(emax, idxCSDA)/(*theDensityFactor)[i]; }
2138  }
2139  theRangeAtMaxEnergy[i] = rmax;
2140  //G4cout << "i= " << i << " Emax(MeV)= " << emax/MeV << " Rmax= "
2141  //<< rmax<< G4endl;
2142  }
2143  }
2144 }
2145 
2146 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2147 
2149 {
2151  if(1 < verboseLevel) {
2152  G4cout << "### Set Range table " << p
2153  << " for " << particle->GetParticleName()
2154  << " and process " << GetProcessName() << G4endl;
2155  }
2156 }
2157 
2158 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2159 
2161 {
2163  if(1 < verboseLevel) {
2164  G4cout << "### Set SecondaryRange table " << p
2165  << " for " << particle->GetParticleName()
2166  << " and process " << GetProcessName() << G4endl;
2167  }
2168 }
2169 
2170 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2171 
2173 {
2175  if(1 < verboseLevel) {
2176  G4cout << "### Set InverseRange table " << p
2177  << " for " << particle->GetParticleName()
2178  << " and process " << GetProcessName() << G4endl;
2179  }
2180 }
2181 
2182 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2183 
2185 {
2186  if(1 < verboseLevel) {
2187  G4cout << "### Set Lambda table " << p
2188  << " for " << particle->GetParticleName()
2189  << " and process " << GetProcessName() << G4endl;
2190  //G4cout << *p << G4endl;
2191  }
2192  theLambdaTable = p;
2193  tablesAreBuilt = true;
2194 
2198 
2199  if(theLambdaTable) {
2200  size_t n = theLambdaTable->length();
2201  G4PhysicsVector* pv = (*theLambdaTable)[0];
2202  G4double e, ss, smax, emax;
2203 
2204  size_t i;
2205 
2206  // first loop on existing vectors
2207  for (i=0; i<n; ++i) {
2208  pv = (*theLambdaTable)[i];
2209  if(pv) {
2210  size_t nb = pv->GetVectorLength();
2211  emax = DBL_MAX;
2212  smax = 0.0;
2213  if(nb > 0) {
2214  for (size_t j=0; j<nb; ++j) {
2215  e = pv->Energy(j);
2216  ss = (*pv)(j);
2217  if(ss > smax) {
2218  smax = ss;
2219  emax = e;
2220  }
2221  }
2222  }
2224  theCrossSectionMax[i] = smax;
2225  if(1 < verboseLevel) {
2226  G4cout << "For " << particle->GetParticleName()
2227  << " Max CS at i= " << i << " emax(MeV)= " << emax/MeV
2228  << " lambda= " << smax << G4endl;
2229  }
2230  }
2231  }
2232  // second loop using base materials
2233  for (i=0; i<n; ++i) {
2234  pv = (*theLambdaTable)[i];
2235  if(!pv){
2236  G4int j = (*theDensityIdx)[i];
2238  theCrossSectionMax[i] = (*theDensityFactor)[i]*theCrossSectionMax[j];
2239  }
2240  }
2241  }
2242 }
2243 
2244 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2245 
2247 {
2248  theSubLambdaTable = p;
2249  if(1 < verboseLevel) {
2250  G4cout << "### Set SebLambda table " << p
2251  << " for " << particle->GetParticleName()
2252  << " and process " << GetProcessName() << G4endl;
2253  }
2254 }
2255 
2256 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2257 
2259 {
2260  const G4Element* elm = nullptr;
2261  if(currentModel) { elm = currentModel->GetCurrentElement(); }
2262  return elm;
2263 }
2264 
2265 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2266 
2268  G4bool flag)
2269 {
2270  if(f > 0.0) {
2271  biasFactor = f;
2272  weightFlag = flag;
2273  if(1 < verboseLevel) {
2274  G4cout << "### SetCrossSectionBiasingFactor: for "
2275  << " process " << GetProcessName()
2276  << " biasFactor= " << f << " weightFlag= " << flag
2277  << G4endl;
2278  }
2279  }
2280 }
2281 
2282 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2283 
2284 void
2286  const G4String& region,
2287  G4bool flag)
2288 {
2289  if(!biasManager) { biasManager = new G4EmBiasingManager(); }
2290  if(1 < verboseLevel) {
2291  G4cout << "### ActivateForcedInteraction: for "
2292  << " process " << GetProcessName()
2293  << " length(mm)= " << length/mm
2294  << " in G4Region <" << region
2295  << "> weightFlag= " << flag
2296  << G4endl;
2297  }
2298  weightFlag = flag;
2299  biasManager->ActivateForcedInteraction(length, region);
2300 }
2301 
2302 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2303 
2304 void
2306  G4double factor,
2307  G4double energyLimit)
2308 {
2309  if (0.0 <= factor) {
2310 
2311  // Range cut can be applied only for e-
2312  if(0.0 == factor && secondaryParticle != G4Electron::Electron())
2313  { return; }
2314 
2315  if(!biasManager) { biasManager = new G4EmBiasingManager(); }
2316  biasManager->ActivateSecondaryBiasing(region, factor, energyLimit);
2317  if(1 < verboseLevel) {
2318  G4cout << "### ActivateSecondaryBiasing: for "
2319  << " process " << GetProcessName()
2320  << " factor= " << factor
2321  << " in G4Region <" << region
2322  << "> energyLimit(MeV)= " << energyLimit/MeV
2323  << G4endl;
2324  }
2325  }
2326 }
2327 
2328 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2329 
2331 {
2332  isIonisation = val;
2333  if(val) { aGPILSelection = CandidateForSelection; }
2335 }
2336 
2337 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2338 
2340 {
2341  if(0.0 < val && val < 1.0) {
2342  linLossLimit = val;
2343  actLinLossLimit = true;
2344  } else { PrintWarning("SetLinearLossLimit", val); }
2345 }
2346 
2347 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2348 
2350 {
2351  if(0.0 < v1 && 0.0 < v2 && v2 < 1.e+50) {
2352  dRoverRange = std::min(1.0, v1);
2353  finalRange = v2;
2354  } else if(v1 <= 0.0) {
2355  PrintWarning("SetStepFunction", v1);
2356  } else {
2357  PrintWarning("SetStepFunction", v2);
2358  }
2359 }
2360 
2361 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2362 
2364 {
2365  if(1.e-18 < val && val < 1.e+50) { lowestKinEnergy = val; }
2366  else { PrintWarning("SetLowestEnergyLimit", val); }
2367 }
2368 
2369 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2370 
2372 {
2373  if(2 < n && n < 1000000000) {
2374  nBins = n;
2375  actBinning = true;
2376  } else {
2377  G4double e = (G4double)n;
2378  PrintWarning("SetDEDXBinning", e);
2379  }
2380 }
2381 
2382 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2383 
2385 {
2386  if(1.e-18 < e && e < maxKinEnergy) {
2387  minKinEnergy = e;
2388  actMinKinEnergy = true;
2389  } else { PrintWarning("SetMinKinEnergy", e); }
2390 }
2391 
2392 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2393 
2395 {
2396  if(minKinEnergy < e && e < 1.e+50) {
2397  maxKinEnergy = e;
2398  actMaxKinEnergy = true;
2399  if(e < maxKinEnergyCSDA) { maxKinEnergyCSDA = e; }
2400  } else { PrintWarning("SetMaxKinEnergy", e); }
2401 }
2402 
2403 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2404 
2406 {
2407  G4String ss = "G4VEnergyLossProcess::" + tit;
2409  ed << "Parameter is out of range: " << val
2410  << " it will have no effect!\n" << " Process "
2411  << GetProcessName() << " nbins= " << nBins
2412  << " Emin(keV)= " << minKinEnergy/keV
2413  << " Emax(GeV)= " << maxKinEnergy/GeV;
2414  G4Exception(ss, "em0044", JustWarning, ed);
2415 }
2416 
2417 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2418 
2419 void G4VEnergyLossProcess::ProcessDescription(std::ostream& outFile) const
2420 {
2421  outFile << "Energy loss process <" << GetProcessName() << ">" << G4endl;
2422 }
2423 
2424 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
G4VEnergyLossProcess(const G4String &name="EnergyLoss", G4ProcessType type=fElectromagnetic)
static const double cm
Definition: G4SIunits.hh:118
G4double condition(const G4ErrorSymMatrix &m)
G4bool UseCutAsFinalRange() const
G4VSubCutProducer * SubCutProducer()
G4double MaxSecondaryKinEnergy(const G4DynamicParticle *dynParticle)
Definition: G4VEmModel.hh:482
G4int NumberOfBinsPerDecade() const
G4ParticleDefinition * GetDefinition() const
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
G4LossTableManager * lManager
std::vector< G4double > theEnergyOfCrossSectionMax
void SetDEDXTable(G4PhysicsTable *p, G4EmTableType tType)
G4int GetParentID() const
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:258
static const double MeV
Definition: G4SIunits.hh:211
G4int NumberOfBins() const
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
void InitializeForPostStep(const G4Track &)
const G4String & GetName() const
G4SafetyHelper * GetSafetyHelper() const
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
const std::vector< G4double > * GetDensityFactors()
static void SetPhysicsVector(G4PhysicsTable *physTable, size_t idx, G4PhysicsVector *vec)
void InitialiseHelper()
G4bool Spline() const
G4int WorkerVerbose() const
G4double MaxKinEnergy() const
static G4LossTableManager * Instance()
void ActivateForcedInteraction(G4double length=0.0, const G4String &r="")
G4GPILSelection aGPILSelection
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)
G4double ApplySecondaryBiasing(std::vector< G4DynamicParticle * > &, const G4Track &track, G4VEmModel *currentModel, G4ParticleChangeForGamma *pParticleChange, G4double &eloss, G4int coupleIdx, G4double tcut, G4double safety=0.0)
G4int verboseLevel
Definition: G4VProcess.hh:368
virtual void CorrectionsAlongStep(const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
Definition: G4VEmModel.cc:362
void UpdateEmModel(const G4String &, G4double, G4double)
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4double GetKineticEnergy() const
void PrintInfoDefinition(const G4ParticleDefinition &part)
CLHEP::Hep3Vector G4ThreeVector
void SetIonisation(G4bool val)
G4double GetStepLength() const
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:634
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
const G4DynamicParticle * GetDynamicParticle() const
G4bool SecondaryBiasingRegion(G4int coupleIdx)
G4bool IsPIXEActive() const
G4PhysicsTable * SubLambdaTable() const
G4double z
Definition: TRTMaterials.hh:39
G4double LowestElectronEnergy() const
G4VEmModel * GetModel(G4int, G4bool ver=false)
G4double GetProductionCut(G4int index) const
void DeRegister(G4VEnergyLossProcess *p)
G4PhysicsTable * RangeTableForLoss() const
virtual G4double GetContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety)
G4String name
Definition: TRTMaterials.hh:40
const std::vector< G4double > * theDensityFactor
G4StepStatus GetStepStatus() const
G4bool BuildCSDARange() const
const G4String & GetName() const
Definition: G4Material.hh:178
G4PhysicsTable * SecondaryRangeTable() const
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
void AddCollaborativeProcess(G4VEnergyLossProcess *)
const G4ThreeVector & GetPosition() const
const std::vector< G4int > * theDensityIdx
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel)
Definition: G4VEmModel.cc:219
void ActivateForcedInteraction(G4double length=0.0, const G4String &region="", G4bool flag=true)
void SetLinearLossLimit(G4double val)
const G4ParticleDefinition * secondaryParticle
virtual G4double GetMeanFreePath(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
void SetStepFunction(G4double v1, G4double v2)
G4bool ForcedInteractionRegion(G4int coupleIdx)
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4PhysicsTable * CSDARangeTable() const
G4PhysicsTable * theSecondaryRangeTable
virtual G4bool StorePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false)
G4VEmFluctuationModel * fluctModel
G4ParticleChangeForLoss fParticleChange
G4VSubCutProducer * subcutProducer
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
void SelectModel(G4double kinEnergy)
G4bool GetFlag(size_t idx) const
G4PhysicsTable * theIonisationTable
void SetTouchableHandle(const G4TouchableHandle &apValue)
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const
virtual G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
G4PhysicsTable * IonisationTableForSubsec() const
G4double GetParentWeight() const
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
G4VEmFluctuationModel * GetModelOfFluctuations()
Definition: G4VEmModel.hh:610
const G4DataVector * SubCutoff() const
size_t GetVectorLength() const
G4PhysicsVector * LambdaPhysicsVector(const G4MaterialCutsCouple *, G4double cut)
void InitializeForAlongStep(const G4Track &)
G4double GetProposedKineticEnergy() const
G4PhysicsTable * IonisationTable() const
void FillLambdaVector(G4PhysicsVector *, const G4MaterialCutsCouple *, G4bool startFromNull=true, G4EmTableType t=fRestricted)
G4ProcessManager * GetProcessManager() const
G4VAtomDeexcitation * atomDeexcitation
int G4int
Definition: G4Types.hh:78
void FillSecondDerivatives()
G4PhysicsTable * BuildLambdaTable(G4EmTableType tType=fRestricted)
const G4String & GetPhysicsTableFileName(const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
Definition: G4VProcess.cc:186
static G4PhysicsTable * PreparePhysicsTable(G4PhysicsTable *physTable)
const G4int smax
G4double MinSubRange() const
const G4String & GetParticleName() const
G4LossTableBuilder * GetTableBuilder()
void SetWeight(G4double aValue)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void SetCreatorModelIndex(G4int idx)
G4bool ExistPhysicsTable(const G4String &fileName) const
const G4ParticleDefinition * theGamma
static const G4double reg
const G4MaterialCutsCouple * currentCouple
void SetSpline(G4bool)
void PrintWarning(G4String, G4double val)
void Initialise(const G4ParticleDefinition &part, const G4String &procName, G4int verbose)
const G4DataVector * theCuts
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:725
virtual G4double Dispersion(const G4Material *, const G4DynamicParticle *, G4double tmax, G4double length)=0
static G4RegionStore * GetInstance()
G4PhysicsTable * theSubLambdaTable
G4StepPoint * GetPreStepPoint() const
G4PhysicsTable * LambdaTable() const
void ProposeWeight(G4double finalWeight)
void SetAngularGeneratorFlag(G4bool)
Definition: G4VEmModel.hh:704
G4double LinearLossLimit() const
void SetSecondaryWeightByProcess(G4bool)
void SetInverseRangeTable(G4PhysicsTable *p)
G4double LambdaFactor() const
void SetCrossSectionBiasingFactor(G4double f, G4bool flag=true)
G4ProcessVector * GetAtRestProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4double GetKineticEnergy() const
G4double SampleSubCutSecondaries(std::vector< G4Track * > &, const G4Step &, G4VEmModel *model, G4int matIdx)
const G4ParticleDefinition * baseParticle
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
void FillDEDXVector(G4PhysicsVector *, const G4MaterialCutsCouple *, G4EmTableType t=fRestricted)
G4int Verbose() const
void FillSecondariesAlongStep(G4double &eloss, G4double &weight)
const G4DataVector * theSubCuts
G4PhysicsTable * DEDXTable() const
G4double ScaledKinEnergyForLoss(G4double range)
G4double GetDEDXForScaledEnergy(G4double scaledKinEnergy)
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=0)
Definition: G4VEmModel.cc:410
const G4ThreeVector & GetPosition() const
const G4Element * GetCurrentElement() const
G4double LowestMuHadEnergy() const
const G4ParticleDefinition * theGenericIon
const G4ParticleDefinition * thePositron
void ComputeLambdaForScaledEnergy(G4double scaledKinEnergy)
bool G4bool
Definition: G4Types.hh:79
std::vector< G4Track * > scTracks
G4double GetStepLimit(G4int coupleIdx, G4double previousStep)
G4double currentInteractionLength
Definition: G4VProcess.hh:297
void SetProposedKineticEnergy(G4double proposedKinEnergy)
G4double GetSubDEDXForScaledEnergy(G4double scaledKinEnergy)
const G4ParticleDefinition * GetParticleDefinition() const
G4PhysicsTable * theCSDARangeTable
static const double GeV
Definition: G4SIunits.hh:214
const G4String & GetParticleType() const
void SetMaxKinEnergy(G4double e)
void Register(G4VEnergyLossProcess *p)
std::vector< G4double > theDEDXAtMaxEnergy
virtual G4VParticleChange * AlongStepDoIt(const G4Track &, const G4Step &)
Definition: G4Step.hh:76
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *fluc=0, const G4Region *region=0)
G4double MinKinEnergy() const
const G4ParticleDefinition * particle
virtual void ProcessDescription(std::ostream &outFile) const
const G4int n
G4double GetGlobalTime() const
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
G4PhysicsTable * DEDXTableForSubsec() const
G4double Energy(size_t index) const
G4double Value(G4double theEnergy, size_t &lastidx) const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
virtual G4double AlongStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection)
void AddSecondary(G4Track *aSecondary)
virtual G4double ChargeSquareRatio(const G4Track &)
Definition: G4VEmModel.cc:337
const G4TouchableHandle & GetTouchableHandle() const
G4Material * GetMaterial() const
void RegisterExtraParticle(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:753
void DumpModelList(G4int verb)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
std::vector< const G4Region * > scoffRegions
const G4DataVector * Initialise(const G4ParticleDefinition *part, const G4ParticleDefinition *secPart, G4double minSubRange, G4int verb)
G4double ComputeSafety(const G4ThreeVector &pGlobalPoint, G4double maxRadius=DBL_MAX)
G4PhysicsTable * theLambdaTable
static G4TransportationManager * GetTransportationManager()
virtual G4bool RetrievePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii)
static const G4double emax
void SetMasterThread(G4bool val)
Definition: G4VEmModel.hh:711
virtual void PreparePhysicsTable(const G4ParticleDefinition &)
size_t length() const
G4double G4Log(G4double x)
Definition: G4Log.hh:230
const G4ParticleDefinition * theElectron
G4int size() const
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
static const double perMillion
Definition: G4SIunits.hh:331
virtual G4double GetParticleCharge(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:354
static G4ProductionCutsTable * GetProductionCutsTable()
const G4Material * currentMaterial
void SetLambdaTable(G4PhysicsTable *p)
G4PhysicsTable * theInverseRangeTable
static G4Positron * Positron()
Definition: G4Positron.cc:94
G4int NumberOfModels() const
G4SafetyHelper * safetyHelper
std::vector< G4double > theCrossSectionMax
G4bool LossFluctuation() const
static const G4double factor
virtual void PrintInfo()=0
G4PhysicsTable * InverseRangeTable() const
void BuildPhysicsTable(const G4ParticleDefinition *aParticle)
G4double GetPDGMass() const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
static G4ParticleTable * GetParticleTable()
virtual void StartTracking(G4Track *)
int G4lrint(double ad)
Definition: templates.hh:163
T max(const T t1, const T t2)
brief Return the largest of the two arguments
void InitialiseBaseMaterials(G4PhysicsTable *table)
void SetNumberOfSecondaries(G4int totSecondaries)
G4bool UseAngularGeneratorForIonisation() const
static G4int Register(const G4String &)
G4double GetLocalEnergyDeposit() const
G4double MeanFreePath(const G4Track &track)
G4PhysicsTable * theDEDXunRestrictedTable
G4StepPoint * GetPostStepPoint() const
void UpdateEmModel(const G4String &, G4double, G4double)
G4bool StoreTable(const G4ParticleDefinition *p, G4PhysicsTable *, G4bool ascii, const G4String &directory, const G4String &tname)
const G4double x[NPOINTSGL]
static G4EmParameters * Instance()
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
G4PhysicsTable * theRangeTableForLoss
void SetEmModel(G4VEmModel *, G4int index=1)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double ContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety)
std::vector< G4VEmModel * > emModels
G4double GetSafety() const
void SetProposedCharge(G4double theCharge)
void PreparePhysicsTable(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p, G4bool theMaster)
void SetCSDARangeTable(G4PhysicsTable *pRange)
G4double GetDEDXDispersion(const G4MaterialCutsCouple *couple, const G4DynamicParticle *dp, G4double length)
void LocalPhysicsTables(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
G4EmTableType
G4double GetLambdaForScaledEnergy(G4double scaledKinEnergy)
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
G4double GetGlobalTime() const
void SetSecondaryRangeTable(G4PhysicsTable *p)
static const double TeV
Definition: G4SIunits.hh:215
G4bool StorePhysicsTable(const G4String &filename, G4bool ascii=false)
G4TrackStatus GetTrackStatus() const
G4VAtomDeexcitation * AtomDeexcitation()
static const double keV
Definition: G4SIunits.hh:213
G4PhysicsTable * DEDXunRestrictedTable() const
std::vector< G4double > theRangeAtMaxEnergy
G4bool RetrieveTable(const G4ParticleDefinition *p, G4PhysicsTable *, G4bool ascii, const G4String &directory, const G4String &tname, G4bool mandatory)
G4PhysicsTable * theDEDXTable
G4double theInitialNumberOfInteractionLength
Definition: G4VProcess.hh:300
G4bool GetFlag(size_t i) const
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
void SetLowestEnergyLimit(G4double)
void SetSubLambdaTable(G4PhysicsTable *p)
virtual G4double MinPrimaryEnergy(const G4ParticleDefinition *, const G4Material *, G4double cut)
G4ProcessVector * GetAlongStepProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4ForceCondition
G4Track * GetTrack() const
G4double MaxEnergyForCSDARange() const
G4double GetPDGCharge() const
G4EmBiasingManager * biasManager
void ActivateSubCutoff(G4bool val, const G4Region *region=0)
void SetRangeTableForLoss(G4PhysicsTable *p)
static G4bool RetrievePhysicsTable(G4PhysicsTable *physTable, const G4String &fileName, G4bool ascii)
G4ProductionCuts * GetProductionCuts() const
G4PhysicsTable * theDEDXSubTable
#define DBL_MAX
Definition: templates.hh:83
G4EmParameters * theParameters
static const double mm
Definition: G4SIunits.hh:114
G4int GetProcessSubType() const
Definition: G4VProcess.hh:426
G4double GetScaledRangeForScaledEnergy(G4double scaledKinEnergy)
std::vector< G4VEnergyLossProcess * > scProcesses
G4PhysicsTable * BuildDEDXTable(G4EmTableType tType=fRestricted)
G4VEmModel * EmModel(G4int index=1) const
void DefineMaterial(const G4MaterialCutsCouple *couple)
void clearAndDestroy()
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:437
void SetDEDXBinning(G4int nbins)
virtual void InitialiseEnergyLossProcess(const G4ParticleDefinition *, const G4ParticleDefinition *)=0
G4GPILSelection
const G4Material * GetMaterial() const
G4PhysicsTable * theIonisationSubTable
G4double CrossSectionPerVolume(G4double kineticEnergy, const G4MaterialCutsCouple *couple)
std::vector< G4DynamicParticle * > secParticles
const std::vector< G4int > * GetCoupleIndexes()
G4ProcessType
const G4Element * GetCurrentElement() const
Definition: G4VEmModel.hh:467
G4bool IsIonisationProcess() const
void SetMinKinEnergy(G4double e)
G4EmModelManager * modelManager