Radar Models API

Radar equation and detection probability calculations.

Overview

from phased_array_systems.models.radar import (
    RadarModel,
    compute_detection_threshold,
    compute_pd_from_snr,
    compute_snr_for_pd,
    albersheim_snr,
    coherent_integration_gain,
    noncoherent_integration_gain,
    integration_loss,
)

Classes

RadarModel

Radar range equation calculator.

Implements the monostatic radar range equation:

P_r = (P_t * G^2 * λ^2 * σ) / ((4π)^3 * R^4 * L_sys)

Or in dB form

SNR = P_t + 2G + 2λ_dB + σ_dBsm - 4*R_dB - L_sys - (4π)^3_dB - N_dB

Where

P_t = Peak transmit power (W) G = Antenna gain (same for Tx/Rx in monostatic) λ = Wavelength (m) σ = Target radar cross section (m^2) R = Range to target (m) L_sys = System losses N = Noise power = kTB

ATTRIBUTE	DESCRIPTION
name	Model block name for identification TYPE: str

evaluate

evaluate(arch: Architecture, scenario: RadarDetectionScenario, context: dict[str, Any]) -> MetricsDict

Evaluate radar detection performance.

PARAMETER	DESCRIPTION
arch	Architecture configuration TYPE: Architecture
scenario	Radar detection scenario TYPE: RadarDetectionScenario
context	Additional context (may include antenna metrics): - g_peak_db: Antenna gain (uses this if provided) - scan_loss_db: Scan loss (uses this if provided) - beamwidth_az_deg: Azimuth beamwidth (for clutter cell) - beamwidth_el_deg: Elevation beamwidth (for clutter cell) TYPE: dict[str, Any]

RETURNS

DESCRIPTION

MetricsDict

Dictionary with radar metrics: - peak_power_w: Peak transmit power (W) - peak_power_dbw: Peak transmit power (dBW) - g_ant_db: Antenna gain (dB) - wavelength_m: Wavelength (m) - target_rcs_dbsm: Target RCS (dBsm) - target_rcs_m2: Target RCS (m^2) - range_m: Target range (m) - noise_power_dbw: Noise power (dBW) - snr_single_pulse_db: Single-pulse SNR (dB) - integration_gain_db: Integration gain (dB) - snr_integrated_db: Integrated SNR (dB) - snr_required_db: Required SNR for Pd/Pfa (dB) - snr_margin_db: SNR margin (dB) - pd_achieved: Achieved probability of detection - detection_range_m: Max detection range for required Pd (m) - clutter_rcs_dbsm: Clutter RCS if applicable (dBsm) - scr_db: Signal-to-clutter ratio if applicable (dB) - scnr_db: Signal-to-clutter-plus-noise ratio (dB) - atmos_loss_db: Two-way atmospheric loss (dB) - rain_loss_db: Two-way rain attenuation (dB) - cfar_loss_db: CFAR processing loss (dB)

Source code in src/phased_array_systems/models/radar/equation.py

