Commit 9dc564c7b58bf4dc2d17e442a72c520c7c879e55
1 parent
734b3482
Move set_obs_statistic() method from ctlikelihood to ctobservation and derive cs…
…cripts now from csobservation instead of cslikelihood since they do not use proper optimizers.
Showing
14 changed files
with
293 additions
and
275 deletions
NEWS
| ... | ... | @@ -164,7 +164,7 @@ None |
| 164 | 164 | |
| 165 | 165 | cslightcrv - Generate light curve |
| 166 | 166 | --------------------------------- |
| 167 | -The script now derived from the cslikelihood base class. A "statistic" | |
| 167 | +The script now derived from the csobservation base class. A "statistic" | |
| 168 | 168 | parameter was added to allow setting the optimizer statistic for the light |
| 169 | 169 | curve generation (#2237). |
| 170 | 170 | |
| ... | ... | @@ -220,14 +220,14 @@ The script was added to allow the generation of PHA, ARF and RMF datasets for cl |
| 220 | 220 | |
| 221 | 221 | csphasecrv - Generate phase curve |
| 222 | 222 | --------------------------------- |
| 223 | -The script now derived from the cslikelihood base class. A "statistic" | |
| 223 | +The script now derived from the csobservation base class. A "statistic" | |
| 224 | 224 | parameter was added to allow setting the optimizer statistic for the light |
| 225 | 225 | curve generation (#2237). |
| 226 | 226 | |
| 227 | 227 | |
| 228 | 228 | cspull - Generate pull distributions |
| 229 | 229 | ------------------------------------ |
| 230 | -The script now derived from the cslikelihood base class. A "statistic" | |
| 230 | +The script now derived from the csobservation base class. A "statistic" | |
| 231 | 231 | parameter was added to allow setting the optimizer statistic for the pull |
| 232 | 232 | distribution generation (#2237). |
| 233 | 233 | |
| ... | ... | @@ -247,14 +247,14 @@ None |
| 247 | 247 | |
| 248 | 248 | cssens - Computes sensitivity curve |
| 249 | 249 | ----------------------------------- |
| 250 | -The script now derived from the cslikelihood base class. A "statistic" | |
| 250 | +The script now derived from the csobservation base class. A "statistic" | |
| 251 | 251 | parameter was added to allow setting the optimizer statistic for the sensitivity |
| 252 | 252 | computation (#2237). |
| 253 | 253 | |
| 254 | 254 | |
| 255 | 255 | csspec - Generates Spectral Energy Distribution |
| 256 | 256 | ----------------------------------------------- |
| 257 | -The script now derived from the cslikelihood base class. A "statistic" | |
| 257 | +The script now derived from the csobservation base class. A "statistic" | |
| 258 | 258 | parameter was added to allow setting the optimizer statistic for the spectrum |
| 259 | 259 | generation (#2237). |
| 260 | 260 | |
| ... | ... | @@ -275,7 +275,7 @@ None |
| 275 | 275 | |
| 276 | 276 | cstsdist - Computes Test Statistic distribution |
| 277 | 277 | ----------------------------------------------- |
| 278 | -The script now derived from the cslikelihood base class. A "statistic" | |
| 278 | +The script now derived from the csobservation base class. A "statistic" | |
| 279 | 279 | parameter was added to allow setting the optimizer statistic for the Test |
| 280 | 280 | Statistic distribution generation (#2237). |
| 281 | 281 | ... | ... |
cscripts/cslightcrv.py
| ... | ... | @@ -28,7 +28,7 @@ from cscripts import ioutils |
| 28 | 28 | # ================ # |
| 29 | 29 | # cslightcrv class # |
| 30 | 30 | # ================ # |
| 31 | -class cslightcrv(ctools.cslikelihood): | |
| 31 | +class cslightcrv(ctools.csobservation): | |
| 32 | 32 | """ |
| 33 | 33 | Generates a lightcurve |
| 34 | 34 | |
| ... | ... | @@ -69,7 +69,7 @@ class cslightcrv(ctools.cslikelihood): |
| 69 | 69 | Constructor |
| 70 | 70 | """ |
| 71 | 71 | # Initialise application by calling the appropriate class constructor |
| 72 | - self._init_cslikelihood('cslightcrv', ctools.__version__, argv) | |
| 72 | + self._init_csobservation('cslightcrv', ctools.__version__, argv) | |
| 73 | 73 | |
| 74 | 74 | # Initialise some members |
| 75 | 75 | self._srcname = '' | ... | ... |
cscripts/csphasecrv.py
| ... | ... | @@ -28,7 +28,7 @@ from cscripts import ioutils |
| 28 | 28 | # ================ # |
| 29 | 29 | # csphasecrv class # |
| 30 | 30 | # ================ # |
| 31 | -class csphasecrv(ctools.cslikelihood): | |
| 31 | +class csphasecrv(ctools.csobservation): | |
| 32 | 32 | """ |
| 33 | 33 | Generates spectra in phase bins |
| 34 | 34 | |
| ... | ... | @@ -57,7 +57,7 @@ class csphasecrv(ctools.cslikelihood): |
| 57 | 57 | Constructor |
| 58 | 58 | """ |
| 59 | 59 | # Initialise application by calling the appropriate class constructor |
| 60 | - self._init_cslikelihood('csphasecrv', ctools.__version__, argv) | |
| 60 | + self._init_csobservation('csphasecrv', ctools.__version__, argv) | |
| 61 | 61 | |
| 62 | 62 | # Initialise some members. Phases are stored in a nested list |
| 63 | 63 | # [[ph1min,ph1max], [ph2min,ph2max],..] | ... | ... |
cscripts/cspull.py
| ... | ... | @@ -28,7 +28,7 @@ from cscripts import ioutils |
| 28 | 28 | # ============ # |
| 29 | 29 | # cspull class # |
| 30 | 30 | # ============ # |
| 31 | -class cspull(ctools.cslikelihood): | |
| 31 | +class cspull(ctools.csobservation): | |
| 32 | 32 | """ |
| 33 | 33 | Generates pull distributions for a model |
| 34 | 34 | """ |
| ... | ... | @@ -38,7 +38,7 @@ class cspull(ctools.cslikelihood): |
| 38 | 38 | Constructor |
| 39 | 39 | """ |
| 40 | 40 | # Initialise application by calling the appropriate class constructor |
| 41 | - self._init_cslikelihood('cspull', ctools.__version__, argv) | |
| 41 | + self._init_csobservation('cspull', ctools.__version__, argv) | |
| 42 | 42 | |
| 43 | 43 | # Initialise some members |
| 44 | 44 | self._edisp = False | ... | ... |
cscripts/cssens.py
| ... | ... | @@ -30,7 +30,7 @@ from cscripts import ioutils |
| 30 | 30 | # ============ # |
| 31 | 31 | # cssens class # |
| 32 | 32 | # ============ # |
| 33 | -class cssens(ctools.cslikelihood): | |
| 33 | +class cssens(ctools.csobservation): | |
| 34 | 34 | """ |
| 35 | 35 | Computes the CTA sensitivity |
| 36 | 36 | |
| ... | ... | @@ -48,7 +48,7 @@ class cssens(ctools.cslikelihood): |
| 48 | 48 | Constructor. |
| 49 | 49 | """ |
| 50 | 50 | # Initialise application by calling the appropriate class constructor |
| 51 | - self._init_cslikelihood('cssens', ctools.__version__, argv) | |
| 51 | + self._init_csobservation('cssens', ctools.__version__, argv) | |
| 52 | 52 | |
| 53 | 53 | # Initialise class members |
| 54 | 54 | self._ebounds = gammalib.GEbounds() | ... | ... |
cscripts/csspec.py
| ... | ... | @@ -26,7 +26,7 @@ import ctools |
| 26 | 26 | # ============ # |
| 27 | 27 | # csspec class # |
| 28 | 28 | # ============ # |
| 29 | -class csspec(ctools.cslikelihood): | |
| 29 | +class csspec(ctools.csobservation): | |
| 30 | 30 | """ |
| 31 | 31 | Generates a spectrum |
| 32 | 32 | |
| ... | ... | @@ -40,7 +40,7 @@ class csspec(ctools.cslikelihood): |
| 40 | 40 | Constructor |
| 41 | 41 | """ |
| 42 | 42 | # Initialise application by calling the appropriate class constructor |
| 43 | - self._init_cslikelihood('csspec', ctools.__version__, argv) | |
| 43 | + self._init_csobservation('csspec', ctools.__version__, argv) | |
| 44 | 44 | |
| 45 | 45 | # Initialise data members |
| 46 | 46 | self._ebounds = gammalib.GEbounds() | ... | ... |
cscripts/cstsdist.py
| ... | ... | @@ -29,7 +29,7 @@ from cscripts import ioutils |
| 29 | 29 | # ============== # |
| 30 | 30 | # cstsdist class # |
| 31 | 31 | # ============== # |
| 32 | -class cstsdist(ctools.cslikelihood): | |
| 32 | +class cstsdist(ctools.csobservation): | |
| 33 | 33 | """ |
| 34 | 34 | Generates Test Statistic distribution for a model |
| 35 | 35 | """ |
| ... | ... | @@ -40,7 +40,7 @@ class cstsdist(ctools.cslikelihood): |
| 40 | 40 | Constructor. |
| 41 | 41 | """ |
| 42 | 42 | # Initialise application by calling the appropriate class constructor |
| 43 | - self._init_cslikelihood('cstsdist', ctools.__version__, argv) | |
| 43 | + self._init_csobservation('cstsdist', ctools.__version__, argv) | |
| 44 | 44 | |
| 45 | 45 | # Initialise some members |
| 46 | 46 | self._srcname = '' | ... | ... |
pyext/ctlikelihood.i
| ... | ... | @@ -91,12 +91,10 @@ public: |
| 91 | 91 | const GOptimizer* opt(void) const; |
| 92 | 92 | |
| 93 | 93 | // Make methods protected in Python by prepending an underscore |
| 94 | - %rename(_evaluate) evaluate; | |
| 95 | - %rename(_set_obs_statistic) set_obs_statistic; | |
| 94 | + %rename(_evaluate) evaluate; | |
| 96 | 95 | |
| 97 | 96 | // Protected methods |
| 98 | 97 | double evaluate(GModelPar& par, const double& value); |
| 99 | - void set_obs_statistic(const std::string& statistic); | |
| 100 | 98 | }; |
| 101 | 99 | |
| 102 | 100 | ... | ... |
pyext/ctobservation.i
| ... | ... | @@ -90,15 +90,22 @@ public: |
| 90 | 90 | // Public methods |
| 91 | 91 | void obs(const GObservations& obs); |
| 92 | 92 | const GObservations& obs(void) const; |
| 93 | - void write_ogip_keywords(GFitsHDU* hdu) const; | |
| 94 | 93 | |
| 95 | 94 | // Make methods private in Python by prepending an underscore |
| 96 | 95 | %rename(_first_unbinned_observation) first_unbinned_observation; |
| 97 | 96 | %rename(_next_unbinned_observation) next_unbinned_observation; |
| 97 | + %rename(_write_ogip_keywords) write_ogip_keywords; | |
| 98 | + %rename(_set_obs_statistic) set_obs_statistic; | |
| 99 | + %rename(_save_events_fits) save_events_fits; | |
| 100 | + %rename(_save_events_xml) save_events_xml; | |
| 98 | 101 | |
| 99 | 102 | // Protected methods |
| 100 | 103 | GCTAObservation* first_unbinned_observation(void); |
| 101 | 104 | GCTAObservation* next_unbinned_observation(void); |
| 105 | + void write_ogip_keywords(GFitsHDU* hdu) const; | |
| 106 | + void set_obs_statistic(const std::string& statistic); | |
| 107 | + void save_events_fits(void); | |
| 108 | + void save_events_xml(void); | |
| 102 | 109 | }; |
| 103 | 110 | |
| 104 | 111 | ... | ... |
src/support/ctlikelihood.cpp
| ... | ... | @@ -196,49 +196,11 @@ ctlikelihood& ctlikelihood::operator=(const ctlikelihood& app) |
| 196 | 196 | |
| 197 | 197 | /*========================================================================== |
| 198 | 198 | = = |
| 199 | - = Private methods = | |
| 199 | + = Protected methods exposed in Python = | |
| 200 | 200 | = = |
| 201 | 201 | ==========================================================================*/ |
| 202 | 202 | |
| 203 | 203 | /***********************************************************************//** |
| 204 | - * @brief Initialise class members | |
| 205 | - ***************************************************************************/ | |
| 206 | -void ctlikelihood::init_members(void) | |
| 207 | -{ | |
| 208 | - // Initialise members | |
| 209 | - m_opt.clear(); | |
| 210 | - | |
| 211 | - // Return | |
| 212 | - return; | |
| 213 | -} | |
| 214 | - | |
| 215 | - | |
| 216 | -/***********************************************************************//** | |
| 217 | - * @brief Copy class members | |
| 218 | - * | |
| 219 | - * @param[in] app Likelihood tool. | |
| 220 | - ***************************************************************************/ | |
| 221 | -void ctlikelihood::copy_members(const ctlikelihood& app) | |
| 222 | -{ | |
| 223 | - // Copy members | |
| 224 | - m_opt = app.m_opt; | |
| 225 | - | |
| 226 | - // Return | |
| 227 | - return; | |
| 228 | -} | |
| 229 | - | |
| 230 | - | |
| 231 | -/***********************************************************************//** | |
| 232 | - * @brief Delete class members | |
| 233 | - ***************************************************************************/ | |
| 234 | -void ctlikelihood::free_members(void) | |
| 235 | -{ | |
| 236 | - // Return | |
| 237 | - return; | |
| 238 | -} | |
| 239 | - | |
| 240 | - | |
| 241 | -/***********************************************************************//** | |
| 242 | 204 | * @brief Evaluates the log-likelihood function |
| 243 | 205 | * |
| 244 | 206 | * @param[in] par Model parameter |
| ... | ... | @@ -254,8 +216,8 @@ double ctlikelihood::evaluate(GModelPar& par, const double& value) |
| 254 | 216 | |
| 255 | 217 | // Throw an exception if the parameter is below the minimum boundary |
| 256 | 218 | if (par.has_min() && value < par.factor_min()) { |
| 257 | - std::string msg = "Value "+gammalib::str(value)+" of parameter \""+ | |
| 258 | - par.name()+"\" is below its minimum boundary "+ | |
| 219 | + std::string msg = "Parameter \""+par.name()+"\" value"+ | |
| 220 | + gammalib::str(value)+" is below its minimum boundary "+ | |
| 259 | 221 | gammalib::str(par.factor_min())+". To omit this " |
| 260 | 222 | "error please lower the minimum parameter boundary."; |
| 261 | 223 | throw GException::invalid_value(G_EVALUATE, msg); |
| ... | ... | @@ -263,8 +225,8 @@ double ctlikelihood::evaluate(GModelPar& par, const double& value) |
| 263 | 225 | |
| 264 | 226 | // Throw an exception if the parameter is above the maximum boundary |
| 265 | 227 | if (par.has_max() && value > par.factor_max()) { |
| 266 | - std::string msg = "Value "+gammalib::str(value)+" of parameter \""+ | |
| 267 | - par.name()+"\" is above its maximum boundary "+ | |
| 228 | + std::string msg = "Parameter \""+par.name()+"\" value"+ | |
| 229 | + gammalib::str(value)+" is above its maximum boundary "+ | |
| 268 | 230 | gammalib::str(par.factor_max())+". To omit this " |
| 269 | 231 | "error please raise the maximum parameter boundary."; |
| 270 | 232 | throw GException::invalid_value(G_EVALUATE, msg); |
| ... | ... | @@ -290,57 +252,45 @@ double ctlikelihood::evaluate(GModelPar& par, const double& value) |
| 290 | 252 | } |
| 291 | 253 | |
| 292 | 254 | |
| 255 | +/*========================================================================== | |
| 256 | + = = | |
| 257 | + = Protected methods = | |
| 258 | + = = | |
| 259 | + ==========================================================================*/ | |
| 260 | + | |
| 293 | 261 | /***********************************************************************//** |
| 294 | - * @brief Set fit statistic for CTA observations | |
| 295 | - * | |
| 296 | - * @param[in] statistic Requested fit statistic. | |
| 297 | - * | |
| 298 | - * Sets the fit statistic for all CTA observations. The method handles | |
| 299 | - * regular CTA observations as well as On/Off observations. | |
| 300 | - * | |
| 301 | - * For regular CTA observations of type GCTAObservation the fit statistic is | |
| 302 | - * only set for binned or stacked observations. Possible values are | |
| 303 | - * @c POISSON, @c GAUSSIAN or @c CHI2 (case insensitive). | |
| 304 | - * | |
| 305 | - * For On/Off CTA observations of type GCTAOnOffObservation the fit statistic | |
| 306 | - * is always set. Possible values are @c POISSON, @c CSTAT or @c WSTAT (case | |
| 307 | - * insensitive). | |
| 262 | + * @brief Initialise class members | |
| 263 | + ***************************************************************************/ | |
| 264 | +void ctlikelihood::init_members(void) | |
| 265 | +{ | |
| 266 | + // Initialise members | |
| 267 | + m_opt.clear(); | |
| 268 | + | |
| 269 | + // Return | |
| 270 | + return; | |
| 271 | +} | |
| 272 | + | |
| 273 | + | |
| 274 | +/***********************************************************************//** | |
| 275 | + * @brief Copy class members | |
| 308 | 276 | * |
| 309 | - * If @p statistic is any other string the fit statistic of the observations | |
| 310 | - * will not be modified. | |
| 277 | + * @param[in] app Likelihood tool. | |
| 311 | 278 | ***************************************************************************/ |
| 312 | -void ctlikelihood::set_obs_statistic(const std::string& statistic) | |
| 279 | +void ctlikelihood::copy_members(const ctlikelihood& app) | |
| 313 | 280 | { |
| 314 | - // Convert statistic to upper case | |
| 315 | - std::string ustatistic = gammalib::toupper(statistic); | |
| 316 | - | |
| 317 | - // Flag which observation types are handled | |
| 318 | - bool handle_cta = ((ustatistic == "POISSON") || | |
| 319 | - (ustatistic == "GAUSSIAN") || | |
| 320 | - (ustatistic == "CHI2")); | |
| 321 | - bool handle_onoff = ((ustatistic == "POISSON") || | |
| 322 | - (ustatistic == "CSTAT") || | |
| 323 | - (ustatistic == "WSTAT")); | |
| 324 | - | |
| 325 | - // Loop over all observation in container | |
| 326 | - for (int i = 0; i < m_obs.size(); ++i) { | |
| 327 | - | |
| 328 | - // For regular CTA observations only set statistic for binned or | |
| 329 | - // stacked observations | |
| 330 | - if ((m_obs[i]->classname() == "GCTAObservation") && handle_cta) { | |
| 331 | - if (m_obs[i]->events()->classname() == "GCTAEventCube") { | |
| 332 | - m_obs[i]->statistic(statistic); | |
| 333 | - } | |
| 334 | - } | |
| 335 | - | |
| 336 | - // For On/Off CTA observations always set statistic | |
| 337 | - else if ((m_obs[i]->classname() == "GCTAOnOffObservation") && | |
| 338 | - handle_onoff) { | |
| 339 | - m_obs[i]->statistic(statistic); | |
| 340 | - } | |
| 341 | - | |
| 342 | - } // endfor: looped over all observations | |
| 281 | + // Copy members | |
| 282 | + m_opt = app.m_opt; | |
| 283 | + | |
| 284 | + // Return | |
| 285 | + return; | |
| 286 | +} | |
| 343 | 287 | |
| 288 | + | |
| 289 | +/***********************************************************************//** | |
| 290 | + * @brief Delete class members | |
| 291 | + ***************************************************************************/ | |
| 292 | +void ctlikelihood::free_members(void) | |
| 293 | +{ | |
| 344 | 294 | // Return |
| 345 | 295 | return; |
| 346 | 296 | } | ... | ... |
src/support/ctlikelihood.hpp
src/support/ctobservation.cpp
| ... | ... | @@ -192,119 +192,14 @@ ctobservation& ctobservation::operator=(const ctobservation& app) |
| 192 | 192 | = = |
| 193 | 193 | ==========================================================================*/ |
| 194 | 194 | |
| 195 | -/***********************************************************************//** | |
| 196 | - * @brief Write OGIP keywords in FITS HDU | |
| 197 | - * | |
| 198 | - * @param[in,out] hdu Pointer to FITS HDU. | |
| 199 | - * | |
| 200 | - * Writes OGIP keywords in FITS HDU. | |
| 201 | - ***************************************************************************/ | |
| 202 | -void ctobservation::write_ogip_keywords(GFitsHDU* hdu) const | |
| 203 | -{ | |
| 204 | - // Continue only if pointer is valid | |
| 205 | - if (hdu != NULL) { | |
| 206 | - | |
| 207 | - // Set creator | |
| 208 | - std::string creator = this->name() + " v" + this->version(); | |
| 209 | - | |
| 210 | - // Compute times | |
| 211 | - std::string utc_obs = m_ogip_tstart.utc(); | |
| 212 | - std::string utc_end = m_ogip_tstop.utc(); | |
| 213 | - std::string date_obs = utc_obs.substr(0, 10); | |
| 214 | - std::string time_obs = utc_obs.substr(11, 8); | |
| 215 | - std::string date_end = utc_end.substr(0, 10); | |
| 216 | - std::string time_end = utc_end.substr(11, 8); | |
| 217 | - | |
| 218 | - // Compute deadtime correction | |
| 219 | - double deadc = (m_ogip_ontime > 0.0) ? m_ogip_livetime / m_ogip_ontime : 1.0; | |
| 220 | - | |
| 221 | - // Set observation information | |
| 222 | - hdu->card("CREATOR", creator, "Program which created the file"); | |
| 223 | - hdu->card("TELESCOP", m_ogip_telescope, "Telescope"); | |
| 224 | - | |
| 225 | - // Set observation time information | |
| 226 | - hdu->card("DATE_OBS", date_obs, "Observation start date"); | |
| 227 | - hdu->card("TIME_OBS", time_obs, "Observation start time"); | |
| 228 | - hdu->card("DATE_END", date_end, "Observation end date"); | |
| 229 | - hdu->card("TIME_END", time_end, "Observation end time"); | |
| 230 | - hdu->card("TELAPSE", m_ogip_telapse, "[s] Elapsed time"); | |
| 231 | - hdu->card("ONTIME", m_ogip_ontime, "[s] Total good time including deadtime"); | |
| 232 | - hdu->card("LIVETIME", m_ogip_livetime, "[s] Total livetime"); | |
| 233 | - hdu->card("EXPOSURE", m_ogip_exposure, "[s] Exposure time"); | |
| 234 | - hdu->card("DEADC", deadc, "Deadtime correction factor"); | |
| 235 | - | |
| 236 | - } // endif: pointer was valid | |
| 237 | - | |
| 238 | - // Return | |
| 239 | - return; | |
| 240 | -} | |
| 241 | - | |
| 242 | 195 | |
| 243 | 196 | /*========================================================================== |
| 244 | 197 | = = |
| 245 | - = Private methods = | |
| 198 | + = Protected methods exposed in Python = | |
| 246 | 199 | = = |
| 247 | 200 | ==========================================================================*/ |
| 248 | 201 | |
| 249 | 202 | /***********************************************************************//** |
| 250 | - * @brief Initialise class members | |
| 251 | - ***************************************************************************/ | |
| 252 | -void ctobservation::init_members(void) | |
| 253 | -{ | |
| 254 | - // Initialise protected members | |
| 255 | - m_obs.clear(); | |
| 256 | - | |
| 257 | - // Initialise private members | |
| 258 | - m_ogip_telescope.clear(); | |
| 259 | - m_ogip_tstart.clear(); | |
| 260 | - m_ogip_tstop.clear(); | |
| 261 | - m_ogip_telapse = 0.0; | |
| 262 | - m_ogip_exposure = 0.0; | |
| 263 | - m_ogip_ontime = 0.0; | |
| 264 | - m_ogip_livetime = 0.0; | |
| 265 | - m_index_unbinned = 0; | |
| 266 | - | |
| 267 | - // Return | |
| 268 | - return; | |
| 269 | -} | |
| 270 | - | |
| 271 | - | |
| 272 | -/***********************************************************************//** | |
| 273 | - * @brief Copy class members | |
| 274 | - * | |
| 275 | - * @param[in] app Observation tool. | |
| 276 | - ***************************************************************************/ | |
| 277 | -void ctobservation::copy_members(const ctobservation& app) | |
| 278 | -{ | |
| 279 | - // Copy protected members | |
| 280 | - m_obs = app.m_obs; | |
| 281 | - | |
| 282 | - // Copy private members | |
| 283 | - m_ogip_telescope = app.m_ogip_telescope; | |
| 284 | - m_ogip_tstart = app.m_ogip_tstart; | |
| 285 | - m_ogip_tstop = app.m_ogip_tstop; | |
| 286 | - m_ogip_telapse = app.m_ogip_telapse; | |
| 287 | - m_ogip_exposure = app.m_ogip_exposure; | |
| 288 | - m_ogip_ontime = app.m_ogip_ontime; | |
| 289 | - m_ogip_livetime = app.m_ogip_livetime; | |
| 290 | - m_index_unbinned = app.m_index_unbinned; | |
| 291 | - | |
| 292 | - // Return | |
| 293 | - return; | |
| 294 | -} | |
| 295 | - | |
| 296 | - | |
| 297 | -/***********************************************************************//** | |
| 298 | - * @brief Delete class members | |
| 299 | - ***************************************************************************/ | |
| 300 | -void ctobservation::free_members(void) | |
| 301 | -{ | |
| 302 | - // Return | |
| 303 | - return; | |
| 304 | -} | |
| 305 | - | |
| 306 | - | |
| 307 | -/***********************************************************************//** | |
| 308 | 203 | * @brief Return first unbinned CTA observation (const version) |
| 309 | 204 | * |
| 310 | 205 | * @return Const pointer to first unbinned CTA observation |
| ... | ... | @@ -424,6 +319,110 @@ const GCTAObservation* ctobservation::next_unbinned_observation(void) const |
| 424 | 319 | |
| 425 | 320 | |
| 426 | 321 | /***********************************************************************//** |
| 322 | + * @brief Write OGIP keywords in FITS HDU | |
| 323 | + * | |
| 324 | + * @param[in,out] hdu Pointer to FITS HDU. | |
| 325 | + * | |
| 326 | + * Writes OGIP keywords in FITS HDU. | |
| 327 | + ***************************************************************************/ | |
| 328 | +void ctobservation::write_ogip_keywords(GFitsHDU* hdu) const | |
| 329 | +{ | |
| 330 | + // Continue only if pointer is valid | |
| 331 | + if (hdu != NULL) { | |
| 332 | + | |
| 333 | + // Set creator | |
| 334 | + std::string creator = this->name() + " v" + this->version(); | |
| 335 | + | |
| 336 | + // Compute times | |
| 337 | + std::string utc_obs = m_ogip_tstart.utc(); | |
| 338 | + std::string utc_end = m_ogip_tstop.utc(); | |
| 339 | + std::string date_obs = utc_obs.substr(0, 10); | |
| 340 | + std::string time_obs = utc_obs.substr(11, 8); | |
| 341 | + std::string date_end = utc_end.substr(0, 10); | |
| 342 | + std::string time_end = utc_end.substr(11, 8); | |
| 343 | + | |
| 344 | + // Compute deadtime correction | |
| 345 | + double deadc = (m_ogip_ontime > 0.0) ? m_ogip_livetime / m_ogip_ontime : 1.0; | |
| 346 | + | |
| 347 | + // Set observation information | |
| 348 | + hdu->card("CREATOR", creator, "Program which created the file"); | |
| 349 | + hdu->card("TELESCOP", m_ogip_telescope, "Telescope"); | |
| 350 | + | |
| 351 | + // Set observation time information | |
| 352 | + hdu->card("DATE_OBS", date_obs, "Observation start date"); | |
| 353 | + hdu->card("TIME_OBS", time_obs, "Observation start time"); | |
| 354 | + hdu->card("DATE_END", date_end, "Observation end date"); | |
| 355 | + hdu->card("TIME_END", time_end, "Observation end time"); | |
| 356 | + hdu->card("TELAPSE", m_ogip_telapse, "[s] Elapsed time"); | |
| 357 | + hdu->card("ONTIME", m_ogip_ontime, "[s] Total good time including deadtime"); | |
| 358 | + hdu->card("LIVETIME", m_ogip_livetime, "[s] Total livetime"); | |
| 359 | + hdu->card("EXPOSURE", m_ogip_exposure, "[s] Exposure time"); | |
| 360 | + hdu->card("DEADC", deadc, "Deadtime correction factor"); | |
| 361 | + | |
| 362 | + } // endif: pointer was valid | |
| 363 | + | |
| 364 | + // Return | |
| 365 | + return; | |
| 366 | +} | |
| 367 | + | |
| 368 | + | |
| 369 | +/***********************************************************************//** | |
| 370 | + * @brief Set fit statistic for CTA observations | |
| 371 | + * | |
| 372 | + * @param[in] statistic Requested fit statistic. | |
| 373 | + * | |
| 374 | + * Sets the fit statistic for all CTA observations. The method handles | |
| 375 | + * regular CTA observations as well as On/Off observations. | |
| 376 | + * | |
| 377 | + * For regular CTA observations of type GCTAObservation the fit statistic is | |
| 378 | + * only set for binned or stacked observations. Possible values are | |
| 379 | + * @c POISSON, @c GAUSSIAN or @c CHI2 (case insensitive). | |
| 380 | + * | |
| 381 | + * For On/Off CTA observations of type GCTAOnOffObservation the fit statistic | |
| 382 | + * is always set. Possible values are @c POISSON, @c CSTAT or @c WSTAT (case | |
| 383 | + * insensitive). | |
| 384 | + * | |
| 385 | + * If @p statistic is any other string the fit statistic of the observations | |
| 386 | + * will not be modified. | |
| 387 | + ***************************************************************************/ | |
| 388 | +void ctobservation::set_obs_statistic(const std::string& statistic) | |
| 389 | +{ | |
| 390 | + // Convert statistic to upper case | |
| 391 | + std::string ustatistic = gammalib::toupper(statistic); | |
| 392 | + | |
| 393 | + // Flag which observation types are handled | |
| 394 | + bool handle_cta = ((ustatistic == "POISSON") || | |
| 395 | + (ustatistic == "GAUSSIAN") || | |
| 396 | + (ustatistic == "CHI2")); | |
| 397 | + bool handle_onoff = ((ustatistic == "POISSON") || | |
| 398 | + (ustatistic == "CSTAT") || | |
| 399 | + (ustatistic == "WSTAT")); | |
| 400 | + | |
| 401 | + // Loop over all observation in container | |
| 402 | + for (int i = 0; i < m_obs.size(); ++i) { | |
| 403 | + | |
| 404 | + // For regular CTA observations only set statistic for binned or | |
| 405 | + // stacked observations | |
| 406 | + if ((m_obs[i]->classname() == "GCTAObservation") && handle_cta) { | |
| 407 | + if (m_obs[i]->events()->classname() == "GCTAEventCube") { | |
| 408 | + m_obs[i]->statistic(statistic); | |
| 409 | + } | |
| 410 | + } | |
| 411 | + | |
| 412 | + // For On/Off CTA observations always set statistic | |
| 413 | + else if ((m_obs[i]->classname() == "GCTAOnOffObservation") && | |
| 414 | + handle_onoff) { | |
| 415 | + m_obs[i]->statistic(statistic); | |
| 416 | + } | |
| 417 | + | |
| 418 | + } // endfor: looped over all observations | |
| 419 | + | |
| 420 | + // Return | |
| 421 | + return; | |
| 422 | +} | |
| 423 | + | |
| 424 | + | |
| 425 | +/***********************************************************************//** | |
| 427 | 426 | * @brief Save event list in FITS format. |
| 428 | 427 | * |
| 429 | 428 | * Save the event list as a FITS file. The file name of the FITS file is |
| ... | ... | @@ -541,3 +540,67 @@ void ctobservation::save_events_xml(void) |
| 541 | 540 | // Return |
| 542 | 541 | return; |
| 543 | 542 | } |
| 543 | + | |
| 544 | + | |
| 545 | +/*========================================================================== | |
| 546 | + = = | |
| 547 | + = Protected methods = | |
| 548 | + = = | |
| 549 | + ==========================================================================*/ | |
| 550 | + | |
| 551 | +/***********************************************************************//** | |
| 552 | + * @brief Initialise class members | |
| 553 | + ***************************************************************************/ | |
| 554 | +void ctobservation::init_members(void) | |
| 555 | +{ | |
| 556 | + // Initialise protected members | |
| 557 | + m_obs.clear(); | |
| 558 | + | |
| 559 | + // Initialise private members | |
| 560 | + m_ogip_telescope.clear(); | |
| 561 | + m_ogip_tstart.clear(); | |
| 562 | + m_ogip_tstop.clear(); | |
| 563 | + m_ogip_telapse = 0.0; | |
| 564 | + m_ogip_exposure = 0.0; | |
| 565 | + m_ogip_ontime = 0.0; | |
| 566 | + m_ogip_livetime = 0.0; | |
| 567 | + m_index_unbinned = 0; | |
| 568 | + | |
| 569 | + // Return | |
| 570 | + return; | |
| 571 | +} | |
| 572 | + | |
| 573 | + | |
| 574 | +/***********************************************************************//** | |
| 575 | + * @brief Copy class members | |
| 576 | + * | |
| 577 | + * @param[in] app Observation tool. | |
| 578 | + ***************************************************************************/ | |
| 579 | +void ctobservation::copy_members(const ctobservation& app) | |
| 580 | +{ | |
| 581 | + // Copy protected members | |
| 582 | + m_obs = app.m_obs; | |
| 583 | + | |
| 584 | + // Copy private members | |
| 585 | + m_ogip_telescope = app.m_ogip_telescope; | |
| 586 | + m_ogip_tstart = app.m_ogip_tstart; | |
| 587 | + m_ogip_tstop = app.m_ogip_tstop; | |
| 588 | + m_ogip_telapse = app.m_ogip_telapse; | |
| 589 | + m_ogip_exposure = app.m_ogip_exposure; | |
| 590 | + m_ogip_ontime = app.m_ogip_ontime; | |
| 591 | + m_ogip_livetime = app.m_ogip_livetime; | |
| 592 | + m_index_unbinned = app.m_index_unbinned; | |
| 593 | + | |
| 594 | + // Return | |
| 595 | + return; | |
| 596 | +} | |
| 597 | + | |
| 598 | + | |
| 599 | +/***********************************************************************//** | |
| 600 | + * @brief Delete class members | |
| 601 | + ***************************************************************************/ | |
| 602 | +void ctobservation::free_members(void) | |
| 603 | +{ | |
| 604 | + // Return | |
| 605 | + return; | |
| 606 | +} | ... | ... |
src/support/ctobservation.hpp
| ... | ... | @@ -64,24 +64,25 @@ public: |
| 64 | 64 | // Methods |
| 65 | 65 | void obs(const GObservations& obs); |
| 66 | 66 | const GObservations& obs(void) const; |
| 67 | - void write_ogip_keywords(GFitsHDU* hdu) const; | |
| 68 | 67 | |
| 69 | 68 | #ifndef SWIG |
| 70 | 69 | protected: |
| 71 | 70 | #endif |
| 72 | - // Iterator methods and members | |
| 71 | + // Protected methods | |
| 73 | 72 | GCTAObservation* first_unbinned_observation(void); |
| 74 | 73 | GCTAObservation* next_unbinned_observation(void); |
| 75 | 74 | const GCTAObservation* first_unbinned_observation(void) const; |
| 76 | 75 | const GCTAObservation* next_unbinned_observation(void) const; |
| 76 | + void write_ogip_keywords(GFitsHDU* hdu) const; | |
| 77 | + void set_obs_statistic(const std::string& statistic); | |
| 78 | + void save_events_fits(void); | |
| 79 | + void save_events_xml(void); | |
| 77 | 80 | |
| 78 | 81 | protected: |
| 79 | 82 | // Protected methods |
| 80 | 83 | void init_members(void); |
| 81 | 84 | void copy_members(const ctobservation& app); |
| 82 | 85 | void free_members(void); |
| 83 | - void save_events_fits(void); | |
| 84 | - void save_events_xml(void); | |
| 85 | 86 | |
| 86 | 87 | // Protected members |
| 87 | 88 | GObservations m_obs; //!< Observation container | ... | ... |
test/test_cscript.py
| ... | ... | @@ -166,32 +166,6 @@ class ctobservation_test(ctools.csobservation): |
| 166 | 166 | def check_next_unbinned_observation(self): |
| 167 | 167 | return self._next_unbinned_observation() |
| 168 | 168 | |
| 169 | - | |
| 170 | -# ======================= # | |
| 171 | -# ctlikelihood_test class # | |
| 172 | -# ======================= # | |
| 173 | -class ctlikelihood_test(ctools.cslikelihood): | |
| 174 | - """ | |
| 175 | - Test class for the ctlikelihood base class | |
| 176 | - """ | |
| 177 | - # Constructor | |
| 178 | - def __init__(self, *argv): | |
| 179 | - """ | |
| 180 | - Constructor | |
| 181 | - """ | |
| 182 | - # Initialise application by calling the appropriate class constructor | |
| 183 | - self._init_cslikelihood('ctlikelihood_test', ctools.__version__, argv) | |
| 184 | - | |
| 185 | - # Return | |
| 186 | - return | |
| 187 | - | |
| 188 | - # Check evaluate() method | |
| 189 | - def check_evaluate(self, value): | |
| 190 | - par = gammalib.GModelPar() | |
| 191 | - par.value(2.0) | |
| 192 | - par.range(1.0, 3.0) | |
| 193 | - return self._evaluate(par, value) | |
| 194 | - | |
| 195 | 169 | # Check set_obs_statistic() method for unbinned observation |
| 196 | 170 | def check_set_obs_statistic_unbinned(self, statistic): |
| 197 | 171 | obs = gammalib.GObservations() |
| ... | ... | @@ -224,6 +198,32 @@ class ctlikelihood_test(ctools.cslikelihood): |
| 224 | 198 | return gammalib.toupper(self.obs()[0].statistic()) |
| 225 | 199 | |
| 226 | 200 | |
| 201 | +# ======================= # | |
| 202 | +# ctlikelihood_test class # | |
| 203 | +# ======================= # | |
| 204 | +class ctlikelihood_test(ctools.cslikelihood): | |
| 205 | + """ | |
| 206 | + Test class for the ctlikelihood base class | |
| 207 | + """ | |
| 208 | + # Constructor | |
| 209 | + def __init__(self, *argv): | |
| 210 | + """ | |
| 211 | + Constructor | |
| 212 | + """ | |
| 213 | + # Initialise application by calling the appropriate class constructor | |
| 214 | + self._init_cslikelihood('ctlikelihood_test', ctools.__version__, argv) | |
| 215 | + | |
| 216 | + # Return | |
| 217 | + return | |
| 218 | + | |
| 219 | + # Check evaluate() method | |
| 220 | + def check_evaluate(self, value): | |
| 221 | + par = gammalib.GModelPar() | |
| 222 | + par.value(2.0) | |
| 223 | + par.range(1.0, 3.0) | |
| 224 | + return self._evaluate(par, value) | |
| 225 | + | |
| 226 | + | |
| 227 | 227 | # ============================ # |
| 228 | 228 | # Test class for base classes # |
| 229 | 229 | # ============================ # |
| ... | ... | @@ -616,33 +616,6 @@ class Test(test): |
| 616 | 616 | self.test_assert(empty.check_next_unbinned_observation() is None, |
| 617 | 617 | 'Check next_unbinned_observation() method') |
| 618 | 618 | |
| 619 | - # Return | |
| 620 | - return | |
| 621 | - | |
| 622 | - # Test ctlikelihood base class | |
| 623 | - def _test_ctlikelihood(self): | |
| 624 | - """ | |
| 625 | - Test ctlikelihood base class | |
| 626 | - """ | |
| 627 | - # Allocate empty | |
| 628 | - empty = ctlikelihood_test() | |
| 629 | - | |
| 630 | - # Test evaluate() method for minimum value violation | |
| 631 | - self.test_try('Check evaluate() for minimum value violation') | |
| 632 | - try: | |
| 633 | - empty.check_evaluate(0.0) | |
| 634 | - self.test_try_failure('Exception not thrown') | |
| 635 | - except ValueError: | |
| 636 | - self.test_try_success() | |
| 637 | - | |
| 638 | - # Test evaluate() method for maximum value violation | |
| 639 | - self.test_try('Check evaluate() for maximum value violation') | |
| 640 | - try: | |
| 641 | - empty.check_evaluate(4.0) | |
| 642 | - self.test_try_failure('Exception not thrown') | |
| 643 | - except ValueError: | |
| 644 | - self.test_try_success() | |
| 645 | - | |
| 646 | 619 | # Test set_obs_statistic() method for unbinned observation |
| 647 | 620 | self.test_value(empty.check_set_obs_statistic_unbinned('CSTAT'), |
| 648 | 621 | 'CSTAT', |
| ... | ... | @@ -687,3 +660,30 @@ class Test(test): |
| 687 | 660 | |
| 688 | 661 | # Return |
| 689 | 662 | return |
| 663 | + | |
| 664 | + # Test ctlikelihood base class | |
| 665 | + def _test_ctlikelihood(self): | |
| 666 | + """ | |
| 667 | + Test ctlikelihood base class | |
| 668 | + """ | |
| 669 | + # Allocate empty | |
| 670 | + empty = ctlikelihood_test() | |
| 671 | + | |
| 672 | + # Test evaluate() method for minimum value violation | |
| 673 | + self.test_try('Check evaluate() for minimum value violation') | |
| 674 | + try: | |
| 675 | + empty.check_evaluate(0.0) | |
| 676 | + self.test_try_failure('Exception not thrown') | |
| 677 | + except ValueError: | |
| 678 | + self.test_try_success() | |
| 679 | + | |
| 680 | + # Test evaluate() method for maximum value violation | |
| 681 | + self.test_try('Check evaluate() for maximum value violation') | |
| 682 | + try: | |
| 683 | + empty.check_evaluate(4.0) | |
| 684 | + self.test_try_failure('Exception not thrown') | |
| 685 | + except ValueError: | |
| 686 | + self.test_try_success() | |
| 687 | + | |
| 688 | + # Return | |
| 689 | + return | ... | ... |