import numpy as np
import matplotlib.pyplot as plt
import astroapers as aapCustom weights and composition
This tutorial shows how to combine bbox-tight weight arrays for more complex shapes, instead of using a simple ndarray-typed mask.
A small image and two component apertures
ny, nx = 48, 54
y, x = np.mgrid[:ny, :nx]
x0, y0 = 26.0, 24.0
data = 8.0 + 0.02 * x - 0.01 * y
data += 120.0 * np.exp(-0.5 * ((x - x0) ** 2 + (y - y0) ** 2) / 3.0**2)
outer = aap.CircAp((x0, y0), r=8.0)
hole = aap.EllipAp((x0 + 2.0, y0), a=3.0, b=1.4, theta=0.5)
outer_weights = outer.weights_exact()[0]
outer_bbox = outer.bboxes()[0]
hole_weights = hole.weights_exact()[0]
hole_bbox = hole.bboxes()[0]
print(f"Outer bbox: {outer_bbox}")
print(f"Outer weights shape: {outer_weights.shape}")
print(f"Hole bbox: {hole_bbox}")
print(f"Hole weights shape: {hole_weights.shape}")Outer bbox: BoundingBox(ixmin=18, ixmax=35, iymin=16, iymax=33)
Outer weights shape: (17, 17)
Hole bbox: BoundingBox(ixmin=25, ixmax=32, iymin=22, iymax=27)
Hole weights shape: (5, 7)The two weight arrays are bbox-tight. They usually have different origins and shapes, so arithmetic only makes sense after placing them on a common bounding box.
Align bbox-tight weights before arithmetic
def union_bbox(boxes):
return aap.BoundingBox(
ixmin=min(box.ixmin for box in boxes),
ixmax=max(box.ixmax for box in boxes),
iymin=min(box.iymin for box in boxes),
iymax=max(box.iymax for box in boxes),
)
def place_on_bbox(weights, box, bbox):
out = np.zeros(bbox.shape, dtype=weights.dtype)
y0 = box.iymin - bbox.iymin
y1 = y0 + weights.shape[0]
x0 = box.ixmin - bbox.ixmin
x1 = x0 + weights.shape[1]
out[y0:y1, x0:x1] = weights
return out
bbox = union_bbox([outer_bbox, hole_bbox])
outer_on_bbox = place_on_bbox(outer_weights, outer_bbox, bbox)
hole_on_bbox = place_on_bbox(hole_weights, hole_bbox, bbox)
weights = np.clip(outer_on_bbox - hole_on_bbox, 0.0, 1.0)Now weights and bbox represent the union of the outer circle and the hole. The weights are clipped to [0, 1] to preserve aperture-mask semantics, but you could also create a non-binary mask by skipping the clipping step.
fig, axes = plt.subplots(1, 3, figsize=(9, 3))
for ax, image, title in [
(axes[0], outer_on_bbox, "outer circle"),
(axes[1], hole_on_bbox, "elliptical exclusion"),
(axes[2], weights, "composed weights"),
]:
ax.imshow(image, origin="lower", vmin=0, vmax=1)
ax.set_title(title)
plt.show()
Use the composed mask
Use BoundingBox methods with the composed weights. Bad-pixel masks are separate full-image boolean arrays where True means excluded.
bad_pixels = np.zeros(data.shape, dtype=bool)
bad_pixels[26:30, 31:36] = True
apsum, npix = bbox.apsum(weights, data, mask=bad_pixels)
mean_value = apsum / npix
apsum, npix, mean_value(6812.401187401858, 177.56867142279435, 38.364882345610404)fig, ax = plt.subplots(figsize=(5, 4))
ax.imshow(data, origin="lower", cmap="viridis")
ax.imshow(bbox.to_image(weights, data.shape) > 0, origin="lower", cmap="gray", alpha=0.25)
ax.imshow(bad_pixels, origin="lower", cmap="gray_r", alpha=0.35)
ax.set_title("Composed aperture and bad pixels")
plt.show()
For unions, use np.maximum.reduce([...]) on aligned component weights. For exclusions or annuli, subtract inner weights from outer weights and clip to [0, 1]. Use the same mask method for every component unless mixed exact/center semantics are intentional.
Raw weights from another source
Of course the user can give hand-written weights and a bounding box:
raw_weights = np.array(
[
[0.00, 0.25, 0.75, 1.00],
[0.50, 1.00, 1.00, 0.50],
[0.00, 0.50, 0.25, 0.00],
],
dtype=np.float32,
)
raw_bbox = aap.BoundingBox(ixmin=20, ixmax=24, iymin=18, iymax=21)
raw_apsum, raw_npix = raw_bbox.apsum(raw_weights, data)
raw_apsum, raw_npix(105.37785765301186, 5.75)When doing this, weights.shape must match bbox.shape, which is (iymax - iymin, ixmax - ixmin). Raw bbox-tight weights are sampled weights, not analytic geometry, so they have no automatic .area, plot path, or polygon engine.
BoundingBox methods preserve user-supplied float32 and float64 weight arrays where practical. Other numeric or boolean weight dtypes are converted to float64, including extended-precision dtypes such as float128/longdouble on platforms where NumPy exposes them. Array-producing helpers such as .to_image(), .weighted_cutout(), and .weighted_values() preserve float32 when both data and weights are float32, but BoundingBox.apsum() and BoundingBox.npix() always accumulate in float64.
Maximum Performance
For custom composed masks, the raw Rust extension can generate the component bbox-tight weights without constructing aperture objects. Import it as aapr. The composition itself is still your Python logic; there is no raw fused Rust function for an arbitrary user-composed shape. For more raw-call patterns, inspect astroapers.kernels.
import astroapers._rust as aapr
raw_data = np.ascontiguousarray(data, dtype=np.float64)
x_outer = np.ascontiguousarray([x0], dtype=np.float64)
y_outer = np.ascontiguousarray([y0], dtype=np.float64)
x_hole = np.ascontiguousarray([x0 + 2.0], dtype=np.float64)
y_hole = np.ascontiguousarray([y0], dtype=np.float64)
cw, c_xmin, c_xmax, c_ymin, c_ymax = aapr.weights_circ_exact(
x_outer,
y_outer,
8.0,
)
ew, e_xmin, e_xmax, e_ymin, e_ymax = aapr.weights_ellip_exact(
x_hole,
y_hole,
3.0,
1.4,
0.5,
)
c_bounds = (c_xmin[0], c_xmax[0], c_ymin[0], c_ymax[0])
e_bounds = (e_xmin[0], e_xmax[0], e_ymin[0], e_ymax[0])
u_bounds = (
min(c_bounds[0], e_bounds[0]),
max(c_bounds[1], e_bounds[1]),
min(c_bounds[2], e_bounds[2]),
max(c_bounds[3], e_bounds[3]),
)
def reshape_raw(flat, bounds):
ixmin, ixmax, iymin, iymax = bounds
return np.asarray(flat, dtype=np.float64).reshape(iymax - iymin, ixmax - ixmin)
def place_raw(weight, bounds, out_bounds):
ixmin, ixmax, iymin, iymax = bounds
uxmin, _, uymin, _ = out_bounds
out = np.zeros((out_bounds[3] - out_bounds[2], out_bounds[1] - out_bounds[0]))
out[iymin - uymin:iymax - uymin, ixmin - uxmin:ixmax - uxmin] = weight
return out
raw_outer = place_raw(reshape_raw(cw[0], c_bounds), c_bounds, u_bounds)
raw_hole = place_raw(reshape_raw(ew[0], e_bounds), e_bounds, u_bounds)
raw_composed = np.clip(raw_outer - raw_hole, 0.0, 1.0)
uxmin, uxmax, uymin, uymax = u_bounds
raw_cutout = raw_data[uymin:uymax, uxmin:uxmax]
raw_apsum = np.sum(raw_composed * raw_cutout, dtype=np.float64)
raw_npix = np.sum(raw_composed, dtype=np.float64)
expected_apsum, expected_npix = bbox.apsum(weights, data)
assert np.allclose(raw_apsum, expected_apsum)
assert np.allclose(raw_npix, expected_npix)
raw_apsum, raw_npix(np.float64(7005.424459504775), np.float64(187.86724068466964))