diff --git a/Paper/GlitchMCNoiseOnly/plot_data.py b/Paper/GlitchMCNoiseOnly/plot_data.py index b04d689036920957b0b601dfdb08fb75473a5034..75a9a5b5153f7abc00179af9b5d50030791c8654 100644 --- a/Paper/GlitchMCNoiseOnly/plot_data.py +++ b/Paper/GlitchMCNoiseOnly/plot_data.py @@ -7,8 +7,35 @@ from tqdm import tqdm from oct2py import octave import glob from scipy.stats import rv_continuous, chi2 +from scipy.special import gammaincc from latex_macro_generator import write_to_macro + +def CF_twoFmax_integrand(theta, twoFmax, Nt): + Fmax = twoFmax/2.0 + return np.exp(1j*theta*twoFmax)*Nt/2.0*Fmax*np.exp(-Fmax)*(1-(1+Fmax)*np.exp(-Fmax))**(Nt-1.) + + +def pdf_twoFhat(twoFhat, Ng, Nts, twoFtildemax=100, dtwoF=0.05): + if np.ndim(Nts) == 0: + Nts = np.zeros(Ng+1) + Nts + + twoFtildemax_int = np.arange(0, twoFtildemax, dtwoF) + theta_int = np.arange(-1./dtwoF, 1./dtwoF, 1./twoFtildemax) + + CF_twoFtildemax_theta = np.array( + [[np.trapz(CF_twoFmax_integrand(t, twoFtildemax_int, Nt), twoFtildemax_int) + for t in theta_int] + for Nt in Nts]) + + CF_twoFhat_theta = np.prod(CF_twoFtildemax_theta, axis=0) + print CF_twoFhat_theta.shape, theta_int.shape + pdf = (1/(2*np.pi)) * np.array( + [np.trapz(np.exp(-1j*theta_int*twoFhat_val)*CF_twoFhat_theta, + theta_int) for twoFhat_val in twoFhat]) + print np.trapz(pdf.real, x=twoFhat) + return pdf + filenames = glob.glob("CollectedOutput/*.txt") plt.style.use('paper') @@ -23,25 +50,41 @@ for fn in filenames: df_list.append(df) df = pd.concat(df_list) +colors = ['C0', 'C1'] fig, ax = plt.subplots() -for ng in np.unique(df.nglitches.values): +handles = [] +labels = [] +for ng, c in zip(np.unique(df.nglitches.values), colors): print 'ng={}'.format(ng) - Nsamples = len(df[df.nglitches == ng]) - MaxtwoF = df[df.nglitches == ng].twoF.max() + df_temp = df[df.nglitches == ng] + #df_temp = df_temp[[str(x).isalpha() for x in df_temp.CLUSTER_ID.values]] + print df_temp.tail() + Nsamples = len(df_temp) + MaxtwoF = df_temp.twoF.max() print 'Number of samples = ', Nsamples print 'Max twoF', MaxtwoF - - ax.hist(df[df.nglitches == ng].twoF, bins=40, histtype='step', normed=True, - linewidth=1, label='$N_\mathrm{{glitches}}={}$'.format(ng)) + print np.any(np.isnan(df_temp.twoF.values)) + ax.hist(df_temp.twoF, bins=40, histtype='stepfilled', + normed=True, align='mid', alpha=0.5, + linewidth=1, label='$N_\mathrm{{glitches}}={}$'.format(ng), + color=c) write_to_macro('DirectedMC{}GlitchNoiseOnlyMaximum'.format(ng), '{:1.1f}'.format(MaxtwoF), '../macros.tex') write_to_macro('DirectedMC{}GlitchNoiseN'.format(ng), '{:1.0f}'.format(Nsamples), '../macros.tex') + twoFmax = np.linspace(0, 100, 200) + ax.plot(twoFmax, pdf_twoFhat(twoFmax, ng, Nsamples, + twoFtildemax=2*MaxtwoF, dtwoF=0.1), + color=c, label='$N_\mathrm{{glitches}}={}$ predicted'.format(ng)) + ax.set_xlabel('$\widehat{2\mathcal{F}}$') ax.set_xlim(0, 90) -ax.legend(frameon=False, fontsize=6) +handles, labels = ax.get_legend_handles_labels() +idxs = np.argsort(labels) +ax.legend(np.array(handles)[idxs], np.array(labels)[idxs], frameon=False, + fontsize=6) fig.tight_layout() fig.savefig('glitch_noise_twoF_histogram.png') diff --git a/Paper/GlitchMCNoiseOnly/submitfile b/Paper/GlitchMCNoiseOnly/submitfile index cb3ef0770ed09e7abe626e26200584c76ab722a4..99f6657352f9552bb5415ebf4ee4764436aed140 100644 --- a/Paper/GlitchMCNoiseOnly/submitfile +++ b/Paper/GlitchMCNoiseOnly/submitfile @@ -9,4 +9,4 @@ Log=CollectedOutput/log.$(Cluster).$(Process) request_cpus = 1 request_memory = 16 GB -Queue 170 +Queue 1806 diff --git a/Paper/glitch_noise_twoF_histogram.png b/Paper/glitch_noise_twoF_histogram.png index 3068a22984eecb41fcfccc6d968bfe98801192a3..16ea1846e98658bf2f0330cf1031edd490a8b2ed 100644 Binary files a/Paper/glitch_noise_twoF_histogram.png and b/Paper/glitch_noise_twoF_histogram.png differ diff --git a/Paper/macros.tex b/Paper/macros.tex index b4477a9b8e5005871cc220e4dde664eba131bd45..3a2ee919d769b6da783584122fff66b0b3071805 100644 --- a/Paper/macros.tex +++ b/Paper/macros.tex @@ -17,8 +17,8 @@ \def\DirectedMCNoiseOnlyMaximum{52.4} \def\DirectedMConeGlitchNoiseN{10000} \def\DirectedMConeGlitchNoiseOnlyMaximum{82.6} -\def\DirectedMCtwoGlitchNoiseN{970} -\def\DirectedMCtwoGlitchNoiseOnlyMaximum{80.8} +\def\DirectedMCtwoGlitchNoiseN{9625} +\def\DirectedMCtwoGlitchNoiseOnlyMaximum{87.8} \def\SingleGlitchDepth{10.0} \def\SingleGlitchFCMismatch{0.7} \def\SingleGlitchFCtwoF{718.5}