G4Physics2DVector Class Reference

#include <G4Physics2DVector.hh>


Public Member Functions

 G4Physics2DVector ()
 G4Physics2DVector (size_t nx, size_t ny)
 G4Physics2DVector (const G4Physics2DVector &)
G4Physics2DVectoroperator= (const G4Physics2DVector &)
 ~G4Physics2DVector ()
G4double Value (G4double x, G4double y)
void PutX (size_t idx, G4double value)
void PutY (size_t idy, G4double value)
void PutValue (size_t idx, size_t idy, G4double value)
void PutVectors (const std::vector< G4double > &vecX, const std::vector< G4double > &vecY)
void ScaleVector (G4double factor)
G4double GetX (size_t index) const
G4double GetY (size_t index) const
G4double GetValue (size_t idx, size_t idy) const
size_t GetLengthX () const
size_t GetLengthY () const
G4PhysicsVectorType GetType () const
void SetBicubicInterpolation (G4bool)
void Store (std::ofstream &fOut)
G4bool Retrieve (std::ifstream &fIn)
G4double GetLastX () const
G4double GetLastY () const
G4double GetLastValue () const
size_t GetLastBinX () const
size_t GetLastBinY () const
void SetVerboseLevel (G4int value)
G4int GetVerboseLevel () const

Protected Member Functions

void PrepareVectors ()
void ClearVectors ()
void CopyData (const G4Physics2DVector &vec)
void ComputeValue (G4double x, G4double y)
void BicubicInterpolation (size_t idx, size_t idy)
size_t FindBinLocation (G4double z, const G4PV2DDataVector &)
void FindBin (G4double z, const G4PV2DDataVector &, size_t &lastidx)
void FindBinLocationX (G4double x)
void FindBinLocationY (G4double y)


Detailed Description

Definition at line 60 of file G4Physics2DVector.hh.


Constructor & Destructor Documentation

G4Physics2DVector::G4Physics2DVector (  ) 

Definition at line 47 of file G4Physics2DVector.cc.

00048  : type(T_G4PhysicsFreeVector),
00049    numberOfXNodes(0), numberOfYNodes(0),
00050    verboseLevel(0), useBicubic(false)
00051 {
00052   cache = new G4Physics2DVectorCache();
00053 }

G4Physics2DVector::G4Physics2DVector ( size_t  nx,
size_t  ny 
)

Definition at line 57 of file G4Physics2DVector.cc.

References PrepareVectors().

00058  : type(T_G4PhysicsFreeVector),
00059    numberOfXNodes(nx), numberOfYNodes(ny),
00060    verboseLevel(0), useBicubic(false)
00061 {
00062   cache = new G4Physics2DVectorCache();
00063   PrepareVectors();
00064 }

G4Physics2DVector::G4Physics2DVector ( const G4Physics2DVector  ) 

Definition at line 76 of file G4Physics2DVector.cc.

References CopyData(), numberOfXNodes, numberOfYNodes, PrepareVectors(), type, useBicubic, verboseLevel, xVector, and yVector.

00077 {
00078   type         = right.type;
00079 
00080   numberOfXNodes = right.numberOfXNodes;
00081   numberOfYNodes = right.numberOfYNodes;
00082 
00083   verboseLevel = right.verboseLevel;
00084   useBicubic   = right.useBicubic;
00085 
00086   xVector      = right.xVector;
00087   yVector      = right.yVector;
00088 
00089   cache = new G4Physics2DVectorCache();
00090   PrepareVectors();
00091   CopyData(right);
00092 }

G4Physics2DVector::~G4Physics2DVector (  ) 

Definition at line 68 of file G4Physics2DVector.cc.

References ClearVectors().

00069 {
00070   delete cache;
00071   ClearVectors();
00072 }


Member Function Documentation

void G4Physics2DVector::BicubicInterpolation ( size_t  idx,
size_t  idy 
) [protected]

