diff --git a/Paper/AllSkyMC/generate_table.py b/Paper/AllSkyMC/generate_table.py index 9fa71a16cc184996ade28e8430bd0473eb5e3118..0cf34001a0895a36aec32785a6ef51cfcf0611a5 100644 --- a/Paper/AllSkyMC/generate_table.py +++ b/Paper/AllSkyMC/generate_table.py @@ -82,7 +82,7 @@ mcmc = pyfstat.MCMCFollowUpSearch( sftfilepath='{}/*{}*sft'.format(outdir, data_label), theta_prior=theta_prior, tref=tref, minStartTime=tstart, maxStartTime=tend, - nwalkers=nwalkers, ntemps=ntemps, + nwalkers=nwalkers, ntemps=ntemps, nsteps=[nsteps, nsteps] log10temperature_min=log10temperature_min) -mcmc.run(Nsegs0=20, R=10) -#mcmc.run(run_setup) +#mcmc.run(Nsegs0=20, R=10) +mcmc.run(run_setup) diff --git a/Paper/AllSkyMCNoiseOnly/plot_data.py b/Paper/AllSkyMCNoiseOnly/plot_data.py index 2f3c878865338a516ca550b48c4c4a6399821e74..a0a3b1200579014fff982f2015e92577cb0d0225 100644 --- a/Paper/AllSkyMCNoiseOnly/plot_data.py +++ b/Paper/AllSkyMCNoiseOnly/plot_data.py @@ -5,6 +5,7 @@ import os from tqdm import tqdm from oct2py import octave import glob +from scipy.stats import rv_continuous, chi2 filenames = glob.glob("CollectedOutput/*.txt") @@ -13,12 +14,13 @@ plt.style.use('paper') Tspan = 100 * 86400 -def maxtwoFinNoise(twoF, Ntrials): - F = twoF/2.0 - alpha = (1 + F)*np.exp(-F) - a = Ntrials/2.0*F*np.exp(-F) - b = (1 - alpha)**(Ntrials-1) - return a*b +class maxtwoFinNoise_gen(rv_continuous): + def _pdf(self, twoF, Ntrials): + F = twoF/2.0 + alpha = (1 + F)*np.exp(-F) + a = Ntrials/2.0*F*np.exp(-F) + b = (1 - alpha)**(Ntrials-1) + return a*b df_list = [] @@ -31,6 +33,21 @@ df = pd.concat(df_list) print 'Number of samples = ', len(df) fig, ax = plt.subplots() + +maxtwoFinNoise = maxtwoFinNoise_gen(a=0) +Ntrials_effective, loc, scale = maxtwoFinNoise.fit(df.twoF.values, floc=0, fscale=1) +print 'Ntrials effective = {:1.2e}'.format(Ntrials_effective) +twoFsmooth = np.linspace(0, df.twoF.max(), 1000) +best_fit_pdf = maxtwoFinNoise.pdf(twoFsmooth, Ntrials_effective) +ax.plot(twoFsmooth, best_fit_pdf, '-r') + +pval = 1e-6 +twoFsmooth_HD = np.linspace( + twoFsmooth[np.argmax(best_fit_pdf)], df.twoF.max(), 100000) +best_fit_pdf_HD = maxtwoFinNoise.pdf(twoFsmooth_HD, Ntrials_effective) +spacing = twoFsmooth_HD[1]-twoFsmooth_HD[0] +print twoFsmooth_HD[np.argmin(np.abs(best_fit_pdf_HD - pval))], spacing + ax.hist(df.twoF, bins=50, histtype='step', color='k', normed=True, linewidth=1) twoFsmooth = np.linspace(0, df.twoF.max(), 100) # ax.plot(twoFsmooth, maxtwoFinNoise(twoFsmooth, 8e5), '-r') diff --git a/Paper/DirectedMC/generate_table.py b/Paper/DirectedMC/generate_table.py index 3877cf5c14401f297f58f08319d8fe2b5234bf75..63cf7015739d0962f4017c3942204d567ea2c032 100644 --- a/Paper/DirectedMC/generate_table.py +++ b/Paper/DirectedMC/generate_table.py @@ -27,7 +27,7 @@ DeltaF1 = VF1 * np.sqrt(45/4.)/(np.pi*Tspan**2) depth = 100 -nsteps = 20 +nsteps = 25 run_setup = [((nsteps, 0), 20, False), ((nsteps, 0), 7, False), ((nsteps, 0), 2, False), @@ -71,6 +71,6 @@ mcmc = pyfstat.MCMCFollowUpSearch( sftfilepath='{}/*{}*sft'.format(outdir, data_label), theta_prior=theta_prior, tref=tref, minStartTime=tstart, maxStartTime=tend, - nwalkers=nwalkers, ntemps=ntemps, + nwalkers=nwalkers, ntemps=ntemps, nsteps=[nsteps, nsteps], log10temperature_min=log10temperature_min) mcmc.run(run_setup) diff --git a/Paper/DirectedMCNoiseOnly/plot_data.py b/Paper/DirectedMCNoiseOnly/plot_data.py index 43ba2c41361190599f16e61380c571b084ed34e3..01b954f1a8623b7bd9350caecee6cb4cf3baf12f 100644 --- a/Paper/DirectedMCNoiseOnly/plot_data.py +++ b/Paper/DirectedMCNoiseOnly/plot_data.py @@ -5,6 +5,7 @@ import os from tqdm import tqdm from oct2py import octave import glob +from scipy.stats import rv_continuous, chi2 filenames = glob.glob("CollectedOutput/*.txt") @@ -13,12 +14,13 @@ plt.style.use('paper') Tspan = 100 * 86400 -def maxtwoFinNoise(twoF, Ntrials): - F = twoF/2.0 - alpha = (1 + F)*np.exp(-F) - a = Ntrials/2.0*F*np.exp(-F) - b = (1 - alpha)**(Ntrials-1) - return a*b +class maxtwoFinNoise_gen(rv_continuous): + def _pdf(self, twoF, Ntrials): + F = twoF/2.0 + alpha = (1 + F)*np.exp(-F) + a = Ntrials/2.0*F*np.exp(-F) + b = (1 - alpha)**(Ntrials-1) + return a*b df_list = [] @@ -32,8 +34,21 @@ print 'Number of samples = ', len(df) fig, ax = plt.subplots() ax.hist(df.twoF, bins=50, histtype='step', color='k', normed=True, linewidth=1) -twoFsmooth = np.linspace(0, df.twoF.max(), 100) -# ax.plot(twoFsmooth, maxtwoFinNoise(twoFsmooth, 2e3), '-r') + +maxtwoFinNoise = maxtwoFinNoise_gen(a=0) +Ntrials_effective, loc, scale = maxtwoFinNoise.fit(df.twoF.values, floc=0, fscale=1) +print 'Ntrials effective = {:1.2e}'.format(Ntrials_effective) +twoFsmooth = np.linspace(0, df.twoF.max(), 1000) +best_fit_pdf = maxtwoFinNoise.pdf(twoFsmooth, Ntrials_effective) +ax.plot(twoFsmooth, best_fit_pdf, '-r') + +pval = 1e-6 +twoFsmooth_HD = np.linspace( + twoFsmooth[np.argmax(best_fit_pdf)], df.twoF.max(), 100000) +best_fit_pdf_HD = maxtwoFinNoise.pdf(twoFsmooth_HD, Ntrials_effective) +spacing = twoFsmooth_HD[1]-twoFsmooth_HD[0] +print twoFsmooth_HD[np.argmin(np.abs(best_fit_pdf_HD - pval))], spacing + ax.set_xlabel('$\widetilde{2\mathcal{F}}$') ax.set_xlim(0, 60) fig.tight_layout() diff --git a/Paper/allsky_noise_twoF_histogram.png b/Paper/allsky_noise_twoF_histogram.png index 9623605b5845b4664f726cfdd5a397ff18b21218..9a32b952d881d57e2a416065164586a2b5469e84 100644 Binary files a/Paper/allsky_noise_twoF_histogram.png and b/Paper/allsky_noise_twoF_histogram.png differ diff --git a/Paper/allsky_setup_run_setup.tex b/Paper/allsky_setup_run_setup.tex index d98ab7f945e231e2516c301f4356b95a29502856..8c4d90d186e3b3b8739ec6c66e555c0b7b58ae91 100644 --- a/Paper/allsky_setup_run_setup.tex +++ b/Paper/allsky_setup_run_setup.tex @@ -1,8 +1,8 @@ \begin{tabular}{c|cccccc} Stage & $\Nseg$ & $\Tcoh^{\rm days}$ &$\Nsteps$ & $\V$ & $\Vsky$ & $\Vpe$ \\ \hline -0 & 20 & 5.0 & 100 & $4{\times}10^{2}$ & 8.0 & 60.0 \\ -1 & 11 & 9.1 & 100 & $4{\times}10^{3}$ & 20.0 & $2{\times}10^{2}$ \\ -2 & 6 & 16.7 & 100 & $4{\times}10^{4}$ & 70.0 & $6{\times}10^{2}$ \\ -3 & 3 & 33.3 & 100 & $3{\times}10^{5}$ & $1{\times}10^{2}$ & $2{\times}10^{3}$ \\ -4 & 1 & 100.0 & 100,100 & $2{\times}10^{6}$ & $2{\times}10^{2}$ & $1{\times}10^{4}$ \\ +0 & 20 & 5.0 & 50 & $4{\times}10^{2}$ & 8.0 & 60.0 \\ +1 & 11 & 9.1 & 50 & $4{\times}10^{3}$ & 20.0 & $2{\times}10^{2}$ \\ +2 & 6 & 16.7 & 50 & $4{\times}10^{4}$ & 70.0 & $6{\times}10^{2}$ \\ +3 & 3 & 33.3 & 50 & $3{\times}10^{5}$ & $1{\times}10^{2}$ & $2{\times}10^{3}$ \\ +4 & 1 & 100.0 & 50,50 & $2{\times}10^{6}$ & $2{\times}10^{2}$ & $1{\times}10^{4}$ \\ \end{tabular} diff --git a/Paper/directed_noise_twoF_histogram.png b/Paper/directed_noise_twoF_histogram.png index 9dfd21c9f0af3f517363d7f6333470d7effab6d3..9d8db4a694ecac2e3f1321a993a36203bdf2295b 100644 Binary files a/Paper/directed_noise_twoF_histogram.png and b/Paper/directed_noise_twoF_histogram.png differ diff --git a/Paper/directed_setup_run_setup.tex b/Paper/directed_setup_run_setup.tex index bf81a4a6451a8a16f1de11fb11edb47edf81b0f2..dbb27e8eb2d1a9310228ff4d720e32734115ee87 100644 --- a/Paper/directed_setup_run_setup.tex +++ b/Paper/directed_setup_run_setup.tex @@ -1,7 +1,7 @@ \begin{tabular}{c|cccc} Stage & $\Nseg$ & $\Tcoh^{\rm days}$ &$\Nsteps$ & $\Vpe$ \\ \hline -0 & 20 & 5.0 & 20 & 10.0 \\ -1 & 7 & 14.3 & 20 & $1{\times}10^{2}$ \\ -2 & 2 & 50.0 & 20 & $1{\times}10^{3}$ \\ -3 & 1 & 100.0 & 20,20 & $2{\times}10^{3}$ \\ +0 & 20 & 5.0 & 25 & 10.0 \\ +1 & 7 & 14.3 & 25 & $1{\times}10^{2}$ \\ +2 & 2 & 50.0 & 25 & $1{\times}10^{3}$ \\ +3 & 1 & 100.0 & 25,25 & $2{\times}10^{3}$ \\ \end{tabular} diff --git a/Paper/macros.tex b/Paper/macros.tex index db1e556bd3b91252e54fd2d3b33b00fad51ee214..df6253acb74f0ed5c54e1bb231391e61af42f057 100644 --- a/Paper/macros.tex +++ b/Paper/macros.tex @@ -1,4 +1,4 @@ -\def\AllSkyMCNoiseN{10000} +\def\AllSkyMCNoiseN{1.0{\times}10^{4}} \def\AllSkyMCNoiseOnlyMaximum{59.8} \def\BasicExampleDeltaFone{1.0{\times}10^{-13}} \def\BasicExampleDeltaFzero{1.0{\times}10^{-7}} @@ -12,5 +12,5 @@ \def\BasicExamplehzero{1.0{\times}10^{-24}} \def\BasicExamplenburn{50.0} \def\BasicExamplenprod{50.0} -\def\DirectedMCNoiseN{10000} +\def\DirectedMCNoiseN{1.0{\times}10^{4}} \def\DirectedMCNoiseOnlyMaximum{52.4}