66 const G4int n_kmpel=36;
68 const G4double kmp_el[n_kmpel]={5.2,.0557,3.5,2.23,.7,.075,.004,.39,.000156,.15,1.,.0156,5.,
69 74.,3.,3.4,.2,.17,.001,8.,.055,3.64,5.e-5,4000.,1500.,.46,
70 1.2e6,3.5e6,5.e-5,1.e10,8.5e8,1.e10,1.1,3.4e6,6.8e6,0.};
73 const G4double HGeVSQ=gigaelectronvolt*gigaelectronvolt/2.;
122 std::vector<G4double*>::iterator
pos;
123 for (pos=
CST.begin(); pos<
CST.end(); pos++)
126 for (pos=
PAR.begin(); pos<
PAR.end(); pos++)
129 for (pos=
SST.begin(); pos<
SST.end(); pos++)
132 for (pos=
S1T.begin(); pos<
S1T.end(); pos++)
135 for (pos=
B1T.begin(); pos<
B1T.end(); pos++)
138 for (pos=
S2T.begin(); pos<
S2T.end(); pos++)
141 for (pos=
B2T.begin(); pos<
B2T.end(); pos++)
144 for (pos=
S3T.begin(); pos<
S3T.end(); pos++)
147 for (pos=
B3T.begin(); pos<
B3T.end(); pos++)
150 for (pos=
S4T.begin(); pos<
S4T.end(); pos++)
153 for (pos=
B4T.begin(); pos<
B4T.end(); pos++)
182 static G4ThreadLocal std::vector <G4int> *colN_G4MT_TLS_ = 0 ;
if (!colN_G4MT_TLS_) colN_G4MT_TLS_ =
new std::vector <G4int> ; std::vector <G4int> &colN = *colN_G4MT_TLS_;
183 static G4ThreadLocal std::vector <G4int> *colZ_G4MT_TLS_ = 0 ;
if (!colZ_G4MT_TLS_) colZ_G4MT_TLS_ =
new std::vector <G4int> ; std::vector <G4int> &colZ = *colZ_G4MT_TLS_;
184 static G4ThreadLocal std::vector <G4double> *colP_G4MT_TLS_ = 0 ;
if (!colP_G4MT_TLS_) colP_G4MT_TLS_ =
new std::vector <G4double> ; std::vector <G4double> &colP = *colP_G4MT_TLS_;
185 static G4ThreadLocal std::vector <G4double> *colTH_G4MT_TLS_ = 0 ;
if (!colTH_G4MT_TLS_) colTH_G4MT_TLS_ =
new std::vector <G4double> ; std::vector <G4double> &colTH = *colTH_G4MT_TLS_;
186 static G4ThreadLocal std::vector <G4double> *colCS_G4MT_TLS_ = 0 ;
if (!colCS_G4MT_TLS_) colCS_G4MT_TLS_ =
new std::vector <G4double> ; std::vector <G4double> &colCS = *colCS_G4MT_TLS_;
201 if(colN[i]==tgN && colZ[i]==tgZ)
220 if(lastCS<=0. && pEn>
lastTH)
241 colP.push_back(pMom);
260 static G4ThreadLocal std::vector <G4double> *PIN_G4MT_TLS_ = 0 ;
if (!PIN_G4MT_TLS_) PIN_G4MT_TLS_ =
new std::vector <G4double> ; std::vector <G4double> &PIN = *PIN_G4MT_TLS_;
327 G4int blast=
static_cast<int>(shift);
346 G4int blast=
static_cast<int>(shift);
401 if ( tgZ == 1 && tgN == 0 )
403 for (
G4int ip=0; ip<n_kmpel; ip++)
lastPAR[ip]=kmp_el[ip];
427 lastPAR[0]=.06*asa/(1.+a*(.01+.1/ssa));
429 lastPAR[2]=.1*a2*ssa/(1.+.0015*a2/ssa);
443 lastPAR[11]=.7/(1.+4.e-12*a16);
444 lastPAR[12]=2.5/a8/(a4+1.e-16*a32);
452 lastPAR[19]=(42.+2.7e-11*a16)/(1.+.14*a);
456 lastPAR[22]=2.4e-3*a8/(1.+2.6e-4*a7);
457 lastPAR[23]=3.5e-36*a32*a8/(1.+5.e-15*a32/
a);
459 lastPAR[24]=1.e5/(a8+2.5e12/a16);
460 lastPAR[25]=8.e7/(a12+1.e-27*a28*a28);
469 lastPAR[32]=1.5e-4/(1.+5.e-6*a12);
477 lastPAR[38]=1.e8*std::exp(.32*asa);
478 lastPAR[39]=20.*std::exp(.45*asa);
480 lastPAR[41]=2.5e5*std::exp(.085*a3);
494 lastPAR[ 9]=4.5*std::pow(a,1.15);
495 lastPAR[10]=.06*std::pow(a,.6);
496 lastPAR[11]=.6*a/(1.+2.e15/a16);
497 lastPAR[12]=.17/(a+9.e5/a3+1.5e33/a32);
498 lastPAR[13]=(.001+7.e-11*
a5)/(1.+4.4e-11*a5);
499 lastPAR[14]=(p1a10*p1a10+2.e-29)/(1.+2.e-22*a12);
501 lastPAR[15]=400./a12+2.e-22*a9;
502 lastPAR[16]=1.e-32*a12/(1.+5.e22/a14);
503 lastPAR[17]=1000./a2+9.5*sa*ssa;
504 lastPAR[18]=4.e-6*a*asa+1.e11/a16;
505 lastPAR[19]=(120./a+.002*
a2)/(1.+2.e14/a16);
513 lastPAR[25]=1.e-5*a2+2.e14/a16;
514 lastPAR[26]=1.2e-11*a2/(1.+1.5e19/a12);
515 lastPAR[27]=.016*asa/(1.+5.e16/a16);
517 lastPAR[28]=.002*a4/(1.+7.e7/std::pow(a-6.83,14));
518 lastPAR[29]=2.e6/a6+7.2/std::pow(a,.11);
519 lastPAR[30]=11.*a3/(1.+7.e23/a16/a8);
523 lastPAR[33]=3.5e-4*a2/(1.+1.e8/a8);
526 lastPAR[36]=1.e-9/a+s4a16*s4a16;
546 lastPAR[51]=1.e15+2.e27/a4/(1.+2.e-18*a16);
576 for(
G4int ip=ini; ip<=fin; ip++)
595 else G4cout<<
"*Warning*G4ChipsKaonMinusElasticXS::GetPTables: PDG="<<PDG
596 <<
", Z="<<tgZ<<
", N="<<tgN<<
", i="<<ini<<
" > fin="<<fin<<
", LP="<<LP
597 <<
" > ILP="<<ILP<<
" nothing is done!"<<
G4endl;
599 else G4cout<<
"*Warning*G4ChipsKaonMinusElasticXS::GetPTables: PDG="<<PDG
600 <<
", Z="<<tgZ<<
", N="<<tgN<<
", i="<<ini<<
">= max="<<
nPoints<<
", LP="<<LP
601 <<
" > ILP="<<ILP<<
", lPMax="<<
lPMax<<
" nothing is done!"<<
G4endl;
610 ed <<
"PDG = " << PDG <<
", Z = " << tgZ <<
", N = " << tgN
611 <<
", while it is defined only for PDG=-321 (K-) " <<
G4endl;
612 G4Exception(
"G4ChipsKaonMinusElasticXS::GetPTables()",
"HAD_CHPS_0000",
621 if(PDG==310 || PDG==130) PDG=-321;
622 if(PDG!=-321)
G4cout<<
"*Warning*G4ChipsKaonMinusElasticXS::GetET:PDG="<<PDG<<
G4endl;
631 G4double R2=(1.-std::exp(-E2*E2*E2));
643 q2=-std::log(1.-ran)/
theB1;
649 q2=-std::log(1.-ran);
651 q2=std::pow(q2,third)/
theB2;
657 q2=-std::log(1.-ran)/
theB3;
671 if(a>6.5)E3*=tm2*tm2*tm2;
686 q2=-std::log(1.-ran)/
theB1;
687 if(std::fabs(tss)>1.e-7) q2=(std::sqrt(
theB1*(
theB1+(tss+tss)*q2))-
theB1)/tss;
693 q2=-std::log(1.-ran)/
theB2;
695 if(a<6.5) q2=std::pow(q2,third);
696 else q2=std::pow(q2,fifth);
702 q2=-std::log(1.-ran)/
theB3;
704 if(a>6.5) q2=std::pow(q2,sevth);
710 q2=-std::log(1.-ran)/
theB4;
711 if(a<6.5) q2=lastTM-q2;
715 if(!(q2>=-1.||q2<=1.))
G4cout<<
"*NAN*G4QKaonMinusElasticCS::GetExchT: -t="<<q2<<
G4endl;
727 if(
lastLP<-4.3)
return 0.;
734 ed <<
"PDG = " << PDG <<
", Z = " << tgZ <<
", N = " << tgN
735 <<
", while it is defined only for PDG=-321 (K-)" <<
G4endl;
752 if(PDG!=-321)
G4cout<<
"*Warning*G4ChipsKaonMinusElasticXS::GetTV:PDG="<<PDG<<
G4endl;
755 G4cout<<
"*Warning*G4QKaonMinusElasticCS::GetTabV:(1-92)NoIsotopes for Z="<<tgZ<<
G4endl;
771 if ( tgZ == 1 && tgN == 0 )
775 theS1=(lastPAR[13]+lastPAR[14]*dl2*dl2)/(1.+lastPAR[15]/p4/p)+
776 (lastPAR[16]/p2+lastPAR[17]*p)/(p4+lastPAR[18]*sp);
777 theB1=lastPAR[19]*std::pow(p,lastPAR[20])/(1.+lastPAR[21]/p3);
778 theS2=lastPAR[22]+lastPAR[23]/(p4+lastPAR[24]*p);
779 theB2=lastPAR[25]+lastPAR[26]/(p4+lastPAR[27]/
sp);
780 theS3=lastPAR[28]+lastPAR[29]/(p4*p4+lastPAR[30]*p2+lastPAR[31]);
781 theB3=lastPAR[32]+lastPAR[33]/(p4+lastPAR[34]);
786 return lastPAR[0]/psp+(lastPAR[1]*dp*dp+lastPAR[3])/(1.-lastPAR[4]/sp+lastPAR[5]/p4)+
787 lastPAR[6]/(
sqr(p-lastPAR[7])+lastPAR[8])+lastPAR[9]/(
sqr(p-lastPAR[10])+lastPAR[11]);
814 lastPAR[40]/(1.+lastPAR[41]*std::pow(p,lastPAR[42])));
815 theB4=lastPAR[43]*pa/p2/(1.+pa*lastPAR[44]);
836 return (lastPAR[0]*dlp*dlp+lastPAR[1]+lastPAR[2]/p3)/(1.+lastPAR[3]/p2/
sp);
851 G4double mds=dmt*std::sqrt(pP2+mK2)+mK2+mt*mt;
852 return dmt*dmt*pP2/mds;
857 ed <<
"PDG = " << PDG <<
", Z = " << tgZ <<
", N = " << tgN
858 <<
", while it is defined only for p projectiles & Z_target>0" <<
G4endl;
859 G4Exception(
"G4ChipsKaonMinusElasticXS::GetQ2max()",
"HAD_CHPS_0000",
std::vector< G4double * > PAR
~G4ChipsKaonMinusElasticXS()
std::vector< G4double * > S1T
std::vector< G4double * > S2T
G4double GetPTables(G4double lpP, G4double lPm, G4int PDG, G4int tZ, G4int tN)
std::ostringstream G4ExceptionDescription
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
std::vector< G4double * > SST
std::vector< G4double * > B2T
G4ChipsKaonMinusElasticXS()
G4ParticleDefinition * GetDefinition() const
#define G4MUTEX_INITIALIZER
static G4KaonMinus * KaonMinus()
G4double GetTotalMomentum() const
std::vector< G4double * > S4T
virtual G4bool IsIsoApplicable(const G4DynamicParticle *Pt, G4int Z, G4int A, const G4Element *elm, const G4Material *mat)
G4IonTable * GetIonTable() const
G4_DECLARE_XS_FACTORY(G4ChipsKaonMinusElasticXS)
G4GLOB_DLL std::ostream G4cout
G4double CalculateCrossSection(G4bool CS, G4int F, G4int I, G4int pPDG, G4int Z, G4int N, G4double pP)
virtual G4double GetChipsCrossSection(G4double momentum, G4int Z, G4int N, G4int pdg)
G4double GetTabValues(G4double lp, G4int pPDG, G4int tgZ, G4int tgN)
static const G4double A[nN]
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
std::vector< G4double * > B3T
G4double GetPDGMass() const
virtual G4double GetIsoCrossSection(const G4DynamicParticle *, G4int tgZ, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
static G4ParticleTable * GetParticleTable()
std::vector< G4double * > CST
static const double gigaelectronvolt
static const double millibarn
G4double GetQ2max(G4int pPDG, G4int tgZ, G4int tgN, G4double pP)
G4double GetSlope(G4int tZ, G4int tN, G4int pPDG)
std::vector< G4double * > B4T
static const G4double pos
std::vector< G4double * > B1T
G4double GetExchangeT(G4int tZ, G4int tN, G4int pPDG)
std::vector< G4double * > S3T