108 std::vector<G4double*>::iterator
pos;
109 for (pos=
CST.begin(); pos<
CST.end(); pos++)
112 for (pos=
PAR.begin(); pos<
PAR.end(); pos++)
115 for (pos=
SST.begin(); pos<
SST.end(); pos++)
118 for (pos=
S1T.begin(); pos<
S1T.end(); pos++)
121 for (pos=
B1T.begin(); pos<
B1T.end(); pos++)
124 for (pos=
S2T.begin(); pos<
S2T.end(); pos++)
127 for (pos=
B2T.begin(); pos<
B2T.end(); pos++)
130 for (pos=
S3T.begin(); pos<
S3T.end(); pos++)
133 for (pos=
B3T.begin(); pos<
B3T.end(); pos++)
136 for (pos=
S4T.begin(); pos<
S4T.end(); pos++)
139 for (pos=
B4T.begin(); pos<
B4T.end(); pos++)
148 outFile <<
"G4ChipsProtonElasticXS provides the elastic cross\n"
149 <<
"section for proton nucleus scattering as a function of incident\n"
150 <<
"momentum. The cross section is calculated using M. Kossov's\n"
151 <<
"CHIPS parameterization of cross section data.\n";
207 if(lastCS<=0. && pEn>
lastTH)
228 colP.push_back(pMom);
311 G4int blast=
static_cast<int>(shift);
330 G4int blast=
static_cast<int>(shift);
371 const G4int n_npel=24;
372 const G4int n_ppel=32;
374 G4double np_el[n_npel]={12.,.05,.0001,5.,.35,6.75,.14,19.,.6,6.75,.14,13.,.14,.6,.00013,
375 75.,.001,7.2,4.32,.012,2.5,0.0,12.,.34};
378 G4double pp_el[n_ppel]={2.865,18.9,.6461,3.,9.,.425,.4276,.0022,5.,74.,3.,3.4,.2,.17,
379 .001,8.,.055,3.64,5.e-5,4000.,1500.,.46,1.2e6,3.5e6,5.e-5,1.e10,
380 8.5e8,1.e10,1.1,3.4e6,6.8e6,0.};
398 if ( tgZ == 0 && tgN == 1 )
400 for (
G4int ip=0; ip<n_npel; ip++)
lastPAR[ip]=np_el[ip];
403 else if ( tgZ == 1 && tgN == 0 )
405 for (
G4int ip=0; ip<n_ppel; ip++)
lastPAR[ip]=pp_el[ip];
430 lastPAR[1]=4.8*std::pow(a,1.14)/(1.+3.6/a3);
431 lastPAR[2]=1./(1.+4.E-3*a4)+2.
E-6*a3/(1.+1.3
E-6*a3);
433 lastPAR[4]=3.E-8*a3/(1.+4.E-7*a4);
434 lastPAR[5]=.07*asa/(1.+.009*a2);
435 lastPAR[6]=(3.+3.E-16*a20)/(1.+a20*(2.
E-16/a+3.
E-19*a));
436 lastPAR[7]=(5.E-9*a4*sa+.27/
a)/(1.+5.E16/a20)/(1.+6.E-9*a4)+.015/a2;
437 lastPAR[8]=(.001*a+.07/
a)/(1.+5.E13/a16+5.
E-7*a3)+.0003/sa;
445 lastPAR[11]=.7/(1.+4.e-12*a16);
446 lastPAR[12]=2.5/a8/(a4+1.e-16*a32);
453 lastPAR[18]=9.e-7/(1.+.035*a5);
454 lastPAR[19]=(42.+2.7e-11*a16)/(1.+.14*a);
458 lastPAR[22]=2.4e-3*a8/(1.+2.6e-4*a7);
459 lastPAR[23]=3.5e-36*a32*a8/(1.+5.e-15*a32/
a);
461 lastPAR[24]=1.e5/(a8+2.5e12/a16);
462 lastPAR[25]=8.e7/(a12+1.e-27*a28*a28);
471 lastPAR[32]=1.5e-4/(1.+5.e-6*a12);
479 lastPAR[38]=1.e8*std::exp(.32*asa);
480 lastPAR[39]=20.*std::exp(.45*asa);
482 lastPAR[41]=2.5e5*std::exp(.085*a3);
496 lastPAR[ 9]=4.5*std::pow(a,1.15);
497 lastPAR[10]=.06*std::pow(a,.6);
498 lastPAR[11]=.6*a/(1.+2.e15/a16);
499 lastPAR[12]=.17/(a+9.e5/a3+1.5e33/a32);
500 lastPAR[13]=(.001+7.e-11*a5)/(1.+4.4
e-11*a5);
501 lastPAR[14]=(p1a10*p1a10+2.e-29)/(1.+2.
e-22*a12);
503 lastPAR[15]=400./a12+2.e-22*a9;
504 lastPAR[16]=1.e-32*a12/(1.+5.e22/a14);
505 lastPAR[17]=1000./a2+9.5*sa*ssa;
506 lastPAR[18]=4.e-6*a*asa+1.e11/a16;
507 lastPAR[19]=(120./a+.002*a2)/(1.+2.e14/a16);
514 lastPAR[24]=.0011*asa/(1.+3.e34/a32/a4);
515 lastPAR[25]=1.e-5*a2+2.e14/a16;
516 lastPAR[26]=1.2e-11*a2/(1.+1.5e19/a12);
517 lastPAR[27]=.016*asa/(1.+5.e16/a16);
519 lastPAR[28]=.002*a4/(1.+7.e7/std::pow(a-6.83,14));
520 lastPAR[29]=2.e6/a6+7.2/std::pow(a,.11);
521 lastPAR[30]=11.*a3/(1.+7.e23/a16/a8);
525 lastPAR[33]=3.5e-4*a2/(1.+1.e8/a8);
528 lastPAR[36]=1.e-9/a+s4a16*s4a16;
545 lastPAR[50]=900.*sa/(1.+500./a3);
548 lastPAR[51]=1.e15+2.e27/a4/(1.+2.e-18*a16);
597 else G4cout<<
"*Warning*G4ChipsProtonElasticXS::GetPTables: PDG="<<PDG<<
", Z="
598 <<tgZ<<
", N="<<tgN<<
", i="<<ini<<
" > fin="<<fin<<
", LP="<<LP<<
" > ILP="
599 <<ILP<<
" nothing is done!"<<
G4endl;
601 else G4cout<<
"*Warning*G4ChipsProtonElasticXS::GetPTables: PDG="<<PDG<<
", Z="
602 <<tgZ<<
", N="<<tgN<<
", i="<<ini<<
">= max="<<
nPoints<<
", LP="<<LP
603 <<
" > ILP="<<ILP<<
", lPMax="<<
lPMax<<
" nothing is done!"<<
G4endl;
612 ed <<
"PDG = " << PDG <<
", Z = " << tgZ <<
", N = " << tgN
613 <<
", while it is defined only for PDG=2212 (p)" <<
G4endl;
614 G4Exception(
"G4ChipsProtonElasticXS::GetPTables()",
"HAD_CHPS_0000",
627 if(PDG!=2212)
G4cout<<
"**Warning*G4ChipsProtonElasticXS::GetExT:PDG="<<PDG<<
G4endl;
636 G4double R2=(1.-std::exp(-E2*E2*E2));
648 q2=-std::log(1.-ran)/
theB1;
654 q2=-std::log(1.-ran);
656 q2=std::pow(q2,third)/
theB2;
662 q2=-std::log(1.-ran)/
theB3;
676 if(a>6.5)E3*=tm2*tm2*tm2;
691 q2=-std::log(1.-ran)/
theB1;
692 if(std::fabs(tss)>1.
e-7) q2=(std::sqrt(
theB1*(
theB1+(tss+tss)*q2))-
theB1)/tss;
698 q2=-std::log(1.-ran)/
theB2;
700 if(a<6.5) q2=std::pow(q2,third);
701 else q2=std::pow(q2,fifth);
707 q2=-std::log(1.-ran)/
theB3;
709 if(a>6.5) q2=std::pow(q2,sevth);
715 q2=-std::log(1.-ran)/
theB4;
716 if(a<6.5) q2=lastTM-q2;
720 if(!(q2>=-1.||q2<=1.))
G4cout<<
"*NAN*G4QElasticCrossSect::GetExchangeT: -t="<<q2<<
G4endl;
733 if(
lastLP<-4.3)
return 0.;
740 ed <<
"PDG = " << PDG <<
", Z = " << tgZ <<
", N = " << tgN
741 <<
", while it is defined only for PDG=2212 (p)" <<
G4endl;
742 G4Exception(
"G4ChipsProtonElasticXS::GetSlope()",
"HAD_CHPS_0000",
761 if(PDG!=2212)
G4cout<<
"*Warning*G4ChipsProtonElasticXS::GetTabV:PDG="<<PDG<<
G4endl;
764 const G4bool isHeavyElementAllowed =
true;
765 if(tgZ<0 || ( !isHeavyElementAllowed && tgZ>92))
767 G4cout<<
"*Warning*G4QProtonElCS::GetTabValue: (1-92) No isotopes for Z="<<tgZ<<
G4endl;
782 if ( tgZ == 1 && tgN == 0 )
787 theS1=(lastPAR[9]+lastPAR[10]*dl2*dl2)/(1.+lastPAR[11]/p4/p)+
788 (lastPAR[12]/p2+lastPAR[13]*
p)/(p4+lastPAR[14]*sp);
789 theB1=lastPAR[15]*std::pow(p,lastPAR[16])/(1.+lastPAR[17]/p3);
790 theS2=lastPAR[18]+lastPAR[19]/(p4+lastPAR[20]*
p);
791 theB2=lastPAR[21]+lastPAR[22]/(p4+lastPAR[23]/
sp);
792 theS3=lastPAR[24]+lastPAR[25]/(p4*p4+lastPAR[26]*p2+lastPAR[27]);
793 theB3=lastPAR[28]+lastPAR[29]/(p4+lastPAR[30]);
798 return lastPAR[0]/p2s/(1.+lastPAR[7]/p2s)+(lastPAR[1]+lastPAR[2]*dl1*dl1+lastPAR[4]/p)
799 /(1.+lastPAR[5]*lp)/(1.+lastPAR[6]/p4);
827 lastPAR[40]/(1.+lastPAR[41]*std::pow(p,lastPAR[42])));
828 theB4=lastPAR[43]*pa/p2/(1.+pa*lastPAR[44]);
863 G4double tMid=std::sqrt(pP2+mProt2)*mProt-mProt2;
870 G4double mds=dmt*std::sqrt(pP2+mProt2)+mProt2+mt*mt;
871 return dmt*dmt*pP2/mds;
879 ed <<
"PDG = " << PDG <<
", Z = " << tgZ <<
", N = " << tgN
880 <<
", while it is defined only for p projectiles & Z_target>0" <<
G4endl;
881 G4Exception(
"G4ChipsProtonElasticXS::GetQ2max()",
"HAD_CHPS_0000",