Triangulating stars on the night sky
In a previous notebook, I have shown properties of the distribution of stars in the sky. Here, I would like to use the existing database of stars' positions and display them as a triangulation
Let's first initialize the notebook:
import numpy as np
np.set_printoptions(precision=6, suppress=True)
%matplotlib inline
%config InlineBackend.figure_format = 'svg'
import matplotlib.pyplot as plt
phi = (np.sqrt(5)+1)/2
fig_width = 10
figsize = (fig_width, fig_width/phi)
importing data¶
Importing all stars' position is as simple as invoking the HYG database:
import pandas as pd
url = "https://github.com/astronexus/HYG-Database/raw/master/hygdata_v3.csv"
space = pd.read_csv(url, index_col=0)
print(f'Columns in the database: {space.columns=}')
print(f'Number of stars in the catalog = {len(space)=}')
extracting ra, dec and mag¶
For which we may extract what interests us: position (right ascension and declination) and visual magnitude:
space_pos = space[['ra', 'dec', 'mag']]
space_pos
First, right ascension is "the angular distance of a particular point measured eastward along the celestial equator from the Sun at the March equinox to the (hour circle of the) point in question above the earth." It is given in hours:
ra_min, ra_max = 0, 24
print(f"RA: {space_pos['ra'].min()=}, {space_pos['ra'].max()=}")
normalizing ra into az¶
Let's convert this in visual angle, that is in an azimuth:
space_norm = space_pos[['dec', 'mag']].copy()
space_norm
az_min, az_max = 0, 360
def ra2az(ra):
return az_max - ra / ra_max * az_max
space_norm["az"] = ra2az(space_pos['ra'])
print(f"AZ: {space_norm['az'].min()=}, {space_norm['az'].max()=}")
Then, declination is "comparable to geographic latitude, projected onto the celestial sphere" and is given here in degrees:
dec_min, dec_max = -90, 90
print(f"DEC: {space_norm['dec'].min()=}, {space_norm['dec'].max()=}")
The magnitude varies a lot (the less the value, the more it is visible - "The apparent magnitudes of known objects range from the Sun at −26.7 to objects in deep Hubble Space Telescope images of magnitude +31.5"):
print(f"mag: {space_norm['mag'].min()=}, {space_norm['mag'].max()=}")
Let's normalize the lower bound:
space_norm['mag'] = space_pos['mag'] - space_pos['mag'].min()
print(f"mag: {space_norm['mag'].min()=}, {space_norm['mag'].max()=}")
Let's define a threshold using the quantile function:
print(f"{space_norm['mag'].quantile(q=0.01)=}")
Which leaves only a limited number of stars
space_bright = space_norm[space_norm['mag']<space_norm['mag'].quantile(q=0.05)]
#space_bright = space_pos[space_pos['mag']<6]
print(f'Number of bright stars = {len(space_bright)=}')
print(f"MAG: {space_bright['mag'].min()=}, {space_bright['mag'].max()=}")
scatter plots¶
From these elements, we may plot the stars on these coordinates:
stars_color = np.array([255., 235., 220.])
stars_color /= 255.
fig, ax = plt.subplots(figsize=(fig_width, fig_width/phi))
ax.set_facecolor('black')
ax.scatter(space_bright['az'], space_bright['dec'], c=stars_color[None, :], ec='none', s=.1 * (space_bright['mag'].max()-space_bright['mag']))
ax.set_xlabel('Azimuth')
ax.set_ylabel('Declination')
ax.set_xlim(az_min, az_max)
ax.set_ylim(dec_min, dec_max);