From 3590ca91a1a215aa2a0840b3871db9c543401581 Mon Sep 17 00:00:00 2001 From: Evgeniy Unzhakov Date: Fri, 2 Oct 2026 08:42:51 +0300 Subject: [PATCH] fix(nuclear_data): use measured transition endpoints per branch paceENSDF stores the measured beta endpoint energy of each transition (ENSDF column 7). The parser instead derived branch endpoints as base-endpoint minus level energy with the base taken from an arbitrary first-listed branch, so Co60 produced negative branch endpoints and a hard ValueError (bs_pnpi Co60 exited 1 regardless of cutoff). - BranchInfo carries the measured endpoint_MeV - overall endpoint = max measured endpoint over all branches - decay_info_to_config prefers measured endpoints; the subtraction is retained only as a legacy fallback for hand-built DecayInfo objects - Co60/Tc99/Cs137 endpoint tests + CLI dry-run/export tests (D3) --- beta_spectrum/nuclear_data.py | 53 +++++++++++++++++++---------- tests/integration/test_cli_flows.py | 17 +++++++++ tests/quality/test_nuclear_data.py | 38 +++++++++++++++++++++ 3 files changed, 90 insertions(+), 18 deletions(-) diff --git a/beta_spectrum/nuclear_data.py b/beta_spectrum/nuclear_data.py index fe69537..aad8929 100644 --- a/beta_spectrum/nuclear_data.py +++ b/beta_spectrum/nuclear_data.py @@ -96,6 +96,10 @@ class BranchInfo: log(ft) value for this branch. transition_type : str Forbiddenness classification for this branch (A, F1, F1U, F2, ...). + endpoint_MeV : float + Measured beta endpoint energy of this transition in MeV (ENSDF + column 7). 0.0 means "not provided": consumers fall back to + base-endpoint minus level-energy subtraction. """ level_index: int @@ -103,6 +107,7 @@ class BranchInfo: intensity: float log_ft: Optional[float] transition_type: str = "A" + endpoint_MeV: float = 0.0 def _parse_nuclide_symbol(symbol: str) -> Tuple[str, int, int]: @@ -384,24 +389,19 @@ def get_decay_info_from_paceENSDF( # Extract Q-value from ground state branch (state[7] = endpoint energy in keV) # and half-life from halflife_data - endpoint_keV = 0.0 half_life_str = "unknown" + # The spectrum endpoint is the maximum measured transition endpoint + # across all branches (ENSDF column 7, keV). Selecting the ground + # state branch or the first branch is wrong whenever those branches + # are absent or non-maximal (e.g. Co60, whose first listed branch + # endpoint is below higher-level transition endpoints). + endpoint_keV = 0.0 if beta_data: - for k, v in beta_data.items(): + for _k, v in beta_data.items(): for state in v: - if state[0] == 0: # ground state branch - endpoint_keV = float(state[7]) - break - if endpoint_keV > 0: - break - - # No ground state branch ā€” use first available branch - if endpoint_keV == 0.0: - for k, v in beta_data.items(): - if v and len(v) > 0: - endpoint_keV = float(v[0][7]) - break + if len(state) > 7 and state[7] is not None: + endpoint_keV = max(endpoint_keV, float(state[7])) if halflife_data: for k, v in halflife_data.items(): @@ -447,6 +447,11 @@ def get_decay_info_from_paceENSDF( intensity=float(state[9]) / 100.0, log_ft=float(state[11]) if state[11] is not None else None, transition_type=branch_transition_type, + endpoint_MeV=( + float(state[7]) / 1000.0 + if len(state) > 7 and state[7] is not None + else 0.0 + ), ) branches.append(branch) @@ -497,12 +502,19 @@ def decay_info_to_config( SpectrumConfig Configuration ready for BetaSpectrum.from_config(). """ - # Convert BranchInfo -> BranchConfig for SpectrumConfig + # Convert BranchInfo -> BranchConfig for SpectrumConfig. Per-branch + # endpoints are the measured transition endpoint energies when + # available; the base-endpoint minus level-energy subtraction is only + # a legacy fallback for hand-built DecayInfo objects without them. branches: Optional[List[BranchConfig]] = None if info.branches: branches = [ BranchConfig( - endpoint_MeV=info.endpoint_MeV - b.level_energy_keV / 1000.0, + endpoint_MeV=( + b.endpoint_MeV + if b.endpoint_MeV > 0.0 + else info.endpoint_MeV - b.level_energy_keV / 1000.0 + ), transition_type=b.transition_type, intensity=b.intensity, ) @@ -615,8 +627,13 @@ def load_json_input(filepath: str) -> Dict[str, Any]: if not path.exists(): raise FileNotFoundError(f"Input file not found: {filepath}") - with open(path, "r") as f: - data = json.load(f) + try: + with open(path, "r", encoding="utf-8") as f: + data = json.load(f) + except OSError as exc: + raise FileNotFoundError(f"Input file not readable: {filepath}") from exc + except json.JSONDecodeError as exc: + raise ValueError(f"Invalid JSON in input file {filepath}: {exc}") from exc # Fill in defaults for missing optional fields for key, default in DEFAULT_JSON_SCHEMA.items(): diff --git a/tests/integration/test_cli_flows.py b/tests/integration/test_cli_flows.py index 8b08a66..57b044d 100644 --- a/tests/integration/test_cli_flows.py +++ b/tests/integration/test_cli_flows.py @@ -163,6 +163,23 @@ class TestCLIPaceENSDF: assert exc.value.code == 1 assert not (tmp_path / "tc99deg.csv").exists() + def test_paceensdf_co60_dry_run_multi_branch(self, capsys: Any) -> None: + """Co60 parses into a healthy multi-branch config (D3): no crash, + all branch endpoints positive.""" + main(["Co60", "--dry-run"]) + out = capsys.readouterr().out + assert "Multi-branch mode: 3 branches" in out + + def test_paceensdf_co60_export(self, tmp_path: Any) -> None: + out = str(tmp_path / "co60.csv") + main(["Co60", "--output", out]) + df = read_export(out) + area = np.trapezoid(df["total_spectrum"], df["energy_MeV"]) + assert np.isclose(area, 1.0, rtol=1e-3) + eps = header_lines_of(out) + assert any("Eā‚€=1492.0" in h for h in eps) + assert not any("excluded" in h for h in eps) + def test_cli_builds_paceensdf_config_via_library( self, monkeypatch: Any, capsys: Any ) -> None: diff --git a/tests/quality/test_nuclear_data.py b/tests/quality/test_nuclear_data.py index 2e54089..c85cd45 100644 --- a/tests/quality/test_nuclear_data.py +++ b/tests/quality/test_nuclear_data.py @@ -294,6 +294,44 @@ class TestDecayInfoToConfig: # --------------------------------------------------------------------------- +@pytest.mark.slow +@pytest.mark.requires_ensdf +@pytest.mark.skipif(not _has_paceENSDF(), reason="paceENSDF not installed") +class TestPaceENSDFBranchEndpoints: + """Per-branch endpoints are the measured transition endpoint energies + (ENSDF column 7), not base-endpoint minus level energy subtractions + with a base taken from an arbitrary first branch (D3).""" + + def _branch_endpoints(self, nuclide: str) -> list[float]: + info = get_decay_info_from_paceENSDF(nuclide, "beta_minus") + cfg = decay_info_to_config(info) + assert cfg.branches is not None + return sorted( + (b.endpoint_MeV for b in cfg.branches), + reverse=True, + ) + + def test_co60_branch_endpoints_positive(self) -> None: + # Measured endpoints: 1492.0, 670.0, 317.88 keV (levels above the + # first branch's endpoint previously produced negative values and + # a hard ValueError from BranchConfig validation). + eps = self._branch_endpoints("Co60") + assert eps == pytest.approx([1.492, 0.670, 0.31788], abs=1e-6) + assert all(e > 0 for e in eps) + + def test_tc99_branch_endpoints(self) -> None: + eps = self._branch_endpoints("Tc99") + assert eps == pytest.approx([0.2975, 0.2080], abs=1e-6) + + def test_cs137_branch_endpoints(self) -> None: + eps = self._branch_endpoints("Cs137") + assert eps == pytest.approx([1.176, 0.89213, 0.51403], abs=1e-6) + + def test_overall_endpoint_is_max_branch_endpoint(self) -> None: + info = get_decay_info_from_paceENSDF("Co60", "beta_minus") + assert info.endpoint_MeV == pytest.approx(1.492, abs=1e-6) + + @pytest.mark.requires_ensdf @pytest.mark.skipif(not _has_paceENSDF(), reason="paceENSDF not installed") class TestPaceENSDFIntegration: