fermipy package¶
Submodules¶
fermipy.config module¶
-
class
fermipy.config.
ConfigManager
[source]¶ Bases:
object
-
static
create
(configfile)[source]¶ Create a configuration dictionary from a yaml config file. This function will first populate the dictionary with defaults taken from pre-defined configuration files. The configuration dictionary is then updated with the user-defined configuration file. Any settings defined by the user will take precedence over the default settings.
-
static
-
class
fermipy.config.
Configurable
(config, **kwargs)[source]¶ Bases:
object
The base class provides common facilities like loading and saving configuration state.
-
config
¶ Return the configuration dictionary of this class.
-
configdir
¶
-
-
fermipy.config.
create_default_config
(defaults)[source]¶ Create a configuration dictionary from a defaults dictionary. The defaults dictionary defines valid configuration keys with default values and docstrings. Each dictionary element should be a tuple or list containing (default value,docstring,type).
fermipy.gtanalysis module¶
-
class
fermipy.gtanalysis.
GTAnalysis
(config, **kwargs)[source]¶ Bases:
fermipy.config.Configurable
,fermipy.sed.SEDGenerator
,fermipy.residmap.ResidMapGenerator
,fermipy.tsmap.TSMapGenerator
,fermipy.tsmap.TSCubeGenerator
,fermipy.sourcefind.SourceFinder
High-level analysis interface that manages a set of analysis component objects. Most of the functionality of the Fermipy package is provided through the methods of this class. The class constructor accepts a dictionary that defines the configuration for the analysis. Keyword arguments to the constructor can be used to override parameters in the configuration dictionary.
-
__delattr__
¶ x.__delattr__(‘name’) <==> del x.name
-
__format__
()¶ default object formatter
-
__getattribute__
¶ x.__getattribute__(‘name’) <==> x.name
-
__hash__
¶
-
__reduce__
()¶ helper for pickle
-
__reduce_ex__
()¶ helper for pickle
-
__repr__
¶
-
__setattr__
¶ x.__setattr__(‘name’, value) <==> x.name = value
-
__sizeof__
() → int¶ size of object in memory, in bytes
-
__str__
¶
-
add_source
(name, src_dict, free=False, init_source=True, save_source_maps=True, **kwargs)[source]¶ Add a source to the ROI model. This function may be called either before or after
setup
.Parameters:
-
add_sources_from_roi
(names, roi, free=False, **kwargs)[source]¶ Add multiple sources to the current ROI model copied from another ROI model.
Parameters:
-
bowtie
(name, fd=None, loge=None)[source]¶ Generate a spectral uncertainty band (bowtie) for the given source. This will create an uncertainty band on the differential flux as a function of energy by propagating the errors on the global fit parameters. Note that this band only reflects the uncertainty for parameters that are currently free in the model.
Parameters: - name (str) – Source name.
- fd (FluxDensity) – Flux density object. If this parameter is None then one will be created.
- loge (array-like) – Sequence of energies in log10(E/MeV) at which the flux band will be evaluated.
-
components
¶ Return the list of analysis components.
-
config
¶ Return the configuration dictionary of this class.
-
configdir
¶
-
configure
(config, **kwargs)¶
-
constrain_norms
(srcNames, cov_scale=1.0)[source]¶ Constrain the normalizations of one or more sources by adding gaussian priors with sigma equal to the parameter error times a scaling factor.
-
static
create
(infile, config=None)[source]¶ Create a new instance of GTAnalysis from an analysis output file generated with
write_roi
. By default the new instance will inherit the configuration of the saved analysis instance. The configuration may be overriden by passing a configuration file path with theconfig
argument.Parameters:
-
defaults
= {'sourcefind': {'max_iter': (3, 'Set the number of search iterations.', <type 'int'>), 'min_separation': (1.0, 'Set the minimum separation in deg for sources added in each iteration.', <type 'float'>), 'tsmap_fitter': ('tsmap', 'Set the method for generating the TS map.', <type 'str'>), 'sqrt_ts_threshold': (5.0, 'Set the threshold on sqrt(TS).', <type 'float'>), 'model': (None, 'Set the source model dictionary. By default the test source will be a PointSource with an Index 2 power-law specturm.', <type 'dict'>), 'sources_per_iter': (3, '', <type 'int'>)}, 'roiopt': {'npred_frac': (0.95, '', <type 'float'>), 'shape_ts_threshold': (25.0, 'Threshold on source TS used for determining the sources that will be fit in the third optimization step.', <type 'float'>), 'npred_threshold': (1.0, '', <type 'float'>), 'skip': (None, 'List of str source names to skip while optimizing.', <type 'list'>), 'max_free_sources': (5, 'Maximum number of sources that will be fit simultaneously in the first optimization step.', <type 'int'>)}, 'selection': {'radius': (None, 'Radius of data selection. If none this will be automatically set from the ROI size.', <type 'float'>), 'tmin': (None, 'Minimum time (MET).', <type 'int'>), 'target': (None, 'Choose an object on which to center the ROI. This option takes precendence over ra/dec or glon/glat.', <type 'str'>), 'glon': (None, '', <type 'float'>), 'emin': (None, 'Minimum Energy (MeV)', <type 'float'>), 'emax': (None, 'Maximum Energy (MeV)', <type 'float'>), 'tmax': (None, 'Maximum time (MET).', <type 'int'>), 'glat': (None, '', <type 'float'>), 'filter': (None, 'Filter string for ``gtmktime`` selection.', <type 'str'>), 'logemax': (None, 'Maximum Energy (log10(MeV))', <type 'float'>), 'ra': (None, '', <type 'float'>), 'evtype': (None, 'Event type selection.', <type 'int'>), 'evclass': (None, 'Event class selection.', <type 'int'>), 'zmax': (None, 'Maximum zenith angle.', <type 'float'>), 'logemin': (None, 'Minimum Energy (log10(MeV))', <type 'float'>), 'dec': (None, '', <type 'float'>), 'roicut': ('no', '', <type 'str'>), 'convtype': (None, 'Conversion type selection.', <type 'int'>)}, 'logging': {'verbosity': (3, '', <type 'int'>), 'chatter': (3, 'Set the chatter parameter of the STs.', <type 'int'>)}, 'tsmap': {'multithread': (False, '', <type 'bool'>), 'model': (None, 'Dictionary defining the properties of the test source.', <type 'dict'>), 'loge_bounds': (None, 'Lower and upper energy bounds in log10(E/MeV). By default the calculation will be performed over the full analysis energy range.', <type 'list'>), 'max_kernel_radius': (3.0, '', <type 'float'>)}, 'mc': {'seed': (None, '', <type 'int'>)}, 'components': (None, '', <type 'list'>), 'localize': {'dtheta_max': (0.3, 'Half-width of the search region in degrees used for the first pass of the localization search.', <type 'float'>), 'nstep': (5, 'Number of steps along each spatial dimension in the refined likelihood scan.', <type 'int'>), 'fix_background': (True, 'Fix background parameters when fitting the source flux in each energy bin.', <type 'bool'>), 'update': (False, 'Update the source model with the best-fit position.', <type 'bool'>)}, 'binning': {'projtype': ('WCS', 'Projection mode (WCS or HPX).', <type 'str'>), 'binsperdec': (8, 'Number of energy bins per decade.', <type 'float'>), 'enumbins': (None, 'Number of energy bins. If none this will be inferred from energy range and ``binsperdec`` parameter.', <type 'int'>), 'roiwidth': (10.0, 'Width of the ROI in degrees. The number of pixels in each spatial dimension will be set from ``roiwidth`` / ``binsz`` (rounded up).', <type 'float'>), 'hpx_ebin': (True, 'Include energy binning', <type 'bool'>), 'binsz': (0.1, 'Spatial bin size in degrees.', <type 'float'>), 'npix': (None, 'Number of pixels. If none then this will be set from ``roiwidth`` and ``binsz``.', <type 'int'>), 'hpx_order': (10, 'Order of the map (int between 0 and 12, included)', <type 'int'>), 'proj': ('AIT', 'Spatial projection for WCS mode.', <type 'str'>), 'coordsys': ('CEL', 'Coordinate system of the spatial projection (CEL or GAL).', <type 'str'>), 'hpx_ordering_scheme': ('RING', 'HEALPix Ordering Scheme', <type 'str'>)}, 'extension': {'width': (None, 'Parameter vector for scan over spatial extent. If none then the parameter vector will be set from ``width_min``, ``width_max``, and ``width_nstep``.', <type 'str'>), 'fix_background': (False, 'Fix any background parameters that are currently free in the model when performing the likelihood scan over extension.', <type 'bool'>), 'width_max': (1.0, 'Maximum value in degrees for the likelihood scan over spatial extent.', <type 'float'>), 'sqrt_ts_threshold': (None, 'Threshold on sqrt(TS_ext) that will be applied when ``update`` is True. If None then nothreshold is applied.', <type 'float'>), 'width_min': (0.01, 'Minimum value in degrees for the likelihood scan over spatial extent.', <type 'float'>), 'spatial_model': ('RadialGaussian', 'Spatial model use for extension test.', <type 'str'>), 'update': (False, 'Update the source model with the best-fit spatial extension.', <type 'bool'>), 'width_nstep': (21, 'Number of steps for the spatial likelihood scan.', <type 'int'>)}, 'sed': {'use_local_index': (False, 'Use a power-law approximation to the shape of the global spectrum in each bin. If this is false then a constant index set to `bin_index` will be used.', <type 'bool'>), 'bin_index': (2.0, 'Spectral index that will be use when fitting the energy distribution within an energy bin.', <type 'float'>), 'cov_scale': (3.0, 'Scale factor that sets the strength of the prior on nuisance parameters when ``fix_background``=True. Setting this to None disables the prior.', <type 'float'>), 'fix_background': (True, 'Fix background normalization parameters when fitting the source flux in each energy bin. If True background normalizations will be profiled with a prior on their value with strength set by ``cov_scale``.', <type 'bool'>), 'ul_confidence': (0.95, 'Confidence level for upper limit calculation.', <type 'float'>)}, 'fileio': {'usescratch': (False, 'Run analysis in a temporary working directory under ``scratchdir``.', <type 'bool'>), 'scratchdir': ('/scratch', 'Path to the scratch directory. If ``usescratch`` is True then a temporary working directory will be created under this directory.', <type 'str'>), 'savefits': (True, 'Save intermediate FITS files.', <type 'bool'>), 'workdir': (None, 'Path to the working directory.', <type 'str'>), 'outdir_regex': (['\\.fits$|\\.fit$|\\.xml$|\\.npy$|\\.png$|\\.pdf$|\\.yaml$'], 'Stage files to the output directory that match at least one of the regular expressions in this list. This option only takes effect when ``usescratch`` is True.', <type 'list'>), 'workdir_regex': (['\\.fits$|\\.fit$|\\.xml$|\\.npy$'], 'Stage files to the working directory that match at least one of the regular expressions in this list. This option only takes effect when ``usescratch`` is True.', <type 'list'>), 'logfile': (None, 'Path to log file. If None then log will be written to fermipy.log.', <type 'str'>), 'outdir': (None, 'Path of the output directory. If none this will default to the directory containing the configuration file.', <type 'str'>)}, 'gtlike': {'irfs': (None, 'Set the IRF string.', <type 'str'>), 'minbinsz': (0.05, 'Set the minimum bin size used for resampling diffuse maps.', <type 'float'>), 'bexpmap': (None, '', <type 'str'>), 'edisp': (True, 'Enable the correction for energy dispersion.', <type 'bool'>), 'srcmap': (None, '', <type 'str'>), 'resample': (True, '', <type 'bool'>), 'llscan_npts': (20, 'Number of evaluation points to use when performing a likelihood scan.', <type 'int'>), 'convolve': (True, '', <type 'bool'>), 'rfactor': (2, '', <type 'int'>), 'edisp_disable': (None, 'Provide a list of sources for which the edisp correction should be disabled.', <type 'list'>)}, 'residmap': {'model': (None, 'Dictionary defining the properties of the test source. By default the test source will be a PointSource with an Index 2 power-law specturm.', <type 'dict'>), 'loge_bounds': (None, 'Lower and upper energy bounds in log10(E/MeV). By default the calculation will be performed over the full analysis energy range.', <type 'list'>)}, 'optimizer': {'retries': (3, 'Set the number of times to retry the fit when the fit quality is less than ``min_fit_quality``.', <type 'int'>), 'optimizer': ('MINUIT', 'Set the optimization algorithm to use when maximizing the likelihood function.', <type 'str'>), 'verbosity': (0, '', <type 'int'>), 'max_iter': (100, 'Maximum number of iterations for the Newtons method fitter.', <type 'int'>), 'min_fit_quality': (2, 'Set the minimum fit quality.', <type 'int'>), 'tol': (0.001, 'Set the optimizer tolerance.', <type 'float'>), 'init_lambda': (0.0001, 'Initial value of damping parameter for step size calculation when using the NEWTON fitter. A value of zero disables damping.', <type 'float'>)}, 'model': {'catalogs': (None, '', <type 'list'>), 'limbdiff': (None, '', <type 'list'>), 'src_radius_roi': (None, 'Half-width of ``src_roiwidth`` selection. This parameter can be used in lieu of ``src_roiwidth``.', <type 'float'>), 'extdir': (None, 'Set a directory that will be searched for extended source FITS templates. Template files in this directory will take precendence over catalog source templates with the same name.', <type 'str'>), 'sources': (None, '', <type 'list'>), 'assoc_xmatch_columns': (['3FGL_Name'], 'Choose a set of association columns on which to cross-match catalogs.', <type 'list'>), 'diffuse': (None, '', <type 'list'>), 'src_roiwidth': (None, 'Width of square selection cut for inclusion of catalog sources in the model. Includes sources within a square region with side ``src_roiwidth`` centered on the ROI. If this parameter is none then no selection is applied. This selection will be ORed with the ``src_radius`` selection.', <type 'float'>), 'isodiff': (None, 'Set the isotropic template.', <type 'list'>), 'merge_sources': (True, 'Merge properties of sources that appear in multiple source catalogs. If merge_sources=false then subsequent sources with the same name will be ignored.', <type 'bool'>), 'extract_diffuse': (False, 'Extract a copy of all mapcube components centered on the ROI.', <type 'bool'>), 'src_radius': (None, 'Radius of circular selection cut for inclusion of catalog sources in the model. Includes sources within a circle of this radius centered on the ROI. If this parameter is none then no selection is applied. This selection will be ORed with the ``src_roiwidth`` selection.', <type 'float'>), 'galdiff': (None, 'Set the galactic IEM mapcube.', <type 'list'>)}, 'data': {'evfile': (None, 'Path to FT1 file or list of FT1 files.', <type 'str'>), 'cacheft1': (True, 'Cache FT1 files when performing binned analysis. If false then only the counts cube is retained.', <type 'bool'>), 'scfile': (None, 'Path to FT2 (spacecraft) file.', <type 'str'>), 'ltcube': (None, 'Path to livetime cube. If none a livetime cube will be generated with ``gtmktime``.', <type 'str'>)}, 'plotting': {'catalogs': (None, '', <type 'list'>), 'format': ('png', '', <type 'str'>), 'loge_bounds': (None, '', <type 'list'>), 'graticule_radii': (None, 'Define a list of radii at which circular graticules will be drawn.', <type 'list'>), 'cmap': ('ds9_b', 'Set the colormap for 2D plots.', <type 'str'>), 'label_ts_threshold': (0.0, 'TS threshold for labeling sources in sky maps. If None then no sources will be labeled.', <type 'float'>)}, 'tscube': {'do_sed': (True, 'Compute the energy bin-by-bin fits', <type 'bool'>), 'remake_test_source': (False, 'If true, recomputes the test source image (otherwise just shifts it)', <type 'bool'>), 'st_scan_level': (0, 'Level to which to do ST-based fitting (for testing)', <type 'int'>), 'cov_scale': (-1.0, 'Scale factor to apply to broadband fitting cov. matrix in bin-by-bin fits ( < 0 -> fixed ) ', <type 'float'>), 'max_iter': (30, 'Maximum number of iterations for the Newtons method fitter.', <type 'int'>), 'nnorm': (10, 'Number of points in the likelihood v. normalization scan', <type 'int'>), 'norm_sigma': (5.0, 'Number of sigma to use for the scan range ', <type 'float'>), 'tol_type': (0, 'Absoulte (0) or relative (1) criteria for convergence.', <type 'int'>), 'cov_scale_bb': (-1.0, 'Scale factor to apply to global fitting cov. matrix in broadband fits. ( < 0 -> no prior ) ', <type 'float'>), 'tol': (0.001, 'Critetia for fit convergence (estimated vertical distance to min < tol )', <type 'float'>), 'model': (None, 'Dictionary defining the properties of the test source. By default the test source will be a PointSource with an Index 2 power-law specturm.', <type 'dict'>), 'init_lambda': (0, 'Initial value of damping parameter for newton step size calculation.', <type 'float'>)}}¶
-
delete_source
(name, save_template=True, delete_source_map=False, build_fixed_wts=True, **kwargs)[source]¶ Delete a source from the ROI model.
Parameters: Returns: src – The deleted source object.
Return type:
-
delete_sources
(cuts=None, distance=None, skydir=None, minmax_ts=None, minmax_npred=None, square=False, exclude_diffuse=True)[source]¶ Delete sources in the ROI model satisfying the given selection criteria.
- cuts : dict
- Dictionary of [min,max] selections on source properties.
- distance : float
- Cut on angular distance from
skydir
. If None then no selection will be applied. - skydir :
SkyCoord
- Reference sky coordinate for
distance
selection. If None then the distance selection will be applied with respect to the ROI center. - minmax_ts : list
- Free sources that have TS in the range [min,max]. If either min or max are None then only a lower (upper) bound will be applied. If this parameter is none no selection will be applied.
- minmax_npred : list
- Free sources that have npred in the range [min,max]. If either min or max are None then only a lower (upper) bound will be applied. If this parameter is none no selection will be applied.
- square : bool
- Switch between applying a circular or square (ROI-like) selection on the maximum projected distance from the ROI center.
Returns: srcs – A list of Model
objects.Return type: list
-
energies
¶ Return the energy bin edges in MeV.
-
enumbins
¶ Return the number of energy bins.
-
extension
(name, **kwargs)[source]¶ Test this source for spatial extension with the likelihood ratio method (TS_ext). This method will substitute an extended spatial model for the given source and perform a one-dimensional scan of the spatial extension parameter over the range specified with the width parameters. The 1-D profile likelihood is then used to compute the best-fit value, upper limit, and TS for extension. Any background parameters that are free will also be simultaneously profiled in the likelihood scan.
Parameters: - name (str) – Source name.
- spatial_model (str) –
Spatial model that will be used to test the source extension. The spatial scale parameter of the respective model will be set such that the 68% containment radius of the model is equal to the width parameter. The following spatial models are supported:
- RadialDisk : Azimuthally symmetric 2D disk.
- RadialGaussian : Azimuthally symmetric 2D gaussian.
- width_min (float) – Minimum value in degrees for the spatial extension scan.
- width_max (float) – Maximum value in degrees for the spatial extension scan.
- width_nstep (int) – Number of scan points between width_min and width_max. Scan points will be spaced evenly on a logarithmic scale between log(width_min) and log(width_max).
- width (array-like) – Sequence of values in degrees for the spatial extension scan. If this argument is None then the scan points will be determined from width_min/width_max/width_nstep.
- fix_background (bool) – Fix all background sources when performing the extension fit.
- update (bool) – Update this source with the best-fit model for spatial
extension if TS_ext >
tsext_threshold
. - sqrt_ts_threshold (float) – Threshold on sqrt(TS_ext) that will be applied when
update
is true. If None then no threshold will be applied. - optimizer (dict) – Dictionary that overrides the default optimizer settings.
Returns: extension – Dictionary containing results of the extension analysis. The same dictionary is also saved to the dictionary of this source under ‘extension’.
Return type:
-
find_sources
(prefix='', **kwargs)¶ An iterative source-finding algorithm.
Parameters: - model (dict) – Dictionary defining the properties of the test source. This is the model that will be used for generating TS maps.
- sqrt_ts_threshold (float) – Source threshold in sqrt(TS). Only peaks with sqrt(TS) exceeding this threshold will be used as seeds for new sources.
- min_separation (float) – Minimum separation in degrees of sources detected in each iteration. The source finder will look for the maximum peak in the TS map within a circular region of this radius.
- max_iter (int) – Maximum number of source finding iterations. The source finder will continue adding sources until no additional peaks are found or the number of iterations exceeds this number.
- sources_per_iter (int) – Maximum number of sources that will be added in each iteration. If the number of detected peaks in a given iteration is larger than this number, only the N peaks with the largest TS will be used as seeds for the current iteration.
- tsmap_fitter (str) –
Set the method used internally for generating TS maps. Valid options:
- tsmap
- tscube
- tsmap (dict) – Keyword arguments dictionary for tsmap method.
- tscube (dict) – Keyword arguments dictionary for tscube method.
Returns: - peaks (list) – List of peak objects.
- sources (list) – List of source objects.
-
fit
(update=True, **kwargs)[source]¶ Run the likelihood optimization. This will execute a fit of all parameters that are currently free in the model and update the charateristics of the corresponding model components (TS, npred, etc.). The fit will be repeated N times (set with the
retries
parameter) until a fit quality greater than or equal tomin_fit_quality
and a fit status code of 0 is obtained. If the fit does not succeed after N retries then all parameter values will be reverted to their state prior to the execution of the fit.Parameters: - update (bool) – Update the model dictionary for all sources with free parameters.
- tol (float) – Set the optimizer tolerance.
- verbosity (int) – Set the optimizer output level.
- optimizer (str) – Set the likelihood optimizer (e.g. MINUIT or NEWMINUIT).
- retries (int) – Set the number of times to rerun the fit when the fit quality is < 3.
- min_fit_quality (int) – Set the minimum fit quality. If the fit quality is smaller than this value then all model parameters will be restored to their values prior to the fit.
- reoptimize (bool) – Refit background sources when updating source properties (TS and likelihood profiles).
Returns: fit – Dictionary containing diagnostic information from the fit (fit quality, parameter covariances, etc.).
Return type:
-
free_source
(name, free=True, pars=None, **kwargs)[source]¶ Free/Fix parameters of a source.
Parameters: - name (str) – Source name.
- free (bool) – Choose whether to free (free=True) or fix (free=False) source parameters.
- pars (list) – Set a list of parameters to be freed/fixed for this source. If none then all source parameters will be freed/fixed with the exception of those defined in the skip_pars list.
-
free_sources
(free=True, pars=None, cuts=None, distance=None, skydir=None, minmax_ts=None, minmax_npred=None, square=False, exclude_diffuse=False, **kwargs)[source]¶ Free or fix sources in the ROI model satisfying the given selection. When multiple selections are defined, the selected sources will be those satisfying the logical AND of all selections (e.g. distance < X && minmax_ts[0] < ts < minmax_ts[1] && ...).
Parameters: - free (bool) – Choose whether to free (free=True) or fix (free=False) source parameters.
- pars (list) – Set a list of parameters to be freed/fixed for each source. If none then all source parameters will be freed/fixed. If pars=’norm’ then only normalization parameters will be freed.
- cuts (dict) – Dictionary of [min,max] selections on source properties.
- distance (float) – Cut on angular distance from
skydir
. If None then no selection will be applied. - skydir (
SkyCoord
) – Reference sky coordinate fordistance
selection. If None then the distance selection will be applied with respect to the ROI center. - minmax_ts (list) – Free sources that have TS in the range [min,max]. If either min or max are None then only a lower (upper) bound will be applied. If this parameter is none no selection will be applied.
- minmax_npred (list) – Free sources that have npred in the range [min,max]. If either min or max are None then only a lower (upper) bound will be applied. If this parameter is none no selection will be applied.
- square (bool) – Switch between applying a circular or square (ROI-like) selection on the maximum projected distance from the ROI center.
- exclude_diffuse (bool) – Exclude diffuse sources.
Returns: srcs – A list of
Model
objects.Return type:
-
generate_model
(model_name=None)[source]¶ Generate model maps for all components. model_name should be a unique identifier for the model. If model_name is None then the model maps will be generated using the current parameters of the ROI.
-
get_config
()¶ Return a default configuration dictionary for this class.
-
get_source_dfde
(name)[source]¶ Return differential flux distribution of a source. For sources with FileFunction spectral type this returns the internal differential flux array.
Returns:
-
get_source_name
(name)[source]¶ Return the name of a source as it is defined in the pyLikelihood model object.
-
get_sources
(cuts=None, distance=None, skydir=None, minmax_ts=None, minmax_npred=None, square=False)[source]¶ Retrieve list of sources in the ROI satisfying the given selections.
Returns: srcs – A list of Model
objects.Return type: list
-
get_src_model
(name, paramsonly=False, reoptimize=False, npts=None, **kwargs)[source]¶ Compose a dictionary for a source with the current best-fit parameters.
Parameters: Returns: src_dict
Return type:
-
like
¶ Return the global likelihood object.
-
load_roi
(infile, reload_sources=False)[source]¶ This function reloads the analysis state from a previously saved instance generated with
write_roi
.Parameters:
-
load_xml
(xmlfile)[source]¶ Load model definition from XML.
Parameters: xmlfile (str) – Name of the input XML file.
-
localize
(name, **kwargs)¶ Find the best-fit position of a source. Localization is performed in two steps. First a TS map is computed centered on the source with half-width set by
dtheta_max
. A fit is then performed to the maximum TS peak in this map. The source position is then further refined by scanning the likelihood in the vicinity of the peak found in the first step. The size of the scan region is set to encompass the 99% positional uncertainty contour as determined from the peak fit.Parameters: - name (str) – Source name.
- dtheta_max (float) – Maximum offset in RA/DEC in deg from the nominal source position that will be used to define the boundaries of the TS map search region.
- nstep (int) – Number of steps in longitude/latitude that will be taken when refining the source position. The bounds of the scan range are set to the 99% positional uncertainty as determined from the TS map peak fit. The total number of sampling points will be nstep**2.
- fix_background (bool) – Fix background parameters when fitting the source position.
- update (bool) – Update the model for this source with the best-fit position. If newname=None this will overwrite the existing source map of this source with one corresponding to its new location.
- newname (str) – Name that will be assigned to the relocalized source when update=True. If newname is None then the existing source name will be used.
- optimizer (dict) – Dictionary that overrides the default optimizer settings.
Returns: localize – Dictionary containing results of the localization analysis. This dictionary is also saved to the dictionary of this source in ‘localize’.
Return type:
-
log_energies
¶ Return the energy bin edges in log10(E/MeV).
-
loge_bounds
¶ Current analysis energy bounds in log10(E/MeV).
-
make_plots
(prefix, mcube_map=None, **kwargs)[source]¶ Make diagnostic plots using the current ROI model.
-
model_counts_map
(name=None, exclude=None)[source]¶ Return the model counts map for a single source, a list of sources, or for the sum of all sources in the ROI. The exclude parameter can be used to exclude one or more components when generating the model map.
Parameters: - name (str or list of str) – Parameter controlling the set of sources for which the model counts map will be calculated. If name=None the model map will be generated for all sources in the ROI.
- exclude (str or list of str) – List of sources that will be excluded when calculating the model map.
Returns: map
Return type:
-
model_counts_spectrum
(name, logemin=None, logemax=None, summed=False)[source]¶ Return the predicted number of model counts versus energy for a given source and energy range. If summed=True return the counts spectrum summed over all components otherwise return a list of model spectra.
-
npix
¶ Return the number of energy bins.
-
optimize
(**kwargs)[source]¶ Iteratively optimize the ROI model. The optimization is performed in three sequential steps:
- Free the normalization of the N largest components (as
determined from NPred) that contain a fraction
npred_frac
of the total predicted counts in the model and perform a simultaneous fit of the normalization parameters of these components. - Individually fit the normalizations of all sources that were
not included in the first step in order of their npred
values. Skip any sources that have NPred <
npred_threshold
. - Individually fit the shape and normalization parameters of
all sources with TS >
shape_ts_threshold
where TS is determined from the first two steps of the ROI optimization.
To ensure that the model is fully optimized this method can be run multiple times.
Parameters: - npred_frac (float) – Threshold on the fractional number of counts in the N largest components in the ROI. This parameter determines the set of sources that are fit in the first optimization step.
- npred_threshold (float) – Threshold on the minimum number of counts of individual sources. This parameter determines the sources that are fit in the second optimization step.
- shape_ts_threshold (float) – Threshold on source TS used for determining the sources that will be fit in the third optimization step.
- max_free_sources (int) – Maximum number of sources that will be fit simultaneously in the first optimization step.
- skip (list) – List of str source names to skip while optimizing.
- optimizer (dict) – Dictionary that overrides the default optimizer settings.
- Free the normalization of the N largest components (as
determined from NPred) that contain a fraction
-
outdir
¶ Return the analysis output directory.
-
print_config
(logger, loglevel=None)¶
-
print_params
(allpars=False, loglevel=20)[source]¶ Print information about the model parameters (values, errors, bounds, scale).
-
print_roi
(loglevel=20)[source]¶ Print information about the spectral and spatial properties of the ROI (sources, diffuse components).
-
profile
(name, parName, logemin=None, logemax=None, reoptimize=False, xvals=None, npts=None, savestate=True, **kwargs)[source]¶ Profile the likelihood for the given source and parameter.
Parameters: Returns: lnlprofile – Dictionary containing results of likelihood scan.
Return type:
-
profile_norm
(name, logemin=None, logemax=None, reoptimize=False, xvals=None, npts=None, fix_shape=True, savestate=True, **kwargs)[source]¶ Profile the normalization of a source.
Parameters:
-
projtype
¶ Return the type of projection to use
-
reload_source
(name, init_source=True)[source]¶ Delete and reload a source in the model. This will refresh the spatial model of this source to the one defined in the XML model.
-
residmap
(prefix='', **kwargs)¶ Generate 2-D spatial residual maps using the current ROI model and the convolution kernel defined with the
model
argument.Parameters: - prefix (str) – String that will be prefixed to the output residual map files.
- model (dict) – Dictionary defining the properties of the convolution kernel.
- exclude (str or list of str) – Source or sources that will be removed from the model when computing the residual map.
- loge_bounds (list) – Restrict the analysis to an energy range (emin,emax) in log10(E/MeV) that is a subset of the analysis energy range. By default the full analysis energy range will be used. If either emin/emax are None then only an upper/lower bound on the energy range wil be applied.
- make_plots (bool) – Write image files.
- write_fits (bool) – Write FITS files.
Returns: maps – A dictionary containing the
Map
objects for the residual significance and amplitude.Return type:
-
roi
¶ Return the ROI object.
-
sed
(name, **kwargs)¶ Generate a spectral energy distribution (SED) for a source. This function will fit the normalization of the source in each energy bin. By default the SED will be generated with the analysis energy bins but a custom binning can be defined with the
loge_bins
parameter.Parameters: - name (str) – Source name.
- prefix (str) – Optional string that will be prepended to all output files (FITS and rendered images).
- loge_bins (
ndarray
) – Sequence of energies in log10(E/MeV) defining the edges of the energy bins. If this argument is None then the analysis energy bins will be used. The energies in this sequence must align with the bin edges of the underyling analysis instance. - bin_index (float) – Spectral index that will be use when fitting the energy distribution within an energy bin.
- use_local_index (bool) – Use a power-law approximation to the shape of the global
spectrum in each bin. If this is false then a constant
index set to
bin_index
will be used. - fix_background (bool) – Fix background components when fitting the flux normalization in each energy bin. If fix_background=False then all background parameters that are currently free in the fit will be profiled. By default fix_background=True.
- ul_confidence (float) – Set the confidence level that will be used for the calculation of flux upper limits in each energy bin.
- cov_scale (float) – Scaling factor that will be applied when setting the gaussian prior on the normalization of free background sources. If this parameter is None then no gaussian prior will be applied.
- write_fits (bool) – Write a FITS file containing the SED analysis results.
- write_npy (bool) – Write a numpy file with the contents of the output dictionary.
- optimizer (dict) – Dictionary that overrides the default optimizer settings.
Returns: sed – Dictionary containing output of the SED analysis. This dictionary is also saved to the ‘sed’ dictionary of the
Source
instance.Return type:
-
set_edisp_flag
(name, flag=True)[source]¶ Enable or disable the energy dispersion correction for the given source.
-
set_energy_range
(logemin, logemax)[source]¶ Set the energy bounds of the analysis. This restricts the evaluation of the likelihood to the data that falls in this range. Input values will be rounded to the closest bin edge value. If either argument is None then the lower or upper bound of the analysis instance will be used.
Parameters: Returns: eminmax – Minimum and maximum energy in log10(E/MeV).
Return type:
-
set_parameter
(name, par, value, true_value=True, scale=None, bounds=None, update_source=True)[source]¶ Update the value of a parameter. Parameter bounds will automatically be adjusted to encompass the new parameter value.
Parameters: - name (str) – Source name.
- par (str) – Parameter name.
- value (float) – Parameter value. By default this argument should be the unscaled (True) parameter value.
- scale (float) – Parameter scale (optional). Value argument is interpreted with respect to the scale parameter if it is provided.
- update_source (bool) – Update the source dictionary for the object.
-
set_parameter_scale
(name, par, scale)[source]¶ Update the scale of a parameter while keeping its value constant.
-
set_source_dfde
(name, dfde, update_source=True)[source]¶ Set the differential flux distribution of a source with the FileFunction spectral type.
Parameters:
-
set_source_spectrum
(name, spectrum_type='PowerLaw', spectrum_pars=None, update_source=True)[source]¶ Set the spectral model of a source. This function can be used to change the spectral type of a source or modify its spectral parameters. If called with spectrum_type=’FileFunction’ and spectrum_pars=None, the source spectrum will be replaced with a FileFunction with the same differential flux distribution as the original spectrum.
Parameters:
-
setup
(init_sources=True, overwrite=False)[source]¶ Run pre-processing for each analysis component and construct a joint likelihood object. This function performs the following tasks: data selection (gtselect, gtmktime), data binning (gtbin), and model generation (gtexpcube2,gtsrcmaps).
Parameters: - init_sources (bool) – Choose whether to compute properties (flux, TS, etc.) for individual sources.
- overwrite (bool) – Run all pre-processing steps even if the output file of that step is present in the working directory. By default this function will skip any steps for which the output file already exists.
-
simulate_roi
(name=None, randomize=True, restore=False)[source]¶ Generate a simulation of the ROI using the current best-fit model and replace the data counts cube with this simulation. The simulation is created by generating an array of Poisson random numbers with expectation values drawn from the model cube of the binned analysis instance. This function will update the counts cube both in memory and in the source map file. The counts cube can be restored to its original state by calling this method with
restore
= True.Parameters:
-
simulate_source
(src_dict=None)[source]¶ Inject simulated source counts into the data.
Parameters: src_dict (dict) – Dictionary defining the spatial and spectral properties of the source that will be injected.
-
tscube
(prefix='', **kwargs)¶ Generate a spatial TS map for a source component with properties defined by the
model
argument. This method uses thegttscube
ST application for source fitting and will simultaneously fit the test source normalization as well as the normalizations of any background components that are currently free. The output of this method is a dictionary containingMap
objects with the TS and amplitude of the best-fit test source. By default this method will also save maps to FITS files and render them as image files.Parameters: - prefix (str) – Optional string that will be prepended to all output files (FITS and rendered images).
- model (dict) – Dictionary defining the properties of the test source.
- do_sed (bool) – Compute the energy bin-by-bin fits.
- nnorm (int) – Number of points in the likelihood v. normalization scan.
- norm_sigma (float) – Number of sigma to use for the scan range.
- tol (float) – Critetia for fit convergence (estimated vertical distance to min < tol ).
- tol_type (int) – Absoulte (0) or relative (1) criteria for convergence.
- max_iter (int) – Maximum number of iterations for the Newton’s method fitter
- remake_test_source (bool) – If true, recomputes the test source image (otherwise just shifts it)
- st_scan_level (int) –
- make_plots (bool) – Write image files.
- write_fits (bool) – Write a FITS file with the results of the analysis.
Returns: maps – A dictionary containing the
Map
objects for TS and source amplitude.Return type:
-
tsmap
(prefix='', **kwargs)¶ Generate a spatial TS map for a source component with properties defined by the
model
argument. The TS map will have the same geometry as the ROI. The output of this method is a dictionary containingMap
objects with the TS and amplitude of the best-fit test source. By default this method will also save maps to FITS files and render them as image files.This method uses a simplified likelihood fitting implementation that only fits for the normalization of the test source. Before running this method it is recommended to first optimize the ROI model (e.g. by running
optimize()
).Parameters: - prefix (str) – Optional string that will be prepended to all output files (FITS and rendered images).
- model (dict) – Dictionary defining the properties of the test source.
- exclude (str or list of str) – Source or sources that will be removed from the model when computing the TS map.
- loge_bounds (list) – Restrict the analysis to an energy range (emin,emax) in log10(E/MeV) that is a subset of the analysis energy range. By default the full analysis energy range will be used. If either emin/emax are None then only an upper/lower bound on the energy range wil be applied.
- max_kernel_radius (float) – Set the maximum radius of the test source kernel. Using a smaller value will speed up the TS calculation at the loss of accuracy. The default value is 3 degrees.
- make_plots (bool) – Write image files.
- write_fits (bool) – Write a FITS file.
- write_npy (bool) – Write a numpy file.
Returns: maps – A dictionary containing the
Map
objects for TS and source amplitude.Return type:
-
update_source
(name, paramsonly=False, reoptimize=False, **kwargs)[source]¶ Update the dictionary for this source.
Parameters:
-
workdir
¶ Return the analysis working directory.
-
write_config
(outfile)¶ Write the configuration dictionary to an output file.
-
write_model_map
(model_name, name=None)[source]¶ Save the counts model map to a FITS file.
Parameters:
-
write_roi
(outfile=None, save_model_map=False, fmt='npy', **kwargs)[source]¶ Write current state of the analysis to a file. This method writes an XML model definition, a ROI dictionary, and a FITS source catalog file. A previously saved analysis state can be reloaded from the ROI dictionary file with the
load_roi
method.Parameters: - outfile (str) – String prefix of the output files. The extension of this string will be stripped when generating the XML, YAML and npy filenames.
- make_plots (bool) – Generate diagnostic plots.
- save_model_map (bool) – Save the current counts model to a FITS file.
- fmt (str) – Set the output file format (yaml or npy).
-
fermipy.logger module¶
-
class
fermipy.logger.
Logger
[source]¶ Bases:
object
This class provides helper functions which facilitate creating instances of the built-in logger class.
fermipy.roi_model module¶
-
class
fermipy.roi_model.
IsoSource
(name, data)[source]¶ Bases:
fermipy.roi_model.Model
-
diffuse
¶
-
filefunction
¶
-
-
class
fermipy.roi_model.
MapCubeSource
(name, data)[source]¶ Bases:
fermipy.roi_model.Model
-
diffuse
¶
-
mapcube
¶
-
-
class
fermipy.roi_model.
Model
(name, data=None)[source]¶ Bases:
object
Base class for source objects. This class is a container for both spectral and spatial parameters as well as other source properties such as TS, Npred, and location within the ROI.
-
assoc
¶
-
data
¶
-
name
¶
-
names
¶
-
params
¶
-
spatial_pars
¶
-
spectral_pars
¶
-
-
class
fermipy.roi_model.
ROIModel
(config=None, **kwargs)[source]¶ Bases:
fermipy.config.Configurable
This class is responsible for managing the ROI model (both sources and diffuse components). Source catalogs can be read from either FITS or XML files. Individual components are represented by instances of
Model
and can be accessed by name using the bracket operator.- Create an ROI with all 3FGL sources and print a summary of its contents:
>>> skydir = astropy.coordinates.SkyCoord(0.0,0.0,unit='deg') >>> roi = ROIModel({'catalogs' : ['3FGL'],'src_roiwidth' : 10.0},skydir=skydir) >>> print roi name SpatialModel SpectrumType offset ts npred -------------------------------------------------------------------------------- 3FGL J2357.3-0150 PointSource PowerLaw 1.956 nan 0.0 3FGL J0006.2+0135 PointSource PowerLaw 2.232 nan 0.0 3FGL J0016.3-0013 PointSource PowerLaw 4.084 nan 0.0 3FGL J0014.3-0455 PointSource PowerLaw 6.085 nan 0.0
- Print a summary of an individual source
>>> print roi['3FGL J0006.2+0135']
- Get the SkyCoord for a source
>>> dir = roi['SourceA'].skydir
- Loop over all sources and print their names
>>> for s in roi.sources: print s.name
-
static
create_from_position
(skydir, config, **kwargs)[source]¶ Create an ROIModel instance centered on a sky direction.
Parameters:
-
static
create_from_source
(name, config, **kwargs)[source]¶ Create an ROI centered on the given source.
-
static
create_roi_from_ft1
(ft1file, config)[source]¶ Create an ROI model by extracting the sources coordinates form an FT1 file.
-
create_source
(name, src_dict, build_index=True, merge_sources=True)[source]¶ Add a new source to the ROI model from a dictionary or an existing source object.
Parameters: Returns: src
Return type:
-
defaults
= {'logfile': (None, '', <type 'str'>), 'catalogs': (None, '', <type 'list'>), 'src_roiwidth': (None, 'Width of square selection cut for inclusion of catalog sources in the model. Includes sources within a square region with side ``src_roiwidth`` centered on the ROI. If this parameter is none then no selection is applied. This selection will be ORed with the ``src_radius`` selection.', <type 'float'>), 'limbdiff': (None, '', <type 'list'>), 'src_radius_roi': (None, 'Half-width of ``src_roiwidth`` selection. This parameter can be used in lieu of ``src_roiwidth``.', <type 'float'>), 'extract_diffuse': (False, 'Extract a copy of all mapcube components centered on the ROI.', <type 'bool'>), 'galdiff': (None, 'Set the galactic IEM mapcube.', <type 'list'>), 'extdir': (None, 'Set a directory that will be searched for extended source FITS templates. Template files in this directory will take precendence over catalog source templates with the same name.', <type 'str'>), 'sources': (None, '', <type 'list'>), 'fileio': {'usescratch': (False, 'Run analysis in a temporary working directory under ``scratchdir``.', <type 'bool'>), 'scratchdir': ('/scratch', 'Path to the scratch directory. If ``usescratch`` is True then a temporary working directory will be created under this directory.', <type 'str'>), 'savefits': (True, 'Save intermediate FITS files.', <type 'bool'>), 'workdir': (None, 'Path to the working directory.', <type 'str'>), 'outdir_regex': (['\\.fits$|\\.fit$|\\.xml$|\\.npy$|\\.png$|\\.pdf$|\\.yaml$'], 'Stage files to the output directory that match at least one of the regular expressions in this list. This option only takes effect when ``usescratch`` is True.', <type 'list'>), 'workdir_regex': (['\\.fits$|\\.fit$|\\.xml$|\\.npy$'], 'Stage files to the working directory that match at least one of the regular expressions in this list. This option only takes effect when ``usescratch`` is True.', <type 'list'>), 'logfile': (None, 'Path to log file. If None then log will be written to fermipy.log.', <type 'str'>), 'outdir': (None, 'Path of the output directory. If none this will default to the directory containing the configuration file.', <type 'str'>)}, 'logging': {'verbosity': (3, '', <type 'int'>), 'chatter': (3, 'Set the chatter parameter of the STs.', <type 'int'>)}, 'isodiff': (None, 'Set the isotropic template.', <type 'list'>), 'merge_sources': (True, 'Merge properties of sources that appear in multiple source catalogs. If merge_sources=false then subsequent sources with the same name will be ignored.', <type 'bool'>), 'assoc_xmatch_columns': (['3FGL_Name'], 'Choose a set of association columns on which to cross-match catalogs.', <type 'list'>), 'diffuse': (None, '', <type 'list'>), 'src_radius': (None, 'Radius of circular selection cut for inclusion of catalog sources in the model. Includes sources within a circle of this radius centered on the ROI. If this parameter is none then no selection is applied. This selection will be ORed with the ``src_roiwidth`` selection.', <type 'float'>)}¶
-
diffuse_sources
¶
-
get_source_by_name
(name)[source]¶ Return a single source in the ROI with the given name. The input name string can match any of the strings in the names property of the source object. Case and whitespace are ignored when matching name strings. If no sources are found or multiple sources then an exception is thrown.
Parameters: name (str) – Name string. Returns: srcs – A source object. Return type: Model
-
get_sources
(skydir=None, distance=None, cuts=None, minmax_ts=None, minmax_npred=None, square=False, exclude_diffuse=False, coordsys='CEL')[source]¶ Retrieve list of sources satisfying the given selections.
Returns: srcs – List of source objects. Return type: list
-
get_sources_by_name
(name)[source]¶ Return a list of sources in the ROI matching the given name. The input name string can match any of the strings in the names property of the source object. Case and whitespace are ignored when matching name strings.
Parameters: name (str) – Returns: srcs – A list of Model
objects.Return type: list
-
get_sources_by_position
(skydir, dist, min_dist=None, square=False, coordsys='CEL')[source]¶ Retrieve sources within a certain angular distance of a sky coordinate. This function supports two types of geometric selections: circular (square=False) and square (square=True). The circular selection finds all sources with a given angular distance of the target position. The square selection finds sources within an ROI-like region of size R x R where R = 2 x dist.
Parameters: - skydir (
SkyCoord
) – Sky direction with respect to which the selection will be applied. - dist (float) – Maximum distance in degrees from the sky coordinate.
- square (bool) – Choose whether to apply a circular or square selection.
- coordsys (str) – Coordinate system to use when applying a selection with square=True.
- skydir (
-
load_fits_catalog
(name, **kwargs)[source]¶ Load sources from a FITS catalog file.
Parameters: name (str) – Catalog name or path to a catalog FITS file.
-
load_source
(src, build_index=True, merge_sources=True, **kwargs)[source]¶ Load a single source.
Parameters:
-
match_source
(src)[source]¶ Look for source or sources in the model that match the given source. Sources are matched by name and any association columns defined in the assoc_xmatch_columns parameter.
-
point_sources
¶
-
projection
¶
-
skydir
¶ Return the sky direction corresponding to the center of the ROI.
-
sources
¶
-
src_name_cols
= ['Source_Name', 'ASSOC', 'ASSOC1', 'ASSOC2', 'ASSOC_GAM', '1FHL_Name', '2FGL_Name', '3FGL_Name', 'ASSOC_GAM1', 'ASSOC_GAM2', 'ASSOC_TEV']¶
-
class
fermipy.roi_model.
Source
(name, data, radec=None)[source]¶ Bases:
fermipy.roi_model.Model
Class representation of a source (non-diffuse) model component. A source object serves as a container for the properties of that source (position, spatial/spectral parameters, TS, etc.) as derived in the current analysis. Most properties of a source object can be accessed with the bracket operator:
# Return the TS of this source >>> print src[‘ts’]
# Get a skycoord representation of the source position >>> print src.skydir
-
associations
¶
-
static
create_from_dict
(src_dict, roi_skydir=None)[source]¶ Create a source object from a python dictionary.
Parameters: src_dict (dict) – Dictionary defining the properties of the source.
-
data
¶
-
diffuse
¶
-
extended
¶
-
radec
¶
-
-
fermipy.roi_model.
create_source_table
(scan_shape)[source]¶ Create an empty source table.
Returns: tab Return type: Table
-
fermipy.roi_model.
get_skydir_distance_mask
(src_skydir, skydir, dist, min_dist=None, square=False, coordsys='CEL')[source]¶ Retrieve sources within a certain angular distance of an (ra,dec) coordinate. This function supports two types of geometric selections: circular (square=False) and square (square=True). The circular selection finds all sources with a given angular distance of the target position. The square selection finds sources within an ROI-like region of size R x R where R = 2 x dist.
Parameters: - src_skydir (
SkyCoord
) – Array of sky directions. - skydir (
SkyCoord
) – Sky direction with respect to which the selection will be applied. - dist (float) – Maximum distance in degrees from the sky coordinate.
- square (bool) – Choose whether to apply a circular or square selection.
- coordsys (str) – Coordinate system to use when applying a selection with square=True.
- src_skydir (
fermipy.utils module¶
-
fermipy.utils.
cl_to_dlnl
(cl)[source]¶ Compute the delta-log-likehood corresponding to an upper limit of the given confidence level.
-
fermipy.utils.
collect_dirs
(path, max_depth=1, followlinks=True)[source]¶ Recursively find directories under the given path.
-
fermipy.utils.
convolve2d_disk
(fn, r, sig, nstep=200)[source]¶ Evaluate the convolution f’(r) = f(r) * g(r) where f(r) is azimuthally symmetric function in two dimensions and g is a step function given by:
g(r) = H(1-r/s)
Parameters:
-
fermipy.utils.
convolve2d_gauss
(fn, r, sig, nstep=200)[source]¶ Evaluate the convolution f’(r) = f(r) * g(r) where f(r) is azimuthally symmetric function in two dimensions and g is a gaussian given by:
g(r) = 1/(2*pi*s^2) Exp[-r^2/(2*s^2)]
Parameters:
-
fermipy.utils.
create_model_name
(src)[source]¶ Generate a name for a source object given its spatial/spectral properties.
Parameters: src ( Source
) – A source object.Returns: name – A source name. Return type: str
-
fermipy.utils.
extend_array
(edges, binsz, lo, hi)[source]¶ Extend an array to encompass lo and hi values.
-
fermipy.utils.
find_function_root
(fn, x0, xb, delta=0.0)[source]¶ Find the root of a function: f(x)+delta in the interval encompassed by x0 and xb.
Parameters: - fn (function) – Python function.
- x0 (float) – Fixed bound for the root search. This will either be used as the lower or upper bound depending on the relative value of xb.
- xb (float) – Upper or lower bound for the root search. If a root is not found in the interval [x0,xb]/[xb,x0] this value will be increased/decreased until a change in sign is found.
-
fermipy.utils.
fit_parabola
(z, ix, iy, dpix=2, zmin=None)[source]¶ Fit a parabola to a 2D numpy array. This function will fit a parabola with the functional form described in
parabola
to a 2D slice of the input arrayz
. The boundaries of the fit region within z are set with the pixel centroid (ix
andiy
) and region size (dpix
).Parameters:
-
fermipy.utils.
get_parameter_limits
(xval, loglike, ul_confidence=0.95, tol=0.001)[source]¶ Compute upper/lower limits, peak position, and 1-sigma errors from a 1-D likelihood function. This function uses the delta-loglikelihood method to evaluate parameter limits by searching for the point at which the change in the log-likelihood value with respect to the maximum equals a specific value. A parabolic spline fit to the log-likelihood values is used to improve the accuracy of the calculation.
Parameters:
-
fermipy.utils.
init_matplotlib_backend
(backend=None)[source]¶ This function initializes the matplotlib backend. When no DISPLAY is available the backend is automatically set to ‘Agg’.
Parameters: backend (str) – matplotlib backend name.
-
fermipy.utils.
isstr
(s)[source]¶ String instance testing method that works under both Python 2.X and 3.X. Returns true if the input is a string.
-
fermipy.utils.
load_data
(infile, workdir=None)[source]¶ Load python data structure from either a YAML or numpy file.
-
fermipy.utils.
make_cdisk_kernel
(psf, sigma, npix, cdelt, xpix, ypix, normalize=False)[source]¶ Make a kernel for a PSF-convolved 2D disk.
Parameters: - psf (
PSFModel
) – - sigma (float) – 68% containment radius in degrees.
- psf (
-
fermipy.utils.
make_cgauss_kernel
(psf, sigma, npix, cdelt, xpix, ypix, normalize=False)[source]¶ Make a kernel for a PSF-convolved 2D gaussian.
Parameters: - psf (
PSFModel
) – - sigma (float) – 68% containment radius in degrees.
- psf (
-
fermipy.utils.
make_disk_kernel
(sigma, npix=501, cdelt=0.01, xpix=0.0, ypix=0.0)[source]¶ Make kernel for a 2D disk.
Parameters: sigma (float) – Disk radius in deg.
-
fermipy.utils.
make_gaussian_kernel
(sigma, npix=501, cdelt=0.01, xpix=0.0, ypix=0.0)[source]¶ Make kernel for a 2D gaussian.
Parameters: sigma (float) – 68% containment radius in degrees.
-
fermipy.utils.
make_pixel_offset
(npix, xpix=0.0, ypix=0.0)[source]¶ Make a 2D array with the distance of each pixel from a reference direction in pixel coordinates. Pixel coordinates are defined such that (0,0) is located at the center of the coordinate grid.
-
fermipy.utils.
make_psf_kernel
(psf, npix, cdelt, xpix, ypix, normalize=False)[source]¶ Generate a kernel for a point-source.
Parameters:
-
fermipy.utils.
match_regex_list
(patterns, string)[source]¶ Perform a regex match of a string against a list of patterns. Returns true if the string matches at least one pattern in the list.
-
fermipy.utils.
merge_dict
(d0, d1, add_new_keys=False, append_arrays=False)[source]¶ Recursively merge the contents of python dictionary d0 with the contents of another python dictionary, d1.
Parameters:
-
fermipy.utils.
parabola
(xy, amplitude, x0, y0, sx, sy, theta)[source]¶ Evaluate a 2D parabola given by:
f(x,y) = f_0 - (1/2) * delta^T * R * Sigma * R^T * delta
where
delta = [(x - x_0), (y - y_0)]
and R is the matrix for a 2D rotation by angle heta and Sigma is the covariance matrix:
- Sigma = [[1/sigma_x^2, 0 ],
- [0 , 1/sigma_y^2]]
Parameters: - xy (tuple) – Tuple containing x and y arrays for the values at which the parabola will be evaluated.
- amplitude (float) – Constant offset value.
- x0 (float) – Centroid in x coordinate.
- y0 (float) – Centroid in y coordinate.
- sx (float) – Standard deviation along first axis (x-axis when theta=0).
- sy (float) – Standard deviation along second axis (y-axis when theta=0).
- theta (float) – Rotation angle in radians.
Returns: vals – Values of the parabola evaluated at the points defined in the
xy
input tuple.Return type:
-
fermipy.utils.
project
(lon0, lat0, lon1, lat1)[source]¶ This function performs a stereographic projection on the unit vector (lon1,lat1) with the pole defined at the reference unit vector (lon0,lat0).
-
fermipy.utils.
tolist
(x)[source]¶ convenience function that takes in a nested structure of lists and dictionaries and converts everything to its base objects. This is useful for dupming a file to yaml.
numpy arrays into python lists
>>> type(tolist(np.asarray(123))) == int True >>> tolist(np.asarray([1,2,3])) == [1,2,3] True
numpy strings into python strings.
>>> tolist([np.asarray('cat')])==['cat'] True
an ordered dict to a dict
>>> ordered=OrderedDict(a=1, b=2) >>> type(tolist(ordered)) == dict True
converts unicode to regular strings
>>> type(u'a') == str False >>> type(tolist(u'a')) == str True
converts numbers & bools in strings to real represntation, (i.e. ‘123’ -> 123)
>>> type(tolist(np.asarray('123'))) == int True >>> type(tolist('123')) == int True >>> tolist('False') == False True
fermipy.sed module¶
Utilities for dealing with SEDs
- Many parts of this code are taken from dsphs/like/lnlfn.py by
- Matthew Wood <mdwood@slac.stanford.edu> Alex Drlica-Wagner <kadrlica@slac.stanford.edu>
-
class
fermipy.sed.
SEDGenerator
[source]¶ Bases:
object
Mixin class that provides SED functionality to
GTAnalysis
.-
sed
(name, **kwargs)[source]¶ Generate a spectral energy distribution (SED) for a source. This function will fit the normalization of the source in each energy bin. By default the SED will be generated with the analysis energy bins but a custom binning can be defined with the
loge_bins
parameter.Parameters: - name (str) – Source name.
- prefix (str) – Optional string that will be prepended to all output files (FITS and rendered images).
- loge_bins (
ndarray
) – Sequence of energies in log10(E/MeV) defining the edges of the energy bins. If this argument is None then the analysis energy bins will be used. The energies in this sequence must align with the bin edges of the underyling analysis instance. - bin_index (float) – Spectral index that will be use when fitting the energy distribution within an energy bin.
- use_local_index (bool) – Use a power-law approximation to the shape of the global
spectrum in each bin. If this is false then a constant
index set to
bin_index
will be used. - fix_background (bool) – Fix background components when fitting the flux normalization in each energy bin. If fix_background=False then all background parameters that are currently free in the fit will be profiled. By default fix_background=True.
- ul_confidence (float) – Set the confidence level that will be used for the calculation of flux upper limits in each energy bin.
- cov_scale (float) – Scaling factor that will be applied when setting the gaussian prior on the normalization of free background sources. If this parameter is None then no gaussian prior will be applied.
- write_fits (bool) – Write a FITS file containing the SED analysis results.
- write_npy (bool) – Write a numpy file with the contents of the output dictionary.
- optimizer (dict) – Dictionary that overrides the default optimizer settings.
Returns: sed – Dictionary containing output of the SED analysis. This dictionary is also saved to the ‘sed’ dictionary of the
Source
instance.Return type:
-
fermipy.sourcefind module¶
-
class
fermipy.sourcefind.
SourceFinder
[source]¶ Bases:
object
Mixin class which provides source-finding functionality to
GTAnalysis
.-
find_sources
(prefix='', **kwargs)[source]¶ An iterative source-finding algorithm.
Parameters: - model (dict) – Dictionary defining the properties of the test source. This is the model that will be used for generating TS maps.
- sqrt_ts_threshold (float) – Source threshold in sqrt(TS). Only peaks with sqrt(TS) exceeding this threshold will be used as seeds for new sources.
- min_separation (float) – Minimum separation in degrees of sources detected in each iteration. The source finder will look for the maximum peak in the TS map within a circular region of this radius.
- max_iter (int) – Maximum number of source finding iterations. The source finder will continue adding sources until no additional peaks are found or the number of iterations exceeds this number.
- sources_per_iter (int) – Maximum number of sources that will be added in each iteration. If the number of detected peaks in a given iteration is larger than this number, only the N peaks with the largest TS will be used as seeds for the current iteration.
- tsmap_fitter (str) –
Set the method used internally for generating TS maps. Valid options:
- tsmap
- tscube
- tsmap (dict) – Keyword arguments dictionary for tsmap method.
- tscube (dict) – Keyword arguments dictionary for tscube method.
Returns: - peaks (list) – List of peak objects.
- sources (list) – List of source objects.
-
localize
(name, **kwargs)[source]¶ Find the best-fit position of a source. Localization is performed in two steps. First a TS map is computed centered on the source with half-width set by
dtheta_max
. A fit is then performed to the maximum TS peak in this map. The source position is then further refined by scanning the likelihood in the vicinity of the peak found in the first step. The size of the scan region is set to encompass the 99% positional uncertainty contour as determined from the peak fit.Parameters: - name (str) – Source name.
- dtheta_max (float) – Maximum offset in RA/DEC in deg from the nominal source position that will be used to define the boundaries of the TS map search region.
- nstep (int) – Number of steps in longitude/latitude that will be taken when refining the source position. The bounds of the scan range are set to the 99% positional uncertainty as determined from the TS map peak fit. The total number of sampling points will be nstep**2.
- fix_background (bool) – Fix background parameters when fitting the source position.
- update (bool) – Update the model for this source with the best-fit position. If newname=None this will overwrite the existing source map of this source with one corresponding to its new location.
- newname (str) – Name that will be assigned to the relocalized source when update=True. If newname is None then the existing source name will be used.
- optimizer (dict) – Dictionary that overrides the default optimizer settings.
Returns: localize – Dictionary containing results of the localization analysis. This dictionary is also saved to the dictionary of this source in ‘localize’.
Return type:
-
fermipy.spectrum module¶
-
class
fermipy.spectrum.
PowerLaw
(params, scale=1.0)[source]¶ Bases:
fermipy.spectrum.SpectralFunction
Class that evaluates a power-law function with the parameterization:
F(x) = p_0 * (x/x_s)^p_1
where x_s is the scale parameter. The `params ` array should be defined with:
- params[0] : Prefactor (p_0)
- params[1] : Index (p_1)
-
class
fermipy.spectrum.
SEDEFluxFunctor
(spectrum, scale, emin, emax)[source]¶ Bases:
fermipy.spectrum.SEDFunctor
Functor that computes the energy flux of a source in a sequence of SED energy bins.
-
class
fermipy.spectrum.
SEDFluxFunctor
(spectrum, scale, emin, emax)[source]¶ Bases:
fermipy.spectrum.SEDFunctor
Functor that computes the flux of a source in a sequence of SED energy bins.
-
class
fermipy.spectrum.
SEDFunctor
(spectrum, scale, emin, emax)[source]¶ Bases:
object
Functor that accepts a model parameter vector and computes the normalization of a spectral model in a sequence of SED energy bins.
-
emax
¶
-
emin
¶
-
scale
¶
-
-
class
fermipy.spectrum.
SpectralFunction
(params, scale=1.0)[source]¶ Bases:
object
Base class for spectral model classes.
-
e2dfde_deriv
(x, params=None)[source]¶ Evaluate derivative of E^2 times differential flux with respect to E.
-
edfde_deriv
(x, params=None)[source]¶ Evaluate derivative of E times differential flux with respect to E.
-
log_params
¶
-
params
¶
-
scale
¶
-
fermipy.skymap module¶
-
class
fermipy.skymap.
HpxMap
(counts, hpx)[source]¶ Bases:
fermipy.skymap.Map_Base
Representation of a 2D or 3D counts map using HEALPix.
-
convert_to_cached_wcs
(hpx_in, sum_ebins=False, normalize=True)[source]¶ Make a WCS object and convert HEALPix data into WCS projection
Parameters:
-
static
create_from_hdu
(hdu, ebins)[source]¶ Creates and returns an HpxMap object from a FITS HDU.
hdu : The FITS ebins : Energy bin edges [optional]
-
static
create_from_hdulist
(hdulist, extname='SKYMAP', ebounds='EBOUNDS')[source]¶ Creates and returns an HpxMap object from a FITS HDUList
extname : The name of the HDU with the map data ebounds : The name of the HDU with the energy bin data
-
hpx
¶
-
-
class
fermipy.skymap.
Map
(counts, wcs)[source]¶ Bases:
fermipy.skymap.Map_Base
Representation of a 2D or 3D counts map using WCS.
-
get_map_values
(lons, lats)[source]¶ Return the indices in the flat array corresponding to a set of coordinates
Parameters: - lons (array-like) – ‘Longitudes’ (RA or GLON)
- lats (array-like) – ‘Latitidues’ (DEC or GLAT)
Returns: vals – Values of pixels in the flattened map, np.nan used to flag coords outside of map
Return type: numpy.ndarray((n))
-
get_pixel_indices
(lons, lats)[source]¶ Return the indices in the flat array corresponding to a set of coordinates
Parameters: - lons (array-like) – ‘Longitudes’ (RA or GLON)
- lats (array-like) – ‘Latitidues’ (DEC or GLAT)
Returns: idxs – Indices of pixels in the flattened map, -1 used to flag coords outside of map
Return type: numpy.ndarray((n),’i’)
-
ipix_swap_axes
(ipix, colwise=False)[source]¶ Return the transposed pixel index from the pixel xy coordinates
if colwise is True (False) this assumes the original index was in column wise scheme
-
ipix_to_xypix
(ipix, colwise=False)[source]¶ Return the pixel xy coordinates from the pixel index
if colwise is True (False) this uses columnwise (rowwise) indexing
-
pix_center
¶ Return the ROI center in pixel coordinates.
-
pix_size
¶ Return the pixel size along the two image dimensions.
-
skydir
¶ Return the sky coordinate of the image center.
-
wcs
¶
-
width
¶ Return the dimensions of the image.
-
fermipy.castro module¶
Utilities for dealing with ‘castro data’, i.e., 2D table of likelihood values.
Castro data can be tabluated in terms of a variety of variables. The most common example is probably a simple SED, where we have the likelihood as a function of Energy and Energy Flux.
However, we could easily convert to the likelihood as a function of other variables, such as the Flux normalization and the spectral index, or the mass and cross-section of a putative dark matter particle.
-
class
fermipy.castro.
CastroData
(norm_vals, nll_vals, refSpec, norm_type)[source]¶ Bases:
fermipy.castro.CastroData_Base
This class wraps the data needed to make a “Castro” plot, namely the log-likelihood as a function of normalization for a series of energy bins.
-
static
create_from_fits
(fitsfile, norm_type='EFLUX', hdu_scan='SCANDATA', hdu_energies='EBOUNDS', irow=None)[source]¶ Create a CastroData object from a fits file.
Parameters: - fitsfile (str) – Name of the fits file
- norm_type (str) –
Type of normalization to use. Valid options are:
- NORM : Normalization w.r.t. to test source
- FLUX : Flux of the test source ( ph cm^-2 s^-1 )
- EFLUX: Energy Flux of the test source ( MeV cm^-2 s^-1 )
- NPRED: Number of predicted photons (Not implemented)
- DFDE : Differential flux of the test source ( ph cm^-2 s^-1 MeV^-1 )
- EDFDE: Differential energy flux of the test source ( MeV cm^-2 s^-1 MeV^-1 ) (Not Implemented)
- hdu_scan (str) – Name of the FITS HDU with the scan data
- hdu_energies (str) – Name of the FITS HDU with the energy binning and normalization data
- irow (int or None) – If none, then this assumes that there is a single row in the scan data table Otherwise, this specifies which row of the table to use
Returns: castro
Return type:
-
static
create_from_flux_points
(txtfile)[source]¶ Create a Castro data object from a text file containing a sequence of differential flux points.
-
static
create_from_sedfile
(fitsfile, norm_type='EFLUX')[source]¶ Create a CastroData object from an SED fits file
Parameters: - fitsfile (str) – Name of the fits file
- norm_type (str) –
Type of normalization to use, options are:
- NORM : Normalization w.r.t. to test source
- FLUX : Flux of the test source ( ph cm^-2 s^-1 )
- EFLUX: Energy Flux of the test source ( MeV cm^-2 s^-1 )
- NPRED: Number of predicted photons (Not implemented)
- DFDE : Differential flux of the test source ( ph cm^-2 s^-1 MeV^-1 )
- EDFDE: Differential energy flux of the test source ( MeV cm^-2 s^-1 MeV^-1 ) (Not Implemented)
Returns: castro
Return type:
-
static
create_from_stack
(shape, components, ylims, weights=None)[source]¶ Combine the log-likelihoods from a number of components.
Parameters: - shape (tuple) – The shape of the return array
- components ([~fermipy.castro.CastroData_Base]) – The components to be stacked
- weights (array-like) –
Returns: castro
Return type:
-
static
create_from_tables
(norm_type='EFLUX', tab_s='SCANDATA', tab_e='EBOUNDS')[source]¶ Create a CastroData object from two tables
Parameters: - norm_type (str) –
Type of normalization to use. Valid options are:
- NORM : Normalization w.r.t. to test source
- FLUX : Flux of the test source ( ph cm^-2 s^-1 )
- EFLUX: Energy Flux of the test source ( MeV cm^-2 s^-1 )
- NPRED: Number of predicted photons (Not implemented)
- DFDE : Differential flux of the test source ( ph cm^-2 s^-1 MeV^-1 )
- EDFDE: Differential energy flux of the test source ( MeV cm^-2 s^-1 MeV^-1 ) (Not Implemented)
- tab_s (str) – table scan data
- tab_e (str) – table energy binning and normalization data
Returns: castro
Return type: - norm_type (str) –
-
create_functor
(specType, scale=1000.0)[source]¶ Create a functor object that computes normalizations in a sequence of energy bins for a given spectral model.
Parameters: Returns: - fn (
SpectralFunction
) – The functor - initPars (
ndarray
) – Default set of initial parameter for this spectral type - scale (float) – Energy scale (same as input)
- fn (
-
nE
¶ Return the number of energy bins. This is also the number of x-axis bins.
-
refSpec
¶ Return a
ReferenceSpec
with the spectral data
-
spectrum_loglike
(specType, params, scale=1000.0)[source]¶ return the log-likelihood for a particular spectrum
Parameters:
-
test_spectra
(spec_types=['PowerLaw', 'LogParabola', 'PLExpCutoff'])[source]¶ Test different spectral types against the SED represented by this CastroData.
Parameters: spec_types ([str,...]) – List of spectral types to try Returns: retDict – A dictionary of dictionaries. The top level dictionary is keyed by spec_type. The sub-dictionaries each contain: - “Function” :
SpectralFunction
- “Result” : tuple with the output of scipy.optimize.fmin
- “Spectrum” :
ndarray
with best-fit spectral values - “ScaleEnergy” : float, the ‘pivot energy’ value
- “TS” : float, the TS for the best-fit spectrum
Return type: dict - “Function” :
-
static
-
class
fermipy.castro.
CastroData_Base
(norm_vals, nll_vals, norm_type)[source]¶ Bases:
object
This class wraps the data needed to make a “Castro” plot, namely the log-likelihood as a function of normalization.
In this case the x-axes and y-axes are generic Sub-classes can implement particul axes choices (e.g., EFlux v. Energy)
-
__call__
(x)[source]¶ Return the negative log-likelihood for an array of values, summed over the energy bins
Parameters: x ( ndarray
) – Array of N x M valuesReturns: nll_val – Array of negative log-likelihood values. Return type: ndarray
-
derivative
(x, der=1)[source]¶ Return the derivate of the log-like summed over the energy bins
Parameters: Returns: der_val – Array of negative log-likelihood values.
Return type:
-
fitNorm_v2
(specVals)[source]¶ Fit the normalization given a set of spectral values that define a spectral shape.
This version uses
scipy.optimize.fmin
.Parameters: - specVals (an array of (nebin values that define a spectral shape) –
- xlims (fit limits) –
Returns: norm – Best-fit normalization value
Return type:
-
fitNormalization
(specVals, xlims)[source]¶ Fit the normalization given a set of spectral values that define a spectral shape
This version is faster, and solves for the root of the derivatvie
Parameters: - specVals (an array of (nebin values that define a spectral shape) –
- xlims (fit limits) –
- the best-fit normalization value (returns) –
-
fit_spectrum
(specFunc, initPars)[source]¶ Fit for the free parameters of a spectral function
Parameters: - specFunc (
SpectralFunction
) – The Spectral Function - initPars (
ndarray
) – The initial values of the parameters
Returns: - result (tuple) – The output of scipy.optimize.fmin
- spec_out (
ndarray
) – The best-fit spectral values - TS_spec (float) – The TS of the best-fit spectrum
- specFunc (
-
fn_mles
()[source]¶ returns the summed likelihood at the maximum likelihood estimate
Note that simply sums the maximum likelihood values at each bin, and does not impose any sort of constrain between bins
-
getLimits
(alpha, upper=True)[source]¶ Evaluate the limits corresponding to a C.L. of (1-alpha)%.
Parameters:
-
nll_null
¶ Return the negative log-likelihood for the null-hypothesis
-
norm_type
¶ Return the normalization type flag
-
nx
¶ Return the number of profiles
-
ny
¶ Return the number of profiles
-
static
stack_nll
(shape, components, ylims, weights=None)[source]¶ Combine the log-likelihoods from a number of components.
Parameters: - shape (tuple) – The shape of the return array
- components (
CastroData_Base
) – The components to be stacked - weights (array-like) –
Returns: - norm_vals (‘numpy.ndarray’) – N X M array of Normalization values
- nll_vals (‘numpy.ndarray’) – N X M array of log-likelihood values
-
-
class
fermipy.castro.
Interpolator
(x, y)[source]¶ Bases:
object
Helper class for interpolating a 1-D function from a set of tabulated values.
Safely deals with overflows and underflows
-
__call__
(x)[source]¶ Return the interpolated values for an array of inputs
x : the inputs
Note that if any x value is outside the interpolation ranges this will return a linear extrapolation based on the slope at the endpoint
-
derivative
(x, der=1)[source]¶ return the derivative a an array of input values
x : the inputs der : the order of derivative
-
x
¶ return the x values used to construct the split
-
xmax
¶ return the maximum value over which the spline is defined
-
xmin
¶ return the minimum value over which the spline is defined
-
y
¶ return the y values used to construct the split
-
-
class
fermipy.castro.
LnLFn
(x, y, norm_type=0)[source]¶ Bases:
object
Helper class for interpolating a 1-D log-likelihood function from a set of tabulated values.
-
getInterval
(alpha)[source]¶ Evaluate the interval corresponding to a C.L. of (1-alpha)%.
Parameters: alpha (limit confidence level.) –
-
getLimit
(alpha, upper=True)[source]¶ Evaluate the limits corresponding to a C.L. of (1-alpha)%.
Parameters: - alpha (limit confidence level.) –
- upper (upper or lower limits.) –
-
interp
¶ return the underlying Interpolator object
-
mle
()[source]¶ return the maximum likelihood estimate
This will return the cached value, if it exists
-
norm_type
¶ Return a string specifying the quantity used for the normalization. This isn’t actually used in this class, but it is carried so that the class is self-describing. The possible values are open-ended.
-
-
class
fermipy.castro.
ReferenceSpec
(emin, emax, ref_dfde, ref_flux, ref_eflux, ref_npred, eref=None)[source]¶ Bases:
object
This class wraps data for a reference spectrum.
Parameters: - ne (
int
) – Number of energy bins - ebins (
ndarray
) – Array of bin edges. - emin (
ndarray
) – Array of lower bin edges. - emax (
ndarray
) – Array of upper bin edges. - bin_widths (
ndarray
) – Array of energy bin widths. - eref (
ndarray
) – Array of reference energies. Typically these are the geometric mean of the energy bins - ref_dfde (
ndarray
) – Array of differential photon flux values. - ref_flux (
ndarray
) – Array of integral photon flux values. - ref_eflux (
ndarray
) – Array of integral energy flux values. - ref_npred (
ndarray
) – Array of predicted number of photons in each energy bin.
-
bin_widths
¶
-
ebins
¶
-
emax
¶
-
emin
¶
-
eref
¶
-
log_ebins
¶
-
nE
¶
-
ref_dfde
¶
-
ref_eflux
¶ return the energy flux values
-
ref_flux
¶ return the flux values
-
ref_npred
¶ return the number of predicted events
- ne (
-
class
fermipy.castro.
SpecData
(ref_spec, norm_dfde, norm_flux, norm_eflux, norm_dfde_err=None)[source]¶ Bases:
fermipy.castro.ReferenceSpec
This class wraps spectral data, e.g., energy bin definitions, flux values and number of predicted photons
-
dfde
¶
-
dfde_err
¶
-
e2dfde
¶
-
e2dfde_err
¶
-
eflux
¶
-
flux
¶
-
norm_dfde
¶
-
norm_eflux
¶
-
norm_flux
¶
-
-
class
fermipy.castro.
TSCube
(tsmap, normmap, tscube, normcube, norm_vals, nll_vals, refSpec, norm_type)[source]¶ Bases:
object
A class wrapping a TSCube, which is a collection of CastroData objects for a set of directions. This class wraps a combination of:
- Pixel data,
- Pixel x Energy bin data,
- Pixel x Energy Bin x Normalization scan point data
-
static
create_from_fits
(fitsfile, norm_type='FLUX')[source]¶ Build a TSCube object from a fits file created by gttscube :param fitsfile: Path to the tscube FITS file.
Parameters: norm_type (str) – String specifying the quantity used for the normalization
-
find_and_refine_peaks
(threshold, min_separation=1.0, use_cumul=False)[source]¶ Run a simple peak-finding algorithm, and fit the peaks to paraboloids to extract their positions and error ellipses.
Parameters: - threshold (float) – Peak threshold in TS.
- min_separation (float) – Radius of region size in degrees. Sets the minimum allowable separation between peaks.
- use_cumul (bool) – If true, used the cumulative TS map (i.e., the TS summed over the energy bins) instead of the TS Map from the fit to and index=2 powerlaw.
Returns: peaks – List of dictionaries containing the location and amplitude of each peak. Output of
find_peaks
Return type:
-
find_sources
(threshold, min_separation=1.0, use_cumul=False, output_peaks=False, output_castro=False, output_specInfo=False, output_src_dicts=False, output_srcs=False)[source]¶
-
nE
¶ return the number of energy bins
-
nN
¶ return the number of sample points in each energy bin
-
normcube
¶ return the Cube of the normalization value per pixel / energy bin
-
normmap
¶ return the Map of the Best-fit normalization value
-
nvals
¶ Return the number of values in the tscube
-
refSpec
¶ Return the Spectral Data object
-
test_spectra_of_peak
(peak, spec_types=None)[source]¶ Test different spectral types against the SED represented by the CastroData corresponding to a single pixel in this TSCube
Parameters: spec_types ([str,...]) – List of spectral types to try Returns: - castro (
CastroData
) – The castro data object for the pixel corresponding to the peak - test_dict (dict) –
The dictionary returned by
test_spectra
- castro (
-
ts_cumul
¶ return the Map of the cumulative TestStatistic value per pixel (summed over energy bin)
-
tscube
¶ return the Cube of the TestStatistic value per pixel / energy bin
-
tsmap
¶ return the Map of the TestStatistic value
fermipy.tsmap module¶
-
class
fermipy.tsmap.
TSCubeGenerator
[source]¶ Bases:
object
-
tscube
(prefix='', **kwargs)[source]¶ Generate a spatial TS map for a source component with properties defined by the
model
argument. This method uses thegttscube
ST application for source fitting and will simultaneously fit the test source normalization as well as the normalizations of any background components that are currently free. The output of this method is a dictionary containingMap
objects with the TS and amplitude of the best-fit test source. By default this method will also save maps to FITS files and render them as image files.Parameters: - prefix (str) – Optional string that will be prepended to all output files (FITS and rendered images).
- model (dict) – Dictionary defining the properties of the test source.
- do_sed (bool) – Compute the energy bin-by-bin fits.
- nnorm (int) – Number of points in the likelihood v. normalization scan.
- norm_sigma (float) – Number of sigma to use for the scan range.
- tol (float) – Critetia for fit convergence (estimated vertical distance to min < tol ).
- tol_type (int) – Absoulte (0) or relative (1) criteria for convergence.
- max_iter (int) – Maximum number of iterations for the Newton’s method fitter
- remake_test_source (bool) – If true, recomputes the test source image (otherwise just shifts it)
- st_scan_level (int) –
- make_plots (bool) – Write image files.
- write_fits (bool) – Write a FITS file with the results of the analysis.
Returns: maps – A dictionary containing the
Map
objects for TS and source amplitude.Return type:
-
-
class
fermipy.tsmap.
TSMapGenerator
[source]¶ Bases:
object
Mixin class for
GTAnalysis
that generates TS maps.-
tsmap
(prefix='', **kwargs)[source]¶ Generate a spatial TS map for a source component with properties defined by the
model
argument. The TS map will have the same geometry as the ROI. The output of this method is a dictionary containingMap
objects with the TS and amplitude of the best-fit test source. By default this method will also save maps to FITS files and render them as image files.This method uses a simplified likelihood fitting implementation that only fits for the normalization of the test source. Before running this method it is recommended to first optimize the ROI model (e.g. by running
optimize()
).Parameters: - prefix (str) – Optional string that will be prepended to all output files (FITS and rendered images).
- model (dict) – Dictionary defining the properties of the test source.
- exclude (str or list of str) – Source or sources that will be removed from the model when computing the TS map.
- loge_bounds (list) – Restrict the analysis to an energy range (emin,emax) in log10(E/MeV) that is a subset of the analysis energy range. By default the full analysis energy range will be used. If either emin/emax are None then only an upper/lower bound on the energy range wil be applied.
- max_kernel_radius (float) – Set the maximum radius of the test source kernel. Using a smaller value will speed up the TS calculation at the loss of accuracy. The default value is 3 degrees.
- make_plots (bool) – Write image files.
- write_fits (bool) – Write a FITS file.
- write_npy (bool) – Write a numpy file.
Returns: maps – A dictionary containing the
Map
objects for TS and source amplitude.Return type:
-
-
fermipy.tsmap.
extract_images_from_tscube
(infile, outfile)[source]¶ Extract data from table HDUs in TSCube file and convert them to FITS images
-
fermipy.tsmap.
f_cash
(x, counts, background, model)[source]¶ Wrapper for cash statistics, that defines the model function.
Parameters:
-
fermipy.tsmap.
overlap_slices
(large_array_shape, small_array_shape, position)[source]¶ Modified version of
overlap_slices
.Get slices for the overlapping part of a small and a large array.
Given a certain position of the center of the small array, with respect to the large array, tuples of slices are returned which can be used to extract, add or subtract the small array at the given position. This function takes care of the correct behavior at the boundaries, where the small array is cut of appropriately.
Parameters: Returns: - slices_large (tuple of slices) –
Slices in all directions for the large array, such that
large_array[slices_large]
extracts the region of the large array that overlaps with the small array. - slices_small (slice) –
Slices in all directions for the small array, such that
small_array[slices_small]
extracts the region that is inside the large array.
- slices_large (tuple of slices) –
Slices in all directions for the large array, such that
fermipy.residmap module¶
-
class
fermipy.residmap.
ResidMapGenerator
[source]¶ Bases:
object
Mixin class for
GTAnalysis
that generates spatial residual maps from the difference of data and model maps smoothed with a user-defined spatial/spectral template. The map of residual significance can be interpreted in the same way as a TS map (the likelihood of a source at the given location).-
residmap
(prefix='', **kwargs)[source]¶ Generate 2-D spatial residual maps using the current ROI model and the convolution kernel defined with the
model
argument.Parameters: - prefix (str) – String that will be prefixed to the output residual map files.
- model (dict) – Dictionary defining the properties of the convolution kernel.
- exclude (str or list of str) – Source or sources that will be removed from the model when computing the residual map.
- loge_bounds (list) – Restrict the analysis to an energy range (emin,emax) in log10(E/MeV) that is a subset of the analysis energy range. By default the full analysis energy range will be used. If either emin/emax are None then only an upper/lower bound on the energy range wil be applied.
- make_plots (bool) – Write image files.
- write_fits (bool) – Write FITS files.
Returns: maps – A dictionary containing the
Map
objects for the residual significance and amplitude.Return type:
-
-
fermipy.residmap.
convolve_map
(m, k, cpix, threshold=0.001, imin=0, imax=None)[source]¶ Perform an energy-dependent convolution on a sequence of 2-D spatial maps.
Parameters: - m (
ndarray
) – 3-D map containing a sequence of 2-D spatial maps. First dimension should be energy. - k (
ndarray
) – 3-D map containing a sequence of convolution kernels (PSF) for each slice in m. This map should have the same dimension as m. - cpix (list) – Indices of kernel reference pixel in the two spatial dimensions.
- threshold (float) – Kernel amplitude
- imin (int) – Minimum index in energy dimension.
- imax (int) – Maximum index in energy dimension.
- m (