Arbitrarily high-order IIR filters

@joslloand Have you ever played or wrote a “n-th order” IIR filter? It would be perfectly possible to write a) an array of size “order” for feedback coefff ; b) an array of size “order+1” for feedforward coeffs, and another array of the same size to store the internal state of the filter. It probably (for sure) will need to be written as a template class (c++) . With one code, we could reuse the code for any IIR filter.

Who wants to try this?

Nick Collins wrote a linear time-invariant filter a long time ago: LTI | SuperCollider 3.12.2 Help – this may be a good start.

hjh

1 Like

As it turns out, I’m unsure if a C++ template is the best idea for all cases. I thought, “What if I want to change the order of the filter at runtime”? Well, in this case, I don’t think that would work. My temporary solution was to define an unrealistic high-limit order (20 or maybe 10) and fill the remaining places with zeros. The JACK client, which can change the order of the filter at runtime, requires one to take care of concurrency; that’s why all the new stuff. Reviews warmly accepted.

filter-nth-order.cpp

Interesting resource: Filter Playground

It’s been a long time since I made this thread and I don’t currently need a multiband compressor, funnily enough I think if I did now, I would probably just read a few papers and code up something quite similar to what Nathan initially did in his (now deleted) initial post. It was a very helpful response from him, but I think the best advice I would want to give myself two years ago is that I should probably just give it a shot. Most of the DSP learning progress I’ve made in the past 8 months - considerably more than in the 3-4 years beforehand - came from just reading stuff and trying to implement it until I eventually got it right. Just thought this might be a helpful point for anyone who’s at the same overcautious point I was at back then.

[Edit: this was from the above linked Multiband Compressor thread]

@smoge I did some thinking about how to make an SC UGen for an nth order IIR filter, and I don’t think this is really possible, because AFAIK you need to have a fixed number of UGen arguments. However, IIUC, there isn’t a DSP difference between, eg. cascading two biquad filters and having a fourth order IIR filter. So I think for a lot of cases, you could just code up eg. biquad or fourth order and then cascade them.

Here is a first order IIR filter that I wrote for myself a few weeks ago. I think it’s correct, although I haven’t had anyone check it - I did a few tests and the results seem fine. I am pretty new to writing UGens so any feedback is very welcome. Adapting it for biquad or whatever shouldn’t be particularly difficult.

2 Likes

That’s a limitation of scsynth, but easily solved using Buffers, maybe?

Calculating second-order coefficients is not hard, but it becomes more involved with higher orders. Scipy can do this (or some linear algebra libraries), but I think it can be computationally intensive. That’s the main problem, I think. In this way, you’re right about using cascading biquads.

I don’t quite understand what you mean sorry. I just meant that you need arguments for all the coefficients.

… which is exactly how LTI does it. As far as I can see, the only thing that needs to be adjusted (perhaps) is to make it modulateable rather than “time invariant.” Otherwise, someone implemented this, many years ago.

hjh

2 Likes

My little application allows one to change coefficients at runtime Using shared pointers and atomic operations, including changing the order of the filter. It’s slightly different than a sc-plugin, with 64-bit signal processing and single-sample operation, but that’s easy to change. I’m not sure there is interest in turning it into a plugin. Maybe not.

Well, as I read the documentation for LTI, I see that it takes one buffer for feed-forward coefficients, and a second buffer for feedback coefficients.

So LTI is an IIR.

At this point, I don’t quite understand what is the disconnect here. The LTI documentation is very clear as to the formula being implemented. “Different type of filter” just isn’t true.

This isn’t the case, if buffers are used to supply the coefficients – as LTI already does.

It’s really quite strange to me – “here’s a thing that exists” and then the conversation goes on as if it doesn’t exist. Anyway, I won’t belabor the point; just read the help and the source code.

hjh

James, I think the documentation uses feedback terms involving previous output values. This feedback means the filter’s impulse response is theoretically infinite. At the same time, it’s a specific type of LTI system because it satisfies both linearity and time invariance. However, the documentation is not precise enough to explain not all LTI systems are IIR filters.

my implementation is different; coeffs can change, including their number, and processing is 64bit, single-sample processing. Also, it uses TDF-II for memory efficiency. For example, maybe a different precision might make a difference for that kind of processing.

and it’s ok it something similar already exists, right? it’s not a problem, even if I fail

I’m pretty sure James is right. The transfer equation in the documentation is the same the transfer equation of an IIR, unless I’m missing something. What is the difference?

Also - and I’m happy to be proven wrong - I doubt that whatever minuscule differences can be found in memory efficiency or whatever will actually make a difference for the sound or usability of a multiband compressor, right?

My bad - I skimmed over your post and didn’t actually check out the LTI documentation before commenting, taking Bernardo’s word for it that the filters are different.

Yes, the transfer (the equation precisely, not in the text) is recursive (The SC documentation equation is an IIR filter). Still, as I said before, the documentation could be more complete.

Again, the IIR filter is a specific type of LTI.

TDF-II is significantly more memory efficient, requiring half the delay elements than TDF-I.

