Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VEmProcess.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: G4VEmProcess.cc 98778 2016-08-09 14:41:08Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4VEmProcess
34 //
35 // Author: Vladimir Ivanchenko on base of Laszlo Urban code
36 //
37 // Creation date: 01.10.2003
38 //
39 // Modifications:
40 // 30-06-04 make it to be pure discrete process (V.Ivanchenko)
41 // 30-09-08 optimise integral option (V.Ivanchenko)
42 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivanchenko)
43 // 11-03-05 Shift verbose level by 1, add applyCuts and killPrimary flags (VI)
44 // 14-03-05 Update logic PostStepDoIt (V.Ivanchenko)
45 // 08-04-05 Major optimisation of internal interfaces (V.Ivanchenko)
46 // 18-04-05 Use G4ParticleChangeForGamma (V.Ivanchenko)
47 // 25-07-05 Add protection: integral mode only for charged particles (VI)
48 // 04-09-05 default lambdaFactor 0.8 (V.Ivanchenko)
49 // 11-01-06 add A to parameters of ComputeCrossSectionPerAtom (VI)
50 // 12-09-06 add SetModel() (mma)
51 // 12-04-07 remove double call to Clear model manager (V.Ivanchenko)
52 // 27-10-07 Virtual functions moved to source (V.Ivanchenko)
53 // 24-06-09 Removed hidden bin in G4PhysicsVector (V.Ivanchenko)
54 // 17-02-10 Added pointer currentParticle (VI)
55 // 30-05-12 allow Russian roulette, brem splitting (D. Sawkey)
56 //
57 // Class Description:
58 //
59 
60 // -------------------------------------------------------------------
61 //
62 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
63 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
64 
65 #include "G4VEmProcess.hh"
66 #include "G4PhysicalConstants.hh"
67 #include "G4SystemOfUnits.hh"
68 #include "G4ProcessManager.hh"
69 #include "G4LossTableManager.hh"
70 #include "G4LossTableBuilder.hh"
71 #include "G4Step.hh"
72 #include "G4ParticleDefinition.hh"
73 #include "G4VEmModel.hh"
74 #include "G4DataVector.hh"
75 #include "G4PhysicsTable.hh"
76 #include "G4PhysicsLogVector.hh"
77 #include "G4VParticleChange.hh"
78 #include "G4ProductionCutsTable.hh"
79 #include "G4Region.hh"
80 #include "G4Gamma.hh"
81 #include "G4Electron.hh"
82 #include "G4Positron.hh"
83 #include "G4PhysicsTableHelper.hh"
84 #include "G4EmBiasingManager.hh"
85 #include "G4GenericIon.hh"
86 #include "G4Log.hh"
87 #include <iostream>
88 
89 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
90 
92  G4VDiscreteProcess(name, type),
93  secondaryParticle(nullptr),
94  buildLambdaTable(true),
95  numberOfModels(0),
96  theLambdaTable(nullptr),
97  theLambdaTablePrim(nullptr),
98  theDensityFactor(nullptr),
99  theDensityIdx(nullptr),
100  integral(false),
101  applyCuts(false),
102  startFromNull(false),
103  splineFlag(true),
104  currentModel(nullptr),
105  particle(nullptr),
106  currentParticle(nullptr),
107  currentCouple(nullptr)
108 {
109  theParameters = G4EmParameters::Instance();
110  SetVerboseLevel(1);
111 
112  // Size of tables assuming spline
113  minKinEnergy = 0.1*keV;
114  maxKinEnergy = 100.0*TeV;
115  nLambdaBins = 77;
116  minKinEnergyPrim = DBL_MAX;
117  actBinning = actSpline = actMinKinEnergy = actMaxKinEnergy = false;
118 
119  // default lambda factor
120  lambdaFactor = 0.8;
121 
122  // default limit on polar angle
123  biasFactor = fFactor = 1.0;
124 
125  // particle types
126  theGamma = G4Gamma::Gamma();
127  theElectron = G4Electron::Electron();
128  thePositron = G4Positron::Positron();
129 
130  theCuts = theCutsGamma = theCutsElectron = theCutsPositron = nullptr;
131 
134  secParticles.reserve(5);
135 
136  baseMaterial = currentMaterial = nullptr;
137 
138  preStepLambda = preStepKinEnergy = 0.0;
139  mfpKinEnergy = DBL_MAX;
140 
141  idxLambda = idxLambdaPrim = currentCoupleIndex
142  = basedCoupleIndex = 0;
143 
144  modelManager = new G4EmModelManager();
145  biasManager = nullptr;
146  biasFlag = false;
147  weightFlag = false;
148  lManager = G4LossTableManager::Instance();
149  lManager->Register(this);
150  secID = fluoID = augerID = biasID = -1;
151  mainSecondaries = 100;
152  if("phot" == GetProcessName() || "compt" == GetProcessName()) {
153  mainSecondaries = 1;
154  }
155 }
156 
157 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
158 
160 {
161  /*
162  if(1 < verboseLevel) {
163  G4cout << "G4VEmProcess destruct " << GetProcessName()
164  << " " << this << " " << theLambdaTable <<G4endl;
165  }
166  */
167  if(lManager->IsMaster()) {
168  if(theLambdaTable) {
169  theLambdaTable->clearAndDestroy();
170  delete theLambdaTable;
171  }
172  if(theLambdaTablePrim) {
173  theLambdaTablePrim->clearAndDestroy();
174  delete theLambdaTablePrim;
175  }
176  }
177  delete modelManager;
178  delete biasManager;
179  lManager->DeRegister(this);
180  //G4cout << "G4VEmProcess removed " << G4endl;
181 }
182 
183 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
184 
185 void G4VEmProcess::Clear()
186 {
187  currentCouple = nullptr;
188  preStepLambda = 0.0;
189  mfpKinEnergy = DBL_MAX;
190  idxLambda = idxLambdaPrim = 0;
191 }
192 
193 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
194 
196  const G4Material*)
197 {
198  return 0.0;
199 }
200 
201 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
202 
204  const G4Region* region)
205 {
206  G4VEmFluctuationModel* fm = nullptr;
207  modelManager->AddEmModel(order, p, fm, region);
208  if(p) { p->SetParticleChange(pParticleChange); }
209 }
210 
211 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
212 
214 {
215  G4int n = emModels.size();
216  if(index >= n) {
217  for(G4int i=n; i<=index; ++i) {emModels.push_back(nullptr);}
218  }
219  emModels[index] = p;
220 }
221 
222 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
223 
225 {
226  G4VEmModel* p = nullptr;
227  if(index >= 0 && index < G4int(emModels.size())) { p = emModels[index]; }
228  return p;
229 }
230 
231 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
232 
234  G4double emin, G4double emax)
235 {
236  modelManager->UpdateEmModel(nam, emin, emax);
237 }
238 
239 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
240 
242 {
243  return modelManager->NumberOfModels();
244 }
245 
246 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
247 
249 {
250  return modelManager->NumberOfRegionModels(couple_index);
251 }
252 
253 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
254 
255 G4VEmModel* G4VEmProcess::GetRegionModel(G4int idx, size_t couple_index) const
256 {
257  return modelManager->GetRegionModel(idx, couple_index);
258 }
259 
260 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
261 
263 {
264  return modelManager->GetModel(idx, ver);
265 }
266 
267 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
268 
270 {
271  G4bool isMaster = true;
272  const G4VEmProcess* masterProcess =
273  static_cast<const G4VEmProcess*>(GetMasterProcess());
274  if(masterProcess && masterProcess != this) { isMaster = false; }
275 
276  if(!particle) { SetParticle(&part); }
277 
278  if(part.GetParticleType() == "nucleus" &&
279  part.GetParticleSubType() == "generic") {
280 
281  G4String pname = part.GetParticleName();
282  if(pname != "deuteron" && pname != "triton" &&
283  pname != "alpha" && pname != "He3" &&
284  pname != "alpha+" && pname != "helium" &&
285  pname != "hydrogen") {
286 
287  particle = G4GenericIon::GenericIon();
288  }
289  }
290 
291  if(1 < verboseLevel) {
292  G4cout << "G4VEmProcess::PreparePhysicsTable() for "
293  << GetProcessName()
294  << " and particle " << part.GetParticleName()
295  << " local particle " << particle->GetParticleName()
296  << G4endl;
297  }
298 
299  if(particle != &part) { return; }
300 
301  G4LossTableBuilder* bld = lManager->GetTableBuilder();
302 
303  lManager->PreparePhysicsTable(&part, this, isMaster);
304 
305  Clear();
306  InitialiseProcess(particle);
307 
308  const G4ProductionCutsTable* theCoupleTable=
310  size_t n = theCoupleTable->GetTableSize();
311 
312  theEnergyOfCrossSectionMax.resize(n, 0.0);
313  theCrossSectionMax.resize(n, DBL_MAX);
314 
315  // initialisation of the process
316  if(!actMinKinEnergy) { minKinEnergy = theParameters->MinKinEnergy(); }
317  if(!actMaxKinEnergy) { maxKinEnergy = theParameters->MaxKinEnergy(); }
318  if(!actSpline) { splineFlag = theParameters->Spline(); }
319 
320  if(isMaster) { SetVerboseLevel(theParameters->Verbose()); }
321  else { SetVerboseLevel(theParameters->WorkerVerbose()); }
322  applyCuts = theParameters->ApplyCuts();
323  lambdaFactor = theParameters->LambdaFactor();
324  theParameters->DefineRegParamForEM(this);
325 
326  // initialisation of models
327  numberOfModels = modelManager->NumberOfModels();
328  for(G4int i=0; i<numberOfModels; ++i) {
329  G4VEmModel* mod = modelManager->GetModel(i);
330  if(0 == i) { currentModel = mod; }
331  mod->SetPolarAngleLimit(theParameters->MscThetaLimit());
332  mod->SetMasterThread(isMaster);
333  if(mod->HighEnergyLimit() > maxKinEnergy) {
334  mod->SetHighEnergyLimit(maxKinEnergy);
335  }
336  }
337 
338  if(lManager->AtomDeexcitation()) { modelManager->SetFluoFlag(true); }
339  theCuts = modelManager->Initialise(particle,secondaryParticle,
340  2.,verboseLevel);
341  theCutsGamma = theCoupleTable->GetEnergyCutsVector(idxG4GammaCut);
342  theCutsElectron = theCoupleTable->GetEnergyCutsVector(idxG4ElectronCut);
343  theCutsPositron = theCoupleTable->GetEnergyCutsVector(idxG4PositronCut);
344 
345  // prepare tables
346  if(buildLambdaTable && isMaster){
347  theLambdaTable =
349  bld->InitialiseBaseMaterials(theLambdaTable);
350  }
351  // high energy table
352  if(isMaster && minKinEnergyPrim < maxKinEnergy){
353  theLambdaTablePrim =
354  G4PhysicsTableHelper::PreparePhysicsTable(theLambdaTablePrim);
355  bld->InitialiseBaseMaterials(theLambdaTablePrim);
356  }
357  // forced biasing
358  if(biasManager) {
359  biasManager->Initialise(part,GetProcessName(),verboseLevel);
360  biasFlag = false;
361  }
362  // defined ID of secondary particles
363  G4String nam1 = GetProcessName();
364  G4String nam2 = nam1 + "_fluo" ;
365  G4String nam3 = nam1 + "_auger";
366  G4String nam4 = nam1 + "_split";
367  secID = G4PhysicsModelCatalog::Register(nam1);
368  fluoID = G4PhysicsModelCatalog::Register(nam2);
369  augerID = G4PhysicsModelCatalog::Register(nam3);
370  biasID = G4PhysicsModelCatalog::Register(nam4);
371 }
372 
373 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
374 
376 {
377  G4bool isMaster = true;
378  const G4VEmProcess* masterProc =
379  static_cast<const G4VEmProcess*>(GetMasterProcess());
380  if(masterProc && masterProc != this) { isMaster = false; }
381 
382  G4String num = part.GetParticleName();
383  if(1 < verboseLevel) {
384  G4cout << "### G4VEmProcess::BuildPhysicsTable() for "
385  << GetProcessName()
386  << " and particle " << num
387  << " buildLambdaTable= " << buildLambdaTable
388  << " isMaster= " << isMaster
389  << G4endl;
390  }
391 
392  if(particle == &part) {
393 
394  G4LossTableBuilder* bld = lManager->GetTableBuilder();
395 
396  // worker initialisation
397  if(!isMaster) {
398  theLambdaTable = masterProc->LambdaTable();
399  theLambdaTablePrim = masterProc->LambdaTablePrim();
400 
401  if(theLambdaTable) {
402  bld->InitialiseBaseMaterials(theLambdaTable);
403  } else if(theLambdaTablePrim) {
404  bld->InitialiseBaseMaterials(theLambdaTablePrim);
405  }
406  theDensityFactor = bld->GetDensityFactors();
407  theDensityIdx = bld->GetCoupleIndexes();
408  if(theLambdaTable) { FindLambdaMax(); }
409 
410  // local initialisation of models
411  G4bool printing = true;
412  numberOfModels = modelManager->NumberOfModels();
413  for(G4int i=0; i<numberOfModels; ++i) {
414  G4VEmModel* mod = GetModelByIndex(i, printing);
415  G4VEmModel* mod0= masterProc->GetModelByIndex(i, printing);
416  mod->InitialiseLocal(particle, mod0);
417  }
418  // master thread
419  } else {
420  theDensityFactor = bld->GetDensityFactors();
421  theDensityIdx = bld->GetCoupleIndexes();
422  if(buildLambdaTable || minKinEnergyPrim < maxKinEnergy) {
423  BuildLambdaTable();
424  }
425  }
426  }
427 
428  // explicitly defined printout by particle name
429  if(1 < verboseLevel ||
430  (0 < verboseLevel && (num == "gamma" || num == "e-" ||
431  num == "e+" || num == "mu+" ||
432  num == "mu-" || num == "proton"||
433  num == "pi+" || num == "pi-" ||
434  num == "kaon+" || num == "kaon-" ||
435  num == "alpha" || num == "anti_proton" ||
436  num == "GenericIon")))
437  {
438  PrintInfoProcess(part);
439  }
440 
441  if(1 < verboseLevel) {
442  G4cout << "### G4VEmProcess::BuildPhysicsTable() done for "
443  << GetProcessName()
444  << " and particle " << num
445  << G4endl;
446  }
447 }
448 
449 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
450 
451 void G4VEmProcess::BuildLambdaTable()
452 {
453  if(1 < verboseLevel) {
454  G4cout << "G4EmProcess::BuildLambdaTable() for process "
455  << GetProcessName() << " and particle "
456  << particle->GetParticleName() << " " << this
457  << G4endl;
458  }
459 
460  // Access to materials
461  const G4ProductionCutsTable* theCoupleTable=
463  size_t numOfCouples = theCoupleTable->GetTableSize();
464 
465  G4LossTableBuilder* bld = lManager->GetTableBuilder();
466 
467  G4PhysicsLogVector* aVector = nullptr;
468  G4PhysicsLogVector* aVectorPrim = nullptr;
469  G4PhysicsLogVector* bVectorPrim = nullptr;
470 
471  G4double scale = theParameters->MaxKinEnergy()/theParameters->MinKinEnergy();
472  G4int nbin = theParameters->NumberOfBinsPerDecade()
473  *G4lrint(std::log10(scale));
474  scale = G4Log(scale);
475  if(actBinning) { nbin = std::max(nbin, nLambdaBins); }
476  G4double emax1 = std::min(maxKinEnergy, minKinEnergyPrim);
477 
478  for(size_t i=0; i<numOfCouples; ++i) {
479 
480  if (bld->GetFlag(i)) {
481 
482  // create physics vector and fill it
483  const G4MaterialCutsCouple* couple =
484  theCoupleTable->GetMaterialCutsCouple(i);
485 
486  // build main table
487  if(buildLambdaTable) {
488  delete (*theLambdaTable)[i];
489 
490  // if start from zero then change the scale
491  G4double emin = minKinEnergy;
492  G4bool startNull = false;
493  if(startFromNull) {
494  G4double e = MinPrimaryEnergy(particle,couple->GetMaterial());
495  if(e >= emin) {
496  emin = e;
497  startNull = true;
498  }
499  }
500  G4double emax = emax1;
501  if(emax <= emin) { emax = 2*emin; }
502  G4int bin = G4lrint(nbin*G4Log(emax/emin)/scale);
503  if(bin < 3) { bin = 3; }
504  aVector = new G4PhysicsLogVector(emin, emax, bin);
505  aVector->SetSpline(splineFlag);
506  modelManager->FillLambdaVector(aVector, couple, startNull);
507  if(splineFlag) { aVector->FillSecondDerivatives(); }
508  G4PhysicsTableHelper::SetPhysicsVector(theLambdaTable, i, aVector);
509  }
510  // build high energy table
511  if(minKinEnergyPrim < maxKinEnergy) {
512  delete (*theLambdaTablePrim)[i];
513 
514  // start not from zero
515  if(!bVectorPrim) {
516  G4int bin = G4lrint(nbin*G4Log(maxKinEnergy/minKinEnergyPrim)/scale);
517  if(bin < 3) { bin = 3; }
518  aVectorPrim =
519  new G4PhysicsLogVector(minKinEnergyPrim, maxKinEnergy, bin);
520  bVectorPrim = aVectorPrim;
521  } else {
522  aVectorPrim = new G4PhysicsLogVector(*bVectorPrim);
523  }
524  // always use spline
525  aVectorPrim->SetSpline(splineFlag);
526  modelManager->FillLambdaVector(aVectorPrim, couple, false,
528  aVectorPrim->FillSecondDerivatives();
529  G4PhysicsTableHelper::SetPhysicsVector(theLambdaTablePrim, i,
530  aVectorPrim);
531  }
532  }
533  }
534 
535  if(buildLambdaTable) { FindLambdaMax(); }
536 
537  if(1 < verboseLevel) {
538  G4cout << "Lambda table is built for "
539  << particle->GetParticleName()
540  << G4endl;
541  }
542 }
543 
544 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
545 
546 void G4VEmProcess::PrintInfoProcess(const G4ParticleDefinition& part)
547 {
548  if(verboseLevel > 0) {
549  G4cout << std::setprecision(6);
550  G4cout << G4endl << GetProcessName() << ": for "
551  << part.GetParticleName();
552  if(integral) { G4cout << ", integral: 1 "; }
553  if(applyCuts) { G4cout << ", applyCuts: 1 "; }
554  G4cout << " SubType= " << GetProcessSubType();;
555  if(biasFactor != 1.0) { G4cout << " BiasingFactor= " << biasFactor; }
556  G4cout << " BuildTable= " << buildLambdaTable;
557  G4cout << G4endl;
558  if(buildLambdaTable) {
559  if(particle == &part) {
560  size_t length = theLambdaTable->length();
561  for(size_t i=0; i<length; ++i) {
562  G4PhysicsVector* v = (*theLambdaTable)[i];
563  if(v) {
564  G4cout << " Lambda table from ";
565  G4double emin = v->Energy(0);
566  G4double emax = v->GetMaxEnergy();
567  G4int nbin = v->GetVectorLength() - 1;
568  if(emin > minKinEnergy) { G4cout << "threshold "; }
569  else { G4cout << G4BestUnit(emin,"Energy"); }
570  G4cout << " to "
571  << G4BestUnit(emax,"Energy")
572  << ", " << G4lrint(nbin/std::log10(emax/emin))
573  << " bins per decade, spline: "
574  << splineFlag
575  << G4endl;
576  break;
577  }
578  }
579  } else {
580  G4cout << " Used Lambda table of "
581  << particle->GetParticleName() << G4endl;;
582  }
583  }
584  if(minKinEnergyPrim < maxKinEnergy) {
585  if(particle == &part) {
586  size_t length = theLambdaTablePrim->length();
587  for(size_t i=0; i<length; ++i) {
588  G4PhysicsVector* v = (*theLambdaTablePrim)[i];
589  if(v) {
590  G4cout << " LambdaPrime table from "
591  << G4BestUnit(v->Energy(0),"Energy")
592  << " to "
593  << G4BestUnit(v->GetMaxEnergy(),"Energy")
594  << " in " << v->GetVectorLength()-1
595  << " bins "
596  << G4endl;
597  break;
598  }
599  }
600  } else {
601  G4cout << " Used LambdaPrime table of "
602  << particle->GetParticleName() << G4endl;;
603  }
604  }
605  PrintInfo();
606  modelManager->DumpModelList(verboseLevel);
607  }
608 
609  if(verboseLevel > 2 && buildLambdaTable) {
610  G4cout << " LambdaTable address= " << theLambdaTable << G4endl;
611  if(theLambdaTable && particle == &part) {
612  G4cout << (*theLambdaTable) << G4endl;
613  }
614  }
615 }
616 
617 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
618 
620 {
621  // reset parameters for the new track
622  currentParticle = track->GetParticleDefinition();
624  mfpKinEnergy = DBL_MAX;
625 
626  // forced biasing only for primary particles
627  if(biasManager) {
628  if(0 == track->GetParentID()) {
629  // primary particle
630  biasFlag = true;
631  biasManager->ResetForcedInteraction();
632  }
633  }
634 }
635 
636 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
637 
639  const G4Track& track,
640  G4double previousStepSize,
642 {
643  *condition = NotForced;
644  G4double x = DBL_MAX;
645 
646  preStepKinEnergy = track.GetKineticEnergy();
647  DefineMaterial(track.GetMaterialCutsCouple());
648  SelectModel(preStepKinEnergy, currentCoupleIndex);
649 
650  if(!currentModel->IsActive(preStepKinEnergy)) {
653  return x;
654  }
655 
656  // forced biasing only for primary particles
657  if(biasManager) {
658  if(0 == track.GetParentID()) {
659  if(biasFlag &&
660  biasManager->ForcedInteractionRegion(currentCoupleIndex)) {
661  return biasManager->GetStepLimit(currentCoupleIndex, previousStepSize);
662  }
663  }
664  }
665 
666  // compute mean free path
667  if(preStepKinEnergy < mfpKinEnergy) {
668  if (integral) { ComputeIntegralLambda(preStepKinEnergy); }
669  else { preStepLambda = GetCurrentLambda(preStepKinEnergy); }
670 
671  // zero cross section
672  if(preStepLambda <= 0.0) {
675  }
676  }
677 
678  // non-zero cross section
679  if(preStepLambda > 0.0) {
680 
682 
683  // beggining of tracking (or just after DoIt of this process)
686 
687  } else if(currentInteractionLength < DBL_MAX) {
688 
690  previousStepSize/currentInteractionLength;
693  }
694 
695  // new mean free path and step limit for the next step
696  currentInteractionLength = 1.0/preStepLambda;
698  /*
699 #ifdef G4VERBOSE
700  if (verboseLevel>2){
701  G4cout << "G4VEmProcess::PostStepGetPhysicalInteractionLength ";
702  G4cout << "[ " << GetProcessName() << "]" << G4endl;
703  G4cout << " for " << currentParticle->GetParticleName()
704  << " in Material " << currentMaterial->GetName()
705  << " Ekin(MeV)= " << preStepKinEnergy/MeV
706  <<G4endl;
707  G4cout << " MeanFreePath = " << currentInteractionLength/cm << "[cm]"
708  << " InteractionLength= " << x/cm <<"[cm] " <<G4endl;
709  }
710 #endif
711  */
712  }
713  return x;
714 }
715 
716 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
717 
719  const G4Step& step)
720 {
721  // In all cases clear number of interaction lengths
723  mfpKinEnergy = DBL_MAX;
724 
726 
727  // Do not make anything if particle is stopped, the annihilation then
728  // should be performed by the AtRestDoIt!
729  if (track.GetTrackStatus() == fStopButAlive) { return &fParticleChange; }
730 
731  G4double finalT = track.GetKineticEnergy();
732 
733  // forced process - should happen only once per track
734  if(biasFlag) {
735  if(biasManager->ForcedInteractionRegion(currentCoupleIndex)) {
736  biasFlag = false;
737  }
738  }
739 
740  // Integral approach
741  if (integral) {
742  G4double lx = GetLambda(finalT, currentCouple);
743  if(preStepLambda<lx && 1 < verboseLevel) {
744  G4cout << "WARNING: for " << currentParticle->GetParticleName()
745  << " and " << GetProcessName()
746  << " E(MeV)= " << finalT/MeV
747  << " preLambda= " << preStepLambda << " < "
748  << lx << " (postLambda) "
749  << G4endl;
750  }
751 
752  if(preStepLambda*G4UniformRand() > lx) {
754  return &fParticleChange;
755  }
756  }
757 
758  SelectModel(finalT, currentCoupleIndex);
759  if(!currentModel->IsActive(finalT)) { return &fParticleChange; }
760 
761  // define new weight for primary and secondaries
763  if(weightFlag) {
764  weight /= biasFactor;
766  }
767 
768  /*
769  if(0 < verboseLevel) {
770  G4cout << "G4VEmProcess::PostStepDoIt: Sample secondary; E= "
771  << finalT/MeV
772  << " MeV; model= (" << currentModel->LowEnergyLimit()
773  << ", " << currentModel->HighEnergyLimit() << ")"
774  << G4endl;
775  }
776  */
777 
778  // sample secondaries
779  secParticles.clear();
780  currentModel->SampleSecondaries(&secParticles,
781  currentCouple,
782  track.GetDynamicParticle(),
783  (*theCuts)[currentCoupleIndex]);
784 
785  G4int num0 = secParticles.size();
786 
787  // splitting or Russian roulette
788  if(biasManager) {
789  if(biasManager->SecondaryBiasingRegion(currentCoupleIndex)) {
790  G4double eloss = 0.0;
791  weight *= biasManager->ApplySecondaryBiasing(
792  secParticles, track, currentModel, &fParticleChange, eloss,
793  currentCoupleIndex, (*theCuts)[currentCoupleIndex],
794  step.GetPostStepPoint()->GetSafety());
795  if(eloss > 0.0) {
798  }
799  }
800  }
801 
802  // save secondaries
803  G4int num = secParticles.size();
804  if(num > 0) {
805 
808  G4double time = track.GetGlobalTime();
809 
810  for (G4int i=0; i<num; ++i) {
811  if (secParticles[i]) {
812  G4DynamicParticle* dp = secParticles[i];
814  G4double e = dp->GetKineticEnergy();
815  G4bool good = true;
816  if(applyCuts) {
817  if (p == theGamma) {
818  if (e < (*theCutsGamma)[currentCoupleIndex]) { good = false; }
819 
820  } else if (p == theElectron) {
821  if (e < (*theCutsElectron)[currentCoupleIndex]) { good = false; }
822 
823  } else if (p == thePositron) {
824  if (electron_mass_c2 < (*theCutsGamma)[currentCoupleIndex] &&
825  e < (*theCutsPositron)[currentCoupleIndex]) {
826  good = false;
827  e += 2.0*electron_mass_c2;
828  }
829  }
830  // added secondary if it is good
831  }
832  if (good) {
833  G4Track* t = new G4Track(dp, time, track.GetPosition());
835  t->SetWeight(weight);
837 
838  // define type of secondary
839  if(i < mainSecondaries) { t->SetCreatorModelIndex(secID); }
840  else if(i < num0) {
841  if(p == theGamma) {
842  t->SetCreatorModelIndex(fluoID);
843  } else {
844  t->SetCreatorModelIndex(augerID);
845  }
846  } else {
847  t->SetCreatorModelIndex(biasID);
848  }
849 
850  //G4cout << "Secondary(post step) has weight " << t->GetWeight()
851  // << ", Ekin= " << t->GetKineticEnergy()/MeV << " MeV" <<G4endl;
852  } else {
853  delete dp;
854  edep += e;
855  }
856  }
857  }
859  }
860 
863  if(particle->GetProcessManager()->GetAtRestProcessVector()->size() > 0)
866  }
867 
868  return &fParticleChange;
869 }
870 
871 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
872 
874  const G4String& directory,
875  G4bool ascii)
876 {
877  G4bool yes = true;
878  const G4VEmProcess* masterProc =
879  static_cast<const G4VEmProcess*>(GetMasterProcess());
880  if(masterProc && masterProc != this) { return yes; }
881 
882  if ( theLambdaTable && part == particle) {
883  const G4String name =
884  GetPhysicsTableFileName(part,directory,"Lambda",ascii);
885  yes = theLambdaTable->StorePhysicsTable(name,ascii);
886 
887  if ( yes ) {
888  G4cout << "Physics table is stored for " << particle->GetParticleName()
889  << " and process " << GetProcessName()
890  << " in the directory <" << directory
891  << "> " << G4endl;
892  } else {
893  G4cout << "Fail to store Physics Table for "
894  << particle->GetParticleName()
895  << " and process " << GetProcessName()
896  << " in the directory <" << directory
897  << "> " << G4endl;
898  }
899  }
900  if ( theLambdaTablePrim && part == particle) {
901  const G4String name =
902  GetPhysicsTableFileName(part,directory,"LambdaPrim",ascii);
903  yes = theLambdaTablePrim->StorePhysicsTable(name,ascii);
904 
905  if ( yes ) {
906  G4cout << "Physics table prim is stored for "
907  << particle->GetParticleName()
908  << " and process " << GetProcessName()
909  << " in the directory <" << directory
910  << "> " << G4endl;
911  } else {
912  G4cout << "Fail to store Physics Table Prim for "
913  << particle->GetParticleName()
914  << " and process " << GetProcessName()
915  << " in the directory <" << directory
916  << "> " << G4endl;
917  }
918  }
919  return yes;
920 }
921 
922 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
923 
925  const G4String& directory,
926  G4bool ascii)
927 {
928  if(1 < verboseLevel) {
929  G4cout << "G4VEmProcess::RetrievePhysicsTable() for "
930  << part->GetParticleName() << " and process "
931  << GetProcessName() << G4endl;
932  }
933  G4bool yes = true;
934 
935  if((!buildLambdaTable && minKinEnergyPrim > maxKinEnergy)
936  || particle != part) { return yes; }
937 
938  const G4String particleName = part->GetParticleName();
939  G4String filename;
940 
941  if(buildLambdaTable) {
942  filename = GetPhysicsTableFileName(part,directory,"Lambda",ascii);
943  yes = G4PhysicsTableHelper::RetrievePhysicsTable(theLambdaTable,
944  filename,ascii);
945  if ( yes ) {
946  if (0 < verboseLevel) {
947  G4cout << "Lambda table for " << particleName
948  << " is Retrieved from <"
949  << filename << ">"
950  << G4endl;
951  }
952  if(theParameters->Spline()) {
953  size_t n = theLambdaTable->length();
954  for(size_t i=0; i<n; ++i) {
955  if((* theLambdaTable)[i]) {
956  (* theLambdaTable)[i]->SetSpline(true);
957  }
958  }
959  }
960  } else {
961  if (1 < verboseLevel) {
962  G4cout << "Lambda table for " << particleName << " in file <"
963  << filename << "> is not exist"
964  << G4endl;
965  }
966  }
967  }
968  if(minKinEnergyPrim < maxKinEnergy) {
969  filename = GetPhysicsTableFileName(part,directory,"LambdaPrim",ascii);
970  yes = G4PhysicsTableHelper::RetrievePhysicsTable(theLambdaTablePrim,
971  filename,ascii);
972  if ( yes ) {
973  if (0 < verboseLevel) {
974  G4cout << "Lambda table prim for " << particleName
975  << " is Retrieved from <"
976  << filename << ">"
977  << G4endl;
978  }
979  if(theParameters->Spline()) {
980  size_t n = theLambdaTablePrim->length();
981  for(size_t i=0; i<n; ++i) {
982  if((* theLambdaTablePrim)[i]) {
983  (* theLambdaTablePrim)[i]->SetSpline(true);
984  }
985  }
986  }
987  } else {
988  if (1 < verboseLevel) {
989  G4cout << "Lambda table prim for " << particleName << " in file <"
990  << filename << "> is not exist"
991  << G4endl;
992  }
993  }
994  }
995 
996  return yes;
997 }
998 
999 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1000 
1001 G4double
1003  const G4MaterialCutsCouple* couple)
1004 {
1005  // Cross section per atom is calculated
1006  DefineMaterial(couple);
1007  G4double cross = 0.0;
1008  if(buildLambdaTable && theLambdaTable) {
1009  cross = GetCurrentLambda(kineticEnergy);
1010  } else {
1011  SelectModel(kineticEnergy, currentCoupleIndex);
1012  cross = fFactor*currentModel->CrossSectionPerVolume(currentMaterial,
1013  currentParticle,
1014  kineticEnergy);
1015  }
1016 
1017  if(cross < 0.0) { cross = 0.0; }
1018  return cross;
1019 }
1020 
1021 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1022 
1024  G4double,
1026 {
1027  *condition = NotForced;
1028  return G4VEmProcess::MeanFreePath(track);
1029 }
1030 
1031 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1032 
1034 {
1035  DefineMaterial(track.GetMaterialCutsCouple());
1036  preStepLambda = GetCurrentLambda(track.GetKineticEnergy());
1037  G4double x = DBL_MAX;
1038  if(0.0 < preStepLambda) { x = 1.0/preStepLambda; }
1039  return x;
1040 }
1041 
1042 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1043 
1044 G4double
1046  G4double Z, G4double A, G4double cut)
1047 {
1048  SelectModel(kineticEnergy, currentCoupleIndex);
1049  G4double x = 0.0;
1050  if(currentModel) {
1051  x = currentModel->ComputeCrossSectionPerAtom(currentParticle,kineticEnergy,
1052  Z,A,cut);
1053  }
1054  return x;
1055 }
1056 
1057 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1058 
1059 void G4VEmProcess::FindLambdaMax()
1060 {
1061  if(1 < verboseLevel) {
1062  G4cout << "### G4VEmProcess::FindLambdaMax: "
1063  << particle->GetParticleName()
1064  << " and process " << GetProcessName() << " " << G4endl;
1065  }
1066  size_t n = theLambdaTable->length();
1067  G4PhysicsVector* pv;
1068  G4double e, ss, emax, smax;
1069 
1070  size_t i;
1071 
1072  // first loop on existing vectors
1073  for (i=0; i<n; ++i) {
1074  pv = (*theLambdaTable)[i];
1075  if(pv) {
1076  size_t nb = pv->GetVectorLength();
1077  emax = DBL_MAX;
1078  smax = 0.0;
1079  if(nb > 0) {
1080  for (size_t j=0; j<nb; ++j) {
1081  e = pv->Energy(j);
1082  ss = (*pv)(j);
1083  if(ss > smax) {
1084  smax = ss;
1085  emax = e;
1086  }
1087  }
1088  }
1089  theEnergyOfCrossSectionMax[i] = emax;
1090  theCrossSectionMax[i] = smax;
1091  if(1 < verboseLevel) {
1092  G4cout << "For " << particle->GetParticleName()
1093  << " Max CS at i= " << i << " emax(MeV)= " << emax/MeV
1094  << " lambda= " << smax << G4endl;
1095  }
1096  }
1097  }
1098  // second loop using base materials
1099  for (i=0; i<n; ++i) {
1100  pv = (*theLambdaTable)[i];
1101  if(!pv){
1102  G4int j = (*theDensityIdx)[i];
1103  theEnergyOfCrossSectionMax[i] = theEnergyOfCrossSectionMax[j];
1104  theCrossSectionMax[i] = (*theDensityFactor)[i]*theCrossSectionMax[j];
1105  }
1106  }
1107 }
1108 
1109 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1110 
1113 {
1114  G4PhysicsVector* v =
1115  new G4PhysicsLogVector(minKinEnergy, maxKinEnergy, nLambdaBins);
1116  v->SetSpline(theParameters->Spline());
1117  return v;
1118 }
1119 
1120 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1121 
1123 {
1124  const G4Element* elm = nullptr;
1125  if(currentModel) {elm = currentModel->GetCurrentElement(); }
1126  return elm;
1127 }
1128 
1129 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1130 
1132 {
1133  if(f > 0.0) {
1134  biasFactor = f;
1135  weightFlag = flag;
1136  if(1 < verboseLevel) {
1137  G4cout << "### SetCrossSectionBiasingFactor: for "
1138  << particle->GetParticleName()
1139  << " and process " << GetProcessName()
1140  << " biasFactor= " << f << " weightFlag= " << flag
1141  << G4endl;
1142  }
1143  }
1144 }
1145 
1146 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1147 
1148 void
1150  G4bool flag)
1151 {
1152  if(!biasManager) { biasManager = new G4EmBiasingManager(); }
1153  if(1 < verboseLevel) {
1154  G4cout << "### ActivateForcedInteraction: for "
1155  << particle->GetParticleName()
1156  << " and process " << GetProcessName()
1157  << " length(mm)= " << length/mm
1158  << " in G4Region <" << r
1159  << "> weightFlag= " << flag
1160  << G4endl;
1161  }
1162  weightFlag = flag;
1163  biasManager->ActivateForcedInteraction(length, r);
1164 }
1165 
1166 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1167 
1168 void
1170  G4double factor,
1171  G4double energyLimit)
1172 {
1173  if (0.0 <= factor) {
1174 
1175  // Range cut can be applied only for e-
1176  if(0.0 == factor && secondaryParticle != G4Electron::Electron())
1177  { return; }
1178 
1179  if(!biasManager) { biasManager = new G4EmBiasingManager(); }
1180  biasManager->ActivateSecondaryBiasing(region, factor, energyLimit);
1181  if(1 < verboseLevel) {
1182  G4cout << "### ActivateSecondaryBiasing: for "
1183  << " process " << GetProcessName()
1184  << " factor= " << factor
1185  << " in G4Region <" << region
1186  << "> energyLimit(MeV)= " << energyLimit/MeV
1187  << G4endl;
1188  }
1189  }
1190 }
1191 
1192 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1193 
1195 {
1196  if(5 < n && n < 10000000) {
1197  nLambdaBins = n;
1198  actBinning = true;
1199  } else {
1200  G4double e = (G4double)n;
1201  PrintWarning("SetLambdaBinning", e);
1202  }
1203 }
1204 
1205 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1206 
1208 {
1209  if(1.e-3*eV < e && e < maxKinEnergy) {
1210  nLambdaBins = G4lrint(nLambdaBins*G4Log(maxKinEnergy/e)
1211  /G4Log(maxKinEnergy/minKinEnergy));
1212  minKinEnergy = e;
1213  actMinKinEnergy = true;
1214  } else { PrintWarning("SetMinKinEnergy", e); }
1215 }
1216 
1217 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1218 
1220 {
1221  if(minKinEnergy < e && e < 1.e+6*TeV) {
1222  nLambdaBins = G4lrint(nLambdaBins*G4Log(e/minKinEnergy)
1223  /G4Log(maxKinEnergy/minKinEnergy));
1224  maxKinEnergy = e;
1225  actMaxKinEnergy = true;
1226  } else { PrintWarning("SetMaxKinEnergy", e); }
1227 }
1228 
1229 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1230 
1232 {
1233  if(theParameters->MinKinEnergy() <= e &&
1234  e <= theParameters->MaxKinEnergy()) { minKinEnergyPrim = e; }
1235  else { PrintWarning("SetMinKinEnergyPrim", e); }
1236 }
1237 
1238 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1239 
1240 void G4VEmProcess::PrintWarning(G4String tit, G4double val)
1241 {
1242  G4String ss = "G4VEmProcess::" + tit;
1244  ed << "Parameter is out of range: " << val
1245  << " it will have no effect!\n" << " Process "
1246  << GetProcessName() << " nbins= " << theParameters->NumberOfBins()
1247  << " Emin(keV)= " << theParameters->MinKinEnergy()/keV
1248  << " Emax(GeV)= " << theParameters->MaxKinEnergy()/GeV;
1249  G4Exception(ss, "em0044", JustWarning, ed);
1250 }
1251 
1252 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1253 
1254 void G4VEmProcess::ProcessDescription(std::ostream& outFile) const
1255 {
1256  outFile << "EM process <" << GetProcessName() << ">" << G4endl;
1257 }
1258 
1259 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
void ActivateForcedInteraction(G4double length=0.0, const G4String &r="", G4bool flag=true)
G4int NumberOfRegionModels(size_t index_couple) const
G4double condition(const G4ErrorSymMatrix &m)
G4PhysicsTable * LambdaTable() const
G4int NumberOfBinsPerDecade() const
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
const XML_Char * name
Definition: expat.h:151
virtual void PrintInfo()=0
G4int GetParentID() const
G4double GetLambda(G4double &kinEnergy, const G4MaterialCutsCouple *couple)
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:257
G4int NumberOfBins() const
virtual void StartTracking(G4Track *) override
G4PhysicsTable * LambdaTablePrim() const
const std::vector< G4double > * GetEnergyCutsVector(size_t pcIdx) const
G4VEmModel * SelectModel(G4double &kinEnergy, size_t index)
const std::vector< G4double > * GetDensityFactors()
static void SetPhysicsVector(G4PhysicsTable *physTable, size_t idx, G4PhysicsVector *vec)
G4bool Spline() const
virtual void BuildPhysicsTable(const G4ParticleDefinition &) override
G4int WorkerVerbose() const
G4double MaxKinEnergy() const
static G4LossTableManager * Instance()
G4int GetNumberOfRegionModels(size_t couple_index) const
void ActivateForcedInteraction(G4double length=0.0, const G4String &r="")
G4double GetMaxEnergy() const
static constexpr double mm
Definition: G4SIunits.hh:115
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
tuple bin
Definition: plottest35.py:22
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
G4VEmModel * GetModel(G4int idx, G4bool ver=false)
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:636
G4int GetNumberOfModels() const
G4VEmProcess(const G4String &name, G4ProcessType type=fElectromagnetic)
Definition: G4VEmProcess.cc:91
const G4DynamicParticle * GetDynamicParticle() const
G4bool SecondaryBiasingRegion(G4int coupleIdx)
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
Definition: G4VEmModel.cc:418
G4double ComputeCrossSectionPerAtom(G4double kineticEnergy, G4double Z, G4double A=0., G4double cut=0.0)
G4VEmModel * EmModel(G4int index=1) const
void DeRegister(G4VEnergyLossProcess *p)
const char * p
Definition: xmltok.h:285
void UpdateEmModel(const G4String &model_name, G4double emin, G4double emax)
G4double MscThetaLimit() const
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &) override
const G4ThreeVector & GetPosition() const
G4ParticleChangeForGamma fParticleChange
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel)
Definition: G4VEmModel.cc:218
G4TrackStatus GetTrackStatus() const
G4bool ForcedInteractionRegion(G4int coupleIdx)
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
virtual G4bool StorePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false) override
virtual void InitialiseProcess(const G4ParticleDefinition *)=0
virtual ~G4VEmProcess()
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
G4bool GetFlag(size_t idx) const
void SetTouchableHandle(const G4TouchableHandle &apValue)
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
G4double GetParentWeight() const
tuple x
Definition: test.py:50
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
G4bool ApplyCuts() const
size_t GetVectorLength() const
const G4String & GetParticleSubType() const
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:447
void FillLambdaVector(G4PhysicsVector *, const G4MaterialCutsCouple *, G4bool startFromNull=true, G4EmTableType t=fRestricted)
virtual G4double GetMeanFreePath(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
int G4int
Definition: G4Types.hh:78
void FillSecondDerivatives()
void UpdateEmModel(const G4String &, G4double, G4double)
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
const G4String & GetParticleName() const
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
G4LossTableBuilder * GetTableBuilder()
void SetWeight(G4double aValue)
virtual G4bool RetrievePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii) override
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void SetCreatorModelIndex(G4int idx)
virtual void PreparePhysicsTable(const G4ParticleDefinition &) override
void SetSpline(G4bool)
void Initialise(const G4ParticleDefinition &part, const G4String &procName, G4int verbose)
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:732
void SetMinKinEnergyPrim(G4double e)
void SetCrossSectionBiasingFactor(G4double f, G4bool flag=true)
static constexpr double TeV
Definition: G4SIunits.hh:218
virtual G4double MinPrimaryEnergy(const G4ParticleDefinition *, const G4Material *)
G4double MeanFreePath(const G4Track &track)
void SetLambdaBinning(G4int nbins)
void ProposeWeight(G4double finalWeight)
void SetSecondaryWeightByProcess(G4bool)
G4double LambdaFactor() const
void SetEmModel(G4VEmModel *, G4int index=1)
G4ProcessVector * GetAtRestProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4double GetKineticEnergy() const
G4VEmModel * GetRegionModel(G4int idx, size_t index_couple)
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
G4int Verbose() const
double A(double temperature)
G4VEmModel * GetRegionModel(G4int idx, size_t couple_index) const
bool G4bool
Definition: G4Types.hh:79
void SetParticle(const G4ParticleDefinition *p)
G4double GetStepLimit(G4int coupleIdx, G4double previousStep)
G4double currentInteractionLength
Definition: G4VProcess.hh:297
const G4ParticleDefinition * GetParticleDefinition() const
static constexpr double eV
Definition: G4SIunits.hh:215
G4bool IsMaster() const
const G4String & GetParticleType() const
void Register(G4VEnergyLossProcess *p)
float electron_mass_c2
Definition: hepunit.py:274
Definition: G4Step.hh:76
G4double MinKinEnergy() const
const G4ParticleDefinition * GetParticleDefinition() const
const G4int n
G4double GetGlobalTime() const
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
G4double Energy(size_t index) const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void AddSecondary(G4Track *aSecondary)
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:321
tuple v
Definition: test.py:18
const G4TouchableHandle & GetTouchableHandle() const
string pname
Definition: eplot.py:33
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:760
void DumpModelList(G4int verb)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
const G4DataVector * Initialise(const G4ParticleDefinition *part, const G4ParticleDefinition *secPart, G4double minSubRange, G4int verb)
static const G4double emax
void SetMasterThread(G4bool val)
Definition: G4VEmModel.hh:718
size_t length() const
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4int size() const
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
void AddEmModel(G4int, G4VEmModel *, const G4Region *region=nullptr)
static G4ProductionCutsTable * GetProductionCutsTable()
G4PhysicsVector * LambdaPhysicsVector(const G4MaterialCutsCouple *)
static G4Positron * Positron()
Definition: G4Positron.cc:94
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:93
G4int NumberOfModels() const
void SetMaxKinEnergy(G4double e)
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
int G4lrint(double ad)
Definition: templates.hh:163
virtual void ProcessDescription(std::ostream &outFile) const
G4double GetProposedKineticEnergy() const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4ProcessManager * GetProcessManager() const
void InitialiseBaseMaterials(G4PhysicsTable *table)
void SetNumberOfSecondaries(G4int totSecondaries)
static G4int Register(const G4String &)
G4double GetLocalEnergyDeposit() const
G4StepPoint * GetPostStepPoint() const
static G4EmParameters * Instance()
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
const G4Element * GetCurrentElement() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double GetSafety() const
void PreparePhysicsTable(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p, G4bool theMaster)
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
void SetFluoFlag(G4bool val)
static constexpr double MeV
Definition: G4SIunits.hh:214
G4bool StorePhysicsTable(const G4String &filename, G4bool ascii=false)
G4TrackStatus GetTrackStatus() const
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const
G4VAtomDeexcitation * AtomDeexcitation()
void SetMinKinEnergy(G4double e)
G4double theInitialNumberOfInteractionLength
Definition: G4VProcess.hh:300
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
void DefineRegParamForEM(G4VEmProcess *) const
void InitializeForPostStep(const G4Track &)
G4ForceCondition
static G4bool RetrievePhysicsTable(G4PhysicsTable *physTable, const G4String &fileName, G4bool ascii)
#define DBL_MAX
Definition: templates.hh:83
static constexpr double keV
Definition: G4SIunits.hh:216
G4double CrossSectionPerVolume(G4double kineticEnergy, const G4MaterialCutsCouple *couple)
G4int GetProcessSubType() const
Definition: G4VProcess.hh:426
void SetPolarAngleLimit(G4double)
Definition: G4VEmModel.hh:767
void clearAndDestroy()
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:437
const G4Material * GetMaterial() const
virtual G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
const std::vector< G4int > * GetCoupleIndexes()
G4ProcessType
const G4Element * GetCurrentElement() const
Definition: G4VEmModel.hh:469