G4INCL::ParticleTable Class Reference

#include <G4INCLParticleTable.hh>


Public Types

typedef G4double(*) NuclearMassFn (const G4int, const G4int)
typedef G4double(*) ParticleMassFn (const ParticleType)
typedef G4double(*) SeparationEnergyFn (const ParticleType, const G4int, const G4int)
 StableCluster
 NeutronDecay
 ProtonDecay
 AlphaDecay
 TwoProtonDecay
 TwoNeutronDecay
 ProtonUnbound
 NeutronUnbound
enum  ClusterDecayType {
  StableCluster, NeutronDecay, ProtonDecay, AlphaDecay,
  TwoProtonDecay, TwoNeutronDecay, ProtonUnbound, NeutronUnbound
}

Static Public Member Functions

static void initialize (Config const *const theConfig=0)
 Initialize the particle table.
static G4int getIsospin (const ParticleType t)
 Get the isospin of a particle.
static std::string getName (const ParticleType t)
 Get the native INCL name of the particle.
static std::string getShortName (const ParticleType t)
 Get the short INCL name of the particle.
static std::string getName (const ParticleSpecies s)
 Get the native INCL name of the particle.
static std::string getShortName (const ParticleSpecies s)
 Get the short INCL name of the particle.
static std::string getName (const G4int A, const G4int Z)
 Get the native INCL name of the ion.
static std::string getShortName (const G4int A, const G4int Z)
 Get the short INCL name of the ion.
static G4double getINCLMass (const G4int A, const G4int Z)
 Get INCL nuclear mass (in MeV/c^2).
static G4double getINCLMass (const ParticleType t)
 Get INCL particle mass (in MeV/c^2).
static G4double getRealMass (const G4INCL::ParticleType t)
 Get particle mass (in MeV/c^2).
static G4double getRealMass (const G4int A, const G4int Z)
 Get nuclear mass (in MeV/c^2).
static G4double getTableQValue (const G4int A1, const G4int Z1, const G4int A2, const G4int Z2)
 Get Q-value (in MeV/c^2).
static G4double getTableQValue (const G4int A1, const G4int Z1, const G4int A2, const G4int Z2, const G4int A3, const G4int Z3)
 Get Q-value (in MeV/c^2).
static G4double getTableSpeciesMass (const ParticleSpecies &p)
static G4int getMassNumber (const ParticleType t)
 Get mass number from particle type.
static G4int getChargeNumber (const ParticleType t)
 Get charge number from particle type.
static G4double getNuclearRadius (const G4int A, const G4int Z)
static G4double getRadiusParameter (const G4int A, const G4int Z)
static G4double getMaximumNuclearRadius (const G4int A, const G4int Z)
static G4double getSurfaceDiffuseness (const G4int A, const G4int Z)
static G4double getMomentumRMS (const G4int A, const G4int Z)
 Return the RMS of the momentum distribution (light clusters).
static G4double getSeparationEnergyINCL (const ParticleType t, const G4int, const G4int)
 Return INCL's default separation energy.
static G4double getSeparationEnergyReal (const ParticleType t, const G4int A, const G4int Z)
 Return the real separation energy.
static G4double getSeparationEnergyRealForLight (const ParticleType t, const G4int A, const G4int Z)
 Return the real separation energy only for light nuclei.
static G4double getProtonSeparationEnergy ()
 Getter for protonSeparationEnergy.
static G4double getNeutronSeparationEnergy ()
 Getter for neutronSeparationEnergy.
static void setProtonSeparationEnergy (const G4double s)
 Setter for protonSeparationEnergy.
static void setNeutronSeparationEnergy (const G4double s)
 Setter for protonSeparationEnergy.
static std::string getElementName (const G4int Z)
 Get the name of the element from the atomic number.
static std::string getIUPACElementName (const G4int Z)
 Get the name of an unnamed element from the IUPAC convention.
static G4int parseIUPACElement (std::string const &pS)
 Parse a IUPAC element name.
static IsotopicDistribution
const & 
getNaturalIsotopicDistribution (const G4int Z)
static G4int drawRandomNaturalIsotope (const G4int Z)

Static Public Attributes

static NuclearMassFn getTableMass
 Static pointer to the mass function for nuclei.
static ParticleMassFn getTableParticleMass
 Static pointer to the mass function for particles.
static SeparationEnergyFn getSeparationEnergy
 Static pointer to the separation-energy function.
static const G4int elementTableSize = 113
static const G4double effectiveNucleonMass = 938.2796
static const G4double effectiveNucleonMass2 = 8.8036860777616e5
static const G4double effectiveDeltaMass = 1232.0
static const G4double effectivePionMass = 138.0
static const G4double effectiveDeltaDecayThreshold
static const G4int maxClusterMass = 12
static const G4int maxClusterCharge = 8
static const G4int clusterTableZSize = ParticleTable::maxClusterCharge+1
static const G4int clusterTableASize = ParticleTable::maxClusterMass+1
static const G4double clusterPosFact [maxClusterMass+1]
static const G4double clusterPosFact2 [maxClusterMass+1]
static const G4int clusterZMin [maxClusterMass+1] = {0, 0, 1, 1, 1, 1, 1, 1, 2, 2, 2, 3, 3}
static const G4int clusterZMax [maxClusterMass+1] = {0, 0, 1, 2, 3, 3, 5, 5, 6, 6, 7, 7, 8}
static const G4double clusterPhaseSpaceCut [maxClusterMass+1]
static G4IonTabletheG4IonTable
static const ClusterDecayType clusterDecayMode [clusterTableZSize][clusterTableASize]
static const G4double eSquared = 1.439964
 Coulomb conversion factor, in MeV*fm.

Protected Member Functions

 ParticleTable ()
 ~ParticleTable ()


Detailed Description

Definition at line 57 of file G4INCLParticleTable.hh.


Member Typedef Documentation

typedef G4double(*) G4INCL::ParticleTable::NuclearMassFn(const G4int, const G4int)

Definition at line 148 of file G4INCLParticleTable.hh.

typedef G4double(*) G4INCL::ParticleTable::ParticleMassFn(const ParticleType)

