Geant4_10
G4LowEnergyPolarizedCompton.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 //
27 // $Id$
28 // GEANT4 tag $Name: $
29 //
30 // ------------------------------------------------------------
31 // GEANT 4 class implementation file
32 // CERN Geneva Switzerland
33 //
34 
35 // --------- G4LowEnergyPolarizedCompton class -----
36 //
37 // by G.Depaola & F.Longo (21 may 2001)
38 //
39 // 21 May 2001 - MGP Modified to inherit from G4VDiscreteProcess
40 // Applies same algorithm as LowEnergyCompton
41 // if the incoming photon is not polarised
42 // Temporary protection to avoid crash in the case
43 // of polarisation || incident photon direction
44 //
45 // 17 October 2001 - F.Longo - Revised according to a design iteration
46 //
47 // 21 February 2002 - F.Longo Revisions with A.Zoglauer and G.Depaola
48 // - better description of parallelism
49 // - system of ref change method improved
50 // 22 January 2003 - V.Ivanchenko Cut per region
51 // 24 April 2003 - V.Ivanchenko Cut per region mfpt
52 //
53 //
54 // ************************************************************
55 //
56 // Corrections by Rui Curado da Silva (2000)
57 // New Implementation by G.Depaola & F.Longo
58 //
59 // - sampling of phi
60 // - polarization of scattered photon
61 //
62 // --------------------------------------------------------------
63 
65 #include "Randomize.hh"
66 #include "G4PhysicalConstants.hh"
67 #include "G4SystemOfUnits.hh"
68 #include "G4ParticleDefinition.hh"
69 #include "G4Track.hh"
70 #include "G4Step.hh"
71 #include "G4ForceCondition.hh"
72 #include "G4Gamma.hh"
73 #include "G4Electron.hh"
74 #include "G4DynamicParticle.hh"
75 #include "G4VParticleChange.hh"
76 #include "G4ThreeVector.hh"
79 #include "G4RDVEMDataSet.hh"
81 #include "G4RDVDataSetAlgorithm.hh"
83 #include "G4RDVRangeTest.hh"
84 #include "G4RDRangeTest.hh"
85 #include "G4MaterialCutsCouple.hh"
86 
87 // constructor
88 
90  : G4VDiscreteProcess(processName),
91  lowEnergyLimit (250*eV), // initialization
92  highEnergyLimit(100*GeV),
93  intrinsicLowEnergyLimit(10*eV),
94  intrinsicHighEnergyLimit(100*GeV)
95 {
96  if (lowEnergyLimit < intrinsicLowEnergyLimit ||
97  highEnergyLimit > intrinsicHighEnergyLimit)
98  {
99  G4Exception("G4LowEnergyPolarizedCompton::G4LowEnergyPolarizedCompton()",
100  "OutOfRange", FatalException,
101  "Energy outside intrinsic process validity range!");
102  }
103 
104  crossSectionHandler = new G4RDCrossSectionHandler;
105 
106 
107  G4RDVDataSetAlgorithm* scatterInterpolation = new G4RDLogLogInterpolation;
108  G4String scatterFile = "comp/ce-sf-";
109  scatterFunctionData = new G4RDCompositeEMDataSet(scatterInterpolation,1.,1.);
110  scatterFunctionData->LoadData(scatterFile);
111 
112  meanFreePathTable = 0;
113 
114  rangeTest = new G4RDRangeTest;
115 
116  // For Doppler broadening
117  shellData.SetOccupancyData();
118 
119 
120  if (verboseLevel > 0)
121  {
122  G4cout << GetProcessName() << " is created " << G4endl
123  << "Energy range: "
124  << lowEnergyLimit / keV << " keV - "
125  << highEnergyLimit / GeV << " GeV"
126  << G4endl;
127  }
128 }
129 
130 // destructor
131 
133 {
134  delete meanFreePathTable;
135  delete crossSectionHandler;
136  delete scatterFunctionData;
137  delete rangeTest;
138 }
139 
140 
142 {
143 
144  crossSectionHandler->Clear();
145  G4String crossSectionFile = "comp/ce-cs-";
146  crossSectionHandler->LoadData(crossSectionFile);
147  delete meanFreePathTable;
148  meanFreePathTable = crossSectionHandler->BuildMeanFreePathForMaterials();
149 
150  // For Doppler broadening
151  G4String file = "/doppler/shell-doppler";
152  shellData.LoadData(file);
153 
154 }
155 
157  const G4Step& aStep)
158 {
159  // The scattered gamma energy is sampled according to Klein - Nishina formula.
160  // The random number techniques of Butcher & Messel are used (Nuc Phys 20(1960),15).
161  // GEANT4 internal units
162  //
163  // Note : Effects due to binding of atomic electrons are negliged.
164 
165  aParticleChange.Initialize(aTrack);
166 
167  // Dynamic particle quantities
168  const G4DynamicParticle* incidentPhoton = aTrack.GetDynamicParticle();
169  G4double gammaEnergy0 = incidentPhoton->GetKineticEnergy();
170  G4ThreeVector gammaPolarization0 = incidentPhoton->GetPolarization();
171 
172  // gammaPolarization0 = gammaPolarization0.unit(); //
173 
174  // Protection: a polarisation parallel to the
175  // direction causes problems;
176  // in that case find a random polarization
177 
178  G4ThreeVector gammaDirection0 = incidentPhoton->GetMomentumDirection();
179  // ---- MGP ---- Next two lines commented out to remove compilation warnings
180  // G4double scalarproduct = gammaPolarization0.dot(gammaDirection0);
181  // G4double angle = gammaPolarization0.angle(gammaDirection0);
182 
183  // Make sure that the polarization vector is perpendicular to the
184  // gamma direction. If not
185 
186  if(!(gammaPolarization0.isOrthogonal(gammaDirection0, 1e-6))||(gammaPolarization0.mag()==0))
187  { // only for testing now
188  gammaPolarization0 = GetRandomPolarization(gammaDirection0);
189  }
190  else
191  {
192  if ( gammaPolarization0.howOrthogonal(gammaDirection0) != 0)
193  {
194  gammaPolarization0 = GetPerpendicularPolarization(gammaDirection0, gammaPolarization0);
195  }
196  }
197 
198  // End of Protection
199 
200  // Within energy limit?
201 
202  if(gammaEnergy0 <= lowEnergyLimit)
203  {
207  return G4VDiscreteProcess::PostStepDoIt(aTrack,aStep);
208  }
209 
210  G4double E0_m = gammaEnergy0 / electron_mass_c2 ;
211 
212  // Select randomly one element in the current material
213 
214  const G4MaterialCutsCouple* couple = aTrack.GetMaterialCutsCouple();
215  G4int Z = crossSectionHandler->SelectRandomAtom(couple,gammaEnergy0);
216 
217  // Sample the energy and the polarization of the scattered photon
218 
219  G4double epsilon, epsilonSq, onecost, sinThetaSqr, greject ;
220 
221  G4double epsilon0 = 1./(1. + 2*E0_m);
222  G4double epsilon0Sq = epsilon0*epsilon0;
223  G4double alpha1 = - std::log(epsilon0);
224  G4double alpha2 = 0.5*(1.- epsilon0Sq);
225 
226  G4double wlGamma = h_Planck*c_light/gammaEnergy0;
227  G4double gammaEnergy1;
228  G4ThreeVector gammaDirection1;
229 
230  do {
231  if ( alpha1/(alpha1+alpha2) > G4UniformRand() )
232  {
233  epsilon = std::exp(-alpha1*G4UniformRand());
234  epsilonSq = epsilon*epsilon;
235  }
236  else
237  {
238  epsilonSq = epsilon0Sq + (1.- epsilon0Sq)*G4UniformRand();
239  epsilon = std::sqrt(epsilonSq);
240  }
241 
242  onecost = (1.- epsilon)/(epsilon*E0_m);
243  sinThetaSqr = onecost*(2.-onecost);
244 
245  // Protection
246  if (sinThetaSqr > 1.)
247  {
248  if (verboseLevel>0) G4cout
249  << " -- Warning -- G4LowEnergyPolarizedCompton::PostStepDoIt "
250  << "sin(theta)**2 = "
251  << sinThetaSqr
252  << "; set to 1"
253  << G4endl;
254  sinThetaSqr = 1.;
255  }
256  if (sinThetaSqr < 0.)
257  {
258  if (verboseLevel>0) G4cout
259  << " -- Warning -- G4LowEnergyPolarizedCompton::PostStepDoIt "
260  << "sin(theta)**2 = "
261  << sinThetaSqr
262  << "; set to 0"
263  << G4endl;
264  sinThetaSqr = 0.;
265  }
266  // End protection
267 
268  G4double x = std::sqrt(onecost/2.) / (wlGamma/cm);;
269  G4double scatteringFunction = scatterFunctionData->FindValue(x,Z-1);
270  greject = (1. - epsilon*sinThetaSqr/(1.+ epsilonSq))*scatteringFunction;
271 
272  } while(greject < G4UniformRand()*Z);
273 
274 
275  // ****************************************************
276  // Phi determination
277  // ****************************************************
278 
279  G4double phi = SetPhi(epsilon,sinThetaSqr);
280 
281  //
282  // scattered gamma angles. ( Z - axis along the parent gamma)
283  //
284 
285  G4double cosTheta = 1. - onecost;
286 
287  // Protection
288 
289  if (cosTheta > 1.)
290  {
291  if (verboseLevel>0) G4cout
292  << " -- Warning -- G4LowEnergyPolarizedCompton::PostStepDoIt "
293  << "cosTheta = "
294  << cosTheta
295  << "; set to 1"
296  << G4endl;
297  cosTheta = 1.;
298  }
299  if (cosTheta < -1.)
300  {
301  if (verboseLevel>0) G4cout
302  << " -- Warning -- G4LowEnergyPolarizedCompton::PostStepDoIt "
303  << "cosTheta = "
304  << cosTheta
305  << "; set to -1"
306  << G4endl;
307  cosTheta = -1.;
308  }
309  // End protection
310 
311 
312  G4double sinTheta = std::sqrt (sinThetaSqr);
313 
314  // Protection
315  if (sinTheta > 1.)
316  {
317  if (verboseLevel>0) G4cout
318  << " -- Warning -- G4LowEnergyPolarizedCompton::PostStepDoIt "
319  << "sinTheta = "
320  << sinTheta
321  << "; set to 1"
322  << G4endl;
323  sinTheta = 1.;
324  }
325  if (sinTheta < -1.)
326  {
327  if (verboseLevel>0) G4cout
328  << " -- Warning -- G4LowEnergyPolarizedCompton::PostStepDoIt "
329  << "sinTheta = "
330  << sinTheta
331  << "; set to -1"
332  << G4endl;
333  sinTheta = -1.;
334  }
335  // End protection
336 
337 
338  G4double dirx = sinTheta*std::cos(phi);
339  G4double diry = sinTheta*std::sin(phi);
340  G4double dirz = cosTheta ;
341 
342 
343  // oneCosT , eom
344 
345 
346 
347  // Doppler broadening - Method based on:
348  // Y. Namito, S. Ban and H. Hirayama,
349  // "Implementation of the Doppler Broadening of a Compton-Scattered Photon Into the EGS4 Code"
350  // NIM A 349, pp. 489-494, 1994
351 
352  // Maximum number of sampling iterations
353 
354  G4int maxDopplerIterations = 1000;
355  G4double bindingE = 0.;
356  G4double photonEoriginal = epsilon * gammaEnergy0;
357  G4double photonE = -1.;
358  G4int iteration = 0;
359  G4double eMax = gammaEnergy0;
360 
361  do
362  {
363  iteration++;
364  // Select shell based on shell occupancy
365  G4int shell = shellData.SelectRandomShell(Z);
366  bindingE = shellData.BindingEnergy(Z,shell);
367 
368  eMax = gammaEnergy0 - bindingE;
369 
370  // Randomly sample bound electron momentum (memento: the data set is in Atomic Units)
371  G4double pSample = profileData.RandomSelectMomentum(Z,shell);
372  // Rescale from atomic units
373  G4double pDoppler = pSample * fine_structure_const;
374  G4double pDoppler2 = pDoppler * pDoppler;
375  G4double var2 = 1. + onecost * E0_m;
376  G4double var3 = var2*var2 - pDoppler2;
377  G4double var4 = var2 - pDoppler2 * cosTheta;
378  G4double var = var4*var4 - var3 + pDoppler2 * var3;
379  if (var > 0.)
380  {
381  G4double varSqrt = std::sqrt(var);
382  G4double scale = gammaEnergy0 / var3;
383  // Random select either root
384  if (G4UniformRand() < 0.5) photonE = (var4 - varSqrt) * scale;
385  else photonE = (var4 + varSqrt) * scale;
386  }
387  else
388  {
389  photonE = -1.;
390  }
391  } while ( iteration <= maxDopplerIterations &&
392  (photonE < 0. || photonE > eMax || photonE < eMax*G4UniformRand()) );
393 
394  // End of recalculation of photon energy with Doppler broadening
395  // Revert to original if maximum number of iterations threshold has been reached
396  if (iteration >= maxDopplerIterations)
397  {
398  photonE = photonEoriginal;
399  bindingE = 0.;
400  }
401 
402  gammaEnergy1 = photonE;
403 
404  // G4cout << "--> PHOTONENERGY1 = " << photonE/keV << G4endl;
405 
406 
408 
409 
410 
411 
412  //
413  // update G4VParticleChange for the scattered photon
414  //
415 
416  // gammaEnergy1 = epsilon*gammaEnergy0;
417 
418 
419  // New polarization
420 
421  G4ThreeVector gammaPolarization1 = SetNewPolarization(epsilon,
422  sinThetaSqr,
423  phi,
424  cosTheta);
425 
426  // Set new direction
427  G4ThreeVector tmpDirection1( dirx,diry,dirz );
428  gammaDirection1 = tmpDirection1;
429 
430  // Change reference frame.
431 
432  SystemOfRefChange(gammaDirection0,gammaDirection1,
433  gammaPolarization0,gammaPolarization1);
434 
435  if (gammaEnergy1 > 0.)
436  {
437  aParticleChange.ProposeEnergy( gammaEnergy1 ) ;
438  aParticleChange.ProposeMomentumDirection( gammaDirection1 );
439  aParticleChange.ProposePolarization( gammaPolarization1 );
440  }
441  else
442  {
445  }
446 
447  //
448  // kinematic of the scattered electron
449  //
450 
451  G4double ElecKineEnergy = gammaEnergy0 - gammaEnergy1 -bindingE;
452 
453 
454  // Generate the electron only if with large enough range w.r.t. cuts and safety
455 
456  G4double safety = aStep.GetPostStepPoint()->GetSafety();
457 
458 
459  if (rangeTest->Escape(G4Electron::Electron(),couple,ElecKineEnergy,safety))
460  {
461  G4double ElecMomentum = std::sqrt(ElecKineEnergy*(ElecKineEnergy+2.*electron_mass_c2));
462  G4ThreeVector ElecDirection((gammaEnergy0 * gammaDirection0 -
463  gammaEnergy1 * gammaDirection1) * (1./ElecMomentum));
464  G4DynamicParticle* electron = new G4DynamicParticle (G4Electron::Electron(),ElecDirection.unit(),ElecKineEnergy) ;
466  aParticleChange.AddSecondary(electron);
467  // aParticleChange.ProposeLocalEnergyDeposit(0.);
469  }
470  else
471  {
473  aParticleChange.ProposeLocalEnergyDeposit(ElecKineEnergy+bindingE);
474  }
475 
476  return G4VDiscreteProcess::PostStepDoIt( aTrack, aStep);
477 
478 }
479 
480 
481 G4double G4LowEnergyPolarizedCompton::SetPhi(G4double energyRate,
482  G4double sinSqrTh)
483 {
484  G4double rand1;
485  G4double rand2;
486  G4double phiProbability;
487  G4double phi;
488  G4double a, b;
489 
490  do
491  {
492  rand1 = G4UniformRand();
493  rand2 = G4UniformRand();
494  phiProbability=0.;
495  phi = twopi*rand1;
496 
497  a = 2*sinSqrTh;
498  b = energyRate + 1/energyRate;
499 
500  phiProbability = 1 - (a/b)*(std::cos(phi)*std::cos(phi));
501 
502 
503 
504  }
505  while ( rand2 > phiProbability );
506  return phi;
507 }
508 
509 
510 G4ThreeVector G4LowEnergyPolarizedCompton::SetPerpendicularVector(G4ThreeVector& a)
511 {
512  G4double dx = a.x();
513  G4double dy = a.y();
514  G4double dz = a.z();
515  G4double x = dx < 0.0 ? -dx : dx;
516  G4double y = dy < 0.0 ? -dy : dy;
517  G4double z = dz < 0.0 ? -dz : dz;
518  if (x < y) {
519  return x < z ? G4ThreeVector(-dy,dx,0) : G4ThreeVector(0,-dz,dy);
520  }else{
521  return y < z ? G4ThreeVector(dz,0,-dx) : G4ThreeVector(-dy,dx,0);
522  }
523 }
524 
525 G4ThreeVector G4LowEnergyPolarizedCompton::GetRandomPolarization(G4ThreeVector& direction0)
526 {
527  G4ThreeVector d0 = direction0.unit();
528  G4ThreeVector a1 = SetPerpendicularVector(d0); //different orthogonal
529  G4ThreeVector a0 = a1.unit(); // unit vector
530 
531  G4double rand1 = G4UniformRand();
532 
533  G4double angle = twopi*rand1; // random polar angle
534  G4ThreeVector b0 = d0.cross(a0); // cross product
535 
537 
538  c.setX(std::cos(angle)*(a0.x())+std::sin(angle)*b0.x());
539  c.setY(std::cos(angle)*(a0.y())+std::sin(angle)*b0.y());
540  c.setZ(std::cos(angle)*(a0.z())+std::sin(angle)*b0.z());
541 
542  G4ThreeVector c0 = c.unit();
543 
544  return c0;
545 
546 }
547 
548 
549 G4ThreeVector G4LowEnergyPolarizedCompton::GetPerpendicularPolarization
550 (const G4ThreeVector& gammaDirection, const G4ThreeVector& gammaPolarization) const
551 {
552 
553  //
554  // The polarization of a photon is always perpendicular to its momentum direction.
555  // Therefore this function removes those vector component of gammaPolarization, which
556  // points in direction of gammaDirection
557  //
558  // Mathematically we search the projection of the vector a on the plane E, where n is the
559  // plains normal vector.
560  // The basic equation can be found in each geometry book (e.g. Bronstein):
561  // p = a - (a o n)/(n o n)*n
562 
563  return gammaPolarization - gammaPolarization.dot(gammaDirection)/gammaDirection.dot(gammaDirection) * gammaDirection;
564 }
565 
566 
567 G4ThreeVector G4LowEnergyPolarizedCompton::SetNewPolarization(G4double epsilon,
568  G4double sinSqrTh,
569  G4double phi,
570  G4double costheta)
571 {
572  G4double rand1;
573  G4double rand2;
574  G4double cosPhi = std::cos(phi);
575  G4double sinPhi = std::sin(phi);
576  G4double sinTheta = std::sqrt(sinSqrTh);
577  G4double cosSqrPhi = cosPhi*cosPhi;
578  // G4double cossqrth = 1.-sinSqrTh;
579  // G4double sinsqrphi = sinPhi*sinPhi;
580  G4double normalisation = std::sqrt(1. - cosSqrPhi*sinSqrTh);
581 
582 
583  // Determination of Theta
584 
585  // ---- MGP ---- Commented out the following 3 lines to avoid compilation
586  // warnings (unused variables)
587  // G4double thetaProbability;
588  G4double theta;
589  // G4double a, b;
590  // G4double cosTheta;
591 
592  /*
593 
594  depaola method
595 
596  do
597  {
598  rand1 = G4UniformRand();
599  rand2 = G4UniformRand();
600  thetaProbability=0.;
601  theta = twopi*rand1;
602  a = 4*normalisation*normalisation;
603  b = (epsilon + 1/epsilon) - 2;
604  thetaProbability = (b + a*std::cos(theta)*std::cos(theta))/(a+b);
605  cosTheta = std::cos(theta);
606  }
607  while ( rand2 > thetaProbability );
608 
609  G4double cosBeta = cosTheta;
610 
611  */
612 
613 
614  // Dan Xu method (IEEE TNS, 52, 1160 (2005))
615 
616  rand1 = G4UniformRand();
617  rand2 = G4UniformRand();
618 
619  if (rand1<(epsilon+1.0/epsilon-2)/(2.0*(epsilon+1.0/epsilon)-4.0*sinSqrTh*cosSqrPhi))
620  {
621  if (rand2<0.5)
622  theta = pi/2.0;
623  else
624  theta = 3.0*pi/2.0;
625  }
626  else
627  {
628  if (rand2<0.5)
629  theta = 0;
630  else
631  theta = pi;
632  }
633  G4double cosBeta = std::cos(theta);
634  G4double sinBeta = std::sqrt(1-cosBeta*cosBeta);
635 
636  G4ThreeVector gammaPolarization1;
637 
638  G4double xParallel = normalisation*cosBeta;
639  G4double yParallel = -(sinSqrTh*cosPhi*sinPhi)*cosBeta/normalisation;
640  G4double zParallel = -(costheta*sinTheta*cosPhi)*cosBeta/normalisation;
641  G4double xPerpendicular = 0.;
642  G4double yPerpendicular = (costheta)*sinBeta/normalisation;
643  G4double zPerpendicular = -(sinTheta*sinPhi)*sinBeta/normalisation;
644 
645  G4double xTotal = (xParallel + xPerpendicular);
646  G4double yTotal = (yParallel + yPerpendicular);
647  G4double zTotal = (zParallel + zPerpendicular);
648 
649  gammaPolarization1.setX(xTotal);
650  gammaPolarization1.setY(yTotal);
651  gammaPolarization1.setZ(zTotal);
652 
653  return gammaPolarization1;
654 
655 }
656 
657 
658 void G4LowEnergyPolarizedCompton::SystemOfRefChange(G4ThreeVector& direction0,
659  G4ThreeVector& direction1,
660  G4ThreeVector& polarization0,
661  G4ThreeVector& polarization1)
662 {
663  // direction0 is the original photon direction ---> z
664  // polarization0 is the original photon polarization ---> x
665  // need to specify y axis in the real reference frame ---> y
666  G4ThreeVector Axis_Z0 = direction0.unit();
667  G4ThreeVector Axis_X0 = polarization0.unit();
668  G4ThreeVector Axis_Y0 = (Axis_Z0.cross(Axis_X0)).unit(); // to be confirmed;
669 
670  G4double direction_x = direction1.getX();
671  G4double direction_y = direction1.getY();
672  G4double direction_z = direction1.getZ();
673 
674  direction1 = (direction_x*Axis_X0 + direction_y*Axis_Y0 + direction_z*Axis_Z0).unit();
675  G4double polarization_x = polarization1.getX();
676  G4double polarization_y = polarization1.getY();
677  G4double polarization_z = polarization1.getZ();
678 
679  polarization1 = (polarization_x*Axis_X0 + polarization_y*Axis_Y0 + polarization_z*Axis_Z0).unit();
680 
681 }
682 
683 
685 {
686  return ( &particle == G4Gamma::Gamma() );
687 }
688 
689 
691  G4double,
693 {
695  G4double energy = photon->GetKineticEnergy();
696  const G4MaterialCutsCouple* couple = track.GetMaterialCutsCouple();
697  size_t materialIndex = couple->GetIndex();
698  G4double meanFreePath;
699  if (energy > highEnergyLimit) meanFreePath = meanFreePathTable->FindValue(highEnergyLimit,materialIndex);
700  else if (energy < lowEnergyLimit) meanFreePath = DBL_MAX;
701  else meanFreePath = meanFreePathTable->FindValue(energy,materialIndex);
702  return meanFreePath;
703 }
704 
705 
706 
707 
708 
709 
710 
711 
712 
713 
714 
715 
716 
717 
void BuildPhysicsTable(const G4ParticleDefinition &photon)
void SetOccupancyData()
G4int verboseLevel
Definition: G4VProcess.hh:368
tuple a
Definition: test.py:11
G4double GetKineticEnergy() const
CLHEP::Hep3Vector G4ThreeVector
double x() const
G4LowEnergyPolarizedCompton(const G4String &processName="polarLowEnCompt")
const G4DynamicParticle * GetDynamicParticle() const
double dot(const Hep3Vector &) const
float h_Planck
Definition: hepunit.py:263
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4VParticleChange * PostStepDoIt(const G4Track &aTrack, const G4Step &aStep)
G4double BindingEnergy(G4int Z, G4int shellIndex) const
void ProposePolarization(G4double Px, G4double Py, G4double Pz)
double getY() const
tuple x
Definition: test.py:50
bool isOrthogonal(const Hep3Vector &v, double epsilon=tolerance) const
Definition: SpaceVector.cc:237
G4RDVEMDataSet * BuildMeanFreePathForMaterials(const G4DataVector *energyCuts=0)
TFile * file
int G4int
Definition: G4Types.hh:78
void setY(double)
double howOrthogonal(const Hep3Vector &v) const
Definition: SpaceVector.cc:219
double z() const
void setZ(double)
void setX(double)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void LoadData(const G4String &fileName)
Double_t y
Definition: plot.C:279
REAL *8 function var(A, B, C, D)
Definition: dpm25nuc1.f:4649
double getX() const
tuple b
Definition: test.py:12
#define G4UniformRand()
Definition: Randomize.hh:87
G4GLOB_DLL std::ostream G4cout
double energy
Definition: plottest35.C:25
G4double GetMeanFreePath(const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)
Float_t Z
Definition: plot.C:39
bool G4bool
Definition: G4Types.hh:79
const G4ThreeVector & GetMomentumDirection() const
float electron_mass_c2
Definition: hepunit.py:274
Definition: G4Step.hh:76
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
Double_t scale
Definition: plot.C:11
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4double RandomSelectMomentum(G4int Z, G4int shellIndex) const
virtual G4bool Escape(const G4ParticleDefinition *particle, const G4MaterialCutsCouple *couple, G4double energy, G4double safety) const =0
virtual void Initialize(const G4Track &)
void SetNumberOfSecondaries(G4int totSecondaries)
Hep3Vector unit() const
G4StepPoint * GetPostStepPoint() const
double getZ() const
double y() const
G4int SelectRandomAtom(const G4MaterialCutsCouple *couple, G4double e) const
const G4ThreeVector & GetPolarization() const
void ProposeEnergy(G4double finalEnergy)
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:289
G4double GetSafety() const
tuple z
Definition: test.py:28
void AddSecondary(G4Track *aSecondary)
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
void LoadData(const G4String &dataFile)
Hep3Vector cross(const Hep3Vector &) const
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
tuple c
Definition: test.py:13
G4ForceCondition
double mag() const
#define DBL_MAX
Definition: templates.hh:83
G4int SelectRandomShell(G4int Z) const
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
float c_light
Definition: hepunit.py:257
G4bool IsApplicable(const G4ParticleDefinition &definition)
virtual G4bool LoadData(const G4String &fileName)=0