Skip to content

Commit

Permalink
Feature: Added command line option to copy the photometry catalogs fi…
Browse files Browse the repository at this point in the history 
…les to the SPEEDYFIT_MODELS directory so the user can adapt them. Adapted photometry_query.py to first try to load the photometry catalogs from the SPEEDYFIT_MODELS directory, and only when they don't exist there, use the catalogs shipped with the package.
  • Loading branch information
Joris Vos committed Feb 19, 2022
1 parent bfa3f4d commit b428cb0
Show file tree
Hide file tree
Showing 3 changed files with 110 additions and 37 deletions.
16 changes: 16 additions & 0 deletions speedyfit/main.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -435,10 +435,20 @@ def perform_fit(args):


def check_grids(args):
"""
Run the check model grids function and report the results to the user.
"""
print_bands = args.bands

model.check_grids(print_bands=print_bands)

def copy_catalogs(args):
"""
Run the copy_photometry_catalogs function so that the user can adapt them to their preferences.
"""
photometry_query.copy_photometry_catalogs(overwrite=args.overwrite)


def main():
parser = argparse.ArgumentParser(description="Speedyfit: obtaining and fitting photometric SEDs")

Expand Down Expand Up @@ -485,6 +495,12 @@ def main():
help="List the photomtric bands included in the integrated lists.")
grid_parser.set_defaults(func=check_grids)

# --copy photometry catalogs--
grid_parser = subparsers.add_parser('copy_catalogs', help='Copy the photometry catalogs files to the directory set '
'in the SPEEDYFIT_MODELS environment variable.')
grid_parser.add_argument('--overwrite', dest='overwrite', action='store_true',
help="When set already existing catalogs files in SPEEDYFIT_MODELS will be overwritten.")
grid_parser.set_defaults(func=copy_catalogs)

args = parser.parse_args()
args.func(args)
Expand Down
126 changes: 89 additions & 37 deletions speedyfit/photometry_query.py
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@

import os

import shutil
import logging
import configparser
import warnings

Expand All @@ -23,20 +24,71 @@
from zero_point import zpt #gaiadr3-zeropoint
zpt.load_tables()

filedir = os.path.dirname(os.path.abspath(__file__))
FILE_DIR = os.path.dirname(os.path.abspath(__file__))
VIZIER_CAT = 'vizier_cats_phot.cfg'
TAP_CAT = 'tap_cats_phot.cfg'

# always request the distance of the observations to the search coordinates and sort on increasing distance
v = Vizier(columns=["*", '+_r'])

#-- read in catalog information
viz_info = configparser.ConfigParser()
viz_info.optionxform = str # make sure the options are case sensitive
viz_info.readfp(open(filedir+'/vizier_cats_phot.cfg'))
def load_photometry_catalogs():
"""
Function to read the photometry catalogs: vizier_cats_phot.cfg and tap_cats_phot.cfg.
The function will first look in the directory set with the SPEEDYFIT_MODELS environment variable. If no catalog
files can be found in that directory, it will use the catalogs shipped with the speedyfit package.
"""
user_dir = os.environ.get('SPEEDYFIT_MODELS', None)

if user_dir is not None and os.path.exists(os.path.join(user_dir, VIZIER_CAT)):
vizier_file = os.path.join(user_dir, VIZIER_CAT)
else:
vizier_file = os.path.join(FILE_DIR, VIZIER_CAT)

if user_dir is not None and os.path.exists(os.path.join(user_dir, TAP_CAT)):
tap_file = os.path.join(user_dir, TAP_CAT)
else:
tap_file = os.path.join(FILE_DIR, TAP_CAT)

print(f"Going to read Vizier photometry catalog from {vizier_file}.")
viz_info = configparser.ConfigParser()
viz_info.optionxform = str # make sure the options are case sensitive
with open(vizier_file) as ifile:
viz_info.readfp(ifile)

print(f"Going to read TAP photometry catalog from {tap_file}.")
tap_info = configparser.ConfigParser()
tap_info.optionxform = str # make sure the options are case sensitive
with open(tap_file) as ifile:
tap_info.readfp(ifile)

return viz_info, tap_info
VIZ_INFO, TAP_INFO = load_photometry_catalogs()

def copy_photometry_catalogs(overwrite=False):
"""
Copy the two photometry catalogs: vizier_cats_phot.cfg and tap_cats_phot.cfg from the package source location
to the SPEEDYFIT_MODELS directory, so that the user can edit them.
"""
if os.environ.get('SPEEDYFIT_MODELS', None) is None:
print("SPEEDYFIT_MODELS env variable is not set, can not copy the photometry catalogs.")
return
else:
destination = os.environ.get('SPEEDYFIT_MODELS', None)

vis_cat = os.path.join(FILE_DIR, VIZIER_CAT)
if os.path.exists(vis_cat) and not overwrite:
print(f"The file {vis_cat} already exists. To overwrite existing files use option --overwrite")
elif (os.path.exists(vis_cat) and overwrite) or not os.path.exists(vis_cat):
shutil.copy(vis_cat, destination)
print(f"Copied {vis_cat} to {destination}")

tap_info = configparser.ConfigParser()
tap_info.optionxform = str # make sure the options are case sensitive
tap_info.readfp(open(filedir+'/tap_cats_phot.cfg'))
tap_cat = os.path.join(FILE_DIR, TAP_CAT)
if os.path.exists(tap_cat) and not overwrite:
print(f"The file {tap_cat} already exists. To overwrite existing files use option --overwrite")
elif (os.path.exists(tap_cat) and overwrite) or not os.path.exists(tap_cat):
shutil.copy(tap_cat, destination)
print(f"Copied {tap_cat} to {destination}")


def get_coordinate(objectname):
Expand Down Expand Up @@ -66,7 +118,7 @@ def get_coordinate(objectname):

def get_vizier_photometry(objectname, radius=5):

catalogs = viz_info.sections()
catalogs = VIZ_INFO.sections()

data = v.query_object(objectname, catalog=catalogs, radius=radius*u.arcsec)

Expand All @@ -76,31 +128,31 @@ def get_vizier_photometry(objectname, radius=5):

distance = (data[catalog]['_r'][0] * u.arcmin).to(u.arcsec).value

if viz_info.has_option(catalog, 'bibcode'):
bibcode = viz_info.get(catalog, 'bibcode')
if VIZ_INFO.has_option(catalog, 'bibcode'):
bibcode = VIZ_INFO.get(catalog, 'bibcode')
else:
bibcode = '-'

print(catalog)

for band in viz_info.options(catalog):
for band in VIZ_INFO.options(catalog):

if '_unit' in band: continue
if '_err' in band: continue
if 'bibcode' in band: continue
if '-' in band: continue

bandname = viz_info.get(catalog, band)
bandname = VIZ_INFO.get(catalog, band)
value = data[catalog][band][0]

errkw = viz_info.get(catalog, band+'_err') if viz_info.has_option(catalog, band+'_err') else 'e_'+band
errkw = VIZ_INFO.get(catalog, band + '_err') if VIZ_INFO.has_option(catalog, band + '_err') else 'e_' + band
try:
err = data[catalog][errkw][0]
except:
err = 0.05

if viz_info.has_option(catalog, band+'_unit'):
unit = viz_info.get(catalog, band+'_unit')
if VIZ_INFO.has_option(catalog, band + '_unit'):
unit = VIZ_INFO.get(catalog, band + '_unit')
else:
unit = data[catalog][band].unit

Expand All @@ -120,29 +172,29 @@ def tap_query_vo(ra, dec, catalog):

service = pyvo.dal.TAPService(catalog)

table = tap_info.get(catalog, 'table')
table = TAP_INFO.get(catalog, 'table')
print(table)

keywords = ""
bands = []
for band in tap_info.options(catalog):
for band in TAP_INFO.options(catalog):
if 'table' in band: continue
if 'rakw' in band: continue
if 'deckw' in band: continue
if '_unit' in band: continue
if '_err' in band: continue
if 'bibcode' in band: continue

errkw = tap_info.get(catalog, band + '_err') if tap_info.has_option(catalog, band + '_err') else 'e_' + band
errkw = TAP_INFO.get(catalog, band + '_err') if TAP_INFO.has_option(catalog, band + '_err') else 'e_' + band

keywords += band + ", " + errkw + ", "
bands.append(band)

if keywords[-2:] == ", ":
keywords = keywords[0:-2]

rakw = tap_info.get(catalog, 'rakw') if tap_info.has_option(catalog, 'rakw') else 'raj2000'
deckw = tap_info.get(catalog, 'deckw') if tap_info.has_option(catalog, 'deckw') else 'dej2000'
rakw = TAP_INFO.get(catalog, 'rakw') if TAP_INFO.has_option(catalog, 'rakw') else 'raj2000'
deckw = TAP_INFO.get(catalog, 'deckw') if TAP_INFO.has_option(catalog, 'deckw') else 'dej2000'

