From be55196be8e694f85b863ca7ebb9602c8ef84e2f Mon Sep 17 00:00:00 2001 From: ApsisRS contributors Date: Sat, 26 Sep 2026 15:29:41 +0800 Subject: [PATCH] workspace: add the S.465 cross-check and document the ITU-R suite --- crates/rf/apsis-antenna/src/reference.rs | 7 +- validation/README.md | 6 ++ validation/python/link/check_s465.py | 130 +++++++++++++++++++++++ 3 files changed, 142 insertions(+), 1 deletion(-) create mode 100644 validation/python/link/check_s465.py diff --git a/crates/rf/apsis-antenna/src/reference.rs b/crates/rf/apsis-antenna/src/reference.rs index 0919221..1238a0b 100644 --- a/crates/rf/apsis-antenna/src/reference.rs +++ b/crates/rf/apsis-antenna/src/reference.rs @@ -115,6 +115,12 @@ pub fn main_lobe_break_angle(diameter_over_wavelength: f64) -> Option { /// head is 29 for a large aperture and 32 for a small one, and the floor /// angle is 36 or 48 degrees accordingly. /// +/// The falling sidelobe expression does not reach exactly minus ten decibels +/// at the floor angle, so the pattern has a step of a few hundredths of a +/// decibel there. This is a property of the recommendation, which defines +/// the regions separately, and it is reproduced here rather than smoothed +/// over. +/// /// Returns `None` when the ratio or the angle is not usable. pub fn reference_pattern_gain(diameter_over_wavelength: f64, off_axis_angle: f64) -> Option { if !off_axis_angle.is_finite() || off_axis_angle < 0.0 { @@ -133,7 +139,6 @@ pub fn reference_pattern_gain(diameter_over_wavelength: f64, off_axis_angle: f64 Some(-10.0) } } - /// Returns the gain of the S.465 reference pattern as a decibel value. /// /// Returns `None` when the ratio or the angle is not usable. diff --git a/validation/README.md b/validation/README.md index a605431..a1938cd 100644 --- a/validation/README.md +++ b/validation/README.md @@ -70,6 +70,10 @@ cargo run -p apsis-rs --features rf --example dump_rain \ # ITU-R P.618-13 雨衰斜路径 cargo run -p apsis-rs --features rf --example dump_rain_path \ > validation/data/expected/rust_rain_path.csv + +# ITU-R S.465 参考天线方向图 +cargo run -p apsis-rs --features rf --example dump_s465 \ + > validation/data/expected/rust_s465.csv ``` Python 侧读取同一批输入，用参考实现重算并输出误差表： @@ -82,6 +86,7 @@ pixi run python validation/python/link/check_link.py pixi run python validation/python/atmosphere/check_gas.py pixi run python validation/python/atmosphere/check_rain.py pixi run python validation/python/atmosphere/check_rain_path.py +pixi run python validation/python/link/check_s465.py ``` 两个脚本都以退出码 0 表示通过、1 表示超差、2 表示缺少输入文件。 @@ -110,6 +115,7 @@ Earth orientation 数据自动下载；本 crate 涉及的换算不需要该数 | `atmosphere/check_gas.py` | 气体吸收，18 个频率，涵盖 60 GHz 氧气复体、22/183 GHz 水汽线 | `itur` 的 ITU-R P.676-12 实现 | 干空气 6.1e-11，水汽 6.4e-9（相对） | | `atmosphere/check_rain.py` | 雨衰系数，17 个频率 × 4 组仰角/极化 | `itur` 的 ITU-R P.838-3 实现 | k 7.5e-9，alpha 6.5e-13（相对） | | `atmosphere/check_rain_path.py` | 雨衰斜路径，8 个场景 | ITU-R P.618-13 方程（Python 转写，系数来自 `itur`） | 3.4e-13 km / 4.3e-13 dB | +| `link/check_s465.py` | 参考天线方向图，3 种口径 × 8 个离轴角 | ITU-R S.465 方程（Python 转写） | 4.0e-10 dB | ### 关于射频层对拍的说明 diff --git a/validation/python/link/check_s465.py b/validation/python/link/check_s465.py new file mode 100644 index 0000000..a8f8e37 --- /dev/null +++ b/validation/python/link/check_s465.py @@ -0,0 +1,130 @@ +"""Cross-check the ApsisRS ITU-R S.465 reference pattern. + +Reads the values written by ``cargo run -p apsis-rs --features rf --example +dump_s465`` and recomputes each one from the equations of the recommendation. + +The script needs no network access. + +Usage: + + pixi run python validation/python/link/check_s465.py [path] + +Exit status is zero when every value agrees within tolerance and one +otherwise. +""" + +from __future__ import annotations + +import csv +import math +import sys +from pathlib import Path + +# The peak gain offset and the branch constants of the recommendation. +PEAK_GAIN_OFFSET_DB = 7.7 +LARGE_APERTURE_RATIO = 100.0 + +# Tolerance in decibels. +TOLERANCE_DB = 1e-6 + +DEFAULT_INPUT = Path("validation/data/expected/rust_s465.csv") + + +def read_rows(path: Path) -> list[dict[str, str]]: + """Reads the CSV rows.""" + with path.open(newline="") as handle: + reader = csv.DictReader(handle) + return [row for row in reader if row.get("ratio")] + + +def break_angle(ratio: float) -> tuple[float, float]: + """Returns the branch constant and the break angle, in degrees.""" + head, floor_angle = (29.0, 36.0) if ratio >= LARGE_APERTURE_RATIO else (32.0, 48.0) + peak = 20.0 * math.log10(ratio) + PEAK_GAIN_OFFSET_DB + + def difference(angle_degrees: float) -> float: + main_lobe = peak - 2.5e-3 * (ratio * angle_degrees) ** 2 + sidelobe = head - 25.0 * math.log10(angle_degrees) + return main_lobe - sidelobe + + low, high = 1.0e-9, floor_angle + if difference(high) > 0.0: + return head, high + for _ in range(200): + middle = 0.5 * (low + high) + if difference(middle) > 0.0: + low = middle + else: + high = middle + return head, 0.5 * (low + high) + + +def reference(ratio: float, angle_degrees: float) -> float: + """Recomputes the pattern gain in decibels. + + The API takes the angle in radians, so the degrees recorded in the file + are converted to radians and back before the branch decision. That + reproduces the rounding the implementation performs, which matters at + the branch boundaries where the recommendation has a small step. + """ + peak = 20.0 * math.log10(ratio) + PEAK_GAIN_OFFSET_DB + head, break_at = break_angle(ratio) + head, floor_angle = ( + (29.0, 36.0) if ratio >= LARGE_APERTURE_RATIO else (32.0, 48.0) + ) + angle_degrees = math.degrees(math.radians(angle_degrees)) + if angle_degrees < break_at: + return peak - 2.5e-3 * (ratio * angle_degrees) ** 2 + if angle_degrees <= floor_angle: + return head - 25.0 * math.log10(angle_degrees) + return -10.0 + + +def main(argv: list[str]) -> int: + path = Path(argv[1]) if len(argv) > 1 else DEFAULT_INPUT + if not path.exists(): + print(f"missing input file: {path}", file=sys.stderr) + print( + "generate it with: cargo run -p apsis-rs --features rf " + "--example dump_s465 > " + str(path), + file=sys.stderr, + ) + return 2 + + rows = read_rows(path) + if not rows: + print("input file has no data rows", file=sys.stderr) + return 2 + + worst_gain = 0.0 + worst_break = 0.0 + worst_at = None + + for row in rows: + ratio = float(row["ratio"]) + angle_degrees = float(row["angle_deg"]) + _, expected_break = break_angle(ratio) + expected_gain = reference(ratio, angle_degrees) + + break_error = abs(float(row["break_deg"]) - expected_break) + gain_error = abs(float(row["gain_dbi"]) - expected_gain) + if break_error > worst_break or gain_error > worst_gain: + worst_at = (ratio, angle_degrees) + worst_break = max(worst_break, break_error) + worst_gain = max(worst_gain, gain_error) + + print(f"points compared: {len(rows)}") + print(f"worst break: {worst_break:.3e} deg") + print(f"worst gain: {worst_gain:.3e} dB") + print(f"worst point: {worst_at}") + print(f"tolerance: {TOLERANCE_DB:.3e} dB") + + if worst_break <= TOLERANCE_DB and worst_gain <= TOLERANCE_DB: + print("result: pass") + return 0 + print("result: FAIL") + return 1 + + +if __name__ == "__main__": + raise SystemExit(main(sys.argv))