Python Sed.manyMagCalc方法代码示例

# 需要导入模块: from lsst.sims.photUtils import Sed [as 别名]
# 或者: from lsst.sims.photUtils.Sed import manyMagCalc [as 别名]
def get_TotalMags(result, bandpasses=('u','g','r','i','z','y')):
    datadir = os.environ.get("SIMS_SED_LIBRARY_DIR")
    tpath = os.getenv('LSST_THROUGHPUTS_DEFAULT')
    bands = {"u":None, "g":None, "r":None, "i":None, "z":None, "y":None}
    for k in bands:
        bands[k] = Bandpass()
        bands[k].readThroughput(os.path.join(tpath, "total_%s.dat"%k))
    # Set up phi, the wavelength-normalized system response for each filter,
    # for each bandpass for manyMagCalc method.
    bplist = []
    for f in ['u','g','r','i','z','y']:
        bands[f].sbTophi()
        bplist.append(bands[f])
    ids = result['galid']
    diskfile = result['sedFilenameDisk']
    bulgefile = result['sedFilenameBulge']
    agnfile = result['sedFilenameAgn']

    diskmn = result['magNormDisk']
    bulgemn = result['magNormBulge']
    agnmn = result['magNormAgn']

    bulgeAv = result['internalAvBulge']
    diskAv = result['internalAvDisk']

    redshift = result['redshift']

    imsimband = Bandpass()
    imsimband.imsimBandpass()
    sedDict = {}
    retMags = dict([(k, []) for k in bands])
    a_int = None
    b_int = None
    tmpwavelen = None
    for id, df, dm, dav, bf, bm, bav, af, am, z in zip(ids, diskfile, diskmn, diskAv, 
            bulgefile, bulgemn, bulgeAv, agnfile, agnmn, redshift):
        tmpflux = None
        for comp in ((df, dm, dav, 'galaxySED', False), (bf, bm, bav, 'galaxySED', False), (af, am, None, 'agnSED', True)):
        #Zero out the AGN contribution
        #for comp in ((df, dm, dav, 'galaxySED', False), (bf, bm, bav, 'galaxySED', False), (af, 99.99, None, 'agnSED', True)):
            if not comp[0] == u'None':
                if sedDict.has_key(comp[0]):
                    sed = copy.deepcopy(sedDict[comp[0]])
                else:
                    sed = Sed()
                    print os.path.join(datadir,comp[3],comp[0])
                    sed.readSED_flambda(os.path.join(datadir,comp[3],comp[0]))
		    if comp[4]:
		        sed.resampleSED(wavelen_match=tmpwavelen)
                    sedDict[comp[0]] = sed
                if a_int is None:
                    phiarray, dlambda = sed.setupPhiArray(bplist)
                    a_int, b_int = sed.setupCCMab()
		    #Careful, this assumes that a disk or bulge sed is read
		    #before any agn sed
		    tmpwavelen = sed.wavelen
                fNorm = sed.calcFluxNorm(comp[1], imsimband)
                sed.multiplyFluxNorm(fNorm)
                #I guess this assumes rv=3.1??
                if comp[2]:
                    sed.addCCMDust(a_int, b_int, A_v=comp[2])
		wavelenArr=sed.wavelen
		if tmpflux is None:
		    tmpflux = sed.flambda
		else:
	            tmpflux += sed.flambda
	newgal = Sed(wavelen=wavelenArr, flambda=tmpflux)
        #a_mw, b_mw = sed.setupCCMab()
        #sed.addCCMDust(a_mw, b_mw, A_v=mwav)
        newgal.redshiftSED(z, dimming=True)
	newgal.resampleSED(wavelen_match=bplist[0].wavelen)
	newgal.flambdaTofnu()
        mags = newgal.manyMagCalc(phiarray, dlambda)
        for i,k in enumerate(['u','g','r','i','z','y']):
            retMags[k].append(mags[i])
    return retMags

# 需要导入模块: from lsst.sims.photUtils import Sed [as 别名]
# 或者: from lsst.sims.photUtils.Sed import manyMagCalc [as 别名]
bplist = []
for f in filterlist:
    lsstbp[f].sbTophi()
    bplist.append(lsstbp[f])
phiarray, dlambda = stars[starlist[0]].setupPhiArray(bplist)
# Set up dictionary + arrays to hold calculated magnitude information.
mags2 = {}
for f in filterlist:
    mags2[f] = numpy.zeros(num_star, dtype="float")
