Update IACT How To's

Jürgen Knödlseder
1 parent 9159e38d
Showing 11 changed files with 537 additions and 1740 deletions
cscripts/csiactobs.par
doc/source/users/tutorials/hess_dr1/index.rst
doc/source/users/tutorials/hess_dr1/jointfermi.rst
doc/source/users/tutorials/howto/iact/analysis.rst
doc/source/users/tutorials/howto/iact/copy.rst
doc/source/users/tutorials/howto/iact/high_level.rst
doc/source/users/tutorials/howto/iact/index.rst
doc/source/users/tutorials/howto/iact/select.rst
doc/source/users/tutorials/howto/iact/simulation.rst
doc/source/users/tutorials/howto/iact/tips.rst
doc/source/users/tutorials/howto/index.rst
@@ -22,7 +22,7 @@ infile,      f, a, runlist.lis,,, &quot;Input runlist file&quot;
 inmodel,     f, h, NONE,,, "Input model definition XML file (optional)"
 master_indx, s, h, master.json,,, "Input master index file"
 outobs,      f, a, obs.xml,,, "Output observation definition XML file"
-outmodel,    f, a, bgmodels.xml,,, "Output model definition XML file"
+outmodel,    f, a, bkgmodels.xml,,, "Output model definition XML file"
 #
 # Script parameters
@@ -31,7 +31,7 @@ bkgpars,          i, a, 1,,, &quot;Number of free parameters per background model&quot;
 bkg_scale,        b, h, yes,,, "Apply background scaling factor from observation index file if available?"
 ev_hiera,         s, h, events,,, "Hierarchy of event formats"
 aeff_hiera,       s, h, aeff_2d,,, "Hierarchy of effective area formats"
-psf_hiera,        s, h, psf_king|psf_3gauss,,, "Hierarchy of psf formats"
+psf_hiera,        s, h, psf_king|psf_3gauss|psf_table,,, "Hierarchy of psf formats"
 edisp_hiera,      s, h, edisp_2d,,, "Hierarchy of energy dispersion formats"
 bkg_hiera,        s, h, bkg_3d,,, "Hierarchy of background formats"
 bkg_mod_hiera,    s, h, irf|aeff|gauss,,, "Hierarchy of background models"
-.. _1dc:
+.. _hess_dr1:
 Analysing public H.E.S.S. data
 ==============================
-.. _hess_dr1_data:
+.. _hess_dr1_fermi:
 Joint analysis of H.E.S.S. and Fermi data
 =========================================
 .. _sec_iact_analysis:
-Analysing IACT data
-===================
-
-In order to analyse data from current IACTs, you first have to retrieve the
-data. Information on how to do that can be found in the section :ref:`sec_iact_copy`.
-
-Check available FITS production
--------------------------------
-Before you start an analysis, it is important that you know what kind of data
-need to be analysed. In the current scheme of
-`storing IACT data <http://gamma-astro-data-formats.readthedocs.org/en/latest/data_storage/index.html>`_,
-every FITS data set is represented by a unique name, e.g.
-``"hess-fits-pa-release-1.0-Prod26-MppStd"``.
-You have to pass this string to some set of scripts to find and use the data
-for analysis.
-
-To check the available names of FITS data sets, you should run the
-:ref:`csiactdata` script:
-
-.. code-block:: bash
-
-  $ csiactdata
-  Path were data is located [] /path/to/fits/data/
-  
-This script will list the names of available FITS data sets into the console.
-In the following it is assumed that you have selected a data set named
-``fits-data-name``.
-You can omit the query for the path where the data is located by setting the
-environment variable ``VHEFITS`` to this location. Also other scripts will
-benefit from having this environment variable set, and this avoids that you
-have to type in the path name every time you use one of these scripts.
-
-
-Find observations
------------------
-
-Now that you know which data are available you have to assemble a list of
-observations (a.k.a. runs) that you want to use in an analysis. You will
-use the :ref:`csfindobs` for that task, which allows you to select the
-observations according to a number of user criteria.
-
-
-Select observations according to pointing direction
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-The most common search for observations is by pointing direction on the sky,
-and you do this by simply invoking the :ref:`csfindobs` script:
-
-.. code-block:: bash
-
-  $ csfindobs
-  Name of FITS production (Run csiactdata to view your options) [fits-data-name] 
-  Right ascension of selection region centre (deg) [83.63]
-  Declination of selection region centre (deg) [22.01]
-  Search radius of selection region (deg) [2.5]
-  Output runlist file [runlist.lis]
-
-The script will create an ASCII file named ``runlist.lis`` that contains the
-identifiers for all observations with pointing directions located within 2.5°
-around the position of the Crab nebula.
-
-
-Select observations using additional criteria
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-You may add further selection criteria using the hidden ``expression`` parameter.
-You may use any FITS selection expression that is
-`supported by the cfitsio row filtering specification <https://heasarc.gsfc.nasa.gov/docs/software/fitsio/c/c_user/node97.html>`_.
-Column names that can be used in the expression can be found
-`here <http://gamma-astro-data-formats.readthedocs.org/en/latest/data_storage/obs_index/index.html>`_
-or by simply browsing the observation index file with a FITS viewer. For example
-
-.. code-block:: bash
-
-  $ csfindobs expression="ZEN_PNT<30&&LIVETIME>1800"
-
-will select all observations that have a zenith angle less than 30° and a
-livetime larger than 1800 seconds.
-
-
-Select observations independent of pointing direction
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-You can also omit the pointing selection, for example in case that you do not
-know the sky coordinates of your object of interest. In that case, you may
-directly use the object name in the expression and pass ``NONE`` as coordinate
-to :ref:`csfindobs`:
-
-.. code-block:: bash
-
-  $ csfindobs expression="OBJECT=='Crab Nebula'"
-  Name of FITS production (Run csiactdata to view your options) [fits-data-name] 
-  Right ascension of selection region centre [NONE]
-  Output runlist file [runlist.lis]
+Performing a maximum likelihood fit
+===================================
+
+  .. admonition:: What you will learn
+
+     You will learn how to perform a maximum likelihood fit of IACT data.
+
+Following a :ref:`csiactobs` run with the ``inmodel=$CTOOLS/share/models/crab.xml``
+you have a model that is ready for maximum likelihood fitting.
+You do such a fit using
+
+.. code-block:: bash
+
+   $ ctlike
+   Input event list, counts cube or observation definition XML file [obs_selected.xml]
+   Input model definition XML file [models.xml]
+   Output model definition XML file [crab_results.xml]
+
+which will produce the following output in the log file ``ctlike.log``:
+
+.. code-block:: none
+
+   2019-04-05T14:18:17: +=================================+
+   2019-04-05T14:18:17: | Maximum likelihood optimisation |
+   2019-04-05T14:18:17: +=================================+
+   2019-04-05T14:18:17:    Parameter "Prefactor" has zero curvature. Fix parameter.
+   2019-04-05T14:18:17:    Parameter "Index" has zero curvature. Fix parameter.
+   2019-04-05T14:18:17:  >Iteration   0: -logL=203790.772, Lambda=1.0e-03
+   2019-04-05T14:18:17:  >Iteration   1: -logL=192855.522, Lambda=1.0e-03, delta=10935.250, step=1.0e+00, max(|grad|)=-2040.622705 [Prefactor:0]
+   2019-04-05T14:18:17:  >Iteration   2: -logL=182999.110, Lambda=1.0e-04, delta=9856.412, step=1.0e+00, max(|grad|)=-981.683624 [Prefactor:0]
+   2019-04-05T14:18:17:  >Iteration   3: -logL=174691.964, Lambda=1.0e-05, delta=8307.146, step=1.0e+00, max(|grad|)=620.662474 [Index:9]
+   2019-04-05T14:18:17:  >Iteration   4: -logL=168496.493, Lambda=1.0e-06, delta=6195.471, step=1.0e+00, max(|grad|)=587.098865 [Index:9]
+   2019-04-05T14:18:17:  >Iteration   5: -logL=164898.487, Lambda=1.0e-07, delta=3598.006, step=1.0e+00, max(|grad|)=423.800306 [Index:5]
+   2019-04-05T14:18:18:  >Iteration   6: -logL=163621.925, Lambda=1.0e-08, delta=1276.562, step=1.0e+00, max(|grad|)=237.546032 [Index:1]
+   2019-04-05T14:18:18:  >Iteration   7: -logL=163423.407, Lambda=1.0e-09, delta=198.518, step=1.0e+00, max(|grad|)=102.646722 [Index:1]
+   2019-04-05T14:18:18:  >Iteration   8: -logL=163409.717, Lambda=1.0e-10, delta=13.690, step=1.0e+00, max(|grad|)=35.371246 [Index:1]
+   2019-04-05T14:18:18:  >Iteration   9: -logL=163408.542, Lambda=1.0e-11, delta=1.175, step=1.0e+00, max(|grad|)=11.708402 [Index:1]
+   2019-04-05T14:18:18:  >Iteration  10: -logL=163408.416, Lambda=1.0e-12, delta=0.126, step=1.0e+00, max(|grad|)=3.907383 [Index:1]
+   2019-04-05T14:18:18:  >Iteration  11: -logL=163408.402, Lambda=1.0e-13, delta=0.014, step=1.0e+00, max(|grad|)=1.309476 [Index:1]
+   2019-04-05T14:18:18:  >Iteration  12: -logL=163408.400, Lambda=1.0e-14, delta=0.002, step=1.0e+00, max(|grad|)=0.439462 [Index:1]
+   2019-04-05T14:18:18:    Free parameter "Prefactor" after convergence was reached with frozen parameter.
+   2019-04-05T14:18:18:    Free parameter "Index" after convergence was reached with frozen parameter.
+   2019-04-05T14:18:18:
+   2019-04-05T14:18:18: +=========================================+
+   2019-04-05T14:18:18: | Maximum likelihood optimisation results |
+   2019-04-05T14:18:18: +=========================================+
+   2019-04-05T14:18:18: === GOptimizerLM ===
+   2019-04-05T14:18:18:  Optimized function value ..: 163408.400
+   2019-04-05T14:18:18:  Absolute precision ........: 0.005
+   2019-04-05T14:18:18:  Acceptable value decrease .: 2
+   2019-04-05T14:18:18:  Optimization status .......: converged
+   2019-04-05T14:18:18:  Number of parameters ......: 26
+   2019-04-05T14:18:18:  Number of free parameters .: 12
+   2019-04-05T14:18:18:  Number of iterations ......: 12
+   2019-04-05T14:18:18:  Lambda ....................: 1e-15
+   2019-04-05T14:18:18:  Maximum log likelihood ....: -163408.400
+   2019-04-05T14:18:18:  Observed events  (Nobs) ...: 18641.000
+   2019-04-05T14:18:18:  Predicted events (Npred) ..: 18640.998 (Nobs - Npred = 0.00164345248049358)
+   2019-04-05T14:18:18: === GModels ===
+   2019-04-05T14:18:18:  Number of models ..........: 6
+   2019-04-05T14:18:18:  Number of parameters ......: 26
+   2019-04-05T14:18:18: === GCTAModelAeffBackground ===
+   2019-04-05T14:18:18:  Name ......................: bkg_23523
+   2019-04-05T14:18:18:  Instruments ...............: HESS
+   2019-04-05T14:18:18:  Instrument scale factors ..: unity
+   2019-04-05T14:18:18:  Observation identifiers ...: 23523
+   2019-04-05T14:18:18:  Model type ................: "PowerLaw" * "Constant"
+   2019-04-05T14:18:18:  Number of parameters ......: 4
+   2019-04-05T14:18:18:  Number of spectral par's ..: 3
+   2019-04-05T14:18:18:   Prefactor ................: 2.98790819794319e-13 +/- 4.43400659819813e-15 [1e-16,1e-12] ph/cm2/s/MeV (free,scale=1e-14,gradient)
+   2019-04-05T14:18:18:   Index ....................: -3.22702359177827 +/- 0.0221521174108544 [-5,5]  (free,scale=1,gradient)
+   2019-04-05T14:18:18:   PivotEnergy ..............: 1000000 MeV (fixed,scale=1000000,gradient)
+   2019-04-05T14:18:18:  Number of temporal par's ..: 1
+   2019-04-05T14:18:18:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+   2019-04-05T14:18:18: === GCTAModelAeffBackground ===
+   2019-04-05T14:18:18:  Name ......................: bkg_23526
+   2019-04-05T14:18:18:  Instruments ...............: HESS
+   2019-04-05T14:18:18:  Instrument scale factors ..: unity
+   2019-04-05T14:18:18:  Observation identifiers ...: 23526
+   2019-04-05T14:18:18:  Model type ................: "PowerLaw" * "Constant"
+   2019-04-05T14:18:18:  Number of parameters ......: 4
+   2019-04-05T14:18:18:  Number of spectral par's ..: 3
+   2019-04-05T14:18:18:   Prefactor ................: 2.44289795381698e-13 +/- 3.75318255643637e-15 [1e-16,1e-12] ph/cm2/s/MeV (free,scale=1e-14,gradient)
+   2019-04-05T14:18:18:   Index ....................: -3.21741632072522 +/- 0.0232805723886611 [-5,5]  (free,scale=1,gradient)
+   2019-04-05T14:18:18:   PivotEnergy ..............: 1000000 MeV (fixed,scale=1000000,gradient)
+   2019-04-05T14:18:18:  Number of temporal par's ..: 1
+   2019-04-05T14:18:18:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+   2019-04-05T14:18:18: === GCTAModelAeffBackground ===
+   2019-04-05T14:18:18:  Name ......................: bkg_23559
+   2019-04-05T14:18:18:  Instruments ...............: HESS
+   2019-04-05T14:18:18:  Instrument scale factors ..: unity
+   2019-04-05T14:18:18:  Observation identifiers ...: 23559
+   2019-04-05T14:18:18:  Model type ................: "PowerLaw" * "Constant"
+   2019-04-05T14:18:18:  Number of parameters ......: 4
+   2019-04-05T14:18:18:  Number of spectral par's ..: 3
+   2019-04-05T14:18:18:   Prefactor ................: 2.26357632227583e-13 +/- 3.56983854227473e-15 [1e-16,1e-12] ph/cm2/s/MeV (free,scale=1e-14,gradient)
+   2019-04-05T14:18:18:   Index ....................: -3.25328282552122 +/- 0.0239266148104449 [-5,5]  (free,scale=1,gradient)
+   2019-04-05T14:18:18:   PivotEnergy ..............: 1000000 MeV (fixed,scale=1000000,gradient)
+   2019-04-05T14:18:18:  Number of temporal par's ..: 1
+   2019-04-05T14:18:18:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+   2019-04-05T14:18:18: === GCTAModelAeffBackground ===
+   2019-04-05T14:18:18:  Name ......................: bkg_23592
+   2019-04-05T14:18:18:  Instruments ...............: HESS
+   2019-04-05T14:18:18:  Instrument scale factors ..: unity
+   2019-04-05T14:18:18:  Observation identifiers ...: 23592
+   2019-04-05T14:18:18:  Model type ................: "PowerLaw" * "Constant"
+   2019-04-05T14:18:18:  Number of parameters ......: 4
+   2019-04-05T14:18:18:  Number of spectral par's ..: 3
+   2019-04-05T14:18:18:   Prefactor ................: 2.87721724170504e-13 +/- 4.31204359584396e-15 [1e-16,1e-12] ph/cm2/s/MeV (free,scale=1e-14,gradient)
+   2019-04-05T14:18:18:   Index ....................: -3.25393054001656 +/- 0.0227802144778143 [-5,5]  (free,scale=1,gradient)
+   2019-04-05T14:18:18:   PivotEnergy ..............: 1000000 MeV (fixed,scale=1000000,gradient)
+   2019-04-05T14:18:18:  Number of temporal par's ..: 1
+   2019-04-05T14:18:18:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+   2019-04-05T14:18:18: === GModelSky ===
+   2019-04-05T14:18:18:  Name ......................: Crab
+   2019-04-05T14:18:18:  Instruments ...............: all
+   2019-04-05T14:18:18:  Instrument scale factors ..: unity
+   2019-04-05T14:18:18:  Observation identifiers ...: all
+   2019-04-05T14:18:18:  Model type ................: PointSource
+   2019-04-05T14:18:18:  Model components ..........: "PointSource" * "PowerLaw" * "Constant"
+   2019-04-05T14:18:18:  Number of parameters ......: 6
+   2019-04-05T14:18:18:  Number of spatial par's ...: 2
+   2019-04-05T14:18:18:   RA .......................: 83.6331 [-360,360] deg (fixed,scale=1)
+   2019-04-05T14:18:18:   DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1)
+   2019-04-05T14:18:18:  Number of spectral par's ..: 3
+   2019-04-05T14:18:18:   Prefactor ................: 1.205595290581e-15 +/- 1.02399444013702e-16 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient)
+   2019-04-05T14:18:18:   Index ....................: -2.64310065863216 +/- 0.0532926238637032 [-0,-5]  (free,scale=-1,gradient)
+   2019-04-05T14:18:18:   PivotEnergy ..............: 300000 [10000,1000000000] MeV (fixed,scale=1000000,gradient)
+   2019-04-05T14:18:18:  Number of temporal par's ..: 1
+   2019-04-05T14:18:18:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
+   2019-04-05T14:18:18: === GCTAModelIrfBackground ===
+   2019-04-05T14:18:18:  Name ......................: CTABackgroundModel
+   2019-04-05T14:18:18:  Instruments ...............: CTA
+   2019-04-05T14:18:18:  Instrument scale factors ..: unity
+   2019-04-05T14:18:18:  Observation identifiers ...: all
+   2019-04-05T14:18:18:  Model type ................: "PowerLaw" * "Constant"
+   2019-04-05T14:18:18:  Number of parameters ......: 4
+   2019-04-05T14:18:18:  Number of spectral par's ..: 3
+   2019-04-05T14:18:18:   Prefactor ................: 1 +/- 0 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
+   2019-04-05T14:18:18:   Index ....................: 0 +/- 0 [-5,5]  (free,scale=1,gradient)
+   2019-04-05T14:18:18:   PivotEnergy ..............: 1000000 [10000,1000000000] MeV (fixed,scale=1000000,gradient)
+   2019-04-05T14:18:18:  Number of temporal par's ..: 1
+   2019-04-05T14:18:18:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
 .. note::
+   Using ``inmodel=$CTOOLS/share/models/crab.xml`` will also append the CTA
+   background model to the ``models.xml`` file, yet since this model only
+   applies to CTA data it is ignored in the fit. This is signalled by
-   By default, :ref:`csfindobs` only selects data of highest quality (i.e. QUALITY=0).
-   You may overwrite this default by specifying the hidden parameter 
-   ``min_qual``. For example, ``min_qual=1`` selects all data with a 
-   looser quality criteria.
+   .. code-block:: none
+      2019-04-05T14:18:17:    Parameter "Prefactor" has zero curvature. Fix parameter.
+      2019-04-05T14:18:17:    Parameter "Index" has zero curvature. Fix parameter.
-Create an observation list
---------------------------
-
-As next step, you must convert the runlist ASCII file into an observation
-definition XML file. The observation definition XML file will contain the
-file names of all files that are needed for the analysis. You do this
-conversion with the :ref:`csiactobs` script:
-
-.. code-block:: bash
-
-  $ csiactobs
-  Data storage name [] fits-data-name
-  Input runlist file [runlist.lis]
-  Number of free parameters per background model [1] 
-  Output model definition XML file [bgmodels.xml]
-  Output observation definition XML file [obs.xml]
-
-The :ref:`csiactobs` script will create two output files: the observations
-definition XML file ``obs.xml`` and an output model definition XML file
-``bgmodels.xml``. To generate ``obs.xml``, :ref:`csiactobs` has used the
-IACT data storage and extracted the relevant file names. ``bgmodels.xml`` is
-a file that is used for background modeling, where each observation will have
-its own independent background model. In the example above, you have set the
-number of free parameters per background model to one, hence the normalisation
-of the background model for each observation will be a free parameter that is
-later adjusted by the maximum likelihood fit.
-
-There are some further hidden parameters to steer the start parameters for the
-effective area and the Gaussian background models. You may have a look at
-:ref:`csiactobs` to see the full list of available parameters.
-
-You may also note that the :ref:`csiactobs` script has created a ``csiactobs.log``
-file that logs the complete energy range of the observations that have been
-selected. These values may be important for later, in particular if you
-intend to do a stacked analysis.
-
-In case that you have already a model definition XML file that describes a model
-of the celestial source distribution (a.k.a. a sky model), you may provide this
-sky model to :ref:`csiactobs` using the hidden ``inmodel`` parameter. The
-output model definition XML file will then contain the sky model and the
-background model, and you can use the model directly in a maximum likelihood
-fit.
-
-.. code-block:: bash
-
-  $ csiactobs inmodel="mymodel.xml"
-  
-Alternatively, you can merged models using the :ref:`csmodelmerge` script:
-
-.. code-block:: bash
-
-  $ csmodelmerge
-  Input model definition XML files [mymodel.xml bgmodels.xml]
-  Output model definition XML file [combined_model.xml]
-
-.. note::
-
-   The number of files that can be merged with :ref:`csmodelmerge` is not
-   limited to two. The list of input file names may be either separated by
-   whitespace or semi-colons, can be specified using wildcards, or can be given
-   in an ASCII file (see :ref:`csmodelmerge`).
+   during the fit, and the corresponding model parameter have an error of
+   zero, although they are signalled as free.
 .. warning::
