Add a new feature: take into account of the condition of the background systematics
@@ -27,6 +27,8 @@ class SensitivityEstimator(Estimator):
Minimum significance. Default is 5.
gamma_min : float, optional
Minimum number of gamma-rays. Default is 10.
+ bkg_syst_fraction : float, optional
+ Fraction of background counts above which the number of gamma-rays is. Default is 0.05
Examples
--------
@@ -37,7 +39,7 @@ class SensitivityEstimator(Estimator):
tag = "SensitivityEstimator"
def __init__(
- self, spectrum=None, n_sigma=5.0, gamma_min=10,
+ self, spectrum=None, n_sigma=5.0, gamma_min=10, bkg_syst_fraction=0.05
):
if spectrum is None:
@@ -46,6 +48,7 @@ class SensitivityEstimator(Estimator):
self.spectrum = spectrum
self.n_sigma = n_sigma
self.gamma_min = gamma_min
+ self.bkg_syst_fraction = bkg_syst_fraction
def estimate_min_excess(self, dataset):
"""Estimate minimum excess to reach the given significance.
@@ -68,6 +71,8 @@ class SensitivityEstimator(Estimator):
excess_counts = stat.excess_matching_significance(self.n_sigma)
is_gamma_limited = excess_counts < self.gamma_min
excess_counts[is_gamma_limited] = self.gamma_min
+ bkg_syst_limited = excess_counts < self.bkg_syst_fraction * dataset.counts_off_normalised.data
+ excess_counts[bkg_syst_limited] = self.bkg_syst_fraction * dataset.counts_off_normalised.data[bkg_syst_limited]
excess = Map.from_geom(geom=dataset._geom, data=excess_counts)
return excess
@@ -97,11 +102,13 @@ class SensitivityEstimator(Estimator):
dnde = phi_0.data[:, 0, 0] * dnde_model * energy ** 2
return dnde.to("erg / (cm2 s)")
- def _get_criterion(self, excess):
+ def _get_criterion(self, excess, bkg):
is_gamma_limited = excess == self.gamma_min
+ is_bkg_syst_limited = excess == bkg * self.bkg_syst_fraction
criterion = np.chararray(excess.shape, itemsize=12)
criterion[is_gamma_limited] = "gamma"
- criterion[~is_gamma_limited] = "significance"
+ criterion[is_bkg_syst_limited] = "bkg"
+ criterion[~np.logical_or(is_gamma_limited, is_bkg_syst_limited)] = "significance"
return criterion
def run(self, dataset):
@@ -120,7 +127,7 @@ class SensitivityEstimator(Estimator):
energy = dataset._geom.axes["energy"].center
excess = self.estimate_min_excess(dataset)
e2dnde = self.estimate_min_e2dnde(excess, dataset)
- criterion = self._get_criterion(excess.data.squeeze())
+ criterion = self._get_criterion(excess.data.squeeze(), dataset.counts_off_normalised.data.squeeze())
return Table(
[
