Merge branch 'master' into devel

Jürgen Knödlseder
2 parents 82ae545c c36b288f
Showing 83 changed files with 1524 additions and 1444 deletions
ChangeLog
NEWS
README
caldb/Makefile.am
caldb/caldb.indx
caldb/cta_dummy_irf.dat
caldb/data/cta/prod2/bcf/North_0.5h/irf_file.fits.gz
caldb/data/cta/prod2/bcf/North_50h/irf_file.fits.gz
caldb/data/cta/prod2/bcf/North_5h/irf_file.fits.gz
caldb/data/cta/prod2/bcf/South_0.5h/irf_file.fits.gz
caldb/data/cta/prod2/bcf/South_50h/irf_file.fits.gz
caldb/data/cta/prod2/bcf/South_5h/irf_file.fits.gz
caldb/data/cta/prod2/caldb.indx
configure.ac
ctools.pc.in
doc/Doxyfile
doc/source/conf.py
doc/source/download.rst
doc/source/reference_manual/cslightcrv.rst
doc/source/reference_manual/csobsdef.rst
-2015-05-20  Juergen Knoedlseder  <jurgen.knodlseder@irap.omp.eu>
+2015-05-22  Juergen Knoedlseder  <jurgen.knodlseder@irap.omp.eu>
-        * Version 1.0.0 released
+        * Version 0.9.0 released
         ========================
-        Adapt to GammaLib 1.0.0 interface
+        Adapt to GammaLib 0.10.0 interface
         Do not import automatically Python scripts
         Support INDEF/NONE/UNDEF/UNDEFINED parameters in ctselect
         Modify cssens parameter interface (#1426)
-New Features and Important Changes in ctools 1.0.0
+New Features and Important Changes in ctools 0.9.0
 Author(s) : Juergen Knoedlseder
-20 May 2015
+22 May 2015
 1. Introduction
@@ -141,7 +141,8 @@ actually returns the differential flux for a given set of energies.
 20. caldb
 ---------
-None
+Official Prod2 instrument response functions for the Southern and Northern
+array have been added.
 21. models
 ctools information
 ==================
-Version:             00-08-01 (8 January 2015)
+Version:             0.9.0 (20 May 2015)
 Author:              Juergen Knodlseder (jurgen.knodlseder@irap.omp.eu)
-GammaLib dependency: 00-09-00
+GammaLib dependency: 0.10.0
 License information
@@ -34,20 +34,27 @@ for the high-level analysis of astronomical gamma-ray data.
 The following tools are available:
-    ctbin    - event binning
-    ctlike   - maximum likelihood model fitting
-    ctmodel  - generation of model counts map
-    ctobssim - simulation of CTA observations
-    ctselect - event selection
-    ctskymap - CTA sky mapping tool
+    ctbin       - event binning
+    ctbkgcube   - generate a background cube
+    ctbutterfly - create a butterfly
+    ctcubemask  - mask bins in binned analysis
+    ctexpcube   - generate an exposure cube
+    ctlike      - maximum likelihood model fitting
+    ctmodel     - generation of model counts map
+    ctobssim    - simulation of CTA observations
+    ctpsfcube   - generate a PSF cube
+    ctselect    - event selection
+    ctskymap    - CTA sky mapping tool
+    cttsmap     - generate a TS map
+    ctulimit    - compute upper limits
 Web sites
 =========
-http://cta.irap.omp.eu/ctools                   - ctools documentation
-https://cta-redmine.irap.omp.eu/projects/ctools - for ctools development
-http://gammalib.sourceforge.net                 - for GammaLib documentation
-https://sourceforge.net/projects/gammalib/      - for GammaLib downloads
+http://cta.irap.omp.eu/ctools                     - ctools documentation
+https://cta-redmine.irap.omp.eu/projects/ctools   - for ctools development
+https://cta-redmine.irap.omp.eu/projects/gammalib - for GammaLib documentation
+https://sourceforge.net/projects/gammalib/        - for GammaLib downloads
 Prerequisites
 ## Process this file with automake to produce Makefile.in
 # Calibration files to distribute
-caldb_indx = caldb.indx
-caldb_irf  = cta_dummy_irf.dat
-
-# Distribute calibration database files
-EXTRA_DIST     = $(caldb_indx) $(caldb_irf)
-caldbindxdir   = $(datadir)/caldb/data/cta/dummy
-caldbindx_DATA = $(caldb_indx)
-caldbirfdir    = $(datadir)/caldb/data/cta/dummy/bcf
-caldbirf_DATA  = $(caldb_irf)
+caldb_prod2_indxdir      = $(datadir)/caldb/data/cta/prod2
+caldb_prod2_indx_DATA    = data/cta/prod2/caldb.indx
+caldb_prod2_north1dir    = $(datadir)/caldb/data/cta/prod2/bcf/North_0.5h
+caldb_prod2_north1_DATA  = data/cta/prod2/bcf/North_0.5h/irf_file.fits.gz
+caldb_prod2_north2dir    = $(datadir)/caldb/data/cta/prod2/bcf/North_5h
+caldb_prod2_north2_DATA  = data/cta/prod2/bcf/North_5h/irf_file.fits.gz
+caldb_prod2_north3dir    = $(datadir)/caldb/data/cta/prod2/bcf/North_50h
+caldb_prod2_north3_DATA  = data/cta/prod2/bcf/North_50h/irf_file.fits.gz
+#
+caldb_prod2_south1dir    = $(datadir)/caldb/data/cta/prod2/bcf/South_0.5h
+caldb_prod2_south1_DATA  = data/cta/prod2/bcf/South_0.5h/irf_file.fits.gz
+caldb_prod2_south2dir    = $(datadir)/caldb/data/cta/prod2/bcf/South_5h
+caldb_prod2_south2_DATA  = data/cta/prod2/bcf/South_5h/irf_file.fits.gz
+caldb_prod2_south3dir    = $(datadir)/caldb/data/cta/prod2/bcf/South_50h
+caldb_prod2_south3_DATA  = data/cta/prod2/bcf/South_50h/irf_file.fits.gz
-log(E) Area  r68  r80  ERes. BG Rate  Diff Sens
--1.7  261.6  0.3621  0.4908 0.5134 1.89924e-02  6.88237e-11
--1.5  5458.2  0.2712  0.3685 0.4129 1.00972e-01  1.72717e-11
--1.3  15590.0  0.1662  0.2103 0.2721 5.75623e-02  6.16963e-12
--1.1  26554.1  0.1253  0.1567 0.2611 2.13008e-02  2.89932e-12
--0.9  52100.5  0.1048  0.1305 0.1987 8.87292e-03  1.39764e-12
--0.7  66132.1  0.0827  0.1024 0.1698 1.09756e-03  6.03531e-13
--0.5  108656.8  0.0703  0.0867 0.1506 4.84287e-04  3.98147e-13
--0.3  129833.0  0.0585  0.0722 0.1338 1.57546e-04  3.23090e-13
--0.1  284604.3  0.0531  0.0656 0.1008 1.36703e-04  2.20178e-13
-0.1  263175.3  0.0410  0.0506 0.0831 2.09694e-05  1.87452e-13
-0.3  778048.6  0.0470  0.0591 0.0842 6.92374e-05  1.53976e-13
-0.5  929818.8  0.0391  0.0492 0.0650 1.45844e-05  1.18947e-13
-0.7  1078450.0  0.0335  0.0415 0.0541 1.15959e-05  1.51927e-13
-0.9  1448579.1  0.0317  0.0397 0.0516 4.71231e-06  1.42439e-13
-1.1  1899905.0  0.0290  0.0372 0.0501 8.14997e-06  1.96670e-13
-1.3  2476403.8  0.0285  0.0367 0.0538 5.91940e-06  2.20695e-13
-1.5  2832570.6  0.0284  0.0372 0.0636 7.33847e-06  3.22523e-13
-1.7  3534065.3  0.0290  0.0386 0.0731 1.34549e-05  4.84153e-13
-1.9  3250103.4  0.0238  0.0308 0.0729 4.42228e-06  6.26265e-13
-2.1  3916071.6  0.0260  0.0354 0.0908 2.26648e-06  7.69921e-13
----------------------------------------------
-Notes
- 1) log(E) = log10(E/TeV) - bin centre
- 2) Eff Area - in square metres after background cut (no theta cut)
- 3) Ang. Res - 68% containment radius of gamma-ray PSF post cuts - in degrees
- 4) Ang. Res - 80% containment radius of gamma-ray PSF post cuts - in degrees
- 5) Fractional Energy Resolution (rms)
- 6) BG Rate  - inside point-source selection region - post call cuts - in Hz
- 7) Diff Sens - differential sensitivity for this bin expressed as E^2 dN/dE
- - in erg cm^-2 s^-1 - for a 50 hours exposure - 5 sigma significance including
- systematics and statistics and at least 10 photons.
@@ -29,7 +29,7 @@
 #############################################################################
 # Initialisation                                                            #
 #############################################################################
-AC_INIT([ctools], [1.0.0], [jurgen.knodlseder@irap.omp.eu], [ctools])
+AC_INIT([ctools], [0.9.0], [jurgen.knodlseder@irap.omp.eu], [ctools])
 AC_CONFIG_SRCDIR([src/ctselect/ctselect.cpp])
 AC_CONFIG_MACRO_DIR([m4])
 #AM_INIT_AUTOMAKE([gnu no-dependencies])
@@ -6,9 +6,9 @@ includedir=${prefix}/include
 Name: ctools
 Description: CTA science analysis tools
 URL: http://cta.irap.omp.eu/ctools
-Version: 1.0.0
+Version: 0.9.0
 Requires: 
-Requires.private: gammalib >= 1.0.0
+Requires.private: gammalib >= 0.10.0
 Conflicts: 
 Cflags: -I${includedir}/ctools @OPENMP_CXXFLAGS@
 Libs: -L${libdir} -lctools
@@ -32,7 +32,7 @@ PROJECT_NAME           = ctools
 # This could be handy for archiving the generated documentation or
 # if some version control system is used.
-PROJECT_NUMBER         = 1.0.0
+PROJECT_NUMBER         = 0.9.0
 # Using the PROJECT_BRIEF tag one can provide an optional one line description
 # for a project that appears at the top of each page and should give viewer
@@ -66,10 +66,10 @@ copyright = u&#39;2014-2015, ctools developers&#39;
 # built documents.
 #
 # The short X.Y version.
-version = '1.0'
+version = '0.9'
 # The full version, including alpha/beta/rc tags.
-release = '1.0.0'
+release = '0.9.0'
 # The language for content autogenerated by Sphinx. Refer to documentation
 # for a list of supported languages.
@@ -13,7 +13,7 @@ bugs.
 Releases
 --------
-The latest ctools release is ``ctools-00-08-01`` (8 January 2015).
+The latest ctools release is ``ctools-0.9.0`` (22 May 2015).
 Below a table of ctools releases. Please note that at this stage of the
 project there is a strict link between the ctools and gammalib versions.
@@ -28,6 +28,9 @@ gammalib.
    * - `ctools <http://cta.irap.omp.eu/ctools/releases/ctools/ChangeLog>`_
      - `gammalib <http://cta.irap.omp.eu/ctools/releases/gammalib/ChangeLog>`_
      - Mac OS X package
+   * - `0.9.0 <http://cta.irap.omp.eu/ctools/releases/ctools/ctools-0.9.0.tar.gz>`_
+     - `0.10.0 <http://cta.irap.omp.eu/ctools/releases/gammalib/gammalib-0.10.0.tar.gz>`_
+     - `ctools-0.9.0-macosx10.3.dmg <http://cta.irap.omp.eu/ctools/releases/ctools/ctools-0.9.0-macosx10.3.dmg>`_
    * - `0.8.1 <http://cta.irap.omp.eu/ctools/releases/ctools/ctools-00-08-01.tar.gz>`_
      - `0.9.1 <http://cta.irap.omp.eu/ctools/releases/gammalib/gammalib-00-09-01.tar.gz>`_
      - `ctools-00-08-01-macosx10.3.dmg <http://cta.irap.omp.eu/ctools/releases/ctools/ctools-00-08-01-macosx10.3.dmg>`_
+.. _cslightcrv:
+
+cslightcrv
+==========
+
+Computes a lightcurve for a given source.
+
+
+Synopsis
+--------
+
+Computes a lightcurve for a given source.
+
+
+General parameters
+------------------
+
+``(inobs = NONE) [file]``
+    Input event list, counts cube or observation definition XML file.
+
+``inmodel [file]``
+    Input source model XML file.
+
+``srcname [string]``
+    Name of the source in the source model XML file which should be used
+    for sensitivity computation.
+
+``outfile [file]``
+    Name of the source spectrum output file.
+
+``expcube [file]``
+    Exposure cube file (only needed for stacked analysis).
+
+``psfcube [file]``
+    PSF cube file (only needed for stacked analysis).
+
+``bkgcube [file]``
+    Background cube file (only needed for stacked analysis).
+
+``caldb [string]``
+    Calibration database.
+ 	 	 
+``irf [string]``
+    Instrumental response function.
+
+``tmin [real]``
+    Lightcurve start time (in MJD).
+
+``tmax [real]``
+    Lightcurve stop time (in MJD).
+
+``tbins [integer]``
+    Number of time bins.
+
+``(tbinalg = GTI) <FILE|LIN|GTI> [string]``
+    Algorithm for defining time bins.
+
+``tbinfile [file]``
+    File defining the time binning.
+
+``binned [boolean]``
+    Use binned analysis in each time bin.
+
+``enumbins [integer]``
+    Number of energy bins per light curve bin (0=unbinned).
+
+``emin [real]``
+    Lower energy limit of events (in TeV).
+ 	 	 
+``emax [real]``
+    Upper energy limit of events (in TeV).
+ 	 	 
+``coordsys <CEL|GAL> [string]``
+    Coordinate system (CEL - celestial, GAL - galactic).
+ 	 	 
+``proj <AIT|AZP|CAR|MER|MOL|STG|TAN> [string]``
+    Projection method.
+
+``xref [real]``
+    Right Ascension / Galactic longitude of image centre (J2000, in degrees).
+ 	 	 
+``yref [real]``
+    Declination / Galactic latitude of image centre (J2000, in degrees).
+
+``nxpix [integer]``
+    Size of the Right Ascension / Galactic longitude axis (in pixels).
+ 	 	 
+``nypix [integer]``
+    Size of the Declination / Galactic latitude axis (in pixels).
+ 	 	 
+``binsz [real]``
+    Pixel size (in degrees/pixel).
+
+``(calc_ts = yes) [boolean]``
+    Compute TS value for each time bin?
+
+``(calc_ulim = yes) [boolean]``
+    Compute upper limit for each time bin?
+
+``(fix_srcs = yes) [boolean]``
+    Fix other sky model parameters?
+
+``(fix_bkg = no) [boolean]``
+    Fix background model parameters?
+
+
+Standard parameters
+-------------------
+
+``(chatter = 2) [integer]``
+    Verbosity of the executable:
+     chatter = 0: no information will be logged
+     
+     chatter = 1: only errors will be logged
+     
+     chatter = 2: errors and actions will be logged
+     
+     chatter = 3: report about the task execution
+     
+     chatter = 4: detailed report about the task execution
+ 	 	 
+``(clobber = yes) [boolean]``
+    Specifies whether an existing output counts cube should be overwritten.
+ 	 	 
+``(debug = no) [boolean]``
+    Enables debug mode. In debug mode the executable will dump any log file output to the console.
+ 	 	 
+``(mode = ql) [string]``
+    Mode of automatic parameters (default is "ql", i.e. "query and learn").
+
+``(logfile = cslightcrv.log) [filename]``
+    Log filename.
+
+
+Related tools
+-------------
+
+None
@@ -45,9 +45,9 @@ Here some usage examples:
 ``csobsdef rad=5``
       Creates observation definition file with a ROI radius of 5 deg.
-``csobsdef caldb=dummy irf=cta_dummy_irf``
-      Creates observation definition file using the "cta_dummy_irf" IRF in the
-      "dummy" calibration database.
+``csobsdef caldb=prod2 irf=South_50h``
+      Creates observation definition file using the "South_50h" IRF in the
+      "prod2" calibration database.
 General parameters
@@ -41,6 +41,9 @@ General parameters
 ``irf [string]``
     Instrumental response function.
+``(deadc = 0.95) [real]``
+    Average deadtime correction factor.
+
 ``(edisp = no) [boolean]``
     Apply energy dispersion to response computation.
@@ -62,12 +65,18 @@ General parameters
 ``enumbins [integer]``
     Number of energy bins (0=unbinned).
-``duration [real]``
-    Duration of observation (in seconds).
+``(tmin = 0.0) [real]``
+    Start time (in seconds).
-``(deadc = 0.95) [real]``
-    Average deadtime correction factor.
+``tmax [real]``
+    Stop time (in seconds).
+ 	 	 
+``npix [integer]``
+    Number of pixels for binned analysis.
+``binsz [real]``
+    Pixel size for binned analysis.
+
 ``(rad = 5.0) [real]``
     ROI radius (in degrees).
@@ -77,12 +86,6 @@ General parameters
 ``(offset = 1.5) [real]``
     Observation pattern offset (in degrees).
-``(npix = 200) [integer]``
-    Number of pixels for binned analysis.
- 	 	 
-``(binsz = 0.05) [real]``
-    Pixel size for binned analysis.
-
 Standard parameters
 -------------------
@@ -108,6 +111,9 @@ Standard parameters
 ``(mode = ql) [string]``
     Mode of automatic parameters (default is "ql", i.e. "query and learn").
+``(logfile = cspull.log) [string]``
+    Log filename.
+
 Related tools
 -------------
@@ -83,16 +83,16 @@ General parameters
     Pixel size (in degrees/pixel).
 ``(calc_ts = yes) [boolean]``
-    Compute TS for each spectral point.
+    Compute TS for each spectral point?
 ``(calc_ulim = yes) [boolean]``
-    Compute upper limit for each spectral point.
+    Compute upper limit for each spectral point?
 ``(fix_srcs = yes) [boolean]``
-    Fix other sky model parameters.
+    Fix other sky model parameters?
 ``(fix_bkg = no) [boolean]``
-    Fix background parameters.
+    Fix background model parameters?
 Standard parameters
@@ -16,30 +16,49 @@ Monte-Carlo simulations.
 General parameters
 ------------------
-``(outfile = ts.dat) [file]``
+``(inobs = NONE) [file]``
+    Input event list, counts cube or observation definition XML file.
+
+``inmodel [file]``
+    Input source model XML file.
+
+``srcname [string]``
+    Name of the source in the source model XML file which should be used
+    for sensitivity computation.
+
+``outfile [file]``
     ASCII file containing the TS distribution values.
- 	 	 
-``ntrials [integer]``
-    Number of Monte Carlo samples.
- 	 	 
+
+``expcube [file]``
+    Exposure cube file (only needed for stacked analysis).
+
+``psfcube [file]``
+    PSF cube file (only needed for stacked analysis).
+
+``bkgcube [file]``
+    Background cube file (only needed for stacked analysis).
+
 ``caldb [string]``
     Calibration database.
 ``irf [string]``
     Instrumental response function.
-``type <point|gauss|shell|disk> [string]``
-    Source model type.
- 	 	 
-``(index = -2.48) [real]``
-    Spectral index for source model.
- 	 	 
-``offset [real]``
-    Source offset angle (in degrees).
- 	 	 
-``bkg [string]``
-    Background model file function (none=power law for subarray E).
+``(deadc = 0.95) [real]``
+    Average deadtime correction factor.
+
+``(edisp = no) [boolean]``
+    Apply energy dispersion to response computation.
+
+``ntrials [integer]``
+    Number of Monte Carlo samples.
+
+``ra [real]``
+    Right Ascension of CTA pointing (J2000, in degrees).
+``dec [real]``
+    Declination of CTA pointing (J2000, in degrees).
+
 ``emin [real]``
     Lower energy limit (in TeV).
@@ -49,18 +68,27 @@ General parameters
 ``enumbins [integer]``
     Number of energy bins (0=unbinned).
-``duration [real]``
-    Effective exposure time (in seconds).
+``(tmin = 0.0) [real]``
+    Start time (in seconds).
-``(rad = 5.0) [real]``
-    Radius of ROI (in degrees).
+``tmax [real]``
+    Stop time (in seconds).
-``(npix = 200) [integer]``
+``npix [integer]``
     Number of pixels for binned analysis.
-``(binsz = 0.05) [real]``
+``binsz [real]``
     Pixel size for binned analysis.
+``(rad = 5.0) [real]``
+    Radius of ROI (in degrees).
+
+``(pattern = single) [string]``
+    Pattern for pointing simulation (single/four).
+
+``(offset = 1.5) [real]``
+    Observation pattern offset (in degrees).
+
 Standard parameters
 -------------------
@@ -86,6 +114,9 @@ Standard parameters
 ``(mode = ql) [string]``
     Mode of automatic parameters (default is "ql", i.e. "query and learn").
+``(logfile = cstsdist.log) [string]``
+    Log filename.
+
 Related tools
 -------------
@@ -38,6 +38,7 @@ cscripts
    cscaldb --- Lists available instrument response functions <cscaldb>
    csobsdef --- Generates observation definition file <csobsdef>
+   cslightcrv --- Computes lightcurve <cslightcrv>
    cspull --- Generates pull distribution <cspull>
    cssens --- Computes CTA sensitivity <cssens>
    csspec --- Computes spectral points <csspec>
@@ -26,7 +26,7 @@ Examples
 ``infile [file]``
     Describes an automatic (queried) file-type parameter with no default value.
-``irf = kb_E_50h_v3 <kb_E_50h_v3|kb_H_50h_v3> [string]``
+``irf = South_50h <North_50h|South_50h> [string]``
     Describes an automatic (queried) string-type parameter with two value
     options.
@@ -39,4 +39,4 @@ Examples
   To set hidden parameters at runtime, specify them explicitly at 
   the command line; for example:
-  ``gtobssim start_time=900``
+  ``ctobssim start_time=900``
@@ -14,5 +14,10 @@ The following C++ and Python examples can help get you started:
 Building your own cscripts
 --------------------------
+To be written ...
+
+
 Building your own ctools
 ------------------------
+
+To be written ...
@@ -4,8 +4,16 @@ Behind the scenes
 CTA data structure
 ------------------
+To be written ...
+
+
 Instrument response functions
 -----------------------------
+To be written ...
+
+
 Maximum likelihood fitting
---------------------------
 \ No newline at end of file
+--------------------------
+
+To be written ...
@@ -31,7 +31,7 @@ Installing gammalib
 ctools is built on top of gammalib, hence gammalib needs to be built,
 installed and configured before ctools can be installed. In case that
 ctools is already installed on your system you may skip reading this
-section and continue with Section TBD that explains how to install ctools.
+section and continue with the :ref:`sec_building_ctools` section.
 You will need the following to build the gammalib:
@@ -59,7 +59,7 @@ in an appropriate location (for example ``$HOME/builds``), and type
 .. code-block:: bash
-  $ tar xvfz gammalib-1.0.0.tar.gz
+  $ tar xvfz gammalib-0.10.0.tar.gz
 (the ``$`` symbol indicates the console prompt and is not part of the
 command that you should type in).
@@ -69,7 +69,7 @@ typing
 .. code-block:: bash
-  $ cd gammalib-1.0.0
+  $ cd gammalib-0.10.0
   $ ./configure
   $ make
   $ make check
@@ -117,6 +117,8 @@ If you use C shell or a variant then add the following to your
   source $GAMMALIB/bin/gammalib-init.csh
+.. _sec_building_ctools:
+
 Building ctools
 ---------------
@@ -125,7 +127,7 @@ an appropriate location (for example ``$HOME/builds``), and type
 .. code-block:: bash
-  $ tar xvfz ctools-1.0.0.tar.gz
+  $ tar xvfz ctools-0.9.0.tar.gz
 (the ``$`` symbol indicates the console prompt and is not part of the
 command that you should type in).
@@ -134,7 +136,7 @@ Step in the directory and build the ctools by typing
 .. code-block:: bash
-  $ cd ctools-1.0.0
+  $ cd ctools-0.9.0
   $ ./configure
   $ make
