Iftach Sadeh / ctools

Commit 1c0d6ef39cee1ac38ae88d4df11a7e3585427430

Authored by Jürgen Knödlseder
1 parent 0e2f0c05

Make example pipeline script compliant with Python coding conventions

Inline Side-by-side
Showing 6 changed files with 586 additions and 425 deletions
examples/pipeline_binned_disk.py
  View file @1c0d6ef
1 #! /usr/bin/env python 1 #! /usr/bin/env python
2 # ========================================================================== 2 # ==========================================================================
3 -# Perform binned CTA analysis based of simulated CTA data. 3 +# Perform binned analysis of simulated CTA data.
4 # 4 #
5 # Copyright (C) 2015-2016 Juergen Knoedlseder 5 # Copyright (C) 2015-2016 Juergen Knoedlseder
6 # 6 #
@@ -30,26 +30,26 @@ def setup_observations(pattern='four', ra=83.63, dec=22.01, offset=1.5, @@ -30,26 +30,26 @@ def setup_observations(pattern='four', ra=83.63, dec=22.01, offset=1.5,
30 emin=0.1, emax=100.0, rad=5.0, duration=180.0, 30 emin=0.1, emax=100.0, rad=5.0, duration=180.0,
31 deadc=0.95, caldb='prod2', irf='South_0.5h'): 31 deadc=0.95, caldb='prod2', irf='South_0.5h'):
32 """ 32 """
33 - Returns an observation container. 33 + Returns an observation container
34 34
35 Parameters 35 Parameters
36 ---------- 36 ----------
37 pattern : str, optional 37 pattern : str, optional
38 Pointing pattern, either 'single' or 'four' 38 Pointing pattern, either 'single' or 'four'
39 ra : float, optional 39 ra : float, optional
40 - Right Ascension of pattern centre [deg] 40 + Right Ascension of pattern centre (deg)
41 dec : float, optional 41 dec : float, optional
42 - Delication of pattern centre [deg] 42 + Declination of pattern centre (deg)
43 offset : float, optional 43 offset : float, optional
44 - Offset between observations of pattern [deg] 44 + Offset between observations of pattern (deg)
45 emin : float, optional 45 emin : float, optional
46 - Minimum energy [TeV] 46 + Minimum energy (TeV)
47 emax : float, optional 47 emax : float, optional
48 - Maximum energy [TeV] 48 + Maximum energy (TeV)
49 rad : float, optional 49 rad : float, optional
50 - ROI radius used for analysis [deg] 50 + ROI radius used for analysis (deg)
51 duration : float, optional 51 duration : float, optional
52 - Duration of one CTA observation [seconds] 52 + Duration of one CTA observation (s)
53 deadc : float, optional 53 deadc : float, optional
54 Deadtime correction factor 54 Deadtime correction factor
55 caldb : str, optional 55 caldb : str, optional
@@ -86,7 +86,7 @@ def setup_observations(pattern='four', ra=83.63, dec=22.01, offset=1.5, @@ -86,7 +86,7 @@ def setup_observations(pattern='four', ra=83.63, dec=22.01, offset=1.5,
86 # =========== # 86 # =========== #
87 def setup_model(obs, model='data/crab.xml'): 87 def setup_model(obs, model='data/crab.xml'):
88 """ 88 """
89 - Setup model for analysis. 89 + Setup model for analysis
90 90
91 Parameters 91 Parameters
92 ---------- 92 ----------
@@ -117,16 +117,16 @@ def run_pipeline(obs, emin=0.1, emax=100.0, @@ -117,16 +117,16 @@ def run_pipeline(obs, emin=0.1, emax=100.0,
117 caldb='prod2', irf='South_0.5h', 117 caldb='prod2', irf='South_0.5h',
118 debug=False): 118 debug=False):
119 """ 119 """
120 - Simulation and binned analysis pipeline. 120 + Simulation and binned analysis pipeline
121 121
122 Parameters 122 Parameters
123 ---------- 123 ----------
124 obs : `~gammalib.GObservations` 124 obs : `~gammalib.GObservations`
125 Observation container 125 Observation container
126 emin : float, optional 126 emin : float, optional
127 - Minimum energy [TeV] 127 + Minimum energy (TeV)
128 emax : float, optional 128 emax : float, optional
129 - Maximum energy [TeV] 129 + Maximum energy (TeV)
130 enumbins : int, optional 130 enumbins : int, optional
131 Number of energy bins 131 Number of energy bins
132 nxpix : int, optional 132 nxpix : int, optional
@@ -134,7 +134,7 @@ def run_pipeline(obs, emin=0.1, emax=100.0, @@ -134,7 +134,7 @@ def run_pipeline(obs, emin=0.1, emax=100.0,
134 nypix : int, optional 134 nypix : int, optional
135 Number of pixels in Y axis 135 Number of pixels in Y axis
136 binsz : float, optional 136 binsz : float, optional
137 - Pixel size [deg/pixel] 137 + Pixel size (deg)
138 coordsys : str, optional 138 coordsys : str, optional
139 Coordinate system 139 Coordinate system
140 proj : str, optional 140 proj : str, optional
examples/pipeline_binned_mem.py
  View file @1c0d6ef
1 #! /usr/bin/env python 1 #! /usr/bin/env python
2 # ========================================================================== 2 # ==========================================================================
3 -# This script illustrates how to perform a binned CTA analysis based on  
4 -# simulated CTA data. You may use and adapt this script to implement your  
5 -# own pipeline. 3 +# Perform binned in-memory analysis of simulated CTA data.
6 # 4 #
7 -# Usage:  
8 -# ./pipeline_binned_mem.py  
9 -#  
10 -# ==========================================================================  
11 -#  
12 -# Copyright (C) 2015 Juergen Knoedlseder 5 +# Copyright (C) 2015-2016 Juergen Knoedlseder
13 # 6 #
14 # This program is free software: you can redistribute it and/or modify 7 # This program is free software: you can redistribute it and/or modify
15 # it under the terms of the GNU General Public License as published by 8 # it under the terms of the GNU General Public License as published by
@@ -33,25 +26,41 @@ from cscripts import obsutils @@ -33,25 +26,41 @@ from cscripts import obsutils
33 # ================== # 26 # ================== #
34 # Setup observations # 27 # Setup observations #
35 # ================== # 28 # ================== #
36 -def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, 29 +def setup_observations(pattern='four', ra=83.63, dec=22.01, offset=1.5,
37 emin=0.1, emax=100.0, rad=5.0, duration=180.0, 30 emin=0.1, emax=100.0, rad=5.0, duration=180.0,
38 - deadc=0.95,  
39 - caldb="prod2", irf="South_50h"): 31 + deadc=0.95, caldb='prod2', irf='South_0.5h'):
40 """ 32 """
41 - Returns an observation container.  
42 -  
43 - Keywords:  
44 - pattern - Pointing pattern, either "single" or "four"  
45 - ra - RA of pattern centre [deg] (default: 83.6331)  
46 - dec - DEC of pattern centre [deg] (default: 22.0145)  
47 - offset - Offset between observations of pattern [deg] (default: 1.5)  
48 - emin - Minimum energy [TeV] (default: 0.1)  
49 - emax - Maximum energy [TeV] (default: 100.0)  
50 - rad - ROI radius used for analysis [deg] (default: 5.0)  
51 - duration - Duration of one CTA observation [seconds] (default: 1800.0)  
52 - deadc - Deadtime correction factor (default: 0.95)  
53 - caldb - Calibration database path (default: "dummy")  
54 - irf - Instrument response function (default: cta_dummy_irf) 33 + Returns an observation container
  34 +
  35 + Parameters
  36 + ----------
  37 + pattern : str, optional
  38 + Pointing pattern, either 'single' or 'four'
  39 + ra : float, optional
  40 + Right Ascension of pattern centre (deg)
  41 + dec : float, optional
  42 + Declination of pattern centre (deg)
  43 + offset : float, optional
  44 + Offset between observations of pattern (deg)
  45 + emin : float, optional
  46 + Minimum energy (TeV)
  47 + emax : float, optional
  48 + Maximum energy (TeV)
  49 + rad : float, optional
  50 + ROI radius used for analysis (deg)
  51 + duration : float, optional
  52 + Duration of one CTA observation (s)
  53 + deadc : float, optional
  54 + Deadtime correction factor
  55 + caldb : str, optional
  56 + Calibration database path
  57 + irf : str, optional
  58 + Instrument response function
  59 +
  60 + Returns
  61 + -------
  62 + obs : `~gammalib.GObservations`
  63 + Observation container
55 """ 64 """
56 # Set list of observations 65 # Set list of observations
57 obs_def_list = obsutils.set_obs_patterns(pattern, 66 obs_def_list = obsutils.set_obs_patterns(pattern,
@@ -77,10 +86,19 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, @@ -77,10 +86,19 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5,
77 # =========== # 86 # =========== #
78 def setup_model(obs, model='data/crab.xml'): 87 def setup_model(obs, model='data/crab.xml'):
79 """ 88 """
80 - Setup model for analysis.  
81 -  
82 - Keywords:  
83 - model - Model Xml file 89 + Setup model for analysis
  90 +
  91 + Parameters
  92 + ----------
  93 + obs : `~gammalib.GObservations`
  94 + Observation container
  95 + model : str, optional
  96 + Model definition XML file
  97 +
  98 + Returns
  99 + -------
  100 + obs : `~gammalib.GObservations`
  101 + Observation container
84 """ 102 """
85 # Append model from file to observation container 103 # Append model from file to observation container
86 obs.models(gammalib.GModels(model)) 104 obs.models(gammalib.GModels(model))
@@ -94,24 +112,36 @@ def setup_model(obs, model='data/crab.xml'): @@ -94,24 +112,36 @@ def setup_model(obs, model='data/crab.xml'):
94 # ================================ # 112 # ================================ #
95 def run_pipeline(obs, emin=0.1, emax=100.0, 113 def run_pipeline(obs, emin=0.1, emax=100.0,
96 enumbins=20, nxpix=200, nypix=200, binsz=0.02, 114 enumbins=20, nxpix=200, nypix=200, binsz=0.02,
97 - coordsys="CEL", proj="CAR", debug=False): 115 + coordsys='CEL', proj='CAR', debug=False):
98 """ 116 """
99 - Simulation and binned analysis pipeline.  
100 -  
101 - Keywords:  
102 - emin - Minimum energy of cube [TeV] (default: 0.1)  
103 - emax - Maximum energy of cube [TeV] (default: 100.0)  
104 - enumbins - Number of energy bins in cube (default: 20)  
105 - nxpix - Number of RA pixels in cube (default: 200)  
106 - nypix - Number of DEC pixels in cube (default: 200)  
107 - binsz - Spatial cube bin size [deg] (default: 0.02)  
108 - coordsys - Cube coordinate system (CEL or GAL)  
109 - proj - Cube World Coordinate System (WCS) projection  
110 - debug - Enable debugging (default: False) 117 + Simulation and binned analysis pipeline
  118 +
  119 + Parameters
  120 + ----------
  121 + obs : `~gammalib.GObservations`
  122 + Observation container
  123 + emin : float, optional
  124 + Minimum energy (TeV)
  125 + emax : float, optional
  126 + Maximum energy (TeV)
  127 + enumbins : int, optional
  128 + Number of energy bins
  129 + nxpix : int, optional
  130 + Number of pixels in X axis
  131 + nypix : int, optional
  132 + Number of pixels in Y axis
  133 + binsz : float, optional
  134 + Pixel size (deg)
  135 + coordsys : str, optional
  136 + Coordinate system
  137 + proj : str, optional
  138 + Coordinate projection
  139 + debug : bool, optional
  140 + Debug function
111 """ 141 """
112 # Simulate events 142 # Simulate events
113 sim = ctools.ctobssim(obs) 143 sim = ctools.ctobssim(obs)
114 - sim["debug"] = debug 144 + sim['debug'] = debug
115 sim.run() 145 sim.run()
116 146
117 # Bin events by looping over all observations in the container 147 # Bin events by looping over all observations in the container
@@ -125,16 +155,16 @@ def run_pipeline(obs, emin=0.1, emax=100.0, @@ -125,16 +155,16 @@ def run_pipeline(obs, emin=0.1, emax=100.0,
125 155
126 # Bin events for that observation 156 # Bin events for that observation
127 bin = ctools.ctbin(container) 157 bin = ctools.ctbin(container)
128 - bin["ebinalg"] = "LOG"  
129 - bin["emin"] = emin  
130 - bin["emax"] = emax  
131 - bin["enumbins"] = enumbins  
132 - bin["nxpix"] = nxpix  
133 - bin["nypix"] = nypix  
134 - bin["binsz"] = binsz  
135 - bin["coordsys"] = coordsys  
136 - bin["usepnt"] = True  
137 - bin["proj"] = proj 158 + bin['ebinalg'] = 'LOG'
  159 + bin['emin'] = emin
  160 + bin['emax'] = emax
  161 + bin['enumbins'] = enumbins
  162 + bin['nxpix'] = nxpix
  163 + bin['nypix'] = nypix
  164 + bin['binsz'] = binsz
  165 + bin['coordsys'] = coordsys
  166 + bin['usepnt'] = True
  167 + bin['proj'] = proj
138 bin.run() 168 bin.run()
139 169
140 # Append result to observations 170 # Append result to observations
@@ -142,7 +172,7 @@ def run_pipeline(obs, emin=0.1, emax=100.0, @@ -142,7 +172,7 @@ def run_pipeline(obs, emin=0.1, emax=100.0,
142 172
143 # Perform maximum likelihood fitting 173 # Perform maximum likelihood fitting
144 like = ctools.ctlike(obs) 174 like = ctools.ctlike(obs)
145 - like["debug"] = True # Switch this always on for results in console 175 + like['debug'] = True # Switch this always on for results in console
146 like.run() 176 like.run()
147 177
148 # Return 178 # Return
@@ -155,9 +185,9 @@ def run_pipeline(obs, emin=0.1, emax=100.0, @@ -155,9 +185,9 @@ def run_pipeline(obs, emin=0.1, emax=100.0,
155 if __name__ == '__main__': 185 if __name__ == '__main__':
156 186
157 # Dump header 187 # Dump header
158 - print("********************************************")  
159 - print("* CTA binned analysis pipeline *")  
160 - print("********************************************") 188 + print('********************************************')
  189 + print('* CTA binned analysis pipeline *')
  190 + print('********************************************')
161 191
162 # Setup observations 192 # Setup observations
163 obs = setup_observations() 193 obs = setup_observations()
examples/pipeline_stacked_disk.py
  View file @1c0d6ef
1 #! /usr/bin/env python 1 #! /usr/bin/env python
2 # ========================================================================== 2 # ==========================================================================
3 -# This script illustrates how to perform a stacked CTA analysis based on  
4 -# simulated CTA data. You may use and adapt this script to implement your  
5 -# own pipeline. 3 +# Perform stacked analysis of simulated CTA data.
6 # 4 #
7 -# Usage:  
8 -# ./pipeline_stacked_disk.py  
9 -#  
10 -# ==========================================================================  
11 -#  
12 -# Copyright (C) 2015 Juergen Knoedlseder 5 +# Copyright (C) 2015-2016 Juergen Knoedlseder
13 # 6 #
14 # This program is free software: you can redistribute it and/or modify 7 # This program is free software: you can redistribute it and/or modify
15 # it under the terms of the GNU General Public License as published by 8 # it under the terms of the GNU General Public License as published by
@@ -33,25 +26,41 @@ from cscripts import obsutils @@ -33,25 +26,41 @@ from cscripts import obsutils
33 # ================== # 26 # ================== #
34 # Setup observations # 27 # Setup observations #
35 # ================== # 28 # ================== #
36 -def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, 29 +def setup_observations(pattern='four', ra=83.63, dec=22.01, offset=1.5,
37 emin=0.1, emax=100.0, rad=5.0, duration=180.0, 30 emin=0.1, emax=100.0, rad=5.0, duration=180.0,
38 - deadc=0.95,  
39 - caldb="prod2", irf="South_50h"): 31 + deadc=0.95, caldb='prod2', irf='South_0.5h'):
40 """ 32 """
41 - Returns an observation container.  
42 -  
43 - Keywords:  
44 - pattern - Pointing pattern, either "single" or "four"  
45 - ra - RA of pattern centre [deg] (default: 83.6331)  
46 - dec - DEC of pattern centre [deg] (default: 22.0145)  
47 - offset - Offset between observations of pattern [deg] (default: 1.5)  
48 - emin - Minimum energy [TeV] (default: 0.1)  
49 - emax - Maximum energy [TeV] (default: 100.0)  
50 - rad - ROI radius used for analysis [deg] (default: 5.0)  
51 - duration - Duration of one CTA observation [seconds] (default: 1800.0)  
52 - deadc - Deadtime correction factor (default: 0.95)  
53 - caldb - Calibration database path (default: "dummy")  
54 - irf - Instrument response function (default: cta_dummy_irf) 33 + Returns an observation container
  34 +
  35 + Parameters
  36 + ----------
  37 + pattern : str, optional
  38 + Pointing pattern, either 'single' or 'four'
  39 + ra : float, optional
  40 + Right Ascension of pattern centre (deg)
  41 + dec : float, optional
  42 + Declination of pattern centre (deg)
  43 + offset : float, optional
  44 + Offset between observations of pattern (deg)
  45 + emin : float, optional
  46 + Minimum energy (TeV)
  47 + emax : float, optional
  48 + Maximum energy (TeV)
  49 + rad : float, optional
  50 + ROI radius used for analysis (deg)
  51 + duration : float, optional
  52 + Duration of one CTA observation (s)
  53 + deadc : float, optional
  54 + Deadtime correction factor
  55 + caldb : str, optional
  56 + Calibration database path
  57 + irf : str, optional
  58 + Instrument response function
  59 +
  60 + Returns
  61 + -------
  62 + obs : `~gammalib.GObservations`
  63 + Observation container
55 """ 64 """
56 # Set list of observations 65 # Set list of observations
57 obs_def_list = obsutils.set_obs_patterns(pattern, 66 obs_def_list = obsutils.set_obs_patterns(pattern,
@@ -75,12 +84,21 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, @@ -75,12 +84,21 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5,
75 # =========== # 84 # =========== #
76 # Setup model # 85 # Setup model #
77 # =========== # 86 # =========== #
78 -def setup_model(obs, model="data/crab.xml"): 87 +def setup_model(obs, model='data/crab.xml'):
79 """ 88 """
80 - Setup model for analysis. 89 + Setup model for analysis
  90 +
  91 + Parameters
  92 + ----------
  93 + obs : `~gammalib.GObservations`
  94 + Observation container
  95 + model : str, optional
  96 + Model definition XML file
81 97
82 - Keywords:  
83 - model - Model Xml file 98 + Returns
  99 + -------
  100 + obs : `~gammalib.GObservations`
  101 + Observation container
84 """ 102 """
85 # Append model from file to observation container 103 # Append model from file to observation container
86 obs.models(gammalib.GModels(model)) 104 obs.models(gammalib.GModels(model))
@@ -94,126 +112,146 @@ def setup_model(obs, model="data/crab.xml"): @@ -94,126 +112,146 @@ def setup_model(obs, model="data/crab.xml"):
94 # ================================ # 112 # ================================ #
95 def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0, 113 def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0,
96 enumbins=20, nxpix=200, nypix=200, binsz=0.02, 114 enumbins=20, nxpix=200, nypix=200, binsz=0.02,
97 - coordsys="CEL", proj="CAR",  
98 - model="data/crab.xml",  
99 - caldb="prod2", irf="South_50h", 115 + coordsys='CEL', proj='CAR',
  116 + model='data/crab.xml',
  117 + caldb='prod2', irf='South_0.5h',
100 debug=False): 118 debug=False):
101 """ 119 """
102 - Simulation and stacked analysis pipeline.  
103 -  
104 - Keywords:  
105 - ra - RA of cube centre [deg] (default: 83.6331)  
106 - dec - DEC of cube centre [deg] (default: 22.0145)  
107 - emin - Minimum energy of cube [TeV] (default: 0.1)  
108 - emax - Maximum energy of cube [TeV] (default: 100.0)  
109 - enumbins - Number of energy bins in cube (default: 20)  
110 - nxpix - Number of RA pixels in cube (default: 200)  
111 - nypix - Number of DEC pixels in cube (default: 200)  
112 - binsz - Spatial cube bin size [deg] (default: 0.02)  
113 - coordsys - Cube coordinate system (CEL or GAL)  
114 - proj - Cube World Coordinate System (WCS) projection  
115 - debug - Enable debugging (default: False) 120 + Simulation and stacked analysis pipeline
  121 +
  122 + Parameters
  123 + ----------
  124 + obs : `~gammalib.GObservations`
  125 + Observation container
  126 + ra : float, optional
  127 + Right Ascension of counts cube centre (deg)
  128 + dec : float, optional
  129 + Declination of Region of counts cube centre (deg)
  130 + emin : float, optional
  131 + Minimum energy (TeV)
  132 + emax : float, optional
  133 + Maximum energy (TeV)
  134 + enumbins : int, optional
  135 + Number of energy bins
  136 + nxpix : int, optional
  137 + Number of pixels in X axis
  138 + nypix : int, optional
  139 + Number of pixels in Y axis
  140 + binsz : float, optional
  141 + Pixel size (deg)
  142 + coordsys : str, optional
  143 + Coordinate system
  144 + proj : str, optional
  145 + Coordinate projection
  146 + model : str, optional
  147 + Model definition XML file
  148 + caldb : str, optional
  149 + Calibration database path
  150 + irf : str, optional
  151 + Instrument response function
  152 + debug : bool, optional
  153 + Debug function
116 """ 154 """
117 # Simulate events 155 # Simulate events
118 sim = ctools.ctobssim(obs) 156 sim = ctools.ctobssim(obs)
119 - sim["debug"] = debug  
120 - sim["outevents"] = "obs.xml" 157 + sim['debug'] = debug
  158 + sim['outevents'] = 'obs.xml'
121 sim.execute() 159 sim.execute()
122 160
123 # Bin events into counts map 161 # Bin events into counts map
124 bin = ctools.ctbin() 162 bin = ctools.ctbin()
125 - bin["inobs"] = "obs.xml"  
126 - bin["outcube"] = "cntcube.fits"  
127 - bin["ebinalg"] = "LOG"  
128 - bin["emin"] = emin  
129 - bin["emax"] = emax  
130 - bin["enumbins"] = enumbins  
131 - bin["nxpix"] = nxpix  
132 - bin["nypix"] = nypix  
133 - bin["binsz"] = binsz  
134 - bin["coordsys"] = coordsys  
135 - bin["proj"] = proj  
136 - bin["xref"] = ra  
137 - bin["yref"] = dec  
138 - bin["debug"] = debug 163 + bin['inobs'] = 'obs.xml'
  164 + bin['outcube'] = 'cntcube.fits'
  165 + bin['ebinalg'] = 'LOG'
  166 + bin['emin'] = emin
  167 + bin['emax'] = emax
  168 + bin['enumbins'] = enumbins
  169 + bin['nxpix'] = nxpix
  170 + bin['nypix'] = nypix
  171 + bin['binsz'] = binsz
  172 + bin['coordsys'] = coordsys
  173 + bin['proj'] = proj
  174 + bin['xref'] = ra
  175 + bin['yref'] = dec
  176 + bin['debug'] = debug
139 bin.execute() 177 bin.execute()
140 178
141 # Create exposure cube 179 # Create exposure cube
142 expcube = ctools.ctexpcube() 180 expcube = ctools.ctexpcube()
143 - expcube["inobs"] = "obs.xml"  
144 - expcube["incube"] = "cntcube.fits"  
145 - expcube["outcube"] = "expcube.fits"  
146 - expcube["caldb"] = caldb  
147 - expcube["irf"] = irf  
148 - expcube["ebinalg"] = "LOG"  
149 - expcube["emin"] = emin  
150 - expcube["emax"] = emax  
151 - expcube["enumbins"] = enumbins  
152 - expcube["nxpix"] = nxpix  
153 - expcube["nypix"] = nypix  
154 - expcube["binsz"] = binsz  
155 - expcube["coordsys"] = coordsys  
156 - expcube["proj"] = proj  
157 - expcube["xref"] = ra  
158 - expcube["yref"] = dec  
159 - expcube["debug"] = debug 181 + expcube['inobs'] = 'obs.xml'
  182 + expcube['incube'] = 'cntcube.fits'
  183 + expcube['outcube'] = 'expcube.fits'
  184 + expcube['caldb'] = caldb
  185 + expcube['irf'] = irf
  186 + expcube['ebinalg'] = 'LOG'
  187 + expcube['emin'] = emin
  188 + expcube['emax'] = emax
  189 + expcube['enumbins'] = enumbins
  190 + expcube['nxpix'] = nxpix
  191 + expcube['nypix'] = nypix
  192 + expcube['binsz'] = binsz
  193 + expcube['coordsys'] = coordsys
  194 + expcube['proj'] = proj
  195 + expcube['xref'] = ra
  196 + expcube['yref'] = dec
  197 + expcube['debug'] = debug
160 expcube.execute() 198 expcube.execute()
161 199
162 # Create PSF cube 200 # Create PSF cube
163 psfcube = ctools.ctpsfcube() 201 psfcube = ctools.ctpsfcube()
164 - psfcube["inobs"] = "obs.xml"  
165 - psfcube["incube"] = "NONE"  
166 - psfcube["outcube"] = "psfcube.fits"  
167 - psfcube["caldb"] = caldb  
168 - psfcube["irf"] = irf  
169 - psfcube["ebinalg"] = "LOG"  
170 - psfcube["emin"] = emin  
171 - psfcube["emax"] = emax  
172 - psfcube["enumbins"] = enumbins  
173 - psfcube["nxpix"] = 10  
174 - psfcube["nypix"] = 10  
175 - psfcube["binsz"] = 1.0  
176 - psfcube["coordsys"] = coordsys  
177 - psfcube["proj"] = proj  
178 - psfcube["xref"] = ra  
179 - psfcube["yref"] = dec  
180 - psfcube["debug"] = debug 202 + psfcube['inobs'] = 'obs.xml'
  203 + psfcube['incube'] = 'NONE'
  204 + psfcube['outcube'] = 'psfcube.fits'
  205 + psfcube['caldb'] = caldb
  206 + psfcube['irf'] = irf
  207 + psfcube['ebinalg'] = 'LOG'
  208 + psfcube['emin'] = emin
  209 + psfcube['emax'] = emax
  210 + psfcube['enumbins'] = enumbins
  211 + psfcube['nxpix'] = 10
  212 + psfcube['nypix'] = 10
  213 + psfcube['binsz'] = 1.0
  214 + psfcube['coordsys'] = coordsys
  215 + psfcube['proj'] = proj
  216 + psfcube['xref'] = ra
  217 + psfcube['yref'] = dec
  218 + psfcube['debug'] = debug
181 psfcube.execute() 219 psfcube.execute()
182 220
183 # Create background cube 221 # Create background cube
184 bkgcube = ctools.ctbkgcube() 222 bkgcube = ctools.ctbkgcube()
185 - bkgcube["inobs"] = "obs.xml"  
186 - bkgcube["inmodel"] = model  
187 - bkgcube["incube"] = "cntcube.fits"  
188 - bkgcube["outcube"] = "bkgcube.fits"  
189 - bkgcube["outmodel"] = "model_bkg.xml"  
190 - bkgcube["caldb"] = caldb  
191 - bkgcube["irf"] = irf  
192 - bkgcube["ebinalg"] = "LOG"  
193 - bkgcube["emin"] = emin  
194 - bkgcube["emax"] = emax  
195 - bkgcube["enumbins"] = enumbins  
196 - bkgcube["nxpix"] = 10  
197 - bkgcube["nypix"] = 10  
198 - bkgcube["binsz"] = 1.0  
199 - bkgcube["coordsys"] = coordsys  
200 - bkgcube["proj"] = proj  
201 - bkgcube["xref"] = ra  
202 - bkgcube["yref"] = dec  
203 - bkgcube["debug"] = debug 223 + bkgcube['inobs'] = 'obs.xml'
  224 + bkgcube['inmodel'] = model
  225 + bkgcube['incube'] = 'NONE'
  226 + bkgcube['outcube'] = 'bkgcube.fits'
  227 + bkgcube['outmodel'] = 'model_bkg.xml'
  228 + bkgcube['caldb'] = caldb
  229 + bkgcube['irf'] = irf
  230 + bkgcube['ebinalg'] = 'LOG'
  231 + bkgcube['emin'] = emin
  232 + bkgcube['emax'] = emax
  233 + bkgcube['enumbins'] = enumbins
  234 + bkgcube['nxpix'] = 10
  235 + bkgcube['nypix'] = 10
  236 + bkgcube['binsz'] = 1.0
  237 + bkgcube['coordsys'] = coordsys
  238 + bkgcube['proj'] = proj
  239 + bkgcube['xref'] = ra
  240 + bkgcube['yref'] = dec
  241 + bkgcube['debug'] = debug
204 bkgcube.execute() 242 bkgcube.execute()
205 243
206 # Perform maximum likelihood fitting 244 # Perform maximum likelihood fitting
207 like = ctools.ctlike() 245 like = ctools.ctlike()
208 - like["inobs"] = "cntcube.fits"  
209 - like["inmodel"] = "model_bkg.xml"  
210 - like["outmodel"] = "fit_results.xml"  
211 - like["expcube"] = "expcube.fits"  
212 - like["psfcube"] = "psfcube.fits"  
213 - like["bkgcube"] = "bkgcube.fits"  
214 - like["caldb"] = caldb  
215 - like["irf"] = irf  
216 - like["debug"] = True # Switch this always on for results in console 246 + like['inobs'] = 'cntcube.fits'
  247 + like['inmodel'] = 'model_bkg.xml'
  248 + like['outmodel'] = 'fit_results.xml'
  249 + like['expcube'] = 'expcube.fits'
  250 + like['psfcube'] = 'psfcube.fits'
  251 + like['bkgcube'] = 'bkgcube.fits'
  252 + like['caldb'] = caldb
  253 + like['irf'] = irf
  254 + like['debug'] = True # Switch this always on for results in console
217 like.execute() 255 like.execute()
218 256
219 # Return 257 # Return
@@ -226,9 +264,9 @@ def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0, @@ -226,9 +264,9 @@ def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0,
226 if __name__ == '__main__': 264 if __name__ == '__main__':
227 265
228 # Dump header 266 # Dump header
229 - print("********************************************")  
230 - print("* CTA stacked analysis pipeline *")  
231 - print("********************************************") 267 + print('********************************************')
  268 + print('* CTA stacked analysis pipeline *')
  269 + print('********************************************')
232 270
233 # Setup observations 271 # Setup observations
234 obs = setup_observations() 272 obs = setup_observations()
examples/pipeline_stacked_mem.py
  View file @1c0d6ef
1 #! /usr/bin/env python 1 #! /usr/bin/env python
2 # ========================================================================== 2 # ==========================================================================
3 -# This script illustrates how to perform a stacked CTA analysis based on  
4 -# simulated CTA data. You may use and adapt this script to implement your  
5 -# own pipeline. 3 +# Perform stacked in-memory analysis of simulated CTA data
6 # 4 #
7 -# Usage:  
8 -# ./pipeline_stacked_mem.py  
9 -#  
10 -# ==========================================================================  
11 -#  
12 -# Copyright (C) 2014-2015 Juergen Knoedlseder 5 +# Copyright (C) 2014-2016 Juergen Knoedlseder
13 # 6 #
14 # This program is free software: you can redistribute it and/or modify 7 # This program is free software: you can redistribute it and/or modify
15 # it under the terms of the GNU General Public License as published by 8 # it under the terms of the GNU General Public License as published by
@@ -33,25 +26,41 @@ from cscripts import obsutils @@ -33,25 +26,41 @@ from cscripts import obsutils
33 # ================== # 26 # ================== #
34 # Setup observations # 27 # Setup observations #
35 # ================== # 28 # ================== #
36 -def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, 29 +def setup_observations(pattern='four', ra=83.63, dec=22.01, offset=1.5,
37 emin=0.1, emax=100.0, rad=5.0, duration=180.0, 30 emin=0.1, emax=100.0, rad=5.0, duration=180.0,
38 - deadc=0.95,  
39 - caldb="prod2", irf="South_50h"): 31 + deadc=0.95, caldb='prod2', irf='South_0.5h'):
40 """ 32 """
41 - Returns an observation container.  
42 -  
43 - Keywords:  
44 - pattern - Pointing pattern, either "single" or "four"  
45 - ra - RA of pattern centre [deg] (default: 83.6331)  
46 - dec - DEC of pattern centre [deg] (default: 22.0145)  
47 - offset - Offset between observations of pattern [deg] (default: 1.5)  
48 - emin - Minimum energy [TeV] (default: 0.1)  
49 - emax - Maximum energy [TeV] (default: 100.0)  
50 - rad - ROI radius used for analysis [deg] (default: 5.0)  
51 - duration - Duration of one CTA observation [seconds] (default: 1800.0)  
52 - deadc - Deadtime correction factor (default: 0.95)  
53 - caldb - Calibration database path (default: "dummy")  
54 - irf - Instrument response function (default: cta_dummy_irf) 33 + Returns an observation container
  34 +
  35 + Parameters
  36 + ----------
  37 + pattern : str, optional
  38 + Pointing pattern, either 'single' or 'four'
  39 + ra : float, optional
  40 + Right Ascension of pattern centre (deg)
  41 + dec : float, optional
  42 + Declination of pattern centre (deg)
  43 + offset : float, optional
  44 + Offset between observations of pattern (deg)
  45 + emin : float, optional
  46 + Minimum energy (TeV)
  47 + emax : float, optional
  48 + Maximum energy (TeV)
  49 + rad : float, optional
  50 + ROI radius used for analysis (deg)
  51 + duration : float, optional
  52 + Duration of one CTA observation (s)
  53 + deadc : float, optional
  54 + Deadtime correction factor
  55 + caldb : str, optional
  56 + Calibration database path
  57 + irf : str, optional
  58 + Instrument response function
  59 +
  60 + Returns
  61 + -------
  62 + obs : `~gammalib.GObservations`
  63 + Observation container
55 """ 64 """
56 # Set list of observations 65 # Set list of observations
57 obs_def_list = obsutils.set_obs_patterns(pattern, 66 obs_def_list = obsutils.set_obs_patterns(pattern,
@@ -75,12 +84,21 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, @@ -75,12 +84,21 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5,
75 # =========== # 84 # =========== #
76 # Setup model # 85 # Setup model #
77 # =========== # 86 # =========== #
78 -def setup_model(obs, model="data/crab.xml"): 87 +def setup_model(obs, model='data/crab.xml'):
79 """ 88 """
80 - Setup model for analysis.  
81 -  
82 - Keywords:  
83 - model - Model Xml file 89 + Setup model for analysis
  90 +
  91 + Parameters
  92 + ----------
  93 + obs : `~gammalib.GObservations`
  94 + Observation container
  95 + model : str, optional
  96 + Model definition XML file
  97 +
  98 + Returns
  99 + -------
  100 + obs : `~gammalib.GObservations`
  101 + Observation container
84 """ 102 """
85 # Append model from file to observation container 103 # Append model from file to observation container
86 obs.models(gammalib.GModels(model)) 104 obs.models(gammalib.GModels(model))
@@ -94,93 +112,107 @@ def setup_model(obs, model="data/crab.xml"): @@ -94,93 +112,107 @@ def setup_model(obs, model="data/crab.xml"):
94 # ================================ # 112 # ================================ #
95 def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0, 113 def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0,
96 enumbins=20, nxpix=200, nypix=200, binsz=0.02, 114 enumbins=20, nxpix=200, nypix=200, binsz=0.02,
97 - coordsys="CEL", proj="CAR", debug=False): 115 + coordsys='CEL', proj='CAR', debug=False):
98 """ 116 """
99 - Simulation and stacked analysis pipeline.  
100 -  
101 - Keywords:  
102 - ra - RA of cube centre [deg] (default: 83.6331)  
103 - dec - DEC of cube centre [deg] (default: 22.0145)  
104 - emin - Minimum energy of cube [TeV] (default: 0.1)  
105 - emax - Maximum energy of cube [TeV] (default: 100.0)  
106 - enumbins - Number of energy bins in cube (default: 20)  
107 - nxpix - Number of RA pixels in cube (default: 200)  
108 - nypix - Number of DEC pixels in cube (default: 200)  
109 - binsz - Spatial cube bin size [deg] (default: 0.02)  
110 - coordsys - Cube coordinate system (CEL or GAL)  
111 - proj - Cube World Coordinate System (WCS) projection  
112 - debug - Enable debugging (default: False) 117 + Simulation and stacked analysis pipeline
  118 +
  119 + Parameters
  120 + ----------
  121 + obs : `~gammalib.GObservations`
  122 + Observation container
  123 + ra : float, optional
  124 + Right Ascension of counts cube centre (deg)
  125 + dec : float, optional
  126 + Declination of Region of counts cube centre (deg)
  127 + emin : float, optional
  128 + Minimum energy (TeV)
  129 + emax : float, optional
  130 + Maximum energy (TeV)
  131 + enumbins : int, optional
  132 + Number of energy bins
  133 + nxpix : int, optional
  134 + Number of pixels in X axis
  135 + nypix : int, optional
  136 + Number of pixels in Y axis
  137 + binsz : float, optional
  138 + Pixel size (deg)
  139 + coordsys : str, optional
  140 + Coordinate system
  141 + proj : str, optional
  142 + Coordinate projection
  143 + debug : bool, optional
  144 + Debug function
113 """ 145 """
114 # Simulate events 146 # Simulate events
115 sim = ctools.ctobssim(obs) 147 sim = ctools.ctobssim(obs)
116 - sim["debug"] = debug 148 + sim['debug'] = debug
117 sim.run() 149 sim.run()
118 150
119 # Bin events into counts map 151 # Bin events into counts map
120 bin = ctools.ctbin(sim.obs()) 152 bin = ctools.ctbin(sim.obs())
121 - bin["ebinalg"] = "LOG"  
122 - bin["emin"] = emin  
123 - bin["emax"] = emax  
124 - bin["enumbins"] = enumbins  
125 - bin["nxpix"] = nxpix  
126 - bin["nypix"] = nypix  
127 - bin["binsz"] = binsz  
128 - bin["coordsys"] = coordsys  
129 - bin["proj"] = proj  
130 - bin["xref"] = ra  
131 - bin["yref"] = dec  
132 - bin["debug"] = debug 153 + bin['ebinalg'] = 'LOG'
  154 + bin['emin'] = emin
  155 + bin['emax'] = emax
  156 + bin['enumbins'] = enumbins
  157 + bin['nxpix'] = nxpix
  158 + bin['nypix'] = nypix
  159 + bin['binsz'] = binsz
  160 + bin['coordsys'] = coordsys
  161 + bin['proj'] = proj
  162 + bin['xref'] = ra
  163 + bin['yref'] = dec
  164 + bin['debug'] = debug
133 bin.run() 165 bin.run()
134 166
135 # Create exposure cube 167 # Create exposure cube
136 expcube = ctools.ctexpcube(sim.obs()) 168 expcube = ctools.ctexpcube(sim.obs())
137 - expcube["incube"] = "NONE"  
138 - expcube["ebinalg"] = "LOG"  
139 - expcube["emin"] = emin  
140 - expcube["emax"] = emax  
141 - expcube["enumbins"] = enumbins  
142 - expcube["nxpix"] = nxpix  
143 - expcube["nypix"] = nypix  
144 - expcube["binsz"] = binsz  
145 - expcube["coordsys"] = coordsys  
146 - expcube["proj"] = proj  
147 - expcube["xref"] = ra  
148 - expcube["yref"] = dec  
149 - expcube["debug"] = debug 169 + expcube['incube'] = 'NONE'
  170 + expcube['ebinalg'] = 'LOG'
  171 + expcube['emin'] = emin
  172 + expcube['emax'] = emax
  173 + expcube['enumbins'] = enumbins
  174 + expcube['nxpix'] = nxpix
  175 + expcube['nypix'] = nypix
  176 + expcube['binsz'] = binsz
  177 + expcube['coordsys'] = coordsys
  178 + expcube['proj'] = proj
  179 + expcube['xref'] = ra
  180 + expcube['yref'] = dec
  181 + expcube['debug'] = debug
150 expcube.run() 182 expcube.run()
151 183
152 # Create PSF cube 184 # Create PSF cube
153 psfcube = ctools.ctpsfcube(sim.obs()) 185 psfcube = ctools.ctpsfcube(sim.obs())
154 - psfcube["incube"] = "NONE"  
155 - psfcube["ebinalg"] = "LOG"  
156 - psfcube["emin"] = emin  
157 - psfcube["emax"] = emax  
158 - psfcube["enumbins"] = enumbins  
159 - psfcube["nxpix"] = 10  
160 - psfcube["nypix"] = 10  
161 - psfcube["binsz"] = 1.0  
162 - psfcube["coordsys"] = coordsys  
163 - psfcube["proj"] = proj  
164 - psfcube["xref"] = ra  
165 - psfcube["yref"] = dec  
166 - psfcube["debug"] = debug 186 + psfcube['incube'] = 'NONE'
  187 + psfcube['ebinalg'] = 'LOG'
  188 + psfcube['emin'] = emin
  189 + psfcube['emax'] = emax
  190 + psfcube['enumbins'] = enumbins
  191 + psfcube['nxpix'] = 10
  192 + psfcube['nypix'] = 10
  193 + psfcube['binsz'] = 1.0
  194 + psfcube['coordsys'] = coordsys
  195 + psfcube['proj'] = proj
  196 + psfcube['xref'] = ra
  197 + psfcube['yref'] = dec
  198 + psfcube['debug'] = debug
167 psfcube.run() 199 psfcube.run()
168 200
169 # Create background cube 201 # Create background cube
170 bkgcube = ctools.ctbkgcube(sim.obs()) 202 bkgcube = ctools.ctbkgcube(sim.obs())
171 - bkgcube["incube"] = "NONE"  
172 - bkgcube["ebinalg"] = "LOG"  
173 - bkgcube["emin"] = emin  
174 - bkgcube["emax"] = emax  
175 - bkgcube["enumbins"] = enumbins  
176 - bkgcube["nxpix"] = 10  
177 - bkgcube["nypix"] = 10  
178 - bkgcube["binsz"] = 1.0  
179 - bkgcube["coordsys"] = coordsys  
180 - bkgcube["proj"] = proj  
181 - bkgcube["xref"] = ra  
182 - bkgcube["yref"] = dec  
183 - bkgcube["debug"] = debug 203 + bkgcube['incube'] = 'NONE'
  204 + bkgcube['ebinalg'] = 'LOG'
  205 + bkgcube['emin'] = emin
  206 + bkgcube['emax'] = emax
  207 + bkgcube['enumbins'] = enumbins
  208 + bkgcube['nxpix'] = 10
  209 + bkgcube['nypix'] = 10
  210 + bkgcube['binsz'] = 1.0
  211 + bkgcube['coordsys'] = coordsys
  212 + bkgcube['proj'] = proj
  213 + bkgcube['xref'] = ra
  214 + bkgcube['yref'] = dec
  215 + bkgcube['debug'] = debug
184 bkgcube.run() 216 bkgcube.run()
185 217
186 # Attach background model to observation container 218 # Attach background model to observation container
@@ -192,7 +224,7 @@ def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0, @@ -192,7 +224,7 @@ def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0,
192 224
193 # Perform maximum likelihood fitting 225 # Perform maximum likelihood fitting
194 like = ctools.ctlike(bin.obs()) 226 like = ctools.ctlike(bin.obs())
195 - like["debug"] = True # Switch this always on for results in console 227 + like['debug'] = True # Switch this always on for results in console
196 like.run() 228 like.run()
197 229
198 # Return 230 # Return
@@ -205,9 +237,9 @@ def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0, @@ -205,9 +237,9 @@ def run_pipeline(obs, ra=83.63, dec=22.01, emin=0.1, emax=100.0,
205 if __name__ == '__main__': 237 if __name__ == '__main__':
206 238
207 # Dump header 239 # Dump header
208 - print("********************************************")  
209 - print("* CTA stacked analysis pipeline *")  
210 - print("********************************************") 240 + print('********************************************')
  241 + print('* CTA stacked analysis pipeline *')
  242 + print('********************************************')
211 243
212 # Setup observations 244 # Setup observations
213 obs = setup_observations() 245 obs = setup_observations()
examples/pipeline_unbinned_disk.py
  View file @1c0d6ef
1 #! /usr/bin/env python 1 #! /usr/bin/env python
2 # ========================================================================== 2 # ==========================================================================
3 -# This script illustrates how to perform an unbinned CTA analysis based on  
4 -# simulated CTA data. You may use and adapt this script to implement your  
5 -# own pipeline. 3 +# Perform unbinned analysis of simulated CTA data.
6 # 4 #
7 -# Usage:  
8 -# ./pipeline_unbinned_disk.py  
9 -#  
10 -# ==========================================================================  
11 -#  
12 -# Copyright (C) 2015 Juergen Knoedlseder 5 +# Copyright (C) 2015-2016 Juergen Knoedlseder
13 # 6 #
14 # This program is free software: you can redistribute it and/or modify 7 # This program is free software: you can redistribute it and/or modify
15 # it under the terms of the GNU General Public License as published by 8 # it under the terms of the GNU General Public License as published by
@@ -33,25 +26,41 @@ from cscripts import obsutils @@ -33,25 +26,41 @@ from cscripts import obsutils
33 # ================== # 26 # ================== #
34 # Setup observations # 27 # Setup observations #
35 # ================== # 28 # ================== #
36 -def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, 29 +def setup_observations(pattern='four', ra=83.63, dec=22.01, offset=1.5,
37 emin=0.1, emax=100.0, rad=5.0, duration=180.0, 30 emin=0.1, emax=100.0, rad=5.0, duration=180.0,
38 - deadc=0.95,  
39 - caldb="prod2", irf="South_50h"): 31 + deadc=0.95, caldb='prod2', irf='South_0.5h'):
40 """ 32 """
41 - Returns an observation container.  
42 -  
43 - Keywords:  
44 - pattern - Pointing pattern, either "single" or "four"  
45 - ra - RA of pattern centre [deg] (default: 83.6331)  
46 - dec - DEC of pattern centre [deg] (default: 22.0145)  
47 - offset - Offset between observations of pattern [deg] (default: 1.5)  
48 - emin - Minimum energy [TeV] (default: 0.1)  
49 - emax - Maximum energy [TeV] (default: 100.0)  
50 - rad - ROI radius used for analysis [deg] (default: 5.0)  
51 - duration - Duration of one CTA observation [seconds] (default: 1800.0)  
52 - deadc - Deadtime correction factor (default: 0.95)  
53 - caldb - Calibration database path (default: "dummy")  
54 - irf - Instrument response function (default: cta_dummy_irf) 33 + Returns an observation container
  34 +
  35 + Parameters
  36 + ----------
  37 + pattern : str, optional
  38 + Pointing pattern, either 'single' or 'four'
  39 + ra : float, optional
  40 + Right Ascension of pattern centre (deg)
  41 + dec : float, optional
  42 + Declination of pattern centre (deg)
  43 + offset : float, optional
  44 + Offset between observations of pattern (deg)
  45 + emin : float, optional
  46 + Minimum energy (TeV)
  47 + emax : float, optional
  48 + Maximum energy (TeV)
  49 + rad : float, optional
  50 + ROI radius used for analysis (deg)
  51 + duration : float, optional
  52 + Duration of one CTA observation (s)
  53 + deadc : float, optional
  54 + Deadtime correction factor
  55 + caldb : str, optional
  56 + Calibration database path
  57 + irf : str, optional
  58 + Instrument response function
  59 +
  60 + Returns
  61 + -------
  62 + obs : `~gammalib.GObservations`
  63 + Observation container
55 """ 64 """
56 # Set list of observations 65 # Set list of observations
57 obs_def_list = obsutils.set_obs_patterns(pattern, 66 obs_def_list = obsutils.set_obs_patterns(pattern,
@@ -75,12 +84,21 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, @@ -75,12 +84,21 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5,
75 # =========== # 84 # =========== #
76 # Setup model # 85 # Setup model #
77 # =========== # 86 # =========== #
78 -def setup_model(obs, model="data/crab.xml"): 87 +def setup_model(obs, model='data/crab.xml'):
79 """ 88 """
80 - Setup model for analysis.  
81 -  
82 - Keywords:  
83 - model - Model Xml file 89 + Setup model for analysis
  90 +
  91 + Parameters
  92 + ----------
  93 + obs : `~gammalib.GObservations`
  94 + Observation container
  95 + model : str, optional
  96 + Model definition XML file
  97 +
  98 + Returns
  99 + -------
  100 + obs : `~gammalib.GObservations`
  101 + Observation container
84 """ 102 """
85 # Append model from file to observation container 103 # Append model from file to observation container
86 obs.models(gammalib.GModels(model)) 104 obs.models(gammalib.GModels(model))
@@ -95,55 +113,69 @@ def setup_model(obs, model="data/crab.xml"): @@ -95,55 +113,69 @@ def setup_model(obs, model="data/crab.xml"):
95 def run_pipeline(obs, ra=83.63, dec=22.01, rad=3.0, 113 def run_pipeline(obs, ra=83.63, dec=22.01, rad=3.0,
96 emin=0.1, emax=100.0, 114 emin=0.1, emax=100.0,
97 tmin=0.0, tmax=0.0, 115 tmin=0.0, tmax=0.0,
98 - model="data/crab.xml",  
99 - caldb="prod2", irf="South_50h", 116 + model='data/crab.xml',
  117 + caldb='prod2', irf='South_0.5h',
100 debug=False): 118 debug=False):
101 """ 119 """
102 - Simulation and unbinned analysis pipeline.  
103 -  
104 - Keywords:  
105 - ra - RA of cube centre [deg] (default: 83.63)  
106 - dec - DEC of cube centre [deg] (default: 22.01)  
107 - rad - Selection radius [deg] (default: 3.0)  
108 - emin - Minimum energy of cube [TeV] (default: 0.1)  
109 - emax - Maximum energy of cube [TeV] (default: 100.0)  
110 - tmin - Start time [MET] (default: 0.0)  
111 - tmax - Stop time [MET] (default: 0.0)  
112 - model - Model Xml file  
113 - caldb - Calibration database path (default: "dummy")  
114 - irf - Instrument response function (default: cta_dummy_irf)  
115 - debug - Enable debugging (default: False) 120 + Simulation and unbinned analysis pipeline
  121 +
  122 + Parameters
  123 + ----------
  124 + obs : `~gammalib.GObservations`
  125 + Observation container
  126 + ra : float, optional
  127 + Right Ascension of Region of Interest centre (deg)
  128 + dec : float, optional
  129 + Declination of Region of Interest centre (deg)
  130 + rad : float, optional
  131 + Radius of Region of Interest (deg)
  132 + emin : float, optional
  133 + Minimum energy (TeV)
  134 + emax : float, optional
  135 + Maximum energy (TeV)
  136 + tmin : float, optional
  137 + Start time (s)
  138 + tmax : float, optional
  139 + Stop time (s)
  140 + model : str, optional
  141 + Model definition XML file
  142 + caldb : str, optional
  143 + Calibration database path
  144 + irf : str, optional
  145 + Instrument response function
  146 + debug : bool, optional
  147 + Debug function
116 """ 148 """
117 # Get model 149 # Get model
118 150
119 # Simulate events 151 # Simulate events
120 sim = ctools.ctobssim(obs) 152 sim = ctools.ctobssim(obs)
121 - sim["debug"] = debug  
122 - sim["outevents"] = "obs.xml" 153 + sim['debug'] = debug
  154 + sim['outevents'] = 'obs.xml'
123 sim.execute() 155 sim.execute()
124 156
125 # Select events 157 # Select events
126 select = ctools.ctselect() 158 select = ctools.ctselect()
127 - select["inobs"] = "obs.xml"  
128 - select["outobs"] = "obs_selected.xml"  
129 - select["ra"] = ra  
130 - select["dec"] = dec  
131 - select["rad"] = rad  
132 - select["emin"] = emin  
133 - select["emax"] = emax  
134 - select["tmin"] = tmin  
135 - select["tmax"] = tmax  
136 - select["debug"] = debug 159 + select['inobs'] = 'obs.xml'
  160 + select['outobs'] = 'obs_selected.xml'
  161 + select['ra'] = ra
  162 + select['dec'] = dec
  163 + select['rad'] = rad
  164 + select['emin'] = emin
  165 + select['emax'] = emax
  166 + select['tmin'] = tmin
  167 + select['tmax'] = tmax
  168 + select['debug'] = debug
137 select.execute() 169 select.execute()
138 170
139 # Perform maximum likelihood fitting 171 # Perform maximum likelihood fitting
140 like = ctools.ctlike() 172 like = ctools.ctlike()
141 - like["inobs"] = "obs_selected.xml"  
142 - like["inmodel"] = model  
143 - like["outmodel"] = "fit_results.xml"  
144 - like["caldb"] = caldb  
145 - like["irf"] = irf  
146 - like["debug"] = True # Switch this always on for results in console 173 + like['inobs'] = 'obs_selected.xml'
  174 + like['inmodel'] = model
  175 + like['outmodel'] = 'fit_results.xml'
  176 + like['caldb'] = caldb
  177 + like['irf'] = irf
  178 + like['debug'] = True # Switch this always on for results in console
147 like.execute() 179 like.execute()
148 180
149 # Return 181 # Return
@@ -156,9 +188,9 @@ def run_pipeline(obs, ra=83.63, dec=22.01, rad=3.0, @@ -156,9 +188,9 @@ def run_pipeline(obs, ra=83.63, dec=22.01, rad=3.0,
156 if __name__ == '__main__': 188 if __name__ == '__main__':
157 189
158 # Dump header 190 # Dump header
159 - print("********************************************")  
160 - print("* CTA unbinned analysis pipeline *")  
161 - print("********************************************") 191 + print('********************************************')
  192 + print('* CTA unbinned analysis pipeline *')
  193 + print('********************************************')
162 194
163 # Setup observations 195 # Setup observations
164 obs = setup_observations() 196 obs = setup_observations()
examples/pipeline_unbinned_mem.py
  View file @1c0d6ef
1 #! /usr/bin/env python 1 #! /usr/bin/env python
2 # ========================================================================== 2 # ==========================================================================
3 -# This script illustrates how to perform an unbinned CTA analysis based on  
4 -# simulated CTA data. You may use and adapt this script to implement your  
5 -# own pipeline. 3 +# Perform unbinned in-memory analysis of simulated CTA data
6 # 4 #
7 -# Usage:  
8 -# ./pipeline_unbinned_mem.py  
9 -#  
10 -# ==========================================================================  
11 -#  
12 -# Copyright (C) 2015 Juergen Knoedlseder 5 +# Copyright (C) 2015-2016 Juergen Knoedlseder
13 # 6 #
14 # This program is free software: you can redistribute it and/or modify 7 # This program is free software: you can redistribute it and/or modify
15 # it under the terms of the GNU General Public License as published by 8 # it under the terms of the GNU General Public License as published by
@@ -33,25 +26,41 @@ from cscripts import obsutils @@ -33,25 +26,41 @@ from cscripts import obsutils
33 # ================== # 26 # ================== #
34 # Setup observations # 27 # Setup observations #
35 # ================== # 28 # ================== #
36 -def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, 29 +def setup_observations(pattern='four', ra=83.63, dec=22.01, offset=1.5,
37 emin=0.1, emax=100.0, rad=5.0, duration=180.0, 30 emin=0.1, emax=100.0, rad=5.0, duration=180.0,
38 - deadc=0.95,  
39 - caldb="prod2", irf="South_50h"): 31 + deadc=0.95, caldb='prod2', irf='South_0.5h'):
40 """ 32 """
41 - Returns an observation container.  
42 -  
43 - Keywords:  
44 - pattern - Pointing pattern, either "single" or "four"  
45 - ra - RA of pattern centre [deg] (default: 83.6331)  
46 - dec - DEC of pattern centre [deg] (default: 22.0145)  
47 - offset - Offset between observations of pattern [deg] (default: 1.5)  
48 - emin - Minimum energy [TeV] (default: 0.1)  
49 - emax - Maximum energy [TeV] (default: 100.0)  
50 - rad - ROI radius used for analysis [deg] (default: 5.0)  
51 - duration - Duration of one CTA observation [seconds] (default: 1800.0)  
52 - deadc - Deadtime correction factor (default: 0.95)  
53 - caldb - Calibration database path (default: "dummy")  
54 - irf - Instrument response function (default: cta_dummy_irf) 33 + Returns an observation container
  34 +
  35 + Parameters
  36 + ----------
  37 + pattern : str, optional
  38 + Pointing pattern, either 'single' or 'four'
  39 + ra : float, optional
  40 + Right Ascension of pattern centre (deg)
  41 + dec : float, optional
  42 + Declination of pattern centre (deg)
  43 + offset : float, optional
  44 + Offset between observations of pattern (deg)
  45 + emin : float, optional
  46 + Minimum energy (TeV)
  47 + emax : float, optional
  48 + Maximum energy (TeV)
  49 + rad : float, optional
  50 + ROI radius used for analysis (deg)
  51 + duration : float, optional
  52 + Duration of one CTA observation (s)
  53 + deadc : float, optional
  54 + Deadtime correction factor
  55 + caldb : str, optional
  56 + Calibration database path
  57 + irf : str, optional
  58 + Instrument response function
  59 +
  60 + Returns
  61 + -------
  62 + obs : `~gammalib.GObservations`
  63 + Observation container
55 """ 64 """
56 # Set list of observations 65 # Set list of observations
57 obs_def_list = obsutils.set_obs_patterns(pattern, 66 obs_def_list = obsutils.set_obs_patterns(pattern,
@@ -75,12 +84,21 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5, @@ -75,12 +84,21 @@ def setup_observations(pattern="four", ra=83.63, dec=22.01, offset=1.5,
75 # =========== # 84 # =========== #
76 # Setup model # 85 # Setup model #
77 # =========== # 86 # =========== #
78 -def setup_model(obs, model="data/crab.xml"): 87 +def setup_model(obs, model='data/crab.xml'):
79 """ 88 """
80 - Setup model for analysis.  
81 -  
82 - Keywords:  
83 - model - Model Xml file 89 + Setup model for analysis
  90 +
  91 + Parameters
  92 + ----------
  93 + obs : `~gammalib.GObservations`
  94 + Observation container
  95 + model : str, optional
  96 + Model definition XML file
  97 +
  98 + Returns
  99 + -------
  100 + obs : `~gammalib.GObservations`
  101 + Observation container
84 """ 102 """
85 # Append model from file to observation container 103 # Append model from file to observation container
86 obs.models(gammalib.GModels(model)) 104 obs.models(gammalib.GModels(model))
@@ -97,38 +115,49 @@ def run_pipeline(obs, ra=83.63, dec=22.01, rad=3.0, @@ -97,38 +115,49 @@ def run_pipeline(obs, ra=83.63, dec=22.01, rad=3.0,
97 tmin=0.0, tmax=0.0, 115 tmin=0.0, tmax=0.0,
98 debug=False): 116 debug=False):
99 """ 117 """
100 - Simulation and unbinned analysis pipeline.  
101 -  
102 - Keywords:  
103 - ra - RA of cube centre [deg] (default: 83.63)  
104 - dec - DEC of cube centre [deg] (default: 22.01)  
105 - rad - Selection radius [deg] (default: 3.0)  
106 - emin - Minimum energy of cube [TeV] (default: 0.1)  
107 - emax - Maximum energy of cube [TeV] (default: 100.0)  
108 - tmin - Start time [MET] (default: 0.0)  
109 - tmax - Stop time [MET] (default: 0.0)  
110 - debug - Enable debugging (default: False) 118 + Simulation and binned analysis pipeline
  119 +
  120 + Parameters
  121 + ----------
  122 + obs : `~gammalib.GObservations`
  123 + Observation container
  124 + ra : float, optional
  125 + Right Ascension of Region of Interest centre (deg)
  126 + dec : float, optional
  127 + Declination of Region of Interest centre (deg)
  128 + rad : float, optional
  129 + Radius of Region of Interest (deg)
  130 + emin : float, optional
  131 + Minimum energy (TeV)
  132 + emax : float, optional
  133 + Maximum energy (TeV)
  134 + tmin : float, optional
  135 + Start time (s)
  136 + tmax : float, optional
  137 + Stop time (s)
  138 + debug : bool, optional
  139 + Debug function
111 """ 140 """
112 # Simulate events 141 # Simulate events
113 sim = ctools.ctobssim(obs) 142 sim = ctools.ctobssim(obs)
114 - sim["debug"] = debug 143 + sim['debug'] = debug
115 sim.run() 144 sim.run()
116 145
117 # Select events 146 # Select events
118 select = ctools.ctselect(sim.obs()) 147 select = ctools.ctselect(sim.obs())
119 - select["ra"] = ra  
120 - select["dec"] = dec  
121 - select["rad"] = rad  
122 - select["emin"] = emin  
123 - select["emax"] = emax  
124 - select["tmin"] = tmin  
125 - select["tmax"] = tmax  
126 - select["debug"] = debug 148 + select['ra'] = ra
  149 + select['dec'] = dec
  150 + select['rad'] = rad
  151 + select['emin'] = emin
  152 + select['emax'] = emax
  153 + select['tmin'] = tmin
  154 + select['tmax'] = tmax
  155 + select['debug'] = debug
127 select.run() 156 select.run()
128 157
129 # Perform maximum likelihood fitting 158 # Perform maximum likelihood fitting
130 like = ctools.ctlike(select.obs()) 159 like = ctools.ctlike(select.obs())
131 - like["debug"] = True # Switch this always on for results in console 160 + like['debug'] = True # Switch this always on for results in console
132 like.run() 161 like.run()
133 162
134 # Return 163 # Return
@@ -141,9 +170,9 @@ def run_pipeline(obs, ra=83.63, dec=22.01, rad=3.0, @@ -141,9 +170,9 @@ def run_pipeline(obs, ra=83.63, dec=22.01, rad=3.0,
141 if __name__ == '__main__': 170 if __name__ == '__main__':
142 171
143 # Dump header 172 # Dump header
144 - print("********************************************")  
145 - print("* CTA unbinned analysis pipeline *")  
146 - print("********************************************") 173 + print('********************************************')
  174 + print('* CTA unbinned analysis pipeline *')
  175 + print('********************************************')
147 176
148 # Setup observations 177 # Setup observations
149 obs = setup_observations() 178 obs = setup_observations()