function [F] = gen_seasonal_indicator_F(s,omega)
% GEN_SEASONAL_INDICATOR_F - Gen(erate) Seasonal Indicator F matrix
%
% Input: 
%   s     : total number of seasons = seasonal indicators.
%   omega : relaxation coefficient
% 
% Our assumed model is a locally constant LINEAR model with "s" seasonal indicators: 
% 
% z_{n+j} = \beta_0 + \beta_1 j + \sum_{i=1}^{s-1} \delta_i \mathrm{IND}_{ji} + \varepsilon_{n+j} 
% 
% See page 159 of the book by Abraham
%   
% Written by:
% -- 
% John L. Weatherwax                2009-06-01
% 
% email: wax@alum.mit.edu
% 
% Please send comments and especially bug reports to the
% above email address.
% 
%-----

n = s+1; % = 2 + ( s - 1 ); 

F = zeros(n); 

oneOverAlpha = 1/(1-omega); 

F(1,1) = oneOverAlpha;
F(2,2) = omega*( 1 + omega )*oneOverAlpha^3;
F(1,2) = -omega*oneOverAlpha^2;
F(2,1) = F(1,2); 

omegaFrac = 1/(1 - omega^s); 
for j=3:(s+1), 
  F(1,j) = ( omega^(s+2-j) ) * omegaFrac;
  F(j,1) = F(1,j); 
  F(j,j) = F(1,j); 
end

omegaFrac = omegaFrac * omegaFrac;
for j=3:(s+1), 
  F(2,j) = -( omega^(s+2-j) ) * ( (s+2-j) + (j-2)*omega^s ) * omegaFrac;
  F(j,2) = F(2,j); 
end