function [] = add_shadow_economy(orig_modfile, ramsey_type,instrument1,instrument2,discount_factor_name,extra_lag_flag)
% extra_lag_flag defaults to 0 if not provided.

global M_ oo_


if strncmp(orig_modfile(length(orig_modfile)-3:length(orig_modfile)),'.mod',4)
    orig_modfile = orig_modfile(1:length(orig_modfile)-4);
end

ramsey_modfile = [orig_modfile,'_',ramsey_type,'_',instrument1,'_',instrument2];
aux_modfile = [orig_modfile,'_aux'];

model_block = 0;
delimiters = char(' ','(',')',',','*','/','+','-','=','^',';');
policy_vars = {'Util';'Welf';'Util1';'Util2';'Welf1';'Welf2';'instr1';'instr2'};

if ~exist('extra_lag_flag','var')
    extra_lag_flag=0;
end
get_aux(orig_modfile,ramsey_type,instrument1,instrument2,extra_lag_flag,discount_factor_name);


% editing _aux.mod to have all policy vars
f_aux = fopen([aux_modfile,'.mod'],'r');
f_tempaux = fopen(['temp_',aux_modfile,'.mod'],'w');

% strategy is to copy file under f_aux (open for reading) to f_tempaux
% (open for writing).
if strcmp(ramsey_type,'ramsey')
    while ~feof(f_aux)  
        s = fgets(f_aux);
        
        if ~model_block && strncmp(s,'model;',6)
            model_block = 1;
            fprintf(f_tempaux,'%s \n','var Welf1, Welf2;')
            fprintf(f_tempaux,'%s \n','parameters omega_welf1, omega_welf2, nbeta;')
            
        elseif model_block
            if strncmp(s,'//Welf1 =',9) || strncmp(s,'//Welf1=',8) || strncmp(s,'//Welf2 =',9) || strncmp(s,'//Welf2=',8)
                s = strrep(s,'//','');
            end;
        end;
        
        fprintf(f_tempaux,'%s',s);
    end;
else
    while ~feof(f_aux)
        s = fgets(f_aux);
        
        if ~model_block && strncmp(s,'model;',6)
            model_block = 1;
            fprintf(f_tempaux,'%s \n','var Util, Welf;')
        elseif model_block
            if strncmp(s,'//Util =',8) || strncmp(s,'//Util=',7) || strncmp(s,'//Welf =',8) || strncmp(s,'//Welf=',7)
                s = strrep(s,'//','');
            end;
        end;
        
        fprintf(f_tempaux,'%s',s);
    end;
    
end;

fclose('all');
%copyfile(['temp_',aux_modfile,'.mod '],[aux_modfile,'.mod '])
if isunix
    eval(['!cp temp_',aux_modfile,'.mod ',aux_modfile,'.mod']);
else
    eval(['!copy temp_',aux_modfile,'.mod ',aux_modfile,'.mod']);
end
delete(['temp_',aux_modfile,'.mod']);
 

eval(['dynare ',orig_modfile,'_aux noclearall']);
M_orig_aux_ = M_;
oo_orig_aux_ = oo_;

% Making sure we will put shadow_ prefix before all policy variables
for i = 1:size(policy_vars,1)
    in_endo_names = 0;
    while in_endo_names == 0
        for j = 1:M_orig_aux_.endo_nbr
            if strcmp(policy_vars{i,:},deblank(M_orig_aux_.endo_names(j,:)))
                in_endo_names = 1;
            end;
            if in_endo_names == 0 && j == M_orig_aux_.endo_nbr
                M_orig_aux_.endo_names = strvcat(M_orig_aux_.endo_names, policy_vars{i,:});
                M_orig_aux_.endo_nbr = M_orig_aux_.endo_nbr + 1;
                in_endo_names = 1;
            end;            
        end;        
    end;    
end;




% copy model block from _aux.mod into a temporary file that will be read
% into the final shadow_... .mod file.
model_block = 0;
f_aux = fopen([aux_modfile,'.mod'],'r');
f_shadowtemp = fopen(['shadowtemp_',aux_modfile,'.mod'],'w');

while ~feof(f_aux)
    s = fgets(f_aux);
    
    if ~model_block && strncmp(s,'model;',6)
        model_block = 1;
        fprintf(f_shadowtemp,'%s',s);        
    elseif model_block
        if strncmp(s,'end;',4)
            model_block = 0;
        else
            fprintf(f_shadowtemp,'%s',s);
        end;            
    end;
end;

fclose('all');

% add shadow_ variables to var block in ramsey file, save as shadow_... .mod
f_ramsey = fopen([ramsey_modfile,'.mod'],'r');
f_shadowramsey = fopen(['shadow_',ramsey_modfile,'.mod'],'w');
f_shadowtemp = fopen(['shadowtemp_',aux_modfile,'.mod'],'r');

while ~feof(f_ramsey);
    s = fgets(f_ramsey);
    
    if strncmp(s,'model;',6)
        
        % first print out shadow variables in a var block        
        fprintf(f_shadowramsey,'%s','var ');
        for i = 1:M_orig_aux_.endo_nbr
            if i ~= M_orig_aux_.endo_nbr
                fprintf(f_shadowramsey,'%s \n',['shadow_',M_orig_aux_.endo_names{i},', ']);
            else
                fprintf(f_shadowramsey,'%s \n',['shadow_',M_orig_aux_.endo_names{i},'; ']);
            end;
        end;
%         

        
        % delimit shadow model block
        fprintf(f_shadowramsey,'%s \n','');
        fprintf(f_shadowramsey,'%s \n','// SHADOW BLOCK');
        fprintf(f_shadowramsey,'%s \n','');
        
        % add shadow var equations from _aux.mod
        while ~feof(f_shadowtemp)
            t = fgets(f_shadowtemp);
            token = tokenize_string(t,delimiters);
            
            for i = 1:size(token,2)
                for j = 1:M_orig_aux_.endo_nbr
                    if strcmp(token{1,i},deblank(M_orig_aux_.endo_names(j,:)))
                        token{1,i} = ['shadow_',token{1,i}];
                    end;
                end;

            end;
            
            shadow_t = strmerge(token);
            fprintf(f_shadowramsey,'%s \n',shadow_t);
            
        end;
        
        % delimit regular model block
        fprintf(f_shadowramsey,'%s \n','');
        fprintf(f_shadowramsey,'%s \n','// REGULAR BLOCK');
        fprintf(f_shadowramsey,'%s \n','');        
    elseif strcmp(deblank(s),'//steady;') || strcmp (deblank(s),'//check;')
        s = strrep(s,'//','');
        fprintf(f_shadowramsey,'%s',s);
    else
        fprintf(f_shadowramsey,'%s',s);
    end;
    
end;


%% _steadystate.m file

% collect original steady states
eval(['dynare ',ramsey_modfile,' noclearall nostrict']);
vec1 = oo_.dr.ys;

eval(['dynare ',aux_modfile,' noclearall nostrict;']);
vec2 = oo_.dr.ys;

ys = [vec1; vec2];

% start to write file
f_shadowramsey_ss = fopen(['shadow_',ramsey_modfile,'_steadystate.m'],'w');

fprintf(f_shadowramsey_ss,'%s \n',['function [optys,params,check] = shadow_',...
                                    ramsey_modfile,'_steadystate(ys,exo,M_,options_)']);
fprintf(f_shadowramsey_ss,'%s \n','');
fprintf(f_shadowramsey_ss,'%s \n','global params paramlist param_switch sensitivity_switch paramfile_name');
fprintf(f_shadowramsey_ss,'%s \n','check = 0;');
fprintf(f_shadowramsey_ss,'%s \n','');
fprintf(f_shadowramsey_ss,'%s \n',[ramsey_modfile,'_paramfile']);
fprintf(f_shadowramsey_ss,'%s \n','');
fprintf(f_shadowramsey_ss,'%s \n','nparams = size(M_.param_names,1);');
fprintf(f_shadowramsey_ss,'%s \n','for icount = 1:nparams');
fprintf(f_shadowramsey_ss,'%s \n','    eval([''M_.params(icount) = '',M_.param_names{icount},'';''])');
fprintf(f_shadowramsey_ss,'%s \n','end');
fprintf(f_shadowramsey_ss,'%s \n','');
fprintf(f_shadowramsey_ss,'%s \n','params=M_.params;');

fprintf(f_shadowramsey_ss,'%s \n','');

% write initial steady-state guess
format = '%12.14f';
fprintf(f_shadowramsey_ss,'%s \n',['ys = [',num2str(ys(1), format),';']);

for i = 2:length(ys)-1
    fprintf(f_shadowramsey_ss,'%s \n',[num2str(ys(i), format),';']);
end;

fprintf(f_shadowramsey_ss,'%s \n',[num2str(ys(length(ys))),'];']);
fprintf(f_shadowramsey_ss,'%s \n','');

% use fsolve get steady-state to solve .mod file

fprintf(f_shadowramsey_ss,'%s \n','nexog = M_.exo_nbr;');
fprintf(f_shadowramsey_ss,'%s \n','exog = zeros(nexog,1);');
fprintf(f_shadowramsey_ss,'%s \n','options = optimset(''Display'',''Iter'',''Jacobian'',''on'',''MaxFunEvals'',1e4,''MaxIter'',1e4,''TolFun'',1e-9,''TolX'',1e-9,''Algorithm'',''trust-region-reflective'');');
fprintf(f_shadowramsey_ss,'%s \n','algo = char(''active-set'', ''trust-region-reflective'', ''interior-point'', ''interior-point-convex'', ''levenberg-marquardt'', ''trust-region-dogleg'', ''lm-line-search'',''sqp'');');
fprintf(f_shadowramsey_ss,'%s \n','options = optimset(''Display'',''Iter'',''Jacobian'',''on'',''MaxFunEvals'',1e4,''MaxIter'',1e4,''TolFun'',1e-9,''TolX'',1e-9,''Algorithm'',algo(5,:));');
fprintf(f_shadowramsey_ss,'%s \n',['optys=fsolve(@(ys) shadow_',ramsey_modfile,'.static(ys,exog,M_.params),ys,options);']);

fclose('all');
rmdir([orig_modfile,'_aux'],'s');
delete([orig_modfile,'_aux*']);
delete('shadowtemp_*');
