This is the simple code to find transfer function between sigout and sigin signals and then are the filter coefficients estimated by invfreqz function.

%% load simlated data
load CFDvsABAQUS_TC.mat

%% simulated data preprocessing (grid orientation B)
time = TCb(:,1);
sigin = TCb(:,2);           % inlet temperature (TC signal by Abaqus)
sigout = TCb(:,3);          % outlet temperature (by cfd FENIX)
Fs = 1/mean(diff(time));    % sampling frequency

%% Transfer Function Estimate outlet->inlet

%% filter design by known transfer function
w = linspace(0,pi,length(fio)); % angular frequency f=[0,fs/2] -> w=[0,pi]
m = 500; % b coeffs
n = 500; % a coeffs
[boi,aoi] = invfreqz(toi,w,n,m);

%% final Digital filtering
esigin = filter(boi,aoi,sigout);

The input data are available at this dropbox link.

The final filter function produce completely wrong results (the result should be similar to original sigin signal). I have no idea what is wrong.

Any help will be very useful for me.

  • 1
    $\begingroup$ Probably all problems comes from extremely high order m=n=500 at invfreqz function... Some numerical instabilities??? But I need high order estimation to fit properly estimated transfer function. $\endgroup$ – michalkvasnicka Apr 28 '17 at 15:04

Your values for m and n are way too big. Matlab starts getting numerical problems for those values above 30, sometimes.

Try reducing them to more reasonable numbers (e.g. 10) and slowly increase them until invfreqz generates bogus outputs.

  • $\begingroup$ Yes, you are right! function "invfreqz" is highly numerical unstable for higher m and/or n values of filter order. I am afraid that this problem has no simple solution... $\endgroup$ – michalkvasnicka May 2 '17 at 6:22
  • $\begingroup$ @michalkvasnicka Just as a matter of engineering practice, lower-order models that satisfy your requirements are usually better than higher-order models because you end up just using the extra degrees of freedom to fit the noise, not the signal $\endgroup$ – Robert L. Nov 1 '18 at 18:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.