Geant4-11
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes | Friends
G4PhysicsVector Class Reference

#include <G4PhysicsVector.hh>

Inheritance diagram for G4PhysicsVector:
G4PhysicsFreeVector G4PhysicsLinearVector G4PhysicsLogVector

Public Member Functions

std::size_t ComputeLogVectorBin (const G4double logenergy) const
 
void DumpValues (G4double unitE=1.0, G4double unitV=1.0) const
 
G4double Energy (const std::size_t index) const
 
void FillSecondDerivatives (const G4SplineType=G4SplineType::Base, const G4double dir1=0.0, const G4double dir2=0.0)
 
std::size_t FindBin (const G4double energy, std::size_t idx) const
 
G4double FindLinearEnergy (const G4double rand) const
 
 G4PhysicsVector (const G4PhysicsVector &&)=delete
 
 G4PhysicsVector (const G4PhysicsVector &)=default
 
 G4PhysicsVector (G4bool spline=false)
 
G4double GetEnergy (const G4double value) const
 
G4double GetLowEdgeEnergy (const std::size_t index) const
 
G4double GetMaxEnergy () const
 
G4double GetMaxValue () const
 
G4double GetMinEnergy () const
 
G4double GetMinValue () const
 
G4bool GetSpline () const
 
G4PhysicsVectorType GetType () const
 
G4double GetValue (const G4double energy, G4bool &isOutRange) const
 
std::size_t GetVectorLength () const
 
G4double LogVectorValue (const G4double energy, const G4double theLogEnergy) const
 
G4bool operator!= (const G4PhysicsVector &right) const =delete
 
G4double operator() (const std::size_t index) const
 
G4PhysicsVectoroperator= (const G4PhysicsVector &&)=delete
 
G4PhysicsVectoroperator= (const G4PhysicsVector &)=default
 
G4bool operator== (const G4PhysicsVector &right) const =delete
 
G4double operator[] (const std::size_t index) const
 
void PutValue (const std::size_t index, const G4double value)
 
G4bool Retrieve (std::ifstream &fIn, G4bool ascii=false)
 
void ScaleVector (const G4double factorE, const G4double factorV)
 
void SetVerboseLevel (G4int value)
 
G4bool Store (std::ofstream &fOut, G4bool ascii=false) const
 
G4double Value (const G4double energy) const
 
G4double Value (const G4double energy, std::size_t &lastidx) const
 
virtual ~G4PhysicsVector ()=default
 

Protected Member Functions

virtual void Initialise ()
 
void PrintPutValueError (std::size_t index, G4double value, const G4String &text)
 

Protected Attributes

std::vector< G4doublebinVector
 
std::vector< G4doubledataVector
 
G4double edgeMax = 0.0
 
G4double edgeMin = 0.0
 
G4int idxmax = 0
 
G4double invdBin = 0.0
 
G4double logemin = 0.0
 
std::size_t numberOfNodes = 0
 
std::vector< G4doublesecDerivative
 
G4PhysicsVectorType type = T_G4PhysicsFreeVector
 
G4int verboseLevel = 0
 

Private Member Functions

void ComputeSecDerivative0 ()
 
void ComputeSecDerivative1 ()
 
void ComputeSecDerivative2 (const G4double firstPointDerivative, const G4double endPointDerivative)
 
std::size_t GetBin (const G4double energy) const
 
G4double Interpolation (const std::size_t idx, const G4double energy) const
 

Private Attributes

G4bool useSpline = false
 

Friends

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

Detailed Description

Definition at line 54 of file G4PhysicsVector.hh.

Constructor & Destructor Documentation

◆ G4PhysicsVector() [1/3]

G4PhysicsVector::G4PhysicsVector ( G4bool  spline = false)
explicit

Definition at line 39 of file G4PhysicsVector.cc.

40 : useSpline(val)
41{}

◆ G4PhysicsVector() [2/3]

G4PhysicsVector::G4PhysicsVector ( const G4PhysicsVector )
default

◆ G4PhysicsVector() [3/3]

G4PhysicsVector::G4PhysicsVector ( const G4PhysicsVector &&  )
delete

◆ ~G4PhysicsVector()

virtual G4PhysicsVector::~G4PhysicsVector ( )
virtualdefault

Member Function Documentation

◆ ComputeLogVectorBin()

std::size_t G4PhysicsVector::ComputeLogVectorBin ( const G4double  logenergy) const
inline

◆ ComputeSecDerivative0()

void G4PhysicsVector::ComputeSecDerivative0 ( )
private

Definition at line 272 of file G4PhysicsVector.cc.

274{
275 std::size_t n = numberOfNodes - 1;
276
277 for(std::size_t i = 1; i < n; ++i)
278 {
279 secDerivative[i] = 3.0 *
280 ((dataVector[i + 1] - dataVector[i]) / (binVector[i + 1] - binVector[i]) -
281 (dataVector[i] - dataVector[i - 1]) /
282 (binVector[i] - binVector[i - 1])) /
283 (binVector[i + 1] - binVector[i - 1]);
284 }
287}
std::size_t numberOfNodes
std::vector< G4double > secDerivative
std::vector< G4double > dataVector
std::vector< G4double > binVector

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

Referenced by FillSecondDerivatives().

◆ ComputeSecDerivative1()

void G4PhysicsVector::ComputeSecDerivative1 ( )
private

Definition at line 290 of file G4PhysicsVector.cc.

294{
295 G4int n = numberOfNodes - 1;
296 G4double* u = new G4double[n];
297 G4double p, sig;
298
299 u[1] = ((dataVector[2] - dataVector[1]) / (binVector[2] - binVector[1]) -
300 (dataVector[1] - dataVector[0]) / (binVector[1] - binVector[0]));
301 u[1] = 6.0 * u[1] * (binVector[2] - binVector[1]) /
302 ((binVector[2] - binVector[0]) * (binVector[2] - binVector[0]));
303
304 // Decomposition loop for tridiagonal algorithm. secDerivative[i]
305 // and u[i] are used for temporary storage of the decomposed factors.
306
307 secDerivative[1] = (2.0 * binVector[1] - binVector[0] - binVector[2]) /
308 (2.0 * binVector[2] - binVector[0] - binVector[1]);
309
310 for(G4int i = 2; i < n - 1; ++i)
311 {
312 sig =
313 (binVector[i] - binVector[i - 1]) / (binVector[i + 1] - binVector[i - 1]);
314 p = sig * secDerivative[i - 1] + 2.0;
315 secDerivative[i] = (sig - 1.0) / p;
316 u[i] =
317 (dataVector[i + 1] - dataVector[i]) / (binVector[i + 1] - binVector[i]) -
318 (dataVector[i] - dataVector[i - 1]) / (binVector[i] - binVector[i - 1]);
319 u[i] =
320 (6.0 * u[i] / (binVector[i + 1] - binVector[i - 1])) - sig * u[i - 1] / p;
321 }
322
323 sig =
324 (binVector[n - 1] - binVector[n - 2]) / (binVector[n] - binVector[n - 2]);
325 p = sig * secDerivative[n - 3] + 2.0;
326 u[n - 1] =
327 (dataVector[n] - dataVector[n - 1]) / (binVector[n] - binVector[n - 1]) -
328 (dataVector[n - 1] - dataVector[n - 2]) /
329 (binVector[n - 1] - binVector[n - 2]);
330 u[n - 1] = 6.0 * sig * u[n - 1] / (binVector[n] - binVector[n - 2]) -
331 (2.0 * sig - 1.0) * u[n - 2] / p;
332
333 p = (1.0 + sig) + (2.0 * sig - 1.0) * secDerivative[n - 2];
334 secDerivative[n - 1] = u[n - 1] / p;
335
336 // The back-substitution loop for the triagonal algorithm of solving
337 // a linear system of equations.
338
339 for(G4int k = n - 2; k > 1; --k)
340 {
341 secDerivative[k] *=
342 (secDerivative[k + 1] - u[k] * (binVector[k + 1] - binVector[k - 1]) /
343 (binVector[k + 1] - binVector[k]));
344 }
346 (secDerivative[n - 1] - (1.0 - sig) * secDerivative[n - 2]) / sig;
347 sig = 1.0 - ((binVector[2] - binVector[1]) / (binVector[2] - binVector[0]));
348 secDerivative[1] *= (secDerivative[2] - u[1] / (1.0 - sig));
349 secDerivative[0] = (secDerivative[1] - sig * secDerivative[2]) / (1.0 - sig);
350
351 delete[] u;
352}
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85

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

Referenced by FillSecondDerivatives().

◆ ComputeSecDerivative2()

void G4PhysicsVector::ComputeSecDerivative2 ( const G4double  firstPointDerivative,
const G4double  endPointDerivative 
)
private

Definition at line 355 of file G4PhysicsVector.cc.

361{
362 G4int n = numberOfNodes - 1;
363 G4double* u = new G4double[n];
364 G4double p, sig, un;
365
366 u[0] = (6.0 / (binVector[1] - binVector[0])) *
367 ((dataVector[1] - dataVector[0]) / (binVector[1] - binVector[0]) -
368 firstPointDerivative);
369
370 secDerivative[0] = -0.5;
371
372 // Decomposition loop for tridiagonal algorithm. secDerivative[i]
373 // and u[i] are used for temporary storage of the decomposed factors.
374
375 for(G4int i = 1; i < n; ++i)
376 {
377 sig =
378 (binVector[i] - binVector[i - 1]) / (binVector[i + 1] - binVector[i - 1]);
379 p = sig * (secDerivative[i - 1]) + 2.0;
380 secDerivative[i] = (sig - 1.0) / p;
381 u[i] =
382 (dataVector[i + 1] - dataVector[i]) / (binVector[i + 1] - binVector[i]) -
383 (dataVector[i] - dataVector[i - 1]) / (binVector[i] - binVector[i - 1]);
384 u[i] =
385 6.0 * u[i] / (binVector[i + 1] - binVector[i - 1]) - sig * u[i - 1] / p;
386 }
387
388 sig =
389 (binVector[n - 1] - binVector[n - 2]) / (binVector[n] - binVector[n - 2]);
390 p = sig * secDerivative[n - 2] + 2.0;
391 un = (6.0 / (binVector[n] - binVector[n - 1])) *
392 (endPointDerivative - (dataVector[n] - dataVector[n - 1]) /
393 (binVector[n] - binVector[n - 1])) -
394 u[n - 1] / p;
395 secDerivative[n] = un / (secDerivative[n - 1] + 2.0);
396
397 // The back-substitution loop for the triagonal algorithm of solving
398 // a linear system of equations.
399
400 for(G4int k = n - 1; k > 0; --k)
401 {
402 secDerivative[k] *=
403 (secDerivative[k + 1] - u[k] * (binVector[k + 1] - binVector[k - 1]) /
404 (binVector[k + 1] - binVector[k]));
405 }
406 secDerivative[0] = 0.5 * (u[0] - secDerivative[1]);
407
408 delete[] u;
409}

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

Referenced by FillSecondDerivatives().

◆ DumpValues()

void G4PhysicsVector::DumpValues ( G4double  unitE = 1.0,
G4double  unitV = 1.0 
) const

Definition at line 163 of file G4PhysicsVector.cc.

164{
165 for(std::size_t i = 0; i < numberOfNodes; ++i)
166 {
167 G4cout << binVector[i] / unitE << " " << dataVector[i] / unitV
168 << G4endl;
169 }
170}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout

References binVector, dataVector, G4cout, G4endl, and numberOfNodes.

Referenced by G4OpWLS::DumpPhysicsTable(), G4OpWLS2::DumpPhysicsTable(), FillSecondDerivatives(), G4SPSEneDistribution::GenUserHistEnergies(), and G4SPSEneDistribution::LinearInterpolation().

◆ Energy()

G4double G4PhysicsVector::Energy ( const std::size_t  index) const
inline

