diff --git a/pyfstat/mcmc_based_searches.py b/pyfstat/mcmc_based_searches.py index f5a9f246a79580f161e646204a6ca0f6b45f216a..42b6139b066b9a7169841d5741fa937fd8ba06e8 100644 --- a/pyfstat/mcmc_based_searches.py +++ b/pyfstat/mcmc_based_searches.py @@ -1872,9 +1872,9 @@ class MCMCFollowUpSearch(MCMCSemiCoherentSearch): d = pickle.load(f) return d - def write_setup_input_file(self, run_setup_input_file, R, Nsegs0, + def write_setup_input_file(self, run_setup_input_file, NstarMax, Nsegs0, nsegs_vals, Nstar_vals, theta_prior): - d = dict(R=R, Nsegs0=Nsegs0, nsegs_vals=nsegs_vals, + d = dict(NstarMax=NstarMax, Nsegs0=Nsegs0, nsegs_vals=nsegs_vals, theta_prior=theta_prior, Nstar_vals=Nstar_vals) with open(run_setup_input_file, 'w+') as f: pickle.dump(d, f) @@ -1886,19 +1886,19 @@ class MCMCFollowUpSearch(MCMCSemiCoherentSearch): return True else: logging.info( - 'Old setup does not match one of R, Nsegs0 or prior') + 'Old setup does not match one of NstarMax, Nsegs0 or prior') except KeyError as e: logging.info( 'Error found when comparing with old setup: {}'.format(e)) return False - def init_run_setup(self, run_setup=None, R=10, Nsegs0=None, log_table=True, - gen_tex_table=True): + def init_run_setup(self, run_setup=None, NstarMax=1000, Nsegs0=None, + log_table=True, gen_tex_table=True): if run_setup is None and Nsegs0 is None: raise ValueError( - 'You must either specify the run_setup, or Nsegs0 from which ' - 'the optimal run_setup given R can be estimated') + 'You must either specify the run_setup, or Nsegs0 and NStarMax' + ' from which the optimal run_setup can be estimated') if run_setup is None: logging.info('No run_setup provided') @@ -1909,11 +1909,12 @@ class MCMCFollowUpSearch(MCMCSemiCoherentSearch): logging.info('Checking old setup input file {}'.format( run_setup_input_file)) old_setup = self.read_setup_input_file(run_setup_input_file) - if self.check_old_run_setup(old_setup, R=R, + if self.check_old_run_setup(old_setup, NstarMax=NstarMax, Nsegs0=Nsegs0, theta_prior=self.theta_prior): - logging.info('Using old setup with R={}, Nsegs0={}'.format( - R, Nsegs0)) + logging.info( + 'Using old setup with NstarMax={}, Nsegs0={}'.format( + NstarMax, Nsegs0)) nsegs_vals = old_setup['nsegs_vals'] Nstar_vals = old_setup['Nstar_vals'] generate_setup = False @@ -1924,11 +1925,11 @@ class MCMCFollowUpSearch(MCMCSemiCoherentSearch): if generate_setup: nsegs_vals, Nstar_vals = get_optimal_setup( - R, Nsegs0, self.tref, self.minStartTime, + NstarMax, Nsegs0, self.tref, self.minStartTime, self.maxStartTime, self.theta_prior, self.search.detector_names, self.earth_ephem, self.sun_ephem) - self.write_setup_input_file(run_setup_input_file, R, Nsegs0, + self.write_setup_input_file(run_setup_input_file, NstarMax, Nsegs0, nsegs_vals, Nstar_vals, self.theta_prior) @@ -2002,9 +2003,9 @@ class MCMCFollowUpSearch(MCMCSemiCoherentSearch): else: return run_setup - def run(self, run_setup=None, proposal_scale_factor=2, R=10, Nsegs0=None, - create_plots=True, log_table=True, gen_tex_table=True, fig=None, - axes=None, return_fig=False, **kwargs): + def run(self, run_setup=None, proposal_scale_factor=2, NstarMax=10, + Nsegs0=None, create_plots=True, log_table=True, gen_tex_table=True, + fig=None, axes=None, return_fig=False, **kwargs): """ Run the follow-up with the given run_setup Parameters @@ -2029,7 +2030,7 @@ class MCMCFollowUpSearch(MCMCSemiCoherentSearch): self.nsegs = 1 self._initiate_search_object() run_setup = self.init_run_setup( - run_setup, R=R, Nsegs0=Nsegs0, log_table=log_table, + run_setup, NstarMax=NstarMax, Nsegs0=Nsegs0, log_table=log_table, gen_tex_table=gen_tex_table) self.run_setup = run_setup diff --git a/pyfstat/optimal_setup_functions.py b/pyfstat/optimal_setup_functions.py index 80a7072d57e4c25ee690573ceb999007b234e818..dda36c233bebb4cb9d95daee0819c34ca83d19d2 100644 --- a/pyfstat/optimal_setup_functions.py +++ b/pyfstat/optimal_setup_functions.py @@ -14,13 +14,13 @@ import pyfstat.helper_functions as helper_functions def get_optimal_setup( - R, Nsegs0, tref, minStartTime, maxStartTime, prior, + NstarMax, Nsegs0, tref, minStartTime, maxStartTime, prior, detector_names, earth_ephem, sun_ephem): """ Using an optimisation step, calculate the optimal setup ladder Parameters ---------- - R : float + NstarMax : float Nsegs0 : int The number of segments for the initial step of the ladder minStartTime, maxStartTime : int @@ -38,8 +38,8 @@ def get_optimal_setup( """ - logging.info('Calculating optimal setup for R={}, Nsegs0={}'.format( - R, Nsegs0)) + logging.info('Calculating optimal setup for NstarMax={}, Nsegs0={}'.format( + NstarMax, Nsegs0)) Nstar_0 = get_Nstar_estimate( Nsegs0, tref, minStartTime, maxStartTime, prior, @@ -54,7 +54,7 @@ def get_optimal_setup( nsegs_i = Nsegs0 while nsegs_i > 1: nsegs_i, Nstar_i = _get_nsegs_ip1( - nsegs_i, R, tref, minStartTime, maxStartTime, prior, + nsegs_i, NstarMax, tref, minStartTime, maxStartTime, prior, detector_names, earth_ephem, sun_ephem) nsegs_vals.append(nsegs_i) Nstar_vals.append(Nstar_i) @@ -65,11 +65,11 @@ def get_optimal_setup( return nsegs_vals, Nstar_vals -def _get_nsegs_ip1(nsegs_i, R, tref, minStartTime, maxStartTime, prior, +def _get_nsegs_ip1(nsegs_i, NstarMax, tref, minStartTime, maxStartTime, prior, detector_names, earth_ephem, sun_ephem): """ Calculate Nsegs_{i+1} given Nsegs_{i} """ - log10R = np.log10(R) + log10NstarMax = np.log10(NstarMax) log10Nstari = np.log10(get_Nstar_estimate( nsegs_i, tref, minStartTime, maxStartTime, prior, detector_names, earth_ephem, sun_ephem)) @@ -89,9 +89,10 @@ def _get_nsegs_ip1(nsegs_i, R, tref, minStartTime, maxStartTime, prior, return 1e6 else: log10Nstarip1 = np.log10(Nstarip1) - return np.abs(log10Nstari + log10R - log10Nstarip1) - res = scipy.optimize.minimize(f, .5*nsegs_i, method='Powell', tol=0.1, + return np.abs(log10Nstari + log10NstarMax - log10Nstarip1) + res = scipy.optimize.minimize(f, .4*nsegs_i, method='Powell', tol=1, options={'maxiter': 10}) + logging.info('{} with {} evaluations'.format(res['message'], res['nfev'])) nsegs_ip1 = int(res.x) if nsegs_ip1 == 0: nsegs_ip1 = 1