0
$\begingroup$

I've discovered that such convolution can be further sped because of symmetric coefficients but I'm not able to get it right.

The filter is a low pass filter at 11025Hz for 44100Hz, with 461 taps.

The convolution that skips half the taps, that works:

public float ProcessHalf(float input)
{
    var length = H.Length;
    var center = H.Length / 2;

    Z[State] = Z[State + length] = input;

    var output = 0.0f;

    for (var i = 1; i < length; i += 2)
    {
        output += H[i] * Z[State + length - i];
    }

    output += H[center] * Z[State + length - center];

    State--;

    if (State < 0)
    {
        State += length;
    }

    return output;
}

The convolution that attempts to further half iterations, does not work:

public float ProcessQuarter(float input)
{
    var length = H.Length;
    var center = H.Length / 2;

    Z[State] = Z[State + length] = input;

    var output = 0.0f;

    for (var i = 1; i < center; i += 2)
    {
        var j = State + length - i;
        var k = State + length + i - center; // this is not right

        var a = Z[j];
        var b = Z[k];

        output += H[i] * (a + b);
    }

    output += H[center] * Z[State + length - center];

    State--;

    if (State < 0)
    {
        State += length;
    }

    return output;
}

When I check the result, not only it's not filtered, but its volume is lower.

Can you tell where the problem is in the algorithm?

$\endgroup$
2
  • 3
    $\begingroup$ With 461 taps, and assuming that H.length =461, what do you expect var center to be? Besides that, you could try to take a smaller example like a filter with 11 taps and write down with a pen what the indexes of a symmetric filter should be and compare this with what you calculate for j and k. This will help you to debug your own code $\endgroup$ Commented Jul 19, 2023 at 11:49
  • $\begingroup$ Thanks for pointing that one out as well, I've added an answer explaining how to deal with both :) $\endgroup$
    – aybe
    Commented Jul 20, 2023 at 3:27

1 Answer 1

1
$\begingroup$

Thanks to @Irreducible tips, here's the answer on what the problem was but also on why center may be a problem in some cases.

The right offset to subtract from state is (center + (center - t)):

public float ProcessQuarter(float input)
{
    var length = H.Length;
    var center = length / 2;
    var offset = State + length;

    Z[State] = Z[offset] = input;

    var output = 0.0f;

    for (var t = 1; t < center; t += 2)
    {
        var l = offset - t;
        var r = offset - (center + (center - t));

        var a = Z[l];
        var b = Z[r];

        output += H[t] * (a + b);
    }

    output += H[center] * Z[offset - center];

    State--;

    if (State < 0)
    {
        State += length;
    }

    return output;
}

Now the center may be a problem depending it's odd or even, this little test shows how it either gets skipped or included when processing half of the taps:

var list = new[] { 11, 13, 15, 17, 19, 21, 23, 459, 461, 463 };

foreach (var taps in list)
{
    var center = taps / 2;
    var mirror = center % 2 == 0;

    Console.WriteLine($"{nameof(taps)}: {taps}, {nameof(center)}: {center}, {nameof(mirror)}: {mirror}");

    Console.WriteLine("\titerating 1 in 2 taps NOPE:");

    for (var i = 1; i < taps; i += 2)
    {
        Console.WriteLine($"\t\t{(i == center ? "\tcenter " : "")}{i}");
    }

    Console.WriteLine("\titerating 1 in 4 taps GOOD:");

    for (var i = 1; i < center; i += 2)
    {
        var j = i;
        var k = center + (center - i);
        Console.WriteLine($"\t\t{j} {k}");

        Assert.AreEqual(taps, j + k + 1);
    }

    Console.WriteLine($"\t\tTODO convolve center: {center}");
    Console.WriteLine();
}

For 11 taps:

taps: 11, center: 5, mirror: False
iterating 1 in 2 taps NOPE:
    1
    3
        center 5
    7
    9
iterating 1 in 4 taps GOOD:
    1 9
    3 7
    TODO convolve center: 5

For 13 taps:

taps: 13, center: 6, mirror: True
    iterating 1 in 2 taps NOPE:
        1
        3
        5
        7
        9
        11
    iterating 1 in 4 taps GOOD:
        1 11
        3 9
        5 7
        TODO convolve center: 6

For 15 taps:

taps: 15, center: 7, mirror: False
    iterating 1 in 2 taps NOPE:
        1
        3
        5
            center 7
        9
        11
        13
    iterating 1 in 4 taps GOOD:
        1 13
        3 11
        5 9
        TODO convolve center: 7

For 17 taps:

taps: 17, center: 8, mirror: True
    iterating 1 in 2 taps NOPE:
        1
        3
        5
        7
        9
        11
        13
        15
    iterating 1 in 4 taps GOOD:
        1 15
        3 13
        5 11
        7 9
        TODO convolve center: 8

Note however that isn't a problem when using symmetry as the loops evaluates t < center, i.e. it is never included and must be computed after the loop.

Here's the corrected code for center when processing 1 taps in 2:

public float ProcessHalf(float input)
{
    var length = H.Length;
    var center = length / 2;
    var offset = State + length;

    Z[State] = Z[offset] = input;

    var output = 0.0f;

    for (var t = 1; t < length; t += 2)
    {
        var l = offset - t;

        output += H[t] * Z[l];
    }

    if ((center & 1) == 0)
    {
        output += H[center] * Z[offset - center];
    }

    State--;

    if (State < 0)
    {
        State += length;
    }

    return output;
}
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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