Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/lowenergy/src/G4PenelopeCrossSection.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /processes/electromagnetic/lowenergy/src/G4PenelopeCrossSection.cc (Version 11.3.0) and /processes/electromagnetic/lowenergy/src/G4PenelopeCrossSection.cc (Version 5.0)


  1 //                                                  1 
  2 // *******************************************    
  3 // * License and Disclaimer                       
  4 // *                                              
  5 // * The  Geant4 software  is  copyright of th    
  6 // * the Geant4 Collaboration.  It is provided    
  7 // * conditions of the Geant4 Software License    
  8 // * LICENSE and available at  http://cern.ch/    
  9 // * include a list of copyright holders.         
 10 // *                                              
 11 // * Neither the authors of this software syst    
 12 // * institutes,nor the agencies providing fin    
 13 // * work  make  any representation or  warran    
 14 // * regarding  this  software system or assum    
 15 // * use.  Please see the license in the file     
 16 // * for the full disclaimer and the limitatio    
 17 // *                                              
 18 // * This  code  implementation is the result     
 19 // * technical work of the GEANT4 collaboratio    
 20 // * By using,  copying,  modifying or  distri    
 21 // * any work based  on the software)  you  ag    
 22 // * use  in  resulting  scientific  publicati    
 23 // * acceptance of all terms of the Geant4 Sof    
 24 // *******************************************    
 25 //                                                
 26 //                                                
 27 // Author: Luciano Pandola                        
 28 //                                                
 29 // History:                                       
 30 // --------                                       
 31 // 18 Mar 2010   L Pandola    First implementa    
 32 // 09 Mar 2012   L. Pandola   Add public metho    
 33 //                            the absolute and    
 34 //                            sections indepen    
 35 //                                                
 36 #include "G4PenelopeCrossSection.hh"              
 37 #include "G4SystemOfUnits.hh"                     
 38 #include "G4PhysicsTable.hh"                      
 39 #include "G4PhysicsFreeVector.hh"                 
 40 #include "G4DataVector.hh"                        
 41 #include "G4Exp.hh"                               
 42 #include "G4Log.hh"                               
 43                                                   
 44 //....oooOO0OOooo........oooOO0OOooo........oo    
 45 G4PenelopeCrossSection::G4PenelopeCrossSection    
 46   fSoftCrossSections(nullptr),                    
 47   fHardCrossSections(nullptr),fShellCrossSecti    
 48   fShellNormalizedCrossSections(nullptr),         
 49   fNumberOfEnergyPoints(nPointsE),fNumberOfShe    
 50 {                                                 
 51   //check the number of points is not zero        
 52   if (!fNumberOfEnergyPoints)                     
 53     {                                             
 54       G4ExceptionDescription ed;                  
 55       ed << "G4PenelopeCrossSection: invalid n    
 56       G4Exception("G4PenelopeCrossSection::G4P    
 57       "em2017",FatalException,ed);                
 58     }                                             
 59                                                   
 60   fIsNormalized = false;                          
 61                                                   
 62   // 1) soft XS table                             
 63   fSoftCrossSections = new G4PhysicsTable();      
 64   //the table contains 3 G4PhysicsFreeVectors,    
 65   //(fSoftCrossSections)[0] -->  log XS0 vs. l    
 66   //(fSoftCrossSections)[1] -->  log XS1 vs. l    
 67   //(fSoftCrossSections)[2] -->  log XS2 vs. l    
 68                                                   
 69   //everything is log-log                         
 70   for (size_t i=0;i<3;i++)                        
 71     fSoftCrossSections->push_back(new G4Physic    
 72                                                   
 73   //2) hard XS table                              
 74   fHardCrossSections = new G4PhysicsTable();      
 75   //the table contains 3 G4PhysicsFreeVectors,    
 76   //(fHardCrossSections)[0] -->  log XH0 vs. l    
 77   //(fHardCrossSections)[1] -->  log XH1 vs. l    
 78   //(fHardCrossSections)[2] -->  log XH2 vs. l    
 79                                                   
 80   //everything is log-log                         
 81   for (size_t i=0;i<3;i++)                        
 82     fHardCrossSections->push_back(new G4Physic    
 83                                                   
 84   //3) shell XS table, if it is the case          
 85   if (fNumberOfShells)                            
 86     {                                             
 87       fShellCrossSections = new G4PhysicsTable    
 88       fShellNormalizedCrossSections = new G4Ph    
 89       //the table has to contain numberofShell    
 90       //(theTable)[ishell] --> cross section f    
 91       for (size_t i=0;i<fNumberOfShells;i++)      
 92   {                                               
 93     fShellCrossSections->push_back(new G4Physi    
 94     fShellNormalizedCrossSections->push_back(n    
 95   }                                               
 96     }                                             
 97 }                                                 
 98                                                   
 99 //....oooOO0OOooo........oooOO0OOooo........oo    
100 G4PenelopeCrossSection::~G4PenelopeCrossSectio    
101 {                                                 
102   //clean up tables                               
103   if (fShellCrossSections)                        
104     {                                             
105       fShellCrossSections->clearAndDestroy();     
106       delete fShellCrossSections;                 
107     }                                             
108   if (fShellNormalizedCrossSections)              
109     {                                             
110       fShellNormalizedCrossSections->clearAndD    
111       delete fShellNormalizedCrossSections;       
112     }                                             
113   if (fSoftCrossSections)                         
114     {                                             
115       fSoftCrossSections->clearAndDestroy();      
116       delete fSoftCrossSections;                  
117     }                                             
118   if (fHardCrossSections)                         
119     {                                             
120       fHardCrossSections->clearAndDestroy();      
121       delete fHardCrossSections;                  
122     }                                             
123 }                                                 
124                                                   
125 //....oooOO0OOooo........oooOO0OOooo........oo    
126 void G4PenelopeCrossSection::AddCrossSectionPo    
127               G4double XH0,                       
128               G4double XH1, G4double XH2,         
129               G4double XS0, G4double XS1,         
130               G4double XS2)                       
131 {                                                 
132   if (!fSoftCrossSections || !fHardCrossSectio    
133     {                                             
134       G4cout << "Something wrong in G4Penelope    
135   G4endl;                                         
136       G4cout << "Trying to fill un-initialized    
137       return;                                     
138     }                                             
139                                                   
140   //fill vectors                                  
141   G4PhysicsFreeVector* theVector = (G4PhysicsF    
142                                                   
143   if (binNumber >= fNumberOfEnergyPoints)         
144     {                                             
145       G4cout << "Something wrong in G4Penelope    
146   G4endl;                                         
147       G4cout << "Trying to register more point    
148       return;                                     
149     }                                             
150    G4double logEne = G4Log(energy);               
151                                                   
152    //XS0                                          
153    G4double val = G4Log(std::max(XS0,1e-42*cm2    
154    theVector->PutValues(binNumber,logEne,val);    
155                                                   
156    //XS1                                          
157    theVector = (G4PhysicsFreeVector*) (*fSoftC    
158    val =  G4Log(std::max(XS1,1e-42*eV*cm2)); /    
159    theVector->PutValues(binNumber,logEne,val);    
160                                                   
161    //XS2                                          
162    theVector = (G4PhysicsFreeVector*) (*fSoftC    
163    val =  G4Log(std::max(XS2,1e-42*eV*eV*cm2))    
164    theVector->PutValues(binNumber,logEne,val);    
165                                                   
166    //XH0                                          
167    theVector = (G4PhysicsFreeVector*) (*fHardC    
168    val =  G4Log(std::max(XH0,1e-42*cm2)); //av    
169    theVector->PutValues(binNumber,logEne,val);    
170                                                   
171    //XH1                                          
172    theVector = (G4PhysicsFreeVector*) (*fHardC    
173    val =  G4Log(std::max(XH1,1e-42*eV*cm2)); /    
174    theVector->PutValues(binNumber,logEne,val);    
175                                                   
176     //XH2                                         
177    theVector = (G4PhysicsFreeVector*) (*fHardC    
178    val =  G4Log(std::max(XH2,1e-42*eV*eV*cm2))    
179    theVector->PutValues(binNumber,logEne,val);    
180                                                   
181    return;                                        
182 }                                                 
183                                                   
184 //....oooOO0OOooo........oooOO0OOooo........oo    
185                                                   
186 void G4PenelopeCrossSection::AddShellCrossSect    
187                    size_t shellID,                
188                    G4double energy,               
189                    G4double xs)                   
190 {                                                 
191   if (!fShellCrossSections)                       
192     {                                             
193       G4cout << "Something wrong in G4Penelope    
194   G4endl;                                         
195       G4cout << "Trying to fill un-initialized    
196       return;                                     
197     }                                             
198                                                   
199   if (shellID >= fNumberOfShells)                 
200     {                                             
201       G4cout << "Something wrong in G4Penelope    
202   G4endl;                                         
203       G4cout << "Trying to fill shell #" << sh    
204        <<  fNumberOfShells-1 << G4endl;           
205       return;                                     
206     }                                             
207                                                   
208   //fill vector                                   
209   G4PhysicsFreeVector* theVector = (G4PhysicsF    
210                                                   
211   if (binNumber >= fNumberOfEnergyPoints)         
212     {                                             
213       G4cout << "Something wrong in G4Penelope    
214   G4endl;                                         
215       G4cout << "Trying to register more point    
216       return;                                     
217     }                                             
218    G4double logEne = G4Log(energy);               
219    G4double val = G4Log(std::max(xs,1e-42*cm2)    
220    theVector->PutValues(binNumber,logEne,val);    
221                                                   
222    return;                                        
223 }                                                 
224                                                   
225 //....oooOO0OOooo........oooOO0OOooo........oo    
226                                                   
227 G4double G4PenelopeCrossSection::GetTotalCross    
228 {                                                 
229   G4double result = 0;                            
230   //take here XS0 + XH0                           
231   if (!fSoftCrossSections || !fHardCrossSectio    
232     {                                             
233       G4cout << "Something wrong in G4Penelope    
234   G4endl;                                         
235       G4cout << "Trying to retrieve from un-in    
236       return result;                              
237     }                                             
238                                                   
239   // 1) soft part                                 
240   G4PhysicsFreeVector* theVector = (G4PhysicsF    
241   if (theVector->GetVectorLength() < fNumberOf    
242     {                                             
243       G4cout << "Something wrong in G4Penelope    
244   G4endl;                                         
245       G4cout << "Soft cross section table look    
246       return result;                              
247     }                                             
248   G4double logene = G4Log(energy);                
249   G4double logXS = theVector->Value(logene);      
250   G4double softXS = G4Exp(logXS);                 
251                                                   
252    // 2) hard part                                
253   theVector = (G4PhysicsFreeVector*) (*fHardCr    
254   if (theVector->GetVectorLength() < fNumberOf    
255     {                                             
256       G4cout << "Something wrong in G4Penelope    
257   G4endl;                                         
258       G4cout << "Hard cross section table look    
259       return result;                              
260     }                                             
261   logXS = theVector->Value(logene);               
262   G4double hardXS = G4Exp(logXS);                 
263                                                   
264   result = hardXS + softXS;                       
265   return result;                                  
266 }                                                 
267                                                   
268 //....oooOO0OOooo........oooOO0OOooo........oo    
269                                                   
270 G4double G4PenelopeCrossSection::GetHardCrossS    
271 {                                                 
272   G4double result = 0;                            
273   //take here XH0                                 
274   if (!fHardCrossSections)                        
275     {                                             
276       G4cout << "Something wrong in G4Penelope    
277   G4endl;                                         
278       G4cout << "Trying to retrieve from un-in    
279       return result;                              
280     }                                             
281                                                   
282   G4PhysicsFreeVector* theVector = (G4PhysicsF    
283   if (theVector->GetVectorLength() < fNumberOf    
284     {                                             
285       G4cout << "Something wrong in G4Penelope    
286   G4endl;                                         
287       G4cout << "Hard cross section table look    
288       return result;                              
289     }                                             
290   G4double logene = G4Log(energy);                
291   G4double logXS = theVector->Value(logene);      
292   result = G4Exp(logXS);                          
293                                                   
294   return result;                                  
295 }                                                 
296                                                   
297 //....oooOO0OOooo........oooOO0OOooo........oo    
298                                                   
299 G4double G4PenelopeCrossSection::GetSoftStoppi    
300 {                                                 
301   G4double result = 0;                            
302   //take here XH0                                 
303   if (!fSoftCrossSections)                        
304     {                                             
305       G4cout << "Something wrong in G4Penelope    
306   G4endl;                                         
307       G4cout << "Trying to retrieve from un-in    
308       return result;                              
309     }                                             
310                                                   
311   G4PhysicsFreeVector* theVector = (G4PhysicsF    
312   if (theVector->GetVectorLength() < fNumberOf    
313     {                                             
314       G4cout << "Something wrong in G4Penelope    
315   G4endl;                                         
316       G4cout << "Soft cross section table look    
317       return result;                              
318     }                                             
319   G4double logene = G4Log(energy);                
320   G4double logXS = theVector->Value(logene);      
321   result = G4Exp(logXS);                          
322                                                   
323   return result;                                  
324 }                                                 
325                                                   
326 //....oooOO0OOooo........oooOO0OOooo........oo    
327                                                   
328 G4double G4PenelopeCrossSection::GetShellCross    
329 {                                                 
330   G4double result = 0;                            
331   if (!fShellCrossSections)                       
332     {                                             
333       G4cout << "Something wrong in G4Penelope    
334   G4endl;                                         
335       G4cout << "Trying to retrieve from un-in    
336       return result;                              
337     }                                             
338   if (shellID >= fNumberOfShells)                 
339     {                                             
340       G4cout << "Something wrong in G4Penelope    
341   G4endl;                                         
342       G4cout << "Trying to retrieve shell #" <    
343        <<  fNumberOfShells-1 << G4endl;           
344       return result;                              
345     }                                             
346                                                   
347   G4PhysicsFreeVector* theVector = (G4PhysicsF    
348                                                   
349   if (theVector->GetVectorLength() < fNumberOf    
350     {                                             
351       G4cout << "Something wrong in G4Penelope    
352   G4endl;                                         
353       G4cout << "Shell cross section table loo    
354       return result;                              
355     }                                             
356   G4double logene = G4Log(energy);                
357   G4double logXS = theVector->Value(logene);      
358   result = G4Exp(logXS);                          
359                                                   
360   return result;                                  
361 }                                                 
362 //....oooOO0OOooo........oooOO0OOooo........oo    
363                                                   
364 G4double G4PenelopeCrossSection::GetNormalized    
365 {                                                 
366   G4double result = 0;                            
367   if (!fShellNormalizedCrossSections)             
368     {                                             
369       G4cout << "Something wrong in G4Penelope    
370   G4endl;                                         
371       G4cout << "Trying to retrieve from un-in    
372       return result;                              
373     }                                             
374                                                   
375   if (!fIsNormalized)                             
376     {                                             
377       G4cout << "Something wrong in G4Penelope    
378       G4cout << "The table of normalized cross    
379     }                                             
380                                                   
381   if (shellID >= fNumberOfShells)                 
382     {                                             
383       G4cout << "Something wrong in G4Penelope    
384   G4endl;                                         
385       G4cout << "Trying to retrieve shell #" <    
386        <<  fNumberOfShells-1 << G4endl;           
387       return result;                              
388     }                                             
389                                                   
390   const G4PhysicsFreeVector* theVector =          
391     (G4PhysicsFreeVector*) (*fShellNormalizedC    
392                                                   
393   if (theVector->GetVectorLength() < fNumberOf    
394     {                                             
395       G4cout << "Something wrong in G4Penelope    
396   G4endl;                                         
397       G4cout << "Shell cross section table loo    
398       return result;                              
399     }                                             
400   G4double logene = G4Log(energy);                
401   G4double logXS = theVector->Value(logene);      
402   result = G4Exp(logXS);                          
403                                                   
404   return result;                                  
405 }                                                 
406                                                   
407 //....oooOO0OOooo........oooOO0OOooo........oo    
408                                                   
409 void G4PenelopeCrossSection::NormalizeShellCro    
410 {                                                 
411   if (fIsNormalized) //already done!              
412     {                                             
413       G4cout << "G4PenelopeCrossSection::Norma    
414       G4cout << "already invoked. Ignore it" <    
415       return;                                     
416     }                                             
417                                                   
418   if (!fShellNormalizedCrossSections)             
419     {                                             
420       G4cout << "Something wrong in G4Penelope    
421   G4endl;                                         
422       G4cout << "Trying to retrieve from un-in    
423       return;                                     
424     }                                             
425                                                   
426   for (size_t i=0;i<fNumberOfEnergyPoints;i++)    
427     {                                             
428       //energy grid is the same for all shells    
429                                                   
430       //Recalculate manually the XS factor, to    
431       //underflows                                
432       G4double normFactor = 0.;                   
433       for (size_t shellID=0;shellID<fNumberOfS    
434   {                                               
435     G4PhysicsFreeVector* theVec =                 
436       (G4PhysicsFreeVector*) (*fShellCrossSect    
437                                                   
438     normFactor += G4Exp((*theVec)[i]);            
439   }                                               
440       G4double logNormFactor = G4Log(normFacto    
441       //Normalize                                 
442       for (size_t shellID=0;shellID<fNumberOfS    
443   {                                               
444    G4PhysicsFreeVector* theVec =                  
445       (G4PhysicsFreeVector*) (*fShellNormalize    
446    G4PhysicsFreeVector* theFullVec =              
447      (G4PhysicsFreeVector*) (*fShellCrossSecti    
448    G4double previousValue = (*theFullVec)[i];     
449    G4double logEnergy = theFullVec->GetLowEdge    
450    //log(XS/normFactor) = log(XS) - log(normFa    
451    theVec->PutValues(i,logEnergy,previousValue    
452   }                                               
453     }                                             
454   fIsNormalized = true;                           
455                                                   
456   return;                                         
457 }                                                 
458