/ sherpa4.17 / ahelp / sample_photon_flux.html

Jump to: Description · Examples · PARAMETERS · Notes · Changes in CIAO · Bugs · See Also

AHELP for CIAO 4.17 Sherpa

sample_photon_flux

Synopsis

Return the photon flux distribution of a model.

Syntax

sample_photon_flux(lo=None, hi=None, id=None, num=1, scales=None,
correlated=False, numcores=None, bkg_id=None, model=None, otherids=(),
clip='hard')

Description

For each iteration, draw the parameter values of the model from a normal distribution, evaluate the model, and sum the model over the given range (the flux). The return array contains the flux and parameter values for each iteration. The units for the flux are as returned by `calc_photon_flux` .

Examples

Example 1

Calculate the photon flux distribution for the range 0.5 to 7, and plot up the resulting flux distribution (as a cumulative distribution):

>>> vals = sample_photon_flux(0.5, 7, num=1000)
>>> plot_cdf(vals[:, 0], name='flux')

Example 2

Repeat the above, but allowing the parameters to be correlated, and then calculate the 5, 50, and 95 percent quantiles of the photon flux distribution:

>>> cvals = sample_photon_flux(0.5, 7, num=1000, correlated=True)
>>> np.percentile(cvals[:, 0], [5, 50, 95])

Example 3

The photon flux of a component (or sub-set of components) can be calculated using the model argument. For the following case, an absorbed power-law was used to fit the data - `xsphabs.gal * powerlaw.pl` - and then the flux of just the power-law component is calculated. Note that the returned array has columns 'flux', 'gal.nh', 'pl.gamma', and 'pl.ampl' (that is flux and then the free parameters in the full model).

>>> vals = sample_photon_flux(0.5, 7, model=pl, num=1000, correlated=True)

Example 4

Calculate the 2-10 keV flux for the pl model using a joint fit to the datasets 1, 2, 3, and 4:

>>> vals = sample_photon_flux(2, 10, model=pl, id=1, otherids=(2,3,4),
...                           num=1000)

Example 5

Use the given parameter errors for sampling the parameter distribution. The fit must have three free parameters, and each parameter is sampled independently (in this case parerrs gives the sigma values for each parameter):

>>> parerrs = [0.25, 1.22, 1.04e-4]
>>> vals = sample_photon_flux(2, 10, num=5000, scales=parerrs)

Example 6

In this case the parameter errors are taken from the covariance analysis, using the `parmaxes` field since these are positive.

>>> covar()
>>> parerrs = get_covar_results().parmaxes
>>> vals = sample_photon_flux(0.5, 2, num=1000, scales=parerrs)

Example 7

Run covariance to estimate the parameter errors and then extract the covariance matrix from the results (as the `cmat` variable). This matrix is then used to define the parameter widths - including correlated terms - in the flux sampling, after being increased by ten percent. This is used to calculate both the absorbed ( `vals1` ) and unabsorbed ( `vals2` ) fluxes. Both arrays have columns: flux, gal.nh, pl.gamma, and pl.ampl.

>>> set_source(xsphabs.gal * powlaw1d.pl)
>>> fit()
>>> covar()
>>> cmat = get_covar_results().extra_output
>>> vals1 = sample_photon_flux(2, 10, num=5000, correlated=True,
...                            scales=1.1 * cmat)
>>> vals2 = sample_photon_flux(2, 10, num=5000, correlated=True,
...                            model=pl, scales=1.1 * cmat)

Example 8

Calculate the flux and error distribution using fits to all datasets:

>>> set_source(xsphabs.gal * xsapec.clus)
>>> set_source(2, gal * clus)
>>> set_source(3, gal * clus)
... fit the data
>>> vals = sample_photon_flux(0.5, 10, model=clus, num=10000)

Example 9

Calculate the flux and error distribution using fits to an explicit set of datasets (in this case datasets 1 and 2):

>>> vals = sample_photon_flux(0.5, 10, id=1, otherids=[2],
...                           model=clus, num=10000)

