MATLAB: Instability of identified model in time domain with the step function

Control System ToolboxidentificationstepSystem Identification Toolbox

I'm trying to fit the frequency response of a system, with a high order model. I need good precision in a wide frequency range. I'm getting good results with high order models (using ARX or tfest) in the frequency domain, but the time response of the system, using the step function is unstable. My physical system is stable.
If I reduce the order of the fit, the step response became stable but the frequency fit is not good enough for what I have to do.
Does it comes from a limit of the step response ? Is my frequency fit wrong ?
t1 = linspace(0,length(out_sweep),length(out_sweep));
in_sweep = cell2mat(struct2cell(in_sweep));
out_sweep = cell2mat(struct2cell(out_sweep));
response_sweep=out_sweep*2e-5./in_sweep;
f_ax = linspace(1,fs,length(in_sweep)+1).';
f_ax(end)=[];
% Get a frd model from the raw frequency data
Gfrd = idfrd(permute(response_sweep,[2 3 1]),f_ax,0,'FrequencyUnit','Hz');
% First fit of the transfer function with a 32 poles 32 zeros model
Gfit = tfest(Gfrd,32,32);
bode(Gfrd, Gfit)
xlim([20,20e3])
legend('Measured','Estimated')
w = Gfrd.Frequency; %#ok<NASGU>
% Crop the frequency response
GfrdProcessed = Gfrd;
GfrdProcessed = fselect(GfrdProcessed,20,20e3);
% create a weight vector
w = GfrdProcessed.Frequency;
f = squeeze(GfrdProcessed.ResponseData);
idx = w<40;
f(idx) = movmean(f(idx),3);
% Low weight for high frequencies
Weight = ones(size(f));
idx = w>100;
Weight(idx) = Weight(idx)/10;
% High weight for low frequencies
idx = w<1e3;
Weight(idx) = Weight(idx)*30;
% Use the weighting option for the tfest identification
tfestOpt = tfestOptions('WeightingFilter',Weight);
% GfitNew is the sys estimation of my frequency data
GfitNew = tfest(GfrdProcessed,10,tfestOpt);
% Convert the sys to plotable data
[mag_gfrd,pha_gfrd,wout]=bode(Gfrd,2*pi*f_ax); %#ok<ASGLU>
mag_gfrd=squeeze(mag_gfrd);
pha_gfrd=squeeze(pha_gfrd);
[mag_gfitnew,pha_gfitnew,wout]=bode(GfitNew,2*pi*f_ax);
mag_gfitnew=squeeze(mag_gfitnew);
pha_gfitnew=squeeze(pha_gfitnew);
Bad_time_response.jpg
Good_frequency_fit.jpg

Best Answer

Can you post your data out_sweep to play with your code and come up with a solution to your problem?