Adapt cslightcrv, csphasecrv and csscs for the usage of rectangular On/Off regions (#3248)

Jürgen Knödlseder
1 parent fba90429
Showing 12 changed files with 116 additions and 23 deletions
cscripts/cslightcrv.par
cscripts/csphasecrv.par
cscripts/csscs.par
cscripts/obsutils.py
doc/source/users/reference_manual/cslightcrv.rst
doc/source/users/reference_manual/csphasecrv.rst
doc/source/users/reference_manual/csscs.rst
test/test_cslightcrv.py
test/test_csphagen.py
test/test_csphasecrv.py
test/test_csscs.py
test/test_obsutils.py
@@ -45,8 +45,11 @@ yref,          r, a, 22.51,-90,90, &quot;Second coordinate of image/source region cen
 nxpix,         i, a, 200,,, "Size of the X axis in pixels"
 nypix,         i, a, 200,,, "Size of the Y axis in pixels"
 binsz,         r, a, 0.02,,, "Pixel size (deg/pixel)"
-srcshape,      s, h, CIRCLE,CIRCLE,, "Shape of the source region"
+srcshape,      s, a, CIRCLE,CIRCLE|RECT,, "Shape of the source region"
 rad,           r, a, 0.2,0,180, "Radius of source region circle (deg)"
+width,         r, a, 0.2,0,180, "Width of source region rectangle (deg)"
+height,        r, a, 0.2,0,180, "Height of source region rectangle (deg)"
+posang,        r, a, 0.0,-360,360, "Position angle of source region rectangle (deg)"
 bkgmethod,     s, h, REFLECTED,REFLECTED,, "Method for background estimation"
 bkgregmin,     i, h, 2,1,1000, "Minimum number of background regions"
 use_model_bkg, b, h, yes,,, "Use background model?"
@@ -42,8 +42,11 @@ yref,      r, a, 22.51,-90,90, &quot;Second coordinate of image/source region center 
 nxpix,     i, a, 200,,, "Size of the X axis in pixels"
 nypix,     i, a, 200,,, "Size of the Y axis in pixels"
 binsz,     r, a, 0.02,,, "Pixel size (deg/pixel)"
-srcshape,  s, h, CIRCLE,CIRCLE,, "Shape of the source region"
+srcshape,  s, a, CIRCLE,CIRCLE|RECT,, "Shape of the source region"
 rad,       r, a, 0.2,0,180, "Radius of source region circle (deg)"
+width,     r, a, 0.2,0,180, "Width of source region rectangle (deg)"
+height,    r, a, 0.2,0,180, "Height of source region rectangle (deg)"
+posang,    r, a, 0.0,-360,360, "Position angle of source region rectangle (deg)"
 bkgmethod, s, h, REFLECTED,REFLECTED,, "Method for background estimation"
 bkgregmin, i, h, 2,1,1000, "Minimum number of background regions"
 maxoffset, r, h, 4.0,,, "Maximum offset of source from camera center to accept observation"
@@ -45,7 +45,10 @@ emax,          r, a, 100.0,,, &quot;Upper energy limit (TeV)&quot;
 method,        s, a, UNBINNED, UNBINNED|ONOFF,, "Analysis method"
 enumbins,      i, a, 20,,, "Number of energy bins"
 bkgmethod,     s, h, REFLECTED,REFLECTED,, "Method for background estimation"
-srcshape,      s, h, CIRCLE,CIRCLE,, "Shape of the source region"
+srcshape,      s, a, CIRCLE,CIRCLE|RECT,, "Shape of the source region"
+width,         r, a, 0.2,0,180, "Width of source region rectangle (deg)"
+height,        r, a, 0.2,0,180, "Height of source region rectangle (deg)"
+posang,        r, a, 0.0,-360,360, "Position angle of source region rectangle (deg)"
 bkgregmin,     i, h, 2,1,1000, "Minimum number of background regions"
 bkgregskip,    i, h, 1,0,10, "Number of background regions to skip next to On region"
 use_model_bkg, b, h, yes,,, "Use background model?"
@@ -129,6 +129,7 @@ def sim(obs, log=False, debug=False, chatter=2, edisp=False, seed=0,
             phagen['emin']        = emin
             phagen['emax']        = emax
             phagen['enumbins']    = nbins
+            phagen['srcshape']    = 'CIRCLE'
             phagen['coordsys']    = 'CEL'
             phagen['ra']          = ra
             phagen['dec']         = dec
@@ -798,8 +799,8 @@ def get_onoff_obs(cls, obs, nthreads=0, ra = None, dec = None, srcname = &#39;&#39;):
     # User has specified custom centre for the source region
     if ra != None and dec != None:
         phagen['coordsys'] = 'CEL'
-        phagen['ra'] = ra
-        phagen['dec'] = dec
+        phagen['ra']       = ra
+        phagen['dec']      = dec
     # Otherwise use default in class
     else:
         phagen['coordsys']    = cls['coordsys'].string()
@@ -811,6 +812,10 @@ def get_onoff_obs(cls, obs, nthreads=0, ra = None, dec = None, srcname = &#39;&#39;):
             phagen['glat'] = cls['yref'].real()
     if cls['srcshape'].string() == 'CIRCLE':
         phagen['rad'] = cls['rad'].real()
+    elif cls['srcshape'].string() == 'RECT':
+        phagen['width']  = cls['width'].real()
+        phagen['height'] = cls['height'].real()
+        phagen['posang'] = cls['posang'].real()
     phagen['bkgmethod'] = cls['bkgmethod'].string()
     if cls['bkgmethod'].string() == 'REFLECTED':
         phagen['bkgregmin'] = cls['bkgregmin'].integer()
@@ -129,13 +129,23 @@ General parameters
 ``binsz [real]``
     Pixel size for ``3D`` analysis (in degrees/pixel).
-``(srcshape = CIRCLE) [string]``
-    Shape of the source region for On/Off analysis.
-    ``CIRCLE``: circular region around given position.
+``srcshape <CIRCLE|RECT> [string]``
+    Shape of the source region. ``CIRCLE`` defines a circular region around given
+    position, ``RECT`` defines a rectangle centred on a given position.
 ``rad [real]``
     Radius of source region circle for On/Off analysis (deg)
+``width [real]``
+    Width of source region rectangle for On/Off analysis (deg).
+
+``height [real]``
+    Height of source region rectangle for On/Off analysis (deg).
+
+``posang [real]``
+    Position angle of source region rectangle for On/Off analysis, counted
+    counterclockwise from celestial North (deg).
+
 ``(bkgmethod = REFLECTED) [string]``
     Method for background estimation in On/Off analysis.
     ``REFLECTED:`` background evaluated in regions with the same shape as
@@ -118,13 +118,23 @@ General parameters
 ``binsz [real]``
     Pixel size for ``3D`` analysis (in degrees/pixel).
-``(srcshape = CIRCLE) [string]``
-    Shape of the source region for On/Off analysis.
-    ``CIRCLE``: circular region around given position.
+``srcshape <CIRCLE|RECT> [string]``
+    Shape of the source region. ``CIRCLE`` defines a circular region around given
+    position, ``RECT`` defines a rectangle centred on a given position.
 ``rad [real]``
     Radius of source region circle for On/Off analysis (deg)
+``width [real]``
+    Width of source region rectangle for On/Off analysis (deg).
+
+``height [real]``
+    Height of source region rectangle for On/Off analysis (deg).
+
+``posang [real]``
+    Position angle of source region rectangle for On/Off analysis, counted
+    counterclockwise from celestial North (deg).
+
 ``(bkgmethod = REFLECTED) [string]``
     Method for background estimation in On/Off analysis.
     ``REFLECTED:`` background evaluated in regions with the same shape as
@@ -159,9 +159,19 @@ On/Off parameters
     with the same shape as source region reflected w.r.t. pointing direction for
     each observation.
-``(srcshape = CIRCLE) [string]``
-    Shape of the source region.
-    ``CIRCLE``: circular region around given position.
+``srcshape <CIRCLE|RECT> [string]``
+    Shape of the source region. ``CIRCLE`` defines a circular region around given
+    position, ``RECT`` defines a rectangle centred on a given position.
+
+``width [real]``
+    Width of source region rectangle for On/Off analysis (deg).
+
+``height [real]``
+    Height of source region rectangle for On/Off analysis (deg).
+
+``posang [real]``
+    Position angle of source region rectangle for On/Off analysis, counted
+    counterclockwise from celestial North (deg).
 ``(bkgregmin = 2) [integer]``
     Minimum number of background regions that are required for an
@@ -2,7 +2,7 @@
 # ==========================================================================
 # This scripts performs unit tests for the cslightcrv script.
 #
-# Copyright (C) 2016-2018 Juergen Knoedlseder
+# Copyright (C) 2016-2020 Juergen Knoedlseder
 #
 # This program is free software: you can redistribute it and/or modify
 # it under the terms of the GNU General Public License as published by
@@ -251,6 +251,7 @@ class Test(test):
         lcrv['tmax']          = '2020-01-01T00:05:00'
         lcrv['tbins']         = 2
         lcrv['method']        = 'ONOFF'
+        lcrv['srcshape']      = 'CIRCLE'
         lcrv['use_model_bkg'] = False  # Needed even if WSTAT
         lcrv['emin']          = 1.0
         lcrv['emax']          = 100.0
@@ -88,7 +88,7 @@ class Test(test):
                          '" ebinalg="LOG" emin=0.1 emax=100. enumbins=' + \
                          str(self._nbins) + \
                          ' etruemin=0.05 etruemax=150 etruebins=5'+ \
-                         ' coordsys="CEL" ra=83.633 dec=22.0145' + \
+                         ' srcshape="CIRCLE" coordsys="CEL" ra=83.633 dec=22.0145' + \
                          ' rad=0.2 stack="no" inexclusion="' + \
                          self._exclusion + \
                          '" bkgmethod="REFLECTED" '+ \
@@ -111,7 +111,7 @@ class Test(test):
                          '" ebinalg="LOG" emin=0.1 emax=100. enumbins=' + \
                          str(self._nbins) + \
                          ' etruemin=0.05 etruemax=150 etruebins=5'+ \
-                         ' coordsys="CEL" ra=83.633 dec=22.0145' + \
+                         ' srcshape="CIRCLE" coordsys="CEL" ra=83.633 dec=22.0145' + \
                          ' rad=0.2 stack="no" inexclusion="' + \
                          self._exclusion + \
                          '" bkgmethod="REFLECTED" '+ \
@@ -143,6 +143,7 @@ class Test(test):
         phagen['emin']        = 0.1
         phagen['emax']        = 100.0
         phagen['enumbins']    = self._nbins
+        phagen['srcshape']    = 'CIRCLE'
         phagen['coordsys']    = 'CEL'
         phagen['ra']          = 83.633
         phagen['dec']         = 22.0145
@@ -176,6 +177,7 @@ class Test(test):
         phagen['emin']      = 0.1
         phagen['emax']      = 100.0
         phagen['enumbins']  = self._nbins
+        phagen['srcshape']  = 'CIRCLE'
         phagen['coordsys']  = 'CEL'
         phagen['ra']        = 83.633
         phagen['dec']       = 22.0145
@@ -215,6 +217,7 @@ class Test(test):
         phagen['emax']        = 100.0
         phagen['enumbins']    = self._nbins
         phagen['ebingamma']   = 1.0
+        phagen['srcshape']    = 'CIRCLE'
         phagen['coordsys']    = 'CEL'
         phagen['ra']          = 83.633
         phagen['dec']         = 22.0145
@@ -247,6 +250,7 @@ class Test(test):
         phagen['emin']        = 0.1
         phagen['emax']        = 100.0
         phagen['enumbins']    = self._nbins
+        phagen['srcshape']    = 'CIRCLE'
         phagen['coordsys']    = 'CEL'
         phagen['ra']          = 83.633
         phagen['dec']         = 22.0145
@@ -346,6 +350,7 @@ class Test(test):
         phagen['emin']        = 0.1
         phagen['emax']        = 100.0
         phagen['enumbins']    = self._nbins
+        phagen['srcshape']    = 'CIRCLE'
         phagen['coordsys']    = 'CEL'
         phagen['ra']          = 83.633
         phagen['dec']         = 22.0145
@@ -372,6 +377,43 @@ class Test(test):
         self._check_output('csphagen_py7', self._nbins, self._nreg_with_excl)
         self._check_outobs('csphagen_py7_obs.xml', 1)
+        # Test with rectangular region
+        phagen = cscripts.csphagen()
+        phagen['inobs']       = self._myevents1
+        phagen['inmodel']     = 'NONE'
+        phagen['caldb']       = self._caldb
+        phagen['irf']         = self._irf
+        phagen['ebinalg']     = 'LOG'
+        phagen['emin']        = 0.1
+        phagen['emax']        = 100.0
+        phagen['enumbins']    = self._nbins
+        phagen['srcshape']    = 'RECT'
+        phagen['width']       = 0.15
+        phagen['height']      = 0.25
+        phagen['posang']      = 45.0
+        phagen['coordsys']    = 'CEL'
+        phagen['ra']          = 83.633
+        phagen['dec']         = 22.0145
+        phagen['rad']         = 0.2
+        phagen['stack']       = False
+        phagen['inexclusion'] = self._exclusion
+        phagen['bkgmethod']   = 'REFLECTED'
+        phagen['etruemin']    = 0.05
+        phagen['etruemax']    = 150.0
+        phagen['etruebins']   = 5
+        phagen['outobs']      = 'csphagen_py8_obs.xml'
+        phagen['outmodel']    = 'csphagen_py8_model.xml'
+        phagen['prefix']      = 'csphagen_py8'
+        phagen['logfile']     = 'csphagen_py8.log'
+        phagen['chatter']     = 2
+
+        # Execute script
+        phagen.execute()
+
+        # Check output
+        self._check_output('csphagen_py8', self._nbins, self._nreg_with_excl)
+        self._check_outobs('csphagen_py8_obs.xml', 1)
+
         # Return
         return
@@ -393,6 +435,7 @@ class Test(test):
         phagen['emin']        = 0.1
         phagen['emax']        = 100.0
         phagen['enumbins']    = self._nbins
+        phagen['srcshape']    = 'CIRCLE'
         phagen['coordsys']    = 'CEL'
         phagen['ra']          = 83.633
         phagen['dec']         = 22.0145
@@ -2,7 +2,7 @@
 # ==========================================================================
 # This scripts performs unit tests for the csphasecrv script.
 #
-# Copyright (C) 2017-2018 Juergen Knoedlseder
+# Copyright (C) 2017-2020 Juergen Knoedlseder
 #
 # This program is free software: you can redistribute it and/or modify
 # it under the terms of the GNU General Public License as published by
@@ -247,6 +247,7 @@ class Test(test):
         pcrv['phbinalg']  = 'LIN'
         pcrv['phbins']    = 2
         pcrv['method']    = 'ONOFF'
+        pcrv['srcshape']  = 'CIRCLE'
         pcrv['statistic'] = 'WSTAT'
         pcrv['emin']      = 1.0
         pcrv['emax']      = 100.0
@@ -299,8 +299,8 @@ class Test(test):
         script['irf']         = self._irf
         script['inmodel']     = self._model
         script['srcnames']    = 'Crab'
-        script['emin']        = 1.
-        script['emax']        = 100.
+        script['emin']        =   1.0
+        script['emax']        = 100.0
         script['nxpix']       = 2
         script['nypix']       = 2
         script['binsz']       = 0.1
@@ -310,6 +310,7 @@ class Test(test):
         script['yref']        = 22.01
         script['rad']         = 0.2
         script['method']      = 'ONOFF'
+        script['srcshape']    = 'CIRCLE'
         script['inexclusion'] = self._exclusion
         script['enumbins']    = 2
         script['outfile']     = 'csscs_py6.fits'
@@ -330,8 +331,8 @@ class Test(test):
         script['irf']      = self._irf
         script['inmodel']  = self._model
         script['srcnames'] = 'Crab'
-        script['emin']     = 1.
-        script['emax']     = 100.
+        script['emin']     =   1.0
+        script['emax']     = 100.0
         script['nxpix']    = 2
         script['nypix']    = 2
         script['binsz']    = 0.1
@@ -341,6 +342,7 @@ class Test(test):
         script['yref']     = 22.01
         script['rad']      = 0.2
         script['method']   = 'ONOFF'
+        script['srcshape'] = 'CIRCLE'
         script['enumbins'] = 2
         script['outfile']  = 'csscs_py7.fits'
         script['chatter']  = 2
@@ -375,6 +377,7 @@ class Test(test):
         script['yref']          = 22.01
         script['rad']           = 0.2
         script['method']        = 'ONOFF'
+        script['srcshape']      = 'CIRCLE'
         script['use_model_bkg'] = False
         script['inexclusion']   = self._exclusion
         script['enumbins']      = 2
@@ -2,7 +2,7 @@
 # ==========================================================================
 # This scripts performs unit tests for the obsutils module
 #
-# Copyright (C) 2017-2018 Juergen Knoedlseder
+# Copyright (C) 2017-2020 Juergen Knoedlseder
 #
 # This program is free software: you can redistribute it and/or modify
 # it under the terms of the GNU General Public License as published by
@@ -670,6 +670,7 @@ class Test(test):
         lcrv['tmax']     = '2020-01-01T00:05:00'
         lcrv['tbins']    = 2
         lcrv['method']   = 'ONOFF'
+        lcrv['srcshape'] = 'CIRCLE'
         lcrv['emin']     = 1.0
         lcrv['emax']     = 100.0
         lcrv['enumbins'] = 2