Example 10

Generate two sets of parameter values, where the parameter values in v1 are generated from a random distribution and then clipped to the hard limits of the parameters, and the values in v2 use the soft limits of the parameters. The last column in both v1 and v2 indicates whether the row had any clipped parameters. The flux1_filt and flux2_filt arrays indicate the photon-flux distribution after it has been filtered to remove any row with clipped parameters:

>>> v1 = sample_photon_flux(0.5, 2, num=1000)
>>> v2 = sample_photon_flux(0.5, 2, num=1000, clip='soft')
>>> flux1 = v1[:, 0]
>>> flux2 = v2[:, 0]
>>> flux1_filt = flux1[v1[:, -1] == 0]
>>> flux2_filt = flux2[v2[:, -1] == 0]

PARAMETERS

The parameters for this function are:

Parameter	Type information	Definition
lo	number, optional	The lower limit to use when summing up the signal. If not given then the lower value of the data grid is used.
hi	optional	The upper limit to use when summing up the signal. If not given then the upper value of the data grid is used.
id	int, str, or None, optional	The identifier of the data set to use. If `None` , the default value, then all datasets with associated models are used to calculate the errors and the model evaluation is done using the default dataset.
num	int, optional	The number of samples to create. The default is 1.
scales	array, optional	The scales used to define the normal distributions for the parameters. The size and shape of the array depends on the number of free parameters in the fit (n) and the value of the `correlated` parameter. When the parameter is `True` , scales must be given the covariance matrix for the free parameters (a n by n matrix that matches the parameter ordering used by Sherpa). For un-correlated parameters the covariance matrix can be used, or a one-dimensional array of n elements can be used, giving the width (specified as the sigma value of a normal distribution) for each parameter (e.g. the square root of the diagonal elements of the covariance matrix). If the scales parameter is not given then the covariance matrix is evaluated for the current model and best-fit parameters.
correlated	bool, optional	Should the correlation between the parameters be included when sampling the parameters? If not, then each parameter is sampled from independent distributions. In both cases a normal distribution is used.
numcores	optional	The number of CPU cores to use. The default is to use all the cores on the machine.
bkg_id	int, str, or None, optional	The identifier of the background component to use. This should only be set when the line to be measured is in the background model.
model	model, optional	The model to integrate. If left as `None` then the source model for the dataset will be used. This can be used to calculate the unabsorbed flux, as shown in the examples. The model must be part of the source expression.
otherids	sequence of integer and string ids, optional	The list of other datasets that should be included when calculating the errors to draw values from.
clip	{'hard', 'soft', 'none'}, optional	What clipping strategy should be applied to the sampled parameters. The default ('hard') is to fix values at their hard limits if they exceed them. A value of 'soft' uses the soft limits instead, and 'none' applies no clipping. The last column in the returned arrays indicates if the row had any clipped parameters (even when clip is set to 'none').

Return value

The return value from this function is:

The return array has the shape (num, n+2) , where n is the number of free parameters in the fit and num is the `num` parameter. The rows of this array contain the flux value, as calculated by `calc_photon_flux` , followed by the values of the thawed parameters used for that iteration, and then a flag column indicating if the parameters were clipped (1) or not (0). The order of the parameters matches the data returned by `get_fit_results` .

Notes

There are two ways to use this function to calculate fluxes from multiple sources. The first is to leave the `id` argument as `None` , in which case all available datasets will be used. Alternatively, the `id` and `otherids` arguments can be set to list the exact datasets to use, such as `id=1, otherids=(2,3,4)` .

The returned value contains all free parameters in the fit, even if they are not included in the model argument (e.g. when calculating an unabsorbed flux).

Changes in CIAO

Changed in CIAO 4.16

The random number generation is now controlled by the `set_rng` routine.

Changed in CIAO 4.13

The model, otherids, and clip parameters were added and the return value has an extra column.

Bugs

See the bugs pages on the Sherpa website for an up-to-date listing of known bugs.