Updates to the paper: transient and follow-up

Paper/paper_cw_mcmc.tex

@@ -792,7 +792,7 @@ We now provide an illustrative example of the follow-up method. We consider a
 directed search over the sky position and frequency in 100 days of data from a
 single detector, with $\sqrt{\Sn}=10^{-23}$~Hz$^{-1/2}$ (at the fiducial
 frequency of the signal). The simulated signal has an amplitude
-$h_0=1.4\times10^{25}$ such that the signal has a depth of $\sqrt{\Sn}/h_0=70$
+$h_0=2\times10^{-25}$ such that the signal has a depth of $\sqrt{\Sn}/h_0=50$
 in the noise.
 
 First, we must define the setup for the run. Using $\mathcal{R}=10$ and
@@ -834,14 +834,28 @@ are listed in Table~\ref{tab_weak_signal_follow_up}.}
 \end{figure}
 
 The key advantage to note here is that all walkers succefully convereged to the
-signal peak, which occupies $\approx 10^{-6}$ of the initial volume. While it
-is possible for this to occur during an ordinary MCMC simulation (with $\Tcoh$
-fixed at $\Tspan$), it would take much longer to converge as the chains explore
-the other `noise peaks' in the data.
+signal peak, which occupies $\sim 10^{-6}$ of the initial volume. While it is
+possible for this to occur during an ordinary MCMC simulation (with $\Tcoh$
+fixed at $\Tspan$), it would take substantially longer to converge as the
+chains explore the other `noise peaks' in the data.
 
 \section{Alternative waveform models: transients}
 \label{sec_transients}
 
+\begin{figure}[htb]
+\centering
+\includegraphics[width=0.5\textwidth]{transient_search_initial_stage_twoFcumulative}
+\caption{}
+\label{fig:}
+\end{figure}
+
+\begin{figure}[htb]
+\centering
+\includegraphics[width=0.5\textwidth]{transient_search_corner}
+\caption{}
+\label{fig:}
+\end{figure}
+
 \section{Alternative waveform models: glitches}
 \label{sec_glitches}
 

Paper/weak_signal_follow_up_run_setup.tex

 \begin{tabular}{c|cccccc}
 Stage & $\Nseg$ & $\Tcoh^{\rm days}$ &$\Nsteps$ & $\V$ & $\Vsky$ & $\Vpe$ \\ \hline
-0 & 93 & 1.1 & 100 & 10.0 & 1.0 & 10.0 \\
-1 & 43 & 2.3 & 100 & $1{\times}10^{2}$ & 6.0 & 20.0 \\
-2 & 20 & 5.0 & 100 & $1{\times}10^{3}$ & 30.0 & 50.0 \\
-3 & 9 & 11.1 & 100 & $1{\times}10^{4}$ & $1{\times}10^{2}$ & $1{\times}10^{2}$ \\
-4 & 4 & 25.0 & 100 & $1{\times}10^{5}$ & $6{\times}10^{2}$ & $2{\times}10^{2}$ \\
-5 & 1 & 100.0 & 100,100 & $1{\times}10^{6}$ & $1{\times}10^{3}$ & $9{\times}10^{2}$ \\
+0 & 93 & 1.1 & 100 & 20.0 & 2.0 & 10.0 \\
+1 & 43 & 2.3 & 100 & $2{\times}10^{2}$ & 10.0 & 20.0 \\
+2 & 20 & 5.0 & 100 & $2{\times}10^{3}$ & 50.0 & 50.0 \\
+3 & 9 & 11.1 & 100 & $2{\times}10^{4}$ & $2{\times}10^{2}$ & $1{\times}10^{2}$ \\
+4 & 4 & 25.0 & 100 & $2{\times}10^{5}$ & $1{\times}10^{3}$ & $2{\times}10^{2}$ \\
+5 & 1 & 100.0 & 100,100 & $2{\times}10^{6}$ & $3{\times}10^{3}$ & $9{\times}10^{2}$ \\
 \end{tabular}

examples/make_fake_data.py

@@ -53,14 +53,3 @@ two_glitch_data = Writer(
     dtglitch=dtglitch, delta_phi=delta_phi, delta_F0=delta_F0,
     delta_F1=delta_F1, delta_F2=delta_F2)
 two_glitch_data.make_data()
-
-
-# Making transient data in the middle third
-data_tstart = tstart - duration
-data_duration = 3 * duration
-transient = Writer(
-    label='transient', outdir='data', tref=tref, tstart=tstart, F0=F0, F1=F1,
-    F2=F2, duration=duration, Alpha=Alpha, Delta=Delta, h0=h0, sqrtSX=sqrtSX,
-    data_tstart=data_tstart, data_duration=data_duration)
-transient.make_data()
-

examples/transient_search_using_MCMC.py

-from pyfstat import MCMCTransientSearch
+import pyfstat
+import numpy as np
 
 F0 = 30.0
 F1 = -1e-10
 F2 = 0
 Alpha = 5e-3
 Delta = 6e-2
-tref = 362750407.0
 
 tstart = 1000000000
 duration = 100*86400
-tstart = 1000000000 - duration
-tend = tstart + 3*duration
+data_tstart = tstart - duration
+data_tend = data_tstart + 3*duration
+tref = .5*(data_tstart+data_tend)
 
-theta_prior = {'F0': {'type': 'unif', 'lower': F0*(1-1e-6), 'upper': F0*(1+1e-6)},
-               'F1': {'type': 'unif', 'lower': F1*(1+1e-2), 'upper': F1*(1-1e-2)},
+h0 = 1e-23
+sqrtSX = 1e-22
+
+transient = pyfstat.Writer(
+    label='transient', outdir='data', tref=tref, tstart=tstart, F0=F0, F1=F1,
+    F2=F2, duration=duration, Alpha=Alpha, Delta=Delta, h0=h0, sqrtSX=sqrtSX,
+    minStartTime=data_tstart, maxStartTime=data_tend)
+transient.make_data()
+
+DeltaF0 = 6e-7
+DeltaF1 = 1e-13
+VF0 = (np.pi * duration * DeltaF0)**2 / 3.0
+VF1 = (np.pi * duration**2 * DeltaF1)**2 * 4/45.
+print '\nV={:1.2e}, VF0={:1.2e}, VF1={:1.2e}\n'.format(VF0*VF1, VF0, VF1)
+
+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,
-               'transient_tstart': {'type': 'unif', 'lower': tstart, 'upper': tend},
-               'transient_duration': {'type': 'halfnorm', 'loc':0, 'scale': duration}
+               'Delta': Delta
                }
 
-ntemps = 4
+ntemps = 3
 log10temperature_min = -1
 nwalkers = 100
-nsteps = [1000, 1000]
+nsteps = [750, 250]
+
+mcmc = pyfstat.MCMCSearch(
+    label='transient_search_initial_stage', outdir='data',
+    sftfilepath='data/*transient*sft', theta_prior=theta_prior, tref=tref,
+    minStartTime=data_tstart, maxStartTime=data_tend, nsteps=nsteps,
+    nwalkers=nwalkers, ntemps=ntemps,
+    log10temperature_min=log10temperature_min)
+mcmc.run()
+mcmc.plot_cumulative_max()
+
+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,
+               'transient_tstart': {'type': 'unif',
+                                    'lower': data_tstart,
+                                    'upper': data_tend},
+               'transient_duration': {'type': 'halfnorm',
+                                      'loc': 0,
+                                      'scale': 0.5*duration}
+               }
+
+nwalkers = 500
+nsteps = [200, 200]
 
-mcmc = MCMCTransientSearch(
-    label='transient_search_using_MCMC', outdir='data',
+mcmc = pyfstat.MCMCTransientSearch(
+    label='transient_search', outdir='data',
     sftfilepath='data/*transient*sft', theta_prior=theta_prior, tref=tref,
-    tstart=tstart, tend=tend, nsteps=nsteps, nwalkers=nwalkers, ntemps=ntemps,
+    minStartTime=data_tstart, maxStartTime=data_tend, nsteps=nsteps,
+    nwalkers=nwalkers, ntemps=ntemps,
     log10temperature_min=log10temperature_min)
 mcmc.run()
 mcmc.plot_corner(add_prior=True)

examples/weak_signal_follow_up.py

@@ -13,7 +13,7 @@ duration = 100*86400
 tend = tstart+duration
 tref = .5*(tstart+tend)
 
-depth = 70
+depth = 50
 data_label = 'weak_signal_follow_up_depth_{:1.0f}'.format(depth)
 
 h0 = sqrtSX / depth
@@ -41,8 +41,8 @@ theta_prior = {'F0': {'type': 'unif', 'lower': F0*(1-1e-6),
                }
 
 ntemps = 3
-log10temperature_min = -1
-nwalkers = 200
+log10temperature_min = -0.5
+nwalkers = 100
 scatter_val = 1e-10
 nsteps = [100, 100]
 
@@ -52,7 +52,6 @@ mcmc = pyfstat.MCMCFollowUpSearch(
     minStartTime=tstart, maxStartTime=tend, nwalkers=nwalkers, nsteps=nsteps,
     ntemps=ntemps, log10temperature_min=log10temperature_min,
     scatter_val=scatter_val)
-mcmc.run(R0=10, Vmin=100)
+mcmc.run(R0=10, Vmin=100, subtractions=[F0, Alpha, Delta], context='paper')
 mcmc.plot_corner(add_prior=True)
 mcmc.print_summary()
-#mcmc.generate_loudest()

pyfstat.py

@@ -655,6 +655,7 @@ class ComputeFstat(object):
         else:
             ax.set_ylabel(r'$\widetilde{2\mathcal{F}}_{\rm cumulative}$')
         ax.set_xlim(0, taus[-1]/86400)
+        if title:
             ax.set_title(title)
         if savefig:
             plt.savefig('{}/{}_twoFcumulative.png'.format(outdir, label))