Object detection and measurement

This page covers object detection and flux measurement (aperture, optimal extraction, SEP-based, and PSF photometry). For photometric calibration (zero points, color terms) see Photometric calibration, for PSF model construction see Point Spread Function (PSF) models, and for detection flag definitions see Object flags.

Detecting objects on the image

STDPipe provides three functions for detecting objects on the image:

The first two have mostly identical signatures and arguments, but differ in minor details like meaning of returned detection flags etc.

The detection in both cases is based on building the noise model through (grid-based) background and background rms estimation, and then extracting the groups of connected pixels above some pre-defined threshold. Optionally, before the thresholding, the image may be smoothed with some small kernel in order to improve the detection of fainter objects. We recommend checking the SExtractor documentation to better understand the concepts of it.

Also, the routines may automatically reject the objects detected too close to frame edge, and probably truncated or improperly measured - this is controlled by edge argument. Less significant detections - with \(\frac1{\rm magerr} < {\rm S/N}\) - may also be automatically rejected by providing the sn argument.

Both routines return the results as a standard Astropy Table, ordered by the object brightness. The table contains at least the following columns:

  • x, y, xerr, yerr - positions of the objects in image coordinates

  • ra, dec - positions of the objects in sky coordinates according to provided WCS

  • a, b, theta - object ellipse semimajor and semiminor axes, as well as its position angle in image coordinates

  • flux, fluxerr - measured flux and its error in a circular aperture

  • mag, magerr - instrumental magnitude computed as -2.5*log10(flux), and its error computed as 2.5/log(10)*fluxerr/flux

  • flags - detection flags

  • fwhm - per-object FWHM estimation

For get_objects_sextractor(), the output may also contain columns related to PSF photometry (x_psf, y_psf, flux_psf, fluxerr_psf, mag_psf, magerr_psf, chi2_psf, spread_model, spreaderr_model), as well as any additional measuremet parameter requested through extra_params argument.

Attention

STDPipe (as well as SEP library) uses pixel coordinate convention with (0, 0) as the origin - that differs from SExtractor that uses (1, 1) as origin of coordinates! So the routine transparently converts x and y (as well as x_psf and y_psf) in the output from SExtractor to proper origin, so that the coordinates of objects detected by both routines are in the same system. However, if you manually add some pixel coordinate parameter to the output through extra_params argument, they will not be appropriately adjusted!

Below are some examples of object detection.

# We will detect objects using SExtractor and get their measurements
# in apertures with 3 pixels radius.
# We will also ignore anything closer than 10 pixels to image edge
obj = photometry.get_objects_sextractor(image, mask=mask, aper=3.0, gain=gain, edge=10)
print(len(obj), 'objects detected')

# Rough estimation of average FWHM of detected objects, taking into account
# only unflagged (e.g. not saturated) ones
fwhm = np.median(obj['fwhm'][obj['flags'] == 0])
print('Average FWHM is %.1f pixels' % fwhm)

Next example shows how to receive extra columns, as well as checkimages, from SExtractor

# Now we will also get segmentation map and a column with object numbers correspoding to this map
obj,segm = photometry.get_objects_sextractor(image, mask=mask,
             aper=3.0, gain=gain, edge=10, extra_params=['NUMBER'],
             checkimages=['SEGMENTATION'], verbose=True)

# We may now e.g. completely mask the footprints of objects having masked pixels
fmask = np.zeros_like(mask)
for _ in obj[(obj['flags'] & 0x100) > 0]:
    fmask |= segm == _['NUMBER']

Finally, using SExtractor star/galaxy separators - CLASS_STAR and SPREAD_MODEL. The latter require PSF model - so we will build it using stdpipe.psf.run_psfex() documented elsewhere

# Get PSF model and store it to temporary file
psf_model = psf.run_psfex(image, mask=mask, order=0, gain=gain, psffile='/tmp/psf.psf', verbose=True)

# Run SExtractor
# You should provide correct path to default.nnw file from SExtractor installation on your system!
obj = photometry.get_objects_sextractor(image, mask=mask, edge=10, wcs=wcs,
             aper=3.0, extra_params=['CLASS_STAR', 'NUMBER'],
             extra={'SEEING_FWHM':fwhm, 'STARNNW_NAME':'/Users/karpov/opt/miniconda3/envs/stdpipe//share/sextractor/default.nnw'},
             psf='/tmp/psf.psf', verbose=True)

for i,cand in enumerate(obj):
    print('Candidate %d at x/y = %.1f %.1d and RA/Dec = %.4f %.4f' %
             (i, cand['x'], cand['y'], cand['ra'], cand['dec']))

    print('SPREAD_MODEL = %.3f +/- %.3f, CLASS_STAR = %.2f' %
             (cand['spread_model'], cand['spreaderr_model'], cand['CLASS_STAR']))

Note

The most important problem with object detection using these routines is handling of blended objects, as the codes we are using can’t properly deblend close groups, except for simplest cases.

stdpipe.photometry.get_objects_sextractor(image, header=None, mask=None, mask_detect=None, err=None, thresh=2.0, aper=3.0, r0=0.0, gain=1, edge=0, minarea=5, wcs=None, sn=3.0, bg_size=None, sort=True, reject_negative=True, mask_to_nans=True, checkimages=[], extra_params=[], extra={}, psf=None, catfile=None, _workdir=None, _tmpdir=None, _exe=None, verbose=False)[source]

Thin wrapper around the SExtractor binary.

Processes the image with optional mask and noise map and returns the list of detected objects, optionally with SExtractor checkimages.

See the SExtractor documentation for details. Detection flags are documented at https://sextractor.readthedocs.io/en/latest/Flagging.html#extraction-flags-flags. Objects with masked pixels in their footprint additionally get flag 0x100.

Parameters:
imagenumpy.ndarray

Input 2-D image.

headerastropy.io.fits.Header, optional

Image header.

masknumpy.ndarray of bool, optional

Pixel mask (True = masked). Used for SExtractor flagging (FLAG_IMAGE / IMAFLAGS_ISO). Objects overlapping this mask get flag 0x100.

mask_detectnumpy.ndarray of bool, optional

Detection mask (True = exclude). When mask_to_nans is True, these pixels (plus non-finite pixels) are set to NaN before running SExtractor.

errnumpy.ndarray, optional

Noise map.

threshfloat

Detection threshold in sigmas (DETECT_THRESH).

aperfloat or list of float

Aperture radius in pixels. If a list, flux is measured for each aperture.

r0float

Smoothing kernel sigma (pixels, i.e. FWHM / 2.355).

gainfloat

Detector gain, e-/ADU.

edgeint

Reject objects closer to edge than this many pixels.

minareaint

Minimum object area (DETECT_MINAREA).

wcsastropy.wcs.WCS, optional

WCS for sky coordinate assignment.

snfloat

Minimum S/N for output.

bg_sizeint, optional

Background grid size (BACK_SIZE).

sortbool

Sort detections by decreasing brightness.

reject_negativebool

Reject objects with negative flux.

mask_to_nansbool

Replace detection-masked pixels with NaN before running SExtractor. Default True.

checkimageslist of str

SExtractor checkimage types to return (e.g. 'BACKGROUND', 'SEGMENTATION').

extra_paramslist of str

Extra SExtractor output parameters (see sex -dp).

extradict

Extra SExtractor configuration parameters (see sex -dd).

psfstr, optional

Path to PSFEx PSF model file for PSF photometry.

catfilestr, optional

Copy the output catalogue to this path.

_workdirstr, optional

Keep temporary files in this directory (for debugging).

_tmpdirstr, optional

Create temp directory inside this path.

_exestr, optional

Full path to SExtractor executable.

verbosebool or callable

Print progress messages.

Returns:
astropy.table.Table or list

Table of detected objects. If checkimages is non-empty, returns [table, checkimage1, checkimage2, ...].

stdpipe.photometry.get_objects_sep(image, header=None, mask=None, mask_detect=None, err=None, thresh=4.0, aper=3.0, bkgann=None, r0=0.5, gain=1, edge=0, minnthresh=2, minarea=5, relfluxradius=2.0, wcs=None, bg_size=64, fwhm=False, optimal=False, group_sources=True, group_radius_factor=1.0, group_halo_factor=1.2, centroid=False, centroid_sigma_scale=1.0, centroid_maxstep_scale=1.0, centroid_maxshift_scale=0.5, use_mask_large=False, subtract_bg=True, npix_large=100, sn=10.0, get_segmentation=False, deblend_fwhm=0, deblend_method='watershed', clip_sigma=3.0, clip_iters=5, fwhm_spatial_order=0, verbose=True, **kwargs)[source]

Object detection and aperture/optimal photometry using SEP.

Signature is as similar as possible to get_objects_sextractor().

Algorithm: background estimation → object detection with optional smoothing and deblending (watershed by default) → optional windowed centroiding → edge rejection → per-object FWHM estimation → aperture or optimal extraction photometry → S/N quality cut.

Detection flags are documented at https://sep.readthedocs.io/en/v1.1.x/reference.html — they differ from SExtractor flags.

Parameters:
imagenumpy.ndarray

Input 2-D image.

headerastropy.io.fits.Header, optional

Image header.

masknumpy.ndarray of bool, optional

Pixel mask (True = masked). Used for background estimation and photometry. Objects whose detection footprint overlaps this mask are flagged with 0x100 but not excluded from detection.

mask_detectnumpy.ndarray of bool, optional

Detection mask (True = exclude from detection). Combined (OR) with non-finite pixel mask. Use this to exclude entire regions (e.g. bad columns, satellite trails). If None, only non-finite pixels are masked during detection.

errnumpy.ndarray, optional

Noise map.

threshfloat

Detection threshold in sigmas above local background.

aperfloat

Aperture radius in pixels, or in FWHM units when fwhm is provided/estimated.

bkganntuple of float, optional

Background annulus (r_in, r_out) in pixels (or FWHM units). Uses arithmetic mean of unmasked pixels inside the annulus. If None, the global background model is used.

r0float

Smoothing kernel sigma (pixels) for detection.

gainfloat

Detector gain, e-/ADU.

edgeint

Reject objects closer to image edge than this many pixels.

minnthreshint

Minimum number of pixels above threshold per detection.

minareaint

Minimum number of pixels in an object footprint.

relfluxradiusfloat

Multiplier for rough FWHM when computing flux_radius search radius (internal use).

wcsastropy.wcs.WCS, optional

WCS for sky coordinate assignment.

bg_sizeint

Background grid cell size in pixels.

fwhmbool or float

FWHM handling. False (default) — per-object FWHM for output only. True — estimate global FWHM and scale aper/bkgann (interpreted as FWHM units). Numeric — use this value directly and scale aper/bkgann.

optimalbool

Use optimal extraction (sep.sum_circle_optimal) instead of aperture photometry. Requires FWHM. SEP 1.4+ only.

group_sourcesbool

Grouped optimal extraction for overlapping sources. Default True (recommended). SEP 1.4+ only.

group_radius_factorfloat

Connectivity scale for grouped optimal extraction. The default 1.0 groups sources whose apertures overlap.

group_halo_factorfloat or None

Local context halo scale for grouped optimal extraction. The default 1.2 keeps a modest expanded solve context without widening source connectivity. None lets SEP use group_radius_factor.

centroidbool

Refine positions via sep.winpos. Default False (windowed positions can be biased in crowded fields). With SEP 1.4+, uses maxstep = 0.2 * FWHM.

centroid_sigma_scalefloat

Multiplicative scale applied to the Gaussian sigma passed to sep.winpos. The default 1.0 keeps the historical behavior sigma = FWHM / 2.355. Values below 1 use a tighter window; values above 1 use a broader window.

centroid_maxstep_scalefloat

Multiplicative scale applied to the SEP 1.4+ maxstep limit used by sep.winpos. The default 1.0 keeps the historical cap maxstep = 0.2 * FWHM. Values below 1 restrict the allowed centroid motion more strongly; values above 1 loosen it.

centroid_maxshift_scalefloat

