Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4PhysicsVector Class Reference

#include <G4PhysicsVector.hh>

Inheritance diagram for G4PhysicsVector:

Public Member Functions

 G4PhysicsVector (G4bool spline=false)
 
 G4PhysicsVector (const G4PhysicsVector &)
 
G4PhysicsVectoroperator= (const G4PhysicsVector &)
 
virtual ~G4PhysicsVector ()
 
G4double Value (G4double theEnergy, size_t &lastidx) const
 
G4double Value (G4double theEnergy) const
 
G4double GetValue (G4double theEnergy, G4bool &isOutRange) const
 
G4int operator== (const G4PhysicsVector &right) const
 
G4int operator!= (const G4PhysicsVector &right) const
 
G4double operator[] (const size_t index) const
 
G4double operator() (const size_t index) const
 
void PutValue (size_t index, G4double theValue)
 
virtual void ScaleVector (G4double factorE, G4double factorV)
 
G4double Energy (size_t index) const
 
G4double GetMaxEnergy () const
 
G4double GetLowEdgeEnergy (size_t binNumber) const
 
size_t GetVectorLength () const
 
size_t FindBin (G4double energy, size_t idx) const
 
void FillSecondDerivatives ()
 
void ComputeSecDerivatives ()
 
void ComputeSecondDerivatives (G4double firstPointDerivative, G4double endPointDerivative)
 
G4double FindLinearEnergy (G4double rand) const
 
G4bool IsFilledVectorExist () const
 
G4PhysicsVectorType GetType () const
 
void SetSpline (G4bool)
 
G4bool Store (std::ofstream &fOut, G4bool ascii=false) const
 
virtual G4bool Retrieve (std::ifstream &fIn, G4bool ascii=false)
 
void DumpValues (G4double unitE=1.0, G4double unitV=1.0) const
 
void SetVerboseLevel (G4int value)
 

Protected Member Functions

void DeleteData ()
 
void CopyData (const G4PhysicsVector &vec)
 
void PrintPutValueError (size_t index)
 

Protected Attributes

G4PhysicsVectorType type
 
G4double edgeMin
 
G4double edgeMax
 
size_t numberOfNodes
 
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)
explicit

Definition at line 61 of file G4PhysicsVector.cc.

63  edgeMin(0.), edgeMax(0.), numberOfNodes(0),
64  useSpline(val),
65  dBin(0.), baseBin(0.),
66  verboseLevel(0)
67 {}
G4PhysicsVectorType type
G4PhysicsVector::G4PhysicsVector ( const G4PhysicsVector right)

Definition at line 76 of file G4PhysicsVector.cc.

77 {
78  dBin = right.dBin;
79  baseBin = right.baseBin;
80  verboseLevel = right.verboseLevel;
81 
82  DeleteData();
83  CopyData(right);
84 }
void CopyData(const G4PhysicsVector &vec)

Here is the call graph for this function:

G4PhysicsVector::~G4PhysicsVector ( )
virtual

Definition at line 71 of file G4PhysicsVector.cc.

72 {}

Member Function Documentation

void G4PhysicsVector::ComputeSecDerivatives ( )

Definition at line 431 of file G4PhysicsVector.cc.

433 {
434  if(3 > numberOfNodes) // cannot compute derivatives for less than 4 bins
435  {
436  useSpline = false;
437  return;
438  }
439 
440  if(!SplinePossible()) { return; }
441 
442  useSpline = true;
443 
444  size_t n = numberOfNodes-1;
445 
446  for(size_t i=1; i<n; ++i)
447  {
448  secDerivative[i] =
449  3.0*((dataVector[i+1]-dataVector[i])/(binVector[i+1]-binVector[i]) -
450  (dataVector[i]-dataVector[i-1])/(binVector[i]-binVector[i-1]))
451  /(binVector[i+1]-binVector[i-1]);
452  }
453  secDerivative[n] = secDerivative[n-1];
455 }
G4PVDataVector dataVector
G4PVDataVector binVector
G4PVDataVector secDerivative
const G4int n

Here is the caller graph for this function:

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

Definition at line 294 of file G4PhysicsVector.cc.

300 {
301  if(4 > numberOfNodes) // cannot compute derivatives for less than 4 bins
302  {
304  return;
305  }
306 
307  if(!SplinePossible()) { return; }
308 
309  useSpline = true;
310 
311  G4int n = numberOfNodes-1;
312 
313  G4double* u = new G4double [n];
314 
315  G4double p, sig, un;
316 
317  u[0] = (6.0/(binVector[1]-binVector[0]))
318  * ((dataVector[1]-dataVector[0])/(binVector[1]-binVector[0])
319  - firstPointDerivative);
320 
321  secDerivative[0] = - 0.5;
322 
323  // Decomposition loop for tridiagonal algorithm. secDerivative[i]
324  // and u[i] are used for temporary storage of the decomposed factors.
325 
326  for(G4int i=1; i<n; ++i)
327  {
328  sig = (binVector[i]-binVector[i-1]) / (binVector[i+1]-binVector[i-1]);
329  p = sig*(secDerivative[i-1]) + 2.0;
330  secDerivative[i] = (sig - 1.0)/p;
331  u[i] = (dataVector[i+1]-dataVector[i])/(binVector[i+1]-binVector[i])
332  - (dataVector[i]-dataVector[i-1])/(binVector[i]-binVector[i-1]);
333  u[i] = 6.0*u[i]/(binVector[i+1]-binVector[i-1]) - sig*u[i-1]/p;
334  }
335 
336  sig = (binVector[n-1]-binVector[n-2]) / (binVector[n]-binVector[n-2]);
337  p = sig*secDerivative[n-2] + 2.0;
338  un = (6.0/(binVector[n]-binVector[n-1]))
339  *(endPointDerivative -
340  (dataVector[n]-dataVector[n-1])/(binVector[n]-binVector[n-1])) - u[n-1]/p;
341  secDerivative[n] = un/(secDerivative[n-1] + 2.0);
342 
343  // The back-substitution loop for the triagonal algorithm of solving
344  // a linear system of equations.
345 
346  for(G4int k=n-1; k>0; --k)
347  {
348  secDerivative[k] *=
349  (secDerivative[k+1] -
350  u[k]*(binVector[k+1]-binVector[k-1])/(binVector[k+1]-binVector[k]));
351  }
352  secDerivative[0] = 0.5*(u[0] - secDerivative[1]);
353 
354  delete [] u;
355 }
G4PVDataVector dataVector
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78
G4PVDataVector binVector
G4PVDataVector secDerivative
const G4int n
double G4double
Definition: G4Types.hh:76
void ComputeSecDerivatives()

Here is the call graph for this function:

void G4PhysicsVector::CopyData ( const G4PhysicsVector vec)
protected

Definition at line 124 of file G4PhysicsVector.cc.

125 {
126  type = vec.type;
127  edgeMin = vec.edgeMin;
128  edgeMax = vec.edgeMax;
130  useSpline = vec.useSpline;
131 
132  size_t i;
133  dataVector.resize(numberOfNodes);
134  for(i=0; i<numberOfNodes; ++i) {
135  dataVector[i] = (vec.dataVector)[i];
136  }
137  binVector.resize(numberOfNodes);
138  for(i=0; i<numberOfNodes; ++i) {
139  binVector[i] = (vec.binVector)[i];
140  }
141  if(0 < (vec.secDerivative).size()) {
142  secDerivative.resize(numberOfNodes);
143  for(i=0; i<numberOfNodes; ++i){
144  secDerivative[i] = (vec.secDerivative)[i];
145  }
146  }
147 }
G4PVDataVector dataVector
G4PVDataVector binVector
G4PhysicsVectorType type
G4PVDataVector secDerivative

Here is the caller graph for this function:

void G4PhysicsVector::DeleteData ( )
protected

Definition at line 116 of file G4PhysicsVector.cc.

117 {
118  useSpline = false;
119  secDerivative.clear();
120 }
G4PVDataVector secDerivative

Here is the caller graph for this function:

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

Definition at line 266 of file G4PhysicsVector.cc.

267 {
268  for (size_t i = 0; i < numberOfNodes; ++i)
269  {
270  G4cout << binVector[i]/unitE << " " << dataVector[i]/unitV << G4endl;
271  }
272 }
G4PVDataVector dataVector
G4PVDataVector binVector
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61

Here is the caller graph for this function:

G4double G4PhysicsVector::Energy ( size_t  index) const
inline

Here is the caller graph for this function:

void G4PhysicsVector::FillSecondDerivatives ( )

Definition at line 359 of file G4PhysicsVector.cc.

