imcascade

imcascade: Fitting astronomical images using a ‘cascade’ of Gaussians

Submodules

Package Contents

Classes

Fitter

A Class used fit images with MultiGaussModel

ImcascadeResults

A class used for collating imcascade results and performing analysis

Attributes

__version__

imcascade.__version__ = '1.0'
class imcascade.Fitter(img, sig, psf_sig, psf_a, weight=None, mask=None, sky_model=True, sky_type='tilted-plane', render_mode='hybrid', log_weight_scale=True, verbose=True, psf_shape=None, init_dict={}, bounds_dict={}, log_file=None)

Bases: imcascade.mgm.MultiGaussModel

A Class used fit images with MultiGaussModel

This is the main class used to fit imcascade models

Parameters

  • img (2D Array) – Data to be fit, it is assumed to be a cutout with the object of interest in the center of the image

  • sig (1D Array) – Widths of Gaussians to be used in MultiGaussModel

  • psf_sig (1D array, None) – Width of Gaussians used to approximate psf

  • psf_a (1D array, None) – Weights of Gaussians used to approximate psf If both psf_sig and psf_a are None then will run in Non-psf mode

  • weight (2D Array, optional) – Array of pixel by pixel weights to be used in fitting. Must be same shape as ‘img’ If None, all the weights will be set to 1.

  • mask (2D Array, optional) – Array with the same shape as ‘img’ denoting which, if any, pixels to mask during fitting process. Values of ‘1’ or ‘True’ values for the pixels to be masked. If set to ‘None’ then will not mask any pixels. In practice, the weights of masked pixels is set to ‘0’.

  • sky_model (bool, optional) – If True will incorperate a tilted plane sky model. Reccomended to be set to True

  • sky_type (str, 'tilted-plane' or 'flat') – Function used to model sky. Default is tilted plane with 3 parameters, const bkg and slopes in each directin. ‘flat’ uses constant background model with 1 parameter.

  • render_mode ('hybrid', 'erf' or 'gauss') – Option to decide how to render models. ‘erf’ analytically computes the integral over the pixel of each profile therefore is more accurate but more computationally intensive. ‘gauss’ assumes the center of a pixel provides a reasonble estimate of the average flux in that pixel. ‘gauss’ is faster but far less accurate for objects which vary on O(pixel size), so use with caution. ‘hybrid’ is the defualt, uses ‘erf’ for components with width < 5 to ensure accuracy and uses ‘gauss’ otherwise as it is accurate enough and faster. Also assumes all flux > 5 sigma for components is 0.

  • log_weight_scale (bool, optional) – Wether to treat weights as log scale, Default True

  • verbose (bool, optional) – If true will log and print out errors

  • psf_shape (dict, Optional) – Dictionary containg at ‘q’ and ‘phi’ that define the shape of the PSF. Note that this slows down model rendering significantly so only reccomended if neccesary.

  • init_dict (dict, Optional) – Dictionary specifying initial guesses for least_squares fitting. The code is desigined to make ‘intelligent’ guesses if none are provided

  • bounds_dict (dict, Optional) – Dictionary specifying boundss for least_squares fitting and priors. The code is desigined to make ‘intelligent’ guesses if none are provided

resid_1d(params)

Given a set of parameters returns the 1-D flattened residuals when compared to the Data, to be used in run_ls_min Function

Parameters

params (Array) – List of parameters to define model

Returns

resid_flatten – 1-D array of the flattened residuals

Return type

array

run_ls_min(ls_kwargs={})

Function to run a least_squares minimization routine using pre-determined inital guesses and bounds.

Utilizes the scipy least_squares routine (https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.least_squares.html)

Parameters

ls_kwargs (dict, optional) – Optional list of arguments to be passes to least_squares routine

Returns

min_param – Returns a 1D array containing the optimized parameters that describe the best fit model.

Return type

1D array

set_up_express_run(set_params=None)

Function to set up ‘express’ run using pre-rendered images with a fixed x0,y0, phi and q. Sets class attribute ‘express_gauss_arr’ which is needed to run dynesty or emcee in express mode

Parameters

set_params (len(4) array-like, optional) – Parameters (x0,y0,q,phi) to use to per-render images. If None will call run_ls_min() or use stored min_res to find parameters.

Returns

express_gauss_arr – Returns a 3-D with pre-rendered images based on input parameters

Return type

array (shape[0],shape[1], Ndof_gauss)

make_express_model(exp_params)

Function to generate a model for a given set of paramters, specifically using the pre-renedered model for the ‘express’ mode

Parameters

exo_params (Array) – List of parameters to define model. Length is Ndof_gauss + Ndof_sky since the structural parameters (x0,y0,q, PA) are set

Returns

model – Model image based on input parameters

Return type

2D-array

chi_sq(params)

Function to calculate chi_sq for a given set of paramters

Parameters

params (Array) – List of parameters to define model

Returns

chi^2 – Chi squared statistic for the given set of parameters

Return type

float

log_like(params)

Function to calculate the log likeliehood for a given set of paramters

Parameters

params (Array) – List of parameters to define model

Returns

log likeliehood – log likeliehood for a given set of paramters, defined as -0.5*chi^2

Return type

float

ptform(u)

Prior transformation function to be used in dynesty ‘full’ mode

Parameters

u (array) – array of random numbers from 0 to 1

Returns

x – array containing distribution of parameters from prior

Return type

array

log_like_exp(exp_params)

Function to calculate the log likeliehood for a given set of paramters, specifically using the pre-renedered model for the ‘express’ mode

Parameters

exo_params (Array) – List of parameters to define model. Length is Ndof_gauss + Ndof_sky since the structural parameters (x0,y0,q, PA) are set

Returns

log likeliehood – log likeliehood for a given set of paramters, defined as -0.5*chi^2

Return type

float

ptform_exp_ls(u)

Prior transformation function to be used in dynesty ‘express’ mode using gaussian priors defined by the results of the least_squares minimization

Parameters

u (array) – array of random numbers from 0 to 1

Returns

x – array containing distribution of parameters from prior

Return type

array

ptform_exp_lin_cauchy(u)

Prior transformation function to be used in dynesty ‘express’ mode using gaussian priors defined by the results of the least_squares minimization

Parameters

u (array) – array of random numbers from 0 to 1

Returns

x – array containing distribution of parameters from prior

Return type

array

ptform_exp_unif(u)

Prior transformation function to be used in dynesty ‘express’ mode using unifrom priors defined by self.lb and self.ub

Parameters

u (array) – array of random numbers from 0 to 1

Returns

x – array containing distribution of parameters from prior

Return type

array

run_dynesty(method='full', sampler_kwargs={}, run_nested_kwargs={}, prior='min_results')

Function to run dynesty to sample the posterior distribution using either the ‘full’ methods which explores all paramters, or the ‘express’ method which sets the structural parameters.

Parameters
Returns

Posterior – posterior distribution derrived. If method is ‘express’, the first 4 columns, containg x0, y0, PA and q, are all the same and equal to values used to pre-render the images

Return type

Array

save_results(file_name)

Function to save results after run_ls_min, run_dynesty and/or run_emcee is performed. Will be saved as an ASDF file.

Parameters

file_name (str) – Str defining location of where to save data

Returns

dict_saved – Dictionary containing all the save quantities

Return type

dict

class imcascade.ImcascadeResults(Obj, thin_posterior=1)

A class used for collating imcascade results and performing analysis

Parameters

  • Obj (imcascade.fitter.Fitter class, dictionary or str) – Object which contains the data to be analyzed. Can be a Fitter object once the run_(ls_min,dynesty, emcee) has been ran. If it is a dictionay needs to contain, at bare minmum the variables sig, Ndof, Ndof_sky, Ndof_gauss, log_weight_scale and either min_param or posterior. If a string is passed it will be interreted as a file locations with an ASDF file containing the neccesary information.

  • thin_posterior (int (optional)) – Factor by which to thin the posterior distribution by. While one wants to ensure the posterior is large enough, some of this analysis can take time if you have >10^6 samples so this is one way to speed up this task but use with caution.

