Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions
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 ( )

Definition at line 79 of file G4PAIySection.cc.

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

Definition at line 108 of file G4PAIySection.cc.

109 {}

Member Function Documentation

void G4PAIySection::ComputeLowEnergyCof ( const G4Material *  material)

Definition at line 239 of file G4PAIySection.cc.

References G4Material::GetElement(), G4Material::GetNumberOfElements(), and G4Element::GetZ().

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

Definition at line 578 of file G4PAIySection.cc.

References python.hepunit::electron_mass_c2, python.hepunit::fine_structure_const, python.hepunit::hbarc, and python.hepunit::pi.

580 {
581  G4double beta, be2,cof,x1,x2,x3,x4,x5,x6,x7,x8,result;
582  //G4double beta, be4;
583  //G4double be4;
584  // G4double betaBohr2 = fine_structure_const*fine_structure_const;
585  // G4double betaBohr4 = betaBohr2*betaBohr2*4.0;
586  be2 = betaGammaSq/(1 + betaGammaSq);
587  //be4 = be2*be2;
588  beta = sqrt(be2);
589  cof = 1;
590  x1 = log(2*electron_mass_c2/fSplineEnergy[i]);
591 
592  if( betaGammaSq < 0.01 ) x2 = log(be2);
593  else
594  {
595  x2 = -log( (1/betaGammaSq - fRePartDielectricConst[i])*
596  (1/betaGammaSq - fRePartDielectricConst[i]) +
597  fImPartDielectricConst[i]*fImPartDielectricConst[i] )/2;
598  }
599  if( fImPartDielectricConst[i] == 0.0 ||betaGammaSq < 0.01 )
600  {
601  x6=0;
602  }
603  else
604  {
605  x3 = -fRePartDielectricConst[i] + 1/betaGammaSq;
606  x5 = -1 - fRePartDielectricConst[i] +
607  be2*((1 +fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
608  fImPartDielectricConst[i]*fImPartDielectricConst[i]);
609 
610  x7 = atan2(fImPartDielectricConst[i],x3);
611  x6 = x5 * x7;
612  }
613  // if(fImPartDielectricConst[i] == 0) x6 = 0;
614 
615  x4 = ((x1 + x2)*fImPartDielectricConst[i] + x6)/hbarc;
616  // if( x4 < 0.0 ) x4 = 0.0;
617  x8 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
618  fImPartDielectricConst[i]*fImPartDielectricConst[i];
619 
620  result = (x4 + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i]);
621  if(result < 1.0e-8) result = 1.0e-8;
622  result *= fine_structure_const/be2/pi;
623  // low energy correction
624 
625  G4double lowCof = fLowEnergyCof; // 6.0 ; // Ar ~ 4.; -> fLowCof as f(Z1,Z2)?
626 
627  result *= (1 - exp(-beta/betaBohr/lowCof));
628 
629  // result *= (1-exp(-beta/betaBohr))*(1-exp(-beta/betaBohr));
630  // result *= (1-exp(-be2/betaBohr2));
631  // result *= (1-exp(-be4/betaBohr4));
632  // if(fDensity >= 0.1)
633  if(x8 > 0.)
634  {
635  result /= x8;
636  }
637  return result;
638 
639 } // end of DifPAIySection
python.hepunit.electron_mass_c2
float electron_mass_c2
Definition: hepunit.py:274
G4double
double G4double
Definition: G4Types.hh:76
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 229 of file G4PAIySection.hh.

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

Definition at line 223 of file G4PAIySection.hh.

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

Definition at line 217 of file G4PAIySection.hh.

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

Definition at line 235 of file G4PAIySection.hh.

236 {
237  if(i < 1 || i > fSplineNumber) { CallError(i, "GetIntegralPlasmon"); }
238  return fIntegralPlasmon[i];
239 }
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 206 of file G4PAIySection.hh.

