Francois Brun / ctools

Commit 75e6edfb8a5045f05d7fdd19ecb872434ace7dd9

Authored by Jürgen Knödlseder
1 parent 1399aa6b

Adapt benchmark to new interface

Inline Side-by-side
Showing 1 changed file with 38 additions and 21 deletions
examples/benchmark_cta_analysis.py
  View file @75e6edf
... ... @@ -167,10 +167,10 @@ def binned_pipeline(duration):
167 167 return telapsed, tctlike
168 168  
169 169  
170   -# ============================ #
171   -# Cube-style analysis pipeline #
172   -# ============================ #
173   -def cube_pipeline(duration):
  170 +# ========================= #
  171 +# Stacked analysis pipeline #
  172 +# ========================= #
  173 +def stacked_pipeline(duration):
174 174 """
175 175 Cube-style analysis pipeline.
176 176 """
... ... @@ -260,17 +260,34 @@ def cube_pipeline(duration):
260 260 psfcube["yref"].real(dec)
261 261 psfcube.run()
262 262  
263   - # Set exposure and PSF cube for first CTA observation
264   - obs = bin.obs()
265   - cta = gammalib.GCTAObservation(obs[0])
266   - cta.response(expcube.expcube(), psfcube.psfcube())
267   - obs[0] = cta
  263 + # Create background cube
  264 + bkgcube = ctools.ctbkgcube(sim.obs())
  265 + bkgcube["incube"].filename("NONE")
  266 + bkgcube["ebinalg"].string("LOG")
  267 + bkgcube["emin"].real(emin)
  268 + bkgcube["emax"].real(emax)
  269 + bkgcube["enumbins"].integer(enumbins)
  270 + bkgcube["nxpix"].integer(10)
  271 + bkgcube["nypix"].integer(10)
  272 + bkgcube["binsz"].real(1.0)
  273 + bkgcube["coordsys"].string(coordsys)
  274 + bkgcube["proj"].string(proj)
  275 + bkgcube["xref"].real(ra)
  276 + bkgcube["yref"].real(dec)
  277 + bkgcube.run()
  278 +
  279 + # Attach background model to observation container
  280 + bin.obs().models(bkgcube.models())
  281 +
  282 + # Set Exposure and Psf cube for first CTA observation
  283 + # (ctbin will create an observation with a single container)
  284 + bin.obs()[0].response(expcube.expcube(), psfcube.psfcube(), bkgcube.bkgcube())
268 285  
269 286 # Get ctlike start CPU time
270 287 tctlike = time.clock()
271 288  
272 289 # Perform maximum likelihood fitting
273   - like = ctools.ctlike(obs)
  290 + like = ctools.ctlike(bin.obs())
274 291 like.run()
275 292  
276 293 # Get stop CPU time
... ... @@ -292,16 +309,16 @@ if __name__ == '__main__':
292 309 print("*************************************")
293 310 print("* CTA analysis benchmark *")
294 311 print("*************************************")
295   - print(" Duration Unbinned Binned Cube-style ctlike-UB ctlike-B ctlike-CS")
  312 + print(" Duration Unbinned Binned Stacked ctlike-UB ctlike-B ctlike-S")
296 313  
297 314 # Initialize arrays
298 315 a_duration = []
299 316 a_unbinned = []
300 317 a_binned = []
301   - a_cube = []
  318 + a_stacked = []
302 319 c_unbinned = []
303 320 c_binned = []
304   - c_cube = []
  321 + c_stacked = []
305 322  
306 323 # Perform benchmarks
307 324 duration = 1800.0
... ... @@ -309,22 +326,22 @@ if __name__ == '__main__':
309 326  
310 327 # Perform analyses
311 328 t_unbinned, ct_unbinned = unbinned_pipeline(duration)
312   - t_binned, ct_binned = binned_pipeline(duration)
313   - t_cube, ct_cube = cube_pipeline(duration)
  329 + t_binned, ct_binned = binned_pipeline(duration)
  330 + t_stacked, ct_stacked = stacked_pipeline(duration)
314 331  
315 332 # Collect results
316 333 a_duration.append(duration/3600.0)
317 334 a_unbinned.append(t_unbinned)
318 335 a_binned.append(t_binned)
319   - a_cube.append(t_cube)
  336 + a_stacked.append(t_stacked)
320 337 c_unbinned.append(ct_unbinned)
321 338 c_binned.append(ct_binned)
322   - c_cube.append(ct_cube)
  339 + c_stacked.append(ct_stacked)
323 340  
324 341 # Print results
325 342 print("%10.0f %10.4f %10.4f %10.4f %10.4f %10.4f %10.4f" % \
326   - (duration, t_unbinned, t_binned, t_cube, \
327   - ct_unbinned, ct_binned, ct_cube))
  343 + (duration, t_unbinned, t_binned, t_stacked, \
  344 + ct_unbinned, ct_binned, ct_stacked))
328 345  
329 346 # Increment duration
330 347 duration *= 2.0
... ... @@ -335,10 +352,10 @@ if __name__ == '__main__':
335 352 plt.title("CTA analysis benchmark")
336 353 plt.loglog(a_duration, a_unbinned, 'ro-', label='unbinned')
337 354 plt.loglog(a_duration, a_binned, 'bo-', label='binned')
338   - plt.loglog(a_duration, a_cube, 'go-', label='cube-style')
  355 + plt.loglog(a_duration, a_stacked, 'go-', label='cube-style')
339 356 plt.loglog(a_duration, c_unbinned, 'ro--')
340 357 plt.loglog(a_duration, c_binned, 'bo--')
341   - plt.loglog(a_duration, c_cube, 'go--')
  358 + plt.loglog(a_duration, c_stacked, 'go--')
342 359 plt.xlabel("Duration (hours)")
343 360 plt.ylabel("CPU time (seconds)")
344 361 plt.legend(loc="lower right")
... ...