My implementation can change coeffs and order at runtime without glitches or data races. Since it’s recursive, I thought 64-bit could be something to TRY.

I almost feel I need to apologize for coding this. It’s a simple filter, friends, why the troubles?

Things can be different and not be a dichotomy. A sc help file is not the canonical reference.

(I’m silencing this thread)

@smoge I apologise, I did not mean to disparage you in any way; my comment was purely a clarification.

I certainly don’t have a problem with you coding filters. I do think that, where there is potential misinformation on the forum, it’s important to clear up misunderstandings to help people learn.

1 Like

Well, so let’s read an actual book. In theory, those words are not the same.

However, if you see the code, you can do similar things with them (it seems that is the case).

However, the implementation is not quite the same regarding real-time safety, seamless coefficients without interrupting the audio stream, and TDF-II, which reduces memory usage and improves numerical stability. One implementation is more optimized for real-time audio processing.

Where is that SC community always at the edge? No more?

Have you considered unrolling the loop to test performance? I think biquad is a good candidate, I would try myself.

// IIR.sc
IIR : UGen {
    *ar { |input, a1, a2, b0, b1, b2|
        ^this.multiNew('audio', input, a1, a2, b0, b1, b2);
    }

    checkInputs {
        ^this.checkValidInputs;
    }
}





// PluginIIR.hpp
// Jordan White (jordanwhitede@gmail.com)

#pragma once

#include "SC_PlugIn.hpp"

namespace IIR {

class IIR : public SCUnit {
public:
    IIR();

private:
    // Calc function
    void next(int nSamples);

    // Member variables 
    float prevSample1;
    float prevSample2;
};

} // namespace IIR




// PluginIIR.cpp
// Jordan White (jordanwhitede@gmail.com)

#include "SC_PlugIn.hpp"
#include "IIR.hpp"

static InterfaceTable* ft;

namespace IIR {

IIR::IIR() {
    mCalcFunc = make_calc_function<IIR, &IIR::next>();
    prevSample1 = 0.0;
    prevSample2 = 0.0;
    next(1);
}

void IIR::next(int nSamples) {

    // args
    const float* input = in(0); 
    const float a1 = in0(1);
    const float a2 = in0(2);
    const float b0 = in0(3);
    const float b1 = in0(4);
    const float b2 = in0(5);

    float* output = out(0);

    int i = 0;


    for (; i < nSamples - 4; i += 4) {
        float filtered0 = ((input[i] - a1 * prevSample1 - a2 * prevSample2) * b0) + (b1 * prevSample1) + (b2 * prevSample2);
        output[i] = zapgremlins(filtered0);
        prevSample2 = prevSample1;
        prevSample1 = output[i];

        float filtered1 = ((input[i + 1] - a1 * prevSample1 - a2 * prevSample2) * b0) + (b1 * prevSample1) + (b2 * prevSample2);
        output[i + 1] = zapgremlins(filtered1);
        prevSample2 = prevSample1;
        prevSample1 = output[i + 1];

        float filtered2 = ((input[i + 2] - a1 * prevSample1 - a2 * prevSample2) * b0) + (b1 * prevSample1) + (b2 * prevSample2);
        output[i + 2] = zapgremlins(filtered2);
        prevSample2 = prevSample1;
        prevSample1 = output[i + 2];

        float filtered3 = ((input[i + 3] - a1 * prevSample1 - a2 * prevSample2) * b0) + (b1 * prevSample1) + (b2 * prevSample2);
        output[i + 3] = zapgremlins(filtered3);
        prevSample2 = prevSample1;
        prevSample1 = output[i + 3];
    }

    // Process remaining samples
    for (; i < nSamples; ++i) {
        float filtered = ((input[i] - a1 * prevSample1 - a2 * prevSample2) * b0) + (b1 * prevSample1) + (b2 * prevSample2);
        output[i] = zapgremlins(filtered);
        prevSample2 = prevSample1;
        prevSample1 = output[i];
    }
}

} // namespace IIR

PluginLoad(IIRUGens) {
    // Plugin magic
    ft = inTable;
    registerUnit<IIR::IIR>(ft, "IIR", false);
}

1 Like

A new implementation, especially if more efficient or precise (since higher-order filters depend on precision), would be very welcome!

It sounded to me like you were saying that LTI is a “different type of filter” from an IIR, hence not an IIR. Probably some confusion over terminology,

In any case, you’re now saying there would be benefit from a different (TDF-II) implementation of higher-order filtering. No dispute there! I was never saying “don’t write a new one” – I was only saying that the LTI UGen as it now exists can do IIR filtering, as it’s a general-case recursive filter. That’s all. I don’t see any room for disagreement about that. When you were saying that LTI as it is today is not an optimal implementation, that is a/ true and b/ a rather different statement from what it seemed like you were saying.

hjh

2 Likes

No Worries. I got confused because LTI is a broader term (often emphasizing one aspect over another, although they are all the same kind), but after seeing the equation and the code, I realized one could do the same things with either of them. I should have checked it out first.

I’m not quite sure what this means sorry. I’m pretty new to c++ and ugen programming…

10 posts were split to a new topic: Implementing various filters, using scipy