#!/usr/bin/env python
# -*- coding: utf-8 -*-

# 08/21/2014

# ***************************************************************************
# *   Copyright (C) 2014, Paul Lutus                                        *
# *                                                                         *
# *   This program is free software; you can redistribute it and/or modify  *
# *   it under the terms of the GNU General Public License as published by  *
# *   the Free Software Foundation; either version 2 of the License, or     *
# *   (at your option) any later version.                                   *
# *                                                                         *
# *   This program is distributed in the hope that it will be useful,       *
# *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
# *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
# *   GNU General Public License for more details.                          *
# *                                                                         *
# *   You should have received a copy of the GNU General Public License     *
# *   along with this program; if not, write to the                         *
# *   Free Software Foundation, Inc.,                                       *
# *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             *
# ***************************************************************************

from math import *

# BEGIN user-defined values

# tank radius
tr = 30

# top of sloped section
st = 30

# top of tank
tt = 70

# slope section integration steps
ssi = 100

# data table step size
tss = 1

# END user-defined values

# classic interpolation function

def ntrp(x,xa,xb,ya,yb):
  return (x-xa) * (yb-ya) / (xb-xa) + ya

# liquid area equation
# from http://arachnoid.com/TankCalc/index.html figure 3

def la(y,r):
  return r*r * acos(-y/r)+y*sqrt(r*r-y*y)

# compute partial volume for a provided y argument

def partial(y,r):
  global st,ssi
  # if y value is above the sloped section,
  # provide 1/2 the volume of the sloped section
  # plus the partial volume of the cylindrical section
  if(y >= st):
    # volume at y in cylindrical section
    # includes 1/2 volume of sloped section
    return pi * r*r * (y-st/2.0)
  # else 0 <= y < st, need to
  # integrate partial volume
  # of sloped section
  else:
    v = 0
    fssi = float(ssi)
    for n in range(ssi+1):
      x = n/fssi
      q = y * x/st
      yy = ntrp(q,0,1,-r,r)
      v += la(yy,r)
    return v * y/(fssi+1)
  
# print table column titles
  
s = '%12s %12s %12s' % ('Sensor','Volume','%')
print(s)
print('-' * len(s))

# old and new table line strings
ols = ""
ls = ""

# table accumulator

def accum(y):
  global ols,ls
  v = partial(y,tr)
  ls = '%12.4f %12.4f %12.4f' % (y,v, 100.0*y/tt)
  # don't print duplicate lines
  if(ols != ls):
    print(ls)
  ols = ls

# generate table

y = 0
while(y <= tt):
  accum(y)
  y += tss

# accumulate a full height reading, added
# to the table only if it was missed
# by the preceding loop
accum(tt)