Definition at line 149 of file G4INCLParticleTable.hh.

typedef G4double(*) G4INCL::ParticleTable::SeparationEnergyFn(const ParticleType, const G4int, const G4int)

Definition at line 161 of file G4INCLParticleTable.hh.


Member Enumeration Documentation

enum G4INCL::ParticleTable::ClusterDecayType

Enumerator:
StableCluster 
NeutronDecay 
ProtonDecay 
AlphaDecay 
TwoProtonDecay 
TwoNeutronDecay 
ProtonUnbound 
NeutronUnbound 

Definition at line 313 of file G4INCLParticleTable.hh.

00313                           {
00314       StableCluster,
00315       NeutronDecay,
00316       ProtonDecay,
00317       AlphaDecay,
00318       TwoProtonDecay,
00319       TwoNeutronDecay,
00320       ProtonUnbound,
00321       NeutronUnbound
00322     };


Constructor & Destructor Documentation

G4INCL::ParticleTable::ParticleTable (  )  [inline, protected]

Definition at line 340 of file G4INCLParticleTable.hh.

00340 {};

G4INCL::ParticleTable::~ParticleTable (  )  [inline, protected]

Definition at line 341 of file G4INCLParticleTable.hh.

00341 {};


Member Function Documentation

static G4int G4INCL::ParticleTable::drawRandomNaturalIsotope ( const G4int  Z  )  [inline, static]

Definition at line 335 of file G4INCLParticleTable.hh.

Referenced by G4INCL::INCL::prepareReaction().

00335                                                          {
00336       return getNaturalIsotopicDistributions()->drawRandomIsotope(Z);
00337     }

static G4int G4INCL::ParticleTable::getChargeNumber ( const ParticleType  t  )  [inline, static]

Get charge number from particle type.

Definition at line 187 of file G4INCLParticleTable.hh.

References G4INCL::DeltaMinus, G4INCL::DeltaPlus, G4INCL::DeltaPlusPlus, G4INCL::DeltaZero, G4INCL::Neutron, G4INCL::PiMinus, G4INCL::PiPlus, G4INCL::PiZero, and G4INCL::Proton.

00187                                                        {
00188       switch(t) {
00189         case DeltaPlusPlus:
00190           return 2;
00191           break;
00192         case Proton:
00193         case DeltaPlus:
00194         case PiPlus:
00195           return 1;
00196           break;
00197         case Neutron:
00198         case DeltaZero:
00199         case PiZero:
00200           return 0;
00201           break;
00202         case DeltaMinus:
00203         case PiMinus:
00204           return -1;
00205           break;
00206         default:
00207           return 0;
00208           break;
00209       }
00210     }

std::string G4INCL::ParticleTable::getElementName ( const G4int  Z  )  [static]

Get the name of the element from the atomic number.

Definition at line 624 of file G4INCLParticleTable.cc.

References elementTableSize, getIUPACElementName(), and WARN.

Referenced by getName(), and getShortName().

00624                                                        {
00625     if(Z<1) {
00626       WARN("getElementName called with Z<1" << std::endl);
00627       return elementTable[0];
00628     } else if(Z<elementTableSize)
00629       return elementTable[Z];
00630     else
00631       return getIUPACElementName(Z);
00632   }

G4double G4INCL::ParticleTable::getINCLMass ( const ParticleType  t  )  [static]

Get INCL particle mass (in MeV/c^2).

Definition at line 451 of file G4INCLParticleTable.cc.

References ERROR, G4INCL::Neutron, G4INCL::PiMinus, G4INCL::PiPlus, G4INCL::PiZero, and G4INCL::Proton.

00451                                                            {
00452     if(pt == Proton) {
00453       return protonMass;
00454     } else if(pt == Neutron) {
00455       return neutronMass;
00456     } else if(pt == PiPlus) {
00457       return piPlusMass;
00458     } else if(pt == PiMinus) {
00459       return piMinusMass;
00460     } else if(pt == PiZero) {
00461       return piZeroMass;
00462     } else {
00463       ERROR("ParticleTable::getMass : Unknown particle type." << std::endl);
00464       return 0.0;
00465     }
00466   }

G4double G4INCL::ParticleTable::getINCLMass ( const G4int  A,
const G4int  Z 
) [static]

Get INCL nuclear mass (in MeV/c^2).

Definition at line 518 of file G4INCLParticleTable.cc.

References G4INCL::Neutron, G4INCL::PiMinus, G4INCL::PiPlus, and G4INCL::Proton.

Referenced by G4INCL::Particle::getEmissionQValueCorrection(), G4INCL::DeltaDecayChannel::getFinalState(), G4INCL::Particle::getINCLMass(), G4INCL::Particle::getTransferQValueCorrection(), initialize(), and G4INCL::CDPP::isBlocked().

00518                                                                   {
00519 // assert(A>=0);
00520     // For nuclei with Z<0 or Z>A, assume that the exotic charge state is due to pions
00521     if(Z<0)
00522       return A*neutronMass - Z*getINCLMass(PiMinus);
00523     else if(Z>A)
00524       return A*protonMass + (A-Z)*getINCLMass(PiPlus);
00525     else if(A>1)
00526       return Z*(protonMass - protonSeparationEnergy) + (A-Z)*(neutronMass - neutronSeparationEnergy);
00527     else if(A==1 && Z==0)
00528       return getINCLMass(Neutron);
00529     else if(A==1 && Z==1)
00530       return getINCLMass(Proton);
00531     else
00532       return 0.;
00533   }

G4int G4INCL::ParticleTable::getIsospin ( const ParticleType  t  )  [static]

Get the isospin of a particle.

Definition at line 347 of file G4INCLParticleTable.cc.

References G4INCL::DeltaMinus, G4INCL::DeltaPlus, G4INCL::DeltaPlusPlus, G4INCL::DeltaZero, ERROR, G4INCL::Neutron, G4INCL::PiMinus, G4INCL::PiPlus, G4INCL::PiZero, and G4INCL::Proton.

Referenced by G4INCL::CrossSections::deltaProduction(), G4INCL::ElasticChannel::getFinalState(), G4INCL::DeltaProductionChannel::getFinalState(), G4INCL::Nucleus::insertParticle(), G4INCL::CrossSections::pionNucleon(), and G4INCL::CrossSections::recombination().

00347                                                       {
00348     // Actually this is the 3rd component of isospin (I_z) multiplied by 2!
00349     if(t == Proton) {
00350       return 1;
00351     } else if(t == Neutron) {
00352       return -1;
00353     } else if(t == PiPlus) {
00354       return 2;
00355     } else if(t == PiMinus) {
00356       return -2;
00357     } else if(t == PiZero) {
00358       return 0;
00359     } else if(t == DeltaPlusPlus) {
00360       return 3;
00361     } else if(t == DeltaPlus) {
00362       return 1;
00363     } else if(t == DeltaZero) {
00364       return -1;
00365     } else if(t == DeltaMinus) {
00366       return -3;
00367     }
00368 
00369     ERROR("Requested isospin of an unknown particle!");
00370     return -10; // Unknown
00371   }

std::string G4INCL::ParticleTable::getIUPACElementName ( const G4int  Z  )  [static]

Get the name of an unnamed element from the IUPAC convention.

Definition at line 677 of file G4INCLParticleTable.cc.

Referenced by getElementName().

00677                                                             {
00678     std::stringstream elementStream;
00679     elementStream << Z;
00680     std::string elementName = elementStream.str();
00681     std::transform(elementName.begin(), elementName.end(), elementName.begin(), ParticleTable::intToIUPAC);
00682     elementName[0] = std::toupper(elementName.at(0));
00683     return elementName;
00684   }

static G4int G4INCL::ParticleTable::getMassNumber ( const ParticleType  t  )  [inline, static]

Get mass number from particle type.

Definition at line 165 of file G4INCLParticleTable.hh.

References G4INCL::DeltaMinus, G4INCL::DeltaPlus, G4INCL::DeltaPlusPlus, G4INCL::DeltaZero, G4INCL::Neutron, G4INCL::PiMinus, G4INCL::PiPlus, G4INCL::PiZero, and G4INCL::Proton.

00165                                                      {
00166       switch(t) {
00167         case Proton:
00168         case Neutron:
00169         case DeltaPlusPlus:
00170         case DeltaPlus:
00171         case DeltaZero:
00172         case DeltaMinus:
00173           return 1;
00174           break;
00175         case PiPlus:
00176         case PiMinus:
00177         case PiZero:
00178           return 0;
00179           break;
00180         default:
00181           return 0;
00182           break;
00183       }
00184     }

G4double G4INCL::ParticleTable::getMaximumNuclearRadius ( const G4int  A,
const G4int  Z 
) [static]

Definition at line 589 of file G4INCLParticleTable.cc.

References ERROR, getNuclearRadius(), and getSurfaceDiffuseness().

Referenced by G4INCL::NuclearDensityFactory::createRCDFTable(), and G4INCL::NuclearDensityFactory::createRPCorrelationTable().

00589                                                                               {
00590     const G4double XFOISA = 8.0;
00591     if(A >= 19) {
00592       return getNuclearRadius(A,Z) + XFOISA * getSurfaceDiffuseness(A,Z);
00593     } else if(A < 19 && A >= 6) {
00594       return 5.5 + 0.3 * (G4double(A) - 6.0)/12.0;
00595     } else if(A >= 2) {
00596       return ParticleTable::getNuclearRadius(A, Z) + 4.5;
00597     } else {
00598       ERROR("ParticleTable::getMaximumNuclearRadius : No maximum radius for nucleus A = " << A << " Z = " << Z << std::endl);
00599       return 0.0;
00600     }
00601   }

static G4double G4INCL::ParticleTable::getMomentumRMS ( const G4int  A,
const G4int  Z 
) [inline, static]

Return the RMS of the momentum distribution (light clusters).

Definition at line 218 of file G4INCLParticleTable.hh.

References G4INCL::PhysicalConstants::Pf, and G4INCL::Math::sqrtThreeFifths.

Referenced by G4INCL::NuclearDensityFactory::createPCDFTable().

00218                                                                  {
00219 // assert(Z>=0 && A>=0 && Z<=A);
00220       if(Z<clusterTableZSize && A<clusterTableASize)
00221         return momentumRMS[Z][A];
00222       else
00223         return Math::sqrtThreeFifths * PhysicalConstants::Pf;
00224     }

std::string G4INCL::ParticleTable::getName ( const G4int  A,
const G4int  Z 
) [static]

Get the native INCL name of the ion.

Definition at line 387 of file G4INCLParticleTable.cc.

References getElementName().

00387                                                                {
00388     std::stringstream stream;
00389     stream << getElementName(Z) << "-" << A;
00390     return stream.str();
00391   }

std::string G4INCL::ParticleTable::getName ( const ParticleSpecies  s  )  [static]

Get the native INCL name of the particle.

Definition at line 380 of file G4INCLParticleTable.cc.

References G4INCL::Composite, getName(), G4INCL::ParticleSpecies::theA, G4INCL::ParticleSpecies::theType, and G4INCL::ParticleSpecies::theZ.

00380                                                           {
00381     if(s.theType==Composite)
00382       return getName(s.theA,s.theZ);
00383     else
00384       return getName(s.theType);
00385   }

std::string G4INCL::ParticleTable::getName ( const ParticleType  t  )  [static]

Get the native INCL name of the particle.

Definition at line 401 of file G4INCLParticleTable.cc.

References G4INCL::Composite, G4INCL::DeltaMinus, G4INCL::DeltaPlus, G4INCL::DeltaPlusPlus, G4INCL::DeltaZero, G4INCL::Neutron, G4INCL::PiMinus, G4INCL::PiPlus, G4INCL::PiZero, and G4INCL::Proton.

Referenced by G4INCL::Particle::dump(), getName(), G4INCL::Particle::print(), G4INCL::Cluster::print(), and G4INCL::Config::summary().

00401                                                        {
00402     if(p == G4INCL::Proton) {
00403       return std::string("proton");
00404     } else if(p == G4INCL::Neutron) {
00405       return std::string("neutron");
00406     } else if(p == G4INCL::DeltaPlusPlus) {
00407       return std::string("delta++");
00408     } else if(p == G4INCL::DeltaPlus) {
00409       return std::string("delta+");
00410     } else if(p == G4INCL::DeltaZero) {
00411       return std::string("delta0");
00412     } else if(p == G4INCL::DeltaMinus) {
00413       return std::string("delta-");
00414     } else if(p == G4INCL::PiPlus) {
00415       return std::string("pi+");
00416     } else if(p == G4INCL::PiZero) {
00417       return std::string("pi0");
00418     } else if(p == G4INCL::PiMinus) {
00419       return std::string("pi-");
00420     } else if(p == G4INCL::Composite) {
00421       return std::string("composite");
00422     }
00423     return std::string("unknown");
00424   }

static IsotopicDistribution const& G4INCL::ParticleTable::getNaturalIsotopicDistribution ( const G4int  Z  )  [inline, static]

Definition at line 331 of file G4INCLParticleTable.hh.

00331                                                                                      {
00332       return getNaturalIsotopicDistributions()->getIsotopicDistribution(Z);
00333     }

static G4double G4INCL::ParticleTable::getNeutronSeparationEnergy (  )  [inline, static]

Getter for neutronSeparationEnergy.

Definition at line 264 of file G4INCLParticleTable.hh.

00264 { return neutronSeparationEnergy; }

G4double G4INCL::ParticleTable::getNuclearRadius ( const G4int  A,
const G4int  Z 
) [static]

Definition at line 535 of file G4INCLParticleTable.cc.

References clusterTableASize, ERROR, getRadiusParameter(), and getSurfaceDiffuseness().

Referenced by getMaximumNuclearRadius(), G4INCL::NuclearDensity::getTransmissionRadius(), G4INCL::CoulombNonRelativistic::maxImpactParameter(), G4INCL::CoulombNone::maxImpactParameter(), and G4INCL::StandardPropagationModel::shootComposite().

00535                                                                        {
00536 // assert(A>0 && Z>=0 && Z<=A);
00537     if(A >= 19 || (A < 6 && A >= 2)) {
00538       // For large (Woods-Saxon or Modified Harmonic Oscillator) or small
00539       // (Gaussian) nuclei, the radius parameter is just the nuclear radius
00540       return getRadiusParameter(A,Z);
00541     } else if(A < clusterTableASize && Z < clusterTableZSize && A >= 6) {
00542       const G4double thisRMS = positionRMS[Z][A];
00543       if(thisRMS>0.0)
00544         return thisRMS;
00545       else {
00546         ERROR("ParticleTable::getRadiusParameter : Radius for nucleus A = " << A << " Z = " << Z << " is ==0.0" << std::endl);
00547         return 0.0;
00548       }
00549     } else if(A < 19) {
00550       const G4double theRadiusParameter = getRadiusParameter(A, Z);
00551       const G4double theDiffusenessParameter = getSurfaceDiffuseness(A, Z);
00552       // The formula yields the nuclear RMS radius based on the parameters of
00553       // the nuclear-density function
00554       return 1.581*theDiffusenessParameter*
00555         (2.+5.*theRadiusParameter)/(2.+3.*theRadiusParameter);
00556     } else {
00557       ERROR("ParticleTable::getNuclearRadius: No radius for nucleus A = " << A << " Z = " << Z << std::endl);
00558       return 0.0;
00559     }
00560   }

static G4double G4INCL::ParticleTable::getProtonSeparationEnergy (  )  [inline, static]

Getter for protonSeparationEnergy.

Definition at line 261 of file G4INCLParticleTable.hh.

00261 { return protonSeparationEnergy; }

G4double G4INCL::ParticleTable::getRadiusParameter ( const G4int  A,
const G4int  Z 
) [static]

Definition at line 562 of file G4INCLParticleTable.cc.

References clusterTableZSize, and ERROR.

Referenced by G4INCL::NuclearDensityFactory::createRCDFTable(), G4INCL::NuclearDensityFactory::createRPCorrelationTable(), and getNuclearRadius().

00562                                                                          {
00563 // assert(A>0 && Z>=0 && Z<=A);
00564     if(A >= 28) {
00565       // phenomenological radius fit
00566       return (2.745e-4 * A + 1.063) * std::pow(A, 1.0/3.0);
00567     } else if(A < 6 && A >= 2) {
00568       if(Z<clusterTableZSize) {
00569         const G4double thisRMS = positionRMS[Z][A];
00570         if(thisRMS>0.0)
00571           return thisRMS;
00572         else {
00573           ERROR("ParticleTable::getRadiusParameter : Radius for nucleus A = " << A << " Z = " << Z << " is ==0.0" << std::endl);
00574           return 0.0;
00575         }
00576       } else {
00577         ERROR("ParticleTable::getRadiusParameter : No radius for nucleus A = " << A << " Z = " << Z << std::endl);
00578         return 0.0;
00579       }
00580     } else if(A < 28 && A >= 6) {
00581       return mediumRadius[A-1];
00582       //      return 1.581*mediumDiffuseness[A-1]*(2.+5.*mediumRadius[A-1])/(2.+3.*mediumRadius[A-1]);
00583     } else {
00584       ERROR("ParticleTable::getRadiusParameter: No radius for nucleus A = " << A << " Z = " << Z << std::endl);
00585       return 0.0;
00586     }
00587   }

G4double G4INCL::ParticleTable::getRealMass ( const G4int  A,
const G4int  Z 
) [static]

Get nuclear mass (in MeV/c^2).

Definition at line 490 of file G4INCLParticleTable.cc.

References DEBUG, G4IonTable::GetNucleusMass(), getRealMass(), G4INCL::Neutron, G4INCL::PiMinus, G4INCL::PiPlus, G4INCL::Proton, and theG4IonTable.

