Cosmic Microwave Background v1

All 24-bit rgb colors mapped to actual Planck cmb temperature data. Made using Chatgpt and Google Colab.

#!/usr/bin/env python3
"""
plot_cmb_map.py

Downloads the Planck 2018 SMICA CMB temperature map (nosz version) from IRSA,
reads it with Astropy + Healpy, projects it, and produces a 4096×4096 PNG
using every 24-bit RGB color exactly once—mapping pixel intensities
by sorting the full RGB spectrum by luminance to create a meaningful heat
scale and assigning colors by rank order of CMB temperature.
The projected map (~2048×4096) is resampled (each row repeated) to
exactly fill a square 4096×4096 canvas before color mapping.
"""

import os
import requests
from astropy.io import fits
import healpy as hp
import numpy as np
from PIL import Image

# IRSA URL for the SMICA “nosz” full-mission CMB map
URL = (
    "https://irsa.ipac.caltech.edu/data/Planck/release_3/"
    "all-sky-maps/maps/component-maps/cmb/"
    "COM_CMB_IQU-smica-nosz_2048_R3.00_full.fits"
)
FILENAME = os.path.basename(URL)
OUTFILE = "cmb_rgb_4096.png"


def download_map(url: str, filename: str) -> None:
    """Download the FITS file if not present."""
    if not os.path.exists(filename):
        r = requests.get(url, stream=True)
        r.raise_for_status()
        with open(filename, "wb") as f:
            for chunk in r.iter_content(4096):
                f.write(chunk)


def read_cmb_map(filename: str):
    """Read I_STOKES from the FITS file."""
    with fits.open(filename, memmap=True) as hdul:
        return hdul[1].data["I_STOKES"]


def project_map(cmb_data, xsize=4096):
    """Project CMB data to a plate carrée array (~2048×xsize)."""
    return hp.cartview(
        cmb_data,
        return_projected_map=True,
        xsize=xsize,
        cbar=False,
        title=None,
        notext=True
    )


def resample_to_square(proj):
    """Resample projected map by repeating each row to reach square shape."""
    h, w = proj.shape
    # Compute repeat factor (should be integer, here 2)
    factor = w // h
    return np.repeat(proj, factor, axis=0)


def map_to_rgb(proj_square):
    """
    Map resampled data to every 24-bit RGB color exactly once.
    Flatten the data, rank-order by temperature,
    generate all RGB triplets sorted by luminance,
    and assign each pixel a unique color.
    """
    h, w = proj_square.shape
    N = h * w
    flat = proj_square.ravel()
    order = np.argsort(flat)

    idx = np.arange(N, dtype=np.uint32)
    r = (idx >> 16) & 0xFF
    g = (idx >> 8) & 0xFF
    b = idx & 0xFF
    lum = 0.2126 * r + 0.7152 * g + 0.0722 * b
    color_order = np.argsort(lum)

    rgb_flat = np.empty((N, 3), dtype=np.uint8)
    colors = np.stack((r, g, b), axis=1)[color_order]
    rgb_flat[order] = colors

    return rgb_flat.reshape((h, w, 3))


def save_image(rgb_img, output: str = OUTFILE):
    """Save the RGB image at exact resolution."""
    img = Image.fromarray(rgb_img, mode='RGB')
    img.save(output)


def main():
    download_map(URL, FILENAME)
    cmb = read_cmb_map(FILENAME)
    proj = project_map(cmb, xsize=4096)
    proj_sq = resample_to_square(proj)
    rgb = map_to_rgb(proj_sq)
    save_image(rgb)
    print(f"Saved RGB image: {OUTFILE} (4096×4096, full 24-bit palette)")


if __name__ == "__main__":
    main()