def evaluate(
    self,
    arch: Architecture,
    scenario: RadarDetectionScenario,
    context: dict[str, Any],
) -> MetricsDict:
    """Evaluate radar detection performance.

    Args:
        arch: Architecture configuration
        scenario: Radar detection scenario
        context: Additional context (may include antenna metrics):
            - g_peak_db: Antenna gain (uses this if provided)
            - scan_loss_db: Scan loss (uses this if provided)
            - beamwidth_az_deg: Azimuth beamwidth (for clutter cell)
            - beamwidth_el_deg: Elevation beamwidth (for clutter cell)

    Returns:
        Dictionary with radar metrics:
            - peak_power_w: Peak transmit power (W)
            - peak_power_dbw: Peak transmit power (dBW)
            - g_ant_db: Antenna gain (dB)
            - wavelength_m: Wavelength (m)
            - target_rcs_dbsm: Target RCS (dBsm)
            - target_rcs_m2: Target RCS (m^2)
            - range_m: Target range (m)
            - noise_power_dbw: Noise power (dBW)
            - snr_single_pulse_db: Single-pulse SNR (dB)
            - integration_gain_db: Integration gain (dB)
            - snr_integrated_db: Integrated SNR (dB)
            - snr_required_db: Required SNR for Pd/Pfa (dB)
            - snr_margin_db: SNR margin (dB)
            - pd_achieved: Achieved probability of detection
            - detection_range_m: Max detection range for required Pd (m)
            - clutter_rcs_dbsm: Clutter RCS if applicable (dBsm)
            - scr_db: Signal-to-clutter ratio if applicable (dB)
            - scnr_db: Signal-to-clutter-plus-noise ratio (dB)
            - atmos_loss_db: Two-way atmospheric loss (dB)
            - rain_loss_db: Two-way rain attenuation (dB)
            - cfar_loss_db: CFAR processing loss (dB)
    """
    # Get antenna gain from context or compute approximate
    if "g_peak_db" in context:
        g_ant_db = context["g_peak_db"]
        # Apply scan loss if provided
        if "scan_loss_db" in context:
            g_ant_db -= context["scan_loss_db"]
    else:
        # Approximate gain for uniform rectangular array
        # G ≈ 4*pi*A/λ^2 = 4*pi * (nx*dx) * (ny*dy) when spacing in wavelengths
        aperture_lambda_sq = (
            arch.array.nx * arch.array.dx_lambda * arch.array.ny * arch.array.dy_lambda
        )
        g_ant_linear = 4 * math.pi * aperture_lambda_sq
        g_ant_db = 10 * math.log10(g_ant_linear)

    # Get beamwidths from context or approximate
    beamwidth_az_deg = context.get("beamwidth_az_deg", 5.0)
    beamwidth_el_deg = context.get("beamwidth_el_deg", 5.0)

    # Transmit power (peak)
    n_elements = arch.array.n_elements
    peak_power_w = arch.rf.tx_power_w_per_elem * n_elements
    peak_power_dbw = W_TO_DBW(peak_power_w)

    # Wavelength
    wavelength_m = C_LIGHT / scenario.freq_hz
    wavelength_db = 10 * math.log10(wavelength_m)

    # Range resolution
    range_resolution_m = scenario.range_resolution_m

    # System losses (feed network + additional system losses)
    system_loss_db = arch.rf.feed_loss_db + arch.rf.system_loss_db

    # Target RCS
    rcs_dbsm = scenario.target_rcs_dbsm
    rcs_m2 = 10 ** (rcs_dbsm / 10)

    # Range
    range_m = scenario.range_m
    range_db = 10 * math.log10(range_m)

    # Compute grazing angle if not specified
    if scenario.grazing_angle_deg is not None:
        grazing_angle = scenario.grazing_angle_deg
    else:
        grazing_angle = compute_grazing_angle(
            range_m,
            scenario.antenna_height_m,
            scenario.target_height_m,
        )
        grazing_angle = max(0.5, min(90.0, grazing_angle))

    # Propagation losses
    atmos_loss = 0.0
    rain_loss = 0.0

    if scenario.include_atmos_loss:
        atmos_loss = atmospheric_loss_db(
            scenario.freq_hz,
            range_m,
            elevation_deg=grazing_angle,
            temperature_c=scenario.temperature_c,
            humidity_pct=scenario.humidity_pct,
        )

    if scenario.rain_rate_mm_hr > 0:
        rain_loss = rain_attenuation_db(
            scenario.freq_hz,
            range_m,
            scenario.rain_rate_mm_hr,
        )

    # Total propagation loss
    propagation_loss_db = atmos_loss + rain_loss

    # Noise convention: rx_noise_temp_k is the ANTENNA temperature.
    # T_sys = T_ant + T0*(F-1), N = k*T_sys*B. Cascaded NF from context
    # (RF cascade model) wins over the flat arch.rf value, matching the
    # comms link budget.
    nf_raw = context.get("cascade_nf_db", arch.rf.noise_figure_db)
    nf_db = float(nf_raw) if isinstance(nf_raw, (int, float)) else arch.rf.noise_figure_db
    noise_factor = 10.0 ** (nf_db / 10.0)
    t_sys_k = scenario.rx_noise_temp_k + 290.0 * (noise_factor - 1.0)
    noise_power_dbw = W_TO_DBW(K_B * t_sys_k * scenario.bandwidth_hz)

    # Radar equation constant: (4π)^3 in dB
    radar_constant_db = 30 * math.log10(4 * math.pi)  # ≈ 32.98 dB

    # Single-pulse SNR (monostatic radar equation in dB)
    # SNR = Pt + 2*G + 2*λ_dB + σ - 4*R_dB - L - (4π)^3_dB - N - L_prop
    snr_single_db = (
        peak_power_dbw
        + 2 * g_ant_db
        + 2 * wavelength_db
        + rcs_dbsm
        - 4 * range_db
        - system_loss_db
        - radar_constant_db
        - noise_power_dbw
        - propagation_loss_db
    )

    # Clutter calculations
    clutter_rcs_dbsm = -100.0  # Default: no clutter
    scr_db = 100.0  # Default: no clutter (infinite SCR)

    if scenario.clutter_type != "none":
        # Compute resolution cell area/volume
        cell_area = compute_resolution_cell_area(range_m, range_resolution_m, beamwidth_az_deg)
        cell_volume = compute_resolution_volume(
            range_m, range_resolution_m, beamwidth_az_deg, beamwidth_el_deg
        )

        if scenario.clutter_type == "sea":
            clutter_rcs_dbsm = sea_clutter_rcs(
                scenario.sea_state,
                grazing_angle,
                scenario.freq_hz,
                cell_area,
                scenario.polarization,
            )
        elif scenario.clutter_type == "ground":
            clutter_rcs_dbsm = ground_clutter_rcs(
                scenario.terrain_type,
                grazing_angle,
                scenario.freq_hz,
                cell_area,
            )
        elif scenario.clutter_type == "rain":
            clutter_rcs_dbsm = rain_clutter_rcs(
                scenario.rain_rate_mm_hr,
                scenario.freq_hz,
                cell_volume,
            )

        scr_db = compute_scr(rcs_dbsm, clutter_rcs_dbsm)

    # MTI clutter suppression, when a canceller is configured. Without it a
    # ground-based radar looking at clutter is judged undetectable, which
    # misrepresents every real MTI system. Improvement factor rather than
    # clutter attenuation is applied, because I = G*CA carries both the
    # filter's gain on the target and its rejection of clutter, and it is
    # the SCR that the detection budget consumes.
    mti_improvement_db = 0.0
    if scenario.mti_n_pulse is not None and scenario.clutter_type != "none":
        from phased_array_systems.models.radar.mti import (
            clutter_spectral_std_hz,
            mti_improvement_factor,
            normalized_clutter_spread_rad,
        )

        if scenario.prf_hz is None:
            raise ValueError("mti_n_pulse requires prf_hz to be set")
        sigma_omega = normalized_clutter_spread_rad(
            clutter_spectral_std_hz(scenario.clutter_velocity_std_ms, scenario.wavelength_m),
            scenario.prf_hz,
        )
        mti_improvement_db = 10.0 * math.log10(
            mti_improvement_factor(scenario.mti_n_pulse, sigma_omega)
        )
        scr_db += mti_improvement_db

    # Compute SCNR (signal-to-clutter-plus-noise ratio)
    scnr_db = compute_scnr(snr_single_db, scr_db)

    # CFAR loss
    cfar_loss = 0.0
    if scenario.cfar_type != "none":
        cfar_loss = cfar_loss_db(
            scenario.cfar_type,
            scenario.cfar_ref_cells,
            scenario.pfa,
        )

    # Integration gain and required SNR must come from the same law to
    # avoid double-counting. Both use the exact detection statistics
    # (noncentral chi-square / gamma mixtures), so target fluctuation
    # affects the margin consistently; the implied noncoherent gain is
    # the drop in required single-pulse SNR vs n=1.
    n_pulses = scenario.n_pulses
    swerling = scenario.swerling
    snr_required_db = _required_snr(
        pd=scenario.pd_required,
        pfa=scenario.pfa,
        swerling=swerling,
        n_pulses=1,
    )
    if scenario.integration_type == "coherent":
        integration_gain_db = coherent_integration_gain(n_pulses)
    else:
        snr_required_single_db = _required_snr(
            pd=scenario.pd_required,
            pfa=scenario.pfa,
            swerling=swerling,
            n_pulses=n_pulses,
        )
        integration_gain_db = snr_required_db - snr_required_single_db

    # Integrated SCNR (use SCNR when clutter is present, SNR otherwise)
    effective_snr_single = scnr_db if scenario.clutter_type != "none" else snr_single_db

    snr_integrated_db = effective_snr_single + integration_gain_db - cfar_loss

    # SNR margin (snr_required_db is referenced to the integrated SNR)
    snr_margin_db = snr_integrated_db - snr_required_db

    # Achieved Pd from per-pulse SNR using exact n-pulse statistics
    pd_achieved = compute_pd_from_snr(
        effective_snr_single - cfar_loss,
        scenario.pfa,
        swerling=swerling,
        n_pulses=n_pulses,
        integration=scenario.integration_type,
    )

    # Detection range (range where margin = 0)
    # From radar equation: R^4 proportional to SNR
    # R_det / R = (SNR_integrated / SNR_required)^(1/4)
    # In dB: R_det = R * 10^(margin_dB / 40)
    detection_range_m = range_m * 10 ** (snr_margin_db / 40) if snr_margin_db > -40 else 0.0

    metrics: MetricsDict = {
        # Power
        "peak_power_w": peak_power_w,
        "peak_power_dbw": peak_power_dbw,
        # Antenna
        "g_ant_db": g_ant_db,
        # Target/Environment
        "wavelength_m": wavelength_m,
        "target_rcs_dbsm": rcs_dbsm,
        "target_rcs_m2": rcs_m2,
        "range_m": range_m,
        "grazing_angle_deg": grazing_angle,
        # Noise
        "noise_power_dbw": noise_power_dbw,
        "noise_temp_system_k": t_sys_k,
        "noise_figure_used_db": nf_db,
        "system_loss_db": system_loss_db,
        # Propagation losses
        "atmos_loss_db": atmos_loss,
        "rain_loss_db": rain_loss,
        "propagation_loss_db": propagation_loss_db,
        # Clutter
        "clutter_type": scenario.clutter_type,
        "clutter_rcs_dbsm": clutter_rcs_dbsm,
        "scr_db": scr_db,
        "scnr_db": scnr_db,
        # CFAR
        "cfar_type": scenario.cfar_type,
        "cfar_loss_db": cfar_loss,
        # SNR
        "snr_single_pulse_db": snr_single_db,
        "integration_gain_db": integration_gain_db,
        "snr_integrated_db": snr_integrated_db,
        "snr_required_db": snr_required_db,
        "snr_margin_db": snr_margin_db,
        # Detection
        "pd_achieved": pd_achieved,
        "pd_required": scenario.pd_required,
        "pfa": scenario.pfa,
        "swerling": swerling,
        "n_pulses": n_pulses,
        "integration_type": scenario.integration_type,
        "detection_range_m": detection_range_m,
    }

    # Emitted only when a canceller is configured, so a run without MTI
    # produces exactly the keys it did before.
    if scenario.mti_n_pulse is not None and scenario.clutter_type != "none":
        metrics["mti_improvement_db"] = mti_improvement_db
        metrics["mti_n_pulse"] = scenario.mti_n_pulse

    return metrics

