Add binary operators to GSkymap class.

Jürgen Knödlseder
1 parent bca7af9e
Showing 3 changed files with 96 additions and 0 deletions
include/GSkymap.hpp
pyext/GSkymap.i
test/test_GSky.cpp
@@ -112,6 +112,10 @@ public:
     GSkymap&      operator*=(const double& factor);
     GSkymap&      operator/=(const GSkymap& map);
     GSkymap&      operator/=(const double& factor);
+    GSkymap       operator+(const GSkymap& map) const;
+    GSkymap       operator-(const GSkymap& map) const;
+    GSkymap       operator*(const GSkymap& map) const;
+    GSkymap       operator/(const GSkymap& map) const;
     double&       operator()(const int& index, const int& map = 0);
     const double& operator()(const int& index, const int& map = 0) const;
     double&       operator()(const GSkyPixel& pixel, const int& map = 0);
@@ -188,6 +192,74 @@ private:
  
  
 /***********************************************************************//**
+ * @brief Binary sky map addition
+ *
+ * @param[in] map Sky map.
+ * @return Sky map to which @p map was added.
+ *
+ * Returns the sum of two sky maps.
+ ***************************************************************************/
+inline
+GSkymap GSkymap::operator+(const GSkymap& map) const
+{
+    GSkymap result = *this;
+    result        += map;
+    return result;
+}
+
+
+/***********************************************************************//**
+ * @brief Binary sky map subtraction
+ *
+ * @param[in] map Sky map.
+ * @return Sky map to which @p map was added.
+ *
+ * Returns the difference of two sky maps.
+ ***************************************************************************/
+inline
+GSkymap GSkymap::operator-(const GSkymap& map) const
+{
+    GSkymap result = *this;
+    result        -= map;
+    return result;
+}
+
+
+/***********************************************************************//**
+ * @brief Binary sky map multiplication
+ *
+ * @param[in] map Sky map.
+ * @return Sky map multiplied by @p map.
+ *
+ * Returns the product of two sky maps.
+ ***************************************************************************/
+inline
+GSkymap GSkymap::operator*(const GSkymap& map) const
+{
+    GSkymap result = *this;
+    result        *= map;
+    return result;
+}
+
+
+/***********************************************************************//**
+ * @brief Binary sky map division
+ *
+ * @param[in] map Sky map.
+ * @return Sky map divided by @p map.
+ *
+ * Returns the ratio of two sky maps.
+ ***************************************************************************/
+inline
+GSkymap GSkymap::operator/(const GSkymap& map) const
+{
+    GSkymap result = *this;
+    result        /= map;
+    return result;
+}
+
+
+/***********************************************************************//**
  * @brief Return class name
  *
  * @return String containing the class name ("GSkymap").
@@ -134,6 +134,10 @@ public:
     GSkymap&      operator-=(const double& value);
     GSkymap&      operator*=(const GSkymap& map);
     GSkymap&      operator*=(const double& factor);
+    GSkymap       operator+(const GSkymap& map) const;
+    GSkymap       operator-(const GSkymap& map) const;
+    GSkymap       operator*(const GSkymap& map) const;
+    GSkymap       operator/(const GSkymap& map) const;
     double&       operator()(const int& index, const int& map = 0);
     double&       operator()(const GSkyPixel& pixel, const int& map = 0);
     double        operator()(const GSkyDir& dir, const int& map = 0) const;
@@ -831,12 +831,15 @@ void TestGSky::test_GSkymap_wcs_io(void)
  
 /***************************************************************************
  * @brief GSkymap
+ *
+ * This method tests the GSkymap class.
  ***************************************************************************/
 void TestGSky::test_GSkymap(void)
 {
     // Define several maps for comparison
     GSkymap map_src("CAR", "GAL", 0.0, 0.0, -1.0, 1.0, 10, 10, 2);
     GSkymap map_dst("CAR", "GAL", 0.0, 0.0, -0.1, 0.1, 100, 100, 2);
+    GSkymap map_new("CAR", "GAL", 0.0, 0.0, -0.1, 0.1, 100, 100, 2);
  
     // Fill map pixels
     double total_src = 0.0;
@@ -848,61 +851,78 @@ void TestGSky::test_GSkymap(void)
     }
  
     // Add pixels to destination map
+    map_new  = map_dst + map_src;
     map_dst += map_src;
  
     // Check total in destination map. Note that the total in the destination
     // map should be 100 times the total in the source map as the operator
     // is expected to perform strict bi-linear interpolation
     double total_dst = 0.0;
+    double total_new = 0.0;
     double total_ref = 0.0;
     for (int pix = 0; pix < map_dst.npix(); ++pix) {
         for (int k = 0; k < map_dst.nmaps(); ++k) {
             total_dst += map_dst(pix,k);
+            total_new += map_new(pix,k);
             total_ref += map_src(map_dst.pix2dir(pix),k);
         }
     }
     total_dst /= 100.0;
+    total_new /= 100.0;
     total_ref /= 100.0;
 	test_value(total_dst, total_ref, 1.0e-3, "Test operator+=(GSkymap)");
+	test_value(total_new, total_ref, 1.0e-3, "Test operator+(GSkymap)");
  
     // Subtract pixels from destination map
+    map_new  = map_dst - map_src;
     map_dst -= map_src;
  
     // Check total in destination map. Note that the total in the destination
     // map should be zero.
     total_dst = 0.0;
+    total_new = 0.0;
     for (int pix = 0; pix < map_dst.npix(); ++pix) {
         for (int k = 0; k < map_dst.nmaps(); ++k) {
             total_dst += map_dst(pix,k);
+            total_new += map_dst(pix,k);
         }
     }
 	test_value(total_dst, 0.0, 1.0e-3, "Test operator-=(GSkymap)");
+	test_value(total_new, 0.0, 1.0e-3, "Test operator-(GSkymap)");
  
     // Check map multiplication
     GSkymap test_map  = map_src;
     test_map         *= map_src;
+    map_new           = map_src * map_src;
     double total_test = 0.0;
     total_ref         = 0.0;
+    total_new         = 0.0;
     for (int pix = 0; pix < test_map.npix(); ++pix) {
         for (int k = 0; k < test_map.nmaps(); ++k) {
             total_test += test_map(pix,k);
+            total_new  += map_new(pix,k);
             total_ref  += map_src(pix,k) * map_src(pix,k);
         }
     }
 	test_value(total_test, total_ref, 1.0e-3, "Test operator*=(GSkymap)");
+	test_value(total_new,  total_ref, 1.0e-3, "Test operator*(GSkymap)");
  
     // Check map division
     test_map   = map_src;
     test_map  /= map_src;
+    map_new    = map_src / map_src;
     total_test = 0.0;
     total_ref  = 0.0;
+    total_new  = 0.0;
     for (int pix = 0; pix < test_map.npix(); ++pix) {
         for (int k = 0; k < test_map.nmaps(); ++k) {
             total_test += test_map(pix,k);
+            total_new  += map_new(pix,k);
             total_ref  += map_src(pix,k) / map_src(pix,k);
         }
     }
 	test_value(total_test, total_ref, 1.0e-3, "Test operator/=(GSkymap)");
+	test_value(total_new,  total_ref, 1.0e-3, "Test operator/(GSkymap)");
  
     // Check map scaling
     test_map   = map_src;