Commit b0375e86b171595a7aa7f38dab7b2ba9be18f65f

Authored by Jürgen Knödlseder
1 parent 4602a860

Add GEventList::has_mc_id() and GEventAtom::mc_id() methods (#1483)

ChangeLog
1   -2017-04-11
  1 +2017-04-13
2 2  
3 3 * Version 1.3.0 released
4 4 =========================
5 5  
  6 + Add GEventList::has_mc_id() and GEventAtom::mc_id() methods (#1483)
6 7 Add GEventList::has_phase() and GEventList::has_detxy() methods (#1982)
7 8 Add GEventList::append_column() method (#1991)
8 9 Add scalar addition and subtraction operators to matrix classes (#1990)
... ...
1 1 New Features and Important Changes in GammaLib 1.3.0.dev1
2 2  
3   -11 April 2017
  3 +13 April 2017
4 4  
5 5  
6 6 1. Introduction
... ... @@ -52,6 +52,10 @@ The following methods have been added:
52 52 - GCTAEventList::has_phase()
53 53 - GCTAEventList::has_detxy(bool&)
54 54 - GCTAEventList::has_detxy()
  55 +- GCTAEventList::has_mc_id(bool&)
  56 +- GCTAEventList::has_mc_id()
  57 +- GCTAEventAtom::mc_id(int&)
  58 +- GCTAEventAtom::mc_id()
55 59  
56 60 The following methods have been removed:
57 61 - ...
... ... @@ -225,6 +229,10 @@ event list (#1991).
225 229 Add GCTAEventList::has_phase() and GCTAEventList::has_detxy() getter and setter
226 230 methods to manipulate the event list information flags (#1982).
227 231  
  232 +Add GCTAEventList::has_mc_id() and GCTAEventAtom::mc_id() getter and setter methods
  233 +to set the Monte Carlo identifier of a given event. This allows to attribute a
  234 +Monte Carlo simulated event to a specific model (#1483).
  235 +
228 236  
229 237 21. COMPTEL interface
230 238 ---------------------
... ...
inst/cta/include/GCTAEventAtom.hpp
1 1 /***************************************************************************
2 2 * GCTAEventAtom.hpp - CTA event atom class *
3 3 * ----------------------------------------------------------------------- *
4   - * copyright (C) 2010-2015 by Juergen Knoedlseder *
  4 + * copyright (C) 2010-2017 by Juergen Knoedlseder *
5 5 * ----------------------------------------------------------------------- *
6 6 * *
7 7 * This program is free software: you can redistribute it and/or modify *
... ... @@ -72,8 +72,10 @@ public:
72 72 // Other methods
73 73 const int& index(void) const;
74 74 const unsigned long& event_id(void) const;
  75 + const int& mc_id(void) const;
75 76 const float& phase(void) const;
76 77 void event_id(const unsigned long& id);
  78 + void mc_id(const int& id);
77 79 void phase(const float& phase);
78 80  
79 81 protected:
... ... @@ -88,6 +90,7 @@ protected:
88 90 GEnergy m_energy; //!< Event energy
89 91 GTime m_time; //!< Event time
90 92 unsigned long m_event_id; //!< Event identifier
  93 + int m_mc_id; //!< Monte Carlo identifier
91 94 float m_phase; //!< Optional phase
92 95 };
93 96  
... ... @@ -206,6 +209,18 @@ const unsigned long&amp; GCTAEventAtom::event_id(void) const
206 209  
207 210  
208 211 /***********************************************************************//**
  212 + * @brief Return Monte Carlo identifier
  213 + *
  214 + * @return Monte Carlo identifier.
  215 + ***************************************************************************/
  216 +inline
  217 +const int& GCTAEventAtom::mc_id(void) const
  218 +{
  219 + return (m_mc_id);
  220 +}
  221 +
  222 +
  223 +/***********************************************************************//**
209 224 * @brief Return event phase
210 225 *
211 226 * @return Event phase.
... ... @@ -231,6 +246,19 @@ void GCTAEventAtom::event_id(const unsigned long&amp; id)
231 246  
232 247  
233 248 /***********************************************************************//**
  249 + * @brief Set Monte Carlo identifier
  250 + *
  251 + * @param[in] id Monte Carlo identifier.
  252 + ***************************************************************************/
  253 +inline
  254 +void GCTAEventAtom::mc_id(const int& id)
  255 +{
  256 + m_mc_id = id;
  257 + return;
  258 +}
  259 +
  260 +
  261 +/***********************************************************************//**
234 262 * @brief Set event phase
235 263 *
236 264 * @param[in] phase Event phase.
... ...
inst/cta/include/GCTAEventList.hpp
... ... @@ -99,8 +99,10 @@ public:
99 99 const std::string& gtiname(void) const;
100 100 void has_phase(const bool& has_phase);
101 101 void has_detxy(const bool& has_detxy);
  102 + void has_mc_id(const bool& has_mc_id);
102 103 const bool& has_phase() const;
103 104 const bool& has_detxy() const;
  105 + const bool& has_mc_id() const;
104 106  
105 107 protected:
106 108 // Protected methods
... ... @@ -122,6 +124,7 @@ protected:
122 124 std::string m_gti_extname; //!< GTI extension name
123 125 bool m_has_phase; //!< Signal presence of phase
124 126 bool m_has_detxy; //!< Signal presence of detector coordinates
  127 + bool m_has_mc_id; //!< Signal presence of MC identifier
125 128  
126 129 // Event list data
127 130 mutable std::vector<GCTAEventAtom> m_events; //!< Events
... ... @@ -239,6 +242,21 @@ void GCTAEventList::has_detxy(const bool&amp; has_detxy)
239 242 }
240 243  
241 244  
  245 +/***********************************************************************//**
  246 + * @brief Set the Monte Carlo identifier flag
  247 + *
  248 + * @param[in] has_mc_id Monte Carlo identifier flag.
  249 + *
  250 + * Set the Monte Carlo identifier flag.
  251 + ***************************************************************************/
  252 +inline
  253 +void GCTAEventList::has_mc_id(const bool& has_mc_id)
  254 +{
  255 + m_has_mc_id = has_mc_id;
  256 + return;
  257 +}
  258 +
  259 +
242 260 /*************************************************************************
243 261 * @brief Return phase information flag
244 262 *
... ... @@ -262,4 +280,16 @@ const bool&amp; GCTAEventList::has_detxy() const
262 280 return m_has_detxy;
263 281 }
264 282  
  283 +
  284 +/*************************************************************************
  285 + * @brief Return Monte Carlo identifier flag
  286 + *
  287 + * @return True if a Monte Carlo identifier is available, false otherwise.
  288 + ***************************************************************************/
  289 +inline
  290 +const bool& GCTAEventList::has_mc_id() const
  291 +{
  292 + return m_has_mc_id;
  293 +}
  294 +
265 295 #endif /* GCTAEVENTLIST_HPP */
... ...
inst/cta/pyext/GCTAEventAtom.i
1 1 /***************************************************************************
2 2 * GCTAEventAtom.i - CTA event atom class *
3 3 * ----------------------------------------------------------------------- *
4   - * copyright (C) 2011-2015 by Juergen Knoedlseder *
  4 + * copyright (C) 2011-2017 by Juergen Knoedlseder *
5 5 * ----------------------------------------------------------------------- *
6 6 * *
7 7 * This program is free software: you can redistribute it and/or modify *
... ... @@ -55,8 +55,10 @@ public:
55 55 // Other methods
56 56 const int& index(void) const;
57 57 const unsigned long& event_id(void) const;
  58 + const int& mc_id(void) const;
58 59 const float& phase(void) const;
59 60 void event_id(const unsigned long& id);
  61 + void mc_id(const int& id);
60 62 void phase(const float& phase);
61 63 };
62 64  
... ...
inst/cta/pyext/GCTAEventList.i
... ... @@ -72,8 +72,10 @@ public:
72 72 const std::string& gtiname(void) const;
73 73 void has_phase(const bool& has_phase);
74 74 void has_detxy(const bool& has_detxy);
  75 + void has_mc_id(const bool& has_mc_id);
75 76 const bool& has_phase() const;
76 77 const bool& has_detxy() const;
  78 + const bool& has_mc_id() const;
77 79 };
78 80  
79 81  
... ...
inst/cta/src/GCTAEventList.cpp
... ... @@ -38,6 +38,7 @@
38 38 #include "GFitsTableCol.hpp"
39 39 #include "GFitsTableFloatCol.hpp"
40 40 #include "GFitsTableDoubleCol.hpp"
  41 +#include "GFitsTableLongCol.hpp"
41 42 #include "GFitsTableULongCol.hpp"
42 43 #include "GTime.hpp"
43 44 #include "GTimeReference.hpp"
... ... @@ -382,6 +383,14 @@ void GCTAEventList::read(const GFits&amp; fits)
382 383 m_has_phase = false;
383 384 }
384 385  
  386 + // Signal presence of "MC_ID" column
  387 + if (events.contains("MC_ID")) {
  388 + m_has_mc_id = true;
  389 + }
  390 + else {
  391 + m_has_mc_id = false;
  392 + }
  393 +
385 394 // Read GTI extension name from data sub-space keyword
386 395 std::string gti_extname = gammalib::read_ds_gti_extname(events);
387 396  
... ... @@ -913,6 +922,7 @@ void GCTAEventList::init_members(void)
913 922 m_gti_extname = gammalib::extname_gti; //!< Default GTI extension name
914 923 m_has_phase = false;
915 924 m_has_detxy = true;
  925 + m_has_mc_id = false;
916 926  
917 927 // Initialise cache
918 928 m_irf_names.clear();
... ... @@ -938,6 +948,7 @@ void GCTAEventList::copy_members(const GCTAEventList&amp; list)
938 948 m_gti_extname = list.m_gti_extname;
939 949 m_has_phase = list.m_has_phase;
940 950 m_has_detxy = list.m_has_detxy;
  951 + m_has_mc_id = list.m_has_mc_id;
941 952  
942 953 // Copy cache
943 954 m_irf_names = list.m_irf_names;
... ... @@ -1022,6 +1033,9 @@ void GCTAEventList::read_events(const GFitsTable&amp; table) const
1022 1033 // Optionally get pointer to PHASE column
1023 1034 const GFitsTableCol* ptr_phase = (m_has_phase) ? table["PHASE"] : NULL;
1024 1035  
  1036 + // Optionally get pointer to MC_ID column
  1037 + const GFitsTableCol* ptr_mc_id = (m_has_mc_id) ? table["MC_ID"] : NULL;
  1038 +
1025 1039 // Copy data from columns into GCTAEventAtom objects
1026 1040 for (int i = 0; i < num; ++i) {
1027 1041  
... ... @@ -1046,6 +1060,11 @@ void GCTAEventList::read_events(const GFitsTable&amp; table) const
1046 1060 event.m_phase = ptr_phase->real(i);
1047 1061 }
1048 1062  
  1063 + // If available, set Monte Carlo identifier
  1064 + if (m_has_mc_id) {
  1065 + event.m_mc_id = ptr_mc_id->integer(i);
  1066 + }
  1067 +
1049 1068 // Append event
1050 1069 m_events.push_back(event);
1051 1070  
... ... @@ -1068,7 +1087,8 @@ void GCTAEventList::read_events(const GFitsTable&amp; table) const
1068 1087 name == "ENERGY" ||
1069 1088 name == "DETX" ||
1070 1089 name == "DETY" ||
1071   - name == "PHASE") {
  1090 + name == "PHASE" ||
  1091 + name == "MC_ID") {
1072 1092 continue;
1073 1093 }
1074 1094  
... ... @@ -1122,11 +1142,11 @@ void GCTAEventList::write_events(GFitsBinTable&amp; hdu) const
1122 1142 if (size() > 0) {
1123 1143  
1124 1144 // Allocate mandatory columns
1125   - GFitsTableULongCol col_eid = GFitsTableULongCol("EVENT_ID", size());
1126   - GFitsTableDoubleCol col_time = GFitsTableDoubleCol("TIME", size());
1127   - GFitsTableFloatCol col_ra = GFitsTableFloatCol("RA", size());
1128   - GFitsTableFloatCol col_dec = GFitsTableFloatCol("DEC", size());
1129   - GFitsTableFloatCol col_energy = GFitsTableFloatCol("ENERGY", size());
  1145 + GFitsTableULongCol col_eid("EVENT_ID", size());
  1146 + GFitsTableDoubleCol col_time("TIME", size());
  1147 + GFitsTableFloatCol col_ra("RA", size());
  1148 + GFitsTableFloatCol col_dec("DEC", size());
  1149 + GFitsTableFloatCol col_energy("ENERGY", size());
1130 1150  
1131 1151 // Set units of columns
1132 1152 // (see http://fits.gsfc.nasa.gov/standard30/fits_standard30aa.pdf)
... ... @@ -1155,8 +1175,8 @@ void GCTAEventList::write_events(GFitsBinTable&amp; hdu) const
1155 1175 if (m_has_detxy) {
1156 1176  
1157 1177 // Allocate columns
1158   - GFitsTableFloatCol col_detx = GFitsTableFloatCol("DETX", size());
1159   - GFitsTableFloatCol col_dety = GFitsTableFloatCol("DETY", size());
  1178 + GFitsTableFloatCol col_detx("DETX", size());
  1179 + GFitsTableFloatCol col_dety("DETY", size());
1160 1180  
1161 1181 // Set units of columns
1162 1182 // (see http://fits.gsfc.nasa.gov/standard30/fits_standard30aa.pdf)
... ... @@ -1179,18 +1199,34 @@ void GCTAEventList::write_events(GFitsBinTable&amp; hdu) const
1179 1199 if (m_has_phase) {
1180 1200  
1181 1201 // Allocate columns
1182   - GFitsTableFloatCol col_phase = GFitsTableFloatCol("PHASE", size());
  1202 + GFitsTableFloatCol col_phase("PHASE", size());
1183 1203  
1184 1204 // Fill columns
1185 1205 for (int i = 0; i < size(); ++i) {
1186 1206 col_phase(i) = m_events[i].m_phase;
1187 1207 }
1188 1208  
1189   - // Append columns to table
  1209 + // Append column to table
1190 1210 hdu.append(col_phase);
1191 1211  
1192 1212 }
1193 1213  
  1214 + // If available, add Monte Carlo identifier
  1215 + if (m_has_mc_id) {
  1216 +
  1217 + // Allocate columns
  1218 + GFitsTableLongCol col_mc_id("MC_ID", size());
  1219 +
  1220 + // Fill columns
  1221 + for (int i = 0; i < size(); ++i) {
  1222 + col_mc_id(i) = m_events[i].m_mc_id;
  1223 + }
  1224 +
  1225 + // Append column to table
  1226 + hdu.append(col_mc_id);
  1227 +
  1228 + }
  1229 +
1194 1230 // Append other columns to table
1195 1231 for (int i = 0; i < m_columns.size(); ++i) {
1196 1232 hdu.append(*m_columns[i]);
... ...