diff --git a/doc/source/contents.rst b/doc/source/contents.rst index 10c9ae9..e9e1af5 100644 --- a/doc/source/contents.rst +++ b/doc/source/contents.rst @@ -10,9 +10,9 @@ Documentation :maxdepth: 2 about + download user_manual/index reference_manual/index faq issues - download develop \ No newline at end of file diff --git a/doc/source/user_manual/getting_started/python.rst b/doc/source/user_manual/getting_started/python.rst index cf6bead..49b1638 100644 --- a/doc/source/user_manual/getting_started/python.rst +++ b/doc/source/user_manual/getting_started/python.rst @@ -28,20 +28,20 @@ way as if it were executed from the command line. .. code-block:: python - import ctools - sim = ctools.ctobssim() - sim["inmodel"] = "${CTOOLS}/share/models/crab.xml" - sim["outevents"] = "events.fits" - sim["caldb"] = "prod2" - sim["irf"] = "South_50h" - sim["ra"] = 83.63 - sim["dec"] = 22.01 - sim["rad"] = 5.0 - sim["tmin"] = 0 - sim["tmax"] = 1800 - sim["emin"] = 0.1 - sim["emax"] = 100 - sim.execute() + >>> import ctools + >>> sim = ctools.ctobssim() + >>> sim["inmodel"] = "${CTOOLS}/share/models/crab.xml" + >>> sim["outevents"] = "events.fits" + >>> sim["caldb"] = "prod2" + >>> sim["irf"] = "South_0.5h" + >>> sim["ra"] = 83.63 + >>> sim["dec"] = 22.01 + >>> sim["rad"] = 5.0 + >>> sim["tmin"] = 0.0 + >>> sim["tmax"] = 1800.0 + >>> sim["emin"] = 0.1 + >>> sim["emax"] = 100.0 + >>> sim.execute() Alternatively, you may "run" the :ref:`ctobssim` tool using @@ -76,7 +76,14 @@ To visualise this container, type: Livetime ..................: 1710 s Deadtime correction .......: 0.95 User energy range .........: undefined - ... + === GCTAPointing === + Pointing direction ........: (RA,Dec)=(83.63,22.01) + === GCTAResponseIrf === + Caldb mission .............: cta + Caldb instrument ..........: prod2 + Response name .............: South_0.5h + Energy dispersion .........: Not used + Save energy range .........: undefined There is one CTA observation in the container and to visualise the events in that observation you may type: @@ -85,7 +92,7 @@ in that observation you may type: >>> print(sim.obs()[0].events()) === GCTAEventList === - Number of events ..........: 31928 + Number of events ..........: 45370 Time interval .............: 51544.5 - 51544.5 days === GEbounds === Number of intervals .......: 1 @@ -118,9 +125,9 @@ terrestrial times (TT) of all events. >>> for event in events: ... print(event) ... - Dir=RA=83.5759, DEC=21.9757 [-0.000597839,-0.000875483] Energy=328.153 GeV Time=-3.15576e+08 s (TT) - Dir=RA=83.5775, DEC=22.0128 [4.88943e-05,-0.000849504] Energy=353.264 GeV Time=-3.15576e+08 s (TT) - Dir=RA=83.4302, DEC=21.931 [-0.00137616,-0.0032343] Energy=117.999 GeV Time=-3.15576e+08 s (TT) + Dir=RA=83.5743, DEC=21.9746 [-0.000617357,-0.000902194] Energy=328.153 GeV Time=-3.15576e+08 s (TT) + Dir=RA=83.5756, DEC=22.0127 [4.78931e-05,-0.000879886] Energy=353.264 GeV Time=-3.15576e+08 s (TT) + Dir=RA=83.6811, DEC=22.0034 [-0.000114612,0.000826825] Energy=550.35 GeV Time=-3.15576e+08 s (TT) ... .. note:: @@ -141,7 +148,7 @@ terrestrial times (TT) of all events. >>> sim.run() >>> print(sim.obs()[0].events()) === GCTAEventList === - Number of events ..........: 15964 + Number of events ..........: 22685 Time interval .............: 51544.5 - 51544.5 days === GEbounds === Number of intervals .......: 1 @@ -180,7 +187,7 @@ To check how the fit went you can inspect the optimiser used by >>> print(like.opt()) === GOptimizerLM === - Optimized function value ..: 101048.382 + Optimized function value ..: 152299.205 Absolute precision ........: 0.005 Acceptable value decrease .: 2 Optimization status .......: converged @@ -213,8 +220,8 @@ member of the observation container: RA .......................: 83.6331 [-360,360] deg (fixed,scale=1) DEC ......................: 22.0145 [-90,90] deg (fixed,scale=1) Number of spectral par's ..: 3 - Prefactor ................: 5.80565e-16 +/- 1.11147e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient) - Index ....................: -2.49266 +/- 0.0163476 [-0,-5] (free,scale=-1,gradient) + Prefactor ................: 5.81591e-16 +/- 1.01725e-17 [1e-23,1e-13] ph/cm2/s/MeV (free,scale=1e-16,gradient) + Index ....................: -2.50618 +/- 0.0153538 [-0,-5] (free,scale=-1,gradient) PivotEnergy ..............: 300000 [10000,1e+09] MeV (fixed,scale=1e+06,gradient) Number of temporal par's ..: 1 Normalization ............: 1 (relative value) (fixed,scale=1,gradient) @@ -226,8 +233,8 @@ member of the observation container: Model type ................: "PowerLaw" * "Constant" Number of parameters ......: 4 Number of spectral par's ..: 3 - Prefactor ................: 1.01701 +/- 0.0216563 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient) - Index ....................: 0.0230864 +/- 0.0123006 [-5,5] (free,scale=1,gradient) + Prefactor ................: 0.998788 +/- 0.0118567 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient) + Index ....................: 0.00990034 +/- 0.00735768 [-5,5] (free,scale=1,gradient) PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient) Number of temporal par's ..: 1 Normalization ............: 1 (relative value) (fixed,scale=1,gradient) @@ -243,8 +250,8 @@ This is easy from Python: >>> like.run() >>> print(like.obs().models()) ... - RA .......................: 83.6341 +/- 0.000958947 [-360,360] deg (free,scale=1) - DEC ......................: 22.0139 +/- 0.000870912 [-90,90] deg (free,scale=1) + RA .......................: 83.6334 +/- 0.000917669 [-360,360] deg (free,scale=1) + DEC ......................: 22.0144 +/- 0.000836554 [-90,90] deg (free,scale=1) The ``like.obs().models()`` method provides the model container, using the ``["Crab"]`` operator we access the Crab model in that container and using @@ -286,7 +293,7 @@ the ``fix()`` method. >>> print(like.obs()) === GObservations === Number of observations ....: 1 - Number of predicted events : 15964 + Number of predicted events : 22685 is okay as the ``ctools.ctlike(sim.obs())`` constructor will create a copy of the observation container that lives within the :ref:`ctlike` @@ -312,15 +319,22 @@ ctools provides the Python module ``obsutils`` that may further simplify your analysis efforts. ``obsutils`` is a Python script that makes use of the GammaLib and ctools modules to create standard analysis steps. +As all Python scripts, ``obsutils`` is part of the cscripts module that +is imported using + +.. code-block:: python + + >>> import cscripts + Here an example of how to use ``obsutils``: .. code-block:: python >>> import gammalib >>> import ctools - >>> from ctools import obsutils + >>> from cscripts import obsutils >>> pattern = obsutils.set_obs_patterns("four", ra=83.63, dec=22.01, offset=1.0) - >>> obs = obsutils.set_obs_list(pattern, duration=1800, emin=0.1, emax=100.0, rad=5.0, caldb="prod2", irf="South_50h") + >>> obs = obsutils.set_obs_list(pattern, duration=1800, emin=0.1, emax=100.0, rad=5.0, caldb="prod2", irf="South_0.5h") >>> print(obs) === GObservations === Number of observations ....: 4 @@ -335,7 +349,7 @@ Here an example of how to use ``obsutils``: Number of parameters ......: 10 ... -The module is imported using the ``from ctools import obsutils`` directive. +The module is imported using the ``from cscripts import obsutils`` directive. The ``obsutils.set_obs_patterns()`` function will create a pointing pattern of four observations located at offset angles of 1 degree from the nominal location of the Crab nebula. @@ -354,10 +368,13 @@ The ``obsutils.sim()`` function then simulates the event data, the functions will be added in the future. -Access common analysis results -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -Here are some examples how to access your analysis results in python. In the following it is -assumed that you have a ctlike object called 'like' which is setup and runs: +Access analysis results +~~~~~~~~~~~~~~~~~~~~~~~ + +Here are some examples that show how to access your analysis results in +python. +In the following it is assumed that you have a ctlike object called +``like`` which is setup and runs: .. code-block:: python @@ -372,9 +389,9 @@ assumed that you have a ctlike object called 'like' which is setup and runs: .. code-block:: python - obs = like.obs() # Get observations object - obs.models()["Crab"]["Prefactor"].value() # This returns the actual fit value - obs.models()["Crab"]["Prefactor"].error() # This returns the actual fit error + obs = like.obs() # Get observations object + obs.models()["Crab"]["Prefactor"].value() # This returns the actual fit value + obs.models()["Crab"]["Prefactor"].error() # This returns the actual fit error * Open XML model file from python and print models on screen: @@ -385,16 +402,25 @@ assumed that you have a ctlike object called 'like' which is setup and runs: * Likelihood value - .. code-block:: python + .. code-block:: python like.opt().value() # Returns likelihood value -* Covariance Matrix +* Curvature Matrix (aka Hessian) - .. code-block:: python + .. code-block:: python - covariance = like.opt().curvature().invert() # Return GMatrix object - print(covariance) + curvature = like.opt().curvature().invert() # Return GMatrix object + print(curvature) + + To get the covariance matrix, the curvature matrix needs to be + inverted: + + .. code-block:: python + + covariance = curvature.invert() + print(covariance) + Modify and work on XML models ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -415,7 +441,8 @@ Modify and work on XML models .. code-block:: python - models["Crab"]["Prefactor"].is_free() + if models["Crab"]["Prefactor"].is_free(): + ... * Set parameter value @@ -430,6 +457,12 @@ Modify and work on XML models models["Crab"]["Index"].min(-5.0) models["Crab"]["Index"].max(-1.0) + or quicker + + .. code-block:: python + + models["Crab"]["Index"].range(-5.0,-1.0) + * Loop over models and parameters (e.g. fix all parameters) .. code-block:: python