Skip to content

Commit

Permalink
Update based on latest comments
Browse files Browse the repository at this point in the history
  • Loading branch information
Christoph Toennis committed Oct 22, 2024
1 parent 0d6c5a4 commit 89f4e37
Show file tree
Hide file tree
Showing 2 changed files with 74 additions and 117 deletions.
3 changes: 2 additions & 1 deletion docs/changes/2625.features.rst
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
Add ''get_bright_stars_with_motion'' to ctapipe.utils.astro
Modifies ''get_bright_stars'' in ctapipe.utils.astro to include proper motion
Add ''get_hipparcos_stars'' to use the Hippoarcos catalog to get star positions with proper motion
this function allows to look up the position of stars
around a specified location including proper motion.
different catalogues can be used.
188 changes: 72 additions & 116 deletions src/ctapipe/utils/astro.py
Original file line number Diff line number Diff line change
Expand Up @@ -7,24 +7,19 @@

import logging
from enum import Enum
from pathlib import Path

import astropy.units as u
from astropy.coordinates import Angle, SkyCoord
from astropy.table import Table
from astropy.time import Time

log = logging.getLogger("main")

__all__ = ["get_bright_stars"]


CACHE_FILE = Path("~/.psf_stars.ecsv").expanduser()
__all__ = ["get_bright_stars", "get_hipparcos_stars"]


class StarCatalogues(Enum):
Yale = "V/50/catalog" # Yale bright star catalogue
Hippoarcos = "I/239/hip_main" # hipparcos catalogue
Yale = "V/50/catalog" #: Yale bright star catalogue
Hippoarcos = "I/239/hip_main" #: hipparcos catalogue


def select_stars(stars, pointing=None, radius=None, Bmag_cut=None, Vmag_cut=None):
Expand Down Expand Up @@ -72,23 +67,28 @@ def select_stars(stars, pointing=None, radius=None, Bmag_cut=None, Vmag_cut=None
return stars


def get_star_catalog(catalog, filename):
def get_star_catalog(catalog):
"""
Utility function to download a star catalog for the get_bright_stars function
Parameters
----------
catalog: string
Name of the catalog to be used. Usable names are found in the Enum StarCatalogues. Default: Yale
filename: string
Name and location of the catalog file to be produced
Returns
----------
Astropy table:
List of all stars after cuts with catalog numbers, magnitudes,
and coordinates as SkyCoord objects including proper motion
"""
from astroquery.vizier import Vizier

vizier = Vizier(
catalog=StarCatalogues[catalog].value,
columns=[
"recno",
"HIP", #: HIP is the Hippoarcos ID available for that catalog
"HR", #: HR is the Harvard Revised Number available for the Yale catalog
"RAJ2000",
"DEJ2000",
"RAICRS",
Expand All @@ -97,138 +97,94 @@ def get_star_catalog(catalog, filename):
"pmDE",
"Vmag",
"Bmag",
"BTmag",
"VTmag",
],
row_limit=1000000,
)

stars = vizier.query_constraints(Vmag="0.0..10.0")[0]

if "Bmag" not in stars.keys():
log.exception("The chosen catalog does not have Bmag data")
if "Vmag" not in stars.keys():
log.exception("The chosen catalog does not have Vmag data")
stars.meta["Catalog"] = StarCatalogues[catalog].value

stars.write(CACHE_FILE, overwrite=True)
return stars


def get_bright_stars(pointing=None, radius=None, magnitude_cut=None):
"""
Get an astropy table of bright stars from the Yale bright star catalog
Parameters
----------
pointing: astropy Skycoord
pointing direction in the sky (if none is given, full sky is returned)
radius: astropy angular units
Radius of the sky region around pointing position. Default: full sky
magnitude_cut: float
Return only stars above a given absolute magnitude. Default: None (all entries)
def get_bright_stars(
pointing=None, radius=None, Bmag_cut=None, Vmag_cut=None, catalog="Yale"
): # max_magnitude):
Returns
-------
Astropy table:
List of all stars after cuts with catalog numbers, magnitudes,
and coordinates as SkyCoord objects including proper motion
"""
Get an astropy table of bright stars from a VizieR catalog

You can browse the catalogs at https://vizier.cds.unistra.fr/viz-bin/VizieR
from ctapipe.utils import get_table_dataset

stars = get_table_dataset("yale_bright_star_catalog5", role="bright star catalog")

stars["ra_dec"] = SkyCoord(
ra=Angle(stars["RAJ2000"], unit=u.deg),
dec=Angle(stars["DEJ2000"], unit=u.deg),
pm_ra_cosdec=stars["pmRA"].quantity, # yes, pmRA is already pm_ra_cosdec
pm_dec=stars["pmDE"].quantity,
frame="icrs",
obstime=Time("J2000.0"),
)
stars.remove_columns(["RAJ2000", "DEJ2000"])

stars = select_stars(
stars, pointing=pointing, radius=radius, Vmag_cut=magnitude_cut
)

return stars


def get_hipparcos_stars(pointing=None, radius=None, magnitude_cut=None):
"""
Get an astropy table of bright stars from the Hippoarcos star catalog
Parameters
----------
pointing: astropy Skycoord
pointing direction in the sky (if none is given, full sky is returned)
radius: astropy angular units
Radius of the sky region around pointing position. Default: full sky
Vmag_cut: float
Return only stars above a given apparent magnitude. Default: None (all entries)
Bmag_cut: float
magnitude_cut: float
Return only stars above a given absolute magnitude. Default: None (all entries)
catalog: string
Name of the catalog to be used. Usable names are found in the Enum StarCatalogues. Default: Yale
Returns
-------
Astropy table:
List of all stars after cuts with names, catalog numbers, magnitudes,
and coordinates
List of all stars after cuts with catalog numbers, magnitudes,
and coordinates as SkyCoord objects including proper motion
"""

stars = None

if CACHE_FILE.exists():
log.info("Loading stars from cached table")
try:
stars = Table.read(CACHE_FILE)
if stars.meta["Catalog"] == StarCatalogues[catalog].value:
if Bmag_cut is not None and "Bmag_cut" in stars.meta.keys():
if stars.meta["Bmag_cut"] >= Bmag_cut:
log.debug(f"Loaded table is valid for { Bmag_cut= }")
else:
log.debug(
"Loaded cache table has smaller magnitude_cut, reloading"
)
stars = None
if Vmag_cut is not None and "Vmag_cut" in stars.meta.keys():
if stars.meta["Vmag_cut"] >= Vmag_cut:
log.debug(f"Loaded table is valid for {Vmag_cut= }")
else:
log.debug(
"Loaded cache table has smaller magnitude_cut, reloading"
)
stars = None
else:
stars = None
except Exception:
log.exception("Cache file exists but reading failed. Recreating")

if stars is None:
from astroquery.vizier import Vizier

log.info("Querying Vizier for bright stars catalog")
# query vizier for stars with 0 <= Vmag <= max_magnitude
vizier = Vizier(
catalog=StarCatalogues[catalog].value,
columns=[
"recno",
"RAJ2000",
"DEJ2000",
"RAICRS",
"DEICRS",
"pmRA",
"pmDE",
"Vmag",
"Bmag",
"BTmag",
"VTmag",
],
row_limit=1000000,
)

stars = vizier.query_constraints(Vmag="0.0..10.0")[0]
if "Bmag" in stars.keys():
if Bmag_cut is not None:
stars.meta["Bmag_cut"] = Bmag_cut
elif Bmag_cut is not None:
log.exception("The chosen catalog does not have Bmag data")
if "Vmag" in stars.keys():
if Vmag_cut is not None:
stars.meta["Vmag_cut"] = Vmag_cut
elif Vmag_cut is not None:
log.exception("The chosen catalog does not have Vmag data")
stars.meta["Catalog"] = StarCatalogues[catalog].value

stars.write(CACHE_FILE, overwrite=True)

if "RAJ2000" in stars.keys():
stars["ra_dec"] = SkyCoord(
ra=Angle(stars["RAJ2000"], unit=u.deg),
dec=Angle(stars["DEJ2000"], unit=u.deg),
pm_ra_cosdec=stars["pmRA"].quantity, # yes, pmRA is already pm_ra_cosdec
pm_dec=stars["pmDE"].quantity,
frame="icrs",
obstime=Time("J2000.0"),
)
elif "RAICRS" in stars.keys():
stars["ra_dec"] = SkyCoord(
ra=Angle(stars["RAICRS"], unit=u.deg),
dec=Angle(stars["DEICRS"], unit=u.deg),
pm_ra_cosdec=stars["pmRA"].quantity, # yes, pmRA is already pm_ra_cosdec
pm_dec=stars["pmDE"].quantity,
frame="icrs",
obstime=Time("J1991.25"),
)
from ctapipe.utils import get_table_dataset

stars = get_table_dataset("hippoarcos_star_catalog5", role="bright star catalog")

stars["ra_dec"] = SkyCoord(
ra=Angle(stars["RAICRS"], unit=u.deg),
dec=Angle(stars["DEICRS"], unit=u.deg),
pm_ra_cosdec=stars["pmRA"].quantity, # yes, pmRA is already pm_ra_cosdec
pm_dec=stars["pmDE"].quantity,
frame="icrs",
obstime=Time("J1991.25"),
)
stars.remove_columns(["RAICRS", "DEICRS"])

stars = select_stars(
stars, pointing=pointing, radius=radius, Bmag_cut=Bmag_cut, Vmag_cut=Vmag_cut
stars, pointing=pointing, radius=radius, Vmag_cut=magnitude_cut
)

return stars

0 comments on commit 89f4e37

Please sign in to comment.