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

#include <G4PAIySection.hh>

Public Member Functions

 G4PAIySection ()
 
 ~G4PAIySection ()
 
void Initialize (const G4Material *material, G4double maxEnergyTransfer, G4double betaGammaSq, G4SandiaTable *)
 
void ComputeLowEnergyCof (const G4Material *material)
 
void InitPAI ()
 
void NormShift (G4double betaGammaSq)
 
void SplainPAI (G4double betaGammaSq)
 
G4double RutherfordIntegral (G4int intervalNumber, G4double limitLow, G4double limitHigh)
 
G4double ImPartDielectricConst (G4int intervalNumber, G4double energy)
 
G4double RePartDielectricConst (G4double energy)
 
G4double DifPAIySection (G4int intervalNumber, G4double betaGammaSq)
 
G4double PAIdNdxCerenkov (G4int intervalNumber, G4double betaGammaSq)
 
G4double PAIdNdxPlasmon (G4int intervalNumber, G4double betaGammaSq)
 
void IntegralPAIySection ()
 
void IntegralCerenkov ()
 
void IntegralPlasmon ()
 
G4double SumOverInterval (G4int intervalNumber)
 
G4double SumOverIntervaldEdx (G4int intervalNumber)
 
G4double SumOverInterCerenkov (G4int intervalNumber)
 
G4double SumOverInterPlasmon (G4int intervalNumber)
 
G4double SumOverBorder (G4int intervalNumber, G4double energy)
 
G4double SumOverBorderdEdx (G4int intervalNumber, G4double energy)
 
G4double SumOverBordCerenkov (G4int intervalNumber, G4double energy)
 
G4double SumOverBordPlasmon (G4int intervalNumber, G4double energy)
 
G4double GetStepEnergyLoss (G4double step)
 
G4double GetStepCerenkovLoss (G4double step)
 
G4double GetStepPlasmonLoss (G4double step)
 
G4int GetNumberOfGammas () const
 
G4int GetSplineSize () const
 
G4int GetIntervalNumber () const
 
G4double GetEnergyInterval (G4int i)
 
G4double GetDifPAIySection (G4int i)
 
G4double GetPAIdNdxCrenkov (G4int i)
 
G4double GetPAIdNdxPlasmon (G4int i)
 
G4double GetMeanEnergyLoss () const
 
G4double GetMeanCerenkovLoss () const
 
G4double GetMeanPlasmonLoss () const
 
G4double GetNormalizationCof () const
 
G4double GetPAItable (G4int i, G4int j) const
 
G4double GetLorentzFactor (G4int i) const
 
G4double GetSplineEnergy (G4int i) const
 
G4double GetIntegralPAIySection (G4int i) const
 
G4double GetIntegralPAIdEdx (G4int i) const
 
G4double GetIntegralCerenkov (G4int i) const
 
G4double GetIntegralPlasmon (G4int i) const
 
void SetVerbose (G4int v)
 

Detailed Description

Definition at line 52 of file G4PAIySection.hh.

Constructor & Destructor Documentation

G4PAIySection::G4PAIySection ( )
explicit

Definition at line 76 of file G4PAIySection.cc.

77 {
78  fSandia = 0;
79  fDensity = fElectronDensity = fNormalizationCof = fLowEnergyCof = 0.0;
80  fIntervalNumber = fSplineNumber = 0;
81  fVerbose = 0;
82 
83  betaBohr = fine_structure_const;
84  G4double cofBetaBohr = 4.0;
86  betaBohr4 = betaBohr2*betaBohr2*cofBetaBohr;
87 
88  fSplineEnergy = G4DataVector(fMaxSplineSize,0.0);
89  fRePartDielectricConst = G4DataVector(fMaxSplineSize,0.0);
90  fImPartDielectricConst = G4DataVector(fMaxSplineSize,0.0);
91  fIntegralTerm = G4DataVector(fMaxSplineSize,0.0);
92  fDifPAIySection = G4DataVector(fMaxSplineSize,0.0);
93  fdNdxCerenkov = G4DataVector(fMaxSplineSize,0.0);
94  fdNdxPlasmon = G4DataVector(fMaxSplineSize,0.0);
95  fIntegralPAIySection = G4DataVector(fMaxSplineSize,0.0);
96  fIntegralPAIdEdx = G4DataVector(fMaxSplineSize,0.0);
97  fIntegralCerenkov = G4DataVector(fMaxSplineSize,0.0);
98  fIntegralPlasmon = G4DataVector(fMaxSplineSize,0.0);
99 
100  for( G4int i = 0; i < 500; ++i )
101  {
102  for( G4int j = 0; j < 112; ++j ) { fPAItable[i][j] = 0.0; }
103  }
104 }
int G4int
Definition: G4Types.hh:78
double G4double
Definition: G4Types.hh:76
G4PAIySection::~G4PAIySection ( )

Definition at line 110 of file G4PAIySection.cc.

111 {}

Member Function Documentation

void G4PAIySection::ComputeLowEnergyCof ( const G4Material material)

Definition at line 244 of file G4PAIySection.cc.

245 {
246  G4int i, numberOfElements = material->GetNumberOfElements();
247  G4double sumZ = 0., sumCof = 0.;
248 
249  static const G4double p0 = 1.20923e+00;
250  static const G4double p1 = 3.53256e-01;
251  static const G4double p2 = -1.45052e-03;
252 
253  G4double* thisMaterialZ = new G4double[numberOfElements];
254  G4double* thisMaterialCof = new G4double[numberOfElements];
255 
256  for( i = 0; i < numberOfElements; i++ )
257  {
258  thisMaterialZ[i] = material->GetElement(i)->GetZ();
259  sumZ += thisMaterialZ[i];
260  thisMaterialCof[i] = p0+p1*thisMaterialZ[i]+p2*thisMaterialZ[i]*thisMaterialZ[i];
261  }
262  for( i = 0; i < numberOfElements; i++ )
263  {
264  sumCof += thisMaterialCof[i]*thisMaterialZ[i]/sumZ;
265  }
266  fLowEnergyCof = sumCof;
267  delete [] thisMaterialZ;
268  delete [] thisMaterialCof;
269  // G4cout<<"fLowEnergyCof = "<<fLowEnergyCof<<G4endl;
270 }
G4double GetZ() const
Definition: G4Element.hh:131
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:202
int G4int
Definition: G4Types.hh:78
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4PAIySection::DifPAIySection ( G4int  intervalNumber,
G4double  betaGammaSq 
)

Definition at line 588 of file G4PAIySection.cc.

590 {
591  G4double beta, be2,cof,x1,x2,x3,x4,x5,x6,x7,x8,result;
592  //G4double beta, be4;
593  //G4double be4;
594  // G4double betaBohr2 = fine_structure_const*fine_structure_const;
595  // G4double betaBohr4 = betaBohr2*betaBohr2*4.0;
596  be2 = betaGammaSq/(1 + betaGammaSq);
597  //be4 = be2*be2;
598  beta = sqrt(be2);
599  cof = 1;
600  x1 = log(2*electron_mass_c2/fSplineEnergy[i]);
601 
602  if( betaGammaSq < 0.01 ) x2 = log(be2);
603  else
604  {
605  x2 = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
606  (1/betaGammaSq - fRePartDielectricConst[i]) +
607  fImPartDielectricConst[i]*fImPartDielectricConst[i] )/2;
608  }
609  if( fImPartDielectricConst[i] == 0.0 ||betaGammaSq < 0.01 )
610  {
611  x6=0;
612  }
613  else
614  {
615  x3 = -fRePartDielectricConst[i] + 1/betaGammaSq;
616  x5 = -1 - fRePartDielectricConst[i] +
617  be2*((1 +fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
618  fImPartDielectricConst[i]*fImPartDielectricConst[i]);
619 
620  x7 = atan2(fImPartDielectricConst[i],x3);
621  x6 = x5 * x7;
622  }
623  // if(fImPartDielectricConst[i] == 0) x6 = 0;
624 
625  x4 = ((x1 + x2)*fImPartDielectricConst[i] + x6)/hbarc;
626  // if( x4 < 0.0 ) x4 = 0.0;
627  x8 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
628  fImPartDielectricConst[i]*fImPartDielectricConst[i];
629 
630  result = (x4 + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i]);
631  result = std::max(result, 1.0e-8);
632  result *= fine_structure_const/be2/pi;
633  // low energy correction
634 
635  G4double lowCof = fLowEnergyCof; // 6.0 ; // Ar ~ 4.; -> fLowCof as f(Z1,Z2)?
636 
637  result *= (1 - exp(-beta/betaBohr/lowCof));
638 
639  // result *= (1-exp(-beta/betaBohr))*(1-exp(-beta/betaBohr));
640  // result *= (1-exp(-be2/betaBohr2));
641  // result *= (1-exp(-be4/betaBohr4));
642  // if(fDensity >= 0.1)
643  if(x8 > 0.)
644  {
645  result /= x8;
646  }
647  return result;
648 
649 } // end of DifPAIySection
G4double G4ParticleHPJENDLHEData::G4double result
float electron_mass_c2
Definition: hepunit.py:274
T max(const T t1, const T t2)
brief Return the largest of the two arguments
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4PAIySection::GetDifPAIySection ( G4int  i)
inline

Definition at line 122 of file G4PAIySection.hh.

122 { return fDifPAIySection[i]; }
G4double G4PAIySection::GetEnergyInterval ( G4int  i)
inline

Definition at line 120 of file G4PAIySection.hh.

120 { return fEnergyInterval[i]; }
G4double G4PAIySection::GetIntegralCerenkov ( G4int  i) const
inline

Definition at line 227 of file G4PAIySection.hh.

228 {
229  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralCerenkov"); }
230  return fIntegralCerenkov[i];
231 }
G4double G4PAIySection::GetIntegralPAIdEdx ( G4int  i) const
inline

Definition at line 221 of file G4PAIySection.hh.

222 {
223  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIdEdx"); }
224  return fIntegralPAIdEdx[i];
225 }
G4double G4PAIySection::GetIntegralPAIySection ( G4int  i) const
inline

Definition at line 215 of file G4PAIySection.hh.

216 {
217  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPAIySection"); }
218  return fIntegralPAIySection[i];
219 }
G4double G4PAIySection::GetIntegralPlasmon ( G4int  i) const
inline

Definition at line 233 of file G4PAIySection.hh.

234 {
235  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPlasmon"); }
236  return fIntegralPlasmon[i];
237 }
G4int G4PAIySection::GetIntervalNumber ( ) const
inline

Definition at line 118 of file G4PAIySection.hh.

118 { return fIntervalNumber; }
G4double G4PAIySection::GetLorentzFactor ( G4int  i) const
inline

Definition at line 204 of file G4PAIySection.hh.

205 {
206  return fLorentzFactor[j];
207 }
G4double G4PAIySection::GetMeanCerenkovLoss ( ) const
inline

Definition at line 127 of file G4PAIySection.hh.

127 {return fIntegralCerenkov[0]; }
G4double G4PAIySection::GetMeanEnergyLoss ( ) const
inline

Definition at line 126 of file G4PAIySection.hh.

126 {return fIntegralPAIySection[0]; }
G4double G4PAIySection::GetMeanPlasmonLoss ( ) const
inline

Definition at line 128 of file G4PAIySection.hh.

128 {return fIntegralPlasmon[0]; }
G4double G4PAIySection::GetNormalizationCof ( ) const
inline

Definition at line 130 of file G4PAIySection.hh.

130 { return fNormalizationCof; }
G4int G4PAIySection::GetNumberOfGammas ( ) const
inline

Definition at line 114 of file G4PAIySection.hh.

114 { return fNumberOfGammas; }
G4double G4PAIySection::GetPAIdNdxCrenkov ( G4int  i)
inline

Definition at line 123 of file G4PAIySection.hh.

123 { return fdNdxCerenkov[i]; }
G4double G4PAIySection::GetPAIdNdxPlasmon ( G4int  i)
inline

Definition at line 124 of file G4PAIySection.hh.

124 { return fdNdxPlasmon[i]; }
G4double G4PAIySection::GetPAItable ( G4int  i,
G4int  j 
) const
inline

Definition at line 199 of file G4PAIySection.hh.

200 {
201  return fPAItable[i][j];
202 }
G4double G4PAIySection::GetSplineEnergy ( G4int  i) const
inline

Definition at line 209 of file G4PAIySection.hh.

210 {
211  if(i < 1 || i > fSplineNumber) { CallError(i, "GetSplineEnergy"); }
212  return fSplineEnergy[i];
213 }
G4int G4PAIySection::GetSplineSize ( ) const
inline

Definition at line 116 of file G4PAIySection.hh.

116 { return fSplineNumber; }
G4double G4PAIySection::GetStepCerenkovLoss ( G4double  step)

Definition at line 1274 of file G4PAIySection.cc.

1275 {
1276  G4int iTransfer ;
1277  G4long numOfCollisions;
1278  G4double loss = 0.0;
1279  G4double meanNumber, position;
1280 
1281  // G4cout<<" G4PAIySection::GetStepCreLosnkovs "<<G4endl;
1282 
1283 
1284 
1285  meanNumber = fIntegralCerenkov[1]*step;
1286  numOfCollisions = G4Poisson(meanNumber);
1287 
1288  // G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
1289 
1290  while(numOfCollisions)
1291  {
1292  position = fIntegralCerenkov[1]*G4UniformRand();
1293 
1294  for( iTransfer=1; iTransfer<=fSplineNumber; iTransfer++ )
1295  {
1296  if( position >= fIntegralCerenkov[iTransfer] ) break;
1297  }
1298  loss += fSplineEnergy[iTransfer] ;
1299  numOfCollisions--;
1300  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
1301  }
1302  // G4cout<<"PAI Cerenkov loss = "<<loss/keV<<" keV"<<G4endl;
1303 
1304  return loss;
1305 }
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:51
long G4long
Definition: G4Types.hh:80
int G4int
Definition: G4Types.hh:78
#define position
Definition: xmlparse.cc:622
#define G4UniformRand()
Definition: Randomize.hh:97
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4PAIySection::GetStepEnergyLoss ( G4double  step)

Definition at line 1237 of file G4PAIySection.cc.

1238 {
1239  G4int iTransfer ;
1240  G4long numOfCollisions;
1241  G4double loss = 0.0;
1242  G4double meanNumber, position;
1243 
1244  // G4cout<<" G4PAIySection::GetStepEnergyLoss "<<G4endl;
1245 
1246 
1247 
1248  meanNumber = fIntegralPAIySection[1]*step;
1249  numOfCollisions = G4Poisson(meanNumber);
1250 
1251  // G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
1252 
1253  while(numOfCollisions)
1254  {
1255  position = fIntegralPAIySection[1]*G4UniformRand();
1256 
1257  for( iTransfer=1; iTransfer<=fSplineNumber; iTransfer++ )
1258  {
1259  if( position >= fIntegralPAIySection[iTransfer] ) break;
1260  }
1261  loss += fSplineEnergy[iTransfer] ;
1262  numOfCollisions--;
1263  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
1264  }
1265  // G4cout<<"PAI energy loss = "<<loss/keV<<" keV"<<G4endl;
1266 
1267  return loss;
1268 }
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:51
long G4long
Definition: G4Types.hh:80
int G4int
Definition: G4Types.hh:78
#define position
Definition: xmlparse.cc:622
#define G4UniformRand()
Definition: Randomize.hh:97
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4PAIySection::GetStepPlasmonLoss ( G4double  step)

Definition at line 1311 of file G4PAIySection.cc.

1312 {
1313  G4int iTransfer ;
1314  G4long numOfCollisions;
1315  G4double loss = 0.0;
1316  G4double meanNumber, position;
1317 
1318  // G4cout<<" G4PAIySection::GetStepCreLosnkovs "<<G4endl;
1319 
1320 
1321 
1322  meanNumber = fIntegralPlasmon[1]*step;
1323  numOfCollisions = G4Poisson(meanNumber);
1324 
1325  // G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
1326 
1327  while(numOfCollisions)
1328  {
1329  position = fIntegralPlasmon[1]*G4UniformRand();
1330 
1331  for( iTransfer=1; iTransfer<=fSplineNumber; iTransfer++ )
1332  {
1333  if( position >= fIntegralPlasmon[iTransfer] ) break;
1334  }
1335  loss += fSplineEnergy[iTransfer] ;
1336  numOfCollisions--;
1337  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
1338  }
1339  // G4cout<<"PAI Plasmon loss = "<<loss/keV<<" keV"<<G4endl;
1340 
1341  return loss;
1342 }
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:51
long G4long
Definition: G4Types.hh:80
int G4int
Definition: G4Types.hh:78
#define position
Definition: xmlparse.cc:622
#define G4UniformRand()
Definition: Randomize.hh:97
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4PAIySection::ImPartDielectricConst ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 510 of file G4PAIySection.cc.

512 {
513  G4double energy2,energy3,energy4,result;
514 
515  energy2 = energy1*energy1;
516  energy3 = energy2*energy1;
517  energy4 = energy3*energy1;
518 
519  result = fA1[k]/energy1+fA2[k]/energy2+fA3[k]/energy3+fA4[k]/energy4;
520  result *=hbarc/energy1;
521 
522  return result;
523 
524 } // end of ImPartDielectricConst
G4double G4ParticleHPJENDLHEData::G4double result
double G4double
Definition: G4Types.hh:76
void G4PAIySection::Initialize ( const G4Material material,
G4double  maxEnergyTransfer,
G4double  betaGammaSq,
G4SandiaTable sandia 
)

Definition at line 117 of file G4PAIySection.cc.

121 {
122  if(fVerbose > 0)
123  {
124  G4cout<<G4endl;
125  G4cout<<"G4PAIySection::Initialize(...,G4SandiaTable* sandia)"<<G4endl;
126  G4cout<<G4endl;
127  }
128  G4int i, j;
129 
130  fSandia = sandia;
131  fIntervalNumber = sandia->GetMaxInterval();
132  fDensity = material->GetDensity();
133  fElectronDensity = material->GetElectronDensity();
134 
135  // fIntervalNumber--;
136 
137  if( fVerbose > 0 )
138  {
139  G4cout<<"fDensity = "<<fDensity<<"\t"<<fElectronDensity<<"\t fIntervalNumber = "
140  <<fIntervalNumber<< " (beta*gamma)^2= " << betaGammaSq << G4endl;
141  }
142  fEnergyInterval = G4DataVector(fIntervalNumber+2,0.0);
143  fA1 = G4DataVector(fIntervalNumber+2,0.0);
144  fA2 = G4DataVector(fIntervalNumber+2,0.0);
145  fA3 = G4DataVector(fIntervalNumber+2,0.0);
146  fA4 = G4DataVector(fIntervalNumber+2,0.0);
147 
148  for( i = 1; i <= fIntervalNumber; i++ )
149  {
150  if ( sandia->GetSandiaMatTablePAI(i-1,0) < 1.*eV )
151  {
152  fIntervalNumber--;
153  continue;
154  }
155  if( ( sandia->GetSandiaMatTablePAI(i-1,0) >= maxEnergyTransfer )
156  || i >= fIntervalNumber )
157  {
158  fEnergyInterval[i] = maxEnergyTransfer;
159  fIntervalNumber = i;
160  break;
161  }
162  fEnergyInterval[i] = sandia->GetSandiaMatTablePAI(i-1,0);
163  fA1[i] = sandia->GetSandiaMatTablePAI(i-1,1);
164  fA2[i] = sandia->GetSandiaMatTablePAI(i-1,2);
165  fA3[i] = sandia->GetSandiaMatTablePAI(i-1,3);
166  fA4[i] = sandia->GetSandiaMatTablePAI(i-1,4);
167 
168  if( fVerbose > 0 ) {
169  G4cout<<i<<"\t"<<fEnergyInterval[i]/keV<<"\t"<<fA1[i]<<"\t"<<fA2[i]<<"\t"
170  <<fA3[i]<<"\t"<<fA4[i]<<"\t"<<G4endl;
171  }
172  }
173  if( fVerbose > 0 ) {
174  G4cout<<"last i = "<<i<<"; "<<"fIntervalNumber = "
175  <<fIntervalNumber<<G4endl;
176  }
177  if( fEnergyInterval[fIntervalNumber] != maxEnergyTransfer )
178  {
179  fIntervalNumber++;
180  fEnergyInterval[fIntervalNumber] = maxEnergyTransfer;
181  }
182  if( fVerbose > 0 )
183  {
184  for( i = 1; i <= fIntervalNumber; i++ )
185  {
186  G4cout<<i<<"\t"<<fEnergyInterval[i]/keV<<"\t"<<fA1[i]<<"\t"<<fA2[i]<<"\t"
187  <<fA3[i]<<"\t"<<fA4[i]<<"\t"<<G4endl;
188  }
189  }
190  if( fVerbose > 0 ) {
191  G4cout<<"Now checking, if two borders are too close together"<<G4endl;
192  }
193  for( i = 1; i < fIntervalNumber; i++ )
194  {
195  if( fEnergyInterval[i+1]-fEnergyInterval[i] >
196  1.5*fDelta*(fEnergyInterval[i+1]+fEnergyInterval[i]) ) continue;
197  else
198  {
199  for( j = i; j < fIntervalNumber; j++ )
200  {
201  fEnergyInterval[j] = fEnergyInterval[j+1];
202  fA1[j] = fA1[j+1];
203  fA2[j] = fA2[j+1];
204  fA3[j] = fA3[j+1];
205  fA4[j] = fA4[j+1];
206  }
207  fIntervalNumber--;
208  }
209  }
210  if( fVerbose > 0 )
211  {
212  for( i = 1; i <= fIntervalNumber; i++ )
213  {
214  G4cout<<i<<"\t"<<fEnergyInterval[i]/keV<<"\t"<<fA1[i]<<"\t"<<fA2[i]<<"\t"
215  <<fA3[i]<<"\t"<<fA4[i]<<"\t"<<G4endl;
216  }
217  }
218  // Preparation of fSplineEnergy array corresponding to min ionisation, G~4
219 
220  ComputeLowEnergyCof(material);
221 
222  G4double betaGammaSqRef =
223  fLorentzFactor[fRefGammaNumber]*fLorentzFactor[fRefGammaNumber] - 1;
224 
225  NormShift(betaGammaSqRef);
226  SplainPAI(betaGammaSqRef);
227 
228  // Preparation of integral PAI cross section for input betaGammaSq
229 
230  for( i = 1; i <= fSplineNumber; i++ )
231  {
232  fDifPAIySection[i] = DifPAIySection(i,betaGammaSq);
233 
234  if( fVerbose > 0 ) G4cout<<i<<"; dNdxPAI = "<<fDifPAIySection[i]<<G4endl;
235  }
237 }
void NormShift(G4double betaGammaSq)
G4double GetDensity() const
Definition: G4Material.hh:180
G4int GetMaxInterval() const
G4double GetSandiaMatTablePAI(G4int, G4int) const
void ComputeLowEnergyCof(const G4Material *material)
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
G4double GetElectronDensity() const
Definition: G4Material.hh:217
static constexpr double eV
Definition: G4SIunits.hh:215
void SplainPAI(G4double betaGammaSq)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4double DifPAIySection(G4int intervalNumber, G4double betaGammaSq)
static constexpr double keV
Definition: G4SIunits.hh:216
void IntegralPAIySection()

