Skip to content

Adding to the output information about the intermediary output to better understand computation performed by BAccMod #90

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
mdebony wants to merge 2 commits into
base: main
Choose a base branch
from
Draft

Adding to the output information about the intermediary output to better understand computation performed by BAccMod #90

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
e14118b
Intermediary output for models
mdebony Mar 30, 2026
288c6f8
Fix
mdebony Apr 21, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
69 changes: 61 additions & 8 deletions baccmod/base_acceptance_map_creator.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -180,7 +180,7 @@ def __init__(self,
self.mini_irf_time_resolution = mini_irf_time_resolution

@abstractmethod
def create_model(self, observations: Observations) -> BackgroundIRF:
def create_model(self, observations: Observations, add_intermediary_output: bool = False) -> BackgroundIRF:
"""
Abstract method to calculate an acceptance map from a list of observations.

Expand All @@ -190,6 +190,8 @@ def create_model(self, observations: Observations) -> BackgroundIRF:
----------
observations : gammapy.data.observations.Observations
The collection of observations used to create the acceptance map.
add_intermediary_output: bool
Add for each model a baccmod field in meta with all the intermediary state of the computation

Returns
-------
Expand Down Expand Up @@ -732,7 +734,8 @@ def _cluster_observations(self,
return observations_dict, split_obs

def _create_model_cos_zenith_binned(self,
observations: Observations
observations: Observations,
add_intermediary_output: bool = False,
) -> BackgroundCollectionZenith:
"""
Calculate a model for each cos zenith bin
Expand All @@ -741,6 +744,8 @@ def _create_model_cos_zenith_binned(self,
----------
observations : gammapy.data.observations.Observations
The collection of observations used to make the background model
add_intermediary_output: bool
Add for each model a baccmod field in meta with all the intermediary state of the computation

Returns
-------
Expand Down Expand Up @@ -877,7 +882,7 @@ def _create_model_cos_zenith_binned(self,
binned_model = []
for i, binned_obs in enumerate(binned_observations):
logger.info(f"Creating model for the bin at cos zenith = {np.round(bin_center[i], 2)}°")
binned_model.append(self.create_model(binned_obs))
binned_model.append(self.create_model(binned_obs, add_intermediary_output))

dict_binned_model = {}
for i in range(len(binned_model)):
Expand All @@ -894,7 +899,8 @@ def _create_model_cos_zenith_binned(self,
def _create_acceptance_model(self,
off_observations: dict[str, Observations],
zenith_binning: bool = False,
zenith_interpolation: bool = False
zenith_interpolation: bool = False,
add_intermediary_output: bool = False,
) -> Union[BackgroundIRF,BackgroundCollectionZenith]:
"""
A function to create a combined background model for observations previously separated in subsets.
Expand All @@ -903,6 +909,8 @@ def _create_acceptance_model(self,
off_observations: gammapy.data.observations.Observations
zenith_binning: bool
zenith_interpolation: bool
add_intermediary_output: bool
Add for each model a baccmod field in meta with all the intermediary state of the computation

Returns
-------
Expand All @@ -911,9 +919,9 @@ def _create_acceptance_model(self,
models={}
for key in off_observations.keys():
if zenith_interpolation or zenith_binning:
models[key] = self._create_model_cos_zenith_binned(observations=off_observations[key])
models[key] = self._create_model_cos_zenith_binned(observations=off_observations[key], add_intermediary_output=add_intermediary_output)
else:
models[key] = self.create_model(observations=off_observations[key])
models[key] = self.create_model(observations=off_observations[key], add_intermediary_output=add_intermediary_output)
if self.azimuth_east_west_splitting:
model=BackgroundCollectionZenithSplitAzimuth(bkg_east=models['east'], bkg_west=models['west'],
interpolation_type=self.interpolation_type,
Expand Down Expand Up @@ -1258,22 +1266,27 @@ def _normalise_model_per_run(self,

def create_acceptance_model(self,
off_observations: Observations,
zenith_binning: bool = False
zenith_binning: bool = False,
add_intermediary_output: bool = False,
) -> Union[BackgroundIRF,BackgroundCollectionZenith]:
"""
A high level function to create a background model from observations.
Parameters
----------
off_observations: gammapy.data.observations.Observations
zenith_binning: bool
add_intermediary_output: bool
Add for each model a baccmod field in meta with all the intermediary state of the computation

Returns
-------
model : BackgroundIRF or BackgroundCollection
"""

off_observations_sets, _ = self._cluster_observations(off_observations)
return self._create_acceptance_model(off_observations_sets, zenith_binning)
return self._create_acceptance_model(off_observations_sets,
zenith_binning,
add_intermediary_output=add_intermediary_output)


def create_acceptance_map_per_observation(self,
Expand Down Expand Up @@ -1338,3 +1351,43 @@ def create_acceptance_map_per_observation(self,
acceptance_map = self._normalise_model_per_run(observations, acceptance_map)

return acceptance_map

@staticmethod
def meta_observations(observations: Observations) -> List[Dict]:
"""
Return the meta block with the information detail of which observation where used to create the model

Parameters
----------
observations : gammapy.data.observations.Observations
The collection of observations used to make the acceptance map

Returns
-------
meta : List of Dict
"""

meta = []

for obs in observations:
meta.append({
"obs_id": obs.obs_id,
"obs_tstart": obs.tstart,
"obs_tstop": obs.tstop,
"obs_mid_alt_az": obs.get_pointing_altaz(obs.tmid),
"obs_pointing": obs.get_pointing_icrs(obs.tmid)
})
return meta

def meta_exclusion_region(self) -> List[Dict]:
"""
Return the meta block with the information detail of which exclusion regions where used to create the model

Returns
-------
meta : List of Dict
"""

if self.exclude_regions is None:
return []
return copy.deepcopy(self.exclude_regions)
19 changes: 16 additions & 3 deletions baccmod/grid3d_acceptance_map_creator.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -165,14 +165,16 @@ def f(*args):

return fnc(x, y, **m.values.to_dict())

def create_model(self, observations: Observations) -> Background3D:
def create_model(self, observations: Observations, add_intermediary_output: bool = False) -> Background3D:
"""
Calculate a 3D grid acceptance map

Parameters
----------
observations : gammapy.data.observations.Observations
The collection of observations used to make the acceptance map
add_intermediary_output: bool
Add for each model a meta_baccmod field with all the intermediary state of the computation

Returns
-------
Expand Down Expand Up @@ -221,12 +223,23 @@ def create_model(self, observations: Observations) -> Background3D:
data_background = corrected_counts / solid_angle[np.newaxis, :, :] / energy_axis_computation.bin_width[:, np.newaxis,
np.newaxis] / livetime

data_background = self._interpolate_bkg_to_energy_axis(data_background, energy_axis_computation)
data_background_interpolated = self._interpolate_bkg_to_energy_axis(data_background, energy_axis_computation)

acceptance_map = Background3D(axes=[self.energy_axis, extended_offset_axis_x, extended_offset_axis_y],
data=data_background.to(u.Unit('s-1 MeV-1 sr-1')),
data=data_background_interpolated.to(u.Unit('s-1 MeV-1 sr-1')),
fov_alignment=FoVAlignment.ALTAZ)


if add_intermediary_output:
meta = {'observations': self.meta_observations(observations),
'exclusion_region': self.meta_exclusion_region(),
'count_map': count_background,
'corrected_exposure': exp_map_background,
'total_exposure': exp_map_background_total,
'energy_axis_computation': energy_axis_computation,
'background_model_computation': data_background}
acceptance_map.meta['baccmod'] = meta

return acceptance_map

def _create_base_computation_map(self, observations: Observations) -> Tuple[np.ndarray, WcsNDMap, WcsNDMap, u.Quantity, MapAxis]:
Expand Down
14 changes: 13 additions & 1 deletion baccmod/radial_acceptance_map_creator.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -65,14 +65,16 @@ def __init__(self,
exclude_regions=exclude_regions,
**kwargs)

def create_model(self, observations: Observations) -> Background2D:
def create_model(self, observations: Observations, add_intermediary_output: bool = False) -> Background2D:
"""
Calculate a radial acceptance map

Parameters
----------
observations : Observations
The collection of observations used to make the acceptance map
add_intermediary_output: bool
Add for each model a baccmod field in meta with all the intermediary state of the computation

Returns
-------
Expand Down Expand Up @@ -104,6 +106,16 @@ def create_model(self, observations: Observations) -> Background2D:

acceptance_map = Background2D(axes=[self.energy_axis, self.offset_axis], data=self._interpolate_bkg_to_energy_axis(data_background, energy_axis_computation))

if add_intermediary_output:
meta = {'observations': self.meta_observations(observations),
'exclusion_region': self.meta_exclusion_region(),
'count_map': count_map_background,
'corrected_exposure': exp_map_background,
'total_exposure': exp_map_background_total,
'energy_axis_computation': energy_axis_computation,
'background_model_computation': data_background}
acceptance_map.meta['baccmod'] = meta

return acceptance_map

def _create_base_computation_map(self, observations: Observations) -> Tuple[np.ndarray, WcsNDMap, WcsNDMap, u.Quantity, MapAxis]:
Expand Down
18 changes: 9 additions & 9 deletions test/integration_test.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@ def test_integration_3D(self):
offset_axis=self.offset_axis,
oversample_map=5,
exclude_regions=self.exclude_region_PKS_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155, add_intermediary_output=True)
assert type(background_model) is Background3D
reference = Background3D.read('ressource/test_data/reference_model/pks_2155_3D.fits')
self._print_model_precision(background_model, reference)
Expand All @@ -77,7 +77,7 @@ def test_integration_spatial_fit(self):
oversample_map=5,
exclude_regions=self.exclude_region_PKS_2155,
method='fit')
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155, add_intermediary_output=True)
assert type(background_model) is Background3D
reference = Background3D.read('ressource/test_data/reference_model/pks_2155_spatial_fit_bkg.fits')
self._print_model_precision(background_model, reference)
Expand All @@ -90,7 +90,7 @@ def test_integration_2D(self):
offset_axis=self.offset_axis,
oversample_map=5,
exclude_regions=self.exclude_region_PKS_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155, add_intermediary_output=True)
assert type(background_model) is Background2D
reference = Background2D.read('ressource/test_data/reference_model/pks_2155_2D.fits')
self._print_model_precision(background_model, reference)
Expand All @@ -104,7 +104,7 @@ def test_integration_3D_irregular_computation_axis(self):
offset_axis=self.offset_axis,
oversample_map=5,
exclude_regions=self.exclude_region_PKS_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155, add_intermediary_output=True)
assert type(background_model) is Background3D
reference = Background3D.read('ressource/test_data/reference_model/pks_2155_3D_bkg_irregular_energy.fits')
self._print_model_precision(background_model, reference)
Expand All @@ -119,7 +119,7 @@ def test_integration_spatial_fit_irregular_computation_axis(self):
oversample_map=5,
exclude_regions=self.exclude_region_PKS_2155,
method='fit')
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155, add_intermediary_output=True)
assert type(background_model) is Background3D
reference = Background3D.read('ressource/test_data/reference_model/pks_2155_spatial_fit_bkg_irregular_energy.fits')
self._print_model_precision(background_model, reference)
Expand All @@ -133,7 +133,7 @@ def test_integration_2D_irregular_computation_axis(self):
offset_axis=self.offset_axis,
oversample_map=5,
exclude_regions=self.exclude_region_PKS_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155, add_intermediary_output=True)
assert type(background_model) is Background2D
reference = Background2D.read('ressource/test_data/reference_model/pks_2155_2D_bkg_irregular_energy.fits')
self._print_model_precision(background_model, reference)
Expand All @@ -151,7 +151,7 @@ def test_integration_3D_dynamic_irregular_computation_axis(self):
dynamic_energy_axis_target_statistics=100,
oversample_map=5,
exclude_regions=self.exclude_region_PKS_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155, add_intermediary_output=True)
assert type(background_model) is Background3D
reference = Background3D.read('ressource/test_data/reference_model/pks_2155_3D_bkg_dynamic_irregular_energy.fits')
self._print_model_precision(background_model, reference)
Expand All @@ -168,7 +168,7 @@ def test_integration_spatial_fit_dynamic_irregular_computation_axis(self):
oversample_map=5,
exclude_regions=self.exclude_region_PKS_2155,
method='fit')
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155, add_intermediary_output=True)
assert type(background_model) is Background3D
reference = Background3D.read('ressource/test_data/reference_model/pks_2155_spatial_fit_bkg_dynamic_irregular_energy.fits')
self._print_model_precision(background_model, reference)
Expand All @@ -184,7 +184,7 @@ def test_integration_2D_dynamic_irregular_computation_axis(self):
dynamic_energy_axis_target_statistics=500,
oversample_map=5,
exclude_regions=self.exclude_region_PKS_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155)
background_model = bkg_maker.create_model(observations=self.obs_collection_pks_2155, add_intermediary_output=True)
assert type(background_model) is Background2D
reference = Background2D.read('ressource/test_data/reference_model/pks_2155_2D_bkg_dynamic_irregular_energy.fits')
self._print_model_precision(background_model, reference)
Expand Down
Loading