function outvec = mpc_simPostprocess(SimData, MPCobj, SimOptions, isexplicit, Tf, r, v, d, ynoise, unoise, mv_signal, simBar)
% Perform closed-loop simulation

% Authors:  A. Bemporad, N.L. Ricker, R. Chen
% Copyright 1994-2020 The MathWorks, Inc.

if simBar
    h = waitbar(0,'MPC simulation');
end
logical_ref_preview = logical(SimData.ref_preview);
logical_md_preview = logical(SimData.md_preview);
logical_simBar = logical(simBar);
% u and y contain offsets but they are scaled, xp and xmpc contain offsets
% lastu is offset-free and scaled, lastx is offset-free

if ~isexplicit
    SimData.H=SimData.Hessian; % his is because explicit MPC .H has overwritten the Hessian
end
Uscale = SimData.Uscale;
Yscale = SimData.Yscale;
mvindex = SimData.mvindex;
mdindex = SimData.mdindex;
MVscale = Uscale(mvindex);

[u,y,xp,xmpc,status,region] = mpcloop_engine(isexplicit, SimData.xp0, SimData.nxp, ...
    SimData.lastx, SimData.lastu,...
    SimData.nx, SimData.nxQP , SimData.nv, SimData.ndp, SimData.nmy, SimData.nmv, SimData.ny,...
    SimData.ypoff, SimData.yoff, SimData.xpoff, SimData.xoff,...
    SimData.dxpoff, SimData.mvoff, SimData.mdoff, SimData.vpoff, SimData.upoff,...
    SimData.myindex, r', v', d', ynoise', unoise', mv_signal', SimData.openloop,...
    SimData.Kx, SimData.Kr, SimData.Ku1, SimData.Kv, SimData.Kut, SimData.utarget,...
    SimData.degrees, SimData.isQP, logical_ref_preview,logical_md_preview,...
    SimData.CustomSolver, false, SimData.UseSuboptimalSolution, ...
    SimData.UseActiveSetSolver, SimData.ASOptions, SimData.IPOptions, SimData.ADOptions, ...
    SimData.MIQPOptions, SimData.isHyb, SimData.Mdis, SimData.Ndis, SimData.Vdis, SimData.numdis, SimData.maxdis, ...
    SimData.Mrows, SimData.Hinv, SimData.Linv, ...
    SimData.Mlim, SimData.Mx, SimData.Mu1, SimData.Mv,...
    SimData.Ac, SimData.p, SimData.A, SimData.Bu, SimData.Bv, SimData.Cm, SimData.Dvm,...
    SimData.L, SimData.M,...
    SimData.Ap, SimData.Bup, SimData.Bvp, SimData.Bdp, SimData.Cp, ...
    SimData.Dvp, SimData.Ddp, Tf, logical_simBar,...
    SimData.F, SimData.G, SimData.H, SimData.K, SimData.i1, SimData.i2, SimData.nr, SimData.np, ...
    SimData.xunfindex ,SimData.refunfindex, ...
    SimData.uunfindex, SimData.mdunfindex,...
    SimData.thmin, SimData.thmax,...
    MVscale,Yscale,Uscale(mdindex),...
    SimData.Jm,SimData.SuJm,SimData.Su1,SimData.Sx,SimData.Hv);
if simBar
    close(h);
end

%% Post-processing
ny = SimData.ny;
nmv = SimData.nmv;
nmd = SimData.nmd;
nud = SimData.nud;
ndp = SimData.ndp;
yoff = SimData.yoff;
mvoff = SimData.mvoff;
mdoff = SimData.mdoff;
% Scale to engineering units
y = y*diag(Yscale);
u = u*diag(MVscale);
% Possible plotting
Ts = MPCobj.Ts;
t = (0:Tf-1)'*Ts;
%% showOutputNoise option
if SimOptions.ShowOutputNoise
    y(:,SimData.myindex) = y(:,SimData.myindex) + ynoise*diag(Yscale(SimData.myindex));
end
%% Plot or return outputs
if nargout==0
    % The number of unmeasured disturbances may be different than nominal.
    if ndp~=nud
        % Change in number of UMDs.  Create a dummy MPC object for plotting.
        model = ss(zeros(ny,nmv+nmd+ndp),'TimeUnit',MPCobj.Model.Plant.TimeUnit);
        model = setmpcsignals(model,'MV',1:nmv,'MD',nmv+1:nmv+nmd,'UD',nmv+nmd+1:nmv+nmd+ndp);
        model.outputname = MPCobj.Model.Plant.OutputName;
        model.inputname(1:nmv) = MPCobj.Model.Plant.InputName(mvindex);
        model.inputname(nmv+1:nmv+nmd) = MPCobj.Model.Plant.InputName(mdindex);
        for i=1:ndp
            umdname = sprintf('UMD%d',i);
            model.inputname{nmv+nmd+i} = umdname;
        end
        MPCobj= mpc(model, Ts);
    end
    % Add offsets and scale.  Note:  d is NOT scaled.
    r = (r(1:Tf,:)+ones(Tf,1)*yoff')*diag(Yscale);
    v=(v(1:Tf,:)+ones(Tf,1)*mdoff')*diag(Uscale(mdindex));
    % plot
    if length(t)>1
        try
            plot(MPCobj,t,y,r,u,v,d);
        catch ME
            throw(ME);
        end
    else
        ctrlMsgUtils.error('MPC:computation:InvalidTfSim');
    end
else
    % return meta data in engineering units
    nxp = SimData.nxp;
    nxdist = SimData.nxdist;
    nxnoise = SimData.nxnoise;
    plantstates = 1:nxp;
    diststates = nxp+1:nxp+nxdist;
    noisestates = nxp+nxdist+1:nxp+nxdist+nxnoise;
    if SimData.openloop
        xmpc = [];
    else
        xmpc = struct('Plant',xmpc(:,plantstates),'Disturbance',xmpc(:,diststates),...
        'Noise',xmpc(:,noisestates),'LastMove',[((SimData.lastu(:)+mvoff).*MVscale)';u(1:end-1,:)]);
    end
    outvec = {y,t,u,xp,xmpc,SimOptions,status,region};
end
