From 7ea034c7befae2fbc18209bb6b16118d77d581f2 Mon Sep 17 00:00:00 2001 From: Jürgen Knödlseder Date: Wed, 12 Apr 2017 10:39:29 +0200 Subject: [PATCH] Enhance unit test of ctprob and ctphase tools --- src/ctphase/ctphase.cpp | 37 +++++++++++++++++++++++-------------- test/test_ctbutterfly.py | 3 ++- test/test_ctphase.py | 66 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ test/test_ctprob.py | 18 ++++++++++++++++++ 4 files changed, 109 insertions(+), 15 deletions(-) diff --git a/src/ctphase/ctphase.cpp b/src/ctphase/ctphase.cpp index ac9a807..a5c2cc3 100644 --- a/src/ctphase/ctphase.cpp +++ b/src/ctphase/ctphase.cpp @@ -504,8 +504,8 @@ void ctphase::get_parameters(void) * @exception GException::invalid_value * No events found in observation. * - * Append phase information to events from a FITS file by adding a "PHASE" - * column to the FITS file. + * Computes the phase for all events in an observation based on a temporal + * phase curve model. ***************************************************************************/ void ctphase::phase_events(GCTAObservation* obs) { @@ -525,20 +525,29 @@ void ctphase::phase_events(GCTAObservation* obs) // Continue only if there are events if (events->size() > 0) { - // Print a warning into the log file if there are some precision concerns + // Get MJD information of observation GCTAEventAtom* event = (*events)[0]; double dt = std::abs(event->time() - m_phase.mjd()); - if (((m_phase.f1() * dt) > 1.0e8) || - (m_phase.f2() * dt * dt) > 1.0e8) { - log_string(NORMAL,"*** WARNING ***"); - log_string(NORMAL," Values supplied for reference MJD and/or f1 and/or f2"); - log_string(NORMAL," are large and may result in numerical precision issues."); - log_string(NORMAL," Consider using an ephemeris derived from a date closer"); - log_string(NORMAL," to the data being analyzed."); - log_value(NORMAL, " Observation ID", obs->obs_id()); - log_value(NORMAL, " Reference MJD", m_phase.mjd().mjd()); - log_value(NORMAL, " First event MJD", event->time().mjd()); - log_string(NORMAL,"*** WARNING ***"); + double f1 = m_phase.f1() * dt; + double f2 = m_phase.f2() * dt * dt; + + // Log information about MJD of the observation + log_value(NORMAL, "Reference MJD", m_phase.mjd().mjd()); + log_value(NORMAL, "First event MJD", event->time().mjd()); + log_value(NORMAL, "First event offset", gammalib::str(dt)+" s"); + log_value(NORMAL, "F1 * dt", gammalib::str(f1)+" Hz"); + log_value(NORMAL, "F2 * dt * dt", gammalib::str(f2)+" Hz"); + + // Print a warning into the log file if there are some precision + // concerns + if ((f1 > 1.0e8) || (f2 > 1.0e8)) { + std::string msg = "\nWARNING: Values supplied for reference MJD " + "and/or F1 and/or F2 are large and may\n" + " result in numerical precision issues. " + "Consider using an ephemeris derived\n" + " from a date closer to the data being " + "analyzed.\n"; + log_string(TERSE, msg); } // Compute the phase for all events diff --git a/test/test_ctbutterfly.py b/test/test_ctbutterfly.py index 0e7562a..30f9116 100755 --- a/test/test_ctbutterfly.py +++ b/test/test_ctbutterfly.py @@ -2,7 +2,7 @@ # ========================================================================== # This scripts performs unit tests for the ctbutterfly tool. # -# Copyright (C) 2014-2016 Michal Mayer +# Copyright (C) 2014-2017 Michal Mayer # # 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 @@ -144,6 +144,7 @@ class Test(test): butterfly['srcname'] = 'Crab' butterfly['caldb'] = self._caldb butterfly['irf'] = self._irf + butterfly['method'] = 'ENVELOPE' butterfly['emin'] = 0.1 butterfly['emax'] = 100.0 butterfly['fit'] = True diff --git a/test/test_ctphase.py b/test/test_ctphase.py index 8fd012c..e4764e2 100644 --- a/test/test_ctphase.py +++ b/test/test_ctphase.py @@ -191,6 +191,72 @@ class Test(test): # Check result file self._check_result_file('ctphase_py3.fits') + # Now run ctphase tool without a model file and a significant + # MJD offset + phase['inobs'] = self._events + phase['inmodel'] = 'NONE' + phase['srcname'] = 'NONE' + phase['mjd'] = 51500.0 + phase['phase'] = 0.0 + phase['f0'] = 1.0 + phase['f1'] = 0.1 + phase['f2'] = 0.0001 + phase['outobs'] = 'ctphase_py4.fits' + phase['logfile'] = 'ctphase_py4.log' + phase['chatter'] = 3 + + # Execute ctphase tool + phase.logFileOpen() # Make sure we get a log file + phase.execute() + + # Check result file + self._check_result_file('ctphase_py4.fits') + + # Test invalid source name + self.test_try('Test ctphase with invalid source name') + try: + test = ctools.ctphase() + test['inobs'] = self._events + test['inmodel'] = self._model_file + test['srcname'] = 'Venus' + test['outobs'] = 'ctphase_py5.fits' + test['logfile'] = 'ctphase_py5.log' + test.logFileOpen() + test.execute() + self.test_try_failure('Exception not thrown') + except (ValueError): + self.test_try_success() + + # Test invalid model type + self.test_try('Test ctphase with invalid model') + try: + test = ctools.ctphase() + test['inobs'] = self._events + test['inmodel'] = self._model + test['srcname'] = 'Background' + test['outobs'] = 'ctphase_py6.fits' + test['logfile'] = 'ctphase_py6.log' + test.logFileOpen() + test.execute() + self.test_try_failure('Exception not thrown') + except (ValueError): + self.test_try_success() + + # Test invalid model without temporal phase curve model + self.test_try('Test ctphase with invalid model') + try: + test = ctools.ctphase() + test['inobs'] = self._events + test['inmodel'] = self._model + test['srcname'] = 'Crab' + test['outobs'] = 'ctphase_py7.fits' + test['logfile'] = 'ctphase_py7.log' + test.logFileOpen() + test.execute() + self.test_try_failure('Exception not thrown') + except (ValueError): + self.test_try_success() + # Return return diff --git a/test/test_ctprob.py b/test/test_ctprob.py index 0ad8306..dd8d9aa 100755 --- a/test/test_ctprob.py +++ b/test/test_ctprob.py @@ -139,6 +139,7 @@ class Test(test): prob['caldb'] = self._caldb prob['irf'] = self._irf prob['edisp'] = False + prob['publish'] = True prob['logfile'] = 'ctprob_py1.log' prob['chatter'] = 2 @@ -189,6 +190,23 @@ class Test(test): # Check result file self._check_result_file('ctprob_py3_crab_events.fits') + # Copy ctprob tool + cpy_prob = prob.copy() + + # Clear phase tool + prob.clear() + + # Execute copy + cpy_prob['prefix'] = 'ctprob_py4_' + cpy_prob['outobs'] = 'ctprob_py4_obs.xml' + cpy_prob['logfile'] = 'ctprob_py4.log' + cpy_prob['chatter'] = 4 + cpy_prob.logFileOpen() # Needed to get a new log file + cpy_prob.execute() + + # Check result file + self._check_result_file('ctprob_py4_ctprob_py3_crab_events.fits') + # Return return -- libgit2 0.21.4