Introduction

GENIE [Andreopoulos:2009rq] is a neutrino event generator created to be the “universal event generator” requested during the NuINT conference series. It was designed to provide an accessible, extensible framework with many convenience tools leveraging existing HEP software like Pythia [Sjostrand:2006za], ROOT [Brun:1997pa], and LHAPDF [Whalley:2005nh].

How to Use this Guide

This guide is meant to bring a GENIE model developer up to speed on contributing to GENIE quickly. The most common ways to contribute are:

  • Modify the parameters of an existing model without changing the code.

  • Modify the source code of an existing model.

  • Debug an error in an existing model.

  • Author an entirely new model.

Developers may safely skip [Building-GENIE] if they have already built the code, but may find it useful if they did not build all of the third-party packages from source or lack debugging symbols in some of the packages.

Please read [Running-GENIE-Basics] first, and then do the exercises in [CCQE-Walkthrough]. Then read [Configuration-XML], although it is not quite complete. The other sections are still very rough and should not be taken seriously quite yet.

We strongly recommend skimming the NIM paper [Andreopoulos:2009rq] and reading section five closely.

The complete source code for this document and some supporting files for the exercises within (cross section spline files, the build_pythia6.sh script, some run scripts, gdb and valgrind templates for the .gdbinit and .valgrindrc files, etc.) may be checked out from Fermilab Redmine at any time with:

git clone http://cdcvs.fnal.gov/projects/developer_manual

Building GENIE

Building from Scratch on Linux

Complete instructions for installing GENIE may be found in the Physics and User Manual (PUM) published on the GENIE web page [geniemanual]. Users who have either already installed GENIE or read the PUM may safely skip this section.

This document logs the experience of working with multiple installations in parallel. We will provide examples of installing GENIE and all third-party software by hand. Developers should use centrally installed software when appropriate, but it is important to be able to run code with debugging symbols present and central installations often have compiler optimizations turned on.

We will eventually build a directory tree that looks something like this:

<some path>/GENIE                           # GENIE Work Directory
           /GENIE/support                   # for 3rd party code
                         /gsl               # install
                         /gsl-<version>     # source
                         /lhapdf            # install / data
                         /lhapdf-<version>  # source
                         /log4cpp
                         /pythia6
                         /root
           /GENIE/R-2_8_0                   # frozen release
           /GENIE/Dev-R-2_8_0               # development

Your directory structure may look different if you install some software via a package manager (for example, you don’t really need to build your own copy of log4cpp). In the examples below, we will sometimes repalce “some path” with /minerva/app/users/perdue.

First, we will check out a couple of copies of the GENIE source code and build all the third-party software (Pythia, ROOT, log4cpp, LHAPDF). It is useful to have a development version of the GENIE code as well as a “frozen release” version for comparison and to help answer questions from users.

$ #-- Get the GENIE code.
$ cd <Your Development Area>
$ mkdir GENIE  # for GENIE and local copies of 3rd party software
$ cd GENIE
$ svn list http://genie.hepforge.org/svn/branches
$ # svn co http://svn.hepforge.org/genie/branches/<genie_tag> <local_dir>
$ svn co http://genie.hepforge.org/svn/branches/R-2_8_0/ R-2_8_0/
$ # Choose a development branch that matches what you are working on.
$ # If there is no appropriate branch, don't worry about checking one out.
$ # Instead, feel free to check out another copy of R-2_8_0 in a ``dev'' area.
$ export GDVTAG="_devel_branch_not_for_users_v280_06jun13_hugh_nd9832nxll3a1/"
$ svn co http://genie.hepforge.org/svn/branches/${GDVTAG} Dev-R-2_8_0
$ unset GDVTAG

With the GENIE raw source in hand, we’ll checkout and build all the required third-part libraries, starting with Pythia (and working in the same directory as above when we checked out GENEI, in our <development area>):

$ mkdir support && cd support  # for 3rd-party
$ mkdir pythia6 && cd pythia6  # get Pythia first (for ROOT)
$ # Get a copy of R. Hatcher's build_pythia6.sh script.
$ . build_pythia6.sh

You may find a copy of build_pythia6.sh here (and via Google):

Next we’ll get GSL, but first check to see if it is already installed by executing locate libgsl. We’ll include an example installing from source below, but GSL is a library you will essentially never need to step through with a debugger, so using a package manager and installing optimized libraries should be fine.

$ cd ..    # go back to top-level support
$ tar -xvzf gsl-1.16.tgz
$ mkdir gsl && cd gsl-1.16

We should author a configure script (use your support area):

do_configure.sh
#!/bin/bash
INSTDIR=/minerva/app/users/perdue/GENIE/support/gsl
./configure --prefix=\$INSTDIR
unset INSTDIR

In this example, we use four cores in the build:

$ ./do_configure
$ time nice make -j4
$ make check
$ make install

After GSL, we will build ROOT next (any recent version 5 tag should work, version 6 has not yet been tested):

$ cd ..    # go back to top-level support
$ git clone http://root.cern.ch/git/root.git
$ cd root/
$ git tag -l         # choose your preferred tag
$ git checkout -b v5-34-08 v5-34-08

We should author a configure script, do_configure.sh:

#!/bin/sh
export GDVS=/minerva/app/users/perdue/GENIE/support
./configure linuxx8664gcc \
 --build=debug \
 --enable-pythia6 \
 --with-pythia6-libdir=$GDVS/pythia6/v6_424/lib \
 --enable-gsl-shared \
 --enable-mathmore \
 --with-gsl-incdir=$GDVS/gsl/include \
 --with-gsl-libdir=$GDVS/gsl/lib
unset GDVS

Note, you may choose to drop the “debug” flag if you plan on running optimized code for performance. In the example below, we do not specify for the make process to use additional cores (sometimes ROOT has touble if you add cores to the build):

$ make
$ cd ..    # go back to top-level support
$ # Get log4cpp source.
$ cd support/log4cpp
$ ./autogen.sh
$ gmake
$ gmake install

For the last third-party library, we will build LHAPDF. We will use LHAPDF 5.9 at this stage because version 6 and later requires the Boost libraries. A good project would be to test LHAPDF 6 with the system version of Boost and, if necessary, install a newer version of Boost.

Caution
Don’t lose track of PDF download below - it is not clear we have chosen the right PDF sets! Do not rely on this guide without verifying the correct choice of PDFs!
$ cd ..    # go back to top-level support
$ # Get source for LHAPDF (e.g., 5.9.0)
$ # LHAPDF 6+ is C++ and requires Boost.
$ cd lhapdf-5.9.0/
$ mkdir ../lhapdf  # Install into this directory.
$ ./configure --prefix=<Your Dev>/GENIE/support/lhapdf
$ gmake
$ gmake install
$ # Now we have to actually download a few PDFs.
$ # Do we choose the right ones? Who knows? Not much guidance...
$ cd bin  # bin dir in lhapdf-5.9.0
$ ./lhapdf-getdata GRV98lo.LHgrid --dest=<Your Dev>/GENIE/support/lhapdf
$ ./lhapdf-getdata GRV98nlo.LHgrid --dest=<Your Dev>/GENIE/support/lhapdf

Next, we will build GENIE. We will keep a pair of scripts to configure our environments, first for the frozen release:

environment_genie.sh
#!/bin/bash
echo "Setting GENIE environment variables..."
export GENIEBASE=/minerva/app/users/perdue/GENIE
export GENIE=\$GENIEBASE/R-2_8_0
export PYTHIA6=\$GENIEBASE/support/pythia6/v6_424/lib
export ROOTSYS=\$GENIEBASE/support/root
export LOG4CPP_INC=\$GENIEBASE/support/log4cpp/include
export LOG4CPP_LIB=\$GENIEBASE/support/log4cpp/lib
export LHAPATH=\$GENIEBASE/support/lhapdf
export LHAPDF_INC=\$GENIEBASE/support/lhapdf/include
export LHAPDF_LIB=\$GENIEBASE/support/lhapdf/lib
export XSECSPLINEDIR=\$GENIEBASE/data
export LD_LIBRARY_PATH=\$LHAPDF_LIB:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=/usr/lib64:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=\$LOG4CPP_LIB:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=\$PYTHIA6:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=\$ROOTSYS/lib:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=\$GENIE/lib:\$LD_LIBRARY_PATH
export PATH=\$GENIE/bin:\$ROOTSYS/bin:\$PATH
unset GENIEBASE

And for the development installation:

environment_devgenie.sh
#!/bin/bash
echo "Setting GENIE environment variables..."
export GENIEBASE=/minerva/app/users/perdue/GENIE
export GENIE=\$GENIEBASE/Dev-R-2_8_0
export PYTHIA6=\$GENIEBASE/support/pythia6/v6_424/lib
export ROOTSYS=\$GENIEBASE/support/root
export LOG4CPP_INC=\$GENIEBASE/support/log4cpp/include
export LOG4CPP_LIB=\$GENIEBASE/support/log4cpp/lib
export LHAPATH=\$GENIEBASE/support/lhapdf
export LHAPDF_INC=\$GENIEBASE/support/lhapdf/include
export LHAPDF_LIB=\$GENIEBASE/support/lhapdf/lib
export XSECSPLINEDIR=\$GENIEBASE/data
export LD_LIBRARY_PATH=\$LHAPDF_LIB:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=/usr/lib64:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=\$LOG4CPP_LIB:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=\$PYTHIA6:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=\$ROOTSYS/lib:\$LD_LIBRARY_PATH
export LD_LIBRARY_PATH=\$GENIE/lib:\$LD_LIBRARY_PATH
export PATH=\$GENIE/bin:\$ROOTSYS/bin:\$PATH
unset GENIEBASE

Note that \$GENIEBASE must be correctly re-defined in the above scripts for your work area.

First, we’ll build the frozen release. Go to the source directory we checked out from the SVN repository and add a configuration script:

do_configure.sh
#!/bin/bash
./configure \
  --enable-doxygen-doc \
  --enable-debug \
  --enable-test \
  --enable-numi \
  --enable-gsl \
  --with-optimiz-level=O0 \
  --with-doxygen-path=/usr/bin/doxygen \
  --with-log4cpp-inc=\$LOG4CPP_INC \
  --with-log4cpp-lib=\$LOG4CPP_LIB \
  --with-libxml2-inc=your path to the xml lib include dir
  --with-libxml2-lib=your path to the xml lib lib dir

Note that we use some of the environment variables from our environment script, so be sure to source that script first. Run the configure script (or type the options by hand). Then execute gmake. (If you choose to enable a prefix or if you would like to make a system-wide installation, you may then gmake install, but this is not necessary.)

Building from Scratch on OSX

The process is very similar to that on Linux. First, check out a copy of the GENIE source code in the same way as described above. On newer versions of OSX, Apple has stopped supporting the Gnu toolchain - this means that g++ is actually clang, the LLVM compiler. This is a problem for GENIE.

You can install the real GNU gcc fairly easily with either Fink, MacPorts, or Homebrew. I recommend Homebrew: http://brew.sh. First, get Xcode from the Mac App Store (it is free). Then, launch Xcode and open the Preferences pane. Click on “Downloads” and install the Command Line Tools. Next, install Homebrew according to the instructions on their webpage. Finally, install GCC with brew install apple-gcc42.

Then check out and build Pythia6. Note, to build Pythia6, you need a Fortran compiler and modern Apple computers may not have one. To install gfortran with homebrew, simply run brew install gfortran.

Next, install GSL (first, check to see if it is installed with mdfind libgsl). It can be installed from source, but it is fine to use a package manager like Homebrew: brew install gsl.

Then, install ROOT as above, with an appropriately modified setup script, e.g.:

#!/bin/sh
./configure macosx64 \
  --build=debug \
  --enable-pythia6 \
  --with-pythia6-libdir=<some path>/GENIE/support/pythia6/v6_424/lib \
  --enable-gsl-shared \
  --enable-mathmore \
  --with-gsl-incdir=/usr/local/include \
  --with-gsl-libdir=/usr/local/lib

Note by default ROOT will build with clang instead of the g++-4.2 we plan on building GENIE with. This may or may not be a problem (not clear yet). lldb appears to be able to step through code built with both. After ROOT, install log4cpp. You may build from source as above or use a package manager, e.g., brew install log4cpp.

For the last third-party library, install LHAPDF. Sadly, there is no Homebrew formula for it.

With everything installed, we are ready to install GENIE. Create environment and setup scripts patterned after those above. If your version of OSX did not support GCC, then we need to edit the GENIE Makefile to use the version we installed with Homebrew. If you installed apple-gcc42, then edit \$GENIE/src/make/Make.include so that it sets its compiler and linker like so:

CXX = g++-4.2
LD  = g++-4.2

Working with the Repository

GENIE is managed with a Subversion (SVN) repository on HEPForge. The repository is browsable on the web: https://genie.hepforge.org/trac/browser

$ svn list http://genie.hepforge.org/svn/branches
$ export GDVTAG="_devel_branch_not_for_users_v280_06jun13_hugh_nd9832nxll3a1/"
$ svn co http://genie.hepforge.org/svn/branches/$GDVTAG Dev-R-2_8_0
$ svn info
Path: .
URL: http://genie.hepforge.org/svn/branches/_devel_...
Repository Root: http://genie.hepforge.org/svn
Repository UUID: cc9776de-3512-45ca-aafc-e2d9ed43c22c
Revision: 3966
Node Kind: directory
Schedule: normal
Last Changed Author: sboyd11
Last Changed Rev: 3966
Last Changed Date: 2013-09-24 06:19:43 -0500 (Tue, 24 Sep 2013)

We can examine the remote repository:

$ svn info http://genie.hepforge.org/svn/branches/${GDVTAG}
Path: _devel_...
URL: http://genie.hepforge.org/svn/branches/_devel_...
Repository Root: http://genie.hepforge.org/svn
Repository UUID: cc9776de-3512-45ca-aafc-e2d9ed43c22c
Revision: 3984
Node Kind: directory
Last Changed Author: sboyd11
Last Changed Rev: 3981
Last Changed Date: 2013-10-04 03:26:08 -0500 (Fri, 04 Oct 2013)

Here is an example update with conflict resolution (choose the repository version):

$ svn up
U    config/master_config.xml
U    config/Messenger.xml
U    config/UserPhysicsOptions.xml
A    config/AlvarezRusoCOHXSec.xml
A    config/COHXSecAR.xml
U    src/scripts/setup/genie-config
U    src/CrossSections/GSLXSecFunc.cxx
U    src/CrossSections/LinkDef.h
U    src/CrossSections/COHXSecAR.h
U    src/CrossSections/GSLXSecFunc.h
U    src/CrossSections/COHXSecAR.cxx
U    src/Conventions/KinePhaseSpace.h
D    src/Coherent/COHHadronicSystemGeneratorAR.cxx
D    src/Coherent/COHPrimaryLeptonGeneratorAR.cxx
D    src/Coherent/COHKinematicsGeneratorAR.h
D    src/Coherent/COHHadronicSystemGeneratorAR.h
D    src/Coherent/COHKinematicsGeneratorAR.cxx
D    src/Coherent/COHPrimaryLeptonGeneratorAR.h
U    src/Coherent/LinkDef.h
D    src/Coherent/COHPrimaryLeptonGenerator.cxx
A    src/Coherent/COHPrimaryLeptonGenerator.cxx
Conflict discovered in 'src/Coherent/COHElKinematicsGenerator.cxx'.
Select: (p) postpone, (df) diff-full, (e) edit,
        (mc) mine-conflict, (tc) theirs-conflict,
        (s) show all options: tc
G    src/Coherent/COHElKinematicsGenerator.cxx
U    src/Coherent/COHKinematicsGenerator.h
D    src/Coherent/COHHadronicSystemGenerator.h
A    src/Coherent/COHHadronicSystemGenerator.h
U    src/Coherent/COHKinematicsGenerator.cxx
D    src/Coherent/COHPrimaryLeptonGenerator.h
A    src/Coherent/COHPrimaryLeptonGenerator.h
D    src/Coherent/COHHadronicSystemGenerator.cxx
A    src/Coherent/COHHadronicSystemGenerator.cxx
U    src/Coherent/COHElKinematicsGenerator.h
U    src/AlvarezRuso/AlvarezRusoCOHXSec.cxx
U    src/AlvarezRuso/coh_multidiff.cxx
A    src/AlvarezRuso/coh_multidiff_current.cxx
U    src/AlvarezRuso/LinkDef.h
A    src/AlvarezRuso/coh_multidiff_current.h
A    src/Numerical/RootNDIntegrator.cxx
A    src/Numerical/MultiDimIntegrator.h
A    src/Numerical/RootNDIntegrator.h
A    src/Numerical/Simpson2DWrap.cxx
A    src/Numerical/MultiDimIntegrator.cxx
A    src/Numerical/Simpson2DWrap.h
U    src/stdapp/Makefile
A    src/stdapp/gCohARLookup.cxx
U    Makefile
Updated to revision 3984.

Now when we run info we see the most recent change has updated.

$ svn info
Path: .
URL: http://genie.hepforge.org/svn/branches/_devel_...
Repository Root: http://genie.hepforge.org/svn
Repository UUID: cc9776de-3512-45ca-aafc-e2d9ed43c22c
Revision: 3984
Node Kind: directory
Schedule: normal
Last Changed Author: sboyd11
Last Changed Rev: 3981
Last Changed Date: 2013-10-04 03:26:08 -0500 (Fri, 04 Oct 2013)

Once we have update we need only gmake in our project directory with the environment configured (use your GENIE setup script).

Running the Code: The Basics

In order to run GENIE as a developer it is important to understand how to run specific interaction channels on specific targets. It is also important to understand how to run the code inside of a debugger. Adding logging commands and recompiling is a very painful way to make progress and in a loosely-coupled framework like GENIE, where so much of the configuration is stored in XML files and enforced only at run-time, it is difficult to rely on the compiler to find some classes of bugs.

Getting Set-up

It is useful to add a “set up” function to your .bash_profile (or the appropriate file if you use a different shell):

function gwork () {
  export GENIEWORKDIR=&lt;Your working area&gt;/GENIE/
  echo "Setting working directory to \$GENIEWORKDIR"
  cd \$GENIEWORKDIR
  source environment_genie.sh
}

Generating an Interaction Spline

These splines are used by the flux driver to calculate whether an interaction takes place. It additionally specifies which interactions are possible: processes that don’t exist in the spline file or have zero cross section at the neutrino’s energy will not be selected. We can use that to our advantage when developing a new model by building a spline that only includes our interaction channel. These splines only provide the total cross section though - to get the kinematics of the event we will need to evaluate the differential cross section during event generation.

Let’s generate cross-section splines for muon-flavor neutrinos on Carbon for CCQE events only:

nice gmkspl -p 14,-14 -t 1000060120 \
  -o ccqe_carbon_splines.xml --event-generator-list CCQE

If we run this for one event generator only, we will find a short XML file with contents like (slightly reformatted for presentation):

<?xml version="1.0" encoding="ISO-8859-1"?>

<!-- generated by genie::XSecSplineList::SaveSplineList() -->

<genie_xsec_spline_list version="2.00" uselog="1">

<spline
    name="genie::LwlynSmithQELCCPXSec/Default/nu:-14;
    tgt:1000060120;
    N:2212;
    proc:Weak[CC],QES;"
    nknots="44">
  <knot> <E>    0.01000 </E> <xsec>          0 </xsec> </knot>
  <knot> <E>   0.030466 </E> <xsec>          0 </xsec> </knot>
  <knot> <E>   0.050932 </E> <xsec>          0 </xsec> </knot>
  <knot> <E>   0.071399 </E> <xsec>          0 </xsec> </knot>
  <knot> <E>   0.091865 </E> <xsec>          0 </xsec> </knot>
  <knot> <E>    0.11233 </E> <xsec>          0 </xsec> </knot>
  <knot> <E>     0.1427 </E> <xsec> 8.487371629e-13 </xsec> </knot>
  <knot> <E>    0.18129 </E> <xsec> 4.103713443e-12 </xsec> </knot>
  <knot> <E>     0.2303 </E> <xsec> 8.367696679e-12 </xsec> </knot>
  <knot> <E>    0.29257 </E> <xsec> 1.290439962e-11 </xsec> </knot>
  <knot> <E>    0.37168 </E> <xsec> 1.807792243e-11 </xsec> </knot>
  ...

