Add drm gen for CTIME (daily) data

gbm_drm_gen/__init__.py

@@ -1,24 +1,24 @@
-#from .drmgen import DRMGen
+# from .drmgen import DRMGen
 from .drmgen_tte import DRMGenTTE
 from .drmgen_trig import DRMGenTrig
+from .drmgen_ctime import DRMGenCTIME
 
 from .input_edges import NaiTTEEdges, BgoTTEEdges
 
 from gbm_drm_gen.io.palantir import Palantir
+
 try:
-    
+
     from gbm_drm_gen.io.balrog_like import BALROGLike
     from gbm_drm_gen.io.balrog_drm import BALROG_DRM
-    from gbm_drm_gen.io.balrog_healpix_map import BALROGHealpixMap    
+    from gbm_drm_gen.io.balrog_healpix_map import BALROGHealpixMap
 
 except:
     BALROGLike = None
     BALROG_DRM = None
     BALROGHealpixMap = None
 
 
-
-
 __all__ = [
     "DRMGenTTE",
     "DRMGenTrig",
@@ -31,5 +31,6 @@
 ]
 
 from ._version import get_versions
-__version__ = get_versions()['version']
+
+__version__ = get_versions()["version"]
 del get_versions

gbm_drm_gen/drmgen_ctime.py

@@ -0,0 +1,153 @@
+__author__ = "fkunzweiler"
+import numpy as np
+import astropy.io.fits as fits
+
+import gbmgeometry
+
+import astropy.units as u
+
+from gbm_drm_gen.drmgen import DRMGen
+from gbm_drm_gen.input_edges import ctime_edges
+
+det_name_lookup = {
+    "NAI_00": 0,
+    "NAI_01": 1,
+    "NAI_02": 2,
+    "NAI_03": 3,
+    "NAI_04": 4,
+    "NAI_05": 5,
+    "NAI_06": 6,
+    "NAI_07": 7,
+    "NAI_08": 8,
+    "NAI_09": 9,
+    "NAI_10": 10,
+    "NAI_11": 11,
+    "BGO_00": 12,
+    "BGO_01": 13,
+}
+
+
+class DRMGenCTIME(DRMGen):
+    """
+    A CTIME specific drmgen already incorporating the standard input edges. Output edges are obtained
+    from the input CTIME file. Spacecraft position is read from the poshist file. For further details see the
+    generic reader (DRMGen).
+
+    :param trigdat: the path to a trigdat file
+    :param det: the number (0-13) of the detector to be used
+    :param mat_type: 0=direct 1=scattered 2=direct+scattered
+    :param time: time relative to trigger to pull spacecraft position or MET if using a poshist file
+    :param cspecfile: the cspecfile to pull energy output edges from
+    :param poshist: read a poshist file
+    """
+
+    def __init__(
+        self,
+        ctime_file=None,
+        det_name=None,
+        time=0.0,
+        trigdat=None,
+        poshist=None,
+        T0=None,
+        mat_type=0,
+        occult=False,
+    ):
+
+        self._occult = occult
+
+        self._time = time
+
+        self._matrix_type = mat_type
+
+        if det_name is None:
+
+            assert ctime_file is not None
+
+            with fits.open(ctime_file) as f:
+
+                det_name = f["PRIMARY"].header["DETNAM"]
+
+        else:
+            assert det_name in det_name_lookup, "must use a valid detector name"
+
+        self._det_number = det_name_lookup[det_name]
+
+        if self._det_number > 11:
+            self._in_edge = ctime_edges["bgo"]
+
+        else:
+            self._in_edge = ctime_edges["nai"]
+
+        # Create the out edge energies
+        with fits.open(ctime_file) as f:
+            out_edge_start = f["EBOUNDS"].data["E_MIN"]
+            out_edge_stop = f["EBOUNDS"].data["E_MAX"][-1]
+
+        out_edge = np.append(out_edge_start, out_edge_stop)
+
+        self._out_edge = out_edge
+
+        if trigdat is not None:
+
+            try:
+                self._position_interpolator = (
+                    gbmgeometry.PositionInterpolator.from_trigdat(trigdat_file=trigdat)
+                )
+
+            except:
+
+                self._position_interpolator = (
+                    gbmgeometry.PositionInterpolator.from_trigdat_hdf5(
+                        trigdat_file=trigdat
+                    )
+                )
+
+            self._gbm = gbmgeometry.GBM(
+                self._position_interpolator.quaternion(time),
+                self._position_interpolator.sc_pos(time) * u.km,
+            )
+
+        elif poshist is not None:
+
+            try:
+
+                self._position_interpolator = (
+                    gbmgeometry.PositionInterpolator.from_poshist(
+                        poshist_file=poshist, T0=T0
+                    )
+                )
+
+            except:
+
+                self._position_interpolator = (
+                    gbmgeometry.PositionInterpolator.from_poshist_hdf5(
+                        poshist_file=poshist, T0=T0
+                    )
+                )
+
+            self._gbm = gbmgeometry.GBM(
+                self._position_interpolator.quaternion(time),
+                self._position_interpolator.sc_pos(time) * u.m,
+            )
+
+        else:
+
+            raise RuntimeError("No trigdat or posthist file used!")
+
+        self._sc_quaternions_updater()
+
+    def _sc_quaternions_updater(self):
+
+        quaternions = self._position_interpolator.quaternion(self._time)
+
+        sc_pos = self._position_interpolator.sc_pos(self._time)
+
+        super(DRMGenCTIME, self).__init__(
+            quaternions=quaternions,
+            sc_pos=sc_pos,
+            det_number=self._det_number,
+            ebin_edge_in=self._in_edge,
+            mat_type=self._matrix_type,
+            ebin_edge_out=self._out_edge,
+            occult=self._occult,
+        )

gbm_drm_gen/input_edges.py

@@ -12,15 +12,15 @@ def __init__(self, edges: np.array):
 
         diff = np.diff(edges)
 
-        if not np.all(diff > 0.):
+        if not np.all(diff > 0.0):
 
-            raise RuntimeError(
-                "The input edges are not monotonically increasing")
+            raise RuntimeError("The input edges are not monotonically increasing")
 
         if len(edges) % 2 == 0:
 
             raise RuntimeError(
-                "The edges must have an odd length with the end points closed")
+                "The edges must have an odd length with the end points closed"
+            )
 
         self._edges = edges
 
@@ -368,16 +368,16 @@ def from_custom_array(cls, edges: np.array):
 
 class NaiTTEEdges(InputEdges):
     _emin = 5.0
-    _emax = 50000.
+    _emax = 50000.0
 
     def __init__(self, edges):
 
         super(NaiTTEEdges, self).__init__(edges)
 
 
 class BgoTTEEdges(InputEdges):
-    _emin = 100.
-    _emax = 200000.
+    _emin = 100.0
+    _emax = 200000.0
 
     def __init__(self, edges):
 
@@ -684,11 +684,15 @@ def __init__(self, edges):
 
 trigdat_out_edge = {}
 trigdat_out_edge["nai"] = np.array(
-    [3.4, 10.0, 22.0, 44.0, 95.0, 300.0, 500.0, 800.0, 2000.0], dtype=np.float32,
+    [3.4, 10.0, 22.0, 44.0, 95.0, 300.0, 500.0, 800.0, 2000.0],
+    dtype=np.float32,
 )
 
 
 trigdat_out_edge["bgo"] = np.array(
-    [150.0, 400.0, 850.0, 1500.0, 3000.0, 5500.0, 10000.0, 20000.0, 50000.0, ],
+    [150.0, 400.0, 850.0, 1500.0, 3000.0, 5500.0, 10000.0, 20000.0, 50000.0],
     dtype=np.float32,
 )
+
+# For CTIME data we use the same in edges as for trigdat data
+ctime_edges = trigdat_edges