stdpipe.photometry_background module

Background estimation routines for astronomical images.

This module provides various methods for estimating and subtracting background from astronomical images, including grid-based methods (SEP, photutils) and non-grid methods (percentile filtering, morphological opening).

stdpipe.photometry_background.get_background_percentile(image, mask=None, size=51, percentile=25.0, maxiters=3, sigma=3.0, use_mode=False)

Background estimation via percentile filtering (non-grid method).

Uses a sliding window to compute low percentile value at each pixel. Stars (positive outliers) are naturally above the percentile and thus rejected.

This is a truly non-grid method: operates on all pixels with sliding window, not limited by grid cell resolution.

Parameters:

imagendarray

Input image

maskndarray, optional

Boolean mask (True = masked pixels to exclude)

sizeint

Filter kernel size (should be larger than typical star size, e.g., 51-101 pixels). Larger values are more robust but may miss fine background structure.

percentilefloat

Percentile value to use (0-100). Typical values: - 25: 25th percentile (default, good for moderate star density) - 10: More aggressive star rejection (for crowded fields) - 50: Median (less star rejection but more robust to noise)

maxitersint

Number of sigma-clipping iterations to refine percentile estimate

sigmafloat

Sigma threshold for clipping

use_modebool

If True, estimate mode instead of percentile (experimental)

Returns:

backgroundndarray

Estimated background map

Notes

• Computationally intensive for large kernels (scales as O(N * kernel^2))

• Consider using smaller kernels (31-51) for large images

• For very crowded fields, use lower percentiles (10-15)

• Handles arbitrary background shapes (not limited to polynomial forms)

stdpipe.photometry_background.get_background_morphology(image, mask=None, size=25, iterations=1, smooth=True, smooth_size=None)

Background estimation via morphological opening (non-grid method).

Performs morphological opening (erosion followed by dilation) which removes structures smaller than the structuring element (stars). Fast and simple.

This is a truly non-grid method: operates on all pixels with sliding window.

Parameters:

imagendarray

Input image

maskndarray, optional

Boolean mask (True = masked pixels to exclude)

sizeint

Size of structuring element (should be larger than typical star size). Typical values: 15-35 pixels depending on seeing.

iterationsint

Number of opening iterations (default: 1, more iterations = more aggressive)

smoothbool

If True, apply additional smoothing to reduce kernel-sized artifacts. This significantly improves the smoothness of the result.

smooth_sizeint, optional

Size of smoothing kernel in pixels. If None, automatically set based on the structuring element size. Use larger values for smoother results.

Returns:

backgroundndarray

Estimated background map

Notes

• Very fast (much faster than percentile filtering)

• Simple and robust

• Can create kernel-sized artifacts without smoothing

• May suppress legitimate background features similar in size to stars

• Best for images with clear point sources on smooth backgrounds

• Use smooth=True (default) to eliminate spotty kernel artifacts

stdpipe.photometry_background.estimate_background_rms_percentile(image, mask=None, size=51)

Estimate background RMS using percentile range.

Uses interquartile range (IQR) as robust RMS estimator.

Parameters:

imagendarray

Input image

maskndarray, optional

Boolean mask (True = masked pixels)

sizeint

Kernel size for local percentile estimation

Returns:

rmsndarray

Background RMS map

stdpipe.photometry_background.estimate_background_rms_local(image, background, mask=None, size=51)

Estimate background RMS from residuals.

Computes local standard deviation of (image - background).

Parameters:

imagendarray

Input image

backgroundndarray

Estimated background

maskndarray, optional

Boolean mask (True = masked pixels)

sizeint

Kernel size for local RMS estimation

Returns:

rmsndarray

Background RMS map

stdpipe.photometry_background.get_background_polynomial(image, mask=None, order=1, sigma=3.0, maxiters=5, directional=False, direction=None)

Estimate background as a low-order 2D polynomial.

Fits a smooth polynomial surface to the image using iterative sigma-clipping to reject stars and other outliers. By construction, low-order polynomials cannot represent localised features such as nebulae, so only the smooth gradient component is captured.

Parameters:

image2D ndarray

Input image.

mask2D bool ndarray, optional

Pixels to exclude (True = masked).

orderint

Polynomial order: 1 = plane, 2 = quadratic, 3 = cubic, etc. Default 1.

sigmafloat

Sigma threshold for iterative outlier clipping. Default 3.0.

maxitersint

Maximum number of sigma-clipping iterations. Default 5.

directionalbool or float

Controls directional constraint for higher-order (≥ 2) terms.

• False (default) — use all 2D polynomial terms.

• True — auto-detect the gradient direction from a linear fit, then restrict order ≥ 2 terms to powers of the directional coordinate t = x cos θ + y sin θ. This prevents high-order polynomials from fitting localised features while still capturing curvature along the gradient.

• A float value — same as True but uses the given angle (degrees, 0 = x-axis, 90 = y-axis) instead of auto-detecting.

For order ≤ 1 this parameter has no effect (a plane already has only 3 coefficients).

directionfloat or None

Deprecated synonym for directional=<angle>. If both are given, directional takes precedence.

Returns:

background2D ndarray

Fitted polynomial background, same shape as image.

rmsfloat

Scalar RMS of the residuals (clipped pixels excluded).

Notes

Full 2D mode (directional=False):

Number of polynomial coefficients = (order+1)(order+2)/2: order 1 → 3, order 2 → 6, order 3 → 10, order 4 → 15.

Directional mode (directional=True or angle):

Number of coefficients = order + 2 (full linear + higher-order directional terms): order 2 → 4, order 3 → 5, order 4 → 6.

Coordinates are normalised to [-1, 1] internally for numerical stability. The fit uses random sub-sampling (100 000 points) for images larger than ~300×300 to keep memory and CPU bounded.

Examples

Remove a linear gradient:

>>> bg, rms = get_background_polynomial(image, order=1)

Fit a cubic curve along the gradient direction (5 parameters):

>>> bg, rms = get_background_polynomial(image, order=3, directional=True)

Specify the gradient direction explicitly (e.g. 45°):

>>> bg, rms = get_background_polynomial(image, order=3, directional=45.0)

stdpipe.photometry_background.get_background_gp(image, mask=None, max_points=3000, grid_step=None, clip_sigma=3.5, n_clip_iter=2, length_scale=64.0, length_scale_bounds=(8.0, 512.0), matern_nu=1.5, white_noise=1.0, white_noise_bounds=(0.001, 1000.0), normalize_y=True, n_restarts_optimizer=1, random_state=0, get_uncertainty=False, chunk_rows=256)

Background estimation via Gaussian Process regression (non-grid method).

Uses a Gaussian Process with a Matern kernel to model the background as a smooth function of spatial coordinates (x, y). This method is very flexible and can handle complex background variations including gradients, vignetting, and other large-scale structures.

Workflow:

1. Select candidate background pixels (masked + sigma clipping to remove stars)

2. Subsample them randomly or on a grid (to keep computation tractable)

3. Fit a GP to the (x, y) -> background mapping

4. Predict background map everywhere (with optional uncertainty)

Parameters:

imagendarray

Input image

maskndarray, optional

Boolean mask (True = masked pixels to exclude)

max_pointsint

Maximum number of training points for the GP. More points = better fit but much slower (O(N^3) cost). Typical values: 1000-5000.

grid_stepint, optional

If set, sample every N pixels on a grid (after masking). If None, random sampling is used.

clip_sigmafloat

Sigma threshold for iterative clipping of outliers (stars) in intensity space

n_clip_iterint

Number of sigma-clipping iterations to remove stars

length_scalefloat

Initial guess for GP length scale in pixels. Should be larger than typical star size. Typical values: 32-128 pixels.

length_scale_boundstuple

Bounds for length scale optimization (min, max)

matern_nufloat

Smoothness parameter for Matern kernel. Typical values: - 0.5: Exponential kernel (rough) - 1.5: Once differentiable (default, good balance) - 2.5: Twice differentiable (very smooth)

white_noisefloat

Initial guess for white noise level (in image units, e.g., ADU)

white_noise_boundstuple

Bounds for white noise optimization (min, max)

normalize_ybool

If True, normalize target values during GP training (recommended)

n_restarts_optimizerint

Number of random restarts for hyperparameter optimization (more = better but slower)

random_stateint

Random seed for reproducibility

get_uncertaintybool

If True, return uncertainty map. If False, return None for RMS.

chunk_rowsint

Process prediction in chunks of this many rows to reduce memory usage

Returns:

backgroundndarray

Estimated background map

uncertaintyndarray or None

If get_uncertainty=True: 2D array of GP posterior standard deviation If get_uncertainty=False: None

Notes

• Requires scikit-learn to be installed

• Computationally expensive for large images (use smaller max_points)

• Very flexible: can handle arbitrary background shapes

• Length scale should be >> PSF size (otherwise may fit stars)

• For very large images, consider using grid_step to reduce training points

Examples

Basic usage with default parameters:

>>> bg = get_background_gp(image, mask=mask)

Custom length scale for wide-field imaging:

>>> bg = get_background_gp(image, mask=mask, length_scale=128.0)

Get uncertainty map:

>>> bg, bg_std = get_background_gp(image, mask=mask, get_uncertainty=True)

Faster computation with grid sampling:

>>> bg = get_background_gp(image, mask=mask, max_points=1000, grid_step=10)

stdpipe.photometry_background.get_background(image, mask=None, method='sep', size=128, get_rms=False, **kwargs)

Estimate background using various methods.

This is the main interface for background estimation, supporting multiple algorithms optimized for different scenarios.

Parameters:

imagendarray

Input image

maskndarray, optional

Boolean mask (True = masked pixels)

methodstr

Background estimation method:

• ‘sep’: SEP grid-based (default) - Fast, accurate for flat/linear backgrounds

• ‘photutils’: photutils grid-based - Similar to SEP

• ‘polynomial’: Low-order polynomial fit - Captures smooth gradients, preserves nebulae

• ‘percentile’: Percentile filtering (non-grid, robust to stars) - Slow but very robust

• ‘morphology’: Morphological opening (non-grid, fast) - Good for moderate gradients

• ‘gp’: Gaussian Process regression (non-grid, very flexible) - Best for complex backgrounds

sizeint

Grid size for sep/photutils, or kernel size for percentile/morphology. For ‘gp’ method, this parameter is used as the initial length_scale. Typical values:

• SEP/photutils: 64-128 pixels

• Percentile: 31-51 pixels (smaller for gradients)

• Morphology: 15-35 pixels (should be > star size)

• GP: 32-128 pixels (length scale, should be >> PSF size)

get_rmsbool

If True, return (background, background_rms). For GP method with get_rms=True, returns the full uncertainty map instead of a scalar RMS.

**kwargsdict

Additional parameters passed to the method:

• For polynomial: order, sigma, maxiters, directional, direction (see get_background_polynomial for full list)

• For percentile: percentile, maxiters, sigma

• For morphology: iterations, smooth, smooth_size

• For gp: max_points, grid_step, clip_sigma, n_clip_iter, matern_nu, etc. (see get_background_gp for full list)

Returns:

backgroundndarray

Estimated background

background_rmsndarray or float, optional

Background RMS (if get_rms=True). For GP method, this is a 2D uncertainty map.

See also

get_background_polynomial

Polynomial fit details

get_background_percentile

Percentile filtering details

get_background_morphology

Morphological opening details

get_background_gp

Gaussian Process regression details

Notes

Method Selection Guide:

• Flat backgrounds: Use ‘sep’ (default) - fastest and most accurate

• Linear gradients: Use ‘sep’ with size=64-128

• Smooth gradients with nebulae: Use ‘polynomial’ — only captures the gradient, preserves localised background features by construction

• Directional gradients: Use ‘polynomial’ with directional=True — higher-order terms confined to the gradient direction

• Quadratic gradients: Use ‘morphology’ or local gradient fitting (bkg_order)

• Complex backgrounds (vignetting, strong curvature): Use ‘gp’

• Very crowded fields: Use ‘percentile’ with percentile=10-15

• Large images: Avoid ‘percentile’ (very slow), use ‘sep’ or ‘morphology’

Performance:

• SEP: ~1-3 ms (512x512 image)

• Morphology: ~300 ms (512x512 image)

• Percentile: ~2-5 seconds (512x512 image)

• GP: ~1-10 seconds (512x512 image, depends on max_points)

For complex backgrounds (strong curvature, vignetting), consider using local gradient fitting instead (bkg_order parameter in measure_objects), or the GP method for maximum flexibility.

Examples

Default grid-based background estimation:

>>> bg = get_background(image, method='sep', size=128)

Non-grid morphological method for moderate gradients:

>>> bg = get_background(image, method='morphology', size=25)

Morphology with extra smoothing to eliminate kernel artifacts:

>>> bg = get_background(image, method='morphology', size=25, smooth_size=15)

Percentile method for very crowded fields:

>>> bg = get_background(image, method='percentile', size=31, percentile=15.0)

Gaussian Process for complex backgrounds with vignetting:

>>> bg = get_background(image, method='gp', size=64, max_points=2000)

Get both background and RMS:

>>> bg, bg_rms = get_background(image, method='sep', get_rms=True)

Get GP background with uncertainty map:

>>> bg, bg_std = get_background(image, method='gp', get_rms=True)

Polynomial fit for gradient removal (preserves nebulae):

>>> bg = get_background(image, method='polynomial', order=2)

Directional gradient with cubic curvature:

>>> bg = get_background(image, method='polynomial', order=3, directional=True)