From 3130ffb6c0b393d0dcd9971d2e34285ee77f6146 Mon Sep 17 00:00:00 2001
From: Karl Kosack <kosack@s69dyn7.ifh.de>
Date: Fri, 31 Jan 2014 14:32:12 +0100
Subject: [PATCH] fixed bug in pointing interpolation (now gives reasonable answers in test case)

---
 include/GHorizDir.hpp             |  1 +
 inst/cta/include/GCTAPointing.hpp |  7 +++++--
 inst/cta/src/GCTAPointing.cpp     |  9 ++++++---
 inst/cta/test/test_CTA.cpp        | 45 ++++++++++++++++++++++++++++++++++++++++++---
 inst/cta/test/test_CTA.hpp        |  1 +
 5 files changed, 55 insertions(+), 8 deletions(-)

diff --git a/include/GHorizDir.hpp b/include/GHorizDir.hpp
index 8f7dbf9..693e90f 100644
--- a/include/GHorizDir.hpp
+++ b/include/GHorizDir.hpp
@@ -31,6 +31,7 @@
 #include <string>
 #include "GBase.hpp"
 #include "GVector.hpp"
+#include "GMath.hpp"
 
 /* __ Compile options ____________________________________________________ */
 #define G_SINCOS_CACHE
diff --git a/inst/cta/include/GCTAPointing.hpp b/inst/cta/include/GCTAPointing.hpp
index ebfcd41..bc1e087 100644
--- a/inst/cta/include/GCTAPointing.hpp
+++ b/inst/cta/include/GCTAPointing.hpp
@@ -33,6 +33,7 @@
 #include "GTime.hpp"
 #include "GMatrix.hpp"
 #include "GNodeArray.hpp"
+#include "GHorizDir.hpp"
 
 /***********************************************************************//**
  * @class GCTAPointing
@@ -91,8 +92,10 @@ protected:
 
     bool                m_has_table; //!< table is loaded
     GNodeArray          m_table_nodes;
-    std::vector<double> m_table_az;
-    std::vector<double> m_table_alt;
+    std::vector<double> m_table_az; //!< table of azimuths (rad)
+    std::vector<double> m_table_alt; //!< table of altitudes (rad)
+    GTime               m_table_tmin; //!< min time bound in table
+    GTime               m_table_tmax; //!<max time bound in table
 
     // Cached members
     mutable bool    m_has_cache;  //!< Has transformation cache
diff --git a/inst/cta/src/GCTAPointing.cpp b/inst/cta/src/GCTAPointing.cpp
index c3632f1..329a322 100644
--- a/inst/cta/src/GCTAPointing.cpp
+++ b/inst/cta/src/GCTAPointing.cpp
@@ -68,8 +68,11 @@ GCTAPointing::load_pointing_table(std::string filename)
     GFitsTableCol *stop  = (*table)["STOP"];
     GFitsTableCol *alt    = (*table)["ALT_PNT"];
     GFitsTableCol *az     = (*table)["AZ_PNT"];
-    //  GFitsTableCol *ra     = (*table)["RA_PNT"];
-    //GFitsTableCol *dec    = (*table)["DEC_PNT"];
+
+    // later on, may also interpolate RA/Dec, in the case of a drift
+    //  scan or other tracking mode:
+    //     GFitsTableCol *ra     = (*table)["RA_PNT"];
+    //     GFitsTableCol *dec    = (*table)["DEC_PNT"];
 
 
     // Loop over the elements and build a lookup table:
@@ -79,7 +82,7 @@ GCTAPointing::load_pointing_table(std::string filename)
 
     for (size_t ii=0; ii < table->nrows(); ii++) {
 
-      double midtime = 0.5*(stop->real(ii) - start->real(ii));
+      double midtime = 0.5*(stop->real(ii) + start->real(ii));
       GTime tt(midtime, "sec" );
       m_table_nodes.append( tt.secs() );
       m_table_az[ii] = az->real(ii) * gammalib::deg2rad;
diff --git a/inst/cta/test/test_CTA.cpp b/inst/cta/test/test_CTA.cpp
index b24b47b..8413a4a 100644
--- a/inst/cta/test/test_CTA.cpp
+++ b/inst/cta/test/test_CTA.cpp
@@ -187,7 +187,11 @@ void TestGCTAPointing::set(void)
     name("GCTAPointing");
 
     // Append tests to test suite
-    append(static_cast<pfunction>(&TestGCTAPointing::test_load_table), "Test load pointing from table");
+    append(static_cast<pfunction>(&TestGCTAPointing::test_load_table), 
+           "Test load pointing from table");
+
+    append(static_cast<pfunction>(&TestGCTAPointing::test_interpolate_altaz),
+           "Test alt/az interpolation given a time");
 
     return;
 }
@@ -206,9 +210,8 @@ TestGCTAPointing* TestGCTAPointing::clone(void) const
 
 
 /***********************************************************************//**
- * @brief Test CTA Npred computation
+ * @brief Test ability to load a CTA pointing table
  *
-
  ***************************************************************************/
 void TestGCTAPointing::test_load_table(void)
 {
@@ -221,6 +224,42 @@ void TestGCTAPointing::test_load_table(void)
 
 
 /***********************************************************************//**
+ * @brief Test interpolation of alt/az pointing dir as a function of
+ * time
+ *
+ ***************************************************************************/
+void TestGCTAPointing::test_interpolate_altaz(void)
+{
+
+    GCTAObservation run;
+
+    GCTAPointing pnt;
+    pnt.load_pointing_table( cta_point_table );
+
+    // test an out-of bounds time
+    
+
+    try {
+      GTime time(155470378.7, "sec"); 
+      GHorizDir dir = pnt.dir_horiz( time );
+    }
+    catch (GException::out_of_range &exc) {
+      std::cout << exc.what() << std::endl;
+    }
+
+    // Get a time that is somewhere in the run:
+
+    GTime time(155470378.7, "sec"); 
+    GHorizDir dir = pnt.dir_horiz( time );
+    
+    std::cout << "time="<<time<<" dir="<<dir<<std::endl;
+
+    return;
+}
+
+
+
+/***********************************************************************//**
  * @brief Test CTA Aeff computation
  ***************************************************************************/
 void TestGCTAResponse::test_response_aeff(void)
diff --git a/inst/cta/test/test_CTA.hpp b/inst/cta/test/test_CTA.hpp
index 3e39340..61e5ea3 100644
--- a/inst/cta/test/test_CTA.hpp
+++ b/inst/cta/test/test_CTA.hpp
@@ -137,6 +137,7 @@ public:
     virtual void              set(void);
     virtual TestGCTAPointing* clone(void) const;
     void                      test_load_table(void);
+    void                      test_interpolate_altaz(void);
 };
 
 
--
libgit2 0.21.4