From c84e552fd89bc41d2c80a49e3c4e2d8a1e977960 Mon Sep 17 00:00:00 2001 From: Ste O'Brien Date: Mon, 12 Aug 2024 11:00:56 -0400 Subject: [PATCH] feat(background_models.py): allowing for user defined default exclusion the default exclusion can be either set by the user or set as the default on region exclusion --- .../make_background/background_models.py | 6 ++++- .../make_background/make_background.py | 16 +++++++++++++- gammapy_tools/utils/exclusion_finder.py | 22 ++++++++++++++----- 3 files changed, 37 insertions(+), 7 deletions(-) diff --git a/gammapy_tools/make_background/background_models.py b/gammapy_tools/make_background/background_models.py index 0d7afd8..6d7c39c 100644 --- a/gammapy_tools/make_background/background_models.py +++ b/gammapy_tools/make_background/background_models.py @@ -17,6 +17,7 @@ # from here from ..utils import ExclusionFinder +from typing import Optional class BackgroundModelEstimator: """ @@ -33,6 +34,7 @@ def __init__( excluded_sources: list = [], smooth_sigma: float = 1, njobs: int = 0, + default_exclusion: Optional[float] = None, ) -> None: """BackgroundModelEstimator class for estimating the background from runs @@ -45,6 +47,8 @@ def __init__( excluded_sources (list) - list of sources to be excluded (astropy.coordinates.Skycoords) njobs (int) - Number of jobs to run + default_exclusion (float) - Default exclusion radius for sources. Defaults to None + which will use the exclusion finder to determine the exclusion radius Returns ---------- @@ -85,7 +89,7 @@ def __init__( # ) # self.hawc = SourceCatalog3HWC("$GAMMAPY_DATA/catalogs/3HWC.ecsv") - self.exclusion = ExclusionFinder() + self.exclusion = ExclusionFinder(default_exclusion = default_exclusion) @staticmethod def _make_bkg2d(energy: MapAxis, offset: MapAxis, unit: str) -> Background2D: diff --git a/gammapy_tools/make_background/make_background.py b/gammapy_tools/make_background/make_background.py index d621c65..07470bd 100644 --- a/gammapy_tools/make_background/make_background.py +++ b/gammapy_tools/make_background/make_background.py @@ -131,12 +131,18 @@ def get_background_for_run(parms: tuple[float, dict]) -> tuple[str, list]: if "smooth_sigma" in config["background_selection"]: smooth_sigma = config["background_selection"]["smooth_sigma"] + default_exclusion = None + if "sky_map" in config: + if "on_exclusion_region" in config["sky_map"]: + default_exclusion = config["sky_map"]["on_exclusion_region"] + estimator = BackgroundModelEstimator( energy, offset, smooth=config["background_selection"]["smooth"], smooth_sigma=smooth_sigma, njobs=config["config"]["njobs"], + default_exclusion=default_exclusion, ) estimator.run(observations) @@ -213,8 +219,16 @@ def generate_background_from_run(parms: tuple[int, dict]) -> str: name="offset", ) + default_exclusion = None + if "sky_map" in config: + if "on_exclusion_region" in config["sky_map"]: + default_exclusion = config["sky_map"]["on_exclusion_region"] + + estimator = BackgroundModelEstimator( - energy, offset, smooth=config["background_selection"]["smooth"] + energy, offset, + smooth=config["background_selection"]["smooth"], + default_exclusion=default_exclusion, ) estimator.run(observations) if "BACKGROUND" in hdul: diff --git a/gammapy_tools/utils/exclusion_finder.py b/gammapy_tools/utils/exclusion_finder.py index e5f9a51..c551b5c 100644 --- a/gammapy_tools/utils/exclusion_finder.py +++ b/gammapy_tools/utils/exclusion_finder.py @@ -4,17 +4,24 @@ from gammapy.catalog import SourceCatalog3HWC, SourceCatalogGammaCat +from typing import Optional + class ExclusionFinder: """ Class for finding exclusion regions for the analysis """ - def __init__(self): + def __init__(self, default_exclusion: Optional[float] = 0.35): """Initilization function Load in the star catalog, GammaCat and 3HWC + + Parameters + ---------- + default_exclusion (float) - Default theta cut to use for exclusion region """ + self.default_exclusion = default_exclusion # ToDo wrap all this into package info this_dir, _ = path.split(__file__) @@ -141,8 +148,9 @@ def find_exclusion( ra: float, dec: float, theta: float = 2.0, - theta_cut: float = 0.35, + theta_cut: Optional[float] = None, mag_cut: float = 7, + star_theta_cut: float = 0.35, ) -> tuple[list[tuple[float, float, float]], list[str]]: """Find sources to exclude @@ -153,10 +161,11 @@ def find_exclusion( theta (float) - Search radius (in degrees) Defaults to 2.0 theta_cut (float) - Default theta cut to use for exclusion region - Defaults to 0.35 + Defaults to 0.35 or the default_exclusion mag_cut (float) - B Magnitude cut on star brightness Defaults to mag 7 - + star_theta_cut (float) - Theta cut for stars + Defaults to 0.35 Returns ---------- @@ -167,10 +176,13 @@ def find_exclusion( regions = [] sources_excluded = [] + if theta_cut is None: + theta_cut = self.default_exclusion + # Get the stars in the region stars = self.find_stars(ra, dec, theta, mag_cut) for star in stars: - regions.append((star["ra"], star["dec"], theta_cut)) + regions.append((star["ra"], star["dec"], star_theta_cut)) sources_excluded.append(f"star_{star['id']:0.0f}") # Find the HAWC sources