tfest() failed with higher order

Im using tfest() function to fit my input and output data. Input is an impulse signal, and ouput is an curve with a one frequency resonance + a random noise.

I using tfest(data, poles, zeros, NaN) with poles and zeros increase. I thought that with higher order the model would fit better, hower it is not.

Here is my fit table. From 1,0 --> 5,4 seems as expected, but after that the fit increase and descrease without a reasons and when poles = 10 and zeros = 9 the fit is only 3.03. I dont understand the reason for this. Im really appriciate any help to explain this behaviours

P/s: Here is my fit with remove the delay time. tfest(data, poles, zeros, 0). The same happens

Answer 1

[Star Strider]

One way to estimate the system order with reasonably accuracy from the output signal itself is to calculate the Fourier transform of the signal using the fft function. Plot the imaginary part of the first half of the fft output as a funciton of frequency. (Zero-padding the signal using the second 'nfft' input to fft, ideally setting that to an integer power of 2 using the nextpow2 function, will increase the frequency resolution and may improve the accuracy of the result.)

The result should resolve itself into a series of curves that look like tangent functions. The poles are the frequencies where the signal approaches infinity, and the zeros are the frequency locations where the signal crosses zero. That should give you a reasonably good idea of the acutal system order. There may be an additional zero or pole at the origin, and an additional zero or pole at infinity, so setting the order to the observed number of poles and zeros and adding 1 or 2 to that number may be necessary. You will have to experiment.

This example does not show the effect well since it does not depict a transfer function, (it would help to have your signal), however it explains the general idea.

Clarification -- As a general rule, this approach is to calculate the Fourier transform of the transfer function, that being the Fourier transform of the output divided element-wise by the Fourier transform of the input. Here, the input is an impulse function and its Fourier transform is identically unity for all frequencies, so here the output is essentially the transfer fnction.

Try something like this --

Fs = 1000;
L = 50;
t = linspace(0, Fs*L, Fs*L+1).'/Fs;
s = sum(exp(-0.1*t).*sin(t*[1 10 100 200 300 450]*2*pi), 2);

figure
plot(t, s)
grid
xlabel('Time')
ylabel('Amplitude')

[FTs1,Fv] = FFT1(s,t);
Sizes = [size(s); size(FTs1)]

Sizes = 2×2

50001 1 32769 1

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

figure
plot(Fv, imag(FTs1))
grid
xlabel('Frequency')
ylabel('Imaginary Part Of Fourier Transform')

figure
plot(Fv, imag(FTs1))
grid
xlabel('Frequency')
ylabel('Imaginary Part Of Fourier Transform')
xlim([0 20])

function [FTs1,Fv] = FFT1(s,t)
% One-Sided Numerical Fourier Transform
% Arguments:
%    s: Signal Vector Or Matrix
%    t: Associated Time Vector 
    t = t(:);
    L = numel(t);
    if size(s,2) == L
        s = s.';
    end
    Fs = 1/mean(diff(t));
    Fn = Fs/2;
    NFFT = 2^nextpow2(L);
    FTs = fft((s - mean(s)) .* hann(L).*ones(1,size(s,2)), NFFT)/sum(hann(L));
    Fv = Fs*(0:(NFFT/2))/NFFT;
%     Fv = linspace(0, 1, NFFT/2+1)*Fn;
    Iv = 1:numel(Fv);
    Fv = Fv(:);
    FTs1 = FTs(Iv,:);
end

EDIT -- (21 Jul 2025 at 22:27)

Added 'Clarification' above. The rest is unchanged.

.