hot_wire.py

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from scipy.fftpack import fft, ifft
from scipy import signal
import math

# GLOBAL CONSTANTS #
g = 9.80011

temp_wire = 533.440029556275
temp_cal = 294.95

max_corrected_voltage = 0.0


def density_water(T):
    a = 0.9998395
    b = 6.7982999E-05
    c = -9.1060255E-06
    d = 1.0052729E-07
    e = -1.1267135E-09
    f = 6.5917956E-12
    return(a+(b*T)+(c*(T**2))+(d*(T**3))+(e*(T**4))+(f*(T**5)))
    c = 0.020582
    return(((a*p)-(h*((b*T)-(c))))/(273.15+T))



rho_water = 1000 * density_water(21.8)



def density_air(T, p, h):
    a = 0.348444
    b = 0.00252
    c = 0.020582
    return(((a*p)-(h*((b*T)-(c))))/(273.15+T))


rho_air = density_air(21.8, 1021.80, 0.36)
print(rho_air)


def u_fluc10(deg10_velocities):
    sds = np.zeros(len(deg10_velocities))
    sums = np.zeros(len(v_10deg_1))
    i = 0
    for v in v_10deg_1:
        sums[i] = (u(float(v)) - deg10_velocities[0])**2
        i += 1
    sds[0] = math.sqrt(np.sum(sums) / (len(v_10deg_1) - 1))
    sums = np.zeros(len(v_10deg_2))
    i = 0
    for v in v_10deg_2:
        sums[i] = (u(float(v)) - deg10_velocities[1])**2
        i += 1
    sds[1] = math.sqrt(np.sum(sums) / (len(v_10deg_1) - 1))
    sums = np.zeros(len(v_10deg_3))
    i = 0
    for v in v_10deg_3:
        sums[i] = (u(float(v)) - deg10_velocities[2])**2
        i += 1
    sds[2] = math.sqrt(np.sum(sums) / (len(v_10deg_2) - 1))
    sums = np.zeros(len(v_10deg_4))
    i = 0
    for v in v_10deg_4:
        sums[i] = (u(float(v)) - deg10_velocities[3])**2
        i += 1
    sds[3] = math.sqrt(np.sum(sums) / (len(v_10deg_3) - 1))
    sums = np.zeros(len(v_10deg_5))
    i = 0
    for v in v_10deg_5:
        sums[i] = (u(float(v)) - deg10_velocities[4])**2
        i += 1
    sds[4] = math.sqrt(np.sum(sums) / (len(v_10deg_4) - 1))
    sums = np.zeros(len(v_10deg_6))
    i = 0
    for v in v_10deg_6:
        sums[i] = (u(float(v)) - deg10_velocities[5])**2
        i += 1
    sds[5] = math.sqrt(np.sum(sums) / (len(v_10deg_5) - 1))
    sums = np.zeros(len(v_10deg_7))
    i = 0
    for v in v_10deg_7:
        sums[i] = (u(float(v)) - deg10_velocities[6])**2
        i += 1
    sds[6] = math.sqrt(np.sum(sums) / (len(v_10deg_6) - 1))
    sums = np.zeros(len(v_10deg_8))
    i = 0
    for v in v_10deg_8:
        sums[i] = (u(float(v)) - deg10_velocities[7])**2
        i += 1
    sds[7] = math.sqrt(np.sum(sums) / (len(v_10deg_7) - 1))
    sums = np.zeros(len(v_10deg_9))
    i = 0
    for v in v_10deg_9:
        sums[i] = (u(float(v)) - deg10_velocities[8])**2
        i += 1
    sds[8] = math.sqrt(np.sum(sums) / (len(v_10deg_8) - 1))
    sums = np.zeros(len(v_10deg_10))
    i = 0
    for v in v_10deg_10:
        sums[i] = (u(float(v)) - deg10_velocities[9])**2
        i += 1
    sds[9] = math.sqrt(np.sum(sums) / (len(v_10deg_9) - 1))
    sums = np.zeros(len(v_10deg_11))
    i = 0
    for v in v_10deg_11:
        sums[i] = (u(float(v)) - deg10_velocities[10])**2
        i += 1
    sds[10] = math.sqrt(np.sum(sums) / (len(v_10deg_10) - 1))
    sums = np.zeros(len(v_10deg_12))
    i = 0
    for v in v_10deg_12:
        sums[i] = (u(float(v)) - deg10_velocities[11])**2
        i += 1
    sds[11] = math.sqrt(np.sum(sums) / (len(v_10deg_11) - 1))
    sums = np.zeros(len(v_10deg_13))
    i = 0
    for v in v_10deg_13:
        sums[i] = (u(float(v)) - deg10_velocities[12])**2
        i += 1
    sds[12] = math.sqrt(np.sum(sums) / (len(v_10deg_13) - 1))
    
    return(sds)    

def u_fluc20(deg20_velocities):
    sds = np.zeros(len(deg20_velocities))
    sums = np.zeros(len(v_20deg_1))
    i = 0
    for v in v_20deg_1:
        sums[i] = (u(float(v)) - deg20_velocities[0])**2
        i += 1
    sds[0] = math.sqrt(np.sum(sums) / (len(v_20deg_1) - 1))
    sums = np.zeros(len(v_20deg_2))
    i = 0
    for v in v_20deg_2:
        sums[i] = (u(float(v)) - deg20_velocities[1])**2
        i += 1
    sds[1] = math.sqrt(np.sum(sums) / (len(v_20deg_2) - 1))
    sums = np.zeros(len(v_20deg_3))
    i = 0
    for v in v_20deg_3:
        sums[i] = (u(float(v)) - deg20_velocities[2])**2
        i += 1
    sds[2] = math.sqrt(np.sum(sums) / (len(v_20deg_3) - 1))
    sums = np.zeros(len(v_20deg_4))
    i = 0
    for v in v_20deg_4:
        sums[i] = (u(float(v)) - deg20_velocities[3])**2
        i += 1
    sds[3] = math.sqrt(np.sum(sums) / (len(v_20deg_4) - 1))
    sums = np.zeros(len(v_20deg_5))
    i = 0
    for v in v_20deg_5:
        sums[i] = (u(float(v)) - deg20_velocities[4])**2
        i += 1
    sds[4] = math.sqrt(np.sum(sums) / (len(v_20deg_5) - 1))
    sums = np.zeros(len(v_20deg_6))
    i = 0
    for v in v_20deg_6:
        sums[i] = (u(float(v)) - deg20_velocities[5])**2
        i += 1
    sds[5] = math.sqrt(np.sum(sums) / (len(v_20deg_6) - 1))
    sums = np.zeros(len(v_20deg_7))
    i = 0
    for v in v_20deg_7:
        sums[i] = (u(float(v)) - deg20_velocities[6])**2
        i += 1
    sds[6] = math.sqrt(np.sum(sums) / (len(v_20deg_7) - 1))
    sums = np.zeros(len(v_20deg_8))
    i = 0
    for v in v_20deg_8:
        sums[i] = (u(float(v)) - deg20_velocities[7])**2
        i += 1
    sds[7] = math.sqrt(np.sum(sums) / (len(v_20deg_8) - 1))
    sums = np.zeros(len(v_20deg_9))
    i = 0
    for v in v_20deg_9:
        sums[i] = (u(float(v)) - deg20_velocities[8])**2
        i += 1
    sds[8] = math.sqrt(np.sum(sums) / (len(v_20deg_9) - 1))
    sums = np.zeros(len(v_20deg_10))
    i = 0
    for v in v_20deg_10:
        sums[i] = (u(float(v)) - deg20_velocities[9])**2
        i += 1
    sds[9] = math.sqrt(np.sum(sums) / (len(v_20deg_10) - 1))
    sums = np.zeros(len(v_20deg_11))
    i = 0
    for v in v_20deg_11:
        sums[i] = (u(float(v)) - deg20_velocities[10])**2
        i += 1
    sds[10] = math.sqrt(np.sum(sums) / (len(v_20deg_11) - 1))
    sums = np.zeros(len(v_20deg_12))
    i = 0
    for v in v_20deg_12:
        sums[i] = (u(float(v)) - deg20_velocities[11])**2
        i += 1
    sds[11] = math.sqrt(np.sum(sums) / (len(v_20deg_12) - 1))
    sums = np.zeros(len(v_20deg_13))
    i = 0
    for v in v_20deg_13:
        sums[i] = (u(float(v)) - deg20_velocities[12])**2
        i += 1
    sds[12] = math.sqrt(np.sum(sums) / (len(v_20deg_13) - 1))
    sums = np.zeros(len(v_20deg_14))
    i = 0
    for v in v_20deg_14:
        sums[i] = (u(float(v)) - deg20_velocities[13])**2
        i += 1
    sds[13] = math.sqrt(np.sum(sums) / (len(v_20deg_14) - 1))
    sums = np.zeros(len(v_20deg_15))
    i = 0
    for v in v_20deg_15:
        sums[i] = (u(float(v)) - deg20_velocities[14])**2
        i += 1
    sds[14] = math.sqrt(np.sum(sums) / (len(v_20deg_15) - 1))
    sums = np.zeros(len(v_20deg_16))
    i = 0
    for v in v_20deg_16:
        sums[i] = (u(float(v)) - deg20_velocities[15])**2
        i += 1
    sds[15] = math.sqrt(np.sum(sums) / (len(v_20deg_16) - 1))
    sums = np.zeros(len(v_20deg_17))
    i = 0
    for v in v_20deg_17:
        sums[i] = (u(float(v)) - deg20_velocities[16])**2
        i += 1
    sds[16] = math.sqrt(np.sum(sums) / (len(v_20deg_17) - 1))
    sums = np.zeros(len(v_20deg_18))
    i = 0
    for v in v_20deg_18:
        sums[i] = (u(float(v)) - deg20_velocities[17])**2
        i += 1
    sds[17] = math.sqrt(np.sum(sums) / (len(v_20deg_18) - 1))
    sums = np.zeros(len(v_20deg_19))
    i = 0
    for v in v_20deg_19:
        sums[i] = (u(float(v)) - deg20_velocities[18])**2
        i += 1
    sds[18] = math.sqrt(np.sum(sums) / (len(v_20deg_19) - 1))
    sums = np.zeros(len(v_20deg_20))
    i = 0
    for v in v_20deg_20:
        sums[i] = (u(float(v)) - deg20_velocities[19])**2
        i += 1
    sds[19] = math.sqrt(np.sum(sums) / (len(v_20deg_20) - 1))
    sums = np.zeros(len(v_20deg_21))
    i = 0
    for v in v_20deg_21:
        sums[i] = (u(float(v)) - deg20_velocities[20])**2
        i += 1
    sds[20] = math.sqrt(np.sum(sums) / (len(v_20deg_21) - 1))
    
    return(sds)    


def av_voltages(voltages):
    global max_corrected_voltage
    total = 0.0
    for v in voltages:
        if v > max_corrected_voltage:
            max_corrected_voltage = v
        total += float(v)
    return(total/len(voltages))

def read_data(name, cali):
    initial_done = False
    fob = open(name, "r")
    a = fob.readlines()
    fob.close()
    b = []
    p_diff = 0
    for line in a:
        b.append(line.strip())
    t = np.zeros(len(b[25:]))
    v = np.zeros(len(b[25:]))
    i = 0
    for item in b[25:]:
        if not initial_done:
            p_diff = (rho_water * g * 0.5 * 1E-03) * float(item.split()[2])
            temp = float(item.split()[3])
            pressure = float(item.split()[4])
            initial_done = True
        t[i] = (item.split()[0])
        v[i] = (item.split()[1])
        i += 1
    if cali:
        return(t, v, p_diff, temp, pressure)
    else:
        return(t, v)


################ IMPORT DATA #####################
t_Cal_1, v_Cal_1, p_diff_1, temp_Cal_1, pres_Cal_1 = read_data("Calibration_1.lvm", True)
t_Cal_2, v_Cal_2, p_diff_2, temp_Cal_2, pres_Cal_2 = read_data("Calibration_2.lvm", True)
t_Cal_3, v_Cal_3, p_diff_3, temp_Cal_3, pres_Cal_3 = read_data("Calibration_3.lvm", True)
t_Cal_4, v_Cal_4, p_diff_4, temp_Cal_4, pres_Cal_4 = read_data("Calibration_4.lvm", True)
t_Cal_5, v_Cal_5, p_diff_5, temp_Cal_5, pres_Cal_5 = read_data("Calibration_5.lvm", True)
t_Cal_6, v_Cal_6, p_diff_6, temp_Cal_6, pres_Cal_6 = read_data("Calibration_6.lvm", True)
t_Cal_7, v_Cal_7, p_diff_7, temp_Cal_7, pres_Cal_7 = read_data("Calibration_7.lvm", True)
t_Cal_8, v_Cal_8, p_diff_8, temp_Cal_8, pres_Cal_8 = read_data("Calibration_8.lvm", True)
t_Cal_9, v_Cal_9, p_diff_9, temp_Cal_9, pres_Cal_9 = read_data("Calibration_9.lvm", True)
t_Cal_10, v_Cal_10, p_diff_10, temp_Cal_10, pres_Cal_10 = read_data("Calibration_10.lvm", True)
t_Cal_11, v_Cal_11, p_diff_11, temp_Cal_11, pres_Cal_11 = read_data("Calibration_11.lvm", True)
t_Cal_12, v_Cal_12, p_diff_12, temp_Cal_12, pres_Cal_12 = read_data("Calibration_12.lvm", True)
t_Cal_13, v_Cal_13, p_diff_13, temp_Cal_13, pres_Cal_13 = read_data("Calibration_13.lvm", True)
t_Cal_14, v_Cal_14, p_diff_14, temp_Cal_14, pres_Cal_14 = read_data("Calibration_14.lvm", True)
t_Cal_15, v_Cal_15, p_diff_15, temp_Cal_15, pres_Cal_15 = read_data("Calibration_15.lvm", True)

# x = (y - b) / A
v_Cal_1_corr = v_Cal_1 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_1+273.15)))
v_Cal_2_corr = v_Cal_2 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_2+273.15)))
v_Cal_3_corr = v_Cal_3 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_3+273.15)))
v_Cal_4_corr = v_Cal_4 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_4+273.15)))
v_Cal_5_corr = v_Cal_5 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_5+273.15)))
v_Cal_6_corr = v_Cal_6 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_6+273.15)))
v_Cal_7_corr = v_Cal_7 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_7+273.15)))
v_Cal_8_corr = v_Cal_8 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_8+273.15)))
v_Cal_9_corr = v_Cal_9 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_9+273.15)))
v_Cal_10_corr = v_Cal_10 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_10+273.15)))
v_Cal_11_corr = v_Cal_11 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_11+273.15)))
v_Cal_12_corr = v_Cal_12 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_12+273.15)))
v_Cal_13_corr = v_Cal_13 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_13+273.15)))
v_Cal_14_corr = v_Cal_14 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_14+273.15)))
v_Cal_15_corr = v_Cal_15 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_Cal_15+273.15)))


