GVOTable skeleton and use through GSkymap on the way

Thierry Louge
1 parent 6afd29a0
Showing 12 changed files with 6503 additions and 7 deletions
include/GSkymap.hpp
include/GVOHub.hpp
include/GVOTable.hpp
src/sky/GSkymap.cpp
src/sky/GSkymap.cpp~
src/vo/GVOApp.cpp~
src/vo/GVOClient.cpp~
src/vo/GVOHub.cpp
src/vo/GVOHub.cpp~
src/vo/GVOTable.cpp
src/vo/GVOTable.cpp~
src/vo/Makefile.am
@@ -41,7 +41,8 @@
 #include "GMatrix.hpp"
 #include "GVector.hpp"
 #include "GBilinear.hpp"
-
+#include "GVOTable.hpp"
+#include "GVOClient.hpp"
  
 /***********************************************************************//**
  * @class GSkymap
@@ -164,6 +165,8 @@ private:
     void              alloc_wcs(const GFitsImage& image);
     GFitsBinTable*    create_healpix_hdu(void) const;
     GFitsImageDouble* create_wcs_hdu(void) const;
+    void 	      publish(void);
+    void 	      make_votable(void);
  
     // Private data area
     int               m_num_pixels; //!< Number of pixels (used for pixel allocation)
@@ -172,6 +175,8 @@ private:
     int               m_num_y;      //!< Number of pixels in y direction (only 2D)
     GSkyProjection*   m_proj;       //!< Pointer to sky projection
     double*           m_pixels;     //!< Pointer to skymap pixels
+    GVOTable*	      m_skytable;   //!< Pointer to skymap VOTable
+    GVOClient*	      m_skyclient;  //!< Pointer to VOClient opened
  
     // Computation cache
     mutable bool      m_hascache;   //!< Cache is valid
@@ -105,6 +105,10 @@ protected:
     std::string              get_hub_lockfile(void) const;
     void 		             send_notifications(const std::string& method,
                                                 const std::string& client);
+    void                     request_call(const GXml& xml,
+                                                  const socklen_t& sock);
+    std::string              get_mtype(const GXml& xml) const;
+    std::string              get_destination(const GXml& xml) const;
  
     // Low-level methods
     void        create_samp_file(void) const;
+/***************************************************************************
+ *                      GVOTable.hpp - VOTable class                     *
+ * ----------------------------------------------------------------------- *
+ *  copyright (C) 2014-2015 by Thierry Louge                               *
+ * ----------------------------------------------------------------------- *
+ *                                                                         *
+ *  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   *
+ *  the Free Software Foundation, either version 3 of the License, or      *
+ *  (at your option) any later version.                                    *
+ *                                                                         *
+ *  This program is distributed in the hope that it will be useful,        *
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of         *
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          *
+ *  GNU General Public License for more details.                           *
+ *                                                                         *
+ *  You should have received a copy of the GNU General Public License      *
+ *  along with this program.  If not, see <http://www.gnu.org/licenses/>.  *
+ *                                                                         *
+ ***************************************************************************/
+/**
+ * @file GVOTable.hpp
+ * @brief VOTable class definition
+ * @author Thierry Louge
+ */
+
+#ifndef GVOTable_HPP
+#define GVOTable_HPP
+
+/* __ Definitions ________________________________________________________ */
+
+/* __ Includes ___________________________________________________________ */
+#include <string>
+#include <vector>
+#include <sys/socket.h>
+#include "GBase.hpp"
+#include "GXml.hpp"
+#include "GXmlNode.hpp"
+
+
+/***********************************************************************//**
+ * @class GVOTable
+ *
+ * @brief VOTable class
+ *
+ * This class implements a VOTable for exchanges through VO-compatible
+ * applications.
+ ***************************************************************************/
+class GVOTable : public GBase {
+
+public:
+    // Constructors and destructors
+    GVOTable(void);
+    GVOTable(const GVOTable& hub);
+    virtual ~GVOTable(void);
+
+    // Operators
+    GVOTable& operator=(const GVOTable& hub);
+
+    // Methods
+    void        	clear(void);
+    GVOTable*     	clone(void) const;
+    std::string 	classname(void) const;
+    std::string		print(const GChatter& chatter = NORMAL) const;
+    void 		create_votable(void);
+    void 		fill_content(const std::string& data) const;
+
+protected:
+    // Protected methods
+    void                init_members(void);
+    void                copy_members(const GVOTable& client);
+    void                free_members(void);
+   
+    // Low-level methods
+private:
+    // Private data area
+    GXml		m_tablexml;   //!< GXml containing the VOTable
+};
+
+
+/***********************************************************************//**
+ * @brief Return class name
+ *
+ * @return String containing the class name ("GVOTable").
+ ***************************************************************************/
+inline
+std::string GVOTable::classname(void) const
+{
+    return ("GVOTable");
+}
+
+#endif /* GVOTable_HPP */
@@ -2261,7 +2261,38 @@ GFitsImageDouble* GSkymap::create_wcs_hdu(void) const
     return hdu;
 }
  
+/***********************************************************************//**
+ * @brief publishes the skymap in a VO Hub
+ *
+ * @param[in] 
+ *
+ * @exception GException::
+ *            
+ ***************************************************************************/
+void GSkymap::publish(void)
+{
+    make_votable();
+    // Return
+    return;
+}
  
+/***********************************************************************//**
+ * @brief makes GVOTable object and fills ot with current skymap(s)
+ *
+ * @param[in] 
+ *
+ * @exception GException::
+ *            
+ ***************************************************************************/
+void GSkymap::make_votable(void)
+{
+    GVOTable* m_skytable;
+    m_skytable->create_votable();
+    //For each pixel set the value of the pixel inside the data part
+    //m_skytable->fill_content();
+    // Return
+    return;
+}
 /*==========================================================================
  =                                                                         =
  =                                 Friends                                 =
+/***************************************************************************
+ *                       GSkymap.cpp - Sky map class                       *
+ * ----------------------------------------------------------------------- *
+ *  copyright (C) 2010-2015 by 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   *
+ *  the Free Software Foundation, either version 3 of the License, or      *
+ *  (at your option) any later version.                                    *
+ *                                                                         *
+ *  This program is distributed in the hope that it will be useful,        *
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of         *
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          *
+ *  GNU General Public License for more details.                           *
+ *                                                                         *
+ *  You should have received a copy of the GNU General Public License      *
+ *  along with this program.  If not, see <http://www.gnu.org/licenses/>.  *
+ *                                                                         *
+ ***************************************************************************/
+/**
+ * @file GSkymap.cpp
+ * @brief Sky map class implementation
+ * @author Juergen Knoedlseder
+ */
+
+/* __ Includes ___________________________________________________________ */
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+#include "GException.hpp"
+#include "GTools.hpp"
+#include "GSkymap.hpp"
+#include "GHealpix.hpp"
+#include "GWcsRegistry.hpp"
+#include "GWcs.hpp"
+#include "GFits.hpp"
+#include "GFitsTableDoubleCol.hpp"
+#include "GFitsImageDouble.hpp"
+
+/* __ Method name definitions ____________________________________________ */
+#define G_CONSTRUCT_HPX                "GSkymap::GSkymap(std::string&, int&,"\
+                                                       " std::string&, int&)"
+#define G_CONSTRUCT_MAP        "GSkymap::GSkymap(std::string&, std::string&,"\
+                      " double&, double&, double& double&, int&, int&, int&)"
+#define G_NMAPS                                        "GSkymap::nmaps(int&)"
+#define G_OP_UNARY_ADD                        "GSkymap::operator+=(GSkymap&)"
+#define G_OP_UNARY_SUB                        "GSkymap::operator-=(GSkymap&)"
+#define G_OP_UNARY_MUL                        "GSkymap::operator-=(GSkymap&)"
+#define G_OP_UNARY_DIV                        "GSkymap::operator-=(GSkymap&)"
+#define G_OP_ACCESS_1D                        "GSkymap::operator(int&, int&)"
+#define G_OP_ACCESS_2D                  "GSkymap::operator(GSkyPixel&, int&)"
+#define G_OP_VALUE                        "GSkymap::operator(GSkyDir&, int&)"
+#define G_INX2DIR                                    "GSkymap::inx2dir(int&)"
+#define G_PIX2DIR                              "GSkymap::pix2dir(GSkyPixel&)"
+#define G_DIR2INX                                "GSkymap::dir2inx(GSkyDir&)"
+#define G_DIR2PIX                                "GSkymap::dir2pix(GSkyDir&)"
+#define G_SOLIDANGLE1                             "GSkymap::solidangle(int&)"
+#define G_SOLIDANGLE2                       "GSkymap::solidangle(GSkyPixel&)"
+#define G_EXTRACT                              "GSkymap::extract(int&, int&)"
+#define G_READ                               "GSkymap::read(const GFitsHDU&)"
+#define G_SET_WCS     "GSkymap::set_wcs(std::string&, std::string&, double&,"\
+                              " double&, double&, double&, double&, double&,"\
+                                                       " GMatrix&, GVector&)"
+#define G_READ_HEALPIX                   "GSkymap::read_healpix(GFitsTable*)"
+#define G_READ_WCS                           "GSkymap::read_wcs(GFitsImage*)"
+#define G_ALLOC_WCS                         "GSkymap::alloc_wcs(GFitsImage*)"
+
+/* __ Macros _____________________________________________________________ */
+
+/* __ Coding definitions _________________________________________________ */
+
+/* __ Debug definitions __________________________________________________ */
+//#define G_READ_HEALPIX_DEBUG                          // Debug read_healpix
+//#define G_READ_WCS_DEBUG                                  // Debug read_wcs
+
+/* __ Prototype __________________________________________________________ */
+
+
+/*==========================================================================
+ =                                                                         =
+ =                          Constructors/destructors                       =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Void constructor
+ ***************************************************************************/
+GSkymap::GSkymap(void)
+{
+    // Initialise class members for clean destruction
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief FITS file constructor
+ *
+ * @param[in] filename FITS file name.
+ ***************************************************************************/
+GSkymap::GSkymap(const std::string& filename)
+{
+    // Initialise class members for clean destruction
+    init_members();
+
+    // Load skymap
+    load(filename);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Healpix sky map constructor
+ *
+ * @param[in] coords Coordinate System (CEL or GAL).
+ * @param[in] nside Nside parameter.
+ * @param[in] order Pixel ordering (RING or NEST).
+ * @param[in] nmaps Number of maps in set (default=1).
+ *
+ * @exception GException::skymap_bad_par
+ *            Invalid sky map parameter.
+ *
+ * Constructs sky map in Healpix pixelisation.
+ ***************************************************************************/
+GSkymap::GSkymap(const std::string& coords,
+                 const int&         nside,
+                 const std::string& order,
+                 const int&         nmaps)
+{
+    // Initialise class members for clean destruction
+    init_members();
+
+    // Check if nmaps parameter is >0
+    if (nmaps < 1) {
+        throw GException::skymap_bad_par(G_CONSTRUCT_HPX, nmaps,
+                                         "nmaps parameter must be >0.");
+    }
+
+    // Allocate Healpix projection
+    GHealpix* projection = new GHealpix(nside, order, coords);
+    m_proj               = projection;
+
+    // Set number of pixels and number of maps
+    m_num_pixels = projection->npix();
+    m_num_maps   = nmaps;
+
+    // Allocate pixels
+    alloc_pixels();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief WCS sky map constructor
+ *
+ * @param[in] wcs World Coordinate System.
+ * @param[in] coords Coordinate System (CEL or GAL).
+ * @param[in] x X coordinate of sky map centre (deg).
+ * @param[in] y Y coordinate of sky map centre (deg).
+ * @param[in] dx Pixel size in x direction at centre (deg/pixel).
+ * @param[in] dy Pixel size in y direction at centre (deg/pixel).
+ * @param[in] nx Number of pixels in x direction.
+ * @param[in] ny Number of pixels in y direction.
+ * @param[in] nmaps Number of maps in set (default=1).
+ *
+ * @exception GException::skymap_bad_par
+ *            Invalid sky map parameter.
+ *
+ * Constructs sky map in World Coordinate System projection.
+ ***************************************************************************/
+GSkymap::GSkymap(const std::string& wcs,
+                 const std::string& coords,
+                 const double&      x,
+                 const double&      y,
+                 const double&      dx,
+                 const double&      dy,
+                 const int&         nx,
+                 const int&         ny,
+                 const int&         nmaps)
+{
+    // Initialise class members for clean destruction
+    init_members();
+
+    // Check parameters
+    if (nx < 1) {
+        throw GException::skymap_bad_par(G_CONSTRUCT_MAP, nx,
+                                         "nx parameter must be >0.");
+    }
+    if (ny < 1) {
+        throw GException::skymap_bad_par(G_CONSTRUCT_MAP, ny,
+                                         "ny parameter must be >0.");
+    }
+    if (nmaps < 1) {
+        throw GException::skymap_bad_par(G_CONSTRUCT_MAP, nmaps,
+                                         "nmaps parameter must be >0.");
+    }
+
+    // Set WCS
+    double  crval1 = x;
+    double  crval2 = y;
+    double  crpix1 = double(nx+1)/2.0;
+    double  crpix2 = double(ny+1)/2.0;
+    double  cdelt1 = dx;
+    double  cdelt2 = dy;
+    GMatrix cd(2,2);
+    GVector pv2(21);
+    set_wcs(wcs, coords, crval1, crval2, crpix1, crpix2, cdelt1, cdelt2, 
+            cd, pv2);
+
+    // Set number of pixels and number of maps
+    m_num_x      = nx;
+    m_num_y      = ny;
+    m_num_pixels = m_num_x * m_num_y;
+    m_num_maps   = nmaps;
+
+    // Allocate pixels
+    alloc_pixels();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy constructor
+ *
+ * @param[in] map Sky map.
+ ***************************************************************************/
+GSkymap::GSkymap(const GSkymap& map)
+{
+    // Initialise class members for clean destruction
+    init_members();
+
+    // Copy members
+    copy_members(map);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Destructor
+ ***************************************************************************/
+GSkymap::~GSkymap(void)
+{
+    // Free members
+    free_members();
+
+    // Return
+    return;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                               Operators                                 =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Assignment operator
+ *
+ * @param[in] map Sky map.
+ * @return Sky map.
+ ***************************************************************************/
+GSkymap& GSkymap::operator=(const GSkymap& map)
+{
+    // Execute only if object is not identical
+    if (this != &map) {
+
+        // Free members
+        free_members();
+
+        // Initialise private members for clean destruction
+        init_members();
+
+        // Copy members
+        copy_members(map);
+
+    } // endif: object was not identical
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Value setting operator
+ *
+ * @param[in] value Value.
+ * @return Sky map.
+ *
+ * Sets all pixels to the specified @p value.
+ ***************************************************************************/
+GSkymap& GSkymap::operator=(const double& value)
+{
+    // Get number of pixels
+    int num = m_num_pixels * m_num_maps;
+
+    // Loop over all pixels
+    for (int i = 0; i < num; ++i) {
+        m_pixels[i] = value;
+    }
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Map addition operator
+ *
+ * @param[in] map Sky map.
+ * @return Sky map.
+ *
+ * @exception GException::invalid_value
+ *            Mismatch between number of maps in skymap object.
+ *
+ * Adds the content of @p map to the skymap. The operator only works on sky
+ * maps with an identical number of layers. The content is added by
+ * bi-linearily interpolating the values in the source sky map, allowing thus
+ * for a reprojection of sky map values.
+ *
+ * @todo The method is not optimized for speed as the transformation is done
+ * for each layer separately. A private method should be introduced that
+ * does the transformation, allowing the to loop more effectively over the
+ * layers.
+ ***************************************************************************/
+GSkymap& GSkymap::operator+=(const GSkymap& map)
+{
+    // Check if number of layers are identical
+    if (map.nmaps() != nmaps()) {
+        std::string msg = "Mismatch of number of maps in skymap object"
+                          " ("+gammalib::str(nmaps())+" maps in destination"
+                          " map, "+gammalib::str(map.nmaps())+" in source"
+                          " map.";
+        throw GException::invalid_value(G_OP_UNARY_ADD, msg);
+    }
+
+    // Loop over all pixels of sky map
+    for (int index = 0; index < npix(); ++index) {
+
+        // Get sky direction of actual pixel
+        GSkyDir dir = inx2dir(index);
+
+        // Loop over all layers
+        for (int layer = 0; layer < nmaps(); ++layer) {
+        
+            // Add value
+            (*this)(index, layer) += map(dir, layer);
+
+        } // endfor: looped over all layers
+        
+    } // endfor: looped over all pixels
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Value addition operator
+ *
+ * @param[in] value Value.
+ * @return Sky map.
+ *
+ * Add @p value to all sky map pixels.
+ ***************************************************************************/
+GSkymap& GSkymap::operator+=(const double& value)
+{
+    // Set total number of sky map pixels
+    int num = m_num_pixels * m_num_maps;
+
+    // Loop over all pixels of sky map
+    for (int i = 0; i < num; ++i) {
+        m_pixels[i] += value;
+    }
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Map subtraction operator
+ *
+ * @param[in] map Sky map.
+ * @return Sky map.
+ *
+ * @exception GException::invalid_value
+ *            Mismatch between number of maps in skymap object.
+ *
+ * Subtracts the content of @p map from the skymap. The operator only works
+ * on sky maps with an identical number of layers. The content is subtracted
+ * by bi-linearily interpolating the values in the source sky map, allowing
+ * thus for a reprojection of sky map values.
+ *
+ * @todo The method is not optimized for speed as the transformation is done
+ * for each layer separately. A private method should be introduced that
+ * does the transformation, allowing the to loop more effectively over the
+ * layers.
+ ***************************************************************************/
+GSkymap& GSkymap::operator-=(const GSkymap& map)
+{
+    // Check if number of layers are identical
+    if (map.nmaps() != nmaps()) {
+        std::string msg = "Mismatch of number of maps in skymap object"
+                          " ("+gammalib::str(nmaps())+" maps in destination"
+                          " map, "+gammalib::str(map.nmaps())+" in source"
+                          " map.";
+        throw GException::invalid_value(G_OP_UNARY_SUB, msg);
+    }
+
+    // Loop over all pixels of sky map
+    for (int index = 0; index < npix(); ++index) {
+
+        // Get sky direction of actual pixel
+        GSkyDir dir = inx2dir(index);
+
+        // Loop over all layers
+        for (int layer = 0; layer < nmaps(); ++layer) {
+        
+            // Subtract value
+            (*this)(index, layer) -= map(dir, layer);
+
+        } // endfor: looped over all layers
+        
+    } // endfor: looped over all pixels
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Value subtraction operator
+ *
+ * @param[in] value Value.
+ * @return Sky map.
+ *
+ * Subtracts @p value from all sky map pixels.
+ ***************************************************************************/
+GSkymap& GSkymap::operator-=(const double& value)
+{
+    // Set total number of sky map pixels
+    int num = m_num_pixels * m_num_maps;
+
+    // Loop over all pixels of sky map
+    for (int i = 0; i < num; ++i) {
+        m_pixels[i] -= value;
+    }
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Multiplication operator
+ *
+ * @param[in] map Sky map.
+ * @return Sky map.
+ *
+ * @exception GException::invalid_value
+ *            Mismatch between number of maps in skymap object.
+ *
+ * Multiplies the content of @p map from the skymap. The operator only works
+ * on sky maps with an identical number of layers. The content is multiplied
+ * by bi-linearily interpolating the values in the source sky map, allowing
+ * thus for a reprojection of sky map values.
+ *
+ * @todo The method is not optimized for speed as the transformation is done
+ * for each layer separately. A private method should be introduced that
+ * does the transformation, allowing to loop more effectively over the
+ * layers.
+ ***************************************************************************/
+GSkymap& GSkymap::operator*=(const GSkymap& map)
+{
+    // Check if number of layers are identical
+    if (map.nmaps() != nmaps()) {
+        std::string msg = "Mismatch of number of maps in skymap object"
+                          " ("+gammalib::str(nmaps())+" maps in destination"
+                          " map, "+gammalib::str(map.nmaps())+" in source"
+                          " map.";
+        throw GException::invalid_value(G_OP_UNARY_MUL, msg);
+    }
+
+    // Loop over all pixels of sky map
+    for (int index = 0; index < npix(); ++index) {
+
+        // Get sky direction of actual pixel
+        GSkyDir dir = inx2dir(index);
+
+        // Loop over all layers
+        for (int layer = 0; layer < nmaps(); ++layer) {
+
+            // Subtract value
+            (*this)(index, layer) *= map(dir, layer);
+
+        } // endfor: looped over all layers
+
+    } // endfor: looped over all pixels
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Sky map scaling operator
+ *
+ * @param[in] factor Scale factor.
+ * @return Sky map.
+ *
+ * Multiplies all pixels of the sky map by the given scale @p factor.
+ ***************************************************************************/
+GSkymap& GSkymap::operator*=(const double& factor)
+{
+    // Compute total number of pixels
+    int n = npix() * nmaps();
+
+    // Loop over all pixels
+    double* pixel = m_pixels;
+    for (int i = 0; i < n; ++i) {
+        *pixel++ *= factor;
+    }
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Division operator
+ *
+ * @param[in] map Sky map.
+ * @return Sky map.
+ *
+ * @exception GException::invalid_value
+ *            Mismatch between number of maps in skymap object.
+ *
+ * Divides the content of @p map from the skymap. The operator only works
+ * on sky maps with an identical number of layers. The content is divided
+ * by bi-linearily interpolating the values in the source sky map, allowing
+ * thus for a reprojection of sky map values.
+ *
+ * @todo The method is not optimized for speed as the transformation is done
+ * for each layer separately. A private method should be introduced that
+ * does the transformation, allowing to loop more effectively over the
+ * layers.
+ ***************************************************************************/
+GSkymap& GSkymap::operator/=(const GSkymap& map)
+{
+    // Check if number of layers are identical
+    if (map.nmaps() != nmaps()) {
+        std::string msg = "Mismatch of number of maps in skymap object"
+                          " ("+gammalib::str(nmaps())+" maps in destination"
+                          " map, "+gammalib::str(map.nmaps())+" in source"
+                          " map.";
+        throw GException::invalid_value(G_OP_UNARY_DIV, msg);
+    }
+
+    // Loop over all pixels of sky map
+    for (int index = 0; index < npix(); ++index) {
+
+        // Get sky direction of actual pixel
+        GSkyDir dir = inx2dir(index);
+
+        // Loop over all layers
+        for (int layer = 0; layer < nmaps(); ++layer) {
+
+            // Check for division by zero
+            if (map(dir,layer) == 0.0) {
+                std::string msg = "Trying to divide by zero."
+                                  " Map entries have to be strictly non-zero"
+                                  " for division";
+                throw GException::invalid_value(G_OP_UNARY_DIV, msg);
+            }
+            
+            // Subtract value
+            (*this)(index, layer) /= map(dir, layer);
+
+        } // endfor: looped over all layers
+
+    } // endfor: looped over all pixels
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Sky map division operator
+ *
+ * @param[in] factor Scale factor.
+ * @return Sky map.
+ *
+ * @exception GException::invalid_argument
+ *            Division by zero error.
+ *
+ * Divides all pixels of the sky map by the given @p factor.
+ ***************************************************************************/
+#define G_OP_UNARY_DIV2                        "GSkymap::operator/=(double&)"
+GSkymap& GSkymap::operator/=(const double& factor)
+{
+    // Check for division by zero
+    if (factor == 0.0) {
+        std::string msg = "Trying to divide sky map pixels by zero.";
+        throw GException::invalid_argument(G_OP_UNARY_DIV2, msg);
+    }
+
+    // Compute total number of pixels
+    int n = npix() * nmaps();
+
+    // Loop over all pixels
+    double* pixel = m_pixels;
+    for (int i = 0; i < n; ++i) {
+        *pixel++ /= factor;
+    }
+
+    // Return this object
+    return *this;
+}
+
+
+/***********************************************************************//**
+ * @brief Pixel index access operator
+ *
+ * @param[in] index Pixel index [0,...,npix()-1].
+ * @param[in] map Map index [0,...,nmaps()-1].
+ * @return Sky map pixel value.
+ *
+ * @exception GException::out_of_range
+ *            Pixel index and/or map index are outside valid range.
+ *
+ * Access sky map pixel by its index, where the most quickly varying axis is
+ * the x axis of the map.
+ ***************************************************************************/
+double& GSkymap::operator()(const int& index, const int& map)
+{
+    // Throw an error if pixel index or map index is not in valid range
+    #if defined(G_RANGE_CHECK)
+    if (index < 0 || index >= m_num_pixels) {
+        throw GException::out_of_range(G_OP_ACCESS_1D,
+                                       "Sky map pixel index",
+                                       index, m_num_pixels);
+    }
+    if (map < 0 || map >= m_num_maps) {
+        throw GException::out_of_range(G_OP_ACCESS_1D,
+                                       "Sky map map index",
+                                       map, m_num_maps);
+    }
+    #endif
+
+    // Return reference to pixel value
+    return m_pixels[index+m_num_pixels*map];
+}
+
+
+/***********************************************************************//**
+ * @brief Pixel index access operator (const variant)
+ *
+ * @param[in] index Pixel index [0,...,npix()-1].
+ * @param[in] map Map index [0,...,nmaps()-1].
+ *
+ * @exception GException::out_of_range
+ *            Pixel index and/or map index are outside valid range.
+ *
+ * Access sky map pixel by its index, where the most quickly varying axis is
+ * the x axis of the map.
+ ***************************************************************************/
+const double& GSkymap::operator()(const int& index, const int& map) const
+{
+    // Throw an error if pixel index or map index is not in valid range
+    #if defined(G_RANGE_CHECK)
+    if (index < 0 || index >= m_num_pixels) {
+        throw GException::out_of_range(G_OP_ACCESS_1D,
+                                       "Sky map pixel index",
+                                       index, m_num_pixels);
+    }
+    if (map < 0 || map >= m_num_maps) {
+        throw GException::out_of_range(G_OP_ACCESS_1D,
+                                       "Sky map map index",
+                                       map, m_num_maps);
+    }
+    #endif
+
+    // Return reference to pixel value
+    return m_pixels[index+m_num_pixels*map];
+}
+
+
+/***********************************************************************//**
+ * @brief Sky map pixel access operator
+ *
+ * @param[in] pixel Sky map pixel.
+ * @param[in] map Map index [0,...,nmaps()-1].
+ *
+ * @exception GException::out_of_range
+ *            Sky pixel and/or map index are outside valid range.
+ *
+ * Access sky map pixel by its 2D index (x,y) that is implemented by the
+ * GSkyPixel class.
+ *
+ * @todo Implement proper skymap exception (actual is for matrix elements)
+ ***************************************************************************/
+double& GSkymap::operator()(const GSkyPixel& pixel, const int& map)
+{
+    // Throw an error if pixel index or map index is not in valid range
+    #if defined(G_RANGE_CHECK)
+    if (!contains(pixel)) {
+        throw GException::out_of_range(G_OP_ACCESS_2D,
+                                       int(pixel.x()), int(pixel.y()),
+                                       m_num_x-1, m_num_y-1);
+    }
+    if (map < 0 || map >= m_num_maps) {
+        throw GException::out_of_range(G_OP_ACCESS_2D,
+                                       "Sky map map index",
+                                       map, m_num_maps);
+    }
+    #endif
+
+    // Get pixel index
+    int index = pix2inx(pixel);
+
+    // Return reference to pixel value
+    return m_pixels[index+m_num_pixels*map];
+}
+
+
+/***********************************************************************//**
+ * @brief Sky map pixel access operator
+ *
+ * @param[in] pixel Sky map pixel.
+ * @param[in] map Map index [0,...,nmaps()-1].
+ *
+ * @exception GException::out_of_range
+ *            Sky pixel and/or map index are outside valid range.
+ *
+ * Access sky map pixel by its 2D index (x,y) that is implemented by the
+ * GSkyPixel class.
+ *
+ * @todo Implement proper skymap exception (actual is for matrix elements)
+ ***************************************************************************/
+const double& GSkymap::operator()(const GSkyPixel& pixel, const int& map) const
+{
+    // Throw an error if pixel index or map index is not in valid range
+    #if defined(G_RANGE_CHECK)
+    if (!contains(pixel)) {
+        throw GException::out_of_range(G_OP_ACCESS_2D,
+                                       int(pixel.x()), int(pixel.y()),
+                                       m_num_x-1, m_num_y-1);
+    }
+    if (map < 0 || map >= m_num_maps) {
+        throw GException::out_of_range(G_OP_ACCESS_2D,
+                                       "Sky map map index",
+                                       map, m_num_maps);
+    }
+    #endif
+
+    // Get pixel index
+    int index = pix2inx(pixel);
+
+    // Return reference to pixel value
+    return m_pixels[index+m_num_pixels*map];
+}
+
+
+/***********************************************************************//**
+ * @brief Return interpolated skymap value for sky direction
+ *
+ * @param[in] dir Sky direction.
+ * @param[in] map Map index [0,...,nmaps()-1].
+ * @return Sky intensity.
+ *
+ * @exception GException::out_of_range
+ *            Map index lies outside valid range.
+ *
+ * Returns the skymap value for a given sky direction, obtained by bi-linear
+ * interpolation of the neighbouring pixels. If the sky direction falls
+ * outside the area covered by the skymap, a value of 0 is returned.
+ *
+ * The method implements a computation cache that avoids recomputation of
+ * interpolation indices and weights in case that the same sky direction
+ * is requested several times. This speeds up use cases where skymap values
+ * for various map indices have to be returned for the same sky direction.
+ ***************************************************************************/
+double GSkymap::operator()(const GSkyDir& dir, const int& map) const
+{
+    // Throw an error if the map index is not in valid range
+    #if defined(G_RANGE_CHECK)
+    if (map < 0 || map >= m_num_maps) {
+        throw GException::out_of_range(G_OP_VALUE,
+                                       "Sky map map index",
+                                       map, m_num_maps);
+    }
+    #endif
+
+    // Initialise intensity
+    double intensity = 0.0;
+
+    // If the pre-computation cache has not been set or if the sky direction
+    // has changed then compute the bilinear interpolator
+    if (!m_hascache || (dir != m_last_dir)) {
+
+        // Perform computation for HealPix map
+        if (m_proj->size() == 1)  {
+
+            // Get Healpix interpolator
+            m_interpol = static_cast<GHealpix*>(m_proj)->interpolator(dir);
+
+            // Set containment flag and store actual sky direction as last
+            // direction
+            m_hascache  = true;
+            m_contained = true;
+            m_last_dir  = dir;
+
+        } // endif: we had a Healpix projection
+
+        // ... otherwise perform computation for WCS map
+        else {
+
+            // Determine sky pixel. At this point an exception may occur
+            // in case that the pixel cannot be represented by the
+            // relevant projection. We catch this exception here
+            try {
+
+                // Determine sky pixel
+                GSkyPixel pixel = dir2pix(dir);
+
+                // Set containment flag and store actual sky direction
+                // as last direction
+                m_hascache  = true;
+                m_contained = contains(pixel);
+                m_last_dir  = dir;
+
+                // Continue only if pixel is within the map
+                if (m_contained) {
+
+                    // Set left indices for interpolation. The left index is
+                    // comprised between 0 and npixels-2. By definition, the
+                    // right index is then the left index + 1
+                    int inx_x = int(pixel.x());
+                    int inx_y = int(pixel.y());
+                    if (inx_x < 0) {
+                        inx_x = 0;
+                    }
+                    else if (inx_x > m_num_x-2) {
+                        inx_x = m_num_x - 2;
+                    }
+                    if (inx_y < 0) {
+                        inx_y = 0;
+                    }
+                    else if (inx_y > m_num_y-2) {
+                        inx_y = m_num_y - 2;
+                    }
+
+                    // Set weighting factors for interpolation
+                    double wgt_x_right = (pixel.x() - inx_x);
+                    double wgt_x_left  = 1.0 - wgt_x_right;
+                    double wgt_y_right = (pixel.y() - inx_y);
+                    double wgt_y_left  = 1.0 - wgt_y_right;
+
+                    // Compute skymap pixel indices for bi-linear interpolation
+                    m_interpol.index1() = inx_x + inx_y * m_num_x;
+                    m_interpol.index2() = m_interpol.index1() + m_num_x;
+                    m_interpol.index3() = m_interpol.index1() + 1;
+                    m_interpol.index4() = m_interpol.index2() + 1;
+
+                    // Compute weighting factors for bi-linear interpolation
+                    m_interpol.weight1() = wgt_x_left  * wgt_y_left;
+                    m_interpol.weight2() = wgt_x_left  * wgt_y_right;
+                    m_interpol.weight3() = wgt_x_right * wgt_y_left;
+                    m_interpol.weight4() = wgt_x_right * wgt_y_right;
+
+                } // endif: pixel was contained in map
+                
+            } // endtry: pixel computation was successful
+            catch (GException::wcs_invalid_phi_theta) {
+                m_contained = false;
+            }
+    
+        } // endelse: we had a WCS map
+
+    } // endif: computation of the bi-linear interpolater was required
+
+    // Compute the interpolated intensity if the pixel is contained in
+    // the map
+    if (m_contained) {
+
+        // Compute map offset
+        int offset = m_num_pixels * map;
+
+        // Compute interpolated skymap value
+        intensity = m_interpol(m_pixels+offset);
+
+    } // endif: direction was contained in map
+
+    // Return intensity
+    return intensity;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                             Public methods                              =
+ =                                                                         =
+ ==========================================================================*/
+
+ /***********************************************************************//**
+ * @brief Clear instance.
+ *
+ * Resets the sky map to the initial state.
+ ***************************************************************************/
+void GSkymap::clear(void)
+{
+    // Free class members
+    free_members();
+
+    // Initialise members
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Clone sky map
+ *
+ * @return Pointer to deep copy of sky map.
+ ***************************************************************************/
+GSkymap* GSkymap::clone(void) const
+{
+    return new GSkymap(*this);
+}
+
+
+/***********************************************************************//**
+ * @brief Set number of maps
+ *
+ * @param[in] nmaps Number of maps.
+ *
+ * @exception GException::invalid_argument
+ *            Invalid number of maps specified.
+ *
+ * Redefines the number of maps in an GSkymap object. If the number of maps
+ * is increased with respect to the existing number, additional maps with
+ * pixel values of zero are append to the object. Existing map pixel values
+ * are kept. If the number of maps is decreased with respect to the existing
+ * number, the excedent maps are dropped. The remaining map pixel values are
+ * kept.
+ ***************************************************************************/
+void GSkymap::nmaps(const int& nmaps)
+{
+    // Throw an exception if less than 1 map is required
+    if (nmaps < 1) {
+        std::string msg = "At least one map is required in an GSkymap object."
+                          " Please specify an argument >=1.";
+        throw GException::invalid_argument(G_NMAPS, msg);
+    }
+
+    // If the map has pixels and if the number of maps changes then set a new
+    // number of maps. Copy over any existing information.
+    if (m_num_pixels > 0 && nmaps != m_num_maps) {
+
+        // Compute new skymap size
+        int new_size = m_num_pixels * nmaps;
+
+        // Allocate memory for new map
+        double* pixels = new double[new_size];
+
+        // Copy over existing pixels
+        int num_copy = (nmaps > m_num_maps) ? m_num_maps : nmaps;
+        num_copy    *= m_num_pixels;
+        for (int i = 0; i < num_copy; ++i) {
+            pixels[i] = m_pixels[i];
+        }
+
+        // Set any additional pixels to zero
+        if (nmaps > m_num_maps) {
+            for (int i = num_copy; i < new_size; ++i) {
+                pixels[i] = 0.0;
+            }
+        }
+
+        // Free existing pixels
+        if (m_pixels != NULL) delete [] m_pixels;
+
+        // Set pointer to stacked pixels
+        m_pixels = pixels;
+
+        // Set number of maps
+        m_num_maps = nmaps;
+
+    } // endif: map had pixels
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Converts pixel index into sky map pixel
+ *
+ * @param[in] index Pixel index [0,...,npix()-1].
+ * @return Sky map pixel.
+ *
+ * Converts the pixel @p index into a sky map pixel (GSkyPixel).The dimension
+ * of GSkyPixel will be identical to the dimension of the sky map (i.e. a 1D
+ * sky map leads to a 1D GSkyPixel object, a 2D sky map leads to a 2D
+ * GSkyPixel object.
+ ***************************************************************************/
+GSkyPixel GSkymap::inx2pix(const int& index) const
+{
+    // Initialise sky map pixel
+    GSkyPixel pixel;
+
+    // Get x and y indices
+    if (m_num_x != 0) { //!< 2D sky map
+        pixel.x(double(index % m_num_x));
+        pixel.y(double(index / m_num_x));
+    }
+    else {              //!< 1D sky map
+        pixel.index(index);
+    }
+
+    // Return pixel
+    return pixel;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns sky direction of pixel
+ *
+ * @param[in] index Pixel index [0,...,npix()-1].
+ * @return Sky direction.
+ *
+ * @exception GException::invalid_value
+ *            No valid sky projection found.
+ *
+ * Returns sky direction for a given pixel index.
+ ***************************************************************************/
+GSkyDir GSkymap::inx2dir(const int& index) const
+{
+    // Throw error if sky projection is not valid
+    if (m_proj == NULL) {
+        std::string msg = "Sky projection has not been defined.";
+        throw GException::invalid_value(G_INX2DIR, msg);
+    }
+
+    // Determine sky direction from pixel index.
+    GSkyDir dir = m_proj->pix2dir(inx2pix(index));
+
+    // Return sky direction
+    return dir;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns sky direction of pixel
+ *
+ * @param[in] pixel Sky map pixel.
+ * @return Sky direction.
+ *
+ * @exception GException::invalid_value
+ *            No valid sky projection found.
+ * @exception GException::invalid_argument
+ *            2D sky map pixel used to access 1D projection.
+ *
+ * Returns sky direction for a given sky map @p pixel.
+ ***************************************************************************/
+GSkyDir GSkymap::pix2dir(const GSkyPixel& pixel) const
+{
+    // Throw error if WCS is not valid
+    if (m_proj == NULL) {
+        std::string msg = "Sky projection has not been defined.";
+        throw GException::invalid_value(G_PIX2DIR, msg);
+    }
+
+    // Initialise sky direction
+    GSkyDir dir;
+
+    // If pixel size matches the projection size then perform a straight
+    // forward conversion
+    if (m_proj->size() == pixel.size()) {
+        dir = m_proj->pix2dir(pixel);
+    }
+
+    // ... otherwise, if we have a 2D projection but a 1D pixel then
+    // interpret the pixel as the linear index in the pixel array
+    else if (m_proj->size() == 2) {
+        dir = m_proj->pix2dir(GSkyPixel(inx2pix(int(pixel))));
+    }
+
+    // ... otherwise we have a 1D projection but a 2D pixel. There is
+    // no unambiguous way to handle this case, hence we throw an exception
+    else {
+        std::string msg = "A 2-dimensional sky map pixel "+pixel.print()+
+                          " is used to determine the sky direction for"
+                          " the 1-dimensional sky projection \""+
+                          m_proj->name()+"\"\n"
+                          "Please specify a 1-dimensional sky map pixel.";
+        throw GException::invalid_argument(G_PIX2DIR, msg);
+    }
+
+    // Return sky direction
+    return dir;
+}
+
+
+/***********************************************************************//**
+ * @brief Converts sky map pixel into pixel index
+ *
+ * @param[in] pixel Sky map pixel.
+ * @return Pixel index [0,...,npix()-1].
+ *
+ * Converts a sky map @p pixel into the pixel index.
+ ***************************************************************************/
+int GSkymap::pix2inx(const GSkyPixel& pixel) const
+{
+    // Initialise pixel index
+    int index = 0;
+
+    // Handle 1D sky map pixel
+    if (pixel.is_1D()) {
+        index = int(pixel);
+    }
+
+    // Handle 2D sky map pixel
+    else if (pixel.is_2D()) {
+
+        // Get x and y indices by rounding the (x,y) values
+        int ix = int(pixel.x()+0.5);
+        int iy = int(pixel.y()+0.5);
+
+        // Set index
+        index = ix + iy * m_num_x;
+
+    }
+
+    // Return index
+    return index;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns pixel index for a given sky direction
+ *
+ * @param[in] dir Sky direction.
+ * @return Pixel index [0,...,npix()-1].
+ *
+ * @exception GException::invalid_value
+ *            No valid sky projection found.
+ *
+ * Returns sky map pixel index for a given sky direction.
+ ***************************************************************************/
+int GSkymap::dir2inx(const GSkyDir& dir) const
+{
+    // Throw error if WCS is not valid
+    if (m_proj == NULL) {
+        std::string msg = "Sky projection has not been defined.";
+        throw GException::invalid_value(G_DIR2INX, msg);
+    }
+
+    // Determine pixel index for a given sky direction
+    int index = pix2inx(m_proj->dir2pix(dir));
+
+    // Return pixel index
+    return index;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns sky map pixel for a given sky direction
+ *
+ * @param[in] dir Sky direction.
+ * @return Sky map pixel.
+ *
+ * @exception GException::invalid_value
+ *            No valid sky projection found.
+ *
+ * Returns sky map pixel for a given sky direction.
+ ***************************************************************************/
+GSkyPixel GSkymap::dir2pix(const GSkyDir& dir) const
+{
+    // Throw error if WCS is not valid
+    if (m_proj == NULL) {
+        std::string msg = "Sky projection has not been defined.";
+        throw GException::invalid_value(G_DIR2PIX, msg);
+    }
+
+    // Determine pixel for a given sky direction
+    GSkyPixel pixel = m_proj->dir2pix(dir);
+
+    // Return pixel
+    return pixel;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns solid angle of pixel
+ *
+ * @param[in] index Pixel index [0,...,npix()-1].
+ * @return Solid angle (steradians)
+ *
+ * @exception GException::invalid_value
+ *            No valid sky projection found.
+ *
+ * Returns the solid angle of the pixel with the specified @p index.
+ ***************************************************************************/
+double GSkymap::solidangle(const int& index) const
+{
+    // Throw error if WCS is not valid
+    if (m_proj == NULL) {
+        std::string msg = "Sky projection has not been defined.";
+        throw GException::invalid_value(G_SOLIDANGLE1, msg);
+    }
+
+    // Determine solid angle from pixel index.
+    double solidangle = m_proj->solidangle(inx2pix(index));
+
+    // Return solid angle
+    return solidangle;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns solid angle of pixel
+ *
+ * @param[in] pixel Sky map pixel.
+ * @return Solid angle (steradians)
+ *
+ * @exception GException::invalid_value
+ *            No valid sky projection found.
+ *
+ * Returns the solid angle of the specified sky map @p pixel.
+ ***************************************************************************/
+double GSkymap::solidangle(const GSkyPixel& pixel) const
+{
+    // Throw error if WCS is not valid
+    if (m_proj == NULL) {
+        std::string msg = "Sky projection has not been defined.";
+        throw GException::invalid_value(G_SOLIDANGLE2, msg);
+    }
+
+    // Initialise solid angle
+    double solidangle = 0.0;
+
+    // If pixel size matches the projection size then perform a straight
+    // forward solid angle determination
+    if (m_proj->size() == pixel.size()) {
+        solidangle = m_proj->solidangle(pixel);
+    }
+
+    // ... otherwise, if we have a 2D projection but a 1D pixel then
+    // interpret the pixel as the linear index in the pixel array
+    else if (m_proj->size() == 2) {
+        solidangle = m_proj->solidangle(GSkyPixel(inx2pix(int(pixel))));
+    }
+
+    // ... otherwise we have a 1D projection but a 2D pixel. There is
+    // no unambiguous way to handle this case, hence we throw an exception
+    else {
+        std::string msg = "A 2-dimensional sky map pixel "+pixel.print()+
+                          " is used to determine the solid angle for"
+                          " the 1-dimensional sky projection \""+
+                          m_proj->name()+"\"\n"
+                          "Please specify a 1-dimensional sky map pixel.";
+        throw GException::invalid_argument(G_SOLIDANGLE2, msg);
+    }
+
+    // Return solid angle
+    return solidangle;
+}
+
+
+/***********************************************************************//**
+ * @brief Set sky projection
+ *
+ * @param[in] proj Sky projection.
+ *
+ * Sets the projection from celestial to pixel coordinates. The method
+ * performs a deep copy of @p proj, allowing to destroy the argument after
+ * using the method.
+ *
+ * Warning: this method may corrupt the GSkymap object as it allows assigning
+ * for example a 1D projection to a 2D skymap. Please use this method only
+ * when you know what you're doing.
+ *
+ * @todo We may restrict this method to not allow changing the projection
+ * dimension.
+ ***************************************************************************/
+void GSkymap::projection(const GSkyProjection& proj)
+{
+    // Free any existing WCS
+    if (m_proj != NULL) delete m_proj;
+
+    // Clone input WCS
+    m_proj = proj.clone();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Verifies if sky direction falls in map
+ *
+ * @param[in] dir Sky direction.
+ *
+ * This method checks if the specified sky direction falls within the pixels
+ * covered by the skymap. The method uses the dir2xy method to convert the
+ * sky direction into 2D pixel indices, and then checks whether the pixel
+ * indices fall in the skymap.
+ ***************************************************************************/
+bool GSkymap::contains(const GSkyDir& dir) const
+{
+    // Convert sky direction into sky pixel
+    GSkyPixel pixel = dir2pix(dir);
+    
+    // Return location flag
+    return (contains(pixel));
+}
+
+
+/***********************************************************************//**
+ * @brief Checks if sky map pixel falls in map
+ *
+ * @param[in] pixel Sky map pixel.
+ * @return Trus if pixels is within map, false otherwise.
+ *
+ * Checks whether the specified sky map @p pixel falls within the skymap
+ * or not. 
+ ***************************************************************************/
+bool GSkymap::contains(const GSkyPixel& pixel) const
+{
+    // Initialise containment flag
+    bool inmap = false;
+
+    // Test 1D pixels
+    if (pixel.is_1D()) {
+        if (pixel.index()+0.5 >= 0.0 && pixel.index()-0.5 < m_num_pixels) {
+            inmap = true;
+        }
+    }
+
+    // Test 2D pixels
+    else if (pixel.is_2D()) {
+
+        // If pixel is in range then set containment flag to true
+        if ((pixel.x()+0.5 >= 0.0 && pixel.x()-0.5 < m_num_x) &&
+            (pixel.y()+0.5 >= 0.0 && pixel.y()-0.5 < m_num_y)) {
+            inmap = true;
+        }
+
+    }
+
+    // Return containment flag
+    return inmap;
+}
+
+
+/***********************************************************************//**
+ * @brief Extract maps into a new sky map object
+ *
+ * @param[in] map First map to extract
+ * @param[in] nmaps Number of maps to extract (default: 1)
+ * @return Extracted map(s).
+ *
+ * @exception GException::out_of_range
+ *            First map index outside valid range
+ * @exception GException::invalid_argument
+ *            Requested number of maps are not available.
+ *
+ * Extracts @p nmaps sky maps starting from @p map from the sky map.
+ ***************************************************************************/
+GSkymap GSkymap::extract(const int& map, const int& nmaps) const
+{
+    // Throw an exception if the first map index is invalid
+    if (map < 0 || map >= m_num_maps) {
+        throw GException::out_of_range(G_EXTRACT, "Sky map index", map,
+                                       m_num_maps);
+    }
+ 
+    // Throw an exception if the number of maps is invalid
+    if (nmaps < 0) {
+        std::string msg = "The number of maps to extract cannot be negative "
+                          "(nmaps="+gammalib::str(nmaps)+").";
+        throw GException::invalid_argument(G_EXTRACT, msg);
+    }
+    if (nmaps > m_num_maps-map) {
+        std::string msg = "The number of maps to extract ("+
+                          gammalib::str(nmaps)+") exceeds the number of maps "
+                          "that are available ("+
+                          gammalib::str(m_num_maps-map)+").";
+        throw GException::invalid_argument(G_EXTRACT, msg);
+    }
+
+    // Compute memory size for extracted map
+    int n_size = m_num_pixels * nmaps;
+
+    // Allocate memory for extracted maps (handle the case that the
+    // extracted map can be empty)
+    double* pixels = NULL;
+    if (n_size > 0) {
+        pixels = new double[n_size];
+    }
+
+    // Extract pixels
+    double *src = m_pixels + map*m_num_pixels;
+    double *dst = pixels;
+    for (int i = 0; i < n_size; ++i) {
+        *dst++ = *src++;
+    }
+
+    // Create a copy of the map
+    GSkymap result = *this;
+
+    // Delete pixels from that map
+    if (result.m_pixels != NULL) delete [] result.m_pixels;
+
+    // Attach copied pixels to the map
+    result.m_pixels = pixels;
+
+    // Set number of maps
+    result.m_num_maps = nmaps;
+
+    // Return map
+    return result;
+}
+
+
+/***********************************************************************//**
+ * @brief Stack all maps into a single map
+ *
+ * @return Stacked map.
+ *
+ * The methods replaces the sky map by a version with only a single map
+ * by summing over the pixel values for all maps. If the sky map has no
+ * pixels or there is only a single map in the object, the method does
+ * nothing.
+ ***************************************************************************/
+void GSkymap::stack_maps(void)
+{
+    // Continue only if the map has pixels and if there is more than 1 map
+    if (m_num_pixels > 0 && m_num_maps > 1) {
+
+        // Allocate memory for stacked map
+        double* pixels = new double[m_num_pixels];
+
+        // Stack map and save in memory
+        for (int i = 0; i < m_num_pixels; ++i) {
+            double sum = 0.0;
+            for (int k = 0; k < m_num_maps; ++k) {
+                sum += (*this)(i,k);
+            }
+            pixels[i] = sum;
+        }
+
+        // Free existing pixels
+        if (m_pixels != NULL) delete [] m_pixels;
+
+        // Set pointer to stacked pixels
+        m_pixels = pixels;
+
+        // Set number of maps to 1
+        m_num_maps = 1;
+    
+    } // endif: map had pixels
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Load skymap from FITS file.
+ *
+ * @param[in] filename FITS file name..
+ *
+ * Loads HEALPix and non HEALPix skymaps. First searches for HEALPix map in
+ * FITS file by scanning all HDUs for PIXTYPE=HEALPIX. If no HEALPix map has
+ * been found then search load first non-empty image.
+ *
+ * @todo Do we have to restrict a HEALPix map to a BinTable and a WCS map
+ * to a Double precision image???
+ ***************************************************************************/
+void GSkymap::load(const std::string& filename)
+{
+    // Free memory and initialise members
+    free_members();
+    init_members();
+
+    // Open FITS file
+    GFits fits(filename);
+
+    // Get number of HDUs
+    int num = fits.size();
+
+    // Initialize load flag
+    bool loaded = false;
+
+    // First search for HEALPix extension. We can skip the first extension
+    // since this is always an image and a HEALPix map is stored in a
+    // binary table
+    for (int extno = 1; extno < num; ++extno) {
+
+        // Get reference to HDU
+        const GFitsHDU& hdu = *fits.at(extno);
+        
+        // If PIXTYPE keyword equals "HEALPIX" then load map
+        if (hdu.has_card("PIXTYPE") && hdu.string("PIXTYPE") == "HEALPIX") {
+            read_healpix(static_cast<const GFitsTable&>(hdu));
+            loaded = true;
+            break;
+        }
+
+    } // endfor: looped over HDUs
+
+    // If we have not found a HEALPIX map then search now for image.
+    // Skip empty images
+    if (!loaded) {
+        for (int extno = 0; extno < num; ++extno) {
+
+            // Get referene to HDU
+            const GFitsHDU& hdu = *fits.at(extno);
+
+            // Skip if extension is not an image
+            if (extno > 0) {
+                if (hdu.string("XTENSION") != "IMAGE")
+                    continue;
+            }
+
+            // Load WCS map
+            read_wcs(static_cast<const GFitsImage&>(hdu));
+            //loaded = true;
+            break;
+
+        } // endfor: looped over HDUs
+    } // endif: no HEALPix map found
+
+    // Close FITS file
+    fits.close();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Save skymap into FITS file.
+ *
+ * @param[in] filename FITS file name.
+ * @param[in] clobber Overwrite existing file? (true=yes)
+ *
+ * The method does nothing if the skymap holds no valid WCS.
+ ***************************************************************************/
+void GSkymap::save(const std::string& filename, bool clobber) const
+{
+    // Continue only if we have data to save
+    if (m_proj != NULL) {
+
+        // Initialise HDU pointer
+        GFitsHDU* hdu = NULL;
+
+        // Case A: Skymap is Healpix
+        if (m_proj->code() == "HPX") {
+            hdu = create_healpix_hdu();
+        }
+
+        // Case B: Skymap is not Healpix
+        else {
+            hdu = create_wcs_hdu();
+        }
+
+        // Create FITS file and save it to disk
+        if (hdu != NULL) {
+            GFits fits;
+            fits.append(*hdu);
+            fits.saveto(filename, clobber);
+        }
+
+        // Delete HDU
+        if (hdu != NULL) delete hdu;
+
+    } // endif: we had data to save
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Read skymap from FITS HDU
+ *
+ * @param[in] hdu FITS HDU.
+ ***************************************************************************/
+void GSkymap::read(const GFitsHDU& hdu)
+{
+    // Free memory and initialise members
+    free_members();
+    init_members();
+
+    // Initialize load flag
+    bool loaded = false;
+
+    // If PIXTYPE keyword equals "HEALPIX" then load map
+    if (hdu.has_card("PIXTYPE") && hdu.string("PIXTYPE") == "HEALPIX") {
+        read_healpix(static_cast<const GFitsTable&>(hdu));
+        loaded = true;
+    }
+
+    // ... otherwise try loading as non HEALPix map
+    if (!loaded) {
+
+        // Load only if HDU contains an image
+        if (hdu.exttype() == 0) {
+            read_wcs(static_cast<const GFitsImage&>(hdu));
+            //loaded = true;
+        }
+
+    } // endif
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Write skymap into FITS file
+ *
+ * @param[in] file FITS file pointer.
+ ***************************************************************************/
+void GSkymap::write(GFits& file) const
+{
+    // Continue only if we have data to save
+    if (m_proj != NULL) {
+
+        // Initialise HDU pointer
+        GFitsHDU* hdu = NULL;
+
+        // Case A: Skymap is Healpix
+        if (m_proj->code() == "HPX") {
+            hdu = create_healpix_hdu();
+        }
+
+        // Case B: Skymap is not Healpix
+        else {
+            hdu = create_wcs_hdu();
+        }
+
+        // Append HDU to FITS file.
+        if (hdu != NULL) {
+            file.append(*hdu);
+        }
+
+        // Delete HDU
+        if (hdu != NULL) delete hdu;
+
+    } // endif: we had data to save
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Print sky map
+ *
+ * @param[in] chatter Chattiness (defaults to NORMAL).
+ * @return String containing sky map information.
+ ***************************************************************************/
+std::string GSkymap::print(const GChatter& chatter) const
+{
+    // Initialise result string
+    std::string result;
+
+    // Continue only if chatter is not silent
+    if (chatter != SILENT) {
+
+        // Append header
+        result.append("=== GSkymap ===");
+
+        // Append information
+        result.append("\n"+gammalib::parformat("Number of pixels"));
+        result.append(gammalib::str(m_num_pixels));
+        result.append("\n"+gammalib::parformat("Number of maps"));
+        result.append(gammalib::str(m_num_maps));
+        if (m_proj != NULL && m_proj->size() == 2) {
+            result.append("\n"+gammalib::parformat("X axis dimension"));
+            result.append(gammalib::str(m_num_x));
+            result.append("\n"+gammalib::parformat("Y axis dimension"));
+            result.append(gammalib::str(m_num_y));
+        }
+
+        // Append sky projection information
+        if (gammalib::reduce(chatter) > SILENT) {
+            if (m_proj != NULL) {
+                result.append("\n"+m_proj->print(gammalib::reduce(chatter)));
+            }
+            else {
+                result.append("\n"+gammalib::parformat("Sky projection"));
+                result.append("not defined");
+            }
+        }
+
+    } // endif: chatter was not silent
+
+    // Return result
+    return result;
+}
+
+/*==========================================================================
+ =                                                                         =
+ =                             Private methods                             =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Initialise class members
+ ***************************************************************************/
+void GSkymap::init_members(void)
+{
+    // Initialise members
+    m_num_pixels = 0;
+    m_num_maps   = 0;
+    m_num_x      = 0;
+    m_num_y      = 0;
+    m_proj       = NULL;
+    m_pixels     = NULL;
+
+    // Initialise computation cache
+    m_hascache  = false;
+    m_contained = false;
+    m_last_dir.clear();
+    m_interpol.clear();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Allocate skymap pixels
+ ***************************************************************************/
+void GSkymap::alloc_pixels(void)
+{
+    // Compute data size
+    int size = m_num_pixels * m_num_maps;
+
+    // Continue only if there are pixels
+    if (size > 0) {
+
+        // Allocate pixels and initialize them to 0
+        m_pixels = new double[size];
+        for (int i = 0; i < size; ++i) {
+            m_pixels[i] = 0.0;
+        }
+
+    } // endif: there were pixels
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy class members
+ *
+ * @param[in] map Sky map.
+ ***************************************************************************/
+void GSkymap::copy_members(const GSkymap& map)
+{
+    // Copy attributes
+    m_num_pixels = map.m_num_pixels;
+    m_num_maps   = map.m_num_maps;
+    m_num_x      = map.m_num_x;
+    m_num_y      = map.m_num_y;
+
+    // Copy computation cache
+    m_hascache  = map.m_hascache;
+    m_contained = map.m_contained;
+    m_last_dir  = map.m_last_dir;
+    m_interpol  = map.m_interpol;
+
+    // Clone sky projection if it is valid
+    if (map.m_proj != NULL) m_proj = map.m_proj->clone();
+
+    // Compute data size
+    int size = m_num_pixels * m_num_maps;
+
+    // Copy pixels
+    if (size > 0 && map.m_pixels != NULL) {
+        alloc_pixels();
+        for (int i = 0; i < size; ++i) {
+            m_pixels[i] = map.m_pixels[i];
+        }
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Delete class members
+ ***************************************************************************/
+void GSkymap::free_members(void)
+{
+    // Free memory
+    if (m_proj   != NULL) delete m_proj;
+    if (m_pixels != NULL) delete [] m_pixels;
+
+    // Signal free pointers
+    m_proj       = NULL;
+    m_pixels     = NULL;
+
+    // Reset number of pixels
+    m_num_pixels = 0;
+    m_num_maps   = 0;
+    m_num_x      = 0;
+    m_num_y      = 0;
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Set World Coordinate System
+ *
+ * @param[in] wcs World Coordinate System code.
+ * @param[in] coords Coordinate system.
+ * @param[in] crval1 X value of reference pixel.
+ * @param[in] crval2 Y value of reference pixel.
+ * @param[in] crpix1 X index of reference pixel.
+ * @param[in] crpix2 Y index of reference pixel.
+ * @param[in] cdelt1 Increment in x direction at reference pixel (deg).
+ * @param[in] cdelt2 Increment in y direction at reference pixel (deg).
+ * @param[in] cd Astrometry parameters (2x2 matrix, deg/pixel).
+ * @param[in] pv2 Projection parameters (length WCS type dependent).
+ *
+ * @exception GException::wcs_invalid
+ *            Invalid wcs parameter (World Coordinate System not known).
+ *
+ * This method sets the WCS projection pointer based on the WCS code and
+ * sky map parameters. It makes use of the GWcsRegistry class to allocate
+ * the correct derived class. Note that this method does not support the
+ * HPX projection.
+ *
+ * @todo Remove cd and pv2 parameters.
+ ***************************************************************************/
+void GSkymap::set_wcs(const std::string& wcs, const std::string& coords,
+                      const double& crval1, const double& crval2,
+                      const double& crpix1, const double& crpix2,
+                      const double& cdelt1, const double& cdelt2,
+                      const GMatrix& cd, const GVector& pv2)
+{
+    // Convert WCS to upper case
+    std::string uwcs = gammalib::toupper(wcs);
+    
+    // Check if HPX was requested (since this is not allowed)
+    if (uwcs == "HPX") {
+       throw GException::wcs_invalid(G_SET_WCS, uwcs,
+                                     "Method not valid for HPX projection.");
+    }
+    
+    // ... otherwise get projection from registry
+    else {
+        // Allocate WCS registry
+        GWcsRegistry registry;
+        
+        // Allocate projection from registry
+        GWcs* projection = registry.alloc(uwcs);
+        m_proj           = projection;
+
+        // Signal if projection type is not known
+        if (projection == NULL) {
+            std::string message = "Projection code not known. "
+                                  "Should be one of "+registry.list()+".";
+            throw GException::wcs_invalid(G_SET_WCS, uwcs, message);
+        }
+
+        // Setup WCS
+        projection->set(coords, crval1, crval2, crpix1, crpix2,
+                        cdelt1, cdelt2);
+    
+    } // endelse: got projection from registry
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Read Healpix data from FITS table.
+ *
+ * @param[in] table FITS table.
+ *
+ * HEALPix data may be stored in various formats depending on the 
+ * application that has writted the data. HEALPix IDL, for example, may
+ * store the data in vectors of length 1024 if the number of pixels is
+ * a multiple of 1024. On the other hand, vectors may also be used to store
+ * several HEALPix maps into a single column. Alternatively, multiple maps
+ * may be stored in multiple columns.
+ ***************************************************************************/
+void GSkymap::read_healpix(const GFitsTable& table)
+{
+    // Determine number of rows and columns in table
+    int nrows = table.nrows();
+    int ncols = table.ncols();
+    #if defined(G_READ_HEALPIX_DEBUG)
+    std::cout << "nrows=" << nrows << " ncols=" << ncols << std::endl;
+    #endif
+
+    // Allocate Healpix projection
+    m_proj = new GHealpix;
+
+    // Read WCS information from FITS header
+    m_proj->read(table);
+
+    // Set number of pixels based on NSIDE parameter
+    m_num_pixels = static_cast<GHealpix*>(m_proj)->npix();
+    #if defined(G_READ_HEALPIX_DEBUG)
+    std::cout << "m_num_pixels=" << m_num_pixels << std::endl;
+    #endif
+
+    // Number of map pixels has to be a multiple of the number of
+    // rows in column
+    if (m_num_pixels % nrows != 0) {
+        throw GException::skymap_bad_size(G_READ_HEALPIX, nrows,
+                                          m_num_pixels);
+    }
+
+    // Determine vector length for HEALPix data storage
+    int nentry = m_num_pixels / nrows;
+    #if defined(G_READ_HEALPIX_DEBUG)
+    std::cout << "nentry=" << nentry << std::endl;
+    #endif
+
+    // Determine number of maps from the number of maps that fit into
+    // all columns. Only count columns that can fully hold the map.
+    m_num_maps = 0;
+    for (int icol = 0; icol < ncols; ++icol) {
+        const GFitsTableCol* col = table[icol];
+        if (col->number() % nentry == 0) {
+            m_num_maps += col->number() / nentry;
+        }
+    }
+    #if defined(G_READ_HEALPIX_DEBUG)
+    std::cout << "m_num_maps=" << m_num_maps << std::endl;
+    #endif
+
+    // Allocate pixels to hold the map
+    alloc_pixels();
+
+    // Initialise map counter
+    int imap = 0;
+
+    // Loop over all columns
+    for (int icol = 0; icol < ncols; ++icol) {
+
+        // Get next column
+        const GFitsTableCol* col = table[icol];
+
+        // Only consider columns that can fully hold maps
+        if (col->number() % nentry == 0) {
+
+            // Determine number of maps in column
+            int num = col->number() / nentry;
+
+            // Loop over all maps in column
+            int inx_start = 0;
+            int inx_end   = nentry;
+            for (int i = 0; i < num; ++i) {
+
+                // Load map
+                double *ptr = m_pixels + m_num_pixels*imap;
+                for (int row = 0; row < col->length(); ++row) {
+                    for (int inx = inx_start; inx < inx_end; ++inx) {
+                        *ptr++ = col->real(row,inx);
+                    }
+                }
+                #if defined(G_READ_HEALPIX_DEBUG)
+                std::cout << "Load map=" << imap << " index="
+                          << inx_start << "-" << inx_end << std::endl;
+                #endif
+
+                // Increment index range
+                inx_start  = inx_end;
+                inx_end   += nentry;
+
+                // Increment map counter
+                imap++;
+
+                // Break if we have loaded all maps
+                if (imap >= m_num_maps) {
+                    break;
+                }
+
+            } // endfor: looped over all maps in column
+        } // endif: column could fully hold maps
+
+        // Break if we have loaded all maps
+        if (imap >= m_num_maps) {
+            break;
+        }
+
+    } // endfor: looped over all columns
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Read WCS image from FITS HDU
+ *
+ * @param[in] image FITS image.
+ *
+ * @exception GException::skymap_bad_image_dim
+ *            WCS image has invalid dimension (naxis=2 or 3).
+ ***************************************************************************/
+void GSkymap::read_wcs(const GFitsImage& image)
+{
+    // Allocate WCS
+    alloc_wcs(image);
+
+    // Read projection information from FITS header
+    m_proj->read(image);
+
+    // Extract map dimension and number of maps from image
+    if (image.naxis() == 2) {
+        m_num_x    = image.naxes(0);
+        m_num_y    = image.naxes(1);
+        m_num_maps = 1;
+    }
+    else if (image.naxis() >= 3) {
+        m_num_x    = image.naxes(0);
+        m_num_y    = image.naxes(1);
+        m_num_maps = image.naxes(2);
+    }
+    else {
+        throw GException::skymap_bad_image_dim(G_READ_WCS, image.naxis());
+    }
+    #if defined(G_READ_WCS_DEBUG)
+    std::cout << "m_num_x=" << m_num_x << std::endl;
+    std::cout << "m_num_y=" << m_num_y << std::endl;
+    std::cout << "m_num_maps=" << m_num_maps << std::endl;
+    #endif
+
+    // Compute number of pixels
+    m_num_pixels = m_num_x * m_num_y;
+    #if defined(G_READ_WCS_DEBUG)
+    std::cout << "m_num_pixels=" << m_num_pixels << std::endl;
+    #endif
+
+    // Allocate pixels to hold the map
+    alloc_pixels();
+
+    // Read image
+    if (image.naxis() == 2) {
+        double* ptr = m_pixels;
+        for (int iy = 0; iy < m_num_y; ++iy) {
+            for (int ix = 0; ix < m_num_x; ++ix) {
+                *ptr++ = image.pixel(ix,iy);
+            }
+        }
+    }
+    else {
+        double* ptr = m_pixels;
+        for (int imap = 0; imap < m_num_maps; ++imap) {
+            for (int iy = 0; iy < m_num_y; ++iy) {
+                for (int ix = 0; ix < m_num_x; ++ix) {
+                    *ptr++ = image.pixel(ix,iy,imap);
+                }
+            }
+        }
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Allocate WCS class
+ *
+ * @param[in] image FITS image.
+ *
+ * @exception GException::fits_key_not_found
+ *            Unable to find required FITS header keyword.
+ * @exception GException::skymap_bad_ctype
+ *            CTYPE1 and CTYPE2 keywords are incompatible.
+ * @exception GException::wcs_invalid
+ *            WCS projection of FITS file not supported by GammaLib.
+ ***************************************************************************/
+void GSkymap::alloc_wcs(const GFitsImage& image)
+{
+    // Get standard keywords
+    std::string ctype1 = image.string("CTYPE1");
+    std::string ctype2 = image.string("CTYPE2");
+
+    // Extract projection type
+    std::string xproj = ctype1.substr(5,3);
+    std::string yproj = ctype2.substr(5,3);
+
+    // Check that projection type is identical on both axes
+    if (xproj != yproj) {
+        throw GException::skymap_bad_ctype(G_ALLOC_WCS,
+                                           ctype1, ctype2);
+    }
+
+    // Allocate WCS registry
+    GWcsRegistry registry;
+    
+    // Allocate projection from registry
+    m_proj = registry.alloc(xproj);
+    
+    // Signal if projection type is not known
+    if (m_proj == NULL) {
+        std::string message = "Projection code not known. "
+                              "Should be one of "+registry.list()+".";
+        throw GException::wcs_invalid(G_ALLOC_WCS, xproj, message);
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Create FITS HDU containing Healpix data
+ *
+ * This method allocates a binary table HDU that contains the Healpix data.
+ * Deallocation of the table has to be done by the client.
+ ***************************************************************************/
+GFitsBinTable* GSkymap::create_healpix_hdu(void) const
+{
+    // Initialise result to NULL pointer
+    GFitsBinTable* hdu = NULL;
+
+    // Compute size of Healpix data
+    int size = m_num_pixels * m_num_maps;
+
+    // Continue only if we have pixels
+    if (size > 0) {
+
+        // Set number of rows and columns
+        int rows   = m_num_pixels;
+        int number = m_num_maps;
+
+        // Create column to hold Healpix data
+        GFitsTableDoubleCol column = GFitsTableDoubleCol("DATA", rows, number);
+
+        // Fill data into column
+        double* ptr = m_pixels;
+        for (int inx = 0; inx < number; ++inx) {
+            for (int row = 0; row < rows; ++row) {
+                column(row,inx) = *ptr++;
+            }
+        }
+
+        // Create HDU that contains Healpix map in a binary table
+        hdu = new GFitsBinTable(rows);
+        hdu->append(column);
+
+    } // endif: there were pixels
+
+    // ... otherwise create an empty header
+    else {
+        hdu = new GFitsBinTable;
+    }
+
+    // Set extension name
+    hdu->extname("HEALPIX");
+
+    // If we have WCS information then write into FITS header
+    if (m_proj != NULL) m_proj->write(*hdu);
+
+    // Set additional keywords
+    hdu->card("NBRBINS", m_num_maps, "Number of HEALPix maps");
+
+    // Return HDU
+    return hdu;
+}
+
+
+/***********************************************************************//**
+ * @brief Create FITS HDU containing WCS image
+ *
+ * This method allocates an image HDU that contains the WCS image data.
+ * Deallocation of the image has to be done by the client.
+ *
+ * @todo Set additional keywords.
+ ***************************************************************************/
+GFitsImageDouble* GSkymap::create_wcs_hdu(void) const
+{
+    // Initialise result to NULL pointer
+    GFitsImageDouble* hdu = NULL;
+
+    // Compute size of Healpix data
+    int size = m_num_pixels * m_num_maps;
+
+    // Continue only if we have pixels
+    if (size > 0) {
+
+        // Set axis parameters for image construction
+        int naxis   = (m_num_maps == 1) ? 2 : 3;
+        int naxes[] = {m_num_x, m_num_y, m_num_maps};
+
+        // Allocate image
+        hdu = new GFitsImageDouble(naxis, naxes);
+
+        // Store data in image
+        if (naxis == 2) {
+            double* ptr = m_pixels;
+            for (int iy = 0; iy < m_num_y; ++iy) {
+                for (int ix = 0; ix < m_num_x; ++ix) {
+                    (*hdu)(ix,iy) = *ptr++;
+                }
+            }
+        }
+        else {
+            double* ptr = m_pixels;
+            for (int imap = 0; imap < m_num_maps; ++imap) {
+                for (int iy = 0; iy < m_num_y; ++iy) {
+                    for (int ix = 0; ix < m_num_x; ++ix) {
+                        (*hdu)(ix,iy,imap) = *ptr++;
+                    }
+                }
+            }
+        }
+
+    } // endif: there were pixels
+
+    // ... otherwise create an empty header
+    else {
+        hdu = new GFitsImageDouble;
+    }
+
+    // Set extension name
+    hdu->extname("IMAGE");
+
+    // If we have sky projection information then write into FITS header
+    if (m_proj != NULL) m_proj->write(*hdu);
+
+    // Set additional keywords
+    //TODO
+
+    // Return HDU
+    return hdu;
+}
+
+/***********************************************************************//**
+ * @brief publishes the skymap in a VO Hub
+ *
+ * @param[in] 
+ *
+ * @exception GException::
+ *            
+ ***************************************************************************/
+void GSkymap::publish(void)
+{
+    make_votable();
+    // Return
+    return;
+}
+
+/***********************************************************************//**
+ * @brief makes GVOTable object and fills ot with current skymap(s)
+ *
+ * @param[in] 
+ *
+ * @exception GException::
+ *            
+ ***************************************************************************/
+void GSkymap::make_votable(void)
+{
+    GVOTable* m_skytable;
+    m_skytable->create_votable();
+    //For each pixel set the value of the pixel inside the data part
+    m_skytable->fill_content("Test value");
+    // Return
+    return;
+}
+/*==========================================================================
+ =                                                                         =
+ =                                 Friends                                 =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Computes square root of skymap elements
+ *
+ * @param[in] map Skymap.
+ * @return Skymap containing the square root of every element.
+ ***************************************************************************/
+GSkymap sqrt(const GSkymap& map)
+{
+    // Initialise result vector
+    GSkymap result(map);
+
+    // Loop over all maps
+    for (int i = 0; i < map.nmaps(); ++i) {
+
+        // Loop over all bins
+        for (int j = 0; j < map.npix(); ++j) {
+
+            // Get the content from the bin
+            double content = map(j,i);
+
+    		// Check if content is not negative
+    		if (content < 0.0) {
+                std::string msg = "Negative value envountered."
+                                  " Cannot take the sqrt from a negative value";
+                throw GException::invalid_value(G_OP_UNARY_DIV, msg);
+            }
+
+            // Set content of the result map
+            result(j,i) = std::sqrt(content);
+
+        } // endfor: Loop over all bins
+
+    } // endfor: Loop over all maps
+
+    // Return vector
+    return result;
+}
+/***************************************************************************
+ *                     GVOApp.hpp - VO SAMP Hub class                      *
+ * ----------------------------------------------------------------------- *
+ *  copyright (C) 2014 by Thierry Louge                                    *
+ * ----------------------------------------------------------------------- *
+ *                                                                         *
+ *  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   *
+ *  the Free Software Foundation, either version 3 of the License, or      *
+ *  (at your option) any later version.                                    *
+ *                                                                         *
+ *  This program is distributed in the hope that it will be useful,        *
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of         *
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          *
+ *  GNU General Public License for more details.                           *
+ *                                                                         *
+ *  You should have received a copy of the GNU General Public License      *
+ *  along with this program.  If not, see <http://www.gnu.org/licenses/>.  *
+ *                                                                         *
+ ***************************************************************************/
+/**
+ * @file GVOApp.cpp
+ * @brief VO  Application connected to SAMP HUB class implementation
+ * @author Thierry Louge
+ */
+
+/* __ Includes ___________________________________________________________ */
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#include <errno.h>
+#endif
+#include <cstdlib>         // std::getenv() function
+#include <cstdio>          // std::fopen(), etc. functions
+#include <cstring>         // std::memset() function
+#include <unistd.h>        // close() function
+#include "GVOApp.hpp"
+#include "GException.hpp"
+#include "GTools.hpp"
+
+/* __ Method name definitions ____________________________________________ */
+
+/* __ Macros _____________________________________________________________ */
+
+/* __ Coding definitions _________________________________________________ */
+
+/* __ Debug definitions __________________________________________________ */
+
+
+/*==========================================================================
+ =                                                                         =
+ =                         Constructors/destructors                        =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Void constructor
+ ***************************************************************************/
+GVOApp::GVOApp(void)
+{
+    // Initialise members
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy constructor
+ *
+ * @param[in] hub VO Hub.
+ ***************************************************************************/
+GVOApp::GVOApp(const GVOApp& hub)
+{
+    // Initialise members
+    init_members();
+
+    // Copy members
+    copy_members(hub);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Destructor
+ ***************************************************************************/
+GVOApp::~GVOApp(void)
+{
+    
+    // Return
+    return;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                                Operators                                =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Assignment operator
+ *
+ * @param[in] hub VO hub.
+ * @return VO hub.
+ ***************************************************************************/
+GVOApp& GVOApp::operator=(const GVOApp& hub)
+{
+    // Execute only if object is not identical
+    if (this != &hub) {
+
+        // Free members
+        free_members();
+
+        // Initialise members
+        init_members();
+
+        // Copy members
+        copy_members(hub);
+
+    } // endif: object was not identical
+
+    // Return
+    return *this;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                             Public methods                              =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Clear object.
+ *
+ * Reset object to a clean initial state.
+ ***************************************************************************/
+void GVOApp::clear(void)
+{
+    // Free memory and initialise members
+    free_members();
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Clone object
+ ***************************************************************************/
+GVOApp* GVOApp::clone(void) const
+{
+    // Clone client
+    return new GVOApp(*this);
+}
+
+
+/***********************************************************************//**
+ * @brief TBD
+ *
+ ***************************************************************************/
+void GVOApp::start(void)
+{
+    // Does nothing
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Print VO hub information
+ *
+ * @param[in] chatter Chattiness (defaults to NORMAL).
+ * @return String containing VO hub information
+ ***************************************************************************/
+std::string GVOApp::print(const GChatter& chatter) const
+{
+    // Initialise result string
+    std::string result;
+    std::list<std::string> temp;
+    temp = registered_calls;
+    
+    // Continue only if chatter is not silent
+    if (chatter != SILENT) {
+
+        // Append header
+        result.append("=== GVOApp ===");
+	result.append(m_id);
+	result.append("%d",temp.size());
+	for (std::list<std::string>::iterator i = temp.begin(); i != temp.end(); ++i) {
+	  result.append("Method:\n");
+	  result.append(*i);
+	}
+	
+      
+
+    } // endif: chatter was not silent
+
+    // Return result
+    return result;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                             Private methods                             =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Initialise class members
+ ***************************************************************************/
+void GVOApp::init_members(void)
+{
+    // Initialise members
+    std::list<std::string> registeredMethods;
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy class members
+ *
+ * @param[in] app registered in Hub.
+ ***************************************************************************/
+void GVOApp::copy_members(const GVOApp& app)
+{
+    // Copy members
+    
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Delete class members
+ ***************************************************************************/
+void GVOApp::free_members(void)
+{
+    // Return
+    return;
+}
+
+/***********************************************************************//**
+ * @brief Creates nex application container
+ ***************************************************************************/
+void GVOApp::start_app(void)
+{
+    // Creates Hub
+    printf("Creating new application Container\n");
+
+    // Return
+    return;
+}
+/***************************************************************************
+ *                     GVOClient.hpp - VO client class                     *
+ * ----------------------------------------------------------------------- *
+ *  copyright (C) 2013 by 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   *
+ *  the Free Software Foundation, either version 3 of the License, or      *
+ *  (at your option) any later version.                                    *
+ *                                                                         *
+ *  This program is distributed in the hope that it will be useful,        *
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of         *
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          *
+ *  GNU General Public License for more details.                           *
+ *                                                                         *
+ *  You should have received a copy of the GNU General Public License      *
+ *  along with this program.  If not, see <http://www.gnu.org/licenses/>.  *
+ *                                                                         *
+ ***************************************************************************/
+/**
+ * @file GVOClient.cpp
+ * @brief VO client class implementation
+ * @author Juergen Knoedlseder
+ */
+
+/* __ Includes ___________________________________________________________ */
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+#include <cstdlib>         // std::getenv() function
+#include <cstdio>          // std::fopen(), etc. functions
+#include <cstring>         // std::memset() function
+#include <unistd.h>        // close() function
+#include <netdb.h>         // getaddrinfo() function
+#include <sys/socket.h>    // socket(), connect() functions
+#include "GVOClient.hpp"
+#include "GException.hpp"
+#include "GTools.hpp"
+
+/* __ Method name definitions ____________________________________________ */
+
+/* __ Macros _____________________________________________________________ */
+
+/* __ Coding definitions _________________________________________________ */
+
+/* __ Debug definitions __________________________________________________ */
+
+
+/*==========================================================================
+ =                                                                         =
+ =                         Constructors/destructors                        =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Void constructor
+ ***************************************************************************/
+GVOClient::GVOClient(void)
+{
+    // Initialise members
+    init_members();
+
+    // Find Hub
+    find_hub();
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy constructor
+ *
+ * @param[in] client VO client.
+ ***************************************************************************/
+GVOClient::GVOClient(const GVOClient& client)
+{
+    // Initialise members
+    init_members();
+
+    // Copy members
+    copy_members(client);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Destructor
+ ***************************************************************************/
+GVOClient::~GVOClient(void)
+{
+    // Free members
+    free_members();
+
+    // Return
+    return;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                                Operators                                =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Assignment operator
+ *
+ * @param[in] client VO client.
+ * @return VO client.
+ ***************************************************************************/
+GVOClient& GVOClient::operator= (const GVOClient& client)
+{
+    // Execute only if object is not identical
+    if (this != &client) {
+
+        // Free members
+        free_members();
+
+        // Initialise members
+        init_members();
+
+        // Copy members
+        copy_members(client);
+
+    } // endif: object was not identical
+
+    // Return
+    return *this;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                             Public methods                              =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Clear object.
+ *
+ * Reset object to a clean initial state.
+ ***************************************************************************/
+void GVOClient::clear(void)
+{
+    // Free memory and initialise members
+    free_members();
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Clone object
+ ***************************************************************************/
+GVOClient* GVOClient::clone(void) const
+{
+    // Clone client
+    return new GVOClient(*this);
+}
+
+
+/***********************************************************************//**
+ * @brief Register client to SAMP Hub
+ *
+ * Connects the VO client to the Hub. A socket to the Hub is opened and the
+ * client is registered at the Hub. Metadata are sent.
+ ***************************************************************************/
+void GVOClient::connect(void)
+{
+    // Register to Hub
+    printf("Sending registration to hub\n");
+    fflush(stdout);
+    register_to_hub();
+    printf("Registration to hub sent\n");
+    fflush(stdout);
+    // Send meta data
+    send_metadata();
+    printf("Metadata sent\n");
+    fflush(stdout);
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Unregister client from SAMP Hub
+ *
+ * Disconnects the VO client from the Hub. The client is unregistered and the
+ * socket to the Hub is closed.
+ ***************************************************************************/
+void GVOClient::disconnect(void)
+{
+    // Unregister from Hub
+    unregister_from_hub();
+
+    // Close socket
+    if (m_socket != -1) {
+        close(m_socket);
+        m_socket = -1;
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Signals if client has Hub information
+ *
+ * @return True if client has Hub information, false otherwise.
+ *
+ * Checks if all mandatory Hub tokens are non-empty.
+ ***************************************************************************/
+bool GVOClient::has_hub(void) const
+{
+    // Return Hub information status
+    return (!m_secret.empty() && !m_hub_url.empty() && !m_version.empty());
+}
+
+
+/***********************************************************************//**
+ * @brief Signals if client is connected to Hub
+ *
+ * @return True if client is connected to Hub, false otherwise.
+ *
+ * Checks if a client is connected to Hub.
+ ***************************************************************************/
+bool GVOClient::is_connected(void) const
+{
+    // Return
+    return false;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns Hub response a XML object
+ *
+ * @return Hub response.
+ ***************************************************************************/
+GXml GVOClient::response(void) const
+{
+    // Declare empty XML document
+    GXml xml;
+
+    // Continue only if connection has been established
+    if (m_socket != -1) {
+	printf("Socket identified, receive string\n");
+	fflush(stdout);
+        // Receive a string
+        std::string response = receive_string();
+
+        // Find start of XML text
+        size_t start = response.find("<?xml");
+
+        // If found then convert text into XML document
+        if (start != std::string::npos) {
+            xml = GXml(response.substr(start, std::string::npos));
+        }
+
+    } // endif: connection has been established
+
+    // Return XML document
+    printf("End of response function\n");
+    fflush(stdout);
+    return xml;
+}
+
+
+/***********************************************************************//**
+ * @brief Print VO client information
+ *
+ * @param[in] chatter Chattiness (defaults to NORMAL).
+ * @return String containing VO client information
+ ***************************************************************************/
+std::string GVOClient::print(const GChatter& chatter) const
+{
+    // Initialise result string
+    std::string result;
+
+    // Continue only if chatter is not silent
+    if (chatter != SILENT) {
+
+        // Append header
+        result.append("=== GVOClient ===");
+
+        // Append client information
+        result.append("\n"+gammalib::parformat("Name")+m_name);
+
+        // Append Hub information
+        if (has_hub()) {
+            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 host (port)"));
+            result.append(m_hub_host+" ("+m_hub_port+")");
+            result.append("\n"+gammalib::parformat("SAMP protocol version"));
+            result.append(m_version);
+            result.append("\n"+gammalib::parformat("Hub connection"));
+            if (m_socket == -1) {
+                result.append("no");
+            }
+            else {
+                if (!m_client_key.empty()) {
+                    result.append("registered as \""+m_client_id);
+                    result.append("\" on Hub \""+m_hub_id);
+                    result.append("\" via socket "+gammalib::str(m_socket));
+                }
+                else {
+                    result.append("established on socket ");
+                    result.append(gammalib::str(m_socket));
+                }
+            }
+        }
+        else {
+            result.append("\n"+gammalib::parformat("Hub key"));
+            result.append("No Hub has been found");
+        }
+
+    } // endif: chatter was not silent
+
+    // Return result
+    return result;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                             Private methods                             =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Initialise class members
+ ***************************************************************************/
+void GVOClient::init_members(void)
+{
+    // Initialise members
+    m_name = "GammaLib";
+    m_secret.clear();
+    m_hub_url.clear();
+    m_hub_host.clear();
+    m_hub_port.clear();
+    m_version.clear();
+    m_client_key.clear();
+    m_hub_id.clear();
+    m_client_id.clear();
+    m_socket = -1;         // Signals no socket
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy class members
+ *
+ * @param[in] client VO client.
+ ***************************************************************************/
+void GVOClient::copy_members(const GVOClient& client)
+{
+    // Copy members
+    m_name       = client.m_name;
+    m_secret     = client.m_secret;
+    m_hub_url    = client.m_hub_url;
+    m_hub_host   = client.m_hub_host;
+    m_hub_port   = client.m_hub_port;
+    m_version    = client.m_version;
+    m_client_key = client.m_client_key;
+    m_hub_id     = client.m_hub_id;
+    m_client_id  = client.m_client_id;
+    m_socket     = client.m_socket;
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Delete class members
+ ***************************************************************************/
+void GVOClient::free_members(void)
+{
+    // Unregister from Hub
+    unregister_from_hub();
+
+    // Close socket
+    if (m_socket != -1) {
+        close(m_socket);
+        m_socket = -1;
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Find SAMP Hub
+ *
+ * @return True of SAMP Hub has been found, false otherwise.
+ *
+ * Search a valid SAMP Hub and retrieve all mandatory token for this Hub.
+ * The manadtory tokens are
+ *
+ *     samp.secret           Opaque text string required for Hub registration
+ *     samp.hub.xmlrpc.url   XML-RPC endpoint for communication
+ *     samp.profile.version  Version of SAMP profile
+ *
+ * Implements IVOA standard REC-SAMP-1.3-20120411.
+ *
+ * @todo If no valid Hub was found the method should attempt to start its
+ * own Hub.
+ ***************************************************************************/
+bool GVOClient::find_hub(void)
+{
+    // Allocate line buffer
+    const int n = 1000; 
+    char      line[n];
+
+    // Initialise find flag to false
+    bool found = false;
+
+    // Get lockfile URL
+    std::string lockurl = get_hub_lockfile();
+
+    // Continue only if a URL has been found
+    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
+        // of URL
+        if (lockurl.compare(0, 7, "file://") == 0) {
+            lockurl = lockurl.substr(7, std::string::npos);
+        }
+    
+        // Open SAMP lockfile. Continue only if opening was successful
+        FILE* fptr = fopen(lockurl.c_str(), "r");
+        if (fptr != NULL) {
+
+            // Parse lockfile and search for mandatory tokens
+            while (fgets(line, n, fptr) != NULL) {
+
+                // Convert line to C++ string
+                std::string cline = std::string(line);
+
+                // Check for secret key
+                if (cline.compare(0, 12, "samp.secret=") == 0) {
+                    m_secret = gammalib::strip_chars(cline.substr(12, std::string::npos), "\r\n");
+                }
+
+                // Check for Hub URL
+                else if (cline.compare(0, 20, "samp.hub.xmlrpc.url=") == 0) {
+                    m_hub_url = gammalib::strip_chars(cline.substr(20, std::string::npos), "\r\n");
+                }
+
+                // Check for profile version
+                else if (cline.compare(0, 21, "samp.profile.version=") == 0) {
+                    m_version = gammalib::strip_chars(cline.substr(21, std::string::npos), "\r\n");
+                }
+
+            }
+
+            // Close SAMP lockfile
+            fclose(fptr);
+
+            // Extract host and port from Hub URL
+            if (m_hub_url.compare(0, 7, "http://") == 0) {
+                size_t length;
+                size_t start = 7;
+                size_t stop  = m_hub_url.find(":", start);
+                if (stop != std::string::npos) {
+                    length = stop - start;
+                }
+                else {
+                    length = std::string::npos;
+                }
+                m_hub_host = m_hub_url.substr(start, length);
+                if (stop != std::string::npos) {
+                    stop = stop + 1;
+                    size_t end = m_hub_url.find("/", stop);
+                    if (end != std::string::npos) {
+                        length = end - stop;
+                    }
+                    else {
+                        length = std::string::npos;
+                    }
+                    m_hub_port = m_hub_url.substr(stop, length);
+                }
+            }
+
+            // Check for existence of mandatory tokens
+            found = has_hub();
+
+        } // endif: SAMP lockfile opened
+
+    } // endif: URL has been found
+
+    //TODO: Create GammaLib own Hub in case that no Hub has been found
+    //if (!found) {
+    //}
+
+    // If no Hub has been found, clear all Hub related members
+    if (!found) {
+        m_secret.clear();
+        m_hub_url.clear();
+        m_version.clear();
+    }
+
+    // Return find flag
+    return found;
+}
+
+
+/***********************************************************************//**
+ * @brief Connect to SAMP Hub
+ *
+ * Connects to Hub by creating a socket and connecting to this socket. The
+ * method expects that a Hub has already been found. If no Hub has been found
+ * (i.e. has_hub() returns false), the method does nothing.
+ *
+ * The method extracts the hostname and the port from the Hub endpoint URL.
+ * If Hub connection fails, any created socket is closed.
+ *
+ * If the method is successful, m_socket will contain on exit a non-negative
+ * number. If any failure occurs, m_socket will be set to -1.
+ ***************************************************************************/
+void GVOClient::connect_to_hub(void)
+{
+    // Close any existing socket
+    if (m_socket != -1) {
+        close(m_socket);
+        m_socket = -1;
+    }
+
+    // Continue only if host and port information was found
+    if (!m_hub_host.empty() && !m_hub_port.empty()) {
+
+        // Set hints
+        struct addrinfo hints;
+        std::memset(&hints, 0, sizeof(hints));
+        hints.ai_family   = AF_INET;
+        hints.ai_socktype = SOCK_STREAM;
+
+        // Get server information
+        struct addrinfo* servinfo;
+        if (getaddrinfo(m_hub_host.c_str(), m_hub_port.c_str(),
+                        &hints, &servinfo) == 0) {
+
+            // Loop through all the results and connect to the first
+            // we can
+            for (struct addrinfo* ptr = servinfo; ptr != NULL; ptr = ptr->ai_next) {
+
+                // Create socket
+                m_socket = socket(ptr->ai_family,
+                                  ptr->ai_socktype,
+                                  ptr->ai_protocol);
+
+                // Connect to socket if socket is valid
+                if (m_socket != -1) {
+                    if (::connect(m_socket,
+                                  ptr->ai_addr,
+                                  ptr->ai_addrlen) == -1) {
+                        close(m_socket);
+                        m_socket = -1;
+                    }
+                }
+
+            } // endfor: looped through all results
+
+        } // endif: server information was valid
+
+    } // endif: host and port information found
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Register client at SAMP Hub
+ ***************************************************************************/
+void GVOClient::register_to_hub(void)
+{
+    // Connect to Hub
+    connect_to_hub();
+
+    // Continue only if Hub connection has been established
+    if (m_socket != -1) {
+
+        // Declare message
+        std::string msg = "";
+
+        // Set metadata header
+        msg.append("<?xml version=\"1.0\"?>\n");
+        msg.append("<methodCall>\n");
+        msg.append("<methodName>samp.hub.register</methodName>\n");
+        msg.append("<params>\n");
+        msg.append("<param><value><string>"+m_secret+"</string></value></param>\n");
+        msg.append("</params>\n");
+        msg.append("</methodCall>\n");
+
+        // Post message
+        post_string(msg);
+	printf("Waiting for Hub response\n");
+	fflush(stdout);
+        // Get Hub response
+        GXml xml = response();
+
+	printf("Hub response received");
+	fflush(stdout);
+	
+        // Extract Hub and client identifiers
+        m_client_key = get_response_value(xml, "samp.private-key");
+        m_hub_id     = get_response_value(xml, "samp.hub-id");
+        m_client_id  = get_response_value(xml, "samp.self-id");
+
+    } // endif: Hub connection has been established
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Unregister client from SAMP Hub
+ ***************************************************************************/
+void GVOClient::unregister_from_hub(void)
+{
+    // Connect to Hub
+    connect_to_hub();
+
+    // Continue only if Hub connection has been established
+    if (m_socket != -1) {
+
+        // Declare message
+        std::string msg = "";
+
+        // Set metadata header
+        msg.append("<?xml version=\"1.0\"?>\n");
+        msg.append("<methodCall>\n");
+        msg.append("<methodName>samp.hub.unregister</methodName>\n");
+        msg.append("<params>\n");
+        msg.append("<param><value><string>"+m_client_key+"</string></value></param>\n");
+        msg.append("</params>\n");
+        msg.append("</methodCall>\n");
+
+        // Post message
+        post_string(msg);
+
+        // Get Hub response
+        GXml xml = response();
+
+        // Reset Hub and client identifiers
+        m_client_key.clear();
+        m_hub_id.clear();
+        m_client_id.clear();
+
+    } // endif: Hub connection has been established
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Send client metadata to Hub
+ ***************************************************************************/
+void GVOClient::send_metadata(void)
+{
+    // Connect to Hub. This was needed to post the metadata to the Hub.
+    // Apparently, the connection is lost after each write/read cycle.
+    // Maybe this is socket standard?
+    connect_to_hub();
+
+    printf("send metadata function \n");
+    fflush(stdout);
+    // Continue only if Hub connection has been established
+    if (m_socket != -1) {
+
+        // Declare message
+        std::string msg = "";
+
+        // Set metadata header
+        msg.append("<?xml version=\"1.0\"?>\n");
+        msg.append("<methodCall>\n");
+        msg.append("<methodName>samp.hub.declareMetadata</methodName>\n");
+        msg.append("<params>\n");
+        msg.append("<param><value><string>"+m_client_key+"</string></value></param>\n");
+        msg.append("<param><value><struct>\n");
+
+        // Set SAMP name
+        msg.append("<member>\n");
+        msg.append("<name>samp.name</name>\n");
+        msg.append("<value><string>"+m_name+"</string></value>\n");
+        msg.append("</member>\n");
+
+        // Set SAMP description text
+        msg.append("<member>\n");
+        msg.append("<name>samp.description.text</name>\n");
+        msg.append("<value><string>GammaLib client</string></value>\n");
+        msg.append("</member>\n");
+
+        // Set SAMP icon URL
+        msg.append("<member>\n");
+        msg.append("<name>samp.icon.url</name>\n");
+        msg.append("<value><string>http://a.fsdn.com/allura/p/gammalib/icon</string></value>\n");
+        msg.append("</member>\n");
+
+        // Set author affiliation
+        msg.append("<member>\n");
+        msg.append("<name>author.affiliation</name>\n");
+        msg.append("<value><string>IRAP, Toulouse, France</string></value>\n");
+        msg.append("</member>\n");
+
+        // Set author e-mail
+        msg.append("<member>\n");
+        msg.append("<name>author.email</name>\n");
+        msg.append("<value><string>jurgen.knodlseder@irap.omp.eu</string></value>\n");
+        msg.append("</member>\n");
+
+        // Set author name
+        msg.append("<member>\n");
+        msg.append("<name>author.name</name>\n");
+        msg.append("<value><string>Juergen Knoedlseder</string></value>\n");
+        msg.append("</member>\n");
+
+        // Set metadata trailer
+        msg.append("</struct></value></param>\n");
+        msg.append("</params>\n");
+        msg.append("</methodCall>\n");
+
+        // Post message
+        post_string(msg);
+
+        // Get Hub response
+        GXml xml = response();
+
+    } // endif: Hub connection has been established
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns value for a SAMP Hub response parameter
+ *
+ * @param[in] xml Hub response XML document.
+ * @param[in] name Parameter name.
+ * @return Parameter value.
+ *
+ * Returns value for a SAMP Hub response parameter. If the specified
+ * parameter was not found or if the response structure is not compliant,
+ * an empty string is returned.
+ ***************************************************************************/
+std::string GVOClient::get_response_value(const GXml&        xml,
+                                          const std::string& name) const
+{
+    // Declare value
+    std::string value = "";
+
+    // Search for value of specified member
+    const GXmlNode* node = xml.element("methodResponse", 0);
+    if (node != NULL) {
+        node = node->element("params", 0);
+        if (node != NULL) {
+            node = node->element("param", 0);
+            if (node != NULL) {
+                node = node->element("value", 0);
+                if (node != NULL) {
+                    node = node->element("struct", 0);
+                    if (node != NULL) {
+                        int num = node->elements("member");
+                        for (int i = 0; i < num; ++i) {
+                            const GXmlNode* member = node->element("member", i);
+                            std::string one_name;
+                            std::string one_value;
+                            get_name_value_pair(member, one_name, one_value);
+                            if (one_name == name) {
+                                value = one_value;
+                                break;
+                            }
+                        }
+                    }
+                }
+            }
+        }
+    }
+
+    // Return value
+    return value;
+}
+
+
+/***********************************************************************//**
+ * @brief Extract name / value pair from XML node
+ *
+ * @param[in] node Pointer to XML node.
+ * @param[out] name Name string.
+ * @param[out] value Value string.
+ *
+ * Extracts a name / value pair from a XML node. If the XML node pointer is
+ * NULL, the name and value strings will be empty.
+ ***************************************************************************/
+void GVOClient::get_name_value_pair(const GXmlNode* node,
+                                    std::string&    name,
+                                    std::string&    value) const
+{
+    // Clear name and value strings
+    name.clear();
+    value.clear();
+
+    // Continue only if node is valid
+    if (node != NULL) {
+
+        // Get name node and extract text content
+        const GXmlNode* ptr = node->element("name", 0);
+        if (ptr != NULL) {
+            const GXmlText* text = static_cast<const GXmlText*>((*ptr)[0]);
+            if (text != NULL) {
+                name = text->text();
+            }
+        }
+
+        // Get value node and extract text content
+        ptr = node->element("value", 0);
+        if (ptr != NULL) {
+            const GXmlText* text = static_cast<const GXmlText*>((*ptr)[0]);
+            if (text != NULL) {
+                value = text->text();
+            }
+        }
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Post string content to Hub
+ *
+ * @param[in] content String content to post
+ *
+ * Posts the content of a string to the Hub.
+ *
+ * The method does nothing if no Hub connection has been established.
+ ***************************************************************************/
+void GVOClient::post_string(const std::string& content) const
+{
+    // Continue only if Hub connection has been established
+    if (m_socket != -1) {
+
+        // Determine content length
+        int length = content.length();
+
+        // Set prefix
+        std::string prefix = "POST /xmlrpc-1 HTTP/1.0\n"
+                             "User-Agent: GammaLib\n"
+                             "Content-Type: text/xml\n"
+                             "Content-Length: "+gammalib::str(length)+"\n\n";
+
+        // Build post string
+        std::string post = prefix + content;
+
+        // Send content to socket
+        bool done = false;
+        do {
+            int length      = post.length();
+            int sent_length = send(m_socket, post.c_str(), length, 0);
+            if (sent_length < length) {
+                post = post.substr(sent_length, std::string::npos);
+            }
+            else {
+                done = true;
+            }
+        } while (!done);
+
+    } // endif: Hub connection had been established
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Receive string content from Hub
+ *
+ * @return String received from Hub.
+ *
+ * Reads information sent by Hub into a string.
+ * 
+ * The method does nothing if no Hub connection has been established.
+ ***************************************************************************/
+std::string GVOClient::receive_string(void) const
+{
+    // Initialise empty string
+    std::string result = "";
+
+    // Continue only if Hub connection has been established
+    if (m_socket != -1) {
+
+        // Define buffer
+        char buffer[1001];
+
+        // Read buffer until it is empty
+        int n = 0;
+        //do {
+	    printf("receive socket instruction\n");
+	    fflush(stdout);
+            n = recv(m_socket, buffer, 1000, 0);
+	    printf("data received: %s\n",buffer);
+	    fflush(stdout);
+            if (n > 0) {
+                buffer[n+1] = '\0';
+                result.append(std::string(buffer));
+            }
+        //} while (n > 0);
+
+    } // endif: Hub connection had been established
+
+    // Return result
+    return result;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns SAMP Hub lockfile URL
+ *
+ * @return SAMP Hub lockfile URL (empty if no lockfile was found).
+ *
+ * Implements IVOA standard REC-SAMP-1.3-20120411.
+ ***************************************************************************/
+std::string GVOClient::get_hub_lockfile(void) const
+{
+    // Initialise result
+    std::string url = "";
+
+    // Check for existence of the SAMP_HUB environment variable first
+    char* hub_ptr = std::getenv("SAMP_HUB");
+    if (hub_ptr != NULL) {
+
+        // Check for mandatory std-lockurl: prefix (no other prefixe is
+        // supported so far)
+        std::string lockurl = std::string(hub_ptr);
+        if (lockurl.compare(0, 12, "std-lockurl:") == 0) {
+
+            // Extract URL
+            url = lockurl.substr(12, std::string::npos);
+
+        } // endif: std-lockurl: prefix found
+  
+    }
+
+    // ... otherwise the lockfile should be $HOME/.samp
+    else {
+
+        // Get user's HOME directory path as the prefix of the full
+        // path. If the HOME environment variable is not set we
+        // expect that .samp is in the local directory. This is non
+        // standard, but prevents for creating an exception here.
+        std::string prefix = "";
+        char* home_ptr = std::getenv("HOME");
+        if (home_ptr != NULL) {
+            prefix = std::string(home_ptr) + "/";
+        }
+
+        // Set filename
+        url = prefix + ".samp";
+
+    } // endelse: no SAMP_HUB environment variable found
+
+    // Return URL
+    return url;
+}
@@ -488,7 +488,7 @@ void GVOHub::handle_request(const socklen_t&amp; sock)
         request_ping(sock);
     }
     else if (method_called.compare("samp.hub.call") == 0) {
-        request_ping(sock);
+        request_call(xml,sock);
     }
     else if (method_called.compare("samp.hub.callAll") == 0) {
         request_ping(sock);
@@ -904,9 +904,8 @@ void GVOHub::request_get_registered_clients(const GXml&amp; xml, const socklen_t&amp; so
  
     // Declare message
     std::string msg = "";
-
     // Set response
-    msg.append("<?xml version=\'1.0\' encoding=\"UTF-8\"?>\n");
+    msg.append("<?xml version=\"1.0\" encoding=\"UTF-8\"?>\n");
     msg.append("<methodResponse>\n");
     msg.append("<params>\n");
     msg.append("  <param><value><array><data>\n");
@@ -925,7 +924,6 @@ void GVOHub::request_get_registered_clients(const GXml&amp; xml, const socklen_t&amp; so
     msg.append("  </data></array></value></param>\n");
     msg.append("</params>\n");
     msg.append("</methodResponse>\n");
-
     // Post response
     post_string(msg, sock);
  
@@ -950,6 +948,14 @@ void GVOHub::request_get_subscribed_clients(const GXml&amp; xml, const socklen_t&amp; so
     std::cout << "GVOHub::request_get_subscribed_clients" << std::endl;
     #endif
  
+    // Declare message
+    std::string msg = "";    
+
+    msg.append("<?xml version=\'1.0\' encoding=\"UTF-8\"?>\n");
+    msg.append("<methodResponse>\n");
+    msg.append("<params>\n");
+    msg.append("  <param><value><struct><member><name>c0</name><value><struct></struct></value></member><member><name>gammalib_hub</name><value><struct></struct></value></member></struct></value></param></params></methodResponse>");
+
 /*
     std::string buffer;
     const int n = 1000; 
@@ -971,8 +977,10 @@ void GVOHub::request_get_subscribed_clients(const GXml&amp; xml, const socklen_t&amp; so
     fclose(fptr);
 */
  
+    // Post response
+    post_string(msg, sock);
     // Post void message
-    post_samp_void(sock);
+    //post_samp_void(sock);
  
     // Return
     return;
@@ -1114,6 +1122,77 @@ void GVOHub::request_get_metadata(const GXml&amp; xml, const socklen_t&amp; sock)
     return;
 }
  
+/***********************************************************************//**
+ * @brief Handles calls to specific client
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming cals for a specific client
+ ***************************************************************************/
+void GVOHub::request_call(const GXml& xml, const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_get_registered_clients" << std::endl;
+    #endif
+
+    // Get client key
+    std::string key = get_client_key(xml);
+    // Get message mtype
+    std::string mtype = get_mtype(xml);
+    // Get the key of client who is to be contacted with the message
+    // Declare message
+    std::string msg = "";
+    // Determine which call has been invoqued
+    // Loop over all clients and send notification
+    for (int i = 0; i < m_clients.size(); ++i) {
+	    msg = "";
+	    
+	    // Set response for table.load.votable
+	    msg.append("<?xml version=\"1.0\"?>\n");
+	    msg.append("<methodCall><methodName>samp.client.receiveCall</methodName>\n");
+	    msg.append("<params>\n");
+	    msg.append("  <param><value>");
+	    msg.append(m_clients[i].private_key);
+	    msg.append("</value></param>\n");
+	    msg.append("  <param><value>"+m_hub_id+"</value></param>\n");
+	    msg.append("  <param><value><struct>\n");
+	    msg.append("    <member>\n");
+	    msg.append("      <name>samp.mtype</name>\n");
+	    msg.append("      <value>table.load.votable</value>\n");
+	    msg.append("    </member>\n");
+	    msg.append("    <member>\n");
+	    msg.append("      <name>samp.params</name>\n");
+	    msg.append("      <value><struct>\n");
+	    msg.append("        <member>\n");
+	    msg.append("          <name>name</name>\n");
+	    msg.append("          <value>imagename</value>\n");
+	    msg.append("        </member>\n");
+	    msg.append("        <member>\n");
+	    msg.append("          <name>table-id</name>\n");
+	    msg.append("          <value>id of table</value>\n");
+	    msg.append("        </member>\n");
+	    msg.append("        <member>\n");
+	    msg.append("           <name>url</name>\n");
+	    msg.append("           <value>url to the table</value>\n");
+	    msg.append("         </member>\n");
+	    msg.append("       </struct></value>\n");
+	    msg.append("      </member>\n");				
+	    msg.append("      </struct></value>\n");
+	    msg.append("    </param>\n");
+	    msg.append("</params>\n");
+	    msg.append("</methodCall>\n");
+              
+	
+            // Post notification
+            notify(m_clients[i].url, msg);
+	}
+    
+    // Return
+    return;
+  
+}
  
 /***********************************************************************//**
  * @brief Handles Hub shutdown requests
@@ -1947,4 +2026,61 @@ std::string GVOHub::random_string(const size_t&amp; length) const
     // Return string
     return str;
 }
+/***********************************************************************//**
+ * @brief Extract mtype from XML request
+ *
+ * @param[in] xml XML message sent by client.
+ *
+ * Extracts the mtype identifying calling message from the XML request.
+ ***************************************************************************/
+std::string GVOHub::get_mtype(const GXml& xml) const
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::get_mtype" << std::endl;
+    #endif
+
+    // Initialise response
+    std::string client_key = "";
+
+    // Get the client's private key
+    const GXmlNode* node = xml.element("methodCall > params > param > value");
+    if (node != NULL) {
+        const GXmlText* text = static_cast<const GXmlText*>((*node)[0]);
+        if (text != NULL) {
+            client_key = text->text();
+        }
+    }
+
+    // Return key
+    return client_key;
+}
+/***********************************************************************//**
+ * @brief Extract client destination from XML request
+ *
+ * @param[in] xml XML message sent by client.
+ *
+ * Extracts identifier of destination client from the XML request.
+ ***************************************************************************/
+std::string GVOHub::get_destination(const GXml& xml) const
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::get_mtype" << std::endl;
+    #endif
  
+    // Initialise response
+    std::string client_key = "";
+
+    // Get the client's private key
+    const GXmlNode* node = xml.element("methodCall > params > param > value");
+    if (node != NULL) {
+        const GXmlText* text = static_cast<const GXmlText*>((*node)[0]);
+        if (text != NULL) {
+            client_key = text->text();
+        }
+    }
+
+    // Return key
+    return client_key;
+}
+/***************************************************************************
+ *                     GVOHub.cpp - VO SAMP Hub class                      *
+ * ----------------------------------------------------------------------- *
+ *  copyright (C) 2014-2015 by Thierry Louge                               *
+ * ----------------------------------------------------------------------- *
+ *                                                                         *
+ *  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   *
+ *  the Free Software Foundation, either version 3 of the License, or      *
+ *  (at your option) any later version.                                    *
+ *                                                                         *
+ *  This program is distributed in the hope that it will be useful,        *
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of         *
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          *
+ *  GNU General Public License for more details.                           *
+ *                                                                         *
+ *  You should have received a copy of the GNU General Public License      *
+ *  along with this program.  If not, see <http://www.gnu.org/licenses/>.  *
+ *                                                                         *
+ ***************************************************************************/
+/**
+ * @file GVOHub.cpp
+ * @brief VO SAMP Hub class implementation
+ * @author Thierry Louge
+ */
+
+/* __ Includes ___________________________________________________________ */
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+#include <cstdlib>         // std::getenv() function
+#include <cstdio>          // std::fopen(), etc. functions
+#include <cstring>         // std::memset() function
+#include <csignal>         // signal() function
+#include <cerrno>          // errno
+#include <unistd.h>        // close() function
+#include <netdb.h>         // getaddrinfo() function
+#include <netinet/in.h>    // sockaddr_in, INADDR_ANY, htons
+#include <fstream>
+#include <sys/shm.h>
+#include <sys/socket.h>    // socket(), connect() functions
+#include <sys/wait.h>      // waitpid() function
+#include <arpa/inet.h>     // inet_addr() function
+#include "GVOHub.hpp"
+#include "GException.hpp"
+#include "GTools.hpp"
+
+/* __ Method name definitions ____________________________________________ */
+#define G_START_HUB                                     "GVOHub::start_hub()"
+
+/* __ Macros _____________________________________________________________ */
+
+/* __ Coding definitions _________________________________________________ */
+
+/* __ Debug definitions __________________________________________________ */
+//#define G_CONSOLE_DUMP
+#define G_SHOW_MESSAGE
+
+
+/*==========================================================================
+ =                                                                         =
+ =                         Constructors/destructors                        =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Void constructor
+ ***************************************************************************/
+GVOHub::GVOHub(void)
+{
+    // Initialise members
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy constructor
+ *
+ * @param[in] hub VO Hub.
+ ***************************************************************************/
+GVOHub::GVOHub(const GVOHub& hub)
+{
+    // Initialise members
+    init_members();
+
+    // Copy members
+    copy_members(hub);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Destructor
+ ***************************************************************************/
+GVOHub::~GVOHub(void)
+{
+    // Free members
+    free_members();
+
+    // Return
+    return;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                                Operators                                =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Assignment operator
+ *
+ * @param[in] hub VO hub.
+ * @return VO hub.
+ ***************************************************************************/
+GVOHub& GVOHub::operator=(const GVOHub& hub)
+{
+    // Execute only if object is not identical
+    if (this != &hub) {
+
+        // Free members
+        free_members();
+
+        // Initialise members
+        init_members();
+
+        // Copy members
+        copy_members(hub);
+
+    } // endif: object was not identical
+
+    // Return
+    return *this;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                             Public methods                              =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Clear object.
+ *
+ * Reset object to a clean initial state.
+ ***************************************************************************/
+void GVOHub::clear(void)
+{
+    // Free members
+    free_members();
+
+    // Initialise members
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Clone object
+ ***************************************************************************/
+GVOHub* GVOHub::clone(void) const
+{
+    // Clone client
+    return new GVOHub(*this);
+}
+
+
+/***********************************************************************//**
+ * @brief Start Hub
+ ***************************************************************************/
+void GVOHub::start(void)
+{
+    // Create SAMP file
+    create_samp_file();
+
+    // Start Hub
+    start_hub();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Print VO hub information
+ *
+ * @param[in] chatter Chattiness (defaults to NORMAL).
+ * @return String containing VO hub information
+ ***************************************************************************/
+std::string GVOHub::print(const GChatter& chatter) const
+{
+    // Initialise result string
+    std::string result;
+
+    // Continue only if chatter is not silent
+    if (chatter != SILENT) {
+
+        // Append header
+        result.append("=== GVOHub ===");
+
+        // 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 host")+m_hub_host);
+        result.append("\n"+gammalib::parformat("Hub port")+m_hub_port);
+        result.append("\n"+gammalib::parformat("Hub path")+m_hub_path);
+        result.append("\n"+gammalib::parformat("SAMP protocol version"));
+        result.append(m_version);
+
+    } // endif: chatter was not silent
+
+    // Return result
+    return result;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                             Private methods                             =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Initialise class members
+ ***************************************************************************/
+void GVOHub::init_members(void)
+{
+    // Initialise members
+    m_secret        = random_string(15);
+    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
+    m_shutdown      = false;
+    m_clients.clear();
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy class members
+ *
+ * @param[in] hub VO hub.
+ ***************************************************************************/
+void GVOHub::copy_members(const GVOHub& 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;
+    m_version  = hub.m_version;
+    m_hub_id   = hub.m_hub_id;
+    m_socket   = hub.m_socket;
+    m_shutdown = hub.m_shutdown;
+    m_clients  = hub.m_clients;
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Delete class members
+ ***************************************************************************/
+void GVOHub::free_members(void)
+{ 
+    // Close sockets
+    if (m_socket != -1) {
+        close(m_socket);
+        m_socket = -1;
+    }
+    
+    // Remove lockfile
+    std::string lockurl = get_hub_lockfile();
+    std::remove(lockurl.c_str());
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Starts the SAMP hub socket and listens on it
+ *
+ * @exception GException::runtime_error
+ *            Problem with creating of or listening on socket
+ *
+ * This is the main Hub event loop.
+ ***************************************************************************/
+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);
+    }
+
+    // Now 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="+
+                          gammalib::str(errno);
+        throw GException::runtime_error(G_START_HUB, msg);
+    }
+
+    // ...
+    struct sockaddr_in cli_addr;
+    socklen_t clilen = sizeof(cli_addr);
+
+    // Main event handling loop
+    while (1) {
+
+        // Accept connection from the client 
+    	socklen_t socket = accept(m_socket, (struct sockaddr *)&cli_addr, &clilen);
+    	if (socket < 0) {
+            std::string msg = "Client connection to socket not accepted.";
+            throw GException::runtime_error(G_START_HUB, msg);
+    	}
+
+        // Handle request
+        handle_request(socket);
+
+        // Close socket
+        close(socket);
+
+        // If we have received a shutdown request then exit the event
+        // loop
+        if (m_shutdown) {
+            break;
+        }
+
+    } // endwhile: main event loop
+
+    // Close socket
+    close(m_socket);
+    m_socket = -1;
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Reads the client message and runs appropriate function
+ *
+ * @param[in] sock Socket to client.
+ ***************************************************************************/
+void GVOHub::handle_request(const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::handle_request" << std::endl;
+    #endif
+
+    // Initialize buffer
+    char buffer[1001];
+
+    // Initialise message
+    std::string message = "";
+
+    // Read from socket until nothing is received anymore.
+    int timeout = 2000; // Initial timeout is 2 sec
+    int n       = 0;
+    do {
+        n = gammalib::recv(sock, buffer, 1000, 0, timeout);
+        if (n > 0) {
+            buffer[n] = '\0';
+            message.append(std::string(buffer));
+        }
+        timeout = 10; // The timeout now is 0.01 sec 
+    } while (n > 0);
+    
+    // Dump the buffer
+    #if defined(G_SHOW_MESSAGE)
+    std::cout << std::endl;
+    std::cout << "GVOHub has received the following message:" << std::endl;
+    std::cout << "==========================================" << std::endl;
+    std::cout << message << std::endl;
+    #endif
+    
+    // Extract response into an XML object
+    GXml   xml;
+    size_t start = message.find("<?xml");
+    if (start != std::string::npos) {
+        xml = GXml(message.substr(start, std::string::npos));
+    }
+
+    // Get methodName value
+    std::string method_called;
+    const GXmlNode* node = xml.element("methodCall > methodName");
+    if (node != NULL) {
+        const GXmlText* text = static_cast<const GXmlText*>((*node)[0]);
+        if (text != NULL) {
+            method_called = text->text();
+        }
+    }
+
+    // Dump the buffer
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "Method called: " << method_called << std::endl;
+    #endif
+
+    // Dispatch according to method
+    if (method_called.compare("samp.hub.ping") == 0) {
+        request_ping(sock);
+    }
+    else if (method_called.compare("samp.hub.register") == 0) {
+        request_register(xml, sock);
+    }
+    else if (method_called.compare("samp.hub.unregister") == 0) {
+        request_unregister(xml, sock);
+    }
+    else if (method_called.compare("samp.hub.declareMetadata") == 0) {
+        request_declare_metadata(xml, sock);
+    }
+    else if (method_called.compare("samp.hub.declareSubscriptions") == 0) {
+        request_declare_subscriptions(xml, sock);
+    }
+    else if (method_called.compare("samp.hub.setXmlrpcCallback") == 0) {
+        request_set_xml_rpc_callback(xml, sock);
+    }
+    else if (method_called.compare("samp.hub.getSubscriptions") == 0) {
+        request_get_subscriptions(xml, sock);
+    }
+    else if (method_called.compare("samp.hub.getRegisteredClients") == 0) {
+        request_get_registered_clients(xml, sock);
+    }
+    else if (method_called.compare("samp.hub.getSubscribedClients") == 0) {
+        request_get_subscribed_clients(xml, sock);
+    }
+    else if (method_called.compare("samp.hub.getMetadata") == 0) {
+        request_get_metadata(xml, sock);
+    }
+    else if (method_called.compare("samp.hub.notify") == 0) {
+        request_ping(sock);
+    }
+    else if (method_called.compare("samp.hub.notifyAll") == 0) {
+        request_ping(sock);
+    }
+    else if (method_called.compare("samp.hub.call") == 0) {
+        request_call(xml,sock);
+    }
+    else if (method_called.compare("samp.hub.callAll") == 0) {
+        request_ping(sock);
+    }
+    else if (method_called.compare("samp.hub.callAndWait") == 0) {
+        request_ping(sock);
+    }
+    else if (method_called.compare("samp.hub.reply") == 0) {
+        request_ping(sock);
+    }
+    else if (method_called.compare("samp.hub.shutdown") == 0) {
+        request_shutdown(sock);
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Handles ping requests
+ *
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming ping requests.
+ ***************************************************************************/
+void GVOHub::request_ping(const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_ping" << std::endl;
+    #endif
+    
+    // Post void message
+    post_samp_void(sock);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Handles registration requests
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming registration requests.
+ ***************************************************************************/
+void GVOHub::request_register(const GXml& xml, const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_register" << std::endl;
+    #endif
+
+    // Set client reference
+    int         counter = 0;
+    std::string reference;
+    do {
+        reference = "c"+gammalib::str(counter);
+        int i = 0;
+        for (; i < m_clients.size(); ++i) {
+            if (m_clients[i].reference == reference) {
+                counter++;
+                break;
+            }
+        }
+        if (i >= m_clients.size()) {
+            break;
+        }
+    } while (true);
+
+    // Create a new client
+    struct client voclient;
+    voclient.reference   = reference;
+    voclient.private_key = random_string(15);
+    
+    // Attached client
+    m_clients.push_back(voclient);
+
+    // Declare response
+    std::string response = "";
+
+    // Set response
+    response.append("<?xml version='1.0' encoding='UTF-8'?>\n");
+    response.append("<methodResponse>\n");
+    response.append("<params>\n");
+    response.append("  <param><value><struct>\n");
+    response.append("    <member>\n");
+    response.append("      <name>samp.private-key</name>\n");
+    response.append("      <value>"+voclient.private_key+"</value>\n");
+    response.append("    </member>\n");
+    response.append("    <member>\n");
+    response.append("      <name>samp.hub-id</name>\n");
+    response.append("      <value>"+m_hub_id+"</value>\n");
+    response.append("    </member>\n");
+    response.append("    <member>\n");
+    response.append("      <name>samp.self-id</name>\n");
+    response.append("      <value>"+voclient.reference+"</value>\n");
+    response.append("    </member>\n");
+    response.append("    <member>\n");
+    response.append("      <name>samp.status</name>\n");
+    response.append("      <value>samp.ok</value>\n");
+    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("    </member>\n");
+    response.append("  </struct></value></param>\n");
+    response.append("</params>\n");
+    response.append("</methodResponse>\n");
+
+    // Post response
+    post_string(response, sock);
+
+    // Send notifications
+    send_notifications("samp.hub.event.register", voclient.reference);	
+
+    //Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Handles unregistration requests
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming unregistration requests.
+ ***************************************************************************/
+void GVOHub::request_unregister(const GXml& xml, const socklen_t& sock)
+{
+    // Optionally dump header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_unregister" << std::endl;
+    #endif 
+
+    // Get client index
+    int i = get_client_index(xml);
+
+    // Continue only if index is valid
+    if (i != -1) {
+
+        // Post void message
+        post_samp_void(sock);
+
+        // Send notifications
+        send_notifications("samp.hub.event.unregister",
+                           m_clients[i].reference);
+
+        // Remove client
+        m_clients.erase(m_clients.begin()+i);
+        
+    } // endif: valid client found
+
+    // Signal if the client was unknown
+    #if defined(G_CONSOLE_DUMP)
+    else {
+        std::cout << " *** ERROR: Unable to find client " << get_client_key(xml) << std::endl;
+    }
+    #endif
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Handles metadata declaration requests
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming metadata declaration requests.
+ ***************************************************************************/
+void GVOHub::request_declare_metadata(const GXml& xml, const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_declare_metadata" << std::endl;
+    #endif
+
+    // Get client index
+    int i = get_client_index(xml);
+
+    // Continue only if index is valid
+    if (i != -1) {
+
+        // Store metadata
+        m_clients[i].name          = get_response_value(xml, "samp.name");;
+        m_clients[i].description   = get_response_value(xml, "samp.description.text");
+        m_clients[i].icon          = get_response_value(xml, "samp.icon.url");
+        m_clients[i].documentation = get_response_value(xml, "samp.documentation.url");
+        m_clients[i].affiliation   = get_response_value(xml, "author.affiliation");
+        m_clients[i].email         = get_response_value(xml, "author.email");
+        m_clients[i].author_name   = get_response_value(xml, "author.name");
+        m_clients[i].homepage      = get_response_value(xml, "home.page");
+
+        // Post void message
+        post_samp_void(sock);
+            
+        // Send notifications
+        send_notifications("samp.hub.event.metadata", m_clients[i].reference);
+
+    } // endif: client found
+        
+    // Signal if the client was unknown
+    #if defined(G_CONSOLE_DUMP)
+    else {
+        std::cout << " *** ERROR: Unable to find client " << get_client_key(xml) << std::endl;
+    }
+    #endif
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Handles subscriptions declaration requests
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming subscription declaration requests.
+ ***************************************************************************/
+void GVOHub::request_declare_subscriptions(const GXml& xml,const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_declare_subscriptions" << std::endl;
+    #endif
+
+    // Get client index
+    int i = get_client_index(xml);
+
+    // Continue only if index is valid
+    if (i != -1) {
+
+        // Get subscriptions
+        std::vector<std::string> subscriptions = get_subscriptions(xml);
+
+        // Store subscriptions
+        for (int k = 0; k < subscriptions.size(); ++k) {
+            m_clients[i].subscriptions.push_back(subscriptions[k]);
+        }
+
+        // Post SAMP ok
+        post_samp_ok(sock);
+        
+    } // endif: valid client found
+
+    // Signal if the client was unknown
+    #if defined(G_CONSOLE_DUMP)
+    else {
+        std::cout << " *** ERROR: Unable to find client " << get_client_key(xml) << std::endl;
+    }
+    #endif
+    
+    //Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Handles XML-RPC callback setting requests
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming XML-RPC callback setting requests.
+ ***************************************************************************/
+void GVOHub::request_set_xml_rpc_callback(const GXml&      xml,
+                                          const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_set_xml_rpc_callback" << std::endl;
+    #endif
+
+    // Get client index
+    int i = get_client_index(xml);
+
+    // Continue only if index is valid
+    if (i != -1) {
+
+        // Get callback URL
+        std::string client_url = get_callback_url(xml);
+
+        // Set callback URL
+        m_clients[i].url = get_callback_url(xml);
+
+        // Post SAMP ok
+        post_samp_ok(sock);
+
+    } // endif: client found
+
+    // Signal if the client was unknown
+    #if defined(G_CONSOLE_DUMP)
+    else {
+        std::cout << " *** ERROR: Unable to find client " << get_client_key(xml) << std::endl;
+    }
+    #endif
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Handles subscriptions getting requests
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming subscriptions getting requests.
+ ***************************************************************************/
+void GVOHub::request_get_subscriptions(const GXml& xml, const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_get_subscriptions" << std::endl;
+    #endif
+
+    // Initialise client name
+    std::string client_name = "";
+
+    // Get the client name
+    const GXmlNode* node = xml.element("methodCall > params > param[1] > value");
+    if (node != NULL) {
+        const GXmlText* text = static_cast<const GXmlText*>((*node)[0]);
+        if (text != NULL) {
+            client_name = text->text();
+        }
+    }
+
+    // Declare response
+    std::string response = "";
+
+    // Set response
+    response.append("<?xml version='1.0' encoding='UTF-8'?>\n");
+    response.append("<methodResponse>\n");
+    response.append("<params>\n");
+    response.append("  <param><value><struct>\n");
+
+    // If key is "gammalib_hub" then return some special information
+    if (client_name == "gammalib_hub") {
+        // Do nothing
+    }
+    
+    // ... handle all other clients
+    else {
+
+        // Get client index
+        int i = get_client_index(client_name);
+
+        // Continue only if index is valid
+        if (i != -1) {
+
+            // Append all subscriptions
+            for (int k = 0; k < m_clients[i].subscriptions.size(); ++k) {
+                response.append("    <member>\n");
+                response.append("      <name>");
+                response.append(m_clients[i].subscriptions[k]);
+                response.append("</name>\n");
+                response.append("      <value></value>\n");
+                response.append("    </member>\n");
+            }
+            
+        }
+
+        // Signal if the client was unknown
+        #if defined(G_CONSOLE_DUMP)
+        else {
+            std::cout << " *** ERROR: Unable to find client " << get_client_key(xml) << std::endl;
+        }
+        #endif
+        
+    }
+
+    // Finish response
+    response.append("  </struct></value></param>\n");
+    response.append("</params>\n");
+    response.append("</methodResponse>\n");
+        
+    // Post response
+    post_string(response, sock);
+   
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Handles registered client information requests
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming registered client information requests.
+ ***************************************************************************/
+void GVOHub::request_get_registered_clients(const GXml& xml, const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_get_registered_clients" << std::endl;
+    #endif
+
+    // Get client key
+    std::string key = get_client_key(xml);
+
+    // Declare message
+    std::string msg = "";
+    // Set response
+    msg.append("<?xml version=\"1.0\" encoding=\"UTF-8\"?>\n");
+    msg.append("<methodResponse>\n");
+    msg.append("<params>\n");
+    msg.append("  <param><value><array><data>\n");
+    msg.append("    <value>gammalib_hub</value>\n");
+
+    // Loop over all clients. Do not send back current client's registration
+    for (int i = 0; i < m_clients.size(); ++i) {
+        if (m_clients[i].private_key != key) {
+            msg.append("    <value>");
+            msg.append(m_clients[i].reference);
+            msg.append("</value>\n");
+        }
+    }
+    
+    // Finish response
+    msg.append("  </data></array></value></param>\n");
+    msg.append("</params>\n");
+    msg.append("</methodResponse>\n");
+    // Post response
+    post_string(msg, sock);
+    
+    // Return
+    return;
+  
+}
+
+
+/***********************************************************************//**
+ * @brief Handles subscribed client information requests
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming subscribed client information requests.
+ ***************************************************************************/
+void GVOHub::request_get_subscribed_clients(const GXml& xml, const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_get_subscribed_clients" << std::endl;
+    #endif
+
+    // Declare message
+    std::string msg = "";    
+
+    msg.append("<?xml version=\'1.0\' encoding=\"UTF-8\"?>\n");
+    msg.append("<methodResponse>\n");
+    msg.append("<params>\n");
+    msg.append("  <param><value><struct><member><name>c0</name><value><struct></struct></value></member><member><name>gammalib_hub</name><value><struct></struct></value></member></struct></value></param></params></methodResponse>");
+
+/*
+    std::string buffer;
+    const int n = 1000; 
+    char      line[n];
+    FILE* fptr = fopen("SAMP-registered.txt", "r");
+    if (fptr != NULL) {
+        while (fgets(line, n, fptr) != NULL) {
+
+                // Convert line to C++ string
+                std::string cline = std::string(line);
+
+                // Check for secret key
+                if (cline.compare(0, 10, "Registered") == 0) {
+                    buffer = gammalib::strip_chars(cline.substr(10, std::string::npos), "\r\n");
+                }
+    }
+    
+    // Closes file
+    fclose(fptr);
+*/
+
+    // Post response
+    post_string(msg, sock);
+    // Post void message
+    //post_samp_void(sock);
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Handles a metadata requests
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all metadata requests.
+ ***************************************************************************/
+void GVOHub::request_get_metadata(const GXml& xml, const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_get_metadata" << std::endl;
+    #endif
+
+    // Initialise client name
+    std::string client_name = "";
+
+    // Get the client name
+    const GXmlNode* node = xml.element("methodCall > params > param[1] > value");
+    if (node != NULL) {
+        const GXmlText* text = static_cast<const GXmlText*>((*node)[0]);
+        if (text != NULL) {
+            client_name = text->text();
+        }
+    }
+
+    // Declare response
+    std::string response = "";
+
+    // Set response
+    response.append("<?xml version='1.0' encoding='UTF-8'?>\n");
+    response.append("<methodResponse>\n");
+    response.append("<params>\n");
+    response.append("  <param><value><struct>\n");
+    response.append("    <member>\n");
+    response.append("      <name>samp.status</name>\n");
+    response.append("      <value>samp.ok</value>\n");
+    response.append("    </member>\n");
+
+    // If key is "gammalib_hub" then return some special information
+    if (client_name == "gammalib_hub") {
+		response.append("    <member>\n");
+        response.append("      <name>samp.name</name>\n");
+        response.append("      <value>gammalib_hub</value>\n");
+        response.append("    </member>\n");
+   		response.append("    <member>\n");
+        response.append("      <name>samp.description.text</name>\n");
+        response.append("      <value>GammaLib VO Hub</value>\n");
+        response.append("    </member>\n");
+		response.append("    <member>\n");
+        response.append("      <name>samp.icon.url</name>\n");
+        response.append("      <value>http://a.fsdn.com/allura/p/gammalib/icon</value>\n");
+        response.append("    </member>\n");
+		response.append("    <member>\n");
+        response.append("      <name>samp.documentation.url</name>\n");
+        response.append("      <value>http://cta.irap.omp.eu/gammalib/user_manual/modules/vo.html</value>\n");
+        response.append("    </member>\n");
+		response.append("    <member>\n");
+        response.append("      <name>author.affiliation</name>\n");
+        response.append("      <value>IRAP, Toulouse, France</value>\n");
+        response.append("    </member>\n");
+		response.append("    <member>\n");
+        response.append("      <name>author.email</name>\n");
+        response.append("      <value>jurgen.knodlseder@irap.omp.eu</value>\n");
+        response.append("    </member>\n");
+		response.append("    <member>\n");
+        response.append("      <name>author.name</name>\n");
+        response.append("      <value>J. Knoedlseder, T. Louge</value>\n");
+        response.append("    </member>\n");
+		response.append("    <member>\n");
+        response.append("      <name>home.page</name>\n");
+        response.append("      <value>http://cta.irap.omp.eu/gammalib/</value>\n");
+        response.append("    </member>\n");
+    }
+
+    // ... otherwise get the index
+    else {
+ 
+        // Get client index
+        int i = get_client_index(client_name);
+
+        // Continue only if index is valid
+        if (i != -1) {
+
+            // Append response
+            response.append("    <member>\n");
+            response.append("      <name>samp.name</name>\n");
+            response.append("      <value>"+m_clients[i].name+"</value>\n");
+            response.append("    </member>\n");
+            response.append("    <member>\n");
+            response.append("      <name>samp.description.text</name>\n");
+            response.append("      <value>"+m_clients[i].description+"</value>\n");
+            response.append("    </member>\n");
+            response.append("    <member>\n");
+            response.append("      <name>samp.icon.url</name>\n");
+            response.append("      <value>"+m_clients[i].icon+"</value>\n");
+            response.append("    </member>\n");
+            response.append("    <member>\n");
+            response.append("      <name>samp.documentation.url</name>\n");
+            response.append("      <value>"+m_clients[i].documentation+"</value>\n");
+            response.append("    </member>\n");
+            response.append("    <member>\n");
+            response.append("      <name>author.affiliation</name>\n");
+            response.append("      <value>"+m_clients[i].affiliation+"</value>\n");
+            response.append("    </member>\n");
+            response.append("    <member>\n");
+            response.append("      <name>author.email</name>\n");
+            response.append("      <value>"+m_clients[i].email+"</value>\n");
+            response.append("    </member>\n");
+            response.append("    <member>\n");
+            response.append("      <name>author.name</name>\n");
+            response.append("      <value>"+m_clients[i].author_name+"</value>\n");
+            response.append("    </member>\n");
+            response.append("    <member>\n");
+            response.append("      <name>home.page</name>\n");
+            response.append("      <value>"+m_clients[i].homepage+"</value>\n");
+            response.append("    </member>\n");
+
+        } // endif: index was valid
+        
+    } // endelse: client was not the hub
+
+    // Finish response
+    response.append("  </struct></value></param>\n");
+    response.append("</params>\n");
+    response.append("</methodResponse>\n");
+    
+    // Post response
+    post_string(response, sock);
+   
+    // Return
+    return;
+}
+
+/***********************************************************************//**
+ * @brief Handles calls to specific client
+ *
+ * @param[in] xml XML message sent by client.
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming cals for a specific client
+ ***************************************************************************/
+void GVOHub::request_call(const GXml& xml, const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_get_registered_clients" << std::endl;
+    #endif
+
+    // Get client key
+    std::string key = get_client_key(xml);
+    // Get message mtype
+    std::string mtype = get_mtype(xml);
+    // Get the key of client who is to be contacted with the message
+    // Declare message
+    std::string msg = "";
+    // Determine which call has been invoqued
+    // Loop over all clients and send notification
+    for (int i = 0; i < m_clients.size(); ++i) {
+	    msg = "";
+	    
+	    // Set response for table.load.votable
+	    msg.append("<?xml version=\"1.0\"?>\n");
+	    msg.append("<methodCall><methodName>samp.client.receiveCall</methodName>\n");
+	    msg.append("<params>\n");
+	    msg.append("  <param><value>");
+	    msg.append(m_clients[i].private_key);
+	    msg.append("</value></param>\n");
+	    msg.append("  <param><value>"+m_hub_id+"</value></param>\n");
+	    msg.append("  <param><value><struct>\n");
+	    msg.append("    <member>\n");
+	    msg.append("      <name>samp.mtype</name>\n");
+	    msg.append("      <value>table.load.votable</value>\n");
+	    msg.append("    </member>\n");
+	    msg.append("    <member>\n");
+	    msg.append("      <name>samp.params</name>\n");
+	    msg.append("      <value><struct>\n");
+	    msg.append("        <member>\n");
+	    msg.append("          <name>name</name>\n");
+	    msg.append("          <value>imagename</value>\n");
+	    msg.append("        </member>\n");
+	    msg.append("        <member>\n");
+	    msg.append("          <name>table-id</name>\n");
+	    msg.append("          <value>id of table</value>\n");
+	    msg.append("        </member>\n");
+	    msg.append("        <member>\n");
+	    msg.append("           <name>url</name>\n");
+	    msg.append("           <value>url to the table</value>\n");
+	    msg.append("         </member>\n");
+	    msg.append("       </struct></value>\n");
+	    msg.append("      </member>\n");				
+	    msg.append("      </struct></value>\n");
+	    msg.append("    </param>\n");
+	    msg.append("</params>\n");
+	    msg.append("</methodCall>\n");
+                        
+	
+            // Post notification
+            notify(m_clients[i].url, msg);
+	}
+    
+    // Return
+    return;
+  
+}
+
+/***********************************************************************//**
+ * @brief Handles Hub shutdown requests
+ *
+ * @param[in] sock Socket.
+ *
+ * Handles all incoming Hub shutdown requests.
+ ***************************************************************************/
+void GVOHub::request_shutdown(const socklen_t& sock)
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::request_shutdown" << std::endl;
+    #endif
+
+    // Set shutdown flag
+    m_shutdown = true;
+    
+    // Post SAMP ok
+    post_samp_ok(sock);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Extract client key from XML request
+ *
+ * @param[in] xml XML message sent by client.
+ *
+ * Extracts the clients key from the XML request.
+ ***************************************************************************/
+std::string GVOHub::get_client_key(const GXml& xml) const
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::get_client_key" << std::endl;
+    #endif
+
+    // Initialise response
+    std::string client_key = "";
+
+    // Get the client's private key
+    const GXmlNode* node = xml.element("methodCall > params > param > value");
+    if (node != NULL) {
+        const GXmlText* text = static_cast<const GXmlText*>((*node)[0]);
+        if (text != NULL) {
+            client_key = text->text();
+        }
+    }
+
+    // Return key
+    return client_key;
+}
+
+
+/***********************************************************************//**
+ * @brief Extract client index in shared memory from XML request
+ *
+ * @param[in] xml XML message sent by client.
+ * @return Client index (-1 if not found).
+ *
+ * Extracts the client index from the XML request. The method returns -1 if
+ * no client was found.
+ ***************************************************************************/
+int GVOHub::get_client_index(const GXml& xml) const
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::get_client_index" << std::endl;
+    #endif
+
+    // Initialise index
+    int index = -1;
+
+    // Get the client's private key
+    std::string key = get_client_key(xml);
+
+    // Continue only if string is not empty
+    if (!key.empty()) {
+
+        // Searches for key in clients
+        for (int i = 0; i < m_clients.size(); ++i) {
+            if (m_clients[i].private_key == key) {
+                index = i;
+                break;
+            }
+        }
+
+    } // endif: key was not empty
+
+    // Return index
+    return index;
+}
+
+
+/***********************************************************************//**
+ * @brief Extract client index in shared memory from XML request
+ *
+ * @param[in] reference Client reference.
+ * @return Client index (-1 if not found).
+ *
+ * Extracts the client index from the client reference. The method returns -1
+ * if no client was found.
+ ***************************************************************************/
+int GVOHub::get_client_index(const std::string& reference) const
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::get_client_index" << std::endl;
+    #endif
+
+    // Initialise index
+    int index = -1;
+
+    // Search for index
+    for (int i = 0; i < m_clients.size(); ++i) {
+        if (m_clients[i].reference == reference) {
+            index = i;
+            break;
+        }
+    }
+
+    // Return index
+    return index;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns value for a SAMP client query parameter
+ *
+ * @param[in] xml client query XML document.
+ * @param[in] name Parameter name.
+ * @return Parameter value.
+ *
+ * Returns value for a SAMP client query parameter. If the specified
+ * parameter was not found or if the response structure is not compliant,
+ * an empty string is returned.
+ ***************************************************************************/
+std::string GVOHub::get_response_value(const GXml&        xml,
+                                       const std::string& name) const
+{
+    // Declare value
+    std::string value = "";
+
+    // Get the client's private key
+    const GXmlNode* node = xml.element("methodCall > params > param[1] > value > struct");
+    if (node != NULL) {
+        int num = node->elements("member");            
+        for (int i = 0; i < num; ++i) {
+            const GXmlNode* member = node->element("member", i);
+            std::string one_name;
+            std::string one_value;
+            get_name_value_pair(member, one_name, one_value);
+            if (one_name == name) {
+                value = one_value;
+                break;
+            }
+        }
+        
+    }
+
+    // Return value
+    return value;
+}
+
+
+/***********************************************************************//**
+ * @brief Extract name / value pair from XML node
+ *
+ * @param[in] node Pointer to XML node.
+ * @param[out] name Name string.
+ * @param[out] value Value string.
+ *
+ * Extracts a name / value pair from a XML node. If the XML node pointer is
+ * NULL, the name and value strings will be empty.
+ ***************************************************************************/
+void GVOHub::get_name_value_pair(const GXmlNode* node,
+                                 std::string&    name,
+                                 std::string&    value) const
+{
+    // Clear name and value strings
+    name.clear();
+    value.clear();
+
+    // Continue only if node is valid
+    if (node != NULL) {
+
+        // Get name node and extract text content
+        const GXmlNode* ptr = node->element("name", 0);
+        if (ptr != NULL) {
+            const GXmlText* text = static_cast<const GXmlText*>((*ptr)[0]);
+            if (text != NULL) {
+                name = text->text();
+            }
+        }
+
+        // Get value node and extract text content
+        ptr = node->element("value", 0);
+        if (ptr != NULL) {
+            const GXmlText* text = static_cast<const GXmlText*>((*ptr)[0]);
+            if (text != NULL) {
+                value = text->text();
+            }
+        }
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns subscriptions from XML document
+ *
+ * @param[in] xml client query XML document.
+ * @return List of subscriptions.
+ ***************************************************************************/
+std::vector<std::string> GVOHub::get_subscriptions(const GXml& xml) const                                         
+{
+    // Declare subscriptions
+    std::vector<std::string> subscriptions;
+
+    // Get the client's private key
+    const GXmlNode* node = xml.element("methodCall > params > param[1] > value > struct");
+    if (node != NULL) {
+        int num = node->elements("member");
+        for (int i = 0; i < num; ++i) {
+            const GXmlNode* member = node->element("member", i);
+            if (member != NULL) {
+                const GXmlNode* name = member->element("name", 0);
+                if (name != NULL) {
+                    const GXmlText* text = static_cast<const GXmlText*>((*name)[0]);
+                    if (text != NULL) {
+                        subscriptions.push_back(text->text());
+                    }
+                }
+            }
+        }
+    }
+
+    // Return subscriptions
+    return subscriptions;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns callback URL of client
+ *
+ * @param[in] xml client query XML document.
+ * @return Callback URL of client.
+ ***************************************************************************/
+std::string GVOHub::get_callback_url(const GXml& xml) const                                         
+{
+    // Declare URL
+    std::string url;
+
+    // Get the client's private key
+    const GXmlNode* node = xml.element("methodCall > params > param[1] > value");
+    if (node != NULL) {
+		const GXmlText* text = static_cast<const GXmlText*>((*node)[0]);
+        if (text != NULL) {
+		    url.append(text->text());
+        }
+    }
+
+    // Return URL
+    return url;
+}
+
+
+/***********************************************************************//**
+ * @brief Returns SAMP Hub lockfile URL
+ *
+ * @return SAMP Hub lockfile URL (empty if no lockfile was found).
+ *
+ * Implements IVOA standard REC-SAMP-1.3-20120411.
+ ***************************************************************************/
+std::string GVOHub::get_hub_lockfile(void) const
+{
+    // Initialise result
+    std::string url = "";
+
+    // Check for existence of the SAMP_HUB environment variable first
+    char* hub_ptr = std::getenv("SAMP_HUB");
+    if (hub_ptr != NULL) {
+
+        // Check for mandatory std-lockurl: prefix (no other prefixe is
+        // supported so far)
+        std::string lockurl = std::string(hub_ptr);
+        if (lockurl.compare(0, 12, "std-lockurl:") == 0) {
+
+            // Extract URL
+            url = lockurl.substr(12, std::string::npos);
+
+        } // endif: std-lockurl: prefix found
+  
+    }
+
+    // ... otherwise the lockfile should be $HOME/.samp
+    else {
+
+        // Get user's HOME directory path as the prefix of the full
+        // path. If the HOME environment variable is not set we
+        // expect that .samp is in the local directory. This is non
+        // standard, but prevents for creating an exception here.
+        std::string prefix = "";
+        char* home_ptr = std::getenv("HOME");
+        if (home_ptr != NULL) {
+            prefix = std::string(home_ptr) + "/";
+        }
+
+        // Set filename
+        url = prefix + ".samp";
+
+    } // endelse: no SAMP_HUB environment variable found
+
+    // Return URL
+    return url;
+}
+
+
+/***********************************************************************//**
+ * @brief Send notifications
+ *
+ * @param[in] method Name of the method
+ * @param[in] client Reference of the client
+ *
+ * The method send notifications to all registered clients that request so.
+ ***************************************************************************/
+void GVOHub::send_notifications(const std::string& method,
+                                const std::string& client)
+{
+    // Loop over all clients
+    for (int i = 0; i < m_clients.size(); ++i) {
+    
+        // Loop over all methods
+        for (int j = 0; j < m_clients[i].subscriptions.size(); ++j) {
+
+            // Continue only if method matches
+            if (method == m_clients[i].subscriptions[j]) {
+    
+                // Declare notification
+                std::string msg = "";
+
+                // Register
+                if (method == "samp.hub.event.register") {
+
+                    // Set response for samp.hub.event.register
+                    msg.append("<?xml version=\"1.0\"?>\n");
+                    msg.append("<methodCall>\n");
+                    msg.append("<methodName>samp.client.receiveNotification</methodName>\n");
+                    msg.append("<params>\n");
+                    msg.append("  <param><value>");
+                    msg.append(m_clients[i].private_key);
+                    msg.append("</value></param>\n");
+                    msg.append("  <param><value>"+m_hub_id+"</value></param>\n");
+                    msg.append("  <param><value><struct>\n");
+                    msg.append("    <member>\n");
+                    msg.append("      <name>samp.mtype</name>\n");
+                    msg.append("      <value>samp.hub.event.register</value>\n");
+                    msg.append("    </member>\n");
+                    msg.append("    <member>\n");
+                    msg.append("      <name>samp.params</name>\n");
+                    msg.append("      <value><struct>\n");
+                    msg.append("        <member>\n");
+                    msg.append("          <name>id</name>\n");
+                    msg.append("          <value>"+client+"</value>\n");
+                    msg.append("        </member>\n");
+                    msg.append("      </struct></value>\n");
+                    msg.append("    </member>\n");
+                    msg.append("  </struct></value></param>\n");
+                    msg.append("</params>\n");
+                    msg.append("</methodCall>\n");
+                    
+                    // Post notification
+                    notify(m_clients[i].url, msg);
+                    
+                }
+                    
+                // Unregister
+                else if (method == "samp.hub.event.unregister") {
+
+                    // Set response for samp.hub.event.unregister
+                    msg.append("<?xml version=\"1.0\"?>\n");
+                    msg.append("<methodCall>\n");
+                    msg.append("<methodName>samp.client.receiveNotification</methodName>\n");
+                    msg.append("<params>\n");
+                    msg.append("  <param><value>");
+                    msg.append(m_clients[i].private_key);
+                    msg.append("</value></param>\n");
+                    msg.append("  <param><value>");
+                    msg.append(m_hub_id);
+                    msg.append("</value></param>\n");
+                    msg.append("  <param><value><struct>\n");
+                    msg.append("    <member>\n");
+                    msg.append("      <name>samp.mtype</name>\n");
+                    msg.append("      <value>samp.hub.event.unregister</value>\n");
+                    msg.append("    </member>\n");
+                    msg.append("    <member>\n");
+                    msg.append("      <name>samp.params</name>\n");
+                    msg.append("      <value><struct>\n");
+                    msg.append("        <member>\n");
+                    msg.append("          <name>id</name>\n");
+                    msg.append("          <value>"+client+"</value>\n");
+                    msg.append("        </member>\n");
+                    msg.append("      </struct></value>\n");
+                    msg.append("    </member>\n");
+                    msg.append("  </struct></value></param>\n");
+                    msg.append("</params>\n");
+                    msg.append("</methodCall>\n");
+                        
+                    // Post notification
+                    notify(m_clients[i].url, msg);
+					
+                }
+                    
+                // Metadata
+			    else if (method == "samp.hub.event.metadata") {
+
+                    // Get index of client
+                    int k;
+                    for (k = 0; k < m_clients.size(); ++k) {
+                        if (m_clients[k].reference == client) {
+                            break;
+                        }
+                    }
+
+                    // Fall through if client index is not valid
+                    if (k >= m_clients.size()) {
+                        continue;
+                    }
+
+                    // Set response for samp.hub.event.metadata
+                    msg.append("<?xml version=\"1.0\"?>\n");
+                    msg.append("<methodCall><methodName>samp.client.receiveNotification</methodName>\n");
+                    msg.append("<params>\n");
+                    msg.append("  <param><value>");
+                    msg.append(m_clients[i].private_key);
+                    msg.append("</value></param>\n");
+                    msg.append("  <param><value>"+m_hub_id+"</value></param>\n");
+                    msg.append("  <param><value><struct>\n");
+                    msg.append("    <member>\n");
+                    msg.append("      <name>samp.mtype</name>\n");
+                    msg.append("      <value>samp.hub.event.metadata</value>\n");
+                    msg.append("    </member>\n");
+                    msg.append("    <member>\n");
+                    msg.append("      <name>samp.params</name>\n");
+                    msg.append("      <value><struct>\n");
+                    msg.append("        <member>\n");
+                    msg.append("          <name>id</name>\n");
+                    msg.append("          <value>"+m_clients[k].reference+"</value>\n");
+                    msg.append("        </member>\n");
+                    msg.append("        <member>\n");
+                    msg.append("          <name>metadata</name>\n");
+                    msg.append("          <value><struct>\n");
+                    msg.append("            <member>\n");
+                    msg.append("              <name>samp.name</name>\n");
+                    msg.append("              <value>"+m_clients[k].name+"</value>\n");
+                    msg.append("            </member>\n");
+                    msg.append("            <member>\n");
+                    msg.append("              <name>samp.description.text</name>\n");
+                    msg.append("              <value>"+m_clients[k].description+"</value>\n");
+                    msg.append("            </member>\n");
+                    msg.append("            <member>\n");
+                    msg.append("              <name>samp.icon.url</name>\n");
+                    msg.append("              <value>"+m_clients[k].icon+"</value>\n");
+                    msg.append("            </member>\n");
+                    msg.append("            <member>\n");
+                    msg.append("              <name>samp.documentation.url</name>\n");
+                    msg.append("              <value>"+m_clients[k].documentation+"</value>\n");
+                    msg.append("            </member>\n");
+                    msg.append("            <member>\n");
+                    msg.append("              <name>author.affiliation</name>\n");
+                    msg.append("              <value>"+m_clients[k].affiliation+"</value>\n");
+                    msg.append("            </member>\n");
+                    msg.append("            <member>\n");
+                    msg.append("              <name>author.name</name>\n");
+                    msg.append("              <value>"+m_clients[k].author_name+"</value>\n");
+                    msg.append("            </member>\n");
+                    msg.append("            <member>\n");
+                    msg.append("              <name>author.email</name>\n");
+                    msg.append("              <value>"+m_clients[k].email+"</value>\n");
+                    msg.append("            </member>\n");
+                    msg.append("            <member>\n");
+                    msg.append("              <name>home.page</name>\n");
+                    msg.append("              <value>"+m_clients[k].homepage+"</value>\n");
+                    msg.append("            </member>\n");
+                    msg.append("          </struct></value>\n");
+                    msg.append("        </member>\n");				
+                    msg.append("      </struct></value>\n");
+                    msg.append("    </member>\n");
+                    msg.append("  </struct></value></param>\n");
+                    msg.append("</params>\n");
+                    msg.append("</methodCall>\n");
+                        
+                    // Post notification
+                    notify(m_clients[i].url, msg);
+					
+                } // endelse: metadata
+
+            } // endif: method was requested one
+
+        } // endfor: looped over all methods
+
+    } // endfor: looped over all clients
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Creates the lockfile, fill it
+ *
+ * Implements IVOA standard REC-SAMP-1.3-20120411.
+ ***************************************************************************/
+void GVOHub::create_samp_file(void) const
+{
+    // Get lockfile URL
+    std::string lockurl = get_hub_lockfile();
+
+    // Open SAMP lockfile. Continue only if opening was successful
+    FILE* fptr = fopen(lockurl.c_str(), "w");
+
+    // Successful?
+    if (fptr != NULL) {
+
+        // Write lockfile
+        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.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());
+        
+        // Close lockfile
+        fclose(fptr);
+
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Post string content to client
+ *
+ * @param[in] content String content to post.
+ * @param[in] sock Socket.
+ *
+ * Posts the content of a string to a client.
+ ***************************************************************************/
+void GVOHub::post_string(const std::string& content, const socklen_t& sock) const
+{
+    // Continue only if socket is valid
+    if (sock != -1) {
+
+        // Determine content length
+        int length = content.length();
+
+        // Set prefix (note that the response has no POST!)
+        std::string prefix = "HTTP/1.1 200 OK\n"
+                             "Connection: close\n"
+                             "Content-Type: text/xml\n"
+                             "Content-Length: "+gammalib::str(length)+"\n\n";
+
+
+        // Build post string
+        std::string post = prefix + content;
+
+        // Dump message to post
+        #if defined(G_SHOW_MESSAGE)
+        std::cout << std::endl;
+        std::cout << "GVOHub response:" << std::endl;
+        std::cout << "================" << std::endl;
+        std::cout << post << std::endl;
+        #endif
+
+        // Send content to socket
+        bool done = false;
+        do {
+            int length      = post.length();
+            int sent_length = send(sock, post.c_str(), length, 0);
+            if (sent_length < length) {
+                post = post.substr(sent_length, std::string::npos);
+            }
+            else {
+                done = true;
+            }
+        } while (!done);
+	
+    } // endif: socket was valid
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Post SAMP ok massage to client
+ *
+ * @param[in] sock Socket.
+ *
+ * Posts a SAMP OK message to a client.
+ ***************************************************************************/
+void GVOHub::post_samp_ok(const socklen_t& sock) const
+{
+    // Declare response
+    std::string response = "";
+
+    // Compose response
+    response.append("<?xml version='1.0' encoding='UTF-8'?>\n");
+    response.append("<methodResponse>\n");
+    response.append("<params>\n");
+    response.append("  <param><value><struct>\n");
+    response.append("    <member>\n");
+    response.append("      <name>samp.status</name>\n");
+    response.append("      <value>samp.ok</value>\n");
+    response.append("    </member>\n");
+    response.append("  </struct></value></param>\n");
+    response.append("</params>\n");
+    response.append("</methodResponse>\n");
+    
+    // Post string
+    post_string(response, sock);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Post SAMP void massage to client
+ *
+ * @param[in] sock Socket.
+ *
+ * Posts a void SAMP message to a client. Void messages can typically be
+ * discared by the client.
+ ***************************************************************************/
+void GVOHub::post_samp_void(const socklen_t& sock) const
+{
+    // Declare response
+    std::string response = "";
+
+    // Compose response
+    response.append("<?xml version='1.0' encoding='UTF-8'?>\n");
+    response.append("<methodResponse>\n");
+    response.append("  <params>\n");
+    response.append("    <param>\n");
+    response.append("      <value/>\n");
+    response.append("    </param>\n");
+    response.append("  </params>\n");
+    response.append("</methodResponse>\n");
+    
+    // Post string
+    post_string(response, sock);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Send notification to client
+ *
+ * @param[in] url URL.
+ * @param[in] notification Notification.
+ *
+ * Send a @p notification to a @p url.
+ ***************************************************************************/
+void GVOHub::notify(const std::string& url,
+                    const std::string& notification) const
+{
+    // Initalise socket
+    int socket = -1;
+
+    // Extract host, port and path from URL
+    std::string host = "127.0.0.1";
+    std::string port = "2525";
+    std::string path = "xmlrpc";
+    if (url.compare(0, 7, "http://") == 0) {
+        size_t length;
+        size_t start = 7;
+        size_t stop  = url.find(":", start);
+        size_t end   = std::string::npos;
+        if (stop != std::string::npos) {
+            length = stop - start;
+        }
+        else {
+            length = std::string::npos;
+        }
+        host = url.substr(start, length);
+        if (stop != std::string::npos) {
+            stop++;
+            end  = url.find("/", stop);
+            if (end != std::string::npos) {
+                length = end - stop;
+            }
+            else {
+                length = std::string::npos;
+            }
+            port = url.substr(stop, length);
+        }
+        if (end != std::string::npos) {
+            end++;
+            length = url.length() - end;
+            path   = url.substr(end, length);
+        }
+    }
+
+    // Set hints
+    struct addrinfo hints;
+    std::memset(&hints, 0, sizeof(hints));
+    hints.ai_family   = AF_INET;
+    hints.ai_socktype = SOCK_STREAM;
+
+    // Get server information
+    struct addrinfo* servinfo;
+    if (getaddrinfo(host.c_str(), port.c_str(), &hints, &servinfo) == 0) {
+
+        // Loop through all the results and connect to the first we can
+        for (struct addrinfo* ptr = servinfo; ptr != NULL; ptr = ptr->ai_next) {
+
+            // Create socket
+            socket = ::socket(ptr->ai_family,
+                              ptr->ai_socktype,
+                              ptr->ai_protocol);
+
+            // Connect to socket if socket is valid
+            if (socket != -1) {
+                if (::connect(socket,
+                              ptr->ai_addr,
+                              ptr->ai_addrlen) == -1) {
+                    close(socket);
+                    socket = -1;
+                }
+                else {
+                    // Socket connection successful, break now
+                    break;
+                }
+            }
+
+        } // endfor: looped through all results
+
+    } // endif: server information was valid
+
+    // Continue only if socket is valid
+    if (socket != -1) {
+
+        // Determine notification length
+        int length = notification.length();
+
+        // Set prefix
+        std::string prefix = "POST /"+path+" HTTP/1.0\n"
+                             "Connection: close\n"
+                             "User-Agent: GammaLib\n"
+                             "Content-Type: text/xml\n"
+                             "Content-Length: "+gammalib::str(length)+"\n\n";
+
+        // Build post string
+        std::string post = prefix + notification;
+
+        // Dump message to post
+        #if defined(G_SHOW_MESSAGE)
+        std::cout << std::endl;
+        std::cout << "GVOHub sends following notification:" << std::endl;
+        std::cout << "====================================" << std::endl;
+        std::cout << post << std::endl;
+        #endif
+
+        // Send content to socket
+        bool done = false;
+        do {
+            int length      = post.length() + 1; // +1 for terminating 0
+            int sent_length = send(socket, post.c_str(), length, 0);
+            if (sent_length < length) {
+                post = post.substr(sent_length, std::string::npos);
+            }
+            else {
+                done = true;
+            }
+        } while (!done);
+
+        /*
+        // Declare response
+        std::string response = "";
+
+        // Read buffer until it is empty
+        char buffer[1001];
+        int timeout = 1000; // Initial timeout is 1 sec
+        int n       = 0;
+        do {
+            n = gammalib::recv(socket, buffer, 1000, 0, timeout);
+            if (n > 0) {
+                buffer[n+1] = '\0';
+                response.append(std::string(buffer));
+            }
+            timeout = 10; // The timeout now is 0.01 sec 
+        } while (n > 0);
+        */
+
+        // Close socket
+        close(socket);
+
+    }
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Generates random string of characters
+ *
+ * @param[in] length Length of random string
+ ***************************************************************************/
+std::string GVOHub::random_string(const size_t& length) const
+{
+    srand(time(0));
+    std::string str = "0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz";
+    int pos;
+    while(str.size() != length) {
+        pos = ((rand() % (str.size() - 1)));
+        str.erase (pos, 1);
+    }
+   
+    // Return string
+    return str;
+}
+/***********************************************************************//**
+ * @brief Extract mtype from XML request
+ *
+ * @param[in] xml XML message sent by client.
+ *
+ * Extracts the mtype identifying calling message from the XML request.
+ ***************************************************************************/
+std::string GVOHub::get_mtype(const GXml& xml) const
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::get_mtype" << std::endl;
+    #endif
+
+    // Initialise response
+    std::string client_key = "";
+
+    // Get the client's private key
+    const GXmlNode* node = xml.element("methodCall > params > param > value");
+    if (node != NULL) {
+        const GXmlText* text = static_cast<const GXmlText*>((*node)[0]);
+        if (text != NULL) {
+            client_key = text->text();
+        }
+    }
+
+    // Return key
+    return client_key;
+}
+/***********************************************************************//**
+ * @brief Extract client destination from XML request
+ *
+ * @param[in] xml XML message sent by client.
+ *
+ * Extracts identifier of destination client from the XML request.
+ ***************************************************************************/
+std::string GVOHub::get_destination(const GXml& xml) const
+{
+    // Header
+    #if defined(G_CONSOLE_DUMP)
+    std::cout << "GVOHub::get_mtype" << std::endl;
+    #endif
+
+    // Initialise response
+    std::string client_key = "";
+
+    // Get the client's private key
+    const GXmlNode* node = xml.element("methodCall > params > param > value");
+    if (node != NULL) {
+        const GXmlText* text = static_cast<const GXmlText*>((*node)[0]);
+        if (text != NULL) {
+            client_key = text->text();
+        }
+    }
+
+    // Return key
+    return client_key;
+}
+/***************************************************************************
+ *                     GVOTable.cpp - VO Table composition class           *
+ * ----------------------------------------------------------------------- *
+ *  copyright (C) 2015 by Thierry Louge                               *
+ * ----------------------------------------------------------------------- *
+ *                                                                         *
+ *  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   *
+ *  the Free Software Foundation, either version 3 of the License, or      *
+ *  (at your option) any later version.                                    *
+ *                                                                         *
+ *  This program is distributed in the hope that it will be useful,        *
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of         *
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          *
+ *  GNU General Public License for more details.                           *
+ *                                                                         *
+ *  You should have received a copy of the GNU General Public License      *
+ *  along with this program.  If not, see <http://www.gnu.org/licenses/>.  *
+ *                                                                         *
+ ***************************************************************************/
+/**
+ * @file GVOTable.cpp
+ * @brief VO Table composition class
+ * @author Thierry Louge
+ */
+
+/* __ Includes ___________________________________________________________ */
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+#include <cstdlib>         // std::getenv() function
+#include <cstdio>          // std::fopen(), etc. functions
+#include <cstring>         // std::memset() function
+#include <csignal>         // signal() function
+#include <cerrno>          // errno
+#include <unistd.h>        // close() function
+#include <netdb.h>         // getaddrinfo() function
+#include <netinet/in.h>    // sockaddr_in, INADDR_ANY, htons
+#include <fstream>
+#include <sys/shm.h>
+#include <sys/socket.h>    // socket(), connect() functions
+#include <sys/wait.h>      // waitpid() function
+#include <arpa/inet.h>     // inet_addr() function
+#include "GVOTable.hpp"
+#include "GException.hpp"
+#include "GTools.hpp"
+
+/* __ Method name definitions ____________________________________________ */
+
+/* __ Macros _____________________________________________________________ */
+
+/* __ Coding definitions _________________________________________________ */
+
+/* __ Debug definitions __________________________________________________ */
+//#define G_CONSOLE_DUMP
+#define G_SHOW_MESSAGE
+
+
+/*==========================================================================
+ =                                                                         =
+ =                         Constructors/destructors                        =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Void constructor
+ ***************************************************************************/
+GVOTable::GVOTable(void)
+{
+    // Initialise members
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy constructor
+ *
+ * @param[in] VOTable
+ ***************************************************************************/
+GVOTable::GVOTable(const GVOTable& votable)
+{
+    // Initialise members
+    init_members();
+
+    // Copy members
+    copy_members(votable);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Destructor
+ ***************************************************************************/
+GVOTable::~GVOTable(void)
+{
+    // Free members
+    free_members();
+
+    // Return
+    return;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                                Operators                                =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Assignment operator
+ *
+ * @param[in] VO table.
+ * @return VO table.
+ ***************************************************************************/
+GVOTable& GVOTable::operator=(const GVOTable& votable)
+{
+    // Execute only if object is not identical
+    if (this != &votable) {
+
+        // Free members
+        free_members();
+
+        // Initialise members
+        init_members();
+
+        // Copy members
+        copy_members(votable);
+
+    } // endif: object was not identical
+
+    // Return
+    return *this;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                             Public methods                              =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Clear object.
+ *
+ * Reset object to a clean initial state.
+ ***************************************************************************/
+void GVOTable::clear(void)
+{
+    // Free members
+    free_members();
+
+    // Initialise members
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Clone object
+ ***************************************************************************/
+GVOTable* GVOTable::clone(void) const
+{
+    // Clone client
+    return new GVOTable(*this);
+}
+
+
+/***********************************************************************//**
+ * @brief Print VO hub information
+ *
+ * @param[in] chatter Chattiness (defaults to NORMAL).
+ * @return String containing VO table information
+ ***************************************************************************/
+std::string GVOTable::print(const GChatter& chatter) const
+{
+    // Initialise result string
+    std::string result;
+
+    // Continue only if chatter is not silent
+    if (chatter != SILENT) {
+
+        // Append header
+        result.append("=== GVOTable ===");
+
+        // Append table information
+
+    } // endif: chatter was not silent
+
+    // Return result
+    return result;
+}
+/***********************************************************************//**
+ * @brief Creates the VO table, fill it
+ *
+ * Implements IVOA standard <standard id>
+ ***************************************************************************/
+void GVOTable::create_votable(void)
+{
+    std::string header = "";
+    header.append("<?xml version=\"1.0\"?>\n");
+    header.append("<VOTABLE version=\"1.3\" xmlns:xsi=\"http://www.w3.org/2001/XMLSchema-instance\"\n");
+    header.append("xmlns=\"http://www.ivoa.net/xml/VOTable/v1.3\"\n");
+    //Change this with relevant datamodel
+    header.append("xmlns:stc=\"http://www.ivoa.net/xml/STC/v1.30\">\n");
+    header.append("\t<RESOURCE name=\"GammalibSkymap\">\n");
+    header.append("\t\t<TABLE name=\"skymap\">\n");
+    header.append("\t\t\t<DESCRIPTION>Skymap coming from gammalib</DESCRIPTION>\n");
+    //Only for testing, we close the document.
+    header.append("\t\t</TABLE>\n");
+    header.append("\t</RESOURCE>\n");
+    header.append("</VOTABLE>\n");
+    
+    m_tablexml = GXml(header);
+    m_tablexml.print(NORMAL,0);
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Fills content of the VO table
+ *
+ * Implements IVOA standard <standard id>
+ ***************************************************************************/
+void GVOTable::fill_content(const std::string& data) const
+{
+
+    // Return
+    return;
+}
+
+/*==========================================================================
+ =                                                                         =
+ =                             Private methods                             =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Initialise class members
+ ***************************************************************************/
+void GVOTable::init_members(void)
+{
+    // Initialise members
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy class members
+ *
+ * @param[in] hub VO hub.
+ ***************************************************************************/
+void GVOTable::copy_members(const GVOTable& hub)
+{
+    // Copy members
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Delete class members
+ ***************************************************************************/
+void GVOTable::free_members(void)
+{ 
+    
+    // Return
+    return;
+}
+
+/***************************************************************************
+ *                     GVOTable.cpp - VO Table composition class           *
+ * ----------------------------------------------------------------------- *
+ *  copyright (C) 2015 by Thierry Louge                               *
+ * ----------------------------------------------------------------------- *
+ *                                                                         *
+ *  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   *
+ *  the Free Software Foundation, either version 3 of the License, or      *
+ *  (at your option) any later version.                                    *
+ *                                                                         *
+ *  This program is distributed in the hope that it will be useful,        *
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of         *
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          *
+ *  GNU General Public License for more details.                           *
+ *                                                                         *
+ *  You should have received a copy of the GNU General Public License      *
+ *  along with this program.  If not, see <http://www.gnu.org/licenses/>.  *
+ *                                                                         *
+ ***************************************************************************/
+/**
+ * @file GVOTable.cpp
+ * @brief VO Table composition class
+ * @author Thierry Louge
+ */
+
+/* __ Includes ___________________________________________________________ */
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+#include <cstdlib>         // std::getenv() function
+#include <cstdio>          // std::fopen(), etc. functions
+#include <cstring>         // std::memset() function
+#include <csignal>         // signal() function
+#include <cerrno>          // errno
+#include <unistd.h>        // close() function
+#include <netdb.h>         // getaddrinfo() function
+#include <netinet/in.h>    // sockaddr_in, INADDR_ANY, htons
+#include <fstream>
+#include <sys/shm.h>
+#include <sys/socket.h>    // socket(), connect() functions
+#include <sys/wait.h>      // waitpid() function
+#include <arpa/inet.h>     // inet_addr() function
+#include "GVOTable.hpp"
+#include "GException.hpp"
+#include "GTools.hpp"
+
+/* __ Method name definitions ____________________________________________ */
+
+/* __ Macros _____________________________________________________________ */
+
+/* __ Coding definitions _________________________________________________ */
+
+/* __ Debug definitions __________________________________________________ */
+//#define G_CONSOLE_DUMP
+#define G_SHOW_MESSAGE
+
+
+/*==========================================================================
+ =                                                                         =
+ =                         Constructors/destructors                        =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Void constructor
+ ***************************************************************************/
+GVOTable::GVOTable(void)
+{
+    // Initialise members
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy constructor
+ *
+ * @param[in] VOTable
+ ***************************************************************************/
+GVOTable::GVOTable(const GVOTable& votable)
+{
+    // Initialise members
+    init_members();
+
+    // Copy members
+    copy_members(votable);
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Destructor
+ ***************************************************************************/
+GVOTable::~GVOTable(void)
+{
+    // Free members
+    free_members();
+
+    // Return
+    return;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                                Operators                                =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Assignment operator
+ *
+ * @param[in] VO table.
+ * @return VO table.
+ ***************************************************************************/
+GVOTable& GVOTable::operator=(const GVOTable& votable)
+{
+    // Execute only if object is not identical
+    if (this != &votable) {
+
+        // Free members
+        free_members();
+
+        // Initialise members
+        init_members();
+
+        // Copy members
+        copy_members(votable);
+
+    } // endif: object was not identical
+
+    // Return
+    return *this;
+}
+
+
+/*==========================================================================
+ =                                                                         =
+ =                             Public methods                              =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Clear object.
+ *
+ * Reset object to a clean initial state.
+ ***************************************************************************/
+void GVOTable::clear(void)
+{
+    // Free members
+    free_members();
+
+    // Initialise members
+    init_members();
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Clone object
+ ***************************************************************************/
+GVOTable* GVOTable::clone(void) const
+{
+    // Clone client
+    return new GVOTable(*this);
+}
+
+
+/***********************************************************************//**
+ * @brief Print VO hub information
+ *
+ * @param[in] chatter Chattiness (defaults to NORMAL).
+ * @return String containing VO table information
+ ***************************************************************************/
+std::string GVOTable::print(const GChatter& chatter) const
+{
+    // Initialise result string
+    std::string result;
+
+    // Continue only if chatter is not silent
+    if (chatter != SILENT) {
+
+        // Append header
+        result.append("=== GVOTable ===");
+
+        // Append table information
+
+    } // endif: chatter was not silent
+
+    // Return result
+    return result;
+}
+/***********************************************************************//**
+ * @brief Creates the VO table, fill it
+ *
+ * Implements IVOA standard <standard id>
+ ***************************************************************************/
+void GVOTable::create_votable(void)
+{
+    std::string header = "";
+    header.append("<?xml version=\"1.0\"?>\n");
+    header.append("<VOTABLE version=\"1.3\" xmlns:xsi=\"http://www.w3.org/2001/XMLSchema-instance\"\n");
+    header.append("xmlns=\"http://www.ivoa.net/xml/VOTable/v1.3\"\n");
+    //Change this with relevant datamodel
+    header.append("xmlns:stc=\"http://www.ivoa.net/xml/STC/v1.30\">\n");
+    header.append("\t<RESOURCE name=\"GammalibSkymap\">\n");
+    header.append("\t\t<TABLE name=\"skymap\">\n");
+    header.append("\t\t\t<DESCRIPTION>Sloaded_map = skymap.load("kymap coming from gammalib</DESCRIPTION>\n");
+    //Only for testing, we close the document.
+    header.append("\t\t</TABLE>\n");
+    header.append("\t</RESOURCE>\n");
+    header.append("</VOTABLE>\n");
+    
+    m_tablexml = GXml(header);
+    m_tablexml.print(NORMAL,0);
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Fills content of the VO table
+ *
+ * Implements IVOA standard <standard id>
+ ***************************************************************************/
+void GVOTable::fill_content(const std::string& data) const
+{
+
+    // Return
+    return;
+}
+
+/*==========================================================================
+ =                                                                         =
+ =                             Private methods                             =
+ =                                                                         =
+ ==========================================================================*/
+
+/***********************************************************************//**
+ * @brief Initialise class members
+ ***************************************************************************/
+void GVOTable::init_members(void)
+{
+    // Initialise members
+    
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Copy class members
+ *
+ * @param[in] hub VO hub.
+ ***************************************************************************/
+void GVOTable::copy_members(const GVOTable& hub)
+{
+    // Copy members
+
+    // Return
+    return;
+}
+
+
+/***********************************************************************//**
+ * @brief Delete class members
+ ***************************************************************************/
+void GVOTable::free_members(void)
+{ 
+    
+    // Return
+    return;
+}
+
@@ -7,7 +7,8 @@ INCLUDES = -I$(top_srcdir)/include \
  
 # Define sources for this directory
 sources = GVOClient.cpp \
-          GVOHub.cpp
+          GVOHub.cpp \
+	  GVOTable.cpp
  
 # Build libtool library
 noinst_LTLIBRARIES = libvo.la