calc_flux(cutoff=None)

Calculate flux of given results

Parameters

cutoff (float (optional)) – Radius out to which to consider the profile. Generally this should be around the half-width of the image or the largest gaussian width use

Returns

Flux – Total flux of best fit model

Return type

float or Array

_min_calc_rX(X, cutoff=None)

Old and slow Function to calculate the radius containing X percent of the light

Parameters

  • X (float) – Fractional radius of intrest to calculate. if X < 1 will take as a fraction, else will interpret as percent and divide X by 100. i.e. to calculate the radius containing 20% of the light, once can either pass X = 20 or 0.2

  • cutoff (float (optional)) – Radius out to which to consider the profile. Generally this should be around the half-width of the image or the largest gaussian width used

Returns

r_X – The radius containg X percent of the light

Return type

float or Array

calc_rX(X, cutoff=None)

Function to calculate the radius containing X percent of the light

Parameters

  • X (float) – Fractional radius of intrest to calculate. if X < 1 will take as a fraction, else will interpret as percent and divide X by 100. i.e. to calculate the radius containing 20% of the light, once can either pass X = 20 or 0.2

  • cutoff (float (optional)) – Radius out to which to consider the profile. Generally this should be around the half-width of the image or the largest gaussian width used

Returns

r_X – The radius containg X percent of the light

Return type

float or Array

calc_r90(cutoff=None)

Wrapper function to calculate the radius containing 90% of the light

Parameters

cutoff (float (optional)) – Radius out to which to consider the profile. Generally this should be around the half-width of the image or the largest gaussian width use

Returns

r_90 – The radius containg 90 percent of the light

Return type

float or Array

calc_r80(cutoff=None)

Wrapper function to calculate the radius containing 80% of the light

Parameters

cutoff (float (optional)) – Radius out to which to consider the profile. Generally this should be around the half-width of the image or the largest gaussian width use

Returns

r_80 – The radius containg 80 percent of the light

Return type

float or Array

calc_r50(cutoff=None)

Wrapper function to calculate the radius containing 50% of the light, or the effective radius

Parameters

cutoff (float (optional)) – Radius out to which to consider the profile. Generally this should be around the half-width of the image or the largest gaussian width use

Returns

r_50 – The radius containg 50 percent of the light

Return type

float or Array

calc_r20(cutoff=None)

Wrapper function to calculate the radius containing 20% of the light

Parameters

cutoff (float (optional)) – Radius out to which to consider the profile. Generally this should be around the half-width of the image or the largest gaussian width use

Returns

r_20 – The radius containg 20 percent of the light

Return type

float or Array

calc_sbp(r, return_ind=False)

Function to calculate surface brightness profiles for the given results

Parameters

  • r (float or array) – Radii (in pixels) at which to evaluate the surface brightness profile

  • return_ind (bool (optional)) – If False will only return the sum of all gaussian, i.e. the best fit profile. If true will return an array with +1 dimensions containing the profiles of each individual gaussian component

Returns

SBP – Surface brightness profiles evaluated at ‘r’. If ‘return_ind = True’, returns the profile of each individual gaussian component

Return type

array

calc_obs_sbp(r, return_ind=False)

Function to calculate the observed surface brightness profiles, i.e. convolved with the PSF for the given results

Parameters

  • r (float or array) – Radii (in pixels) at which to evaluate the surface brightness profile

  • return_ind (bool (optional)) – If False will only return the sum of all gaussian, i.e. the best fit profile. If true will return an array with +1 dimensions containing the profiles of each individual gaussian component

Returns

obsereved SBP – Observed surface brightness profiles evaluated at ‘r’. If ‘return_ind = True’, returns the profile of each individual gaussian component

Return type

array

calc_cog(r, return_ind=False, norm=False, cutoff=None)

Function to calculate curves-of-growth for the given results

Parameters

  • r (float or array) – Radii (in pixels) at which to evaluate the surface brightness profile

  • return_ind (bool (optional)) – If False will only return the sum of all gaussian, i.e. the best fit profile. If true will return an array with +1 dimensions containing the profiles of each individual gaussian component

  • norm (Bool (optional)) – Wether to normalize curves-of-growth to total flux, calculated using ‘self.calc_flux’. Does nothing if ‘return_ind = True’

  • cutoff (Float (optional)) – Cutoff radius used in ‘self.calc_flux’, only is used if ‘norm’ is True

Returns

COG – curves-of-growth evaluated at ‘r’. If ‘return_ind = True’, returns the profile of each individual gaussian component

Return type

array

calc_obs_cog(r, return_ind=False, norm=False, cutoff=None)

Function to calculate the observed curve of growth, i.e. convolved with the PSF for the given results

Parameters

  • r (float or array) – Radii (in pixels) at which to evaluate the surface brightness profile

  • return_ind (bool (optional)) – If False will only return the sum of all gaussian, i.e. the best fit profile. If true will return an array with +1 dimensions containing the profiles of each individual gaussian component

  • norm (Bool (optional)) – Wether to normalize curves-of-growth to total flux, calculated using ‘self.calc_flux’. Does nothing if ‘return_ind = True’

  • cutoff (Float (optional)) – Cutoff radius used in ‘self.calc_flux’, only is used if ‘norm’ is True

Returns

observed COG – curves-of-growth evaluated at ‘r’. If ‘return_ind = True’, returns the profile of each individual gaussian component

Return type

array

run_basic_analysis(zpt=None, cutoff=None, errp_lo=16, errp_hi=84, save_results=False, save_file='./imcascade_results.asdf')

Function to calculate a set of common variables and save the save the results

Parameters

  • zpt (float (optional)) – photometric zeropoint for the data. if not ‘None’, will also calculate magnitude

  • cutoff (float (optional)) – Radius out to which to consider the profile. Generally this should be around the half-width of the image or the largest gaussian width use

  • errp_(lo,hi) (float (optional)) – percentiles to be used to calculate the lower and upper error bars from the posterior distribution. Default is 16 and 84, corresponding to 1-sigma for a guassian distribtuion

  • save_results (bool (optional)) – If true will save results to file. If input is a file, will add to given file, else will save to file denoted by ‘save_file’ (see below)

  • save_file (str) – String to describe where to save file, only applicaple if the input is not a file.

Returns

res_dict – Dictionary contining the results of the analysis

Return type

dictionary

calc_iso_r(I, zpt=None, pix_scale=None)

Function to calculate the isophotal radius

Parameters

  • I (float) – Surface brightness target to define the isophotal radii. By defualt this shoud be in image units unless both zpt and pix_scale are given, then I is interpreted as mag per arcsec^2.

  • zpt (float (optional)) – zeropoint magnitude of image, used to convert I to mag per arcsec^2

  • pix_scale (float (optional)) – pixel scale in units of arcseconds/pixel, used to convert I to mag per arcsec^2

Returns

r_I – The radius, in pixel units, where the surface brightness profile matches I

Return type

float or Array

calc_petro_r(P_ratio=0.2, r_fac_min=0.8, r_fac_max=1.25)

Function to calculate the petrosian radii of a galaxy

Parameters

  • P_ratio (float (optional)) – The Petrosian ratio which defines the Petrosian radii, default is 0.2

  • r_fac_min (float (optional)) – lower multiplicative factor which is used to integrate flux, default 0.8

  • r_fac_max (float (optional)) – higher multiplicative factor which is used to inegrate flux, default 1.25

Returns

r_I – The radius, in pixel units, where the surface brightness profile matches I

Return type

float or Array

make_diagnostic_fig()

Function which generates a diagnostic figure to assess fit

Returns

fig – matplotlib figure object

Return type

matplotlib figure