Here is the call graph for this function:

void G4PAIySection::InitPAI ( )

Definition at line 277 of file G4PAIySection.cc.

278 {
279  G4int i;
280  G4double betaGammaSq = fLorentzFactor[fRefGammaNumber]*
281  fLorentzFactor[fRefGammaNumber] - 1;
282 
283  // Preparation of integral PAI cross section for reference gamma
284 
285  NormShift(betaGammaSq);
286  SplainPAI(betaGammaSq);
287 
290  IntegralPlasmon();
291 
292  for( i = 0; i<= fSplineNumber; i++)
293  {
294  fPAItable[i][fRefGammaNumber] = fIntegralPAIySection[i];
295 
296  if(i != 0) fPAItable[i][0] = fSplineEnergy[i];
297  }
298  fPAItable[0][0] = fSplineNumber;
299 
300  for( G4int j = 1; j < 112; j++) // for other gammas
301  {
302  if( j == fRefGammaNumber ) continue;
303 
304  betaGammaSq = fLorentzFactor[j]*fLorentzFactor[j] - 1;
305 
306  for(i = 1; i <= fSplineNumber; i++)
307  {
308  fDifPAIySection[i] = DifPAIySection(i,betaGammaSq);
309  fdNdxCerenkov[i] = PAIdNdxCerenkov(i,betaGammaSq);
310  fdNdxPlasmon[i] = PAIdNdxPlasmon(i,betaGammaSq);
311  }
314  IntegralPlasmon();
315 
316  for(i = 0; i <= fSplineNumber; i++)
317  {
318  fPAItable[i][j] = fIntegralPAIySection[i];
319  }
320  }
321 }
void NormShift(G4double betaGammaSq)
G4double PAIdNdxCerenkov(G4int intervalNumber, G4double betaGammaSq)
void IntegralCerenkov()
int G4int
Definition: G4Types.hh:78
void IntegralPlasmon()
void SplainPAI(G4double betaGammaSq)
double G4double
Definition: G4Types.hh:76
G4double DifPAIySection(G4int intervalNumber, G4double betaGammaSq)
G4double PAIdNdxPlasmon(G4int intervalNumber, G4double betaGammaSq)
void IntegralPAIySection()
void G4PAIySection::IntegralCerenkov ( )

Definition at line 784 of file G4PAIySection.cc.

785 {
786  G4int i, k;
787  fIntegralCerenkov[fSplineNumber] = 0;
788  fIntegralCerenkov[0] = 0;
789  k = fIntervalNumber -1;
790 
791  for( i = fSplineNumber-1; i >= 1; i-- )
792  {
793  if(fSplineEnergy[i] >= fEnergyInterval[k])
794  {
795  fIntegralCerenkov[i] = fIntegralCerenkov[i+1] + SumOverInterCerenkov(i);
796  // G4cout<<"int: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
797  }
798  else
799  {
800  fIntegralCerenkov[i] = fIntegralCerenkov[i+1] +
801  SumOverBordCerenkov(i+1,fEnergyInterval[k]);
802  k--;
803  // G4cout<<"bord: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
804  }
805  }
806 
807 } // end of IntegralCerenkov
G4double SumOverInterCerenkov(G4int intervalNumber)
int G4int
Definition: G4Types.hh:78
G4double SumOverBordCerenkov(G4int intervalNumber, G4double energy)
void G4PAIySection::IntegralPAIySection ( )

Definition at line 753 of file G4PAIySection.cc.

754 {
755  fIntegralPAIySection[fSplineNumber] = 0;
756  fIntegralPAIdEdx[fSplineNumber] = 0;
757  fIntegralPAIySection[0] = 0;
758  G4int k = fIntervalNumber -1;
759 
760  for(G4int i = fSplineNumber-1; i >= 1; i--)
761  {
762  if(fSplineEnergy[i] >= fEnergyInterval[k])
763  {
764  fIntegralPAIySection[i] = fIntegralPAIySection[i+1] + SumOverInterval(i);
765  fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] + SumOverIntervaldEdx(i);
766  }
767  else
768  {
769  fIntegralPAIySection[i] = fIntegralPAIySection[i+1] +
770  SumOverBorder(i+1,fEnergyInterval[k]);
771  fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] +
772  SumOverBorderdEdx(i+1,fEnergyInterval[k]);
773  k--;
774  }
775  }
776 } // end of IntegralPAIySection
G4double SumOverIntervaldEdx(G4int intervalNumber)
G4double SumOverBorderdEdx(G4int intervalNumber, G4double energy)
int G4int
Definition: G4Types.hh:78
G4double SumOverBorder(G4int intervalNumber, G4double energy)
G4double SumOverInterval(G4int intervalNumber)
void G4PAIySection::IntegralPlasmon ( )

Definition at line 815 of file G4PAIySection.cc.

816 {
817  fIntegralPlasmon[fSplineNumber] = 0;
818  fIntegralPlasmon[0] = 0;
819  G4int k = fIntervalNumber -1;
820  for(G4int i=fSplineNumber-1;i>=1;i--)
821  {
822  if(fSplineEnergy[i] >= fEnergyInterval[k])
823  {
824  fIntegralPlasmon[i] = fIntegralPlasmon[i+1] + SumOverInterPlasmon(i);
825  }
826  else
827  {
828  fIntegralPlasmon[i] = fIntegralPlasmon[i+1] +
829  SumOverBordPlasmon(i+1,fEnergyInterval[k]);
830  k--;
831  }
832  }
833 
834 } // end of IntegralPlasmon
int G4int
Definition: G4Types.hh:78
G4double SumOverBordPlasmon(G4int intervalNumber, G4double energy)
G4double SumOverInterPlasmon(G4int intervalNumber)
void G4PAIySection::NormShift ( G4double  betaGammaSq)

Definition at line 328 of file G4PAIySection.cc.

329 {
330  G4int i, j;
331 
332  for( i = 1; i <= fIntervalNumber-1; i++ )
333  {
334  for( j = 1; j <= 2; j++ )
335  {
336  fSplineNumber = (i-1)*2 + j;
337 
338  if( j == 1 ) fSplineEnergy[fSplineNumber] = fEnergyInterval[i ]*(1+fDelta);
339  else fSplineEnergy[fSplineNumber] = fEnergyInterval[i+1]*(1-fDelta);
340  // G4cout<<"cn = "<<fSplineNumber<<"; "<<"energy = "
341  // <<fSplineEnergy[fSplineNumber]<<G4endl;
342  }
343  }
344  fIntegralTerm[1]=RutherfordIntegral(1,fEnergyInterval[1],fSplineEnergy[1]);
345 
346  j = 1;
347 
348  for(i=2;i<=fSplineNumber;i++)
349  {
350  if(fSplineEnergy[i]<fEnergyInterval[j+1])
351  {
352  fIntegralTerm[i] = fIntegralTerm[i-1] +
353  RutherfordIntegral(j,fSplineEnergy[i-1],
354  fSplineEnergy[i] );
355  }
356  else
357  {
358  G4double x = RutherfordIntegral(j,fSplineEnergy[i-1],
359  fEnergyInterval[j+1] );
360  j++;
361  fIntegralTerm[i] = fIntegralTerm[i-1] + x +
362  RutherfordIntegral(j,fEnergyInterval[j],
363  fSplineEnergy[i] );
364  }
365  // G4cout<<i<<"\t"<<fSplineEnergy[i]<<"\t"<<fIntegralTerm[i]<<"\n"<<G4endl;
366  }
367  static const G4double nfactor =
369  fNormalizationCof = nfactor*fElectronDensity/fIntegralTerm[fSplineNumber];
370 
371  // G4cout<<"fNormalizationCof = "<<fNormalizationCof<<G4endl;
372 
373  // Calculation of PAI differrential cross-section (1/(keV*cm))
374  // in the energy points near borders of energy intervals
375 
376  for(G4int k=1;k<=fIntervalNumber-1;k++)
377  {
378  for(j=1;j<=2;j++)
379  {
380  i = (k-1)*2 + j;
381  fImPartDielectricConst[i] = fNormalizationCof*
382  ImPartDielectricConst(k,fSplineEnergy[i]);
383  fRePartDielectricConst[i] = fNormalizationCof*
384  RePartDielectricConst(fSplineEnergy[i]);
385  fIntegralTerm[i] *= fNormalizationCof;
386 
387  fDifPAIySection[i] = DifPAIySection(i,betaGammaSq);
388  fdNdxCerenkov[i] = PAIdNdxCerenkov(i,betaGammaSq);
389  fdNdxPlasmon[i] = PAIdNdxPlasmon(i,betaGammaSq);
390  }
391  }
392 
393 } // end of NormShift
G4double RutherfordIntegral(G4int intervalNumber, G4double limitLow, G4double limitHigh)
G4double ImPartDielectricConst(G4int intervalNumber, G4double energy)
G4double PAIdNdxCerenkov(G4int intervalNumber, G4double betaGammaSq)
tuple x
Definition: test.py:50
int G4int
Definition: G4Types.hh:78
float electron_mass_c2
Definition: hepunit.py:274
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
G4double RePartDielectricConst(G4double energy)
G4double DifPAIySection(G4int intervalNumber, G4double betaGammaSq)
G4double PAIdNdxPlasmon(G4int intervalNumber, G4double betaGammaSq)
G4double G4PAIySection::PAIdNdxCerenkov ( G4int  intervalNumber,
G4double  betaGammaSq 
)

Definition at line 655 of file G4PAIySection.cc.

657 {
658  G4double logarithm, x3, x5, argument, modul2, dNdxC;
659  G4double be2, be4;
660 
661  //G4double cof = 1.0;
662 
663  be2 = betaGammaSq/(1 + betaGammaSq);
664  be4 = be2*be2;
665 
666  if( betaGammaSq < 0.01 ) logarithm = log(1.0+betaGammaSq); // 0.0;
667  else
668  {
669  logarithm = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
670  (1/betaGammaSq - fRePartDielectricConst[i]) +
671  fImPartDielectricConst[i]*fImPartDielectricConst[i] )*0.5;
672  logarithm += log(1+1.0/betaGammaSq);
673  }
674 
675  if( fImPartDielectricConst[i] == 0.0 || betaGammaSq < 0.01 )
676  {
677  argument = 0.0;
678  }
679  else
680  {
681  x3 = -fRePartDielectricConst[i] + 1.0/betaGammaSq;
682  x5 = -1.0 - fRePartDielectricConst[i] +
683  be2*((1.0 +fRePartDielectricConst[i])*(1.0 + fRePartDielectricConst[i]) +
684  fImPartDielectricConst[i]*fImPartDielectricConst[i]);
685  if( x3 == 0.0 ) argument = 0.5*pi;
686  else argument = atan2(fImPartDielectricConst[i],x3);
687  argument *= x5 ;
688  }
689  dNdxC = ( logarithm*fImPartDielectricConst[i] + argument )/hbarc;
690 
691  if(dNdxC < 1.0e-8) dNdxC = 1.0e-8;
692 
693  dNdxC *= fine_structure_const/be2/pi;
694 
695  dNdxC *= (1-exp(-be4/betaBohr4));
696 
697  // if(fDensity >= 0.1)
698  // {
699  modul2 = (1.0 + fRePartDielectricConst[i])*(1.0 + fRePartDielectricConst[i]) +
700  fImPartDielectricConst[i]*fImPartDielectricConst[i];
701  if(modul2 > 0.)
702  {
703  dNdxC /= modul2;
704  }
705  return dNdxC;
706 
707 } // end of PAIdNdxCerenkov
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::PAIdNdxPlasmon ( G4int  intervalNumber,
G4double  betaGammaSq 
)

Definition at line 714 of file G4PAIySection.cc.

716 {
717  G4double cof, resonance, modul2, dNdxP;
718  G4double be2, be4;
719 
720  cof = 1;
721 
722  be2 = betaGammaSq/(1 + betaGammaSq);
723  be4 = be2*be2;
724 
725  resonance = log(2*electron_mass_c2*be2/fSplineEnergy[i]);
726  resonance *= fImPartDielectricConst[i]/hbarc;
727 
728  dNdxP = ( resonance + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i] );
729 
730  dNdxP = std::max(dNdxP, 1.0e-8);
731 
732  dNdxP *= fine_structure_const/be2/pi;
733  dNdxP *= (1-exp(-be4/betaBohr4));
734 
735 // if( fDensity >= 0.1 )
736 // {
737  modul2 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
738  fImPartDielectricConst[i]*fImPartDielectricConst[i];
739  if(modul2 > 0.)
740  {
741  dNdxP /= modul2;
742  }
743  return dNdxP;
744 
745 } // end of PAIdNdxPlasmon
float electron_mass_c2
Definition: hepunit.py:274
T max(const T t1, const T t2)
brief Return the largest of the two arguments
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4PAIySection::RePartDielectricConst ( G4double  energy)

Definition at line 533 of file G4PAIySection.cc.

534 {
535  G4double x0, x02, x03, x04, x05, x1, x2, xx1 ,xx2 , xx12,
536  c1, c2, c3, cof1, cof2, xln1, xln2, xln3, result;
537 
538  x0 = enb;
539  result = 0;
540 
541  for(G4int i=1;i<=fIntervalNumber-1;i++)
542  {
543  x1 = fEnergyInterval[i];
544  x2 = fEnergyInterval[i+1];
545  xx1 = x1 - x0;
546  xx2 = x2 - x0;
547  xx12 = xx2/xx1;
548 
549  if(xx12<0)
550  {
551  xx12 = -xx12;
552  }
553  xln1 = log(x2/x1);
554  xln2 = log(xx12);
555  xln3 = log((x2 + x0)/(x1 + x0));
556  x02 = x0*x0;
557  x03 = x02*x0;
558  x04 = x03*x0;
559  x05 = x04*x0;
560  G4double x12 = x1*x2;
561  c1 = (x2 - x1)/x12;
562  c2 = (x2 - x1)*(x2 +x1)/(x12*x12);
563  c3 = (x2 -x1)*(x1*x1 + x1*x2 + x2*x2)/(x12*x12*x12);
564 
565  result -= (fA1[i]/x02 + fA3[i]/x04)*xln1;
566  result -= (fA2[i]/x02 + fA4[i]/x04)*c1;
567  result -= fA3[i]*c2/2/x02;
568  result -= fA4[i]*c3/3/x02;
569 
570  cof1 = fA1[i]/x02 + fA3[i]/x04;
571  cof2 = fA2[i]/x03 + fA4[i]/x05;
572 
573  result += 0.5*(cof1 +cof2)*xln2;
574  result += 0.5*(cof1 - cof2)*xln3;
575  }
576  result *= 2*hbarc/pi;
577 
578  return result;
579 
580 } // end of RePartDielectricConst
G4double G4ParticleHPJENDLHEData::G4double result
static c2_factory< G4double > c2
int G4int
Definition: G4Types.hh:78
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
tuple c1
Definition: plottest35.py:14
G4double G4PAIySection::RutherfordIntegral ( G4int  intervalNumber,
G4double  limitLow,
G4double  limitHigh 
)

Definition at line 488 of file G4PAIySection.cc.

491 {
492  G4double c1, c2, c3;
493  // G4cout<<"RI: x1 = "<<x1<<"; "<<"x2 = "<<x2<<G4endl;
494  G4double x12 = x1*x2;
495  c1 = (x2 - x1)/x12;
496  c2 = (x2 - x1)*(x2 + x1)/(x12*x12);
497  c3 = (x2 - x1)*(x1*x1 + x1*x2 + x2*x2)/(x12*x12*x12);
498  // G4cout<<" RI: c1 = "<<c1<<"; "<<"c2 = "<<c2<<"; "<<"c3 = "<<c3<<G4endl;
499 
500  return fA1[k]*log(x2/x1) + fA2[k]*c1 + fA3[k]*c2/2 + fA4[k]*c3/3;
501 
502 } // end of RutherfordIntegral
static c2_factory< G4double > c2
double G4double
Definition: G4Types.hh:76
tuple c1
Definition: plottest35.py:14
void G4PAIySection::SetVerbose ( G4int  v)
inline

Definition at line 143 of file G4PAIySection.hh.

143 { fVerbose = v; };
tuple v
Definition: test.py:18
void G4PAIySection::SplainPAI ( G4double  betaGammaSq)

Definition at line 401 of file G4PAIySection.cc.

402 {
403  G4int k = 1;
404  G4int i = 1;
405 
406  while ( (i < fSplineNumber) && (fSplineNumber < fMaxSplineSize-1) )
407  {
408  if(fSplineEnergy[i+1] > fEnergyInterval[k+1])
409  {
410  k++; // Here next energy point is in next energy interval
411  i++;
412  continue;
413  }
414  // Shifting of arrayes for inserting the geometrical
415  // average of 'i' and 'i+1' energy points to 'i+1' place
416  fSplineNumber++;
417 
418  for(G4int j = fSplineNumber; j >= i+2; j-- )
419  {
420  fSplineEnergy[j] = fSplineEnergy[j-1];
421  fImPartDielectricConst[j] = fImPartDielectricConst[j-1];
422  fRePartDielectricConst[j] = fRePartDielectricConst[j-1];
423  fIntegralTerm[j] = fIntegralTerm[j-1];
424 
425  fDifPAIySection[j] = fDifPAIySection[j-1];
426  fdNdxCerenkov[j] = fdNdxCerenkov[j-1];
427  fdNdxPlasmon[j] = fdNdxPlasmon[j-1];
428  }
429  G4double x1 = fSplineEnergy[i];
430  G4double x2 = fSplineEnergy[i+1];
431  G4double yy1 = fDifPAIySection[i];
432  G4double y2 = fDifPAIySection[i+1];
433 
434  G4double en1 = sqrt(x1*x2);
435  fSplineEnergy[i+1] = en1;
436 
437  // Calculation of logarithmic linear approximation
438  // in this (enr) energy point, which number is 'i+1' now
439 
440  G4double a = log10(y2/yy1)/log10(x2/x1);
441  G4double b = log10(yy1) - a*log10(x1);
442  G4double y = a*log10(en1) + b;
443  y = pow(10.,y);
444 
445  // Calculation of the PAI dif. cross-section at this point
446 
447  fImPartDielectricConst[i+1] = fNormalizationCof*
448  ImPartDielectricConst(k,fSplineEnergy[i+1]);
449  fRePartDielectricConst[i+1] = fNormalizationCof*
450  RePartDielectricConst(fSplineEnergy[i+1]);
451  fIntegralTerm[i+1] = fIntegralTerm[i] + fNormalizationCof*
452  RutherfordIntegral(k,fSplineEnergy[i],
453  fSplineEnergy[i+1]);
454 
455  fDifPAIySection[i+1] = DifPAIySection(i+1,betaGammaSq);
456  fdNdxCerenkov[i+1] = PAIdNdxCerenkov(i+1,betaGammaSq);
457  fdNdxPlasmon[i+1] = PAIdNdxPlasmon(i+1,betaGammaSq);
458 
459  // Condition for next division of this segment or to pass
460  // to higher energies
461 
462  G4double x = 2*(fDifPAIySection[i+1] - y)/(fDifPAIySection[i+1] + y);
463 
464  G4double delta = 2.*(fSplineEnergy[i+1]-fSplineEnergy[i])
465  /(fSplineEnergy[i+1]+fSplineEnergy[i]);
466 
467  if( x < 0 )
468  {
469  x = -x;
470  }
471  if( x > fError && fSplineNumber < fMaxSplineSize-1 && delta > 2.*fDelta )
472  {
473  continue; // next division
474  }
475  i += 2; // pass to next segment
476 
477  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
478  } // close 'while'
479 
480 } // end of SplainPAI
G4double RutherfordIntegral(G4int intervalNumber, G4double limitLow, G4double limitHigh)
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
G4double ImPartDielectricConst(G4int intervalNumber, G4double energy)
G4double PAIdNdxCerenkov(G4int intervalNumber, G4double betaGammaSq)
tuple x
Definition: test.py:50
int G4int
Definition: G4Types.hh:78
tuple b
Definition: test.py:12
double G4double
Definition: G4Types.hh:76
G4double RePartDielectricConst(G4double energy)
G4double DifPAIySection(G4int intervalNumber, G4double betaGammaSq)
G4double PAIdNdxPlasmon(G4int intervalNumber, G4double betaGammaSq)
G4double G4PAIySection::SumOverBordCerenkov ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1115 of file G4PAIySection.cc.

