From a0eee809ebe8a4f5daec839d39d2f1619f997241 Mon Sep 17 00:00:00 2001
From: Gregory Ashton <gregory.ashton@ligo.org>
Date: Fri, 29 Sep 2017 18:16:27 +0200
Subject: [PATCH] Renaming functions for clarity of twoF and fix issue of F ->
 likelihood

Fixes #4 .
---
 pyfstat/core.py                | 28 ++++++++----------------
 pyfstat/grid_based_searches.py |  6 +++---
 pyfstat/mcmc_based_searches.py | 26 +++++++++++------------
 tests.py                       | 39 ++++++++++++----------------------
 4 files changed, 38 insertions(+), 61 deletions(-)

diff --git a/pyfstat/core.py b/pyfstat/core.py
index 88bac00..0ff84cb 100755
--- a/pyfstat/core.py
+++ b/pyfstat/core.py
@@ -576,19 +576,9 @@ class ComputeFstat(BaseSearchClass):
             self.windowRange.tauBand = 0
             self.windowRange.dtau = 1
 
-    def compute_fullycoherent_det_stat_single_point(
-            self, F0, F1, F2, Alpha, Delta, asini=None, period=None, ecc=None,
-            tp=None, argp=None):
-        """ Compute the fully-coherent det. statistic at a single point """
-
-        return self.run_computefstatistic_single_point(
-            self.minStartTime, self.maxStartTime, F0, F1, F2, Alpha, Delta,
-            asini, period, ecc, tp, argp)
-
-    def run_computefstatistic_single_point(self, tstart, tend, F0, F1,
-                                           F2, Alpha, Delta, asini=None,
-                                           period=None, ecc=None, tp=None,
-                                           argp=None):
+    def get_fullycoherent_twoF(self, tstart, tend, F0, F1, F2, Alpha, Delta,
+                               asini=None, period=None, ecc=None, tp=None,
+                               argp=None):
         """ Returns twoF or ln(BSGL) fully-coherently at a single point """
 
         self.PulsarDopplerParams.fkdot = np.array([F0, F1, F2, 0, 0, 0, 0])
@@ -684,7 +674,7 @@ class ComputeFstat(BaseSearchClass):
             self.transient = True
             self.init_computefstatistic_single_point()
         for tau in taus:
-            twoFs.append(self.run_computefstatistic_single_point(
+            twoFs.append(self.get_fullycoherent_twoF(
                 tstart=tstart, tend=tstart+tau, F0=F0, F1=F1, F2=F2,
                 Alpha=Alpha, Delta=Delta, asini=asini, period=period, ecc=ecc,
                 tp=tp, argp=argp))
@@ -875,7 +865,7 @@ class SemiCoherentSearch(ComputeFstat):
                     'Semi-coherent end time {} after last SFT timestamp {}'
                     .format(self.tboundaries[-1], self.SFT_timestamps[-1]))
 
-    def run_semi_coherent_computefstatistic_single_point(
+    def get_semicoherent_twoF(
             self, F0, F1, F2, Alpha, Delta, asini=None,
             period=None, ecc=None, tp=None, argp=None,
             record_segments=False):
@@ -992,7 +982,7 @@ class SemiCoherentGlitchSearch(ComputeFstat):
         self.binary = False
         self.init_computefstatistic_single_point()
 
-    def compute_nglitch_fstat(self, F0, F1, F2, Alpha, Delta, *args):
+    def get_semicoherent_nglitch_twoF(self, F0, F1, F2, Alpha, Delta, *args):
         """ Returns the semi-coherent glitch summed twoF """
 
         args = list(args)
@@ -1013,7 +1003,7 @@ class SemiCoherentGlitchSearch(ComputeFstat):
         for i, theta_i_at_tref in enumerate(thetas):
             ts, te = tboundaries[i], tboundaries[i+1]
 
-            twoFVal = self.run_computefstatistic_single_point(
+            twoFVal = self.get_fullycoherent_twoF(
                 ts, te, theta_i_at_tref[1], theta_i_at_tref[2],
                 theta_i_at_tref[3], Alpha, Delta)
             twoFSum += twoFVal
@@ -1039,14 +1029,14 @@ class SemiCoherentGlitchSearch(ComputeFstat):
         theta_post_glitch = self._shift_coefficients(
             theta_post_glitch_at_glitch, tref - tglitch)
 
-        twoFsegA = self.run_computefstatistic_single_point(
+        twoFsegA = self.get_fullycoherent_twoF(
             self.minStartTime, tglitch, theta[0], theta[1], theta[2], Alpha,
             Delta)
 
         if tglitch == self.maxStartTime:
             return twoFsegA
 
-        twoFsegB = self.run_computefstatistic_single_point(
+        twoFsegB = self.get_fullycoherent_twoF(
             tglitch, self.maxStartTime, theta_post_glitch[0],
             theta_post_glitch[1], theta_post_glitch[2], Alpha,
             Delta)
diff --git a/pyfstat/grid_based_searches.py b/pyfstat/grid_based_searches.py
index 99e16e8..42c84fa 100644
--- a/pyfstat/grid_based_searches.py
+++ b/pyfstat/grid_based_searches.py
@@ -60,7 +60,7 @@ class GridSearch(BaseSearchClass):
                 BSGL=self.BSGL, SSBprec=self.SSBprec,
                 injectSources=self.injectSources,
                 assumeSqrtSX=self.assumeSqrtSX)
-            self.search.get_det_stat = self.search.run_computefstatistic_single_point
+            self.search.get_det_stat = self.search.get_fullycoherent_twoF
         else:
             self.search = SemiCoherentSearch(
                 label=self.label, outdir=self.outdir, tref=self.tref,
@@ -71,7 +71,7 @@ class GridSearch(BaseSearchClass):
                 injectSources=self.injectSources)
 
             def cut_out_tstart_tend(*vals):
-                return self.search.run_semi_coherent_computefstatistic_single_point(*vals[2:])
+                return self.search.get_semicoherent_twoF(*vals[2:])
             self.search.get_det_stat = cut_out_tstart_tend
 
     def get_array_from_tuple(self, x):
@@ -397,7 +397,7 @@ class FrequencySlidingWindow(GridSearch):
             BSGL=self.BSGL, SSBprec=self.SSBprec,
             injectSources=self.injectSources)
         self.search.get_det_stat = (
-            self.search.run_computefstatistic_single_point)
+            self.search.get_fullycoherent_twoF)
 
     def get_input_data_array(self):
         arrays = []
diff --git a/pyfstat/mcmc_based_searches.py b/pyfstat/mcmc_based_searches.py
index cc3ebaf..9e12524 100644
--- a/pyfstat/mcmc_based_searches.py
+++ b/pyfstat/mcmc_based_searches.py
@@ -160,9 +160,9 @@ class MCMCSearch(core.BaseSearchClass):
     def logl(self, theta, search):
         for j, theta_i in enumerate(self.theta_idxs):
             self.fixed_theta[theta_i] = theta[j]
-        FS = search.compute_fullycoherent_det_stat_single_point(
-            *self.fixed_theta)
-        return FS + self.likelihoodcoef
+        twoF = search.get_fullycoherent_twoF(
+            self.minStartTime, self.maxStartTime, *self.fixed_theta)
+        return twoF/2.0 + self.likelihoodcoef
 
     def _unpack_input_theta(self):
         full_theta_keys = ['F0', 'F1', 'F2', 'Alpha', 'Delta']
@@ -1244,11 +1244,11 @@ class MCMCSearch(core.BaseSearchClass):
 
         """
         if any(np.isposinf(self.lnlikes)):
-            logging.info('twoF values contain positive infinite values')
+            logging.info('lnlike values contain positive infinite values')
         if any(np.isneginf(self.lnlikes)):
-            logging.info('twoF values contain negative infinite values')
+            logging.info('lnlike values contain negative infinite values')
         if any(np.isnan(self.lnlikes)):
-            logging.info('twoF values contain nan')
+            logging.info('lnlike values contain nan')
         idxs = np.isfinite(self.lnlikes)
         jmax = np.nanargmax(self.lnlikes[idxs])
         maxlogl = self.lnlikes[jmax]
@@ -1262,7 +1262,7 @@ class MCMCSearch(core.BaseSearchClass):
             maxtwoF = self.logl(p, self.search)
             self.search.BSGL = self.BSGL
         else:
-            maxtwoF = maxlogl - self.likelihoodcoef
+            maxtwoF = (maxlogl - self.likelihoodcoef)*2
 
         repeats = []
         for i, k in enumerate(self.theta_keys):
@@ -1633,8 +1633,8 @@ class MCMCGlitchSearch(MCMCSearch):
 
         for j, theta_i in enumerate(self.theta_idxs):
             self.fixed_theta[theta_i] = theta[j]
-        FS = search.compute_nglitch_fstat(*self.fixed_theta)
-        return FS + self.likelihoodcoef
+        twoF = search.get_semicoherent_nglitch_twoF(*self.fixed_theta)
+        return twoF/2.0 + self.likelihoodcoef
 
     def _unpack_input_theta(self):
         glitch_keys = ['delta_F0', 'delta_F1', 'tglitch']
@@ -1839,9 +1839,9 @@ class MCMCSemiCoherentSearch(MCMCSearch):
     def logl(self, theta, search):
         for j, theta_i in enumerate(self.theta_idxs):
             self.fixed_theta[theta_i] = theta[j]
-        FS = search.run_semi_coherent_computefstatistic_single_point(
+        twoF = search.get_semicoherent_twoF(
             *self.fixed_theta)
-        return FS + self.likelihoodcoef
+        return twoF/2.0 + self.likelihoodcoef
 
 
 class MCMCFollowUpSearch(MCMCSemiCoherentSearch):
@@ -2152,8 +2152,8 @@ class MCMCTransientSearch(MCMCSearch):
         in_theta[1] = in_theta[0] + in_theta[1]
         if in_theta[1] > self.maxStartTime:
             return -np.inf
-        FS = search.run_computefstatistic_single_point(*in_theta)
-        return FS + self.likelihoodcoef
+        twoF = search.get_fullycoherent_twoF(*in_theta)
+        return twoF/2.0 + self.likelihoodcoef
 
     def _unpack_input_theta(self):
         full_theta_keys = ['transient_tstart',
diff --git a/tests.py b/tests.py
index c6408bb..4ece9d9 100644
--- a/tests.py
+++ b/tests.py
@@ -97,13 +97,9 @@ class TestComputeFstat(Test):
         search = pyfstat.ComputeFstat(
             tref=Writer.tref,
             sftfilepattern='{}/*{}*sft'.format(Writer.outdir, Writer.label))
-        FS = search.run_computefstatistic_single_point(Writer.tstart,
-                                                       Writer.tend,
-                                                       Writer.F0,
-                                                       Writer.F1,
-                                                       Writer.F2,
-                                                       Writer.Alpha,
-                                                       Writer.Delta)
+        FS = search.get_fullycoherent_twoF(
+            Writer.tstart, Writer.tend, Writer.F0, Writer.F1, Writer.F2,
+            Writer.Alpha, Writer.Delta)
         print predicted_FS, FS
         self.assertTrue(np.abs(predicted_FS-FS)/FS < 0.2)
 
@@ -116,13 +112,9 @@ class TestComputeFstat(Test):
         search = pyfstat.ComputeFstat(
             tref=Writer.tref, assumeSqrtSX=1,
             sftfilepattern='{}/*{}*sft'.format(Writer.outdir, Writer.label))
-        FS = search.run_computefstatistic_single_point(Writer.tstart,
-                                                       Writer.tend,
-                                                       Writer.F0,
-                                                       Writer.F1,
-                                                       Writer.F2,
-                                                       Writer.Alpha,
-                                                       Writer.Delta)
+        FS = search.get_fullycoherent_twoF(
+            Writer.tstart, Writer.tend, Writer.F0, Writer.F1, Writer.F2,
+            Writer.Alpha, Writer.Delta)
         print predicted_FS, FS
         self.assertTrue(np.abs(predicted_FS-FS)/FS < 0.2)
 
@@ -137,13 +129,9 @@ class TestComputeFstat(Test):
             minCoverFreq=28, maxCoverFreq=32, minStartTime=Writer.tstart,
             maxStartTime=Writer.tstart+Writer.duration,
             detectors=Writer.detectors)
-        FS = search.run_computefstatistic_single_point(Writer.tstart,
-                                                       Writer.tend,
-                                                       Writer.F0,
-                                                       Writer.F1,
-                                                       Writer.F2,
-                                                       Writer.Alpha,
-                                                       Writer.Delta)
+        FS = search.get_fullycoherent_twoF(
+            Writer.tstart, Writer.tend, Writer.F0, Writer.F1, Writer.F2,
+            Writer.Alpha, Writer.Delta)
         Writer.make_data()
         predicted_FS = Writer.predict_fstat()
         print predicted_FS, FS
@@ -153,7 +141,7 @@ class TestComputeFstat(Test):
 class TestSemiCoherentGlitchSearch(Test):
     label = "TestSemiCoherentGlitchSearch"
 
-    def test_compute_nglitch_fstat(self):
+    def test_get_semicoherent_nglitch_twoF(self):
         duration = 10*86400
         dtglitch = .5*duration
         delta_F0 = 0
@@ -171,10 +159,9 @@ class TestSemiCoherentGlitchSearch(Test):
             tref=Writer.tref, minStartTime=Writer.tstart,
             maxStartTime=Writer.tend, nglitch=1)
 
-        FS = search.compute_nglitch_fstat(Writer.F0, Writer.F1, Writer.F2,
-                                          Writer.Alpha, Writer.Delta,
-                                          Writer.delta_F0, Writer.delta_F1,
-                                          search.minStartTime+dtglitch)
+        FS = search.get_semicoherent_nglitch_twoF(
+            Writer.F0, Writer.F1, Writer.F2, Writer.Alpha, Writer.Delta,
+            Writer.delta_F0, Writer.delta_F1, search.minStartTime+dtglitch)
 
         # Compute the predicted semi-coherent glitch Fstat
         minStartTime = Writer.tstart
-- 
GitLab