Merge pull request #168 from kbarbary/datacov

add data flux covariance

docs/conf.py

@@ -55,12 +55,12 @@
     'sphinx.ext.intersphinx',
     'sphinx.ext.inheritance_diagram',
     'sphinx.ext.mathjax',
-    'sphinx.ext.napoleon',
     'sphinx_gallery.gen_gallery',
+    'numpydoc',
     matplotlib.sphinxext.plot_directive.__name__
 ]
 
-#numpydoc_show_class_members = False
+numpydoc_show_class_members = False
 autosummary_generate = True
 autoclass_content = "class"
 autodoc_default_flags = ["members", "no-special-members"]

docs/history.rst

@@ -11,6 +11,16 @@ example, due to changes in integration method.)
 v1.5.0 (unreleased)
 ===================
 
+- Added support for covariance in photometric data measurements. Covariance
+  is stored as a ``'fluxcov'`` column in the table of measurements.
+
+- Added support for reading snfit-format "covmat" files into a table of
+  photometry::
+
+    >>> data = read_lc('filename', format='salt2', read_covmat=True)
+    >>> data['Fluxcov'].shape == (len(data), len(data))
+    True
+
 v1.4.0 (2016-11-16)
 ===================
 

docs/reference.rst

@@ -84,6 +84,7 @@ Fitting Photometric Data
 .. autosummary::
    :toctree: api
 
+   select_data
    chisq
    flatten_result
 

sncosmo/__init__.py

@@ -162,6 +162,7 @@ class Conf(ConfigNamespace):
     from .fitting import *
     from .simulation import *
     from .plotting import *
+    from .photdata import *
     from .registry import *
 
     from . import registry  # deprecated in v1.2; use previous import.

sncosmo/io.py

@@ -8,7 +8,7 @@
 import sys
 import re
 import json
-from collections import OrderedDict as odict
+from collections import OrderedDict
 
 import numpy as np
 from astropy.table import Table
@@ -207,13 +207,9 @@ def write_griddata_fits(x0, x1, y, name_or_obj):
 
 # -----------------------------------------------------------------------------
 # Reader: ascii
-def _read_ascii(f, **kwargs):
+def _read_ascii(f, delim=None, metachar='@', commentchar='#'):
 
-    delim = kwargs.get('delim', None)
-    metachar = kwargs.get('metachar', '@')
-    commentchar = kwargs.get('commentchar', '#')
-
-    meta = odict()
+    meta = OrderedDict()
     colnames = []
     cols = []
     readingdata = False
@@ -249,22 +245,27 @@ def _read_ascii(f, **kwargs):
         for col, item in zip(cols, items):
             col.append(_cast_str(item))
 
-    data = odict(zip(colnames, cols))
+    data = OrderedDict(zip(colnames, cols))
     return meta, data
 
 
 # -----------------------------------------------------------------------------
 # Reader: salt2
 
-def _read_salt2(f, **kwargs):
+def _read_salt2(name_or_obj, read_covmat=False):
     """Read a new-style SALT2 file.
 
     Such a file has metadata on lines starting with '@' and column names
     on lines starting with '#' and containing a ':' after the column name.
     There is optionally a line containing '#end' before the start of data.
     """
 
-    meta = odict()
+    if isinstance(name_or_obj, six.string_types):
+        f = open(name_or_obj, 'r')
+    else:
+        f = name_or_obj
+
+    meta = OrderedDict()
     colnames = []
     cols = []
     readingdata = False
@@ -313,22 +314,35 @@ def _read_salt2(f, **kwargs):
         for col, item in zip(cols, items):
             col.append(_cast_str(item))
 
-    data = odict(zip(colnames, cols))
+    if isinstance(name_or_obj, six.string_types):
+        f.close()
+
+    # read covariance matrix file, if requested and present
+    if read_covmat and 'COVMAT' in meta:
+        fname = os.path.join(os.path.dirname(f.name), meta['COVMAT'])
+
+        # use skiprows=1 because first row has array dimensions
+        fluxcov = np.loadtxt(fname, skiprows=1)
+
+        # asethetics: capitalize 'Fluxcov' to match salt2 colnames
+        # such as 'Fluxerr'
+        colnames.append('Fluxcov')
+        cols.append(fluxcov)
+
+    data = OrderedDict(zip(colnames, cols))
 
     return meta, data
 
 
 # -----------------------------------------------------------------------------
 # Reader: salt2-old
 
-def _read_salt2_old(dirname, **kwargs):
+def _read_salt2_old(dirname, filenames=None):
     """Read old-style SALT2 files from a directory.
 
     A file named 'lightfile' must exist in the directory.
     """
 
-    filenames = kwargs.get('filenames', None)
-
     # Get list of files in directory.
     if not (os.path.exists(dirname) and os.path.isdir(dirname)):
         raise IOError("Not a directory: '{0}'".format(dirname))
@@ -338,7 +352,7 @@ def _read_salt2_old(dirname, **kwargs):
     if 'lightfile' not in dirfilenames:
         raise IOError("no lightfile in directory: '{0}'".format(dirname))
     with open(os.path.join(dirname, 'lightfile'), 'r') as lightfile:
-        meta = odict()
+        meta = OrderedDict()
         for line in lightfile.readlines():
             line = line.strip()
             if len(line) == 0:
@@ -401,7 +415,7 @@ def _read_salt2_old(dirname, **kwargs):
 
 # -----------------------------------------------------------------------------
 # Reader: json
-def _read_json(f, **kwargs):
+def _read_json(f):
     t = json.load(f)
 
     # Encode data keys as ascii rather than UTF-8 so that they can be
@@ -434,6 +448,15 @@ def read_lc(file_or_dir, format='ascii', **kwargs):
     format : {'ascii', 'json', 'salt2', 'salt2-old'}, optional
         Format of file. Default is 'ascii'. 'salt2' is the new format available
         in snfit version >= 2.3.0.
+    read_covmat : bool, optional
+        **[salt2 only]** If True, and if a ``COVMAT`` keyword is present in
+        header, read the covariance matrix from the filename specified
+        by ``COVMAT`` (assumed to be in the same directory as the lightcurve
+        file) and include it as a column named ``Fluxcov`` in the returned
+        table. Default is False.
+
+        *New in version 1.5.0*
+
     delim : str, optional
         **[ascii only]** Used to split entries on a line. Default is `None`.
         Extra whitespace is ignored.
@@ -658,8 +681,8 @@ def _write_snana(f, data, meta, **kwargs):
 def _write_json(f, data, meta, **kwargs):
 
     # Build a dictionary of pure-python objects
-    output = odict([('meta', meta),
-                    ('data', odict())])
+    output = OrderedDict([('meta', meta),
+                          ('data', OrderedDict())])
     for key in data.dtype.names:
         output['data'][key] = data[key].tolist()
     json.dump(output, f)
@@ -710,7 +733,7 @@ def write_lc(data, fname, format='ascii', **kwargs):
         meta = data.meta
         data = np.asarray(data)
     else:
-        meta = odict()
+        meta = OrderedDict()
         if not isinstance(data, np.ndarray):
             data = dict_to_array(data)
     with open(fname, 'w') as f: