83 using std::ostringstream;
85 using namespace genie;
116 <<
"************ Running hA2018 MODE INTRANUKE ************";
117 GHepParticle * nuclearTarget = evrec -> TargetNucleus();
118 nuclA = nuclearTarget ->
A();
122 LOG(
"HAIntranuke2018",
pINFO) <<
"Done with this event";
132 LOG(
"HAIntranuke2018",
pERROR) <<
"** Null input!";
142 bool is_handled = (is_baryon || is_pion || is_kaon || is_gamma);
144 LOG(
"HAIntranuke2018",
pERROR) <<
"** Can not handle particle: " << p->
Name();
155 LOG(
"HAIntranuke2018",
pERROR) <<
"** Couldn't select a fate";
176 <<
"Generating kinematics for " << p->
Name()
199 LOG(
"HAIntranuke2018",
pWARN) <<
"Running PreEquilibrium for kIHAFtCmp";
209 <<
"Failed attempt to generate kinematics for " 215 <<
"Failed attempt to generate kinematics for " 218 <<
" attempts. Trying a new fate...";
238 <<
"Selecting hA fate for " << p->
Name() <<
" with KE = " << ke <<
" MeV";
241 unsigned int iter = 0;
266 double frac_rescale = 1./(frac_cex + frac_inel + frac_abs + frac_piprod);
268 frac_cex *= frac_rescale;
269 frac_inel *= frac_rescale;
270 frac_abs *= frac_rescale;
271 frac_piprod *= frac_rescale;
275 double tf = frac_cex +
281 double r = tf * rnd->
RndFsi().Rndm();
282 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 283 LOG(
"HAIntranuke2018",
pDEBUG) <<
"r = " << r <<
" (max = " << tf <<
")";
286 if(r < (cf += frac_cex ))
return kIHAFtCEx;
289 if(r < (cf += frac_abs ))
return kIHAFtAbs;
293 <<
"No selection after going through all fates! " 294 <<
"Total fraction = " << tf <<
" (r = " << r <<
")";
320 double frac_rescale = 1./(frac_cex + frac_inel + frac_abs + frac_pipro);
322 frac_cex *= frac_rescale;
323 frac_inel *= frac_rescale;
324 frac_abs *= frac_rescale;
325 frac_pipro *= frac_rescale;
328 double tf = frac_cex +
335 double r = tf * rnd->
RndFsi().Rndm();
336 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 337 LOG(
"HAIntranuke2018",
pDEBUG) <<
"r = " << r <<
" (max = " << tf <<
")";
340 if(r < (cf += frac_cex ))
return kIHAFtCEx;
343 if(r < (cf += frac_abs ))
return kIHAFtAbs;
345 if(r < (cf += frac_cmp ))
return kIHAFtCmp;
348 <<
"No selection after going through all fates! " 349 <<
"Total fraction = " << tf <<
" (r = " << r <<
")";
360 double tf = frac_inel +
362 double r = tf * rnd->
RndFsi().Rndm();
363 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 364 LOG(
"HAIntranuke2018",
pDEBUG) <<
"r = " << r <<
" (max = " << tf <<
")";
368 if(r < (cf += frac_abs ))
return kIHAFtAbs;
384 const int nprob = 25;
385 double dintor = 0.0174533;
386 double denom = 47979.453;
387 double rprob[nprob] = {
388 5000., 4200., 3000., 2600., 2100., 1800., 1200., 750., 500., 230., 120.,
389 35., 9., 3., 11., 18., 29., 27., 20., 14., 10., 6., 2., 0.14, 0.19 };
391 double angles[nprob];
392 for(
int i=0;
i<nprob;
i++) angles[
i] = 2.5*
i;
395 double r = rnd->
RndFsi().Rndm();
402 for(
int i=0;
i<60;
i++) {
404 for(
int j=0;
j < nprob-1;
j++) {
408 if(binl<=theta && binh>=theta)
break;
412 double tfract = (theta-binl)/2.5;
413 double delp = rprob[itj+1] - rprob[itj];
414 xsum += (rprob[itj] + tfract*delp)/denom;
422 <<
"Generated pi+A elastic scattering angle = " << theta <<
" radians";
437 const int nprob = 20;
438 double dintor = 0.0174533;
439 double denom = 11967.0;
440 double rprob[nprob] = {
441 2400., 2350., 2200., 2000., 1728., 1261., 713., 312., 106., 35.,
442 6., 5., 10., 12., 11., 9., 6., 1., 1., 1. };
444 double angles[nprob];
445 for(
int i=0;
i<nprob;
i++) angles[
i] = 1.0*
i;
448 double r = rnd->
RndFsi().Rndm();
455 for(
int i=0;
i< nprob;
i++) {
457 for(
int j=0;
j < nprob-1;
j++) {
461 if(binl<=theta && binh>=theta)
break;
465 double tfract = (theta-binl)/2.5;
466 double delp = rprob[itj+1] - rprob[itj];
467 xsum += (rprob[itj] + tfract*delp)/denom;
475 <<
"Generated N+A elastic scattering angle = " << theta <<
" radians";
486 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 488 <<
"ElasHA() is invoked for a : " << p->
Name()
509 int pcode = p->
Pdg();
510 double Mp = p->
Mass();
518 TLorentzVector t4PpL = *p->
P4();
519 TLorentzVector t4PtL =
fRemnP4;
524 else C3CM = TMath::Cos(this->
PiBounce());
527 TLorentzVector t4P3L, t4P4L;
531 LOG(
"HAIntranuke2018",
pNOTICE) <<
"ElasHA() failed";
533 exception.
SetReason(
"TwoBodyKinematics failed in ElasHA");
544 <<
"C3cm = " << C3CM;
546 <<
"|p3| = " << t4P3L.Vect().Mag() <<
", E3 = " << t4P3L.E() <<
",Mp = " << Mp;
548 <<
"|p4| = " <<
fRemnP4.Vect().Mag() <<
", E4 = " <<
fRemnP4.E() <<
",Mt = " << Mt;
560 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 562 <<
"InelasticHA() is invoked for a : " << p->
Name()
579 int pcode = p->
Pdg();
580 int tcode, scode, s2code;
606 {
LOG(
"HAIntranuke2018",
pWARN) <<
"InelasticHA() cannot handle fate: " 608 <<
" for particle " << p->
Name();
623 LOG(
"HAIntranuke2018",
pNOTICE) <<
"InelasticHA() stops : not enough nucleons";
632 LOG(
"HAIntranuke2018",
pWARN) <<
"InelasticHA() failed : too few protons in nucleus";
643 double tM = t.
Mass();
648 target.SetHitNucPdg(tcode);
661 double pM = p->
Mass();
662 double E_p = ((*p->
P4() + *t.
P4()).
Mag2() - tM*tM - pM*pM)/(2.0*tM);
670 LOG(
"HAIntranuke2018",
pWARN) <<
"unphysical angle chosen in InelasicHA - put particle outside nucleus";
675 double KE1L = p->
KinE();
676 double KE2L = t.
KinE();
678 <<
" KE1L = " << KE1L <<
" " << KE1L <<
" KE2L = " << KE2L;
687 double E3L = cl1.
KinE();
688 double E4L = cl2.
KinE();
689 LOG (
"HAIntranuke2018",
pINFO) <<
"Successful quasielastic scattering or charge exchange";
691 <<
"C3CM = " << C3CM <<
"\n P3L, E3L = " 692 << P3L <<
" " << E3L <<
" P4L, E4L = "<< P4L <<
" " << E4L ;
694 <<
"P4L = " << P4L <<
" ;E4L= " << E4L <<
"\n probe KE = " << ev->
Probe()->
KinE() <<
"\n";
696 TParticlePDG * remn = 0;
703 <<
"NO Particle with [A = " <<
fRemnA <<
", Z = " << fRemnZ
704 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
708 MassRem = remn->Mass();
710 <<
"Particle with [A = " <<
fRemnA <<
", Z = " << fRemnZ
711 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
715 double MRemn =
TMath::Sqrt(ERemn*ERemn - PRemn*PRemn);
716 LOG(
"HAIntranuke2018",
pINFO) <<
"PRemn = " << PRemn <<
" ;ERemn= " << ERemn;
717 LOG(
"HAIntranuke2018",
pINFO) <<
"MRemn= " << MRemn <<
" ;true Mass= " << MassRem <<
" ; excitation energy= " << (MRemn-MassRem)*1000. <<
" MeV";
724 exception.
SetReason(
"TwoBodyCollison gives KE> probe KE in hA simulation");
734 exception.
SetReason(
"TwoBodyCollison failed in hA simulation");
768 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 770 <<
"Inelastic() is invoked for a : " << p->
Name()
774 bool allow_dup =
true;
785 ev,p,&s1,&s2,&s3,
fRemnA,
fRemnZ,
fRemnP4,
fDoFermi,
fFermiFac,
fFermiMomentum,
fNuclmodel);
788 LOG (
"HAIntranuke2018",
pINFO) <<
" successful pion production fate";
803 LOG(
"HAIntranuke2018",
pNOTICE) <<
"Error: could not create pion production final state";
805 exception.
SetReason(
"PionProduction kinematics failed - retry kinematics");
823 LOG(
"HAIntranuke2018",
pNOTICE) <<
"stop propagation - could not create absorption final state: too few particles";
830 LOG(
"HAIntranuke2018",
pNOTICE) <<
"stop propagation - could not create absorption final state: Pi- or K- cannot be absorbed by only neutrons";
837 LOG(
"HAIntranuke2018",
pINFO) <<
"stop propagation - could not create absorption final state: Pi+ or K+ cannot be absorbed by only protons";
851 int t1code,t2code,scode,s2code;
857 double Prob_pipd_pp=2.*ppcnt*(1.-ppcnt);
858 double Prob_pipnn_pn=.083*(1.-ppcnt)*(1.-ppcnt);
859 if (rnd->
RndFsi().Rndm()*(Prob_pipd_pp+Prob_pipnn_pn)<Prob_pipd_pp){
867 double Prob_pimd_nn=2.*ppcnt*(1.-ppcnt);
868 double Prob_pimpp_pn=.083*ppcnt*ppcnt;
869 if (rnd->
RndFsi().Rndm()*(Prob_pimd_nn+Prob_pimpp_pn)<Prob_pimd_nn){
877 double Prob_pi0d_pn=0.88*ppcnt*(1.-ppcnt);
878 double Prob_pi0pp_pp=.14*ppcnt*ppcnt;
879 double Prob_pi0nn_nn=.14*(1.-ppcnt)*(1.-ppcnt);
880 if (rnd->
RndFsi().Rndm()*(Prob_pi0d_pn+Prob_pi0pp_pp+Prob_pi0nn_nn)<Prob_pi0d_pn){
883 else if (rnd->
RndFsi().Rndm()*(Prob_pi0d_pn+Prob_pi0pp_pp+Prob_pi0nn_nn)<(Prob_pi0d_pn+Prob_pi0pp_pp)){
890 LOG(
"HAIntranuke2018",
pINFO) <<
"choose 2 body absorption, probe, fs = " << pdgc <<
" "<< scode <<
" "<<s2code;
893 double M2_1 = pLib->
Find(t1code)->Mass();
894 double M2_2 = pLib->
Find(t2code)->Mass();
896 double M3 = pLib->
Find(scode) ->Mass();
897 double M4 = pLib->
Find(s2code)->Mass();
901 TVector3 tP2_1L, tP2_2L;
906 target.SetHitNucPdg(t1code);
912 target.SetHitNucPdg(t2code);
921 tP2_1L.SetXYZ(0.0, 0.0, 0.0);
924 tP2_2L.SetXYZ(0.0, 0.0, 0.0);
927 TLorentzVector dNucl_P4=TLorentzVector(tP2_1L+tP2_2L,E2_1L+E2_2L);
929 double E2L = E2_1L + E2_2L;
936 LOG(
"HAIntranuke2018",
pWARN) <<
"Inelastic() failed: IntBounce returned bad angle - try again";
938 exception.
SetReason(
"unphysical angle for hN scattering");
943 TLorentzVector t4P1L,t4P2L,t4P3L,t4P4L;
945 t4P2L=TLorentzVector(TVector3(tP2_1L+tP2_2L),E2L);
960 TParticlePDG * remn = 0;
967 <<
"NO Particle with [A = " <<
fRemnA <<
", Z = " <<
fRemnZ 968 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
972 MassRem = remn->Mass();
974 <<
"Particle with [A = " <<
fRemnA <<
", Z = " <<
fRemnZ 975 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
979 double MRemn =
TMath::Sqrt(ERemn*ERemn - PRemn*PRemn);
980 LOG(
"HAIntranuke2018",
pINFO) <<
"PRemn = " << PRemn <<
" ;ERemn= " << ERemn;
981 LOG(
"HAIntranuke2018",
pINFO) <<
"expt MRemn= " << MRemn <<
" ;true Mass= " << MassRem <<
" ; excitation energy (>0 good)= " << (MRemn-MassRem)*1000. <<
" MeV";
993 t1->SetPdgCode(scode);
994 t1->SetMomentum(t4P3L);
1011 LOG(
"HAIntranuke2018",
pNOTICE) <<
"Inelastic in hA failed calling TwoBodyKineamtics";
1013 exception.
SetReason(
"Pion absorption kinematics through TwoBodyKinematics failed");
1037 nd0 = -135.227 * TMath::Exp(-7.124*
fRemnZ /
double(
fRemnA)) + 4.914;
1039 Sig_nd = 2.034 +
fRemnA * 0.007846;
1041 double c1 = 0.041 + ke * 0.0001525;
1042 double c2 = -0.003444 - ke * 0.00002324;
1045 double c3 = 0.064 - ke * 0.000015;
1046 gam_ns = c1 * TMath::Exp(c2*
fRemnA) +
c3;
1047 if(gam_ns<0.002) gam_ns = 0.002;
1049 LOG(
"HAIntranuke2018",
pINFO) <<
"nucleon absorption";
1050 LOG(
"HAIntranuke2018",
pINFO) <<
"--> mean diff distr = " << nd0 <<
", stand dev = " << Sig_nd;
1052 LOG(
"HAIntranuke2018",
pINFO) <<
"--> gam_ns = " << gam_ns;
1056 ns0 = .0001*(1.+ke/250.) * (
fRemnA-10)*(
fRemnA-10) + 3.5;
1057 nd0 = (1.+ke/250.) - ((
fRemnA/200.)*(1. + 2.*ke/250.));
1058 Sig_ns = (10. + 4. * ke/250.)*TMath::Power(
fRemnA/250.,0.9);
1059 Sig_nd = 4*(1 - TMath::Exp(-0.03*ke));
1060 LOG(
"HAIntranuke2018",
pINFO) <<
"pion absorption";
1061 LOG(
"HAIntranuke2018",
pINFO) <<
"--> mean diff distr = " << nd0 <<
", stand dev = " << Sig_nd;
1062 LOG(
"HAIntranuke2018",
pINFO) <<
"--> mean sum distr = " << ns0 <<
", Stand dev = " << Sig_ns;
1066 ns0 = (rnd->
RndFsi().Rndm()>0.5?3:2);
1070 LOG(
"HAIntranuke2018",
pINFO) <<
"kaon absorption - set ns, nd later";
1076 LOG(
"HAIntranuke2018",
pWARN) <<
"Inelastic() cannot handle absorption reaction for " << p->
Name();
1086 double u1 = 0, u2 = 0;
1092 LOG(
"HAIntranuke2018",
pNOTICE) <<
"Error: could not choose absorption final state";
1093 LOG(
"HAIntranuke2018",
pNOTICE) <<
"--> mean diff distr = " << nd0 <<
", stand dev = " << Sig_nd;
1094 LOG(
"HAIntranuke2018",
pNOTICE) <<
"--> mean sum distr = " << ns0 <<
", Stand dev = " << Sig_ns;
1095 LOG(
"HAIntranuke2018",
pNOTICE) <<
"--> gam_ns = " << gam_ns;
1098 exception.
SetReason(
"Absorption choice of # of p,n failed");
1107 u1 = rnd->
RndFsi().Rndm();
1108 u2 = rnd->
RndFsi().Rndm();
1109 if (u1==0) u1 = rnd->
RndFsi().Rndm();
1110 if (u2==0) u2 = rnd->
RndFsi().Rndm();
1119 ns = -TMath::Log(rnd->
RndFsi().Rndm())/gam_ns;
1123 ns = (rnd->
RndFsi().Rndm()<0.5?2:3);
1134 double max = ns0 + Sig_ns * 10;
1137 bool not_found =
true;
1145 LOG(
"HAIntranuke2018",
pNOTICE) <<
"Error: stuck in random variable loop for ns";
1146 LOG(
"HAIntranuke2018",
pNOTICE) <<
"--> mean of sum parent distr = " << ns0 <<
", Stand dev = " << Sig_ns;
1150 exception.
SetReason(
"Random number generator for choice of #p,n final state failed - unusual - redo kinematics");
1155 u1 = rnd->
RndFsi().Rndm();
1156 u2 = rnd->
RndFsi().Rndm();
1157 if (u1==0) u1 = rnd->
RndFsi().Rndm();
1158 if (u2==0) u2 = rnd->
RndFsi().Rndm();
1161 ns = ns0 + Sig_ns *
x1;
1162 if ( ns>max || ns<0 ) {iter2++;
continue;}
1163 else if ( rnd->
RndFsi().Rndm() > (ns/
max) ) {iter2++;
continue;}
1171 double nd = nd0 + Sig_nd *
x2;
1176 np =
int((ns+nd)/2.+.5);
1177 nn =
int((ns-nd)/2.+.5);
1179 LOG(
"HAIntranuke2018",
pINFO) <<
"ns = "<<ns<<
", nd = "<<nd<<
", np = "<<np<<
", nn = "<<nn;
1185 if (np < 0 || nn < 0 ) {iter++;
continue;}
1186 else if (np + nn < 2. ) {iter++;
continue;}
1189 - ((pdgc==
kPdgPiM || pdgc==
kPdgKM)?1:0)) {iter++;
continue;}
1194 LOG(
"HAIntranuke2018",
pINFO) <<
"success, iter = " << iter <<
" np, nn = " << np <<
" " << nn;
1197 double frac = 85./double(np+nn);
1205 if (rnd->
RndFsi().Rndm()<np/(double)(np+nn)) np--;
1209 LOG(
"HAIntranuke2018",
pNOTICE) <<
"Final state chosen; # protons : " 1210 << np <<
", # neutrons : " << nn;
1238 PDGCodeList* listar[5] = {&list0, &list1, &list2, &list3, &list4};
1243 double probM = pLib->
Find(pdgc) ->Mass();
1245 TVector3 pP3 = p->
P4()->Vect() * (1./5.);
1246 double probKE = p->
P4()->E() -probM;
1247 double clusKE = probKE * (1./5.);
1248 TLorentzVector clusP4(pP3,clusKE);
1249 LOG(
"HAIntranuke2018",
pINFO) <<
"probM = " << probM <<
" ;clusKE= " << clusKE;
1250 TLorentzVector X4(*p->
X4());
1265 for (
int i=0;
i<(np+nn);
i++)
1275 for (
int i=0;
i<5;
i++)
1277 LOG(
"HAIntranuke2018",
pINFO) <<
"List" <<
i <<
" size: " << listar[
i]->size();
1278 if (listar[
i]->
size() <2)
1281 exception.
SetReason(
"too few particles for Phase Space decay - try again");
1321 if(success1 && success2 && success3 && success4 && success5)
1323 LOG(
"HAIntranuke2018",
pINFO)<<
"Successful many-body absorption - n>=18";
1325 TParticlePDG * remn = 0;
1326 double MassRem = 0.;
1332 <<
"NO Particle with [A = " <<
fRemnA <<
", Z = " << fRemnZ
1333 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
1337 MassRem = remn->Mass();
1339 <<
"Particle with [A = " <<
fRemnA <<
", Z = " << fRemnZ
1340 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
1344 double MRemn =
TMath::Sqrt(ERemn*ERemn - PRemn*PRemn);
1345 LOG(
"HAIntranuke2018",
pINFO) <<
"PRemn = " << PRemn <<
" ;ERemn= " << ERemn;
1346 LOG(
"HAIntranuke2018",
pINFO) <<
"MRemn= " << MRemn <<
" ;true Mass= " << MassRem <<
" ; excitation energy (>0 good)= " << (MRemn-MassRem)*1000. <<
" MeV";
1352 LOG(
"HAIntranuke2018",
pWARN) <<
"PhaseSpace decay fails for HadrCluster- recovery likely incorrect - rethrow event";
1391 double probM = pLib->
Find(pdgc) ->Mass();
1392 double probBE = (np+nn)*.005;
1393 TVector3 pP3 = p->
P4()->Vect();
1394 double probKE = p->
P4()->E() - (probM - probBE);
1395 double clusKE = probKE;
1396 TLorentzVector clusP4(pP3,clusKE);
1397 LOG(
"HAIntranuke2018",
pINFO) <<
"probM = " << probM <<
" ;clusKE= " << clusKE;
1398 TLorentzVector X4(*p->
X4());
1407 for (
int i=0;
i<np;
i++)
1413 for (
int i=0;
i<nn;
i++)
1449 <<
"Remnant nucleus (A,Z) = (" <<
fRemnA <<
", " <<
fRemnZ <<
")";
1450 LOG(
"HAIntranuke2018",
pINFO) <<
" list size: " << np+nn;
1454 exception.
SetReason(
"too few particles for Phase Space decay - try again");
1463 LOG (
"HAIntranuke2018",
pINFO) <<
"Successful many-body absorption, n<=18";
1464 LOG(
"HAIntranuke2018",
pDEBUG) <<
"Nucleus : (A,Z) = ("<<
fRemnA<<
','<<fRemnZ<<
')';
1465 TParticlePDG * remn = 0;
1466 double MassRem = 0.;
1472 <<
"NO Particle with [A = " <<
fRemnA <<
", Z = " << fRemnZ
1473 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
1477 MassRem = remn->Mass();
1479 <<
"Particle with [A = " <<
fRemnA <<
", Z = " << fRemnZ
1480 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
1484 double MRemn =
TMath::Sqrt(ERemn*ERemn - PRemn*PRemn);
1485 LOG(
"HAIntranuke2018",
pINFO) <<
"PRemn = " << PRemn <<
" ;ERemn= " << ERemn;
1486 LOG(
"HAIntranuke2018",
pINFO) <<
"expt MRemn= " << MRemn <<
" ;true Mass= " << MassRem <<
" ; excitation energy (>0 good)= " << (MRemn-MassRem)*1000. <<
" MeV";
1495 if ( pdgc==
kPdgPiM ) fRemnZ++;
1498 exception.
SetReason(
"Phase space generation of absorption final state failed");
1565 LOG(
"HAIntranuke2018",
pINFO) <<
"R0 = " <<
fR0 <<
" fermi";
1575 LOG(
"HAIntranuke2018",
pINFO) <<
"DoFermi? = " << ((
fDoFermi)?(
true):(
false));
static string AsString(INukeMode_t mode)
An exception thrown by SimulateHadronState for kinematics problems. TwoBodyCollision/Kinematics used ...
void SetFirstMother(int m)
void PreEquilibrium(GHepRecord *ev, GHepParticle *p, int &RemnA, int &RemnZ, TLorentzVector &RemnP4, bool DoFermi, double FermiFac, const NuclearModelI *Nuclmodel, double NucRmvE, EINukeMode mode=kIMdHN)
virtual GHepParticle * Particle(int position) const
TRandom3 & RndFsi(void) const
rnd number generator used by intranuclear cascade monte carlos
double fNucleonFracAbsScale
double fFermiMomentum
whether or not particle collision is pauli blocked
THE MAIN GENIE PROJECT NAMESPACE
void SimulateHadronicFinalState(GHepRecord *ev, GHepParticle *p) const
int HandleCompoundNucleus(GHepRecord *ev, GHepParticle *p, int mom) const
double fPionFracPiProdScale
static RandomGen * Instance()
Access instance.
const TLorentzVector * P4(void) const
double fDelRPion
factor by which Pion Compton wavelength gets multiplied to become nuclear size enhancement ...
bool PhaseSpaceDecay(GHepRecord *ev, GHepParticle *p, const PDGCodeList &pdgv, TLorentzVector &RemnP4, double NucRmvE, EINukeMode mode=kIMdHA)
general phase space decay method
double fEPreEq
threshold for pre-equilibrium reaction
double fNR
param multiplying the nuclear radius, determining how far to track hadrons beyond the "nuclear bounda...
Float_t x1[n_points_granero]
enum genie::EGHepStatus GHepStatus_t
double fFermiFac
testing parameter to modify fermi momentum
unsigned int fNumIterations
int fRemnA
remnant nucleus A
double fNucleonFracCExScale
bool PionProduction(GHepRecord *ev, GHepParticle *p, GHepParticle *s1, GHepParticle *s2, GHepParticle *s3, int &RemnA, int &RemnZ, TLorentzVector &RemnP4, bool DoFermi, double FermiFac, double FermiMomentum, const NuclearModelI *Nuclmodel)
double PiBounce(void) const
void Inelastic(GHepRecord *ev, GHepParticle *p, INukeFateHA_t fate) const
Pure abstract base class. Defines the NuclearModelI interface to be implemented by any physics model ...
::xsd::cxx::tree::exception< char > exception
double Mass(Resonance_t res)
resonance mass (GeV)
bool fXsecNNCorr
use nuclear medium correction for NN cross section
void SetMomentum(const TLorentzVector &p4)
double Mass(void) const
Mass that corresponds to the PDG code.
void ProcessEventRecord(GHepRecord *event_rec) const
void InelasticHA(GHepRecord *ev, GHepParticle *p, INukeFateHA_t fate) const
double fNucAbsFac
absorption xsec correction factor (hN Mode)
A singleton holding random number generator classes. All random number generation in GENIE should tak...
double Pz(void) const
Get Pz.
const NuclearModelI * fNuclmodel
nuclear model used to generate fermi momentum
const TVector3 & Momentum3(void) const
double Px(void) const
Get Px.
virtual GHepParticle * Probe(void) const
int LastMother(void) const
double FracADep(int hpdgc, INukeFateHA_t fate, double ke, int targA) const
int FirstMother(void) const
double fPionFracInelScale
string Name(void) const
Name that corresponds to the PDG code.
const int kPdgCompNuclCluster
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
enum genie::EINukeFateHN_t INukeFateHN_t
double fPionMFPScale
tweaking factors for tuning
double fNucleonFracPiProdScale
void SetReason(string reason)
A Neutrino Interaction Target. Is a transparent encapsulation of quite different physical systems suc...
bool TwoBodyKinematics(double M3, double M4, TLorentzVector tP1L, TLorentzVector tP2L, TLorentzVector &tP3L, TLorentzVector &tP4L, double C3CM, TLorentzVector &RemnP4, double bindE=0)
bool TwoBodyCollision(GHepRecord *ev, int pcode, int tcode, int scode, int s2code, double C3CM, GHepParticle *p, GHepParticle *t, int &RemnA, int &RemnZ, TLorentzVector &RemnP4, EINukeMode mode=kIMdHA)
Intranuke utility functions.
double FracAIndep(int hpdgc, INukeFateHA_t fate, double ke) const
virtual GHepParticle * TargetNucleus(void) const
bool fDoFermi
whether or not to do fermi mom.
double frac(double x)
Fractional part.
double fNucCEXFac
charge exchange xsec correction factor (hN Mode)
static const double ns
Module that plots metrics from reconstructed cosmic ray data.
Misc GENIE control constants.
double fR0
effective nuclear size param
void ElasHA(GHepRecord *ev, GHepParticle *p, INukeFateHA_t fate) const
bool fDoCompoundNucleus
whether or not to do compound nucleus considerations
virtual void ProcessEventRecord(GHepRecord *event_rec) const
double KinE(bool mass_from_pdg=false) const
Get kinetic energy.
static PDGLibrary * Instance(void)
double PnBounce(void) const
static INukeHadroData2018 * Instance(void)
void SetLastMother(int m)
Singleton class to load & serve a TDatabasePDG.
Var Sqrt(const Var &v)
Use to take sqrt of a var.
const TLorentzVector * X4(void) const
void SetStatus(GHepStatus_t s)
static string AsString(INukeFateHN_t fate)
bool fAltOset
NuWro's table-based implementation (not recommended)
double fNucleonFracInelScale
INukeHadroData2018 * fHadroData2018
a collection of h+N,h+A data & calculations
double fDelRNucleon
factor by which Nucleon Compton wavelength gets multiplied to become nuclear size enhancement ...
bool IsNeutronOrProton(int pdgc)
TParticlePDG * Find(int pdgc)
int IonPdgCode(int A, int Z)
double fHadStep
step size for intranuclear hadron transport
void SimulateHadronicFinalStateKinematics(GHepRecord *ev, GHepParticle *p) const
virtual bool GenerateNucleon(const Target &) const =0
virtual void AddParticle(const GHepParticle &p)
int nuclA
value of A for the target nucleus in hA mode
bool GetParamDef(const RgKey &name, T &p, const T &def) const
bool GetParam(const RgKey &name, T &p, bool is_top_call=true) const
double fNucRmvE
binding energy to subtract from cascade nucleons
static const unsigned int kRjMaxIterations
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
GENIE's GHEP MC event record.
bool fUseOset
Oset model for low energy pion in hN.
Most commonly used PDG codes. A set of utility functions to handle PDG codes is provided in PDGUtils...
STDHEP-like event record entry that can fit a particle or a nucleus.
INukeFateHA_t HadronFateHA(const GHepParticle *p) const
enum genie::EINukeFateHA_t INukeFateHA_t
Root of GENIE utility namespaces.
double IntBounce(const GHepParticle *p, int target, int s1, INukeFateHN_t fate)
void push_back(int pdg_code)
int fRemnZ
remnant nucleus Z
TLorentzVector fRemnP4
P4 of remnant system.
double Py(void) const
Get Py.
const Algorithm * SubAlg(const RgKey ®istry_key) const