Lili Yang / ctools

Commit 84c6461fdbdafb5bd14278c098f7092ff145901e

Authored by Jürgen Knödlseder
1 parent 44315163

Enhance ctselect branch coverage (#1828)

Inline Side-by-side
Showing 2 changed files with 177 additions and 54 deletions
src/ctselect/ctselect.cpp
  View file @84c6461
... ... @@ -728,6 +728,9 @@ void ctselect::select_events(GCTAObservation* obs,
728 728 // Make RoI selection
729 729 if (m_select_roi) {
730 730  
  731 + // Store the original radius
  732 + double original_rad = rad;
  733 +
731 734 // Log the requested RoI
732 735 log_value(NORMAL, "Requested RoI",
733 736 "Centre(RA,DEC)=("+gammalib::str(ra)+", "+
... ... @@ -751,17 +754,20 @@ void ctselect::select_events(GCTAObservation* obs,
751 754 }
752 755 }
753 756  
754   - // Log the selected RoI
755   - if (logNormal()) {
756   - log << gammalib::parformat("Selected RoI");
757   - if (rad <= 0.0) {
758   - log << "None. There is no overlap between existing ";
759   - log << "and requested ROI." << std::endl;
760   - }
761   - else {
762   - log << "Centre(RA,DEC)=(" << ra << ", " << dec << ") deg, ";
763   - log << "Radius=" << rad << " deg" << std::endl;
764   - }
  757 + // If the RoI radius is negative then there is no overlap between
  758 + // the existing and the requested RoI and the radius will be restored
  759 + // to the original value. We signal that all events should be
  760 + // removed.
  761 + if (rad <= 0.0) {
  762 + rad = original_rad;
  763 + remove_all = true;
  764 + log_value(NORMAL, "Selected RoI", "None. There is no overlap "
  765 + "between existing and requested RoI.");
  766 + }
  767 + else {
  768 + log_value(NORMAL, "Selected RoI", "Centre(RA,DEC)=("+
  769 + gammalib::str(ra)+", "+gammalib::str(dec)+") deg, "
  770 + "Radius="+gammalib::str(rad)+" deg");
765 771 }
766 772  
767 773 // Format RoI selection
... ... @@ -813,12 +819,13 @@ void ctselect::select_events(GCTAObservation* obs,
813 819 throw GException::invalid_value(G_SELECT_EVENTS, msg);
814 820 }
815 821  
816   - // Determine number of events in the events HDU
817   - int nevents = file.table(evtname)->nrows();
  822 + // Determine number of events in the events HDU. If removal of all events
  823 + // has been requested then set the number of events to zero.
  824 + int nevents = (remove_all) ? 0 : file.table(evtname)->nrows();
818 825  
819   - // If the selected event list is empty or if removal of all events
820   - // has been requested then append an empty event list to the observation.
821   - if ((nevents < 1) || (remove_all)) {
  826 + // If the selected event list is empty then append an empty event list
  827 + // to the observation
  828 + if (nevents < 1) {
822 829  
823 830 // Create empty event list
824 831 GCTAEventList eventlist;
... ...
test/test_ctselect.py
  View file @84c6461
... ... @@ -72,8 +72,8 @@ class Test(test):
72 72 cmd = ctselect+' inobs="'+self._events+'"'+ \
73 73 ' outobs="ctselect_cmd1.fits"'+ \
74 74 ' ra=83.63 dec=22.01 rad=3.0'+ \
75   - ' tmin=0.0 tmax=1800.0'+ \
76   - ' emin=0.1 emax=100.0'+ \
  75 + ' tmin=500.0 tmax=1000.0'+ \
  76 + ' emin=0.2 emax=80.0'+ \
77 77 ' logfile="ctselect_cmd1.log" chatter=1'
78 78  
79 79 # Check if execution of wrong command fails
... ... @@ -91,8 +91,8 @@ class Test(test):
91 91 cmd = ctselect+' inobs="event_file_that_does_not_exist.fits"'+ \
92 92 ' outobs="ctselect_cmd2.fits"'+ \
93 93 ' ra=83.63 dec=22.01 rad=3.0'+ \
94   - ' tmin=0.0 tmax=1800.0'+ \
95   - ' emin=0.1 emax=100.0'+ \
  94 + ' tmin=500.0 tmax=1000.0'+ \
  95 + ' emin=0.2 emax=80.0'+ \
96 96 ' logfile="ctselect_cmd2.log" chatter=1'
97 97  
98 98 # Check if execution failed
... ... @@ -143,10 +143,10 @@ class Test(test):
143 143 select['ra'] = 83.63
144 144 select['dec'] = 22.01
145 145 select['rad'] = 3
146   - select['tmin'] = 0
147   - select['tmax'] = 1800
148   - select['emin'] = 0.1
149   - select['emax'] = 100.0
  146 + select['tmin'] = 500.0
  147 + select['tmax'] = 1000.0
  148 + select['emin'] = 0.2
  149 + select['emax'] = 80.0
150 150 select['outobs'] = 'ctselect_py1.fits'
151 151 select['logfile'] = 'ctselect_py1.log'
152 152 select['chatter'] = 2
... ... @@ -166,7 +166,9 @@ class Test(test):
166 166 self._check_obs(cpy_select.obs())
167 167  
168 168 # Execute copy of ctselect tool again, now with a higher chatter
169   - # level than before and without any selection
  169 + # level than before and without any selection. Since the tools
  170 + # still holds the same selected observation container from before
  171 + # the number of events will be identical.
170 172 cpy_select['ra'] = 'NONE'
171 173 cpy_select['dec'] = 'NONE'
172 174 cpy_select['rad'] = 'NONE'
... ... @@ -184,11 +186,16 @@ class Test(test):
184 186 self._check_result_file('ctselect_py2.fits')
185 187  
186 188 # Execute again the copy of ctselect tool again, now using the pointing
187   - # information in the input observation
188   - cpy_select['usepnt'] = True
189   - cpy_select['outobs'] = 'ctselect_py3.fits'
190   - cpy_select['logfile'] = 'ctselect_py3.log'
191   - cpy_select['chatter'] = 4
  189 + # information in the input observation. We also set the minimum values
  190 + # to valid event selections, but since the maximum values are still
  191 + # 'NONE', no event selections should occur.
  192 + cpy_select['usepnt'] = True
  193 + cpy_select['tmin'] = 500.0
  194 + cpy_select['emin'] = 0.2
  195 + cpy_select['usethres'] = 'USER'
  196 + cpy_select['outobs'] = 'ctselect_py3.fits'
  197 + cpy_select['logfile'] = 'ctselect_py3.log'
  198 + cpy_select['chatter'] = 4
192 199 cpy_select.logFileOpen() # Needed to get a new log file
193 200 cpy_select.execute()
194 201  
... ... @@ -204,20 +211,27 @@ class Test(test):
204 211 # Get mixed observation container
205 212 obs = self._obs_mixed()
206 213  
207   - # Now setup ctselect tool from observation container. We make an
208   - # event selection that will result in zero events in the file.
  214 + # Attach response function to first observation which is the
  215 + # event list. This is necessary to run the ctselect tool with
  216 + # "DEFAULT" thresholds.
  217 + obs[0].response('South_0.5h', gammalib.GCaldb('cta', 'prod2'))
  218 +
  219 + # Setup ctselect tool from observation container. An energy range
  220 + # beyond the energies covered in the event file is specified, hence
  221 + # an empty event list will be saved.
209 222 select = ctools.ctselect(obs)
210   - select['ra'] = 83.63
211   - select['dec'] = 22.01
212   - select['rad'] = 3.0
213   - select['tmin'] = 0
214   - select['tmax'] = 1800
215   - select['emin'] = 120.0
216   - select['emax'] = 130.0
217   - select['expr'] = 'DETX == 0'
218   - select['outobs'] = 'ctselect_py4.fits'
219   - select['logfile'] = 'ctselect_py4.log'
220   - select['chatter'] = 3
  223 + select['ra'] = 83.63
  224 + select['dec'] = 22.01
  225 + select['rad'] = 3.0
  226 + select['tmin'] = 500.0
  227 + select['tmax'] = 1000.0
  228 + select['emin'] = 120.0
  229 + select['emax'] = 130.0
  230 + select['expr'] = 'DETX == 0'
  231 + select['usethres'] = 'DEFAULT' # Has no impact as IRF has no keywords
  232 + select['outobs'] = 'ctselect_py4.fits'
  233 + select['logfile'] = 'ctselect_py4.log'
  234 + select['chatter'] = 3
221 235  
222 236 # Execute tool
223 237 select.logFileOpen() # Needed to get a new log file
... ... @@ -232,10 +246,10 @@ class Test(test):
232 246 select['ra'] = 83.63
233 247 select['dec'] = 22.01
234 248 select['rad'] = 3
235   - select['tmin'] = 0
236   - select['tmax'] = 1800
237   - select['emin'] = 0.1
238   - select['emax'] = 100.0
  249 + select['tmin'] = 500.0
  250 + select['tmax'] = 1000.0
  251 + select['emin'] = 0.2
  252 + select['emax'] = 80.0
239 253 select['outobs'] = 'ctselect_py5.fits'
240 254 select['logfile'] = 'ctselect_py5.log'
241 255 select['chatter'] = 3
... ... @@ -249,12 +263,106 @@ class Test(test):
249 263 except ValueError:
250 264 self.test_try_success()
251 265  
  266 + # Setup ctselect tool from event list with extension name. The "emax"
  267 + # value should be ignored.
  268 + select = ctools.ctselect()
  269 + select['inobs'] = self._events+'[EVENTS]'
  270 + select['ra'] = 83.63
  271 + select['dec'] = 22.01
  272 + select['rad'] = 3.0
  273 + select['tmin'] = 500.0
  274 + select['tmax'] = 1000.0
  275 + select['emin'] = 0.2
  276 + select['emax'] = 0.0 # Signals that "emax" should be ignored
  277 + select['outobs'] = 'ctselect_py6.fits'
  278 + select['logfile'] = 'ctselect_py6.log'
  279 + select['chatter'] = 3
  280 +
  281 + # Execute tool
  282 + select.logFileOpen() # Needed to get a new log file
  283 + select.execute()
  284 +
  285 + # Check result file
  286 + self._check_result_file('ctselect_py6.fits', nevents=595)
  287 +
  288 + # Now ignore the "emin" value.
  289 + select = ctools.ctselect()
  290 + select['inobs'] = self._events+'[EVENTS]'
  291 + select['ra'] = 83.63
  292 + select['dec'] = 22.01
  293 + select['rad'] = 3.0
  294 + select['tmin'] = 500.0
  295 + select['tmax'] = 1000.0
  296 + select['emin'] = 0.0 # Signals that "emin" should be ignored
  297 + select['emax'] = 80.0
  298 + select['outobs'] = 'ctselect_py7.fits'
  299 + select['logfile'] = 'ctselect_py7.log'
  300 + select['chatter'] = 3
  301 +
  302 + # Execute tool
  303 + select.logFileOpen() # Needed to get a new log file
  304 + select.execute()
  305 +
  306 + # Check result file
  307 + self._check_result_file('ctselect_py7.fits', nevents=1684)
  308 +
  309 + # Now set "emin > emax"
  310 + select['emin'] = 150.0
  311 + select['emax'] = 80.0
  312 + select['outobs'] = 'ctselect_py8.fits'
  313 + select['logfile'] = 'ctselect_py8.log'
  314 + select['chatter'] = 3
  315 +
  316 + # Execute tool
  317 + select.logFileOpen() # Needed to get a new log file
  318 + select.execute()
  319 +
  320 + # Check result file
  321 + self._check_result_file('ctselect_py8.fits', nevents=0)
  322 +
  323 + # Setup ctselect tool with an RoI that is displaced and of the same
  324 + # size as the original RoI. This should reduce the RoI radius
  325 + # automatically.
  326 + select = ctools.ctselect()
  327 + select['inobs'] = self._events
  328 + select['ra'] = 83.63
  329 + select['dec'] = 24.01
  330 + select['rad'] = 5.0
  331 + select['tmin'] = 500.0
  332 + select['tmax'] = 1000.0
  333 + select['emin'] = 0.2
  334 + select['emax'] = 80.0
  335 + select['outobs'] = 'ctselect_py9.fits'
  336 + select['logfile'] = 'ctselect_py9.log'
  337 + select['chatter'] = 3
  338 +
  339 + # Execute tool
  340 + select.logFileOpen() # Needed to get a new log file
  341 + select.execute()
  342 +
  343 + # Check result file
  344 + self._check_result_file('ctselect_py9.fits', nevents=466, dec=24.01)
  345 +
  346 + # Now put the RoI outside the existing RoI. This should leave to
  347 + # an empty event list.
  348 + select['ra'] = 83.63
  349 + select['dec'] = -22.01
  350 + select['outobs'] = 'ctselect_py10.fits[EVENTS;GTI]'
  351 + select['logfile'] = 'ctselect_py10.log'
  352 + select['chatter'] = 3
  353 +
  354 + # Execute tool
  355 + select.logFileOpen() # Needed to get a new log file
  356 + select.execute()
  357 +
  358 + # Check result file
  359 + self._check_result_file('ctselect_py10.fits', nevents=0, dec=-22.01, rad=5.0)
252 360  
253 361 # Return
254 362 return
255 363  
256 364 # Check result file
257   - def _check_result_file(self, filename, nevents=6127):
  365 + def _check_result_file(self, filename, nevents=591, dec=22.01, rad=3.0):
258 366 """
259 367 Check result file
260 368  
... ... @@ -264,18 +372,22 @@ class Test(test):
264 372 Event list file name
265 373 nevents : int, optional
266 374 Expected number of events
  375 + dec : float, optional
  376 + Expected Declination (deg)
  377 + rad : float, optional
  378 + Expected radius (deg)
267 379 """
268 380 # Open result file
269 381 events = gammalib.GCTAEventList(filename)
270 382  
271 383 # Check event list
272   - self._check_events(events, nevents=nevents)
  384 + self._check_events(events, nevents=nevents, dec=dec, rad=rad)
273 385  
274 386 # Return
275 387 return
276 388  
277 389 # Check observation and event list
278   - def _check_obs(self, obs, nobs=1, nevents=6127):
  390 + def _check_obs(self, obs, nobs=1, nevents=591):
279 391 """
280 392 Check observation and event list
281 393  
... ... @@ -300,7 +412,7 @@ class Test(test):
300 412 return
301 413  
302 414 # Check events
303   - def _check_events(self, events, nevents=6127):
  415 + def _check_events(self, events, nevents=591, dec=22.01, rad=3.0):
304 416 """
305 417 Check event list
306 418  
... ... @@ -310,14 +422,18 @@ class Test(test):
310 422 Event list
311 423 nevents : int, optional
312 424 Expected number of events
  425 + dec : float, optional
  426 + Expected Declination (deg)
  427 + rad : float, optional
  428 + Expected radius (deg)
313 429 """
314 430 # Check event list
315 431 self.test_value(events.size(), nevents, 'Check number of events')
316   - self.test_value(events.roi().centre().dir().ra_deg(), 83.63, 1.0e-6,
  432 + self.test_value(events.roi().centre().dir().ra_deg(), 83.63,
317 433 'Check for RoI centre Right Ascension')
318   - self.test_value(events.roi().centre().dir().dec_deg(), 22.01, 1.0e-6,
  434 + self.test_value(events.roi().centre().dir().dec_deg(), dec,
319 435 'Check for RoI centre Declination')
320   - self.test_value(events.roi().radius(), 3.0, 1.0e-6,
  436 + self.test_value(events.roi().radius(), rad,
321 437 'Check for RoI radius')
322 438  
323 439 # Return
... ...