diff --git a/matlab/mcforecast3.m b/matlab/+conditional_forecasts/get_shock_paths.m
similarity index 96%
rename from matlab/mcforecast3.m
rename to matlab/+conditional_forecasts/get_shock_paths.m
index 78ff228763fbe85fddc0deb2a59d3a65c8880f37..d8434f3c212393316170a3bed82cbf12691d75cd 100644
--- a/matlab/mcforecast3.m
+++ b/matlab/+conditional_forecasts/get_shock_paths.m
@@ -1,5 +1,5 @@
-function [forcs, e] = mcforecast3(cL, H, mcValue, shocks, forcs, T, R, mv, mu)
-% [forcs, e] = mcforecast3(cL, H, mcValue, shocks, forcs, T, R, mv, mu)
+function [forcs, e] = get_shock_paths(cL, H, mcValue, shocks, forcs, T, R, mv, mu)
+% [forcs, e] = get_shock_paths(cL, H, mcValue, shocks, forcs, T, R, mv, mu)
 % Computes the shock values for constrained forecasts necessary to keep
 % endogenous variables at their constrained paths
 %
diff --git a/matlab/plot_icforecast.m b/matlab/+conditional_forecasts/plot.m
similarity index 89%
rename from matlab/plot_icforecast.m
rename to matlab/+conditional_forecasts/plot.m
index da4d88ebbd9eb576597a16231ec5c5883f9afb91..562d60a1868e6c175bd51d7a6ad297176c87d157 100644
--- a/matlab/plot_icforecast.m
+++ b/matlab/+conditional_forecasts/plot.m
@@ -1,4 +1,5 @@
-function plot_icforecast(Variables,periods,options_,oo_)
+function plot(Variables,periods,options_,oo_)
+% plot(Variables,periods,options_,oo_)
 % Build plots for the conditional forecasts.
 %
 % INPUTS
@@ -11,7 +12,7 @@ function plot_icforecast(Variables,periods,options_,oo_)
 %  None.
 %
 % SPECIAL REQUIREMENTS
-%  This routine has to be called after imcforecast.m.
+%  This routine has to be called after conditional_forecasts.run
 
 % Copyright © 2006-2023 Dynare Team
 %
@@ -44,7 +45,7 @@ else
 end
 
 if periods>forecast_periods
-    fprintf('\nplot_icforecast:: Number of periods for plotting exceeds forecast horizon. Setting periods to %d.\n',forecast_periods-1)
+    fprintf('\nconditional_forecasts.plot:: Number of periods for plotting exceeds forecast horizon. Setting periods to %d.\n',forecast_periods-1)
     periods=forecast_periods;
 end
 
diff --git a/matlab/imcforecast.m b/matlab/+conditional_forecasts/run.m
similarity index 90%
rename from matlab/imcforecast.m
rename to matlab/+conditional_forecasts/run.m
index e202943dc475d8bd487907b267a076db7f3820b0..9737147557f31084ae9952ba9578fdc8bc84864d 100644
--- a/matlab/imcforecast.m
+++ b/matlab/+conditional_forecasts/run.m
@@ -1,11 +1,16 @@
-function imcforecast(constrained_paths, constrained_vars, options_cond_fcst)
-
+function forecasts=run(M_,options_,oo_,bayestopt_,estim_params_,constrained_paths, constrained_vars, options_cond_fcst)
+% run(constrained_paths, constrained_vars, options_cond_fcst)
 % Computes conditional forecasts.
 %
 % INPUTS
-% - constrained_paths  [double]        m*p array, where m is the number of constrained endogenous variables and p is the number of constrained periods.
-% - constrained_vars     [integer]     m*1 array, indices in M_.endo_names of the constrained variables.
-% - options_cond_fcst    [structure]   containing the options. The fields are:
+%   M_                  [structure]     describing the model
+%   options_            [structure]     describing the options
+%   oo_:                [structure]     storing the results
+%   bayestopt_          [structure]     describing the priors
+%   estim_params_       [structure]     characterizing parameters to be estimated
+% - constrained_paths   [double]        m*p array, where m is the number of constrained endogenous variables and p is the number of constrained periods.
+% - constrained_vars    [integer]       m*1 array, indices in M_.endo_names of the constrained variables.
+% - options_cond_fcst   [structure]     containing the options. The fields are:
 %
 %                                      + replic              [integer]   scalar, number of monte carlo simulations.
 %                                      + parameter_set       [char]      values of the estimated parameters:
@@ -19,16 +24,16 @@ function imcforecast(constrained_paths, constrained_vars, options_cond_fcst)
 %                                      + conf_sig            [double]     scalar in [0,1], probability mass covered by the confidence bands.
 %
 % OUTPUTS
-% None.
+%  - forecasts          [structure]     results of the conditional forecast exercise
 %
 % SPECIAL REQUIREMENTS
 % This routine has to be called after an estimation statement or an estimated_params block.
 %
 % REMARKS
 % [1] Results are stored in oo_.conditional_forecast.
-% [2] Use the function plot_icforecast to plot the results.
+% [2] Use the function conditional_forecasts.plot to plot the results.
 
