1d_separation.py 2.94 KB
 Jakob Knollmueller committed Feb 14, 2018 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 from point_separation import build_problem, problem_iteration from nifty2go import * import numpy as np from matplotlib import rc rc('font',**{'family':'serif','serif':['Palatino']}) rc('text', usetex=True) from matplotlib import pyplot as plt np.random.seed(42) if __name__ == '__main__': s_space = RGSpace([1024]) FFT = FFTOperator(s_space) h_space = FFT.target[0] p_space = PowerSpace(h_space) sp = Field(p_space, val=1./(1+p_space.k_lengths)**2.5 ) sh = power_synthesize(sp) s = FFT.adjoint_times(sh) u = Field(s_space, val = -12) u.val[200] = 1 u.val[300] = 3 u.val[500] = 4 u.val[700] = 5 u.val[900] = 2 u.val[154] = 0.5 u.val[421] = 0.25 u.val[652] = 1 u.val[1002] = 2.5 d = exp(s) + exp(u) data = d.val energy1 = build_problem(data,1.25) energy2 = build_problem(data,1.5) energy3 = build_problem(data,1.75) for i in range(20): energy1 = problem_iteration(energy1) energy2 = problem_iteration(energy2) energy3 = problem_iteration(energy3)  Jakob Knollmueller committed Feb 15, 2018 42  size = 15  Jakob Knollmueller committed Feb 14, 2018 43 44 45  plt.figure() # plt.plot(data, 'k-') f, (ax0, ax1,ax2) = plt.subplots(3, sharex=True, sharey=True)  Jakob Knollmueller committed Feb 15, 2018 46  plt.suptitle('diffuse components', size=size)  Jakob Knollmueller committed Feb 14, 2018 47 48 49  ax0.plot(exp(energy1.s).val, 'k-') ax0.yaxis.set_label_position("right")  Jakob Knollmueller committed Feb 15, 2018 50  ax0.set_ylabel(r'$\alpha = 1.25$', size=size)  Jakob Knollmueller committed Feb 14, 2018 51 52 53 54 55  ax0.set_ylim(1e-1,1e3) ax0.set_yscale("log") ax1.plot(exp(energy2.s).val, 'k-') ax1.yaxis.set_label_position("right")  Jakob Knollmueller committed Feb 15, 2018 56  ax1.set_ylabel(r'$\alpha = 1.5$', size=size)  Jakob Knollmueller committed Feb 14, 2018 57 58 59  ax2.plot(exp(energy3.s).val, 'k-') ax2.yaxis.set_label_position("right")  Jakob Knollmueller committed Feb 15, 2018 60  ax2.set_ylabel(r'$\alpha = 1.75$', size=size)  Jakob Knollmueller committed Feb 14, 2018 61 62 63 64 65 66  plt.savefig('1d_diffuse.pdf') plt.figure() f, (ax0, ax1,ax2) = plt.subplots(3, sharex=True, sharey=True)  Jakob Knollmueller committed Feb 15, 2018 67  plt.suptitle('point-like components', size=size)  Jakob Knollmueller committed Feb 14, 2018 68 69 70  ax0.plot(exp(energy1.u).val, 'k-') ax0.yaxis.set_label_position("right")  Jakob Knollmueller committed Feb 15, 2018 71  ax0.set_ylabel(r'$\alpha = 1.25$', size=size)  Jakob Knollmueller committed Feb 14, 2018 72 73 74 75 76  ax0.set_ylim(1e-1,1e3) ax0.set_yscale("log") ax1.plot(exp(energy2.u).val, 'k-') ax1.yaxis.set_label_position("right")  Jakob Knollmueller committed Feb 15, 2018 77  ax1.set_ylabel(r'$\alpha = 1.5$', size=size)  Jakob Knollmueller committed Feb 14, 2018 78 79 80  ax2.plot(exp(energy3.u).val, 'k-') ax2.yaxis.set_label_position("right")  Jakob Knollmueller committed Feb 15, 2018 81  ax2.set_ylabel(r'$\alpha = 1.75$', size=size)  Jakob Knollmueller committed Feb 14, 2018 82 83 84 85 86 87 88 89 90  ax0.set_yscale("log") ax0.set_ylim(1e-1,1e3) # plt.ylim(1e-0) plt.savefig('1d_points.pdf') plt.figure() f, (ax0, ax1,ax2) = plt.subplots(3, sharex=True, sharey=True)  Jakob Knollmueller committed Feb 15, 2018 91  plt.suptitle('data and true components', size=size)  Jakob Knollmueller committed Feb 14, 2018 92 93 94 95 96  ax0.plot(data, 'k-') ax0.set_yscale("log") ax0.set_ylim(1e-1,1e3) ax0.yaxis.set_label_position("right")  Jakob Knollmueller committed Feb 15, 2018 97  ax0.set_ylabel(r'data', size=size)  Jakob Knollmueller committed Feb 14, 2018 98 99 100 101  ax1.plot(exp(s).val, 'k-') ax1.yaxis.set_label_position("right")  Jakob Knollmueller committed Feb 15, 2018 102  ax1.set_ylabel(r'diffuse', size=size)  Jakob Knollmueller committed Feb 14, 2018 103 104  ax2.plot(exp(u).val, 'k-') ax2.yaxis.set_label_position("right")  Jakob Knollmueller committed Feb 15, 2018 105  ax2.set_ylabel(r'point-like', size=size)  Jakob Knollmueller committed Feb 14, 2018 106 107 108 109  # plt.ylim(1e-0) plt.savefig('1d_data.pdf')