%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Get IRFS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function s = makeirfsecondorder_stochastic(M_,oo_,burn_in,shock_label,shock_size,nperiods,nreps,std_vec)

order = 2;

shock = zeros(1,M_.exo_nbr);

shock_pos = strmatch(shock_label,M_.exo_names,'exact');
shock(shock_pos)=shock_size;


% first get the stochastic fixed point -- that way we avoid mixing the
% transtition to a stochastic environment.


[~, ss_fixed_pt] = mymkirf2(oo_.dr,oo_.dr.ys,1,zeros(1,M_.exo_nbr),0,order,1,zeros(M_.endo_nbr,1),burn_in);  % returns matrix containing IRFs for all variables
                    
history = zeros(M_.endo_nbr,nperiods);

init = ss_fixed_pt;

for i = 1:nreps
    [history_noshock] = mymkirf2_stoch(oo_.dr,oo_.dr.ys,nperiods+1,std_vec'.*randn(1,M_.exo_nbr),init);
    init = history_noshock(:,1);
    history_shock = mymkirf2_stoch(oo_.dr,oo_.dr.ys,nperiods,shock,init);
                                  
    history = history+history_shock-history_noshock(:,2:end);
end

history = history/nreps;

reordered = M_.endo_names(oo_.dr.order_var,:);          % the vector reordered contains names of variables
                                            % corresponding to rows of
                                            % history

% this loop assigns each row of history to its variable
% it also assigns steady state values
% IRFs are stored under varname_irf  where varname is the variable name
% SS values are stored under varname
% NB: SS values are scalar, IRFs are column vectors.

% history = history(dr_.order_var,:);
history = history';

for indxi = 1:M_.endo_nbr
    eval(['s','.',deblank(M_.endo_names(indxi,:)),'_ss = oo_.dr.ys(indxi);']);
    eval(['s','.',deblank(reordered(indxi,:)),'_irf=history(:,indxi);']);
end
