89 vT =
new std::vector<G4double*>;
90 vL =
new std::vector<G4double*>;
91 vX =
new std::vector<std::pair<G4double,G4double>*>;
117 std::vector<G4double*>::iterator
pos;
118 for(
pos=
vT->begin(); pos<vT->end();
pos++)
123 for(
pos=
vL->begin(); pos<vL->end();
pos++)
128 std::vector<std::pair<G4double,G4double>*>::iterator pos2;
129 for(pos2=
vX->begin(); pos2<vX->end(); pos2++)
142 if(tgA<2)
return std::make_pair(QF2In,R);
143 std::pair<G4double,G4double> ElTot=
GetElTot(pIU, pPDG, tgZ, tgN);
145 if(pPDG>999 && pIU<227. && tgZ+tgN>1) R=1.;
146 else if(ElTot.second>0.)
148 R=ElTot.first/ElTot.second;
151 return std::make_pair(QF2In,R);
158 if(m_s<
toler ||
A<2)
return 1.;
159 if(m_s>
min_s)
return 0.;
170 if(nDB)
for (i=0; i<nDB; i++)
if(
A==
vARatio[i])
314 G4double ss=std::sqrt(std::sqrt(m_s));
315 G4double P=7.48e-5*s2/(1.+8.77e12/s4/s4/s2);
330 G4cout<<
"-Warning-G4QuasiElRatios::CalcElTot: p="<<p<<
" is zero or negative"<<
G4endl;
331 return std::make_pair(El,To);
338 El=1./(.00012+p2*.2);
355 El=
LE+(
pbe*lp2+6.72+32.6/p)/(1.+rp2/p);
356 To=
LE+(
pbt*lp2+38.2+52.7*rp2)/(1.+2.72*rp2*rp2);
364 El=1./(.00012+p2*(.051+.1*p2));
381 El=
LE+(
pbe*lp2+6.72+30./p)/(1.+.49*rp2/p);
382 To=
LE+(
pbt*lp2+38.2)/(1.+.54*rp2*rp2);
391 El=1.53/(lr*lr+.0676);
399 El=
pbe*ld2+2.4+7./
sp;
400 To=
pbt*ld2+22.3+12./
sp;
415 El=
LE+(
pbe*ld2+2.4+7./
sp)/(1.+.7/p4)+.6/md+.05/hd;
416 To=
LE*3+(
pbt*ld2+22.3+12./
sp)/(1.+.4/p4)+1./md+.06/hd;
426 El=13./(lr2+lr2*lr2+.0676);
434 El=
pbe*ld2+2.4+6./
sp;
435 To=
pbt*ld2+22.3+5./
sp;
449 El=
LE+(
pbe*ld2+2.4+6./
sp)/(1.+3./p4)+.7/md;
450 To=
LE+(
pbt*ld2+22.3+5./
sp)/(1.+1./p4)+.8/md;
481 El=5.2/psp+(
pbe*ld2+2.23)/(1.-.7/
sp+.075/p4)+.004/md+.15/hd;
482 To=14./psp+(
pbt*ld2+19.5)/(1.-.21/
sp+.52/p4)+.006/md+.30/hd;
492 El=.7/(lr*lr+.0676)+2./md;
513 El=
LE+(
pbe*ld2+2.23)/(1.-.7/
sp+.1/p4)+2./md;
514 To=
LE+(
pbt*ld2+19.5)/(1.+.46/
sp+1.6/p4)+2.6/md;
522 El=1./(.002+p2*(.12+p2));
530 El=(
pbe*lp2+6.72)/(1.+2./
sp);
531 To=(
pbt*lp2+38.2+900./
sp)/(1.+27./
sp);
541 El=
LE+(
pbe*lp2+6.72+99./p2)/(1.+2./
sp+2./p4);
542 To=
LE+(
pbt*lp2+38.2+900./
sp)/(1.+27./
sp+3./p4);
560 El=80./(ye+1.)+
pbe*lp2+6.72;
561 To=(80./yt+.3)/yt+
pbt*lp2+38.2;
566 G4cout<<
"*Error*G4QuasiElRatios::CalcElTot:ind="<<I<<
" is not defined (0-7)"<<
G4endl;
570 return std::make_pair(El,To);
579 if(PDG==130||PDG==310)
584 if ( (PDG == 2212 && F) || (PDG == 2112 && !F) ) ind=0;
585 else if ( (PDG == 2112 && F) || (PDG == 2212 && !F) ) ind=1;
586 else if ( (PDG == -211 && F) || (PDG == 211 && !F) ) ind=2;
587 else if ( (PDG == 211 && F) || (PDG == -211 && !F) ) ind=3;
593 else if ( PDG == -321 || PDG == -311 || (kf && !kfl) ||
594 PDG == 411 || PDG == 421 || PDG == 431 ||
595 PDG == -521 || PDG == -511 || PDG == -531 || PDG == -541 ) ind=4;
596 else if ( PDG == 321 || PDG == 311 || (kf && kfl) ||
597 PDG == -411 || PDG == -421 || PDG == -431 ||
598 PDG == 521 || PDG == 511 || PDG == 531 || PDG == 541 ) ind=5;
599 else if ( PDG > 3000 && PDG < 5333 ) ind=6;
600 else if ( PDG > -5333 && PDG < -2000 ) ind=7;
602 G4cout<<
"*Error*G4QuasiElRatios::CalcElTotXS: PDG="<<PDG
603 <<
", while it is defined only for p,n,hyperons,anti-baryons,pi,K/antiK"<<
G4endl;
623 if(PDG==130||PDG==310)
628 if ( (PDG == 2212 && F) || (PDG == 2112 && !F) ) ind=0;
629 else if ( (PDG == 2112 && F) || (PDG == 2212 && !F) ) ind=1;
630 else if ( (PDG == -211 && F) || (PDG == 211 && !F) ) ind=2;
631 else if ( (PDG == 211 && F) || (PDG == -211 && !F) ) ind=3;
637 else if ( PDG == -321 || PDG == -311 || (kf && !kfl) ||
638 PDG == 411 || PDG == 421 || PDG == 431 ||
639 PDG == -521 || PDG == -511 || PDG == -531 || PDG == -541 ) ind=4;
640 else if ( PDG == 321 || PDG == 311 || (kf && kfl) ||
641 PDG == -411 || PDG == -421 || PDG == -431 ||
642 PDG == 521 || PDG == 511 || PDG == 531 || PDG == 541 ) ind=5;
643 else if ( PDG > 3000 && PDG < 5333 ) ind=6;
644 else if ( PDG > -5333 && PDG < -2000 ) ind=7;
646 G4cout<<
"*Error*G4QuasiElRatios::FetchElTot: PDG="<<PDG
647 <<
", while it is defined only for p,n,hyperons,anti-baryons,pi,K/antiK"<<
G4endl;
655 if(nDB)
for (i=0; i<nDB; i++)
if(ind==
vItot[i])
735 G4cout<<
"-Warning-G4QuasiElRatio::GetElTot:Z="<<
Z<<
",N="<<N<<
", return zero"<<
G4endl;
736 return std::make_pair(0.,0.);
738 std::pair<G4double,G4double> hp=
FetchElTot(pGeV, hPDG,
true);
739 std::pair<G4double,G4double> hn=
FetchElTot(pGeV, hPDG,
false);
741 return std::make_pair((
Z*hp.first+N*hn.first)/
A,(
Z*hp.second+N*hn.second)/
A);
753 G4cout<<
"-Warning-G4QuasiElRatio::GetChExF:Z="<<
Z<<
",N="<<N<<
", return zero"<<
G4endl;
754 return std::make_pair(resP,resN);
759 if (hPDG==-211||hPDG==-321||hPDG==3112||hPDG==3212||hPDG==3312) pf=
Z/(
A+N);
760 else if(hPDG==211||hPDG==321||hPDG==3222||hPDG==3212||hPDG==3322) nf=N/(
A+
Z);
761 else if(hPDG==-311||hPDG==311||hPDG==130||hPDG==310)
770 mult=1./(1.+
G4Log(pGeV+pGeV))/pGeV;
775 std::pair<G4double,G4double> hp=
FetchElTot(pGeV, hPDG,
true);
776 resP=pf*(hp.second/hp.first-1.)*mult;
780 std::pair<G4double,G4double> hn=
FetchElTot(pGeV, hPDG,
false);
781 resN=nf*(hn.second/hn.first-1.)*mult;
783 return std::make_pair(resP,resN);
804 if(NPDG==2212||NPDG==90001000)
810 else if(NPDG==90001001)
816 else if(NPDG==90002001)
822 else if(NPDG==90001002)
828 else if(NPDG==90002002)
834 else if(NPDG!=2112&&NPDG!=90000001)
836 G4cout<<
"Error:G4QuasiElRatios::Scatter:NPDG="<<NPDG<<
" is not 2212 or 2112"<<
G4endl;
852 if(pPDG==2112||pPDG==-211||pPDG==-321) PDG=2112;
857 if (PDG==2212) PDG=2112;
858 else if(PDG==2112) PDG=2212;
875 if(cost>1. || cost<-1. || !(cost>-1. || cost<=1.))
877 if (cost>1.) cost=1.;
878 else if(cost<-1.) cost=-1.;
884 G4cerr<<
"G4QuasiFreeRatio::Scat:*NAN* cost="<<cost<<
",-t="<<mint<<
",tm="<<
tm<<
G4endl;
890 if(!
RelDecayIn2(tot4M, pr4M, reco4M, dir4M, cost, cost))
892 G4cerr<<
"G4QFR::Scat:t="<<tot4M<<tot4M.
m()<<
",mT="<<mT<<
",mP="<<std::sqrt(mP2)<<
G4endl;
922 if(pPDG==-211) sPDG=111;
928 else if(pPDG==-311||pPDG==311||pPDG==130||pPDG==310) sPDG=321;
929 else if(pPDG==3112) sPDG=3212;
930 else if(pPDG==3212) sPDG=3222;
931 else if(pPDG==3312) sPDG=3322;
935 if(pPDG==211) sPDG=111;
941 else if(pPDG==-311||pPDG==311||pPDG==130||pPDG==310) sPDG=-321;
942 else if(pPDG==3222) sPDG=3212;
943 else if(pPDG==3212) sPDG=3112;
944 else if(pPDG==3322) sPDG=3312;
948 G4cout<<
"Error:G4QuasiElRatios::ChExer: NPDG="<<NPDG<<
" is not 2212 or 2112"<<
G4endl;
955 G4cout<<
"Error:G4QuasiElRatios::ChExer: BAD pPDG="<<pPDG<<
", NPDG="<<NPDG<<
G4endl;
970 if(pPDG==2112||pPDG==-211||pPDG==-321) PDG=2112;
975 if (PDG==2212) PDG=2112;
976 else if(PDG==2112) PDG=2212;
993 if(cost>1. || cost<-1. || !(cost>-1. || cost<=1.))
995 if (cost>1.) cost=1.;
996 else if(cost<-1.) cost=-1.;
999 G4cerr<<
"G4QuasiFreeRatio::ChExer:*NAN* c="<<cost<<
",t="<<mint<<
",tm="<<maxt<<
G4endl;
1006 if(!
RelDecayIn2(tot4M, pr4M, reco4M, dir4M, cost, cost))
1008 G4cerr<<
"G4QFR::ChEx:t="<<tot4M<<tot4M.
m()<<
",mT="<<mT<<
",mP="<<mS<<
G4endl;
1021 if(
A<1.5)
return 0.;
1023 if (pPDG==2212) Cex=N/(
A+
Z);
1024 else if(pPDG==2112) Cex=
Z/(
A+N);
1025 else G4cout<<
"*Warning*G4CohChrgExchange::ChExElCoef: wrong PDG="<<pPDG<<
G4endl;
1032 G4double T=(6.75+.14*dl1*dl1+13./p)/(1.+.14/p4)+.6/(p4+.00013);
1033 G4double U=(6.25+8.33e-5/p4/p)*(p*
sp+.34)/p2/p;
1052 G4cerr<<
"***G4QHad::RelDecIn2: Tachionic 4-mom="<<theMomentum<<
", E-p="<<
dE-vP<<
G4endl;
1057 G4cerr<<
"G4QHadron::RelDecIn2: *Boost* E-p shift is corrected to "<<emodif<<
G4endl;
1058 theMomentum.
setE(vP+emodif+.01*accuracy);
1069 if(vdir.
mag2() > 0.)
1076 if(maxCost> 1.) maxCost= 1.;
1077 if(minCost<-1.) minCost=-1.;
1078 if(maxCost<-1.) maxCost=-1.;
1079 if(minCost> 1.) minCost= 1.;
1080 if(minCost> maxCost) minCost=maxCost;
1081 if(std::fabs(iM-fM-sM)<.00000001)
1085 f4Mom=fR*theMomentum;
1086 s4Mom=sR*theMomentum;
1089 else if (iM+.001<fM+sM || iM==0.)
1091 G4cerr<<
"***G4QH::RelDecIn2: fM="<<fM<<
"+sM="<<sM<<
">iM="<<iM<<
",d="<<iM-fM-sM<<
G4endl;
1109 if(std::fabs(ct)<1.) pss = p * std::sqrt(1.-ct*ct);
1115 G4ThreeVector pVect=(pss*std::sin(phi))*vz+(pss*std::cos(phi))*vy+p*ct*vx;
1118 f4Mom.
setE(std::sqrt(fM2+p2));
1120 s4Mom.
setE(std::sqrt(sM2+p2));
1122 if(f4Mom.
e()+.001<f4Mom.
rho())
G4cerr<<
"*G4QH::RDIn2:*Boost* f4M="<<f4Mom<<
",e-p="
1125 if(s4Mom.
e()+.001<s4Mom.
rho())
G4cerr<<
"*G4QH::RDIn2:*Boost* s4M="<<s4Mom<<
",e-p="
G4double C(G4double temp)
static const G4double pos
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
G4double G4Log(G4double x)
CLHEP::HepLorentzVector G4LorentzVector
static constexpr double megaelectronvolt
static constexpr double gigaelectronvolt
static constexpr double twopi
static constexpr double millibarn
static constexpr double MeV
G4GLOB_DLL std::ostream G4cerr
G4GLOB_DLL std::ostream G4cout
Hep3Vector orthogonal() const
Hep3Vector cross(const Hep3Vector &) const
Hep3Vector boostVector() const
HepLorentzVector & boost(double, double, double)
void setVect(const Hep3Vector &)
G4double GetExchangeT(G4int tZ, G4int tN, G4int pPDG)
static const char * Default_Name()
virtual G4double GetChipsCrossSection(G4double momentum, G4int Z, G4int N, G4int pdg)
G4double GetExchangeT(G4int tZ, G4int tN, G4int pPDG)
virtual G4double GetChipsCrossSection(G4double momentum, G4int Z, G4int N, G4int pdg)
static const char * Default_Name()
G4VCrossSectionDataSet * GetCrossSectionDataSet(const G4String &name, G4bool warning=false)
static G4CrossSectionDataSetRegistry * Instance()
static G4Deuteron * Deuteron()
static G4Neutron * Neutron()
G4double GetPDGMass() const
static G4Pow * GetInstance()
G4double powA(G4double A, G4double y) const
static G4Proton * Proton()
std::vector< G4double * > * vL
std::pair< G4double, G4double > GetChExFactor(G4double pIU, G4int pPDG, G4int Z, G4int N)
std::vector< G4double > vMtot
G4double ChExElCoef(G4double p, G4int Z, G4int N, G4int pPDG)
std::vector< G4int > vARatio
std::vector< G4double * > * vT
std::pair< G4double, G4double > * lastXtot
G4bool RelDecayIn2(G4LorentzVector &theMomentum, G4LorentzVector &f4Mom, G4LorentzVector &s4Mom, G4LorentzVector &dir, G4double maxCost=1., G4double minCost=-1.)
G4double GetQF2IN_Ratio(G4double TotCS_mb, G4int A)
std::pair< G4double, G4double > CalcElTot(G4double pGeV, G4int Index)
std::vector< G4int > vKtot
std::pair< G4double, G4double > FetchElTot(G4double pGeV, G4int PDG, G4bool F)
std::vector< std::pair< G4double, G4double > * > * vX
std::vector< G4double > vHRatio
std::vector< G4int > vKRatio
std::pair< G4LorentzVector, G4LorentzVector > ChExer(G4int NPDG, G4LorentzVector N4M, G4int pPDG, G4LorentzVector p4M)
std::vector< G4double > vMRatio
std::pair< G4double, G4double > GetElTot(G4double pIU, G4int hPDG, G4int Z, G4int N)
G4double CalcQF2IN_Ratio(G4double TCSmb, G4int A)
G4ChipsNeutronElasticXS * NCSmanager
std::pair< G4double, G4double > lastRtot
G4ChipsProtonElasticXS * PCSmanager
std::pair< G4double, G4double > GetRatios(G4double pIU, G4int prPDG, G4int tgZ, G4int tgN)
std::vector< G4int > vNRatio
std::pair< G4double, G4double > GetElTotXS(G4double Mom, G4int PDG, G4bool F)
std::pair< G4LorentzVector, G4LorentzVector > Scatter(G4int NPDG, G4LorentzVector N4M, G4int pPDG, G4LorentzVector p4M)
std::vector< G4int > vItot
static G4Triton * Triton()