Copy this file to your \$XSECSPLINEDIR. Note that for single-model work, the value of the cross section (in milibarns) need not be exactly correct and sometimes it is convenient to “hack” a spline file (for example, when working on a differential cross section formula for a new model in order to be able to generate a spline file).

See config/EventGeneratorListAssembler.xml for a list of valid event generator lists. If we check there, we find entries like:

<param_set name="CCQE">
   <param type="int" name="NGenerators"> 1                            </param>
   <param type="alg" name="Generator-0"> genie::EventGenerator/QEL-CC </param>
</param_set>

We can see that if we specify a list containing only CCQE, that the number of generators is set to one and the algorithm “category” is identified. We can see what that means by investigating \$GENIE/config/UserPhysicsOptions.xml:

 <!--
  ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  Specify which cross section model is to be used by each GENIE
  event generation thread. The parameter name is build as:
  "XSecModel@[name of thread]"
  -->
<param type="alg" name="XSecModel@genie::EventGenerator/QEL-CC">
  genie::LwlynSmithQELCCPXSec/Default
</param>

Here we find a useful comment and the actual algorithm used. It is a GENIE convention to name files after classes, so we can quickly find the code with (in the GENIE source directory):

find . -name "*LwlynSmithQELCCPXSec*" | grep -v svn | grep -v ROOT
./LlewellynSmith/LwlynSmithQELCCPXSec.cxx
./LlewellynSmith/LwlynSmithQELCCPXSec.h

Generating Events on Simple Targets

For physics work, we obviously want to use the rich detector geometry and flux facilities offered by GENIE, but for initial coding development work,

Create a shell script to run your simple job, do_a_ccqe_run.sh:

#!/bin/sh
gevgen -n 5 -p 14 -t 1000060120 -e 1,2 --run 100 \
  -f 'x*exp(-x)' \
  --seed 2989819 --cross-sections \$XSECSPLINEDIR/ccqe_carbon_splines.xml \
  --event-generator-list CCQE

(Be careful of trailing spaces after the line-continuation characters if copying-and-pasting this code.) It is also useful to specify mono-energetic beams. You can do this by dropping the -f flag and giving just one number to the -e flag.

Running Inside gdb

gdb is a large and full-featured program and it is impossible to go over all the features here. Instead, we will present a short list of common commands. These will be demonstrated in the walkthrough.

Important
For Mac users, Apple’s recent decision to stop supporting the GNU toolchain makes using gdb painful on a Mac. There are comments online suggesting ways to get gdb working again, but I have had limited success. The alternative is to use lldb as a debugger instead. Google searches will turn up tutorials and gdb to lldb cheatsheets. What lldb lacks (at the moment) is a “graphical” mode that shows large sections of the source code concurrent with the stepper. We are investigating an extension that may solve this problem. In the interim, please try to run this tutorial on Linux or an older version of OSX that supports gdb. We are workging on an lldb version of the instructions and walkthrough.

Note that gdb will automatically extend unambiguous abbreviations of words, so b, br, bre, brea, and break are all equivalent. However, f is ambiguous because there are commands file, finish, flushregs, and so on.

gdb Commands
  1. n - continue to the next line.

  2. s - step into the next function call.

  3. fin - step out of a function call.

  4. c - continue to the next break point.

  5. u - finish the current loop. (Note: this can be a bit slow. You may wish to use a break after the loop instead.)

  6. b - set a breakpoint at the current line.

  7. info b - look at the existing set of break points.

  8. disable N - disable breakpoint N. disable will disable all breakpoints.

  9. enable N - disable breakpoint N. enable will enable all breakpoints.

  10. p var - print a variable’s value. This has surprising power. print i to look at a loop counter, but also print *p to de-reference a pointer, and even print p->method(). You can use tab-completion with p.

  11. Control-l - Clean the screen (very useful after large blocks of text go to stdout).

  12. Control-x o - Switch focus in TUI mode (e.g., up and down arrow keys will scroll code or move through gdb history).

  13. Control-x Control a - Switch TUI mode on and off.

Put the set args... command on one line. Environment variable resolution (e.g., \$XSECSPLINEDIR) is temperamental in gdb, so it will be useful to save a few fully-formatted commands in local text files.

gdb -tui gevgen

set args -n 5 -p 14 -t 1000060120 \
  -e 0,10 --run 100 -f 'x*exp(-x)' --seed 2989819 \
  --cross-sections <Your Dev>/data/ccqe_carbon_splines.xml \
  --event-generator-list CCQE

Note: you will often want to turn the default logging verbosity down. Copy the file Messenger_laconic.xml over the default Messenger.xml (save the default if you’d like). You may even prefer to create a Messenger_fatal.xml file with all output levels set to FATAL only. (The logging output is not needed while inside gdb since you may inspect variables directly.)

gdb: Printing Variables

We can print (or just p) variable names and even de-reference pointers:

(gdb) p i
(gdb) p fXSecModel
(gdb) p *fXSecModel

Printing a std::string is sometimes tricky. If you’re having trouble, try (for string “inpfile”):

(gdb) p *(char**)inpfile._M_dataplus._M_p
$3 = 0xa682038 "/minerva/app/users/perdue/GENIE/data/coh_carbon_splines.xml"

gdb: Breakpoint Commands

We can break on a line in the file we’re stepping through:

(gdb) break 158
Breakpoint 13 at 0x2b8e0f0c4781: file COHKinematicsGenerator.cxx, line 158.

Then, we can add a set of commands for gdb to execute every time we hit that breakpoint:

(gdb) commands 13
Type commands for when breakpoint 13 is hit, one per line.
End with a line saying just "end".
>silent
>printf"xsec =%g ; Q2 = %f ; y= %f ; t = %f \n", xsec, gQ2, gy, gt
>end
(gdb) c
Continuing.
xsec =5.17161e-16 ; Q2 = 0.602666 ; y= 0.730704;   t = 0.191961
Continuing.
xsec =7.96193e-15 ; Q2 = 0.415240 ; y= 0.625730; t = 0.148100
Continuing.

gdb: Conditional Breakpoints

It is also often useful to only break if a certain condition is met:

(gdb) break 119 if tpi < 0
(gdb) info b
...
10      breakpoint     keep y   ...
        stop only if tpi<0

We can turn conditions off with condition number where number is the breakpoint number and then reset them if we choose:

(gdb) condition 11
Breakpoint 11 now unconditional.
(gdb) condition 11 xsec      > 5e-10
(gdb) info b
Num     Type           Disp EnbAddress               What
...
11      breakpoint     keep y   0x00002ad85a8c3fa4 in ...
        stop only if xsec > 5e-10
        breakpoint already hit 254 times
        printf "xsec = %g\n", xsec

Running valgrind

Our goal is naturally that there be no resource leaks, but that is impossible given the resource leaks that invariably appear in third-party code. ROOT [Brun:1997pa], for example, contains some leaks by design. The most we can expect is that new work on model development will not introduce any new resource leaks.

