# Configure Jupyter so figures appear in the notebook
%matplotlib inline

# Configure Jupyter to display the assigned value after an assignment
%config InteractiveShell.ast_node_interactivity='last_expr_or_assign'

# import functions from the modsim.py module
from modsim import *

km = UNITS.kilometer
m = UNITS.meter
s = UNITS.second
kg = UNITS.kilogram
days = UNITS.day
degree = UNITS.degree

params = Params(distance = 149.6*10**9 *m,
                v_init = 0 * m/ s,
                G = 6.67 * 10**-11 * (m**3/(s**2 * kg)),
                mass_earth = 5.972 * 10**24* kg,
                mass_sun = 1.989 * 10**30 * kg)
                

params2 = Params(x_distance = 149.6*10**9 *m,
                 y_distance = 0 *m,
                v_init = 0 * m/ s,
                angle = 90 * degree,
                G = 6.67 * 10**-11 * (m**3/(s**2 * kg)),
                mass_earth = 5.972 * 10**24* kg,
                mass_sun = 1.989 * 10**30 * kg)
                

def make_system2(params2):
    """Makes a System object for the given conditions.
    
    params: Params object
    
    returns: System object
    """
    angle, v_init = params2.angle, params2.v_init
    x_distance = params2.x_distance
    y_distance = params2.y_distance
    
    theta = np.deg2rad(angle)
    vx, vy = pol2cart(theta, v_init)
    
    R = Vector(x_distance, y_distance)
    V = Vector(vx, vy)
    
    init = State(R=R, V=V)
    
    t_end =  3.154 *10**7 *s
    dt = t_end / 10**4
    
    return System(params2, init=init, t_end=t_end, dt=dt)

def make_system(params):
    """Makes a System object for the given conditions.
    
    params: Params object
    
    returns: System object
    """
    
    distance = params.distance
    v_init = params.v_init
    
    init = State(r=distance, v=v_init)
    t_end =  9.072* 10**6 *s
    dt = t_end / 10**4
    
    return System(params, init=init, t_end=t_end, dt=dt)

system = make_system2(params2)

def slope_func2(state, t, system):
    
    

    R,V = state
    mass_sun = system.mass_sun
    G = system.G
    r_mag = R.mag
    r_unit_vector = R.hat()
    
    drdt = V
    dvdt = ((-G* mass_sun)/(r_mag**2))*r_unit_vector
    
    return drdt, dvdt

def slope_func(state, t, system):
    """Compute derivatives of the state.
    
    state: position, velocity
    t: time
    system: System object
    
    returns: derivatives of y and v
    """
    r, v = state
    mass_sun = system.mass_sun
    G = system.G
    
    drdt = v
    dvdt = (-G* mass_sun)/(r**2)
    
    return drdt, dvdt

slope_func2(system.init, 0, system)

def event_func(state, t, system):

    R, V = state
    surfaces_touching = R.mag - ((6.371 *10**6 * m) + 695.51 *10**6 * m)
    
    return surfaces_touching

results, details = run_ode_solver(system, slope_func2, events=event_func)
details

results.head()

results.tail()

x = results.R.extract('x')
y = results.R.extract('y')

plot(x,y, label='trajectory')