Hello
I am trying to implement an odd symmetric FIR filter (using C code) but i
am struggling with the algorithm .
Here is the piece of code that I am trying to understand. Can anyone give
me some steps on how the array index are working here. I�am still reading
my C book on arrays but still I cannot get the whole idea.Please if
someone out there can give me some hints I will be grateful to learn and
to save some time too. I hope my question is understandable!!!
Thanks for your time and help!
//Lowpass Filter coefficients
const float coeff[FIR_LENGHT ] =
{
0.0008 0.0054 0.9908 0.0054 -0.0008};
//A/D converted values
const unsigned int input[] =
{
0x0400, 0x0800, 0x0C00, 0x1000,0x2400, 0x2000, 0x1C00, 0x1800};
void main(void)
{
int y,i; /*loop counter*/
for (y = 0; y<3; y++)
{
for (i = 0; < FIR_LENGHT/2; i++)
{
output[y] =output[y]+ coeff[i]*(input [ y + 4-i]+ input [y +i]);
};
output[y] = output[y] + (input [ y + 4 -i]* coeff[i] );
};
}
This message was sent using the Comp.DSP web interface on
www.DSPRelated.com
Odd symmetric FIR filter
Started by ●September 23, 2005
Reply by ●September 23, 20052005-09-23
Umutesi Faith wrote:> Hello > I am trying to implement an odd symmetric FIR filter (using C code) but i > am struggling with the algorithm . > Here is the piece of code that I am trying to understand. Can anyone give > me some steps on how the array index are working here. I�am still reading > my C book on arrays but still I cannot get the whole idea.Please if > someone out there can give me some hints I will be grateful to learn and > to save some time too. I hope my question is understandable!!! > Thanks for your time and help! > > //Lowpass Filter coefficients > const float coeff[FIR_LENGHT ] = > { > 0.0008 0.0054 0.9908 0.0054 -0.0008}; > > //A/D converted values > const unsigned int input[] = > { > 0x0400, 0x0800, 0x0C00, 0x1000,0x2400, 0x2000, 0x1C00, 0x1800}; > > void main(void) > { > int y,i; /*loop counter*/ > > for (y = 0; y<3; y++) > { > for (i = 0; < FIR_LENGHT/2; i++) > { > output[y] =output[y]+ coeff[i]*(input [ y + 4-i]+ input [y +i]); > }; > output[y] = output[y] + (input [ y + 4 -i]* coeff[i] ); > }; > } >The code is executing an FIR filter, which convolves the filter impulse response function with the input. It is trying to execute w_n = sum from {k = 0} to N {x_{n+k} * h_k}, which is your basic convolution where x_n is the n'th input, h_k is the k'th filter tap, and w_n is the n'th output. Since the writer knew that the filter was odd symmetrical, he or she knew that h_k = k_{FIR_LENGHT - 1 - k}, and saved two steps in the for loop. C arrays are zero referenced, which means that if you have an array int bob[] = {1, 2, 3, 4}; then bob[0] is 1, bob[1] is 2, etc. The code is, IMHO, poorly written. The writer goes to the trouble to define the length of the FIR filter and use it in a couple of places in the code, but not others (and they misspell "LENGTH"). So the inner FOR loop is counted with FIR_LENGHT, but the outer FOR loop which must traverse the length of the input array minus FIR_LENGHT just hard-codes a 3; similarly the input[y + 4 - i] is hardcoding FIR_LENGHT - 1 as 4. Finally, the trick with the FOR loops will play well on some compilers and machines, particularly CISC machines, but on a real DSP chip you want to set up a repeated MAC in hardware and let it rip. The TI compiler is reported to be able to recognize when this can be done for some chips (I can attest to older Code Composter versions not being able to do it for the '2812). I have yet to actually use a compiler that does this or that supplies a working library function to do it to my satisfaction, so I always end up writing the ten lines or so of assembly (well, 25 if it's a matrix multiply) myself. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com
Reply by ●September 24, 20052005-09-24
"Umutesi Faith" <ma_nz1@yahoo.com> wrote in message news:raCdndq6_oHc-KneRVn-pQ@giganews.com...> Hello > I am trying to implement an odd symmetric FIR filter (using C code) but i > am struggling with the algorithm . > Here is the piece of code that I am trying to understand. Can anyone give > me some steps on how the array index are working here. I�am still reading > my C book on arrays but still I cannot get the whole idea.Please if > someone out there can give me some hints I will be grateful to learn and > to save some time too. I hope my question is understandable!!! > Thanks for your time and help! > > //Lowpass Filter coefficients > const float coeff[FIR_LENGHT ] = > { > 0.0008 0.0054 0.9908 0.0054 -0.0008};ummmmmm.... this isn't a symmetric filter because the end coefficients have opposite signs. Maybe that's a typo? Fred
Reply by ●September 24, 20052005-09-24
It was a typing mistake the FIR filter was intended to be odd symmetric!! Thank you for your replies, at least i can have a small idea on how this filter was implemented! I'm myself trying to implement a FIR filter,i cannot use assembler as i don't know much about it, but i can try to do it in a simple C code, but so far i have found the following steps which are not obvious to understand : Basic Algorithm for implementing FIR filters -Put the input signal in the delay line.// what do they mean by this ?? -Multiply each sample in the delay line by the corresponding coefficient and accumulate the result. -Shift the delay line by one sample to make room for the next input sample. //this is how i see this algorithm For instance if the length of the filter is 4 and the number of input samples is 16, putting the input in the delay does it mean starting from input[0]?? . For the second step the multiply and accumulate is to sum up the the multiplication of the input sample and fir coefficients! So this shifting step , the filtering start at input[1] and the oldest sample is dropped?? I am not sure if i do understand it in the correct way, why i am so struggling implementing these filters , are there some tricks to follow??? Thanks Umutesi This message was sent using the Comp.DSP web interface on www.DSPRelated.com
Reply by ●September 24, 20052005-09-24
Umutesi Faith wrote:> It was a typing mistake the FIR filter was intended to be odd symmetric!! > Thank you for your replies, at least i can have a small idea on how this > filter was implemented! I'm myself trying to implement a FIR filter,i > cannot use assembler as i don't know much about it, but i can try to do it > in a simple C code, but so far i have found the following steps which are > not obvious to understand : > Basic Algorithm for implementing FIR filters > > -Put the input signal in the delay line.// what do they mean by this ?? > -Multiply each sample in the delay line by the corresponding coefficient > and accumulate the result. > -Shift the delay line by one sample to make room for the next input > sample. > > //this is how i see this algorithm > For instance if the length of the filter is 4 and the number of input > samples is 16, putting the input in the delay does it mean starting from > input[0]?? . > For the second step the multiply and accumulate is to sum up the the > multiplication of the input sample and fir coefficients! > So this shifting step , the filtering start at input[1] and the oldest > sample is dropped?? > I am not sure if i do understand it in the correct way, why i am so > struggling implementing these filters , are there some tricks to > follow??? > Thanks > > Umutesi > > This message was sent using the Comp.DSP web interface on > www.DSPRelated.comYou don't need tricks, but you do need to understand the process. You have most of it, with some missing details. Read Chapter 6 in http://www.dspguide.com/ and more starting with Chapter 14. To be sure you understand the details, follow the code examples. When confused, ask more questions here. You will also find good material at http://www.dspguru.com/ and the free on-line courses at http://www.bores.com/. Jerry -- Engineering is the art of making what you want from things you can get. �����������������������������������������������������������������������
Reply by ●September 24, 20052005-09-24
"Umutesi Faith" <ma_nz1@yahoo.com> wrote in message news:XYSdnQCIAq5JiKjeRVn-uw@giganews.com...> It was a typing mistake the FIR filter was intended to be odd symmetric!!***Actually, I'm not used to the term "odd symmetric" and I'm a little concerned that you could mean: - antisymmetric as in 1 2 3 0 -3 -2 -1 or in 1 2 3 -3 -2 -1 instead of: - odd length *and* symmetric as in 1 2 3 4 3 2 1 That's because the latter is an Even sequence of Odd length (if you define time zero at the center anyway - for discussion purposes) and The former examples are both Odd sequences of Odd or Even length (again if you define time zero at the center).> Thank you for your replies, at least i can have a small idea on how this > filter was implemented! I'm myself trying to implement a FIR filter,i > cannot use assembler as i don't know much about it, but i can try to do it > in a simple C code, but so far i have found the following steps which are > not obvious to understand : > Basic Algorithm for implementing FIR filters > > -Put the input signal in the delay line.// what do they mean by this ??***In the olden days the filters operated in continuous time using "analog" delay lines. The delay lines had taps on them (that is, interim output points) The outputs (which include the input and the final output as well as the interim outputs) were each weighted and then all the weighted results were summed to get the filter output.> -Multiply each sample in the delay line by the corresponding coefficient > and accumulate the result.***Note that a discrete version of a delay line is a memory. In the not quite so olden days we built FIR filters out of bucket-brigade delay lines (continuous values with discrete "buckets") and we built them out of shift-register types of memories (discrete values ... bits ... with discrete memory registers that were coupled together in a serial manner). Both versions were clocked to accomplish the shifts. So, if the shifting clock rate was the same as the sample rate the data would move along the chain as each new sample entered and as the oldest sample disappeared as it was replaced by the next oldest sample. Note that this ties the sample rate and the delay values together .. which makes sense, eh? ***So, the software version of one of these would have a number of memory locations to represent the "buckets" or the locations in the shift register. A brute-force method will move the data in the same fashion like this: b(n)=b(n-1) b(n-1)=b(n-2) . . b(1)=b(0) b(0)=xnew then multiply the locations by the same coefficients each time A more sophisticated method (look up "circular buffer") avoids moving all the data and moves a reference pointer /offset instead: offset=k; 0<=k<=(n-1) i=k to (n-1)+k modulo n b(k)=xnew then multiply the *offset* locations by the coefficients> -Shift the delay line by one sample to make room for the next input > sample. > > //this is how i see this algorithm > For instance if the length of the filter is 4 and the number of input > samples is 16, putting the input in the delay does it mean starting from > input[0]?? .***yes. If the filter is a0 a1 a2 a3 and the data is d0 d1 d2 .... d15 then the multiplies start: y0 = a0*d0 y1 = a0*d1 + a1*d0 y2 = a0*d2 + a1*d1 + a2*d0 y3 = a0*d3 + a1*d2 + a2*d1 + a3*d0 y4 = a0*d4 + a1*d3 + a2*d2 + a3*d1 . . . . . . . . . . y15 = a0*d15 + a1*d14 + a2*d13 + a3*d12 y16 = a1*d15 + a2*d14 + a3*d13 y17= a2*d15 + a3*d14 y18= a3*d15 Note that the filter fills up in N-1, that is, 3 sample times and it "unloads" over 3 samples times. These are startup and ending transients that are of length 1 less than the length of the filter. You might well disregard these outputs.> For the second step the multiply and accumulate is to sum up the the > multiplication of the input sample and fir coefficients! > So this shifting step , the filtering start at input[1] and the oldest > sample is dropped??Yes, as above.> I am not sure if i do understand it in the correct way, why i am so > struggling implementing these filters , are there some tricks to > follow???The address change vs. actually shifting the data is one. If the filter is symmetrical then you might add two samples and multiply once - assuming that makes any difference in your machine. Then there might be things about scaling and affects in roundoff errors. Then you might modify the coefficients to be sums of powers of two to eliminate multiplies altogether in favor of shifts. Again, it depends on your machine. This works well in FPGA implementations, for example. That is .. a coefficient of 0.50222 might be rounded to 0.5 and a simple shift will replace a multiply. Or, 0.627 would be 0.5 + 0.125 which would be two shifts and and add to replace the multiply. Do read the things that Jerry pointed to. Fred
Reply by ●September 25, 20052005-09-25
Thank you all for your help. I did read quite a lot tutorial about these FIR filters, but i guess there have been some missing details for me to understand them properly. With all the help i got from here i hope finally to be able to do my implementation! I will keep posting once i get stuck again!! Umutesi This message was sent using the Comp.DSP web interface on www.DSPRelated.com
Reply by ●October 1, 20052005-10-01
Hello
I have tried to understand this shifting method you discribed below, if
understood it correct ,xnew is the new sample. So if for instance
/*input samples*/
Input[7] ={s0,s1,s2,s3,s4,s5,s6,s7};
this would be like Input[0]=xnew;
How is this "xnew" accessing the data in the Input[7] array?
Thanks!
Umutesi
******************************************************************
Re: Odd symmetric FIR filter - Fred Marshall - 12:57 24-09-05
***So, the software version of one of these would have a number of memory
locations to represent the "buckets" or the locations in the shift
register.
A brute-force method will move the data in the same fashion like this:
b(n)=b(n-1)
b(n-1)=b(n-2)
Reply by ●October 1, 20052005-10-01
"Umutesi Faith" <ma_nz1@yahoo.com> wrote in message news:quKdnVpLQtqq8aPeRVn-vA@giganews.com...> Hello > > I have tried to understand this shifting method you discribed below, if > understood it correct ,xnew is the new sample. So if for instance > /*input samples*/ > Input[7] ={s0,s1,s2,s3,s4,s5,s6,s7}; > > this would be like Input[0]=xnew; > How is this "xnew" accessing the data in the Input[7] array? > Thanks! > > UmutesiI don't understand your notation. You seem to be using Input as both a vector and a scalar. That is, an array of multiple elements and an array of 1 element. Note that I'm not trying to write code for you, so this is all pseudo code gibberish that you have to put properly into the language of your choice. I said there is a very brute force and inefficient method that may be instructive: b(n)=b(n-1) b(n-1)=b(n-2) . . b(1)=b(0) b(0)=xnew where "b()" is the buffer. And, I was assuming this was a streaming implementation - so you get xnew at each time step and don't know what it is beforehand. If you are doing this in something like Matlab and already have all the inputs in memory then it might look like this: inputs=[1 2 3 4 3 4 3 5 6 5 4 3 2 1 2 ...... 8 5 9 3 4 6 9] (and of course the number format can be anything you like - I just showed positive integers because I'm being lazy) Now you have inputs(1) inputs(2), etc. and I will state that inputs(1) is the oldest / first data point / sample. Then convolve the inputs vector with the filter unit sample response. loop k=1 y(k)=inputs(k)*filter(1) + inputs(1+k)*filter(2) .... + inputs(N-1+k)*filter(N) k=k+1 continue until (N-1+k) is the maximum index in the array "inputs" where "filter" is the array of filter coefficients (reversed I believe- that is, filter(N) is the first unit sample output value) and I have left out the startup and ending transient computations altogether for simplicity. In this implementation, the filter slides along the data and the data of interest is whatever lines up with the filter coefficients at each sample time. Since the data is all contained in the vector "inputs", the implementation just changes "k" to move the filter along the data. Fred
Reply by ●October 1, 20052005-10-01
Well, this is the most common form of array initialization in C language. This array has 8 elements running from s0 to s7! My intention was just to make sure if i understood well "the pseudo code" so that i can continue with fixed ideas in my studying about the FIR filter implementation!>> /*input samples*/int Input[7] ={s0,s1,s2,s3,s4,s5,s6,s7};>> >> this would be like Input[0]=xnew; >> How is this "xnew" accessing the data in the Input[7] array? >> Thanks! >> >> Umutesi > >I don't understand your notation. You seem to be using Input as both a >vector and a scalar. That is, an array of multiple elements and an arrayof>1 element. Note that I'm not trying to write code for you, so this isall>pseudo code gibberish that you have to put properly into the language of>your choice. > >I said there is a very brute force and inefficient method that may be >instructive: > >b(n)=b(n-1) >b(n-1)=b(n-2) >. >Fred > > >This message was sent using the Comp.DSP web interface on www.DSPRelated.com
