In [8]:
import astropy.table as table
import numpy as np
from defcuts import *
from defflags import *
from def_get_mags import *
bands=['g', 'r', 'i','z', 'y']
daperture=[1.01,1.51,2.02,3.02,4.03,5.71,8.40,11.8,16.8,23.5]
aperture=[x*0.5 for x in daperture]
ty='mean'
tag='outcut'
txtdist= 'Figure2'
txtslope='Figure1'
outdir='/Users/amandanewmark/repositories/galaxy_dark_matter/lumprofplots/clumps/'+ty+tag
doutdir='/Users/amandanewmark/repositories/galaxy_dark_matter/lumprofplots/distribution/'+ty+tag
Flags=['flags_pixel_bright_object_center', 'brobj_cen_flag-', 'No Bright Ojbect Centers', 'Only Bright Object Centers', 'brobj_cen_flag']
indir='/Users/amandanewmark/repositories/galaxy_dark_matter/GAH/'
bigdata = table.Table.read(indir+ 'LOWZ_HSCGAMA15_apmgs+cmodmag.fits')
def do_cuts(datatab):
parm=['flags_pixel_saturated_center','flags_pixel_edge','flags_pixel_interpolated_center','flags_pixel_cr_center','flags_pixel_suspect_center', 'flags_pixel_clipped_any','flags_pixel_bad']
ne=[99.99, 199.99, 0.0]
mincut=0.1
maxcut=''
cutdata=not_cut(datatab, bands, 'mag_forced_cmodel', ne)
for b in range(0, len(bands)-1):
newdata=many_flags(cutdata, parm, bands[b]) #flags not in y?
cutdata=newdata
return newdata
def get_TF(data):
bandi=['i']
Flag, Not,lab= TFflag(bandi,Flags, data)
return Flag, Not
newdata=do_cuts(bigdata)
Flagdat, Notdat=get_TF(newdata)
In [13]:
def get_ind_lums(newdata, bands, aperture, scale=''):
import numpy as np
from def_get_mags import get_zdistmod, get_kcorrect2, aper_and_comov, abs2lum, lumdensity, abs_mag
import math
from defclump import meanlum2
from my_def_plots import halflight_plot, scatter_fit
from scipy import interpolate
import matplotlib.pyplot as plt
from def_mymath import halflight
Naps=len(aperture)
Ndat=len(newdata)
try:
redshifts=newdata['Z']
DM= get_zdistmod(newdata, 'Z')
except:
redshifts=newdata['Z_2']
DM= get_zdistmod(newdata, 'Z_2')
kcorrect=get_kcorrect2(newdata,'mag_forced_cmodel', '_err', bands, '','hsc_filters.dat',redshifts)
fig=plt.figure()
bigLI=[]
bigrad=[]
bigden=[]
for n in range(0, Ndat):
LI=[]
LI2=[]
lumdi=[]
string=str(n)
radkpc=aper_and_comov(aperture, redshifts[n])
#print('redshifts is ', redshifts[n])
for a in range(0, Naps): #this goes through every aperture
ns=str(a)
#print('aperture0',ns)
absg, absr, absi, absz, absy= abs_mag(newdata[n], 'mag_aperture0', kcorrect, DM[n], bands, ns, n)
Lumg, Lumr, Lumi, Lumz, Lumy=abs2lum(absg, absr, absi, absz, absy)
Lg, Lr, Li, Lz, Ly=lumdensity(Lumg, Lumr, Lumi, Lumz, Lumy, radkpc[a])
if scale== 'log':
#print('getting logs')
logLumi=math.log10(Lumi)
logLi=math.log10(Li)
LI.append(logLumi)
lumdi.append(logLi)
else:
LI.append(Lumi)
lumdi.append(Li)
#print('LI for ',n,' galaxy is ', LI)
bigLI.append(LI)
bigden.append(lumdi)
if scale== 'log':
lograd=[math.log10(radkpc[n]) for n in range(len(radkpc))]
bigrad.append(lograd)
else:
bigrad.append(radkpc)
bigLIs=np.array(bigLI)
bigrads=np.array(bigrad)
lumdensi=np.array(bigden)
return bigLIs, bigrads, lumdensi
def my_halflight2(dat1):
loglum, lograd, loglumd= get_ind_lums(dat1, bands, aperture, scale='log')
return loglum, lograd, loglumd
In [14]:
loglum, lograd, loglumd=my_halflight2(Flagdat)
In [17]:
def get_halfrad(lograds, loglums):
from scipy import interpolate
import math
import numpy as np
print('from halflight_math')
maxlums=10**np.max(loglums)
halfL=maxlums/2
loghalfL=np.log10(halfL)
f=interpolate.interp1d(loglums,lograds, kind='linear', axis=-1)
logr12=f(loghalfL)
return logr12
def upper_rad_cut1(loglum, lograd, logden, m): #this should get rid of galaxies outside 4r1/2
from def_halflight_math import get_halfrad
nloglum=[]
nlograd=[]
nlogden=[]
mult=m
print(len(loglum), len(lograd))
N=len(lograd)
for n in range(N):
loglums=loglum[n]
lograds=lograd[n]
logdens=logden[n]
#print(loglums, lograds)
logr12=get_halfrad(lograds,loglums)
print(n)
r12=10**logr12
print(logr12)
r412=mult*r12
logr412=np.log10(r412)
print(logr412)
if np.max(lograds) >= logr412:
lograd=lograds[(lograds>=logr412)&(lograds<=logr412)]
if len(lograd)>=4:
nloglum.append(loglums)
nlograd.append(lograds)
nlogden.append(logdens)
else:
print('not enough data points')
else:
print('Upper limit out of range')
nloglum=np.array(nloglum)
nlograd=np.array(nlograd)
nlogden=np.array(nlogden)
return nloglum, nlograd, nlogden
loglum, lograd, loglumd= upper_rad_cut(loglum, lograd, loglumd, 4)
In [ ]:
In [ ]:
In [ ]:
In [ ]: