Dk2nu GENIE FluxDriver


This class presents itself as a concrete instance of the genie::flux::GFluxI interface. The purpose of this case is to read and interpret flux (really, hadron/muon decay) files generated by the various beam simulations. These simulations could come from gnumi (the initial GEANT3 based simulation), g4numi (a Geant4-based simulation) and flugg (a combination of a modern fluka and geant4), g4lbne or a reinterpretation of the FNAL Booster beam simulation.

The most comprehensive document on the basics of these beam simulations ntuples can be found in MINOS-DocDB-6316: The NuMI Beam Simulation with Flugg. It is important to note that while entries in the ntuple represent decays of hadrons or muons, they don't reflect every decay. Generally there is an importance weight attached to each entry that must be accounted for. Also while each ntuple entry has information about the neutrino from a randomly oriented decay and possibly decays that have the neutrino ray go through two fixed points in space, these are not particularly useful for any individual experiment. The random decays generally miss any detector and the two fixed points (one "near", one "far") fail to illuminate the entirety of the detector (if they are even relevant for that detector, i.e. NOvA vs. MINOS FarDet).

The GDk2NuFlux code uses information found in the ntuple to generate neutrino rays over a flux window by evaluating new weights and energies for a ntuple entry based on the position of the random point on the flux window in the beam coordinate system. This approach was first described in MINOS-DocDB-109. The flux window should normally be chosen to shadow the entirety of the detector accounting for the orientation of the beam with respect to the detector and any divergence in the beam.

Using the code

The basics of using the code is to create and configure an object:

  genie::flux::GDk2NuFlux* gdk2nu = new genie::flux::GDk2NuFlux();
  GFluxI* fdriver = dynamic_cast<GFluxI*>(gnumi);

and then make repeated calls to generate neutrino rays:

The two strings passed to LoadBeamSimData() are a pattern that specifies flux files to use and a parameter set name in a XML file.

specifying flux file names

The input files are expected to be ROOT ntuple files of the dk2nu format. Wildcards (normal shell globs, not regexp style) can be used.

XML parameter sets

Experiment specific configurations are encapsulated in <param_set> sections of an XML file; multiple <param_set name="myname"> instances can exist in a single file as long as each has a unique "myname". By default GDk2NuFlux will look for the file GDk2Nu.xml, but this can be overridden by use of the environment variable GDK2NUFLUXXML. The GDk2NuFlux.xml file must be found in a location defined by the environment variable GXMLPATH; include "." if you want the current location searched before the default.

Configuring a XML param_set

An example parameter set might look like:

  <param_set name="MINOS-FarDet">
    <units> m </units>
    <beamdir type="series">
      <rotation axis="x" units="rad"> +0.057184957 </rotation>
    <beampos> (0, 0, 0 ) = ( 0, 0, 735340.0 ) </beampos>
      <point coord="det"> -7.000, -9.474, -13.4796 </point>
      <point coord="det">  7.000, -9.474, -13.4796 </point>
      <point coord="det"> -7.000,  7.453, -14.4487 </point>
    <enumax> 120. 1.05 1.05 2500000 </enumax>
    <reuse> 10 </reuse>
    <upstreamz> -400.0 </upstreamz>

The string "MINOS-FarDet" identifies the configuration; this is the string that would get passed as the second argument to LoadBeamSimData().

The top level entities in a <param_set> are described below. If an entity is specified more than once in a parameter set, the last version trumps any prior settings for that parameter set. One will find in the default file the <param_set name="MINOS-NearDet"> specifies the <beamdir> and <beampos> tags multiple times: this is for pedagogical purposes. All instances give (to within precision) the same result but are there to demonstrate the various alternatives. Other parameter sets do not need to duplicate this; pick one specification and stick with it.

The <beamdir> and <beampos> tags are used together to define the coordinate transformation between the user (i.e. detector) system and the beam system. The <beamdir> tag must come before the <beampos> tag (n.b. to self: code should rigorously enforce this). These two tags and the <window> are the most critical for an individual detector setup.


This controls how noisy GDk2NuFlux is when interpreting the <param_set> and is mostly for debugging purposes.


The <units> tag identifies the user coordinate system length unit (e.g. "cm", "m" or "meter", etc). This default is meter (as is the GENIE standard units); the beam flux files are assumed to be in cm. This should match the units used by the detector geometry. This feature has not been extensively tested and it is probably best to stick with meter if at all possible.


The <beamdir> defines the rotation matrix between the beam and user coordinates. There are three ways of specifying this; use the one that is most convenient for the particular setup in question. If in doubt about the result use the PrintConfig() method to output the final result and verify that it is as you expect.

The first approach is to specify the transformation as a series of rotations around various axes. Getting this right is easiest if the transformation is a single rotation about one of the primary axes; going beyond this often leads to confusion about order and the correct sign. An example might look like:

    <beamdir type="series"> 
       <rotation axis="x" units="rad"> -0.0582977560 </rotation>
       <!-- additional rotations might follow -->

