From c40f23c20eb64f68069fe5461d5076612ab8c262 Mon Sep 17 00:00:00 2001 From: jesse-ux Date: Mon, 5 Oct 2026 00:58:00 +0800 Subject: [PATCH] research(rectification): nadi second-level study helpers and N0 checks (BUG-1240) Five-level Vimshottari and equal-division D150 helpers, plus a runner that refuses to score until the preregistration file exists. First three lords match the production chain on the v5 events. No production files changed. --- docs/BUG_HISTORY.md | 16 + docs/tasks/README.md | 2 +- scripts/research/nadi_seconds_lib.py | 569 ++++++++++++++++++++++ scripts/research/nadi_seconds_run.py | 676 +++++++++++++++++++++++++++ tests/test_nadi_seconds_research.py | 179 +++++++ 5 files changed, 1441 insertions(+), 1 deletion(-) create mode 100644 scripts/research/nadi_seconds_lib.py create mode 100644 scripts/research/nadi_seconds_run.py create mode 100644 tests/test_nadi_seconds_research.py diff --git a/docs/BUG_HISTORY.md b/docs/BUG_HISTORY.md index afc1b700..9822d5c5 100644 --- a/docs/BUG_HISTORY.md +++ b/docs/BUG_HISTORY.md @@ -16657,3 +16657,19 @@ - 相关记录：10-04 咨询耗时审计（任务书事故实证第 2 条）。 - 复发自：无 - 修复版本：待发布 + +## BUG-1240 | 秒级 / 纳迪校准从未经过检验 + +- 状态：investigating +- 首次发现：2026-10-05 +- 最近更新：2026-10-05 +- 影响面：生时校正对外精度口径；`scripts/active_rectification_event_engine.py` 的 Vimshottari 计分只使用前三层。Sookshma、Prana 与 D150 不在评价集内。 +- 用户现象：校正结果不能被说成秒级。记录时间多数取整到 5 分钟，问前事对上某一天也不能单独证明某一秒。 +- 触发条件：无用户路径。离线研究，任务书 `docs/tasks/TASK-rectification-nadi-seconds-research-20261005.md`。 +- 根因：未知。已知的是主链不评第 4、第 5 层，且公开评价集的记录时钟没有秒。纳迪式规则能否胜过安慰剂日期，要等预先登记的检验跑完才能写。 +- 修复：尚未改生产代码。研究代码在 `scripts/research/nadi_seconds_lib.py` 与 `scripts/research/nadi_seconds_run.py`。 +- 验证：N0 对账测试已通过（前三层与主链 0 差）。N1–N5 尚未在 v5 上跑。 +- 防复发：规则、阈值和种子先写入预登记再跑；看过结果后改的规则不能用来判过门。不得重调 BUG-1091 已关闭的权重。 +- 相关记录：BUG-560、BUG-1090、BUG-1091、BUG-1105。 +- 复发自：无 +- 修复版本：研究分支 `codex/rectification-nadi-seconds-research-20261005`（未推送、未部署） diff --git a/docs/tasks/README.md b/docs/tasks/README.md index 776e38b8..e2ad943b 100644 --- a/docs/tasks/README.md +++ b/docs/tasks/README.md @@ -422,4 +422,4 @@ | `TASK-astrologer-rulings-batch6-20261004.md` | `PROGRESS-astrologer-rulings-batch6-20261004.md` | 第七轮裁定（共享仓书面回复，产品采用；问 1 选 A）：Rath 版双主星按 p.43 (a)–(e)（BUG-1227 解除 blocked，推翻第五批红线 2 的 Table 17 年数底线）；第 5 级宫主度数只倒算计都（p.71 脚注 42）；罗计旺陷 ±1 年；同宫两主比经度；BUG 从 1228 起 | 待验收 | 分支 `codex/astrologer-rulings-batch6-20261004`（BUG-1227～1229；校正分数文件未改、未升版本） | | `TASK-gate-log-volume-20261004.md` | `PROGRESS-gate-log-volume-20261004.md` | 门禁 validate 单步日志 4.4 万行 / 2 MB 网页打不开（run 3170）：快速门只打摘要（失败给末 200 行 + 日志文件）、前端测试门禁上 dot + 失败汇总（本机仍 TAP）、拆分 validate 为 7 个 step（产品授权改 workflow，只限拆分与重定向）；检查一项不少 | **已实现待验收**（BUG-1230；未推送、未部署；Gitea 各 step 页面是否打得开留待推送后由产品确认） | 分支 `codex/gate-log-volume-20261004` | | `TASK-consult-latency-quickwins-20261005.md` | `PROGRESS-consult-latency-quickwins-20261005.md` | 普通对话耗时两项快修：咨询链校正闸不再同步白等 VedAstro 官方快照子进程（每域约 4 s，复用顶层缓存 + 负缓存，不推翻 BUG-301）；补分段计时（第 0/1 步耗时、推理 token、分类耗时进日志与 usage）。推理强度/精简说明待模型 key 另单（BUG-1231、BUG-1232） | 已验收（Claude 10-05 Linux：Python 62 / 前端 24 失败与基线逐名相同，Static、gzip 0%；装 SDK 实测每领域 4.47→1.66 s，剩余 1.33 s 是保留的顶层前台等待；/admin/usage 真机欠） | 分支 `codex/consult-latency-quickwins-20261005` | -| `TASK-rectification-nadi-seconds-research-20261005.md` | `PROGRESS-rectification-nadi-seconds-research-20261005.md` | **「纳迪秒级校准」可证伪检验（离线）**：竞品宣传「问前事到天 → 秒级」。本仓主链只用三层小运，Sookshma / Prana 与 D150 从未进评价集；v5 真值 52/77 是整 5 分钟（秒级无真值可对）。N0 五层小运 + D150 底座（前三层与主链对账 0 差）、N1 拟合率 vs 安慰剂日期（核心）、N2 留一件预测、N3 六题后区间内再细分能否提头名、N4 岁差 / 坐标 / 时间扰动的噪声地板、N5 D150 结构层（原文比对 blocked）、N6 结论 + 对外口径草稿。规则先登记再跑；不改生产代码；不得重调 BUG-1091 已关的权重 | 待领取 | 分支 `codex/rectification-nadi-seconds-research-20261005`（BUG-1240） | +| `TASK-rectification-nadi-seconds-research-20261005.md` | `PROGRESS-rectification-nadi-seconds-research-20261005.md` | **「纳迪秒级校准」可证伪检验（离线）**：竞品宣传「问前事到天 → 秒级」。本仓主链只用三层小运，Sookshma / Prana 与 D150 从未进评价集；v5 真值 52/77 是整 5 分钟（秒级无真值可对）。N0 五层小运 + D150 底座（前三层与主链对账 0 差）、N1 拟合率 vs 安慰剂日期（核心）、N2 留一件预测、N3 六题后区间内再细分能否提头名、N4 岁差 / 坐标 / 时间扰动的噪声地板、N5 D150 结构层（原文比对 blocked）、N6 结论 + 对外口径草稿。规则先登记再跑；不改生产代码；不得重调 BUG-1091 已关的权重 | 执行中 | 分支 `codex/rectification-nadi-seconds-research-20261005`（BUG-1240；未推送、未部署） | diff --git a/scripts/research/nadi_seconds_lib.py b/scripts/research/nadi_seconds_lib.py new file mode 100644 index 00000000..62f3d641 --- /dev/null +++ b/scripts/research/nadi_seconds_lib.py @@ -0,0 +1,569 @@ +#!/usr/bin/env python3 +"""Offline falsification helpers for second-scale nadi / sookshma / prana claims. + +Research only. This module imports production chart and Vimshottari functions +and does not change them. Fitters and rankers never read ``true_minute``. + +The five-level lords are the production chain +(``dasha_analyzer.build_antardasha`` applied twice more past pratyantar). +That is the same recursion ``_active_vimshottari`` uses for the first three +levels, including its midnight date anchor and 365.25-day year. +``calculate_five_level_dasha`` stays available as a side check; it is not the +scoring chain, because its 365.25636-day year does not reconcile to zero +against the production lords. + +Classical D150 order (movable direct, fixed reverse, dual from the middle) +is marked ``variant_unverified``. No Chandra Kala Nadi verse text is stored. +""" + +from __future__ import annotations + +import hashlib +import json +import math +from collections.abc import Mapping, Sequence +from datetime import date, datetime, timedelta +from pathlib import Path +from random import Random +from typing import Any + +import dasha_analyzer +import divisional_charts_extended +import domain_calculation_service +import narayana_dasha +from scripts.active_rectification_event_engine import ( + AYANAMSA, + DOMAIN_CONFIG, + NODE_MODE, + _active_vimshottari, + _event_datetime, + _house_lords, +) + +ROOT = Path(__file__).resolve().parents[2] +HOLDOUT_V5 = ROOT / "references" / "real_case_calibration" / "minute_rectification_holdout_v5.json" +PREREG_PATH = ROOT / "docs" / "research" / "nadi_seconds_preregistration_2026_10_05.json" +EVENT_CUTOFF = date(2026, 9, 14) +LEVELS = ("md", "ad", "pd", "sookshma", "prana") +# Structural reading only. See classical_nadi_index. +CLASSICAL_NADI_STATUS = "variant_unverified" +KM_PER_DEGREE_LAT = 111.32 + +_VARGA = divisional_charts_extended.DivisionalChartsCalculator() +_SKY_CACHE: dict[tuple, dict[str, Any]] = {} + + +def load_cases(path: Path = HOLDOUT_V5) -> list[dict[str, Any]]: + return list(json.loads(path.read_text(encoding="utf-8"))["cases"]) + + +def canonical_bytes(payload: Any) -> bytes: + return ( + json.dumps(payload, ensure_ascii=False, sort_keys=True, separators=(",", ":")) + + "\n" + ).encode("utf-8") + + +def sha256_bytes(raw: bytes) -> str: + return hashlib.sha256(raw).hexdigest() + + +def file_sha256(path: Path) -> str: + return hashlib.sha256(path.read_bytes()).hexdigest() + + +def round6(value: float) -> float: + return round(float(value), 6) + + +def subdivide(period: Mapping[str, Any]) -> list[dict[str, Any]]: + """One production antardasha split. Adjacent rows abut by construction.""" + return list(dasha_analyzer.build_antardasha(dict(period))) + + +def chain_at(birth_date: str, moon_longitude: float, event_at: datetime) -> list[dict[str, Any]]: + """Five nested production periods active at ``event_at``. + + ``birth_date`` is the candidate's local calendar date (YYYY-MM-DD), the + same anchor ``_active_vimshottari`` passes through. + """ + return _chain_from_timeline(_timeline(birth_date, moon_longitude), event_at) + + +def five_lords(birth_date: str, moon_longitude: float, event_at: datetime) -> tuple[str, str, str, str, str]: + return tuple(str(period["lord"]) for period in chain_at(birth_date, moon_longitude, event_at)) # type: ignore[return-value] + + +def production_three(birth_date: str, moon_longitude: float, event_at: datetime) -> tuple[str, str, str]: + return _active_vimshottari(birth_date, float(moon_longitude) % 360.0, event_at) + + +def d150_equal(longitude: float) -> dict[str, Any]: + """Equal-slice D150 via the repository ``calc_custom_varga`` (N=150).""" + row = _VARGA.calc_custom_varga(float(longitude) % 360.0, 150) + sign = str(row["sign"]) + return { + "part_index": int(row["part_index"]), + "sign": sign, + "sign_idx": int(row["sign_idx"]), + "sign_lord": narayana_dasha.SIGN_LORDS[sign], + "source": "calc_custom_varga", + } + + +def classical_nadi_index(longitude: float) -> dict[str, Any]: + """Reorder the 150 equal slices. Not a gate input. + + Movable / fixed / dual follow ``DivisionalChartsCalculator`` sign classes. + Movable: part 0..149 in longitude order. Fixed: 149 down to 0. + Dual: start at the middle index 75 and wrap. The dual start is a + structural reading of "从中段起排", not a quotation, so the status is + ``variant_unverified``. Named nadi verses are not used. + """ + longitude = float(longitude) % 360.0 + sign_index = int(longitude // 30.0) % 12 + part = int(d150_equal(longitude)["part_index"]) + if sign_index in _VARGA.MOVABLE_SIGNS: + index = part + order = "direct" + elif sign_index in _VARGA.FIXED_SIGNS: + index = 149 - part + order = "reverse" + else: + index = (part + 75) % 150 + order = "from_middle" + return { + "index": index, + "equal_part": part, + "sign_index": sign_index, + "order": order, + "status": CLASSICAL_NADI_STATUS, + } + + +def domain_target_lords(ascendant_index: int, domain: str) -> frozenset[str]: + houses = DOMAIN_CONFIG[domain][1] + return frozenset(_house_lords(int(ascendant_index), houses)) + + +def explained(lords: Sequence[str], targets: frozenset[str], rule: str, segment_lord: str | None = None) -> bool: + """Binary 'this second explains this event' for one pre-registered rule. + + A4 / A5: level 4 or 5 lord is a domain house lord. + B1..B5: at least that many of the five levels are domain house lords. + N5A: equal-division D150 lagna sign lord is a domain house lord. + N5B: that same sign lord appears in any of the five levels. + """ + lords = tuple(lords) + if rule == "A4": + return lords[3] in targets + if rule == "A5": + return lords[4] in targets + if rule.startswith("B") and rule[1:].isdigit(): + need = int(rule[1:]) + return sum(lord in targets for lord in lords) >= need + if rule == "N5A": + return segment_lord is not None and segment_lord in targets + if rule == "N5B": + return segment_lord is not None and segment_lord in set(lords) + raise ValueError(f"unknown rule {rule}") + + +def fit_count( + lord_rows: Sequence[Sequence[str]], + domains: Sequence[str], + ascendant_index: int, + rule: str, + segment_lord: str | None = None, +) -> int: + """How many events ``rule`` explains. No clock and no truth label.""" + total = 0 + for lords, domain in zip(lord_rows, domains, strict=True): + targets = domain_target_lords(ascendant_index, domain) + if explained(lords, targets, rule, segment_lord): + total += 1 + return total + + +def rank_offsets(counts: Mapping[int, int]) -> list[int]: + """Offsets tied for the highest count, ascending. Ties are all kept.""" + if not counts: + return [] + best = max(counts.values()) + return sorted(offset for offset, count in counts.items() if count == best) + + +def candidate_moment(birth: Mapping[str, Any], offset_seconds: int) -> datetime: + base = datetime.strptime(f"{birth['date']} {str(birth['time'])[:5]}", "%Y-%m-%d %H:%M") + return base + timedelta(seconds=int(offset_seconds)) + + +def sky_at( + birth: Mapping[str, Any], + moment: datetime, + *, + ayanamsa: str = AYANAMSA, + node_mode: str = NODE_MODE, + longitude: float | None = None, +) -> dict[str, Any]: + """Moon and ascendant from the production chart. Seconds are passed through.""" + lat = round(float(birth["latitude"]), 6) + lon = round(float(birth["longitude"] if longitude is None else longitude), 6) + tz = round(float(birth["timezone_offset"]), 4) + key = (moment.year, moment.month, moment.day, moment.hour, moment.minute, moment.second, lat, lon, tz, str(ayanamsa), str(node_mode)) + cached = _SKY_CACHE.get(key) + if cached is not None: + return cached + chart = domain_calculation_service.compute_chart({ + "year": moment.year, + "month": moment.month, + "day": moment.day, + "hour": moment.hour, + "minute": moment.minute, + "second": moment.second, + "lat": lat, + "lon": lon, + "tz": tz, + "ayanamsa": ayanamsa, + "node_mode": node_mode, + }) + moon = float(chart["planets"]["Moon"]["lon"]) % 360.0 + asc = float(chart["ascendant"]["lon"]) % 360.0 + rahu = chart.get("planets", {}).get("Rahu") or {} + rahu_lon = float(rahu["lon"]) % 360.0 if isinstance(rahu.get("lon"), (int, float)) else None + sky = { + "birth_date": moment.date().isoformat(), + "moon_lon": moon, + "asc_lon": asc, + "asc_index": int(asc // 30.0) % 12, + "rahu_lon": rahu_lon, + } + _SKY_CACHE[key] = sky + return sky + + +def shift_longitude_km(longitude: float, latitude: float, kilometers_east: float) -> float: + """Move the birth longitude by an east-west ground distance. + + One degree of longitude is ``111.32 * cos(latitude)`` kilometres. + Near a pole the east-west degree collapses and the shift is refused. + """ + scale = KM_PER_DEGREE_LAT * math.cos(math.radians(float(latitude))) + if abs(scale) < 1e-3: + raise ValueError("east-west shift is undefined this close to a pole") + moved = float(longitude) + float(kilometers_east) / scale + return ((moved + 180.0) % 360.0) - 180.0 + + +def day_events(case: Mapping[str, Any]) -> list[dict[str, Any]]: + """Day-precision events copied without truth labels or narrative text.""" + rows = [] + for event in case["events"]: + if event.get("precision") != "day": + continue + rows.append({ + "id": str(event["id"]), + "domain": str(event["domain"]), + "date": str(event["date"])[:10], + "precision": "day", + }) + return rows + + +def _clamp_life_date(day: date, birth: date, cutoff: date) -> date: + guard = 0 + while day <= birth and guard < 160: + day = day + timedelta(days=365) + guard += 1 + guard = 0 + while day > cutoff and guard < 160: + day = day - timedelta(days=365) + guard += 1 + if birth < day <= cutoff: + return day + fallback = birth + timedelta(days=400) + if fallback > cutoff: + fallback = cutoff + if fallback <= birth: + fallback = birth + timedelta(days=1) + return fallback + + +def placebo_shift_dates( + events: Sequence[Mapping[str, Any]], + birth_date: str, + *, + seed: str, + low_days: int, + high_days: int, + cutoff: date = EVENT_CUTOFF, +) -> list[dict[str, Any]]: + """Shift each day-event date by a seeded ±[low, high] day offset. + + Domain stays. A result outside (birth, cutoff] is flipped, then repaired + by 365-day steps. The seed includes the case id so case order cannot + change one case's offsets. + """ + rng = Random(seed) + birth = date.fromisoformat(birth_date) + shifted = [] + for event in events: + original = date.fromisoformat(str(event["date"])[:10]) + magnitude = rng.randint(int(low_days), int(high_days)) + sign = rng.choice((-1, 1)) + candidate = original + timedelta(days=sign * magnitude) + if not (birth < candidate <= cutoff): + candidate = original + timedelta(days=-sign * magnitude) + if not (birth < candidate <= cutoff): + candidate = _clamp_life_date(original + timedelta(days=magnitude), birth, cutoff) + shifted.append({**event, "date": candidate.isoformat()}) + return shifted + + +def placebo_swap_dates( + cases: Sequence[Mapping[str, Any]], + *, + seed: int, + cutoff: date = EVENT_CUTOFF, +) -> dict[str, list[dict[str, Any]]]: + """Permute day-event dates across cases. Domains stay on the recipient event.""" + slots: list[tuple[str, dict[str, Any]]] = [] + for case in cases: + for event in day_events(case): + slots.append((str(case["case_id"]), event)) + slots.sort(key=lambda item: (item[0], item[1]["id"])) + dates = [item[1]["date"] for item in slots] + order = list(range(len(dates))) + Random(int(seed)).shuffle(order) + by_case: dict[str, list[dict[str, Any]]] = {} + births = {str(case["case_id"]): str(case["birth"]["date"]) for case in cases} + for index, (case_id, event) in enumerate(slots): + raw = date.fromisoformat(dates[order[index]]) + clamped = _clamp_life_date(raw, date.fromisoformat(births[case_id]), cutoff) + by_case.setdefault(case_id, []).append({**event, "date": clamped.isoformat()}) + return by_case + + +def event_moments(events: Sequence[Mapping[str, Any]]) -> list[datetime]: + return [_event_datetime(event) for event in events] # type: ignore[arg-type] + + +def _timeline(birth_date: str, moon_longitude: float) -> list[dict[str, Any]]: + nakshatra, progress, _pada = dasha_analyzer.lon_to_nakshatra(float(moon_longitude) % 360.0) + timeline, _elapsed, _remaining, _lord = dasha_analyzer.build_dasha_timeline( + birth_date, nakshatra, progress, + ) + return timeline + + +def _chain_from_timeline(timeline: Sequence[Mapping[str, Any]], event_at: datetime) -> list[dict[str, Any]]: + _index, major = dasha_analyzer.find_current(list(timeline), event_at) + levels = [major] + for _depth in range(4): + levels.append(dasha_analyzer.find_current_sub(subdivide(levels[-1]), event_at)) + return levels + + +def lords_for_events( + birth_date: str, + moon_longitude: float, + events: Sequence[Mapping[str, Any]], +) -> list[tuple[str, str, str, str, str]]: + """Five lords for every event. The mahadasha timeline is built once.""" + timeline = _timeline(birth_date, moon_longitude) + rows = [] + for moment in event_moments(events): + rows.append(tuple(str(period["lord"]) for period in _chain_from_timeline(timeline, moment))) # type: ignore[arg-type] + return rows + + +def counts_by_offset( + skies: Mapping[int, Mapping[str, Any]], + events: Sequence[Mapping[str, Any]], + rule: str, + *, + use_segment_lord: bool, +) -> dict[int, int]: + domains = [str(event["domain"]) for event in events] + counts: dict[int, int] = {} + for offset, sky in skies.items(): + rows = lords_for_events(str(sky["birth_date"]), float(sky["moon_lon"]), events) + segment_lord = str(d150_equal(float(sky["asc_lon"]))["sign_lord"]) if use_segment_lord else None + counts[int(offset)] = fit_count(rows, domains, int(sky["asc_index"]), rule, segment_lord) + return counts + + +def mean_rate(counts: Mapping[int, int], offsets: Sequence[int], event_count: int) -> float | None: + if event_count <= 0 or not offsets: + return None + total = sum(counts[offset] for offset in offsets) + return total / (len(list(offsets)) * event_count) + + +def fraction_explaining_all(counts: Mapping[int, int], offsets: Sequence[int], event_count: int) -> float | None: + if event_count <= 0 or not offsets: + return None + hits = sum(1 for offset in offsets if counts[offset] == event_count) + return hits / len(list(offsets)) + + +def paired_sign_flip_p( + differences: Sequence[float], + *, + seed: int, + permutations: int, +) -> dict[str, Any]: + """One-sided paired permutation: share of sign-flips with mean >= observed. + + p = (count + 1) / (permutations + 1), counting the observed assignment. + """ + values = [float(item) for item in differences] + n = len(values) + if n == 0: + return {"n": 0, "mean": None, "p": None, "extreme": None, "permutations": permutations} + observed = sum(values) / n + rng = Random(int(seed)) + extreme = 0 + for _ in range(int(permutations)): + total = 0.0 + for value in values: + total += value if rng.randrange(2) == 0 else -value + if total / n >= observed - 1e-15: + extreme += 1 + return { + "n": n, + "mean": round6(observed), + "p": round6((extreme + 1) / (int(permutations) + 1)), + "extreme": extreme, + "permutations": int(permutations), + } + + +def bootstrap_mean_ci( + values: Sequence[float], + *, + seed: int, + resamples: int, +) -> dict[str, Any]: + """Percentile 95% interval of the mean. Passes only when the lower bound is > 0.""" + pool = [float(item) for item in values] + n = len(pool) + if n == 0: + return {"n": 0, "mean": None, "low": None, "high": None, "excludes_zero_positive": False} + rng = Random(int(seed)) + means = [] + for _ in range(int(resamples)): + draw = 0.0 + for _item in range(n): + draw += pool[rng.randrange(n)] + means.append(draw / n) + means.sort() + low_index = int(0.025 * (len(means) - 1)) + high_index = int(0.975 * (len(means) - 1)) + low = means[low_index] + high = means[high_index] + return { + "n": n, + "mean": round6(sum(pool) / n), + "low": round6(low), + "high": round6(high), + "excludes_zero_positive": bool(low > 0.0), + } + + +def minute_of(clock: str) -> int: + return int(str(clock)[3:5]) + + +def rounded_to_five_minutes(clock: str) -> bool: + return minute_of(clock) % 5 == 0 + + +def clock_delta_minutes(stamp: str, recorded: str) -> int: + def _clock(value: str) -> int: + hours, minutes = str(value)[:5].split(":") + return int(hours) * 60 + int(minutes) + + delta = (_clock(stamp) - _clock(recorded)) % 1440 + return delta - 1440 if delta > 720 else delta + + +def leaders_within(offsets: Sequence[int], limit_seconds: int) -> bool: + """True only when every tied leader is inside the limit. An empty tie is a miss.""" + return bool(offsets) and all(abs(int(offset)) <= int(limit_seconds) for offset in offsets) + + +def truth_audit(cases: Sequence[Mapping[str, Any]]) -> dict[str, Any]: + """Evaluator-side audit. This function is allowed to read the recorded clock. + + It does not read ``true_minute``; the recorded minute is ``birth.time``. + """ + rounded: list[str] = [] + unrounded: list[str] = [] + minute_hist: dict[str, int] = {} + per_case = [] + for case in cases: + clock = str(case["birth"]["time"])[:5] + minute = f"{minute_of(clock):02d}" + minute_hist[minute] = minute_hist.get(minute, 0) + 1 + case_id = str(case["case_id"]) + (rounded if rounded_to_five_minutes(clock) else unrounded).append(case_id) + day_count = sum(1 for event in case["events"] if event.get("precision") == "day") + per_case.append({"case_id": case_id, "minute": minute, "day_events": day_count}) + day_counts = [row["day_events"] for row in per_case] + return { + "cases": len(per_case), + "rounded_5min_count": len(rounded), + "unrounded_count": len(unrounded), + "rounded_5min_case_ids": rounded, + "unrounded_case_ids": unrounded, + "minute_histogram": minute_hist, + "day_event_total": sum(day_counts), + "day_events_ge_3": sum(1 for count in day_counts if count >= 3), + "day_events_ge_4": sum(1 for count in day_counts if count >= 4), + "day_events_ge_5": sum(1 for count in day_counts if count >= 5), + "day_events_zero": sum(1 for count in day_counts if count == 0), + "per_case": per_case, + } + + +def reconcile_case(case: Mapping[str, Any]) -> dict[str, Any]: + """First three lords versus ``_active_vimshottari`` for every event at the recorded second.""" + moment = candidate_moment(case["birth"], 0) + sky = sky_at(case["birth"], moment) + mismatches = [] + for event in case["events"]: + event_at = _event_datetime(event) # type: ignore[arg-type] + ours = five_lords(sky["birth_date"], sky["moon_lon"], event_at) + theirs = production_three(sky["birth_date"], sky["moon_lon"], event_at) + if ours[:3] != theirs: + mismatches.append(str(event["id"])) + return { + "case_id": str(case["case_id"]), + "events": len(case["events"]), + "mismatches": mismatches, + } + + +def public_fields(case: Mapping[str, Any]) -> dict[str, Any]: + """Copy the birth clock, place, and events a fitter may see.""" + birth = case["birth"] + return { + "case_id": str(case["case_id"]), + "birth": { + "date": str(birth["date"]), + "time": str(birth["time"])[:5], + "latitude": float(birth["latitude"]), + "longitude": float(birth["longitude"]), + "timezone_offset": float(birth["timezone_offset"]), + }, + "events": [ + { + "id": str(event["id"]), + "domain": str(event["domain"]), + "date": str(event["date"]), + "precision": str(event["precision"]), + } + for event in case["events"] + ], + } diff --git a/scripts/research/nadi_seconds_run.py b/scripts/research/nadi_seconds_run.py new file mode 100644 index 00000000..fa77e89b --- /dev/null +++ b/scripts/research/nadi_seconds_run.py @@ -0,0 +1,676 @@ +#!/usr/bin/env python3 +"""Run the pre-registered nadi / sookshma / prana falsification study. + +Reads ``docs/research/nadi_seconds_preregistration_2026_10_05.json`` and +refuses to score without it. The output JSON has no timestamps, so two runs +with ``PYTHONHASHSEED=0`` are byte-identical. This script does not import +``true_minute``. +""" + +from __future__ import annotations + +import argparse +import json +import sys +import traceback +from datetime import date +from pathlib import Path +from random import Random +from typing import Any, Mapping, Sequence + +ROOT = Path(__file__).resolve().parents[2] +if str(ROOT) not in sys.path: + sys.path.insert(0, str(ROOT)) + +from scripts.research import nadi_seconds_lib as lib # noqa: E402 + +REPORT = ROOT / "docs" / "research" / "nadi_seconds_results_2026_10_05.json" + + +def load_prereg() -> dict[str, Any]: + if not lib.PREREG_PATH.is_file(): + raise SystemExit(f"missing preregistration: {lib.PREREG_PATH}") + return json.loads(lib.PREREG_PATH.read_text(encoding="utf-8")) + + +def offsets_of(radius: int, step: int) -> list[int]: + return list(range(-int(radius), int(radius) + 1, int(step))) + + +def truth_minute_offsets() -> list[int]: + return list(range(0, 60)) + + +def collect_skies(birth: Mapping[str, Any], offsets: Sequence[int], **kwargs: Any) -> dict[int, dict[str, Any]]: + return {int(offset): lib.sky_at(birth, lib.candidate_moment(birth, int(offset)), **kwargs) for offset in offsets} + + +def lord_bundle(skies: Mapping[int, Mapping[str, Any]], events: Sequence[Mapping[str, Any]]) -> dict[str, Any]: + lords: dict[int, list] = {} + asc: dict[int, int] = {} + segment: dict[int, str] = {} + part: dict[int, int] = {} + for offset, sky in skies.items(): + lords[int(offset)] = lib.lords_for_events(str(sky["birth_date"]), float(sky["moon_lon"]), events) + asc[int(offset)] = int(sky["asc_index"]) + row = lib.d150_equal(float(sky["asc_lon"])) + segment[int(offset)] = str(row["sign_lord"]) + part[int(offset)] = int(row["part_index"]) + return {"lords": lords, "asc": asc, "segment": segment, "part": part} + + +def rule_counts(bundle: Mapping[str, Any], events: Sequence[Mapping[str, Any]], rule: str) -> dict[int, int]: + domains = [str(event["domain"]) for event in events] + use_segment = rule.startswith("N5") + counts: dict[int, int] = {} + for offset, rows in bundle["lords"].items(): + targets = [lib.domain_target_lords(bundle["asc"][offset], domain) for domain in domains] + segment = bundle["segment"][offset] if use_segment else None + counts[offset] = sum( + lib.explained(lords, target, rule, segment) + for lords, target in zip(rows, targets, strict=True) + ) + return counts + + +def rule_bits(bundle: Mapping[str, Any], events: Sequence[Mapping[str, Any]], rule: str) -> dict[int, list[bool]]: + domains = [str(event["domain"]) for event in events] + use_segment = rule.startswith("N5") + bits: dict[int, list[bool]] = {} + for offset, rows in bundle["lords"].items(): + targets = [lib.domain_target_lords(bundle["asc"][offset], domain) for domain in domains] + segment = bundle["segment"][offset] if use_segment else None + bits[offset] = [ + lib.explained(lords, target, rule, segment) + for lords, target in zip(rows, targets, strict=True) + ] + return bits + + +def mean_defined(values: Sequence[float | None]) -> float | None: + nums = [float(value) for value in values if value is not None] + if not nums: + return None + return lib.round6(sum(nums) / len(nums)) + + +def rate(flags: Sequence[bool]) -> float | None: + if not flags: + return None + return lib.round6(sum(bool(flag) for flag in flags) / len(flags)) + + +def split_rates(rows: Sequence[Mapping[str, Any]], key: str, rounded: set[str]) -> dict[str, float | None]: + return { + "all": rate([bool(row[key]) for row in rows]), + "rounded_5min": rate([bool(row[key]) for row in rows if row["case_id"] in rounded]), + "unrounded": rate([bool(row[key]) for row in rows if row["case_id"] not in rounded]), + "n_all": len(rows), + "n_rounded_5min": sum(1 for row in rows if row["case_id"] in rounded), + "n_unrounded": sum(1 for row in rows if row["case_id"] not in rounded), + } + + +def delivery_offsets(start: str | None, end: str | None, recorded: str, radius_minutes: int) -> list[int]: + """One-second offsets covering each delivered clock minute, clipped to the search radius.""" + if not start or not end: + return [] + cap = int(radius_minutes) * 60 + start_min = lib.clock_delta_minutes(start, recorded) + end_min = lib.clock_delta_minutes(end, recorded) + + def clip(lo: int, hi: int) -> list[int]: + lo = max(lo, -cap) + hi = min(hi, cap) + if lo > hi: + return [] + return list(range(lo, hi + 1)) + + if start_min <= end_min: + return clip(start_min * 60, end_min * 60 + 59) + return sorted(set(clip(start_min * 60, cap) + clip(-cap, end_min * 60 + 59))) + + +def six_probe(case: Mapping[str, Any], radius: int) -> dict[str, Any]: + """futile_collect_stop_replay d1: six probes, no guided-window injection.""" + from scripts.active_rectification_event_engine import AYANAMSA, NODE_MODE, compute_candidate_static_contexts + from scripts.rectification.event_probes import discriminating_event_probes + from scripts.rectification.refinement_packet import window_scan + from scripts.rectification.scoring_service import build_event_contribution_matrix, score_from_matrix + from scripts.research.fewer_probes_card_replay import _outcome + from scripts.research.guided_collect_holdout_replay import _hhmm, posterior_state + from scripts.research.minute_resolution_sweep import MINUTE_STEP, scoring_request_for + from scripts.research.probe_supply_after_six import ASK_COUNT, TODAY + + recorded = str(case["birth"]["time"])[:5] + request = scoring_request_for(dict(case), radius) + request["ayanamsa"] = AYANAMSA + request["node_mode"] = NODE_MODE + request["minute_step"] = MINUTE_STEP + static_contexts = compute_candidate_static_contexts(request) + built = build_event_contribution_matrix(request, static_contexts=static_contexts) + rows = score_from_matrix(request, built) + times = [stamp for row in rows if (stamp := _hhmm(row.get("time")))] + probes = discriminating_event_probes( + {**request, "refresh_probes": False, "asked_probe_keys": []}, + built, + scan=window_scan(built), + candidate_times=times, + representative_time=recorded, + today=TODAY, + )[:ASK_COUNT] + state = posterior_state(rows=rows, contexts=static_contexts, probes=probes, true_time=recorded) + outcome = _outcome(state, recorded) + scores = state["scores"] + valid = list(state["valid"]) + if not valid: + leaders: list[str] = [] + else: + def _score(row: Mapping[str, Any]) -> float: + stamp = _hhmm(row.get("time")) or "" + return float(scores.get(stamp, row.get("score") or 0)) + + best = max(_score(row) for row in valid) + leaders = sorted({ + str(row.get("time"))[:5] + for row in valid + if abs(_score(row) - best) <= 1e-9 and str(row.get("time") or "")[:5] + }) + deltas = [lib.clock_delta_minutes(stamp, recorded) for stamp in leaders] + return { + "start": outcome["start"], + "end": outcome["end"], + "width": outcome["width"], + "truth_in_range": bool(outcome["truth_in_range"]), + "leader_times": leaders, + "all_within_1min": bool(deltas) and all(abs(delta) <= 1 for delta in deltas), + "any_exact": recorded in leaders, + "error": None, + } + + +def gate_p(result: Mapping[str, Any]) -> bool: + if result.get("mean") is None or result.get("extreme") is None: + return False + exact = (int(result["extreme"]) + 1) / (int(result["permutations"]) + 1) + return bool(result["mean"] > 0 and exact < 0.05) + + +def summarize_rule( + rows: Sequence[Mapping[str, Any]], + *, + seed: int, + permutations: int, +) -> dict[str, Any]: + diffs = [float(row["diff"]) for row in rows] + test = lib.paired_sign_flip_p(diffs, seed=seed, permutations=permutations) + return { + **test, + "pass": gate_p(test), + "truth_mean": mean_defined([row["truth"] for row in rows]), + "placebo_mean": mean_defined([row["placebo"] for row in rows]), + "window_mean": mean_defined([row["window"] for row in rows]), + "fraction_grid_explains_all_mean": mean_defined([row["fraction_all"] for row in rows]), + } + + +def n2_case( + bits: Mapping[int, Sequence[bool]], + coarse: Sequence[int], + events: Sequence[Mapping[str, Any]], + placebo_bits: Mapping[int, Sequence[bool]], + *, + case_id: str, + seed: str, +) -> dict[str, Any]: + deltas_random = [] + deltas_placebo = [] + held_rates = [] + for index, event in enumerate(events): + counts = { + offset: sum(flag for event_index, flag in enumerate(flags) if event_index != index) + for offset, flags in bits.items() + if offset in coarse + } + leaders = [offset for offset in lib.rank_offsets(counts) if offset in set(coarse)] + if not leaders: + continue + held = sum(1 for offset in leaders if bits[offset][index]) / len(leaders) + rng = Random(f"{seed}:{case_id}:{event['id']}") + if len(leaders) >= len(coarse): + sample = list(coarse) + else: + sample = rng.sample(list(coarse), len(leaders)) + random_rate = sum(1 for offset in sample if bits[offset][index]) / len(sample) + placebo_rate = sum(1 for offset in leaders if placebo_bits[offset][index]) / len(leaders) + held_rates.append(held) + deltas_random.append(held - random_rate) + deltas_placebo.append(held - placebo_rate) + return { + "held_rate": mean_defined(held_rates), + "delta_random": mean_defined(deltas_random), + "delta_placebo": mean_defined(deltas_placebo), + } + + +def level_change(base: Sequence[Sequence[str]], other: Sequence[Sequence[str]], index: int) -> float | None: + if not base or len(base) != len(other): + return None + changed = sum(left[index] != right[index] for left, right in zip(base, other, strict=True)) + return lib.round6(changed / len(base)) + + +def build(prereg: Mapping[str, Any], *, limit: int, radii: Sequence[int]) -> dict[str, Any]: + cases = [lib.public_fields(case) for case in lib.load_cases()] + if limit: + cases = cases[:limit] + grids = prereg["grids"] + coarse = offsets_of(grids["coarse_radius_seconds"], grids["coarse_step_seconds"]) + fine = offsets_of(grids["fine_radius_seconds"], grids["fine_step_seconds"]) + truth_offsets = truth_minute_offsets() + union = sorted(set(coarse) | set(fine) | set(truth_offsets)) + rules = list(prereg["rules"]["computed"]) + primary = list(prereg["rules"]["n1_primary"]) + n1_min = int(prereg["samples"]["n1_min_day_events"]) + n2_min = int(prereg["samples"]["n2_min_day_events"]) + permutations = int(prereg["tests"]["permutations"]) + resamples = int(prereg["tests"]["bootstrap_resamples"]) + shift_low = int(prereg["placebo"]["shift_low_days"]) + shift_high = int(prereg["placebo"]["shift_high_days"]) + swap = lib.placebo_swap_dates(cases, seed=int(prereg["placebo"]["swap_seed"])) + rounded = {case_id for case_id in lib.truth_audit(cases)["rounded_5min_case_ids"]} + audit = lib.truth_audit(cases) + recon = [lib.reconcile_case(case) for case in cases] + + n1_rows: dict[str, list[dict[str, Any]]] = {rule: [] for rule in primary} + n1_placebos = ("P1", "P2") + n5_rows: dict[str, list[dict[str, Any]]] = {rule: [] for rule in prereg["rules"]["n5_primary"]} + n2_rows = [] + n3a_rows = [] + n4_rows = [] + fraction_rows: dict[str, list[float]] = {rule: [] for rule in ("A4", "A5", "B1", "B2")} + case_fit: list[dict[str, Any]] = [] + + for index, case in enumerate(cases, start=1): + events = lib.day_events(case) + birth = case["birth"] + print(f"nadi {index}/{len(cases)} {case['case_id']} day_events={len(events)}", flush=True) + skies = collect_skies(birth, union) + real = lord_bundle(skies, events) if events else None + shifted = lib.placebo_shift_dates( + events, + birth["date"], + seed=prereg["placebo"]["shift_seed_template"].format(case_id=case["case_id"]), + low_days=shift_low, + high_days=shift_high, + ) if events else [] + p1 = lord_bundle(skies, shifted) if shifted else None + swapped = swap.get(case["case_id"], []) + p2 = lord_bundle(skies, swapped) if swapped else None + fit_row: dict[str, Any] = { + "case_id": case["case_id"], + "day_events": len(events), + "rounded_5min": case["case_id"] in rounded, + } + if real is not None: + for rule in rules: + if rule.startswith("N5") and rule not in prereg["rules"]["n5_primary"] and rule not in primary: + continue + counts = rule_counts(real, events, rule) + truth = lib.mean_rate(counts, truth_offsets, len(events)) + window = lib.mean_rate(counts, coarse, len(events)) + fraction = lib.fraction_explaining_all(counts, coarse, len(events)) + fit_row[rule] = { + "truth": None if truth is None else lib.round6(truth), + "window": None if window is None else lib.round6(window), + "fraction_all": None if fraction is None else lib.round6(fraction), + } + if rule in fraction_rows and fraction is not None and len(events) >= n1_min: + fraction_rows[rule].append(fraction) + if rule in primary and len(events) >= n1_min and p1 is not None and p2 is not None: + for label, bundle in (("P1", p1), ("P2", p2)): + placebo_rate = lib.mean_rate(rule_counts(bundle, shifted if label == "P1" else swapped, rule), truth_offsets, len(events)) + n1_rows[rule].append({ + "case_id": case["case_id"], + "placebo": label, + "truth": lib.round6(truth or 0.0), + "placebo_rate": None if placebo_rate is None else lib.round6(placebo_rate), + "diff": None if placebo_rate is None else lib.round6((truth or 0.0) - placebo_rate), + "window": None if window is None else lib.round6(window), + "fraction_all": None if fraction is None else lib.round6(fraction), + }) + if rule in n5_rows and len(events) >= n1_min and p1 is not None and p2 is not None: + for label, bundle in (("P1", p1), ("P2", p2)): + placebo_rate = lib.mean_rate(rule_counts(bundle, shifted if label == "P1" else swapped, rule), truth_offsets, len(events)) + n5_rows[rule].append({ + "case_id": case["case_id"], + "placebo": label, + "truth": lib.round6(truth or 0.0), + "placebo_rate": None if placebo_rate is None else lib.round6(placebo_rate), + "diff": None if placebo_rate is None else lib.round6((truth or 0.0) - placebo_rate), + "window": None if window is None else lib.round6(window), + "fraction_all": None if fraction is None else lib.round6(fraction), + }) + if len(events) >= n2_min and p1 is not None: + real_bits = rule_bits(real, events, prereg["rules"]["n2_primary"]) + p1_bits = rule_bits(p1, shifted, prereg["rules"]["n2_primary"]) + # P1 bits are aligned to shifted events, same index as real events. + n2_rows.append({ + "case_id": case["case_id"], + "rounded_5min": case["case_id"] in rounded, + **n2_case( + real_bits, coarse, events, p1_bits, + case_id=case["case_id"], seed=str(prereg["tests"]["n2_random_second_seed"]), + ), + }) + a5_counts = rule_counts(real, events, prereg["rules"]["n3_primary"]) if events else {} + leaders = lib.rank_offsets({offset: a5_counts.get(offset, 0) for offset in coarse}) if events else [] + n3a_rows.append({ + "case_id": case["case_id"], + "rounded_5min": case["case_id"] in rounded, + "day_events": len(events), + "hit_60s": lib.leaders_within(leaders, 60) if events else False, + "hit_120s": lib.leaders_within(leaders, 120) if events else False, + "tied_leaders": len(leaders), + }) + else: + n3a_rows.append({ + "case_id": case["case_id"], + "rounded_5min": case["case_id"] in rounded, + "day_events": 0, + "hit_60s": False, + "hit_120s": False, + "tied_leaders": len(coarse), + }) + case_fit.append(fit_row) + + # N4 at the recorded second only. + if events: + base_sky = skies[0] + base_lords = lib.lords_for_events(base_sky["birth_date"], base_sky["moon_lon"], events) + base_part = int(lib.d150_equal(base_sky["asc_lon"])["part_index"]) + perturbations: dict[str, Any] = {} + for name in prereg["n4"]["ayanamsas"]: + if name == prereg["ayanamsa"]: + continue + sky = lib.sky_at(birth, lib.candidate_moment(birth, 0), ayanamsa=name) + other = lib.lords_for_events(sky["birth_date"], sky["moon_lon"], events) + perturbations[f"ayanamsa:{name}"] = { + "level4": level_change(base_lords, other, 3), + "level5": level_change(base_lords, other, 4), + "d150_changed": int(lib.d150_equal(sky["asc_lon"])["part_index"]) != base_part, + } + for seconds in prereg["n4"]["time_offsets_seconds"]: + sky = skies.get(int(seconds)) or lib.sky_at(birth, lib.candidate_moment(birth, int(seconds))) + other = lib.lords_for_events(sky["birth_date"], sky["moon_lon"], events) + perturbations[f"time:{int(seconds)}"] = { + "level4": level_change(base_lords, other, 3), + "level5": level_change(base_lords, other, 4), + "d150_changed": int(lib.d150_equal(sky["asc_lon"])["part_index"]) != base_part, + } + for km in prereg["n4"]["longitude_shifts_km"]: + try: + moved = lib.shift_longitude_km(float(birth["longitude"]), float(birth["latitude"]), float(km)) + sky = lib.sky_at(birth, lib.candidate_moment(birth, 0), longitude=moved) + other = lib.lords_for_events(sky["birth_date"], sky["moon_lon"], events) + perturbations[f"east_km:{km}"] = { + "level4": level_change(base_lords, other, 3), + "level5": level_change(base_lords, other, 4), + "d150_changed": int(lib.d150_equal(sky["asc_lon"])["part_index"]) != base_part, + "blocked": False, + } + except ValueError: + perturbations[f"east_km:{km}"] = {"blocked": True} + node_sky = lib.sky_at(birth, lib.candidate_moment(birth, 0), node_mode="true") + perturbations["node:true"] = { + "moon_changed": abs(float(node_sky["moon_lon"]) - float(base_sky["moon_lon"])) > 1e-6, + "asc_changed": abs(float(node_sky["asc_lon"]) - float(base_sky["asc_lon"])) > 1e-6, + "rahu_changed": node_sky["rahu_lon"] != base_sky["rahu_lon"], + "level5": level_change(base_lords, lib.lords_for_events(node_sky["birth_date"], node_sky["moon_lon"], events), 4), + } + n4_rows.append({"case_id": case["case_id"], "day_events": len(events), "perturbations": perturbations}) + + def pack_comparisons(bucket: Mapping[str, list[dict[str, Any]]], seed_base: int) -> dict[str, Any]: + packed = {} + cursor = 0 + for rule, rows in bucket.items(): + packed[rule] = {} + for label in n1_placebos: + subset = [ + { + "diff": row["diff"], + "truth": row["truth"], + "placebo": row["placebo_rate"], + "window": row["window"], + "fraction_all": row["fraction_all"], + } + for row in rows + if row["placebo"] == label and row["diff"] is not None + ] + packed[rule][label] = summarize_rule(subset, seed=seed_base + cursor, permutations=permutations) + cursor += 1 + return packed + + n1 = pack_comparisons(n1_rows, int(prereg["tests"]["n1_permutation_seed_base"])) + n5 = pack_comparisons(n5_rows, int(prereg["tests"]["n5_permutation_seed_base"])) + n1_pass = all(item["pass"] for rule in n1.values() for item in rule.values()) + n5_pass = all(item["pass"] for rule in n5.values() for item in rule.values()) + n2_values = [float(row["delta_random"]) for row in n2_rows if row["delta_random"] is not None] + n2_placebo_values = [float(row["delta_placebo"]) for row in n2_rows if row["delta_placebo"] is not None] + n2_ci = lib.bootstrap_mean_ci(n2_values, seed=int(prereg["tests"]["n2_bootstrap_seed"]), resamples=resamples) + n2_placebo_ci = lib.bootstrap_mean_ci( + n2_placebo_values, seed=int(prereg["tests"]["n2_bootstrap_seed"]) + 1, resamples=resamples, + ) + uniform_60 = lib.round6(sum(1 for offset in coarse if abs(offset) <= 60) / len(coarse)) + uniform_120 = lib.round6(sum(1 for offset in coarse if abs(offset) <= 120) / len(coarse)) + + # N3b. Separate from the fit loop so a replay failure does not drop N1. + n3b_rows = [] + gate_radius = int(prereg["n3_protocol"]["radius_gate"]) + for index, case in enumerate(cases, start=1): + print(f"replay {index}/{len(cases)} {case['case_id']}", flush=True) + by_radius = {} + for radius in radii: + try: + snapshot = six_probe(case, radius) + except Exception as exc: # noqa: BLE001 — recorded as a miss, message dropped + snapshot = { + "start": None, + "end": None, + "width": None, + "truth_in_range": False, + "leader_times": [], + "all_within_1min": False, + "any_exact": False, + "error": type(exc).__name__, + } + traceback.print_exc(limit=2) + events = lib.day_events(case) + span = delivery_offsets(snapshot["start"], snapshot["end"], case["birth"]["time"], radius) + if events and span: + skies = collect_skies(case["birth"], span) + bundle = lord_bundle(skies, events) + counts = rule_counts(bundle, events, prereg["rules"]["n3_primary"]) + leaders = lib.rank_offsets(counts) + nadi_hit = lib.leaders_within(leaders, int(prereg["n3_protocol"]["nadi_limit_seconds"])) + tied = len(leaders) + else: + nadi_hit = False + tied = 0 + by_radius[str(radius)] = { + "truth_in_range": snapshot["truth_in_range"], + "width": snapshot["width"], + "baseline_all_within_1min": snapshot["all_within_1min"], + "baseline_any_exact": snapshot["any_exact"], + "nadi_all_within_60s": nadi_hit, + "nadi_tied_leaders": tied, + "improvement": int(nadi_hit) - int(snapshot["all_within_1min"]), + "error": snapshot["error"], + } + n3b_rows.append({ + "case_id": case["case_id"], + "rounded_5min": case["case_id"] in rounded, + "radii": by_radius, + }) + + def n3b_summary(radius: int) -> dict[str, Any]: + key = str(radius) + improvements = [int(row["radii"][key]["improvement"]) for row in n3b_rows] + inside = sum(1 for row in n3b_rows if row["radii"][key]["truth_in_range"]) + ci = lib.bootstrap_mean_ci(improvements, seed=int(prereg["tests"]["n3_bootstrap_seed"]) + radius, resamples=resamples) + return { + "radius": radius, + "truth_in_range": f"{inside}/{len(n3b_rows)}", + "truth_in_range_count": inside, + "truth_in_range_pass": inside >= 76 if len(n3b_rows) == 77 else None, + "baseline_all_within_1min": split_rates( + [{"case_id": row["case_id"], "hit": row["radii"][key]["baseline_all_within_1min"]} for row in n3b_rows], + "hit", rounded, + ), + "baseline_any_exact": rate([row["radii"][key]["baseline_any_exact"] for row in n3b_rows]), + "nadi_all_within_60s": split_rates( + [{"case_id": row["case_id"], "hit": row["radii"][key]["nadi_all_within_60s"]} for row in n3b_rows], + "hit", rounded, + ), + "improvement_ci": ci, + "pass": bool(ci["excludes_zero_positive"] and len(n3b_rows) == 77 and inside >= 76), + } + + def n4_aggregate(label: str, field: str) -> dict[str, Any]: + values = [] + d150_flags = [] + weights = [] + for row in n4_rows: + item = row["perturbations"].get(label) or {} + if item.get("blocked") or item.get(field) is None: + continue + values.append(float(item[field])) + weights.append(int(row["day_events"])) + if "d150_changed" in item: + d150_flags.append(bool(item["d150_changed"])) + event_rate = None + if weights and sum(weights): + event_rate = lib.round6(sum(value * weight for value, weight in zip(values, weights, strict=True)) / sum(weights)) + return { + "event_rate": event_rate, + "case_rate": rate(d150_flags) if field == "level5" else None, + "d150_case_rate": rate(d150_flags), + "cases": len(values), + } + + wording_hits = [] + for name in prereg["n4"]["ayanamsas"]: + if name == prereg["ayanamsa"]: + continue + label = f"ayanamsa:{name}" + level5 = n4_aggregate(label, "level5") + wording_hits.append({ + "switch": f"{prereg['ayanamsa']}->{name}", + "level5_event_rate": level5["event_rate"], + "d150_case_rate": level5["d150_case_rate"], + "over_20pct": bool( + (level5["event_rate"] is not None and level5["event_rate"] > 0.20) + or (level5["d150_case_rate"] is not None and level5["d150_case_rate"] > 0.20) + ), + }) + wording_ban = any(item["over_20pct"] for item in wording_hits) + + n3a = { + "grid": "coarse", + "uniform_within_60s": uniform_60, + "uniform_within_120s": uniform_120, + "published_engine_prior_top1_pm10": 0.18, + "hit_60s": split_rates( + [{"case_id": row["case_id"], "hit": row["hit_60s"]} for row in n3a_rows], "hit", rounded, + ), + "hit_120s": split_rates( + [{"case_id": row["case_id"], "hit": row["hit_120s"]} for row in n3a_rows], "hit", rounded, + ), + } + n3b = {str(radius): n3b_summary(radius) for radius in radii} + n3_pass = bool(n3b.get(str(gate_radius), {}).get("pass")) + + return { + "prereg_sha256": lib.file_sha256(lib.PREREG_PATH), + "library_sha256": lib.file_sha256(Path(__file__).with_name("nadi_seconds_lib.py")), + "runner_sha256": lib.file_sha256(Path(__file__)), + "dataset_sha256": lib.file_sha256(lib.HOLDOUT_V5), + "swisseph": __import__("swisseph").version, + "limit": limit, + "n0": { + "audit": {key: value for key, value in audit.items() if key != "per_case"}, + "per_case_day_events": audit["per_case"], + "reconciliation_mismatch_cases": [row["case_id"] for row in recon if row["mismatches"]], + "reconciliation_events": sum(row["events"] for row in recon), + }, + "n1": { + "pass": n1_pass and not limit, + "comparisons": n1, + "fraction_grid_explains_all_mean": {rule: mean_defined(values) for rule, values in fraction_rows.items()}, + "cases": case_fit, + }, + "n2": { + "pass": bool(n2_ci["excludes_zero_positive"]) and not limit, + "versus_random": n2_ci, + "versus_placebo_date": n2_placebo_ci, + "cases": n2_rows, + }, + "n3": { + "pass": n3_pass and not limit, + "n3a": n3a, + "n3b": n3b, + "cases_n3a": n3a_rows, + "cases_n3b": n3b_rows, + }, + "n4": { + "wording_ban_second_scale": wording_ban, + "ayanamsa_switches": wording_hits, + "time": {str(seconds): n4_aggregate(f"time:{int(seconds)}", "level5") for seconds in prereg["n4"]["time_offsets_seconds"]}, + "time_level4": {str(seconds): n4_aggregate(f"time:{int(seconds)}", "level4") for seconds in prereg["n4"]["time_offsets_seconds"]}, + "longitude_km": {str(km): n4_aggregate(f"east_km:{km}", "level5") for km in prereg["n4"]["longitude_shifts_km"]}, + "longitude_km_d150": {str(km): n4_aggregate(f"east_km:{km}", "level5")["d150_case_rate"] for km in prereg["n4"]["longitude_shifts_km"]}, + "node_mean_to_true": { + "moon_changed_cases": sum(1 for row in n4_rows if row["perturbations"]["node:true"]["moon_changed"]), + "asc_changed_cases": sum(1 for row in n4_rows if row["perturbations"]["node:true"]["asc_changed"]), + "rahu_changed_cases": sum(1 for row in n4_rows if row["perturbations"]["node:true"]["rahu_changed"]), + "level5_event_rate": n4_aggregate("node:true", "level5")["event_rate"], + }, + "cases": n4_rows, + }, + "n5": { + "pass": n5_pass and not limit, + "named_text": "blocked", + "named_text_reason": "no legal Chandra Kala Nadi source in the repository", + "classical_index_status": lib.CLASSICAL_NADI_STATUS, + "comparisons": n5, + }, + "gates": { + "n1": n1_pass and not limit, + "n2": bool(n2_ci["excludes_zero_positive"]) and not limit, + "n3": n3_pass and not limit, + "n4_blocks_second_scale_wording": wording_ban, + "n5": n5_pass and not limit, + }, + } + + +def main() -> int: + parser = argparse.ArgumentParser() + parser.add_argument("--limit", type=int, default=0) + parser.add_argument("--radii", default="10,30,60") + parser.add_argument("--out", default=str(REPORT)) + args = parser.parse_args() + prereg = load_prereg() + radii = tuple(int(item) for item in str(args.radii).split(",") if item.strip()) + payload = build(prereg, limit=int(args.limit), radii=radii) + text = json.dumps(payload, ensure_ascii=False, sort_keys=True, indent=2) + "\n" + out = Path(args.out) + out.parent.mkdir(parents=True, exist_ok=True) + out.write_text(text, encoding="utf-8") + print(json.dumps(payload["gates"], ensure_ascii=False, sort_keys=True), flush=True) + print(f"wrote {out}", flush=True) + return 0 + + +if __name__ == "__main__": + raise SystemExit(main()) diff --git a/tests/test_nadi_seconds_research.py b/tests/test_nadi_seconds_research.py new file mode 100644 index 00000000..b6b66e08 --- /dev/null +++ b/tests/test_nadi_seconds_research.py @@ -0,0 +1,179 @@ +"""Nadi / five-level dasha research helpers. Public v5 cases only; no production edits.""" + +from __future__ import annotations + +import inspect +import math +from datetime import datetime, timedelta + +import pytest + +from scripts.research import nadi_seconds_lib as lib + + +class _Guard(dict): + """Raises if a fitter touches the truth label.""" + + def __getitem__(self, key): + if key == "true_minute": + raise AssertionError("read true_minute") + return dict.__getitem__(self, key) + + def get(self, key, default=None): + if key == "true_minute": + raise AssertionError("read true_minute") + return dict.get(self, key, default) + + +def _packet() -> dict: + moment = datetime(2010, 1, 1) + lords = lib.five_lords("1990-06-15", 120.0, moment) + return { + "lords": lords, + "again": lib.five_lords("1990-06-15", 120.0, moment), + "d150": lib.d150_equal(10.25), + "classical": lib.classical_nadi_index(40.1), + "p": lib.paired_sign_flip_p([0.1, -0.2, 0.3, 0.0], seed=20261005, permutations=200), + "ci": lib.bootstrap_mean_ci([1.0, 0.0, 1.0, 0.0], seed=20261008, resamples=300), + } + + +def test_two_runs_are_byte_identical() -> None: + assert lib.canonical_bytes(_packet()) == lib.canonical_bytes(_packet()) + + +def test_five_levels_abut_nest_and_flip_on_the_boundary() -> None: + birth = "1990-06-15" + moon = 0.0 + event = datetime(1990, 6, 25) + chain = lib.chain_at(birth, moon, event) + assert len(chain) == 5 + for period in chain: + assert period["start"] <= event < period["end"] + for period in chain[:4]: + rows = lib.subdivide(period) + assert rows[0]["start"] == period["start"] + for left, right in zip(rows, rows[1:], strict=False): + assert left["end"] == right["start"] + rows = lib.subdivide(chain[0]) + boundary = rows[0]["end"] + before = lib.five_lords(birth, moon, boundary - timedelta(seconds=1)) + after = lib.five_lords(birth, moon, boundary) + assert before[1] == rows[0]["lord"] + assert after[1] == rows[1]["lord"] + assert before[1] != after[1] + + +def test_d150_part_flips_across_slice_and_sign_boundaries() -> None: + assert lib.d150_equal(0.199)["part_index"] == 0 + assert lib.d150_equal(0.201)["part_index"] == 1 + assert lib.d150_equal(29.999)["part_index"] == 149 + assert lib.d150_equal(30.001)["part_index"] == 0 + assert lib.classical_nadi_index(0.1)["status"] == "variant_unverified" + assert lib.classical_nadi_index(0.1)["index"] == 0 + assert lib.classical_nadi_index(30.1)["order"] == "reverse" + assert lib.classical_nadi_index(30.1)["index"] == 149 + assert lib.classical_nadi_index(60.0)["order"] == "from_middle" + assert lib.classical_nadi_index(60.0)["index"] == 75 + + +def test_rule_families_use_domain_house_lords_only() -> None: + assert lib.domain_target_lords(0, "career") == frozenset({"Saturn"}) + lords = ("Sun", "Moon", "Mars", "Mercury", "Jupiter") + targets = frozenset({"Mercury"}) + assert lib.explained(lords, targets, "A4") is True + assert lib.explained(lords, targets, "A5") is False + assert lib.explained(lords, targets, "B1") is True + assert lib.explained(lords, targets, "B2") is False + assert lib.explained(lords, targets, "N5A", segment_lord="Saturn") is False + assert lib.explained(lords, targets, "N5B", segment_lord="Moon") is True + assert lib.fit_count([lords, lords], ["career", "career"], 0, "A5") == 0 + assert lib.rank_offsets({-5: 1, 0: 3, 5: 3, 10: 2}) == [0, 5] + + +def test_fitter_does_not_read_true_minute() -> None: + for function in ( + lib.fit_count, + lib.rank_offsets, + lib.explained, + lib.counts_by_offset, + lib.public_fields, + lib.lords_for_events, + ): + assert "true_minute" not in inspect.getsource(function) + raw = lib.load_cases()[0] + guarded = _Guard(raw) + guarded["birth"] = _Guard(raw["birth"]) + guarded["events"] = [_Guard(event) for event in raw["events"]] + public = lib.public_fields(guarded) + assert "true_minute" not in public + sky = lib.sky_at(public["birth"], lib.candidate_moment(public["birth"], 0)) + events = lib.day_events(public)[:1] + assert events + rows = lib.lords_for_events(sky["birth_date"], sky["moon_lon"], events) + count = lib.fit_count(rows, [events[0]["domain"]], sky["asc_index"], "A5") + assert count in (0, 1) + assert lib.rank_offsets({0: 2, 5: 2, -5: 1}) == [0, 5] + + +def test_midnight_candidate_uses_the_previous_local_date_and_matches_production() -> None: + birth = { + "date": "1990-06-15", + "time": "00:00", + "latitude": 31.2, + "longitude": 121.5, + "timezone_offset": 8.0, + } + moment = lib.candidate_moment(birth, -30) + assert moment == datetime(1990, 6, 14, 23, 59, 30) + sky = lib.sky_at(birth, moment) + assert sky["birth_date"] == "1990-06-14" + event = datetime(2001, 3, 4) + assert lib.five_lords(sky["birth_date"], sky["moon_lon"], event)[:3] == lib.production_three( + sky["birth_date"], sky["moon_lon"], event, + ) + + +def test_v5_truth_audit_matches_the_task_counts() -> None: + audit = lib.truth_audit(lib.load_cases()) + assert audit["cases"] == 77 + assert audit["rounded_5min_count"] == 52 + assert audit["unrounded_count"] == 25 + assert audit["day_event_total"] == 244 + assert audit["day_events_ge_3"] == 43 + assert audit["day_events_ge_5"] == 19 + assert audit["day_events_zero"] == 4 + assert sum(audit["minute_histogram"].values()) == 77 + + +def test_first_three_levels_match_production_on_every_v5_event() -> None: + mismatches = [] + for case in lib.load_cases(): + row = lib.reconcile_case(case) + if row["mismatches"]: + mismatches.append(row["case_id"]) + assert mismatches == [] + + +def test_placebo_shift_is_seeded_and_stays_after_birth() -> None: + case = lib.public_fields(lib.load_cases()[0]) + events = lib.day_events(case) + first = lib.placebo_shift_dates( + events, case["birth"]["date"], seed=f"20261005:{case['case_id']}", low_days=30, high_days=180, + ) + second = lib.placebo_shift_dates( + events, case["birth"]["date"], seed=f"20261005:{case['case_id']}", low_days=30, high_days=180, + ) + assert first == second + birth = case["birth"]["date"] + assert all(birth < row["date"] <= "2026-09-14" for row in first) + assert [row["domain"] for row in first] == [row["domain"] for row in events] + assert [row["date"] for row in first] != [row["date"] for row in events] + + +def test_longitude_shift_is_about_ten_kilometres() -> None: + moved = lib.shift_longitude_km(121.5, 31.2, 10.0) + scale = 111.32 * math.cos(math.radians(31.2)) + assert moved == pytest.approx(121.5 + 10.0 / scale) + with pytest.raises(ValueError): + lib.shift_longitude_km(0.0, 90.0, 10.0)