G4PhysicsVector Class Reference

#include <G4PhysicsVector.hh>

Inheritance diagram for G4PhysicsVector:

G4LPhysicsFreeVector G4PhysicsFreeVector G4PhysicsLinearVector G4PhysicsLnVector G4PhysicsLogVector G4PhysicsOrderedFreeVector

Public Member Functions

 G4PhysicsVector (G4bool spline=false)
 G4PhysicsVector (const G4PhysicsVector &)
G4PhysicsVectoroperator= (const G4PhysicsVector &)
virtual ~G4PhysicsVector ()
void * operator new (size_t)
void operator delete (void *)
G4double Value (G4double theEnergy)
G4double GetValue (G4double theEnergy, G4bool &isOutRange)
G4int operator== (const G4PhysicsVector &right) const
G4int operator!= (const G4PhysicsVector &right) const
G4double operator[] (const size_t binNumber) const
G4double operator() (const size_t binNumber) const
void PutValue (size_t index, G4double theValue)
virtual void ScaleVector (G4double factorE, G4double factorV)
G4double Energy (size_t index) const
G4double GetMaxEnergy () const
virtual G4double GetLowEdgeEnergy (size_t binNumber) const
size_t GetVectorLength () const
void FillSecondDerivatives ()
void ComputeSecDerivatives ()
void ComputeSecondDerivatives (G4double firstPointDerivative, G4double endPointDerivative)
G4bool IsFilledVectorExist () const
G4PhysicsVectorType GetType () const
void SetSpline (G4bool)
virtual G4bool Store (std::ofstream &fOut, G4bool ascii=false)
virtual G4bool Retrieve (std::ifstream &fIn, G4bool ascii=false)
G4double GetLastEnergy () const
G4double GetLastValue () const
size_t GetLastBin () const
void SetVerboseLevel (G4int value)
G4int GetVerboseLevel (G4int)

Protected Member Functions

virtual size_t FindBinLocation (G4double theEnergy) const =0
void DeleteData ()
void CopyData (const G4PhysicsVector &vec)

Protected Attributes

G4PhysicsVectorType type
G4double edgeMin
G4double edgeMax
size_t numberOfNodes
G4PhysicsVectorCachecache
G4PVDataVector dataVector
G4PVDataVector binVector
G4PVDataVector secDerivative
G4double dBin
G4double baseBin
G4int verboseLevel

Friends

std::ostream & operator<< (std::ostream &, const G4PhysicsVector &)

Detailed Description

Definition at line 76 of file G4PhysicsVector.hh.


Constructor & Destructor Documentation

G4PhysicsVector::G4PhysicsVector ( G4bool  spline = false  ) 

Definition at line 63 of file G4PhysicsVector.cc.

References cache.

00064  : type(T_G4PhysicsVector),
00065    edgeMin(0.), edgeMax(0.), numberOfNodes(0),
00066    useSpline(spline), 
00067    dBin(0.), baseBin(0.),
00068    verboseLevel(0)
00069 {
00070   cache      = new G4PhysicsVectorCache();
00071 }

G4PhysicsVector::G4PhysicsVector ( const G4PhysicsVector  ) 

Definition at line 82 of file G4PhysicsVector.cc.

References baseBin, cache, CopyData(), dBin, DeleteData(), and verboseLevel.

00083 {
00084   cache        = new G4PhysicsVectorCache();
00085   dBin         = right.dBin;
00086   baseBin      = right.baseBin;
00087   verboseLevel = right.verboseLevel;
00088 
00089   DeleteData();
00090   CopyData(right);
00091 }

G4PhysicsVector::~G4PhysicsVector (  )  [virtual]

Definition at line 75 of file G4PhysicsVector.cc.

References cache.

00076 {
00077   delete cache; cache =0;
00078 }


Member Function Documentation

void G4PhysicsVector::ComputeSecDerivatives (  ) 

Definition at line 441 of file G4PhysicsVector.cc.

References binVector, dataVector, CLHEP::detail::n, numberOfNodes, and secDerivative.

Referenced by ComputeSecondDerivatives(), and FillSecondDerivatives().

00443 {
00444   if(!SplinePossible())  { return; }
00445 
00446   if(3 > numberOfNodes)  // cannot compute derivatives for less than 4 bins
00447   {
00448     useSpline = false;
00449     return;
00450   }
00451 
00452   size_t n = numberOfNodes-1;
00453 
00454   for(size_t i=1; i<n; ++i)
00455   {
00456     secDerivative[i] =
00457       3.0*((dataVector[i+1]-dataVector[i])/(binVector[i+1]-binVector[i]) -
00458            (dataVector[i]-dataVector[i-1])/(binVector[i]-binVector[i-1]))
00459       /(binVector[i+1]-binVector[i-1]);
00460   }
00461   secDerivative[n] = secDerivative[n-1];
00462   secDerivative[0] = secDerivative[1];
00463 }

void G4PhysicsVector::ComputeSecondDerivatives ( G4double  firstPointDerivative,
G4double  endPointDerivative 
)

Definition at line 308 of file G4PhysicsVector.cc.

References binVector, ComputeSecDerivatives(), dataVector, CLHEP::detail::n, numberOfNodes, and secDerivative.