Referenced by G4SPSEneDistribution::ArbInterpolate(), G4EmCorrections::BarkasCorrection(), G4DiffuseElasticV2::BuildAngleTable(), G4PolarizedCompton::BuildAsymmetryTable(), G4PolarizedAnnihilation::BuildAsymmetryTables(), G4PolarizedIonisation::BuildAsymmetryTables(), G4EmCorrections::BuildCorrectionVector(), G4KokoulinMuonNuclearXS::BuildCrossSectionTable(), G4LossTableBuilder::BuildInverseRangeTable(), G4VLEPTSModel::BuildMeanFreePathTable(), G4Cerenkov::BuildPhysicsTable(), G4Scintillation::BuildPhysicsTable(), G4OpWLS::BuildPhysicsTable(), G4OpWLS2::BuildPhysicsTable(), G4LossTableBuilder::BuildRangeTable(), G4PenelopeBremsstrahlungFS::BuildScaledXSTable(), G4LossTableBuilder::BuildTableForModel(), G4AdjointCSManager::BuildTotalSigmaTables(), G4MaterialPropertiesTable::CalculateGROUPVEL(), G4OpRayleigh::CalculateRayleighMeanFreePaths(), G4PenelopeRayleighModelMI::CalculateThetaAndAngFun(), G4eeToHadronsModel::ComputeCMCrossSectionPerElectron(), G4JAEAElasticScatteringModel::ComputeCrossSectionPerAtom(), G4JAEAPolarizedElasticScatteringModel::ComputeCrossSectionPerAtom(), G4LivermoreComptonModel::ComputeCrossSectionPerAtom(), G4LivermorePolarizedComptonModel::ComputeCrossSectionPerAtom(), G4LivermorePolarizedRayleighModel::ComputeCrossSectionPerAtom(), G4LivermoreRayleighModel::ComputeCrossSectionPerAtom(), G4LowEPComptonModel::ComputeCrossSectionPerAtom(), G4LowEPPolarizedComptonModel::ComputeCrossSectionPerAtom(), G4TablesForExtrapolator::ComputeElectronDEDX(), G4TablesForExtrapolator::ComputeMuonDEDX(), G4TablesForExtrapolator::ComputeProtonDEDX(), G4TablesForExtrapolator::ComputeTrasportXS(), G4PenelopeRayleighModelMI::CrossSectionPerVolume(), G4PAIModelData::CrossSectionPerVolume(), G4PAIPhotData::CrossSectionPerVolume(), G4PAIModelData::DEDXPerVolume(), G4PAIPhotData::DEDXPerVolume(), G4EmModelManager::DumpModelList(), G4EmModelManager::FillDEDXVector(), G4EmModelManager::FillLambdaVector(), G4HadElementSelector::G4HadElementSelector(), G4LivermorePhotoElectricModel::GetBindingEnergy(), G4ICRU90StoppingData::GetDEDX(), G4PAIPhotData::GetEnergyPhotonTransfer(), G4PAIPhotData::GetEnergyPlasmonTransfer(), G4PAIModelData::GetEnergyTransfer(), G4PAIPhotData::GetEnergyTransfer(), G4PAIPhotData::GetPlasmonRatio(), G4PAIPhotData::Initialise(), G4PAIModelData::Initialise(), G4eeToHadronsModel::Initialise(), G4eplusTo2GammaOKVIModel::Initialise(), G4WentzelVIModel::Initialise(), G4DNABornExcitationModel2::Initialise(), G4HadronXSDataTable::Initialise(), G4NeutronCaptureXS::IsoCrossSection(), G4Cerenkov::PostStepDoIt(), G4GDMLWriteMaterials::PropertyVectorWrite(), G4IonICRU73Data::RetrieveVector(), G4PAIPhotData::SampleAlongStepPhotonTransfer(), G4PAIPhotData::SampleAlongStepPlasmonTransfer(), G4PAIPhotData::SampleAlongStepTransfer(), G4PAIModelData::SampleAlongStepTransfer(), G4PAIPhotData::SamplePostStepPhotonTransfer(), G4PAIPhotData::SamplePostStepPlasmonTransfer(), G4PAIPhotData::SamplePostStepTransfer(), G4PAIModelData::SamplePostStepTransfer(), G4DiffuseElasticV2::SampleTableThetaCMS(), G4VEnergyLossProcess::ScaledKinEnergyForLoss(), G4VEnergyLossProcess::SetLambdaTable(), and G4VEmProcess::StreamInfo().

◆ FillSecondDerivatives()

void G4PhysicsVector::FillSecondDerivatives ( const G4SplineType  stype = G4SplineType::Base,
const G4double  dir1 = 0.0,
const G4double  dir2 = 0.0 
)

Definition at line 205 of file G4PhysicsVector.cc.

208{
209 if(!useSpline) { return; }
210 // cannot compute derivatives for less than 5 points
211 const std::size_t nmin = (stype == G4SplineType::Base) ? 5 : 4;
212 if(nmin > numberOfNodes)
213 {
214 if(0 < verboseLevel)
215 {
216 G4cout << "### G4PhysicsVector: spline cannot be used for "
217 << numberOfNodes << " points - spline disabled"
218 << G4endl;
219 DumpValues();
220 }
221 useSpline = false;
222 return;
223 }
224 // check energies of free vector
226 {
227 for(G4int i=0; i<=idxmax; ++i)
228 {
229 if(binVector[i + 1] <= binVector[i])
230 {
231 if(0 < verboseLevel)
232 {
233 G4cout << "### G4PhysicsVector: spline cannot be used, because "
234 << " E[" << i << "]=" << binVector[i]
235 << " >= E[" << i+1 << "]=" << binVector[i + 1]
236 << G4endl;
237 DumpValues();
238 }
239 useSpline = false;
240 return;
241 }
242 }
243 }
244
245 // spline is possible
246 Initialise();
248
249 if(1 < verboseLevel)
250 {
251 G4cout << "### G4PhysicsVector:: FillSecondDerivatives N="
252 << numberOfNodes << G4endl;
253 DumpValues();
254 }
255
256 switch(stype)
257 {
260 break;
261
263 ComputeSecDerivative2(dir1, dir2);
264 break;
265
266 default:
268 }
269}
@ T_G4PhysicsFreeVector
G4PhysicsVectorType type
void ComputeSecDerivative0()
void ComputeSecDerivative2(const G4double firstPointDerivative, const G4double endPointDerivative)
void ComputeSecDerivative1()
virtual void Initialise()
void DumpValues(G4double unitE=1.0, G4double unitV=1.0) const

