I3GenieTracReader¶
This module reads GENIE Generator output, samples primary neutrino direction and interaction vertex position and then writes event information to I3File. Requires GENIE’s native GHEP format root files to be converted to both GST and gRooTracker formats.
Event geometry sampling¶
It is expected that input GENIE events are produced with GENIE’s generic event
generation app gevgen
, which uses point geometry (defined in
genie::geometry::PointGeomAnalyzer
).
Caution
Using GENIE simulation with non-point geometry as genie-reader input will cause incorrect results!
The fact that GENIE event was simulated with point geometry means that it’s interaction vertex is located in the center of coordinates and it has maximum density-weighted path length of \(1 kg/m^2\) for the given composition of target materials.
genie-reader samples event geometry by moving and rotating all particles in the event record and by rescaling interaction probability for the primary neutrino. More precisely, it is done in the following way:
The code first identifies any daughter muons of the initial neutrino interaction. The energy of the most energetic muon is used to determine the generation cylinder length extension using the muon_volume_scaling configuration parameter. If no muons are found amongst the daughter particles, the length extension is 0.
\[L' = muon\_volume\_scaling * muon.energy\]Interaction vertex position is sampled uniformly inside the cylinder according to (in cylindrical coordinates):
\[ \begin{align}\begin{aligned}\rho = R * \sqrt{Uniform(0, 1)}\\\theta = Uniform(0, 2\pi)\\h = Uniform(-L/2 - L', , +L/2)\end{aligned}\end{align} \]here \(R\) is radius of the cylinder, \(L\) is cylinder length (\(R\) and \(L\) are input parameters for the module), and L’ is the volume extension from any daughter muons. \(Uniform(x, y)\) is a random number from a uniform distribution between \(x\) and \(y\).
Cylinder axis align with incoming neutrino direction.
Neutrino direction is sampled uniformly from the surface of the shere:
\[ \begin{align}\begin{aligned}\theta_{zen} = \arccos(Uniform(\cos(\theta_{zen,min}), \cos(\theta_{zen,max})))\\\phi_{az} = Uniform(\phi_{az,min}, \phi_{az,max})\end{aligned}\end{align} \]where \(\theta_{zen}\) and \(\phi_{az}\) are zenith and azimuth angles respectively; \(\theta_{zen,min}\), \(\theta_{zen,max}\), \(\phi_{az,min}\) and \(\phi_{az,max}\) are limits on zenith and azimuth, which are taken as input parameters.
After new neutrino direction is sampled, all simulated particles from GENIE event record are rotated according to it.
Neutrino interaction probability (and related event weight) is rescaled to account for larger generation volume. Interaction probability in GENIE is defined as:
\[P = \frac{N_{A}\cdot \sigma\cdot p_{L}}{A}\]where \(N_{A}\) is Avogadro number; \(\sigma\) is total cross section; \(p_{L}\) is density-weighted neutrino path length; \(A\) is target material molar mass.
Because in the point geometry from
gevgen
\(p_{L} = 1 g/cm^2\) for a given mix of target materials, interaction probability can be rescaled to match different path length and density:\[P_{new} = P_{old} * L[cm] * \rho_{ice}[g/cm^3]\]where \(L\) is total cylinder length (including the muon extension) and \(\rho_{ice}\) is ice density.
In total two parameters are rescaled:
GlobalProbabilityScale
(used inOneWeight
calculation) andInteractionProbability
.
Example usage¶
tray.AddModule(I3GenieTracReader,
GenieFileGST = 'gntp.NuMu.A.1400002.gst.root',
GenieFileGRooTracker = 'gntp.NuMu.A.1400002.gtrac.root',
GenieInfo = info,
IceDensity = 0.93*I3Units.g/I3Units.cm3,
RandomService = RandomService,
ResultDictName = "I3GenieResult",
MCTreeName = "I3MCTree",
MCWeightDictName = "I3MCWeightDict",
OutputResultDict = True,
GenVolCenter = I3Position(46.29,-34.88,-330.0))
Signature¶
- class icecube.genie_reader.I3GenieTracReader(context)¶
I3GenieTracReader module.
- Parameters:
GenieFileGST (str, default='') – Name of the GENIE gst file to read events from.
GenieFileGRooTracker (str, default='') – Name of the GENIE gRooTracker file to read events from. Assumes that it is converted from the same GHEP file, as gst input file!
GenieInfo (I3GenieInfo (struct), default=None) – I3GenieInfo object with generation info (defined in simclasses).
IceDensity (float, default=0.93*I3Units.g/I3Units.cm3) – Ice density.
RandomService (default=None) – The I3RandomService we are going to use.
ResultDictName (str, default="I3GenieResult") – Name of the output I3GenieResult frame object.
MCTreeName (str, default="I3MCTree") – Name of the output I3MCTree frame object.
MCWeightDictName (str, default=I3MCWeightDict") – Name of the output I3MCWeightDict frame object.
OutputResultDict (bool, default=True) – Choose whether or not to output I3GenieResult.
GenVolCenter (dataclasses.I3Position, default=I3Position(0.,0.,0.)*I3Units.m) – I3Position of generation volume center.
Methods:
- Configure()¶
Configures I3GenieTracReader:
First opens GenieFileGST and GenieFileGRooTracker, then sets ice density and probability scaling factor, then gets random service from context, then sets volume center, then sets frame index to zero.
- static calc_norm(px, py, pz)¶
Calculates Euclidean norm for 3D vector.
- static rotate(dir_pos, zenith, azimuth)¶
Rotates I3Direction/I3Position by provided zenith and azimuth angles.
- set_i3direction(px, py, pz, zenith, azimuth)¶
Sets I3Direction according to provided projections px, py and pz. Then rotates created I3Direction by provided zenith and azimuth angles using the I3GenieTracReader.rotate method.
- sample_primary(gst)¶
Samples primary neutrino direction and vertex position.
- get_resultdict(gst, gtrac, primary)¶
Fills the I3GenieResult object (resultdict) with information from GST and gRooTracker files.
Sampled vertex position is used and all particles directions are rotated according to sampled direction. Interaction probability and global probability scale are scaled to match new maximum path length.
- static get_mctree(resultdict, primary)¶
Fills the I3MCTree object (mctree) using information from resultdict.
- get_weightdict(resultdict, primary)¶
Fills the I3MCWeightDict object (weightdict).
- DAQ(frame)¶
If this is the first frame to process (frame index=0), writes I3GENIEInfo to I3Frame. Samples primary neutrino direction and vertex position, gets I3GenieResult, I3MCTree and I3MCWeightDict and writes them to DAQ frame.
- Finish()¶
Finish, check if number of processed frames matches to the number of events in I3GenieInfo.
Output objects¶
I3GenieTracReader produces four types of objects:
I3GenieInfo
(S-Frame)Structure defined in simcalsses:
icetray::simclasses::I3GenieInfo
. Contains file-wise information: input parameters (e.g. neutrino type, energy and directional limits, number of events in the file) andGlobalProbabilityScale
. Example:I3GenieInfo [I3GenieInfo]: [ I3GenieInfo run_id : 140000 n_events : 10 n_flux__events : 12 cylinder_height : 500 cylinder_radius : 250 min_zenith : 0 max_zenith : 3.14159 min_azimuth : 0 max_azimuth : 6.28319 min_energy : 1 max_energy : 5 power_law_index : 2 oxygen_fraction : 0.888102 ice_density : 5.8046e+30 global_probability_scale : 8.03895e-11 ]
I3MCTree
(Q-Frame)I3MCTree. Contains true primary neutrino, outgoing lepton and particles from hadronic shower (after intranuclear rescattering).
Caution
Energy of I3Particles, saved in the I3MCTree is their kinetic energy.
Example:
I3MCTree [TreeBase::Tree<I3Particle, I3ParticleID, i3hash<I3ParticleID> >]: [I3MCTree: 5 NuMu (193.923m, -185.174m, -397.057m) (39.1559deg, 289.836deg) 0ns 1.04671GeV nanm 6 MuMinus (193.923m, -185.174m, -397.057m) (40.0252deg, 29.5117deg) 0ns 0.341934GeV nanm 7 PPlus (193.923m, -185.174m, -397.057m) (66.4657deg, 263.303deg) 0ns 1.32416GeV nanm 8 PiPlus (193.923m, -185.174m, -397.057m) (49.1233deg, 48.5163deg) 0ns 0.213221GeV nanm ]
I3MCWeightDict
(Q-Frame)Dictionary (
I3MapStringDouble
) which contains essential information about event (e.g. neutrino energy, interaction type, cross section, OneWeight and other weighting related variables). Format and variable names are the same as in genie-icetray. Example:I3MCWeightDict [I3Map<string, double>]: [Crosssection => 5.1486e-12, EnergyLost => 0, GENIEWeight => 0.261366, ProbabilityScalingFactor => 88357320, GeneratorVolume => 9.81748e+07, GlobalProbabilityScale => 8.03895e-11, InjectionSurfaceR => 250, InteractionProbabilityWeight => 2.10832e-13, InteractionType => 1, LengthInVolume => 250.655, MaxAzimuth => 6.28319, MaxEnergyLog => 0.69897, MaxZenith => 3.14159, MinAzimuth => 0, MinEnergyLog => 0, MinZenith => 0, NEvents => 15, NEventsRequested => 10, OneWeight => 4.5439e-05, OneWeightPerType => 4.5439e-05, PowerLawIndex => 2, PrimaryNeutrinoEnergy => 1.04671, TargetPDGCode => 2212, TotalDetectionLength => 500, TotalInteractionProbabilityWeight => 2.10111e-11]
I3GenieResult
(Q-Frame)Structure defined in simcalsses:
icetray::simclasses::I3GenieResult
. Contains combined information from GST and gRooTracker event records.