Rotation units for this tag can be either deg or rad and the rotation specifies an axis about which to rotate. A note about the rotation sign: the MINOS NearDet sees a downward going beam (mixing y and z) which is specified by a negative value as shown here to give the beam a negative y component in the user frame.

The second approach is to specify the beam coordinate's axes in the user system as a set of (θ,φ) pairs, ala a GEANT3 approach:

    <beamdir type="thetaphi3" units="deg">  
      ( 90,        0 )  ( 86.65978462612, 90 ) ( 3.34021537388, 270 ) 

Punctuation (i.e. parentheses and commas) aren't relevant and are only given for clarity; tags of this type need six values.

The third and most unambiguous specification is to give the complete 3x3 matrix that would transform a direction vector in beam coordinates into one in user coordinates, ala a Geant4 approach:

    <beamdir type="newxyz"> 
     [ 1  0               0              ]
     [ 0  0.998301167046 -0.058264739543 ]
     [ 0  0.058264739543  0.998301167046 ]

A column vector of directions in beam coordinates, when multiplied by this rotation matrix from the left should give the direction in the user system, i.e. dir user = Rot 3x3 * dir beam . Again punctuation and layout is irrelevant and shown for only for clarity; tags of this type need nine values.


The <beampos> specification is the second piece needed to complete the coordinate transformations. There are two variants. If specified with three values, these are the position of the beam origin in user coordinates:

    <!-- position of beam origin in user coords -->
    <beampos> 1.4828, 60.629225, -1034.72755 </beampos>

Positions are transformed between coordinate systems by: xyz user = ( scale Rot 3x3 * xyz beam ) + xyz beampos where scale is the ratio of user units (e.g. m) to beam units (cm).

Alternatively it is sometimes easier to know a single point in space in both systems (6 values):

    <!-- position given as ( user coords ) = ( beam coords )                   -->
    <!-- both should be in user coord units, though                            -->
    <!-- convenient for MINOS-NearDet because it is how the survey is reported -->
    <!-- XML config must have already set the rotation matrix                  -->
    <beampos> ( 1.4828, 0.2385, 0.0 ) = ( 0, 0, 1036.48837 ) </beampos>

As with <beamdir> punctuation is irrelevant but useful for clarity.


The flux window on which neutrino rays will originate is given by specifying three points in user coordinate space.

      <point coord="det"> -4.358, -2.903, -0.1835 </point>
      <point coord="det">  4.796, -2.903, -0.1835 </point>
      <point coord="det"> -4.358, 10.377,  0.5923 </point>

The order in which they are specified is important; the 2nd and 3rd value are used to define vectors relative to the first, so give the common point first (unless you do actually intend to sweep out a rhombus) ala:

flux-window specification

In principle is should not matter whether the window is square to the beam direction or the detector as long as it is large enough that rays passing through it completely shadow the volume of interest (including any beam divergence). (n.b. to self: still need to empirically verify that window orientation doesn't matter)


The rejection method used by the de-weighting (importance and x-y over the window) needs an upper estimate of the neutrino energy. By default it will attempt to guess a conservative value, but CPU cycles might be saved by helping this process along. (n.b. to self: explain other values on this line)


Most entries are rejected (due to importance and x-y weights). The evaluation process can be sped up some if ntuple entries are reused by choosing a new position on the flux window a number of times. This increases the very low probability that the same decay entry might actually be accepted twice in a row. This is very rare and due to the choice of a new position on the flux window results in a different ray direction and energy, but is also undesirable. This can be tuned based on statistics given after a run by the PrintConfig() method. Typically a value of 10 is a quite safe choice.


By default neutrino rays originate on the specified flux window and GENIE in generating events only considers material in the geometry downstream of this starting position. Sometimes one would like to leave a flux window in an understood position and orientation and simply start neutrino rays elsewhere. If abs(upstreamz) < 10^30^ then each ray will be adjusted along its trajectory to start at that z position.


A convenient way to build off other configurations is to use the <using_param_set> element:

  <param_set name="MINOS-FarRock">
    <using_param_set> MINOS-FarDet </using_param_set>
      <point coord="det"> -15.0000, -17.4612, -13.0224 </point>
      <point coord="det">  15.0000, -17.4612, -13.0224 </point>
      <point coord="det"> -15.0000,  15.4402, -14.9059 </point>

This set takes as its base the "MINOS-FarDet" configuration and modifies the window size.


Sometimes you just want to put a "floor" on the estimated max weight (perhaps the distribution has a long tail that the estimator doesn't "find")

    <!-- force maxwgt to have an initial miminum of -->
    <!-- <minmaxwgt> 4.2e-07 </minmaxwgt> -->
    <minmaxwgt> 2.5e-07 </minmaxwgt>


What to do if the estimated max weight is exceeded

    <!-- what to do if exceeded [ "bump" ] "frozen" "abort" -->
    <maxwgtfail> bump </maxwgtfail>

Generating single events after a "bump" is probably fine. But if generating (using the art nutools GENIEHelper) "pileup" (multiple events in a "spill") probably it incorrect to allow it to "bump". Better if that occurs regularly, to set a reasonable minimum (see minmaxweight above) and use "frozen".