ApertureStats

class photutils.aperture.ApertureStats(data, aperture, *, error=None, mask=None, wcs=None, sigma_clip=None, sum_method='exact', subpixels=5, local_bkg=None)[source]

Bases: object

Class to create a catalog of statistics for pixels within an aperture.

Note that this class returns the statistics of the input data values within the aperture. It does not convert data in surface brightness units to flux or counts. Conversion from surface-brightness units should be performed before using this function.

Parameters:
data2D ndarray, Quantity, NDData

The 2D array from which to calculate the source properties. For accurate source properties, data should be background-subtracted. Non-finite data values (NaN and inf) are automatically masked.

apertureAperture

The aperture to apply to the data. The aperture object may contain more than one position. If aperture is a SkyAperture object, then a WCS must be input using the wcs keyword.

error2D ndarray or Quantity, optional

The total error array corresponding to the input data array. error is assumed to include all sources of error, including the Poisson error of the sources (see calc_total_error) . error must have the same shape as the input data. If data is a Quantity array then error must be a Quantity array (and vice versa) with identical units. Non-finite error values (NaN and +/- inf) are not automatically masked, unless they are at the same position of non-finite values in the input data array. Such pixels can be masked using the mask keyword.

mask2D ndarray (bool), optional

A boolean mask with the same shape as data where a True value indicates the corresponding element of data is masked. Masked data are excluded from all calculations. Non-finite values (NaN and inf) in the input data are automatically masked.

wcsWCS object or None, optional

A world coordinate system (WCS) transformation that supports the astropy shared interface for WCS (e.g., astropy.wcs.WCS, gwcs.wcs.WCS). wcs is required if the input aperture is a SkyAperture. If None, then all sky-based properties will be set to None.

sigma_clipNone or astropy.stats.SigmaClip instance, optional

A SigmaClip object that defines the sigma clipping parameters. If None then no sigma clipping will be performed.

sum_method{‘exact’, ‘center’, ‘subpixel’}, optional

The method used to determine the overlap of the aperture on the pixel grid. This method is used only for calculating the sum, sum_error, sum_aper_area, data_sumcutout, and error_sumcutout properties. All other properties use the “center” aperture mask method. Not all options are available for all aperture types. The following methods are available:

  • 'exact' (default): The the exact fractional overlap of the aperture and each pixel is calculated. The aperture weights will contain values between 0 and 1.

  • 'center': A pixel is considered to be entirely in or out of the aperture depending on whether its center is in or out of the aperture. The aperture weights will contain values only of 0 (out) and 1 (in).

  • 'subpixel': A pixel is divided into subpixels (see the subpixels keyword), each of which are considered to be entirely in or out of the aperture depending on whether its center is in or out of the aperture. If subpixels=1, this method is equivalent to 'center'. The aperture weights will contain values between 0 and 1.

subpixelsint, optional

For the 'subpixel' method, resample pixels by this factor in each dimension. That is, each pixel is divided into subpixels**2 subpixels. This keyword is ignored unless sum_method='subpixel'.

local_bkgfloat, ndarray, Quantity, or None

The per-pixel local background values to subtract from the data before performing measurements. If input as any array, the order of local_bkg values corresponds to the order of the input aperture positions. local_bkg must have the same length as the the input aperture or must be a scalar value, which will be broadcast to all apertures. If None, then no local background subtraction is performed. If the input data has units, then local_bkg must be a Quantity with the same units.

Notes

data should be background-subtracted for accurate source properties. In addition to global background subtraction, local background subtraction can be performed using the local_bkg keyword values.

Most source properties are calculated using the “center” aperture-mask method, which gives aperture weights of 0 or 1. This avoids the need to compute weighted statistics — the data pixel values are directly used.

The input sum_method and subpixels keywords are used to determine the aperture-mask method when calculating the sum-related properties: sum, sum_error, sum_aper_area, data_sumcutout, and error_sumcutout. The default is sum_method='exact', which produces exact aperture-weighted photometry.

Examples

>>> from photutils.datasets import make_4gaussians_image
>>> from photutils.aperture import CircularAperture, ApertureStats
>>> data = make_4gaussians_image()
>>> aper = CircularAperture((150, 25), 8)
>>> aperstats = ApertureStats(data, aper)
>>> print(aperstats.xcentroid)  
149.98737072209013
>>> print(aperstats.ycentroid)  
24.99729176183652
>>> print(aperstats.centroid)  
[149.98737072  24.99729176]
>>> print(aperstats.mean, aperstats.median, aperstats.std) 
46.861845146453526 33.743501730319 38.25291812758177
>>> print(aperstats.sum)  
9118.129697119366
>>> print(aperstats.sum_aper_area) 
201.0619298297468 pix2
>>> # more than one aperture position
>>> aper2 = CircularAperture(((150, 25), (90, 60)), 10)
>>> aperstats2 = ApertureStats(data, aper2)
>>> print(aperstats2.xcentroid)  
[149.97230436  90.00833613]
>>> print(aperstats2.sum)  
[ 9863.56195844 36629.52906175]

Attributes Summary

bbox

The BoundingBox of the aperture.

bbox_xmax

The maximum x pixel index of the bounding box.

bbox_xmin

The minimum x pixel index of the bounding box.

bbox_ymax

The maximum y pixel index of the bounding box.

bbox_ymin

The minimum y pixel index of the bounding box.

biweight_location

The biweight location of the unmasked pixel values within the aperture.

biweight_midvariance

The biweight midvariance of the unmasked pixel values within the aperture.

center_aper_area

The total area of the unmasked pixels within the aperture using the "center" aperture mask method.

centroid

The (x, y) coordinate of the centroid.

covar_sigx2

The (0, 0) element of the covariance matrix, representing \(\sigma_x^2\), in units of pixel**2.

covar_sigxy

The (0, 1) and (1, 0) elements of the covariance matrix, representing \(\sigma_x \sigma_y\), in units of pixel**2.

covar_sigy2

The (1, 1) element of the covariance matrix, representing \(\sigma_y^2\), in units of pixel**2.

covariance

The covariance matrix of the 2D Gaussian function that has the same second-order moments as the source.

covariance_eigvals

The two eigenvalues of the covariance matrix in decreasing order.

cutout_centroid

The (x, y) coordinate, relative to the cutout data, of the centroid within the aperture.

cxx

SourceExtractor's CXX ellipse parameter in units of pixel**(-2).

cxy

SourceExtractor's CXY ellipse parameter in units of pixel**(-2).

cyy

SourceExtractor's CYY ellipse parameter in units of pixel**(-2).

data_cutout

A 2D aperture-weighted cutout from the data using the aperture mask with the "center" method as a MaskedArray.

data_sumcutout

A 2D aperture-weighted cutout from the data using the aperture mask with the input sum_method method as a MaskedArray.

eccentricity

The eccentricity of the 2D Gaussian function that has the same second-order moments as the source.

ellipticity

1.0 minus the ratio of the lengths of the semimajor and semiminor axes (or 1.0 minus the elongation).

elongation

The ratio of the lengths of the semimajor and semiminor axes.

error_sumcutout

A 2D aperture-weighted error cutout using the aperture mask with the input sum_method method as a MaskedArray.

fwhm

The circularized full width at half maximum (FWHM) of the 2D Gaussian function that has the same second-order central moments as the source.

gini

The Gini coefficient of the unmasked pixel values within the aperture.

id

The aperture identification number(s).

ids

The aperture identification number(s), always as an iterable ndarray.

inertia_tensor

The inertia tensor of the source for the rotation around its center of mass.

isscalar

Whether the instance is scalar (e.g., a single aperture position).

mad_std

The standard deviation calculated using the median absolute deviation (MAD).

max

The maximum of the unmasked pixel values within the aperture.

mean

The mean of the unmasked pixel values within the aperture.

median

The median of the unmasked pixel values within the aperture.

min

The minimum of the unmasked pixel values within the aperture.

mode

The mode of the unmasked pixel values within the aperture.

moments

