Commit 6f41050678d12b9c8572d4f94abddf32bff3329a

Authored by L. Tibaldo
1 parent 9cc4a742

logging and fixed bug

cscripts/csresspec.par
... ... @@ -61,4 +61,4 @@ chatter, i, h, 2,0,4, "Chattiness of output"
61 61 clobber, b, h, yes,,, "Overwrite existing output files with new output files?"
62 62 debug, b, h, no,,, "Debugging mode activated"
63 63 mode, s, h, ql,,, "Mode of automatic parameters"
64   -logfile, f, h, csresmap.log,,, "Log filename"
  64 +logfile, f, h, csresspec.log,,, "Log filename"
... ...
cscripts/csresspec.py 100644 → 100755
... ... @@ -79,10 +79,11 @@ class csresspec(ctools.csobservation):
79 79 n_other = self.obs().size() - n_cta
80 80  
81 81 # Log census of input observations
82   - self._log_value(gammalib.TERSE, 'Unbinned CTA observations', n_unbinned)
83   - self._log_value(gammalib.TERSE, 'Binned CTA observations', n_binned)
84   - self._log_value(gammalib.TERSE, 'On/off CTA observations', n_onoff)
85   - self._log_value(gammalib.TERSE, 'Other observations', n_other)
  82 + self._log_value(gammalib.NORMAL, 'Unbinned CTA observations',
  83 + n_unbinned)
  84 + self._log_value(gammalib.NORMAL, 'Binned CTA observations', n_binned)
  85 + self._log_value(gammalib.NORMAL, 'On/off CTA observations', n_onoff)
  86 + self._log_value(gammalib.NORMAL, 'Other observations', n_other)
86 87 if n_other > 0:
87 88 msg = 'WARNING: Only CTA observation can be handled, all non-CTA ' \
88 89 + 'observations will be ignored.\n'
... ... @@ -182,13 +183,16 @@ class csresspec(ctools.csobservation):
182 183 # If On/Off use the GCTAOnOffObservation constructor
183 184 if self.obs()[0].classname() == 'GCTAOnOffObservation':
184 185 msg = 'Stacking %s On/Off observations.' % (self.obs().size())
185   - self._log_string(gammalib.TERSE, msg)
186   - stacked_obs = gammalib.GCTAOnOffObservation(self.obs())
  186 + self._log_string(gammalib.NORMAL, msg)
  187 + stacked = gammalib.GCTAOnOffObservation(self.obs())
  188 + stacked_obs = gammalib.GObservations()
  189 + stacked_obs.append(stacked)
  190 + stacked_obs.models(self.obs().models())
187 191  
188 192 # If event list then bin observations
189 193 elif self.obs()[0].classname() == 'GCTAObservation':
190 194 msg = 'Stacking %s event lists.' % (self.obs().size())
191   - self._log_string(gammalib.TERSE, msg)
  195 + self._log_string(gammalib.NORMAL, msg)
192 196 stacked_obs = obsutils.get_stacked_obs(self, self.obs())
193 197  
194 198 # Return
... ... @@ -212,6 +216,9 @@ class csresspec(ctools.csobservation):
212 216 npix = int(2 * rad / 0.02) + 1
213 217  
214 218 # Bin events
  219 + msg = 'Binning events'
  220 + self._log_string(gammalib.EXPLICIT, msg)
  221 + #
215 222 cntcube = ctools.ctbin(obs)
216 223 cntcube['xref'] = ra
217 224 cntcube['yref'] = dec
... ... @@ -233,6 +240,8 @@ class csresspec(ctools.csobservation):
233 240 binned_obs = cntcube.obs().copy()
234 241  
235 242 # Compute binned response
  243 + msg = 'Computing binned response'
  244 + self._log_string(gammalib.EXPLICIT, msg)