av_Cal_1 = av_voltages(v_Cal_1_corr)
av_Cal_2 = av_voltages(v_Cal_2_corr)
av_Cal_3 = av_voltages(v_Cal_3_corr)
av_Cal_4 = av_voltages(v_Cal_4_corr)
av_Cal_5 = av_voltages(v_Cal_5_corr)
av_Cal_6 = av_voltages(v_Cal_6_corr)
av_Cal_7 = av_voltages(v_Cal_7_corr)
av_Cal_8 = av_voltages(v_Cal_8_corr)
av_Cal_9 = av_voltages(v_Cal_9_corr)
av_Cal_10 = av_voltages(v_Cal_10_corr)
av_Cal_11 = av_voltages(v_Cal_11_corr)
av_Cal_12 = av_voltages(v_Cal_12_corr)
av_Cal_13 = av_voltages(v_Cal_13_corr)
av_Cal_14 = av_voltages(v_Cal_14_corr)
av_Cal_15 = av_voltages(v_Cal_15_corr)


#av_Cal_1 = av_voltages(v_Cal_1)
#av_Cal_2 = av_voltages(v_Cal_2)
#av_Cal_3 = av_voltages(v_Cal_3)
#av_Cal_4 = av_voltages(v_Cal_4)
#av_Cal_5 = av_voltages(v_Cal_5)
#av_Cal_6 = av_voltages(v_Cal_6)
#av_Cal_7 = av_voltages(v_Cal_7)
#av_Cal_8 = av_voltages(v_Cal_8)
#av_Cal_9 = av_voltages(v_Cal_9)
#av_Cal_10 = av_voltages(v_Cal_10)
#av_Cal_11 = av_voltages(v_Cal_11)
#av_Cal_12 = av_voltages(v_Cal_12)
#av_Cal_13 = av_voltages(v_Cal_13)
#av_Cal_14 = av_voltages(v_Cal_14)
#av_Cal_15 = av_voltages(v_Cal_15)




t_10deg_1, v_10deg_1, p_diff_10deg_1, temp_10deg_1, pres_10deg_1 = read_data("10deg_1.lvm", True)
t_10deg_2, v_10deg_2, p_diff_10deg_2, temp_10deg_2, pres_10deg_2 = read_data("10deg_2.lvm", True)
t_10deg_3, v_10deg_3, p_diff_10deg_3, temp_10deg_3, pres_10deg_3 = read_data("10deg_3.lvm", True)
t_10deg_4, v_10deg_4, p_diff_10deg_4, temp_10deg_4, pres_10deg_4 = read_data("10deg_4.lvm", True)
t_10deg_5, v_10deg_5, p_diff_10deg_5, temp_10deg_5, pres_10deg_5 = read_data("10deg_5.lvm", True)
t_10deg_6, v_10deg_6, p_diff_10deg_6, temp_10deg_6, pres_10deg_6 = read_data("10deg_6.lvm", True)
t_10deg_7, v_10deg_7, p_diff_10deg_7, temp_10deg_7, pres_10deg_7 = read_data("10deg_7.lvm", True)
t_10deg_8, v_10deg_8, p_diff_10deg_8, temp_10deg_8, pres_10deg_8 = read_data("10deg_8.lvm", True)
t_10deg_9, v_10deg_9, p_diff_10deg_9, temp_10deg_9, pres_10deg_9 = read_data("10deg_9.lvm", True)
t_10deg_10, v_10deg_10, p_diff_10deg_10, temp_10deg_10, pres_10deg_10 = read_data("10deg_10.lvm", True)
t_10deg_11, v_10deg_11, p_diff_10deg_11, temp_10deg_11, pres_10deg_11 = read_data("10deg_11.lvm", True)
t_10deg_12, v_10deg_12, p_diff_10deg_12, temp_10deg_12, pres_10deg_12 = read_data("10deg_12.lvm", True)
t_10deg_13, v_10deg_13, p_diff_10deg_13, temp_10deg_13, pres_10deg_13 = read_data("10deg_13.lvm", True)

v_10deg_1_corr = v_10deg_1 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_1+273.15)))
v_10deg_2_corr = v_10deg_2 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_2+273.15)))
v_10deg_3_corr = v_10deg_3 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_3+273.15)))
v_10deg_4_corr = v_10deg_4 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_4+273.15)))
v_10deg_5_corr = v_10deg_5 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_5+273.15)))
v_10deg_6_corr = v_10deg_6 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_6+273.15)))
v_10deg_7_corr = v_10deg_7 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_7+273.15)))
v_10deg_8_corr = v_10deg_8 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_8+273.15)))
v_10deg_9_corr = v_10deg_9 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_9+273.15)))
v_10deg_10_corr = v_10deg_10 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_10+273.15)))
v_10deg_11_corr = v_10deg_11 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_11+273.15)))
v_10deg_12_corr = v_10deg_12 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_12+273.15)))
v_10deg_13_corr = v_10deg_13 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_10deg_13+273.15)))

av_10deg_1 = av_voltages(v_10deg_1_corr)
av_10deg_2 = av_voltages(v_10deg_2_corr)
av_10deg_3 = av_voltages(v_10deg_3_corr)
av_10deg_4 = av_voltages(v_10deg_4_corr)
av_10deg_5 = av_voltages(v_10deg_5_corr)
av_10deg_6 = av_voltages(v_10deg_6_corr)
av_10deg_7 = av_voltages(v_10deg_7_corr)
av_10deg_8 = av_voltages(v_10deg_8_corr)
av_10deg_9 = av_voltages(v_10deg_9_corr)
av_10deg_10 = av_voltages(v_10deg_10_corr)
av_10deg_11 = av_voltages(v_10deg_11_corr)
av_10deg_12 = av_voltages(v_10deg_12_corr)
av_10deg_13 = av_voltages(v_10deg_13_corr)

#av_10deg_1 = av_voltages(v_10deg_1)
#av_10deg_2 = av_voltages(v_10deg_2)
#av_10deg_3 = av_voltages(v_10deg_3)
#av_10deg_4 = av_voltages(v_10deg_4)
#av_10deg_5 = av_voltages(v_10deg_5)
#av_10deg_6 = av_voltages(v_10deg_6)
#av_10deg_7 = av_voltages(v_10deg_7)
#av_10deg_8 = av_voltages(v_10deg_8)
#av_10deg_9 = av_voltages(v_10deg_9)
#av_10deg_10 = av_voltages(v_10deg_10)
#av_10deg_11 = av_voltages(v_10deg_11)
#av_10deg_12 = av_voltages(v_10deg_12)
#av_10deg_13 = av_voltages(v_10deg_13)


voltages_av_10deg = np.array([av_10deg_1,av_10deg_2,av_10deg_3,av_10deg_4,av_10deg_5,av_10deg_6,av_10deg_7,av_10deg_8,av_10deg_9,av_10deg_10,av_10deg_11,av_10deg_12,av_10deg_13])

t_20deg_1, v_20deg_1, p_diff_20deg_1, temp_20deg_1, pres_20deg_1 = read_data("20deg_1.lvm", True)
t_20deg_2, v_20deg_2, p_diff_20deg_2, temp_20deg_2, pres_20deg_2 = read_data("20deg_2.lvm", True)
t_20deg_3, v_20deg_3, p_diff_20deg_3, temp_20deg_3, pres_20deg_3 = read_data("20deg_3.lvm", True)
t_20deg_4, v_20deg_4, p_diff_20deg_4, temp_20deg_4, pres_20deg_4 = read_data("20deg_4.lvm", True)
t_20deg_5, v_20deg_5, p_diff_20deg_5, temp_20deg_5, pres_20deg_5 = read_data("20deg_5.lvm", True)
t_20deg_6, v_20deg_6, p_diff_20deg_6, temp_20deg_6, pres_20deg_6 = read_data("20deg_6.lvm", True)
t_20deg_7, v_20deg_7, p_diff_20deg_7, temp_20deg_7, pres_20deg_7 = read_data("20deg_7.lvm", True)
t_20deg_8, v_20deg_8, p_diff_20deg_8, temp_20deg_8, pres_20deg_8 = read_data("20deg_8.lvm", True)
t_20deg_9, v_20deg_9, p_diff_20deg_9, temp_20deg_9, pres_20deg_9 = read_data("20deg_9.lvm", True)
t_20deg_10, v_20deg_10, p_diff_20deg_10, temp_20deg_10, pres_20deg_10 = read_data("20deg_10.lvm", True)
t_20deg_11, v_20deg_11, p_diff_20deg_11, temp_20deg_11, pres_20deg_11 = read_data("20deg_11.lvm", True)
t_20deg_12, v_20deg_12, p_diff_20deg_12, temp_20deg_12, pres_20deg_12 = read_data("20deg_12.lvm", True)
t_20deg_13, v_20deg_13, p_diff_20deg_13, temp_20deg_13, pres_20deg_13 = read_data("20deg_13.lvm", True)
t_20deg_14, v_20deg_14, p_diff_20deg_14, temp_20deg_14, pres_20deg_14 = read_data("20deg_14.lvm", True)
t_20deg_15, v_20deg_15, p_diff_20deg_15, temp_20deg_15, pres_20deg_15 = read_data("20deg_15.lvm", True)
t_20deg_16, v_20deg_16, p_diff_20deg_16, temp_20deg_16, pres_20deg_16 = read_data("20deg_16.lvm", True)
t_20deg_17, v_20deg_17, p_diff_20deg_17, temp_20deg_17, pres_20deg_17 = read_data("20deg_17.lvm", True)
t_20deg_18, v_20deg_18, p_diff_20deg_18, temp_20deg_18, pres_20deg_18 = read_data("20deg_18.lvm", True)
t_20deg_19, v_20deg_19, p_diff_20deg_19, temp_20deg_19, pres_20deg_19 = read_data("20deg_19.lvm", True)
t_20deg_20, v_20deg_20, p_diff_20deg_20, temp_20deg_20, pres_20deg_20 = read_data("20deg_20.lvm", True)
t_20deg_21, v_20deg_21, p_diff_20deg_21, temp_20deg_21, pres_20deg_21 = read_data("20deg_21.lvm", True)

v_20deg_1_corr = v_20deg_1 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_1+273.15)))
v_20deg_2_corr = v_20deg_2 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_2+273.15)))
v_20deg_3_corr = v_20deg_3 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_3+273.15)))
v_20deg_4_corr = v_20deg_4 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_4+273.15)))
v_20deg_5_corr = v_20deg_5 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_5+273.15)))
v_20deg_6_corr = v_20deg_6 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_6+273.15)))
v_20deg_7_corr = v_20deg_7 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_7+273.15)))
v_20deg_8_corr = v_20deg_8 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_8+273.15)))
v_20deg_9_corr = v_20deg_9 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_9+273.15)))
v_20deg_10_corr = v_20deg_10 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_10+273.15)))
v_20deg_11_corr = v_20deg_11 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_11+273.15)))
v_20deg_12_corr = v_20deg_12 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_12+273.15)))
v_20deg_13_corr = v_20deg_13 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_13+273.15)))
v_20deg_14_corr = v_20deg_14 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_14+273.15)))
v_20deg_15_corr = v_20deg_15 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_15+273.15)))
v_20deg_16_corr = v_20deg_16 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_16+273.15)))
v_20deg_17_corr = v_20deg_17 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_17+273.15)))
v_20deg_18_corr = v_20deg_18 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_18+273.15)))
v_20deg_19_corr = v_20deg_19 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_19+273.15)))
v_20deg_20_corr = v_20deg_20 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_20+273.15)))
v_20deg_21_corr = v_20deg_21 * math.sqrt((temp_wire - temp_cal)/(temp_wire - (temp_20deg_21+273.15)))


#av_20deg_1 = av_voltages(v_20deg_1)
#av_20deg_2 = av_voltages(v_20deg_2)
#av_20deg_3 = av_voltages(v_20deg_3)
#av_20deg_4 = av_voltages(v_20deg_4)
#av_20deg_5 = av_voltages(v_20deg_5)
#av_20deg_6 = av_voltages(v_20deg_6)
#av_20deg_7 = av_voltages(v_20deg_7)
#av_20deg_8 = av_voltages(v_20deg_8)
#av_20deg_9 = av_voltages(v_20deg_9)
#av_20deg_10 = av_voltages(v_20deg_10)
#av_20deg_11 = av_voltages(v_20deg_11)
#av_20deg_12 = av_voltages(v_20deg_12)
#av_20deg_13 = av_voltages(v_20deg_13)
#av_20deg_14 = av_voltages(v_20deg_14)
#av_20deg_15 = av_voltages(v_20deg_15)
#av_20deg_16 = av_voltages(v_20deg_16)
#av_20deg_17 = av_voltages(v_20deg_17)
#av_20deg_18 = av_voltages(v_20deg_18)
#av_20deg_19 = av_voltages(v_20deg_19)
#av_20deg_20 = av_voltages(v_20deg_20)
#av_20deg_21 = av_voltages(v_20deg_21)
av_20deg_1 = av_voltages(v_20deg_1_corr)
av_20deg_2 = av_voltages(v_20deg_2_corr)
av_20deg_3 = av_voltages(v_20deg_3_corr)
av_20deg_4 = av_voltages(v_20deg_4_corr)
av_20deg_5 = av_voltages(v_20deg_5_corr)
av_20deg_6 = av_voltages(v_20deg_6_corr)
av_20deg_7 = av_voltages(v_20deg_7_corr)
av_20deg_8 = av_voltages(v_20deg_8_corr)
av_20deg_9 = av_voltages(v_20deg_9_corr)
av_20deg_10 = av_voltages(v_20deg_10_corr)
av_20deg_11 = av_voltages(v_20deg_11_corr)
av_20deg_12 = av_voltages(v_20deg_12_corr)
av_20deg_13 = av_voltages(v_20deg_13_corr)
av_20deg_14 = av_voltages(v_20deg_14_corr)
av_20deg_15 = av_voltages(v_20deg_15_corr)
av_20deg_16 = av_voltages(v_20deg_16_corr)
av_20deg_17 = av_voltages(v_20deg_17_corr)
av_20deg_18 = av_voltages(v_20deg_18_corr)
av_20deg_19 = av_voltages(v_20deg_19_corr)
av_20deg_20 = av_voltages(v_20deg_20_corr)
av_20deg_21 = av_voltages(v_20deg_21_corr)


##################################################

# y=0.8889x+2.8162
# 559.12

wake_velocities = np.zeros(15)


def f(x, A, B):
    return(A*x + B)

def fourth_order(x, A, B, C, D, E):
    return(A*x**4 + B*x**3 + C*x**2 + D*x + E)

def pressure_callibration():
    p_wake_mm = np.array([238, 234, 221, 200, 164, 239, 220, 179, 139, 116, 238, 231, 215, 180, 150])
    p_top_mm = np.array([240, 234, 221, 201, 164, 240, 220, 180, 140, 115, 239, 233, 215, 181, 151])
    p_test_mm = np.array([125, 125, 124, 122, 118, 125, 123, 119, 115, 112, 125, 125, 123, 119, 116])
    p_amb_mm = np.array([110, 111, 112, 112, 112, 110, 111, 112, 112, 112, 110, 111, 112, 112 ,112])
    p_wake = (rho_water * g * 0.5 * (p_wake_mm - p_amb_mm) * 1E-03) + 102180
    p_top = (rho_water * g * 0.5 * (p_top_mm - p_amb_mm) * 1E-03) + 102180
    p_test = (rho_water * g * 0.5 * (p_test_mm - p_amb_mm) * 1E-03) + 102180
    y_mm = p_top_mm - p_test_mm
    x_mm = p_wake_mm - p_amb_mm
    y = p_top - p_test
    x = p_wake - 102180
    return(x, y)

x, y = pressure_callibration()


A, B = curve_fit(f, x, y)[0]

print(type(rho_air))

v1 = math.sqrt(((2*((p_diff_1 - B) / A)))/rho_air)
v2 = math.sqrt(((2*((p_diff_2 - B) / A)))/rho_air)
v3 = math.sqrt(((2*((p_diff_3 - B) / A)))/rho_air)
v4 = math.sqrt(((2*((p_diff_4 - B) / A)))/rho_air)
v5 = math.sqrt(((2*((p_diff_5 - B) / A)))/rho_air)
v6 = math.sqrt(((2*((p_diff_6 - B) / A)))/rho_air)
v7 = math.sqrt(((2*((p_diff_7 - B) / A)))/rho_air)
v8 = math.sqrt(((2*((p_diff_8 - B) / A)))/rho_air)
v9 = math.sqrt(((2*((p_diff_9 - B) / A)))/rho_air)
v10 = math.sqrt(((2*((p_diff_10 - B) / A)))/rho_air)
v11 = math.sqrt(((2*((p_diff_11 - B) / A)))/rho_air)
v12 = math.sqrt(((2*((p_diff_12 - B) / A)))/rho_air)
v13 = math.sqrt(((2*((p_diff_13 - B) / A)))/rho_air)
v14 = math.sqrt(((2*((p_diff_14 - B) / A)))/rho_air)
v15 = math.sqrt(((2*((p_diff_15 - B) / A)))/rho_air)
print(" ")

wake_vels = np.array([v1,v2,v3,v4,v5,v6,v7,v8,v9,v10,v11,v12,v13,v14,v15])
wake_avg_volts = np.array([av_Cal_1,av_Cal_2,av_Cal_3,av_Cal_4,av_Cal_5,av_Cal_6,av_Cal_7,av_Cal_8,av_Cal_9,av_Cal_10,av_Cal_11,av_Cal_12,av_Cal_13,av_Cal_14,av_Cal_15 ])

###############
fs = []
for point in x:
    fs.append(f(point, A, B))
###############

iss = np.linspace(1.8, 2.15, 50)
wsss = []
wsss_uncorrected = []
wA, wB, wC, wD, wE = curve_fit(fourth_order, wake_avg_volts, wake_vels)[0]
print(wA,wB,wC,wD,wE)
coeffs_original = np.array([-481.256547108,3796.04214114,-11116.3797762,14386.1659252,-6954.53632154])

#wA = coeffs_original[0]
#wB = coeffs_original[1]
#wC = coeffs_original[2]
#wD = coeffs_original[3]
#wE = coeffs_original[4]

def u(vol):
    return (wA*vol**4 + wB*vol**3 + wC*vol**2 + wD*vol + wE)


wake_pos_10deg = np.array([12, 10, 8, 6, 4, 2, 0, -2, -4, -6, -8, -10, -12])
wake_vel_10deg = np.zeros(len(wake_pos_10deg))

wake_pos_20deg = np.linspace(20, -20, 21)
wake_vel_20deg = np.zeros(len(wake_pos_20deg))

wake_vel_10deg[0] = u(av_10deg_1)
wake_vel_10deg[1] = u(av_10deg_2)
wake_vel_10deg[2] = u(av_10deg_3)
wake_vel_10deg[3] = u(av_10deg_4)
wake_vel_10deg[4] = u(av_10deg_5)
wake_vel_10deg[5] = u(av_10deg_6)
wake_vel_10deg[6] = u(av_10deg_7)
wake_vel_10deg[7] = u(av_10deg_8)
wake_vel_10deg[8] = u(av_10deg_9)
wake_vel_10deg[9] = u(av_10deg_10)
wake_vel_10deg[10] = u(av_10deg_11)
wake_vel_10deg[11] = u(av_10deg_12)
wake_vel_10deg[12] = u(av_10deg_13)

wake_vel_20deg[0] = u(av_20deg_1)
wake_vel_20deg[1] = u(av_20deg_2)
wake_vel_20deg[2] = u(av_20deg_3)
wake_vel_20deg[3] = u(av_20deg_4)
wake_vel_20deg[4] = u(av_20deg_5)
wake_vel_20deg[5] = u(av_20deg_6)
wake_vel_20deg[6] = u(av_20deg_7)
wake_vel_20deg[7] = u(av_20deg_8)
wake_vel_20deg[8] = u(av_20deg_9)
wake_vel_20deg[9] = u(av_20deg_10)
wake_vel_20deg[10] = u(av_20deg_11)
wake_vel_20deg[11] = u(av_20deg_12)
wake_vel_20deg[12] = u(av_20deg_13)
wake_vel_20deg[13] = u(av_20deg_14)
wake_vel_20deg[14] = u(av_20deg_15)
wake_vel_20deg[15] = u(av_20deg_16)
wake_vel_20deg[16] = u(av_20deg_17)
wake_vel_20deg[17] = u(av_20deg_18)
wake_vel_20deg[18] = u(av_20deg_19)
wake_vel_20deg[19] = u(av_20deg_20)
wake_vel_20deg[20] = u(av_20deg_21)

#for i in iss:
#    wsss.append(fourth_order(i, wA, wB, wC, wD, wE))
#    wsss_uncorrected.append(fourth_order(i, coeffs_original[0],coeffs_original[1],coeffs_original[2],coeffs_original[3],coeffs_original[4]))
#
#plt.plot(iss, wsss)
#plt.plot(iss, wsss_uncorrected)
#plt.show()

uinf10deg = np.zeros(len(wake_vel_10deg))
uinf20deg = np.zeros(len(wake_vel_20deg))

p_diff_10degs =  np.array([p_diff_10deg_1,p_diff_10deg_2,p_diff_10deg_3,p_diff_10deg_4,p_diff_10deg_5,p_diff_10deg_6,p_diff_10deg_7,p_diff_10deg_8,p_diff_10deg_9,p_diff_10deg_10,p_diff_10deg_11,p_diff_10deg_12,p_diff_10deg_13])
p_diff_20degs =  np.array([p_diff_20deg_1,p_diff_20deg_2,p_diff_20deg_3,p_diff_20deg_4,p_diff_20deg_5,p_diff_20deg_6,p_diff_20deg_7,p_diff_20deg_8,p_diff_20deg_9,p_diff_20deg_10,p_diff_20deg_11,p_diff_20deg_12,p_diff_20deg_13,p_diff_20deg_14,p_diff_20deg_15,p_diff_20deg_16,p_diff_20deg_17,p_diff_20deg_18,p_diff_20deg_19,p_diff_20deg_20,p_diff_20deg_21])


AR = (0.005133308 / 0.05)

for i in range(len(wake_vel_10deg)):
    uinf10deg[i] = math.sqrt(  p_diff_10degs[i]   /   (  0.5 * rho_air * (1-(AR**2)) )  )
for i in range(len(wake_vel_20deg)):
    uinf20deg[i] = math.sqrt(  p_diff_20degs[i]   /   (  0.5 * rho_air * (1-(AR**2)) )  )

drag_integrand_10deg = wake_vel_10deg * (uinf10deg - wake_vel_10deg)
drag_integrand_20deg = wake_vel_20deg * (uinf20deg - wake_vel_20deg)
D10 = (rho_air*0.04875*np.trapz(drag_integrand_10deg,wake_pos_10deg/1000))
D20 = (rho_air*0.04875*np.trapz(drag_integrand_20deg,wake_pos_20deg/1000))
print(D10)
print(D20)
print(max_corrected_voltage)

N_fft = len(v_10deg_4)
Fs = 10000
freqsig = abs(np.fft.fft(v_10deg_4, n=N_fft))
freqsig[0:300] = 0.0
tot = 0.0
for f in freqsig:
    tot += f
av = tot / len(freqsig)

freqsig -= av

#freqaxis = np.arange(0, Fs, Fs/N_fft)
#axes = plt.gca()
#axes.set_ylim([0, 100])
#axes.set_xlim([0,Fs/2])
#plt.plot(freqaxis, freqsig)
#plt.show()

#u_flucs_10deg = np.zeros(len(wake_pos_20deg))
#sds = u_fluc(wake_vel_20deg)
#
#plt.plot(wake_pos_20deg, sds)

u_flucs_10deg = u_fluc10(wake_vel_10deg)
u_flucs_20deg = u_fluc20(wake_vel_20deg)