Luigi Tibaldo / gammalib

Commit d26b95a8e5ebed9999f73b51493b1c29fcfc1937

Authored by Josh Cardenzana
1 parent f3c8e222

Updated tests for GModelSpectralSmoothBrokenPlaw class (#1948) 

Added file "GModelSpectralSmoothBrokenPlaw.i" for python wrappers.
Added file "model_point_smoothbplaw.xml" for testing the construction of the class from xml file.
Updated "test_GModel.cpp" to include test of the 'GModelSpectralSmoothBrokenPlaw' class.
Some work is still needed to ensure the test passes when running 'make check'.
Inline Side-by-side
Showing 5 changed files with 236 additions and 19 deletions
pyext/GModelSpectralSmoothBrokenPlaw.i 0 → 100644
  View file @d26b95a
  1 +/***************************************************************************
  2 + * GModelSpectralSmoothBrokenPlaw.i - *
  3 + * Smoothly broken power law spectrum class *
  4 + * ----------------------------------------------------------------------- *
  5 + * copyright (C) 2017 by Josh Cardenzana *
  6 + * ----------------------------------------------------------------------- *
  7 + * *
  8 + * This program is free software: you can redistribute it and/or modify *
  9 + * it under the terms of the GNU General Public License as published by *
  10 + * the Free Software Foundation, either version 3 of the License, or *
  11 + * (at your option) any later version. *
  12 + * *
  13 + * This program is distributed in the hope that it will be useful, *
  14 + * but WITHOUT ANY WARRANTY; without even the implied warranty of *
  15 + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
  16 + * GNU General Public License for more details. *
  17 + * *
  18 + * You should have received a copy of the GNU General Public License *
  19 + * along with this program. If not, see <http://www.gnu.org/licenses/>. *
  20 + * *
  21 + ***************************************************************************/
  22 +/**
  23 + * @file GModelSpectralSmoothBrokenPlaw.i
  24 + * @brief Smoothly broken power law spectrum class definition
  25 + * @author Josh Cardenzana
  26 + */
  27 +%{
  28 + /* Put headers and other declarations here that are needed for compilation */
  29 +#include "GModelSpectralSmoothBrokenPlaw.hpp"
  30 + %}
  31 +
  32 +
  33 +/***********************************************************************//**
  34 + * @class GModelSpectralSmoothBrokenPlaw
  35 + *
  36 + * @brief Smoothly broken Power law spectral model class
  37 + ***************************************************************************/
  38 +class GModelSpectralSmoothBrokenPlaw : public GModelSpectral {
  39 +public:
  40 + // Constructors and destructors
  41 + GModelSpectralSmoothBrokenPlaw(void);
  42 + GModelSpectralSmoothBrokenPlaw(const std::string& type,
  43 + const std::string& prefactor,
  44 + const std::string& index1,
  45 + const std::string& pivot,
  46 + const std::string& index2,
  47 + const std::string& breakenergy,
  48 + const std::string& beta);
  49 + GModelSpectralSmoothBrokenPlaw(const double& prefactor,
  50 + const double& index1,
  51 + const GEnergy& pivot,
  52 + const double& index2,
  53 + const GEnergy& breakenergy,
  54 + const double& beta);
  55 + explicit GModelSpectralSmoothBrokenPlaw(const GXmlElement& xml);
  56 + GModelSpectralSmoothBrokenPlaw(const GModelSpectralSmoothBrokenPlaw& model);
  57 + virtual ~GModelSpectralSmoothBrokenPlaw(void);
  58 +
  59 + // Implemented pure virtual methods
  60 + virtual void clear(void);
  61 + virtual GModelSpectralSmoothBrokenPlaw* clone(void) const;
  62 + virtual std::string classname(void) const;
  63 + virtual std::string type(void) const;
  64 + virtual double eval(const GEnergy& srcEng,
  65 + const GTime& srcTime = GTime(),
  66 + const bool& gradients = false) const;
  67 + virtual double flux(const GEnergy& emin,
  68 + const GEnergy& emax) const;
  69 + virtual double eflux(const GEnergy& emin,
  70 + const GEnergy& emax) const;
  71 + virtual GEnergy mc(const GEnergy& emin,
  72 + const GEnergy& emax,
  73 + const GTime& time,
  74 + GRan& ran) const;
  75 + virtual void read(const GXmlElement& xml);
  76 + virtual void write(GXmlElement& xml) const;
  77 +
  78 + // Other methods
  79 + double prefactor(void) const;
  80 + double index1(void) const;
  81 + double index2(void) const;
  82 + GEnergy pivot(void) const;
  83 + GEnergy breakenergy(void) const;
  84 + void prefactor(const double& prefactor);
  85 + void index1(const double& index1);
  86 + void index2(const double& index2);
  87 + void pivot(const GEnergy& pivot);
  88 + void breakenergy(const GEnergy& breakenergy);
  89 + void beta(const double& beta);
  90 +};
  91 +
  92 +
  93 +/***********************************************************************//**
  94 + * @brief GModelSpectralSmoothBrokenPlaw class extension
  95 + ***************************************************************************/
  96 +%extend GModelSpectralSmoothBrokenPlaw {
  97 + GModelSpectralSmoothBrokenPlaw copy() {
  98 + return (*self);
  99 + }
  100 +};
... ...
src/model/GModelSpectralSmoothBrokenPlaw.cpp
  View file @d26b95a
... ... @@ -40,24 +40,23 @@
40 40 /* __ Constants __________________________________________________________ */
41 41  
42 42 /* __ Globals ____________________________________________________________ */
43   -const GModelSpectralSmoothBrokenPlaw g_spectral_blaw_seed1("SmoothBrokenPowerLaw",
44   - "Prefactor",
45   - "Index1",
46   - "PivotEnergy",
47   - "Index2",
48   - "BreakEnergy",
49   - "BreakSmoothness");
50   -const GModelSpectralRegistry g_spectral_blaw_registry1(&g_spectral_blaw_seed1);
51   -#if defined(G_LEGACY_XML_FORMAT)
52   -const GModelSpectralSmoothBrokenPlaw g_spectral_blaw_seed2("SmoothBrokenPowerLaw",
53   - "Prefactor",
54   - "Index1",
55   - "Scale",
56   - "Index2",
57   - "BreakValue",
58   - "Beta");
59   -const GModelSpectralRegistry g_spectral_blaw_registry2(&g_spectral_blaw_seed2);
60   -#endif
  43 +const GModelSpectralSmoothBrokenPlaw g_spectral_sblaw_seed1("SmoothBrokenPowerLaw",
  44 + "Prefactor",
  45 + "Index1",
  46 + "PivotEnergy",
  47 + "Index2",
  48 + "BreakEnergy",
  49 + "BreakSmoothness");
  50 +const GModelSpectralSmoothBrokenPlaw g_spectral_sblaw_seed2("SmoothBrokenPowerLaw",
  51 + "Prefactor",
  52 + "Index1",
  53 + "Scale",
  54 + "Index2",
  55 + "BreakValue",
  56 + "Beta");
  57 +const GModelSpectralRegistry g_spectral_sblaw_registry1(&g_spectral_sblaw_seed1);
  58 +const GModelSpectralRegistry g_spectral_sblaw_registry2(&g_spectral_sblaw_seed2);
  59 +
61 60  
62 61 /* __ Method name definitions ____________________________________________ */
63 62 #define G_MC "GModelSpectralSmoothBrokenPlaw::mc(GEnergy&, GEnergy&, GTime&,"\
... ... @@ -851,7 +850,7 @@ void GModelSpectralSmoothBrokenPlaw::init_members(void)
851 850  
852 851 // Initialize beta (break smoothness parameter)
853 852 m_beta.clear();
854   - m_beta.name("Beta");
  853 + m_beta.name("BreakSmoothness");
855 854 m_beta.scale(1.0);
856 855 m_beta.value(1.0);
857 856 m_beta.range(0.01,10.0);
... ...
test/data/model_point_smoothbplaw.xml 0 → 100644
  View file @d26b95a
  1 +<?xml version="1.0" standalone="no"?>
  2 +<source_library title="source library">
  3 + <source name="Crab1" type="PointSource">
  4 + <spectrum type="SmoothBrokenPowerLaw">
  5 + <parameter name="Prefactor" scale="1e-16" value="5.7" min="1e-07" max="1000.0" free="1"/>
  6 + <parameter name="Index1" scale="-1" value="2.48" min="0.0" max="+5.0" free="1"/>
  7 + <parameter name="Scale" scale="1e6" value="1.0" min="0.01" max="1000.0" free="0"/>
  8 + <parameter name="Index2" scale="-1" value="2.70" min="0.01" max="+5.0" free="1"/>
  9 + <parameter name="BreakValue" scale="1e6" value="0.3" min="0.01" max="1000.0" free="1"/>
  10 + <parameter name="Beta" scale="1.0" value="0.2" min="0.01" max="10" free="1"/>
  11 + </spectrum>
  12 + <spatialModel type="PointSource">
  13 + <parameter free="0" max="360" min="-360" name="RA" scale="1" value="83.6331" />
  14 + <parameter free="0" max="90" min="-90" name="DEC" scale="1" value="22.0145" />
  15 + </spatialModel>
  16 + </source>
  17 + <source name="Crab2" type="PointSource">
  18 + <spectrum type="SmoothBrokenPowerLaw">
  19 + <parameter name="Prefactor" scale="1e-16" value="5.7" min="1e-07" max="1000.0" free="1"/>
  20 + <parameter name="Index1" scale="-1" value="2.48" min="0.0" max="+5.0" free="1"/>
  21 + <parameter name="PivotEnergy" scale="1e6" value="1.0" min="0.01" max="1000.0" free="0"/>
  22 + <parameter name="Index2" scale="-1" value="2.70" min="0.01" max="+5.0" free="1"/>
  23 + <parameter name="BreakEnergy" scale="1e6" value="0.3" min="0.01" max="1000.0" free="1"/>
  24 + <parameter name="BreakSmoothness" scale="1.0" value="0.2" min="0.01" max="10" free="1"/>
  25 + </spectrum>
  26 + <spatialModel type="PointSource">
  27 + <parameter free="0" max="360" min="-360" name="RA" scale="1" value="83.6331" />
  28 + <parameter free="0" max="90" min="-90" name="DEC" scale="1" value="22.0145" />
  29 + </spatialModel>
  30 + </source>
  31 +</source_library>
... ...
test/test_GModel.cpp
  View file @d26b95a
... ... @@ -64,6 +64,7 @@ void TestGModel::set(void)
64 64 m_xml_model_point_eplaw = datadir + "/model_point_eplaw.xml";
65 65 m_xml_model_point_einvplaw = datadir + "/model_point_einvplaw.xml";
66 66 m_xml_model_point_bplaw = datadir + "/model_point_bplaw.xml";
  67 + m_xml_model_point_smoothbplaw = datadir + "/model_point_smoothbplaw.xml";
67 68 m_xml_model_point_supeplaw = datadir + "/model_point_supeplaw.xml";
68 69 m_xml_model_point_logparabola = datadir + "/model_point_logparabola.xml";
69 70 m_xml_point_multiplicative = datadir + "/model_point_multiplicative.xml";
... ... @@ -148,6 +149,8 @@ void TestGModel::set(void)
148 149 "Test GModelSpectralSuperExpPlaw");
149 150 append(static_cast<pfunction>(&TestGModel::test_bplaw),
150 151 "Test GModelSpectralBrokenPlaw");
  152 + append(static_cast<pfunction>(&TestGModel::test_smoothbplaw),
  153 + "Test GModelSpectralSmoothBrokenPlaw");
