Merge pull request #479 from desihub/same_speedoflight

import speed of light from scipy

py/desisim/lya_mock_p1d.py

@@ -1,5 +1,11 @@
 import numpy as np
 
+try:
+    from scipy import constants
+    C_LIGHT = constants.c/1000.0
+except TypeError: # This can happen during documentation builds.
+    C_LIGHT = 299792458.0/1000.0
+
 # code to make mock Lya spectra following McDonald et al. (2006)
 # copied from c++ code in Cosmology/LNLyaF
 
@@ -64,7 +70,7 @@ def get_redshifts(self):
         """Get redshifts for each cell in the array (centered at z_c)."""
         N = self.N
         L_kms = N * self.dv_kms
-        c_kms = 2.998e5
+        c_kms = C_LIGHT
         if (L_kms > 4 * c_kms):
             print('Array is too long, approximations break down.')
             raise SystemExit

py/desisim/qso_template/desi_qso_templ.py

@@ -17,6 +17,12 @@
 
 from astropy.io import fits
 
+try:
+    from scipy import constants
+    C_LIGHT = constants.c/1000.0
+except TypeError: # This can happen during documentation builds.
+    C_LIGHT = 299792458.0/1000.0
+
 from desisim.qso_template import fit_boss_qsos as fbq
 from desiutil.stats import perc
 import desisim.io
@@ -368,7 +374,7 @@ def desi_qso_templates(z_wind=0.2, zmnx=(0.4,4.), outfil=None, N_perz=500,
 
     # Rebin
     if rebin_wave is None:
-        light = 2.99792458e5        # [km/s]
+        light = C_LIGHT        # [km/s]
         velpixsize = 10.            # [km/s]
         pixsize = velpixsize/light/np.log(10) # [pixel size in log-10 A]
         minwave = np.log10(wvmnx[0])          # minimum wavelength [log10-A]

py/desisim/spec_qa/redshifts.py

@@ -20,6 +20,12 @@
 from astropy.io import fits
 from astropy.table import Table, vstack, hstack, MaskedColumn, join
 
+try:
+    from scipy import constants
+    C_LIGHT = constants.c/1000.0
+except TypeError: # This can happen during documentation builds.
+    C_LIGHT = 299792458.0/1000.0
+
 import desispec.io
 from .utils import elg_flux_lim, get_sty_otype, catastrophic_dv, match_otype
 
@@ -530,7 +536,7 @@ def criteria(simz_tab, objtype=None, dvlimit=None):
         else:
             dv = dvlimit
         dz = calc_dz(simz_tab[omask]) # dz/1+z
-        cat = np.where(np.abs(dz)*3e5 > dv)[0]
+        cat = np.where(np.abs(dz)*C_LIGHT > dv)[0]
         dv_mask[omask[cat]] = False
     # Return
     return objtype_mask, z_mask, survey_mask, dv_mask, zwarn_mask
@@ -716,7 +722,7 @@ def obj_fig(simz_tab, objtype, summ_stats, outfile=None):
             ax.set_xlabel(lbl)
             ax.set_ylabel(ylbl)
             ax.set_xlim(xmin,xmax)
-            v_ylim = ylim * 3e5  # redshift to km/s
+            v_ylim = ylim * C_LIGHT  # redshift to km/s
             ax.set_ylim(-v_ylim+yoff, v_ylim+yoff)
 
             # Points
@@ -730,7 +736,7 @@ def obj_fig(simz_tab, objtype, summ_stats, outfile=None):
                 xbins = np.linspace(xmin, xmax, 20)
             ybins = np.linspace(-v_ylim+yoff, v_ylim+yoff, 40) # km/s
             #import pdb; pdb.set_trace()
-            counts, xedges, yedges = np.histogram2d(xval, yval * 3e5, bins=(xbins, ybins))
+            counts, xedges, yedges = np.histogram2d(xval, yval * C_LIGHT, bins=(xbins, ybins))
             max_c = np.max(counts)
             #if kk == 3:
             ax.pcolormesh(xedges, yedges, counts.transpose(), cmap=cm, vmin=0, vmax=max_c/5.)
@@ -1061,7 +1067,7 @@ def dz_summ(simz_tab, outfile=None, pdict=None, min_count=20):
 
             # Simple stats
             ok = survey['ZWARN'] == 0
-            dv = calc_dz(survey)*3e5 # dz/1+z
+            dv = calc_dz(survey)*C_LIGHT # dz/1+z
             bad = dv > catastrophic_dv(otype)
             #if i==2:
             #    pdb.set_trace()

py/desisim/templates.py

@@ -14,7 +14,11 @@
 from desiutil.log import get_logger, DEBUG
 from desisim.io import empty_metatable
 
-LIGHT = 2.99792458E5  #- speed of light in km/s
+try:
+    from scipy import constants
+    C_LIGHT = constants.c/1000.0
+except TypeError: # This can happen during documentation builds.
+    C_LIGHT = 299792458.0/1000.0
 
 def _check_input_meta(input_meta, ignore_templateid=False):
     log = get_logger()
@@ -98,7 +102,7 @@ def __init__(self, minwave=3650.0, maxwave=7075.0, cdelt_kms=20.0, log10wave=Non
 
         # Build a wavelength array if one is not given.
         if log10wave is None:
-            cdelt = cdelt_kms/LIGHT/np.log(10) # pixel size [log-10 A]
+            cdelt = cdelt_kms/C_LIGHT/np.log(10) # pixel size [log-10 A]
             npix = int(round((np.log10(maxwave)-np.log10(minwave))/cdelt))+1
             self.log10wave = np.linspace(np.log10(minwave), np.log10(maxwave), npix)
         else:
@@ -298,7 +302,7 @@ def spectrum(self, oiiihbeta=None, oiihbeta=None, niihbeta=None,
                 line['flux'][ii] = hbetaflux*line['ratio'][ii]
 
         # Finally build the emission-line spectrum
-        log10sigma = linesigma /LIGHT / np.log(10) # line-width [log-10 Angstrom]
+        log10sigma = linesigma /C_LIGHT / np.log(10) # line-width [log-10 Angstrom]
         emspec = np.zeros_like(self.log10wave)
 
         loglinewave = np.log10(line['wave'])
@@ -472,7 +476,7 @@ def _blurmatrix(self, vdisp, log=None):
 
         blurmatrix = dict()
         for uvv in uvdisp:
-            sigma = 1.0 + (self.basewave * uvv / LIGHT)
+            sigma = 1.0 + (self.basewave * uvv / C_LIGHT)
             blurmatrix[uvv] = pxs.gauss_blur_matrix(self.pixbound, sigma).astype('f4')
 
         return blurmatrix
@@ -1415,7 +1419,7 @@ def make_star_templates(self, nmodel=100, vrad_meansig=(0.0, 200.0),
                     else:
                         vrad = np.repeat(vrad_meansig[0], nmodel)
 
-                    redshift = np.array(vrad) / LIGHT
+                    redshift = np.array(vrad) / C_LIGHT
 
                 if mag is None:
                     mag = rand.uniform(magrange[0], magrange[1], nmodel).astype('f4')