NuRadioMC.EvtGen.generator module
- NuRadioMC.EvtGen.generator.load_input_hdf5(filename)[source]
reads input file into memory
- Parameters:
- filename: string
Name of the file
- Returns:
- fin: dictionary
Dictionary containing the elements in filename
- NuRadioMC.EvtGen.generator.get_tau_95_length(energies)[source]
Returns a fit to the 95% percentile of the tau track length calculated with PROPOSAL. We calculate the 95% percentile for the largest energy.
- NuRadioMC.EvtGen.generator.write_events_to_hdf5(filename, data_sets, attributes, n_events_per_file=None, start_file_id=0)[source]
writes NuRadioMC input parameters to hdf5 file
this function can automatically split the dataset up into multiple files for easy multiprocessing
- Parameters:
- filename: string
the desired output filename (if multiple files are generated, a ‘part000x’ is appended to the filename
- data_sets: dict
a dictionary with the data sets
- attributes: dict
a dictionary containing the meta attributes
- n_events_per_file: int (optional, default None)
the number of events per file
- additional_interactions: dict or None (default)
a dictionary containing potential additional interactions, such as the second tau interaction vertex.
- NuRadioMC.EvtGen.generator.primary_energy_from_deposited(Edep, ccnc, flavor, inelasticity)[source]
Calculates the primary energy of the neutrino from the deposited energy in the medium.
- Parameters:
- Edep: float
deposited energy
- ccnc: string
indicates ‘nc’, neutral current; ‘cc’, charged current
- flavor: int
neutrino flavor
- inelasticity: float
inelasticity of the interaction
- NuRadioMC.EvtGen.generator.get_energy_from_flux(Emin, Emax, n_events, flux, rnd=None)[source]
returns randomly distribution of energy according to a flux
- Parameters:
- Emin: float
minumum energy
- Emax: float
maximum energy
- n_event: int
number of events to generate
- flux: function
must return flux as function of energy in units of events per energy, time, solid angle and area
- rnd: random generator object
if None is provided, a new default random generator object is initialized
- Returns:
- energies: array of floats
array of energies
- NuRadioMC.EvtGen.generator.get_product_position_time(data_sets, product, iE)[source]
Calculates the position of a product particle given by the NuRadioProposal module that has been created by a PROPOSAL lepton propagation.
- Parameters:
- data_sets: dictionary
Dictionary with the data sets from the generating functions (generate_eventlist_cylinder and generate_surface_muons)
- product: secondary_properties class from NuRadioPROPOSAL
Contains the properties of the shower-inducing particle given by PROPOSAL
- iE: int
Number of the event in data_sets corresponding to the product particle
- Returns:
- x, y, z, time: tuple
The 3-D position of the shower-inducing product particle and the time elapsed since the first interaction to the present interaction
- NuRadioMC.EvtGen.generator.get_energies(n_events, Emin, Emax, spectrum_type, rnd=None)[source]
generates a random distribution of enrgies following a certain spectrum
- Parameters:
- n_events: int
the total number of events
- Emin: float
the minimal energy
- Emax: float
the maximum energy
- spectrum_type: string
defines the probability distribution for which the neutrino energies are generated
‘log_uniform’: uniformly distributed in the logarithm of energy
‘E-?’: E to the -? spectrum where ? can be any float
‘IceCube-nu-2017’: astrophysical neutrino flux measured with IceCube muon sample (https://doi.org/10.22323/1.301.1005)
‘GZK-1’: GZK neutrino flux model from van Vliet et al., 2019, https://arxiv.org/abs/1901.01899v1 for 10% proton fraction (see get_proton_10 in examples/Sensitivities/E2_fluxes3.py for details)
‘GZK-2’: GZK neutrino flux model from fit to TA data (https://pos.sissa.it/395/338/)
‘GZK-1+IceCube-nu-2017’: a combination of the cosmogenic (GZK-1) and astrophysical (IceCube nu 2017) flux
‘GZK-1+IceCube-nu-2022’: a combination of the cosmogenic (GZK-1) and astrophysical (IceCube nu 2022) flux
‘GZK-2+IceCube-nu-2022’: a combination of the cosmogenic (GZK-2) and astrophysical (IceCube nu 2022) flux
- rnd: random generator object
if None is provided, a new default random generator object is initialized
- NuRadioMC.EvtGen.generator.set_volume_attributes(volume, proposal, attributes)[source]
Helper function that interprets the volume settings and sets the relevant quantities as attributes dictinary. The relevant quantities to define the volume in which neutrinos are generated are (rmin, rmax, zmin, zmax) or (xmin, xmax, ymin, ymax, zmin, zmax) based on the user settings (see explanation blow).
The user can specify a “fiducial” (mandatory) and “full” (optional) volume. In the end only events are considered (= passed to radio simulation) which produce a shower (over a configurable energy threshold) within the fiducial volume. However, neutrino vertices might be generated outside this fiducial volume but the generated charge lepton (mu or tau) reaches the fiducial volume and produces a trigger (this obviously only makes sense when running proposal). Hence, a “full” volume can be defined in which neutrinos interactions are generated. If the “full” volume is not specified in the case you run proposal, the volume is determind using the energy dependent mean-free-path length of the tau.
TL;DR - The fiducial volume needs to be chosen large enough such that no events outside of it will trigger. - The full volume is optional and most of the times it is not necessary to specify it (because it is automatically determined).
- Parameters:
- volume: dict
A dictionary specifying the simulation volume. Can be either a cylinder (specified with min/max radius and z-coordinate) or a cube (specified with min/max x-, y-, z-coordinates).
Cylinder:
- fiducial_{rmin,rmax,zmin,zmax}: float
lower r / upper r / lower z / upper z coordinate of fiducial volume
- full_{rmin,rmax,zmin,zmax}: float (optional)
lower r / upper r / lower z / upper z coordinate of simulated volume
Cube:
- fiducial_{xmin,xmax,ymin,ymax,zmin,zmax}: float
lower / upper x-, y-, z-coordinate of fiducial volume
- full_{xmin,xmax,ymin,ymax,zmin,zmax}: float (optional)
lower / upper x-, y-, z-coordinate of simulated volume
Optinal you can also define the horizontal center of the volume (if you want to displace it from the origin)
- x0, y0: float
x-, y-coordinate this shift the center of the simulation volume
- proposal: bool
specifies if secondary interaction via proposal are calculated
- attributes: dictionary
dict storing hdf5 attributes
- NuRadioMC.EvtGen.generator.generate_vertex_positions(attributes, n_events, rnd=None)[source]
Helper function that generates the vertex position randomly distributed in simulation volume. The relevant quantities are also saved into the hdf5 attributes.
- Parameters:
- attributes: dicitionary
dict storing hdf5 attributes
- rnd: random generator object
if None is provided, a new default random generator object is initialized
- NuRadioMC.EvtGen.generator.intersection_box_ray(bounds, ray)[source]
this function calculates the intersection between a ray and an axis-aligned box code adapted from https://www.scratchapixel.com/lessons/3d-basic-rendering/minimal-ray-tracer-rendering-simple-shapes/ray-box-intersection
- Parameters:
- box: array with shape (2,3)
definition of box with two points
- ray: array with shape (2,3)
definiton of ray using origin and direction 3-dim vectors
- NuRadioMC.EvtGen.generator.generate_surface_muons(filename, n_events, Emin, Emax, volume, thetamin=0.0, thetamax=3.141592653589793, phimin=0.0, phimax=6.283185307179586, start_event_id=1, plus_minus='mix', n_events_per_file=None, spectrum='log_uniform', start_file_id=0, config_file='SouthPole', tables_path=None, proposal_kwargs={}, log_level=None, max_n_events_batch=100000.0, seed=None)[source]
Event generator for surface muons
Generates muons at the surface for the atmospheric muon acceptance studies. All events are saved in an hdf5 file.
- Parameters:
- filename: string
the output filename of the hdf5 file
- n_events: int
number of events to generate
- Emin: float
the minimum neutrino energy (energies are randomly chosen assuming a uniform distribution in the logarithm of the energy)
- Emax: float
the maximum neutrino energy (energies are randomly chosen assuming a uniform distribution in the logarithm of the energy)
- volume: dict
A dictionary specifying the simulation volume. Can be either a cylinder (specified with min/max radius and z-coordinate) or a cube (specified with min/max x-, y-, z-coordinates). For more information see set_volume_attributes
Cylinder:
- fiducial_{rmin,rmax,zmin,zmax}: float
lower r / upper r / lower z / upper z coordinate of fiducial volume
- full_{rmin,rmax,zmin,zmax}: float (optional)
lower r / upper r / lower z / upper z coordinate of simulated volume
Cube:
- fiducial_{xmin,xmax,ymin,ymax,zmin,zmax}: float
lower / upper x-, y-, z-coordinate of fiducial volume
- full_{xmin,xmax,ymin,ymax,zmin,zmax}: float (optional)
lower / upper x-, y-, z-coordinate of simulated volume
Optinal you can also define the horizontal center of the volume (if you want to displace it from the origin)
- x0, y0: float
x-, y-coordinate this shift the center of the simulation volume
- thetamin: float
lower zenith angle for neutrino arrival direction
- thetamax: float
upper zenith angle for neutrino arrival direction
- phimin: float
lower azimuth angle for neutrino arrival direction
- phimax: float
upper azimuth angle for neutrino arrival direction
- start_event: int
default: 1 event number of first event
- plus_minus: string
if ‘plus’: generates only positive muons if ‘minus’: generates only negative muons else generates positive and negative muons randomly
- n_events_per_file: int or None
the maximum number of events per output files. Default is None, which means that all events are saved in one file. If ‘n_events_per_file’ is smaller than ‘n_events’ the event list is split up into multiple files. This is useful to split up the computing on multiple cores.
- spectrum: string
defines the probability distribution for which the neutrino energies are generated
‘log_uniform’: uniformly distributed in the logarithm of energy
‘E-?’: E to the -? spectrum where ? can be any float
‘IceCube-nu-2017’: astrophysical neutrino flux measured with IceCube muon sample (https://doi.org/10.22323/1.301.1005)
‘GZK-1’: GZK neutrino flux model from van Vliet et al., 2019, https://arxiv.org/abs/1901.01899v1 for 10% proton fraction (see get_proton_10 in examples/Sensitivities/E2_fluxes3.py for details)
‘GZK-1+IceCube-nu-2017’: a combination of the cosmogenic (GZK-1) and astrophysical (IceCube nu 2017) flux
- start_file_id: int (default 0)
in case the data set is distributed over several files, this number specifies the id of the first file (useful if an existing data set is extended) if True, generate deposited energies instead of primary neutrino energies
- config_file: string
The user can specify the path to their own config file or choose among the available options:
‘SouthPole’, a config file for the South Pole (spherical Earth). It consists of a 2.7 km deep layer of ice, bedrock below and air above.
‘MooresBay’, a config file for Moore’s Bay (spherical Earth). It consists of a 576 m deep ice layer with a 2234 m deep water layer below, and bedrock below that.
‘InfIce’, a config file with a medium of infinite ice
‘Greenland’, a config file for Summit Station, Greenland (spherical Earth), same as SouthPole but with a 3 km deep ice layer.
- tables_path: path
path where the proposal cross section tables are stored or should be generated
- proposal_kwargs: dict
additional kwargs that are passed to the get_secondaries_array function of the NuRadioProposal class
- log_level: logging log level or None
sets the log level in the event generation. None means system default.
- max_n_events_batch: int (default 1e6)
the maximum numbe of events that get generated per batch. Relevant if a fiducial volume cut is applied)
- seed: None of int
seed of the random state
- NuRadioMC.EvtGen.generator.generate_eventlist_cylinder(filename, n_events, Emin, Emax, volume, thetamin=0.0, thetamax=3.141592653589793, phimin=0.0, phimax=6.283185307179586, start_event_id=1, flavor=[12, -12, 14, -14, 16, -16], n_events_per_file=None, spectrum='log_uniform', deposited=False, proposal=False, proposal_config='SouthPole', proposal_tables_path=None, start_file_id=0, log_level=None, proposal_kwargs={}, max_n_events_batch=100000.0, write_events=True, seed=None, interaction_type='ccnc')[source]
Event generator
Generates neutrino interactions, i.e., vertex positions, neutrino directions, neutrino flavor, charged currend/neutral current and inelastiviy distributions. All events are saved in an hdf5 file.
- Parameters:
- filename: string
the output filename of the hdf5 file
- n_events: int
number of events to generate
- Emin: float
the minimum neutrino energy
- Emax: float
the maximum neutrino energy
- volume: dict
A dictionary specifying the simulation volume. Can be either a cylinder (specified with min/max radius and z-coordinate) or a cube (specified with min/max x-, y-, z-coordinates). For more information see set_volume_attributes
Cylinder:
- fiducial_{rmin,rmax,zmin,zmax}: float
lower r / upper r / lower z / upper z coordinate of fiducial volume
- full_{rmin,rmax,zmin,zmax}: float (optional)
lower r / upper r / lower z / upper z coordinate of simulated volume
Cube:
- fiducial_{xmin,xmax,ymin,ymax,zmin,zmax}: float
lower / upper x-, y-, z-coordinate of fiducial volume
- full_{xmin,xmax,ymin,ymax,zmin,zmax}: float (optional)
lower / upper x-, y-, z-coordinate of simulated volume
Optinal you can also define the horizontal center of the volume (if you want to displace it from the origin)
- x0, y0: float
x-, y-coordinate this shift the center of the simulation volume
- thetamin: float
lower zenith angle for neutrino arrival direction (default 0deg)
- thetamax: float
upper zenith angle for neutrino arrival direction
- phimin: float
lower azimuth angle for neutrino arrival direction
- phimax: float
upper azimuth angle for neutrino arrival direction
- start_event_id: int
default: 1 event number of first event
- flavor: array of ints
default: [12, -12, 14, -14, 16, -16] specify which neutrino flavors to generate. A uniform distribution of all specified flavors is assumed.
The neutrino flavor (integer) encoded as using PDF numbering scheme, particles have positive sign, anti-particles have negative sign, relevant for us are:
12: electron neutrino
14: muon neutrino
16: tau neutrino
- n_events_per_file: int or None
the maximum number of events per output files. Default is None, which means that all events are saved in one file. If ‘n_events_per_file’ is smaller than ‘n_events’ the event list is split up into multiple files. This is useful to split up the computing on multiple cores.
- spectrum: string
defines the probability distribution for which the neutrino energies are generated
‘log_uniform’: uniformly distributed in the logarithm of energy
‘E-?’: E to the -? spectrum where ? can be any float
‘IceCube-nu-2017’: astrophysical neutrino flux measured with IceCube muon sample (https://doi.org/10.22323/1.301.1005)
‘IceCube-nu-2022’: astrophysical neutrino flux measured with IceCube muon sample (https://doi.org/10.48550/arXiv.2111.10299)
‘GZK-1’: GZK neutrino flux model from van Vliet et al., 2019, https://arxiv.org/abs/1901.01899v1 for 10% proton fraction (see get_proton_10 in examples/Sensitivities/E2_fluxes3.py for details)
‘GZK-2’: GZK neutrino flux model from fit to TA data (https://pos.sissa.it/395/338/)
‘GZK-1+IceCube-nu-2017’: a combination of the cosmogenic (GZK-1) and astrophysical (IceCube nu 2017) flux
‘GZK-1+IceCube-nu-2022’: a combination of the cosmogenic (GZK-1) and astrophysical (IceCube nu 2022) flux
‘GZK-2+IceCube-nu-2022’: a combination of the cosmogenic (GZK-2) and astrophysical (IceCube nu 2022) flux
- deposited: bool
if True, generate deposited energies instead of primary neutrino energies
- proposal: bool
if True, the tau and muon secondaries are calculated using PROPOSAL
- proposal_config: string or path
The user can specify the path to their own config file or choose among the available options:
‘SouthPole’, a config file for the South Pole (spherical Earth). It consists of a 2.7 km deep layer of ice, bedrock below and air above.
‘MooresBay’, a config file for Moore’s Bay (spherical Earth). It consists of a 576 m deep ice layer with a 2234 m deep water layer below, and bedrock below that.
‘InfIce’, a config file with a medium of infinite ice
‘Greenland’, a config file for Summit Station, Greenland (spherical Earth), same as SouthPole but with a 3 km deep ice layer.
- proposal_tables_path: path
path where the proposal cross section tables are stored or should be generated
- start_file_id: int (default 0)
in case the data set is distributed over several files, this number specifies the id of the first file (useful if an existing data set is extended)
- log_level: logging log level or None
sets the log level in the event generation. None means system default.
- proposal_kwargs: dict
additional kwargs that are passed to the get_secondaries_array function of the NuRadioProposal class
- max_n_events_batch: int (default 1e6)
the maximum numbe of events that get generated per batch. Relevant if a fiducial volume cut is applied)
- write_events: bool
if True (default) the event list is written to an output file if False the event datasets + atrributes are returned
- seed: None of int
seed of the random state
- interaction_type: string
the interaction type. default is “ccnc” which randomly choses neutral current (NC) or charged-current (CC) interactions. The use can also specify “nc” or “cc” to exclusively simulate NC or CC interactions