Repository: gammapy/gammapy
Pull request #5655 (2025-01-02 10:26:29)

Support irf exposure map creation from dl4 irf with MapDatasetMaker
Support irf exposure map creation from dl4 irf with MapDatasetMaker
Related to #5649


Patch:
diff --git a/gammapy/makers/map.py b/gammapy/makers/map.py
index d8d39cb99..aea3e79aa 100644
--- a/gammapy/makers/map.py
+++ b/gammapy/makers/map.py
@@ -4,7 +4,7 @@ import astropy.units as u
 from astropy.table import Table
 from regions import PointSkyRegion
 from gammapy.datasets import MapDatasetMetaData
-from gammapy.irf import EDispKernelMap, PSFMap, RecoPSFMap
+from gammapy.irf import EDispKernelMap, PSFMap
 from gammapy.maps import Map
 from .core import Maker
 from .utils import (
@@ -285,9 +285,15 @@ class MapDatasetMaker(Maker):
         edisp : `~gammapy.irf.EDispKernelMap`
             Energy dispersion kernel map.
         """
-        if isinstance(observation.edisp, EDispKernelMap):
-            exposure = None
-            interp_map = observation.edisp.edisp_map.interp_to_geom(geom)
+        edisp = observation.edisp
+        if isinstance(edisp, EDispKernelMap):
+            if edisp.exposure_map:
+                exposure = edisp.exposure_map.interp_to_geom(
+                    geom.squash(axis_name="energy")
+                )
+            else:
+                exposure = None
+            interp_map = edisp.edisp_map.interp_to_geom(geom)
             return EDispKernelMap(edisp_kernel_map=interp_map, exposure_map=exposure)
 
         exposure = self.make_exposure_irf(geom.squash(axis_name="energy"), observation)
@@ -319,10 +325,15 @@ class MapDatasetMaker(Maker):
         """
         psf = observation.psf
 
-        if isinstance(psf, RecoPSFMap):
-            return RecoPSFMap(psf.psf_map.interp_to_geom(geom))
-        elif isinstance(psf, PSFMap):
-            return PSFMap(psf.psf_map.interp_to_geom(geom))
+        if isinstance(psf, PSFMap):
+            if psf.exposure_map:
+                exposure_map = psf.exposure_map.interp_to_geom(
+                    geom.squash(axis_name="rad")
+                )
+            else:
+                exposure_map = None
+            return psf.__class__(psf.psf_map.interp_to_geom(geom), exposure_map)
+
         exposure = self.make_exposure_irf(geom.squash(axis_name="rad"), observation)
 
         return make_psf_map(


Test patch:
diff --git a/gammapy/makers/tests/test_map.py b/gammapy/makers/tests/test_map.py
index a51363372..07b5088cf 100644
--- a/gammapy/makers/tests/test_map.py
+++ b/gammapy/makers/tests/test_map.py
@@ -367,6 +367,8 @@ def test_interpolate_map_dataset():
         geom=geom_target, energy_axis_true=energy_true, rad_axis=rad_axis, name="test"
     )
     dataset = maker.run(dataset, obs)
+    assert dataset.psf.exposure_map and np.all(dataset.psf.exposure_map.data == 1.0)
+    assert dataset.edisp.exposure_map and np.all(dataset.edisp.exposure_map.data == 0.0)
 
     # test counts
     assert dataset.counts.data.sum() == nr_ev
@@ -483,7 +485,6 @@ def test_dataset_hawc():
     results["NN"] = [6.57154247837e16, 62, 0.76743538]
 
     for which in ["GP", "NN"]:
-
         # paths and file names
         data_path = "$GAMMAPY_DATA/hawc/crab_events_pass4/"
         hdu_filename = "hdu-index-table-" + which + "-Crab.fits.gz"