00314 {
00315   if(4 > numberOfNodes)   // cannot compute derivatives for less than 4 bins
00316   {
00317     ComputeSecDerivatives();
00318     return;
00319   }
00320 
00321   if(!SplinePossible()) { return; }
00322 
00323   G4int n = numberOfNodes-1;
00324 
00325   G4double* u = new G4double [n];
00326   
00327   G4double p, sig, un;
00328 
00329   u[0] = (6.0/(binVector[1]-binVector[0]))
00330     * ((dataVector[1]-dataVector[0])/(binVector[1]-binVector[0])
00331        - firstPointDerivative);
00332  
00333   secDerivative[0] = - 0.5;
00334 
00335   // Decomposition loop for tridiagonal algorithm. secDerivative[i]
00336   // and u[i] are used for temporary storage of the decomposed factors.
00337 
00338   for(G4int i=1; i<n; ++i)
00339   {
00340     sig = (binVector[i]-binVector[i-1]) / (binVector[i+1]-binVector[i-1]);
00341     p = sig*(secDerivative[i-1]) + 2.0;
00342     secDerivative[i] = (sig - 1.0)/p;
00343     u[i] = (dataVector[i+1]-dataVector[i])/(binVector[i+1]-binVector[i])
00344          - (dataVector[i]-dataVector[i-1])/(binVector[i]-binVector[i-1]);
00345     u[i] = 6.0*u[i]/(binVector[i+1]-binVector[i-1]) - sig*u[i-1]/p;
00346   }
00347 
00348   sig = (binVector[n-1]-binVector[n-2]) / (binVector[n]-binVector[n-2]);
00349   p = sig*secDerivative[n-2] + 2.0;
00350   un = (6.0/(binVector[n]-binVector[n-1]))
00351     *(endPointDerivative - 
00352       (dataVector[n]-dataVector[n-1])/(binVector[n]-binVector[n-1])) - u[n-1]/p;
00353   secDerivative[n] = un/(secDerivative[n-1] + 2.0);
00354 
00355   // The back-substitution loop for the triagonal algorithm of solving
00356   // a linear system of equations.
00357    
00358   for(G4int k=n-1; k>0; --k)
00359   {
00360     secDerivative[k] *= 
00361       (secDerivative[k+1] - 
00362        u[k]*(binVector[k+1]-binVector[k-1])/(binVector[k+1]-binVector[k]));
00363   }
00364   secDerivative[0] = 0.5*(u[0] - secDerivative[1]);
00365 
00366   delete [] u;
00367 }

void G4PhysicsVector::CopyData ( const G4PhysicsVector vec  )  [protected]

Definition at line 130 of file G4PhysicsVector.cc.

References binVector, cache, dataVector, edgeMax, edgeMin, GetLastBin(), GetLastEnergy(), GetLastValue(), G4PhysicsVectorCache::lastBin, G4PhysicsVectorCache::lastEnergy, G4PhysicsVectorCache::lastValue, numberOfNodes, secDerivative, type, and useSpline.

Referenced by G4PhysicsVector(), and operator=().

00131 {
00132   type = vec.type;
00133   edgeMin = vec.edgeMin;
00134   edgeMax = vec.edgeMax;
00135   numberOfNodes = vec.numberOfNodes;
00136   cache->lastEnergy = vec.GetLastEnergy();
00137   cache->lastValue = vec.GetLastValue();
00138   cache->lastBin = vec.GetLastBin();
00139   useSpline = vec.useSpline;
00140 
00141   size_t i;
00142   dataVector.clear();
00143   for(i=0; i<(vec.dataVector).size(); i++){ 
00144     dataVector.push_back( (vec.dataVector)[i] );
00145   }
00146   binVector.clear();
00147   for(i=0; i<(vec.binVector).size(); i++){ 
00148     binVector.push_back( (vec.binVector)[i] );
00149   }
00150   secDerivative.clear();
00151   for(i=0; i<(vec.secDerivative).size(); i++){ 
00152     secDerivative.push_back( (vec.secDerivative)[i] );
00153   }
00154 }

void G4PhysicsVector::DeleteData (  )  [protected]

Definition at line 123 of file G4PhysicsVector.cc.

References secDerivative.

Referenced by G4PhysicsVector(), and operator=().

00124 {
00125   secDerivative.clear();
00126 }

G4double G4PhysicsVector::Energy ( size_t  index  )  const [inline]

Definition at line 106 of file G4PhysicsVector.icc.

References binVector.

Referenced by G4EmCorrections::BarkasCorrection(), G4KokoulinMuonNuclearXS::BuildCrossSectionTable(), G4LossTableBuilder::BuildInverseRangeTable(), G4ePolarizedIonisation::BuildPhysicsTable(), G4LossTableBuilder::BuildRangeTable(), G4LossTableBuilder::BuildTableForModel(), G4Scintillation::BuildThePhysicsTable(), G4LivermoreRayleighModel::ComputeCrossSectionPerAtom(), G4EmModelManager::FillDEDXVector(), G4EmModelManager::FillLambdaVector(), G4NeutronInelasticXS::GetElementCrossSection(), G4NeutronElasticXS::GetElementCrossSection(), G4NeutronCaptureXS::GetElementCrossSection(), G4NeutronCaptureXS::GetIsoCrossSection(), and G4GDMLWriteMaterials::PropertyVectorWrite().

00107 {
00108   return binVector[binNumber];
00109 }

void G4PhysicsVector::FillSecondDerivatives (  ) 

Definition at line 371 of file G4PhysicsVector.cc.

References binVector, ComputeSecDerivatives(), dataVector, CLHEP::detail::n, numberOfNodes, and secDerivative.

Referenced by G4VEnergyLossProcess::BuildDEDXTable(), G4LossTableBuilder::BuildDEDXTable(), G4LossTableBuilder::BuildInverseRangeTable(), G4VEnergyLossProcess::BuildLambdaTable(), G4LossTableBuilder::BuildRangeTable(), and G4LossTableBuilder::BuildTableForModel().