151 154 append(static_cast<pfunction>(&TestGModel::test_logparabola),
152 155 "Test GModelSpectralLogParabola");
153 156 append(static_cast<pfunction>(&TestGModel::test_multiplicative),
... ... @@ -1839,6 +1842,87 @@ void TestGModel::test_bplaw(void)
1839 1842  
1840 1843  
1841 1844 /***********************************************************************//**
  1845 + * @brief Test GModelSpectralSmoothBrokenPlaw class
  1846 + ***************************************************************************/
  1847 +void TestGModel::test_smoothbplaw(void)
  1848 +{
  1849 + // Test void constructor
  1850 + GModelSpectralSmoothBrokenPlaw model1;
  1851 + test_value(model1.type(), "SmoothBrokenPowerLaw", "Check void model type");
  1852 +
  1853 + // Test value constructor
  1854 + GModelSpectralSmoothBrokenPlaw model2(2.0, -2.1, GEnergy(100.0, "MeV"), -2.8,
  1855 + GEnergy(1000.0, "MeV"), 2.0);
  1856 + test_value(model2.prefactor(), 2.0);
  1857 + test_value(model2.index1(), -2.1);
  1858 + test_value(model2.pivot().MeV(), 100.0);
  1859 + test_value(model2.index2(), -2.8);
  1860 + test_value(model2.breakenergy().MeV(), 1000.0);
  1861 + test_value(model2.beta(), 4.0);
  1862 +
  1863 + // Test XML constructor - 1
  1864 + GXml xml(m_xml_model_point_smoothbplaw);
  1865 + GXmlElement* element = xml.element(0)->element(0)->element("spectrum", 0);
  1866 + GModelSpectralSmoothBrokenPlaw model3(*element);
  1867 + test_value(model3.size(), 4);
  1868 + test_value(model3.type(), "SmoothBrokenPowerLaw", "Check model type");
  1869 + test_value(model3.prefactor(), 5.7e-16);
  1870 + test_value(model3.index1(), -2.48);
  1871 + test_value(model3.pivot().TeV(), 1.0);
  1872 + test_value(model3.index2(), -2.70);
  1873 + test_value(model3.breakenergy().TeV(), 0.3);
  1874 +
  1875 + // Test XML constructor - 2
  1876 + GXmlElement* element2 = xml.element(0)->element(1)->element("spectrum", 0);
  1877 + GModelSpectralSmoothBrokenPlaw model4(*element2);
  1878 + test_value(model4.size(), 4);
  1879 + test_value(model4.type(), "SmoothBrokenPowerLaw", "Check model type");
  1880 + test_value(model4.prefactor(), 5.7e-16);
  1881 + test_value(model4.index1(), -2.48);
  1882 + test_value(model4.pivot().TeV(), 1.0);
  1883 + test_value(model4.index2(), -2.70);
  1884 + test_value(model4.breakenergy().TeV(), 0.3);
  1885 +
  1886 + // Test prefactor method
  1887 + model3.prefactor(2.3e-16);
  1888 + test_value(model3.prefactor(), 2.3e-16);
  1889 +
  1890 + // Test index1 method
  1891 + model3.index1(-2.6);
  1892 + test_value(model3.index1(), -2.6);
  1893 +
  1894 + // Test pivot method
  1895 + model3.pivot(GEnergy(0.5, "TeV"));
  1896 + test_value(model3.pivot().TeV(), 0.5);
  1897 +
  1898 + // Test index2 method
  1899 + model3.index2(-3.6);
  1900 + test_value(model3.index2(), -3.6);
  1901 +
  1902 + // Test breakenergy method
  1903 + model3.breakenergy(GEnergy(0.5, "TeV"));
  1904 + test_value(model3.breakenergy().TeV(), 0.5);
  1905 +
  1906 + // Test operator access
  1907 + const char* strarray[] = {"Prefactor", "Index1", "Pivot", "Index2",
  1908 + "BreakEnergy", "Beta"};
  1909 + for (int i = 0; i < 4; ++i) {
  1910 + std::string keyname(strarray[i]);
  1911 + model3[keyname].remove_range(); // To allow setting of any value
  1912 + model3[keyname].value(2.1);
  1913 + model3[keyname].error(1.9);
  1914 + model3[keyname].gradient(0.8);
  1915 + test_value(model3[keyname].value(), 2.1);
  1916 + test_value(model3[keyname].error(), 1.9);
  1917 + test_value(model3[keyname].gradient(), 0.8);
  1918 + }
  1919 +
  1920 + // Exit test
  1921 + return;
  1922 +}
  1923 +
  1924 +
  1925 +/***********************************************************************//**
1842 1926 * @brief Test GModelSpectralLogParabola class
1843 1927 ***************************************************************************/
1844 1928 void TestGModel::test_logparabola(void)
... ... @@ -2346,6 +2430,7 @@ void TestGModel::test_spectral_model(void)
2346 2430 test_xml_model("GModelSpectralExpPaw", m_xml_model_point_eplaw);
2347 2431 test_xml_model("GModelSpectralExpInvPaw", m_xml_model_point_einvplaw);
2348 2432 test_xml_model("GModelSpectralBrokenPlaw", m_xml_model_point_bplaw);
  2433 + test_xml_model("GModelSpectralSmoothBrokenPlaw",m_xml_model_point_smoothbplaw);
2349 2434 test_xml_model("GModelSpectralSuperExpPlaw", m_xml_model_point_supeplaw);
2350 2435 test_xml_model("GModelSpectralLogParabola", m_xml_model_point_logparabola);
2351 2436 test_xml_model("GModelSpectralNodes", m_xml_model_point_nodes);
... ...
test/test_GModel.hpp
  View file @d26b95a
... ... @@ -68,6 +68,7 @@ public:
68 68 void test_eplaw(void);
69 69 void test_einvplaw(void);
70 70 void test_bplaw(void);
  71 + void test_smoothbplaw(void);
71 72 void test_supeplaw(void);
72 73 void test_logparabola(void);
73 74 void test_multiplicative(void);
... ... @@ -116,6 +117,7 @@ private:
116 117 std::string m_xml_model_point_eplaw;
117 118 std::string m_xml_model_point_einvplaw;
118 119 std::string m_xml_model_point_bplaw;
  120 + std::string m_xml_model_point_smoothbplaw;
119 121 std::string m_xml_model_point_supeplaw;
120 122 std::string m_xml_model_point_logparabola;
121 123 std::string m_xml_point_multiplicative;
... ...