function cpav = CpAavg(T)
a = 6.8132;
b = 278.6;
c = -156.28;
d = 34.76;
Tref = 25+273.15;
Trefk = Tref/1000;
Tk = T/1000;
CpT = a*Tk + b/2*Tk^2 + c/3*Tk^3 + d/4*Tk^4;
CpTref = a*Trefk + b/2*Trefk^2 + c/3*Trefk^3 + d/4*Trefk^4;
cpav = (CpT-CpTref)/(Tk-Trefk);