00375 {  
00376   if(5 > numberOfNodes)  // cannot compute derivatives for less than 4 points
00377   {
00378     ComputeSecDerivatives();
00379     return;
00380   }
00381 
00382   if(!SplinePossible()) { return; }
00383  
00384   G4int n = numberOfNodes-1;
00385 
00386   G4double* u = new G4double [n];
00387   
00388   G4double p, sig;
00389 
00390   u[1] = ((dataVector[2]-dataVector[1])/(binVector[2]-binVector[1]) -
00391           (dataVector[1]-dataVector[0])/(binVector[1]-binVector[0]));
00392   u[1] = 6.0*u[1]*(binVector[2]-binVector[1])
00393     / ((binVector[2]-binVector[0])*(binVector[2]-binVector[0]));
00394  
00395   // Decomposition loop for tridiagonal algorithm. secDerivative[i]
00396   // and u[i] are used for temporary storage of the decomposed factors.
00397 
00398   secDerivative[1] = (2.0*binVector[1]-binVector[0]-binVector[2])
00399     / (2.0*binVector[2]-binVector[0]-binVector[1]);
00400 
00401   for(G4int i=2; i<n-1; ++i)
00402   {
00403     sig = (binVector[i]-binVector[i-1]) / (binVector[i+1]-binVector[i-1]);
00404     p = sig*secDerivative[i-1] + 2.0;
00405     secDerivative[i] = (sig - 1.0)/p;
00406     u[i] = (dataVector[i+1]-dataVector[i])/(binVector[i+1]-binVector[i])
00407          - (dataVector[i]-dataVector[i-1])/(binVector[i]-binVector[i-1]);
00408     u[i] = (6.0*u[i]/(binVector[i+1]-binVector[i-1])) - sig*u[i-1]/p;
00409   }
00410 
00411   sig = (binVector[n-1]-binVector[n-2]) / (binVector[n]-binVector[n-2]);
00412   p = sig*secDerivative[n-3] + 2.0;
00413   u[n-1] = (dataVector[n]-dataVector[n-1])/(binVector[n]-binVector[n-1])
00414     - (dataVector[n-1]-dataVector[n-2])/(binVector[n-1]-binVector[n-2]);
00415   u[n-1] = 6.0*sig*u[n-1]/(binVector[n]-binVector[n-2])
00416     - (2.0*sig - 1.0)*u[n-2]/p;  
00417 
00418   p = (1.0+sig) + (2.0*sig-1.0)*secDerivative[n-2];
00419   secDerivative[n-1] = u[n-1]/p;
00420 
00421   // The back-substitution loop for the triagonal algorithm of solving
00422   // a linear system of equations.
00423    
00424   for(G4int k=n-2; k>1; --k)
00425   {
00426     secDerivative[k] *= 
00427       (secDerivative[k+1] - 
00428        u[k]*(binVector[k+1]-binVector[k-1])/(binVector[k+1]-binVector[k]));
00429   }
00430   secDerivative[n] = (secDerivative[n-1] - (1.0-sig)*secDerivative[n-2])/sig;
00431   sig = 1.0 - ((binVector[2]-binVector[1])/(binVector[2]-binVector[0]));
00432   secDerivative[1] *= (secDerivative[2] - u[1]/(1.0-sig));
00433   secDerivative[0]  = (secDerivative[1] - sig*secDerivative[2])/(1.0-sig);
00434 
00435   delete [] u;
00436 }

virtual size_t G4PhysicsVector::FindBinLocation ( G4double  theEnergy  )  const [protected, pure virtual]

Implemented in G4PhysicsFreeVector, G4PhysicsLinearVector, G4PhysicsLnVector, and G4PhysicsLogVector.

size_t G4PhysicsVector::GetLastBin (  )  const [inline]

Definition at line 84 of file G4PhysicsVector.icc.

References cache, and G4PhysicsVectorCache::lastBin.

Referenced by CopyData().

00085 {
00086    return cache->lastBin;
00087 }

G4double G4PhysicsVector::GetLastEnergy (  )  const [inline]

Definition at line 72 of file G4PhysicsVector.icc.

References cache, and G4PhysicsVectorCache::lastEnergy.

Referenced by CopyData().

00073 {
00074   return cache->lastEnergy;     
00075 }

G4double G4PhysicsVector::GetLastValue (  )  const [inline]

Definition at line 78 of file G4PhysicsVector.icc.

References cache, and G4PhysicsVectorCache::lastValue.

Referenced by CopyData().

00079 {
00080    return cache->lastValue;
00081 }

G4double G4PhysicsVector::GetLowEdgeEnergy ( size_t  binNumber  )  const [virtual]

Reimplemented in G4PhysicsOrderedFreeVector.

Definition at line 158 of file G4PhysicsVector.cc.

References binVector.

Referenced by G4VXTRenergyLoss::BuildAngleForEnergyBank(), G4VXTRenergyLoss::BuildAngleTable(), G4NuclNuclDiffuseElastic::BuildAngleTable(), G4DiffuseElastic::BuildAngleTable(), G4PolarizedCompton::BuildAsymmetryTable(), G4eplusPolarizedAnnihilation::BuildAsymmetryTable(), G4VXTRenergyLoss::BuildGlobalAngleTable(), G4VeLowEnergyLoss::BuildInverseRangeTable(), G4VRangeToEnergyConverter::BuildLossTable(), G4PAIPhotonModel::BuildPAIonisationTable(), G4PAIModel::BuildPAIonisationTable(), G4VRangeToEnergyConverter::BuildRangeVector(), G4AdjointCSManager::BuildTotalSigmaTables(), G4ForwardXrayTR::BuildXrayTRtables(), G4VRangeToEnergyConverter::ConvertCutToKineticEnergy(), G4PenelopeRayleighModel::DumpFormFactorTable(), G4PAIPhotonModel::GetAlongStepTransfer(), G4ForwardXrayTR::GetEnergyTR(), G4VXTRenergyLoss::GetMeanFreePath(), G4PAIPhotonModel::GetPostStepTransfer(), G4PAIModel::GetPostStepTransfer(), G4VXTRenergyLoss::GetXTRrandomEnergy(), G4eeToHadronsModel::Initialise(), G4InitXscPAI::IntegralCherenkov(), G4InitXscPAI::IntegralPAIdEdx(), G4InitXscPAI::IntegralPAIxSection(), G4InitXscPAI::IntegralPlasmon(), G4ForwardXrayTR::PostStepDoIt(), G4XnpTotalLowE::Print(), G4XnpElasticLowE::Print(), G4XNNElasticLowE::Print(), G4PAIModel::SampleFluctuations(), G4NuclNuclDiffuseElastic::SampleTableThetaCMS(), and G4DiffuseElastic::SampleTableThetaCMS().

