Login / Register  0 items | $0.00 NewWhat is KVR? Submit News Advertise
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Wed Apr 29, 2009 2:24 pm A Collection of Useful C++ Classes for Signal Processing

PROJECT IS HERE:
https://github.com/vinniefalco/DSPFilters
************************************************

It implements an arbitrary-order lowpass and highpass filter with Butterworth response characteristics. There is no Q factor. I will be adding the rest of the filter types including band pass and band reject, and all variations of responses (chebyshev, inverse chebyshev, elliptic, etc).

It uses templates so you can have your choice of processing floats, doubles, or other real valued type. The base class supports arbitrary placement of poles. The Process() function that does the work supports N-channel interleaved or N-channel data using the 'skip' template parameter. Calculations can be performed using a different underlying type for accuracy.

Given an instance of a filter, a method is provided to determine the complex-valued filter response at a given angular frequency, using a closed-form expression.

EXAMPLE CODE REMOVED

UPDATE: The distribution now includes source code for a JUCE test application that demonstrates use of the filters.

All of this code has been tested and listened to using various audio source material. It compiles under Visual Studio 2008, Visual Studio 2010, and XCode.

Comments, constructive criticisms welcome.
Last edited by thevinn on Thu Mar 08, 2012 1:04 pm, edited 5 times in total.
asomers
KVRian
 
646 posts since 17 Feb, 2006, from California

Postby asomers; Thu Apr 30, 2009 3:26 am

:love:

thanks!
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Thu Apr 30, 2009 3:41 am

You're welcome.

I have the band pass working, currently debugging the band stop. I had to make some changes to the class interfaces to make the derived templates interact correctly. I'm unhappy the the GetSPole function is duplicated now but I will hopefully find a way to fix it. Once I have the band stop working and eliminate the duplicate GetSPole() implementation I will add Chebyshev, Inverse Chebyshev, and Elliptic response characteristics.

Oh yeah and there is still the open problem of determining a suitable angular frequency in the passband of the band stop for sampling the response and computing a scale factor.

The final version will also explicitly support N-channel interleaved data by calculating and storing the coefficients only once for all channels, while maintaining a separate history buffer for each channel. Since the computation of elliptic curve coefficients use an iterative numeric approach, this will be much faster.

Channels will be specified through the use of a template paramter. If your data is planar (i.e. N-channel non-interleaved) you will have to create a separate object per channel. Can't win 'em all!
duncanparsons
KVRAF
 
8366 posts since 11 Apr, 2003, from now on the flat

Postby duncanparsons; Thu Apr 30, 2009 5:47 am

looks a good bit of code there vincent :) thx for sharing!
Image
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Thu Apr 30, 2009 5:55 am

I updated the file without documenting the changes unfortunately...still working on the remaining filters. The latest file has band stop and band pass, and I have also done some tests to determine the limit of poles.

From the header:

Butterworth Low pass stable to 53 poles (27 second order stages)
Butterworth High pass stable to 110 poles (55 second order stages)
Butterworth Band pass stable to 80 poles (80 second order stages)
Butterworth Band stop stable to 109 poles (109 second order stages)

It is no longer necessary to call Prepare(). Additionally, there is support for N-channel interleaved data. For this case set skip=0. If your data is not interleaved you will need to create on object per channel and just set skip=0. If you want to filter just one channel of an N-channel interleaved dataset set skip=N-1 and pass the pointer to the beginning of the channel of interest.

Code: Select all
// Create a 6 pole low pass filter to
// process stereo interleaved data
Dsp::ButterLowPass<6,2> butter;
butter.Setup( frequency/sampleRate );
butter.Process<0>( count, stereoBuffer );


Note that Process() takes the number of frames, where each frame contains N channels of interleaved data.

This handling of channel data is optimized over the previous method, which required a separate object per channel. So for N-channels of data, the coefficients are calculated only ones and reused for each channel. There is a separate history buffer for each channel, however.
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Thu Apr 30, 2009 7:24 am

Added a Biquad class that maps coefficients to Cascade Filter

Provided Biquad formula implementations for the following:

Low Pass
High Pass
Band Pass (constant skirt dB)
Band Pass (constant peak dB)
Band Stop
All Pass
Peak Equalizer
Low Shelf
High Shelf

Formula are from RBJ Audio Cookbook so now you can have your "Q" parameter for resonance (although not with orders greater than 2).

All biquads were implemented and tested using the frequency response formula. Allpass phase response was checked for correctness.

They sound great!
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Thu Apr 30, 2009 2:48 pm

File has been updated.

Added the following filters:

Cheby1LowPass
Cheby1HighPass
Cheby1BandPass
Cheby1BandStop
Cheby2LowPass
Cheby2HighPass
Cheby2BandPass
Cheby2BandStop

I changed the method used to normalize the gain. For most cases there is a "passband hint" where we sample at a certain place in the passband. For the Chebyshev Type I filters, Brent's method is used to find the local maximum.

A bug was fixed in the computation for the number of poles required for band pass and band stop, it was using twice the amount of needed stages.

All of these filters have been tested at their extremes.

If anyone has any advice on how the parameters should be passed I would love to hear it. How is the width of the band pass or band stop region defined? Is it typically in angular frequency? Or normalized frequency? Or does it matter? What about the center frequency for the band pass and band stop?

How can I turn my low/hi/band/stop filters into shelving filters? Would there be utility in an Allpass cascade? What about a peaking Eq? How do I get those from these building blocks?

I've renamed this to

"A Collection of Useful C++ Classes for Digital Signal Processing"

Work remaining includes elliptic calculations, a few more utility classes, and I would love some feedback before I publish this to musicdsp.
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Thu Apr 30, 2009 5:36 pm

Updated the file again.

Some improvements were made to the calculation of the gain adjustment. It is more accurate, and for some cases we sample the passband once instead of using Brent's method (thanks, Robin from www.rs-met.com!).
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Fri May 01, 2009 4:21 am

File is updated again.

The original implementation used Direct Form II for the cascade stages. This is a little bit faster than Direct Form I but as I discovered, causes discontinuities in the output when filter parameters are changed during processing.

Therefore, I provide both implementations.

Use ProcessI() for Direct Form I. Its a little bit slower but does not cause discontinuities (i.e. clicks, blurps, and chirps) when filter parameters are changed mid stream.

If your filter parameters do not change you can use ProcessII().
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Fri May 01, 2009 2:45 pm

Comments? Criticisms? Can anyone tell me if it compiles, or what I need to change, under environments other than Visual Studio 2008? I would like to groom this as much as possible before posting it to musicdsp
asomers
KVRian
 
646 posts since 17 Feb, 2006, from California

Postby asomers; Fri May 01, 2009 6:07 pm

This is the most complete self-contained filter implementation I have ever seen, and I intend to use it. First step will be to create an open source multi-format/platform plugin as a test harness, with scope and zplane views like the dfilter applet. I imagine one of the other juce nerds may get to it first, and if so, great! But I'll try to have it done by the end of the next week.

side note:
Your effort here has been truly inspirational! Thanks!
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Fri May 01, 2009 6:28 pm

asomers wrote:This is the most complete self-contained filter implementation I have ever seen, and I intend to use it.


Wow thats high praise, thanks a lot!

I was starting to get discouraged but you've given me renewed vigor to get those ellipticals in there (they are a bit scary). You will be very pleased with the results of the filter output. Let me know if I should change the parameters I provide other ways of specifying filter parameters.

I still feel that the library is incomplete. There are the following pending issues:

- Shelving, peak, all-pass for Butterworth, Chebyshev, and Elliptic.

The Biquads have these versions and I would like the others to have them
as well. It would also be super awesome if higher order filters could have a "Q" parameter for resonance but I'm not expecting miracles, it would require redistributing the poles and zeroes. But if there's a research paper or code out there...I could incorporate it.

- Denormal testing and fixing

I'd like to know if denormals are a problem. And if so, it would be nice to have a small function that can reproduce the denormal problem. This way I can test the fix under all conditions. I will include the function as a "unit test" object in the header file so anyone can verify its correctness. But I'm a little lost.

- Optimized template specializations for stages=1, channels={1,2}

There are some pretty obvious optimizations I am saving for "the end". I don't want to do them until the code is finalized.

- Optimized template specializations for SSE, other special instructions

- Optimized trigonometric functions for fast parameter changes

- Elliptic curve based filter coefficients

- Some utility functions for manipulating sample buffers

So I will still be adding more!
MackTuesday
KVRist
 
450 posts since 11 Jul, 2004, from Southern California, USA

Postby MackTuesday; Sat May 02, 2009 12:59 am

I probably will use this too if I ever make another plugin. Great work.
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Sat May 02, 2009 6:18 am

I wrote a simple loop to time the performance of the filter. I generate random noise in a buffer and then run the low pass on the data over and over again.

After just 3000 iterations, I get denormals. The debugger shows the samples as -1.#IND00000000

Is running a filter on the same block of data over and over again supposed to do that?
User avatar
thevinn
KVRian
 
775 posts since 30 Nov, 2008

Postby thevinn; Sat May 02, 2009 6:41 am

File is updated:

- Added optimized approximations for trignometric functions
- Optimized Biquad coefficient calculations for speed
- Optimized Biquad special case for setting CascadeFilter with 1 stage
- Added SetupFast() to Biquads for high speed parameter changes using
trigonometric approximations.
- Added utility function templates for manipulating sample buffers

I tried to specialize the CascadeFilter processing function but amazingly, the MSVC compiler is so good that the unspecialized version still comes within 2% of the speed. Further improvements will require CPU-specific extensions.
Next

Moderator: Moderators (Main)

Return to DSP and Plug-in Development

Who is online

Users browsing this forum: CCBot (commoncrawl)