1117 {
1118  G4double x0,x1,y0,yy1,a,e0,c,d,result;
1119 
1120  e0 = en0;
1121  x0 = fSplineEnergy[i];
1122  x1 = fSplineEnergy[i+1];
1123  y0 = fdNdxCerenkov[i];
1124  yy1 = fdNdxCerenkov[i+1];
1125 
1126  // G4cout<<G4endl;
1127  //G4cout<<"SumBordC, i = "<<i<<"; en0 = "<<en0<<"; x0 ="<<x0<<"; x1 = "<<x1
1128  // <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
1129  c = x1/x0;
1130  d = e0/x0;
1131  a = log10(yy1/y0)/log10(c);
1132 
1133  G4double b = 0.0;
1134  if(a < 20.) b = y0/pow(x0,a);
1135 
1136  a += 1.0;
1137  if( a == 0 ) result = b*log(x0/e0);
1138  else result = y0*(x0 - e0*pow(d,a-1))/a;
1139  a += 1.0;
1140 
1141  if( a == 0 ) fIntegralCerenkov[0] += b*log(x0/e0);
1142  else fIntegralCerenkov[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
1143 
1144  //G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
1145 
1146  x0 = fSplineEnergy[i - 1];
1147  x1 = fSplineEnergy[i - 2];
1148  y0 = fdNdxCerenkov[i - 1];
1149  yy1 = fdNdxCerenkov[i - 2];
1150 
1151  //G4cout<<"x0 ="<<x0<<"; x1 = "<<x1
1152  // <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
1153 
1154  c = x1/x0;
1155  d = e0/x0;
1156  a = log10(yy1/y0)/log10(x1/x0);
1157 
1158  // G4cout << "a= " << a << G4endl;
1159  if(a > 20.0) b = 0.0;
1160  else b = y0/pow(x0,a);
1161 
1162  //G4cout << "b= " << b << G4endl;
1163 
1164  a += 1.0;
1165  if( a == 0 ) result += b*log(e0/x0);
1166  else result += y0*(e0*pow(d,a-1) - x0 )/a;
1167  a += 1.0;
1168  //G4cout << "result= " << result << G4endl;
1169 
1170  if( a == 0 ) fIntegralCerenkov[0] += b*log(e0/x0);
1171  else fIntegralCerenkov[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
1172 
1173  //G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
1174 
1175  return result;
1176 
1177 }
G4double G4ParticleHPJENDLHEData::G4double result
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
tuple b
Definition: test.py:12
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13
G4double G4PAIySection::SumOverBorder ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 989 of file G4PAIySection.cc.

991 {
992  G4double x0,x1,y0,yy1,a,/*c,*/d,e0,result;
993 
994  e0 = en0;
995  x0 = fSplineEnergy[i];
996  x1 = fSplineEnergy[i+1];
997  y0 = fDifPAIySection[i];
998  yy1 = fDifPAIySection[i+1];
999 
1000  //c = x1/x0;
1001  d = e0/x0;
1002  a = log10(yy1/y0)/log10(x1/x0);
1003 
1004  G4double b = 0.0;
1005  if(a < 20.) b = y0/pow(x0,a);
1006 
1007  a += 1;
1008  if(a == 0)
1009  {
1010  result = b*log(x0/e0);
1011  }
1012  else
1013  {
1014  result = y0*(x0 - e0*pow(d,a-1))/a;
1015  }
1016  a++;
1017  if(a == 0)
1018  {
1019  fIntegralPAIySection[0] += b*log(x0/e0);
1020  }
1021  else
1022  {
1023  fIntegralPAIySection[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
1024  }
1025  x0 = fSplineEnergy[i - 1];
1026  x1 = fSplineEnergy[i - 2];
1027  y0 = fDifPAIySection[i - 1];
1028  yy1 = fDifPAIySection[i - 2];
1029 
1030  //c = x1/x0;
1031  d = e0/x0;
1032  a = log10(yy1/y0)/log10(x1/x0);
1033  // b0 = log10(y0) - a*log10(x0);
1034  b = y0/pow(x0,a);
1035  a += 1;
1036  if(a == 0)
1037  {
1038  result += b*log(e0/x0);
1039  }
1040  else
1041  {
1042  result += y0*(e0*pow(d,a-1) - x0)/a;
1043  }
1044  a++;
1045  if(a == 0)
1046  {
1047  fIntegralPAIySection[0] += b*log(e0/x0);
1048  }
1049  else
1050  {
1051  fIntegralPAIySection[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
1052  }
1053  return result;
1054 
1055 }
G4double G4ParticleHPJENDLHEData::G4double result
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
tuple b
Definition: test.py:12
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::SumOverBorderdEdx ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1059 of file G4PAIySection.cc.

1061 {
1062  G4double x0,x1,y0,yy1,a,/*c,*/d,e0,result;
1063 
1064  e0 = en0;
1065  x0 = fSplineEnergy[i];
1066  x1 = fSplineEnergy[i+1];
1067  y0 = fDifPAIySection[i];
1068  yy1 = fDifPAIySection[i+1];
1069 
1070  //c = x1/x0;
1071  d = e0/x0;
1072  a = log10(yy1/y0)/log10(x1/x0);
1073 
1074  G4double b = 0.0;
1075  if(a < 20.) b = y0/pow(x0,a);
1076 
1077  a += 2;
1078  if(a == 0)
1079  {
1080  result = b*log(x0/e0);
1081  }
1082  else
1083  {
1084  result = y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
1085  }
1086  x0 = fSplineEnergy[i - 1];
1087  x1 = fSplineEnergy[i - 2];
1088  y0 = fDifPAIySection[i - 1];
1089  yy1 = fDifPAIySection[i - 2];
1090 
1091  //c = x1/x0;
1092  d = e0/x0;
1093  a = log10(yy1/y0)/log10(x1/x0);
1094 
1095  if(a < 20.) b = y0/pow(x0,a);
1096 
1097  a += 2;
1098  if(a == 0)
1099  {
1100  result += b*log(e0/x0);
1101  }
1102  else
1103  {
1104  result += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
1105  }
1106  return result;
1107 
1108 }
G4double G4ParticleHPJENDLHEData::G4double result
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
tuple b
Definition: test.py:12
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::SumOverBordPlasmon ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1184 of file G4PAIySection.cc.

1186 {
1187  G4double x0,x1,y0,yy1,a,c,d,e0,result;
1188 
1189  e0 = en0;
1190  x0 = fSplineEnergy[i];
1191  x1 = fSplineEnergy[i+1];
1192  y0 = fdNdxPlasmon[i];
1193  yy1 = fdNdxPlasmon[i+1];
1194 
1195  c = x1/x0;
1196  d = e0/x0;
1197  a = log10(yy1/y0)/log10(c);
1198 
1199  G4double b = 0.0;
1200  if(a < 20.) b = y0/pow(x0,a);
1201 
1202  a += 1.0;
1203  if( a == 0 ) result = b*log(x0/e0);
1204  else result = y0*(x0 - e0*pow(d,a-1))/a;
1205  a += 1.0;
1206 
1207  if( a == 0 ) fIntegralPlasmon[0] += b*log(x0/e0);
1208  else fIntegralPlasmon[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
1209 
1210  x0 = fSplineEnergy[i - 1];
1211  x1 = fSplineEnergy[i - 2];
1212  y0 = fdNdxPlasmon[i - 1];
1213  yy1 = fdNdxPlasmon[i - 2];
1214 
1215  c = x1/x0;
1216  d = e0/x0;
1217  a = log10(yy1/y0)/log10(c);
1218 
1219  if(a < 20.) b = y0/pow(x0,a);
1220 
1221  a += 1.0;
1222  if( a == 0 ) result += b*log(e0/x0);
1223  else result += y0*(e0*pow(d,a-1) - x0)/a;
1224  a += 1.0;
1225 
1226  if( a == 0 ) fIntegralPlasmon[0] += b*log(e0/x0);
1227  else fIntegralPlasmon[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
1228 
1229  return result;
1230 
1231 }
G4double G4ParticleHPJENDLHEData::G4double result
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
tuple b
Definition: test.py:12
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13
G4double G4PAIySection::SumOverInterCerenkov ( G4int  intervalNumber)

Definition at line 918 of file G4PAIySection.cc.

919 {
920  G4double x0,x1,y0,yy1,a,c,result;
921 
922  x0 = fSplineEnergy[i];
923  x1 = fSplineEnergy[i+1];
924 
925  if( std::abs( 2.*(x1-x0)/(x1+x0) ) < 1.e-6) return 0.;
926 
927  y0 = fdNdxCerenkov[i];
928  yy1 = fdNdxCerenkov[i+1];
929  // G4cout<<"SumC, i = "<<i<<"; x0 ="<<x0<<"; x1 = "<<x1
930  // <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
931 
932  c = x1/x0;
933  a = log10(yy1/y0)/log10(c);
934  G4double b = 0.0;
935  if(a < 20.) b = y0/pow(x0,a);
936 
937  a += 1.0;
938  if(a == 0) result = b*log(c);
939  else result = y0*(x1*pow(c,a-1) - x0)/a;
940  a += 1.0;
941 
942  if( a == 0 ) fIntegralCerenkov[0] += b*log(x1/x0);
943  else fIntegralCerenkov[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
944  // G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
945  return result;
946 
947 } // end of SumOverInterCerenkov
G4double G4ParticleHPJENDLHEData::G4double result
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
tuple b
Definition: test.py:12
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13
G4double G4PAIySection::SumOverInterPlasmon ( G4int  intervalNumber)

Definition at line 955 of file G4PAIySection.cc.

956 {
957  G4double x0,x1,y0,yy1,a,c,result;
958 
959  x0 = fSplineEnergy[i];
960  x1 = fSplineEnergy[i+1];
961 
962  if( std::abs( 2.*(x1-x0)/(x1+x0) ) < 1.e-6) return 0.;
963 
964  y0 = fdNdxPlasmon[i];
965  yy1 = fdNdxPlasmon[i+1];
966  c = x1/x0;
967  a = log10(yy1/y0)/log10(c);
968 
969  G4double b = 0.0;
970  if(a < 20.) b = y0/pow(x0,a);
971 
972  a += 1.0;
973  if(a == 0) result = b*log(x1/x0);
974  else result = y0*(x1*pow(c,a-1) - x0)/a;
975  a += 1.0;
976 
977  if( a == 0 ) fIntegralPlasmon[0] += b*log(x1/x0);
978  else fIntegralPlasmon[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
979 
980  return result;
981 
982 } // end of SumOverInterPlasmon
G4double G4ParticleHPJENDLHEData::G4double result
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
tuple b
Definition: test.py:12
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13
G4double G4PAIySection::SumOverInterval ( G4int  intervalNumber)

Definition at line 842 of file G4PAIySection.cc.

843 {
844  G4double x0,x1,y0,yy1,a,b,c,result;
845 
846  x0 = fSplineEnergy[i];
847  x1 = fSplineEnergy[i+1];
848 
849  if( std::abs( 2.*(x1-x0)/(x1+x0) ) < 1.e-6) return 0.;
850 
851  y0 = fDifPAIySection[i];
852  yy1 = fDifPAIySection[i+1];
853  //G4cout << "## x0= " << x0 << " x1= " << x1 << G4endl;
854  c = x1/x0;
855  //G4cout << "c= " << c << " y0= " << y0 << " yy1= " << yy1 << G4endl;
856  a = log10(yy1/y0)/log10(c);
857  //G4cout << "a= " << a << G4endl;
858  // b = log10(y0) - a*log10(x0);
859  b = y0/pow(x0,a);
860  a += 1;
861  if(a == 0)
862  {
863  result = b*log(x1/x0);
864  }
865  else
866  {
867  result = y0*(x1*pow(c,a-1) - x0)/a;
868  }
869  a++;
870  if(a == 0)
871  {
872  fIntegralPAIySection[0] += b*log(x1/x0);
873  }
874  else
875  {
876  fIntegralPAIySection[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
877  }
878  return result;
879 
880 } // end of SumOverInterval
G4double G4ParticleHPJENDLHEData::G4double result
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
tuple b
Definition: test.py:12
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13
G4double G4PAIySection::SumOverIntervaldEdx ( G4int  intervalNumber)

Definition at line 884 of file G4PAIySection.cc.

885 {
886  G4double x0,x1,y0,yy1,a,b,c,result;
887 
888  x0 = fSplineEnergy[i];
889  x1 = fSplineEnergy[i+1];
890 
891  if( std::abs( 2.*(x1-x0)/(x1+x0) ) < 1.e-6) return 0.;
892 
893  y0 = fDifPAIySection[i];
894  yy1 = fDifPAIySection[i+1];
895  c = x1/x0;
896  a = log10(yy1/y0)/log10(c);
897  // b = log10(y0) - a*log10(x0);
898  b = y0/pow(x0,a);
899  a += 2;
900  if(a == 0)
901  {
902  result = b*log(x1/x0);
903  }
904  else
905  {
906  result = y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
907  }
908  return result;
909 
910 } // end of SumOverInterval
G4double G4ParticleHPJENDLHEData::G4double result
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
tuple b
Definition: test.py:12
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13

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