From b9314a28008eee2646f3c18ee104285fd5acebc7 Mon Sep 17 00:00:00 2001 From: Eric Teunis de Boone Date: Tue, 14 Nov 2023 16:25:49 +0100 Subject: [PATCH] PDFs: revamped as used in Thesis --- simulations/12_noise_phase.py | 63 ++++++++++++++++++++++++----------- 1 file changed, 43 insertions(+), 20 deletions(-) diff --git a/simulations/12_noise_phase.py b/simulations/12_noise_phase.py index 0d2e709..389be0a 100644 --- a/simulations/12_noise_phase.py +++ b/simulations/12_noise_phase.py @@ -13,6 +13,7 @@ import scipy.stats as stat from scipy import special from lib.util import MethodMappingProxy as MethodProxy +from lib.ft_plot import axis_pi_ticker def expectation(x,pdfx): dx = x[1]-x[0] @@ -48,6 +49,23 @@ def amplitude_distribution_gauss(a,sigma,s): k=s/sigma return (2*np.pi)**-0.5*np.exp(-((a-s)/sigma)**2/2) + +figsize = (8,6) +if True: + from matplotlib import rcParams + #rcParams["text.usetex"] = True + rcParams["font.family"] = "serif" + plt.rc('lines',lw=2) + if True:# small + figsize = (6, 4) + rcParams["font.size"] = "15" # 15 at 6,4 looks fine + elif True: # large + figsize = (9, 6) + rcParams["font.size"] = "16" # 15 at 9,6 looks fine + rcParams["grid.linestyle"] = 'dotted' + rcParams["figure.figsize"] = figsize + + signal_max= 4 amp_max = signal_max*2 thetas = np.linspace(-np.pi,np.pi,500) @@ -58,14 +76,14 @@ sigma = 1 ## figure 1 if True: - fig, ax = plt.subplots(1,2,figsize=(2*8,1*8)) - _fig1, _ax0 =plt.subplots(1,1, figsize=(1*8, 1*8)) - _fig2, _ax1 =plt.subplots(1,1, figsize=(1*8, 1*8)) + fig, ax = plt.subplots(1,2,figsize=(2*figsize[0],figsize[1])) + _fig1, _ax0 =plt.subplots(1,1) + _fig2, _ax1 =plt.subplots(1,1) ax0 = MethodProxy(ax[0], _ax0) ax1 = MethodProxy(ax[1], _ax1) for s in signals: - pdfs_label='s/$\sigma$ ='+str(s) + pdfs_label='s = '+str(int(s)) phase_vals= phase_distribution(thetas,sigma,s) amp_vals= amplitude_distribution(amplitudes,sigma,s) phase_vals_g = phase_distribution_gauss(thetas,sigma,s) @@ -73,13 +91,16 @@ if True: ax0.plot(amplitudes,amp_vals, label=pdfs_label) ax1.plot(thetas,phase_vals, label=pdfs_label) ax0.legend() + ax0.set_xlabel(r'Amplitude $a$') + ax0.set_ylabel(r'$p(a)$') + + ax1.set_xlabel(r'Phase $\varphi$') + ax1.set_ylabel(r'$p(\varphi)$') _ax1.legend()# only in the separate figure - ax0.set_xlabel(r'$a$') - ax0.set_xlabel(r'$\theta$') - ax1.set_ylabel(r'$p(a)$') - ax1.set_ylabel(r'$p(\theta)$') - # ax[0].grid() - # ax[1].grid() + [ axis_pi_ticker(ax.xaxis, major_divider=3) for ax in ax1.elements ] + + for a in [ax0, ax1]: + a.grid() MethodProxy(fig, _fig1, _fig2).tight_layout() fig.savefig('pdfs.pdf') @@ -92,7 +113,7 @@ if True: ## figure 2 amplitudes = np.linspace(0,amp_max*5,500) signals = np.linspace(0.1,signal_max*5,101) -if False: +if True: V_theta = [variance(thetas,phase_distribution(thetas,sigma,s)) for s in signals ] E_theta=[expectation(thetas,phase_distribution(thetas,sigma,s)) for s in signals ] V_theta_g = [variance(thetas,phase_distribution_gauss(thetas,sigma,s)) for s in signals ] @@ -102,17 +123,17 @@ if False: V_a_g = [variance(amplitudes,amplitude_distribution_gauss(amplitudes,sigma,s)) for s in signals ] E_a_g=[expectation(amplitudes,amplitude_distribution_gauss(amplitudes,sigma,s)) for s in signals ] - fig2, _ax2 = plt.subplots(2,2,figsize=(2*8,2*8)) + fig2, _ax2 = plt.subplots(2,2,figsize=(2*figsize[0],2*figsize[1])) ax2 = fig2.get_axes() if True: _figs = [] _axs = [] - for i, ax in enumerate(_ax2): - _f, _a = plt.subplots(1,1, figsize=(1*8, 1*8)) + for i, ax in enumerate(ax2): + _f, _a = plt.subplots(1,1) _figs.append(_f) _axs.append(_a) - ax2[i] = MethodProxy(ax2[0], _a) + ax2[i] = MethodProxy(ax, _a) ax2[0].plot(signals,E_a,label='$p(a)$') ax2[0].plot(signals,E_a_g,ls='dashed',label='Gaussian approx.') @@ -125,17 +146,17 @@ if False: ax2[1].set_xscale('log') ax2[1].set_ylabel('$\sigma_a^2$') - ax2[2].plot(signals,E_theta,label=r'$p(\theta)$') + ax2[2].plot(signals,E_theta,label=r'$p(\varphi)$') ax2[2].plot(signals,E_theta_g,ls='dashed',label='Gaussian approx.') ax2[2].set_xscale('log') ax2[2].set_ylim(-1.1,1.1) - ax2[2].set_ylabel(r'$\mu_\theta$') + ax2[2].set_ylabel(r'$\mu_\varphi$') - ax2[3].plot(signals,V_theta,label=r'$p(\theta)$') + ax2[3].plot(signals,V_theta,label=r'$p(\varphi)$') ax2[3].plot(signals,V_theta_g,ls='dashed',label='Gaussian approx.') ax2[3].set_xscale('log') ax2[3].set_yscale('log') - ax2[3].set_ylabel(r'$\sigma_\theta^2$') + ax2[3].set_ylabel(r'$\sigma_\varphi^2$') for a in ax2: a.grid(which='both') a.set_xlabel(r'$s/\sigma$') @@ -150,12 +171,14 @@ if False: 'phase_mean', 'phase_sigma', ][i] + + _f.tight_layout() _f.savefig(fnames+'.pdf') plt.close(_f) ## figure 3, beacon timing accuracy if True: - fig3, ax3 = plt.subplots(1,1,figsize=(1*8,1*8)) + fig3, ax3 = plt.subplots(1,1) ax3 = fig3.get_axes() sigma_t = [variance(thetas,phase_distribution(thetas,sigma,s)) for s in signals ] lfs=np.linspace(np.log10(50.),4,1)