Definition at line 206 of file G4Physics2DVector.cc.

References GetValue(), G4Physics2DVectorCache::lastValue, G4Physics2DVectorCache::lastX, and G4Physics2DVectorCache::lastY.

Referenced by ComputeValue().

00207 {
00208     // Bicubic interpolation according to 
00209     // 1. H.M. Antia, "Numerical Methods for Scientists and Engineers",
00210     //    MGH, 1991. 
00211     // 2. W.H. Press et al., "Numerical recipes. The Art of Scientific 
00212     //    Computing", Cambridge University Press, 2007. 
00213     G4double x1 = xVector[idx];
00214     G4double x2 = xVector[idx+1];
00215     G4double y1 = yVector[idy];
00216     G4double y2 = yVector[idy+1];
00217     G4double x  = cache->lastX;
00218     G4double y  = cache->lastY;
00219     G4double f1 = GetValue(idx,   idy);
00220     G4double f2 = GetValue(idx+1, idy);
00221     G4double f3 = GetValue(idx+1, idy+1);
00222     G4double f4 = GetValue(idx,   idy+1);
00223 
00224     G4double dx = x2 - x1;
00225     G4double dy = y2 - y1;
00226 
00227     G4double h1 = (x - x1)/dx;
00228     G4double h2 = (y - y1)/dy;   
00229 
00230     G4double h12 = h1*h1;
00231     G4double h13 = h12*h1;
00232     G4double h22 = h2*h2;
00233     G4double h23 = h22*h2;
00234 
00235     // Three derivatives at each of four points (1-4) defining the 
00236     // subregion are computed by numerical centered differencing from 
00237     // the functional values already tabulated on the grid. 
00238 
00239     G4double f1x = DerivativeX(idx, idy, dx);
00240     G4double f2x = DerivativeX(idx+1, idy, dx);
00241     G4double f3x = DerivativeX(idx+1, idy+1, dx);
00242     G4double f4x = DerivativeX(idx, idy+1, dx);
00243 
00244     G4double f1y = DerivativeY(idx, idy, dy);
00245     G4double f2y = DerivativeY(idx+1, idy, dy);
00246     G4double f3y = DerivativeY(idx+1, idy+1, dy);
00247     G4double f4y = DerivativeY(idx, idy+1, dy);
00248 
00249     G4double dxy = dx*dy;
00250     G4double f1xy = DerivativeXY(idx, idy, dxy);
00251     G4double f2xy = DerivativeXY(idx+1, idy, dxy);
00252     G4double f3xy = DerivativeXY(idx+1, idy+1, dxy);
00253     G4double f4xy = DerivativeXY(idx, idy+1, dxy);
00254 
00255     cache->lastValue = 
00256       f1 + f1y*h2 + (3*(f4-f1) - 2*f1y - f4y)*h22 + (2*(f1 - f4) + f1y + f4y)*h23
00257       + f1x*h1 + f1xy*h1*h2 +(3*(f4x - f1x) - 2*f1xy - f4xy)*h1*h22
00258       + (2*(f1x - f4x) + f1xy + f4xy)*h1*h23
00259       + (3*(f2 - f1) - 2*f1x - f2x)*h12 + (3*f2y - 3*f1y - 2*f1xy - f2xy)*h12*h2
00260       + (9*(f1 - f2 + f3 - f4) + 6*f1x + 3*f2x - 3*f3x - 6*f4x + 6*f1y - 6*f2y
00261          - 3*f3y + 3*f4y + 4*f1xy + 2*f2xy + f3xy + 2*f4xy)*h12*h22
00262       + (6*(-f1 + f2 - f3 + f4) - 4*f1x - 2*f2x + 2*f3x + 4*f4x - 3*f1y
00263          + 3*f2y + 3*f3y - 3*f4y - 2*f1xy - f2xy - f3xy - 2*f4xy)*h12*h23
00264       + (2*(f1 - f2) + f1x + f2x)*h13 + (2*(f1y - f2y) + f1xy + f2xy)*h13*h2
00265       + (6*(-f1 + f2 -f3 + f4) + 3*(-f1x - f2x + f3x + f4x) - 4*f1y
00266          + 4*f2y + 2*f3y - 2*f4y - 2*f1xy - 2*f2xy - f3xy - f4xy)*h13*h22
00267       + (4*(f1 - f2 + f3 - f4) + 2*(f1x + f2x - f3x - f4x) 
00268          + 2*(f1y - f2y - f3y + f4y) + f1xy + f2xy + f3xy + f4xy)*h13*h23;
00269 }

