mga revised this gist 1 year ago. Go to revision
5 files changed, 184 insertions
README.md(file created)
| @@ -0,0 +1,3 @@ | |||
| 1 | + | ## Dataset | |
| 2 | + | ||
| 3 | + | https://www.kaggle.com/datasets/konivat/hipparcos-star-catalog | |
distribution_plot.py(file created)
| @@ -0,0 +1,50 @@ | |||
| 1 | + | import os | |
| 2 | + | import pandas as pd | |
| 3 | + | import numpy as np | |
| 4 | + | import matplotlib.pyplot as plt | |
| 5 | + | ||
| 6 | + | # Open the file "./hipparcos-voidmain.csv" for reading | |
| 7 | + | f = open("./hipparcos-voidmain.csv", "r") | |
| 8 | + | ||
| 9 | + | dataset = pd.read_csv(f, sep=",", header=0) | |
| 10 | + | ||
| 11 | + | # Get the data we want to plot | |
| 12 | + | # Apparent magnitude in the V band (visible light) | |
| 13 | + | apparentMag = np.array(dataset["Vmag"], dtype=float) | |
| 14 | + | # B-V color index | |
| 15 | + | bmv = np.array(dataset["B-V"], dtype=float) | |
| 16 | + | # Paralax | |
| 17 | + | paralax = np.array(dataset["Plx"].notna().notnull(), dtype=float) | |
| 18 | + | paralax = paralax[paralax > 0.1] | |
| 19 | + | # Compute distance only if parallax is not zero | |
| 20 | + | distance = 1000 / paralax # unit: parsec | |
| 21 | + | # Absolute magnitude | |
| 22 | + | absoluteMag = apparentMag - 5 * np.log10(distance / 10) | |
| 23 | + | ||
| 24 | + | # Plot the KDE of apparent magnitude, B-V color index, distance and absolute magnitude | |
| 25 | + | figure, axis = plt.subplots(2, 2, figsize=(10, 10)) | |
| 26 | + | axis[0, 0].set_title("KDE of apparent magnitude") | |
| 27 | + | axis[0, 0].set_xlabel("Apparent magnitude") | |
| 28 | + | axis[0, 0].set_ylabel("Density") | |
| 29 | + | ||
| 30 | + | axis[0, 1].set_title("KDE of B-V color index") | |
| 31 | + | axis[0, 1].set_xlabel("B-V color index") | |
| 32 | + | axis[0, 1].set_ylabel("Density") | |
| 33 | + | ||
| 34 | + | axis[1, 0].set_title("KDE of distance") | |
| 35 | + | axis[1, 0].set_xlabel("Distance") | |
| 36 | + | axis[1, 0].set_ylabel("Density") | |
| 37 | + | ||
| 38 | + | axis[1, 1].set_title("KDE of absolute magnitude") | |
| 39 | + | axis[1, 1].set_xlabel("Absolute magnitude") | |
| 40 | + | axis[1, 1].set_ylabel("Density") | |
| 41 | + | ||
| 42 | + | axis[0, 0].hist(apparentMag, bins=100, density=True) | |
| 43 | + | axis[0, 1].hist(bmv, bins=100, density=True) | |
| 44 | + | axis[1, 0].hist(paralax, bins=100, density=True) | |
| 45 | + | axis[1, 1].hist(absoluteMag, bins=100, density=True) | |
| 46 | + | ||
| 47 | + | plt.show() | |
| 48 | + | ||
| 49 | + | # Close the file | |
| 50 | + | f.close() | |
hr_plot.py(file created)
| @@ -0,0 +1,47 @@ | |||
| 1 | + | import os | |
| 2 | + | import pandas as pd | |
| 3 | + | import numpy as np | |
| 4 | + | import matplotlib.pyplot as plt | |
| 5 | + | ||
| 6 | + | import utilities | |
| 7 | + | ||
| 8 | + | # Open the file "/hipparcos-voidmain.csv" for reading | |
| 9 | + | f = open("./hipparcos-voidmain.csv", "r") | |
| 10 | + | ||
| 11 | + | dataset = pd.read_csv(f, sep=",", header=0) | |
| 12 | + | ||
| 13 | + | # Isolate the data we want to plot | |
| 14 | + | inputData = np.array(dataset[["B-V", "e_B-V", "Vmag", "Plx", "e_Plx"]], dtype=float) | |
| 15 | + | ||
| 16 | + | # Filter the data | |
| 17 | + | # Get rid of the NaN values | |
| 18 | + | inputData = inputData[~np.isnan(inputData).any(axis=1)] | |
| 19 | + | # Apply the following conditions: | |
| 20 | + | # - e_B-V should be smaller than 0.1 | |
| 21 | + | # - Plx should be larger than 0.1 | |
| 22 | + | # - e_Plx should be smaller than 5 | |
| 23 | + | inputData = inputData[(inputData[:, 1] < 0.1) & (inputData[:, 3] > 0.1) & (inputData[:, 4] < 5)] | |
| 24 | + | ||
| 25 | + | # Compute the distance | |
| 26 | + | distance = 1000 / inputData[:, 3] # unit: parsec | |
| 27 | + | # Compute the absolute magnitude | |
| 28 | + | absoluteMag = inputData[:, 2] - 5 * np.log10(distance / 10) | |
| 29 | + | bmv = inputData[:, 0] | |
| 30 | + | ||
| 31 | + | # Map B-V to RGB | |
| 32 | + | rgb = np.array([utilities.bv2rgb(item) for item in bmv]) | |
| 33 | + | ||
| 34 | + | # Plot the Hertzsprung-Russell diagram | |
| 35 | + | plt.figure(figsize=(10, 10)) | |
| 36 | + | plt.scatter(bmv, absoluteMag, s=0.05, c=rgb, alpha=0.75) | |
| 37 | + | plt.xlim(-0.5, 2.5) | |
| 38 | + | plt.ylim(-15, 15) | |
| 39 | + | plt.xlabel("B-V") | |
| 40 | + | plt.ylabel("Vmag") | |
| 41 | + | plt.title("Hertzsprung-Russell diagram") | |
| 42 | + | plt.gca().set_facecolor((0, 0, 0)) | |
| 43 | + | plt.gca().invert_yaxis() | |
| 44 | + | plt.show() | |
| 45 | + | ||
| 46 | + | # Close the file | |
| 47 | + | f.close() | |
position_plot.py(file created)
| @@ -0,0 +1,37 @@ | |||
| 1 | + | import os | |
| 2 | + | import pandas as pd | |
| 3 | + | import numpy as np | |
| 4 | + | import matplotlib.pyplot as plt | |
| 5 | + | ||
| 6 | + | import utilities | |
| 7 | + | ||
| 8 | + | # Open the file "./hipparcos-voidmain.csv" for reading | |
| 9 | + | f = open("./hipparcos-voidmain.csv", "r") | |
| 10 | + | ||
| 11 | + | dataset = pd.read_csv(f, sep=",", header=0) | |
| 12 | + | ||
| 13 | + | # Get the position data | |
| 14 | + | x = np.array(dataset["RAdeg"], dtype=float) | |
| 15 | + | y = np.array(dataset["DEdeg"], dtype=float) | |
| 16 | + | ||
| 17 | + | # Get the magnitude data (will be use to determine the size of the point) | |
| 18 | + | apparentMag = np.array(dataset["Vmag"], dtype=float) | |
| 19 | + | # Convert apparentMag to size | |
| 20 | + | size = 10 ** (apparentMag / -2.5) | |
| 21 | + | ||
| 22 | + | # Get the color data (will be use to determine the color of the point) | |
| 23 | + | bmv = np.array(dataset["B-V"], dtype=float) | |
| 24 | + | # Map B-V to RGB | |
| 25 | + | rgb = np.array([utilities.bv2rgb(bv) for bv in bmv]) | |
| 26 | + | ||
| 27 | + | # Plot the stars position distribution | |
| 28 | + | plt.scatter(x, y, s=size, c=rgb, alpha=0.75) | |
| 29 | + | plt.title("Stars position") | |
| 30 | + | plt.xlabel("RAdeg") | |
| 31 | + | plt.ylabel("DEdeg") | |
| 32 | + | plt.gca().set_facecolor('black') | |
| 33 | + | plt.gca().invert_xaxis() | |
| 34 | + | plt.show() | |
| 35 | + | ||
| 36 | + | # Close the file | |
| 37 | + | f.close() | |
utilities.py(file created)
| @@ -0,0 +1,47 @@ | |||
| 1 | + | # Get the rgb color from B-V color | |
| 2 | + | # More info: https://wikipedia.org/wiki/Color_index | |
| 3 | + | # Function from https://stackoverflow.com/questions/21977786/star-b-v-color-index-to-apparent-rgb-color | |
| 4 | + | # By AymericG | |
| 5 | + | def bv2rgb(bv): | |
| 6 | + | if bv < -0.40: bv = -0.40 | |
| 7 | + | if bv > 2.00: bv = 2.00 | |
| 8 | + | r = 0.0 | |
| 9 | + | g = 0.0 | |
| 10 | + | b = 0.0 | |
| 11 | + | if -0.40 <= bv<0.00: | |
| 12 | + | t=(bv+0.40)/(0.00+0.40) | |
| 13 | + | r=0.61+(0.11*t)+(0.1*t*t) | |
| 14 | + | elif 0.00 <= bv<0.40: | |
| 15 | + | t=(bv-0.00)/(0.40-0.00) | |
| 16 | + | r=0.83+(0.17*t) | |
| 17 | + | elif 0.40 <= bv<2.10: | |
| 18 | + | t=(bv-0.40)/(2.10-0.40) | |
| 19 | + | r=1.00 | |
| 20 | + | if -0.40 <= bv<0.00: | |
| 21 | + | t=(bv+0.40)/(0.00+0.40) | |
| 22 | + | g=0.70+(0.07*t)+(0.1*t*t) | |
| 23 | + | elif 0.00 <= bv<0.40: | |
| 24 | + | t=(bv-0.00)/(0.40-0.00) | |
| 25 | + | g=0.87+(0.11*t) | |
| 26 | + | elif 0.40 <= bv<1.60: | |
| 27 | + | t=(bv-0.40)/(1.60-0.40) | |
| 28 | + | g=0.98-(0.16*t) | |
| 29 | + | elif 1.60 <= bv<2.00: | |
| 30 | + | t=(bv-1.60)/(2.00-1.60) | |
| 31 | + | g=0.82-(0.5*t*t) | |
| 32 | + | if -0.40 <= bv<0.40: | |
| 33 | + | t=(bv+0.40)/(0.40+0.40) | |
| 34 | + | b=1.00 | |
| 35 | + | elif 0.40 <= bv<1.50: | |
| 36 | + | t=(bv-0.40)/(1.50-0.40) | |
| 37 | + | b=1.00-(0.47*t)+(0.1*t*t) | |
| 38 | + | elif 1.50 <= bv<1.94: | |
| 39 | + | t=(bv-1.50)/(1.94-1.50) | |
| 40 | + | b=0.63-(0.6*t*t) | |
| 41 | + | return (r, g, b) | |
| 42 | + | ||
| 43 | + | # Get the temperature from B-V color | |
| 44 | + | # More info: https://wikipedia.org/wiki/Color_index | |
| 45 | + | def bv2temp(bv): | |
| 46 | + | temp=4600.0*(1.0/(0.92*bv+1.7)+1.0/(0.92*bv+0.62)) | |
| 47 | + | return temp | |