Plot hubble constant

Hubble_plot.py

import matplotlib.pyplot as plt

# Define data: (label, value, uncertainty, category)
# https://lambda.gsfc.nasa.gov/education/graphic_history/hubb_const.html
data = [
    ("Cepheids+SNIa 2021", 73.04, 1.04, "Late"),
    ("LIGO+Virgo+KAGRA 2021", 68, 9.5, "Late"),
    ("SBF_distances 2021", 73.3, 2.5, "Late"),
    ("eBOSS_BBN+BAO 2020", 67.35, 0.97, "Early"),
    ("Megamasers 2020", 73.9, 3.0, "Late"),
    ("TRGB_Dist_Ladder 2019", 69.8, 1.9, "Late"),
    ("GravLens_Time_Delay 2019", 73.3, 1.75, "Late"),
    ("XMM+Planck_tSZ 2018", 67, 3.0, "Early"),
    ("Inv_Dist_Ladder 2015", 67.3, 1.1, "Late"),
    ("Cepheids+SNIa 2011", 73.8, 2.4, "Late"),
    ("CHANDRA+tSZ 2006", 73.7, 9.5, "Early"),
    ("HST_Key_Project 2001", 72, 8.0, "Late"),
    ("Planck+ACTPol+SPTpol_EE 2021", 68.7, 1.3, "Early"),
    ("SPT-3G_TE+EE 2021", 68.8, 1.5, "Early"),
    ("ACTPol_DR4 2020", 67.9, 1.5, "Early"),
    ("SPTPol 2017", 71.2, 2.1, "Early"),
    ("SPT 2013", 75.0, 3.5, "Early"),
]

# Split into early and late
early = [(label, v, e) for label, v, e, cat in data if cat == "Early"]
late = [(label, v, e) for label, v, e, cat in data if cat == "Late"]

# Plot setup
fig, ax = plt.subplots(figsize=(9, 10))
y_pos = list(range(len(data)))[::-1]  # top to bottom

# Plot each point
for i, (label, value, err, cat) in enumerate(data):
    color = "tab:blue" if cat == "Early" else "tab:orange"
    ax.errorbar(value, y_pos[i], xerr=err, fmt='o', color='black', ecolor=color, elinewidth=3, capsize=5)

# Labels and formatting
ax.set_yticks(y_pos)
ax.set_yticklabels([label for label, _, _, _ in data])
ax.set_xlabel("Hubble Constant (km/s/Mpc)")
ax.set_title("Hubble Constant Measurements (Early vs Late Universe)")
ax.invert_yaxis()
ax.grid(axis='x', linestyle='--', alpha=0.6)

# Weighted average function
def weighted_avg(values):
    weights = [1 / (err ** 2) for _, _, err in values]
    weighted_values = [v * w for (_, v, err), w in zip(values, weights)]
    avg = sum(weighted_values) / sum(weights)
    std = (1 / sum(weights)) ** 0.5
    return avg, std

# Compute averages
early_avg, early_std = weighted_avg(early)
late_avg, late_std = weighted_avg(late)

# Plot horizontal average lines
ax.axvline(early_avg, color='tab:blue', linestyle='--', label=f"Early avg: {early_avg:.2f} ± {early_std:.2f}")
ax.axvline(late_avg, color='tab:orange', linestyle='--', label=f"Late avg: {late_avg:.2f} ± {late_std:.2f}")
ax.legend()

plt.tight_layout()
plt.show()