From e6d349bfe027c02228dd04f63f1cff22bf2eddfc Mon Sep 17 00:00:00 2001 From: Jürgen Knödlseder Date: Tue, 28 Mar 2017 11:37:01 +0200 Subject: [PATCH] Support of map cubes lying outside RoI in ctobssim --- ChangeLog | 3 ++- NEWS | 7 +++++-- src/ctobssim/ctobssim.cpp | 53 ++++++++++++++++++++++++++++++++++++++++++++++++----- 3 files changed, 55 insertions(+), 8 deletions(-) diff --git a/ChangeLog b/ChangeLog index abb74f4..1f435ad 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,8 +1,9 @@ -2017-03-24 +2017-03-28 * Version 1.3.0 released ======================== + Support of map cubes lying outside RoI in ctobssim Support energy boundary FITS file in ctmapcube (#1978) Add csebins script (#1976) Consider RoI in ctskymap and use numerical integration for IRF (#1975) diff --git a/NEWS b/NEWS index be31e60..0ac2b8c 100644 --- a/NEWS +++ b/NEWS @@ -1,6 +1,6 @@ New Features and Important Changes in ctools 1.3.0.dev1 -24 March 2017 +28 March 2017 Introduction @@ -27,7 +27,10 @@ so that writing commas instead of dots is avoided (#1370). ctobssim - CTA observation simulation ------------------------------------- -None +The tool now supports map cubes that are not overlapping with the RoI. +In this case, no spectral nodes exist and before the spectral flux() +method called in ctobssim::get_model_flux() was throwing an exception. +Now, for map cubes outside the RoI the flux is set to zero. ctselect - CTA event selection diff --git a/src/ctobssim/ctobssim.cpp b/src/ctobssim/ctobssim.cpp index 108eb00..d6bd2b0 100644 --- a/src/ctobssim/ctobssim.cpp +++ b/src/ctobssim/ctobssim.cpp @@ -722,6 +722,11 @@ void ctobssim::simulate_source(GCTAObservation* obs, GRan& ran, GLog* wrklog) { + // Debug code: signal that we step into the simulate_source method + #if defined(G_SOURCE_DEBUG) + std::cout << "ctobssim::simulate_source: in" << std::endl; + #endif + // Continue only if observation pointer is valid if (obs != NULL) { @@ -898,6 +903,11 @@ void ctobssim::simulate_source(GCTAObservation* obs, } // endif: observation pointer was valid + // Debug code: signal that we are down with the simulate_source method + #if defined(G_SOURCE_DEBUG) + std::cout << "ctobssim::simulate_source: out" << std::endl; + #endif + // Return return; } @@ -946,6 +956,14 @@ void ctobssim::simulate_interval(GCTAObservation* obs, std::vector& nphotons, std::vector& nevents) { + // Debug code: signal that we step into the simulate_interval method + #if defined(G_SOURCE_DEBUG) + std::cout << "ctobssim::simulate_interval: in"; + std::cout << " Etrue=[" << etrue_min << "," << etrue_max << "]"; + std::cout << " Ereco=[" << ereco_min << "," << ereco_max << "]"; + std::cout << std::endl; + #endif + // Loop over all models for (int i = 0; i < models.size(); ++i) { @@ -1084,6 +1102,11 @@ void ctobssim::simulate_interval(GCTAObservation* obs, } // endfor: looped over models + // Debug code: signal that we are down with the simulate_interval method + #if defined(G_SOURCE_DEBUG) + std::cout << "ctobssim::simulate_interval: out" << std::endl; + #endif + // Return return; } @@ -1289,6 +1312,13 @@ double ctobssim::get_model_flux(const GModelSky* model, const int& indent, GLog* wrklog) { + // Debug code: signal that we step into the get_model_flux method + #if defined(G_SOURCE_DEBUG) + std::cout << "ctobssim::get_model_flux: in"; + std::cout << " Etrue=[" << emin << "," << emax << "]"; + std::cout << std::endl; + #endif + // Initialise flux double flux = 0.0; @@ -1329,19 +1359,20 @@ double ctobssim::get_model_flux(const GModelSky* model, // function spectral model that is the product of the diffuse // cube node function and the spectral model evaluated at the // energies of the node function - GModelSpatialDiffuseCube* cube = dynamic_cast(model->spatial()); + GModelSpatialDiffuseCube* cube = + dynamic_cast(model->spatial()); if (cube != NULL) { // Set MC cone cube->set_mc_cone(centre, radius); // Allocate node function to replace the spectral component - GModelSpectralNodes* nodes = new GModelSpectralNodes(cube->spectrum()); + GModelSpectralNodes* nodes = + new GModelSpectralNodes(cube->spectrum()); for (int i = 0; i < nodes->nodes(); ++i) { GEnergy energy = nodes->energy(i); - GTime time; // Dummy time double intensity = nodes->intensity(i); - double value = spectral->eval(energy, time); + double value = spectral->eval(energy); nodes->intensity(i, value * intensity); } @@ -1351,11 +1382,18 @@ double ctobssim::get_model_flux(const GModelSky* model, // Set the spectral model pointer to the node function spectral = nodes; + // Kluge: if there are no nodes then the spectral->flux method + // will throw an exception. We therefore set here the use_model + // flag to false in case that there are no spectral nodes + if (nodes->nodes() == 0) { + use_model = false; + } + } // endif: spatial model was a diffuse cube // Compute flux within [emin, emax] in model from spectral // component (units: ph/cm2/s) - double flux0 = spectral->flux(emin, emax); + double flux0 = (use_model) ? spectral->flux(emin, emax) : 0.0; flux = flux0 * norm; // Dump flux @@ -1379,6 +1417,11 @@ double ctobssim::get_model_flux(const GModelSky* model, } // endif: model overlaps with simulation region + // Debug code: signal that we step out of the get_model_flux method + #if defined(G_SOURCE_DEBUG) + std::cout << "ctobssim::get_model_flux: in" << std::endl; + #endif + // Return model flux return flux; } -- libgit2 0.21.4