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:

1. 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/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). \(Uniform(x, y)\) is a random number from a uniform distribution between \(x\) and \(y\).

Cylinder axis align with incoming neutrino direction.

  1. 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.

  1. Neutrino interaction probability (and related event weight) is rescaled to account for larger generation volume. Interaction probablility 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 cylinder length and \(\rho_{ice}\) is ice density.

In total two parameters are rescaled: GlobalProbabilityScale (used in OneWeight calculation) and InteractionProbability.

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()

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:

  1. 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) and GlobalProbabilityScale. 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
    ]
    
  2. I3MCTree (Q-Frame)

    I3MCTree. Contains true primary neutrino, outgoing lepton and pariticles 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
    ]
    
  3. 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,
    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 => 10,
    OneWeight => 4.5439e-05,
    PowerLawIndex => 2,
    PrimaryNeutrinoEnergy => 1.04671,
    TargetPDGCode => 2212,
    TotalDetectionLength => 500,
    TotalInteractionProbabilityWeight => 2.10111e-11]
    
  4. I3GenieResult (Q-Frame)

    Structure defined in simcalsses: icetray::simclasses::I3GenieResult. Contains combined information from GST and gRooTracker event records.