Forecasting

bb4cast/noise_calc.py

@@ -0,0 +1,98 @@
+from __future__ import print_function
+import numpy as np
+from scipy.interpolate import interp1d
+
+
+def so_SAT_beams(ls):
+    fwhm = np.array([91., 63., 30., 17., 11., 9.])
+    sig = fwhm*np.pi/(np.sqrt(8.*np.log(2))*180.*60.)
+    return np.array([np.exp(-0.5*ls*(ls+1)*sig**2)
+                     for s in sig])
+
+
+def so_SAT_Nl(info, ell_max, include_kludge=True):
+    """
+    sensitivity_mode
+         1: baseline,
+         2: goal
+    one_over_f_mode
+         0: pessimistic
+         1: optimistic
+    """
+    sens_dict = {'baseline': 1, 'goal': 2}
+    oof_dict = {'pessimistic': 0, 'optimistic': 1}
+    sensitivity_mode = sens_dict[info.get('sensitivity', 'baseline')]
+    ydet_LF = info.get('ydet_LF', 1.)
+    ydet_MF = info.get('ydet_MF', 9.)
+    ydet_UHF = info.get('ydet_UHF', 5.)
+    f_sky = info.get('f_sky', 0.1)
+    S_SA_27 = np.array([1.e9, 21, 15]) * np.sqrt(1./ydet_LF)
+    S_SA_39 = np.array([1.e9, 13, 10]) * np.sqrt(1./ydet_LF)
+    # Factor 2 because numbers assumed a default 2 MF tubes
+    S_SA_93 = np.array([1.e9, 3.4, 2.4]) * np.sqrt(2./ydet_MF)
+    S_SA_145 = np.array([1.e9, 4.3, 2.7]) * np.sqrt(2./ydet_MF)
+    S_SA_225 = np.array([1.e9, 8.6, 5.7]) * np.sqrt(1./ydet_UHF)
+    S_SA_280 = np.array([1.e9, 22, 14]) * np.sqrt(1./ydet_UHF)
+
+    t = 365. * 24. * 3600  # 1Y in s
+    t = t * 0.2  # retention after observing efficiency and cuts
+    if include_kludge:
+        t = t * 0.85  # a kludge to account for map edges
+    A_SR = 4 * np.pi * f_sky  # sky area in steradians
+
+    ell = np.arange(ell_max+1)
+    # Avoid division by zero
+    ell[0] = 1
+
+    # White
+    W_T_27 = S_SA_27[sensitivity_mode] / np.sqrt(t)
+    W_T_39 = S_SA_39[sensitivity_mode] / np.sqrt(t)
+    W_T_93 = S_SA_93[sensitivity_mode] / np.sqrt(t)
+    W_T_145 = S_SA_145[sensitivity_mode] / np.sqrt(t)
+    W_T_225 = S_SA_225[sensitivity_mode] / np.sqrt(t)
+    W_T_280 = S_SA_280[sensitivity_mode] / np.sqrt(t)
+
+    # 1/f
+    oof_key = info.get('one_over_f', 'optimistic')
+    if oof_key == 'custom':
+        # If custom 1/f, read transfer function
+        fname_transfer = info.get('transfer', None)
+        if fname_transfer is None:
+            tell = np.ones_like(ell)
+        else:
+            ls, tls = np.loadtxt(fname_transfer, unpack=True)
+            tlf = interp1d(ls, tls,
+                           fill_value=(tls[0], tls[-1]),
+                           bounds_error=False)
+            tell = tlf(ell)
+        # One-over-f modulation is 1/(transfer function)
+        AN_P_27 = AN_P_39 = AN_P_93 = AN_P_145 = AN_P_225 = AN_P_280 = 1./tell
+    else:
+        alpha_pol = np.array([-2.4, -2.4, -2.5, -3, -3, -3])
+        one_over_f_mode = oof_dict[info.get('one_over_f', 'optimistic')]
+        f_knee_pol_SA_27 = np.array([30., 15.])
+        f_knee_pol_SA_39 = np.array([30., 15.])  # from QUIET
+        f_knee_pol_SA_93 = np.array([50., 25.])
+        f_knee_pol_SA_145 = np.array([50., 25.])  # from ABS
+        f_knee_pol_SA_225 = np.array([70., 35.])
+        f_knee_pol_SA_280 = np.array([100., 40.])
+        AN_P_27 = (ell / f_knee_pol_SA_27[one_over_f_mode])**alpha_pol[0] + 1.
+        AN_P_39 = (ell / f_knee_pol_SA_39[one_over_f_mode])**alpha_pol[1] + 1.
+        AN_P_93 = (ell / f_knee_pol_SA_93[one_over_f_mode])**alpha_pol[2] + 1.
+        AN_P_145 = (ell / f_knee_pol_SA_145[one_over_f_mode])**alpha_pol[3] + 1.
+        AN_P_225 = (ell / f_knee_pol_SA_225[one_over_f_mode])**alpha_pol[4] + 1.
+        AN_P_280 = (ell / f_knee_pol_SA_280[one_over_f_mode])**alpha_pol[5] + 1.
+
+    # combined
+    N_ell_P_27 = (W_T_27 * np.sqrt(2))**2. * A_SR * AN_P_27
+    N_ell_P_39 = (W_T_39 * np.sqrt(2))**2. * A_SR * AN_P_39
+    N_ell_P_93 = (W_T_93 * np.sqrt(2))**2. * A_SR * AN_P_93
+    N_ell_P_145 = (W_T_145 * np.sqrt(2))**2. * A_SR * AN_P_145
+    N_ell_P_225 = (W_T_225 * np.sqrt(2))**2. * A_SR * AN_P_225
+    N_ell_P_280 = (W_T_280 * np.sqrt(2))**2. * A_SR * AN_P_280
+    N_ell_P_SA = np.array([N_ell_P_27, N_ell_P_39, N_ell_P_93,
+                           N_ell_P_145, N_ell_P_225, N_ell_P_280])
+    N_ell_P_SA[:, :2] = 0
+    ell[0] = 0
+
+    return ell, N_ell_P_SA

bb4cast/params/fgs_fiducial.yml

@@ -0,0 +1,13 @@
+fg:
+  A_sync_BB: 1.6
+  EB_sync: 2.0
+  alpha_sync_EE: -0.9
+  alpha_sync_BB: -0.9
+  beta_sync: -3.1
+  A_dust_BB: 28.0
+  EB_dust: 2.0
+  alpha_dust_EE: -0.2
+  alpha_dust_BB: -0.2
+  Alens: 1.
+  r: 0.
+  use_moments: False

bb4cast/params/noise_config.yml

@@ -0,0 +1,19 @@
+# Allowed values: baseline/goal
+sensitivity: baseline
+# Allowed values:
+#  - pessimistic/optimistic
+#  - custom: in this case, include also a transfer function (below)
+one_over_f: optimistic
+# If `one_over_f == 'custom'`, pass a transfer function below.
+# The format of this file should be two columns: ell, T_ell
+transfer: bb4cast/params/transfer_dummy.txt
+# Number of detector-years in each band
+# Values below are the nominal SO 5-year survey
+ydet_LF: 1.
+ydet_MF: 9.
+ydet_UHF: 5.
+recompute_data: False
+# How do you call python?
+pyexec: "python3"
+# Sampler to use (Fisher matrix)
+sampler: fisher

bb4cast/params/so_nominal.yml

@@ -0,0 +1,30 @@
+bands:
+  LF1:
+    freq: 27.
+    beam: 91.
+    bandpass: examples/data/bandpasses/LF1.txt
+  LF2:
+    freq: 39.
+    beam: 63.
+    bandpass: examples/data/bandpasses/LF2.txt
+  MF1:
+    freq: 93.
+    beam: 30.
+    bandpass: examples/data/bandpasses/MF1.txt
+  MF2:
+    freq: 145.
+    beam: 17.
+    bandpass: examples/data/bandpasses/MF2.txt
+  UHF1:
+    freq: 225.
+    beam: 11.
+    bandpass: examples/data/bandpasses/UHF1.txt
+  UHF2:
+    freq: 280.
+    beam: 9.
+    bandpass: examples/data/bandpasses/UHF2.txt
+
+f_sky: 0.184
+f_sky_ss: 0.114
+f_sky_sn: 0.152
+f_sky_nn: 0.184

bb4cast/params/so_nominal_beampert.yml

@@ -0,0 +1,36 @@
+bands:
+  LF1:
+    freq: 27.
+    beam: 91.
+    bandpass: examples/data/bandpasses/LF1.txt
+    delta_beam: bb4cast/example_data/delta_beam_LF1_FWHM0p01.txt
+  LF2:
+    freq: 39.
+    beam: 63.
+    bandpass: examples/data/bandpasses/LF2.txt
+    delta_beam: bb4cast/example_data/delta_beam_LF2_FWHM0p01.txt
+  MF1:
+    freq: 93.
+    beam: 30.
+    bandpass: examples/data/bandpasses/MF1.txt
+    delta_beam: bb4cast/example_data/delta_beam_MF1_FWHM0p01.txt
+  MF2:
+    freq: 145.
+    beam: 17.
+    bandpass: examples/data/bandpasses/MF2.txt
+    delta_beam: bb4cast/example_data/delta_beam_MF2_FWHM0p01.txt
+  UHF1:
+    freq: 225.
+    beam: 11.
+    bandpass: examples/data/bandpasses/UHF1.txt
+    delta_beam: bb4cast/example_data/delta_beam_UHF1_FWHM0p01.txt
+  UHF2:
+    freq: 280.
+    beam: 9.
+    bandpass: examples/data/bandpasses/UHF2.txt
+    delta_beam: bb4cast/example_data/delta_beam_UHF2_FWHM0p01.txt
+
+f_sky: 0.1
+f_sky_ss: 0.1
+f_sky_sn: 0.1
+f_sky_nn: 0.1

bb4cast/params/st_baseline_optimistic_1_9_5.yml

@@ -0,0 +1,11 @@
+prefix_out: so_nominal_baseline_optimistic
+sensitivity: baseline
+one_over_f: optimistic
+ydet_LF: 1.
+ydet_MF: 9.
+ydet_UHF: 5.
+recompute_data: False
+pyexec: "addqueue -s -c baseline_optimistic_1_9_5 -n 1x12 -q cmb -m 2 /usr/bin/python3"
+sampler: emcee
+nwalkers: 128
+nsamples: 10000