From 0ad65eddf5abda46d9a1d05ea3f3d8a5c160c4d8 Mon Sep 17 00:00:00 2001
From: Aaron LI <aaronly.me@outlook.com>
Date: Thu, 16 Feb 2017 15:57:00 +0800
Subject: Add acis.py and spectrum.py modules

Also add astropy to the runtime dependencies.
---
 README.md           |  1 +
 scripts/acis.py     | 34 ++++++++++++++++++++++++++++++++++
 scripts/spectrum.py | 44 ++++++++++++++++++++++++++++++++++++++++++++
 3 files changed, 79 insertions(+)
 create mode 100644 scripts/acis.py
 create mode 100644 scripts/spectrum.py

diff --git a/README.md b/README.md
index 386b664..dbf0d3b 100644
--- a/README.md
+++ b/README.md
@@ -85,6 +85,7 @@ Installation
    * [ruamel.yaml](https://bitbucket.org/ruamel/yaml)
    * [numpy](http://numpy.org/)
    * [scipy](https://scipy.org/)
+   * [astropy](http://www.astropy.org/)
 
 
 Settings
diff --git a/scripts/acis.py b/scripts/acis.py
new file mode 100644
index 0000000..5200f4e
--- /dev/null
+++ b/scripts/acis.py
@@ -0,0 +1,34 @@
+# Copyright (c) 2017 Weitian LI <liweitianux@live.com>
+# MIT license
+
+"""
+Chandra ACIS utilities
+"""
+
+import math
+
+
+class ACIS:
+    """
+    Chandra ACIS detector properties and utilities.
+
+    References
+    ----------
+    [1] CIAO - dictionary - ACIS (advanced camera for imaging and spectroscopy)
+        http://cxc.harvard.edu/ciao/dictionary/acis.html
+    [2] CIAO - Dictionary - PI (pulse invariant)
+        http://cxc.harvard.edu/ciao/dictionary/pi.html
+    """
+    # Pixel size
+    pixel2arcsec = 0.492  # [arcsec]
+    # Number of channels
+    nchannel = 1024
+    # Channel energy width
+    echannel = 14.6  # [eV]
+
+    @classmethod
+    def energy2channel(self, energy):
+        """
+        Convert energy [eV] to channel number.
+        """
+        return math.floor(energy/self.echannel + 1)
diff --git a/scripts/spectrum.py b/scripts/spectrum.py
new file mode 100644
index 0000000..ad38b0d
--- /dev/null
+++ b/scripts/spectrum.py
@@ -0,0 +1,44 @@
+# Copyright (c) 2017 Weitian LI <liweitianux@live.com>
+# MIT license
+
+"""
+Chandra ACIS spectrum.
+"""
+
+
+from astropy.io import fits
+
+from acis import ACIS
+
+
+class Spectrum:
+    """
+    Chandra ACIS spectrum
+    """
+    def __init__(self, filepath):
+        self.filepath = filepath
+        self.fitsobj = fits.open(filepath)
+        ext_spec = self.fitsobj["SPECTRUM"]
+        self.header = ext_spec.header
+        # spectral data
+        self.channel = ext_spec.data.columns["CHANNEL"].array
+        self.counts = ext_spec.data.columns["COUNTS"].array
+        # spectral keywords
+        self.EXPOSURE = self.header.get("EXPOSURE")
+        self.BACKSCAL = self.header.get("BACKSCAL")
+
+    def calc_pb_flux(self, elow=9500, ehigh=12000):
+        """
+        Calculate the particle background flux:
+            flux = counts / exposure / area
+
+        Parameters
+        ----------
+        elow, ehigh : float, optional
+            Lower and upper energy limit for the particle background.
+        """
+        chlow = ACIS.energy2channel(elow)
+        chhigh = ACIS.energy2channel(ehigh)
+        counts = self.counts[(chlow-1):chhigh].sum()
+        flux = counts / self.EXPOSURE / self.BACKSCAL
+        return flux
-- 
cgit v1.2.2