207 {
208  return fLorentzFactor[j];
209 }
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 201 of file G4PAIySection.hh.

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

Definition at line 211 of file G4PAIySection.hh.

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

Definition at line 116 of file G4PAIySection.hh.

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

Definition at line 1265 of file G4PAIySection.cc.

References G4Poisson(), G4UniformRand, and position.

1266 {
1267  G4int iTransfer ;
1268  G4long numOfCollisions;
1269  G4double loss = 0.0;
1270  G4double meanNumber, position;
1271 
1272  // G4cout<<" G4PAIySection::GetStepCreLosnkovs "<<G4endl;
1273 
1274 
1275 
1276  meanNumber = fIntegralCerenkov[1]*step;
1277  numOfCollisions = G4Poisson(meanNumber);
1278 
1279  // G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
1280 
1281  while(numOfCollisions)
1282  {
1283  position = fIntegralCerenkov[1]*G4UniformRand();
1284 
1285  for( iTransfer=1; iTransfer<=fSplineNumber; iTransfer++ )
1286  {
1287  if( position >= fIntegralCerenkov[iTransfer] ) break;
1288  }
1289  loss += fSplineEnergy[iTransfer] ;
1290  numOfCollisions--;
1291  }
1292  // G4cout<<"PAI Cerenkov loss = "<<loss/keV<<" keV"<<G4endl;
1293 
1294  return loss;
1295 }
G4Poisson
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:51
G4long
long G4long
Definition: G4Types.hh:80
G4int
int G4int
Definition: G4Types.hh:78
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:87
position
int position
Definition: filter.cc:7
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::GetStepEnergyLoss ( G4double  step)

Definition at line 1229 of file G4PAIySection.cc.

References G4Poisson(), G4UniformRand, and position.

1230 {
1231  G4int iTransfer ;
1232  G4long numOfCollisions;
1233  G4double loss = 0.0;
1234  G4double meanNumber, position;
1235 
1236  // G4cout<<" G4PAIySection::GetStepEnergyLoss "<<G4endl;
1237 
1238 
1239 
1240  meanNumber = fIntegralPAIySection[1]*step;
1241  numOfCollisions = G4Poisson(meanNumber);
1242 
1243  // G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
1244 
1245  while(numOfCollisions)
1246  {
1247  position = fIntegralPAIySection[1]*G4UniformRand();
1248 
1249  for( iTransfer=1; iTransfer<=fSplineNumber; iTransfer++ )
1250  {
1251  if( position >= fIntegralPAIySection[iTransfer] ) break;
1252  }
1253  loss += fSplineEnergy[iTransfer] ;
1254  numOfCollisions--;
1255  }
1256  // G4cout<<"PAI energy loss = "<<loss/keV<<" keV"<<G4endl;
1257 
1258  return loss;
1259 }
G4Poisson
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:51
G4long
long G4long
Definition: G4Types.hh:80
G4int
int G4int
Definition: G4Types.hh:78
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:87
position
int position
Definition: filter.cc:7
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::GetStepPlasmonLoss ( G4double  step)

Definition at line 1301 of file G4PAIySection.cc.

References G4Poisson(), G4UniformRand, and position.

1302 {
1303  G4int iTransfer ;
1304  G4long numOfCollisions;
1305  G4double loss = 0.0;
1306  G4double meanNumber, position;
1307 
1308  // G4cout<<" G4PAIySection::GetStepCreLosnkovs "<<G4endl;
1309 
1310 
1311 
1312  meanNumber = fIntegralPlasmon[1]*step;
1313  numOfCollisions = G4Poisson(meanNumber);
1314 
1315  // G4cout<<"numOfCollisions = "<<numOfCollisions<<G4endl;
1316 
1317  while(numOfCollisions)
1318  {
1319  position = fIntegralPlasmon[1]*G4UniformRand();
1320 
1321  for( iTransfer=1; iTransfer<=fSplineNumber; iTransfer++ )
1322  {
1323  if( position >= fIntegralPlasmon[iTransfer] ) break;
1324  }
1325  loss += fSplineEnergy[iTransfer] ;
1326  numOfCollisions--;
1327  }
1328  // G4cout<<"PAI Plasmon loss = "<<loss/keV<<" keV"<<G4endl;
1329 
1330  return loss;
1331 }
G4Poisson
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:51
G4long
long G4long
Definition: G4Types.hh:80
G4int
int G4int
Definition: G4Types.hh:78
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:87
position
int position
Definition: filter.cc:7
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::ImPartDielectricConst ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 501 of file G4PAIySection.cc.

References python.hepunit::hbarc.

503 {
504  G4double energy2,energy3,energy4,result;
505 
506  energy2 = energy1*energy1;
507  energy3 = energy2*energy1;
508  energy4 = energy3*energy1;
509 
510  result = fA1[k]/energy1+fA2[k]/energy2+fA3[k]/energy3+fA4[k]/energy4;
511  result *=hbarc/energy1;
512 
513  return result;
514 
515 } // end of ImPartDielectricConst
G4double
double G4double
Definition: G4Types.hh:76
void G4PAIySection::Initialize ( const G4Material *  material,
G4double  maxEnergyTransfer,
G4double  betaGammaSq,
G4SandiaTable *  sandia 
)

Definition at line 115 of file G4PAIySection.cc.

References python.hepunit::eV, G4cout, G4endl, G4Material::GetDensity(), G4Material::GetElectronDensity(), G4SandiaTable::GetMaxInterval(), G4SandiaTable::GetSandiaMatTablePAI(), and python.hepunit::keV.

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

Definition at line 272 of file G4PAIySection.cc.

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

Definition at line 774 of file G4PAIySection.cc.

775 {
776  G4int i, k;
777  fIntegralCerenkov[fSplineNumber] = 0;
778  fIntegralCerenkov[0] = 0;
779  k = fIntervalNumber -1;
780 
781  for( i = fSplineNumber-1; i >= 1; i-- )
782  {
783  if(fSplineEnergy[i] >= fEnergyInterval[k])
784  {
785  fIntegralCerenkov[i] = fIntegralCerenkov[i+1] + SumOverInterCerenkov(i);
786  // G4cout<<"int: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
787  }
788  else
789  {
790  fIntegralCerenkov[i] = fIntegralCerenkov[i+1] +
791  SumOverBordCerenkov(i+1,fEnergyInterval[k]);
792  k--;
793  // G4cout<<"bord: i = "<<i<<"; sumC = "<<fIntegralCerenkov[i]<<G4endl;
794  }
795  }
796 
797 } // end of IntegralCerenkov
G4PAIySection::SumOverInterCerenkov
G4double SumOverInterCerenkov(G4int intervalNumber)
Definition: G4PAIySection.cc:908
G4int
int G4int
Definition: G4Types.hh:78
G4PAIySection::SumOverBordCerenkov
G4double SumOverBordCerenkov(G4int intervalNumber, G4double energy)
Definition: G4PAIySection.cc:1105
void G4PAIySection::IntegralPAIySection ( )

Definition at line 743 of file G4PAIySection.cc.

744 {
745  fIntegralPAIySection[fSplineNumber] = 0;
746  fIntegralPAIdEdx[fSplineNumber] = 0;
747  fIntegralPAIySection[0] = 0;
748  G4int k = fIntervalNumber -1;
749 
750  for(G4int i = fSplineNumber-1; i >= 1; i--)
751  {
752  if(fSplineEnergy[i] >= fEnergyInterval[k])
753  {
754  fIntegralPAIySection[i] = fIntegralPAIySection[i+1] + SumOverInterval(i);
755  fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] + SumOverIntervaldEdx(i);
756  }
757  else
758  {
759  fIntegralPAIySection[i] = fIntegralPAIySection[i+1] +
760  SumOverBorder(i+1,fEnergyInterval[k]);
761  fIntegralPAIdEdx[i] = fIntegralPAIdEdx[i+1] +
762  SumOverBorderdEdx(i+1,fEnergyInterval[k]);
763  k--;
764  }
765  }
766 } // end of IntegralPAIySection
G4PAIySection::SumOverIntervaldEdx
G4double SumOverIntervaldEdx(G4int intervalNumber)
Definition: G4PAIySection.cc:874
G4PAIySection::SumOverBorderdEdx
G4double SumOverBorderdEdx(G4int intervalNumber, G4double energy)
Definition: G4PAIySection.cc:1049
G4int
int G4int
Definition: G4Types.hh:78
G4PAIySection::SumOverBorder
G4double SumOverBorder(G4int intervalNumber, G4double energy)
Definition: G4PAIySection.cc:979
G4PAIySection::SumOverInterval
G4double SumOverInterval(G4int intervalNumber)
Definition: G4PAIySection.cc:832
void G4PAIySection::IntegralPlasmon ( )

Definition at line 805 of file G4PAIySection.cc.

806 {
807  fIntegralPlasmon[fSplineNumber] = 0;
808  fIntegralPlasmon[0] = 0;
809  G4int k = fIntervalNumber -1;
810  for(G4int i=fSplineNumber-1;i>=1;i--)
811  {
812  if(fSplineEnergy[i] >= fEnergyInterval[k])
813  {
814  fIntegralPlasmon[i] = fIntegralPlasmon[i+1] + SumOverInterPlasmon(i);
815  }
816  else
817  {
818  fIntegralPlasmon[i] = fIntegralPlasmon[i+1] +
819  SumOverBordPlasmon(i+1,fEnergyInterval[k]);
820  k--;
821  }
822  }
823 
824 } // end of IntegralPlasmon
G4int
int G4int
Definition: G4Types.hh:78
G4PAIySection::SumOverBordPlasmon
G4double SumOverBordPlasmon(G4int intervalNumber, G4double energy)
Definition: G4PAIySection.cc:1176
G4PAIySection::SumOverInterPlasmon
G4double SumOverInterPlasmon(G4int intervalNumber)
Definition: G4PAIySection.cc:945
void G4PAIySection::NormShift ( G4double  betaGammaSq)

Definition at line 323 of file G4PAIySection.cc.

References python.hepunit::electron_mass_c2, python.hepunit::fine_structure_const, python.hepunit::hbarc, python.hepunit::pi, and test::x.

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

Definition at line 645 of file G4PAIySection.cc.

References python.hepunit::fine_structure_const, python.hepunit::hbarc, and python.hepunit::pi.

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

Definition at line 704 of file G4PAIySection.cc.

References python.hepunit::electron_mass_c2, python.hepunit::fine_structure_const, python.hepunit::hbarc, and python.hepunit::pi.

706 {
707  G4double cof, resonance, modul2, dNdxP;
708  G4double be2, be4;
709 
710  cof = 1;
711 
712  be2 = betaGammaSq/(1 + betaGammaSq);
713  be4 = be2*be2;
714 
715  resonance = log(2*electron_mass_c2*be2/fSplineEnergy[i]);
716  resonance *= fImPartDielectricConst[i]/hbarc;
717 
718  dNdxP = ( resonance + cof*fIntegralTerm[i]/fSplineEnergy[i]/fSplineEnergy[i] );
719 
720  if( dNdxP < 1.0e-8 ) dNdxP = 1.0e-8;
721 
722  dNdxP *= fine_structure_const/be2/pi;
723  dNdxP *= (1-exp(-be4/betaBohr4));
724 
725 // if( fDensity >= 0.1 )
726 // {
727  modul2 = (1 + fRePartDielectricConst[i])*(1 + fRePartDielectricConst[i]) +
728  fImPartDielectricConst[i]*fImPartDielectricConst[i];
729  if(modul2 > 0.)
730  {
731  dNdxP /= modul2;
732  }
733  return dNdxP;
734 
735 } // end of PAIdNdxPlasmon
python.hepunit.electron_mass_c2
float electron_mass_c2
Definition: hepunit.py:274
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::RePartDielectricConst ( G4double  energy)

Definition at line 524 of file G4PAIySection.cc.

References plottest35::c1, python.hepunit::hbarc, and python.hepunit::pi.

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

Definition at line 480 of file G4PAIySection.cc.

References plottest35::c1.

483 {
484  G4double c1, c2, c3;
485  // G4cout<<"RI: x1 = "<<x1<<"; "<<"x2 = "<<x2<<G4endl;
486  c1 = (x2 - x1)/x1/x2;
487  c2 = (x2 - x1)*(x2 + x1)/x1/x1/x2/x2;
488  c3 = (x2 - x1)*(x1*x1 + x1*x2 + x2*x2)/x1/x1/x1/x2/x2/x2;
489  // G4cout<<" RI: c1 = "<<c1<<"; "<<"c2 = "<<c2<<"; "<<"c3 = "<<c3<<G4endl;
490 
491  return fA1[k]*log(x2/x1) + fA2[k]*c1 + fA3[k]*c2/2 + fA4[k]*c3/3;
492 
493 } // end of RutherfordIntegral
G4double
double G4double
Definition: G4Types.hh:76
plottest35.c1
tuple c1
Definition: plottest35.py:14
void G4PAIySection::SetVerbose ( G4int  v)
inline

Definition at line 143 of file G4PAIySection.hh.

References test::v.

143 { fVerbose = v; };
void G4PAIySection::SplainPAI ( G4double  betaGammaSq)

Definition at line 395 of file G4PAIySection.cc.

References test::a, test::b, and test::x.

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

Definition at line 1105 of file G4PAIySection.cc.

References test::a, test::b, and test::c.

1107 {
1108  G4double x0,x1,y0,yy1,a,e0,c,d,result;
1109 
1110  e0 = en0;
1111  x0 = fSplineEnergy[i];
1112  x1 = fSplineEnergy[i+1];
1113  y0 = fdNdxCerenkov[i];
1114  yy1 = fdNdxCerenkov[i+1];
1115 
1116  // G4cout<<G4endl;
1117  //G4cout<<"SumBordC, i = "<<i<<"; en0 = "<<en0<<"; x0 ="<<x0<<"; x1 = "<<x1
1118  // <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
1119  c = x1/x0;
1120  d = e0/x0;
1121  a = log10(yy1/y0)/log10(c);
1122 
1123  G4double b = 0.0;
1124  if(a < 20.) b = y0/pow(x0,a);
1125 
1126  a += 1.0;
1127  if( a == 0 ) result = b*log(x0/e0);
1128  else result = y0*(x0 - e0*pow(d,a-1))/a;
1129  a += 1.0;
1130 
1131  if( a == 0 ) fIntegralCerenkov[0] += b*log(x0/e0);
1132  else fIntegralCerenkov[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
1133 
1134  //G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
1135 
1136  x0 = fSplineEnergy[i - 1];
1137  x1 = fSplineEnergy[i - 2];
1138  y0 = fdNdxCerenkov[i - 1];
1139  yy1 = fdNdxCerenkov[i - 2];
1140 
1141  //G4cout<<"x0 ="<<x0<<"; x1 = "<<x1
1142  // <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
1143 
1144  c = x1/x0;
1145  d = e0/x0;
1146  a = log10(yy1/y0)/log10(x1/x0);
1147 
1148  // G4cout << "a= " << a << G4endl;
1149  if(a < 20.) b = y0/pow(x0,a);
1150 
1151  if(a > 20.0) b = 0.0;
1152  else b = y0/pow(x0,a); // pow(10.,b0);
1153 
1154  //G4cout << "b= " << b << G4endl;
1155 
1156  a += 1.0;
1157  if( a == 0 ) result += b*log(e0/x0);
1158  else result += y0*(e0*pow(d,a-1) - x0 )/a;
1159  a += 1.0;
1160  //G4cout << "result= " << result << G4endl;
1161 
1162  if( a == 0 ) fIntegralCerenkov[0] += b*log(e0/x0);
1163  else fIntegralCerenkov[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
1164 
1165  //G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
1166 
1167  return result;
1168 
1169 }
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::SumOverBorder ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 979 of file G4PAIySection.cc.

References test::a, and test::b.

981 {
982  G4double x0,x1,y0,yy1,a,/*c,*/d,e0,result;
983 
984  e0 = en0;
985  x0 = fSplineEnergy[i];
986  x1 = fSplineEnergy[i+1];
987  y0 = fDifPAIySection[i];
988  yy1 = fDifPAIySection[i+1];
989 
990  //c = x1/x0;
991  d = e0/x0;
992  a = log10(yy1/y0)/log10(x1/x0);
993 
994  G4double b = 0.0;
995  if(a < 20.) b = y0/pow(x0,a);
996 
997  a += 1;
998  if(a == 0)
999  {
1000  result = b*log(x0/e0);
1001  }
1002  else
1003  {
1004  result = y0*(x0 - e0*pow(d,a-1))/a;
1005  }
1006  a++;
1007  if(a == 0)
1008  {
1009  fIntegralPAIySection[0] += b*log(x0/e0);
1010  }
1011  else
1012  {
1013  fIntegralPAIySection[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
1014  }
1015  x0 = fSplineEnergy[i - 1];
1016  x1 = fSplineEnergy[i - 2];
1017  y0 = fDifPAIySection[i - 1];
1018  yy1 = fDifPAIySection[i - 2];
1019 
1020  //c = x1/x0;
1021  d = e0/x0;
1022  a = log10(yy1/y0)/log10(x1/x0);
1023  // b0 = log10(y0) - a*log10(x0);
1024  b = y0/pow(x0,a);
1025  a += 1;
1026  if(a == 0)
1027  {
1028  result += b*log(e0/x0);
1029  }
1030  else
1031  {
1032  result += y0*(e0*pow(d,a-1) - x0)/a;
1033  }
1034  a++;
1035  if(a == 0)
1036  {
1037  fIntegralPAIySection[0] += b*log(e0/x0);
1038  }
1039  else
1040  {
1041  fIntegralPAIySection[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
1042  }
1043  return result;
1044 
1045 }
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::SumOverBorderdEdx ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1049 of file G4PAIySection.cc.

References test::a, and test::b.

1051 {
1052  G4double x0,x1,y0,yy1,a,/*c,*/d,e0,result;
1053 
1054  e0 = en0;
1055  x0 = fSplineEnergy[i];
1056  x1 = fSplineEnergy[i+1];
1057  y0 = fDifPAIySection[i];
1058  yy1 = fDifPAIySection[i+1];
1059 
1060  //c = x1/x0;
1061  d = e0/x0;
1062  a = log10(yy1/y0)/log10(x1/x0);
1063 
1064  G4double b = 0.0;
1065  if(a < 20.) b = y0/pow(x0,a);
1066 
1067  a += 2;
1068  if(a == 0)
1069  {
1070  result = b*log(x0/e0);
1071  }
1072  else
1073  {
1074  result = y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
1075  }
1076  x0 = fSplineEnergy[i - 1];
1077  x1 = fSplineEnergy[i - 2];
1078  y0 = fDifPAIySection[i - 1];
1079  yy1 = fDifPAIySection[i - 2];
1080 
1081  //c = x1/x0;
1082  d = e0/x0;
1083  a = log10(yy1/y0)/log10(x1/x0);
1084 
1085  if(a < 20.) b = y0/pow(x0,a);
1086 
1087  a += 2;
1088  if(a == 0)
1089  {
1090  result += b*log(e0/x0);
1091  }
1092  else
1093  {
1094  result += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
1095  }
1096  return result;
1097 
1098 }
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::SumOverBordPlasmon ( G4int  intervalNumber,
G4double  energy 
)

Definition at line 1176 of file G4PAIySection.cc.

References test::a, test::b, and test::c.

1178 {
1179  G4double x0,x1,y0,yy1,a,c,d,e0,result;
1180 
1181  e0 = en0;
1182  x0 = fSplineEnergy[i];
1183  x1 = fSplineEnergy[i+1];
1184  y0 = fdNdxPlasmon[i];
1185  yy1 = fdNdxPlasmon[i+1];
1186 
1187  c = x1/x0;
1188  d = e0/x0;
1189  a = log10(yy1/y0)/log10(c);
1190 
1191  G4double b = 0.0;
1192  if(a < 20.) b = y0/pow(x0,a);
1193 
1194  a += 1.0;
1195  if( a == 0 ) result = b*log(x0/e0);
1196  else result = y0*(x0 - e0*pow(d,a-1))/a;
1197  a += 1.0;
1198 
1199  if( a == 0 ) fIntegralPlasmon[0] += b*log(x0/e0);
1200  else fIntegralPlasmon[0] += y0*(x0*x0 - e0*e0*pow(d,a-2))/a;
1201 
1202  x0 = fSplineEnergy[i - 1];
1203  x1 = fSplineEnergy[i - 2];
1204  y0 = fdNdxPlasmon[i - 1];
1205  yy1 = fdNdxPlasmon[i - 2];
1206 
1207  c = x1/x0;
1208  d = e0/x0;
1209  a = log10(yy1/y0)/log10(c);
1210 
1211  if(a < 20.) b = y0/pow(x0,a);
1212 
1213  a += 1.0;
1214  if( a == 0 ) result += b*log(e0/x0);
1215  else result += y0*(e0*pow(d,a-1) - x0)/a;
1216  a += 1.0;
1217 
1218  if( a == 0 ) fIntegralPlasmon[0] += b*log(e0/x0);
1219  else fIntegralPlasmon[0] += y0*(e0*e0*pow(d,a-2) - x0*x0)/a;
1220 
1221  return result;
1222 
1223 }
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::SumOverInterCerenkov ( G4int  intervalNumber)

Definition at line 908 of file G4PAIySection.cc.

References test::a, test::b, and test::c.

909 {
910  G4double x0,x1,y0,yy1,a,c,result;
911 
912  x0 = fSplineEnergy[i];
913  x1 = fSplineEnergy[i+1];
914 
915  if( std::abs( 2.*(x1-x0)/(x1+x0) ) < 1.e-6) return 0.;
916 
917  y0 = fdNdxCerenkov[i];
918  yy1 = fdNdxCerenkov[i+1];
919  // G4cout<<"SumC, i = "<<i<<"; x0 ="<<x0<<"; x1 = "<<x1
920  // <<"; y0 = "<<y0<<"; yy1 = "<<yy1<<G4endl;
921 
922  c = x1/x0;
923  a = log10(yy1/y0)/log10(c);
924  G4double b = 0.0;
925  if(a < 20.) b = y0/pow(x0,a);
926 
927  a += 1.0;
928  if(a == 0) result = b*log(c);
929  else result = y0*(x1*pow(c,a-1) - x0)/a;
930  a += 1.0;
931 
932  if( a == 0 ) fIntegralCerenkov[0] += b*log(x1/x0);
933  else fIntegralCerenkov[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
934  // G4cout<<"a = "<<a<<"; b = "<<b<<"; result = "<<result<<G4endl;
935  return result;
936 
937 } // end of SumOverInterCerenkov
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::SumOverInterPlasmon ( G4int  intervalNumber)

Definition at line 945 of file G4PAIySection.cc.

References test::a, test::b, and test::c.

946 {
947  G4double x0,x1,y0,yy1,a,c,result;
948 
949  x0 = fSplineEnergy[i];
950  x1 = fSplineEnergy[i+1];
951 
952  if( std::abs( 2.*(x1-x0)/(x1+x0) ) < 1.e-6) return 0.;
953 
954  y0 = fdNdxPlasmon[i];
955  yy1 = fdNdxPlasmon[i+1];
956  c = x1/x0;
957  a = log10(yy1/y0)/log10(c);
958 
959  G4double b = 0.0;
960  if(a < 20.) b = y0/pow(x0,a);
961 
962  a += 1.0;
963  if(a == 0) result = b*log(x1/x0);
964  else result = y0*(x1*pow(c,a-1) - x0)/a;
965  a += 1.0;
966 
967  if( a == 0 ) fIntegralPlasmon[0] += b*log(x1/x0);
968  else fIntegralPlasmon[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
969 
970  return result;
971 
972 } // end of SumOverInterPlasmon
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::SumOverInterval ( G4int  intervalNumber)

Definition at line 832 of file G4PAIySection.cc.

References test::a, test::b, and test::c.

833 {
834  G4double x0,x1,y0,yy1,a,b,c,result;
835 
836  x0 = fSplineEnergy[i];
837  x1 = fSplineEnergy[i+1];
838 
839  if( std::abs( 2.*(x1-x0)/(x1+x0) ) < 1.e-6) return 0.;
840 
841  y0 = fDifPAIySection[i];
842  yy1 = fDifPAIySection[i+1];
843  //G4cout << "## x0= " << x0 << " x1= " << x1 << G4endl;
844  c = x1/x0;
845  //G4cout << "c= " << c << " y0= " << y0 << " yy1= " << yy1 << G4endl;
846  a = log10(yy1/y0)/log10(c);
847  //G4cout << "a= " << a << G4endl;
848  // b = log10(y0) - a*log10(x0);
849  b = y0/pow(x0,a);
850  a += 1;
851  if(a == 0)
852  {
853  result = b*log(x1/x0);
854  }
855  else
856  {
857  result = y0*(x1*pow(c,a-1) - x0)/a;
858  }
859  a++;
860  if(a == 0)
861  {
862  fIntegralPAIySection[0] += b*log(x1/x0);
863  }
864  else
865  {
866  fIntegralPAIySection[0] += y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
867  }
868  return result;
869 
870 } // end of SumOverInterval
G4double
double G4double
Definition: G4Types.hh:76
G4double G4PAIySection::SumOverIntervaldEdx ( G4int  intervalNumber)

Definition at line 874 of file G4PAIySection.cc.

References test::a, test::b, and test::c.

875 {
876  G4double x0,x1,y0,yy1,a,b,c,result;
877 
878  x0 = fSplineEnergy[i];
879  x1 = fSplineEnergy[i+1];
880 
881  if( std::abs( 2.*(x1-x0)/(x1+x0) ) < 1.e-6) return 0.;
882 
883  y0 = fDifPAIySection[i];
884  yy1 = fDifPAIySection[i+1];
885  c = x1/x0;
886  a = log10(yy1/y0)/log10(c);
887  // b = log10(y0) - a*log10(x0);
888  b = y0/pow(x0,a);
889  a += 2;
890  if(a == 0)
891  {
892  result = b*log(x1/x0);
893  }
894  else
895  {
896  result = y0*(x1*x1*pow(c,a-2) - x0*x0)/a;
897  }
898  return result;
899 
900 } // end of SumOverInterval
G4double
double G4double
Definition: G4Types.hh:76
The documentation for this class was generated from the following files:
Generated on Wed Apr 30 2014 15:56:51 for Geant4.10 by   doxygen 1.8.7