References Base, binVector, ComputeSecDerivative0(), ComputeSecDerivative1(), ComputeSecDerivative2(), DumpValues(), FixedEdges, G4cout, G4endl, idxmax, Initialise(), numberOfNodes, secDerivative, T_G4PhysicsFreeVector, type, useSpline, and verboseLevel.

Referenced by G4WaterStopping::AddData(), G4PSTARStopping::AddData(), G4ASTARStopping::AddData(), G4ICRU90StoppingData::AddData(), G4UPiNuclearCrossSection::AddDataSet(), G4MaterialPropertiesTable::AddProperty(), G4PolarizedCompton::BuildAsymmetryTable(), G4PolarizedAnnihilation::BuildAsymmetryTables(), G4VEnergyLossProcess::BuildDEDXTable(), G4LossTableBuilder::BuildDEDXTable(), G4PenelopeRayleighModel::BuildFormFactorTable(), G4PenelopeRayleighModelMI::BuildFormFactorTable(), G4LossTableBuilder::BuildInverseRangeTable(), G4VEmProcess::BuildLambdaTable(), G4VEnergyLossProcess::BuildLambdaTable(), G4LossTableBuilder::BuildRangeTable(), G4LossTableBuilder::BuildTableForModel(), G4TablesForExtrapolator::ComputeElectronDEDX(), G4TablesForExtrapolator::ComputeMuonDEDX(), G4TablesForExtrapolator::ComputeProtonDEDX(), G4TablesForExtrapolator::ComputeTrasportXS(), G4eplusTo2GammaOKVIModel::Initialise(), G4WentzelVIModel::Initialise(), G4HadronXSDataTable::Initialise(), G4LindhardSorensenData::InitialiseData(), G4PenelopeBremsstrahlungAngular::PrepareTables(), G4ComponentSAIDTotalXS::ReadData(), G4LivermorePhotoElectricModel::ReadData(), G4JAEAElasticScatteringModel::ReadData(), G4JAEAPolarizedElasticScatteringModel::ReadData(), G4LivermoreNuclearGammaConversionModel::ReadData(), G4LivermorePolarizedGammaConversionModel::ReadData(), G4BoldyshevTripletModel::ReadData(), G4IonICRU73Data::ReadElementData(), and G4IonICRU73Data::RetrieveVector().

◆ FindBin()

std::size_t G4PhysicsVector::FindBin ( const G4double  energy,
std::size_t  idx 
) const

◆ FindLinearEnergy()

G4double G4PhysicsVector::FindLinearEnergy ( const G4double  rand) const
inline

◆ GetBin()

std::size_t G4PhysicsVector::GetBin ( const G4double  energy) const
inlineprivate

Referenced by FindBin().

◆ GetEnergy()

G4double G4PhysicsVector::GetEnergy ( const G4double  value) const

Definition at line 431 of file G4PhysicsVector.cc.

432{
433 if(0 == numberOfNodes)
434 {
435 return 0.0;
436 }
437 if(1 == numberOfNodes || val <= dataVector[0])
438 {
439 return edgeMin;
440 }
441 if(val >= dataVector[numberOfNodes - 1])
442 {
443 return edgeMax;
444 }
445 std::size_t bin =
446 std::lower_bound(dataVector.begin(), dataVector.end(), val) -
447 dataVector.begin() - 1;
448 if(static_cast<G4int>(bin) > idxmax) { bin = idxmax; }
449 G4double res = binVector[bin];
450 G4double del = dataVector[bin + 1] - dataVector[bin];
451 if(del > 0.0)
452 {
453 res += (val - dataVector[bin]) * (binVector[bin + 1] - res) / del;
454 }
455 return res;
456}

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

Referenced by G4SPSEneDistribution::GenEpnHistEnergies(), G4SPSAngDistribution::GenerateUserDefPhi(), G4SPSAngDistribution::GenerateUserDefTheta(), G4SPSRandomGenerator::GenRandEnergy(), G4SPSRandomGenerator::GenRandPhi(), G4SPSRandomGenerator::GenRandPosPhi(), G4SPSRandomGenerator::GenRandPosTheta(), G4SPSRandomGenerator::GenRandTheta(), G4SPSRandomGenerator::GenRandX(), G4SPSRandomGenerator::GenRandY(), G4SPSRandomGenerator::GenRandZ(), G4SPSEneDistribution::GenUserHistEnergies(), G4Scintillation::PostStepDoIt(), G4OpWLS::PostStepDoIt(), and G4OpWLS2::PostStepDoIt().

◆ GetLowEdgeEnergy()

G4double G4PhysicsVector::GetLowEdgeEnergy ( const std::size_t  index) const
inline

