Take into account RoI in ctskymap

Jürgen Knödlseder
1 parent 052e74fa
Showing 4 changed files with 27 additions and 5 deletions
ChangeLog
NEWS
src/ctskymap/ctskymap.cpp
src/ctskymap/ctskymap.hpp
-2017-03-13
+2017-03-18
         * Version 1.3.0 released
         ========================
+        Take into account RoI in ctskymap
         Set local language to english (#1370)
         Enhance logfile output of csroot2caldb script
         Set default deadc value to 0.98 (CTA has <2% deadtime)
 New Features and Important Changes in ctools 1.3.0.dev1
-13 March 2017
+18 March 2017
 Introduction
@@ -52,7 +52,8 @@ None
 ctskymap - CTA sky mapping
 --------------------------
-None
+Any RoI is now taken into account in ctskymap so that if IRF background
+subtraction is requested the counts and background maps are coherent.
 ctmodel - CTA counts cube modelling
@@ -489,10 +489,14 @@ void ctskymap::map_events(GCTAObservation* obs)
     emax.TeV(m_emax);
     // Initialise binning statistics
+    int num_outside_roi    = 0;
     int num_outside_map    = 0;
     int num_outside_erange = 0;
     int num_in_map         = 0;
+    // Extract region of interest from observation
+    GCTARoi roi = obs->roi();
+
     // Fill sky map
     for (int i = 0; i < events->size(); ++i) {
@@ -517,6 +521,12 @@ void ctskymap::map_events(GCTAObservation* obs)
             continue;
         }
+        // If RoI is valid then skip if  instrument direction is not within RoI
+        if (roi.is_valid() && !roi.contains(*inst)) {
+            continue;
+            num_outside_roi++;
+        }
+
         // Fill event in skymap
         m_skymap(pixel, 0) += 1.0;
         num_in_map++;
@@ -526,6 +536,7 @@ void ctskymap::map_events(GCTAObservation* obs)
     // Log binning results
     log_value(NORMAL, "Events in list", obs->events()->size());
     log_value(NORMAL, "Events in map", num_in_map);
+    log_value(NORMAL, "Events outside RoI", num_outside_roi);
     log_value(NORMAL, "Events outside map area", num_outside_map);
     log_value(NORMAL, "Events outside energies", num_outside_erange);
@@ -614,6 +625,9 @@ void ctskymap::map_background_irf(GCTAObservation* obs)
         nebin = 10;
     }
+    // Extract region of interest from observation
+    GCTARoi roi = obs->roi();
+
     // Setup energy nodes for computation
     GEnergies energies(nebin, emin, emax); // More complicated nbins later
@@ -626,6 +640,12 @@ void ctskymap::map_background_irf(GCTAObservation* obs)
         // Convert sky direction in instrument direction
         GCTAInstDir instdir = obs->pointing().instdir(skydir);
+        // If RoI is valid and instrument direction is not within RoI then
+        // skip pixel
+        if (roi.is_valid() && !roi.contains(instdir)) {
+            continue;
+        }
+
         // Compute background value for each energy node
         std::vector<double> values(energies.size(), 0.0);
         for (int k = 0; k < energies.size(); ++k) {
 /***************************************************************************
  *                       ctskymap - Sky mapping tool                       *
  * ----------------------------------------------------------------------- *
- *  copyright (C) 2010-2016 by Juergen Knoedlseder                         *
+ *  copyright (C) 2010-2017 by Juergen Knoedlseder                         *
  * ----------------------------------------------------------------------- *
  *                                                                         *
  *  This program is free software: you can redistribute it and/or modify   *
@@ -33,7 +33,7 @@
 /* __Definitions _________________________________________________________ */
 #define CTSKYMAP_NAME    "ctskymap"
-#define CTSKYMAP_VERSION "1.2.0"
+#define CTSKYMAP_VERSION "1.3.0"
 /***********************************************************************//**