marley::ElectronReaction Class Reference

Inheritance diagram for marley::ElectronReaction:

Public Member Functions
	ElectronReaction (int pdg_a, int target_atom_pdg, const std::string &source_file)
virtual marley::TargetAtom	atomic_target () const override final
	Returns the target atom involved in this reaction.
virtual std::shared_ptr< HepMC3::GenEvent >	create_event (int particle_id_a, double KEa, marley::Generator &gen) const override
	Create an event object for this reaction.
virtual double	diff_xs (int pdg_a, double KEa, double cos_theta_c_cm) const
double	g1 () const
double	g2 () const
virtual double	threshold_kinetic_energy () const override
	Get the minimum lab-frame kinetic energy (MeV) of the projectile that allows this reaction to proceed via a transition to the residue's ground state.
virtual double	total_xs (int pdg_a, double KEa) const override
	Compute the reaction's total cross section (MeV -2)
Public Member Functions inherited from marley::Reaction
	Reaction (const std::string &source_file)
	Construct a Reaction with the resolved path of the data file.
const std::string &	get_description () const
	Get a string that contains the formula for this reaction.
int	pdg_a () const
	Get the projectile PDG code.
int	pdg_b () const
	Get the target PDG code.
ProcessType	process_type () const
	Get the process type for this reaction.
const std::string &	source_file () const
	Get the resolved path of the reaction data file used to construct this Reaction.

Additional Inherited Members
Public Types inherited from marley::Reaction
enum	DataFormat { DiscreteStrengths = 0 , MultipoleResponses = 1 }
	Enumerated type describing the file format for reaction data. More...
enum	ProcessType {   Unknown = -1 , NeutrinoCC_Discrete = 0 , AntiNeutrinoCC_Discrete = 1 , NC_Discrete = 2 ,   NuElectronElastic = 3 , NeutrinoCC_Continuum = 4 , AntiNeutrinoCC_Continuum = 5 , NC_Continuum = 6 ,   StandaloneDecay = 7 }
	Enumerated type describing the kind of scattering process represented by a Reaction. More...
Static Public Member Functions inherited from marley::Reaction
static int	get_ejectile_pdg (int pdg_a, ProcessType proc_type)
static void	get_residue_pdg_and_charge (ProcessType proc_type, int pdg_b, int &pdg_d, int &q_d)
static std::vector< std::unique_ptr< Reaction > >	load_from_file (const std::string &filename, StructureDatabase &db, CoulombCorrector::CoulombMode coulomb_mode, const JSON &ff_config)
static std::string	proc_type_to_string (const ProcessType &pt)
Protected Member Functions inherited from marley::Reaction
virtual std::shared_ptr< HepMC3::GenEvent >	make_event_object (double KEa, const std::shared_ptr< HepMC3::GenParticle > &ejectile, const std::shared_ptr< HepMC3::GenParticle > &residue) const
	Helper function that makes an event object.
virtual std::shared_ptr< HepMC3::GenEvent >	make_event_object (double KEa, double pc_cm, double cos_theta_c_cm, double phi_c_cm, double Ec_cm, double Ed_cm, int residue_status) const
	Helper function that makes an event object.
void	two_two_scatter (double KEa, double &s, double &Ec_cm, double &pc_cm, double &Ed_cm) const
	Helper function that handles CM frame kinematics for the reaction.
Static Protected Member Functions inherited from marley::Reaction
static const std::vector< int > &	get_projectiles (ProcessType proc_type)
Protected Attributes inherited from marley::Reaction
std::string	description_
	String that contains a formula describing the reaction.
double	ma_
	Projectile mass (MeV)
double	mb_
	Target mass (MeV)
double	mc_
	Ejectile mass (MeV)
double	md_
	Residue mass (MeV)
int	pdg_a_
	PDG code for the projectile.
int	pdg_b_
	PDG code for the target.
int	pdg_c_
	PDG code for the ejectile.
int	pdg_d_
	PDG code for the residue.
ProcessType	process_type_
	Type of scattering process (CC, NC) represented by this reaction.
const std::string	source_file_
	Resolved path of the reaction data file.

Detailed Description

Definition at line 35 of file ElectronReaction.hh.

Constructor & Destructor Documentation

ElectronReaction()

marley::ElectronReaction::ElectronReaction	(	int	pdg_a,
		int	target_atom_pdg,
		const std::string &	source_file )

Definition at line 31 of file ElectronReaction.cc.

  : Reaction( source_file ), atom_( target_atom_pdg )
{
  process_type_ = ProcessType::NuElectronElastic;
 
  pdg_a_ = pdg_a;
  pdg_b_ = marley_utils::ELECTRON;
  pdg_c_ = get_ejectile_pdg(pdg_a_, process_type_);
  pdg_d_ = marley_utils::ELECTRON;
 
  // Set the description string based on the particle PDG codes
  description_ = marley_utils::get_particle_symbol( pdg_a_ )
    + " + " + marley_utils::get_particle_symbol( pdg_b_ );
  description_ += " --> " + marley_utils::get_particle_symbol( pdg_c_ );
  description_ += " + " + marley_utils::get_particle_symbol( pdg_d_ );
 
  const auto& mt = marley::MassTable::Instance();
  ma_ = mt.get_particle_mass( pdg_a_ );
  mb_ = mt.get_particle_mass( pdg_b_ );
  mc_ = mt.get_particle_mass( pdg_c_ );
  md_ = mt.get_particle_mass( pdg_d_ );
 
  this->set_coupling_constants();
 
  // Compute the kinetic energy of the projectile needed
  // for this reaction to proceed at threshold
  KEa_threshold_ = ( std::pow(mc_ + md_, 2)
    - std::pow(ma_ + mb_, 2) ) / ( 2.*mb_ );
 
  MARLEY_LOG( DEBUG, "physics.reaction" ) << "ElectronReaction: "
    << description_ << ", g1 = " << g1_ << ", g2 = " << g2_
    << ", threshold KE = " << KEa_threshold_ << " MeV";
}

Member Function Documentation

atomic_target()

virtual marley::TargetAtom marley::ElectronReaction::atomic_target ( ) const

inlinefinaloverridevirtual

Returns the target atom involved in this reaction.

For nuclear reactions, this is identical to the pdg_b_ member variable. For electron reactions, it is distinct (since particle b is the initial struck electron).

Implements marley::Reaction.

Definition at line 42 of file ElectronReaction.hh.

        { return atom_; }

create_event()

std::shared_ptr< HepMC3::GenEvent > marley::ElectronReaction::create_event ( int pdg_a, double KEa, marley::Generator & gen ) const

overridevirtual

Create an event object for this reaction.

Parameters

pdg_a	PDG code for the incident projectile
KEa	Lab-frame kinetic energy of the projectile
gen	Reference to the Generator to use for random sampling

Note

Functions that override create_event() should throw an Error if pdg_a != pdg_a_.

Implements marley::Reaction.

Definition at line 166 of file ElectronReaction.cc.

{
  // If the projectile's PDG code doesn't match that stored in this object,
  // complain and refuse to create an event.
  if ( pdg_a != pdg_a_ ) throw marley::Error("Could"
    " not create this event. The requested projectile PDG code "
    + std::to_string(pdg_a) + " does not match the value "
    + std::to_string(pdg_a_) + " stored in the Reaction object.");
 
  // Also complain if we're below threshold
  if ( KEa < KEa_threshold_ ) throw marley::Error("Could not create"
    " this event. The kinetic energy (" + std::to_string(KEa) + " MeV)"
    " of the projectile is below the reaction threshold of "
    + std::to_string(KEa_threshold_) + " MeV.");
 
  double s, Ec_cm, pc_cm, Ed_cm;
  two_two_scatter(KEa, s, Ec_cm, pc_cm, Ed_cm);
 
  // Compute the maximum differential cross section to use for rejection
  // sampling. To do this, we analytically solve for the value of
  // cos_theta_c_cm (labeled cth below) for which the derivative of the
  // differential cross section vanishes. This is an extremum of the function
  // and might correspond to the maximum. If cth is within the allowed angular
  // range, then it is considered alongside the two endpoints (COS_MIN and
  // COS_MAX), and the largest of the three values (or just the endpoint values
  // if cth is outside the allowed range) is chosen as the maximum.
  double me2_over_s = md_*md_ / s;
  double B = marley_utils::ONE_HALF * std::pow(g2_*(1. - me2_over_s), 2);
  double A = g1_*g2_*me2_over_s + g2_*g2_*(1. - me2_over_s) - B;
  double cth = -A / B;
 
  // Set the differential cross section to a huge negative value
  // at cth. This value will be compared to those at the angular
  // endpoints if cth does not lie in the allowed range. Otherwise,
  // the correct value will replace this one.
  double dxs_at_cth = std::numeric_limits<double>::lowest();
 
  double max = 0.;
  if ( cth >= COS_MIN && cth <= COS_MAX ) {
    dxs_at_cth = this->diff_xs(pdg_a, KEa, cth);
  }
 
  double dxs_at_min = this->diff_xs(pdg_a, KEa, COS_MIN);
  double dxs_at_max = this->diff_xs(pdg_a, KEa, COS_MAX);
 
  // Find the maximum value of the differential cross section
  max = std::max( { dxs_at_min, dxs_at_max, dxs_at_cth } );
 
  // Sample a CM frame scattering cosine for the ejectile.
  double cos_theta_c_cm = gen.rejection_sample(
    [this, pdg_a, KEa](double ctheta) -> double
    { return this->diff_xs(pdg_a, KEa, ctheta); }, COS_MIN, COS_MAX, max);
 
  // Sample a CM frame azimuthal scattering angle (phi) uniformly on [0, 2*pi).
  // We can do this because the differential cross section is independent of
  // the azimuthal angle.
  double phi_c_cm = gen.uniform_random_double(0., marley_utils::two_pi, false);
 
  MARLEY_LOG( DEBUG, "physics.reaction" ) << "ElectronReaction::create_event:"
    " KEa = " << KEa << " MeV, sampled cos_theta_cm = " << cos_theta_c_cm
    << ", phi_cm = " << phi_c_cm;
 
  // Create and return the completed event object
  return make_event_object( KEa, pc_cm, cos_theta_c_cm, phi_c_cm, Ec_cm, Ed_cm,
    marley_hepmc3::NUHEPMC_FINAL_STATE_STATUS );
}

References diff_xs(), marley::Reaction::make_event_object(), marley::Reaction::md_, marley::Reaction::pdg_a(), marley::Reaction::pdg_a_, marley::Generator::rejection_sample(), marley::Reaction::two_two_scatter(), and marley::Generator::uniform_random_double().

diff_xs()

double marley::ElectronReaction::diff_xs ( int pdg_a, double KEa, double cos_theta_c_cm ) const

virtual

Todo

Take into account effects of electron binding energy

Definition at line 127 of file ElectronReaction.cc.

