00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042 #include "G4VelocityTable.hh"
00043 #include "G4PhysicalConstants.hh"
00044 #include "G4StateManager.hh"
00045 #include "G4ApplicationState.hh"
00046
00047 #include "G4ios.hh"
00048
00049 G4VelocityTable* G4VelocityTable::theInstance = new G4VelocityTable();
00050
00052 G4VelocityTable::G4VelocityTable()
00054 : edgeMin(0.), edgeMax(0.), numberOfNodes(0),
00055 dBin(0.), baseBin(0.),
00056 lastEnergy(-DBL_MAX), lastValue(0.), lastBin(0),
00057 maxT( 1000.0 ), minT( 0.0001 ), NbinT( 500 )
00058 {
00059 PrepareVelocityTable();
00060 }
00061
00063 G4VelocityTable::~G4VelocityTable()
00065 {
00066 dataVector.clear();
00067 binVector.clear();
00068 }
00069
00070
00072 void G4VelocityTable::PrepareVelocityTable()
00074 {
00075 dataVector.clear();
00076 binVector.clear();
00077
00078 dBin = std::log10(maxT/minT)/NbinT;
00079 baseBin = std::log10(minT)/dBin;
00080
00081 numberOfNodes = NbinT + 1;
00082 dataVector.reserve(numberOfNodes);
00083 binVector.reserve(numberOfNodes);
00084
00085 const G4double g4log10 = std::log(10.);
00086
00087 binVector.push_back(minT);
00088 dataVector.push_back(0.0);
00089
00090 for (size_t i=1; i<numberOfNodes-1; i++){
00091 binVector.push_back(std::exp(g4log10*(baseBin+i)*dBin));
00092 dataVector.push_back(0.0);
00093 }
00094 binVector.push_back(maxT);
00095 dataVector.push_back(0.0);
00096
00097 edgeMin = binVector[0];
00098 edgeMax = binVector[numberOfNodes-1];
00099
00100 for (G4int i=0; i<=NbinT; i++){
00101 G4double T = binVector[i];
00102 dataVector[i]= c_light*std::sqrt(T*(T+2.))/(T+1.0);
00103 }
00104
00105 return;
00106 }
00107
00108 size_t G4VelocityTable::FindBinLocation(G4double theEnergy) const
00109 {
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120 return size_t( std::log10(theEnergy)/dBin - baseBin );
00121 }
00122
00123 G4double G4VelocityTable::Value(G4double theEnergy)
00124 {
00125
00126
00127
00128
00129
00130 if( theEnergy == lastEnergy ) {
00131
00132 } else if( theEnergy < lastEnergy
00133 && theEnergy >= binVector[lastBin]) {
00134 lastEnergy = theEnergy;
00135 lastValue = Interpolation();
00136
00137 } else if( theEnergy <= edgeMin ) {
00138 lastBin = 0;
00139 lastEnergy = edgeMin;
00140 lastValue = dataVector[0];
00141
00142 } else if( theEnergy >= edgeMax ){
00143 lastBin = numberOfNodes-1;
00144 lastEnergy = edgeMax;
00145 lastValue = dataVector[lastBin];
00146
00147 } else {
00148 lastBin = FindBinLocation(theEnergy);
00149 lastEnergy = theEnergy;
00150 lastValue = Interpolation();
00151
00152 }
00153 return lastValue;
00154 }
00155
00157
00158
00160 G4VelocityTable* G4VelocityTable::GetVelocityTable()
00162 {
00163 return theInstance;
00164 }
00165
00167 void G4VelocityTable::SetVelocityTableProperties(G4double t_max, G4double t_min, G4int nbin)
00169 {
00170 G4StateManager* stateManager = G4StateManager::GetStateManager();
00171 G4ApplicationState currentState = stateManager->GetCurrentState();
00172
00173
00174 if(!(currentState==G4State_Idle||currentState==G4State_PreInit)){
00175 G4Exception("G4VelocityTable::SetVelocityTableProperties",
00176 "Track101", JustWarning,
00177 "Can modify only in PreInit or Idle state : Method ignored.");
00178 return;
00179 }
00180
00181 if (nbin > 100 ) theInstance->NbinT = nbin;
00182 if ((t_min < t_max)&&(t_min>0.)) {
00183 theInstance->minT = t_min;
00184 theInstance->maxT = t_max;
00185 }
00186 theInstance->PrepareVelocityTable();
00187 }
00188
00190 G4double G4VelocityTable::GetMaxTOfVelocityTable()
00192 { return theInstance->maxT; }
00193
00195 G4double G4VelocityTable::GetMinTOfVelocityTable()
00196
00197 { return theInstance->minT; }
00198
00200 G4int G4VelocityTable::GetNbinOfVelocityTable()
00201
00202 { return theInstance->NbinT; }
00203