Add ctbutterfly::butterfly() method (#1512)

Jürgen Knödlseder
1 parent 9ecb561a
Showing 7 changed files with 72 additions and 36 deletions
ChangeLog
NEWS
doc/source/admin/release_history/1.5.rst
pyext/ctbutterfly.i
src/ctbutterfly/ctbutterfly.cpp
src/ctbutterfly/ctbutterfly.hpp
src/ctlike/ctlike.cpp
-2018-01-22
+2018-01-23
         * Version 1.5.0 released
         ========================
@@ -7,7 +7,8 @@
           DATE-END, TIME-OBS and TIME-END
         Add CUSTOM background method to csphagen script (#2284)
         Add csresspec script (#2289)
-        Add ctools::time_reference constant and remove ctools::time_reference() method (#2295)
+        Add ctools::time_reference constant and remove ctools::time_reference()
+          method (#2295)
         Correctly strip .gz suffix from output filenames (#2287)
         Support On/Off analysis in cslightcrv and csphasecrv (#2273)
         ctselect informs about non-overlapping time interval (#2272)
@@ -25,6 +26,7 @@
         Add "NODES" method to csspec to support non-CTA analysis (#1927)
         Use t-type parameters for all "tmin" and "tmax" parameters (#1864)
         Set OBS_ID in ctobssim (#1836)
+        Add ctbutterfly::butterfly() method (#1512)
 2017-07-28
-22 January 2018
+23 January 2018
 Introduction
@@ -42,6 +42,8 @@ ctbutterfly - Generate butterfly diagrams
 -----------------------------------------
 Add "statistic" parameter to allow setting of optimizer statistic (#2237).
+Add butterfly() method to return the butterfly diagram as CSV object (#1512).
+
 ctcubemask - CTA cube filtering
 -------------------------------
@@ -58,6 +58,8 @@ Improvements
   Add ``NODES`` method to :ref:`csspec` to support non-CTA analysis
 * [`1836 <https://cta-redmine.irap.omp.eu/issues/1836>`_] -
   Set ``OBS_ID`` in :ref:`ctobssim`
+* [`1512 <https://cta-redmine.irap.omp.eu/issues/1512>`_] -
+  Add ``ctbutterfly::butterfly()`` method to return butterfly diagram
 New features
 /***************************************************************************
  *                  ctbutterfly - butterfly calculation tool               *
  * ----------------------------------------------------------------------- *
- *  copyright (C) 2014 by Michael Mayer                                    *
+ *  copyright (C) 2014-2018 by Michael Mayer                               *
  * ----------------------------------------------------------------------- *
  *                                                                         *
  *  This program is free software: you can redistribute it and/or modify   *
@@ -44,9 +44,10 @@ public:
     virtual ~ctbutterfly(void);
     // Methods
-    void clear(void);
-    void run(void);
-    void save(void);
+    void        clear(void);
+    void        run(void);
+    void        save(void);
+    const GCsv& butterfly(void) const;
 };
 /***************************************************************************
  *                 ctbutterfly - butterfly calculation tool                *
  * ----------------------------------------------------------------------- *
- *  copyright (C) 2014-2017 by Michael Mayer                               *
+ *  copyright (C) 2014-2018 by Michael Mayer                               *
  * ----------------------------------------------------------------------- *
  *                                                                         *
  *  This program is free software: you can redistribute it and/or modify   *
@@ -299,6 +299,9 @@ void ctbutterfly::run(void)
         ellipsoid_boundary(models);
     }
+    // Set CSV table
+    set_csv();
+
     // Write header into logger
     log_header1(TERSE, "Results");
@@ -311,7 +314,7 @@ void ctbutterfly::run(void)
                 gammalib::str(m_min_intensities[k])+", "+
                 gammalib::str(m_max_intensities[k])+") ph/cm2/s/MeV");
-    } //endfor: loop over energy bin  
+    } // endfor: loop over energy bin
     // Restore energy dispersion flag for all CTA observations
     for (int i = 0; i < m_obs.size(); ++i) {
@@ -329,7 +332,7 @@ void ctbutterfly::run(void)
 /***********************************************************************//**
  * @brief Save butterfly diagram
  *
- * Saves the butterfly diagram into an ascii file using a column separated
+ * Saves the butterfly diagram into an ASCII file using a column separated
  * value (CSV) format with blanks as separators.
  ***************************************************************************/
 void ctbutterfly::save(void)
@@ -346,19 +349,8 @@ void ctbutterfly::save(void)
         // Log butterfly diagram file name
         log_value(NORMAL, "Butterfly file", m_outfile.url());
-        // Create CSV table with 4 columns
-        GCsv table(m_energies.size(), 4);
-
-        // Fill CSV table
-        for (int i = 0; i < m_energies.size(); ++i) {
-            table.real(i, 0, m_energies[i]);
-            table.real(i, 1, m_intensities[i]);
-            table.real(i, 2, m_min_intensities[i]);
-            table.real(i, 3, m_max_intensities[i]);
-        }
-
         // Save CSV table
-        table.save(m_outfile.url(), " ", clobber());
+        m_csv.save(m_outfile.url(), " ", clobber());
     }
@@ -395,6 +387,7 @@ void ctbutterfly::init_members(void)
     m_intensities.clear();
     m_min_intensities.clear();
     m_max_intensities.clear();
+    m_csv.clear();
     // Set optimizer parameters
     m_opt.max_iter(m_max_iter);
@@ -429,6 +422,7 @@ void ctbutterfly::copy_members(const ctbutterfly&amp; app)
     m_intensities     = app.m_intensities;
     m_min_intensities = app.m_min_intensities;
     m_max_intensities = app.m_max_intensities;
+    m_csv             = app.m_csv;
     // Return
     return;
@@ -851,13 +845,13 @@ void ctbutterfly::check_model(void)
  * http://www.math.harvard.edu/archive/21b_fall_04/exhibits/2dmatrices/index.html
  ***************************************************************************/
 void ctbutterfly::eigenvectors(const double& a,
-        const double& b,
-        const double& c,
-        const double& d,
-        double*       lambda1,
-        double*       lambda2,
-        GVector*      vector1,
-        GVector*      vector2)
+                               const double& b,
+                               const double& c,
+                               const double& d,
+                               double*       lambda1,
+                               double*       lambda2,
+                               GVector*      vector1,
+                               GVector*      vector2)
 {
     // Compute trace and determinant
     double trace = a + d;
@@ -896,3 +890,26 @@ void ctbutterfly::eigenvectors(const double&amp; a,
     // Return
     return;
 }
+
+
+/***********************************************************************//**
+ * @brief Set CSV object
+ *
+ * Set CSV object with results.
+ ***************************************************************************/
+void ctbutterfly::set_csv(void)
+{
+    // Create CSV table with 4 columns
+    m_csv = GCsv(m_energies.size(), 4);
+
+    // Fill CSV table
+    for (int i = 0; i < m_energies.size(); ++i) {
+        m_csv.real(i, 0, m_energies[i]);
+        m_csv.real(i, 1, m_intensities[i]);
+        m_csv.real(i, 2, m_min_intensities[i]);
+        m_csv.real(i, 3, m_max_intensities[i]);
+    }
+
+    // Return
+    return;
+}
 /***************************************************************************
  *                ctbutterfly - butterfly calculation tool                 *
  * ----------------------------------------------------------------------- *
- *  copyright (C) 2014-2017 by Michael Mayer                               *
+ *  copyright (C) 2014-2018 by Michael Mayer                               *
  * ----------------------------------------------------------------------- *
  *                                                                         *
  *  This program is free software: you can redistribute it and/or modify   *
@@ -66,9 +66,10 @@ public:
     ctbutterfly& operator=(const ctbutterfly& app);
     // Methods
-    void clear(void);
-    void run(void);
-    void save(void);
+    void        clear(void);
+    void        run(void);
+    void        save(void);
+    const GCsv& butterfly(void) const;
 protected:
     // Protected methods
@@ -78,8 +79,6 @@ protected:
     void get_parameters(void);
     void gaussian_error_propagation(GModels& models);
     void ellipsoid_boundary(GModels& models);
-
-
     void check_model(void);
     void eigenvectors(const double& a,
                       const double& b,
@@ -89,6 +88,7 @@ protected:
                       double*       lambda2,
                       GVector*      vector1,
                       GVector*      vector2);
+    void set_csv(void);
     // User parameters
     std::string m_srcname;      //!< Name of source to compute butterfly
@@ -107,7 +107,19 @@ protected:
     std::vector<double> m_intensities;     //!< Power law intensity
     std::vector<double> m_min_intensities; //!< Minimum intensities
     std::vector<double> m_max_intensities; //!< Maximum intensities
+    GCsv                m_csv;             //!< CSV table filled by run() method
 };
+/***********************************************************************//**
+ * @brief Return butterfly
+ *
+ * @return Butterfly as CSV table.
+ ***************************************************************************/
+inline
+const GCsv& ctbutterfly::butterfly(void) const
+{
+    return (m_csv);
+}
+
 #endif /* CTBUTTERFLY_HPP */
@@ -360,7 +360,7 @@ void ctlike::save(void)
         // Log filename
         log_value(NORMAL, "Covariance matrix file", m_outcovmat.url());
-        // Save covarianve matrix
+        // Save covariance matrix
         m_obs.function().save(m_outcovmat);
     }