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:
stdpipe.photometry.get_objects_sextractor()- based on external SExtractor binarystdpipe.photometry.get_objects_sep()- based on Python SEP librarystdpipe.photometry.get_objects_photutils()- based on photutils library (pure Python, no compiled dependencies)
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 coordinatesra,dec- positions of the objects in sky coordinates according to provided WCSa,b,theta- object ellipse semimajor and semiminor axes, as well as its position angle in image coordinatesflux,fluxerr- measured flux and its error in a circular aperturemag,magerr- instrumental magnitude computed as-2.5*log10(flux), and its error computed as2.5/log(10)*fluxerr/fluxflags- detection flagsfwhm- 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_radiussearch 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.
Nonelets SEP usegroup_radius_factor.- centroidbool
Refine positions via
sep.winpos. Default False (windowed positions can be biased in crowded fields). With SEP 1.4+, usesmaxstep = 0.2 * FWHM.- centroid_sigma_scalefloat
Multiplicative scale applied to the Gaussian sigma passed to
sep.winpos. The default 1.0 keeps the historical behaviorsigma = 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+
maxsteplimit used bysep.winpos. The default 1.0 keeps the historical capmaxstep = 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+
maxshiftcap on the total centroid drift across all winpos iterations (in units of FWHM). The default 0.5 means a centroid can move at most0.5 * FWHMfrom its initial position, complementing the per-iterationmaxstepcap. 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_idcolumn.- 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=Trueor via the optimal path), fit a 2-D polynomialFWHM(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 numericfwhmis supplied directly. SeeFWHMMap.- 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:
estimate_fwhm_from_objects()— higher-level, works directly on detection tables with automatic quality filtering (S/N, flags, ellipticity)estimate_fwhm()— lower-level, operates on a raw array of FWHM values
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_radiuscolumn is present (half-light radius, as produced byget_objects_sep()), it is converted to an equivalent Gaussian FWHM viaFWHM = 2 * flux_radiusand 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_radiusis not available, the function falls back to thefwhm/FWHM_IMAGEcolumn or moment-derived FWHM, and applies ana*b-weighted median correction if the mode appears biased (weighted median exceeds mode by > 10 %).- Parameters:
- obj
~astropy.table.Tableor structured ndarray Detection results. Recognized columns (checked in priority order):
flux_radius(half-light radius),fwhm,FWHM_IMAGE, or derived froma/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.>=1attempts to fit a 2-D polynomialFWHM(x, y)to the per-object values and returns aFWHMMapcallable. 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]whenspatial_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.
- obj
- Returns:
- float or FWHMMap
Scalar FWHM in pixels when
spatial_order == 0. AFWHMMapcallable otherwise, or a scalar on fallback. Returnsnp.nanif 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()whenspatial_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,ymay be scalars or broadcastable arrays; the result has the broadcast shape). Instances are also convertible viafloat(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)withi + j <= orderandjvaried 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.nanon 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
snis set and positive. It will also filter out events with negative flux unlesskeep_negativeis 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
bkgannis set.- fwhmfloat, callable, or None
If provided,
aperandbkgannare measured in units of FWHM. Also used to define Gaussian PSF for optimal extraction ifpsfis not provided. A callable (e.g. astdpipe.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 summaryfloat(fwhm).- psfdict or None
PSF model for optimal extraction and PSF-weighted centroiding. Can be a dict from
psf.run_psfex(),psf.load_psf(), orpsf.create_psf_model(). If None, a Gaussian PSF is created from thefwhmparameter.- optimalbool
If True, use optimal extraction instead of aperture photometry. Requires either
psforfwhmto 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
psforfwhmparameter.- 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,magandmagerrcolumns 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
flagscolumn: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 (viasep.psf_fit), and iterative centroiding (viawinpos).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,
aperandbkgannare measured in units of FWHM. Also used for Gaussian PSF in optimal extraction, PSF fitting, and centroiding. A callable (e.g. astdpipe.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 feedssep.sum_circle_optimal. The scalar summaryfloat(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 asep.PSFobject.- optimalbool
If True, use optimal extraction via
sep.sum_circle_optimal(). Requiresfwhm. Ignored whenpsfis 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.
Nonelets SEP usegroup_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
psfif 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
bkgannis provided.- clip_itersint
Maximum number of clipping iterations when
bkgannis 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,magandmagerrcolumns from SEP measurements. Whenpsfis provided, also includesx_psf,y_psf(fitted positions),chi2_psf,niter_psf, andflags_psfcolumns.- 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
flagscolumn: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:
stdpipe.photometry_psf.measure_objects_psf()— photutils backend (pure Python, flexible PSF models)stdpipe.photometry_iraf.measure_objects_psf()— DAOPHOT backend (classic IRAF, requires PyRAF)
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 onimage, 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 betweenmodel_fluxandself_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_seedis 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
fwhmwhenfwhmis given (matching the convention ofmeasure_objects()).- psfdict
PSF model as returned by
stdpipe.psf.run_psfex(),stdpipe.psf.load_psf(), orstdpipe.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
targetstill appear in the neighbour model. Default: every detection.- fwhmfloat, optional
Source FWHM in pixels. When given,
aperis interpreted as a multiple of FWHM.- err2-D ndarray, optional
Per-pixel 1-sigma error map for the input image. Used by
sep.sum_circleto compute the aperture noise.- mask2-D bool ndarray, optional
Bad-pixel mask.
Truepixels are excluded from aperture sums.- gainfloat, optional
Detector gain in e-/ADU; forwarded to
sep.sum_circlefor Poisson noise on the aperture sum.- bg2-D ndarray, optional
Background map. If given, it is subtracted from
imagebefore the aperture sums (soimageitself need not be background- subtracted).- aperture_correction{‘ratio_field’, ‘fixed’, ‘none’}
How to map
flux_seedto “total flux” when scaling neighbour PSF stamps.'ratio_field'fits a smooth 2-D field of the per-source PSF/seed ratio (default; seecorrection_field_kwargsandmin_for_field_fit);'fixed'uses a single MAD-clipped median ratio;'none'usesflux_seedas the total flux.- correction_field_kwargsdict, optional
Forwarded to
stdpipe.smoothing.fit_vector_field_2d(). The default uses the LOESS backend withk=100and 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
flagsfor measurements that could not be computed.- propagate_neighbour_errorsbool
If True (default), inflate
fluxerrby the quadrature sum of each near neighbour’s PSF-flux error scaled by its overlap with the target aperture. Requiresflux_psf_err.- overlap_subpixelint
Sub-pixel sampling factor for the per-neighbour aperture overlap integrals.
1for 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
printlogging, 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 whendiagnostics=True.