Soy Alfonso Navas (@navasa) y en este cuaderno presentaré cómo generé un diagrama de Hertzprung-Russell, basándome en el diseño que se encuentra en esta referencia. Luego, mostraré como hice una animación de la gráfica creada.

Ejercicio No. 1

Introducción

Un digrama de diagrama de Hertzprung-Russel, también llamado simplemente como diagrama H-R, es un gráfico de dispersión que indica la relación entre la temperatura y la luminosidad de las estrellas. En estos diagramas se pueden distinguir 5 grupos de estrellas. Por una lado, están las estrellas de la secuacia principal. A este primer grupo corresponden la gran mayoría de la estrellas, cuya luminosidad es directamente proporcional a su temperatura. Por otra parte, las gigantes y supergigantes rojas son estrellas de gran tamaño pero con una masa pequeña o mediana. Estas estrellas son muy luminosas a pesar de no ser muy calientes. En contraposición, tenemos a las enanas blancas, las cuales son estrellas muy calientes y de pequeño tamaño, lo cual implica que éstas no sean muy luminosas a pesar de su temperatura. Por último, las gigantes azules son estrellas que sobresalen notoriamente de la secuacia principal debido a su gran tamaño y luminosidad.

Solución

En primer lugar, importamos las bibliotecas que vamos a estar utilizando en el desarrollo de este ejercicio.

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

# Para modificar el colormap
from matplotlib.colors import ListedColormap,LinearSegmentedColormap

# Para personalizar más los ejes de las gráficas (e.g. ajustar el grosor)
from pylab import *

# para "desordenar" los datos
import random

# para crear la animación necesitamos:
import glob
import os 
os.environ["IMAGEIO_FFMPEG_EXE"] = "/usr/bin/ffmpeg" #soluciona los problemas al importar moviepy
import moviepy.editor as mpy

Ahora, cargamos los datos.

dwarfs = np.loadtxt('./data/dwarfs.csv', delimiter = ',', skiprows = 1 )
giants = np.loadtxt('./data/giants.txt', delimiter = ' ', skiprows = 1 )
ms = np.loadtxt('./data/ms.csv', delimiter = ',', skiprows = 1)
supergiants = np.loadtxt('./data/supergiants.txt', delimiter = ' ', skiprows = 1)

Creamos un arreglo que contenga todos los datos y creamos le asignamos un vector a cada variable.

datos = np.concatenate((dwarfs, giants, ms, supergiants))

lum = datos[:, 0] # luminosidad
temp = datos [:, 1] # temperatura
radius = datos[:, 2] # radio

Antes de graficar, configuraremos la paleta de colores (colormap) que vamos a utilizar, para que nuestro gráfico parezca más al modelo.

RdYlBu = plt.get_cmap('RdYlBu') #guardamos el colormap de referencia (RdYlBu) en una variable
RdYlBu_m = ListedColormap(RdYlBu(np.linspace(0.15, 0.95, 256)))# selecionamos un rango solamente de este colormap

A continuación, procedemos a graficar.

# Ajustamos el tamaño de nuestro gráfica
plt.figure(figsize=(11,7))
ax = plt.gca()

# Ajustamos el grosor de los ejes
rc('axes', linewidth= 1.5)

# Graficamos
plt.scatter(temp, lum, # graficamos los valores de temperatura vs. luminosidad
            s = radius*11, # hacemos que el tamaño de cada burbuja sea proporcional al radio de la estrella 
            c = np.log2(temp), # ponemos al color de cada burbuja como función de la temperatura
            cmap = RdYlBu_m, # definimos el "conjunto de llegada" de colores 
            edgecolor = 'grey', # dibujamos bordes grises a cada burbuja
            linewidths = 1.25) # ajustamos el grosor de los bordes de cada burbuja

# Colocamos las etiquetas a los ejes
plt.xlabel('Temperature (K)', fontsize = 12, fontweight='bold')                  
plt.ylabel('Luminosity($\mathregular{L_{sun}}$)', fontsize = 12, fontweight='bold')

# Modificamos el eje x 
ax.invert_xaxis() # invertimos el eje x
plt.xscale('log', base = 2) # pasamos a escala logarítmica
plt.xticks([3500, 7000, 14000], ['3,500', '7,000', '14,000']) # colocamos los ticks
plt.xlim(14000, 3400)

# Ajustamos la escala y los límites del eje y
plt.yscale('log')
plt.ylim(0.2*10**-4, 10**7)

# Ponemos en negrita y ajustamos el tamaño de fuente de los ticks sobre los ejes
for tick in ax.xaxis.get_major_ticks():
    tick.label1.set_fontsize(9)
    tick.label1.set_fontweight('bold')
for tick in ax.yaxis.get_major_ticks():
    tick.label1.set_fontsize(9)
    tick.label1.set_fontweight('bold')
    
# Eliminamos spines superior y derecho
ax.spines['right'].set_color('none')    
ax.spines['top'].set_color('none')

