Stellar Alignment

Register one image onto another's pixel grid — using stars, plate solutions, or both — the engine behind Stacking Suite, Dither Analysis, Exoplanet Hunter, and every other SASpro tool that needs frames in the same coordinate system.

1. Picking a source and target

The two top boxes each hold a selector:

• Source Image (Reference) — the image whose pixel grid everything else gets aligned onto. The output image inherits its dimensions and WCS.
• Target Image (To be Aligned) — the image that gets warped. The output is a warped copy of this image, same shape as the source.

Each box offers From File (open a FITS or image file) or From View (pick any open SASpro document from the dropdown). Click Refresh next to a view dropdown if you've opened new images since the dialog launched. Source and target can mix freely — file-against-view, view-against-view, file-against-file.

2. The three alignment methods

The Method dropdown at the top of the Transform / Distortion box picks the alignment strategy:

3. Check WCS and Trial Detect

Two small buttons let you verify the inputs before committing to a long alignment job:

Both buttons populate inline status labels — they never pop a modal — so you can tweak settings and re-check freely.

4. Transform models

Below the method, the Model dropdown picks the mathematical form of the final transform:

Polynomial models always build on top of a base affine or homography — the polynomial fits only the residual correction after the linear model is applied. This two-pass structure gives much better stability than trying to fit a pure polynomial from scratch.

5. Detection and matching parameters

6. Running and inspecting

Run Alignment kicks the job off. Status updates stream to the left panel. When it finishes, the aligned result appears in the right-panel preview. Use AutoStretch: ON/OFF to toggle a simple percentile stretch on the preview — doesn't modify the data, just helps you see the low-contrast detail.

For Affine or Homography alignments, a Transformation Matrix Details popup shows the full matrix and derived parameters: translation (t_x, t_y), scale (x, y), rotation angle, and shear. When the alignment ran as WCS→Stellar refinement, the dialog notes that the matrix shown is only the residual star-refinement delta on top of the (much larger) WCS reprojection.

7. Delivering the result

8. Workflow recipes

9. Headless / shortcut use

The tool integrates with SASpro's drag-grip shortcut system via a command ID of star_alignment. A preset payload has this shape:

Reference types supported:

• {"type": "active"} — use whatever document is active when the drop happens.
• {"type": "view_name", "name": "..."} — match by visible document title.
• {"type": "view_ptr", "doc_ptr": 1234} — exact document pointer (for round-trippable shortcuts).
• {"type": "file", "path": "/abs/path.fits"} — load the reference from disk.

Dropping the preset onto a target view runs the alignment headlessly — detect stars, solve the transform, warp, deliver — with no dialog. When overwrite is true the target doc is replaced in place; otherwise a new view is created. The default behavior across SASpro tools that depend on alignment (Stacking Suite, Dither, Exoplanet Hunter) is to call this engine directly rather than through the dialog, so the same math and the same tuning parameters apply everywhere.

10. Technical details

The star detection, triangle matching, RANSAC geometry fitting, downsample lift math, multi-pass refinement, and the anti-ring warp.

10.1 Star detection via SEP

Both the standalone dialog and the batch registration engine use SEP (Barbary 2016, a Python port of Source Extractor by Bertin & Arnouts 1996) for detection. The pipeline:

1. Compute a 64×64-cell background mesh and global RMS: sep.Background(img, bw=64, bh=64).
2. Threshold at σ · RMS_bg:
   \text{threshold} = \sigma \cdot \mathrm{RMS}_{\mathrm{bg}}
3. Extract sources via connected-component labeling with minimum area minarea.
4. For each detected object, SEP reports flux, centroid (x, y), semi-major axis a, semi-minor axis b, and pixel count.

SEP is explicitly not thread-safe at the C level, so every call runs inside a module-level lock (_AA_LOCK). Pixel-stack size is grown lazily if the "internal pixel buffer full" error appears — on dense fields the default 1.5M pixels isn't enough.

10.2 Hot pixel rejection

Raw SEP detections include hot pixels, cosmic ray hits, and other spurious sources that would corrupt the geometric fit. Three filters are applied:

FWHM proxy — from the semi-major axis a, estimate the full-width at half-maximum assuming a Gaussian PSF:

Reject detections where FWHM < 1.2 px. Real stars, even undersampled, have a ≥ 1.0 px; hot pixels tend to be a ~ 0.5 px.

Ellipticity — from the axis ratio:

Reject detections where e > 0.6. Real stars are nearly round (e < 0.3 typically). Cosmic ray streaks and satellite trails have e > 0.8.

Pixel count — require npix ≥ minarea (default 10). Belt and suspenders with the morphological filters above; sometimes SEP's connected-component bookkeeping includes very small objects that slip past the FWHM filter.

If all detections are rejected (a legitimate scenario on heavily under-sampled frames), the pipeline falls back to a size-only filter of npix ≥ minarea/2 — the alignment should never fail silently just because the hot-pixel filter was aggressive.

10.3 Uniform-grid control point selection

After hot-pixel rejection, the surviving detections are sorted by flux (brightest first) and allocated into a 4×4 or 6×6 grid of cells over the image. Each cell can hold up to max_per_cell stars (default 25); the total is capped at max_total (default 500). This gives a spatially uniform sample of bright stars rather than a flux-sorted cluster all in the brightest region of the frame.

Two quality gates on the resulting set:

Coverage — fraction of grid cells with at least one star:

Below 0.5 triggers a log warning — the alignment may still succeed but the solution is probably biased toward one side of the frame.

Spread — require the 5th-95th percentile span of x (and y) to cover at least 35% of the image width (height):

and at least 6 of 9 cells in a 3×3 grid occupied. Points clustered in a small region produce ill-conditioned transforms that fit the local area perfectly but fail wildly at the edges.

10.4 Astroalign triangle matching

With control-point sets in hand, astroalign (Beroiz et al. 2020) builds Delaunay triangulations of each, then matches triangles by their shape-invariant features. For a triangle with sorted side lengths d₁ ≤ d₂ ≤ d₃:

The invariant is a 2D feature that's unchanged by translation, rotation, and uniform scale. Two triangles with the same invariant (within a tolerance) are candidate matches; multiple candidate matches generate a RANSAC vote for the global transform.

If initial detection fails (not enough triangles agree), the backoff retries with progressively stricter detection σ: 15, 25, 50, 80, 120. Each stricter pass yields fewer but higher-quality stars, which can succeed where a noise-contaminated dense set didn't.

10.5 Center-crop solve with coordinate lift

The initial astroalign solve runs not on the full reference image but on a center crop at ~1.2× the target image's dimensions. This has two benefits: the crop is faster to detect in, and it avoids the "target footprint sticks outside the reference" edge case that produces zero valid matches.

After matches are found in crop coordinates, the target-side x and y are shifted by the crop's top-left offset (x₀, y₀) to lift them into full-reference coordinates. The transform is then re-estimated with these lifted correspondences for maximum accuracy. If the first solve returns too few spatially-spread matches, the pipeline retries with a 3×3 grid of tiles across the reference — this gives much better coverage on frames where distortion makes a single-crop solve under-sample one quadrant.

10.6 Downsample affine lift

For performance, the solve runs on a factor-d downsampled grid. The transform returned in the downsampled coordinate system must be lifted back to full resolution. If S is the downsample scaling matrix:

then an affine A measured in downsampled space maps to full-resolution space via conjugation:

and the inverse (used when applying an already-known full-resolution transform to a downsampled frame) is:

The same conjugation applies to homographies (3×3). Doing the solve at d=3 is roughly 9× faster than at full resolution and loses almost no accuracy for ordinary astronomical sampling — the lifted transform is still sub-pixel accurate.

10.7 Similarity projection

A similarity transform is an affine with no shear and uniform scale (rotation + uniform scale + translation only — 4 degrees of freedom). It's useful for cases where you know the geometry is physically a rigid rotation plus zoom and want to prevent spurious shear from fitting detection noise. The nearest similarity to a given 2×3 affine is found via SVD of its 2×2 linear part M:

where U, V come from M = U Σ V^T, R is the closest rotation (determinant corrected if needed), and s is the uniform scale (mean of the two singular values). The resulting similarity is sR with the original translation preserved. This is the operation that cleans up shear from a star-matched affine when the user has picked the Similarity model.

10.8 Anti-ring warp (Lanczos + bilinear blend)

Lanczos interpolation produces sharp, high-quality warps everywhere except around saturated stars, where its negative lobes create a visible dark ring (the "Lanczos undershoot" halo). Clipping Lanczos to a bilinear-interpolation envelope creates rosette-pattern facets around the stars from quantizing the kernel's lobes. The right solution is a blend:

1. Compute the full-frame bilinear warp (soft, no ring).
2. Compute the full-frame Lanczos warp (sharp, ring on bright stars).
3. Build a bright-pixel mask from the bilinear warp using robust statistics:
   \tau = \mathrm{median}(L) + k \cdot 1.4826 \cdot \mathrm{MAD}(L)
   with k = 3.0 by default (1.4826 is the Gaussian MAD → sigma conversion).
4. Dilate the mask by ~6 px (the Lanczos radius plus safety margin) so it covers the ring zone, not just the star core.
5. Feather the mask edges with a 2-px-sigma Gaussian so the seam is invisible.
6. Blend:
   I_{\mathrm{out}} = I_{\mathrm{Lanczos}} \cdot (1 - m) + I_{\mathrm{bilinear}} \cdot m

Where there's no bright star (m ≈ 0), the output is pure Lanczos — all the sharpness benefits. On and around bright stars (m ≈ 1), it's pure bilinear — no ring. The transition is smooth, imperceptible on inspection, and preserves the Lanczos quality on 99%+ of the frame.

10.9 Multi-pass refinement (batch registration)

The batch Star Registration engine (used by Stacking Suite) runs multiple refinement passes. Pass 1 solves on the original frame; subsequent passes solve on the already-partially-aligned frame, producing a small residual delta each time. The accumulated transform is updated by matrix composition:

Convergence is declared when the maximum delta shift across all frames falls below a tolerance (default 0.2 px):

Frames that have already converged are skipped in later passes, concentrating work on the stragglers. Default is 3 passes maximum; in practice most frames converge after pass 2.

10.10 WCS reprojection pipeline

For the WCS and WCS→Stellar methods, the target image is reprojected onto the source grid via a two-step projection: for every output pixel (x, y) in the source grid, the source WCS forward-projects to a sky position (α, δ), then the target WCS inverse-projects back to a target pixel location (x', y'), where the image is sampled. This uses astropy's all_pix2world and all_world2pix, which honor the full SIP distortion polynomial on both ends.

Tiled processing: the output grid is swept in 512×512 tiles rather than one monolithic allocation, giving progress feedback and bounded memory. A round-trip check — forward through the target WCS, back through the source WCS, compare to the original coordinate — rejects pixels where the SIP inverse iteration diverged, replacing them with the no-data fill value. Rejection count is logged.

10.11 Polynomial residual warp

When Polynomial (order 3 or 4) is selected, two warps are chained. The first applies the base affine or homography to get the target image onto an approximately-correct grid. The second is a polynomial correction applied via skimage.transform.warp with an inverse_map callable: for each output pixel, compute its normalized location, evaluate the polynomial to get its source location in the base-warped image, and sample. The polynomial is estimated from inlier matches (those within RANSAC threshold of the base model's prediction), normalized to the [0, 1] domain for numerical conditioning. Points per channel are matched between the color stripes to keep RGB alignment identical.

10.12 Parallelism

Batch registration dispatches per-frame solves across a process pool (ProcessPoolExecutor) on native (non-frozen) builds and spawn-aware processes on PyInstaller-frozen builds (guarded by freeze_support() in the entry point to catch worker re-execs). Each worker caps BLAS and OpenCV thread counts to 1 so N frames × M internal cores doesn't oversubscribe the CPU. The finalize pass (full-resolution warp + save) uses a smaller worker pool but allows each worker's OpenCV to use multiple threads — the trade-off between number of concurrent workers and threads-per-worker is scaled by core count.

When process pool initialization fails (common on mixed-env Windows installations where child workers die on NumPy/SciPy import), the engine falls back to a thread pool so registration still completes. The thread fallback serializes astroalign calls through the global _AA_LOCK rather than parallelizing them.

10.13 Drizzle self-description stamp

Each aligned frame's FITS header carries a set of SASDZ* keywords that let the drizzle integrator reconstruct the frame's transform without needing the .sasd sidecar file. SASDZKND stores the transform kind ("affine", "homography", "poly3", "poly4"); SASDZRFH and SASDZRFW give the reference grid dimensions; SASDZORG gives the source basename; and for affine/homography, SASDZNR, SASDZNC, and SASDZ00 through SASDZ11 flatten the matrix row-major into individual FITS cards (which can only hold 8-character keyword names). Polynomial frames carry only the kind and reference dimensions, since their transforms are dense coordinate maps rather than small matrices — the drizzle integrator treats these as already-warped and deposits them on an identity canvas.