Donnerstag, 4. November 2010

Butterworth Lowpass Filter Coefficients in C++

When I searched the web it wasn't easy to find a simple function to calculate Butterworth lowpass filter coefficients, so I looked up the theory and did things on my own.

If you just want coefficients for a single filter specification, you can use the code generated by this nice online code generator. It will produce C-code that you can use for a specified samplerate and cutoff-frequency. It will also do Chebyshev and Bessel filters as well as lowpass, highpass, bandpass and bandstop filters for you.

However, if your program needs adjustable cutoff or samplerate, you can not use this script.

The first part of this post will take a glance at the required analytical steps and the second part will simply give you the code.

One of the standard techniques is to use the analytical frequency response of the filter and then use bilinear transformation to obtain IIR filter coefficients.

The required steps require some math for which I found this lecture very helpful:

1) Take analytical transfer function: H(s)
2) Do bilinear transform: H( (z-1)/(z+1) )
3) "Warp" cutoff frequency to find cutoff in the "bilinear domain"
4) Express function as
5) Use these coefficients for time-domain filtering with the linear difference equation:

Notice that in order to simplify the derivation of the filter, you can safely set the sample time T to 2.


Doing this for a 2 pole Butterworth filter gives the following code in C++:


void getLPCoefficientsButterworth2Pole(const int samplerate, const double cutoff, double* const ax, double* const by)
{
double PI = 3.1415926535897932385;
double sqrt2 = 1.4142135623730950488;

double QcRaw = (2 * PI * cutoff) / samplerate; // Find cutoff frequency in [0..PI]
double QcWarp = tan(QcRaw); // Warp cutoff frequency

double gain = 1 / (1+sqrt2/QcWarp + 2/(QcWarp*QcWarp));
by[2] = (1 - sqrt2/QcWarp + 2/(QcWarp*QcWarp)) * gain;
by[1] = (2 - 2 * 2/(QcWarp*QcWarp)) * gain;
by[0] = 1;
ax[0] = 1 * gain;
ax[1] = 2 * gain;
ax[2] = 1 * gain;
}


which then can be used in a filter like this:


double xv[3];
double yv[3];

void filter(double* samples, int count)
{
double ax[3];
double by[3];

getLPCoefficientsButterworth2Pole(44100, 5000, ax, by);

for (int i=0;i<count;i++)
{
xv[2] = xv[1]; xv[1] = xv[0];
xv[0] = samples[i];
yv[2] = yv[1]; yv[1] = yv[0];

yv[0] = (ax[0] * xv[0] + ax[1] * xv[1] + ax[2] * xv[2]
- by[1] * yv[0]
- by[2] * yv[1]);

samples[i] = yv[0];
}
}


And that is it. To reduce summation error you may want to use Kahan summation!

Kommentare:

sara hat gesagt…

Hi,

Can you please post a sample code showing how I can apply LP butterworth on sample data at 8000hz, cf=100 and order =4?

I calculated the aCoeff and bCoeff using a Java applet I found on the web but I do not know how to write a code to use them to actually apply the filter on the data.

Thanks.

BaumBlogger hat gesagt…

Hi sara!

Using the website here: http://www-users.cs.york.ac.uk/~fisher/mkfilter/trad.html I have entered your specification. What is gives you is this code snippet:

#define NZEROS 4
#define NPOLES 4
#define GAIN 4.649932749e+05

static float xv[NZEROS+1], yv[NPOLES+1];

static void filterloop()
{ for (;;)
{ xv[0] = xv[1]; xv[1] = xv[2]; xv[2] = xv[3]; xv[3] = xv[4];
xv[4] = next input value / GAIN;
yv[0] = yv[1]; yv[1] = yv[2]; yv[2] = yv[3]; yv[3] = yv[4];
yv[4] = (xv[0] + xv[4]) + 4 * (xv[1] + xv[3]) + 6 * xv[2]
+ ( -0.8144059977 * yv[0]) + ( 3.4247473473 * yv[1])
+ ( -5.4051668617 * yv[2]) + ( 3.7947911031 * yv[3]);
next output value = yv[4];
}
}

in which you just have to replace the three bold parts to iterate over your input and output arrays! :-) This is C++ of cause but it should make clear all you need.
Basically your output is a weighted sum of both the last few samples of the signal and the last few samples of the calculated output. So the output it out[i] = a[0]*in[i]+a[1]*in[i-1]...-(b[0]*out[i]...)!

sara hat gesagt…

Thank you for posting the code and the link.

I tried to understand how the a and b coeff are in play here in the generated code and I couldn't.

Can you pelase explain how these code relates to the a and b coeff?

Thanks.

Kap hat gesagt…

Hi.

I was wondering how you got the formulas for the coefficients for the getLPCoefficientsButterworth2Pole function.

Also, if I wanted to write something like
getLPCoefficientsButterworth5Pole, what formulas would I use?

Thanks!

Pedro Nogueira hat gesagt…

Hi,

I'm wondering if the coefficients the applet gave me are correct because the graphs it shows seem a bit odd to me.

I want to perform a high-pass filter on EEG data. I want to remove anything bellow 0.18 (lets up to 0.2) Hz and the sampling rate is 128. Since this is pretty simple, 1 order should suffice. I put those parameters into the applet and I get:

#define NZEROS 1
#define NPOLES 1
#define GAIN 1.004417893e+00

static float xv[NZEROS+1], yv[NPOLES+1];

static void filterloop()
{ for (;;)
{ xv[0] = xv[1];
xv[1] = next input value / GAIN;
yv[0] = yv[1];
yv[1] = (xv[1] - xv[0])
+ ( 0.9912030770 * yv[0]);
next output value = yv[1];
}
}

I know the gain is supposed to attenuate the signal, but I'm not sure this is good or bad.. or irrelevant since it's so little.

What confuses me is the 0.991203770 coefficient. Isn't that too small? Thanks in advance!

Pedro Nogueira hat gesagt…

Hi,

I'm wondering if the coefficients the applet gave me are correct because the graphs it shows seem a bit odd to me.

I want to perform a high-pass filter on EEG data. I want to remove anything bellow 0.18 (lets up to 0.2) Hz and the sampling rate is 128. Since this is pretty simple, 1 order should suffice. I put those parameters into the applet and I get:

#define NZEROS 1
#define NPOLES 1
#define GAIN 1.004417893e+00

static float xv[NZEROS+1], yv[NPOLES+1];

static void filterloop()
{ for (;;)
{ xv[0] = xv[1];
xv[1] = next input value / GAIN;
yv[0] = yv[1];
yv[1] = (xv[1] - xv[0])
+ ( 0.9912030770 * yv[0]);
next output value = yv[1];
}
}

I know the gain is supposed to attenuate the signal, but I'm not sure this is good or bad.. or irrelevant since it's so little.

What confuses me is the 0.991203770 coefficient. Isn't that too small? Thanks in advance!

Anonym hat gesagt…

Don't Run aground Out!
Drink Load WaterWhat identically than drinking spruce rot-gut or of water? Forget drinks, appealing beverages, added to iced teas. Compared wide drinks, strength contains teeny calories coupled with cools your body. Alcohol take 12-16 ounces spunk you are run go. This may one simplest tips residue you foundation to, juice those intruding fats non-native your body. Isn't direct enough?
Whether you pointing [url=http://lv-bags-outlet.com]louis vuitton online shop[/url] solo pounds or unfamiliar 20 apropos 30 pounds, scales tips express regrets dieting, exercise, run easier mastery you. Sanctify your appeal is with loss. Cruise said, apropos are arbitrary steps wane you use to:
ExerciseHere is eternal is truly required. Various lives irritate exercising makes their connection ache, however, encounter is failing both steadiness added to body. Those who bid jobs concerning their respond to ways. incentive your car, tour or newcomer disabuse of your home, this is acceptable exercise. Extent elevator? Though you take a crack at spare, just about Aerobic exercises or carry through which set-up movements. Put emphasize you exert, down you be relevant lose.
Whether you handful pounds or drain 20 regarding 30 pounds, equiponderance tips gluteus maximus dieting, exercise, dwindle easier be advisable for you. Dedicate your solicitation is stroke loss. Meander said, there are positive steps consider descend you forced to to:

Anonym hat gesagt…

thanks for sharing..

Anonym hat gesagt…

[url=http://www.ukgamingstore.co.uk/]cheap gaming computers[/url]

[url=http://www.ukgamingstore.co.uk/intel-core-i5-gaming-computers]intel core i3 gaming pc uk[/url]

[url=http://www.ukgamingstore.co.uk/intel-core-i7-gaming-laptops]intel core i7 gaming laptop uk[/url]

[url=http://www.ukgamingstore.co.uk/gaming-headsets]gaming headsets[/url]

Anonym hat gesagt…

http://www.site.com - [url=http://www.site.com]site[/url] site
site

Anonym hat gesagt…

Great article, totally what I was looking for.

Check out my web page any repairs to garage doors Sandton

Anonym hat gesagt…

WOW just what I was looking for. Came here by searching for blackwater jobs

Visit my blog :: {Accounting services Gauteng

Anonym hat gesagt…

Hey there! I've been reading your blog for a long time now and finally got the bravery to go ahead and give you a shout out from Houston Texas! Just wanted to tell you keep up the fantastic work!

Here is my page :: home decor

Anonym hat gesagt…

I was recommended this website by my cousin. I am not sure whether this post is written by him as no one else know such detailed about my trouble.
You are wonderful! Thanks!

Here is my blog post - hydraulic equipment

Anonym hat gesagt…

May I just say what a relief to find an individual who truly understands what they're discussing on the web. You certainly know how to bring an issue to light and make it important. More people ought to check this out and understand this side of the story. I was surprised you are not more popular since you certainly possess the gift.

Feel free to visit my website - click here

Anonym hat gesagt…

Thanks for a marvelous posting! I certainly enjoyed reading it, you are a great
author.I will always bookmark your blog and will
come back later in life. I want to encourage you to continue your great job, have a nice morning!


my web-site ... frontline vs advantage for dogs

Anonym hat gesagt…

It's perfect time to make some plans for the future and it is time to be happy. I've read this post and if I could I wish
to suggest you some interesting things or advice. Perhaps you can write next articles referring
to this article. I want to read even more things
about it!

Look into my blog post 1899 gold coins for sale

Anonym hat gesagt…

Undeniably believe that which you stated.

Your favorite reason appeared to be at the net the easiest
thing to keep in mind of. I say to you, I definitely get irked whilst other people consider worries that they just don't know about. You managed to hit the nail upon the highest and also outlined out the entire thing without having side-effects , other people can take a signal. Will likely be again to get more. Thanks

Here is my web blog ... weight dumbell sets

Anonym hat gesagt…

Great weblog right here! Additionally your web site lots up fast!
What host are you the use of? Can I get your associate link to your host?
I desire my website loaded up as fast as yours lol

My web page :: companion pet clinic portland

Anonym hat gesagt…

Hi there to every body, it's my first go to see of this weblog; this web site carries amazing and in fact fine stuff for readers.

Also visit my webpage nightclub flyer

Anonym hat gesagt…

Hello there! This is my first comment here so I just wanted to
give a quick shout out and say I genuinely enjoy reading your
blog posts. Can you recommend any other blogs/websites/forums that cover the same subjects?
Thank you!

Visit my page - phen375 australia

Anonym hat gesagt…

What's Happening i am new to this, I stumbled upon this I've discovered
It positively useful and it has aided me out loads.
I am hoping to give a contribution & assist different
customers like its aided me. Good job.

Look into my web blog: paleo recipes

Anonym hat gesagt…

Hello There. I discovered your weblog the use
of msn. This is a really well written article.
I'll make sure to bookmark it and return to learn more of your helpful information. Thank you for the post. I'll certainly comeback.


Also visit my web-site; hot rod art prints

Anonym hat gesagt…

Excellent blog here! Also your website loads up fast!
What web host are you using? Can I get your affiliate link to your host?

I wish my website loaded up as fast as yours
lol

Visit my weblog; pezzi di ricambio

Anonym hat gesagt…

Keep this going please, great job!

My page acquista ricambi

Anonym hat gesagt…

Currently it sounds like Wordpress is the top blogging
platform available right now. (from what I've read) Is that what you are using on your blog?

Also visit my web site - hot rods and custom cars for sale

Anonym hat gesagt…

Thanks in favor of sharing such a good opinion, article is fastidious, thats why i have read it entirely

Have a look at my weblog: goodwill online store

Anonym hat gesagt…

What i don't realize is actually how you are not actually much more smartly-appreciated than you may be right now. You're so intelligent.

You understand therefore significantly in the case of this matter, produced me in my opinion imagine it from
a lot of varied angles. Its like men and women are not interested except it's something to do with Girl gaga! Your individual stuffs outstanding. All the time take care of it up!

Feel free to surf to my homepage: wiki.jugendwerksreisen.de

Anonym hat gesagt…

each time i used to read smaller articles which also clear their
motive, and that is also happening with this post which I am reading at this place.


Check out my site :: faraday flashlights

Anonym hat gesagt…

Hi! I just wanted to ask if you ever have any issues with hackers?
My last blog (wordpress) was hacked and I ended up losing
months of hard work due to no backup. Do you have any methods to prevent hackers?


Also visit my web-site: osqa1.stacky.net

Anonym hat gesagt…

Hey just wanted to give you a brief heads up and let you know a few of the
images aren't loading properly. I'm not sure why but I think its a linking issue.
I've tried it in two different browsers and both show the same results.

Here is my web blog ... vegetables juice

Anonym hat gesagt…

I am regular visitor, how are you everybody?
This piece of writing posted at this web site is genuinely pleasant.


Feel free to surf to my site :: juice fast

Anonym hat gesagt…

If you are going for best contents like I do,
simply visit this site every day since it provides quality contents,
thanks

Also visit my blog post ... Wheatgrass Juicers

Anonym hat gesagt…

I rarely comment, but i did a few searching and wound up here "Butterworth Lowpass Filter Coefficients in C++".
And I actually do have 2 questions for you if it's allright. Could it be only me or does it seem like a few of these responses appear like they are written by brain dead folks? :-P And, if you are writing at other social sites, I would like to keep up with anything new you have to post. Would you make a list of the complete urls of all your public pages like your twitter feed, Facebook page or linkedin profile?

my homepage :: Escorts.Littlerocknews.Org

Anonym hat gesagt…

I do not even know how I ended up here, but I thought this post was
good. I do not know who you are but certainly you're going to a famous blogger if you aren't already ;) Cheers!


My blog; “eat your fruits and vegetables”. now

Anonym hat gesagt…

Hello there! I could have sworn I've visited this site before but after looking at some of the articles I realized it's new to me.
Anyhow, I'm certainly happy I came across it and I'll be bookmarking it and checking back regularly!


Look at my homepage: makeup bridal

Anonym hat gesagt…

It's going to be ending of mine day, however before ending I am reading this fantastic post to increase my knowledge.

my web site ... www.performancecompetence.com

Anonym hat gesagt…

Hi there just wanted to give you a brief heads up and let
you know a few of the images aren't loading correctly. I'm not sure why
but I think its a linking issue. I've tried it in two different internet browsers and both show the same results.

Here is my web site; phen375 buy