From 62360f3e5088a7e14ed2efb2698a32262f6dd426 Mon Sep 17 00:00:00 2001 From: Jan Weldert Date: Thu, 8 Aug 2024 06:07:01 -0500 Subject: [PATCH 1/6] Revive parametric Aeff stage --- pisa/stages/aeff/param.py | 139 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 139 insertions(+) create mode 100644 pisa/stages/aeff/param.py diff --git a/pisa/stages/aeff/param.py b/pisa/stages/aeff/param.py new file mode 100644 index 000000000..e2ffa8ae1 --- /dev/null +++ b/pisa/stages/aeff/param.py @@ -0,0 +1,139 @@ +""" +This is an effective area stage designed for quick studies of how effective +areas affect experimental observables and sensitivities. In addition, it is +supposed to be easily reproducible as it may rely on (phenomenological) +functions or interpolated discrete data points, dependent on energy +(and optionally cosine zenith), and which can thus be used as reference or +benchmark scenarios. +""" + + +from __future__ import absolute_import, division + +from collections import Mapping, OrderedDict + +import numpy as np +from scipy.interpolate import interp1d + +from pisa.core.stage import Stage +from pisa.utils.fileio import from_file +from pisa.utils.hash import hash_obj +from pisa.utils.log import logging +from pisa.utils import vectorizer +from pisa.utils.profiler import profile + + +__all__ = ['load_aeff_param', 'param'] + +__author__ = 'T.C. Arlen, T. Ehrhardt, S. Wren' + +__license__ = '''Copyright (c) 2014-2017, The IceCube Collaboration + Licensed under the Apache License, Version 2.0 (the "License"); + you may not use this file except in compliance with the License. + You may obtain a copy of the License at + http://www.apache.org/licenses/LICENSE-2.0 + Unless required by applicable law or agreed to in writing, software + distributed under the License is distributed on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + See the License for the specific language governing permissions and + limitations under the License.''' + + +def load_aeff_param(source): + """Load aeff parameterisation (energy- or coszen-dependent) from file + or dictionary. + Parameters + ---------- + source : string or mapping + Source of the parameterization. If string, treat as file path or + resource location and load from the file; this must yield a mapping. If + `source` is a mapping, it is used directly. See notes below on format. + Returns + ------- + aeff_params : OrderedDict + Keys are stringified flavintgroups and values are the callables that + produce aeff when called with energy or coszen values. + Notes + ----- + The mapping passed via `source` or loaded therefrom must have the format: + { + : val, + : val, + ... + } + Note that the `transform_groups` (container.name) defined + in a pipeline config file using this must match the groupings defined + above. + `val`s can be one of the following: + - Callable with one argument + - String such that `eval(val)` yields a callable with one argument + """ + if not (source is None or isinstance(source, (str, Mapping))): + raise TypeError('`source` must be string, mapping, or None') + + if isinstance(source, str): + orig_dict = from_file(source) + elif isinstance(source, Mapping): + orig_dict = source + else: + raise TypeError('Cannot load aeff parameterizations from a %s' + % type(source)) + + return orig_dict + + +class param(Stage): + """Effective area service based on parameterisation functions stored in a + .json file. + Transforms an input map of a flux of a given flavour into maps of + event rates for the two types of weak current (charged or neutral), + according to energy and cosine zenith dependent effective areas specified + by parameterisation functions. + Parameters + ---------- + params : ParamSet + Must exclusively have parameters: + aeff_energy_paramfile + aeff_coszen_paramfile + livetime + aeff_scale + """ + def __init__( + self, + **std_kwargs, + ): + expected_params = ( + 'aeff_energy_paramfile', + 'aeff_coszen_paramfile', + 'livetime', + 'aeff_scale' + ) + + super(self.__class__, self).__init__( + expected_params=expected_params, + **std_kwargs, + ) + + energy_param_source = self.params.aeff_energy_paramfile.value + coszen_param_source = self.params.aeff_coszen_paramfile.value + + self.energy_param = load_aeff_param(energy_param_source) + self.coszen_param = load_aeff_param(coszen_param_source) + + @profile + def apply_function(self): + # read out + aeff_scale = self.params.aeff_scale.m_as('dimensionless') + livetime_s = self.params.livetime.m_as('sec') + + for container in self.data: + scale = aeff_scale * livetime_s * np.ones(container.size) + if container.name in self.energy_param.keys(): + func = eval(self.energy_param[container.name]) + scale *= func(container['true_energy']) + if container.name in self.coszen_param.keys(): + func = eval(self.coszen_param[container.name]) + scale *= func(container['true_coszen']) + + container['weights'] *= scale + container.mark_changed('weights') From 6ba004191bb39dc7ff0e4ed436e3b384eea67440 Mon Sep 17 00:00:00 2001 From: Jan Weldert Date: Thu, 8 Aug 2024 06:34:10 -0500 Subject: [PATCH 2/6] Make python 3.10 happy --- pisa/stages/aeff/param.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/pisa/stages/aeff/param.py b/pisa/stages/aeff/param.py index e2ffa8ae1..c12ea7d4c 100644 --- a/pisa/stages/aeff/param.py +++ b/pisa/stages/aeff/param.py @@ -10,7 +10,8 @@ from __future__ import absolute_import, division -from collections import Mapping, OrderedDict +from collections import OrderedDict +from collections.abc import Mapping import numpy as np from scipy.interpolate import interp1d From e602aef80f4e2e07c9e6440c2345619b75444dc1 Mon Sep 17 00:00:00 2001 From: Jan Weldert Date: Fri, 9 Aug 2024 04:40:18 -0500 Subject: [PATCH 3/6] Make Thomas happy --- pisa/stages/aeff/param.py | 37 ++++++++++++++----------------------- 1 file changed, 14 insertions(+), 23 deletions(-) diff --git a/pisa/stages/aeff/param.py b/pisa/stages/aeff/param.py index c12ea7d4c..e56cedaca 100644 --- a/pisa/stages/aeff/param.py +++ b/pisa/stages/aeff/param.py @@ -10,23 +10,18 @@ from __future__ import absolute_import, division -from collections import OrderedDict from collections.abc import Mapping import numpy as np -from scipy.interpolate import interp1d from pisa.core.stage import Stage from pisa.utils.fileio import from_file -from pisa.utils.hash import hash_obj -from pisa.utils.log import logging -from pisa.utils import vectorizer from pisa.utils.profiler import profile __all__ = ['load_aeff_param', 'param'] -__author__ = 'T.C. Arlen, T. Ehrhardt, S. Wren' +__author__ = 'T.C. Arlen, T. Ehrhardt, S. Wren, J. Weldert' __license__ = '''Copyright (c) 2014-2017, The IceCube Collaboration Licensed under the Apache License, Version 2.0 (the "License"); @@ -69,8 +64,8 @@ def load_aeff_param(source): - Callable with one argument - String such that `eval(val)` yields a callable with one argument """ - if not (source is None or isinstance(source, (str, Mapping))): - raise TypeError('`source` must be string, mapping, or None') + if not isinstance(source, (str, Mapping)): + raise TypeError('`source` must be string or mapping') if isinstance(source, str): orig_dict = from_file(source) @@ -83,13 +78,13 @@ def load_aeff_param(source): return orig_dict -class param(Stage): +class param(Stage): # pylint: disable=invalid-name """Effective area service based on parameterisation functions stored in a .json file. - Transforms an input map of a flux of a given flavour into maps of - event rates for the two types of weak current (charged or neutral), - according to energy and cosine zenith dependent effective areas specified - by parameterisation functions. + Transforms an input map of a flux of a given flavour (and interaction) + into maps of event rates, according to energy and cosine zenith dependent + effective areas specified by parameterisation functions. + Requires true_energy, true_coszen, and weights to be present in the container. Parameters ---------- params : ParamSet @@ -106,27 +101,23 @@ def __init__( expected_params = ( 'aeff_energy_paramfile', 'aeff_coszen_paramfile', - 'livetime', + 'livetime', 'aeff_scale' ) - super(self.__class__, self).__init__( + super().__init__( expected_params=expected_params, **std_kwargs, ) - - energy_param_source = self.params.aeff_energy_paramfile.value - coszen_param_source = self.params.aeff_coszen_paramfile.value - - self.energy_param = load_aeff_param(energy_param_source) - self.coszen_param = load_aeff_param(coszen_param_source) + + self.energy_param = load_aeff_param(self.params.aeff_energy_paramfile.value) + self.coszen_param = load_aeff_param(self.params.aeff_coszen_paramfile.value) @profile def apply_function(self): - # read out aeff_scale = self.params.aeff_scale.m_as('dimensionless') livetime_s = self.params.livetime.m_as('sec') - + for container in self.data: scale = aeff_scale * livetime_s * np.ones(container.size) if container.name in self.energy_param.keys(): From 488240aa5c3630bb175aa801303033d1308fdf13 Mon Sep 17 00:00:00 2001 From: Jan Weldert Date: Fri, 9 Aug 2024 07:08:10 -0500 Subject: [PATCH 4/6] Include 1d interpolation in aeff param stage --- pisa/stages/aeff/param.py | 54 ++++++++++++++++++++++++++++++++++----- 1 file changed, 48 insertions(+), 6 deletions(-) diff --git a/pisa/stages/aeff/param.py b/pisa/stages/aeff/param.py index e56cedaca..d946b91a1 100644 --- a/pisa/stages/aeff/param.py +++ b/pisa/stages/aeff/param.py @@ -13,6 +13,7 @@ from collections.abc import Mapping import numpy as np +from scipy.interpolate import interp1d from pisa.core.stage import Stage from pisa.utils.fileio import from_file @@ -63,19 +64,61 @@ def load_aeff_param(source): `val`s can be one of the following: - Callable with one argument - String such that `eval(val)` yields a callable with one argument + - Mapping with the format: + { + <"energy" or "coszen">: [sequence of values], + "aeff": [sequence of values] + } + the two sequences are used to form a linear interpolant callable that + maps energy or coszen values to aeff values. """ if not isinstance(source, (str, Mapping)): raise TypeError('`source` must be string or mapping') if isinstance(source, str): - orig_dict = from_file(source) + aeff_dict = from_file(source) elif isinstance(source, Mapping): - orig_dict = source + aeff_dict = source else: raise TypeError('Cannot load aeff parameterizations from a %s' % type(source)) - return orig_dict + for k in aeff_dict: + func = aeff_dict[k] + + if isinstance(func, str): + param_func = eval(func) + + elif callable(func): + param_func = func + + elif isinstance(func, Mapping): + is_energy = 'energy' in func + is_coszen = 'coszen' in func + + if 'aeff' not in func: + raise ValueError('No effective area values are provided for %s'%(k)) + if not (is_energy or is_coszen): + raise ValueError('No energy or coszen values are provided for %s'%(k)) + + var = 'energy' if is_energy else 'coszen' + x_vals = func[var] + aeff_vals = func['aeff'] + + param_func = interp1d(x_vals, aeff_vals, kind='linear', + bounds_error=False, fill_value=0) + + else: + raise TypeError( + 'Expected parameteriation spec to be either a string that' + ' can be interpreted by eval or as a mapping of values' + ' from which to construct a spline. Got "%s".' + % type(param_spec) + ) + + aeff_dict[k] = param_func + + return aeff_dict class param(Stage): # pylint: disable=invalid-name @@ -113,7 +156,6 @@ def __init__( self.energy_param = load_aeff_param(self.params.aeff_energy_paramfile.value) self.coszen_param = load_aeff_param(self.params.aeff_coszen_paramfile.value) - @profile def apply_function(self): aeff_scale = self.params.aeff_scale.m_as('dimensionless') livetime_s = self.params.livetime.m_as('sec') @@ -121,10 +163,10 @@ def apply_function(self): for container in self.data: scale = aeff_scale * livetime_s * np.ones(container.size) if container.name in self.energy_param.keys(): - func = eval(self.energy_param[container.name]) + func = self.energy_param[container.name] scale *= func(container['true_energy']) if container.name in self.coszen_param.keys(): - func = eval(self.coszen_param[container.name]) + func = self.coszen_param[container.name] scale *= func(container['true_coszen']) container['weights'] *= scale From 0f6fca7952ab47a350cdf10a09d6fc4e76ab06de Mon Sep 17 00:00:00 2001 From: Jan Weldert Date: Fri, 9 Aug 2024 07:48:10 -0500 Subject: [PATCH 5/6] Forgot to rename --- pisa/stages/aeff/param.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/pisa/stages/aeff/param.py b/pisa/stages/aeff/param.py index d946b91a1..e2a9c733b 100644 --- a/pisa/stages/aeff/param.py +++ b/pisa/stages/aeff/param.py @@ -17,7 +17,6 @@ from pisa.core.stage import Stage from pisa.utils.fileio import from_file -from pisa.utils.profiler import profile __all__ = ['load_aeff_param', 'param'] @@ -110,10 +109,10 @@ def load_aeff_param(source): else: raise TypeError( - 'Expected parameteriation spec to be either a string that' + 'Expected parameteriation to be either a string that' ' can be interpreted by eval or as a mapping of values' ' from which to construct a spline. Got "%s".' - % type(param_spec) + % type(func) ) aeff_dict[k] = param_func From 731f73d822145743efb54e504772a94c32ab35b5 Mon Sep 17 00:00:00 2001 From: Jan Weldert Date: Mon, 12 Aug 2024 07:24:01 -0500 Subject: [PATCH 6/6] Expand interpolation description --- pisa/stages/aeff/param.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/pisa/stages/aeff/param.py b/pisa/stages/aeff/param.py index e2a9c733b..3b38480b6 100644 --- a/pisa/stages/aeff/param.py +++ b/pisa/stages/aeff/param.py @@ -69,7 +69,9 @@ def load_aeff_param(source): "aeff": [sequence of values] } the two sequences are used to form a linear interpolant callable that - maps energy or coszen values to aeff values. + maps energy or coszen values to aeff values. The effective area for any + energy or coszen outside the bounds of the corresponding sequence is + assumed to be 0. """ if not isinstance(source, (str, Mapping)): raise TypeError('`source` must be string or mapping')