#!/usr/bin/env python3

import sys
import numpy as np
import phoebe
from astropy import units

#logger = phoebe.logger(clevel='INFO')
#logger = phoebe.logger(clevel='DEBUG')

b = phoebe.default_binary()

times, u, v, wavelengths, vises, sigmas = np.loadtxt("Vis.dat", usecols=[0, 1, 2, 3, 5, 6], unpack=True)

b.add_dataset('vis', times=times, u=u, v=v, wavelengths=wavelengths, vises=vises, sigmas=sigmas, if_method='integrate', passband='Johnson:R')

b.set_value('distance', context = 'system', value=100*units.pc)

f = open('twigs.txt0', 'w')
for twig in b.twigs:
  f.write("%s\n" % (twig))
f.close()

#b.run_compute()
b.run_compute(if_method='simple')

f = open('twigs.txt', 'w')
for twig in b.twigs:
  f.write("%s\n" % (twig))
f.close()

# Note: Model (synthetic) baselines and wavelengths must be copied ex-post!
#b.set_value('u@latest@model', value=b.get_value('u@dataset'), ignore_readonly=True)
#b.set_value('v@latest@model', value=b.get_value('v@dataset'), ignore_readonly=True)
#b.set_value('wavelengths@latest@model', value=b.get_value('wavelengths@dataset'), ignore_readonly=True)

#print("b['times@vis01@vis@dataset'] = ", b['times@vis01@vis@dataset'])
#print("b['u@vis01@vis@dataset'] = ", b['u@vis01@vis@dataset'])
#print("b['v@vis01@vis@dataset'] = ", b['v@vis01@vis@dataset'])
#print("b['wavelengths@vis01@vis@dataset'] = ", b['wavelengths@vis01@vis@dataset'])
#print("b['vises@vis01@vis@dataset'] = ", b['vises@vis01@vis@dataset'])
#print("")
#print("b['times@vis01@phoebe01@latest@vis@model'] = ", b['times@vis01@phoebe01@latest@vis@model'])
#print("b['u@vis01@phoebe01@latest@vis@model'] = ", b['u@vis01@phoebe01@latest@vis@model'])
#print("b['v@vis01@phoebe01@latest@vis@model'] = ", b['v@vis01@phoebe01@latest@vis@model'])
#print("b['wavelengths@vis01@phoebe01@latest@vis@model'] = ", b['wavelengths@vis01@phoebe01@latest@vis@model'])
#print("b['vises@vis01@phoebe01@latest@vis@model'] = ", b['vises@vis01@phoebe01@latest@vis@model'])

times = b['times@vis01@phoebe01@latest@vis@model'].value
u = b['u@vis01@phoebe01@latest@vis@model'].value
v = b['v@vis01@phoebe01@latest@vis@model'].value
wavelengths = b['wavelengths@vis01@phoebe01@latest@vis@model'].value
vises = b['vises@vis01@phoebe01@latest@vis@model'].value
#vises = b.get_value('vises', dataset='vis01', context='model')

np.savetxt('test_phoebe3.out', np.c_[times, u, v, wavelengths, vises], header='times u v wavelenghts vises')

#b.plot(show=True)
b.plot(x='u', show=True)
#b.plot(x='v', show=True)