# -- first try to query with distance
query = """SELECT
Expand All @@ -167,19 +219,19 @@ def tap_query_vo(ra, dec, catalog):
c1 = SkyCoord(ra * u.degree, dec * u.degree)
distance = SkyCoord(results['ra'][0] * u.degree, results['de'][0] * u.degree).separation(c1).to(u.arcsec).value

bibcode = tap_info.get(catalog, 'bibcode')
bibcode = TAP_INFO.get(catalog, 'bibcode')

# -- get the photometric measurements, errors and units
photometry = []
for band in bands:
bandname = tap_info.get(catalog, band)
bandname = TAP_INFO.get(catalog, band)
value = results[band][0]

errkw = tap_info.get(catalog, band + '_err') if tap_info.has_option(catalog, band + '_err') else 'e_' + band
errkw = TAP_INFO.get(catalog, band + '_err') if TAP_INFO.has_option(catalog, band + '_err') else 'e_' + band
err = results[errkw][0]

if tap_info.has_option(catalog, band + '_unit'):
unit = tap_info.get(catalog, band + '_unit')
if TAP_INFO.has_option(catalog, band + '_unit'):
unit = TAP_INFO.get(catalog, band + '_unit')
else:
unit = results[band].unit

Expand All @@ -195,29 +247,29 @@ def tap_query(ra, dec, catalog):

tp = TapPlus(url=catalog)

table = tap_info.get(catalog, 'table')
table = TAP_INFO.get(catalog, 'table')
print(table)

keywords = ""
bands = []
for band in tap_info.options(catalog):
for band in TAP_INFO.options(catalog):
if 'table' in band: continue
if 'rakw' in band: continue
if 'deckw' in band: continue
if '_unit' in band: continue
if '_err' in band: continue
if 'bibcode' in band: continue

errkw = tap_info.get(catalog, band+'_err') if tap_info.has_option(catalog, band+'_err') else 'e_'+band
errkw = TAP_INFO.get(catalog, band + '_err') if TAP_INFO.has_option(catalog, band + '_err') else 'e_' + band

keywords += band + ", " + errkw + ", "
bands.append(band)

if keywords[-2:] == ", ":
keywords = keywords[0:-2]

rakw = tap_info.get(catalog, 'rakw') if tap_info.has_option(catalog, 'rakw') else 'raj2000'
deckw = tap_info.get(catalog, 'deckw') if tap_info.has_option(catalog, 'deckw') else 'dej2000'
rakw = TAP_INFO.get(catalog, 'rakw') if TAP_INFO.has_option(catalog, 'rakw') else 'raj2000'
deckw = TAP_INFO.get(catalog, 'deckw') if TAP_INFO.has_option(catalog, 'deckw') else 'dej2000'

try:
#-- first try to query with distance
Expand Down Expand Up @@ -262,19 +314,19 @@ def tap_query(ra, dec, catalog):
distance = SkyCoord(results['ra'][0] * u.degree, results['de'][0] * u.degree).separation(c1).to(u.arcsec).value


bibcode = tap_info.get(catalog, 'bibcode')
bibcode = TAP_INFO.get(catalog, 'bibcode')

#-- get the photometric measurements, errors and units
photometry = []
for band in bands:
bandname = tap_info.get(catalog, band)
bandname = TAP_INFO.get(catalog, band)
value = results[band][0]

errkw = tap_info.get(catalog, band+'_err') if tap_info.has_option(catalog, band+'_err') else 'e_'+band
errkw = TAP_INFO.get(catalog, band + '_err') if TAP_INFO.has_option(catalog, band + '_err') else 'e_' + band
err = results[errkw][0]

if tap_info.has_option(catalog, band+'_unit'):
unit = tap_info.get(catalog, band+'_unit')
if TAP_INFO.has_option(catalog, band + '_unit'):
unit = TAP_INFO.get(catalog, band + '_unit')
else:
unit = results[band].unit

Expand All @@ -287,7 +339,7 @@ def tap_query(ra, dec, catalog):


def get_tap_photometry(ra, dec):
catalogs = tap_info.sections()
catalogs = TAP_INFO.sections()

dtypes = [('band', '<U20'), ('meas', 'f8'), ('emeas', 'f8'), ('unit', '<U10'), ('distance', 'f8'), ('bibcode', '<U20')]
photometry = np.array([], dtype=dtypes)
Expand Down
5 changes: 5 additions & 0 deletions speedyfit/tests/test_model.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,13 +17,18 @@ class TestCheckGrids:

@pytest.fixture(scope='class')
def make_modeldir(self, tmpdir_factory):
org_dir = os.environ['SPEEDYFIT_MODELS']

models_dir = tmpdir_factory.mktemp('SED_models')
models_dir.join('grid_description.yaml').write(grid_description_ex)
models_dir.join('kurucz93_z0.0_k2odfnew_sed.fits').write('test')
models_dir.join('ikurucz93_z0.0_k2odfnew_sed_lawfitzpatrick2004_Rv3.10.fits').write('test')

yield models_dir

os.environ['SPEEDYFIT_MODELS'] = org_dir
reload(model)

def test_check_grids__model_dir_without_dash(self, make_modeldir):
models_dir = make_modeldir

Expand Down

0 comments on commit b428cb0

Please sign in to comment.