Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4WentzelVIRelModel.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: G4WentzelVIRelModel.cc 74581 2013-10-15 12:03:25Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4WentzelVIRelModel
34 //
35 // Author: V.Ivanchenko
36 //
37 // Creation date: 08.06.2012 from G4WentzelVIRelModel
38 //
39 // Modifications:
40 //
41 // Class Description:
42 //
43 // Implementation of the model of multiple scattering based on
44 // G.Wentzel, Z. Phys. 40 (1927) 590.
45 // H.W.Lewis, Phys Rev 78 (1950) 526.
46 // J.M. Fernandez-Varea et al., NIM B73 (1993) 447.
47 // L.Urban, CERN-OPEN-2006-077.
48 
49 // -------------------------------------------------------------------
50 //
51 
52 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
53 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
54 
55 #include "G4WentzelVIRelModel.hh"
56 #include "G4PhysicalConstants.hh"
57 #include "G4SystemOfUnits.hh"
58 #include "Randomize.hh"
60 #include "G4PhysicsTableHelper.hh"
61 #include "G4ElementVector.hh"
62 #include "G4ProductionCutsTable.hh"
63 #include "G4LossTableManager.hh"
64 #include "G4Pow.hh"
65 #include "G4NistManager.hh"
66 #include "G4Log.hh"
67 #include "G4Exp.hh"
68 
69 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
70 
71 using namespace std;
72 
74  G4VMscModel(nam),
75  numlimit(0.1),
76  currentCouple(0),
77  cosThetaMin(1.0),
78  inside(false),
79  singleScatteringMode(false)
80 {
81  invsqrt12 = 1./sqrt(12.);
82  tlimitminfix = 1.e-6*mm;
83  lowEnergyLimit = 1.0*eV;
84  particle = 0;
85  nelments = 5;
86  xsecn.resize(nelments);
87  prob.resize(nelments);
88  theManager = G4LossTableManager::Instance();
89  fNistManager = G4NistManager::Instance();
90  fG4pow = G4Pow::GetInstance();
91  wokvi = new G4WentzelVIRelXSection();
92 
93  preKinEnergy = tPathLength = zPathLength = lambdaeff = currentRange = xtsec = 0;
94  currentMaterialIndex = 0;
95  cosThetaMax = cosTetMaxNuc = 1.0;
96 
97  fParticleChange = 0;
98  currentCuts = 0;
99  currentMaterial = 0;
100 }
101 
102 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
103 
105 {
106  delete wokvi;
107 }
108 
109 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
110 
112  const G4DataVector& cuts)
113 {
114  // reset parameters
115  SetupParticle(p);
116  currentRange = 0.0;
117 
118  cosThetaMax = cos(PolarAngleLimit());
119  wokvi->Initialise(p, cosThetaMax);
120  /*
121  G4cout << "G4WentzelVIRelModel: " << particle->GetParticleName()
122  << " 1-cos(ThetaLimit)= " << 1 - cosThetaMax
123  << G4endl;
124  */
125  currentCuts = &cuts;
126 
127  // set values of some data members
128  fParticleChange = GetParticleChangeForMSC(p);
129 }
130 
131 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
132 
134  const G4ParticleDefinition* p,
135  G4double kinEnergy,
136  G4double Z, G4double,
137  G4double cutEnergy, G4double)
138 {
139  G4double cross = 0.0;
140  if(p != particle) { SetupParticle(p); }
141  if(kinEnergy < lowEnergyLimit) { return cross; }
142  if(!CurrentCouple()) {
143  G4Exception("G4WentzelVIRelModel::ComputeCrossSectionPerAtom", "em0011",
144  FatalException, " G4MaterialCutsCouple is not defined");
145  return 0.0;
146  }
147  DefineMaterial(CurrentCouple());
148  cosTetMaxNuc = wokvi->SetupKinematic(kinEnergy, currentMaterial);
149  if(cosTetMaxNuc < 1.0) {
150  cosTetMaxNuc = wokvi->SetupTarget(G4lrint(Z), cutEnergy);
151  cross = wokvi->ComputeTransportCrossSectionPerAtom(cosTetMaxNuc);
152  /*
153  if(p->GetParticleName() == "e-")
154  G4cout << "G4WentzelVIRelModel::CS: Z= " << G4int(Z) << " e(MeV)= " << kinEnergy
155  << " 1-cosN= " << 1 - cosTetMaxNuc << " cross(bn)= " << cross/barn
156  << " " << particle->GetParticleName() << G4endl;
157  */
158  }
159  return cross;
160 }
161 
162 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
163 
165 {
166  SetupParticle(track->GetDynamicParticle()->GetDefinition());
167  inside = false;
168 }
169 
170 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
171 
173  const G4Track& track,
174  G4double& currentMinimalStep)
175 {
176  G4double tlimit = currentMinimalStep;
177  const G4DynamicParticle* dp = track.GetDynamicParticle();
178  G4StepPoint* sp = track.GetStep()->GetPreStepPoint();
179  G4StepStatus stepStatus = sp->GetStepStatus();
180  singleScatteringMode = false;
181  //G4cout << "G4WentzelVIRelModel::ComputeTruePathLengthLimit stepStatus= "
182  // << stepStatus << G4endl;
183 
184 
185  // initialisation for each step, lambda may be computed from scratch
186  preKinEnergy = dp->GetKineticEnergy();
187  DefineMaterial(track.GetMaterialCutsCouple());
188  lambdaeff = GetTransportMeanFreePath(particle,preKinEnergy);
189  currentRange = GetRange(particle,preKinEnergy,currentCouple);
190  cosTetMaxNuc = wokvi->SetupKinematic(preKinEnergy, currentMaterial);
191 
192  // extra check for abnormal situation
193  // this check needed to run MSC with eIoni and eBrem inactivated
194  if(tlimit > currentRange) { tlimit = currentRange; }
195 
196  // stop here if small range particle
197  if(inside || tlimit < tlimitminfix) {
198  return ConvertTrueToGeom(tlimit, currentMinimalStep);
199  }
200 
201  // pre step
202  G4double presafety = sp->GetSafety();
203  // far from geometry boundary
204  if(currentRange < presafety) {
205  inside = true;
206  return ConvertTrueToGeom(tlimit, currentMinimalStep);
207  }
208 
209  // compute presafety again if presafety <= 0 and no boundary
210  // i.e. when it is needed for optimization purposes
211  if(stepStatus != fGeomBoundary && presafety < tlimitminfix) {
212  presafety = ComputeSafety(sp->GetPosition(), tlimit);
213  if(currentRange < presafety) {
214  inside = true;
215  return ConvertTrueToGeom(tlimit, currentMinimalStep);
216  }
217  }
218  /*
219  G4cout << "e(MeV)= " << preKinEnergy/MeV
220  << " " << particle->GetParticleName()
221  << " CurLimit(mm)= " << tlimit/mm <<" safety(mm)= " << presafety/mm
222  << " R(mm)= " <<currentRange/mm
223  << " L0(mm^-1)= " << lambdaeff*mm
224  <<G4endl;
225  */
226 
227  // natural limit for high energy
228  G4double rlimit = std::max(facrange*currentRange,
229  0.7*(1.0 - cosTetMaxNuc)*lambdaeff);
230 
231  // low-energy e-
232  if(cosThetaMax > cosTetMaxNuc) {
233  rlimit = std::min(rlimit, facsafety*presafety);
234  }
235 
236  // cut correction
237  G4double rcut = currentCouple->GetProductionCuts()->GetProductionCut(1);
238  //G4cout << "rcut= " << rcut << " rlimit= " << rlimit << " presafety= " << presafety
239  // << " 1-cosThetaMax= " <<1-cosThetaMax << " 1-cosTetMaxNuc= " << 1-cosTetMaxNuc
240  // << G4endl;
241  if(rcut > rlimit) { rlimit = std::min(rlimit, rcut*sqrt(rlimit/rcut)); }
242 
243  if(rlimit < tlimit) { tlimit = rlimit; }
244 
245  tlimit = std::max(tlimit, tlimitminfix);
246 
247  // step limit in infinite media
248  tlimit = std::min(tlimit, 50*currentMaterial->GetRadlen()/facgeom);
249 
250  //compute geomlimit and force few steps within a volume
252  {
253  G4double geomlimit = ComputeGeomLimit(track, presafety, currentRange);
254  tlimit = std::min(tlimit, geomlimit/facgeom);
255  }
256 
257  /*
258  G4cout << particle->GetParticleName() << " e= " << preKinEnergy
259  << " L0= " << lambdaeff << " R= " << currentRange
260  << "tlimit= " << tlimit
261  << " currentMinimalStep= " << currentMinimalStep << G4endl;
262  */
263  return ConvertTrueToGeom(tlimit, currentMinimalStep);
264 }
265 
266 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
267 
269 {
270  tPathLength = truelength;
271  zPathLength = tPathLength;
272 
273  if(lambdaeff > 0.0 && lambdaeff != DBL_MAX) {
274  G4double tau = tPathLength/lambdaeff;
275  //G4cout << "ComputeGeomPathLength: tLength= " << tPathLength
276  // << " Leff= " << lambdaeff << " tau= " << tau << G4endl;
277  // small step
278  if(tau < numlimit) {
279  zPathLength *= (1.0 - 0.5*tau + tau*tau/6.0);
280 
281  // medium step
282  } else {
283  G4double e1 = 0.0;
284  if(currentRange > tPathLength) {
285  e1 = GetEnergy(particle,currentRange-tPathLength,currentCouple);
286  }
287  e1 = 0.5*(e1 + preKinEnergy);
288  cosTetMaxNuc = wokvi->SetupKinematic(e1, currentMaterial);
289  lambdaeff = GetTransportMeanFreePath(particle,e1);
290  zPathLength = lambdaeff*(1.0 - G4Exp(-tPathLength/lambdaeff));
291  }
292  } else { lambdaeff = DBL_MAX; }
293  //G4cout<<"Comp.geom: zLength= "<<zPathLength<<" tLength= "<<tPathLength<<G4endl;
294  return zPathLength;
295 }
296 
297 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
298 
300 {
301  // initialisation of single scattering x-section
302  xtsec = 0.0;
303  cosThetaMin = cosTetMaxNuc;
304 
305  //G4cout << "Step= " << geomStepLength << " Lambda= " << lambdaeff
306  // << " 1-cosThetaMaxNuc= " << 1 - cosTetMaxNuc << G4endl;
307  // pathalogical case
308  if(lambdaeff == DBL_MAX) {
309  singleScatteringMode = true;
310  zPathLength = geomStepLength;
311  tPathLength = geomStepLength;
312  cosThetaMin = 1.0;
313 
314  // normal case
315  } else {
316 
317  // small step use only single scattering
318  static const G4double singleScatLimit = 1.0e-7;
319  if(geomStepLength < lambdaeff*singleScatLimit*(1.0 - cosTetMaxNuc)) {
320  singleScatteringMode = true;
321  cosThetaMin = 1.0;
322  zPathLength = geomStepLength;
323  tPathLength = geomStepLength;
324 
325  // step defined by transportation
326  } else if(geomStepLength != zPathLength) {
327 
328  // step defined by transportation
329  zPathLength = geomStepLength;
330  G4double tau = geomStepLength/lambdaeff;
331  tPathLength = zPathLength*(1.0 + 0.5*tau + tau*tau/3.0);
332 
333  // energy correction for a big step
334  if(tau > numlimit) {
335  G4double e1 = 0.0;
336  if(currentRange > tPathLength) {
337  e1 = GetEnergy(particle,currentRange-tPathLength,currentCouple);
338  }
339  e1 = 0.5*(e1 + preKinEnergy);
340  cosTetMaxNuc = wokvi->SetupKinematic(e1, currentMaterial);
341  lambdaeff = GetTransportMeanFreePath(particle,e1);
342  tau = zPathLength/lambdaeff;
343 
344  if(tau < 0.999999) { tPathLength = -lambdaeff*G4Log(1.0 - tau); }
345  else { tPathLength = currentRange; }
346  }
347  }
348  }
349 
350  // check of step length
351  // define threshold angle between single and multiple scattering
352  if(!singleScatteringMode) { cosThetaMin = 1.0 - 1.5*tPathLength/lambdaeff; }
353 
354  // recompute transport cross section - do not change energy
355  // anymore - cannot be applied for big steps
356  if(cosThetaMin > cosTetMaxNuc) {
357 
358  // new computation
359  G4double cross = ComputeXSectionPerVolume();
360  //G4cout << "%%%% cross= " << cross << " xtsec= " << xtsec << G4endl;
361  if(cross <= 0.0) {
362  singleScatteringMode = true;
363  tPathLength = zPathLength;
364  lambdaeff = DBL_MAX;
365  } else if(xtsec > 0.0) {
366 
367  lambdaeff = 1./cross;
368  G4double tau = zPathLength*cross;
369  if(tau < numlimit) { tPathLength = zPathLength*(1.0 + 0.5*tau + tau*tau/3.0); }
370  else if(tau < 0.999999) { tPathLength = -lambdaeff*G4Log(1.0 - tau); }
371  else { tPathLength = currentRange; }
372 
373  if(tPathLength > currentRange) { tPathLength = currentRange; }
374  }
375  }
376 
377  /*
378  G4cout <<"Comp.true: zLength= "<<zPathLength<<" tLength= "<<tPathLength
379  <<" Leff(mm)= "<<lambdaeff/mm<<" sig0(1/mm)= " << xtsec <<G4endl;
380  G4cout << particle->GetParticleName() << " 1-cosThetaMin= " << 1-cosThetaMin
381  << " 1-cosTetMaxNuc= " << 1-cosTetMaxNuc
382  << " e(MeV)= " << preKinEnergy/MeV << " " << singleScatteringMode << G4endl;
383  */
384  return tPathLength;
385 }
386 
387 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
388 
391  G4double safety)
392 {
393  fDisplacement.set(0.0,0.0,0.0);
394  //G4cout << "!##! G4WentzelVIRelModel::SampleScattering for "
395  // << particle->GetParticleName() << G4endl;
396 
397  // ignore scattering for zero step length and energy below the limit
398  if(preKinEnergy < lowEnergyLimit || tPathLength <= 0.0)
399  { return fDisplacement; }
400 
401  G4double invlambda = 0.0;
402  if(lambdaeff < DBL_MAX) { invlambda = 0.5/lambdaeff; }
403 
404  // use average kinetic energy over the step
405  G4double cut = (*currentCuts)[currentMaterialIndex];
406  /*
407  G4cout <<"SampleScat: E0(MeV)= "<< preKinEnergy/MeV
408  << " Leff= " << lambdaeff <<" sig0(1/mm)= " << xtsec
409  << " xmsc= " << tPathLength*invlambda << " safety= " << safety << G4endl;
410  */
411 
412  // step limit due msc
413  G4double x0 = tPathLength;
414  // const G4double thinlimit = 0.5;
415  static const G4double thinlimit = 0.1;
416  G4int nMscSteps = 1;
417  // large scattering angle case - two step approach
418  if(tPathLength*invlambda > thinlimit && safety > tlimitminfix) {
419  x0 *= 0.5;
420  nMscSteps = 2;
421  }
422 
423  // step limit due to single scattering
424  G4double x1 = 2*tPathLength;
425  if(0.0 < xtsec) { x1 = -G4Log(G4UniformRand())/xtsec; }
426 
427  const G4ElementVector* theElementVector =
428  currentMaterial->GetElementVector();
429  G4int nelm = currentMaterial->GetNumberOfElements();
430 
431  // geometry
432  G4double sint, cost, phi;
433  G4ThreeVector temp(0.0,0.0,1.0);
434 
435  // current position and direction relative to the end point
436  // because of magnetic field geometry is computed relatively to the
437  // end point of the step
438  G4ThreeVector dir(0.0,0.0,1.0);
439  fDisplacement.set(0.0,0.0,-zPathLength);
440  G4double mscfac = zPathLength/tPathLength;
441 
442  // start a loop
443  G4double x2 = x0;
444  G4double step, z;
445  G4bool singleScat;
446  /*
447  G4cout << "Start of the loop x1(mm)= " << x1 << " x2(mm)= " << x2
448  << " 1-cost1= " << 1 - cosThetaMin << " " << singleScatteringMode
449  << " xtsec= " << xtsec << G4endl;
450  */
451  do {
452 
453  // single scattering case
454  if(x1 < x2) {
455  step = x1;
456  singleScat = true;
457  } else {
458  step = x2;
459  singleScat = false;
460  }
461 
462  // new position
463  fDisplacement += step*mscfac*dir;
464 
465  if(singleScat) {
466 
467  // select element
468  G4int i = 0;
469  if(nelm > 1) {
470  G4double qsec = G4UniformRand()*xtsec;
471  for (; i<nelm; ++i) { if(xsecn[i] >= qsec) { break; } }
472  }
473  G4double cosTetM =
474  wokvi->SetupTarget(G4lrint((*theElementVector)[i]->GetZ()), cut);
475  //G4cout << "!!! " << cosThetaMin << " " << cosTetM << " " << prob[i] << G4endl;
476  temp = wokvi->SampleSingleScattering(cosThetaMin, cosTetM, prob[i]);
477 
478  // direction is changed
479  temp.rotateUz(dir);
480  dir = temp;
481 
482  // new proposed step length
483  x1 = -G4Log(G4UniformRand())/xtsec;
484  x2 -= step;
485  if(x2 <= 0.0) { --nMscSteps; }
486 
487  // multiple scattering
488  } else {
489  --nMscSteps;
490  x1 -= step;
491  x2 = x0;
492 
493  // for multiple scattering x0 should be used as a step size
494  if(!singleScatteringMode) {
495  G4double z0 = x0*invlambda;
496 
497  // correction to keep first moment
498 
499  // sample z in interval 0 - 1
500  if(z0 > 5.0) { z = G4UniformRand(); }
501  else {
502  G4double zzz = 0.0;
503  if(z0 > 0.01) { zzz = G4Exp(-1.0/z0); }
504  z = -z0*G4Log(1.0 - (1.0 - zzz)*G4UniformRand());
505  // /(1.0 - (1.0/z0 + 1.0)*zzz);
506  }
507 
508  cost = 1.0 - 2.0*z/*factCM*/;
509  if(cost > 1.0) { cost = 1.0; }
510  else if(cost < -1.0) { cost =-1.0; }
511  sint = sqrt((1.0 - cost)*(1.0 + cost));
512  phi = twopi*G4UniformRand();
513  G4double vx1 = sint*cos(phi);
514  G4double vy1 = sint*sin(phi);
515 
516  // change direction
517  temp.set(vx1,vy1,cost);
518  temp.rotateUz(dir);
519  dir = temp;
520 
521  // lateral displacement
522  if (latDisplasment && x0 > tlimitminfix) {
523  G4double rms = invsqrt12*sqrt(2*z0);
524  G4double r = x0*mscfac;
525  G4double dx = r*(0.5*vx1 + rms*G4RandGauss::shoot(0.0,1.0));
526  G4double dy = r*(0.5*vy1 + rms*G4RandGauss::shoot(0.0,1.0));
527  G4double dz;
528  G4double d = r*r - dx*dx - dy*dy;
529  if(d >= 0.0) { dz = sqrt(d) - r; }
530  else { dx = dy = dz = 0.0; }
531 
532  // change position
533  temp.set(dx,dy,dz);
534  temp.rotateUz(dir);
535  fDisplacement += temp;
536  }
537  }
538  }
539  } while (0 < nMscSteps);
540 
541  dir.rotateUz(oldDirection);
542 
543  //G4cout << "G4WentzelVIRelModel sampling of scattering is done" << G4endl;
544  // end of sampling -------------------------------
545 
546  fParticleChange->ProposeMomentumDirection(dir);
547 
548  // lateral displacement
549  fDisplacement.rotateUz(oldDirection);
550 
551  /*
552  G4cout << " r(mm)= " << fDisplacement.mag()
553  << " safety= " << safety
554  << " trueStep(mm)= " << tPathLength
555  << " geomStep(mm)= " << zPathLength
556  << " x= " << fDisplacement.x()
557  << " y= " << fDisplacement.y()
558  << " z= " << fDisplacement.z()
559  << G4endl;
560  */
561 
562  //G4cout<< "G4WentzelVIRelModel::SampleScattering end NewDir= " << dir<< G4endl;
563  return fDisplacement;
564 }
565 
566 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
567 
568 G4double G4WentzelVIRelModel::ComputeXSectionPerVolume()
569 {
570  // prepare recomputation of x-sections
571  const G4ElementVector* theElementVector = currentMaterial->GetElementVector();
572  const G4double* theAtomNumDensityVector =
573  currentMaterial->GetVecNbOfAtomsPerVolume();
574  G4int nelm = currentMaterial->GetNumberOfElements();
575  if(nelm > nelments) {
576  nelments = nelm;
577  xsecn.resize(nelm);
578  prob.resize(nelm);
579  }
580  G4double cut = (*currentCuts)[currentMaterialIndex];
581  // cosTetMaxNuc = wokvi->GetCosThetaNuc();
582 
583  // check consistency
584  xtsec = 0.0;
585  if(cosTetMaxNuc > cosThetaMin) { return 0.0; }
586 
587  // loop over elements
588  G4double xs = 0.0;
589  for (G4int i=0; i<nelm; ++i) {
590  G4double costm =
591  wokvi->SetupTarget(G4lrint((*theElementVector)[i]->GetZ()), cut);
592  G4double density = theAtomNumDensityVector[i];
593 
594  G4double esec = 0.0;
595  if(costm < cosThetaMin) {
596 
597  // recompute the transport x-section
598  if(1.0 > cosThetaMin) {
599  xs += density*wokvi->ComputeTransportCrossSectionPerAtom(cosThetaMin);
600  }
601  // recompute the total x-section
602  G4double nucsec = wokvi->ComputeNuclearCrossSection(cosThetaMin, costm);
603  esec = wokvi->ComputeElectronCrossSection(cosThetaMin, costm);
604  nucsec += esec;
605  if(nucsec > 0.0) { esec /= nucsec; }
606  xtsec += nucsec*density;
607  }
608  xsecn[i] = xtsec;
609  prob[i] = esec;
610  //G4cout << i << " xs= " << xs << " xtsec= " << xtsec
611  // << " 1-cosThetaMin= " << 1-cosThetaMin
612  // << " 1-cosTetMaxNuc2= " <<1-cosTetMaxNuc2<< G4endl;
613  }
614 
615  //G4cout << "ComputeXS result: xsec(1/mm)= " << xs
616  // << " txsec(1/mm)= " << xtsec <<G4endl;
617  return xs;
618 }
619 
620 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
void set(double x, double y, double z)
static G4Pow * GetInstance()
Definition: G4Pow.cc:53
G4double facgeom
Definition: G4VMscModel.hh:177
void StartTracking(G4Track *)
ThreeVector shoot(const G4int Ap, const G4int Af)
G4double ComputeNuclearCrossSection(G4double CosThetaMin, G4double CosThetaMax)
virtual G4double ComputeTruePathLengthLimit(const G4Track &track, G4double &currentMinimalStep)
G4ThreeVector SampleSingleScattering(G4double CosThetaMin, G4double CosThetaMax, G4double elecRatio=0.0)
static G4LossTableManager * Instance()
virtual G4ThreeVector & SampleScattering(const G4ThreeVector &, G4double safety)
std::vector< G4Element * > G4ElementVector
G4double GetKineticEnergy() const
const G4DynamicParticle * GetDynamicParticle() const
G4double facrange
Definition: G4VMscModel.hh:176
G4double z
Definition: TRTMaterials.hh:39
G4double GetProductionCut(G4int index) const
const char * p
Definition: xmltok.h:285
G4StepStatus GetStepStatus() const
G4double ConvertTrueToGeom(G4double &tLength, G4double &gLength)
Definition: G4VMscModel.hh:246
G4WentzelVIRelModel(const G4String &nam="WentzelVIUni")
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4bool latDisplasment
Definition: G4VMscModel.hh:189
virtual G4double ComputeGeomPathLength(G4double truePathLength)
G4ParticleDefinition * GetDefinition() const
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double KineticEnergy, G4double AtomicNumber, G4double AtomicWeight=0., G4double cut=DBL_MAX, G4double emax=DBL_MAX)
G4double ComputeTransportCrossSectionPerAtom(G4double CosThetaMax)
const G4Step * GetStep() const
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:188
int G4int
Definition: G4Types.hh:78
static G4NistManager * Instance()
G4double SetupKinematic(G4double kinEnergy, const G4Material *mat)
G4double SetupTarget(G4int Z, G4double cut=DBL_MAX)
G4StepStatus
Definition: G4StepStatus.hh:49
G4StepPoint * GetPreStepPoint() const
G4double GetEnergy(const G4ParticleDefinition *part, G4double range, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.hh:308
G4double density
Definition: TRTMaterials.hh:39
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:204
const G4MaterialCutsCouple * CurrentCouple() const
Definition: G4VEmModel.hh:426
#define G4UniformRand()
Definition: Randomize.hh:87
const G4ThreeVector & GetPosition() const
G4double GetRange(const G4ParticleDefinition *part, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.hh:288
G4ParticleChangeForMSC * GetParticleChangeForMSC(const G4ParticleDefinition *p=0)
Definition: G4VMscModel.cc:89
bool G4bool
Definition: G4Types.hh:79
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:72
G4double ComputeGeomLimit(const G4Track &, G4double &presafety, G4double limit)
Definition: G4VMscModel.hh:256
G4ThreeVector fDisplacement
Definition: G4VMscModel.hh:185
G4double GetTransportMeanFreePath(const G4ParticleDefinition *part, G4double kinEnergy)
Definition: G4VMscModel.hh:345
G4double GetRadlen() const
Definition: G4Material.hh:218
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4double ComputeSafety(const G4ThreeVector &position, G4double limit)
Definition: G4VMscModel.hh:238
G4double G4Log(G4double x)
Definition: G4Log.hh:227
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:180
void ProposeMomentumDirection(const G4ThreeVector &Pfinal)
int G4lrint(double ad)
Definition: templates.hh:163
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4double facsafety
Definition: G4VMscModel.hh:178
virtual G4double ComputeTrueStepLength(G4double geomStepLength)
G4double PolarAngleLimit() const
Definition: G4VEmModel.hh:620
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double GetSafety() const
G4MscStepLimitType steppingAlgorithm
Definition: G4VMscModel.hh:186
size_t GetNumberOfElements() const
Definition: G4Material.hh:184
G4double ComputeElectronCrossSection(G4double CosThetaMin, G4double CosThetaMax)
void Initialise(const G4ParticleDefinition *, G4double CosThetaLim)
double G4double
Definition: G4Types.hh:76
G4ProductionCuts * GetProductionCuts() const
#define DBL_MAX
Definition: templates.hh:83