Spatial moments up to 3rd order of the source.

moments_central

Central moments (translation invariant) of the source up to 3rd order.

n_apertures

The number of positions in the input aperture.

orientation

The angle between the x axis and the major axis of the 2D Gaussian function that has the same second-order moments as the source.

properties

A sorted list of built-in source properties.

semimajor_sigma

The 1-sigma standard deviation along the semimajor axis of the 2D Gaussian function that has the same second-order central moments as the source.

semiminor_sigma

The 1-sigma standard deviation along the semiminor axis of the 2D Gaussian function that has the same second-order central moments as the source.

sky_centroid

The sky coordinate of the centroid of the unmasked pixels within the aperture, returned as a SkyCoord object.

sky_centroid_icrs

The sky coordinate in the International Celestial Reference System (ICRS) frame of the centroid of the unmasked pixels within the aperture, returned as a SkyCoord object.

std

The standard deviation of the unmasked pixel values within the aperture.

sum

The sum of the unmasked data values within the aperture.

sum_aper_area

The total area of the unmasked pixels within the aperture using the input sum_method aperture mask method.

sum_err

The uncertainty of sum , propagated from the input error array.

var

The variance of the unmasked pixel values within the aperture.

xcentroid

The x coordinate of the centroid.

ycentroid

The y coordinate of the centroid.

Methods Summary

copy()

Return a deep copy of this object.

get_id(id_num)

Return a new ApertureStats object for the input ID number only.

get_ids(id_nums)

Return a new ApertureStats object for the input ID numbers only.

to_table([columns])

Create a QTable of source properties.

Attributes Documentation

bbox

The BoundingBox of the aperture.

Note that the aperture bounding box is calculated using the exact size of the aperture, which may be slightly larger than the aperture mask calculated using the “center” mode.

bbox_xmax

The maximum x pixel index of the bounding box.

Note that this value is inclusive, unlike numpy slice indices.

bbox_xmin

The minimum x pixel index of the bounding box.

bbox_ymax

The maximum y pixel index of the bounding box.

Note that this value is inclusive, unlike numpy slice indices.

bbox_ymin

The minimum y pixel index of the bounding box.

biweight_location

The biweight location of the unmasked pixel values within the aperture.

See astropy.stats.biweight_location.

biweight_midvariance

The biweight midvariance of the unmasked pixel values within the aperture.

See astropy.stats.biweight_midvariance

center_aper_area

The total area of the unmasked pixels within the aperture using the “center” aperture mask method.

centroid

The (x, y) coordinate of the centroid.

The centroid is computed as the center of mass of the unmasked pixels within the aperture.

covar_sigx2

The (0, 0) element of the covariance matrix, representing \(\sigma_x^2\), in units of pixel**2.

covar_sigxy

The (0, 1) and (1, 0) elements of the covariance matrix, representing \(\sigma_x \sigma_y\), in units of pixel**2.

covar_sigy2

The (1, 1) element of the covariance matrix, representing \(\sigma_y^2\), in units of pixel**2.

covariance

The covariance matrix of the 2D Gaussian function that has the same second-order moments as the source.

covariance_eigvals

The two eigenvalues of the covariance matrix in decreasing order.

cutout_centroid

The (x, y) coordinate, relative to the cutout data, of the centroid within the aperture.

The centroid is computed as the center of mass of the unmasked pixels within the aperture.

cxx

SourceExtractor’s CXX ellipse parameter in units of pixel**(-2).

The ellipse is defined as

\[cxx (x - \bar{x})^2 + cxy (x - \bar{x}) (y - \bar{y}) + cyy (y - \bar{y})^2 = R^2\]

where \(R\) is a parameter which scales the ellipse (in units of the axes lengths). SourceExtractor reports that the isophotal limit of a source is well represented by \(R \approx 3\).

cxy

SourceExtractor’s CXY ellipse parameter in units of pixel**(-2).

The ellipse is defined as

\[cxx (x - \bar{x})^2 + cxy (x - \bar{x}) (y - \bar{y}) + cyy (y - \bar{y})^2 = R^2\]

where \(R\) is a parameter which scales the ellipse (in units of the axes lengths). SourceExtractor reports that the isophotal limit of a source is well represented by \(R \approx 3\).

cyy

SourceExtractor’s CYY ellipse parameter in units of pixel**(-2).

The ellipse is defined as

\[cxx (x - \bar{x})^2 + cxy (x - \bar{x}) (y - \bar{y}) + cyy (y - \bar{y})^2 = R^2\]

where \(R\) is a parameter which scales the ellipse (in units of the axes lengths). SourceExtractor reports that the isophotal limit of a source is well represented by \(R \approx 3\).

data_cutout

A 2D aperture-weighted cutout from the data using the aperture mask with the “center” method as a MaskedArray.

The cutout does not have units due to current limitations of masked quantity arrays.

The mask is True for pixels from the input mask, non-finite data values (NaN and inf), sigma-clipped pixels within the aperture, and pixels where the aperture mask has zero weight.

data_sumcutout

A 2D aperture-weighted cutout from the data using the aperture mask with the input sum_method method as a MaskedArray.

The cutout does not have units due to current limitations of masked quantity arrays.

The mask is True for pixels from the input mask, non-finite data values (NaN and inf), sigma-clipped pixels within the aperture, and pixels where the aperture mask has zero weight.

eccentricity

The eccentricity of the 2D Gaussian function that has the same second-order moments as the source.

The eccentricity is the fraction of the distance along the semimajor axis at which the focus lies.

\[e = \sqrt{1 - \frac{b^2}{a^2}}\]

where \(a\) and \(b\) are the lengths of the semimajor and semiminor axes, respectively.

ellipticity

1.0 minus the ratio of the lengths of the semimajor and semiminor axes (or 1.0 minus the elongation).

\[\mathrm{ellipticity} = 1 - \frac{b}{a}\]

where \(a\) and \(b\) are the lengths of the semimajor and semiminor axes, respectively.

elongation

The ratio of the lengths of the semimajor and semiminor axes.

\[\mathrm{elongation} = \frac{a}{b}\]

where \(a\) and \(b\) are the lengths of the semimajor and semiminor axes, respectively.

error_sumcutout

A 2D aperture-weighted error cutout using the aperture mask with the input sum_method method as a MaskedArray.

The cutout does not have units due to current limitations of masked quantity arrays.

The mask is True for pixels from the input mask, non-finite data values (NaN and inf), sigma-clipped pixels within the aperture, and pixels where the aperture mask has zero weight.

fwhm

The circularized full width at half maximum (FWHM) of the 2D Gaussian function that has the same second-order central moments as the source.

\[\begin{split}\mathrm{FWHM} & = 2 \sqrt{2 \ln(2)} \sqrt{0.5 (a^2 + b^2)} \\ & = 2 \sqrt{\ln(2) \ (a^2 + b^2)}\end{split}\]

where \(a\) and \(b\) are the 1-sigma lengths of the semimajor (semimajor_sigma) and semiminor (semiminor_sigma) axes, respectively.

gini

The Gini coefficient of the unmasked pixel values within the aperture.

The Gini coefficient is calculated using the prescription from Lotz et al. 2004 as:

\[G = \frac{1}{\left | \bar{x} \right | n (n - 1)} \sum^{n}_{i} (2i - n - 1) \left | x_i \right |\]

where \(\bar{x}\) is the mean over pixel values \(x_i\) within the aperture.

The Gini coefficient is a way of measuring the inequality in a given set of values. In the context of galaxy morphology, it measures how the light of a galaxy image is distributed among its pixels. A Gini coefficient value of 0 corresponds to a galaxy image with the light evenly distributed over all pixels while a Gini coefficient value of 1 represents a galaxy image with all its light concentrated in just one pixel.

id

The aperture identification number(s).

ids

The aperture identification number(s), always as an iterable ndarray.

inertia_tensor

The inertia tensor of the source for the rotation around its center of mass.

isscalar

Whether the instance is scalar (e.g., a single aperture position).

mad_std

The standard deviation calculated using the median absolute deviation (MAD).

The standard deviation estimator is given by:

\[\sigma \approx \frac{\textrm{MAD}}{\Phi^{-1}(3/4)} \approx 1.4826 \ \textrm{MAD}\]

where \(\Phi^{-1}(P)\) is the normal inverse cumulative distribution function evaluated at probability \(P = 3/4\).

max

The maximum of the unmasked pixel values within the aperture.

mean

The mean of the unmasked pixel values within the aperture.

median

The median of the unmasked pixel values within the aperture.

min

The minimum of the unmasked pixel values within the aperture.

mode

The mode of the unmasked pixel values within the aperture.

The mode is estimated as (3 * median) - (2 * mean).

moments

Spatial moments up to 3rd order of the source.

moments_central

Central moments (translation invariant) of the source up to 3rd order.

n_apertures

The number of positions in the input aperture.

orientation

The angle between the x axis and the major axis of the 2D Gaussian function that has the same second-order moments as the source.

The angle increases in the counter-clockwise direction.

properties

A sorted list of built-in source properties.

semimajor_sigma

The 1-sigma standard deviation along the semimajor axis of the 2D Gaussian function that has the same second-order central moments as the source.

semiminor_sigma

The 1-sigma standard deviation along the semiminor axis of the 2D Gaussian function that has the same second-order central moments as the source.

sky_centroid

The sky coordinate of the centroid of the unmasked pixels within the aperture, returned as a SkyCoord object.

The output coordinate frame is the same as the input wcs.

None if wcs is not input.

sky_centroid_icrs

The sky coordinate in the International Celestial Reference System (ICRS) frame of the centroid of the unmasked pixels within the aperture, returned as a SkyCoord object.

None if wcs is not input.

std

The standard deviation of the unmasked pixel values within the aperture.

sum

The sum of the unmasked data values within the aperture.

\[F = \sum_{i \in A} I_i\]

where \(F\) is sum, \(I_i\) is the background-subtracted data, and \(A\) are the unmasked pixels in the aperture.

Non-finite pixel values (NaN and inf) are excluded (automatically masked).

sum_aper_area

The total area of the unmasked pixels within the aperture using the input sum_method aperture mask method.

sum_err

The uncertainty of sum , propagated from the input error array.

sum_err is the quadrature sum of the total errors over the unmasked pixels within the aperture:

\[\Delta F = \sqrt{\sum_{i \in A} \sigma_{\mathrm{tot}, i}^2}\]

where \(\Delta F\) is the sum, \(\sigma_{\mathrm{tot, i}}\) are the pixel-wise total errors (error), and \(A\) are the unmasked pixels in the aperture.

Pixel values that are masked in the input data, including any non-finite pixel values (NaN and inf) that are automatically masked, are also masked in the error array.

var

The variance of the unmasked pixel values within the aperture.

xcentroid

The x coordinate of the centroid.

The centroid is computed as the center of mass of the unmasked pixels within the aperture.

ycentroid

The y coordinate of the centroid.

The centroid is computed as the center of mass of the unmasked pixels within the aperture.

Methods Documentation

copy()[source]

Return a deep copy of this object.

get_id(id_num)[source]

Return a new ApertureStats object for the input ID number only.

Parameters:
id_numint

The aperture ID number.

Returns:
resultApertureStats

A new ApertureStats object containing only the source with the input ID number.

get_ids(id_nums)[source]

Return a new ApertureStats object for the input ID numbers only.

Parameters:
id_numslist, tuple, or ndarray of int

The aperture ID number(s).

Returns:
resultApertureStats

A new ApertureStats object containing only the sources with the input ID numbers.

to_table(columns=None)[source]

Create a QTable of source properties.

Parameters:
columnsstr, list of str, None, optional

Names of columns, in order, to include in the output QTable. The allowed column names are any of the ApertureStats properties. If columns is None, then a default list of scalar-valued properties (as defined by the default_columns attribute) will be used.

Returns:
tableQTable

A table of sources properties with one row per source.