Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Protected Member Functions
G4SeltzerBergerModel Class Reference

#include <G4SeltzerBergerModel.hh>

Inheritance diagram for G4SeltzerBergerModel:
G4eBremsstrahlungRelModel G4VEmModel G4ePolarizedBremsstrahlungModel

Public Member Functions

 G4SeltzerBergerModel (const G4ParticleDefinition *p=0, const G4String &nam="eBremSB")
 
virtual ~G4SeltzerBergerModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
 
virtual void InitialiseForElement (const G4ParticleDefinition *, G4int Z)
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double cutEnergy, G4double maxEnergy)
 
void SetBicubicInterpolationFlag (G4bool)
 
- Public Member Functions inherited from G4eBremsstrahlungRelModel
 G4eBremsstrahlungRelModel (const G4ParticleDefinition *p=0, const G4String &nam="eBremLPM")
 
virtual ~G4eBremsstrahlungRelModel ()
 
virtual void InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel)
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy)
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double tkin, G4double Z, G4double, G4double cutEnergy, G4double maxEnergy=DBL_MAX)
 
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double)
 
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *, G4double cut)
 
void SetLPMconstant (G4double val)
 
G4double LPMconstant () const
 
void SetLowestKinEnergy (G4double)
 
G4double LowestKinEnergy () const
 
- Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
 
virtual ~G4VEmModel ()
 
virtual void InitialiseForMaterial (const G4ParticleDefinition *, const G4Material *)
 
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double ChargeSquareRatio (const G4Track &)
 
virtual G4double GetChargeSquareRatio (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual G4double GetParticleCharge (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void StartTracking (G4Track *)
 
virtual void CorrectionsAlongStep (const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
 
virtual G4double Value (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *)
 
virtual void DefineForRegion (const G4Region *)
 
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
 
std::vector
< G4EmElementSelector * > * 
GetElementSelectors ()
 
void SetElementSelectors (std::vector< G4EmElementSelector * > *)
 
G4double ComputeDEDX (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
G4double CrossSection (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeMeanFreePath (const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, const G4Element *, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectIsotopeNumber (const G4Element *)
 
const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectRandomAtom (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectRandomAtomNumber (const G4Material *)
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=0)
 
void SetCrossSectionTable (G4PhysicsTable *, G4bool isLocal)
 
G4ElementDataGetElementData ()
 
G4PhysicsTableGetCrossSectionTable ()
 
G4VEmFluctuationModelGetModelOfFluctuations ()
 
G4VEmAngularDistributionGetAngularDistribution ()
 
void SetAngularDistribution (G4VEmAngularDistribution *)
 
G4double HighEnergyLimit () const
 
G4double LowEnergyLimit () const
 
G4double HighEnergyActivationLimit () const
 
G4double LowEnergyActivationLimit () const
 
G4double PolarAngleLimit () const
 
G4double SecondaryThreshold () const
 
G4bool LPMFlag () const
 
G4bool DeexcitationFlag () const
 
G4bool ForceBuildTableFlag () const
 
G4bool UseAngularGeneratorFlag () const
 
void SetAngularGeneratorFlag (G4bool)
 
void SetHighEnergyLimit (G4double)
 
void SetLowEnergyLimit (G4double)
 
void SetActivationHighEnergyLimit (G4double)
 
void SetActivationLowEnergyLimit (G4double)
 
G4bool IsActive (G4double kinEnergy)
 
void SetPolarAngleLimit (G4double)
 
void SetSecondaryThreshold (G4double)
 
void SetLPMFlag (G4bool val)
 
void SetDeexcitationFlag (G4bool val)
 
void SetForceBuildTable (G4bool val)
 
void SetMasterThread (G4bool val)
 
G4bool IsMaster () const
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
const G4StringGetName () const
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
const G4ElementGetCurrentElement () const
 

Protected Member Functions

virtual G4double ComputeDXSectionPerAtom (G4double gammaEnergy)
 
virtual G4String DirectoryPath () const
 
- Protected Member Functions inherited from G4eBremsstrahlungRelModel
void SetCurrentElement (const G4double)
 
- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLossGetParticleChangeForLoss ()
 
G4ParticleChangeForGammaGetParticleChangeForGamma ()
 
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 

Additional Inherited Members

- Protected Attributes inherited from G4eBremsstrahlungRelModel
G4NistManagernist
 
const G4ParticleDefinitionparticle
 
G4ParticleDefinitiontheGamma
 
G4ParticleChangeForLossfParticleChange
 
G4double bremFactor
 
G4double particleMass
 
G4double kinEnergy
 
G4double totalEnergy
 
G4double currentZ
 
G4double densityFactor
 
G4double densityCorr
 
G4bool isElectron
 
- Protected Attributes inherited from G4VEmModel
G4ElementDatafElementData
 
G4VParticleChangepParticleChange
 
G4PhysicsTablexSectionTable
 
const std::vector< G4double > * theDensityFactor
 
const std::vector< G4int > * theDensityIdx
 
size_t idxTable
 

Detailed Description

Definition at line 61 of file G4SeltzerBergerModel.hh.

Constructor & Destructor Documentation

G4SeltzerBergerModel::G4SeltzerBergerModel ( const G4ParticleDefinition p = 0,
const G4String nam = "eBremSB" 
)

Definition at line 83 of file G4SeltzerBergerModel.cc.

References python.hepunit::keV, G4eBremsstrahlungRelModel::LowestKinEnergy(), G4VEmModel::SetLowEnergyLimit(), G4eBremsstrahlungRelModel::SetLowestKinEnergy(), and G4VEmModel::SetLPMFlag().

85  : G4eBremsstrahlungRelModel(p,nam),useBicubicInterpolation(false)
86 {
89  SetLPMFlag(false);
90  nwarn = 0;
91  idx = idy = 0;
92 }
G4eBremsstrahlungRelModel(const G4ParticleDefinition *p=0, const G4String &nam="eBremLPM")
void SetLPMFlag(G4bool val)
Definition: G4VEmModel.hh:732
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:690
G4SeltzerBergerModel::~G4SeltzerBergerModel ( )
virtual

Definition at line 96 of file G4SeltzerBergerModel.cc.

References G4VEmModel::IsMaster().

97 {
98  if(IsMaster()) {
99  for(size_t i=0; i<101; ++i) {
100  if(dataSB[i]) {
101  delete dataSB[i];
102  dataSB[i] = 0;
103  }
104  }
105  }
106 }
G4bool IsMaster() const
Definition: G4VEmModel.hh:676

Member Function Documentation

G4double G4SeltzerBergerModel::ComputeDXSectionPerAtom ( G4double  gammaEnergy)
protectedvirtual

Reimplemented from G4eBremsstrahlungRelModel.

Definition at line 193 of file G4SeltzerBergerModel.cc.

References G4eBremsstrahlungRelModel::bremFactor, G4eBremsstrahlungRelModel::currentZ, G4Exp(), G4Log(), G4lrint(), InitialiseForElement(), G4eBremsstrahlungRelModel::isElectron, G4eBremsstrahlungRelModel::kinEnergy, python.hepunit::MeV, python.hepunit::millibarn, G4eBremsstrahlungRelModel::particleMass, G4eBremsstrahlungRelModel::totalEnergy, G4Physics2DVector::Value(), and test::x.

194 {
195 
196  if(gammaEnergy < 0.0 || kinEnergy <= 0.0) { return 0.0; }
197  G4double x = gammaEnergy/kinEnergy;
199  G4int Z = G4lrint(currentZ);
200 
201  //G4cout << "G4SeltzerBergerModel::ComputeDXSectionPerAtom Z= " << Z
202  // << " x= " << x << " y= " << y << " " << dataSB[Z] << G4endl;
203  if(!dataSB[Z]) { InitialiseForElement(0, Z); }
204  /*
205  G4ExceptionDescription ed;
206  ed << "Bremsstrahlung data for Z= " << Z
207  << " are not initialized!";
208  G4Exception("G4SeltzerBergerModel::ComputeDXSectionPerAtom()","em0005",
209  FatalException, ed,
210  "G4LEDATA version should be G4EMLOW6.23 or later.");
211  }
212  */
213  G4double invb2 =
215  G4double cross = dataSB[Z]->Value(x,y,idx,idy)*invb2*millibarn/bremFactor;
216 
217  if(!isElectron) {
218  G4double invbeta1 = sqrt(invb2);
219  G4double e2 = kinEnergy - gammaEnergy;
220  if(e2 > 0.0) {
221  G4double invbeta2 = (e2 + particleMass)/sqrt(e2*(e2 + 2*particleMass));
222  G4double xxx = alpha*currentZ*(invbeta1 - invbeta2);
223  if(xxx < expnumlim) { cross = 0.0; }
224  else { cross *= G4Exp(xxx); }
225  } else {
226  cross = 0.0;
227  }
228  }
229 
230  return cross;
231 }
int G4int
Definition: G4Types.hh:78
int millibarn
Definition: hepunit.py:40
G4double Value(G4double x, G4double y, size_t &lastidx, size_t &lastidy) const
G4double G4Log(G4double x)
Definition: G4Log.hh:227
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:180
virtual void InitialiseForElement(const G4ParticleDefinition *, G4int Z)
int G4lrint(double ad)
Definition: templates.hh:163
double G4double
Definition: G4Types.hh:76
G4String G4SeltzerBergerModel::DirectoryPath ( ) const
protectedvirtual

Definition at line 139 of file G4SeltzerBergerModel.cc.

140 {
141  return "/brem_SB/br";
142 }
void G4SeltzerBergerModel::Initialise ( const G4ParticleDefinition p,
const G4DataVector cuts 
)
virtual

Reimplemented from G4eBremsstrahlungRelModel.

Reimplemented in G4ePolarizedBremsstrahlungModel.

Definition at line 110 of file G4SeltzerBergerModel.cc.

References G4Element::GetElementTable(), G4Element::GetNumberOfElements(), G4eBremsstrahlungRelModel::Initialise(), and G4VEmModel::IsMaster().

Referenced by G4ePolarizedBremsstrahlungModel::Initialise().

112 {
113  // Access to elements
114  if(IsMaster()) {
115 
116  // check environment variable
117  // Build the complete string identifying the file with the data set
118  char* path = getenv("G4LEDATA");
119 
120  const G4ElementTable* theElmTable = G4Element::GetElementTable();
121  size_t numOfElm = G4Element::GetNumberOfElements();
122  if(numOfElm > 0) {
123  for(size_t i=0; i<numOfElm; ++i) {
124  G4int Z = G4int(((*theElmTable)[i])->GetZ());
125  if(Z < 1) { Z = 1; }
126  else if(Z > 100) { Z = 100; }
127  //G4cout << "Z= " << Z << G4endl;
128  // Initialisation
129  if(!dataSB[Z]) { ReadData(Z, path); }
130  }
131  }
132  }
133 
135 }
G4bool IsMaster() const
Definition: G4VEmModel.hh:676
int G4int
Definition: G4Types.hh:78
static size_t GetNumberOfElements()
Definition: G4Element.cc:402
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
std::vector< G4Element * > G4ElementTable
static G4ElementTable * GetElementTable()
Definition: G4Element.cc:395
void G4SeltzerBergerModel::InitialiseForElement ( const G4ParticleDefinition ,
G4int  Z 
)
virtual

Reimplemented from G4VEmModel.

Definition at line 370 of file G4SeltzerBergerModel.cc.

Referenced by ComputeDXSectionPerAtom().

372 {
373  G4AutoLock l(&SeltzerBergerModelMutex);
374  // G4cout << "G4SeltzerBergerModel::InitialiseForElement Z= " << Z << G4endl;
375  if(!dataSB[Z]) { ReadData(Z); }
376 }
void G4SeltzerBergerModel::SampleSecondaries ( std::vector< G4DynamicParticle * > *  vdp,
const G4MaterialCutsCouple couple,
const G4DynamicParticle dp,
G4double  cutEnergy,
G4double  maxEnergy 
)
virtual

Reimplemented from G4eBremsstrahlungRelModel.

Reimplemented in G4ePolarizedBremsstrahlungModel.

Definition at line 236 of file G4SeltzerBergerModel.cc.

References G4eBremsstrahlungRelModel::currentZ, G4eBremsstrahlungRelModel::densityCorr, G4eBremsstrahlungRelModel::densityFactor, python.hepunit::electron_mass_c2, python.hepunit::fine_structure_const, G4eBremsstrahlungRelModel::fParticleChange, fStopAndKill, G4Exception(), G4Exp(), G4Log(), G4UniformRand, G4VEmModel::GetAngularDistribution(), G4DynamicParticle::GetKineticEnergy(), G4MaterialCutsCouple::GetMaterial(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetParticleName(), G4eBremsstrahlungRelModel::isElectron, JustWarning, python.hepunit::MeV, G4INCL::Math::min(), G4eBremsstrahlungRelModel::particle, G4eBremsstrahlungRelModel::particleMass, G4VParticleChange::ProposeTrackStatus(), G4VEmAngularDistribution::SampleDirection(), G4VEmModel::SecondaryThreshold(), G4VEmModel::SelectRandomAtom(), G4eBremsstrahlungRelModel::SetCurrentElement(), G4ParticleChangeForLoss::SetProposedKineticEnergy(), G4ParticleChangeForLoss::SetProposedMomentumDirection(), G4eBremsstrahlungRelModel::SetupForMaterial(), G4eBremsstrahlungRelModel::theGamma, G4eBremsstrahlungRelModel::totalEnergy, python.hepunit::twopi, test::v, G4Physics2DVector::Value(), G4VEmModel::Value(), and test::x.

Referenced by G4ePolarizedBremsstrahlungModel::SampleSecondaries().

241 {
242  G4double kineticEnergy = dp->GetKineticEnergy();
243  G4double cut = std::min(cutEnergy, kineticEnergy);
244  G4double emax = std::min(maxEnergy, kineticEnergy);
245  if(cut >= emax) { return; }
246 
247  SetupForMaterial(particle, couple->GetMaterial(), kineticEnergy);
248 
249  const G4Element* elm =
250  SelectRandomAtom(couple,particle,kineticEnergy,cut,emax);
251  SetCurrentElement(elm->GetZ());
252  G4int Z = G4int(currentZ);
253 
254  totalEnergy = kineticEnergy + particleMass;
256  G4double totMomentum = sqrt(kineticEnergy*(totalEnergy + electron_mass_c2));
257  /*
258  G4cout << "G4SeltzerBergerModel::SampleSecondaries E(MeV)= "
259  << kineticEnergy/MeV
260  << " Z= " << Z << " cut(MeV)= " << cut/MeV
261  << " emax(MeV)= " << emax/MeV << " corr= " << densityCorr << G4endl;
262  */
263  G4double xmin = G4Log(cut*cut + densityCorr);
264  G4double xmax = G4Log(emax*emax + densityCorr);
265  G4double y = G4Log(kineticEnergy/MeV);
266 
267  G4double gammaEnergy, v;
268 
269  // majoranta
270  G4double x0 = cut/kineticEnergy;
271  G4double vmax = dataSB[Z]->Value(x0, y, idx, idy)*1.02;
272  // G4double invbeta1 = 0;
273 
274  // majoranta corrected for e-
275  if(isElectron && x0 < 0.97 &&
276  ((kineticEnergy > epeaklimit) || (kineticEnergy < elowlimit))) {
277  G4double ylim = std::min(ylimit[Z],1.1*dataSB[Z]->Value(0.97,y,idx,idy));
278  if(ylim > vmax) { vmax = ylim; }
279  }
280  if(x0 < 0.05) { vmax *= 1.2; }
281 
282  //G4cout<<"y= "<<y<<" xmin= "<<xmin<<" xmax= "<<xmax
283  //<<" vmax= "<<vmax<<G4endl;
284  // G4int ncount = 0;
285  do {
286  //++ncount;
287  G4double x = G4Exp(xmin + G4UniformRand()*(xmax - xmin)) - densityCorr;
288  if(x < 0.0) { x = 0.0; }
289  gammaEnergy = sqrt(x);
290  G4double x1 = gammaEnergy/kineticEnergy;
291  v = dataSB[Z]->Value(x1, y, idx, idy);
292 
293  // correction for positrons
294  if(!isElectron) {
295  G4double e1 = kineticEnergy - cut;
296  G4double invbeta1 = (e1 + particleMass)/sqrt(e1*(e1 + 2*particleMass));
297  G4double e2 = kineticEnergy - gammaEnergy;
298  G4double invbeta2 = (e2 + particleMass)/sqrt(e2*(e2 + 2*particleMass));
299  G4double xxx = twopi*fine_structure_const*currentZ*(invbeta1 - invbeta2);
300 
301  if(xxx < expnumlim) { v = 0.0; }
302  else { v *= G4Exp(xxx); }
303  }
304 
305  if (v > 1.05*vmax && nwarn < 5) {
306  ++nwarn;
308  ed << "### G4SeltzerBergerModel Warning: Majoranta exceeded! "
309  << v << " > " << vmax << " by " << v/vmax
310  << " Egamma(MeV)= " << gammaEnergy
311  << " Ee(MeV)= " << kineticEnergy
312  << " Z= " << Z << " " << particle->GetParticleName();
313 
314  if ( 20 == nwarn ) {
315  ed << "\n ### G4SeltzerBergerModel Warnings stopped";
316  }
317  G4Exception("G4SeltzerBergerModel::SampleScattering","em0044",
318  JustWarning, ed,"");
319 
320  }
321  } while (v < vmax*G4UniformRand());
322 
323  //
324  // angles of the emitted gamma. ( Z - axis along the parent particle)
325  // use general interface
326  //
327 
328  G4ThreeVector gammaDirection =
329  GetAngularDistribution()->SampleDirection(dp, totalEnergy-gammaEnergy,
330  Z, couple->GetMaterial());
331 
332  // create G4DynamicParticle object for the Gamma
333  G4DynamicParticle* gamma =
334  new G4DynamicParticle(theGamma,gammaDirection,gammaEnergy);
335  vdp->push_back(gamma);
336 
337  G4ThreeVector direction = (totMomentum*dp->GetMomentumDirection()
338  - gammaEnergy*gammaDirection).unit();
339 
340  /*
341  G4cout << "### G4SBModel: v= "
342  << " Eg(MeV)= " << gammaEnergy
343  << " Ee(MeV)= " << kineticEnergy
344  << " DirE " << direction << " DirG " << gammaDirection
345  << G4endl;
346  */
347  // energy of primary
348  G4double finalE = kineticEnergy - gammaEnergy;
349 
350  // stop tracking and create new secondary instead of primary
351  if(gammaEnergy > SecondaryThreshold()) {
354  G4DynamicParticle* el =
355  new G4DynamicParticle(const_cast<G4ParticleDefinition*>(particle),
356  direction, finalE);
357  vdp->push_back(el);
358 
359  // continue tracking
360  } else {
363  }
364 }
G4double SecondaryThreshold() const
Definition: G4VEmModel.hh:627
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4double GetKineticEnergy() const
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:578
int G4int
Definition: G4Types.hh:78
virtual void SetupForMaterial(const G4ParticleDefinition *, const G4Material *, G4double)
const G4String & GetParticleName() const
const G4ParticleDefinition * particle
#define G4UniformRand()
Definition: Randomize.hh:87
G4double Value(G4double x, G4double y, size_t &lastidx, size_t &lastidy) const
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
const G4ThreeVector & GetMomentumDirection() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
float electron_mass_c2
Definition: hepunit.py:274
void SetProposedMomentumDirection(const G4ThreeVector &dir)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4double G4Log(G4double x)
Definition: G4Log.hh:227
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:180
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4ParticleChangeForLoss * fParticleChange
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:510
virtual G4double Value(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
Definition: G4VEmModel.cc:346
const G4Material * GetMaterial() const
void G4SeltzerBergerModel::SetBicubicInterpolationFlag ( G4bool  val)
inline

Definition at line 106 of file G4SeltzerBergerModel.hh.

107 {
108  useBicubicInterpolation = val;
109 }

The documentation for this class was generated from the following files: