from math import *

def gamma(y, T, alpha):
  l = 2.65
  M = 10.142
  m = 0.020
  r = 0.0075
  a = 0.024
  b = 0.1015
  y0 = 1.265 
  y = y - y0
  phi = 0.5 * atan(y/l)
  alpha = 2*pi/360.0 * alpha
  return 4.0*(pi**2.0)*phi*(2.0/5.0*r**2 + a**2)*(sqrt(a**2 + b**2 - 2*a*b*cos(alpha))**3)/(T**2 * M * a * b* sin(alpha))

def tors(y, alpha):
  l = 2.65
  M = 10.142
  m = 0.020
  r = 0.0075
  a = 0.024
  b = 0.1015
  gam = 6.67384E-11
  phi = 0.5 * atan(y/l)
  alpha = 2*pi/360.0 * alpha
  R = 0.00001
  L = 0.72
  return 8.0*L*a*b*gam*M*m*sin(alpha)/(pi*R**4*2*phi*sqrt(a**2 + b**2 - 2*a*b*cos(alpha))**3)

def torserror(y, alpha):
  l = 2.65
  M = 10.142
  m = 0.020
  r = 0.0075
  a = 0.024
  b = 0.1015
  gam = 6.67384E-11
  phi = 0.5 * atan(y/l)
  alpha = 2*pi/360.0 * alpha
  R = 0.00001
  L = 0.72
  G = 8.0*L*a*b*gam*M*m*sin(alpha)/(pi*R**4*2*phi*sqrt(a**2 + b**2 - 2*a*b*cos(alpha))**3)
  yerr = 0.0016/sqrt(10)
  alphaerr = pi/90
  Gerr = G*sqrt((yerr*L/((L**2 + y**2)*phi))**2 + (alphaerr*(2*(a**2 + b**2)*cos(alpha) + a*b*(cos(2*alpha)-5))/(2*(a**2 + b**2 -2*a*b*cos(alpha))))**2)
  return Gerr

y = 1.1394
T = 548.0
alpha = -45
gamma1 = gamma(y, T, alpha)
print "Gravitationskonstante: gamma1 = ", gamma1

y = 1.3703
T = 539.9
alpha = +45
gamma2 = gamma(y, T, alpha)
print "Gravitationskonstante: gamma2 = ", gamma2
print "Mittelwert: gamma = ", (gamma1+gamma2)/2.0
print "Literaturwert: gamma = ", 6.67384E-11
print "relative Abweichung: ", 100*(6.67384E-11-(gamma1+gamma2)/2.0)/6.67384E-11, "%"

alpha = 45
y0 = 1.265 
y1 = 1.1394
y2 = 1.3703
y  = 0.5 * (fabs(y1-y0) + fabs(y2-y0))
torsmom = tors(y, alpha)
terror = torserror(y, alpha)
print ""
print "Torsionsmodul des Fadens: G = ", torsmom, " plusminus ", terror

