Python FlatLambdaCDM.luminosity_distance方法代码示例

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
def test_approx_DL():
    for z in np.linspace(0.01, 4, num=10):
        z = 2.
        v1 = approx_DL()(z)
        cosmo = FlatLambdaCDM(H0=70, Om0=0.3, Ob0=None)
        v2 = cosmo.luminosity_distance(z).value
        assert abs(v1/v2 - 1) < 0.01

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
def absolute_magnitude_lim(z, app_mag_lim, cosmo=None):
    """
    give the absolute magnitude limit as a function of redshift for a flux-limited survey.
    
    Parameters
    ----------
    z: array_like
        redshift
    
    app_mag_lim: float
       apparent magnitude limit
    
    cosmo: cosmology object
    
    Returns
    -------
    M,z: np.array, np.array
        absolute magnitude in mag+5loh(h) units
    """
    if cosmo==None:
        from astropy.cosmology import FlatLambdaCDM
        cosmo = FlatLambdaCDM(H0=70, Om0=0.3)
        print('Warning, no cosmology specified, using default:',cosmo)
    
    d_L = cosmo.luminosity_distance(z).value
    M = apparent_to_absolute_magnitude(app_mag_lim, d_L)
    
    return M-5.0*np.log10(cosmo.h)

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
def get_lum(rmid, zfinal_dict):
    flux = get_raw_lum(rmid)
    z = zfinal_dict[int(rmid)]
    cosmo = FlatLambdaCDM(H0=70.0, Om0=0.3)
    dl_MPC = cosmo.luminosity_distance(z)
    dl_cm = float(dl_MPC.value * 3.085677581 * (10.0 ** 24.0))
    lum = 4.0 * 3.1415926 * dl_cm * dl_cm * flux[0] * 10. ** (0. - 17.) * 5100.0
    return [lum, flux[1] * lum]

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
def angular_diameter(size,redshift): ## size must be in kpc

    cosmo = FlatLambdaCDM(H0=70, Om0=0.3)  ## defining the cosmology  => see the documentation and details of the parameters and their values
    lumdist=cosmo.luminosity_distance(redshift)  ## value in kpc

    zang=size*(1+redshift)**2/(1000.*lumdist)*180./np.pi*3600. ## arcsec

    return zang

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
def simple_flux_from_greybody(lambdavector, Trf = None, b = None, Lrf = None, zin = None, ngal = None):
	''' 
	Return flux densities at any wavelength of interest (in the range 1-10000 micron),
	assuming a galaxy (at given redshift) graybody spectral energy distribution (SED),
	with a power law replacing the Wien part of the spectrum to account for the
	variability of dust temperatures within the galaxy. The two different functional
	forms are stitched together by imposing that the two functions and their first
	derivatives coincide. The code contains the nitty-gritty details explicitly.
	Cosmology assumed: H0=70.5, Omega_M=0.274, Omega_L=0.726 (Hinshaw et al. 2009)

	Inputs:
	alphain = spectral index of the power law replacing the Wien part of the spectrum, to account for the variability of dust temperatures within a galaxy [default = 2; see Blain 1999 and Blain et al. 2003]
	betain = spectral index of the emissivity law for the graybody [default = 2; see Hildebrand 1985]
	Trf = rest-frame temperature [in K; default = 20K]
	Lrf = rest-frame FIR bolometric luminosity [in L_sun; default = 10^10]
	zin = galaxy redshift [default = 0.001]
	lambdavector = array of wavelengths of interest [in microns; default = (24, 70, 160, 250, 350, 500)];
	
	AUTHOR:
	Lorenzo Moncelsi [[email protected]]
	
	HISTORY:
	20June2012: created in IDL
	November2015: converted to Python
	'''

	nwv = len(lambdavector)
	nuvector = c * 1.e6 / lambdavector # Hz

	nsed = 1e4
	lambda_mod = loggen(1e3, 8.0, nsed) # microns
	nu_mod = c * 1.e6/lambda_mod # Hz

	#Lorenzo's version had: H0=70.5, Omega_M=0.274, Omega_L=0.726 (Hinshaw et al. 2009)
	cosmo = FlatLambdaCDM(H0 = 70.5 * u.km / u.s / u.Mpc, Om0 = 0.273)
	conversion = 4.0 * np.pi *(1.0E-13 * cosmo.luminosity_distance(zin) * 3.08568025E22)**2.0 / L_sun # 4 * pi * D_L^2    units are L_sun/(Jy x Hz)

	Lir = Lrf / conversion # Jy x Hz

	Ain = np.zeros(ngal) + 1.0e-36 #good starting parameter
	betain =  np.zeros(ngal) + b 
	alphain=  np.zeros(ngal) + 2.0

	fit_params = Parameters()
	fit_params.add('Ain', value= Ain)
	#fit_params.add('Tin', value= Trf/(1.+zin), vary = False)
	#fit_params.add('betain', value= b, vary = False)
	#fit_params.add('alphain', value= alphain, vary = False)

	#pdb.set_trace()
	#THE LM FIT IS HERE
	#Pfin = minimize(sedint, fit_params, args=(nu_mod,Lir.value,ngal))
	Pfin = minimize(sedint, fit_params, args=(nu_mod,Lir.value,ngal,Trf/(1.+zin),b,alphain))

	#pdb.set_trace()
	flux_mJy=sed(Pfin.params,nuvector,ngal,Trf/(1.+zin),b,alphain)

	return flux_mJy

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
def ag_mods(model,ti,tf,nui,nuf,j_ang,E_k,n,eps_b,eps_e,p,z,Gamma):
	""" Creates GRB afterglow light curves or spectra """

	# define physical constants in cgs units
	c = 2.998e10
	me = 9.109e-28
	mp = 1.673e-24
	e = 4.803e-10
	sigma_t = 6.652e-25

	# Calculate luminosity distance
	cosmo = FlatLambdaCDM(H0=70,Tcmb0=2.725,Om0=0.3)
	d_l = (cosmo.luminosity_distance(z).value*3.086e+24)

	t_dec,B,P_max,R,Ne = bw_props(Gamma,E_k,n,eps_b,p,mp,me,c,e,sigma_t)

	if model == 'SPECTRUM':
		freq = np.logspace(np.log10(nui),np.log10(nuf),1000)
		ymod = np.zeros(len(freq))
		for i in range(len(freq)):
			ymod[i] = model_flux(t_dec,B,P_max,R,Ne,d_l,z,mp,me,e,c,sigma_t,
				ti,freq[i],Gamma,E_k,n,eps_b,eps_e,p,j_ang)
					  
		plt.figure()
		plt.plot(freq,ymod,color='black',marker=' ')
		plt.xlabel('$\\nu$ (Hz)')
		plt.ylabel('Flux Density (mJy)')
		plt.xscale('log')
		plt.yscale('log')
		#plt.savefig('spectrum_mod.png')
		plt.show()
		
	if model == 'LC':
		times = np.logspace(np.log10(ti),np.log10(tf),1000)
		ymod = np.zeros(len(times))
		for i in range(len(times)):
			ymod[i] = model_flux(t_dec,B,P_max,R,Ne,d_l,z,mp,me,e,c,sigma_t,
				times[i],nui,Gamma,E_k,n,eps_b,eps_e,p,j_ang)		

		plt.figure()
		plt.plot(times,ymod,color='black',marker=' ')
		plt.xlabel('Time since trigger (s)')
		plt.ylabel('Flux Density (mJy)')
		plt.xscale('log')
		plt.yscale('log')
		#plt.savefig('lightcurve_mod.png')
		plt.show()
	return 0

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
def get_sdss_sample():
    
    from scipy import interpolate
    
    filepath = cu.get_output_path() + 'processed_data/NYU_VAGC/'
    galaxy_catalogue = 'nyu_lss_mpa_vagc_dr7'
    f =  h5py.File(filepath+galaxy_catalogue+'.hdf5', 'r')
    GC = f.get(galaxy_catalogue) #halo catalogue
    GC = np.array(GC)
    #for name in GC.dtype.names: print(name)
    
    #trim the catalogue a bit
    zmin = 0.01
    zmax = 0.2
    selection = (GC['M']<=17.6) & (GC['Z']<zmax) & (GC['Z']>zmin) &\
                (GC['ZTYPE']==1) & (GC['FGOTMAIN']>0)
    GC = GC[selection]
    
    sm_key = 'sm_MEDIAN'
    GC[sm_key] = GC[sm_key]+np.log10(0.7**2.0)
    
    #make cuts on data to get a clean and complete sample
    cosmo = FlatLambdaCDM(H0=100, Om0=0.3)
    z = np.linspace(0.001,1,1000)
    dL = cosmo.luminosity_distance(z).value
    
    #make cheater fit to dl(z) function
    dL_z = interpolate.interp1d(z,dL,kind='linear')
    Mstar_lim = (4.852 + 2.246*np.log10(dL) + 1.123*np.log10(1+z) - 1.186*z)/(1.0-0.067*z)
    
    #dL = cosmo.luminosity_distance(GC['Z']).value
    dL = dL_z(GC['Z'])
    z = GC['Z']
    LHS = (4.852 + 2.246*np.log10(dL) + 1.123*np.log10(1.0+z) - 1.186*z)/(1.0-0.067*z)
    keep = (GC[sm_key]>LHS)
    GC = GC[keep]
    
    return GC

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
def convert_Fline2Lbol(lineflux,linefluxerr,redshift,verbose=True):
    """
    Converting an observed integrated line flux [erg/s/cm2] to bolometric luminoisity [erg/s]

    --- EXAMPLE OF USE ---
    import NEOGALmodels as nm
    LbolLsun, LbolLsunerr = nm.convert_Fline2Lbol(2000,100,4.1)

    """
    cosmo       = FlatLambdaCDM(H0=70, Om0=0.3)
    if verbose: print ' - Estimating bolometric luminoisity for flat standard cosmology (H0=70, Om0=0.3, OL0=0.7)'
    DL          = cosmo.luminosity_distance(redshift).value
    # DLplus    = cosmo.luminosity_distance(redshift+redshifterr).value
    # DLminus   = cosmo.luminosity_distance(redshift-redshifterr).value
    Mpc2cm      = 3.086                          # 10**24 cm/Mpc
    Asphere     = 4*np.pi*(DL*Mpc2cm)**2         # 10**48 cm2
    Lbol        = lineflux*Asphere               # 10**28 erg/s ; assuming line fluxes are in 10**-20 erg/s/cm2
    Lbolerr     = linefluxerr*Asphere            # 10**28 erg/s ; assuming line fluxes are in 10**-20 erg/s/cm2
    LbolLsun    = Lbol/3.826*10**-5              # in units of Lbol_sun = 3.826*10**33 erg/s
    LbolLsunerr = Lbolerr/3.826*10**-5           # in units of Lbol_sun = 3.826*10**33 erg/s
    if verbose: print ' - Retunring luminoisity in units of Lbol_sun = 3.826e33 erg/s'
    if verbose: print ' - Result is: '+str(LbolLsun)+' +/- '+str(LbolLsunerr)+' [3.826e33 erg/s]'
    return LbolLsun, LbolLsunerr

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
class Data(object):
	"""docstring for Data"""
	def __init__(self, N_sn, seed, pop2=False, asifIa = False):
		super(Data, self).__init__()
		self.N_sn = N_sn
		self.seed = seed

		self.omega_M=0.28
		self.cosmo=FlatLambdaCDM(70,self.omega_M)

		self.zmin=0.1
		self.zmax=1.4

		self.sigma_snIa=0.1
		self.sigma_nonIa=1.
		self.sigma_nonIa_2=0.25

		self.alpha_snIa=2.
		self.alpha_nonIa=self.alpha_snIa*10**(-2./2.5)
		self.alpha_nonIa_2=self.alpha_snIa*10**(-0.5/2.5)

		self.frac_Ia_0=.95
		self.frac_Ia_1=.2

		self.asifIa = asifIa


		if pop2:
			self.frac_nonIa_0=1.
			self.frac_nonIa_1=0.2
		else:
			self.frac_nonIa_0=1.
			self.frac_nonIa_1=1.


		numpy.random.seed(seed)
		self.initialize_()

	def initialize_(self):
		self.redshifts_()
		self.types_()
		self.adus_()
		self.hosts_()

	def redshifts_(self):
		# zs : the true redshifts
		# volume_zero = self.cosmo.comoving_volume(self.zmin).value
		# volume_norm = self.cosmo.comoving_volume(self.zmax).value- volume_zero

		# self.zs = numpy.sort(numpy.random.uniform(size=self.N_sn))

		#self.zs=numpy.sort(numpy.random.uniform((1+self.zmin)**3,(1+self.zmax)**3,size=self.N_sn)**(1./3)-1)
		#self.zs=numpy.sort(numpy.random.uniform((1+self.zmin)**3,(1+self.zmax)**3,size=self.N_sn)**(1./3)-1)
		# for i in xrange(len(self.zs)):
		# 	self.zs[i]=scipy.optimize.newton(lambda z: (self.cosmo.comoving_volume(z).value - volume_zero)/volume_norm-self.zs[i], 0.5)		
		self.zs = numpy.sort(numpy.random.uniform(self.zmin**3, self.zmax**3,size=self.N_sn)**(1./3))


	def types_(self):
		self.snIa = numpy.random.uniform(size=self.N_sn)
		self.snIa = numpy.less_equal(self.snIa, self.frac_Ia_0 + (self.frac_Ia_1-self.frac_Ia_0)*self.zs/self.zmax/1.1).astype(int)
		self.nonIa = numpy.random.uniform(size=self.N_sn)
		self.nonIa = numpy.less_equal(self.nonIa, self.frac_nonIa_0 + (self.frac_nonIa_1-self.frac_nonIa_0)*(self.zs**2)/((self.zmax/1.1))**2).astype(int)

	def adus_(self):
		adu = 1/(self.cosmo.luminosity_distance(self.zs).value/self.cosmo.hubble_distance.value)**2
		adu_random = numpy.random.normal(size=self.N_sn)
		wsnIa = numpy.where(self.snIa)[0]
		wnonIa = numpy.where(numpy.logical_and(self.snIa ==0, self.nonIa ==1))[0]
		adu[wsnIa] = self.alpha_snIa*adu[wsnIa]*10**(adu_random[wsnIa]*self.sigma_snIa/2.5)
		adu[wnonIa] = self.alpha_nonIa*adu[wnonIa]*10**(adu_random[wnonIa]*self.sigma_nonIa/2.5)
		wpop2 = numpy.where(numpy.logical_and(self.snIa ==0, self.nonIa ==0))[0]
		adu[wpop2] = self.alpha_nonIa_2*adu[wpop2]*10**(adu_random[wpop2]*self.sigma_nonIa_2/2.5)
		self.adu = adu

	def hosts_(self):
		self.host_choice = numpy.random.binomial(1,0.98,size=self.N_sn)

		# zs : the true redshifts