363 {
364  if(5 > numberOfNodes) // cannot compute derivatives for less than 4 points
365  {
367  return;
368  }
369 
370  if(!SplinePossible()) { return; }
371 
372  useSpline = true;
373 
374  G4int n = numberOfNodes-1;
375 
376  G4double* u = new G4double [n];
377 
378  G4double p, sig;
379 
380  u[1] = ((dataVector[2]-dataVector[1])/(binVector[2]-binVector[1]) -
381  (dataVector[1]-dataVector[0])/(binVector[1]-binVector[0]));
382  u[1] = 6.0*u[1]*(binVector[2]-binVector[1])
383  / ((binVector[2]-binVector[0])*(binVector[2]-binVector[0]));
384 
385  // Decomposition loop for tridiagonal algorithm. secDerivative[i]
386  // and u[i] are used for temporary storage of the decomposed factors.
387 
388  secDerivative[1] = (2.0*binVector[1]-binVector[0]-binVector[2])
389  / (2.0*binVector[2]-binVector[0]-binVector[1]);
390 
391  for(G4int i=2; i<n-1; ++i)
392  {
393  sig = (binVector[i]-binVector[i-1]) / (binVector[i+1]-binVector[i-1]);
394  p = sig*secDerivative[i-1] + 2.0;
395  secDerivative[i] = (sig - 1.0)/p;
396  u[i] = (dataVector[i+1]-dataVector[i])/(binVector[i+1]-binVector[i])
397  - (dataVector[i]-dataVector[i-1])/(binVector[i]-binVector[i-1]);
398  u[i] = (6.0*u[i]/(binVector[i+1]-binVector[i-1])) - sig*u[i-1]/p;
399  }
400 
401  sig = (binVector[n-1]-binVector[n-2]) / (binVector[n]-binVector[n-2]);
402  p = sig*secDerivative[n-3] + 2.0;
403  u[n-1] = (dataVector[n]-dataVector[n-1])/(binVector[n]-binVector[n-1])
404  - (dataVector[n-1]-dataVector[n-2])/(binVector[n-1]-binVector[n-2]);
405  u[n-1] = 6.0*sig*u[n-1]/(binVector[n]-binVector[n-2])
406  - (2.0*sig - 1.0)*u[n-2]/p;
407 
408  p = (1.0+sig) + (2.0*sig-1.0)*secDerivative[n-2];
409  secDerivative[n-1] = u[n-1]/p;
410 
411  // The back-substitution loop for the triagonal algorithm of solving
412  // a linear system of equations.
413 
414  for(G4int k=n-2; k>1; --k)
415  {
416  secDerivative[k] *=
417  (secDerivative[k+1] -
418  u[k]*(binVector[k+1]-binVector[k-1])/(binVector[k+1]-binVector[k]));
419  }
420  secDerivative[n] = (secDerivative[n-1] - (1.0-sig)*secDerivative[n-2])/sig;
421  sig = 1.0 - ((binVector[2]-binVector[1])/(binVector[2]-binVector[0]));
422  secDerivative[1] *= (secDerivative[2] - u[1]/(1.0-sig));
423  secDerivative[0] = (secDerivative[1] - sig*secDerivative[2])/(1.0-sig);
424 
425  delete [] u;
426 }
G4PVDataVector dataVector
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78
G4PVDataVector binVector
G4PVDataVector secDerivative
const G4int n
double G4double
Definition: G4Types.hh:76
void ComputeSecDerivatives()

Here is the call graph for this function:

Here is the caller graph for this function:

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

Here is the caller graph for this function:

G4double G4PhysicsVector::FindLinearEnergy ( G4double  rand) const

Definition at line 516 of file G4PhysicsVector.cc.

517 {
518  if(1 >= numberOfNodes) { return 0.0; }
519  G4double y = rand*dataVector[numberOfNodes-1];
520  size_t bin = std::lower_bound(dataVector.begin(), dataVector.end(), y)
521  - dataVector.begin() - 1;
522  bin = std::min(bin, numberOfNodes-2);
523  G4double res = binVector[bin];
524  G4double del = dataVector[bin+1] - dataVector[bin];
525  if(del > 0.0) {
526  res += (y - dataVector[bin])*(binVector[bin+1] - res)/del;
527  }
528  return res;
529 }
G4PVDataVector dataVector
tuple bin
Definition: plottest35.py:22
G4PVDataVector binVector
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4PhysicsVector::GetLowEdgeEnergy ( size_t  binNumber) const

Definition at line 151 of file G4PhysicsVector.cc.

152 {
153  return binVector[binNumber];
154 }
G4PVDataVector binVector

Here is the caller graph for this function:

G4double G4PhysicsVector::GetMaxEnergy ( ) const
inline

Here is the caller graph for this function:

G4PhysicsVectorType G4PhysicsVector::GetType ( ) const
inline
G4double G4PhysicsVector::GetValue ( G4double  theEnergy,
G4bool isOutRange 
) const
inline

Here is the caller graph for this function:

size_t G4PhysicsVector::GetVectorLength ( ) const
inline

Here is the caller graph for this function:

G4bool G4PhysicsVector::IsFilledVectorExist ( ) const
inline
G4int G4PhysicsVector::operator!= ( const G4PhysicsVector right) const

Definition at line 109 of file G4PhysicsVector.cc.

110 {
111  return (this != &right);
112 }
G4double G4PhysicsVector::operator() ( const size_t  index) const
inline
G4PhysicsVector & G4PhysicsVector::operator= ( const G4PhysicsVector right)

Definition at line 88 of file G4PhysicsVector.cc.

89 {
90  if (&right==this) { return *this; }
91  dBin = right.dBin;
92  baseBin = right.baseBin;
93  verboseLevel = right.verboseLevel;
94 
95  DeleteData();
96  CopyData(right);
97  return *this;
98 }
void CopyData(const G4PhysicsVector &vec)

Here is the call graph for this function:

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

Definition at line 102 of file G4PhysicsVector.cc.

103 {
104  return (this == &right);
105 }
G4double G4PhysicsVector::operator[] ( const size_t  index) const
inline
void G4PhysicsVector::PrintPutValueError ( size_t  index)
protected

Definition at line 533 of file G4PhysicsVector.cc.

534 {
536  ed << "Vector type " << type << " length= " << numberOfNodes
537  << " an attempt to put data at index= " << index;
538  G4Exception("G4PhysicsVector::PutValue()","gl0005",FatalException,
539  ed,"Memory overwritten");
540 
541 }
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4PhysicsVectorType type
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

Here is the call graph for this function:

Here is the caller graph for this function:

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

Here is the caller graph for this function:

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

Reimplemented in G4PhysicsLogVector, and G4PhysicsLinearVector.

Definition at line 191 of file G4PhysicsVector.cc.

192 {
193  // clear properties;
194  dataVector.clear();
195  binVector.clear();
196  secDerivative.clear();
197 
198  // retrieve in ascii mode
199  if (ascii){
200  // binning
201  fIn >> edgeMin >> edgeMax >> numberOfNodes;
202  if (fIn.fail()) { return false; }
203  // contents
204  G4int siz=0;
205  fIn >> siz;
206  if (fIn.fail() || siz<=0) { return false; }
207 
208  binVector.reserve(siz);
209  dataVector.reserve(siz);
210  G4double vBin, vData;
211 
212  for(G4int i = 0; i < siz ; i++)
213  {
214  vBin = 0.;
215  vData= 0.;
216  fIn >> vBin >> vData;
217  if (fIn.fail()) { return false; }
218  binVector.push_back(vBin);
219  dataVector.push_back(vData);
220  }
221 
222  // to remove any inconsistency
223  numberOfNodes = siz;
224  edgeMin = binVector[0];
225  edgeMax = binVector[numberOfNodes-1];
226  return true ;
227  }
228 
229  // retrieve in binary mode
230  // binning
231  fIn.read((char*)(&edgeMin), sizeof edgeMin);
232  fIn.read((char*)(&edgeMax), sizeof edgeMax);
233  fIn.read((char*)(&numberOfNodes), sizeof numberOfNodes );
234 
235  // contents
236  size_t size;
237  fIn.read((char*)(&size), sizeof size);
238 
239  G4double* value = new G4double[2*size];
240  fIn.read((char*)(value), 2*size*(sizeof(G4double)) );
241  if (G4int(fIn.gcount()) != G4int(2*size*(sizeof(G4double))) )
242  {
243  delete [] value;
244  return false;
245  }
246 
247  binVector.reserve(size);
248  dataVector.reserve(size);
249  for(size_t i = 0; i < size; ++i)
250  {
251  binVector.push_back(value[2*i]);
252  dataVector.push_back(value[2*i+1]);
253  }
254  delete [] value;
255 
256  // to remove any inconsistency
257  numberOfNodes = size;
258  edgeMin = binVector[0];
259  edgeMax = binVector[numberOfNodes-1];
260 
261  return true;
262 }
G4PVDataVector dataVector
int G4int
Definition: G4Types.hh:78
G4PVDataVector binVector
const XML_Char int const XML_Char * value
Definition: expat.h:331
G4PVDataVector secDerivative
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

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

Reimplemented in G4PhysicsLogVector, and G4PhysicsLinearVector.

Definition at line 277 of file G4PhysicsVector.cc.

278 {
279  size_t n = dataVector.size();
280  size_t i;
281  for(i=0; i<n; ++i) {
282  binVector[i] *= factorE;
283  dataVector[i] *= factorV;
284  }
285  secDerivative.clear();
286 
287  edgeMin = binVector[0];
288  edgeMax = binVector[n-1];
289 }
G4PVDataVector dataVector
G4PVDataVector binVector
G4PVDataVector secDerivative
const G4int n

Here is the caller graph for this function:

void G4PhysicsVector::SetSpline ( G4bool  )
inline

Here is the caller graph for this function:

void G4PhysicsVector::SetVerboseLevel ( G4int  value)
inline
G4bool G4PhysicsVector::Store ( std::ofstream &  fOut,
G4bool  ascii = false 
) const

Definition at line 158 of file G4PhysicsVector.cc.

159 {
160  // Ascii mode
161  if (ascii)
162  {
163  fOut << *this;
164  return true;
165  }
166  // Binary Mode
167 
168  // binning
169  fOut.write((char*)(&edgeMin), sizeof edgeMin);
170  fOut.write((char*)(&edgeMax), sizeof edgeMax);
171  fOut.write((char*)(&numberOfNodes), sizeof numberOfNodes);
172 
173  // contents
174  size_t size = dataVector.size();
175  fOut.write((char*)(&size), sizeof size);
176 
177  G4double* value = new G4double[2*size];
178  for(size_t i = 0; i < size; ++i)
179  {
180  value[2*i] = binVector[i];
181  value[2*i+1]= dataVector[i];
182  }
183  fOut.write((char*)(value), 2*size*(sizeof (G4double)));
184  delete [] value;
185 
186  return true;
187 }
G4PVDataVector dataVector
G4PVDataVector binVector
const XML_Char int const XML_Char * value
Definition: expat.h:331
double G4double
Definition: G4Types.hh:76
G4double G4PhysicsVector::Value ( G4double  theEnergy,
size_t &  lastidx 
) const

Definition at line 498 of file G4PhysicsVector.cc.

499 {
500  G4double y;
501  if(theEnergy <= edgeMin) {
502  lastIdx = 0;
503  y = dataVector[0];
504  } else if(theEnergy >= edgeMax) {
505  lastIdx = numberOfNodes-1;
506  y = dataVector[lastIdx];
507  } else {
508  lastIdx = FindBin(theEnergy, lastIdx);
509  y = Interpolation(lastIdx, theEnergy);
510  }
511  return y;
512 }
G4PVDataVector dataVector
size_t FindBin(G4double energy, size_t idx) const
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4PhysicsVector::Value ( G4double  theEnergy) const
inline

Friends And Related Function Documentation

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

Definition at line 479 of file G4PhysicsVector.cc.

480 {
481  // binning
482  out << std::setprecision(12) << pv.edgeMin << " "
483  << pv.edgeMax << " " << pv.numberOfNodes << G4endl;
484 
485  // contents
486  out << pv.dataVector.size() << G4endl;
487  for(size_t i = 0; i < pv.dataVector.size(); i++)
488  {
489  out << pv.binVector[i] << " " << pv.dataVector[i] << G4endl;
490  }
491  out << std::setprecision(6);
492 
493  return out;
494 }
G4PVDataVector dataVector
G4PVDataVector binVector
#define G4endl
Definition: G4ios.hh:61

Member Data Documentation

G4double G4PhysicsVector::baseBin
protected

Definition at line 238 of file G4PhysicsVector.hh.

G4PVDataVector G4PhysicsVector::binVector
protected

Definition at line 215 of file G4PhysicsVector.hh.

G4PVDataVector G4PhysicsVector::dataVector
protected

Definition at line 214 of file G4PhysicsVector.hh.

G4double G4PhysicsVector::dBin
protected

Definition at line 237 of file G4PhysicsVector.hh.

G4double G4PhysicsVector::edgeMax
protected

Definition at line 210 of file G4PhysicsVector.hh.

G4double G4PhysicsVector::edgeMin
protected

Definition at line 209 of file G4PhysicsVector.hh.

size_t G4PhysicsVector::numberOfNodes
protected

Definition at line 212 of file G4PhysicsVector.hh.

G4PVDataVector G4PhysicsVector::secDerivative
protected

Definition at line 216 of file G4PhysicsVector.hh.

G4PhysicsVectorType G4PhysicsVector::type
protected

Definition at line 207 of file G4PhysicsVector.hh.

G4int G4PhysicsVector::verboseLevel
protected

Definition at line 240 of file G4PhysicsVector.hh.


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