Browse Source

Initial commit - Playing with star catalogs :)

main
Jennie Rose 6 months ago
commit
3afd7a6d96
6 changed files with 248 additions and 0 deletions
  1. +0
    -0
      galsim/__init__.py
  2. +89
    -0
      galsim/catalog.py
  3. +0
    -0
      galsim/model.py
  4. +67
    -0
      galsim/planet.py
  5. +44
    -0
      requirements.txt
  6. +48
    -0
      sandbox.py

+ 0
- 0
galsim/__init__.py View File


+ 89
- 0
galsim/catalog.py View File

@ -0,0 +1,89 @@
import astropy
import json
import numpy
import pyvo
from enum import Enum
class FieldSet(Enum):
ID = ('pl_name', 'hostname', 'pl_letter', 'gaia_id', 'hip_name', 'hd_name')
SYSTEM = ('sy_snum', 'sy_pnum', 'sy_mnum', 'cb_flag')
LOCATION = ('sy_dist', 'ra', 'dec', 'glat', 'glon')
STELLAR_DATA = ('st_mass', 'st_rad', 'st_dens', 'st_teff', 'st_age', 'st_met', 'st_lum')
PLANETARY_DATA = ('pl_bmasse', 'pl_rade', 'pl_bmassj', 'pl_radj', 'pl_dens')
ORBIT = ('pl_orbper', 'pl_orbsmax', 'pl_orbeccen', 'pl_orbincl')
MAGNITUDE = ('sy_bmag', 'sy_vmag')
class Catalog:
archive_url = 'https://exoplanetarchive.ipac.caltech.edu/TAP'
def __init__(self, max_distance, spacing=12.0):
'''
max_distance: maximum distance, in parsecs, of stars to include in query
spacing: distance, in parsecs, between concentric groups
'''
self.max_distance = max_distance
self.spacing = spacing
def create_query(self, *field_sets):
fields = []
for field_set in field_sets:
fields.extend([field for field in field_set.value])
query = [
'select',
', '.join(fields),
'from PSCompPars',
f'where sy_dist between 0 and {self.max_distance}',
'order by',
', '.join(FieldSet.LOCATION.value)
]
return ' '.join(query)
def fetch(self, *field_sets):
service = pyvo.dal.TAPService(self.archive_url)
query = service.create_query(self.create_query(*field_sets))
self.data = query.execute().to_table()
def lookup_index(self, planet_name):
if hasattr(self, 'data'):
name = self.data['pl_name']
return numpy.where(name == planet_name)
def lookup_row(self, planet_name):
if hasattr(self, 'data'):
return self.data[self.lookup_index(planet_name)]
@staticmethod
def distance_group(data, spacing=12.0):
return numpy.trunc(data['sy_dist'] / spacing)
@staticmethod
def latitude_subgroup(data, spacing=12.0):
ring = Catalog.distance_group(data, spacing)
latitude = numpy.pi * data['glat'] / 180.0
subgroup = numpy.zeros(ring.shape)
for i in range(len(subgroup)):
crit_angle = numpy.arcsin(ring[i] / (ring[i] + 1))
if crit_angle > 0.0:
bins = [numpy.pi, crit_angle, 0.0, -crit_angle, -numpy.pi]
subgroup[i] = numpy.digitize(latitude[i], bins)
else:
subgroup[i] = 0.0
return subgroup
def group_regions(self):
if hasattr(self, 'data'):
rings = self.distance_group(self.data, self.spacing)
return self.data.group_by(rings)
def subgroup_regions(self):
if hasattr(self, 'data'):
rings = self.distance_group(self.data, self.spacing)
tiers = self.latitude_subgroup(self.data, self.spacing)
key_table = astropy.table.Table([rings, tiers], names=('ring_group', 'tier_subgroup'))
return self.data.group_by(key_table)

+ 0
- 0
galsim/model.py View File


+ 67
- 0
galsim/planet.py View File

