Deproject¶
- class deproject.deproject.Deproject(radii, theta=<Quantity 360. deg>, angdist=None, cosmology=None)¶
Bases:
SpecStack
Support deprojecting a set of spectra (2-d concentric circular annuli).
- Parameters:
radii (AstroPy Quantity representing an angle on the sky) – The edges of each annulus, which must be circular, concentric, in ascending order, and >= 0. If there are n annuli then there are n+1 radii, since the start and end of the sequence must be given. The units are expected to be arcsec, arcminute, or degree.
theta (AstroPy Quantity (scalar or array) representing an angle) – The “fill factor” of each annulus, given by the azimuthal coverage of the shell in degrees. The value can be a scalar, so the same value is used for all annuli, or a sequence with a length matching the number of annuli. Since the annulus assumes circular symmetry there is no need to define the starting point of the measurement, for cases when the value is less than 360 degrees.
angdist (None or AstroPy.Quantity, optional) – The angular-diameter distance to the source. If not given then it is calculated using the source redshift along with the cosmology attribute.
cosmology (None or astropy.cosmology object, optional) – The cosmology used to convert redshift to an angular-diameter distance. This is used when angdist is None. If cosmology is None then the astropy.cosmology.Planck15 Cosmology object is used.
Examples
The following highly-simplified example fits a deprojected model to data from three annuli - ann1.pi, ann2.pi, and ann3.pi - and also calculates errors on the parameters using the confidence method:
>>> dep = Deproject([0, 10, 40, 100] * u.arcsec) >>> dep.load_pha('ann1.pi', 0) >>> dep.load_pha('ann2.pi', 1) >>> dep.load_pha('ann3.pi', 2) >>> dep.subtract() >>> dep.notice(0.5, 7.0) >>> dep.set_source('xsphabs * xsapec') >>> dep.set_par('xsapec.redshift', 0.23) >>> dep.thaw('xsapec.abundanc') >>> dep.set_par('xsphabs.nh', 0.087) >>> dep.freeze('xsphabs.nh') >>> dep.fit() >>> dep.fit_plot('rstat') >>> errs = dep.conf() >>> dep.conf_plot('density')
Attributes Summary
Angular size distance (an AstroPy quantity)
Return the cosmology object (only used if angdist not set)
Information (dataset identifier, annulus, name) about the loaded data.
How many datasets are registered?
Source redshift
Methods Summary
conf
()Estimate errors using confidence, using the "onion-peeling" method.
conf_plot
(field[, results, units, xlog, ...])Plot up the confidence errors as a function of radius.
covar
()Estimate errors using covariance, using the "onion-peeling" method.
covar_plot
(field[, results, units, xlog, ...])Plot up the covariance errors as a function of radius.
density_plot
([units, xlog, ylog, overplot, ...])Plot up the electron density as a function of radius.
dummyfunc
(*args, **kwargs)find_norm
(shell)Return the normalization value for the given shell.
find_parval
(parname)Return the value of the first parameter matching the name.
fit
()Fit the data using the "onion-peeling" method.
fit_plot
(field[, results, units, xlog, ...])Plot up the fit results as a function of radius.
freeze
(par)Freeze the given parameter in each shell.
What are the conf results, per annulus?
What are the covar results, per annulus?
Calculate the electron density for each shell.
What are the fit results, per annulus?
get_par
(par)Return the parameter value for each shell.
get_radii
([units])What are the radii of the shells?
How are the annuli grouped?
group
(*args)Apply the grouping for each data set.
guess
()Guess the starting point by fitting the projected data.
ignore
(*args)Apply Sherpa ignore command to each dataset.
load_pha
(specfile, annulus)Load a pha file and add to the datasets for stacked analysis.
notice
(*args)Apply Sherpa notice command to each dataset.
par_plot
(par[, units, xlog, ylog, overplot, ...])Plot up the parameter as a function of radius.
plot_arf
(*args, **kwargs)plot_bkg
(*args, **kwargs)plot_bkg_chisqr
(*args, **kwargs)plot_bkg_delchi
(*args, **kwargs)plot_bkg_fit
(*args, **kwargs)plot_bkg_fit_delchi
(*args, **kwargs)plot_bkg_fit_resid
(*args, **kwargs)plot_bkg_model
(*args, **kwargs)plot_bkg_ratio
(*args, **kwargs)plot_bkg_resid
(*args, **kwargs)plot_bkg_source
(*args, **kwargs)plot_bkg_unconvolved
(*args, **kwargs)plot_chisqr
(*args, **kwargs)plot_data
(*args, **kwargs)plot_delchi
(*args, **kwargs)plot_fit
(*args, **kwargs)plot_fit_delchi
(*args, **kwargs)plot_fit_resid
(*args, **kwargs)plot_model
(*args, **kwargs)plot_order
(*args, **kwargs)plot_psf
(*args, **kwargs)plot_ratio
(*args, **kwargs)plot_resid
(*args, **kwargs)plot_source
(*args, **kwargs)print_window
(*args, **kwargs)Create a hardcopy version of each plot window.
set_bkg_model
(bkgmodel)Create a background model for each annulus.
set_par
(par, val)Set the parameter value in each shell.
set_source
([srcmodel])Create a source model for each annulus.
subtract
(*args)Subtract the background from each dataset.
thaw
(par)Thaw the given parameter in each shell.
tie_par
(par, base, *others)Tie parameters in one or more shells to the base shell.
ungroup
(*args)Turn off the grouping for each data set.
unsubtract
(*args)Un-subtract the background from each dataset.
untie_par
(par, *others)Remove the parameter tie/link in the shell.
Attributes Documentation
- angdist¶
Angular size distance (an AstroPy quantity)
- cosmology¶
Return the cosmology object (only used if angdist not set)
- datasets = None¶
Information (dataset identifier, annulus, name) about the loaded data.
- n_datasets¶
How many datasets are registered?
This is not the same as the number of annuli.
- Returns:
ndata – The number of datasets.
- Return type:
int
- redshift¶
Source redshift
Methods Documentation
- conf()¶
Estimate errors using confidence, using the “onion-peeling” method.
It is assumed that
fit
has been called. The results can also be retrieved withget_conf_results
.- Returns:
errors – This records per-annulus data, such as the inner and outer radius (rlo_ang, rhi_ang, rlo_phys, rhi_phys), the sigma and percent values, and parameter results (accessed using <model name>.<par name>, <model name>.<par name>_lo, and <model name>.<par name>_hi syntax, where the match is case sensitive).
- Return type:
astropy.table.Table instance
See also
Examples
Run a fit and then error analysis, then plot up the abundance against temperature values including the error bars. Note that the Matplotlib errorbar routine requires “positive” error values whereas the <param>_lo values are negative, hence they are negated in the creation of
dkt
anddabund
:>>> dep.fit() >>> errs = dep.conf() >>> kt, abund = errs['xsapec.kT'], errs['xsapec.Abundanc'] >>> ktlo, kthi = errs['xsapec.kT_lo'], errs['xsapec.kT_hi'] >>> ablo, abhi = errs['xsapec.Abundanc_lo'], errs['xsapec.Abundanc_hi'] >>> dkt = np.vstack((-ktlo, kthi)) >>> dabund = np.vstack((-ablo, abhi)) >>> plt.clf() >>> plt.errorbar(kt, abund, xerr=dkt, yerr=dabund, fmt='.')
Plot up the temperature distibution as a function of radius, including the error bars calculated by the conf routine:
>>> dep.fit() >>> dep.conf() >>> dep.conf_plot('xsmekal.kt')
- conf_plot(field, results=None, units='kpc', xlog=True, ylog=False, overplot=False, clearwindow=True)¶
Plot up the confidence errors as a function of radius.
This method can be used to plot up the last conf results or a previously-stored set. Any error bars are shown at the scale they were calculated (as given by the
sigma
andpercent
columns of the results).- Parameters:
field (str) – The column to plot from the fit results (the match is case insensitive).
results (None or astropy.table.Table instance) – The return value from the
conf
orget_conf_results
methods.units (str or astropy.units.Unit, optional) – The X-axis units (a length or angle, such as ‘Mpc’ or ‘arcsec’, where the case is important).
xlog (bool, optional) – Should the x axis be drawn with a log scale (default True)?
ylog (bool, optional) – Should the y axis be drawn with a log scale (default False)?
overplot (bool, optional) – Clear the plot or add to existing plot?
clearwindow (bool, optional) – How does this interact with overplot?
See also
fit
,get_conf_results
,fit_plot
,covar_plot
,density_plot
,par_plot
Notes
Error bars are included on the dependent axis if the results contain columns that match the requested field with suffixes of ‘_lo’ and ‘_hi’. These error bars are asymmetric, which is different to
covar_plot
.If a limit is missing (i.e. it is a NaN) then no error bar is drawn. This can make it look like the error is very small.
Examples
Plot the temperature as a function of radius from the last fit, including error bars:
>>> dep.conf_plot('xsapec.kt')
Plot the density with the radii labelled in arcminutes and the density shown on a log scale:
>>> dep.conf_plot('density', units='arcmin', ylog=True)
Overplot the current conf results on those from a previous fit, where
conf1
was returned from theconf
orget_conf_results
methods:>>> dep.conf_plot('xsapec.abundanc', results=conf1) >>> dep.conf_plot('xsapec.abundanc', overplot=True)
- covar()¶
Estimate errors using covariance, using the “onion-peeling” method.
It is assumed that
fit
has been called. The results can also be retrieved withget_covar_results
.- Returns:
errors – This records per-annulus data, such as the inner and outer radius (rlo_ang, rhi_ang, rlo_phys, rhi_phys), the sigma and percent values, and parameter results (accessed using <model name>.<par name>, <model name>.<par name>_lo, and <model name>.<par name>_hi syntax, where the match is case sensitive). The _lo and _hi values are symmetric for covar, that is the _lo value will be the negative of the _hi value.
- Return type:
astropy.table.Table instance
See also
Examples
Run a fit and then error analysis, then plot up the abundance against temperature values including the error bars. Since the covariance routine returns symmetric error bars, the <param>_hi values are used in the plot:
>>> dep.fit() >>> errs = dep.covar() >>> kt, abund = errs['xsapec.kT'], errs['xsapec.Abundanc'] >>> dkt = errs['xsapec.kT_hi'] >>> dabund = errs['xsapec.Abundanc_hi'] >>> plt.clf() >>> plt.errorbar(kt, abund, xerr=dkt, yerr=dabund, fmt='.')
Plot up the temperature distibution as a function of radius, including the error bars calculated by the covar routine:
>>> dep.fit() >>> dep.covar() >>> dep.covar_plot('xsmekal.kt')
- covar_plot(field, results=None, units='kpc', xlog=True, ylog=False, overplot=False, clearwindow=True)¶
Plot up the covariance errors as a function of radius.
This method can be used to plot up the last covar results or a previously-stored set. Any error bars are shown at the scale they were calculated (as given by the
sigma
andpercent
columns of the results).- Parameters:
field (str) – The column to plot from the fit results (the match is case insensitive).
results (None or astropy.table.Table instance) – The return value from the
covar
orget_covar_results
methods.units (str or astropy.units.Unit, optional) – The X-axis units (a length or angle, such as ‘Mpc’ or ‘arcsec’, where the case is important).
xlog (bool, optional) – Should the x axis be drawn with a log scale (default True)?
ylog (bool, optional) – Should the y axis be drawn with a log scale (default False)?
overplot (bool, optional) – Clear the plot or add to existing plot?
clearwindow (bool, optional) – How does this interact with overplot?
See also
fit
,get_covar_results
,fit_plot
,conf_plot
,density_plot
,par_plot
Notes
Error bars are included on the dependent axis if the results contain columns that match the requested field with the suffix ‘_hi’. The error bars are therefore symmetric, which is different to
conf_plot
.If a limit is missing (i.e. it is a NaN) then no error bar is drawn. This can make it look like the error is very small.
Examples
Plot the temperature as a function of radius from the last fit, including error bars:
>>> dep.covar_plot('xsapec.kt')
Plot the density with the radii labelled in arcminutes and the density shown on a log scale:
>>> dep.covar_plot('density', units='arcmin', ylog=True)
Overplot the current covar results on those from a previous fit, where
covar1
was returned from thecovar
orget_covar_results
methods:>>> dep.covar_plot('xsapec.abundanc', results=covar1) >>> dep.covar_plot('xsapec.abundanc', overplot=True)
- density_plot(units='kpc', xlog=True, ylog=True, overplot=False, clearwindow=True)¶
Plot up the electron density as a function of radius.
The density is displayed with units of cm^-3. This plots up the density calculated using the current normalization parameter values. The
fit_plot
,conf_plot
, andcovar_plot
routines display the fit and error results for these parameters.- Parameters:
units (str or astropy.units.Unit, optional) – The X-axis units (a length or angle, such as ‘Mpc’ or ‘arcsec’, where the case is important).
xlog (bool, optional) – Should the x axis be drawn with a log scale (default True)?
ylog (bool, optional) – Should the y axis be drawn with a log scale (default False)?
overplot (bool, optional) – Clear the plot or add to existing plot?
clearwindow (bool, optional) – How does this interact with overplot?
See also
Examples
Plot the density as a function of radius.
>>> dep.density_plot()
Label the radii with units of arcminutes:
>>> dep.density_plot(units='arcmin')
- dummyfunc(*args, **kwargs)¶
- find_norm(shell)¶
Return the normalization value for the given shell.
This is limited to XSPEC-style models, where the parameter is called “norm”.
- Parameters:
shell (int) – The shell number.
- Returns:
norm – The normalization of the shell.
- Return type:
float
- Raises:
ValueError – If there is not one norm parameter for the shell.
See also
- find_parval(parname)¶
Return the value of the first parameter matching the name.
- Parameters:
parname (str) – The parameter name. The case is ignored in the match, and the first match is returned.
- Returns:
parval – The parameter value
- Return type:
float
- Raises:
ValueError – There is no match for the parameter.
Examples
>>> kt = dep.find_parval('kt')
- fit()¶
Fit the data using the “onion-peeling” method.
Unlike the normal Sherpa fit, this does not fit all the data simultaneously, but instead fits the outermost annulus first, then freezes its parameters and fits the annulus inside it, repeating this until all annuli have been fit. At the end of the fit all the parameters that were frozen are freed. The results can also be retrieved with
get_fit_results
.- Returns:
fits – This records per-annulus data, such as the inner and outer radius (rlo_ang, rhi_ang, rlo_phys, rhi_phys), the final fit statistic and change in fit statistic (statval and dstatval), the reduced statistic and q value (as rstat and qval) if appropriate, and the thawed parameter values (accessed using <model name>.<par name> syntax, where the match is case sensitive).
- Return type:
astropy.table.Table instance
See also
Notes
If there are any tied parameters between annuli then these annuli are combined together (fit simultaneously). The results from the fits to each annulus can be retrieved after
fit
has been called with theget_fit_results
method.The results have been separated out per annulus, even if several annuli were combined in a fit due to tied parameters, and there is no information in the returned structure to note this.
Examples
Fit the annuli using the onion-peeling approach, and then plot up the reduced statistic for each dataset:
>>> res = dep.fit() >>> plt.clf() >>> rmid = 0.5 * (res['rlo_phys'] + res['rhi_phys']) >>> plt.plot(rmid, res['rstat'])
Plot the temperature-abundance values per shell, color-coded by annulus:
>>> plt.clf() >>> plt.plot(res['xsapec.kT'], res['xsapec.Abundanc'], ... c=res['annulus']) >>> plt.colorbar() >>> plt.xlabel('kT') >>> plt.ylabel('Abundance')
Plot up the temperature distibution as a function of radius from the fit:
>>> dep.fit() >>> dep.fit_plot('xsmekal.kt')
- fit_plot(field, results=None, units='kpc', xlog=True, ylog=False, overplot=False, clearwindow=True)¶
Plot up the fit results as a function of radius.
This method can be used to plot up the last fit results or a previously-stored set. To include error bars on the dependent values use the conf_plot or covar_plot methods.
- Parameters:
field (str) – The column to plot from the fit results (the match is case insensitive).
results (None or astropy.table.Table instance) – The return value from the
fit
orget_fit_results
methods.units (str or astropy.units.Unit, optional) – The X-axis units (a length or angle, such as ‘Mpc’ or ‘arcsec’, where the case is important).
xlog (bool, optional) – Should the x axis be drawn with a log scale (default True)?
ylog (bool, optional) – Should the y axis be drawn with a log scale (default False)?
overplot (bool, optional) – Clear the plot or add to existing plot?
clearwindow (bool, optional) – How does this interact with overplot?
See also
fit
,get_fit_results
,conf_plot
,covar_plot
,density_plot
,par_plot
Examples
Plot the temperature as a function of radius from the last fit:
>>> dep.fit_plot('xsapec.kt')
Plot the reduced fit statistic from the last fit:
>>> dep.fit_plot('rstat')
Plot the density with the radii labelled in arcminutes and the density shown on a log scale:
>>> dep.fit_plot('density', units='arcmin', ylog=True)
Overplot the current fit results on those from a previous fit, where
fit1
was returned from thefit
orget_fit_results
methods:>>> dep.fit_plot('xsapec.abundanc', results=fit1) >>> dep.fit_plot('xsapec.abundanc', overplot=True)
- freeze(par)¶
Freeze the given parameter in each shell.
- Parameters:
par (str) – The parameter name, specified as <model_type>.<par_name>
Examples
>>> dep.freeze('clus.abundanc')
- get_conf_results()¶
What are the conf results, per annulus?
This returns the fit result for each annulus from the last time that the
conf
method was called. It does not check to see if anything has changed since the lastconf
call (e.g. parameters being tied together or untied, or a manual fit to a shell). Note thatget_shells
should be used to find out if the shells were grouped together (although this can be reconstructed from the datasets field of each ErrorEstResults instance).- Returns:
errors – This records per-annulus data, such as the inner and outer radius (rlo_ang, rhi_ang, rlo_phys, rhi_phys), the sigma and percent values, and parameter results (accessed using <model name>.<par name>, <model name>.<par name>_lo, and <model name>.<par name>_hi syntax, where the match is case sensitive).
- Return type:
astropy.table.Table instance
See also
fit
,get_covar_results
,get_fit_results
,get_radii
,get_shells
,conf_plot
- get_covar_results()¶
What are the covar results, per annulus?
This returns the fit result for each annulus from the last time that the
covar
method was called. It does not check to see if anything has changed since the lastcovar
call (e.g. parameters being tied together or untied, or a manual fit to a shell). Note thatget_shells
should be used to find out if the shells were grouped together.- Returns:
errors – This records per-annulus data, such as the inner and outer radius (rlo_ang, rhi_ang, rlo_phys, rhi_phys), the sigma and percent values, and parameter results (accessed using <model name>.<par name>, <model name>.<par name>_lo, and <model name>.<par name>_hi syntax, where the match is case sensitive).
- Return type:
astropy.table.Table instance
See also
fit
,get_conf_results
,get_fit_results
,get_radii
,get_shells
,covar_plot
- get_density()¶
Calculate the electron density for each shell.
Convert the model normalzations (assumed to match the standard definition for XSPEC thermal-plasma models) for each shell.
- Returns:
dens – The densities calculated for each shell, in units of cm^-3.
- Return type:
astropy.units.quantity.Quantity instance
See also
Notes
The electron density is taken to be:
n_e^2 = norm * 4*pi * DA^2 * 1e14 * (1+z)^2 / volume * ne_nh_ratio
where:
norm = model normalization from sherpa fit DA = angular size distance (cm) volume = volume (cm^3) ne_nh_ratio = 1.18
The model components for each volume element (the intersection of the annular cylinder
a
with the spherical shells
) are multiplied by a volume normalization:vol_norm[s,a] = volume[s,a] / v_sphere v_sphere = volume of sphere enclosing outer annulus
With this convention the
volume
used in calculating the electron density is simplyv_sphere
.
- get_fit_results()¶
What are the fit results, per annulus?
This returns the fit result for each annulus from the last time that the
fit
method was called. It does not check to see if anything has changed since the lastfit
call (e.g. parameters being tied together or untied, or a manual fit to a shell). Note thatget_shells
should be used to find out if the shells were grouped together.- Returns:
fits – This records per-annulus data, such as the inner and outer radius (rlo_ang, rhi_ang, rlo_phys, rhi_phys), the final fit statistic and change in fit statistic (statval and dstatval), the reduced statistic and q value (as rstat and qval) if appropriate, and the thawed parameter values (accessed using <model name>.<par name> syntax, where the match is case sensitive).
- Return type:
astropy.table.Table instance
See also
fit
,get_conf_results
,get_covar_results
,get_radii
,get_shells
,fit_plot
- get_par(par)¶
Return the parameter value for each shell.
- Parameters:
par (str) – The parameter name, specified as <model_type>.<par_name>
- Returns:
vals – The parameter values, in shell order.
- Return type:
ndarray
See also
Examples
>>> kts = dep.get_par('xsapec.kt')
- get_radii(units='arcsec')¶
What are the radii of the shells?
Return the inner and outer edge of each annulus, in the given units. Physical units (e.g. ‘kpc’) can only be used if a redshift or angular-diameter distance has been set. This does not apply the grouping that get_shells does.
- Parameters:
units (str or astropy.units.Unit, optional) – The name of the units to use for the returned radii. They must be an angle - such as ‘arcsec’ - or a length - such as ‘kpc’ or ‘Mpc’ (case is important).
See also
- Returns:
rlo, rhi – The inner and outer radius for each annulus.
- Return type:
astropy.units.Quantity, astropy.units.Quantity
- get_shells()¶
How are the annuli grouped?
An annulus may have multiple data sets associated with it, but it may also be linked to other annuli due to tied parameters. The return value is per group, in the ordering needed for the outside-to-inside onion skin fit, where the keys for the dictionary are ‘annuli’ and ‘dataids’.
- Returns:
groups – Each dictionary has the keys ‘annuli’ and ‘dataids’, and lists the annuli and data identifiers that are fit together. The ordering matches that of the onion-skin approach, so the outermost group first.
- Return type:
list of dicts
Examples
For a 3-annulus deprojection where there are no parameter ties to combine annului:
>>> dep.get_shells() [{'annuli': [2], 'dataids': [2]}, {'annuli': [1], 'dataids': [1]}, {'annuli': [0], 'dataids': [0]}]
After tie-ing the abundance parameter for the outer two shells, there are now two groups of annuli:
>>> dep.tie_par('xsapec.abundanc', 1, 2) Tying xsapec_2.Abundanc to xsapec_1.Abundanc >>> dep.get_shells() [{'annuli': [1, 2], 'dataids': [1, 2]}, {'annuli': [0], 'dataids': [0]}]
- guess()¶
Guess the starting point by fitting the projected data.
Use a fitting scheme - based on the suggestion in the XSPEC projct documention - to estimate the starting position of the fit (the initial fit parameters). This can be useful since it can reduce the time taken to fit the deprojected data and help avoid the deprojection from getting stuck in a local minimum.
See also
Notes
Each annulus, from outer to inner, is fit individually, ignoring the contribution from any outer annulus. After the fit, the model normalisation is corrected for the volume-filling factor of the annulus. If there are any tied parameters between annuli then these annuli are combined together (fit simultaneously).
Unlike the Sherpa guess function, this does not change the limits of any parameter.
- Possible improvements include:
re-normalize each spectrum before fitting.
transfer the model parameters of the inner-most shell in a group to the next set of shells to fit.
- ignore(*args)¶
Apply Sherpa ignore command to each dataset.
The filter is applied to each data set separately.
See also
Examples
Restrict the analysis to those bins which fall in the range 0.5 to 7.0 keV, where the limits are not included in the noticed range. The call to notice is used to clear any existing filter.
>>> dep.notice(None, None) >>> dep.ignore(None, 0.5) >>> dep.ignore(7.0, None)
- load_pha(specfile, annulus)¶
Load a pha file and add to the datasets for stacked analysis.
It is required that datasets for all annuli are loaded before the source model is created (to ensure that components are created for each annulus).
- Parameters:
specfile (str or sherpa.astro.data.DataPHA object) – If a string, the name of the file containing the source spectrum, which must be in PHA format (the data is expected to be extracted on the PI column). If a DataPHA object, then this is used (and is assumed to contain any needed background data).
annulus (int) – The annulus number for the data.
- Returns:
dataid – The Sherpa dataset identifier used for this spectrum.
- Return type:
int
Examples
Load the data for four annuli from the files ‘ann1.pi’ to ‘ann4.pi’.
>>> dep.load_pha('ann1.pi', 0) >>> dep.load_pha('ann2.pi', 1) >>> dep.load_pha('ann3.pi', 2) >>> dep.load_pha('ann4.pi', 3)
Load in the PHA files into Sherpa DataPHA objects, and then use these objects:
>>> s1 = ui.unpack_pha('src1.pi') >>> s2 = ui.unpack_pha('src2.pi') >>> s3 = ui.unpack_pha('src3.pi') >>> dep.load_pha(s1, 0) >>> dep.load_pha(s2, 1) >>> dep.load_pha(s3, 2)
- notice(*args)¶
Apply Sherpa notice command to each dataset.
The filter is applied to each data set separately.
See also
Examples
Restrict the analysis to those bins which fall in the range 0.5 to 7.0 keV, where the limits are included in the noticed range. The first call to notice is used to clear any existing filter.
>>> dep.notice(None, None) >>> dep.notice(0.5, 7.0)
- par_plot(par, units='kpc', xlog=True, ylog=False, overplot=False, clearwindow=True)¶
Plot up the parameter as a function of radius.
This plots up the current parameter values. The
fit_plot
,conf_plot
, andcovar_plot
routines display the fit and error results for these parameters.- Parameters:
par (str) – The parameter name, specified as <model_type>.<par_name>.
units (str or astropy.units.Unit, optional) – The X-axis units (a length or angle, such as ‘Mpc’ or ‘arcsec’, where the case is important).
xlog (bool, optional) – Should the x axis be drawn with a log scale (default True)?
ylog (bool, optional) – Should the y axis be drawn with a log scale (default False)?
overplot (bool, optional) – Clear the plot or add to existing plot?
clearwindow (bool, optional) – How does this interact with overplot?
See also
Examples
Plot the temperature as a function of radius.
>>> dep.par_plot('xsapec.kt')
Label the radii with units of arcminutes for the abundanc parameter of the xsapec model:
>>> dep.par_plot('xsapec.abundanc', units='arcmin')
- plot_arf(*args, **kwargs)¶
- plot_bkg(*args, **kwargs)¶
- plot_bkg_chisqr(*args, **kwargs)¶
- plot_bkg_delchi(*args, **kwargs)¶
- plot_bkg_fit(*args, **kwargs)¶
- plot_bkg_fit_delchi(*args, **kwargs)¶
- plot_bkg_fit_resid(*args, **kwargs)¶
- plot_bkg_model(*args, **kwargs)¶
- plot_bkg_ratio(*args, **kwargs)¶
- plot_bkg_resid(*args, **kwargs)¶
- plot_bkg_source(*args, **kwargs)¶
- plot_bkg_unconvolved(*args, **kwargs)¶
- plot_chisqr(*args, **kwargs)¶
- plot_data(*args, **kwargs)¶
- plot_delchi(*args, **kwargs)¶
- plot_fit(*args, **kwargs)¶
- plot_fit_delchi(*args, **kwargs)¶
- plot_fit_resid(*args, **kwargs)¶
- plot_model(*args, **kwargs)¶
- plot_order(*args, **kwargs)¶
- plot_psf(*args, **kwargs)¶
- plot_ratio(*args, **kwargs)¶
- plot_resid(*args, **kwargs)¶
- plot_source(*args, **kwargs)¶
- print_window(*args, **kwargs)¶
Create a hardcopy version of each plot window.
- Parameters:
args – The arguments to be sent to the “create a hardcopy” routine. The first argument, if given, is assumed to be the file name and so will have the shell number added to it.
kwargs – Keyword arguments for the call.
Notes
This is not guaranteed to work properly for Matplotlib.
Examples
Create hardcopy versions of the data plots, called “data0”, “data1”, …
>>> dep.plot_data() >>> dep.print_window('data')
- set_bkg_model(bkgmodel)¶
Create a background model for each annulus.
The background model is the same between the annuli, except that a scaling factor is added for each annulus (to allow for normalization uncertainities). The scaling factors are labelled ‘bkg_norm_<obsid>’, and at least one of these must be frozen (otherwise it is likely to be degenerate with the background normalization, causing difficulties for the optimiser).
- Parameters:
bkgmodel (model instance) – The background model expression applied to each annulus. Unlike set_source this should be the actual model instance, and not a string.
See also
Examples
Model the background with a single power-law component:
>>> dep.set_bkg_model(xspowerlaw.bpl)
- set_par(par, val)¶
Set the parameter value in each shell.
- Parameters:
par (str) – The parameter name, specified as <model_type>.<par_name>
val (float) – The parameter value.
Examples
>>> dep.set_par('xsapec.abundanc', 0.25)
- set_source(srcmodel='xsphabs*xsapec')¶
Create a source model for each annulus.
Unlike the standard set_source command, this version just uses the <model name>, not <model name>.<username>, since the <username> is automatically created for users by appending the annulus number to <model name>.
- Parameters:
srcmodel (str, optional) – The source model expression applied to each annulus.
See also
Notes
The data must have been read in for all the data before calling this method (this matches Sherpa, where you can not call set_source unless you have already loaded the data to fit).
Examples
The following two calls have the same result: model instances called ‘xsphabs<annulus>’ and ‘xsapec<annulus>’ are created for each annulus, and the source expression for the annulus set to their multiplication:
>>> dep.set_source() >>> dep.set_source('xsphabs * xsapec')
Use the XSPEC vapec model rather than the apec model to represent the plasma emission:
>>> dep.set_source('xsphabs * xsvapec')
- subtract(*args)¶
Subtract the background from each dataset.
See also
Examples
>>> dep.substract()
- thaw(par)¶
Thaw the given parameter in each shell.
- Parameters:
par (str) – The parameter name, specified as <model_type>.<par_name>
Examples
>>> dep.thaw('clus.abundanc')
- tie_par(par, base, *others)¶
Tie parameters in one or more shells to the base shell.
This is a limited form of the Sherpa ability to link parameters, since it sets the parameter in the other shells to the same value as the parameter in the base shell. More complex situations will require direct calls to sherpa.astro.ui.link.
- Parameters:
par (str) – The parameter specifier, as <model_type>.<par_name>.
base (int) – The base shell number.
*others (scalar) – The shell, or shells to link to the base shell.
Examples
Tie the temperature and abundance parameters in shell 9 to that in shell 8, so that any fits will set the shell 9 values to those used in shell 8 (so reducing the number of free parameters in the fit).
>>> dep.tie_par('xsapec.kt', 8, 9) Tying xsapec_9.kT to xsapec_8.kT >>> dep.tie_par('xsapec.abundanc', 8, 9) Tying xsapec_9.Abundanc to xsapec_8.Abundanc
Tie three annuli together:
>>> dep.tie_par('xsapec.kt', 12, 13, 14) Tying xsapec_13.kT to xsapec_12.kT Tying xsapec_14.kT to xsapec_12.kT
- unsubtract(*args)¶
Un-subtract the background from each dataset.
This can be useful when you want to compare the results to the “wstat” stat (a Poisson-based stat which includes the background data as a component and provides a goodness-of-fit estimate).
See also
Examples
>>> dep.unsubstract()
- untie_par(par, *others)¶
Remove the parameter tie/link in the shell.
This is intended to remove links between shells created by tie_par, but will remove any links created by sherpa.astro.ui.link.
- Parameters:
par (str) – The parameter specifier, as <model_type>.<par_name>.
*others (scalar) – The shell, or shells to un-tie/unlink.
See also
Notes
It is safe to call on a parameter that is not tied or linked to another parameter.
Examples
Untie the abundance parameter in shell 9; that is, it is now free to vary independently in a fit.
>>> dep.untie_par('xsapec.abundanc', 9) Untying xsapec_9.Abundanc
Untie multiple annuli:
>>> dep.untie_par('xsmekal.kt', 13, 14) Untying xsmekal_13.kT Untying xsmekal_14.kT