% Perform Stolz's method to estimate heat flux
%   -dimensionless data must be supplied
%   -geometry is the insulated finite plate
%   -uses function fdX22B10T0.m
% dt is the dimensionless time step in the data
%  x is a scalar dimensionless measurement location
%  Y is the vector of measurements beginning at time=dt
% function [qhat,Tsur] = Stolz( dt, x, Y )
function [qhat] = Stolz( dt, x, Y )
if( nargin == 0 )
   dt = 1;
   x = 0.5;
   Y = [1:1:10]';
end

N = size(Y,1);   % the total number of observations
if( N == 1 )
   error('Y must be a column vector');
end
time = [dt:dt:N*dt]';
phi = fdX22B10T0( x, time ); % get the sensitivity values
qhat = zeros(N, 1); 
for M = 1 : N
   % compute the influence of all the previous q's
   That = 0;
   for n = 1:M-1
       That = That + (phi(M-n+1)-phi(M-n)) * qhat(n);
   end
   % now use the Stolz formula to compute the next qhat
   qhat(M) = ( Y(M) - That ) / phi(1);
end