diff --git a/test/data/model_point_gauss.xml b/test/data/model_point_gauss.xml
new file mode 100644
index 0000000..c9a01fb
--- /dev/null
+++ b/test/data/model_point_gauss.xml
@@ -0,0 +1,14 @@
+<?xml version="1.0" standalone="no"?>
+<source_library title="source library">
+	<source name="Crab" type="PointSource">
+        <spectrum type="Gaussian">
+           <parameter name="Prefactor" scale="1e-16" value="5.7"  min="1e-07" max="1000.0" free="1"/>
+           <parameter name="Mean"      scale="1e6"   value="2.48" min="0.0"   max="+5.0"   free="1"/>
+           <parameter name="Sigma"     scale="1e6"   value="1.0"  min="0.0001"  max="1000.0" free="1"/>
+        </spectrum>
+		<spatialModel type="SkyDirFunction">
+			<parameter free="0" max="360" min="-360" name="RA" scale="1" value="83.6331" />
+			<parameter free="0" max="90" min="-90" name="DEC" scale="1" value="22.0145" />
+		</spatialModel>
+	</source>
+</source_library>
diff --git a/test/test_GModel.cpp b/test/test_GModel.cpp
index 7624408..a7c3bcb 100644
--- a/test/test_GModel.cpp
+++ b/test/test_GModel.cpp
@@ -49,6 +49,7 @@ void TestGModel::set(void)
     m_cube_file                   = "data/test_cube.fits";
     m_xml_file                    = "data/crab.xml";
     m_xml_model_point_const       = "data/model_point_const.xml";
+    m_xml_model_point_gauss       = "data/model_point_gauss.xml";
     m_xml_model_point_plaw        = "data/model_point_plaw.xml";
     m_xml_model_point_plaw2       = "data/model_point_plaw2.xml";
     m_xml_model_point_eplaw       = "data/model_point_eplaw.xml";
@@ -81,6 +82,7 @@ void TestGModel::set(void)
 
     // Append spectral model tests
     append(static_cast<pfunction>(&TestGModel::test_const), "Test GModelSpectralConst");
+    append(static_cast<pfunction>(&TestGModel::test_gauss), "Test GModelSpectralGauss");
     append(static_cast<pfunction>(&TestGModel::test_plaw), "Test GModelSpectralPlaw");
     append(static_cast<pfunction>(&TestGModel::test_plaw2), "Test GModelSpectralPlaw2");
     append(static_cast<pfunction>(&TestGModel::test_eplaw), "Test GModelSpectralExpPlaw");
@@ -1154,6 +1156,86 @@ void TestGModel::test_const(void)
 
 
 /***********************************************************************//**
+ * @brief Test GModelSpectralGauss class
+ ***************************************************************************/
+void TestGModel::test_gauss(void)
+{
+    // Test void constructor
+    test_try("Test void constructor");
+    try {
+        GModelSpectralGauss model;
+        test_assert(model.type() == "Gaussian",
+                                    "Model type \"Gaussian\" expected.");
+        test_try_success();
+    }
+    catch (std::exception &e) {
+        test_try_failure(e);
+    }
+
+    // Test value constructor
+    test_try("Test value constructor");
+    try {
+        GModelSpectralGauss model(1.0, GEnergy(42.0, "MeV"), GEnergy(0.5, "MeV"));
+        test_value(model.prefactor(), 1.0);
+        test_value(model.mean().MeV(), 42.0);
+        test_value(model.sigma().MeV(), 0.5);
+        test_try_success();
+    }
+    catch (std::exception &e) {
+        test_try_failure(e);
+    }
+
+    // Test XML constructor and value
+    test_try("Test XML constructor, value and gradients");
+    try {
+        // Test XML constructor
+        GXml                  xml(m_xml_model_point_gauss);
+        GXmlElement*          element = xml.element(0)->element(0)->element("spectrum", 0);
+        GModelSpectralGauss model(*element);
+        test_value(model.size(), 3);
+        test_assert(model.type() == "Gaussian", "Expected \"Gaussian\"");
+        test_value(model.prefactor(), 5.7e-16);
+        test_value(model.mean().MeV(), 2.48e6);
+        test_value(model.sigma().MeV(), 1e6);
+
+        // Test prefactor method
+        model.prefactor(2.3e-16);
+        test_value(model.prefactor(), 2.3e-16);
+
+        // Test mean method
+        model.mean(GEnergy(96.0, "MeV"));
+        test_value(model.mean().MeV(), 96.0);
+
+        // Test sigma method
+        model.sigma(GEnergy(0.5e7, "keV"));
+        test_value(model.sigma().keV(), 0.5e7);
+
+        // Test operator access
+        const char* strarray[] = {"Prefactor", "Mean", "Sigma"};
+        for (int i = 0; i < 3; ++i) {
+            std::string keyname(strarray[i]);
+            model[keyname].remove_range(); // To allow setting of any value
+            model[keyname].value(2.1);
+            model[keyname].error(1.9);
+            model[keyname].gradient(0.8);
+            test_value(model[keyname].value(), 2.1);
+            test_value(model[keyname].error(), 1.9);
+            test_value(model[keyname].gradient(), 0.8);
+        }
+
+        // Success if we reached this point
+        test_try_success();
+    }
+    catch (std::exception &e) {
+        test_try_failure(e);
+    }
+
+    // Exit test
+    return;
+}
+
+
+/***********************************************************************//**
  * @brief Test GModelSpectralPlaw class
  ***************************************************************************/
 void TestGModel::test_plaw(void)
diff --git a/test/test_GModel.hpp b/test/test_GModel.hpp
index baac6ac..9c0836b 100644
--- a/test/test_GModel.hpp
+++ b/test/test_GModel.hpp
@@ -58,6 +58,7 @@ public:
     void                test_elliptical_disk(void);
     void                test_spatial_model(void);
     void                test_const(void);
+    void                test_gauss(void);
     void                test_plaw(void);
     void                test_plaw2(void);
     void                test_eplaw(void);
@@ -79,6 +80,7 @@ private:
     std::string m_cube_file;
     std::string m_xml_file;
     std::string m_xml_model_point_const;
+    std::string m_xml_model_point_gauss;
     std::string m_xml_model_point_plaw;
     std::string m_xml_model_point_plaw2;
     std::string m_xml_model_point_eplaw;