diff --git a/pyfstat/core.py b/pyfstat/core.py index 2b57fa2e3c19b7a11d97a752f284f6c9645a85bc..0fd551788f93140cd11852f7f5c4f4a6f4975578 100755 --- a/pyfstat/core.py +++ b/pyfstat/core.py @@ -151,8 +151,8 @@ class ComputeFstat(object): def __init__(self, tref, sftfilepath=None, minStartTime=None, maxStartTime=None, binary=False, transient=True, BSGL=False, detector=None, minCoverFreq=None, maxCoverFreq=None, - earth_ephem=None, sun_ephem=None, injectSources=None - ): + earth_ephem=None, sun_ephem=None, injectSources=None, + assumeSqrtSX=None): """ Parameters ---------- @@ -180,6 +180,10 @@ class ComputeFstat(object): Paths of the two files containing positions of Earth and Sun, respectively at evenly spaced times, as passed to CreateFstatInput. If None defaults defined in BaseSearchClass will be used. + assumeSqrtSX: float + Don't estimate noise-floors but assume (stationary) per-IFO + sqrt{SX} (if single value: use for all IFOs). If signal only, + set sqrtSX=1 """ @@ -250,7 +254,14 @@ class ComputeFstat(object): FstatOAs.runningMedianWindow = lalpulsar.FstatOptionalArgsDefaults.runningMedianWindow FstatOAs.FstatMethod = lalpulsar.FstatOptionalArgsDefaults.FstatMethod FstatOAs.InjectSqrtSX = lalpulsar.FstatOptionalArgsDefaults.injectSqrtSX - FstatOAs.assumeSqrtSX = lalpulsar.FstatOptionalArgsDefaults.assumeSqrtSX + if self.assumeSqrtSX is None: + FstatOAs.assumeSqrtSX = lalpulsar.FstatOptionalArgsDefaults.assumeSqrtSX + else: + mnf = lalpulsar.MultiNoiseFloor() + assumeSqrtSX = np.atleast_1d(self.assumeSqrtSX) + mnf.sqrtSn[:len(assumeSqrtSX)] = assumeSqrtSX + mnf.length = len(assumeSqrtSX) + FstatOAs.assumeSqrtSX = mnf FstatOAs.prevInput = lalpulsar.FstatOptionalArgsDefaults.prevInput FstatOAs.collectTiming = lalpulsar.FstatOptionalArgsDefaults.collectTiming @@ -460,7 +471,7 @@ class SemiCoherentSearch(BaseSearchClass, ComputeFstat): binary=False, BSGL=False, minStartTime=None, maxStartTime=None, minCoverFreq=None, maxCoverFreq=None, detector=None, earth_ephem=None, sun_ephem=None, - injectSources=None): + injectSources=None, assumeSqrtSX=None): """ Parameters ---------- @@ -570,7 +581,7 @@ class SemiCoherentGlitchSearch(BaseSearchClass, ComputeFstat): @helper_functions.initializer def __init__(self, label, outdir, tref, minStartTime, maxStartTime, nglitch=0, sftfilepath=None, theta0_idx=0, BSGL=False, - minCoverFreq=None, maxCoverFreq=None, + minCoverFreq=None, maxCoverFreq=None, assumeSqrtSX=None, detector=None, earth_ephem=None, sun_ephem=None): """ Parameters diff --git a/pyfstat/mcmc_based_searches.py b/pyfstat/mcmc_based_searches.py index f7ce22c84f313a20ab40dff31221505356180826..8b3ddab4d6115a8fb2b13632c12af15c6da18c76 100644 --- a/pyfstat/mcmc_based_searches.py +++ b/pyfstat/mcmc_based_searches.py @@ -13,7 +13,8 @@ import emcee import corner import dill as pickle -from core import BaseSearchClass, ComputeFstat +from core import BaseSearchClass, ComputeFstat, SemiCoherentSearch +from optimal_setup_functions import get_V_estimate from core import tqdm, args, earth_ephem, sun_ephem from optimal_setup_functions import get_optimal_setup import helper_functions @@ -28,7 +29,7 @@ class MCMCSearch(BaseSearchClass): theta_initial=None, scatter_val=1e-10, binary=False, BSGL=False, minCoverFreq=None, maxCoverFreq=None, detector=None, earth_ephem=None, - sun_ephem=None, injectSources=None): + sun_ephem=None, injectSources=None, assumeSqrtSX=None): """ Parameters label, outdir: str @@ -113,7 +114,8 @@ class MCMCSearch(BaseSearchClass): earth_ephem=self.earth_ephem, sun_ephem=self.sun_ephem, detector=self.detector, BSGL=self.BSGL, transient=False, minStartTime=self.minStartTime, maxStartTime=self.maxStartTime, - binary=self.binary, injectSources=self.injectSources) + binary=self.binary, injectSources=self.injectSources, + assumeSqrtSX=self.assumeSqrtSX) def logp(self, theta_vals, theta_prior, theta_keys, search): H = [self.generic_lnprior(**theta_prior[key])(p) for p, key in @@ -1247,7 +1249,7 @@ class MCMCSemiCoherentSearch(MCMCSearch): scatter_val=1e-10, detector=None, BSGL=False, minStartTime=None, maxStartTime=None, minCoverFreq=None, maxCoverFreq=None, earth_ephem=None, sun_ephem=None, - injectSources=None): + injectSources=None, assumeSqrtSX=None): """ """ @@ -1284,7 +1286,7 @@ class MCMCSemiCoherentSearch(MCMCSearch): maxStartTime=self.maxStartTime, minCoverFreq=self.minCoverFreq, maxCoverFreq=self.maxCoverFreq, detector=self.detector, earth_ephem=self.earth_ephem, sun_ephem=self.sun_ephem, - injectSources=self.injectSources) + injectSources=self.injectSources, assumeSqrtSX=self.assumeSqrtSX) def logp(self, theta_vals, theta_prior, theta_keys, search): H = [self.generic_lnprior(**theta_prior[key])(p) for p, key in diff --git a/tests.py b/tests.py index 04ed03987dca454f6e388d779a7ae88fc5a8fbba..4c1f5651e4901bc5e9257ccf1222dd2926a4be02 100644 --- a/tests.py +++ b/tests.py @@ -105,6 +105,23 @@ class TestComputeFstat(Test): print predicted_FS, FS self.assertTrue(np.abs(predicted_FS-FS)/FS < 0.2) + def run_computefstatistic_single_point_no_noise(self): + Writer = pyfstat.Writer(self.label, outdir=outdir, add_noise=False) + Writer.make_data() + predicted_FS = Writer.predict_fstat() + + search = pyfstat.ComputeFstat( + tref=Writer.tref, assumeSqrtSX=1, + sftfilepath='{}/*{}*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) + print predicted_FS, FS + self.assertTrue(np.abs(predicted_FS-FS)/FS < 0.2) class TestSemiCoherentGlitchSearch(Test): label = "Test"