Reputation: 278
I am trying to do a moving average filter in C language, I've adapted a matlab program that works correctly, the input of my filter is a .pcm archive (a sweep audio signal), the problem to me is the output archive of the moving average in C language, the output goes wrong, the signal only decreases along the time (don't filtter).
Below My C code:
#include <stdio.h>
#include <stdlib.h>
#define m 16 // MEDIA LENGTH
int main()
{
short *x; // X VECTOR = 0
double coef[m]; // COEF VECTOR = 0
int i; // FOR COUNT VARIABLE
int n; // FOR COUNT VARIABLE
int itera; // LENGTH OF THE INPUT FILE
double aux; // AUX VARIABLE TO MAKE THE SUM
double *y; // AUX VECTOR TO RECEIVE X VECTOR
double *saida; // OUTPUT VECTOR
FILE *arq; // POINTER TO THE FILE
for (i = 0; i < m; i++) {
coef[i] = 1.0/m;
}
arq = fopen("sweep_100_3400.pcm", "rb");
if (arq == NULL)
{
printf("ERROR WHILE OPENING THE FILE\n");
return;
}
// compute size of file
fseek(arq,0,SEEK_END);
itera = ftell(arq)/sizeof(short);
rewind(arq);
// alloc mem for x, read the vector from input file
x = malloc(itera*sizeof(short));
fread(x,sizeof(short),itera,arq);
fclose(arq);
// alloc mem for y
y = malloc(itera*sizeof(double));
saida = malloc(itera*sizeof(double));
for (i=0; i < itera; i++) {
y[0] = x[i];
aux=0;
for (n=0; n<m; n++){
aux += coef[n] * y[n];
}
saida[i]=aux;
for (n=m; n <2; n--){
x[n] = x[n-1];
}
}
arq=fopen("saida_medial_movel_c.pcm","wb");
fwrite(saida, sizeof(double),itera,arq);
fclose(arq);
free(saida);
free(y);
free(x);
}
The image below is the output of the matlab program to moving average with length 16:
This image is the output in C language with moving average with length 16:
Someone knows what may be?
Below the code in matlab, that I've adapted:
%MOVING AVERAGE EXAMPLE SCRIPT
clear all;
close all;
clc;
% DEFINES MEDIA LENGTH
m = 16;
%VECTOS EQUAL ZERO
x = zeros (m,1);
coef = zeros (m,1);
%INITIALIZE VECTOR
for j = 1 : m,
coef (j,1) = 1/m;
end
%READ INPUT FILE
fid = fopen ('sweep_100_3400.pcm','rb');
s = fread (fid, 'int16');
fclose(fid);
subplot(2,1,1);
plot(s);
grid on;
title('ENTRADA DO FILTRO');
%PROCESS
itera = length(s);
sav_y = zeros (itera, 1);
%EXECUTE PROCESS
for j = 1 : itera,
x(1,1) = s (j,1);
%PRODUCTS SUM
y=0;
for n = 1 : m,
y = y + coef(n,1) * x(n,1);
end
sav_y(j,1) = y;
%SHIFT THE VECTOR
for n = m: -1 : 2,
x (n,1) = x(n-1,1);
end
end
%PLOT OUTPUT
subplot (2,1,2);
plot (sav_y);
grid on;
title('SAÍDA DO FLITRO');
%SAVE THE OUTPUT IN ANOTHER FILE
fid = fopen('saida_mm_manual.pcm','wb');
fwrite(fid,sav_y,'int16');
fclose(fid);
Update 1 (Using the answer above):
The begining of the signal still with interference, but from the middle to the end of the signal is correctly.
Upvotes: 2
Views: 7513
Reputation: 138
There are a few issues with your main signal processing loop.
for (i=0; i < itera; i++) {
y[0] = x[i];
aux=0;
for (n=0; n<m; n++){
aux += coef[n] * y[n];
}
saida[i]=aux;
for (n=m; n <2; n--){
x[n] = x[n-1];
}
}
Every element of your coef array is the same, so you may as well make it a single constant value of 1.0 / m. For the y array, you are not setting any of its elements apart from the first. Thus your loop adds a constant multiplied by an uninitialized value 16 times over to aux. This is what produces the garbage output. I'm not sure what your (n=m; n <2; n--) loop is supposed to do, but it will never even run since m > 2.
A (inefficient) simple moving average would look something more like this:
for (i = 0; i < itera - m; i++) {
aux = 0;
for (n = 0; n < m; n++){
aux += x[i+n];
}
saida[i] = aux * (1.0 / m);
}
A more efficient version could avoid repeatedly processing intermediate elements, and instead just add the new number entering the sliding window and subtract the old number leaving the window each iteration. If dealing with floating point numbers, care would have to be taken with numerical stability when handling ill conditioned numbers etc, but that's a whole different problem that you don't have to concern yourself about for now.
Upvotes: 2