Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/inclxx/incl_physics/src/G4INCLCrossSectionsMultiPions.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/hadronic/models/inclxx/incl_physics/src/G4INCLCrossSectionsMultiPions.cc (Version 11.3.0) and /processes/hadronic/models/inclxx/incl_physics/src/G4INCLCrossSectionsMultiPions.cc (Version 3.2)


  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 // INCL++ intra-nuclear cascade model             
 27 // Alain Boudard, CEA-Saclay, France              
 28 // Joseph Cugnon, University of Liege, Belgium    
 29 // Jean-Christophe David, CEA-Saclay, France      
 30 // Pekka Kaitaniemi, CEA-Saclay, France, and H    
 31 // Sylvie Leray, CEA-Saclay, France               
 32 // Davide Mancusi, CEA-Saclay, France             
 33 //                                                
 34 #define INCLXX_IN_GEANT4_MODE 1                   
 35                                                   
 36 #include "globals.hh"                             
 37                                                   
 38 #include "G4INCLCrossSectionsMultiPions.hh"       
 39 #include "G4INCLKinematicsUtils.hh"               
 40 #include "G4INCLParticleTable.hh"                 
 41 #include "G4INCLLogger.hh"                        
 42 // #include <cassert>                             
 43                                                   
 44 namespace G4INCL {                                
 45                                                   
 46   template<G4int N>                               
 47     struct BystrickyEvaluator {                   
 48       static G4double eval(const G4double pLab    
 49         const G4double pMeV = pLab*1E3;           
 50         const G4double ekin=std::sqrt(Particle    
 51         const G4double xrat=ekin*oneOverThresh    
 52         const G4double x=std::log(xrat);          
 53         return HornerEvaluator<N>::eval(x, coe    
 54       }                                           
 55     };                                            
 56                                                   
 57   const G4int CrossSectionsMultiPions::nMaxPiN    
 58   const G4int CrossSectionsMultiPions::nMaxPiP    
 59                                                   
 60   const G4double CrossSectionsMultiPions::s11p    
 61   const G4double CrossSectionsMultiPions::s01p    
 62   const G4double CrossSectionsMultiPions::s01p    
 63   const G4double CrossSectionsMultiPions::s11p    
 64   const G4double CrossSectionsMultiPions::s12p    
 65   const G4double CrossSectionsMultiPions::s12p    
 66   const G4double CrossSectionsMultiPions::s12z    
 67   const G4double CrossSectionsMultiPions::s02p    
 68   const G4double CrossSectionsMultiPions::s02p    
 69   const G4double CrossSectionsMultiPions::s12m    
 70                                                   
 71   CrossSectionsMultiPions::CrossSectionsMultiP    
 72     s11pzHC(-2.228000000000294018,8.7560000000    
 73     s01ppHC(2.0570000000126518344,-6.029000000    
 74     s01pzHC(0.18030000000000441851,7.870099999    
 75     s11pmHC(0.20590000000000031866,3.345099999    
 76     s12pmHC(-0.77235999999999901328,4.26265999    
 77     s12ppHC(-0.75724999999999975664,2.09343999    
 78     s12zzHC(-0.89599999999996965072,7.88299999    
 79     s02pzHC(-1.0579999999999967036,11.11399999    
 80     s02pmHC(2.4009000000012553286,-7.768000000    
 81     s12mzHC(-0.21858699999999976269,1.91489999    
 82   {                                               
 83   }                                               
 84                                                   
 85   G4double CrossSectionsMultiPions::NNElastic(    
 86                                                   
 87     /* The NN cross section is parametrised as    
 88      * of one of the nucleons. For NDelta or D    
 89      * assumption is that the cross section is    
 90      * total CM energy*. Thus, we calculate s     
 91      * we convert this value to the lab moment    
 92      * an NN collision*.                          
 93      */                                           
 94     const G4double s = KinematicsUtils::square    
 95                                                   
 96     if(part1->isNucleon() && part2->isNucleon(    
 97       const G4int i = ParticleTable::getIsospi    
 98         + ParticleTable::getIsospin(part2->get    
 99       return NNElasticFixed(s, i);                
100     }                                             
101     else {  // Nucleon-Delta and Delta-Delta      
102       const G4double plab = 0.001*KinematicsUt    
103       if (plab < 0.440) {                         
104         return 34.*std::pow(plab/0.4, (-2.104)    
105       }                                           
106       else if (plab < 0.800) {                    
107         return 23.5+1000.*std::pow(plab-0.7, 4    
108       }                                           
109       else if (plab <= 2.0) {                     
110         return 1250./(50.+plab)-4.*std::pow(pl    
111       }                                           
112       else {                                      
113         return 77./(plab+1.5);                    
114       }                                           
115     }                                             
116   }                                               
117                                                   
118     G4double CrossSectionsMultiPions::NNElasti    
119                                                   
120       /* From NNElastic, with isospin fixed an    
121       */                                          
122                                                   
123       G4double plab = 0.001*KinematicsUtils::m    
124       G4double sigma = 0.;                        
125                                                   
126       if (i == 0) {  // pn                        
127         if (plab < 0.446) {                       
128           G4double alp=std::log(plab);            
129           sigma = 6.3555*std::exp(-3.2481*alp-    
130         }                                         
131         else if (plab < 0.851) {                  
132           sigma = 33.+196.*std::pow(std::fabs(    
133         }                                         
134         else if (plab <= 2.0) {                   
135           sigma = 31./std::sqrt(plab);            
136         }                                         
137         else {                                    
138           sigma = 77./(plab+1.5);                 
139         }                                         
140         //if(plab < 0.9 && plab > 0.802) sigma    
141         //if(plab < 1.4 && plab > 1.31) sigma     
142         return sigma;                             
143       }                                           
144       else {  // pp and nn                        
145         if (plab < 0.440) {                       
146           return 34.*std::pow(plab/0.4, (-2.10    
147         }                                         
148         else if (plab < 0.8067) {                 
149           return 23.5+1000.*std::pow(plab-0.7,    
150         }                                         
151         else if (plab <= 2.0) {                   
152           return 1250./(50.+plab)-4.*std::pow(    
153         }                                         
154         else if (plab <= 3.0956) {                
155           return 77./(plab+1.5);                  
156         }                                         
157         else {                                    
158           G4double alp=std::log(plab);            
159           return 11.2+25.5*std::pow(plab, -1.1    
160         }                                         
161       }                                           
162     }                                             
163                                                   
164     G4double CrossSectionsMultiPions::NNTot(Pa    
165                                                   
166         G4int i = ParticleTable::getIsospin(pa    
167         + ParticleTable::getIsospin(part2->get    
168                                                   
169         if(part1->isNucleon() && part2->isNucl    
170           const G4double s = KinematicsUtils::    
171           return NNTotFixed(s, i);                
172         }                                         
173         else if (part1->isDelta() && part2->is    
174             return elastic(part1, part2);         
175         }                                         
176         else {  // Nucleon-Delta                  
177             return NDeltaToNN(part1, part2) +     
178         }                                         
179     }                                             
180                                                   
181     G4double CrossSectionsMultiPions::NNTotFix    
182                                                   
183       /* From NNTot, with isospin fixed and fo    
184       */                                          
185                                                   
186       G4double plab = 0.001*KinematicsUtils::m    
187                                                   
188       if (i == 0) {  // pn                        
189         if (plab < 0.446) {                       
190           G4double alp=std::log(plab);            
191           return 6.3555*std::exp(-3.2481*alp-0    
192         }                                         
193         else if (plab < 1.0) {                    
194           return 33.+196.*std::sqrt(std::pow(s    
195         }                                         
196         else if (plab < 1.924) {                  
197           return 24.2+8.9*plab;                   
198         }                                         
199         else {                                    
200           G4double alp=std::log(plab);            
201           return 48.9-33.7*std::pow(plab, -3.0    
202         }                                         
203       }                                           
204       else {  // pp and nn                        
205         if (plab < 0.440) {                       
206           return 34.*std::pow(plab/0.4, (-2.10    
207         }                                         
208         else if (plab < 0.8734) {                 
209           return 23.5+1000.*std::pow(plab-0.7,    
210         }                                         
211         else if (plab < 1.5) {                    
212           return 23.5+24.6/(1.+std::exp(-10.*(    
213         }                                         
214         else if (plab < 3.0044) {                 
215           return 41.+60.*(plab-0.9)*std::exp(-    
216         }                                         
217         else {                                    
218           G4double alp=std::log(plab);            
219           return 45.6+219.*std::pow(plab, -4.2    
220         }                                         
221       }                                           
222     }                                             
223                                                   
224     G4double CrossSectionsMultiPions::NNInelas    
225                                                   
226       const G4double s = ener*ener;               
227       G4double sincl;                             
228                                                   
229       if (iso != 0) {                             
230         if(s>=4074595.287720512986) { // plab>    
231           sincl = NNTotFixed(s, 2)-NNElasticFi    
232         }                                         
233         else {                                    
234           sincl =  0. ;                           
235         }                                         
236       } else {                                    
237         if(s>=4074595.287720512986) { // plab>    
238           sincl = 2*(NNTotFixed(s, 0)-NNElasti    
239         }                                         
240         else {                                    
241           return 0. ;                             
242         }                                         
243       }                                           
244       if (sincl < 0.) sincl = 0.;                 
245       return sincl;                               
246     }                                             
247                                                   
248     G4double CrossSectionsMultiPions::NNOnePiO    
249                                                   
250         /* Article J. Physique 48 (1987)1901-1    
251          nucleon-cucleon inelastic total cross    
252          J. Bystricky, P. La France, F. Lehar,    
253          S11PZ= section pp->pp pi0                
254          S01PP= section pp->pn pi+                
255          S01PZ= section pn->pn pi0                
256          S11PM= section pn->pp pi-                
257          S= X-Section, 1st number : 1 if pp an    
258          2nd number = number of pions, PP= pi+    
259          */                                       
260                                                   
261         const G4double s = ener*ener;             
262         G4double plab = 0.001*KinematicsUtils:    
263                                                   
264         G4double snnpit1=0.;                      
265         G4double snnpit=0.;                       
266         G4double s11pz=0.;                        
267         G4double s01pp=0.;                        
268         G4double s01pz=0.;                        
269         G4double s11pm=0.;                        
270                                                   
271         if ((iso != 0) && (plab < 2.1989)) {      
272             snnpit = xsiso - NNTwoPi(ener, iso    
273             if (snnpit < 1.e-8) snnpit=0.;        
274             return snnpit;                        
275         }                                         
276         else if ((iso == 0) && (plab < 1.7369)    
277             snnpit = xsiso;                       
278             if (snnpit < 1.e-8) snnpit=0.;        
279             return snnpit;                        
280         }                                         
281                                                   
282 //s11pz                                           
283         if (plab > 18.) {                         
284             s11pz=55.185/std::pow((0.1412*plab    
285         }                                         
286         else if (plab > 13.9) {                   
287             G4double alp=std::log(plab);          
288             s11pz=6.67-13.3*std::pow(plab, -6.    
289         }                                         
290         else if (plab >= 0.7765) {                
291             const G4double b=BystrickyEvaluato    
292             s11pz=b*b;                            
293         }                                         
294 //s01pp                                           
295         if (plab >= 0.79624) {                    
296             const G4double b=BystrickyEvaluato    
297             s01pp=b*b;                            
298         }                                         
299                                                   
300 // channel T=1                                    
301         snnpit1=s11pz+s01pp;                      
302         if (snnpit1 < 1.e-8) snnpit1=0.;          
303         if (iso != 0) {                           
304             return snnpit1;                       
305         }                                         
306                                                   
307 //s01pz                                           
308         if (plab > 4.5) {                         
309             s01pz=15289.4/std::pow((11.573*pla    
310         }                                         
311         else if (plab >= 0.777) {                 
312             const G4double b=BystrickyEvaluato    
313             s01pz=b*b;                            
314         }                                         
315 //s11pm                                           
316         if (plab > 14.) {                         
317             s11pm=46.68/std::pow((0.2231*plab+    
318         }                                         
319         else if (plab >= 0.788) {                 
320             const G4double b=BystrickyEvaluato    
321             s11pm=b*b;                            
322         }                                         
323                                                   
324 // channel T=0                                    
325 //        snnpit=s01pz+2*s11pm-snnpit1; //modi    
326         snnpit = 2*(s01pz+2*s11pm)-snnpit1;       
327         if (snnpit < 1.e-8) snnpit=0.;            
328         return snnpit;                            
329     }                                             
330                                                   
331     G4double CrossSectionsMultiPions::NNTwoPi(    
332                                                   
333         /* Article J. Physique 48 (1987)1901-1    
334            J. Bystricky, P. La France, F. Leha    
335            S12PM : pp -> pp Pi+ Pi-               
336            S12ZZ : pp -> pp Pi0 Pi0               
337            S12PP : pp -> nn Pi+ Pi+               
338            S02PZ : pp -> pn Pi+ Pi0               
339            S02PM : pn -> pn Pi+ Pi-               
340            S12MZ : pn -> pp Pi- Pi0               
341         */                                        
342                                                   
343         const G4double s = ener*ener;             
344         G4double plab = 0.001*KinematicsUtils:    
345                                                   
346         G4double snn2pit=0.;                      
347         G4double s12pm=0.;                        
348         G4double s12pp=0.;                        
349         G4double s12zz=0.;                        
350         G4double s02pz=0.;                        
351         G4double s02pm=0.;                        
352         G4double s12mz=0.;                        
353                                                   
354         if (iso==0 && plab<3.33) {                
355             snn2pit = xsiso - NNOnePiOrDelta(e    
356             if (snn2pit < 1.e-8) snn2pit=0.;      
357             return snn2pit;                       
358         }                                         
359                                                   
360         if (iso != 0) {                           
361 //s12pm                                           
362          if (plab > 15.) {                        
363             s12pm=25.977/plab;                    
364          }                                        
365          else if (plab >= 1.3817) {               
366             const G4double b=BystrickyEvaluato    
367             s12pm=b*b;                            
368          }                                        
369 //s12pp                                           
370          if (plab > 10.) {                        
371             s12pp=141.505/std::pow((-0.1016*pl    
372          }                                        
373          else if (plab >= 1.5739) {               
374             const G4double b=BystrickyEvaluato    
375             s12pp=b*b;                            
376          }                                        
377         }                                         
378 //s12zz                                           
379         if (plab > 4.) {                          
380             s12zz=97.355/std::pow((1.1579*plab    
381         }                                         
382         else if (plab >= 1.72207) {               
383             const G4double b=BystrickyEvaluato    
384             s12zz=b*b;                            
385         }                                         
386 //s02pz                                           
387         if (plab > 4.5) {                         
388             s02pz=178.082/std::pow((0.2014*pla    
389         }                                         
390         else if (plab >= 1.5656) {                
391             const G4double b=BystrickyEvaluato    
392             s02pz=b*b;                            
393         }                                         
394                                                   
395 // channel T=1                                    
396         if (iso != 0) {                           
397             snn2pit=s12pm+s12pp+s12zz+s02pz;      
398             if (snn2pit < 1.e-8) snn2pit=0.;      
399             return snn2pit;                       
400         }                                         
401                                                   
402 //s02pm                                           
403         if (plab > 5.) {                          
404             s02pm=135.826/std::pow(plab,2);       
405         }                                         
406         else if (plab >= 1.21925) {               
407             const G4double b=BystrickyEvaluato    
408             s02pm=b*b;                            
409         }                                         
410 //s12mz                                           
411         if (plab >= 1.29269) {                    
412             const G4double b=BystrickyEvaluato    
413             s12mz=b*b;                            
414         }                                         
415                                                   
416 // channel T=0                                    
417 //        snn2pit=3*(0.5*s02pm+0.5*s12mz-0.5*s    
418         snn2pit=3*(s02pm+0.5*s12mz-0.5*s02pz-s    
419         if (snn2pit < 1.e-8) snn2pit=0.;          
420         return snn2pit;                           
421     }                                             
422                                                   
423     G4double CrossSectionsMultiPions::NNThreeP    
424                                                   
425         const G4double s = ener*ener;             
426         G4double plab = 0.001*KinematicsUtils:    
427                                                   
428         G4double snn3pit=0.;                      
429                                                   
430         if (iso == 0) {                           
431 // channel T=0                                    
432             if (plab > 7.2355) {                  
433                 return 46.72/std::pow((plab -     
434             }                                     
435             else {                                
436                 snn3pit=xsiso-xs1pi-xs2pi;        
437                 if (snn3pit < 1.e-8) snn3pit=0    
438                 return snn3pit;                   
439             }                                     
440         }                                         
441         else {                                    
442 // channel T=1                                    
443             if (plab > 7.206) {                   
444                 return 5592.92/std::pow((plab+    
445             }                                     
446             else if (plab > 2.1989){              
447                 snn3pit=xsiso-xs1pi-xs2pi;        
448                 if (snn3pit < 1.e-8) snn3pit=0    
449                 return snn3pit;                   
450             }                                     
451             else return snn3pit;                  
452         }                                         
453     }                                             
454                                                   
455     G4double CrossSectionsMultiPions::NNOnePi(    
456         // Cross section for nucleon-nucleon d    
457                                                   
458         const G4int iso=ParticleTable::getIsos    
459         if (iso!=0) // If pp or nn we choose t    
460           return 0.;                              
461                                                   
462         const G4double ener=KinematicsUtils::t    
463                                                   
464         const G4double xsiso2=NNInelasticIso(e    
465         const G4double xsiso0=NNInelasticIso(e    
466         return 0.25*(NNOnePiOrDelta(ener, 0, x    
467     }                                             
468                                                   
469     G4double CrossSectionsMultiPions::NNOnePiO    
470         // Cross section for nucleon-nucleon d    
471         const G4double ener=KinematicsUtils::t    
472         const G4int iso=ParticleTable::getIsos    
473                                                   
474         const G4double xsiso2=NNInelasticIso(e    
475         if (iso != 0)                             
476           return NNOnePiOrDelta(ener, iso, xsi    
477         else {                                    
478           const G4double xsiso0=NNInelasticIso    
479           return 0.5*(NNOnePiOrDelta(ener, 0,     
480         }                                         
481     }                                             
482                                                   
483     G4double CrossSectionsMultiPions::NNTwoPi(    
484         //                                        
485         //     Nucleon-Nucleon producing one p    
486         //                                        
487         const G4double ener=KinematicsUtils::t    
488         const G4int iso=ParticleTable::getIsos    
489                                                   
490                                                   
491         const G4double xsiso2=NNInelasticIso(e    
492         if (iso != 0) {                           
493             return NNTwoPi(ener, 2, xsiso2);      
494         }                                         
495         else {                                    
496             const G4double xsiso0=NNInelasticI    
497             return 0.5*(NNTwoPi(ener, 0, xsiso    
498         }                                         
499         return 0.0; // Should never reach this    
500     }                                             
501                                                   
502     G4double CrossSectionsMultiPions::NNThreeP    
503         //                                        
504         //     Nucleon-Nucleon producing one p    
505         //                                        
506                                                   
507         const G4double ener=KinematicsUtils::t    
508         const G4int iso=ParticleTable::getIsos    
509                                                   
510                                                   
511         const G4double xsiso2=NNInelasticIso(e    
512         const G4double xs1pi2=NNOnePiOrDelta(e    
513         const G4double xs2pi2=NNTwoPi(ener, 2,    
514         if (iso != 0)                             
515           return NNThreePi(ener, 2, xsiso2, xs    
516         else {                                    
517           const G4double xsiso0=NNInelasticIso    
518           const G4double xs1pi0=NNOnePiOrDelta    
519           const G4double xs2pi0=NNTwoPi(ener,     
520           return 0.5*(NNThreePi(ener, 0, xsiso    
521         }                                         
522     }                                             
523                                                   
524     G4double CrossSectionsMultiPions::NNFourPi    
525       const G4double s = KinematicsUtils::squa    
526       if(s<6.25E6)                                
527         return 0.;                                
528       const G4double sigma = NNTot(particle1,     
529       return ((sigma>1.e-9) ? sigma : 0.);        
530     }                                             
531                                                   
532     G4double CrossSectionsMultiPions::NNToxPiN    
533       //                                          
534       //     Nucleon-Nucleon producing xpi pio    
535       //                                          
536 // assert(xpi>0 && xpi<=nMaxPiNN);                
537 // assert(particle1->isNucleon() && particle2-    
538                                                   
539       if (xpi == 1)                               
540         return NNOnePi(particle1, particle2);     
541       else if (xpi == 2)                          
542         return NNTwoPi(particle1, particle2);     
543       else if (xpi == 3)                          
544         return NNThreePi(particle1, particle2)    
545       else if (xpi == 4)                          
546         return NNFourPi(particle1, particle2);    
547       else // should never reach this point       
548         return 0.;                                
549     }                                             
550                                                   
551                                                   
552   G4double CrossSectionsMultiPions::spnPiPlusP    
553     // HE and LE pi- p and pi+ n                  
554     G4double ramass = 0.0;                        
555                                                   
556     if(x <= 1306.78) {                            
557        G4double y = x*x;                          
558        G4double q2;                               
559        q2=(y-std::pow(1076.0, 2))*(y-std::pow(    
560        if (q2 > 0.) {                             
561           G4double q3=std::pow(q2, 3./2.);        
562           G4double f3=q3/(q3+std::pow(180.0, 3    
563     G4double sdel;                                
564     sdel=326.5/(std::pow((x-1215.0-ramass)*2.0    
565     return sdel*f3*(1.0-5.0*ramass/1215.0);       
566        }                                          
567        else {                                     
568           return 0;                               
569        }                                          
570     }                                             
571     if(x <= 1754.0) {                             
572       return -2.33730e-06*std::pow(x, 3)+1.138    
573         -1.83993e+01*x+9893.4;                    
574     } else if (x <= 2150.0) {                     
575       return 1.13531e-06*std::pow(x, 3)-6.9169    
576         +1.39907e+01*x-9360.76;                   
577     } else {                                      
578       return -3.18087*std::log(x)+52.9784;        
579     }                                             
580   }                                               
581                                                   
582   G4double CrossSectionsMultiPions::spnPiMinus    
583     // HE pi- p and pi+ n                         
584     G4double ramass = 0.0;                        
585                                                   
586     if(x <= 1275.8) {                             
587        G4double y = x*x;                          
588        G4double q2;                               
589        q2=(y-std::pow(1076.0, 2))*(y-std::pow(    
590        if (q2 > 0.) {                             
591           G4double q3=std::pow(q2, 3./2.);        
592           G4double f3=q3/(q3+std::pow(180.0, 3    
593     G4double sdel;                                
594     sdel=326.5/(std::pow((x-1215.0-ramass)*2.0    
595     return sdel*f3*(1.0-5.0*ramass/1215.0)/3.;    
596        }                                          
597        else {                                     
598           return 0;                               
599        }                                          
600     }                                             
601     if(x <= 1495.0) {                             
602       return 0.00120683*(x-1372.52)*(x-1372.52    
603     } else if(x <= 1578.0) {                      
604       return 1.15873e-05*x*x+49965.6/((x-1519.    
605     } else if(x <= 2028.4) {                      
606       return 34.0248+43262.2/((x-1681.65)*(x-1    
607     } else if(x <= 7500.0) {                      
608       return 3.3e-7*(x-7500.0)*(x-7500.0)+24.5    
609     } else {                                      
610       return 24.5;                                
611     }                                             
612   }                                               
613                                                   
614   G4double CrossSectionsMultiPions::total(Part    
615     G4double inelastic;                           
616     if(p1->isNucleon() && p2->isNucleon()) {      
617       return NNTot(p1, p2);                       
618     } else if((p1->isNucleon() && p2->isDelta(    
619               (p1->isDelta() && p2->isNucleon(    
620       inelastic = NDeltaToNN(p1, p2);             
621     } else if((p1->isNucleon() && p2->isPion()    
622               (p1->isPion() && p2->isNucleon()    
623       return piNTot(p1,p2);                       
624     } else {                                      
625       inelastic = 0.;                             
626     }                                             
627                                                   
628     return inelastic + elastic(p1, p2);           
629   }                                               
630                                                   
631                                                   
632   G4double CrossSectionsMultiPions::piNIne(Par    
633     //      piN inelastic cross section (Delta    
634                                                   
635     const Particle *pion;                         
636     const Particle *nucleon;                      
637     if(particle1->isNucleon()) {                  
638       nucleon = particle1;                        
639       pion = particle2;                           
640     } else {                                      
641       pion = particle1;                           
642       nucleon = particle2;                        
643     }                                             
644 // assert(pion->isPion());                        
645                                                   
646     const G4double pLab = KinematicsUtils::mom    
647                                                   
648     // these limits correspond to sqrt(s)=1230    
649     if(pLab>212677. || pLab<296.367)              
650       return 0.0;                                 
651                                                   
652     const G4int ipit3 = ParticleTable::getIsos    
653     const G4int ind2t3 = ParticleTable::getIso    
654     const G4int cg = 4 + ind2t3*ipit3;            
655 // assert(cg==2 || cg==4 || cg==6);               
656                                                   
657 //    const G4double p1=1e-3*pLab;                
658 //    const G4double p2=std::log(p1);             
659     G4double xpipp = 0.0;                         
660     G4double xpimp = 0.0;                         
661                                                   
662     if(cg!=2) {                                   
663       // x-section pi+ p inelastique :            
664       xpipp=piPluspIne(pion,nucleon);             
665                                                   
666       if(cg==6) // cas pi+ p et pi- n             
667         return xpipp;                             
668     }                                             
669                                                   
670     // x-section pi- p inelastique :              
671     xpimp=piMinuspIne(pion,nucleon);              
672                                                   
673     if(cg==2) // cas pi- p et pi+ n               
674       return xpimp;                               
675     else      // cas pi0 p et pi0 n               
676       return 0.5*(xpipp+xpimp);                   
677   }                                               
678                                                   
679   G4double CrossSectionsMultiPions::piNToDelta    
680     //      piN Delta production                  
681                                                   
682     G4double x = KinematicsUtils::totalEnergyI    
683     if(x>20000.) return 0.0; // no cross secti    
684                                                   
685     G4int ipit3 = 0;                              
686     G4int ind2t3 = 0;                             
687     const G4double ramass = 0.0;                  
688                                                   
689     if(particle1->isPion()) {                     
690       ipit3 = ParticleTable::getIsospin(partic    
691       ind2t3 = ParticleTable::getIsospin(parti    
692     } else if(particle2->isPion()) {              
693       ipit3 = ParticleTable::getIsospin(partic    
694       ind2t3 = ParticleTable::getIsospin(parti    
695     }                                             
696                                                   
697     const G4double y=x*x;                         
698     const G4double q2=(y-1076.0*1076.0)*(y-800    
699     if (q2 <= 0.) {                               
700       return 0.0;                                 
701     }                                             
702     const G4double q3 = std::pow(std::sqrt(q2)    
703     const G4double f3 = q3/(q3 + 5832000.); //    
704     G4double sdelResult = 326.5/(std::pow((x-1    
705     sdelResult = sdelResult*(1.0-5.0*ramass/12    
706     const G4int cg = 4 + ind2t3*ipit3;            
707     sdelResult = sdelResult*f3*cg/6.0;            
708                                                   
709     return sdelResult;                            
710   }                                               
711                                                   
712   G4double CrossSectionsMultiPions::piNTot(Par    
713     //      FUNCTION SPN(X,IND2T3,IPIT3,f17)      
714     // SIGMA(PI+ + P) IN THE (3,3) REGION         
715     // NEW FIT BY J.VANDERMEULEN  + FIT BY Th     
716     //                              CONST AT L    
717     //      COMMON/BL8/RATHR,RAMASS               
718     //      integer f17                           
719     // RATHR and RAMASS are always 0.0!!!         
720                                                   
721     G4double x = KinematicsUtils::totalEnergyI    
722                                                   
723     G4int ipit3 = 0;                              
724     G4int ind2t3 = 0;                             
725                                                   
726     if(particle1->isPion()) {                     
727       ipit3 = ParticleTable::getIsospin(partic    
728       ind2t3 = ParticleTable::getIsospin(parti    
729     } else if(particle2->isPion()) {              
730       ipit3 = ParticleTable::getIsospin(partic    
731       ind2t3 = ParticleTable::getIsospin(parti    
732     }                                             
733                                                   
734     G4double spnResult=0.0;                       
735                                                   
736     // HE pi+ p and pi- n                         
737       if((ind2t3 == 1 && ipit3 == 2) || (ind2t    
738         spnResult=spnPiPlusPHE(x);                
739       else if((ind2t3 == 1 && ipit3 == -2) ||     
740         spnResult=spnPiMinusPHE(x);               
741       else if(ipit3 == 0) spnResult = (spnPiPl    
742       else {                                      
743         INCL_ERROR("Unknown configuration!\n"     
744       }                                           
745                                                   
746     return spnResult;                             
747   }                                               
748                                                   
749   G4double CrossSectionsMultiPions::NDeltaToNN    
750     const G4int isospin = ParticleTable::getIs    
751     if(isospin==4 || isospin==-4) return 0.0;     
752                                                   
753     G4double s = KinematicsUtils::squareTotalE    
754     G4double Ecm = std::sqrt(s);                  
755     G4int deltaIsospin;                           
756     G4double deltaMass;                           
757     if(p1->isDelta()) {                           
758       deltaIsospin = ParticleTable::getIsospin    
759       deltaMass = p1->getMass();                  
760     } else {                                      
761       deltaIsospin = ParticleTable::getIsospin    
762       deltaMass = p2->getMass();                  
763     }                                             
764                                                   
765     if(Ecm <= 938.3 + deltaMass) {                
766       return 0.0;                                 
767     }                                             
768                                                   
769     if(Ecm < 938.3 + deltaMass + 2.0) {           
770       Ecm = 938.3 + deltaMass + 2.0;              
771       s = Ecm*Ecm;                                
772     }                                             
773                                                   
774     const G4double x = (s - 4.*ParticleTable::    
775       (s - std::pow(ParticleTable::effectiveNu    
776     const G4double y = s/(s - std::pow(deltaMa    
777     /* Concerning the way we calculate the lab    
778      * in CrossSections::elasticNNLegacy().       
779      */                                           
780     G4double sDelta;                              
781     const G4double xsiso2=NNInelasticIso(Ecm,     
782     if (isospin != 0)                             
783       sDelta = NNOnePiOrDelta(Ecm, isospin, xs    
784     else {                                        
785       const G4double xsiso0=NNInelasticIso(Ecm    
786       sDelta = 0.25*(NNOnePiOrDelta(Ecm, 0, xs    
787     }                                             
788     G4double result = 0.5 * x * y * sDelta;       
789     /* modification for pion-induced cascade (    
790      * result=3.*result                           
791      * pi absorption increased also for intern    
792      */                                           
793     result *= 3.*(32.0 + isospin * isospin * (    
794     result /= 1.0 + 0.25 * (isospin * isospin)    
795     return result;                                
796   }                                               
797                                                   
798   G4double CrossSectionsMultiPions::NNToNDelta    
799 // assert(p1->isNucleon() && p2->isNucleon());    
800     const G4int isospin = ParticleTable::getIs    
801     G4double sigma = NNOnePiOrDelta(p1, p2);      
802     if(isospin==0)                                
803       sigma *= 0.5;                               
804     return sigma;                                 
805   }                                               
806                                                   
807   G4double CrossSectionsMultiPions::elastic(Pa    
808 //    if(!p1->isPion() && !p2->isPion()){         
809     if((p1->isNucleon()||p1->isDelta()) && (p2    
810       return NNElastic(p1, p2);                   
811       }                                           
812 //    else if (p1->isNucleon() || p2->isNucleo    
813   else if ((p1->isNucleon() && p2->isPion()) |    
814       G4double pielas = piNTot(p1,p2) - piNIne    
815         if (pielas < 0.){                         
816             pielas = 0.;                          
817         }                                         
818 //        return piNTot(p1,p2) - piNIne(p1,p2)    
819         return pielas;                            
820       }                                           
821     else {                                        
822        return 0.0;                                
823       }                                           
824   }                                               
825                                                   
826   G4double CrossSectionsMultiPions::calculateN    
827     G4double x = 0.001 * pl; // Change to GeV     
828     if(iso != 0) {                                
829       if(pl <= 2000.0) {                          
830         x = std::pow(x, 8);                       
831         return 5.5e-6 * x/(7.7 + x);              
832       } else {                                    
833         return (5.34 + 0.67*(x - 2.0)) * 1.0e-    
834       }                                           
835     } else {                                      
836       if(pl < 800.0) {                            
837         G4double b = (7.16 - 1.63*x) * 1.0e-6;    
838         return b/(1.0 + std::exp(-(x - 0.45)/0    
839       } else if(pl < 1100.0) {                    
840         return (9.87 - 4.88 * x) * 1.0e-6;        
841       } else {                                    
842         return (3.68 + 0.76*x) * 1.0e-6;          
843       }                                           
844     }                                             
845     return 0.0; // Should never reach this poi    
846   }                                               
847                                                   
848                                                   
849     G4double CrossSectionsMultiPions::piNToxPi    
850         //                                        
851         //     pion-Nucleon producing xpi pion    
852         //                                        
853     const Particle *pion;                         
854     const Particle *nucleon;                      
855     if(particle1->isNucleon()) {                  
856       nucleon = particle1;                        
857       pion = particle2;                           
858     } else {                                      
859       pion = particle1;                           
860       nucleon = particle2;                        
861     }                                             
862 // assert(xpi>1 && xpi<=nMaxPiPiN);               
863 // assert((particle1->isNucleon() && particle2    
864         const G4double plab = KinematicsUtils:    
865     if (xpi == 2) {                               
866       G4double OnePi=piNOnePi(particle1,partic    
867       if (OnePi < 1.e-09) OnePi = 0.;             
868             return OnePi;                         
869         }                                         
870         else if (xpi == 3){                       
871       G4double TwoPi=piNTwoPi(particle1,partic    
872       if (TwoPi < 1.e-09) TwoPi = 0.;             
873             return TwoPi;                         
874         }                                         
875         else if (xpi == 4) {                      
876             G4double piNThreePi = piNIne(parti    
877             if (piNThreePi < 1.e-09 || plab <     
878             return piNThreePi;                    
879         } else // should never reach this poin    
880           return 0.0;                             
881     }                                             
882                                                   
883   G4double CrossSectionsMultiPions::piNOnePi(P    
884     const Particle *pion;                         
885     const Particle *nucleon;                      
886     if(particle1->isNucleon()) {                  
887       nucleon = particle1;                        
888       pion = particle2;                           
889     } else {                                      
890       pion = particle1;                           
891       nucleon = particle2;                        
892     }                                             
893 // assert(pion->isPion());                        
894                                                   
895     const G4double pLab = KinematicsUtils::mom    
896                                                   
897     // this limit corresponds to sqrt(s)=1230     
898     if(pLab<296.367)                              
899       return 0.0;                                 
900                                                   
901     const G4int ipi = ParticleTable::getIsospi    
902     const G4int ind2 = ParticleTable::getIsosp    
903     const G4int cg = 4 + ind2*ipi;                
904 // assert(cg==2 || cg==4 || cg==6);               
905                                                   
906     //  const G4double p1=1e-3*pLab;              
907     G4double tamp6=0.;                            
908     G4double tamp2=0.;                            
909     const G4double elas = elastic(particle1, p    
910                                                   
911     //   X-SECTION PI+ P INELASTIQUE :            
912     if(cg != 2) {                                 
913       tamp6=piPluspOnePi(particle1,particle2);    
914       if (cg == 6){ //   CAS PI+ P ET PI- N       
915         if(tamp6 >= elas && pLab < 410.) tamp6    
916         return tamp6;                             
917       }                                           
918     }                                             
919                                                   
920     //   X-SECTION PI- P INELASTIQUE :            
921     tamp2=piMinuspOnePi(particle1,particle2);     
922     if (tamp2 < 0.0) tamp2=0;                     
923                                                   
924     if (cg == 2) //   CAS PI- P ET PI+ N          
925       return tamp2;                               
926     else {       //   CAS PI0 P ET PI0 N          
927       G4double s1pin = 0.5*(tamp6+tamp2);         
928       const G4double inelastic = piNIne(partic    
929       if(s1pin >= elas && pLab < 410.) s1pin =    
930       if (s1pin > inelastic)                      
931         s1pin = inelastic;                        
932       return s1pin;                               
933     }                                             
934   }                                               
935                                                   
936   G4double CrossSectionsMultiPions::piNTwoPi(P    
937     //                                            
938     //     pion-nucleon interaction, producing    
939     //     fit from Landolt-Bornstein multipli    
940     //                                            
941                                                   
942     const Particle *pion;                         
943     const Particle *nucleon;                      
944     if(particle1->isNucleon()) {                  
945       nucleon = particle1;                        
946       pion = particle2;                           
947     } else {                                      
948       pion = particle1;                           
949       nucleon = particle2;                        
950     }                                             
951 // assert(pion->isPion());                        
952                                                   
953     const G4double pLab = KinematicsUtils::mom    
954     const G4double elas = elastic(pion, nucleo    
955                                                   
956     // this limit corresponds to sqrt(s)=1230     
957     if(pLab<296.367)                              
958       return 0.0;                                 
959                                                   
960     const G4int ipi = ParticleTable::getIsospi    
961     const G4int ind2 = ParticleTable::getIsosp    
962     const G4int cg = 4 + ind2*ipi;                
963 // assert(cg==2 || cg==4 || cg==6);               
964                                                   
965     G4double tamp6=0.;                            
966     G4double tamp2=0.;                            
967                                                   
968     //   X-SECTION PI+ P INELASTIQUE :            
969     if(cg!=2) {                                   
970       tamp6=piPluspTwoPi(particle1,particle2);    
971       if(cg==6){ //   CAS PI+ P ET PI- N          
972         if(tamp6 >= elas && pLab < 410.) tamp6    
973         return tamp6;}                            
974     }                                             
975                                                   
976     //   X-SECTION PI- P INELASTIQUE :            
977     tamp2=piMinuspTwoPi(particle1,particle2);     
978                                                   
979     if(cg==2) //   CAS PI- P ET PI+ N             
980       return tamp2;                               
981     else {    //   CAS PI0 P ET PI0 N             
982       const G4double s2pin=0.5*(tamp6+tamp2);     
983       return s2pin;                               
984     }                                             
985   }                                               
986                                                   
987   G4double CrossSectionsMultiPions::piPluspIne    
988     //      piPlusP inelastic cross section (D    
989                                                   
990     const Particle *pion;                         
991     const Particle *nucleon;                      
992     if(particle1->isNucleon()) {                  
993       nucleon = particle1;                        
994       pion = particle2;                           
995     } else {                                      
996       pion = particle1;                           
997       nucleon = particle2;                        
998     }                                             
999 // assert(pion->isPion());                        
1000                                                  
1001     const G4double pLab = KinematicsUtils::mo    
1002                                                  
1003     // these limits correspond to sqrt(s)=123    
1004     if(pLab>212677. || pLab<296.367)             
1005       return 0.0;                                
1006                                                  
1007 //    const G4int ipit3 = ParticleTable::getI    
1008 //    const G4int ind2t3 = ParticleTable::get    
1009 //    const G4int cg = 4 + ind2t3*ipit3;         
1010 //    assert(cg==2 || cg==4 || cg==6);           
1011                                                  
1012     const G4double p1=1e-3*pLab;                 
1013     const G4double p2=std::log(p1);              
1014     G4double xpipp = 0.0;                        
1015                                                  
1016     // x-section pi+ p inelastique :             
1017     if(p1<=0.75)                                 
1018       xpipp=17.965*std::pow(p1, 5.4606);         
1019     else                                         
1020       xpipp=24.3-12.3*std::pow(p1, -1.91)+0.3    
1021     // cas pi+ p et pi- n                        
1022     return xpipp;                                
1023                                                  
1024   }                                              
1025                                                  
1026   G4double CrossSectionsMultiPions::piMinuspI    
1027     //      piMinusp inelastic cross section     
1028                                                  
1029     const Particle *pion;                        
1030     const Particle *nucleon;                     
1031     if(particle1->isNucleon()) {                 
1032       nucleon = particle1;                       
1033       pion = particle2;                          
1034     } else {                                     
1035       pion = particle1;                          
1036       nucleon = particle2;                       
1037     }                                            
1038 // assert(pion->isPion());                       
1039                                                  
1040     const G4double pLab = KinematicsUtils::mo    
1041                                                  
1042     // these limits correspond to sqrt(s)=123    
1043     if(pLab>212677. || pLab<296.367)             
1044       return 0.0;                                
1045                                                  
1046 //    const G4int ipit3 = ParticleTable::getI    
1047 //    const G4int ind2t3 = ParticleTable::get    
1048 //    const G4int cg = 4 + ind2t3*ipit3;         
1049 //    assert(cg==2 || cg==4 || cg==6);           
1050                                                  
1051     const G4double p1=1e-3*pLab;                 
1052     const G4double p2=std::log(p1);              
1053     G4double xpimp = 0.0;                        
1054                                                  
1055     // x-section pi- p inelastique :             
1056     if(p1 <= 0.4731)                             
1057       xpimp=0;                                   
1058     else                                         
1059       xpimp=26.6-7.18*std::pow(p1, -1.86)+0.3    
1060     if(xpimp<0.)                                 
1061       xpimp=0;                                   
1062                                                  
1063     // cas pi- p et pi+ n                        
1064     return xpimp;                                
1065                                                  
1066   }                                              
1067                                                  
1068   G4double CrossSectionsMultiPions::piPluspOn    
1069     const Particle *pion;                        
1070     const Particle *nucleon;                     
1071     if(particle1->isNucleon()) {                 
1072       nucleon = particle1;                       
1073       pion = particle2;                          
1074     } else {                                     
1075       pion = particle1;                          
1076       nucleon = particle2;                       
1077     }                                            
1078 // assert(pion->isPion());                       
1079                                                  
1080     const G4double pLab = KinematicsUtils::mo    
1081                                                  
1082     // this limit corresponds to sqrt(s)=1230    
1083     if(pLab<296.367)                             
1084       return 0.0;                                
1085                                                  
1086     //  const G4int ipi = ParticleTable::getI    
1087     //  const G4int ind2 = ParticleTable::get    
1088     //  const G4int cg = 4 + ind2*ipi;           
1089     //  assert(cg==2 || cg==4 || cg==6);         
1090                                                  
1091     const G4double p1=1e-3*pLab;                 
1092     G4double tamp6=0.;                           
1093                                                  
1094     //   X-SECTION PI+ P INELASTIQUE :           
1095     if(pLab < 1532.52) // corresponds to sqrt    
1096       tamp6=piPluspIne(particle1, particle2);    
1097     else                                         
1098       tamp6=0.204+18.2*std::pow(p1, -1.72)+6.    
1099                                                  
1100     //   CAS PI+ P ET PI- N                      
1101     return tamp6;                                
1102                                                  
1103   }                                              
1104                                                  
1105   G4double CrossSectionsMultiPions::piMinuspO    
1106     const Particle *pion;                        
1107     const Particle *nucleon;                     
1108     if(particle1->isNucleon()) {                 
1109       nucleon = particle1;                       
1110       pion = particle2;                          
1111     } else {                                     
1112       pion = particle1;                          
1113       nucleon = particle2;                       
1114     }                                            
1115 // assert(pion->isPion());                       
1116                                                  
1117     const G4double pLab = KinematicsUtils::mo    
1118                                                  
1119     // this limit corresponds to sqrt(s)=1230    
1120     if(pLab<296.367)                             
1121       return 0.0;                                
1122                                                  
1123     //  const G4int ipi = ParticleTable::getI    
1124     //  const G4int ind2 = ParticleTable::get    
1125     //  const G4int cg = 4 + ind2*ipi;           
1126     //  assert(cg==2 || cg==4 || cg==6);         
1127                                                  
1128     const G4double p1=1e-3*pLab;                 
1129     G4double tamp2=0.;                           
1130                                                  
1131     //   X-SECTION PI- P INELASTIQUE :           
1132     if (pLab < 1228.06) // corresponds to sqr    
1133       tamp2=piMinuspIne(particle1, particle2)    
1134     else                                         
1135       tamp2=9.04*std::pow(p1, -1.17)+18.*std:    
1136     if (tamp2 < 0.0) tamp2=0;                    
1137                                                  
1138     //   CAS PI- P ET PI+ N                      
1139     return tamp2;                                
1140   }                                              
1141                                                  
1142   G4double CrossSectionsMultiPions::piPluspTw    
1143     //                                           
1144     //     pion-nucleon interaction, producin    
1145     //     fit from Landolt-Bornstein multipl    
1146     //                                           
1147                                                  
1148     const Particle *pion;                        
1149     const Particle *nucleon;                     
1150     if(particle1->isNucleon()) {                 
1151       nucleon = particle1;                       
1152       pion = particle2;                          
1153     } else {                                     
1154       pion = particle1;                          
1155       nucleon = particle2;                       
1156     }                                            
1157 // assert(pion->isPion());                       
1158                                                  
1159     const G4double pLab = KinematicsUtils::mo    
1160                                                  
1161     // this limit corresponds to sqrt(s)=1230    
1162     if(pLab<296.367)                             
1163       return 0.0;                                
1164                                                  
1165     //  const G4int ipi = ParticleTable::getI    
1166     //  const G4int ind2 = ParticleTable::get    
1167     //  const G4int cg = 4 + ind2*ipi;           
1168     //  assert(cg==2 || cg==4 || cg==6);         
1169                                                  
1170     const G4double p1=1e-3*pLab;                 
1171     G4double tamp6=0.;                           
1172                                                  
1173     //   X-SECTION PI+ P INELASTIQUE :           
1174     if(pLab < 2444.7) // corresponds to sqrt(    
1175       tamp6=piPluspIne(particle1, particle2)-    
1176     else                                         
1177       tamp6=1.59+25.5*std::pow(p1, -1.04); //    
1178                                                  
1179     //   CAS PI+ P ET PI- N                      
1180     return tamp6;                                
1181   }                                              
1182                                                  
1183     G4double CrossSectionsMultiPions::piMinus    
1184   //                                             
1185   //     pion-nucleon interaction, producing     
1186   //     fit from Landolt-Bornstein multiplie    
1187   //                                             
1188                                                  
1189   const Particle *pion;                          
1190   const Particle *nucleon;                       
1191   if(particle1->isNucleon()) {                   
1192     nucleon = particle1;                         
1193     pion = particle2;                            
1194   } else {                                       
1195     pion = particle1;                            
1196     nucleon = particle2;                         
1197   }                                              
1198 // assert(pion->isPion());                       
1199                                                  
1200   const G4double pLab = KinematicsUtils::mome    
1201                                                  
1202   // this limit corresponds to sqrt(s)=1230 M    
1203   if(pLab<296.367)                               
1204     return 0.0;                                  
1205                                                  
1206   //  const G4int ipi = ParticleTable::getIso    
1207   //  const G4int ind2 = ParticleTable::getIs    
1208   //  const G4int cg = 4 + ind2*ipi;             
1209   //  assert(cg==2 || cg==4 || cg==6);           
1210                                                  
1211   const G4double p1=1e-3*pLab;                   
1212   G4double tamp2=0.;                             
1213                                                  
1214   //   X-SECTION PI- P INELASTIQUE :             
1215   if(pLab<2083.63) // corresponds to sqrt(s)=    
1216     tamp2=piMinuspIne(particle1, particle2)-p    
1217   else                                           
1218     tamp2=2.457794117647+18.066176470588*std:    
1219                                                  
1220   //   CAS PI- P ET PI+ N                        
1221   return tamp2;                                  
1222 }                                                
1223                                                  
1224                                                  
1225                                                  
1226                                                  
1227     G4double CrossSectionsMultiPions::piNToEt    
1228     //                                           
1229     //     Pion-Nucleon producing Eta cross s    
1230     //                                           
1231         return 0.;                               
1232     }                                            
1233                                                  
1234     G4double CrossSectionsMultiPions::piNToOm    
1235     //                                           
1236     //     Pion-Nucleon producing Omega cross    
1237     //                                           
1238         return 0.;                               
1239     }                                            
1240                                                  
1241     G4double CrossSectionsMultiPions::piNToEt    
1242     //                                           
1243     //     Pion-Nucleon producing EtaPrime cr    
1244     //                                           
1245         return 0.;                               
1246     }                                            
1247                                                  
1248     G4double CrossSectionsMultiPions::etaNToP    
1249     //                                           
1250     //     Eta-Nucleon producing Pion cross s    
1251     //                                           
1252           return 0.;                             
1253     }                                            
1254                                                  
1255                                                  
1256      G4double CrossSectionsMultiPions::etaNTo    
1257     //                                           
1258     //     Eta-Nucleon producing Two Pions cr    
1259     //                                           
1260           return 0.;                             
1261      }                                           
1262                                                  
1263                                                  
1264     G4double CrossSectionsMultiPions::omegaNT    
1265     //                                           
1266     //     Omega-Nucleon producing Pion cross    
1267     //                                           
1268         return 0.;                               
1269     }                                            
1270                                                  
1271     G4double CrossSectionsMultiPions::omegaNT    
1272     //                                           
1273     //     Omega-Nucleon producing Two Pions     
1274     //                                           
1275         return 0.;                               
1276     }                                            
1277                                                  
1278     G4double CrossSectionsMultiPions::etaPrim    
1279     //                                           
1280     //     EtaPrime-Nucleon producing Pion cr    
1281     //                                           
1282         return 0.;                               
1283     }                                            
1284                                                  
1285     G4double CrossSectionsMultiPions::NNToNNE    
1286     //                                           
1287     //     Nucleon-Nucleon producing Eta cros    
1288     //                                           
1289         return 0.;                               
1290     }                                            
1291                                                  
1292   G4double CrossSectionsMultiPions::NNToNNEta    
1293     //                                           
1294     //     Nucleon-Nucleon producing Eta cros    
1295     //                                           
1296       return 0.;                                 
1297      }                                           
1298                                                  
1299   G4double CrossSectionsMultiPions::NNToNNEta    
1300       return 0.;                                 
1301      }                                           
1302                                                  
1303     G4double CrossSectionsMultiPions::NNToNDe    
1304     //                                           
1305     //     Nucleon-Nucleon producing N-Delta-    
1306     //                                           
1307     return 0.;                                   
1308     }                                            
1309                                                  
1310     G4double CrossSectionsMultiPions::NNToNNO    
1311     //                                           
1312     //     Nucleon-Nucleon producing Omega cr    
1313     //                                           
1314      return 0.;                                  
1315     }                                            
1316                                                  
1317     G4double CrossSectionsMultiPions::NNToNNO    
1318     //                                           
1319     //     Nucleon-Nucleon producing Omega cr    
1320     //                                           
1321      return 0.;                                  
1322     }                                            
1323                                                  
1324     G4double CrossSectionsMultiPions::NNToNNO    
1325      return 0.;                                  
1326     }                                            
1327                                                  
1328     G4double CrossSectionsMultiPions::NNToNDe    
1329   //                                             
1330   //     Nucleon-Nucleon producing N-Delta-Om    
1331   //                                             
1332      return 0.;                                  
1333     }                                            
1334                                                  
1335                                                  
1336                                                  
1337                                                  
1338     G4double CrossSectionsMultiPions::NYelast    
1339         //                                       
1340         //      Hyperon-Nucleon elastic cross    
1341         //                                       
1342     return 0.;                                   
1343     }                                            
1344                                                  
1345     G4double CrossSectionsMultiPions::NKelast    
1346         //                                       
1347         //      Kaon-Nucleon elastic cross se    
1348         //                                       
1349     return 0.;                                   
1350   }                                              
1351                                                  
1352     G4double CrossSectionsMultiPions::NKbelas    
1353         //                                       
1354         //      antiKaon-Nucleon elastic cros    
1355         //                                       
1356     return 0.;                                   
1357   }                                              
1358                                                  
1359                                                  
1360   G4double CrossSectionsMultiPions::NNToNLK(P    
1361         //                                       
1362         //      Nucleon-Nucleon producing N-L    
1363         //                                       
1364         return 0.;                               
1365     }                                            
1366                                                  
1367     G4double CrossSectionsMultiPions::NNToNSK    
1368         //                                       
1369         //      Nucleon-Nucleon producing N-S    
1370         //                                       
1371         return 0.;                               
1372     }                                            
1373                                                  
1374     G4double CrossSectionsMultiPions::NNToNLK    
1375         //                                       
1376         //      Nucleon-Nucleon producing N-L    
1377         //                                       
1378         return 0.;                               
1379     }                                            
1380                                                  
1381     G4double CrossSectionsMultiPions::NNToNSK    
1382         //                                       
1383         //      Nucleon-Nucleon producing N-S    
1384         //                                       
1385         return 0.;                               
1386     }                                            
1387                                                  
1388     G4double CrossSectionsMultiPions::NNToNLK    
1389         //                                       
1390         //     Nucleon-Nucleon producing N-La    
1391         //                                       
1392         return 0.;                               
1393     }                                            
1394                                                  
1395     G4double CrossSectionsMultiPions::NNToNSK    
1396         //                                       
1397         //      Nucleon-Nucleon producing N-S    
1398         //                                       
1399         return 0.;                               
1400     }                                            
1401                                                  
1402     G4double CrossSectionsMultiPions::NNToNNK    
1403         //                                       
1404         //      Nucleon-Nucleon producing Nuc    
1405         //                                       
1406         return 0.;                               
1407     }                                            
1408                                                  
1409     G4double CrossSectionsMultiPions::NNToMis    
1410         //                                       
1411         //      Nucleon-Nucleon missing stran    
1412         //                                       
1413         return 0.;                               
1414     }                                            
1415                                                  
1416     G4double CrossSectionsMultiPions::NDeltaT    
1417         // Nucleon-Delta producing Nucleon La    
1418         return 0;                                
1419     }                                            
1420     G4double CrossSectionsMultiPions::NDeltaT    
1421         // Nucleon-Delta producing Nucleon Si    
1422         return 0;                                
1423     }                                            
1424     G4double CrossSectionsMultiPions::NDeltaT    
1425         // Nucleon-Delta producing Delta Lamb    
1426         return 0;                                
1427     }                                            
1428     G4double CrossSectionsMultiPions::NDeltaT    
1429         // Nucleon-Delta producing Delta Sigm    
1430         return 0;                                
1431     }                                            
1432                                                  
1433     G4double CrossSectionsMultiPions::NDeltaT    
1434         // Nucleon-Delta producing Nucleon-Nu    
1435         return 0;                                
1436     }                                            
1437                                                  
1438                                                  
1439     G4double CrossSectionsMultiPions::NpiToLK    
1440         //                                       
1441         //      Pion-Nucleon producing Lambda    
1442         //                                       
1443         return 0.;                               
1444     }                                            
1445                                                  
1446     G4double CrossSectionsMultiPions::NpiToSK    
1447         //                                       
1448         //      Pion-Nucleon producing Sigma-    
1449         //                                       
1450         return 0.;                               
1451     }                                            
1452     G4double CrossSectionsMultiPions::p_pimTo    
1453         return 0.;                               
1454     }                                            
1455     G4double CrossSectionsMultiPions::p_pimTo    
1456         return 0.;                               
1457     }                                            
1458     G4double CrossSectionsMultiPions::p_pizTo    
1459         return 0.;                               
1460     }                                            
1461                                                  
1462     G4double CrossSectionsMultiPions::NpiToLK    
1463         //                                       
1464         //      Pion-Nucleon producing Lambda    
1465         //                                       
1466         return 0.;                               
1467     }                                            
1468                                                  
1469     G4double CrossSectionsMultiPions::NpiToSK    
1470         //                                       
1471         //      Pion-Nucleon producing Sigma-    
1472         //                                       
1473         return 0.;                               
1474     }                                            
1475                                                  
1476     G4double CrossSectionsMultiPions::NpiToLK    
1477         //                                       
1478         //      Pion-Nucleon producing Lambda    
1479         //                                       
1480         return 0.;                               
1481     }                                            
1482                                                  
1483     G4double CrossSectionsMultiPions::NpiToSK    
1484         //                                       
1485         //      Pion-Nucleon producing Lambda    
1486         //                                       
1487         return 0.;                               
1488     }                                            
1489                                                  
1490     G4double CrossSectionsMultiPions::NpiToNK    
1491         //                                       
1492         //      Pion-Nucleon producing Nucleo    
1493         //                                       
1494         return 0.;                               
1495     }                                            
1496                                                  
1497     G4double CrossSectionsMultiPions::NpiToMi    
1498         //                                       
1499         //      Pion-Nucleon missing strangen    
1500         //                                       
1501         return 0.;                               
1502     }                                            
1503                                                  
1504     G4double CrossSectionsMultiPions::NLToNS(    
1505         //                                       
1506         //      Nucleon-Hyperon multiplet cha    
1507         //                                       
1508         return 0.;                               
1509     }                                            
1510                                                  
1511     G4double CrossSectionsMultiPions::NSToNL(    
1512         //                                       
1513         //      Nucleon-Sigma quasi-elastic c    
1514         //                                       
1515         return 0.;                               
1516     }                                            
1517                                                  
1518     G4double CrossSectionsMultiPions::NSToNS(    
1519         //                                       
1520         //      Nucleon-Sigma quasi-elastic c    
1521         //                                       
1522         return 0.;                               
1523     }                                            
1524                                                  
1525     G4double CrossSectionsMultiPions::NKToNK(    
1526         //                                       
1527         //      Nucleon-Kaon quasi-elastic cr    
1528         //                                       
1529         return 0.;                               
1530     }                                            
1531                                                  
1532     G4double CrossSectionsMultiPions::NKToNKp    
1533         //                                       
1534         //      Nucleon-Kaon producing Nucleo    
1535         //                                       
1536         return 0.;                               
1537     }                                            
1538                                                  
1539     G4double CrossSectionsMultiPions::NKToNK2    
1540         //                                       
1541         //      Nucleon-Kaon producing Nucleo    
1542         //                                       
1543         return 0.;                               
1544     }                                            
1545                                                  
1546     G4double CrossSectionsMultiPions::NKbToNK    
1547         //                                       
1548         //      Nucleon-antiKaon quasi-elasti    
1549         //                                       
1550         return 0.;                               
1551     }                                            
1552                                                  
1553     G4double CrossSectionsMultiPions::NKbToSp    
1554         //                                       
1555         //      Nucleon-antiKaon producing Si    
1556         //                                       
1557         return 0.;                               
1558     }                                            
1559                                                  
1560     G4double CrossSectionsMultiPions::NKbToLp    
1561         //                                       
1562         //      Nucleon-antiKaon producing La    
1563         //                                       
1564         return 0.;                               
1565     }                                            
1566                                                  
1567     G4double CrossSectionsMultiPions::NKbToS2    
1568         //                                       
1569         //      Nucleon-antiKaon producing Si    
1570         //                                       
1571         return 0.;                               
1572     }                                            
1573                                                  
1574     G4double CrossSectionsMultiPions::NKbToL2    
1575         //                                       
1576         //      Nucleon-antiKaon producing La    
1577         //                                       
1578         return 0.;                               
1579     }                                            
1580                                                  
1581     G4double CrossSectionsMultiPions::NKbToNK    
1582         //                                       
1583         //      Nucleon-antiKaon producing Nu    
1584         //                                       
1585         return 0.;                               
1586     }                                            
1587                                                  
1588     G4double CrossSectionsMultiPions::NKbToNK    
1589         //                                       
1590         //      Nucleon-antiKaon producing Nu    
1591         //                                       
1592         return 0.;                               
1593     }                                            
1594                                                  
1595     G4double CrossSectionsMultiPions::NNbarEl    
1596         //                                       
1597         //     Nucleon-AntiNucleon to Nucleon    
1598         //                                       
1599       return 0.;                                 
1600     }                                            
1601                                                  
1602     G4double CrossSectionsMultiPions::NNbarCE    
1603         //                                       
1604         //     Nucleon-AntiNucleon charge exc    
1605         //                                       
1606       return 0.;                                 
1607     }                                            
1608                                                  
1609     G4double CrossSectionsMultiPions::NNbarTo    
1610         //                                       
1611         //     Nucleon-AntiNucleon to Lambda-    
1612         //                                       
1613       return 0.;                                 
1614     }                                            
1615                                                  
1616     G4double CrossSectionsMultiPions::NNbarTo    
1617         //                                       
1618         //     Nucleon-AntiNucleon to Nucleon    
1619         //                                       
1620       return 0.;                                 
1621     }                                            
1622                                                  
1623     G4double CrossSectionsMultiPions::NNbarTo    
1624         //                                       
1625         //     Nucleon-AntiNucleon to Nucleon    
1626         //                                       
1627       return 0.;                                 
1628     }                                            
1629                                                  
1630     G4double CrossSectionsMultiPions::NNbarTo    
1631         //                                       
1632         //     Nucleon-AntiNucleon to Nucleon    
1633         //                                       
1634       return 0.;                                 
1635     }                                            
1636                                                  
1637     G4double CrossSectionsMultiPions::NNbarTo    
1638         //                                       
1639         //     Nucleon-AntiNucleon total anni    
1640         //                                       
1641       return 0.;                                 
1642     }                                            
1643 } // namespace G4INCL                            
1644                                                  
1645