Digital filtering « LeafLabs Garden

LeafLabs Garden » Support » Maple Hardware Support

Digital filtering

(10 posts) (3 voices)

Started 2 years ago by pyrohaz
Latest reply from mlundinse

pyrohaz
Member

Hey again guys (sorry for the mass questions!)

I'm wondering how I can make a variable low pass filter in code? It needs to be a chebyshev lowpass with a passband ripple of minimum 6dB (to give that resonant peak), cutoff frequency needs to be variable from about 100 up to about 5k or so. Sample rate is 28014Hz and my wavetable consists of values between 0 and 1024.

I had a go with the website: http://www-users.cs.york.ac.uk/~fisher/mkfilter/trad.html in creating a couple of filters and these work in code but I can't figure out how to vary them! I know its dependent on the output gain, and the two values that preceed yv[0] and yv[1] but I can't seem to find any correlation in how these affect the cutoff frequency.

As you can probably tell, i'm not very versed in filter knowledge!

Cheers,

Posted 2 years ago #
pyrohaz
Member

Hmm after searcing the MusicDSP database, I implemented the 1-RC and C filter. So far it seems to be working well though its being used in a realtime situation and unfortunately seems to be running quite slow. Since its dependent on float math, how can I convert it to use integer math instead?

filterval = 1+analogRead(filterpot) //Between 1 and 4096
c = 2*3.14*filterval/samplerate //Small angles, tanx = x, samplerate is 28014hz
r = pow(0.5, resonance+24/16)) //Resonance value lies in the range of 1 - 32ish
v0 = (1-r*c)*v0 - (c)*v1 + (c)*inputval;
v1 = (1-r*c)*v1 + (c)*v0;

Output is v1, written to pwm, inputval is a sample in the range of 1 and 1024.

Could anyone aid me in converting this to integer math?

Thanks, Harris

Posted 2 years ago #
mlundinse
Member

The slow part here is probably the calculation of the filter coefficients. So if you read the filterpot at a lower frequency (10Hz ?) and reclaculate r and c only when the pot values or the resonance have changed, then the main filter loop should run real fast.

Use 16.16 or 8.24 fixed point arithmetic. That means the low order 16 or 8 bits of a 32 bit integer are fraction, and the high bits are integer part.
- Add numbers as usual
- Multiply numbers and then normalize by shifting 16 or 8 bits to the right.
- Flot to pixed point, just multiply with the power of 2, 65536 or 256.

Regards, Magnus

Posted 2 years ago #
pyrohaz
Member

Hey, cheers for the reply! Well the v0 and v1 section is done within an interrupt to write the next sample to the PWM register.

Am I right in assuming:
filterval = (100+analogRead(filterpot));
c = pow(2, 16)*(6.28*(100+filterval)/SampleRate);
r = pow(2, 16)*pow(0.5, ((filterval/128)+30)/16);

v0 = ((65535-r*c)*v0 - (c)*v1 + (c)*arrayu[IncrementV1])>>16;
v1 = ((65536-r*c)*v1 + (c)*v0)>>16;

Or do I need to shift every multiplication by 16?

By multiplying the c by 2^16, gives a range from 63114 to 1483 (within range of unsigned int) and multiplying r by 2^16 gives a range of 129818 to 89750 (within the range of a long).

Am I correct in assuming that (1-r*c) needs to become (2^16-r*c)?

The only reason I really find this hard to grasp is for example
0.1*0.3 in integer math would be

((0.1*2^16)*(0.3*2^16))>>32

Where the floats (0.1 and 0.3) are being multiplied to make them large integers then divided later again to turn them into floats if you get me.

Thanks, Harris

Posted 2 years ago #
gbulmer
Moderator

pyrohaz - are you comfortable with the command-line tools (used under the covers by Maple IDE?)

I ask because the ARM Launchpad embedded gcc could either replace the parts in the IDE, or used from the command line.

The important point about the Launchpad gcc is it's version 4.7, which supports fixed point, and saturating data types using integer operations http://gcc.gnu.org/onlinedocs/gcc/Fixed_002dPoint.html

There is a thread on the forum where someone has used it.

I think this is non-trivial, and would take time and effort to get comfortable. However if you're intending to do a lot of coding, it might make the arithmetic expressions easy.

Posted 2 years ago #
mlundinse
Member

as long as there is only two factors in the product you can sum them before shifting.
r and c are float, dont multiply them with 2^16 until you convert them to fixed point integers.

fixed point versions of the filter coefficients:

int oneminusrc = 1<<16-(int)(r*c*65536);
int cfp = (int)(c*65536)

v0 = (oneminusrc*v0 - cfp*v1 + cfp*arrayu[IncrementV1])>>16;

Your example, in fixed point arithmetic the real, floating point, value is multiplied by 2^16. When you multily two such numbers the real walue will be scaled by 2^32, but then we shift by 16 bits to get the scaling back to 2^16, as a fixed point integer used in your filter loop. Then afterwards the fixed point value value is divided by floating point 65536.0 to get the real float value back.

((0.1*2^16)*(0.3*2^16))>>16 == 0.03*2^16

Posted 2 years ago #
pyrohaz
Member

gbulmer - unfortunately, I know absolutely nothing about the command line tools! Would they be hard to implement and would they increase the speed of the filter command i'm trying to use?

mlundinse - cheers for the help! I tried the code you gave me and it didn't seem to work properly, the output just seemed to produce white noise. Since the code is dependent on both v0 and v1, do I need to do any bitshifts for v1?

v0 = (omrc*v0 - (cfp)*v1 + (c)*arrayu[IncrementV2]) >> 16;
v1 = (omrc*v1 + (cfp)*v0) >> 16;

Would the above work? Sorry for my lack of ingenuity, this is the final step of my project and it would be great if it worked haha.

I just realised that as you say scaling by 2^32 because you can rearrange it to (0.1*0.3) * (2^16*2^16) which using simple indices maths turns out as (0.1*0.3*2^32).

Cheers so far!

Posted 2 years ago #
mlundinse
Member

yes I was a bit lazy, as you point out, the v1 must also be bitshifted after the multiplications.

The values in the arrayu must also be scaled as 16.16 bit fixed point values (is this a sine table ?), and the output v1 scaled to fit the number of bits in the pwm.

Posted 2 years ago #
pyrohaz
Member

No problem, thank you for the help so far anyway!

The values in arrayu are just the values for a square wave between 1 and 1024:
int arrayu[512] = {1024, 1024, 1023, 1024, 1024, 1023, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1023, 1024, 1024, 1023, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1024, 1023, 1024, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1023, 1024, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1023, 1024, 1023, 1023, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1023, 1023, 1024, 1024, 1024, 1023, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1023, 1024, 1023, 1024, 1024, 1023, 1023, 1024, 1024, 1024, 1024, 1024, 1023, 1024, 1024, 1023, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1023, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1023, 1024, 1023, 1023, 1024, 1024, 1024, 1024, 1024, 1023, 1024, 1023, 1023, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1023, 1023, 1024, 1024, 1023, 1023, 1024, 1024, 1024, 1023, 1024, 1024, 1023, 1023, 1023, 1023, 1024, 1024, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1023, 1023, 1023, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1023, 1024, 1024, 1024, 1023, 1024, 1024, 1024, 1023, 1024, 1023, 1024, 1023, 1024, 1024, 1024, 1024, 1023, 1023, 1024, 1024, 1023, 1024, 1024, 1023, 1024, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1023, 1023, 1024, 1024, 1024, 1024, 1023, 1024, 1023, 1024, 1024, 1024, 1024, 1024, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1};

I created them on matlab, hence why a few are 1023 instead but thats not a problem.

So sorry to be a pain but would you be able to write me the final piece of code that would be required, all relating v0, v1, r and c?

Thanks a bunch!

Posted 2 years ago #
mlundinse
Member

Check the pow call, perhaps it should be pow( resonance + 24/16, 0.5 ) to calculate a square root.
Check the values of r and c with a debug println so the filter is stable
With array values between 0 and 1024 this corresponds to a very small amplitude, the Maple pwm resolution is 16 bits so scale it up by a factor of 32.
As long as the filter is stable the calculated values of the filtered signal, v1, are in the same range as the inputs, so the filtered values can be directly fed to the pwm.

int v0,v1;
int ci, omrc;

float c = 2*3.14*filterval/samplerate; //Small angles, tanx = x, samplerate is 28014hz
float r = pow(resonance+24/16, 0.5); //Resonance value lies in the range of 1 - 32ish

/* convert filter coefficients to 16.16 fixed points
/* Both int and long are 32 bits on arm architecture */
ci = (int)(c*65536);
omrc = 65536-(int)(r*c*65536);

/* In sample rate loop, Update the phase value, IncrementV2, then read wave table and update filter */

v0 = (omrc*v0 - ci*v1 + ci*32*arrayu[IncrementV2]) >> 16;
v1 = (omrc*v1 + ci*v0) >> 16;

/* Write output */
pwmWrite(BOARD_ANALOG_OUT, v1);

Good luck, I have not tested it or even tried to compile it, but hope it helps.
As you have Matlab you can use Matlab to debug the sanity of the algorithms used.
Magnus

Posted 2 years ago #

RSS feed for this topic

Reply

You must log in to post.