We will want to run valgrind in a deep diagnostic mode. Edit a file in your GENIE run directory called .valgrindrc and include these lines:

--num-callers=50
--leak-check=full
--verbose
--show-reachable=yes
--suppressions=valgrind-genie.supp

The first step towards sensible output with ROOT is to create a suppressions file. ROOT is distributed with a default file in \$ROOTSYS/etc, but the default version may not work with older releases of valgrind. In that case, first create a suppressions file with this script:

#!/bin/sh
valgrind --log-file-exactly=tmp.txt --gen-suppressions=all \
 gevgen -n 5 -p 14 -t 1000060120 -e 0,10 --run 100 \
 -f x*exp\(-x\) --seed 2989819 \
 --cross-sections \$XSECSPLINEDIR/gxspl-vA-v2.8.0.xml

Next, we will want to scrub the log file to get only the suppressions. Use the Perl script trim_suppressions_log.pl:

#!/usr/bin/env perl

# Run on the captured output of a valgrind job with suppression
# generation turned on to remove everything but the suppression
# messages.

if ($ARGV < 2 ) {
  print " USAGE:\n";
  print "   $0 [-i valgrind log file] [-o output suppressions file]\n";
  exit 1;
}

foreach $argnum (0 .. $ARGV) {
  print $argnum."\n";
  if (lc($ARGV[$argnum]) eq "-i") {
    $logfile = $ARGV[$argnum+1];
  }
  elsif (lc($ARGV[$argnum]) eq "-o") {
    $suppfile = $ARGV[$argnum+1];
  }
}
print $logfile."\n";
print $suppfile."\n";

open INPUT, $logfile or die "Can't open raw valgrind log file!";
open OUTPUT, ">$suppfile" or die "Can't open suppressions file!";

while (<INPUT>) {
  chomp;
  if ( !(/==\d+==/ or /--\d+--/) ) { print OUTPUT $_."\n"; }
}

With the trimmed suppressions file in hand, we can run valgrind and check for resource leaks using only our development model:

#!/bin/sh

NUMEVT=1
if [ $# -gt 0 ]; then
  NUMEVT=$1
fi

EXE=gevgen
A1=$(echo -e "-n $NUMEVT -p 14 -t 1000060120 -e 1,2")
A2=$(echo -e "--run 100 -f x*exp(-x) --seed 2989819")
A3=$(echo -e "--cross-sections $XSECSPLINEDIR/ccqe_carbon_splines.xml")
A4=$(echo -e "--event-generator-list CCQE")
DAT=`date -u +%s`
VALLOG="${DAT}_CCQE_grindlog.txt"

valgrind --log-file-exactly=${VALLOG} ${EXE} ${A1} ${A2} ${A3} ${A4}

cp $VALLOG leaksum.txt
perl -i -e 'while(<>) { chomp; if (/definitely/) { print $_,"\n"; } }' leaksum.txt
cat leaksum.txt
echo

(The different argument pieces may all be combined into one line, they are only broken up here for presentation purposes. Only the line producing the function requires the -e flag on the echo.)

It is a good idea to periodically run the memory test and look for changes. Once a resource leak enters the code, it can be challenging to track down.

Running the Code: The CCQE Walkthrough

This section assumes you have a cross section spline and have tested the steps outlined in [Running-GENIE-Basics].

Setting up Breakpoints in Advance

The following breaks can be very useful for understanding the flow control, and are valid in the frozen release v2-8-0. Put these commands in a file named .gdbinit in the directory where you will run GENIE.

Warning
You may also place .gdbinit in your \$HOME directory, but this sets values for all instances of gdb globally.
# we can't set breakpoints for .so's we haven't loaded yet without this
set breakpoint pending on

# can often get away with typing: pv or pvis, etc.
define pvisit
  p visitor-&gt;Id().Key()
end

# print C++ standard strings
define pcppstr
 p *(char**)(\$arg0)._M_dataplus._M_p
end

#
# Breakpoints ----------
#

# 1
break main

# 2
break gEvGen.cxx:274

# 3
# watch the loop through event record visitors!
# top of the loop. `p visitor-&gt;Id().Key()` (pvis above)
# check the processor and see if you want to step in
break EventGenerator.cxx:107 if ffwd==false
commands 3
pvisit
end

# 4 actual call to ProcessEventRecord(event_rec)
break EventGenerator.cxx:118 if ffwd==false

# 5 once an event is accepted
break QELKinematicsGenerator.cxx:193
commands 5
silent
printf "xsec = %g ; E = %f ; Q^2 = %f\n", xsec, E, gQ2
end

# 6
break QELKinematicsGenerator.cxx:116

# 7
break KineGeneratorWithCache.cxx:78

# 8
break KineGeneratorWithCache.cxx:83

# 9
break QELKinematicsGenerator.cxx:499

# 10
break LwlynSmithQELCCPXSec.cxx:82

# 11
break LwlynSmithQELCCPXSec.cxx:114
commands 11
printf "xsec = %g; ml = %f; Q2 = %f", xsec, ml, -1.0*q2
end

# 12 - After algorithm loop
break EventGenerator.cxx:195

Of course, these breakpoints are fragile if you are changing the source code because they rely on line numbers. You may also use function names to make a more generic “starting” .gdbinit file (they may not always work though, while line numbers are reliable).

It is important to understand some of GENIE’s design patterns [wikidesignpatterns], [wikigof] when stepping through framework code. GENIE makes heavy use of the Singleton Pattern, Visitor Pattern, and Chain of Responsibility Pattern [gof]. The Singleton Pattern restricts object creation for a given class to one and only one instance. The Visitor Pattern separates an algorithm from an object structure by moving the hierarchy of methods into one object. The Chain of Responsibility Pattern delegates commands to a chain of processing objects.

These patterns are how you see pointers on the stack accessing global data, and why sometimes loops over vectors of generic classes calling the same method over and over are able to process such a wide variety of behaviors with very compact code.

Don’t worry overly much about the patterns in the main framework. Stepping through them with gdb will make their operation clear in any case.

Feel free to explore more along the way and step into whatever you like. It is always okay to set a break point just after the function you step into so you leap out with a continue at any time.

First,

gdb -tui gevgen

Then, at the gdb prompt,

set args -n 5 -p 14 -t 1000060120 -e 3 --run 100 \
--seed 2989819 \
--cross-sections /minerva/app/users/perdue/GENIE/data/ccqe_carbon_splines.xml \
--event-generator-list CCQE

It is convenient to simply copy and paste these lines. Note that we use “\” to break the line (the whole command needs to be on one line). Then,

run

Or just r. Try to solve all of the exercises below. Note that you may not have the complete set of tools the first time through, or you may come to a conclusion and then change your mind later after you learn more about GENIE and how to investigate with gdb. Now, for the walkthrough:

n

We were automatically set to break on main. Move past GetCommandLineArgs

s

Step into Initialize. You may need to ctrl-l in order to clear logging text - keep this in mind going forward! If the screen gets cluttered with logging text (gdb does not redraw cleverly), use ctrl-l to get an understandable screen.

b 237

Set a breakpoint at the end of this function. We’ll test it shortly.

inf b

Note the number of the breakpoint you just set.

Exercise

What is the RunOpt::Instance()?

n

Step down to XSecTable (line 233).

s

Step into XSecTable.

Exercise

Try to figure out what is going on in this function. Step into anything you need to and experiment with fin to step out (and back up into XSecTable).

c

“Continue” once you are satisfied in XSecTable. You should end up on line 237 of gEvGen.cxx back where we just set our breakpoint.

del #

Delete the breakpoint we set (we won’t need it any longer).

inf b

Verify the breakpoint is gone.

n

Step back into main. You only need to press n once; after that you may simply keep pressing enter.

n

Step down to GenerateEventsAtFixedInitState (line 221).

s

Step into GenerateEventsAtFixedInitState.

n

Step down over these functions down to the break at line 274. There are many important things happening along the way here for setting up the run. For the purposes of understanding model development, we’ll skip them for now. (They are important for understanding the GENIE framework though.) Depending on your interests, you may want to go ahead and step through the instantiation of the InitialState to get a quick sense for what is happening there.

Note

From this point, is is understood you need to use n to step through from one instruction to the next where appropriate. We will be causual with its use.

s

Step into GenerateEvent (274). We have a break at that line, so if you are “far away,” you may press c to continue to that line.

b 276

We don’t want to accidently go past line 276.

s

Walk down through GenerateEvent and step into SelectInteraction at the breakpoint we just set (GEVGDriver.cxx:276).

b 178

We have now stepped into PhysInteractionSelector::SelectInteraction. Set a breakpoint on line 178.

b 214

Set a breakpoint on line 214.

Exercise

Walk through the code in this method and find out the names of the cross section algorithms as they are loaded. You will need to use the gdb print command and may need to dereference the variable if it is a pointer.

Exercise

When you get to line 178, print the size of the xseclist. Remember you can tab-complete variable names and even methods! Then you may delete the breakpoint. (But you aren’t obligated to.)

Exercise

When you get to line 214, try the following commands (first, recall that ctrl-x o will switch the focus between the code and command windows of gdb while in TUI mode - this will enable us to press “up” to go back through command history): p evrec, p *evrec, p *(evrec->fInteraction), p *(evrec->fInteraction->fProcInfo), and p evrec->fInteraction->fProcInfo->fScatteringType. Be careful about when we are and aren’t dereferencing a pointer. (Tab completion may not work once you’re several methods deep in that chain.) What do these commands show you? Are you processing a QuasiElastic event? If not, make sure you followed the steps at the beginning of this section and in [Running-GENIE-Basics] correctly. Finally, you may delete the breakpoint.

Note

We just learned some nice tricks for inspecting objects on the heap!

n

Continue stepping through the program until you reach line 322 of GEVGDriver.cxx.

s

Step into ProcessEventRecord on line 322.

n

Walk down to our breakpoint on line 117.

Exercise

Hitting the breakpoint should trigger a print statement. What is the information it is telling us?

c

Continue to the breakpoint on line 118. Step into the function if you would like to see what is happening, otherwise, continue.

Exercise

Later, return to and study these algorithms we’re skipping now. What do they do?

c

Continute again to the breakpoint on line 107. What is the name of the algorithm being processed now.

c

Keep pressing c until the break on line 107 tells us the active algorithm is QELKinematicsGenerator. Then press c again to get to line 118.

s

Step into the Kinematics generator.

n

Walk down to the breakpoint on line 116.

Exercise

What are the Q-squared limits for this event?

s

Step in at the breakpoint on line 116.

Exercise

We end up in KineGeneratorWithCache::MaxXSec. Based purely on the class name, what do you think this class does. What does the method do?

c

Continue to the breakpoint on line 78.

Exercise

Step into FindMaxXSec (on line 78). What does this function do? When you’re satisfied, type fin to exit the funciton.

inf b

What is the number of the breakpoint we just crossed?

dis #

Disable the breakpoint we just crossed by number.

c

Continue to the breakpoint on line 83.

s

Step into ComputeMaxXSec. Then, walk down to the breakpoint on line 499.

Exercise

What is going on in this function? How are we computing the maximum cross section?

s

Step into the breakpoint on line 499. We are now in the XSec funciton. Walk down to the breakpoint on line 82.

s

Step into the breakpoint on line 82.

Exercise

We are now in QELFormFactors::Calculate. What do the functions here do? How are the form factors organized? Step into whichever functions look interesting to you. It might be good to set a breakpoint on line 68 so you can hop out quickly if you dig deep. Alternatively, you might want to use the fin command to hop out functions.

Note

If you use fin, gdb may set temporary helper breakpoints. Don’t worry, it will clean them up.

n

When you get back into LwlynSmithQELCCPXSec::XSec, keep walking through the function and watch what is happening. When we hit the breakpoint on line 14, note the cross section. Walk the rest of the way out of the function.

disable

When we get back into QELKinematicsGenerator::ComputeMaxXSec, disable all breakpoints.

u

Walk through the rest of the for loop, but when you reach line 496, use until, or just u to finish the loop rather than stepping through it again. (Although you should feel to step through it again to watch what happens and make sure you understand what is going on.)

enable

Now turn all breakpoints back on.

n

Keep walking and note what happens next in KineGeneratorWithCache::MaxXSec.

n

Once we get back to QELKinematicsGenerator::ProcessEventRecord, keep walking.

Exercise

By what mechanism do we finally accept a cross section value? Note that while walking through, we will end up stopping at our breakpoints inside LwlynSmithQELCCPXSec::XSec. Use info b to figure out their number, or look at the number gdb reports and disable them (unless you want to watch them).

c

Once you are satisfied with mechanics of how we accept a cross section, continue to the breakpoint on line 193.

n

Walk through the rest of the function and note what is going on. What happens between lines 230 and 234? fin when you are satisfied.

c

We are now in EventGenerator::ProcessEventRecord again. Continue to the breakpoint on line 107. What is the next algorithm we will process? Step through it if you are interested.

c

Continue until we are alerted the next algorithm is QELHadronicSystemGenerator. Continue again and step into the breakpoint on line 118.

s

Inside QELHadronicSystemGenerator::ProcessEventRecord, step into AddRecoilBaryon. Walk through that function and back out to EventGenerator::ProcessEventRecord. Continue again to the breakpoint on line 107. Keep continuing and noting the names of the algorithms we’re walking past. Step into any that look interesting.

Exercise

Check out what is going on in NucDeExcitationSim. Can you figure out what the most important element is in the far detector of the oscillation experiment that Costas happens to belong to?

n

After the algorithm loop, continue to the breakpoint on line 195 and walk out the method, or use fin to get out. (The breakpoint on 195 is mostly to guard against over-zealous continuing in the algorithm loop.)

n

Back in GEVGDriver::GenerateEvent, keep walking and observing. Walk all the way until we end up in GenerateEventsAtFixedInitState in gEvGen.cxx. Keep walking until we return to the top of the while loop on line 269. Feel free to step into anything interesting along the way.

Congratulations!

We finished a walthrough of one event.

disable

Disable all breakpoints.

c

Continue to the end. There are no breakpoints, so this will end the program.

r

You may press r at any time to re-run. You may even do this while the program is running. In fact, a common (if somewhat unstable) situation is to find a bug while stepping through code, fix it, recompile, and then press r to restart. You do not need to exit gdb in this case! This way your current breakpoints, etc. are preserved.

Where did the long list of algorithms we stepped through come from? For CCQE, the classes that process the event can be found and specified in config/EventGenerator.xml:

<param_set name="QEL-CC">
   <param type="string" name="VldContext"> </param>
   <param type="int"    name="NModules"> 12 </param>
   <param type="alg"    name="Module-0">
    genie::InitialStateAppender/Default
   </param>
   <param type="alg"    name="Module-1">
    genie::VertexGenerator/Default
   </param>
   <param type="alg"    name="Module-2">
    genie::FermiMover/Default
   </param>
   <param type="alg"    name="Module-3">
    genie::QELKinematicsGenerator/CC-Default
   </param>
   <param type="alg"    name="Module-4">
    genie::QELPrimaryLeptonGenerator/Default
   </param>
   <param type="alg"    name="Module-5">
    genie::QELHadronicSystemGenerator/Default
   </param>
   <param type="alg"    name="Module-6">
    genie::PauliBlocker/Default
   </param>
   <param type="alg"    name="Module-7">
    genie::UnstableParticleDecayer/BeforeHadronTransport
   </param>
   <param type="alg"    name="Module-8">
    genie::NucDeExcitationSim/Default
   </param>
   <param type="alg"    name="Module-9">
    genie::HadronTransporter/Default
   </param>
   <param type="alg"    name="Module-10">
    genie::NucBindEnergyAggregator/Default
   </param>
   <param type="alg"    name="Module-11">
    genie::UnstableParticleDecayer/AfterHadronTransport
   </param>
   <param type="alg"    name="ILstGen">
    genie::QELInteractionListGenerator/CC-Default
   </param>
</param_set>

One very important feature of this approach is that these classes can be swapped or re-ordered. Of course, re-ordering these specific classes might break the event, but nothing in GENIE forces us to handle the target first or the final-state lepton second, etc. We can sequence the appropriate set of computations as required by the physics and we can even make some changes dynamically (without recompiling).

Configuration XML

GENIE manages much of its configuration via XML text files and specific conventions. The compiler cannot help you find problems if you are using a non-sensical configuration file or have typo in your specification XML. On the other hand, the wide use of XML enables a very loose coupling between GENIE components, making maintenance easier, and allows for dynamic run-time flexibility and decision-making. The secret is to understand and use the XML configuration correctly and make heavy use of grep and other tools to look for typos inside the strings.

The AlgConfigPool class controls much of the parameter state in GENIE. Many of these files are loaded when the class is created in methods like AlgConfigPool::LoadAlgConfig().

Important Files

Physics Parameters
  1. config/UserPhysicsOptions.xml: Contains important physical parameters (e.g. CKM matrix parameters, Axial Mass for quasi-elastic scattering, etc.). This is a good place to start if you want to change something about the physics in GENIE, but you aren’t sure where to look first. It also contains some model execution options (e.g., include lepton mass correction in charged-current coherent pion production or not, etc.).

We also often set the algorithm names for implementing cross section models here, e.g.

<param type="alg" name="XSecModel@genie::EventGenerator/COH-CC">
  genie::ReinSeghalCOHPiPXSec/Default
</param>
<param type="alg" name="XSecModel@genie::EventGenerator/COH-NC">
  genie::ReinSeghalCOHPiPXSec/Default
</param>

To replace the “COH-CC” model with a different one, we supply a different string here, e.g.:

<param type="alg" name="XSecModel@genie::EventGenerator/COH-CC">
  genie::BergerSehgalCOHPiPXSec/Default
</param>
<param type="alg" name="XSecModel@genie::EventGenerator/COH-NC">
  genie::BergerSehgalCOHPiPXSec/Default
</param>
Event Generator Processes
  1. config/EventGeneratorListAssembler.xml: Contains physics lists. In general, you will want to use a comprehensive list for physics work, but for development, it can be very useful to specify a target list like:

<param_set name="CCMEC">
  <param type="int" name="NGenerators"> 1                            </param>
  <param type="alg" name="Generator-0"> genie::EventGenerator/MEC-CC </param>
</param_set>

Here, NGenerators is accessed in EventGeneratorList::AssembleGeneratorList() to count the number of event generation algorithms to load. . config/EventGenerator.xml: For each process specified in config/EventGeneratorListAssembler.xml, the list of event record visitors is specified. An event record visitor implements the interface EventRecordVisitorI and the execution process will loop over all of these in order and call the ProcessEventRecord method. Because this list may specified here in XML, we can set the list of computations dynamically and flexibly. These lists look like:

<param_set name="MEC-CC">
   <param type="string" name="VldContext"> </param>
   <param type="int"    name="NModules">   5                                                   </param>
   <param type="alg"    name="Module-0">
      genie::InitialStateAppender/Default                 </param>
   <param type="alg"    name="Module-1">
      genie::VertexGenerator/Default                      </param>
   <param type="alg"    name="Module-2">
      genie::MECGenerator/Default                         </param>
   <param type="alg"    name="Module-3">
      genie::HadronTransporter/Default                    </param>
   <param type="alg"    name="Module-4">
      genie::UnstableParticleDecayer/AfterHadronTransport </param>
   <param type="alg"    name="ILstGen">
      genie::MECInteractionListGenerator/CC-Default       </param>
  </param_set>

The algorithms are specified as “name/configuration” - this informaiton is parsed when the algorithm is instantiated in AlgFactory::GetAlgorithm(string name, string config). The strings are not arbitrary - check, for example, EventGenerator::LoadConfig() to see how GENIE looks for algorithms keyed by “Module-”.

General Issues in Model Development

There are a number of conventions for organizing details of the model. Developers should strive to implement their code in accordance with these conventions to promote maintainability.

Conventions for Organization
  1. Developers should be aware of the Utils/KineUtils class. This class collects kinematic limits and methods, for example, channel specific (QEL vs Coherent, etc.) limits on Q^2, W, and others as well as methods for variable transformations, etc.

  2. The Interaction/KPhaseSpace class collects kinematic limits in from KineUtils and offers them through a more generic interface in addition to making calculations as to whether the required phase space is available for a reaction to occur (for example, is there enough energy to produce all the final state particles, etc.).

Implementing the Berger-Sehgal Model

Priorities
  1. We need to fix the kinematic term in the triple-differential cross-section that is currently valid only at Q^2 = 0.

  2. We need to fix the kinematics to reflect the muon mass.

  3. The pion scattering cross-section makes a large contribution and needs to be fixed.

  4. Event kinematics are problematic.

How to get the code to build?

How to best test? Is it possible to run coherent only (I know it is - I did some study a long time ago in this fashion - but what is the most effective)?

Elastic
  1. COHElasticPXSec.h

  2. COHElHadronicSystemGenerator.h

  3. COHElKinematicsGenerator.h

General Coherent
  1. COHHadronicSystemGenerator.h

  2. COHInteractionListGenerator.h

  3. COHKinematicsGenerator.h

  4. COHPrimaryLeptonGenerator.h

Alvarez-Russo
  1. COHHadronicSystemGeneratorAR.h

  2. COHKinematicsGeneratorAR.h

  3. COHPrimaryLeptonGeneratorAR.h

Where to start for each piece?

First, the code in BergerSehgalCOHPiPXSec.{h,cxx} is identical to the ReinSeghalCOHPiPXSec.{h,cxx} files, except for a global replace of "ReinSeghal" with "BergerSehgal" (even refs are still 1983 paper). This did at least fix the spelling of “Sehgal.”

Does BergerSehgal have code already in existence, just not committed to the repository?

If I poke around in the source, I find these:

$ ls -l BergerSehgal/*.cxx
-rw-r--r-- 1 perdue e938 8373 Sep 30 10:26 BergerSehgal/BergerSehgalCOHPiPXSec.cxx

$ ls -l Coherent/*.cxx | grep -v El
-rw-r--r-- 1 perdue e938  8547 Oct  8 13:53 Coherent/COHHadronicSystemGenerator.cxx
-rw-r--r-- 1 perdue e938  4082 Sep 30 10:25 Coherent/COHInteractionListGenerator.cxx
-rw-r--r-- 1 perdue e938 25105 Oct  8 13:53 Coherent/COHKinematicsGenerator.cxx
-rw-r--r-- 1 perdue e938  4697 Oct  8 13:53 Coherent/COHPrimaryLeptonGenerator.cxx

During event generation we go through classes in the second group and they all have clauses like:

  if (fXSecModel->Id().Name() == "genie::ReinSeghalCOHPiPXSec") {
       CalculateHadronicSystem_ReinSeghal(evrec);
  } else if ((fXSecModel->Id().Name() == "genie::AlvarezRusoCOHXSec")) {
       CalculateHadronicSystem_AlvarezRuso(evrec);
  }

So they need to get the lines:

  } else if ((fXSecModel->Id().Name() == "genie::BergerSehgalCOHPiXSec")) {
       CalculateHadronicSystem_BergerSehgal(evrec);

along with the appropriate function implementations. And they need to know how to call BergerSehgal instead of ReinSehgal when a cross-section is needed.

In implementing this, we go from choosing a max cross section based on neutrino energy and a two-dimensional binning in x and y to a two-dimensional binning in x and $Q^2$. Is this reasonable? How do we choose the maximum $Q^2$? How do we choose the right step size and step spacing?

Conventions - add 3q (magnitude of q three-vector) to KineVar enums?

Update config/master_config.xml:

Add line:

<config alg="genie::BergerSehgalCOHPiPXSec"> BergerSehgalCOHPiPXSec.xml </config>

to <!--  ****** CONFIGURATION FOR XSEC ALGORITHMS ****** --> block.

In config directory, make sure we have a configuration file matching the name of the file specified in master_config.xml:

cp ReinSeghalCOHPiPXSec.xml BergerSehgalCOHPiPXSec.xml
Update config/UserPhysicsOptions.xml:

Add lines:

<param type="alg" name="XSecModel@genie::EventGenerator/COH-CC"> genie::BergerSehgalCOHPiPXSec/Default </param>
<param type="alg" name="XSecModel@genie::EventGenerator/COH-NC"> genie::BergerSehgalCOHPiPXSec/Default </param>

Comment out lines that set AlvarezRuso or ReinSeghal [sic] to the COH default.

Need to be sure BergerSehgal/LinkDef.h has the right class names defined.

Investigate double genie::utils::hadxs::TotalPionNucleonXSec(double Epion)

Also, utils::hadxs::InelasticPionNucleonXSec(Epi)
Add BergerSehgal to the top-level Makefile (follow, for example, AlvarezRuso pattern).

Update libs var in src/scripts/genie-config


Update with BergerSehgal functions in:

Coherent/COHHadronicSystemGenerator.h
Coherent/COHKinematicsGenerator.h
Coherent/COHPrimaryLeptonGenerator.h

Coherent/COHKinematicsGenerator.cxx
Coherent/COHHadronicSystemGenerator.cxx
Coherent/COHPrimaryLeptonGenerator.cxx

Be sure before running the code that you have the model represented in the cross section spline. You may need to edit a file by hand for quick-and-dirty testing, but you will also need to eventually be able to run gmkspl to produce the splines correctly. This may require you to update and/or add an integrator, so it is reasonable to put that step off at the very beginning and use a “by hand” file with fake numbers and a monochromatic neutrino beam for early-stage testing.

Bibliography

  • [gof] E. Gamma, R. Helm, R. Johnson, and J. Vlissides, "Design Patterns: Elements of Reusable Object-Oriented Software", Addison-Wesley, 1994.

  • [Brun:1997pa] Brun, R. and Rademakers, F. "ROOT: An object oriented data analysis framework". Nucl.Instrum.Meth. A389. 81-86. 1997.

  • [Whalley:2005nh] Whalley, M.R. and Bourilkov, D. and Group, R.C. "The Les Houches accord PDFs (LHAPDF) and LHAGLUE". 2005. hep-ph/0508110.

  • [Sjostrand:2006za] Sjostrand, Torbjorn and Mrenna, Stephen and Skands, Peter Z. "PYTHIA 6.4 Physics and Manual". JHEP 0605. 2006

  • [Andreopoulos:2009rq] Andreopoulos, C. and Bell, A. and Bhattacharya, D. and Cavanna, F. and Dobson, J. and others. "The GENIE Neutrino Monte Carlo Generator". Nucl.Instrum.Meth. A614. 87-104. 2010.

  • [geniemanual] The GENIE Neutrino Monte Carlo Generator. http://genie.hepforge.org/manuals/GENIE_PhysicsAndUserManual_20130615.pdf. 2013.

  • [wikidesignpatterns] Software design pattern. http://en.wikipedia.org/wiki/Design_pattern_(computer_science). 2013

  • [wikigof] Design Patterns. http://en.wikipedia.org/wiki/Design_Patterns. 2013

  • [atfnal] GENIE at Fermilab. https://cdcvs.fnal.gov/redmine/projects/genie/wiki