Multiplicative scale applied to the SEP 1.4+ maxshift cap on the total centroid drift across all winpos iterations (in units of FWHM). The default 0.5 means a centroid can move at most 0.5 * FWHM from its initial position, complementing the per-iteration maxstep cap. Set to a large value to effectively disable. Ignored when SEP is older than 1.4 or when FWHM is not available.

use_mask_largebool

Filter out objects with footprints larger than npix_large.

npix_largeint

Pixel-count threshold for large-object rejection.

subtract_bgbool

Subtract background before detection (default True).

snfloat

Minimum S/N (magerr < 1/sn) for output.

get_segmentationbool

Also return the segmentation map and add seg_id column.

deblend_fwhmfloat

Fixed Gaussian FWHM (pixels) for deterministic deblending. 0 = adaptive shapes with stochastic assignment. SEP 1.4+ only.

deblend_method{‘watershed’, ‘threshold’}

Deblending algorithm. SEP 1.4+ only.

clip_sigmafloat

Sigma-clipping threshold for annulus background. SEP 1.4+ only.

clip_itersint

Max sigma-clipping iterations for annulus background. 0 = disable. SEP 1.4+ only.

fwhm_spatial_orderint

When FWHM is auto-estimated (fwhm=True or via the optimal path), fit a 2-D polynomial FWHM(x, y) of this order and apply it per source. 0 (default) keeps the historical scalar behaviour. Per-source values are used for aperture/background-annulus scaling and for the optimal-extraction PSF width; the windowed centroider still uses the scalar median. Ignored when a numeric fwhm is supplied directly. See FWHMMap.

verbosebool or callable

Print progress messages.

Returns:
astropy.table.Table or tuple

Table of detected objects. If get_segmentation is True, returns (table, segmentation_map).

Columns: x, y, xerr, yerr, flux, fluxerr, mag, magerr, flags, ra, dec, bg, fwhm, a, b, theta, and optionally seg_id.

Metadata: aper, bkgann, optimal, group_sources, fwhm_phot.

Notes

  • FWHM estimation prefers SEP’s built-in method (ported from SExtractor). Falls back to 2nd moments for older SEP versions.

  • Optimal extraction assumes a Gaussian PSF.

  • Grouped extraction fits nearby sources simultaneously — dramatically improves accuracy for close pairs.

  • Watershed deblending often gives better completeness for pairs closer than ~2 FWHM.

Detection using photutils

The get_objects_photutils() function provides a pure-Python alternative using the photutils library. It offers three detection algorithms:

  • segmentation (default) - Image segmentation with connected pixels above threshold. Supports optional multi-level deblending for crowded fields.

  • dao - DAOStarFinder algorithm using Gaussian convolution, optimized for stellar fields with known FWHM.

  • iraf - IRAFStarFinder, similar to DAOStarFinder but with IRAF-compatible behavior.

The segmentation method with deblending is particularly useful for crowded fields where sources overlap:

# Detect objects using photutils segmentation with deblending
obj = photometry.get_objects_photutils(image, mask=mask, thresh=3.0,
             method='segmentation', deblend=True, aper=3.0, edge=10)
print(len(obj), 'objects detected')

# For stellar fields, DAOStarFinder may be faster
obj = photometry.get_objects_photutils(image, mask=mask, thresh=5.0,
             method='dao', fwhm=3.0, aper=3.0)
stdpipe.photometry.get_objects_photutils(image, header=None, mask=None, mask_detect=None, err=None, thresh=2.0, method='segmentation', deblend=True, minarea=5, saturation=None, npixels=5, nlevels=32, contrast=0.001, connectivity=8, fwhm=3.0, sharplo=0.2, sharphi=1.0, roundlo=-1.0, roundhi=1.0, aper=3.0, bkgann=None, bg_size=64, subtract_bg=True, edge=0, sn=3.0, wcs=None, get_segmentation=False, verbose=False, **kwargs)[source]

Detect sources in image using photutils detection algorithms.

This function provides an alternative to get_objects_sep() using photutils instead of SEP. It supports multiple detection methods including segmentation-based detection with deblending and point source detection using DAOStarFinder or IRAFStarFinder.

Parameters:
imagenumpy.ndarray

2D image array

headerastropy.io.fits.Header, optional

FITS header for WCS information

masknumpy.ndarray, optional

Boolean mask (True = masked pixel). Used for background estimation, photometry, and flagging (0x100). Not used for detection itself.

mask_detectnumpy.ndarray, optional

Boolean mask for detection (True = exclude from detection). OR-ed with non-finite pixel mask. If None, only non-finite pixels are excluded from detection.

errnumpy.ndarray, optional

Error/uncertainty map. If None, uses background RMS

threshfloat, optional

Detection threshold in sigma units (default: 2.0)

methodstr, optional

Detection method: ‘segmentation’, ‘dao’, or ‘iraf’ (default: ‘segmentation’)

deblendbool, optional

Apply deblending for segmentation method (default: True)

minareaint, optional

Minimum source area in pixels (default: 5)

saturationfloat, optional

Saturation threshold for flagging saturated sources. Sources with peak pixel values exceeding this threshold will have flag 0x004 set. If None (default), saturation detection is disabled.

npixelsint, optional

Minimum connected pixels for segmentation (default: 5)

nlevelsint, optional

Number of deblending levels (default: 32)

contrastfloat, optional

Deblending contrast threshold (default: 0.001)

connectivityint, optional

Pixel connectivity (4 or 8) for segmentation (default: 8)

fwhmfloat, optional

FWHM for StarFinder methods in pixels (default: 3.0)

sharplofloat, optional

Lower sharpness bound for StarFinder (default: 0.2)

sharphifloat, optional

Upper sharpness bound for StarFinder (default: 1.0)

roundlofloat, optional

Lower roundness bound for StarFinder (default: -1.0)

roundhifloat, optional

Upper roundness bound for StarFinder (default: 1.0)

aperfloat, optional

Aperture radius in pixels (default: 3.0)

bkganntuple, optional

Background annulus (inner, outer) radii in pixels

bg_sizeint, optional

Background estimation box size (default: 64)

subtract_bgbool, optional

Subtract background before detection (default: True)

edgeint, optional

Edge exclusion in pixels (default: 0)

snfloat, optional

Minimum S/N ratio (default: 3.0)

wcsastropy.wcs.WCS, optional

WCS for coordinate conversion. If None, extracted from header

get_segmentationbool, optional

Return segmentation map (only for method=’segmentation’) (default: False)

verbosebool or callable, optional

Verbose output (default: False)

**kwargs

Additional parameters passed to photutils functions

Returns:
tableastropy.table.Table

Table of detected sources with columns: - x, y: pixel coordinates (0-indexed) - xerr, yerr: position errors - flux, fluxerr: aperture flux and error - mag, magerr: instrumental magnitude and error - fwhm: FWHM estimate - a, b: semi-major and semi-minor axes - theta: position angle - bg: local background - flags: detection flags - ra, dec: WCS coordinates (if WCS available)

segmphotutils.segmentation.SegmentationImage, optional

Segmentation map (only if get_segmentation=True and method=’segmentation’)

Notes

The function follows the same workflow as get_objects_sep(): 1. Validate inputs and create mask from NaNs 2. Estimate background using Background2D 3. Optionally subtract background 4. Detect sources using selected method 5. Measure aperture photometry 6. Apply quality filters 7. Add WCS coordinates if available 8. Return table with standard columns

Detection methods: - ‘segmentation’: Uses detect_sources() with optional deblending.

Best for extended sources and crowded fields.

  • ‘dao’: Uses DAOStarFinder for point sources. Best for stellar fields with well-defined PSF.

  • ‘iraf’: Uses IRAFStarFinder for IRAF-compatible detection. Similar to ‘dao’ but with IRAF-style parameters.

Examples

Basic segmentation detection:

>>> obj = get_objects_photutils(image, thresh=3.0, aper=5.0)

Segmentation with deblending:

>>> obj = get_objects_photutils(
...     image, method='segmentation', deblend=True,
...     nlevels=32, contrast=0.001
... )

Point source detection:

>>> obj = get_objects_photutils(
...     image, method='dao', fwhm=3.0, thresh=5.0
... )

With background annulus:

>>> obj = get_objects_photutils(
...     image, aper=5.0, bkgann=(10.0, 15.0)
... )

Get segmentation map:

>>> obj, segm = get_objects_photutils(
...     image, method='segmentation', get_segmentation=True
... )

FWHM estimation

After detection, you typically need a robust estimate of the image FWHM for aperture sizing, optimal extraction, PSF photometry, and matching radius selection. STDPipe provides two functions for this:

Both use a sliding-window mode estimator that is resistant to outliers from galaxies and blends. When flux_radius is available (from get_objects_sep()), it is preferred as it is 2-3x more stable than Gaussian-core FWHM.

# Higher-level: automatic quality filtering
fwhm = photometry.estimate_fwhm_from_objects(obj, verbose=True)
print('FWHM: %.2f pixels' % fwhm)

# Lower-level: manual filtering + raw array
good = (obj['flags'] == 0) & (obj['magerr'] < 0.1)
fwhm = photometry.estimate_fwhm(obj['fwhm'], good=good)

Position-dependent FWHM

Wide-field images often show FWHM variation across the field from focus gradients, optical aberrations, or tracking. Passing spatial_order >= 1 fits a 2-D polynomial FWHM(x, y) to the per-object values (with a*b weights and sigma-clipping) and returns a FWHMMap callable instead of a scalar. The callable evaluates the model at any position, and float(fwhm_map) gives a scalar median summary for APIs that require one.

# Fit FWHM(x, y) as a 2-D polynomial
fwhm_map = photometry.estimate_fwhm_from_objects(
    obj, spatial_order=2, image_shape=image.shape, verbose=True,
)
print('median FWHM:', float(fwhm_map), 'pix')
print('FWHM at image centre:', fwhm_map(image.shape[1] / 2, image.shape[0] / 2))

If too few high-quality candidates survive for the requested order, the function transparently falls back to the scalar path.

get_objects_sep() exposes the same feature via the fwhm_spatial_order parameter. When set, aperture and background-annulus radii are scaled per source and the per-source width is passed to sep.sum_circle_optimal (the windowed centroider still uses the scalar summary). The fitted model is stored in obj.meta['fwhm_phot_model']:

obj = photometry.get_objects_sep(image, fwhm=True, fwhm_spatial_order=2,
                                 optimal=True, aper=1.5, bkgann=(3, 5))
fmap = obj.meta['fwhm_phot_model']  # FWHMMap, or absent for scalar order=0

The downstream measurement routines measure_objects() and measure_objects_sep() also accept a callable fwhm. The SEP backend broadcasts a per-source width and per-source aperture/bkgann arrays directly to SEP; the pure-Python backend uses the per-source width in ungrouped optimal extraction and falls back to the scalar median where the underlying photutils APIs require a single radius.

stdpipe.photometry.estimate_fwhm_from_objects(obj, *, snr_min=10.0, max_ellipticity=0.3, use_flags=True, fwhm_range=(1.0, 20.0), min_candidates=5, spatial_order=0, xy_cols=('x', 'y'), image_shape=None, verbose=False)[source]

Estimate global FWHM from an object table (SEP or SExtractor).

Extracts per-object FWHM values and builds a quality mask from the available catalog columns, then delegates to estimate_fwhm().

When a flux_radius column is present (half-light radius, as produced by get_objects_sep()), it is converted to an equivalent Gaussian FWHM via FWHM = 2 * flux_radius and used as the primary estimator. Half-light radii are 2–3× more stable than SEP’s Gaussian-core FWHM (relative IQR ~0.2 vs ~0.7) and far less susceptible to bias from sub-stellar contaminants, because the measurement integrates flux over many pixels rather than fitting a slope to a handful of core pixels.

When flux_radius is not available, the function falls back to the fwhm / FWHM_IMAGE column or moment-derived FWHM, and applies an a*b-weighted median correction if the mode appears biased (weighted median exceeds mode by > 10 %).

Parameters:
obj~astropy.table.Table or structured ndarray

Detection results. Recognized columns (checked in priority order): flux_radius (half-light radius), fwhm, FWHM_IMAGE, or derived from a/b (A_IMAGE/B_IMAGE) moments.

snr_minfloat or None

Minimum S/N for candidate selection (via magerr). Set to None to disable.

max_ellipticityfloat or None

Maximum ellipticity 1 - b/a. Set to None to disable.

use_flagsbool

If True, reject objects with nonzero flags.

fwhm_rangetuple of float

Passed through to estimate_fwhm().

min_candidatesint

Passed through to estimate_fwhm().

spatial_orderint

Polynomial order for a position-dependent FWHM model. 0 (default) preserves the historical behaviour and returns a scalar. >=1 attempts to fit a 2-D polynomial FWHM(x, y) to the per-object values and returns a FWHMMap callable. If too few candidates survive the quality cuts for the requested order, the function falls back to the scalar path.

xy_colstuple of str

Column names used for pixel coordinates when spatial_order >= 1. Defaults to ('x', 'y') (SEP-style). Use ('X_IMAGE', 'Y_IMAGE') for SExtractor output. Falls back to the scalar path if either column is missing.

image_shapetuple of int, optional

(H, W) image shape used to normalise pixel coordinates to roughly [-1, 1] when spatial_order >= 1. If not provided, the min/max of the candidate coordinates are used.

verbosebool or callable

If True, log diagnostic information about the estimation.

Returns:
float or FWHMMap

Scalar FWHM in pixels when spatial_order == 0. A FWHMMap callable otherwise, or a scalar on fallback. Returns np.nan if no usable candidates survive.

class stdpipe.photometry.FWHMMap(coeffs, order, x0, y0, sx, sy, median, n_used, fwhm_range)[source]

Position-dependent FWHM estimator backed by a 2-D polynomial.

Produced by estimate_fwhm_from_objects() when spatial_order >= 1. Coordinates are normalised internally to roughly [-1, 1] for numerical stability.

Instances are callable: fwhm_map(x, y) returns FWHM in pixels at the requested pixel coordinates (x, y may be scalars or broadcastable arrays; the result has the broadcast shape). Instances are also convertible via float(fwhm_map) to the scalar median summary, so they can be dropped into APIs that still expect a single number.

Attributes:
coeffsndarray

Flat array of polynomial coefficients, ordered by (i, j) with i + j <= order and j varied fastest.

orderint

Polynomial order.

x0, y0, sx, syfloat

Coordinate normalisation: xn = (x - x0) / sx.

medianfloat

Representative scalar FWHM (median of the fitted surface over the good candidates). Used by float(self).

n_usedint

Number of candidates retained after sigma-clipping.

fwhm_rangetuple of float

Output clipping range in pixels.

Methods

__call__(x, y)

Call self as a function.

__call__(x, y)[source]

Call self as a function.

stdpipe.photometry.estimate_fwhm(values, *, good=None, fwhm_range=(1.0, 20.0), min_candidates=5)[source]

Robust global FWHM estimate from per-object FWHM measurements.

Uses a sliding-window mode estimator (narrowest quartile interval) that is resistant to outliers from galaxies, blends, and non-Gaussian tails.

Parameters:
valuesarray-like

1-D array of per-object FWHM values in pixels.

goodarray-like of bool, optional

Boolean mask, True = use this object. The caller is responsible for applying any S/N, flag, or ellipticity cuts before calling, or encoding them in this mask.

fwhm_rangetuple of float

Hard clip range (lo, hi) in pixels. Values outside this range are excluded before the mode calculation.

min_candidatesint

Minimum number of surviving candidates required.

Returns:
float

Scalar FWHM estimate in pixels, or np.nan on failure.

More accurate photometric measurements

The photometric measurements returned by detection routines are sometimes not optimal. The stdpipe.photometry.measure_objects() function provides several improvements:

  • Local background estimation using annuli around each object instead of global background model

  • Centroiding to refine object positions before measurement

  • Optimal extraction (Naylor 1998) for ~10% S/N improvement on point sources

  • Grouped fitting for crowded fields with overlapping PSFs

Local background estimation

The detection routines use a globally estimated background model (built as a low-resolution spatial map and then interpolated to original pixels, see here). On rapidly varying backgrounds, you may expect better results from locally estimated background using annuli around every object:

# We will pass FWHM so that aperture and background radii are relative to it
# We will also reject all objects with measured S/N < 5
obj = photometry.measure_objects(obj, image, mask=mask, fwhm=fwhm, gain=gain,
             aper=1.0, bkgann=[5, 7], sn=5, verbose=True)
print(len(obj), 'objects properly measured')

Centroiding

You can refine object positions before photometry using iterative centroiding. This improves aperture placement, especially when initial positions from detection are approximate:

# Refine positions with 3 centroiding iterations before aperture photometry
obj = photometry.measure_objects(obj, image, mask=mask, fwhm=fwhm,
             aper=1.0, centroid_iter=3, verbose=True)

# Original positions are stored in x_orig, y_orig columns
print('Position shift:', np.median(np.hypot(obj['x'] - obj['x_orig'],
                                             obj['y'] - obj['y_orig'])))

Optimal extraction

For point sources, optimal extraction (Naylor 1998) provides ~10% improvement in signal-to-noise ratio by using PSF-weighted photometry instead of simple aperture sums:

# Use optimal extraction with Gaussian PSF based on FWHM
obj = photometry.measure_objects(obj, image, mask=mask,
             fwhm=fwhm, aper=1.5, optimal=True, verbose=True)

# Or provide a PSF model from PSFEx for better accuracy
psf_model = psf.run_psfex(image, mask=mask, gain=gain)
obj = photometry.measure_objects(obj, image, mask=mask,
             psf=psf_model, aper=1.5, optimal=True, verbose=True)

Optimal extraction adds quality metrics to the output:

  • chi2_optimal - chi-squared of the fit (lower is better)

  • norm_optimal - PSF normalization factor (should be ~1 for point sources)

  • npix_optimal - number of unmasked pixels used in the fit

Grouped optimal extraction for crowded fields

In crowded fields where source PSFs overlap, grouped optimal extraction simultaneously fits fluxes for nearby sources using weighted least squares. This properly accounts for flux sharing between overlapping PSFs:

# Grouped extraction - sources within 2*FWHM are fitted together
obj = photometry.measure_objects(obj, image, mask=mask,
             fwhm=fwhm, aper=1.5, optimal=True,
             group_sources=True, grouper_radius=2*fwhm,
             verbose=True)

# Check group information
print('Sources in groups of 2+:', len(obj[obj['group_size'] > 1]))

This adds group_id and group_size columns to identify which sources were fitted together.

stdpipe.photometry.measure_objects(obj, image, aper=3, bkgann=None, bkg_order=0, fwhm=None, psf=None, optimal=False, group_sources=True, grouper_radius=None, mask=None, bg=None, err=None, gain=None, bg_size=64, sn=None, centroid_iter=0, centroid_method='com', keep_negative=True, get_bg=False, verbose=False)[source]

Photometry at the positions of already detected objects.

Supports both standard aperture photometry and optimal extraction that provides ~10% S/N improvement for point sources (Naylor 1998).

It will estimate and subtract the background unless external background estimation (bg) is provided, and use user-provided noise map (err) if requested.

The results may optionally be filtered to drop detections with low signal to noise ratio if sn is set and positive. It will also filter out events with negative flux unless keep_negative is True.

Parameters:
objastropy.table.Table

Table with initial object detections to be measured.

imagendarray

Input image as a NumPy array.

aperfloat

Circular aperture radius in pixels for flux measurement. For optimal extraction, this is the clipping radius.

bkganntuple of float or None

Background annulus (inner, outer radii) for local background estimation. If None, global background model is used instead.

bkg_orderint

Polynomial order for local background fitting. 0 = constant (mean), 1 = plane (linear gradient, recommended), 2 = quadratic surface. Only used if bkgann is set.

fwhmfloat, callable, or None

If provided, aper and bkgann are measured in units of FWHM. Also used to define Gaussian PSF for optimal extraction if psf is not provided. A callable (e.g. a stdpipe.photometry.FWHMMap) is accepted and evaluated at each source position for the optimal-extraction PSF width; aperture/bkgann scaling, PSF-weighted centroiding, and local background annuli still use the scalar summary float(fwhm).

psfdict or None

PSF model for optimal extraction and PSF-weighted centroiding. Can be a dict from psf.run_psfex(), psf.load_psf(), or psf.create_psf_model(). If None, a Gaussian PSF is created from the fwhm parameter.

optimalbool

If True, use optimal extraction instead of aperture photometry. Requires either psf or fwhm to define the PSF profile.

group_sourcesbool

If True and optimal=True, use grouped optimal extraction for overlapping sources. Fits nearby sources simultaneously for better accuracy in crowded fields.

grouper_radiusfloat or None

Radius in pixels for grouping nearby sources. Sources within this distance are fitted simultaneously. If None, defaults to 2*aper. Only used if group_sources=True.

maskndarray of bool or None

Image mask (True values are masked).

bgndarray or None

If provided, use this background instead of automatically computed one. Must have the same shape as input image.

errndarray or None

Image noise map to use instead of automatically computed one.

gainfloat or None

Image gain in e-/ADU, used to build image noise model.

bg_sizeint

Background grid size in pixels.

snfloat or None

Minimal S/N ratio. If set, measurements with magnitude errors exceeding 1/sn will be discarded.

centroid_iterint

Number of centroiding iterations to run before photometry.

centroid_methodstr

Centroiding method: ‘com’ (center-of-mass) or ‘psf’ (PSF-weighted). PSF-weighted is more accurate but degrades with heavy random masking (>20%). Requires psf or fwhm parameter.

keep_negativebool

If False, measurements with negative fluxes will be discarded.

get_bgbool

If True, also return estimated background and background noise images.

verbosebool or callable

Whether to show verbose messages. May be boolean or a print-like function.

Returns:
objastropy.table.Table

Copy of original table with flux, fluxerr, mag and magerr columns replaced with measured values.

bg_imagendarray

Background image (only returned if get_bg=True).

err_imagendarray

Error image (only returned if get_bg=True).

Notes

Quality flags set in the flags column:

  • 0x200: At least one aperture pixel is masked

  • 0x400: Invalid position or local background estimation failed

  • 0x800: Optimal extraction failed (NaN result)

  • 0x1000: Poor fit quality (chi2 > 1000, numerical instability)

SEP-based photometry

The stdpipe.photometry_measure.measure_objects_sep() function provides an alternative measurement backend using SEP 1.4+ features. It offers sigma-clipped local background estimation and grouped optimal extraction implemented in C for better performance.

from stdpipe.photometry_measure import measure_objects_sep, _HAS_SEP_OPTIMAL

if _HAS_SEP_OPTIMAL:
    # Aperture photometry with robust local background
    obj = measure_objects_sep(obj, image, aper=1.5, fwhm=fwhm,
                 bkgann=(3.0, 5.0), sn=5, verbose=True)

    # Grouped optimal extraction for crowded fields
    obj = measure_objects_sep(obj, image, aper=1.5, fwhm=fwhm,
                 optimal=True, group_sources=True, sn=5, verbose=True)
stdpipe.photometry_measure.measure_objects_sep(obj, image, aper=3, bkgann=None, fwhm=None, psf=None, optimal=False, group_sources=True, group_factor=2.0, group_radius_factor=1.0, group_halo_factor=1.2, maxiter=20, fit_positions=True, fit_radius=0.0, damp_snthresh=0.0, mask=None, bg=None, err=None, gain=None, bg_size=64, sn=None, centroid_iter=0, centroid_psf=None, keep_negative=True, get_bg=False, clip_sigma=3.0, clip_iters=5, verbose=False)[source]

Photometry at the positions of already detected objects using SEP routines.

Uses SEP’s built-in features for optimal extraction with sigma-clipped background (via sum_circle_optimal), PSF fitting photometry (via sep.psf_fit), and iterative centroiding (via winpos).

