From bb9410c7c406e17d9246fc02016e65b46518d4a8 Mon Sep 17 00:00:00 2001 From: Josh Cardenzana Date: Wed, 5 Apr 2017 11:29:58 +0200 Subject: [PATCH] Updated test cases for ctphase (#1982) --- pyext/Makefile.am | 1 + pyext/ctools/tools.i | 1 + pyext/ctphase.i | 61 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ src/ctphase/ctphase.cpp | 179 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++-------------------------------------------------------------------------------------------------- src/ctphase/ctphase.hpp | 8 +++----- src/ctphase/ctphase.par | 10 +++++----- src/ctphase/main.cpp | 4 ++-- test/test_ctphase.py | 223 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ test/test_python_ctools.py | 2 ++ 9 files changed, 379 insertions(+), 110 deletions(-) create mode 100644 pyext/ctphase.i create mode 100644 test/test_ctphase.py diff --git a/pyext/Makefile.am b/pyext/Makefile.am index 1326dab..e38a57a 100755 --- a/pyext/Makefile.am +++ b/pyext/Makefile.am @@ -44,6 +44,7 @@ testscripts = ctools/tests/__init__.py \ $(top_srcdir)/test/test_python_ctools.py \ $(top_srcdir)/test/test_ctobssim.py \ $(top_srcdir)/test/test_ctselect.py \ + $(top_srcdir)/test/test_ctphase.py \ $(top_srcdir)/test/test_ctbin.py \ $(top_srcdir)/test/test_ctlike.py \ $(top_srcdir)/test/test_cttsmap.py \ diff --git a/pyext/ctools/tools.i b/pyext/ctools/tools.i index d5350be..e15ffba 100755 --- a/pyext/ctools/tools.i +++ b/pyext/ctools/tools.i @@ -61,6 +61,7 @@ namespace std { %include "ctlike.i" %include "ctmodel.i" %include "ctselect.i" +%include "ctphase.i" %include "ctskymap.i" %include "ctexpcube.i" %include "ctpsfcube.i" diff --git a/pyext/ctphase.i b/pyext/ctphase.i new file mode 100644 index 0000000..13e681b --- /dev/null +++ b/pyext/ctphase.i @@ -0,0 +1,61 @@ +/*************************************************************************** + * ctphase - Append phase information to CTA events file * + * ----------------------------------------------------------------------- * + * copyright (C) 2017 by Leonardo Di Venere * + * ----------------------------------------------------------------------- * + * * + * 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 . * + * * + ***************************************************************************/ +/** + * @file ctphase.i + * @brief Append phase information to CTA events file + * @author Leonardo Di Venere + */ +%{ +/* Put headers and other declarations here that are needed for compilation */ +#include "ctphase.hpp" +%} + + +/***********************************************************************//** + * @class ctphase + * + * @brief Data selection tool + ***************************************************************************/ +class ctphase : public ctobservation { +public: + // Constructors and destructors + ctphase(void); + explicit ctphase(const GObservations& obs); + ctphase(int argc, char *argv[]); + ctphase(const ctphase& app); + virtual ~ctphase(void); + + // Methods + void clear(void); + void run(void); + void save(void); + void publish(const std::string& name = ""); +}; + + +/***********************************************************************//** + * @brief Data selection tool Python extension + ***************************************************************************/ +%extend ctphase { + ctphase copy() { + return (*self); + } +} diff --git a/src/ctphase/ctphase.cpp b/src/ctphase/ctphase.cpp index e71641e..34e1bd2 100644 --- a/src/ctphase/ctphase.cpp +++ b/src/ctphase/ctphase.cpp @@ -1,7 +1,7 @@ /*************************************************************************** - * ctphase - Data selection tool * + * ctphase - Append phase information to CTA events file * * ----------------------------------------------------------------------- * - * copyright (C) 2010-2017 by Juergen Knoedlseder * + * copyright (C) 2017 by Leonardo Di Venere * * ----------------------------------------------------------------------- * * * * This program is free software: you can redistribute it and/or modify * @@ -20,8 +20,8 @@ ***************************************************************************/ /** * @file ctphase.cpp - * @brief Data selection tool definition - * @author Juergen Knoedlseder + * @brief Append phase information to CTA events file + * @author Leonardo Di Venere */ /* __ Includes ___________________________________________________________ */ @@ -34,7 +34,7 @@ /* __ Method name definitions ____________________________________________ */ #define G_RUN "ctphase::run()" -#define G_SELECT_EVENTS "ctphase::select_events(GCTAObservation*, "\ +#define G_PHASE_EVENTS "ctphase::phase_events(GCTAObservation*, "\ "std::string&, std::string&, std::string&)" #define G_SET_EBOUNDS "ctphase::set_ebounds(GCTAObservation*, GEbounds&)" @@ -203,14 +203,10 @@ void ctphase::clear(void) /***********************************************************************//** - * @brief Select event data + * @brief Append phase information to event data * * This method reads in the application parameters and loops over all - * observations that were found to perform an event selection. Event - * selection is done by writing each observation to a temporary file and - * re-opening the temporary file using the cfitsio event filter syntax. - * The temporary file is deleted after this action so that no disk overflow - * will occur. + * observations that were found to append phase information to the data. ***************************************************************************/ void ctphase::run(void) { @@ -226,7 +222,7 @@ void ctphase::run(void) log_observations(NORMAL, m_obs, "Input observation"); // Write header into logger - log_header1(TERSE, "Event selection"); + log_header1(TERSE, "Appending phase informaiton"); // Initialise counters int n_observations = 0; @@ -283,7 +279,7 @@ void ctphase::run(void) // Fall through in case that the event file is empty if (obs->events()->size() == 0) { log_string(NORMAL, " Warning: No events in event file \""+ - m_infiles[i]+"\". Event selection skipped."); + m_infiles[i]+"\". Event phasing skipped."); continue; } @@ -326,7 +322,7 @@ void ctphase::run(void) } } - // Load observation from temporary file, including event selection + // Load observation from temporary file and compute phase phase_events(obs, filename, m_evtname[i], m_gtiname[i]); // Close temporary file @@ -359,22 +355,22 @@ void ctphase::run(void) /***********************************************************************//** - * @brief Save the selected event list(s) + * @brief Save the output event list(s) * - * This method saves the selected event list(s) into FITS file(s). There are + * This method saves the updated event list(s) into FITS file(s). There are * two modes, depending on the m_use_xml flag. * - * If m_use_xml is true, all selected event list(s) will be saved into FITS + * If m_use_xml is true, all output event list(s) will be saved into FITS * files, where the output filenames are constructued from the input * filenames by prepending the m_prefix string to name. Any path information * will be stripped form the input name, hence event files will be written * into the local working directory (unless some path information is present * in the prefix). In addition, an XML file will be created that gathers - * the filename information for the selected event list(s). If an XML file + * the filename information for the output event list(s). If an XML file * was present on input, all metadata information will be copied from this * input file. * - * If m_use_xml is false, the selected event list will be saved into a FITS + * If m_use_xml is false, the output event list will be saved into a FITS * file. ***************************************************************************/ void ctphase::save(void) @@ -506,6 +502,8 @@ void ctphase::copy_members(const ctphase& app) ***************************************************************************/ void ctphase::free_members(void) { + m_phase.clear(); + // Return return; } @@ -528,19 +526,70 @@ void ctphase::get_parameters(void) // information and do not accept counts cubes. setup_observations(m_obs, false, true, false); - // Check for sanity of energy selection parameters + // Check for sanity of frequency and phase information if ((*this)["p0"].is_valid() && (*this)["mjd"].is_valid() && (*this)["f0"].is_valid()) { - m_phase = GModelTemporalPhaseCurve(); - m_phase.norm(1.0); - m_phase.mjd((*this)["mjd"].real()); - m_phase.phase((*this)["p0"].real()); - m_phase.f0((*this)["f0"].real()); - m_phase.f1((*this)["f1"].real()); - m_phase.f2((*this)["f2"].real()); + // The following creates an XML object responsible for initializing + // 'm_phase'. This was done as the default values for the parameters in + // 'GModelTemporalPhaseCurve' are such that the ranges are predefined. + // Initializing 'm_phase' from an XML object circumvents this restriction. + + GXmlElement norm; + norm.name("parameter"); + norm.attribute("name", "Normalization"); + norm.attribute("value", "1"); + norm.attribute("min", "0"); + norm.attribute("max", "1000"); + + GXmlElement mjd; + mjd.name("parameter"); + mjd.attribute("name", "MJD"); + mjd.attribute("value", (*this)["mjd"].value()); + mjd.attribute("min", (*this)["mjd"].min()); + mjd.attribute("max", (*this)["mjd"].max()); + + GXmlElement p0; + p0.name("parameter"); + p0.attribute("name", "Phase"); + p0.attribute("value", (*this)["p0"].value()); + p0.attribute("min", (*this)["p0"].min()); + p0.attribute("max", (*this)["p0"].max()); + + GXmlElement f0; + f0.name("parameter"); + f0.attribute("name", "F0"); + f0.attribute("value", (*this)["f0"].value()); + f0.attribute("min", (*this)["f0"].min()); + f0.attribute("max", (*this)["f0"].max()); + + GXmlElement f1; + f1.name("parameter"); + f1.attribute("name", "F1"); + f1.attribute("value", (*this)["f1"].value()); + f1.attribute("min", (*this)["f1"].min()); + f1.attribute("max", (*this)["f1"].max()); + + GXmlElement f2; + f2.name("parameter"); + f2.attribute("name", "F2"); + f2.attribute("value", (*this)["f2"].value()); + f2.attribute("min", (*this)["f2"].min()); + f2.attribute("max", (*this)["f2"].max()); + GXmlElement info; + info.append(norm); + info.append(mjd); + info.append(p0); + info.append(f0); + info.append(f1); + info.append(f2); + + info.attribute("file",""); + + m_phase = GModelTemporalPhaseCurve(); + m_phase.read(info); } else { GException::invalid_value("ctphase::get_parameters()", @@ -576,12 +625,8 @@ void ctphase::get_parameters(void) * @exception GException::invalid_value * No events extension found in FITS file. * - * Select events from a FITS file by making use of the selection possibility - * of the cfitsio library on loading a file. A selection string is created - * from the specified criteria that is appended to the filename so that - * cfitsio will automatically filter the event data. This selection string - * is then applied when opening the FITS file. The event list in the current - * observation is replaced by selected event list read from the FITS file. + * Append phase information to events from a FITS file by adding a "PHASE" + * column to the FITS file. ***************************************************************************/ void ctphase::phase_events(GCTAObservation* obs, const std::string& filename, @@ -591,13 +636,6 @@ void ctphase::phase_events(GCTAObservation* obs, // Write header into logger log_header3(NORMAL, "Events phasing"); - // Initialise selection and addition strings - std::string selection; - std::string add; - - // Initialise selection flags - bool remove_all = false; - // Get CTA event list pointer GCTAEventList* list = static_cast(const_cast(obs->events())); @@ -617,62 +655,7 @@ void ctphase::phase_events(GCTAObservation* obs, if (!list->has_phase()) { list->has_phase(true); } -/* - // Save GTI for later usage - GGti gti = list->gti(); - - // Build input filename including selection expression - std::string expression = filename; - // Dump FITS filename including selection expression - log_value(NORMAL, "FITS filename", expression); - - // Open FITS file - GFits file(expression); - - // Log selected FITS file - log_header3(EXPLICIT, "FITS file content after selection"); - log_string(EXPLICIT, file.print(m_chatter)); - - // Check if we have an events HDU - if (!file.contains(evtname)) { - std::string msg = "No events extension \""+evtname+"\" found in " - "FITS file. The expression \""+expression+"\" " - "was used to open the FITS file."; - throw GException::invalid_value(G_SELECT_EVENTS, msg); - } - - // Determine number of events in the events HDU. If removal of all events - // has been requested then set the number of events to zero. - int nevents = (remove_all) ? 0 : file.table(evtname)->nrows(); - - // If the selected event list is empty then append an empty event list - // to the observation - if (nevents < 1) { - - // Create empty event list - GCTAEventList eventlist; - - // Append list to observation - obs->events(eventlist); - - } - - // ... otherwise load the data from the temporary file - else { - obs->read(file); - } - - // Get CTA event list pointer - list = static_cast(const_cast(obs->events())); - - // Make sure that events are fetched since the temporary file will be - // closed later - list->fetch(); - - // Set event list GTI (in any case as any event list has a GTI) - list->gti(gti); -*/ // Return return; } @@ -687,7 +670,7 @@ void ctphase::phase_events(GCTAObservation* obs, * This method checks if the input FITS file is correct. ***************************************************************************/ std::string ctphase::check_infile(const std::string& filename, - const std::string& evtname) const + const std::string& evtname) const { // Initialise message string std::string message = ""; @@ -805,7 +788,7 @@ std::string ctphase::get_gtiname(const std::string& filename, // Get GTI name gtiname = events.gtiname(); - } + } // endif: filename was not empty // Return GTI name return (gtiname); @@ -957,7 +940,7 @@ void ctphase::save_xml(void) * * myfits.fits[EVENTS1;GTI1] * - * will write the selected events into the "EVENTS1" extension and the + * will write the updated events into the "EVENTS1" extension and the * Good Time Intervals into the "GTI1" extension of the "myfits.fits" FITS * file. If the Good Time Intervals extension name is skipped, e.g. * diff --git a/src/ctphase/ctphase.hpp b/src/ctphase/ctphase.hpp index dac9f8a..0a58c5b 100644 --- a/src/ctphase/ctphase.hpp +++ b/src/ctphase/ctphase.hpp @@ -1,5 +1,5 @@ /*************************************************************************** - * ctphase - Data selection tool * + * ctphase - Append phase information to CTA events file * * ----------------------------------------------------------------------- * * copyright (C) 2017 by Leonardo Di Venere * * ----------------------------------------------------------------------- * @@ -20,8 +20,8 @@ ***************************************************************************/ /** * @file ctphase.hpp - * @brief Data selection tool definition - * @author Juergen Knoedlseder + * @brief Append phase information to CTA events file + * @author Leonardo Di Venere */ #ifndef CTPHASE_HPP @@ -63,8 +63,6 @@ public: void save(void); void publish(const std::string& name = ""); - // TODO: Add Methods for updating model parameters - protected: // Protected methods void init_members(void); diff --git a/src/ctphase/ctphase.par b/src/ctphase/ctphase.par index 9ee197b..c3bfd80 100644 --- a/src/ctphase/ctphase.par +++ b/src/ctphase/ctphase.par @@ -22,12 +22,12 @@ outobs, f, a, "phased_events.fits",,, "Output event list or observation definiti prefix, s, h, "phased_",,, "Prefix for event lists in observation definition XML file" # -# Selection parameters +# Phase information parameters #===================== mjd, r, a, 52400,,, "Reference MJD for phase information" -f0, r, a, 1,,, "Frequency at reference MJD" -f1, r, a, 0,,, "First freq. derivative at reference MJD" -f2, r, a, 0,,, "Second freq. derivative at reference MJD" +f0, r, a, 1,,, "Frequency at reference MJD [Hz]" +f1, r, a, -0.2,,, "First freq. derivative at reference MJD [s^-1]" +f2, r, a, 0,,, "Second freq. derivative at reference MJD [s^-2]" p0, r, a, 0.1,0.0,1.0, "Phase at reference MJD" # @@ -37,5 +37,5 @@ publish, b, h, no,,, "Publish event list on VO Hub?" chatter, i, h, 2,0,4, "Chattiness of output" clobber, b, h, yes,,, "Overwrite existing output files with new output files?" debug, b, h, no,,, "Debugging mode activated" -mode, s, h, "h",,, "Mode of automatic parameters" +mode, s, h, "ql",,, "Mode of automatic parameters" logfile, f, h, "ctphase.log",,, "Log filename" diff --git a/src/ctphase/main.cpp b/src/ctphase/main.cpp index 3afb1ab..c0761b6 100644 --- a/src/ctphase/main.cpp +++ b/src/ctphase/main.cpp @@ -1,7 +1,7 @@ /*************************************************************************** - * ctphase - Append phase information to CTA events file * + * ctphase - Append phase information to CTA events file * * ----------------------------------------------------------------------- * - * copyright (C) 2017 by Leonardo Di Venere * + * copyright (C) 2017 by Leonardo Di Venere * * ----------------------------------------------------------------------- * * * * This program is free software: you can redistribute it and/or modify * diff --git a/test/test_ctphase.py b/test/test_ctphase.py new file mode 100644 index 0000000..d02646c --- /dev/null +++ b/test/test_ctphase.py @@ -0,0 +1,223 @@ +#! /usr/bin/env python +# ========================================================================== +# This scripts performs unit tests for the ctselect tool. +# +# Copyright (C) 2014-2016 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 . +# +# ========================================================================== +import os +import gammalib +import ctools +from testing import test + + +# ============================ # +# Test class for ctselect tool # +# ============================ # +class Test(test): + """ + Test class for ctselect tool + """ + # Constructor + def __init__(self): + """ + Constructor + """ + # Call base class constructor + test.__init__(self) + + # Set test data + self._invalid_events = self._datadir + '/invalid_event_list.fits' + + # Return + return + + # Set test functions + def set(self): + """ + Set all test functions + """ + # Set test name + self.name('ctphase') + + # Append tests + self.append(self._test_cmd, 'Test ctphase on command line') + self.append(self._test_python, 'Test ctphase from Python') + + # Return + return + + # Test ctphase on command line + def _test_cmd(self): + """ + Test ctphase on the command line + """ + # Set tool name + ctphase = self._tool('ctphase') + + # Setup ctphase command + cmd = ctphase+' inobs="'+self._events+'"'+ \ + ' outobs="ctphase_cmd1.fits"'+ \ + ' p0=0.1 mjd=52400'+ \ + ' f0=0.1 f1=0. f2=0.'+ \ + ' logfile="ctphase_cmd1.log" chatter=1' + + # Check if execution of wrong command fails + self.test_assert(self._execute('command_that_does_not_exist') != 0, + 'Self test of test script') + + # Check if execution was successful + self.test_assert(self._execute(cmd) == 0, + 'Check successful execution from command line') + + # Check result file + self._check_result_file('ctphase_cmd1.fits') + + # Setup ctphase command + cmd = ctphase+' inobs="event_file_that_does_not_exist.fits"'+ \ + ' outobs="ctphase_cmd2.fits"'+ \ + ' p0=0.1 mjd=52400'+ \ + ' f0=0.1 f1=0. f2=0.'+ \ + ' logfile="ctphase_cmd2.log" chatter=1' + + # Check if execution failed + self.test_assert(self._execute(cmd) != 0, + 'Check invalid input file when executed from command line') + + # Setup ctphase --help option + cmd = ctphase+' --help' + + # Check if execution was successful in case that the CTOOLS + # environment variable was set or failed otherwise + if 'CTOOLS' in os.environ: + self.test_value(self._execute(cmd), 0, + 'Check successful execution with --help option') + else: + self.test_assert(self._execute(cmd) != 0, + 'Check execution failure with --help option') + + # Return + return + + # Test ctphase from Python + def _test_python(self): + """ + Test ctphase from Python + """ + # Allocate empty ctphase tool + phase = ctools.ctphase() + + # Check that empty ctphase tool holds an empty observation + self._check_obs(phase.obs(), nobs=0) + + # Check that saving does nothing + phase['outobs'] = 'ctphase_py0.fits' + phase['logfile'] = 'ctphase_py0.log' + phase.logFileOpen() + phase.save() + self.test_assert(not os.path.isfile('ctphase_py0.fits'), + 'Check that no event list has been created') + + # Check that clearing does not lead to an exception or segfault + phase.clear() + + # Now set ctphase parameters + phase['inobs'] = self._events + phase['p0'] = 0. + phase['mjd'] = 52400 + phase['f0'] = 1. + phase['f1'] = -0.2 + phase['f2'] = 0. + phase['outobs'] = 'ctphase_py1.fits' + phase['logfile'] = 'ctphase_py1.log' + phase['chatter'] = 2 + + # Run ctphase tool + phase.logFileOpen() # Make sure we get a log file + phase.run() + phase.save() + + # Check result file + self._check_result_file('ctphase_py1.fits', nevents=6141) + + # Copy ctphase tool + cpy_phase = phase.copy() + + # Execute copy of ctphase tool again, now with a negative value + # of initial phase to test for invalid parameter value + try: + self.test_try('Test invalid parameter value') + cpy_phase['p0'] = -0.5 + cpy_phase['mjd'] = 52400 + cpy_phase['f0'] = 1. + cpy_phase['f1'] = -0.2 + cpy_phase['f2'] = 0. + cpy_phase['outobs'] = 'ctphase_py2.fits' + cpy_phase['logfile'] = 'ctphase_py2.log' + cpy_phase['chatter'] = 3 + cpy_phase['publish'] = True + cpy_phase.logFileOpen() # Needed to get a new log file + cpy_phase.execute() + self.test_try_failure('Exception not thrown') + except ValueError: + self.test_try_success() + + return + + # Check result file + def _check_result_file(self, filename, nevents=591, dec=22.01, rad=3.0): + """ + Check result file + + Parameters + ---------- + filename : str + Event list file name + nevents : int, optional + Expected number of events + dec : float, optional + Expected Declination (deg) + rad : float, optional + Expected radius (deg) + """ + # Open result file + events = gammalib.GCTAEventList(filename) + phases = [ events[i].phase() for i in range(events.size())] + test_result = int(sum(phases)>0) + self.test_value(test_result, 1, 'Check number of events') + + # Return + return + + # Check observation and event list + def _check_obs(self, obs, nobs=1): + """ + Check observation and event list + + Parameters + ---------- + obs : `~gammalib.GObservations` + Models + nobs : int, optional + Expected number of observations + nevents : int, optional + Expected number of events + """ + # Check number of observations + self.test_value(obs.size(), nobs, 'Check number of observations') + + # Return + return diff --git a/test/test_python_ctools.py b/test/test_python_ctools.py index 9f546da..0ba05e0 100755 --- a/test/test_python_ctools.py +++ b/test/test_python_ctools.py @@ -25,6 +25,7 @@ import os sys.path.append(os.path.dirname(os.path.abspath(__file__))) import test_ctobssim import test_ctselect +import test_ctphase import test_ctbin import test_ctlike import test_cttsmap @@ -102,6 +103,7 @@ def test(installed=False): # Define list of test suites tests = [test_ctobssim.Test(), test_ctselect.Test(), + test_ctphase.Test(), test_ctbin.Test(), test_ctlike.Test(), test_cttsmap.Test(), -- libgit2 0.21.4