@@ -171,7 +173,7 @@ following output at the end of the unit testing:
   PASS: test_python.py
   make[4]: Nothing to be done for `all'.
   ============================================================================
-  Testsuite summary for ctools 1.0.0
+  Testsuite summary for ctools 0.9.0
   ============================================================================
   # TOTAL: 3
   # PASS:  3
@@ -189,26 +191,33 @@ For older automake version, you should see
   ***************
   * Test ctools *
   ***************
-  Test ctobssim: .... ok
+  Test ctobssim: .. ok
   Test ctskymap: .. ok
   Test ctbin: .... ok
-  Test ctmodel: ...... ok
-  Test ctselect: .. ok
-  Test ctlike: ........ ok
-  Test cttsmap: .. ok
   Test ctexpcube: .... ok
   Test ctpsfcube: .... ok
   Test ctbkgcube: ...... ok
+  Test ctmodel: ...... ok
+  Test ctselect: .. ok
+  Test ctlike: ...... ok
+  Test cttsmap: .. ok
   Test ctcubemask: .. ok
   Test ctbutterfly: .. ok
+  Test ctulimit: .. ok
   PASS: test_ctools.sh
   *****************
   * Test cscripts *
   *****************
+  Test csobsdef: .. ok
   Test cspull: .. ok
+  Test cssens: .. ok
   Test cstsdist: .. ok
   Test csresmap: .. ok
+  Test csspec (unbinned): .. ok
+  Test csspec (binned): .. ok
+  Test cslightcrv (unbinned): .. ok
+  Test cslightcrv (binned): .. ok
   PASS: test_cscripts.sh
   ***********************
@@ -227,6 +236,7 @@ For older automake version, you should see
   Test ctbkgcube functionality: ... ok
   Test ctcubemask functionality: ... ok
   Test ctbutterfly functionality: ... ok
+  Test ctulimit functionality: ... ok
   Test unbinned pipeline with FITS file saving: .... ok
   Test unbinned in-memory pipeline: .... ok
   PASS: test_python.py
@@ -307,9 +317,8 @@ GammaLib unit tests fail
 Some users have reported failure of a large fraction of the GammaLib unit
 tests after after typing ``make check``. In all cases, this was related to
-the absence of the directory where the shared ``libcfitsio`` library 
-resides in the library load path. To solve the issue, locate the directory
-where the shared ``libcfitsio`` library resides and then type
+not finding the ``libcfitsio`` library. To solve the issue, locate the 
+directory where the shared ``libcfitsio`` library resides and then type
 .. code-block:: bash
@@ -20,12 +20,12 @@ nebula, each offset by 0.5 deg from the nebula in opposite directions:
   RA of pointing (degrees) (0-360) [83.63] 
   Dec of pointing (degrees) (-90-90) [22.01] 21.51
   Radius of FOV (degrees) (0-180) [5.0] 
-  Start time (MET in s) (0) [0.0] 
-  End time (MET in s) (0) [1800.0] 
-  Lower energy limit (TeV) (0) [0.1] 
-  Upper energy limit (TeV) (0) [100.0] 
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
+  Start time (MET in s) [0.0] 
+  End time (MET in s) [1800.0] 
+  Lower energy limit (TeV) [0.1] 
+  Upper energy limit (TeV) [100.0] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Model [$CTOOLS/share/models/crab.xml] 
   Output event data file or observation definition file [events.fits] events1.fits
@@ -35,12 +35,12 @@ nebula, each offset by 0.5 deg from the nebula in opposite directions:
   RA of pointing (degrees) (0-360) [83.63] 
   Dec of pointing (degrees) (-90-90) [21.51] 22.51
   Radius of FOV (degrees) (0-180) [5.0] 
-  Start time (MET in s) (0) [0.0] 
-  End time (MET in s) (0) [1800.0] 
-  Lower energy limit (TeV) (0) [0.1] 
-  Upper energy limit (TeV) (0) [100.0] 
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
+  Start time (MET in s) [0.0] 
+  End time (MET in s) [1800.0] 
+  Lower energy limit (TeV) [0.1] 
+  Upper energy limit (TeV) [100.0] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Model [$CTOOLS/share/models/crab.xml] 
   Output event data file or observation definition file [events1.fits] events2.fits
@@ -92,11 +92,11 @@ Now we are ready to do a joint maximum likelihood analysis using
 .. code-block:: bash
   $ ctlike
-  Event list, counts cube or observation definition file [events.fits] obs.xml
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
+  Event list, counts cube or observation definition file [selected_events.fits] obs.xml
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Source model [$CTOOLS/share/models/crab.xml] 
-  Source model output file [crab_results.xml]
+  Source model output file [crab_results.xml] 
 Instead of providing an event list or a counts cube, we now provided the
 filename of the observation definition XML file (here ``obs.xml``) as input
@@ -106,102 +106,96 @@ joint analysis, as indicated in the log file:
 .. code-block:: xml
-  2015-02-04T20:33:26: +==============+
-  2015-02-04T20:33:26: | Observations |
-  2015-02-04T20:33:26: +==============+
-  2015-02-04T20:33:26: === GObservations ===
-  2015-02-04T20:33:26:  Number of observations ....: 2
-  2015-02-04T20:33:26:  Number of predicted events : 0
-  2015-02-04T20:33:26: === GCTAObservation ===
-  2015-02-04T20:33:26:  Name ......................: Crab
-  2015-02-04T20:33:26:  Identifier ................: 00001
-  2015-02-04T20:33:26:  Instrument ................: CTA
-  2015-02-04T20:33:26:  Event file ................: events1.fits
-  2015-02-04T20:33:26:  Event type ................: EventList
-  2015-02-04T20:33:26:  Statistics ................: Poisson
-  2015-02-04T20:33:26:  Ontime ....................: 1800 s
-  2015-02-04T20:33:26:  Livetime ..................: 1710 s
-  2015-02-04T20:33:26:  Deadtime correction .......: 0.95
-  2015-02-04T20:33:26:  User energy range .........: undefined
-  2015-02-04T20:33:26: === GCTAPointing ===
-  2015-02-04T20:33:26:  Pointing direction ........: (RA,Dec)=(83.63,21.51)
+  2015-05-22T20:17:07: +==============+
+  2015-05-22T20:17:07: | Observations |
+  2015-05-22T20:17:07: +==============+
+  2015-05-22T20:17:07: === GObservations ===
+  2015-05-22T20:17:07:  Number of observations ....: 2
+  2015-05-22T20:17:07:  Number of predicted events : 0
+  2015-05-22T20:17:07: === GCTAObservation ===
+  2015-05-22T20:17:07:  Name ......................: 
+  2015-05-22T20:17:07:  Identifier ................: 00001
+  2015-05-22T20:17:07:  Instrument ................: CTA
+  2015-05-22T20:17:07:  Event file ................: events1.fits
+  2015-05-22T20:17:07:  Event type ................: EventList
+  2015-05-22T20:17:07:  Statistics ................: Poisson
+  2015-05-22T20:17:07:  Ontime ....................: 1800 s
+  2015-05-22T20:17:07:  Livetime ..................: 1710 s
+  2015-05-22T20:17:07:  Deadtime correction .......: 0.95
+  2015-05-22T20:17:07:  User energy range .........: undefined
+  2015-05-22T20:17:07: === GCTAPointing ===
+  2015-05-22T20:17:07:  Pointing direction ........: (RA,Dec)=(83.63,21.51)
   ...
-  2015-02-04T20:33:26: === GCTAObservation ===
-  2015-02-04T20:33:26:  Name ......................: Crab
-  2015-02-04T20:33:26:  Identifier ................: 00002
-  2015-02-04T20:33:26:  Instrument ................: CTA
-  2015-02-04T20:33:26:  Event file ................: events2.fits
-  2015-02-04T20:33:26:  Event type ................: EventList
-  2015-02-04T20:33:26:  Statistics ................: Poisson
-  2015-02-04T20:33:26:  Ontime ....................: 1800 s
-  2015-02-04T20:33:26:  Livetime ..................: 1710 s
-  2015-02-04T20:33:26:  Deadtime correction .......: 0.95
-  2015-02-04T20:33:26:  User energy range .........: undefined
-  2015-02-04T20:33:26: === GCTAPointing ===
-  2015-02-04T20:33:26:  Pointing direction ........: (RA,Dec)=(83.63,22.51)
+  2015-05-22T20:17:07: === GCTAObservation ===
+  2015-05-22T20:17:07:  Name ......................: 
+  2015-05-22T20:17:07:  Identifier ................: 00002
+  2015-05-22T20:17:07:  Instrument ................: CTA
+  2015-05-22T20:17:07:  Event file ................: events2.fits
+  2015-05-22T20:17:07:  Event type ................: EventList
+  2015-05-22T20:17:07:  Statistics ................: Poisson
+  2015-05-22T20:17:07:  Ontime ....................: 1800 s
+  2015-05-22T20:17:07:  Livetime ..................: 1710 s
+  2015-05-22T20:17:07:  Deadtime correction .......: 0.95
+  2015-05-22T20:17:07:  User energy range .........: undefined
+  2015-05-22T20:17:07: === GCTAPointing ===
+  2015-05-22T20:17:07:  Pointing direction ........: (RA,Dec)=(83.63,22.51)
   ...
-  2015-02-04T20:33:26: +=================================+
-  2015-02-04T20:33:26: | Maximum likelihood optimisation |
-  2015-02-04T20:33:26: +=================================+
-  2015-02-04T20:33:26:  >Iteration   0: -logL=89099.788, Lambda=1.0e-03
-  2015-02-04T20:33:26:  >Iteration   1: -logL=89094.435, Lambda=1.0e-03, delta=5.353, max(|grad|)=-4.388872 [Normalization:7]
-  2015-02-04T20:33:26:  >Iteration   2: -logL=89094.420, Lambda=1.0e-04, delta=0.015, max(|grad|)=-0.116432 [Index:3]
-  2015-02-04T20:33:26:  >Iteration   3: -logL=89094.420, Lambda=1.0e-05, delta=0.000, max(|grad|)=-0.003989 [Index:3]
-  2015-02-04T20:33:26: 
-  2015-02-04T20:33:26: +=========================================+
-  2015-02-04T20:33:26: | Maximum likelihood optimization results |
-  2015-02-04T20:33:26: +=========================================+
-  2015-02-04T20:33:26: === GOptimizerLM ===
-  2015-02-04T20:33:26:  Optimized function value ..: 89094.420
-  2015-02-04T20:33:26:  Absolute precision ........: 0.005
-  2015-02-04T20:33:26:  Acceptable value decrease .: 2
-  2015-02-04T20:33:26:  Optimization status .......: converged
-  2015-02-04T20:33:26:  Number of parameters ......: 9
-  2015-02-04T20:33:26:  Number of free parameters .: 4
-  2015-02-04T20:33:26:  Number of iterations ......: 3
-  2015-02-04T20:33:26:  Lambda ....................: 1e-06
-  2015-02-04T20:33:26:  Maximum log likelihood ....: -89094.420
-  2015-02-04T20:33:26:  Observed events  (Nobs) ...: 12260.000
-  2015-02-04T20:33:26:  Predicted events (Npred) ..: 12260.000 (Nobs - Npred = 6.95764e-06)
-  2015-02-04T20:33:26: === GModels ===
-  2015-02-04T20:33:26:  Number of models ..........: 2
-  2015-02-04T20:33:26:  Number of parameters ......: 9
-  2015-02-04T20:33:26: === GModelSky ===
-  2015-02-04T20:33:26:  Name ......................: Crab
-  2015-02-04T20:33:26:  Instruments ...............: all
-  2015-02-04T20:33:26:  Instrument scale factors ..: unity
-  2015-02-04T20:33:26:  Observation identifiers ...: all
-  2015-02-04T20:33:26:  Model type ................: PointSource
-  2015-02-04T20:33:26:  Model components ..........: "SkyDirFunction" * "PowerLaw" * "Constant"
-  2015-02-04T20:33:26:  Number of parameters ......: 6
-  2015-02-04T20:33:26:  Number of spatial par's ...: 2
-  2015-02-04T20:33:26:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
-  2015-02-04T20:33:26:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
-  2015-02-04T20:33:26:  Number of spectral par's ..: 3
-  2015-02-04T20:33:26:   Prefactor ................: 6.06322e-16 +/- 1.44948e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
-  2015-02-04T20:33:26:   Index ....................: -2.49814 +/- 0.0178524 [-0,-5]  (free,scale=-1,gradient)
-  2015-02-04T20:33:26:   PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
-  2015-02-04T20:33:26:  Number of temporal par's ..: 1
-  2015-02-04T20:33:26:   Constant .................: 1 (relative value) (fixed,scale=1,gradient)
-  2015-02-04T20:33:26: === GCTAModelRadialAcceptance ===
-  2015-02-04T20:33:26:  Name ......................: Background
-  2015-02-04T20:33:26:  Instruments ...............: CTA
-  2015-02-04T20:33:26:  Instrument scale factors ..: unity
-  2015-02-04T20:33:26:  Observation identifiers ...: all
-  2015-02-04T20:33:26:  Model type ................: "Gaussian" * "FileFunction" * "Constant"
-  2015-02-04T20:33:26:  Number of parameters ......: 3
-  2015-02-04T20:33:26:  Number of radial par's ....: 1
-  2015-02-04T20:33:26:   Sigma ....................: 3.0356 +/- 0.0215259 [0.01,10] deg2 (free,scale=1,gradient)
-  2015-02-04T20:33:26:  Number of spectral par's ..: 1
-  2015-02-04T20:33:26:   Normalization ............: 0.99936 +/- 0.0122024 [0,1000]  (free,scale=1,gradient)
-  2015-02-04T20:33:26:  Number of temporal par's ..: 1
-  2015-02-04T20:33:26:   Constant .................: 1 (relative value) (fixed,scale=1,gradient)
-  2015-02-04T20:33:26: 
-  2015-02-04T20:33:26: +==============+
-  2015-02-04T20:33:26: | Save results |
-  2015-02-04T20:33:26: +==============+
-  2015-02-04T20:33:26: 
-  2015-02-04T20:33:26: Application "ctlike" terminated after 10 wall clock seconds, consuming 0.385432 seconds of CPU time.
+  2015-05-22T20:17:07: +=================================+
+  2015-05-22T20:17:07: | Maximum likelihood optimisation |
+  2015-05-22T20:17:07: +=================================+
+  2015-05-22T20:17:07:  >Iteration   0: -logL=201836.657, Lambda=1.0e-03
+  2015-05-22T20:17:07:  >Iteration   1: -logL=201829.438, Lambda=1.0e-03, delta=7.220, max(|grad|)=22.756507 [Index:7]
+  2015-05-22T20:17:07:  >Iteration   2: -logL=201829.432, Lambda=1.0e-04, delta=0.005, max(|grad|)=0.086395 [Index:3]
+  2015-05-22T20:17:07:  >Iteration   3: -logL=201829.432, Lambda=1.0e-05, delta=0.000, max(|grad|)=0.001287 [Index:3]
+  2015-05-22T20:17:07: 
+  2015-05-22T20:17:07: +=========================================+
+  2015-05-22T20:17:07: | Maximum likelihood optimization results |
+  2015-05-22T20:17:07: +=========================================+
+  2015-05-22T20:17:07: === GOptimizerLM ===
+  2015-05-22T20:17:07:  Optimized function value ..: 201829.432
+  2015-05-22T20:17:07:  Absolute precision ........: 0.005
+  2015-05-22T20:17:07:  Acceptable value decrease .: 2
+  2015-05-22T20:17:07:  Optimization status .......: converged
+  2015-05-22T20:17:07:  Number of parameters ......: 10
+  2015-05-22T20:17:07:  Number of free parameters .: 4
+  2015-05-22T20:17:07:  Number of iterations ......: 3
+  2015-05-22T20:17:07:  Lambda ....................: 1e-06
+  2015-05-22T20:17:07:  Maximum log likelihood ....: -201829.432
+  2015-05-22T20:17:07:  Observed events  (Nobs) ...: 31799.000
+  2015-05-22T20:17:07:  Predicted events (Npred) ..: 31799.000 (Nobs - Npred = 2.67544e-06)
+  2015-05-22T20:17:07: === GModels ===
+  2015-05-22T20:17:07:  Number of models ..........: 2
+  2015-05-22T20:17:07:  Number of parameters ......: 10
+  2015-05-22T20:17:07: === GModelSky ===
+  2015-05-22T20:17:07:  Name ......................: Crab
+  2015-05-22T20:17:07:  Instruments ...............: all
+  2015-05-22T20:17:07:  Instrument scale factors ..: unity
+  2015-05-22T20:17:07:  Observation identifiers ...: all
+  2015-05-22T20:17:07:  Model type ................: PointSource
+  2015-05-22T20:17:07:  Model components ..........: "SkyDirFunction" * "PowerLaw" * "Constant"
+  2015-05-22T20:17:07:  Number of parameters ......: 6
+  2015-05-22T20:17:07:  Number of spatial par's ...: 2
+  2015-05-22T20:17:07:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
+  2015-05-22T20:17:07:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
+  2015-05-22T20:17:07:  Number of spectral par's ..: 3
+  2015-05-22T20:17:07:   Prefactor ................: 5.8727e-16 +/- 8.04314e-18 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
+  2015-05-22T20:17:07:   Index ....................: -2.49825 +/- 0.0118418 [-0,-5]  (free,scale=-1,gradient)
+  2015-05-22T20:17:07:   PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T20:17:07:  Number of temporal par's ..: 1
+  2015-05-22T20:17:07:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+  2015-05-22T20:17:07: === GCTAModelIrfBackground ===
+  2015-05-22T20:17:07:  Name ......................: CTABackgroundModel
+  2015-05-22T20:17:07:  Instruments ...............: CTA
+  2015-05-22T20:17:07:  Instrument scale factors ..: unity
+  2015-05-22T20:17:07:  Observation identifiers ...: all
+  2015-05-22T20:17:07:  Model type ................: "PowerLaw" * "Constant"
+  2015-05-22T20:17:07:  Number of parameters ......: 4
+  2015-05-22T20:17:07:  Number of spectral par's ..: 3
+  2015-05-22T20:17:07:   Prefactor ................: 1.01608 +/- 0.0152751 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
+  2015-05-22T20:17:07:   Index ....................: 0.0226023 +/- 0.0086794 [-5,5]  (free,scale=1,gradient)
+  2015-05-22T20:17:07:   PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T20:17:07:  Number of temporal par's ..: 1
+  2015-05-22T20:17:07:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
 The log file indicates that the fit converged quickly, the spectral
 parameters of the Crab nebula have now been constrained using the events
@@ -3,33 +3,50 @@
 Connecting observations to dedicated instrument response functions
 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-.. note::
-
-   To execute the following examples on your machine you need to
-   install a specific calibration database. Download the file
-   :download:`caldb.tar.gz <caldb.tar.gz>` into a local working
-   directory, decompress and untar the folder, and set the calibration
-   database to this folder by typing:
-
-   .. code-block:: bash
-
-      $ export CALDB=$PWD/caldb
-
-   Note that the response functions in the calibration database
-   are all identical to the ``cta_dummy_irf`` response, yet there
-   a multiple database and response entries so that the usage of
-   specific response functions can be emulated.
-
 So far the same instrument response function has been used for all
 observations in our analysis examples.
 You may have recognised that we always specified the calibration database
-``dummy`` and the instrument response function ``cta_dummy_irf`` that
-are shipped with the ctools release.
+``prod2`` and the instrument response function ``South_50h`` that
+are shipped with ctools.
 When combining multiple observations in a joint analysis we may however
 need to specify dedicated instrument response functions for the
 different observations.
-This can be done by adding the calibration information to the observation
-definition XML file:
+
+Let's simulate two observations of the Crab nebula, one done with
+the northern array and the other with the southern array:
+
+.. code-block:: bash
+
+  $ ctobssim
+  RA of pointing (degrees) (0-360) [83.63] 
+  Dec of pointing (degrees) (-90-90) [22.01] 21.51
+  Radius of FOV (degrees) (0-180) [5.0] 
+  Start time (MET in s) [0.0] 
+  End time (MET in s) [1800.0] 
+  Lower energy limit (TeV) [0.1] 
+  Upper energy limit (TeV) [100.0] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] North_50h
+  Model [$CTOOLS/share/models/crab.xml] 
+  Output event data file or observation definition file [events.fits] events1.fits
+
+.. code-block:: bash
+
+  $ ctobssim
+  RA of pointing (degrees) (0-360) [83.63] 
+  Dec of pointing (degrees) (-90-90) [22.01] 22.51
+  Radius of FOV (degrees) (0-180) [5.0] 
+  Start time (MET in s) [0.0] 
+  End time (MET in s) [1800.0] 
+  Lower energy limit (TeV) [0.1] 
+  Upper energy limit (TeV) [100.0] 
+  Calibration database [prod2] 
+  Instrument response function [North_50h] South_50h
+  Model [$CTOOLS/share/models/crab.xml] 
+  Output event data file or observation definition file [events1.fits] events2.fits
+
+Now we combine the information about these two observations in an
+observation definition XML file:
 .. code-block:: xml
@@ -37,11 +54,11 @@ definition XML file:
   <observation_list title="observation library">
     <observation name="Crab" id="00001" instrument="CTA">
       <parameter name="EventList" file="events1.fits"/>
-      <parameter name="Calibration" database="north" response="cta_north_irf"/>
+      <parameter name="Calibration" database="prod2" response="North_50h"/>
     </observation>
     <observation name="Crab" id="00002" instrument="CTA">
       <parameter name="EventList" file="events2.fits"/>
-      <parameter name="Calibration" database="south" response="cta_south_irf"/>
+      <parameter name="Calibration" database="prod2" response="South_50h"/>
     </observation>
   </observation_list>
@@ -55,41 +72,39 @@ option to get the output of the script on the console):
    $ cscaldb debug=yes
    Parfile cscaldb.par not found. Create default parfile.
-   2015-02-05T14:26:50: +============+
-   2015-02-05T14:26:50: | Parameters |
-   2015-02-05T14:26:50: +============+
-   2015-02-05T14:26:50:  chatter ...................: 2
-   2015-02-05T14:26:50:  clobber ...................: yes
-   2015-02-05T14:26:50:  debug .....................: yes
-   2015-02-05T14:26:50:  mode ......................: ql
-   2015-02-05T14:26:50: 
-   2015-02-05T14:26:50: +==============+
-   2015-02-05T14:26:50: | Mission: cta |
-   2015-02-05T14:26:50: +==============+
-   2015-02-05T14:26:50: === Response functions in database "dummy" ===
-   2015-02-05T14:26:50: cta_dummy_irf
-   2015-02-05T14:26:50: cta_prod1_irf
-   2015-02-05T14:26:50: cta_prod2_irf
-   2015-02-05T14:26:50: 
-   2015-02-05T14:26:50: === Response functions in database "north" ===
-   2015-02-05T14:26:50: cta_north_irf
-   2015-02-05T14:26:50: 
-   2015-02-05T14:26:50: === Response functions in database "south" ===
-   2015-02-05T14:26:50: cta_south_irf
-
-You will see that for the CTA observatory three databases exist.
-The ``dummy`` database contains three response functions
-(``cta_dummy_irf``, ``cta_prod1_irf``, ``cta_prod2_irf``) while the
-``north`` and ``south`` databases contain a single response function.
-You may now do a :ref:`ctlike` fit using the XML file with the specific
-response definitions:
+   2015-05-22T21:35:39: +============+
+   2015-05-22T21:35:39: | Parameters |
+   2015-05-22T21:35:39: +============+
+   2015-05-22T21:35:39:  chatter ...................: 2
+   2015-05-22T21:35:39:  clobber ...................: yes
+   2015-05-22T21:35:39:  debug .....................: yes
+   2015-05-22T21:35:39:  mode ......................: ql
+   2015-05-22T21:35:39: 
+   2015-05-22T21:35:39: +==============+
+   2015-05-22T21:35:39: | Mission: cta |
+   2015-05-22T21:35:39: +==============+
+   2015-05-22T21:35:40: === Response functions in database "prod2" ===
+   2015-05-22T21:35:40: North_0.5h
+   2015-05-22T21:35:40: North_50h
+   2015-05-22T21:35:40: North_5h
+   2015-05-22T21:35:40: South_0.5h
+   2015-05-22T21:35:40: South_50h
+   2015-05-22T21:35:40: South_5h
+
+You will see that ctools ships with one database for the CTA observatory.
+This is the ``prod2`` database.
+Within this database there are six response functions:
+``North_0.5h``, ``North_5h``, ``North_50h``,
+``South_0.5h``, ``South_5h``, and ``South_50h``.
+
+We now jointly fit both observations using :ref:`ctlike`:
 .. code-block:: bash
    $ ctlike
-   Event list, counts cube or observation definition file [events.fits] obs_irf.xml
+   Event list, counts cube or observation definition file [obs.xml] obs_irf.xml
    Source model [$CTOOLS/share/models/crab.xml] 
-   Source model output file [crab_results.xml]
+   Source model output file [crab_results.xml] 
 To see the usage of the individual response functions you may inspect the
 log file. You will notice that each observation now has a specific response
@@ -98,69 +113,110 @@ both observations.
 .. code-block:: xml
-   2015-02-05T14:48:46: +==============+
-   2015-02-05T14:48:46: | Observations |
-   2015-02-05T14:48:46: +==============+
-   2015-02-05T14:48:46: === GObservations ===
-   2015-02-05T14:48:46:  Number of observations ....: 2
-   2015-02-05T14:48:46:  Number of predicted events : 0
-   2015-02-05T14:48:46: === GCTAObservation ===
-   2015-02-05T14:48:46:  Name ......................: Crab
-   2015-02-05T14:48:46:  Identifier ................: 00001
-   2015-02-05T14:48:46:  Instrument ................: CTA
-   2015-02-05T14:48:46:  Event file ................: events1.fits
-   2015-02-05T14:48:46:  Event type ................: EventList
-   2015-02-05T14:48:46:  Statistics ................: Poisson
-   2015-02-05T14:48:46:  Ontime ....................: 1800 s
-   2015-02-05T14:48:46:  Livetime ..................: 1710 s
-   2015-02-05T14:48:46:  Deadtime correction .......: 0.95
-   2015-02-05T14:48:46:  User energy range .........: undefined
-   2015-02-05T14:48:46: === GCTAPointing ===
-   2015-02-05T14:48:46:  Pointing direction ........: (RA,Dec)=(83.63,21.51)
-   2015-02-05T14:48:46: === GCTAResponseIrf ===
-   2015-02-05T14:48:46:  Response name .............: cta_north_irf
-   2015-02-05T14:48:46:  Energy dispersion .........: Not used
-   2015-02-05T14:48:46:  Save energy range .........: undefined
-   2015-02-05T14:48:46: === GCaldb ===
-   2015-02-05T14:48:46:  Database root .............: /Users/jurgen/Desktop/tmp/beyond/caldb
-   2015-02-05T14:48:46:  Selected Mission ..........: CTA
-   2015-02-05T14:48:46:  Selected Instrument .......: NORTH
-   2015-02-05T14:48:46:  Calibration Index File ....: /Users/jurgen/Desktop/tmp/beyond/caldb/data/cta/north/caldb.indx
-   2015-02-05T14:48:46:  Number of entries .........: 4
-   2015-02-05T14:48:46: === GCTAAeffPerfTable ===
-   2015-02-05T14:48:46:  Filename ..................: /Users/jurgen/Desktop/tmp/beyond/caldb/data/cta/north/bcf/cta_north_irf.dat
-   2015-02-05T14:48:46:  Number of energy bins .....: 20
-   2015-02-05T14:48:46:  Log10(Energy) range .......: 0.0199526 - 125.893 TeV
-   2015-02-05T14:48:46:  Offset angle dependence ...: Fixed sigma=3
-   ...
-   2015-02-05T14:48:46: === GCTAObservation ===
-   2015-02-05T14:48:46:  Name ......................: Crab
-   2015-02-05T14:48:46:  Identifier ................: 00002
-   2015-02-05T14:48:46:  Instrument ................: CTA
-   2015-02-05T14:48:46:  Event file ................: events2.fits
-   2015-02-05T14:48:46:  Event type ................: EventList
-   2015-02-05T14:48:46:  Statistics ................: Poisson
-   2015-02-05T14:48:46:  Ontime ....................: 1800 s
-   2015-02-05T14:48:46:  Livetime ..................: 1710 s
-   2015-02-05T14:48:46:  Deadtime correction .......: 0.95
-   2015-02-05T14:48:46:  User energy range .........: undefined
-   2015-02-05T14:48:46: === GCTAPointing ===
-   2015-02-05T14:48:46:  Pointing direction ........: (RA,Dec)=(83.63,22.51)
-   2015-02-05T14:48:46: === GCTAResponseIrf ===
-   2015-02-05T14:48:46:  Response name .............: cta_south_irf
-   2015-02-05T14:48:46:  Energy dispersion .........: Not used
-   2015-02-05T14:48:46:  Save energy range .........: undefined
-   2015-02-05T14:48:46: === GCaldb ===
-   2015-02-05T14:48:46:  Database root .............: /Users/jurgen/Desktop/tmp/beyond/caldb
-   2015-02-05T14:48:46:  Selected Mission ..........: CTA
-   2015-02-05T14:48:46:  Selected Instrument .......: SOUTH
-   2015-02-05T14:48:46:  Calibration Index File ....: /Users/jurgen/Desktop/tmp/beyond/caldb/data/cta/south/caldb.indx
-   2015-02-05T14:48:46:  Number of entries .........: 4
-   2015-02-05T14:48:46: === GCTAAeffPerfTable ===
-   2015-02-05T14:48:46:  Filename ..................: /Users/jurgen/Desktop/tmp/beyond/caldb/data/cta/south/bcf/cta_south_irf.dat
-   2015-02-05T14:48:46:  Number of energy bins .....: 20
-   2015-02-05T14:48:46:  Log10(Energy) range .......: 0.0199526 - 125.893 TeV
-   2015-02-05T14:48:46:  Offset angle dependence ...: Fixed sigma=3
+  2015-05-22T22:06:43: +==============+
+  2015-05-22T22:06:43: | Observations |
+  2015-05-22T22:06:43: +==============+
+  2015-05-22T22:06:43: === GObservations ===
+  2015-05-22T22:06:43:  Number of observations ....: 2
+  2015-05-22T22:06:43:  Number of predicted events : 0
+  2015-05-22T22:06:43: === GCTAObservation ===
+  2015-05-22T22:06:43:  Name ......................: 
+  2015-05-22T22:06:43:  Identifier ................: 00001
+  2015-05-22T22:06:43:  Instrument ................: CTA
+  2015-05-22T22:06:43:  Event file ................: events1.fits
+  2015-05-22T22:06:43:  Event type ................: EventList
+  2015-05-22T22:06:43:  Statistics ................: Poisson
+  2015-05-22T22:06:43:  Ontime ....................: 1800 s
+  2015-05-22T22:06:43:  Livetime ..................: 1710 s
+  2015-05-22T22:06:43:  Deadtime correction .......: 0.95
+  2015-05-22T22:06:43:  User energy range .........: undefined
+  2015-05-22T22:06:43: === GCTAPointing ===
+  2015-05-22T22:06:43:  Pointing direction ........: (RA,Dec)=(83.63,21.51)
+  2015-05-22T22:06:43: === GCTAResponseIrf ===
+  2015-05-22T22:06:43:  Caldb mission .............: cta
+  2015-05-22T22:06:43:  Caldb instrument ..........: prod2
+  2015-05-22T22:06:43:  Response name .............: North_50h
+  2015-05-22T22:06:43:  Energy dispersion .........: Not used
+  2015-05-22T22:06:43:  Save energy range .........: undefined
+  ...
+  2015-05-22T22:06:43: === GCTAObservation ===
+  2015-05-22T22:06:43:  Name ......................: 
+  2015-05-22T22:06:43:  Identifier ................: 00002
+  2015-05-22T22:06:43:  Instrument ................: CTA
+  2015-05-22T22:06:43:  Event file ................: events2.fits
+  2015-05-22T22:06:43:  Event type ................: EventList
+  2015-05-22T22:06:43:  Statistics ................: Poisson
+  2015-05-22T22:06:43:  Ontime ....................: 1800 s
+  2015-05-22T22:06:43:  Livetime ..................: 1710 s
+  2015-05-22T22:06:43:  Deadtime correction .......: 0.95
+  2015-05-22T22:06:43:  User energy range .........: undefined
+  2015-05-22T22:06:43: === GCTAPointing ===
+  2015-05-22T22:06:43:  Pointing direction ........: (RA,Dec)=(83.63,22.51)
+  2015-05-22T22:06:43: === GCTAResponseIrf ===
+  2015-05-22T22:06:43:  Caldb mission .............: cta
+  2015-05-22T22:06:43:  Caldb instrument ..........: prod2
+  2015-05-22T22:06:43:  Response name .............: South_50h
+  2015-05-22T22:06:43:  Energy dispersion .........: Not used
+  2015-05-22T22:06:43:  Save energy range .........: undefined
+  ...
+  2015-05-22T22:06:43: +=================================+
+  2015-05-22T22:06:43: | Maximum likelihood optimisation |
+  2015-05-22T22:06:43: +=================================+
+  2015-05-22T22:06:43:  >Iteration   0: -logL=138039.472, Lambda=1.0e-03
+  2015-05-22T22:06:43:  >Iteration   1: -logL=138035.501, Lambda=1.0e-03, delta=3.971, max(|grad|)=16.203761 [Index:7]
+  2015-05-22T22:06:43:  >Iteration   2: -logL=138035.496, Lambda=1.0e-04, delta=0.005, max(|grad|)=-0.053284 [Index:7]
+  2015-05-22T22:06:43:  >Iteration   3: -logL=138035.496, Lambda=1.0e-05, delta=0.000, max(|grad|)=0.001386 [Index:7]
+  2015-05-22T22:06:44: 
+  2015-05-22T22:06:44: +=========================================+
+  2015-05-22T22:06:44: | Maximum likelihood optimization results |
+  2015-05-22T22:06:44: +=========================================+
+  2015-05-22T22:06:44: === GOptimizerLM ===
+  2015-05-22T22:06:44:  Optimized function value ..: 138035.496
+  2015-05-22T22:06:44:  Absolute precision ........: 0.005
+  2015-05-22T22:06:44:  Acceptable value decrease .: 2
+  2015-05-22T22:06:44:  Optimization status .......: converged
+  2015-05-22T22:06:44:  Number of parameters ......: 10
+  2015-05-22T22:06:44:  Number of free parameters .: 4
+  2015-05-22T22:06:44:  Number of iterations ......: 3
+  2015-05-22T22:06:44:  Lambda ....................: 1e-06
+  2015-05-22T22:06:44:  Maximum log likelihood ....: -138035.496
+  2015-05-22T22:06:44:  Observed events  (Nobs) ...: 20814.000
+  2015-05-22T22:06:44:  Predicted events (Npred) ..: 20814.000 (Nobs - Npred = 2.67698e-06)
+  2015-05-22T22:06:44: === GModels ===
+  2015-05-22T22:06:44:  Number of models ..........: 2
+  2015-05-22T22:06:44:  Number of parameters ......: 10
+  2015-05-22T22:06:44: === GModelSky ===
+  2015-05-22T22:06:44:  Name ......................: Crab
+  2015-05-22T22:06:44:  Instruments ...............: all
+  2015-05-22T22:06:44:  Instrument scale factors ..: unity
+  2015-05-22T22:06:44:  Observation identifiers ...: all
+  2015-05-22T22:06:44:  Model type ................: PointSource
+  2015-05-22T22:06:44:  Model components ..........: "SkyDirFunction" * "PowerLaw" * "Constant"
+  2015-05-22T22:06:44:  Number of parameters ......: 6
+  2015-05-22T22:06:44:  Number of spatial par's ...: 2
+  2015-05-22T22:06:44:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
+  2015-05-22T22:06:44:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
+  2015-05-22T22:06:44:  Number of spectral par's ..: 3
+  2015-05-22T22:06:44:   Prefactor ................: 5.9184e-16 +/- 9.68924e-18 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
+  2015-05-22T22:06:44:   Index ....................: -2.49569 +/- 0.0145376 [-0,-5]  (free,scale=-1,gradient)
+  2015-05-22T22:06:44:   PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T22:06:44:  Number of temporal par's ..: 1
+  2015-05-22T22:06:44:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+  2015-05-22T22:06:44: === GCTAModelIrfBackground ===
+  2015-05-22T22:06:44:  Name ......................: CTABackgroundModel
+  2015-05-22T22:06:44:  Instruments ...............: CTA
+  2015-05-22T22:06:44:  Instrument scale factors ..: unity
+  2015-05-22T22:06:44:  Observation identifiers ...: all
+  2015-05-22T22:06:44:  Model type ................: "PowerLaw" * "Constant"
+  2015-05-22T22:06:44:  Number of parameters ......: 4
+  2015-05-22T22:06:44:  Number of spectral par's ..: 3
+  2015-05-22T22:06:44:   Prefactor ................: 1.02524 +/- 0.0182054 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
+  2015-05-22T22:06:44:   Index ....................: 0.0229728 +/- 0.0104713 [-5,5]  (free,scale=1,gradient)
+  2015-05-22T22:06:44:   PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T22:06:44:  Number of temporal par's ..: 1
+  2015-05-22T22:06:44:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+
+
 A more fine grained control over the response function can be achieved by 
 specifying individual filenames for the various response components.
@@ -173,17 +229,17 @@ This is definitely expert mode, to be used with utmost care.
   <observation_list title="observation library">
     <observation name="Crab" id="00001" instrument="CTA">
       <parameter name="EventList"           file="events1.fits"/>
-      <parameter name="EffectiveArea"       file="$(CALDB)/data/cta/dummy/bcf/cta_dummy_irf.dat"/>
-      <parameter name="PointSpreadFunction" file="$(CALDB)/data/cta/dummy/bcf/cta_prod1_irf.dat"/>
-      <parameter name="EnergyDispersion"    file="$(CALDB)/data/cta/dummy/bcf/cta_prod2_irf.dat"/>
-      <parameter name="Background"          file="$(CALDB)/data/cta/north/bcf/cta_north_irf.dat"/>
+      <parameter name="EffectiveArea"       file="$(CALDB)/data/cta/prod2/bcf/North_50h/irf_file.fits.gz"/>
+      <parameter name="PointSpreadFunction" file="$(CALDB)/data/cta/prod2/bcf/North_50h/irf_file.fits.gz"/>
+      <parameter name="EnergyDispersion"    file="$(CALDB)/data/cta/prod2/bcf/North_50h/irf_file.fits.gz"/>
+      <parameter name="Background"          file="$(CALDB)/data/cta/prod2/bcf/North_50h/irf_file.fits.gz"/>
     </observation>
     <observation name="Crab" id="00002" instrument="CTA">
       <parameter name="EventList"           file="events2.fits"/>
-      <parameter name="EffectiveArea"       file="$(CALDB)/data/cta/dummy/bcf/cta_dummy_irf.dat"/>
-      <parameter name="PointSpreadFunction" file="$(CALDB)/data/cta/dummy/bcf/cta_prod1_irf.dat"/>
-      <parameter name="EnergyDispersion"    file="$(CALDB)/data/cta/dummy/bcf/cta_prod2_irf.dat"/>
-      <parameter name="Background"          file="$(CALDB)/data/cta/south/bcf/cta_south_irf.dat"/>
+      <parameter name="EffectiveArea"       file="$(CALDB)/data/cta/prod2/bcf/South_50h/irf_file.fits.gz"/>
+      <parameter name="PointSpreadFunction" file="$(CALDB)/data/cta/prod2/bcf/South_50h/irf_file.fits.gz"/>
+      <parameter name="EnergyDispersion"    file="$(CALDB)/data/cta/prod2/bcf/South_50h/irf_file.fits.gz"/>
+      <parameter name="Background"          file="$(CALDB)/data/cta/prod2/bcf/South_50h/irf_file.fits.gz"/>
     </observation>
   </observation_list>
@@ -199,11 +255,13 @@ file is given below:
       <parameter name="CountsCube"   file="cntcube1.fits"/>
       <parameter name="ExposureCube" file="expcube1.fits"/>
       <parameter name="PsfCube"      file="psfcube1.fits"/>
+      <parameter name="BkgCube"      file="bkgcube1.fits"/>
     </observation>
     <observation name="Crab" id="00002" instrument="CTA">
       <parameter name="CountsCube"   file="cntcube2.fits"/>
       <parameter name="ExposureCube" file="expcube2.fits"/>
       <parameter name="PsfCube"      file="psfcube2.fits"/>
+      <parameter name="BkgCube"      file="bkgcube2.fits"/>
     </observation>
   </observation_list>
@@ -56,23 +56,21 @@ event list.
         <parameter name="DEC" scale="1.0" value="22.0145" min="-90"  max="90"  free="0"/>
       </spatialModel>
     </source>
-    <source name="Background_00001" type="RadialAcceptance" instrument="CTA" id="00001">
-      <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-        <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="1.0" free="1"/>
+    <source name="Background_00001" type="CTAIrfBackground" instrument="CTA" id="00001">
+      <spectrum type="PowerLaw">
+        <parameter name="Prefactor" scale="1.0"  value="1.0"  min="1e-3" max="1e+3"   free="1"/>
+        <parameter name="Index"     scale="1.0"  value="0.0"  min="-5.0" max="+5.0"   free="1"/>
+        <parameter name="Scale"     scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>
       </spectrum>
-      <radialModel type="Gaussian">
-         <parameter name="Sigma" scale="1.0" value="3.0" min="0.01" max="10.0" free="1"/>
-      </radialModel>
     </source>
-    <source name="Background_00002" type="RadialAcceptance" instrument="CTA" id="00002">
-      <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-        <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="10.0" free="1"/>
+    <source name="Background_00002" type="CTAIrfBackground" instrument="CTA" id="00002">
+      <spectrum type="PowerLaw">
+        <parameter name="Prefactor" scale="1.0"  value="10.0"  min="1e-3" max="1e+3"   free="1"/>
+        <parameter name="Index"     scale="1.0"  value="0.0"  min="-5.0" max="+5.0"   free="1"/>
+        <parameter name="Scale"     scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>
       </spectrum>
-      <radialModel type="Gaussian">
-         <parameter name="Sigma" scale="1.0" value="3.0" min="0.01" max="10.0" free="1"/>
-      </radialModel>
     </source>
-    </source_library>
+  </source_library>
 Note that both background model components distinguish in their event rate;
 component ``Background_00002`` has a ten times larger background rate than
@@ -112,21 +110,19 @@ This is illustrated in the model definition file below where component
         <parameter name="DEC" scale="1.0" value="22.0145" min="-90"  max="90"  free="0"/>
       </spatialModel>
     </source>
-    <source name="Background_00001" type="RadialAcceptance" instrument="CTA" id="00001,00003,00004">
-      <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-        <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="1.0" free="1"/>
+    <source name="Background_00001" type="CTAIrfBackground" instrument="CTA" id="00001,00003,00004">
+      <spectrum type="PowerLaw">
+        <parameter name="Prefactor" scale="1.0"  value="1.0"  min="1e-3" max="1e+3"   free="1"/>
+        <parameter name="Index"     scale="1.0"  value="0.0"  min="-5.0" max="+5.0"   free="1"/>
+        <parameter name="Scale"     scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>
       </spectrum>
-      <radialModel type="Gaussian">
-         <parameter name="Sigma" scale="1.0" value="3.0" min="0.01" max="10.0" free="1"/>
-      </radialModel>
     </source>
-    <source name="Background_00002" type="RadialAcceptance" instrument="CTA" id="00002">
-      <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-        <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="10.0" free="1"/>
+    <source name="Background_00002" type="CTAIrfBackground" instrument="CTA" id="00002">
+      <spectrum type="PowerLaw">
+        <parameter name="Prefactor" scale="1.0"  value="10.0"  min="1e-3" max="1e+3"   free="1"/>
+        <parameter name="Index"     scale="1.0"  value="0.0"  min="-5.0" max="+5.0"   free="1"/>
+        <parameter name="Scale"     scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>
       </spectrum>
-      <radialModel type="Gaussian">
-         <parameter name="Sigma" scale="1.0" value="3.0" min="0.01" max="10.0" free="1"/>
-      </radialModel>
     </source>
   </source_library>
@@ -16,15 +16,27 @@ the model value to obtained the fractional residuals
    Residual = \frac{Counts-Model}{Model}
 for all map pixels.
-The :ref:`csresmap` script is used as follows:
+
+Suppose that you have done the following unbinned ``ctlike`` analysis:
+
+.. code-block:: bash
+
+  $ ctlike
+  Event list, counts cube or observation definition file [events.fits] obs.xml
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
+  Source model [$CTOOLS/share/models/crab.xml] 
+  Source model output file [crab_results.xml] 
+
+The :ref:`csresmap` script is then used as follows:
 .. code-block:: bash
   $ csresmap
   Parfile csresmap.par not found. Create default parfile.
   Event list, counts cube, or observation definition file [events.fits] obs.xml
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Source model [$CTOOLS/share/models/crab.xml] crab_results.xml
   Output residual map [resmap.fits] 
   First coordinate of image center in degrees (RA or galactic l) [83.63] 
@@ -15,10 +15,10 @@ into a single counts cube:
 .. code-block:: bash
   $ ctbin
-  Event list or observation definition file [events2.fits] obs.xml
-  First coordinate of image center in degrees (RA or galactic l) [83.63] 
-  Second coordinate of image center in degrees (DEC or galactic b) [22.01] 
-  Projection method e.g. AIT|AZP|CAR|MER|STG|TAN (AIT|AZP|CAR|MER|STG|TAN) [CAR] 
+  Event list or observation definition file [events.fits] obs.xml
+  First coordinate of image center in degrees (RA or galactic l) (0-360) [83.63] 
+  Second coordinate of image center in degrees (DEC or galactic b) (-90-90) [22.01] 
+  Projection method (AIT|AZP|CAR|MER|MOL|STG|TAN) [CAR] 
   Coordinate system (CEL - celestial, GAL - galactic) (CEL|GAL) [CEL] 
   Image scale (in degrees/pixel) [0.02] 
   Size of the X axis in pixels [200] 
@@ -27,15 +27,16 @@ into a single counts cube:
   Start value for first energy bin in TeV [0.1] 
   Stop value for last energy bin in TeV [100.0] 
   Number of energy bins [20] 
-  Output counts cube [cntcube2.fits] cntcube.fits
+  Output counts cube [cntcube.fits] 
 We now have a stacked counts cube ``cntcube.fits`` on disk.
 Before we can use that counts cube in a maximum likelihood
 analysis, we have to compute the instrument response and the
 background model that are needed to describe the stacked data.
+
 For the former, we have to compute the total exposure for the stacked
 cube (i.e. the sum of the effective areas for each observation multiplied
-by the corresponding lifetimes) and an effective point spread function
+by the corresponding livetimes) and an effective point spread function
 (i.e. the point spread function of the different observations weighted by
 the corresponding exposures).
 To get both informations we use the :ref:`ctexpcube` and 
@@ -45,29 +46,29 @@ To get both informations we use the :ref:`ctexpcube` and
   $ ctexpcube
   Event list or observation definition file [NONE] obs.xml
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Counts cube for exposure cube definition [NONE] cntcube.fits
-  Output exposure cube file [expcube.fits]
+  Output exposure cube file [expcube.fits] 
 .. code-block:: bash
   $ ctpsfcube
   Event list or observation definition file [NONE] obs.xml
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Counts cube for psf cube definition [NONE] 
-  First coordinate of image center in degrees (RA or galactic l) [83.63] 
-  Second coordinate of image center in degrees (DEC or galactic b) [22.01] 
-  Projection method e.g. AIT|AZP|CAR|MER|MOL|STG|TAN (AIT|AZP|CAR|MER|MOL|STG|TAN) [CAR] 
+  First coordinate of image center in degrees (RA or galactic l) (0-360) [83.63] 
+  Second coordinate of image center in degrees (DEC or galactic b) (-90-90) [22.01] 
+  Projection method (AIT|AZP|CAR|MER|MOL|STG|TAN) [CAR] 
   Coordinate system (CEL - celestial, GAL - galactic) (CEL|GAL) [CEL] 
   Image scale (in degrees/pixel) [1.0] 
   Size of the X axis in pixels [10] 
   Size of the Y axis in pixels [10] 
-  Start value for first energy bin in TeV [0.1] 
-  Stop value for last energy bin in TeV [100.0] 
+  Lower energy limit (TeV) [0.1] 
+  Upper energy limit (TeV) [100.0] 
   Number of energy bins [20] 
-  Output psf cube file [psfcube.fits]
+  Output psf cube file [psfcube.fits] 
 We provide the ``obs.xml`` file on input to inform both tools which
 observations have been combined.
@@ -93,17 +94,17 @@ background cube using the :ref:`ctbkgcube` tool:
   $ ctbkgcube
   Input event list or observation definition file [NONE] obs.xml
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
-  Input (background) model XML file [NONE] $CTOOLS/share/models/crab.xml
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Counts cube for background cube definition [NONE] cntcube.fits
+  Input model XML file [NONE] $CTOOLS/share/models/crab.xml
   Output background cube file [bkgcube.fits] 
-  Output (background) model XML file [NONE] model.xml
+  Output model XML file [NONE] model.xml
 The usage of :ref:`ctbkgcube` is very similar to that of :ref:`ctexpcube`,
 yet it takes the model XML file as an additional input parameter.
-We here use the usual ``$CTOOLS/share/models/crab.xml`` Crab plus
-background model file that is shipped with the ctools.
+We here use the usual ``$CTOOLS/share/models/crab.xml`` model
+file that is shipped with the ctools.
 :ref:`ctbkgcube` provides on output the background cube file
 ``bkgcube.fits`` and the model XML file ``model.xml`` that can
 be used for further analysis.
@@ -125,25 +126,21 @@ modelling works:
         <parameter name="DEC" value="22.0145" scale="1" min="-90" max="90" free="0" />
       </spatialModel>
     </source>
-    <source name="ctbkgcube default background model" type="CTACubeBackground">
+    <source name="BackgroundModel" type="CTACubeBackground" instrument="CTA,HESS,MAGIC,VERITAS">
       <spectrum type="PowerLaw">
         <parameter name="Prefactor" value="1" error="0" scale="1" min="0" free="1" />
         <parameter name="Index" value="0" error="0" scale="1" min="-10" max="10" free="1" />
         <parameter name="Scale" value="1" scale="1e+06" free="0" />
       </spectrum>
-      <spatialModel type="MapCubeFunction" file="bkgcube.fits">
-        <parameter name="Normalization" value="1" scale="1" min="0.001" max="1000" free="0" />
-      </spatialModel>
     </source>
   </source_library>
 The Crab source component is the same that is also present in
 ``$CTOOLS/share/models/crab.xml`` and is not modified.
-The background component, however, has been replaced and now is
-the ``ctbkgcube default background model``.
-This model is of type ``CTACubeBackground`` which is a 3-dimensional
-data cube that describes the expected background rate as function
-of spatial position and energy.
+The background component, however, has been replaced by a model of
+type ``CTACubeBackground``.
+This model is a 3-dimensional data cube that describes the expected 
+background rate as function of spatial position and energy.
 The data cube is multiplied by a power law spectrum that allows to adjust
 the normalization and slope of the background spectrum in the fit.
 This power law could be replaced by any spectral model that is found
@@ -166,87 +163,82 @@ analysis using the :ref:`ctlike` tool:
 .. code-block:: bash
   $ ctlike
-  Event list, counts cube or observation definition file [events.fits] cntcube.fits
-  Exposure cube file [NONE] expcube.fits
-  PSF cube file [NONE] psfcube.fits
+  Event list, counts cube or observation definition file [obs.xml] cntcube.fits
+  Exposure cube file (only needed for stacked analysis) [NONE] expcube.fits
+  PSF cube file (only needed for stacked analysis) [NONE] psfcube.fits
+  Background cube file (only needed for stacked analysis) [NONE] bkgcube.fits
   Source model [$CTOOLS/share/models/crab.xml] model.xml
-  Source model output file [crab_results.xml]
+  Source model output file [crab_results.xml] 
 :ref:`ctlike` recognises that a counts cube should be analysed and queries
-for the exposure cube and PSF cube file names.
-We specified the names of the files produced by the :ref:`ctexpcube` and
-:ref:`ctpsfcube` tools and furthermore provided the ``model.xml`` file
-that is generated by the :ref:`ctbkgcube` tool as source model.
+for the exposure cube, the PSF cube, and the background cube file names.
+We specified the names of the files produced by the :ref:`ctexpcube`,
+the :ref:`ctpsfcube` and the :ref:`ctbkgcube` tools.
+Furthermore we provided as model the ``model.xml`` file that has been
+generated by the :ref:`ctbkgcube` tool.
+
 The log file of the :ref:`ctlike` run is shown below.
 Note that the spectral model that is multiplied with the background
-cube has a Prefactor of 1.06 +/- 0.02 and an Index of 0.004 +/- 0.009,
+cube has a Prefactor of 1.09 +/- 0.02 and an Index of 0.03 +/- 0.01,
 indicating a very small correction to the actual spectrum of the background
 cube.
 Real life situations may of course require larger correction factors.
 .. code-block:: xml
-  2015-02-04T23:00:25: +=================================+
-  2015-02-04T23:00:25: | Maximum likelihood optimisation |
-  2015-02-04T23:00:25: +=================================+
-  2015-02-04T23:00:29:  >Iteration   0: -logL=36748.112, Lambda=1.0e-03
-  2015-02-04T23:00:33:  >Iteration   1: -logL=36734.170, Lambda=1.0e-03, delta=13.942, max(|grad|)=34.291911 [Index:8]
-  2015-02-04T23:00:38:  >Iteration   2: -logL=36734.127, Lambda=1.0e-04, delta=0.044, max(|grad|)=0.178088 [Index:8]
-  2015-02-04T23:00:42:  >Iteration   3: -logL=36734.127, Lambda=1.0e-05, delta=0.000, max(|grad|)=-0.001993 [Index:8]
-  2015-02-04T23:00:46: 
-  2015-02-04T23:00:46: +=========================================+
-  2015-02-04T23:00:46: | Maximum likelihood optimization results |
-  2015-02-04T23:00:46: +=========================================+
-  2015-02-04T23:00:46: === GOptimizerLM ===
-  2015-02-04T23:00:46:  Optimized function value ..: 36734.127
-  2015-02-04T23:00:46:  Absolute precision ........: 0.005
-  2015-02-04T23:00:46:  Acceptable value decrease .: 2
-  2015-02-04T23:00:46:  Optimization status .......: converged
-  2015-02-04T23:00:46:  Number of parameters ......: 11
-  2015-02-04T23:00:46:  Number of free parameters .: 4
-  2015-02-04T23:00:46:  Number of iterations ......: 3
-  2015-02-04T23:00:46:  Lambda ....................: 1e-06
-  2015-02-04T23:00:46:  Maximum log likelihood ....: -36734.127
-  2015-02-04T23:00:46:  Observed events  (Nobs) ...: 10685.000
-  2015-02-04T23:00:46:  Predicted events (Npred) ..: 10685.000 (Nobs - Npred = 2.06522e-06)
-  2015-02-04T23:00:46: === GModels ===
-  2015-02-04T23:00:46:  Number of models ..........: 2
-  2015-02-04T23:00:46:  Number of parameters ......: 11
-  2015-02-04T23:00:46: === GModelSky ===
-  2015-02-04T23:00:46:  Name ......................: Crab
-  2015-02-04T23:00:46:  Instruments ...............: all
-  2015-02-04T23:00:46:  Instrument scale factors ..: unity
-  2015-02-04T23:00:46:  Observation identifiers ...: all
-  2015-02-04T23:00:46:  Model type ................: PointSource
-  2015-02-04T23:00:46:  Model components ..........: "SkyDirFunction" * "PowerLaw" * "Constant"
-  2015-02-04T23:00:46:  Number of parameters ......: 6
-  2015-02-04T23:00:46:  Number of spatial par's ...: 2
-  2015-02-04T23:00:46:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
-  2015-02-04T23:00:46:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
-  2015-02-04T23:00:46:  Number of spectral par's ..: 3
-  2015-02-04T23:00:46:   Prefactor ................: 6.01471e-16 +/- 1.44183e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
-  2015-02-04T23:00:46:   Index ....................: -2.49533 +/- 0.0179769 [-0,-5]  (free,scale=-1,gradient)
-  2015-02-04T23:00:46:   PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
-  2015-02-04T23:00:46:  Number of temporal par's ..: 1
-  2015-02-04T23:00:46:   Constant .................: 1 (relative value) (fixed,scale=1,gradient)
-  2015-02-04T23:00:46: === GCTAModelCubeBackground ===
-  2015-02-04T23:00:46:  Name ......................: ctbkgcube default background model
-  2015-02-04T23:00:46:  Instruments ...............: all
-  2015-02-04T23:00:46:  Instrument scale factors ..: unity
-  2015-02-04T23:00:46:  Observation identifiers ...: all
-  2015-02-04T23:00:46:  Model type ................: "MapCubeFunction" * "PowerLaw" * "Constant"
-  2015-02-04T23:00:46:  Number of parameters ......: 5
-  2015-02-04T23:00:46:  Number of spatial par's ...: 1
-  2015-02-04T23:00:46:   Normalization ............: 1 [0.001,1000]  (fixed,scale=1,gradient)
-  2015-02-04T23:00:46:  Number of spectral par's ..: 3
-  2015-02-04T23:00:46:   Prefactor ................: 1.05876 +/- 0.0178132 [0,infty[ ph/cm2/s/MeV (free,scale=1,gradient)
-  2015-02-04T23:00:46:   Index ....................: 0.00396962 +/- 0.00881842 [-10,10]  (free,scale=1,gradient)
-  2015-02-04T23:00:46:   PivotEnergy ..............: 1e+06 MeV (fixed,scale=1e+06,gradient)
-  2015-02-04T23:00:46:  Number of temporal par's ..: 1
-  2015-02-04T23:00:46:   Constant .................: 1 (relative value) (fixed,scale=1,gradient)
-  2015-02-04T23:00:46: 
-  2015-02-04T23:00:46: +==============+
-  2015-02-04T23:00:46: | Save results |
-  2015-02-04T23:00:46: +==============+
-  2015-02-04T23:00:46: 
-  2015-02-04T23:00:46: Application "ctlike" terminated after 43 wall clock seconds, consuming 21.5624 seconds of CPU time.
+  2015-05-22T21:10:09: +=================================+
+  2015-05-22T21:10:09: | Maximum likelihood optimisation |
+  2015-05-22T21:10:09: +=================================+
+  2015-05-22T21:10:10:  >Iteration   0: -logL=60663.798, Lambda=1.0e-03
+  2015-05-22T21:10:11:  >Iteration   1: -logL=60642.413, Lambda=1.0e-03, delta=21.385, max(|grad|)=120.031102 [Index:7]
+  2015-05-22T21:10:13:  >Iteration   2: -logL=60642.233, Lambda=1.0e-04, delta=0.180, max(|grad|)=1.139187 [Index:7]
+  2015-05-22T21:10:14:  >Iteration   3: -logL=60642.233, Lambda=1.0e-05, delta=0.000, max(|grad|)=-0.001247 [Index:7]
+  2015-05-22T21:10:16: 
+  2015-05-22T21:10:16: +=========================================+
+  2015-05-22T21:10:16: | Maximum likelihood optimization results |
+  2015-05-22T21:10:16: +=========================================+
+  2015-05-22T21:10:16: === GOptimizerLM ===
+  2015-05-22T21:10:16:  Optimized function value ..: 60642.233 
+  2015-05-22T21:10:16:  Absolute precision ........: 0.005
+  2015-05-22T21:10:16:  Acceptable value decrease .: 2
+  2015-05-22T21:10:16:  Optimization status .......: converged
+  2015-05-22T21:10:16:  Number of parameters ......: 10
+  2015-05-22T21:10:16:  Number of free parameters .: 4
+  2015-05-22T21:10:16:  Number of iterations ......: 3
+  2015-05-22T21:10:16:  Lambda ....................: 1e-06
+  2015-05-22T21:10:16:  Maximum log likelihood ....: -60642.233
+  2015-05-22T21:10:16:  Observed events  (Nobs) ...: 25963.000
+  2015-05-22T21:10:16:  Predicted events (Npred) ..: 25963.000 (Nobs - Npred = 1.5749e-05)
+  2015-05-22T21:10:16: === GModels ===
+  2015-05-22T21:10:16:  Number of models ..........: 2
+  2015-05-22T21:10:16:  Number of parameters ......: 10
+  2015-05-22T21:10:16: === GModelSky ===
+  2015-05-22T21:10:16:  Name ......................: Crab
+  2015-05-22T21:10:16:  Instruments ...............: all
+  2015-05-22T21:10:16:  Instrument scale factors ..: unity
+  2015-05-22T21:10:16:  Observation identifiers ...: all
+  2015-05-22T21:10:16:  Model type ................: PointSource
+  2015-05-22T21:10:16:  Model components ..........: "SkyDirFunction" * "PowerLaw" * "Constant"
+  2015-05-22T21:10:16:  Number of parameters ......: 6
+  2015-05-22T21:10:16:  Number of spatial par's ...: 2
+  2015-05-22T21:10:16:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
+  2015-05-22T21:10:16:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
+  2015-05-22T21:10:16:  Number of spectral par's ..: 3
+  2015-05-22T21:10:16:   Prefactor ................: 5.8499e-16 +/- 8.02289e-18 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
+  2015-05-22T21:10:16:   Index ....................: -2.49909 +/- 0.0118976 [-0,-5]  (free,scale=-1,gradient)
+  2015-05-22T21:10:16:   PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T21:10:16:  Number of temporal par's ..: 1
+  2015-05-22T21:10:16:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+  2015-05-22T21:10:16: === GCTAModelCubeBackground ===
+  2015-05-22T21:10:16:  Name ......................: BackgroundModel
+  2015-05-22T21:10:16:  Instruments ...............: CTA, HESS, MAGIC, VERITAS
+  2015-05-22T21:10:16:  Instrument scale factors ..: unity
+  2015-05-22T21:10:16:  Observation identifiers ...: all
+  2015-05-22T21:10:16:  Model type ................: "PowerLaw" * "Constant"
+  2015-05-22T21:10:16:  Number of parameters ......: 4
+  2015-05-22T21:10:16:  Number of spectral par's ..: 3
+  2015-05-22T21:10:16:   Prefactor ................: 1.08781 +/- 0.0212678 [0,infty[ ph/cm2/s/MeV (free,scale=1,gradient)
+  2015-05-22T21:10:16:   Index ....................: 0.0269681 +/- 0.0111532 [-10,10]  (free,scale=1,gradient)
+  2015-05-22T21:10:16:   PivotEnergy ..............: 1e+06 MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T21:10:16:  Number of temporal par's ..: 1
+  2015-05-22T21:10:16:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
@@ -6,7 +6,7 @@ Generating a Test Statistic map
 In the above examples we were always considering that we knew exactly
 where the source is located, but in real life situations this may not
 necessarily be the case.
-We thus need a mean to search for significant excess emission in our data.
+We thus need a means to search for significant excess emission in our data.
 One possibility to accomplish this task is to compute the Test Statistic
 for a specific source at a grid of trial positions.
 The :ref:`cttsmap` tool will exactly do that job.
@@ -32,25 +32,25 @@ nominal position of the Crab nebula with a grid spacing of 0.02 degrees.
   $ cttsmap
   Input event list, counts cube or observation definition file [events.fits] obs.xml
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Source model [$CTOOLS/share/models/crab.xml] 
   Test source [Crab] 
-  First coordinate of image center in degrees (RA or galactic l) [83.63] 
-  Second coordinate of image center in degrees (DEC or galactic b) [22.01] 
-  Projection method e.g. AIT|AZP|CAR|MER|STG|TAN (AIT|AZP|CAR|MER|STG|TAN) [CAR] 
+  First coordinate of image center in degrees (RA or galactic l) (0-360) [83.63] 
+  Second coordinate of image center in degrees (DEC or galactic b) (-90-90) [22.01] 
+  Projection method (AIT|AZP|CAR|MER|MOL|STG|TAN) [CAR] 
   Coordinate system (CEL - celestial, GAL - galactic) (CEL|GAL) [CEL] 
   Image scale (in degrees/pixel) [0.02] 
   Size of the X axis in pixels [200] 10
   Size of the Y axis in pixels [200] 10
-  Output Test Statistic map [tsmap.fits]
+  Output Test Statistic map [tsmap.fits] 
 :ref:`cttsmap` writes the Test Statistic map in the ``tsmap.fits`` file
 that contains one extension for the Test Statistic value and further
 extensions for the spectral parameters that have been fitted for the
 source at each position of the grid.
 The figure below show the Test Statistic map which reaches a maximum
-value of 13782 near the centre of the map.
+value of 37484 near the centre of the map.
 .. figure:: tsmap-crab.png
    :height: 400px
@@ -319,22 +319,6 @@ CTA cube background
   .. code-block:: xml
     <source name="Background" type="CTACubeBackground" instrument="CTA">
-      <spatialModel type="GaussFunction">
-        <parameter name="RA"    scale="1.0" value="83.6221" min="0.0"   max="360.0" free="0"/>
-        <parameter name="DEC"   scale="1.0" value="22.01"   min="-90.0" max="90.0"  free="0"/>
-        <parameter name="Sigma" scale="1.0" value="1.0"     min="0.0"   max="10.0"  free="0"/>
-      </spatialModel>
-      <spectrum type="...">
-        ...
-      </spectrum>
-    </source>
-
-  .. code-block:: xml
-
-    <source name="Background" type="CTACubeBackground" instrument="CTA">
-      <spatialModel type="MapCubeFunction" file="map_cube.fits">
-        <parameter name="Normalization" scale="1" value="1" min="0.001" max="1000.0" free="0"/>
-      </spatialModel>
       <spectrum type="...">
         ...
       </spectrum>
@@ -20,7 +20,7 @@ generally not very useful.
 .. warning::
-   The GammaLib module needs to be imported before the ctools 
+   The GammaLib module needs to be imported **before** the ctools 
    module so that GammaLib class types are handled correctly.
    So make sure that you **always import the GammaLib module before the ctools
    module.**
@@ -42,8 +42,8 @@ way as if it were executed from the command line.
    >>> sim = ctools.ctobssim()
    >>> sim["inmodel"].filename("${CTOOLS}/share/models/crab.xml")
    >>> sim["outevents"].filename("events.fits")
-   >>> sim["caldb"].string("dummy")
-   >>> sim["irf"].string("cta_dummy_irf")
+   >>> sim["caldb"].string("prod2")
+   >>> sim["irf"].string("South_50h")
    >>> sim["ra"].real(83.63)
    >>> sim["dec"].real(22.01)
    >>> sim["rad"].real(5.0)
@@ -95,13 +95,13 @@ in that observation you may type:
    >>> print(sim.obs()[0].events())
    === GCTAEventList ===
-    Number of events ..........: 12152
+    Number of events ..........: 31928
     Time interval .............: 51544.5 - 51544.5 days
    === GEbounds ===
     Number of intervals .......: 1
     Energy range ..............: 100 GeV - 100 TeV
    === GCTARoi ===
-    ROI centre ................: RA=83.63, DEC=22.01 [0, 0]
+    ROI centre ................: RA=83.63, DEC=22.01 [0,0]
     ROI radius ................: 5 deg
 The ``obs()[0]`` operator returns the first observation in the observation 
@@ -114,11 +114,12 @@ To see what kind of object you actually got, use:
    >>> type(sim.obs()[0].events())
    <class 'gammalib.cta.GCTAEventList'>
-The CTA event list is implement as the ``GCTAEventList`` class in the 
+The CTA event list is implement by the ``GCTAEventList`` class in the 
 ``cta`` module of GammaLib.
 To visualise the individual events you can iterate over the events using a 
 for loop.
-This will show the simulated celestial coordinates (RA, DEC), energies and 
+This will show the simulated celestial coordinates (RA, DEC), the 
+coordinate in the camera system [DETX, DETY], energies and 
 terrestrial times (TT) of all events. 
 .. code-block:: python
@@ -127,9 +128,9 @@ terrestrial times (TT) of all events.
    >>> for event in events:
    ...     print(event)
    ...
-   Dir=RA=83.6477, DEC=22.0202 [0.000178879, 0.000286637] Energy=3.235 TeV Time=-3.15576e+08 s (TT)
-   Dir=RA=83.482, DEC=22.0189 [0.000155643, -0.00239471] Energy=141.949 GeV Time=-3.15576e+08 s (TT)
-   Dir=RA=83.6058, DEC=22.1586 [0.00259306, -0.000391919] Energy=316.376 GeV Time=-3.15576e+08 s (TT)
+   Dir=RA=83.5759, DEC=21.9757 [-0.000597839,-0.000875483] Energy=328.153 GeV Time=-3.15576e+08 s (TT)
+   Dir=RA=83.5775, DEC=22.0128 [4.88943e-05,-0.000849504] Energy=353.264 GeV Time=-3.15576e+08 s (TT)
+   Dir=RA=83.4302, DEC=21.931 [-0.00137616,-0.0032343] Energy=117.999 GeV Time=-3.15576e+08 s (TT)
    ...
 .. note::
@@ -150,13 +151,13 @@ terrestrial times (TT) of all events.
      >>> sim.run()                   
      >>> print(sim.obs()[0].events())
      === GCTAEventList ===
-      Number of events ..........: 6141
+      Number of events ..........: 15964
       Time interval .............: 51544.5 - 51544.5 days
      === GEbounds ===
       Number of intervals .......: 1
       Energy range ..............: 100 GeV - 100 TeV
      === GCTARoi ===
-      ROI centre ................: RA=83.63, DEC=22.01 [0, 0]
+      ROI centre ................: RA=83.63, DEC=22.01 [0,0]
       ROI radius ................: 5 deg
@@ -178,28 +179,30 @@ Where are the user parameters?
 :ref:`ctlike` doesn't in fact need any as all the relevant information is 
 already contained in the observation container produced by the 
 :ref:`ctobssim` class.
-And you make have recognised that we constructed the :ref:`ctlike` 
+And you may have recognised that we constructed the :ref:`ctlike` 
 instance by using the :ref:`ctobssim` observation container as
 constructor argument.
-To check how the fit went you may inspect the optimiser class used by
+To check how the fit went you can inspect the optimiser used by
 :ref:`ctlike`:
 .. code-block:: python
    >>> print(like.opt())         
    === GOptimizerLM ===
-    Optimized function value ..: 44578.761
+    Optimized function value ..: 101048.382
     Absolute precision ........: 0.005
     Acceptable value decrease .: 2
     Optimization status .......: converged
-    Number of parameters ......: 9
+    Number of parameters ......: 10
     Number of free parameters .: 4
-    Number of iterations ......: 3
-    Lambda ....................: 1e-06
+    Number of iterations ......: 2
+    Lambda ....................: 1e-05
-Apparently, the fit converged after fitting 4 parameters in 3 iterations.
-To inspect the fit results you may print the model container that is a 
+Apparently, the fit converged after 2 iterations.
+Out of 10 parameters in the model 4 have been fitted (the others were kept
+fixed).
+To inspect the fit results you can print the model container that is a 
 member of the observation container:
 .. code-block:: python
@@ -207,7 +210,7 @@ member of the observation container:
    >>> print(like.obs().models())
    === GModels ===
     Number of models ..........: 2
-    Number of parameters ......: 9
+    Number of parameters ......: 10
    === GModelSky ===
     Name ......................: Crab
     Instruments ...............: all
@@ -220,24 +223,24 @@ member of the observation container:
      RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
      DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
     Number of spectral par's ..: 3
-     Prefactor ................: 6.13265e-16 +/- 2.05734e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
-     Index ....................: -2.50565 +/- 0.0250831 [-0,-5]  (free,scale=-1,gradient)
+     Prefactor ................: 5.80565e-16 +/- 1.11147e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
+     Index ....................: -2.49266 +/- 0.0163476 [-0,-5]  (free,scale=-1,gradient)
      PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
     Number of temporal par's ..: 1
-     Constant .................: 1 (relative value) (fixed,scale=1,gradient)
-   === GCTAModelRadialAcceptance ===
-    Name ......................: Background
+     Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+   === GCTAModelIrfBackground ===
+    Name ......................: CTABackgroundModel
     Instruments ...............: CTA
     Instrument scale factors ..: unity
     Observation identifiers ...: all
-    Model type ................: "Gaussian" * "FileFunction" * "Constant"
-    Number of parameters ......: 3
-    Number of radial par's ....: 1
-     Sigma ....................: 3.03693 +/- 0.0304896 [0.01,10] deg2 (free,scale=1,gradient)
-    Number of spectral par's ..: 1
-     Normalization ............: 0.998667 +/- 0.0172585 [0,1000]  (free,scale=1,gradient)
+    Model type ................: "PowerLaw" * "Constant"
+    Number of parameters ......: 4
+    Number of spectral par's ..: 3
+     Prefactor ................: 1.01701 +/- 0.0216563 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
+     Index ....................: 0.0230864 +/- 0.0123006 [-5,5]  (free,scale=1,gradient)
+     PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
     Number of temporal par's ..: 1
-     Constant .................: 1 (relative value) (fixed,scale=1,gradient)
+     Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
 Suppose you want to repeat the fit by optimising also the position of the 
 point source.
@@ -250,8 +253,8 @@ This is easy from Python:
    >>> like.run()
    >>> print(like.obs().models())
    ...
-     RA .......................: 83.633 +/- 0.00137216 [-360,360] deg (free,scale=1)
-     DEC ......................: 22.0144 +/- 0.00127247 [-90,90] deg (free,scale=1)
+     RA .......................: 83.6341 +/- 0.000958947 [-360,360] deg (free,scale=1)
+     DEC ......................: 22.0139 +/- 0.000870912 [-90,90] deg (free,scale=1)
 The ``like.obs().models()`` method provides the model container, using the 
 ``["Crab"]`` operator we access the Crab model in that container and using 
@@ -293,12 +296,12 @@ the ``fix()`` method.
       >>> print(like.obs())
       === GObservations ===
        Number of observations ....: 1
-       Number of predicted events : 6141
+       Number of predicted events : 15964
    is okay as the ``ctools.ctlike(sim.obs())`` constructor will create
    a copy of the observation container that lives within the :ref:`ctlike`
    instance.
-   To preserve an observation container after a ctools object went out 
+   To preserve an observation container after a ctools object goes out 
    of scope you have to create a local copy of the container using the
    ``copy()`` method:
@@ -327,8 +330,8 @@ Here an example of how to use ``obsutils``:
    >>> import ctools
    >>> from ctools import obsutils
    >>> pattern = obsutils.set_obs_patterns("four", ra=83.63, dec=22.01, offset=1.0)
-   >>> obs = obsutils.set_obs_list(pattern, duration=1800, emin=0.1, emax=100.0, rad=5.0, caldb="dummy", irf="cta_dummy_irf")
-   >>> print(like.obs())   
+   >>> obs = obsutils.set_obs_list(pattern, duration=1800, emin=0.1, emax=100.0, rad=5.0, caldb="prod2", irf="South_50h")
+   >>> print(obs)   
    === GObservations ===
     Number of observations ....: 4
     Number of predicted events : 0
@@ -339,7 +342,7 @@ Here an example of how to use ``obsutils``:
    >>> print(like.obs().models())   
    === GModels ===
     Number of models ..........: 2
-    Number of parameters ......: 9
+    Number of parameters ......: 10
    ...
 The module is imported using the ``from ctools import obsutils`` directive.
@@ -39,14 +39,14 @@ type :ref:`ctobssim` at the console prompt (which is denoted by ``$``).
   RA of pointing (degrees) (0-360) [83.63] 
   Dec of pointing (degrees) (-90-90) [22.01] 
   Radius of FOV (degrees) (0-180) [5.0] 
-  Start time (MET in s) (0) [0.0] 
-  End time (MET in s) (0) [1800.0] 
-  Lower energy limit (TeV) (0) [0.1] 
-  Upper energy limit (TeV) (0) [100.0] 
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
+  Start time (MET in s) [0.0] 
+  End time (MET in s) [1800.0] 
+  Lower energy limit (TeV) [0.1] 
+  Upper energy limit (TeV) [100.0] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Model [$CTOOLS/share/models/crab.xml] 
-  Output event data file or observation definition file [events.fits] 
+  Output event data file or observation definition file [events.fits]
 Each line represents a query for one parameter value.
 The line starts with a short description of the parameter, followed by 
@@ -62,8 +62,18 @@ You may have recognised that the environment variable ``$CTOOLS`` has
 been used in the path name of the model. ctools will automatically expand
 the environment variables in parameter inputs.
-The CTA instrument properties (effective area, PSF width) are taken for
-the moment from a dummy calibration database.
+The CTA instrument properties (effective area, PSF width) are taken from
+the prod2 database. The response for the southern array using the cuts
+optimised for 50 hours of observing time are used.
+
+.. note::
+
+   ctools comes bundled with CTA response functions for the northern and
+   the southern array. The response functions are based on a Prod2
+   analysis with cuts optimised for 0.5 hours, 5 hours and 50 hours of
+   observing time. The following Instrument Response Functions
+   (IRFs) are available: ``North_0.5h``, ``North_5h``, ``North_50h``,
+   ``South_0.5h``, ``South_5h``, and ``South_50h``.
 Events are simulated based on the instrument properties and based on a
 source and background model. Only events that fall within the specified
@@ -90,14 +100,13 @@ The source and background model is defined by the XML file
         <parameter name="DEC" scale="1.0" value="22.0145" min="-90"  max="90"  free="0"/>
       </spatialModel>
     </source>
-    <source name="Background" type="RadialAcceptance" instrument="CTA">
-      <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-        <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="1.0" free="1"/>
+    <source name="CTABackgroundModel" type="CTAIrfBackground" instrument="CTA">
+      <spectrum type="PowerLaw">  
+        <parameter name="Prefactor" scale="1.0"  value="1.0"  min="1e-3" max="1e+3"   free="1"/>  
+        <parameter name="Index"     scale="1.0"  value="0.0"  min="-5.0" max="+5.0"   free="1"/>  
+        <parameter name="Scale"     scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>  
       </spectrum>
-      <radialModel type="Gaussian">
-         <parameter name="Sigma" scale="1.0" value="3.0" min="0.01" max="10.0" free="1"/>
-      </radialModel>
-    </source>
+    </source>     
   </source_library>
 The model consists of a source library that contains two components:
@@ -107,14 +116,15 @@ The Crab nebula is modelled by a factorized sky model that has a spectral
 and a spatial component (tags ``<spectrum>`` and ``<spatialModel>``,
 respectively). The spectrum is modelled by a power law, which is defined by 
 three parameters: the ``Prefactor``, the ``Index`` and the ``Scale``.
-The spatial model has two parameters: Right Ascension in degrees (RA), and 
-Declination in degrees (DEC). Each parameter has a value and a scale factor, 
+The spatial model has two parameters: Right Ascension in degrees (``RA``), and 
+Declination in degrees (``DEC``). Each parameter has a value and a scale factor, 
 the real value of the parameter being the product ``value * scale``. Typically,
-scale is chosen so that value is of the order of 1 (this is relevant for 
-model fitting). In addition, value is bound by a minimum (min) and 
-maximum (max) value, and a parameter may be free (``free="1"``) or fixed
-(``free="0"``). The min, max, and free attributes are not relevant here for
-the simulations, but they will be important for the model fitting later.
+``scale`` is chosen so that ``value`` is of the order of 1 (this is relevant for 
+model fitting). In addition, ``value`` is bound by a minimum (``min``) and 
+maximum (``max``) value, and a parameter may be free (``free="1"``) or fixed
+(``free="0"``). The ``min``, ``max``, and ``free`` attributes are not
+relevant here for the simulations, but they will be important for the model 
+fitting later.
 The spectral intensity I(E) (in units of photons/cm2/s/MeV) of the power
 law is given by 
@@ -128,19 +138,20 @@ where the parameters in the XML definition have the following mappings:
 * :math:`\gamma` = ``Index``
 * :math:`E_0` = ``Scale``
-Note that energies are given in MeV. This is a GammaLib convention that
-can not be modified. **So make sure you always use MeV as energy unit in
-an XML file.**
+.. warning::
+
+   Energies are given in the XML file in MeV units. This is a GammaLib
+   convention that can not be modified. **So make sure you always use 
+   MeV as energy unit in an XML file.**
-The instrumental background of CTA is modelled by a factorized data model
-that has a spectral and a radial component (tags ``<spectrum>`` and
-``<radialModel>``, respectively). The spectral component describes the 
-on-axis background counting rate of CTA as function of energy in units of 
-counts/s/sr/MeV. The radial component describes the variation of the 
-background rate with offset angle squared, (i.e. square of the offset angle 
-with respect to the pointing direction) which is modelled here by a Gaussian. 
-The only parameter of the radial component is the width of the Gaussian Sigma, 
-which is here set to 3 degrees squared.
+The instrumental background of CTA is modelled using the background
+information provided in the IRF (``CTAIrfBackground``) multipled
+by a power law. As it is defined here, the power law represents a
+constant of 1, hence the background IRF will be used without any
+modification. The power law will become active when fitting the data
+later and allows a spectral adjustment of the background model that
+may account for uncertainties in the background information provided
+in the IRF.
 :ref:`ctobssim` has a couple of hidden parameters, the most important one being
 certainly ``seed``. ``seed`` is an integer that specifies the seed value
@@ -167,25 +178,25 @@ from this file are shown here:
 .. code-block:: xml
-  2015-02-04T15:04:56: +======================+
-  2015-02-04T15:04:56: | Simulate observation |
-  2015-02-04T15:04:56: +======================+
-  2015-02-04T15:04:56: === Observation ===
-  2015-02-04T15:04:56:  Simulation area ...........: 1.9635e+11 cm2
-  2015-02-04T15:04:56:  Simulation cone ...........: RA=83.63 deg, Dec=22.01 deg, r=5.5 deg
-  2015-02-04T15:04:56:  Time interval .............: 0 - 1800 s
-  2015-02-04T15:04:56:  Photon energy range .......: 100 GeV - 100 TeV
-  2015-02-04T15:04:56:  Event energy range ........: 100 GeV - 100 TeV
-  2015-02-04T15:04:56:  MC source photons .........: 207547 [Crab]
-  2015-02-04T15:04:56:  MC source events ..........: 995 [Crab]
-  2015-02-04T15:04:56:  MC source events ..........: 995 (all source models)
-  2015-02-04T15:04:56:  MC background events ......: 5146
-  2015-02-04T15:04:56:  MC events .................: 6141 (all models)
+  2015-05-22T19:24:53: +======================+
+  2015-05-22T19:24:53: | Simulate observation |
+  2015-05-22T19:24:53: +======================+
+  2015-05-22T19:24:53: === Observation ===
+  2015-05-22T19:24:53:  Simulation area ...........: 1.9635e+11 cm2
+  2015-05-22T19:24:53:  Simulation cone ...........: RA=83.63 deg, Dec=22.01 deg, r=5.5 deg
+  2015-05-22T19:24:53:  Time interval .............: 0 - 1800 s
+  2015-05-22T19:24:53:  Photon energy range .......: 100 GeV - 100 TeV
+  2015-05-22T19:24:53:  Event energy range ........: 100 GeV - 100 TeV
+  2015-05-22T19:24:53:  MC source photons .........: 207547 [Crab]
+  2015-05-22T19:24:53:  MC source events ..........: 2961 [Crab]
+  2015-05-22T19:24:53:  MC source events ..........: 2961 (all source models)
+  2015-05-22T19:24:55:  MC background events ......: 13003
+  2015-05-22T19:24:55:  MC events .................: 15964 (all models)
 Each line starts with the UTC time at which the line has been written. In
 this run, 207547 Crab photons have been thrown over an area of 19.6 square
-kilometres during a time interval of 1800 seconds. 995 of these photons have
-been registered by CTA as events. In the same time interval, 5146 background
+kilometres during a time interval of 1800 seconds. 2961 of these photons have
+been registered by CTA as events. In the same time interval, 13003 background
 events have been registred by CTA.
 You may change the name of the log file using the hidden parameter 
@@ -233,9 +244,9 @@ control).
   $ ctbin
   Event list or observation definition file [events.fits] 
-  First coordinate of image center in degrees (RA or galactic l) [83.63] 
-  Second coordinate of image center in degrees (DEC or galactic b) [22.01] 
-  Projection method e.g. AIT|AZP|CAR|MER|STG|TAN (AIT|AZP|CAR|MER|STG|TAN) [CAR] 
+  First coordinate of image center in degrees (RA or galactic l) (0-360) [83.63] 
+  Second coordinate of image center in degrees (DEC or galactic b) (-90-90) [22.01] 
+  Projection method (AIT|AZP|CAR|MER|MOL|STG|TAN) [CAR] 
   Coordinate system (CEL - celestial, GAL - galactic) (CEL|GAL) [CEL] 
   Image scale (in degrees/pixel) [0.02] 
   Size of the X axis in pixels [200] 
@@ -244,7 +255,7 @@ control).
   Start value for first energy bin in TeV [0.1] 
   Stop value for last energy bin in TeV [100.0] 
   Number of energy bins [20] 
-  Output counts cube [cntcube.fits]
+  Output counts cube [cntcube.fits] 
 The counts cube will be centred on the location of the Crab (Right Ascension
 83.63 degrees, Declination 22.01 degrees) and will be aligned in celestial
@@ -285,41 +296,18 @@ reproduced below:
 .. code-block:: xml
-  2015-02-04T15:33:34: +====================+
-  2015-02-04T15:33:34: | Binned observation |
-  2015-02-04T15:33:34: +====================+
-  2015-02-04T15:33:34: === GObservations ===
-  2015-02-04T15:33:34:  Number of observations ....: 1
-  2015-02-04T15:33:34:  Number of predicted events : 0
-  2015-02-04T15:33:34: === GCTAObservation ===
-  2015-02-04T15:33:34:  Name ......................: 
-  2015-02-04T15:33:34:  Identifier ................: 
-  2015-02-04T15:33:34:  Instrument ................: CTA
-  2015-02-04T15:33:34:  Event file ................: 
-  2015-02-04T15:33:34:  Event type ................: CountsCube
-  2015-02-04T15:33:34:  Statistics ................: Poisson
-  2015-02-04T15:33:34:  Ontime ....................: 1800 s
-  2015-02-04T15:33:34:  Livetime ..................: 1710 s
-  2015-02-04T15:33:34:  Deadtime correction .......: 0.95
-  2015-02-04T15:33:34:  User energy range .........: undefined
-  2015-02-04T15:33:34: === GCTAPointing ===
-  2015-02-04T15:33:34:  Pointing direction ........: (RA,Dec)=(83.63,22.01)
-  2015-02-04T15:33:34:  Response function .........: undefined
-  2015-02-04T15:33:34: === GCTAEventCube ===
-  2015-02-04T15:33:34:  Number of events ..........: 5542
-  2015-02-04T15:33:34:  Number of elements ........: 800000
-  2015-02-04T15:33:34:  Number of pixels ..........: 40000
-  2015-02-04T15:33:34:  Number of energy bins .....: 20
-  2015-02-04T15:33:34:  Time interval .............: 0 - 1800 sec
-  2015-02-04T15:33:34: 
-  2015-02-04T15:33:34: +==================+
-  2015-02-04T15:33:34: | Save observation |
-  2015-02-04T15:33:34: +==================+
-  2015-02-04T15:33:34: 
-  2015-02-04T15:33:34: Application "ctbin" terminated after 14 wall clock seconds, consuming 0.332675 seconds of CPU time.
-
-From the 6141 events that have been simulated and stored in the 
-``events.fits`` file, 5542 lie within the cube boundaries and are thus put
+  2015-05-22T19:40:30: +=================+
+  2015-05-22T19:40:30: | Bin observation |
+  2015-05-22T19:40:30: +=================+
+  2015-05-22T19:40:30: === Observation ===
+  2015-05-22T19:40:30:  Events in list ............: 15964
+  2015-05-22T19:40:30:  Events in cube ............: 13360
+  2015-05-22T19:40:30:  Event bins outside RoI ....: 0
+  2015-05-22T19:40:30:  Events outside cube area ..: 2604
+  2015-05-22T19:40:30:  Events outside energy bins : 0
+
+From the 15964 events that have been simulated and stored in the 
+``events.fits`` file, 13360 lie within the cube boundaries and are thus put
 into the resulting counts cube. The counts cube is stored in a cartesian
 projection in a World Coordinate System (WCS) compliant format.
@@ -340,12 +328,13 @@ and we do the fit by typing:
   Event list, counts cube or observation definition file [events.fits] cntcube.fits
   Exposure cube file (only needed for stacked analysis) [NONE] 
   PSF cube file (only needed for stacked analysis) [NONE] 
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
+  Background cube file (only needed for stacked analysis) [NONE] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
   Source model [$CTOOLS/share/models/crab.xml] 
   Source model output file [crab_results.xml]
-Fitting of the data is done in *binned* mode, which means that the events
+The data are fitted in *binned* mode, which means that the events
 have been binned into a counts cube and the fit computes the log-likelihood
 function by summing over all 200 x 200 x 20 bins of the counts cube. There is
 an alternative method, the so called *unbinned* mode, where the events are
@@ -356,8 +345,8 @@ One of the parameters of :ref:`ctlike` is a source model output file
 (we specified ``crab_results.xml`` in the example), and this file will be
 a copy of the source model input XML file where the parameter values have
 been replaced by the fit results. In addition, the statistical uncertainties
-are added for each fitted parameter using the attribute error. Below we show 
-the XML result file that has been produced by the run:
+are added for each fitted parameter using the attribute ``error``.
+Below we show the XML result file that has been produced by the run:
 .. code-block:: xml
@@ -365,8 +354,8 @@ the XML result file that has been produced by the run:
   <source_library title="source library">
     <source name="Crab" type="PointSource">
       <spectrum type="PowerLaw">
-        <parameter name="Prefactor" value="6.07928" error="0.204582" scale="1e-16" min="1e-07" max="1000" free="1" />
-        <parameter name="Index" value="2.5009" error="0.0252057" scale="-1" min="0" max="5" free="1" />
+        <parameter name="Prefactor" value="5.78197" error="0.110889" scale="1e-16" min="1e-07" max="1000" free="1" />
+        <parameter name="Index" value="2.49174" error="0.0164369" scale="-1" min="0" max="5" free="1" />
         <parameter name="Scale" value="0.3" scale="1e+06" min="0.01" max="1000" free="0" />
       </spectrum>
       <spatialModel type="SkyDirFunction">
@@ -374,20 +363,31 @@ the XML result file that has been produced by the run:
         <parameter name="DEC" value="22.0145" scale="1" min="-90" max="90" free="0" />
       </spatialModel>
     </source>
-    <source name="Background" type="RadialAcceptance" instrument="CTA">
-      <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-        <parameter name="Normalization" value="0.990708" error="0.0192887" scale="1" min="0" max="1000" free="1" />
+    <source name="CTABackgroundModel" type="CTAIrfBackground" instrument="CTA">
+      <spectrum type="PowerLaw">
+        <parameter name="Prefactor" value="1.02048" error="0.0280684" scale="1" min="0.001" max="1000" free="1" />
+        <parameter name="Index" value="0.0197856" error="0.0156686" scale="1" min="-5" max="5" free="1" />
+        <parameter name="Scale" value="1" scale="1e+06" min="0.01" max="1000" free="0" />
       </spectrum>
-      <radialModel type="Gaussian">
-        <parameter name="Sigma" value="3.02523" error="0.0594117" scale="1" min="0.01" max="10" free="1" />
-      </radialModel>
     </source>
   </source_library>
 In this example, the ``Prefactor`` and ``Index`` of the spectral model for the
-Crab as well as the ``Normalization`` and ``Sigma`` parameter of the radial
-acceptance model have been fitted (all parameters having the attribute 
-``free="1"`` are fitted).
+Crab as well as the ``Prefactor`` and ``Index`` of the background spectral
+model have been fitted (all parameters having the attribute ``free="1"`` are
+fitted). Now we see also the effect of having multiplied the background 
+model with a power law. In that way, the amplitude of the background as
+well as it's spectral slope is adjusted by the fit. Obviously, in this
+example the adjustment compensates only for the statistical fluctuations
+of the background, but with real data, the adjustment may account also for
+some of the systematic uncertainties.
+
+.. warning::
+
+   As good practice, the amplitude of the background model should always be
+   left as a free parameter of the fit. Otherwise, any uncertainty in the
+   background rate will immediately propagate into the flux estimate of the 
+   source. 
 To get more details about the model fitting you can inspect the log file.
 Below the last lines of the ctlike.log log file that has been produced by
@@ -395,70 +395,63 @@ this run:
 .. code-block:: xml
-  2015-02-04T15:50:43: +=================================+
-  2015-02-04T15:50:43: | Maximum likelihood optimisation |
-  2015-02-04T15:50:43: +=================================+
-  2015-02-04T15:50:44:  >Iteration   0: -logL=22752.553, Lambda=1.0e-03
-  2015-02-04T15:50:46:  >Iteration   1: -logL=22750.640, Lambda=1.0e-03, delta=1.913, max(|grad|)=-0.941594 [Index:3]
-  2015-02-04T15:50:47:  >Iteration   2: -logL=22750.632, Lambda=1.0e-04, delta=0.008, max(|grad|)=-0.079270 [Index:3]
-  2015-02-04T15:50:48:  >Iteration   3: -logL=22750.632, Lambda=1.0e-05, delta=0.000, max(|grad|)=-0.003073 [Index:3]
-  2015-02-04T15:50:49: 
-  2015-02-04T15:50:49: +=========================================+
-  2015-02-04T15:50:49: | Maximum likelihood optimization results |
-  2015-02-04T15:50:49: +=========================================+
-  2015-02-04T15:50:49: === GOptimizerLM ===
-  2015-02-04T15:50:49:  Optimized function value ..: 22750.632
-  2015-02-04T15:50:49:  Absolute precision ........: 0.005
-  2015-02-04T15:50:49:  Acceptable value decrease .: 2
-  2015-02-04T15:50:49:  Optimization status .......: converged
-  2015-02-04T15:50:49:  Number of parameters ......: 9
-  2015-02-04T15:50:49:  Number of free parameters .: 4
-  2015-02-04T15:50:49:  Number of iterations ......: 3
-  2015-02-04T15:50:49:  Lambda ....................: 1e-06
-  2015-02-04T15:50:49:  Maximum log likelihood ....: -22750.632
-  2015-02-04T15:50:49:  Observed events  (Nobs) ...: 5542.000
-  2015-02-04T15:50:49:  Predicted events (Npred) ..: 5542.000 (Nobs - Npred = 3.3498e-06)
-  2015-02-04T15:50:49: === GModels ===
-  2015-02-04T15:50:49:  Number of models ..........: 2
-  2015-02-04T15:50:49:  Number of parameters ......: 9
-  2015-02-04T15:50:49: === GModelSky ===
-  2015-02-04T15:50:49:  Name ......................: Crab
-  2015-02-04T15:50:49:  Instruments ...............: all
-  2015-02-04T15:50:49:  Instrument scale factors ..: unity
-  2015-02-04T15:50:49:  Observation identifiers ...: all
-  2015-02-04T15:50:49:  Model type ................: PointSource
-  2015-02-04T15:50:49:  Model components ..........: "SkyDirFunction" * "PowerLaw" * "Constant"
-  2015-02-04T15:50:49:  Number of parameters ......: 6
-  2015-02-04T15:50:49:  Number of spatial par's ...: 2
-  2015-02-04T15:50:49:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
-  2015-02-04T15:50:49:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
-  2015-02-04T15:50:49:  Number of spectral par's ..: 3
-  2015-02-04T15:50:49:   Prefactor ................: 6.07928e-16 +/- 2.04582e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
-  2015-02-04T15:50:49:   Index ....................: -2.5009 +/- 0.0252057 [-0,-5]  (free,scale=-1,gradient)
-  2015-02-04T15:50:49:   PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
-  2015-02-04T15:50:49:  Number of temporal par's ..: 1
-  2015-02-04T15:50:49:   Constant .................: 1 (relative value) (fixed,scale=1,gradient)
-  2015-02-04T15:50:49: === GCTAModelRadialAcceptance ===
-  2015-02-04T15:50:49:  Name ......................: Background
-  2015-02-04T15:50:49:  Instruments ...............: CTA
-  2015-02-04T15:50:49:  Instrument scale factors ..: unity
-  2015-02-04T15:50:49:  Observation identifiers ...: all
-  2015-02-04T15:50:49:  Model type ................: "Gaussian" * "FileFunction" * "Constant"
-  2015-02-04T15:50:49:  Number of parameters ......: 3
-  2015-02-04T15:50:49:  Number of radial par's ....: 1
-  2015-02-04T15:50:49:   Sigma ....................: 3.02523 +/- 0.0594117 [0.01,10] deg2 (free,scale=1,gradient)
-  2015-02-04T15:50:49:  Number of spectral par's ..: 1
-  2015-02-04T15:50:49:   Normalization ............: 0.990708 +/- 0.0192887 [0,1000]  (free,scale=1,gradient)
-  2015-02-04T15:50:49:  Number of temporal par's ..: 1
-  2015-02-04T15:50:49:   Constant .................: 1 (relative value) (fixed,scale=1,gradient)
-  2015-02-04T15:50:49: 
-  2015-02-04T15:50:49: +==============+
-  2015-02-04T15:50:49: | Save results |
-  2015-02-04T15:50:49: +==============+
-  2015-02-04T15:50:49: 
-  2015-02-04T15:50:49: Application "ctlike" terminated after 15 wall clock seconds, consuming 5.92686 seconds of CPU time.
-
-The maximum likelihood optimizer required 3 iterations to converge. This
+  2015-05-22T19:45:15: +=================================+
+  2015-05-22T19:45:15: | Maximum likelihood optimisation |
+  2015-05-22T19:45:15: +=================================+
+  2015-05-22T19:45:16:  >Iteration   0: -logL=39910.748, Lambda=1.0e-03
+  2015-05-22T19:45:18:  >Iteration   1: -logL=39908.788, Lambda=1.0e-03, delta=1.960, max(|grad|)=6.701856 [Index:7]
+  2015-05-22T19:45:19:  >Iteration   2: -logL=39908.787, Lambda=1.0e-04, delta=0.001, max(|grad|)=-0.057041 [Index:3]
+  2015-05-22T19:45:21: 
+  2015-05-22T19:45:21: +=========================================+
+  2015-05-22T19:45:21: | Maximum likelihood optimization results |
+  2015-05-22T19:45:21: +=========================================+
+  2015-05-22T19:45:21: === GOptimizerLM ===
+  2015-05-22T19:45:21:  Optimized function value ..: 39908.787
+  2015-05-22T19:45:21:  Absolute precision ........: 0.005
+  2015-05-22T19:45:21:  Acceptable value decrease .: 2
+  2015-05-22T19:45:21:  Optimization status .......: converged
+  2015-05-22T19:45:21:  Number of parameters ......: 10
+  2015-05-22T19:45:21:  Number of free parameters .: 4
+  2015-05-22T19:45:21:  Number of iterations ......: 2
+  2015-05-22T19:45:21:  Lambda ....................: 1e-05
+  2015-05-22T19:45:21:  Maximum log likelihood ....: -39908.787
+  2015-05-22T19:45:21:  Observed events  (Nobs) ...: 13360.000
+  2015-05-22T19:45:21:  Predicted events (Npred) ..: 13359.998 (Nobs - Npred = 0.00180799)
+  2015-05-22T19:45:21: === GModels ===
+  2015-05-22T19:45:21:  Number of models ..........: 2
+  2015-05-22T19:45:21:  Number of parameters ......: 10
+  2015-05-22T19:45:21: === GModelSky ===
+  2015-05-22T19:45:21:  Name ......................: Crab
+  2015-05-22T19:45:21:  Instruments ...............: all
+  2015-05-22T19:45:21:  Instrument scale factors ..: unity
+  2015-05-22T19:45:21:  Observation identifiers ...: all
+  2015-05-22T19:45:21:  Model type ................: PointSource
+  2015-05-22T19:45:21:  Model components ..........: "SkyDirFunction" * "PowerLaw" * "Constant"
+  2015-05-22T19:45:21:  Number of parameters ......: 6
+  2015-05-22T19:45:21:  Number of spatial par's ...: 2
+  2015-05-22T19:45:21:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
+  2015-05-22T19:45:21:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
+  2015-05-22T19:45:21:  Number of spectral par's ..: 3
+  2015-05-22T19:45:21:   Prefactor ................: 5.78197e-16 +/- 1.10889e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
+  2015-05-22T19:45:21:   Index ....................: -2.49174 +/- 0.0164369 [-0,-5]  (free,scale=-1,gradient)
+  2015-05-22T19:45:21:   PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T19:45:21:  Number of temporal par's ..: 1
+  2015-05-22T19:45:21:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+  2015-05-22T19:45:21: === GCTAModelIrfBackground ===
+  2015-05-22T19:45:21:  Name ......................: CTABackgroundModel
+  2015-05-22T19:45:21:  Instruments ...............: CTA
+  2015-05-22T19:45:21:  Instrument scale factors ..: unity
+  2015-05-22T19:45:21:  Observation identifiers ...: all
+  2015-05-22T19:45:21:  Model type ................: "PowerLaw" * "Constant"
+  2015-05-22T19:45:21:  Number of parameters ......: 4
+  2015-05-22T19:45:21:  Number of spectral par's ..: 3
+  2015-05-22T19:45:21:   Prefactor ................: 1.02048 +/- 0.0280684 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
+  2015-05-22T19:45:21:   Index ....................: 0.0197856 +/- 0.0156686 [-5,5]  (free,scale=1,gradient)
+  2015-05-22T19:45:21:   PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T19:45:21:  Number of temporal par's ..: 1
+  2015-05-22T19:45:21:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+
+The maximum likelihood optimizer required 2 iterations to converge. This
 is pretty fast, but recall that we used the same model file for the simulation
 and for fitting, hence the initial parameter values were already very close
 to the best fitting values. To see the impact of the initial parameters on
@@ -466,7 +459,7 @@ the fit result, you may re-run :ref:`ctlike` using another copy of the model
 XML file where you change the value attributes of the parameters that should be 
 fitted. You will see that the optimizer requires a couple of more iterations,
 but it should converge to the same solution (provided that the initial values
-are not too far of the best fitting values).
+are not too far from the best fitting values).
 .. note::
@@ -502,17 +495,17 @@ are not too far of the best fitting values).
    .. code-block:: xml
-      2015-02-05T08:29:07: === GModelSky ===
-      2015-02-05T08:29:07:  Name ......................: Crab
-      2015-02-05T08:29:07:  Instruments ...............: all
-      2015-02-05T08:29:07:  Test Statistic ............: 6875.35
+    2015-05-22T19:58:43: === GModelSky ===
+    2015-05-22T19:58:43:  Name ......................: Crab
+    2015-05-22T19:58:43:  Instruments ...............: all
+    2015-05-22T19:58:43:  Test Statistic ............: 18662.6
    The Test Statistic value will also be added as new attribute
    ``ts`` to the XML result file:
    .. code-block:: xml
-      <source name="Crab" type="PointSource" ts="6875.350" tscalc="1">
+    <source name="Crab" type="PointSource" ts="18662.576" tscalc="1">
 .. _sec_unbinned_cta:
@@ -521,8 +514,8 @@ Doing an unbinned analysis
 ~~~~~~~~~~~~~~~~~~~~~~~~~~
 As gamma-ray events are rare, the counts cubes generated by :ref:`ctbin`
-will in general be sparse, having many empty pixels, in particular at
-high energies.
+will in general be sparsly populated, having many empty pixels, in
+particular at high energies.
 An alternative analysis technique consists of working directly on the event
 list without binning the events in a counts cube. We will see the benefit of
 such an analysis later once you re-run :ref:`ctlike` in unbinned mode.
@@ -553,10 +546,10 @@ following example shows how to run ctselect:
   RA for ROI centre (degrees) (0-360) [83.63] 
   Dec for ROI centre (degrees) (-90-90) [22.01] 
   Radius of ROI (degrees) (0-180) [3.0] 
-  Start time (CTA MET in seconds) (0) [0.0] 
-  End time (CTA MET in seconds) (0) [0.0] 
-  Lower energy limit (TeV) (0) [0.1] 
-  Upper energy limit (TeV) (0) [100.0] 
+  Start time (CTA MET in seconds) [0.0] 
+  End time (CTA MET in seconds) [0.0] 
+  Lower energy limit (TeV) [0.1] 
+  Upper energy limit (TeV) [100.0] 
   Output event list or observation definition file [selected_events.fits]
 :ref:`ctselect` takes the input event list ``events.fits``, performs an
@@ -585,16 +578,16 @@ part:
 .. code-block:: xml
-  2015-02-04T17:21:41: +=================+
-  2015-02-04T17:21:41: | Event selection |
-  2015-02-04T17:21:41: +=================+
-  2015-02-04T17:21:41: === Observation ===
-  2015-02-04T17:21:41:  Selected energy range .....: 0.1 - 100 TeV
-  2015-02-04T17:21:41:  Requested ROI .............: Centre(RA,DEC)=(83.63, 22.01) deg, Radius=3 deg
-  2015-02-04T17:21:41:  ROI of data ...............: Centre(RA,DEC)=(83.63, 22.01) deg, Radius=5 deg
-  2015-02-04T17:21:41:  Selected ROI ..............: Centre(RA,DEC)=(83.63, 22.01) deg, Radius=3 deg
-  2015-02-04T17:21:41:  cfitsio selection .........: ENERGY >= 0.10000000 && ENERGY <= 100.00000000 && ANGSEP(83.630000,22.010000,RA,DEC) <= 3.000000
-  2015-02-04T17:21:41:  FITS filename .............: /var/tmp/tmp.0.pTIaTL[EVENTS][ENERGY >= 0.10000000 && ENERGY <= 100.00000000 && ANGSEP(83.630000,22.010000,RA,DEC) <= 3.000000]
+  2015-05-22T20:01:46: +=================+
+  2015-05-22T20:01:46: | Event selection |
+  2015-05-22T20:01:46: +=================+
+  2015-05-22T20:01:46: === Observation ===
+  2015-05-22T20:01:47:  Selected energy range .....: 0.1 - 100 TeV
+  2015-05-22T20:01:47:  Requested ROI .............: Centre(RA,DEC)=(83.63, 22.01) deg, Radius=3 deg
+  2015-05-22T20:01:47:  ROI of data ...............: Centre(RA,DEC)=(83.63, 22.01) deg, Radius=5 deg
+  2015-05-22T20:01:47:  Selected ROI ..............: Centre(RA,DEC)=(83.63, 22.01) deg, Radius=3 deg
+  2015-05-22T20:01:47:  cfitsio selection .........: ENERGY >= 0.10000000 && ENERGY <= 100.00000000 && ANGSEP(83.630000,22.010000,RA,DEC) <= 3.000000
+  2015-05-22T20:01:47:  FITS filename .............: /var/tmp/tmp.0.rRgLpO[EVENTS][ENERGY >= 0.10000000 && ENERGY <= 100.00000000 && ANGSEP(83.630000,22.010000,RA,DEC) <= 3.000000]
 .. note::
@@ -608,10 +601,10 @@ part:
    :ref:`ctselect` will automatically adjust
    the selection parameters to guarantee full enclosure. To keep track of
    this adjustment, the :ref:`ctselect` log file quotes the requested
-   selection, and existing selections, and the selection that was finally 
+   selection, any existing selections, and the selection that was finally 
    applied.
-.. note::
+.. warning::
    :ref:`ctselect` may of course also be used for event selection prior to
    binned analysis, for example to select events for a given period in
@@ -627,24 +620,25 @@ instead of the counts cube:
 .. code-block:: bash
-  $ ctlike
-  Event list, counts cube or observation definition file [cntmap.fits] selected_events.fits
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
-  Source model [$CTOOLS/share/models/crab.xml] 
-  Source model output file [crab_results.xml] 
+  $ ctlike 
+  Event list, counts cube or observation definition file [cntcube.fits] selected_events.fits
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
+  Source model [$CTOOLS/share/models/crab.xml]
+  Source model output file [crab_results.xml]
 You will recognise that :ref:`ctlike` runs much faster in unbinned mode
 compared to binned mode.
 This is understandable as the selected event list contains
-only 6127 events, while the binned counts cube we used before had 
+only 15312 events, while the binned counts cube we used before had 
 200 x 200 x 20 = 800000 pixels. As unbinned maximum likelihood fitting loops
 over the events (while binned maximum likelihood loops over the pixels),
 there are much less operations to perform in unbinned than in binned mode
 (there is some additional overhead in unbinned mode that comes from
 integrating the models over the region of interest, yet this is negligible
-compared to the operations needed when looping over the pixels). So as long
-as you work with short event lists, unbinned mode is faster.
+compared to the operations needed when looping over all pixels). So as long
+as you work with small event lists, unbinned mode is faster (this
+typically holds up to 100 hours of observing time).
 Unbinned :ref:`ctlike` should also be more precise as no binning is performed,
 hence there is no loss of information due to histogramming.
@@ -658,68 +652,62 @@ the same events for both analyses.
 .. code-block:: xml
-  2015-02-04T17:31:58: +=================================+
-  2015-02-04T17:31:58: | Maximum likelihood optimisation |
-  2015-02-04T17:31:58: +=================================+
-  2015-02-04T17:31:58:  >Iteration   0: -logL=44395.244, Lambda=1.0e-03
-  2015-02-04T17:31:58:   Iteration   1: -logL=44396.223, Lambda=1.0e-03, delta=-0.979, max(|grad|)=9.644743 [Sigma:6] (stalled)
-  2015-02-04T17:31:58:  >Iteration   2: -logL=44396.162, Lambda=1.0e-02, delta=0.061, max(|grad|)=-0.309683 [Normalization:7]
-  2015-02-04T17:31:59:  >Iteration   3: -logL=44396.162, Lambda=1.0e-03, delta=0.000, max(|grad|)=-0.008931 [Index:3]
-  2015-02-04T17:31:59: 
-  2015-02-04T17:31:59: +=========================================+
-  2015-02-04T17:31:59: | Maximum likelihood optimization results |
-  2015-02-04T17:31:59: +=========================================+
-  2015-02-04T17:31:59: === GOptimizerLM ===
-  2015-02-04T17:31:59:  Optimized function value ..: 44396.162
-  2015-02-04T17:31:59:  Absolute precision ........: 0.005
-  2015-02-04T17:31:59:  Acceptable value decrease .: 2
-  2015-02-04T17:31:59:  Optimization status .......: converged
-  2015-02-04T17:31:59:  Number of parameters ......: 9
-  2015-02-04T17:31:59:  Number of free parameters .: 4
-  2015-02-04T17:31:59:  Number of iterations ......: 3
-  2015-02-04T17:31:59:  Lambda ....................: 0.0001
-  2015-02-04T17:31:59:  Maximum log likelihood ....: -44396.162
-  2015-02-04T17:31:59:  Observed events  (Nobs) ...: 6127.000
-  2015-02-04T17:31:59:  Predicted events (Npred) ..: 6127.000 (Nobs - Npred = 0.000303534)
-  2015-02-04T17:31:59: === GModels ===
-  2015-02-04T17:31:59:  Number of models ..........: 2
-  2015-02-04T17:31:59:  Number of parameters ......: 9
-  2015-02-04T17:31:59: === GModelSky ===
-  2015-02-04T17:31:59:  Name ......................: Crab
-  2015-02-04T17:31:59:  Instruments ...............: all
-  2015-02-04T17:31:59:  Instrument scale factors ..: unity
-  2015-02-04T17:31:59:  Observation identifiers ...: all
-  2015-02-04T17:31:59:  Model type ................: PointSource
-  2015-02-04T17:31:59:  Model components ..........: "SkyDirFunction" * "PowerLaw" * "Constant"
-  2015-02-04T17:31:59:  Number of parameters ......: 6
-  2015-02-04T17:31:59:  Number of spatial par's ...: 2
-  2015-02-04T17:31:59:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
-  2015-02-04T17:31:59:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
-  2015-02-04T17:31:59:  Number of spectral par's ..: 3
-  2015-02-04T17:31:59:   Prefactor ................: 6.13325e-16 +/- 2.05738e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
-  2015-02-04T17:31:59:   Index ....................: -2.5057 +/- 0.0250818 [-0,-5]  (free,scale=-1,gradient)
-  2015-02-04T17:31:59:   PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
-  2015-02-04T17:31:59:  Number of temporal par's ..: 1
-  2015-02-04T17:31:59:   Constant .................: 1 (relative value) (fixed,scale=1,gradient)
-  2015-02-04T17:31:59: === GCTAModelRadialAcceptance ===
-  2015-02-04T17:31:59:  Name ......................: Background
-  2015-02-04T17:31:59:  Instruments ...............: CTA
-  2015-02-04T17:31:59:  Instrument scale factors ..: unity
-  2015-02-04T17:31:59:  Observation identifiers ...: all
-  2015-02-04T17:31:59:  Model type ................: "Gaussian" * "FileFunction" * "Constant"
-  2015-02-04T17:31:59:  Number of parameters ......: 3
-  2015-02-04T17:31:59:  Number of radial par's ....: 1
-  2015-02-04T17:31:59:   Sigma ....................: 3.04429 +/- 0.0330227 [0.01,10] deg2 (free,scale=1,gradient)
-  2015-02-04T17:31:59:  Number of spectral par's ..: 1
-  2015-02-04T17:31:59:   Normalization ............: 0.996767 +/- 0.0175227 [0,1000]  (free,scale=1,gradient)
-  2015-02-04T17:31:59:  Number of temporal par's ..: 1
-  2015-02-04T17:31:59:   Constant .................: 1 (relative value) (fixed,scale=1,gradient)
-  2015-02-04T17:31:59: 
-  2015-02-04T17:31:59: +==============+
-  2015-02-04T17:31:59: | Save results |
-  2015-02-04T17:31:59: +==============+
-  2015-02-04T17:31:59: 
-  2015-02-04T17:31:59: Application "ctlike" terminated after 8 wall clock seconds, consuming 0.21275 seconds of CPU time.
+  2015-05-22T20:05:54: +=================================+
+  2015-05-22T20:05:54: | Maximum likelihood optimisation |
+  2015-05-22T20:05:54: +=================================+
+  2015-05-22T20:05:54:  >Iteration   0: -logL=93854.543, Lambda=1.0e-03
+  2015-05-22T20:05:54:  >Iteration   1: -logL=93848.595, Lambda=1.0e-03, delta=5.948, max(|grad|)=22.480581 [Index:7]
+  2015-05-22T20:05:54:  >Iteration   2: -logL=93848.587, Lambda=1.0e-04, delta=0.008, max(|grad|)=-0.081846 [Index:3]
+  2015-05-22T20:05:54:  >Iteration   3: -logL=93848.587, Lambda=1.0e-05, delta=0.000, max(|grad|)=-0.003485 [Index:3]
+  2015-05-22T20:05:54: 
+  2015-05-22T20:05:54: +=========================================+
+  2015-05-22T20:05:54: | Maximum likelihood optimization results |
+  2015-05-22T20:05:54: +=========================================+
+  2015-05-22T20:05:54: === GOptimizerLM ===
+  2015-05-22T20:05:54:  Optimized function value ..: 93848.587
+  2015-05-22T20:05:54:  Absolute precision ........: 0.005
+  2015-05-22T20:05:54:  Acceptable value decrease .: 2
+  2015-05-22T20:05:54:  Optimization status .......: converged
+  2015-05-22T20:05:54:  Number of parameters ......: 10
+  2015-05-22T20:05:54:  Number of free parameters .: 4
+  2015-05-22T20:05:54:  Number of iterations ......: 3
+  2015-05-22T20:05:54:  Lambda ....................: 1e-06
+  2015-05-22T20:05:54:  Maximum log likelihood ....: -93848.587
+  2015-05-22T20:05:54:  Observed events  (Nobs) ...: 15312.000
+  2015-05-22T20:05:54:  Predicted events (Npred) ..: 15312.000 (Nobs - Npred = 3.03823e-06)
+  2015-05-22T20:05:54: === GModels ===
+  2015-05-22T20:05:54:  Number of models ..........: 2
+  2015-05-22T20:05:54:  Number of parameters ......: 10
+  2015-05-22T20:05:54: === GModelSky ===
+  2015-05-22T20:05:54:  Name ......................: Crab
+  2015-05-22T20:05:54:  Instruments ...............: all
+  2015-05-22T20:05:54:  Instrument scale factors ..: unity
+  2015-05-22T20:05:54:  Observation identifiers ...: all
+  2015-05-22T20:05:54:  Model type ................: PointSource
+  2015-05-22T20:05:54:  Model components ..........: "SkyDirFunction" * "PowerLaw" * "Constant"
+  2015-05-22T20:05:54:  Number of parameters ......: 6
+  2015-05-22T20:05:54:  Number of spatial par's ...: 2
+  2015-05-22T20:05:54:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
+  2015-05-22T20:05:54:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
+  2015-05-22T20:05:54:  Number of spectral par's ..: 3
+  2015-05-22T20:05:54:   Prefactor ................: 5.80619e-16 +/- 1.11149e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
+  2015-05-22T20:05:54:   Index ....................: -2.49282 +/- 0.0163467 [-0,-5]  (free,scale=-1,gradient)
+  2015-05-22T20:05:54:   PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T20:05:54:  Number of temporal par's ..: 1
+  2015-05-22T20:05:54:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+  2015-05-22T20:05:54: === GCTAModelIrfBackground ===
+  2015-05-22T20:05:54:  Name ......................: CTABackgroundModel
+  2015-05-22T20:05:54:  Instruments ...............: CTA
+  2015-05-22T20:05:54:  Instrument scale factors ..: unity
+  2015-05-22T20:05:54:  Observation identifiers ...: all
+  2015-05-22T20:05:54:  Model type ................: "PowerLaw" * "Constant"
+  2015-05-22T20:05:54:  Number of parameters ......: 4
+  2015-05-22T20:05:54:  Number of spectral par's ..: 3
+  2015-05-22T20:05:54:   Prefactor ................: 1.03511 +/- 0.024311 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
+  2015-05-22T20:05:54:   Index ....................: 0.0333176 +/- 0.0134386 [-5,5]  (free,scale=1,gradient)
+  2015-05-22T20:05:54:   PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
+  2015-05-22T20:05:54:  Number of temporal par's ..: 1
+  2015-05-22T20:05:54:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
 .. _sec_butterfly_cta:
@@ -739,14 +727,14 @@ line.
 .. code-block:: bash
   $ ctbutterfly
-  Input event list, cube or observation definition file [events.fits] 
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
-  Source model [$CTOOLS/share/models/crab.xml] crab_results.xml
+  Input event list, counts cube or observation definition file [events.fits] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
+  Source model [$CTOOLS/share/models/crab.xml] 
   Source of interest [Crab] 
   Start value for first energy bin in TeV [0.1] 
   Stop value for last energy bin in TeV [100.0] 
-  Output ascii file [butterfly.txt]
+  Output ascii file [butterfly.txt] 
 Now that you have computed the confidence band of the spectral fit and 
 that you will have an ASCII file named ``butterfly.txt`` on disk you can
@@ -3,20 +3,6 @@
 CTA Instrument Response Functions
 ---------------------------------
-.. note ::
-
-   Instrument response functions are formally not part of ctools as they
-   should be provided by the instrument teams. The GammaLib framework on
-   which ctools are built comes with a set of basic response functions, but
-   for getting the latest instrument response functions you should contact
-   the relevant instrument teams. CTA consortium members should 
-   `download 
-   <https://portal.cta-observatory.org/WG/DM/DM_wiki/DATA_Access/Pages/Science%20Tools.aspx>`_
-   the latest calibration database containing Prod1 and Prod2 instrument
-   response functions from the consortium Sharepoint site (requires 
-   password).
-
-
 What are instrument response functions?
 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
@@ -47,98 +33,64 @@ atmospheric conditions, etc.). All these quantities and hence the
 instrument response function may depend on time.
-Data formats
-~~~~~~~~~~~~
-
-So far the data format of the CTA instrument response functions has not
-yet been settled, and information about instrument response functions is
-provided in various ways. Three types of data formats are supported:
-
--  :ref:`sec_cta_rsp_perftable`
-
--  :ref:`sec_cta_rsp_xspec`
-
--  :ref:`sec_cta_rsp_rsptable`
-
-All formats provide the instrument response as function of true (and 
-sometimes also measured) photon energy, typically from about 20 GeV to
-about 125 TeV.
-Performance tables provide the instrument response for on-axis sources.
-ARF, RMF and PSF files provide the instrument response for a specific
-source position (and region) within the field of view.
-Response tables provide the instrument response as function of position in 
-the field of view. Response tables are thus the most universal form of 
-instrument response functions available, and we recommend using this form 
-for any analysis.
-
-.. note ::
-
-   Instrument response functions are so far only provided for a fixed
-   zenith angle of 20 deg. Therefore, no zenith or azimuth angle dependence has 
-   been implemented so far.
-
-
-Installing the CTA calibration database
-~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
-After `downloading 
-<https://portal.cta-observatory.org/WG/DM/DM_wiki/DATA_Access/Pages/Science%20Tools.aspx>`_
-the latest calibration database (only possible for CTA consortium members),
-the database is installed using
-
-.. code-block:: bash
-
-  $ [sudo] tar -C $CALDB -zxvf cta-caldb-20150220.tar.gz
-
-or
-
-.. code-block:: bash
-
-  $ [sudo] tar -C $CALDB -xvf cta-caldb-20150220.tar
-
-depending on whether the database has been retrieved as a gzipped file or 
-not. Depending on your access rights to the ``$CALDB`` directory,
-installation of the calibration database may require root privileges
-(type ``sudo`` in this case, otherwise omit this part of the command).
-
-The calibration database must be installed into the directory to which you 
-``CALDB`` environment variable points. By default, ctools sets this 
-enviroment variable to
-
-.. code-block:: bash
-
-  $CTOOLS/share/caldb
-
-but it may be that due to some other analysis software on your system the
-``CALDB`` environment variable points to another directory. If in doubt, type
-
-.. code-block:: bash
-
-  echo $CALDB
-
-to find out to which directory your ``CALDB`` environment variable points. 
-You may always overwrite this setting using
-
-.. code-block:: bash
-
-  export CALDB=/my/preferred/caldb/directory
-
-which you can add to your ``.bashrc`` file to set the directory 
-permanently.
-ctools will use the ``CALDB`` environment variable to find out where the 
-calibration database is located. If needed, this mechanism can be 
-circumvented (see :ref:`sec_cta_rsp_abspath`).
-
-
-Specifying the CTA Instrument Response Functions
-~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+CTA response functions
+~~~~~~~~~~~~~~~~~~~~~~
+
+The instrument response functions for CTA are factorised into 
+the effective area :math:`A_{\rm eff}(d, p, E, t)` (units :math:`cm^2`),
+the point spread function :math:`PSF(p' | d, p, E, t)`,
+and the energy dispersion :math:`E_{\rm disp}(E' | d, p, E, t)`
+following:
+
+.. math::
+    R(p', E', t' | d, p, E, t) =
+    A_{\rm eff}(d, p, E, t) \times
+    PSF(p' | d, p, E, t) \times
+    E_{\rm disp}(E' | d, p, E, t)
+
+ctools are shipped with response functions for the northern and southern
+arrays, and variants are available that have been optimised for exposure
+times of 0.5 hours, 5 hours and 50 hours.
+In total, the following six instrument response functions are available:
+``North_0.5h``, ``North_5h``, ``North_50h``, ``South_0.5h``,
+``South_5h``, and ``South_50h``.
+
+Each response is stored in a single FITS file, and each component of
+the response factorisation is stored in a binary table of that FITS
+file.
+In addition, the response files contain an additional table that
+describes the background rate as function of energy and position in
+the field of view.
+An example of a CTA response file is shown below:
+
+.. figure:: irf-file.png
+   :width: 100%
+
+Each table in the response file is in a standardised format that is
+the one that is also used for the Fermi/LAT telescope.
+As an example, the effective area component of the response file
+is shown below.
+Response information is stored in a n-dimensional cube, and each axis
+of this cube is described by the lower and upper edges of the axis bins.
+In this example the effective area is stored as a 2D matrix with the
+first axis being energy and the second axis being offaxis angle.
+Effective area information is stored for true (``EFFAREA``) and
+reconstructed (``EFFAREA_RECO``) energy.
+Vector columns are used to store all information.
+
+.. figure:: irf-aeff.png
+   :width: 100%
+
+
+Specifying CTA response functions
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 The specification of the CTA Instrument Response Functions depends on the 
 way how ctools are used. Common to all methods is that the IRFs are 
 defined by a response name and a calibration database name.
-
-There are different means to specify the CTA Instrument Response Functions 
-when using ctools, and the following section describe the
+ctools makes use of HEASARC's CALDB format to index and store
+IRFs, and specification of the database and response names is
+sufficient to access the response.
 Specifying the response function as input parameters
@@ -146,23 +98,22 @@ Specifying the response function as input parameters
 ctools that require instrument response functions have two parameters
 to specify the calibration database name and the response function name.
-The following example shows a ``ctobssim`` run using the ``dummy``
-calibration database and the ``cta_dummy_irf`` response function:
+The following example shows a ``ctobssim`` run using the ``prod2``
+calibration database and the ``South_50h`` response function:
 .. code-block:: bash
   $ ctobssim
-  Model [$CTOOLS/share/models/crab.xml] 
-  Calibration database [dummy] 
-  Instrument response function [cta_dummy_irf] 
   RA of pointing (degrees) (0-360) [83.63] 
-  Dec of pointing (degrees) (-90-90) [22.01] 
+  Dec of pointing (degrees) (-90-90) [22.51] 
   Radius of FOV (degrees) (0-180) [5.0] 
-  Start time (MET in s) (0) [0.0] 
-  End time (MET in s) (0) [1800.0] 
-  Lower energy limit (TeV) (0) [0.1] 
-  Upper energy limit (TeV) (0) [100.0] 
-  Output event data file or observation definition file [events.fits]
+  Start time (MET in s) [0.0] 
+  End time (MET in s) [1800.0] 
+  Lower energy limit (TeV) [0.1] 
+  Upper energy limit (TeV) [100.0] 
+  Calibration database [prod2] 
+  Instrument response function [South_50h] 
+  Model [$CTOOLS/share/models/crab.xml] 
 Running the other tools is equivalent.
@@ -180,7 +131,7 @@ definition file:
   <observation_list title="observation library">
     <observation name="Crab" id="00001" instrument="CTA">
       <parameter name="EventList"   file="events.fits"/>
-      <parameter name="Calibration" database="dummy" response="cta_dummy_irf"/>
+      <parameter name="Calibration" database="prod2" response="South_50h"/>
     </observation>
   </observation_list>
@@ -190,9 +141,9 @@ response name. You can then pass this file directly to, e.g., ``ctlike``:
 .. code-block:: bash
   $ ctlike
-  Event list, counts map or observation definition file [events.fits] obs_rsp.xml
+  Event list, counts cube or observation definition file [events.fits] obs_irf.xml
   Source model [$CTOOLS/share/models/crab.xml] 
-  Source model output file [crab_results.xml]
+  Source model output file [crab_results.xml] 
 Note that ``ctlike`` does not ask for the calibration database and
 response name as it found the relevant information in the XML file.
@@ -211,10 +162,10 @@ specify them individually in the XML observation file as follows:
   <observation_list title="observation library">
     <observation name="Crab" id="00001" instrument="CTA">
       <parameter name="EventList"           file="events.fits"/>
-      <parameter name="EffectiveArea"       file="$CALDB/data/cta/dummy/bcf/cta_dummy_irf.dat"/>
-      <parameter name="PointSpreadFunction" file="$CALDB/data/cta/dummy/bcf/cta_dummy_irf.dat"/>
-      <parameter name="EnergyDispersion"    file="$CALDB/data/cta/dummy/bcf/cta_dummy_irf.dat"/>
-      <parameter name="Background"          file="$CALDB/data/cta/dummy/bcf/cta_dummy_irf.dat"/>
+      <parameter name="EffectiveArea"       file="$CALDB/data/cta/prod2/bcf/North_0.5h/irf_file.fits.gz"/>
+      <parameter name="PointSpreadFunction" file="$CALDB/data/cta/prod2/bcf/North_0.5h/irf_file.fits.gz"/>
+      <parameter name="EnergyDispersion"    file="$CALDB/data/cta/prod2/bcf/North_0.5h/irf_file.fits.gz"/>
+      <parameter name="Background"          file="$CALDB/data/cta/prod2/bcf/North_0.5h/irf_file.fits.gz"/>
     </observation>
   </observation_list>
@@ -229,169 +180,12 @@ and response name from within Python:
   import gammalib
   obs   = gammalib.GCTAObservation()
-  caldb = gammalib.GCaldb("cta", "dummy")
-  irf   = "cta_dummy_irf"
+  caldb = gammalib.GCaldb("cta", "prod2")
+  irf   = "South_50h"
   obs.response(irf, caldb)
 The calibration database is set by creating a ``GCaldb`` object. The
 constructor takes as argument the mission (always ``cta``) and the 
-database name, in our case ``dummy``. The response function is then set
-by passing the response name (here ``cta_dummy_irf``) and the calibration
+database name, in our case ``prod2``. The response function is then set
+by passing the response name (here ``South_50h``) and the calibration
 database object to the ``response`` method.
-
-Alternatively, you can specify a specific response file by using
-
-.. code-block:: python
-
-  import gammalib
-  obs   = gammalib.GCTAObservation()
-  caldb = gammalib.GCaldb("$CALDB/data/cta/dummy/bcf")
-  irf   = "cta_dummy_irf.dat"
-  obs.response(irf, caldb)
-
-Here, the directory in which the calibration file resides is specified
-as argument to the ``GCaldb`` constructor, and the filename is passed
-as response function name to the ``response`` method.
-
-
-What response functions are available?
-~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
-To see what response functions are available on your system, type
-
-.. code-block:: bash
-
-  $ cscaldb
-
-This creates the file ``cscaldb.log`` in your current working directory
-that contains information about the available response functions. The 
-format of the file is:
-
-.. code-block:: bash
-
-  2014-02-20T17:15:21: +==============+
-  2014-02-20T17:15:21: | Mission: cta |
-  2014-02-20T17:15:21: +==============+
-  2014-02-20T17:15:21: === Calibration: aar ===
-  2014-02-20T17:15:21: DESY20140105_50h
-  2014-02-20T17:15:21: DESY20140105_50h_0deg
-  2014-02-20T17:15:21: DESY20140105_50h_180deg
-  2014-02-20T17:15:21: 
-  2014-02-20T17:15:21: === Calibration: aar500 ===
-  2014-02-20T17:15:21: DESY20140105_50h
-  2014-02-20T17:15:21: DESY20140105_50h_0deg
-  2014-02-20T17:15:21: DESY20140105_50h_180deg
-
-The name after ``Calibration`` is what we called earlier the calibration
-database, and the names that are following each 
-``=== Calibration: xxx ===`` header are the response file names.
-
-To get the same information also onto the console, type
-
-.. code-block:: bash
-
-  $ cscaldb debug=yes
-
-
-Data format details
-~~~~~~~~~~~~~~~~~~~
-
-.. _sec_cta_rsp_perftable:
-
-Performance tables
-^^^^^^^^^^^^^^^^^^
-
-In the early days, the instrument performances derived from Monte-Carlo
-simulations have been summarised in what we call here Performance 
-Tables, which are ASCII files that contain as function of energy the
-on-axis performance parameters of CTA, such as effective area, point spread
-function containment radius, energy resolution, background count rate and
-differential sensitivity.
-
-Below an example of a CTA performance table::
-
-  log(E)     Area     r68     r80  ERes. BG Rate    Diff Sens
-  -1.7      261.6  0.3621  0.4908 0.5134 0.0189924  6.88237e-11
-  -1.5     5458.2  0.2712  0.3685 0.4129 0.1009715  1.72717e-11
-  -1.3    15590.0  0.1662  0.2103 0.2721 0.0575623  6.16963e-12
-  -1.1    26554.1  0.1253  0.1567 0.2611 0.0213008  2.89932e-12
-  -0.9    52100.5  0.1048  0.1305 0.1987 0.0088729  1.39764e-12
-  -0.7    66132.1  0.0827  0.1024 0.1698 0.0010976  6.03531e-13
-  -0.5   108656.8  0.0703  0.0867 0.1506 0.0004843  3.98147e-13
-  -0.3   129833.0  0.0585  0.0722 0.1338 0.0001575  3.23090e-13
-  -0.1   284604.3  0.0531  0.0656 0.1008 0.0001367  2.20178e-13
-   0.1   263175.3  0.0410  0.0506 0.0831 0.0000210  1.87452e-13
-   0.3   778048.6  0.0470  0.0591 0.0842 0.0000692  1.53976e-13
-   0.5   929818.8  0.0391  0.0492 0.0650 0.0000146  1.18947e-13
-   0.7  1078450.0  0.0335  0.0415 0.0541 0.0000116  1.51927e-13
-   0.9  1448579.1  0.0317  0.0397 0.0516 0.0000047  1.42439e-13
-   1.1  1899905.0  0.0290  0.0372 0.0501 0.0000081  1.96670e-13
-   1.3  2476403.8  0.0285  0.0367 0.0538 0.0000059  2.20695e-13
-   1.5  2832570.6  0.0284  0.0372 0.0636 0.0000073  3.22523e-13
-   1.7  3534065.3  0.0290  0.0386 0.0731 0.0000135  4.84153e-13
-   1.9  3250103.4  0.0238  0.0308 0.0729 0.0000044  6.26265e-13
-   2.1  3916071.6  0.0260  0.0354 0.0908 0.0000023  7.69921e-13
-   ---------------------------------------------
-   1) log(E) = log10(E/TeV) - bin centre
-   2) Eff Area - in square metres after background cut (no theta cut)
-   3) Ang. Res - 68% containment radius of gamma-ray PSF post cuts - in degrees
-   4) Ang. Res - 80% containment radius of gamma-ray PSF post cuts - in degrees
-   5) Fractional Energy Resolution (rms)
-   6) BG Rate  - inside point-source selection region - post call cuts - in Hz
-   7) Diff Sens - differential sensitivity for this bin expressed as E^2 dN/dE
-      - in erg cm^-2 s^-1 - for a 50 hours exposure - 5 sigma significance including
-      systematics and statistics and at least 10 photons.
-
-
-.. _sec_cta_rsp_xspec:
-
-ARF, RMF and PSF files
-^^^^^^^^^^^^^^^^^^^^^^
-
-Instrument response information for the first CTA Data Challenge (1DC) has been
-provided in a format that was heavily inspired by the 
-`RMF and ARF file formats 
-<http://heasarc.gsfc.nasa.gov/docs/heasarc/caldb/docs/memos/cal_gen_92_002/cal_gen_92_002.html>`_
-that have been introduced by
-`HEASARC <http://heasarc.gsfc.nasa.gov/>`_
-for spectral analysis of X-ray data.
-
-The Redistribution Matrix File (RMF) describes how an incoming photon with 
-a given true energy is redistributed in measured energy. In other words, it
-describes the energy dispersion of the instrument. The RMF is organised as
-a two-dimensional matrix, with the first axis being in true energies while
-the second axis represents the measured energies.
-
-The Ancillary Response File (ARF) describes the sensitivity of the 
-instrument to photons of a given true energy. The ARF gives the effective
-area of the detector system after applying an event selection cut
-(theta cut). For computing the ARF, some knowledge of the Point Spread
-Function (PSF) is needed, so that the fraction of photons that falls into
-the event selection region for a given source can be properly estimated.
-Once this is done, no further PSF information is needed for the analysis.
-
-ctools however relies on the full modelling of the instrument, including 
-the Point Spread Function, and hence, PSF information has also been 
-provided for 1DC. Two data formats have been used for this: a first that 
-is based on a simple one-dimensional vector, providing the width of a
-2-dimensional Gaussian as function of energy; and a second that is based
-on a three-component 2D Gaussian as function of energy and offset angle
-in the camera system. While the former was implement by a simply FITS 
-table with two columns, the latter was implemented by
-:ref:`sec_cta_rsp_rsptable`.
-
-
-.. _sec_cta_rsp_rsptable:
-
-Response tables
-^^^^^^^^^^^^^^^
-
-The CTA response table class ``GCTAResponseTable`` provides a generic 
-handle for multi-dimensional response information. It is based on the 
-response format used for storing response information for the
-*Fermi*/LAT telescope. In this format, all information is stored in
-a single row of a FITS binary table. Each element of the row contains
-a vector column, that describes the axes of the  multi-dimensional response
-cube and the response information. Note that this class may in the future
-be promoted to the GammaLib core, as a similar class has been implemented
-in the *Fermi*/LAT interface. 
@@ -9,6 +9,13 @@ benchmark_cta_analysis.py
 check_models.py
   Checks a variety of models using an unbinned analysis.
+make_pointings.py
+  This script creates pointings for CTA surveys. For the moment the
+  Galactic Plane Survey (gps), the Extragalactic Survey (extgal),
+  the Galactic Centre Survey (gc) and the LMC survey (lmc) are
+  implemented. The script outputs and ASCII file that can be used by
+  csobsdef.
+
 make_pull_at_sensitivity_limit.py
   This script creates pull distributions as function of energy
   for a test source that is at the 5 sigma sensitivity limit for
@@ -39,8 +39,8 @@ def unbinned_pipeline(duration):
     """
     # Set script parameters
     model_name  = "${CTOOLS}/share/models/crab.xml"
-    caldb       = "dummy"
-    irf         = "cta_dummy_irf"
+    caldb       = "prod2"
+    irf         = "South_50h"
     ra          =   83.63
     dec         =   22.01
     rad_sim     =   10.0
@@ -103,8 +103,8 @@ def binned_pipeline(duration):
     """
     # Set script parameters
     model_name  = "${CTOOLS}/share/models/crab.xml"
-    caldb       = "dummy"
-    irf         = "cta_dummy_irf"
+    caldb       = "prod2"
+    irf         = "South_50h"
     ra          =   83.63
     dec         =   22.01
     rad_sim     =   10.0
@@ -167,17 +167,17 @@ def binned_pipeline(duration):
     return telapsed, tctlike
-# ============================ #
-# Cube-style analysis pipeline #
-# ============================ #
-def cube_pipeline(duration):
+# ========================= #
+# Stacked analysis pipeline #
+# ========================= #
+def stacked_pipeline(duration):
     """
     Cube-style analysis pipeline.
     """
     # Set script parameters
     model_name  = "${CTOOLS}/share/models/crab.xml"
-    caldb       = "dummy"
-    irf         = "cta_dummy_irf"
+    caldb       = "prod2"
+    irf         = "South_50h"
     ra          =   83.63
     dec         =   22.01
     rad_sim     =   10.0
@@ -260,17 +260,34 @@ def cube_pipeline(duration):
     psfcube["yref"].real(dec)
     psfcube.run()
-    # Set exposure and PSF cube for first CTA observation
-    obs = bin.obs()
-    cta = gammalib.GCTAObservation(obs[0])
-    cta.response(expcube.expcube(), psfcube.psfcube())
-    obs[0] = cta
+    # Create background cube
+    bkgcube = ctools.ctbkgcube(sim.obs())
+    bkgcube["incube"].filename("NONE")
+    bkgcube["ebinalg"].string("LOG")
+    bkgcube["emin"].real(emin)
+    bkgcube["emax"].real(emax)
+    bkgcube["enumbins"].integer(enumbins)
+    bkgcube["nxpix"].integer(10)
+    bkgcube["nypix"].integer(10)
+    bkgcube["binsz"].real(1.0)
+    bkgcube["coordsys"].string(coordsys)
+    bkgcube["proj"].string(proj)
+    bkgcube["xref"].real(ra)
+    bkgcube["yref"].real(dec)
+    bkgcube.run()
+
+    # Attach background model to observation container
+    bin.obs().models(bkgcube.models())
+
+    # Set Exposure and Psf cube for first CTA observation
+    # (ctbin will create an observation with a single container)
+    bin.obs()[0].response(expcube.expcube(), psfcube.psfcube(), bkgcube.bkgcube())
     # Get ctlike start CPU time
     tctlike = time.clock()
     # Perform maximum likelihood fitting
