% prediction error filter - laguardia airport
clear all;

% load laguardia airport temperature data
T = load('lgameanT.txt');
N = length(T);
dt = 1.0;
tmin = 0.0;
t = tmin + dt*[1:N]';
tmax = tmin + dt*(N-1);

M=100;
G = zeros(N+1,M);
d = zeros(N+1,1);

% matrix formulation of temperature
for p = [1:M]
    G(p:N,p) = T(1:N-p+1);
    d(p)=0;
end

% single constrait with large epsilon forcing f(1) to be (-1)
G(N+1,1)=1e6;
d(N+1)=-1e6;

% gh
f = inv(G'*G)*G'*d;

y = conv(f,T);

figure(1);
clf;
plot( t(1:M), f );
hold on;
xlabel( 'time, days' );
ylabel( 'filter' );

figure(2);
clf;
Np = N;
plot( t(1:Np), T(1:Np), 'b' );
hold on;
plot( t(1:Np), y(1:Np), 'r' );
axis( [tmin, dt*Np, -50, 150]' );
xlabel( 'time, days' );
ylabel( 'temp, degF' );

E = y'*y;
E