From 3749220efc8b76b33971d5a279360c9f3b024c62 Mon Sep 17 00:00:00 2001 From: Evgeniy Unzhakov Date: Tue, 5 May 2026 12:43:43 +0300 Subject: [PATCH] Fix pyright type errors and correct paceENSDF Q-value extraction MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit - nuclear_data.py: Extract Q-value from ground state branch state[7] instead of halflife_data (was reading half-life seconds as Q-value) - nuclear_data.py: Use ground state branch for forbiddenness code, with fallback - nuclear_data.py: Add None guards for beta_data and pairs - spectrum.py: Pass Z=config.Z_daughter, A=config.A_number to RadiativeCorrection - phase_space.py: Add __init__ debug logging (W0, transition_type, m_e, m_nu) - phase_space.py: Add output range logging in __call__ - phase_space.py: Fix possibly unbound forbid_factor variable - logging_utils.py: Fix broken get_git_short_hash indentation - cli.py: Add explicit SpectrumConfig annotation to fix pyright type identity - cw_extractor.py: Fix T_to_W ndarray→float conversions (6 locations) - cw_extractor.py: Fix FermiFunction() calls — was passing Z as __call__ arg - fitter.py: Add explicit dict type annotation for least_squares kwargs --- beta_spectrum/cli.py | 2 ++ beta_spectrum/components/phase_space.py | 20 +++++++---- beta_spectrum/cw_extractor.py | 26 +++++++------- beta_spectrum/fitter.py | 26 +++++++------- beta_spectrum/logging_utils.py | 4 +-- beta_spectrum/nuclear_data.py | 48 +++++++++++++++++++------ beta_spectrum/spectrum.py | 8 ++++- 7 files changed, 89 insertions(+), 45 deletions(-) diff --git a/beta_spectrum/cli.py b/beta_spectrum/cli.py index 0564cb9..928b72b 100644 --- a/beta_spectrum/cli.py +++ b/beta_spectrum/cli.py @@ -28,6 +28,7 @@ import sys from typing import Any, Dict from beta_spectrum.logging_utils import LoggingConfig, setup_logging +from beta_spectrum.spectrum import SpectrumConfig def _build_parser() -> argparse.ArgumentParser: @@ -210,6 +211,7 @@ def _run(args: argparse.Namespace) -> None: detector_kwargs = _parse_detector_args(args) source_type: str + config: SpectrumConfig if args.nuclide: logger.info("Source: paceENSDF nuclide=%s", args.nuclide) config = create_config_from_source( diff --git a/beta_spectrum/components/phase_space.py b/beta_spectrum/components/phase_space.py index 669cd42..bc7a837 100644 --- a/beta_spectrum/components/phase_space.py +++ b/beta_spectrum/components/phase_space.py @@ -30,6 +30,11 @@ class PhaseSpace(SpectrumComponent): """ super().__init__(logger=logger) self.tr_type = transition_type + if self._logger: + self._logger.debug( + "PhaseSpace: W0=%.5f, transition_type=%s, m_e=%.6f, m_nu=%.6f", + W0, transition_type, m_e, m_nu, + ) if self.tr_type not in ["A", "F1", "F1U", "F2", "F2U", "F3", "F3U", "F4"]: raise ValueError( f"Cannot calculate phase space for forbiddance '{self.tr_type}'.\nShould be [A, F1, F1U, F2, F2U, F3, F3U, F4].\n" @@ -47,17 +52,20 @@ class PhaseSpace(SpectrumComponent): phase_space = p_e * w_e * p_nu * w_nu if self.tr_type in ["A", "F1"]: - forbid_factor = 1 - + forbid_factor: float | np.ndarray = 1 elif self.tr_type in ["F1U", "F2"]: forbid_factor = p_nu**2 + p_e**2 - elif self.tr_type in ["F2U", "F3"]: forbid_factor = p_nu**4 + 3 / 10 * p_nu**2 * p_e**2 + p_e**4 - - elif self.tr_type in ["F3U", "F4"]: + else: # F3U, F4 forbid_factor = ( p_nu**6 + 7 * p_nu**4 * p_e**2 + 7 * p_nu**2 * p_e**4 + p_e**6 ) - return np.asarray(phase_space * forbid_factor, dtype=np.float64) + result = np.asarray(phase_space * forbid_factor, dtype=np.float64) + if self._logger: + self._logger.debug( + "PhaseSpace: output range=[%.6e, %.6e]", + result.min(), result.max(), + ) + return result diff --git a/beta_spectrum/cw_extractor.py b/beta_spectrum/cw_extractor.py index 3e048af..abdbd40 100644 --- a/beta_spectrum/cw_extractor.py +++ b/beta_spectrum/cw_extractor.py @@ -256,9 +256,9 @@ class CWExtractor: if self.config.use_radiative: from .components.radiative import RadiativeCorrection - W0 = T_to_W(self.config.endpoint_MeV) + W0_float = float(T_to_W(self.config.endpoint_MeV)) rad = RadiativeCorrection( - W0=W0, + W0=W0_float, Z=self.config.Z_parent, A=self.config.A_number, use_endpoint_resummation=True, @@ -291,9 +291,9 @@ class CWExtractor: Extracted C(W) values with uncertainties. """ if endpoint_keV is not None: - W0 = T_to_W(endpoint_keV / 1000.0) + W0_float = float(T_to_W(endpoint_keV / 1000.0)) else: - W0 = self._determine_endpoint() + W0_float = float(self._determine_endpoint()) # Compute theoretical factors (excluding C(W)) factors = self._theoretical_factors(self.energies_W) @@ -310,12 +310,12 @@ class CWExtractor: product *= factors[key] # Endpoint factor: (W0 - W)^2 - endpoint_factor = (W0 - self.energies_W) ** 2 + endpoint_factor = (W0_float - self.energies_W) ** 2 # Overall normalization if normalize_method == "endpoint": # Use flux at lowest energy where C(W) ≈ 1 for allowed decay - low_E_mask = self.energies_W < W0 - 0.02 + low_E_mask = self.energies_W < W0_float - 0.02 if np.any(low_E_mask): C0 = np.median( self.flux[low_E_mask] @@ -325,7 +325,7 @@ class CWExtractor: C0 = 1.0 else: # Fit normalization - C0 = self._fit_normalization(W0, product, endpoint_factor) + C0 = self._fit_normalization(W0_float, product, endpoint_factor) # Extract C(W) cw = self.flux / (C0 * product * endpoint_factor) @@ -347,8 +347,8 @@ class CWExtractor: cw_errors=cw_errors, flux=self.flux.copy(), flux_errors=self.flux_errors.copy(), - endpoint_W=W0, - endpoint_keV=W_to_T(W0) * 1000.0, + endpoint_W=W0_float, + endpoint_keV=float(W_to_T(W0_float)) * 1000.0, fit_result=None, ) @@ -585,7 +585,7 @@ class CWExtractor: p = momentum(W) from .components.fermi import FermiFunction - F = FermiFunction()(W, self.config.Z_parent) + F = FermiFunction(Z=self.config.Z_parent, A=self.config.A_number)(W) kurie = np.sqrt(flux / (p * W * F)) @@ -605,7 +605,7 @@ class CWExtractor: except np.linalg.LinAlgError: continue - return T_to_W(best_endpoint / 1000.0) + return float(T_to_W(best_endpoint / 1000.0)) def _fit_kurie_endpoint(self, cw_result: CWExtractionResult) -> FitResult: """Fit Kurie plot to determine endpoint energy.""" @@ -620,7 +620,7 @@ class CWExtractor: p = momentum(W) from .components.fermi import FermiFunction - F = FermiFunction()(W, self.config.Z_parent) + F = FermiFunction(Z=self.config.Z_parent, A=self.config.A_number)(W) kurie = np.sqrt(flux / (p * W * F)) E_keV = self.energies_keV[mask] @@ -639,7 +639,7 @@ class CWExtractor: ) return fitter.fit( - x0=[np.mean(kurie), -np.mean(kurie) / E_keV[-1]], + x0=[float(np.mean(kurie)), float(-np.mean(kurie) / E_keV[-1])], param_names=["a", "b"], ) diff --git a/beta_spectrum/fitter.py b/beta_spectrum/fitter.py index 804c3cf..43cfd6a 100644 --- a/beta_spectrum/fitter.py +++ b/beta_spectrum/fitter.py @@ -390,22 +390,22 @@ class CurveFitter: np.asarray(bounds[1], dtype=np.float64), ) - ls_kwargs = dict( - fun=self._residuals, - x0=x0_array, - method=self.config.method, - max_nfev=self.config.max_nfev, - x_scale=self.config.x_scale, - loss=self.config.loss, - f_scale=self.config.f_scale, - ftol=self.config.ftol, - xtol=self.config.xtol, - gtol=self.config.gtol, - ) + ls_kwargs: dict[str, object] = { + "fun": self._residuals, + "x0": x0_array, + "method": self.config.method, + "max_nfev": self.config.max_nfev, + "x_scale": self.config.x_scale, + "loss": self.config.loss, + "f_scale": self.config.f_scale, + "ftol": self.config.ftol, + "xtol": self.config.xtol, + "gtol": self.config.gtol, + } if bounds is not None: ls_kwargs["bounds"] = bounds_array - result = least_squares(**ls_kwargs) + result = least_squares(**ls_kwargs) # type: ignore[arg-type] # Compute covariance matrix from Jacobian # Cov = (J^T J)^(-1) * χ² / (n - p) diff --git a/beta_spectrum/logging_utils.py b/beta_spectrum/logging_utils.py index e4e6e8e..98b412d 100644 --- a/beta_spectrum/logging_utils.py +++ b/beta_spectrum/logging_utils.py @@ -137,9 +137,9 @@ def get_git_short_hash(length: int = 7) -> str: str Short git hash, or "unknown" if not in a git repo. """ - try: - import subprocess + import subprocess + try: result = subprocess.run( ["git", "rev-parse", f"--short={length}", "HEAD"], capture_output=True, diff --git a/beta_spectrum/nuclear_data.py b/beta_spectrum/nuclear_data.py index 36b0df9..5fb0a0e 100644 --- a/beta_spectrum/nuclear_data.py +++ b/beta_spectrum/nuclear_data.py @@ -309,6 +309,7 @@ def _resolve_decay_index( e = pe.ENSDF() pairs = e.ensdf_pairs(edata, mode) + assert pairs is not None, "paceENSDF returned None for ensdf_pairs" matching = [(k, v) for k, v in pairs.items() if str(k[0]) == parent_symbol] @@ -378,13 +379,31 @@ def get_decay_info_from_paceENSDF( jpi_data = e.get_parent_jpi(edata, nuclide, decay_idx, mode=mode) beta_data = e.get_beta_minus(edata, nuclide, decay_idx, units="best") - # Extract Q-value and half-life - q_value_keV = 0.0 + # 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" - for k, v in halflife_data.items(): - q_value_keV = v[0] - half_life_str = str(v[2]) if len(v) > 2 else "unknown" - break + + if beta_data: + 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 halflife_data: + for k, v in halflife_data.items(): + half_life_str = str(v[2]) if len(v) > 2 else "unknown" + break # Extract parent spin/parity parent_J = jpi_data[0][0] if jpi_data else 0.0 @@ -398,10 +417,19 @@ def get_decay_info_from_paceENSDF( forbidden_code = "0A" # default transition_type = "A" # default + if beta_data is None: + beta_data = {} + for k, v in beta_data.items(): for state in v: - # state[15] contains forbiddenness code (e.g. '2F', '1UF') - if len(state) > 15 and state[15] is not None: + # Use ground state branch (level_index == 0) for forbiddenness + if state[0] == 0 and len(state) > 15 and state[15] is not None: + code = str(state[15]) + if code in FORBIDDENNESS_MAP: + forbidden_code = code + transition_type = FORBIDDENNESS_MAP[code] + # Fallback: use first branch's forbiddenness if no ground state + elif transition_type == "A" and len(state) > 15 and state[15] is not None: code = str(state[15]) if code in FORBIDDENNESS_MAP: forbidden_code = code @@ -420,8 +448,8 @@ def get_decay_info_from_paceENSDF( Z_parent=Z_parent, Z_daughter=Z_daughter, A_number=A_number, - endpoint_keV=q_value_keV, - endpoint_MeV=q_value_keV / 1000.0, + endpoint_keV=endpoint_keV, + endpoint_MeV=endpoint_keV / 1000.0, transition_type=transition_type, forbiddenness_code=forbidden_code, parent_J=parent_J, diff --git a/beta_spectrum/spectrum.py b/beta_spectrum/spectrum.py index 22adb19..3b8c755 100644 --- a/beta_spectrum/spectrum.py +++ b/beta_spectrum/spectrum.py @@ -169,7 +169,13 @@ class BetaSpectrum: if config.use_radiative: components.append( - RadiativeCorrection(W0=W0, use_endpoint_resummation=True, logger=logger) + RadiativeCorrection( + W0=W0, + Z=config.Z_daughter, + A=config.A_number, + use_endpoint_resummation=True, + logger=logger, + ) ) return cls(components, logger=logger)