Observations | Baryon Acoustic Oscillations

Baryon Acoustic Oscillations (BAO)

Sound waves frozen into the large-scale distribution of galaxies: a cosmic standard ruler for mapping the expansion history of the universe and constraining dark energy

1. Introduction: Sound Waves in the Primordial Plasma

The Origin of Baryon Acoustic Oscillations

Before recombination at $z \approx 1100$, the universe consisted of a tightly coupled plasma of baryons (protons and electrons) and photons. Tiny primordial density perturbations — seeded during cosmic inflation — drove acoustic oscillations in this photon-baryon fluid. Each initial overdensity launched a spherical sound wave that propagated outward at the sound speed of the relativistic plasma.

At recombination, photons decoupled from baryons and the sound waves froze. The baryon perturbation was left as a spherical shell of enhanced density at the sound horizon radius $r_s \approx 147$ Mpc from each initial overdensity. This shell became imprinted in the distribution of all matter, and subsequently in the galaxy distribution observed today — a faint but measurable excess of galaxy pairs separated by approximately 150 Mpc.

The BAO signal provides one of the cleanest geometric probes in cosmology. Unlike supernovae (which require empirical standardization) or the CMB (which probes a single redshift), BAO can be measured at multiple redshifts using galaxy surveys. The characteristic scale $r_s$ is set by well-understood pre-recombination physics, making it a robust standard ruler. By measuring the apparent angular and radial extent of this ruler at different redshifts, we map out the expansion history $H(z)$ and the angular diameter distance $D_A(z)$, thereby constraining dark energy and spatial curvature.

Physical Picture

Consider an initial point-like overdensity in the early universe. In the tightly coupled photon-baryon fluid, an outgoing spherical pressure wave propagates at the sound speed. Dark matter, being collisionless, stays near the original overdensity. At decoupling:

The dark matter perturbation remains at the center
The baryon shell has traveled to radius $r_s$ (the sound horizon)
Photons free-stream away and no longer contribute to clustering
Baryons fall into the combined dark matter + baryon gravitational potential
The final matter distribution has a small excess at $r \approx r_s$

This excess is tiny — only about 1% enhancement in the correlation function — but it is statistically robust when measured over large cosmological volumes containing millions of galaxies.

2. Derivation 1: Sound Horizon at the Drag Epoch

The fundamental scale of BAO is the comoving sound horizon $r_s$ at the drag epoch — the moment when baryons decouple from the photon drag force. This occurs at $z_{\rm drag} \approx 1060$, slightly after photon decoupling because the baryon-to-photon momentum transfer rate drops below the expansion rate slightly later than the photon last scattering.

2.1 Sound Speed in the Photon-Baryon Fluid

The photon-baryon fluid has pressure provided by photons and inertia provided by both baryons and photons. The photon equation of state is $P_\gamma = \rho_\gamma c^2 / 3$. The total energy density contributing to inertia is $\rho_\gamma + \rho_b$ (in the relativistic limit where we include rest-mass energy). The adiabatic sound speed is:

$$c_s^2 = \frac{\partial P}{\partial \rho}\bigg|_s = \frac{\partial P_\gamma / \partial T}{\partial (\rho_\gamma + \rho_b) / \partial T} = \frac{\frac{4}{3}\rho_\gamma c^2 / (3T)}{\frac{4}{3}\rho_\gamma / T + \rho_b / T}$$

Since $\rho_\gamma \propto T^4$ implies $\partial \rho_\gamma / \partial T = 4\rho_\gamma/T$, and for non-relativistic baryons $\rho_b \propto T^3$ (through the adiabatic relation during tight coupling), we can simplify by defining the baryon loading parameter:

$$R \equiv \frac{3\rho_b}{4\rho_\gamma}$$

This gives the photon-baryon sound speed:

$$\boxed{c_s = \frac{c}{\sqrt{3(1 + R)}}}$$

Note that $R$ evolves with redshift. Since $\rho_b \propto a^{-3}$and $\rho_\gamma \propto a^{-4}$, we have $R \propto a \propto 1/(1+z)$:

$$R(z) = \frac{3\rho_{b,0}}{4\rho_{\gamma,0}} \cdot \frac{1}{1+z} = \frac{3\Omega_b}{4\Omega_\gamma} \cdot \frac{1}{1+z}$$

For Planck 2018 parameters ($\Omega_b h^2 = 0.02237$, $T_{\rm CMB} = 2.7255$ K), we find $R(z_{\rm drag}) \approx 0.63$, giving $c_s(z_{\rm drag}) \approx 0.454\,c$. At early times ($z \gg z_{\rm drag}$), $R \to 0$ and $c_s \to c/\sqrt{3} \approx 0.577\,c$, the pure radiation sound speed.

2.2 Computing the Sound Horizon

The comoving sound horizon is the total comoving distance a sound wave can travel from the Big Bang to the drag epoch:

$$r_s = \int_0^{t_{\rm drag}} \frac{c_s(t)}{a(t)}\,dt$$

Converting to a redshift integral using $dt = -dz / [H(z)(1+z)]$:

$$\boxed{r_s = \int_{z_{\rm drag}}^{\infty} \frac{c_s(z)}{H(z)}\,dz = \int_{z_{\rm drag}}^{\infty} \frac{c}{\sqrt{3(1 + R(z))}}\,\frac{dz}{H(z)}}$$

In the radiation-matter era (relevant for $z > 1000$), the Hubble parameter is:

$$H(z) = H_0 \sqrt{\Omega_r (1+z)^4 + \Omega_m (1+z)^3}$$

where $\Omega_r \approx 9.15 \times 10^{-5}$ includes photons and three species of massless neutrinos. For the Planck 2018 best-fit $\Lambda$CDM cosmology:

Sound Horizon: Numerical Result

Performing the integral numerically with $z_{\rm drag} = 1059.94$,$\Omega_b h^2 = 0.02237$, $\Omega_m h^2 = 0.1430$, and $H_0 = 67.36$ km/s/Mpc:

$$r_s(z_{\rm drag}) \approx 147.09 \;\text{Mpc}$$

This value is determined to sub-percent precision by Planck and is the foundation of all BAO distance measurements. The sound horizon is effectively calibrated by the CMB, making BAO an “absolute” distance indicator when combined with CMB data.

3. Derivation 2: BAO in the Galaxy Correlation Function

3.1 The Two-Point Correlation Function

The galaxy two-point correlation function $\xi(r)$ measures the excess probability of finding two galaxies separated by comoving distance $r$, relative to a random distribution:

$$\boxed{\xi(r) = \langle \delta(\mathbf{x})\,\delta(\mathbf{x} + \mathbf{r}) \rangle}$$

where $\delta(\mathbf{x}) = [\rho(\mathbf{x}) - \bar{\rho}] / \bar{\rho}$ is the density contrast. By statistical isotropy, $\xi$ depends only on $|\mathbf{r}| = r$. The BAO signal appears as a broad bump (not a sharp spike) in $\xi(r)$ at $r \approx r_s \approx 150$ Mpc.

3.2 Connection to the Power Spectrum

The correlation function and the matter power spectrum $P(k)$ form a Fourier transform pair. In three dimensions, for an isotropic field:

$$\boxed{\xi(r) = \int_0^\infty P(k)\,j_0(kr)\,\frac{k^2\,dk}{2\pi^2}}$$

where $j_0(x) = \sin(x)/x$ is the zeroth-order spherical Bessel function. The inverse relation is:

$$P(k) = 4\pi \int_0^\infty \xi(r)\,j_0(kr)\,r^2\,dr$$

3.3 Oscillations in the Power Spectrum

The BAO feature in the power spectrum appears as quasi-sinusoidal oscillations (wiggles) superimposed on the smooth broadband shape. The characteristic wavenumber of these oscillations is:

$$\boxed{k_{\rm BAO} = \frac{2\pi}{r_s} \approx \frac{2\pi}{147\;\text{Mpc}} \approx 0.043\;\text{Mpc}^{-1}}$$

The oscillatory component of the transfer function can be written schematically as:

$$T(k) = T_{\rm smooth}(k)\left[1 + A_{\rm BAO}\,\frac{\sin(k\,r_s)}{k\,r_s}\,e^{-k^2\Sigma^2/2}\right]$$

where $T_{\rm smooth}(k)$ is the smooth (no-wiggle) transfer function,$A_{\rm BAO} \approx \Omega_b / \Omega_m$ sets the oscillation amplitude (proportional to the baryon fraction), and $\Sigma$ accounts for damping due to Silk diffusion and nonlinear structure formation. The matter power spectrum is then:

$$P(k) = A_s\,k^{n_s}\,T^2(k)\,D^2(z)$$

where $D(z)$ is the linear growth factor. The sinusoidal oscillations in $P(k)$Fourier-transform into the localized bump in $\xi(r)$ at $r = r_s$. This Fourier-pair relationship — oscillations in $k$-space mapping to a peak in real space — is the mathematical heart of BAO.

