NuWro / nuwro

The main distribution of NuWro
GNU General Public License v3.0
13 stars 15 forks source link

Changes to Include QEL Hyperon Production and Generic Fluxes/Geometries #20

Closed cthorpe123 closed 3 years ago

cthorpe123 commented 3 years ago

QEL Hyperon Production, Hyperon FSI and Flexible Flux/Geometry

These additions affect two areas:

  1. A new quasielastic hyperon production channel, simulating the antineutrino-nucleon interaction and final state interactions for hyperons. A related change to the DIS channel preserves hyperons produced this way.
  2. Changes to flux/geometry handling to enable users to use their own unit conventions for lengths, density and window positioning.

Below is a detailed list of changes:

QEL Hyperon Production

This simulates the quasielastic production of hyperons (Lambda, Sigma0 and Sigma- baryons) from nucleons and nuclei. The model includes pseudotensor second class current and SU(3) symmetry breaking corrections.

  1. New event generator code hypevent.cc generates the primary interaction for the three sub channels and decides which to use. The cross sections themselves are calculated inside hyperon_interaction.cc. This file has several functions that calculate the hyperon production differential cross sections, the one currently in use is Singh_Model2. The corresponding form factor calculators have been added to ff.cc.

  2. nuwro.cc is modified to include this new channel.

  3. params.xml includes several new parameters: A switch for the hyperon production channel: dyn_hyp_cc. Switches for the three sub channels: hyp_lambda, hyp_sigma_zero and hyp_sigma_minus. Setting of the second class current strength: hyp_g2_Re_part, hyp_g2_Im_part. Switch for the SU(3) symmetry corrections: hyp_su3_sym_breaking. Separate axial mass parameter: hyp_axial_mass. Strengths of the Lambda-nucleus and Sigma-nucleus potentials: hyp_Lambda_Eb and hyp_Sigma_Eb.

These have all been added to data/params.txt.

  1. A new function PDG::hyperon was added to pdg.h that returns true if given the pdg code of a Lambda or Sigma.

  2. Several new functions are added to scatter.cc: scatter2_with_BE and scatter2_with_BE_SC generate kinematics for a 2->2 process where one of the final state particles is off shell. These are not currently used. rescale_momenta used to converting kinematics for a Lambda production event to those of a Sigma0 event, keeping all the angles in the CMS the same. Used for comparing their cross sections in hypevent.cc.

  3. New functions to retrieve binding energy for Bodek-Ritchie nuclear model. Not currently used.

DIS Hyperons

It appears that by default Pythia will decay hyperons. Changes made to DIS code to stabilize hyperons produced this way.

  1. dis/dis_pdg.h has pdg codes for hyperons added.

  2. dis/disevent.cc,dis/dishadr.cc and dis/reshadr.cc were changed to instruct Pythia not to decay hyperons.

  3. dis/disevent.cc must now be passed the nucleus object to calculate the hyperon potential.

Hyperon FSI

Changes include final state interactions for hyperons.

  1. Changes to kaskada7.cc :

Now adds hyperons to the particle queue alongside nucleons and pions inside prepare_particles and modifies the hyperon's kinematics to propagate it in the potential (energy is added if the potential is attractive, subtracted if the potential is repulsive).

move_particle modified to check if hyperon should be jailed at the end of every step. Similarly in leave_nucleus and finalize_interaction, check if the hyperon should be jailed when trying to exit the nucleus of after an interaction.

make_interaction checks if a secondary interaction involved a hyperon and if that hyperon has switched families (Lambda <-> Sigma) it will change the value of the potential that hyperon is propagating through accordingly. If the hyperon cannot be moved to the new potential this secondary interaction is ignored.

  1. Changes to Interaction.cc:

total_cross_sections was modified to check if it has a hyperon, There is a new function, get_hyp_xsec that retrieves the cross sections for hyperon FSI which is invoked.

Similarly in particle_scattering, if the particle given is a hyperon, calls a new function hyperon_scattering to generate the interaction.

New function get_hyp_xsec calculates some kinematic variables for the hyperon-nucleon cross section formulae, then calls functions in hyperon_cascade.cc that gets the cross sections themselves.

New function hyperon_scattering performs the hyperon-nucleon scatter itself. Calls the hyperon_state function from hyperon_cascade.cc to choose a final state by Monte Carlo, sets the process ID and then generates the kinematics using functions in scatter.cc.

New function hyperon_process_id to get hyperon process ID. enum at the top of Interaction.h modified to accommodate hyperon process codes.

  1. proctable.h was expanded to include three categories of hyperon secondary interactions, elastic, Lambda->Sigma and Sigma->Lambda. The corresponding counters are nod[15],nod[16] and nod[17].

  2. hyperon_cascade.cc contains a collection of functions for calculating the hyperon-nucleon cross sections used in the cascade, calculating the phase space factors appearing in the formulae and selecting a hyperon-nucleon final state by Monte Carlo.

Flux/Geometry Changes

Unit conventions used by ND280 were previously hard coded in. Several new parameters have been introduced to perform unit conversions so other fluxes and geometries can be used.

  1. Changes to params.xml:

The following new parameters are used:

geom_length_units, the length units used by the detector geometry. geom_density_convert takes the density used by the geometry and converts to g/cm3. The required value has been added to data/target/ND280.txt. beam_length_units, the length units used for the beam window placement (units of the xnu,ynu and znu branches of the tree containing the flux rays). beam_pot_per_file the number of protons on target for each flux file used. Currently assumes this is the same for every file. If this is not supplied in the parameters then NuWro assumes ND280 POT conventions.

These use ND280 conventions by default.

  1. Changes to beamRF.h:

The beamRF constructor now uses beam_length_units to calculate any conversion factors. It also reads in the POT per file from params.

nu_from_event now reads the z coordinate from the event and applies and required unit conversions.

  1. Changes to geomy.h:

geomy constructor sets up length and density unit conversions using the values of geom_length_units and geom_density_convert.

These conversions are used in calculating masses of detector volumes and number densities.