diff --git a/.gitignore b/.gitignore index 255f396..a8defc0 100644 --- a/.gitignore +++ b/.gitignore @@ -12,6 +12,8 @@ *.toc *.Pi *.Tpi +*.gcno +*.gcda .DS_Store .cproject .project @@ -41,6 +43,7 @@ depcomp test-driver testconf.sh coverage-report-python.xml +coverage-report-cxx.xml doc/html doc/man diff --git a/Makefile.am b/Makefile.am index 41b9ad7..8054c73 100644 --- a/Makefile.am +++ b/Makefile.am @@ -93,7 +93,7 @@ cppncss: # Add coverage rule coverage: @if type gcovr >/dev/null 2>&1; then \ - gcovr -r . -x -o coverage-report-cxx.xml; \ + gcovr -r . -x --exclude-unreachable-branches -o coverage-report-cxx.xml; \ else \ echo "gcovr script not found in PATH."; \ fi; \ diff --git a/configure.ac b/configure.ac index 4b03958..405b508 100644 --- a/configure.ac +++ b/configure.ac @@ -87,6 +87,22 @@ AC_CANONICAL_HOST ############################################################################# +# Check for --enable-coverage option +############################################################################# +AC_MSG_CHECKING(whether to build ctools for code coverage measure) +AC_ARG_ENABLE([coverage], + AS_HELP_STRING([--enable-coverage], + [turn on code coverage measure [default=no]]), + [enable_coverage="$enableval"], + [enable_coverage="no"]) +if test "x$enable_coverage" = "xyes"; then + AC_MSG_RESULT(yes) +else + AC_MSG_RESULT(no) +fi + + +############################################################################# # Check for valid C++ compiler. This has to come before the check for # # header files, otherwise we won't detect if the C++ compiler is missing. # # # @@ -140,11 +156,23 @@ else fi AM_CONDITIONAL(IS_ELCAPITAN, test "x$osx_el_capitan" = "xyes") +# If coverage is selected then save current C++ options +if test "x$enable_coverage" = "xyes"; then + SAVED_CXXFLAGS=$CXXFLAGS + SAVED_CFLAGS=$CFLAGS +fi + # Our language is C++ AC_LANG(C++) AC_PROG_CXX($cxx_list) AC_PROG_CC($cc_list) +# If coverage is selected then restore former C++ options +if test "x$enable_coverage" = "xyes"; then + CXXFLAGS=$SAVED_CXXFLAGS + CFLAGS=$SAVED_CFLAGS +fi + ############################################################################# # Darwin option: Determine architecture of universal library to be built # @@ -758,25 +786,22 @@ fi ############################################################################# # Initialise Python coverage with empty string PYTHON_COVERAGE= +COVERAGE_CXXFLAGS= -# Check for coverage option -AC_MSG_CHECKING(whether to build ctools for code coverage measure) -AC_ARG_ENABLE([coverage], - AS_HELP_STRING([--enable-coverage], - [turn on code coverage measure [default=no]]), - [enable_coverage="$enableval"], - [enable_coverage="no"]) +# If coverage option has been selected then set some variables if test "x$enable_coverage" = "xyes"; then CXXFLAGS="--coverage $CXXFLAGS" + COVERAGE_CXXFLAGS="--coverage" LDFLAGS="--coverage $LDFLAGS" if test "x$has_python" = "xyes"; then - PYTHON_COVERAGE="coverage run --append " + PYTHON_COVERAGE="coverage run --branch --append " fi AC_MSG_RESULT(yes) else AC_MSG_RESULT(no) fi AC_SUBST(PYTHON_COVERAGE) +AC_SUBST(COVERAGE_CXXFLAGS) ############################################################################# @@ -937,7 +962,7 @@ OPENMP_LIB_LDFLAGS= #fi AC_SUBST(OPENMP_LIB_LDFLAGS) -# Signal if ctools cannot be compiled as gammalibe requires OpenMP but ctools +# Signal if ctools cannot be compiled as gammalib requires OpenMP but ctools # cannot be compiled with OpenMP if test "x$enable_openmp" = "xno" && test "x$gammalib_with_openmp" = "xyes"; then AC_MSG_ERROR([The GammaLib library has been built with openMP support but OpenMP is not available on this system. Please install OpenMP or re-build GammaLib without OpenMP support.]) diff --git a/cscripts/Makefile.am b/cscripts/Makefile.am index 7fbfd42..5400d1c 100644 --- a/cscripts/Makefile.am +++ b/cscripts/Makefile.am @@ -26,13 +26,17 @@ linkscripts = cspull cssens cstsdist cscaldb csresmap \ # Test scripts to be distributed in $(pyexecdir)/cscripts/tests testscripts = tests/__init__.py \ ../test/test_python_cscripts.py \ - ../test/test_cscaldb.py + ../test/test_cscaldb.py \ + ../test/test_cslightcrv.py # Test data to be distributed in $(pyexecdir)/cscripts/tests/data -testdata = +testdata = ../test/data/crab.xml \ + ../test/data/crab_events.fits.gz \ + ../test/data/crab_cntmap.fits.gz \ + ../test/data/exclusion.reg # Test IRFs to be distributed in $(pyexecdir)/cscripts/tests/irf -testirf = +#testirf = # Distribute scripts and parfiles EXTRA_DIST = $(binscripts) $(pythonscripts) @@ -46,8 +50,8 @@ testscriptsdir = $(pyexecdir)/cscripts/tests testscripts_DATA = $(testscripts) testdatadir = $(pyexecdir)/cscripts/tests/data testdata_DATA = $(testdata) -testirfdir = $(pyexecdir)/cscripts/tests/irf -testirf_DATA = $(testirf) +#testirfdir = $(pyexecdir)/cscripts/tests/irf +#testirf_DATA = $(testirf) # Add hook rule to create symbolic links for cscripts install-data-hook: diff --git a/cscripts/csiactobs.py b/cscripts/csiactobs.py index 3c74c5c..26dccd9 100755 --- a/cscripts/csiactobs.py +++ b/cscripts/csiactobs.py @@ -501,14 +501,14 @@ class csiactobs(ctools.cscript): self._log("\n") # Append observation to XML and set attributes - obs = lib.append("observation"); - obs.attribute("name", object_name); - obs.attribute("id", str(obs_id)); - obs.attribute("instrument", telescope); + obs = lib.append("observation") + obs.attribute("name", object_name) + obs.attribute("id", str(obs_id)) + obs.attribute("instrument", telescope) # Append event file ev = gammalib.GXmlElement("parameter name=\"EventList\"") - ev.attribute("file", eventfile); + ev.attribute("file", eventfile) # Append effective area aeff = gammalib.GXmlElement("parameter name=\"EffectiveArea\"") diff --git a/cscripts/cslightcrv.py b/cscripts/cslightcrv.py index 313bc55..76ff963 100755 --- a/cscripts/cslightcrv.py +++ b/cscripts/cslightcrv.py @@ -76,12 +76,15 @@ class cslightcrv(ctools.cscript): self._stacked = False # Initialise observation container from constructor arguments. - self._obs = self._set_input_obs(argv) + self._obs, argv = self._set_input_obs(argv) - # Initialise application by calling the appropriate class - # constructor. + # Initialise script by calling the appropriate class constructor. self._init_cscript(argv) + # Set logger properties + self._log_header() + self._log.date(True) + # Return return @@ -91,12 +94,10 @@ class cslightcrv(ctools.cscript): """ Get parameters from parfile. """ - # Set observation if not done before - if self._obs == None or self._obs.size() == 0: - self._require_inobs("cslightcrv::_get_parameters()") - self._obs = self._get_observations() + # Setup observations + self._setup_observations(self._obs) - # Set models if we have none + # Set models if there are none in the container if self._obs.models().size() == 0: self._obs.models(self["inmodel"].filename()) @@ -106,11 +107,18 @@ class cslightcrv(ctools.cscript): # Get time boundaries self._tbins = self._create_tbounds() - # Set stacked analysis flag if the requested number of energy bins - # is not zero; otherwise an unbinned analysis will be done. - self._stacked = False if self["enumbins"].integer() == 0 else True + # Set stacked analysis flag to True if the requested number of + # energy bins is positive. Otherwise an unbinned analysis will + # be done and the stacked analysis flag will be set to False. + if self["enumbins"].integer() > 0: + self._stacked = True + else: + self._stacked = False - # Make sure that other parameters are queried now and not later + # Make sure that remaining user parameters are queried now. We + # do not store the actual parameter values as we do not want + # too many instance attributes with enhances the maintenance + # costs. self["emin"].real() self["emax"].real() if self._stacked: @@ -429,16 +437,26 @@ class cslightcrv(ctools.cscript): if self._logExplicit(): self._log.header3("Creating point spread function cube") - # Create psf cube + # Compute spatial binning for point spread function and + # energy dispersion cubes + binsz = 0.1 * self["binsz"].real() + nxpix = self["nxpix"].integer() / 10 + nypix = self["nypix"].integer() / 10 + if nxpix < 2: + nxpix = 2 + if nypix < 2: + nypix = 2 + + # Create point spread function cube psfcube = ctools.ctpsfcube(obs) psfcube["incube"] = "NONE" psfcube["usepnt"] = False psfcube["ebinalg"] = "LOG" psfcube["xref"] = self["xref"].real() psfcube["yref"] = self["yref"].real() - psfcube["binsz"] = self["binsz"].real() - psfcube["nxpix"] = self["nxpix"].integer() - psfcube["nypix"] = self["nypix"].integer() + psfcube["binsz"] = binsz + psfcube["nxpix"] = nxpix + psfcube["nypix"] = nypix psfcube["enumbins"] = self["enumbins"].integer() psfcube["emin"] = self["emin"].real() psfcube["emax"] = self["emax"].real() @@ -453,16 +471,16 @@ class cslightcrv(ctools.cscript): if self._logExplicit(): self._log.header3("Creating energy dispersion cube") - # Create edisp cube + # Create energy dispersion cube edispcube = ctools.ctedispcube(obs) edispcube["incube"] = "NONE" edispcube["usepnt"] = False edispcube["ebinalg"] = "LOG" edispcube["xref"] = self["xref"].real() edispcube["yref"] = self["yref"].real() - edispcube["binsz"] = self["binsz"].real() - edispcube["nxpix"] = self["nxpix"].integer() - edispcube["nypix"] = self["nypix"].integer() + edispcube["binsz"] = binsz + edispcube["nxpix"] = nxpix + edispcube["nypix"] = nypix edispcube["enumbins"] = self["enumbins"].integer() edispcube["emin"] = self["emin"].real() edispcube["emax"] = self["emax"].real() diff --git a/cscripts/csobsinfo.py b/cscripts/csobsinfo.py index 10ef2cc..30a0b32 100755 --- a/cscripts/csobsinfo.py +++ b/cscripts/csobsinfo.py @@ -52,7 +52,7 @@ class csobsinfo(ctools.cscript): self._gti = gammalib.GGti() # Initialise observation container from constructor arguments. - self._obs = self._set_input_obs(argv) + self._obs, argv = self._set_input_obs(argv) # Initialise application by calling the appropriate class # constructor. diff --git a/cscripts/cspull.py b/cscripts/cspull.py index ec659ef..0f1e47d 100755 --- a/cscripts/cspull.py +++ b/cscripts/cspull.py @@ -18,7 +18,7 @@ # along with this program. If not, see . # # ========================================================================== -import gammalib +#import gammalib import ctools from cscripts import obsutils import sys @@ -62,7 +62,7 @@ class cspull(ctools.cscript): self._log_clients = False # Initialise observation container from constructor arguments. - self._obs = self._set_input_obs(argv) + self._obs, argv = self._set_input_obs(argv) # Initialise application by calling the appropriate class # constructor. @@ -410,17 +410,17 @@ class cspull(ctools.cscript): # Write out result immediately if seed == 0: - file = open(self._outfile.url(), 'w') - writer = csv.DictWriter(file, result['colnames']) + f = open(self._outfile.url(), 'w') + writer = csv.DictWriter(f, result['colnames']) headers = {} for n in result['colnames']: headers[n] = n writer.writerow(headers) else: - file = open(self._outfile.url(), 'a') - writer = csv.DictWriter(file, result['colnames']) + f = open(self._outfile.url(), 'a') + writer = csv.DictWriter(f, result['colnames']) writer.writerow(result['values']) - file.close() + f.close() # Return return diff --git a/cscripts/csresmap.py b/cscripts/csresmap.py index 8559c87..4197d23 100755 --- a/cscripts/csresmap.py +++ b/cscripts/csresmap.py @@ -74,7 +74,7 @@ class csresmap(ctools.cscript): self._log_clients = False # Initialise observation container from constructor arguments. - self._obs = self._set_input_obs(argv) + self._obs, argv = self._set_input_obs(argv) # Initialise application by calling the appropriate class # constructor. diff --git a/cscripts/cssens.py b/cscripts/cssens.py index e526b67..89e3e89 100755 --- a/cscripts/cssens.py +++ b/cscripts/cssens.py @@ -347,18 +347,18 @@ class cssens(ctools.cscript): photon_flux_value = [] energy_flux_value = [] sensitivity_value = [] - iter = 0 + iterations = 0 test_crab_flux = 0.1 # Initial test flux in Crab units (100 mCrab) # Loop until we break while True: # Update iteration counter - iter += 1 + iterations += 1 # Write header if self._logExplicit(): - self._log.header2("Iteration "+str(iter)) + self._log.header2("Iteration "+str(iterations)) # Set source model. crab_prefactor is the Prefactor that # corresponds to 1 Crab @@ -370,7 +370,7 @@ class cssens(ctools.cscript): self._obs.models(src_model) # Simulate events - sim = obsutils.sim(self._obs, nbins=self._enumbins, seed=iter, + sim = obsutils.sim(self._obs, nbins=self._enumbins, seed=iterations, binsz=self._binsz, npix=self._npix, log=self._log_clients, debug=self._debug, edisp=self._edisp) @@ -464,8 +464,8 @@ class cssens(ctools.cscript): # Get fitted Crab, photon and energy fluxes crab_flux = result_all[self._srcname]['Prefactor'].value() / \ crab_prefactor - crab_flux_err = result_all[self._srcname]['Prefactor'].error() / \ - crab_prefactor + #crab_flux_err = result_all[self._srcname]['Prefactor'].error() / \ + # crab_prefactor photon_flux = result_all[self._srcname].spectral().flux(emin, emax) energy_flux = result_all[self._srcname].spectral().eflux(emin, emax) @@ -511,7 +511,7 @@ class cssens(ctools.cscript): for par in model: self._log(str(par)+"\n") elif self._logTerse(): - self._log.parformat("Iteration "+str(iter)) + self._log.parformat("Iteration "+str(iterations)) self._log("TS=") self._log(ts) self._log(" ") @@ -549,7 +549,7 @@ class cssens(ctools.cscript): sensitivity /= num # Compare average flux to last average - if iter > self._num_avg: + if iterations > self._num_avg: if test_crab_flux > 0: ratio = crab_flux/test_crab_flux @@ -570,7 +570,7 @@ class cssens(ctools.cscript): test_crab_flux = crab_flux # Exit loop if number of trials exhausted - if (iter >= self._max_iter): + if (iterations >= self._max_iter): if self._logTerse(): self._log(" Test ended after "+str(self._max_iter)+ " iterations.\n") diff --git a/cscripts/csspec.py b/cscripts/csspec.py index f011e20..dc4d131 100755 --- a/cscripts/csspec.py +++ b/cscripts/csspec.py @@ -60,7 +60,7 @@ class csspec(ctools.cscript): self._debug = False # Initialise observation container from constructor arguments. - self._obs = self._set_input_obs(argv) + self._obs, argv = self._set_input_obs(argv) # Initialise application by calling the appropriate class # constructor. @@ -158,9 +158,10 @@ class csspec(ctools.cscript): # If no layers overlap then throw an exception if (len(use_layers) == 0): - msg = "Energy range ["+str(cube_ebounds.emin(i))+ \ - ", "+str(cube_ebounds.emax(i))+"] of counts "+ \ - "cube does not overlap with specified energy range ["+ \ + msg = "Energy range ["+str(cube_ebounds.emin())+ \ + ", "+str(cube_ebounds.emax())+"] of counts "+ \ + "cube does not overlap with specified energy "+ \ + "range ["+ \ str(emin)+", "+str(emax)+"]. Specify "+ \ "overlapping energy range." raise RuntimeError(msg) diff --git a/cscripts/cstsdist.py b/cscripts/cstsdist.py index 137671b..070b681 100755 --- a/cscripts/cstsdist.py +++ b/cscripts/cstsdist.py @@ -419,17 +419,17 @@ class cstsdist(ctools.cscript): # Write out result immediately if seed == 0: - file = open(self._outfile.url(), 'w') - writer = csv.DictWriter(file, result['colnames']) + f = open(self._outfile.url(), 'w') + writer = csv.DictWriter(f, result['colnames']) headers = {} for n in result['colnames']: headers[n] = n writer.writerow(headers) else: - file = open(self._outfile.url(), 'a') - writer = csv.DictWriter(file, result['colnames']) + f = open(self._outfile.url(), 'a') + writer = csv.DictWriter(f, result['colnames']) writer.writerow(result['values']) - file.close() + f.close() # Return return diff --git a/cscripts/cstsmapmerge.py b/cscripts/cstsmapmerge.py index 4eccec7..a19a0ad 100755 --- a/cscripts/cstsmapmerge.py +++ b/cscripts/cstsmapmerge.py @@ -414,8 +414,8 @@ class cstsmapmerge(ctools.cscript): # Delete TS input maps if requested if self._delete: for filename in self._merged_files: - os.remove(filename) - if self._logTerse(): + os.remove(filename) + if self._logTerse(): self._log(gammalib.parformat("Deleted input file")) self._log(filename) self._log("\n") diff --git a/cscripts/cta_kb_bkgrates.py b/cscripts/cta_kb_bkgrates.py index adebc0a..ca16f81 100755 --- a/cscripts/cta_kb_bkgrates.py +++ b/cscripts/cta_kb_bkgrates.py @@ -4,7 +4,7 @@ # into a column separated file function table. The table is automatically # named "bkg_XXX.txt" where "XXX" is for example "kb_E_50h_v3". # -# Copyright (C) 2011-2014 Juergen Knoedlseder +# Copyright (C) 2011-2016 Juergen Knoedlseder # # 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 @@ -38,8 +38,8 @@ def make_file_function(irfname, filename): filename - File function filename """ # Open performance file and file function - irf = open(irfname, "r") - file = open(filename, "w") + irf = open(irfname, "r") + f = open(filename, "w") # Read data for row in irf: @@ -53,12 +53,12 @@ def make_file_function(irfname, filename): # Extract data fields = row.replace(" ",",").replace(",,",",").rstrip("\n").split(",") logE = float(fields[0]) - aeff = float(fields[1]) - r68 = float(fields[2]) + #aeff = float(fields[1]) + #r68 = float(fields[2]) r80 = float(fields[3]) - eres = float(fields[4]) + #eres = float(fields[4]) bgd = float(fields[5]) - sens = float(fields[6]) + #sens = float(fields[6]) # Compute energy (in MeV) energy = math.pow(10.0, logE)*1.0e6 @@ -76,11 +76,11 @@ def make_file_function(irfname, filename): bkg_rate = 0.0 # Write - file.write(str(energy)+" "+str(bkg_rate)+"\n") + f.write(str(energy)+" "+str(bkg_rate)+"\n") # Close files irf.close() - file.close() + f.close() # Return return @@ -90,9 +90,7 @@ def make_file_function(irfname, filename): # Main entry point # # ================ # if __name__ == '__main__': - """ - Main entry point - """ + # Set parameters caldb = "/usr/local/gamma/share/caldb/data/cta/kb/" diff --git a/cscripts/cta_kb_root2caldb.py b/cscripts/cta_kb_root2caldb.py index 2ae96ab..7e3bdac 100755 --- a/cscripts/cta_kb_root2caldb.py +++ b/cscripts/cta_kb_root2caldb.py @@ -23,7 +23,8 @@ # along with this program. If not, see . # # ========================================================================== -from ROOT import TFile, TH1F, TH2F +#from ROOT import TFile, TH1F, TH2F +from ROOT import TFile, TH1F from datetime import datetime import gammalib import math @@ -506,11 +507,11 @@ if __name__ == '__main__': """ """ # Get ROOT filename from command line argument - if (len(sys.argv) > 1): - filename = sys.argv[1] - else: - db = "/Users/jurgen/Documents/Travail/Projects/CTA/WP-MC/root/Kb_all_20deg_root" - #db = "/Users/jurgen/Documents/Travail/Projects/CTA/WP-MC/root/Kb_all_50deg_root" + #if (len(sys.argv) > 1): + # filename = sys.argv[1] + #else: + db = "/Users/jurgen/Documents/Travail/Projects/CTA/WP-MC/root/Kb_all_20deg_root" + #db = "/Users/jurgen/Documents/Travail/Projects/CTA/WP-MC/root/Kb_all_50deg_root" # Get list of all performance files files = glob.glob(db+"/*.root") @@ -520,26 +521,26 @@ if __name__ == '__main__': irf = caldb("kb") # Loop over all files - for file in files: + for f in files: # Get filename - head, tail = os.path.split(file) + head, tail = os.path.split(f) # Build name name = tail.strip(".root") if name.find("kb_") != -1: name = name.replace("_20deg", "") elif name.find("IFAE_") != -1: - name = name.replace("Subarray", "") - name = name.replace("IFAE_", "") - name = "IFAE_"+name - i = name.find("hours") - if i != -1: - name = name[0:i+1] + name = name.replace("Subarray", "") + name = name.replace("IFAE_", "") + name = "IFAE_"+name + index = name.find("hours") + if index != -1: + name = name[0:index+1] print(tail) # Open calibration file irf.add(name) # Translate ROOT to CALDB information - irf.root2caldb(file) + irf.root2caldb(f) diff --git a/cscripts/cta_make_nodes.py b/cscripts/cta_make_nodes.py index 425fa0b..9d8a3dd 100755 --- a/cscripts/cta_make_nodes.py +++ b/cscripts/cta_make_nodes.py @@ -20,7 +20,6 @@ # along with this program. If not, see . # # ========================================================================== -#import gammalib import sys import math @@ -28,7 +27,7 @@ import math # ============== # # Write one node # # ============== # -def write_node(file, energy, scale): +def write_node(f, energy, scale): """ Writes one node to XML file. """ @@ -37,22 +36,22 @@ def write_node(file, energy, scale): str_scale = str(scale) # Write start tag - file.write(' \n') + f.write(' \n') # Write energy - file.write(' \n') + f.write(' \n') # Write intensity - file.write(' \n') + f.write(' \n') # Write end tag - file.write(' \n') + f.write(' \n') # Return return @@ -66,7 +65,7 @@ def create_nodes(emin, emax, enumbins): Create nodes (energies in TeV). """ # Open file - file = open("nodes.xml", "w") + f = open("nodes.xml", "w") # Set node boundaries elogmin = math.log10(float(emin)) @@ -83,14 +82,14 @@ def create_nodes(emin, emax, enumbins): scale = 3.45e-17 * math.pow(energy, -2.63) # Write node - write_node(file, energy, scale) + write_node(f, energy, scale) # Debug #sys.stdout.write(energy) #sys.stdout.write(math.pow(10.0, int(math.log10(scale))-1.0)±"\n") # Close file - file.close() + f.close() # Return return @@ -100,9 +99,7 @@ def create_nodes(emin, emax, enumbins): # Main script entry point # # ======================= # if __name__ == '__main__': - """ - Compute nodes. - """ + # Check command line usage = "Usage: cta_make_nodes emin emax enumbins" if len(sys.argv) < 3: diff --git a/cscripts/cta_root2caldb.py b/cscripts/cta_root2caldb.py index aa2a6d7..4b88cf4 100755 --- a/cscripts/cta_root2caldb.py +++ b/cscripts/cta_root2caldb.py @@ -3,7 +3,7 @@ # This script generates IRFs in the CALDB format of HEASARC using ROOT # offaxis performance files. # ------------------------------------------------------------------------- -# Copyright (C) 2011-2015 Juergen Knoedlseder +# Copyright (C) 2011-2016 Juergen Knoedlseder # ------------------------------------------------------------------------- # # This program is free software: you can redistribute it and/or modify @@ -24,12 +24,13 @@ # - gammalib # - pyroot # ========================================================================== -from ROOT import TFile, TH1F, TH2F, TH3F +#from ROOT import TFile, TH1F, TH2F, TH3F +from ROOT import TFile, TH2F, TH3F import gammalib from datetime import datetime import math -import sys -import glob +#import sys +#import glob import os diff --git a/pyext/ctool.i b/pyext/ctool.i index 89234ec..d16b956 100644 --- a/pyext/ctool.i +++ b/pyext/ctool.i @@ -87,6 +87,7 @@ public: %rename(_read_ahead) read_ahead() const; %rename(_time_reference) time_reference() const; %rename(_get_observations) get_observations(const bool& get_response = true); + %rename(_setup_observations) setup_observations(GObservations& obs); %rename(_create_ebounds) create_ebounds(); %rename(_create_map) create_map(const GObservations& obs); %rename(_create_cube) create_cube(const GObservations& obs); @@ -105,6 +106,7 @@ public: const bool& read_ahead(void) const; const GTimeReference& time_reference(void) const; GObservations get_observations(const bool& get_response = true); + void setup_observations(GObservations& obs); // Protected methods that create objects from user parameters GEbounds create_ebounds(void); @@ -170,14 +172,15 @@ public: # Initialise observation container from constructor arguments. In case that # an observation container is provided as argument, this container will be -# used to initialise the class. +# used to initialise the class. The function returns a tuple containing the +# observation container and the eventually reduced argument list. def _set_input_obs(self, argv): if len(argv) > 0 and isinstance(argv[0],gammalib.GObservations): obs = argv[0] argv = argv[1:] else: obs = gammalib.GObservations() - return obs + return (obs, argv) cscript._set_input_obs = _set_input_obs # Initialise application by calling the appropriate class constructor. diff --git a/pyext/setup.py.in b/pyext/setup.py.in index 1d45b41..30ddad3 100644 --- a/pyext/setup.py.in +++ b/pyext/setup.py.in @@ -1,7 +1,7 @@ ############################################################################# # Disutils setup file for ctools Python extension module # # ------------------------------------------------------------------------- # -# Copyright (C) 2011-2014 Juergen Knoedlseder # +# Copyright (C) 2011-2016 Juergen Knoedlseder # # ------------------------------------------------------------------------- # # # # This program is free software: you can redistribute it and/or modify # @@ -27,35 +27,35 @@ import sys # Detect execution mode is_build = True if (len(sys.argv) >= 2): - if (sys.argv[1] == 'install'): - is_build = False + if (sys.argv[1] == 'install'): + is_build = False if is_build: - sys.stdout.write("Python extension module building:\n") + sys.stdout.write("Python extension module building:\n") else: - sys.stdout.write("Python extension module installation:\n") + sys.stdout.write("Python extension module installation:\n") # Darwin MacPorts kluge: Remove -L/opt/local/lib path from LDSHARED as this # overrides the GammaLib settings. Otherwise, MacPorts libraries would be # used instead of the libraries selected by GammaLib if sys.platform.lower().startswith("darwin"): - - # Get configuration variables - cvars = sysconfig.get_config_vars() - - # Extract LDSHARED - ldshared_old = cvars.get('LDSHARED') - ldshared_new = '' - ldshared_rm = '' - for option in ldshared_old.split(' '): - if option.lstrip().find('-L/opt/local/lib') != 0: - ldshared_new += option + ' ' - else: - ldshared_rm += option + ' ' - if (len(ldshared_rm) > 0): - sys.stdout.write("Removed following option(s) from LDSHARED: "+ldshared_rm+"\n"); - - # Set LDSHARED - cvars['LDSHARED'] = ldshared_new + + # Get configuration variables + cvars = sysconfig.get_config_vars() + + # Extract LDSHARED + ldshared_old = cvars.get('LDSHARED') + ldshared_new = '' + ldshared_rm = '' + for option in ldshared_old.split(' '): + if option.lstrip().find('-L/opt/local/lib') != 0: + ldshared_new += option + ' ' + else: + ldshared_rm += option + ' ' + if (len(ldshared_rm) > 0): + sys.stdout.write("Removed following option(s) from LDSHARED: "+ldshared_rm+"\n"); + + # Set LDSHARED + cvars['LDSHARED'] = ldshared_new # Set default include directories and data files include_dirs = ['@prefix@/include/gammalib'] @@ -66,8 +66,8 @@ libraries = ['ctools','gamma'] # Add install directory to runtime library path if os.path.exists('@prefix@/lib'): - library_dirs.append('@prefix@/lib') - runtime_library_dirs.append('@prefix@/lib') + library_dirs.append('@prefix@/lib') + runtime_library_dirs.append('@prefix@/lib') # Add application include files. The method we use here is not very # intelligent as it simple will add all names that start with 'ct' @@ -75,44 +75,49 @@ if os.path.exists('@prefix@/lib'): # add an explicit filter to avoid appending those. dirs = glob.glob('../src/ct*') for dir in dirs: - if len(dir) > 0 and 'ctools' not in dir: - path = os.path.abspath(dir) - include_dirs.append(path) + if len(dir) > 0 and 'ctools' not in dir: + path = os.path.abspath(dir) + include_dirs.append(path) include_dirs.append(os.path.abspath('../src/support')) # Add extra libraries extra_libs = ('@PYTHON_EXTRA_LIBS@').split(',') for extra_lib in extra_libs: - if len(extra_lib) > 0: - libraries.append(extra_lib) + if len(extra_lib) > 0: + libraries.append(extra_lib) # Get extra library paths extra_lib_paths = ('@LDFLAGS@').split(' ') for extra_lib_path in extra_lib_paths: - if extra_lib_path.lstrip().find('-L') == 0: - dir = extra_lib_path.lstrip().lstrip('-L') - if len(dir) > 0: - library_dirs.append(dir) - runtime_library_dirs.append(dir) + if extra_lib_path.lstrip().find('-L') == 0: + dir = extra_lib_path.lstrip().lstrip('-L') + if len(dir) > 0: + library_dirs.append(dir) + runtime_library_dirs.append(dir) # Get extra include directories extra_include_dirs = ('@CPPFLAGS@').split(' ') for extra_include_dir in extra_include_dirs: - if extra_include_dir.lstrip().find('-I') == 0: - dir = extra_include_dir.lstrip().lstrip('-I') - if len(dir) > 0: - include_dirs.append(dir) + if extra_include_dir.lstrip().find('-I') == 0: + dir = extra_include_dir.lstrip().lstrip('-I') + if len(dir) > 0: + include_dirs.append(dir) # Set data files to install. data_files=[('@prefix@/share/ctools/swig', install_files)] # Optionally enable OpenMP if '@OPENMP_TRUE@' != '#': - extra_compile_args = ['@OPENMP_CXXFLAGS@'] - extra_link_args = ['@OPENMP_CXXFLAGS@'] + extra_compile_args = ['@OPENMP_CXXFLAGS@'] + extra_link_args = ['@OPENMP_CXXFLAGS@'] else: - extra_compile_args = [] - extra_link_args = [] + extra_compile_args = [] + extra_link_args = [] + +# Append coverage CXXFLAGS +if len('@COVERAGE_CXXFLAGS@') > 0: + extra_compile_args.append('@COVERAGE_CXXFLAGS@') + extra_link_args.append('@COVERAGE_CXXFLAGS@') # Optionally add link flag for El Capitan post-processing of modules if '@IS_ELCAPITAN_TRUE@' != '#': @@ -133,25 +138,25 @@ setup(name='@PACKAGE_NAME@', runtime_library_dirs=runtime_library_dirs, extra_compile_args=extra_compile_args, extra_link_args=extra_link_args)], - data_files=data_files, - license="GNU GPL v3", - long_description="CTA science analysis tools", - platforms=["Linux", "Solaris", "Mac OS-X", "Unix", "BSD"], + data_files=data_files, + license="GNU GPL v3", + long_description="CTA science analysis tools", + platforms=["Linux", "Solaris", "Mac OS-X", "Unix", "BSD"], ) # Post-processing. A symbolic link is created to the directory containing the # shared wrapper libraries so that we can add this directory to the PYTHONPATH # for testing. if is_build: - try: - os.chdir('build') - os.remove('ctools') - sys.stdout.write("Removed symbolic link build/ctools.\n") - except: - pass - try: - path = glob.glob('lib.*') - os.symlink(path[0], 'ctools') - sys.stdout.write("Created symbolic link build/ctools.\n") - except: - sys.stdout.write("WARNING: Failed to create symbolic link build/ctools.\n") + try: + os.chdir('build') + os.remove('ctools') + sys.stdout.write("Removed symbolic link build/ctools.\n") + except: + pass + try: + path = glob.glob('lib.*') + os.symlink(path[0], 'ctools') + sys.stdout.write("Created symbolic link build/ctools.\n") + except: + sys.stdout.write("WARNING: Failed to create symbolic link build/ctools.\n") diff --git a/src/ctbin/ctbin.cpp b/src/ctbin/ctbin.cpp index 4da9f0e..cf03319 100644 --- a/src/ctbin/ctbin.cpp +++ b/src/ctbin/ctbin.cpp @@ -62,7 +62,7 @@ ctbin::ctbin(void) : ctool(CTBIN_NAME, CTBIN_VERSION) // Return return; -} +} // GCOV_EXCL_LINE /***********************************************************************//** diff --git a/src/support/ctool.cpp b/src/support/ctool.cpp index dd9608a..ec19f93 100644 --- a/src/support/ctool.cpp +++ b/src/support/ctool.cpp @@ -49,6 +49,7 @@ /* __ Method name definitions ____________________________________________ */ +#define G_SETUP_OBSERVATION "ctool::setup_observations(GObservations&)" #define G_GET_MEAN_POINTING "ctool::get_mean_pointing(GObservations&)" #define G_CREATE_EBOUNDS "ctool::create_ebounds()" #define G_PROVIDE_HELP "ctool::provide_help()" @@ -378,7 +379,98 @@ GObservations ctool::get_observations(const bool& get_response) // Return observation container return obs; +} + + +/***********************************************************************//** + * @brief Setup observation container + * + * @param[in] obs Observation container + * + * @exception GException::invalid_value + * Invalid "inobs" parameter encountered. + * + * Setup an observation container by appending all required information that + * is missing using user parameters. + * + * If the observation container @p obs is empty, the method will extract the + * filename from the "inobs" parameter, and if the filename is neither empty + * nor "NONE", will load the observation container from the specified file. + * If the file is a FITS file, and dependent on the information that is found + * in the FITS file, the method will append a single CTA observation to the + * observation container, containing either an event list or a counts cube. + * Otherwise, the method assumes that the file is an observation definition + * XML file and will load that file into the observation container. + * + * If the filename is empty or "NONE", the method will throw an exception. + * + * The method will then setup the response for all CTA observations in the + * container (see the set_response() method for more information). + * + * The method will furthermore setup the observation boundaries for all CTA + * observations in the container that contain event lists (see the + * set_obs_bounds() method for more information). + ***************************************************************************/ +void ctool::setup_observations(GObservations& obs) +{ + // Load observations if there are none in the container + if (obs.size() == 0) { + + // Get the filename from the "inobs" parameter + GFilename filename((*this)["inobs"].filename()); + + // Load observation container from file if the filename is neither + // "NONE" nor empty ... + if ((gammalib::toupper(filename.url()) != "NONE") && + !filename.is_empty()) { + + // If file does not exist then throw an exceptions + if (!filename.exists()) { + std::string msg = "The \"inobs\" parameter \""+filename+ + "\" specifies a file that does not " + "exist. Please specify a valid event " + "list, counts cube or observation " + "definition XML file."; + throw GException::invalid_value(G_SETUP_OBSERVATION, msg); + } + + // If the file is a FITS file then load it into a CTA + // observation and append that observation to the container ... + else if (filename.is_fits()) { + GCTAObservation cta(filename); + obs.append(cta); + } + + // ... otherwise it is assumed that the file is an observation + // definition XML file and the file is loaded into the + // observation container + else { + obs.load(filename); + } + + } + + // ... otherwise, if the filename is "NONE" or empty then throw + // an exception. + else { + std::string msg = "The \"inobs\" parameter \""+filename+ + "\" is not a valid filename. Please specify " + "a valid event list, counts cube or " + "observation definition XML file."; + throw GException::invalid_value(G_SETUP_OBSERVATION, msg); + } + + } // endif: there were no observations in the container + + // Setup the response for all CTA observations in the observation + // container + set_response(obs); + // Setup boundaries for all CTA observations containing event lists + set_obs_bounds(obs); + + // Return + return; } @@ -595,15 +687,16 @@ GCTAEventCube ctool::create_cube(const GObservations& obs) * method furthermore sets the ontime, livetime and deadtime correction * factor. * - * The following parameters are read: - * ra - Right Ascension of pointing and RoI centre (deg) - * dec - Declination of pointing and RoI centre (deg) - * rad - Radius of RoI (deg) - * deadc - Deadtime correction factor - * tmin - Start time - * tmax - Stop time - * emin - Minimum energy (TeV) - * emax - Maximum energy (TeV) + * The following parameters are read + * + * ra: Right Ascension of pointing and RoI centre (deg) + * dec: Declination of pointing and RoI centre (deg) + * rad: Radius of RoI (deg) + * deadc: Deadtime correction factor + * tmin: Start time + * tmax: Stop time + * emin: Minimum energy (TeV) + * emax: Maximum energy (TeV) ***************************************************************************/ GCTAObservation ctool::create_cta_obs(void) { @@ -804,12 +897,9 @@ void ctool::require_inobs_nocube(const std::string& method) * * @param[in,out] obs Observation container * - * Set the response for a CTA observations in the container that so far have - * no response using the "caldb" and "irf" task parameters. - * - * The following parameters are read: - * caldb - Calibration database - * irf - Instrument response function + * Set the response for all CTA observations in the container that so far + * have no response information available. See the set_obs_response() + * method for more information about the user parameters that are read. ***************************************************************************/ void ctool::set_response(GObservations& obs) { @@ -841,69 +931,83 @@ void ctool::set_response(GObservations& obs) * * @param[in,out] obs CTA observation * - * Set the response for a CTA observation. If the CTA observation contains + * Set the response for one CTA observation. If the CTA observation contains * a counts cube the method first attempts to set the response information - * using the @a expcube and @a psfcube user parameters. If this is not - * successful (for example because the parameters do not exist or their - * value is NONE or blank), the method will set the response information - * from the @a caldb and @a irf user parameters. + * using the following user parameters + * + * expcube: Exposure cube file + * psfcube: Point spread function cube file + * edispcube: Energy dispersion cube file (optional) + * bkgcube: Background cube file + * + * If none of the @p expcube, @p psfcube and @p bkgcube parameters is + * "NONE" or empty, the method will allocate a stacked response for the + * observation. In case that an energy dispersion cube is provided using + * the @p edispcube parameter, this energy dispersion cube is also loaded + * and attached to the stacked response function. + * + * In case that any of the @p expcube, @p psfcube and @p bkgcube parameters + * is "NONE" or empty, the method will assume that a IRF response should + * be used, and will attempty to set the response information using the + * following user parameters + * + * caldb: Calibration database + * irf: Instrument response function * - * The following parameters are read: - * expcube - Exposure cube file (optional) - * psfcube - PSF cube file (optional) - * caldb - Calibration database (optional) - * irf - Instrument response function (optional) ***************************************************************************/ void ctool::set_obs_response(GCTAObservation* obs) { // Initialise response flag bool has_response = false; - // If the observation contains a counts cube, then first check whether - // the expcube and psfcube parameters exist and are not NONE + // If the observation contains a counts cube, then check whether + // response information for a stacked analysis is provided, and if so, + // set that response information for the CTA observation if (dynamic_cast(obs->events()) != NULL) { if (has_par("expcube") && has_par("psfcube") && has_par("edispcube") && has_par("bkgcube")) { // Get filenames - std::string expcube = (*this)["expcube"].filename(); - std::string psfcube = (*this)["psfcube"].filename(); - std::string edispcube = (*this)["edispcube"].filename(); - std::string bkgcube = (*this)["bkgcube"].filename(); - - // Extract response information if available - if ((expcube != "NONE") && (psfcube != "NONE") && - (bkgcube != "NONE") && - (gammalib::strip_whitespace(expcube) != "") && - (gammalib::strip_whitespace(psfcube) != "") && - (gammalib::strip_whitespace(bkgcube) != "")) { + GFilename expcube((*this)["expcube"].filename()); + GFilename psfcube((*this)["psfcube"].filename()); + GFilename edispcube((*this)["edispcube"].filename()); + GFilename bkgcube((*this)["bkgcube"].filename()); + + // Extract stacked response information if available + if ((gammalib::toupper(expcube.url()) != "NONE") && + (gammalib::toupper(psfcube.url()) != "NONE") && + (gammalib::toupper(bkgcube.url()) != "NONE") && + !expcube.is_empty() && + !psfcube.is_empty() && + !bkgcube.is_empty()) { // Get exposure, PSF and background cubes GCTACubeExposure exposure(expcube); GCTACubePsf psf(psfcube); GCTACubeBackground background(bkgcube); - // Optionally load energy dispersion cube - if ((edispcube != "NONE") && - (gammalib::strip_whitespace(edispcube) != "")) { + // If an energy dispersion cube is available, get also + // this cube, and set all information as stacked response + // for this CTA observation + if ((gammalib::toupper(edispcube.url()) != "NONE") && + !edispcube.is_empty()) { GCTACubeEdisp edisp(edispcube); obs->response(exposure, psf, edisp, background); } else { - // Set reponse obs->response(exposure, psf, background); } - // Signal response availability + // Signal that response is available has_response = true; } // endif: filenames were available - } // endif: expcube and psfcube parameters exist + } // endif: expcube, psfcube and bkdcube parameters valid } // endif: observation contains a counts cube - // If we have not yet response information then get it now - // from the caldb and irf parameters + // If we have not yet response information then get it now from the + // caldb and irf user parameters if (!has_response) { // Load response information @@ -938,59 +1042,70 @@ void ctool::set_obs_response(GCTAObservation* obs) /***********************************************************************//** - * @brief Set observation boundaries + * @brief Set observation boundaries for all CTA observations * * @param[in,out] obs Observation container + * + * Sets the observation boundaries for all CTA observations in the container + * that contain event lists. If the energy boundaries of an event list are + * empty, the method gets the + * + * emin: Minimum energy (TeV) + * emax: Maximum energy (TeV) + * + * user parameters and sets a single energy boundary element @p [emin,emax]. + * + * If the RoI radius of the event list is invalid, the method gets the + * + * rad: Region of Interest radius (deg) + * + * user parameter to set the radius of the RoI, and uses the CTA pointing + * direction to set the centre of the RoI. ***************************************************************************/ void ctool::set_obs_bounds(GObservations& obs) { - // Setup response for all observations + // Loop over all observations in the container for (int i = 0; i < obs.size(); ++i) { - // Is this observation a CTA observation? + // Get pointer to CTA observation GCTAObservation* cta = dynamic_cast(obs[i]); - // Continue only if observation is CTA and has events - if ((cta != NULL) && (cta->has_events())) { - - // Get pointer on event list - GCTAEventList* list = const_cast(dynamic_cast(cta->events())); - - // Continue only if it's valid - if (list != NULL) { - - // If there are no energy boundaries then read the user - // parameters and add them - if (list->ebounds().is_empty()) { - - // If there are "emin" and "emax" parameters then use - // them - if (has_par("emin") && has_par("emax")) { - double emin = (*this)["emin"].real(); - double emax = (*this)["emax"].real(); - GEbounds ebounds(GEnergy(emin, "TeV"), - GEnergy(emax, "TeV")); - list->ebounds(ebounds); - } - - } + // Fall through if observation is not a CTA observation or if it + // does not contain events + if ((cta == NULL) || (!cta->has_events())) { + continue; + } - // If there is no ROI then read the user parameters and add - // them - if (list->roi().radius() == 0) { - - // If there is a "rad" parameter then use it - if (has_par("rad")) { - double rad = (*this)["rad"].real(); - GCTARoi roi(GCTAInstDir(cta->pointing().dir()), rad); - list->roi(roi); - } + // Get pointer on CTA event list + GCTAEventList* list = const_cast + (dynamic_cast(cta->events())); - } + // Fall through if CTA observation does not contain an event list + if (list == NULL) { + continue; + } - } // endif: list was valid + // If there are no energy boundaries then read the "emin" and "emax" + // user parameters and add them + if (list->ebounds().is_empty()) { + if (has_par("emin") && has_par("emax")) { + double emin((*this)["emin"].real()); + double emax((*this)["emax"].real()); + GEbounds ebounds(GEnergy(emin, "TeV"), + GEnergy(emax, "TeV")); + list->ebounds(ebounds); + } + } - } // endif: observation was CTA + // If there is no RoI then read the "rad" user parameters and use + // the pointing direction to set the RoI + if (list->roi().radius() == 0) { + if (has_par("rad")) { + double rad = (*this)["rad"].real(); + GCTARoi roi(GCTAInstDir(cta->pointing().dir()), rad); + list->roi(roi); + } + } } // endfor: looped over observations diff --git a/src/support/ctool.hpp b/src/support/ctool.hpp index e26a008..dfb92ee 100644 --- a/src/support/ctool.hpp +++ b/src/support/ctool.hpp @@ -72,6 +72,7 @@ protected: const bool& read_ahead(void) const; const GTimeReference& time_reference(void) const; GObservations get_observations(const bool& get_response = true); + void setup_observations(GObservations& obs); // Protected methods that create objects from user parameters GEbounds create_ebounds(void); diff --git a/test/test_cscaldb.py b/test/test_cscaldb.py index c9de3d9..b2d8984 100755 --- a/test/test_cscaldb.py +++ b/test/test_cscaldb.py @@ -18,6 +18,7 @@ # along with this program. If not, see . # # ========================================================================== +import os import gammalib import cscripts @@ -28,7 +29,11 @@ import cscripts class Test(gammalib.GPythonTestSuite): """ Test class for cscaldb script. + + This test class makes unit tests for the cslightcrv script by using it + from the command line and from Python. """ + # Constructor def __init__(self): """ @@ -49,15 +54,15 @@ class Test(gammalib.GPythonTestSuite): self.name("cscaldb") # Append tests - self.append(self._test_functional, "Test cscaldb functionality") + self.append(self._test_python, "Test cscaldb from Python") # Return return - # Test cscaldb functionality - def _test_functional(self): + # Test cscaldb from Python + def _test_python(self): """ - Test cscaldb functionality. + Test cscaldb from Python. """ # Set-up cscaldb caldb = cscripts.cscaldb() diff --git a/test/test_cslightcrv.py b/test/test_cslightcrv.py new file mode 100755 index 0000000..e66bd55 --- /dev/null +++ b/test/test_cslightcrv.py @@ -0,0 +1,271 @@ +#! /usr/bin/env python +# ========================================================================== +# This scripts performs unit tests for the cslightcrv script. +# +# Copyright (C) 2016 Juergen Knoedlseder +# +# 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 . +# +# ========================================================================== +import os +import gammalib +import cscripts + + +# ================================ # +# Test class for cslightcrv script # +# ================================ # +class Test(gammalib.GPythonTestSuite): + """ + Test class for cslightcrv script. + + This test class makes unit tests for the cslightcrv script by using it + from the command line and from Python. + """ + + # Constructor + def __init__(self): + """ + Constructor. + """ + # Call base class constructor + gammalib.GPythonTestSuite.__init__(self) + + # Set members + self._events_name = "data/crab_events.fits" + self._model_name = "data/crab.xml" + + # Return + return + + # Set test functions + def set(self): + """ + Set all test functions. + """ + # Set test name + self.name("cslightcrv") + + # Append tests + self.append(self._test_cmd, "Test cslightcrv on command line") + self.append(self._test_python, "Test cslightcrv from Python") + + # Return + return + + # Test cslightcrv on command line + def _test_cmd(self): + """ + Test cslightcrv on the command line. + """ + # Kluge to set the command (installed version has no README file) + if os.path.isfile("README"): + cslightcrv = "../cscripts/cslightcrv.py" + else: + cslightcrv = "cslightcrv" + + # Setup cslightcrv command + cmd = cslightcrv+' inobs="'+self._events_name+'"'+ \ + ' inmodel="'+self._model_name+'"'+ \ + ' srcname="Crab" caldb="prod2" irf="South_0.5h"'+ \ + ' tbinalg="LIN" tmin=51544.50 tmax=51544.53'+ \ + ' tbins=3 enumbins=0 emin=0.1 emax=100.0'+ \ + ' outfile="lightcurve_cmd1.fits"' + + # Execute cslightcrv, make sure we catch any exception + try: + rc = os.system(cmd+" >/dev/null 2>&1") + except: + pass + + # Check if execution was successful + self.test_assert(rc == 0, + "Successful cslightcrv execution on command line") + + # Check light curve + self._check_light_curve("lightcurve_cmd1.fits", 3) + + # Setup cslightcrv command + cmd = cslightcrv+' inobs="events_that_do_not_exist.fits"'+ \ + ' inmodel="'+self._model_name+'"'+ \ + ' srcname="Crab" caldb="prod2" irf="South_0.5h"'+ \ + ' tbinalg="LIN" tmin=51544.50 tmax=51544.53'+ \ + ' tbins=3 enumbins=0 emin=0.1 emax=100.0'+ \ + ' outfile="lightcurve_cmd1.fits"' + + # Execute cslightcrv, make sure we catch any exception + try: + rc = os.system(cmd+" >/dev/null 2>&1") + except: + pass + + # Check if execution failed + self.test_assert(rc != 0, + "Failure of cslightcrv execution on command line") + + # Return + return + + # Test cslightcrv from Python + def _test_python(self): + """ + Test cslightcrv from Python. + """ + # Set-up unbinned cslightcrv + lcrv = cscripts.cslightcrv() + lcrv["inobs"] = self._events_name + lcrv["inmodel"] = self._model_name + lcrv["srcname"] = "Crab" + lcrv["caldb"] = "prod2" + lcrv["irf"] = "South_0.5h" + lcrv["tbinalg"] = "LIN" + lcrv["tmin"] = 51544.50 + lcrv["tmax"] = 51544.53 + lcrv["tbins"] = 3 + lcrv["enumbins"] = 0 + lcrv["emin"] = 0.1 + lcrv["emax"] = 100.0 + lcrv["outfile"] = "lightcurve_py1.fits" + + # Run cslightcrv script and save light curve + lcrv.run() + lcrv.save() + + # Check light curve + self._check_light_curve("lightcurve_py1.fits", 3) + + # Now use FILE as time bin algorithm. For this we need first to + # create an ASCII file. We use now 6 time bins. The ASCII file + # is saved into the file "lightcurve_py2.dat". + csv = gammalib.GCsv(2,2) + tmin = 51544.50 + tdelta = 0.01 + for i in range(csv.nrows()): + csv[i,0] = "%.5f" % (tmin + i * tdelta) + csv[i,1] = "%.5f" % (tmin + (i+1) * tdelta) + csv.save("lightcurve_py2.dat", " ", True) + + # Set-up unbinned cslightcrv + lcrv = cscripts.cslightcrv() + lcrv["inobs"] = self._events_name + lcrv["inmodel"] = self._model_name + lcrv["srcname"] = "Crab" + lcrv["caldb"] = "prod2" + lcrv["irf"] = "South_0.5h" + lcrv["tbinalg"] = "FILE" + lcrv["tbinfile"] = "lightcurve_py2.dat" + lcrv["enumbins"] = 0 + lcrv["emin"] = 0.1 + lcrv["emax"] = 100.0 + lcrv["outfile"] = "lightcurve_py2.fits" + + # Execute cslightcrv script + lcrv.execute() + + # Check light curve + self._check_light_curve("lightcurve_py2.fits", 2) + + # Now we setup an observation container on input. We attached the + # model to the observation container so that cslightcrv should + # no longer query for the parameter. + cta = gammalib.GCTAObservation(self._events_name) + obs = gammalib.GObservations() + obs.append(cta) + obs.models(self._model_name) + + # Set-up unbinned cslightcrv from observation container. Now use + # the GTI algorithm so that we test all timing algorithms. + lcrv = cscripts.cslightcrv(obs) + lcrv["srcname"] = "Crab" + lcrv["caldb"] = "prod2" + lcrv["irf"] = "South_0.5h" + lcrv["tbinalg"] = "GTI" + lcrv["enumbins"] = 0 + lcrv["emin"] = 0.1 + lcrv["emax"] = 100.0 + lcrv["outfile"] = "lightcurve_py3.fits" + + # Execute cslightcrv script + lcrv.execute() + + # Check light curve + self._check_light_curve("lightcurve_py3.fits", 1) + + # Finally we set-up a binned cslightcrv + lcrv = cscripts.cslightcrv() + lcrv["inobs"] = self._events_name + lcrv["inmodel"] = self._model_name + lcrv["srcname"] = "Crab" + lcrv["caldb"] = "prod2" + lcrv["irf"] = "South_0.5h" + lcrv["tbinalg"] = "LIN" + lcrv["tmin"] = 51544.50 + lcrv["tmax"] = 51544.53 + lcrv["tbins"] = 2 + lcrv["emin"] = 0.1 + lcrv["emax"] = 100.0 + lcrv["enumbins"] = 10 + lcrv["coordsys"] = "CEL" + lcrv["proj"] = "TAN" + lcrv["xref"] = 83.63 + lcrv["yref"] = 22.01 + lcrv["nxpix"] = 20 + lcrv["nypix"] = 20 + lcrv["binsz"] = 0.02 + lcrv["outfile"] = "lightcurve_py4.fits" + + # Execute cslightcrv script + lcrv.execute() + + # Check light curve + self._check_light_curve("lightcurve_py4.fits", 2) + + # Return + return + + # Check light curve result file + def _check_light_curve(self, filename, bins): + """ + Check light curve file. + """ + # Expected column names + cols = ["MJD", "e_MJD", "Prefactor", "e_Prefactor", + "Index", "e_Index", "TS", "UpperLimit"] + + # Open FITS file + fits = gammalib.GFits(filename) + + # Check FITS file structure + self.test_value(fits.size(), 2, + 'Check for 2 extensions in light curve FITS file') + self.test_assert(fits.contains("LIGHTCURVE"), + 'FITS file contains "LIGHTCURVE" extension') + + # Get LIGHTCURVE table + table = fits["LIGHTCURVE"] + + # Check FITS table structure + self.test_value(table.ncols(), len(cols), + 'Check for %d columns in light curve FITS table' % len(cols)) + self.test_value(table.nrows(), bins, + 'Check for %d rows in light curve FITS table' % bins) + for col in cols: + self.test_assert(table.contains(col), + 'FITS file contains "'+col+'" column') + + # Close FITS file + fits.close() + + # Return + return diff --git a/test/test_ctbin.py b/test/test_ctbin.py index 8884f0a..eab2c70 100755 --- a/test/test_ctbin.py +++ b/test/test_ctbin.py @@ -1,6 +1,6 @@ #! /usr/bin/env python # ========================================================================== -# This scripts performs unit tests for the ctbin module. +# This scripts performs unit tests for the ctbin tool. # # Copyright (C) 2014-2016 Juergen Knoedlseder # @@ -18,16 +18,20 @@ # along with this program. If not, see . # # ========================================================================== +import os import gammalib import ctools -# =========================== # -# Test class for ctbin module # -# =========================== # +# ========================= # +# Test class for ctbin tool # +# ========================= # class Test(gammalib.GPythonTestSuite): """ - Test class for ctbin module. + Test class for ctbin tool. + + This test class makes unit tests for the ctbin tool by using it from + the command line and from Python. """ # Constructor @@ -53,15 +57,68 @@ class Test(gammalib.GPythonTestSuite): self.name("ctbin") # Append tests - self.append(self._test_ctbin, "Test ctbin") + self.append(self._test_cmd, "Test ctbin on command line") + self.append(self._test_python, "Test ctbin from Python") # Return return - # Test ctbin - def _test_ctbin(self): + # Test ctbin on command line + def _test_cmd(self): """ - Test ctbin. + Test ctbin on the command line. + """ + # Kluge to set the command (installed version has no README file) + if os.path.isfile("README"): + ctbin = "../src/ctbin/ctbin" + else: + ctbin = "ctbin" + + # Setup ctbin command + cmd = ctbin+' inobs="'+self._events_name+'"'+ \ + ' outcube="cntmap_cmd1.fits"'+\ + ' emin=0.1 emax=100.0 enumbins=20 ebinalg="LOG"'+ \ + ' nxpix=200 nypix=200 binsz=0.02 coordsys="CEL"'+ \ + ' xref=83.63 yref=22.01 proj="CAR"' + + # Execute ctbin, make sure we catch any exception + try: + rc = os.system(cmd+" >/dev/null 2>&1") + except: + pass + + # Check if execution was successful + self.test_assert(rc == 0, + "Successful ctbin execution on command line") + + # Load counts cube and check content. + evt = gammalib.GCTAEventCube("cntmap_cmd1.fits") + self._check_cube(evt, 5542) + + # Setup ctbin command + cmd = ctbin+' inobs="events_that_do_not_exist.fits"'+ \ + ' outcube="cntmap_cmd2.fits"'+\ + ' emin=0.1 emax=100.0 enumbins=20 ebinalg="LOG"'+ \ + ' nxpix=200 nypix=200 binsz=0.02 coordsys="CEL"'+ \ + ' xref=83.63 yref=22.01 proj="CAR"' + + # Execute ctbin, make sure we catch any exception + try: + rc = os.system(cmd+" >/dev/null 2>&1") + except: + pass + + # Check if execution failed + self.test_assert(rc != 0, + "Failure of ctbin execution on command line") + + # Return + return + + # Test ctbin from Python + def _test_python(self): + """ + Test ctbin from Python. """ # Set-up ctbin bin = ctools.ctbin() @@ -83,15 +140,15 @@ class Test(gammalib.GPythonTestSuite): bin.run() # Check content of observation and cube - self._test_observation(bin, 5542) - self._test_cube(bin.cube(), 5542) + self._check_observation(bin, 5542) + self._check_cube(bin.cube(), 5542) # Test copy constructor cpy_bin = bin.copy() # Check content of observation and cube - self._test_observation(cpy_bin, 5542) - self._test_cube(cpy_bin.cube(), 5542) + self._check_observation(cpy_bin, 5542) + self._check_cube(cpy_bin.cube(), 5542) # Run copy of ctbin tool again cpy_bin.run() @@ -100,8 +157,8 @@ class Test(gammalib.GPythonTestSuite): # event cube as on input the observation is binned, and any binned # observation will be skipped, hence the counts cube should be # empty. - self._test_observation(cpy_bin, 0) - self._test_cube(cpy_bin.cube(), 0) + self._check_observation(cpy_bin, 0) + self._check_cube(cpy_bin.cube(), 0) # Clear and run copy of ctbin again #cpy_bin.clear() @@ -115,7 +172,7 @@ class Test(gammalib.GPythonTestSuite): # Load counts cube and check content. evt = gammalib.GCTAEventCube("cntmap.fits") - self._test_cube(evt, 5542) + self._check_cube(evt, 5542) # Prepare observation container for stacked analysis cta = gammalib.GCTAObservation(self._events_name) @@ -147,13 +204,12 @@ class Test(gammalib.GPythonTestSuite): # Check content of observation and cube (need multiplier=3 since # three identical observations have been appended) - self._test_observation(bin, 5542, multiplier=3) - self._test_cube(bin.cube(), 5542, multiplier=3) + self._check_observation(bin, 5542, multiplier=3) + self._check_cube(bin.cube(), 5542, multiplier=3) - # Set-up ctbin with invalid event file - bin = ctools.ctbin() - bin["inobs"] = "event_file_that_does_not_exist.fits" - bin["outcube"] = "cntmap.fits" + # Set-up ctbin using an observation container + bin = ctools.ctbin(obs) + bin["outcube"] = "cntmap2.fits" bin["ebinalg"] = "LOG" bin["emin"] = 0.1 bin["emax"] = 100.0 @@ -166,21 +222,29 @@ class Test(gammalib.GPythonTestSuite): bin["xref"] = 83.63 bin["yref"] = 22.01 - # Run ctbin tool - self.test_try("Run ctbin with invalid event file") - try: - bin.run() - self.test_try_failure() - except: - self.test_try_success() + # Execute ctbin tool + bin.execute() + + # Load counts cube and check content. + evt = gammalib.GCTAEventCube("cntmap2.fits") + self._check_cube(evt, 5542, multiplier=3) + + # Clear tool + #bin.clear() +#what(): *** ERROR in GApplicationPars::write(std::string&): Unable to open parameter file "pfiles/". + + # Check content of observation and cube + #self.test_value(bin.obs().size(), 0, "There is no observation") + #self.test_value(bin.cube().size(), 0, "0 event bins") + #self.test_value(bin.cube().number(), 0, "0 events") # Return return # Check observation - def _test_observation(self, ctbin, nevents, multiplier=1): + def _check_observation(self, ctbin, nevents, multiplier=1): """ - Test content of an observation. + Check content of an observation. Args: ctbin: ctbin instance @@ -206,15 +270,15 @@ class Test(gammalib.GPythonTestSuite): "Pointing Declination is 22.01 deg") # Test event cube - self._test_cube(obs.events(), nevents, multiplier=multiplier) + self._check_cube(obs.events(), nevents, multiplier=multiplier) # Return return # Check event cube - def _test_cube(self, cube, nevents, multiplier=1): + def _check_cube(self, cube, nevents, multiplier=1): """ - Test content of event cube." + Check content of event cube." Args: cube: Event cube diff --git a/test/test_ctobssim.py b/test/test_ctobssim.py index 5816c06..4413eee 100755 --- a/test/test_ctobssim.py +++ b/test/test_ctobssim.py @@ -1,6 +1,6 @@ #! /usr/bin/env python # ========================================================================== -# This scripts performs unit tests for the ctobssim module. +# This scripts performs unit tests for the ctobssim tool. # # Copyright (C) 2014-2016 Juergen Knoedlseder # @@ -18,17 +18,21 @@ # along with this program. If not, see . # # ========================================================================== +import os import gammalib import ctools from cscripts import obsutils -# ============================== # -# Test class for ctobssim module # -# ============================== # +# ============================ # +# Test class for ctobssim tool # +# ============================ # class Test(gammalib.GPythonTestSuite): """ - Test class for ctobssim module. + Test class for ctobssim tools. + + This test class makes unit tests for the ctobssim tool by using it from + the command line and from Python. """ # Constructor @@ -56,16 +60,67 @@ class Test(gammalib.GPythonTestSuite): self.name("ctobssim") # Append tests - self.append(self._test_ctobssim, "Test ctobssim") - #self.append(self._test_container, "Test ctobssim on observation container") + self.append(self._test_ctobssim_cmd, "Test ctobssim on command line") + self.append(self._test_ctobssim_python, "Test ctobssim from Python") # Return return - # Test ctobssim - def _test_ctobssim(self): + # Test ctobssim on command line + def _test_ctobssim_cmd(self): """ - Test ctobssim. + Test ctobssim on the command line. + """ + # Kluge to set the command (installed version has no README file) + if os.path.isfile("README"): + ctobssim = "../src/ctobssim/ctobssim" + else: + ctobssim = "ctobssim" + + # Setup ctobssim command + cmd = ctobssim+' inmodel="data/crab.xml" outevents="events.fits"'+ \ + ' caldb="irf" irf="cta_dummy_irf" '+ \ + ' ra=83.63 dec=22.01 rad=10.0'+ \ + ' tmin=0.0 tmax=1800.0 emin=0.1 emax=100.0' + + # Execute ctobssim, make sure we catch any exception + try: + rc = os.system(cmd+" >/dev/null 2>&1") + except: + pass + + # Check if execution was successful + self.test_assert(rc == 0, + "Successful ctobssim execution on command line") + + # Load counts cube and check content. + evt = gammalib.GCTAEventList("events.fits") + self._test_list(evt, 4119) + + # Setup ctobssim command + cmd = ctobssim+' inmodel="model_that_does_not_exist.xml"'+ \ + ' outevents="events.fits"'+ \ + ' caldb="irf" irf="cta_dummy_irf" '+ \ + ' ra=83.63 dec=22.01 rad=10.0'+ \ + ' tmin=0.0 tmax=1800.0 emin=0.1 emax=100.0' + + # Execute ctobssim, make sure we catch any exception + try: + rc = os.system(cmd+" >/dev/null 2>&1") + except: + pass + + # Check if execution failed + self.test_assert(rc != 0, + "Failure of ctobssim execution on command line") + + # Return + return + + # Test ctobssim from Python + def _test_ctobssim_python(self): + """ + Test ctobssim from Python. """ # Set-up ctobssim sim = ctools.ctobssim() diff --git a/test/test_python_cscripts.py b/test/test_python_cscripts.py index c3266af..5067331 100755 --- a/test/test_python_cscripts.py +++ b/test/test_python_cscripts.py @@ -18,12 +18,13 @@ # along with this program. If not, see . # # ========================================================================== +import os +import sys import gammalib import cscripts -import sys -import os sys.path.append(os.path.dirname(os.path.abspath(__file__))) import test_cscaldb +import test_cslightcrv # ================== # @@ -37,13 +38,16 @@ def test(installed=False): suites = gammalib.GTestSuites("cscripts Python module unit testing") # Allocate test suites and append them to the container - suite_cscaldb = test_cscaldb.Test() + suite_cscaldb = test_cscaldb.Test() + suite_cslightcrv = test_cslightcrv.Test() # Setup unit tests suite_cscaldb.set() + suite_cslightcrv.set() # Append tests to container suites.append(suite_cscaldb) + suites.append(suite_cslightcrv) # If we have an installed version then create a temporary directory and # copy over all information that is needed @@ -60,9 +64,13 @@ def test(installed=False): dirname = os.path.dirname(testdir) # Copy over test data and irf - #os.system("cp -r %s %s" % (dirname+"/data", "data")) + os.system("cp -r %s %s" % (dirname+"/data", "data")) #os.system("cp -r %s %s" % (dirname+"/irf", "irf")) + # ... otherwise set the calibration database + else: + os.environ['CALDB'] = "../caldb" + # Set PFILES environment variable try: os.mkdir("pfiles") @@ -97,8 +105,6 @@ def test(installed=False): # Main routine entry point # # ======================== # if __name__ == '__main__': - """ - Perform unit testing for cscripts used from Python. - """ + # Run tests test() diff --git a/test/test_python_ctools.py b/test/test_python_ctools.py index 1f72500..0a1e579 100755 --- a/test/test_python_ctools.py +++ b/test/test_python_ctools.py @@ -40,15 +40,15 @@ import test_cterror import test_pipelines -# ================== # -# Perform unit tests # -# ================== # +# ========================= # +# Perform ctools unit tests # +# ========================= # def test(installed=False): """ - Perform unit testing for Python interface. + Perform unit testing for ctools. """ # Allocate test suite container - suites = gammalib.GTestSuites("ctools Python module unit testing") + suites = gammalib.GTestSuites("ctools unit testing") # Allocate test suites and append them to the container suite_ctobssim = test_ctobssim.Test() @@ -153,8 +153,6 @@ def test(installed=False): # Main routine entry point # # ======================== # if __name__ == '__main__': - """ - Perform unit testing for Python interface. - """ + # Run tests test()