Commit 42aae29758cf375a4709aa5d69646bb0ccd024b7

Authored by Jürgen Knödlseder
1 parent 809bc8f6

Integrate code

ChangeLog
1 -2015-02-23 Juergen Knoedlseder <jurgen.knodlseder@irap.omp.eu> 1 +2015-03-28 Juergen Knoedlseder <jurgen.knodlseder@irap.omp.eu>
2 2
3 * Version 1.0.0 released 3 * Version 1.0.0 released
4 ======================== 4 ========================
@@ -9,6 +9,7 @@ @@ -9,6 +9,7 @@
9 Modify cssens parameter interface (#1426) 9 Modify cssens parameter interface (#1426)
10 Add csobsdef script (#1429) 10 Add csobsdef script (#1429)
11 Add ctulimit tool (#1363) 11 Add ctulimit tool (#1363)
  12 + Handle background cubes for stacked analysis (#1371)
12 13
13 14
14 2015-01-08 Juergen Knoedlseder <jurgen.knodlseder@irap.omp.eu> 15 2015-01-08 Juergen Knoedlseder <jurgen.knodlseder@irap.omp.eu>
@@ -2,7 +2,7 @@ New Features and Important Changes in ctools 1.0.0 @@ -2,7 +2,7 @@ New Features and Important Changes in ctools 1.0.0
2 2
3 Author(s) : Juergen Knoedlseder 3 Author(s) : Juergen Knoedlseder
4 4
5 -23 February 2015 5 +28 March 2015
6 6
7 7
8 1. Introduction 8 1. Introduction
@@ -49,12 +49,14 @@ None @@ -49,12 +49,14 @@ None
49 49
50 8. ctlike - Maximum likelihood fitting 50 8. ctlike - Maximum likelihood fitting
51 -------------------------------------- 51 --------------------------------------
52 -None 52 +The bkgcube parameter has been added to allow specification of a background
  53 +cube for stacked analysis.
53 54
54 55
55 9. cttsmap - Generation of Test Statistics map 56 9. cttsmap - Generation of Test Statistics map
56 ---------------------------------------------- 57 ----------------------------------------------
57 -None 58 +The bkgcube parameter has been added to allow specification of a background
  59 +cube for stacked analysis.
58 60
59 61
60 10. ctskymap - CTA sky mapping 62 10. ctskymap - CTA sky mapping
@@ -64,7 +66,8 @@ None @@ -64,7 +66,8 @@ None
64 66
65 11. ctmodel - CTA counts cube modelling 67 11. ctmodel - CTA counts cube modelling
66 --------------------------------------- 68 ---------------------------------------
67 -None 69 +The bkgcube parameter has been added to allow specification of a background
  70 +cube for stacked analysis.
68 71
69 72
70 12. ctexpcube - CTA exposure cube generation 73 12. ctexpcube - CTA exposure cube generation
@@ -89,13 +92,17 @@ None @@ -89,13 +92,17 @@ None
89 92
90 16. ctbutterfly - Generate butterfly diagrams 93 16. ctbutterfly - Generate butterfly diagrams
91 --------------------------------------------- 94 ---------------------------------------------
92 -None 95 +The bkgcube parameter has been added to allow specification of a background
  96 +cube for stacked analysis.
93 97
94 98
95 17. ctulimit - Compute upper limit 99 17. ctulimit - Compute upper limit
96 ---------------------------------- 100 ----------------------------------
97 Added an initial version of the upper limit computation tool. 101 Added an initial version of the upper limit computation tool.
98 102
  103 +The bkgcube parameter has been added to allow specification of a background
  104 +cube for stacked analysis.
  105 +
99 106
100 18. Examples 107 18. Examples
101 ------------ 108 ------------
pyext/ctbkgcube.i
1 /*************************************************************************** 1 /***************************************************************************
2 * ctbkgcube - Background cube generation tool * 2 * ctbkgcube - Background cube generation tool *
3 * ----------------------------------------------------------------------- * 3 * ----------------------------------------------------------------------- *
4 - * copyright (C) 2014 by Chia-Chun Lu * 4 + * copyright (C) 2014-2015 by Chia-Chun Lu *
5 * ----------------------------------------------------------------------- * 5 * ----------------------------------------------------------------------- *
6 * * 6 * *
7 * 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 *
@@ -44,11 +44,11 @@ public: @@ -44,11 +44,11 @@ public:
44 virtual ~ctbkgcube(void); 44 virtual ~ctbkgcube(void);
45 45
46 // Methods 46 // Methods
47 - void clear(void);  
48 - void run(void);  
49 - void save(void); 47 + void clear(void);
  48 + void run(void);
  49 + void save(void);
50 const GCTACubeBackground& cube(void) const; 50 const GCTACubeBackground& cube(void) const;
51 - const GModels& models(void) const; 51 + const GModels& models(void) const;
52 }; 52 };
53 53
54 /***********************************************************************//** 54 /***********************************************************************//**
src/ctbkgcube/ctbkgcube.cpp
@@ -34,7 +34,6 @@ @@ -34,7 +34,6 @@
34 34
35 /* __ Method name definitions ____________________________________________ */ 35 /* __ Method name definitions ____________________________________________ */
36 #define G_GET_PARAMETERS "ctbkgcube::get_parameters()" 36 #define G_GET_PARAMETERS "ctbkgcube::get_parameters()"
37 -#define G_FILL_CUBE "ctbkgcube::fill_cube(GCTAObservation*)"  
38 37
39 /* __ Debug definitions __________________________________________________ */ 38 /* __ Debug definitions __________________________________________________ */
40 39
@@ -246,14 +245,13 @@ void ctbkgcube::run(void) @@ -246,14 +245,13 @@ void ctbkgcube::run(void)
246 // Copy models from observation container and reset output model 245 // Copy models from observation container and reset output model
247 // container 246 // container
248 GModels models_orig = m_obs.models(); 247 GModels models_orig = m_obs.models();
249 - m_bkgmdl = m_obs.models(); 248 + m_bkgmdl = m_obs.models();
250 m_outmdl.clear(); 249 m_outmdl.clear();
251 250
252 // Remove all models that are not CTA background models from the 251 // Remove all models that are not CTA background models from the
253 // container and put all removed components in the output 252 // container and put all removed components in the output
254 // container 253 // container
255 int num = m_bkgmdl.size(); 254 int num = m_bkgmdl.size();
256 -  
257 for (int i = num-1; i >= 0; --i) { 255 for (int i = num-1; i >= 0; --i) {
258 256
259 // Flag removal 257 // Flag removal
@@ -315,21 +313,18 @@ void ctbkgcube::run(void) @@ -315,21 +313,18 @@ void ctbkgcube::run(void)
315 // Create a background model for the output background cube and append 313 // Create a background model for the output background cube and append
316 // that model to the input model in place of the original 314 // that model to the input model in place of the original
317 // background models 315 // background models
318 - // Todo: We might think of creating the spectral model via user parameter  
319 - GCTAModelCubeBackground model(GModelSpectralConst(1.0)); 316 + // TODO: We might think of creating the spectral model via user parameter
  317 + GModelSpectralPlaw spectral(1.0, 0.0, GEnergy(1.0, "TeV"));
  318 + GCTAModelCubeBackground model(spectral);
320 319
321 // Set model name 320 // Set model name
322 model.name("BackgroundModel"); 321 model.name("BackgroundModel");
323 322
324 // Set model instrument 323 // Set model instrument
325 - // Todo: Account for possibility to have observations from different 324 + // TODO: Account for possibility to have observations from different
326 // IACTs in the same container 325 // IACTs in the same container
327 model.instruments("CTA,HESS,MAGIC,VERITAS"); 326 model.instruments("CTA,HESS,MAGIC,VERITAS");
328 327
329 - // Set model identifier to "0". This is coupled to ctbin, where the binned output  
330 - // observation id is set to 0.  
331 - model.ids("0");  
332 -  
333 // Append model to output container 328 // Append model to output container
334 m_outmdl.append(model); 329 m_outmdl.append(model);
335 330
@@ -361,7 +356,6 @@ void ctbkgcube::save(void) @@ -361,7 +356,6 @@ void ctbkgcube::save(void)
361 log.header1("Save background cube"); 356 log.header1("Save background cube");
362 } 357 }
363 358
364 -  
365 // Get output filenames 359 // Get output filenames
366 std::string outfile = (*this)["outcube"].filename(); 360 std::string outfile = (*this)["outcube"].filename();
367 std::string outmodel = (*this)["outmodel"].filename(); 361 std::string outmodel = (*this)["outmodel"].filename();
@@ -374,6 +368,7 @@ void ctbkgcube::save(void) @@ -374,6 +368,7 @@ void ctbkgcube::save(void)
374 368
375 // Save output model for binned analyses 369 // Save output model for binned analyses
376 m_outmdl.save(outmodel); 370 m_outmdl.save(outmodel);
  371 +
377 } 372 }
378 373
379 // Return 374 // Return
@@ -399,7 +394,6 @@ void ctbkgcube::init_members(void) @@ -399,7 +394,6 @@ void ctbkgcube::init_members(void)
399 m_background.clear(); 394 m_background.clear();
400 m_bkgmdl.clear(); 395 m_bkgmdl.clear();
401 m_outmdl.clear(); 396 m_outmdl.clear();
402 - m_ebounds.clear();  
403 397
404 // Return 398 // Return
405 return; 399 return;
@@ -420,7 +414,6 @@ void ctbkgcube::copy_members(const ctbkgcube&amp; app) @@ -420,7 +414,6 @@ void ctbkgcube::copy_members(const ctbkgcube&amp; app)
420 m_background = app.m_background; 414 m_background = app.m_background;
421 m_bkgmdl = app.m_bkgmdl; 415 m_bkgmdl = app.m_bkgmdl;
422 m_outmdl = app.m_outmdl; 416 m_outmdl = app.m_outmdl;
423 - m_ebounds = app.m_ebounds;  
424 417
425 // Return 418 // Return
426 return; 419 return;
@@ -498,9 +491,6 @@ void ctbkgcube::get_parameters(void) @@ -498,9 +491,6 @@ void ctbkgcube::get_parameters(void)
498 m_obs.models(models); 491 m_obs.models(models);
499 } 492 }
500 493
501 -// // Get energy definition  
502 -// m_ebounds = m_background.ebounds();  
503 -  
504 // Read output filenames (if needed) 494 // Read output filenames (if needed)
505 if (read_ahead()) { 495 if (read_ahead()) {
506 m_outcube = (*this)["outcube"].filename(); 496 m_outcube = (*this)["outcube"].filename();
src/ctbkgcube/ctbkgcube.hpp
@@ -56,11 +56,11 @@ public: @@ -56,11 +56,11 @@ public:
56 ctbkgcube& operator=(const ctbkgcube& app); 56 ctbkgcube& operator=(const ctbkgcube& app);
57 57
58 // Methods 58 // Methods
59 - void clear(void);  
60 - void run(void);  
61 - void save(void); 59 + void clear(void);
  60 + void run(void);
  61 + void save(void);
62 const GCTACubeBackground& cube(void) const; 62 const GCTACubeBackground& cube(void) const;
63 - const GModels& models(void) const; 63 + const GModels& models(void) const;
64 64
65 protected: 65 protected:
66 // Protected methods 66 // Protected methods
@@ -71,14 +71,12 @@ protected: @@ -71,14 +71,12 @@ protected:
71 void fill_cube(GCTAObservation* obs); 71 void fill_cube(GCTAObservation* obs);
72 72
73 // Protected members 73 // Protected members
74 - std::string m_outcube; //!< Filename of output cube  
75 - std::string m_outmodel; //!< Filename of output XML model  
76 - GObservations m_obs; //!< Observation container  
77 - GCTACubeBackground m_background; //!< Background cube response  
78 - GModels m_bkgmdl; //!< CTA background models  
79 - GModels m_outmdl; //!< Output models  
80 - GEbounds m_ebounds; //!< Energy boundaries  
81 - 74 + std::string m_outcube; //!< Filename of output cube
  75 + std::string m_outmodel; //!< Filename of output XML model
  76 + GObservations m_obs; //!< Observation container
  77 + GCTACubeBackground m_background; //!< Background cube response
  78 + GModels m_bkgmdl; //!< CTA background models
  79 + GModels m_outmdl; //!< Output models
82 }; 80 };
83 81
84 82
test/test_ctools.sh
@@ -179,20 +179,20 @@ $ECHO -n &quot;Test ctexpcube: &quot; @@ -179,20 +179,20 @@ $ECHO -n &quot;Test ctexpcube: &quot;
179 # 179 #
180 # Run 1 180 # Run 1
181 $ctexpcube inobs="data/crab_events.fits" \ 181 $ctexpcube inobs="data/crab_events.fits" \
182 -incube="NONE" \  
183 -outcube="expcube1.fits" \  
184 -caldb="irf" \  
185 -irf="cta_dummy_irf" \  
186 -emin=0.1 \  
187 -emax=100.0 \  
188 -enumbins=20 \  
189 -nxpix=200 \  
190 -nypix=200 \  
191 -binsz=0.02 \  
192 -coordsys="CEL" \  
193 -xref=83.63 \  
194 -yref=22.01 \  
195 -proj="CAR" 182 + incube="NONE" \
  183 + outcube="expcube1.fits" \
  184 + caldb="irf" \
  185 + irf="cta_dummy_irf" \
  186 + emin=0.1 \
  187 + emax=100.0 \
  188 + enumbins=20 \
  189 + nxpix=200 \
  190 + nypix=200 \
  191 + binsz=0.02 \
  192 + coordsys="CEL" \
  193 + xref=83.63 \
  194 + yref=22.01 \
  195 + proj="CAR"
196 $ECHO -n "." 196 $ECHO -n "."
197 if [ -s "expcube1.fits" ] 197 if [ -s "expcube1.fits" ]
198 then 198 then
@@ -204,10 +204,10 @@ fi @@ -204,10 +204,10 @@ fi
204 # 204 #
205 # Run 2 205 # Run 2
206 $ctexpcube inobs="events.fits" \ 206 $ctexpcube inobs="events.fits" \
207 -incube="cntmap2.fits" \  
208 -outcube="expcube2.fits" \  
209 -caldb="irf" \  
210 -irf="cta_dummy_irf" 207 + incube="cntmap2.fits" \
  208 + outcube="expcube2.fits" \
  209 + caldb="irf" \
  210 + irf="cta_dummy_irf"
211 $ECHO -n "." 211 $ECHO -n "."
212 if [ -s "expcube2.fits" ] 212 if [ -s "expcube2.fits" ]
213 then 213 then
@@ -226,22 +226,22 @@ $ECHO -n &quot;Test ctpsfcube: &quot; @@ -226,22 +226,22 @@ $ECHO -n &quot;Test ctpsfcube: &quot;
226 # 226 #
227 # Run 1 227 # Run 1
228 $ctpsfcube inobs="data/crab_events.fits" \ 228 $ctpsfcube inobs="data/crab_events.fits" \
229 -incube="NONE" \  
230 -outcube="psfcube1.fits" \  
231 -caldb="irf" \  
232 -irf="cta_dummy_irf" \  
233 -emin=0.1 \  
234 -emax=100.0 \  
235 -enumbins=20 \  
236 -nxpix=10 \  
237 -nypix=10 \  
238 -binsz=0.4 \  
239 -coordsys="CEL" \  
240 -xref=83.63 \  
241 -yref=22.01 \  
242 -proj="CAR" \  
243 -amax=0.3 \  
244 -anumbins=10 229 + incube="NONE" \
  230 + outcube="psfcube1.fits" \
  231 + caldb="irf" \
  232 + irf="cta_dummy_irf" \
  233 + emin=0.1 \
  234 + emax=100.0 \
  235 + enumbins=20 \
  236 + nxpix=10 \
  237 + nypix=10 \
  238 + binsz=0.4 \
  239 + coordsys="CEL" \
  240 + xref=83.63 \
  241 + yref=22.01 \
  242 + proj="CAR" \
  243 + amax=0.3 \
  244 + anumbins=10
245 $ECHO -n "." 245 $ECHO -n "."
246 if [ -s "psfcube1.fits" ] 246 if [ -s "psfcube1.fits" ]
247 then 247 then
@@ -253,12 +253,12 @@ fi @@ -253,12 +253,12 @@ fi
253 # 253 #
254 # Run 2 254 # Run 2
255 $ctpsfcube inobs="events.fits" \ 255 $ctpsfcube inobs="events.fits" \
256 -incube="cntmap2.fits" \  
257 -outcube="psfcube2.fits" \  
258 -caldb="irf" \  
259 -irf="cta_dummy_irf" \  
260 -amax=0.3 \  
261 -anumbins=2 256 + incube="cntmap2.fits" \
  257 + outcube="psfcube2.fits" \
  258 + caldb="irf" \
  259 + irf="cta_dummy_irf" \
  260 + amax=0.3 \
  261 + anumbins=2
262 $ECHO -n "." 262 $ECHO -n "."
263 if [ -s "psfcube2.fits" ] 263 if [ -s "psfcube2.fits" ]
264 then 264 then
@@ -277,22 +277,22 @@ $ECHO -n &quot;Test ctbkgcube: &quot; @@ -277,22 +277,22 @@ $ECHO -n &quot;Test ctbkgcube: &quot;
277 # 277 #
278 # Run 1 278 # Run 1
279 $ctbkgcube inobs="data/crab_events.fits" \ 279 $ctbkgcube inobs="data/crab_events.fits" \
280 -inmodel="data/crab.xml" \  
281 -incube="NONE" \  
282 -outcube="bkgcube1.fits" \  
283 -outmodel="bkgcube1.xml" \  
284 -caldb="irf" \  
285 -irf="cta_dummy_irf" \  
286 -emin=0.1 \  
287 -emax=100.0 \  
288 -enumbins=20 \  
289 -nxpix=10 \  
290 -nypix=10 \  
291 -binsz=0.4 \  
292 -coordsys="CEL" \  
293 -xref=83.63 \  
294 -yref=22.01 \  
295 -proj="CAR" 280 + inmodel="data/crab.xml" \
  281 + incube="NONE" \
  282 + outcube="bkgcube1.fits" \
  283 + outmodel="bkgcube1.xml" \
  284 + caldb="irf" \
  285 + irf="cta_dummy_irf" \
  286 + emin=0.1 \
  287 + emax=100.0 \
  288 + enumbins=20 \
  289 + nxpix=10 \
  290 + nypix=10 \
  291 + binsz=0.4 \
  292 + coordsys="CEL" \
  293 + xref=83.63 \
  294 + yref=22.01 \
  295 + proj="CAR"
