import numpy as np import fitsio import sys import os def tselectparams(ddphoto): """ Computes intermediate quantities needed for target selection from the input collate files. """ ngal = len(ddphoto) psfmag = np.zeros([ngal,5]) r_dev = np.zeros([ngal,5]) fib2mag = np.zeros([ngal,5]) modlmag = np.zeros([ngal,5]) cmodmag = np.zeros([ngal,5]) cmodflux = np.zeros(ngal) for i in range(5): cmodflux[:] = ddphoto['DEVFLUX'][:,i] * ddphoto['FRACPSF'][:,i] + ddphoto['EXPFLUX'][:,i] * (1. - ddphoto['FRACPSF'][:,i]) r_dev[:,i] = ddphoto['R_DEV'][:,i] ## clip doesn't work for this version of numpy, for some reason.. i guess clip is new? if 0==1: psfmag[:,i] = 22.5 - 2.5 * np.log10(ddphoto['PSFFLUX'][:,i].clip(min=0.001,max=None)) - ddphoto['EXTINCTION'][:,i] fib2mag[:,i] = 22.5 - 2.5 * np.log10(ddphoto['FIBER2FLUX'][:,i].clip(min=0.001,max=None)) - ddphoto['EXTINCTION'][:,i] modlmag[:,i] = 22.5 - 2.5 * np.log10(ddphoto['MODELFLUX'][:,i].clip(min=0.001,max=None)) - ddphoto['EXTINCTION'][:,i] cmodmag[:,i] = 22.5 - 2.5 * np.log10(cmodflux.clip(min=0.001,max=None)) - ddphoto['EXTINCTION'][:,i] else: psfmag[:,i] = 22.5 - 2.5 * np.log10(ddphoto['PSFFLUX'][:,i]) - ddphoto['EXTINCTION'][:,i] xx = np.where(ddphoto['PSFFLUX'][:,i] <= 0.001)[0]; psfmag[xx,i] = 22.5 - 2.5 * np.log10(0.001) - ddphoto['EXTINCTION'][xx,i] xx = np.where(np.isnan(psfmag))[0]; assert len(xx) == 0 fib2mag[:,i] = 22.5 - 2.5 * np.log10(ddphoto['FIBER2FLUX'][:,i]) - ddphoto['EXTINCTION'][:,i] xx = np.where(ddphoto['FIBER2FLUX'][:,i] <= 0.001)[0]; fib2mag[xx,i] = 22.5 - 2.5 * np.log10(0.001) - ddphoto['EXTINCTION'][xx,i] xx = np.where(np.isnan(fib2mag))[0]; assert len(xx) == 0 modlmag[:,i] = 22.5 - 2.5 * np.log10(ddphoto['MODELFLUX'][:,i]) - ddphoto['EXTINCTION'][:,i] xx = np.where(ddphoto['MODELFLUX'][:,i] <= 0.001)[0]; modlmag[xx,i] = 22.5 - 2.5 * np.log10(0.001) - ddphoto['EXTINCTION'][xx,i] xx = np.where(np.isnan(modlmag))[0]; assert len(xx) == 0 cmodmag[:,i] = 22.5 - 2.5 * np.log10(cmodflux) - ddphoto['EXTINCTION'][:,i] xx = np.where(cmodflux <= 0.001)[0]; cmodmag[xx,i] = 22.5 - 2.5 * np.log10(0.001) - ddphoto['EXTINCTION'][xx,i] xx = np.where(np.isnan(cmodmag))[0]; assert len(xx) == 0 cperp = (modlmag[:,2] - modlmag[:,3]) - (modlmag[:,1] - modlmag[:,2])/4.-0.18 cpar = 0.7 * (modlmag[:,1] - modlmag[:,2]) + 1.2 * (modlmag[:,2] - modlmag[:,3] - 0.18) dperp = (modlmag[:,2] - modlmag[:,3]) - (modlmag[:,1] - modlmag[:,2])/8. return psfmag, r_dev, fib2mag, modlmag, cmodmag, cperp, cpar, dperp def lrg_select(dd, ddphoto, cmassorlowz=0): """ cmassorlowz = 0 for CMASS, 1 for LOWZ. copied from make_samples_OLD.pro. Based on my work on Oct 31 2014 (see email to boss-gal-clustering), this is my best guess at the chunk 2 and chunk 3-6 selection. We can run this selection on chunks >= 7 for comparison if necessary. If you're using the new dtype created below, submit same structure for dd and ddphoto. ddphoto used for mag/color cut info dd used for BOSS_TARGET1_009 entry. """ psfmag, r_dev, fib2mag, modlmag, cmodmag, cperp, cpar, dperp = tselectparams(ddphoto) if cmassorlowz == 0: # cmass xx = np.where((cmodmag[:,3] > 17.5) & \ (cmodmag[:,3] < 19.9) & \ (fib2mag[:,3] < 21.5) & \ ## closer inspection of hte targeting code implies this is not actually used (??) Double check with Nikhil. ## Never mind, it looks like it WAS applied, since if not some LOWZ objects pass CMASS cuts. (r_dev[:,3] < 20.0) & \ (dperp > 0.55) & \ (modlmag[:,2] - modlmag[:,3] < 2.0) & \ (cmodmag[:,3] < 19.86 + 1.6*(dperp - 0.8)) & \ (psfmag[:,4]-modlmag[:,4] > 9.125 - 0.46*modlmag[:,4]) & \ (psfmag[:,3]-modlmag[:,3] > 0.2 + 0.2*(20.0-modlmag[:,3])))[0] if cmassorlowz == 1: # lowz xx = np.where((cmodmag[:,2] > 16.0) & \ (cmodmag[:,2] < 19.6) & \ (np.fabs(cperp) < 0.2) & \ (cmodmag[:,2] < 13.5 + cpar/0.3) & \ (psfmag[:,2]-modlmag[:,2] > 0.3))[0] ## just to safety, ALSO require LOWZ flag to be set. I verified that this doesn't remove any targets ## that otherwise pass the cuts. if cmassorlowz == 2: ## chunk 2. xx = np.where((dd['BOSS_TARGET1_009'] & 2**0 == 2**0) & (dd['BOSS_TARGET1_009'] > 0) & \ (cmodmag[:,2] > 16.0) & \ (cmodmag[:,2] < 19.5) & \ (np.fabs(cperp) < 0.2) & \ (cmodmag[:,2] < 13.4 + cpar/0.3) & \ (psfmag[:,2]-modlmag[:,2] > 0.3) & \ (psfmag[:,3]-modlmag[:,3] > 0.2 + 0.2*(20.0-modlmag[:,3])))[0] if cmassorlowz == 3: ## chunk 3. xx = np.where((dd['BOSS_TARGET1_009'] & 2**0 == 2**0) & (dd['BOSS_TARGET1_009'] > 0) & \ (cmodmag[:,2] > 17.0) & \ (cmodmag[:,2] < 19.5) & \ (np.fabs(cperp) < 0.2) & \ (cmodmag[:,2] < 13.4 + cpar/0.3) & \ (psfmag[:,2]-modlmag[:,2] > 0.3) & \ (psfmag[:,4]-modlmag[:,4] > 9.125 - 0.46*modlmag[:,4]) & \ (psfmag[:,3]-modlmag[:,3] > 0.2 + 0.2*(20.0-modlmag[:,3])))[0] return xx def get_collate_targets(collatefile,collatefilephoto,cmass_output=None,lowz_output=None, lowzearly2_output=None,lowzearly3_output=None): """ If you want a target output, fill in cmass_output, lowz_output, lowzearly_output. If you want to skip those outputs, don't fill those fields (default is "None") """ dd = fitsio.read(collatefile) ddphoto = fitsio.read(collatefilephoto) assert len(dd) == len(ddphoto) #get_collate_targets() ## note to self, ichunk comes from ddphoto! But they're equal for dr12v4, all good. ## >>> (ddphoto['ICHUNK'] == dd['CHUNK']).all() ##True targdtype = np.dtype([('RA', '>f8'), ('DEC', '>f8'), ('BOSS_TARGET1', '>i4'), ('BOSS_TARGET1_009', '>i4'), ('RUN', '>i4'), ('RERUN', '|S12'), ('CAMCOL', '>i4'), ('FIELD', '>i4'), ('ID', '>i4'), ('ICHUNK', '>i4'), ('IPOLY', '>i4'), ('ISECT', '>i4'), ('FINALN', '>i4'), ('FINALID', '>i4'), ('NASSIGNED', '>i4'), ('WHICHFINAL', '>i4'), ('TILE', '>i2', (3L,)), ('COLLATEINDX', '>i4'), ('ICOLLIDED', '>i4'), ('COLL_GROUP_FINALID', '>i4', (10L,)), ('COLL_GROUP_BOSSTARGET1', '>i4', (10L,)), ('INGROUP', '>i4'), ('FIRSTGROUP', '>i4'), ('MULTGROUP', '>i4'), ('NEXTGROUP', '>i4'), ('FRACPSF', '>f4', (5L,)), ('EXPFLUX', '>f4', (5L,)), ('DEVFLUX', '>f4', (5L,)), ('PSFFLUX', '>f4', (5L,)), ('MODELFLUX', '>f4', (5L,)), ('FIBER2FLUX', '>f4', (5L,)), ('R_DEV', '>f4', (5L,)), ('PSF_FWHM', '>f4', (5L,)), ('MODELFLUX_IVAR', '>f4', (5L,)), ('AIRMASS', '>f4'), ('SKYFLUX', '>f4', (5L,)), ('IMAGE_DEPTH', '>f4', (5L,)), ('EB_MINUS_V', '>f4'), ('EXTINCTION', '>f4', (5L,))]) eetmp = np.recarray(len(dd),dtype=targdtype) for ff in targdtype.fields: if ff in dd.dtype.fields: # print 'in1',ff eetmp[ff] = dd[ff] else: if ff in ddphoto.dtype.fields: # print 'in2',ff eetmp[ff] = ddphoto[ff] else: print 'none',ff ## write to a temporary file since I don't know how to do veto stuff in python. eetmp['COLLATEINDX'][:] = np.arange(0,len(eetmp),1,dtype='>i4') fitsio.write('tmp.fits',eetmp,clobber=True) norig = len(eetmp) del eetmp del dd del ddphoto ## veto with all the veto masks. os.system("""idl -e 'vetofile,"tmp.fits","tmp.fits"'""") dd = fitsio.read('tmp.fits') ## now has been vetoed, just need to do target selection and then print final files out. ntot = len(dd) print 'non-veto fraction = ',float(ntot)/float(norig) cmassi = np.where((dd['BOSS_TARGET1_009'] & 2**1 == 2**1) & (dd['BOSS_TARGET1_009'] > 0))[0] lowzi = np.where((dd['BOSS_TARGET1_009'] & 2**0 == 2**0) & (dd['BOSS_TARGET1_009'] > 0))[0] print 'Number of CMASS targets after trim:',len(cmassi) print 'Number of LOWZ targets after trim:',len(lowzi) if cmass_output is not None: ddcmass = dd[cmassi] print 'writing data to ',cmass_output fitsio.write(cmass_output,ddcmass,clobber=True) if lowz_output is not None: ddlowz = dd[lowzi] print 'writing data to ',lowz_output fitsio.write(lowz_output,ddlowz,clobber=True) if lowzearly2_output is not None: ## see message 3917 on boss-gal-clustering list. ## chunk 2 and chunks 3-6 appear to be selected differently. lowzearlyi2 = lrg_select(dd, dd, cmassorlowz=2) ddlowzearly2 = dd[lowzearlyi2] fitsio.write(lowzearly2_output,ddlowzearly2,clobber=True) if lowzearly3_output is not None: lowzearlyi3 = lrg_select(dd, dd, cmassorlowz=3) ddlowzearly3 = dd[lowzearlyi3] fitsio.write(lowzearly3_output,ddlowzearly3,clobber=True) def get_collate_targets_vetotmp(collatefile,collatefilephoto,cmassorlowz=0): """ this is for figuring out veto fraction involved in seeing and extinction. """ dd = fitsio.read(collatefile) ddphoto = fitsio.read(collatefilephoto) assert len(dd) == len(ddphoto) #get_collate_targets() ## note to self, ichunk comes from ddphoto! But they're equal for dr12v4, all good. ## >>> (ddphoto['ICHUNK'] == dd['CHUNK']).all() ##True targdtype = np.dtype([('RA', '>f8'), ('DEC', '>f8'), ('BOSS_TARGET1', '>i4'), ('BOSS_TARGET1_009', '>i4'), ('RUN', '>i4'), ('RERUN', '|S12'), ('CAMCOL', '>i4'), ('FIELD', '>i4'), ('ID', '>i4'), ('ICHUNK', '>i4'), ('IPOLY', '>i4'), ('ISECT', '>i4'), ('FINALN', '>i4'), ('FINALID', '>i4'), ('NASSIGNED', '>i4'), ('WHICHFINAL', '>i4'), ('TILE', '>i2', (3L,)), ('COLLATEINDX', '>i4'), ('ICOLLIDED', '>i4'), ('COLL_GROUP_FINALID', '>i4', (10L,)), ('COLL_GROUP_BOSSTARGET1', '>i4', (10L,)), ('INGROUP', '>i4'), ('FIRSTGROUP', '>i4'), ('MULTGROUP', '>i4'), ('NEXTGROUP', '>i4'), ('FRACPSF', '>f4', (5L,)), ('EXPFLUX', '>f4', (5L,)), ('DEVFLUX', '>f4', (5L,)), ('PSFFLUX', '>f4', (5L,)), ('MODELFLUX', '>f4', (5L,)), ('FIBER2FLUX', '>f4', (5L,)), ('R_DEV', '>f4', (5L,)), ('PSF_FWHM', '>f4', (5L,)), ('MODELFLUX_IVAR', '>f4', (5L,)), ('AIRMASS', '>f4'), ('SKYFLUX', '>f4', (5L,)), ('IMAGE_DEPTH', '>f4', (5L,)), ('EB_MINUS_V', '>f4'), ('EXTINCTION', '>f4', (5L,))]) eetmp = np.recarray(len(dd),dtype=targdtype) for ff in targdtype.fields: if ff in dd.dtype.fields: # print 'in1',ff eetmp[ff] = dd[ff] else: if ff in ddphoto.dtype.fields: # print 'in2',ff eetmp[ff] = ddphoto[ff] else: print 'none',ff ## write to a temporary file since I don't know how to do veto stuff in python. eetmp['COLLATEINDX'][:] = np.arange(0,len(eetmp),1,dtype='>i4') fitsio.write('tmpvetotest.fits',eetmp,clobber=True) norig = len(eetmp) del eetmp del dd del ddphoto ## veto with all the veto masks. os.system("""idl -e 'vetofiletmp,"tmpvetotest.fits","tmpvetotest.fits"'""") dd = fitsio.read('tmpvetotest.fits') ## now has been vetoed, just need to do target selection and then print final files out. ntot = len(dd) print 'non-veto fraction = ',float(ntot)/float(norig) del dd if __name__ == '__main__': get_collate_targets("inputfiles/bosstile-final-collated-boss2-boss37.fits","inputfiles/bosstile-final-collated-boss2-boss37-galimaging.fits","../data/trimmed-collate-cmasspy-dr12v0.fits","../data/trimmed-collate-lowzzpy-dr12v0.fits","../data/trimmed-collate-lowzztbugpy-dr12v0.fits") ## get this to work. #print """idl -e 'vetofile,"tmp.fits","tmp2.fits"'"""