00490                                                                   {
00491 // assert(A>=0);
00492     // For nuclei with Z<0 or Z>A, assume that the exotic charge state is due to pions
00493     if(Z<0)
00494       return A*neutronMass - Z*getRealMass(PiMinus);
00495     else if(Z>A)
00496       return A*protonMass + (A-Z)*getRealMass(PiPlus);
00497     else if(Z==0)
00498       return A*getRealMass(Neutron);
00499     else if(A==Z)
00500       return A*getRealMass(Proton);
00501     else if(A>1) {
00502 #ifndef INCLXX_IN_GEANT4_MODE
00503       if(ParticleTable::hasMassTable(A,Z))
00504         return ParticleTable::massTable.at(Z).at(A);
00505       else {
00506         DEBUG("Real mass unavailable for isotope A=" << A << ", Z=" << Z
00507             << ", using Weizsaecker's formula"
00508             << std::endl);
00509         return getWeizsaeckerMass(A,Z);
00510       }
00511 #else
00512       return theG4IonTable->GetNucleusMass(Z,A) / MeV;
00513 #endif
00514     } else
00515       return 0.;
00516   }

G4double G4INCL::ParticleTable::getRealMass ( const G4INCL::ParticleType  t  )  [static]

Get particle mass (in MeV/c^2).

Definition at line 468 of file G4INCLParticleTable.cc.

References ERROR, G4INCL::Neutron, G4INCL::PiMinus, G4INCL::PiPlus, G4INCL::PiZero, and G4INCL::Proton.

Referenced by getRealMass(), G4INCL::Particle::getRealMass(), and initialize().

00468                                                           {
00469     switch(t) {
00470       case Proton:
00471         return theRealProtonMass;
00472         break;
00473       case Neutron:
00474         return theRealNeutronMass;
00475         break;
00476       case PiPlus:
00477       case PiMinus:
00478         return theRealChargedPiMass;
00479         break;
00480       case PiZero:
00481         return theRealPiZeroMass;
00482         break;
00483       default:
00484         ERROR("Particle::getRealMass : Unknown particle type." << std::endl);
00485         return 0.0;
00486         break;
00487     }
00488   }

static G4double G4INCL::ParticleTable::getSeparationEnergyINCL ( const ParticleType  t,
const   G4int,
const   G4int 
) [inline, static]

Return INCL's default separation energy.

Definition at line 227 of file G4INCLParticleTable.hh.

References ERROR, G4INCL::Neutron, and G4INCL::Proton.

Referenced by initialize().

00227                                                                                                         {
00228       if(t==Proton)
00229         return theINCLProtonSeparationEnergy;
00230       else if(t==Neutron)
00231         return theINCLNeutronSeparationEnergy;
00232       else {
00233         ERROR("ParticleTable::getSeparationEnergyINCL : Unknown particle type." << std::endl);
00234         return 0.0;
00235       }
00236     }

static G4double G4INCL::ParticleTable::getSeparationEnergyReal ( const ParticleType  t,
const G4int  A,
const G4int  Z 
) [inline, static]

Return the real separation energy.

Definition at line 239 of file G4INCLParticleTable.hh.

References ERROR, G4INCL::Neutron, and G4INCL::Proton.

Referenced by initialize().

00239                                                                                                 {
00240       // Real separation energies for all nuclei
00241       if(t==Proton)
00242         return (*getTableParticleMass)(Proton) + (*getTableMass)(A-1,Z-1) - (*getTableMass)(A,Z);
00243       else if(t==Neutron)
00244         return (*getTableParticleMass)(Neutron) + (*getTableMass)(A-1,Z) - (*getTableMass)(A,Z);
00245       else {
00246         ERROR("ParticleTable::getSeparationEnergyReal : Unknown particle type." << std::endl);
00247         return 0.0;
00248       }
00249     }

static G4double G4INCL::ParticleTable::getSeparationEnergyRealForLight ( const ParticleType  t,
const G4int  A,
const G4int  Z 
) [inline, static]

Return the real separation energy only for light nuclei.

Definition at line 252 of file G4INCLParticleTable.hh.

Referenced by initialize().

00252                                                                                                         {
00253       // Real separation energies for light nuclei, fixed values for heavy nuclei
00254       if(Z<clusterTableZSize && A<clusterTableASize)
00255         return getSeparationEnergyReal(t, A, Z);
00256       else
00257         return getSeparationEnergyINCL(t, A, Z);
00258     }

std::string G4INCL::ParticleTable::getShortName ( const G4int  A,
const G4int  Z 
) [static]

Get the short INCL name of the ion.

Definition at line 393 of file G4INCLParticleTable.cc.

References getElementName().

00393                                                                     {
00394     std::stringstream stream;
00395     stream << getElementName(Z);
00396     if(A>0)
00397       stream << A;
00398     return stream.str();
00399   }

std::string G4INCL::ParticleTable::getShortName ( const ParticleSpecies  s  )  [static]

Get the short INCL name of the particle.

Definition at line 373 of file G4INCLParticleTable.cc.

References G4INCL::Composite, getShortName(), G4INCL::ParticleSpecies::theA, G4INCL::ParticleSpecies::theType, and G4INCL::ParticleSpecies::theZ.

00373                                                                {
00374     if(s.theType==Composite)
00375       return getShortName(s.theA,s.theZ);
00376     else
00377       return getShortName(s.theType);
00378   }

std::string G4INCL::ParticleTable::getShortName ( const ParticleType  t  )  [static]

Get the short INCL name of the particle.

Definition at line 426 of file G4INCLParticleTable.cc.

References G4INCL::Composite, G4INCL::DeltaMinus, G4INCL::DeltaPlus, G4INCL::DeltaPlusPlus, G4INCL::DeltaZero, G4INCL::Neutron, G4INCL::PiMinus, G4INCL::PiPlus, G4INCL::PiZero, and G4INCL::Proton.

Referenced by getShortName().

00426                                                             {
00427     if(p == G4INCL::Proton) {
00428       return std::string("p");
00429     } else if(p == G4INCL::Neutron) {
00430       return std::string("n");
00431     } else if(p == G4INCL::DeltaPlusPlus) {
00432       return std::string("d++");
00433     } else if(p == G4INCL::DeltaPlus) {
00434       return std::string("d+");
00435     } else if(p == G4INCL::DeltaZero) {
00436       return std::string("d0");
00437     } else if(p == G4INCL::DeltaMinus) {
00438       return std::string("d-");
00439     } else if(p == G4INCL::PiPlus) {
00440       return std::string("pi+");
00441     } else if(p == G4INCL::PiZero) {
00442       return std::string("pi0");
00443     } else if(p == G4INCL::PiMinus) {
00444       return std::string("pi-");
00445     } else if(p == G4INCL::Composite) {
00446       return std::string("comp");
00447     }
00448     return std::string("unknown");
00449   }

