ivs.catalogs.gator

1 # -*- coding: utf-8 -*- 2 """ 3 Interface to the GATOR search engine 4 """ 5 import os 6 import urllib 7 import logging 8 import ConfigParser 9 import numpy as np 10 11 from ivs.aux import loggers 12 from ivs.aux import numpy_ext 13 from ivs.inout import ascii 14 from ivs.sed import filters 15 from ivs.units import conversions 16 17 18 logger = logging.getLogger("CAT.GATOR") 19 logger.addHandler(loggers.NullHandler()) 20 21 22 basedir = os.path.dirname(os.path.abspath(__file__)) 23 24 #-- read in catalog information 25 cat_info = ConfigParser.ConfigParser() 26 cat_info.optionxform = str # make sure the options are case sensitive 27 cat_info.readfp(open(os.path.join(basedir,'gator_cats_phot.cfg'))) 28 29 30 #{ Basic interfaces 31

32 -def search(catalog,**kwargs):

33 """ 34 Search and retrieve information from a Gator catalog. 35 36 Two ways to search for data within a catalog C{name}: 37 38 1. You're looking for info on B{one target}, then give the target's 39 C{ID} or coordinates (C{ra} and C{dec}), and a search C{radius}. 40 41 2. You're looking for information of B{a whole field}, then give the 42 field's coordinates (C{ra} and C{dec}), and C{radius}. 43 44 If you have a list of targets, you need to loop this function. 45 46 If you supply a filename, the results will be saved to that path, and you 47 will get the filename back as received from urllib.URLopener (should be the 48 same as the input name, unless something went wrong). 49 50 If you don't supply a filename, you should leave C{filetype} to the default 51 C{tsv}, and the results will be saved to a temporary 52 file and deleted after the function is finished. The content of the file 53 will be read into a dictionary, as well as the units (two separate 54 dictionaries with the same keys, depending on the colum names in the 55 catalog). The entries in the dictionary are of type C{ndarray}, and will 56 be converted to a float-array if possible. If not, the array will consist 57 of strings. The comments are also returned as a list of strings. 58 59 60 @param catalog: name of a GATOR catalog (e.g. 'II/246/out') 61 @type catalog: str 62 @keyword filename: name of the file to write the results to (no extension) 63 @type filename: str 64 @return: filename / catalog data columns, units, comments 65 @rtype: str/ record array, dict, list of str 66 """ 67 filename = kwargs.pop('filename',None) # remove filename from kwargs 68 filetype = kwargs.setdefault('filetype','1') 69 if filename is not None and '.' in os.path.basename(filename): 70 filetype = os.path.splitext(filename)[1][1:] 71 elif filename is not None: 72 filename = '%s.%s'%(filename,filetype) 73 74 #-- gradually build URI 75 base_url = _get_URI(catalog,**kwargs) 76 #-- prepare to open URI 77 url = urllib.URLopener() 78 filen,msg = url.retrieve(base_url,filename=filename) 79 # maybe we are just interest in the file, not immediately in the content 80 if filename is not None: 81 logger.info('Querying GATOR source %s and downloading to %s'%(catalog,filen)) 82 url.close() 83 return filen 84 85 # otherwise, we read everything into a dictionary 86 if filetype=='1': 87 try: 88 results,units,comms = txt2recarray(filen) 89 #-- raise an exception when multiple catalogs were specified 90 except ValueError: 91 raise ValueError, "failed to read %s, perhaps multiple catalogs specified (e.g. III/168 instead of III/168/catalog)"%(catalog) 92 url.close() 93 logger.info('Querying GATOR source %s (%d)'%(catalog,(results is not None and len(results) or 0))) 94 return results,units,comms

95 96 97 98

99 -def list_catalogs():

100 """ 101 Return a list of all availabe GATOR catalogues. 102 103 @return: list of gator catalogs and discriptions 104 @rtype: list of string tuples 105 """ 106 url = urllib.URLopener() 107 filen,msg = url.retrieve('http://irsa.ipac.caltech.edu/cgi-bin/Gator/nph-scan?mode=ascii') 108 results,units,comms = txt2recarray(filen) 109 cats = [] 110 for i in range(len(results)): 111 cats += [(results['catname'][i],results['description'][i])] 112 logger.info('%-20s: %s'%cats[-1]) 113 url.close() 114 return cats

115 116 117 118 119

120 -def get_photometry(ID=None,extra_fields=['dist','ra','dec'],**kwargs):

121 """ 122 Download all available photometry from a star to a record array. 123 124 For extra kwargs, see L{_get_URI} and L{gator2phot} 125 126 Example usage: 127 128 >>> import pylab 129 >>> import vizier 130 >>> name = 'kr cam' 131 >>> master = vizier.get_photometry(name,to_units='erg/s/cm2/AA',extra_fields=[]) 132 >>> master = get_photometry(name,to_units='erg/s/cm2/AA',extra_fields=[],master=master) 133 >>> p = pylab.figure() 134 >>> wise = np.array(['WISE' in photband and True or False for photband in master['photband']]) 135 >>> p = pylab.errorbar(master['cwave'],master['cmeas'],yerr=master['e_cmeas'],fmt='ko') 136 >>> p = pylab.errorbar(master['cwave'][wise],master['cmeas'][wise],yerr=master['e_cmeas'][wise],fmt='ro',ms=8) 137 >>> p = pylab.gca().set_xscale('log') 138 >>> p = pylab.gca().set_yscale('log') 139 >>> p = pylab.show() 140 141 Other examples: 142 >>> master = get_photometry(ra=71.239527,dec=-70.589427,to_units='erg/s/cm2/AA',extra_fields=[],radius=1.) 143 >>> master = get_photometry(ID='J044458.39-703522.6',to_units='W/m2',extra_fields=[],radius=1.) 144 """ 145 kwargs['ID'] = ID 146 to_units = kwargs.pop('to_units','erg/s/cm2/AA') 147 master_ = kwargs.get('master',None) 148 master = None 149 #-- retrieve all measurements 150 for source in cat_info.sections(): 151 results,units,comms = search(source,**kwargs) 152 if results is not None: 153 master = gator2phot(source,results,units,master,extra_fields=extra_fields) 154 155 #-- convert the measurement to a common unit. 156 if to_units and master is not None: 157 #-- prepare columns to extend to basic master 158 dtypes = [('cwave','f8'),('cmeas','f8'),('e_cmeas','f8'),('cunit','a50')] 159 cols = [[],[],[],[]] 160 #-- forget about 'nan' errors for the moment 161 no_errors = np.isnan(master['e_meas']) 162 master['e_meas'][no_errors] = 0. 163 #-- extend basic master 164 zp = filters.get_info(master['photband']) 165 for i in range(len(master)): 166 try: 167 value,e_value = conversions.convert(master['unit'][i],to_units,master['meas'][i],master['e_meas'][i],photband=master['photband'][i]) 168 except ValueError: # calibrations not available 169 value,e_value = np.nan,np.nan 170 except AssertionError: # the error or flux must be positive number 171 value,e_value = np.nan,np.nan 172 try: 173 eff_wave = filters.eff_wave(master['photband'][i]) 174 except IOError: 175 eff_wave = np.nan 176 cols[0].append(eff_wave) 177 cols[1].append(value) 178 cols[2].append(e_value) 179 cols[3].append(to_units) 180 master = numpy_ext.recarr_addcols(master,cols,dtypes) 181 #-- reset errors 182 master['e_meas'][no_errors] = np.nan 183 master['e_cmeas'][no_errors] = np.nan 184 185 if master_ is not None and master is not None: 186 master = numpy_ext.recarr_addrows(master_,master.tolist()) 187 elif master_ is not None: 188 master = master_ 189 190 #-- and return the results 191 return master

192 193 194 #} 195 196 #{ Convenience functions 197 198

199 -def txt2recarray(filename):

200 """ 201 Read a GATOR outfmt=1__ (ASCII) file into a record array. 202 203 @param filename: name of the TSV file 204 @type filename: str 205 @return: catalog data columns, units, comments 206 @rtype: record array, dict, list of str 207 """ 208 ff = open(filename,'r') 209 data = [] 210 comms = [] 211 indices = None 212 while 1: # might call read several times for a file 213 line = ff.readline() 214 if not line: break # end of file 215 216 #-- strip return character from line 217 if line.isspace(): 218 continue # empty line 219 220 #-- remove return characters 221 line = line.replace('\n','') 222 #-- when reading a comment line 223 if line[0] =='\\': 224 continue # treat next line 225 if line[0] == '|': 226 comms.append(line) 227 indices = [] 228 for i in range(len(line)): 229 if line[i]=='|': indices.append(i) 230 continue 231 if indices is None: break 232 data.append([]) 233 for i in range(len(indices)-1): 234 data[-1].append(line[indices[i]:indices[i+1]].strip()) 235 results = None 236 units = {} 237 #-- retrieve the data and put it into a record array 238 if comms: 239 #-- now convert this thing into a nice dictionary 240 data = np.array(data) 241 #-- retrieve the format of the columns. They are given in the 242 # Fortran format. In rare cases, columns contain multiple values 243 # themselves (so called vectors). In those cases, we interpret 244 # the contents as a long string 245 names = [head.strip() for head in comms[0].split('|')[1:-1]] 246 formats = comms[1] 247 formats = formats.replace('double','f8') 248 formats = formats.replace('char','a100') 249 formats = formats.replace('int','f8') 250 formats = formats.replace('long','f8') 251 formats = [head.strip() for head in formats.split('|')[1:-1]] 252 units_ = [head.strip() for head in comms[2].split('|')[1:-1]] 253 #-- define dtypes for record array 254 dtypes = np.dtype([(i,j) for i,j in zip(names,formats)]) 255 #-- remove spaces or empty values 256 cols = [] 257 for i,key in enumerate(names): 258 units[key] = units_[i] 259 if len(data)==0: break 260 col = data[:,i] 261 #-- fill empty or null values with nan 262 col = [(row.isspace() or not row or row=='null') and 'nan' or row for row in col] 263 cols.append(col) 264 265 #-- define columns for record array and construct record array 266 if len(data)==0: 267 results = None 268 else: 269 cols = np.array(cols) 270 cols = [np.cast[dtypes[i]](cols[i]) for i in range(len(cols))] 271 results = np.rec.array(cols, dtype=dtypes) 272 return results,units,comms

273 274 275

276 -def gator2phot(source,results,units,master=None,extra_fields=['_r','_RAJ2000','_DEJ2000']):

277 """ 278 Convert/combine Gator record arrays to measurement record arrays. 279 280 Every line in the combined array represents a measurement in a certain band. 281 282 This is probably only useful if C{results} contains only information on 283 one target (or you have to give 'ID' as an extra field, maybe). 284 285 The standard columns are: 286 287 1. C{meas}: containing the photometric measurement 288 2. C{e_meas}: the error on the photometric measurement 289 3. C{flag}: an optional quality flag 290 4. C{unit}: the unit of the measurement 291 5. C{photband}: the photometric passband (FILTER.BAND) 292 6. C{source}: name of the source catalog 293 294 You can add extra information from the Gator catalog via the list of keys 295 C{extra_fields}. 296 297 If you give a C{master}, the information will be added to a previous 298 record array. If not, a new master will be created. 299 300 The result is a record array with each row a measurement. 301 302 @param source: name of the Gator source 303 @type source: str 304 @param results: results from Gator C{search} 305 @type results: record array 306 @param units: header of Gator catalog with key name the column name and 307 key value the units of that column 308 @type units: dict 309 @param master: master record array to add information to 310 @type master: record array 311 @param extra_fields: any extra columns you want to add information from 312 @type extra_fields: list of str 313 @return: array with each row a measurement 314 @rtype: record array 315 """ 316 e_flag = 'e_' 317 q_flag = 'q_' 318 #-- basic dtypes 319 dtypes = [('meas','f8'),('e_meas','f8'),('flag','a20'), 320 ('unit','a30'),('photband','a30'),('source','a50')] 321 322 #-- extra can be added: 323 names = list(results.dtype.names) 324 translation = {'_r':'dist','_RAJ2000':'ra','_DEJ2000':'dec'} 325 if extra_fields is not None: 326 for e_dtype in extra_fields: 327 dtypes.append((e_dtype,results.dtype[names.index(translation[e_dtype])].str)) 328 329 #-- create empty master if not given 330 newmaster = False 331 if master is None or len(master)==0: 332 master = np.rec.array([tuple([('f' in dt[1]) and np.nan or 'none' for dt in dtypes])],dtype=dtypes) 333 newmaster = True 334 335 #-- add fluxes and magnitudes to the record array 336 cols_added = 0 337 for key in cat_info.options(source): 338 if key[:2] in [e_flag,q_flag] or key =='bibcode': 339 continue 340 photband = cat_info.get(source,key) 341 #-- contains measurement, error, quality, units, photometric band and 342 # source 343 cols = [results[key][:1], 344 (e_flag+key in cat_info.options(source) and results[cat_info.get(source,e_flag+key)][:1] or np.ones(len(results[:1]))*np.nan), 345 (q_flag+key in cat_info.options(source) and results[cat_info.get(source,q_flag+key)][:1] or np.ones(len(results[:1]))*np.nan), 346 np.array(len(results[:1])*[units[key]],str), 347 np.array(len(results[:1])*[photband],str), 348 np.array(len(results[:1])*[source],str)] 349 #-- add any extra fields if desired. 350 if extra_fields is not None: 351 for e_dtype in extra_fields: 352 cols += [translation[e_dtype] in results.dtype.names and results[translation[e_dtype]][:1] or np.ones(len(results[:1]))*np.nan] 353 #-- add to the master 354 rows = [] 355 for i in range(len(cols[0])): 356 rows.append(tuple([col[i] for col in cols])) 357 master = np.core.records.fromrecords(master.tolist()+rows,dtype=dtypes) 358 cols_added += 1 359 360 #-- skip first line from building 361 if newmaster: master = master[1:] 362 return master

363 364 #} 365 366 #{ Internal helper functions 367

368 -def _get_URI(name,ID=None,ra=None,dec=None,radius=1.,filetype='1',spatial='cone',**kwargs):

369 """ 370 Build GATOR URI from available options. 371 372 Filetype should be one of: 373 374 6___ returns a program interface in XML 375 3___ returns a VO Table (XML) 376 2___ returns SVC (Software handshaking structure) message 377 1___ returns an ASCII table 378 0___ returns Gator Status Page in HTML (default) 379 380 kwargs are to catch unused arguments. 381 382 @param name: name of a GATOR catalog (e.g. 'wise_prelim_p3as_psd') 383 @type name: str 384 @keyword ID: target name 385 @type ID: str 386 @param filetype: type of the retrieved file 387 @type filetype: str (one of '0','1','2','3','6'... see GATOR site) 388 @keyword radius: search radius (arcseconds) 389 @type radius: float 390 @param ra: target's right ascension 391 @type ra: float 392 @param dec: target's declination 393 @type dec: float 394 @param spatial: type of spatial query 395 @type spatial: str (one of Cone, Box (NIY), Polygon (NIY), NONE (NIY)) 396 @return: url 397 @rtype: str 398 """ 399 base_url = 'http://irsa.ipac.caltech.edu/cgi-bin/Gator/nph-query?' 400 base_url += 'catalog=%s'%(name) 401 #base_url += '&spatial=cone' 402 403 #-- in GATOR, many things should be given via the keyword 'objstr': right 404 # ascension, declination, target name... 405 objstr = None 406 407 #-- transform input to the 'objstr' paradigm 408 if ID is not None: 409 #-- if the ID is given in the form 'J??????+??????', derive the 410 # coordinates of the target from the name. 411 if ID[0]=='J': 412 ra = int(ID[1:3]),int(ID[3:5]),float(ID[5:10]) 413 dec = int(ID[10:13]),int(ID[13:15]),float(ID[15:]) 414 ra = '%02d+%02d+%.2f'%ra 415 dec = '+%+02d+%02d+%.2f'%dec 416 objstr = ra+dec 417 else: 418 objstr = ID 419 420 #-- this is when ra and dec are given 421 if ra is not None and dec is not None: 422 ra = str(ra) 423 dec = str(dec) 424 ra = ra.replace(' ','+').replace(':','+') 425 dec = dec.replace(' ','+').replace(':','+') 426 objstr = '+'.join([ra,dec]) 427 428 #-- build up the URI 429 if 'objstr' is not None: base_url += '&objstr=%s'%(objstr) 430 if 'filetype' is not None: base_url += '&outfmt=%s'%(filetype) 431 if 'radius' is not None: base_url += '&radius=%s'%(radius) 432 433 logger.debug(base_url) 434 435 return base_url

