diff --git a/beta_spectrum/visualize/__init__.py b/beta_spectrum/visualize/__init__.py index 2b32e9a..149f3e8 100644 --- a/beta_spectrum/visualize/__init__.py +++ b/beta_spectrum/visualize/__init__.py @@ -206,7 +206,20 @@ class BetaSpectrumAnalyzer: n = self.config.A_number title = f"Beta-decay: {n}{parent} -> {n}{daughter}" if self._n_branches > 1: - title += f" ({self._n_branches} branches)" + # Count active branches (not excluded by neutrino mass or intensity cutoff) + active_branches = 0 + for i, branch in enumerate(self._unified_branches): + branch_endpoint_keV = branch.endpoint_MeV * 1000 + neutrino_mass_keV = self.config.neutrino_mass_MeV * 1000 + is_excluded_by_mass = branch_endpoint_keV <= neutrino_mass_keV + is_excluded_by_intensity = branch.intensity < self.config.intensity_cutoff + if not is_excluded_by_mass and not is_excluded_by_intensity: + active_branches += 1 + + if active_branches < self._n_branches: + title += f" ({active_branches} of {self._n_branches} branches)" + else: + title += f" ({self._n_branches} branches)" return title def _format_branch_label(self, i: int) -> str: @@ -214,7 +227,19 @@ class BetaSpectrumAnalyzer: branch = self._unified_branches[i] e0 = branch.endpoint_MeV * 1000 intensity_pct = branch.intensity * 100 - return f"Br. {i + 1}: {e0:.1f} keV ({intensity_pct:.1f}%)" + + # Check if branch is excluded + branch_endpoint_keV = branch.endpoint_MeV * 1000 + neutrino_mass_keV = self.config.neutrino_mass_MeV * 1000 + is_excluded_by_mass = branch_endpoint_keV <= neutrino_mass_keV + is_excluded_by_intensity = branch.intensity < self.config.intensity_cutoff + + if is_excluded_by_mass: + return f"Br. {i + 1}: {e0:.1f} keV ({intensity_pct:.1f}%) - excluded (endpoint < m_nu)" + elif is_excluded_by_intensity: + return f"Br. {i + 1}: {e0:.1f} keV ({intensity_pct:.1f}%) - excluded (intensity < cutoff)" + else: + return f"Br. {i + 1}: {e0:.1f} keV ({intensity_pct:.1f}%)" # --------------------------------------------------------------------------- # Plotting — unified multi-branch approach @@ -302,7 +327,18 @@ class BetaSpectrumAnalyzer: zorder=10, ) - for i, contrib in enumerate(branch_contribs): + # Only plot active branches (not excluded by neutrino mass or intensity cutoff) + active_branch_indices = [] + for i in range(self._n_branches): + branch_endpoint_keV = self._unified_branches[i].endpoint_MeV * 1000 + neutrino_mass_keV = self.config.neutrino_mass_MeV * 1000 + is_excluded_by_mass = branch_endpoint_keV <= neutrino_mass_keV + is_excluded_by_intensity = self._unified_branches[i].intensity < self.config.intensity_cutoff + if not is_excluded_by_mass and not is_excluded_by_intensity: + active_branch_indices.append(i) + + for i in active_branch_indices: + contrib = branch_contribs[i] color = branch_colors[i % len(branch_colors)] ax.plot( self.energies_MeV, @@ -321,7 +357,7 @@ class BetaSpectrumAnalyzer: ax.legend( loc="best", fontsize=9, - ncol=min(2, self._n_branches + 1), + ncol=min(2, len(active_branch_indices) + 1), ) self._build_figure_header(fig) @@ -390,7 +426,19 @@ class BetaSpectrumAnalyzer: label="Total spectrum", zorder=10, ) - for i, contrib in enumerate(branch_contribs): + + # Only plot active branches (not excluded by neutrino mass or intensity cutoff) + active_branch_indices = [] + for i in range(n_branches): + branch_endpoint_keV = self._unified_branches[i].endpoint_MeV * 1000 + neutrino_mass_keV = self.config.neutrino_mass_MeV * 1000 + is_excluded_by_mass = branch_endpoint_keV <= neutrino_mass_keV + is_excluded_by_intensity = self._unified_branches[i].intensity < self.config.intensity_cutoff + if not is_excluded_by_mass and not is_excluded_by_intensity: + active_branch_indices.append(i) + + for i in active_branch_indices: + contrib = branch_contribs[i] color = branch_colors[i % len(branch_colors)] ax.plot( self.energies_MeV, @@ -404,7 +452,7 @@ class BetaSpectrumAnalyzer: ax.set_ylabel(f"Normalized Counts / {e_step_keV:.1f} keV", fontsize=10) ax.set_title("Spectrum", fontsize=11) ax.grid(True, alpha=0.3) - ax.legend(fontsize=8, loc="best", ncol=min(2, n_branches + 1)) + ax.legend(fontsize=8, loc="best", ncol=min(2, len(active_branch_indices) + 1)) # Panel 2: Fermi function + all phase spaces ax = axes[1] @@ -419,7 +467,7 @@ class BetaSpectrumAnalyzer: label="Fermi", alpha=0.8, ) - for i in range(n_branches): + for i in active_branch_indices: color = branch_colors[i % len(branch_colors)] branch = self._unified_branches[i] key = self._get_component_key(i, "PhaseSpace") @@ -474,7 +522,7 @@ class BetaSpectrumAnalyzer: # Panel 4: W0-dependent corrections (Radiative for each branch) ax = axes[3] - for i in range(n_branches): + for i in active_branch_indices: color = branch_colors[i % len(branch_colors)] key = self._get_component_key(i, "Radiative") if key in components: @@ -491,7 +539,7 @@ class BetaSpectrumAnalyzer: ax.set_xlabel(r"Electron kinetic energy $E$ [MeV]", fontsize=10) ax.set_title("Radiative Corrections", fontsize=11) ax.grid(True, alpha=0.3) - ax.legend(fontsize=8, loc="best", ncol=min(2, n_branches + 1)) + ax.legend(fontsize=8, loc="best", ncol=min(2, len(active_branch_indices) + 1)) self._build_figure_header(fig)