from math import *
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
import copy
[docs]def parse_xy_data(filename):
"""Read the given file,
parse the xy-data and
return the x and y coordinates as two numpy-arrays.
Assume the file has xy-data so that each row has the x and y
coordinates of a single data point, separated by a comma and a space.
That is, the file format is
.. code-block:: python
x1, y1
x2, y2
...
.. note ::
This function is incomplete!
Args:
filename (str): name of the file to be read
Returns:
array, array: xdata, ydata
"""
xdata = None
ydata = None
return xdata, ydata
[docs]def linear(x, a, b):
"""Calculates the linear function :math:`y=ax+b` and returns the result.
Args:
x (float): the variable
a (float): the slope
b (float): the constant term
Returns:
float: the result
"""
return a * x + b
[docs]def linear_fit(xdata, ydata, name="linear", printout=True):
"""Fit a linear function :math:`y = ax+b` to the given xy data.
Also return the optimal fit values obtained for slope a and constant b.
If printout is True, the fit values will also be printed on screen.
To identify this information, you may also name the operation.
Args:
xdata (array): x values
ydata (array): y values
name (str): identifier for printing the result
printout (bool): if True, results are printed on screen
Returns:
float, float: slope, constant
"""
params, covariance = curve_fit(linear,
np.array(xdata),
np.array(ydata))
if printout:
print()
print("slope for "+name+" fit: "+str(params[0]))
print("constant for "+name+" fit: "+str(params[1]))
print()
return params[0], params[1]
[docs]def plot_line(xdata, slope, constant):
"""Plot the linear function y = slope * x + constant
from the first x value in xdata to the last.
The plot is created in memory. You still need to use show()
or savefig() to actually see the plot.
Args:
xdata (array): x values as a list or numpy array
slope (float): slope a in a*x+b
constant (float): constant b in a*x+b
"""
# We can plot the linear function as a line between two points.
# Calculate the y-values of the start and end point of this line.
x_start = xdata[0]
x_end = xdata[-1]
y_start = linear(x_start, slope, constant)
y_end = linear(x_end, slope, constant)
plt.plot([x_start,x_end],
[y_start,y_end],"-")
[docs]def plot_xy_data(xdata, ydata, style='o'):
"""Plot the given xy data in the given style.
The plot is created in memory. You still need to use show()
or savefig() to actually see the plot.
Args:
xdata (array): x values as a list or numpy array
ydata (array): y values as a list or numpy array
style (str): matplotlib plotting style string
"""
# plot the numerical data
plt.plot(xdata, ydata, style)
[docs]def save_plot(filename):
"""Save all the plots currently in memory in a file and show them on screen.
Afterwards, plots are cleared from memory.
Args:
filename (str): name of the file to write
"""
plt.savefig(filename)
plt.show()
plt.clf()
[docs]def write_file(text, filename):
"""Write the given text to a new file.
Args:
text (str): text to be written
filename (str): name of the file to write
"""
f = open(filename,'w')
f.write(text)
f.close()
[docs]def write_data(xdata, ydata, filename):
"""Write data in a file.
The format is the same as in :meth:`parse_xy_data`,
.. code-block:: python
x1, y1
x2, y2
...
.. note ::
This function is incomplete!
Args:
xdata (array): x values as a list or numpy array
ydata (array): y values as a list or numpy array
filename (str): name of the file to write
"""
writelines = ""
write_file(writelines, filename)
[docs]def linlin_analysis(x_data, y_data, printout=True):
"""
Performs a linear fit for x versus y.
Finds the values for constants :math:`a` and :math:`b`
that minimize the sum of squares error between the given data and
the function
.. math::
y = a x + b
Args:
x_data (array): x values
y_data (array): y values, must be the same size as x_data
printout (bool): if True, results are printed on screen and a plot is saved
Returns:
float, float: slope a, constant b
"""
#
# Do a y vs. x fit and plot the results.
#
linlin_slope, linlin_constant = linear_fit(x_data, y_data, "y - x", printout)
if printout:
plot_xy_data(x_data, y_data)
plot_line(x_data, linlin_slope, linlin_constant)
save_plot("linlin.png")
return linlin_slope, linlin_constant
[docs]def linlog_analysis(x_data, y_data, printout=True):
"""
Performs a linear fit for x versus ln(y).
Finds the values for constants :math:`a` and :math:`b`
that minimize the sum of squares error between the given data and
the function
.. math::
\\ln(y) = a x + b
.. note ::
This function is incomplete!
Args:
x_data (array): x values
y_data (array): y values, must be the same size as x_data
printout (bool): if True, results are printed on screen and a plot is saved
Returns:
float, float: slope a, constant b
"""
#
# Create new array storing the values of log(y).
#
log_y_data = y_data
#
# Do a log(y) vs. x fit and plot the results
#
linlog_slope = 0
linlog_constant = 0
if printout:
plot_xy_data(x_data, log_y_data)
plot_line(x_data, linlog_slope, linlog_constant)
save_plot("linlog.png")
return linlog_slope, linlog_constant
[docs]def loglog_analysis(x_data, y_data, printout=True):
"""
Performs a linear fit for ln(x) versus ln(y).
Finds the values for constants :math:`a` and :math:`b`
that minimize the sum of squares error between the given data and
the function
.. math::
\\ln(y) = a \\ln(x) + b
.. note ::
This function is incomplete!
Args:
x_data (array): x values
y_data (array): y values, must be the same size as x_data
printout (bool): if True, results are printed on screen and a plot is saved
Returns:
float, float: slope a, constant b
"""
#
# Create new arrays storing the values of
# log(x) and log(y) for the x and y.
#
log_x_data = x_data
log_y_data = y_data
#
# Do a log(y) vs. log(x) fit and plot the results.
#
loglog_slope = 0
loglog_constant = 0
if printout:
plot_xy_data(log_x_data, log_y_data)
plot_line(log_x_data, loglog_slope, loglog_constant)
save_plot("loglog.png")
return loglog_slope, loglog_constant
# the main program
[docs]def main(filename):
"""
Reads xy data from the given file and performs linear and logarithmic analysis.
Args:
filename (str): name of the file to read
"""
# read and parse the data
x_data, y_data = parse_xy_data(filename)
#
# Do a y vs. x fit and plot the results
#
linlin_analysis(x_data, y_data)
#
# Do a log(y) vs. log(x) fit and plot the results.
#
loglog_analysis(x_data, y_data)
#
# Do a log(y) vs. x fit and plot the results
#
linlog_analysis(x_data, y_data)
# execute main() if the program is launched as a script,
# but not if it is imported as a module
if __name__ == "__main__":
main('A1_data.txt')