void G4Physics2DVector::ClearVectors (  )  [protected]

Definition at line 132 of file G4Physics2DVector.cc.

Referenced by operator=(), PutVectors(), Retrieve(), and ~G4Physics2DVector().

00133 {
00134   for(size_t j=0; j<numberOfYNodes; ++j) {
00135     delete value[j];
00136   }
00137 }

void G4Physics2DVector::ComputeValue ( G4double  x,
G4double  y 
) [protected]

Definition at line 157 of file G4Physics2DVector.cc.

References BicubicInterpolation(), FindBinLocationX(), FindBinLocationY(), GetValue(), G4Physics2DVectorCache::lastBinX, G4Physics2DVectorCache::lastBinY, G4Physics2DVectorCache::lastValue, G4Physics2DVectorCache::lastX, and G4Physics2DVectorCache::lastY.

Referenced by Value().

00158 {
00159   if(xx != cache->lastBinX) {
00160     if(xx <= xVector[0]) {
00161       cache->lastX = xVector[0];
00162       cache->lastBinX = 0;
00163     } else if(xx >= xVector[numberOfXNodes-1]) {
00164       cache->lastX = xVector[numberOfXNodes-1];
00165       cache->lastBinX = numberOfXNodes-2;
00166     } else {
00167       cache->lastX = xx;
00168       FindBinLocationX(xx);
00169     }
00170   }
00171   if(yy != cache->lastBinY) {
00172     if(yy <= yVector[0]) {
00173       cache->lastY = yVector[0];
00174       cache->lastBinY = 0;
00175     } else if(yy >= yVector[numberOfYNodes-1]) {
00176       cache->lastY = yVector[numberOfYNodes-1];
00177       cache->lastBinY = numberOfYNodes-2;
00178     } else {
00179       cache->lastY = yy;
00180       FindBinLocationY(yy);
00181     }
00182   }
00183   size_t idx  = cache->lastBinX;
00184   size_t idy  = cache->lastBinY;
00185   if(useBicubic) {
00186     BicubicInterpolation(idx, idy);
00187   } else {
00188     G4double x1 = xVector[idx];
00189     G4double x2 = xVector[idx+1];
00190     G4double y1 = yVector[idy];
00191     G4double y2 = yVector[idy+1];
00192     G4double x  = cache->lastX;
00193     G4double y  = cache->lastY;
00194     G4double v11= GetValue(idx,   idy);
00195     G4double v12= GetValue(idx+1, idy);
00196     G4double v21= GetValue(idx,   idy+1);
00197     G4double v22= GetValue(idx+1, idy+1);
00198     cache->lastValue = 
00199       ((y2 - y)*(v11*(x2 - x) + v12*(x - x1)) + 
00200        ((y - y1)*(v21*(x2 - x) + v22*(x - x1))))/((x2 - x1)*(y2 - y1)); 
00201   }
00202 }

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

Definition at line 141 of file G4Physics2DVector.cc.

References PutValue(), value, xVector, and yVector.

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

00142 {
00143   for(size_t i=0; i<numberOfXNodes; ++i) {
00144     xVector[i] = right.xVector[i];
00145   }
00146   for(size_t j=0; j<numberOfYNodes; ++j) {
00147     yVector[j] = right.yVector[j];
00148     G4PV2DDataVector* v0 = right.value[j];
00149     for(size_t i=0; i<numberOfXNodes; ++i) { 
00150       PutValue(i,j,(*v0)[i]); 
00151     }
00152   }
00153 }

