diff --git a/cscripts/csinfo.py b/cscripts/csinfo.py index ff73255..5c8852e 100755 --- a/cscripts/csinfo.py +++ b/cscripts/csinfo.py @@ -1,6 +1,6 @@ #!/usr/bin/env python # ========================================================================== -# Print info about Gammalib / ctools to the console. +# Print info about Gammalib / ctools to the console # # Copyright (C) 2015-2016 Christoph Deil # @@ -87,8 +87,8 @@ def get_command_output(cmd): output : str Command output """ - # Try to import getoutput from subprocess, alternatively get it from - # commands + # Try to import getoutput from the "subprocess" module, alternatively get + # it from the "commands" module try: from subprocess import getoutput except ImportError: @@ -287,8 +287,7 @@ def get_pkg_config_info(info, library): # Execute command out = get_command_output(cmd) - # If the command returns 'was not found' then pkg-config is not - # available + # If the command returns "not found" then pkg-config is not available if 'not found' in out: out = 'Not available' @@ -322,8 +321,8 @@ def csinfo(argv): """ Print information about Gammalib and ctools to the console """ - # If there are no command line arguments then print help and exit with - # success + # If there are no command line arguments then print help text and exit + # with success if len(argv) <= 1: csinfo_print_help() sys.exit(0) diff --git a/cscripts/csobsdef.py b/cscripts/csobsdef.py index ea5d87e..e715957 100755 --- a/cscripts/csobsdef.py +++ b/cscripts/csobsdef.py @@ -28,7 +28,7 @@ import ctools # ============== # class csobsdef(ctools.cscript): """ - Creates an observation definition file from a pointing list. + Creates an observation definition file from a pointing list The pointing list is a comma-separated value ASCII file with header keywords in the first row followed by a list of pointings (one @@ -92,8 +92,7 @@ class csobsdef(ctools.cscript): self._pntdef = gammalib.GCsv() self._tmin = 0.0 - # Initialise application by calling the appropriate class - # constructor. + # Initialise application by calling the appropriate class constructor self._init_cscript(argv) # Return @@ -128,6 +127,20 @@ class csobsdef(ctools.cscript): The method creates an XML element so that that the response XML writer will write the database and response name into the observation definition file. + + Parameters + ---------- + obs : `~gammalib.GCTAObservation` + CTA observation + caldb : str + Calibration database + irf : str + Instrument response function + + Returns + obs : `~gammalib.GCTAObservation` + CTA observation with response attached + ------- """ # Create XML element xml = gammalib.GXmlElement() @@ -156,7 +169,7 @@ class csobsdef(ctools.cscript): Raises ------ RuntimeError - Invalid pointing definition file format. + Invalid pointing definition file format """ # Switch screen logging on in debug mode if self._logDebug(): @@ -170,21 +183,25 @@ class csobsdef(ctools.cscript): self._log('\n') self._log.header1('Creating observation definition XML file') - # Load pointing definition file if it is not already set + # Load pointing definition file if it is not already set. Extract + # the number of columns and pointings if self._pntdef.size() == 0: self._pntdef = gammalib.GCsv(self['inpnt'].filename(), ',') ncols = self._pntdef.ncols() npnt = self._pntdef.nrows()-1 - # Throw an exception is there is no header information + # Raise an exception if there is no header information if self._pntdef.nrows() < 1: raise RuntimeError('No header found in pointing definition file.') # Clear observation container self._obs.clear() + + # Initialise observation identifier counter identifier = 1 - # Extract header from pointing definition file + # Extract header columns from pointing definition file and put them + # into a list header = [] for col in range(ncols): header.append(self._pntdef[0,col]) @@ -192,28 +209,31 @@ class csobsdef(ctools.cscript): # Loop over all pointings for pnt in range(npnt): - # Set row index + # Set pointing definition CSV file row index row = pnt + 1 - # Create CTA observation + # Create empty CTA observation obs = gammalib.GCTAObservation() - # Set observation name + # Set observation name. If no observation name was given then + # use "None". if 'name' in header: name = self._pntdef[row, header.index('name')] else: name = 'None' obs.name(name) - # Set identifier + # Set observation identifier. If no observation identified was + # given the use the internal counter. if 'id' in header: - id_ = self._pntdef[row, header.index('id')] + obsid = self._pntdef[row, header.index('id')] else: - id_ = '%6.6d' % identifier + obsid = '%6.6d' % identifier identifier += 1 - obs.id(id_) + obs.id(obsid) - # Set pointing + # Set pointing. Either use "ra" and "dec" or "lon" and "lat". + # If none of these pairs are given then raise an exception. if 'ra' in header and 'dec' in header: ra = float(self._pntdef[row, header.index('ra')]) dec = float(self._pntdef[row, header.index('dec')]) @@ -229,7 +249,8 @@ class csobsdef(ctools.cscript): 'found in pointing definition file.') obs.pointing(gammalib.GCTAPointing(pntdir)) - # Set response function + # Set response function. If no "caldb" or "irf" information is + # provided then use the user parameter values. if 'caldb' in header: caldb = self._pntdef[row, header.index('caldb')] else: @@ -241,14 +262,16 @@ class csobsdef(ctools.cscript): if caldb != '' and irf != '': obs = self._set_response(obs, caldb, irf) - # Set deadtime correction factor + # Set deadtime correction factor. If no information is provided + # then use the user parameter value "deadc". if 'deadc' in header: deadc = float(self._pntdef[row, header.index('deadc')]) else: deadc = self['deadc'].real() obs.deadc(deadc) - # Set Good Time Interval + # Set Good Time Interval. If no information is provided then use + # the user parameter value "duration". if 'duration' in header: duration = float(self._pntdef[row, header.index('duration')]) else: @@ -263,7 +286,9 @@ class csobsdef(ctools.cscript): obs.ontime(gti.ontime()) obs.livetime(gti.ontime()*deadc) - # Set Energy Boundaries + # Set Energy Boundaries. If no "emin" or "emax" information is + # provided then use the user parameter values in case they are + # valid. has_emin = False has_emax = False if 'emin' in header: @@ -285,7 +310,8 @@ class csobsdef(ctools.cscript): ebounds = gammalib.GEbounds(gammalib.GEnergy(emin, 'TeV'), gammalib.GEnergy(emax, 'TeV')) - # Set ROI + # Set ROI. If no ROI radius is provided then use the user + # parameters "rad". has_roi = False if 'rad' in header: rad = float(self._pntdef[row, header.index('rad')]) @@ -298,17 +324,17 @@ class csobsdef(ctools.cscript): roi = gammalib.GCTARoi(gammalib.GCTAInstDir(pntdir), rad) # Create an empty event list - list_ = gammalib.GCTAEventList() - list_.gti(gti) + event_list = gammalib.GCTAEventList() + event_list.gti(gti) - # Set optional information + # If available, set the energy boundaries and the ROI if has_ebounds: - list_.ebounds(ebounds) + event_list.ebounds(ebounds) if has_roi: - list_.roi(roi) + event_list.roi(roi) # Attach event list to CTA observation - obs.events(list_) + obs.events(event_list) # Write observation into logger if self._logExplicit(): diff --git a/cscripts/obsutils.py b/cscripts/obsutils.py index af2c9c4..adc3f19 100755 --- a/cscripts/obsutils.py +++ b/cscripts/obsutils.py @@ -136,15 +136,11 @@ def sim(obs, log=False, debug=False, chatter=2, edisp=False, seed=0, # cubes and append them to the observation if len(sim.obs()) > 1: - # First get the cube. We need this because ctbin builds a cube on - # the fly. We should avoid this and make the cube part of ctbin - cube = bin.cube() - # Get xref, yref, coordsys and proj from counts cube - _xref = cube.counts().projection().crval(0) - _yref = cube.counts().projection().crval(1) - _coordsys = cube.counts().projection().coordsys() - _proj = cube.counts().projection().code() + _xref = bin.cube().counts().projection().crval(0) + _yref = bin.cube().counts().projection().crval(1) + _coordsys = bin.cube().counts().projection().coordsys() + _proj = bin.cube().counts().projection().code() # Kluge: map "EQU" to "CEL" coordinate system (not sure why the # GSkyProjection::coordsys returns "EQU" instead of "CEL" diff --git a/examples/show_butterfly.py b/examples/show_butterfly.py index 09af4d5..006a4c7 100755 --- a/examples/show_butterfly.py +++ b/examples/show_butterfly.py @@ -25,7 +25,7 @@ try: import matplotlib.pyplot as plt plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() diff --git a/examples/show_model.py b/examples/show_model.py index b3f0e70..d2b31c4 100755 --- a/examples/show_model.py +++ b/examples/show_model.py @@ -25,7 +25,7 @@ try: import matplotlib.pyplot as plt plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() diff --git a/examples/show_obs.py b/examples/show_obs.py index 4fbcbac..daaa160 100755 --- a/examples/show_obs.py +++ b/examples/show_obs.py @@ -25,7 +25,7 @@ try: import matplotlib.pyplot as plt plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() diff --git a/examples/show_pha.py b/examples/show_pha.py index afb34ae..2c42866 100755 --- a/examples/show_pha.py +++ b/examples/show_pha.py @@ -26,7 +26,7 @@ try: import matplotlib.pyplot as plt plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() diff --git a/examples/show_pull_evolution.py b/examples/show_pull_evolution.py index 6b5d008..8892b10 100755 --- a/examples/show_pull_evolution.py +++ b/examples/show_pull_evolution.py @@ -25,7 +25,7 @@ try: import matplotlib.pyplot as plt plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() try: diff --git a/examples/show_pull_histogram.py b/examples/show_pull_histogram.py index ebe3aee..7c777f4 100755 --- a/examples/show_pull_histogram.py +++ b/examples/show_pull_histogram.py @@ -26,7 +26,7 @@ try: import matplotlib.mlab as mlab plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() try: diff --git a/examples/show_response.py b/examples/show_response.py index eb7deda..533c910 100755 --- a/examples/show_response.py +++ b/examples/show_response.py @@ -26,7 +26,7 @@ try: import matplotlib.pyplot as plt plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() diff --git a/examples/show_sensitivity.py b/examples/show_sensitivity.py index 58b29f0..5913347 100755 --- a/examples/show_sensitivity.py +++ b/examples/show_sensitivity.py @@ -26,7 +26,7 @@ try: import matplotlib.pyplot as plt plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() diff --git a/examples/show_spectrum.py b/examples/show_spectrum.py index 0003db6..bf08c24 100755 --- a/examples/show_spectrum.py +++ b/examples/show_spectrum.py @@ -25,7 +25,7 @@ try: import matplotlib.pyplot as plt plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() diff --git a/examples/show_ts_distribution.py b/examples/show_ts_distribution.py index 8c23c11..0a97f36 100755 --- a/examples/show_ts_distribution.py +++ b/examples/show_ts_distribution.py @@ -26,7 +26,7 @@ try: import matplotlib.pyplot as plt plt.figure() plt.close() -except ImportError, RuntimeError: +except (ImportError, RuntimeError): print('This script needs the "matplotlib" module') sys.exit() try: diff --git a/pyext/ctbin.i b/pyext/ctbin.i index 8fcea66..22c3384 100644 --- a/pyext/ctbin.i +++ b/pyext/ctbin.i @@ -49,7 +49,7 @@ public: void save(void); void publish(const std::string& name = ""); const GObservations& obs(void) const; - GCTAEventCube cube(void) const; + const GCTAEventCube& cube(void) const; }; diff --git a/src/ctbin/ctbin.cpp b/src/ctbin/ctbin.cpp index aa56393..7f0f632 100644 --- a/src/ctbin/ctbin.cpp +++ b/src/ctbin/ctbin.cpp @@ -281,6 +281,10 @@ void ctbin::run(void) } // endfor: looped over observations + // Build event cube (needs to come before obs_cube() since this method + // relies on correct setting of m_cube) + m_cube = GCTAEventCube(m_counts, m_weights, m_ebounds, m_gti); + // Set a single cube in the observation container obs_cube(); @@ -400,6 +404,7 @@ void ctbin::init_members(void) m_weights.clear(); m_ebounds.clear(); m_gti.clear(); + m_cube.clear(); m_ontime = 0.0; m_livetime = 0.0; @@ -427,6 +432,7 @@ void ctbin::copy_members(const ctbin& app) m_weights = app.m_weights; m_ebounds = app.m_ebounds; m_gti = app.m_gti; + m_cube = app.m_cube; m_ontime = app.m_ontime; m_livetime = app.m_livetime; diff --git a/src/ctbin/ctbin.hpp b/src/ctbin/ctbin.hpp index 58098f0..138f6c0 100644 --- a/src/ctbin/ctbin.hpp +++ b/src/ctbin/ctbin.hpp @@ -81,7 +81,7 @@ public: void save(void); void publish(const std::string& name = ""); const GObservations& obs(void) const; - GCTAEventCube cube(void) const; + const GCTAEventCube& cube(void) const; protected: // Protected methods @@ -95,25 +95,28 @@ protected: std::vector cube_layer_usage(const GEbounds& ebounds) const; // User parameters - GFilename m_outcube; //!< Output counts map file name - bool m_usepnt; //!< Use pointing instead of xref/yref parameters - bool m_publish; //!< Publish counts cube? - GChatter m_chatter; //!< Chattiness + GFilename m_outcube; //!< Output counts map file name + bool m_usepnt; //!< Use pointing instead of xref/yref parameters + bool m_publish; //!< Publish counts cube? + GChatter m_chatter; //!< Chattiness // Protected members - GObservations m_obs; //!< Observation container - GSkyMap m_counts; //!< Event cube counts - GSkyMap m_weights; //!< Event cube weights - GEbounds m_ebounds; //!< Energy boundaries - GGti m_gti; //!< Good time intervals - double m_ontime; //!< Total ontime - double m_livetime; //!< Total livetime + GObservations m_obs; //!< Observation container + GSkyMap m_counts; //!< Event cube counts + GSkyMap m_weights; //!< Event cube weights + GEbounds m_ebounds; //!< Energy boundaries + GGti m_gti; //!< Good time intervals + GCTAEventCube m_cube; //!< Events cube (for cube() method) + double m_ontime; //!< Total ontime + double m_livetime; //!< Total livetime }; /***********************************************************************//** * @brief Return observation container * * @return Reference to observation container + * + * Returns a reference to an observation container. ***************************************************************************/ inline const GObservations& ctbin::obs(void) const @@ -125,18 +128,14 @@ const GObservations& ctbin::obs(void) const /***********************************************************************//** * @brief Return event cube * - * @return Event cube + * @return Reference to event cube * - * Returns a event cube built from the ctbin member information. + * Returns a reference to the event cube. ***************************************************************************/ inline -GCTAEventCube ctbin::cube(void) const +const GCTAEventCube& ctbin::cube(void) const { - // Build event cube - GCTAEventCube cube(m_counts, m_weights, m_ebounds, m_gti); - - // Return cube - return cube; + return m_cube; } #endif /* CTBIN_HPP */