Functions

compute_snr_for_pd

compute_snr_for_pd(pd: float, pfa: float, swerling: SwerlingModel = 0, n_pulses: int = 1, integration: Literal['coherent', 'noncoherent'] = 'noncoherent') -> float

Compute required SNR for given Pd and Pfa.

Inverse of compute_pd_from_snr using numerical root finding.

PARAMETER	DESCRIPTION
pd	Required probability of detection (0 < pd < 1) TYPE: float
pfa	Probability of false alarm (0 < pfa < 1) TYPE: float
swerling	Swerling target model (0-4) TYPE: SwerlingModel DEFAULT: 0
n_pulses	Number of pulses integrated TYPE: int DEFAULT: 1
integration	Integration type TYPE: Literal['coherent', 'noncoherent'] DEFAULT: 'noncoherent'

RETURNS	DESCRIPTION
float	Required single-pulse SNR in dB

Source code in src/phased_array_systems/models/radar/detection.py

def compute_snr_for_pd(
    pd: float,
    pfa: float,
    swerling: SwerlingModel = 0,
    n_pulses: int = 1,
    integration: Literal["coherent", "noncoherent"] = "noncoherent",
) -> float:
    """Compute required SNR for given Pd and Pfa.

    Inverse of compute_pd_from_snr using numerical root finding.

    Args:
        pd: Required probability of detection (0 < pd < 1)
        pfa: Probability of false alarm (0 < pfa < 1)
        swerling: Swerling target model (0-4)
        n_pulses: Number of pulses integrated
        integration: Integration type

    Returns:
        Required single-pulse SNR in dB
    """
    if not 0 < pd < 1:
        raise ValueError("pd must be between 0 and 1")
    if not 0 < pfa < 1:
        raise ValueError("pfa must be between 0 and 1")

    def objective(snr_db: float) -> float:
        pd_calc = compute_pd_from_snr(snr_db, pfa, swerling, n_pulses, integration)
        return pd_calc - pd

    # Use Albersheim as initial guess
    snr_guess = albersheim_snr(pd, pfa, n_pulses)

    try:
        result = optimize.brentq(objective, snr_guess - 20, snr_guess + 20)
        return float(result)
    except ValueError:
        # If brentq fails, return Albersheim estimate
        return snr_guess

compute_pd_from_snr

compute_pd_from_snr(snr_db: float, pfa: float, swerling: SwerlingModel = 0, n_pulses: int = 1, integration: Literal['coherent', 'noncoherent'] = 'noncoherent') -> float

Compute probability of detection for given per-pulse SNR.

Square-law detector statistics. Conditioned on the total received signal power s (in noise-power units), the normalized detector output follows a noncentral chi-square distribution with 2n degrees of freedom and noncentrality 2s, so Pd = Q_chi2'(2T; 2n, 2s) where T is the normalized threshold. Swerling fluctuation is the gamma-distributed mixture of s:

Swerling 0: s = n*SNR (deterministic; Pd is the Marcum Q result)
Swerling 1: s ~ Gamma(1, n*SNR)   (scan-to-scan Rayleigh)
Swerling 2: s ~ Gamma(n, SNR)     (pulse-to-pulse Rayleigh; closed form)
Swerling 3: s ~ Gamma(2, n*SNR/2) (scan-to-scan chi-4)
Swerling 4: s ~ Gamma(2n, SNR/2)  (pulse-to-pulse chi-4)

Coherent integration multiplies SNR by n and detects on a single sample; noncoherent integration uses the n-sample statistics directly (no separate empirical gain factor).

PARAMETER	DESCRIPTION
snr_db	Signal-to-noise ratio per pulse (dB) TYPE: float
pfa	Probability of false alarm TYPE: float
swerling	Swerling target model (0 = non-fluctuating) TYPE: SwerlingModel DEFAULT: 0
n_pulses	Number of pulses integrated TYPE: int DEFAULT: 1
integration	Integration type ("coherent" or "noncoherent") TYPE: Literal['coherent', 'noncoherent'] DEFAULT: 'noncoherent'

RETURNS	DESCRIPTION
float	Probability of detection (0-1)

Source code in src/phased_array_systems/models/radar/detection.py

def compute_pd_from_snr(
    snr_db: float,
    pfa: float,
    swerling: SwerlingModel = 0,
    n_pulses: int = 1,
    integration: Literal["coherent", "noncoherent"] = "noncoherent",
) -> float:
    """Compute probability of detection for given per-pulse SNR.

    Square-law detector statistics. Conditioned on the total received signal
    power s (in noise-power units), the normalized detector output follows a
    noncentral chi-square distribution with 2n degrees of freedom and
    noncentrality 2s, so Pd = Q_chi2'(2T; 2n, 2s) where T is the normalized
    threshold. Swerling fluctuation is the gamma-distributed mixture of s:

        Swerling 0: s = n*SNR (deterministic; Pd is the Marcum Q result)
        Swerling 1: s ~ Gamma(1, n*SNR)   (scan-to-scan Rayleigh)
        Swerling 2: s ~ Gamma(n, SNR)     (pulse-to-pulse Rayleigh; closed form)
        Swerling 3: s ~ Gamma(2, n*SNR/2) (scan-to-scan chi-4)
        Swerling 4: s ~ Gamma(2n, SNR/2)  (pulse-to-pulse chi-4)

    Coherent integration multiplies SNR by n and detects on a single sample;
    noncoherent integration uses the n-sample statistics directly (no separate
    empirical gain factor).

    Args:
        snr_db: Signal-to-noise ratio per pulse (dB)
        pfa: Probability of false alarm
        swerling: Swerling target model (0 = non-fluctuating)
        n_pulses: Number of pulses integrated
        integration: Integration type ("coherent" or "noncoherent")

    Returns:
        Probability of detection (0-1)
    """
    if not 0 < pfa < 1:
        raise ValueError("pfa must be between 0 and 1")
    if n_pulses < 1:
        raise ValueError("n_pulses must be >= 1")

    snr_linear = 10 ** (snr_db / 10)

    if integration == "coherent":
        # Coherent integration: full n-times SNR gain, single detection sample
        snr_linear *= n_pulses
        n = 1
    else:
        n = n_pulses

    threshold = compute_detection_threshold(pfa, n_samples=n)

    if swerling == 0:
        s = n * snr_linear
        pd = stats.ncx2.sf(2 * threshold, 2 * n, 2 * s)
    elif swerling == 2:
        # Sum of n independent exponential pulses with mean (1 + SNR)
        pd = special.gammaincc(n, threshold / (1 + snr_linear))
    elif swerling in (1, 3, 4):
        if swerling == 1:
            shape, scale = 1.0, n * snr_linear
        elif swerling == 3:
            shape, scale = 2.0, n * snr_linear / 2
        else:  # swerling == 4
            shape, scale = 2.0 * n, snr_linear / 2

        def integrand(s: float) -> float:
            return float(
                stats.ncx2.sf(2 * threshold, 2 * n, 2 * s) * stats.gamma.pdf(s, shape, scale=scale)
            )

        pd, _ = integrate.quad(integrand, 0, stats.gamma.ppf(1 - 1e-10, shape, scale=scale))
    else:
        raise ValueError(f"Unknown Swerling model: {swerling}")

    return max(0.0, min(1.0, float(pd)))

albersheim_snr

albersheim_snr(pd: float, pfa: float, n_pulses: int = 1) -> float

Albersheim's equation for required SNR (Swerling 0).

Empirical approximation valid for: - 0.1 <= Pd <= 0.99 - 1e-9 <= Pfa <= 1e-3 - 1 <= n_pulses <= 8096

