32 #include <TLorentzVector.h> 34 #include <TParticlePDG.h> 61 using namespace genie;
100 InitialState * init_state = interaction -> InitStatePtr();
101 Target * tgt = init_state -> TgtPtr();
110 if(p4->Px()>0 || p4->Py()>0 || p4->Pz()>0)
return;
120 double rad = nucleon->
X4()->Vect().Mag();
127 <<
"Generated nucleon momentum: (" 128 << p3.Px() <<
", " << p3.Py() <<
", " << p3.Pz() <<
"), " 129 <<
"|p| = " << p3.Mag();
131 <<
"Generated nucleon removal energy: w = " <<
w;
133 double pF2 = p3.Mag2();
143 int eject_nucleon_pdg = 0;
147 EN = nucleon->
Mass() - w -
148 pF2 / (2 * (nucleus->
Mass() - nucleon->
Mass()));
151 TParticlePDG * other_nucleon;
159 EN = deuteron->Mass() - 2 * w -
160 TMath::Sqrt(pF2 + other_nucleon->Mass() * other_nucleon->Mass());
162 deuteron->PdgCode() != nucleus->
Pdg()) {
166 eject_nucleon_pdg = other_nucleon->PdgCode();
172 int nucleon_pdgc = nucleon->
Pdg();
174 int Z = (is_p) ? nucleus->
Z()-1 : nucleus->
Z();
175 int A = nucleus->
A() - 1;
177 TParticlePDG * fnucleus = 0;
182 <<
"No particle with [A = " << A <<
", Z = " << Z
183 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
188 double Mf = fnucleus ->
Mass();
189 double Mi = nucleus ->
Mass();
193 double MN = nucleon->
Mass();
194 double MN2 = TMath::Power(MN,2);
198 const int nucleus_pdgc = nucleus->
Pdg();
199 const int nucleon_pdgc = nucleon->
Pdg();
208 double numNuc = (is_p) ? (
double)tgt->
Z():(double)tgt->
N();
209 double radius = nucleon->
X4()->Vect().Mag();
211 kF= TMath::Power(3*
kPi2*numNuc*
221 double prob = rnd->
RndGen().Rndm();
229 p4->SetPx( p3.Px() );
230 p4->SetPy( p3.Py() );
231 p4->SetPz( p3.Pz() );
243 <<
"Event below threshold after generating Fermi momentum";
248 <<
"Ev (@ nucleon rest frame) = " << Ev <<
", Ethr = " << Ethr;
252 exception.
SetReason(
"E < Ethr after generating nucleon Fermi momentum");
256 if(eject_nucleon_pdg != 0) {
262 const int eject_pdg_code)
const 264 LOG(
"FermiMover",
pINFO) <<
"Adding a recoil nucleon " 265 "due to short range corelation";
273 double vx = nucleon->
Vx();
274 double vy = nucleon->
Vy();
275 double vz = nucleon->
Vz();
276 double px = -1.* nucleon->
Px();
277 double py = -1.* nucleon->
Py();
278 double pz = -1.* nucleon->
Pz();
283 eject_pdg_code, status, imom, -1, -1, -1, px, py, pz, E, vx, vy, vz, 0);
290 LOG(
"FermiMover",
pINFO) <<
"Adding final state nucleus";
298 int A = nucleus->
A();
299 int Z = nucleus->
Z();
304 for(
int id = fd;
id <= ld;
id++) {
309 int pdgc = particle->
Pdg();
314 if (is_p || is_n) A--;
316 Px += particle->
Px();
317 Py += particle->
Py();
318 Pz += particle->
Pz();
323 TParticlePDG * remn = 0;
328 <<
"No particle with [A = " << A <<
", Z = " << Z
329 <<
", pdgc = " << ipdgc <<
"] in PDGLibrary!";
333 double Mi = nucleus->
Mass();
341 <<
"Adding nucleus [A = " << A <<
", Z = " << Z
342 <<
", pdgc = " << ipdgc <<
"]";
346 ipdgc,
kIStStableFinalState, imom,-1,-1,-1, Px,Py,Pz,E, 0,0,0,0);
363 RgKey nuclkey =
"NuclearModel";
bool fSRCRecoilNucleon
simulate recoil nucleon due to short range corellation?
const KPhaseSpace & PhaseSpace(void) const
virtual GHepParticle * Particle(int position) const
void Configure(const Registry &config)
double RemovalEnergy(void) const
THE MAIN GENIE PROJECT NAMESPACE
double E(void) const
Get energy.
virtual Interaction * Summary(void) const
static RandomGen * Instance()
Access instance.
Defines the EventRecordVisitorI interface. Concrete implementations of this interface use the 'Visito...
double Density(double r, int A, double ring=0.)
int FirstDaughter(void) const
static const double fermi
double Threshold(void) const
Energy threshold.
enum genie::EGHepStatus GHepStatus_t
bool IsNucleus(void) const
void ProcessEventRecord(GHepRecord *event_rec) const
static FermiMomentumTablePool * Instance(void)
void KickHitNucleon(GHepRecord *evrec) const
give hit nucleon a momentum
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)
A table of Fermi momentum constants.
void SetMomentum(const TLorentzVector &p4)
double Mass(void) const
Mass that corresponds to the PDG code.
A singleton holding random number generator classes. All random number generation in GENIE should tak...
double Pz(void) const
Get Pz.
const TVector3 & Momentum3(void) const
FermiMoverInteractionType_t GetFermiMoverInteractionType(void) const
double Px(void) const
Get Px.
virtual NuclearModel_t ModelType(const Target &) const =0
Summary information for an interaction.
int LastDaughter(void) const
An exception thrown by EventRecordVisitorI when the normal processing sequence has to be disrupted (f...
Singleton class to load & serve tables of Fermi momentum constants.
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
const FermiMomentumTable * GetTable(string name)
A Neutrino Interaction Target. Is a transparent encapsulation of quite different physical systems suc...
virtual void Configure(const Registry &config)
const NuclearModelI * fNuclModel
nuclear model
virtual GHepParticle * TargetNucleus(void) const
static constexpr Double_t rad
static const double kLightSpeed
enum genie::EFermiMoverInteractionType FermiMoverInteractionType_t
void SetRemovalEnergy(double Erm)
void SwitchOnFastForward(void)
TLorentzVector * HitNucP4Ptr(void) const
static PDGLibrary * Instance(void)
void SetReason(string reason)
virtual TBits * EventFlags(void) const
Var Sqrt(const Var &v)
Use to take sqrt of a var.
double Vz(void) const
Get production z.
const TLorentzVector * X4(void) const
A registry. Provides the container for algorithm configuration parameters.
virtual GHepParticle * HitNucleon(void) const
TRandom3 & RndGen(void) const
rnd number generator for generic usage
TParticlePDG * Find(int pdgc)
int IonPdgCode(int A, int Z)
bool IsAboveThreshold(void) const
Checks whether the interaction is above the energy threshold.
assert(nhit_max >=nhit_nbins)
virtual bool GenerateNucleon(const Target &) const =0
virtual void AddParticle(const GHepParticle &p)
const InitialState & InitState(void) const
void AddTargetNucleusRemnant(GHepRecord *evrec) const
^ emit a 2nd nucleon due to short range corellations
double Vy(void) const
Get production y.
double FindClosestKF(int target_pdgc, int nucleon_pdgc) const
bool GetParamDef(const RgKey &name, T &p, const T &def) const
const Target & Tgt(void) const
bool fKeepNuclOnMassShell
keep hit bound nucleon on the mass shell?
static const double kPlankConstant
void Emit2ndNucleonFromSRC(GHepRecord *evrec, const int eject_nucleon_pdg) const
double ProbeE(RefFrame_t rf) const
GENIE's GHEP MC event record.
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.
virtual int TargetNucleusPosition(void) const
double Vx(void) const
Get production x.
Initial State information.
double Py(void) const
Get Py.
const Algorithm * SubAlg(const RgKey ®istry_key) const