diff --git a/Paper/AllSkyMCNoiseOnly/allsky_noise_twoF_histogram.png b/Paper/AllSkyMCNoiseOnly/allsky_noise_twoF_histogram.png
new file mode 100644
index 0000000000000000000000000000000000000000..2e4a80ec22349299deb6933f833912fb7db325d9
Binary files /dev/null and b/Paper/AllSkyMCNoiseOnly/allsky_noise_twoF_histogram.png differ
diff --git a/Paper/AllSkyMCNoiseOnly/generate_data.py b/Paper/AllSkyMCNoiseOnly/generate_data.py
new file mode 100644
index 0000000000000000000000000000000000000000..cc26b72bdc649338d211dd25eb64825c116f8979
--- /dev/null
+++ b/Paper/AllSkyMCNoiseOnly/generate_data.py
@@ -0,0 +1,101 @@
+import pyfstat
+import numpy as np
+import os
+import sys
+import time
+
+ID = sys.argv[1]
+outdir = sys.argv[2]
+
+label = 'run_{}'.format(ID)
+data_label = '{}_data'.format(label)
+results_file_name = '{}/MCResults_{}.txt'.format(outdir, ID)
+
+# Properties of the GW data
+sqrtSX = 2e-23
+tstart = 1000000000
+Tspan = 100*86400
+tend = tstart + Tspan
+
+# Fixed properties of the signal
+F0_center = 30
+F1_center = 1e-10
+F2 = 0
+tref = .5*(tstart+tend)
+
+
+VF0 = VF1 = 100
+DeltaF0 = VF0 * np.sqrt(3)/(np.pi*Tspan)
+DeltaF1 = VF1 * np.sqrt(45/4.)/(np.pi*Tspan**2)
+
+nsteps = 50
+run_setup = [((nsteps, 0), 20, False),
+ ((nsteps, 0), 11, False),
+ ((nsteps, 0), 6, False),
+ ((nsteps, 0), 3, False),
+ ((nsteps, nsteps), 1, False)]
+
+DeltaAlpha = 0.05
+DeltaDelta = 0.05
+
+h0 = 0
+F0 = F0_center + np.random.uniform(-0.5, 0.5)*DeltaF0
+F1 = F1_center + np.random.uniform(-0.5, 0.5)*DeltaF1
+Alpha = np.random.uniform(0, 2*np.pi)
+Delta = np.arccos(2*np.random.uniform(0, 1)-1)-np.pi/2
+fAlpha = np.random.uniform(0, 1)
+Alpha_min = Alpha - DeltaAlpha*(1-fAlpha)
+Alpha_max = Alpha + DeltaAlpha*fAlpha
+fDelta = np.random.uniform(0, 1)
+Delta_min = Delta - DeltaDelta*(1-fDelta)
+Delta_max = Delta + DeltaDelta*fDelta
+
+psi = np.random.uniform(-np.pi/4, np.pi/4)
+phi = np.random.uniform(0, 2*np.pi)
+cosi = np.random.uniform(-1, 1)
+
+data = pyfstat.Writer(
+ label=data_label, outdir=outdir, tref=tref,
+ tstart=tstart, F0=F0, F1=F1, F2=F2, duration=Tspan, Alpha=Alpha,
+ Delta=Delta, h0=h0, sqrtSX=sqrtSX, psi=psi, phi=phi, cosi=cosi,
+ detector='H1,L1')
+data.make_data()
+predicted_twoF = data.predict_fstat()
+
+startTime = time.time()
+theta_prior = {'F0': {'type': 'unif',
+ 'lower': F0-DeltaF0/2.,
+ 'upper': F0+DeltaF0/2.},
+ 'F1': {'type': 'unif',
+ 'lower': F1-DeltaF1/2.,
+ 'upper': F1+DeltaF1/2.},
+ 'F2': F2,
+ 'Alpha': {'type': 'unif',
+ 'lower': Alpha_min,
+ 'upper': Alpha_max},
+ 'Delta': {'type': 'unif',
+ 'lower': Delta_min,
+ 'upper': Delta_max},
+ }
+
+ntemps = 1
+log10temperature_min = -1
+nwalkers = 100
+
+mcmc = pyfstat.MCMCFollowUpSearch(
+ label=label, outdir=outdir,
+ sftfilepath='{}/*{}*sft'.format(outdir, data_label),
+ theta_prior=theta_prior,
+ tref=tref, minStartTime=tstart, maxStartTime=tend,
+ nwalkers=nwalkers, ntemps=ntemps,
+ log10temperature_min=log10temperature_min)
+mcmc.run(run_setup=run_setup, create_plots=False, log_table=False,
+ gen_tex_table=False)
+d, maxtwoF = mcmc.get_max_twoF()
+dF0 = F0 - d['F0']
+dF1 = F1 - d['F1']
+runTime = time.time() - startTime
+with open(results_file_name, 'a') as f:
+ f.write('{:1.8e} {:1.8e} {:1.8e} {:1.8e} {}\n'
+ .format(dF0, dF1, predicted_twoF, maxtwoF, runTime))
+os.system('rm {}/*{}*'.format(outdir, label))
diff --git a/Paper/AllSkyMCNoiseOnly/generate_table.py b/Paper/AllSkyMCNoiseOnly/generate_table.py
new file mode 100644
index 0000000000000000000000000000000000000000..e304f107748cf58b9c78625f9c1f816cb70bf8ac
--- /dev/null
+++ b/Paper/AllSkyMCNoiseOnly/generate_table.py
@@ -0,0 +1,85 @@
+import pyfstat
+import numpy as np
+
+outdir = 'data'
+
+label = 'allsky_setup'
+data_label = '{}_data'.format(label)
+
+# Properties of the GW data
+sqrtSX = 2e-23
+tstart = 1000000000
+Tspan = 100*86400
+tend = tstart + Tspan
+
+# Fixed properties of the signal
+F0_center = 30
+F1_center = 1e-10
+F2 = 0
+tref = .5*(tstart+tend)
+
+
+VF0 = VF1 = 100
+DeltaF0 = VF0 * np.sqrt(3)/(np.pi*Tspan)
+DeltaF1 = VF1 * np.sqrt(45/4.)/(np.pi*Tspan**2)
+
+DeltaAlpha = 0.05
+DeltaDelta = 0.05
+
+depth = 100
+
+nsteps = 50
+run_setup = [((nsteps, 0), 20, False),
+ ((nsteps, 0), 11, False),
+ ((nsteps, 0), 6, False),
+ ((nsteps, 0), 3, False),
+ ((nsteps, nsteps), 1, False)]
+
+h0 = sqrtSX / float(depth)
+r = np.random.uniform(0, 1)
+theta = np.random.uniform(0, 2*np.pi)
+F0 = F0_center + 3*np.sqrt(r)*np.cos(theta)/(np.pi**2 * Tspan**2)
+F1 = F1_center + 45*np.sqrt(r)*np.sin(theta)/(4*np.pi**2 * Tspan**4)
+
+Alpha = 0
+Delta = 0
+
+psi = np.random.uniform(-np.pi/4, np.pi/4)
+phi = np.random.uniform(0, 2*np.pi)
+cosi = np.random.uniform(-1, 1)
+
+data = pyfstat.Writer(
+ label=data_label, outdir=outdir, tref=tref,
+ tstart=tstart, F0=F0, F1=F1, F2=F2, duration=Tspan, Alpha=Alpha,
+ Delta=Delta, h0=h0, sqrtSX=sqrtSX, psi=psi, phi=phi, cosi=cosi,
+ detector='H1,L1')
+data.make_data()
+predicted_twoF = data.predict_fstat()
+
+theta_prior = {'F0': {'type': 'unif',
+ 'lower': F0-DeltaF0/2.,
+ 'upper': F0+DeltaF0/2.},
+ 'F1': {'type': 'unif',
+ 'lower': F1-DeltaF1/2.,
+ 'upper': F1+DeltaF1/2.},
+ 'F2': F2,
+ 'Alpha': {'type': 'unif',
+ 'lower': Alpha-DeltaAlpha/2.,
+ 'upper': Alpha+DeltaAlpha/2.},
+ 'Delta': {'type': 'unif',
+ 'lower': Delta-DeltaDelta/2.,
+ 'upper': Delta+DeltaDelta/2.},
+ }
+
+ntemps = 1
+log10temperature_min = -1
+nwalkers = 100
+
+mcmc = pyfstat.MCMCFollowUpSearch(
+ label=label, outdir=outdir,
+ sftfilepath='{}/*{}*sft'.format(outdir, data_label),
+ theta_prior=theta_prior,
+ tref=tref, minStartTime=tstart, maxStartTime=tend,
+ nwalkers=nwalkers, ntemps=ntemps,
+ log10temperature_min=log10temperature_min)
+mcmc.run(run_setup)
diff --git a/Paper/AllSkyMCNoiseOnly/plot_data.py b/Paper/AllSkyMCNoiseOnly/plot_data.py
new file mode 100644
index 0000000000000000000000000000000000000000..8165be5dd1ad2355967e278aed4a78887b9cf062
--- /dev/null
+++ b/Paper/AllSkyMCNoiseOnly/plot_data.py
@@ -0,0 +1,40 @@
+import matplotlib.pyplot as plt
+import pandas as pd
+import numpy as np
+import os
+from tqdm import tqdm
+from oct2py import octave
+import glob
+
+filenames = glob.glob("CollectedOutput/*.txt")
+
+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
+
+
+df_list = []
+for fn in filenames:
+ df = pd.read_csv(
+ fn, sep=' ', names=['dF0', 'dF1', 'twoF_predicted',
+ 'twoF', 'runTime'])
+ df['CLUSTER_ID'] = fn.split('_')[1]
+ df_list.append(df)
+df = pd.concat(df_list)
+
+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, 8e5), '-r')
+ax.set_xlabel('$\widetilde{2\mathcal{F}}$')
+ax.set_xlim(0, 60)
+fig.tight_layout()
+fig.savefig('allsky_noise_twoF_histogram.png')
diff --git a/Paper/AllSkyMCNoiseOnly/repeat.sh b/Paper/AllSkyMCNoiseOnly/repeat.sh
new file mode 100755
index 0000000000000000000000000000000000000000..5f43c35213b7e51665a94cdbbc5d67fc9b7b1e3c
--- /dev/null
+++ b/Paper/AllSkyMCNoiseOnly/repeat.sh
@@ -0,0 +1,12 @@
+#!/bin/bash
+
+. /home/gregory.ashton/lalsuite-install/etc/lalapps-user-env.sh
+export PATH="/home/gregory.ashton/anaconda2/bin:$PATH"
+export MPLCONFIGDIR=/home/gregory.ashton/.config/matplotlib
+
+rm /local/user/gregory.ashton/MCResults*txt
+for ((n=0;n<100;n++))
+do
+/home/gregory.ashton/anaconda2/bin/python generate_data.py "$1" /local/user/gregory.ashton --quite --no-template-counting
+done
+cp /local/user/gregory.ashton/MCResults*txt /home/gregory.ashton/PyFstat/Paper/AllSkyMCNoiseOnly/CollectedOutput
diff --git a/Paper/AllSkyMCNoiseOnly/submitfile b/Paper/AllSkyMCNoiseOnly/submitfile
new file mode 100644
index 0000000000000000000000000000000000000000..191a502764a6e221a87d0ccc89d49be609c8587c
--- /dev/null
+++ b/Paper/AllSkyMCNoiseOnly/submitfile
@@ -0,0 +1,12 @@
+Executable= repeat.sh
+Arguments= $(Cluster)_$(Process)
+Universe=vanilla
+Input=/dev/null
+accounting_group = ligo.dev.o2.cw.explore.test
+Output=CollectedOutput/out.$(Process)
+Error=CollectedOutput/err.$(Process)
+Log=CollectedOutput/log.$(Process)
+request_cpus = 1
+request_memory = 16 GB
+
+Queue 100
diff --git a/Paper/DirectedMCNoiseOnly/generate_data.py b/Paper/DirectedMCNoiseOnly/generate_data.py
new file mode 100644
index 0000000000000000000000000000000000000000..952ccd4cfac6276b1c3251d66bbdeddfbb64209d
--- /dev/null
+++ b/Paper/DirectedMCNoiseOnly/generate_data.py
@@ -0,0 +1,87 @@
+import pyfstat
+import numpy as np
+import os
+import sys
+import time
+
+
+ID = sys.argv[1]
+outdir = sys.argv[2]
+
+label = 'run_{}'.format(ID)
+data_label = '{}_data'.format(label)
+results_file_name = '{}/MCResults_{}.txt'.format(outdir, ID)
+
+# Properties of the GW data
+sqrtSX = 2e-23
+tstart = 1000000000
+Tspan = 100*86400
+tend = tstart + Tspan
+
+# Fixed properties of the signal
+F0_center = 30
+F1_center = 1e-10
+F2 = 0
+Alpha = 5e-3
+Delta = 6e-2
+tref = .5*(tstart+tend)
+
+VF0 = VF1 = 100
+DeltaF0 = VF0 * np.sqrt(3)/(np.pi*Tspan)
+DeltaF1 = VF1 * np.sqrt(45/4.)/(np.pi*Tspan**2)
+
+nsteps = 20
+run_setup = [((nsteps, 0), 20, False),
+ ((nsteps, 0), 7, False),
+ ((nsteps, 0), 2, False),
+ ((nsteps, nsteps), 1, False)]
+
+h0 = 0
+F0 = F0_center + np.random.uniform(-0.5, 0.5)*DeltaF0
+F1 = F1_center + np.random.uniform(-0.5, 0.5)*DeltaF1
+
+psi = np.random.uniform(-np.pi/4, np.pi/4)
+phi = np.random.uniform(0, 2*np.pi)
+cosi = np.random.uniform(-1, 1)
+
+data = pyfstat.Writer(
+ label=data_label, outdir=outdir, tref=tref,
+ tstart=tstart, F0=F0, F1=F1, F2=F2, duration=Tspan, Alpha=Alpha,
+ Delta=Delta, h0=h0, sqrtSX=sqrtSX, psi=psi, phi=phi, cosi=cosi,
+ detector='H1,L1')
+data.make_data()
+predicted_twoF = data.predict_fstat()
+
+startTime = time.time()
+theta_prior = {'F0': {'type': 'unif',
+ 'lower': F0-DeltaF0/2.,
+ 'upper': F0+DeltaF0/2.},
+ 'F1': {'type': 'unif',
+ 'lower': F1-DeltaF1/2.,
+ 'upper': F1+DeltaF1/2.},
+ 'F2': F2,
+ 'Alpha': Alpha,
+ 'Delta': Delta
+ }
+
+ntemps = 1
+log10temperature_min = -1
+nwalkers = 100
+
+mcmc = pyfstat.MCMCFollowUpSearch(
+ label=label, outdir=outdir,
+ sftfilepath='{}/*{}*sft'.format(outdir, data_label),
+ theta_prior=theta_prior,
+ tref=tref, minStartTime=tstart, maxStartTime=tend,
+ nwalkers=nwalkers, ntemps=ntemps,
+ log10temperature_min=log10temperature_min)
+mcmc.run(run_setup=run_setup, create_plots=False, log_table=False,
+ gen_tex_table=False)
+d, maxtwoF = mcmc.get_max_twoF()
+dF0 = F0 - d['F0']
+dF1 = F1 - d['F1']
+runTime = time.time() - startTime
+with open(results_file_name, 'a') as f:
+ f.write('{:1.8e} {:1.8e} {:1.8e} {:1.8e} {}\n'
+ .format(dF0, dF1, predicted_twoF, maxtwoF, runTime))
+os.system('rm {}/*{}*'.format(outdir, label))
diff --git a/Paper/DirectedMCNoiseOnly/repeat.sh b/Paper/DirectedMCNoiseOnly/repeat.sh
new file mode 100755
index 0000000000000000000000000000000000000000..0142257850dd0ab734558a72caa5ffae365cedfa
--- /dev/null
+++ b/Paper/DirectedMCNoiseOnly/repeat.sh
@@ -0,0 +1,12 @@
+#!/bin/bash
+
+. /home/gregory.ashton/lalsuite-install/etc/lalapps-user-env.sh
+export PATH="/home/gregory.ashton/anaconda2/bin:$PATH"
+export MPLCONFIGDIR=/home/gregory.ashton/.config/matplotlib
+
+rm /local/user/gregory.ashton/MCResults*txt
+for ((n=0;n<100;n++))
+do
+/home/gregory.ashton/anaconda2/bin/python generate_data.py "$1" /local/user/gregory.ashton --quite --no-template-counting
+done
+cp /local/user/gregory.ashton/MCResults*txt /home/gregory.ashton/PyFstat/Paper/DirectedMCNoiseOnly/CollectedOutput
diff --git a/Paper/DirectedMCNoiseOnly/submitfile b/Paper/DirectedMCNoiseOnly/submitfile
new file mode 100644
index 0000000000000000000000000000000000000000..7f057b75593f21a6a3a04b6c577f164c038340e0
--- /dev/null
+++ b/Paper/DirectedMCNoiseOnly/submitfile
@@ -0,0 +1,12 @@
+Executable= repeat.sh
+Arguments= $(Cluster)_$(Process)
+Universe=vanilla
+Input=/dev/null
+accounting_group = ligo.dev.o2.cw.explore.test
+Output=CollectedOutput/out.$(Process)
+Error=CollectedOutput/err.$(Process)
+Log=CollectedOutput/log.$(Process)
+request_cpus = 1
+request_memory = 8 GB
+
+Queue 100
diff --git a/Paper/directed_noise_twoF_histogram.png b/Paper/directed_noise_twoF_histogram.png
new file mode 100644
index 0000000000000000000000000000000000000000..e8b23a5f85c29004d657dc6b578ae161fcdb22f2
Binary files /dev/null and b/Paper/directed_noise_twoF_histogram.png differ