Referenced by G4EMDissociation::ApplyYourself(), G4VXTRenergyLoss::BuildAngleForEnergyBank(), G4VXTRenergyLoss::BuildAngleTable(), G4DiffuseElastic::BuildAngleTable(), G4NuclNuclDiffuseElastic::BuildAngleTable(), G4hRDEnergyLoss::BuildDEDXTable(), G4PenelopeIonisationXSHandler::BuildDeltaTable(), G4VXTRenergyLoss::BuildEnergyTable(), G4VXTRenergyLoss::BuildGlobalAngleTable(), G4hRDEnergyLoss::BuildInverseRangeTable(), G4hRDEnergyLoss::BuildLabTimeVector(), G4hImpactIonisation::BuildLambdaTable(), G4hImpactIonisation::BuildLossTable(), G4hRDEnergyLoss::BuildProperTimeVector(), G4hRDEnergyLoss::BuildRangeVector(), G4hhElastic::BuildTableT(), G4ForwardXrayTR::BuildXrayTRtables(), G4PenelopeIonisationXSHandler::BuildXSTable(), G4PenelopeBremsstrahlungModel::BuildXSTable(), G4SPSEneDistribution::ConvertEPNToEnergy(), G4PenelopeRayleighModel::DumpFormFactorTable(), G4PenelopeRayleighModelMI::DumpFormFactorTable(), G4SPSEneDistribution::ExpInterpolation(), G4SPSEneDistribution::GenArbPointEnergies(), G4SPSEneDistribution::GenEpnHistEnergies(), G4SPSAngDistribution::GenerateUserDefPhi(), G4SPSAngDistribution::GenerateUserDefTheta(), G4SPSRandomGenerator::GenRandEnergy(), G4SPSRandomGenerator::GenRandPhi(), G4SPSRandomGenerator::GenRandPosPhi(), G4SPSRandomGenerator::GenRandPosTheta(), G4SPSRandomGenerator::GenRandTheta(), G4SPSRandomGenerator::GenRandX(), G4SPSRandomGenerator::GenRandY(), G4SPSRandomGenerator::GenRandZ(), G4SPSEneDistribution::GenUserHistEnergies(), G4ForwardXrayTR::GetEnergyTR(), G4VXTRenergyLoss::GetMeanFreePath(), G4VXTRenergyLoss::GetNumberOfPhotons(), G4VXTRenergyLoss::GetXTRrandomEnergy(), G4NeutronElectronElXsc::Initialise(), G4NeutronElectronElModel::Initialise(), G4InitXscPAI::IntegralCherenkov(), G4InitXscPAI::IntegralPAIdEdx(), G4InitXscPAI::IntegralPAIxSection(), G4InitXscPAI::IntegralPlasmon(), G4hRDEnergyLoss::InvertRangeVector(), G4SPSEneDistribution::LinearInterpolation(), G4SPSEneDistribution::LogInterpolation(), G4PenelopeCrossSection::NormalizeShellCrossSections(), G4ForwardXrayTR::PostStepDoIt(), G4XNNElasticLowE::Print(), G4XnpElasticLowE::Print(), G4XnpTotalLowE::Print(), G4DiffuseElastic::SampleTableThetaCMS(), G4NuclNuclDiffuseElastic::SampleTableThetaCMS(), and G4SPSEneDistribution::SplineInterpolation().

◆ GetMaxEnergy()

G4double G4PhysicsVector::GetMaxEnergy ( ) const
inline

◆ GetMaxValue()

G4double G4PhysicsVector::GetMaxValue ( ) const
inline

◆ GetMinEnergy()

G4double G4PhysicsVector::GetMinEnergy ( ) const
inline

◆ GetMinValue()

G4double G4PhysicsVector::GetMinValue ( ) const
inline

◆ GetSpline()

G4bool G4PhysicsVector::GetSpline ( ) const
inline

◆ GetType()

G4PhysicsVectorType G4PhysicsVector::GetType ( ) const
inline

◆ GetValue()

G4double G4PhysicsVector::GetValue ( const G4double  energy,
G4bool isOutRange 
) const
inline

◆ GetVectorLength()

std::size_t G4PhysicsVector::GetVectorLength ( ) const
inline

Referenced by G4PolarizedAnnihilation::BuildAsymmetryTables(), G4PolarizedIonisation::BuildAsymmetryTables(), G4LossTableBuilder::BuildDEDXTable(), G4PenelopeIonisationXSHandler::BuildDeltaTable(), G4hRDEnergyLoss::BuildInverseRangeTable(), G4LossTableBuilder::BuildInverseRangeTable(), G4Cerenkov::BuildPhysicsTable(), G4Scintillation::BuildPhysicsTable(), G4OpWLS::BuildPhysicsTable(), G4OpWLS2::BuildPhysicsTable(), G4LossTableBuilder::BuildRangeTable(), G4PenelopeIonisationXSHandler::BuildXSTable(), G4PenelopeBremsstrahlungModel::BuildXSTable(), G4MaterialPropertiesTable::CalculateGROUPVEL(), G4OpRayleigh::CalculateRayleighMeanFreePaths(), G4JAEAElasticScatteringModel::ComputeCrossSectionPerAtom(), G4JAEAPolarizedElasticScatteringModel::ComputeCrossSectionPerAtom(), G4LivermoreComptonModel::ComputeCrossSectionPerAtom(), G4LivermoreNuclearGammaConversionModel::ComputeCrossSectionPerAtom(), G4LivermorePolarizedComptonModel::ComputeCrossSectionPerAtom(), G4LivermorePolarizedGammaConversionModel::ComputeCrossSectionPerAtom(), G4LivermorePolarizedRayleighModel::ComputeCrossSectionPerAtom(), G4LivermoreRayleighModel::ComputeCrossSectionPerAtom(), G4LowEPComptonModel::ComputeCrossSectionPerAtom(), G4LowEPPolarizedComptonModel::ComputeCrossSectionPerAtom(), G4SPSEneDistribution::ConvertEPNToEnergy(), G4PAIModelData::CrossSectionPerVolume(), G4PAIPhotData::CrossSectionPerVolume(), G4PAIModelData::DEDXPerVolume(), G4PAIPhotData::DEDXPerVolume(), G4PenelopeRayleighModel::DumpFormFactorTable(), G4PenelopeRayleighModelMI::DumpFormFactorTable(), G4EmModelManager::DumpModelList(), G4SPSEneDistribution::ExpInterpolation(), G4EmModelManager::FillDEDXVector(), G4EmModelManager::FillLambdaVector(), G4SPSEneDistribution::GenArbPointEnergies(), G4SPSEneDistribution::GenEpnHistEnergies(), G4SPSAngDistribution::GenerateUserDefPhi(), G4SPSAngDistribution::GenerateUserDefTheta(), G4SPSRandomGenerator::GenRandEnergy(), G4SPSRandomGenerator::GenRandPhi(), G4SPSRandomGenerator::GenRandPosPhi(), G4SPSRandomGenerator::GenRandPosTheta(), G4SPSRandomGenerator::GenRandTheta(), G4SPSRandomGenerator::GenRandX(), G4SPSRandomGenerator::GenRandY(), G4SPSRandomGenerator::GenRandZ(), G4SPSEneDistribution::GenUserHistEnergies(), G4SPSEneDistribution::GetArbEneWeight(), G4Cerenkov::GetAverageNumberOfPhotons(), G4PAIPhotData::GetEnergyPhotonTransfer(), G4PAIPhotData::GetEnergyPlasmonTransfer(), G4PAIModelData::GetEnergyTransfer(), G4PAIPhotData::GetEnergyTransfer(), G4PenelopeCrossSection::GetHardCrossSection(), G4PenelopeCrossSection::GetNormalizedShellCrossSection(), G4PAIPhotData::GetPlasmonRatio(), G4PenelopeCrossSection::GetShellCrossSection(), G4PenelopeCrossSection::GetSoftStoppingPower(), G4PenelopeCrossSection::GetTotalCrossSection(), G4eeToHadronsModel::Initialise(), G4NeutronElasticXS::Initialise(), G4NeutronInelasticXS::Initialise(), G4ParticleInelasticXS::Initialise(), G4SPSEneDistribution::LinearInterpolation(), G4SPSEneDistribution::LogInterpolation(), G4GDMLWriteMaterials::PropertyVectorWrite(), G4PAIPhotData::SampleAlongStepPhotonTransfer(), G4PAIPhotData::SampleAlongStepPlasmonTransfer(), G4PAIPhotData::SampleAlongStepTransfer(), G4PAIModelData::SampleAlongStepTransfer(), G4PAIPhotData::SamplePostStepPhotonTransfer(), G4PAIPhotData::SamplePostStepPlasmonTransfer(), G4PAIPhotData::SamplePostStepTransfer(), G4PAIModelData::SamplePostStepTransfer(), G4VEnergyLossProcess::SetLambdaTable(), G4SPSEneDistribution::SplineInterpolation(), and G4VEmProcess::StreamInfo().

