由于工作内容涉及信号系统、信号处理相关知识,本人本硕均为计算机相关专业,专业、研究方向均未涉及信号相关知识,因此需进行系统地学习。之前已将《信号与系统》快速过了一遍,但感觉较抽象且理解较浅显。在此系统地学习如何使用Matlab进行信号处理,以加深对信号相关知识的理解。
一个简单的示例:
设计Matlab程序,去掉左图信号中大于20Hz的频率,使其变为右图。
% 1.生成原信号并绘图
% 生成10s的信号,每隔0.0001s生成1次
t1 = 0:0.0001:10;
% 原信号
y1 = 3*sin(2*pi*10*t1)+sin(2*pi*40*t1)+sin(2*pi*200*t1);
% 原信号时域图
subplot(321);
plot(t1,y1);
% 设置横纵轴范围
xlim([0 1]);
ylim([-6 6]);
% 设置图名、横纵轴标签
title('原信号时域');
xlabel('时间(s)');
ylabel('幅值');
% 设置横纵轴分度值
set(gca,'XTick',0:0.25:1);
set(gca,'YTick',-6:3:6);
% 2.模拟采样
% 采样时长1s,采样频率100Hz
t2 = 0:0.01:1;
y2 = 3*sin(2*pi*10*t2)+sin(2*pi*40*t2)+sin(2*pi*200*t2);
% 采样信号图
subplot(322);
plot(t2,y2);
% 设置横纵轴范围
xlim([0 1]);
ylim([-6 6]);
% 设置图名、横纵轴标签
title('采样信号波形');
xlabel('时间(s)');
ylabel('幅值');
% 设置横纵轴分度值
set(gca,'XTick',0:0.25:1);
set(gca,'YTick',-6:3:6);
% 3.fft
% 采样后信号的数据长度
dataLength = length(y2);
% 对采样后的信号进行fft
y2FFT = fft(y2,dataLength);
% 求幅值
mag = abs(y2FFT);
% 幅值归一化
mag = mag*2/dataLength;
pha = angle(y2FFT)*180/pi;
% 无效相位置0
for i = 1:dataLength
if (mag(1,i)<0.3)
pha(1,i) = 0;
end
end
% fft后的序列坐标
n = 0:dataLength-1;
% 采样频率
fs = 1/0.01;
% 序列频率
f = (0:dataLength-1)*fs/dataLength;
% fft频域图
% 绘制幅度谱,stem:绘制离散信号图
subplot(323);
stem(f(1:dataLength/2),mag(1:dataLength/2));
xlim([0 50]);
ylim([0 4]);
set(gca,"XTick",0:10:50);
xlabel("频率(Hz)");
title('幅度谱');
% 绘制相位谱
subplot(324);
stem(f(1:dataLength/2),pha(1:dataLength/2));
xlim([0 50]);
set(gca,'XTick',0:10:50);
xlabel("频率(Hz)");
ylabel("相位");
title('相位谱');
% 低通滤波器滤波
% 载入设计的滤波器
load("my_filter.mat");
% 滤波
% 滤波器系数存放于filtercoe数组中
filtercoe = my_filter;
y_Filtered = filter(filtercoe,1,y2);
% 滤波后的时域图
subplot(326);
plot(t2,y_Filtered);
xlim([0 1]);
ylim([-6 6]);
set(gca,'XTick',0:0.25:1);
set(gca,'YTick',-6:3:6);
title('滤波后时域');
xlabel('时间(s)');
ylabel('幅值');
注:
1. 第95行低通滤波器的设计需要在命令行窗口输入filterDesigner以打开滤波器设计窗口,如下图:
在窗口中设置如下红框中主要参数,然后点击设计滤波器:
生成成后,会出现粉框中的波形。
点击文件-导出,点击弹出窗口中的导出。
此时滤波器已经导出到工作区中。
右击工作区中导出的滤波器,另存为,之后便可通过该滤波器的名称导入到程序中使用。
运行效果:
遗留问题:
1. 图3中存在频谱泄露,后续学习如何解决;
2. 程序中的相关设计仍在学习;