import numpy as np
import matplotlib.pyplot as pl
import pandas as pd
from scipy.optimize import curve_fit # automated curve fitting
df_alex = pd.read_csv('alex1.csv')
df_alex.head()
| timefloat64 | anglefloat64 | ||
|---|---|---|---|
0 | 0 | 0 | |
1 | 0.001 | 0 | |
2 | 0.002 | 0 | |
3 | 0.003 | 0 | |
4 | 0.003 | 0 | |
| 5 rows × 2 columns | |||
pl.plot(df_alex.time, df_alex.angle, 'b-')
good_time = (df_alex.time > 8.5) & (df_alex.time<16)
time = df_alex[good_time].time
angle = df_alex[good_time].angle
pl.plot(time, angle ,'b.')
pl.grid()
T0 = 2
A0 = 45
def model(t, A, omega, phi, off):
return A*np.sin(omega*t + phi) + off
par, cov = curve_fit(model, time, angle, p0 = (A0, 2*np.pi/T0, 0, 0))
A, omega, phi, off = par
dA, dOm, dPhi, dOff = np.sqrt(np.diag(cov))
print("A = {0:3.4f}+/-{1:3.4f} degrees, omega = {2:3.4f} +/- {3:3.4f} rad/sec".format(A, dA, omega, dOm))
pl.plot(time, angle ,'b.', label='data')
pl.plot(time, model(time, A, omega, phi, off), 'r-', label="fit")
pl.grid()
pl.title("Fit at around 50 degrees")
pl.xlabel('time (sec)')
pl.ylabel('angle (deg)')
pl.legend()
A = 49.6569+/-0.0853 degrees, omega = 3.4957 +/- 0.0008 rad/sec
df_dimi1 = pd.read_csv('dimitris1.csv')
df_dimi1.head()
| timefloat64 | anglefloat64 | ||
|---|---|---|---|
0 | 0 | 0 | |
1 | 0.001 | 0 | |
2 | 0.002 | 0 | |
3 | 0.002 | 0 | |
4 | 0.003 | 0 | |
| 5 rows × 2 columns | |||
pl.plot(df_dimi1.time, df_dimi1.angle, 'b-')
good_time = (df_dimi1.time > 8.5) & (df_dimi1.time<16)
time = df_dimi1[good_time].time
angle = df_dimi1[good_time].angle
pl.plot(time, angle ,'b.')
pl.grid()
T0 = 2
A0 = 45
def model(t, A, omega, phi, off):
return A*np.sin(omega*t + phi) + off
par, cov = curve_fit(model, time, angle, p0 = (A0, 2*np.pi/T0, 0, 0))
A, omega, phi, off = par
dA, dOm, dPhi, dOff = np.sqrt(np.diag(cov))
print("A = {0:3.4f}+/-{1:3.4f} degrees, omega = {2:3.4f} +/- {3:3.4f} rad/sec".format(A, dA, omega, dOm))
pl.plot(time, angle ,'b.', label='data')
pl.plot(time, model(time, A, omega, phi, off), 'r-', label="fit")
pl.grid()
pl.xlabel('time (sec)')
pl.ylabel('angle (deg)')
pl.legend()
A = 23.9288+/-0.0583 degrees, omega = 3.6381 +/- 0.0011 rad/sec
df_dimi2 = pd.read_csv('dimitris2.csv')
df_dimi2.head()
| timefloat64 | anglefloat64 | ||
|---|---|---|---|
0 | 0 | 0 | |
1 | 0.001 | 0 | |
2 | 0.002 | 0 | |
3 | 0.003 | 0 | |
4 | 0.004 | 0 | |
| 5 rows × 2 columns | |||
pl.plot(df_dimi2.time, df_dimi2.angle, 'b-')
good_time = (df_dimi2.time > 8.5) & (df_dimi2.time<16)
time = df_dimi2[good_time].time
angle = df_dimi2[good_time].angle
pl.plot(time, angle ,'b.')
pl.grid()
T0 = 2
A0 = 45
def model(t, A, omega, phi, off):
return A*np.sin(omega*t + phi) + off
par, cov = curve_fit(model, time, angle, p0 = (A0, 2*np.pi/T0, 0, 0))
A, omega, phi, off = par
dA, dOm, dPhi, dOff = np.sqrt(np.diag(cov))
print("A = {0:3.4f}+/-{1:3.4f} degrees, omega = {2:3.4f} +/- {3:3.4f} rad/sec".format(A, dA, omega, dOm))
pl.plot(time, angle ,'b.', label='data')
pl.plot(time, model(time, A, omega, phi, off), 'r-', label="fit")
pl.grid()
pl.title("Fit at around 50 degrees")
pl.xlabel('time (sec)')
pl.ylabel('angle (deg)')
pl.legend()
A = 36.5456+/-1.6959 degrees, omega = 3.4562 +/- 0.0209 rad/sec
df_mitch1 = pd.read_csv('mitchell1.csv')
df_mitch1.head()
| timefloat64 | anglefloat64 | ||
|---|---|---|---|
0 | 0 | 0 | |
1 | 0 | 0 | |
2 | 0.001 | 0 | |
3 | 0.002 | 0 | |
4 | 0.003 | 0 | |
| 5 rows × 2 columns | |||
pl.plot(df_mitch1.time, df_mitch1.angle, 'b-')
good_time = (df_mitch1.time > 8.5) & (df_mitch1.time<16)
time = df_mitch1[good_time].time
angle = df_mitch1[good_time].angle
pl.plot(time, angle ,'b.')
pl.grid()
T0 = 2
A0 = 45
def model(t, A, omega, phi, off):
return A*np.sin(omega*t + phi) + off
par, cov = curve_fit(model, time, angle, p0 = (A0, 2*np.pi/T0, 0, 0))
A, omega, phi, off = par
dA, dOm, dPhi, dOff = np.sqrt(np.diag(cov))
print("A = {0:3.4f}+/-{1:3.4f} degrees, omega = {2:3.4f} +/- {3:3.4f} rad/sec".format(A, dA, omega, dOm))
pl.plot(time, angle ,'b.', label='data')
pl.plot(time, model(time, A, omega, phi, off), 'r-', label="fit")
pl.grid()
pl.title("Fit at around 50 degrees")
pl.xlabel('time (sec)')
pl.ylabel('angle (deg)')
pl.legend()
A = 31.3807+/-0.0897 degrees, omega = 3.6081 +/- 0.0021 rad/sec
df_soren1 = pd.read_csv('soren1.csv')
df_soren1.head()
| timefloat64 | anglefloat64 | ||
|---|---|---|---|
0 | 0 | 0 | |
1 | 0.001 | 0 | |
2 | 0.002 | 0 | |
3 | 0.002 | 0 | |
4 | 0.003 | 0 | |
| 5 rows × 2 columns | |||
pl.plot(df_soren1.time, df_soren1.angle, 'b-')
good_time = (df_soren1.time > 8.5) & (df_soren1.time<16)
time = df_soren1[good_time].time
angle = df_soren1[good_time].angle
pl.plot(time, angle ,'b.')
pl.grid()
T0 = 2
A0 = 45
def model(t, A, omega, phi, off):
return A*np.sin(omega*t + phi) + off
par, cov = curve_fit(model, time, angle, p0 = (A0, 2*np.pi/T0, 0, 0))
A, omega, phi, off = par
dA, dOm, dPhi, dOff = np.sqrt(np.diag(cov))
print("A = {0:3.4f}+/-{1:3.4f} degrees, omega = {2:3.4f} +/- {3:3.4f} rad/sec".format(A, dA, omega, dOm))
pl.plot(time, angle ,'b.', label='data')
pl.plot(time, model(time, A, omega, phi, off), 'r-', label="fit")
pl.grid()
pl.title("Fit at around 50 degrees")
pl.xlabel('time (sec)')
pl.ylabel('angle (deg)')
pl.legend()
A = 9.2986+/-0.0468 degrees, omega = 3.6744 +/- 0.0024 rad/sec
df_soren2 = pd.read_csv('soren2.csv')
df_soren2.head()
| timefloat64 | anglefloat64 | ||
|---|---|---|---|
0 | 0 | 0 | |
1 | 0.001 | 0 | |
2 | 0.002 | 0 | |
3 | 0.003 | 0 | |
4 | 0.004 | 0 | |
| 5 rows × 2 columns | |||
pl.plot(df_soren2.time, df_soren2.angle, 'b-')
bad_time = (df_soren2.time > 8.5) & (df_soren2.time<16)
time = df_soren2[bad_time].time
angle = df_soren2[bad_time].angle
pl.plot(time, angle ,'b.')
pl.grid()
T0 = 2
A0 = 45
def model(t, A, omega, phi, off):
return A*np.sin(omega*t + phi) + off
par, cov = curve_fit(model, time, angle, p0 = (A0, 2*np.pi/T0, 0, 0))
A, omega, phi, off = par
dA, dOm, dPhi, dOff = np.sqrt(np.diag(cov))
print("A = {0:3.4f}+/-{1:3.4f} degrees, omega = {2:3.4f} +/- {3:3.4f} rad/sec".format(A, dA, omega, dOm))
pl.plot(time, angle ,'b.', label='data')
pl.plot(time, model(time, A, omega, phi, off), 'r-', label="fit")
pl.grid()
pl.title("Fit at around 50 degrees")
pl.xlabel('time (sec)')
pl.ylabel('angle (deg)')
pl.legend()
A = 3.3849+/-0.7718 degrees, omega = 3.3386 +/- 0.1057 rad/sec
import os
import csv
fname = 'rawdata.csv'
dataFile = open(os.path.join('rawdata',fname))
reader = csv.reader(dataFile, delimiter=',')
lineCnt=0
t1_list=[]
t2_list=[]
theta1_list=[]
theta2_list=[]
for line in reader:
if lineCnt in [0]:
pass
else:
try:
t1_list.append(float(line[1]))
theta1_list.append(float(line[2])*np.pi/180.0) # convert to radians
t2_list.append(float(line[3]))
theta2_list.append(float(line[4])*np.pi/180.0) # convert to radians
except:
print("Ack!", line)
lineCnt+=1
N=len(t1_list)
t1arr=np.array(t1_list)
t2arr=np.array(t2_list)
th1arr=np.array(theta1_list)+off # theta 1 corrected for offset
th2arr=np.array(theta2_list)+off # theta 2 corrected for offset
ratio = (t2arr-t1arr)/T0 # array of experimental periods devided by zero amplitude period
thAvg = (th1arr+th2arr)/2.0
vsaa = 1.0*np.ones(N)
saa = 1.0+thAvg*thAvg/16
def simpson(f, a, b, N, args=None):
if args is None:
args=[]
x=np.linspace(a,b,2*N) # sample the function
h=(b-a)/(2*N-1)
y=f(*([x] +list(args)))
evens = y[2:-2:2]
odds = y[1:-1:2]
return (y[0] + y[-1] + 2*odds.sum() + 4*evens.sum())*h/3.0
def K_int(phi, m):
return 1.0/np.sqrt(1-m*np.sin(phi)**2)
def K(m, N=1000):
return simpson(K_int, 0, np.pi/2, N=N, args=(m,))
sraList=[]
for theta in thAvg:
sraList.append(2*K(np.sin(theta/2.0)**2)/np.pi)
thDeg=thAvg*180/np.pi
# array of average amplitudes (in radians)
pl.title("Amplitude of a simple pendulum")
pl.xlabel("$\\theta_{\\rm max}$ (degrees)")
pl.ylabel("T/To")
pl.axis([30,70,0.99,1.1])
pl.plot(thDeg, ratio,'b.', label="experiment")
pl.plot(thDeg, saa, 'r-',label="small angle approx.")
pl.plot(thDeg, vsaa, 'g-', label="very small angle approx.")
pl.plot(thDeg, sraList, 'c-', label="exact theory")
pl.legend(loc=2)