Commit c7ac3443 by Jakob Knollmueller

 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) plt.figure() # plt.plot(data, 'k-') f, (ax0, ax1,ax2) = plt.subplots(3, sharex=True, sharey=True) plt.suptitle('diffuse components') ax0.plot(exp(energy1.s).val, 'k-') ax0.yaxis.set_label_position("right") ax0.set_ylabel(r'$\alpha = 1.25$') ax0.set_ylim(1e-1,1e3) ax0.set_yscale("log") ax1.plot(exp(energy2.s).val, 'k-') ax1.yaxis.set_label_position("right") ax1.set_ylabel(r'$\alpha = 1.5$') ax2.plot(exp(energy3.s).val, 'k-') ax2.yaxis.set_label_position("right") ax2.set_ylabel(r'$\alpha = 1.75$') plt.savefig('1d_diffuse.pdf') plt.figure() f, (ax0, ax1,ax2) = plt.subplots(3, sharex=True, sharey=True) plt.suptitle('point-like components') ax0.plot(exp(energy1.u).val, 'k-') ax0.yaxis.set_label_position("right") ax0.set_ylabel(r'$\alpha = 1.25$') ax0.set_ylim(1e-1,1e3) ax0.set_yscale("log") ax1.plot(exp(energy2.u).val, 'k-') ax1.yaxis.set_label_position("right") ax1.set_ylabel(r'$\alpha = 1.5$') ax2.plot(exp(energy3.u).val, 'k-') ax2.yaxis.set_label_position("right") ax2.set_ylabel(r'$\alpha = 1.75$') 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) plt.suptitle('data and true components') ax0.plot(data, 'k-') ax0.set_yscale("log") ax0.set_ylim(1e-1,1e3) ax0.yaxis.set_label_position("right") ax0.set_ylabel(r'data') ax1.plot(exp(s).val, 'k-') ax1.yaxis.set_label_position("right") ax1.set_ylabel(r'diffuse') ax2.plot(exp(u).val, 'k-') ax2.yaxis.set_label_position("right") ax2.set_ylabel(r'point-like') # plt.ylim(1e-0) plt.savefig('1d_data.pdf')