滤波器设计和使用

实验目的

利用C语言集成开发环境设计FIR滤波器软件,并进行测试:

  1. 练习随机数、静态局部变量、滤波器等知识的综合运用;
  2. 培养学生综合运用C程序语言进行程序设计的能力;
  3. 培养结构化程序设计能力;
  4. 分析FIR低通滤波器的特点;
  5. 分析wav文件格式;
  6. 设计fir低通滤波器
  7. 用C或C++编程实现滤波器;
  8. 调式程序,测试滤波器效果

实验内容

上机完成以下编程实验,调试运行程序并完成报告。
1、 完成随机数产生函数:double RandDouble(double Max),即生产0到Max之间的随机数,且返回的随机数带2位小数,如6.28等;
2、 编写一函数,实现FIR滤波器功能:输入xi,输出yi,其中xi的定义如5-1所示,即i大于0时,为输入信号,小于0时为0,i为整数;xi和yi的关系如5-2所示。(函数原型自行确定,可以为分帧处理,也可以单点处理) 在这里插入图片描述
在这里插入图片描述
c0,c1,c2,c3,c4,…,cn的为FIR滤波器系数;
3、 写一主函数,在主函数中调用随机数生成函数,产生序列x,设计不小于5阶的滤波器系数,并求出序列y的值,输出x和y的值。
4、 设Cn为某低通滤波器系数,Cn的值见后(可利用Matlab,自行设计滤波器系数)。改写滤波器函数,实现低通滤波功能。编写主函数,读取一个被噪声污染的语音文件,并进行低通滤波处理,把处理后的文件保存下来,并进行试听,观察处理效果。如某一低通FIR滤波器系数如下:
CN[82] = {
0.0002088566718795,0.0005019666320485, 0.00101659692763, 0.001761517078465,
0.002717496869779, 0.003805742549565, 0.004881532449878, 0.005744013542788,
0.006161971401866, 0.005917629852047, 0.004859494981977, 0.002955831081991,
0.0003356045898205,-0.002697639896313,-0.005687071486917,-0.008084032953745,
-0.009344469124396,-0.009045419365927,-0.006999963923105,-0.003344508688152,
0.001426457367602, 0.006496401543416, 0.01083826348381, 0.01339119282049,
0.01327693896666, 0.0100166598574, 0.003703629355468,-0.004909743327555,
-0.01444030683204, -0.02302828304832, -0.02860345497317, -0.02922565654424,
-0.02344372245441, -0.01060796387519, 0.008925525047616, 0.03374481923033,
0.06150720551817, 0.08922868601273, 0.1136966478108, 0.1319402460404,
0.1416800797754, 0.1416800797754, 0.1319402460404, 0.1136966478108,
0.08922868601273, 0.06150720551817, 0.03374481923033, 0.008925525047616,
-0.01060796387519, -0.02344372245441, -0.02922565654424, -0.02860345497317,
-0.02302828304832, -0.01444030683204,-0.004909743327555, 0.003703629355468,
0.0100166598574, 0.01327693896666, 0.01339119282049, 0.01083826348381,
0.006496401543416, 0.001426457367602,-0.003344508688152,-0.006999963923105,
-0.009045419365927,-0.009344469124396,-0.008084032953745,-0.005687071486917,
-0.002697639896313,0.0003356045898205, 0.002955831081991, 0.004859494981977,
0.005917629852047, 0.006161971401866, 0.005744013542788, 0.004881532449878,
0.003805742549565, 0.002717496869779, 0.001761517078465, 0.00101659692763,
0.0005019666320485,0.0002088566718795};
该低通滤波器的设计基本参数如下:

 * Discrete-Time FIR Filter (real)
 * -------------------------------
 * Filter Structure  : Direct-Form FIR
 * Filter Length     : 82
 * Stable            : Yes
 * Linear Phase      : Yes (Type 2)
 * Sample Frequency  : 16e3
 * Fpass             : 1000
 * Fstop             : 1500
 * Astop             : -80dB
 * Design Method     : Equiripple

对于采样率为16000Hz的语音信号,大于1500Hz的信号被抑制,小于1000Hz的信号被完整保留;如对于一个500Hz和2.5KHz的叠加正弦波信号,滤波处理后几乎仅有低频信号;滤波前后的音频文件见附件。
提示:
wav文件格式的为文件头+数据,其中文件头的定义如下:

typedef struct _tagMsWavPcmHeader44{
	unsigned char ChunkID[4]; // "RIFF"; The "RIFF" the mainchunk;
	unsigned long ChunkSize; // FileSize - 8; The size following this data
	unsigned char Format[4]; 
// "WAVE"; The "WAVE" format consists of two subchunks: "fmt " and "data"
	unsigned char SubChunk1ID[4]; // "fmt "
	unsigned long SubChunk1Size; 
// 16 for PCM. This is the size of the rest of the subchunk which follows this data.
	unsigned short AudioFormat; // 1 for PCM. Linear quantization
	unsigned short NumChannels; // 1->Mono, 2->stereo, etc..
	unsigned long SampleRate; // 8000, 11025, 16000, 44100, 48000, etc..
	unsigned long ByteRate; // = SampleRate * NumChannels * BitsPerSample/8
	unsigned short BlockAlign; // = NumChannels * BitsPerSample / 8
	unsigned short BitsPerSample; // 8->8bits, 16->16bits, etc..
	unsigned char SubChunk2ID[4]; // "data"
	unsigned long SubChun2Size; 
// = NumSamples * NumChannels * BitsPerSample / 8. The size of data
} wav_pcm_header44;

实际测试结果

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
以上为滤波前后wav文件的波形利用matlab辅助工具画出滤波前wav文件的频谱
利用matlab辅助工具画出滤波前wav文件的波形
以上为使用matlab画出的滤波前的频谱和波形图
利用matlab辅助工具画出滤波后wav文件的频谱
用matlab辅助工具画出滤波后wav文件的波形
以上为使用matlab画出的滤波后的频谱和波形图

代码结构

在这里插入图片描述

附件(具体代码)

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>
#include<graphics.h>
#pragma warning(disable:4996)
#define FRAME_READ_LENGTH  16
#define FIR_FILTER_LENGTH 82
double max = 32.25;
double c[5]{ 0.01,0.1,0.78,0.1,0.01 };//filter2函数5阶滤波的系数
double CN[82] = {
  0.0002088566718795,0.0005019666320485,  0.00101659692763, 0.001761517078465,
   0.002717496869779, 0.003805742549565, 0.004881532449878, 0.005744013542788,
   0.006161971401866, 0.005917629852047, 0.004859494981977, 0.002955831081991,
  0.0003356045898205,-0.002697639896313,-0.005687071486917,-0.008084032953745,
  -0.009344469124396,-0.009045419365927,-0.006999963923105,-0.003344508688152,
   0.001426457367602, 0.006496401543416,  0.01083826348381,  0.01339119282049,
	0.01327693896666,   0.0100166598574, 0.003703629355468,-0.004909743327555,
   -0.01444030683204, -0.02302828304832, -0.02860345497317, -0.02922565654424,
   -0.02344372245441, -0.01060796387519, 0.008925525047616,  0.03374481923033,
	0.06150720551817,  0.08922868601273,   0.1136966478108,   0.1319402460404,
	 0.1416800797754,   0.1416800797754,   0.1319402460404,   0.1136966478108,
	0.08922868601273,  0.06150720551817,  0.03374481923033, 0.008925525047616,
   -0.01060796387519, -0.02344372245441, -0.02922565654424, -0.02860345497317,
   -0.02302828304832, -0.01444030683204,-0.004909743327555, 0.003703629355468,
	 0.0100166598574,  0.01327693896666,  0.01339119282049,  0.01083826348381,
   0.006496401543416, 0.001426457367602,-0.003344508688152,-0.006999963923105,
  -0.009045419365927,-0.009344469124396,-0.008084032953745,-0.005687071486917,
  -0.002697639896313,0.0003356045898205, 0.002955831081991, 0.004859494981977,
   0.005917629852047, 0.006161971401866, 0.005744013542788, 0.004881532449878,
   0.003805742549565, 0.002717496869779, 0.001761517078465,  0.00101659692763,
  0.0005019666320485,0.0002088566718795 };//filter函数82阶低通滤波的系数,注意其数据可能存在的对称性

typedef struct _tagMsWavPcmHeader44 {
	unsigned char ChunkID[4]; // "RIFF"; The "RIFF" the mainchunk;
	unsigned long ChunkSize; // FileSize - 8; The size following this data
	unsigned char Format[4];// "WAVE"; The "WAVE" format consists of two subchunks: "fmt " and "data"
	unsigned char SubChunk1ID[4]; // "fmt "
	unsigned long SubChunk1Size;// 16 for PCM. This is the size of the rest of the subchunk which follows this data.
	unsigned short AudioFormat; // 1 for PCM. Linear quantization
	unsigned short NumChannels; // 1->Mono, 2->stereo, etc..
	unsigned long SampleRate; // 8000, 11025, 16000, 44100, 48000, etc..
	unsigned long ByteRate; // = SampleRate * NumChannels * BitsPerSample/8
	unsigned short BlockAlign; // = NumChannels * BitsPerSample / 8
	unsigned short BitsPerSample; // 8->8bits, 16->16bits, etc..
	unsigned char SubChunk2ID[4]; // "data"
	unsigned long SubChunk2Size;// = NumSamples * NumChannels * BitsPerSample / 8. The size of data
} wav_pcm_header44;

