diff --git a/CHANGELOG-unreleased.md b/CHANGELOG-unreleased.md index e3df9c826..bd2cb7039 100644 --- a/CHANGELOG-unreleased.md +++ b/CHANGELOG-unreleased.md @@ -18,6 +18,7 @@ the released changes. - `Parameter.as_latex` method for latex representation of a parameter. - `pint.output.publish` module and `pintpublish` script for generating publication (LaTeX) output. - Added radial velocity methods for binary models +- Support for wideband data in `pint.bayesian` (no correlated noise). - Added `DMWaveX` model (Fourier representation of DM noise) ### Fixed - Wave model `validate()` can correctly use PEPOCH to assign WAVEEPOCH parameter diff --git a/docs/conf.py b/docs/conf.py index 6a0476968..5516fbf3e 100755 --- a/docs/conf.py +++ b/docs/conf.py @@ -119,6 +119,8 @@ "examples/compare_tempo2_B1855.py", "examples/example_dmx_ranges.py", "examples/example_pulse_numbers.py", + # "examples/bayesian-example-NGC6440E.py", + # "examples/bayesian-wideband-example.py", "conf.py", "_ext", ] diff --git a/docs/examples/bayesian-example-NGC6440E.py b/docs/examples/bayesian-example-NGC6440E.py index daf571e8f..7b1952736 100644 --- a/docs/examples/bayesian-example-NGC6440E.py +++ b/docs/examples/bayesian-example-NGC6440E.py @@ -1,3 +1,18 @@ +# --- +# jupyter: +# jupytext: +# formats: ipynb,py:percent +# text_representation: +# extension: .py +# format_name: percent +# format_version: '1.3' +# jupytext_version: 1.14.7 +# kernelspec: +# display_name: Python 3 (ipykernel) +# language: python +# name: python3 +# --- + # %% [markdown] # # PINT Bayesian Interface Examples @@ -6,6 +21,7 @@ from pint.bayesian import BayesianTiming from pint.config import examplefile from pint.models.priors import Prior +from pint.logging import setup as setup_log from scipy.stats import uniform # %% @@ -16,6 +32,10 @@ import io import matplotlib.pyplot as plt +# %% +# Turn off log messages. They can slow down the processing. +setup_log(level="WARNING") + # %% # Read the par and tim files parfile = examplefile("NGC6440E.par.good") @@ -70,36 +90,48 @@ + np.random.randn(nwalkers, bt.nparams) * maxlike_errors ) +# %% +# ** IMPORTANT!!! ** +# This is used to exclude some of the following time-consuming steps from the readthedocs build. +# Set this to False while actually using this example. +rtd = True + # %% # Use longer chain_length for real runs. It is kept small here so that # the sampling finishes quickly (and because I know the burn in is short # because of the cheating priors above). -print("Running emcee...") -chain_length = 1000 -sampler.run_mcmc( - start_points, - chain_length, - progress=True, -) +if not rtd: + print("Running emcee...") -# %% -# Merge all the chains together after discarding the first 100 samples as 'burn-in'. -# The burn-in should be decided after looking at the chains in the real world. -samples_emcee = sampler.get_chain(flat=True, discard=100) + chain_length = 1000 + + sampler.run_mcmc( + start_points, + chain_length, + progress=True, + ) # %% -# Plot the MCMC chains to make sure that the burn-in has been removed properly. -# Otherwise, go back and discard more points. -for idx, param_chain in enumerate(samples_emcee.T): - plt.subplot(bt.nparams, 1, idx + 1) - plt.plot(param_chain, label=bt.param_labels[idx]) - plt.legend() -plt.show() +if not rtd: + # Merge all the chains together after discarding the first 100 samples as 'burn-in'. + # The burn-in should be decided after looking at the chains in the real world. + samples_emcee = sampler.get_chain(flat=True, discard=100) + + # %% +if not rtd: + # Plot the MCMC chains to make sure that the burn-in has been removed properly. + # Otherwise, go back and discard more points. + for idx, param_chain in enumerate(samples_emcee.T): + plt.subplot(bt.nparams, 1, idx + 1) + plt.plot(param_chain, label=bt.param_labels[idx]) + plt.legend() + plt.show() # %% # Plot the posterior distribution. -fig = corner.corner(samples_emcee, labels=bt.param_labels) -plt.show() +if not rtd: + fig = corner.corner(samples_emcee, labels=bt.param_labels) + plt.show() # %% [markdown] # ## Nested sampling with nestle @@ -118,28 +150,30 @@ # `dlogz` is the target accuracy in the computed Bayesian evidence. # Increasing `npoints` or decreasing `dlogz` gives more accurate results, # but at the cost of time. -print("Running nestle...") -result_nestle_1 = nestle.sample( - bt.lnlikelihood, - bt.prior_transform, - bt.nparams, - method="multi", - npoints=150, - dlogz=0.5, - callback=nestle.print_progress, -) +if not rtd: + print("Running nestle...") + result_nestle_1 = nestle.sample( + bt.lnlikelihood, + bt.prior_transform, + bt.nparams, + method="multi", + npoints=150, + dlogz=0.5, + callback=nestle.print_progress, + ) # %% # Plot the posterior # The nested samples come with weights, which must be taken into account # while plotting. -fig = corner.corner( - result_nestle_1.samples, - weights=result_nestle_1.weights, - labels=bt.param_labels, - range=[0.999] * bt.nparams, -) -plt.show() +if not rtd: + fig = corner.corner( + result_nestle_1.samples, + weights=result_nestle_1.weights, + labels=bt.param_labels, + range=[0.999] * bt.nparams, + ) + plt.show() # %% [markdown] # Let us create a new model with an EFAC applied to all toas (all @@ -170,36 +204,41 @@ print(bt2.likelihood_method) # %% -result_nestle_2 = nestle.sample( - bt2.lnlikelihood, - bt2.prior_transform, - bt2.nparams, - method="multi", - npoints=150, - dlogz=0.5, - callback=nestle.print_progress, -) +if not rtd: + result_nestle_2 = nestle.sample( + bt2.lnlikelihood, + bt2.prior_transform, + bt2.nparams, + method="multi", + npoints=150, + dlogz=0.5, + callback=nestle.print_progress, + ) # %% # Plot the posterior. # The EFAC looks consistent with 1. -fig2 = corner.corner( - result_nestle_2.samples, - weights=result_nestle_2.weights, - labels=bt2.param_labels, - range=[0.999] * bt2.nparams, -) -plt.show() +if not rtd: + fig2 = corner.corner( + result_nestle_2.samples, + weights=result_nestle_2.weights, + labels=bt2.param_labels, + range=[0.999] * bt2.nparams, + ) + plt.show() # %% [markdown] # Now let us look at the evidences and compute the Bayes factor. # %% -print(f"Evidence without EFAC : {result_nestle_1.logz} +/- {result_nestle_1.logzerr}") -print(f"Evidence with EFAC : {result_nestle_2.logz} +/- {result_nestle_2.logzerr}") +if not rtd: + print( + f"Evidence without EFAC : {result_nestle_1.logz} +/- {result_nestle_1.logzerr}" + ) + print(f"Evidence with EFAC : {result_nestle_2.logz} +/- {result_nestle_2.logzerr}") -bf = np.exp(result_nestle_1.logz - result_nestle_2.logz) -print(f"Bayes factor : {bf} (in favor of no EFAC)") + bf = np.exp(result_nestle_1.logz - result_nestle_2.logz) + print(f"Bayes factor : {bf} (in favor of no EFAC)") # %% [markdown] # The Bayes factor tells us that the EFAC is unnecessary for this dataset. diff --git a/docs/examples/bayesian-wideband-example.py b/docs/examples/bayesian-wideband-example.py new file mode 100644 index 000000000..6ab04bb32 --- /dev/null +++ b/docs/examples/bayesian-wideband-example.py @@ -0,0 +1,110 @@ +# %% [markdown] +# # PINT Bayesian Interface Example (Wideband) + +# %% +import corner +import emcee +import matplotlib.pyplot as plt +import numpy as np + +from pint.bayesian import BayesianTiming +from pint.config import examplefile +from pint.fitter import WidebandDownhillFitter +from pint.logging import setup as setup_log +from pint.models import get_model_and_toas + +# %% +# Turn off log messages. They can slow down the processing. +setup_log(level="WARNING") + +# %% +# This is a simulated dataset. +m, t = get_model_and_toas(examplefile("test-wb-0.par"), examplefile("test-wb-0.tim")) + +# %% +# Fit the model to the data to get the parameter uncertainties. +ftr = WidebandDownhillFitter(t, m) +ftr.fit_toas() +m = ftr.model + +# %% +# Now set the priors. +# I am cheating here by setting the priors around the maximum likelihood estimates. +# This is a bad idea for real datasets and can bias the estimates. I am doing this +# here just to make everything finish faster. In the real world, these priors should +# be informed by, e.g. previous (independent) timing solutions, pulsar search results, +# VLBI localization etc. Note that unbounded uniform priors don't work here. +prior_info = {} +for par in m.free_params: + param = getattr(m, par) + param_min = float(param.value - 10 * param.uncertainty_value) + param_max = float(param.value + 10 * param.uncertainty_value) + prior_info[par] = {"distr": "uniform", "pmin": param_min, "pmax": param_max} + +# %% +# Set the EFAC and DMEFAC priors and unfreeze them. +# Don't do this before the fitting step. The fitter doesn't know +# how to deal with noise parameters. +prior_info["EFAC1"] = {"distr": "normal", "mu": 1, "sigma": 0.1} +prior_info["DMEFAC1"] = {"distr": "normal", "mu": 1, "sigma": 0.1} + +m.EFAC1.frozen = False +m.EFAC1.uncertainty_value = 0.01 +m.DMEFAC1.frozen = False +m.DMEFAC1.uncertainty_value = 0.01 + +# %% +# The likelihood function behaves better if `use_pulse_numbers==True`. +bt = BayesianTiming(m, t, use_pulse_numbers=True, prior_info=prior_info) + +# %% +print("Number of parameters = ", bt.nparams) +print("Likelihood method = ", bt.likelihood_method) + +# %% +nwalkers = 25 +sampler = emcee.EnsembleSampler(nwalkers, bt.nparams, bt.lnposterior) + +# %% +# Start the sampler close to the maximul likelihood estimate. +maxlike_params = np.array([param.value for param in bt.params], dtype=float) +maxlike_errors = [param.uncertainty_value for param in bt.params] +start_points = ( + np.repeat([maxlike_params], nwalkers).reshape(bt.nparams, nwalkers).T + + np.random.randn(nwalkers, bt.nparams) * maxlike_errors +) + +# %% +# ** IMPORTANT!!! ** +# This is used to exclude the following time-consuming steps from the readthedocs build. +# Set this to False while actually using this example. +rtd = True + +# %% +if not rtd: + print("Running emcee...") + chain_length = 1000 + sampler.run_mcmc( + start_points, + chain_length, + progress=True, + ) + + samples_emcee = sampler.get_chain(flat=True, discard=100) + +# %% +# Plot the chains to make sure they have converged and the burn-in has been removed properly. +if not rtd: + for idx, param_chain in enumerate(samples_emcee.T): + plt.subplot(bt.nparams, 1, idx + 1) + plt.plot(param_chain) + plt.ylabel(bt.param_labels[idx]) + plt.autoscale() + plt.show() + +# %% +if not rtd: + fig = corner.corner( + samples_emcee, labels=bt.param_labels, quantiles=[0.5], truths=maxlike_params + ) + plt.show() diff --git a/docs/tutorials.rst b/docs/tutorials.rst index fea519176..06f7237f5 100644 --- a/docs/tutorials.rst +++ b/docs/tutorials.rst @@ -67,6 +67,8 @@ are not included in the default build because they take too long, but you can do examples/check_phase_connection.ipynb examples/PINT_observatories.ipynb examples/solar_wind.ipynb + examples/bayesian-example-NGC6440E.py + examples/bayesian-wideband-example.py examples-rendered/paper_validation_example.ipynb .. _`Time a Pulsar`: examples/time_a_pulsar.html diff --git a/src/pint/bayesian.py b/src/pint/bayesian.py index 0e6c99e3d..c83a3aeca 100644 --- a/src/pint/bayesian.py +++ b/src/pint/bayesian.py @@ -6,7 +6,7 @@ from scipy.stats import norm, uniform from pint.models.priors import Prior, UniformUnboundedRV -from pint.residuals import Residuals +from pint.residuals import Residuals, WidebandTOAResiduals class BayesianTiming: @@ -34,19 +34,32 @@ class BayesianTiming: 1. The `prior` attribute of each free parameter in the `model` object should be set to an instance of :class:`pint.models.priors.Prior`. - 2. The parameters of BayesianTiming.model will change for every likelihood function call. + 2. The parameters of `BayesianTiming.model` will change for every likelihood function call. These parameters in general will not be the best-fit values. Hence, it is NOT a good idea to save it as a par file. - 3. Only narow-band TOAs are supported at present. + 3. Both narrow-band and wide-band TOAs are supported. 4. Currently, only uniform and normal distributions are supported in prior_info. More general priors should be set directly in the TimingModel object before creating the - BayesianTiming object. Here is an example prior_info object: - - `prior_info = { "F0" : {"distr" : "normal", "mu" : 1, "sigma" : 0.00001}, "EFAC1" : {"distr" : "uniform", "pmin" : 0.5, "pmax" : 2.0} }` - - See examples/bayesian-example-NGC6440E.py for detailed example. + BayesianTiming object. Here is an example prior_info object:: + + ``` + prior_info = { + "F0" : { + "distr" : "normal", + "mu" : 1, + "sigma" : 0.00001 + }, + "EFAC1" : { + "distr" : "uniform", + "pmin" : 0.5, + "pmax" : 2.0 + } + } + ``` + + See `examples/bayesian-example-NGC6440E.py` and `examples/bayesian-wideband-example` for detailed examples. """ def __init__(self, model, toas, use_pulse_numbers=False, prior_info=None): @@ -54,8 +67,12 @@ def __init__(self, model, toas, use_pulse_numbers=False, prior_info=None): self.model = deepcopy(model) self.toas = toas - if toas.is_wideband(): - raise NotImplementedError("Wideband TOAs are not yet supported.") + if use_pulse_numbers: + self.toas.compute_pulse_numbers(self.model) + + self.track_mode = "use_pulse_numbers" if use_pulse_numbers else "nearest" + + self.is_wideband = toas.is_wideband() self.param_labels = self.model.free_params self.params = [getattr(self.model, par) for par in self.param_labels] @@ -79,8 +96,6 @@ def __init__(self, model, toas, use_pulse_numbers=False, prior_info=None): self.likelihood_method = self._decide_likelihood_method() - self.track_mode = "use_pulse_numbers" if use_pulse_numbers else "nearest" - def _validate_priors(self): for param in self.params: if not hasattr(param, "prior") or param.prior is None: @@ -91,8 +106,9 @@ def _validate_priors(self): ) def _decide_likelihood_method(self): - """Weighted least squares with normalization term (wls), or Generalized - least-squares with normalization term (gls).""" + """Weighted least squares with normalization term (wls), or Generalized least + squares with normalization term (gls), for narrow-band (nb) or wide-band (wb) + dataset.""" if "NoiseComponent" not in self.model.component_types: return "wls" @@ -104,7 +120,6 @@ def _decide_likelihood_method(self): ) else: return "wls" - # return "gls" def lnprior(self, params): """Basic implementation of a factorized log prior. @@ -118,7 +133,8 @@ def lnprior(self, params): """ if len(params) != self.nparams: raise IndexError( - f"The number of input parameters ({len(params)}) should be the same as the number of free parameters ({self.nparams})." + f"The number of input parameters ({len(params)}) should be the same " + f"as the number of free parameters ({self.nparams})." ) lnsum = 0.0 @@ -158,7 +174,11 @@ def lnlikelihood(self, params): float: The value of the log-likelihood at params """ if self.likelihood_method == "wls": - return self._wls_lnlikelihood(params) + return ( + self._wls_wb_lnlikelihood(params) + if self.is_wideband + else self._wls_nb_lnlikelihood(params) + ) elif self.likelihood_method == "gls": raise NotImplementedError( "GLS likelihood for correlated noise is not yet implemented." @@ -179,16 +199,19 @@ def lnposterior(self, params): lnpr = self.lnprior(params) return lnpr + self.lnlikelihood(params) if np.isfinite(lnpr) else -np.inf - def _wls_lnlikelihood(self, params): - """Implementation of Log-Likelihood function for uncorrelated noise only. - `wls' stands for weighted least squares. Also includes the normalization - term to enable sampling over white noise parameters (EFAC and EQUAD). + def _wls_nb_lnlikelihood(self, params): + """Implementation of Log-Likelihood function for uncorrelated noise only for + narrow-band TOAs. `wls' stands for weighted least squares. Also includes the + normalization term to enable sampling over white noise parameters (EFAC and + EQUAD). Args: - params (array-like): Parameters + params : (array-like) + Parameters Returns: - float: The value of the log-likelihood at params + float : + The value of the log-likelihood at params """ params_dict = dict(zip(self.param_labels, params)) self.model.set_param_values(params_dict) @@ -196,3 +219,34 @@ def _wls_lnlikelihood(self, params): chi2 = res.calc_chi2() sigmas = self.model.scaled_toa_uncertainty(self.toas).si.value return -chi2 / 2 - np.sum(np.log(sigmas)) + + def _wls_wb_lnlikelihood(self, params): + """Implementation of Log-Likelihood function for uncorrelated noise only for + wide-band TOAs. `wls' stands for weighted least squares. Also includes the + normalization terms to enable sampling over white noise parameters (EFAC, EQUAD, + DMEFAC and DMEQUAD). + + Args: + params : (array-like) + Parameters + + Returns: + float : + The value of the log-likelihood at params + """ + params_dict = dict(zip(self.param_labels, params)) + self.model.set_param_values(params_dict) + + res = WidebandTOAResiduals( + self.toas, self.model, toa_resid_args={"track_mode": self.track_mode} + ) + + chi2_toa = res.toa.calc_chi2() + sigmas_toa = self.model.scaled_toa_uncertainty(self.toas).si.value + lnL_toa = -chi2_toa / 2 - np.sum(np.log(sigmas_toa)) + + chi2_dm = res.dm.calc_chi2() + sigmas_dm = self.model.scaled_dm_uncertainty(self.toas).si.value + lnL_dm = -chi2_dm / 2 - np.sum(np.log(sigmas_dm)) + + return lnL_toa + lnL_dm diff --git a/src/pint/data/examples/test-wb-0.par b/src/pint/data/examples/test-wb-0.par new file mode 100644 index 000000000..568e9b5b1 --- /dev/null +++ b/src/pint/data/examples/test-wb-0.par @@ -0,0 +1,31 @@ +# Created: 2023-08-18T12:08:17.544471 +# PINT_version: 0.9.6+119.g22716c49 +# User: Abhimanyu Susobhanan (abhimanyu) +# Host: abhimanyu-VirtualBox +# OS: Linux-6.2.0-26-generic-x86_64-with-glibc2.35 +# Python: 3.10.12 (main, Jul 5 2023, 18:54:27) [GCC 11.2.0] +# Format: pint +PSRJ WBTEST +EPHEM DE440 +DILATEFREQ N +DMDATA N +NTOA 0 +CHI2 0.0 +RAJ 4:00:00.00000000 1 0.00000000000000000000 +DECJ 15:00:00.00000000 1 0.00000000000000000000 +PMRA 0.0 +PMDEC 0.0 +PX 0.0 +POSEPOCH 55000.0000000000000000 +F0 100.0 1 0.0 +F1 1e-15 1 0.0 +PEPOCH 55000.0000000000000000 +PLANET_SHAPIRO N +DM 10.0 1 0.0 +DM1 0.0001 1 0.0 +DMEPOCH 55000.0000000000000000 +TZRMJD 55005.0251256281408097 +TZRSITE ssb +TZRFRQ inf +EFAC tel gmrt 1 +DMEFAC tel gmrt 1 \ No newline at end of file diff --git a/src/pint/data/examples/test-wb-0.tim b/src/pint/data/examples/test-wb-0.tim new file mode 100644 index 000000000..a7d3d52c1 --- /dev/null +++ b/src/pint/data/examples/test-wb-0.tim @@ -0,0 +1,207 @@ +FORMAT 1 +C Created: 2023-08-18T12:08:45.036835 +C PINT_version: 0.9.6+119.g22716c49 +C User: Abhimanyu Susobhanan (abhimanyu) +C Host: abhimanyu-VirtualBox +C OS: Linux-6.2.0-26-generic-x86_64-with-glibc2.35 +C Python: 3.10.12 (main, Jul 5 2023, 18:54:27) [GCC 11.2.0] +fake 400.000000 54000.0000000338171065 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999689348312622012 -pn -8683387108.0 +fake 800.000000 54010.0502512459059028 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999817394902294414 -pn -8596545796.0 +fake 1200.000000 54020.1005025362641782 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999891183572727027 -pn -8509705405.0 +fake 400.000000 54030.1507537448405324 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999631939252680306 -pn -8422866134.0 +fake 800.000000 54040.2010049845390741 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9995842407658120776 -pn -8336028086.0 +fake 1200.000000 54050.2512562274484143 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999814299821878541 -pn -8249191442.0 +fake 400.000000 54060.3015075897839931 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9996197994101480066 -pn -8162356296.0 +fake 800.000000 54070.3517587435774421 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999725184871209322 -pn -8075522630.0 +fake 1200.000000 54080.4020101060183797 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999767790770363772 -pn -7988690482.0 +fake 400.000000 54090.4522613409235186 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999686831665311229 -pn -7901859813.0 +fake 800.000000 54100.5025125235372685 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999841888231202267 -pn -7815030452.0 +fake 1200.000000 54110.5527637856061109 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999855240286045288 -pn -7728202298.0 +fake 400.000000 54120.6030150412243056 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999827872380426933 -pn -7641375182.0 +fake 800.000000 54130.6532663832630324 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999810559785563071 -pn -7554548810.0 +fake 1200.000000 54140.7035175807754050 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999944845353083878 -pn -7467722991.0 +fake 400.000000 54150.7537688665358102 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999636388273509543 -pn -7380897481.0 +fake 800.000000 54160.8040201505214583 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999806942922156322 -pn -7294071936.0 +fake 1200.000000 54170.8542713282116436 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999794849102732609 -pn -7207246151.0 +fake 400.000000 54180.9045226137458797 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999755115559108437 -pn -7120419888.0 +fake 800.000000 54190.9547738557419675 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999909233122702556 -pn -7033592836.0 +fake 1200.000000 54201.0050251104141204 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999546290756704471 -pn -6946764852.0 +fake 400.000000 54211.0552763242439815 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999844631460155727 -pn -6859935773.0 +fake 800.000000 54221.1055276941583796 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999672173798245581 -pn -6773105381.0 +fake 1200.000000 54231.1557789516169097 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999734640037218725 -pn -6686273645.0 +fake 400.000000 54241.2060301601370255 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999652091265824406 -pn -6599440510.0 +fake 800.000000 54251.2562814636363426 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999855249592944891 -pn -6512605883.0 +fake 1200.000000 54261.3065327059755903 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9996220393044169695 -pn -6425769852.0 +fake 400.000000 54271.3567838821766782 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999763661348765612 -pn -6338932478.0 +fake 800.000000 54281.4070351975145023 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999840281158982123 -pn -6252093789.0 +fake 1200.000000 54291.4572863853987268 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999763278947545926 -pn -6165253979.0 +fake 400.000000 54301.5075376360473032 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999901200977491657 -pn -6078413208.0 +fake 800.000000 54311.5577889072436342 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999719461962584742 -pn -5991571602.0 +fake 1200.000000 54321.6080401637861922 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999673262590097732 -pn -5904729427.0 +fake 400.000000 54331.6582915012324653 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999966387002557102 -pn -5817886914.0 +fake 800.000000 54341.7085426773096991 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999787029487756375 -pn -5731044244.0 +fake 1200.000000 54351.7587939677584142 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999932454784072644 -pn -5644201711.0 +fake 400.000000 54361.8090452441937731 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999987919231125403 -pn -5557359572.0 +fake 800.000000 54371.8592964550520139 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999755049704442288 -pn -5470518003.0 +fake 1200.000000 54381.9095477115447685 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999640428002742632 -pn -5383677276.0 +fake 400.000000 54391.9597989894186458 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9998219673397774525 -pn -5296837611.0 +fake 800.000000 54402.0100502489624421 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999871983091432774 -pn -5209999116.0 +fake 1200.000000 54412.0603015108758565 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.00000703329270763 -pn -5123161982.0 +fake 400.000000 54422.1105527431482870 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999884527228331438 -pn -5036326330.0 +fake 800.000000 54432.1608040676221759 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999816662373855615 -pn -4949492142.0 +fake 1200.000000 54442.2110552311124421 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.99964302668641098 -pn -4862659484.0 +fake 400.000000 54452.2613065788043287 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9998872578856580476 -pn -4775828332.0 +fake 800.000000 54462.3115577385225695 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999944859944271937 -pn -4688998523.0 +fake 1200.000000 54472.3618090102764815 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999879509734557669 -pn -4602169981.0 +fake 400.000000 54482.4120603123242361 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999995082747862556 -pn -4515342546.0 +fake 800.000000 54492.4623115789709374 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999816373048808085 -pn -4428515929.0 +fake 1200.000000 54502.5125628608401852 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999774228300115434 -pn -4341689957.0 +fake 400.000000 54512.5628140679162617 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9997085408988969114 -pn -4254864383.0 +fake 800.000000 54522.6130652806756481 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9999465770119952785 -pn -4168038864.0 +fake 1200.000000 54532.6633165776520833 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999716794891414485 -pn -4081213200.0 +fake 400.000000 54542.7135678561986690 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999874148552329637 -pn -3994387138.0 +fake 800.000000 54552.7638191202778934 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999771394414748284 -pn -3907560367.0 +fake 1200.000000 54562.8140703119680556 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999811388630607196 -pn -3820732736.0 +fake 400.000000 54572.8643216100541089 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9998608980168683345 -pn -3733904058.0 +fake 800.000000 54582.9145728791742246 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999906682211421074 -pn -3647074118.0 +fake 1200.000000 54592.9648241367796759 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999895041223835885 -pn -3560242862.0 +fake 400.000000 54603.0150753694817014 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9998716615977875505 -pn -3473410218.0 +fake 800.000000 54613.0653266786111111 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000010231048680057 -pn -3386576091.0 +fake 1200.000000 54623.1155778719024421 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000004776792867532 -pn -3299740544.0 +fake 400.000000 54633.1658291538335648 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.99967172207249051 -pn -3212903628.0 +fake 800.000000 54643.2160804013887037 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000108122754010403 -pn -3126065367.0 +fake 1200.000000 54653.2663316774607060 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999980552841559337 -pn -3039225929.0 +fake 400.000000 54663.3165828598465394 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000068670936557959 -pn -2952385478.0 +fake 800.000000 54673.3668341711507523 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999679371789044481 -pn -2865544126.0 +fake 1200.000000 54683.4170853739764931 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999990687399299406 -pn -2778702126.0 +fake 400.000000 54693.4673367320890972 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9999720094319615475 -pn -2691859715.0 +fake 800.000000 54703.5175879366010185 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.00005208609530464 -pn -2605017061.0 +fake 1200.000000 54713.5678392190929050 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999987574852184185 -pn -2518174457.0 +fake 400.000000 54723.6180904375442015 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.99989236262914892 -pn -2431332170.0 +fake 800.000000 54733.6683416922952431 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9997956601147193045 -pn -2344490365.0 +fake 1200.000000 54743.7185929525367246 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999866982796133246 -pn -2257649329.0 +fake 400.000000 54753.7688441728829630 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999859586492998137 -pn -2170809292.0 +fake 800.000000 54763.8190955198084374 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9998250224034890485 -pn -2083970360.0 +fake 1200.000000 54773.8693467226388774 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0000284739417204605 -pn -1997132752.0 +fake 400.000000 54783.9195980059268866 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0002369037720472405 -pn -1910296595.0 +fake 800.000000 54793.9698491886990046 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000065942886162499 -pn -1823461885.0 +fake 1200.000000 54804.0201005511857755 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999944891336582828 -pn -1736628712.0 +fake 400.000000 54814.0703517722123264 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999927984918873445 -pn -1649797063.0 +fake 800.000000 54824.1206030679627083 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999787109292992464 -pn -1562966790.0 +fake 1200.000000 54834.1708543070402546 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000085109166870542 -pn -1476137742.0 +fake 400.000000 54844.2211055407940278 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999892553417434964 -pn -1389309959.0 +fake 800.000000 54854.2713567972613657 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000030425752240226 -pn -1302483073.0 +fake 1200.000000 54864.3216080450358565 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000126721528663949 -pn -1215656918.0 +fake 400.000000 54874.3718593442969676 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999967224511925204 -pn -1128831243.0 +fake 800.000000 54884.4221106018054630 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999973596921750957 -pn -1042005722.0 +fake 1200.000000 54894.4723617696319560 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999859830658437138 -pn -955180146.0 +fake 400.000000 54904.5226130309290625 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999938186615571616 -pn -868354255.0 +fake 800.000000 54914.5728643355897106 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9998099117997843585 -pn -781527742.0 +fake 1200.000000 54924.6231155570151967 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999933195919728759 -pn -694700435.0 +fake 400.000000 54934.6733667954781017 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999879310638236722 -pn -607872144.0 +fake 800.000000 54944.7236180665839236 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.00000920829822565 -pn -521042643.0 +fake 1200.000000 54954.7738693193937037 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999937597913816637 -pn -434211854.0 +fake 400.000000 54964.8241205821242014 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999776063534331207 -pn -347379703.0 +fake 800.000000 54974.8743718712418403 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000173321152803697 -pn -260546076.0 +fake 1200.000000 54984.9246230913911458 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000124638396050736 -pn -173711016.0 +fake 400.000000 54994.9748743358192478 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000005764303560627 -pn -86874574.0 +fake 800.000000 55005.0251256172320255 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999930822646455541 -pn -36749.0 +fake 1200.000000 55015.0753769322097917 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999909343956333098 -pn 86802299.0 +fake 400.000000 55025.1256281353382986 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999997018102579935 -pn 173642409.0 +fake 800.000000 55035.1758794046910995 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999960011834497916 -pn 260483490.0 +fake 1200.000000 55045.2261306872286574 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000095822229545901 -pn 347325288.0 +fake 400.000000 55055.2763818888063657 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999909720239615588 -pn 434167570.0 +fake 800.000000 55065.3266331769281134 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000019166610654427 -pn 521010185.0 +fake 1200.000000 55075.3768844260846644 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999921504496390928 -pn 607852826.0 +fake 400.000000 55085.4271357309605325 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9999044840334824925 -pn 694695237.0 +fake 800.000000 55095.4773869083597221 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999911289996185105 -pn 781537250.0 +fake 1200.000000 55105.5276381787532987 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000066185921840471 -pn 868378565.0 +fake 400.000000 55115.5778894215778819 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999946260456715445 -pn 955218957.0 +fake 800.000000 55125.6281406752153588 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.99994596279539534 -pn 1042058303.0 +fake 1200.000000 55135.6783919489981597 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000072974981106686 -pn 1128896371.0 +fake 400.000000 55145.7286431653392593 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000047762124557193 -pn 1215733029.0 +fake 800.000000 55155.7788944329273264 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9999371452968319925 -pn 1302568256.0 +fake 1200.000000 55165.8291457748755208 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999967032042070764 -pn 1389401948.0 +fake 400.000000 55175.8793970056304398 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000127460559598744 -pn 1476234108.0 +fake 800.000000 55185.9296482781195716 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999949353445512128 -pn 1563064854.0 +fake 1200.000000 55195.9798995514579283 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000038469599295368 -pn 1649894233.0 +fake 400.000000 55206.0301507896968171 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000184488525271191 -pn 1736722388.0 +fake 800.000000 55216.0804020083592593 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000020397349509981 -pn 1823549567.0 +fake 1200.000000 55226.1306532466018171 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999832160375013049 -pn 1910375940.0 +fake 400.000000 55236.1809045240011805 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000220882169742772 -pn 1997201742.0 +fake 800.000000 55246.2311557515021991 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0000712212124736955 -pn 2084027296.0 +fake 1200.000000 55256.2814070796484143 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.99988297755272494 -pn 2170852822.0 +fake 400.000000 55266.3316582414325116 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000027401418097545 -pn 2257678566.0 +fake 800.000000 55276.3819095680059607 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000059497362002379 -pn 2344504852.0 +fake 1200.000000 55286.4321607618225463 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000242234020140315 -pn 2431331861.0 +fake 400.000000 55296.4824120603521874 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000334406591437854 -pn 2518159783.0 +fake 800.000000 55306.5326633376824653 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999886245295953048 -pn 2604988867.0 +fake 1200.000000 55316.5829146155144906 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000059986473515738 -pn 2691819199.0 +fake 400.000000 55326.6331658011917592 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000066715297739709 -pn 2778650861.0 +fake 800.000000 55336.6834170336900000 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000088876729311907 -pn 2865483994.0 +fake 1200.000000 55346.7336683820114005 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0001555801232671445 -pn 2952318559.0 +fake 400.000000 55356.7839196332011459 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999993792149525078 -pn 3039154520.0 +fake 800.000000 55366.8341708451521991 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0000175345260493365 -pn 3125991898.0 +fake 1200.000000 55376.8844221036844907 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999927402643881202 -pn 3212830533.0 +fake 400.000000 55386.9346733113793634 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000126990980531342 -pn 3299670286.0 +fake 800.000000 55396.9849246155898148 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.9999658517649233766 -pn 3386511074.0 +fake 1200.000000 55407.0351759348490509 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999947721079837201 -pn 3473352642.0 +fake 400.000000 55417.0854271475039583 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000089478169504019 -pn 3560194776.0 +fake 800.000000 55427.1356784344858333 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000100905303029587 -pn 3647037321.0 +fake 1200.000000 55437.1859296868953009 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0001193380051912945 -pn 3733879972.0 +fake 400.000000 55447.2361809303057870 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000020504533308162 -pn 3820722483.0 +fake 800.000000 55457.2864321273571759 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000068139775406002 -pn 3907564673.0 +fake 1200.000000 55467.3366834195775347 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000128256706893354 -pn 3994406245.0 +fake 400.000000 55477.3869346802707523 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000173487726014024 -pn 4081246971.0 +fake 800.000000 55487.4371859654377315 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999922519302911486 -pn 4168086708.0 +fake 1200.000000 55497.4874372176980671 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000077837822342276 -pn 4254925225.0 +fake 400.000000 55507.5376884554523611 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000207002326607559 -pn 4341762373.0 +fake 800.000000 55517.5879397130096992 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000239330535125347 -pn 4428598110.0 +fake 1200.000000 55527.6381908988099653 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000162814231407506 -pn 4515432327.0 +fake 400.000000 55537.6884421805568172 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000120149866405688 -pn 4602265002.0 +fake 800.000000 55547.7386934888747339 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000237508928458069 -pn 4689096237.0 +fake 1200.000000 55557.7889446945634028 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000141186570696267 -pn 4775926069.0 +fake 400.000000 55567.8391959928024653 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0002246170119403965 -pn 4862754615.0 +fake 800.000000 55577.8894472530937152 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000219027074768849 -pn 4949582119.0 +fake 1200.000000 55587.9396984947352199 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000268042086678121 -pn 5036408740.0 +fake 400.000000 55597.9899497214191552 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000040338675025169 -pn 5123234696.0 +fake 800.000000 55608.0402009677436111 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000263667756289297 -pn 5210060320.0 +fake 1200.000000 55618.0904522558979745 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000257447578972783 -pn 5296885822.0 +fake 400.000000 55628.1407035277368403 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000263034605704899 -pn 5383711450.0 +fake 800.000000 55638.1909548204844676 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000185262083079626 -pn 5470537541.0 +fake 1200.000000 55648.2412060724572569 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000217004341520001 -pn 5557364273.0 +fake 400.000000 55658.2914572935849074 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000134242018416416 -pn 5644191850.0 +fake 800.000000 55668.3417085787174190 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000200539617497755 -pn 5731020538.0 +fake 1200.000000 55678.3919597921981829 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000234690861819577 -pn 5817850422.0 +fake 400.000000 55688.4422110333989121 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000117595413135451 -pn 5904681611.0 +fake 800.000000 55698.4924623667830324 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000242095632078261 -pn 5991514255.0 +fake 1200.000000 55708.5427135788948612 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000156527178478722 -pn 6078348320.0 +fake 400.000000 55718.5929648814338888 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000298488753902373 -pn 6165183800.0 +fake 800.000000 55728.6432161318090741 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.00020249294367174 -pn 6252020717.0 +fake 1200.000000 55738.6934673512211690 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0004272154687253536 -pn 6338858925.0 +fake 400.000000 55748.7437186170883681 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000058013579587481 -pn 6425698306.0 +fake 800.000000 55758.7939698085431249 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000378942298283491 -pn 6512538772.0 +fake 1200.000000 55768.8442211370403588 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000301078007239742 -pn 6599380089.0 +fake 400.000000 55778.8944723171171297 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000222538928678034 -pn 6686222049.0 +fake 800.000000 55788.9447236245216435 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000197672546348252 -pn 6773064493.0 +fake 1200.000000 55798.9949748948978704 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000208636614789991 -pn 6859907132.0 +fake 400.000000 55809.0452261702578357 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000151785477645361 -pn 6946749714.0 +fake 800.000000 55819.0954773940804167 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000344632926116159 -pn 7033592054.0 +fake 1200.000000 55829.1457285995028009 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0003719383921038395 -pn 7120433862.0 +fake 400.000000 55839.1959798529702083 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000278484473008056 -pn 7207274894.0 +fake 800.000000 55849.2462311071718403 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000126392667082668 -pn 7294115004.0 +fake 1200.000000 55859.2964824622073033 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000133314170880341 -pn 7380953957.0 +fake 400.000000 55869.3467336315645370 1.000 gmrt -pp_dme 0.0001 -pp_dm 9.999932149006183948 -pn 7467791577.0 +fake 800.000000 55879.3969849065807871 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000301782830795302 -pn 7554627822.0 +fake 1200.000000 55889.4472361936712269 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000220044521693582 -pn 7641462564.0 +fake 400.000000 55899.4974873990938426 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0000454418055698845 -pn 7728295755.0 +fake 800.000000 55909.5477386928395254 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000062125494163085 -pn 7815127495.0 +fake 1200.000000 55919.5979899230785996 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000270080431560062 -pn 7901957794.0 +fake 400.000000 55929.6482412396903356 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0002914779061594 -pn 7988786754.0 +fake 800.000000 55939.6984924249845602 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.00028655900684835 -pn 8075614613.0 +fake 1200.000000 55949.7487436704994676 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0003792051184662885 -pn 8162441507.0 +fake 400.000000 55959.7989949440956827 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000325634328277777 -pn 8249267654.0 +fake 800.000000 55969.8492462052852083 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000326622727975257 -pn 8336093382.0 +fake 1200.000000 55979.8994974572198611 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000160297267337141 -pn 8422918888.0 +fake 400.000000 55989.9497487285055439 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.000299306219806401 -pn 8509744435.0 +fake 800.000000 56000.0000000023089352 1.000 gmrt -pp_dme 0.0001 -pp_dm 10.0003150571158367504 -pn 8596570357.0 diff --git a/src/pint/utils.py b/src/pint/utils.py index 785ae36be..c4410618a 100644 --- a/src/pint/utils.py +++ b/src/pint/utils.py @@ -2354,8 +2354,8 @@ def compute_hash(filename): cryptographically robust. It uses the SHA256 algorithm, which is known to be vulnerable to a length-extension attack. - Parameter - --------- + Parameters + ---------- f : str or Path or file-like The source of input. If file-like, it should return ``bytes`` not ``str`` - that is, the file should be opened in binary mode. diff --git a/tests/test_bayesian.py b/tests/test_bayesian.py index 816723295..4487984ed 100644 --- a/tests/test_bayesian.py +++ b/tests/test_bayesian.py @@ -178,10 +178,24 @@ def test_prior_dict(data_NGC6440E_efac): bt = BayesianTiming(model, toas, use_pulse_numbers=True, prior_info=prior_info) -def test_wideband_exception(data_J0740p6620_wb): +def test_wideband_data(data_J0740p6620_wb): model, toas = data_J0740p6620_wb - with pytest.raises(NotImplementedError): - bt = BayesianTiming(model, toas) + bt = BayesianTiming(model, toas) + + assert bt.is_wideband and bt.likelihood_method == "wls" + + test_cube = 0.5 * np.ones(bt.nparams) + test_params = bt.prior_transform(test_cube) + assert np.all(np.isfinite(test_params)) + + lnpr = bt.lnprior(test_params) + assert np.isfinite(lnpr) + + lnl = bt.lnlikelihood(test_params) + assert np.isfinite(lnl) + + lnp = bt.lnposterior(test_params) + assert np.isfinite(lnp) and np.isclose(lnp, lnpr + lnl) def test_gls_exception(data_NGC6440E, data_NGC6440E_red): @@ -189,18 +203,6 @@ def test_gls_exception(data_NGC6440E, data_NGC6440E_red): with pytest.raises(NotImplementedError): bt = BayesianTiming(model, toas) - model, toas = data_NGC6440E - bt = BayesianTiming(model, toas) - bt.likelihood_method = "gls" - test_cube = 0.5 * np.ones(bt.nparams) - test_params = bt.prior_transform(test_cube) - with pytest.raises(NotImplementedError): - bt.lnlikelihood(test_params) - - bt.likelihood_method = "bla" - with pytest.raises(ValueError): - bt.lnlikelihood(test_params) - def test_badprior_exception(data_NGC6440E_default_priors): model, toas = data_NGC6440E_default_priors