Interpolating between two filters

[KG_is_back]

Here's the thing... I have two filters (sets of coefficients) and I need A->B morphing. Something like wet/dry mixing, except this would be wet1/wet2. Problem is, I need to have only one peace of code to preform the filter computations. How do I interpolate the coefficients?

[trogluddite]

I'm assuming that these are regular bi-quads?

You could just try linear interpolating all of the co-efficents, but I would do it with your speakers turned down! The relationship between co-efficients is pretty complex, so there's a chance that the filter becomes unstable part way through the fade.
MyCo's Filter Construction Kit might be worth a look. I don't understand a lot of it - but he has transformations there for the poles and zeros stuff on the 'Z-plane'. Interpolating the positions of the poles and zeros should ensure that the poles don't stray outside the "stable" zone. Nice sound too from what I remember of the old 'Z-plane morphing' synths of years gone by.

[KG_is_back]

I'm assuming that these are regular bi-quads?

That's the problem... they aren't. They are filters of order upto 6 (downto 0 which is basically a gain) and they may be of different orders too.
I have just checked myco's thing and the thing is truly complex. Basically I have to add the two Hilbert transforms, which would take more CPU then just having the two filters running separate and scaling the gain.
I thought I would save CPU by having only one filter with slowly changing coefs - the thing will be running in poly, so bypass options are very limited.
The thing is a physical modeling of a string. I have a set of 3 filter coefficients (one for the attack, one for the decay and one for sort of transition between the two). I want the filter to slowly change from response 1 through response 2 to response 3. It seems most CPU optimal solution will be to load response 1 and 2 morph from 1 to 2 and when filter is fully 2, change the response 1 to response 3 and morph slowly back. Hope it makes cents.
Working memin would be very helpful in this case... I will probably use your shared mem system to have direct access to those coefs...

[martinvicanek]

If they are biquads you can use linear interpolation for the coefficients without the filter going unstable. I have done this in my vowel filter over at SM. Anyway, Trog's recommendation to turn the speakers off during the test phase is well justified.

[martinvicanek]

I'm assuming that these are regular bi-quads?

That's the problem... they aren't. They are filters of order upto 6 (downto 0 which is basically a gain) and they may be of different orders too.

You can always factorize them (offline) into biquads, though.

[Exo]

I'm assuming that these are regular bi-quads?

That's the problem... they aren't. They are filters of order upto 6 (downto 0 which is basically a gain) and they may be of different orders too.
I have just checked myco's thing and the thing is truly complex. Basically I have to add the two Hilbert transforms, which would take more CPU then just having the two filters running separate and scaling the gain.
I thought I would save CPU by having only one filter with slowly changing coefs - the thing will be running in poly, so bypass options are very limited.
The thing is a physical modeling of a string. I have a set of 3 filter coefficients (one for the attack, one for the decay and one for sort of transition between the two). I want the filter to slowly change from response 1 through response 2 to response 3. It seems most CPU optimal solution will be to load response 1 and 2 morph from 1 to 2 and when filter is fully 2, change the response 1 to response 3 and morph slowly back. Hope it makes cents.
Working memin would be very helpful in this case... I will probably use your shared mem system to have direct access to those coefs...

The next release has support for Mem input in code and ASM. Although the Mem is copied on every voice. which makes it quite inefficient for large mems. Malc said he would look into working with the mem address directly instead of copying in the release after this.

I managed to get him to add a Mem to address primitive though . So we now have efficient Mem input to ASM.

[KG_is_back]

If they are biquads you can use linear interpolation for the coefficients without the filter going unstable. I have done this in my vowel filter over at SM. Anyway, Trog's recommendation to turn the speakers off during the test phase is well justified.

that's why I first test this sort of stuff connected to scope or readout module.

You can always factorize them (offline) into biquads, though.

I was always wondering how that is done.

[martinvicanek]

You can always factorize them (offline) into biquads, though.

I was always wondering how that is done.

Consider the Formant Filter from musicdsp.org (ported directly to FS):

Code: Select all

// Formant Filter
// http://www.musicdsp.org/showone.php?id=110

streamin in;
streamout out;

float out1, out2, out3, out4, out5, out6, out7, out8, out9, out10;

// vowel "A" coefficients:
float b0 = 8.11044e-06,
a1 = -8.943665402, a2 = 36.83889529, a3 = -92.01697887, a4 = 154.337906, a5 = -181.6233289,
a6 = 151.8651235, a7 = -89.09614114, a8 = 35.10298511, a9 = -8.388101016, a10 = 0.923313471;

out = b0*in - a1*out1 - a2*out2 - a3*out3 - a4*out4 - a5*out5 -
              a6*out6 - a7*out7 - a8*out8 - a9*out9 - a10*out10;

out10 = out9;
out9 = out8;
out8 = out7;
out7 = out6;
out6 = out5;
out5 = out4;
out4 = out3;
out3 = out2;
out2 = out1;
out1 = out;

It is a 10th order allpole filter. Its direct implementation in FS' single precision is unstable. However, we can fix it by breaking it up in five 2nd order sections.

The transfer function of the formant filter is:

Code: Select all

                       b0
H(z) = ---------------------------------------  (1)
       1 + a1*z^-1 + a2*z^-2 + ... + a10*z^-10

We want to write it as:

Code: Select all

                       b0
H(z) = ---------------------------------------,  (2)
       (1 - z1/z)*(1 + z2/z)* ... *(1 - z10/z)

where the z1, z2, etc. are the poles of the filter. They are the roots of the denominator polynomial in equation (1). They can be real or complex conjugate pairs. There is a very useful root finder at http://www.akiti.ca/PolyRootRe.html from which we get the following result:

zeroes.png (61.83 KiB) Viewed 32373 times

Observe that the ten roots come in five complex conjugate pairs:

Code: Select all

z1 = x + iy, z2 = x - iy, etc.

in equation (2) each such pair yields a section with real coefficients p and q:

Code: Select all

(1 - z1/z)*(1 - z2/z) = 1 + p*z^-1 + q*z^-2

with

Code: Select all

p = -2*x and q = x^2 + y^2.

The transfer function can be written as a product of 2nd order sections:

Code: Select all

H(z) = H1(z)*H2(z)* ... *H5(z)

where

Code: Select all

                b
H1(z) = ------------------- etc.
        1 + p*z^-1 + q*z^-2

All that remains is to implement the 2nd order sections:

Code: Select all

// 2nd Order Section
streamin in;
streamin b;
streamin p;
streamin q;
streamout out;

float out1, out2;

out = b*in - p*out1 - q*out2;

out2 = out1;
out1 = out;

Below is a schematic for demonstration.

[KG_is_back]

Thanks martin, Now I finally understand the thing... and here comes the hard part... implementation in FS.

My schematic uses Least-squares method to find the filter coefficients from known input and output wave
(the equation of a filter is linear: y=x0*b0+x1*b1+...-y1*a1-y2*a2...). It is something like deconvolution to create IR for FIR filter, however this creates IIR filter.

One more question. I what to implement nonlinearity in the algorithm too by adding c0*htan(x0)+c1*htan(x1)+... and various other functions like sqrt(abs(x)) x^2 etc. Will this factorization still work?

[martinvicanek]

I want to implement nonlinearity in the algorithm too by adding c0*htan(x0)+c1*htan(x1)+... and various other functions like sqrt(abs(x)) x^2 etc. Will this factorization still work?

Probably not in a rigorous sense. If you are lucky it might work approximately, but it really depends on the details.