-
-   Each model component in the input model definition XML files that are provided
-   to :ref:`csmodelmerge` must have a different and unique name. Merging for
-   example the same XML model twice will raise an exception.
-
-A list of available sky models can be found
-`here <http://gammalib.sourceforge.net/user_manual/modules/model.html>`_.
-If you are not familiar with creating sky models you should read the section
-about :ref:`modelling CTA data <um_models>`.
-
-
-Example XML files
------------------
-To get familiar with the XML syntax and format, example files for an observation
-definition XML file and a model definition XML file are shown in the following.
-
-Observation definition XML file
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-.. code-block:: xml
-
-	<?xml version="1.0" encoding="UTF-8" standalone="no"?>
-	<observation_list title="observation list">
-	  <observation name="Crab Nebula" id="11111" instrument="HESS">
-	    <parameter name="EventList" file="/path/to/fits/file/events_11111.fits.gz" />
-	    <parameter name="EffectiveArea" file="/path/to/fits/file/aeff_11111.fits.gz" />
-	    <parameter name="PointSpreadFunction" file="/path/to/fits/file/psf_11111.fits.gz" />
-	    <parameter name="EnergyDispersion" file="/path/to/fits/file/edisp_11111.fits.gz" />
-	    <parameter name="Background" file="/path/to/fits/file/bgmodel_11111.fits.gz" />
-	  </observation>
-	  <observation name="Crab Nebula" id="11112" instrument="HESS">
-	    <parameter name="EventList" file="/path/to/fits/file/events_11112.fits.gz" />
-	    <parameter name="EffectiveArea" file="/path/to/fits/file/aeff_11112.fits.gz" />
-	    <parameter name="PointSpreadFunction" file="/path/to/fits/file/psf_11112.fits.gz" />
-	    <parameter name="EnergyDispersion" file="/path/to/fits/file/edisp_11112.fits.gz" />
-	    <parameter name="Background" file="/path/to/fits/file/bgmodel_11112.fits.gz" />
-	  </observation>
-	  <observation name="Crab Nebula" id="11113" instrument="HESS">
-	    <parameter name="EventList" file="/path/to/fits/file/events_11113.fits.gz" />
-	    <parameter name="EffectiveArea" file="/path/to/fits/file/aeff_11113.fits.gz" />
-	    <parameter name="PointSpreadFunction" file="/path/to/fits/file/psf_11113.fits.gz" />
-	    <parameter name="EnergyDispersion" file="/path/to/fits/file/edisp_11113.fits.gz" />
-	    <parameter name="Background" file="/path/to/fits/file/bgmodel_11113.fits.gz" />
-	  </observation>
-	  <observation name="Crab Nebula" id="11114" instrument="HESS">
-	    <parameter name="EventList" file="/path/to/fits/file/events_11114.fits.gz" />
-	    <parameter name="EffectiveArea" file="/path/to/fits/file/aeff_11114.fits.gz" />
-	    <parameter name="PointSpreadFunction" file="/path/to/fits/file/psf_11114.fits.gz" />
-	    <parameter name="EnergyDispersion" file="/path/to/fits/file/edisp_11114.fits.gz" />
-	    <parameter name="Background" file="/path/to/fits/file/bgmodel_11114.fits.gz" />
-	  </observation>
-	</observation_list>
-
-
-Model definition XML file
-^^^^^^^^^^^^^^^^^^^^^^^^^
-
-.. code-block:: xml
-
-	<?xml version="1.0" encoding="UTF-8" standalone="no"?>
-	<source_library title="source library">
-	  <source name="bkg_11111" type="CTAIrfBackground" instrument="HESS" id="11111">
-	    <spectrum type="ConstantValue">
-	      <parameter name="Normalization" value="1" error="0" scale="1" min="0.01" max="100" free="1" />
-	    </spectrum>
-	  </source>
-	  <source name="bkg_11112" type="CTAIrfBackground" instrument="HESS" id="11112">
-	    <spectrum type="ConstantValue">
-	      <parameter name="Normalization" value="1" error="0" scale="1" min="0.01" max="100" free="1" />
-	    </spectrum>
-	  </source>
-	  <source name="bkg_11113" type="CTAIrfBackground" instrument="HESS" id="11113">
-	    <spectrum type="ConstantValue">
-	      <parameter name="Normalization" value="1" error="0" scale="1" min="0.01" max="100" free="1" />
-	    </spectrum>
-	  </source>
-	  <source name="bkg_11114" type="CTAIrfBackground" instrument="HESS" id="11114">
-	    <spectrum type="ConstantValue">
-	      <parameter name="Normalization" value="1" error="0" scale="1" min="0.01" max="100" free="1" />
-	    </spectrum>
-	  </source>
-	    <source name="Crab" type="PointSource">
-	    <spectrum type="PowerLaw">
-	       <parameter name="Prefactor"   scale="1e-17" value="3.0"  min="1e-07" max="1000.0" free="1"/>
-	       <parameter name="Index"       scale="-1"    value="2.5" min="0.0"   max="+5.0"   free="1"/>
-	       <parameter name="PivotEnergy" scale="1e6"   value="1.0"  min="0.01"  max="1000.0" free="0"/>
-	    </spectrum>
-	    <spatialModel type="PointSource">
-	      <parameter name="RA"  scale="1.0" value="83.6331" min="-360" max="360" free="1"/>
-	      <parameter name="DEC" scale="1.0" value="22.0145" min="-90"  max="90"  free="1"/>
-	    </spatialModel>
-	  </source>
-	</source_library>
-
-.. note::
-
-   It is important to ensure background models are properly linked to their
-   respective observation. Therefore it is required to keep the attributes
-   ``instrument`` and ``id`` the same for the observation and the corresponding
-   background model. The tool :ref:`csiactobs` assures this automatically.
-
-
-Run ctselect
-------------
-To prepare the data for analysis, cuts have to be applied to the event data.
-The selection is performed by :ref:`ctselect`.
-This tool writes out selected event lists into the local directory. If the
-observation definition XML file contains several runs, it is recommended
-to first create a separate folder and specify this folder in the hidden
-``prefix`` argument.
-
-.. code-block:: bash
-
-  $ mkdir selected
-  
-.. code-block:: bash
-
-  $ ctselect usethres=DEFAULT usepnt=yes prefix=selected/selected_
-  Input event list or observation definition file [events.fits] obs.xml
-  Lower energy limit (TeV) [0.1]
-  Upper energy limit (TeV) [100.0]
-  Radius of ROI (degrees) (0-180) [3.0] 2.5
-  Start time (CTA MET in seconds) [0.0]
-  End time (CTA MET in seconds) [0.0]
-  Output event list or observation definition file [selected_events.fits] selected_obs.xml 
-  
-For IACT analysis, it is recommended to use the hidden parameter
-``usethres="DEFAULT"``. This instructs :ref:`ctselect` to extract the safe
-energy range from the instrument response functions and apply them to the
-data. This safe energy range is thus superior to the energy limit passed via
-the user parameters. In addition, to analyse the complete field of view, the
-parameter ``usepnt=yes`` uses, for each observation, the pointing position as
-centre for the selection radius. The radius parameter is dependent on the
-intrument, for an instrument with a 5° field of view, a radius of 2.5°
-seems reasonable. The time selection is not applied in the above example; 
-specifying 0 as start and end time skips the time selection. For time-resolved
-analysis, it is important to know the MET time that is required to extract.
-The result of the selection step is written into the observation XML file
-``selected_obs.xml``, which now contains references to the new selected event
-FITS files.
-
-
-Unbinned analysis
------------------
-Once the data is selected, the easiest way to analyse is an unbinned analysis.
-Note that the input model definition XML file must now contain the background
-models and source components to describe the field of view.
-
-.. code-block:: bash
-
-  $ ctlike
-  Event list, counts cube or observation definition file [events.fits] selected_obs.xml
-  Source model [$CTOOLS/share/models/crab.xml] crab_models.xml
-  Source model output file [crab_results.xml]
-  
-The result of the fit was stored in ``crab_results.xml``. Note that fitted
-parameters, ``Prefactors`` in particular, typically use MeV as energy unit.
-To monitor the progress of the fit on the screen, one can simply run with the
-option ``debug=yes``. Alternatively, the logfile ``ctlike.log`` can be inspected
-after the fit.
-
-On default, energy dispersion is not considered in the fit. To switch on the
-usage of the energy migration matrix, the hidden parameter ``edisp=yes`` can
-be provided. Note that this will cause a significant reduction of the
-computing speed.
-
-
-Stacked analysis
-----------------
-In a stacked analysis the events of all observations are stacked into a single
-counts cube. This analysis mode is much faster than unbinned analysis when
-having a large dataset (e.g. > 100 hours). For this type of analysis, some
-intermediate data products have to be produced. The products are a binned counts
-cube, an exposure cube, a point spead function cube, optionally an energy
-dispersion cube, and a background cube.
-
-Bin events
-^^^^^^^^^^
-
-.. code-block:: bash
-
-  $ ctbin
-  Event list or observation definition file [events.fits] selected_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] 
-  Size of the Y axis in pixels [200]
-  Algorithm for defining energy bins (FILE|LIN|LOG) [LOG]
-  Start value for first energy bin in TeV [0.1] 0.5
-  Stop value for last energy bin in TeV [100.0] 50
-  Number of energy bins [20]
-  Output counts cube [cntcube.fits]
-
-
-Create exposure cube
-^^^^^^^^^^^^^^^^^^^^
-After binning the events into a three-dimensional cube, an exposure cube has to
-be computed. The exposure is defined as the effective area times the dead-time
-corrected observation time. Each observation from the input container gets
-stacked in the resulting cube. The exposure is stored in units of :math:`cm^2 s`.
-The exposure cube does not have to contain the same binning as the event cube
-but for simplicity, the event cube can be passed to adopt the binning
-parameters. Note, however, that the exposure cube is defined in true sky
-coordinates and energy while the counts cube is defined in reconstructed sky
-coordinates and energy. Consequently, the sky area and energy range covered by
-the exposure cube should be slightly larger than that of the counts cube to
-accommodate for spill over of events due to the point spread function and energy
-dispersion. By computing the exposure cube on the same grid as the counts cube,
-the spill over of events from sources at the edge of cube will not be handled
-correctly. In the example below, however, no source at the edge of the field of
-view is present. Therefore, for simplicity, the count cube is used as input cube
-to extract the binning.
-
-This task of computing the exposure cube is is performed by :ref:`ctexpcube`. 
-
-.. code-block:: bash
-
-  $ ctexpcube
-  Event list or observation definition file [NONE] selected_obs.xml
-  Input counts cube file to extract exposure cube definition [NONE] cntcube.fits
-  Output exposure cube file [expcube.fits]
-  
-Alternatively, the exposure cube can be created with different binning than the
-counts cube:
-
-.. code-block:: bash
-
-  $ ctexpcube
-  Input event list or observation definition XML file [NONE] selected_obs.xml 
-  Input counts cube file to extract exposure cube definition [NONE] 
-  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] 0.04
-  Size of the X axis in pixels [200] 100
-  Size of the Y axis in pixels [200] 100
-  Lower energy limit (TeV) [0.5] 
-  Upper energy limit (TeV) [50.0] 
-  Number of energy bins [20] 30
-  Output exposure cube file [expcube.fits] 
-
-Create point spread function cube
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-As a next step for the stacked analysis, a cube containing the point spread
-function (PSF) must be computed. Since the PSF cannot be stored by one single
-parameter, the PSF cube computed by :ref:`ctpsfcube` has a fourth dimension.
-In each bin of the cube, the PSF is stored as a function of offset from source.
-The granularity of the PSF histogram is determined by the hidden parameter
-``anumbins`` (default: 200). Therefore, when passing the event cube to adopt
-the sky binning for the PSF cube, the resulting FITS file can become quite large
-due to the fourth dimension. Usually in IACT analysis, the PSF doesn't change
-too dramatically across the field of view. Therefore the user can think about
-reducing the spatial binning of the PSF cube:
-
-.. code-block:: bash
-
-  $ ctpsfcube
-  Input event list or observation definition XML file [NONE] selected_obs.xml 
-  Input counts cube file to extract PSF cube definition [NONE] 
-  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] 0.2
-  Size of the X axis in pixels [10] 20
-  Size of the Y axis in pixels [10] 20
-  Lower energy limit (TeV) [0.1] 0.5
-  Upper energy limit (TeV) [100.0] 50
-  Number of energy bins [20] 
-  Output PSF cube file [psfcube.fits]
-  
-Depending on the required PSF precision, one could reduce the number of offset
-bins via the hidden parameter ``anumbins``:
-
-.. code-block:: bash
-
-  $ ctpsfcube anumbins=100
-
-Create energy dispersion cube (optional)
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-In case energy dispersion should be considered used in the analysis, the tool
-:ref:`ctedispcube` computes the energy migration response for the given
-observations. Analoguous to the PSF cube, the energy dispersion cube is stored
-in a four-dimensional map. In each bin of the cube, the migration value
-``E_reco/E_true`` is stored. The granularity of the migration histogram can be
-steered via the hidden parameter ``migrabins`` (default: 100). In addition, the
-maximum migration value can be set via the hidden parameter ``migramax``
-(default: 2.0). Similar to the PSF cube, the energy dispersion cube FITS file
-can become quite large. Therefore, since the energy dispersion shouldn't vary
-too much across the sky region, one could reduce the spatial binning:
-
-.. code-block:: bash
-
-  $ ctedispcube
-  Input event list or observation definition XML file [NONE] selected_obs.xml 
-  Input counts cube file to extract energy dispersion cube definition [NONE] 
-  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] 0.2
-  Size of the X axis in pixels [10] 20
-  Size of the Y axis in pixels [10] 20
-  Lower energy limit (TeV) [0.1] 0.5
-  Upper energy limit (TeV) [100.0] 50
-  Number of energy bins [20] 
-  Output energy dispersion cube file [edispcube.fits]
-
-Create background cube
-^^^^^^^^^^^^^^^^^^^^^^
-Last but not least, for a stacked IACT analysis a cube containing the background
-rate in sky coordinates and reconstructed energy has to be computed. This task
-is performed by :ref:`ctbkgcube`. The binning here can also differ from the counts
-cube. For simplicity, however, the example below uses the counts cube to adopt
-the binning.
-
-.. code-block:: bash
-
-  $ ctbkgcube debug=yes
-  Input event list or observation definition XML file [NONE] selected_obs.xml 
-  Input counts cube file to extract background cube definition [NONE] cntcube.fits 
-  Input model XML file [NONE] crab_models.xml 
-  Output background cube file [bkgcube.fits] 
-  Output model XML file [NONE] binned_models.xml
-
-Note that this tool also requires the parameters of an input and output model.
-In the model XML file that came out of :ref:`csiactobs`, one background model
-per observation is included. This models get merged and averaged in the background
-sky cube. In the output model the background models per obseervation will be
-removed. Instead, a global background model for the newly created background
-cube is included. Sky models present in the input model XML file will also be
-included in the new XML file, which subsequently can be used for stacked
-:ref:`ctlike`.
-
-Example for stacked model XML file
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-The output model of :ref:`ctbkgcube` looks the following:
-
-.. code-block:: xml
-
-	<?xml version="1.0" encoding="UTF-8" standalone="no"?>
-	<source_library title="source library">
-	  <source name="BackgroundModel" type="CTACubeBackground" instrument="CTA,HESS,MAGIC,VERITAS">
-	    <spectrum type="PowerLaw">
-	      <parameter name="Prefactor"   value="1" error="0" scale="1" min="0.01" max="100" free="1" />
-	      <parameter name="Index"       value="0" error="0" scale="1" min="-5" max="5" free="1" />
-	      <parameter name="PivotEnergy" value="1" scale="1e+06" free="0" />
-	    </spectrum>
-	  </source>
-	  <source name="Crab" type="PointSource">
-	    <spectrum type="PowerLaw">
-	       <parameter name="Prefactor"   scale="1e-17" value="3.0"  min="1e-07" max="1000.0" free="1"/>
-	       <parameter name="Index"       scale="-1"    value="2.48" min="0.0"   max="+5.0"   free="1"/>
-	       <parameter name="PivotEnergy" scale="1e6"   value="1."  min="0.01"  max="1000.0" free="0"/>
-	    </spectrum>
-	    <spatialModel type="PointSource">
-	      <parameter name="RA"  scale="1.0" value="83.6331" min="-360" max="360" free="0"/>
-	      <parameter name="DEC" scale="1.0" value="22.0145" min="-90"  max="90"  free="0"/>
-	    </spatialModel>
-	  </source>
-	</source_library>
-	
-The background model with ``type=CTACubeBackground`` is used to scale the
-background cube stored in the FITS file created by :ref:`ctbkgcube`.
-
-Run ctlike
-^^^^^^^^^^
-Having all the intermediate data products ready, a stacked analysis can be
-conducted using :ref:`ctlike`.
-
-.. code-block:: bash
-
-  $ ctlike
-  Input event list, counts cube or observation definition XML file [events.fits] cntcube.fits 
-  Input exposure cube file (only needed for stacked analysis) [NONE] expcube.fits 
-  Input PSF cube file (only needed for stacked analysis) [NONE] psfcube.fits 
-  Input background cube file (only needed for stacked analysis) [NONE] bkgcube.fits 
-  Input model XML file [binned_models.xml]
-  Output model XML file [binned_results.xml]
-	
-Note that when passing an counts cube to :ref:`ctlike`, the tool behaves
-differently than in unbinned mode. It queries directly for the additional
-ingredients for the stacked analysis. It is important to pass the background
-model generated by :ref:`ctbkgcube` here to ensure the proper modelling of
-the background in the fit.
-
-To consider also the energy dispersion, you have to execute the :ref:`ctlike`
-tool with the hidden parameter ``edisp=yes``. The tool will then also query for
-the energy dispersion cube:
-
-.. code-block:: bash
-  
-  $ ctlike debug=yes
-  Input event list, counts cube or observation definition XML file [events.fits] cntcube.fits 
-  Input exposure cube file (only needed for stacked analysis) [NONE] expcube.fits 
-  Input PSF cube file (only needed for stacked analysis) [NONE] psfcube.fits 
-  Input energy dispersion cube file (only needed for stacked analysis) [NONE] edispcube.fits 
-  Input background cube file (only needed for stacked analysis) [NONE] bkgcube.fits 
-  Input model XML file [binned_models.xml]
-  Output model XML file [binned_results.xml] 
+   The effective area background model generated by :ref:`csiactobs` is very
+   simplistic, and should not be used for a serious data analysis of H.E.S.S.
+   data. Please follow the procedure :ref:`hess_dr1` for the generation of
+   a more reliable background model.
@@ -3,66 +3,162 @@
 Copying IACT data
 ==================
+  .. admonition:: What you will learn
+
+     You will learn how copy data from an IACT data store to your disk.
+
 Before you start
 ----------------
-Currently the way how IACT data is stored is quite prototypical. The ctools can work
-with IACT data that is stored following the `specifications here <http://gamma-astro-data-formats.readthedocs.org/en/latest/data_storage/index.html>`_.
-The script :ref:`csiactcopy` can be used to transfer these data from one system to another.
-IACT data is the property of collaborations (such as HESS, MAGIC or VERITAS) and not publicly available. If you are member of one of these
-collaborations, you should have access to a central place where the FITS data is stored and maintained.
-
-Mount remote file system
-------------------------
-Note that :ref:`csiactcopy` does not include functionality to establish a conection to a remote file system.
-For the moment, it is capable of copying data from one location to another. Therefore it is necessary to
-mount the remote file system where the data is located on the local machine. This can be achieved using e.g.
-`sshfs <https://github.com/libfuse/sshfs>`_:
+
+Currently the way how IACT data is stored is quite prototypical.
+The ctools can work with IACT data that is stored following the specifications
+provided in the
+`Data formats for gamma-ray astronomy <http://gamma-astro-data-formats.readthedocs.org/en/latest/data_storage/index.html>`_.
+The script :ref:`csiactcopy` can be used to transfer IACT data from one system
+to another.
+IACT data is the property of collaborations (such as HESS, MAGIC or VERITAS)
+and not publicly available, except for the recent
+`public release of some H.E.S.S. data <https://www.mpi-hd.mpg.de/hfm/HESS/pages/dl3-dr1>`_.
+In this tutorial we use these public data.
+After downloading the data and uncompression the tarball you should have a
+folder named ``hess_dl3_dr1`` in your current workingcdirectory with the
+following structure:
 .. code-block:: bash
-  $ sshfs user@remote.server.com:/path/to/remote/fits/data/ /path/to/mountpoint/
+   hess_dl3_dr1/
+   hess_dl3_dr1/data
+   hess_dl3_dr1/hdu-index.fits.gz
+   hess_dl3_dr1/obs-index.fits.gz
+   hess_dl3_dr1/README.txt
+
+You have to add the file ``master.json`` to this folder with the following
+content:
+
+.. code-block:: json
+
+   {
+     "datasets": [
+       {
+         "name": "hess_dl3_dr1",
+         "hduindx": "hdu-index.fits.gz",
+         "obsindx": "obs-index.fits.gz",
+         "Total observations": 106,
+         "Total time": "48.6 hours",
+         "Total events": 662826,
+         "Disk usage": "43 MB",
+         "Date created": "2018-09-10"
+       }
+     ]
+   }
+
+.. note::
+   If you are member of an IACT collaboration and if the collaboration has
+   setup a datastore in the
+   `compliant format <http://gamma-astro-data-formats.readthedocs.org/en/latest/data_storage/index.html>`_
+   you can also establish a remote connection to this datastore using for
+   example
+
+   .. code-block:: bash
+
+      $ sshfs user@remote.server.com:/path/to/remote/fits/data/ /project-data/hess
+
+   The success of this operation can be verified by checking if there is a file
+   called ``master.json`` in the mounted directory.
+
+
+Check what data are available
+-----------------------------
+
+Before you start an analysis, it is important that you know what kind of data
+are available for analysis.
+In the `current scheme <http://gamma-astro-data-formats.readthedocs.org/en/latest/data_storage/index.html>`_
+every dataset is represented by a unique name, in our case ``hess_dl3_dr1``.
+You have to pass this string to some set of scripts to find and use the data
+for analysis.
+
+To check what datasets are available you should run the :ref:`csiactdata`
+script:
+
+.. code-block:: bash
+
+   $ csiactdata debug=yes
+   Path were data is located [] /project-data/hess/hess_dl3_dr1
+
+This script will list the names of available datasets. Here at excerpt of the
+log file ``csiactdata.log``:
+
+.. code-block:: none
+
+   2019-04-05T13:10:35: +======================+
+   2019-04-05T13:10:35: | Data storage entries |
+   2019-04-05T13:10:35: +======================+
+   2019-04-05T13:10:35:  Master index file .........: /project-data/hess/hess_dl3_dr1/master.json
+   2019-04-05T13:10:35:
+   2019-04-05T13:10:35: +------------------------+
+   2019-04-05T13:10:35: | Available data configs |
+   2019-04-05T13:10:35: +------------------------+
+   2019-04-05T13:10:35: === hess_dl3_dr1 ===
+   2019-04-05T13:10:35:  Name ......................: hess_dl3_dr1
+   2019-04-05T13:10:35:  Observation index .........: obs-index.fits.gz
+   2019-04-05T13:10:35:  HDU index .................: hdu-index.fits.gz
+   2019-04-05T13:10:35:  Date created ..............: 2018-09-10
+   2019-04-05T13:10:35:  Total events ..............: 662826
+   2019-04-05T13:10:35:  Total time ................: 48.6 hours
+   2019-04-05T13:10:35:  Disk usage ................: 43 MB
+   2019-04-05T13:10:35:  Total observations ........: 105
-The success of this operation can be verified by checking if there is a file called ``master.json``
-in the mounted directory.
 Copy data 
 ---------
-After mounting the file system to a folder of your choice (``/path/to/mountpoint/`` in this case), the copying of the data can be initiated.
-To monitor the progress on the screen, the hidden parameter ``debug=yes`` can be used. To increase the chattiness of the output,
-the hidden parameter ``chatter`` can e.g. be set to 3.
-.. code-block:: bash
-  
-  $ csiactcopy debug=yes chatter=3
-  Location of remote master file [/path/to/mountpoint/master.json] 
-  Name of FITS production to download [iact-fits]
-  Destination path of FITS data [/path/to/local/fits/data/] 
-  
-The name of the FITS production must be available in the remote 
-``master.json`` file. In case the names of the available productions are unkown
-there are two options to determine an appropriate name:
-
-*  Run :ref:`csiactdata` and set as ``datapath`` the mountpoint
-*  Try any name, :ref:`csiactcopy` will throw a ``RuntimeError`` with the available production names if the input name is not found
-
-Copy only a subset of the data
-------------------------------
-If the user is interested only in a small subset of the remote FITS production,
-it is also possible to download only parts of the data. The user simply has to provide an
-ASCII file containing a list of observation IDs to be downloaded (one ID per line). The hidden
-parameter ``runlist`` then needs to be specified.
+You can now copy the data from the data store to your disk. In this example we
+assume that the H.E.S.S. datastore is at ``/project-data/hess/hess_dl3_dr1``.
+To do the copy you should type
 .. code-block:: bash
-  
-  $ csiactcopy runlist=myrunlist.txt debug=yes
-  
+
+   $ csiactcopy
+   Location of remote master file [/path/to/mountpoint/master.json] /project-data/hess/hess_dl3_dr1/master.json
+   Name of FITS production to download [iact-fits] hess_dl3_dr1
+   Destination path of FITS data [/path/to/local/fits/data/] data
+
+.. note::
+   To monitor the progress on the screen you should set the hidden parameter
+   ``debug=yes``. You may also increase the chattiness of the output by
+   increasing the hidden parameter ``chatter`` from its default value of 2
+   to for example 3.
+
+.. note::
+   If you are only interested in a subset of the dataset you can also download
+   only that subset. For that purpose you have to specify an ASCII file
+   containing a list of observation identifiers to be downloaded, e.g.
+
+   .. code-block:: none
+
+      20136
+      20137
+      20151
+
+   You can then run :ref:`csiactcopy` by specifying the ASCII filename
+   for the hidden ``runlist`` parameter:
+
+   .. code-block:: bash
+
+      $ csiactcopy runlist=runlist.txt
+      Location of remote master file [/project-data/hess/hess_dl3_dr1/master.json]
+      Name of FITS production to download [hess_dl3_dr1]
+      Destination path of FITS data [data]
+
+
 Troubleshooting
 ---------------
-In case a download was interupted (e.g. the remote file system was disconnected), the script can simply be rerun.
-The index files which list the files available on the system are refreshed after the copying process.
-The hidden boolean parameter ``clobber`` (default=yes) specifies if the user wishes that local content should be overwritten with
-remote content. In turn, specifying ``clobber=no`` can speed up the process since remote content does not need to be copied
-if it is already available on the user machine.
-
-To verify the copying process was successful, it is recommended to run :ref:`csiactdata` after the process in order to list available productions
-on the local machine. The output should now include the copied FITS production name.
+In case a download was interupted (e.g. the remote file system was
+disconnected), the script can simply be rerun.
+The index files, which list the files available on the system, are refreshed
+after the copying process.
+The hidden boolean parameter ``clobber`` (by default ``yes``) specifies if the
+user wishes that local content should be overwritten with remote content.
+In turn, specifying ``clobber=no`` can speed up the process since remote
+content does not need to be copied if it is already available on the user
+machine.
-.. _sec_iact_high_level:
-
-High level analysis tools for IACT data
-=======================================
-
-This section describes and explains how to create high level output. The
-following scripts and tools are explained:
-
-* :ref:`csobsinfo` (inspect and observation XML files)
-* :ref:`csmodelinfo` (inspect model XML files)
-* :ref:`ctulimit` (compute upper limits)
-* :ref:`cterror` (compute asymmetric errors)
-* :ref:`csspec` (compute spectral points)
-* :ref:`cslightcrv` (compute a light curve)
-* :ref:`csresmap` (compute a residual map)
-* :ref:`cttsmap` (compute a test statistics map)
-
-Inspecting observation definition files
----------------------------------------
-After having created an observation definition file, the user might want to
-know a summary of its content. Important quantities usually are the number of
-observations, the total livetime, or the mean zenith angle. For this purpose,
-:ref:`csobsinfo` inspects an observation XML container for its content.
-
-.. code-block:: bash
-
-  $ csobsinfo
-  Event list, counts cube, or observation definition file [obs.xml] selected_obs.xml 
-  
-This dumps general information of the observation container in a logfile ``csobsinfo.log``.
-For a direct output on the screen, pass the hidden parameter ``debug=yes`` on the command line:
-
-.. code-block:: bash
-
-  $ csobsinfo debug=yes
-  Event list, counts cube, or observation definition file [obs.xml] selected_obs.xml 
-  ...
-  2016-02-19T15:43:08: +=========+
-  2016-02-19T15:43:08: | Summary |
-  2016-02-19T15:43:08: +=========+
-  2016-02-19T15:43:08:  Unbinned observations .....: 4
-  2016-02-19T15:43:08:  Number of events ..........: 6085
-  2016-02-19T15:43:08:  Binned observations .......: 0
-  2016-02-19T15:43:08:  Number of bins ............: 0
-  2016-02-19T15:43:08: 
-  2016-02-19T15:43:08: === Pointings ===
-  2016-02-19T15:43:08:  Mean zenith angle .........: 46.89
-  2016-02-19T15:43:08:  Mean azimuth angle ........: 13.63
-  2016-02-19T15:43:08: 
-  2016-02-19T15:43:08: === Energy range ===
-  2016-02-19T15:43:08:  Emin ......................: 461.48 GeV
-  2016-02-19T15:43:08:  Emax ......................: 37.8832 TeV
-  2016-02-19T15:43:08: 
-  2016-02-19T15:43:08: === Time range ===
-  2016-02-19T15:43:08:  Start [MJD] ...............: 53343.9212153
-  2016-02-19T15:43:08:  Stop [MJD] ................: 53347.9321065
-  2016-02-19T15:43:08:  Total ontime ..............: 6877.00 s = 114.62 min = 1.91 h
-  2016-02-19T15:43:08:  Total livetime ............: 6427.52 s = 107.13 min = 1.79 h
-  ... 
-
-Further options are described on the reference page :ref:`csobsinfo`.
-
-Inspecting model XML files
---------------------------
-Similar to the observation containers, model XML file can be inspected, too.
-The tool :ref:`csmodelinfo` gives a summary of the model container. In particular,
-the number of free parameters, number of sky and background models and also number of
-parameters that are at their limits might be of interest.
-
-.. code-block:: bash
-
-  $ csmodelinfo debug=yes
-  Input model definiton file XML file [$CTOOLS/share/models/crab.xml]
-  Output DS9 region file [ds9.reg]
-  2016-02-19T15:50:40: +=========+
-  2016-02-19T15:50:40: | Summary |
-  2016-02-19T15:50:40: +=========+
-  2016-02-19T15:50:40: === Instrument specific models ===
-  2016-02-19T15:50:40:  All .......................: 1
-  2016-02-19T15:50:40:  CTA .......................: 1
-  2016-02-19T15:50:40: === Model types ===
-  2016-02-19T15:50:40:  PointSource ...............: 1
-  2016-02-19T15:50:40:  CTAIrfBackground ..........: 1
-  2016-02-19T15:50:40: 
-  2016-02-19T15:50:40: +=======================+
-  2016-02-19T15:50:40: | Parameter information |
-  2016-02-19T15:50:40: +=======================+
-  2016-02-19T15:50:40:  All parameters ............: 10
-  2016-02-19T15:50:40:  Fixed parameters ..........: 6
-  2016-02-19T15:50:40:  Free parameters (total) ...: 4
-  2016-02-19T15:50:40:  Free background parameters : 2
-  2016-02-19T15:50:40:  Free source parameters ....: 2
-  2016-02-19T15:50:40:  Free spectral parameters ..: 2
-  2016-02-19T15:50:40:  Parameters at limit .......: 0
-
-Similar to :ref:`csobsinfo`, the script also generates a ds9 region file
-including the positions of all the sky model components.
-
-Compute upper limit
--------------------
-Very often in gamma-ray astronomy sources are at the verge of detection or even not detectable.
-In such cases, it is useful to derive an upper limit using :ref:`ctulimit`.
-
-.. code-block:: bash
-
-	$ ctulimit
-	Input event list, counts cube or observation definition XML file [selected_obs.xml]
-	Input model XML file [crab_models.xml]  
-	Source of interest [Crab] 
-	
-The upper limit will be stored in the log file. To get the limit printed on screen, use the hidden parameter
-``debug=yes``.
-
-Compute asymmetric errors
--------------------------
-When an analysis approaches its final state, it makes sense to have asymmetric errors on the parameters of the
-source of interest. For this purpose, the tool :ref:`cterror` can be used:
-
-.. code-block:: bash
-
-	$ cterror debug=yes
-	Input event list, counts cube or observation definition XML file [selected_obs.xml] 
-	Input model XML file [crab_models.xml] 
-	Source of interest [Crab] 
-	Output model XML file [cterror_results.xml] 
-
-The output model does not contain asymmetric errors yet. The positive and negative uncertainties can be read from the
-logfile (or from screen if ``debug=yes`` was specified).
-
-Compute spectral points
------------------------
-A very common task in astronomy is to compute spectral data points. To determine a spectral data point, a small energy
-range is considered and the model prefactor and its uncertainty is evaluated. The tool :ref:`csspec` works on both, binned
-and unbinned data. The hidden parameter ``edisp=yes`` can be specified in both cases to consider the energy migration
-matrix in the fit.
-
-.. note::
-  
-  The resulting spectral points are provided as a function of reconstructed energy. 
-
-Unbinned
-^^^^^^^^
-
-.. code-block:: bash
-
-	$ csspec debug=yes
-	Input event list, counts cube, or observation definition XML file [events.fits] selected_obs.xml 
-	Input model XML file [$CTOOLS/share/models/crab.xml] crab_models.xml 
-	Source name [Crab] 
-	Algorithm for defining energy bins (FILE|LIN|LOG) [LOG] 
-	Lower energy limit for spectral points (TeV) [0.1] 0.5
-	Upper energy limit for spectral points (TeV) [100.0] 50.0
-	Number of spectral points (1-10000) [20] 10
-	Output spectrum file [spectrum.fits] 
-	
-Binned
-^^^^^^
-
-.. code-block:: bash
-
-	$ csspec debug=yes
-	Input event list, counts cube, or observation definition XML file [cntcube.fits] 
-	Input exposure cube file (only needed for stacked analysis) [expcube.fits] 
-	Input PSF cube file (only needed for stacked analysis) [psfcube.fits] 
-	Input background cube file (only needed for stacked analysis) [bkgcube.fits] 
-	Input model XML file [binned_models.xml] 
-	Source name [Crab] 
-	Number of spectral points (1-10000) [10] 
-	Lower energy limit for spectral points (TeV) [0.5] 
-	Upper energy limit for spectral points (TeV) [50.0] 
-	Output spectrum file [spectrum.fits]  
-
-Plot spectral points
-^^^^^^^^^^^^^^^^^^^^
-Instead of supporting a plotting library, simple example scripts are available to visualise data products.
-Have a look at ``$CTOOLS/examples/show_spectrum.py`` how a spectrum can be plotted using ``matplotlib``.
-To have first glance at the above computed spectrum one can use this script in the following way:
-
-.. code-block:: bash
-  
-  $ python $CTOOLS/examples/show_spectrum.py spectrum.fits
-
-
-Compute light curves
---------------------
-Euqally common in astronomy are light curves from a time-series analysis. The tool :ref:`cslightcrv` takes as input only
-an unbinned observation container (or a single event list). Counts cube are not possible as inout since the time informstion
-is lost during the binning procedure. Nevertheless, in case of large time bins, it is possible to require a binned analysis
-in each time bin. The tool will accordingly slice the event list and create the data products for the stacked analysis in each
-time span.
-
-.. code-block:: bash
-
-	$ cslightcrv
-	Input event list, counts cube, or observation definition XML file [selected_obs.xml]
-	Input model XML file [crab_models.xml]  
-	Source name [Crab] 
-	Algorithm for defining time bins (FILE|LIN|GTI) [GTI] 
-	Number of energy bins per light curve bin (0=unbinned) [0] 
-	Lower energy limit of events (TeV) [0.1] 0.5
-	Upper energy limit of events (TeV) [100.0] 50
-	Output light curve file [lightcurve.fits]  
-
-Compute a residual map
-----------------------
-A frequent means to visually inspect the fitted model with respect to the input data is to create residual maps.
-Using :ref:`csresmap` such a map can easily be computed. The tool internally bins the data according to user parameters.
-(if data is not already provided in a binned state). Taking into account the corresponding IRFs,
-:reF:`csresmap` computes a model map (running :ref:`ctmodel`). Subsequently, there are several choices how data
-and model should be compared. There are three options:
-
-* ``SUB``: the subtraction of the model from the counts. The resulting map will display differences in absolute counts.
-* ``SUBDIV``: the subtraction and division by the model. The resulting map will display relative differences of the data with respect to the model.
-* ``SUBDIVSQRT``: the subtraction and division by the square root of the model. The resulting map will display an approximation of a residual significance map. In case of sufficient count statistic per bin, the bin value represents the significance. 
-
-Example unbinned observation container
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-.. code-block:: bash
-
-	$ csresmap
-	Input event list, counts cube, or observation definition XML file [selected_obs.xml]
-	Input model XML file [crab_results.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] 
-	Coordinate System (CEL|GAL) [CEL] 
-	Projection method (AIT|AZP|CAR|MER|MOL|STG|TAN) [CAR] 
-	Size of the X axis in pixels [100]
-	Size of the Y axis in pixels [100]
-	Pixel size (deg/pixel) [0.02]
-	Residual map computation algorithm (SUB|SUBDIV|SUBDIVSQRT) [SUBDIV]
-	Output residual map file [resmap.fits] 
-
-Example binned/stacked observation container
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-The tool behaves differently if the input observation parameter is a count cube. It will query for the additional binned reponse components.
-
-.. code-block:: bash
-
-	$ csresmap 
-	Input event list, counts cube, or observation definition XML file [cntcube.fits]
-	Input model cube file (generated with ctmodel) [NONE] 
-	Input exposure cube file (only needed for stacked analysis) [expcube.fits]
-	Input PSF cube file (only needed for stacked analysis) [psfcube.fits]  
-	Input background cube file (only needed for stacked analysis) [bkgcube.fits] 
-	Input model XML file [binned_results.xml] 
-	Residual map computation algorithm (SUB|SUBDIV|SUBDIVSQRT) [SUBDIV] 
-	Output residual map file [resmap.fits] 
-	
-Note that the tool also queries for an input model cube file (which we set to ``NONE`` here). This is very convenient in case the model cube has already been precomputed
-using :ref:`ctmodel`.
-
-Example with input model cube
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-In case a model cube can be provided, :ref:`csresmap` simply collapses both cube into one skymap and applies the specified algorithm.
-
-.. code-block:: bash
-
-	$ csresmap
-	Input event list, counts cube, or observation definition XML file [cntcube.fits]
-	Input model cube file (generated with ctmodel) [modcube.fits] 
-	Residual map computation algorithm (SUB|SUBDIV|SUBDIVSQRT) [SUB] 
-	Output residual map file [resmap.fits] 
-
-
-Compute a test statistics (TS) map
-----------------------------------
-The test statistics quantity is a more precise measurement of statistic than the ``SUBDIVSQRT`` option in :ref:`csresmap`. 
-To get familiar with the concept of this procedure, :ref:`read more about TS maps <howto_ts>`.
-
-.. code-block:: bash
-
-	$ cttsmap
-	Input event list, counts cube or observation definition XML file [selected_obs.xml]  
-	Input model XML file [crab_results.xml] 
-	Test source name [Crab] 
-	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.05]
-	Size of the X axis in pixels [20]
-	Size of the Y axis in pixels [20]
-	Output Test Statistic map file [tsmap.fits] 
-
-The output file will contain a sky map holding the TS map. In complementary extensions, for each free parameter of the source model, a map displaying the parameter
-value is added to the output file. For instance a map of the prefactor can be visualised with `DS9 <http://ds9.si.edu/>`_ as follows:
-
-.. code-block:: bash
-  
-  $ ds9 tsmap.fits[Prefactor] 
-
-Since :ref:`cttsmap` has to run a parameter optimisation in each skymap bin, it is very time comsuming to
-compute a fine-granulated TS sky map. Read here how to split up the :ref:`computation into several jobs <sec_iact_tips>`. 
-Similar as many other tools :ref:`cttsmap` can also work on binned observation input:
-
-.. code-block:: bash
-  
-	$ cttsmap debug=yes
-	Input event list, counts cube or observation definition XML file [cntcube.fits]
-	Input exposure cube file (only needed for stacked analysis) [expcube.fits] 
-	Input PSF cube file (only needed for stacked analysis) [psfcube.fits] 
-	Input background cube file (only needed for stacked analysis) [bkgcube.fits] 
-	Input model XML file [binned_results.xml] 
-	Test source name [Crab] 
-	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.05] 
-	Size of the X axis in pixels [20] 
-	Size of the Y axis in pixels [20] 
-	Output Test Statistic map file [tsmap.fits] 
-
-	
 .. _sec_iact:
-Analysing current IACT data
-===========================
+Analysing data from existing IACTs
+==================================
+
+This tutorial is specifically written for the analysis of H.E.S.S. data but
+it may also be applicable to VERITAS and MAGIC, provided that these instruments
+set up a data store that is compliant with the format described in
+`Data formats for gamma-ray astronomy <http://gamma-astro-data-formats.readthedocs.org/en/latest/data_storage/index.html>`_.
 .. toctree::
    :maxdepth: 1
    copy
+   select
    analysis
-   high_level
-   simulation
-   tips
+.. _sec_iact_selection:
+
+Selecting IACT observations
+===========================
+
+  .. admonition:: What you will learn
+
+     You will learn how to select observations from a dataset.
+
+As next step you have to select observations from a dataset and to build
+an observation definition XML file for the analysis.
+
+
+Selecting observations
+----------------------
+
+First run :ref:`csfindobs` to select the observations of interest.
+The most common search for observations is by pointing direction on the sky,
+and you do this by typing
+
+.. code-block:: bash
+
+   $ csfindobs
+   Path were data are located (NONE uses $VHEFITS environment variable) [NONE] data
+   Name of FITS production (Run csiactdata to view your options) [prod-name] hess_dl3_dr1
+   Right ascension of selection region centre (deg) [83.63]
+   Declination of selection region centre (deg) [22.51]
+   Search radius of selection region (deg) [2.5]
+   Output runlist file [runlist.lis]
+
+This creates an output ASCII file named ``runlist.lis`` that contains the
+identifiers for all observations with pointing directions located within 2.5 deg
+around the position of the Crab nebula. The content of the file is
+
+.. code-block:: none
+
+   23523
+   23526
+   23559
+   23592
+
+hence four observations were found in the production. An excerpt of the log
+file ``csfindobs.log`` of the run is shown below:
+
+.. code-block:: none
+
+   2019-04-05T13:22:08: +===================+
+   2019-04-05T13:22:08: | Find observations |
+   2019-04-05T13:22:08: +===================+
+   2019-04-05T13:22:08:  Expression ................: ANGSEP(83.63,22.51,RA_PNT,DEC_PNT)<=2.5&&QUALITY<=0
+   2019-04-05T13:22:08:  Observations ..............: 4
+   2019-04-05T13:22:08:  Observation 1 .............: 23523
+   2019-04-05T13:22:08:  Observation 2 .............: 23526
+   2019-04-05T13:22:08:  Observation 3 .............: 23559
+   2019-04-05T13:22:08:  Observation 4 .............: 23592
+
+.. note::
+   You may add further selection criteria using the hidden ``expression``
+   parameter.
+   You may use any FITS selection expression that is
+   `supported by the cfitsio row filtering specification <https://heasarc.gsfc.nasa.gov/docs/software/fitsio/c/c_user/node97.html>`_.
+   Column names that can be used in the expression can be found
+   `here <http://gamma-astro-data-formats.readthedocs.org/en/latest/data_storage/obs_index/index.html>`_
+   or by simply browsing the observation index file with a FITS viewer.
+   For example
+
+   .. code-block:: bash
+
+      $ csfindobs expression="ZEN_PNT < 30 && LIVETIME > 1600"
+      Path were data are located (NONE uses $VHEFITS environment variable) [data] 
+      Name of FITS production (Run csiactdata to view your options) [hess_dl3_dr1]
+      Right ascension of selection region centre (deg) [83.63] NONE
+      Output runlist file [runlist.lis]
+
+   will select all observations that have a zenith angle less than 30 deg and a
+   livetime larger than 1600 seconds. Note that pointing selection was omitted
+   by specifying ``NONE`` for the Right Ascension.
+
+.. note::
+   By default, :ref:`csfindobs` only selects data of highest quality
+   (i.e. ``QUALITY=0``).
+   You may overwrite this default by specifying the hidden parameter 
+   ``min_qual``. For example, ``min_qual=1`` selects all data with a 
+   looser quality criteria.
+
+
+Create an observation list
+--------------------------
+
+Now you can convert the runlist into an observation definition XML file.
+You do this conversion with the :ref:`csiactobs` script:
+
+.. code-block:: bash
+
+   $ csiactobs
+   Path were data are located (NONE uses $VHEFITS environment variable) [NONE] data
+   Data storage name [fits-prod-name] hess_dl3_dr1
+   Input runlist file [runlist.lis]
+   Number of free parameters per background model [1] 2
+   Output model definition XML file [bkgmodels.xml]
+   Output observation definition XML file [obs.xml]
+
+The :ref:`csiactobs` script will create two output files: the observation
+definition XML file ``obs.xml`` and an output model definition XML file
+``bkgmodels.xml``. To generate ``obs.xml``, :ref:`csiactobs` has used the
+IACT data storage and extracted the relevant file names. ``bgmodels.xml`` is
+a file that is used for background modeling, where each observation will have
+its own independent background model. In the example above, you have set the
+number of free parameters per background model to one, hence the normalisation
+of the background model for each observation will be a free parameter that is
+later adjusted by the maximum likelihood fit.
+
+Here the observation definition file that was created
+
+.. code-block:: xml
+
+   <?xml version="1.0" encoding="UTF-8" standalone="no"?>
+   <observation_list title="observation list">
+     <observation name="Crab Nebula" id="23523" instrument="HESS">
+       <parameter name="EventList" file="data/data/hess_dl3_dr1_obs_id_023523.fits.gz[events]" />
+       <parameter name="EffectiveArea" file="data/data/hess_dl3_dr1_obs_id_023523.fits.gz[aeff]" />
+       <parameter name="PointSpreadFunction" file="data/data/hess_dl3_dr1_obs_id_023523.fits.gz[psf]" />
+       <parameter name="EnergyDispersion" file="data/data/hess_dl3_dr1_obs_id_023523.fits.gz[edisp]" />
+       <parameter name="Background" file="" />
+     </observation>
+     <observation name="Crab Nebula" id="23526" instrument="HESS">
+       <parameter name="EventList" file="data/data/hess_dl3_dr1_obs_id_023526.fits.gz[events]" />
+       <parameter name="EffectiveArea" file="data/data/hess_dl3_dr1_obs_id_023526.fits.gz[aeff]" />
+       <parameter name="PointSpreadFunction" file="data/data/hess_dl3_dr1_obs_id_023526.fits.gz[psf]" />
+       <parameter name="EnergyDispersion" file="data/data/hess_dl3_dr1_obs_id_023526.fits.gz[edisp]" />
+       <parameter name="Background" file="" />
+     </observation>
+     <observation name="Crab Nebula" id="23559" instrument="HESS">
+       <parameter name="EventList" file="data/data/hess_dl3_dr1_obs_id_023559.fits.gz[events]" />
+       <parameter name="EffectiveArea" file="data/data/hess_dl3_dr1_obs_id_023559.fits.gz[aeff]" />
+       <parameter name="PointSpreadFunction" file="data/data/hess_dl3_dr1_obs_id_023559.fits.gz[psf]" />
+       <parameter name="EnergyDispersion" file="data/data/hess_dl3_dr1_obs_id_023559.fits.gz[edisp]" />
+       <parameter name="Background" file="" />
+     </observation>
+     <observation name="Crab Nebula" id="23592" instrument="HESS">
+       <parameter name="EventList" file="data/data/hess_dl3_dr1_obs_id_023592.fits.gz[events]" />
+       <parameter name="EffectiveArea" file="data/data/hess_dl3_dr1_obs_id_023592.fits.gz[aeff]" />
+       <parameter name="PointSpreadFunction" file="data/data/hess_dl3_dr1_obs_id_023592.fits.gz[psf]" />
+       <parameter name="EnergyDispersion" file="data/data/hess_dl3_dr1_obs_id_023592.fits.gz[edisp]" />
+       <parameter name="Background" file="" />
+     </observation>
+   </observation_list>
+   <observation_list title="observation list" />
+
+and here the model definition file
+
+.. code-block:: xml
+
+  <?xml version="1.0" encoding="UTF-8" standalone="no"?>
+  <source_library title="source library">
+    <source name="bkg_23523" type="CTAAeffBackground" instrument="HESS" id="23523" tscalc="0">
+      <spectrum type="PowerLaw">
+        <parameter name="Prefactor" value="1" error="0" scale="1e-14" min="0.01" max="100" free="1" />
+        <parameter name="Index" value="-2" scale="1" min="-5" max="5" free="1" />
+        <parameter name="PivotEnergy" value="1" scale="1000000" free="0" />
+      </spectrum>
+    </source>
+    <source name="bkg_23526" type="CTAAeffBackground" instrument="HESS" id="23526" tscalc="0">
+      <spectrum type="PowerLaw">
+        <parameter name="Prefactor" value="1" error="0" scale="1e-14" min="0.01" max="100" free="1" />
+        <parameter name="Index" value="-2" scale="1" min="-5" max="5" free="1" />
+        <parameter name="PivotEnergy" value="1" scale="1000000" free="0" />
+      </spectrum>
+    </source>
+    <source name="bkg_23559" type="CTAAeffBackground" instrument="HESS" id="23559" tscalc="0">
+      <spectrum type="PowerLaw">
+        <parameter name="Prefactor" value="1" error="0" scale="1e-14" min="0.01" max="100" free="1" />
+        <parameter name="Index" value="-2" scale="1" min="-5" max="5" free="1" />
+        <parameter name="PivotEnergy" value="1" scale="1000000" free="0" />
+      </spectrum>
+    </source>
+    <source name="bkg_23592" type="CTAAeffBackground" instrument="HESS" id="23592" tscalc="0">
+      <spectrum type="PowerLaw">
+        <parameter name="Prefactor" value="1" error="0" scale="1e-14" min="0.01" max="100" free="1" />
+        <parameter name="Index" value="-2" scale="1" min="-5" max="5" free="1" />
+        <parameter name="PivotEnergy" value="1" scale="1000000" free="0" />
+      </spectrum>
+    </source>
+  </source_library>
+
+.. note::
+   In case that you have already a model definition XML file that describes a
+   model of the celestial source distribution (a.k.a. a sky model), you may
+   provide this sky model to :ref:`csiactobs` using the hidden ``inmodel``
+   parameter. The output model definition XML file will then contain the sky
+   model and the background model, and you can use the model directly in a
+   maximum likelihood fit.
+
+   .. code-block:: bash
+
+      $ csiactobs inmodel=$CTOOLS/share/models/crab.xml
+      Path were data are located (NONE uses $VHEFITS environment variable) [data]
+      Data storage name [hess_dl3_dr1]
+      Input runlist file [runlist.lis]
+      Number of free parameters per background model [2]
+      Output model definition XML file [bkgmodels.xml] models.xml
+      Output observation definition XML file [obs.xml]
+
+
+Select events
+-------------
+
+The final step needed is the selection of events, to define in particular
+the region of interest and energy range for the analysis.
+You do this by typing
+
+.. code-block:: bash
+
+   $ ctselect
+   Input event list or observation definition XML file [events.fits] obs.xml
+   Radius of ROI around pointing or specified RA/DEC (degrees) (0-180) [3.0] 2.0
+   Start time (UTC string, JD, MJD or MET in seconds) [NONE]
+   Lower energy limit (TeV) [0.1] 0.5
+   Upper energy limit (TeV) [100.0] 10.0
+   Output event list or observation definition XML file [selected_events.fits] obs_selected.xml
+
+which will select all events within a radius of 2 degrees around the pointing
+direction with energies between 500 GeV and 10 TeV.
+
-.. _sec_iact_simulation:
-
-Simulating IACT data
-====================
-This section describes how to run simulations with current IACT data. In particular, a guide how to simulate data for a given set of IRFs, or a given set of
-observations is provided. In addition, tools to compute sensitivity curves or pull distributions are presented.
-All tools are based on the event list simulator :ref:`ctobssim`.
-
-Simulations for a specific observation container 
-------------------------------------------------
-On default, :ref:`ctobssim` builds a new observation container by querying relevant user parameters. In case
-an observation XML container is already in place, the hidden parameter ``inobs`` can be specified (:ref:`see here how to create an
-observation container <sec_iact_analysis>`). It might be useful to run :ref:`ctselect` beforehand to ensure that the simulation knows about
-energy thresholds, field of view, etc. If an unselected observation container is used, the simulation tool will query for global energy
-boundaries and a global RoI radius. The following example illustrates the usage of an already selected observation container.
-
-.. code-block:: bash
- 
-	$ ctobssim inobs=selected_obs.xml
-	Input model XML file [crab_models.xml] 
-	Output event data file or observation definition XML file [sim_obs.xml] 
-	
-Note that :ref:`ctobssim` writes a new event list file per observation into the current directory. If the observation container
-is large, it might be convenient to create a new folder and specify the ``prefix`` keyword (similar to :ref:`ctselect`):
-
-.. code-block:: bash
-  
-  $ mkdir simulated
-  $ ctobssim inobs=selected_obs.xml prefix=simulated/sim_events_
-  ...
-  
-The newly created observation container can be used in the same way as a normal observation container, i.e. it can be analysed with
-:ref:`ctlike` etc. 
-
-.. note::
-
-	On default, :ref:`ctobssim` uses the same seed in every execution, i.e. the tool will always return the same simulated event list, regardless
-	when, where and how it was executed. To steer this the user can provide the hidden integer parameter ``seed``:
-	
-	.. code-block:: bash
-	
-		$ ctobssim inobs=selected_obs.xml seed=2
-  
-  This might be important if simulations are run in parallel.
-
-Simulating observation time for a given set of IRFs
----------------------------------------------------
-There are other use cases to run simulations, for instance when there is no observation container available yet. This might be the case e.g
-for an observation proposal: one would want to use a certain set of IRFs to simulate e.g. 50 hours of data to check whether a source
-of interest can be detected. For this purpose it is necessary to first find a set of suitable IRFs to be used for
-simulations. This can be easily done using :ref:`csfindobs` and using the ``expression`` parameter. The following example extracts
-observations with at a certain zenith angle range (20-25 degrees) and azimuth angle range (90-100 degrees):
-
-.. code-block:: bash
-	
-	$ csfindobs expression="AZ_PNT>90.0&&AZ_PNT<100.0&&ZEN_PNT>20.0&&ZEN_PNT<25.0"
-	Name of FITS production (Run csiactdata to view your options) [fits-prod-name]
-	Right ascension [83.6331] NONE
-	Runlist outfile [zen20_az90.lis]
-	
-Of course, further details could be specified, e.g. one might want to check for a certain range of muon efficieny, or any other parameter which impacts IRFs.
-The output runlist subsequently has to be converted to an observation XML file using :ref:`csiactobs` (:ref:`see here how to execute this tool <sec_iact_analysis>`).
-In the following it is assumed that the observation XML file of the selected representative observations is produced and called ``myobs.xml``.
-
-Create caldb entry
-^^^^^^^^^^^^^^^^^^
-Now, one of the observations in this container have to be added to the local calibration database. This database is used
-by ctools to extract IRFs for simulations.
-
-.. code-block:: bash
-
-	$ csobs2caldb
-	Input observation definition file [myobs.xml] 
-	Response output name (e.g. Zenith50) [zen20_az90] 
-	
-Note that the response name has been specified to "zen20_az90". Any name can be chosen which represents the observation which was inserted into the database.
-On default, :ref:`csobs2caldb` uses the first observation in the provided XML file. This can be modifyied by specifying the hidden integer parameter index, e.g.:
-
-.. code-block:: bash
-
-  $ csobs2caldb index=1
-  
-This would use the second entry instead of the first from the observation container to extract the IRFs.
-
-.. note::
-
-  The ``index`` parameter starts counting from zero.  
-
-Check caldb entries
-^^^^^^^^^^^^^^^^^^^
-In order to check if the calibration database has successfully been added, the tool :ref:`cscaldb` can be used to inspect the calibration database:
-
-.. code-block:: bash
-
-	$ cscaldb debug=yes
-	...
-	2016-02-25T15:57:11: +==============+
-	2016-02-25T15:57:11: | Mission: cta |
-	2016-02-25T15:57:11: +==============+
-	2016-02-25T15:57:11: === Response functions in database "hess" ===
-	2016-02-25T15:57:11: zen20_az90
-	2016-02-25T15:57:11: 
-	2016-02-25T15:57:11: === Response functions in database "prod2" ===
-	2016-02-25T15:57:11: North_0.5h
-	2016-02-25T15:57:11: North_50h
-	2016-02-25T15:57:11: North_5h
-	2016-02-25T15:57:11: South_0.5h
-	2016-02-25T15:57:11: South_50h
-	2016-02-25T15:57:11: South_5h
-	...
-
-Note that the observation of the actual instrument ("hess" in this example) was added to the CTA mission. This is mainly for structural and simplicity reasons
-and might be improved eventually. Nevertheless, this caldb entry can now be used to simulate any observation duration.
-
-Simulate
-^^^^^^^^
-
-.. code-block:: bash
-  
-	$ ctobssim
-	RA of pointing (degrees) (0-360) [83.63] 
-	Dec of pointing (degrees) (-90-90) [22.01] 
-	Radius of FOV (degrees) (0-180) [2.5] 
-	Start time (MET in s) [0.0] 
-	End time (MET in s) [18000.0] 
-	Lower energy limit (TeV) [0.5] 
-	Upper energy limit (TeV) [50] 
-	Calibration database [hess] 
-	Instrument response function [zen20_az90] 
-	Input model XML file [input_model.xml] 
-	Output event data file or observation definition XML file [simulated_events.fits]
-	
-Here, a total duration of 5 hours at the crab position was simulated using the newly created caldb entry. Note that the choice of the energy range is up to the user.
-One should bear in mind that simulations only make sense in energy ranges where the IRFs are well-defined. The simulation
-assumes one observation of defined duration here. Accordingly, only one FITS event list is written as output. The input XML model
-must now contain only one background model and (of course) all sky models that are needed. Here is an example of an input model
-e.g. for the Crab Nebula:
-
-.. code-block:: xml
-
-	<?xml version="1.0" standalone="no"?>
-	<source_library title="source library">
-	  <source name="Crab" type="PointSource">
-	    <spectrum type="PowerLaw">
-	       <parameter name="Prefactor"   scale="1e-17" value="3.5"  min="1e-07" max="1000.0" free="1"/>
-	       <parameter name="Index"       scale="-1"    value="2.5"  min="0.0"   max="+5.0"   free="1"/>
-	       <parameter name="PivotEnergy" scale="1e6"   value="1.0"  min="0.01"  max="1000.0" free="0"/>
-	    </spectrum>
-	    <spatialModel type="PointSource">
-	      <parameter name="RA"  scale="1.0" value="83.6331" min="-360" max="360" free="0"/>
-	      <parameter name="DEC" scale="1.0" value="22.0145" min="-90"  max="90"  free="0"/>
-	    </spatialModel>
-	  </source>
-	  <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="PivotEnergy" scale="1e6"  value="1.0"  min="0.01" max="1000.0" free="0"/>
-	    </spectrum>
-	  </source>	
-	</source_library>
-
-Note that since the caldb entry (:ref:`csobs2caldb`) has been produced under the mission "CTA", the event list is now considered a CTA observation (regardless of the input instrument).
-Therefore, it is important to have "CTA" specified in the instrument attribute of the background model XML component. Of course the addition of any sky model can be
-done without a problem.
-
-The output of this kind of simulation is a single FITS event list file that can be used for further analysis.
-For instance one could run :ref:`ctlike` in the following way:
-
-.. code-block:: bash
-
-	$ ctlike
-	Input event list, counts cube or observation definition XML file [simulated_events.fits] 
-	Calibration database [hess] 
-	Instrument response function [zen20_az90] 
-	Input model XML file [input_model.xml] 
-	Output model XML file [simulated_fit_results.xml] 
-	
-Note that :ref:`ctlike` automatically queries for the IRF database in case no full observation container is provided. Any other tool taking
-unbinned observations as input will behave the same way. Now all means to run simulations with current IACTs are at hand.
-The following sections will describe tools for running higher level simulation analysis.
-
-Create pull distributions
--------------------------
-Pull distributions inspect the capability of the science tools to find source parameters in the data that were injected into the simulations.
-The tool :ref:`cspull` will run several times the sequence of simulating and fitting the data. After each sequence, the pull value of the
-fitted parameters is written into an ASCII file. The pull value ``p`` of parameter ``X`` is defined as follows:
-
-.. math::
-  p = \frac{X_{\rm fitted} - X_{\rm true}}{\Delta X_{\rm fitted}}
-  
-Where :math:`\Delta X_{\rm fitted}` is the derived error of the fitted parameter :math:`X_{\rm fitted}`. After several trials, the pull value should follow
-a normal distribution. Knowing this, the tool may help to distinguish if a source is detectable and can be normally reconstructed.
-The tool can either work on an input observation container or on the caldb entry that was created above. The number of trials is queried by the tool, too.
-
-Example: observation container
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-.. code-block:: bash
-
-	$ cspull inobs=selected_obs.xml 
-	Number of energy bins (0=unbinned) [0] 
-	Input model XML file [crab_results.xml] 
-	Output pull distribution file [pull.dat] 
-	Number of trials [100] 
-
-Note that the tool asks for the energy binning. It will use an unbinned analysis if zero bins are provided. Otherwise
-it will query for binning parameters and will convert the simulated observations into a binned/stacked observation by running
-:ref:`ctbin`, :ref:`ctexpcube`, :ref:`ctpsfcube` and :ref:`ctbkgcube` in order before fitting the simulated data with :ref:`ctlike`.
-
-Example: caldb entry
-^^^^^^^^^^^^^^^^^^^^
-
-.. code-block:: bash
-
-	$ cspull
-	RA of pointing (deg) (0-360) [83.6331] 
-	Dec of pointing (deg) (-90-90) [22.0145] 
-	Duration (in s) [18000.0] 
-	Lower energy limit (TeV) [0.5]
-	Upper energy limit (TeV) [50]
-	Calibration database [hess] 
-	Instrument response function [zen20_az90] 
-	Number of energy bins (0=unbinned) [0]
-	Input model XML file [input_models.xml] 
-	Output pull distribution file [pull.dat] 
-	Number of trials [100] 
-
-Again, the caldb entry ``zen20_az90`` which was created above is used to simulate 5 hours of data. The data is simulated and analysed 100
-times in a row and the output is written into the file ``pull.dat``.
-
-Visualise pull distributions
-----------------------------
-The output ASCII file of :ref:`cspull` can be inspected by the user using own scripts. As a starting point, the ``example`` folder
-contains a script that can be used to visualise the pull distribution histogram:
-
-.. code-block:: bash
-
-	$ python $CTOOLS/examples/show_pull_histogram.py pull.dat Pull_Crab_Prefactor 50
-
-Run the script without arguments to see the usage. In case a wrong parameter name was provided, the tool will print available parameters on the screen.
-The last integer parameter defines the binning of the histogram. In addition, the tool shows the standard normal distribution which should
-match the pull distribution for a high number of trials if everything went well.
-
-Create sensitivity curve
-------------------------
-A very important means to illustrate if a source can be detected with a certain significance is a sensitivity curve.
-The tool :ref:`cssens` can be used to compute either an integral, or a differential sensitivity curve for a given observation container.
-Analogous to :ref:`cspull`, the tool works either on an unbinned observation container or using a caldb entry. 
-The script :ref:`cssens` computes the differential or integrated CTA sensitivity using maximum likelihood fitting of a test source.
-The differential sensitivity is determined for a number of energy bins, the integral sensitivity is determined for a number
-of energy thresholds. The test source is fitted to simulated data using ctlike to determine it’s detection significance as a
-function of source flux. The source flux is then varied until the source significance achieves a given level, specified by the 
-(hidden) significance parameter ``sigma``. To damp variations between individual Monte Carlo simulations, a sliding average is
-applied in the significance computation (controlled by hidden ``num_avg parameter``). Note that this procedure might take a while to compute.
-To switch from diffferential to integral sensitivity, the hidden parameter ``type="Integral"`` has to be specified. It should be noted that
-the results are written to a file ``sensitivity.dat``. To change the output file name, the hidden parameter ``outfile=mysensitivity.dat``
-can be provided.
-
-
-Example: observation container
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-.. code-block:: bash
-
-	$ cssens inobs=selected_obs.xml
-	Input model XML file [crab_results.xml] 
-	Source name [Crab] 
-	Lower energy limit (TeV) [0.5] 
-	Upper energy limit (TeV) [50] 
-	Number of energy bins for differential sensitivity computation [21] 
-
-Example: caldb entry
-^^^^^^^^^^^^^^^^^^^^
-
-.. code-block:: bash
-
-	$ cssens 
-	Calibration database [hess] 
-	Instrument response function [zen20_az90] 
-	Effective exposure time (s) [18000.0] 
-	Radius of ROI (deg) [2.5] 
-	Input model XML file [input_models.xml] 
-	Source name [Crab] 
-	Lower energy limit (TeV) [0.5] 
-	Upper energy limit (TeV) [50] 
-	Number of energy bins for differential sensitivity computation [21]
-	
-Note that the input model in the latter case only requires to contain one background model (see section above). Using this setup,
-the user can specify the observation time for the sensitivity curve. In this example, again, 5 hours are used.
-
-Visualise sensitivity curves
-----------------------------
-In order to show a sensitivity curve, there is also a script in the example folder:
-
-.. code-block:: bash
-
-	$ python $CTOOLS/examples/show_sensitivity.py sensitivity.dat
-	
-This will display a sensitivity curve for the applied settings.
-.. _sec_iact_tips:
-
-Tips, Tricks and FAQs
-=====================
-This section shows some of the wrinkles that could ease your life when
-running IACT analyses.
-
-
-Visualise observations
-----------------------
-A very common task is to inspect the observation container and visualise several
-quantities.
-
-Create a skymap of counts
-^^^^^^^^^^^^^^^^^^^^^^^^^
-The tool :ref:`ctskymap` allows to create a single map of counts contained
-inside an observation container. In contrast to :ref:`ctbin`, ctskymap takes
-into account all events that fall into the skymap energy range and spatial
-range. Since :ref:`ctbin` is intended for stacked analysis, it respects energy
-thresholds, RoI selections, etc. In contrast, :ref:`ctskymap` may be used for
-a simple visualisation of the counts.
-
-.. code-block:: bash
-
-   $ ctskymap debug=yes
-   Input event list or observation definition XML file [events.fits] obs.xml
-   Lower energy limit (TeV) [0.1]
-   Upper energy limit (TeV) [100.0]
-   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]
-   Size of the Y axis in pixels [200]
-   Output skymap file [skymap.fits]
-	
-Alternatively, you may start from a counts cube and collapse the maps into one
-single map via Python:
-
-.. code-block:: python
- 
-   >>> import gammalib
-   >>> map = gammalib.GSkyMap("myCountCube.fits") # Load counts cube
-   >>> map.stack_maps()                           # Stack all cube maps into a single map
-   >>> map.save("myCountMap.fits")                # Save map to file
-	
-Create region file of pointings
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-You may want to know how the pointings of the observation container distribute
-across the sky. For this purpose, :ref:`csobsinfo` provides the hidden parameter
-``ds9file`` which instructs the script produce a region file that contains the
-pointing information:
-
-.. code-block:: bash
-
-   $ csobsinfo ds9file=pointings.reg debug=yes
-   Event list, counts cube, or observation definition file [obs.xml]
-	
-Together with the previously produced skymap, you can plot the pointings using
-for example ds9:
-
-.. code-block:: bash
-  
-   $ ds9 skymap.fits -regions pointings.reg
-
-
-Plot zenith angle distribution of observations
-----------------------------------------------
-The script :ref:`csobsinfo` logs information about an observation container
-into a log file. You may access the information through dedicated class methods
-and plot the information using for example matplotlib:
-
-.. code-block:: python
-
-   >>> # Run script
-   >>> import cscripts
-   >>> info = cscripts.csobsinfo()
-   >>> info["inobs"] = "selected_obs.xml"
-   >>> info.run()
-
-   >>> # Plot data
-   >>> import matplotlib.pyplot as plt
-   >>> plt.hist(info.zeniths(), bins=10, range=[0,90])
-   >>> plt.show()
-  
-This example will plot a histogram showing the zenith angle distribution of the
-observation container. In addition there are more methods to access information
-in Python:
-
-.. code-block:: python
-
-   >>> info.azimuths() # list of azimuth values
-   >>> info.offsets()  # list of offset values (only computed if hidden parameter offset=yes was specified)
-   >>> info.ebounds()  # gammalib.GEbounds object of energy ranges
-   >>> info.gti()      # gammalib.GGti object containing good time intervals
-  
-There is also an example script to plot information about an observation
-container which might be useful for a quick visualisation:
-
-.. code-block:: bash
-  
-   $ python $CTOOLS/examples/show_obs.py selected_obs.xml
-
-
-Visualise models
-----------------
-Similar to :ref:`csobsinfo`, there is the :ref:`csmodelinfo` script to find out
-more about the content of a model XML file. In order to show the position and
-sizes of the model on top of a skymap, this tool has the hidden parameter
-``ds9file``, too.
-
-.. code-block:: bash
-
-   $ csmodelinfo ds9file=models.reg debug=yes
-   Input model XML file [$CTOOLS/share/models/crab.xml]
-
-There are several options regarding the color, text and other attributes of the
-region file. To see a full list, :ref:`visit the reference page <csmodelinfo>`.
-Analogous to the pointings of the observation container, the models can be
-visualised using for example ds9:
-
-.. code-block:: bash
-  
-   $ ds9 skymap.fits -regions models.reg
-
-
-Manipulating models in python
------------------------------
-This example will give you some guidance on how to work with model XML files in
-Python and how to manipulate their content:
-
-.. code-block:: python
-
-   >>> import gammalib
-
-   >>> # Open model file
-   >>> models = gammalib.GModels("$CTOOLS/share/models/crab.xml")
-
-   >>> # Access a model component
-   >>> src = models["Crab"]
-   >>> print(src)
-
-   >>> # Retrieve the spectral or spatial component
-   >>> print(src.spectral())
-   >>> print(src.spatial())
-
-   >>> # Print the spectral parameter "Prefactor"
-   >>> prefactor = src.spectral()["Prefactor"]
-   >>> print("Prefactor: "+str(prefactor.value())+" +- "+str(prefactor.error()))
-
-   >>> # Set the prefactor value
-   >>> prefactor.value(3.5e-17)
-
-   >>> # Loop over models and fix all parameters
-   >>> for model in models:
-   >>>     for par in model:   # Loop over model parameters
-   >>>         par.fix()       # Fix parameter
-
-   >>> # Release one specific parameter
-   >>> models["Crab"]["Prefactor"].free()
-
-   >>> # Set the parameter range
-   >>> models["Crab"]["Prefactor"].min(1e-18)
-   >>> models["Crab"]["Prefactor"].max(1e-16)
-
-   >>> # ... or in one step:
-   >>> models["Crab"]["Prefactor"].range(1e-18, 1e-16)
-
-   >>> # Remove model from container
-   >>> models.remove("Crab")
-
-   >>> # Save model to another XML file
-   >>> models.save("mymodels.xml")
-  
-
-Retrieve likelihood values from ctlike
---------------------------------------
-For some purposes, it might be useful to retrieve the fitted likelihood value
-and other results of the fit with :ref:`ctlike`.
-
-.. code-block:: python
-
-   >>> import gammalib
-   >>> import ctools
-
-   >>> # Create and run ctlike
-   >>> like = ctools.ctlike()
-   >>> like["inobs"]   = "selected_obs.xml"
-   >>> like["inmodel"] = "mymodels.xml"
-   >>> like.run()
-
-   >>> # Get optimizer from ctlike
-   >>> opt = like.opt()
-   >>> print(opt)
-
-   >>> # Get likelihood value from optimizer
-   >>> print(opt.value())
-
-   >>> # Get fit status:
-   >>> print(opt.status())
-   >>> # 0: converged
-   >>> # 1: stalled
-   >>> # 2: singular curvature matrix encountered
-   >>> # 3: curvature matrix not positive definite
-   >>> # 4: errors are inaccurate
-
-
-Speed up analysis
------------------
-In some cases it may occur that the analysis takes very long. There are several
-reasons that can slow down the analysis quite severly. Not all issues can be
-tackled. Nevertheless, here is a list of actions that have proven to speed up
-the fit.
-
-* Consider switching to stacked analysis if the observation time exceeds 50-100
-  hours.
-* Diffuse and extended sources in the RoI are computational expensive. However,
-  in unbinned analysis they impact the speed much more.
-* Usage of energy dispersion could also cause a long fit. Try switching it off
-  if speed is important. Keep in mind this will have an effect on the results.
-* Switch to a machine that supports OpenMP. This allows to parallelise the fit
-  onto several cores.
-
-
-Compute excess maps
--------------------
-A very important means to visualise images of sources are excess maps. The script
-:ref:`csresmap` can also be used to create such a map. The user, however, has to
-manually remove the source(s) of interest from the input model XML file.
-
-.. code-block:: python
-
-   >>> import gammalib
-   >>> models = gammalib.GModels("crab_models.xml")
-   >>> models.remove("Crab")
-   >>> models.save("models_without_crab.xml")
-  
-Of course, this can also be done by editing the XML file with an editor.
-
-Subsequently, :ref:`csresmap` can be executed using the ``algorithm=SUB``
-parameter:
-
-.. code-block:: bash
-
-   $ csresmap
-   Input event list, counts cube, or observation definition XML file [selected_obs.xml]
-   Input model XML file [models_without_crab.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]
-   Coordinate System (CEL|GAL) [CEL]
-   Projection method (AIT|AZP|CAR|MER|MOL|STG|TAN) [CAR]
-   Size of the X axis in pixels [100]
-   Size of the Y axis in pixels [100]
-   Pixel size (deg/pixel) [0.02]
-   Residual map computation algorithm (SUB|SUBDIV|SUBDIVSQRT) [SUB]
-   Output residual map file [excessmap.fits]
-
-
-Compute flux maps
------------------
-The creation of a flux map is much more expensive in terms of computing than an
-excess map. In the section about :ref:`high level analysis tools <sec_iact_high_level>`,
-the usage of :ref:`cttsmap` was demonstrated. This tool can also be used to
-generate flux maps. For this, the user has to incorporate a test source in the
-XML model file. This source should consist of a spectral component using a
-``PowerLawPhotonFlux`` model:
-
-.. code-block:: xml
-
-	<source name="TestSource" type="PointSource">
-	 <spectrum type="PowerLawPhotonFlux">
-	  <parameter scale="1e-07" name="PhotonFlux" min="1e-07" max="1000.0"    value="1.0" free="1"/>
-	  <parameter scale="1.0"   name="Index"      min="-5.0"  max="+5.0"      value="-2.0" free="1"/>
-	  <parameter scale="1.0"   name="LowerLimit" min="10.0"  max="1000000.0" value="100.0" free="0"/>
-	  <parameter scale="1.0"   name="UpperLimit" min="10.0"  max="1000000.0" value="500000.0" free="0"/>
-	 </spectrum>
-	 <spatialModel type="PointSource">
-	  <parameter free="0" max="360" min="-360" name="RA"  scale="1" value="83.6331" />
-	  <parameter free="0" max="90"  min="-90"  name="DEC" scale="1" value="22.0145" />
-	 </spatialModel>
-	</source>
-
-The further content of the rest of the XML model file depends on the user
-requirements:
-
-* For an absolute flux map, all other sky models should be removed
-* For a residual flux map, all other sky models should be kept
-
-The tool :ref:`cttsmap` will create one skymap per free spectral parameter in
-the model. In the result FITS file, there will be an extension called ``Integral``
-that contains the flux map. For the above XML example, the source name
-``TestSource`` should be specified to :ref:`cttsmap`.
-
-
-Speed up TS map computation
----------------------------
-Since the computation of a TS map can be extremely time consuming, the option of
-splitting the computation into several jobs is supported. This might be of
-particular interest if the user has access to a batch farm.
-
-Split TS map computation
-^^^^^^^^^^^^^^^^^^^^^^^^
-For the purpose of job splitting, the hidden parameter ``binmin`` and ``binmax``
-were included in the tool. These are integer parameter that specify which bins
-should be computed. For instance, if the map should consist of 30x30(=900) pixels,
-the user could for example execute
-
-.. code-block:: bash
-
-	$ cttsmap binmin=0 binmax=299 outmap=tsmap_0_299.fits
-	...
-	$ cttsmap binmin=300 binmax=599 outmap=tsmap_300_599.fits
-	...
-	$ cttsmap binmin=600 binmax=899 outmap=tsmap_600_899.fits
-	...
-	
-Each command could run on a different machine.
-
-.. note::
-
-  The output file name of individual jobs should be different. Otherwise files
-  could overwrite each other. The naming of the individual slices is up to the
-  user.
-  
-The script :ref:`cstsmapsplit` will take care of this bookkeeping. It creates an
-ASCII file containing all the commands according to the user input. The following
-example will create split the computation of a TS map into 2000 separate task.
-Each task will only compute 5 bins on its own. This is very useful if the
-observation container is large or the fit simply takes a long time.
-
-.. code-block:: bash
-
-  $ cstsmapsplit
-  Input event list, counts cube or observation definition XML file [selected_obs.xml]
-  Input model definition XML file [$CTOOLS/share/models/crab.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 [100]
-  Size of the Y axis in pixels [100]
-  Test source name [Crab] 
-  Output Test Statistic map file [tsmap.fits] 
-  Number of TS map bins per task [5] 
-  Compute null hypothesis first? [yes] 
-  ASCII file containing all commands [commands.dat] 
-  
-.. note::
-  
-  One can decide if we want to compute the null hypothesis first. This way, a
-  parameter optimisation will be performed without the test source. The obtained
-  likelihood value will then be passed to each individual task of :ref:`cttsmap`.
-  If ``compute_null=no``, each task has to compute the null hypothesis itself.
-  
-Merge splitted TS maps
-^^^^^^^^^^^^^^^^^^^^^^
-The script :ref:`cstsmapmerge` is intended to take care of merging the files
-that were produced while splitting the TS map computation. There are several
-options to pass the files to be merged as arguments:
-
-* a space-separated list of file names (e.g. tsmap1.fits tsmap2.fits)
-* a comma-separated list of file names (e.g. tsmap1.fits,tsmap2.fits)
-* a wildcard string (e.g. tsmap*.xml)
-* an ASCII file containing the file names, one file per line (e.g. @mymaps.txt)
-
-In this example the ASCII file method is presented:
-
-.. code-block:: bash 
-
-   # Put slice files into an ascii file
-   $ ls tsmap_*.fits > tsmapfiles.txt
-
-   # Run cstsmapmerge
-   $ cstsmapmerge
-   Input TS map FITS files [@tsmapfiles.txt]
-   Output TS map FITS file [mytsmap.fits]
-
-
-Creating a python analysis pipeline
------------------------------------
-It is easily possible to build an own analysis workflow with a simple python
-script. The following source code example shows a python script running from
-gathering observations until fitting spectral points without storing intermediate
-data products on disk. It assumes that the environment variable ``$VHEFITS`` is
-set to the path where IACT FITS data is located.
-
-.. code-block:: python
-  
-   >>> import gammalib
-   >>> import ctools
-   >>> import cscripts
-
-   >>> # Set debug flag
-   >>> debug = True
-
-   >>> # Set flag to use energy dispersion
-   >>> edisp = False
-
-   >>> # Set inmodel file name
-   >>> inmodel = "$GAMMALIB/test/data/model_point_plaw.xml"
-
-   >>> # Expand environment variable
-   >>> inmodel = gammalib.expand_env(inmodel)
-
-   >>> # Extract coordinates and model properties
-   >>> models  = gammalib.GModels(inmodel)
-   >>> srcname = models[0].name()
-   >>> ra      = models[0]["RA"].value()
-   >>> dec     = models[0]["DEC"].value()
-
-   >>> # Find FITS production name
-   >>> iactdata          = cscripts.csiactdata()
-   >>> iactdata["debug"] = debug
-   >>> iactdata.run()
-
-   >>> # Use first available production
-   >>> prodname = iactdata.names()[0]
-
-   >>> # Run csfindobs
-   >>> findobs             = cscripts.csfindobs()
-   >>> findobs["ra"]       = ra
-   >>> findobs["dec"]      = dec
-   >>> findobs["rad"]      = 2.5
-   >>> findobs["prodname"] = prodname
-   >>> findobs["debug"]    = debug
-   >>> findobs["outfile"]  = "NONE"
-   >>> findobs.run()
-
-   >>> # Retrieve obervation IDs (runlist)
-   >>> obs_ids = findobs.obs_ids()
-
-   >>> # Build observation container
-   >>> iactobs             = cscripts.csiactobs()
-   >>> iactobs["prodname"] = prodname
-   >>> iactobs["inmodel"]  = inmodel
-   >>> iactobs["bkgpars"]  = 1
-   >>> iactobs["outobs"]   = "NONE"
-   >>> iactobs["outmodel"] = "NONE"
-   >>> iactobs["debug"]    = debug
-   >>> iactobs.runlist(obs_ids)
-   >>> iactobs.run()
-
-   >>> # Retrieve observation container and energy boundaries
-   >>> obs     = iactobs.obs()
-   >>> ebounds = iactobs.ebounds()
-
-   >>> # Run ctselect
-   >>> select             = ctools.ctselect(obs)
-   >>> select["usepnt"]   = True
-   >>> select["rad"]      = 2.5
-   >>> select["usethres"] = "DEFAULT"
-   >>> select["tmin"]     = 0.0
-   >>> select["tmax"]     = 0.0
-   >>> select["emin"]     = 0.1
-   >>> select["emax"]     = 100.0
-   >>> select["debug"]    = debug
-   >>> select.run()
-
-   >>> # Pass selected observations to ctlike
-   >>> like          = ctools.ctlike(select.obs())
-   >>> like["debug"] = debug
-   >>> like["edisp"] = edisp
-   >>> like.run()
-
-   >>> # Compute a spectrum and save
-   >>> spec             = cscripts.csspec(like.obs())
-   >>> spec["srcname"]  = srcname
-   >>> spec["emin"]     = ebounds.emin().TeV()
-   >>> spec["emax"]     = ebounds.emax().TeV()
-   >>> spec["enumbins"] = 10
-   >>> spec["edisp"]    = edisp
-   >>> spec["ebinalg"]  = "LOG"
-   >>> spec["debug"]    = debug
-   >>> spec["outfile"]  = "spectrum.fits"
-   >>> spec.execute()
-
-   >>> # Remove model for excess map computation
-   >>> like.obs().models().remove(srcname)
-
-   >>> # Compute an excess map
-   >>> resmap             = cscripts.csresmap(like.obs())
-   >>> resmap["xref"]     = ra
-   >>> resmap["yref"]     = dec
-   >>> resmap["proj"]     = "CAR"
-   >>> resmap["coordsys"] = "CEL"
-   >>> resmap["emin"]     = ebounds.emin().TeV()
-   >>> resmap["emax"]     = ebounds.emax().TeV()
-   >>> resmap["nxpix"]    = 100
-   >>> resmap["nypix"]    = 100
-   >>> resmap["binsz"]    = 0.02
-   >>> resmap["outmap"]   = "excessmap.fits"
-   >>> resmap.execute()
-
-The results of this analysis workflow can be inspected using an example Python
-script and for example ds9:
-
-.. code-block:: bash
-  
-   $ python $CTOOLS/examples/show_spectrum.py spectrum.fits
-   $ ds9 excessmap.fits
-  
@@ -46,5 +46,5 @@ Analysis of data from other instruments
    How to analyse Fermi-LAT data with ctools? <fermi/index>
    How to analyse COMPTEL data with ctools? <comptel/index>
-   How to analyse H.E.S.S. data with ctools? <iact/index>
+   How to analyse data from existing IACTs (H.E.S.S., VERITAS or MAGIC) with ctools? <iact/index>
    howto_combine_instruments