matlab、C语言实现时域卷积运算

背景

某次面试（岗位为音频算法）遇到了c语言实现卷积的编程题。当时不够精通c语言，写的程序比较垃圾。现在重新整理了一下。

原理

卷积公式：

         

matlab有自带的计算卷积的函数conv()

根据公式，现编写实现卷积的函数并与自带函数做对比。

关键点

假设x的长度为a，y的长度为b；在求和时要注意满足条件：m>=0&&m<a,n-m>=0&n-m<b.

具体实现

clc;
clear;
close all;

x=[1,2,3,4,5,6,7,8,9];
y=[11，12，13];

z1=conv(x,y);
z2=myconv(x,y);

z1=roundn(z1,-11);
z2=roundn(z2,-11);

disp(z1);
disp(z2);

if(all(z1(:)== z2(: )))
    fprintf('true\n');
else
    fprintf('error\n');
    for i=1:1:length(z1)
        if(z1(i)~=z2(i))
           fprintf('z1(%d)=%f\n',i,z1(i)); 
           fprintf('z2(%d)=%f\n',i,z2(i)); 
        end
    end
end

function output=myconv(x,y)
%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%x,y------------输入数组
%%%%output---------输出数组

a=length(x);
b=length(y);

if(a==0||b==0)
    fprintf('error:输入数组长度为0!\n');
    return;
end

z=zeros(1,a+b-1);

for i=1:1:a+b-1
    z(i)=0;
    for j=max(i+1-b,1):1:min(i,a)
        z(i)=z(i)+x(j)*y(i-j+1);
    end
end

output=z;
end

for irr=1:1:10
len_x=ceil(100*rand)+1;
x=zeros(1,len_x);
for i=1:1:len_x
    x(i)=100*rand;
end
len_y=ceil(100*rand)+1;
y=zeros(1,len_y);
for i=1:1:len_y    
    y(i)=100*rand;
end

z1=conv(x,y);
z2=myconv(x,y);

z1=roundn(z1,-11);
z2=roundn(z2,-11);

if(all(z1(:)== z2(: )))
    fprintf('true\n');
else
    fprintf('error\n');
    for i=1:1:length(z1)
        if(z1(i)~=z2(i))
           fprintf('z1(%d)=%f\n',i,z1(i)); 
           fprintf('z2(%d)=%f\n',i,z2(i)); 
        end
    end
end
end

经过多次随机验证，效果非常好。

运算结果：

 

 

c语言实现版：

#include<stdio.h>
#include<malloc.h>

int main()
{
	float x[] = { 1,2,3,4,5,6,7,8,9};
	float y[] = { 11，12，13 };
	int a = sizeof(x) / sizeof(x[0]);
	int b = sizeof(y) / sizeof(y[0]);
	int length =a + b - 1;
	float* z = (float*)malloc(length * sizeof(float));
	float* myconv(float x[], float y[],float z[],int a,int b);
	float* p;
	p = (float*)myconv(x,y,z,a,b);
	for (int i = 0; i <length; i++)
	{
		printf("%f\n", *(p+i));
	}

}

int max(int m, int n)
{
	return (m >= n) ? m : n;
}

int min(int m, int n)
{
	return (m <= n) ? m : n;
}

float* myconv(float x[], float y[],float z[],int a,int b)
{
	if (a == 0 || b == 0)
	{
		return 0;
	}
	int len = a + b - 1;
	for (int i = 0; i < len; i++)
	{
		z[i] = 0;
		for (int j = max(i-b+1,0); j <= min(i, a-1);j++)
		{
			z[i] =z[i] + x[j] * y[i - j];
		}
	}
	return z;
}

 结论

 matlab: roundn()控制精度，即小数点后多少位。matlab的库函数conv计算结果的保留位数与自定义函数不完全相同，所以统一使用roundn()保留相同位数，便于比较。

c语言：实现起来非常简单，注意内存分配和函数声明即可。