Add csobs2caldb unit tests and change user parameter interface.

Jürgen Knödlseder
1 parent f662447a
Showing 7 changed files with 339 additions and 80 deletions
cscripts/Makefile.am
cscripts/csobs2caldb.par
cscripts/csobs2caldb.py
doc/source/reference_manual/csobs2caldb.rst
test/test_cscripts.sh.in
test/test_csobs2caldb.py
test/test_python_cscripts.py
@@ -27,10 +27,12 @@ linkscripts = cspull cssens cstsdist cscaldb csresmap \
 testscripts = tests/__init__.py \
               ../test/test_python_cscripts.py \
               ../test/test_cscaldb.py \
+              ../test/test_csobs2caldb.py \
               ../test/test_cslightcrv.py
 # Test data to be distributed in $(pyexecdir)/cscripts/tests/data
 testdata    = ../test/data/crab.xml \
+              ../test/data/obs_unbinned.xml \
               ../test/data/crab_events.fits.gz \
               ../test/data/crab_cntmap.fits.gz \
               ../test/data/exclusion.reg
@@ -15,11 +15,13 @@
 #                                                                           #
 #############################################################################
 #
-# General parameters
-#===================
+# Input and output parameters
+#============================
 inobs,   f, a, obs.xml,,, "Observation definition XML file"
-rspname, s, a, NONE,,, "Response output name (e.g. Zenith50)"
-outirfs, f, h, irf_file.fits,,, "Output IRF file name"
+caldb,   s, a, cta,,, "Calibration database"
+irf,     s, a, Standard,,, "Instrument response function (e.g. Standard)"
+outfile, f, h, irf_file.fits,,, "Output IRF file name"
+rootdir, s, h, ,,, "Optional CALDB root directory"
 index,   i, h, 0,,, "Index of observation to be used from XML file"
 #
 #!/usr/bin/env python
 # ==========================================================================
-# Generation of an caldb entry from general IACT observation.
+# Creates a calibration database entry for an IACT observation.
 #
 # Copyright (C) 2015-2016 Michael Mayer
 #
@@ -30,7 +30,7 @@ import datetime
 # ================= #
 class csobs2caldb(ctools.cscript):
     """
-    Creates a calibration database entry for an observation.
+    Creates a calibration database entry for an IACT observation.
     The creation of a calibration database entry is useful for performing
     simulations for current Imaging Air Cherenkov Telescopes (IACTs).
@@ -50,17 +50,17 @@ class csobs2caldb(ctools.cscript):
         self._version = "1.1.0"
         # Initialise members
-        self._obs         = []     # Observation container
-        self._observation = None
-        self._mission     = "CTA"
-        self._instrument  = ""
-        self._rspname     = ""
-        self._outfile     = ""
+        self._observation = gammalib.GCTAObservation()
+        self._mission     = "cta"
+        self._caldb       = "cta"
+        self._outfile     = gammalib.GFilename("irf_file.fits")
         self._base_dir    = ""
         self._cal_dir     = ""
+
+        # Initialise observation container from constructor arguments.
+        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)
         # Return
@@ -70,36 +70,64 @@ class csobs2caldb(ctools.cscript):
     # Private methods
     def _get_parameters(self):
         """
-        Get parameters from parfile and setup the observation.
+        Get parameters from parfile.
+
+        Raises:
+            ValueError, IndexError, & RuntimeError.
         """
-        # Load observation container
-        self._obs = gammalib.GObservations(self["inobs"].filename())
+        # Load observation container from "inobs" parameter in case that
+        # the observation container is still empty.
+        if self._obs.is_empty():
+            self._obs = gammalib.GObservations(self["inobs"].filename())
+
+        # Raise an exception if the observation container is empty
+        if self._obs.is_empty():
+            raise ValueError('No or empty observation provided for '
+                             'parameter "inobs".')
-        # Get index of observation to be used (0=default)
-        inx = self["index"].integer()
+        # Get index of observation to be used
+        index = self["index"].integer()
+
+        # Raise an exception if the index is not valid
+        if index < 0 or index >= self._obs.size():
+            raise IndexError('Parameter "index=%d" outside the validity '
+                             'range [0,%d].' % (index, self._obs.size()))
         # Get observation
-        self._observation = self._obs[inx]
-        
-        # Get instrument name
-        self._instrument = self._observation.instrument()
-        
-        # Get config and file name
-        self._rspname = self["rspname"].string()
-        self._outfile = self["outirfs"].filename()
-        
-        # Make sure we have a cta observation
+        self._observation = self._obs[index]
+
+        # Make sure we have a CTA observation
         if not self._observation.classname() == "GCTAObservation":
-            raise RuntimeError("Input observation not of type "+
-                               "\"GCTAObservation\".")
+            raise RuntimeError('Input observation not of type '
+                               '"GCTAObservation".')
-        # Make sure we have an irf response
+        # Make sure we have an IRF response associated with the CTA
+        # observation
         if not self._observation.response().classname() == "GCTAResponseIrf":
-            raise RuntimeError("Response of input observation not of "+
-                               "type \"GCTAResponseIrf\".")
+            raise RuntimeError('Response of input observation not of '
+                               'type "GCTAResponseIrf".')
+        
+        # Get calibration database name. If the "caldb" parameter is "NONE"
+        # or empty then use the instrument name from the observation as
+        # calibration database name.
+        self._caldb = self["caldb"].string()
+        if (self._caldb == "NONE" or len(self._caldb) == 0): 
+            self._caldb = self._observation.instrument().lower()
+
+        # Get instrument response function file name
+        self._outfile = self["outfile"].filename()
+
+        # 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["irf"].string()
+        self["rootdir"].string()
-        # Get other parameters
-        self._clobber = self["clobber"].boolean()
+        # Write input parameters and header into logger
+        if self._logTerse():
+            self._log_parameters()
+            self._log("\n")
         # Return
         return
@@ -122,10 +150,11 @@ class csobs2caldb(ctools.cscript):
         fname_edisp = rsp.edisp().filename()
         fname_bkg   = rsp.background().filename()
-        # Log filename instances
-        if self._logTerse():
+        # Log filenames
+        if self._logNormal():
+            self._log.header3("IRF input files")
             self._log(gammalib.parformat("Effective area"))
-            self._log(repr(fname_aeff))
+            self._log(fname_aeff.url())
             self._log("\n")
             self._log(gammalib.parformat("Point spread function"))
             self._log(fname_psf.url())
@@ -142,15 +171,35 @@ class csobs2caldb(ctools.cscript):
         fits_psf   = gammalib.GFits(fname_psf)
         fits_edisp = gammalib.GFits(fname_edisp)
         fits_bkg   = gammalib.GFits(fname_bkg)
+
+        # Get extension names
+        ext_aeff  = fname_aeff.extname("EFFECTIVE AREA")
+        ext_psf   = fname_psf.extname("POINT SPREAD FUNCTION")
+        ext_edisp = fname_edisp.extname("ENERGY DISPERSION")
+        ext_bkg   = fname_bkg.extname("BACKGROUND")
         # Create empty FITS file
         fits = gammalib.GFits()
         # Append IRF component to FITS file
-        fits.append(fits_aeff[fname_aeff.extname("EFFECTIVE AREA")])
-        fits.append(fits_psf[fname_psf.extname("POINT SPREAD FUNCTION")])
-        fits.append(fits_edisp[fname_edisp.extname("ENERGY DISPERSION")])
-        fits.append(fits_bkg[fname_bkg.extname("BACKGROUND")])
+        fits.append(fits_aeff[ext_aeff])
+        fits.append(fits_psf[ext_psf])
+        fits.append(fits_edisp[ext_edisp])
+        fits.append(fits_bkg[ext_bkg])
+
+        # Log resulting FITS file
+        if self._logNormal():
+            self._log(str(fits))
+            self._log("\n")
+        if self._logExplicit():
+            self._log(str(fits[ext_aeff].header()))
+            self._log("\n")
+            self._log(str(fits[ext_psf].header()))
+            self._log("\n")
+            self._log(str(fits[ext_edisp].header()))
+            self._log("\n")
+            self._log(str(fits[ext_bkg].header()))
+            self._log("\n")
         # Return fits file
         return fits
@@ -235,7 +284,7 @@ class csobs2caldb(ctools.cscript):
             # Set element
             table["TELESCOP"][row] = self._mission
-            table["INSTRUME"][row] = self._instrument
+            table["INSTRUME"][row] = self._caldb
             table["DETNAM"][row]   = "NONE"
             table["FILTER"][row]   = "NONE"
             table["CAL_DEV"][row]  = "ONLINE"
@@ -245,7 +294,7 @@ class csobs2caldb(ctools.cscript):
             table["CAL_VST"][row]  = time.split(".")[0]
             table["REF_TIME"][row] = 51544.0
             table["CAL_QUAL"][row] = 0
-            table["CAL_CBD"][row]  = "NAME("+self._rspname+")"
+            table["CAL_CBD"][row]  = "NAME("+self["irf"].string()+")"
             table["CAL_DATE"][row] = today.replace("-","/")[2:]
             table["CAL_DIR"][row]  = self._cal_dir
             table["CAL_FILE"][row] = self._outfile
@@ -253,22 +302,30 @@ class csobs2caldb(ctools.cscript):
         # Add effective area information
         row = row_index-4
         table["CAL_CNAM"][row] = "EFF_AREA"
-        table["CAL_DESC"][row] = self._instrument+" effective area"
+        table["CAL_DESC"][row] = self._caldb+" effective area"
         # Add point spread function information
         row = row_index-3
         table["CAL_CNAM"][row] = "RPSF"
-        table["CAL_DESC"][row] = self._instrument+" point spread function"
+        table["CAL_DESC"][row] = self._caldb+" point spread function"
         # Add energy dispersion information
         row = row_index-2
         table["CAL_CNAM"][row] = "EDISP"
-        table["CAL_DESC"][row] = self._instrument+" energy dispersion"
+        table["CAL_DESC"][row] = self._caldb+" energy dispersion"
         # Add background information
         row = row_index-1
         table["CAL_CNAM"][row] = "BGD"
-        table["CAL_DESC"][row] = self._instrument+" background"
+        table["CAL_DESC"][row] = self._caldb+" background"
+
+        # Log resulting FITS table and header
+        if self._logNormal():
+            self._log(str(table))
+            self._log("\n")
+        if self._logExplicit():
+            self._log(str(table.header()))
+            self._log("\n")
         # Return CIF FITS file
         return cif
@@ -283,34 +340,40 @@ class csobs2caldb(ctools.cscript):
             self._log.header2("Creating directory structure")
         # Create calibration database        
-        caldb = gammalib.GCaldb()
+        caldb = gammalib.GCaldb(self["rootdir"].string())
         # Set calibration directory 
         self._cal_dir  = "data"
         self._cal_dir += "/"+self._mission.lower()
-        self._cal_dir += "/"+self._instrument.lower()
-        self._cal_dir += "/bcf/"+self._rspname
+        self._cal_dir += "/"+self._caldb.lower()
+        self._cal_dir += "/bcf/"+self["irf"].string()
         # Set absolute path
         self._base_dir = caldb.rootdir() +"/data"
         self._base_dir += "/"+self._mission.lower()
-        self._base_dir += "/"+self._instrument.lower()
+        self._base_dir += "/"+self._caldb.lower()
         # Set directory for irf file
         self._rsp_dir = caldb.rootdir() + "/" + self._cal_dir
+
+        # Log resulting FITS table
+        if self._logNormal():
+            self._log(gammalib.parformat("Calibration directory"))
+            self._log(self._cal_dir)
+            self._log("\n")
+            self._log(gammalib.parformat("Base directory"))
+            self._log(self._base_dir)
+            self._log("\n")
+            if not os.path.isdir(self._rsp_dir):
+                self._log(gammalib.parformat("IRF directory"))
+            else:
+                self._log(gammalib.parformat("IRF directory (existing)"))
+            self._log(self._rsp_dir)
+            self._log("\n")
-        # Create directories and log information
+        # Create IRF directory is it does not yet exist
         if not os.path.isdir(self._rsp_dir):
-            if self._logExplicit():
-                self._log(gammalib.parformat("Directory"))
-                self._log(self._rsp_dir)
-                self._log("\n")
             os.makedirs(self._rsp_dir)
-        else:
-            if self._logExplicit():
-                self._log(gammalib.parformat("Directory (existing)"))
-                self._log(self._rsp_dir)
-                self._log("\n")
         # Return
         return
@@ -330,8 +393,7 @@ class csobs2caldb(ctools.cscript):
         # Write input parameters and header into logger
         if self._logTerse():
-            self._log_parameters()
-            self._log("\n\n")
+            self._log("\n")
             self._log.header1("Creating CALDB entry")
         # Create directory structure
@@ -359,7 +421,7 @@ class csobs2caldb(ctools.cscript):
         filename = self._rsp_dir + "/" + self._outfile
         # Write filenames into logger
-        if self._logTerse():
+        if self._logNormal():
             self._log(gammalib.parformat("CALDB index file"))
             self._log(self._caldb_inx.filename().url())
             self._log("\n")
@@ -368,10 +430,10 @@ class csobs2caldb(ctools.cscript):
             self._log("\n")
         # Save caldb index file
-        self._caldb_inx.save(self._clobber)
+        self._caldb_inx.save(self._clobber())
         # Save response file
-        self._irf_fits.saveto(filename, self._clobber)
+        self._irf_fits.saveto(filename, self._clobber())
         # Return
         return
@@ -381,7 +443,7 @@ class csobs2caldb(ctools.cscript):
         Execute the script.
         """
         # Open logfile
-        self._logFileOpen()
+        self.logFileOpen()
         # Run the script
         self.run()
@@ -28,12 +28,18 @@ General parameters
 ``inobs [file]``
     Input event list, counts cube or observation definition XML file.
-``rspname [string]``
-    Response output name (e.g. Zenith50).
+``caldb [string]``
+    Calibration database.
+
+``irf [string]``
+    Instrument response function (e.g. Standard).
 ``(outirfs = irf_file.fits) [file]``
     Output IRF file name.
+``(rootdir = ) [string]``
+    Optional CALDB root directory.
+
 ``(index = 0) [integer]``
     Index of observation to be used from XML file.
@@ -77,14 +77,15 @@ $ECHO -n &quot;Test csobs2caldb: &quot;
 chmod +x @abs_top_srcdir@/cscripts/csobs2caldb.py
 @PYTHON_COVERAGE@@abs_top_srcdir@/cscripts/csobs2caldb.py \
   inobs="data/obs_unbinned.xml" \
-  rspname="test_irfs" \
+  caldb="cta" \
+  irf="Standard" \
   chatter=4 > csobs2caldb.log
 $ECHO -n "."
-if [ -s "$CALDB/data/cta/cta/bcf/test_irfs/irf_file.fits" ]
+if [ -s "$CALDB/data/cta/cta/bcf/Standard/irf_file.fits" ]
 then
   $ECHO -n "."
 else
-  $ECHO " $CALDB/data/cta/cta/bcf/test_irfs/irf_file.fits file is not valid"
+  $ECHO " $CALDB/data/cta/cta/bcf/Standard/irf_file.fits file is not valid"
   exit 1
 fi
 $ECHO " ok"
@@ -282,7 +283,7 @@ chmod +x @abs_top_srcdir@/cscripts/cspull.py
   outfile="pull.dat" \
   ntrials="3" \
   caldb="cta" \
-  irf="test_irfs" \
+  irf="Standard" \
   ra=83.6331 \
   dec=22.0145 \
   emin=0.1 \
+#! /usr/bin/env python
+# ==========================================================================
+# This scripts performs unit tests for the csobs2caldb 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 <http://www.gnu.org/licenses/>.
+#
+# ==========================================================================
+import os
+import gammalib
+import cscripts
+
+
+# ================================= #
+# Test class for csobs2caldb script #
+# ================================= #
+class Test(gammalib.GPythonTestSuite):
+    """
+    Test class for csobs2caldb script.
+
+    This test class makes unit tests for the csobs2caldb 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._inobs   = "data/obs_unbinned.xml"
+
+        # Clean calibration database
+        os.system("rm -rf caldb >/dev/null 2>&1")
+        os.system("mkdir -p caldb >/dev/null 2>&1")
+
+        # Return
+        return
+
+    # Set test functions
+    def set(self):
+        """
+        Set all test functions.
+        """
+        # Set test name
+        self.name("csobs2caldb")
+
+        # Append tests
+        self.append(self._test_cmd, "Test csobs2caldb on command line")
+        self.append(self._test_python, "Test csobs2caldb from Python")
+
+        # Return
+        return
+
+    # Test csobs2caldb on command line
+    def _test_cmd(self):
+        """
+        Test csobs2caldb on the command line.
+        """
+        # Kluge to set the command (installed version has no README file)
+        if os.path.isfile("README"):
+            csobs2caldb = "../cscripts/csobs2caldb.py"
+        else:
+            csobs2caldb = "csobs2caldb"
+
+        # Setup csobs2caldb command
+        cmd = csobs2caldb+' inobs="'+self._inobs+'"'+ \
+                          ' caldb="cta" irf="Cmd1" rootdir="caldb"'+ \
+                          ' logfile="csobs2caldb_cmd1.log" chatter=1'
+
+        # Execute csobs2caldb, 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 csobs2caldb execution on command line")
+
+        # Check response
+        self._check_response("cta", "Cmd1")
+
+        # Setup csobs2caldb command
+        cmd = csobs2caldb+' inobs="xml_file_that_does_not_exist.xml"'+ \
+                          ' caldb="cta" irf="Cmd2" rootdir="caldb"'+ \
+                          ' logfile="csobs2caldb_cmd2.log"'
+
+        # Execute csobs2caldb, 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 csobs2caldb execution on command line")
+
+        # Return
+        return
+
+    # Test csobs2caldb from Python
+    def _test_python(self):
+        """
+        Test csobs2caldb from Python.
+        """
+        # Set-up csobs2caldb using "inobs" parameter
+        irf = cscripts.csobs2caldb()
+        irf["inobs"]   = self._inobs
+        irf["caldb"]   = "cta"
+        irf["irf"]     = "Py1"
+        irf["rootdir"] = "caldb"
+        irf["logfile"] = "csobs2caldb_py1.log"
+        irf["chatter"] = 2
+
+        # Run csobs2caldb script and save IRF
+        irf.logFileOpen()   # Make sure we get a log file
+        irf.run()
+        irf.save()
+
+        # Check response
+        self._check_response("cta", "Py1")
+
+        # Load observation container
+        obs = gammalib.GObservations(self._inobs)
+
+        # Set-up csobs2caldb using observation container
+        irf = cscripts.csobs2caldb(obs)
+        irf["caldb"]   = "cta"
+        irf["irf"]     = "Py2"
+        irf["rootdir"] = "caldb"
+        irf["logfile"] = "csobs2caldb_py2.log"
+        irf["chatter"] = 3
+
+        # Execute csobs2caldb script
+        irf.execute()
+
+        # Check response
+        self._check_response("cta", "Py2")
+
+        # Return
+        return
+
+    # Check response
+    def _check_response(self, caldb, irf):
+        """
+        Check response.
+        """
+        # Open calibration database
+        db = gammalib.GCaldb("caldb")
+        db.open("cta", caldb)
+
+        # Get filenames of response components
+        expr  = "NAME("+irf+")"
+        aeff  = gammalib.GFilename(db.filename("","","EFF_AREA","","",expr))
+        psf   = gammalib.GFilename(db.filename("","","RPSF","","",expr))
+        edisp = gammalib.GFilename(db.filename("","","EDISP","","",expr))
+        bgd   = gammalib.GFilename(db.filename("","","BGD","","",expr))
+
+        # Check whether files exist
+        self.test_assert(aeff.exists(), "Effective area file exists")
+        self.test_assert(psf.exists(), "Point spread function file exists")
+        self.test_assert(edisp.exists(), "Energy dispersion file exists")
+        self.test_assert(bgd.exists(), "Background file exists")
+        
+        # Return
+        return
 #! /usr/bin/env python
 # ==========================================================================
-# This scripts performs unit tests for cscripts used from Python.
+# This scripts performs unit tests for cscripts.
 #
 # Copyright (C) 2016 Juergen Knoedlseder
 #
@@ -24,6 +24,7 @@ import gammalib
 import cscripts
 sys.path.append(os.path.dirname(os.path.abspath(__file__)))
 import test_cscaldb
+import test_csobs2caldb
 import test_cslightcrv
@@ -32,21 +33,24 @@ import test_cslightcrv
 # ================== #
 def test(installed=False):
     """
-    Perform unit testing for Python interface.
+    Perform cscripts unit testing.
     """
     # Allocate test suite container
-    suites = gammalib.GTestSuites("cscripts Python module unit testing")
+    suites = gammalib.GTestSuites("cscripts unit testing")
     # Allocate test suites and append them to the container
-    suite_cscaldb    = test_cscaldb.Test()
-    suite_cslightcrv = test_cslightcrv.Test()
+    suite_cscaldb     = test_cscaldb.Test()
+    suite_csobs2caldb = test_csobs2caldb.Test()
+    suite_cslightcrv  = test_cslightcrv.Test()
     # Setup unit tests
     suite_cscaldb.set()
+    suite_csobs2caldb.set()
     suite_cslightcrv.set()
     # Append tests to container
     suites.append(suite_cscaldb)
+    suites.append(suite_csobs2caldb)
     suites.append(suite_cslightcrv)
     # If we have an installed version then create a temporary directory and