"""These solutions provide a minimalist answer to the question"""
from __future__ import print_function, division
import matplotlib.pyplot as plt
import numpy as np
import emcee

def v_func(a, t):
    """With a list or numpy array of parameters and t a list of times, this function 
    returns the function defined in equation 1. """
    #INSERT CODE HERE!
    
def simulate_v(t_min=10, t_max=30, a_in=[0,1,1,1,0], n_points=100, sigma=1):
    """Some default keyword arguments, for a function to generate random data.
    Not all of the default keyword arguments match the assignment question.
    
    Here is some example code. Does it make sense? 
    """
    times = t_min + (t_max - t_min)*np.random.random(n_points)
    v_with_errors = v_func(a_in, times) + sigma*np.random.normal(size=n_points)
    return times, v_with_errors
        
def lnprob_v_func(a_try, times, vs, sigma=1.0):
    #ADD CODE:
    #Compute chi-squared, and return -chi-squared/2

if __name__=="__main__":  
    #Simulate our data. After this point, the times at which the data
    #were taken, and the data values are fixed.
    times, v_with_errors = simulate_v()
    
    #ADD CODE FOR PLOTTING HERE.

    #Define some key variables.
    ndim = 5
    nwalkers = 20
    n_burnin = int(1000)
    n_chain = int(10000)

    #Lets make a starting point as a bunch of random numbers
    a_init = np.random.normal( size=(nwalkers,ndim) )

    #sampler = INSERT_CODE
    pos, lnprob, state = sampler.run_mcmc(a_init, n_burnin)

    #A check plot, to see if burn in is complete.
    #plt.plot(sampler.lnprobability.T)
    #plt.title("Is burn in complete?")

    #Some helper code to help eliminate the worst chains: 
    best_chain = np.argmax(lnprob)
    poor_chains = np.where(lnprob < np.percentile(lnprob, 33))
    for ix in poor_chains:
        pos[ix] = pos[best_chain]
        
    #Some helper code to eliminate awkward negative numbers
    #YOU MIGHT WANT TO TRY SOMETHING HERE, ESPECIALLY FOR QUESTIONS 1-3
    
    #Some helper code to prevent a[4] being many multiples of 2 \pi away from zero.
    #Lets make sure that a[4] is always in the interval [-pi,pi)
    #Uncomment the following line for question (1-3)
    #pos[:,4] += ((pos[:,4]+np.pi) % (2*np.pi))- np.pi
    
    #Reset and go again!
    sampler.reset()
    sampler.run_mcmc(pos, n_chain)
    
    #Use plt.hist or even http://corner.readthedocs.io/en/latest/
    #to examine the outputs:
    #sampler.chain
    #sampler.flatchain (if you're feeling lazy)
    #e.g. 
    #plt.hist(sampler.flatchain[:,0])