Cosmología - Actividad Guiada I
Jorge Martínez Garrido
January 1, 2024
Introducción
El objetivo del presente trabajo es investigar la relación existente entre el corrimiento al rojo y la distancia de un cuerpo en el Universo. Para ello, se han calculado los datos medios fotométricos de diferentes cúmulos de galaxias. Obtenidas las magnitudes medias en las principales bandas del espectro electromagnético, se ha representado su valor frente al corrimiento al rojo para cada cúmulo de galaxias.
Todos los datos han sido procesados con Python, empleando herramientas open-source. El código se adjunta junto con las explicaciones para asegurar una reproducibilidad de los resultados. Además, reduce el número de errores humanos y permite escalar de forma rápida el número de resultados.
Recopilando datos
Se ha recopilado información de más de 150000 clústers de galaxias a través de NASA/IPAC Extragalactic Database (NED). Estos datos se han transformado a un archivo tipo Apache Parquet para mejorar su almacenamiento y lectura.
import pathlib
import pandas as pd
DATAFILE = pathlib.Path() / "dat" / "NED.parquet"
galaxy_cluster = pd.read_parquet(DATAFILE).drop("No.", axis=1)
Estos datos contienen el nombre del clúster, sus coordenadas en J2000 y su desplazamiento al rojo. La estuctura resulta la siguiente para las primeras cinco filas:
DESIRED_CLUSTERS = 25
galaxy_cluster = galaxy_cluster.head(DESIRED_CLUSTERS)
print(galaxy_cluster.head(5))
Object Name RA DEC Type Redshift
0 RM J080000.5+314126.6 120.00204 31.69072 GClstr 0.356692
1 WHL J080000.7+415337 120.00292 41.89361 GClstr 0.131000
2 WHL J080001.1+304141 120.00458 30.69472 GClstr 0.485200
3 GMBCG J120.00618+53.97479 120.00618 53.97479 GClstr 0.464000
4 GMBCG J120.00897+48.07423 120.00897 48.07423 GClstr 0.381000
El objetivo final es representar la distancia en función del corrimiento al rojo.
Medir distancias no es tarea fácil pero sí recoger la luz que llega de un punto. Además, conocida la relación:
$$ m - M = 5 \log{\frac{d_L}{10}} $$
es posible relacionar la magnitud con la distancia.
Por ello, la primera tarea es recopilar información en lo que a magnitud se refiere para las diferentes galaxias que forman cada uno de los clúster.
Identificando las galaxias que contiene un clúster
Para averiguar qué galaxias contiene cada uno de los clúster podemos asumir que las galaxias se encuentran en un radio de 1Mpc. Este valor será más grande para ciertos clúster pero es un valor bastante aproximado que permite trabajar con los datos.
Los datos se han recopilado de Sloan Digital Sky Survey (SDSS) utilizando el módulo de astroquery.sdss.
Es necesario pasar un radio de búsqueda en el método SDSS.query_region. Dicho radio se calcula mediante la expresión:
$$ \theta = \frac{d}{D_A} $$
donde $d$ es el tamaño del clúster asumido y $D_A$ es la distancia radial comóvil al clúster.
Para calcular la distancia radial $D_A$, se ha utilizado el script de Python proporcionado en CosmoCalc. El script ha sido adaptado para usar Python 3.
Los valores utilizados para las contantes son $H_0=69.6$ km/s/Mpc, $\Omega_m=0.286$ y $\Omega_{\Lambda}=0.714$:
import numpy as np
import cosmocalc
H0, Omega_m, Omega_vac = 69.6, 0.286, 0.714
indices = [
'Age at Redshift', 'Comoving Radial Distance',
'Angular Size Distance', 'Distance Modulus'
]
def calculate_cosmological_values(cluster):
result = cosmocalc.solve(cluster['Redshift'], H0, Omega_m, Omega_vac)
return pd.Series(result, index=indices)
galaxy_cluster[indices] = galaxy_cluster.apply(calculate_cosmological_values, axis=1)
GENERIC_CLUSTER_SIZE = 1
RAD_TO_ARCMIN = 180 / np.pi * 60
galaxy_cluster["Angular Radius"] = (
1 / galaxy_cluster["Comoving Radial Distance"] * RAD_TO_ARCMIN
)
Los datos calculados se muestran a continuación:
print(galaxy_cluster[['Object Name', 'Comoving Radial Distance', 'Angular Size Distance']].head(5))
Object Name Comoving Radial Distance Angular Size Distance
0 RM J080000.5+314126.6 1412.465618 5.047436
1 WHL J080000.7+415337 547.979098 2.348963
2 WHL J080001.1+304141 1859.330082 6.069393
3 GMBCG J120.00618+53.97479 1787.833163 5.920518
4 GMBCG J120.00897+48.07423 1499.486599 5.264087
Una vez calculado el radio de búsqueda, se resuelven los objetos que contengan datos de fotometría. Para ello se definen las siguientes funciones auxiliares:
import astropy.coordinates as coords
from astroquery.sdss import SDSS
def search_cluster_from_position(ra, dec, photoobj_fields=None):
cluster_coord = coords.SkyCoord(ra=ra, dec=dec)
return SDSS.query_region(
cluster_coord, spectro=True, photoobj_fields=photoobj_fields
)
def search_galaxies_within_radius_from_position(ra, dec, radius, photoobj_fields):
cluster_coord = coords.SkyCoord(ra=ra, dec=dec)
return SDSS.query_region(
cluster_coord, radius=radius, photoobj_fields=photoobj_fields
)
Se itera sobre cada uno de los clústers, guardando la información de los objetos encontrados así como el clúster al que pertenecen:
from astropy import units as u
BANDS = "UGRIZ".lower()
BAND_INDICES = [f"modelMag_{band}" for band in BANDS]
COLORS = ["royalblue", "darkgreen", "red", "darkred", "black"]
photoobj_fields = ["objid", "ra", "dec", "extinction_r"] + BAND_INDICES
threshold = 0.05
for index, cluster in galaxy_cluster.iterrows():
ra, dec, redshift, radius = (
cluster["RA"] * u.deg, cluster["DEC"] * u.deg,
cluster["Redshift"], cluster["Angular Radius"] * u.arcmin
)
galaxies_table = search_galaxies_within_radius_from_position(
ra, dec, radius, photoobj_fields
)
if not galaxies_table:
continue
galaxies = galaxies_table.to_pandas()
# Discard failing values
valid_galaxies = (galaxies["modelMag_r"] > 0) & (galaxies["extinction_r"] > 0)
galaxies = galaxies[valid_galaxies]
# Apply redshift threshold
mean_redshift = galaxies["extinction_r"].mean()
lower_limit, upper_limit = mean_redshift - threshold, mean_redshift + threshold
valid_galaxies = (
(galaxies["extinction_r"] > lower_limit) &
(galaxies["extinction_r"] < upper_limit)
)
galaxies = galaxies[valid_galaxies]
if galaxies.empty:
continue
indices = BAND_INDICES + ["extinction_r"]
def calculate_mean_values(galaxy):
mean_values = (galaxy[band].mean() for band in indices)
return pd.Series(mean_values, index=indices)
galaxy_cluster[indices] = galaxies.apply(calculate_mean_values, axis=1)
galaxy_cluster = galaxy_cluster.dropna()
print(galaxy_cluster[["Redshift"] + BAND_INDICES].head(5))
Redshift modelMag_u modelMag_g modelMag_r modelMag_i modelMag_z
0 0.356692 17.77366 16.55025 16.14032 15.98137 15.93708
1 0.131000 17.78428 16.55184 16.14135 15.98272 15.93935
2 0.485200 17.78422 16.55191 16.14106 15.98271 15.93946
3 0.464000 23.96814 21.55711 19.88091 19.27671 19.20989
4 0.381000 16.50165 15.07830 14.53191 14.31796 14.22687
Grabamos los datos en un archivo para los primeros 15 clústers:
galaxy_cluster.head(15).to_excel("plantilla.xlsx")
Representando la relación entre magnitud aparente y corrimiento al rojo
Podemos representar la evolución de la magnitud de cada banda en función del corrimiento al rojo:
import matplotlib.pyplot as plt
fig, ax = plt.subplots(nrows=1, ncols=1)
for mag_band, color in zip(BAND_INDICES, COLORS):
z, mag = (
np.asarray(galaxy_cluster["Redshift"]),
np.asarray(galaxy_cluster[mag_band])
)
ax.scatter(z, mag, c=color, alpha=0.35)
a, b = np.polyfit(np.log10(z), mag, 1)
# Calcular el coeficiente de determinación (R²)
correlation_matrix = np.corrcoef(z, mag)
r_squared = correlation_matrix[0, 1] ** 2
ax.plot(
sorted(z), sorted(a * z + b),
color=color,
label=f"Banda {mag_band[-1].upper()} con $R^2$ = {r_squared:.2f}"
)
ax.legend(shadow=True)
ax.set_title(
"Evolución de la magnitud en función del corrimiento al rojo"
)
ax.set_xlabel("Corrimiento al rojo\n log(z)")
ax.set_ylabel("Magnitud")
ax.set_xlim([0.1, 0.7])
ax.set_ylim([19.5, 23])
ax.set_xscale("log")
ax.legend(shadow=True, loc="upper left")
plt.show()
Acerca de los resultados obtenidos y su precisión
La figura anterior muestra una relación entre la magnitud y el corrimiento al rojo. Sin embargo, el valor de la correlación es muy pequeño. Esto indica que los datos deben de ser mejorados.
El autor ha corregido el valor del “threshold” en varias ocasiones. Los histogramas demuestran que los datos presentan una fuerte dispersión, por lo que tomar la media de las magnitudes no es suficiente. Sería necesario aplicar un mejor filtrado para identificar de forma robusta qué cuerpos devueltos por SDSS pertenecen al cluster observado y cuales no.
Otro factor importante es la incertidumbre en las diferentes bandas. El polvo interestelar tiende a absover las altas frecuencias, dejando pasar las frecuencias rojizas. Esto implica que los valores para la banda R e I son más precisos. Por ello, el análisis debería de centrarse en estas bandas.
Calculando el valor de $H_0$
Es posible utilizar cefeidas y supernovas tipo Ia como puntos de referencia.
Las cefeidas son estrellas variables cuyo período de variación está relacionado con su luminosidad intrínseca. Al medir el período y la magnitud aparente de las cefeidas en los cúmulos de galaxias, puedes establecer una relación entre la magnitud y la distancia basada en la relación periodo-luminosidad
Las supernovas tipo Ia también son objetos astronómicos estándar ampliamente utilizados para medir distancias en el universo. La magnitud aparente de una supernova tipo Ia está relacionada con su luminosidad intrínseca, lo que permite estimar la distancia a la supernova y, por lo tanto, al cúmulo en el que se encuentra.
Calculadas las distancias de cada cúmulo, es posible generar un gráfico de distancia frente a redshift, el valor de $H_0$ es la pendiente de la línea que mejor ajusta a los datos.