Merge branch 'master' into emsystem

docs/abscomp.rst

@@ -149,11 +149,24 @@ Misc
 I/O
 +++
 
-One may generate a *dict* of the key properties of the AbsSystem
+One may generate a *dict* of the key properties of the AbsComponent
 with the to_dict() method::
 
    cdict = component.to_dict()
 
+One may also wish to Voigt profile fit components with one
+of a number of software packages (e.g. ALIS, JoeBVP).  To
+generate an input file for JoeBVP use::
+
+    from linetools.isgm.io import write_joebvp_from_components
+    write_joebvp_from_components(component_list, specfile, 'output_file.ascii')
+
+Similarly, one can generate a list of components from an outputted
+JoeBVP file::
+
+    from linetools.isgm.io import read_joebvp_to_components
+    comp_list = read_joebvp_to_components('joebvp_file.out')
+
 
 Synthesize Components
 +++++++++++++++++++++

linetools/isgm/abscomponent.py

@@ -65,7 +65,7 @@ class AbsComponent(object):
         A comment, default is ``
     """
     @classmethod
-    def from_abslines(cls, abslines, stars=None, reliability='none', **kwargs):
+    def from_abslines(cls, abslines, stars=None, comment="", reliability='none', **kwargs):
         """Instantiate from a list of AbsLine objects
 
         Parameters
@@ -74,6 +74,9 @@ def from_abslines(cls, abslines, stars=None, reliability='none', **kwargs):
             List of AbsLine objects
         stars : str, optional
             Asterisks to append to the ion name (e.g. fine-structure, CII*)
+        comment : str
+            A comment associated to the resulting AbsComponent
+            Default is `""`
         reliability : str, optional
             Reliability of AbsComponent
                 'a' - reliable
@@ -91,7 +94,7 @@ def from_abslines(cls, abslines, stars=None, reliability='none', **kwargs):
         init_line = abslines[0]
         slf = cls( init_line.attrib['coord'], (init_line.data['Z'],init_line.data['ion']),
                    init_line.z, init_line.limits.vlim,
-                   Ej=init_line.data['Ej'], stars=stars, reliability=reliability)
+                   Ej=init_line.data['Ej'], stars=stars, reliability=reliability, comment=comment)
         slf._abslines.append(init_line)
         # Append with component checking
         if len(abslines) > 1:

linetools/isgm/io.py

@@ -13,8 +13,18 @@
 import numpy as np
 import warnings
 
+from astropy.table import Table
+from astropy import constants as const
+from astropy import units as u
+
 from linetools import utils as ltu
+from linetools.analysis.absline import linear_clm
 from linetools.isgm.abssystem import GenericAbsSystem
+from linetools.isgm.abscomponent import AbsComponent
+from linetools.spectralline import AbsLine
+from linetools.lists.linelist import LineList
+
+ckms = const.c.to('km/s').value
 
 def abssys_from_json(filename):
     """
@@ -41,3 +51,130 @@ def abssys_from_json(filename):
 
     # Return
     return abs_sys
+
+
+def read_joebvp_to_components(filename, coord, llist=None, specfile=None, chk_vel=False):
+    """ Generate a list of AbsComponent objects from a JoeB VP output file
+
+    Parameters
+    ----------
+    filename : str
+      joeB VP filename
+    coord : SkyCoord
+      QSO sightline
+    llist : LineList, optional
+      Used to construct AbsLine objects
+    specfile : str, optional
+    chk_vel : bool, optional
+      Demand that the velocities of a given ion all be the same
+
+    Returns
+    -------
+    comps : list
+      list of AbsComponent objects
+    """
+    # init
+    if llist is None:
+        llist = LineList('ISM')
+    comps = []
+    # Read
+    vp_data = Table.read(filename, format='ascii')
+
+    # Subset by zsys + trans
+    lbls = []
+    for izsys, itrans in zip(vp_data['zsys'], vp_data['trans']):
+        lbls.append('{:.6f}_{:s}'.format(izsys, itrans))
+    lbls = np.array(lbls)
+    ulbls = np.unique(lbls)
+
+    # Subset by nflag; Build components
+    for lbl in ulbls:
+        mt_lines = np.where(lbls == lbl)[0]
+        if chk_vel:
+            if len(np.unique(vp_data['vel'][mt_lines])) != 1:
+                pdb.set_trace()
+        z_fit = ltu.z_from_dv(vp_data['vel'][mt_lines[0]]*u.km/u.s, vp_data['zsys'][mt_lines[0]])
+        # Loop on abs lines
+        alines = []
+        for idx in mt_lines:
+            zlim = [vp_data['zsys'][idx] +
+                    vp_data[vkey][idx] * (1 + vp_data['zsys'][idx]) / ckms
+                    for vkey in ['vlim1', 'vlim2']]
+            absline = AbsLine(vp_data['restwave'][idx] * u.AA, z=z_fit,
+                              zlim=zlim, linelist=llist)
+            # Add measurements [JB -- Want to capture anything else??]
+            absline.attrib['coord'] = coord
+            absline.attrib['flag_N'] = 1
+            absline.attrib['logN'] = vp_data['col'][idx]
+            absline.attrib['sig_logN'] = vp_data['sigcol'][idx]
+            absline.attrib['b'] = vp_data['bval'][idx]
+            absline.attrib['sig_b'] = vp_data['sigbval'][idx]
+            absline.attrib['z'] = z_fit
+            absline.attrib['sig_z'] = ltu.dz_from_dv(vp_data['sigvel'][idx]*u.km/u.s, vp_data['z_comp'][idx])
+            if specfile is None:
+                absline.attrib['specfile'] = vp_data['specfile'][idx]
+            else:
+                absline.attrib['specfile'] = specfile
+            # Fill N, sig_N
+            _, _, = linear_clm(absline.attrib)
+            alines.append(absline)
+
+        # AbsComponent
+        stars = '*' * alines[0].ion_name.count('*')
+        if 'comment' in vp_data.keys():
+            comment = vp_data['comment'][mt_lines[0]]
+        else:
+            comment = ''
+        if 'rely' in vp_data.keys():
+            reliability = vp_data['rely'][mt_lines[0]]
+        else:
+            reliability = 'none'
+        abscomp = AbsComponent.from_abslines(alines, stars=stars, comment=comment, reliability=reliability)
+
+        # Add measurements [JB -- Want to capture anything else??]
+        abscomp.attrib = alines[0].attrib.copy()
+        # Remove undesired keys
+        for key in ['EW', 'sig_EW', 'flag_EW', 'N', 'sig_N']:
+            abscomp.attrib.pop(key)
+        # And more required
+        for key in ['flag_N', 'logN', 'sig_logN']:
+            setattr(abscomp, key, abscomp.attrib[key])
+        # Errors must be in first line!
+        assert abscomp.sig_logN > 0.
+        comps.append(abscomp)
+    # Finish
+    return comps
+
+
+def write_joebvp_from_components(comp_list, specfile, outfile):
+    """ From a given component list, it produces an
+    input file for JOEBVP (Voigt profile fitter).
+
+    Parameters
+    ----------
+    comp_list : list of AbsComponent
+        Input list of components to group
+    specfile : str
+        Name of the spectrum file associated to the components
+        in comp_list
+    outfile : str
+        Name of the output file
+
+    """
+    # Open new file to write out
+    f = open(outfile, 'w')
+
+    # Print header
+    s = 'specfile|restwave|zsys|col|bval|vel|nflag|bflag|vflag|vlim1|vlim2|wobs1|wobs2|z_comp|trans|rely|comment\n'
+    f.write(s)
+
+    # Components
+    for ii, comp in enumerate(comp_list):
+        flags = (ii+2,ii+2,ii+2)
+        try:
+            b_val = comp.attrib['b']
+        except KeyError:
+            b_val = 10*u.km/u.s
+        s = comp.repr_joebvp(specfile, flags=flags, b_default=b_val)  # still, b values from abslines take precedence if they exist
+        f.write(s)
+    f.close()

linetools/isgm/setup_package.py

@@ -1,4 +1,4 @@
 def get_package_data():
     # Installs the testing data files. Unable to get package_data
     # to deal with a directory hierarchy of files, so just explicitly list.
-    return {'linetools.isgm.tests': ['files/*.fits', 'files/*.dat', 'files/*.json', 'files/*.all', 'files/*.fits.gz', 'files/*.out']}
+    return {'linetools.isgm.tests': ['files/*.fits', 'files/*.dat', 'files/*.json', 'files/*.all', 'files/*.fits.gz', 'files/*.VP', 'files/*.out']}

linetools/isgm/tests/files/group_9.VP

@@ -0,0 +1,4 @@
+specfile|restwave|zsys|col|sigcol|bval|sigbval|vel|sigvel|nflag|bflag|vflag|vlim1|vlim2|wobs1|wobs2|pix1|pix2|z_comp|trans 
+TON580_cont.fits|1215.67|0.2895|12.789|0.201|25.558|5.802|33.939|8.036|3|3|3|-100.0|100.0|1567.20683282|1568.01784055|13707|13729|0.28965|H I 
+TON580_cont.fits|1025.7222|0.28978|12.738|0.212|8.697|5.126|-16.281|2.0|2|2|2|-100.0|100.0|1322.61207078|1323.29635897|6256|6279|0.28971|H I 
+TON580_cont.fits|1215.67|0.28978|12.738|0.0|8.697|0.0|-16.281|0.0|2|2|2|-100.0|100.0|1567.53925778|1568.35026551|13716|13738|0.28971|H I 

linetools/isgm/tests/test_io.py

@@ -0,0 +1,50 @@
+# Module to run tests on using AbsComponent
+
+from __future__ import print_function, absolute_import, division, unicode_literals
+
+# TEST_UNICODE_LITERALS
+
+import pytest
+import os
+
+from pkg_resources import resource_filename
+
+from astropy.coordinates import SkyCoord
+from astropy import units as u
+
+from linetools.isgm.io import read_joebvp_to_components
+from linetools.isgm.io import write_joebvp_from_components
+
+from .utils import mk_comp, compare_two_files
+
+def data_path(filename):
+    data_dir = os.path.join(os.path.dirname(__file__), 'files')
+    return os.path.join(data_dir, filename)
+
+
+def test_read_joebvp():
+    vp_file = data_path('group_9.VP')
+    # Load
+    icoord = SkyCoord(ra=12., dec=-12, unit='deg')
+    comps = read_joebvp_to_components(vp_file, icoord)
+    # Test
+    assert isinstance(comps, list)
+    assert len(comps) == 2
+    assert len(comps[1]._abslines) == 2
+    assert comps[1]._abslines[0].attrib['N'].value > 0.
+
+
+def test_complist_to_joebvp():
+    # will write a file in directory ./files/
+    abscomp, HIlines = mk_comp('HI', b=15*u.km/u.s, use_rand=False)
+    comp_list = [abscomp, abscomp]
+    write_joebvp_from_components(comp_list, 'test.fits', data_path('test_joebvp_repr.joebvp'))
+    # now read the output and compare to reference
+    compare_two_files(data_path('test_joebvp_repr.joebvp'),
+                      resource_filename('linetools', '/data/tests/test_joebvp_repr_reference.joebvp'))
+    # now add attribute to comp and compare again
+    abscomp.attrib['b'] = 15*u.km/u.s
+    write_joebvp_from_components(comp_list, 'test.fits', data_path('test_joebvp_repr.joebvp'))
+    compare_two_files(data_path('test_joebvp_repr.joebvp'),
+                      resource_filename('linetools', '/data/tests/test_joebvp_repr_reference.joebvp'))
+

linetools/isgm/tests/test_use_abscomp.py

@@ -16,10 +16,10 @@
 
 from linetools.isgm.abscomponent import AbsComponent
 from linetools.spectralline import AbsLine
-from linetools.spectra import io as lsio
-from linetools.analysis import absline as ltaa
 from linetools.isgm import utils as ltiu
 
+from .utils import mk_comp, mk_comptable
+
 import imp, os
 lt_path = imp.find_module('linetools')[1]
 
@@ -28,72 +28,12 @@
 # Set of Input lines
 
 
-def mk_comp(ctype,vlim=[-300.,300]*u.km/u.s,add_spec=False, use_rand=True,
-            add_trans=False, zcomp=2.92939, b=20*u.km/u.s):
-    # Read a spectrum Spec
-    if add_spec:
-        xspec = lsio.readspec(lt_path+'/spectra/tests/files/UM184_nF.fits')
-    else:
-        xspec = None
-    # AbsLines
-    if ctype == 'HI':
-        all_trans = ['HI 1215', 'HI 1025']
-    elif ctype == 'SiII':
-        all_trans = ['SiII 1260', 'SiII 1304', 'SiII 1526', 'SiII 1808']
-        if add_trans:
-            all_trans += ['SiII 1193']
-    elif ctype == 'SiII*':
-        all_trans = ['SiII* 1264', 'SiII* 1533']
-    abslines = []
-    for trans in all_trans:
-        iline = AbsLine(trans, z=zcomp)
-        if use_rand:
-            rnd = np.random.rand()
-        else:
-            rnd = 0.
-        iline.attrib['logN'] = 13.3 + rnd
-        iline.attrib['sig_logN'] = 0.15
-        iline.attrib['flag_N'] = 1
-        iline.attrib['b'] = b
-        iline.analy['spec'] = xspec
-        iline.limits.set(vlim)
-        _,_ = ltaa.linear_clm(iline.attrib)  # Loads N, sig_N
-        abslines.append(iline)
-    # Component
-    abscomp = AbsComponent.from_abslines(abslines)
-    return abscomp, abslines
-
-
-def mk_comptable():
-    tab = Table()
-    tab['ion_name'] = ['HI', 'HI', 'CIV', 'SiII', 'OVI']
-    tab['Z'] = [1,1,4,14,8]
-    tab['ion'] = [1,1,4,2,6]
-    tab['z_comp'] = [0.05, 0.0999, 0.1, 0.1001, 0.6]
-    tab['RA'] = [100.0] * len(tab) * u.deg
-    tab['DEC'] = [-0.8] * len(tab) * u.deg
-    tab['vmin'] = [-50.] * len(tab) * u.km/u.s
-    tab['vmax'] = [100.] * len(tab) * u.km/u.s
-    tab['Ej'] = [0.] *len(tab) / u.cm
-    return tab
-
 
 def data_path(filename):
     data_dir = os.path.join(os.path.dirname(__file__), 'files')
     return os.path.join(data_dir, filename)
 
 
-def compare_two_files(file1, file2):
-    f1 = open(file1, 'r')
-    f2 = open(file2, 'r')
-    lines1 = f1.readlines()
-    lines2 = f2.readlines()
-    for l1,l2 in zip(lines1,lines2):
-        assert l1 == l2
-    f1.close()
-    f2.close()
-
-
 
 
 def test_get_components_at_z():
@@ -210,18 +150,6 @@ def test_repr_joebvp():
                 'test.fits|1025.72220|2.92939000|13.3000|15.0000|0.|2|2|2|-300.0000|300.0000|4026.43132|4034.49783|2.92939000|HI|none|Something\n'
 
 
-def test_complist_to_joebvp():
-    # will write a file in directory ./files/
-    abscomp, HIlines = mk_comp('HI', b=15*u.km/u.s, use_rand=False)
-    comp_list = [abscomp, abscomp]
-    ltiu.joebvp_from_components(comp_list, 'test.fits', data_path('test_joebvp_repr.joebvp'))
-    # now read the output and compare to reference
-    compare_two_files(data_path('test_joebvp_repr.joebvp'), lt_path + '/data/tests/test_joebvp_repr_reference.joebvp')
-    # now add attribute to comp and compare again
-    abscomp.attrib['b'] = 15*u.km/u.s
-    ltiu.joebvp_from_components(comp_list, 'test.fits', data_path('test_joebvp_repr.joebvp'))
-    compare_two_files(data_path('test_joebvp_repr.joebvp'), lt_path + '/data/tests/test_joebvp_repr_reference.joebvp')
-
 
 def test_get_wvobs_chunks():
     abscomp, HIlines = mk_comp('HI', zcomp=0., vlim=[0,10]*u.km/u.s)