plot(beam_strengths, trap_strengths, 'bo', label='Trap Force')
plot(line, fit, 'r-', label='Best Fit Trap Force vs Beam Strength')
title('Trap Force vs Beam Strength')
xlabel('Beam Power (mW)')
ylabel('Trap Force (pN)')
legend()
print("Slope is " + str(p[0]))

mu = 8.9e-4 #Pa*S
R = 2.56e-6/2 #meters

T = .7 #seconds
V = 54 #Volts
dist = (20e-6*54)/75
v_s = dist/T
print("v_s is " + str(v_s))
F_d = 6*pi*mu*R*v_s
print("F_d is " + str(F_d))
k1 = F_d

T = .43 #seconds
V = 54 #Volts
dist = (20e-6*54)/75
v_s = dist/T
print("v_s is " + str(v_s))
F_d = 6*pi*mu*R*v_s
print("F_d is " + str(F_d))
k2 = F_d

T = .36 #seconds
V = 54 #Volts
dist = (20e-6*54)/75
v_s = dist/T
print("v_s is " + str(v_s))
F_d = 6*pi*mu*R*v_s
print("F_d is " + str(F_d))
k3 = F_d

T = .23 #seconds
V = 54 #Volts
dist = (20e-6*54)/75
v_s = dist/T
print("v_s is " + str(v_s))
F_d = 6*pi*mu*R*v_s
print("F_d is " + str(F_d))
k4 = F_d

T = .17 #seconds
V = 54 #Volts
dist = (20e-6*54)/75
v_s = dist/T
print("v_s is " + str(v_s))
F_d = 6*pi*mu*R*v_s
print("F_d is " + str(F_d))
k5 = F_d

T = .13 #seconds
V = 54 #Volts
dist = (20e-6*54)/75
v_s = dist/T
print("v_s is " + str(v_s))
F_d = 6*pi*mu*R*v_s
print("F_d is " + str(F_d))
k6 = F_d

beam_strengths = [5.5, 8.5, 11.7, 15.3, 19.0, 23.0] #mW
trap_strengths = array([k1, k2, k3, k4, k5, k6]) #N
trap_strengths = trap_strengths*1e12 #convert to pN

plot(beam_strengths,trap_strengths, 'bo')
p = polyfit(beam_strengths, trap_strengths, 1)
line = linspace(5,30,2)
fit = p[1] + p[0]*line
plot(line, fit, 'r-')

p[1] + p[0]*67.8