Search All of the Math Forum:

Views expressed in these public forums are not endorsed by NCTM or The Math Forum.

Notice: We are no longer accepting new posts, but the forums will continue to be readable.

Topic: filter.m to C - algorithm problem
Replies: 4   Last Post: Feb 11, 2015 4:41 PM

 Messages: [ Previous | Next ]
 yllerena Posts: 3 From: Portugal Registered: 6/4/10
filter.m to C - algorithm problem
Posted: Jun 7, 2010 3:44 AM

I'm trying to create my own C implementation of matlab filter.m funtion to down-sample a signal around 2 kHz to speed up autocorrelation process. But I have some inconsistences on results. That's way I need some help on this.

To test the veracity of my implementation and under assuntion of:
y(n) = b(1)*x(n) + b(2)*x(n-1) + ... + b(nb+1)*x(n-nb)
- a(2)*y(n-1) - ... - a(na+1)*y(n-na)

In filte.c (attached) file are the test set and result set and also the matlab result set for comparition)
==================================
#include <stdio.h>
#include <stdlib.h>

void filter(double *b, int ord, double a, double *x, int np, double *y) {
y[0]=b[0]*x[0];
int n,k;
for ( n=0; n<ord; n++) {

for (k=0; k<n+1; k++)
y[n] += b[k]*x[n-k];

for (k=0; k<n; k++)
y[n] -= a*y[n-k-1];
}

for (n=ord; n<np; n++) {

for (k=0; k<ord; k++)
y[n] += b[k]*x[n-k];

for (k=0; k<ord; k++)
y[n] -= a*y[n-k-1];
}
}

void main() {
double b[24] = {0.000000000000000, 0.000000333087156, 0.000001496873275, 0.000003958627788, 0.000008120976691, 0.000014289216468, 0.000022641741288, 0.000033205011502,
0.000045834367448, 0.000060201817729, 0.000075791709912, 0.000091904930112, 0.000107671983786, 0.000122074992086, 0.000133978306167, 0.000142167106464,
0.000145393026160, 0.000142425528942, 0.000132107491488, 0.000113413201241, 0.000085506789149, 0.000047798983375, -0.000000000000000,-0.000057833614664};
double x[72] = {0.018035888671875, -0.019561767578125, -0.021911621093750, -0.022125244140625, -0.022583007812500, -0.023925781250000, -0.024719238281250, -0.025207519531250,
-0.026062011718750, -0.027343750000000, -0.028350830078125, -0.030029296875000, -0.030975341796875, -0.031372070312500, -0.032012939453125, -0.032806396484375,
-0.033569335937500, -0.035156250000000, -0.035186767578125, -0.035461425781250, -0.036071777343750, -0.037292480468750, -0.037322998046875, -0.038085937500000,
-0.039215087890625, -0.039245605468750, -0.039672851562500, -0.039916992187500, -0.040374755859375, -0.039184570312500, -0.038635253906250, -0.039123535156250,
-0.040679931640625, -0.040679931640625, -0.039825439453125, -0.039672851562500, -0.039276123046875, -0.039093017578125, -0.038635253906250, -0.038543701171875,
-0.038879394531250, -0.039001464843750, -0.038238525390625, -0.037597656250000, -0.036773681640625, -0.036376953125000, -0.035766601562500, -0.034881591796875,
-0.034667968750000, -0.034088134765625, -0.033111572265625, -0.032562255859375, -0.031585693359375, -0.030548095703125, -0.029327392578125, -0.029144287109375,
-0.028320312500000, -0.027709960937500, -0.026672363281250, -0.026062011718750, -0.025482177734375, -0.025115966796875, -0.024230957031250, -0.023498535156250,
-0.023101806640625, -0.022247314453125, -0.021301269531250, -0.020904541015625, -0.020080566406250, -0.018554687500000, -0.018188476562500, -0.017486572265625};

double *y;
y = (double *) malloc(72 * sizeof(double));

filter(b,24,1,x,72,y);
int i;
for (i=0; i<72; i++) {
printf("%.15f ",y[i]);
}

}

/* Result
0.000000000000000 0.000000006007523 0.000000014474143 0.000000014335737 -0.000000005954683 -0.000000057385410 -0.000000149938358 -0.000000291793276
-0.000000488407714 -0.000000741937453 -0.000001050908125 -0.000001409802937 -0.000001809000273 -0.000002234576123 -0.000002668108635 -0.000003086976954
-0.000003464981861 -0.000003773101643 -0.000003980694477 -0.000004056219216 -0.000003968281081 -0.000003687129048 -0.000003186343097 -0.000002443992560
0.000000815087300 -0.000001254812358 -0.000001366742322 -0.000001307769569 -0.000001235503455 -0.000001253497965 -0.000001256025278 -0.000001244436353
-0.000001286030533 -0.000001418951286 -0.000001591368623 -0.000001904559844 -0.000002231430220 -0.000002538024524 -0.000002870260838 -0.000003215156027
-0.000003540022139 -0.000003926858378 -0.000004115407794 -0.000004192924874 -0.000004156223620 -0.000004035167979 -0.000003625825633 -0.000003061956385
-0.000002331996243 0.000000987735378 -0.000001061928367 -0.000001152228361 -0.000001110064527 -0.000000878691016 -0.000000726927706 -0.000000638802812
-0.000000696999764 -0.000000729513200 -0.000000760205157 -0.000000869949092 -0.000001091533397 -0.000001339949343 -0.000001538558821 -0.000001797326394
-0.000002141919099 -0.000002498693914 -0.000002859401337 -0.000003030116043 -0.000003069301805 -0.000003037576816 -0.000002916745842 -0.000002482098645
*/

/*
Matlab results
0.000000000000000, 0.000038728700701, 0.000212468158825, 0.000669952090429, 0.001603869555242, 0.003239789754726, 0.005821697041210, 0.009595712388848,
0.014788450415421, 0.021584557421138, 0.030107784610860, 0.040399973744600, 0.052400321296955, 0.065929790887386, 0.080678650284441, 0.096198113468023,
0.111898106348139, 0.127052412074315, 0.140811853415091, 0.152223263872654, 0.160255890238859, 0.163831872700465, 0.161860570084496, 0.153278230462050,
0.151654132385857, 0.150283032845994, 0.148819248112176, 0.147502321307337, 0.146569989138610, 0.146029414393046, 0.146376158011450, 0.146483290451214,
0.146723609633208, 0.147979917324537, 0.149170865285480, 0.150130608131086, 0.151282454133807, 0.152053641843408, 0.152245396299540, 0.151768031227638,
0.150816292130522, 0.149814912956444, 0.148571612894952, 0.147797061447995, 0.146939840157426, 0.145917255327909, 0.144987225361455, 0.143622055161194,
0.141875587930011, 0.140309739204180, 0.138772362906510, 0.137599749478392, 0.136113778987935, 0.135812276863310, 0.136298453733806, 0.136624098564743,
0.136895645597126, 0.137684533116548, 0.138179364639825, 0.139324399180288, 0.140602894249594, 0.141566998378014, 0.141782363390037, 0.141445743997878,
0.140566082677064, 0.139567963532260, 0.138825406303969, 0.137300249608882, 0.135515699671258, 0.134191265772710, 0.132725082619978, 0.131163927795268
*/
=====================================
Best regards

Date Subject Author
6/7/10 yllerena
4/10/13 Thang Ngo