Merge pull request #32 from mtakahiro/niriss_spectrum

Niriss spectrum

docs/parameters.rst

@@ -81,7 +81,7 @@ An example config file can be generated by executing the following command (TBD)
      - Comma separated list for spectral files. Path is relative to DIR_EXTR.
    * - MORP
      - str
-     - Profile shape, if grism spectra are provided. Will be used to convolve templates. "moffat", "gauss", or none.
+     - "moffat", "gauss", or none. Profile shape, if grism spectra are provided. Will be used to convolve templates.
    * - MORP_FILE
      - str
      - Ascii file for morphology parameters. 

gsf/function.py

@@ -5,6 +5,7 @@
 import pickle as cPickle
 import os
 import scipy.interpolate as interpolate
+import logging
 
 c = 3.e18 # A/s
 #d = 10**(73.6/2.5) # From [ergs/s/cm2/A] to [ergs/s/cm2/Hz]
@@ -17,6 +18,44 @@
 fLW = np.zeros(len(LW0), dtype='int') # flag.
 
 
+def str2bool(v):
+    '''
+    '''
+    if isinstance(v, bool):
+       return v
+    if v.lower() in ('yes', 'true', 't', 'y', '1'):
+        return True
+    elif v.lower() in ('no', 'false', 'f', 'n', '0'):
+        return False
+    else:
+        raise argparse.ArgumentTypeError('Boolean value expected.')
+
+
+def get_uvbeta(lm, flam, zbes, lam_blue=1650, lam_red=2300):
+    '''
+    Purpose
+    -------
+    get UV beta_lambda slope.
+
+    Parameters
+    ----------
+    lm : float array
+        in lambda
+    flam : float array
+        in flambda 
+    '''
+    con_uv = (lm/(1.+zbes)>lam_blue) & (lm/(1.+zbes)<lam_red)
+    #flam = fnutolam(lm,fl)
+    try:
+        fit_results = np.polyfit(np.log10(lm/(1.+zbes))[con_uv], np.log10(flam)[con_uv], 1) #, w=flam[con_uv])
+        beta = fit_results[0]
+        if np.isnan(beta):
+            beta = -99
+    except:
+        beta = -99
+    return beta
+
+
 def func_tmp(xint,eobs,fmodel):
     '''
     '''
@@ -523,6 +562,22 @@ def fnutolam(lam, fnu, m0set=25.0):
     return flam
 
 
+def delta(x,A):
+    '''
+    '''
+    yy = np.zeros(len(x),float)
+    iix = np.argmin(np.abs(x))
+    if len(x)%2 == 0:
+        logger= logging.getLogger( __name__ )
+        print(logger)
+        print('Input array has an even number.')
+        yy[iix] = A/2.
+        yy[iix+1] = A/2.
+    else:
+        yy[iix] = A
+    return yy
+
+
 def gauss(x,A,sig):
     return A * np.exp(-0.5*x**2/sig**2)
 

gsf/function_class.py

@@ -241,7 +241,7 @@ def open_spec_fits_dir(self, nage:int, nz:int, kk, Av00:float, zgal:float, A00:f
         return A00 * yyd_sort, xxd_sort
 
 
-    def get_template(self, lib, Amp:float = 1.0, T:float = 1.0, Av:float = 0.0, Z:float = 0.0, zgal:float = 1.0, f_bb:bool = False):
+    def get_template(self, lib, Amp:float = 1.0, T:float = 1.0, Av:float = 0.0, Z:float = 0.0, zgal:float = 1.0, f_bb:bool = False, fneb=False):
         '''
         Gets an element template given a set of parameters.
         Not necessarily the most efficient way, but easy to use.
@@ -282,7 +282,10 @@ def get_template(self, lib, Amp:float = 1.0, T:float = 1.0, Av:float = 0.0, Z:fl
         if T - self.age[nmodel] != 0:
             print('T=%.2f is not found in age library. T=%.2f is used.'%(T,self.age[nmodel]))
 
-        coln = int(2 + pp*len(self.ZZ)*len(self.AA) + NZ*len(self.AA) + nmodel)
+        if fneb:
+            coln = int(2 + pp*len(self.ZZ)*1 + NZ*1 + 0)
+        else:
+            coln = int(2 + pp*len(self.ZZ)*len(self.AA) + NZ*len(self.AA) + nmodel)
         nr = lib[:, 0]
         xx = lib[:, 1] # This is OBSERVED wavelength range at z=zgal
         yy = lib[:, coln]
@@ -776,12 +779,11 @@ def open_spec_fits(self, fall:int=0, orig:bool=False):
             hdu0 = self.MB.af['spec_full']
 
         DIR_TMP = self.DIR_TMP
-
         NZ = len(ZZ)
         NT = self.MB.ntau
         NA = self.MB.nage
         for zz,Z in enumerate(ZZ):
-            for tt,TT in enumerate(self.MB.tau):
+            for tt,TT in enumerate(self.MB.tau):                
                 for ss,TA in enumerate(self.MB.ageparam):
                     if zz == 0 and tt == 0 and ss == 0:
                         nr = hdu0['colnum']

gsf/function_igm.py

@@ -77,25 +77,28 @@ def get_sig_lya(lam_o, z_s, T=1e4, c=3e18):
 
 def get_nH(z):
 	'''
+	Purpose
+	-------
+	Get HI density by using Cen & Haiman 2000.
+
 	Returns
 	-------
-	HI density in IGM.
+	HI density in IGM, in cm^-3
 	'''
 	try:
 		nH = np.zeros(len(z),dtype='float')
 	except:
 		nH = 0
 
-	# From Cen & Haiman 2000
-	nH = 8.5e-5 * ((1.+z)/8)**3 # in cm^-3
+	nH = 8.5e-5 * ((1.+z)/8)**3
 	return nH
 
 
 def get_column(zin, cosmo, Mpc_cm=3.08568025e+24, z_r=6.0, delz=0.1):
 	'''
 	Returns
 	-------
-	HI column density of IGM at zin, in cm^-3.
+	HI column density of IGM at zin, in cm^-2.
 	'''
 	z = np.arange(z_r, zin, delz)
 	try:

gsf/gsf.py

@@ -109,7 +109,7 @@ def run_gsf_template(inputs, fplt=0, tau_lim=0.001, idman=None, nthin=1, delwave
 def run_gsf_all(parfile, fplt, cornerplot=True, f_Alog=True, idman=None, zman=None, f_label=True, f_symbol=True, 
     f_SFMS=True, f_fill=True, save_sed=True, figpdf=False, mmax=300, skip_sfh=False, f_fancyplot=False, 
     skip_zhist=False, tau_lim=0.001, tset_SFR_SED=0.1, f_shuffle=False, amp_shuffle=1e-2, Zini=None, 
-    nthin=1, delwave=1):
+    nthin=1, delwave=1, f_plot_resid=False):
     '''
     Purpose
     -------
@@ -234,7 +234,7 @@ def run_gsf_all(parfile, fplt, cornerplot=True, f_Alog=True, idman=None, zman=No
         plot_sed(MB, fil_path=MB.DIR_FILT,
         figpdf=figpdf, save_sed=save_sed, mmax=mmax,
         dust_model=MB.dust_model, DIR_TMP=MB.DIR_TMP, f_label=f_label, f_fill=f_fill, 
-        f_fancyplot=f_fancyplot, f_plot_resid=True)
+        f_fancyplot=f_fancyplot, f_plot_resid=f_plot_resid)
 
     '''
     if fplt == 4: