Invert phase response

[juha_p]

I prepared a 2nd order Peak filter (fs=44.1kHz, fc=20kHz, Q=0.089084, Gain=12dB) using RBJ's cookbook formula (+ additional Q pre-warping) which results these responses:

and coeffs:

Code: Select all

den = [-2.86512631333559 1.915296632059649 4.86512631333559]
num = [-14.38734500433902 1.915296632059649 16.38734500433902]

Phase response is wrong so, to turn the phase I did change the sign of the variable Q and got the right responses:

and coeffs:

Code: Select all

den = [4.86512631333559 1.915296632059649 -2.86512631333559]
num = [16.38734500433902 1.915296632059649 -14.38734500433902]

Is this the right way to turn the phase response?

EDIT: Changed the title

[Miles1981]

Are you sure you have not just inverted the filter because the plynomials were in the wrong order?

[mystran]

Phase response is wrong so, to turn the phase I did change the sign of the variable Q and got the right responses:

Most likely you are calculating the phase of the digital filter wrong (maybe calculating polynomials of z rather than z^-1?) and it comes out with an inverted sign. If the actual filter is stable (and it should be, if you got it from RBJ), then this is almost certainly the problem.

[juha_p]

Are you sure you have not just inverted the filter because the plynomials were in the wrong order?

How can this be checked?
I sure noticed the 'mirrored' coefficients but not sure what it means. Is the Octave showing the magnitude plot wrong for the other (inverted) filter? Mentioned additional Q pre-warping is tan() based so could this have some effect in this?

Most likely you are calculating the phase of the digital filter wrong (maybe calculating polynomials of z rather than z^-1?) and it comes out with an inverted sign. If the actual filter is stable (and it should be, if you got it from RBJ), then this is almost certainly the problem.

I'm using Octave.

Code: Select all

% Octave packages -------------------------------
pkg load signal

% filter parameters -----------------------------
fs = 44100;
fc = 20000;
Q  = 0.089084; 
db = 12.0;

%% Analog prototype (RBJ EQ Cookbook formulas) --
K  = 10^(db/40);
QA = K*Q;
wA0 = 2*pi*fc;
BA = [1 K^2*wA0/QA wA0^2];
AA = [1 wA0/QA wA0^2];

A = tf(BA,AA);

%% Q inversion ----------------------------------
%Q=Q*-1.0;

% Q pre-warping ---------------------------------
% use the one from https://www.google.fi/?gws_rd=ssl#q=RobClarkPhD.pdf

%% Digital filter (RBJ EQ Cookbook formulas) --
w0 = 2*pi*fc/fs;
cw0 = cos(w0);
sw0 = sin(w0);
alpha = sw0/(2*Q);

b0 = 1+alpha*K;
b1 = -2*cw0;
b2 = 1-alpha*K;
a0 = 1+alpha/K;
a1 = -2*cw0;
a2 = 1-alpha/K;

BD = [b2 b1 b0];
AD = [a2 a1 a0];      

RBJQ = tf(BD,AD,1/fs);

% bode plot -------------------------------------
wmin = 2 * pi * 20; % 20Hz
wmax = 2 * pi * ((fs/2.0) - (fs/2 - 20000)); % 20kHz

bode(A,'r', RBJQ, '--',{wmin, wmax});
legend('Analog','RBJ case Q w/ Q pre-warping','location', 'northwest');

[juha_p]

Is the script I'm using OK or is there something wrong with the analog prototype definition (i.e., is the phase response for analog prototype OK)?

[mystran]

While I'm not really an Octave user, I'm pretty sure you have your BD and AD vectors in wrong order.

They should almost certainly be [b0 b1 b2] and [a0 a1 a2], which should give the same results with inverted phase.

[juha_p]

While I'm not really an Octave user, I'm pretty sure you have your BD and AD vectors in wrong order.

They should almost certainly be [b0 b1 b2] and [a0 a1 a2], which should give the same results with inverted phase.

Yes, it does (it can be seen by comparing the coefficients of those two filters in my opening post).

Actually I was wondering if changing the sign of Q is correct (secure) way to do the inversion or is there something more happening in background one can't see in the plot (as like what Miles1981 suggested earlier ... if I sum those two filters and then divide by 2 it results coefficients which give flat response for both mag and phase (but the Octave shows mag plots equal).

[mystran]

Actually I was wondering if changing the sign of Q is correct (secure) way to do the inversion or is there something more happening in background one can't see in the plot

Just reverse the coefficient order... although it's worth noting that if you do this to stable poles it'll make them unstable (and doing it to minimum-phase zeroes will make them maximum-phase). If you want to sweep Q from stability to instability and back, then you might want to work in terms of 1/Q instead (which happens to be the band-pass feedback coefficient for damping a state-variable filter .. what a coincidence), since the latter goes through zero (at self-osc) while the former goes through the point at infinity.

The only "stable" method for inverting the phase is to reverse the input signal, run the filter and then reverse the output signal. The good news is that if you take a stable analog filter and manage to convert it into a stable digital filter and end up with inverted phase plot then that's just a matter of plotting inconsistency.

Actually, now that I think about it, I guess it's probably actually the analog filter that is going the wrong way, since the usual substitutions with s=i*w or z=exp(i*w) would give you negative phase-shift for stable filters in both cases.. but maybe Octave just likes to plot it the other way around, since I'd assume it's supposed to be Matlab compatible, the documentation says your analog filter is in the correct order and Matlab documentation says you're suppose to give the digital filter coefficients in [b0, b1 .. ] order... so like .. whatever.

Ultimately it doesn't really matter anyway, as long as you pick a consistent convention, since stable filters can only ever shift the phase (ie. delay the signal) one direction anyway and this direction happens to be the same for analog and digital filters because continuous and discrete time normally move in the same direction.

[juha_p]

Tänkjuu verimats mystran! Looks clearer now.

[juha_p]

When one prepare filters using coefficients from my opening post the other case is unstable (one pole and one zero are well out from the unit circle limits) and does not give output responses as like how it is shown above (I did try this using SE and VST Analyzer).
So, Octave must do something for the coefficients of the unstable filter before plotting (Bode plot in question). What is needed to do to there?

[mystran]

When one prepare filters using coefficients from my opening post the other case is unstable (one pole and one zero are well out from the unit circle limits) and does not give output responses as like how it is shown above (I did try this using SE and VST Analyzer).
So, Octave must do something for the coefficients of the unstable filter before plotting (Bode plot in question). What is needed to do to there?

You can plot unstable filters just fine if you do it from the transfer function, but tools like VST analyzer are trying to analyze the actual output of the filter instead and this doesn't work so very well if the filter isn't stable.

ps. From the point of view of transfer function analysis, the stability of a filter is essentially an "implementation detail" and unstable poles can even be implemented by reversing "time" (which is not terribly useful for real-time signals, but it would work fine for batch processing of recordings, images, whatever).

[juha_p]

OK, thanks for the info.