Allow using spectral models without a 'Prefactor' parameter in cssens and cstsdist (#1889)

Jürgen Knödlseder
1 parent 0a9552e7
Showing 6 changed files with 88 additions and 22 deletions
ChangeLog
NEWS
cscripts/cssens.py
cscripts/cstsdist.py
cscripts/modutils.py
doc/source/admin/release_history/1.7.rst
-2019-07-23
+2019-10-31
  
         * Version 1.7.0 released
         ========================
  
         Remove response cache before source testing in cttsmap (#2973)
         Enabling running of example scripts from objects in memory (#2925)
+        Allow using spectral models without "Prefactor" parameter in cssens
+          and cstsdist (#1889)
  
  
 2019-07-15
 New Features and Important Changes in ctools 1.7.0
  
-23 July 2019
+31 October 2019
  
  
 Introduction
@@ -254,7 +254,8 @@ None
  
 cssens - Computes sensitivity curve
 -----------------------------------
-None
+Allow using spectral models without "Prefactor" parameter for sensitivity
+computation (#1889).
  
  
 csspec - Generates Spectral Energy Distribution
@@ -269,7 +270,8 @@ None
  
 cstsdist - Computes Test Statistic distribution
 -----------------------------------------------
-None
+Allow using spectral models without "Prefactor" parameter for computation of
+TS distribution (#1889).
  
  
 cstsmapmerge - Merge Test Statistic maps
@@ -304,7 +306,9 @@ None
  
 modutils - Model utilities
 --------------------------
-None
+Allow using spectral models without "Prefactor" parameter as test source in
+test_source() function. Add normalisation_parameter() function to retrieve
+spectral normalisation parameter from any model (#1889).
  
  
 mputils - Multiprocessing utilities
@@ -2,7 +2,7 @@
 # ==========================================================================
 # Computes the array sensitivity using the Test Statistic for a test source
 #
-# Copyright (C) 2011-2018 Juergen Knoedlseder
+# Copyright (C) 2011-2019 Juergen Knoedlseder
 #
 # This program is free software: you can redistribute it and/or modify
 # it under the terms of the GNU General Public License as published by
@@ -377,13 +377,14 @@ class cssens(ctools.csobservation):
             # Write header for iteration into logger
             self._log_header2(gammalib.EXPLICIT, 'Iteration '+str(iterations))
  
-            # Create a copy of the test models, set the prefactor of the test
-            # source in the models, and append the models to the observation.
-            # "crab_prefactor" is the Prefactor that corresponds to a flux of
-            # 1 Crab.
+            # Create a copy of the test models, set the normalisation parameter
+            # of the test source in the models, and append the models to the
+            # observation. "crab_prefactor" is the prefactor that corresponds
+            # to a flux of 1 Crab.
             models         = test_model.copy()
-            crab_prefactor = models[self._srcname]['Prefactor'].value() * crab_unit
-            models[self._srcname]['Prefactor'].value(crab_prefactor * test_crab_flux)
+            prefactor      = modutils.normalisation_parameter(models[self._srcname])
+            crab_prefactor = prefactor.value() * crab_unit
+            prefactor.value(crab_prefactor * test_crab_flux)
             self.obs().models(models)
  
             # Simulate events for the models. "sim" holds an observation
@@ -422,7 +423,8 @@ class cssens(ctools.csobservation):
             ts     = source.ts()
  
             # Get fitted Crab, photon and energy fluxes
-            crab_flux   = source['Prefactor'].value() / crab_prefactor
+            prefactor   = modutils.normalisation_parameter(source)
+            crab_flux   = prefactor.value() / crab_prefactor
             photon_flux = source.spectral().flux(emin, emax)
             energy_flux = source.spectral().eflux(emin, emax)
  
@@ -679,7 +681,7 @@ class cssens(ctools.csobservation):
  
         # Set test source model for this observation
         self._models = modutils.test_source(self.obs().models(), self._srcname,
-                                      ra=self._ra, dec=self._dec)
+                                            ra=self._ra, dec=self._dec)
  
         # Write observation into logger
         self._log_observations(gammalib.NORMAL, self.obs(), 'Input observation')
@@ -232,9 +232,10 @@ class cstsdist(ctools.csobservation):
                 for par in model:
                     self._log_string(gammalib.EXPLICIT, str(par))
         elif self._logNormal():
-            name  = 'Trial %d' % seed
-            value = 'TS=%.3f  Prefactor=%e +/- %e' % \
-                    (ts, model['Prefactor'].value(), model['Prefactor'].error())
+            prefactor = modutils.normalisation_parameter(model)
+            name      = 'Trial %d' % seed
+            value     = 'TS=%.3f  %s=%e +/- %e' % \
+                        (ts, prefactor.name(), prefactor.value(), prefactor.error())
             self._log_value(gammalib.TERSE, name, value)
  
         # Initialise results
 # ==========================================================================
 # Utility functions for model handling
 #
-# Copyright (C) 2016-2017 Juergen Knoedlseder
+# Copyright (C) 2016-2019 Juergen Knoedlseder
 #
 # This program is free software: you can redistribute it and/or modify
 # it under the terms of the GNU General Public License as published by
@@ -61,10 +61,20 @@ def test_source(models, srcname, ra=None, dec=None, fitspat=False,
     model = outmodels[srcname]
     model.tscalc(True)
  
-    # If source model has no "Prefactor" parameter then raise an exception
-    if not model.has_par('Prefactor'):
-        msg = ('Model "%s" has no parameter "Prefactor". Only spectral '
-               'models with a "Prefactor" parameter are supported.' %
+    # Get spectral model component
+    spectral = model.spectral()
+
+    # If source model has no normalisation parameter then raise an exception
+    if not (spectral.has_par('Normalization') or
+            spectral.has_par('Value') or
+            spectral.has_par('Prefactor') or
+            spectral.has_par('Integral') or
+            spectral.has_par('PhotonFlux') or
+            spectral.has_par('EnergyFlux')):
+        msg = ('Model "%s" has no flux normalisation parameter. Only spectral '
+               'models with a flux normalisation parameter (one of '
+               '"Normalization", "Value", "Prefactor", "Integral", '
+               '"PhotonFlux", and "EnergyFlux") are supported.' %
                srcname)
         raise RuntimeError(msg)
  
@@ -100,6 +110,51 @@ def test_source(models, srcname, ra=None, dec=None, fitspat=False,
     return outmodels
  
  
+# ============================================= #
+# Return spectral model normalisation parameter #
+# ============================================= #
+def normalisation_parameter(model):
+    """
+    Return spectral model normalisation parameter
+
+    Parameters
+    ----------
+    model : `~gammalib.GModel`
+        Input model
+
+    Returns
+    -------
+    par : `~gammalib.GModelPar`
+        Spectral model normalisation parameter
+    """
+    # Initialise normalisation parameter as None
+    par = None
+
+    # Set list of possible normalisation parameters
+    norm_pars = ['Normalization', 'Value', 'Prefactor', 'Integral',
+                 'PhotonFlux', 'EnergyFlux']
+
+    # Get spectral model component
+    spectral = model.spectral()
+
+    # Find normalisation parameter
+    for norm_par in norm_pars:
+        if spectral.has_par(norm_par):
+            par = spectral[norm_par]
+
+    # If no normalisation parameter then raise an exception
+    if par == None:
+        msg = ('Model "%s" has no flux normalisation parameter. Only spectral '
+               'models with a flux normalisation parameter (one of '
+               '"Normalization", "Value", "Prefactor", "Integral", '
+               '"PhotonFlux", and "EnergyFlux") are supported.' %
+               model.name())
+        raise RuntimeError(msg)
+
+    # Return normalisation parameter
+    return par
+
+
 # ================================= #
 # Returns attributes for ds9 string #
 # ================================= #
@@ -26,6 +26,8 @@ Improvements
   Remove response cache before source testing in cttsmap
 * [`2925 <https://cta-redmine.irap.omp.eu/issues/2925>`_] -
   Enabling running of example scripts from objects in memory
+* [`1889 <https://cta-redmine.irap.omp.eu/issues/1889>`_] -
+  Allow using spectral models without ``Prefactor`` parameter in :ref:`cssens` and :ref:`cstsdist`
  
  
 New features