EventList#
- class gammapy.data.EventList(table)[source]#
Bases:
object
Event list.
Event list data is stored as
table
(Table
) data member.The most important reconstructed event parameters are available as the following columns:
TIME
- Mission elapsed time (sec)RA
,DEC
- ICRS system position (deg)ENERGY
- Energy (usually MeV for Fermi and TeV for IACTs)
Note that
TIME
is usually sorted, but sometimes it is not. E.g. when simulating data, or processing it in certain ways. So generally any analysis code should assumeTIME
is not sorted.Other optional (columns) that are sometimes useful for high level analysis:
GLON
,GLAT
- Galactic coordinates (deg)DETX
,DETY
- Field of view coordinates (deg)
Note that when reading data for analysis you shouldn’t use those values directly, but access them via properties which create objects of the appropriate class:
- Parameters
- table
Table
Event list table
- table
Examples
>>> from gammapy.data import EventList >>> events = EventList.read("$GAMMAPY_DATA/cta-1dc/data/baseline/gps/gps_baseline_110380.fits") >>> print(events) EventList --------- Instrument : None Telescope : CTA Obs. ID : 110380 Number of events : 106217 Event rate : 59.273 1 / s Time start : 59235.5 Time stop : 59235.52074074074 Min. energy : 3.00e-02 TeV Max. energy : 1.46e+02 TeV Median energy : 1.02e-01 TeV Max. offset : 5.0 deg
Attributes Summary
ALT / AZ position computed from RA / DEC (
SkyCoord
).ALT / AZ frame (
AltAz
).ALT / AZ position from table (
SkyCoord
).Event energies (
Quantity
).Event Galactic sky coordinates (
SkyCoord
).Median position in radec
Whether observation is pointed
Dead-time fraction (float).
Live-time duration in seconds (
Quantity
).Observation time duration in seconds (
Quantity
).Observation start time (
Time
).Observation stop time (
Time
).Observatory location (
EarthLocation
).Event offset from the array pointing position (
Angle
).Event offset from the median position (
Angle
).Pointing RA / DEC sky coordinates (
SkyCoord
).Event RA / DEC sky coordinates (
SkyCoord
).Event times (
Time
).Time reference (
Time
).Methods Summary
check
([checks])Run checks.
copy
()Copy event list (
EventList
)from_stack
(event_lists, **kwargs)Stack (concatenate) list of event lists.
map_coord
(geom)Event map coordinates for a given geometry.
peek
([allsky])Quick look plots.
plot_energy
([ax])Plot counts as a function of energy.
plot_energy_offset
([ax, center])Plot counts histogram with energy and offset axes
plot_image
([ax, allsky])Quick look counts map sky plot.
plot_offset2_distribution
([ax, center])Plot offset^2 distribution of the events.
plot_time
([ax])Plots an event rate time curve.
read
(filename, **kwargs)Read from FITS file.
select_energy
(energy_range)Select events in energy band.
select_mask
(mask)Select events inside a mask (
EventList
).select_offset
(offset_band)Select events in offset band.
select_parameter
(parameter, band)Select events with respect to a specified parameter.
select_rad_max
(rad_max[, position])Select energy dependent offset
select_region
(regions[, wcs])Select events in given region.
select_row_subset
(row_specifier)Select table row subset.
select_time
(time_interval)Select events in time interval.
stack
(other)Stack with another EventList in place.
to_table_hdu
([format])Convert event list to a
BinTableHDU
write
(filename[, gti, overwrite, format])Write the event list to a FITS file.
Attributes Documentation
- galactic_median#
Median position in radec
- is_pointed_observation#
Whether observation is pointed
- observation_dead_time_fraction#
Dead-time fraction (float).
This is a keyword related to IACTs Defined as dead-time over observation time.
Dead-time is defined as the time during the observation where the detector didn’t record events: http://en.wikipedia.org/wiki/Dead_time https://ui.adsabs.harvard.edu/abs/2004APh….22..285F
The dead-time fraction is used in the live-time computation, which in turn is used in the exposure and flux computation.
- observation_live_time_duration#
Live-time duration in seconds (
Quantity
).The dead-time-corrected observation time.
In Fermi-LAT it is automatically provided in the header of the event list.
In IACTs is computed as
t_live = t_observation * (1 - f_dead)
where
f_dead
is the dead-time fraction.
- observation_time_duration#
Observation time duration in seconds (
Quantity
).This is a keyword related to IACTs The wall time, including dead-time.
- observatory_earth_location#
Observatory location (
EarthLocation
).
- time#
Event times (
Time
).Notes
Times are automatically converted to 64-bit floats. With 32-bit floats times will be incorrect by a few seconds when e.g. adding them to the reference time.
Methods Documentation
- classmethod from_stack(event_lists, **kwargs)[source]#
Stack (concatenate) list of event lists.
Calls
vstack
.- Parameters
- event_listslist
list of
EventList
to stack
- peek(allsky=False)[source]#
Quick look plots.
- Parameters
- allskybool
Whether to look at the events allsky
- plot_energy_offset(ax=None, center=None, **kwargs)[source]#
Plot counts histogram with energy and offset axes
- Parameters
- ax
Axis
Plot axis
- center
SkyCoord
Sky coord from which offset is computed
- **kwargsdict
Keyword arguments forwarded to
pcolormesh
- ax
- Returns
- ax
Axis
Plot axis
- ax
- plot_image(ax=None, allsky=False)[source]#
Quick look counts map sky plot.
- Parameters
- ax
Axes
Axes to plot on.
- allskybool,
Whether to plot on an all sky geom
- ax
- plot_offset2_distribution(ax=None, center=None, **kwargs)[source]#
Plot offset^2 distribution of the events.
The distribution shown in this plot is for this quantity:
offset = center.separation(events.radec).deg offset2 = offset ** 2
Note that this method is just for a quicklook plot.
If you want to do computations with the offset or offset^2 values, you can use the line above. As an example, here’s how to compute the 68% event containment radius using
numpy.percentile
:import numpy as np r68 = np.percentile(offset, q=68)
- Parameters
- ax
Axes
(optional) Axes
- center
astropy.coordinates.SkyCoord
Center position for the offset^2 distribution. Default is the observation pointing position.
- **kwargs
Extra keyword arguments are passed to
hist
.
- ax
- Returns
- ax
Axes
Axes
- ax
Examples
Load an example event list:
>>> from gammapy.data import EventList >>> from astropy import units as u >>> filename = "$GAMMAPY_DATA/hess-dl3-dr1/data/hess_dl3_dr1_obs_id_023523.fits.gz" >>> events = EventList.read(filename)
>>> #Plot the offset^2 distribution wrt. the observation pointing position >>> #(this is a commonly used plot to check the background spatial distribution): >>> events.plot_offset2_distribution() Plot the offset^2 distribution wrt. the Crab pulsar position (this is commonly used to check both the gamma-ray signal and the background spatial distribution):
>>> import numpy as np >>> from astropy.coordinates import SkyCoord >>> center = SkyCoord(83.63307, 22.01449, unit='deg') >>> bins = np.linspace(start=0, stop=0.3 ** 2, num=30) * u.deg ** 2 >>> events.plot_offset2_distribution(center=center, bins=bins)
Note how we passed the
bins
option ofmatplotlib.pyplot.hist
to control the histogram binning, in this case 30 bins ranging from 0 to (0.3 deg)^2.
- classmethod read(filename, **kwargs)[source]#
Read from FITS file.
Format specification: EVENTS
- Parameters
- filename
pathlib.Path
, str Filename
- filename
- select_energy(energy_range)[source]#
Select events in energy band.
- Parameters
- energy_range
Quantity
Energy range
[energy_min, energy_max)
- energy_range
- Returns
- event_list
EventList
Copy of event list with selection applied.
- event_list
Examples
>>> from astropy import units as u >>> from gammapy.data import EventList >>> filename = "$GAMMAPY_DATA/fermi_3fhl/fermi_3fhl_events_selected.fits.gz" >>> event_list = EventList.read(filename) >>> energy_range =[1, 20] * u.TeV >>> event_list = event_list.select_energy(energy_range=energy_range)
- select_mask(mask)[source]#
Select events inside a mask (
EventList
).Examples
>>> from gammapy.data import EventList >>> from gammapy.maps import WcsGeom, Map >>> geom = WcsGeom.create(skydir=(0,0), width=(4, 4), frame="galactic") >>> mask = geom.region_mask("galactic;circle(0, 0, 0.5)") >>> filename = "$GAMMAPY_DATA/cta-1dc/data/baseline/gps/gps_baseline_110380.fits" >>> events = EventList.read(filename) >>> masked_event = events.select_mask(mask) >>> len(masked_event.table) 5594
- select_offset(offset_band)[source]#
Select events in offset band.
- Parameters
- offset_band
Angle
offset band
[offset_min, offset_max)
- offset_band
- Returns
- event_list
EventList
Copy of event list with selection applied.
- event_list
Examples
>>> from gammapy.data import EventList >>> import astropy.units as u >>> filename = "$GAMMAPY_DATA/cta-1dc/data/baseline/gps/gps_baseline_110380.fits" >>> events = EventList.read(filename) >>> selected_events = events.select_offset([0.3, 0.9]*u.deg) >>> len(selected_events.table) 12688
- select_parameter(parameter, band)[source]#
Select events with respect to a specified parameter.
- Parameters
- parameterstr
Parameter used for the selection. Must be present in
self.table
.- bandtuple or
astropy.units.Quantity
Min and max value for the parameter to be selected (min <= parameter < max). If parameter is not dimensionless you have to provide a Quantity.
- Returns
- event_list
EventList
Copy of event list with selection applied.
- event_list
Examples
>>> from astropy import units as u >>> from gammapy.data import EventList >>> filename = "$GAMMAPY_DATA/fermi_3fhl/fermi_3fhl_events_selected.fits.gz" >>> event_list = EventList.read(filename) >>> zd = (0, 30) * u.deg >>> event_list = event_list.select_parameter(parameter='ZENITH_ANGLE', band=zd) >>> print(len(event_list.table)) 123944
- select_region(regions, wcs=None)[source]#
Select events in given region.
- Parameters
- regionsstr,
Region
or list ofRegion
Region or list of regions (pixel or sky regions accepted). A region can be defined as a string ind DS9 format as well. See http://ds9.si.edu/doc/ref/region.html for details.
- wcs
WCS
World coordinate system transformation
- regionsstr,
- Returns
- event_list
EventList
Copy of event list with selection applied.
- event_list
- select_row_subset(row_specifier)[source]#
Select table row subset.
- Parameters
- row_specifierslice, int, or array of ints
Specification for rows to select, passed on to
self.table[row_specifier]
.
- Returns
- event_list
EventList
New event list with table row subset selected
- event_list
Examples
>>> from gammapy.data import EventList >>> import numpy as np >>> filename = "$GAMMAPY_DATA/cta-1dc/data/baseline/gps/gps_baseline_110380.fits" >>> events = EventList.read(filename) >>> #Use a boolean mask as ``row_specifier``: >>> mask = events.table['MC_ID'] == 1 >>> events2 = events.select_row_subset(mask) >>> print(len(events2.table)) 97978 >>> #Use row index array as ``row_specifier``: >>> idx = np.where(events.table['MC_ID'] == 1)[0] >>> events2 = events.select_row_subset(idx) >>> print(len(events2.table)) 97978
- select_time(time_interval)[source]#
Select events in time interval.
- Parameters
- time_interval
astropy.time.Time
Start time (inclusive) and stop time (exclusive) for the selection.
- time_interval
- Returns
- events
EventList
Copy of event list with selection applied.
- events
- stack(other)[source]#
Stack with another EventList in place.
Calls
vstack
.- Parameters
- other
EventList
Event list to stack to self
- other
- to_table_hdu(format='gadf')[source]#
Convert event list to a
BinTableHDU
- Parameters
- format: str
Output format, currently only “gadf” is supported
- Returns
- hdu:
astropy.io.fits.BinTableHDU
EventList converted to FITS representation
- hdu:
- write(filename, gti=None, overwrite=False, format='gadf')[source]#
Write the event list to a FITS file.
If a GTI object is provided, it is saved into a second extension in the file.
- Parameters
- filename
pathlib.Path
, str Filename
- gti
GTI
Good Time Intervals object to save to the same file. Default is None.
- overwritebool
Overwrite existing file?
- formatstr, optional
FITS format convention. By default files will be written to the gamma-astro-data-formats (GADF) format.
- filename