% Example 2.11
% modified from solution of Example Problem 2.10, cooling billet.
clear all; format short g
n=17; %Measured temperature in Celsius. Time step is 96 sec.
y=[137.55 129.37 121.96 115.17 108.99 102.91 97.7 92.98 88.69 84.69 ...
    80.91 77.45 73.91 71.05 68.44 65.99 63.71]';
R=y'*y; t=0:1:n-1; t=t'; %Only for compactness is more than one statement on line.
dt = 96;
figure(1); clf;
% plot(t*96, y)
hold on
p = 5;  % from Example 2.4.7
  X=[ones(n,1)];
  for jj=1:p-1
      tj=t.^jj; X=[X tj];
  end
  b=X\y; e=y-X*b; sum=e'*e; DR=R-sum; 
  s2=sum/(n-p); Con=cond(X); [S]=svd(X);
  
dYdt = ( b(2) + 2*b(3)*t + 3*b(4)*( t ).^2 + 4*b(5).*(t).^3)/dt;
Y = X*b;
rho = 8890; % kg/m^3
c = 385;    % J/kg-K
VoverA = 0.00498;  % m
Tinf = 27.5;  % C
h = rho * c *VoverA./(Tinf-Y).*dYdt;
figure(1);
% plot( [ 0 1600], [0 0], 'r-' )
xlabel('time, seconds')
ylabel('Estimated h_i, W/m^2-C')
plot( t*dt, h, 'bo', 'linewidth', 1.0 )
% print the 600 dpi figure to tiff format
fig=figure(1);
print(fig,fullfile(pwd, 'Figure2_5.tif'),'-dtiff','-r600')