ExcessProfileEstimator

class gammapy.estimators.ExcessProfileEstimator(regions, energy_edges=None, spectrum=None, n_sigma=1.0, n_sigma_ul=3.0, selection_optional='all')[source]

Bases: gammapy.estimators.Estimator

Estimate profile from a DataSet.

Parameters
regionslist of regions

regions to use

energy_edgesQuantity

Energy edges of the profiles to be computed.

n_sigmafloat (optional)

Number of sigma to compute errors. By default, it is 1.

n_sigma_ulfloat (optional)

Number of sigma to compute upper limit. By default, it is 3.

spectrumSpectralModel (optional)

Spectral model to compute the fluxes or brightness. Default is power-law with spectral index of 2.

selection_optionallist of str

Additional quantities to be estimated. Possible options are:

  • “errn-errp”: estimate asymmetric errors.

  • “ul”: estimate upper limits.

By default all quantities are estimated.

Examples

This example shows how to compute a counts profile for the Fermi galactic center region:

import matplotlib.pyplot as plt
from astropy import units as u
from astropy.coordinates import SkyCoord
from gammapy.data import GTI
from gammapy.estimators import ExcessProfileEstimator, ImageProfile
from gammapy.utils.regions import make_orthogonal_rectangle_sky_regions
from gammapy.datasets import Datasets

# load example data
datasets = Datasets.read("$GAMMAPY_DATA/fermi-3fhl-crab/",
    "Fermi-LAT-3FHL_datasets.yaml", "Fermi-LAT-3FHL_models.yaml")
# configuration
datasets[0].gti = GTI.create("0s", "1e7s", "2010-01-01")

# creation of the boxes and axis
start_line = SkyCoord(182.5, -5.8, unit='deg', frame='galactic')
end_line = SkyCoord(186.5, -5.8, unit='deg', frame='galactic')
boxes, axis = make_orthogonal_rectangle_sky_regions(start_line,
                                end_line,
                                datasets[0].counts.geom.wcs,
                                1.*u.deg,
                                11)

# set up profile estimator and run
prof_maker = ExcessProfileEstimator(boxes, axis)
fermi_prof = prof_maker.run(datasets[0])

# smooth and plot the data using the ImageProfile class
fermi_prof.peek()
plt.show()

ax = plt.gca()
ax.set_yscale('log')
ax = fermi_prof.plot("flux", ax=ax)

Attributes Summary

config_parameters

Config parameters

selection_optional

tag

Methods Summary

copy()

Copy estimator

get_spectrum_datasets(dataset)

Utility to make the final Datasets

get_sqrt_ts(ts, norm)

Compute sqrt(TS) value.

make_prof(sp_datasets)

Utility to make the profile in each region

run(dataset)

Make the profiles

Attributes Documentation

config_parameters

Config parameters

selection_optional
tag = 'ExcessProfileEstimator'

Methods Documentation

copy()

Copy estimator

get_spectrum_datasets(dataset)[source]

Utility to make the final Datasets

Parameters
datasetMapDataset or MapDatasetOnOff

the dataset to use for profile extraction

Returns
——–
sp_datasetsarray of SpectrumDataset

the list of SpectrumDataset computed in each box

static get_sqrt_ts(ts, norm)

Compute sqrt(TS) value.

Compute sqrt(TS) as defined by:

\[\begin{split}\sqrt{TS} = \left \{ \begin{array}{ll} -\sqrt{TS} & : \text{if} \ norm < 0 \\ \sqrt{TS} & : \text{else} \end{array} \right.\end{split}\]
Parameters
tsndarray

TS value.

normndarray

norm value

Returns
——-
sqrt_tsndarray

Sqrt(TS) value.

make_prof(sp_datasets)[source]

Utility to make the profile in each region

Parameters
sp_datasetsMapDatasets of SpectrumDataset or SpectrumDatasetOnOff

the dataset to use for profile extraction

Returns
resultslist of dictionary

the list of results (list of keys: x_min, x_ref, x_max, alpha, counts, background, excess, ts, sqrt_ts, err, errn, errp, ul, exposure, solid_angle)

run(dataset)[source]

Make the profiles

Parameters
datasetMapDataset or MapDatasetOnOff

the dataset to use for profile extraction

Returns
imageprofileImageProfile

Return an image profile class containing the result