Open kecnry opened 1 year ago
Script to reproduce
import phoebe import numpy as np b = phoebe.default_binary() # set parameter values b.set_value('q', value = 0.6) b.set_value('incl', component='binary', value = 84.5) b.set_value('ecc', 0.2) b.set_value('per0', 63.7) b.set_value('requiv', component='primary', value=1.) b.set_value('requiv', component='secondary', value=0.6) b.set_value('teff', component='secondary', value=5500.) # add an lc dataset b.add_dataset('lc', compute_phases=phoebe.linspace(0,1,101)) b.set_value_all('ld_mode', 'manual') b.set_value_all('ld_mode_bol', 'manual') b.set_value_all('atm', 'blackbody') #compute the model b.run_compute(irrad_method='none') # extract the arrays from the model that we'll use as observables in the next step times = b.get_value('times', context='model', dataset='lc01') # here we're adding noise to the fluxes as well to make the fake data more "realistic" np.random.seed(0) # to ensure reproducibility with added noise fluxes = b.get_value('fluxes', context='model', dataset='lc01') + np.random.normal(size=times.shape) * 0.02 sigmas_lc = np.ones_like(times) * 0.04 b = phoebe.default_binary() b.add_dataset('lc', times=times, fluxes=fluxes, sigmas=sigmas_lc) b.set_value_all('ld_mode', 'manual') b.set_value_all('ld_mode_bol', 'manual') b.set_value_all('atm', 'blackbody') b.add_distribution({'q': phoebe.gaussian_around(0.2), 'incl@binary': phoebe.gaussian_around(7), 'teff@primary': phoebe.gaussian_around(12), 'teff@secondary': phoebe.gaussian_around(200)}, distribution='s1') b.add_solver('sampler.emcee', solver='emcee_solver',overwrite=True) b.set_value('init_from', ['s1'],overwrite=True) b.set_value('compute', solver='emcee_solver', value='phoebe01') b.run_solver('emcee_solver', niters=10, nwalkers=16, solution='emcee_sol',overwrite=True) b.get_parameter(qualifier='q', component='binary', context='component')._latexfmt = 'q' b.get_parameter(qualifier='incl', component='binary', context='component')._latexfmt = 'i' b.get_parameter(qualifier='teff', component='primary', context='component')._latexfmt = 'TT_\mathrm{{1}}' b.get_parameter(qualifier='teff', component='secondary', context='component')._latexfmt = 'TT_\mathrm{{2}}' b.save('test.bundle')
import phoebe b = phoebe.load('test.bundle') afig, mplfig = b.plot(solution='emcee_sol', style='corner', show=True, parameters=['q', 'incl@binary', 'teff@primary', 'teff@secondary'])
See #726
Script to reproduce
See #726