Only available if SEP version 1.4+ with these features is installed.

Parameters:
objastropy.table.Table

Table with initial object detections to be measured.

imagendarray

Input image as a NumPy array.

aperfloat

Circular aperture radius in pixels for flux measurement.

bkganntuple of float or None

Background annulus (inner, outer radii) for local background. For optimal extraction, SEP handles this internally with sigma-clipping. For aperture photometry, sep.stats_circann() is used.

fwhmfloat, callable, or None

If provided, aper and bkgann are measured in units of FWHM. Also used for Gaussian PSF in optimal extraction, PSF fitting, and centroiding. A callable (e.g. a stdpipe.photometry.FWHMMap) is accepted and evaluated at each source position: per-source aperture/bkgann arrays feed the SEP photometry routines and the per-source width feeds sep.sum_circle_optimal. The scalar summary float(fwhm) is used for the SEP windowed centroider (which requires a scalar sigma).

psfdict, sep.PSF, or None

PSF model for PSF fitting photometry. When provided, PSF fitting is used instead of aperture or optimal extraction. Can be a PSFEx dict from stdpipe.psf.run_psfex() or a sep.PSF object.

optimalbool

If True, use optimal extraction via sep.sum_circle_optimal(). Requires fwhm. Ignored when psf is provided.

group_sourcesbool

If True, use grouped fitting for optimal or PSF photometry.

group_factorfloat

Local fitting halo factor for grouped PSF fitting. Only used for PSF fitting when group_sources=True; SEP determines PSF-fit connectivity from direct fit-support overlap.

group_radius_factorfloat

Connectivity scale for grouped optimal extraction. Only used for optimal extraction when group_sources=True. The default 1.0 groups sources whose apertures overlap.

group_halo_factorfloat or None

Local context halo scale for grouped optimal extraction. The default 1.2 keeps a modest expanded solve context without widening source connectivity. None lets SEP use group_radius_factor.

maxiterint

Maximum number of PSF fitting iterations.

fit_positionsbool

If True, fit source positions during PSF fitting.

fit_radiusfloat

If > 0, only pixels within this radius (in pixels) of the source center participate in PSF fitting. 0 means use the full PSF stamp. Values of 2-3x FWHM reduce scatter in crowded fields.

damp_snthreshfloat

S/N threshold for damping PSF fit updates.

maskndarray of bool or None

Image mask (True values are masked).

bgndarray or None

If provided, use this background instead of automatically computed one.

errndarray or None

Image noise map.

gainfloat or None

Image gain in e-/ADU, used to build image noise model.

bg_sizeint

Background grid size in pixels.

snfloat or None

Minimal S/N ratio for filtering.

centroid_iterint

Number of centroiding iterations (uses SEP’s built-in iteration).

centroid_psfdict, sep.PSF, or None

PSF model for centroiding. If None, uses psf if available, otherwise Gaussian windowed centroiding.

keep_negativebool

If False, measurements with negative fluxes will be discarded.

get_bgbool

If True, also return estimated background and noise images.

clip_sigmafloat

Sigma value for clipping when bkgann is provided.

clip_itersint

Maximum number of clipping iterations when bkgann is provided. Set to 0 to disable clipping.

verbosebool or callable

Whether to show verbose messages. May be boolean or print-like function.

Returns:
objastropy.table.Table

Copy of table with flux, fluxerr, mag and magerr columns from SEP measurements. When psf is provided, also includes x_psf, y_psf (fitted positions), chi2_psf, niter_psf, and flags_psf columns.

bg_imagendarray

Background image (only returned if get_bg=True).

err_imagendarray

Error image (only returned if get_bg=True).

Notes

Quality flags set in the flags column:

  • 0x200: At least one aperture pixel is masked (aperture/optimal paths)

  • 0x400: Invalid position

  • 0x800: Optimal extraction failed (NaN result)

  • 0x1000: PSF fit returned non-zero quality flag, or PSF fit failed (NaN result)

  • 0x2000: Large centroid shift during PSF fit (>1 pixel)

PSF photometry

For the most accurate flux measurements, especially in crowded fields, PSF fitting photometry is available through two backends:

Both support grouped fitting for crowded fields and return fitted positions alongside fluxes. See Point Spread Function (PSF) models for detailed documentation on PSF model construction, position-dependent PSF, grouped fitting, and empirical ePSF models.

from stdpipe import photometry_psf

# PSF photometry with automatic Gaussian PSF (photutils backend)
result = photometry_psf.measure_objects_psf(obj, image, fwhm=fwhm)

# With PSFEx model and grouped fitting
result = photometry_psf.measure_objects_psf(obj, image, psf=psf_model,
             group_sources=True, grouper_radius=10.0)

# DAOPHOT backend (requires PyRAF/IRAF)
from stdpipe import photometry_iraf
result = photometry_iraf.measure_objects_psf(obj, image, fwhm=fwhm)

Hybrid aperture photometry with PSF-based deblending

stdpipe.photometry_measure.measure_aperture_deblended() combines aperture photometry on the target with PSF-modelled subtraction of any neighbours that fall inside the aperture. It gives PSF-quality crowding handling while keeping the target flux as an aperture sum, which is robust to PSF-model mismatch on the target itself: any model error feeds only into the small neighbour-subtraction term and partially cancels between the modelled and self-flux contributions.

The function expects each detection to come with both an upstream aperture flux (flux_seed) and a PSF-fit flux (flux_psf); these are produced by get_objects_sep() followed by measure_objects_sep() with psf=.... The per-source aperture-correction ratio flux_psf / flux_seed is then fit as a smooth field across the frame (LOESS by default via stdpipe.smoothing.fit_vector_field_2d()) and used to scale the PSF stamps placed in the neighbour model.

from stdpipe import photometry, photometry_measure, psf

# Detect, build a PSF model and run PSF fitting upstream
obj = photometry.get_objects_sep(image, mask=mask, thresh=3.0, fwhm=True)
fwhm = float(obj.meta['fwhm_phot'])
psf_model = psf.run_psfex(image, mask=mask, gain=gain, verbose=True)
meas = photometry_measure.measure_objects_sep(
    obj, image, psf=psf_model, fwhm=fwhm,
    err=err_map, gain=gain, mask=mask, group_sources=True,
)

