Summary¶
Comparing spectra from
- POLARBEAR's Planck 2015 CAMB input,
- theory spectra from Planck 2018 release
- class command line output
- class Python interface output
class' agrees with first 2 when setting format as camb.
class without setting this format, and the result from the Python interface doesn't agree. Why?
In [1]:
from classy import Class
In [2]:
import os
import io
In [3]:
import pandas as pd
import numpy as np
In [4]:
import matplotlib.pyplot as plt
In [5]:
%matplotlib inline
In [6]:
def read_txt(path):
with open(path, 'r') as f:
# remove beginning `#`
text = f.read()[1:]
with io.StringIO(text) as f:
# input has 6 sig. fig. which fits in float32
return pd.read_csv(f, delim_whitespace=True, index_col=0, dtype=np.float32)
In [7]:
def read_class_txt(path, camb=False):
# read once
with open(path, 'r') as f:
text = f.read()
# get last comment line
comment = None
for line in text.split('\n'):
if line.startswith('#'):
comment = line
# remove beginning '#'
comment = comment[1:]
# parse comment line: get name after ':'
names = [name.split(':')[1] for name in comment.strip().split()]
with io.StringIO(text) as f:
df = pd.read_csv(f, delim_whitespace=True, index_col=0, comment='#', header=None, names=names)
return df if camb else df * 1.e12
In [30]:
def plot_compare(df1, df2, keys, log=False, relative=False):
'''plot to compare columns from `df1` and `df2`,
using columns from `df1`.
`keys`: tuple of str of names of `df1` and `df2`.
`relative`: if True, plot the relative error w.r.t. `df1` instead.
'''
for col in df1.columns:
# pd.merge(df1[col], df2[col], left_index=True, right_index=True).plot()
df_temp = pd.concat((df1[col], df2[col]), axis=1, join='inner', keys=(' '.join((key, col)) for key in keys))
if relative:
temp = df_temp.values
if log:
plt.loglog(df_temp.index, np.abs((temp[:, 0] - temp[:, 1]) / temp[:, 0]), label=col)
else:
plt.plot(df_temp.index, np.abs((temp[:, 0] - temp[:, 1]) / temp[:, 0]), label=col)
plt.legend()
plt.show()
else:
if log:
df_temp.plot(logx=True, logy=True)
else:
df_temp.plot()
Planck 2018 release¶
c.f.
- wiki
- filename:
COM_PowerSpect_CMB-base-plikHM_TTTEEE-lowl-lowE-lensing-minimum-theory_R3.01.txt - archive: http://pla.esac.esa.int/pla/
In [9]:
path = os.path.expanduser('~/temp/COM_PowerSpect_CMB-base-plikHM_TTTEEE-lowl-lowE-lensing-minimum-theory_R3.01.txt')
In [10]:
!head -n 20 $path
In [11]:
df_planck = read_txt(path)
In [12]:
df_planck.head()
Out[12]:
from CAMB¶
In [13]:
df_camb = pd.read_hdf('/scratch/data/largepatch_high/Planck2015_base_plikHM_TT_lowTEB_lensing_r0p025_lensedtotCls.hdf5')
In [14]:
df_camb.head()
Out[14]:
In [15]:
plot_compare(df_camb, df_planck, ('CAMB', 'Planck'))
In [31]:
plot_compare(df_camb, df_planck, ('CAMB', 'Planck'), relative=True, log=True)
Class result from command line¶
In [17]:
path = os.path.expanduser('~/git/fork/class/output/base_2018_plikHM_TTTEEE_lowl_lowE_lensing_cl_lensed.dat')
In [18]:
!head -n 20 $path
In [19]:
df_class = read_class_txt(path, camb=True)
df_class.head()
Out[19]:
In [33]:
plot_compare(df_class[['TT', 'EE', 'BB', 'TE']], df_planck, ('CAMB', 'Planck'), relative=True, log=True)
Class¶
In [21]:
LambdaCDM = Class()
In [22]:
# optional: clear content of LambdaCDM (to reuse it for another model)
LambdaCDM.struct_cleanup()
# optional: reset parameters to default
LambdaCDM.empty()
In [23]:
kwargs = {
# background parameters
'H0': 67.32117,
'omega_b': 0.02238280,
'N_ur': 2.03066666667,
'omega_cdm': 0.1201075,
'N_ncdm': 1,
'omega_ncdm': 0.0006451439,
'YHe': 0.2454006,
'tau_reio': 0.05430842,
'n_s': 0.9660499,
'A_s': 2.100549e-09,
'output': 'tCl,pCl,lCl,mPk',
'lensing': 'yes',
'P_k_max_1/Mpc': 3.,
'l_max_scalars': 3000,
}
LambdaCDM.set(kwargs)
Out[23]:
In [24]:
LambdaCDM.compute()
In [25]:
df = pd.DataFrame(LambdaCDM.lensed_cl())
df.set_index('ell', inplace=True)
ell = df.index.values.astype(np.int32)
df *= ((ell * (ell + 1)) * 0.5 / np.pi)[:, None]
df.columns = df.columns.str.upper()
df *= 1.e12
In [26]:
df.head()
Out[26]:
In [27]:
df_planck.head()
Out[27]:
In [35]:
plot_compare(df_planck, df, ('planck', 'classy'), relative=True)
In [ ]: