X Y XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX

Data File 1.txt
X Y
0.05 0.956
0.11 1.09
0.15 1.332
0.31 0.717
0.46 0.771
0.52 0.539
0.70 0.378
0.74 0.370
0.82 0.306
0.98 0.242
1.17 0.104
HW7b_OOP.py
import numpy as np
import matplotlib.pyplot as pyplot
from scipy.optimize import curve_fit
from math import *
class MyLS_Class():
def __init__(self, xdata, ydata):
self.x = xdata
self.y = ydata
def RSquared(self, a):
'''
To calculate the R**2 value for a set of x,y data and a LeastSquares fit with polynomial having coefficients a
:param x:
:param y:
:param a:
:return:
'''
AvgY = np.mean(self.y) # calculates the average value of y
SSTot = 0
SSRes = 0
for i in range(len(self.y)):
SSTot += (self.y[i] - AvgY) ** 2
SSRes += (self.y[i] - self.Poly(self.x[i], a)) ** 2
RSq = 1 - SSRes / SSTot
return RSq
def Poly(self, xval, a):
'''
calculates the value for a polynomial given a value for x and the coefficients of the polynomial.
f(x)=y=a[0]+a[1]x+a[2]x**2+a[3]x**3+...
:param a:
:return:
'''
fx = 0
for i in range(len(a)):
fx += a[i] * xval ** i
return fx
def LeastSquares(self, power):
'''
Calculates the coefficients for a polynomial of degree power to best
fit a data set (x,y) using the least squares approach.
:param x: the independent variable of the data set
:param y: the value of the function evaluated at each value of x
:param power: the degree of the polynomial
:return: the array of coefficients (i.e., f(x)=a[0]+a[1]x+a[2]x**2+...)
'''
A = np.zeros((power + 1, power + 1)) # make a square array of zeros
b = np.zeros(power + 1) # make a vector of zeros
for i in range(len(self.x)): # fill out the elements of the A matrix and b vector
for r in range(len(A)):
for c in range(len(A[0])):
A[r][c] += self.x[i] ** (r + c)
b[r] += self.y[i] * self.x[i] ** r
a = np.linalg.solve(A, b) # note: I am using numpy to solve the linear equation [A][x]=[b]
return a
def GetEq(self, a):
strEq = ""
lbl = "const"
if len(a) == 2: lbl = "linear"
if len(a) == 3: lbl = "quadratic"
if len(a) == 4: lbl = "cubic"
# if len(a) == 5: lbl = "exponential"
for i in range(len(a)):
if i == 0:
strEq += "{:0.4f}".format(a[i])
else:
if a[i] >= 0.0: strEq += "+"
if i > 1:
strEq += "{:0.4f}*x^{:0d}".format(a[i], i)
else:
strEq += "{:0.4f}*x".format(a[i])
return strEq, lbl
# Function for generating exponential function
def expofun(self, x, a, b):
return a * np.exp(b * x)
def GetPlotInfo(self, power, npoints=500):
# For Exponential function
if power == 4:
Xmin = min(self.x)
Xmax = max(self.x)
Ymin = min(self.y)
Ymax = max(self.y)
dX = 1.0 * (Xmax - Xmin) / npoints
# Finding coefficients of exponential function using "Scipy.optimize" function
param, param_cov = curve_fit(self.expofun, self.x, self.y)
a = [param[0], param[1]] # Coefficients
xvals = []
yvals = []
for i in range(npoints):
xvals.append(Xmin + i * dX)
yvals.append(a[0] * np.exp(a[1] * xvals[i]))
RSq = self.RSquared(a)
eq = '{:0.4f}*e^({:0.4f}x)'.format(a[0], a[1])
lbl = "exponential"
return xvals, yvals, RSq, eq, lbl
# For other functions
else:
Xmin = min(self.x)
Xmax = max(self.x)
Ymin = min(self.y)
Ymax = max(self.y)
dX = 1.0 * (Xmax - Xmin) / npoints
a = self.LeastSquares(power)
eq, lbl = self.GetEq(a)
...