Geant4  10.01.p03
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 88981 2015-03-17 10:14:15Z 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 
164 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
165 
167  G4ProcessType type):
168  G4VContinuousDiscreteProcess(name, type),
169  secondaryParticle(0),
170  nSCoffRegions(0),
171  idxSCoffRegions(0),
172  nProcesses(0),
173  theDEDXTable(0),
174  theDEDXSubTable(0),
175  theDEDXunRestrictedTable(0),
176  theIonisationTable(0),
177  theIonisationSubTable(0),
178  theRangeTableForLoss(0),
179  theCSDARangeTable(0),
180  theSecondaryRangeTable(0),
181  theInverseRangeTable(0),
182  theLambdaTable(0),
183  theSubLambdaTable(0),
184  theDensityFactor(0),
185  theDensityIdx(0),
186  baseParticle(0),
187  lossFluctuationFlag(true),
188  rndmStepFlag(false),
189  tablesAreBuilt(false),
190  integral(true),
191  isIon(false),
192  isIonisation(true),
193  useSubCutoff(false),
194  useDeexcitation(false),
195  particle(0),
196  currentCouple(0),
197  nWarnings(0),
198  mfpKinEnergy(0.0)
199 {
201  SetVerboseLevel(1);
202 
203  // low energy limit
204  lowestKinEnergy = 1.*eV;
205  preStepKinEnergy = 0.0;
206  preStepRangeEnergy = 0.0;
208 
209  // Size of tables assuming spline
210  minKinEnergy = 0.1*keV;
211  maxKinEnergy = 10.0*TeV;
212  nBins = 77;
213  maxKinEnergyCSDA = 1.0*GeV;
214  nBinsCSDA = 35;
216  = actLossFluc = false;
217 
218  // default linear loss limit for spline
219  linLossLimit = 0.01;
220 
221  // default dRoverRange and finalRange
222  SetStepFunction(0.2, 1.0*mm);
223 
224  // default lambda factor
225  lambdaFactor = 0.8;
226 
227  // cross section biasing
228  biasFactor = 1.0;
229 
230  // particle types
234  theGenericIon = 0;
235 
236  // run time objects
241  ->GetSafetyHelper();
243 
244  // initialise model
246  lManager->Register(this);
247  fluctModel = 0;
248  currentModel = 0;
249  atomDeexcitation = 0;
250  subcutProducer = 0;
251 
252  biasManager = 0;
253  biasFlag = false;
254  weightFlag = false;
255  isMaster = true;
256  lastIdx = 0;
257 
261 
262  scTracks.reserve(5);
263  secParticles.reserve(5);
264 
265  theCuts = theSubCuts = 0;
266  currentMaterial = 0;
270 
271  secID = biasID = subsecID = -1;
272 }
273 
274 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
275 
277 {
278  /*
279  G4cout << "** G4VEnergyLossProcess::~G4VEnergyLossProcess() for "
280  << GetProcessName() << " isMaster: " << isMaster
281  << " basePart: " << baseParticle
282  << G4endl;
283  */
284  Clean();
285 
286  // G4cout << " isIonisation " << isIonisation << " "
287  // << theDEDXTable << " " << theIonisationTable << G4endl;
288 
289  if (isMaster && !baseParticle) {
290  if(theDEDXTable) {
291 
292  //G4cout << " theIonisationTable " << theIonisationTable << G4endl;
294  //G4cout << " delete theDEDXTable " << theDEDXTable << G4endl;
296  delete theDEDXTable;
297  theDEDXTable = 0;
298  if(theDEDXSubTable) {
300  { theIonisationSubTable = 0; }
302  delete theDEDXSubTable;
303  theDEDXSubTable = 0;
304  }
305  }
306  //G4cout << " theIonisationTable " << theIonisationTable << G4endl;
307  if(theIonisationTable) {
308  //G4cout << " delete theIonisationTable " << theIonisationTable << G4endl;
310  delete theIonisationTable;
311  theIonisationTable = 0;
312  }
315  delete theIonisationSubTable;
317  }
322  }
325  delete theCSDARangeTable;
326  theCSDARangeTable = 0;
327  }
328  //G4cout << "delete RangeTable: " << theRangeTableForLoss << G4endl;
331  delete theRangeTableForLoss;
333  }
334  //G4cout << "delete InvRangeTable: " << theInverseRangeTable << G4endl;
335  if(theInverseRangeTable && isIonisation /*&& !isIon*/) {
337  delete theInverseRangeTable;
339  }
340  //G4cout << "delete LambdaTable: " << theLambdaTable << G4endl;
341  if(theLambdaTable) {
343  delete theLambdaTable;
344  theLambdaTable = 0;
345  }
346  if(theSubLambdaTable) {
348  delete theSubLambdaTable;
349  theSubLambdaTable = 0;
350  }
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 = 0;
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 = 0;
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  // Tables preparation
542  if (isMaster && !baseParticle) {
543 
544  if(theDEDXTable && isIonisation) {
547  delete theDEDXTable;
549  }
553  delete theDEDXSubTable;
555  }
556  }
557 
560 
561  if(theDEDXSubTable) {
562  theDEDXSubTable =
564  }
565 
566  if (lManager->BuildCSDARange()) {
571  }
572 
574 
575  if(isIonisation) {
580  }
581 
583  theDEDXSubTable =
587  }
588  }
589  /*
590  G4cout << "** G4VEnergyLossProcess::PreparePhysicsTable() for "
591  << GetProcessName() << " and " << particle->GetParticleName()
592  << " isMaster: " << isMaster << " isIonisation: "
593  << isIonisation << G4endl;
594  G4cout << " theDEDX: " << theDEDXTable
595  << " theRange: " << theRangeTableForLoss
596  << " theInverse: " << theInverseRangeTable
597  << " theLambda: " << theLambdaTable << G4endl;
598  */
599  // forced biasing
600  if(biasManager) {
602  biasFlag = false;
603  }
604 
605  G4double initialCharge = particle->GetPDGCharge();
606  G4double initialMass = particle->GetPDGMass();
607 
608  if (baseParticle) {
609  massRatio = (baseParticle->GetPDGMass())/initialMass;
610  G4double q = initialCharge/baseParticle->GetPDGCharge();
611  chargeSqRatio = q*q;
612  if(chargeSqRatio > 0.0) { reduceFactor = 1.0/(chargeSqRatio*massRatio); }
613  }
614 
615  // defined ID of secondary particles
616  if(isMaster) {
617  G4String nam1 = GetProcessName();
618  G4String nam4 = nam1 + "_split";
619  G4String nam5 = nam1 + "_subcut";
623  }
624 
625  // initialisation of models
627  for(G4int i=0; i<nmod; ++i) {
628  G4VEmModel* mod = modelManager->GetModel(i);
632  if(mod->HighEnergyLimit() > maxKinEnergy) {
634  }
635  }
638  verboseLevel);
639 
640  // Sub Cutoff
641  if(nSCoffRegions > 0) {
642  if(theParameters->MinSubRange() < 1.0) { useSubCutoff = true; }
643 
645 
646  idxSCoffRegions = new G4bool[n];
647  for (size_t j=0; j<n; ++j) {
648 
649  const G4MaterialCutsCouple* couple =
650  theCoupleTable->GetMaterialCutsCouple(j);
651  const G4ProductionCuts* pcuts = couple->GetProductionCuts();
652 
653  G4bool reg = false;
654  for(G4int i=0; i<nSCoffRegions; ++i) {
655  if( pcuts == scoffRegions[i]->GetProductionCuts()) {
656  reg = true;
657  break;
658  }
659  }
660  idxSCoffRegions[j] = reg;
661  }
662  }
663 
664  if(1 < verboseLevel) {
665  G4cout << "G4VEnergyLossProcess::PrepearPhysicsTable() is done "
666  << " for local " << particle->GetParticleName()
667  << " isIon= " << isIon;
668  if(baseParticle) { G4cout << "; base: " << baseParticle->GetParticleName(); }
669  G4cout << " chargeSqRatio= " << chargeSqRatio
670  << " massRatio= " << massRatio
671  << " reduceFactor= " << reduceFactor << G4endl;
672  if (nSCoffRegions) {
673  G4cout << " SubCutoff Regime is ON for regions: " << G4endl;
674  for (G4int i=0; i<nSCoffRegions; ++i) {
675  const G4Region* r = scoffRegions[i];
676  G4cout << " " << r->GetName() << G4endl;
677  }
678  }
679  }
680 }
681 
682 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
683 
685 {
686  G4bool verb = false;
687  if(1 < verboseLevel || verb) {
688 
689  G4cout << "### G4VEnergyLossProcess::BuildPhysicsTable() for "
690  << GetProcessName()
691  << " and particle " << part.GetParticleName()
692  << "; local: " << particle->GetParticleName();
693  if(baseParticle) {
694  G4cout << "; base: " << baseParticle->GetParticleName();
695  }
696  G4cout << " TablesAreBuilt= " << tablesAreBuilt
697  << " isIon= " << isIon << " " << this << G4endl;
698  }
699 
700  if(&part == particle) {
701 
703  if(isMaster) {
707 
708  } else {
709 
710  const G4VEnergyLossProcess* masterProcess =
711  static_cast<const G4VEnergyLossProcess*>(GetMasterProcess());
712 
713  // define density factors for worker thread
714  bld->InitialiseBaseMaterials(masterProcess->DEDXTable());
717 
718  // copy table pointers from master thread
719  SetDEDXTable(masterProcess->DEDXTable(),fRestricted);
721  SetDEDXTable(masterProcess->DEDXunRestrictedTable(),fTotal);
722  SetDEDXTable(masterProcess->IonisationTable(),fIsIonisation);
724  SetRangeTableForLoss(masterProcess->RangeTableForLoss());
725  SetCSDARangeTable(masterProcess->CSDARangeTable());
727  SetInverseRangeTable(masterProcess->InverseRangeTable());
728  SetLambdaTable(masterProcess->LambdaTable());
729  SetSubLambdaTable(masterProcess->SubLambdaTable());
730  isIonisation = masterProcess->IsIonisationProcess();
731 
732  tablesAreBuilt = true;
733  // local initialisation of models
734  G4bool printing = true;
735  G4int numberOfModels = modelManager->NumberOfModels();
736  for(G4int i=0; i<numberOfModels; ++i) {
737  G4VEmModel* mod = GetModelByIndex(i, printing);
738  G4VEmModel* mod0= masterProcess->GetModelByIndex(i,printing);
739  mod->InitialiseLocal(particle, mod0);
740  }
741 
743  }
744 
745  // needs to be done only once
747  }
748  // explicitly defined printout by particle name
749  G4String num = part.GetParticleName();
750  if(1 < verboseLevel ||
751  (0 < verboseLevel && (num == "e-" ||
752  num == "e+" || num == "mu+" ||
753  num == "mu-" || num == "proton"||
754  num == "pi+" || num == "pi-" ||
755  num == "kaon+" || num == "kaon-" ||
756  num == "alpha" || num == "anti_proton" ||
757  num == "GenericIon")))
758  {
759  PrintInfoDefinition(part);
760  }
761 
762  // Added tracking cut to avoid tracking artifacts
763  // identify deexcitation flag
764  if(isIonisation) {
768  if(atomDeexcitation) {
769  if(atomDeexcitation->IsPIXEActive()) { useDeexcitation = true; }
770  }
771  }
772  /*
773  G4cout << "** G4VEnergyLossProcess::BuildPhysicsTable() for "
774  << GetProcessName() << " and " << particle->GetParticleName()
775  << " isMaster: " << isMaster << " isIonisation: "
776  << isIonisation << G4endl;
777  G4cout << " theDEDX: " << theDEDXTable
778  << " theRange: " << theRangeTableForLoss
779  << " theInverse: " << theInverseRangeTable
780  << " theLambda: " << theLambdaTable << G4endl;
781  */
782  //if(1 < verboseLevel || verb) {
783  if(1 < verboseLevel) {
784  G4cout << "### G4VEnergyLossProcess::BuildPhysicsTable() done for "
785  << GetProcessName()
786  << " and particle " << part.GetParticleName();
787  if(isIonisation) { G4cout << " isIonisation flag = 1"; }
788  G4cout << G4endl;
789  }
790 }
791 
792 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
793 
795 {
796  G4bool verb = false;
797  if(1 < verboseLevel || verb) {
798  G4cout << "G4VEnergyLossProcess::BuildDEDXTable() of type " << tType
799  << " for " << GetProcessName()
800  << " and particle " << particle->GetParticleName()
801  << G4endl;
802  }
803  G4PhysicsTable* table = 0;
804  G4double emax = maxKinEnergy;
805  G4int bin = nBins;
806 
807  if(fTotal == tType) {
808  emax = maxKinEnergyCSDA;
809  bin = nBinsCSDA;
810  table = theDEDXunRestrictedTable;
811  } else if(fRestricted == tType) {
812  table = theDEDXTable;
813  } else if(fSubRestricted == tType) {
814  table = theDEDXSubTable;
815  } else {
816  G4cout << "G4VEnergyLossProcess::BuildDEDXTable WARNING: wrong type "
817  << tType << G4endl;
818  }
819 
820  // Access to materials
821  const G4ProductionCutsTable* theCoupleTable=
823  size_t numOfCouples = theCoupleTable->GetTableSize();
824 
825  if(1 < verboseLevel || verb) {
826  G4cout << numOfCouples << " materials"
827  << " minKinEnergy= " << minKinEnergy
828  << " maxKinEnergy= " << emax
829  << " nbin= " << bin
830  << " EmTableType= " << tType
831  << " table= " << table << " " << this
832  << G4endl;
833  }
834  if(!table) { return table; }
835 
837  G4bool splineFlag = theParameters->Spline();
838  G4PhysicsLogVector* aVector = 0;
839  G4PhysicsLogVector* bVector = 0;
840 
841  for(size_t i=0; i<numOfCouples; ++i) {
842 
843  if(1 < verboseLevel || verb) {
844  G4cout << "G4VEnergyLossProcess::BuildDEDXVector Idx= " << i
845  << " flagTable= " << table->GetFlag(i)
846  << " Flag= " << bld->GetFlag(i) << G4endl;
847  }
848  if(bld->GetFlag(i)) {
849 
850  // create physics vector and fill it
851  const G4MaterialCutsCouple* couple =
852  theCoupleTable->GetMaterialCutsCouple(i);
853  if((*table)[i]) { delete (*table)[i]; }
854  if(bVector) {
855  aVector = new G4PhysicsLogVector(*bVector);
856  } else {
857  bVector = new G4PhysicsLogVector(minKinEnergy, emax, bin);
858  aVector = bVector;
859  }
860  aVector->SetSpline(splineFlag);
861 
862  modelManager->FillDEDXVector(aVector, couple, tType);
863  if(splineFlag) { aVector->FillSecondDerivatives(); }
864 
865  // Insert vector for this material into the table
866  G4PhysicsTableHelper::SetPhysicsVector(table, i, aVector);
867  }
868  }
869 
870  if(1 < verboseLevel || verb) {
871  G4cout << "G4VEnergyLossProcess::BuildDEDXTable(): table is built for "
873  << " and process " << GetProcessName()
874  << G4endl;
875  //if(2 < verboseLevel) G4cout << (*table) << G4endl;
876  }
877 
878  return table;
879 }
880 
881 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
882 
884 {
885  G4PhysicsTable* table = 0;
886 
887  if(fRestricted == tType) {
888  table = theLambdaTable;
889  } else if(fSubRestricted == tType) {
890  table = theSubLambdaTable;
891  } else {
892  G4cout << "G4VEnergyLossProcess::BuildLambdaTable WARNING: wrong type "
893  << tType << G4endl;
894  }
895 
896  if(1 < verboseLevel) {
897  G4cout << "G4VEnergyLossProcess::BuildLambdaTable() of type "
898  << tType << " for process "
899  << GetProcessName() << " and particle "
901  << " EmTableType= " << tType
902  << " table= " << table
903  << G4endl;
904  }
905  if(!table) {return table;}
906 
907  // Access to materials
908  const G4ProductionCutsTable* theCoupleTable=
910  size_t numOfCouples = theCoupleTable->GetTableSize();
911 
915 
916  G4bool splineFlag = theParameters->Spline();
917  G4PhysicsLogVector* aVector = 0;
919 
920  for(size_t i=0; i<numOfCouples; ++i) {
921 
922  if (bld->GetFlag(i)) {
923 
924  // create physics vector and fill it
925  const G4MaterialCutsCouple* couple =
926  theCoupleTable->GetMaterialCutsCouple(i);
927  delete (*table)[i];
928 
929  G4bool startNull = true;
930  G4double emin =
931  MinPrimaryEnergy(particle,couple->GetMaterial(),(*theCuts)[i]);
932  if(minKinEnergy > emin) {
933  emin = minKinEnergy;
934  startNull = false;
935  }
936 
937  G4double emax = maxKinEnergy;
938  if(emax <= emin) { emax = 2*emin; }
939  G4int bin = G4lrint(nBins*G4Log(emax/emin)/scale);
940  bin = std::max(bin, 3);
941  aVector = new G4PhysicsLogVector(emin, emax, bin);
942  aVector->SetSpline(splineFlag);
943 
944  modelManager->FillLambdaVector(aVector, couple, startNull, tType);
945  if(splineFlag) { aVector->FillSecondDerivatives(); }
946 
947  // Insert vector for this material into the table
948  G4PhysicsTableHelper::SetPhysicsVector(table, i, aVector);
949  }
950  }
951 
952  if(1 < verboseLevel) {
953  G4cout << "Lambda table is built for "
955  << G4endl;
956  }
957 
958  return table;
959 }
960 
961 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
962 
963 void
965 {
966  if(0 < verboseLevel) {
967  G4cout << std::setprecision(6);
968  G4cout << G4endl << GetProcessName() << ": for "
969  << part.GetParticleName()
970  << " SubType= " << GetProcessSubType()
971  << G4endl;
972  G4cout << " dE/dx and range tables from "
973  << G4BestUnit(minKinEnergy,"Energy")
974  << " to " << G4BestUnit(maxKinEnergy,"Energy")
975  << " in " << nBins << " bins" << G4endl
976  << " Lambda tables from threshold to "
977  << G4BestUnit(maxKinEnergy,"Energy")
978  << ", " << theParameters->NumberOfBinsPerDecade()
979  << " bins per decade, spline: "
980  << theParameters->Spline()
981  << G4endl;
983  G4cout << " finalRange(mm)= " << finalRange/mm
984  << ", dRoverRange= " << dRoverRange
985  << ", integral: " << integral
986  << ", fluct: " << lossFluctuationFlag
987  << ", linLossLimit= " << linLossLimit
988  << G4endl;
989  }
990  PrintInfo();
993  G4cout << " CSDA range table up"
994  << " to " << G4BestUnit(maxKinEnergyCSDA,"Energy")
995  << " in " << nBinsCSDA << " bins" << G4endl;
996  }
997  if(nSCoffRegions>0 && isIonisation) {
998  G4cout << " Subcutoff sampling in " << nSCoffRegions
999  << " regions" << G4endl;
1000  }
1001  if(2 < verboseLevel) {
1002  G4cout << " DEDXTable address= " << theDEDXTable << G4endl;
1003  if(theDEDXTable && isIonisation) G4cout << (*theDEDXTable) << G4endl;
1004  G4cout << "non restricted DEDXTable address= "
1007  G4cout << (*theDEDXunRestrictedTable) << G4endl;
1008  }
1009  if(theDEDXSubTable && isIonisation) {
1010  G4cout << (*theDEDXSubTable) << G4endl;
1011  }
1012  G4cout << " CSDARangeTable address= " << theCSDARangeTable
1013  << G4endl;
1015  G4cout << (*theCSDARangeTable) << G4endl;
1016  }
1017  G4cout << " RangeTableForLoss address= " << theRangeTableForLoss
1018  << G4endl;
1020  G4cout << (*theRangeTableForLoss) << G4endl;
1021  }
1022  G4cout << " InverseRangeTable address= " << theInverseRangeTable
1023  << G4endl;
1025  G4cout << (*theInverseRangeTable) << G4endl;
1026  }
1027  G4cout << " LambdaTable address= " << theLambdaTable << G4endl;
1028  if(theLambdaTable && isIonisation) {
1029  G4cout << (*theLambdaTable) << G4endl;
1030  }
1031  G4cout << " SubLambdaTable address= " << theSubLambdaTable
1032  << G4endl;
1034  G4cout << (*theSubLambdaTable) << G4endl;
1035  }
1036  }
1037  }
1038 }
1039 
1040 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1041 
1043 {
1044  G4RegionStore* regionStore = G4RegionStore::GetInstance();
1045  const G4Region* reg = r;
1046  if (!reg) {
1047  reg = regionStore->GetRegion("DefaultRegionForTheWorld", false);
1048  }
1049 
1050  // the region is in the list
1051  if (nSCoffRegions > 0) {
1052  for (G4int i=0; i<nSCoffRegions; ++i) {
1053  if (reg == scoffRegions[i]) {
1054  return;
1055  }
1056  }
1057  }
1058 
1059  // new region
1060  if(val) {
1061  scoffRegions.push_back(reg);
1062  ++nSCoffRegions;
1063  }
1064 }
1065 
1066 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1067 
1069 {
1070  /*
1071  G4cout << track->GetDefinition()->GetParticleName()
1072  << " e(MeV)= " << track->GetKineticEnergy()
1073  << " baseParticle " << baseParticle << " proc " << this;
1074  if(particle) G4cout << " " << particle->GetParticleName();
1075  G4cout << " isIon= " << isIon << " dedx " << theDEDXTable <<G4endl;
1076  */
1077  // reset parameters for the new track
1080  preStepRangeEnergy = 0.0;
1081 
1082  // reset ion
1083  if(isIon) {
1084  chargeSqRatio = 0.5;
1085 
1086  G4double newmass = track->GetDefinition()->GetPDGMass();
1087  if(baseParticle) {
1088  massRatio = baseParticle->GetPDGMass()/newmass;
1089  } else if(theGenericIon) {
1090  massRatio = proton_mass_c2/newmass;
1091  } else {
1092  massRatio = 1.0;
1093  }
1094  }
1095  // forced biasing only for primary particles
1096  if(biasManager) {
1097  if(0 == track->GetParentID()) {
1098  // primary particle
1099  biasFlag = true;
1101  }
1102  }
1103 }
1104 
1105 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1106 
1109  G4GPILSelection* selection)
1110 {
1111  G4double x = DBL_MAX;
1112  *selection = aGPILSelection;
1115  x = fRange;
1116  G4double finR = finalRange;
1117  if(rndmStepFlag) {
1118  finR = std::min(finR,
1120  }
1121  if(fRange > finR) {
1122  x = fRange*dRoverRange + finR*(1.0 - dRoverRange)*(2.0 - finR/fRange);
1123  }
1124 
1125  // if(particle->GetPDGMass() > 0.9*GeV)
1126  /*
1127  G4cout<<GetProcessName()<<": e= "<<preStepKinEnergy
1128  <<" range= "<<fRange << " idx= " << basedCoupleIndex
1129  << " finR= " << finR
1130  << " limit= " << x <<G4endl;
1131  G4cout << "massRatio= " << massRatio << " Q^2= " << chargeSqRatio
1132  << " finR= " << finR << " dRoverRange= " << dRoverRange
1133  << " finalRange= " << finalRange << G4endl;
1134  */
1135  }
1136  //G4cout<<GetProcessName()<<": e= "<<preStepKinEnergy
1137  //<<" stepLimit= "<<x<<G4endl;
1138  return x;
1139 }
1140 
1141 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1142 
1144  const G4Track& track,
1145  G4double previousStepSize,
1147 {
1148  // condition is set to "Not Forced"
1149  *condition = NotForced;
1150  G4double x = DBL_MAX;
1151 
1152  // initialisation of material, mass, charge, model
1153  // at the beginning of the step
1156  preStepScaledEnergy = preStepKinEnergy*massRatio;
1158 
1162  return x;
1163  }
1164 
1165  // change effective charge of an ion on fly
1166  if(isIon) {
1168  if(q2 != chargeSqRatio && q2 > 0.0) {
1169  chargeSqRatio = q2;
1170  fFactor = q2*biasFactor*(*theDensityFactor)[currentCoupleIndex];
1171  reduceFactor = 1.0/(fFactor*massRatio);
1172  }
1173  }
1174  // if(particle->GetPDGMass() > 0.9*GeV)
1175  //G4cout << "q2= "<<chargeSqRatio << " massRatio= " << massRatio << G4endl;
1176 
1177  // forced biasing only for primary particles
1178  if(biasManager) {
1179  if(0 == track.GetParentID()) {
1180  if(biasFlag &&
1182  return biasManager->GetStepLimit(currentCoupleIndex, previousStepSize);
1183  }
1184  }
1185  }
1186 
1187  // compute mean free path
1191 
1192  // zero cross section
1193  if(preStepLambda <= 0.0) {
1196  }
1197  }
1198 
1199  // non-zero cross section
1200  if(preStepLambda > 0.0) {
1202 
1203  // beggining of tracking (or just after DoIt of this process)
1206 
1207  } else if(currentInteractionLength < DBL_MAX) {
1208 
1209  // subtract NumberOfInteractionLengthLeft using previous step
1211  previousStepSize/currentInteractionLength;
1212 
1215  }
1216 
1217  // new mean free path and step limit
1220  }
1221 #ifdef G4VERBOSE
1222  if (verboseLevel>2){
1223  // if(particle->GetPDGMass() > 0.9*GeV){
1224  G4cout << "G4VEnergyLossProcess::PostStepGetPhysicalInteractionLength ";
1225  G4cout << "[ " << GetProcessName() << "]" << G4endl;
1226  G4cout << " for " << track.GetDefinition()->GetParticleName()
1227  << " in Material " << currentMaterial->GetName()
1228  << " Ekin(MeV)= " << preStepKinEnergy/MeV
1229  << " " << track.GetMaterial()->GetName()
1230  <<G4endl;
1231  G4cout << "MeanFreePath = " << currentInteractionLength/cm << "[cm]"
1232  << "InteractionLength= " << x/cm <<"[cm] " <<G4endl;
1233  }
1234 #endif
1235  return x;
1236 }
1237 
1238 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1239 
1241  const G4Step& step)
1242 {
1244  // The process has range table - calculate energy loss
1246  return &fParticleChange;
1247  }
1248 
1249  // Get the actual (true) Step length
1250  G4double length = step.GetStepLength();
1251  if(length <= 0.0) { return &fParticleChange; }
1252  G4double eloss = 0.0;
1253 
1254  /*
1255  if(-1 < verboseLevel) {
1256  const G4ParticleDefinition* d = track.GetParticleDefinition();
1257  G4cout << "AlongStepDoIt for "
1258  << GetProcessName() << " and particle "
1259  << d->GetParticleName()
1260  << " eScaled(MeV)= " << preStepScaledEnergy/MeV
1261  << " range(mm)= " << fRange/mm
1262  << " s(mm)= " << length/mm
1263  << " rf= " << reduceFactor
1264  << " q^2= " << chargeSqRatio
1265  << " md= " << d->GetPDGMass()
1266  << " status= " << track.GetTrackStatus()
1267  << " " << track.GetMaterial()->GetName()
1268  << G4endl;
1269  }
1270  */
1271 
1272  const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
1273 
1274  // define new weight for primary and secondaries
1276  if(weightFlag) {
1277  weight /= biasFactor;
1279  }
1280 
1281  // stopping
1282  if (length >= fRange) {
1283  eloss = preStepKinEnergy;
1284  if (useDeexcitation) {
1286  eloss, currentCoupleIndex);
1287  if(scTracks.size() > 0) { FillSecondariesAlongStep(eloss, weight); }
1288  if(eloss < 0.0) { eloss = 0.0; }
1289  }
1292  return &fParticleChange;
1293  }
1294  //G4cout << theDEDXTable << " idx= " << basedCoupleIndex
1295  // << " " << GetProcessName() << " "<< currentMaterial->GetName()<<G4endl;
1296  //if(particle->GetParticleName() == "e-")G4cout << (*theDEDXTable) <<G4endl;
1297  // Short step
1299 
1300  //G4cout << "eloss= " << eloss << G4endl;
1301 
1302  // Long step
1303  if(eloss > preStepKinEnergy*linLossLimit) {
1304 
1305  G4double x = (fRange - length)/reduceFactor;
1306  //G4cout << "x= " << x << " " << theInverseRangeTable << G4endl;
1308 
1309  /*
1310  if(-1 < verboseLevel)
1311  G4cout << "Long STEP: rPre(mm)= "
1312  << GetScaledRangeForScaledEnergy(preStepScaledEnergy)/mm
1313  << " rPost(mm)= " << x/mm
1314  << " ePre(MeV)= " << preStepScaledEnergy/MeV
1315  << " eloss(MeV)= " << eloss/MeV
1316  << " eloss0(MeV)= "
1317  << GetDEDXForScaledEnergy(preStepScaledEnergy)*length/MeV
1318  << " lim(MeV)= " << preStepKinEnergy*linLossLimit/MeV
1319  << G4endl;
1320  */
1321  }
1322 
1323  /*
1324  G4double eloss0 = eloss;
1325  if(-1 < verboseLevel ) {
1326  G4cout << "Before fluct: eloss(MeV)= " << eloss/MeV
1327  << " e-eloss= " << preStepKinEnergy-eloss
1328  << " step(mm)= " << length/mm
1329  << " range(mm)= " << fRange/mm
1330  << " fluct= " << lossFluctuationFlag
1331  << G4endl;
1332  }
1333  */
1334 
1335  G4double cut = (*theCuts)[currentCoupleIndex];
1336  G4double esec = 0.0;
1337 
1338  //G4cout << "cut= " << cut << " useSubCut= " << useSubCutoff << G4endl;
1339 
1340  // SubCutOff
1341  if(useSubCutoff && !subcutProducer) {
1343 
1344  G4bool yes = false;
1345  G4StepPoint* prePoint = step.GetPreStepPoint();
1346 
1347  // Check boundary
1348  if(prePoint->GetStepStatus() == fGeomBoundary) { yes = true; }
1349 
1350  // Check PrePoint
1351  else {
1352  G4double preSafety = prePoint->GetSafety();
1353  G4double rcut =
1355 
1356  // recompute presafety
1357  if(preSafety < rcut) {
1358  preSafety = safetyHelper->ComputeSafety(prePoint->GetPosition(),
1359  rcut);
1360  }
1361 
1362  if(preSafety < rcut) { yes = true; }
1363 
1364  // Check PostPoint
1365  else {
1366  G4double postSafety = preSafety - length;
1367  if(postSafety < rcut) {
1368  postSafety = safetyHelper->ComputeSafety(
1369  step.GetPostStepPoint()->GetPosition(), rcut);
1370  if(postSafety < rcut) { yes = true; }
1371  }
1372  }
1373  }
1374 
1375  // Decided to start subcut sampling
1376  if(yes) {
1377 
1378  cut = (*theSubCuts)[currentCoupleIndex];
1380  esec = SampleSubCutSecondaries(scTracks, step,
1381  currentModel,currentCoupleIndex);
1382  // add bremsstrahlung sampling
1383  /*
1384  if(nProcesses > 0) {
1385  for(G4int i=0; i<nProcesses; ++i) {
1386  (scProcesses[i])->SampleSubCutSecondaries(
1387  scTracks, step, (scProcesses[i])->
1388  SelectModelForMaterial(preStepKinEnergy, currentCoupleIndex),
1389  currentCoupleIndex);
1390  }
1391  }
1392  */
1393  }
1394  }
1395  }
1396 
1397  // Corrections, which cannot be tabulated
1398  if(isIon) {
1399  G4double eadd = 0.0;
1400  G4double eloss_before = eloss;
1402  eloss, eadd, length);
1403  if(eloss < 0.0) { eloss = 0.5*eloss_before; }
1404  }
1405 
1406  // Sample fluctuations
1407  if (lossFluctuationFlag) {
1409  if(fluc &&
1410  (eloss + esec + lowestKinEnergy) < preStepKinEnergy) {
1411 
1412  G4double tmax =
1413  std::min(currentModel->MaxSecondaryKinEnergy(dynParticle),cut);
1414  eloss = fluc->SampleFluctuations(currentCouple,dynParticle,
1415  tmax,length,eloss);
1416  /*
1417  if(-1 < verboseLevel)
1418  G4cout << "After fluct: eloss(MeV)= " << eloss/MeV
1419  << " fluc= " << (eloss-eloss0)/MeV
1420  << " ChargeSqRatio= " << chargeSqRatio
1421  << " massRatio= " << massRatio
1422  << " tmax= " << tmax
1423  << G4endl;
1424  */
1425  }
1426  }
1427 
1428  // deexcitation
1429  if (useDeexcitation) {
1430  G4double esecfluo = preStepKinEnergy - esec;
1431  G4double de = esecfluo;
1432  //G4double eloss0 = eloss;
1433  /*
1434  G4cout << "### 1: E(keV)= " << preStepKinEnergy/keV
1435  << " Efluomax(keV)= " << de/keV
1436  << " Eloss(keV)= " << eloss/keV << G4endl;
1437  */
1439  de, currentCoupleIndex);
1440 
1441  // sum of de-excitation energies
1442  esecfluo -= de;
1443 
1444  // subtracted from energy loss
1445  if(eloss >= esecfluo) {
1446  esec += esecfluo;
1447  eloss -= esecfluo;
1448  } else {
1449  esec += esecfluo;
1450  eloss = 0.0;
1451  }
1452  /*
1453  if(esecfluo > 0.0) {
1454  G4cout << "### 2: E(keV)= " << preStepKinEnergy/keV
1455  << " Esec(keV)= " << esec/keV
1456  << " Esecf(kV)= " << esecfluo/keV
1457  << " Eloss0(kV)= " << eloss0/keV
1458  << " Eloss(keV)= " << eloss/keV
1459  << G4endl;
1460  }
1461  */
1462  }
1464  subcutProducer->SampleSecondaries(step, scTracks, eloss, cut);
1465  }
1466  if(scTracks.size() > 0) { FillSecondariesAlongStep(eloss, weight); }
1467 
1468  // Energy balanse
1469  G4double finalT = preStepKinEnergy - eloss - esec;
1470  if (finalT <= lowestKinEnergy) {
1471  eloss += finalT;
1472  finalT = 0.0;
1473  } else if(isIon) {
1476  currentMaterial,finalT));
1477  }
1478 
1479  eloss = std::max(eloss, 0.0);
1480 
1483  /*
1484  if(-1 < verboseLevel) {
1485  G4double del = finalT + eloss + esec - preStepKinEnergy;
1486  G4cout << "Final value eloss(MeV)= " << eloss/MeV
1487  << " preStepKinEnergy= " << preStepKinEnergy
1488  << " postStepKinEnergy= " << finalT
1489  << " de(keV)= " << del/keV
1490  << " lossFlag= " << lossFluctuationFlag
1491  << " status= " << track.GetTrackStatus()
1492  << G4endl;
1493  }
1494  */
1495  return &fParticleChange;
1496 }
1497 
1498 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1499 
1500 void
1502 {
1503  G4int n0 = scTracks.size();
1504 
1505  // weight may be changed by biasing manager
1506  if(biasManager) {
1508  weight *=
1510  }
1511  }
1512 
1513  // fill secondaries
1514  G4int n = scTracks.size();
1516 
1517  for(G4int i=0; i<n; ++i) {
1518  G4Track* t = scTracks[i];
1519  if(t) {
1520  t->SetWeight(weight);
1522  if(i >= n0) { t->SetCreatorModelIndex(biasID); }
1523  //G4cout << "Secondary(along step) has weight " << t->GetWeight()
1524  //<< ", kenergy " << t->GetKineticEnergy()/MeV << " MeV" <<G4endl;
1525  }
1526  }
1527  scTracks.clear();
1528 }
1529 
1530 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1531 
1532 G4double
1533 G4VEnergyLossProcess::SampleSubCutSecondaries(std::vector<G4Track*>& tracks,
1534  const G4Step& step,
1535  G4VEmModel* model,
1536  G4int idx)
1537 {
1538  // Fast check weather subcutoff can work
1539  G4double esec = 0.0;
1540  G4double subcut = (*theSubCuts)[idx];
1541  G4double cut = (*theCuts)[idx];
1542  if(cut <= subcut) { return esec; }
1543 
1544  const G4Track* track = step.GetTrack();
1545  const G4DynamicParticle* dp = track->GetDynamicParticle();
1547  G4double cross = (*theDensityFactor)[idx]*chargeSqRatio
1548  *(((*theSubLambdaTable)[(*theDensityIdx)[idx]])->Value(e, idxSubLambda));
1549  G4double length = step.GetStepLength();
1550 
1551  // negligible probability to get any interaction
1552  if(length*cross < perMillion) { return esec; }
1553  /*
1554  if(-1 < verboseLevel)
1555  G4cout << "<<< Subcutoff for " << GetProcessName()
1556  << " cross(1/mm)= " << cross*mm << ">>>"
1557  << " e(MeV)= " << preStepScaledEnergy
1558  << " matIdx= " << currentCoupleIndex
1559  << G4endl;
1560  */
1561 
1562  // Sample subcutoff secondaries
1563  G4StepPoint* preStepPoint = step.GetPreStepPoint();
1564  G4StepPoint* postStepPoint = step.GetPostStepPoint();
1565  G4ThreeVector prepoint = preStepPoint->GetPosition();
1566  G4ThreeVector dr = postStepPoint->GetPosition() - prepoint;
1567  G4double pretime = preStepPoint->GetGlobalTime();
1568  G4double dt = postStepPoint->GetGlobalTime() - pretime;
1569  //G4double dt = length/preStepPoint->GetVelocity();
1570  G4double fragment = 0.0;
1571 
1572  do {
1573  G4double del = -G4Log(G4UniformRand())/cross;
1574  fragment += del/length;
1575  if (fragment > 1.0) { break; }
1576 
1577  // sample secondaries
1578  secParticles.clear();
1580  dp,subcut,cut);
1581 
1582  // position of subcutoff particles
1583  G4ThreeVector r = prepoint + fragment*dr;
1584  std::vector<G4DynamicParticle*>::iterator it;
1585  for(it=secParticles.begin(); it!=secParticles.end(); ++it) {
1586 
1587  G4Track* t = new G4Track((*it), pretime + fragment*dt, r);
1590  tracks.push_back(t);
1591  esec += t->GetKineticEnergy();
1592  if (t->GetParticleDefinition() == thePositron) {
1593  esec += 2.0*electron_mass_c2;
1594  }
1595 
1596  /*
1597  if(-1 < verboseLevel)
1598  G4cout << "New track "
1599  << t->GetParticleDefinition()->GetParticleName()
1600  << " e(keV)= " << t->GetKineticEnergy()/keV
1601  << " fragment= " << fragment
1602  << G4endl;
1603  */
1604  }
1605  } while (fragment <= 1.0);
1606  return esec;
1607 }
1608 
1609 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1610 
1612  const G4Step& step)
1613 {
1614  // In all cases clear number of interaction lengths
1617 
1619  G4double finalT = track.GetKineticEnergy();
1620  if(finalT <= lowestKinEnergy) { return &fParticleChange; }
1621 
1622  G4double postStepScaledEnergy = finalT*massRatio;
1623 
1624  if(!currentModel->IsActive(postStepScaledEnergy)) {
1625  return &fParticleChange;
1626  }
1627  /*
1628  if(-1 < verboseLevel) {
1629  G4cout << GetProcessName()
1630  << "::PostStepDoIt: E(MeV)= " << finalT/MeV
1631  << G4endl;
1632  }
1633  */
1634 
1635  // forced process - should happen only once per track
1636  if(biasFlag) {
1638  biasFlag = false;
1639  }
1640  }
1641 
1642  // Integral approach
1643  if (integral) {
1644  G4double lx = GetLambdaForScaledEnergy(postStepScaledEnergy);
1645  /*
1646  if(preStepLambda<lx && 1 < verboseLevel && nWarnings<200) {
1647  G4cout << "WARNING: for " << particle->GetParticleName()
1648  << " and " << GetProcessName()
1649  << " E(MeV)= " << finalT/MeV
1650  << " preLambda= " << preStepLambda
1651  << " < " << lx << " (postLambda) "
1652  << G4endl;
1653  ++nWarnings;
1654  }
1655  */
1656  if(lx <= 0.0 || preStepLambda*G4UniformRand() > lx) {
1657  return &fParticleChange;
1658  }
1659  }
1660 
1661  SelectModel(postStepScaledEnergy);
1662 
1663  // define new weight for primary and secondaries
1665  if(weightFlag) {
1666  weight /= biasFactor;
1668  }
1669 
1670  const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
1671  G4double tcut = (*theCuts)[currentCoupleIndex];
1672 
1673  // sample secondaries
1674  secParticles.clear();
1675  //G4cout<< "@@@ Eprimary= "<<dynParticle->GetKineticEnergy()/MeV
1676  // << " cut= " << tcut/MeV << G4endl;
1678  dynParticle, tcut);
1679 
1680  G4int num0 = secParticles.size();
1681 
1682  // bremsstrahlung splitting or Russian roulette
1683  if(biasManager) {
1685  G4double eloss = 0.0;
1687  secParticles,
1688  track, currentModel,
1689  &fParticleChange, eloss,
1690  currentCoupleIndex, tcut,
1691  step.GetPostStepPoint()->GetSafety());
1692  if(eloss > 0.0) {
1695  }
1696  }
1697  }
1698 
1699  // save secondaries
1700  G4int num = secParticles.size();
1701  if(num > 0) {
1702 
1704  G4double time = track.GetGlobalTime();
1705 
1706  for (G4int i=0; i<num; ++i) {
1707  if(secParticles[i]) {
1708  G4Track* t = new G4Track(secParticles[i], time, track.GetPosition());
1710  t->SetWeight(weight);
1711  if(i < num0) { t->SetCreatorModelIndex(secID); }
1712  else { t->SetCreatorModelIndex(biasID); }
1713 
1714  //G4cout << "Secondary(post step) has weight " << t->GetWeight()
1715  // << ", kenergy " << t->GetKineticEnergy()/MeV << " MeV"
1716  // << " time= " << time/ns << " ns " << G4endl;
1718  }
1719  }
1720  }
1721 
1727  }
1728 
1729  /*
1730  if(-1 < verboseLevel) {
1731  G4cout << "::PostStepDoIt: Sample secondary; Efin= "
1732  << fParticleChange.GetProposedKineticEnergy()/MeV
1733  << " MeV; model= (" << currentModel->LowEnergyLimit()
1734  << ", " << currentModel->HighEnergyLimit() << ")"
1735  << " preStepLambda= " << preStepLambda
1736  << " dir= " << track.GetMomentumDirection()
1737  << " status= " << track.GetTrackStatus()
1738  << G4endl;
1739  }
1740  */
1741  return &fParticleChange;
1742 }
1743 
1744 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1745 
1747  const G4ParticleDefinition* part, const G4String& directory,
1748  G4bool ascii)
1749 {
1750  G4bool res = true;
1751  //G4cout << "G4VEnergyLossProcess::StorePhysicsTable: " << part->GetParticleName()
1752  // << " " << directory << " " << ascii << G4endl;
1753  if (!isMaster || baseParticle || part != particle ) return res;
1754 
1755  if(!StoreTable(part,theDEDXTable,ascii,directory,"DEDX"))
1756  {res = false;}
1757 
1758  if(!StoreTable(part,theDEDXunRestrictedTable,ascii,directory,"DEDXnr"))
1759  {res = false;}
1760 
1761  if(!StoreTable(part,theDEDXSubTable,ascii,directory,"SubDEDX"))
1762  {res = false;}
1763 
1764  if(!StoreTable(part,theIonisationTable,ascii,directory,"Ionisation"))
1765  {res = false;}
1766 
1767  if(!StoreTable(part,theIonisationSubTable,ascii,directory,"SubIonisation"))
1768  {res = false;}
1769 
1770  if(isIonisation &&
1771  !StoreTable(part,theCSDARangeTable,ascii,directory,"CSDARange"))
1772  {res = false;}
1773 
1774  if(isIonisation &&
1775  !StoreTable(part,theRangeTableForLoss,ascii,directory,"Range"))
1776  {res = false;}
1777 
1778  if(isIonisation &&
1779  !StoreTable(part,theInverseRangeTable,ascii,directory,"InverseRange"))
1780  {res = false;}
1781 
1782  if(!StoreTable(part,theLambdaTable,ascii,directory,"Lambda"))
1783  {res = false;}
1784 
1785  if(!StoreTable(part,theSubLambdaTable,ascii,directory,"SubLambda"))
1786  {res = false;}
1787 
1788  if ( !res ) {
1789  if(1 < verboseLevel) {
1790  G4cout << "Physics tables are stored for "
1792  << " and process " << GetProcessName()
1793  << " in the directory <" << directory
1794  << "> " << G4endl;
1795  }
1796  } else {
1797  G4cout << "Fail to store Physics Tables for "
1799  << " and process " << GetProcessName()
1800  << " in the directory <" << directory
1801  << "> " << G4endl;
1802  }
1803  return res;
1804 }
1805 
1806 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
1807 
1808 G4bool
1810  const G4String& directory,
1811  G4bool ascii)
1812 {
1813  G4bool res = true;
1814  if (!isMaster) return res;
1815  const G4String particleName = part->GetParticleName();
1816 
1817  if(1 < verboseLevel) {
1818  G4cout << "G4VEnergyLossProcess::RetrievePhysicsTable() for "
1819  << particleName << " and process " << GetProcessName()
1820  << "; tables_are_built= " << tablesAreBuilt
1821  << G4endl;
1822  }
1823  if(particle == part) {
1824 
1825  if ( !baseParticle ) {
1826 
1827  G4bool fpi = true;
1828  if(!RetrieveTable(part,theDEDXTable,ascii,directory,"DEDX",fpi))
1829  {fpi = false;}
1830 
1831  // ionisation table keeps individual dEdx and not sum of sub-processes
1832  if(!RetrieveTable(part,theDEDXTable,ascii,directory,"Ionisation",false))
1833  {fpi = false;}
1834 
1835  if(!RetrieveTable(part,theRangeTableForLoss,ascii,directory,"Range",fpi))
1836  {res = false;}
1837 
1838  if(!RetrieveTable(part,theDEDXunRestrictedTable,ascii,directory,
1839  "DEDXnr",false))
1840  {res = false;}
1841 
1842  if(!RetrieveTable(part,theCSDARangeTable,ascii,directory,
1843  "CSDARange",false))
1844  {res = false;}
1845 
1846  if(!RetrieveTable(part,theInverseRangeTable,ascii,directory,
1847  "InverseRange",fpi))
1848  {res = false;}
1849 
1850  if(!RetrieveTable(part,theLambdaTable,ascii,directory,"Lambda",true))
1851  {res = false;}
1852 
1853  G4bool yes = false;
1854  if(nSCoffRegions > 0) {yes = true;}
1855 
1856  if(!RetrieveTable(part,theDEDXSubTable,ascii,directory,"SubDEDX",yes))
1857  {res = false;}
1858 
1859  if(!RetrieveTable(part,theSubLambdaTable,ascii,directory,
1860  "SubLambda",yes))
1861  {res = false;}
1862 
1863  if(!fpi) yes = false;
1864  if(!RetrieveTable(part,theIonisationSubTable,ascii,directory,
1865  "SubIonisation",yes))
1866  {res = false;}
1867  }
1868  }
1869 
1870  return res;
1871 }
1872 
1873 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
1874 
1876  G4PhysicsTable* aTable, G4bool ascii,
1877  const G4String& directory,
1878  const G4String& tname)
1879 {
1880  //G4cout << "G4VEnergyLossProcess::StoreTable: " << aTable
1881  // << " " << directory << " " << tname << G4endl;
1882  G4bool res = true;
1883  if ( aTable ) {
1884  const G4String name = GetPhysicsTableFileName(part,directory,tname,ascii);
1885  G4cout << name << G4endl;
1886  //G4cout << *aTable << G4endl;
1887  if( !aTable->StorePhysicsTable(name,ascii)) res = false;
1888  }
1889  return res;
1890 }
1891 
1892 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
1893 
1894 G4bool
1896  G4PhysicsTable* aTable,
1897  G4bool ascii,
1898  const G4String& directory,
1899  const G4String& tname,
1900  G4bool mandatory)
1901 {
1902  G4bool isRetrieved = false;
1903  G4String filename = GetPhysicsTableFileName(part,directory,tname,ascii);
1904  if(aTable) {
1905  if(aTable->ExistPhysicsTable(filename)) {
1906  if(G4PhysicsTableHelper::RetrievePhysicsTable(aTable,filename,ascii)) {
1907  isRetrieved = true;
1908  if(theParameters->Spline()) {
1909  size_t n = aTable->length();
1910  for(size_t i=0; i<n; ++i) {
1911  if((*aTable)[i]) { (*aTable)[i]->SetSpline(true); }
1912  }
1913  }
1914  if (0 < verboseLevel) {
1915  G4cout << tname << " table for " << part->GetParticleName()
1916  << " is Retrieved from <" << filename << ">"
1917  << G4endl;
1918  }
1919  }
1920  }
1921  }
1922  if(mandatory && !isRetrieved) {
1923  if(0 < verboseLevel) {
1924  G4cout << tname << " table for " << part->GetParticleName()
1925  << " from file <"
1926  << filename << "> is not Retrieved"
1927  << G4endl;
1928  }
1929  return false;
1930  }
1931  return true;
1932 }
1933 
1934 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1935 
1937  const G4MaterialCutsCouple *couple,
1938  const G4DynamicParticle* dp,
1939  G4double length)
1940 {
1941  DefineMaterial(couple);
1942  G4double ekin = dp->GetKineticEnergy();
1943  SelectModel(ekin*massRatio);
1945  tmax = std::min(tmax,(*theCuts)[currentCoupleIndex]);
1946  G4double d = 0.0;
1948  if(fm) { d = fm->Dispersion(currentMaterial,dp,tmax,length); }
1949  return d;
1950 }
1951 
1952 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1953 
1955  G4double kineticEnergy, const G4MaterialCutsCouple* couple)
1956 {
1957  // Cross section per volume is calculated
1958  DefineMaterial(couple);
1959  G4double cross = 0.0;
1960  if(theLambdaTable) {
1961  cross = GetLambdaForScaledEnergy(kineticEnergy*massRatio);
1962  } else {
1963  SelectModel(kineticEnergy*massRatio);
1964  cross = biasFactor*(*theDensityFactor)[currentCoupleIndex]
1966  particle, kineticEnergy,
1968  }
1969  if(cross < 0.0) { cross = 0.0; }
1970  return cross;
1971 }
1972 
1973 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1974 
1976 {
1979  G4double x = DBL_MAX;
1980  if(0.0 < preStepLambda) { x = 1.0/preStepLambda; }
1981  return x;
1982 }
1983 
1984 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1985 
1987  G4double x, G4double y,
1988  G4double& z)
1989 {
1990  G4GPILSelection sel;
1991  return AlongStepGetPhysicalInteractionLength(track, x, y, z, &sel);
1992 }
1993 
1994 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1995 
1997  const G4Track& track,
1998  G4double,
2000 
2001 {
2002  *condition = NotForced;
2003  return MeanFreePath(track);
2004 }
2005 
2006 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2007 
2009  const G4Track&,
2011 {
2012  return DBL_MAX;
2013 }
2014 
2015 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2016 
2019  G4double)
2020 {
2021  G4PhysicsVector* v =
2024  return v;
2025 }
2026 
2027 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2028 
2031 {
2032  G4bool add = true;
2033  if(p->GetProcessName() != "eBrem") { add = false; }
2034  if(add && nProcesses > 0) {
2035  for(G4int i=0; i<nProcesses; ++i) {
2036  if(p == scProcesses[i]) {
2037  add = false;
2038  break;
2039  }
2040  }
2041  }
2042  if(add) {
2043  scProcesses.push_back(p);
2044  ++nProcesses;
2045  if (1 < verboseLevel) {
2046  G4cout << "### The process " << p->GetProcessName()
2047  << " is added to the list of collaborative processes of "
2048  << GetProcessName() << G4endl;
2049  }
2050  }
2051 }
2052 
2053 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2054 
2055 void
2057 {
2058  G4bool verb = false;
2059  if(fTotal == tType) {
2061  if(p) {
2062  size_t n = p->length();
2063  G4PhysicsVector* pv = (*p)[0];
2064  G4double emax = maxKinEnergyCSDA;
2065 
2069 
2070  for (size_t i=0; i<n; ++i) {
2071  G4double dedx = 0.0;
2072  pv = (*p)[i];
2073  if(pv) {
2074  dedx = pv->Value(emax, idxDEDXunRestricted);
2075  } else {
2076  pv = (*p)[(*theDensityIdx)[i]];
2077  if(pv) {
2078  dedx =
2079  pv->Value(emax, idxDEDXunRestricted)*(*theDensityFactor)[i];
2080  }
2081  }
2082  theDEDXAtMaxEnergy[i] = dedx;
2083  //G4cout << "i= " << i << " emax(MeV)= " << emax/MeV<< " dedx= "
2084  // << dedx << G4endl;
2085  }
2086  }
2087 
2088  } else if(fRestricted == tType) {
2089  if(verb) {
2090  G4cout<< "G4VEnergyLossProcess::SetDEDXTable "
2092  << " oldTable " << theDEDXTable << " newTable " << p
2093  << " ion " << theIonisationTable
2094  << " IsMaster " << isMaster
2095  << " " << GetProcessName() << G4endl;
2096  G4cout << (*p) << G4endl;
2097  }
2098  theDEDXTable = p;
2099  } else if(fSubRestricted == tType) {
2100  theDEDXSubTable = p;
2101  } else if(fIsIonisation == tType) {
2102  if(verb) {
2103  G4cout<< "G4VEnergyLossProcess::SetIonisationTable "
2105  << " oldTable " << theDEDXTable << " newTable " << p
2106  << " ion " << theIonisationTable
2107  << " IsMaster " << isMaster
2108  << " " << GetProcessName() << G4endl;
2109  }
2110  theIonisationTable = p;
2111  } else if(fIsSubIonisation == tType) {
2113  }
2114 }
2115 
2116 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2117 
2119 {
2120  theCSDARangeTable = p;
2121 
2122  if(p) {
2123  size_t n = p->length();
2124  G4PhysicsVector* pv;
2125  G4double emax = maxKinEnergyCSDA;
2126 
2127  for (size_t i=0; i<n; ++i) {
2128  pv = (*p)[i];
2129  G4double rmax = 0.0;
2130  if(pv) { rmax = pv->Value(emax, idxCSDA); }
2131  else {
2132  pv = (*p)[(*theDensityIdx)[i]];
2133  if(pv) { rmax = pv->Value(emax, idxCSDA)/(*theDensityFactor)[i]; }
2134  }
2135  theRangeAtMaxEnergy[i] = rmax;
2136  //G4cout << "i= " << i << " Emax(MeV)= " << emax/MeV << " Rmax= "
2137  //<< rmax<< G4endl;
2138  }
2139  }
2140 }
2141 
2142 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2143 
2145 {
2147  if(1 < verboseLevel) {
2148  G4cout << "### Set Range table " << p
2149  << " for " << particle->GetParticleName()
2150  << " and process " << GetProcessName() << G4endl;
2151  }
2152 }
2153 
2154 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2155 
2157 {
2159  if(1 < verboseLevel) {
2160  G4cout << "### Set SecondaryRange table " << p
2161  << " for " << particle->GetParticleName()
2162  << " and process " << GetProcessName() << G4endl;
2163  }
2164 }
2165 
2166 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2167 
2169 {
2171  if(1 < verboseLevel) {
2172  G4cout << "### Set InverseRange table " << p
2173  << " for " << particle->GetParticleName()
2174  << " and process " << GetProcessName() << G4endl;
2175  }
2176 }
2177 
2178 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2179 
2181 {
2182  if(1 < verboseLevel) {
2183  G4cout << "### Set Lambda table " << p
2184  << " for " << particle->GetParticleName()
2185  << " and process " << GetProcessName() << G4endl;
2186  //G4cout << *p << G4endl;
2187  }
2188  theLambdaTable = p;
2189  tablesAreBuilt = true;
2190 
2194 
2195  if(theLambdaTable) {
2196  size_t n = theLambdaTable->length();
2197  G4PhysicsVector* pv = (*theLambdaTable)[0];
2198  G4double e, ss, smax, emax;
2199 
2200  size_t i;
2201 
2202  // first loop on existing vectors
2203  for (i=0; i<n; ++i) {
2204  pv = (*theLambdaTable)[i];
2205  if(pv) {
2206  size_t nb = pv->GetVectorLength();
2207  emax = DBL_MAX;
2208  smax = 0.0;
2209  if(nb > 0) {
2210  for (size_t j=0; j<nb; ++j) {
2211  e = pv->Energy(j);
2212  ss = (*pv)(j);
2213  if(ss > smax) {
2214  smax = ss;
2215  emax = e;
2216  }
2217  }
2218  }
2219  theEnergyOfCrossSectionMax[i] = emax;
2220  theCrossSectionMax[i] = smax;
2221  if(1 < verboseLevel) {
2222  G4cout << "For " << particle->GetParticleName()
2223  << " Max CS at i= " << i << " emax(MeV)= " << emax/MeV
2224  << " lambda= " << smax << G4endl;
2225  }
2226  }
2227  }
2228  // second loop using base materials
2229  for (i=0; i<n; ++i) {
2230  pv = (*theLambdaTable)[i];
2231  if(!pv){
2232  G4int j = (*theDensityIdx)[i];
2234  theCrossSectionMax[i] = (*theDensityFactor)[i]*theCrossSectionMax[j];
2235  }
2236  }
2237  }
2238 }
2239 
2240 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2241 
2243 {
2244  theSubLambdaTable = p;
2245  if(1 < verboseLevel) {
2246  G4cout << "### Set SebLambda table " << p
2247  << " for " << particle->GetParticleName()
2248  << " and process " << GetProcessName() << G4endl;
2249  }
2250 }
2251 
2252 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2253 
2255 {
2256  const G4Element* elm = 0;
2257  if(currentModel) { elm = currentModel->GetCurrentElement(); }
2258  return elm;
2259 }
2260 
2261 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2262 
2264  G4bool flag)
2265 {
2266  if(f > 0.0) {
2267  biasFactor = f;
2268  weightFlag = flag;
2269  if(1 < verboseLevel) {
2270  G4cout << "### SetCrossSectionBiasingFactor: for "
2271  << " process " << GetProcessName()
2272  << " biasFactor= " << f << " weightFlag= " << flag
2273  << G4endl;
2274  }
2275  }
2276 }
2277 
2278 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2279 
2280 void
2282  const G4String& region,
2283  G4bool flag)
2284 {
2285  if(!biasManager) { biasManager = new G4EmBiasingManager(); }
2286  if(1 < verboseLevel) {
2287  G4cout << "### ActivateForcedInteraction: for "
2288  << " process " << GetProcessName()
2289  << " length(mm)= " << length/mm
2290  << " in G4Region <" << region
2291  << "> weightFlag= " << flag
2292  << G4endl;
2293  }
2294  weightFlag = flag;
2295  biasManager->ActivateForcedInteraction(length, region);
2296 }
2297 
2298 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2299 
2300 void
2302  G4double factor,
2303  G4double energyLimit)
2304 {
2305  if (0.0 <= factor) {
2306 
2307  // Range cut can be applied only for e-
2308  if(0.0 == factor && secondaryParticle != G4Electron::Electron())
2309  { return; }
2310 
2311  if(!biasManager) { biasManager = new G4EmBiasingManager(); }
2312  biasManager->ActivateSecondaryBiasing(region, factor, energyLimit);
2313  if(1 < verboseLevel) {
2314  G4cout << "### ActivateSecondaryBiasing: for "
2315  << " process " << GetProcessName()
2316  << " factor= " << factor
2317  << " in G4Region <" << region
2318  << "> energyLimit(MeV)= " << energyLimit/MeV
2319  << G4endl;
2320  }
2321  }
2322 }
2323 
2324 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2325 
2327 {
2328  isIonisation = val;
2329  if(val) { aGPILSelection = CandidateForSelection; }
2331 }
2332 
2333 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2334 
2336 {
2337  if(0.0 < val && val < 1.0) {
2338  linLossLimit = val;
2339  actLinLossLimit = true;
2340  } else { PrintWarning("SetLinearLossLimit", val); }
2341 }
2342 
2343 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2344 
2346 {
2347  if(0.0 < v1 && 0.0 < v2 && v2 < 1.e+50) {
2348  dRoverRange = std::min(1.0, v1);
2349  finalRange = v2;
2350  } else if(v1 <= 0.0) {
2351  PrintWarning("SetStepFunction", v1);
2352  } else {
2353  PrintWarning("SetStepFunction", v2);
2354  }
2355 }
2356 
2357 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2358 
2360 {
2361  if(1.e-18 < val && val < 1.e+50) { lowestKinEnergy = val; }
2362  else { PrintWarning("SetLowestEnergyLimit", val); }
2363 }
2364 
2365 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2366 
2368 {
2369  if(2 < n && n < 1000000000) {
2370  nBins = n;
2371  actBinning = true;
2372  } else {
2373  G4double e = (G4double)n;
2374  PrintWarning("SetDEDXBinning", e);
2375  }
2376 }
2377 
2378 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2379 
2381 {
2382  if(1.e-18 < e && e < maxKinEnergy) {
2383  minKinEnergy = e;
2384  actMinKinEnergy = true;
2385  } else { PrintWarning("SetMinKinEnergy", e); }
2386 }
2387 
2388 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2389 
2391 {
2392  if(minKinEnergy < e && e < 1.e+50) {
2393  maxKinEnergy = e;
2394  actMaxKinEnergy = true;
2395  if(e < maxKinEnergyCSDA) { maxKinEnergyCSDA = e; }
2396  } else { PrintWarning("SetMaxKinEnergy", e); }
2397 }
2398 
2399 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
2400 
2402 {
2403  G4String ss = "G4VEnergyLossProcess::" + tit;
2405  ed << "Parameter is out of range: " << val
2406  << " it will have no effect!\n" << " Process "
2407  << GetProcessName() << " nbins= " << nBins
2408  << " Emin(keV)= " << minKinEnergy/keV
2409  << " Emax(GeV)= " << maxKinEnergy/GeV;
2410  G4Exception(ss, "em0044", JustWarning, ed);
2411 }
2412 
2413 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
G4VEnergyLossProcess(const G4String &name="EnergyLoss", G4ProcessType type=fElectromagnetic)
static const double cm
Definition: G4SIunits.hh:106
G4double condition(const G4ErrorSymMatrix &m)
G4bool UseCutAsFinalRange() const
G4VSubCutProducer * SubCutProducer()
G4double MaxSecondaryKinEnergy(const G4DynamicParticle *dynParticle)
Definition: G4VEmModel.hh:464
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:193
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:354
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:616
void BuildPhysicsTable(const G4ParticleDefinition &)
const G4DynamicParticle * GetDynamicParticle() const
G4bool SecondaryBiasingRegion(G4int coupleIdx)
G4bool IsPIXEActive() const
G4PhysicsTable * SubLambdaTable() const
G4double z
Definition: TRTMaterials.hh:39
void SetLowEnergyLimit(G4double elimit)
G4VEmModel * GetModel(G4int, G4bool ver=false)
G4double GetProductionCut(G4int index) const
void DeRegister(G4VEnergyLossProcess *p)
G4PhysicsTable * RangeTableForLoss() const
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
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
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
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
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:592
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:707
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:686
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:93
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:402
const G4ThreeVector & GetPosition() const
const G4Element * GetCurrentElement() 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:196
const G4String & GetParticleType() const
void SetMaxKinEnergy(G4double e)
void Register(G4VEnergyLossProcess *p)
std::vector< G4double > theDEDXAtMaxEnergy
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
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
G4double AlongStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection)
void AddSecondary(G4Track *aSecondary)
virtual G4double ChargeSquareRatio(const G4Track &)
Definition: G4VEmModel.cc:329
const G4TouchableHandle & GetTouchableHandle() const
G4Material * GetMaterial() const
void RegisterExtraParticle(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p)
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:735
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()
G4bool RetrievePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii)
void SetMasterThread(G4bool val)
Definition: G4VEmModel.hh:693
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:298
virtual G4double GetParticleCharge(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:346
static G4ProductionCutsTable * GetProductionCutsTable()
const G4Material * currentMaterial
void SetLambdaTable(G4PhysicsTable *p)
static const double eV
Definition: G4SIunits.hh:194
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()
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)
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:197
G4bool StorePhysicsTable(const G4String &filename, G4bool ascii=false)
G4TrackStatus GetTrackStatus() const
G4VAtomDeexcitation * AtomDeexcitation()
static const double keV
Definition: G4SIunits.hh:195
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:102
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:449
G4bool IsIonisationProcess() const
void SetMinKinEnergy(G4double e)
G4EmModelManager * modelManager