◆ Initialise()

void G4PhysicsVector::Initialise ( )
protectedvirtual

◆ Interpolation()

G4double G4PhysicsVector::Interpolation ( const std::size_t  idx,
const G4double  energy 
) const
inlineprivate

◆ LogVectorValue()

G4double G4PhysicsVector::LogVectorValue ( const G4double  energy,
const G4double  theLogEnergy 
) const
inline

◆ operator!=()

G4bool G4PhysicsVector::operator!= ( const G4PhysicsVector right) const
delete

◆ operator()()

G4double G4PhysicsVector::operator() ( const std::size_t  index) const
inline

◆ operator=() [1/2]

G4PhysicsVector & G4PhysicsVector::operator= ( const G4PhysicsVector &&  )
delete

◆ operator=() [2/2]

G4PhysicsVector & G4PhysicsVector::operator= ( const G4PhysicsVector )
default

◆ operator==()

G4bool G4PhysicsVector::operator== ( const G4PhysicsVector right) const
delete

◆ operator[]()

G4double G4PhysicsVector::operator[] ( const std::size_t  index) const
inline

◆ PrintPutValueError()

void G4PhysicsVector::PrintPutValueError ( std::size_t  index,
G4double  value,
const G4String text 
)
protected

Definition at line 459 of file G4PhysicsVector.cc.

462{
464 ed << "Vector type: " << type << " length= " << numberOfNodes
465 << "; an attempt to put data at index= " << index
466 << " value= " << val << " in " << text;
467 G4Exception("G4PhysicsVector:", "gl0005",
468 FatalException, ed, "Wrong operation");
469}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40

References FatalException, G4Exception(), numberOfNodes, and type.

Referenced by G4PhysicsFreeVector::PutValues().

◆ PutValue()

void G4PhysicsVector::PutValue ( const std::size_t  index,
const G4double  value 
)
inline

Referenced by G4VXTRenergyLoss::BuildAngleForEnergyBank(), G4PolarizedCompton::BuildAsymmetryTable(), G4PolarizedAnnihilation::BuildAsymmetryTables(), G4PolarizedIonisation::BuildAsymmetryTables(), G4EmCorrections::BuildCorrectionVector(), G4KokoulinMuonNuclearXS::BuildCrossSectionTable(), G4hRDEnergyLoss::BuildDEDXTable(), G4LossTableBuilder::BuildDEDXTable(), G4VXTRenergyLoss::BuildEnergyTable(), G4VXTRenergyLoss::BuildGlobalAngleTable(), G4hRDEnergyLoss::BuildInverseRangeTable(), G4hRDEnergyLoss::BuildLabTimeVector(), G4hImpactIonisation::BuildLambdaTable(), G4hImpactIonisation::BuildLossTable(), G4VLEPTSModel::BuildMeanFreePathTable(), G4VLEPTSModel::BuildPhysicsTable(), G4ChargeExchangeProcess::BuildPhysicsTable(), G4hRDEnergyLoss::BuildProperTimeVector(), G4hRDEnergyLoss::BuildRangeCoeffATable(), G4hRDEnergyLoss::BuildRangeCoeffBTable(), G4hRDEnergyLoss::BuildRangeCoeffCTable(), G4LossTableBuilder::BuildRangeTable(), G4hRDEnergyLoss::BuildRangeVector(), G4LossTableBuilder::BuildTableForModel(), G4AdjointCSManager::BuildTotalSigmaTables(), G4ForwardXrayTR::BuildXrayTRtables(), G4eeToHadronsModel::ComputeCMCrossSectionPerElectron(), G4TablesForExtrapolator::ComputeElectronDEDX(), G4TablesForExtrapolator::ComputeMuonDEDX(), G4TablesForExtrapolator::ComputeProtonDEDX(), G4TablesForExtrapolator::ComputeTrasportXS(), G4EmModelManager::FillDEDXVector(), G4EmModelManager::FillLambdaVector(), G4NeutronElectronElXsc::G4NeutronElectronElXsc(), G4XNNElasticLowE::G4XNNElasticLowE(), G4XnpElasticLowE::G4XnpElasticLowE(), G4XnpTotalLowE::G4XnpTotalLowE(), G4NeutronElectronElXsc::Initialise(), G4PAIPhotData::Initialise(), G4PAIModelData::Initialise(), G4eeToHadronsModel::Initialise(), G4eplusTo2GammaOKVIModel::Initialise(), G4WentzelVIModel::Initialise(), G4DNABornExcitationModel2::Initialise(), G4HadronXSDataTable::Initialise(), G4LindhardSorensenData::InitialiseData(), G4InitXscPAI::IntegralCherenkov(), G4InitXscPAI::IntegralPAIdEdx(), G4InitXscPAI::IntegralPAIxSection(), G4InitXscPAI::IntegralPlasmon(), G4hRDEnergyLoss::InvertRangeVector(), G4IonICRU73Data::ReadElementData(), G4ChannelingECHARM::ReadFromECHARM(), G4IonICRU73Data::RetrieveVector(), and G4ChannelingMaterialData::SetBR().

◆ Retrieve()

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

Definition at line 87 of file G4PhysicsVector.cc.

88{
89 // clear properties;
90 dataVector.clear();
91 binVector.clear();
92 secDerivative.clear();
93
94 // retrieve in ascii mode
95 if(ascii)
96 {
97 // binning
98 fIn >> edgeMin >> edgeMax >> numberOfNodes;
99 if(fIn.fail() || numberOfNodes < 2)
100 {
101 return false;
102 }
103 // contents
104 G4int siz = 0;
105 fIn >> siz;
106 if(fIn.fail() || siz != G4int(numberOfNodes))
107 {
108 return false;
109 }
110
111 binVector.reserve(siz);
112 dataVector.reserve(siz);
113 G4double vBin, vData;
114
115 for(G4int i = 0; i < siz; ++i)
116 {
117 vBin = 0.;
118 vData = 0.;
119 fIn >> vBin >> vData;
120 if(fIn.fail())
121 {
122 return false;
123 }
124 binVector.push_back(vBin);
125 dataVector.push_back(vData);
126 }
127 Initialise();
128 return true;
129 }
130
131 // retrieve in binary mode
132 // binning
133 fIn.read((char*) (&edgeMin), sizeof edgeMin);
134 fIn.read((char*) (&edgeMax), sizeof edgeMax);
135 fIn.read((char*) (&numberOfNodes), sizeof numberOfNodes);
136
137 // contents
138 std::size_t size;
139 fIn.read((char*) (&size), sizeof size);
140
141 G4double* value = new G4double[2 * size];
142 fIn.read((char*) (value), 2 * size * (sizeof(G4double)));
143 if(G4int(fIn.gcount()) != G4int(2 * size * (sizeof(G4double))))
144 {
145 delete[] value;
146 return false;
147 }
148
149 binVector.reserve(size);
150 dataVector.reserve(size);
151 for(std::size_t i = 0; i < size; ++i)
152 {
153 binVector.push_back(value[2 * i]);
154 dataVector.push_back(value[2 * i + 1]);
155 }
156 delete[] value;
157
158 Initialise();
159 return true;
160}

References binVector, dataVector, edgeMax, edgeMin, Initialise(), numberOfNodes, and secDerivative.

Referenced by G4ComponentSAIDTotalXS::ReadData(), G4LivermorePhotoElectricModel::ReadData(), G4LivermoreComptonModel::ReadData(), G4LivermoreNuclearGammaConversionModel::ReadData(), G4LivermorePolarizedComptonModel::ReadData(), G4LivermorePolarizedGammaConversionModel::ReadData(), G4LivermoreRayleighModel::ReadData(), G4LowEPComptonModel::ReadData(), G4LowEPPolarizedComptonModel::ReadData(), G4BoldyshevTripletModel::ReadData(), G4LivermoreGammaConversion5DModel::ReadData(), G4LivermoreGammaConversionModel::ReadData(), G4LivermorePolarizedRayleighModel::ReadData(), G4PhysicsTable::RetrievePhysicsTable(), G4GammaNuclearXS::RetrieveVector(), G4IonICRU73Data::RetrieveVector(), G4NeutronCaptureXS::RetrieveVector(), G4NeutronInelasticXS::RetrieveVector(), and G4ParticleInelasticXS::RetrieveVector().

◆ ScaleVector()

void G4PhysicsVector::ScaleVector ( const G4double  factorE,
const G4double  factorV 
)

◆ SetVerboseLevel()

void G4PhysicsVector::SetVerboseLevel ( G4int  value)
inline

◆ Store()

G4bool G4PhysicsVector::Store ( std::ofstream &  fOut,
G4bool  ascii = false 
) const

Definition at line 55 of file G4PhysicsVector.cc.

56{
57 // Ascii mode
58 if(ascii)
59 {
60 fOut << *this;
61 return true;
62 }
63 // Binary Mode
64
65 // binning
66 fOut.write((char*) (&edgeMin), sizeof edgeMin);
67 fOut.write((char*) (&edgeMax), sizeof edgeMax);
68 fOut.write((char*) (&numberOfNodes), sizeof numberOfNodes);
69
70 // contents
71 std::size_t size = dataVector.size();
72 fOut.write((char*) (&size), sizeof size);
73
74 G4double* value = new G4double[2 * size];
75 for(std::size_t i = 0; i < size; ++i)
76 {
77 value[2 * i] = binVector[i];
78 value[2 * i + 1] = dataVector[i];
79 }
80 fOut.write((char*) (value), 2 * size * (sizeof(G4double)));
81 delete[] value;
82
83 return true;
84}

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

◆ Value() [1/2]

G4double G4PhysicsVector::Value ( const G4double  energy) const
inline

◆ Value() [2/2]

G4double G4PhysicsVector::Value ( const G4double  energy,
std::size_t &  lastidx 
) const
inline

Referenced by G4EmCorrections::BarkasCorrection(), G4EmCorrections::BuildCorrectionVector(), G4PenelopeRayleighModelMI::BuildFormFactorTable(), G4LossTableBuilder::BuildRangeTable(), G4OpBoundaryProcess::CalculateReflectivity(), G4Track::CalculateVelocityForOpticalPhoton(), G4eeToHadronsModel::ComputeCMCrossSectionPerElectron(), G4JAEAElasticScatteringModel::ComputeCrossSectionPerAtom(), G4JAEAPolarizedElasticScatteringModel::ComputeCrossSectionPerAtom(), G4LivermoreComptonModel::ComputeCrossSectionPerAtom(), G4LivermoreNuclearGammaConversionModel::ComputeCrossSectionPerAtom(), G4LivermorePolarizedComptonModel::ComputeCrossSectionPerAtom(), G4LivermorePolarizedGammaConversionModel::ComputeCrossSectionPerAtom(), G4LivermorePolarizedRayleighModel::ComputeCrossSectionPerAtom(), G4LivermoreRayleighModel::ComputeCrossSectionPerAtom(), G4LowEPComptonModel::ComputeCrossSectionPerAtom(), G4LowEPPolarizedComptonModel::ComputeCrossSectionPerAtom(), G4PenelopeGammaConversionModel::ComputeCrossSectionPerAtom(), G4PenelopePhotoElectricModel::ComputeCrossSectionPerAtom(), G4PenelopeRayleighModel::ComputeCrossSectionPerAtom(), G4PenelopeRayleighModelMI::ComputeCrossSectionPerAtom(), G4BoldyshevTripletModel::ComputeCrossSectionPerAtom(), G4eplusTo2GammaOKVIModel::ComputeCrossSectionPerAtom(), G4LivermoreGammaConversion5DModel::ComputeCrossSectionPerAtom(), G4LivermoreGammaConversionModel::ComputeCrossSectionPerAtom(), G4eeToHadronsModel::ComputeCrossSectionPerElectron(), G4PolarizedAnnihilation::ComputeSaturationFactor(), G4PolarizedCompton::ComputeSaturationFactor(), G4PolarizedIonisation::ComputeSaturationFactor(), G4eplusTo2GammaOKVIModel::CrossSectionPerVolume(), G4DNABornExcitationModel2::CrossSectionPerVolume(), G4EmCorrections::EffectiveChargeCorrection(), G4eeToHadronsModel::GenerateCMPhoton(), G4LivermorePolarizedRayleighModel::GenerateCosTheta(), G4Cerenkov::GetAverageNumberOfPhotons(), G4ChannelingMaterialData::GetBR(), G4ICRU90StoppingData::GetDEDX(), G4PenelopeIonisationXSHandler::GetDensityCorrection(), G4ChannelingECHARM::GetEC(), G4NeutronElectronElXsc::GetElementCrossSection(), G4KokoulinMuonNuclearXS::GetElementCrossSection(), G4PAIPhotData::GetEnergyPhotonTransfer(), G4PAIPhotData::GetEnergyPlasmonTransfer(), G4PAIModelData::GetEnergyTransfer(), G4PAIPhotData::GetEnergyTransfer(), G4PenelopeRayleighModel::GetFSquared(), G4PenelopeRayleighModelMI::GetFSquared(), G4PenelopeCrossSection::GetHardCrossSection(), G4OpAbsorption::GetMeanFreePath(), G4OpMieHG::GetMeanFreePath(), G4OpRayleigh::GetMeanFreePath(), G4OpWLS::GetMeanFreePath(), G4OpWLS2::GetMeanFreePath(), G4PenelopeCrossSection::GetNormalizedShellCrossSection(), G4SPSEneDistribution::GetProbability(), G4OpBoundaryProcess::GetReflectivity(), G4Scintillation::GetScintillationYieldByParticleType(), G4PenelopePhotoElectricModel::GetShellCrossSection(), G4PenelopeCrossSection::GetShellCrossSection(), G4PenelopeCrossSection::GetSoftStoppingPower(), G4PenelopeCrossSection::GetTotalCrossSection(), G4ElementData::GetValueForElement(), G4NeutronElectronElXsc::Initialise(), G4PAIPhotData::Initialise(), G4eeToHadronsModel::Initialise(), G4EmCorrections::KShellCorrection(), G4EmCorrections::LShellCorrection(), G4Cerenkov::PostStepDoIt(), G4OpBoundaryProcess::PostStepDoIt(), G4PenelopeBremsstrahlungAngular::PrepareTables(), G4DNABornExcitationModel2::RandomSelect(), G4IonICRU73Data::RetrieveVector(), G4PAIModelData::SampleAlongStepTransfer(), G4PenelopeBremsstrahlungAngular::SampleDirection(), G4PAIModelData::SamplePostStepTransfer(), G4PenelopeRayleighModel::SampleSecondaries(), G4PenelopeRayleighModelMI::SampleSecondaries(), G4eplusTo2GammaOKVIModel::SampleSecondaries(), G4VEnergyLossProcess::ScaledKinEnergyForLoss(), G4PenelopePhotoElectricModel::SelectRandomShell(), and G4EmCorrections::ShellCorrection().

Friends And Related Function Documentation

◆ operator<<

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

Definition at line 412 of file G4PhysicsVector.cc.

413{
414 // binning
415 G4int prec = out.precision();
416 out << std::setprecision(12) << pv.edgeMin << " " << pv.edgeMax << " "
417 << pv.numberOfNodes << G4endl;
418
419 // contents
420 out << pv.dataVector.size() << G4endl;
421 for(std::size_t i = 0; i < pv.dataVector.size(); ++i)
422 {
423 out << pv.binVector[i] << " " << pv.dataVector[i] << G4endl;
424 }
425 out << std::setprecision(prec);
426
427 return out;
428}
static const double prec
Definition: RanecuEngine.cc:61

Field Documentation

◆ binVector

std::vector<G4double> G4PhysicsVector::binVector
protected

◆ dataVector

std::vector<G4double> G4PhysicsVector::dataVector
protected

◆ edgeMax

G4double G4PhysicsVector::edgeMax = 0.0
protected

◆ edgeMin

G4double G4PhysicsVector::edgeMin = 0.0
protected

◆ idxmax

G4int G4PhysicsVector::idxmax = 0
protected

◆ invdBin

G4double G4PhysicsVector::invdBin = 0.0
protected

◆ logemin

G4double G4PhysicsVector::logemin = 0.0
protected

Definition at line 205 of file G4PhysicsVector.hh.

Referenced by G4PhysicsLogVector::Initialise().

◆ numberOfNodes

std::size_t G4PhysicsVector::numberOfNodes = 0
protected

◆ secDerivative

std::vector<G4double> G4PhysicsVector::secDerivative
protected

◆ type

G4PhysicsVectorType G4PhysicsVector::type = T_G4PhysicsFreeVector
protected

◆ useSpline

G4bool G4PhysicsVector::useSpline = false
private

Definition at line 220 of file G4PhysicsVector.hh.

Referenced by FillSecondDerivatives().

◆ verboseLevel

G4int G4PhysicsVector::verboseLevel = 0
protected

Definition at line 207 of file G4PhysicsVector.hh.

Referenced by FillSecondDerivatives().


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