void draw(short* data, int len);//画图函数
double filter(short* x, short* y, int leng);//FIR低通滤波函数,82阶,输入为short类型
double RandDouble(double Max);//产生随机浮点数,保留两位小数
double filter2(double* x, double* y, int leng);//FIR低通滤波函数,5阶,输入为double类型

int main()
{
	int nReadLen = 0, i = 0;
	int choice, num;
	double max, Max;
	double x[25], y[25];
	srand((unsigned)time(0));//用于生成不同的随机数

	FILE* pSrcFile = fopen("C:\\Users\\Jungle Chen\\Desktop\\计算机编程class1\\实验四滤波器\\filter\\sin0K5+2K.wav", "rb");//读取需要处理的音频文件
	FILE* pDstFile = fopen("C:\\Users\\Jungle Chen\\Desktop\\计算机编程class1\\实验四滤波器\\filter\\sin0K5+2K+FIR.wav", "wb");//用于写入存储滤波结果的音频文件
	wav_pcm_header44 hwav = { 0 };//初始化
	short pReadData1[FRAME_READ_LENGTH];//用于暂时存放wav文件部分数据
	printf("请选择:\n");
	printf("1.随机产生N个数据进行滤波\n");
	printf("1.对wav文件进行滤波\n");
	scanf("%d", &choice);

	if (choice==1)//产生随机数进行5阶滤波
	{ 

		printf("请输入Max,在0`Max之间产生随机数:\n");
		scanf("%f", &Max);
		printf("随机数如下:\n");
		Max = RandDouble(32.25);//产生随机浮点数,介于0~32.25
		printf("%3.2f  ", Max);
		printf("\n请输入需要滤波随机数的数量:\n");
		scanf("%d", &num);
		for (i = 0; i < num; i++)//产生需要滤波的随机浮点数,介于0~32.25
		{
			x[i] = RandDouble(32.25);
			printf("%3.2f\n", x[i]);
		}
		filter2(x, y, num);//对生成的数据进行滤波
	}

	else if (choice == 2)//读取wav文件进行82阶低通滤波
	{
		if ((pSrcFile == nullptr) || (pDstFile == nullptr))//判断文件是否打开成功
		{
			printf("can't open audio file\n");
			return 0;
		}
		fread(&hwav, sizeof(wav_pcm_header44), 1, pSrcFile);//读取wav文件头并放入hwav中

		printf("-----Information of The File:-------\n");
		printf("ChunkID: %c%c%c%c\n", hwav.ChunkID[0], hwav.ChunkID[1], hwav.ChunkID[2], hwav.ChunkID[3]);
		printf("ChunkSize: %u\n", hwav.ChunkSize);
		printf("Format: %c%c%c%c\n", hwav.Format[0], hwav.Format[1], hwav.Format[2], hwav.Format[3]);
		printf("SubChunk1Size: %u\n", hwav.SubChunk1Size);
		printf("AudioFormat: %u\n", hwav.AudioFormat);
		printf("NumChannels: %u\n", hwav.NumChannels);
		printf("SampleRate: %u\n", hwav.SampleRate);
		printf("ByteRate: %u\n", hwav.ByteRate);
		printf("BlockAlign: %u\n", hwav.BlockAlign);
		printf("BitsPerSample: %u\n", hwav.BitsPerSample);
		printf("SubChunk2ID: %c%c%c%c\n", hwav.SubChunk2ID[0], hwav.SubChunk2ID[1], hwav.SubChunk2ID[2], hwav.SubChunk2ID[3]);
		printf("SubChunk2Size: %u\n", hwav.SubChunk2Size);
		printf("----------end!--------\n");

		fwrite(&hwav, sizeof(wav_pcm_header44), 1, pDstFile);//将原音频文件的文件头直接写出新音频文件内
		fseek(pSrcFile, 44L, SEEK_SET);
		fseek(pDstFile, 44L, SEEK_SET);//将文件指针移到两个文件,距离文件首位置44的位置,即跳过文件头
		int number = 0;

		while (1)//直到判断读取到文件结尾才跳出循环
		{
			nReadLen = fread(pReadData1, sizeof(short), FRAME_READ_LENGTH, pSrcFile);//选择16为一段,读一段处理一段
			if (nReadLen == FRAME_READ_LENGTH)//当读取文件的数据数目等于16时,将数据打印出来
			{
				for (i = 0; i < FRAME_READ_LENGTH; i++)
				{
					printf("%7d", pReadData1[i]);
				}
				printf("\n");
				number++;//16个数据为一段,number是段数,即16*number+最后一次nReadLen认为是文件的数据长度
			}
			else
			{
				for (i = 0; i < nReadLen; i++)
				{
					printf("%7d", pReadData1[i]);
				}
				printf("\n");
				break;//若最后一次读取数据少于16,则判断读取结束,打印剩余的数据然后直接跳出
			}
		}
		short pReadData2[100000];//用于存储wav文件除文件头以外的所有数据
		short pReadData3[100000];//用于存储wav文件滤波后的数据

		fseek(pSrcFile, 44L, SEEK_SET);//定位到文件的44位置,即跳过文件头
		fseek(pDstFile, 44L, SEEK_SET);//定位到文件的44位置,即跳过文件头
		fread(pReadData2, sizeof(short), FRAME_READ_LENGTH * number, pSrcFile);//将文件数据读取到pReadData2
		filter(pReadData2, pReadData3, FRAME_READ_LENGTH * number);//对pReadData2进行滤波
		fwrite(pReadData3, sizeof(short), FRAME_READ_LENGTH * number, pDstFile);//将pReadData3写入到新wav文件中

		draw(pReadData3, 256);//画出wav文件的部分波形图

		fclose(pSrcFile);
		fclose(pDstFile);//关闭文件
		printf("\n----------finished!--------\n");
	}
	else
	{
		printf("error,请退出程序!\n");
	}
	return 0;
}

void draw(short* data, int len) //画图函数
{
	int i = 0;
	const int High = 100;
	const int Len = 500;
	initgraph(640, 480);
	setlinecolor(RED);
	setlinestyle(PS_DASH);
	line(50, High * 1, 600, High * 1);
	line(50, High * 2, 600, High * 2);
	line(50, High * 3, 600, High * 3);
	for (i = 50; i <= 600; i = i + 50)        line(i, High, i, High * 3);
	setlinecolor(GREEN);    setlinestyle(PS_SOLID, 2);
	for (i = 1; i < len; i++) {
		int x1 = i * 2 - 1 + 50, x2 = x1 + 2;
		int y1 = -data[i - 1] / 32768.0 * High + High * 2;
		int y2 = -data[i] / 32768.0 * High + High * 2;
		line(x1, y1, x2, y2);
	}
	getchar();
	closegraph();
	return;
}

double filter(short* x, short* y, int leng)//FIR低通滤波函数,82阶,输入为short类型
{
	int i = 0, n = 0;
	short sum = 0;
	for (i = 0; i < leng; i++)
	{
		sum = 0;
		for (n = 0; n < FIR_FILTER_LENGTH; n++)
		{
			if (i - n >= 0)
			{
				sum = sum + x[i - n] * CN[n];//卷积
			}
		}
		y[i] = sum;
		//printf("%3.10f,", y[i]);
	}
	return 0;
}

double RandDouble(double Max)//产生随机浮点数,保留两位小数
{
	double randnum[20];
	int i;
	int r;
	int num = 10;
	for (i = 0; i < num; i++)
	{
		randnum[i] = 0 + 1.0 * rand() / RAND_MAX * int(max);//生成十个随机浮点数
		//区别在这一行,RAND_MAX虽然可以取得很好的随机效果,但是很难取到闭区间的值,可通过四舍五入函数解决以适应需求
	}	
	return randnum[4];//返回第5个随机浮点数
}

double filter2(double* x, double* y, int leng)//FIR低通滤波函数,5阶,输入为double类型
{
	int i = 0, n = 0;
	double sum = 0; \
		printf("\n5阶FIR低通滤波器滤波处理后的数据:\n");
	for (i = 0; i < leng; i++)
	{
		sum = 0;
		for (n = 0; n < 5; n++)
		{
			if (i - n >= 0)
			{
				sum = sum + x[i - n] * c[n];//卷积
			}
		}
		y[i] = sum;
		printf("%5.5f\n", y[i]);
	}
	return 0;
}


  • 1
    点赞
  • 2
    收藏
    觉得还不错? 一键收藏
  • 1
    评论
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值