# Untitled

a guest
Sep 20th, 2019
98
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
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.