Merge branch 'master' into llist_addlines

CHANGES.rst

@@ -3,7 +3,9 @@
 
 Updates
 .......
-
+- Added extra attributes to AbsComponents
+- Added script to get HST/COS life-time position from date
+- Added Cashman+2017 catalog in LineList
 
 Bug fixes
 .........

linetools/isgm/abscomponent.py

@@ -30,6 +30,13 @@
 # Global import for speed
 c_kms = const.c.to('km/s').value
 
+# Global initial attrib dict (speeds up performance b/c astropy.Quantity issue?)
+init_attrib = {'N': 0./u.cm**2, 'sig_N': 0./u.cm**2, 'flag_N': 0, # Column    ## NOT ENOUGH SPEED-UP
+              'logN': 0., 'sig_logN': 0.,
+              'b': 0.*u.km/u.s, 'sig_b': 0.*u.km/u.s,  # Doppler
+              'vel': 0*u.km/u.s, 'sig_vel': 0*u.km/u.s
+              }
+
 # Class for Components
 class AbsComponent(object):
     """
@@ -55,6 +62,8 @@ class AbsComponent(object):
         Atomic mass -- used to distinguish isotopes
     Ej : Quantity
         Energy of lower level (1/cm)
+    attrib : dict
+        Absorption properties (e.g. column, Doppler b, centroid)
     reliability : str, optional
         Reliability of AbsComponent
             'a' - reliable
@@ -65,7 +74,8 @@ class AbsComponent(object):
         A comment, default is ``
     """
     @classmethod
-    def from_abslines(cls, abslines, stars=None, comment="", reliability='none', **kwargs):
+    def from_abslines(cls, abslines, stars=None, comment="", reliability='none',
+                      tol=0.01, chk_meas=False, verbose=True, **kwargs):
         """Instantiate from a list of AbsLine objects
 
         Parameters
@@ -83,7 +93,15 @@ def from_abslines(cls, abslines, stars=None, comment="", reliability='none', **k
                 'b' - possible
                 'c' - uncertain
                 'none' - not defined (default)
+        chk_meas : bool, optional
+            If true, require that measurements for lines composing a component
+            have matching values before setting component attribs.
+            Otherwise, throw a warning
+        tol : float, optional
+            Fractional tolerance for line attributes (N,b) to match one another
+
         """
+        from linetools.analysis import absline as ltaa
         # Check
         if not isinstance(abslines, list):
             raise IOError("Need a list of AbsLine objects")
@@ -100,6 +118,58 @@ def from_abslines(cls, abslines, stars=None, comment="", reliability='none', **k
         if len(abslines) > 1:
             for absline in abslines[1:]:
                 slf.add_absline(absline, **kwargs)
+
+        ### Add attribs to component
+        # First check that the measurements for all the lines match
+        # Grab them all
+        cols = np.array([al.attrib['N'].value for al in abslines])
+        colerrs = np.array([al.attrib['sig_N'].value for al in abslines])
+        bs = np.array([al.attrib['b'].value for al in abslines])
+        berrs = np.array([al.attrib['sig_b'].value for al in abslines])
+        vels = np.array([al.attrib['vel'].value for al in abslines])
+        velerrs = np.array([al.attrib['sig_vel'].value for al in abslines])
+        medcol = np.median(cols)
+        medcolerr = np.median(colerrs)
+        medb = np.median(bs)
+        medberr = np.median(berrs)
+        medvel = np.median(vels)
+        medvelerr = np.median(velerrs)
+        # Perform checks on measurements
+        if medcol == 0.:
+            colcrit = np.all((cols) == 0.)
+        else:
+            colcrit = all(np.abs(cols - medcol) / medcol < tol) is True
+        if medb == 0.:
+            bcrit = np.all((bs) == 0.)
+        else:
+            bcrit = all(np.abs(bs - medb) / medb < tol) is True
+        # Set attribs
+        slf.attrib['N'] = medcol / u.cm ** 2
+        slf.attrib['sig_N'] = medcolerr / u.cm ** 2
+        slf.attrib['b'] = medb * u.km / u.s
+        slf.attrib['sig_b'] = medberr * u.km / u.s
+        slf.attrib['vel'] = medvel * u.km / u.s
+        slf.attrib['sig_vel'] = medvelerr * u.km / u.s
+        logN, sig_logN = ltaa.log_clm(slf.attrib)
+        slf.attrib['logN'] = logN
+        slf.attrib['sig_logN'] = sig_logN
+        # Stop or throw warnings
+        if (bcrit&colcrit):
+            pass
+        elif chk_meas:
+            raise ValueError('The line measurements for the lines in this '
+                             'component are not consistent with one another.')
+        else:
+            if verbose:
+                warnings.warn('The line measurements for the lines in this component'
+                          ' may not be consistent with one another.')
+                if bcrit:
+                    print('Problem lies in the column density values')
+                elif colcrit:
+                    print('Problem lies in the b values')
+                else:
+                    print('Problems lie in the column densities and b values')
+
         # Return
         return slf
 
@@ -158,7 +228,14 @@ def from_dict(cls, idict, coord=None, **kwargs):
                   Ej=idict['Ej']/u.cm, A=idict['A'],
                   Ntup = tuple([idict[key] for key in ['flag_N', 'logN', 'sig_logN']]),
                   comment=idict['comment'], name=idict['Name'], reliability=reliability)
-
+        # Attributes
+        if 'attrib' in idict.keys():
+            attrkeys = idict['attrib'].keys()
+            for ak in attrkeys:
+                if isinstance(idict['attrib'][ak],dict):
+                    slf.attrib[ak] = ltu.convert_quantity_in_dict(idict['attrib'][ak])
+                else:
+                    slf.attrib[ak] = idict['attrib'][ak]
         # Add lines
         for key in idict['lines'].keys():
             iline = AbsLine.from_dict(idict['lines'][key], coord=coord, **kwargs)
@@ -251,7 +328,7 @@ def __init__(self, radec, Zion, zcomp, vlim, Ej=0./u.cm, A=None,
         self.reliability = reliability
 
         # Potential for attributes
-        self.attrib = dict()
+        self.attrib = init_attrib.copy()
 
         # Other
         self._abslines = []
@@ -264,6 +341,22 @@ def vlim(self):
     def zcomp(self):
         return self.limits.z
 
+    @property
+    def b(self):
+        return self.attrib['b']
+
+    @property
+    def sig_b(self):
+        return self.attrib['sig_b']
+
+    @property
+    def vel(self):  #velocity centroid
+        return self.attrib['vel']
+
+    @property
+    def sig_vel(self):  #velocity centroid
+        return self.attrib['sig_vel']
+
     def add_absline(self, absline, tol=0.1*u.arcsec, chk_vel=True,
                     chk_sep=True, vtoler=1., **kwargs):
         """Add an AbsLine object to the component if it satisfies
@@ -828,7 +921,8 @@ def to_dict(self):
                      Name=self.name,
                      RA=self.coord.icrs.ra.value, DEC=self.coord.icrs.dec.value,
                      A=self.A, Ej=self.Ej.to('1/cm').value, comment=self.comment,
-                     flag_N=self.flag_N, logN=self.logN, sig_logN=self.sig_logN)
+                     flag_N=self.flag_N, logN=self.logN, sig_logN=self.sig_logN,
+                     attrib=self.attrib)
         cdict['class'] = self.__class__.__name__
         # AbsLines
         cdict['lines'] = {}

linetools/isgm/abssystem.py

@@ -101,7 +101,7 @@ def from_components(cls, components, vlim=None, NHI=None):
           List of AbsComponent objects
         vlim : list, optional
           Velocity limits for the system
-          If not set, the first components sets vlim
+          If not set, the first component sets vlim
         NHI : float, optional
           Set the NHI value of the system.  If not set,
           the method sums the NHI values of all the HI
@@ -494,6 +494,20 @@ def update_component_colm(self, **kwargs):
         for comp in self._components:
             comp.synthesize_colm(**kwargs)
 
+    def update_component_vel(self):
+        """Change the velocities of each component to rest frame of zsys
+        """
+        from linetools.analysis.zlimits import zLimits
+
+        for i,comp in enumerate(self._components):
+            dv = ltu.dv_from_z(comp.zcomp, self.zabs)
+            zmin = ltu.z_from_dv(comp.limits.vlim[0]+dv,self.zabs)
+            zmax = ltu.z_from_dv(comp.limits.vlim[1]+dv,self.zabs)
+            newzlim = zLimits(self.zabs, (zmin, zmax))
+            comp.limits = newzlim
+            comp.attrib['vel'] = comp.attrib['vel'] + dv
+
+
     def stack_plot(self, pvlim=None, **kwargs):
         """Show a stack plot of the system, if spec are loaded
         Assumes the data are normalized.
@@ -564,6 +578,7 @@ def get_vmnx(components):
         else:
             components = self._components
             self.vlim = get_vmnx(components)  # Using system z
+
 
     def write_json(self, outfil=None, overwrite=True):
         """ Generate a JSON file from the system

linetools/isgm/io.py

@@ -107,8 +107,8 @@ def read_joebvp_to_components(filename, coord, llist=None, specfile=None, chk_ve
             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['b'] = vp_data['bval'][idx] * u.km/u.s
+            absline.attrib['sig_b'] = vp_data['sigbval'][idx] * u.km/u.s
             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:
@@ -140,7 +140,8 @@ def read_joebvp_to_components(filename, coord, llist=None, specfile=None, chk_ve
         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.
+        assert abscomp.sig_logN > 0., "AbsComponent has sig_logN=0 {}".format(abscomp)
+
         comps.append(abscomp)
     # Finish
     return comps

linetools/isgm/tests/test_init_abscomp.py

@@ -32,6 +32,15 @@ def test_init_failures():
         AbsComponent.from_abslines(['blah'])
     with pytest.raises(IOError):
         AbsComponent.from_component('blah')
+    # Inconsistent abslines
+    lya = AbsLine(1215.670*u.AA,z=2.92939)
+    lya.attrib['N'] = 1e12 / u.cm**2
+    lya.attrib['b'] = 30 * u.km/u.s
+    lyb = AbsLine('HI 1025',z=2.92939)
+    lyb.attrib['N'] = 3e12 / u.cm**2
+    lyb.attrib['b'] = 30 * u.km/u.s
+    with pytest.raises(ValueError):
+        AbsComponent.from_abslines([lya,lyb], chk_meas=True)
 
 
 def test_init_single_absline():
@@ -42,7 +51,6 @@ def test_init_single_absline():
     # Test
     assert abscomp.Zion[0] == 1
     np.testing.assert_allclose(abscomp.zcomp,2.92939)
-    print(abscomp)
 
 
 def test_copy():

linetools/isgm/utils.py

@@ -361,7 +361,7 @@ def table_from_complist(complist):
     """
     key_order = ['RA', 'DEC', 'name', 'z_comp', 'sig_z', 'Z', 'ion', 'Ej',
                  'vmin', 'vmax','ion_name', 'flag_N', 'logN', 'sig_logN',
-                 'b','sig_b', 'specfile']
+                 'b','sig_b', 'vel', 'sig_vel','specfile']
 
     tab = Table()
 
@@ -376,7 +376,8 @@ def table_from_complist(complist):
     tab['ion'] = [icomp.Zion[1] for icomp in complist]
 
     # . attributes (required ones)
-    for attrib in ['zcomp', 'Ej', 'flag_N', 'logN', 'sig_logN']:
+    for attrib in ['zcomp', 'Ej', 'flag_N', 'logN', 'sig_logN',
+                   'b','sig_b','vel','sig_vel']:
         values = [getattr(icomp,attrib) for icomp in complist]
         if isinstance(values[0], u.Quantity):
             values = u.Quantity(values)
@@ -394,7 +395,7 @@ def table_from_complist(complist):
     tab['ion_name'] = ion_names
 
     # attrib dict
-    for attrib in ['sig_z', 'b', 'sig_b', 'specfile']:
+    for attrib in ['sig_z', 'b', 'sig_b', 'vel', 'sig_vel', 'specfile']:
         try:
             values = [icomp.attrib[attrib] for icomp in complist]
         except KeyError:

linetools/lists/linelist.py

@@ -399,10 +399,11 @@ def make_extra_table(self, abundance_type='solar', ion_correction='none',
         # Set QM strength as MaskedColumn (self._data['f'] is MaskedColumn)
         qm_strength = self._data['f'] * self._data['wrest']
         qm_strength.name = 'qm_strength'
-        self._extra_table['log(w*f)'] = np.log10(qm_strength)
+        self._extra_table['log(w*f)'] = np.zeros_like(qm_strength)
+        gdqm = qm_strength > 0.
+        self._extra_table['log(w*f)'][gdqm] = np.log10(qm_strength[gdqm])
         # mask out potential nans
-        cond = np.isnan(self._extra_table['log(w*f)'])
-        self._extra_table['log(w*f)'].mask = np.where(cond, True, self._extra_table['log(w*f)'].mask)
+        self._extra_table['log(w*f)'].mask = ~gdqm
 
         # Set Abundance
         if abundance_type not in ['none', 'solar']:

linetools/scripts/lt_get_COS_LP.py

@@ -0,0 +1,35 @@
+#!/usr/bin/env python
+from __future__ import (print_function, absolute_import, division, unicode_literals)
+
+
+try:
+    ustr = unicode
+except NameError:
+    ustr = str
+
+"""Obtain the COS life-time position given a date.
+"""
+def parser(options=None):
+    import argparse
+    # Parse arguments
+    parser = argparse.ArgumentParser(
+        description='Obtain the HST/COS life-time position (LP) given a date.')
+    parser.add_argument("date", type=str, default=None, help='Enter a date, e.g. "2017-10-09 00:00:00", "2017-10-09", etc.')
+    if options is None:
+        args = parser.parse_args()
+    else:
+        args = parser.parse_args(options)
+    return args
+
+def main(args=None):
+    from linetools.spectra import utils as lsu
+
+    pargs = parser(options=args)
+    date = pargs.date
+
+    # return
+    lp = lsu.get_COS_LP_from_date(date)
+    print(lp)
+
+if __name__ == '__main__':
+    main()

linetools/scripts/tests/test_scripts.py

@@ -8,6 +8,7 @@
 from linetools.scripts import lt_line
 from linetools.scripts import lt_radec
 from linetools.scripts import lt_solabnd
+from linetools.scripts import lt_get_COS_LP
 
 
 def test_lt_absline():
@@ -33,3 +34,6 @@ def test_lt_solabnd():
     lt_solabnd.main(['Fe'])
     lt_solabnd.main(['-a'])
     lt_solabnd.main(['-a', '--sortZ'])
+
+def test_lt_get_COS_LP():
+    lt_get_COS_LP.main(["2017-10-01"])

linetools/spectralline.py

@@ -38,7 +38,8 @@
 
 abs_attrib = {'N': 0./u.cm**2, 'sig_N': 0./u.cm**2, 'flag_N': 0, # Column    ## NOT ENOUGH SPEED-UP
               'logN': 0., 'sig_logN': 0.,
-                    'b': 0.*u.km/u.s, 'sig_b': 0.*u.km/u.s  # Doppler
+                    'b': 0.*u.km/u.s, 'sig_b': 0.*u.km/u.s,  # Doppler
+                    'vel': 0.*u.km/u.s, 'sig_vel': 0.*u.km/u.s  # Doppler
                     }
 
 emiss_attrib = {'flux': 0.*u.erg/u.s, 'sig_flux': 0.*u.erg/u.s, 'flag_flux': 0,

setup.py

@@ -106,6 +106,8 @@ def recursive_glob(basedir, pattern):
     'lt_plot = linetools.scripts.lt_plot:main',
     'lt_radec = linetools.scripts.lt_radec:main',
     'lt_solabnd = linetools.scripts.lt_solabnd:main',
+    'lt_get_COS_LP = linetools.scripts.lt_get_COS_LP:main',
+
 #    'astropy-package-template-example = packagename.example_mod:main',
 ]