diff --git a/pyext/ctskymap.i b/pyext/ctskymap.i index 158f03e..8fef10f 100755 --- a/pyext/ctskymap.i +++ b/pyext/ctskymap.i @@ -49,6 +49,7 @@ public: void save(void); void publish(const std::string& name = ""); const GSkyMap& skymap(void) const; + const GFits& fits(void) const; void exclusion_map(const GSkyRegionMap& exclusion_map); GSkyRegionMap exclusion_map(void) const; }; diff --git a/src/ctskymap/ctskymap.cpp b/src/ctskymap/ctskymap.cpp index aafcd0a..c784dcb 100755 --- a/src/ctskymap/ctskymap.cpp +++ b/src/ctskymap/ctskymap.cpp @@ -240,6 +240,9 @@ void ctskymap::run(void) publish(); } + // Create the output fits object + construct_fits(); + // Return return; } @@ -262,51 +265,60 @@ void ctskymap::save(void) // Get sky map filename m_outmap = (*this)["outmap"].filename(); - // Create empty FITS file - GFits fits; + // Save FITS file to disk + m_fits.saveto(m_outmap, clobber()); - // Write sky map into FITS file - write_map(fits, m_skymap, "SKYMAP"); + } // endif: filename and map were not empty - // If background subtraction is requested then write the background map, - // the significance map, the counts map and the acceptance map to the - // FITS file - if (m_bkgsubtract != "NONE") { + // Write into logger what has been done + std::string fname = (m_outmap.is_empty()) ? "NONE" : m_outmap.url(); + if (m_skymap.is_empty()) { + fname.append(" (map is empty, no file created)"); + } + log_value(NORMAL, "Sky map file", fname); - // Write background map into FITS file - write_map(fits, m_bkgmap, "BACKGROUND"); + // Return + return; +} - // Write significance map into FITS file - write_map(fits, m_sigmap, "SIGNIFICANCE"); +/***********************************************************************//** + * @brief Construct GFits object consisting of all maps + * + * Assembles all of the skymaps into a GFits object. This object is saved + * when ctskymap::save() is called. the object can also be accessed by + * calling ctskymap::fits(). + ***************************************************************************/ +void ctskymap::construct_fits(void) +{ + // Write sky map into FITS object + write_map(m_fits, m_skymap, "SKYMAP"); - // Write counts map into FITS file - write_map(fits, m_counts, "COUNTS"); + // If background subtraction is requested then write the background map, + // the significance map, the counts map and the acceptance map to the + // FITS file + if (m_bkgsubtract != "NONE") { - // Write acceptance map into FITS file - write_map(fits, m_acceptance, "ACCEPTANCE"); + // Write background map into FITS file + write_map(m_fits, m_bkgmap, "BACKGROUND"); - // If RING background subtraction is requested then write also - // the exclusion map to the FITS file - if (m_bkgsubtract == "RING") { - write_map(fits, m_exclmap, "EXCLUSION"); - } + // Write significance map into FITS file + write_map(m_fits, m_sigmap, "SIGNIFICANCE"); - } // endif: background subtraction was requested + // Write counts map into FITS file + write_map(m_fits, m_counts, "COUNTS"); - // Save FITS file to disk - fits.saveto(m_outmap, clobber()); + // Write acceptance map into FITS file + write_map(m_fits, m_acceptance, "ACCEPTANCE"); - } // endif: filename and map were not empty + // If RING background subtraction is requested then write also + // the exclusion map to the FITS file + if (m_bkgsubtract == "RING") { + write_map(m_fits, m_exclmap, "EXCLUSION"); + } - // Write into logger what has been done - std::string fname = (m_outmap.is_empty()) ? "NONE" : m_outmap.url(); - if (m_skymap.is_empty()) { - fname.append(" (map is empty, no file created)"); - } - log_value(NORMAL, "Sky map file", fname); + } // endif: background subtraction was requested - // Return return; } @@ -372,6 +384,7 @@ void ctskymap::init_members(void) m_counts.clear(); m_acceptance.clear(); m_exclmap.clear(); + m_fits.clear(); // Initialise cache m_solidangle.clear(); @@ -415,6 +428,7 @@ void ctskymap::copy_members(const ctskymap& app) m_counts = app.m_counts; m_acceptance = app.m_acceptance; m_exclmap = app.m_exclmap; + m_fits = app.m_fits; // Copy cache m_solidangle = app.m_solidangle; diff --git a/src/ctskymap/ctskymap.hpp b/src/ctskymap/ctskymap.hpp index e63ed43..a239b23 100755 --- a/src/ctskymap/ctskymap.hpp +++ b/src/ctskymap/ctskymap.hpp @@ -61,6 +61,7 @@ public: void save(void); void publish(const std::string& name = ""); const GSkyMap& skymap(void) const; + const GFits& fits(void) const; void exclusion_map(const GSkyRegionMap& exclusion_map); GSkyRegionMap exclusion_map(void) const; @@ -87,6 +88,7 @@ protected: double& non, double& noff, double& alpha); + void construct_fits(void); void ring_bounding_box(const int& ipixel, int& ix1, int& ix2, int& iy1, int& iy2) const; GSkyMap ring_convolve(const GSkyMap& map, const double& rmin, @@ -123,6 +125,7 @@ protected: GSkyMap m_counts; //!< Counts map GSkyMap m_acceptance; //!< Acceptance map GSkyMap m_exclmap; //!< Exclusion map for RING background + GFits m_fits; //!< Output GFits object // Caching variables to prevent multiple computations of the same thing std::vector m_solidangle; //!< Cached pixel solid angles @@ -145,6 +148,13 @@ const GSkyMap& ctskymap::skymap(void) const } +inline +const GFits& ctskymap::fits(void) const +{ + return m_fits; +} + + /***********************************************************************//** * @brief Set exclusion region map * diff --git a/test/test_ctskymap.py b/test/test_ctskymap.py index cf3ee70..101fb25 100755 --- a/test/test_ctskymap.py +++ b/test/test_ctskymap.py @@ -290,6 +290,11 @@ class Test(test): self._check_result_file('ctskymap_ring_py1.fits[BACKGROUND]', 10, 10) self._check_result_file('ctskymap_ring_py1.fits[SIGNIFICANCE]', 10, 10) + # Check generated fits object + self._check_map(gammalib.GSkyMap(skymap.fits().image('SKYMAP')), nx=10, ny=10) + self._check_map(gammalib.GSkyMap(skymap.fits().image('BACKGROUND')), nx=10, ny=10) + self._check_map(gammalib.GSkyMap(skymap.fits().image('SIGNIFICANCE')), nx=10, ny=10) + # Test iterative RING background subtraction skymap = ctools.ctskymap() skymap['inobs'] = self._events