%--------------------------------------------------------------------------
%  sosfinbeg.m
%
%  Test program for the functions in ComItAlg.m, DiscAlg, and SimpRulT.m.
%  These functions call on qzswitch.m (by Simms) and reorder.m (by Klein).
%  The program also uses Var1SimPs to calculate impulse response functions.
%
%  The program uses the model in
%  Goyal, A. (2008), "The Structure of Inflation, Information and Labour Markets: 
%  Implications for monetary policy", WP-2008-010
%  and "The Natural Interest Rate in Emerging Markets" WP-2008-014
%
%  The MatLab code is modified from an original program written by 
%  Paul Soderlind (nepso@hhs.se), 20 June 2000.
%
%  Ashima Goyal (ashima@igidr.ac.in), 13 October 2008.
%----------------------------------------------------------------------------

diary('sosfindemo.out');            %printing: to screen AND the file
%-----------------------------------------------------------------------

                      
al   = 0.3;              
eta  = 0.4;
phi  = 0.25;
ke   = 0.1;
lam  = 0.24;
bet  = 0.99;
rho  = 0.01;

sigr = 1;
%sig = sigr/eta
w = sigr+(1-al)*(sigr-1);
%siget = (sig*eta*(1-al) +al*sigr)/(eta*(1-al)+w*al)
sigd = sigr/(eta*(1-al)+w*al);
%sigal = sigr/((1-al)+w*al)
%mcs = (sig*eta*(1-al)+al*sigr)/sigr
%d = siget/(phi + siget)
%sum = d*((1-eta)*(sig-siget))

vard = 0.19^2;                   %shocks
vars  = 0.2^2;
vare = 0.2^2;

qy    = 0.7;                      %Loss function
qpi   = 2;                      %at=1; cp=0.08; y*=1
qf    = 1;


[A,B,K,Q,U,R] = sosfindemo(al,eta,phi,ke,lam,bet,rho,vard,vars,vare,qy,qpi,qf);

n1 = 6;
n2 = 3;
k  = 1;

Tbig   = 12;                       %periods to simulate in VAR


%Shock0 = [([sqrt(vard) zeros(11,1)]); Ts ];   %Impulse response wrt ep

%Shock0 = [0;sqrt(vard);zeros(4,1)]; %dd shock

 %Shock0 = [0;0;0;sqrt(vare);zeros(2,1)]; %e shock

 Shock0 = [0;0;sqrt(vars);0;zeros(2,1)];   %Impulse response wrt cost shock

 %Shock0 = [1;0;0;0;zeros(2,1)];            %Impulse response wrt rn shock
%-----------------------------------------------------------------------
 
                             
                                  %DISCRETION

disp('Please wait, discretionary case takes some time to solve');
if exist('OCTAVE_VERSION');
  page_output_immediately = 1;   %seems to make Octave work better
end;
[M_Disc,C_Disc,V,F] = DiscAlg( A,B,Q,R,U,bet,n1,n2,eye(n1),zeros(n2,n1),...
                               [1E-1;1E-5],0,1,0,1,0,1e+4 );
F
x1 = Var1SimPs( M_Disc,Shock0,Tbig );  %VAR of x1(t)
x2 = x1*C_Disc';                        %x2(t)
x  = [x1, x2];                          %x(t) = [x1(t),x2(t)]
uu = - x1*F';                            %u(t)

%picypiie_Disc = [ x(:,5),x(:,7), x(:,8), x(:,9), uu ];

picypiie_Disc = [ ((1-al)*(x(:,8))+ al*(x(:,9))-al*(x(:,4))),x(:,7), x(:,8), x(:,9), uu ]; %fl cpi

tid = linspace(1,Tbig,Tbig)';
if exist('OCTAVE_VERSION');
  title('Discretion: impulse response to a one std of price shock');
  plot( tid,picypiie_Disc(:,1),'-1;Consumer Inflation;',...
        tid,picypiie_Disc(:,2),'-2;Output;',...
        tid,picypiie_Disc(:,2),'-3;Domestic Inflation;',...
        tid,picypiie_Disc(:,2),'-4;Exchange rate;',...
        tid,picypiie_Disc(:,3),'-5;Interest rate;');
else;
  figure(3);
  %plot( tid,picypiie_Disc );
  hold on
  subplot(3,1,1);
  %plot( tid,x(:,5),'-*k');
  plot( tid,((1-al)*(x(:,8))+ al*(x(:,9))-al*(x(:,4))),'-*k'); %current cpi
  title('Figure : Benchmark DIT: cost shock'); %CPI negative shock .01 to natural interest rate 
  hold on 
  plot( tid,x(:,7),'-+r');
  hold on 
  plot( tid,x(:,8),'--dm');
  hold on 
  plot( tid,x(:,9),':hb'); hold on 
  plot( tid,uu,'-^c');
  %title('Full "DIT": impulse response to a calibrated cost shock');
  legend('Consumer inflation','Output','Domestic inflation','Exchange rate','Interest rate',-1);
 end;


disp(' ');disp('Benchmark DIT:  (pi-c,y,pi,e,i)'); %- .01Aa"impulse response to an rn shock:
disp( [linspace(1,Tbig,Tbig)', picypiie_Disc] );
Y = std(picypiie_Disc,1)
disp Y

%-----------------------------------------------------------------------

diary off;
