From cc97b1079697578e26dd99eff5c22de4e88c3f0e Mon Sep 17 00:00:00 2001 From: Jürgen Knödlseder Date: Fri, 14 Apr 2017 15:06:30 +0200 Subject: [PATCH] ctobssim adds Monte Carlo identifiers to simulated events (#1483) --- ChangeLog | 3 ++- NEWS | 6 +++++- src/ctobssim/ctobssim.cpp | 80 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++- src/ctobssim/ctobssim.hpp | 4 ++++ 4 files changed, 90 insertions(+), 3 deletions(-) diff --git a/ChangeLog b/ChangeLog index b66b010..281d036 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,8 +1,9 @@ -2017-04-13 +2017-04-14 * Version 1.3.0 released ======================== + ctobssim adds Monte Carlo identifiers to simulated events (#1483) Add "phase" parameter to ctselect (#1996) Add ctphase tool (#1982) Add ctprob tool (#1991) diff --git a/NEWS b/NEWS index 792ae2e..07b32d3 100644 --- a/NEWS +++ b/NEWS @@ -1,6 +1,6 @@ New Features and Important Changes in ctools 1.3.0.dev1 -13 April 2017 +14 April 2017 Introduction @@ -44,6 +44,10 @@ 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. +The tool now also added Monte Carlo identifiers to the simulated events +so that each event can be traced back to the model from which it originated +(#1483). + ctselect - CTA event selection ------------------------------ diff --git a/src/ctobssim/ctobssim.cpp b/src/ctobssim/ctobssim.cpp index d66c8c3..7a3fc41 100644 --- a/src/ctobssim/ctobssim.cpp +++ b/src/ctobssim/ctobssim.cpp @@ -363,6 +363,9 @@ void ctobssim::run(void) // Simulate source events simulate_background(&obs_clone, models, m_rans[i], &wrklog); + // Set Monte Carlo identifier and name correspondance + set_mc_id_names(&obs_clone, models); + // Dump simulation results if (logTerse()) { wrklog << gammalib::parformat("MC events"); @@ -1063,7 +1066,7 @@ void ctobssim::simulate_interval(GCTAObservation* obs, } // Simulate time slice - simulate_time_slice(obs, rsp, events, model, tstart, tstop, + simulate_time_slice(obs, rsp, events, model, i+1, tstart, tstop, etrue_min, etrue_max, ereco_min, ereco_max, dir, rad, area, ran, wrklog, indent, @@ -1119,6 +1122,7 @@ void ctobssim::simulate_interval(GCTAObservation* obs, * @param[in] rsp Pointer on CTA response. * @param[in,out] events Pointer on CTA event list. * @param[in] model Model. + * @param[in] mc_id Monte Carlo identifier of model. * @param[in] tstart Start time. * @param[in] tstop Stop time. * @param[in] etrue_min Minimum true energy. @@ -1140,6 +1144,7 @@ void ctobssim::simulate_time_slice(GCTAObservation* obs, const GCTAResponseIrf* rsp, GCTAEventList* events, const GModelSky* model, + const int& mc_id, const GTime& tstart, const GTime& tstop, const GEnergy& etrue_min, @@ -1207,8 +1212,17 @@ void ctobssim::simulate_time_slice(GCTAObservation* obs, event->energy() <= ereco_max && event->time() >= tstart && event->time() <= tstop) { + + // Set event identifier event->event_id(m_event_id); + + // Set Monte Carlo identifier + event->mc_id(mc_id); + + // Append event events->append(*event); + + // Increment counters m_event_id++; nevents++; } @@ -1217,6 +1231,9 @@ void ctobssim::simulate_time_slice(GCTAObservation* obs, } // endfor: looped over events + // Signal that event list contains Monte Carlo identifiers + events->has_mc_id(true); + // Return return; } @@ -1509,6 +1526,9 @@ void ctobssim::simulate_background(GCTAObservation* obs, event->event_id(m_event_id); m_event_id++; + // Set Monte Carlo identifier + event->mc_id(i+1); + // Append event events->append(*event); @@ -1541,6 +1561,64 @@ void ctobssim::simulate_background(GCTAObservation* obs, } // endfor: looped over all models + // Signal that event list contains Monte Carlo identifiers + events->has_mc_id(true); + + } // endif: observation pointer was valid + + // Return + return; +} + + +/***********************************************************************//** + * @brief Set correspondance between Monte Carlo identifier and model names + * + * @param[in] obs Pointer on CTA observation. + * @param[in] models Models. + * @param[in] wrklog Pointer to logger. + * + * Sets the correspondance between Monte Carlo identifier and model names. + ***************************************************************************/ +void ctobssim::set_mc_id_names(GCTAObservation* obs, + const GModels& models, + GLog* wrklog) +{ + // Continue only if observation pointer is valid + if (obs != NULL) { + + // If no logger is specified then use the default logger + if (wrklog == NULL) { + wrklog = &log; + } + + // Get pointer on event list (circumvent const correctness) + GCTAEventList* events = + static_cast(const_cast(obs->events())); + + // Initialise lists + std::vector ids; + std::vector names; + + // Loop over all models + for (int i = 0; i < models.size(); ++i) { + + // Append identifier and names + ids.push_back(i+1); + names.push_back(models[i]->name()); + + // Dump identifier name correspondance + if (logNormal()) { + *wrklog << gammalib::parformat("MC identifier " + + gammalib::str(i+1)); + *wrklog << models[i]->name() << std::endl; + } + + } // endfor: looped over all models + + // Set correspondance for the event list + events->set_mc_id_names(ids, names); + } // endif: observation pointer was valid // Return diff --git a/src/ctobssim/ctobssim.hpp b/src/ctobssim/ctobssim.hpp index 9982549..724f205 100644 --- a/src/ctobssim/ctobssim.hpp +++ b/src/ctobssim/ctobssim.hpp @@ -100,6 +100,7 @@ protected: const GCTAResponseIrf* rsp, GCTAEventList* events, const GModelSky* model, + const int& mc_id, const GTime& tstart, const GTime& tstop, const GEnergy& etrue_min, @@ -129,6 +130,9 @@ protected: const GModels& models, GRan& ran, GLog* wrklog = NULL); + void set_mc_id_names(GCTAObservation* obs, + const GModels& models, + GLog* wrklog = NULL); void save_fits(void); void save_xml(void); std::string outfile(const int& index); -- libgit2 0.21.4