The Kinetic Sunyaev–Zeldovich Effect

Let's make some maps of the kinetic Sunyaev-Zeldovich (SZ) effect. First, let's load up the example halo catalog included in this package. These will be automatically downloaded the first time you load them.

using XGPaint, Plots

# example , ra and dec in radians, halo mass in M200c (Msun)
ra, dec, redshift, halo_mass = XGPaint.load_example_halos()
print("Number of halos: ", length(halo_mass))
Number of halos: 2455621

This small catalog is limited to a relatively small patch of the sky. Before we generate some SZ maps, let's take a look at the halo mass distribution.

b = 10.0 .^ (11:0.25:16)
histogram(halo_mass, bins=b, xaxis=(:log10, (1e11, 1e16)),
    yscale=:log10, label="", xlabel="Halo mass (solar masses)", ylabel="counts")
Example block output

To allow for safe threaded painting on a single map, we'll sort the halo catalog by declinations.

ra, dec, redshift, halo_mass = sort_halo_catalog(ra, dec, redshift, halo_mass);
(Float32[0.03305517, -0.052936707, 0.07140066, -0.06599737, -0.05942857, -0.05192846, 0.022953046, 0.01955623, 0.05080256, 0.00839793  …  0.060685053, -0.061119348, -0.040161338, 0.047465798, -0.041020952, -0.048511826, -0.018889412, 0.058788735, -0.050192643, -0.0027171609], Float32[-0.11871151, -0.118711114, -0.118711084, -0.11871107, -0.1187108, -0.118710734, -0.11871068, -0.118710674, -0.11871066, -0.11871049  …  0.10907136, 0.1090714, 0.10907147, 0.109071486, 0.1090715, 0.109071545, 0.10907169, 0.10907171, 0.10907183, 0.10907187], Float32[2.583177, 0.6798465, 1.4767964, 2.4593022, 1.3147808, 1.7026237, 1.313409, 1.2228533, 2.061896, 1.7019778  …  0.4401552, 3.107578, 1.2217613, 2.3513637, 2.2199023, 1.7342392, 0.6818484, 0.8421148, 1.2661853, 2.0941606], Float32[1.8216356f12, 3.4135917f12, 4.958732f12, 3.3517478f12, 2.6270643f12, 2.1461618f12, 1.5084386f12, 4.1571304f12, 1.3568221f12, 2.2222684f12  …  4.438609f12, 1.3447792f12, 1.5955485f12, 1.7920556f12, 8.5669514f12, 2.6231164f12, 2.2672384f12, 2.409821f13, 5.5713213f12, 2.7740707f12])

Now, we'll set up the map to generate. We will construct a small CAR (Clenshaw-Curtis variant) patch on the sky. You have to construct a workspace for each new sky patch.

using Pixell
box = [4.5   -4.5;           # RA
       -3     3] * Pixell.degree  # DEC
shape, wcs = geometry(CarClenshawCurtis{Float64}, box, 0.5 * Pixell.arcminute)
m = Enmap(zeros(shape), wcs)
1080×720 Pixell.Enmap{Float64, 2, Matrix{Float64}, Pixell.CarClenshawCurtis{Float64}}:
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  …  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  …  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 ⋮                        ⋮              ⋱            ⋮                   
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  …  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0

Now let's set up an electron profile.

p = BattagliaTauProfile(Omega_c=0.267, Omega_b=0.0493,  h=0.6712)
BattagliaTauProfile{Float64, Cosmology.FlatLCDM{Float64}}(0.1558646854252292, Cosmology.FlatLCDM{Float64}(0.6712, 0.6836071403864511, 0.3163, 9.285961354891755e-5))