Skip to content
spec1d.py 36.5 KiB
Newer Older
Shuai Feng's avatar
Shuai Feng committed
        lam  = hdulist[1].data['Wave']
        flux = hdulist[2].data
        par  = hdulist[3].data
        
        lam_range_temp = np.array([3500, 12000])
        TemNew, log_lam_temp = log_rebin(lam_range_temp, flux[1, :], velscale = velscale)[:2]
        
        # FWHM of XLS templates
        Temp_wave = np.exp(log_lam_temp)
        Temp_FWHM = np.zeros(len(log_lam_temp))
        Temp_FWHM[(Temp_wave < 5330)] = 13 * 2.35 / 3e5 * Temp_wave[(Temp_wave < 5330)]   # sigma = 13km/s at lambda <5330
        Temp_FWHM[(Temp_wave >= 5330) & (Temp_wave < 9440)] = 11 * 2.35 / 3e5 * Temp_wave[(Temp_wave >= 5330) & (Temp_wave < 9440)]
        # sigma = 13km/s at 5330 < lambda < 9440
        Temp_FWHM[(Temp_wave >= 9440)] = 16 * 2.35 / 3e5 * Temp_wave[(Temp_wave >= 9440)] # sigma=16km/s at lambda > 9440
        
        LSF = Temp_FWHM
        
        FWHM_eff = Temp_FWHM.copy()   # combined FWHM from stellar library and instrument(input)
        if np.isscalar(FWHM_inst):
            FWHM_eff[Temp_FWHM < FWHM_inst] = FWHM_inst
            LSF[Temp_FWHM < FWHM_inst]      = FWHM_inst
        else:
            FWHM_eff[Temp_FWHM < FWHM_inst] = FWHM_inst[Temp_FWHM < FWHM_inst]
            LSF[Temp_FWHM < FWHM_inst]      = FWHM_inst[Temp_FWHM < FWHM_inst]
        FWHM_dif  = np.sqrt(FWHM_eff ** 2 - Temp_FWHM ** 2)
        sigma_dif = FWHM_dif / 2.355 / (lam[1] - lam[0])  # Sigma difference in pixels

        temp = np.empty((TemNew.size, par.size))
        for i in range(par.size):
            temp0 = log_rebin(lam_range_temp, flux[i, :], velscale=velscale)[0]
            if np.isscalar(FWHM_dif):
                temp1 = ndimage.gaussian_filter1d(temp0, sigma_dif)
            else:
                temp1 = gaussian_filter1d(temp0, sigma_dif)             # convolution with variable sigma
            tempNew = temp1 / np.mean(temp1)
            temp[:, i] = tempNew
            

        self.templates    = temp
        self.log_lam_temp = log_lam_temp
        self.teff_grid = par['Teff']
        self.feh_grid  = par['FeH']
        self.logg_grid = par['logg']
        self.LSF       = Temp_FWHM
        self.velscale  = velscale
        
class SingleStar():

    """
    Class of single stelar spectrum

    Parameters
    ----------
    config : class
        Class of configuration
    template : class
        Class of single stellar population template
    mag : float, optional
        Magnitude in SDSS r-band, by default 15.0
    Teff : float, optional
        Effective tempreture, by default 10000.0K
    FeH : float, optional
        Metallicity of stellar, by default 0.0
    vel : float, optional
        Line of sight velocity, by default 100.0km/s
    Ebv : float, optional
        Dust extinction, by default 0.1
    """

    def __init__(self, config, template, mag = 15.0, teff = 10000.0, feh = 0.0, vel = 100.0, ebv = 0.0):
    
        StarTemp = template.templates
        
        # Select metal bins
        idx_FeH = (np.abs(template.feh_grid - feh) < 0.5)
        tpls = StarTemp[:, idx_FeH]
        
        # Select Teff bins
        Teff_FeH = template.teff_grid[idx_FeH]
        minloc   = np.argmin(abs(teff - Teff_FeH))
        starspec = tpls[:, minloc]
        
        wave = np.exp(template.log_lam_temp)
        
        # Dust Reddening
        if np.isscalar(ebv):
            starspec = reddening(wave, starspec, ebv = ebv)
            
        # Redshift
        redshift = vel / 3e5
        wave_r = wave * (1 + redshift)
        
        flux = np.interp(config.wave, wave_r, starspec)
        
        # Calibration
        if np.isscalar(mag):
            flux = calibrate(config.wave, flux, mag, filtername='SLOAN_SDSS.r')
        
        # Convert to input wavelength
        self.wave = config.wave
        self.flux = flux