API Reference

tform = mismatch2affine(mms, thresh, knots)

Return an AffineMap that is a least-squares "best initial guess" for the transformation minimizing the mismatch.

mms is an array of MismatchArrays (one per aperture, in the format returned by RegisterMismatch). thresh is the denominator threshold that determines which aperture regions have sufficient pixel/intensity overlap to be valid. knots specifies the aperture centers (see RegisterDeformation).

The algorithm fits each aperture to a quadratic, then solves a global least-squares problem over all apertures — guaranteeing a global solution at the cost of the quadratic approximation. If thresh is too restrictive, it is halved up to three times before raising an error.

Returns

• tform::AffineMap — affine transformation (linear map + translation)

To refine tform beyond the quadratic approximation, see optimize in the registration pipeline.

tfms = pat_rotation(fixed, moving)
tfms = pat_rotation(fixed, moving, SD)
tfms = pat_rotation(fixedpa, moving)
tfms = pat_rotation(fixedpa, moving, SD)

Compute the Principal Axes Transform (PAT) aligning the low-order intensity moments of two images. fixed is the reference image and moving is the image to align. fixedpa is a (center, cov) tuple from principalaxes, useful when aligning many images to the same reference to avoid recomputing its principal axes.

SD is an optional spatial-dimensions matrix that accounts for non-isotropic sampling (e.g., SD = Diagonal(voxelspacing)). Defaults to the identity.

Because intensity ellipsoids are ambiguous up to 180° rotations (sign-flips of an even number of coordinate axes), the function returns multiple candidate transforms. Evaluate alignment quality for each candidate and select the best.

Returns

• tfms::Vector{AffineMap} — 2 candidates in 2D, 4 candidates in 3D

Examples

julia> fixed = zeros(5, 7); fixed[3, 2:6] .= 1.0;   # horizontal bar

julia> moving = zeros(7, 5); moving[2:6, 3] .= 1.0;  # vertical bar

julia> tfms = pat_rotation(fixed, moving);

julia> length(tfms)
2

julia> tfms[1].linear   # ≈ 90° rotation
2×2 Matrix{Float64}:
  0.0  1.0
 -1.0  0.0

u = optimize_per_aperture(mms, thresh)

Compute the naive per-aperture displacement that minimizes the mismatch, treating each aperture independently. mms is a Vector of MismatchArrays (one per aperture) and thresh is the denominator threshold.

For each aperture, the first shift-dimension component of the argmin is recorded. The returned array has size (1, length(mms)).

See also RegisterCore.argmin_mismatch.

Returns

• u::Matrix{Float64} of size (1, n) — first shift component at the minimum for each of the n apertures

Examples

julia> using RegisterCore

julia> num1 = [(i - 1)^2 + j^2 for i in -5:5, j in -5:5];

julia> num2 = [i^2 + j^2 for i in -5:5, j in -5:5];

julia> mms = [MismatchArray(num1, ones(11, 11)), MismatchArray(num2, ones(11, 11))];

julia> optimize_per_aperture(mms, 0.5)
1×2 Matrix{Float64}:
 1.0  0.0

E0, u0, Q = qfit(mm, thresh; maxsep=size(mm), opt=true)

Perform a quadratic fit of the mismatch data in mm. Returns the mismatch value E0 and shift u0 at the minimum, plus the curvature matrix Q of the best-fit model:

    mm ≈ E0 + (u - u0)' * Q * (u - u0)

Only shift-locations where mm[i].denom > thresh are used. If no valid locations exist, returns (zero(T), zeros(T, d), zeros(T, d, d)).

maxsep restricts the fit to shifts satisfying |u[d] - u0[d]| ≤ maxsep[d]. Setting opt=false uses a fast heuristic for Q instead of a full nonlinear solve, trading accuracy for speed.

Returns

• E0::T — mismatch value at the fitted minimum
• u0::Vector{T} — shift coordinates of the fitted minimum (length ndims(mm))
• Q::Matrix{T} — symmetric positive-semidefinite curvature matrix of size (d, d)

Examples

julia> using RegisterCore