00159 {
00160   return binVector[binNumber];
00161 }

G4double G4PhysicsVector::GetMaxEnergy (  )  const [inline]

Definition at line 114 of file G4PhysicsVector.icc.

References edgeMax.

00115 {
00116   return edgeMax;
00117 }

G4PhysicsVectorType G4PhysicsVector::GetType (  )  const [inline]

Definition at line 216 of file G4PhysicsVector.icc.

References type.

00217 {
00218   return type;
00219 }

G4double G4PhysicsVector::GetValue ( G4double  theEnergy,
G4bool isOutRange 
) [inline]

Definition at line 130 of file G4PhysicsVector.icc.

References Value().

Referenced by G4VEnergyLoss::BuildRangeCoeffATable(), G4VeLowEnergyLoss::BuildRangeCoeffATable(), G4VEnergyLoss::BuildRangeCoeffBTable(), G4VeLowEnergyLoss::BuildRangeCoeffBTable(), G4VEnergyLoss::BuildRangeCoeffCTable(), G4VeLowEnergyLoss::BuildRangeCoeffCTable(), G4eeToHadronsModel::ComputeCrossSectionPerElectron(), G4XResonance::CrossSection(), G4XnpTotalLowE::CrossSection(), G4XnpElasticLowE::CrossSection(), G4XNNElasticLowE::CrossSection(), G4PartialWidthTable::Dump(), G4eeToHadronsModel::GenerateCMPhoton(), G4ChargeExchangeProcess::GetElementCrossSection(), G4PolarizedCompton::GetMeanFreePath(), G4eeToHadronsModel::Initialise(), G4PolarizedCompton::PostStepGetPhysicalInteractionLength(), G4XnpTotalLowE::Print(), G4XnpElasticLowE::Print(), and G4XNNElasticLowE::Print().

00131 {
00132   return Value(theEnergy);
00133 }

size_t G4PhysicsVector::GetVectorLength (  )  const [inline]

Definition at line 122 of file G4PhysicsVector.icc.

References numberOfNodes.

Referenced by G4LossTableBuilder::BuildInverseRangeTable(), G4ePolarizedIonisation::BuildPhysicsTable(), G4LossTableBuilder::BuildRangeTable(), G4Scintillation::BuildThePhysicsTable(), G4AdjointCSManager::BuildTotalSigmaTables(), G4LivermoreRayleighModel::ComputeCrossSectionPerAtom(), G4LivermoreGammaConversionModel::ComputeCrossSectionPerAtom(), G4PenelopeRayleighModel::DumpFormFactorTable(), G4EmModelManager::FillDEDXVector(), G4EmModelManager::FillLambdaVector(), G4SPSRandomGenerator::GenRandEnergy(), G4SPSRandomGenerator::GenRandPhi(), G4SPSRandomGenerator::GenRandPosPhi(), G4SPSRandomGenerator::GenRandPosTheta(), G4SPSRandomGenerator::GenRandTheta(), G4SPSRandomGenerator::GenRandX(), G4SPSRandomGenerator::GenRandY(), G4SPSRandomGenerator::GenRandZ(), G4NeutronInelasticXS::GetElementCrossSection(), G4NeutronElasticXS::GetElementCrossSection(), G4PenelopeCrossSection::GetHardCrossSection(), G4PenelopeCrossSection::GetNormalizedShellCrossSection(), G4PenelopeCrossSection::GetShellCrossSection(), G4PenelopeCrossSection::GetSoftStoppingPower(), G4PenelopeCrossSection::GetTotalCrossSection(), G4eeToHadronsModel::Initialise(), G4GDMLWriteMaterials::PropertyVectorWrite(), and G4VEnergyLossProcess::SetLambdaTable().

00123 {
00124   return numberOfNodes;
00125 }

G4int G4PhysicsVector::GetVerboseLevel ( G4int   )  [inline]

Definition at line 240 of file G4PhysicsVector.icc.

References verboseLevel.

00241 {
00242    return verboseLevel;
00243 }

G4bool G4PhysicsVector::IsFilledVectorExist (  )  const [inline]

Definition at line 205 of file G4PhysicsVector.icc.

References numberOfNodes.

00206 {
00207   G4bool status=false;
00208 
00209   if(numberOfNodes > 0) { status=true; }
00210   return status;
00211 }

void G4PhysicsVector::operator delete ( void *   )  [inline]

Definition at line 57 of file G4PhysicsVector.icc.

References aPVAllocator.

00058 {
00059   aPVAllocator.FreeSingle((G4PhysicsVector*)aVector);
00060 }

void * G4PhysicsVector::operator new ( size_t   )  [inline]

Definition at line 50 of file G4PhysicsVector.icc.

References aPVAllocator.

00051 {
00052   void* aVector;
00053   aVector = (void*)aPVAllocator.MallocSingle();
00054   return aVector;
00055 }

G4int G4PhysicsVector::operator!= ( const G4PhysicsVector right  )  const

Definition at line 116 of file G4PhysicsVector.cc.

00117 {
00118   return (this != &right);
00119 }

G4double G4PhysicsVector::operator() ( const size_t  binNumber  )  const [inline]

Definition at line 98 of file G4PhysicsVector.icc.

References dataVector.

00099 {
00100   return dataVector[binNumber];
00101 }

G4PhysicsVector & G4PhysicsVector::operator= ( const G4PhysicsVector  ) 

Definition at line 95 of file G4PhysicsVector.cc.

References baseBin, CopyData(), dBin, DeleteData(), and verboseLevel.

00096 {
00097   if (&right==this)  { return *this; }
00098   dBin         = right.dBin;
00099   baseBin      = right.baseBin;
00100   verboseLevel = right.verboseLevel;
00101 
00102   DeleteData();
00103   CopyData(right);
00104   return *this;
00105 }

G4int G4PhysicsVector::operator== ( const G4PhysicsVector right  )  const

Definition at line 109 of file G4PhysicsVector.cc.

00110 {
00111   return (this == &right);
00112 }

G4double G4PhysicsVector::operator[] ( const size_t  binNumber  )  const [inline]

Definition at line 90 of file G4PhysicsVector.icc.

References dataVector.

00091 {
00092   return  dataVector[binNumber];
00093 }

void G4PhysicsVector::PutValue ( size_t  index,
G4double  theValue 
) [inline]

Definition at line 197 of file G4PhysicsVector.icc.

References dataVector.

Referenced by G4VXTRenergyLoss::BuildAngleForEnergyBank(), G4PolarizedCompton::BuildAsymmetryTable(), G4eplusPolarizedAnnihilation::BuildAsymmetryTable(), G4KokoulinMuonNuclearXS::BuildCrossSectionTable(), G4LossTableBuilder::BuildDEDXTable(), G4VXTRenergyLoss::BuildGlobalAngleTable(), G4VeLowEnergyLoss::BuildInverseRangeTable(), G4PAIPhotonModel::BuildLambdaVector(), G4PAIModel::BuildLambdaVector(), G4VRangeToEnergyConverter::BuildLossTable(), G4PAIPhotonModel::BuildPAIonisationTable(), G4PAIModel::BuildPAIonisationTable(), G4ePolarizedIonisation::BuildPhysicsTable(), G4ChargeExchangeProcess::BuildPhysicsTable(), G4VEnergyLoss::BuildRangeCoeffATable(), G4VeLowEnergyLoss::BuildRangeCoeffATable(), G4VEnergyLoss::BuildRangeCoeffBTable(), G4VeLowEnergyLoss::BuildRangeCoeffBTable(), G4VEnergyLoss::BuildRangeCoeffCTable(), G4VeLowEnergyLoss::BuildRangeCoeffCTable(), G4LossTableBuilder::BuildRangeTable(), G4VRangeToEnergyConverter::BuildRangeVector(), G4LossTableBuilder::BuildTableForModel(), G4AdjointCSManager::BuildTotalSigmaTables(), G4ForwardXrayTR::BuildXrayTRtables(), G4EmModelManager::FillDEDXVector(), G4EmModelManager::FillLambdaVector(), G4XNNElasticLowE::G4XNNElasticLowE(), G4XnpElasticLowE::G4XnpElasticLowE(), G4XnpTotalLowE::G4XnpTotalLowE(), G4eeToHadronsModel::Initialise(), G4InitXscPAI::IntegralCherenkov(), G4InitXscPAI::IntegralPAIdEdx(), G4InitXscPAI::IntegralPAIxSection(), G4InitXscPAI::IntegralPlasmon(), and G4VRangeToEnergyConverter::operator=().

00198 {
00199   dataVector[binNumber] = theValue;
00200 }

G4bool G4PhysicsVector::Retrieve ( std::ifstream &  fIn,
G4bool  ascii = false 
) [virtual]

Reimplemented in G4PhysicsLinearVector, G4PhysicsLnVector, and G4PhysicsLogVector.

Definition at line 198 of file G4PhysicsVector.cc.

References binVector, cache, dataVector, DBL_MAX, edgeMax, edgeMin, G4cerr, G4endl, G4PhysicsVectorCache::lastBin, G4PhysicsVectorCache::lastEnergy, G4PhysicsVectorCache::lastValue, numberOfNodes, and secDerivative.

Referenced by G4PhysicsLogVector::Retrieve(), G4PhysicsLnVector::Retrieve(), G4PhysicsLinearVector::Retrieve(), and G4PhysicsTable::RetrievePhysicsTable().

00199 {
00200   // clear properties;
00201   cache->lastEnergy=-DBL_MAX;
00202   cache->lastValue =0.;
00203   cache->lastBin   =0;
00204   dataVector.clear();
00205   binVector.clear();
00206   secDerivative.clear();
00207 
00208   // retrieve in ascii mode
00209   if (ascii){
00210     // binning
00211     fIn >> edgeMin >> edgeMax >> numberOfNodes; 
00212     if (fIn.fail())  { return false; }
00213     // contents
00214     G4int siz=0;
00215     fIn >> siz;
00216     if (fIn.fail())  { return false; }
00217     if (siz<=0)
00218     {
00219 #ifdef G4VERBOSE  
00220       G4cerr << "G4PhysicsVector::Retrieve():";
00221       G4cerr << " Invalid vector size: " << siz << G4endl;
00222 #endif
00223       return false;
00224     }
00225 
00226     binVector.reserve(siz);
00227     dataVector.reserve(siz);
00228     G4double vBin, vData;
00229 
00230     for(G4int i = 0; i < siz ; i++)
00231     {
00232       vBin = 0.;
00233       vData= 0.;
00234       fIn >> vBin >> vData;
00235       if (fIn.fail())  { return false; }
00236       binVector.push_back(vBin);
00237       dataVector.push_back(vData);
00238     }
00239 
00240     // to remove any inconsistency 
00241     numberOfNodes = siz;
00242     edgeMin = binVector[0];
00243     edgeMax = binVector[numberOfNodes-1];
00244     return true ;
00245   }
00246 
00247   // retrieve in binary mode
00248   // binning
00249   fIn.read((char*)(&edgeMin), sizeof edgeMin);
00250   fIn.read((char*)(&edgeMax), sizeof edgeMax);
00251   fIn.read((char*)(&numberOfNodes), sizeof numberOfNodes ); 
00252  
00253   // contents
00254   size_t size;
00255   fIn.read((char*)(&size), sizeof size); 
00256  
00257   G4double* value = new G4double[2*size];
00258   fIn.read((char*)(value),  2*size*(sizeof(G4double)) );
00259   if (G4int(fIn.gcount()) != G4int(2*size*(sizeof(G4double))) )
00260   {
00261     delete [] value;
00262     return false;
00263   }
00264 
00265   binVector.reserve(size);
00266   dataVector.reserve(size);
00267   for(size_t i = 0; i < size; ++i)
00268   {
00269     binVector.push_back(value[2*i]);
00270     dataVector.push_back(value[2*i+1]);
00271   }
00272   delete [] value;
00273 
00274   // to remove any inconsistency 
00275   numberOfNodes = size;
00276   edgeMin = binVector[0];
00277   edgeMax = binVector[numberOfNodes-1];
00278 
00279   return true;
00280 }

void G4PhysicsVector::ScaleVector ( G4double  factorE,
G4double  factorV 
) [virtual]

Reimplemented in G4PhysicsLinearVector, G4PhysicsLnVector, and G4PhysicsLogVector.

Definition at line 285 of file G4PhysicsVector.cc.

References binVector, cache, dataVector, edgeMax, edgeMin, G4PhysicsVectorCache::lastEnergy, G4PhysicsVectorCache::lastValue, CLHEP::detail::n, and secDerivative.

Referenced by G4PhysicsLogVector::ScaleVector(), G4PhysicsLnVector::ScaleVector(), and G4PhysicsLinearVector::ScaleVector().

00286 {
00287   size_t n = dataVector.size();
00288   size_t i;
00289   if(n > 0) { 
00290     for(i=0; i<n; ++i) {
00291       binVector[i]  *= factorE;
00292       dataVector[i] *= factorV;
00293     } 
00294   }
00295   //  n = secDerivative.size();
00296   // if(n > 0) { for(i=0; i<n; ++i) { secDerivative[i] *= factorV; } }
00297   secDerivative.clear();
00298 
00299   edgeMin *= factorE;
00300   edgeMax *= factorE;
00301   cache->lastEnergy = factorE*(cache->lastEnergy);
00302   cache->lastValue  = factorV*(cache->lastValue);
00303 }

void G4PhysicsVector::SetSpline ( G4bool   )  [inline]

Definition at line 224 of file G4PhysicsVector.icc.

Referenced by G4VEnergyLossProcess::BuildDEDXTable(), G4LossTableBuilder::BuildDEDXTable(), G4LossTableBuilder::BuildInverseRangeTable(), G4VEnergyLossProcess::BuildLambdaTable(), G4LossTableBuilder::BuildRangeTable(), G4LossTableBuilder::BuildTableForModel(), G4EmElementSelector::G4EmElementSelector(), G4HadronNucleonXsc::InitialiseKaonNucleonTotXsc(), G4VEnergyLossProcess::LambdaPhysicsVector(), G4VEmProcess::LambdaPhysicsVector(), G4eeToTwoPiModel::PhysicsVector(), G4eeToPGammaModel::PhysicsVector(), G4eeTo3PiModel::PhysicsVector(), G4ee2KNeutralModel::PhysicsVector(), and G4ee2KChargedModel::PhysicsVector().

00225 {
00226   useSpline = val;
00227 }

void G4PhysicsVector::SetVerboseLevel ( G4int  value  )  [inline]

Definition at line 232 of file G4PhysicsVector.icc.

References verboseLevel.

00233 {
00234    verboseLevel = value;
00235 }

G4bool G4PhysicsVector::Store ( std::ofstream &  fOut,
G4bool  ascii = false 
) [virtual]

Definition at line 165 of file G4PhysicsVector.cc.

References binVector, dataVector, edgeMax, edgeMin, and numberOfNodes.

00166 {
00167   // Ascii mode
00168   if (ascii)
00169   {
00170     fOut << *this;
00171     return true;
00172   } 
00173   // Binary Mode
00174 
00175   // binning
00176   fOut.write((char*)(&edgeMin), sizeof edgeMin);
00177   fOut.write((char*)(&edgeMax), sizeof edgeMax);
00178   fOut.write((char*)(&numberOfNodes), sizeof numberOfNodes);
00179 
00180   // contents
00181   size_t size = dataVector.size(); 
00182   fOut.write((char*)(&size), sizeof size);
00183 
00184   G4double* value = new G4double[2*size];
00185   for(size_t i = 0; i < size; ++i)
00186   {
00187     value[2*i]  =  binVector[i];
00188     value[2*i+1]=  dataVector[i];
00189   }
00190   fOut.write((char*)(value), 2*size*(sizeof (G4double)));
00191   delete [] value;
00192 
00193   return true;
00194 }

G4double G4PhysicsVector::Value ( G4double  theEnergy  )  [inline]

Definition at line 62 of file G4PhysicsVector.icc.

References cache, G4PhysicsVectorCache::lastEnergy, and G4PhysicsVectorCache::lastValue.

Referenced by G4EmCorrections::BarkasCorrection(), G4LossTableBuilder::BuildRangeTable(), G4Track::CalculateVelocityForOpticalPhoton(), G4PenelopeRayleighModel::ComputeCrossSectionPerAtom(), G4PenelopePhotoElectricModel::ComputeCrossSectionPerAtom(), G4PenelopeGammaConversionModel::ComputeCrossSectionPerAtom(), G4LivermoreRayleighModel::ComputeCrossSectionPerAtom(), G4LivermoreGammaConversionModel::ComputeCrossSectionPerAtom(), G4VRangeToEnergyConverter::ConvertCutToKineticEnergy(), G4EmCorrections::EffectiveChargeCorrection(), G4NeutronInelasticXS::GetElementCrossSection(), G4NeutronElasticXS::GetElementCrossSection(), G4NeutronCaptureXS::GetElementCrossSection(), G4PenelopeCrossSection::GetHardCrossSection(), G4NeutronCaptureXS::GetIsoCrossSection(), G4HadronNucleonXsc::GetKmNeutronTotXscVector(), G4HadronNucleonXsc::GetKmProtonTotXscVector(), G4HadronNucleonXsc::GetKpNeutronTotXscVector(), G4HadronNucleonXsc::GetKpProtonTotXscVector(), G4PenelopeCrossSection::GetNormalizedShellCrossSection(), G4SPSEneDistribution::GetProbability(), G4PenelopePhotoElectricModel::GetShellCrossSection(), G4PenelopeCrossSection::GetShellCrossSection(), G4PenelopeCrossSection::GetSoftStoppingPower(), G4PenelopeCrossSection::GetTotalCrossSection(), GetValue(), G4ElementData::GetValueForElement(), G4EmCorrections::KShellCorrection(), G4EmCorrections::LShellCorrection(), G4Scintillation::PostStepDoIt(), G4OpBoundaryProcess::PostStepDoIt(), G4Cerenkov::PostStepDoIt(), G4PenelopeBremsstrahlungAngular::SampleDirection(), G4PenelopeRayleighModel::SampleSecondaries(), G4NeutronCaptureXS::SelectIsotope(), G4VEnergyLossProcess::SetCSDARangeTable(), G4VEnergyLossProcess::SetDEDXTable(), and G4EmCorrections::ShellCorrection().

00063 {
00064   // Use cache for speed up - check if the value 'theEnergy' is same as the 
00065   // last call. If it is same, then use the last value, if not - recompute
00066 
00067   if( theEnergy != cache->lastEnergy ) { ComputeValue(theEnergy); }
00068   return cache->lastValue;
00069 }


Friends And Related Function Documentation

std::ostream& operator<< ( std::ostream &  out,
const G4PhysicsVector pv 
) [friend]

Definition at line 487 of file G4PhysicsVector.cc.

00488 {
00489   // binning
00490   out << std::setprecision(12) << pv.edgeMin << " "
00491       << pv.edgeMax << " " << pv.numberOfNodes << G4endl; 
00492 
00493   // contents
00494   out << pv.dataVector.size() << G4endl; 
00495   for(size_t i = 0; i < pv.dataVector.size(); i++)
00496   {
00497     out << pv.binVector[i] << "  " << pv.dataVector[i] << G4endl;
00498   }
00499   out << std::setprecision(6);
00500 
00501   return out;
00502 }


Field Documentation

G4double G4PhysicsVector::baseBin [protected]

Definition at line 232 of file G4PhysicsVector.hh.

Referenced by G4PhysicsLogVector::FindBinLocation(), G4PhysicsLnVector::FindBinLocation(), G4PhysicsLinearVector::FindBinLocation(), G4PhysicsLinearVector::G4PhysicsLinearVector(), G4PhysicsLnVector::G4PhysicsLnVector(), G4PhysicsLogVector::G4PhysicsLogVector(), G4PhysicsVector(), operator=(), G4PhysicsLogVector::Retrieve(), G4PhysicsLnVector::Retrieve(), G4PhysicsLinearVector::Retrieve(), G4PhysicsLogVector::ScaleVector(), G4PhysicsLnVector::ScaleVector(), and G4PhysicsLinearVector::ScaleVector().

G4PVDataVector G4PhysicsVector::binVector [protected]

Definition at line 212 of file G4PhysicsVector.hh.

Referenced by ComputeSecDerivatives(), ComputeSecondDerivatives(), CopyData(), G4PhysicsOrderedFreeVector::DumpValues(), G4LPhysicsFreeVector::DumpValues(), Energy(), FillSecondDerivatives(), G4PhysicsFreeVector::FindBinLocation(), G4LPhysicsFreeVector::G4LPhysicsFreeVector(), G4PhysicsFreeVector::G4PhysicsFreeVector(), G4PhysicsLinearVector::G4PhysicsLinearVector(), G4PhysicsLnVector::G4PhysicsLnVector(), G4PhysicsLogVector::G4PhysicsLogVector(), GetLowEdgeEnergy(), G4PhysicsOrderedFreeVector::GetLowEdgeEnergy(), G4PhysicsOrderedFreeVector::GetMaxLowEdgeEnergy(), G4PhysicsOrderedFreeVector::GetMinLowEdgeEnergy(), G4PhysicsOrderedFreeVector::InsertValues(), operator<<(), G4PhysicsFreeVector::PutValue(), G4LPhysicsFreeVector::PutValues(), Retrieve(), G4PhysicsLogVector::Retrieve(), G4PhysicsLnVector::Retrieve(), G4PhysicsLinearVector::Retrieve(), ScaleVector(), G4PhysicsLogVector::ScaleVector(), G4PhysicsLnVector::ScaleVector(), G4PhysicsLinearVector::ScaleVector(), and Store().

G4PhysicsVectorCache* G4PhysicsVector::cache [protected]

Definition at line 209 of file G4PhysicsVector.hh.

Referenced by CopyData(), G4PhysicsVector(), GetLastBin(), GetLastEnergy(), GetLastValue(), Retrieve(), ScaleVector(), Value(), and ~G4PhysicsVector().

G4PVDataVector G4PhysicsVector::dataVector [protected]

Definition at line 211 of file G4PhysicsVector.hh.

Referenced by ComputeSecDerivatives(), ComputeSecondDerivatives(), CopyData(), G4PhysicsOrderedFreeVector::DumpValues(), G4LPhysicsFreeVector::DumpValues(), FillSecondDerivatives(), G4LPhysicsFreeVector::G4LPhysicsFreeVector(), G4PhysicsFreeVector::G4PhysicsFreeVector(), G4PhysicsLinearVector::G4PhysicsLinearVector(), G4PhysicsLnVector::G4PhysicsLnVector(), G4PhysicsLogVector::G4PhysicsLogVector(), G4PhysicsOrderedFreeVector::GetMaxValue(), G4PhysicsOrderedFreeVector::GetMinValue(), G4PhysicsOrderedFreeVector::InsertValues(), operator()(), operator<<(), operator[](), PutValue(), G4PhysicsFreeVector::PutValue(), G4LPhysicsFreeVector::PutValues(), Retrieve(), ScaleVector(), and Store().

G4double G4PhysicsVector::dBin [protected]

Definition at line 231 of file G4PhysicsVector.hh.

Referenced by G4PhysicsLogVector::FindBinLocation(), G4PhysicsLnVector::FindBinLocation(), G4PhysicsLinearVector::FindBinLocation(), G4PhysicsLinearVector::G4PhysicsLinearVector(), G4PhysicsLnVector::G4PhysicsLnVector(), G4PhysicsLogVector::G4PhysicsLogVector(), G4PhysicsVector(), operator=(), G4PhysicsLogVector::Retrieve(), G4PhysicsLnVector::Retrieve(), G4PhysicsLinearVector::Retrieve(), G4PhysicsLogVector::ScaleVector(), G4PhysicsLnVector::ScaleVector(), and G4PhysicsLinearVector::ScaleVector().

G4double G4PhysicsVector::edgeMax [protected]

Definition at line 205 of file G4PhysicsVector.hh.

Referenced by CopyData(), G4LPhysicsFreeVector::G4LPhysicsFreeVector(), G4PhysicsFreeVector::G4PhysicsFreeVector(), G4PhysicsLinearVector::G4PhysicsLinearVector(), G4PhysicsLnVector::G4PhysicsLnVector(), G4PhysicsLogVector::G4PhysicsLogVector(), GetMaxEnergy(), G4PhysicsOrderedFreeVector::InsertValues(), operator<<(), G4PhysicsFreeVector::PutValue(), G4LPhysicsFreeVector::PutValues(), Retrieve(), ScaleVector(), and Store().

G4double G4PhysicsVector::edgeMin [protected]

Definition at line 204 of file G4PhysicsVector.hh.

Referenced by CopyData(), G4LPhysicsFreeVector::G4LPhysicsFreeVector(), G4PhysicsFreeVector::G4PhysicsFreeVector(), G4PhysicsLinearVector::G4PhysicsLinearVector(), G4PhysicsLnVector::G4PhysicsLnVector(), G4PhysicsLogVector::G4PhysicsLogVector(), G4PhysicsOrderedFreeVector::InsertValues(), operator<<(), G4PhysicsFreeVector::PutValue(), G4LPhysicsFreeVector::PutValues(), Retrieve(), ScaleVector(), and Store().

size_t G4PhysicsVector::numberOfNodes [protected]

Definition at line 207 of file G4PhysicsVector.hh.

Referenced by ComputeSecDerivatives(), ComputeSecondDerivatives(), CopyData(), G4PhysicsOrderedFreeVector::DumpValues(), G4LPhysicsFreeVector::DumpValues(), FillSecondDerivatives(), G4PhysicsFreeVector::FindBinLocation(), G4LPhysicsFreeVector::G4LPhysicsFreeVector(), G4PhysicsFreeVector::G4PhysicsFreeVector(), G4PhysicsLinearVector::G4PhysicsLinearVector(), G4PhysicsLnVector::G4PhysicsLnVector(), G4PhysicsLogVector::G4PhysicsLogVector(), GetVectorLength(), G4PhysicsOrderedFreeVector::InsertValues(), IsFilledVectorExist(), operator<<(), G4PhysicsFreeVector::PutValue(), G4LPhysicsFreeVector::PutValues(), Retrieve(), and Store().

G4PVDataVector G4PhysicsVector::secDerivative [protected]

Definition at line 213 of file G4PhysicsVector.hh.

Referenced by ComputeSecDerivatives(), ComputeSecondDerivatives(), CopyData(), DeleteData(), FillSecondDerivatives(), Retrieve(), and ScaleVector().

G4PhysicsVectorType G4PhysicsVector::type [protected]

Definition at line 202 of file G4PhysicsVector.hh.

Referenced by CopyData(), G4LPhysicsFreeVector::G4LPhysicsFreeVector(), G4PhysicsFreeVector::G4PhysicsFreeVector(), G4PhysicsLinearVector::G4PhysicsLinearVector(), G4PhysicsLnVector::G4PhysicsLnVector(), G4PhysicsLogVector::G4PhysicsLogVector(), G4PhysicsOrderedFreeVector::G4PhysicsOrderedFreeVector(), and GetType().

G4int G4PhysicsVector::verboseLevel [protected]

Definition at line 234 of file G4PhysicsVector.hh.

Referenced by G4PhysicsVector(), GetVerboseLevel(), operator=(), and SetVerboseLevel().


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