Python bindings to dataclasses¶
Examples¶
Get an object out of the frame¶
With a frame like this:
>>> from icecube import icetray, dataclasses, dataio
>>> i3f = dataio.I3File("/v/icecube/data/level3.nugen_numu.000768.000111.i3.gz")
>>> frame = i3f.pop_physics()
>>> print frame
[ I3Frame (Physics):
'TrackParaboloidFit32' [Physics] ==> I3Particle
'FirstSubtrackBayesianFit' [Physics] ==> I3Particle
'TrackLlhFit32_e' [Physics] ==> I3Particle
'TrackBayesianFitFitParams' [Physics] ==> I3LogLikelihoodFitParams
'TrackLlhFit32' [Physics] ==> I3Particle
'ToIParams' [Physics] ==> I3TensorOfInertiaFitParams
'TrackLlhFitFitParams' [Physics] ==> I3LogLikelihoodFitParams
'TrackBayesianFit32FitParams' [Physics] ==> I3LogLikelihoodFitParams
'SecondSubtrackLlhFitFitParams' [Physics] ==> I3LogLikelihoodFitParams
'InitialPulseSeriesReco' [Physics] ==> I3Map<OMKey, std::vector<I3RecoPulse, std::allocator<I3RecoPulse> > >
'InIceRawData' [Physics] ==> I3Map<OMKey, std::vector<I3DOMLaunch, std::allocator<I3DOMLaunch> > >
'I3MCTree_MaxInIceTrackIn' [Physics] ==> I3Particle
'IceTopRawData' [Physics] ==> I3Map<OMKey, std::vector<I3DOMLaunch, std::allocator<I3DOMLaunch> > >
'IcLinefit' [Physics] ==> I3Particle
'I3EventHeader' [Physics] ==> I3EventHeader
'TrackLlhFit32Umbrella' [Physics] ==> I3Particle
'I3TriggerHierarchy' [Physics] ==> I3Tree<I3Trigger>
'IcLinefitParams' [Physics] ==> I3LineFitParams
'FirstSubtrackLlhFitFitParams' [Physics] ==> I3LogLikelihoodFitParams
'FilterMask' [Physics] ==> I3Map<std::string, I3FilterResult>
'SecondSubtrackLlhFit' [Physics] ==> I3Particle
'CascadeFirst' [Physics] ==> I3Particle
'TrackLlhFit' [Physics] ==> I3Particle
'TrackLlhFit32_LogL_dof' [Physics] ==> I3Double
'I3MCTree_MaxInIceTrackCenter' [Physics] ==> I3Particle
'TrackParaboloidFit' [Physics] ==> I3Particle
'FirstSubtrackLlhFit' [Physics] ==> I3Particle
'MMCTrackList' [Physics] ==> I3Vector<I3MMCTrack>
'wfllhFitParams' [Physics] ==> I3LogLikelihoodFitParams
'rloglFinitephpnh_rusage' [Physics] ==> I3RUsage
'TrackParaboloidFitFitParams' [Physics] ==> I3ParaboloidFitParams
'CascadeFirstParams' [Physics] ==> I3CFirstFitParams
'DrivingTime' [Physics] ==> I3Time
'TrackParaboloidFit32FitParams' [Physics] ==> I3ParaboloidFitParams
'I3MCWeightDict' [Physics] ==> I3Map<std::string, double>
'SecondSubtrackBayesianFit' [Physics] ==> I3Particle
'TrackLlhFit32FitParams' [Physics] ==> I3LogLikelihoodFitParams
'TrackBayesianFit32' [Physics] ==> I3Particle
'I3MCTree' [Physics] ==> I3Tree<I3Particle>
'MCHitSeriesMap' [Physics] ==> I3Map<OMKey, std::vector<I3MCHit, std::allocator<I3MCHit> > >
'ToI' [Physics] ==> I3Particle
'FirstSubtrackBayesianFitFitParams' [Physics] ==> I3LogLikelihoodFitParams
'I3MCTree_MaxInIceTrackOut' [Physics] ==> I3Particle
'CascadeLinefitParams' [Physics] ==> I3LineFitParams
'wfllh' [Physics] ==> I3Particle
'TrackLlhFit32UmbrellaFitParams' [Physics] ==> I3LogLikelihoodFitParams
'I3MCTree_MaxInIceTrackCoG' [Physics] ==> I3Particle
'CascadeLinefit' [Physics] ==> I3Particle
'TrackLlhFit32_LogL' [Physics] ==> I3Double
'TrackBayesianFit' [Physics] ==> I3Particle
'SecondSubtrackBayesianFitFitParams' [Physics] ==> I3LogLikelihoodFitParams
]
Looking at the mc hits,
>>> hits = frame['MCHitSeriesMap']
>>> hits
<icecube.dataclasses.I3MCHitSeriesMap object at 0x8511534>
nchannel¶
We can get ‘nchannel’:
>>> len(hits)
54
or the list of hit OMs:
>>> [e.key() for e in hits]
[OMKey(21,24), OMKey(21,31),
OMKey(29,41), OMKey(29,60), OMKey(30,47), OMKey(38,41),
OMKey(38,58), OMKey(39,13), OMKey(39,46), OMKey(40,42),
OMKey(40,43), OMKey(40,44), OMKey(40,45), OMKey(40,54),
OMKey(46,58), OMKey(47,42), OMKey(48,3), OMKey(48,38), OMKey(48,41),
OMKey(48,44), OMKey(49,15), OMKey(49,36), OMKey(49,43),
OMKey(49,44), OMKey(49,45), OMKey(49,46), OMKey(50,46),
OMKey(50,57), OMKey(56,32), OMKey(56,38), OMKey(57,5), OMKey(57,41),
OMKey(57,42), OMKey(57,43), OMKey(57,44), OMKey(58,10),
OMKey(58,42), OMKey(58,47), OMKey(59,39), OMKey(59,40), OMKey(65,2),
OMKey(65,5), OMKey(65,40), OMKey(65,44), OMKey(65,50), OMKey(66,18),
OMKey(66,52), OMKey(67,40), OMKey(72,37), OMKey(72,41),
OMKey(72,44), OMKey(73,12), OMKey(73,50), OMKey(74,10)]
nhit¶
or the total number of hits in all oms; Create a list of hits per OM:
>>> [len(e.data()) for e in hits]
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 5, 2, 1, 2, 1, 1, 1, 1, 1, 1, 5, 2, 1, 1, 2, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 2, 1, 1, 1, 1]
and take the sum:
>>> sum([len(e.data()) for e in hits])
70
one can get a list of lists of hit times:
>>> [[hit.Time for hit in e.data()] for e in hits]
[[204604.9742586276], [211559.12329093015], [217228.25888738033],
[202762.55772166929], [203625.54720526567], [223173.56777643564],
[203903.865336891], [196707.33572210337], [206972.02583769718],
[206583.63700918012], [205939.05947181155, 207110.31857003138],
[205962.32537366869, 206270.57986411842], [206487.34792608922,
214417.52757325419], [195490.18964865553], [195203.18451091682],
[211080.56898628554], [199945.66162354519], [202687.33103472061],
[205426.34651168939], [212989.69564489473], [218339.10920387553],
[223291.45251082265], [205808.24168743056], [205496.01178302418,
205496.71124786529, 205525.32546579896, 205562.98169248717,
205609.40689863981], [205537.32113750148, 207372.70384354892],
[205891.7556336984], [199383.30950566352, 211678.89014103252],
[195923.0251990071], [194945.97831539251], [196519.39973339965],
[206981.80230232462], [207692.90224580758], [205368.31951649411],
[204809.47011011466, 204848.77878282889, 204885.80844991628,
205360.19684061655, 210251.45781586654], [204799.03252794495,
205630.88520065596], [198762.71528484017], [205454.05597575437],
[206051.72036208294, 214004.1349712183], [198149.68830273568],
[213105.61270986524], [216064.29977506239], [200573.88347293844],
[205112.0535136415], [204456.76376554792], [205180.55285840927],
[202969.06534386508], [197885.66878645719], [213020.2059352111],
[221589.93938874462], [204374.1205059239, 204678.91635859985],
[204925.51566635544], [194694.86695246131], [197903.13150816062],
[222911.53660383544]]
or with a nested list comprehension get one flat list:
>>> [hit.Time for entry in hits for hit in entry.data()]
[204604.9742586276, 211559.12329093015, 217228.25888738033,
202762.55772166929, 203625.54720526567, 223173.56777643564,
203903.865336891, 196707.33572210337, 206972.02583769718,
206583.63700918012, 205939.05947181155, 207110.31857003138,
205962.32537366869, 206270.57986411842, 206487.34792608922,
214417.52757325419, 195490.18964865553, 195203.18451091682,
211080.56898628554, 199945.66162354519, 202687.33103472061,
205426.34651168939, 212989.69564489473, 218339.10920387553,
223291.45251082265, 205808.24168743056, 205496.01178302418,
205496.71124786529, 205525.32546579896, 205562.98169248717,
205609.40689863981, 205537.32113750148, 207372.70384354892,
205891.7556336984, 199383.30950566352, 211678.89014103252,
195923.0251990071, 194945.97831539251, 196519.39973339965,
206981.80230232462, 207692.90224580758, 205368.31951649411,
204809.47011011466, 204848.77878282889, 204885.80844991628,
205360.19684061655, 210251.45781586654, 204799.03252794495,
205630.88520065596, 198762.71528484017, 205454.05597575437,
206051.72036208294, 214004.1349712183, 198149.68830273568,
213105.61270986524, 216064.29977506239, 200573.88347293844,
205112.0535136415, 204456.76376554792, 205180.55285840927,
202969.06534386508, 197885.66878645719, 213020.2059352111,
221589.93938874462, 204374.1205059239, 204678.91635859985,
204925.51566635544, 194694.86695246131, 197903.13150816062,
222911.53660383544]
Getting the energy of the most energetic primary¶
In [1]: from icecube import icetray, dataclasses, dataio
In [2]: f = dataio.I3File("/v/icecube/data/combinedMC600.i3")
In [3]: frame = f.pop_physics()
In [4]: print frame
[ I3Frame (P):
'BadOMSelection' [P] ==> I3Vector<OMKey>
'CalibratedATWD' [P] ==> I3Map<OMKey, std::vector<I3Waveform, std::allocator<I3Waveform> > >
'CalibratedFADC' [P] ==> I3Map<OMKey, std::vector<I3Waveform, std::allocator<I3Waveform> > >
'DOMLaunch' [P] ==> I3Map<OMKey, std::vector<I3DOMLaunch, std::allocator<I3DOMLaunch> > >
'DOMLaunch_rc' [P] ==> I3Map<OMKey, std::vector<I3DOMLaunch, std::allocator<I3DOMLaunch> > >
'DrivingTime' [P] ==> I3Time
'DrivingTimeUCR' [P] ==> I3Time
'GlobalTriggerHierarchy' [P] ==> I3Tree<I3Trigger>
'I3-mult' [P] ==> I3Tree<I3Trigger>
'I3EventHeader' [P] ==> I3EventHeader
'I3MCTree' [P] ==> I3Tree<I3Particle>
'IIEventHead' [P] ==> I3EventHeader
'IceTopRawData' [P] ==> I3Map<OMKey, std::vector<I3DOMLaunch, std::allocator<I3DOMLaunch> > >
'MCHitSeriesMap' [P] ==> I3Map<OMKey, std::vector<I3MCHit, std::allocator<I3MCHit> > >
'PrimaryTree' [P] ==> I3Tree<I3Particle>
'TWRLaunchSeriesMap' [P] ==> I3Map<OMKey, std::vector<I3TWRLaunch, std::allocator<I3TWRLaunch> > >
'icetop_trig' [P] ==> I3Tree<I3Trigger>
'noisyMCMap' [P] ==> I3Map<OMKey, std::vector<I3MCHit, std::allocator<I3MCHit> > >
'twrSelected' [P] ==> I3Map<OMKey, std::vector<I3TWRLaunch, std::allocator<I3TWRLaunch> > >
'twrSelectedCleanedKeys' [P] ==> I3Vector<OMKey>
]
In [5]: mctree = frame['I3MCTree']
In [6]: prim = mctree.GetMostEnergeticPrimary()
In [7]: prim
Out[7]: <icecube.dataclasses.I3Particle object at 0x8372f7c>
In [8]: prim.GetEnergy()
Out[8]: 497.685
Getting an I3RecoPulseSeriesMap from an I3RecoPulseSeriesMapMask¶
In [1]: frame.keys()
Out[1]:
['CleanTriggerHierarchy_IT',
'SPEFitSingleFitParams',
'FilterMask',
'TWOfflinePulsesFR',
'SPEFitSingle',
'SPEFit2CramerRaoParams',
'MaskedOfflinePulses',
'I3EventHeader',
'SPEFit2MuE',
'OfflinePulsesHLC',
'I3DST11',
'I3SuperDST',
'TWOfflinePulsesHLC',
'SRTOfflinePulses',
'TWOfflinePulsesDC',
'RTTWOfflinePulsesFR',
'I3TriggerHierarchy',
'LineFit',
'OfflinePulses',
'TWSRTOfflinePulses',
'SPEFit2',
'LineFitParams',
'SPEFit2FitParams']
In [2]: pulse_map = dataclasses.I3RecoPulseSeriesMap.from_frame(frame,'TWSRTOfflinePulses')
In [3]: pulse_map.keys()
Out[3]:
[OMKey(4,8,0),
OMKey(4,10,0),
OMKey(10,10,0),
OMKey(11,14,0),
OMKey(11,19,0),
OMKey(11,20,0),
OMKey(11,21,0),
OMKey(11,22,0),
OMKey(11,23,0),
OMKey(18,11,0)]
Efficiently getting pulse properties with numpy¶
A series of I3RecoPulses (the elements in an I3RecoPulseSeriesMap
) can be
efficiently accessed as a numpy array using numpy.asarray(). Note the following
limitations:
The data in the array are read-only
Only the pulse time and charge are in the numpy array. For other metadata, you have to go back to the original pulses.
An example follows:
In [14]: p = d['OfflinePulses'].values()[165]
In [15]: p
Out[15]: [<icecube.dataclasses.I3RecoPulse object at 0x8676345d0>, <icecube.dataclasses.I3RecoPulse object at 0x80aabed50>, <icecube.dataclasses.I3RecoPulse object at 0x80b1a3da0>]
In [16]: numpy.asarray(p)
Out[16]:
array([[1.10534619e+004, 1.07492745e+000],
[1.52605963e+000, 1.22596182e+004],
[2.56995640e+064, 6.47956996e+170]])
In [17]: numpy.asarray(p)[0]
Out[17]: array([1.10534619e+04, 1.07492745e+00])
In [18]: print(p[0].time, p[0].charge)
11053.4619140625 1.0749274492263794
You can also access an entire I3RecoPulseSeriesMap
as a numpy array using numpy.asarray()
:
In [10]: numpy.asarray(d['OfflinePulses'])
Out[10]:
array([[4.00000000e+00, 3.20000000e+01, 0.00000000e+00, 2.36043042e+03,
5.78553498e-01, 8.33333302e+00],
[8.00000000e+00, 4.90000000e+01, 0.00000000e+00, 1.23705576e+04,
9.13645923e-01, 8.33333302e+00],
[8.00000000e+00, 5.00000000e+01, 0.00000000e+00, 1.45267861e+04,
5.65040946e-01, 8.33333302e+00],
...,
[8.60000000e+01, 2.40000000e+01, 0.00000000e+00, 2.54425366e+03,
9.20636594e-01, 8.33333302e+00],
[8.60000000e+01, 5.20000000e+01, 0.00000000e+00, 1.43962598e+04,
7.19962656e-01, 8.33333302e+00],
[8.60000000e+01, 5.90000000e+01, 0.00000000e+00, 1.30112207e+04,
1.60624707e+00, 8.33333302e+00]])
The columns are (String, OM, PMT, Time, Charge, Width), with one row per pulse. You can use this to quickly compute total charge, for example:
In [11]: numpy.asarray(d['OfflinePulses'])[:,4].sum()
Out[11]: 6854.178911853582
In addition, a helper function is provided (I3RecoPulseSeriesMap.pmt_array_offsets()
) that will return a list of indices for the rows when the array moves from one PMT to the next, which can be helpful for slicing the array representation into per-PMT blocks from Python.