Why the BAO Peak is Broad

The BAO peak in $\xi(r)$ is not a delta function because the oscillations in $P(k)$are damped at high $k$ by Silk damping (pre-recombination photon diffusion) and by nonlinear evolution (post-recombination gravitational collapse). This damping corresponds to an exponential envelope $e^{-k^2\Sigma^2/2}$ with $\Sigma \sim 8$ Mpc, which broadens the real-space peak to a width of roughly 10–20 Mpc. BAO reconstruction techniques can partially reverse the nonlinear broadening, sharpening the peak and improving the distance measurement.

4. Derivation 3: The Alcock-Paczy&nacute;ski Effect and BAO Distance Measures

A spherical standard ruler observed at redshift $z$ subtends both an angular size on the sky and a radial extent in redshift space. These two measurements probe different distance-redshift relations, making BAO sensitive to both $D_A(z)$ and $H(z)$ independently.

4.1 Transverse (Angular) Distance

A standard ruler of physical length $r_s$ oriented perpendicular to the line of sight subtends an angle:

$$\Delta\theta = \frac{r_s}{(1+z)\,D_A(z)}$$

where $D_A(z)$ is the angular diameter distance. Equivalently, using the comoving angular diameter distance $D_M(z) = (1+z)\,D_A(z)$:

$$\boxed{D_M(z) = \frac{r_s}{\Delta\theta}}$$

4.2 Radial (Line-of-Sight) Distance

The same ruler oriented along the line of sight spans a redshift interval:

$$\Delta z = \frac{r_s\,H(z)}{c}$$

This defines the Hubble distance:

$$\boxed{D_H(z) = \frac{c}{H(z)} = \frac{r_s}{\Delta z}}$$

4.3 Volume-Averaged Distance

When the data are not sufficient to separately constrain the transverse and radial BAO scales (as in early, lower-statistics surveys), one measures a volume-averaged, angle-averaged combination. The effective volume element at redshift $z$ scales as $D_M^2(z) \cdot D_H(z)$. This motivates the definition:

$$\boxed{D_V(z) = \left[z\,D_M^2(z)\,D_H(z)\right]^{1/3} = \left[\frac{c\,z\,D_M^2(z)}{H(z)}\right]^{1/3}}$$

The factor of $z$ is included by convention so that $D_V$ has dimensions of distance and behaves sensibly at low $z$. The key observables from BAO are the dimensionless ratios:

Transverse

$\frac{D_M(z)}{r_s}$

Radial

$\frac{D_H(z)}{r_s}$

Isotropic

$\frac{D_V(z)}{r_s}$

These ratios are independent of the absolute value of $r_s$ — the BAO measurement itself provides relative distances. Combining with CMB calibration of $r_s$ then yields absolute distances.

The Alcock-Paczy&nacute;ski Test

Alcock & Paczy&nacute;ski (1979) noted that any spherically symmetric feature, when analyzed with the wrong cosmology, will appear distorted — compressed or elongated along the line of sight relative to the transverse direction. The ratio $D_M(z) \cdot H(z) / c$ must be consistent across all redshifts for the correct cosmology. BAO provides a clean implementation of this test because the sound horizon is known to be isotropic to high precision.

5. Derivation 4: BAO as a Dark Energy Probe

5.1 Distance-Redshift Relations and Dark Energy

In a flat $\Lambda$CDM universe, the comoving distance is:

$$D_M(z) = \frac{c}{H_0}\int_0^z \frac{dz'}{E(z')}$$

where $E(z) = H(z)/H_0$ is the dimensionless Hubble parameter. For a dark energy component with equation of state $w(z)$:

