Gonzalo Rodriguez Fernandez / gammalib

Commit a16cd0db5fca9db0f9d00b95d9a92d529fc90ea2

Authored by Jürgen Knödlseder
1 parent 9fb9976d

Remove event fetching and disposal from GCTAObservation as this is now handled a… 

…t GCTAEventList level.
Inline Side-by-side
Showing 5 changed files with 140 additions and 439 deletions
ChangeLog
  View file @a16cd0d
... ... @@ -3,10 +3,11 @@
3 3 * Version 1.1.0 released
4 4 =========================
5 5  
6   - Implement event fetching and disposal in GCTAEventList (#1648)
  6 + Remove event fetching from GCTAObservation class (#1648)
  7 + Implement event fetching and disposal in GCTAEventList class (#1648)
7 8 Add GFilename::fullname() method
8 9 Add GCTAEventList::remove() method
9   - Handle also additional columns in GCTAEventList::append method
  10 + Handle also additional columns in GCTAEventList::append() method
10 11 Correct GFitsTableCol::insert() method for empty columns
11 12 Correct GCTAEdisp2D::compute_ebounds_xxx() methods (#1662)
12 13 Use GFilename in GFits to store filename (#1649)
... ...
... ... @@ -13,9 +13,12 @@ made since the last release of GammaLib.
13 13  
14 14 2. Public interface modifications
15 15 ---------------------------------
16   -The following methods return now GFilename instead of std::string:
17   -- GFits::filename()
  16 +The following methods have been removed from the public interface:
  17 +- GCTAObservation::dispose_events()
18 18  
  19 +The return value has been changed for the following methods:
  20 +- GFits::filename() returns now GFilename instead of std::string
  21 +- GCTAObservation::eventtype() returns now std::string instead of const std::string&
19 22  
20 23 3. Configuration
21 24 -----------------
... ... @@ -112,7 +115,8 @@ None
112 115 -----------------
113 116 Event fetching and disposal has been implemented in the GCTAEventList class by
114 117 additing the fetch() and dispose() methods. Now, events are only loaded upon
115   -request and can be dispatched after user to save memory (#1648).
  118 +request and can be dispatched after user to save memory. The event fetching and
  119 +disposal logic has been removed from GCTAObservation (#1648).
116 120  
117 121 Add GCTAEventList::remove() method to remove events from an event list.
118 122  
... ...
inst/cta/include/GCTAObservation.hpp
  View file @a16cd0d
... ... @@ -87,10 +87,6 @@ public:
87 87 virtual void write(GXmlElement& xml) const;
88 88 virtual std::string print(const GChatter& chatter = NORMAL) const;
89 89  
90   - // Overwrite virtual base class methods
91   - virtual const GEvents* events(void) const;
92   - virtual void events(const GEvents& events);
93   -
94 90 // Other methods
95 91 bool has_response(void) const;
96 92 bool has_events(void) const;
... ... @@ -128,8 +124,7 @@ public:
128 124 void deadc(const double& deadc);
129 125 void eventfile(const std::string& filename);
130 126 const std::string& eventfile(void) const;
131   - const std::string& eventtype(void) const;
132   - void dispose_events(void);
  127 + std::string eventtype(void) const;
133 128 const double& lo_user_thres(void) const;
134 129 const double& hi_user_thres(void) const;
135 130 void n_tels(const int& tels);
... ... @@ -142,12 +137,10 @@ protected:
142 137 void free_members(void);
143 138 void read_attributes(const GFitsHDU& hdu);
144 139 void write_attributes(GFitsHDU& hdu) const;
145   - void set_event_type(void);
146 140  
147 141 // Protected members
148 142 std::string m_instrument; //!< Instrument name
149 143 std::string m_eventfile; //!< Event filename
150   - std::string m_eventtype; //!< Event type (for XML file)
151 144 GCTAResponse* m_response; //!< Pointer to instrument response functions
152 145 GCTAPointing m_pointing; //!< Pointing direction
153 146 int m_obs_id; //!< Observation ID
... ... @@ -449,18 +442,6 @@ const std::string&amp; GCTAObservation::eventfile(void) const
449 442  
450 443  
451 444 /***********************************************************************//**
452   - * @brief Return type of events
453   - *
454   - * @return Type of events.
455   - ***************************************************************************/
456   -inline
457   -const std::string& GCTAObservation::eventtype(void) const
458   -{
459   - return m_eventtype;
460   -}
461   -
462   -
463   -/***********************************************************************//**
464 445 * @brief Return user low energy threshold
465 446 *
466 447 * @return User low energy threshold.
... ...
inst/cta/pyext/GCTAObservation.i
  View file @a16cd0d
... ... @@ -59,10 +59,6 @@ public:
59 59 virtual void read(const GXmlElement& xml);
60 60 virtual void write(GXmlElement& xml) const;
61 61  
62   - // Overwrite virtual base class methods
63   - virtual const GEvents* events(void) const;
64   - virtual void events(const GEvents& events);
65   -
66 62 // Other methods
67 63 bool has_response(void) const;
68 64 bool has_events(void) const;
... ... @@ -100,8 +96,7 @@ public:
100 96 void deadc(const double& deadc);
101 97 void eventfile(const std::string& filename);
102 98 const std::string& eventfile(void) const;
103   - const std::string& eventtype(void) const;
104   - void dispose_events(void);
  99 + std::string eventtype(void) const;
105 100 const double& lo_user_thres(void) const;
106 101 const double& hi_user_thres(void) const;
107 102 void n_tels(const int& tels);
... ...
inst/cta/src/GCTAObservation.cpp
  View file @a16cd0d
... ... @@ -35,7 +35,6 @@
35 35 #include "GGti.hpp"
36 36 #include "GTools.hpp"
37 37 #include "GIntegral.hpp"
38   -#include "GCTAException.hpp"
39 38 #include "GCTASupport.hpp"
40 39 #include "GCTAObservation.hpp"
41 40 #include "GCTAResponseIrf.hpp"
... ... @@ -64,7 +63,6 @@ const GObservationRegistry g_obs_veritas_registry(&amp;g_obs_veritas_seed);
64 63 #define G_WRITE "GCTAObservation::write(GXmlElement&)"
65 64 #define G_LOAD "GCTAObservation::load(std::string&, std::string&, "\
66 65 "std::string&, std::string&)"
67   -#define G_EVENTS "GCTAObservation::events()"
68 66  
69 67 /* __ Macros _____________________________________________________________ */
70 68  
... ... @@ -82,7 +80,7 @@ const GObservationRegistry g_obs_veritas_registry(&amp;g_obs_veritas_seed);
82 80 /***********************************************************************//**
83 81 * @brief Void constructor
84 82 *
85   - * Creates an empty CTA observation.
  83 + * Constructs an empty CTA observation.
86 84 ***************************************************************************/
87 85 GCTAObservation::GCTAObservation(void) : GObservation()
88 86 {
... ... @@ -99,7 +97,7 @@ GCTAObservation::GCTAObservation(void) : GObservation()
99 97 *
100 98 * @param[in] instrument Instrument name.
101 99 *
102   - * Creates an empty CTA observation for a given instrument. This enables
  100 + * Constructs an empty CTA observation for a given instrument. This enables
103 101 * using the CTA specific interface for any other VHE instrument. Note that
104 102 * each other VHE instruments needs a specific registry at the beginning
105 103 * of the GCTAObservation.cpp file. So far the following instruments are
... ... @@ -126,8 +124,8 @@ GCTAObservation::GCTAObservation(const std::string&amp; instrument) : GObservation()
126 124 * @param[in] psfcube Psf cube file name.
127 125 * @param[in] bkgcube Backgorund cube file name.
128 126 *
129   - * Creates a CTA observation from a counts cube, an exposure cube and a Psf
130   - * cube.
  127 + * Constructs a CTA observation from a counts cube, an exposure cube, a Psf
  128 + * cube and a background cube.
131 129 ***************************************************************************/
132 130 GCTAObservation::GCTAObservation(const std::string& cntcube,
133 131 const std::string& expcube,
... ... @@ -150,7 +148,7 @@ GCTAObservation::GCTAObservation(const std::string&amp; cntcube,
150 148 *
151 149 * @param[in] obs CTA observation.
152 150 *
153   - * Creates a CTA observation by copying an existing CTA observation.
  151 + * Constructs a CTA observation by copying an existing CTA observation.
154 152 ***************************************************************************/
155 153 GCTAObservation::GCTAObservation(const GCTAObservation& obs) : GObservation(obs)
156 154 {
... ... @@ -167,6 +165,8 @@ GCTAObservation::GCTAObservation(const GCTAObservation&amp; obs) : GObservation(obs)
167 165  
168 166 /***********************************************************************//**
169 167 * @brief Destructor
  168 + *
  169 + * Destructs CTA observation.
170 170 ***************************************************************************/
171 171 GCTAObservation::~GCTAObservation(void)
172 172 {
... ... @@ -223,6 +223,8 @@ GCTAObservation&amp; GCTAObservation::operator=(const GCTAObservation&amp; obs)
223 223  
224 224 /***********************************************************************//**
225 225 * @brief Clear CTA observation
  226 + *
  227 + * Clear CTA observation.
226 228 ***************************************************************************/
227 229 void GCTAObservation::clear(void)
228 230 {
... ... @@ -243,6 +245,8 @@ void GCTAObservation::clear(void)
243 245 * @brief Clone CTA observation
244 246 *
245 247 * @return Pointer to deep copy of CTA observation.
  248 + *
  249 + * Returns a pointer to a deep copy of a CTA observation.
246 250 ***************************************************************************/
247 251 GCTAObservation* GCTAObservation::clone(void) const
248 252 {
... ... @@ -255,28 +259,30 @@ GCTAObservation* GCTAObservation::clone(void) const
255 259 *
256 260 * @param[in] rsp Response function.
257 261 *
258   - * @exception GException::invalid_value
259   - * Specified response in not of type GCTAResponse.
  262 + * @exception GException::invalid_argument
  263 + * Invalid response class specified.
260 264 *
261   - * Sets the response function for the observation. The argument has to be of
262   - * type GCTAResponse, otherwise an exception is thrown.
  265 + * Sets the response function for the observation.
263 266 ***************************************************************************/
264 267 void GCTAObservation::response(const GResponse& rsp)
265 268 {
  269 + // Cast response dynamically
  270 + const GCTAResponse* ptr = dynamic_cast<const GCTAResponse*>(&rsp);
  271 +
  272 + // Throw exception if response is not of correct type
  273 + if (ptr == NULL) {
  274 + std::string cls = std::string(typeid(&rsp).name());
  275 + std::string msg = "Invalid response type \""+cls+"\" provided on "
  276 + "input. Please specify a \"GCTAResponse\" "
  277 + "object as argument.";
  278 + throw GException::invalid_argument(G_RESPONSE_SET, msg);
  279 + }
  280 +
266 281 // Free response
267 282 if (m_response != NULL) delete m_response;
268   - m_response = NULL;
269   -
270   - // Get pointer on CTA response
271   - const GCTAResponse* cta = dynamic_cast<const GCTAResponse*>(&rsp);
272   - if (cta == NULL) {
273   - std::string msg = "Specified response function is not a CTA "
274   - "response function.\n" + rsp.print();
275   - throw GException::invalid_value(G_RESPONSE_SET, msg);
276   - }
277 283  
278 284 // Clone response function
279   - m_response = cta->clone();
  285 + m_response = ptr->clone();
280 286  
281 287 // Return
282 288 return;
... ... @@ -291,15 +297,16 @@ void GCTAObservation::response(const GResponse&amp; rsp)
291 297 * @exception GException::invalid_value
292 298 * No valid response found in CTA observation.
293 299 *
294   - * Returns a pointer to the CTA response function. The pointer returned is
295   - * never NULL.
  300 + * Returns a pointer to the CTA response function. An exception is thrown if
  301 + * the pointer is not valid, hence the user does not need to verify the
  302 + * validity of the pointer.
296 303 ***************************************************************************/
297 304 const GCTAResponse* GCTAObservation::response(void) const
298 305 {
299 306 // Throw an exception if the response pointer is not valid
300 307 if (m_response == NULL) {
301   - std::string msg = "No valid response function found in CTA"
302   - " observation.\n";
  308 + std::string msg = "No valid response function found in CTA "
  309 + "observation.\n";
303 310 throw GException::invalid_value(G_RESPONSE_GET, msg);
304 311 }
305 312  
... ... @@ -380,287 +387,84 @@ void GCTAObservation::response(const GCTACubeExposure&amp; expcube,
380 387 /***********************************************************************//**
381 388 * @brief Get Region of Interest
382 389 *
  390 + * @return Region of Interest.
  391 + *
383 392 * @exception GException::invalid_value
384   - * Events is not an event list or unable to open an event file.
  393 + * Observation does not contain events.
  394 + * Observation does not contain an event list.
385 395 *
386   - * Extract the Region of Interest from the event list. If events are loaded
387   - * into the observation, the RoI is simply copied over. If no events are
388   - * present, the method attemps to read the RoI from the header keywords of
389   - * the event file. This obviously only works if a valid event file name is
390   - * available.
  396 + * Extracts the Region of Interest from the event list. An exception is
  397 + * thrown if the observation does not contain an event list.
391 398 ***************************************************************************/
392 399 GCTARoi GCTAObservation::roi(void) const
393 400 {
394   - // Initialise ROI
395   - GCTARoi roi;
396   -
397   - // If CTA has events then simply retrieve the ROI from the event list
398   - if (m_events != NULL) {
399   -
400   - // Extract region of interest from CTA observation
401   - const GCTAEventList* list = dynamic_cast<const GCTAEventList*>(m_events);
402   - if (list == NULL) {
403   - std::string msg = "The observation does not contain an event "
404   - "list. An event list is needed to retrieve "
405   - "the Region of Interest.";
406   - throw GException::invalid_value(G_ROI, msg);
407   - }
408   - roi = list->roi();
409   -
  401 + // Throw an exception if no events exist
  402 + if (m_events == NULL) {
  403 + std::string msg = "Region of Interest is not defined since the "
  404 + "observation does not contain events.";
  405 + throw GException::invalid_value(G_ROI, msg);
410 406 }
411 407  
412   - // ... otherwise try reading the ROI from the event file header
413   - else {
414   -
415   - // Try opening event file
416   - GFits fits;
417   - try {
418   - fits.open(m_eventfile);
419   - }
420   - catch (std::exception &e) {
421   - std::string msg;
422   - if (m_eventfile.length() == 0) {
423   - msg = "The observation does not contain any information about "
424   - "the event file, hence the Region of Interest cannot be "
425   - "determined.";
426   - }
427   - else {
428   - msg = "Could no open the event file \""+m_eventfile+"\". "
429   - "Please check that the file name is valid.";
430   - }
431   - throw GException::invalid_value(G_ROI, msg);
432   - }
433   -
434   - GFilename fname(m_eventfile);
435   -
436   - // Get event list HDU
437   - const GFitsTable& hdu = *fits.table(fname.extname("EVENTS"));
438   -
439   - // Read ROI from data selection keywords
440   - roi = gammalib::read_ds_roi(hdu);
441   -
442   - // Close FITS file
443   - fits.close();
444   -
  408 + // Get pointer to event list. Throw an exception if no event list is found
  409 + const GCTAEventList* list = dynamic_cast<const GCTAEventList*>(m_events);
  410 + if (list == NULL) {
  411 + std::string msg = "Region of Interest is not defined since the "
  412 + "observation does not contain an event list.";
  413 + throw GException::invalid_value(G_ROI, msg);
445 414 }
446 415  
447 416 // Return ROI
448   - return roi;
  417 + return (list->roi());
449 418 }
450 419  
451 420  
452 421 /***********************************************************************//**
453 422 * @brief Get Good Time Intervals
454 423 *
455   - * @return Good Time Interval
  424 + * @return Good Time Intervals.
456 425 *
457 426 * @exception GException::invalid_value
458   - * Unable to open the event file or to extract GTI information
459   - * from file.
460   - *
461   - * Extract the Good Time Intervals from the events. If events are loaded into
462   - * the observation, the GTIs are simply copied over. If no events are present
463   - * the method attemps to read the GTIs first from a GTI file, and if no GTI
464   - * file is given, by building a GTI from the event TSTART and TSTOP values.
465   - * The latter obviously only works if a valid event file name is available.
  427 + * Observation does not contain events.
  428 + *
  429 + * Extracts the Good Time Intervals from the events. An exception is thrown
  430 + * if the observation does not contain events.
466 431 ***************************************************************************/
467 432 GGti GCTAObservation::gti(void) const
468 433 {
469   - // Initialise GTIs
470   - GGti gti;
471   -
472   - // If CTA has events then simply retrieve the GTIs from the events
473   - if (m_events != NULL) {
474   -
475   - // Copy over GTIs
476   - gti = m_events->gti();
477   -
478   - }
479   -
480   - // ... otherwise try reading the GTIs from file
481   - else {
482   -
483   - // Extract the GTIs on basis on the DSS keywords, and if these
484   - // keywords are not present, we build the GTIs from the start and
485   - // stop time that is provided in the event header
486   - if (!m_eventfile.empty()) {
487   -
488   - // Initialise filename
489   - GFilename fname(m_eventfile);
490   -
491   - // Open FITS file
492   - GFits fits(fname.filename());
493   -
494   - // Get extension name
495   - std::string extname = fname.extname("EVENTS");
496   -
497   - // If the FITS file contains the specified extension then
498   - // first search for DSS keywords, and if they were not found,
499   - // use TSTART and TSTOP to build the GTI
500   - if (fits.contains(extname)) {
501   -
502   - // Get event list HDU
503   - const GFitsTable& events = *fits.table(extname);
504   -
505   - // Read GTI extension name from data sub-space keyword
506   - std::string gti_extname = gammalib::read_ds_gti_extname(events);
507   -
508   - // If no GTI extension name was found then
509   - if (gti_extname.empty()) {
510   - gti_extname = "GTI";
511   - }
512   -
513   - // If GTI extension is present in FITS file then read Good
514   - // Time Intervals from that extension
515   - if (fits.contains(gti_extname)) {
516   - const GFitsTable& hdu = *fits.table(gti_extname);
517   - gti.read(hdu);
518   - }
519   -
520   - // ... otherwise build GTI from TSTART and TSTOP
521   - else {
522   -
523   - // Read start and stop time
524   - double tstart = events.real("TSTART");
525   - double tstop = events.real("TSTOP");
526   -
527   - // Create time reference from header information
528   - GTimeReference timeref(events);
529   -
530   - // Set GTI time reference
531   - gti.reference(timeref);
532   -
533   - // Set start and stop time
534   - GTime start(tstart, gti.reference());
535   - GTime stop(tstop, gti.reference());
536   -
537   - // Append start and stop time as single time interval to GTI
538   - gti.append(start, stop);
539   -
540   - } // endelse: GTI built from TSTART and TSTOP
541   -
542   - // Close FITS file
543   - fits.close();
544   -
545   - } // endif: FITS files contained event extension
546   -
547   - // ... otherwise throw an exception
548   - else {
549   - std::string msg = "No event extension \""+extname+"\" found "
550   - "in file \""+fname.filename()+"\". Please "
551   - "specify a file name with a valid event "
552   - "extension name.";
553   - throw GException::invalid_value(G_GTI, msg);
554   - }
555   -
556   - } // endif: event file name was not empty
557   -
558   - // ... else throw an exception
559   - else {
560   - std::string msg = "The observation does not contain any "
561   - "information about the event file or the GTI "
562   - "file, hence the Good Time Intervals cannot "
563   - "be determined.";
564   - throw GException::invalid_value(G_GTI, msg);
565   - } // endelse: filenames were empty
566   -
567   - } // endelse: GTI had to be loaded
568   -
569   - // Check for GTI validity
570   - if (gti.size() == 0) {
571   - std::string msg = "Could not find Good Time Interval information "
572   - "in the event file \""+m_eventfile+"\". Please "
573   - "check that Good Time Interval information is "
574   - "properly set in the event file.";
  434 + // Throw an exception if no events exist
  435 + if (m_events == NULL) {
  436 + std::string msg = "Good Time Intervals are not defined since the "
  437 + "observation does not contain events.";
575 438 throw GException::invalid_value(G_GTI, msg);
576 439 }
577 440  
578   - // Return GTIs
579   - return gti;
  441 + // Return GTI
  442 + return (m_events->gti());
580 443 }
581 444  
582 445  
583 446 /***********************************************************************//**
584 447 * @brief Get energy boundaries
585 448 *
  449 + * @return Energy boundaries.
  450 + *
586 451 * @exception GException::invalid_value
587   - * Unable to open the event file or to extract energy boundaries
588   - * from file.
589   - *
590   - * Extract the energy boundaries from the events. If events are loaded into
591   - * the observation, the energy boundaries are simply copied over. If no
592   - * events are present the method attemps to read the energy boundaries from
593   - * the event file. This obviously only works if a valid event file name is
594   - * available.
  452 + * Observation does not contain events.
  453 + *
  454 + * Extract the energy boundaries from the events. An exception is thrown if
  455 + * the observation does not contain events.
595 456 ***************************************************************************/
596 457 GEbounds GCTAObservation::ebounds(void) const
597 458 {
598   - // Initialise energy boundaries
599   - GEbounds ebounds;
600   -
601   - // If CTA has events then simply retrieve the energy boundaries from the
602   - // events
603   - if (m_events != NULL) {
604   -
605   - // Copy over energy boundaries
606   - ebounds = m_events->ebounds();
607   -
608   - }
609   -
610   - // ... otherwise try reading the energy boundaries from the event file
611   - else {
612   -
613   - // Try opening event file
614   - GFits fits;
615   - try {
616   - fits.open(m_eventfile);
617   - }
618   - catch (std::exception &e) {
619   - std::string msg;
620   - if (m_eventfile.length() == 0) {
621   - msg = "The observation does not contain any information about "
622   - "the event file, hence the energy boundaries cannot be "
623   - "determined.";
624   - }
625   - else {
626   - msg = "Could not open the event file \""+m_eventfile+"\". "
627   - "Please check that the file name is valid.";
628   - }
629   - throw GException::invalid_value(G_EBOUNDS, msg);
630   - }
631   -
632   - GFilename fname(m_eventfile);
633   -
634   - // If file contains an EBOUNDS extension then load the energy
635   - // boundaries from that extension
636   - if (fits.contains("EBOUNDS")) {
637   - const GFitsTable& hdu = *fits.table("EBOUNDS");
638   - ebounds.read(hdu);
639   - }
640   -
641   - // ... otherwise, if file contains an EVENTS extension then load the
642   - // energy boundaries from the data selection keywords
643   - else if (fits.contains(fname.extname("EVENTS"))) {
644   - const GFitsTable& hdu = *fits.table(fname.extname("EVENTS"));
645   - ebounds = gammalib::read_ds_ebounds(hdu);
646   - }
647   -
648   - // ... else throw an exception
649   - else {
650   - std::string msg = "Could no find energy boundary information in "
651   - "the event file \""+m_eventfile+"\". Please "
652   - "check that energy boundary information is "
653   - "properly set in the event file.";
654   - throw GException::invalid_value(G_EBOUNDS, msg);
655   - }
656   -
657   - // Close FITS file
658   - fits.close();
659   -
  459 + // Throw an exception if no events exist
  460 + if (m_events == NULL) {
  461 + std::string msg = "Energy boundaries are not defined since the "
  462 + "observation does not contain events.";
  463 + throw GException::invalid_value(G_EBOUNDS, msg);
660 464 }
661 465  
662 466 // Return energy boundaries
663   - return ebounds;
  467 + return (m_events->ebounds());
664 468 }
665 469  
666 470  
... ... @@ -726,6 +530,7 @@ GEbounds GCTAObservation::ebounds(void) const
726 530 * ...
727 531 * <parameter name="ExposureCube" file="..."/>
728 532 * <parameter name="PsfCube" file="..."/>
  533 + * <parameter name="BkgCube" file="..."/>
729 534 * </observation>
730 535 *
731 536 * Optional user energy thresholds can be specified by adding the @a emin
... ... @@ -775,19 +580,22 @@ void GCTAObservation::read(const GXmlElement&amp; xml)
775 580 (par->attribute("name") == "CountsCube")) {
776 581  
777 582 // Store event type
778   - m_eventtype = par->attribute("name");
  583 + //m_eventtype = par->attribute("name");
779 584  
780 585 // Read eventlist file name
781 586 std::string filename = par->attribute("file");
782 587  
  588 + // Load events
  589 + load(filename);
  590 +
783 591 // Initialise filename
784   - GFilename fname(filename);
  592 + //GFilename fname(filename);
785 593  
786 594 // Open FITS file
787   - GFits fits(fname.filename());
  595 + //GFits fits(fname.filename());
788 596  
789 597 // Get extension name
790   - std::string extname = fname.extname("EVENTS");
  598 + //std::string extname = fname.extname("EVENTS");
791 599  
792 600 // Read observation attributes. Note that we do not load the
793 601 // events and GTI here as this will be done on-the-fly when
... ... @@ -796,6 +604,7 @@ void GCTAObservation::read(const GXmlElement&amp; xml)
796 604 // GCTAObservation to access th GTIs without accessing the
797 605 // events. It would be better if the GCTAEventList class handles
798 606 // the loading of the events on request.
  607 + /*
799 608 if (fits.contains(extname)) {
800 609 const GFitsHDU& hdu = *fits.at(extname);
801 610 read_attributes(hdu);
... ... @@ -804,12 +613,13 @@ void GCTAObservation::read(const GXmlElement&amp; xml)
804 613 const GFitsHDU& hdu = *fits.at(0);
805 614 read_attributes(hdu);
806 615 }
  616 + */
807 617  
808 618 // Close FITS file
809   - fits.close();
  619 + //fits.close();
810 620  
811 621 // Store event filename
812   - m_eventfile = filename;
  622 + //m_eventfile = filename;
813 623  
814 624 // Increment parameter counter
815 625 n_evfile++;
... ... @@ -1024,6 +834,7 @@ void GCTAObservation::read(const GXmlElement&amp; xml)
1024 834 * ...
1025 835 * <parameter name="ExposureCube" file="..."/>
1026 836 * <parameter name="PsfCube" file="..."/>
  837 + * <parameter name="BkgCube" file="..."/>
1027 838 * </observation>
1028 839 *
1029 840 * If user energy thresholds are defined (i.e. threshold values are >0) the
... ... @@ -1035,18 +846,17 @@ void GCTAObservation::read(const GXmlElement&amp; xml)
1035 846 ***************************************************************************/
1036 847 void GCTAObservation::write(GXmlElement& xml) const
1037 848 {
1038   - // Throw an exception if m_eventtype is neither "EventList" nor
1039   - // "CountsCube"
1040   - if ((m_eventtype != "EventList") &&
1041   - (m_eventtype != "CountsCube")) {
  849 + // Throw an exception if eventtype() is neither "EventList" nor
  850 + std::string evttype = eventtype();
  851 + if ((evttype != "EventList") && (evttype != "CountsCube")) {
1042 852 std::string msg;
1043   - if (m_eventtype.length() == 0) {
  853 + if (evttype.empty()) {
1044 854 msg = "The observation does not contain any events, hence "
1045 855 "it cannot be written into an XML file.";
1046 856 }
1047 857 else {
1048 858 msg = "The observation contains an unknown event type \""+
1049   - m_eventtype+"\". The event type needs to be either "
  859 + evttype+"\". The event type needs to be either "
1050 860 "\"EventList\" or \"CountsCube\".";
1051 861 }
1052 862 throw GException::invalid_value(G_WRITE, msg);
... ... @@ -1065,7 +875,7 @@ void GCTAObservation::write(GXmlElement&amp; xml) const
1065 875 if (!m_eventfile.empty()) {
1066 876  
1067 877 // Write event file name
1068   - GXmlElement* par = gammalib::xml_need_par(G_WRITE, xml, m_eventtype);
  878 + GXmlElement* par = gammalib::xml_need_par(G_WRITE, xml, evttype);
1069 879 par->attribute("file", m_eventfile);
1070 880  
1071 881 }
... ... @@ -1104,6 +914,15 @@ void GCTAObservation::write(GXmlElement&amp; xml) const
1104 914 * @brief Read data from FITS file
1105 915 *
1106 916 * @param[in] fits FITS file.
  917 + *
  918 + * Reads event data from a FITS file and sets the observation attributes.
  919 + *
  920 + * The method automatically switches between an event list and a counts cube,
  921 + * depending of the information provided in the FITS file. If an extension
  922 + * name is specified, the method checks whether the extension exists and
  923 + * loads the extension as event list. Otherwise, it checks whether the file
  924 + * contains an "EVENTS" extension and loads the extension as event list.
  925 + * If none of the above are satistified, the method loads a counts cube.
1107 926 ***************************************************************************/
1108 927 void GCTAObservation::read(const GFits& fits)
1109 928 {
... ... @@ -1155,9 +974,6 @@ void GCTAObservation::read(const GFits&amp; fits)
1155 974  
1156 975 }
1157 976  
1158   - // Set the event type
1159   - set_event_type();
1160   -
1161 977 // Return
1162 978 return;
1163 979 }
... ... @@ -1352,97 +1168,46 @@ void GCTAObservation::save(const std::string&amp; filename, const bool&amp; clobber) con
1352 1168  
1353 1169  
1354 1170 /***********************************************************************//**
1355   - * @brief Set event container
  1171 + * @brief Return event type string
1356 1172 *
1357   - * @param[in] events Event container.
  1173 + * @return Event type string.
1358 1174 *
1359   - * Set the event container for this observation by cloning the @p events.
  1175 + * Returns "EventList" if the observation contains an event list,
  1176 + * "CountsCube" if it contains a counts cube, "Events" if it container an
  1177 + * unknown type of events (which should never occur), and an empty string if
  1178 + * no events have been allocated.
1360 1179 ***************************************************************************/
1361   -void GCTAObservation::events(const GEvents& events)
  1180 +std::string GCTAObservation::eventtype(void) const
1362 1181 {
1363   - // Remove an existing event container
1364   - if (m_events != NULL) delete m_events;
1365   -
1366   - // Clone events
1367   - m_events = events.clone();
1368   -
1369   - // Set event type
1370   - set_event_type();
1371   -
1372   - // Return
1373   - return;
1374   -}
1375   -
  1182 + // Initialise empty event type string
  1183 + std::string eventtype = "";
1376 1184  
1377   -/***********************************************************************//**
1378   - * @brief Return pointer to event container
1379   - *
1380   - * @return Pointer to event container (always valid).
1381   - *
1382   - * @exception GException::invalid_value
1383   - * Events could not be fetched from disk.
1384   - *
1385   - * Returns a pointer to the event container. If no events exist the method
1386   - * tries to fetch events from the file specified by the m_eventfile member.
1387   - * If this fails, an exception will be thrown.
1388   - ***************************************************************************/
1389   -const GEvents* GCTAObservation::events(void) const
1390   -{
1391   - // If event container is offline then bring it online now
1392   - if (m_events == NULL) {
  1185 + // Continue only if events are allocated
  1186 + if (m_events != NULL) {
1393 1187  
1394   - // Try loading the events from FITS file. Catch any exception. Put
1395   - // the code into a critical zone as it might be called from within
1396   - // a parallelized thread.
1397   - #pragma omp critical
1398   - {
1399   - try {
1400   - const_cast<GCTAObservation*>(this)->load(m_eventfile);
1401   - }
1402   - catch (...) {
1403   - ;
1404   - }
  1188 + // Case A: we have a list
  1189 + GCTAEventList* list = dynamic_cast<GCTAEventList*>(m_events);
  1190 + if (list != NULL) {
  1191 + eventtype = "EventList";
1405 1192 }
1406 1193  
1407   - // Throw an exception if the event container is still not valid
1408   - if (m_events == NULL) {
1409   - std::string msg;
1410   - if (m_eventfile.length() == 0) {
1411   - msg = "The observation does not contain any information about "
1412   - "the event file, hence the events cannot be fetched "
1413   - "from disk for access.";
  1194 + // Case B: we have a cube
  1195 + else {
  1196 + GCTAEventCube* cube = dynamic_cast<GCTAEventCube*>(m_events);
  1197 + if (cube != NULL) {
  1198 + eventtype = "CountsCube";
1414 1199 }
  1200 +
  1201 + // Case C: we don't know what we have
1415 1202 else {
1416   - msg = "Could no load the event file \""+m_eventfile+"\" into "
1417   - "the observation. Please check that the file name is "
1418   - "valid and that the file contains events.";
  1203 + eventtype = "Events";
1419 1204 }
1420   - throw GException::invalid_value(G_EVENTS, msg);
1421 1205 }
1422 1206  
1423   - }
1424   -
1425   - // Return pointer to event container
1426   - return m_events;
1427   -}
1428   -
  1207 + } // endif: events were allocated
1429 1208  
1430   -/***********************************************************************//**
1431   - * @brief Dispose events
1432   - *
1433   - * Drops the events from the observation. Be careful with using this method
1434   - * as the events are not saved before being disposed.
1435   - ***************************************************************************/
1436   -void GCTAObservation::dispose_events(void)
1437   -{
1438   - // Delete any existing event container
1439   - if (m_events != NULL) delete m_events;
1440   -
1441   - // Signal that we disposed the events
1442   - m_events = NULL;
1443   -
1444   - // Return
1445   - return;
  1209 + // Return event type
  1210 + return eventtype;
1446 1211 }
1447 1212  
1448 1213  
... ... @@ -1468,7 +1233,7 @@ std::string GCTAObservation::print(const GChatter&amp; chatter) const
1468 1233 result.append("\n"+gammalib::parformat("Identifier")+id());
1469 1234 result.append("\n"+gammalib::parformat("Instrument")+instrument());
1470 1235 result.append("\n"+gammalib::parformat("Event file")+eventfile());
1471   - result.append("\n"+gammalib::parformat("Event type")+m_eventtype);
  1236 + result.append("\n"+gammalib::parformat("Event type")+eventtype());
1472 1237 result.append("\n"+gammalib::parformat("Statistics")+statistics());
1473 1238 result.append("\n"+gammalib::parformat("Ontime"));
1474 1239 result.append(gammalib::str(ontime())+" s");
... ... @@ -1544,7 +1309,6 @@ void GCTAObservation::init_members(void)
1544 1309 m_instrument = "CTA";
1545 1310 m_object.clear();
1546 1311 m_eventfile.clear();
1547   - m_eventtype.clear();
1548 1312 m_bgdfile.clear();
1549 1313 m_response = NULL;
1550 1314 m_pointing.clear();
... ... @@ -1574,7 +1338,6 @@ void GCTAObservation::copy_members(const GCTAObservation&amp; obs)
1574 1338 m_instrument = obs.m_instrument;
1575 1339 m_object = obs.m_object;
1576 1340 m_eventfile = obs.m_eventfile;
1577   - m_eventtype = obs.m_eventtype;
1578 1341 m_bgdfile = obs.m_bgdfile;
1579 1342 m_pointing = obs.m_pointing;
1580 1343 m_obs_id = obs.m_obs_id;
... ... @@ -1759,46 +1522,3 @@ void GCTAObservation::write_attributes(GFitsHDU&amp; hdu) const
1759 1522 // Return
1760 1523 return;
1761 1524 }
1762   -
1763   -
1764   -/***********************************************************************//**
1765   - * @brief Set event type
1766   - *
1767   - * Depending on the type of events that are stored in the m_events member,
1768   - * set the event type to
1769   - *
1770   - * "EventList" if m_events is of type GCTAEventList
1771   - * "CountsCube" if m_events is of type GCTAEventCube
1772   - * "Events" if m_events is not NULL but neither GCTAEventList nor GCTAEventCube
1773   - * "" if m_events is NULL
1774   - *
1775   - ***************************************************************************/
1776   -void GCTAObservation::set_event_type(void)
1777   -{
1778   - // Continue only if events are loaded
1779   - if (m_events != NULL) {
1780   -
1781   - // Case A: we have a list
1782   - GCTAEventList* list = dynamic_cast<GCTAEventList*>(m_events);
1783   - if (list != NULL) {
1784   - m_eventtype = "EventList";
1785   - }
1786   -
1787   - // Case B: we have a cube
1788   - else {
1789   - GCTAEventCube* cube = dynamic_cast<GCTAEventCube*>(m_events);
1790   - if (cube != NULL) {
1791   - m_eventtype = "CountsCube";
1792   - }
1793   -
1794   - // Case C: we don't know what we have
1795   - else {
1796   - m_eventtype = "Events";
1797   - }
1798   - }
1799   -
1800   - }
1801   -
1802   - // Return
1803   - return;
1804   -}
... ...