@ -0,0 +1,67 @@
import numpy
from galsim.catalog import FieldSet
regions = [
'Andes Line', 'Rocky Mountain Line', 'Kunlun Line',
'Ural Line', 'Atlas Line', 'Himalaya Line',
'Altai Line', 'Carpathian Line', 'Sierra Madre Line',
'Karakoram Line', 'Cascade Line', 'Annamite Line',
]
class Planet:
def __init__(self, row):
self._data = {
field: numpy.array(row[field]).item()
for field in row.colnames
}
self.name = self._data['pl_name']
self.distance = self._data['sy_dist']
self.latitude = self._data['glat']
self.longitude = self._data['glon']
@property
def coordinates(self):
'''
return cartesian coordinates [x, y, z]
x : positive axis toward the galactic core
y : positive axis away from galactic spin direction
z : positive axis "up" toward the galactic north pole
'''
el = numpy.pi * self.latitude / 180.0
az = numpy.pi * self.longitude / 180.0
x = self.distance * numpy.cos(az) * numpy.cos(el)
y = self.distance * numpy.sin(az) * numpy.cos(el)
z = self.distance * numpy.sin(el)
return numpy.array([x, y, z])
@property
def x(self):
return self.coordinates[0]
@property
def y(self):
return self.coordinates[1]
@property
def z(self):
return self.coordinates[2]
def get_data(self, field_set):
return {field: self._data[field] for field in field_set.value}
def to_dict(self, *field_sets):
out = {
'name': self.name,
'distance': self.distance,
'latitude': self.latitude,
'longitude': self.longitude,
'x': self.x,
'y': self.y,
'z': self.z,
}
for field_set in field_sets:
out.update(self.get_data(field_set))
return out

+ 44
- 0
requirements.txt View File

@ -0,0 +1,44 @@
astropy==4.2.1
astroquery==0.4.2
backcall==0.2.0
beautifulsoup4==4.9.3
certifi==2021.5.30
cffi==1.14.5
chardet==4.0.0
cryptography==3.4.7
decorator==5.0.9
html5lib==1.1
idna==2.10
importlib-metadata==4.4.0
ipykernel==5.5.5
ipython==7.24.1
ipython-genutils==0.2.0
jedi==0.18.0
jeepney==0.6.0
jupyter-client==6.1.12
jupyter-core==4.7.1
keyring==23.0.1
matplotlib-inline==0.1.2
mimeparse==0.1.3
numpy==1.20.3
parso==0.8.2
pexpect==4.8.0
pickleshare==0.7.5
prompt-toolkit==3.0.18
ptyprocess==0.7.0
pycparser==2.20
pyerfa==2.0.0
Pygments==2.9.0
python-dateutil==2.8.1
pyvo==1.1
pyzmq==22.1.0
requests==2.25.1
SecretStorage==3.3.1
six==1.16.0
soupsieve==2.2.1
tornado==6.1
traitlets==5.0.5
urllib3==1.26.5
wcwidth==0.2.5
webencodings==0.5.1
zipp==3.4.1

+ 48
- 0
sandbox.py View File

@ -0,0 +1,48 @@
import json
import numpy
import astropy
from galsim.planet import Planet
from galsim.catalog import Catalog, FieldSet
catalog = Catalog(144)
catalog.fetch(*FieldSet)
Catalog.latitude_subgroup(catalog.data)
rings = catalog.subgroup_regions()
for key, group in zip(rings.groups.keys, rings.groups):
print('==========')
print(key)
print(group[FieldSet.LOCATION.value])
print('==========')
minx, miny, minz = 0.0, 0.0, 0.0
maxx, maxy, maxz = 0.0, 0.0, 0.0
for row in catalog.data:
planet = Planet(row)
if planet.x < minx:
minx = planet.x
if planet.x > maxx:
maxx = planet.x
if planet.y < miny:
miny = planet.y
if planet.y > maxy:
maxy = planet.y
if planet.z < minz:
minz = planet.z
if planet.z > maxz:
maxz = planet.z
(minx, miny, minz)
(maxx, maxy, maxz)
row = catalog.lookup_row('tau Cet g')
planet = Planet(row)
planet.to_dict()

Loading…
Cancel
Save