$$\boxed{E^2(z) = \Omega_m(1+z)^3 + \Omega_r(1+z)^4 + \Omega_{\rm DE}\,\exp\!\left[3\int_0^z \frac{1+w(z')}{1+z'}\,dz'\right]}$$

For constant $w$, the dark energy term simplifies to $\Omega_{\rm DE}(1+z)^{3(1+w)}$. For $w = -1$ (cosmological constant), it is simply $\Omega_\Lambda$. The CPL parameterization $w(z) = w_0 + w_a\,z/(1+z)$ captures time-dependent dark energy:

$$\Omega_{\rm DE}(z) = \Omega_{\rm DE}\,(1+z)^{3(1+w_0+w_a)}\,\exp\!\left[-3w_a\,\frac{z}{1+z}\right]$$

5.2 Constraining w from BAO

BAO measurements at multiple redshifts trace out $D_M(z)/r_s$ and $D_H(z)/r_s$, which have different sensitivities to the dark energy parameters. The key insight is:

$D_H(z) = c/H(z)$ directly probes the instantaneous expansion rate, making it sensitive to $w(z)$ at each redshift
$D_M(z)$ is an integral of $1/H(z')$from 0 to $z$, so it accumulates sensitivity and provides a smoother constraint
The combination of both breaks degeneracies: $D_H$ is more sensitive to $w_a$(time evolution), while $D_M$ constrains the average equation of state

5.3 The Inverse Distance Ladder

The traditional distance ladder uses local measurements (Cepheids, tip of the red giant branch) to calibrate Type Ia supernovae, yielding $H_0$. The inverse distance ladder works in the opposite direction:

The CMB calibrates $r_s$ to sub-percent precision (from $\Omega_b h^2$and $\Omega_m h^2$)
BAO measures $D_V(z)/r_s$ at multiple redshifts, giving absolute distances
Extrapolating to $z = 0$ yields $H_0$ independently of any local calibration

$$H_0^{\rm BAO+CMB} = 67.4 \pm 0.5\;\text{km/s/Mpc}$$

This value is in $4$–$6\sigma$ tension with the SH0ES local measurement of $H_0 = 73.0 \pm 1.0$ km/s/Mpc, constituting the “Hubble tension.” Importantly, the BAO+CMB determination is completely independent of the local distance ladder, so the tension cannot be attributed to shared systematics.

DESI 2024 Constraints on Dark Energy

The Dark Energy Spectroscopic Instrument (DESI) released its first-year BAO results in April 2024, measuring $D_M/r_s$ and $D_H/r_s$ at seven effective redshifts from$z = 0.3$ to $z = 2.33$. Combined with CMB and supernova data, the DESI results show intriguing hints of evolving dark energy:

$$w_0 = -0.55^{+0.39}_{-0.21}, \qquad w_a = -1.32^{+0.58}_{-0.82}\qquad\text{(DESI + CMB + SNe)}$$

These values prefer $w_0 > -1$ and $w_a < 0$, suggesting dark energy that was weaker in the past and is becoming more dominant over time. The evidence for $w \neq -1$is at the 2–3$\sigma$ level — tantalizing but not yet conclusive.

6. Observational Applications and Survey Results

Major BAO Surveys

Survey	Years	Tracers	Redshift Range	Key Result
SDSS-II LRG	2005	46,748 LRGs	$z \approx 0.35$	First BAO detection (3.4$\sigma$)
BOSS	2009–2014	1.5M galaxies	$z = 0.2\text{--}0.75$	1% distance measurements
eBOSS	2014–2019	LRGs, ELGs, QSOs	$z = 0.6\text{--}2.2$	Extended to $z > 1$ with multiple tracers
DESI	2021–	40M+ targets	$z = 0.1\text{--}3.5$	Hints of evolving $w(z)$ (2024 DR1)
Euclid	2023–	30M+ $\text{H}\alpha$ emitters	$z = 0.9\text{--}1.8$	Sub-percent BAO from space

6.1 Lyman-$\alpha$ Forest BAO

At high redshift ($z > 2$), galaxies become too sparse for efficient BAO measurements. Instead, the Lyman-$\alpha$ forest — absorption features in quasar spectra caused by intervening neutral hydrogen — traces the matter distribution. BOSS and eBOSS measured the BAO scale in the Ly$\alpha$ forest at $z \approx 2.33$, providing the highest-redshift BAO constraints. DESI has dramatically increased the Ly$\alpha$ sample, improving these measurements further.

The Ly$\alpha$ forest BAO is particularly powerful because it probes the epoch when dark energy was subdominant — this provides a long lever arm for testing dark energy evolution. Combined with lower-redshift galaxy BAO, it constrains the full shape of $D_M(z)$ and$D_H(z)$ from $z = 0.1$ to $z \approx 2.5$.

6.2 BAO Reconstruction

Nonlinear gravitational evolution smears the BAO peak in $\xi(r)$, degrading the precision of the distance measurement. BAO reconstruction (Eisenstein et al. 2007) reverses this degradation by using the observed galaxy density field to estimate the displacement field and move galaxies back toward their original positions. This sharpens the BAO peak and typically improves the distance constraint by a factor of 1.5–2. All modern BAO analyses (BOSS, eBOSS, DESI) apply reconstruction as a standard step.

7. Historical Context

1970: Theoretical Prediction

Peebles & Yu and Sunyaev & Zel'dovich independently predicted that acoustic oscillations in the photon-baryon fluid would leave imprints in both the CMB and the matter distribution. The calculation showed that the coupling between photons and baryons before recombination would produce a characteristic scale in the galaxy correlation function corresponding to the sound horizon.

1998–2003: Pre-detection Forecasts

Eisenstein, Hu, & Tegmark (1998) and Blake & Glazebrook (2003) developed the framework for using BAO as a dark energy probe, establishing it as a primary science driver for large spectroscopic surveys.

2005: First Detection

Eisenstein et al. (2005) detected the BAO peak in the correlation function of 46,748 luminous red galaxies from SDSS at $z \approx 0.35$, at a significance of 3.4$\sigma$. Simultaneously, Cole et al. (2005) detected BAO oscillations in the power spectrum of the 2dF Galaxy Redshift Survey. These landmark papers confirmed a 35-year-old theoretical prediction and opened a new era in precision cosmology.

2012–2020: BOSS and eBOSS Era

The Baryon Oscillation Spectroscopic Survey (BOSS) achieved 1% distance measurements at$z = 0.38$ and $z = 0.61$ using 1.5 million galaxies and 160,000 Ly$\alpha$ quasars. eBOSS extended coverage to $z \approx 2.2$using emission-line galaxies and quasars as tracers.

2024: DESI First Results

The Dark Energy Spectroscopic Instrument released its Data Release 1 BAO measurements in April 2024, based on 5.7 million galaxy and quasar redshifts plus 420,000 Ly$\alpha$ forest sightlines. The DESI results provide the most precise BAO measurements to date and hint at dark energy evolution ($w_0 > -1$, $w_a < 0$) at the 2–3$\sigma$ level.

8. Numerical Exploration: BAO Correlation Function and Distance Measures

The following Python code computes two key BAO visualizations: (1) the galaxy two-point correlation function showing the BAO peak at $r \approx 150$ Mpc, and (2) the volume-averaged distance$D_V(z)/r_s$ compared with SDSS/BOSS/DESI data points.

Python

bao_analysis.py252 lines

#!/usr/bin/env python3
"""
bao_analysis.py
Compute the BAO correlation function and distance-redshift relation.
Uses numpy only (no scipy).
"""
import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

# ============================================
# Cosmological parameters (Planck 2018 LCDM)
# ============================================
H0_val = 67.36            # km/s/Mpc
h = H0_val / 100.0
c_light = 299792.458      # km/s
Omega_b = 0.02237 / h**2
Omega_c = 0.1200 / h**2
Omega_m = Omega_b + Omega_c
Omega_r = 9.15e-5
Omega_L = 1.0 - Omega_m - Omega_r
n_s = 0.9649
r_s_fid = 147.09          # Mpc (Planck sound horizon at drag epoch)

# ============================================
# Hubble parameter E(z) = H(z)/H0
# ============================================
def E_of_z(z):
    return np.sqrt(Omega_r * (1+z)**4 + Omega_m * (1+z)**3 + Omega_L)

def H_of_z(z):
    return H0_val * E_of_z(z)

# ============================================
# Comoving distance via trapezoidal integration
# ============================================
def comoving_distance(z):
    """D_M(z) in Mpc."""
    if np.isscalar(z):
        z = np.array([z])
    results = np.zeros_like(z, dtype=float)
    for i, zz in enumerate(z):
        zp = np.linspace(0, zz, 5000)
        integrand = 1.0 / E_of_z(zp)
        results[i] = np.trapezoid(integrand, zp) * c_light / H0_val
    return results

# ============================================
# Distance measures
# ============================================
def D_M(z):
    return comoving_distance(z)

def D_H(z):
    return c_light / H_of_z(z)

def D_V(z):
    dm = D_M(z)
    dh = D_H(z)
    return (z * dm**2 * dh)**(1.0/3.0)

# ============================================
# FIGURE 1: Galaxy Correlation Function with BAO
# ============================================
# Build a model correlation function using the Eisenstein & Hu
# transfer function approach (simplified).

def transfer_EH_nowiggle(k, Om, Ob, h_val):
    """Eisenstein & Hu no-wiggle transfer function."""
    Omh2 = Om * h_val**2
    Obh2 = Ob * h_val**2
    f_b = Ob / Om
    theta = 2.7255 / 2.7  # T_CMB/2.7
    z_eq = 2.5e4 * Omh2 * theta**(-4)
    k_eq = 7.46e-2 * Omh2 * theta**(-2)  # Mpc^-1
    s = 44.5 * np.log(9.83 / Omh2) / np.sqrt(1.0 + 10.0 * Obh2**0.75)
    alpha_G = 1.0 - 0.328 * np.log(431.0 * Omh2) * f_b + 0.38 * np.log(22.3 * Omh2) * f_b**2
    q = k * theta**2 / (Omh2 * np.exp(-Ob * (1.0 + np.sqrt(2*h_val)/Ob)))
    q = k / (13.41 * k_eq)
    L = np.log(2*np.e + 1.8 * q)
    C = 14.2 / alpha_G + 386.0 / (1.0 + 69.9 * q**1.08)
    T0 = L / (L + C * q**2)
    return T0

def pk_model(k, with_bao=True):
    """Model P(k) with or without BAO wiggles."""
    T_nw = transfer_EH_nowiggle(k, Omega_m, Omega_b, h)
    Pk_nw = k**n_s * T_nw**2

if with_bao:
        # Add BAO oscillations
        f_b = Omega_b / Omega_m
        Sigma_nl = 8.0   # Mpc, nonlinear damping scale
        A_bao = f_b * 0.75
        wiggle = A_bao * np.sin(k * r_s_fid) / (k * r_s_fid) * np.exp(-0.5 * (k * Sigma_nl)**2)
        Pk = Pk_nw * (1.0 + wiggle)
    else:
        Pk = Pk_nw
    return Pk

# Compute correlation function via numerical Hankel transform
# xi(r) = int P(k) j0(kr) k^2 dk / (2 pi^2)
def compute_xi(r_arr, with_bao=True):
    """Compute correlation function from P(k)."""
    k_arr = np.logspace(-4, 1, 8000)
    dk = np.diff(k_arr)
    k_mid = 0.5 * (k_arr[:-1] + k_arr[1:])
    Pk_mid = pk_model(k_mid, with_bao=with_bao)

xi_out = np.zeros_like(r_arr)
    for i, r in enumerate(r_arr):
        kr = k_mid * r
        j0 = np.sin(kr) / kr
        integrand = Pk_mid * j0 * k_mid**2 / (2.0 * np.pi**2)
        xi_out[i] = np.sum(integrand * dk)
    return xi_out

# Compute correlation functions
r_arr = np.linspace(20, 200, 400)

print("Computing correlation function with BAO...")
xi_bao = compute_xi(r_arr, with_bao=True)
print("Computing smooth (no-BAO) correlation function...")
xi_nw = compute_xi(r_arr, with_bao=False)

# Scale for plotting: multiply by r^2
xi_r2_bao = r_arr**2 * xi_bao
xi_r2_nw = r_arr**2 * xi_nw

# Normalize for display
norm = np.max(np.abs(xi_r2_bao))
xi_r2_bao /= norm
xi_r2_nw /= norm

# ============================================
# FIGURE 2: D_V(z)/r_s vs z with data
# ============================================
# BAO data points: (z_eff, D_V/r_s, sigma)
# From SDSS, BOSS, eBOSS, and DESI
bao_data = [
    # SDSS-II (6dFGS, MGS)
    (0.106, 2.98, 0.13, 'SDSS-6dFGS'),
    (0.15,  4.47, 0.17, 'SDSS-MGS'),
    # BOSS DR12
    (0.38, 10.23, 0.17, 'BOSS'),
    (0.51, 13.36, 0.21, 'BOSS'),
    (0.61, 15.45, 0.26, 'BOSS'),
    # eBOSS
    (0.70, 17.86, 0.33, 'eBOSS'),
    (1.48, 30.69, 0.80, 'eBOSS'),
    # DESI DR1 (2024)
    (0.30,  8.09, 0.15, 'DESI'),
    (0.51, 13.62, 0.18, 'DESI'),
    (0.71, 18.33, 0.22, 'DESI'),
    (0.93, 22.43, 0.48, 'DESI'),
    (1.32, 27.79, 0.69, 'DESI'),
    (2.33, 39.71, 0.94, 'DESI'),
]

# Model D_V/r_s
z_model = np.linspace(0.01, 2.8, 500)
DV_model = D_V(z_model) / r_s_fid

# ============================================
# Create the figure
# ============================================
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(16, 6.5))

# --- Panel 1: Correlation function ---
ax1.plot(r_arr, xi_r2_bao, 'b-', linewidth=2.0, label='With BAO wiggles')
ax1.plot(r_arr, xi_r2_nw, 'r--', linewidth=1.5, alpha=0.7, label='Smooth (no BAO)')
ax1.axvline(x=r_s_fid, color='green', linestyle=':', linewidth=1.5, alpha=0.7,
            label=f'Sound horizon $r_s$ = {r_s_fid:.0f} Mpc')
ax1.set_xlabel('Comoving Separation r [Mpc]', fontsize=13)
ax1.set_ylabel(r'$r^2 \xi(r)$ [arbitrary units]', fontsize=13)
ax1.set_title('Galaxy Two-Point Correlation Function', fontsize=14, fontweight='bold')
ax1.legend(fontsize=11, loc='upper right')
ax1.set_xlim(20, 200)
ax1.grid(True, alpha=0.3)

# Annotate the BAO peak
peak_idx = np.argmax(xi_r2_bao[r_arr > 100])
peak_r = r_arr[r_arr > 100][peak_idx]
peak_val = xi_r2_bao[r_arr > 100][peak_idx]
ax1.annotate('BAO peak', xy=(peak_r, peak_val),
             xytext=(peak_r + 15, peak_val + 0.1),
             fontsize=12, color='blue',
             arrowprops=dict(arrowstyle='->', color='blue', lw=1.5))

# --- Panel 2: D_V(z)/r_s vs z ---
ax2.plot(z_model, DV_model, 'k-', linewidth=2.0, label=r'$\Lambda$CDM (Planck 2018)')

# Plot data points by survey
colors_map = {'SDSS-6dFGS': '#e74c3c', 'SDSS-MGS': '#e74c3c',
              'BOSS': '#3498db', 'eBOSS': '#2ecc71', 'DESI': '#9b59b6'}
markers_map = {'SDSS-6dFGS': 's', 'SDSS-MGS': 's',
               'BOSS': 'o', 'eBOSS': 'D', 'DESI': '^'}
plotted = set()
for z_d, dv_d, sig_d, label_d in bao_data:
    show_label = label_d if label_d not in plotted else None
    plotted.add(label_d)
    ax2.errorbar(z_d, dv_d, yerr=sig_d, fmt=markers_map[label_d],
                 color=colors_map[label_d], markersize=7, capsize=3,
                 label=show_label, zorder=5)

ax2.set_xlabel('Redshift z', fontsize=13)
ax2.set_ylabel(r'$D_V(z) / r_s$', fontsize=13)
ax2.set_title('BAO Distance-Redshift Relation', fontsize=14, fontweight='bold')
ax2.legend(fontsize=10, loc='upper left')
ax2.set_xlim(0, 2.8)
ax2.set_ylim(0, 45)
ax2.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('bao_analysis.png', dpi=150, bbox_inches='tight')
plt.show()

# ============================================
# Print summary statistics
# ============================================
print("\n" + "="*60)
print("BARYON ACOUSTIC OSCILLATIONS: KEY PARAMETERS")
print("="*60)

# Sound horizon computation
print(f"\nSound horizon at drag epoch: r_s = {r_s_fid:.2f} Mpc")
print(f"BAO wavenumber: k_BAO = 2pi/r_s = {2*np.pi/r_s_fid:.4f} Mpc^-1")

# Sound speed at drag epoch
z_drag = 1059.94
Omega_g = 2.47e-5 / h**2
R_drag = 3 * Omega_b / (4 * Omega_g) / (1 + z_drag)
c_s_drag = 1.0 / np.sqrt(3 * (1 + R_drag))
print(f"\nBaryon loading at drag epoch: R(z_drag) = {R_drag:.4f}")
print(f"Sound speed at drag epoch: c_s = {c_s_drag:.4f} c")
print(f"Sound speed at early times (R->0): c_s = {1/np.sqrt(3):.4f} c")

# Distance measures at key redshifts
print(f"\nDistance measures (Planck LCDM):")
print(f"{'z':>6s}  {'D_M [Mpc]':>12s}  {'D_H [Mpc]':>12s}  {'D_V [Mpc]':>12s}  {'D_V/r_s':>10s}")
print("-" * 60)
for z_val in [0.15, 0.38, 0.51, 0.61, 0.70, 1.0, 1.5, 2.0, 2.33]:
    z_a = np.array([z_val])
    dm = D_M(z_a)[0]
    dh = D_H(z_a)
    dv = D_V(z_a)[0]
    print(f"{z_val:6.2f}  {dm:12.1f}  {dh:12.1f}  {dv:12.1f}  {dv/r_s_fid:10.2f}")

print(f"\nH0 from inverse distance ladder (BAO+CMB): ~67.4 km/s/Mpc")
print(f"H0 from local measurements (SH0ES): ~73.0 km/s/Mpc")
print(f"Tension: ~5 sigma")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

Code Output

The script generates two panels:

Left: $r^2 \xi(r)$ showing the BAO peak at $r \approx 150$ Mpc (blue) vs smooth reference (red dashed)
Right: $D_V(z)/r_s$ vs $z$ with SDSS/BOSS/eBOSS/DESI data and $\Lambda$CDM prediction

Uses simplified Eisenstein & Hu transfer function; real analyses use CAMB/CLASS.

9. Summary of Key Results

Photon-Baryon Sound Speed

$$c_s = \frac{c}{\sqrt{3(1+R)}}, \qquad R = \frac{3\rho_b}{4\rho_\gamma} \approx 0.63\;\text{at}\;z_{\rm drag}$$

Sound Horizon (Standard Ruler)

$$r_s = \int_{z_{\rm drag}}^{\infty} \frac{c_s(z)}{H(z)}\,dz \approx 147.09\;\text{Mpc}$$

Correlation Function & Power Spectrum

$$\xi(r) = \int_0^\infty P(k)\,j_0(kr)\,\frac{k^2\,dk}{2\pi^2}, \qquad k_{\rm BAO} = \frac{2\pi}{r_s} \approx 0.043\;\text{Mpc}^{-1}$$

BAO Distance Measures

$$D_M(z) = \frac{r_s}{\Delta\theta}, \qquad D_H(z) = \frac{c}{H(z)}, \qquad D_V(z) = \left[\frac{c\,z\,D_M^2(z)}{H(z)}\right]^{1/3}$$

Dark Energy Equation of State

$$E^2(z) = \Omega_m(1+z)^3 + \Omega_{\rm DE}(1+z)^{3(1+w_0+w_a)}\exp\!\left[-3w_a\frac{z}{1+z}\right]$$

Inverse Distance Ladder

$$\text{CMB} \xrightarrow{r_s} \text{BAO}(D_V/r_s) \xrightarrow{z \to 0} H_0 = 67.4 \pm 0.5\;\text{km/s/Mpc}$$

BAO: Key Connections

1. Sound waves in the primordial plasma → frozen at decoupling into a $\sim 150$ Mpc standard ruler

2. BAO peak in $\xi(r)$ → oscillations in $P(k)$ (Fourier duality)

3. Transverse BAO → $D_A(z)/r_s$ (angular diameter distance)

4. Radial BAO → $H(z)\,r_s$ (expansion rate)

5. Alcock-Paczy&nacute;ski test → consistency check on cosmological model

6. CMB + BAO inverse distance ladder → $H_0$ independent of local measurements

7. Multi-redshift BAO → constraints on $w_0$ and $w_a$ (dark energy evolution)

8. DESI 2024 hints → possible $w(z) \neq -1$ at 2–3$\sigma$

Bibliography

Textbooks & Reviews

Dodelson, S. & Schmidt, F. (2020). Modern Cosmology, 2nd ed. Academic Press.
Weinberg, D.H. et al. (2013). “Observational Probes of Cosmic Acceleration,” Phys. Rep. 530, 87. arXiv:1201.2434.
Bassett, B.A. & Hlozek, R. (2010). “Baryon Acoustic Oscillations,” in Dark Energy, CUP. arXiv:0910.5224.

Key Papers

Peebles, P.J.E. & Yu, J.T. (1970). ApJ 162, 815. — Theoretical prediction of BAO.
Sunyaev, R.A. & Zel'dovich, Ya.B. (1970). Ap&SS 7, 3. — Independent prediction of acoustic oscillations.
Eisenstein, D.J. et al. (2005). ApJ 633, 560. arXiv:astro-ph/0501171. — First BAO detection.
Cole, S. et al. (2005). MNRAS 362, 505. — BAO in the 2dFGRS.
Eisenstein, D.J., Seo, H.-J. & White, M. (2007). ApJ 664, 660. — BAO reconstruction.
Alam, S. et al. (BOSS) (2017). MNRAS 470, 2617. — Final BOSS BAO.
DESI Collaboration (2024). arXiv:2404.03000. — DESI DR1 BAO and dark energy constraints.
Alcock, C. & Paczy&nacute;ski, B. (1979). Nature 281, 358. — The AP geometric test.

← Course Overview CMB Physics →