ChebyshevField#

final class lsst.images.fields.ChebyshevField(bounds: Bounds, coefficients: ndarray, *, unit: UnitBase | None = None)#

Bases: BaseField

A 2-d Chebyshev polynomial over a rectangular region.

Parameters#

bounds

The region where this field can be evaluated. The bbox of this region is grown by half a pixel on all sides and then used to remap coordinates to [-1, 1]x[-1, 1], which is the natural domain of a 2-d Chebyshev polynomial.

coefficients

Coefficients for the 2-d Chebyshev polynomial of the first kind, as a 2-d matrix in which element [p, q] corresponds to the coefficient of T_p(y) T_q(x). Will be set to read-only in place.

unit

Units of the field.

Attributes Summary

bounds

The region over which this field can be evaluated (Bounds).

coefficients

Coefficients for the 2-d Chebyshev polynomial of the first kind, as a 2-d matrix in which element [p, q] corresponds to the coefficient of T_p(y) T_q(x).

order

Maximum polynomial order in either dimension (int).

unit

The units of the field (astropy.units.UnitBase or None).

x_order

Maximum polynomial order in the column dimension (int).

y_order

Maximum polynomial order in the row dimension (int).

Methods Summary

__call__(*, x, y[, quantity])

Call self as a function.

deserialize(model, archive)

Deserialize the Chebyshev field from an input archive.

evaluate(*, x, y, quantity)

Evaluate at non-gridded points.

fit(bounds, data[, order, weight, y_order, ...])

Fit a Chebyshev field to data points using linear least squares.

from_legacy(legacy[, unit])

Convert from a legacy lsst.afw.math.ChebyshevBoundedField.

from_legacy_background(legacy_background[, unit])

Convert from a legacy lsst.afw.math.BackgroundMI instance.

multiply_constant(factor)

Multiply by a constant, returning a new field of the same type.

render([bbox, dtype])

Create an image realization of the field.

serialize(archive)

Serialize the Chebyshev field to an output archive.

to_legacy()

Convert to a legacy lsst.afw.math.ChebyshevBoundedField.

Attributes Documentation

bounds#
coefficients#

Coefficients for the 2-d Chebyshev polynomial of the first kind, as a 2-d matrix in which element [p, q] corresponds to the coefficient of T_p(y) T_q(x).

order#

Maximum polynomial order in either dimension (int).

unit#
x_order#

Maximum polynomial order in the column dimension (int).

y_order#

Maximum polynomial order in the row dimension (int).

Methods Documentation

__call__(*, x: ndarray, y: ndarray, quantity: bool = False) ndarray | Quantity#

Call self as a function.

static deserialize(model: ChebyshevFieldSerializationModel, archive: InputArchive[Any]) ChebyshevField#

Deserialize the Chebyshev field from an input archive.

evaluate(*, x: ndarray, y: ndarray, quantity: bool) ndarray | Quantity#

Evaluate at non-gridded points.

Parameters#

x

X coordinates to evaluate at.

y

Y coordinates to evaluate at; must be broadcast-compatible with x.

quantity

If True, return an astropy.units.Quantity instead of a numpy.ndarray. If unit is None, the returned object will be a dimensionless Quantity.

static fit(bounds: Bounds, data: ndarray | Quantity, order: int | None = None, *, y: ndarray, x: ndarray, weight: ndarray | None = None, y_order: int | None = None, x_order: int | None = None, triangular: bool = True, unit: UnitBase | None = None) ChebyshevField#

Fit a Chebyshev field to data points using linear least squares.

Parameters#

data

Data points to fit. If this is an astropy.units.Quantity, this sets the units of the field and the unit argument cannot also be provided.

order

Maximum order for the Chebyshev polynomial in both dimensions.

y

Y coordinates of the data points. Must have either the same shape as data (providing the coordinates for all points directly), or be a 1-d array with the same size as data.shape[0] (when data is a 2-d image and y provides the coordinates of the rows).

x

X coordinates of the data points. Must have either the same shape as data (providing the coordinates for all points directly), or be a 1-d array with the same size as data.shape[1] (when data is a 2-d image and x provides the coordinates of the columns).

weight

Weights to apply to the data points. Must have the same shape as data.

y_order

Maximum order for the Chebyshev polynomial in y. Requires x_order to also be provided. Incompatible with order.

x_order

Maximum order for the Chebyshev polynomial in x. Requires y_order to also be provided. Incompatible with order.

triangular

If True, only fit for coefficients of T_p(y) T_q(x) where p + q <= max(y_order, x_order).

unit

Units of the returned field.

static from_legacy(legacy: LegacyChebyshevBoundedField, unit: astropy.units.UnitBase | None = None) ChebyshevField#

Convert from a legacy lsst.afw.math.ChebyshevBoundedField.

static from_legacy_background(legacy_background: LegacyBackground, unit: astropy.units.UnitBase | None = None) ChebyshevField#

Convert from a legacy lsst.afw.math.BackgroundMI instance.

multiply_constant(factor: float | Quantity | UnitBase) ChebyshevField#

Multiply by a constant, returning a new field of the same type.

Parameters#

factor

Factor to multiply by. When this has units, those should multiply self.unit or set the units of the returned field if self.unit is None.

render(bbox: Box | None = None, *, dtype: type[Any] | dtype[Any] | _SupportsDType[dtype[Any]] | tuple[Any, Any] | list[Any] | _DTypeDict | str | None = None) Image#

Create an image realization of the field.

Parameters#

bbox

Bounding box of the image. If not provided, self.bounds.bbox will be used.

dtype

Pixel data type for the returned image.

serialize(archive: OutputArchive[Any]) ChebyshevFieldSerializationModel#

Serialize the Chebyshev field to an output archive.

to_legacy() LegacyChebyshevBoundedField#

Convert to a legacy lsst.afw.math.ChebyshevBoundedField.