Version of the Control Package, version of Octave
~/workspace/pkg-control > octave --no-gui --quiet --eval "disp(version()); pkg load control; pkg list control"
10.3.0
Package Name | Version | Installation directory
--------------+---------+-----------------------
control *| 3.6.1 | .../share/octave/octave_packages/control-3.6.1What happens?
Matched-Z discretiziaion for c2d returns slightly wrong uniform gain for high-pass filters (or any filter with poles or zeros at the origin).
What is the expected result?
Gains should probably match matlab's accuracy.
Steps to reproduce
I was working on my own filter library where I have an optional matched-z discretization method. I used Octave’s method from this package and didn’t think about it again. An AI reviewer caught an inconsistency with the apparent canonical DC gain matching for Matched-Z in my implementation. I wrote a quick test script to compare Octave’s behavior here against Matlab. Note that you only see this when designing a filter with zeros at s=0 (like a high-pass filter).
% Octave c2d matched gain vs MATLAB.
if exist ('OCTAVE_VERSION', 'builtin')
pkg load control
fprintf ('Running with octave\n');
else
fprintf ('Running with matlab\n');
end
sys = tf ([1 0], [1 1]); % high-pass filtr
sys_d = c2d (sys, 0.1, 'matched');
[num, den] = tfdata (sys_d, 'v');
fprintf ('k_d = %.10f\n', num(1) / den(1));Matlab
>> test_matched_z
Running with matlab
k_d = 0.9516258196
Octave (today)
~/workspace/pkg-control > octave --no-gui test_matched_z_matlab.m
octave: X11 DISPLAY environment variable not set
octave: disabling GUI features
Running with octave
k_d = 0.9524187090
Sure enough, they are different, not by a lot (~0.083% error at the 0.1 sampling period when taking Matlab as the source of truth).
What I think the bug is
The reason it differs, I think, is Octave’s decision to gain match along the real axis when it should be matching along the imaginary axis. This is probably clear to all authors but I had to think about it some. Gain matching in general is only required at all b/c the matched Z directly maps the poles/zeros from the continuous to discrete domains, but completely ignores scaling. So the idea should be to pick a frequency at DC, and scale the whole transfer function in the discrete domain to match the continuous domain’s scaling. If there are zeros at s=0 you can’t gain match at DC because the DC gain is blocked by the transfer function itself (i.e. H_c(0) = 0 aka a high-pass filter). You’d just end up with a DC gain of zero which isn’t exactly useful. That’s why you walk to the next non-zero frequency in order to get a useful scalar.
The implementation here today does:
But w (omega) should be defined in terms of the imaginary axis as it’s the frequency response, s = j*omega and z = e^(j*omega*T).
That’s also why you don’t see this small error for any low-pass filter because you evaluate H_c(0) = some useful number. So w_d = exp (w_c * tsam) = exp(0) = 1 which is the same as w_d = exp (j * w_c * tsam) = exp(j * 0 * T) = 1. No zeros at the origin means you stick with wc=0 so the math works out the same.
So when the implementation here walks the frequency to the next step it's walking along the real axis of the s-plane which AFAIK doesn’t have a physical interpretation.
Version of the Control Package, version of Octave
What happens?
Matched-Z discretiziaion for
c2dreturns slightly wrong uniform gain for high-pass filters (or any filter with poles or zeros at the origin).What is the expected result?
Gains should probably match matlab's accuracy.
Steps to reproduce
I was working on my own filter library where I have an optional matched-z discretization method. I used Octave’s method from this package and didn’t think about it again. An AI reviewer caught an inconsistency with the apparent canonical DC gain matching for Matched-Z in my implementation. I wrote a quick test script to compare Octave’s behavior here against Matlab. Note that you only see this when designing a filter with zeros at s=0 (like a high-pass filter).
Matlab
>> test_matched_z Running with matlab k_d = 0.9516258196Octave (today)
Sure enough, they are different, not by a lot (~0.083% error at the 0.1 sampling period when taking Matlab as the source of truth).
What I think the bug is
The reason it differs, I think, is Octave’s decision to gain match along the real axis when it should be matching along the imaginary axis. This is probably clear to all authors but I had to think about it some. Gain matching in general is only required at all b/c the matched Z directly maps the poles/zeros from the continuous to discrete domains, but completely ignores scaling. So the idea should be to pick a frequency at DC, and scale the whole transfer function in the discrete domain to match the continuous domain’s scaling. If there are zeros at s=0 you can’t gain match at DC because the DC gain is blocked by the transfer function itself (i.e. H_c(0) = 0 aka a high-pass filter). You’d just end up with a DC gain of zero which isn’t exactly useful. That’s why you walk to the next non-zero frequency in order to get a useful scalar.
The implementation here today does:
But w (omega) should be defined in terms of the imaginary axis as it’s the frequency response,
s = j*omegaandz = e^(j*omega*T).That’s also why you don’t see this small error for any low-pass filter because you evaluate
H_c(0) = some useful number. Sow_d = exp (w_c * tsam) = exp(0) = 1which is the same asw_d = exp (j * w_c * tsam) = exp(j * 0 * T) = 1. No zeros at the origin means you stick withwc=0so the math works out the same.So when the implementation here walks the frequency to the next step it's walking along the real axis of the s-plane which AFAIK doesn’t have a physical interpretation.