pyproffit package

Submodules

pyproffit.data module

class pyproffit.data.Data(imglink, explink=None, bkglink=None, voronoi=False, rmsmap=None)[source]

Bases: object

Class containing the data to be loaded and used by other pyproffit routines

Parameters:
  • imglink (str) – Path to input image
  • explink (str , optional) – Path to exposure map. If none, assume a flat exposure of 1s or an input error map provided through rmsmap
  • bkglink (str , optional) – Path to background map. If none, assume zero background
  • voronoi (bool , optional) – Define whether the input image is a Voronoi image or not (default=False)
  • rmsmap (str , optional) – Path to error map if the data is not Poisson distributed
dmfilth(outfile=None, smoothing_scale=8)[source]

Mask the regions provided in a region file and fill in the holes by interpolating the smoothed image into the gaps and generating a Poisson realization

Parameters:
  • outfile (str , optional) – If outfile is not None, file name to output the dmfilth image into a FITS file
  • smoothing_scale (int) – Size of smoothing scale (in pixel) to estimate the surface brightness distribution outside of the masked areas
region(regfile)[source]

Filter out regions provided in an input DS9 region file

Parameters:regfile (str) – Path to region file. Accepted region file formats are fk5 and image.
reset_exposure()[source]

Revert to the original exposure map and ignore the current region file

pyproffit.data.get_extnum(fitsfile)[source]

Find the extension number of the first IMAGE extension in an input FITS file

Parameters:fitsfile (str) – Input FITS file to be read
Returns:extension number
Return type:int

pyproffit.deproject module

class pyproffit.deproject.Deproject(z=None, profile=None, cf=None, f_abund='aspl')[source]

Bases: object

Class to perform all calculations of deprojection, density profile, gas mass, count rate, and luminosity

Parameters:
  • z (float) – Source redshift. If z=None, only the surface brightness reconstruction can be done.
  • profile (class:pyproffit.profextract.Profile) – Object of type pyproffit.profextract.Profile containing the surface brightness profile data
  • cf (float) – Conversion factor from count rate to emissivity. If cf=None, only the surface brightness reconstruction can be done.
  • f_abund (str) – Solar abundance table to compute the electron-to-proton ratio and mean molecular weight. Available tables are ‘aspl’ (Aspling+09, default), ‘angr’ (Anders & Grevesse 89), and ‘grsa’ (Grevesse & Sauval 98)
CSB(rin=40.0, rout=400.0, plot=True, outfile=None, figsize=(13, 10), nbins=30, fontsize=40, yscale='linear', **kwargs)[source]

Compute the surface brightness concentration from a loaded brightness profile reconstruction. The surface brightness concentration is defined as the ratio of fluxes computed within two apertures.

Parameters:
  • rin (float) – Lower aperture value in kpc (default=40)
  • rout (float) – Higher aperture value in kpc (default=400)
  • plot (bool) – Plot the posterior CSB distribution (default=True)
  • outfile (str) – Output file name to save the figure. If outfile=None, plot only to stdout
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • nbins (int , optional) – Number of bins on the X axis to construct the posterior distribution. Defaults to 30
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • yscale (str , optional) – Scale on the Y axis. Defaults to ‘linear’
  • kwargs – Options to be passed to matplotplib.pyplot.hist
Returns:

Median count rate, 16th and 84th percentiles

Return type:

float

CountRate(a, b, plot=True, outfile=None, figsize=(13, 10), nbins=30, fontsize=40, yscale='linear', **kwargs)[source]

Compute the model count rate integrated between radii a and b. Optionally, the count rate distribution can be plotted and saved.

Parameters:
  • a (float) – Inner integration boundary in arcmin
  • b (float) – Outer integration boundary in arcmin
  • plot (bool) – Plot the posterior count rate distribution (default=True)
  • outfile (str) – Output file name to save the figure. If outfile=None, plot only to stdout
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • nbins (int , optional) – Number of bins on the X axis to construct the posterior distribution. Defaults to 30
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • yscale (str , optional) – Scale on the Y axis. Defaults to ‘linear’
  • kwargs – Options to be passed to matplotplib.pyplot.hist
Returns:

Median count rate, 16th and 84th percentiles

Return type:

float

Density(rout=None)[source]

Compute a density profile from a multiscale reconstruction

Parameters:rout (numpy.ndarray) – Radial binning of the density profile. If rout=None, the original binning of the surface brightness profile is used
Luminosity(a, b, plot=True, outfile=None, figsize=(13, 10), nbins=30, fontsize=40, yscale='linear', **kwargs)[source]

Compute the luminosity integrated between radii a and b. Optionally, the luminosity distribution can be plotted and saved. Requires the luminosity factor to be computed using the pyproffit.profextract.Profile.Emissivity() method.

Parameters:
  • a (float) – Inner integration boundary in arcmin
  • b (float) – Outer integration boundary in arcmin
  • plot (bool) – Plot the posterior count rate distribution (default=True)
  • outfile (str) – Output file name to save the figure. If outfile=None, plot only to stdout
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • nbins (int , optional) – Number of bins on the X axis to construct the posterior distribution. Defaults to 30
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • yscale (str , optional) – Scale on the Y axis. Defaults to ‘linear’
  • kwargs – Options to be passed to matplotplib.pyplot.hist
Returns:

Median luminosity, 16th and 84th percentiles

Return type:

float

Mgas(radius, radius_err=None, plot=True, outfile=None, rad_scale='normal', figsize=(13, 10), nbins=30, fontsize=40, yscale='linear', **kwargs)[source]

Compute the posterior cumulative gas mass within a given radius. Optionally, the posterior distribution can be plotted and saved.

Parameters:
  • radius (float) – Gas mass integration radius in kpc
  • radius_err (float , optional) – (Gaussian) error on the input radius to be propagated to the gas mass measurement. To be used in case one wants to evaluate \(M_{gas}\) at an overdensity radius with a known uncertainty
  • plot (bool) – Plot the posterior Mgas distribution (default=True)
  • outfile (str , optional) – Output file name to save the figure. If outfile=None, plot only to stdout
  • rad_scale (str) – If radius_err is not None, specify whether the distribution of radii is drawn from a normal distribution (rad_scale=’normal’) or a log-normal distribution (rad_scale=’lognormal’). Defaults to ‘normal’.
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • nbins (int , optional) – Number of bins on the X axis to construct the posterior distribution. Defaults to 30
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • yscale (str , optional) – Scale on the Y axis. Defaults to ‘linear’
  • kwargs – Options to be passed to matplotplib.pyplot.hist
Returns:

Median \(M_{gas}\), 16th and 84th percentiles

Return type:

float

Multiscale(backend='pymc3', nmcmc=1000, tune=500, bkglim=None, back=None, samplefile=None, nrc=None, nbetas=6, depth=10, min_beta=0.6)[source]

Run Multiscale deprojection using the method described in Eckert+20

Parameters:
  • backend (str) – Backend to run the optimization problem. Available backends are ‘pymc3’ (default) and ‘stan’
  • nmcmc (int) – Number of HMC points in the output sample
  • tune (int) – Number of HMC tuning steps
  • bkglim (float) – Limit beyond which it is assumed that the background dominates, i.e. the source is set to 0. If bkglim=None (default), the entire radial range is used
  • back (float) – Input value for the background, around which a gaussian prior is set. If back=None (default), the input background value will be computed as the average of the source-free region
  • samplefile (str) – Path to output file to write the output samples. If samplefile=None (default), the data are not written to file and only loaded into memory
  • nrc (int) – Number of core radii. If nrc=None (default), the number of core radiis will be defined on-the-fly
  • nbetas (int) – Number of beta values. Default=6
  • depth (int) – Set the max_treedepth parameter for Stan (default=10)
  • min_beta (float) – Minimum value of beta. Default=0.6
Ncounts(plot=True, outfile=None, figsize=(13, 10), nbins=30, fontsize=40, yscale='linear', **kwargs)[source]

Compute the total model number of counts. Optionally, the posterior distribution can be plotted and saved.

Parameters:
  • plot (bool) – Plot the posterior distribution of number of counts (default=True)
  • outfile (str) – Output file name to save the figure. If outfile=None, plot only to stdout
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • nbins (int , optional) – Number of bins on the X axis to construct the posterior distribution. Defaults to 30
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • yscale (str , optional) – Scale on the Y axis. Defaults to ‘linear’
  • kwargs – Options to be passed to matplotplib.pyplot.hist
Returns:

Median number of counts, 16th and 84th percentiles

Return type:

float

OnionPeeling(nmc=1000)[source]

Run standard Onion Peeling deprojection using the Kriss+83 method and the McLaughlin+99 edge correction

Parameters:nmc (int) – Number of Monte Carlo realizations of the profile to compute the uncertainties (default=1000)
PlotDensity(outfile=None, xunit='kpc', figsize=(13, 10), fontsize=40, color='C0', lw=2, **kwargs)[source]

Plot the loaded density profile

Parameters:
  • outfile (str) – Output file name. If outfile=None (default) plot only to stdout
  • xunit (str) – Choose whether the x axis should be in unit of ‘kpc’ (default), ‘arcmin’, or ‘both’, in which case two axes are drawn at the top and the bottom of the plot
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • color (str , optional) – Line color following matplotlib conventions. Defaults to ‘blue’
  • lw (int , optional) – Line width. Defaults to 2
  • kwargs – Additional arguments to be passed to matplotlib.pyplot.plot
PlotMgas(rout=None, outfile=None, xunit='kpc', figsize=(13, 10), color='C0', lw=2, fontsize=40, xscale='log', yscale='log')[source]

Plot the cumulative gas mass profile from the output of a reconstruction

Parameters:
  • rout (numpy.ndarray) – Radial binning of the gas mass profile. If rout=None, the original binning of the surface brightness profile is used
  • outfile (str) – Output file name to save the figure. If outfile=None, plot only to stdout
  • xunit (str) – Choose whether the x axis should be in unit of ‘kpc’ (default), ‘arcmin’, or ‘both’, in which case two axes are drawn at the top and the bottom of the plot
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • color (str , optional) – Line color following matplotlib conventions. Defaults to ‘C0’
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • xscale (str , optional) – Scale of the X axis. Defaults to ‘log’
  • yscale (str , optional) – Scale of the Y axis. Defaults to ‘log’
  • lw (int , optional) – Line width. Defaults to 2
PlotSB(outfile=None, figsize=(13, 10), fontsize=40, xscale='log', yscale='log', lw=2, fmt='d', markersize=7, data_color='red', bkg_color='green', model_color='C0', skybkg_color='k')[source]

Plot the surface brightness profile reconstructed after applying the multiscale deprojection and PSF deconvolution technique, and compare it with the input brightness profile

Parameters:
  • outfile (str) – File name of saving the output figure. If outfile=None (default), plot only to stdout
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • xscale (str , optional) – Scale of the X axis. Defaults to ‘log’
  • yscale (str , optional) – Scale of the Y axis. Defaults to ‘log’
  • lw (int , optional) – Line width. Defaults to 2
  • fmt (str , optional) – Marker type following matplotlib convention. Defaults to ‘d’
  • markersize (int , optional) – Marker size. Defaults to 7
  • data_color (str , optional) – Color of the data points following matplotlib convention. Defaults to ‘red’
  • bkg_color (str , optional) – Color of the particle background following matplotlib convention. Defaults to ‘green’
  • model_color (str , optional) – Color of the surface brightness model following matplotlib convention. Defaults to ‘C0’
  • skybkg_color (str , optional) – Color of the fitted sky background model following matplotlib convention. Defaults to ‘k’
Reload(samplefile, bkglim=None)[source]

Reload the samples stored from a previous reconstruction run

Parameters:
  • samplefile (str) – Path to file containing the saved HMC samples
  • bkglim (float) – Limit beyond which it is assumed that the background dominates, i.e. the source is set to 0. This parameter needs to be the same as the value used to run the reconstruction. If bkglim=None (default), the entire radial range is used
SaveAll(outfile=None)[source]

Save the results of a profile reconstruction into an output FITS file First extension is data Second extension is density Third extension is Mgas Fourth extension is PSF

Parameters:outfile (str) – Output file name
pyproffit.deproject.Deproject_Multiscale_PyMC3(deproj, bkglim=None, nmcmc=1000, tune=500, back=None, samplefile=None, nrc=None, nbetas=6, min_beta=0.6)[source]

Run the multiscale deprojection optimization using the PyMC3 backend

Parameters:
  • deproj (class:pyproffit.deproject.Deproject) – Object of type pyproffit.deproject.Deproject containing the data and parameters
  • bkglim (float) – Limit beyond which it is assumed that the background dominates, i.e. the source is set to 0. If bkglim=None (default), the entire radial range is used
  • nmcmc (int) – Number of HMC points in the output sample
  • back (float) – Input value for the background, around which a gaussian prior is set. If back=None (default), the input background value will be computed as the average of the source-free region
  • samplefile (str) – Path to output file to write the output samples. If samplefile=None (default), the data are not written to file and only loaded into memory
  • nrc (int) – Number of core radii. If nrc=None (default), the number of core radiis will be defined on-the-fly
  • nbetas (int) – Number of beta values. Default=6
  • min_beta (float) – Minimum value of beta. Default=0.6
pyproffit.deproject.Deproject_Multiscale_Stan(deproj, bkglim=None, nmcmc=1000, back=None, samplefile=None, nrc=None, nbetas=6, depth=10, min_beta=0.6)[source]

Run the multiscale deprojection optimization using the Stan backend

Parameters:
  • deproj (class:pyproffit.deproject.Deproject) – Object of type pyproffit.deproject.Deproject containing the data and parameters
  • bkglim (float) – Limit beyond which it is assumed that the background dominates, i.e. the source is set to 0. If bkglim=None (default), the entire radial range is used
  • nmcmc (int) – Number of HMC points in the output sample
  • back (float) – Input value for the background, around which a gaussian prior is set. If back=None (default), the input background value will be computed as the average of the source-free region
  • samplefile (str) – Path to output file to write the output samples. If samplefile=None (default), the data are not written to file and only loaded into memory
  • nrc (int) – Number of core radii. If nrc=None (default), the number of core radiis will be defined on-the-fly
  • nbetas (int) – Number of beta values. Default=6
  • depth (int) – Set the max_treedepth parameter for Stan (default=10)
  • min_beta (float) – Minimum value of beta. Default=0.6
pyproffit.deproject.EdgeCorr(nbin, rin_cm, rout_cm, em0)[source]

For Onion Peeling deprojection, correct for edge effects by estimating the contribution of flux beyond the edge, using McLaughlin+99 method

Parameters:
  • nbin (int) – Number of bins in input profile
  • rin_cm (numpy.ndarray) – Inner radii of the bins in cm
  • rout_cm (numpy.ndarray) – Outer radii of the bins in cm
  • em0 (numpy.ndarray) – Deprojected emissivity profile before edge correction
