Advertisement
Guest User

Untitled

a guest
Sep 20th, 2019
121
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 2.28 KB | None | 0 0
  1. 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):
  2. import camb.correlations
  3. import numpy as np
  4.  
  5. kappa_lmin = int(klmin)
  6. kappa_lmax = int(klmax)
  7.  
  8. reload(np) # necessary to call this function multiple times (?)
  9.  
  10. # want theory to start at l = 0
  11. lpp_th = np.zeros(len(lpp_theo)+int(lpp_theo[0]))
  12. lpp_th[0] = 0.
  13. lpp_th[1] = 1.
  14. lpp_th[int(lpp_theo[0]):] = np.array(lpp_theo)
  15. # fill in unlensed theory and lensing potential
  16. clpp_th = np.zeros(len(lpp_th))
  17. tt_th = np.zeros(len(lpp_th))
  18. te_th = np.zeros(len(lpp_th))
  19. ee_th = np.zeros(len(lpp_th))
  20. idx = int(lpp_theo[0])
  21. clpp_th[idx:] = np.array(clpp_theo)
  22. tt_th[idx:] = np.array(tt_unlen)
  23. te_th[idx:] = np.array(te_unlen)
  24. ee_th[idx:] = np.array(ee_unlen)
  25. # fill in the noise curve so it also starts at L=0 and is same size as theory
  26. NL = np.zeros(lpp_th.shape)
  27. idx = int(L[0])
  28. NL[:idx] = np.inf
  29. if len(N_L) >= len(NL) + idx:
  30. NL[idx:] = np.array(N_L[:len(NL)+idx])
  31. else:
  32. NL[idx:len(N_L)+idx] = np.array(N_L)
  33. NL[len(N_L)+idx:] = np.inf
  34. L = lpp_th.copy()
  35.  
  36. def _weiner1d_filter(signal, noise):
  37. return (signal/(signal+noise))*signal
  38.  
  39. # get residual lensing potential
  40. clkk = clpp_th*(2*np.pi/4)
  41. clkk_filt = _weiner1d_filter(clkk, NL)
  42. clpp_filt = clkk_filt*4/(2*np.pi)
  43. loc = np.where(L > kappa_lmax)
  44. if len(loc[0]) > 0:
  45. if loc[0][0] > len(clpp_filt):
  46. clpp_filt[loc[0][0]:] = 0.
  47. loc = np.where(L < kappa_lmin)
  48. clpp_filt[loc] = 0.
  49.  
  50. clpp_res = clpp_th-clpp_filt
  51.  
  52. # create array to hold delensed spectra
  53. dls = np.zeros((len(lpp_th),4))
  54. dls[:,0] = tt_th
  55. dls[:,1] = ee_th
  56. dls[:,2] = np.zeros(len(tt_th))
  57. dls[:,3] = te_th
  58.  
  59. reload(camb)
  60. reload(camb.correlations)
  61. camb_ret = camb.correlations.lensed_cls(dls, clpp_res, delta_cls=False, sampling_factor=2.)
  62.  
  63. ret = {}
  64. ret['dltt'] = camb_ret[:,0]
  65. ret['dlee'] = camb_ret[:,1]
  66. ret['dlbb'] = camb_ret[:,2]
  67. ret['dlte'] = camb_ret[:,3]
  68.  
  69. return ret['dltt'][1:], ret['dlte'][1:], ret['dlee'][1:]
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement