diff --git a/doc/source/user_manual/modules/cta.rst b/doc/source/user_manual/modules/cta.rst
index 088b812..41b3167 100644
--- a/doc/source/user_manual/modules/cta.rst
+++ b/doc/source/user_manual/modules/cta.rst
@@ -98,6 +98,7 @@ the format:
+
for an unbinned observation and:
@@ -109,6 +110,7 @@ for an unbinned observation and:
+
for a binned observation. Each parameter within the ```` tag
@@ -129,6 +131,8 @@ about the off-axis reponse and event selection (see
So far, energy dispersion is not supported by GammaLib, hence the
``file`` attribute for the ``EnergyDispersion`` parameter can be left
blank.
+The ``Background`` parameter allows specifying of a file that contains
+the instrument background model for the observation.
The ``instrument`` attribute of CTA observations can be one of
``CTA``, ``HESS``, ``MAGIC`` or ``VERITAS``. This allows mixing of
diff --git a/inst/cta/include/GCTAResponse.hpp b/inst/cta/include/GCTAResponse.hpp
index c940b06..9fb1ad3 100644
--- a/inst/cta/include/GCTAResponse.hpp
+++ b/inst/cta/include/GCTAResponse.hpp
@@ -51,6 +51,7 @@ class GCTAInstDir;
class GCTAAeff;
class GCTAPsf;
class GCTAEdisp;
+class GCTABackground;
/***********************************************************************//**
@@ -103,24 +104,27 @@ public:
virtual GEbounds ebounds_src(const GEnergy& obsEnergy) const;
// Other Methods
- GCTAEventAtom* mc(const double& area, const GPhoton& photon,
- const GObservation& obs, GRan& ran) const;
- void caldb(const GCaldb& caldb);
- const GCaldb& caldb(void) const;
- void load(const std::string& rspname);
- void eps(const double& eps);
- const double& eps(void) const;
- void load_aeff(const std::string& filename);
- void load_psf(const std::string& filename);
- void load_edisp(const std::string& filename);
- void offset_sigma(const double& sigma);
- double offset_sigma(void) const;
- const GCTAAeff* aeff(void) const;
- void aeff(GCTAAeff* aeff);
- const GCTAPsf* psf(void) const;
- void psf(GCTAPsf* psf);
- const GCTAEdisp* edisp(void) const;
- void edisp(GCTAEdisp* edisp);
+ GCTAEventAtom* mc(const double& area, const GPhoton& photon,
+ const GObservation& obs, GRan& ran) const;
+ void caldb(const GCaldb& caldb);
+ const GCaldb& caldb(void) const;
+ void load(const std::string& rspname);
+ void eps(const double& eps);
+ const double& eps(void) const;
+ void load_aeff(const std::string& filename);
+ void load_psf(const std::string& filename);
+ void load_edisp(const std::string& filename);
+ void load_background(const std::string& filename);
+ void offset_sigma(const double& sigma);
+ double offset_sigma(void) const;
+ const GCTAAeff* aeff(void) const;
+ void aeff(GCTAAeff* aeff);
+ const GCTAPsf* psf(void) const;
+ void psf(GCTAPsf* psf);
+ const GCTAEdisp* edisp(void) const;
+ void edisp(GCTAEdisp* edisp);
+ const GCTABackground* background(void) const;
+ void background(GCTABackground* background);
// Low-level response methods
double aeff(const double& theta,
@@ -172,13 +176,14 @@ private:
std::string irf_filename(const std::string& filename) const;
// Private data members
- GCaldb m_caldb; //!< Calibration database
- std::string m_rspname; //!< Name of the instrument response
- double m_eps; //!< Integration precision
- GCTAAeff* m_aeff; //!< Effective area
- GCTAPsf* m_psf; //!< Point spread function
- GCTAEdisp* m_edisp; //!< Energy dispersion
- mutable bool m_apply_edisp; //!< Apply energy dispersion
+ GCaldb m_caldb; //!< Calibration database
+ std::string m_rspname; //!< Name of the instrument response
+ double m_eps; //!< Integration precision
+ GCTAAeff* m_aeff; //!< Effective area
+ GCTAPsf* m_psf; //!< Point spread function
+ GCTAEdisp* m_edisp; //!< Energy dispersion
+ GCTABackground* m_background; //!< Energy dispersion
+ mutable bool m_apply_edisp; //!< Apply energy dispersion
// Npred cache
mutable std::vector m_npred_names; //!< Model names
@@ -366,4 +371,31 @@ void GCTAResponse::edisp(GCTAEdisp* edisp)
return;
}
+
+
+
+/***********************************************************************//**
+ * @brief Return pointer to background model
+ *
+ * @return Pointer to background model.
+ ***************************************************************************/
+inline
+const GCTABackground* GCTAResponse::background(void) const
+{
+ return m_background;
+}
+
+
+/***********************************************************************//**
+ * @brief Set pointer to background model
+ *
+ * @param[in] edisp Pointer to background model.
+ ***************************************************************************/
+inline
+void GCTAResponse::background(GCTABackground* background)
+{
+ m_background = background;
+ return;
+}
+
#endif /* GCTARESPONSE_HPP */
diff --git a/inst/cta/pyext/GCTAResponse.i b/inst/cta/pyext/GCTAResponse.i
index 18bd4ee..a558320 100644
--- a/inst/cta/pyext/GCTAResponse.i
+++ b/inst/cta/pyext/GCTAResponse.i
@@ -74,24 +74,27 @@ public:
virtual GEbounds ebounds_src(const GEnergy& obsEnergy) const;
// Other Methods
- GCTAEventAtom* mc(const double& area, const GPhoton& photon,
- const GObservation& obs, GRan& ran) const;
- void caldb(const GCaldb& caldb);
- const GCaldb& caldb(void) const;
- void load(const std::string& rspname);
- void eps(const double& eps);
- const double& eps(void) const;
- void load_aeff(const std::string& filename);
- void load_psf(const std::string& filename);
- void load_edisp(const std::string& filename);
- void offset_sigma(const double& sigma);
- double offset_sigma(void) const;
- const GCTAAeff* aeff(void) const;
- void aeff(GCTAAeff* aeff);
- const GCTAPsf* psf(void) const;
- void psf(GCTAPsf* psf);
- const GCTAEdisp* edisp(void) const;
- void edisp(GCTAEdisp* edisp);
+ GCTAEventAtom* mc(const double& area, const GPhoton& photon,
+ const GObservation& obs, GRan& ran) const;
+ void caldb(const GCaldb& caldb);
+ const GCaldb& caldb(void) const;
+ void load(const std::string& rspname);
+ void eps(const double& eps);
+ const double& eps(void) const;
+ void load_aeff(const std::string& filename);
+ void load_psf(const std::string& filename);
+ void load_edisp(const std::string& filename);
+ void load_background(const std::string& filename);
+ void offset_sigma(const double& sigma);
+ double offset_sigma(void) const;
+ const GCTAAeff* aeff(void) const;
+ void aeff(GCTAAeff* aeff);
+ const GCTAPsf* psf(void) const;
+ void psf(GCTAPsf* psf);
+ const GCTAEdisp* edisp(void) const;
+ void edisp(GCTAEdisp* edisp);
+ const GCTABackground* background(void) const;
+ void background(GCTABackground* background);
// Low-level response methods
double aeff(const double& theta,
diff --git a/inst/cta/src/GCTAObservation.cpp b/inst/cta/src/GCTAObservation.cpp
index 40a3e0f..8ffae43 100644
--- a/inst/cta/src/GCTAObservation.cpp
+++ b/inst/cta/src/GCTAObservation.cpp
@@ -44,6 +44,7 @@
#include "GCTAAeffArf.hpp"
#include "GCTAAeffPerfTable.hpp"
#include "GCTAEdisp.hpp"
+#include "GCTABackground.hpp"
/* __ Globals ____________________________________________________________ */
const GCTAObservation g_obs_cta_seed("CTA");
@@ -288,6 +289,7 @@ void GCTAObservation::response(const std::string& rspname, const GCaldb& caldb)
*
*
*
+ *
*
*
* for an unbinned observation and
@@ -297,12 +299,12 @@ void GCTAObservation::response(const std::string& rspname, const GCaldb& caldb)
*
*
*
+ *
*
*
* for a binned observation.
*
* @todo Still supports old ARF, PSF and RMF parameter names.
- * @todo PSF correction for ARF not yet implemented.
***************************************************************************/
void GCTAObservation::read(const GXmlElement& xml)
{
@@ -315,14 +317,14 @@ void GCTAObservation::read(const GXmlElement& xml)
// Determine number of parameter nodes in XML element
int npars = xml.elements("parameter");
- // Verify that XML element has exactly 4 parameters
- if (xml.elements() != 4 || npars != 4) {
+ // Verify that XML element has exactly 5 parameters
+ if (xml.elements() != 5 || npars != 5) {
throw GException::xml_invalid_parnum(G_READ, xml,
- "CTA observation requires exactly 4 parameters.");
+ "CTA observation requires exactly 5 parameters.");
}
// Extract parameters
- int npar[] = {0, 0, 0, 0};
+ int npar[] = {0, 0, 0, 0, 0};
for (int i = 0; i < npars; ++i) {
// Get parameter element
@@ -443,13 +445,29 @@ void GCTAObservation::read(const GXmlElement& xml)
npar[3]++;
}
+ // Handle background model
+ else if (par->attribute("name") == "Background") {
+
+ // Get filename
+ std::string filename = par->attribute("file");
+
+ // If filename is not empty then load background model
+ if (gammalib::strip_whitespace(filename).length() > 0) {
+ m_response.load_background(filename);
+ }
+
+ // Increase number of parameters
+ npar[4]++;
+ }
+
} // endfor: looped over all parameters
// Verify that all parameters were found
- if (npar[0] != 1 || npar[1] != 1 || npar[2] != 1 || npar[3] != 1) {
+ if (npar[0] != 1 || npar[1] != 1 || npar[2] != 1 || npar[3] != 1 || npar[4] != 1) {
throw GException::xml_invalid_parnames(G_READ, xml,
"Require \"EventList\" or \"CountsMap\" and \"EffectiveArea\""
- ", \"PointSpreadFunction\" and \"EnergyDispersion\" parameters.");
+ ", \"PointSpreadFunction\", \"EnergyDispersion\" and"
+ " \"Background\" parameters.");
}
// If we have an ARF then remove thetacut if necessary
@@ -485,6 +503,7 @@ void GCTAObservation::read(const GXmlElement& xml)
*
*
*
+ *
*
*
* for an unbinned observation and
@@ -494,6 +513,7 @@ void GCTAObservation::read(const GXmlElement& xml)
*
*
*
+ *
*
*
* for a binned observation.
@@ -514,7 +534,7 @@ void GCTAObservation::write(GXmlElement& xml) const
// Set event list flag
bool is_list = (list != NULL);
- // If XML element has 0 nodes then append 4 parameter nodes
+ // If XML element has 0 nodes then append 5 parameter nodes
if (xml.elements() == 0) {
if (is_list) {
xml.append(GXmlElement("parameter name=\"EventList\""));
@@ -525,17 +545,18 @@ void GCTAObservation::write(GXmlElement& xml) const
xml.append(GXmlElement("parameter name=\"EffectiveArea\""));
xml.append(GXmlElement("parameter name=\"PointSpreadFunction\""));
xml.append(GXmlElement("parameter name=\"EnergyDispersion\""));
+ xml.append(GXmlElement("parameter name=\"Background\""));
}
- // Verify that XML element has exactly 4 parameters
- if (xml.elements() != 4 || xml.elements("parameter") != 4) {
+ // Verify that XML element has exactly 5 parameters
+ if (xml.elements() != 5 || xml.elements("parameter") != 5) {
throw GException::xml_invalid_parnum(G_WRITE, xml,
- "CTA observation requires exactly 4 parameters.");
+ "CTA observation requires exactly 5 parameters.");
}
// Set or update parameter attributes
- int npar[] = {0, 0, 0, 0};
- for (int i = 0; i < 4; ++i) {
+ int npar[] = {0, 0, 0, 0, 0};
+ for (int i = 0; i < 5; ++i) {
// Get parameter element
GXmlElement* par = xml.element("parameter", i);
@@ -619,13 +640,24 @@ void GCTAObservation::write(GXmlElement& xml) const
npar[3]++;
}
+ // Handle Background
+ else if (par->attribute("name") == "Background") {
+ std::string filename = "";
+ if (m_response.background() != NULL) {
+ filename = m_response.background()->filename();
+ }
+ par->attribute("file", filename);
+ npar[4]++;
+ }
+
} // endfor: looped over all parameters
// Verify that all required parameters are present
- if (npar[0] != 1 || npar[1] != 1 || npar[2] != 1 || npar[3] != 1) {
+ if (npar[0] != 1 || npar[1] != 1 || npar[2] != 1 || npar[3] != 1 || npar[4] != 1) {
throw GException::xml_invalid_parnames(G_WRITE, xml,
- "Require \"EventList\" or \"CountsMap\" and \"ARF\", \"RMF\""
- " \"PSF\" parameters.");
+ "Require \"EventList\" or \"CountsMap\" and"
+ " \"PointSpreadFunction\", \"EnergyDispersion\""
+ " and \"Background\" parameters.");
}
// Return
diff --git a/inst/cta/src/GCTAResponse.cpp b/inst/cta/src/GCTAResponse.cpp
index 9e5f5fd..3332051 100644
--- a/inst/cta/src/GCTAResponse.cpp
+++ b/inst/cta/src/GCTAResponse.cpp
@@ -58,6 +58,7 @@
#include "GCTAPsf.hpp"
#include "GCTAEdisp.hpp"
#include "GCTAEdispPerfTable.hpp"
+#include "GCTABackground.hpp"
/* __ Method name definitions ____________________________________________ */
#define G_CALDB "GCTAResponse::caldb(GCaldb&)"
@@ -737,7 +738,7 @@ void GCTAResponse::load_psf(const std::string& filename)
***************************************************************************/
void GCTAResponse::load_edisp(const std::string& filename)
{
- // Free any existing point spread function instance
+ // Free any existing energy dispersion instance
if (m_edisp != NULL) delete m_edisp;
m_edisp = NULL;
@@ -752,6 +753,25 @@ void GCTAResponse::load_edisp(const std::string& filename)
m_edisp = new GCTAEdispPerfTable(filename);
}
+ // Return
+ return;
+}
+
+
+/***********************************************************************//**
+ * @brief Load background model
+ *
+ * @param[in] filename Background model file name.
+ *
+ * @todo Method to be implemented.
+ ***************************************************************************/
+void GCTAResponse::load_background(const std::string& filename)
+{
+ // Free any existing background model instance
+ if (m_background != NULL) delete m_background;
+ m_background = NULL;
+
+ // Return
return;
}
@@ -2349,6 +2369,7 @@ void GCTAResponse::init_members(void)
m_aeff = NULL;
m_psf = NULL;
m_edisp = NULL;
+ m_background = NULL;
m_apply_edisp = false; //!< Switched off by default
// Initialise Npred cache
@@ -2382,9 +2403,10 @@ void GCTAResponse::copy_members(const GCTAResponse& rsp)
m_npred_values = rsp.m_npred_values;
// Clone members
- m_aeff = (rsp.m_aeff != NULL) ? rsp.m_aeff->clone() : NULL;
- m_psf = (rsp.m_psf != NULL) ? rsp.m_psf->clone() : NULL;
- m_edisp = (rsp.m_edisp != NULL) ? rsp.m_edisp->clone() : NULL;
+ m_aeff = (rsp.m_aeff != NULL) ? rsp.m_aeff->clone() : NULL;
+ m_psf = (rsp.m_psf != NULL) ? rsp.m_psf->clone() : NULL;
+ m_edisp = (rsp.m_edisp != NULL) ? rsp.m_edisp->clone() : NULL;
+ m_background = (rsp.m_background != NULL) ? rsp.m_background->clone() : NULL;
// Return
return;
@@ -2397,14 +2419,16 @@ void GCTAResponse::copy_members(const GCTAResponse& rsp)
void GCTAResponse::free_members(void)
{
// Free memory
- if (m_aeff != NULL) delete m_aeff;
- if (m_psf != NULL) delete m_psf;
- if (m_edisp != NULL) delete m_edisp;
+ if (m_aeff != NULL) delete m_aeff;
+ if (m_psf != NULL) delete m_psf;
+ if (m_edisp != NULL) delete m_edisp;
+ if (m_background != NULL) delete m_background;
// Initialise pointers
- m_aeff = NULL;
- m_psf = NULL;
- m_edisp = NULL;
+ m_aeff = NULL;
+ m_psf = NULL;
+ m_edisp = NULL;
+ m_background = NULL;
// Return
return;
diff --git a/inst/cta/test/data/irf_binned.xml b/inst/cta/test/data/irf_binned.xml
index 3a6d749..a71f8ef 100644
--- a/inst/cta/test/data/irf_binned.xml
+++ b/inst/cta/test/data/irf_binned.xml
@@ -1,8 +1,9 @@
-
-
-
+
+
+
+
diff --git a/inst/cta/test/data/irf_unbinned.xml b/inst/cta/test/data/irf_unbinned.xml
index 13dcb3c..6d7b3a5 100644
--- a/inst/cta/test/data/irf_unbinned.xml
+++ b/inst/cta/test/data/irf_unbinned.xml
@@ -1,8 +1,9 @@
-
-
-
+
+
+
+
diff --git a/inst/cta/test/data/obs_binned.xml b/inst/cta/test/data/obs_binned.xml
index 1057e7e..b4be423 100644
--- a/inst/cta/test/data/obs_binned.xml
+++ b/inst/cta/test/data/obs_binned.xml
@@ -4,5 +4,6 @@
+
diff --git a/inst/cta/test/data/obs_unbinned.xml b/inst/cta/test/data/obs_unbinned.xml
index 63564a9..d5d07be 100644
--- a/inst/cta/test/data/obs_unbinned.xml
+++ b/inst/cta/test/data/obs_unbinned.xml
@@ -4,5 +4,6 @@
+