############################################################################
## CALCULUS - BIOINFORMATICS, 2018-19 T3
## Authors: Toni Guillamon, Alex Haro
## Sources: [G. Ledder, "Mathematics for the Life Sciences", Springer 2013]
############################################################################

#### Python script to import numpy as np
import numpy as np
import matplotlib.pyplot as plt
from pylab import figure, title, show, savefig

#Try to understand the function of each sentence. Pay also attention to optional features, mostly within the plots.
#Ensure that you are working in the folder you want.
#Copy and paste the commands from this file into the Octave commands' window.

################################
### Lab on Approximation - a ###
################################

#We enter the points of the table.
tp=[0.0,0.5,1.0,1.5,2.0,2.5,3.0,3.5,4.0]
yp=[650.0,552.3,471.5,406.7,347.1,297.6,254.9,219.8,184.2]#JIM
##A second dataset.
#yp1=[650.0,559.5,481.5,414.5,356.7,307.0,264.3,227.5,195.8];#JAN

#Plot in the t-y plane:
figure(1, figsize=(6, 6))

#Print on file
fig, ax = plt.subplots( nrows=1, ncols=1 )
ax.plot(tp,yp,"bo")
title('Exponential decay', bbox={'facecolor': '0.8', 'pad': 5})
fig.savefig('Figure1.png')   # save the figure to file
plt.close(fig)   

#Print on screen
#plt.plot(tp,yp,"bo")
#title('Exponential decay', bbox={'facecolor': '0.8', 'pad': 5})
#plt.show()

################################
### Lab on Approximation - b ###
################################

#Logatirthmic transform of y datasets.
zp=[np.log(650/x) for x in yp]
#zp1=[np.log(650/x) for x in yp1]

################################
### Lab on Approximation - c ###
################################

#Plot in the t-y plane:
figure(2, figsize=(6, 6))

#Print on file
fig, ax = plt.subplots( nrows=1, ncols=1 )
ax.plot(tp,zp,"bo")
title('Exponential decay, logarithmic plot', bbox={'facecolor': '0.8', 'pad': 5})
fig.savefig('FigureLog.png')   # save the figure to file
plt.close(fig)   

################################
### Lab on Approximation - e ###
################################

import numpy as np

#We apply the formulas to fit the parameter kappa for both datasets:
kappa=np.dot(zp,tp)/np.dot(tp,tp)
#kappa1=np.dot(zp1,tp)/np.dot(tp,tp)

#A more refined grid of t's to plot apparently continuous functions.
tf = np.linspace(0, 4, num=101, endpoint=True)

#We put the the results into the t-y coordinates:
S=650*np.exp(-kappa*tf);
#S1=650*exp(-kappa1*tf);

#Plot in the t-y plane:
figure(3, figsize=(6, 6))

#Print on file
fig, ax = plt.subplots( nrows=1, ncols=1 )
ax.plot(tp,yp,"bo")
ax.plot(tf,S,"r")
title('Exponential decay', bbox={'facecolor': '0.8', 'pad': 5})
fig.savefig('Figure1fit.png')   # save the figure to file
plt.close(fig)   

#Plot in the t-log y plane:
figure(4, figsize=(6, 6))

#Print on file
fig, ax = plt.subplots( nrows=1, ncols=1 )
ax.plot(tp,zp,"bo")
ax.plot(tf,kappa*tf,"r")
title('Exponential decay, log plot', bbox={'facecolor': '0.8', 'pad': 5})
fig.savefig('FigureLogFit.png')   # save the figure to file
plt.close(fig)