PARAMETER	DESCRIPTION
pd	Probability of detection TYPE: float
pfa	Probability of false alarm TYPE: float
n_pulses	Number of pulses (non-coherent integration) TYPE: int DEFAULT: 1

RETURNS	DESCRIPTION
float	Required single-pulse SNR in dB

Source code in src/phased_array_systems/models/radar/detection.py

def albersheim_snr(
    pd: float,
    pfa: float,
    n_pulses: int = 1,
) -> float:
    """Albersheim's equation for required SNR (Swerling 0).

    Empirical approximation valid for:
    - 0.1 <= Pd <= 0.99
    - 1e-9 <= Pfa <= 1e-3
    - 1 <= n_pulses <= 8096

    Args:
        pd: Probability of detection
        pfa: Probability of false alarm
        n_pulses: Number of pulses (non-coherent integration)

    Returns:
        Required single-pulse SNR in dB
    """
    if not 0.1 <= pd <= 0.9999:
        raise ValueError("pd must be between 0.1 and 0.9999 for Albersheim")
    if not 1e-10 <= pfa <= 0.1:
        raise ValueError("pfa must be between 1e-10 and 0.1 for Albersheim")
    if n_pulses < 1:
        raise ValueError("n_pulses must be >= 1")

    # Albersheim's equation
    A = math.log(0.62 / pfa)
    B = math.log(pd / (1 - pd))

    # SNR required for n pulses (non-coherent integration)
    snr_n_db = -5 * math.log10(n_pulses) + (6.2 + 4.54 / math.sqrt(n_pulses + 0.44)) * math.log10(
        A + 0.12 * A * B + 1.7 * B
    )

    return snr_n_db

coherent_integration_gain

coherent_integration_gain(n_pulses: int) -> float

Coherent integration gain in dB.

Coherent integration (phase-preserving) provides full N-times improvement in SNR because signals add coherently while noise adds incoherently.

PARAMETER	DESCRIPTION
n_pulses	Number of pulses integrated (must be >= 1) TYPE: int

RETURNS	DESCRIPTION
float	Integration gain in dB: 10 * log10(n_pulses)

RAISES	DESCRIPTION
ValueError	If n_pulses < 1

Source code in src/phased_array_systems/models/radar/integration.py

def coherent_integration_gain(n_pulses: int) -> float:
    """Coherent integration gain in dB.

    Coherent integration (phase-preserving) provides full N-times
    improvement in SNR because signals add coherently while noise
    adds incoherently.

    Args:
        n_pulses: Number of pulses integrated (must be >= 1)

    Returns:
        Integration gain in dB: 10 * log10(n_pulses)

    Raises:
        ValueError: If n_pulses < 1
    """
    if n_pulses < 1:
        raise ValueError("n_pulses must be >= 1")

    if n_pulses == 1:
        return 0.0

    return 10 * math.log10(n_pulses)

noncoherent_integration_gain

noncoherent_integration_gain(n_pulses: int, pd: float = 0.9, pfa: float = 1e-06) -> float

Non-coherent integration gain in dB.

Non-coherent integration (magnitude-only) provides less than full N-times gain because both signal and noise magnitudes are combined. The efficiency depends on SNR and Pd/Pfa.

Uses empirical approximation: gain ≈ 10 * log10(n^efficiency) where efficiency ≈ 0.8 for typical radar parameters.

PARAMETER	DESCRIPTION
n_pulses	Number of pulses integrated (must be >= 1) TYPE: int
pd	Probability of detection (affects efficiency) TYPE: float DEFAULT: 0.9
pfa	Probability of false alarm (affects efficiency) TYPE: float DEFAULT: 1e-06

RETURNS	DESCRIPTION
float	Integration gain in dB (always <= coherent gain)

RAISES	DESCRIPTION
ValueError	If n_pulses < 1

Source code in src/phased_array_systems/models/radar/integration.py

def noncoherent_integration_gain(
    n_pulses: int,
    pd: float = 0.9,
    pfa: float = 1e-6,
) -> float:
    """Non-coherent integration gain in dB.

    Non-coherent integration (magnitude-only) provides less than
    full N-times gain because both signal and noise magnitudes
    are combined. The efficiency depends on SNR and Pd/Pfa.

    Uses empirical approximation: gain ≈ 10 * log10(n^efficiency)
    where efficiency ≈ 0.8 for typical radar parameters.

    Args:
        n_pulses: Number of pulses integrated (must be >= 1)
        pd: Probability of detection (affects efficiency)
        pfa: Probability of false alarm (affects efficiency)

    Returns:
        Integration gain in dB (always <= coherent gain)

    Raises:
        ValueError: If n_pulses < 1
    """
    if n_pulses < 1:
        raise ValueError("n_pulses must be >= 1")

    if n_pulses == 1:
        return 0.0

    # Efficiency factor depends on operating point
    # Higher Pd requires higher SNR, reducing integration efficiency
    if pd >= 0.99:
        efficiency = 0.7
    elif pd >= 0.9:
        efficiency = 0.8
    elif pd >= 0.5:
        efficiency = 0.85
    else:
        efficiency = 0.9

    # Non-coherent gain: approximately n^efficiency
    return 10 * efficiency * math.log10(n_pulses)

Output Metrics

Metric	Units	Description
snr_single_pulse_db	dB	Single-pulse SNR
snr_integrated_db	dB	SNR after integration
snr_required_db	dB	Required SNR for Pd/Pfa
snr_margin_db	dB	Margin above required
detection_range_m	m	Maximum detection range
integration_gain_db	dB	Gain from pulse integration

Usage Examples

Using RadarModel

from phased_array_systems.models.radar import RadarModel
from phased_array_systems.scenarios import RadarDetectionScenario

scenario = RadarDetectionScenario(
    freq_hz=10e9,
    bandwidth_hz=100e3,
    range_m=100e3,
    target_rcs_dbsm=0.0,
    pd_required=0.9,
    pfa=1e-6,
    prf_hz=1000,
    n_pulses=10,
    integration_type="coherent",
    swerling=1,
)

model = RadarModel()
metrics = model.evaluate(arch, scenario, context={})

print(f"Single-Pulse SNR: {metrics['snr_single_pulse_db']:.1f} dB")
print(f"Integrated SNR: {metrics['snr_integrated_db']:.1f} dB")
print(f"SNR Margin: {metrics['snr_margin_db']:.1f} dB")

Computing Required SNR

from phased_array_systems.models.radar import compute_snr_for_pd

snr_req = compute_snr_for_pd(
    pd=0.9,
    pfa=1e-6,
    swerling=1,
)
print(f"Required SNR: {snr_req:.1f} dB")

Computing Detection Probability

from phased_array_systems.models.radar import compute_pd_from_snr

pd = compute_pd_from_snr(
    snr_db=15.0,
    pfa=1e-6,
    swerling=1,
)
print(f"Detection Probability: {pd:.3f}")

Integration Gain

from phased_array_systems.models.radar import coherent_integration_gain, noncoherent_integration_gain

# Coherent integration
gain_coherent = coherent_integration_gain(n_pulses=16)
print(f"Coherent Gain: {gain_coherent:.1f} dB")  # 12.0 dB

# Non-coherent integration
gain_noncoherent = noncoherent_integration_gain(n_pulses=16)
print(f"Non-coherent Gain: {gain_noncoherent:.1f} dB")  # ~6.0 dB

Radar Range Equation

\[ SNR = \frac{P_t G^2 \lambda^2 \sigma}{(4\pi)^3 R^4 k T_s B_n L_s} \]

Where:

• \(P_t\) = Peak transmit power (W)
• \(G\) = Antenna gain (linear)
• \(\lambda\) = Wavelength (m)
• \(\sigma\) = Target RCS (m²)
• \(R\) = Target range (m)
• \(k\) = Boltzmann constant
• \(T_s\) = System noise temperature (K)
• \(B_n\) = Noise bandwidth (Hz)
• \(L_s\) = System losses (linear)

Swerling Models

Model	PDF	Decorrelation
0	Constant	None
1	Rayleigh	Scan-to-scan
2	Rayleigh	Pulse-to-pulse
3	Chi-squared (4 DOF)	Scan-to-scan
4	Chi-squared (4 DOF)	Pulse-to-pulse

MTI Clutter Suppression

Doppler-domain clutter rejection. See Theory: MTI Clutter Suppression.

clutter_spectral_std_hz

clutter_spectral_std_hz(clutter_velocity_std_ms: float, wavelength_m: float) -> float

sigma_c = 2 sigma_v / lambda (Skolnik ch. 15).

The velocity spread is a property of the clutter, not the radar, so the same wooded hillside produces a wider Doppler spectrum at higher frequency.

Source code in src/phased_array_systems/models/radar/mti.py

def clutter_spectral_std_hz(clutter_velocity_std_ms: float, wavelength_m: float) -> float:
    """sigma_c = 2 sigma_v / lambda (Skolnik ch. 15).

    The velocity spread is a property of the clutter, not the radar, so the
    same wooded hillside produces a wider Doppler spectrum at higher frequency.
    """
    if clutter_velocity_std_ms < 0:
        raise ValueError("clutter_velocity_std_ms must be >= 0")
    if wavelength_m <= 0:
        raise ValueError("wavelength_m must be > 0")
    return float(2.0 * clutter_velocity_std_ms / wavelength_m)

normalized_clutter_spread_rad

normalized_clutter_spread_rad(clutter_std_hz: float, prf_hz: float) -> float

sigma_omega = 2 pi sigma_c / PRF, the spread in normalized angular frequency.

This is the only clutter quantity the canceller math needs: everything downstream depends on the spectrum's width relative to the PRF, not on its absolute width.

Source code in src/phased_array_systems/models/radar/mti.py

def normalized_clutter_spread_rad(clutter_std_hz: float, prf_hz: float) -> float:
    """sigma_omega = 2 pi sigma_c / PRF, the spread in normalized angular frequency.

    This is the only clutter quantity the canceller math needs: everything
    downstream depends on the spectrum's width relative to the PRF, not on its
    absolute width.
    """
    if clutter_std_hz < 0:
        raise ValueError("clutter_std_hz must be >= 0")
    if prf_hz <= 0:
        raise ValueError("prf_hz must be > 0")
    return float(2.0 * math.pi * clutter_std_hz / prf_hz)

clutter_autocorrelation

clutter_autocorrelation(sigma_omega_rad: float, lag: int) -> float

rho_c[k] = exp(-(sigma_omega k)^2 / 2), Richards FRSP Eq. (5.53).

The normalized autocorrelation of a Gaussian clutter spectrum, valid for sigma_omega << pi. At the wide-spectrum limit the approximation breaks down along with the premise that the clutter is narrowband relative to the PRF.

Source code in src/phased_array_systems/models/radar/mti.py

def clutter_autocorrelation(sigma_omega_rad: float, lag: int) -> float:
    """rho_c[k] = exp(-(sigma_omega k)^2 / 2), Richards FRSP Eq. (5.53).

    The normalized autocorrelation of a Gaussian clutter spectrum, valid for
    sigma_omega << pi. At the wide-spectrum limit the approximation breaks down
    along with the premise that the clutter is narrowband relative to the PRF.
    """
    if sigma_omega_rad < 0:
        raise ValueError("sigma_omega_rad must be >= 0")
    return float(math.exp(-((sigma_omega_rad * lag) ** 2) / 2.0))

canceller_weights

canceller_weights(n_pulse: int) -> list[int]

Binomial (N-1)th-difference canceller weights w_k = (-1)^k C(N-1, k).

N = 2 gives [1, -1] and N = 3 gives [1, -2, 1], the conventional two- and three-pulse cancellers.

Source code in src/phased_array_systems/models/radar/mti.py

def canceller_weights(n_pulse: int) -> list[int]:
    """Binomial (N-1)th-difference canceller weights w_k = (-1)^k C(N-1, k).

    N = 2 gives [1, -1] and N = 3 gives [1, -2, 1], the conventional two- and
    three-pulse cancellers.
    """
    if n_pulse < 2:
        raise ValueError("n_pulse must be >= 2")
    return [(-1) ** k * comb(n_pulse - 1, k) for k in range(n_pulse)]

mti_signal_gain

mti_signal_gain(n_pulse: int) -> float

Average signal gain over all Doppler shifts, G = sum_k w_k^2.

Richards FRSP p. 246 defines the gain as the mean of |H(F)|^2 over the unambiguous Doppler band, which for an FIR filter is the sum of the squared weights by Parseval. Gives G = 2 (3.0 dB) for the two-pulse canceller and G = 6 (7.8 dB) for the three-pulse, matching FRSP p. 247.

The target velocity is assumed unknown a priori; a radar that knows where to look does better than this average.

Source code in src/phased_array_systems/models/radar/mti.py

def mti_signal_gain(n_pulse: int) -> float:
    """Average signal gain over all Doppler shifts, G = sum_k w_k^2.

    Richards FRSP p. 246 defines the gain as the mean of |H(F)|^2 over the
    unambiguous Doppler band, which for an FIR filter is the sum of the squared
    weights by Parseval. Gives G = 2 (3.0 dB) for the two-pulse canceller and
    G = 6 (7.8 dB) for the three-pulse, matching FRSP p. 247.

    The target velocity is assumed unknown a priori; a radar that knows where
    to look does better than this average.
    """
    return float(sum(w * w for w in canceller_weights(n_pulse)))

mti_improvement_factor

mti_improvement_factor(n_pulse: int, sigma_omega_rad: float) -> float

Improvement factor I = G * CA for an N-pulse binomial canceller.

Computed from the general quadratic form

I = sum_k w_k^2 / sum_i sum_j w_i w_j rho_c[i-j]

rather than from the published closed forms, which it reproduces exactly: Richards FRSP Eq. (5.52) gives 1/(1 - rho[1]) for N = 2 and Eq. (5.54) gives 1/(1 - (4/3) rho[1] + (1/3) rho[2]) for N = 3. Both are asserted against this function in the oracle tests.

Returned as a linear ratio. The value grows without bound as the clutter spectrum narrows, which is physical -- perfectly stationary clutter is perfectly cancellable -- but a design should not lean on figures far beyond the system's phase noise and stability limits, which this model does not represent.

Source code in src/phased_array_systems/models/radar/mti.py

def mti_improvement_factor(n_pulse: int, sigma_omega_rad: float) -> float:
    """Improvement factor I = G * CA for an N-pulse binomial canceller.

    Computed from the general quadratic form

        I = sum_k w_k^2 / sum_i sum_j w_i w_j rho_c[i-j]

    rather than from the published closed forms, which it reproduces exactly:
    Richards FRSP Eq. (5.52) gives 1/(1 - rho[1]) for N = 2 and Eq. (5.54)
    gives 1/(1 - (4/3) rho[1] + (1/3) rho[2]) for N = 3. Both are asserted
    against this function in the oracle tests.

    Returned as a linear ratio. The value grows without bound as the clutter
    spectrum narrows, which is physical -- perfectly stationary clutter is
    perfectly cancellable -- but a design should not lean on figures far beyond
    the system's phase noise and stability limits, which this model does not
    represent.
    """
    weights = canceller_weights(n_pulse)
    gain = sum(w * w for w in weights)
    residue = sum(
        weights[i] * weights[j] * clutter_autocorrelation(sigma_omega_rad, abs(i - j))
        for i in range(n_pulse)
        for j in range(n_pulse)
    )
    if residue <= 0:
        raise ValueError("clutter residue is non-positive; sigma_omega out of valid range")
    return float(gain / residue)

mti_clutter_attenuation

mti_clutter_attenuation(n_pulse: int, sigma_omega_rad: float) -> float

CA = I / G, the clutter power ratio across the filter (FRSP Eq. 5.43).

Distinct from the improvement factor: this is rejection alone, with no credit for the filter's gain on the target.

Source code in src/phased_array_systems/models/radar/mti.py

def mti_clutter_attenuation(n_pulse: int, sigma_omega_rad: float) -> float:
    """CA = I / G, the clutter power ratio across the filter (FRSP Eq. 5.43).

    Distinct from the improvement factor: this is rejection alone, with no
    credit for the filter's gain on the target.
    """
    return float(mti_improvement_factor(n_pulse, sigma_omega_rad) / mti_signal_gain(n_pulse))

required_clutter_attenuation_db

required_clutter_attenuation_db(target_rcs_dbsm: float, clutter_rcs_dbsm: float, required_scr_db: float) -> float

CA_req = (S/C)_required - (sigma_target - sigma_clutter), in dB.

How much clutter power must be removed for the target to clear the required signal-to-clutter ratio. Positive means suppression is needed.

Source code in src/phased_array_systems/models/radar/mti.py

def required_clutter_attenuation_db(
    target_rcs_dbsm: float,
    clutter_rcs_dbsm: float,
    required_scr_db: float,
) -> float:
    """CA_req = (S/C)_required - (sigma_target - sigma_clutter), in dB.

    How much clutter power must be removed for the target to clear the required
    signal-to-clutter ratio. Positive means suppression is needed.
    """
    return float(required_scr_db - (target_rcs_dbsm - clutter_rcs_dbsm))

blind_speed_ms

blind_speed_ms(prf_hz: float, wavelength_m: float, harmonic: int = 1) -> float

v_blind = n PRF lambda / 2: target speeds the canceller nulls along with clutter.

A target at a blind speed produces the same phase advance per pulse as stationary clutter and is cancelled with it. The first blind speed bounds the useful Doppler coverage of a single-PRF MTI.

Source code in src/phased_array_systems/models/radar/mti.py

def blind_speed_ms(prf_hz: float, wavelength_m: float, harmonic: int = 1) -> float:
    """v_blind = n PRF lambda / 2: target speeds the canceller nulls along with clutter.

    A target at a blind speed produces the same phase advance per pulse as
    stationary clutter and is cancelled with it. The first blind speed bounds
    the useful Doppler coverage of a single-PRF MTI.
    """
    if prf_hz <= 0:
        raise ValueError("prf_hz must be > 0")
    if wavelength_m <= 0:
        raise ValueError("wavelength_m must be > 0")
    if harmonic < 1:
        raise ValueError("harmonic must be >= 1")
    return float(harmonic * prf_hz * wavelength_m / 2.0)

doppler_shift_hz

doppler_shift_hz(radial_velocity_ms: float, wavelength_m: float) -> float

f_d = 2 v_r / lambda, the monostatic two-way Doppler shift.

Source code in src/phased_array_systems/models/radar/mti.py

def doppler_shift_hz(radial_velocity_ms: float, wavelength_m: float) -> float:
    """f_d = 2 v_r / lambda, the monostatic two-way Doppler shift."""
    if wavelength_m <= 0:
        raise ValueError("wavelength_m must be > 0")
    return float(2.0 * radial_velocity_ms / wavelength_m)

unambiguous_range_m

unambiguous_range_m(prf_hz: float) -> float

R_ua = c / (2 PRF).

Source code in src/phased_array_systems/models/radar/mti.py

def unambiguous_range_m(prf_hz: float) -> float:
    """R_ua = c / (2 PRF)."""
    if prf_hz <= 0:
        raise ValueError("prf_hz must be > 0")
    return float(C / (2.0 * prf_hz))

Track Accuracy

Closed-form measurement error and steady-state track performance. See Theory: Track Accuracy for the equations and the SNR-convention note.

range_sigma_m

range_sigma_m(snr_db: float, bandwidth_hz: float, resolution_alpha: float = 1.0) -> float

sigma_R = dR / sqrt(2 * SNR), Curry Eq. (8.6) / POMR Eq. (18.33).

Thermal (SNR-driven) term only. Curry Eq. (8.5) adds fixed-random and bias terms in quadrature; those are hardware assertions, not consequences of the design, so they are left to the caller.

Source code in src/phased_array_systems/models/radar/tracking.py

def range_sigma_m(
    snr_db: float,
    bandwidth_hz: float,
    resolution_alpha: float = 1.0,
) -> float:
    """sigma_R = dR / sqrt(2 * SNR), Curry Eq. (8.6) / POMR Eq. (18.33).

    Thermal (SNR-driven) term only. Curry Eq. (8.5) adds fixed-random and bias
    terms in quadrature; those are hardware assertions, not consequences of the
    design, so they are left to the caller.
    """
    snr = _snr_linear(snr_db)
    if snr <= 0:
        raise ValueError("snr_db must give positive linear SNR")
    return float(range_resolution_m(bandwidth_hz, resolution_alpha) / math.sqrt(2.0 * snr))

angle_sigma_deg

angle_sigma_deg(snr_db: float, beamwidth_deg: float, monopulse_slope: float = DEFAULT_MONOPULSE_SLOPE) -> float

sigma_theta = theta_3dB / (k_m * sqrt(2 * SNR)), POMR Eq. (18.63) / Curry Eq. (8.8).

Valid above MONOPULSE_SNR_FLOOR_DB; see the module note.

Thermal (SNR-driven) term only, as with :func:range_sigma_m. Curry Eq. (8.7) adds fixed-random and bias terms in quadrature, and target glint can dominate angular error at short range (Curry p. 170). None of those are consequences of the array design, so they are left to the caller.

Source code in src/phased_array_systems/models/radar/tracking.py

def angle_sigma_deg(
    snr_db: float,
    beamwidth_deg: float,
    monopulse_slope: float = DEFAULT_MONOPULSE_SLOPE,
) -> float:
    """sigma_theta = theta_3dB / (k_m * sqrt(2 * SNR)), POMR Eq. (18.63) / Curry Eq. (8.8).

    Valid above ``MONOPULSE_SNR_FLOOR_DB``; see the module note.

    Thermal (SNR-driven) term only, as with :func:`range_sigma_m`. Curry
    Eq. (8.7) adds fixed-random and bias terms in quadrature, and target glint
    can dominate angular error at short range (Curry p. 170). None of those are
    consequences of the array design, so they are left to the caller.
    """
    if beamwidth_deg <= 0:
        raise ValueError("beamwidth_deg must be > 0")
    if monopulse_slope <= 0:
        raise ValueError("monopulse_slope must be > 0")
    snr = _snr_linear(snr_db)
    return float(beamwidth_deg / (monopulse_slope * math.sqrt(2.0 * snr)))

scan_broadened_beamwidth_deg

scan_broadened_beamwidth_deg(beamwidth_deg: float, scan_angle_deg: float) -> float

theta_phi = theta_B / cos(phi), Curry Eq. (8.9).

A phased array's beam broadens off broadside, so angle accuracy degrades with scan angle even at constant SNR.

Source code in src/phased_array_systems/models/radar/tracking.py

def scan_broadened_beamwidth_deg(beamwidth_deg: float, scan_angle_deg: float) -> float:
    """theta_phi = theta_B / cos(phi), Curry Eq. (8.9).

    A phased array's beam broadens off broadside, so angle accuracy degrades
    with scan angle even at constant SNR.
    """
    if abs(scan_angle_deg) >= 90.0:
        raise ValueError("scan_angle_deg must be within +/-90 degrees")
    return float(beamwidth_deg / math.cos(math.radians(scan_angle_deg)))

crossrange_sigma_m

crossrange_sigma_m(range_m: float, angle_sigma_deg_value: float) -> float

sigma_D = R * sigma_A, Curry Eq. (8.10). Angle in degrees, result in metres.

Source code in src/phased_array_systems/models/radar/tracking.py

def crossrange_sigma_m(range_m: float, angle_sigma_deg_value: float) -> float:
    """sigma_D = R * sigma_A, Curry Eq. (8.10). Angle in degrees, result in metres."""
    if range_m < 0:
        raise ValueError("range_m must be >= 0")
    return float(range_m * math.radians(angle_sigma_deg_value))

velocity_sigma_ms

velocity_sigma_ms(snr_db: float, coherent_dwell_s: float, wavelength_m: float) -> float

sigma_V = lambda / (2 tau sqrt(2 SNR)), Curry Eq. (8.13).

Curry attributes the form to Barton & Ward, Handbook of Radar Measurement, pp. 101-103. Equivalent to POMR Eq. (18.31) under the convention noted in the module docstring.

Source code in src/phased_array_systems/models/radar/tracking.py

def velocity_sigma_ms(
    snr_db: float,
    coherent_dwell_s: float,
    wavelength_m: float,
) -> float:
    """sigma_V = lambda / (2 tau sqrt(2 SNR)), Curry Eq. (8.13).

    Curry attributes the form to Barton & Ward, *Handbook of Radar
    Measurement*, pp. 101-103. Equivalent to POMR Eq. (18.31) under the
    convention noted in the module docstring.
    """
    if coherent_dwell_s <= 0:
        raise ValueError("coherent_dwell_s must be > 0")
    if wavelength_m <= 0:
        raise ValueError("wavelength_m must be > 0")
    snr = _snr_linear(snr_db)
    return float(wavelength_m / (2.0 * coherent_dwell_s * math.sqrt(2.0 * snr)))

combine_angle_errors_deg

combine_angle_errors_deg(*sigma_deg: float) -> float

Root-sum-square of independent angle error terms.

The seam that lets a thermal-noise angle error combine with the hardware pointing error from models/antenna/errors.py: phase-shifter bits and calibration residue then propagate all the way to track accuracy, which is the connection no tracking library can make.

Source code in src/phased_array_systems/models/radar/tracking.py

def combine_angle_errors_deg(*sigma_deg: float) -> float:
    """Root-sum-square of independent angle error terms.

    The seam that lets a thermal-noise angle error combine with the hardware
    pointing error from ``models/antenna/errors.py``: phase-shifter bits and
    calibration residue then propagate all the way to track accuracy, which is
    the connection no tracking library can make.
    """
    total_sq = 0.0
    for term in sigma_deg:
        if term < 0:
            raise ValueError("angle error terms must be >= 0")
        total_sq += term * term
    return float(math.sqrt(total_sq))

tracking_index

tracking_index(sigma_v: float, sigma_w: float, revisit_s: float) -> float

Random tracking index Gamma = sigma_v T^2 / sigma_w, POMR Eq. (19.47).

Kalata's parameter: the ratio of position uncertainty from target maneuverability to that from the sensor measurement. It is the single number that sets the steady-state filter.

Source code in src/phased_array_systems/models/radar/tracking.py

def tracking_index(sigma_v: float, sigma_w: float, revisit_s: float) -> float:
    """Random tracking index Gamma = sigma_v T^2 / sigma_w, POMR Eq. (19.47).

    Kalata's parameter: the ratio of position uncertainty from target
    maneuverability to that from the sensor measurement. It is the single
    number that sets the steady-state filter.
    """
    if sigma_v < 0:
        raise ValueError("sigma_v must be >= 0")
    if sigma_w <= 0:
        raise ValueError("sigma_w must be > 0")
    if revisit_s <= 0:
        raise ValueError("revisit_s must be > 0")
    return float(sigma_v * revisit_s**2 / sigma_w)

deterministic_tracking_index

deterministic_tracking_index(accel_max_ms2: float, revisit_s: float, sigma_w: float) -> float

Gamma_D = A_max T^2 / sigma_w, POMR Eq. (19.59).

Source code in src/phased_array_systems/models/radar/tracking.py

def deterministic_tracking_index(
    accel_max_ms2: float,
    revisit_s: float,
    sigma_w: float,
) -> float:
    """Gamma_D = A_max T^2 / sigma_w, POMR Eq. (19.59)."""
    if accel_max_ms2 < 0:
        raise ValueError("accel_max_ms2 must be >= 0")
    if revisit_s <= 0:
        raise ValueError("revisit_s must be > 0")
    if sigma_w <= 0:
        raise ValueError("sigma_w must be > 0")
    return float(accel_max_ms2 * revisit_s**2 / sigma_w)

process_noise_from_maneuver

process_noise_from_maneuver(gamma_d: float, accel_max_ms2: float) -> float

sigma_v = kappa_1_min(Gamma_D) * A_max, POMR Eqs. (19.63)/(19.66).

kappa_1_min = 0.87 - 0.09 log10(Gamma_D) - 0.02 [log10(Gamma_D)]^2

Lets the caller state a physical maneuver ("the target pulls 4 g") instead of tuning a process-noise variance. POMR fits the curve over 0.01 <= Gamma_D <= 10; outside that band the fit is extrapolated and the caller should treat the result as indicative.

Source code in src/phased_array_systems/models/radar/tracking.py

def process_noise_from_maneuver(gamma_d: float, accel_max_ms2: float) -> float:
    """sigma_v = kappa_1_min(Gamma_D) * A_max, POMR Eqs. (19.63)/(19.66).

    kappa_1_min = 0.87 - 0.09 log10(Gamma_D) - 0.02 [log10(Gamma_D)]^2

    Lets the caller state a physical maneuver ("the target pulls 4 g") instead
    of tuning a process-noise variance. POMR fits the curve over
    0.01 <= Gamma_D <= 10; outside that band the fit is extrapolated and the
    caller should treat the result as indicative.
    """
    if gamma_d <= 0:
        raise ValueError("gamma_d must be > 0")
    if accel_max_ms2 < 0:
        raise ValueError("accel_max_ms2 must be >= 0")
    log_gd = math.log10(gamma_d)
    kappa = 0.87 - 0.09 * log_gd - 0.02 * log_gd**2
    return float(kappa * accel_max_ms2)

alpha_beta_gains

alpha_beta_gains(gamma: float) -> tuple[float, float]

Steady-state alpha-beta gains from the tracking index, POMR Eqs. (19.54)/(19.55).

Satisfies the Kalata relation beta = 2(2 - alpha) - 4 sqrt(1 - alpha) (Eq. 19.56) and inverts exactly through Gamma = beta / sqrt(1 - alpha).

Computed through a rearrangement rather than from Eqs. (19.54)/(19.55) literally, because those forms lose precision as alpha approaches 1: they build alpha from a difference of large terms and the caller then needs 1 - alpha, so the relative error in beta reaches 4e-5 by Gamma = 1e5. Substituting r = sqrt(1 - alpha) into Eq. (19.47) with Eq. (19.56) gives

Gamma = 2 (1 - r)^2 / r   ->   2 r^2 - (4 + Gamma) r + 2 = 0
r = (4 + Gamma - sqrt(Gamma^2 + 8 Gamma)) / 4      (smaller root)
alpha = 1 - r^2,   beta = 2 (1 - r)^2

which is the same solution with 1 - alpha carried exactly as r^2. It holds the identity to ~1e-12 at Gamma = 1e5. Both forms are pinned against each other in the oracle tests; this is presentation, not a different model.

