Samuel Flis / gammalib

Commit 19f117ae2dc878b856b8150da4368d1049938bdc

Authored by Jürgen Knödlseder
1 parent 34b8b420

Make sure that energy boundaries are updated in GCTAEdisp2D

Inline Side-by-side
Showing 1 changed file with 13 additions and 6 deletions
inst/cta/src/GCTAEdisp2D.cpp
  View file @19f117a
... ... @@ -447,6 +447,7 @@ void GCTAEdisp2D::save(const GFilename& filename, const bool& clobber) const
447 447 * @param[in] phi Azimuth angle in camera system (rad).
448 448 * @param[in] zenith Zenith angle in Earth system (rad).
449 449 * @param[in] azimuth Azimuth angle in Earth system (rad).
  450 + * @return Observed energy.
450 451 *
451 452 * @exception GException::invalid_return_value
452 453 * No energy dispersion information found for parameters or
... ... @@ -475,11 +476,12 @@ GEnergy GCTAEdisp2D::mc(GRan& ran,
475 476 // energy (this means that no energy dispersion information is available
476 477 // for that energy)
477 478 if (emin >= emax) {
  479 + double eng = std::pow(10.0, logEsrc);
478 480 std::string msg = "No valid energy dispersion information available "
479   - "for true photon energy "+ebounds.emin().print()+
480   - ", offset angle "+
  481 + "for true photon energy "+gammalib::str(eng)+" TeV,"
  482 + "offset angle "+
481 483 gammalib::str(theta*gammalib::rad2deg)+" deg "
482   - " and azimuth angle "+
  484 + "and azimuth angle "+
483 485 gammalib::str(phi*gammalib::rad2deg)+" deg. Either "
484 486 "provide an energy dispersion matrix covering these "
485 487 "parameters or restrict the parameter space.";
... ... @@ -516,9 +518,10 @@ GEnergy GCTAEdisp2D::mc(GRan& ran,
516 518 if (f == 0.0) {
517 519 zeros++;
518 520 if (zeros > max_subsequent_zeros) {
  521 + double eng = std::pow(10.0, logEsrc);
519 522 std::string msg = "No valid energy dispersion information "
520 523 "available for true photon energy "+
521   - ebounds.emin().print()+", offset angle "+
  524 + gammalib::str(eng)+" TeV, offset angle "+
522 525 gammalib::str(theta*gammalib::rad2deg)+
523 526 " deg and azimuth angle "+
524 527 gammalib::str(phi*gammalib::rad2deg)+
... ... @@ -534,7 +537,7 @@ GEnergy GCTAEdisp2D::mc(GRan& ran,
534 537 // Get uniform random value between zero and the maximum energy
535 538 // dispersion value. The loop is quite if the random number is smaller
536 539 // than the energy dispersion value
537   - ftest = ran.uniform() * m_max_edisp;
  540 + ftest = ran.uniform() * m_max_edisp;
538 541  
539 542 } // endwhile:
540 543  
... ... @@ -555,6 +558,7 @@ GEnergy GCTAEdisp2D::mc(GRan& ran,
555 558 * @param[in] phi Azimuth angle in camera system (rad). Not used.
556 559 * @param[in] zenith Zenith angle in Earth system (rad). Not used.
557 560 * @param[in] azimuth Azimuth angle in Earth system (rad). Not used.
  561 + * @return Observed energy boundaries.
558 562 *
559 563 * Returns the band of observed energies outside of which the energy
560 564 * dispersion becomes negligible for a given true energy @p logEsrc and
... ... @@ -576,6 +580,7 @@ GEbounds GCTAEdisp2D::ebounds_obs(const double& logEsrc,
576 580 // Set computation flag
577 581 m_ebounds_obs_computed = true;
578 582 m_last_theta_obs = theta;
  583 + m_last_logEsrc = -30.0; // force update
579 584  
580 585 // Compute ebounds_obs
581 586 compute_ebounds_obs(theta, phi, zenith, azimuth);
... ... @@ -620,6 +625,7 @@ GEbounds GCTAEdisp2D::ebounds_obs(const double& logEsrc,
620 625 * @param[in] phi Azimuth angle in camera system (rad). Not used.
621 626 * @param[in] zenith Zenith angle in Earth system (rad). Not used.
622 627 * @param[in] azimuth Azimuth angle in Earth system (rad). Not used.
  628 + * @return True energy boundaries.
623 629 *
624 630 * Returns the band of true photon energies outside of which the energy
625 631 * dispersion becomes negligible for a given observed energy @p logEobs and
... ... @@ -635,12 +641,13 @@ GEbounds GCTAEdisp2D::ebounds_src(const double& logEobs,
635 641 // Make sure that energy dispersion is online
636 642 fetch();
637 643  
638   - // Compute only if parameters changed
  644 + // Compute true energy boundaries only if parameters changed
639 645 if (!m_ebounds_src_computed || theta != m_last_theta_src) {
640 646  
641 647 // Set computation flag
642 648 m_ebounds_src_computed = true;
643 649 m_last_theta_src = theta;
  650 + m_last_logEobs = -30.0; // force update
644 651  
645 652 // Compute ebounds_src
646 653 compute_ebounds_src(theta, phi, zenith, azimuth);
... ...