296 $ECHO -n "." 296 $ECHO -n "."
297 if [ -s "bkgcube1.fits" ] 297 if [ -s "bkgcube1.fits" ]
298 then 298 then
@@ -311,12 +311,12 @@ fi @@ -311,12 +311,12 @@ fi
311 # 311 #
312 # Run 2 312 # Run 2
313 $ctbkgcube inobs="events.fits" \ 313 $ctbkgcube inobs="events.fits" \
314 -inmodel="data/crab.xml" \  
315 -incube="cntmap2.fits" \  
316 -caldb="irf" \  
317 -irf="cta_dummy_irf" \  
318 -outcube="bkgcube2.fits" \  
319 -outmodel="bkgcube2.xml" 314 + inmodel="data/crab.xml" \
  315 + incube="cntmap2.fits" \
  316 + caldb="irf" \
  317 + irf="cta_dummy_irf" \
  318 + outcube="bkgcube2.fits" \
  319 + outmodel="bkgcube2.xml"
320 $ECHO -n "." 320 $ECHO -n "."
321 if [ -s "bkgcube2.fits" ] 321 if [ -s "bkgcube2.fits" ]
322 then 322 then
@@ -347,7 +347,7 @@ $ctmodel inobs=&quot;cntmap1.fits&quot; \ @@ -347,7 +347,7 @@ $ctmodel inobs=&quot;cntmap1.fits&quot; \
347 outcube="modmap1.fits" \ 347 outcube="modmap1.fits" \
348 expcube="NONE" \ 348 expcube="NONE" \
349 psfcube="NONE" \ 349 psfcube="NONE" \
350 - bgcube="NONE" \ 350 + bkgcube="NONE" \
351 caldb="irf" \ 351 caldb="irf" \
352 irf="cta_dummy_irf" \ 352 irf="cta_dummy_irf" \
353 inmodel="data/crab.xml" 353 inmodel="data/crab.xml"
@@ -362,14 +362,14 @@ fi @@ -362,14 +362,14 @@ fi
362 # 362 #
363 # Run 2 363 # Run 2
364 $ctmodel inobs="cntmap2.fits" \ 364 $ctmodel inobs="cntmap2.fits" \
365 -incube="cntmap2.fits" \  
366 -outcube="modmap2.fits" \  
367 -expcube="expcube2.fits" \  
368 -psfcube="psfcube2.fits" \  
369 -bgcube="bkgcube2.fits" \  
370 -caldb="irf" \  
371 -irf="cta_dummy_irf" \  
372 -inmodel="bkgcube2.xml" 365 + incube="cntmap2.fits" \
  366 + outcube="modmap2.fits" \
  367 + expcube="expcube2.fits" \
  368 + psfcube="psfcube2.fits" \
  369 + bkgcube="bkgcube2.fits" \
  370 + caldb="irf" \
  371 + irf="cta_dummy_irf" \
  372 + inmodel="bkgcube2.xml"
373 $ECHO -n "." 373 $ECHO -n "."
374 if [ -s "modmap2.fits" ] 374 if [ -s "modmap2.fits" ]
375 then 375 then
@@ -448,7 +448,7 @@ $ctlike inobs=&quot;cntmap1.fits&quot; \ @@ -448,7 +448,7 @@ $ctlike inobs=&quot;cntmap1.fits&quot; \
448 outmodel="results_binned.xml" \ 448 outmodel="results_binned.xml" \
449 expcube="NONE" \ 449 expcube="NONE" \
450 psfcube="NONE" \ 450 psfcube="NONE" \
451 - bgcube="NONE" \ 451 + bkgcube="NONE" \
452 caldb="irf" \ 452 caldb="irf" \
453 irf="cta_dummy_irf" 453 irf="cta_dummy_irf"
454 $ECHO -n "." 454 $ECHO -n "."
@@ -482,7 +482,7 @@ $ctlike inobs=&quot;cntmap2.fits&quot; \ @@ -482,7 +482,7 @@ $ctlike inobs=&quot;cntmap2.fits&quot; \
482 outmodel="results_binned_cube_background.xml" \ 482 outmodel="results_binned_cube_background.xml" \
483 expcube="expcube2.fits" \ 483 expcube="expcube2.fits" \
484 psfcube="psfcube2.fits" \ 484 psfcube="psfcube2.fits" \
485 - bgcube="bkgcube2.fits" \ 485 + bkgcube="bkgcube2.fits" \
486 caldb="irf" \ 486 caldb="irf" \
487 irf="cta_dummy_irf" 487 irf="cta_dummy_irf"
488 $ECHO -n "." 488 $ECHO -n "."
@@ -577,10 +577,10 @@ $ECHO &quot; ok&quot; @@ -577,10 +577,10 @@ $ECHO &quot; ok&quot;
577 # =============== 577 # ===============
578 $ECHO -n "Test ctulimit: " 578 $ECHO -n "Test ctulimit: "
579 $ctulimit inobs="data/crab_events.fits.gz" \ 579 $ctulimit inobs="data/crab_events.fits.gz" \
580 - inmodel="data/crab.xml" \  
581 - srcname="Crab" \  
582 - caldb="irf" \  
583 - irf="cta_dummy_irf" 580 + inmodel="data/crab.xml" \
  581 + srcname="Crab" \
  582 + caldb="irf" \
  583 + irf="cta_dummy_irf"
584 $ECHO -n "." 584 $ECHO -n "."
585 if [ -s "ctulimit.log" ] 585 if [ -s "ctulimit.log" ]
586 then 586 then