实验目的
利用C语言集成开发环境设计FIR滤波器软件,并进行测试:
- 练习随机数、静态局部变量、滤波器等知识的综合运用;
- 培养学生综合运用C程序语言进行程序设计的能力;
- 培养结构化程序设计能力;
- 分析FIR低通滤波器的特点;
- 分析wav文件格式;
- 设计fir低通滤波器
- 用C或C++编程实现滤波器;
- 调式程序,测试滤波器效果
实验内容
上机完成以下编程实验,调试运行程序并完成报告。
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画出的滤波前的频谱和波形图
以上为使用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;
}