-    like = ctools.ctlike(obs)
+    like = ctools.ctlike(bin.obs())
     like.run()
     # Get stop CPU time
@@ -292,16 +309,16 @@ if __name__ == &#39;__main__&#39;:
     print("*************************************")
     print("*      CTA analysis benchmark       *")
     print("*************************************")
-    print("  Duration   Unbinned     Binned Cube-style  ctlike-UB   ctlike-B  ctlike-CS")
+    print("  Duration   Unbinned     Binned    Stacked  ctlike-UB   ctlike-B   ctlike-S")
     # Initialize arrays
     a_duration = []
     a_unbinned = []
     a_binned   = []
-    a_cube     = []
+    a_stacked  = []
     c_unbinned = []
     c_binned   = []
-    c_cube     = []
+    c_stacked  = []
     # Perform benchmarks
     duration = 1800.0
@@ -309,22 +326,22 @@ if __name__ == &#39;__main__&#39;:
         # Perform analyses
         t_unbinned, ct_unbinned = unbinned_pipeline(duration)
-        t_binned, ct_binned     = binned_pipeline(duration)
-        t_cube, ct_cube         = cube_pipeline(duration)
+        t_binned,   ct_binned   = binned_pipeline(duration)
+        t_stacked,  ct_stacked  = stacked_pipeline(duration)
         # Collect results
         a_duration.append(duration/3600.0)
         a_unbinned.append(t_unbinned)
         a_binned.append(t_binned)
-        a_cube.append(t_cube)
+        a_stacked.append(t_stacked)
         c_unbinned.append(ct_unbinned)
         c_binned.append(ct_binned)
-        c_cube.append(ct_cube)
+        c_stacked.append(ct_stacked)
         # Print results
         print("%10.0f %10.4f %10.4f %10.4f %10.4f %10.4f %10.4f" % \
-              (duration, t_unbinned, t_binned, t_cube, \
-                         ct_unbinned, ct_binned, ct_cube))
+              (duration, t_unbinned, t_binned, t_stacked, \
+                         ct_unbinned, ct_binned, ct_stacked))
         # Increment duration
         duration *= 2.0
@@ -335,10 +352,10 @@ if __name__ == &#39;__main__&#39;:
         plt.title("CTA analysis benchmark")
         plt.loglog(a_duration, a_unbinned, 'ro-', label='unbinned')
         plt.loglog(a_duration, a_binned, 'bo-', label='binned')
-        plt.loglog(a_duration, a_cube, 'go-', label='cube-style')
+        plt.loglog(a_duration, a_stacked, 'go-', label='cube-style')
         plt.loglog(a_duration, c_unbinned, 'ro--')
         plt.loglog(a_duration, c_binned, 'bo--')
-        plt.loglog(a_duration, c_cube, 'go--')
+        plt.loglog(a_duration, c_stacked, 'go--')
         plt.xlabel("Duration (hours)")
         plt.ylabel("CPU time (seconds)")
         plt.legend(loc="lower right")
@@ -41,8 +41,8 @@ def pipeline(model_name):
 	FITS files on disk. All data is only hold in memory.
 	"""
 	# Set script parameters
-	caldb       = "dummy"
-	irf         = "cta_dummy_irf"
+	caldb       = "prod2"
+	irf         = "South_50h"
 	ra          =   83.63
 	dec         =   22.01
 	rad_sim     =   10.0
 #! /usr/bin/env python
 # ==========================================================================
-# This script generated pointing patterns for CTA observations.
+# This script generates pointing patterns for CTA observations.
 #
 # Usage:
 #   ./make_pointings.py
@@ -124,19 +124,18 @@ def set_patch(lmin=-30.0, lmax=30.0, bmin=-1.5, bmax=+1.5, \
                 obs_site = site
             # Set site IRF
+            caldb = "prod2"
             if obs_site == "North":
-                caldb = "tenerife"
                 if lst:
-                    irf = "DESY20140105_50h"
+                    irf = "North_50h"
                 else:
-                    irf = "DESY20140105_50h_noLST"
+                    irf = "North_50h"
                 exposure_north += duration
             else:
-                caldb = "aar"
                 if lst:
-                    irf = "DESY20140105_50h"
+                    irf = "South_50h"
                 else:
-                    irf = "DESY20140105_50h_noLST"
+                    irf = "South_50h"
                 exposure_south += duration
             # Set IRF information
@@ -241,6 +240,59 @@ def set_gc(lst=True):
     return obsdef
+# ============= #
+# Setup LMC KSP #
+# ============= #
+def set_lmc(hours=250.0, lst=True):
+    """
+    Setup LMC KSP.
+    """
+    # Initialise observation definition
+    obsdef = []
+
+    # Set LMC centre
+    centre = gammalib.GSkyDir()
+    centre.radec_deg(80.0, -69.0)
+
+    # Set offset angle and number of pointings
+    offset = 1.0 # degrees
+    n_pnt  = 12
+
+    # Prepare computations
+    dphi     = 360.0/n_pnt
+    duration = hours*3600.0/float(n_pnt)
+
+    # Loop over pointings
+    for ipnt in range(n_pnt):
+
+        # Compute pointing direction
+        pnt = centre.copy()
+        pnt.rotate_deg(ipnt*dphi, offset)
+        lon = pnt.l_deg()
+        lat = pnt.b_deg()
+
+        # Set positions and duration
+        obs = {'lon': lon, 'lat': lat, 'duration': duration}
+
+        # Add IRF
+        caldb = "prod2"
+        if lst:
+            irf = "South_50h"
+        else:
+            irf = "South_50h"
+        obs['caldb'] = caldb
+        obs['irf']   = irf
+
+        # Append observation
+        obsdef.append(obs)
+
+    # Dump statistics
+    print("Number of pointings: %d (%.2f sec)" % (n_pnt,duration))
+
+    # Return observation definition
+    return obsdef
+
+
 # ======================================= #
 # Write observation definition dictionary #
 # ======================================= #
@@ -291,7 +343,6 @@ if __name__ == &#39;__main__&#39;:
     need_help = False
     # Test for command line arguments
-    #print(sys.argv[0])
     if (len(sys.argv) > 1):
         if sys.argv[1] == "-h":
             need_help = True
@@ -304,13 +355,14 @@ if __name__ == &#39;__main__&#39;:
     if need_help:
         print("Usage: make_pointing.py [OPTIONS]")
         print("     -h       Display this usage message")
-        print("     gps      Galactic plane scan (2 row scheme)")
-        print("     gps3     Galactic plane scan (3 row scheme)")
+        print("     gps      Galactic plane survey (2 row scheme)")
+        print("     gps3     Galactic plane survey (3 row scheme)")
         print("     extgal   Extragalactic survey")
-        print("     gc       Galactic centre")
+        print("     gc       Galactic centre survey")
+        print("     lmc      LMC survey")
         sys.exit()
-    # Galactic plane scan
+    # Galactic plane survey
     if obsname == "gps":
         obsdef = set_gps(lst=True)
         write_obsdef("gps.dat", obsdef)
@@ -328,6 +380,11 @@ if __name__ == &#39;__main__&#39;:
         obsdef = set_gc(lst=True)
         write_obsdef("gc.dat", obsdef)
+    # LMC
+    elif obsname == "lmc":
+        obsdef = set_lmc(lst=True)
+        write_obsdef("lmc.dat", obsdef)
+
     # Invalid pattern
     else:
         print("Unknown option \""+obsname+"\"")
@@ -61,8 +61,8 @@ def create_pull(loge, emin, emax, models, ntrials=100, duration=180000.0, \
 	pull.models(models)
 	pull["outfile"]  = outfile
 	pull["ntrials"]  = ntrials
-	pull["caldb"]    = "${CALDB}/data/cta/dummy/bcf"
-	pull["irf"]      = "cta_dummy_irf"
+	pull["caldb"]    = "prod2"
+	pull["irf"]      = "South_50h"
 	pull["ra"]       = 83.6331
 	pull["dec"]      = 22.0145
 	pull["emin"]     = emin
@@ -36,8 +36,8 @@ def make_spectrum():
     """
     # Set script parameters
     model_name  = "${CTOOLS}/share/models/crab.xml"
-    caldb       = "dummy"
-    irf         = "cta_dummy_irf"
+    caldb       = "prod2"
+    irf         = "South_50h"
     ra          =   83.63
     dec         =   22.01
     rad_sim     =    3.0
@@ -88,7 +88,7 @@ def plot_spectrum(spectrum):
                spectrum['flux']['value'], 'ro', label='Crab')
     plt.errorbar(spectrum['energy']['value'], \
                  spectrum['flux']['value'], \
-                 spectrum['flux']['ed_value'], fmt=None, ecolor='r')
+                 spectrum['flux']['ed_value'], ecolor='r')
     # Put labels
     plt.xlabel("Energy ("+spectrum['energy']['unit']+")")
@@ -61,21 +61,23 @@ def create_ts(loge, emin, emax, ntrials=100, duration=180000.0, \
 	# Setup cstsdist tool
 	tsdist = cstsdist()
-	tsdist["outfile"]  = outfile
-	tsdist["ntrials"]  = ntrials
-	tsdist["caldb"]    = "${CALDB}/data/cta/dummy/bcf"
-	tsdist["irf"]      = "cta_dummy_irf"
-	tsdist["type"]     = "point"
-	tsdist["index"]    = -2.48
-	tsdist["offset"]   = 0.0
-	tsdist["bkg"]      = "${CTOOLS}/share/models/bkg_dummy.txt"
-	tsdist["emin"]     = float(emin)
-	tsdist["emax"]     = float(emax)
-	tsdist["enumbins"] = int(enumbins)
-	tsdist["duration"] = float(duration)
-	tsdist["rad"]      = 5.0
-	tsdist["npix"]     = 200
-	tsdist["binsz"]    = 0.05
+	tsdist["inmodel"].filename("$CTOOLS/share/models/crab.xml")
+	tsdist["srcname"].string("Crab")
+	tsdist["outfile"].filename(outfile)
+	tsdist["ntrials"].integer(ntrials)
+	tsdist["caldb"].string("prod2")
+	tsdist["irf"].string("South_50h")
+	tsdist["ra"].real(83.63)
+	tsdist["dec"].real(22.01)
+	tsdist["emin"].real(float(emin))
+	tsdist["emax"].real(float(emax))
+	tsdist["enumbins"].integer(int(enumbins))
+	tsdist["tmin"].real(0.0)
+	tsdist["tmax"].real(float(duration))
+	tsdist["rad"].real(5.0)
+	tsdist["npix"].integer(200)
+	tsdist["binsz"].real(0.05)
+	#tsdist["debug"].boolean(True)
 	# Optionally open the log file
 	if log:
@@ -139,14 +141,10 @@ if __name__ == &#39;__main__&#39;:
 	enumbins    = pars[1]['value']
 	duration    = pars[2]['value']
 	max_threads = pars[3]['value']
-	#print ntrials
-	#print enumbins
-	#print duration
-	#print max_threads
 	# Loop over energy bands. The energy bands are those that are also
 	# used for sensitivity computation.
-	for ieng in range(21):
+	for ieng in range(20):
 		# Set energies
 		loge  = -1.7 + ieng * 0.2
@@ -36,7 +36,7 @@ from ctools import obsutils
 def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, \
                        emin=0.1, emax=100.0, rad=5.0, duration=1800.0, \
                        deadc=0.95, \
-                       caldb="dummy", irf="cta_dummy_irf"):
+                       caldb="prod2", irf="South_50h"):
     """
     Returns an observation container.
@@ -96,7 +96,7 @@ def run_pipeline(obs, emin=0.1, emax=100.0, \
                  enumbins=20, nxpix=200, nypix=200, binsz=0.02, \
                  coordsys="CEL", proj="CAR", \
                  model="${CTOOLS}/share/models/crab.xml", \
-                 caldb="dummy", irf="cta_dummy_irf", \
+                 caldb="prod2", irf="South_50h", \
                  debug=False):
     """
     Simulation and binned analysis pipeline.
@@ -165,6 +165,7 @@ def run_pipeline(obs, emin=0.1, emax=100.0, \
     like["outmodel"].filename("fit_results.xml")
     like["expcube"].filename("NONE")
     like["psfcube"].filename("NONE")
+    like["bkgcube"].filename("NONE")
     like["caldb"].string(caldb)
     like["irf"].string(irf)
     like["debug"].boolean(True) # Switch this always on for results in console
@@ -36,7 +36,7 @@ from ctools import obsutils
 def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, \
                        emin=0.1, emax=100.0, rad=5.0, duration=1800.0, \
                        deadc=0.95, \
-                       caldb="dummy", irf="cta_dummy_irf"):
+                       caldb="prod2", irf="South_50h"):
     """
     Returns an observation container.
@@ -36,7 +36,7 @@ from ctools import obsutils
 def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, \
                        emin=0.1, emax=100.0, rad=5.0, duration=1800.0, \
                        deadc=0.95, \
-                       caldb="dummy", irf="cta_dummy_irf"):
+                       caldb="prod2", irf="South_50h"):
     """
     Returns an observation container.
@@ -96,7 +96,7 @@ def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0, \
                  enumbins=20, nxpix=200, nypix=200, binsz=0.02, \
                  coordsys="CEL", proj="CAR", \
                  model="${CTOOLS}/share/models/crab.xml", \
-                 caldb="dummy", irf="cta_dummy_irf", \
+                 caldb="prod2", irf="South_50h", \
                  debug=False):
     """
     Simulation and stacked analysis pipeline.
@@ -210,6 +210,7 @@ def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0, \
     like["outmodel"].filename("fit_results.xml")
     like["expcube"].filename("expcube.fits")
     like["psfcube"].filename("psfcube.fits")
+    like["bkgcube"].filename("bkgcube.fits")
     like["caldb"].string(caldb)
     like["irf"].string(irf)
     like["debug"].boolean(True) # Switch this always on for results in console
@@ -36,7 +36,7 @@ from ctools import obsutils
 def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, \
                        emin=0.1, emax=100.0, rad=5.0, duration=1800.0, \
                        deadc=0.95, \
-                       caldb="dummy", irf="cta_dummy_irf"):
+                       caldb="prod2", irf="South_50h"):
     """
     Returns an observation container.
@@ -188,7 +188,7 @@ def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0, \
     # Set Exposure and Psf cube for first CTA observation
     # (ctbin will create an observation with a single container)
-    bin.obs()[0].response(expcube.expcube(), psfcube.psfcube())
+    bin.obs()[0].response(expcube.expcube(), psfcube.psfcube(), bkgcube.bkgcube())
     # Perform maximum likelihood fitting
     like = ctools.ctlike(bin.obs())
@@ -36,7 +36,7 @@ from ctools import obsutils
 def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, \
                        emin=0.1, emax=100.0, rad=5.0, duration=1800.0, \
                        deadc=0.95, \
-                       caldb="dummy", irf="cta_dummy_irf"):
+                       caldb="prod2", irf="South_50h"):
     """
     Returns an observation container.
@@ -96,7 +96,7 @@ def run_pipeline(obs, ra=83.63, dec=22.01, rad=3.0, \
                  emin=0.1, emax=100.0, \
                  tmin=0.0, tmax=0.0, \
                  model="${CTOOLS}/share/models/crab.xml", \
-                 caldb="dummy", irf="cta_dummy_irf", \
+                 caldb="prod2", irf="South_50h", \
                  debug=False):
     """
     Simulation and unbinned analysis pipeline.
@@ -36,7 +36,7 @@ from ctools import obsutils
 def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, \
                        emin=0.1, emax=100.0, rad=5.0, duration=1800.0, \
                        deadc=0.95, \
-                       caldb="dummy", irf="cta_dummy_irf"):
+                       caldb="prod2", irf="South_50h"):
     """
     Returns an observation container.
@@ -38,7 +38,7 @@ except:
 # ======================== #
 def createobs(ra=86.171648, dec=-1.4774586, rad=5.0, \
               emin=0.1, emax=100.0, duration=360000.0, deadc=0.95, \
-              irf="cta_dummy_irf", caldb="$CTOOLS/share/caldb/data/cta/dummy"):
+              irf="South_50h", caldb="prod2"):
     """
     Create CTA observation.
     """
 ## Process this file with automake to produce Makefile.in
 # Background rate files to distribute
-bkgs = bkg_dummy.txt
-bkgsdir   = $(datadir)/models
-bkgs_DATA = $(bkgs)
+#bkgs = bkg_dummy.txt
+#bkgsdir   = $(datadir)/models
+#bkgs_DATA = $(bkgs)
 # Models to distribute
 models      = crab.xml \
@@ -16,4 +16,4 @@ modelsdir   = $(datadir)/models
 models_DATA = $(models)
 # Extra files to be distributed in the tar file
-EXTRA_DIST = $(bkgs) $(models)
+EXTRA_DIST = $(models)
-19952.6231497 0.00888770093889
-31622.7766017 0.0528862764556
-50118.7233627 0.0584088588708
-79432.8234724 0.0245627698931
-125892.541179 0.00930817022553
-199526.231497 0.00117990856387
-316227.766017 0.000458229904465
-501187.233627 0.000135628482432
-794328.234724 8.99472876298e-05
-1258925.41179 1.46319988679e-05
-1995262.31497 2.23452132224e-05
-3162277.66017 4.28526605833e-06
-5011872.33627 3.02154407423e-06
-7943282.34724 8.46588074886e-07
-12589254.1179 1.05217831413e-06
-19952623.1497 4.95410035794e-07
-31622776.6017 3.77171043065e-07
-50118723.3627 4.05250267135e-07
-79432823.4724 1.31996926243e-07
-125892541.179 3.23120422483e-08
@@ -11,12 +11,11 @@
       <parameter name="DEC" scale="1.0" value="22.0145" min="-90"  max="90"  free="0"/>
     </spatialModel>
   </source>
-  <source name="Background" type="RadialAcceptance" instrument="CTA">
-    <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-      <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="1.0" free="1"/>
+  <source name="CTABackgroundModel" type="CTAIrfBackground" instrument="CTA">
+    <spectrum type="PowerLaw">	
+      <parameter name="Prefactor" scale="1.0"  value="1.0"  min="1e-3" max="1e+3"   free="1"/>	
+      <parameter name="Index"     scale="1.0"  value="0.0"  min="-5.0" max="+5.0"   free="1"/>	
+      <parameter name="Scale"     scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>	
     </spectrum>
-    <radialModel type="Gaussian">
-       <parameter name="Sigma" scale="1.0" value="3.0" min="0.01" max="10.0" free="1"/>
-    </radialModel>
-  </source>
+  </source>	
 </source_library>
@@ -12,12 +12,11 @@
       <parameter name="Radius" scale="1.0" value="0.45"    min="0.01" max="10"  free="1"/>
     </spatialModel>
   </source>
-  <source name="Background" type="RadialAcceptance" instrument="CTA">
-    <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-      <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="1.0" free="1"/>
+  <source name="CTABackgroundModel" type="CTAIrfBackground" instrument="CTA">
+    <spectrum type="PowerLaw">	
+      <parameter name="Prefactor" scale="1.0"  value="1.0"  min="1e-3" max="1e+3"   free="1"/>	
+      <parameter name="Index"     scale="1.0"  value="0.0"  min="-5.0" max="+5.0"   free="1"/>	
+      <parameter name="Scale"     scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>	
     </spectrum>
-    <radialModel type="Gaussian">
-       <parameter name="Sigma" scale="1.0" value="3.0" min="0.01" max="10.0" free="1"/>
-    </radialModel>
-  </source>
+  </source>	
 </source_library>
@@ -12,12 +12,11 @@
       <parameter name="Sigma" scale="1.0" value="0.20"    min="0.01" max="10"  free="1"/>
     </spatialModel>
   </source>
-  <source name="Background" type="RadialAcceptance" instrument="CTA">
-    <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-      <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="1.0" free="1"/>
+  <source name="CTABackgroundModel" type="CTAIrfBackground" instrument="CTA">
+    <spectrum type="PowerLaw">	
+      <parameter name="Prefactor" scale="1.0"  value="1.0"  min="1e-3" max="1e+3"   free="1"/>	
+      <parameter name="Index"     scale="1.0"  value="0.0"  min="-5.0" max="+5.0"   free="1"/>	
+      <parameter name="Scale"     scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>	
     </spectrum>
-    <radialModel type="Gaussian">
-       <parameter name="Sigma" scale="1.0" value="3.0" min="0.01" max="10.0" free="1"/>
-    </radialModel>
-  </source>
+  </source>	
 </source_library>
@@ -13,12 +13,11 @@
       <parameter name="Width"  scale="1.0" value="0.10"    min="0.01" max="10"  free="1"/>
     </spatialModel>
   </source>
-  <source name="Background" type="RadialAcceptance" instrument="CTA">
-    <spectrum type="FileFunction" file="$CTOOLS/share/models/bkg_dummy.txt">
-      <parameter scale="1.0" name="Normalization" min="0.0" max="1000.0" value="1.0" free="1"/>
+  <source name="CTABackgroundModel" type="CTAIrfBackground" instrument="CTA">
+    <spectrum type="PowerLaw">	
+      <parameter name="Prefactor" scale="1.0"  value="1.0"  min="1e-3" max="1e+3"   free="1"/>	
+      <parameter name="Index"     scale="1.0"  value="0.0"  min="-5.0" max="+5.0"   free="1"/>	
+      <parameter name="Scale"     scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>	
     </spectrum>
-    <radialModel type="Gaussian">
-       <parameter name="Sigma" scale="1.0" value="3.0" min="0.01" max="10.0" free="1"/>
-    </radialModel>
-  </source>
+  </source>	
 </source_library>
@@ -47,7 +47,7 @@ public:
     void                      clear(void);
     void                      run(void);
     void                      save(void);
-    const GCTACubeBackground& cube(void) const;
+    const GCTACubeBackground& bkgcube(void) const;
     const GModels&            models(void) const;
 };
@@ -288,7 +288,7 @@ def modmap(obs, eref=0.1, proj=&quot;TAN&quot;, coord=&quot;GAL&quot;, xval=0.0, yval=0.0, \
 # ===================================== #
 def set(pntdir, tstart=0.0, duration=1800.0, deadc=0.95, \
         emin=0.1, emax=100.0, rad=5.0, \
-        irf="cta_dummy_irf", caldb="dummy"):
+        irf="South_50h", caldb="prod2"):
     """
     Obsolete function, use set_obs instead.
     """
@@ -309,7 +309,7 @@ def set(pntdir, tstart=0.0, duration=1800.0, deadc=0.95, \
 # ======================= #
 def set_obs(pntdir, tstart=0.0, duration=1800.0, deadc=0.95, \
             emin=0.1, emax=100.0, rad=5.0, \
-            irf="cta_dummy_irf", caldb="dummy", id="000000", instrument="CTA"):
+            irf="South_50h", caldb="prod2", id="000000", instrument="CTA"):
     """
     Returns a single CTA observation containing an empty CTA event list.
     By looping over this function you can add CTA observations to the
@@ -384,7 +384,7 @@ def set_obs(pntdir, tstart=0.0, duration=1800.0, deadc=0.95, \
 # ============================ #
 def set_obs_list(obsdeflist, tstart=0.0, duration=1800.0, deadc=0.95, \
         emin=0.1, emax=100.0, rad=5.0, \
-        irf="cta_dummy_irf", caldb="dummy"):
+        irf="South_50h", caldb="prod2"):
     """
     Returns an observation container filled with a list of CTA observations.
     The list is defined by the obsdeflist parameter which is a dictionnary
@@ -584,9 +584,9 @@ def spectrum(obs, source, ebounds):
         # Clone observations and reset energy thresholds
         select = ctools.ctselect(obs)
-        select["ra"].real(-1.0)
-        select["dec"].real(-1.0)
-        select["rad"].real(-1.0)
+        select["ra"].value("UNDEF")
+        select["dec"].value("UNDEF")
+        select["rad"].value("UNDEF")
         select["tmin"].real(0.0)
         select["tmax"].real(0.0)
         select["emin"].real(ebounds.emin(i).TeV())
@@ -105,15 +105,17 @@ class cslightcrv(ctools.cscript):
             pars.append(gammalib.GApplicationPar("expcube","f","a","NONE","","","Exposure cube file (only needed for stacked analysis)"))
             pars.append(gammalib.GApplicationPar("psfcube","f","a","NONE","","","PSF cube file (only needed for stacked analysis)"))
             pars.append(gammalib.GApplicationPar("bkgcube","s","a","NONE","","","Background cube file (only needed for stacked analysis)"))
-            pars.append(gammalib.GApplicationPar("caldb","s","a","dummy","","","Calibration database"))
-            pars.append(gammalib.GApplicationPar("irf","s","a","cta_dummy_irf","","","Instrument response function"))
-            pars.append(gammalib.GApplicationPar("tmin","r","a","51544.5","","days","Lightcurve Time start [MJD]"))
-            pars.append(gammalib.GApplicationPar("tmax","r","a","51544.6","","days","Lightcurve Time stop [MJD]"))
-            pars.append(gammalib.GApplicationPar("tbins","i","a","5","","","Number of time bins points"))
-            pars.append(gammalib.GApplicationPar("tbinalg","s", "a","GTI","FILE|LIN|GTI","", "Algorithm for defining energy bins"))
-            pars.append(gammalib.GApplicationPar("tbinfile","s", "a","tbins.fits","","", "File defining the time binning"))
-            pars.append(gammalib.GApplicationPar("binned","b","a","no","yes|no","","Use binned analysis in each energy bin"))
-            pars.append(gammalib.GApplicationPar("nebins","i","h","20","","","Number of energy bins per light curve point"))
+            pars.append(gammalib.GApplicationPar("caldb","s","a","prod2","","","Calibration database"))
+            pars.append(gammalib.GApplicationPar("irf","s","a","South_50h","","","Instrument response function"))
+            pars.append(gammalib.GApplicationPar("tbinalg","s","h","GTI","FILE|LIN|GTI","", "Algorithm for defining time bins"))
+            pars.append(gammalib.GApplicationPar("tmin","r","a","51544.5","","days","Lightcurve start time [MJD]"))
+            pars.append(gammalib.GApplicationPar("tmax","r","a","51544.6","","days","Lightcurve stop time [MJD]"))
+            pars.append(gammalib.GApplicationPar("tbins","i","a","5","","","Number of time bins"))
+            pars.append(gammalib.GApplicationPar("tbinfile","f","a","tbins.fits","","", "File defining the time binning"))
+            #pars.append(gammalib.GApplicationPar("binned","b","a","no","yes|no","","Use binned analysis in each time bin"))
+            pars.append(gammalib.GApplicationPar("enumbins","i","a","0","","","Number of energy bins per light curve bin (0=unbinned)"))
+            pars.append(gammalib.GApplicationPar("emin","r","a","0.1","","","Lower energy limit of events (TeV)"))
+            pars.append(gammalib.GApplicationPar("emax","r","a","100.0","","","Upper energy limit of events (TeV)"))
             pars.append(gammalib.GApplicationPar("coordsys","s","a","CEL","CEL|GAL","","Coordinate System"))
             pars.append(gammalib.GApplicationPar("proj","s","a","CAR","AIT|AZP|CAR|MER|MOL|STG|TAN","","Projection method"))
             pars.append(gammalib.GApplicationPar("xref","r","a","83.63","0","360","First coordinate of image center in degrees (RA or galactic l)"))
@@ -121,12 +123,10 @@ class cslightcrv(ctools.cscript):
             pars.append(gammalib.GApplicationPar("nxpix","i","a","200","","","Size of the X axis in pixels"))
             pars.append(gammalib.GApplicationPar("nypix","i","a","200","","","Size of the Y axis in pixels"))
             pars.append(gammalib.GApplicationPar("binsz","r","a","0.02","","","Pixel size (deg/pixel)"))
-            pars.append(gammalib.GApplicationPar("emin","r","a","0.1","","","Lower energy limit for binned analysis (TeV)"))
-            pars.append(gammalib.GApplicationPar("emax","r","a","100.0","","","Upper energy limit for binned analysis (TeV)"))
-            pars.append(gammalib.GApplicationPar("calc_ts","b","h","yes","yes|no","","Compute TS value in each bin"))
-            pars.append(gammalib.GApplicationPar("calc_ulim","b","h","yes","yes|no","","Compute upper limit in each bin"))
-            pars.append(gammalib.GApplicationPar("fix_srcs","b","h","yes","yes|no","","Fix other skymodel parameters"))
-            pars.append(gammalib.GApplicationPar("fix_bkg","b","h","no","yes|no","","Fix background parameters"))
+            pars.append(gammalib.GApplicationPar("calc_ts","b","h","yes","yes|no","","Compute TS value for each bin?"))
+            pars.append(gammalib.GApplicationPar("calc_ulim","b","h","yes","yes|no","","Compute upper limit for each bin?"))
+            pars.append(gammalib.GApplicationPar("fix_srcs","b","h","yes","yes|no","","Fix other sky model parameters?"))
+            pars.append(gammalib.GApplicationPar("fix_bkg","b","h","no","yes|no","","Fix background model parameters?"))
             pars.append_standard()
             pars.append(gammalib.GApplicationPar("logfile","f","h","cslightcrv.log","","","Log filename"))
             pars.save(parfile)
@@ -190,22 +190,33 @@ class cslightcrv(ctools.cscript):
         # Get source name   
         self.m_srcname = self["srcname"].string()
-        # Get ebounds             
+        # Get time boundaries             
         self.create_tbounds()
+        # Unbinned or binned analysis?
+        self.m_ebins    = self["enumbins"].integer()
+        if self.m_ebins == 0:
+            self.m_binned = False
+        else:
+            self.m_binned = True
+
+        # Get energy range
+        self.m_emin = self["emin"].real()
+        self.m_emax = self["emax"].real()
+
         # Get binning flag
-        self.m_binned = self["binned"].boolean()
+        #self.m_binned = self["binned"].boolean()
         if self.m_binned:
+            self.m_coordsys = self["coordsys"].string()
+            self.m_proj     = self["proj"].string()
             self.m_xref     = self["xref"].real()
             self.m_yref     = self["yref"].real()
             self.m_nxpix    = self["nxpix"].integer()
             self.m_nypix    = self["nypix"].integer()
             self.m_binsz    = self["binsz"].real()
-            self.m_coordsys = self["coordsys"].string()
-            self.m_proj     = self["proj"].string()
-            self.m_ebins    = self["nebins"].integer()
-            self.m_emin     = self["emin"].real()
-            self.m_emax     = self["emax"].real()
+            #self.m_ebins    = self["enumbins"].integer()
+            #self.m_emin     = self["emin"].real()
+            #self.m_emax     = self["emax"].real()
         # Read other parameters
         self.m_outfile = self["outfile"].filename()
@@ -367,8 +378,8 @@ class cslightcrv(ctools.cscript):
             # Select events
             select = ctools.ctselect(self.obs)
-            select["emin"].value("UNDEFINED")    
-            select["emax"].value("UNDEFINED") 
+            select["emin"].real(self.m_emin)    
+            select["emax"].real(self.m_emax) 
             select["tmin"].real(tmin.convert(select.time_reference()))
             select["tmax"].real(tmax.convert(select.time_reference()))
             select["rad"].value("UNDEFINED")
@@ -469,7 +480,7 @@ class cslightcrv(ctools.cscript):
                 obs = bin.obs()
                 # Set precomputed binned response
-                obs[0].response(expcube.expcube(), psfcube.psfcube(), bkgcube.cube())
+                obs[0].response(expcube.expcube(), psfcube.psfcube(), bkgcube.bkgcube())
                 # Get new models
                 models = bkgcube.models()
@@ -99,20 +99,20 @@ class cspull(ctools.cscript):
             pars.append(gammalib.GApplicationPar("expcube","s","a","NONE","","","Exposure cube file (only needed for stacked analysis)"))
             pars.append(gammalib.GApplicationPar("psfcube","s","a","NONE","","","PSF cube file (only needed for stacked analysis)"))
             pars.append(gammalib.GApplicationPar("bkgcube","s","a","NONE","","","Background cube file (only needed for stacked analysis)"))
-            pars.append(gammalib.GApplicationPar("caldb","s","a","dummy","","","Calibration database"))
-            pars.append(gammalib.GApplicationPar("irf","s","a","cta_dummy_irf","","","Instrument response function"))
+            pars.append(gammalib.GApplicationPar("caldb","s","a","prod2","","","Calibration database"))
+            pars.append(gammalib.GApplicationPar("irf","s","a","South_50h","","","Instrument response function"))
+            pars.append(gammalib.GApplicationPar("deadc","r","h","0.95","0","1","Deadtime correction factor"))
             pars.append(gammalib.GApplicationPar("edisp","b","h","no","","","Apply energy dispersion?"))
             pars.append(gammalib.GApplicationPar("ntrials","i","a","10","","","Number of trials"))
             pars.append(gammalib.GApplicationPar("ra","r","a","83.6331","0","360","RA of pointing (deg)"))
             pars.append(gammalib.GApplicationPar("dec","r","a","22.0145","-90","90","Dec of pointing (deg)"))
             pars.append(gammalib.GApplicationPar("emin","r","a","0.1","","","Lower energy limit (TeV)"))
             pars.append(gammalib.GApplicationPar("emax","r","a","100.0","","","Upper energy limit (TeV)"))
+            pars.append(gammalib.GApplicationPar("enumbins","i","a","0","","","Number of energy bins (0=unbinned)"))
             pars.append(gammalib.GApplicationPar("tmin","r","h","0.0","","","Start time (MET in s)"))
             pars.append(gammalib.GApplicationPar("tmax","r","a","1800.0","","","Duration (in s)"))
-            pars.append(gammalib.GApplicationPar("enumbins","i","a","0","","","Number of energy bins (0=unbinned)"))
             pars.append(gammalib.GApplicationPar("npix","i","a","200","","","Number of pixels for binned"))
             pars.append(gammalib.GApplicationPar("binsz","r","a","0.05","","","Pixel size for binned (deg/pixel)"))
-            pars.append(gammalib.GApplicationPar("deadc","r","h","0.95","0","1","Deadtime correction factor"))
             pars.append(gammalib.GApplicationPar("rad","r","h","5.0","0","180","Radius of ROI (deg)"))         
             pars.append(gammalib.GApplicationPar("pattern","s","h","single","","","Observation pattern (single/four)"))
             pars.append(gammalib.GApplicationPar("offset","r","h","1.5","","","Observation pattern offset (deg)"))
@@ -103,8 +103,8 @@ class csresmap(ctools.cscript):
             pars.append(gammalib.GApplicationPar("inobs","f","a","events.fits","","","Event list, counts cube, or observation definition file"))
             pars.append(gammalib.GApplicationPar("inmodel","f","a","$CTOOLS/share/models/crab.xml","","","Source model"))
             pars.append(gammalib.GApplicationPar("outmap","f","a","resmap.fits","","","Output residual map"))
-            pars.append(gammalib.GApplicationPar("caldb","s","a","dummy","","","Calibration database"))
-            pars.append(gammalib.GApplicationPar("irf","s","a","cta_dummy_irf","","","Instrument response function"))
+            pars.append(gammalib.GApplicationPar("caldb","s","a","prod2","","","Calibration database"))
+            pars.append(gammalib.GApplicationPar("irf","s","a","South_50h","","","Instrument response function"))
             pars.append(gammalib.GApplicationPar("ebinalg","s","h","LOG","LIN|LOG|FILE","","Binning algorithm"))
             pars.append(gammalib.GApplicationPar("emin","r","h","0.1","","","Lower energy limit (TeV)"))
             pars.append(gammalib.GApplicationPar("emax","r","h","100.0","","","Upper energy limit (TeV)"))
@@ -101,8 +101,8 @@ class cssens(ctools.cscript):
             pars.append(gammalib.GApplicationPar("inmodel","f","a","$CTOOLS/share/models/crab.xml","","","Source model"))
             pars.append(gammalib.GApplicationPar("srcname","s","a","Crab","","","Source name"))
             pars.append(gammalib.GApplicationPar("offset","r","h","0.0","","","Source offset angle (deg)"))
-            pars.append(gammalib.GApplicationPar("caldb","s","a","dummy","","","Calibration database"))
-            pars.append(gammalib.GApplicationPar("irf","s","a","cta_dummy_irf","","","Instrument response function"))
+            pars.append(gammalib.GApplicationPar("caldb","s","a","prod2","","","Calibration database"))
+            pars.append(gammalib.GApplicationPar("irf","s","a","South_50h","","","Instrument response function"))
             pars.append(gammalib.GApplicationPar("deadc","r","h","0.95","0","1","Deadtime correction factor"))
             pars.append(gammalib.GApplicationPar("outfile","f","h","sensitivity.dat","","","Output file name"))
             pars.append(gammalib.GApplicationPar("duration","r","a","180000.0","","","Effective exposure time (s)"))
@@ -106,8 +106,8 @@ class csspec(ctools.cscript):
             pars.append(gammalib.GApplicationPar("expcube","f","a","NONE","","","Exposure cube file (only needed for stacked analysis)"))
             pars.append(gammalib.GApplicationPar("psfcube","f","a","NONE","","","PSF cube file (only needed for stacked analysis)"))
             pars.append(gammalib.GApplicationPar("bkgcube","s","a","NONE","","","Background cube file (only needed for stacked analysis)"))
-            pars.append(gammalib.GApplicationPar("caldb","s","a","dummy","","","Calibration database"))
-            pars.append(gammalib.GApplicationPar("irf","s","a","cta_dummy_irf","","","Instrument response function"))
+            pars.append(gammalib.GApplicationPar("caldb","s","a","prod2","","","Calibration database"))
+            pars.append(gammalib.GApplicationPar("irf","s","a","South_50h","","","Instrument response function"))
             pars.append(gammalib.GApplicationPar("emin","r","h","0.1","","","Lower energy limit for spectral points(TeV)"))
             pars.append(gammalib.GApplicationPar("emax","r","h","100.0","","","Upper energy limit for spectral points(TeV)"))
             pars.append(gammalib.GApplicationPar("enumbins","i","a","20","","","Number of spectral points"))
@@ -424,7 +424,7 @@ class csspec(ctools.cscript):
                 obs = bin.obs()
                 # Set precomputed binned response
-                obs[0].response(expcube.expcube(), psfcube.psfcube(), bkgcube.cube())
+                obs[0].response(expcube.expcube(), psfcube.psfcube(), bkgcube.bkgcube())
                 # Get new models
                 models = bkgcube.models()
@@ -95,27 +95,29 @@ class cstsdist(ctools.cscript):
             # Create default parfile
             pars = gammalib.GApplicationPars()
             pars.append(gammalib.GApplicationPar("inobs","f","h","NONE","","","Event list, counts cube, or observation definition file"))
-            pars.append(gammalib.GApplicationPar("outfile","f","h","ts.dat","","","Output file name"))
+            pars.append(gammalib.GApplicationPar("inmodel","f","a","$CTOOLS/share/models/crab.xml","","","Source model"))
+            pars.append(gammalib.GApplicationPar("srcname","s","a","Crab","","","Source name"))
+            pars.append(gammalib.GApplicationPar("outfile","f","a","ts.dat","","","Output file name"))
             pars.append(gammalib.GApplicationPar("expcube","s","a","NONE","","","Exposure cube file (only needed for stacked analysis)"))
             pars.append(gammalib.GApplicationPar("psfcube","s","a","NONE","","","PSF cube file (only needed for stacked analysis)"))
-            pars.append(gammalib.GApplicationPar("caldb","s","a","dummy","","","Calibration database"))
-            pars.append(gammalib.GApplicationPar("irf","s","a","cta_dummy_irf","","","Instrument response function"))
-            pars.append(gammalib.GApplicationPar("ra","r","h","83.63","0","360","RA of pointing (deg)"))
-            pars.append(gammalib.GApplicationPar("dec","r","h","22.01","-90","90","Dec of pointing (deg)"))
+            pars.append(gammalib.GApplicationPar("bkgcube","s","a","NONE","","","Background cube file (only needed for stacked analysis)"))
+            pars.append(gammalib.GApplicationPar("caldb","s","a","prod2","","","Calibration database"))
+            pars.append(gammalib.GApplicationPar("irf","s","a","South_50h","","","Instrument response function"))
+            pars.append(gammalib.GApplicationPar("deadc","r","h","0.95","","","Deadtime correction factor"))
+            pars.append(gammalib.GApplicationPar("edisp","b","h","no","","","Apply energy dispersion?"))
+            pars.append(gammalib.GApplicationPar("ntrials","i","a","10","","","Number of trials"))
+            pars.append(gammalib.GApplicationPar("ra","r","a","83.63","0","360","RA of pointing (deg)"))
+            pars.append(gammalib.GApplicationPar("dec","r","a","22.01","-90","90","Dec of pointing (deg)"))
             pars.append(gammalib.GApplicationPar("emin","r","a","0.1","","","Lower energy limit (TeV)"))
             pars.append(gammalib.GApplicationPar("emax","r","a","100.0","","","Upper energy limit (TeV)"))
+            pars.append(gammalib.GApplicationPar("enumbins","i","a","0","","","Number of energy bins (0=unbinned)"))
             pars.append(gammalib.GApplicationPar("tmin","r","h","0.0","","","Start time (MET in s)"))
             pars.append(gammalib.GApplicationPar("tmax","r","a","1800.0","","","Observation duration (in s)"))
-            pars.append(gammalib.GApplicationPar("enumbins","i","a","0","","","Number of energy bins (0=unbinned)"))
-            pars.append(gammalib.GApplicationPar("npix","i","h","200","","","Number of pixels for binned"))
-            pars.append(gammalib.GApplicationPar("binsz","r","h","0.05","","","Pixel size for binned (deg/pixel)"))
-            pars.append(gammalib.GApplicationPar("deadc","r","h","0.95","","","Deadtime correction factor"))
+            pars.append(gammalib.GApplicationPar("npix","i","a","200","","","Number of pixels for binned"))
+            pars.append(gammalib.GApplicationPar("binsz","r","a","0.05","","","Pixel size for binned (deg/pixel)"))
             pars.append(gammalib.GApplicationPar("rad","r","h","5.0","","","Radius of ROI (deg)"))         
-            pars.append(gammalib.GApplicationPar("ntrials","i","a","10","","","Number of trials"))
-            pars.append(gammalib.GApplicationPar("type","s","a","point","","","Source model type (point/gauss/shell/disk)"))
-            pars.append(gammalib.GApplicationPar("index","r","h","-2.48","","","Spectral index"))       
-            pars.append(gammalib.GApplicationPar("offset","r","a","0.0","","","Source offset angle (deg)"))
-            pars.append(gammalib.GApplicationPar("bkg","s","a","$CTOOLS/share/models/bkg_dummy.txt","","","Background model file function"))
+            pars.append(gammalib.GApplicationPar("pattern","s","h","single","","","Observation pattern (single/four)"))
+            pars.append(gammalib.GApplicationPar("offset","r","h","1.5","","","Observation pattern offset (deg)"))
             pars.append_standard()
             pars.append(gammalib.GApplicationPar("logfile","f","h","cstsdist.log","","","Log filename"))
             pars.save(parfile)
@@ -131,6 +133,15 @@ class cstsdist(ctools.cscript):
         # Set observation if not done before
         if self.obs == None or self.obs.size() == 0:
             self.obs = self.get_observations()
+
+            # Check for requested pattern and use above
+            # observation parameters to set wobble pattern
+            self.m_pattern = self["pattern"].string()
+            if self.m_pattern == "four":
+                self.obs = self.set_obs()
+
+        # Get source name
+        self.m_srcname = self["srcname"].string()
         # Get number of energy bins
         self.m_enumbins = self["enumbins"].integer()
@@ -143,36 +154,21 @@ class cstsdist(ctools.cscript):
             # Set dummy values (required by obsutils)
             self.m_npix = 0
             self.m_binsz = 0.0
+
+        # Set models if we have none
+        if self.obs.models().size() == 0:
+            self.obs.models(self["inmodel"].filename())
         # Get other parameters
-        self.m_outfile  = self["outfile"].filename()
-        self.m_ntrials  = self["ntrials"].integer()
-        self.m_type     = self["type"].string()
-        self.m_index    = self["index"].real()
-        self.m_offset   = self["offset"].real()
-        self.m_bkg      = self["bkg"].string()
+        self.m_outfile = self["outfile"].filename()
+        self.m_ntrials = self["ntrials"].integer()
+        self.m_edisp   = self["edisp"].boolean()
+        self.m_offset  = self["offset"].real()
         # Set some fixed parameters
-        self.m_log   = False # Logging in client tools
-        self.m_debug = False # Debugging in client tools
-          
-        pnt = gammalib.GSkyDir()
-        pnt.radec_deg(self["ra"].real(),self["dec"].real())  
-          
-        # Initialise models. Note that we centre the point source at the
-        # center of our observation, so we're onaxis.
-        self.bkg_model  = gammalib.GModels()
-        self.full_model = gammalib.GModels()
-        self.bkg_model.append(self.set_bkg_model())
-        self.full_model.append(self.set_bkg_model())
-        self.full_model.append(self.set_src_model(pnt.l_deg(), pnt.b_deg()+self.m_offset, \
-                                                  flux=0.010, \
-                                                  type=self.m_type, \
-                                                  index=self.m_index))
-
-        # Attach background model to observation container
-        self.obs.models(self.bkg_model)
-        
+        self.m_log   = False                    # Logging in client tools
+        self.m_debug = self["debug"].boolean()  # Debugging in client tools
+
         # Return
         return
@@ -211,6 +207,20 @@ class cstsdist(ctools.cscript):
         if self.logTerse():
             self.log_parameters()
             self.log("\n")
+
+        # Initialise models
+        full_model, bkg_model = self.set_models()
+
+        # Write models into logger
+        if self.logTerse():
+            self.log("\n")
+            self.log.header1("Models")
+            self.log.header2("Background model")
+            self.log(str(bkg_model))
+            self.log("\n\n")
+            self.log.header2("Full model")
+            self.log(str(full_model))
+            self.log("\n")
         # Write observation into logger
         if self.logTerse():
@@ -235,7 +245,7 @@ class cstsdist(ctools.cscript):
         for seed in range(self.m_ntrials):
             # Make a trial
-            result = self.trial(seed)
+            result = self.trial(seed, full_model, bkg_model)
             # Write out result immediately
             if seed == 0:
@@ -251,97 +261,89 @@ class cstsdist(ctools.cscript):
         # Return
         return
-    def set_bkg_model(self, fitsigma=False):
+    def set_models(self, fitpos=False, fitspec=False):
         """
