Merge branch 'devel' into release

# Conflicts: # ChangeLog # NEWS # README # inst/cta/src/GCTAEdisp2D.cpp

Merge branch 'devel' into release
# Conflicts: # ChangeLog # NEWS # README # inst/cta/src/GCTAEdisp2D.cpp
Jürgen Knödlseder
2 parents 89317305 bc79a4f9
Showing 16 changed files with 306 additions and 486 deletions
ChangeLog
NEWS
configure.ac
doc/Doxyfile
doc/source/conf.py
include/GModel.hpp
include/GVOHub.hpp
inst/cta/include/GCTAEventAtom.hpp
inst/cta/include/GCTAEventList.hpp
inst/cta/pyext/GCTAEventAtom.i
inst/cta/src/GCTAEventAtom.cpp
inst/cta/src/GCTAEventList.cpp
pyext/GModel.i
src/vo/GVOClient.cpp
src/vo/GVOHub.cpp
test/test_GVO.cpp
-2015-12-17  Juergen Knoedlseder  <jurgen.knodlseder@irap.omp.eu>
+2015-12-18  Juergen Knoedlseder  <jurgen.knodlseder@irap.omp.eu>
         * Version 1.0.0 released
         =========================
         Ignore errors for El Capitan post processing (#1605)
+        Read only mandatary columns in GCTAEventList (#1600)
+        Add GModel::has_ts() method (#1578)
         Remove extname arguments from classes (#1601)
         Add GFilename class (#1599)
         Check existence of files in CTA response loading (#1586)
         GSkyMap::__getitem__ throw exception if out of bound (#1581)
         Handle ~ in filenames (#1580)
+        Add GModel::has_ts() method (#1578)
         Perform __init__.py creation only in build step (#1569)
         Remove access by name methods from GSkyRegions (#1564)
         Resolve Python test problem on El Capitan (#1563)
@@ -58,7 +58,8 @@ None
 8. Model module
 ---------------
-None
+Added GModel::has_ts() method that signals if a TS value has been computed
+for a given model.
 9. Numerics module
@@ -134,9 +135,12 @@ None
 20. CTA interface
 -----------------
-The IRF loaders now check the existence of the file before loading, and now throw
+The IRF loaders now check the existence of the file before loading, and now throws
 appropriate error messages in case that the files do not exist (#1586).
+Only manadatory columns are read by GCTAEventsList, while any other columns will only
+be stored for back writing in a FITS file.
+
 21. COMPTEL interface
 ---------------------
@@ -29,7 +29,7 @@
 #############################################################################
 # Initialisation                                                            #
 #############################################################################
-AC_INIT([gammalib], [1.0.0], [jurgen.knodlseder@irap.omp.eu], [gammalib])
+AC_INIT([gammalib], [1.1.0.dev], [jurgen.knodlseder@irap.omp.eu], [gammalib])
 AC_CONFIG_SRCDIR([src/support/GTools.cpp])
 AC_CONFIG_MACRO_DIR([m4])
 AM_INIT_AUTOMAKE([gnu no-dependencies -Wno-portability subdir-objects])
@@ -32,7 +32,7 @@ PROJECT_NAME           = GammaLib
 # This could be handy for archiving the generated documentation or
 # if some version control system is used.
-PROJECT_NUMBER         = 1.0.0
+PROJECT_NUMBER         = 1.1.0.dev
 # Using the PROJECT_BRIEF tag one can provide an optional one line description
 # for a project that appears at the top of each page and should give viewer
@@ -55,9 +55,10 @@ copyright = u&#39;2006-2015, GammaLib developers&#39;
 # built documents.
 #
 # The short X.Y version.
-version = '1.0'
+version = '1.1'
+
 # The full version, including alpha/beta/rc tags.
-release = '1.0.0'
+release = '1.1.0.dev'
 # The language for content autogenerated by Sphinx. Refer to documentation
 # for a list of supported languages.
@@ -139,6 +139,7 @@ public:
     void                ts(const double& ts);
     const bool&         tscalc(void) const;
     void                tscalc(const bool& tscalc);
+    const bool&         has_ts(void) const;
     std::string         instruments(void) const;
     void                instruments(const std::string& instruments);
     GModelPar           scale(const std::string& instrument) const;
@@ -314,7 +315,7 @@ void GModel::tscalc(const bool&amp; tscalc)
 /***********************************************************************//**
  * @brief Signals that model has scales
  *
- * @return True is model has scale factors.
+ * @return True if model has scale factors.
  ***************************************************************************/
 inline
 bool GModel::has_scales(void) const
@@ -322,4 +323,16 @@ bool GModel::has_scales(void) const
     return (!m_scales.empty());
 }
+
+/***********************************************************************//**
+ * @brief Signals that model has Test Statistics value
+ *
+ * @return True if model has TS value.
+ ***************************************************************************/
+inline
+const bool& GModel::has_ts(void) const
+{
+    return m_has_ts;
+}
+
 #endif /* GMODEL_HPP */
@@ -108,11 +108,14 @@ protected:
     // Low-level methods
     void        create_samp_file(void) const;
+    void        delete_samp_file(void) const;
+    int         get_socket(void);
     void        post_samp_ok(const socklen_t& sock) const;
     void        post_samp_void(const socklen_t& sock) const;
     void	    post_string(const std::string& content, const socklen_t& sock) const;
     void        notify(const std::string& url, const std::string& notification) const;
     std::string	random_string(const size_t& length) const;
+    std::string hub_url(void) const;
     // Protected structure
     struct client {
@@ -132,7 +135,6 @@ protected:
     // Protected members
     std::string         m_secret;    //!< Secret Hub key
-    std::string         m_hub_url;   //!< The XML-RPC endpoint for communication with the hub
     std::string         m_hub_host;  //!< Hub host
     std::string         m_hub_port;  //!< Hub port
     std::string         m_hub_path;  //!< Hub path
@@ -72,10 +72,8 @@ public:
     // Other methods
     const int&           index(void) const;
     const unsigned long& event_id(void) const;
-    const unsigned long& obs_id(void) const;
     const float&         phase(void) const;
     void                 event_id(const unsigned long& id);
-    void                 obs_id(const unsigned long& id);
     void                 phase(const float& phase);
 protected:
@@ -90,24 +88,6 @@ protected:
     GEnergy       m_energy;         //!< Event energy
     GTime         m_time;           //!< Event time
     unsigned long m_event_id;       //!< Event identifier
-    unsigned long m_obs_id;         //!< Observation identifier
-    int           m_multip;         //!< Multiplicity
-    char          m_telmask;        //!< Telescope mask
-    float         m_dir_err;        //!< Error on event direction
-    float         m_alt;            //!< Event altitude
-    float         m_az;             //!< Event azimuth
-    float         m_corex;          //!< Position on ground (m)
-    float         m_corey;          //!< Position on ground (m)
-    float         m_core_err;       //!< Error on core reconstruction
-    float         m_xmax;           //!< Position of shower max (g/cm2)
-    float         m_xmax_err;       //!< Error on shower max (g/cm2)
-    float         m_shwidth;        //!< Shower width (m)
-    float         m_shlength;       //!< Shower length (m)
-    float         m_energy_err;     //!< Error on event energy (MeV)
-    float         m_hil_msw;        //!< Hillas MSW
-    float         m_hil_msw_err;    //!< Hillas MSW error
-    float         m_hil_msl;        //!< Hillas MSL
-    float         m_hil_msl_err;    //!< Hillas MSL error
     float         m_phase;          //!< Optional phase
 };
@@ -226,18 +206,6 @@ const unsigned long&amp; GCTAEventAtom::event_id(void) const
 /***********************************************************************//**
- * @brief Return observation identifier
- *
- * @return Observation identifier.
- ***************************************************************************/
-inline
-const unsigned long& GCTAEventAtom::obs_id(void) const
-{
-    return (m_obs_id);
-}
-
-
-/***********************************************************************//**
  * @brief Return event phase
  *
  * @return Event phase.
@@ -263,19 +231,6 @@ void GCTAEventAtom::event_id(const unsigned long&amp; id)
 /***********************************************************************//**
- * @brief Set observation identifier
- *
- * @param[in] id Observation identifier.
- ***************************************************************************/
-inline
-void GCTAEventAtom::obs_id(const unsigned long& id)
-{
-    m_obs_id = id;
-    return;
-}
-
-
-/***********************************************************************//**
  * @brief Set event phase
  *
  * @param[in] phase Event phase.
@@ -88,19 +88,18 @@ protected:
     void         free_members(void);
     virtual void set_energies(void) { return; }
     virtual void set_times(void) { return; }
-    void         read_events(const GFitsTable& hdu);
-    void         read_events_v0(const GFitsTable& hdu);
-    void         read_events_v1(const GFitsTable& hdu);
-    void         read_events_hillas(const GFitsTable& hdu);
-    void         write_events(GFitsBinTable& hdu) const;
+    void         read_events(const GFitsTable& table);
+    void         write_events(GFitsBinTable& table) const;
     void         write_ds_keys(GFitsHDU& hdu) const;
     int          irf_cache_init(const std::string& name) const;
     int          irf_cache_index(const std::string& name) const;
     // Protected members
-    GCTARoi                    m_roi;       //!< Region of interest
-    std::vector<GCTAEventAtom> m_events;    //!< Events
-    bool                       m_has_phase; //!< Signal presence of phase
+    GCTARoi                     m_roi;       //!< Region of interest
+    std::vector<GCTAEventAtom>  m_events;    //!< Events
+    std::vector<GFitsTableCol*> m_columns;   //!< Pointers to optional columns
+    bool                        m_has_phase; //!< Signal presence of phase
+    bool                        m_has_detxy; //!< Signal presence of detector coordinates
     // IRF cache for diffuse models
     mutable std::vector<std::string>          m_irf_names;  //!< Model names
@@ -55,10 +55,8 @@ public:
     // Other methods
     const int&           index(void) const;
     const unsigned long& event_id(void) const;
-    const unsigned long& obs_id(void) const;
     const float&         phase(void) const;
     void                 event_id(const unsigned long& id);
-    void                 obs_id(const unsigned long& id);
     void                 phase(const float& phase);
 };
@@ -207,27 +207,9 @@ void GCTAEventAtom::init_members(void)
     m_dir.clear();
     m_time.clear();
     m_energy.clear();
-    m_index       = 0;
-    m_event_id    = 0;
-    m_obs_id      = 0;
-    m_multip      = 0;
-    m_telmask     = 0;
-    m_dir_err     = 0.0;
-    m_alt         = 0.0;
-    m_az          = 0.0;
-    m_corex       = 0.0;
-    m_corey       = 0.0;
-    m_core_err    = 0.0;
-    m_xmax        = 0.0;
-    m_xmax_err    = 0.0;
-    m_shwidth     = 0.0;
-    m_shlength    = 0.0;
-    m_energy_err  = 0.0;
-    m_hil_msw     = 0.0;
-    m_hil_msw_err = 0.0;
-    m_hil_msl     = 0.0;
-    m_hil_msl_err = 0.0;
-    m_phase       = 0.0;
+    m_index    = 0;
+    m_event_id = 0;
+    m_phase    = 0.0;
     // Return
     return;
@@ -242,30 +224,12 @@ void GCTAEventAtom::init_members(void)
 void GCTAEventAtom::copy_members(const GCTAEventAtom& atom)
 {
     // Copy members
-    m_dir         = atom.m_dir;
-    m_time        = atom.m_time;
-    m_energy      = atom.m_energy;
-    m_index       = atom.m_index;
-    m_event_id    = atom.m_event_id;
-    m_obs_id      = atom.m_obs_id;
-    m_multip      = atom.m_multip;
-    m_telmask     = atom.m_telmask; 
-    m_dir_err     = atom.m_dir_err;
-    m_alt         = atom.m_alt;
-    m_az          = atom.m_az;
-    m_corex       = atom.m_corex;
-    m_corey       = atom.m_corey;
-    m_core_err    = atom.m_core_err;
-    m_xmax        = atom.m_xmax;
-    m_xmax_err    = atom.m_xmax_err;
-    m_shwidth     = atom.m_shwidth;
-    m_shlength    = atom.m_shlength;
-    m_energy_err  = atom.m_energy_err;
-    m_hil_msw     = atom.m_hil_msw;
-    m_hil_msw_err = atom.m_hil_msw_err;
-    m_hil_msl     = atom.m_hil_msl;
-    m_hil_msl_err = atom.m_hil_msl_err;
-    m_phase       = atom.m_phase;
+    m_dir      = atom.m_dir;
+    m_time     = atom.m_time;
+    m_energy   = atom.m_energy;
+    m_index    = atom.m_index;
+    m_event_id = atom.m_event_id;
+    m_phase    = atom.m_phase;
     // Return
     return;
@@ -572,7 +572,9 @@ void GCTAEventList::init_members(void)
     // Initialise members
     m_roi.clear();
     m_events.clear();
+    m_columns.clear();
     m_has_phase = false;
+    m_has_detxy = true;
     // Initialise cache
     m_irf_names.clear();
@@ -594,11 +596,18 @@ void GCTAEventList::copy_members(const GCTAEventList&amp; list)
     m_roi       = list.m_roi;
     m_events    = list.m_events;
     m_has_phase = list.m_has_phase;
+    m_has_detxy = list.m_has_detxy;
     // Copy cache
     m_irf_names  = list.m_irf_names;
     m_irf_values = list.m_irf_values;
+    // Copy column pointers
+    m_columns.clear();
+    for (int i = 0; i < list.m_columns.size(); ++i) {
+        m_columns.push_back((list.m_columns[i]->clone()));
+    }
+
     // Return
     return;
 }
@@ -609,6 +618,12 @@ void GCTAEventList::copy_members(const GCTAEventList&amp; list)
  ***************************************************************************/
 void GCTAEventList::free_members(void)
 {
+    // Free columns
+    for (int i = 0; i < m_columns.size(); ++i) {
+        delete m_columns[i];
+        m_columns[i] = NULL;
+    }
+
     // Return
     return;
 }
@@ -619,52 +634,26 @@ void GCTAEventList::free_members(void)
  *
  * @param[in] table FITS table.
  *
- * This method reads the CTA event list from a FITS table HDU into memory.
- * Depending on the columns existing in the file, it either selects v0 or
- * v1 of the event list reader.
- ***************************************************************************/
-void GCTAEventList::read_events(const GFitsTable& table)
-{
-    // Clear existing events
-    m_events.clear();
-
-    // Extract number of events in FITS file
-    int num = table.integer("NAXIS2");
-
-    // Continue only if there are events
-    if (num > 0) {
-
-        // Read events for v1
-        if (table.contains("SHWIDTH") && table.contains("SHLENGTH")) {
-            read_events_v1(table);
-        }
-
-        // ... otherwise read events for v0
-        else {
-            read_events_v0(table);
-        }
-
-        // Read (optional) Hillas parameters
-        read_events_hillas(table);
-
-    } // endif: there were events
-
-    // Return
-    return;
-}
-
-
-/***********************************************************************//**
- * @brief Read CTA events from FITS table (version 0)
+ * This method reads the CTA event list from a FITS table into memory.
+ * The following columns are mandatory in the FITS table:
  *
- * @param[in] table FITS table.
+ *      - EVENT_ID
+ *      - TIME
+ *      - RA
+ *      - DEC
+ *      - ENERGY
+ *
+ * Optionally, the
  *
- * This method reads the CTA event list from a FITS table HDU into memory.
- * It is a minimal event reader that is compliant with the initial data
- * format distributed by Karl Kosack. Information that is not present in
- * that format is set to 0.
+ *      - DETX
+ *      - DETY
+ *      - PHASE
+ *
+ * columns are read. Any other columns present in the FITS table will also be
+ * read, and written into the FITS table when the write_events() method is
+ * called.
  ***************************************************************************/
-void GCTAEventList::read_events_v0(const GFitsTable& table)
+void GCTAEventList::read_events(const GFitsTable& table)
 {
     // Clear existing events
     m_events.clear();
@@ -679,120 +668,25 @@ void GCTAEventList::read_events_v0(const GFitsTable&amp; table)
         m_events.reserve(num);
         // Get column pointers
-        const GFitsTableCol* ptr_eid         = table["EVENT_ID"];
-        const GFitsTableCol* ptr_time        = table["TIME"];
-        const GFitsTableCol* ptr_multip      = table["MULTIP"];
-        const GFitsTableCol* ptr_ra          = table["RA"];
-        const GFitsTableCol* ptr_dec         = table["DEC"];
-        const GFitsTableCol* ptr_dir_err     = table["DIR_ERR"];
-        const GFitsTableCol* ptr_detx        = table["DETX"];
-        const GFitsTableCol* ptr_dety        = table["DETY"];
-        const GFitsTableCol* ptr_alt         = table["ALT"];
-        const GFitsTableCol* ptr_az          = table["AZ"];
-        const GFitsTableCol* ptr_corex       = table["COREX"];
-        const GFitsTableCol* ptr_corey       = table["COREY"];
-        const GFitsTableCol* ptr_core_err    = table["CORE_ERR"];
-        const GFitsTableCol* ptr_xmax        = table["XMAX"];
-        const GFitsTableCol* ptr_xmax_err    = table["XMAX_ERR"];
-        const GFitsTableCol* ptr_energy      = table["ENERGY"];
-        const GFitsTableCol* ptr_energy_err  = table["ENERGY_ERR"];
-
-        // Check for phase column
-        const GFitsTableCol* ptr_phase;
-        if (table.contains("PHASE")) {
-            m_has_phase = true;
-            ptr_phase   = table["PHASE"];
+        const GFitsTableCol* ptr_eid    = table["EVENT_ID"];
+        const GFitsTableCol* ptr_time   = table["TIME"];
+        const GFitsTableCol* ptr_ra     = table["RA"];
+        const GFitsTableCol* ptr_dec    = table["DEC"];
+        const GFitsTableCol* ptr_energy = table["ENERGY"];
+
+        // Check for DETX and DETY columns
+        const GFitsTableCol* ptr_detx;
+        const GFitsTableCol* ptr_dety;
+        if (table.contains("DETX") && table.contains("DETY")) {
+            m_has_detxy = true;
+            ptr_detx    = table["DETX"];
+            ptr_dety    = table["DETY"];
         }
         else {
-            m_has_phase = false;
-        }
-
-        // Copy data from columns into GCTAEventAtom objects
-        GCTAEventAtom event;
-        for (int i = 0; i < num; ++i) {
-            event.m_index     = i;
-            event.m_time.set(ptr_time->real(i), m_gti.reference());
-            event.m_dir.dir().radec_deg(ptr_ra->real(i), ptr_dec->real(i));
-            event.m_dir.detx(ptr_detx->real(i)*gammalib::deg2rad);
-            event.m_dir.dety(ptr_dety->real(i)*gammalib::deg2rad);
-            event.m_energy.TeV(ptr_energy->real(i));
-            event.m_event_id    = ptr_eid->integer(i);
-            event.m_obs_id      = 0;
-            event.m_multip      = ptr_multip->integer(i);
-            event.m_telmask     = 0;
-            event.m_dir_err     = ptr_dir_err->real(i);
-            event.m_alt         = ptr_alt->real(i);
-            event.m_az          = ptr_az->real(i);
-            event.m_corex       = ptr_corex->real(i);
-            event.m_corey       = ptr_corey->real(i);
-            event.m_core_err    = ptr_core_err->real(i);
-            event.m_xmax        = ptr_xmax->real(i);
-            event.m_xmax_err    = ptr_xmax_err->real(i);
-            event.m_shwidth     = 0.0;
-            event.m_shlength    = 0.0;
-            event.m_energy_err  = ptr_energy_err->real(i);
-
-            // Set pulse phase if available
-            if (m_has_phase) {
-                event.m_phase = ptr_phase->real(i);
-            }
-
-            m_events.push_back(event);
+            m_has_detxy = false;
         }
-    } // endif: there were events
-
-    // Return
-    return;
-}
-
-
-/***********************************************************************//**
- * @brief Read CTA events from FITS table (version 1)
- *
- * @param[in] table FITS table.
- *
- * This method reads the CTA event list from a FITS table HDU into memory.
- *
- * @todo Implement agreed column format
- ***************************************************************************/
-void GCTAEventList::read_events_v1(const GFitsTable& table)
-{
-    // Clear existing events
-    m_events.clear();
-
-    // Extract number of events in FITS file
-    int num = table.nrows();
-
-    // If there are events then load them
-    if (num > 0) {
-
-        // Reserve data
-        m_events.reserve(num);
-
-        // Get column pointers
-        const GFitsTableCol* ptr_eid         = table["EVENT_ID"];
-        const GFitsTableCol* ptr_oid         = table["OBS_ID"];
-        const GFitsTableCol* ptr_time        = table["TIME"];
-        const GFitsTableCol* ptr_multip      = table["MULTIP"];
-        const GFitsTableCol* ptr_ra          = table["RA"];
-        const GFitsTableCol* ptr_dec         = table["DEC"];
-        const GFitsTableCol* ptr_dir_err     = table["DIR_ERR"];
-        const GFitsTableCol* ptr_detx        = table["DETX"];
-        const GFitsTableCol* ptr_dety        = table["DETY"];
-        const GFitsTableCol* ptr_alt         = table["ALT"];
-        const GFitsTableCol* ptr_az          = table["AZ"];
-        const GFitsTableCol* ptr_corex       = table["COREX"];
-        const GFitsTableCol* ptr_corey       = table["COREY"];
-        const GFitsTableCol* ptr_core_err    = table["CORE_ERR"];
-        const GFitsTableCol* ptr_xmax        = table["XMAX"];
-        const GFitsTableCol* ptr_xmax_err    = table["XMAX_ERR"];
-        const GFitsTableCol* ptr_shw         = table["SHWIDTH"];
-        const GFitsTableCol* ptr_shl         = table["SHLENGTH"];
-        const GFitsTableCol* ptr_energy      = table["ENERGY"];
-        const GFitsTableCol* ptr_energy_err  = table["ENERGY_ERR"];
-
-        // Check for phase column
+        // Check for PHASE column
         const GFitsTableCol* ptr_phase;
         if (table.contains("PHASE")) {
             m_has_phase = true;
@@ -805,98 +699,57 @@ void GCTAEventList::read_events_v1(const GFitsTable&amp; table)
         // Copy data from columns into GCTAEventAtom objects
         GCTAEventAtom event;
         for (int i = 0; i < num; ++i) {
-            event.m_index     = i;
+
+            // Set mandatory information
             event.m_time.set(ptr_time->real(i), m_gti.reference());
             event.m_dir.dir().radec_deg(ptr_ra->real(i), ptr_dec->real(i));
-            event.m_dir.detx(ptr_detx->real(i)*gammalib::deg2rad);
-            event.m_dir.dety(ptr_dety->real(i)*gammalib::deg2rad);
             event.m_energy.TeV(ptr_energy->real(i));
-            event.m_event_id    = ptr_eid->integer(i);
-            event.m_obs_id      = ptr_oid->integer(i);
-            event.m_multip      = ptr_multip->integer(i);
-            event.m_telmask     = 0;
-            event.m_dir_err     = ptr_dir_err->real(i);
-            event.m_alt         = ptr_alt->real(i);
-            event.m_az          = ptr_az->real(i);
-            event.m_corex       = ptr_corex->real(i);
-            event.m_corey       = ptr_corey->real(i);
-            event.m_core_err    = ptr_core_err->real(i);
-            event.m_xmax        = ptr_xmax->real(i);
-            event.m_xmax_err    = ptr_xmax_err->real(i);
-            event.m_shwidth     = ptr_shw->real(i);
-            event.m_shlength    = ptr_shl->real(i);
-            event.m_energy_err  = ptr_energy_err->real(i);
-
-            // Set phase if available
+            event.m_event_id = ptr_eid->integer(i);
+            event.m_index    = i;
+
+            // If available, set detector coordinates in radians
+            if (m_has_detxy) {
+                event.m_dir.detx(ptr_detx->real(i)*gammalib::deg2rad);
+                event.m_dir.dety(ptr_dety->real(i)*gammalib::deg2rad);
+            }
+
+            // If available, set pulse phase
             if (m_has_phase) {
                 event.m_phase = ptr_phase->real(i);
             }
+            // Append event
             m_events.push_back(event);
-        }
-    } // endif: there were events
+        } // endfor: looped over all events
-    // Return
-    return;
-}
+        // Loop over table and find optional columns
+        for (int i = 0; i < table.ncols(); ++i) {
+            // Get column pointer
+            const GFitsTableCol* col = table[i];
-/***********************************************************************//**
- * @brief Read Hillas information for CTA events from FITS table
- *
- * @param[in] table FITS table.
- *
- * This method reads the Hillas reconstruction information for CTA events
- * from an EVENTS file. It searches for the columns HIL_MSW, HIL_MSW_ERR,
- * HIL_MSL, and HIL_MSL_ERR in the FITS table and extracts the relevant
- * columns from the FITS file. If a column is not found, no action is
- * performed.
- *
- * @todo Verify consistency of event list size
- ***************************************************************************/
-void GCTAEventList::read_events_hillas(const GFitsTable& table)
-{
-    // Extract number of events in FITS file
-    int num = table.nrows();
+            // Get column name
+            const std::string name = col->name();
-    //TODO: Make sure that dimension is consistent with existing event
-    //      list
-
-    // Continue only if there are events
-    if (num > 0) {
-
-        // HIL_MSW
-        if (table.contains("HIL_MSW")) {
-            const GFitsTableCol* ptr = table["HIL_MSW"];
-            for (int i = 0; i < num; ++i) {
-                m_events[i].m_hil_msw = ptr->real(i);
-            }
-        }
-
-        // HIL_MSW_ERR
-        if (table.contains("HIL_MSW_ERR")) {
-            const GFitsTableCol* ptr = table["HIL_MSW_ERR"];
-            for (int i = 0; i < num; ++i) {
-                m_events[i].m_hil_msw_err = ptr->real(i);
+            // If column was mandatory or optional then skip it ...
+            if (name == "EVENT_ID" ||
+                name == "TIME"     ||
+                name == "RA"       ||
+                name == "DEC"      ||
+                name == "ENERGY"   ||
+                name == "DETX"     ||
+                name == "DETY"     ||
+                name == "PHASE") {
+                continue;
             }
-        }
-        // HIL_MSL
-        if (table.contains("HIL_MSL")) {
-            const GFitsTableCol* ptr = table["HIL_MSL"];
-            for (int i = 0; i < num; ++i) {
-                m_events[i].m_hil_msl = ptr->real(i);
+            // ... otherwise keep a copy
+            else {
+                m_columns.push_back(col->clone());
             }
-        }
-        // HIL_MSL_ERR
-        if (table.contains("HIL_MSL_ERR")) {
-            const GFitsTableCol* ptr = table["HIL_MSL_ERR"];
-            for (int i = 0; i < num; ++i) {
-                m_events[i].m_hil_msl_err = ptr->real(i);
-            }
-        }
+        } // endfor: loop over optional columns
     } // endif: there were events
@@ -908,12 +761,26 @@ void GCTAEventList::read_events_hillas(const GFitsTable&amp; table)
 /***********************************************************************//**
  * @brief Write CTA events into FITS table
  *
- * @param[in] hdu FITS table HDU.
+ * @param[in] table FITS binary table.
  *
- * Write the CTA event list into FITS table.
+ * Write the CTA event list into a FITS binary table.
  *
- * @todo The TELMASK column is allocated with a dummy length of 100.
- * @todo Implement agreed column format
+ * The following mandatory columns will be written into the FITS table:
+ *
+ *      - EVENT_ID
+ *      - TIME
+ *      - RA
+ *      - DEC
+ *      - ENERGY
+ *
+ * If available, also the
+ *
+ *      - DETX
+ *      - DETY
+ *      - PHASE
+ *
+ * columns are written. Any other columns that were read by the read_events()
+ * method will be also written into the table.
  ***************************************************************************/
 void GCTAEventList::write_events(GFitsBinTable& hdu) const
 {
@@ -923,109 +790,55 @@ void GCTAEventList::write_events(GFitsBinTable&amp; hdu) const
     // If there are events then write them now
     if (size() > 0) {
-        // Allocate columns
-        GFitsTableULongCol  col_eid         = GFitsTableULongCol("EVENT_ID", size());
-        GFitsTableULongCol  col_oid         = GFitsTableULongCol("OBS_ID", size());
-        GFitsTableDoubleCol col_time        = GFitsTableDoubleCol("TIME", size());
-        GFitsTableDoubleCol col_live        = GFitsTableDoubleCol("TLIVE", size());
-        GFitsTableShortCol  col_multip      = GFitsTableShortCol("MULTIP", size());
-        GFitsTableBitCol    col_telmask     = GFitsTableBitCol("TELMASK", size(), 100);
-        GFitsTableFloatCol  col_ra          = GFitsTableFloatCol("RA", size());
-        GFitsTableFloatCol  col_dec         = GFitsTableFloatCol("DEC", size());
-        GFitsTableFloatCol  col_direrr      = GFitsTableFloatCol("DIR_ERR", size());
-        GFitsTableFloatCol  col_detx        = GFitsTableFloatCol("DETX", size());
-        GFitsTableFloatCol  col_dety        = GFitsTableFloatCol("DETY", size());
-        GFitsTableFloatCol  col_alt         = GFitsTableFloatCol("ALT", size());
-        GFitsTableFloatCol  col_az          = GFitsTableFloatCol("AZ", size());
-        GFitsTableFloatCol  col_corex       = GFitsTableFloatCol("COREX", size());
-        GFitsTableFloatCol  col_corey       = GFitsTableFloatCol("COREY", size());
-        GFitsTableFloatCol  col_core_err    = GFitsTableFloatCol("CORE_ERR", size());
-        GFitsTableFloatCol  col_xmax        = GFitsTableFloatCol("XMAX", size());
-        GFitsTableFloatCol  col_xmax_err    = GFitsTableFloatCol("XMAX_ERR", size());
-        GFitsTableFloatCol  col_shw         = GFitsTableFloatCol("SHWIDTH", size());
-        GFitsTableFloatCol  col_shl         = GFitsTableFloatCol("SHLENGTH", size());
-        GFitsTableFloatCol  col_energy      = GFitsTableFloatCol("ENERGY", size());
-        GFitsTableFloatCol  col_energy_err  = GFitsTableFloatCol("ENERGY_ERR", size());
-        GFitsTableFloatCol  col_hil_msw     = GFitsTableFloatCol("HIL_MSW", size());
-        GFitsTableFloatCol  col_hil_msw_err = GFitsTableFloatCol("HIL_MSW_ERR", size());
-        GFitsTableFloatCol  col_hil_msl     = GFitsTableFloatCol("HIL_MSL", size());
-        GFitsTableFloatCol  col_hil_msl_err = GFitsTableFloatCol("HIL_MSL_ERR", size());
-
-        // Create phase column if required
-        GFitsTableFloatCol  col_phase;
-        if (m_has_phase) {
-            col_phase = GFitsTableFloatCol("PHASE", size());
-        }
+        // Allocate mandatory columns
+        GFitsTableULongCol  col_eid    = GFitsTableULongCol("EVENT_ID", size());
+        GFitsTableDoubleCol col_time   = GFitsTableDoubleCol("TIME", size());
+        GFitsTableFloatCol  col_ra     = GFitsTableFloatCol("RA", size());
+        GFitsTableFloatCol  col_dec    = GFitsTableFloatCol("DEC", size());
+        GFitsTableFloatCol  col_energy = GFitsTableFloatCol("ENERGY", size());
-        // Fill columns
+        // Fill mandatory columns
         for (int i = 0; i < size(); ++i) {
-            col_eid(i)         = m_events[i].m_event_id;
-            col_oid(i)         = m_events[i].m_obs_id;
-            col_time(i)        = m_events[i].time().convert(m_gti.reference());
-            col_live(i)        = 0.0;
-            col_multip(i)      = 0;
-            //col_telmask
-            col_ra(i)          = m_events[i].dir().dir().ra_deg();
-            col_dec(i)         = m_events[i].dir().dir().dec_deg();
-            col_direrr(i)      = m_events[i].m_dir_err;
-            col_detx(i)        = m_events[i].dir().detx() * gammalib::rad2deg;
-            col_dety(i)        = m_events[i].dir().dety() * gammalib::rad2deg;
-            col_alt(i)         = m_events[i].m_alt;
-            col_az(i)          = m_events[i].m_az;
-            col_corex(i)       = m_events[i].m_corex;
-            col_corey(i)       = m_events[i].m_corey;
-            col_core_err(i)    = m_events[i].m_core_err;
-            col_xmax(i)        = m_events[i].m_xmax;
-            col_xmax_err(i)    = m_events[i].m_xmax_err;
-            col_shw(i)         = m_events[i].m_shwidth;
-            col_shl(i)         = m_events[i].m_shlength;
-            col_energy(i)      = m_events[i].energy().TeV();
-            col_energy_err(i)  = m_events[i].m_energy_err;
-            col_hil_msw(i)     = m_events[i].m_hil_msw;
-            col_hil_msw_err(i) = m_events[i].m_hil_msw_err;
-            col_hil_msl(i)     = m_events[i].m_hil_msl;
-            col_hil_msl_err(i) = m_events[i].m_hil_msl_err;
-
-            // Optionally fill pulse phase column
-            if (m_has_phase) {
-                col_phase(i) = m_events[i].m_phase;
-            }
-
-        } // endfor: looped over rows
+            col_eid(i)    = m_events[i].m_event_id;
+            col_time(i)   = m_events[i].time().convert(m_gti.reference());
+            col_ra(i)     = m_events[i].dir().dir().ra_deg();
+            col_dec(i)    = m_events[i].dir().dir().dec_deg();
+            col_energy(i) = m_events[i].energy().TeV();
+        }
-        // Append columns to table
+        // Append mandatory columns to table
         hdu.append(col_eid);
-        hdu.append(col_oid);
         hdu.append(col_time);
-        hdu.append(col_live);
-        hdu.append(col_multip);
-        hdu.append(col_telmask);
         hdu.append(col_ra);
         hdu.append(col_dec);
-        hdu.append(col_direrr);
-        hdu.append(col_detx);
-        hdu.append(col_dety);
-        hdu.append(col_alt);
-        hdu.append(col_az);
-        hdu.append(col_corex);
-        hdu.append(col_corey);
-        hdu.append(col_core_err);
-        hdu.append(col_xmax);
-        hdu.append(col_xmax_err);
-        hdu.append(col_shw);
-        hdu.append(col_shl);
         hdu.append(col_energy);
-        hdu.append(col_energy_err);
-        hdu.append(col_hil_msw);
-        hdu.append(col_hil_msw_err);
-        hdu.append(col_hil_msl);
-        hdu.append(col_hil_msl_err);
-        // Optionally write phase column
+        // If available, add detector coordinates in degrees
+        if (m_has_detxy) {
+            GFitsTableFloatCol col_detx = GFitsTableFloatCol("DETX", size());
+            GFitsTableFloatCol col_dety = GFitsTableFloatCol("DETY", size());
+            for (int i = 0; i < size(); ++i) {
+                col_detx(i) = m_events[i].dir().detx() * gammalib::rad2deg;
+                col_dety(i) = m_events[i].dir().dety() * gammalib::rad2deg;
+            }
+            hdu.append(col_detx);
+            hdu.append(col_dety);
+        }
+
+        // If available, add event phase
         if (m_has_phase) {
+            GFitsTableFloatCol col_phase = GFitsTableFloatCol("PHASE", size());
+            for (int i = 0; i < size(); ++i) {
+                col_phase(i) = m_events[i].m_phase;
+            }
             hdu.append(col_phase);
         }
+        // Append other columns to table
+        for (int i = 0; i < m_columns.size(); ++i) {
+            hdu.append(*m_columns[i]);
+        }
+
     } // endif: there were events to write
     // Return
@@ -72,6 +72,7 @@ public:
     void                ts(const double& ts);
     const bool&         tscalc(void) const;
     void                tscalc(const bool& tscalc);
+    const bool&         has_ts(void) const;
     GModelPar           scale(const std::string& instrument) const;
     void                scale(const GModelPar& par);
     std::string         ids(void) const;
@@ -517,7 +517,7 @@ bool GVOClient::find_hub(void)
     if (!lockurl.empty()) {
         // If we have a file:// prefix then strip it now. This is a kluge
-        // and should be remplaced by a method that allows opening any kind
+        // and should be replaced by a method that allows opening any kind
         // of URL
         if (lockurl.compare(0, 7, "file://") == 0) {
             lockurl = lockurl.substr(7, std::string::npos);
@@ -955,7 +955,6 @@ std::string GVOClient::receive_string(void) const
         int n       = 0;
         do {
             n = gammalib::recv(m_socket, buffer, 1000, 0, timeout);
-            //n = recv(m_socket, buffer, 1000, 0);
             if (n > 0) {
                 buffer[n] = '\0';
                 result.append(std::string(buffer));
@@ -1074,7 +1073,7 @@ std::string GVOClient::get_hub_lockfile(void) const
     char* hub_ptr = std::getenv("SAMP_HUB");
     if (hub_ptr != NULL) {
-        // Check for mandatory std-lockurl: prefix (no other prefixe is
+        // Check for mandatory std-lockurl: prefix (no other prefix is
         // supported so far)
         std::string lockurl = std::string(hub_ptr);
         if (lockurl.compare(0, 12, "std-lockurl:") == 0) {
@@ -1084,7 +1083,7 @@ std::string GVOClient::get_hub_lockfile(void) const
         } // endif: std-lockurl: prefix found
-    }
+    } // endif: SAMP_HUB environment variable found
     // ... otherwise the lockfile should be $HOME/.samp
     else {
@@ -47,6 +47,7 @@
 /* __ Method name definitions ____________________________________________ */
 #define G_START_HUB                                     "GVOHub::start_hub()"
+#define G_START_SOCKET                               "GVOHub::start_socket()"
 /* __ Macros _____________________________________________________________ */
@@ -179,9 +180,6 @@ GVOHub* GVOHub::clone(void) const
  ***************************************************************************/
 void GVOHub::start(void)
 {
-    // Create SAMP file
-    create_samp_file();
-
     // Start Hub
     start_hub();
@@ -210,7 +208,7 @@ std::string GVOHub::print(const GChatter&amp; chatter) const
         // Append Hub information
         result.append("\n"+gammalib::parformat("Hub identifier")+m_hub_id);
         result.append("\n"+gammalib::parformat("Hub key")+m_secret);
-        result.append("\n"+gammalib::parformat("Hub URL")+m_hub_url);
+        result.append("\n"+gammalib::parformat("Hub URL")+hub_url());
         result.append("\n"+gammalib::parformat("Hub host")+m_hub_host);
         result.append("\n"+gammalib::parformat("Hub port")+m_hub_port);
         result.append("\n"+gammalib::parformat("Hub path")+m_hub_path);
@@ -240,7 +238,6 @@ void GVOHub::init_members(void)
     m_hub_host      = "127.0.0.1";
     m_hub_port      = "2526";
     m_hub_path      = "xmlrpc";
-    m_hub_url       = "http://"+m_hub_host+":"+m_hub_port+"/"+m_hub_path;
     m_version       = "1.3";
     m_hub_id        = "gammalib_hub";
     m_socket        = -1;        // Signals no socket
@@ -261,7 +258,6 @@ void GVOHub::copy_members(const GVOHub&amp; hub)
 {
     // Copy members
     m_secret   = hub.m_secret;
-    m_hub_url  = hub.m_hub_url;
     m_hub_host = hub.m_hub_host;
     m_hub_port = hub.m_hub_port;
     m_hub_path = hub.m_hub_path;
@@ -288,8 +284,7 @@ void GVOHub::free_members(void)
     }
     // Remove lockfile
-    std::string lockurl = get_hub_lockfile();
-    std::remove(lockurl.c_str());
+    delete_samp_file();
     // Return
     return;
@@ -306,46 +301,13 @@ void GVOHub::free_members(void)
  ***************************************************************************/
 void GVOHub::start_hub(void)
 {
-    // Prepare TCP/IP structure
-    struct sockaddr_in serv_addr;
-    std::memset(&serv_addr, 0, sizeof(serv_addr));
-    serv_addr.sin_family      = AF_INET;
-    serv_addr.sin_addr.s_addr = inet_addr(m_hub_host.c_str());
-    serv_addr.sin_port        = htons(gammalib::toint(m_hub_port));
-    
-    // Create Hub socket
-    m_socket = socket(AF_INET, SOCK_STREAM, 0);
-    
-    // Creation of hub main socket
-    if (m_socket < 0) {
-        std::string msg = "Unable to create Hub socket. Errno="+
-                          gammalib::str(errno);
-        throw GException::runtime_error(G_START_HUB, msg);
-    }
-    
-    // Set hub main socket to allow multiple connections
-    int opt = 1;
-    if (setsockopt(m_socket, SOL_SOCKET, SO_REUSEADDR, (char*)&opt, sizeof(opt)) < 0) {
-        std::string msg = "Unable to set Hub socket to multiple connections."
-                          " Errno="+gammalib::str(errno);
-        throw GException::runtime_error(G_START_HUB, msg);
-    }
-    #ifdef SO_REUSEPORT
-    if (setsockopt(m_socket, SOL_SOCKET, SO_REUSEPORT, (char*)&opt, sizeof(opt)) < 0) {
-        std::string msg = "Unable to set Hub socket to multiple connections."
-                          " Errno="+gammalib::str(errno);
-        throw GException::runtime_error(G_START_HUB, msg);
-    }
-    #endif
-    
-    // Server socket is opened. Now, bind it to the port, with family etc.
-    if (bind(m_socket, (struct sockaddr*)&serv_addr, sizeof(serv_addr)) < 0) {
-        std::string msg = "Unable to bind Hub socket to server socket. Errno="+
-                          gammalib::str(errno);
-        throw GException::runtime_error(G_START_HUB, msg);
-    }
+    // Get socket
+    m_socket = get_socket();
+
+    // Create SAMP file
+    create_samp_file();
-    // Now start listening for the clients: 5 requests simultaneously pending
+    // Start listening for the clients: 5 requests simultaneously pending
     // at maximum
     if (listen(m_socket, 5) < 0) {
         std::string msg = "Unable to start listening on Hub socket. Errno="+
@@ -385,6 +347,9 @@ void GVOHub::start_hub(void)
     close(m_socket);
     m_socket = -1;
+    // Delete SAMP file
+    delete_samp_file();
+
     // Return
     return;
 }
@@ -596,7 +561,7 @@ void GVOHub::request_register(const GXml&amp; xml, const socklen_t&amp; sock)
     response.append("    </member>\n");
     response.append("    <member>\n");
     response.append("      <name>samp.url-translator</name>\n");
-    response.append("      <value>"+m_hub_url+"</value>\n");
+    response.append("      <value>"+hub_url()+"</value>\n");
     response.append("    </member>\n");
     response.append("  </struct></value></param>\n");
     response.append("</params>\n");
@@ -1628,9 +1593,10 @@ void GVOHub::send_notifications(const std::string&amp; method,
 /***********************************************************************//**
- * @brief Creates the lockfile, fill it
+ * @brief Create the lockfile
  *
- * Implements IVOA standard REC-SAMP-1.3-20120411.
+ * Creates the SAMP lockfile and fill it with information. Implements IVOA
+ * standard REC-SAMP-1.3-20120411.
  ***************************************************************************/
 void GVOHub::create_samp_file(void) const
 {
@@ -1647,7 +1613,7 @@ void GVOHub::create_samp_file(void) const
         fprintf(fptr, "# SAMP lockfile\n");
         fprintf(fptr, "# Required keys:\n");
         fprintf(fptr, "samp.secret=%s\n", m_secret.c_str());
-        fprintf(fptr, "samp.hub.xmlrpc.url=%s\n", m_hub_url.c_str());
+        fprintf(fptr, "samp.hub.xmlrpc.url=%s\n", hub_url().c_str());
         fprintf(fptr, "samp.profile.version=%s\n", m_version.c_str());
         fprintf(fptr, "# Info stored by hub for some private reason:\n");
         fprintf(fptr, "gammalib.hubid=%s\n", m_hub_id.c_str());
@@ -1663,6 +1629,87 @@ void GVOHub::create_samp_file(void) const
 /***********************************************************************//**
+ * @brief Delete the lockfile
+ *
+ * Deletes the SAMP lockfile on disk.
+ ***************************************************************************/
+void GVOHub::delete_samp_file(void) const
+{
+    // Get lockfile URL
+    std::string lockurl = get_hub_lockfile();
+
+    // Delete lockfile
+    std::remove(lockurl.c_str());
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Get Hub socket
+ *
+ * @return Hub socket (-1 if no socket was found).
+ *
+ * Returns Hub socket. Starting from an initial port of 2526 the method
+ * searches for the next free port and binds it to a socket.
+ ***************************************************************************/
+int GVOHub::get_socket(void)
+{
+    // Initialise socket, port and socket address structure
+    int    sock = -1;
+    int    port = 2526;
+    struct sockaddr_in serv_addr;
+
+    // Loop over 100 ports at most
+    for (int i = 0; i < 100; ++i) {
+    
+        // Clean TCP/IP structure
+        std::memset(&serv_addr, 0, sizeof(serv_addr));
+
+        // Set TCP/IP structure
+        serv_addr.sin_family      = AF_INET;
+        serv_addr.sin_addr.s_addr = inet_addr(m_hub_host.c_str());
+        serv_addr.sin_port        = htons(port);
+    
+        // Create Hub socket
+        sock = socket(AF_INET, SOCK_STREAM, 0);
+    
+        // Creation of hub main socket
+        if (sock < 0) {
+            std::string msg = "Unable to create Hub socket. Errno="+
+                              gammalib::str(errno);
+            throw GException::runtime_error(G_START_SOCKET, msg);
+        }
+
+        // Now bind socket to the address and port
+        if (bind(sock, (struct sockaddr*)&serv_addr, sizeof(serv_addr)) < 0) {
+            if (errno == EADDRINUSE) {
+                port++;
+                close(sock);
+            }
+            else {
+                std::string msg = "Unable to bind Hub socket to address "+
+                                  m_hub_host+":"+gammalib::str(port)+
+                                  ". Errno="+gammalib::str(errno);
+                throw GException::runtime_error(G_START_HUB, msg);
+            }
+        }
+        else {
+            break;
+        }
+
+    } // endfor: looped over ports
+
+    // Store port as string
+    m_hub_port = gammalib::str(port);
+
+    // Return socket
+    return sock;
+}
+
+
+/***********************************************************************//**
  * @brief Post string content to client
  *
  * @param[in] content String content to post.
@@ -1948,3 +1995,20 @@ std::string GVOHub::random_string(const size_t&amp; length) const
     return str;
 }
+
+/***********************************************************************//**
+ * @brief Return Hub URL
+ *
+ * @return Hub URL.
+ *
+ * Returns the XML-RPC endpoint for communication with the hub.
+ ***************************************************************************/
+std::string GVOHub::hub_url(void) const
+{
+    // Set Hub URL
+    std::string hub_url = "http://"+m_hub_host+":"+m_hub_port+"/"+m_hub_path;
+
+    // Return Hub URL
+    return (hub_url);
+}
+
@@ -175,6 +175,10 @@ int main(void)
     // Allocate test suit container
     GTestSuites testsuites("VO module");
+    // Unset HOME environment variable (so that the .samp file is written
+    // into the test directory)
+    unsetenv("HOME");
+
     // Initially assume that we pass all tests
     bool success = true;