236 245 response = obsutils.get_stacked_response(obs, ra, dec,
237 246 binsz=0.02,
238 247 nxpix=npix,
... ... @@ -320,6 +329,8 @@ class csresspec(ctools.csobservation):
320 329 # If we started from event list mask the ROI only
321 330 # for residual computation
322 331 if evlist_info['was_list']:
  332 + msg = 'Masking ROI from original event list'
  333 + self._log_string(gammalib.EXPLICIT, msg)
323 334 cube = self._masked_cube(cube, evlist_info['roi_ra'],
324 335 evlist_info['roi_dec'],
325 336 evlist_info['roi_rad'],
... ... @@ -328,6 +339,8 @@ class csresspec(ctools.csobservation):
328 339  
329 340 # Apply user mask
330 341 if self._mask:
  342 + msg = 'Masking ROI requested by user'
  343 + self._log_string(gammalib.EXPLICIT, msg)
331 344 cube = self._masked_cube(cube, self['ra'], self['dec'],
332 345 self['rad'],
333 346 regfile=self['regfile'])
... ... @@ -437,8 +450,10 @@ class csresspec(ctools.csobservation):
437 450  
438 451 # Stack observations if requested
439 452 if self._stack:
  453 + self._log_header1(gammalib.TERSE, 'Stacking Observations')
440 454 self.obs(self._stack_observations())
441 455  
  456 + self._log_header1(gammalib.TERSE, 'Processing Observations')
442 457 # Loop over observations and calculate residuals
443 458 for s, observation in enumerate(self.obs()):
444 459  
... ... @@ -449,6 +464,10 @@ class csresspec(ctools.csobservation):
449 464 if obs_id == '' and self.obs().size() > 1:
450 465 obs_id = str(s)
451 466  
  467 + if self.obs().size() > 1:
  468 + self._log_header2(gammalib.NORMAL,
  469 + 'Processing observation %s' % obs_id)
  470 +
452 471 # Turn into observation container and assign models
453 472 obs = gammalib.GObservations()
454 473 obs.append(observation)
... ... @@ -467,6 +486,8 @@ class csresspec(ctools.csobservation):
467 486 else:
468 487 # we remember if we binned an event list
469 488 # so that we can mask only the ROI for residual calculation
  489 + msg = 'Setting up binned observation'
  490 + self._log_string(gammalib.NORMAL, msg)
470 491 obs, evlist_info = self._bin_evlist(obs)
471 492  
472 493 ## Calculate Model and residuals
... ... @@ -475,6 +496,8 @@ class csresspec(ctools.csobservation):
475 496 modcube = gammalib.GCTAEventCube(self['inmodel'].filename())
476 497 # Otherwise calculate it now
477 498 else:
  499 + msg = 'Computing model cube'
  500 + self._log_string(gammalib.NORMAL, msg)
478 501 modelcube = ctools.ctmodel(obs)
479 502 modelcube.run()
480 503 modcube = modelcube.cube().copy()
... ... @@ -483,6 +506,8 @@ class csresspec(ctools.csobservation):
483 506 cntcube = obs[0].events().copy()
484 507  
485 508 # Derive count spectra from cubes
  509 + msg = 'Computing counts, model, and residual spectra'
  510 + self._log_string(gammalib.NORMAL, msg)
486 511 counts = self._cube_to_spectrum(cntcube, evlist_info)
487 512 model = self._cube_to_spectrum(modcube, evlist_info)
488 513  
... ... @@ -493,12 +518,17 @@ class csresspec(ctools.csobservation):
493 518 ebounds = cntcube.ebounds()
494 519  
495 520 # Fill results table
  521 + msg = 'Filling residual table'
  522 + self._log_string(gammalib.NORMAL, msg)
496 523 table = self._residuals_table(obs_id, ebounds, counts, model,
497 524 residuals)
498 525  
499 526 ## Calculate models of individual components if requested
500 527 if self['components'].boolean():
501 528 for component in obs.models():
  529 + self._log_value(gammalib.NORMAL,
  530 + 'Computing model for component',
  531 + component.name())
502 532 # Set model cube models to
503 533 # individual component
504 534 model_cont = gammalib.GModels()
... ... @@ -522,6 +552,8 @@ class csresspec(ctools.csobservation):
522 552 onoff = obs[0]
523 553  
524 554 ## Calculate Counts, Model and residuals
  555 + msg = 'Computing counts, model, and residual spectra'
  556 + self._log_string(gammalib.NORMAL, msg)
525 557  
526 558 # On spectrum
527 559 counts = onoff.on_spec().counts_spectrum()
... ... @@ -541,10 +573,14 @@ class csresspec(ctools.csobservation):
541 573 ebounds = onoff.on_spec().ebounds()
542 574  
543 575 # Fill results table
  576 + msg = 'Filling residual table'
  577 + self._log_string(gammalib.NORMAL, msg)
544 578 table = self._residuals_table(obs_id, ebounds, counts, model,
545 579 residuals)
546 580  
547 581 # Get Off spectrum and add to table
  582 + msg = 'Computing counts, model, and residual spectra for Off regions'
  583 + self._log_string(gammalib.NORMAL, msg)
548 584 counts_off = onoff.off_spec().counts_spectrum()
549 585 table = self._append_column(table, 'Counts_Off',
550 586 counts_off)
... ... @@ -561,13 +597,17 @@ class csresspec(ctools.csobservation):
561 597 ## Calculate models of individual components if requested
562 598 if self['components'].boolean():
563 599 for component in obs.models():
564   - # If background, append the values already calculated
  600 + # If background pass
  601 + # We always add the background at the end so that
  602 + # we accommodate WSTAT for which the background is not
  603 + # mandatory in the model
565 604 if component.classname() == 'GCTAModelIrfBackground':
566   - background *= onoff.on_spec().backscal_spectrum()
567   - table = self._append_column(table, 'Background',
568   - background)
  605 + pass
569 606 # Otherwise calculate gamma component and append to Table
570 607 else:
  608 + self._log_value(gammalib.NORMAL,
  609 + 'Computing model for component',
  610 + component.name())
571 611 # Create observation container for individual components
572 612 model_cont = gammalib.GModels()
573 613 model_cont.append(component)
... ... @@ -577,6 +617,13 @@ class csresspec(ctools.csobservation):
577 617 # Append to table
578 618 table = self._append_column(table, component.name(),
579 619 model)
  620 + # Add background already calculated
  621 + self._log_value(gammalib.NORMAL,
  622 + 'Computing model for component',
  623 + 'Background')
  624 + background *= onoff.on_spec().backscal_spectrum()
  625 + table = self._append_column(table, 'Background',
  626 + background)
580 627  
581 628 ## Append results table to output file
582 629 self._fits.append(table)
... ...