IceTray Introduction

What is IceTray?

IceTray is the software framework, IceCube is using to process data, perform simulations, do reconstructions etc.
It consists of so-called modules and services.

Connection to WIPAC computation resources

The WIPAC computer resources can be utilized for IceTray computations. The icerec and simulation meta-projects are available in cvmfs. You can connect to the cobalt machines via:

ssh -X username@cobalt.icecube.wisc.edu

The -X enables a graphical interface, so you don't have to work in the shell only, but can use editors such as gedit. However, it is not recommended to use heavy load graphical applications, such as steamshovel or eclipse, as this will be very slow, especially if you are outside of the 222 network.

Layout of an IceTray script

IceTray script are python scripts. We want to make a directory to store our scripts: mkdir /home/icecube/i3/scripts
and switch to that directory: cd /home/icecube/i3/scripts.
Then we create a new file: touch my_first_icetray_script.py

When writing an IceTray script, there is a certain minimum of code we have to set-up to get going:

In []:
#!/usr/bin/env python

from icecube import icetray, dataio

from I3Tray import *

tray = I3Tray()

tray.Add("I3Reader","my_reader",Filename='/home/icecube/i3/data/IC86-2011_nugen_numu/Level2_nugen_numu_IC86.2011.011184.000000.i3.bz2')

tray.Execute()
tray.Finish()

A script starts with a so-called shebang line:
#!/usr/bin/env python

If you have several versions of Python installed, /usr/bin/env will ensure the interpreter used is the first one on your environment's $PATH. The alternative would be to hardcode something line #!/usr/bin/python or the like -- that's OK but less flexible. In Unix, an executable file that's meant to be interpreted can indicate what interpreter to use by having a #! at the start of the first line, followed by the interpreter (and any flags it may need).

For this exercise, we need certain icetray projects, which have to be imported from the topmost instance icecube:
from icecube import icetray, dataio

icetray is the heart and soul of basically any data-handling script you are going to write.
dataio allows you to read-in/write-out data files, e.g. i3.

Next, we have to import the I3Tray:
from I3Tray import *

In general, it is not recommended to use this syntax for most of the projects, as this might result in an increased start-up time of your program as many (unnecessary) things will be loaded.

First, we have to create the tray. This is basically our peon, that will do all the work we tell it to do, as long as our instructions are clear.
tray = I3Tray()

As our peon is eager to work, so we have to give it something to do ("add a module"). IceTray works on frames, so we either need to read an I3File that contains frames or create an empty container for frames (events) that we will fill later using a generator. Here, we will just read in one of the I3Files available in the VM:
tray.Add("I3Reader","my_reader",Filename='/path/to/file.i3')

Our peon shall execute on all frames in the input file:
tray.Execute()
In the parantheses we can also specify the number of frames it should process.

At last, we allow our peon to finish its tasks:
tray.Finish()

Notes:

We can execute this script now: python my_first_icetray_script.py.
Very likely you will run into this error:

In []:
Traceback (most recent call last):
  File "my_first_icetray_script.py", line 3, in <module>
    from icecube import icetray, dataio
ImportError: No module named icecube

In order to use the script, we need to load our IceTray software!

In []:
/home/icecube/i3_software/combo/build/env-shell.sh

Now we can use everything that was built when running cmake some time ago (all modules that show up in the I3Build folder).
Executing the script again should run with the following output:

In []:
INFO (I3Tray): Last module ("my_reader") has a dangling outbox. Adding TrashCan to end of tray (I3Tray.cxx:328 in void I3Tray::Configure())
INFO (I3Module): Opened file /home/icecube/i3/data/IC86-2011_nugen_numu/Level2_nugen_numu_IC86.2011.011184.000000.i3.bz2 (I3Reader.cxx:180 in void I3Reader::OpenNextFile())
NOTICE (I3Tray): I3Tray finishing... (I3Tray.cxx:464 in void I3Tray::Execute())

These are all Info and Notice messages, and the I3Tray actually finished, so the script ran successfully.
The first info is a convenience that was added some time ago. Previously, you had to add the TrashCan to the tray:
tray.Add("TrashCan","trash_it")

Now we actually want to give our script some functionality.
Our script will look like this in the end:

In []:
#!/usr/bin/env python

import sys

from icecube import icetray, dataio
from icecube import dataclasses

from I3Tray import *

outfile = sys.argv[1]
infiles = sys.argv[2:]

tray = I3Tray()

tray.Add("I3Reader","my_reader",FilenameList=infiles)

tray.Add("Dump")

def PrintEventHeaderID(frame):
    # ToDo: print the Event ID to the shell

tray.Add(PrintEventHeaderID, Streams=[icetray.I3Frame.DAQ])

tray.Add("I3Writer", filename=outfile)

tray.Execute(10)
tray.Finish()
First, we want to read the filename from the command line instead of having it hard-coded. Since I am lazy, I am using argv[] for this purpose. To use it, we have to import python standard library modules:
import sys

So we want to write out in the first command line argument and take all of the following arguments as input files:
outfile = sys.argv[1] infiles = sys.argv[1:-1]

As we now want to work with data from the file, we need to import dataclasses.
dataclasses gives access to data objects such as I3Particle, I3EventHeader, I3RecoPulse and allows manipulating them.
We can have a look at the dataclasses project in the repository: http://code.icecube.wisc.edu/svn/projects/dataclasses/releases/V15-04-01/
Most projects contain the following folders:

We now want to add a function that reads the I3EventHeader and extracts the event_id:

In []:
def PrintEventHeaderStartTime(frame): 
    #ToDo print Event ID

We need to pass our function the frame it should work on. Then we check if the frame is actually present in the frame. Otherwise, if it does not exist and we try to access it anyway, we will get an error.
How do we actually know how to call the event ID within the I3EventHeader?
The safest way to check is to look it up in the corresponding pybindings, in this case: http://code.icecube.wisc.edu/projects/icecube/browser/IceCube/projects/dataclasses/releases/V15-04-01/private/pybindings/I3EventHeader.cxx
Now we actually have to call our newly created function:
tray.Add(PrintEventHeaderStartTime, Streams=[icetray.I3Frame.DAQ])
We can specify on which type of frame our function should run. If nothing is specified, such functions will run on P-frames.

We can now modify our function to actually perform a cut. If we add a return statement to our function, it will only keep those P-frames where the return value is True. If none of the P-frames survived our cuts, we actually also want to get rid of the Q-frames, as they don't seem to be relevant for our analysis. This can be achieved by adding DropOrphanStreams=[icetray.I3Frame.DAQ] to the writer module.

In []:
def RemoveIceTop(frame):
    if frame["I3EventHeader"].sub_event_stream == 'ice_top': #ToDo cut all events except for those passing the MuonFilter
        return False
    else:
        return True
tray.Add(RemoveIceTop)

tray.Add("I3Writer", DropOrphanStreams=[icetray.I3Frame.DAQ], filename=outfile)
Since people are generally lazy and want to type as little as possible, you might encounter so-called Lambda Function.

Example:

In []:
g = lambda x: x**2
print g(8)

The above lambda function returns the value it is called with (8) to the power of 2, i.e. 64.

We can now replace the entire RemoveIceTop function with the following line:

In []:
tray.AddModule(lambda frame: frame["I3EventHeader"].sub_event_stream != 'ice_top')

Not all frame objects are actually relevant for each analysis, so we might want to get rid of them:

In []:
tray.AddModule("Delete","Cleanup", Keys = [ 'InIceRawData'])

Even more likely you might want to add new frame objects. As an example, we would like to sum up all charges in OfflinePulses and save them as a double:

In []:
def CountCharges(frame):
    #ToDo: count all charges contained in OfflinePulses and save them as a new frame object frame["TotalCharge"]

tray.Add(CountCharges)

If we want to add icecube modules to the tray, but do not know how to instantiate them, we can use the icetray-inspect feature:
icetray-inspect dataio
Here, all the modules are listed which we can use. Also, their arguments are shown.
We can also modify our function to print out the total charge:

In []:
def CountCharges(frame, do_print=False):
    #ToDo: count all charges contained in OfflinePulses and save them as a new frame object frame["TotalCharge"]
    if (do_print):
        print "Total Charge of OfflinePulses", frame["TotalCharge"].value

tray.Add(CountCharges, do_print=True)

Besides modules, we can also add services to the tray. While modules operate on frames and execute commands in the given order, services are classes that are available throughout the entire tray and can be used by several modules. Probably the most commonly used service is the RandomService, which provides you with a sequence of random numbers.

In []:
from I3Tray import *
from icecube import phys_services
tray = I3Tray()
tray.context["I3RandomService"] = phys_services.I3GSLRandomService(seed=1)
random_service = tray.context["I3RandomService"]
print random_service.uniform(1,5)

The RNGs live in the phys-services project. Note that python can't handle dashes, so we translated it to phys_services for use with python.
The service is stored in a dictionary, called context.
If we want to know which RNGs are available and how we have to call them, we can use icetray-inspect phys-services (here with dash!). A huge list of available modules will pop up and we have to pick an RNG, e.g. the I3GSLRandomService. If we do not specify a context, it will store the RNG in the context["I3RandomService"], which is also the default context where modules that require an RNG will search for it.
If we want to produce different sequences of random numbers, we have to change the seed.
The above example can e.g. be executed in ipython, which also allows you to see which kind of random distributions are available. Here, I chose a uniform distribution between 1 and 5, which should return the same value for everyone using the same seed.