Tutorials on Using the Project

Below is a description on how to use this class to simulate an event. Note that there is an example event located in $I3_TESTDATA/radcube/100000. This event was hand-made to some extent to be used for testing purposes. If you would like real simulated events, they are available on Cobalt (or you could make some yourself).

You can find example scripts to run in $I3_BUILD/radcube/resources/examples/ which include calling typical module sequences as well as custom functions to augment the data for plotting.

• The file ExampleScript_GetEFieldInVxB.py reads in a simulated event and converts the antenna locations to the \(\vec V \times \vec B\) coordinate system, and makes a plot of the log of the peak amplitude of the Hilbert envelope.

• The file ExampleScript_SimulateRadioEvent.py is an example of how to perform a full simulation of an event and save it to a file. You will see that this is actually fairly similar thanks to the RadioInjection segment.

There are two services which are important and unique to this project, I3AntennaResponse and I3ElectronicsResponse (see Services).

Use of Basic RadCube Modules

For a more detailed overview of what the various modules within radcube, see Radcube Modules, Services, and Segments.

In general a simulation if you want to do it “by hand” would follow the code below

#First we import a library to handle python inputs smoothly
import argparse
parser = argparse.ArgumentParser()
parser.add_argument('input', type=str, nargs='+', help='Input data files.')
args = parser.parse_args()

#import the standard libararies for icetray
from icecube.icetray import I3Tray
from icecube import icetray, dataio
from icecube.icetray import I3Units

#import radcube!
from icecube import radcube

tray = I3Tray()

#We make one copy of all the services that we need
#These will give you the "standard" setting for the respective services
#The names of the services are returned so you can pass them to the respective module(s)
antennaResponse = radcube.defaults.CreateDefaultAntennaResponse(tray)
electronicsResponse = radcube.defaults.CreateDefaultElectronicsResponse(tray)
randomService = radcube.defaults.CreateDefaultRandomService(tray)

#Alternatively you could create them yourself with non-default parameters
tray.AddService("I3ElectronicsResponseFactory", "MyNonDefaultER",
                 IncludeCables=True,
                 CableLength=100*I3Units.meter,
                 CableTemperature=274*I3Units.kelvin,
                 IncludeLNA=True,
                 IncludeRadioBoard=True,
                 IncludeTaxi=True,
                 AntennaType=antType,
                 IncludeElectronicNoise = False,
                 IncludeRadioBoardNoise = False,
                 InstallServiceAs="MyNonDefaultER"
                )

#Add the CoREAS reading module
tray.Add("CoreasReader", "CoreasReader",
         DirectoryList=args.input,      #List of CoREAS output directories (one per shower))
         MakeGCDFrames=True,            #Make a GCD frames from simulations
         MakeDAQFrames=True,            #Make the Q frames from simulations
        )

#Inject the electric fields into the antennas
tray.AddModule("ChannelInjector", "_ChannelInjector",
               InputName="CoREASEFieldMap",  #This is the default name that the CoreasReader puts in the frame
               OutputName="RawVoltageMap",
               AntennaResponseName=antennaResponse
              )

#Convolve with the electronics response
#Here we will use the default response we created above, but you can also use
#the non-default one as well. Plug and play!
tray.AddModule("ElectronicResponseAdder", "AddElectronics",
               InputName="RawVoltageMap",
               OutputName="FoldedVoltageMap",
               ElectronicsResponse=electronicsResponse
              )

#Make yourself a P Frame (CoreasReader only makes GCD and Qs)
tray.AddModule("I3NullSplitter","splitter",
               SubEventStreamName="RadioEvent"
              )

#Filter the voltages to only look at a narrower band in frequency space
#Here we will filter between 110 and 190 MHz using the Butterworth algorithm
tray.AddModule("BandpassFilter", "BandpassFilterButter",
               InputName="RadioEvent",
               OutputName="Filtered",
               FilterType=radcube.eButterworth,
               FilterLimits=[110*I3Units.megahertz,190*I3Units.megahertz],
               ButterworthOrder=radcube.constants.butterworthOrder
              )

#Add any other modules you want here for: analyzing, plotting, saving, etc.

tray.Execute()

Likewise one would run this code as

$ python the-script-defined-above.py $I3_TESTDATA/radcube/100000

Full Radio Reconstruction Using Segments

In general, there are two main segments for radcube, the RadioInjection and RadioReconPrep (defined in ./radcube/python/segments/). The first reads in coreas simulations, creates GCD frames, as needed, and makes one Q frame per passed in event. The second takes simulated files (expected to be in a format that comes from the Pole) and prepares them for reconstruction. This includes removing various antennas from the event, estimating the shower properties, etc.

These two segments will be kept up to date such that an “official” (or at least recommended) simulation would instead be the more simple:

import argparse
import os

parser = argparse.ArgumentParser()
parser.add_argument('--output', dest='outfile', type=str, default=str(os.getenv('I3_BUILD'))+'/simulated-radio-event.i3.gz', help='Select to add radio to input gcd file')
parser.add_argument('input', type=str, nargs='+', help='Input data files.')
args = parser.parse_args()

from icecube.icetray import I3Tray
from icecube import icetray, radcube, dataio
from icecube.icetray import I3Units
from icecube.icetray.i3logging import *
from icecube.radcube.defaults import *

tray = I3Tray()

#These lines create basic response services. See radcube.defaults for more
#The functions place them in the frame and give you back the name
antennaResponse = CreateDefaultAntennaResponse(tray)
electronicsResponse = CreateDefaultElectronicsResponse(tray)
randomService = CreateDefaultRandomService(tray)

#This segment reads in the CoREAS event, adds noise, convolves it with the
#antenna/electronics response, "digitizes" it, and puts it in the frame.
tray.AddSegment(radcube.segments.RadioInjection, 'RadioInjection',
                 InputFiles = args.input,
                 AntennaServiceName = antennaResponse,
                 ElectronicServiceName = electronicsResponse,
                 RandomServiceName = randomService,
                 OutputName = "InjectedVoltageMap"
                )

#Make yourself a P Frame (CoreasReader only makes GCD and Qs)
tray.AddModule("I3NullSplitter","splitter",
              SubEventStreamName="RadioEvent"
             )

#This segment takes the simulated data, converts from ADC bins to a
#voltage, removes electronics effects, and estimates shower parameters
tray.AddSegment(radcube.segments.RadioReconPrep, "ReconPrep",
                ElectronicServiceName = "electronicResponseName",
                InputName = "InjectedVoltageMap"
               )

tray.AddModule("I3Writer","i3writer",
               filename=args.outfile,
               # streams=[icetray.I3Frame.DAQ,icetray.I3Frame.Physics]
               # streams=[icetray.I3Frame.Geometry,icetray.I3Frame.Calibration, icetray.I3Frame.DetectorStatus]
               )

tray.Execute()

Again, this sctipt can be run by passing in an event

$ python the-script-defined-above.py $I3_TESTDATA/radcube/100000

You can also pass in an output location to put the file

$ python the-script-defined-above.py --output ./myfile.i3.gz $I3_TESTDATA/radcube/100000