Add "like_accuracy" parameter to likelihood tools and scripts (#2662)

Jürgen Knödlseder
1 parent 9a08012a
Showing 21 changed files with 99 additions and 50 deletions
ChangeLog
NEWS
cscripts/cstsmapsplit.par
doc/source/admin/release_history/1.6.rst
doc/source/users/reference_manual/cstsmapsplit.rst
doc/source/users/reference_manual/ctbutterfly.rst
doc/source/users/reference_manual/cterror.rst
doc/source/users/reference_manual/ctlike.rst
doc/source/users/reference_manual/cttsmap.rst
doc/source/users/reference_manual/ctulimit.rst
src/ctbutterfly/ctbutterfly.cpp
src/ctbutterfly/ctbutterfly.par
src/cterror/cterror.cpp
src/cterror/cterror.par
src/ctlike/ctlike.cpp
src/ctlike/ctlike.hpp
src/ctlike/ctlike.par
src/cttsmap/cttsmap.par
src/ctulimit/ctulimit.cpp
src/ctulimit/ctulimit.par
-2018-09-10
+2018-09-11
         * Version 1.6.0 released
         ========================
+        Add "like_accuracy" parameter to likelihood tools and scripts (#2662)
         Do not skip bins with zero counts in csspec (#2671)
         Use "ENERG" instead of "ETRUE" in energy dispersion axis names in csroot2caldb (#2670)
         Add hidden "use_irf_bkg" parameter to csphagen (#2659)
 New Features and Important Changes in ctools 1.6.0
-10 September 2018
+11 September 2018
 Introduction
@@ -61,7 +61,8 @@ None
 ctbutterfly - Generate butterfly diagrams
 -----------------------------------------
-None
+Add hidden "like_accuracy" parameter to allow control of the accuracy of
+the maximum log-likelikelihood value (#2662).
 ctcubemask - CTA cube filtering
@@ -76,6 +77,9 @@ None
 cterror - Likelihood profile error computation
 ----------------------------------------------
+Add hidden "like_accuracy" parameter to allow control of the accuracy of
+the maximum log-likelikelihood value (#2662).
+
 Add hidden "nthreads" parameter to control parallel computations (#2421).
@@ -86,6 +90,9 @@ None
 ctlike - Maximum likelihood fitting
 -----------------------------------
+Add hidden "like_accuracy" parameter to allow control of the accuracy of
+the maximum log-likelikelihood value (#2662).
+
 Add hidden "nthreads" parameter to control parallel computations (#2421).
@@ -150,11 +157,15 @@ accurately computes the distances between sky map pixels (#2309).
 cttsmap - Generation of Test Statistic map
 ------------------------------------------
-None
+Add hidden "like_accuracy" parameter to allow control of the accuracy of
+the maximum log-likelikelihood value (#2662).
 ctulimit - Compute upper limit
 ------------------------------
+Add hidden "like_accuracy" parameter to allow control of the accuracy of
+the maximum log-likelikelihood value (#2662).
+
 Add "parname" parameter to allow upper limit computations also for a node
 function spectral model. In addition, it is now checked that the parameter
 values actually bracket the upper limit (#2233).
@@ -356,7 +367,8 @@ None
 cstsmapsplit - Prepare split Test Statistic computation
 -------------------------------------------------------
-None
+Add hidden "like_accuracy" parameter to allow control of the accuracy of
+the maximum log-likelikelihood value (#2662).
 csviscube - Generate visibility cube
@@ -28,6 +28,7 @@ irf,           s, a, South_0.5h,,, &quot;Instrument response function&quot;
 edisp,         b, h, no,,, "Apply energy dispersion?"
 errors,        b, h, no,,, "Compute/store parameter errors?"
 statistic,     s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
+like_accuracy, r, h, 0.005,1.0e-30,0.1, "Absolute accuracy of maximum likelihood value"
 outmap,        f, a, tsmap.fits,,, "Output Test Statistic map file"
 bins_per_job,  i, a, 5,,, "Number of TS map bins per task"
 compute_null,  b, a, yes,,,"Compute null hypothesis first?"
@@ -37,6 +37,10 @@ Improvements
 * [`2670 <https://cta-redmine.irap.omp.eu/issues/2670>`_] -
   Use ``ENERG`` instead of ``ETRUE`` in energy dispersion axis names in
   :ref:`csroot2caldb`
+* [`2662 <https://cta-redmine.irap.omp.eu/issues/2662>`_] -
+  Add ``like_accuracy`` parameter to :ref:`ctbutterfly`, :ref:`cterror`,
+  :ref:`ctlike`, :ref:`cttsmap`, :ref:`ctulimit`, and :ref:`cstsmapsplit`
+  to allow controlling the accuracy of the maximum likelihood computation
 * [`2659 <https://cta-redmine.irap.omp.eu/issues/2659>`_] -
   Add hidden ``use_irf_bkg`` parameter to :ref:`csphagen`
 * [`2656 <https://cta-redmine.irap.omp.eu/issues/2656>`_] -
@@ -61,6 +61,12 @@ General parameters
     observations, and ``CHI2`` uses the Chi squared statistic for all
     binned or stacked observations.
+``(like_accuracy = 0.005) [real]``
+    Absolute accuracy of maximum likelihood value. Reducing this value will
+    increase the number of iterations and provide a more accurate maximum
+    log likelihood value. Converserly, decreasing the value will result in less
+    iterations at the expense of a less accurate maximum likelihood value.
+
 ``outmap [file]``
     Output Test Statistic map file.
@@ -94,6 +94,12 @@ General parameters
     observations, and ``CHI2`` uses the Chi squared statistic for all
     binned or stacked observations.
+``(like_accuracy = 0.005) [real]``
+    Absolute accuracy of maximum likelihood value. Reducing this value will
+    increase the number of iterations and provide a more accurate maximum
+    log likelihood value. Converserly, decreasing the value will result in less
+    iterations at the expense of a less accurate maximum likelihood value.
+
 ``(matrix = "NONE") [file]``
     Input covariance matrix file (not used)
@@ -75,9 +75,15 @@ General parameters
     observations, and ``CHI2`` uses the Chi squared statistic for all
     binned or stacked observations.
+``(like_accuracy = 0.005) [real]``
+    Absolute accuracy of maximum likelihood value. Reducing this value will
+    increase the number of iterations and provide a more accurate maximum
+    log likelihood value. Converserly, decreasing the value will result in less
+    iterations at the expense of a less accurate maximum likelihood value.
+
 ``(max_iter = 50) [integer]``
-    Maximum number of iterations before stopping the likelihood
-    profil computations.
+    Maximum number of iterations before stopping the likelihood profile
+    computations.
 Standard parameters
@@ -105,6 +105,12 @@ General parameters
 ``(refit = no) [boolean]``
     Perform refitting of solution after initial fit.
+``(like_accuracy = 0.005) [real]``
+    Absolute accuracy of maximum likelihood value. Reducing this value will
+    increase the number of iterations and provide a more accurate maximum
+    log likelihood value. Converserly, decreasing the value will result in less
+    iterations at the expense of a less accurate maximum likelihood value.
+
 ``(fix_spat_for_ts = no) [boolean]``
     Fix spatial parameters for TS computation.
@@ -74,6 +74,12 @@ General parameters
     observations, and ``CHI2`` uses the Chi squared statistic for all
     binned or stacked observations.
+``(like_accuracy = 0.005) [real]``
+    Absolute accuracy of maximum likelihood value. Reducing this value will
+    increase the number of iterations and provide a more accurate maximum
+    log likelihood value. Converserly, decreasing the value will result in less
+    iterations at the expense of a less accurate maximum likelihood value.
+
 ``(usepnt = no) [boolean]``
     Use CTA pointing direction for map centre instead of xref/yref parameters?
@@ -87,6 +87,12 @@ General parameters
     observations, and ``CHI2`` uses the Chi squared statistic for all
     binned or stacked observations.
+``(like_accuracy = 0.005) [real]``
+    Absolute accuracy of maximum likelihood value. Reducing this value will
+    increase the number of iterations and provide a more accurate maximum
+    log likelihood value. Converserly, decreasing the value will result in less
+    iterations at the expense of a less accurate maximum likelihood value.
+
 ``(max_iter = 50) [integer]``
     Maximum number of iterations before stopping the upper
     limit computations.
@@ -374,7 +374,6 @@ void ctbutterfly::init_members(void)
     m_srcname.clear();
     m_method.clear();
     m_confidence  = 0.68;
-    m_max_iter    = 50;
     m_apply_edisp = false;
     m_fit         = false;
     m_chatter     = static_cast<GChatter>(2);
@@ -389,10 +388,6 @@ void ctbutterfly::init_members(void)
     m_max_intensities.clear();
     m_csv.clear();
-    // Set optimizer parameters
-    m_opt.max_iter(m_max_iter);
-    m_opt.max_stalls(10);
-
     // Return
     return;
 }
@@ -32,10 +32,11 @@ method,     s, h, GAUSSIAN, GAUSSIAN|ENVELOPE,,&quot;Butterfly computation method&quot;
 #
 # Butterfly parameters
 #=====================
-fit,        b, h, no,,, "Perform maximum likelihood fitting?"
-confidence, r, h, 0.68,0.0,1.0, "Confidence level"
-statistic,  s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
-matrix,     f, h, NONE,,, "Input covariance Matrix FITS file"
+fit,           b, h, no,,, "Perform maximum likelihood fitting?"
+confidence,    r, h, 0.68,0.0,1.0, "Confidence level"
+statistic,     s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
+like_accuracy, r, h, 0.005,1.0e-30,0.1, "Absolute accuracy of maximum likelihood value"
+matrix,        f, h, NONE,,, "Input covariance Matrix FITS file"
 #
 # Energy binning parameters
@@ -399,10 +399,6 @@ void cterror::init_members(void)
     m_best_logL   = 0.0;
     m_model_par   = NULL;
-    // Set optimizer parameters
-    m_opt.max_iter(m_max_iter);
-    m_opt.max_stalls(10);
-
     // Return
     return;
 }
@@ -32,10 +32,11 @@ outmodel, f, a, results.xml,,, &quot;Output model definiton XML file&quot;
 #
 # Error calculation parameters
 #=============================
-confidence, r, h, 0.68,0.0,1.0, "Confidence level"
-tol,        r, h, 1e-3,,, "Computation tolerance"
-statistic,  s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
-max_iter,   i, h, 50,1,1000, "Maximum number of iterations"
+confidence,    r, h, 0.68,0.0,1.0, "Confidence level"
+tol,           r, h, 1e-3,,, "Computation tolerance"
+statistic,     s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
+like_accuracy, r, h, 0.005,1.0e-30,0.1, "Absolute accuracy of maximum likelihood value"
+max_iter,      i, h, 50,1,1000, "Maximum number of iterations"
 #
 # Standard parameters
@@ -559,8 +559,6 @@ void ctlike::init_members(void)
     m_outmodel.clear();
     m_outcovmat.clear();
     m_refit           = false;
-    m_max_iter        = 100;   // Set maximum number of iterations
-    m_max_stall       = 10;    // Set maximum number of stalls
     m_logL            = 0.0;
     m_apply_edisp     = false;
     m_fix_spat_for_ts = false;
@@ -569,10 +567,6 @@ void ctlike::init_members(void)
     // Set logger properties
     log.date(true);
-    // Set optimizer parameters
-    m_opt.max_iter(m_max_iter);
-    m_opt.max_stalls(m_max_stall);
-
     // Return
     return;
 }
@@ -589,8 +583,6 @@ void ctlike::copy_members(const ctlike&amp; app)
     m_refit           = app.m_refit;
     m_outmodel        = app.m_outmodel;
     m_outcovmat       = app.m_outcovmat;
-    m_max_iter        = app.m_max_iter;
-    m_max_stall       = app.m_max_stall;
     m_logL            = app.m_logL;
     m_apply_edisp     = app.m_apply_edisp;
     m_fix_spat_for_ts = app.m_fix_spat_for_ts;
 /***************************************************************************
  *                ctlike - Maximum likelihood fitting tool                 *
  * ----------------------------------------------------------------------- *
- *  copyright (C) 2010-2017 by Juergen Knoedlseder                         *
+ *  copyright (C) 2010-2018 by Juergen Knoedlseder                         *
  * ----------------------------------------------------------------------- *
  *                                                                         *
  *  This program is free software: you can redistribute it and/or modify   *
@@ -75,8 +75,6 @@ protected:
     GChatter  m_chatter;         //!< Chattiness
     // Members
-    int       m_max_iter;        //!< Maximum number of iterations
-    int       m_max_stall;       //!< Maximum number of stalls
     double    m_logL;            //!< Maximum log likelihood
 };
@@ -33,6 +33,7 @@ outcovmat, f, h, NONE,,, &quot;Output covariance matrix FITS file&quot;
 #========================
 statistic,       s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
 refit,           b, h, no,,, "Do refitting?"
+like_accuracy,   r, h, 0.005,1.0e-30,0.1, "Absolute accuracy of maximum likelihood value"
 fix_spat_for_ts, b, h, no,,, "Fix spatial parameters for TS computation?"
 #
@@ -31,8 +31,9 @@ outmap,    f, a, tsmap.fits,,, &quot;Output Test Statistic map file&quot;
 #
 # TS map computation parameters
 #==============================
-errors,    b, h, no,,, "Compute and store parameter errors?"
-statistic, s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
+errors,        b, h, no,,, "Compute and store parameter errors?"
+statistic,     s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
+like_accuracy, r, h, 0.005,1.0e-30,0.1, "Absolute accuracy of maximum likelihood value"
 #
 # Spatial binning parameters
@@ -372,10 +372,6 @@ void ctulimit::init_members(void)
     m_diff_ulimit  = 0.0;
     m_eflux_ulimit = 0.0;
-    // Set optimizer parameters
-    m_opt.max_iter(m_max_iter);
-    m_opt.max_stalls(10);
-
     // Return
     return;
 }
@@ -31,15 +31,16 @@ edisp,    b, h, no,,, &quot;Apply energy dispersion?&quot;
 #
 # Upper limit calculation parameters
 #===================================
-confidence, r, h, 0.95,0.0,1.0, "Confidence level"
-sigma_min,  r, h, 0.0,,, "Start minimum value (multiple sigmas above fit value)"
-sigma_max,  r, h, 10.0,,, "Start maximum value (multiple sigmas above fit value)"
-eref,       r, h, 1.0,,, "Reference energy for differential limit (TeV)"
-emin,       r, h, 1.0,,, "Minimum energy for integral flux limit (TeV)"
-emax,       r, h, 100.0,,, "Maximum energy for integral flux limit (TeV)"
-tol,        r, h, 1e-5,,, "Computation tolerance"
-statistic,  s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
-max_iter,   i, h, 50,1,1000, "Maximum number of iterations"
+confidence,    r, h, 0.95,0.0,1.0, "Confidence level"
+sigma_min,     r, h, 0.0,,, "Start minimum value (multiple sigmas above fit value)"
+sigma_max,     r, h, 10.0,,, "Start maximum value (multiple sigmas above fit value)"
+eref,          r, h, 1.0,,, "Reference energy for differential limit (TeV)"
+emin,          r, h, 1.0,,, "Minimum energy for integral flux limit (TeV)"
+emax,          r, h, 100.0,,, "Maximum energy for integral flux limit (TeV)"
+tol,           r, h, 1e-5,,, "Computation tolerance"
+statistic,     s, h, DEFAULT,DEFAULT|CSTAT|WSTAT|CHI2,, "Optimization statistic"
+like_accuracy, r, h, 0.005,1.0e-30,0.1, "Absolute accuracy of maximum likelihood value"
+max_iter,      i, h, 50,1,1000, "Maximum number of iterations"
 #
 # Standard parameters
@@ -290,6 +290,19 @@ void ctlikelihood::init_members(void)
     // Initialise members
     m_opt.clear();
+    // Initialise optimizer characteristics
+    m_opt.max_iter(50);     // Maximum number of iterations
+    m_opt.max_stalls(10);   // Maximum number of stalls
+    m_opt.eps(0.005);       // Accuracy of maximum likelihood value
+
+    // Set optimizer characteristics from optional user parameters
+    if (has_par("like_accuracy")) {
+        m_opt.eps((*this)["like_accuracy"].real());
+    }
+    if (has_par("max_iter")) {
+        m_opt.max_iter((*this)["max_iter"].integer());
+    }
+
     // Return
     return;
 }