From bed2fc7f812ebdf4ac41f604652b0b16ed52b62e Mon Sep 17 00:00:00 2001 From: Jürgen Knödlseder Date: Fri, 28 Jun 2013 12:46:23 +0200 Subject: [PATCH] Add CTA on-off observation handling classes. --- inst/cta/Makefile.am | 4 ++++ inst/cta/include/GCTAOnOffObservation.hpp | 102 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ inst/cta/include/GCTAOnOffObservations.hpp | 209 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ inst/cta/pyext/GCTAOnOffObservation.i | 76 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ inst/cta/pyext/GCTAOnOffObservations.i | 97 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ inst/cta/pyext/cta.i | 2 ++ inst/cta/src/GCTAOnOffObservation.cpp | 642 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ inst/cta/src/GCTAOnOffObservations.cpp | 832 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 8 files changed, 1964 insertions(+), 0 deletions(-) create mode 100644 inst/cta/include/GCTAOnOffObservation.hpp create mode 100644 inst/cta/include/GCTAOnOffObservations.hpp create mode 100644 inst/cta/pyext/GCTAOnOffObservation.i create mode 100644 inst/cta/pyext/GCTAOnOffObservations.i create mode 100644 inst/cta/src/GCTAOnOffObservation.cpp create mode 100644 inst/cta/src/GCTAOnOffObservations.cpp diff --git a/inst/cta/Makefile.am b/inst/cta/Makefile.am index d4b7524..6776831 100644 --- a/inst/cta/Makefile.am +++ b/inst/cta/Makefile.am @@ -11,6 +11,8 @@ INCLUDES = -I$(top_srcdir)/include \ sources = src/GCTASupport.cpp \ src/GCTAException.cpp \ src/GCTAObservation.cpp \ + src/GCTAOnOffObservation.cpp \ + src/GCTAOnOffObservations.cpp \ src/GCTAEventList.cpp \ src/GCTAEventAtom.cpp \ src/GCTAEventCube.cpp \ @@ -41,6 +43,8 @@ sources = src/GCTASupport.cpp \ # Define headers to be installed pkginclude_HEADERS = include/GCTAException.hpp \ include/GCTAObservation.hpp \ + include/GCTAOnOffObservation.hpp \ + include/GCTAOnOffObservations.hpp \ include/GCTAEventList.hpp \ include/GCTAEventAtom.hpp \ include/GCTAEventCube.hpp \ diff --git a/inst/cta/include/GCTAOnOffObservation.hpp b/inst/cta/include/GCTAOnOffObservation.hpp new file mode 100644 index 0000000..0b58f5f --- /dev/null +++ b/inst/cta/include/GCTAOnOffObservation.hpp @@ -0,0 +1,102 @@ +/*************************************************************************** + * GCTAOnOffObservation.hpp - CTA on-off observation class * + * ----------------------------------------------------------------------- * + * copyright (C) 2013 by Chia-Chun Lu & Christoph Deil * + * ----------------------------------------------------------------------- * + * * + * 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 * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see . * + * * + ***************************************************************************/ +/** + * @file GCTAOnOffObservation.hpp + * @brief CTA on-off observation class definition + * @author Chia-Chun Lu & Christoph Deil + */ + +#ifndef GCTAONOFFOBSERVATION_HPP +#define GCTAONOFFOBSERVATION_HPP + +/* __ Includes ___________________________________________________________ */ +#include +#include "GBase.hpp" +#include "GPha.hpp" +#include "GArf.hpp" +#include "GRmf.hpp" +#include "GCTAEventList.hpp" +#include "GCTAEventAtom.hpp" +#include "GCTAObservation.hpp" +#include "GSkyRegions.hpp" + + +/***********************************************************************//** + * @class GCTAOnOffObservation + * + * @brief CTA on-off observation class + ***************************************************************************/ +class GCTAOnOffObservation : public GBase { + +public: + // Constructors and destructors + GCTAOnOffObservation(void); + GCTAOnOffObservation(const GEbounds& ereco, const GSkyRegions& on, + const GSkyRegions& off); + GCTAOnOffObservation(const GCTAOnOffObservation& obs); + virtual ~GCTAOnOffObservation(void); + + // Operators + GCTAOnOffObservation& operator=(const GCTAOnOffObservation& obs); + + // Methods + void clear(void); + GCTAOnOffObservation* clone(void) const; + void name(const std::string& name) { m_name = name; } + void instrument(const std::string& instrument) { m_instrument = instrument; } + void id(const std::string& id) { m_id = id; } + void off_regions(const GSkyRegions& regions) {m_off_regions = regions;} + void on_regions(const GSkyRegions& regions) {m_on_regions = regions;} + const std::string& name(void) const { return m_name; } + const std::string& instrument(void) const { return m_instrument; } + const std::string& id(void) const { return m_id; } + const GPha& on_spec(void) const { return m_on_spec; } + const GPha& off_spec(void) const { return m_off_spec; } + const GArf& arf(void) { return m_arf; }; + const GRmf& rmf(void) {return m_rmf;} + void fill(const GCTAObservation& obs); + void compute_response(const GCTAObservation& obs, + const GEbounds& etrue); + void read(const GXmlElement& xml); + void write(GXmlElement& xml) const; + std::string print(const GChatter& chatter = NORMAL) const; + +protected: + // Protected methods + void init_members(void); + void copy_members(const GCTAOnOffObservation& obs); + void free_members(void); + void compute_arf(const GCTAObservation& obs); + void compute_rmf(const GCTAObservation& obs, const GEbounds& etrue); + + // Protected data members + std::string m_name; //!< Name + std::string m_instrument; //!< Instrument name + std::string m_id; //!< Observation identifier + GPha m_on_spec; + GPha m_off_spec; + GArf m_arf; + GRmf m_rmf; + GSkyRegions m_on_regions; + GSkyRegions m_off_regions; +}; + +#endif /* GCTAONOFFOBSERVATION_HPP */ diff --git a/inst/cta/include/GCTAOnOffObservations.hpp b/inst/cta/include/GCTAOnOffObservations.hpp new file mode 100644 index 0000000..57ef8ac --- /dev/null +++ b/inst/cta/include/GCTAOnOffObservations.hpp @@ -0,0 +1,209 @@ +/*************************************************************************** + * GCTAOnOffObservations.hpp - CTA on-off observation container class * + * ----------------------------------------------------------------------- * + * copyright (C) 2013 by Pierrick Martin * + * ----------------------------------------------------------------------- * + * * + * 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 * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see . * + * * + ***************************************************************************/ +/** + * @file GCTAOnOffObservations.hpp + * @brief CTA on-off observation container class definition + * @author Pierrick Martin + */ + +#ifndef GCTAONOFFOBSERVATIONS_HPP +#define GCTAONOFFOBSERVATIONS_HPP + +/* __ Includes ___________________________________________________________ */ +#include +#include +#include "GContainer.hpp" +#include "GCTAOnOffObservation.hpp" +#include "GModels.hpp" + + +/***********************************************************************//** + * @class GCTAOnOffObservations + * + * @brief ON/OFF Observation container class + * + * This class in a container of GCTAOnOffObservation objects. Still some + * work to do and things to be clarified... + ***************************************************************************/ +class GCTAOnOffObservations : public GContainer { + +public: + // Constructors and destructors + GCTAOnOffObservations(void); + GCTAOnOffObservations(const GCTAOnOffObservations& obs); + explicit GCTAOnOffObservations(const std::string& filename); + virtual ~GCTAOnOffObservations(void); + + // Operators + GCTAOnOffObservations& operator=(const GCTAOnOffObservations& obs); + GCTAOnOffObservation* operator[](const int& index); + const GCTAOnOffObservation* operator[](const int& index) const; + + // Methods + void clear(void); + GCTAOnOffObservations* clone(void) const; + GCTAOnOffObservation* at(const int& index); + const GCTAOnOffObservation* at(const int& index) const; + int size(void) const; + bool isempty(void) const; + GCTAOnOffObservation* set(const int& index, const GCTAOnOffObservation& obs); + GCTAOnOffObservation* append(const GCTAOnOffObservation& obs); + GCTAOnOffObservation* insert(const int& index, const GCTAOnOffObservation& obs); + void remove(const int& index); + void reserve(const int& num); + void extend(const GCTAOnOffObservations& obs); + bool contains(const std::string& instrument, + const std::string& id) const; + void load(const std::string& filename); + void save(const std::string& filename) const; + void read(const GXml& xml); + void write(GXml& xml) const; + void models(const GModels& models); + void models(const std::string& filename); + const GModels& models(void) const; + std::string print(const GChatter& chatter = NORMAL) const; + + /* + // To do + Implement fit functions for both cases summed over runs or not ? + Implement significance or have it at ctool level ? + Implement flux points computation and fit of a spectral function ? + // Do we want that ? If yes, need to retrieve inline functions from GObservations... + double npred(void) const; + // Left aside at the moment + void optimize(GOptimizer& opt); + */ + +protected: + // Protected methods + void init_members(void); + void copy_members(const GCTAOnOffObservations& obs); + void free_members(void); + int get_index(const std::string& instrument, + const std::string& id) const; + + // Protected members + std::vector m_obs; //!< List of observations + GModels m_models; //!< List of models +}; + + +/***********************************************************************//** + * @brief Return pointer to observation + * + * @param[in] index Observation index [0,...,size()-1]. + * @return Observation. + * + * Returns a pointer to the observation with specified @p index. + ***************************************************************************/ +inline +GCTAOnOffObservation* GCTAOnOffObservations::operator[](const int& index) +{ + return m_obs[index]; +} + + +/***********************************************************************//** + * @brief Return pointer to observation (const version) + * + * @param[in] index Observation index [0,...,size()-1]. + * + * Returns a const pointer to the observation with specified @p index. + ***************************************************************************/ +inline +const GCTAOnOffObservation* GCTAOnOffObservations::operator[](const int& index) const +{ + return m_obs[index]; +} + + +/***********************************************************************//** + * @brief Return number of observations in container + * + * @return Number of observations in container. + * + * Returns the number of observations in the observation container. + ***************************************************************************/ +inline +int GCTAOnOffObservations::size(void) const +{ + return (m_obs.size()); +} + + +/***********************************************************************//** + * @brief Signals if there are no observations in container + * + * @return True if container is empty, false otherwise. + * + * Signals if the observation container does not contain any observation. + ***************************************************************************/ +inline +bool GCTAOnOffObservations::isempty(void) const +{ + return (m_obs.empty()); +} + + +/***********************************************************************//** + * @brief Reserves space for observations in container + * + * @param[in] num Number of observations. + * + * Reserves space for @p num observations in the container. + ***************************************************************************/ +inline +void GCTAOnOffObservations::reserve(const int& num) +{ + m_obs.reserve(num); + return; +} + + +/***********************************************************************//** + * @brief Set model container + * + * @param[in] models Model container. + * + * Sets the model container for the observations. + ***************************************************************************/ +inline +void GCTAOnOffObservations::models(const GModels& models) +{ + m_models = models; + return; +} + + +/***********************************************************************//** + * @brief Return model container + * + * @return Model container. + * + * Returns the model container of the observations. + ***************************************************************************/ +inline +const GModels& GCTAOnOffObservations::models(void) const +{ + return m_models; +} + +#endif /* GCTAONOFFOBSERVATIONS_HPP */ diff --git a/inst/cta/pyext/GCTAOnOffObservation.i b/inst/cta/pyext/GCTAOnOffObservation.i new file mode 100644 index 0000000..1b92d3d --- /dev/null +++ b/inst/cta/pyext/GCTAOnOffObservation.i @@ -0,0 +1,76 @@ +/*************************************************************************** + * GCTAOnOffObservation.i - CTA on-off observation class * + * ----------------------------------------------------------------------- * + * copyright (C) 2013 by Michael 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 * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see . * + * * + ***************************************************************************/ +/** + * @file GCTAOnOffObservation.i + * @brief CTA on-off observation class definition + * @author Michael Mayer + */ +%{ +/* Put headers and other declarations here that are needed for compilation */ +#include "GCTAOnOffObservation.hpp" +%} + + +/***********************************************************************//** + * @class GCTAOnOffObservation + * + * @brief CTA on-off observation class + ***************************************************************************/ +class GCTAOnOffObservation : public GBase { +public: + // Constructors and destructors + GCTAOnOffObservation(void); + GCTAOnOffObservation(const GEbounds& ereco, const GSkyRegions& on, + const GSkyRegions& off); + GCTAOnOffObservation(const GCTAOnOffObservation& obs); + virtual ~GCTAOnOffObservation(void); + + // Methods + void clear(void); + GCTAOnOffObservation* clone(void) const; + void name(const std::string& name) { m_name = name; } + void instrument(const std::string& instrument) { m_instrument = instrument; } + void id(const std::string& id) { m_id = id; } + void off_regions(const GSkyRegions& regions) {m_off_regions = regions;} + void on_regions(const GSkyRegions& regions) {m_on_regions = regions;} + const std::string& name(void) const { return m_name; } + const std::string& instrument(void) const { return m_instrument; } + const std::string& id(void) const { return m_id; } + const GPha& on_spec(void) const { return m_on_spec; } + const GPha& off_spec(void) const { return m_off_spec; } + const GArf& arf(void) { return m_arf; }; + const GRmf& rmf(void) {return m_rmf;} + void fill(const GCTAObservation& obs); + void compute_response(const GCTAObservation& obs, + const GEbounds& etrue); + void read(const GXmlElement& xml); + void write(GXmlElement& xml) const; +}; + + +/***********************************************************************//** + * @brief GCTOnOffAObservation class extension + ***************************************************************************/ +%extend GCTAOnOffObservation { + GCTAOnOffObservation copy() { + return (*self); + } +} diff --git a/inst/cta/pyext/GCTAOnOffObservations.i b/inst/cta/pyext/GCTAOnOffObservations.i new file mode 100644 index 0000000..fc7cc8d --- /dev/null +++ b/inst/cta/pyext/GCTAOnOffObservations.i @@ -0,0 +1,97 @@ +/*************************************************************************** + * GCTAOnOffObservations.i - CTA on-off observation container class * + * ----------------------------------------------------------------------- * + * copyright (C) 2013 by Pierrick Martin * + * ----------------------------------------------------------------------- * + * * + * 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 * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see . * + * * + ***************************************************************************/ +/** + * @file GCTAOnOffObservations.i + * @brief CTA on-off observation container class definition + * @author Pierrick Martin + */ +%{ +/* Put headers and other declarations here that are needed for compilation */ +#include "GCTAOnOffObservations.hpp" +%} + + +/***********************************************************************//** + * @class GCTAOnOffObservations + * + * @brief ON/OFF Observation container class + ***************************************************************************/ +class GCTAOnOffObservations : public GContainer { +public: + // Constructors and destructors + GCTAOnOffObservations(void); + GCTAOnOffObservations(const GCTAOnOffObservations& obs); + explicit GCTAOnOffObservations(const std::string& filename); + virtual ~GCTAOnOffObservations(void); + + // Methods + void clear(void); + GCTAOnOffObservations* clone(void) const; + GCTAOnOffObservation* at(const int& index); + const GCTAOnOffObservation* at(const int& index) const; + int size(void) const; + bool isempty(void) const; + GCTAOnOffObservation* set(const int& index, const GCTAOnOffObservation& obs); + GCTAOnOffObservation* append(const GCTAOnOffObservation& obs); + GCTAOnOffObservation* insert(const int& index, const GCTAOnOffObservation& obs); + void remove(const int& index); + void reserve(const int& num); + void extend(const GCTAOnOffObservations& obs); + bool contains(const std::string& instrument, + const std::string& id) const; + void load(const std::string& filename); + void save(const std::string& filename) const; + void read(const GXml& xml); + void write(GXml& xml) const; + void models(const GModels& models); + void models(const std::string& filename); + const GModels& models(void) const; +}; + + +/***********************************************************************//** + * @brief GCTAOnOffObservations class extension + ***************************************************************************/ +%extend GCTAOnOffObservations { + GCTAOnOffObservation* __getitem__(const int& index) { + if (index >= 0 && index < self->size()) { + return (*self)[index]; + } + else { + throw GException::out_of_range("__getitem__(int)", index, self->size()); + } + } + void __setitem__(const int& index, const GCTAOnOffObservation& val) { + if (index>=0 && index < self->size()) { + self->set(index, val); + return; + } + else { + throw GException::out_of_range("__setitem__(int)", index, self->size()); + } + } + int __len__() { + return (self->size()); + } + GCTAOnOffObservations copy() { + return (*self); + } +}; diff --git a/inst/cta/pyext/cta.i b/inst/cta/pyext/cta.i index 0fdafe5..ce16010 100644 --- a/inst/cta/pyext/cta.i +++ b/inst/cta/pyext/cta.i @@ -63,6 +63,8 @@ /* __ CTA ________________________________________________________________ */ %include "GCTAObservation.i" +%include "GCTAOnOffObservation.i" +%include "GCTAOnOffObservations.i" %include "GCTAEventCube.i" %include "GCTAEventList.i" %include "GCTAEventBin.i" diff --git a/inst/cta/src/GCTAOnOffObservation.cpp b/inst/cta/src/GCTAOnOffObservation.cpp new file mode 100644 index 0000000..cc5ff5c --- /dev/null +++ b/inst/cta/src/GCTAOnOffObservation.cpp @@ -0,0 +1,642 @@ +/*************************************************************************** + * GCTAOnOffObservation.cpp - CTA on-off observation class * + * ----------------------------------------------------------------------- * + * copyright (C) 2013 by Chia-Chun Lu & Christoph Deil * + * ----------------------------------------------------------------------- * + * * + * 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 * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see . * + * * + ***************************************************************************/ +/** + * @file GCTAOnOffObservation.cpp + * @brief CTA on-off observation class implementation + * @author Chia-Chun Lu & Christoph Deil + */ + +/* __ Includes ___________________________________________________________ */ +#ifdef HAVE_CONFIG_H +#include +#endif +#include "GCTAOnOffObservation.hpp" +#include "GTools.hpp" + +/* __ Method name definitions ____________________________________________ */ +#define G_WRITE "GCTAOnOffObservation::write(GXmlElement&)" +#define G_READ "GCTAOnOffObservation::read(GXmlElement&)" +#define G_FILL "GCTAOnOffObservation::fill(GCTAObservation&)" + +/* __ Macros _____________________________________________________________ */ + +/* __ Coding definitions _________________________________________________ */ + +/* __ Debug definitions __________________________________________________ */ + + +/*========================================================================== + = = + = Constructors/destructors = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Void constructor + ***************************************************************************/ +GCTAOnOffObservation::GCTAOnOffObservation(void) +{ + // Initialise private members + init_members(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Copy constructor + * + * @param[in] obs CTA on-off observation. + ***************************************************************************/ +GCTAOnOffObservation::GCTAOnOffObservation(const GCTAOnOffObservation& obs) +{ + // Initialise private + init_members(); + + // Copy members + copy_members(obs); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Parameter constructor + * + * @param[in] ereco Reconstructed energy bins. + * @param[in] on On regions. + * @param[in] off Off regions. + ***************************************************************************/ +GCTAOnOffObservation::GCTAOnOffObservation(const GEbounds& ereco, + const GSkyRegions& on, + const GSkyRegions& off) +{ + // Initialise private + init_members(); + + // Initialise spectra + m_on_spec = GPha(ereco); + m_off_spec = GPha(ereco); + + // Store regions + m_on_regions = on; + m_off_regions = off; + + // Return + return; +} + + +/***********************************************************************//** + * @brief Destructor + ***************************************************************************/ +GCTAOnOffObservation::~GCTAOnOffObservation(void) +{ + // Free members + free_members(); + + // Return + return; +} + + +/*========================================================================== + = = + = Operators = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Assignment operator + * + * @param[in] obs CTA on-off observation. + ***************************************************************************/ +GCTAOnOffObservation& GCTAOnOffObservation::operator=(const GCTAOnOffObservation& obs) +{ + // Execute only if object is not identical + if (this != &obs) { + + // Free members + free_members(); + + // Initialise private members for clean destruction + init_members(); + + // Copy members + copy_members(obs); + + } // endif: object was not identical + + // Return + return *this; +} + + +/*========================================================================== + = = + = Public methods = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Clear instance + * + * This method properly resets the object to an initial state. + ***************************************************************************/ +void GCTAOnOffObservation::clear(void) +{ + // Free class members + free_members(); + + // Initialise members + init_members(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Clone instance + **************************************************************************/ +GCTAOnOffObservation* GCTAOnOffObservation::clone(void) const +{ + return new GCTAOnOffObservation(*this); +} + + +/***********************************************************************//** + * @brief read observation from an xml element + * + * @param[in] xml GXmlElement. + ***************************************************************************/ +void GCTAOnOffObservation::read(const GXmlElement& xml) +{ + // clean object + clear(); + + // Extract instrument name + m_instrument = xml.attribute("instrument"); + + // Determine number of parameter nodes in XML element + int npars = xml.elements("parameter"); + + // Verify that XML element has exactly 6 parameters + if (xml.elements() != 6 || npars != 6) { + throw GException::xml_invalid_parnum(G_READ, xml, + "CTA observation requires exactly 6 parameters."); + } + + // Extract parameters + int npar[] = {0, 0, 0, 0, 0, 0}; + for (int i = 0; i < npars; ++i) { + + // Get parameter element + const GXmlElement* par = xml.element("parameter", i); + + // Handle Pha_on + if (par->attribute("name") == "Pha_on") { + + // Read pha_on file name + std::string filename = par->attribute("file"); + + // Load pha_on + m_on_spec.load(filename); + + // Increment parameter counter + npar[0]++; + } + + // Handle Pha_off + else if (par->attribute("name") == "Pha_off") { + + // Read countsmap file name + std::string filename = par->attribute("file"); + + // Load pha_off + m_off_spec.load(filename); + + // Increment parameter counter + npar[1]++; + } + + // Handle on regions + else if(par->attribute("name") == "Regions_on") { + + // Get filename + std::string filename = par->attribute("file"); + + // load on regions + m_on_regions.load(filename); + + // Increase number of parameters + npar[2]++; + } + + // Handle off regions + else if(par->attribute("name") == "Regions_off") { + + // Get filename + std::string filename = par->attribute("file"); + + // load off regions + m_on_regions.load(filename); + + // Increase number of parameters + npar[3]++; + } + + // Handle Arf + else if(par->attribute("name") == "Arf") { + + // Get filename + std::string filename = par->attribute("file"); + + // load arf + m_arf.load(filename); + + // Increase number of parameters + npar[4]++; + } + + // Handle Rmf + else if(par->attribute("name") == "Rmf") { + + // Get filename + std::string filename = par->attribute("file"); + + // load Rmf + m_rmf.load(filename); + + // Increase number of parameters + npar[5]++; + } + + } // endfor: looped over all parameters + + // Verify that all parameters were found + if (npar[0] != 1 || npar[1] != 1 || npar[2] != 1 || npar[3] != 1 || npar[4] != 1 || npar[5] != 1) { + throw GException::xml_invalid_parnames(G_READ, xml, + "Require \"Pha_on\" or \"Pha_off\" and \"Regions_on\"" + ", \"Regions_off\",\"Rmf\" and \"Arf\" parameters."); + } + + // Return + return; +} + + +/***********************************************************************//** + * @brief write observation to an xml element + * + * @param[in] xml GXmlElement. + ***************************************************************************/ +void GCTAOnOffObservation::write(GXmlElement& xml) const +{ + // If XML element has 0 nodes then append 6 parameter nodes + if (xml.elements() == 0) { + xml.append(GXmlElement("parameter name=\"Pha_on\"")); + xml.append(GXmlElement("parameter name=\"Pha_off\"")); + xml.append(GXmlElement("parameter name=\"Regions_on\"")); + xml.append(GXmlElement("parameter name=\"Regions_off\"")); + xml.append(GXmlElement("parameter name=\"Arf\"")); + xml.append(GXmlElement("parameter name=\"Rmf\"")); + } + + // Verify that XML element has exactly 4 parameters + if (xml.elements() != 6 || xml.elements("parameter") != 6) { + throw GException::xml_invalid_parnum(G_WRITE, xml, + "CTAOnOffObservation requires exactly 6 parameters."); + } + + // Set or update parameter attributes + int npar[] = {0, 0, 0, 0, 0, 0}; + for (int i = 0; i < 6; ++i) { + + // Get parameter element + GXmlElement* par = xml.element("parameter", i); + + // Handle on counts + if (par->attribute("name") == "Pha_on") { + par->attribute("file", m_on_spec.filename()); + npar[0]++; + } + + // Handle off counts + else if (par->attribute("name") == "Pha_off") { + par->attribute("file", m_off_spec.filename()); + npar[0]++; + } + + // Handle on regions + else if (par->attribute("name") == "Regions_on") { + par->attribute("file", m_on_regions.filename()); + npar[0]++; + } + + // Handle off regions + else if (par->attribute("name") == "Regions_off") { + par->attribute("file", m_off_regions.filename()); + npar[0]++; + } + + // Handle effective area + else if (par->attribute("name") == "Arf") { + par->attribute("file", m_arf.filename()); + npar[0]++; + } + + // Handle energy resolution + else if (par->attribute("name") == "Rmf") { + par->attribute("file", m_rmf.filename()); + npar[0]++; + } + } + // Verify that all required parameters are present + if (npar[0] != 1 || npar[1] != 1 || npar[2] != 1 || npar[3] != 1 || npar[4] != 1 || npar[5] !=1) { + throw GException::xml_invalid_parnames(G_WRITE, xml, + "Require \"Pha_on\" or \"Pha_off\" and \"Regions_on\", \"Regions_off\"" + " \"Arf\", \"Rmf\" parameters."); + } + + // Return + return; +} + + +/***********************************************************************//** + * @brief Fill events in the on and off spectra + * + * @param[in] obs CTA observation. + * + * @exception GException::invalid_value + * No CTA event list found in CTA observation. + ***************************************************************************/ +void GCTAOnOffObservation::fill(const GCTAObservation& obs) +{ + // Get CTA event list pointer + const GCTAEventList* events = dynamic_cast(obs.events()); + if (events == NULL) { + throw GException::invalid_value(G_FILL, + "GCTAObservation must provide an event list"); + } + + // Loop over all events + for (int i = 0; i < events->size(); ++i) { + + // Get measured event direction + const GCTAEventAtom* atom = (*events)[i]; + GSkyDir dir = atom->dir().dir(); + + // Fill in spectrum according to region containment + if (m_on_regions.contains(dir)) { + m_on_spec.fill(atom->energy()); + } + else if (m_off_regions.contains(dir)) { + m_off_spec.fill(atom->energy()); + } + + } // endfor: looped over all events + + // Return + return; +} + + +/***********************************************************************//** + * @brief Compute response of an on/off observation + * + * @param[in] obs GCTAObservation. + * @param[in] etrue True energy boundaries. + ***************************************************************************/ +void GCTAOnOffObservation::compute_response(const GCTAObservation& obs, + const GEbounds& etrue) +{ + // Compute response components + compute_arf(obs); + compute_rmf(obs, etrue); + + // Return + return; +} + +/***********************************************************************//** + * @brief Compute arf of an on/off observation + * + * @param[in] obs GCTAObservation. + * + * @todo Implement GCTAResponse::npred usage. + ***************************************************************************/ +void GCTAOnOffObservation::compute_arf(const GCTAObservation& obs) +{ + // Set constant response parameters + const double theta = 0.0; + const double phi = 0.0; + const double zenith = 0.0; + const double azimuth = 0.0; + + // Get energy boundaries from on spectrum + GEbounds ereco = m_on_spec.ebounds(); + + // Continue only if there are spectral bins + int nreco = ereco.size(); + if (nreco > 0) { + + // Get CTA response pointer + GCTAResponse* response = obs.response(); + + // Initialize ARF + m_arf = GArf(ereco); + + // Loop over reconstructed energies + for (int i = 0; i < nreco; ++i) { + + // Get mean energy of bin + double logenergy = ereco.elogmean(i).log10TeV(); + + // Set specresp value + m_arf[i] = response->aeff(theta, + phi, + zenith, + azimuth, + logenergy); + + } // endfor: looped over reconstructed energies + + } // endif: there were energy bins + + // Return + return; +} + + +/***********************************************************************//** + * @brief compute rmf of an on/off observation + * + * @param[in] obs GCTAObservation. + * @param[in] etrue True energy boundaries. + ***************************************************************************/ +void GCTAOnOffObservation::compute_rmf(const GCTAObservation& obs, + const GEbounds& etrue) +{ + // Set constant response parameters + const double theta = 0.0; + const double phi = 0.0; + const double zenith = 0.0; + const double azimuth = 0.0; + + // Get reconstructed energy boundaries from on spectrum + GEbounds ereco = m_on_spec.ebounds(); + + // Continue only if there are spectral bins + int ntrue = etrue.size(); + int nreco = ereco.size(); + if (ntrue > 0 && nreco > 0) { + + // Get CTA response pointer + GCTAResponse* response = obs.response(); + + // Initialize RMF + m_rmf = GRmf(etrue, ereco); + + // Loop over reconstructed energy + for (int ireco = 0; ireco < nreco; ++ireco) { + + // Compute reconstructed energy + double eng_reco = ereco.elogmean(ireco).log10TeV(); + + // Loop over true energy + for (int itrue = 0; itrue < ntrue; ++itrue) { + + // Compute true energy + double eng_true = etrue.elogmean(itrue).log10TeV(); + + // Set RMF value + m_rmf(itrue, ireco) = response->edisp(eng_reco, + theta, + phi, + zenith, + azimuth, + eng_true); + } // endfor: looped over true energy + } // endfor: looped over reconstructed energy + } // endif: there were energy bins + + // Return + return; +} + + +/***********************************************************************//** + * @brief Print class information + ***************************************************************************/ +std::string GCTAOnOffObservation::print(const GChatter& chatter) const +{ + // Initialise result string + std::string result; + + // Continue only if chatter is not silent + if (chatter != SILENT) { + + // Append header + result.append("=== GCTAOnOffObservation ==="); + + // Append parameters + result.append("\n"+gammalib::parformat("Name")+m_name); + result.append("\n"+gammalib::parformat("Identifier")+m_id); + + // Append spectra, ARF and RMF + result.append("\n"+m_on_spec.print()); + result.append("\n"+m_off_spec.print()); + result.append("\n"+m_arf.print()); + result.append("\n"+m_rmf.print()); + + // Append regions + result.append("\n"+m_on_regions.print()); + result.append("\n"+m_off_regions.print()); + } + + // Return result + return result; +} + + +/*========================================================================== + = = + = Private methods = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Initialise class members + ***************************************************************************/ +void GCTAOnOffObservation::init_members(void) +{ + // Initialise members + m_name.clear(); + m_id.clear(); + m_instrument = "CTA"; + m_on_spec.clear(); + m_off_spec.clear(); + m_arf.clear(); + m_rmf.clear(); + m_on_regions.clear(); + m_off_regions.clear(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Copy class members + * + * @param[in] obs CTA on-off observation. + ***************************************************************************/ +void GCTAOnOffObservation::copy_members(const GCTAOnOffObservation& obs) +{ + // Copy attributes + m_name = obs.m_name; + m_id = obs.m_id; + m_instrument = obs.m_instrument; + m_on_spec = obs.m_on_spec; + m_off_spec = obs.m_off_spec; + m_arf = obs.m_arf; + m_rmf = obs.m_rmf; + m_on_regions = obs.m_on_regions; + m_off_regions = obs.m_off_regions; + + // Return + return; +} + + + +/***********************************************************************//** + * @brief Delete class members + ***************************************************************************/ +void GCTAOnOffObservation::free_members(void) +{ + // Return + return; +} diff --git a/inst/cta/src/GCTAOnOffObservations.cpp b/inst/cta/src/GCTAOnOffObservations.cpp new file mode 100644 index 0000000..22023a4 --- /dev/null +++ b/inst/cta/src/GCTAOnOffObservations.cpp @@ -0,0 +1,832 @@ +/*************************************************************************** + * GCTAOnOffObservations.cpp - ON/OFF Observation container class * + * ----------------------------------------------------------------------- * + * copyright (C) 2013 by Pierrick Martin * + * ----------------------------------------------------------------------- * + * * + * 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 * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see . * + * * + ***************************************************************************/ +/** + * @file GCTAOnOffObservations.cpp + * @brief ON/OFF Observations container class implementation + * @author Pierrick Martin + */ + +/* __ Includes ___________________________________________________________ */ +#ifdef HAVE_CONFIG_H +#include +#endif +#include "GTools.hpp" +#include "GException.hpp" +#include "GCTAOnOffObservations.hpp" + +/* __ Method name definitions ____________________________________________ */ +#define G_AT "GCTAOnOffObservations::at(int&)" +#define G_SET "GCTAOnOffObservations::set(int&, GCTAOnOffObservation&)" +#define G_APPEND "GCTAOnOffObservations::append(GCTAOnOffObservation&)" +#define G_INSERT "GCTAOnOffObservations::insert(int&, GCTAOnOffObservation&)" +#define G_REMOVE "GCTAOnOffObservations::remove(int&)" +#define G_EXTEND "GCTAOnOffObservations::extend(GCTAOnOffObservations&)" +#define G_READ "GCTAOnOffObservations::read(GXml&)" + +/* __ Macros _____________________________________________________________ */ + +/* __ Coding definitions _________________________________________________ */ + +/* __ Debug definitions __________________________________________________ */ + + +/*========================================================================== + = = + = Constructors/destructors = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Void constructor + ***************************************************************************/ +GCTAOnOffObservations::GCTAOnOffObservations(void) +{ + // Initialise members + init_members(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Copy constructor + * + * @param obs Observation container. + ***************************************************************************/ +GCTAOnOffObservations::GCTAOnOffObservations(const GCTAOnOffObservations& obs) +{ + // Initialise members + init_members(); + + // Copy members + copy_members(obs); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Load constructor + * + * @param[in] filename XML filename. + * + * Construct the observation container by loading all relevant information + * from a XML file. Please refer to the read() method for more information + * about the structure of the XML file. + ***************************************************************************/ +GCTAOnOffObservations::GCTAOnOffObservations(const std::string& filename) +{ + // Initialise members + init_members(); + + // Load XML file + load(filename); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Destructor + ***************************************************************************/ +GCTAOnOffObservations::~GCTAOnOffObservations(void) +{ + // Free members + free_members(); + + // Return + return; +} + + +/*========================================================================== + = = + = Operators = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Assignment operator + * + * @param[in] obs Observation container. + * @return Observation container. + ***************************************************************************/ +GCTAOnOffObservations& GCTAOnOffObservations::operator=(const GCTAOnOffObservations& obs) +{ + // Execute only if object is not identical + if (this != &obs) { + + // Free members + free_members(); + + // Initialise members + init_members(); + + // Copy members + copy_members(obs); + + } // endif: object was not identical + + // Return this object + return *this; +} + + +/*========================================================================== + = = + = Public methods = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Clear observations + ***************************************************************************/ +void GCTAOnOffObservations::clear(void) +{ + // Free members + free_members(); + + // Initialise members + init_members(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Clone observations + * + * @return Pointer to deep copy of observation container. + ***************************************************************************/ +GCTAOnOffObservations* GCTAOnOffObservations::clone(void) const +{ + // Clone observations + return new GCTAOnOffObservations(*this); +} + + +/***********************************************************************//** + * @brief Return pointer to observation + * + * @param[in] index Observation index [0,...,size()-1]. + * @return Observation. + * + * @exception GException::out_of_range + * Operation index is out of range. + * + * Returns a pointer to the observation with specified @p index. + ***************************************************************************/ +GCTAOnOffObservation* GCTAOnOffObservations::at(const int& index) +{ + // Raise exception if index is out of range + if (index < 0 || index >= size()) { + throw GException::out_of_range(G_AT, index, 0, size()-1); + } + + // Return pointer + return m_obs[index]; +} + + +/***********************************************************************//** + * @brief Return pointer to observation (const version) + * + * @param[in] index Observation index [0,...,size()-1]. + * + * @exception GException::out_of_range + * Operation index is out of range. + * + * Returns a const pointer to the observation with specified @p index. + ***************************************************************************/ +const GCTAOnOffObservation* GCTAOnOffObservations::at(const int& index) const +{ + // Raise exception if index is out of range + if (index < 0 || index >= size()) { + throw GException::out_of_range(G_AT, index, 0, size()-1); + } + + // Return pointer + return m_obs[index]; +} + + +/***********************************************************************//** + * @brief Set observation in container + * + * @param[in] index Observation index [0,...,size()-1]. + * @param[in] obs Observation. + * @return Pointer to deep copy of observation. + * + * @exception GException::out_of_range + * Observation index is out of range. + * @exception GException::invalid_value + * Observation with same instrument and identifier already + * exists in container. + * + * Set a deep copy and observation @p obs at the specified @p index in the + * container. + ***************************************************************************/ +GCTAOnOffObservation* GCTAOnOffObservations::set(const int& index, const GCTAOnOffObservation& obs) +{ + // Compile option: raise exception if index is out of range + #if defined(G_RANGE_CHECK) + if (index < 0 || index >= size()) { + throw GException::out_of_range(G_SET, index, 0, size()-1); + } + #endif + + // Raise an exception if an observation with specified instrument and + // identifier already exists + int inx = get_index(obs.instrument(), obs.id()); + if (inx != -1 && inx != index) { + std::string msg = + "Attempt to set \""+obs.instrument()+"\" observation with" + " identifier \""+obs.id()+"\" in observation container at" + " index "+gammalib::str(index)+", but an observation with the" + " same attributes exists already at index "+gammalib::str(inx)+ + " in the container.\n" + "Every observation for a given instrument in the observation" + " container needs a unique identifier."; + throw GException::invalid_value(G_SET, msg); + } + + // Delete existing observation + if (m_obs[index] != NULL) delete m_obs[index]; + + // Store pointer to a deep copy of the observation + m_obs[index] = obs.clone(); + + // Return pointer + return m_obs[index]; +} + + +/***********************************************************************//** + * @brief Append observation to container + * + * @param[in] obs Observation. + * @return Pointer to deep copy of observation. + * + * @exception GException::invalid_value + * Observation with same instrument and identifier already + * exists in container. + * + * Appends a deep copy of an observation to the container. + ***************************************************************************/ +GCTAOnOffObservation* GCTAOnOffObservations::append(const GCTAOnOffObservation& obs) +{ + // Raise an exception if an observation with specified instrument and + // identifier already exists + int inx = get_index(obs.instrument(), obs.id()); + if (inx != -1) { + std::string msg = + "Attempt to append \""+obs.instrument()+"\" observation with" + " identifier \""+obs.id()+"\" to observation container, but an" + " observation with the same attributes exists already at" + " index "+gammalib::str(inx)+" in the container.\n" + "Every observation for a given instrument in the observation" + " container needs a unique identifier."; + throw GException::invalid_value(G_APPEND, msg); + } + + // Clone observation + GCTAOnOffObservation* ptr = obs.clone(); + + // Append to list + m_obs.push_back(ptr); + + // Return pointer + return ptr; +} + + +/***********************************************************************//** + * @brief Insert observation into container + * + * @param[in] index Observation index [0,...,size()-1]. + * @param[in] obs Observation. + * @return Pointer to deep copy of observation. + * + * @exception GException::out_of_range + * Observation index is out of range. + * @exception GException::invalid_value + * Observation with same instrument and identifier already + * exists in container. + * + * Inserts a deep copy of an observation into the container before the + * observation with the specified @p index. + ***************************************************************************/ +GCTAOnOffObservation* GCTAOnOffObservations::insert(const int& index, const GCTAOnOffObservation& obs) +{ + // Compile option: raise exception if index is out of range + #if defined(G_RANGE_CHECK) + if (isempty()) { + if (index > 0) { + throw GException::out_of_range(G_INSERT, index, 0, size()-1); + } + } + else { + if (index < 0 || index >= size()) { + throw GException::out_of_range(G_INSERT, index, 0, size()-1); + } + } + #endif + + // Raise an exception if an observation with specified instrument and + // identifier already exists + int inx = get_index(obs.instrument(), obs.id()); + if (inx != -1 && inx != index) { + std::string msg = + "Attempt to insert \""+obs.instrument()+"\" observation with" + " identifier \""+obs.id()+"\" in observation container before" + " index "+gammalib::str(index)+", but an observation with the" + " same attributes exists already at index "+gammalib::str(inx)+ + " in the container.\n" + "Every observation for a given instrument in the observation" + " container needs a unique identifier."; + throw GException::invalid_value(G_INSERT, msg); + } + + // Clone observation + GCTAOnOffObservation* ptr = obs.clone(); + + // Clone observation and insert into list + m_obs.insert(m_obs.begin()+index, ptr); + + // Return pointer + return ptr; +} + + +/***********************************************************************//** + * @brief Remove observation from container + * + * @param[in] index Observation index [0,...,size()-1]. + * + * @exception GException::out_of_range + * Observation index is out of range. + * + * Removes observation of specified @p index from the container. + ***************************************************************************/ +void GCTAOnOffObservations::remove(const int& index) +{ + // Compile option: If index is outside boundary then raise exception + #if defined(G_RANGE_CHECK) + if (index < 0 || index >= size()) { + throw GException::out_of_range(G_REMOVE, index, 0, size()-1); + } + #endif + + // Erase observation component from container + m_obs.erase(m_obs.begin() + index); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Append observations from observation container + * + * @param[in] obs Observations. + * + * @exception GException::invalid_value + * Observation with same instrument and identifier already + * exists in container. + * + * Appends deep copies of observations to the container. + ***************************************************************************/ +void GCTAOnOffObservations::extend(const GCTAOnOffObservations& obs) +{ + // Do nothing if observation container is empty + if (!obs.isempty()) { + + // Get size. Note that we extract the size first to avoid an + // endless loop that arises when a container is appended to + // itself. + int num = obs.size(); + + // Reserve enough space + reserve(size() + num); + + // Loop over all observations, clone them and append them to the + // list + for (int i = 0; i < num; ++i) { + + // Raise an exception if an observation with specified + // instrument and identifier already exists + int inx = get_index(obs[i]->instrument(), obs[i]->id()); + if (inx != -1) { + std::string msg = + "Attempt to append \""+obs[i]->instrument()+"\"" + " observation with identifier \""+obs[i]->id()+"\"" + " to observation container, but an observation with the" + " same attributes exists already at index "+ + gammalib::str(inx)+" in the container.\n" + "Every observation for a given instrument in the" + " observation container needs a unique identifier."; + throw GException::invalid_value(G_EXTEND, msg); + } + + // Append observation to container + m_obs.push_back(obs[i]->clone()); + } + + } // endif: observation container was not empty + + // Return + return; +} + + +/***********************************************************************//** + * @brief Signals if observation exists + * + * @param[in] instrument Instrument. + * @param[in] id Observation identifier. + * @return True if observation with specified @p instrument and identifier + * @p id exists. + * + * Searches all observations for a match with the specified @p instrument + * and @p identifier. If the specified attributes have been found, true is + * returned. + ***************************************************************************/ +bool GCTAOnOffObservations::contains(const std::string& instrument, + const std::string& id) const +{ + // Get observation index + int index = get_index(instrument, id); + + // Return + return (index != -1); +} + + +/***********************************************************************//** + * @brief Load observations from XML file + * + * @param[in] filename XML filename. + * + * Loads observation from a XML file into the container. Please refer to the + * read() method for more information about the structure of the XML file. + ***************************************************************************/ +void GCTAOnOffObservations::load(const std::string& filename) +{ + // Clear any existing observations + clear(); + + // Load XML document + GXml xml(filename); + + // Read observations from XML document + read(xml); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Save observations into XML file + * + * @param[in] filename XML filename. + * + * Saves observations into a XML file. Please refer to the read() method for + * more information about the structure of the XML file. + ***************************************************************************/ +void GCTAOnOffObservations::save(const std::string& filename) const +{ + // Declare empty XML document + GXml xml; + + // Write observations into XML file + write(xml); + + // Save XML document + xml.save(filename); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Read observations from XML document + * + * @param[in] xml XML document. + * + * @exception GException::invalid_instrument + * Invalid instrument encountered in XML file. + * + * Reads observations from the first observation list that is found in the + * XML document. The decoding of the instrument specific observation + * definition is done within the observation's GCTAOnOffObservation::read() method. + * The following file structure is expected: + * + * + * + * ... + * + * + * ... + * + * ... + * + * + * The @p name and @p id attributes allow for a unique identification of an + * observation within the observation container. The @p instrument + * attributes specifies the instrument to which the observation applies. + * At the moment there is only CTA but more instruments may be implemented + * later and require the need for a registry (see GObservations.cpp). + * + * The structure within the @p observation tag is defined by the instrument + * specific GCTAOnOffObservation class. + * + ***************************************************************************/ +void GCTAOnOffObservations::read(const GXml& xml) +{ + // Get pointer on observation library + const GXmlElement* lib = xml.element("observation_list", 0); + + // Loop over all observations + int n = lib->elements("observation"); + for (int i = 0; i < n; ++i) { + + // Get pointer on observation + const GXmlElement* obs = lib->element("observation", i); + + // Get attributes + std::string name = obs->attribute("name"); + std::string id = obs->attribute("id"); + std::string instrument = obs->attribute("instrument"); + + // Allocate observation + // (only CTA at the moment, implement registry if more) + GCTAOnOffObservation* ptr = new GCTAOnOffObservation; + + // If observation is valid then read its definition from XML file + if (ptr != NULL) { + + // Read definition + ptr->read(*obs); + + // Set attributes + ptr->name(name); + ptr->id(id); + + } // endif: observation was valid + + // ... otherwise throw an exception + else { + std::string msg = + "Problem allocating GCTAOnOffObservation object."; + throw GException::invalid_instrument(G_READ, msg); + } + + // Append observation to container + append(*ptr); + + // Free observation (the append method made a deep copy) + delete ptr; + + } // endfor: looped over all observations + + // Return + return; +} + + +/***********************************************************************//** + * @brief Write observations into XML document + * + * @param[in] xml XML document. + * + * Write observations into the first observation library that is found in the + * XML document. In case that no observation library exists, one is added to + * the document. Please refer to the read() method for more information + * about the structure of the XML document. + ***************************************************************************/ +void GCTAOnOffObservations::write(GXml& xml) const +{ + // If there is no observation library then append one + if (xml.elements("observation_list") == 0) { + xml.append(GXmlElement("observation_list title=\"observation list\"")); + } + + // Get pointer on observation library + GXmlElement* lib = xml.element("observation_list", 0); + + // Write all observations into library + for (int i = 0; i < size(); ++i) { + + // Initialise pointer on observation + GXmlElement* obs = NULL; + + // Search corresponding observation + int n = xml.elements("observation"); + for (int k = 0; k < n; ++k) { + GXmlElement* element = xml.element("observation", k); + if (element->attribute("name") == m_obs[i]->name() && + element->attribute("id") == m_obs[i]->id() && + element->attribute("instrument") == m_obs[i]->instrument()) { + obs = element; + break; + } + } + + // If no observation with corresponding name, ID and instrument was + // found then append one now + if (obs == NULL) { + obs = lib->append("observation"); + obs->attribute("name", m_obs[i]->name()); + obs->attribute("id", m_obs[i]->id()); + obs->attribute("instrument", m_obs[i]->instrument()); + } + + // Write now observation + m_obs[i]->write(*obs); + + } // endfor: looped over all observaitons + + // Return + return; +} + + +/***********************************************************************//** + * @brief Load models from XML file + * + * @param[in] filename XML filename. + * + * Loads the models from a XML file. Please refer to the GModels::read() + * method for more information about the expected structure of the XML + * file. + ***************************************************************************/ +void GCTAOnOffObservations::models(const std::string& filename) +{ + // Load models + m_models.load(filename); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Print observation list information + * + * @param[in] chatter Chattiness (defaults to NORMAL). + * @return String containing observation list information + ***************************************************************************/ +std::string GCTAOnOffObservations::print(const GChatter& chatter) const +{ + // Initialise result string + std::string result; + + // Continue only if chatter is not silent + if (chatter != SILENT) { + + // Append header + result.append("=== GCTAOnOffObservations ==="); + + // Append information + result.append("\n"+gammalib::parformat("Number of observations")); + result.append(gammalib::str(size())); + + // NORMAL: Append observations + if (chatter >= NORMAL) { + for (int i = 0; i < size(); ++i) { + result.append("\n"); + result.append((*this)[i]->print()); + } + } // endif: chatter was normal + + // EXPLICIT: Append models + if (chatter >= EXPLICIT) { + result.append("\n"+m_models.print()); + } // endif: chatter was explicit + + } // endif: chatter was not silent + + // Return result + return result; +} + + +/*========================================================================== + = = + = Private methods = + = = + ==========================================================================*/ + +/***********************************************************************//** + * @brief Initialise class members + ***************************************************************************/ +void GCTAOnOffObservations::init_members(void) +{ + // Initialise members + m_obs.clear(); + m_models.clear(); + + // Return + return; +} + + +/***********************************************************************//** + * @brief Copy class members + * + * @param[in] obs Observation container. + * + * Copies all class members. Deep copies are created for the observations, + * which is what is needed to have a real physical copy of the members. + ***************************************************************************/ +void GCTAOnOffObservations::copy_members(const GCTAOnOffObservations& obs) +{ + // Copy attributes + m_models = obs.m_models; + + // Copy observations + m_obs.clear(); + for (int i = 0; i < obs.m_obs.size(); ++i) { + m_obs.push_back((obs.m_obs[i]->clone())); + } + + // Return + return; +} + + +/***********************************************************************//** + * @brief Delete class members + * + * Deallocates all observations. Since container classes that hold pointers + * need to handle the proper deallocation of memory, we loop here over all + * pointers and make sure that we deallocate all observations in the + * container. + ***************************************************************************/ +void GCTAOnOffObservations::free_members(void) +{ + // Free observations + for (int i = 0; i < m_obs.size(); ++i) { + delete m_obs[i]; + m_obs[i] = NULL; + } + + // Return + return; +} + + +/***********************************************************************//** + * @brief Return observation index by instrument and identifier + * + * @param[in] instrument Instrument. + * @param[in] id Observation identifier. + * @return Observation index (-1 of instrument and observation identifier + * has not been found) + * + * Returns observation index based on the specified @p instrument and + * observation identifier @p id. If no observation with the specified + * attributes has been found, the method returns -1. + ***************************************************************************/ +int GCTAOnOffObservations::get_index(const std::string& instrument, + const std::string& id) const +{ + // Initialise index + int index = -1; + + // Search observation with specified instrument and id + for (int i = 0; i < size(); ++i) { + if ((m_obs[i]->instrument() == instrument) && + (m_obs[i]->id() == id)) { + index = i; + break; + } + } + + // Return index + return index; +} -- libgit2 0.21.4