#		volume_zero = self.cosmo.comoving_volume(self.zmin/1.1).value
#		volume_norm = self.cosmo.comoving_volume(self.zmax*1.1).value- volume_zero
#		self.neighbor_zs = numpy.sort(numpy.random.uniform(size=self.N_sn))
#		for i in xrange(len(self.zs)):
#			self.neighbor_zs[i]=scipy.optimize.newton(lambda z: (self.cosmo.comoving_volume(z).value - volume_zero)/volume_norm-self.zs[i], 0.5)		
		# self.neighbor_zs=numpy.random.uniform((1+self.zmin/1.1)**3,(1+self.zmax*1.1)**3,size=self.N_sn)**(1./3)-1
		self.neighbor_zs=numpy.random.uniform((self.zmin/1.1)**3, (self.zmax*1.1)**3,size=self.N_sn)**(1./3)
		self.host_zs_random = self.zs*self.host_choice + self.neighbor_zs*(1-self.host_choice)
		self.neighbor_zs_random = self.zs*(1-self.host_choice) + self.neighbor_zs*self.host_choice


	def found(self, ADU0):
		self.ADU0 =ADU0
		self.found_ = self.adu>=ADU0
		self.s_obs_random = numpy.where(self.found_)[0]
		numpy.random.shuffle(self.s_obs_random)

	def spectrum(self, frac_obs):
		if self.asifIa:
			frac_obs=1.
		self.N_s_obs = frac_obs*self.found_.sum()
#.........这里部分代码省略.........

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
        for j in range(115, 195, 3):#107 = 50Myr, 69 = 10Myr
     
            wavs, sfd = sp.get_spectrum(tage=ages[j])
            spectrum = np.array([wavs, sfd]).T
            spectrum[:,1] = spectrum[:,1]*(3.83*10**33)*10**(-EBV*coef/2.5) #erg/s/Hz
            spectrum[:,1] = spectrum[:,1]*(3*10**18)/spectrum[:,0]/spectrum[:,0] #erg/s/A
            synmags = np.zeros(13*700, dtype="float")
            synmags.shape = (700, 13)
            synmags[:,0] = np.arange(0.01, 7.01, 0.01)

            for i in range(1, 701):
                z = 0.01*i
                if ages[j]*(10**-9) < 14.00 - cosmo.lookback_time(z).value:
                    print "Calculating colours for age: " + str(ages[j]) + " and redshift: " + str(z) + ", total age: " + str(ages[j]*(10**-9) + cosmo.lookback_time(z).value) + ", tau = " + tauvals[t]
                    zspectrum = np.copy(spectrum)
                    for k in range(len(zspectrum)):
                        if zspectrum[k,0] < 912.:
                            zspectrum[k,1] = 0.
                        elif zspectrum[k,0] > 912. and zspectrum[k,0] < 1026.:
                            zspectrum[k,1] = zspectrum[k,1]*(1-D[i-1,1])
                        elif zspectrum[k,0] > 1026. and zspectrum[k,0] < 1216.:
                            zspectrum[k,1] = zspectrum[k,1]*(1-D[i-1,0])
                    zspectrum[:,1] = zspectrum[:,1]/(4*np.pi*(cosmo.luminosity_distance(z)*3.086*10**24)**2) #convert to observed flux at given redshift in erg/s/A/cm^2
                    zspectrum[:,1] = zspectrum[:,1]/(1+z) #reduce flux by a factor of 1/(1+z)   to account for redshifting
                    zspectrum[:,0] = zspectrum[:,0]*(1+z) #change wavelength values to account for redshifting

                    th_mags = s.AB_mags_UDS(zspectrum)
                    synmags[i-1, 1:] = th_mags

            np.savetxt("../../models/tau_fsps/" + tauvals[t] + "/synmags_age_" + str(ages[j]) + "_EBV_" + str(EBV) + ".txt", synmags)

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]

#.........这里部分代码省略.........
                
    def  initialize(self, plot=False):
        '''
        Will transform the magnitudes to fluxes and use the distance to the object to
        calculate the luminosity at each wavelength.
        '''

        if (not os.path.isdir(self.plotdir)):
            os.makedirs(self.plotdir)

        # generate the data  
        with warnings.catch_warnings():
            warnings.simplefilter("ignore")   
            self.wls, self.fluxes, self.fluxerrs = self._band2flux()

        #Plot the raw fluxes before correcting them.
        '''if (plot):
            plt.figure(figsize=(8,6))
            plt.errorbar(self.wls, self.fluxes, yerr=self.fluxerrs, marker="o", lw=0)
            for i in range(len(self.wls)):
                plt.text(self.wls[i], self.fluxes[i]*1.01, self.bands[i].split(",")[-1], alpha=.4)
            name = self._get_save_path(None, "fluxes_obs")
            plt.yscale("log")
            plt.xlabel("Wavelength [A]")
            plt.ylabel("log (Flux/[erg/cm2/s])")
            plt.tight_layout()
            plt.savefig(name, dpi=200)'''
            
        if not self.distMpc is None and self.distMpc !=0:
            print ("Using distance to the source of %.1e Mpc"%self.distMpc)
            fluxFactor = (4*np.pi*((self.distMpc*u.Mpc).to(u.cm) )**2).value

        elif (self.distMpc is None or self.distMpc==0 )and (not self.z is None and self.z != 0):
            self.distMpc = self.cosmo.luminosity_distance(self.z)
                
            #Compute the flux multiplication factor for the object if it is at distance distMpc
            #We transform that to cm, as the flux is in erg cm-2 s-1
            fluxFactor = (4*np.pi*(self.distMpc.to(u.cm) )**2).value
        
        else: # self.distMpc is None and self.z is None:
            #Here we do not use any multiplication flux factor
            print ("Warning: no redshift or distance provided!")
            fluxFactor = 1
            
        self.fluxes = self.fluxes * fluxFactor
        self.fluxerrs = self.fluxerrs * fluxFactor

        self._initialize_parameters(plot)


    def _initialize_parameters(self, plot=True):
        '''
        Runs the least squares optimiztion routine to find the best initial parameters 
        to start the MCMC with.
        '''

        lam = np.linspace(3000, 25000, 2000)
        
        if self.model == "BlackBody2_Av":
            flux_ini = self._model2_av_r_2(lam, self.initT1, self.initR1, self.av_host, self.initT2, self.initR2)

            p0 = (self.initT1, self.initR1, self.av_host, self.initT2, self.initR2)
            
            print ("Initial ", p0, "Fluxerrs", self.fluxes/self.fluxerrs)
            
            params, covar = curve_fit(self._model2_av_r_2, self.wls , self.fluxes, \

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
class Calculator:
    """
    Calculate various quantities under a specific cosmology, as well
    as some values with respect to Chandra ACIS detector properties.
    """
    def __init__(self, H0=71.0, Om0=0.27, Ob0=0.046):
        self.H0 = H0  # [km/s/Mpc]
        self.Om0 = Om0
        self.Ode0 = 1.0 - Om0
        self.Ob0 = Ob0
        self._cosmo = FlatLambdaCDM(H0=H0, Om0=Om0, Ob0=Ob0)

    def evolution_factor(self, z):
        return self._cosmo.efunc(z)

    def luminosity_distance(self, z, unit="Mpc"):
        dist = self._cosmo.luminosity_distance(z)
        return dist.to(au.Unit(unit)).value

    def angular_diameter_distance(self, z, unit="Mpc"):
        dist = self._cosmo.angular_diameter_distance(z)
        return dist.to(au.Unit(unit)).value

    def kpc_per_arcsec(self, z):
        """
        Calculate the transversal length (unit: kpc) corresponding to
        1 arcsec at the *angular diameter distance* of z.
        """
        dist_kpc = self.angular_diameter_distance(z, unit="kpc")
        return dist_kpc * au.arcsec.to(au.rad)

    def kpc_per_pix(self, z):
        """
        Calculate the transversal length (unit: kpc) corresponding to
        1 ACIS pixel (i.e., 0.492 arcsec) at the *angular diameter distance*
        of z.
        """
        pix = ACIS.pixel2arcsec * au.arcsec
        dist_kpc = self.angular_diameter_distance(z, unit="kpc")
        return dist_kpc * pix.to(au.rad).value

    def cm_per_pix(self, z):
        """
        Calculate the transversal length (unit: cm) corresponding to
        1 ACIS pixel (i.e., 0.492 arcsec) at the *angular diameter distance*
        of z.
        """
        return self.kpc_per_pix(z) * au.kpc.to(au.cm)

    def norm_apec(self, z):
        """
        The normalization factor of the XSPEC APEC model assuming
        EM = 1 (i.e., n_e = n_H = 1 cm^-3, and V = 1 cm^3)

        norm = 1e-14 / (4*pi* (D_A * (1+z))^2) * int(n_e * n_H) dV
        unit: [cm^-5]

        This value will be used to calculate the cooling function values.

        References
        ----------
        * XSPEC: APEC model
          https://heasarc.gsfc.nasa.gov/docs/xanadu/xspec/manual/XSmodelApec.html
        """
        da = self.angular_diameter_distance(z, unit="cm")
        norm = 1e-14 / (4*math.pi * (da * (1+z))**2)
        return norm

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
            if y[2] < redshiftlimit and y[0] < redchilimit:

                sntype = types[sn]
                if sntype[4] =='b':
                    ib+=1
                    redchib.append(y[0])
                else:
                    redchic.append(y[0])
                    ic+=1
                # redchi.append(y[2])
                #if y[0] > redchilimit:
                #   pass
                    #print sn, 'redchi: ', y[2]
                model.set(z=y[1][0],t0=y[1][1], amplitude=y[1][2]) 
                cosmo=FlatLambdaCDM(H0=70,Om0=0.3)
                dispc = cosmo.luminosity_distance(y[2])*(10**(6))  #input redshift, will return in Mpc
                
                timearray = np.arange(y[1][1],y[1][1]+60,1)
                fluxlc = model.bandflux('ptf48r',timearray, zp = 27.0, zpsys='ab')
                
                
                obsmag =[(-21.49 - 2.5*np.log10(x)) for x in fluxlc]
                absmag =[(np.min(point) - 5*(np.log10(dispc.value) - 1)) for point in obsmag]
                
                fluxdatapoints = np.array(y[3][:,4], dtype = float)
                fluxdatapointsall = np.array([(float(i)*10**((-27-21.49)/2.5)) for i in hml[:,1]], dtype = float)
                
                fluxyerr = np.array(y[3][:,5], dtype = float)
                fluxyerrall = np.array([(float(i)*10**((-27-21.49)/2.5)) for i in hml[:,2]], dtype = float)
                
                

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
    def __init__(self, Om0=0.315, H0=67.7, zmin=0.5,
                 zmax=7.0, ninterp=100):
        """ Initializer

        Parameters
        ----------
        Om0: float
          Matter density parameter

        H0: float
          Hubble constant, in km / s / Mpc

        zmin: float
          Minimum redshift supported.

        zmax: float
          Maximum redshift supported.

        ninterp: int
          Number of interpolation points for luminosity distance.

        """

        from astropy.cosmology import FlatLambdaCDM
        from astropy.units import Quantity
        import astropy.io.fits as fits
        from pkg_resources import resource_filename
        from scipy.interpolate import RectBivariateSpline as rbs
        from scipy.interpolate import interp1d

        self._zmin = float(zmin)
        self._zmax = float(zmax)
        self._Om0 = float(Om0)
        self._H0 = float(H0)
        self._ninterp = int(ninterp)

        if self._zmin == self._zmax:
            raise ValueError("No range between zmin and zmax")
        if self._zmin > self._zmax:
            self._zmin, self._zmax = self._zmax, self._zmin
        if self._zmin < 0.0:
            raise ValueError("zmin must be >= 0: {:f}".format(self._zmin))
        if self._Om0 <= 0.0:
            raise ValueError("Om0 must be positive: {:f}".format(self._Om0))
        if self._H0 <= 0.0:
            raise ValueError("H0 must be positive: {:f}".format(self._H0))
        if self._ninterp <= 0:
            raise ValueError("Ninterp must be > 0: {:d}".format(self._ninterp))

        # Set up luminosity distance interpolant.  We actually
        # interpolate log((1+z) / (4 pi d_L^2)) in log(1+z)
        cos = FlatLambdaCDM(H0=self._H0, Om0=self._Om0, Tcmb0=0.0, Neff=0.0)
        zrange = np.linspace(self._zmin, self._zmax, self._ninterp)
        mpc_in_cm = 3.0857e24
        prefac = 1.0 / (4 * math.pi * mpc_in_cm**2)
        lumdist = cos.luminosity_distance(zrange)
        if isinstance(lumdist, Quantity):
            lumdist = lumdist.value
        dlval = prefac * (1.0 + zrange) / lumdist**2
        self._dlfac = interp1d(np.log(1 + zrange), np.log(dlval))

        # Read in the data products, and set up interpolations on them
        sb_tpl = resource_filename(__name__, 'resources/SED_sb.fits')
        hdu = fits.open(sb_tpl)
        dat = hdu['SEDS'].data
        hdu.close()
        self._sblam = dat['LAMBDA'][0]
        self._sbumean = dat['UMEAN'][0]
        arg = np.argsort(self._sbumean)
        self._sbumean = self._sbumean[arg]
        self._sbrange = np.array([self._sbumean[0], self._sbumean[-1]])
        self._sbseds = dat['SEDS'][0, :, :].transpose()[arg, :]
        self._sbinterp = rbs(self._sbumean, self._sblam, self._sbseds,
                             kx=1, ky=1)

        ms_tpl = resource_filename(__name__, 'resources/SED_ms.fits')
        hdu = fits.open(ms_tpl)
        dat = hdu['SEDS'].data
        hdu.close()
        self._mslam = dat['LAMBDA'][0]
        self._msumean = dat['UMEAN'][0]
        arg = np.argsort(self._msumean)
        self._msumean = self._msumean[arg]
        self._msrange = np.array([self._msumean[0], self._msumean[-1]])
        self._msseds = dat['SEDS'][0, :, :].transpose()[arg, :]
        self._msinterp = []
        self._msinterp = rbs(self._msumean, self._mslam, self._msseds,
                             kx=1, ky=1)

# 需要导入模块: from astropy.cosmology import FlatLambdaCDM [as 别名]
# 或者: from astropy.cosmology.FlatLambdaCDM import luminosity_distance [as 别名]
            
            
            
        if w=="FUV":
            continuous_FUV=continuous_UV
            sp.data=sp.data - continuous_FUV 
            

    number,name,mag,group,redshift=np.genfromtxt("quasar_data.txt", dtype="|S10",unpack=True)       
    
    
    
    
    redshift=np.array([float(x) for x in redshift])
    mag=np.array([float(x) for x in mag])
    dl=a=cosmology.luminosity_distance(redshift).to(u.cm)
    
    F_nuv_galex=galex_mAB_Fnuv(mAB_nuv_galex)*u.erg/(u.second*u.cm*u.cm*u.angstrom)
    galex_central_wl=2315.7*u.angstrom/(1+redshift)
    L_wl=((1+redshift)*4*np.pi*F_nuv_galex*dl.to(u.cm)*dl.to(u.cm)).to(u.erg/(u.cm*u.s))
    L_nu=(L_wl*np.power(2315.7*u.angstrom/(1+redshift),2)/(3e10*u.cm/u.second)).to(u.erg)
    




    L_wl=np.multiply( L_wl ,10**(-1*mag_order) )
    
    L_wl=L_wl.value
    galex_central_wl=galex_central_wl.value