% solve whole domain IHCP with no regularization
% assume dimensionless variables
A = 15;
x = [ 1. ];   % measurement locations
t = [ -0.24: 0.06: 1.8 ]';
dt = t(2) - t(1);
rng(129);
noise = randn( size(t) );
sigma_pct = 0.00;
% T = 0.6*T_1pulse(t, 0, 1.2, x )';
% T = round(T,6);      % round to 6 decimal places
T = [ 0.000000 0.000000 0.000000 0.000000 0.000000 0.000007 ...
      0.000374 0.002171 0.006323 0.013381 0.023656 0.037319 ...
      0.054465 0.075145 0.099389 0.127201 0.157880 0.189293 ...
      0.219594 0.247679 0.272930 0.295005 0.313715 0.328954 ...
      0.340665 0.348823 0.353763 0.356545 0.358089 0.358943 ...
      0.359415 0.359677 0.359821 0.359901 0.359945 ]';

sigma = sigma_pct*max(Y);
Y = T + sigma*noise;
% can't have negative times
tzero = t - t(1);
phi = fdX22B10T0( x, tzero, A );
dphi = diff( phi );
n = length( dphi );  % compute also for first data point
X = zeros( n, n );
for ic = 1:n
    X( ic:n, ic ) = dphi(1:n-ic+1);
end
cond(X); 
qhat = (X'*X) \ X'*Y(2:n+1);  % treat first Y(1) as initial condition
figure(99); clf;
plot( t(2:end) - dt/2, qhat, 'ro', t, 0.6*q_1pulse( t, 0, 1.2), 'k--', 'linewidth', 1.0 )
xlabel('time'); ylabel('Heat flux')
legend({'$\hat{\mathbf q}$', 'Exact'}, 'Interpreter', 'latex')
axis( [ -0.5 2.0 -.2 0.8 ] )

%    fig=figure(99);
%    print(fig,fullfile(pwd, 'Figure4_2.tif'),'-dtiff','-r600')