void G4Physics2DVector::FindBin ( G4double  z,
const G4PV2DDataVector ,
size_t &  lastidx 
) [inline, protected]

Definition at line 122 of file G4Physics2DVector.icc.

References FindBinLocation().

Referenced by FindBinLocationX(), and FindBinLocationY().

00123 {
00124   if(z < v[idz] || z >=  v[idz]) { idz = FindBinLocation(z, v); }
00125 }

size_t G4Physics2DVector::FindBinLocation ( G4double  z,
const G4PV2DDataVector  
) [protected]

Definition at line 373 of file G4Physics2DVector.cc.

Referenced by FindBin().

00375 {
00376   size_t lowerBound = 0;
00377   size_t upperBound = v.size() - 2;
00378 
00379   while (lowerBound <= upperBound)
00380   {
00381     size_t midBin = (lowerBound + upperBound)/2;
00382     if( z < v[midBin] ) { upperBound = midBin-1; }
00383     else                { lowerBound = midBin+1; }
00384   }
00385 
00386   return upperBound;
00387 }

void G4Physics2DVector::FindBinLocationX ( G4double  x  )  [inline, protected]

Definition at line 127 of file G4Physics2DVector.icc.

References FindBin(), and G4Physics2DVectorCache::lastBinX.

Referenced by ComputeValue().

00128 {
00129   FindBin(x, xVector, cache->lastBinX);
00130 }

void G4Physics2DVector::FindBinLocationY ( G4double  y  )  [inline, protected]

Definition at line 132 of file G4Physics2DVector.icc.

References FindBin(), and G4Physics2DVectorCache::lastBinY.

Referenced by ComputeValue().

00133 {
00134   FindBin(y, yVector, cache->lastBinY);
00135 }

size_t G4Physics2DVector::GetLastBinX (  )  const [inline]

Definition at line 111 of file G4Physics2DVector.icc.

References G4Physics2DVectorCache::lastBinX.

00112 {
00113    return cache->lastBinX;
00114 }

size_t G4Physics2DVector::GetLastBinY (  )  const [inline]

Definition at line 116 of file G4Physics2DVector.icc.

References G4Physics2DVectorCache::lastBinY.

00117 {
00118    return cache->lastBinY;
00119 }

G4double G4Physics2DVector::GetLastValue (  )  const [inline]

Definition at line 106 of file G4Physics2DVector.icc.

References G4Physics2DVectorCache::lastValue.

00107 {
00108   return cache->lastValue;
00109 }

G4double G4Physics2DVector::GetLastX (  )  const [inline]

Definition at line 96 of file G4Physics2DVector.icc.

References G4Physics2DVectorCache::lastX.

00097 {
00098   return cache->lastX;
00099 }

G4double G4Physics2DVector::GetLastY (  )  const [inline]

Definition at line 101 of file G4Physics2DVector.icc.

References G4Physics2DVectorCache::lastY.

00102 {
00103   return cache->lastY;
00104 }

size_t G4Physics2DVector::GetLengthX (  )  const [inline]

Definition at line 76 of file G4Physics2DVector.icc.

00077 {
00078   return numberOfXNodes;
00079 }

size_t G4Physics2DVector::GetLengthY (  )  const [inline]

Definition at line 81 of file G4Physics2DVector.icc.

00082 {
00083   return numberOfYNodes;
00084 }

G4PhysicsVectorType G4Physics2DVector::GetType (  )  const [inline]

Definition at line 86 of file G4Physics2DVector.icc.

00087 {
00088   return type;
00089 }

G4double G4Physics2DVector::GetValue ( size_t  idx,
size_t  idy 
) const [inline]

Definition at line 71 of file G4Physics2DVector.icc.

Referenced by BicubicInterpolation(), ComputeValue(), ScaleVector(), and Store().

00072 {
00073   return (*(value[idy]))[idx];
00074 }

G4int G4Physics2DVector::GetVerboseLevel (  )  const [inline]

Definition at line 142 of file G4Physics2DVector.icc.

00143 {
00144    return verboseLevel;
00145 }

G4double G4Physics2DVector::GetX ( size_t  index  )  const [inline]

Definition at line 60 of file G4Physics2DVector.icc.

00061 {
00062   return xVector[index];
00063 }

G4double G4Physics2DVector::GetY ( size_t  index  )  const [inline]

Definition at line 65 of file G4Physics2DVector.icc.

00066 {
00067   return yVector[index];
00068 }

G4Physics2DVector & G4Physics2DVector::operator= ( const G4Physics2DVector  ) 

Definition at line 96 of file G4Physics2DVector.cc.

References G4Physics2DVectorCache::Clear(), ClearVectors(), CopyData(), numberOfXNodes, numberOfYNodes, PrepareVectors(), type, useBicubic, and verboseLevel.

00097 {
00098   if (&right==this)  { return *this; }
00099   ClearVectors();
00100 
00101   type         = right.type;
00102 
00103   numberOfXNodes = right.numberOfXNodes;
00104   numberOfYNodes = right.numberOfYNodes;
00105 
00106   verboseLevel = right.verboseLevel;
00107   useBicubic   = right.useBicubic;
00108 
00109   cache->Clear();
00110   PrepareVectors();
00111   CopyData(right);
00112 
00113   return *this;
00114 }

void G4Physics2DVector::PrepareVectors (  )  [protected]

Definition at line 118 of file G4Physics2DVector.cc.

Referenced by G4Physics2DVector(), operator=(), PutVectors(), and Retrieve().

00119 {
00120   xVector.resize(numberOfXNodes,0.);
00121   yVector.resize(numberOfYNodes,0.);
00122   value.resize(numberOfYNodes,0);
00123   for(size_t j=0; j<numberOfYNodes; ++j) {
00124     G4PV2DDataVector* v = new G4PV2DDataVector();
00125     v->resize(numberOfXNodes,0.);
00126     value[j] = v;
00127   }
00128 }

void G4Physics2DVector::PutValue ( size_t  idx,
size_t  idy,
G4double  value 
) [inline]

Definition at line 55 of file G4Physics2DVector.icc.

Referenced by CopyData(), Retrieve(), and ScaleVector().

00056 {
00057   (*(value[idy]))[idx] = val;
00058 }

void G4Physics2DVector::PutVectors ( const std::vector< G4double > &  vecX,
const std::vector< G4double > &  vecY 
)

Definition at line 274 of file G4Physics2DVector.cc.

References G4Physics2DVectorCache::Clear(), ClearVectors(), and PrepareVectors().

00276 {
00277   ClearVectors();
00278   numberOfXNodes = vecX.size();
00279   numberOfYNodes = vecY.size();
00280   PrepareVectors();
00281   if(!cache) { cache = new G4Physics2DVectorCache(); }
00282   cache->Clear();
00283   for(size_t i = 0; i<numberOfXNodes; ++i) {
00284     xVector[i] = vecX[i];
00285   }
00286   for(size_t j = 0; j<numberOfYNodes; ++j) {
00287     yVector[j] = vecY[j];
00288   }
00289 }

void G4Physics2DVector::PutX ( size_t  idx,
G4double  value 
) [inline]

Definition at line 44 of file G4Physics2DVector.icc.

00045 {
00046   xVector[idx] = val;
00047 }

void G4Physics2DVector::PutY ( size_t  idy,
G4double  value 
) [inline]

Definition at line 49 of file G4Physics2DVector.icc.

00050 {
00051   yVector[idy] = val;
00052 }

G4bool G4Physics2DVector::Retrieve ( std::ifstream &  fIn  ) 

Definition at line 322 of file G4Physics2DVector.cc.

References G4Physics2DVectorCache::Clear(), ClearVectors(), PrepareVectors(), and PutValue().

00323 {
00324   // initialisation
00325   cache->Clear();
00326   ClearVectors();
00327 
00328   // binning
00329   G4int k;
00330   in >> k >> numberOfXNodes >> numberOfYNodes;
00331   if (in.fail())  { return false; }
00332   PrepareVectors();
00333   type = G4PhysicsVectorType(k); 
00334 
00335   // contents
00336   G4double val;
00337   for(size_t i = 0; i<numberOfXNodes; ++i) {
00338     in >> xVector[i];
00339     if (in.fail())  { return false; }
00340   }
00341   for(size_t j = 0; j<numberOfYNodes; ++j) {
00342     in >> yVector[j];
00343     if (in.fail())  { return false; }
00344   }
00345   for(size_t j = 0; j<numberOfYNodes; ++j) {
00346     for(size_t i = 0; i<numberOfXNodes; ++i) {
00347       in >> val;
00348       if (in.fail())  { return false; }
00349       PutValue(i, j, val);
00350     }
00351   }
00352   in.close();
00353   return true;
00354 }

void G4Physics2DVector::ScaleVector ( G4double  factor  ) 

Definition at line 359 of file G4Physics2DVector.cc.

References GetValue(), and PutValue().

00360 {
00361   G4double val;
00362   for(size_t j = 0; j<numberOfYNodes; ++j) {
00363     for(size_t i = 0; i<numberOfXNodes; ++i) {
00364       val = GetValue(i, j)*factor;
00365       PutValue(i, j, val);
00366     }
00367   }
00368 }

void G4Physics2DVector::SetBicubicInterpolation ( G4bool   )  [inline]

Definition at line 91 of file G4Physics2DVector.icc.

00092 {
00093   useBicubic = val;
00094 }

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

Definition at line 137 of file G4Physics2DVector.icc.

00138 {
00139    verboseLevel = val;
00140 }

void G4Physics2DVector::Store ( std::ofstream &  fOut  ) 

Definition at line 293 of file G4Physics2DVector.cc.

References G4endl, and GetValue().

00294 {
00295   // binning
00296   G4int prec = out.precision();
00297   out << G4int(type) << " " << numberOfXNodes << " " << numberOfYNodes 
00298       << G4endl; 
00299   out << std::setprecision(5);
00300 
00301   // contents
00302   for(size_t i = 0; i<numberOfXNodes-1; ++i) {
00303     out << xVector[i] << "  ";
00304   }
00305   out << xVector[numberOfXNodes-1] << G4endl;
00306   for(size_t j = 0; j<numberOfYNodes-1; ++j) {
00307     out << yVector[j] << "  ";
00308   }
00309   out << yVector[numberOfYNodes-1] << G4endl;
00310   for(size_t j = 0; j<numberOfYNodes; ++j) {
00311     for(size_t i = 0; i<numberOfXNodes-1; ++i) {
00312       out << GetValue(i, j) << "  ";
00313     }
00314     out << GetValue(numberOfXNodes-1,j) << G4endl;
00315   }
00316   out.precision(prec);
00317   out.close();
00318 }

G4double G4Physics2DVector::Value ( G4double  x,
G4double  y 
) [inline]

Definition at line 38 of file G4Physics2DVector.icc.

References ComputeValue(), G4Physics2DVectorCache::lastValue, G4Physics2DVectorCache::lastX, and G4Physics2DVectorCache::lastY.

Referenced by G4SeltzerBergerModel::ComputeDXSectionPerAtom(), and G4SeltzerBergerModel::SampleSecondaries().

00039 {
00040   if(x != cache->lastX || y != cache->lastY) { ComputeValue(x, y); }
00041   return cache->lastValue;
00042 }


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