# For each star (in num_star's), apply MW dust, fluxnorm, and then calculate mags using manyMagCalc.
for i in range(num_star):
    starname = starlist[star_name[i]]
    tmpstar = Sed(wavelen=stars[starname].wavelen, flambda=stars[starname].flambda)
    tmpstar.addCCMDust(a_mw, b_mw, ebv=ebv_mw[i])
    tmpstar.multiplyFluxNorm(fluxnorm[i])
    tmpmags = tmpstar.manyMagCalc(phiarray, dlambda)
    j = 0
    for f in filterlist:
        mags2[f][i] = tmpmags[j]
        j = j + 1
dt, t = dtime(t)
print "Calculating dust/fluxnorm/%d magnitudes for %d stars optimized with Sed.manyMagCalc took %f s" % (
    len(filterlist),
    num_star,
    dt,
)


"""
This test turned out so badly that I removed these methods from Bandpass.
Basically, it was slower than just the normal calculation above, and took up too much memory.

# 需要导入模块: from lsst.sims.photUtils import Sed [as 别名]
# 或者: from lsst.sims.photUtils.Sed import manyMagCalc [as 别名]
# Now, we can calculate the magnitudes in multiple bandpasses, for each SED.

# make the bandpass list
bplist = []
for filter in lsstfilterlist:
    bplist.append(lsst[filter])
phiArray, dlambda = tmpstar.setupPhiArray(bplist)
# store the values in a 2-d array (could do in a dictionary of arrays too)
mags = n.empty((len(starskeys), len(bplist)), dtype="float")
for i in range(len(starskeys)):
    # make a copy of the original SED *if* you want to 'reuse' the SED for multiple magnitude
    # calculations with various fluxnorms and dust applications (otherwise just use the object
    # you instantiated above)
    tmpstar = Sed(wavelen=stars[starskeys[i]].wavelen, flambda=stars[starskeys[i]].flambda)
    tmpstar.addCCMDust(a, b, ebv=ebv[i])
    tmpstar.multiplyFluxNorm(fluxnorm[i])
    mags[i] = tmpstar.manyMagCalc(phiArray, dlambda)

print "#sedname mag_u  mag_g   mag_r   mag_i   mag_z   mag_y"
for i in range(len(starskeys)):
    print "%s  %.4f %.4f %.4f %.4f %.4f %.4f" % (
        starskeys[i],
        mags[i][0],
        mags[i][1],
        mags[i][2],
        mags[i][3],
        mags[i][4],
        mags[i][5],
    )

# 需要导入模块: from lsst.sims.photUtils import Sed [as 别名]
# 或者: from lsst.sims.photUtils.Sed import manyMagCalc [as 别名]
# Set up dictionary + arrays to hold calculated magnitude information. 
mags2 = {}
for f in filterlist:
    mags2[f] = numpy.zeros(num_gal, dtype='float')
# For each galaxy (in num_gal's), apply internal dust, redshift, resample to 300-1200 nm, apply MW dust on
#   shorter (and standardized) wavelength range, fluxnorm, and then calculate mags using manyMagCalc. 
for i in range(num_gal):
    galname = gallist[gal_name[i]]
    tmpgal = Sed(wavelen=gals[galname].wavelen, flambda=gals[galname].flambda)
    tmpgal.addCCMDust(a_int, b_int, ebv=ebv_int[i])
    tmpgal.redshiftSED(redshifts[i])
    tmpgal.resampleSED(wavelen_min=wavelen_min, wavelen_max=wavelen_max, wavelen_step=wavelen_step)
    tmpgal.addCCMDust(a_mw, b_mw, ebv=ebv_mw[i])
    tmpgal.multiplyFluxNorm(fluxnorm[i])
    #tmpgal.flambdaTofnu() #Not needed because multiplyFluxNorm calculates fnu
    tmpmags = tmpgal.manyMagCalc(phiarray, dlambda)
    j = 0
    for f in filterlist:
        mags2[f][i] = tmpmags[j]
        j = j+1
dt, t = dtime(t)
print "Calculating dust/redshift/dust/fluxnorm/%d magnitudes for %d galaxies optimized way took %f s" \
      %(len(filterlist), num_gal, dt)


# Check for differences in magnitudes.
import pylab
pylab.figure()
colors = ['m', 'g', 'r', 'b', 'k', 'y']
i = 0
diff = {}