Análisis de anisotropías de los rayos cósmicos de ultra alta energía

Introducción

Los rayos cósmicos de ultra alta energía (UHECR, con energías por encima de 10¹⁸ eV) son los proyectiles más energéticos que se observan en la naturaleza, pero su origen sigue siendo un misterio. La dificultad es que, al tratarse de partículas cargadas, los campos magnéticos de nuestra galaxia y del medio intergaláctico desvían sus trayectorias y borran la línea recta que apuntaría directamente a la fuente.

Aab et al. (2017) plantean una pregunta concreta: pese a esas desviaciones, ¿la distribución de las direcciones de llegada de los rayos cósmicos por encima de 8 EeV (1 EeV ≡ 10¹⁸ eV) es isótropa o presenta una anisotropía a gran escala? La respuesta del artículo es contundente: existe una anisotropía de tipo dipolar (un exceso suave de eventos hacia un lado del cielo y un defecto hacia el opuesto) con una significancia estadística superior a 5,2σ. Su dirección apunta a unos 125° del centro galáctico, lo que descarta un origen puramente galáctico y señala un origen extragaláctico.

Aquí se reproduce el análisis sobre un subconjunto de los datos del Observatorio Pierre Auger empleados en el artículo, proporcionado para esta actividad y disponible públicamente a través de Pierre Auger Open Data.

Energías: fluorescencia (FD) frente a superficie (SD)

El Observatorio Pierre Auger (Malargüe, Argentina, latitud −35,2°) es un sistema híbrido con dos técnicas de detección complementarias (The Pierre Auger Collaboration, 2015):

Detectores de superficie (SD): 1600 tanques de agua repartidos en 3000 km² que registran, mediante luz Cherenkov, las partículas de la cascada cuando llegan al suelo. Funcionan casi el 100 % del tiempo, de día y de noche.
Telescopios de fluorescencia (FD): observan la tenue luz ultravioleta que emiten las cascadas al excitar el nitrógeno del aire. Esta medida es casi calorimétrica, porque mide de forma casi directa la energía depositada, con poca dependencia de modelos, pero solo puede realizarse en noches despejadas y sin luna, alrededor del 15 % del tiempo.

Estas diferencias explican lo que se observa en la Figura 1:

El SD registra muchísimos más eventos que el FD. En los datos hay 2465 sucesos con E_SD > 8 EeV frente a solo 82 con medida FD útil, por el distinto ciclo de trabajo (100 % frente a 15 %).
El espectro del SD se extiende hasta aproximadamente 150 EeV y cae de forma muy abrupta (cuanto más energéticos, más raros son los eventos), mientras que las pocas energías FD se concentran en la zona de menor energía.
FD y SD son dos estimadores independientes de la misma cascada. La energía del SD (deducida de la señal a 1000 m del eje) se calibra con eventos híbridos vistos por el FD. Por eso, para un mismo suceso, ambas energías difieren dentro de la incertidumbre del 16 %, y algunos eventos con E_SD > 8 EeV tienen E_FD por debajo de 8 EeV.

Histograma de energías con muchos más eventos del detector de superficie que de fluorescencia — Figura 1. Distribución de energía de los eventos con E_SD > 8 EeV (escala logarítmica). En azul, la energía medida por los detectores de superficie (SD); en naranja, la de los telescopios de fluorescencia (FD).

Distribuciones angulares: ángulo cenital y azimutal

La Figura 2 muestra las dos distribuciones angulares, que tienen formas muy distintas: la del ángulo cenital θ y la del ángulo azimutal φ.

Histogramas del ángulo cenital, con pico hacia 45 grados, y del azimutal, casi plano — Figura 2. Distribuciones angulares de los eventos con E > 8 EeV.

El ángulo cenital θ (medido desde la vertical) crece desde casi cero hasta un máximo en torno a 40–45° y luego decae, con un corte limpio en 60°. La forma es puramente geométrica: para un flujo isótropo, el número de eventos por intervalo de ángulo cenital es proporcional a sen θ · cos θ. El factor sen θ procede del ángulo sólido (hay “más cielo” disponible a ángulos grandes); el factor cos θ, del área proyectada del detector (un detector plano “ve” menos superficie efectiva cuanto más inclinada llega la cascada). El producto sen θ · cos θ = ½ sen(2θ) alcanza su máximo en θ = 45°, justo donde la distribución tiene su pico; el corte en 60° elimina los eventos muy inclinados, que se reconstruyen con otro método.

El ángulo azimutal φ (la orientación en el plano horizontal, como la aguja de una brújula) es, en cambio, aproximadamente plano. La razón es la simetría: el detector no tiene ninguna dirección horizontal privilegiada, de modo que todas las orientaciones son equivalentes y reciben, en promedio, el mismo número de eventos. Las pequeñas fluctuaciones son estadísticas, con ligeras desviaciones debidas al campo geomagnético y a una leve inclinación del terreno. En resumen, θ tiene estructura porque combina la geometría del cielo y la del detector, mientras que φ es plano por simetría.

Mapas de direcciones de llegada y de exposición

Del subconjunto proporcionado (24.319 eventos, 22.751 con energía reconstruida por el detector de superficie), solo 2465 superan los 8 EeV: son los que se emplean en el análisis de direcciones de llegada.

Mapa del cielo en coordenadas galácticas con los eventos concentrados en el hemisferio sur — Figura 3. Mapa de direcciones de llegada (coordenadas galácticas) de los eventos con E > 8 EeV. Cada píxel indica el número de eventos registrados; el cielo austral aparece más poblado, siguiendo el patrón de exposición.

La Figura 3 representa el mapa de direcciones de llegada (número de eventos por píxel, en coordenadas galácticas). Los mapas se construyen dividiendo el cielo en píxeles de igual área mediante el esquema HEALPix (Górski et al., 2005), implementado en la librería healpy.

La Figura 4 muestra el mapa de exposición (cuánto tiempo “mira” el detector cada región del cielo), en ambos sistemas de coordenadas. Se calcula con la expresión analítica de la exposición de un detector de superficie en función de la declinación (Sommers, 2001). Como el observatorio está en el hemisferio sur, ve mejor el cielo austral y no alcanza las declinaciones muy boreales; en coordenadas ecuatoriales esto aparece como franjas horizontales. El punto clave es que la exposición es prácticamente uniforme en ascensión recta: al girar la Tierra sobre su eje, el detector barre todas las ascensiones rectas por igual.

Mapas de exposición en coordenadas ecuatoriales y galácticas, más intensos hacia el sur — Figura 4. Mapa de exposición en coordenadas ecuatoriales (izquierda) y galácticas (derecha). La exposición depende de la declinación (franjas horizontales en ecuatoriales) y es uniforme en ascensión recta.

¿Es lo que se espera? Sí. El mapa de llegadas reproduce el patrón de la exposición: hay más eventos donde el detector observa más tiempo (el sur), y no porque esa zona del cielo sea intrínsecamente “más rica”. Esta es la diferencia esencial entre ambos mapas y su distinto significado: el mapa de exposición es la referencia de lo que se vería si el cielo fuera perfectamente isótropo, mientras que el mapa de llegadas es lo que realmente se observa. A simple vista, el mapa de llegadas está dominado por la exposición, y cualquier anisotropía física queda enmascarada por la forma en que se observa; por eso, para extraer una señal real hay que comparar ambos mapas o, mejor, analizar la ascensión recta, donde la exposición es plana. El cuaderno incluye además el mapa de llegadas suavizado con una función top-hat de 45° de radio (como en el artículo), que realza las estructuras a gran escala.

En el Anexo A se incluyen las versiones suavizadas del mapa de exposición y del mapa de recuento. Su comparación es ilustrativa: el mapa de recuento suavizado reproduce casi exactamente el patrón del mapa de exposición suavizado, lo que confirma que, a simple vista, las direcciones de llegada están dominadas por la exposición y no por una anisotropía intrínseca. Esto justifica recurrir al análisis de Rayleigh en ascensión recta para extraer la señal real.

El análisis de Rayleigh

Para buscar una anisotropía a gran escala se realiza un análisis armónico (de Fourier) del número de eventos en función de la ascensión recta α. Se usa esta coordenada precisamente porque, como se ha visto, la exposición es casi constante en ella: cualquier modulación en α es señal de una anisotropía real, en concreto de la componente del dipolo perpendicular al eje de rotación de la Tierra. El análisis de Rayleigh consiste en ajustar el primer armónico de esa distribución mediante dos coeficientes de Fourier:

a = (2/N) Σ_i=1^N cos α_i b = (2/N) Σ_i=1^N sen α_i

A partir de ellos se obtienen la amplitud r = √(a² + b²) y la fase φ = arctan(b/a), que indica la ascensión recta en la que la tasa de eventos es máxima.

Lo más valioso del método es que cuantifica la confianza del resultado: proporciona la probabilidad de que, partiendo de un cielo isótropo, el azar produjera una amplitud igual o mayor que la observada:

P = exp(−N r² / 4)

Esa probabilidad se traduce después en una significancia estadística (en número de σ). Dicho de forma sencilla, el análisis de Rayleigh responde a tres preguntas a la vez: ¿hay una dirección preferente en ascensión recta?, ¿cuál es? y ¿qué confianza hay en que no sea una casualidad? Conviene notar que la amplitud que cabría esperar “por azar” disminuye al aumentar el número de eventos N (el ruido estadístico va como 1/√N), de modo que, con más estadística, una misma amplitud real resulta cada vez más significativa.

Resultados: máximo de la modulación y significancia

Aplicando el análisis de Rayleigh a los 2465 eventos con E > 8 EeV del subconjunto (Figura 5) se obtiene una amplitud del primer armónico r = 5,8 % ± 2,8 %, un máximo de la modulación en α₀ = 144° ± 28° de ascensión recta y una probabilidad de isotropía P = 0,12, equivalente a una significancia de 1,5σ. Es decir, la tasa de eventos es máxima hacia α ≈ 144°, con una amplitud del orden del 6 %.

Tasa normalizada de eventos frente a la ascensión recta con una sinusoide ajustada y una línea plana — Figura 5. Tasa de eventos normalizada en función de la ascensión recta para E > 8 EeV. Los puntos son los datos con su error 1σ; la curva azul es el primer armónico de Rayleigh (máximo en α₀ ≈ 144°) y la línea naranja discontinua, el ajuste constante (isotropía).

Ahora bien, una significancia de 1,5σ está muy por debajo del umbral habitual de detección (se suele exigir ≥ 3σ para hablar de “evidencia” y 5σ para un “descubrimiento”). Con este subconjunto de datos, el resultado es solo un indicio y sigue siendo compatible con la hipótesis de isotropía. No puede afirmarse que se haya detectado la anisotropía; únicamente, que los datos insinúan una modulación en la dirección esperada.

Comparación con la bibliografía

El artículo de Aab et al. (2017), con 32.187 eventos por encima de 8 EeV, obtiene una amplitud del primer armónico r = 4,7 % (+0,8/−0,7), un máximo en α = 100° ± 10° y una significancia de 5,2σ (un dipolo total del 6,5 % que apunta a un origen extragaláctico). La comparación con el resultado de este trabajo es muy informativa:

Amplitud: el 5,8 % ± 2,8 % obtenido es compatible con el 4,7 % del artículo dentro de las barras de error.
Fase: los 144° ± 28° distan unos 44° de los 100° ± 10° del artículo, una diferencia de aproximadamente 1,5σ: marginalmente compatibles, dada la gran incertidumbre del subconjunto.
Significancia: aquí está la diferencia más llamativa, 1,5σ frente a 5,2σ. Lejos de ser una contradicción, es justo lo que cabe esperar, porque la significancia de Rayleigh crece como √N. Con unas 13 veces menos eventos (2465 frente a 32.187), escalando la significancia obtenida, 1,5σ × √(32.187 / 2465) ≈ 5,4σ, se recupera aproximadamente la significancia que reporta el artículo. En otras palabras, con toda la estadística cabría esperar que el mismo análisis alcanzara el nivel de un descubrimiento.

Por tanto, los resultados son coincidentes, no contradictorios, con los del artículo: se reproduce la misma metodología y se obtienen una amplitud y una dirección compatibles dentro de las incertidumbres. Solo falta estadística: un subconjunto de en torno al 8 % de los datos no basta para alcanzar el nivel de significancia de un descubrimiento, pero sí apunta en la misma dirección, la de una anisotropía dipolar de probable origen extragaláctico. Con más años de datos, el Observatorio Pierre Auger ha reforzado este resultado: con 19 años de observaciones, la significancia del dipolo por encima de 8 EeV alcanza 6,8σ (The Pierre Auger Collaboration, 2025).

Conclusiones

Se ha reproducido el análisis de Aab et al. (2017) sobre un subconjunto de datos del Observatorio Pierre Auger. Se ha explicado por qué difieren las distribuciones de energía SD y FD (distinto ciclo de trabajo y calibración cruzada de dos estimadores independientes), por qué las distribuciones angulares tienen formas distintas (geometría sen θ · cos θ frente a simetría azimutal) y por qué los mapas de llegada están dominados por la exposición. El análisis de Rayleigh en ascensión recta da un máximo en α₀ ≈ 144° con una amplitud del orden del 6 % y una significancia de 1,5σ: un indicio que, dentro de la limitación estadística del subconjunto, es plenamente compatible con la anisotropía dipolar de origen extragaláctico de 5,2σ que reporta el artículo.

Anexo A. Material complementario del análisis de mapas

Para completar el análisis de los mapas se muestran sus versiones suavizadas con una función top-hat de 45° de radio, en coordenadas ecuatoriales y galácticas.

Mapa de exposición suavizado en coordenadas ecuatoriales y galácticas — Figura A1. Mapa de exposición suavizado con una función top-hat de 45° de radio, en coordenadas ecuatoriales (izquierda) y galácticas (derecha). Comparado con la Figura 4, el patrón es prácticamente idéntico: como la exposición es una función suave de la declinación, el suavizado solo reduce ligeramente el valor máximo por píxel.

Mapa de recuento suavizado casi idéntico al mapa de exposición — Figura A2. Mapa de recuento (direcciones de llegada) suavizado con la misma función top-hat de 45°, para E > 8 EeV, en coordenadas ecuatoriales (izquierda) y galácticas (derecha). Su gran parecido con la Figura A1 evidencia que las llegadas están dominadas por la exposición.

Referencias

Aab, A., et al. (2017). Observation of a large-scale anisotropy in the arrival directions of cosmic rays above 8 × 10¹⁸ eV. Science, 357(6357), 1266–1270.
Górski, K. M., Hivon, E., Banday, A. J., Wandelt, B. D., Hansen, F. K., Reinecke, M., & Bartelmann, M. (2005). HEALPix: A framework for high-resolution discretization and fast analysis of data distributed on the sphere. The Astrophysical Journal, 622(2), 759–771.
Sommers, P. (2001). Cosmic ray anisotropy analysis with a full-sky observatory. Astroparticle Physics, 14(4), 271–286.
The Pierre Auger Collaboration (2015). The Pierre Auger Cosmic Ray Observatory. Nuclear Instruments and Methods in Physics Research Section A, 798, 172–213.
The Pierre Auger Collaboration (2025). Large-scale anisotropies of ultra-high-energy cosmic rays measured at the Pierre Auger Observatory. arXiv:2507.19243.