2 // ********************************************************************
3 // * License and Disclaimer *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
27 // $Id: G4PhysicsVector.icc 93409 2015-10-21 13:26:27Z gcosmo $
30 //---------------------------------------------------------------
31 // GEANT 4 class source file
33 // G4PhysicsVector.icc
36 // A physics vector which has values of energy-loss, cross-section,
37 // and other physics values of a particle in matter in a given
38 // range of the energy, momentum, etc.
39 // This class serves as the base class for a vector having various
40 // energy scale, for example like 'log', 'linear', 'free', etc.
42 //---------------------------------------------------------------
45 G4double G4PhysicsVector::operator[](const size_t binNumber) const
47 return dataVector[binNumber];
50 //---------------------------------------------------------------
53 G4double G4PhysicsVector::operator()(const size_t binNumber) const
55 return dataVector[binNumber];
58 //---------------------------------------------------------------
61 G4double G4PhysicsVector::Energy(const size_t binNumber) const
63 return binVector[binNumber];
66 //---------------------------------------------------------------
69 G4double G4PhysicsVector::GetMaxEnergy() const
74 //---------------------------------------------------------------
77 size_t G4PhysicsVector::GetVectorLength() const
82 //---------------------------------------------------------------
85 G4double G4PhysicsVector::GetValue(G4double theEnergy, G4bool&) const
88 return Value(theEnergy, idx);
91 //------------------------------------------------
94 G4double G4PhysicsVector::LinearInterpolation(size_t idx, G4double e) const
96 // Linear interpolation is used to get the value. Before this method
97 // is called it is ensured that the energy is inside the bin
98 // 0 < idx < numberOfNodes-1
100 return dataVector[idx] +
101 ( dataVector[idx + 1]-dataVector[idx] ) * (e - binVector[idx])
102 /( binVector[idx + 1]-binVector[idx] );
105 //---------------------------------------------------------------
108 G4double G4PhysicsVector::SplineInterpolation(size_t idx, G4double e) const
110 // Spline interpolation is used to get the value. Before this method
111 // is called it is ensured that the energy is inside the bin
112 // 0 < idx < numberOfNodes-1
114 static const G4double onesixth = 1.0/6.0;
117 G4double x1 = binVector[idx];
118 G4double x2 = binVector[idx + 1];
119 G4double delta = x2 - x1;
121 G4double a = (x2 - e)/delta;
122 G4double b = (e - x1)/delta;
124 // Final evaluation of cubic spline polynomial for return
125 G4double y1 = dataVector[idx];
126 G4double y2 = dataVector[idx + 1];
128 G4double res = a*y1 + b*y2 +
129 ( (a*a*a - a)*secDerivative[idx] +
130 (b*b*b - b)*secDerivative[idx + 1] )*delta*delta*onesixth;
135 //---------------------------------------------------------------
138 G4double G4PhysicsVector::Interpolation(size_t idx, G4double e) const
141 if(useSpline) { res = SplineInterpolation(idx, e); }
142 else { res = LinearInterpolation(idx, e); }
146 //---------------------------------------------------------------
149 void G4PhysicsVector::PutValue(size_t binNumber, G4double theValue)
151 dataVector[binNumber] = theValue;
154 //---------------------------------------------------------------
157 G4bool G4PhysicsVector::IsFilledVectorExist() const
161 if(numberOfNodes > 0) { status=true; }
165 //---------------------------------------------------------------
168 G4PhysicsVectorType G4PhysicsVector::GetType() const
173 //---------------------------------------------------------------
175 // Flag useSpline is "true" only if second derivatives are filled
177 void G4PhysicsVector::SetSpline(G4bool val)
180 if(0 == secDerivative.size() && 0 < dataVector.size()) {
181 FillSecondDerivatives();
185 secDerivative.clear();
189 //---------------------------------------------------------------
192 void G4PhysicsVector::SetVerboseLevel(G4int value)
194 verboseLevel = value;
197 //---------------------------------------------------------------
200 G4int G4PhysicsVector::GetVerboseLevel(G4int)
205 //---------------------------------------------------------------
208 size_t G4PhysicsVector::FindBinLocation(G4double theEnergy) const
211 if(type == T_G4PhysicsLogVector) {
212 bin = size_t(G4Log(theEnergy)/dBin - baseBin);
213 if(bin > 0 && theEnergy < binVector[bin]) { --bin; }
214 else if(theEnergy > binVector[bin+1]) { ++bin; }
215 } else if(type == T_G4PhysicsLinearVector) {
216 bin = size_t( theEnergy/dBin - baseBin );
217 if(bin > 0 && theEnergy < binVector[bin]) { --bin; }
218 else if(theEnergy > binVector[bin+1]) { ++bin; }
220 // Bin location proposed by K.Genser (FNAL)
221 bin = std::lower_bound(binVector.begin(), binVector.end(), theEnergy)
222 - binVector.begin() - 1;
224 return std::min(bin, numberOfNodes-2);
227 //---------------------------------------------------------------
229 inline size_t G4PhysicsVector::FindBin(G4double e, size_t idx) const
232 if(e < binVector[1]) {
234 } else if(e >= binVector[numberOfNodes-2]) {
235 id = numberOfNodes - 2;
236 } else if(idx >= numberOfNodes || e < binVector[idx]
237 || e > binVector[idx+1]) {
238 id = FindBinLocation(e);
243 //---------------------------------------------------------------
246 G4double G4PhysicsVector::Value(G4double theEnergy) const
249 return Value(theEnergy, idx);
252 //---------------------------------------------------------------