Merge pull request #52 from VERITAS-Observatory/0.5.0-dev3

0.5.0 dev3

Author: GernotMaier

LICENSE

@@ -1,6 +1,6 @@
 BSD 3-Clause License
 
-Copyright (c) 2023, 2024, 2025 Gernot Maier
+Copyright (c) 2023-2025 Gernot Maier
 
 Redistribution and use in source and binary forms, with or without
 modification, are permitted provided that the following conditions are met:

README.md

@@ -53,7 +53,7 @@ Reflected region analysis using gammapy:
 ```console
 python v2dl5/scripts/reflected_region_analysis.py \
     --run_list my_output_dir/runlist.txt \
-    --config examples/reflected_region_analysis.yml \
+    --config examples/reflected_region.yml \
     --output_dir my_output_dir
 ```
 

examples/reflected_region.yml

@@ -34,13 +34,13 @@ datasets:
     background:
         method: reflected
 fit:
-    fit_range: {min: 0.1 TeV, max: 20 TeV}
     model: pl
     index: 2.5
     reference_energy: 1 TeV
 
 flux_points:
     energy: {min: 0.1 TeV, max: 20 TeV, nbins: 10}
+    flux: {min: 1.e-15, max: 1.e-9}
 
 light_curve:
     energy: {min: 1. TeV, max: 100 TeV}

v2dl5/analysis.py

@@ -86,6 +86,8 @@ def plot(self):
             plotter.plot_spectra(
                 flux_points=self.flux_points,
                 model=self.spectral_model,
+                y_min=self.args_dict["flux_points"].get("flux", {}).get("min"),
+                y_max=self.args_dict["flux_points"].get("flux", {}).get("max"),
             )
         else:
             self._logger.warning("Skipping spectral plots because fit failed")

v2dl5/data.py

@@ -14,7 +14,7 @@ class Data:
     Data class holding data store and observations.
 
     Allows to select data from run list or based on
-    target coordinates (and observation cone).
+    target coordinates and observation cone.
 
     Parameters
     ----------
@@ -33,7 +33,7 @@ def __init__(self, args_dict, target=None):
         """
         Initialize Data object.
 
-        Uses run_list if not set to None, otherwise selects data
+        Uses 'run_list' from args_dict if not set to None, otherwise selects data
         according to target coordinates and observation cone.
 
         """
@@ -84,7 +84,7 @@ def get_observations(self, reflected_region=True, skip_missing=False):
 
     def _from_run_list(self, run_list):
         """
-        Read run_list from file and select data.
+        Read run list from file and select data.
 
         Parameters
         ----------
@@ -140,15 +140,15 @@ def get_on_region_radius(self):
         try:
             rad_max = {obs.rad_max.data[0][0] for obs in observations}
         except IndexError:
-            self._logger.error("Rad max not found in observations.")
+            self._logger.error("On region radius not found in observations.")
             raise
 
         if len(rad_max) > 1:
-            self._logger.error("Rad max is not the same for all observations.")
+            self._logger.error("On region radius not the same for all observations.")
             raise ValueError
 
         on_region = rad_max.pop() * u.deg
-        self._logger.info(f"On region size: {on_region}")
+        self._logger.info(f"On region radius: {on_region}")
 
         return on_region
 
@@ -166,7 +166,7 @@ def get_max_wobble_distance(self, fov=3.5 * u.deg):
         Returns
         -------
         max_offset : astropy.units.Quantity
-            Maximum offset.
+            Maximum offset (radius of FoV).
 
         """
         woff = np.array(

v2dl5/plot.py

@@ -44,14 +44,30 @@ def plot_maps(self, exclusion_mask=None):
         self.plot_regions(exclusion_mask=exclusion_mask)
         self.plot_theta2()
 
-    def plot_spectra(self, flux_points=None, model=None):
-        """Spectrum related plots."""
+    def plot_spectra(self, flux_points=None, model=None, y_min=None, y_max=None):
+        """
+        Spectrum related plots.
+
+        Parameters
+        ----------
+        flux_points : `~gammapy.datasets.FluxPointsDataset`
+            Flux points dataset to plot
+        model : `~gammapy.modeling.models.Model`
+            Model to plot
+        y_min : float, optional
+            Minimum y-axis value for the SED plot
+        y_max : float, optional
+            Maximum y-axis value for the SED plot
+
+
+        """
         for dataset in self.data_set:
             self.plot_fit(dataset)
 
         self.plot_flux_points(flux_points)
         self.plot_sed(
             FluxPointsDataset(data=flux_points, models=model.copy()),
+            y_min=y_min, y_max=y_max,
         )
 
     def plot_light_curves(self, light_curves):
@@ -124,6 +140,7 @@ def plot_fit(self, data_set):
             ax_spectrum, _ = data_set.plot_fit()
         except ValueError:
             return
+        # TODO: Adjust axis limits and add any necessary labels or annotations to improve plot clarity.
         ax_spectrum.set_ylim(0.1, 40)
         data_set.plot_masks(ax=ax_spectrum)
         try:
@@ -141,11 +158,13 @@ def plot_flux_points(self, flux_point_dataset):
         flux_point_dataset.plot_ts_profiles(ax=ax, sed_type="dnde")
         self._plot(plot_name="flux_points", output_dir=self.output_dir)
 
-    def plot_sed(self, flux_point_dataset):
+    def plot_sed(self, flux_point_dataset, y_min=None, y_max=None):
         """Plot spectral energy distribution."""
         kwargs_model = {"color": "grey", "ls": "--", "sed_type": "dnde"}
         kwargs_fp = {"color": "black", "marker": "o", "sed_type": "dnde"}
-        flux_point_dataset.plot_spectrum(kwargs_fp=kwargs_fp, kwargs_model=kwargs_model)
+        ax = flux_point_dataset.plot_spectrum(kwargs_fp=kwargs_fp, kwargs_model=kwargs_model)
+        if y_min is not None and y_max is not None:
+            ax.set_ylim(y_min, y_max)
         self._plot(plot_name="spectrum", output_dir=self.output_dir)
         try:
             flux_point_dataset.plot_residuals(method="diff/model")

v2dl5/sky_regions.py

@@ -52,8 +52,7 @@ def get_target(self, sky_coord=None, print_target_info=True):
 
         """
         target = None
-        lon = sky_coord.get("lon", None)
-        lat = sky_coord.get("lat", None)
+        lon, lat = sky_coord.get("lon", None), sky_coord.get("lat", None)
         if lon is not None and lat is not None:
             if sky_coord.get("frame", None) == "icrs":
                 target = SkyCoord(