montepython_public icon indicating copy to clipboard operation
montepython_public copied to clipboard

BK18 likelihood added

Open ThomasTram opened this issue 10 months ago • 2 comments

I implemented the BK18 likelihood based on the previous BK14 likelihood. I initially wanted to do it by inheriting from the BK14 likelihood class, but in the end a lot of small things had changed so I decided to copy-paste it and then edit it. There are actually two likelihoods, one using only B-modes (which is the default) called BK18lf_dust and one that includes E-modes which is called BK18lf_dust_incEE. The two likelihoods are identical and differ only in the datafiles and the maps used inside the .data file.

To validate the likelihood I recreated Figure 4 of BKXIII in 2110.00483. Using this script:

import getdist.plots as gdplt
from getdist import MCSamples, loadMCSamples
from pathlib import Path

# Load chain file
samples = loadMCSamples('chains/BK18_fig4_20/2024-04-15_20000_', settings={"ignore_rows": 0.3})
samples.paramNames.setLabelsFromParamNames('chains/BK18_fig4_20/2024-04-15_20000_.paramnames')
# Plot settings
plot_settings = {
    "smooth_scale_1D": 0.15,
    "smooth_scale_2D": 2,
    "num_bins": 100,
    "num_bins_2D": 40,
    "contours": [0.68, 0.95],
    "triangle_plot": True,
    "plot_2D_param": 0,
    "plot_2D_num": 1,
    "params": ["r", "BBdust", "BBsync"],
    "triangle_params": ["r", "BBdust", "BBsync"],
    "param_limits": {"r": [0, 0.18], "BBdust": [0, 11], "BBsync": [0, 8], "BBdustsynccorr": [-1, 1], "BBbetadust": [1, 2], 
                     "BBbetasync": [-4.2, -1.8], "BBalphadust": [-1, 0], "BBalphasync": [-1, 0]},
}

# Make triangle plot of figure 4 of BKXIII
gdplt.get_subplot_plotter()
gdplot = gdplt.getSubplotPlotter()
gdplot.settings.num_plot_contours = 2
gdplot.settings.alpha_filled_add = 0.6
gdplot.triangle_plot(samples, **plot_settings, filled=True)

# Make 1d posterior plots of figure 4 of BKXIII
plot_settings['params'] = ["BBbetadust", "BBbetasync", "BBdustsynccorr", "BBalphadust", "BBalphasync"]
plot_settings['xlims'] = [plot_settings['param_limits'][p] for p in plot_settings['params']]
gdplot.plots_1d(samples, **plot_settings)

I get these figures: image

image

which at least by eye agree.

ThomasTram avatar Apr 15 '24 13:04 ThomasTram

Hi Thomas,

Thanks a lot for doing this, would you like to go ahead and add it to the private devel branch? I can also do it, but will be quite busy the next three weeks and might not get around to it soon. Then once we have everything there (e.g. DESI BAO and new WMAP Python 3 wrapper from @schoeneberg , NPIPE, ground based CMB, etc) we can put out a new release.

If there are any other fixes or improvements you've done that haven't been included yet (e.g. PR #372 ) you can also add those. Thanks a lot!

Best, Thejs

brinckmann avatar May 17 '24 10:05 brinckmann

@brinckmann I actually forgot that I had access to the private repo. I can add the changes yes.

ThomasTram avatar May 21 '24 09:05 ThomasTram