ivs.catalogs.crossmatch

1 # -*- coding: utf-8 -*- 2 """ 3 Crossmatch catalogs for data on one target, or crossmatch entire catalogs. 4 5 This module is callable from the command line. E.g., simply type:: 6 7 $:> python crossmatch.py get_photometry ID=HD180642 8 9 or ask for help on any function:: 10 11 $:> python crossmatch.py get_photometry --help 12 """ 13 import logging 14 import itertools 15 import pylab as pl 16 import numpy as np 17 18 from ivs.catalogs import vizier 19 from ivs.catalogs import gator 20 from ivs.catalogs import gcpd 21 from ivs.catalogs import mast 22 from ivs.catalogs import sesame 23 from ivs.units import conversions 24 from ivs.aux import numpy_ext 25 from ivs.aux import progressMeter 26 from ivs.aux import loggers 27 from ivs.aux import argkwargparser 28 from ivs.inout import ascii 29 30 from scipy.spatial import KDTree 31 32 logger = logging.getLogger("CAT.XMATCH") 33

34 -def get_photometry(ID=None,to_units='erg/s/cm2/AA',extra_fields=[],include=None, 35 exclude=None,**kwargs):

36 """ 37 Collect photometry from different sources. 38 39 The output consists of a record array containing the following keys: 40 41 'meas': the measurement's value directly from the catalog in original units 42 'e_meas': the error on the measurements 43 'flag': any flags that are associated with the measurement in a catalog 44 'unit': the unit of the original measurement 45 'source' the source catalog's name 46 'photband': photometric pass bands' name 47 'cwave': the effective wavelength of the passband 48 'cmeas': converted measurement (to C{to_units}) 49 'e_cmeas': error on converted measurment 50 'cunit': converted unit 51 52 Be aware that some of the values can be 'nan': e.g. sometimes no error is 53 listed in the catalog, or no flag. Also the 'cwave' column will be nan for 54 all photometric colours (e.g. B-V) 55 56 If you define C{extra_fields}, make sure all the {get_photometry} know how 57 to handle it: probably some default values need to be inserted if these 58 extra columns are not available in some catalog. It is safest just to leave 59 it blank. 60 61 You can include or exclude search sources via C{include} and C{exclude}. 62 When given, these should be a list containing strings. The default is to 63 include C{gator}, C{vizier} and C{gcpd}. 64 65 Extra keyword arguments are passed to each C{get_photometry} functions in 66 this package's modules. 67 68 Example usage: 69 70 1. You want to download all available photometry and write the results to 71 an ASCII file for later reference. 72 73 >>> master = get_photometry('vega') 74 >>> ascii.write_array(master,header=True,auto_width=True) 75 76 2. You want to plot the raw, unmodelled SED of an object: 77 78 >>> master = get_photometry('vega') 79 >>> pl.errorbar(master['cwave'],master['cmeas'],yerr=master['e_cmeas'],fmt='ko') 80 >>> pl.gca().set_xscale('log',nonposx='clip') 81 >>> pl.gca().set_yscale('log',nonposy='clip') 82 83 We made no difference between colors (B-V) and magnitude (V), because the 84 'cwave' for colors is 'nan', so they will not be plotted anyway.The final 85 two lines are just to correct errorbars that go below zero in a logarithmic 86 plot. 87 88 @param ID: the target's name, understandable by SIMBAD 89 @type ID: str 90 @param to_units: units to convert everything to. 91 @type to_units: 92 @param include: sources to include 93 @type include: list of strings (from C{gator}, C{vizier} or C{gcpd}) 94 @param exclude: sources to include 95 @type exclude: list of strings (from C{gator}, C{vizier} or C{gcpd}) 96 @return: record array where eacht entry is a photometric measurement 97 @rtype: record array 98 """ 99 #-- make sure all catalog names are lower case 100 if include is not None: include = [i.lower() for i in include] 101 if exclude is not None: exclude = [i.lower() for i in exclude] 102 103 #-- check which sources to include/exclude 104 searchables = ['gator','vizier','gcpd','mast'] 105 if include is not None: 106 searchables = include 107 if exclude is not None: 108 searchables = list( set(searchables)- set(exclude)) 109 110 #-- and search photometry 111 if 'mast' in searchables: 112 kwargs['master'] = mast.get_photometry(ID=ID,to_units=to_units,extra_fields=extra_fields,**kwargs) 113 if 'gator' in searchables: 114 kwargs['master'] = gator.get_photometry(ID=ID,to_units=to_units,extra_fields=extra_fields,**kwargs) 115 if 'vizier' in searchables: 116 #-- first query catalogs that can only be queried via HD number 117 info = sesame.search(ID=ID,fix=True) 118 if 'alias' in info: 119 HDnumber = [name for name in info['alias'] if name[:2]=='HD'] 120 if HDnumber: 121 kwargs['master'] = vizier.get_photometry(extra_fields=extra_fields,constraints=['HD=%s'%(HDnumber[0][3:])],sources=['II/83/catalog','V/33/phot'],sort=None,**kwargs) 122 #-- then query catalogs that can only be queried via another catalog 123 results,units,comms = vizier.search('J/A+A/380/609/table1',ID=ID) 124 if results is not None: 125 catname = results[0]['Name'].strip() 126 kwargs['master'] = vizier.get_photometry(take_mean=True,extra_fields=extra_fields,constraints=['Name={0}'.format(catname)],sources=['J/A+A/380/609/table{0}'.format(tnr) for tnr in range(2,5)],sort=None,**kwargs) 127 #-- then query normal catalogs 128 kwargs['master'] = vizier.get_photometry(ID=ID,to_units=to_units,extra_fields=extra_fields,**kwargs) 129 130 131 if 'gcpd' in searchables: 132 kwargs['master'] = gcpd.get_photometry(ID=ID,to_units=to_units,extra_fields=extra_fields,**kwargs) 133 master = kwargs['master'] 134 135 #-- now make a summary of the contents: 136 photbands = [phot.split('.')[0] for phot in master['photband']] 137 contents = [(i,photbands.count(i)) for i in sorted(list(set(photbands)))] 138 for phot in contents: 139 logger.info('%10s: found %d measurements'%phot) 140 return master

141

142 -def add_bibcodes(master):

143 """ 144 Add bibcodes to a master record. 145 146 @param master: master record of photometry 147 @type master: numpy record array 148 @return: master record of photometry with additional column 149 @rtype: numpy record array 150 """ 151 bibcodes = [] 152 for source in master['source']: 153 bibcodes.append('None') 154 #-- check if it's a ViZieR catalog 155 if source.count('/')>=2: 156 #-- some catalogs have no ADS links but do have a bibcode, or the 157 # bibcode is manually put in the *cats_phot.cfg file 158 if 'bibcode' in vizier.cat_info.options(source): 159 bibcodes[-1] = vizier.cat_info.get(source,'bibcode') 160 #-- others may be retrieved immediately from the ViZieR site 161 else: 162 bibcodes[-1] = vizier.catalog2bibcode(source) 163 elif source=='GCPD': 164 bibcodes[-1] = '1997A&AS..124..349M' 165 elif source in gator.cat_info.sections(): 166 if 'bibcode' in gator.cat_info.options(source): 167 bibcodes[-1] = gator.cat_info.get(source,'bibcode') 168 elif source in mast.cat_info.sections(): 169 if 'bibcode' in mast.cat_info.options(source): 170 bibcodes[-1] = mast.cat_info.get(source,'bibcode') 171 if bibcodes[-1]=='None' or bibcodes[-1] is None: 172 logger.info('Bibcode of source {0:25s} not found'.format(source)) 173 bibcodes = np.array(bibcodes,str) 174 bibcodes = np.rec.fromarrays([bibcodes],names=['bibcode']) 175 return numpy_ext.recarr_join(master,bibcodes)

176 177

178 -def make_bibtex(master,filename):

179 """ 180 Make a bib file from a master record 181 """ 182 bibcodes = sorted(set(list(master['bibcode']))) 183 with open(filename,'w') as ff: 184 ff.write('% \citet{{{0}}}\n\n\n'.format(','.join(bibcodes))) 185 for bibcode in bibcodes: 186 try: 187 ff.write(vizier.bibcode2bibtex(bibcode)+'\n') 188 except IOError: 189 logger.info('Error retrieving bibtex for {0}'.format(bibcode))

190 191 192 193 194

195 -def xmatch(coords1,coords2):

196 """ 197 Crossmatch two sets of coordinates. 198 """ 199 tree = KDTree(coords1) 200 distance,order = tree.query(coords2) 201 202 match_success = distance<(tol/(60.)) 203 204 205 #-- first tuple contains matches (cat1,cat2) 206 # second tuple contains mismatches (cat1,cat2) 207 # so you if you call the function via 208 # >>> match,mismatch = xmatch(coord1,coord2) 209 # you probably want to do something like 210 # >>> matched_cat1,matched_cat2 = cat1[match[0]],cat2[match[1]] 211 # >>> mismatched_cat1,mismatched_cat2 = cat1[mismatch[0]],cat2[mismatch[1]] 212 # Then, you can just glue together these four parts: 213 # [ matched_cat1 matched_cat2 ] 214 # [ mismatched_cat1 nan ] 215 # [ nan mismatched_cat2 ] 216 return (order[match_success],match_success),(order[-match_success],-match_success)

217

218 -def photometry2str(master):

219 """ 220 String representation of master record array 221 222 @param master: master record array containing photometry 223 @type master: numpy record array 224 """ 225 master = master[np.argsort(master['photband'])] 226 txt = '#%19s %12s %12s %12s %10s %12s %12s %11s %s\n'%('PHOTBAND','MEAS','E_MEAS','UNIT','CWAVE','CMEAS','E_CMEAS','UNIT','SOURCE') 227 txt+= '#=========================================================================================================================\n' 228 for i,j,k,l,m,n,o,p in zip(master['photband'],master['meas'],master['e_meas'],master['unit'],master['cwave'],master['cmeas'],master['e_cmeas'],master['source']): 229 txt+='%20s %12g %12g %12s %10.0f %12g %12g erg/s/cm2/AA %s\n'%(i,j,k,l,m,n,o,p) 230 return txt

231 232 if __name__=="__main__": 233 logger = loggers.get_basic_logger() 234 235 method,args,kwargs = argkwargparser.parse() 236 print_help = '--help' in args or '-h' in args 237 if print_help: 238 help(globals()[method]) 239 else: 240 master = globals()[method](*args,**kwargs) 241 print photometry2str(master) 242 243 244 245 #xmatch_vizier('II/169/main','B/pastel/pastel') 246

Source Code for Module ivs.catalogs.crossmatch