436 437 438 #} 439 440 if __name__=="__main__": 441 import vizier 442 logger = loggers.get_basic_logger("") 443 #-- example 1 444 #master = get_photometry(ra=71.239527,dec=-70.589427,to_units='erg/s/cm2/AA',extra_fields=[],radius=1.) 445 #master = vizier.get_photometry(ra=71.239527,dec=-70.589427,to_units='erg/s/cm2/AA',extra_fields=[],radius=5.,master=master) 446 447 #-- example 2 448 #master = get_photometry(ID='J044458.39-703522.6',to_units='W/m2',extra_fields=[],radius=1.) 449 #master = vizier.get_photometry(ID='J044458.39-703522.6',to_units='W/m2',extra_fields=[],radius=5.,master=master) 450 451 452 #-- example 3 453 #master = get_photometry(ID='HD43317',to_units='W/m2',extra_fields=[],radius=1.) 454 #master = vizier.get_photometry(ID='HD43317',to_units='W/m2',extra_fields=[],radius=5.,master=master) 455 456 #-- example 4 457 #master = get_photometry(ID='HD143454',to_units='erg/s/cm2/AA',extra_fields=[],radius=1.) 458 #master = vizier.get_photometry(ID='HD143454',to_units='erg/s/cm2/AA',extra_fields=[],radius=30.,master=master) 459 #-- example 5 460 master = get_photometry(ID='RR Aql',to_units='erg/s/cm2/AA',extra_fields=[],radius=1.) 461 master = vizier.get_photometry(ID='RR Aql',to_units='erg/s/cm2/AA',extra_fields=[],radius=30.,master=master) 462 463 464 print master 465 466 from pylab import * 467 figure() 468 #loglog(master['cwave'],master['cmeas'],'ko') 469 #errorbar(master['cwave'],master['cmeas'],yerr=master['e_cmeas'],fmt='ro') 470 plot(np.log10(master['cwave']),np.log10(master['cmeas']),'ko') 471 #errorbar(master['cwave'],master['cmeas'],yerr=master['e_cmeas'],fmt='ro') 472 #plot(np.log10(1538.62),np.log10(conversions.convert('muJy','W/m2',1202.758,photband='GALEX.FUV'))+13,'bo') 473 #plot(np.log10(2315.66),np.log10(conversions.convert('muJy','W/m2',2896.989,photband='GALEX.NUV'))+13,'bo') 474 #plot(np.log10(1538.62),np.log10(conversions.convert('muJy','W/m2',5025.582,photband='GALEX.FUV'))+13,'go') 475 #plot(np.log10(2315.66),np.log10(conversions.convert('muJy','W/m2',7106.731,photband='GALEX.NUV'))+13,'go') 476 #ylim(-1.5,2.5) 477 #xlim(3.1,4.7) 478 479 480 show() 481 482 sys.exit() 483 import doctest 484 doctest.testmod() 485 #logger = loggers.get_basic_logger() 486 #catalogs = ['wise_prelim_p3as_psd'] 487 #import vizier 488 #import pylab 489 ##list_gator_catalogs() 490 ##sys.exit() 491 #name = 'kr cam' 492 #master = vizier.get_photometry(name,to_units='erg/s/cm2/AA',extra_fields=[]) 493 #master = get_photometry(name,to_units='erg/s/cm2/AA',extra_fields=[],master=master) 494 495 #pylab.figure() 496 #wise = np.array(['WISE' in photband and True or False for photband in master['photband']]) 497 #pylab.errorbar(master['cwave'],master['cmeas'],yerr=master['e_cmeas'],fmt='ko') 498 #pylab.errorbar(master['cwave'][wise],master['cmeas'][wise],yerr=master['e_cmeas'][wise],fmt='ro',ms=8) 499 #pylab.gca().set_xscale('log') 500 #pylab.gca().set_yscale('log') 501 #pylab.show() 502

Source Code for Module ivs.catalogs.gator