# Hybrid aperture photometry on a subset of the detections
target = meas['flags'] == 0
res = photometry_measure.measure_aperture_deblended(
    image, meas['x'], meas['y'], aper=1.5, fwhm=fwhm, psf=psf_model,
    flux_seed=meas['flux_seed'], flux_psf=meas['flux'],
    flux_psf_err=meas['fluxerr'],
    target=target,
    err=err_map, mask=mask, gain=gain,
    propagate_neighbour_errors=True,
    diagnostics=True,
)

The output table has one row per measured target (x, y, flux, fluxerr, mag, magerr, flags). Setting diagnostics=True adds the per-target flux_ap_raw, flux_model, flux_total_model, flux_ratio_model, self_frac, and flux_neighbour_err columns.

With propagate_neighbour_errors=True (the default), the reported fluxerr includes the per-near-neighbour PSF-flux error scaled by its overlap with the target aperture, summed in quadrature on top of the bare aperture noise. Disable for speed when neighbour errors are negligible.

stdpipe.photometry_measure.measure_aperture_deblended(image, x, y, aper, psf, *, flux_seed, flux_psf, flux_psf_err=None, target=None, fwhm=None, err=None, mask=None, gain=None, bg=None, aperture_correction='ratio_field', correction_field_kwargs=None, ratio_clip=(0.5, 3.0), outlier_clip=5.0, min_for_field_fit=6, bad_flag=4096, propagate_neighbour_errors=True, overlap_subpixel=4, diagnostics=False, verbose=False)[source]

Aperture photometry with PSF-based neighbour subtraction.

For every position in (x[target], y[target]) measure the aperture flux on image, subtract the modelled contribution of every other detection inside that aperture (using each one’s PSF-scaled total flux), and add back the target’s own modelled contribution. The result is an aperture flux that handles crowding the way PSF photometry would, while remaining robust to PSF-model mismatch on the target itself: PSF-model error feeds only into the small neighbour-subtraction term, partially cancelling between model_flux and self_flux.

The “total flux” used to scale each detection’s PSF stamp in the neighbour model is

total_i = flux_seed_i * ratio(x_i, y_i)

where ratio = flux_psf / flux_seed is fit as a smooth field across the frame from the per-source ratios. This provides a position-dependent aperture correction even when the PSF model has smoothly varying shape.

Parameters:
image2-D ndarray

Background-subtracted image (or pass the background separately via bg).

x, y(N,) array-like

Pixel positions of all detections in the region (targets + neighbours used only to model contamination).

aperfloat

Aperture radius. In image pixels by default; interpreted as a multiple of fwhm when fwhm is given (matching the convention of measure_objects()).

psfdict

PSF model as returned by stdpipe.psf.run_psfex(), stdpipe.psf.load_psf(), or stdpipe.psf.create_psf_model(). May be position-dependent.

flux_seed(N,) array-like

Per-source aperture flux from the upstream detection pass. Used to anchor the smooth aperture-correction ratio field.

flux_psf(N,) array-like

Per-source PSF-fitted flux from the upstream PSF-fitting pass. Used to anchor the ratio field and to scale neighbour PSF stamps in the model image.

flux_psf_err(N,) array-like, optional

Per-source PSF-fit flux errors. Required when propagate_neighbour_errors=True.

target(N,) bool mask or array of int indices, optional

Which detections to measure. Detections not in target still appear in the neighbour model. Default: every detection.

fwhmfloat, optional

Source FWHM in pixels. When given, aper is interpreted as a multiple of FWHM.

err2-D ndarray, optional

Per-pixel 1-sigma error map for the input image. Used by sep.sum_circle to compute the aperture noise.

mask2-D bool ndarray, optional

Bad-pixel mask. True pixels are excluded from aperture sums.

gainfloat, optional

Detector gain in e-/ADU; forwarded to sep.sum_circle for Poisson noise on the aperture sum.

bg2-D ndarray, optional

Background map. If given, it is subtracted from image before the aperture sums (so image itself need not be background- subtracted).

aperture_correction{‘ratio_field’, ‘fixed’, ‘none’}

How to map flux_seed to “total flux” when scaling neighbour PSF stamps. 'ratio_field' fits a smooth 2-D field of the per-source PSF/seed ratio (default; see correction_field_kwargs and min_for_field_fit); 'fixed' uses a single MAD-clipped median ratio; 'none' uses flux_seed as the total flux.

correction_field_kwargsdict, optional

Forwarded to stdpipe.smoothing.fit_vector_field_2d(). The default uses the LOESS backend with k=100 and per-axis scales chosen from the spread of the input positions, which is a quality upgrade over the rigid 2-D quadratic used in earlier prototypes.

ratio_clip(lo, hi)

Hard limits applied to the modelled ratio at every position to guard against extrapolation pathologies.

outlier_clipfloat

MAD multiple used to reject outliers from the per-source ratio sample before fitting the field.

min_for_field_fitint

Minimum number of clean ratios required to fit the smooth field; below this, the median ratio is used everywhere.

bad_flagint

Flag bit OR’ed into the output flags for measurements that could not be computed.

propagate_neighbour_errorsbool

If True (default), inflate fluxerr by the quadrature sum of each near neighbour’s PSF-flux error scaled by its overlap with the target aperture. Requires flux_psf_err.

overlap_subpixelint

Sub-pixel sampling factor for the per-neighbour aperture overlap integrals. 1 for a faster pixel-centre mask.

diagnosticsbool

If True, attach extra columns to the output table: flux_ap_raw, flux_model, flux_total_model, flux_ratio_model, self_frac, flux_neighbour_err.

verbosebool or callable

Boolean to enable default print logging, or a print-like callable.

Returns:
astropy.table.Table

One row per measured target, in the order they appeared in the input. Columns: x, y, flux, fluxerr, mag, magerr, flags. Diagnostic columns when diagnostics=True.