diff --git a/NEWS b/NEWS index 94fd687..77a62f6 100644 --- a/NEWS +++ b/NEWS @@ -23,7 +23,8 @@ The following classes have been renamed: The following methods have been added: - GSkyRegion::overlaps(GSkyRegions&) - GSkyMap::overlaps(GSkyRegion&) -- GModelTemporalPhaseCurve::normalise() +- GModelTemporalPhaseCurve::normalize() +- GModelTemporalPhaseCurve::value(double&) The following methods have been removed: - n.a. diff --git a/include/GModelTemporalPhaseCurve.hpp b/include/GModelTemporalPhaseCurve.hpp index c8bdfa8..7f28230 100644 --- a/include/GModelTemporalPhaseCurve.hpp +++ b/include/GModelTemporalPhaseCurve.hpp @@ -116,6 +116,7 @@ public: double phase(void) const; double phase(const GTime& time) const; void phase(const double& phase); + double value(const double& phase) const; double f0(void) const; void f0(const double& f0); double f1(void) const; @@ -130,7 +131,7 @@ protected: void copy_members(const GModelTemporalPhaseCurve& model); void free_members(void); void load_nodes(const GFilename& filename); - void normalise_nodes(void); + void normalize_nodes(void); // Protected parameter members GModelPar m_norm; //!< Normalization factor diff --git a/pyext/GModelTemporalPhaseCurve.i b/pyext/GModelTemporalPhaseCurve.i index d8317ec..c7caf91 100644 --- a/pyext/GModelTemporalPhaseCurve.i +++ b/pyext/GModelTemporalPhaseCurve.i @@ -73,6 +73,7 @@ public: double phase(void) const; double phase(const GTime& time) const; void phase(const double& phase); + double value(const double& phase) const; double f0(void) const; void f0(const double& f0); double f1(void) const; diff --git a/src/model/GModelTemporalPhaseCurve.cpp b/src/model/GModelTemporalPhaseCurve.cpp index 90d0e00..77990aa 100644 --- a/src/model/GModelTemporalPhaseCurve.cpp +++ b/src/model/GModelTemporalPhaseCurve.cpp @@ -47,7 +47,7 @@ const GModelTemporalRegistry g_temporal_phase_registry(&g_temporal_phase_seed) #define G_READ "GModelTemporalPhaseCurve::read(GXmlElement&)" #define G_WRITE "GModelTemporalPhaseCurve::write(GXmlElement&)" #define G_LOAD_NODES "GModelTemporalPhaseCurve::load_nodes(GFilename&)" -#define G_NORMALISE_NODES "GModelTemporalPhaseCurve::normalise_nodes()" +#define G_NORMALISE_NODES "GModelTemporalPhaseCurve::normalize_nodes()" /* __ Macros _____________________________________________________________ */ @@ -522,6 +522,34 @@ double GModelTemporalPhaseCurve::phase(const GTime& time) const /***********************************************************************//** + * @brief Evaluate phase curve value for a given phase + * + * @param[in] phase Phase. + * @return Value of phase curve. + * + * Computes + * + * \f[ + * S_{\rm t}(\Phi) = r(\Phi) \times {\tt m\_norm} + * \f] + * + * where + * + * \f$r(\Phi)\f$ is the phase dependent rate, defined by linear interpolation + * between the nodes in a FITS file, \f$\Phi\f$ is the phase, and + * \f${\tt m\_norm}\f$ is a normalisation constant. + ***************************************************************************/ +double GModelTemporalPhaseCurve::value(const double& phase) const +{ + // Compute function value + double value = m_norm.value() * m_nodes.interpolate(phase, m_values); + + // Return + return value; +} + + +/***********************************************************************//** * @brief Print phase curve information * * @param[in] chatter Chattiness. @@ -833,9 +861,9 @@ void GModelTemporalPhaseCurve::load_nodes(const GFilename& filename) // Close FITS file fits.close(); - // Optionally normalise nodes + // Optionally normalize nodes if (m_normalize) { - normalise_nodes(); + normalize_nodes(); } // Return @@ -852,7 +880,7 @@ void GModelTemporalPhaseCurve::load_nodes(const GFilename& filename) * Normalise the node values so that the integral over the phase interval * [0,1] gives an average normalisation of 1. ***************************************************************************/ -void GModelTemporalPhaseCurve::normalise_nodes(void) +void GModelTemporalPhaseCurve::normalize_nodes(void) { // Initialise node integral double sum = 0.0; diff --git a/test/test_GModel.cpp b/test/test_GModel.cpp index 66f639f..462a93c 100644 --- a/test/test_GModel.cpp +++ b/test/test_GModel.cpp @@ -2322,52 +2322,78 @@ void TestGModel::test_temp_phasecurve(void) test_value(model2.f1(), 0.1); test_value(model2.f2(), 0.01); test_value(model2.norm(), 2.0); - - // Test XML constructor - GXml xml(m_xml_model_point_temp_phasecurve); - GXmlElement* element = xml.element(0)->element(0)->element("temporal", 0); - GModelTemporalPhaseCurve model3(*element); - test_value(model3.size(), 6); - test_value(model3.type(), "PhaseCurve", "Check model type"); + test_value(model2.value(0.1), 0.0); + test_value(model2.value(0.2), 9.09090909090909); + test_value(model2.value(0.6), 8.18181818181818); + test_assert(model2.normalize(), "Check that phase curve is normalised"); + + // Test value constructor without phase curve normalization + GModelTemporalPhaseCurve model3(m_temp_phasecurve, + GTime(3.0, "secs"), + 0.0, + 1.0, + 0.1, + 0.01, + 2.0, + false); test_value(model3.filename().url(), m_temp_phasecurve, "Check phase curve file name"); - test_value(model3.norm(), 1.0); - test_value(model3.mjd().mjd(), 51544.5); + test_value(model3.mjd().secs(), 3.0); test_value(model3.phase(), 0.0); test_value(model3.f0(), 1.0); test_value(model3.f1(), 0.1); test_value(model3.f2(), 0.01); + test_value(model3.norm(), 2.0); + test_value(model3.value(0.1), 0.0); + test_value(model3.value(0.2), 2.0); + test_value(model3.value(0.6), 1.8); + test_assert(!model3.normalize(), "Check that phase curve is not normalised"); + + // Test XML constructor + GXml xml(m_xml_model_point_temp_phasecurve); + GXmlElement* element = xml.element(0)->element(0)->element("temporal", 0); + GModelTemporalPhaseCurve model4(*element); + test_value(model4.size(), 6); + test_value(model4.type(), "PhaseCurve", "Check model type"); + test_value(model4.filename().url(), m_temp_phasecurve, + "Check phase curve file name"); + test_value(model4.norm(), 1.0); + test_value(model4.mjd().mjd(), 51544.5); + test_value(model4.phase(), 0.0); + test_value(model4.f0(), 1.0); + test_value(model4.f1(), 0.1); + test_value(model4.f2(), 0.01); // Test filename method - model3.filename(m_temp_phasecurve); - test_value(model3.filename().url(), m_temp_phasecurve, + model4.filename(m_temp_phasecurve); + test_value(model4.filename().url(), m_temp_phasecurve, "Check phase curve file name"); // Test parameter methods - model3.norm(3.0); - test_value(model3.norm(), 3.0); - model3.mjd(GTime(3.0, "s")); - test_value(model3.mjd().secs(), 3.0); - model3.phase(0.3); - test_value(model3.phase(), 0.3); - model3.f0(3.0); - test_value(model3.f0(), 3.0); - model3.f1(3.0); - test_value(model3.f1(), 3.0); - model3.f2(3.0); - test_value(model3.f2(), 3.0); + model4.norm(3.0); + test_value(model4.norm(), 3.0); + model4.mjd(GTime(3.0, "s")); + test_value(model4.mjd().secs(), 3.0); + model4.phase(0.3); + test_value(model4.phase(), 0.3); + model4.f0(3.0); + test_value(model4.f0(), 3.0); + model4.f1(3.0); + test_value(model4.f1(), 3.0); + model4.f2(3.0); + test_value(model4.f2(), 3.0); // Test operator access const char* strarray[] = {"Normalization", "MJD", "Phase", "F0", "F1", "F2"}; for (int i = 0; i < 6; ++i) { std::string keyname(strarray[i]); - model3[keyname].remove_range(); // To allow setting of any value - model3[keyname].value(2.1); - model3[keyname].error(1.9); - model3[keyname].gradient(0.8); - test_value(model3[keyname].value(), 2.1); - test_value(model3[keyname].error(), 1.9); - test_value(model3[keyname].gradient(), 0.8); + model4[keyname].remove_range(); // To allow setting of any value + model4[keyname].value(2.1); + model4[keyname].error(1.9); + model4[keyname].gradient(0.8); + test_value(model4[keyname].value(), 2.1); + test_value(model4[keyname].error(), 1.9); + test_value(model4[keyname].gradient(), 0.8); } // Return