G4double G4INCL::ParticleTable::getSurfaceDiffuseness ( const G4int  A,
const G4int  Z 
) [static]

Definition at line 604 of file G4INCLParticleTable.cc.

References ERROR.

Referenced by G4INCL::NuclearDensityFactory::createRCDFTable(), G4INCL::NuclearDensityFactory::createRPCorrelationTable(), getMaximumNuclearRadius(), and getNuclearRadius().

00604                                                                                  {
00605 #else
00606   G4double ParticleTable::getSurfaceDiffuseness(const G4int A, const G4int Z) {
00607 #endif // INCLXX_IN_GEANT4_MODE
00608 
00609     if(A >= 28) {
00610       return 1.63e-4 * A + 0.510;
00611     } else if(A < 28 && A >= 19) {
00612       return mediumDiffuseness[A-1];
00613     } else if(A < 19 && A >= 6) {
00614       return mediumDiffuseness[A-1];
00615     } else if(A < 6 && A >= 2) {
00616       ERROR("ParticleTable::getSurfaceDiffuseness: was called for A = " << A << " Z = " << Z << std::endl);
00617       return 0.0;
00618     } else {
00619       ERROR("ParticleTable::getSurfaceDiffuseness: No diffuseness for nucleus A = " << A << " Z = " << Z << std::endl);
00620       return 0.0;
00621     }
00622   }

static G4double G4INCL::ParticleTable::getTableQValue ( const G4int  A1,
const G4int  Z1,
const G4int  A2,
const G4int  Z2,
const G4int  A3,
const G4int  Z3 
) [inline, static]

Get Q-value (in MeV/c^2).

Uses the getTableMass function to compute the Q-value for the following reaction:

\[ (A_1,Z_1) + (A_2, Z_2) --> (A_3,Z_3) + (A1+A2-A3,Z1+Z2-Z3) \]

Definition at line 143 of file G4INCLParticleTable.hh.

00143                                                                                                                                    {
00144       return getTableMass(A1,Z1) + getTableMass(A2,Z2) - getTableMass(A3,Z3) - getTableMass(A1+A2-A3,Z1+Z2-Z3);
00145     }

static G4double G4INCL::ParticleTable::getTableQValue ( const G4int  A1,
const G4int  Z1,
const G4int  A2,
const G4int  Z2 
) [inline, static]

Get Q-value (in MeV/c^2).

Uses the getTableMass function to compute the Q-value for the following reaction:

\[ (A_1,Z_1) + (A_2, Z_2) --> (A_1+A_2,Z_1+Z_2) \]

Definition at line 133 of file G4INCLParticleTable.hh.

Referenced by G4INCL::Particle::getEmissionQValueCorrection(), and G4INCL::Particle::getTransferQValueCorrection().

00133                                                                                                    {
00134       return getTableMass(A1,Z1) + getTableMass(A2,Z2) - getTableMass(A1+A2,Z1+Z2);
00135     }

static G4double G4INCL::ParticleTable::getTableSpeciesMass ( const ParticleSpecies p  )  [inline, static]

Definition at line 153 of file G4INCLParticleTable.hh.

References G4INCL::Composite, G4INCL::ParticleSpecies::theA, G4INCL::ParticleSpecies::theType, and G4INCL::ParticleSpecies::theZ.

00153                                                                   {
00154       if(p.theType == Composite)
00155         return (*getTableMass)(p.theA, p.theZ);
00156       else
00157         return (*getTableParticleMass)(p.theType);
00158     }

void G4INCL::ParticleTable::initialize ( Config const *const   theConfig = 0  )  [static]

Initialize the particle table.

Definition at line 304 of file G4INCLParticleTable.cc.

References FATAL, G4ParticleTable::FindParticle(), getINCLMass(), G4INCL::Config::getINCLXXDataFilePath(), G4ParticleTable::GetIonTable(), G4ParticleTable::GetParticleTable(), G4ParticleDefinition::GetPDGMass(), getRealMass(), getSeparationEnergy, getSeparationEnergyINCL(), getSeparationEnergyReal(), getSeparationEnergyRealForLight(), G4INCL::Config::getSeparationEnergyType(), getTableMass, getTableParticleMass, G4INCL::Config::getUseRealMasses(), G4INCL::INCLSeparationEnergy, G4INCL::RealForLightSeparationEnergy, G4INCL::RealSeparationEnergy, and theG4IonTable.

Referenced by G4INCL::INCL::INCL().

00304                                                                       {
00305     protonMass = theINCLNucleonMass;
00306     neutronMass = theINCLNucleonMass;
00307     piPlusMass = theINCLPionMass;
00308     piMinusMass = theINCLPionMass;
00309     piZeroMass = theINCLPionMass;
00310 
00311 #ifndef INCLXX_IN_GEANT4_MODE
00312     std::string dataFilePath;
00313     if(theConfig)
00314       dataFilePath = theConfig->getINCLXXDataFilePath();
00315     readRealMasses(dataFilePath);
00316 #endif
00317 
00318     if(theConfig && theConfig->getUseRealMasses()) {
00319       getTableMass = ParticleTable::getRealMass;
00320       getTableParticleMass = ParticleTable::getRealMass;
00321     } else {
00322       getTableMass = ParticleTable::getINCLMass;
00323       getTableParticleMass = ParticleTable::getINCLMass;
00324     }
00325 #ifdef INCLXX_IN_GEANT4_MODE
00326     G4ParticleTable *theG4ParticleTable = G4ParticleTable::GetParticleTable();
00327     theG4IonTable = theG4ParticleTable->GetIonTable();
00328     theRealProtonMass = theG4ParticleTable->FindParticle("proton")->GetPDGMass() / MeV;
00329     theRealNeutronMass = theG4ParticleTable->FindParticle("neutron")->GetPDGMass() / MeV;
00330     theRealChargedPiMass = theG4ParticleTable->FindParticle("pi+")->GetPDGMass() / MeV;
00331     theRealPiZeroMass = theG4ParticleTable->FindParticle("pi0")->GetPDGMass() / MeV;
00332 #endif
00333 
00334     // Initialise the separation-energy function
00335     if(!theConfig || theConfig->getSeparationEnergyType()==INCLSeparationEnergy)
00336       getSeparationEnergy = ParticleTable::getSeparationEnergyINCL;
00337     else if(theConfig->getSeparationEnergyType()==RealSeparationEnergy)
00338       getSeparationEnergy = ParticleTable::getSeparationEnergyReal;
00339     else if(theConfig->getSeparationEnergyType()==RealForLightSeparationEnergy)
00340       getSeparationEnergy = ParticleTable::getSeparationEnergyRealForLight;
00341     else {
00342       FATAL("Unrecognized separation-energy type in ParticleTable initialization: " << theConfig->getSeparationEnergyType() << std::endl);
00343     }
00344 
00345   }

G4int G4INCL::ParticleTable::parseIUPACElement ( std::string const &  pS  )  [static]

Parse a IUPAC element name.

Note: this function is UGLY. Look at it at your own peril.

Parameters:
pS a normalised string (lowercase)
Returns:
the charge number of the nuclide, or zero on fail

Definition at line 686 of file G4INCLParticleTable.cc.

00686                                                            {
00687     // Normalise to lower case
00688     std::string elementName(s);
00689     std::transform(elementName.begin(), elementName.end(), elementName.begin(), ::tolower);
00690     // Return 0 if the element name contains anything but IUPAC digits
00691     if(elementName.find_first_not_of(elementIUPACDigits)!=std::string::npos)
00692       return 0;
00693     std::transform(elementName.begin(), elementName.end(), elementName.begin(), ParticleTable::iupacToInt);
00694     std::stringstream elementStream(elementName);
00695     G4int Z;
00696     elementStream >> Z;
00697     return Z;
00698   }

static void G4INCL::ParticleTable::setNeutronSeparationEnergy ( const G4double  s  )  [inline, static]

Setter for protonSeparationEnergy.

Definition at line 270 of file G4INCLParticleTable.hh.

Referenced by G4INCL::Nucleus::Nucleus().

00270 { neutronSeparationEnergy  = s; }

static void G4INCL::ParticleTable::setProtonSeparationEnergy ( const G4double  s  )  [inline, static]

Setter for protonSeparationEnergy.

Definition at line 267 of file G4INCLParticleTable.hh.

Referenced by G4INCL::Nucleus::Nucleus().

00267 { protonSeparationEnergy = s; }


Field Documentation

const ParticleTable::ClusterDecayType G4INCL::ParticleTable::clusterDecayMode [static]

Initial value:

  {

    {StableCluster, StableCluster,   NeutronDecay, NeutronUnbound, NeutronUnbound,  NeutronUnbound, NeutronUnbound,  NeutronUnbound, NeutronUnbound, NeutronUnbound,  NeutronUnbound, NeutronUnbound, NeutronUnbound},
    {StableCluster, StableCluster,  StableCluster,  StableCluster,   NeutronDecay, TwoNeutronDecay,   NeutronDecay, TwoNeutronDecay, NeutronUnbound, NeutronUnbound,  NeutronUnbound, NeutronUnbound, NeutronUnbound},
    {StableCluster, StableCluster,    ProtonDecay,  StableCluster,  StableCluster,    NeutronDecay,  StableCluster,    NeutronDecay,  StableCluster,   NeutronDecay, TwoNeutronDecay, NeutronUnbound, NeutronUnbound},
    {StableCluster, StableCluster,  StableCluster,  ProtonUnbound,    ProtonDecay,     ProtonDecay,  StableCluster,   StableCluster,  StableCluster,  StableCluster,    NeutronDecay,  StableCluster,   NeutronDecay},
    {StableCluster, StableCluster,  StableCluster,  StableCluster,  ProtonUnbound,     ProtonDecay, TwoProtonDecay,   StableCluster,     AlphaDecay,  StableCluster,   StableCluster,  StableCluster,  StableCluster},
    {StableCluster, StableCluster,  StableCluster,  StableCluster,  StableCluster,   ProtonUnbound, TwoProtonDecay,     ProtonDecay,  StableCluster,    ProtonDecay,   StableCluster,  StableCluster,  StableCluster},
    {StableCluster, StableCluster,  StableCluster,  StableCluster,  StableCluster,   StableCluster,  ProtonUnbound,   ProtonUnbound, TwoProtonDecay,  StableCluster,   StableCluster,  StableCluster,  StableCluster},
    {StableCluster, StableCluster,  StableCluster,  StableCluster,  StableCluster,   StableCluster,  StableCluster,   ProtonUnbound,  ProtonUnbound,  ProtonUnbound,     ProtonDecay,    ProtonDecay,  StableCluster},
    {StableCluster, StableCluster,  StableCluster,  StableCluster,  StableCluster,   StableCluster,  StableCluster,   StableCluster,  ProtonUnbound,  ProtonUnbound,   ProtonUnbound,  ProtonUnbound,    ProtonDecay}
  }

Definition at line 323 of file G4INCLParticleTable.hh.

Referenced by G4INCL::ClusterDecay::isStable().

const G4double G4INCL::ParticleTable::clusterPhaseSpaceCut [static]

Initial value:

 {0.0, 70000.0, 180000.0,
                                                            90000.0, 90000.0,
                                                            128941.0 ,145607.0,
                                                            161365.0, 176389.0,
                                                            190798.0, 204681.0,
                                                            218109.0, 231135.0}

Definition at line 303 of file G4INCLParticleTable.hh.

Referenced by G4INCL::ClusteringModelIntercomparison::getCluster().

const G4double G4INCL::ParticleTable::clusterPosFact [static]

Initial value:

 {0.0, 1.0, 0.5,
                                                      0.33333, 0.25,
                                                      0.2, 0.16667,
                                                      0.14286, 0.125,
                                                      0.11111, 0.1,
                                                      0.09091, 0.083333}
Precomputed factor 1.0/A

Definition at line 299 of file G4INCLParticleTable.hh.

Referenced by G4INCL::ClusterUtils::getNewPositionVector().

const G4double G4INCL::ParticleTable::clusterPosFact2 [static]

Initial value:

 {0.0, 1.0, 0.25,
                                                       0.11111, 0.0625,
                                                       0.04, 0.0277778,
                                                       0.020408, 0.015625,
                                                       0.012346, 0.01,
                                                       0.0082645, 0.0069444}
Precomputed factor (1.0/A)^2

Definition at line 300 of file G4INCLParticleTable.hh.

Referenced by G4INCL::ClusteringModelIntercomparison::getCluster(), and G4INCL::ClusterUtils::getPhaseSpace().

const G4int G4INCL::ParticleTable::clusterTableASize = ParticleTable::maxClusterMass+1 [static]

Definition at line 298 of file G4INCLParticleTable.hh.

Referenced by getNuclearRadius().

const G4int G4INCL::ParticleTable::clusterTableZSize = ParticleTable::maxClusterCharge+1 [static]

Definition at line 297 of file G4INCLParticleTable.hh.

Referenced by getRadiusParameter().

const G4int G4INCL::ParticleTable::clusterZMax = {0, 0, 1, 2, 3, 3, 5, 5, 6, 6, 7, 7, 8} [static]

Definition at line 302 of file G4INCLParticleTable.hh.

Referenced by G4INCL::ClusteringModelIntercomparison::ClusteringModelIntercomparison().

const G4int G4INCL::ParticleTable::clusterZMin = {0, 0, 1, 1, 1, 1, 1, 1, 2, 2, 2, 3, 3} [static]

Definition at line 301 of file G4INCLParticleTable.hh.

Referenced by G4INCL::ClusteringModelIntercomparison::ClusteringModelIntercomparison().

const G4double G4INCL::ParticleTable::effectiveDeltaDecayThreshold [static]

Initial value:

    ParticleTable::theRealNeutronMass + ParticleTable::theRealChargedPiMass
    + 0.5

Definition at line 292 of file G4INCLParticleTable.hh.

Referenced by G4INCL::InteractionAvatar::enforceEnergyConservation(), and G4INCL::InteractionAvatar::postInteraction().

const G4double G4INCL::ParticleTable::effectiveDeltaMass = 1232.0 [static]

Definition at line 290 of file G4INCLParticleTable.hh.

const G4double G4INCL::ParticleTable::effectiveNucleonMass = 938.2796 [static]

Definition at line 288 of file G4INCLParticleTable.hh.

Referenced by G4INCL::DeltaDecayChannel::computeDecayTime(), G4INCL::Nucleus::computeTotalEnergy(), G4INCL::CrossSections::deltaProduction(), G4INCL::ElasticChannel::getFinalState(), G4INCL::DeltaProductionChannel::getFinalState(), and G4INCL::CrossSections::recombination().

const G4double G4INCL::ParticleTable::effectiveNucleonMass2 = 8.8036860777616e5 [static]

Definition at line 289 of file G4INCLParticleTable.hh.

Referenced by G4INCL::DeltaProductionChannel::getFinalState(), and G4INCL::CrossSections::recombination().

const G4double G4INCL::ParticleTable::effectivePionMass = 138.0 [static]

Definition at line 291 of file G4INCLParticleTable.hh.

Referenced by G4INCL::DeltaDecayChannel::computeDecayTime(), and G4INCL::CrossSections::deltaProduction().

const G4int G4INCL::ParticleTable::elementTableSize = 113 [static]

Definition at line 286 of file G4INCLParticleTable.hh.

Referenced by getElementName().

const G4double G4INCL::ParticleTable::eSquared = 1.439964 [static]

Coulomb conversion factor, in MeV*fm.

\[ e^2/(4 pi epsilon_0) \]

Definition at line 329 of file G4INCLParticleTable.hh.

Referenced by G4INCL::Nucleus::getCoulombRadius().

ParticleTable::SeparationEnergyFn G4INCL::ParticleTable::getSeparationEnergy [static]

Static pointer to the separation-energy function.

Definition at line 162 of file G4INCLParticleTable.hh.

Referenced by initialize().

ParticleTable::NuclearMassFn G4INCL::ParticleTable::getTableMass [static]

Static pointer to the mass function for nuclei.

Definition at line 150 of file G4INCLParticleTable.hh.

Referenced by G4INCL::ProjectileRemnant::addMostDynamicalSpectators(), G4INCL::Nucleus::computeRecoilKinematics(), G4INCL::Nucleus::finalizeProjectileRemnant(), G4INCL::KinematicsUtils::gammaFromKineticEnergy(), G4INCL::Nucleus::getConservationBalance(), G4INCL::Particle::getEmissionQValueCorrection(), G4INCL::ParticleEntryChannel::getFinalState(), G4INCL::Particle::getTableMass(), initialize(), G4INCL::StandardPropagationModel::shootComposite(), and G4INCL::StandardPropagationModel::shootParticle().

ParticleTable::ParticleMassFn G4INCL::ParticleTable::getTableParticleMass [static]

Static pointer to the mass function for particles.

Definition at line 151 of file G4INCLParticleTable.hh.

Referenced by G4INCL::KinematicsUtils::gammaFromKineticEnergy(), G4INCL::Particle::getTableMass(), initialize(), and G4INCL::StandardPropagationModel::shootParticle().

const G4int G4INCL::ParticleTable::maxClusterCharge = 8 [static]

Definition at line 295 of file G4INCLParticleTable.hh.

const G4int G4INCL::ParticleTable::maxClusterMass = 12 [static]

Definition at line 294 of file G4INCLParticleTable.hh.

Referenced by G4INCL::ClusteringModelIntercomparison::ClusteringModelIntercomparison().

G4IonTable * G4INCL::ParticleTable::theG4IonTable [static]

Definition at line 306 of file G4INCLParticleTable.hh.

Referenced by getRealMass(), and initialize().


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:54:09 2013 for Geant4 by  doxygen 1.4.7