-% Copyright © 2006-2020 Dynare Team
+% Copyright © 2006-2023 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -45,8 +50,6 @@ function imcforecast(constrained_paths, constrained_vars, options_cond_fcst)
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-global options_ oo_ M_ bayestopt_ estim_params_
-
 if ~isfield(options_cond_fcst,'parameter_set') || isempty(options_cond_fcst.parameter_set)
     if isfield(oo_,'posterior_mode')
         options_cond_fcst.parameter_set = 'posterior_mode';
@@ -123,7 +126,6 @@ if estimated_model
     missing_value = dataset_info.missing.state;
 
     %store qz_criterium
-    qz_criterium_old=options_.qz_criterium;
     options_=select_qz_criterium_value(options_);
     options_smoothed_state_uncertainty_old = options_.smoothed_state_uncertainty;
     [atT, ~, ~, ~,ys, ~, ~, ~, ~, ~, ~, ~, ~, ~,oo_,bayestopt_] = ...
@@ -155,7 +157,6 @@ if estimated_model
     trend = constant(oo_.dr.order_var,:);
     InitState(:,1) = atT(:,end);
 else
-    qz_criterium_old=options_.qz_criterium;
     if isempty(options_.qz_criterium)
         options_.qz_criterium = 1+1e-6;
     end
@@ -246,7 +247,7 @@ FORCS1_shocks = zeros(n1,cL,options_cond_fcst.replic);
 for b=1:options_cond_fcst.replic %conditional forecast using cL set to constrained values
     shocks = sQ*randn(ExoSize,options_cond_fcst.periods);
     shocks(controlled_varexo,:) = zeros(n1, options_cond_fcst.periods);
-    [FORCS1(:,:,b), FORCS1_shocks(:,:,b)] = mcforecast3(cL,options_cond_fcst.periods,constrained_paths,shocks,FORCS1(:,:,b),T,R,mv, mu);
+    [FORCS1(:,:,b), FORCS1_shocks(:,:,b)] = conditional_forecasts.get_shock_paths(cL,options_cond_fcst.periods,constrained_paths,shocks,FORCS1(:,:,b),T,R,mv, mu);
     FORCS1(:,:,b)=FORCS1(:,:,b)+trend; %add trend
 end
 if max(max(max(abs(bsxfun(@minus,FORCS1(constrained_vars,1+1:1+cL,:),trend(constrained_vars,1:cL)+constrained_paths)))))>1e-4
@@ -292,7 +293,7 @@ FORCS2(:,1,:) = repmat(InitState,1,options_cond_fcst.replic); %set initial stead
 for b=1:options_cond_fcst.replic %conditional forecast using cL set to 0
     shocks = sQ*randn(ExoSize,options_cond_fcst.periods);
     shocks(controlled_varexo,:) = zeros(n1, options_cond_fcst.periods);
-    FORCS2(:,:,b) = mcforecast3(0,options_cond_fcst.periods,constrained_paths,shocks,FORCS2(:,:,b),T,R,mv, mu)+trend;
+    FORCS2(:,:,b) = conditional_forecasts.get_shock_paths(0,options_cond_fcst.periods,constrained_paths,shocks,FORCS2(:,:,b),T,R,mv, mu)+trend;
 end
 
 mFORCS2 = mean(FORCS2,3);
@@ -307,8 +308,4 @@ for i = 1:EndoSize
     forecasts.uncond.ci.(M_.endo_names{oo_.dr.order_var(i)}) = [tmp(t1,:)' ,tmp(t2,:)' ]';
 end
 forecasts.graph.title = graph_title;
-forecasts.graph.fname = M_.fname;
-
-%reset qz_criterium
-options_.qz_criterium=qz_criterium_old;
-oo_.conditional_forecast = forecasts;
+forecasts.graph.fname = M_.fname;
\ No newline at end of file
diff --git a/preprocessor b/preprocessor
index cf45c77343c50cddf89ab848b4300932e3bdf673..3dadac8f191dfa1dd660442d5bc4526c4c218149 160000
--- a/preprocessor
+++ b/preprocessor
@@ -1 +1 @@
-Subproject commit cf45c77343c50cddf89ab848b4300932e3bdf673
+Subproject commit 3dadac8f191dfa1dd660442d5bc4526c4c218149
diff --git a/tests/conditional_forecasts/3/fs2000_conditional_forecast_initval.mod b/tests/conditional_forecasts/3/fs2000_conditional_forecast_initval.mod
index 4bd3e6ae063906bfde4ea9e7494ec163eb1a37f1..5a061b3e13da6e6c70d05a65f8db7a2903f0e652 100644
--- a/tests/conditional_forecasts/3/fs2000_conditional_forecast_initval.mod
+++ b/tests/conditional_forecasts/3/fs2000_conditional_forecast_initval.mod
@@ -59,7 +59,7 @@ var e_a; stderr 0.014;
 var e_m; stderr 0.005;
 end;
 
-steady;
+steady(tolx=1e-10,tolf=1e-12);
 
 check;