How to recover the response to a step forcing upon system identification using impulseest?

I am unable to recover the response to a step forcing. The recreated signal reaches a level lower than the actual one (phi_asym), and then it is constant in time for the rest of the time. Otherwise, before reacing the plateau, the reconstruction seems very good. I also wonder why the first nonzero value of the step forcing (mystep) is supposed to correspond with a still zero value of the response (output).

t = 0:2^8; 
dt = t(2)-t(1);
tau = 20; % time scale of relaxation in years
phi_asym = 4; % relative increase of the observable phi t -> infinity

sig = 0.; % noise strength corrupting the signal

phi_r2s = phi_asym*(1 - exp(-t/tau)); % uncorrupted

output = [zeros(1,length(phi_r2s)) phi_r2s+sig*randn(1,length(phi_r2s))]; % corrupted

mystep = [zeros(1,length(phi_r2s)) ones(1,length(phi_r2s))];

figure; plot(1:length(mystep),mystep,1:length(output),output)

io_r2s = iddata(output',mystep',dt);

sys = impulseest(io_r2s);

[y, ~, ~, ysd] = step(sys,t);

figure; plot(t, y, 'b', t, y+3*ysd, 'b:', t, y-3*ysd, 'b:', t, output(end-length(t)+1:end))

Answer 1

[Tamas Bodai]

Perhaps the default order of the impulse response model is 70, which is where the plateau suddenly starts. There is no mention of a default value in the Help article for impulseest.

With noise corrupting the response (e.g. sig = 0.1), it appears that an order somewhat smaller than the time series length gives a better result in ways of "smoothing", as opposed to taking the maximally allowed value.