Trace Detection Algorithm

Overview

Spectral trace detection identifies the spatial position of stellar/galactic spectra in 2D longslit images. pyKOSMOS++ uses cross-correlation with Gaussian templates to robustly detect traces even in low signal-to-noise conditions.

Method: Cross-Correlation with Gaussian Templates

Algorithm Steps

  1. Collapse 2D spectrum spatially (median along spectral direction) to create 1D spatial profile

  2. Create Gaussian template with expected FWHM:

    template[y] = exp(-0.5 * ((y - center) / sigma)^2)

    where sigma = FWHM / 2.355

  3. Cross-correlate spatial profile with template using scipy.ndimage.correlate1d

  4. Detect peaks in correlation signal using scipy.signal.find_peaks with:

    • Height threshold: min_snr * background_noise

    • Minimum separation: min_separation pixels

  5. Trace spatial position along spectral direction:

    • For each detected peak position

    • Fit Gaussian centroid in overlapping spectral windows

    • Track center position vs wavelength

Implementation

from scipy.ndimage import correlate1d
from scipy.signal import find_peaks

def detect_traces_cross_correlation(
    data: np.ndarray,
    expected_fwhm: float = 4.0,
    min_snr: float = 3.0,
    min_separation: int = 20
) -> List[Trace]:
    """
    Detect spectral traces via cross-correlation.

    Parameters
    ----------
    data : ndarray
        2D calibrated spectrum (spatial x spectral)
    expected_fwhm : float
        Expected trace FWHM in spatial direction (pixels)
    min_snr : float
        Minimum signal-to-noise ratio for detection
    min_separation : int
        Minimum pixel separation between distinct traces

    Returns
    -------
    traces : List[Trace]
        Detected trace objects with spatial positions
    """
    # Collapse spatial profile
    spatial_profile = np.median(data, axis=1)

    # Create Gaussian template
    sigma = expected_fwhm / 2.355
    y = np.arange(len(spatial_profile))
    template = np.exp(-0.5 * ((y - len(y)//2) / sigma)**2)

    # Cross-correlate
    correlation = correlate1d(spatial_profile, template, mode='constant')

    # Detect peaks
    noise = np.std(correlation)
    height_threshold = min_snr * noise
    peaks, properties = find_peaks(
        correlation,
        height=height_threshold,
        distance=min_separation
    )

    # Create Trace objects for each detection
    traces = []
    for peak in peaks:
        trace = trace_spatial_position(data, peak, expected_fwhm)
        traces.append(trace)

    return traces

Advantages

Robust to Noise:

Cross-correlation averages over multiple pixels, improving SNR compared to simple peak-finding.

Template Matching:

Using expected seeing profile improves detection accuracy for extended sources.

Multi-Trace Support:

Automatically detects multiple traces in multi-slit observations.

Parameters

expected_fwhm

Type: float Default: 4.0 pixels Range: 2.0 - 20.0 pixels

Full-width half-maximum of expected trace profile. Should match observing conditions:

  • Good seeing (~0.8”): FWHM ≈ 3-4 pixels

  • Average seeing (~1.5”): FWHM ≈ 5-6 pixels

  • Poor seeing (~2.5”): FWHM ≈ 8-10 pixels

  • Extended galaxy: FWHM ≈ 10-20 pixels

How to determine:

# Measure from spatial profile
import numpy as np
spatial_profile = np.median(data, axis=1)

# Find FWHM by fitting Gaussian
from scipy.optimize import curve_fit
def gaussian(x, amp, center, sigma):
    return amp * np.exp(-0.5 * ((x - center) / sigma)**2)

popt, _ = curve_fit(gaussian, np.arange(len(spatial_profile)), spatial_profile)
fwhm = 2.355 * popt[2]
print(f"Measured FWHM: {fwhm:.2f} pixels")

min_snr

Type: float Default: 3.0 Range: 1.0 - 10.0

Signal-to-noise threshold for detection. Lower values detect fainter traces but increase false positives.

Guidelines:

  • SNR ≥ 5: Very reliable detections

  • SNR ≥ 3: Standard threshold, good balance

  • SNR ≥ 2: For faint extended sources (more false positives)

  • SNR ≥ 10: Very stringent, only brightest traces

min_separation

Type: int Default: 20 pixels Range: 5 - 100 pixels

Minimum spatial separation between distinct traces. Prevents duplicate detections of the same trace.

Guidelines:

  • Multi-slit: Set to ~0.8 × slit separation

  • Single object: Default (20 pixels) sufficient

  • Crowded field: Reduce to 10-15 pixels

Trace Position Fitting

After initial detection, trace positions are refined by fitting the spatial profile in spectral windows:

def trace_spatial_position(
    data: np.ndarray,
    initial_position: int,
    window_width: int = 50
) -> Trace:
    """Refine trace position along spectral direction."""
    n_spatial, n_spectral = data.shape
    spatial_positions = []
    spectral_pixels = []

    # Slide window along spectral direction
    for i in range(0, n_spectral, window_width // 2):
        window = data[:, i:i+window_width]
        profile = np.median(window, axis=1)

        # Fit Gaussian centroid
        center = fit_gaussian_center(profile, initial_position)

        spatial_positions.append(center)
        spectral_pixels.append(i + window_width // 2)

    return Trace(
        trace_id=...,
        spatial_positions=np.array(spatial_positions),
        spectral_pixels=np.array(spectral_pixels)
    )

Performance

Typical Performance:

  • Detection time: 0.5-2 seconds for 2048×515 image

  • Memory usage: ~50 MB for full-resolution KOSMOS frame

  • False positive rate: <5% with default parameters

Failure Modes:

  1. No traces detected:

    • Spectrum too faint (SNR < 2)

    • Wrong expected_fwhm (off by factor of 2+)

    • Incorrect calibration (flat/bias issues)

  2. Multiple detections of same trace:

    • min_separation too small

    • Trace has brightness variations

  3. Missed faint traces:

    • min_snr too high

    • Need spatial binning to boost SNR

References

  • pyKOSMOS trace detection (Davenport et al. 2023)

  • scipy.ndimage.correlate1d documentation

  • scipy.signal.find_peaks for peak detection

See Also