-        Setup CTA background model.
+        Set full and background model.
         """
-        # Define radial component
-        radial = gammalib.GCTAModelRadialGauss(3.0)
-        if fitsigma:
-            radial["Sigma"].free()
-        else:
-            radial["Sigma"].fix()
-        
-        # Define spectral component
-        spectrum = gammalib.GModelSpectralFunc(self.m_bkg, 1.0)
+        # Retrieve full model from observation container
+        full_model = self.obs.models().copy()
+
+        # Get source model
+        model = full_model[self.m_srcname]
-        # Create background model
-        model = gammalib.GCTAModelRadialAcceptance(radial, spectrum)
-        model.name("Background")
-        model.instruments("CTA")
-    
-        # Return background model
-        return model
-    
-    def set_src_model(self, l, b, flux=1.0, index=-2.48, \
-                      type="point", sigma=1.0, radius=1.0, width=0.1, \
-                      fitpos=False, fitidx=False):
-        """
-        Returns a single source with Crab-like spectrum. The source flux
-        can be scaled in Crab units. The Crab spectrum is based on MAGIC
-        observations (Albert et al. 2008, ApJ, 674, 1037).
+        # Check that model has a Prefactor
+        if not model.has_par("Prefactor"):
+            msg = "Model \""+self.m_srcname+"\" has no parameter \"Prefactor\"."+ \
+                  " Only spectral models with a \"Prefactor\" parameter are supported."
+            raise gammalib.GException.invalid_value("cssens", msg)
-        Parameters:
-         l      - Galactic longitude of source location [deg]
-         b      - Galactic latitude of source location [deg]
-        Keywords:
-         flux   - Source flux [Crabs]
-         index  - Spectral index
-         type   - Source type ("point", "gauss", "disk", "shell")
-         sigma  - Gaussian sigma (for type="gauss")
-         radius - Disk or shell inner radius [deg] (for type="disk" and type="shell")
-         width  - Shell width [deg] (for type="shell")
-         fitpos - Fit position? (default: True)
-         fitidx - Fit index? (default: True)
-        """
-        # Set source location
-        location = gammalib.GSkyDir()
-        location.lb_deg(l, b)
-    
-        # Set source spectrum
-        spectrum = gammalib.GModelSpectralPlaw(flux*5.7e-16, index, gammalib.GEnergy(0.3, "TeV"))
-        if fitidx:
-            spectrum["Index"].free()
+        # Fit or fix spatial parameters
+        if fitpos:
+            if model.has_par("RA"):
+                model["RA"].free()
+            if model.has_par("DEC"):
+                model["DEC"].free()
+            if model.has_par("Sigma"):
+                model["Sigma"].free()
+            if model.has_par("Radius"):
+                model["Radius"].free()
+            if model.has_par("Width"):
+                model["Width"].free()
         else:
-            spectrum["Index"].fix() 
-
-        # Set source
-        if type == "point":
-            spatial = gammalib.GModelSpatialPointSource(location)
-            if fitpos:
-                spatial[0].free()
-                spatial[1].free()
-        elif type == "gauss":
-            spatial = gammalib.GModelSpatialRadialGauss(location, sigma)
-            if fitpos:
-                spatial[0].free()
-                spatial[1].free()
-        elif type == "disk":
-            spatial = gammalib.GModelSpatialRadialDisk(location, radius)
-            if fitpos:
-                spatial[0].free()
-                spatial[1].free()
-        elif type == "shell":
-            spatial = gammalib.GModelSpatialRadialShell(location, radius, width)
-            if fitpos:
-                spatial[0].free()
-                spatial[1].free()
+            if model.has_par("RA"):
+                model["RA"].fix()
+            if model.has_par("DEC"):
+                model["DEC"].fix()
+            if model.has_par("Sigma"):
+                model["Sigma"].fix()
+            if model.has_par("Radius"):
+                model["Radius"].fix()
+            if model.has_par("Width"):
+                model["Width"].fix()
+
+        # Fit or fix spectral parameters
+        if fitspec:
+            if model.has_par("Index"):
+                model["Index"].free()
+            if model.has_par("Cutoff"):
+                model["Cutoff"].free()
         else:
-            self.log("ERROR: Unknown source type '"+type+"'.\n")
-            return None
-        source = gammalib.GModelSky(spatial, spectrum)
+            if model.has_par("Index"):
+                model["Index"].fix()
+            if model.has_par("Cutoff"):
+                model["Cutoff"].fix()
-        # Set source name
-        source.name("Test")
+        # Create background model
+        bkg_model = full_model.copy()
+        bkg_model.remove(self.m_srcname)
+
+        # Return models
+        return full_model, bkg_model
+
+    def set_obs(self):
+        """
+        Returns an observation container with a set of CTA observations.
-        # Turn TS value computation off
-        # since we compute it by hand later on
-        source.tscalc(False)
-    
-        # Return source
-        return source
+        Keywords:
+        """
+        
+        # Setup observation definition list
+        obsdeflist = obsutils.set_obs_patterns(self.m_pattern, \
+                                               ra=self["ra"].real(), dec=self["dec"].real(), \
+                                               offset=self["offset"].real())
+        
+        # Create list of observations
+        obs = obsutils.set_obs_list(obsdeflist, \
+                                    tstart=self["tmin"].real(), duration=self["tmax"].real()-self["tmin"].real(), \
+                                    deadc=self["deadc"].real(), \
+                                    emin=self["emin"].real(), emax=self["emax"].real(), \
+                                    rad=self["rad"].real(), \
+                                    irf=self["irf"].string(), caldb=self["caldb"].string())
-    def trial(self, seed):
+        # Return observation container
+        return obs
+
+    def trial(self, seed, full_model, bkg_model):
         """
         Create the TS for a single trial.
@@ -373,9 +375,9 @@ class cstsdist(ctools.cscript):
             self.log("\n")
         # Fit background only
-        sim.models(self.bkg_model)
+        sim.models(bkg_model)
         like_bgm   = obsutils.fit(sim, log=self.m_log, debug=self.m_debug)
-        result_bgm = like_bgm.obs().models()
+        result_bgm = like_bgm.obs().models().copy()
         LogL_bgm   = like_bgm.opt().value()
         npred_bgm  = like_bgm.obs().npred()
@@ -396,9 +398,9 @@ class cstsdist(ctools.cscript):
                     self.log(str(par)+"\n")
         # Fit background and test source
-        sim.models(self.full_model)
+        sim.models(full_model)
         like_all   = obsutils.fit(sim, log=self.m_log, debug=self.m_debug)
-        result_all = like_all.obs().models()
+        result_all = like_all.obs().models().copy()
         LogL_all   = like_all.opt().value()
         npred_all  = like_all.obs().npred()
         ts         = 2.0*(LogL_bgm-LogL_all)
@@ -428,9 +430,9 @@ class cstsdist(ctools.cscript):
             self.log("TS=")
             self.log(ts)
             self.log("  Prefactor=")
-            self.log(result_all["Test"]["Prefactor"].value())
+            self.log(result_all[self.m_srcname]["Prefactor"].value())
             self.log("+/-")
-            self.log(result_all["Test"]["Prefactor"].error())
+            self.log(result_all[self.m_srcname]["Prefactor"].error())
             self.log("\n")
         # Initialise results
@@ -59,7 +59,7 @@ public:
     void                      clear(void);
     void                      run(void);
     void                      save(void);
-    const GCTACubeBackground& cube(void) const;
+    const GCTACubeBackground& bkgcube(void) const;
     const GModels&            models(void) const;
 protected:
@@ -86,7 +86,7 @@ protected:
  * @return Background response cube.
  ***************************************************************************/
 inline
-const GCTACubeBackground& ctbkgcube::cube(void) const
+const GCTACubeBackground& ctbkgcube::bkgcube(void) const
 {
     // Return background cube
     return (m_background);
@@ -22,8 +22,8 @@ inmodel,  f, a, &quot;NONE&quot;,,, &quot;Input model XML file&quot;
 incube,   f, a, "NONE",,, "Counts cube for background cube definition"
 outcube,  f, a, "bkgcube.fits",,, "Output background cube file"
 outmodel, f, a, "NONE",,, "Output model XML file"
-caldb,    s, a, "dummy",,, "Calibration database"
-irf,      s, a, "cta_dummy_irf",,, "Instrument response function"
+caldb,    s, a, "prod2",,, "Calibration database"
+irf,      s, a, "South_50h",,, "Instrument response function"
 #
 # Energy binning parameters
 #==========================
@@ -23,8 +23,8 @@ srcname,  s, a, &quot;Crab&quot;,,,&quot;Source of interest&quot;
 expcube,  f, a, "NONE",,, "Exposure cube file (only needed for stacked analysis)"
 psfcube,  f, a, "NONE",,, "PSF cube file (only needed for stacked analysis)"
 bkgcube,  f, a, "NONE",,, "Background cube file (only needed for stacked analysis)"
-caldb,    s, a, "dummy",,, "Calibration database"
-irf,      s, a, "cta_dummy_irf",,, "Instrument response function"
+caldb,    s, a, "prod2",,, "Calibration database"
+irf,      s, a, "South_50h",,, "Instrument response function"
 outfile,  f, a, "butterfly.txt",,, "Output ascii file"
 matrix,   f, h, "NONE",,, "Input covariance Matrix FITS file"
@@ -20,8 +20,8 @@
 inobs,    f, a, "NONE",,, "Event list or observation definition file"
 incube,   f, a, "NONE",,, "Counts cube for exposure cube definition"
 outcube,  f, a, "expcube.fits",,, "Output exposure cube file"
-caldb,    s, a, "dummy",,, "Calibration database"
-irf,      s, a, "cta_dummy_irf",,, "Instrument response function"
+caldb,    s, a, "prod2",,, "Calibration database"
+irf,      s, a, "South_50h",,, "Instrument response function"
 edisp,    b, h, no,,, "Apply energy dispersion?"
 #
 # Energy binning parameters
@@ -22,8 +22,8 @@ inmodel,  f, a, &quot;$CTOOLS/share/models/crab.xml&quot;,,, &quot;Source model&quot;
 expcube,  f, a, "NONE",,, "Exposure cube file (only needed for stacked analysis)"
 psfcube,  f, a, "NONE",,, "PSF cube file (only needed for stacked analysis)"
 bkgcube,  f, a, "NONE",,, "Background cube file (only needed for stacked analysis)"
-caldb,    s, a, "dummy",,, "Calibration database"
-irf,      s, a, "cta_dummy_irf",,, "Instrument response function"
+caldb,    s, a, "prod2",,, "Calibration database"
+irf,      s, a, "South_50h",,, "Instrument response function"
 outmodel, f, a, "crab_results.xml",,, "Source model output file"
 stat,     s, h, "POISSON",,, "Optimization statistics"
 edisp,    b, h, no,,, "Apply energy dispersion?"
@@ -28,8 +28,8 @@ inobs,    f, a, &quot;NONE&quot;,,, &quot;Event list, counts cube or observation definition fil
 expcube,  f, a, "NONE",,, "Exposure cube file (only needed for stacked analysis)"
 psfcube,  f, a, "NONE",,, "PSF cube file (only needed for stacked analysis)"
 bkgcube,  f, a, "NONE",,, "Background cube file (only needed for stacked analysis)"
-caldb,    s, a, "dummy",,, "Calibration database"
-irf,      s, a, "cta_dummy_irf",,, "Instrument response function"
+caldb,    s, a, "prod2",,, "Calibration database"
+irf,      s, a, "South_50h",,, "Instrument response function"
 ra,       r, a, 83.63,0,360, "RA of pointing (degrees)"
 dec,      r, a, 22.01,-90,90, "Dec of pointing (degrees)"
 rad,      r, a, 5.0,0,180, "Radius of FOV (degrees)"
@@ -965,7 +965,7 @@ void ctobssim::simulate_source(GCTAObservation* obs, const GModels&amp; models,
  *
  * @param[in] model Sky model.
  * @param[in] emin Minimum energy.
- * @param[in] emin Maximum energy.
+ * @param[in] emax Maximum energy.
  * @param[in] centre Centre of region for photon rate determination.
  * @param[in] radius Radius of region for photon rate determination.
  * @return Model flux (photons/cm2/sec).
@@ -24,8 +24,8 @@ prefix,    s, h, &quot;sim_events_&quot;,,, &quot;Prefix for event lists in observation definit
 #
 # Observation parameters
 #=======================
-caldb, s, a, "dummy",,, "Calibration database"
-irf,   s, a, "cta_dummy_irf",,, "Instrument response function"
+caldb, s, a, "prod2",,, "Calibration database"
+irf,   s, a, "South_50h",,, "Instrument response function"
 #
 # Simulation parameters
 #======================
@@ -20,8 +20,8 @@
 inobs,    f, a, "NONE",,, "Event list or observation definition file"
 incube,   f, a, "NONE",,, "Counts cube for psf cube definition"
 outcube,  f, a, "psfcube.fits",,, "Output psf cube file"
-caldb,    s, a, "dummy",,, "Calibration database"
-irf,      s, a, "cta_dummy_irf",,, "Instrument response function"
+caldb,    s, a, "prod2",,, "Calibration database"
+irf,      s, a, "South_50h",,, "Instrument response function"
 edisp,    b, h, no,,, "Apply energy dispersion?"
 #
 # Energy binning parameters
@@ -23,8 +23,8 @@ srcname, s, a, &quot;Crab&quot;,,, &quot;Test source&quot;
 expcube, f, a, "NONE",,, "Exposure cube file (only needed for stacked analysis)"
 psfcube, f, a, "NONE",,, "PSF cube file (only needed for stacked analysis)"
 bkgcube, f, a, "NONE",,, "Background cube file (only needed for stacked analysis)"
-caldb,   s, a, "dummy",,, "Calibration database"
-irf,     s, a, "cta_dummy_irf",,, "Instrument response function"
+caldb,   s, a, "prod2",,, "Calibration database"
+irf,     s, a, "South_50h",,, "Instrument response function"
 outmap,  f, a, "tsmap.fits",,, "Output Test Statistic map"
 #
@@ -23,8 +23,8 @@ srcname,  s, a, &quot;Crab&quot;,,, &quot;Source of interest&quot;
 expcube,  f, a, "NONE",,, "Exposure cube file (only needed for stacked analysis)"
 psfcube,  f, a, "NONE",,, "PSF cube file (only needed for stacked analysis)"
 bkgcube,  f, a, "NONE",,, "Background cube file (only needed for stacked analysis)"
-caldb,    s, a, "dummy",,, "Calibration database"
-irf,      s, a, "cta_dummy_irf",,, "Instrument response function"
+caldb,    s, a, "prod2",,, "Calibration database"
+irf,      s, a, "South_50h",,, "Instrument response function"
 #
 # Upper limit calculation parameters
@@ -136,13 +136,14 @@ $ECHO &quot; ok&quot;
 $ECHO -n "Test cstsdist: "
 chmod +x @abs_top_srcdir@/scripts/cstsdist.py
 @abs_top_srcdir@/scripts/cstsdist.py \
+  inmodel="data/crab.xml" \
+  srcname="Crab" \
   outfile="ts.dat" \
   ntrials="1" \
   caldb="irf" \
   irf="cta_dummy_irf" \
-  type="point" \
-  offset="0.0" \
-  bkg="data/bkg_dummy.txt" \
+  ra=83.63 \
+  dec=22.01 \
   emin=0.1 \
   emax=100.0 \
   enumbins=0 \
@@ -160,6 +161,7 @@ else
 fi
 $ECHO " ok"
+
 #
 # Test csresmap
 # =============
@@ -189,6 +191,7 @@ else
 fi
 $ECHO " ok"
+
 #
 # Test csspec (unbinned)
 # =====================
@@ -215,6 +218,7 @@ else
 fi
 $ECHO " ok"
+
 #
 # Test csspec (binned)
 # ====================
@@ -249,6 +253,7 @@ else
 fi
 $ECHO " ok"
+
 #
 # Test cslightcrv (unbinned)
 # ==========================
@@ -257,14 +262,16 @@ chmod +x @abs_top_srcdir@/scripts/cslightcrv.py
 @abs_top_srcdir@/scripts/cslightcrv.py \
   inobs="data/crab_events.fits" \
   inmodel="data/crab.xml" \
+  srcname="Crab" \
   caldb="irf" \
   irf="cta_dummy_irf" \
-  tbins=3 \
+  tbinalg="LIN" \
   tmin=51544.5 \
   tmax=51544.53 \
-  srcname="Crab" \
-  tbinalg="LIN" \
-  binned=no \
+  tbins=3 \
+  enumbins=0 \
+  emin=0.1 \
+  emax=100.0 \
   outfile="lightcurve.fits" > cslightcrv.log
 $ECHO -n "."
 if [ -s "lightcurve.fits" ]
@@ -276,6 +283,7 @@ else
 fi
 $ECHO " ok"
+
 #
 # Test cslightcrv (binned)
 # ========================
@@ -284,21 +292,20 @@ chmod +x @abs_top_srcdir@/scripts/cslightcrv.py
 @abs_top_srcdir@/scripts/cslightcrv.py \
   inobs="data/crab_events.fits" \
   inmodel="data/crab.xml" \
+  srcname="Crab" \
   caldb="irf" \
   irf="cta_dummy_irf" \
   tbinalg="GTI" \
-  srcname="Crab" \
+  enumbins=2 \
   emin=0.1 \
   emax=100.0 \
-  binned=yes \
-  nxpix=10 \
-  nypix=10 \
-  binsz=0.1 \
-  nebins=2 \
   coordsys="CEL" \
+  proj="CAR" \
   xref=83.63 \
   yref=22.01 \
-  proj="CAR" \
+  nxpix=10 \
+  nypix=10 \
+  binsz=0.1 \
   outfile="lightcurve.fits" > cslightcrv.log
 $ECHO -n "."
 if [ -s "lightcurve.fits" ]
@@ -33,7 +33,6 @@ cd examples
 # Create local environment
 rm -rf ctools
-#mkdir -p ctools/share
 mkdir -p ctools/share/caldb
 ln -s $base/caldb ctools/share/caldb/cta
 ln -s $base/models ctools/share/models
@@ -52,11 +51,11 @@ export CTOOLS=$PWD/ctools
 # Run makers
 ./make_survey.py
-./make_ts_distributions.py 2 2
+./make_ts_distributions.py -n 2 -e 0 -d 1800
 #./make_pull_at_sensitivity_limit.py
 # Remove local environment
-rm -rf gammalib
+rm -rf ctools
 # Step back to base directory
 cd $base