# -*- coding: utf-8 -*-
"""
Created on Wed Jan 24 13:59:38 2021
@author: neil
"""

# this program illustrates creating a function: that is a stand-alone piece of code that has defined
# input and output variables, and which can be called from your main program to perform some specific
# function. (all the numpy and scipy methods are just fancy 'functions' like this simple function)
# One important note about 'functions': The names of variables inside of a function are completely separate
# from the variable names in the main program. For example if the main program has a variable 'x' this is 
# unrelated to a variable 'x' in the function. This makes the function stand-alone. (technical term is: different 'name_spaces')
# Read the comments below to see how a function is set up and how it returns values
# the program plots a sine curve in radians, and then calculates the derivative, and plots that. And then plots a cosine.

import numpy as np
import matplotlib.pyplot as plt

# function definitions have to be constructed before they are called, so here is the function 'mydiff'
#'mydiff' is a simple numerical derivative calculator. It expects evenly spaced data, and returns
#  2 arrays, the x positions of the derivatives, and the derivative values
# A function is a piece of indented code starting with the word 'def', followed by a name, followed by a list of arguments.
# Arguments in the function are dummy names that replace the names in the argument list from where the function was called 
#  in the main program.  As with all python indented code, the 'def' line ends with a ':'

def mydiff(xin,yin):                     # create a function called 'mydiff', which is a primative derivative calculator
    """ takes evenly spaced data and returns the derivative on 
    a half step shifted grid """         #by convention we add a comment about usage (python will include this in console error msgs)
    diff = np.zeros(len(yin)-1)          # 'len' is a built-in method which returns the length of a countable object
    delx = xin[1]-xin[0]                 # assumes evenly spaced data, (we should really check this with an 'if' statement)
    for i in np.arange(0,len(yin)-1):    # 'np.arange()' defaults to integer (1) steps
        diff[i] = yin[i+1] - yin[i]
    diff = diff/delx                     # divide all the diffs by delx, this constructs 'dy/dx'
    xout = xin[:-1]+delx*.5              # create a x array that is moved to the middle of each old x interval
    return xout,diff                     # return 2 arrays (actually vectors), the x positions and the slopes or dy/dx at each
	                                     # note the names 'xout,diff' are not know outside the function.  What is returned is 
										 # actually just the arrays.

# ********************************************************************************************************    
# this is the main program, uses improved coding practise, eg a header section (only does basic plotting, using 'plt.plot')
# ******************* header information *****************************************************************
xstart  = 0                             # we need to dfine the starting point of our region of the sine curve to differentiate
xstop   = 7                             # the 'span' of our problem, in radians
numberdx = 50                           # number of points at which to calculate dy/dx

x = np.linspace(xstart,xstop,numberdx)  # linspace is a very useful numpy method to create a list (array) of numbers
w = np.sin(x)                           # built-in trig functions, takes an array and returns an array of sine values
plt.plot(x,w,'r-o')                     # minimal plotting (hard to modify or add stuff, see other examples for better plotting)

x2,d2 = mydiff(x,w)                     # 'call' our function, the 'x,w' are called arguements and are sent to the function.
                                        # Since the function returns 2 arrays, we assign the function results to 2 arrays 'x2,d2'
plt.plot(x2,d2,'k-o')                   # by default plots on the previous plot, would need a new plt.figure() to get a new plot
plt.plot(x,np.cos(x),'b+')              # Note you do not need to create unnecessary variables: 'np.cos(x)' IS an array
plt.show()