SHARE
TWEET

Untitled

a guest Sep 20th, 2019 87 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. theo_unlen = cmblens.theory.load_theory_cls('cosmo2017_10K_acc3_unlensed_cmb', unit=unit)
  2. L    = theo_unlen['l']
  3. clkk = theo_unlen['clkk']
  4. clpp = theo_unlen['clpp']
  5.  
  6. def get_delensed_amanda(tt_unlen, te_unlen, ee_unlen, lpp_theo, clpp_theo, L, N_L, klmin, klmax, load_dls=True, unit='uk', lmax_ret=8000):
  7.     import camb.correlations
  8.     import numpy as np
  9.  
  10.     kappa_lmin = int(klmin)
  11.     kappa_lmax = int(klmax)
  12.  
  13.     reload(np)  # necessary to call this function multiple times (?)
  14.  
  15.     # want theory to start at l = 0
  16.     lpp_th = np.zeros(len(lpp_theo)+int(lpp_theo[0]))
  17.     lpp_th[0] = 0.
  18.     lpp_th[1] = 1.
  19.     lpp_th[int(lpp_theo[0]):] = np.array(lpp_theo)
  20.     # fill in unlensed theory and lensing potential
  21.     clpp_th = np.zeros(len(lpp_th))
  22.     tt_th = np.zeros(len(lpp_th))
  23.     te_th = np.zeros(len(lpp_th))
  24.     ee_th = np.zeros(len(lpp_th))
  25.     idx = int(lpp_theo[0])
  26.     clpp_th[idx:] = np.array(clpp_theo)
  27.     tt_th[idx:] = np.array(tt_unlen)
  28.     te_th[idx:] = np.array(te_unlen)
  29.     ee_th[idx:] = np.array(ee_unlen)
  30.     # fill in the noise curve so it also starts at L=0 and is same size as theory
  31.     NL = np.zeros(lpp_th.shape)
  32.     idx = int(L[0])
  33.     NL[:idx] = np.inf
  34.     if len(N_L) >= len(NL) + idx:
  35.         NL[idx:] = np.array(N_L[:len(NL)+idx])
  36.     else:
  37.         NL[idx:len(N_L)+idx] = np.array(N_L)
  38.         NL[len(N_L)+idx:] = np.inf     
  39.     L = lpp_th.copy()  
  40.  
  41.     def _weiner1d_filter(signal, noise):
  42.         return (signal/(signal+noise))*signal
  43.  
  44.     # get residual lensing potential
  45.     clpp_filt      = _weiner1d_filter(clpp_th, NL)
  46.     loc            = np.where(L > kappa_lmax)
  47.     if len(loc[0]) > 0:
  48.         if loc[0][0] > len(clpp_filt):
  49.             clpp_filt[loc[0][0]:] = 0.
  50.     loc            = np.where(L < kappa_lmin)
  51.     clpp_filt[loc] = 0.
  52.  
  53.     clpp_res = clpp_th-clpp_filt
  54.  
  55.     # create array to hold delensed spectra
  56.     dls = np.zeros((len(lpp_th),4))
  57.     dls[:,0] = tt_th
  58.     dls[:,1] = ee_th
  59.     dls[:,2] = np.zeros(len(tt_th))
  60.     dls[:,3] = te_th
  61.  
  62.     reload(camb)
  63.     reload(camb.correlations)
  64.     camb_ret = camb.correlations.lensed_cls(dls, clpp_res, delta_cls=False, sampling_factor=2.)
  65.  
  66.     ret         = {}
  67.     ret['dltt'] = camb_ret[:,0]
  68.     ret['dlee'] = camb_ret[:,1]
  69.     ret['dlbb'] = camb_ret[:,2]
  70.     ret['dlte'] = camb_ret[:,3]
  71.    
  72.     return ret['dltt'][1:], ret['dlte'][1:], ret['dlee'][1:]
  73.  
  74. L_n, cln0_mv    = cln0['l'], cln0[PCOMB.MV].copy()
  75. amanda_dl = get_delensed_amanda(theo_unlen['dltt'],theo_unlen['dlte'],theo_unlen['dlee'],L, clpp,L, cln0_mv, kappa_lmin, kappa_lmax)
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top