#Insertamos el texto
text(5200, 0.4*10**7, 'Red Supergiants', fontsize = 14, fontweight = 'bold') # Red Supergiants
text(4750, 0.25*10**3, 'Red Giants', fontsize = 14, fontweight = 'bold') # Red Giants
text(10000, 0.3*10**7, 'Blue Giants', fontsize = 14, fontweight = 'bold') # Blue Giants
text(12000, 0.8*10**3, 'Main Sequence', fontsize = 14, fontweight = 'bold')  # Main sequence (top)
text(5500, 1.5*10**-2, 'Main Sequence', fontsize = 14, fontweight = 'bold')  # Main sequence (down)
text(7500, 0.7*10**-3, 'White Dwarfs', fontsize = 14, fontweight = 'bold') # White Dwarfs

# Mostramos la gráfica
plt.show()

Ejercicio No. 2

A continuación vamos a generar una animación a partir de la gráfica que hemos creado.

Para esto, primero creamos una función que genere, que a partir de los datos, un gráfico como el que hicimos y lo guarde con un nombre de salida que se recibe por parámetro ( fname ).

def graficar(lum, temp, radius, fname):
    
    # Ajustamos el tamaño de nuestro gráfica
    plt.figure(figsize=(11,7))
    plt.subplot(111)
    
    ax = plt.gca()
    
    # Ajustamos el grosor de los ejes
    rc('axes', linewidth= 1.5)
    
    # Graficamos
    plt.scatter(temp, lum, # graficamos los valores de temperatura vs. luminosidad
                s = radius*11, # hacemos que el tamaño de cada burbuja sea proporcional al radio de la estrella 
                c = np.log2(temp), # ponemos al color de cada burbuja como función de la temperatura
                cmap = RdYlBu_m, # definimos el "conjunto de llegada" de colores 
                edgecolor = 'grey', # dibujamos bordes grises a cada burbuja
                linewidths = 1.25) # ajustamos el grosor de los bordes de cada burbuja
    # Colocamos las etiquetas a los ejes
    plt.xlabel('Temperature (K)', fontsize = 12, fontweight='bold')                  
    plt.ylabel('Luminosity($\mathregular{L_{sun}}$)', fontsize = 12, fontweight='bold')

    # Modificamos el eje x 
    ax.invert_xaxis() # invertimos el eje x
    plt.xscale('log', base = 2) # pasamos a escala logarítmica
    plt.xticks([3500, 7000, 14000], ['3,500', '7,000', '14,000']) # colocamos los ticks
    plt.xlim(14000, 3400)
    
    # Ajustamos la escala y los límites del eje y
    plt.yscale('log')
    plt.ylim(0.2*10**-4, 10**7)

    # Ponemos en negrita y ajustamos el tamaño de fuente de los ticks sobre los ejes
    for tick in ax.xaxis.get_major_ticks():
        tick.label1.set_fontsize(9)
        tick.label1.set_fontweight('bold')
    for tick in ax.yaxis.get_major_ticks():
        tick.label1.set_fontsize(9)
        tick.label1.set_fontweight('bold')
    
    # Eliminamos spines superior y derecho
    ax.spines['right'].set_color('none')    
    ax.spines['top'].set_color('none') 
    
    # Insertamos el texto 
    text(5200, 0.4*10**7, 'Red Supergiants', fontsize = 14, fontweight = 'bold') # Red Supergiants
    text(4750, 0.25*10**3, 'Red Giants', fontsize = 14, fontweight = 'bold') # Red Giants
    text(10000, 0.3*10**7, 'Blue Giants', fontsize = 14, fontweight = 'bold') # Blue Giants
    text(12000, 0.8*10**3, 'Main Sequence', fontsize = 14, fontweight = 'bold')  # Main sequence (top)
    text(5500, 1.5*10**-2, 'Main Sequence', fontsize = 14, fontweight = 'bold')  # Main sequence (down)
    text(7500, 0.7*10**-3, 'White Dwarfs', fontsize = 14, fontweight = 'bold') # White Dwarfs
    
    #Guardamos la gráfica
    plt.savefig('./animacion/' + fname, dpi = 200, bbox_inches='tight')

Ahora, generamos las gráficas para nuestra simulación.

# Primero, desordenamos nuestros datos para la simulación
shuffle(datos)

lum = datos[:, 0] # luminosidad
temp = datos [:, 1] # temperatura
radius = datos[:, 2] # radio

# Empezamos a generar las gráficas
for i in range(len(datos)):
    nsalida = 'HR_' + str(i).zfill(5)+'.png' # generamos el nombre el archivo de la gráfica
    graficar(lum[:i + 1], temp[:i + 1], radius[:i +1], nsalida) # y hacemos la gráfica con la función definida
    plt.close()

ngif = 'HR_animacion' # nombre de la animación
fps = 24 # cuadros por segundo

#Creamos una lista (ordenada) con la dirreccion  de todos las gráficas que hemos generado
lista_graficas = sorted(glob.glob('./animacion/HR_*.png')) 

# Usando la función ImageSequenceClip de moviepy, generamos nuestra animación
animacion = mpy.ImageSequenceClip(lista_graficas[1:], fps=fps)

# Guardamos nuestra animación en forma de gif
animacion.write_gif(ngif + '.gif')

t:   0%|          | 0/105 [00:00<?, ?it/s, now=None]

MoviePy - Building file HR_animacion.gif with imageio.