julia> num = [(i - 1)^2 + (j + 2)^2 for i in -5:5, j in -5:5];

julia> mm = MismatchArray(num, ones(11, 11));

julia> E0, u0, Q = qfit(mm, 0.5);

julia> E0
0.0

julia> u0
2-element Vector{Float64}:
  1.0
 -2.0

julia> Q ≈ [1.0 0.0; 0.0 1.0]
true

cs, Qs, mmis = mms2fit!(mms, thresh)

Compute per-aperture shifts and quadratic-fit matrices for the N-dimensional array-of-MismatchArrays mms, using thresh as the denominator threshold. Also prepares mms for interpolation, modifying it in-place after extracting cs and Qs.

The dimension N of the container mms must equal the number of spatial dimensions of each MismatchArray element. For example, a 2×3 matrix of 2D mismatch arrays is valid; a Vector of 2D mismatch arrays is not.

Returns

• cs: Array{SVector{N,T},N} — per-aperture shift positions
• Qs: Array{SMatrix{N,N,T},N} — per-aperture quadratic curvature matrices
• mmis: array of interpolated MismatchArrays, suitable for input to RegisterPenalty.fixed_λ and RegisterPenalty.auto_λ

Examples

julia> using RegisterCore

julia> num1 = [(i - 1)^2 + (j + 2)^2 for i in -5:5, j in -5:5];

julia> num2 = [(i + 1)^2 + (j - 1)^2 for i in -5:5, j in -5:5];

julia> denom = ones(11, 11);

julia> mms = reshape([MismatchArray(num1, denom), MismatchArray(num2, denom)], 1, 2);

julia> cs, Qs, mmis = mms2fit!(mms, 0.5);

julia> cs[1, 1] ≈ [1.0, -2.0]
true

r = qbuild(E0, umin, Q, maxshift)

Build a CenterIndexedArray representing the quadratic mismatch approximation over the full shift domain [-maxshift[d], maxshift[d]]. The quadratic model is:

    r[u] = E0 + (u - umin)' * Q * (u - umin)

E0, umin, and Q are the outputs of qfit. Useful for debugging and visualizing the quadratic fit.

Returns

• r::CenterIndexedArray — evaluated mismatch, indexed from -maxshift to +maxshift

Examples

julia> using RegisterCore

julia> num = [(i - 1)^2 + (j + 2)^2 for i in -5:5, j in -5:5];

julia> mm = MismatchArray(num, ones(11, 11));

julia> E0, umin, Q = qfit(mm, 0.5);

julia> r = qbuild(E0, umin, Q, (5, 5));

julia> r[1, -2]
0.0

julia> r[0, 0]
5.0

tf = uisvalid(u, maxshift)

Return true if every entry of the displacement array u is within the allowed domain: |u[idim, j]| < maxshift[idim] - 0.5001 for all aperture positions j and displacement dimensions idim. Returns false as soon as any entry violates this condition.

Examples

julia> uisvalid([1.5, 0.5], (3, 3))
true

julia> uisvalid([2.5, 0.5], (3, 3))
false

uclamp!(u, maxshift)

Clamp the entries of the displacement array u in-place so that each satisfies |u[idim, j]| ≤ maxshift[idim] - 0.51. Returns u.

Accepts both numeric arrays of shape (nd, apertures...) and arrays whose elements are StaticVectors (e.g., Array{SVector{N,T}}); both representations are mutated in-place.

Examples

julia> u = [4.0, -5.0];

julia> uclamp!(u, (3, 3))
2-element Vector{Float64}:
  2.49
 -2.49

center, cov = principalaxes(img)

Compute the intensity-weighted centroid and covariance of image img. Coordinates are 1-based array indices. NaN pixels are ignored.

Returns

• center::Vector{T} — intensity-weighted centroid, length ndims(img)
• cov::Matrix{T} — N×N intensity-weighted covariance matrix

Examples

julia> img = zeros(5, 5); img[3, 3] = 1.0;

julia> center, cov = principalaxes(img);

julia> center
2-element Vector{Float64}:
 3.0
 3.0

julia> cov
2×2 Matrix{Float64}:
 0.0  0.0
 0.0  0.0