Source code in src/phased_array_systems/models/radar/tracking.py

def alpha_beta_gains(gamma: float) -> tuple[float, float]:
    """Steady-state alpha-beta gains from the tracking index, POMR Eqs. (19.54)/(19.55).

    Satisfies the Kalata relation beta = 2(2 - alpha) - 4 sqrt(1 - alpha)
    (Eq. 19.56) and inverts exactly through Gamma = beta / sqrt(1 - alpha).

    Computed through a rearrangement rather than from Eqs. (19.54)/(19.55)
    literally, because those forms lose precision as alpha approaches 1: they
    build alpha from a difference of large terms and the caller then needs
    1 - alpha, so the relative error in beta reaches 4e-5 by Gamma = 1e5.
    Substituting r = sqrt(1 - alpha) into Eq. (19.47) with Eq. (19.56) gives

        Gamma = 2 (1 - r)^2 / r   ->   2 r^2 - (4 + Gamma) r + 2 = 0
        r = (4 + Gamma - sqrt(Gamma^2 + 8 Gamma)) / 4      (smaller root)
        alpha = 1 - r^2,   beta = 2 (1 - r)^2

    which is the same solution with 1 - alpha carried exactly as r^2. It holds
    the identity to ~1e-12 at Gamma = 1e5. Both forms are pinned against each
    other in the oracle tests; this is presentation, not a different model.
    """
    if gamma < 0:
        raise ValueError("gamma must be >= 0")
    if gamma == 0:
        return 0.0, 0.0
    r = (4.0 + gamma - math.sqrt(gamma**2 + 8.0 * gamma)) / 4.0
    return float(1.0 - r * r), float(2.0 * (1.0 - r) ** 2)

steady_state_sigmas

steady_state_sigmas(sigma_w: float, alpha: float, beta: float, revisit_s: float) -> tuple[float, float]

Steady-state position and velocity RMS error, POMR Eq. (19.53).

P = sigma_w^2 * [[alpha, beta/T ], [beta/T, beta(2 alpha - beta)/(2(1-alpha)T^2)]]

so sigma_pos = sigma_w sqrt(alpha) and sigma_vel = (sigma_w/T) sqrt(beta(2 alpha - beta)/(2(1 - alpha))).

This is the total steady-state error, process noise included: it matches the fixed-gain covariance recursion iterated to convergence with Q present. For the no-maneuver figure use :func:variance_reduction_position.

Source code in src/phased_array_systems/models/radar/tracking.py

def steady_state_sigmas(
    sigma_w: float,
    alpha: float,
    beta: float,
    revisit_s: float,
) -> tuple[float, float]:
    """Steady-state position and velocity RMS error, POMR Eq. (19.53).

    P = sigma_w^2 * [[alpha,            beta/T                     ],
                     [beta/T,  beta(2 alpha - beta)/(2(1-alpha)T^2)]]

    so sigma_pos = sigma_w sqrt(alpha) and
    sigma_vel = (sigma_w/T) sqrt(beta(2 alpha - beta)/(2(1 - alpha))).

    This is the *total* steady-state error, process noise included: it matches
    the fixed-gain covariance recursion iterated to convergence with Q present.
    For the no-maneuver figure use :func:`variance_reduction_position`.
    """
    if sigma_w <= 0:
        raise ValueError("sigma_w must be > 0")
    if revisit_s <= 0:
        raise ValueError("revisit_s must be > 0")
    if not 0.0 <= alpha < 1.0:
        raise ValueError("alpha must satisfy 0 <= alpha < 1")
    if beta < 0:
        raise ValueError("beta must be >= 0")
    sigma_pos = sigma_w * math.sqrt(alpha)
    vel_var = beta * (2.0 * alpha - beta) / (2.0 * (1.0 - alpha) * revisit_s**2)
    sigma_vel = sigma_w * math.sqrt(max(0.0, vel_var))
    return float(sigma_pos), float(sigma_vel)

variance_reduction_position

variance_reduction_position(alpha: float, beta: float) -> float

Sensor-noise-only VRR, Mahafza Eq. (11.94).

(VRR)_x = (2 alpha^2 - 3 alpha beta + 2 beta) / (alpha (4 - 2 alpha - beta))

Distinct from :func:steady_state_sigmas, and both are correct: this is the variance ratio with the process noise removed (no maneuver), verified against the Q = 0 covariance recursion. It answers "how much does filtering reduce measurement noise", while Eq. (19.53) answers "how well is the target actually located". A third form circulating in the literature, (2 alpha^2 + 2 beta + alpha beta)/(alpha(4 - 2 alpha - beta)), is wrong: it exceeds unity, i.e. claims filtering amplifies noise.

Source code in src/phased_array_systems/models/radar/tracking.py

def variance_reduction_position(alpha: float, beta: float) -> float:
    """Sensor-noise-only VRR, Mahafza Eq. (11.94).

    (VRR)_x = (2 alpha^2 - 3 alpha beta + 2 beta) / (alpha (4 - 2 alpha - beta))

    Distinct from :func:`steady_state_sigmas`, and both are correct: this is
    the variance ratio with the process noise removed (no maneuver), verified
    against the Q = 0 covariance recursion. It answers "how much does filtering
    reduce measurement noise", while Eq. (19.53) answers "how well is the
    target actually located". A third form circulating in the literature,
    (2 alpha^2 + 2 beta + alpha beta)/(alpha(4 - 2 alpha - beta)), is wrong: it
    exceeds unity, i.e. claims filtering amplifies noise.
    """
    denom = alpha * (4.0 - 2.0 * alpha - beta)
    if denom <= 0:
        raise ValueError("alpha, beta outside the stable region")
    return float((2.0 * alpha**2 - 3.0 * alpha * beta + 2.0 * beta) / denom)

maneuver_lag_m

maneuver_lag_m(sigma_w: float, alpha: float, beta: float, gamma_d: float) -> float

Maximum position MSE under a sustained maneuver, POMR Eq. (19.60).

MMSE_p = sigma_w^2 [ (2 alpha^2 + beta(2 - 3 alpha))/(alpha(4 - 2 alpha - beta)) + (1 - alpha)^2 Gamma_D^2 / beta^2 ]

Returned as an RMS distance. The second term is the deterministic lag: a filter tuned for a quiet target falls progressively behind a maneuvering one, and no amount of SNR fixes it.

Source code in src/phased_array_systems/models/radar/tracking.py

def maneuver_lag_m(
    sigma_w: float,
    alpha: float,
    beta: float,
    gamma_d: float,
) -> float:
    """Maximum position MSE under a sustained maneuver, POMR Eq. (19.60).

    MMSE_p = sigma_w^2 [ (2 alpha^2 + beta(2 - 3 alpha))/(alpha(4 - 2 alpha - beta))
                         + (1 - alpha)^2 Gamma_D^2 / beta^2 ]

    Returned as an RMS distance. The second term is the deterministic lag: a
    filter tuned for a quiet target falls progressively behind a maneuvering
    one, and no amount of SNR fixes it.
    """
    if beta <= 0:
        raise ValueError("beta must be > 0")
    denom = alpha * (4.0 - 2.0 * alpha - beta)
    if denom <= 0:
        raise ValueError("alpha, beta outside the stable region")
    noise_term = (2.0 * alpha**2 + beta * (2.0 - 3.0 * alpha)) / denom
    lag_term = (1.0 - alpha) ** 2 * gamma_d**2 / beta**2
    return float(sigma_w * math.sqrt(noise_term + lag_term))

See Also

• Theory: Radar Equation
• Theory: MTI Clutter Suppression
• Theory: Track Accuracy
• User Guide: Radar Detection
• Scenarios API