基于matlab通过多个FIR,IIR形成对violin信号进行均衡

2024-03-17 12:40

本文主要是介绍基于matlab通过多个FIR,IIR形成对violin信号进行均衡,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

基于《数字信号处理》这门课程的课程设计,特意写下这篇博客,仅供有需要的看官参考。

理论知识:

FIR和IIR滤波器是数字信号处理中常用的滤波器类型。设计带通滤波器时,需要考虑滤波器的通带、阻带、过渡带等参数,以及滤波器的阶数、截止频率等设计参数。

FIR滤波器是一种无限脉冲响应滤波器,其基本理论是通过对输入信号的加权和延迟来实现滤波。FIR滤波器的特点是稳定性好、易于设计和实现。

IIR滤波器是一种有限脉冲响应滤波器,其基本理论是通过对输入信号的加权和延迟以及对输出信号的反馈来实现滤波。IIR滤波器的特点是具有较高的性能和效率,但设计和实现相对复杂。

任务要求:设计三种不同增益的均衡器,分别增强低音、中音和高音

首先,我们需要通过matlab获取violin信号时域图和幅频特性:

clc;
clear;
[x0,fs]=audioread('violin.wav');%获取violin信号的波形数据以及采样率
N=length(x0);%获取音频信号长度
t=0:1/fs:(N-1)/fs;%得到时间轴
faxis=fs/N:fs/N:fs/2-fs/N;
%获取波形的时域图以及频谱图
figure
ffx0=fftshift(abs(fft(x0)));ffx=ffx0(N/2+2:N);
subplot(2,1,1);plot(t,x0,'k-','LineWidth',1.5);box off;
xlim([0 max(t)]);%ylim([-3 3]);
xlabel('Time (s)');ylabel('Amplitude');
subplot(2,1,2);plot(faxis,ffx,'k-','LineWidth',1.5);box off;

b08b054a947c4988a2b9bbfed7ff8cff.png

ad86174f5cc24aa8b550ea5d4260243b.png

接着设计多个FIR滤波器实现均衡器:

clc;
clear;
[x0,fs]=audioread('violin.wav');
N=length(x0);
t=0:1/fs:(N-1)/fs;
faxis=fs/N:fs/N:fs/2-fs/N;figure
ffx0=fftshift(abs(fft(x0)));ffx=ffx0(N/2+2:N);
subplot(2,1,1);plot(t,x0,'k-','LineWidth',1.5);box off;
xlim([0 max(t)]);%ylim([-3 3]);
xlabel('Time(s)');ylabel('Amplitude');
subplot(2,1,2);plot(faxis,ffx,'k-','LineWidth',1.5);box off;%%  滤波器指标
rp=0.1;%rp表示通带最大衰减量
rs=0.001;%rs表示阻带最小衰减量Ap=-20*log10(1-rp);
As=-20*log10(rs);%将其转换为dB%该频率是根据上述幅频特性曲线进行选取
wn1=[40 60 700 720];%设置具有过渡带宽为20的通带和阻带,[60,700]在该频率通过,%小于40,大于720的不允许通过,下述9个与之相同
wn2=[740 760 1100 1120]; 
wn3=[1180 1200 1400 1420];
wn4=[1500 1520 1800 1820];  
wn5=[1840 1860 2000 2020]; 
wn6=[2040 2060 2340 2360];  
wn7=[2380 2400 2700 2720]; 
wn8=[2780 2800 3100 3120];  
wn9=[3270 3290 3600 3620];  
%% FIR 均衡器设计 用firpm
%%% 计算滤波器参数 firpmord
mags=[0 1 0];%描述带通
devs=[rs rp rs];%用于计算滤波器参数的值
[n1,f01,ao1,w1]=firpmord(wn1,mags,devs,fs);%计算滤波器参数
[n2,f02,ao2,w2]=firpmord(wn2,mags,devs,fs);
[n3,f03,ao3,w3]=firpmord(wn3,mags,devs,fs);
[n4,f04,ao4,w4]=firpmord(wn4,mags,devs,fs);
[n5,f05,ao5,w5]=firpmord(wn5,mags,devs,fs);
[n6,f06,ao6,w6]=firpmord(wn6,mags,devs,fs);
[n7,f07,ao7,w7]=firpmord(wn7,mags,devs,fs);
[n8,f08,ao8,w8]=firpmord(wn8,mags,devs,fs);
[n9,f09,ao9,w9]=firpmord(wn9,mags,devs,fs);
%%% 计算滤波器系数
hh1=firpm(n1,f01,ao1,w1);%计算滤波器系数
hh2=firpm(n2,f02,ao2,w2);
hh3=firpm(n3,f03,ao3,w3);
hh4=firpm(n4,f04,ao4,w4);
hh5=firpm(n5,f05,ao5,w5);
hh6=firpm(n6,f06,ao6,w6);
hh7=firpm(n7,f07,ao7,w7);
hh8=firpm(n8,f08,ao8,w8);
hh9=firpm(n9,f09,ao9,w9);
%%% 画出滤波器幅频特性曲线
[h1 W1]=freqz(hh1,1,1024);
[h2,W2]=freqz(hh2,1,1024);
[h3,W3]=freqz(hh3,1,1024);
[h4,W4]=freqz(hh4,1,1024);
[h5,W5]=freqz(hh5,1,1024);
[h6,W6]=freqz(hh6,1,1024);
[h7,W7]=freqz(hh7,1,1024);
[h8,W8]=freqz(hh8,1,1024);
[h9,W9]=freqz(hh9,1,1024);
figure; %绘图,均衡器幅频响应
loglog(W1,abs(h1),W2,abs(h2), W3,abs(h3),...
W4,abs(h4),W5,abs(h5),W6,abs(h6),W7,abs(h7),...
W8,abs(h8),W9,abs(h9));
xlabel('Frequency (Hz)r) ;ylabel(filter Gain');
axis([10 10^5   10^(-9) 1]);
%%%  对音乐信号进行均衡并画出频谱
y1=filter(hh1,1,x0);%设计滤波器对x0信号进行滤波
y2=filter(hh2,1,x0);
y3=filter(hh3,1,x0);
y4=filter(hh4,1,x0);
y5=filter(hh5,1,x0);
y6=filter(hh6,1,x0);
y7=filter(hh7,1,x0);
y8=filter(hh8,1,x0);
y9=filter(hh9,1,x0);
g1= 10;g2 = 10;g3 = 10;g4 = 0;
g5 = 0;g6 = 0;g7 = 0;g8 = 0;g9 = 0;
y=y1.*g1+y2.*g2+y3.*g3+y4.*g4+y5.*g5+y6.*g6+y7.*g7+y8.*g8+y9.*g9+x0;
%对前三个信号即低音进行均衡则选取g1,g2,g3设置为10,依次向后推理获得中音和高音均衡
%对均衡信号画出时域和频域波形图
figure
ffy=fftshift(abs(fft(y)));ffy0=ffy(N/2+2:N);
subplot(2,1,1);plot(t,y,'k-','LineWidth',1.5);box off;
xlim([0 max(t)]);%ylim([-3 3]);
xlabel('Time(s)');ylabel('Amplitude');
subplot(2,1,2);plot(faxis,ffy0,'k-','LineWidth',1.5);box off;
%%% 听效果并存储
sound(y);

该均衡器获得的幅频特性曲线如下图:

02a359560ef54b1dabb6ce8175f69743.png

对低音进行均衡后信号的时域频域波形图如下:

44d9b95a932940dca775192c9196d199.png

17f3288000574c01b73dbc61ad805a7b.png

 

对中音进行均衡后信号的时域频域波形图如下:

6bfbd36be1354a0baaa53675c2f20eb4.png

8d18145638a74adbbf927a5790bec49e.png

 

对高音进行均衡后信号的时域频域波形图如下:

f93e4e1c49734fa9a0e6eab421cc02ef.png

967b8dca894848639e3943ffcb2c72ea.png

基于matlab的多个IIR设计实现均衡器:

clc;
clear;
[x0,fs]=audioread('violin.wav');
N=length(x0);
t=0:1/fs:(N-1)/fs;
faxis=fs/N:fs/N:fs/2-fs/N;figure
ffx0=fftshift(abs(fft(x0)));ffx=ffx0(N/2+2:N);
subplot(2,1,1);plot(t,x0,'k-','LineWidth',1.5);box off;
xlim([0 max(t)]);%ylim([-3 3]);
xlabel('Time (s)');ylabel('Amplitude');
subplot(2,1,2);plot(faxis,ffx,'k-','LineWidth',1.5);box off;%% 滤波器指标
Ap=1;
As=20;wp1=[60 700]; ws1=[40 720];
wp2=[760 1100]; ws2=[740 1120];
wp3=[1200 1400]; ws3=[1180 1420];
wp4=[1520 1800]; ws4=[1500 1820];
wp5=[1860 2000]; ws5=[1840 2020];
wp6=[2060 2340]; ws6=[2040 2360];
wp7=[2400 2700]; ws7=[2380 2720];
wp8=[2800 3100]; ws8=[2780 3120];
wp9=[3290 3600]; ws9=[3270 3620];
wp10=[3660 6000];ws10=[3640 6020];%% IIR 均衡器设计 ellip
wpI1=2*wp1/fs;wsI1=2*ws1/fs;%将频率进行归一化
wpI2=2*wp2/fs;wsI2=2*ws2/fs;
wpI3=2*wp3/fs;wsI3=2*ws3/fs;
wpI4=2*wp4/fs;wsI4=2*ws4/fs;
wpI5=2*wp5/fs;wsI5=2*ws5/fs;
wpI6=2*wp6/fs;wsI6=2*ws6/fs;
wpI7=2*wp7/fs;wsI7=2*ws7/fs;
wpI8=2*wp8/fs;wsI8=2*ws8/fs;
wpI9=2*wp9/fs;wsI9=2*ws9/fs;
wpI10=2*wp10/fs;wsI10=2*ws10/fs;%%% 计算滤波器参数 ellipord
[n1,wn1]=ellipord(wpI1,wsI1,Ap,As);Niir1=n1;
[n2,wn2]=ellipord(wpI2,wsI2,Ap,As);Niir2=n2;
[n3,wn3]=ellipord(wpI3,wsI3,Ap,As);Niir3=n3;
[n4,wn4]=ellipord(wpI4,wsI4,Ap,As);Niir4=n4;
[n5,wn5]=ellipord(wpI5,wsI5,Ap,As);Niir5=n5;
[n6,wn6]=ellipord(wpI6,wsI6,Ap,As);Niir6=n6;
[n7,wn7]=ellipord(wpI7,wsI7,Ap,As);Niir7=n7;
[n8,wn8]=ellipord(wpI8,wsI8,Ap,As);Niir8=n8;
[n9,wn9]=ellipord(wpI9,wsI9,Ap,As);Niir9=n9;
[n10,wn10]=ellipord(wpI10,wsI10,Ap,As);Niir10=n10;%%% 计算滤波器系数
[b1,a1]=ellip(n1,Ap,As,wn1);
[b2,a2]=ellip(n2,Ap,As,wn2);
[b3,a3]=ellip(n3,Ap,As,wn3);
[b4,a4]=ellip(n4,Ap,As,wn4);
[b5,a5]=ellip(n5,Ap,As,wn5);
[b6,a6]=ellip(n6,Ap,As,wn6);
[b7,a7]=ellip(n7,Ap,As,wn7);
[b8,a8]=ellip(n8,Ap,As,wn8);
[b9,a9]=ellip(n9,Ap,As,wn9);
[b10,a10]=ellip(n10,Ap,As,wn10);%%% 画出滤波器幅频特性曲线
[h1,W1]=freqz(b1,a1);
[h2,W2]=freqz(b2,a2);
[h3,W3]=freqz(b3,a3);
[h4,W4]=freqz(b4,a4);
[h5,W5]=freqz(b5,a5);
[h6,W6]=freqz(b6,a6);
[h7,W7]=freqz(b7,a7);
[h8,W8]=freqz(b8,a8);
[h9,W9]=freqz(b9,a9);
[h10,W10]=freqz(b10,a10);figure; %绘图,均衡器幅频响应
loglog(W1,abs(h1),W2,abs(h2), W3,abs(h3),...
W4,abs(h4),W5,abs(h5),W6,abs(h6),W7,abs(h7),...
W8,abs(h8),W9,abs(h9),W10,abs(h10));
xlabel('Frequency (Hz)'); ylabel('Filter Gain');
axis([10 10^5 10^-9 1]);%%% 对音乐信号进行均衡并画出频谱
y1=filter(b1,a1,x0);
y2=filter(b2,a2,x0);
y3=filter(b3,a3,x0);
y4=filter(b4,a4,x0);
y5=filter(b5,a5,x0);
y6=filter(b6,a6,x0);
y7=filter(b7,a7,x0);
y8=filter(b8,a8,x0);
y9=filter(b9,a9,x0);
y10=filter(b10,a10,x0);
g1=0; g2=0; g3=0; g4=0;
g5=10; g6=10; g7=10; g8=0; g9=0; g10=0;
y = y1.*g1 + y2.*g2 + y3.*g3 + y4.*g4 + y5.*g5 +...
y6.*g6 + y7.*g7 + y8.*g8 + y9.*g9 + y10.*g10 + x0;figure
ffy=fftshift(abs(fft(y)));ffy0=ffy(N/2+2:N);
subplot(2,1,1);plot(t,y,'k-','LineWidth',1.5);box off;
xlim([0 max(t)]);%ylim([-3 3]);
xlabel('Time(s)');ylabel('Amplitude');
subplot(2,1,2);plot(faxis,ffy0,'k-','LineWidth',1.5);box off;
%%% 听效果并存储
sound(y,fs);

该均衡器获得的幅频特性曲线如下图:

241bf65e6fe84e3bbbd1e0a49e4b1908.png

IIR设计均衡器对低音进行增强后均衡处理后信号的时域,频域图如下:

d111a531f89f43128b09b2ae1240b699.png

3420d689da1a41ed9261aa1ad4f19164.png

IIR设计均衡器对中音进行增强后均衡处理后信号的时域,频域图如下:

dd6e3ecb5a6b46ff8c857dba4c2d4ad5.png

9b836e5d18e34b889705a221c2d11924.png

IIR设计均衡器对高音进行增强后均衡处理后信号的时域,频域图如下:

216ec8d45fb54973a9066a6f66e9bb45.png

c892441676a24c6f9134bb411e2c1dc6.png

上述代码有点呆(很早之前写的代码),各位优化代码可通过数组或for循环进行简化,

另外,若存在错误,请各位多多指正!

 

 

这篇关于基于matlab通过多个FIR,IIR形成对violin信号进行均衡的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/818979

相关文章

Hadoop集群数据均衡之磁盘间数据均衡

生产环境,由于硬盘空间不足,往往需要增加一块硬盘。刚加载的硬盘没有数据时,可以执行磁盘数据均衡命令。(Hadoop3.x新特性) plan后面带的节点的名字必须是已经存在的,并且是需要均衡的节点。 如果节点不存在,会报如下错误: 如果节点只有一个硬盘的话,不会创建均衡计划: (1)生成均衡计划 hdfs diskbalancer -plan hadoop102 (2)执行均衡计划 hd

【Prometheus】PromQL向量匹配实现不同标签的向量数据进行运算

✨✨ 欢迎大家来到景天科技苑✨✨ 🎈🎈 养成好习惯,先赞后看哦~🎈🎈 🏆 作者简介:景天科技苑 🏆《头衔》:大厂架构师,华为云开发者社区专家博主,阿里云开发者社区专家博主,CSDN全栈领域优质创作者,掘金优秀博主,51CTO博客专家等。 🏆《博客》:Python全栈,前后端开发,小程序开发,人工智能,js逆向,App逆向,网络系统安全,数据分析,Django,fastapi

业务中14个需要进行A/B测试的时刻[信息图]

在本指南中,我们将全面了解有关 A/B测试 的所有内容。 我们将介绍不同类型的A/B测试,如何有效地规划和启动测试,如何评估测试是否成功,您应该关注哪些指标,多年来我们发现的常见错误等等。 什么是A/B测试? A/B测试(有时称为“分割测试”)是一种实验类型,其中您创建两种或多种内容变体——如登录页面、电子邮件或广告——并将它们显示给不同的受众群体,以查看哪一种效果最好。 本质上,A/B测

matlab读取NC文件(含group)

matlab读取NC文件(含group): NC文件数据结构: 代码: % 打开 NetCDF 文件filename = 'your_file.nc'; % 替换为你的文件名% 使用 netcdf.open 函数打开文件ncid = netcdf.open(filename, 'NC_NOWRITE');% 查看文件中的组% 假设我们想读取名为 "group1" 的组groupName

遮罩,在指定元素上进行遮罩

废话不多说,直接上代码: ps:依赖 jquer.js 1.首先,定义一个 Overlay.js  代码如下: /*遮罩 Overlay js 对象*/function Overlay(options){//{targetId:'',viewHtml:'',viewWidth:'',viewHeight:''}try{this.state=false;//遮罩状态 true 激活,f

利用matlab bar函数绘制较为复杂的柱状图,并在图中进行适当标注

示例代码和结果如下:小疑问:如何自动选择合适的坐标位置对柱状图的数值大小进行标注?😂 clear; close all;x = 1:3;aa=[28.6321521955954 26.2453660695847 21.69102348512086.93747104431360 6.25442246899816 3.342835958564245.51365061796319 4.87

C# double[] 和Matlab数组MWArray[]转换

C# double[] 转换成MWArray[], 直接赋值就行             MWNumericArray[] ma = new MWNumericArray[4];             double[] dT = new double[] { 0 };             double[] dT1 = new double[] { 0,2 };

一种改进的red5集群方案的应用、基于Red5服务器集群负载均衡调度算法研究

转自: 一种改进的red5集群方案的应用: http://wenku.baidu.com/link?url=jYQ1wNwHVBqJ-5XCYq0PRligp6Y5q6BYXyISUsF56My8DP8dc9CZ4pZvpPz1abxJn8fojMrL0IyfmMHStpvkotqC1RWlRMGnzVL1X4IPOa_  基于Red5服务器集群负载均衡调度算法研究 http://ww

Python脚本:对文件进行批量重命名

字符替换:批量对文件名中指定字符进行替换添加前缀:批量向原文件名添加前缀添加后缀:批量向原文件名添加后缀 import osdef Rename_CharReplace():#对文件名中某字符进行替换(已完结)re_dir = os.getcwd()re_list = os.listdir(re_dir)original_char = input('请输入你要替换的字符:')replace_ch

列举你能想到的UNIX信号,并说明信号用途

信号是一种软中断,是一种处理异步事件的方法。一般来说,操作系统都支持许多信号。尤其是UNIX,比较重要应用程序一般都会处理信号。 UNIX定义了许多信号,比如SIGINT表示中断字符信号,也就是Ctrl+C的信号,SIGBUS表示硬件故障的信号;SIGCHLD表示子进程状态改变信号;SIGKILL表示终止程序运行的信号,等等。信号量编程是UNIX下非常重要的一种技术。 Unix信号量也可以