Returns:

Edge-corrected deprojected profile

Return type:

numpy.ndarray

class pyproffit.deproject.MyDeprojVol(radin, radout)[source]

Bases: object

Class to compute the projection volumes

Parameters:
  • radin (class:numpy.ndarray) – Array of inner radii of the bins
  • radout (class:numpy.ndarray) – Array of outer radii of the bins
deproj_vol()[source]

Compute the projection volumes

Returns:Volume matrix
Return type:numpy.ndarray
pyproffit.deproject.OP(deproj, nmc=1000)[source]

Run standard Onion Peeling deprojection including edge correction

Parameters:
  • deproj (class:pyproffit.deproject.Deproject) – Object of type pyproffit.deproject.Deproject containing the input data
  • nmc (int) – Number of Monte Carlo realizations of the profile to compute the uncertainties (default=1000)
pyproffit.deproject.calc_density_operator(rad, pars, z, cosmo)[source]

Compute linear operator to transform parameters into gas density profiles

Parameters:
  • rad (numpy.ndarray) – Array of input radii in arcmin
  • pars (numpy.ndarray) – List of beta model parameters obtained through list_params
  • z (float) – Source redshift
Returns:

Linear operator for gas density

Return type:

numpy.ndarray

pyproffit.deproject.calc_int_operator(a, b, pars)[source]

Compute a linear operator to integrate analytically the basis functions within some radial range and return count rate and luminosities

Parameters:
  • a (float) – Lower integration boundary
  • b (float) – Upper integration boundary
  • pars (numpy.ndarray) – List of beta model parameters obtained through list_params
Returns:

Linear integration operator

Return type:

numpy.ndarray

pyproffit.deproject.calc_linear_operator(rad, sourcereg, pars, area, expo, psf)[source]

Function to calculate a linear operator transforming parameter vector into predicted model counts

Parameters:
  • rad (numpy.ndarray) – Array of input radii in arcmin
  • sourcereg (numpy.ndarray) – Selection array for the source region
  • pars (numpy.ndarray) – List of beta model parameters obtained through list_params
  • area (numpy.ndarray) – Bin area in arcmin^2
  • expo (numpy.ndarray) – Bin effective exposure in s
  • psf (numpy.ndarray) – PSF mixing matrix
Returns:

Linear projection and PSF mixing operator

Return type:

numpy.ndarray

pyproffit.deproject.calc_sb_operator(rad, sourcereg, pars)[source]

Function to calculate a linear operator transforming parameter vector into surface brightness

Parameters:
  • rad (numpy.ndarray) – Array of input radii in arcmin
  • sourcereg (numpy.ndarray) – Selection array for the source region
  • pars (numpy.ndarray) – List of beta model parameters obtained through list_params
Returns:

Linear projection operator

Return type:

numpy.ndarray

pyproffit.deproject.calc_sb_operator_psf(rad, sourcereg, pars, area, expo, psf)[source]

Same as calc_sb_operator but convolving the model surface brightness with the PSF model

Parameters:
  • rad (numpy.ndarray) – Array of input radii in arcmin
  • sourcereg (numpy.ndarray) – Selection array for the source region
  • pars (numpy.ndarray) – List of beta model parameters obtained through list_params
  • area (numpy.ndarray) – Bin area in arcmin^2
  • expo (numpy.ndarray) – Bin effective exposure in s
  • psf (numpy.ndarray) – PSF mixing matrix
Returns:

Linear projection and PSF mixing operator

Return type:

numpy.ndarray

pyproffit.deproject.fbul19(R, z, cosmo, Runit='kpc')[source]

Compute Mgas from input R500 using Bulbul+19 M-Mgas scaling relation

Parameters:
  • R (float) – Input R500 value
  • z (float) – Input redshift
  • Runit (str) – Unit of input radis, kpc or arcmin (default=’kpc’)
Returns:

Mgas

Return type:

float

pyproffit.deproject.list_params(rad, sourcereg, nrc=None, nbetas=6, min_beta=0.6)[source]

Define a list of parameters to define the dictionary of basis functions

Parameters:
  • rad (numpy.ndarray) – Array of input radii in arcmin
  • sourcereg (numpy.ndarray) – Selection array for the source region
  • nrc (int) – Number of core radii. If nrc=None (default), the number of core radiis will be defined on-the-fly
  • nbetas (int) – Number of beta values. Default=6
  • min_beta (float) – Minimum value of beta. Default=0.6
Returns:

Array containing sets of values to set up the function dictionary

Return type:

numpy.ndarray

pyproffit.deproject.list_params_density(rad, sourcereg, z, cosmo, nrc=None, nbetas=6, min_beta=0.6)[source]

Define a list of parameters to transform the basis functions into gas density profiles

Parameters:
  • rad (numpy.ndarray) – Array of input radii in arcmin
  • sourcereg (numpy.ndarray) – Selection array for the source region
  • z (float) – Source redshift
  • nrc (int) – Number of core radii. If nrc=None (default), the number of core radiis will be defined on-the-fly
  • nbetas (int) – Number of beta values. Default=6
  • min_beta (float) – Minimum value of beta. Default=0.6
Returns:

Array containing sets of values to set up the function dictionary

Return type:

numpy.ndarray

pyproffit.deproject.medsmooth(prof)[source]

Smooth a given profile by taking the median value of surrounding points instead of the initial value

Parameters:prof (numpy.ndarray) – Input profile to be smoothed
Returns:Smoothd profile
Return type:numpy.ndarray
pyproffit.deproject.plot_multi_methods(profs, deps, labels=None, outfile=None, xunit='kpc', figsize=(13, 10), fontsize=40, xscale='log', yscale='log', fmt='.', markersize=7)[source]

Plot multiple gas density profiles (e.g. obtained through several methods, centers or sectors) to compare them

Parameters:
  • profs (tuple) – List of Profile objects to be plotted
  • deps (tuple) – List of Deproject objects to be plotted
  • labels (tuple) – List of labels for the legend (default=None)
  • outfile (str) – If outfile is not None, path to file name to output the plot
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • xscale (str , optional) – Scale of the X axis. Defaults to ‘log’
  • yscale (str , optional) – Scale of the Y axis. Defaults to ‘log’
  • fmt (str , optional) – Marker type following matplotlib convention. Defaults to ‘d’
  • markersize (int , optional) – Marker size. Defaults to 7

pyproffit.emissivity module

pyproffit.emissivity.calc_emissivity(cosmo, z, nh, kt, rmf, Z=0.3, elow=0.5, ehigh=2.0, arf=None, type='cr', lum_elow=0.5, lum_ehigh=2.0, abund='aspl')[source]

Function calc_emissivity. The function computes the scaling factor between count rate and APEC/MEKAL norm using XSPEC, which is needed to extract density profiles. Requires XSPEC to be in PATH

Parameters:
  • cosmo (class: astropy.cosmology) – Astropy cosmology object
  • z (float) – Source redshift
  • nh (float) – Source NH in units of 1e22 cm**(-2)
  • kt (float) – Source temperature in keV
  • rmf (str) – Path to response file (RMF/RSP)
  • Z (float) – Metallicity with respect to solar (default = 0.3)
  • elow (float) – Low-energy bound of the input image in keV (default = 0.5)
  • ehigh (float) – High-energy bound of the input image in keV (default = 2.0)
  • arf (str) – Path to on-axis ARF (optional, in case response file is RMF)
  • type (str) – Specify whether the exposure map is in units of sec (type=’cr’) or photon flux (type=’photon’). By default type=’cr’.
  • lum_elow (float) – Low energy bound (rest frame) for luminosity calculation. Defaults to 0.5
  • lum_ehigh (float) – High energy bound (rest frame) for luminosity calculation. Defaults to 2.0
Returns:

Conversion factor

Return type:

float

pyproffit.emissivity.is_tool(name)[source]

Check whether name is on PATH.

pyproffit.fitting module

class pyproffit.fitting.ChiSquared(model, x, dx, y, dy, psfmat=None, fitlow=None, fithigh=None)[source]

Bases: object

Class defining a chi-square likelihood based on a surface brightness profile and a model. Let \(S_i\) be the measured surface brightness in annulus i and \(\sigma_i\) the corresponding Gaussian error. The likelihood function to be optimized is

\[-2\log \mathcal{L} = \sum_{i=1}^N \frac{(S_i - f(r_i))^2}{\sigma_i^2}\]

This class is called by the Fitter object when using the method=’chi2’ option.

Parameters:
  • model (class:pyproffit.models.Model) – Model definition. A pyproffit.models.Model object defining the model to be used.
  • x (class:numpy.ndarray) – x axis data
  • dx (class:numpy.ndarray) – x bin size data. dx is defined as half of the total bin size.
  • y (class:numpy.ndarray) – y axis data
  • dy (class:numpy.ndarray) – y error data
  • psfmat (class:numpy.ndarray , optional) – PSF convolution matrix
  • fitlow (float , optional) – Lower fitting boundary in arcmin. If fitlow=None the entire radial range is used, default to None
  • fithigh (float , optional) – Upper fitting boundary in arcmin. If fithigh=None the entire radial range is used, default to None
errordef = 1.0
class pyproffit.fitting.Cstat(model, x, dx, counts, area, effexp, bkgc, psfmat=None, fitlow=None, fithigh=None)[source]

Bases: object

Class defining a C-stat likelihood based on a surface brightness profile and a model. Let \(A_i\) , \(T_i\) be the area and the effective exposure time of annulus i. We set \(F_{i} = f(r_{i})A_{i}T_{i}\) the predicted number of counts in the annulus. The Poisson likelihood is then given by

\[-2\log \mathcal{L} = 2 \sum_{i=1}^N F_i - C_i \log F_i - C_i + C_i \log C_i\]

with \(C_i\) the observed number of counts in annulus i.

This class is called by the Fitter object when using the method=’cstat’ option.

Parameters:
  • model (class:pyproffit.models.Model) – Model definition. A pyproffit.models.Model object defining the model to be used.
  • x (numpy.ndarray) – x axis data
  • counts (numpy.ndarray) – counts per bin data
  • area (numpy.ndarray) – bin are in arcmin^2
  • effexp (numpy.ndarray) – bin effective exposure in s
  • bkgc (numpy.ndarray) – number of background counts per bin
  • psfmat (numpy.ndarray) – PSF convolution matrix
  • fitlow (float) – Lower fitting boundary in arcmin. If fitlow=None (default) the entire radial range is used
  • fithigh (float) – Upper fitting boundary in arcmin. If fithigh=None (default) the entire radial range is used
errordef = 1.0
class pyproffit.fitting.Fitter(model, profile, method='chi2', fitlow=None, fithigh=None, **kwargs)[source]

Bases: object

Class containing the tools to fit surface brightness profiles with a model. Sets up the likelihood and optimizes for the parameters.

Parameters:
  • model (class:pyproffit.models.Model) – Object of type pyproffit.models.Model defining the model to be used.
  • profile (class:pyproffit.profextract.Profile) – Object of type pyproffit.profextract.Profile containing the surface brightness profile to be fitted
  • method (str) – Likelihood function to be optimized. Available likelihoods are ‘chi2’ (chi-squared) and ‘cstat’ (C statistic). Defaults to ‘chi2’.
  • fitlow (float) – Lower boundary of the active fitting radial range. If fitlow=None the entire range is used. Defaults to None
  • fithigh (float) – Upper boundary of the active fitting radial range. If fithigh=None the entire range is used. Defaults to None
  • kwargs – List of arguments to be passed to the iminuit library. For instance, setting parameter boundaries, optimization options or fixing parameters. See the iminuit documentation: https://iminuit.readthedocs.io/en/stable/index.html
Corner(labels=None, **kwargs)[source]

Produce a parameter corner plot from a loaded set of samples. Uses the corner library: https://corner.readthedocs.io/en/latest/

Parameters:
Returns:

Output matplotlib figure

Emcee(nmcmc=5000, burnin=100, start=None, prior=None, walkers=32, thin=15)[source]

Run a Markov Chain Monte Carlo optimization using the affine-invariant ensemble sampler emcee. See https://emcee.readthedocs.io/en/stable/ for details.

Parameters:
  • nmcmc (int) – Number of MCMC samples. Defaults to 5000
  • burnin (int) – Size of the burn-in phase that will eventually be ignored. Defaults to 100
  • start (class:numpy.ndarray) – Array of input parameter values. If None, the code will look for the results of a previous Migrad optimization and use the corresponding parameters as starting values. Defaults to None
  • prior (function) – Function defining the priors on the parameters. The function should take the parameter set as input and return the log prior probability. If None, the code will search for the results of a previous Migrad optimization and set up broad Gaussian priors on each parameter with sigma set to 5 times the Migrad errors. Defaults to None.
  • walkers (int) – Number of emcee walkers. Defaults to 32.
  • thin (int) – Thinning number for the output samples. The total number of sample values will be nmcmc*walkers/thin. Defaults to 15.
Migrad(fixed=None)[source]

Perform maximum-likelihood optimization of the model using the MIGRAD routine from the MINUIT library.

Parameters:fixed (class:numpy.ndarray) – A boolean array setting up whether parameters are fixed (True) or left free (False) while fitting. If None, all parameters are fitted. Defaults to None.

pyproffit.hmc module

pyproffit.hmc.BetaModelPM(x, beta, rc, norm, bkg)[source]
pyproffit.hmc.BknPowPM(x, alpha1, alpha2, norm, jump, bkg, rf=3.0)[source]

Broken power law 3D model projected along the line of sight for discontinuity modeling

\[I(r) = I_0 \int F(\omega)^2 d\ell + B\]

with \(\omega^2 = r^2 + \ell^2\) and

\[\begin{split}F(\omega) = \left\{ \begin{array}{ll} \omega^{-\alpha_1}, & \omega<r_f \\ \frac{1}{C}\omega ^{-\alpha_2}, & \omega\geq r_f \end{array} \right.\end{split}\]
Parameters:
  • x (numpy.ndarray) – Radius in arcmin
  • alpha1 (class:theano.tensor) – \(\alpha_1\) parameter
  • alpha2 (class:theano.tensor) – \(\alpha_2\) parameter
  • rf (float) – rf parameter
  • norm (class:theano.tensor) – log of I0 parameter
  • jump (class:theano.tensor) – C parameter
  • bkg (class:theano.tensor) – log of B parameter
Returns:

Calculated model

Return type:

numpy.ndarray

pyproffit.hmc.ConstPM(x, bkg)[source]
pyproffit.hmc.DoubleBetaPM(x, beta, rc1, rc2, ratio, norm, bkg)[source]
class pyproffit.hmc.HMCModel(model, start=None, sd=None, limits=None, fix=None)[source]

Bases: object

Class containing pyproffit model structure for HMC optimization

Parameters:
  • model (function) – Function to be used as surface brightness model
  • vals (numpy.ndarray) – Array containing initial values for the parameters (optional)
SetPriors(start=None, sd=None, limits=None, fix=None)[source]

Set prior definition for the function parameters

Parameters:
  • start
  • sd
  • limits
  • fix
pyproffit.hmc.IntFuncPM(omega, rf, alpha, xmin, xmax)[source]
pyproffit.hmc.PowerLawPM(x, alpha, norm, pivot, bkg)[source]
pyproffit.hmc.VikhlininPM(x, beta, rc, alpha, rs, epsilon, gamma, norm, bkg)[source]
pyproffit.hmc.fit_profile_pymc3(hmcmod, prof, nmcmc=1000, tune=500, find_map=True, fitlow=0.0, fithigh=10000000000.0)[source]

pyproffit.miscellaneous module

pyproffit.miscellaneous.bkg_smooth(data, smoothing_scale=25)[source]

Smooth an input background image using a broad Gaussian

Parameters:
  • data (class:numpy.ndarray) – 2D array
  • smoothing_scale (float) – Size of Gaussian smoothing kernel in pixel. Defaults to 25
Returns:

smoothed image

Return type:

class:numpy.ndarray

pyproffit.miscellaneous.clean_bkg(img, bkg)[source]

Subtract statistically the background from a Poisson image

Parameters:
  • img (class:numpy.ndarray) – Input image
  • bkg (class:numpy.ndarray) – Background map
Returns:

Background subtracted Poisson image

Return type:

class:numpy.ndarray

pyproffit.miscellaneous.dist_eval(coords2d, index=None, x_c=None, y_c=None, metric='euclidean', selected=None)[source]

Computed distance of a set of points to a centroid

Parameters:
  • coords2d (class:numpy.ndarray) – 2D array with shape (N,2), N is the number of points
  • index (class:numpy.ndarray) – index defining input value for the centroid
  • x_c (float) – input centroid X axis coordinate
  • y_c (float) – input centroid Y axis coordinate
  • metric (str , optional) – metric system (defaults to ‘euclidean’)
  • selected (class:numpy.ndarray , optional) – index defining a selected subset of points to be used
Returns:

distance

Return type:

class:numpy.ndarray

pyproffit.miscellaneous.get_bary(x, y, x_c=None, y_c=None, weight=None, wdist=False)[source]

Compute centroid position and ellipse parameters from a set of points using principle component analysis

Parameters:
  • x (class:numpy.ndarray) – Array of positions on the X axis
  • y (class:numpy.ndarray) – Array of positions on the Y axis
  • x_c (float , optional) – Initial guess for the centroid X axis value
  • y_c (float , optional) – Initial guess for the centroid X axis value
  • weight (class:numpy.ndarray , optional) – Weights to be applied to the points when computing the average
  • wdist (bool) – Switch to apply the weights. Defaults to False
Returns:

  • x_c_w (float): Centroid X coordinate
  • y_c_w (float): Centroid Y coordinate
  • sig_x (float): X axis standard deviation
  • sig_y (float): Y axis standard deviation
  • r_cluster (float): characteristic size of cluster
  • semi_major_angle (float): rotation angle of ellipse
  • pos_err (float): positional error on the centroid

pyproffit.miscellaneous.heaviside(x)[source]

Heavyside theta function

Parameters:x (class:numpy.ndarray) – 2D array
Returns:filter
Return type:class:numpy.ndarray
pyproffit.miscellaneous.logbinning(binsize, maxrad)[source]

Set up a logarithmic binning scheme with a minimum bin size

Parameters:
  • binsize (float) – Minimum bin size in arcsec
  • maxrad (float) – Maximum extraction radius in arcmin
Returns:

Bins and bin width

Return type:

class:numpy.ndarray

pyproffit.miscellaneous.median_all_cov(dat, bins, ebins, rads, nsim=1000, fitter=None, thin=10)[source]

Generate Monte Carlo simulations of a Voronoi image and compute the median profile for each of them. The function returns an array of size (nbin, nsim) with nbin the number of bins in the profile and nsim the number of Monte Carlo simulations.

Parameters:
  • dat (class:pyproffit.data.Data) – A pyproffit.data.Data object containing the input Voronoi image and error map
  • bins (class:numpy.ndarray) – Central value of radial binning
  • ebins (class:numpy.ndarray) – Half-width of radial binning
  • rads (class:numpy.ndarray) – Array containing the distance to the center, in arcmin, for each pixel
  • nsim (int) – Number of Monte Carlo simulations to generate
  • fitter (class:pyproffit.fitter.Fitter) – A pyproffit.fitter.Fitter object containing the result of a fit to the background region, for subtraction of the background to the resulting profile
  • thin (int) – Number of blocks into which the calculation of the bootstrap will be divided. Increasing thin reduces memory usage drastically, at the cost of a modest increase in computation time.
Returns:

  • Samples of median profiles
  • Area of each bin

Return type:

class:numpy.ndarray

pyproffit.miscellaneous.medianval(vals, errs, nsim)[source]

Compute the median value of a sample of values and compute its uncertainty using Monte Carlo simulations

Parameters:
  • vals (class:numpy.ndarray) – Array containing the set of values in the sample
  • errs (class:numpy.ndarray) – Array contatining the error on each value
  • nsim (int) – Number of Monte Carlo simulations to be performed
Returns:

  • med (float): Median of sample
  • err (float): Error on median

pyproffit.miscellaneous.model_from_samples(x, model, samples, psfmat=None)[source]

Compute the median model and 1-sigma model envelope from a loaded chain, either from HMC or Emcee

Parameters:
  • x (class:numpy.ndarray) – Vector containing the X axis definition
  • model (class:pyproffit.models.Model) – Fitted model
  • samples (class:numpy.ndarray) – 2-dimensional array containing the parameter samples
Returns:

  • median (class:numpy.ndarray): Median model array
  • model_lo (class:numpy.ndarray): Lower 1-sigma envelope array
  • model_hi (class:numpy.ndarray): Upper 1-sigma envelope array

pyproffit.models module

pyproffit.models.BetaModel(x, beta, rc, norm, bkg)[source]

Single Beta model

\[I(r) = I_0 \left(1 + (x/r_c)^2\right) ^ {-3 \beta + 0.5} + B\]
Parameters:
  • x (numpy.ndarray) – Radius in arcmin
  • beta (float) – \(\beta\) parameter
  • rc (float) – rc parameter
  • norm (float) – log of I0 parameter
  • bkg (float) – log of B parameter
Returns:

Calculated model

Return type:

numpy.ndarray

pyproffit.models.BknPow(x, alpha1, alpha2, rf, norm, jump, bkg)[source]

Broken power law 3D model projected along the line of sight for discontinuity modeling

\[I(r) = I_0 \int F(\omega)^2 d\ell + B\]

with \(\omega^2 = r^2 + \ell^2\) and

\[\begin{split}F(\omega) = \left\{ \begin{array}{ll} \omega^{-\alpha_1}, & \omega<r_f \\ \frac{1}{C}\omega ^{-\alpha_2}, & \omega\geq r_f \end{array} \right.\end{split}\]
Parameters:
  • x (numpy.ndarray) – Radius in arcmin
  • alpha1 (float) – \(\alpha_1\) parameter
  • alpha2 (float) – \(\alpha_2\) parameter
  • rf (float) – rf parameter
  • norm (float) – log of I0 parameter
  • jump (float) – C parameter
  • bkg (float) – log of B parameter
Returns:

Calculated model

Return type:

numpy.ndarray

pyproffit.models.Const(x, bkg)[source]

Constal model for background fitting

Parameters:
  • x (numpy.ndarray) – Radius in arcmin
  • bkg (float) – log of B parameter
Returns:

Calculated model

Return type:

numpy.ndarray

pyproffit.models.DoubleBeta(x, beta, rc1, rc2, ratio, norm, bkg)[source]

Double beta model

\[I(r) = I_0 \left[ (1 + (x/r_{c,1}) ^ 2)^{-3 \beta + 0.5} + R ( 1 + (x/r_{c,2}) ^ 2 )^{-3 \beta+0.5} \right] + B\]
Parameters:
  • x (numpy.ndarray) – Radius in arcmin
  • beta (float) – \(\beta\) parameter
  • rc1 (float) – rc1 parameter
  • rc2 (float) – rc2 parameters
  • ratio (float) – R parameter
  • norm (float) – log of I0 parameter
  • bkg (float) – log of B parameter
Returns:

Calculated model

Return type:

numpy.ndarray

pyproffit.models.IntFunc(omega, rf, alpha, xmin, xmax)[source]

Numerical integration of a power law along the line of sight

\[\int_{x_{min}}^{x_{max}} \left(\frac{\omega^2 + \ell^2}{r_f^2}\right)^{-\alpha} d\ell\]
Parameters:
  • omega (float) – Projected radius
  • rf (float) – rf parameter
  • alpha (float) – \(\alpha\) parameter
  • xmin (float) – xmin parameter
  • xmax (float) – xmax parameter
Returns:

Line-of-sight integral

Return type:

float

class pyproffit.models.Model(model, vals=None)[source]

Bases: object

Class containing pyproffit models

Parameters:
  • model (function) – Function to be used as surface brightness model
  • vals (numpy.ndarray) – Array containing initial values for the parameters (optional)
SetErrors(vals)[source]

Set input values for the errors on the parameters

Parameters:vals (numpy.ndarray) – Array containing initial values for the errors
SetParameters(vals)[source]

Set input values for the model parameters

Parameters:vals (numpy.ndarray) – Array containing initial values for the parameters
pyproffit.models.PowerLaw(x, alpha, norm, pivot, bkg)[source]

Single power law model in projected space

\[I(r)=I_0\left( \frac{x}{x_p} \right)^{-\alpha} + B\]
Parameters:
  • x (numpy.ndarray) – Radius in arcmin
  • alpha (float) – \(\alpha\) parameter
  • norm (float) – log of I0 parameter
  • pivot (float) – \(x_p\) parameter
  • bkg (float) – log of B parameter
Returns:

Calculated model

Return type:

numpy.ndarray

pyproffit.models.Vikhlinin(x, beta, rc, alpha, rs, epsilon, gamma, norm, bkg)[source]

Simplified Vikhlinin+06 surface brightness model used in Ghirardini+19

\[I(r) = I_0\left( \frac{x}{r_c}\right)^{-\alpha} \left(1+(x/r_c)^2 \right)^{-3\beta + \alpha/2} \left( 1 + (x/r_s)^{\gamma} \right)^{-\epsilon/\gamma}\]
Parameters:
  • x (numpy.ndarray) – Radius in arcmin
  • beta (float) – \(\beta\) parameter
  • rc (float) – rc parameter
  • norm (float) – log of I0 parameter
  • alpha (float) – \(\alpha\) parameter
  • rs (float) – rs parameter
  • epsilon (float) – \(\epsilon\) parameter
  • gamma (float) – \(\gamma\) parameter
  • bkg (float) – log of B parameter
Returns:

Calculated model

Return type:

numpy.ndarray

pyproffit.power_spectrum module

class pyproffit.power_spectrum.PowerSpectrum(data, profile)[source]

Bases: object

Class to perform fluctuation power spectrum analysis from Poisson count images. This is the code used in Eckert et al. 2017.

Parameters:
  • data (class:pyproffit.data.Data) – Object of type pyproffit.data.Data including the image, exposure map, background map, and region definition
  • profile (class:pyproffit.profextract.Profile) – Object of type pyproffit.profextract.Profile including the extracted surface brightness profile
MexicanHat(modimg_file, z, region_size=1.0, factshift=1.5)[source]

Convolve the input image and model image with a set of Mexican Hat filters at various scales. The convolved images are automatically stored into FITS images called conv_scale_xx.fits and conv_beta_xx.fits, with xx the scale in kpc.

Parameters:
  • modimg_file (str) – Path to a FITS file including the model image, typically produced with pyproffit.profextract.Profile.SaveModelImage()
  • z (float) – Source redshift
  • region_size (float) – Size of the region of interest in Mpc. Defaults to 1.0
  • factshift (float) – Size of the border around the region, i.e. a region of size factshift * region_size is used for the computation. Defaults to 1.5
PS(z, region_size=1.0, radius_in=0.0, radius_out=1.0)[source]

Function to compute the power spectrum from existing Mexican Hat images in a given circle or annulus

Parameters:
  • z (float) – Source redshift
  • region_size (float) – Size of the region of interest in Mpc. Defaults to 1.0. This value must be equal to the region_size parameter used in pyproffit.power_spectrum.PowerSpectrum.MexicanHat().
  • radius_in (float) – Inner boundary in Mpc of the annulus to be used. Defaults to 0.0
  • radius_out (float) – Outer boundary in Mpc of the annulus to be used. Defaults to 1.0
Plot(save_plots=True, outps='power_spectrum.pdf', outamp='a2d.pdf', plot_3d=False, cfact=None)[source]

Plot the loaded power spectrum

Parameters:
  • save_plots (bool) – Indicate whether the plot should be saved to disk or not. Defaults to True.
  • outps (str) – Name of output file. Defaults to ‘power_spectrum.pdf’
  • outamp (str) – Name of output file to save the 2D amplitude plot. Defaults to ‘a2d.pdf’
  • plot_3d (bool) – Add or not the 3D power spectrum to the plot. Defaults to False
  • cfact (class:numpy.ndarray , optional) – 2D to 3D projection factor
ProjectionFactor(z, betaparams, region_size=1.0)[source]

Compute numerically the 2D to 3D deprojection factor. The routine is simulating 3D fluctuations using the surface brightness model and the region mask, projecting the 3D data along one axis, and computing the ratio of 2D to 3D power as a function of scale.

Caution: this is a memory intensive computation which will run into memory overflow if the image size is too large or the available memory is insufficient.

Parameters:
  • z (float) – Source redshift
  • betaparams (class:numpy.ndarray) – Parameters of the beta model or double beta model
  • region_size (float) – Size of the region of interest in Mpc. Defaults to 1.0. This value must be equal to the region_size parameter used in pyproffit.power_spectrum.PowerSpectrum.MexicanHat().
Returns:

Array of projection factors

Return type:

class:numpy.ndarray

Save(outfile, outcov='covariance.txt')[source]

Save the loaded power spectra to an output ASCII file

Parameters:
  • outfile (str) – Name of output ASCII file
  • outcov (str) – Output covariance matrix. Defaults to ‘covariance.txt’
pyproffit.power_spectrum.betamodel(x, par)[source]
pyproffit.power_spectrum.calc_mexicanhat(sc, img, mask, simmod)[source]

Filter an input image with a Mexican-hat filter

Parameters:
  • sc (float) – Mexican Hat scale in pixel
  • img (class:numpy.ndarray) – Image to be smoothed
  • mask (class:numpy.ndarray) – Mask image
  • simmod (class:numpy.ndarray) – Model surface brightness image
Returns:

Mexican Hat convolved image and SB model

Return type:

class:numpy.ndarray

pyproffit.power_spectrum.calc_projection_factor(nn, mask, betaparams, scale)[source]

Compute numerically the 2D to 3D deprojection factor. The routine is simulating 3D fluctuations using the surface brightness model and the region mask, projecting the 3D data along one axis, and computing the ratio of 2D to 3D power as a function of scale.

Caution: this is a memory intensive computation which will run into memory overflow if the image size is too large or the available memory is insufficient.

Parameters:
  • nn (int) – Image size in pixel
  • mask (class:numpy.ndarray) – Array defining the mask of active pixels, of size (nn , nn)
  • betaparams (class:numpy.ndarray) – Parameters of the beta model or double beta model
  • scale (class:numpy.ndarray) – Array of scales at which the projection factor should be computed
Returns:

Array containing wave number, 2D power, 2D amplitude, and 3D power

Return type:

class:numpy.ndarray

pyproffit.power_spectrum.calc_ps(region, img, mod, kr, nreg)[source]

Function to compute the power at a given scale kr from the Mexican Hat filtered images

Parameters:
  • region (class:numpy.ndarray) – Index defining the region from which the power spectrum will be extracted
  • img (class:numpy.ndarray) – Mexican Hat filtered image
  • mod (class:numpy.ndarray) – Mexican Hat filtered SB model
  • kr (float) – Extraction scale
  • nreg (int) – Number of subregions into which the image should be splitted to perform the bootstrap
Returns:

  • ps (float): Power at scale kr
  • psnoise (float): Noise at scale kr
  • vals (class:numpy.ndarray): set of values for bootstrap error calculation

pyproffit.power_spectrum.do_bootstrap(vals, nsample)[source]

Compute the covariance matrix of power spectra by bootstrapping the image

Parameters:
  • vals (class:numpy.ndarray) – Set of values
  • nsample (int) – Number of bootstrap samples
Returns:

2D covariance matrix

Return type:

class:numpy.ndarray

pyproffit.power_spectrum.doublebeta(x, pars)[source]

pyproffit.profextract module

class pyproffit.profextract.Profile(data=None, center_choice=None, maxrad=None, binsize=None, center_ra=None, center_dec=None, binning='linear', centroid_region=None, bins=None, cosmo=None)[source]

Bases: object

pyproffit.Profile class. The class allows the user to extract surface brightness profiles and use them to fit models, extract density profiles, etc.

Parameters:
  • data (class:pyproffit.Data) – Object of type pyproffit.data.Data containing the data to be used
  • center_choice (str) –

    Choice of the center of the surface brightness profile. Available options are “centroid”, “peak”, “custom_ima” and “custom_fk5”. Args:

    • ’centroid’: Compute image centroid and ellipticity. This is done by performing principle component analysis on the count image. If a dmfilth image exists, it will be used instead of the original count image.
    • ’peak’: Compute the surface brightness peak, The peak is computed as the maximum of the count image after a light smoothing. If a dmfilth image exists, it will be used instead of the original count image.
    • ’custom_fk5’: Use any custom center in FK5 coordinates, provided by the “center_ra” and “center_dec” arguments
    • ’custom_ima’: Similar to ‘custom_fk5’ but with input coordinates in image pixels
  • maxrad (float) – The maximum radius (in arcmin) out to which the surface brightness profile will be computed
  • binsize (float) – Minumum bin size (in arcsec).
  • center_ra (float) – User defined center R.A. If center_choice=’custom_fk5’ this is the right ascension in degrees. If center_choice=’custom_ima’ this is the image pixel on the X axis. If center_choice=’peak’ or ‘centroid’ this is not used.
  • center_dec (float) – User defined center declination. If center_choice=’custom_fk5’ this is the declination in degrees. If center_choice=’custom_ima’ this is the image pixel on the Y axis. If center_choice=’peak’ or ‘centroid’ this is not used.
  • binning (str) –

    Binning type. Available types are ‘linear’, ‘log’ or ‘custom’. Defaults to ‘linear’. Args:

    • ’linear’: Use a linear radial binning with bin size equal to ‘binsize’
    • ’log’: Use logarithmic binning, i.e. bin size increasing logarithmically with radius with a minimum bin size given by ‘binsize’
    • ’custom’: Any user-defined binning in the form of an input numpy array provided through the ‘bins’ option
  • centroid_region (float) – If center_choice=’centroid’, this option defines the radius of the region (in arcmin), centered on the center of the image, within which the centroid will be calculated. If centroid_region=None the entire image is used. Defaults to None.
  • bins (class:numpy.ndarray) – in case binning is set to ‘custom’, a numpy array containing the binning definition. For an input array of length N, the binning will contain N-1 bins with boundaries set as the values of the input array.
  • cosmo (class:astropy.cosmology) – An astropy.cosmology object containing the definition of the cosmological model. If cosmo=None, Planck 2015 cosmology is used.
AzimuthalScatter(nsect=12, model=None)[source]

Compute the azimuthal scatter profile around the loaded profile. The azimuthal scatter is defined as the standard deviation of the surface brightness in equispaced sectors with respect to the azimuthal mean,

\[\Sigma_X(r) = \frac{1}{N} \sum_{i=1}^N \frac{(S_i(r) - \langle S(r) \rangle)^2}{\langle S(r) \rangle^2}\]

with N the number of sectors and \(\langle S(r) \rangle\) the loaded mean surface brightness profile.

Parameters:
  • nsect (int) – Number of sectors from which the azimuthal scatter will be computed. Defaults to nsect=12
  • model (class:pyproffit.models.Model) – A pyproffit.models.Model object containing the background to be subtracted, in case the scatter is to be computed on background-subtracted profiles. Defaults to None (i.e. no background subtraction).
Backsub(fitter)[source]

Subtract a fitted background value from the loaded surface brightness profile. Each pyproffit model contains a ‘bkg’ parameter, which will be fitted and loaded in a Fitter object. The routine reads the value of ‘bkg’, subtracts it from the data, and adds its error in quadrature to the error profile.

Parameters:fitter (class:pyproffit.fitting.Fitter) – Object of type pyproffit.fitting.Fitter containing a model and optimization results
Emissivity(z=None, nh=None, kt=6.0, rmf=None, Z=0.3, elow=0.5, ehigh=2.0, arf=None, type='cr', lum_elow=0.5, lum_ehigh=2.0)[source]

Use XSPEC to compute the conversion from count rate to emissivity using the pyproffit.calc_emissivity routine (see its description)

Parameters:
  • z (float) – Source redshift
  • nh (float) – Source NH in units of 1e22 cm**(-2)
  • kt (float) – Source temperature in keV. Default to 6.0
  • rmf (str) – Path to response file (RMF/RSP)
  • Z (float) – Metallicity with respect to solar. Defaults to 0.3
  • elow (float) – Low-energy bound of the input image in keV. Defaults to 0.5
  • ehigh (float) – High-energy bound of the input image in keV. Defaults to 2.0
  • arf (str , optional) – Path to on-axis ARF in case response file type is RMF)
  • type (str) – Specify whether the exposure map is in units of sec (type=’cr’) or photon flux (type=’photon’). Defaults to ‘cr’
  • lum_elow (float) – Low energy bound (rest frame) for luminosity calculation. Defaults to 0.5
  • lum_ehigh (float) – High energy bound (rest frame) for luminosity calculation. Defaults to 2.0
Returns:

Conversion factor

Return type:

float

MedianSB(ellipse_ratio=1.0, rotation_angle=0.0, nsim=1000, outsamples=None, fitter=None, thin=10)[source]

Extract the median surface brightness profile in circular annuli from a provided Voronoi binned image, following the method outlined in Eckert et al. 2015

Parameters:
  • ellipse_ratio (float) – Ratio a/b of major to minor axis in the case of an elliptical annulus definition. Defaults to 1.0, i.e. circular annuli.
  • rotation_angle (float) – Rotation angle of the ellipse or box respective to the R.A. axis. Defaults 0.
  • nsim (int) – Number of Monte Carlo realizations of the Voronoi image to be performed
  • outsamples (str) – Name of output FITS file to store the bootstrap realizations of the median profile. Defaults to None
  • fitter (class:pyproffit.fitter.Fitter) – A pyproffit.fitter.Fitter object containing the result of a fit to the background region, for subtraction of the background to the resulting profile
  • thin (int) – Number of blocks into which the calculation of the bootstrap will be divided. Increasing thin reduces memory usage drastically, at the cost of a modest increase in computation time.
PSF(psffunc=None, psffile=None, psfimage=None, psfpixsize=None, sourcemodel=None)[source]

Function to calculate a PSF convolution matrix given an input PSF image or function. To compute the PSF mixing matrix, images of each annuli are convolved with the PSF image using FFT and determine the fraction of photons leaking into neighbouring annuli. FFT-convolved images are then used to determine a mixing matrix. See Eckert et al. 2020 for more details.

Parameters:
  • psffunc (function) – Function describing the radial shape of the PSF, with the radius in arcmin
  • psffile (str) – Path to file containing an image of the PSF. The pixel size must be equal to the pixel size of the image.
  • psfimage (class:numpy.ndarray) – Array containing an image of the PSF. The pixel size must be equal to the pixel size of the image.
  • psfpixsize (float) – (currently inactive) Pixel size of the PSF image in arcsec. Currently not implemented.
  • sourcemodel (class:pyproffit.models.Model) – Object of type pyproffit.models.Model including a surface brightness model to account for surface brightness gradients across the bins. If sourcemodel=None a flat distribution is assumed across each bin. Defaults to None
Plot(model=None, samples=None, outfile=None, axes=None, scatter=False, figsize=(13, 10), fontsize=40.0, xscale='log', yscale='log', fmt='o', markersize=7, lw=2, data_color='black', bkg_color='green', model_color='blue', **kwargs)[source]

Plot the loaded surface brightness profile

Parameters:
  • model (class:pyproffit.models.Model , optional) – If model is not None, plot the provided model of type pyproffit.models.Model together with the data. Defaults to None
  • samples (class:numpy.ndarray) – Use parameter samples outputted either by Emcee or HMC optimization to compute the median optimized models and upper and lower envelopes. Defaults to None
  • outfile (str , optional) – If outfile is not None, name of output file to save the plot. Defaults to None
  • axes (list , optional) – List of 4 numbers defining the X and Y axis ranges for the plot. Gives axes=[x1, x2, y1, y2], the X axis will be set between x1 and x2, and the Y axis will be set between y1 and y2.
  • scatter (bool) – Set whether the azimuthal scatter profile will be displayed instead of the surface brightness profile. Defaults to False.
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • xscale (str , optional) – Scale of the X axis. Defaults to ‘log’
  • yscale (str , optional) – Scale of the Y axis. Defaults to ‘log’
  • lw (int , optional) – Line width. Defaults to 2
  • fmt (str , optional) – Marker type following matplotlib convention. Defaults to ‘d’
  • markersize (int , optional) – Marker size. Defaults to 7
  • data_color (str , optional) – Color of the data points following matplotlib convention. Defaults to ‘red’
  • bkg_color (str , optional) – Color of the particle background following matplotlib convention. Defaults to ‘green’
  • model_color (str , optional) – Color of the model curve following matplotlib convention. Defaults to ‘blue’
  • kwargs – Arguments to be passed to matplotlib.pyplot.errorbar
SBprofile(ellipse_ratio=1.0, rotation_angle=0.0, angle_low=0.0, angle_high=360.0, box=False, width=None)[source]

Extract a surface brightness profile and store the results in the input Profile object

Parameters:
  • ellipse_ratio (float) – Ratio a/b of major to minor axis in the case of an elliptical annulus definition. Defaults to 1.0, i.e. circular annuli.
  • rotation_angle (float) – Rotation angle of the ellipse or box respective to the R.A. axis. Defaults 0.
  • angle_low (float) – In case the surface brightness profile should be extracted across a sector instead of the whole azimuth, lower position angle of the sector respective to the R.A. axis. Defaults to 0
  • angle_high (float) – In case the surface brightness profile should be extracted across a sector instead of the whole azimuth, upper position angle of the sector respective to the R.A. axis. Defaults to 360
  • voronoi (bool) – Set whether the input data is a Voronoi binned image (True) or a standard raw count image (False). Defaults to False.
  • box (bool) – Define whether the profile should be extract along an annulus or a box. The parameter definition of the box matches the DS9 definition. Defaults to False.
  • width (float) – In case box=True, set the full width of the box (in arcmin)
Save(outfile=None, model=None)[source]

Save the data loaded in the Profile class into an output FITS file.

Parameters:
  • outfile (str) – Output file name
  • model (class:pyproffit.models.Model , optional) – If model is not None, Object of type pyproffit.models.Model including the fitted model. Defaults to None
SaveModelImage(outfile, model=None, vignetting=True)[source]

Compute a model image and output it to a FITS file

Parameters:
  • model (class:pyproffit.models.Model) – Object of type pyproffit.models.Model including the surface brightness model. If model=None (default), the azimuthally averaged radial profile is interpolated at each point.
  • outfile (str) – Name of output file
  • vignetting (bool) – Choose whether the model will be convolved with the vignetting model (i.e. multiplied by the exposure map) or if the actual surface brightness will be extracted (False). Defaults to True
pyproffit.profextract.plot_multi_profiles(profs, labels=None, outfile=None, axes=None, figsize=(13, 10), fontsize=40, xscale='log', yscale='log', fmt='o', markersize=7)[source]

Plot multiple surface brightness profiles on a single plot. This feature can be useful e.g. to compare profiles across multiple sectors

Parameters:
  • profs (tuple) – List of Profile objects to be plotted
  • labels (tuple) – List of labels for the legend (default=None)
  • outfile (str) – If outfile is not None, path to file name to output the plot
  • axes (list , optional) – List of 4 numbers defining the X and Y axis ranges for the plot. Gives axes=[x1, x2, y1, y2], the X axis will be set between x1 and x2, and the Y axis will be set between y1 and y2.
  • figsize (tuple , optional) – Size of figure. Defaults to (13, 10)
  • fontsize (int , optional) – Font size of the axis labels. Defaults to 40
  • xscale (str , optional) – Scale of the X axis. Defaults to ‘log’
  • yscale (str , optional) – Scale of the Y axis. Defaults to ‘log’
  • lw (int , optional) – Line width. Defaults to 2
  • fmt (str , optional) – Marker type following matplotlib convention. Defaults to ‘d’
  • markersize (int , optional) – Marker size. Defaults to 7

pyproffit.reload module

pyproffit.reload.Reload(infile, model=None)[source]

Reload the results of a previous session into new Profile and Data objects. If model=None and a Model structure is present, it will be read and its values will be stored into the input Model structure.

Parameters:
Returns:

Module contents