65 const G4double HGeVSQ=gigaelectronvolt*gigaelectronvolt/2.;
67 const G4int n_kppel=35;
69 G4double kpp_el[n_kppel]={.7,.38,.0676,.0557,3.5,2.23,.7,.1,2.,1.,.372,5.,74.,3.,3.4,
70 .2,.17,.001,8.,.055,3.64,5.e-5,4000.,1500.,.46,1.2e6,3.5e6,
71 5.e-5,1.e10,8.5e8,1.e10,1.1,3.4e6,6.8e6,0.};
124 std::vector<G4double*>::iterator
pos;
125 for (pos=
CST.begin(); pos<
CST.end(); pos++)
128 for (pos=
PAR.begin(); pos<
PAR.end(); pos++)
131 for (pos=
SST.begin(); pos<
SST.end(); pos++)
134 for (pos=
S1T.begin(); pos<
S1T.end(); pos++)
137 for (pos=
B1T.begin(); pos<
B1T.end(); pos++)
140 for (pos=
S2T.begin(); pos<
S2T.end(); pos++)
143 for (pos=
B2T.begin(); pos<
B2T.end(); pos++)
146 for (pos=
S3T.begin(); pos<
S3T.end(); pos++)
149 for (pos=
B3T.begin(); pos<
B3T.end(); pos++)
152 for (pos=
S4T.begin(); pos<
S4T.end(); pos++)
155 for (pos=
B4T.begin(); pos<
B4T.end(); pos++)
213 if(lastCS<=0. && pEn>
lastTH)
234 colP.push_back(pMom);
317 G4int blast=
static_cast<int>(shift);
336 G4int blast=
static_cast<int>(shift);
390 if ( tgZ == 1 && tgN == 0 )
392 for (
G4int ip=0; ip<n_kppel; ip++)
lastPAR[ip]=kpp_el[ip];
416 lastPAR[0]=.06*asa/(1.+a*(.01+.1/ssa));
432 lastPAR[11]=.7/(1.+4.e-12*a16);
433 lastPAR[12]=2.5/a8/(a4+1.e-16*a32);
441 lastPAR[19]=(42.+2.7e-11*a16)/(1.+.14*a);
445 lastPAR[22]=2.4e-3*a8/(1.+2.6e-4*a7);
446 lastPAR[23]=3.5e-36*a32*a8/(1.+5.e-15*a32/
a);
448 lastPAR[24]=1.e5/(a8+2.5e12/a16);
449 lastPAR[25]=8.e7/(a12+1.e-27*a28*a28);
458 lastPAR[32]=1.5e-4/(1.+5.e-6*a12);
466 lastPAR[38]=1.e8*std::exp(.32*asa);
467 lastPAR[39]=20.*std::exp(.45*asa);
469 lastPAR[41]=2.5e5*std::exp(.085*a3);
483 lastPAR[ 9]=4.5*std::pow(a,1.15);
484 lastPAR[10]=.06*std::pow(a,.6);
485 lastPAR[11]=.6*a/(1.+2.e15/a16);
486 lastPAR[12]=.17/(a+9.e5/a3+1.5e33/a32);
487 lastPAR[13]=(.001+7.e-11*
a5)/(1.+4.4e-11*a5);
488 lastPAR[14]=(p1a10*p1a10+2.e-29)/(1.+2.e-22*a12);
490 lastPAR[15]=400./a12+2.e-22*a9;
491 lastPAR[16]=1.e-32*a12/(1.+5.e22/a14);
492 lastPAR[17]=1000./a2+9.5*sa*ssa;
493 lastPAR[18]=4.e-6*a*asa+1.e11/a16;
494 lastPAR[19]=(120./a+.002*
a2)/(1.+2.e14/a16);
502 lastPAR[25]=1.e-5*a2+2.e14/a16;
503 lastPAR[26]=1.2e-11*a2/(1.+1.5e19/a12);
504 lastPAR[27]=.016*asa/(1.+5.e16/a16);
506 lastPAR[28]=.002*a4/(1.+7.e7/std::pow(a-6.83,14));
507 lastPAR[29]=2.e6/a6+7.2/std::pow(a,.11);
508 lastPAR[30]=11.*a3/(1.+7.e23/a16/a8);
512 lastPAR[33]=3.5e-4*a2/(1.+1.e8/a8);
515 lastPAR[36]=1.e-9/a+s4a16*s4a16;
535 lastPAR[51]=1.e15+2.e27/a4/(1.+2.e-18*a16);
565 for(
G4int ip=ini; ip<=fin; ip++)
584 else G4cout<<
"*Warning*G4ChipsKaonPlusElasticXS::GetPTables: PDG="<<PDG
585 <<
", Z="<<tgZ<<
", N="<<tgN<<
", i="<<ini<<
" > fin="<<fin<<
", LP="<<LP
586 <<
" > ILP="<<ILP<<
" nothing is done!"<<
G4endl;
588 else G4cout<<
"*Warning*G4ChipsKaonPlusElasticXS::GetPTables: PDG="<<PDG
589 <<
", Z="<<tgZ<<
", N="<<tgN<<
", i="<<ini<<
">= max="<<
nPoints<<
", LP="<<LP
590 <<
" > ILP="<<ILP<<
", lPMax="<<
lPMax<<
" nothing is done!"<<
G4endl;
599 ed <<
"PDG = " << PDG <<
", Z = " << tgZ <<
", N = " << tgN
600 <<
", while it is defined only for PDG=321 (K+) " <<
G4endl;
601 G4Exception(
"G4ChipsKaonPlusElasticXS::GetPTables()",
"HAD_CHPS_0000",
610 if(PDG!=321)
G4cout<<
"*Warning*G4ChipsKaonPlusElasticXS::GetExT:PDG="<<PDG<<
G4endl;
619 G4double R2=(1.-std::exp(-E2*E2*E2));
631 q2=-std::log(1.-ran)/
theB1;
637 q2=-std::log(1.-ran);
639 q2=std::pow(q2,third)/
theB2;
645 q2=-std::log(1.-ran)/
theB3;
659 if(a>6.5)E3*=tm2*tm2*tm2;
674 q2=-std::log(1.-ran)/
theB1;
675 if(std::fabs(tss)>1.e-7) q2=(std::sqrt(
theB1*(
theB1+(tss+tss)*q2))-
theB1)/tss;
681 q2=-std::log(1.-ran)/
theB2;
683 if(a<6.5) q2=std::pow(q2,third);
684 else q2=std::pow(q2,fifth);
690 q2=-std::log(1.-ran)/
theB3;
692 if(a>6.5) q2=std::pow(q2,sevth);
698 q2=-std::log(1.-ran)/
theB4;
699 if(a<6.5) q2=lastTM-q2;
703 if(!(q2>=-1.||q2<=1.))
G4cout<<
"*NAN*G4QKaonPlusElasticCS::GetExchT: -t="<<q2<<
G4endl;
715 if(
lastLP<-4.3)
return 0.;
719 ed <<
"PDG = " << PDG <<
", Z = " << tgZ <<
", N = " << tgN
720 <<
", while it is defined only for PDG=321 (K+)" <<
G4endl;
721 G4Exception(
"G4ChipsKaonPlusElasticXS::GetSlope()",
"HAD_CHPS_0000",
739 if(PDG!=321)
G4cout<<
"*Warning*G4ChipsKaonPlusElasticXS::GetTaV:PDG="<<PDG<<
G4endl;
742 G4cout<<
"*Warning*G4QKaonPlusElasticCS::GetTabV:(1-92)NoIsotopes for Z="<<tgZ<<
G4endl;
757 if ( tgZ == 1 && tgN == 0 )
761 theS1=(lastPAR[12]+lastPAR[13]*dl2*dl2)/(1.+lastPAR[14]/p4/p)+
762 (lastPAR[15]/p2+lastPAR[16]*p)/(p4+lastPAR[17]*sp);
763 theB1=lastPAR[18]*std::pow(p,lastPAR[19])/(1.+lastPAR[20]/p3);
764 theS2=lastPAR[21]+lastPAR[22]/(p4+lastPAR[23]*p);
765 theB2=lastPAR[24]+lastPAR[25]/(p4+lastPAR[26]/
sp);
766 theS3=lastPAR[27]+lastPAR[28]/(p4*p4+lastPAR[29]*p2+lastPAR[30]);
767 theB3=lastPAR[31]+lastPAR[32]/(p4+lastPAR[33]);
773 return lastPAR[0]/(lastPAR[2]+
sqr(p-lastPAR[1]))+(lastPAR[3]*dp*dp+lastPAR[5])/
774 (1.-lastPAR[6]/sp+lastPAR[7]/p4)
775 +lastPAR[8]/(
sqr(p-lastPAR[9])+lastPAR[10]);
803 lastPAR[40]/(1.+lastPAR[41]*std::pow(p,lastPAR[42])));
804 theB4=lastPAR[43]*pa/p2/(1.+pa*lastPAR[44]);
825 return (lastPAR[0]*dlp*dlp+lastPAR[1]+lastPAR[2]/p2)/(1.+lastPAR[3]/p2/
sp);
840 G4double mds=dmt*std::sqrt(pP2+mK2)+mK2+mt*mt;
841 return dmt*dmt*pP2/mds;
846 ed <<
"PDG = " << PDG <<
",Z = " << tgZ <<
", N = " << tgN
847 <<
", while it is defined only for p projectiles & Z_target>0" <<
G4endl;
848 G4Exception(
"G4ChipsKaonPlusElasticXS::GetQ2max()",
"HAD_CHPS_0000",
G4ChipsKaonPlusElasticXS()
std::vector< G4double * > SST
std::ostringstream G4ExceptionDescription
virtual G4double GetIsoCrossSection(const G4DynamicParticle *, G4int tgZ, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
G4_DECLARE_XS_FACTORY(G4ChipsKaonPlusElasticXS)
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
G4double GetPTables(G4double lpP, G4double lPm, G4int PDG, G4int tZ, G4int tN)
std::vector< G4double * > S1T
std::vector< G4int > colZ
#define G4MUTEX_INITIALIZER
G4double CalculateCrossSection(G4bool CS, G4int F, G4int I, G4int pPDG, G4int Z, G4int N, G4double pP)
std::vector< G4double * > B2T
G4double GetTotalMomentum() const
std::vector< G4double * > S4T
G4double GetQ2max(G4int pPDG, G4int tgZ, G4int tgN, G4double pP)
G4IonTable * GetIonTable() const
std::vector< G4double * > B1T
G4GLOB_DLL std::ostream G4cout
std::vector< G4double * > CST
std::vector< G4double * > B4T
std::vector< G4double * > B3T
virtual G4bool IsIsoApplicable(const G4DynamicParticle *Pt, G4int Z, G4int A, const G4Element *elm, const G4Material *mat)
static const G4double A[nN]
std::vector< G4double * > S2T
std::vector< G4double > colP
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
G4double GetExchangeT(G4int tZ, G4int tN, G4int pPDG)
G4double GetPDGMass() const
static G4ParticleTable * GetParticleTable()
G4double GetSlope(G4int tZ, G4int tN, G4int pPDG)
virtual G4double GetChipsCrossSection(G4double momentum, G4int Z, G4int N, G4int pdg)
static const double gigaelectronvolt
static const double millibarn
G4double GetTabValues(G4double lp, G4int pPDG, G4int tgZ, G4int tgN)
~G4ChipsKaonPlusElasticXS()
std::vector< G4double > colTH
std::vector< G4double > colCS
std::vector< G4double * > S3T
std::vector< G4int > colN
static G4KaonPlus * KaonPlus()
std::vector< G4double > PIN
std::vector< G4double * > PAR
static const G4double pos