{
  // If we're asked for the wrong projectile, then just return zero
  if ( pdg_a != pdg_a_ ) return 0.;
 
  // If we're below threshold, then just return zero
  if ( KEa < KEa_threshold_ ) return 0.;
 
  // If the kinetic energy is not positive, then just return zero
  if ( KEa <= 0. ) return 0.;
 
  // Mandelstam s (square of the total center of momentum frame energy)
  double s = std::pow(ma_ + mb_, 2) + 2.*mb_*KEa;
 
  // CM frame ejectile total energy
  double Ec_cm = (s + mc_*mc_ - md_*md_) / ( 2. * std::sqrt(s) );
 
  // Helper variables
  double me2_over_s = md_*md_ / s;
  double overall_factor = (2. / marley_utils::pi)
    * std::pow(marley_utils::GF * Ec_cm, 2);
  double terms = std::pow(g1_, 2) + g1_*g2_*me2_over_s*(cos_theta_c_cm - 1.)
    + std::pow(g2_ * (1. + marley_utils::ONE_HALF*(1. - me2_over_s)
    * (cos_theta_c_cm - 1.)), 2);
 
  // Compute and return the cross section
  double diff_xsec = overall_factor * terms;
 
  // Multiply the single-electron cross section by the number of electrons
  // present in this atom
  diff_xsec *= atom_.Z();
 
  return diff_xsec;
}

References marley::Reaction::ma_, marley::Reaction::mb_, marley::Reaction::mc_, marley::Reaction::md_, marley::Reaction::pdg_a(), and marley::Reaction::pdg_a_.

Referenced by create_event().

g1()

double marley::ElectronReaction::g1 ( ) const

inline

Definition at line 60 of file ElectronReaction.hh.

{ return g1_; }

g2()

double marley::ElectronReaction::g2 ( ) const

inline

Definition at line 61 of file ElectronReaction.hh.

{ return g2_; }

threshold_kinetic_energy()

virtual double marley::ElectronReaction::threshold_kinetic_energy ( ) const

inlineoverridevirtual

Get the minimum lab-frame kinetic energy (MeV) of the projectile that allows this reaction to proceed via a transition to the residue's ground state.

Implements marley::Reaction.

Definition at line 57 of file ElectronReaction.hh.

        { return KEa_threshold_; }

total_xs()

double marley::ElectronReaction::total_xs ( int pdg_a, double KEa ) const

overridevirtual

Compute the reaction's total cross section (MeV ^-2)

Parameters

pdg_a	Projectile's PDG code
KEa	Lab-frame kinetic energy of the incident projectile

Returns

Reaction total cross section (MeV ^-2)

Note

Functions that override total_xs() should always return zero if pdg_a != pdg_a_.

Todo

Take into account effects of electron binding energy

Implements marley::Reaction.

Definition at line 93 of file ElectronReaction.cc.

                                                                   {
 
  // If the cross section was requested for a different projectile,
  // then just return zero.
  if ( pdg_a != pdg_a_ ) return 0.;
 
  // If we're below threshold, then just return zero
  if ( KEa < KEa_threshold_ ) return 0.;
 
  // If the kinetic energy is not positive, then just return zero
  if ( KEa <= 0. ) return 0.;
 
  // Mandelstam s (square of the total center of momentum frame energy)
  double s = std::pow(ma_ + mb_, 2) + 2.*mb_*KEa;
 
  // CM frame ejectile total energy
  double Ec_cm = (s + mc_*mc_ - md_*md_) / ( 2. * std::sqrt(s) );
 
  // Helper variables
  double me2_over_s = md_*md_ / s;
  double g2_squared_over_three = g2_*g2_ / 3.;
 
  // Total cross section in natural units (MeV^(-2))
  double xs = (4. / marley_utils::pi) * std::pow(marley_utils::GF * Ec_cm, 2)
    * (std::pow(g1_, 2) + (g2_squared_over_three - g1_*g2_)*me2_over_s
    + g2_squared_over_three*(1. + std::pow(me2_over_s, 2)));
 
  // Multiply the single-electron cross section by the number of electrons
  // present in this atom
  xs *= atom_.Z();
  return xs;
}

References marley::Reaction::ma_, marley::Reaction::mb_, marley::Reaction::mc_, marley::Reaction::md_, marley::Reaction::pdg_a(), and marley::Reaction::pdg_a_.

---

The documentation for this class was generated from the following files:

• /home/sjg/Desktop/marley/include/marley/ElectronReaction.hh
• /home/sjg/Desktop/marley/src/ElectronReaction.cc