186
社区成员
发帖
与我相关
我的任务
分享clc,clear,close all;
%matlab滑动平均滤波器基础例程:
x = linspace(-pi,pi,100);
y = sin(2*pi*x);
n = randn(size(x));
t = y+n;
y1 = movmean(t,8);
figure
subplot(211);
plot(x,t,'LineWidth',1.1);
grid on
title('带噪声输入信号');
subplot(212);
plot(x,y1,'LineWidth',1.1);
grid on
title('平滑滤波后的信号');
% 滑动平均滤波器的幅频响应和相频响应
% 定义滤波器的阶数
Ns = [3,11,31];
for N = Ns
% 定义采样频率和截止频率
fs = 1000; % 采样频率为1000Hz
fc = 50; % 截止频率为50Hz
% 计算滤波器系数
b = ones(1,N)/N; % 滑动平均滤波器的系数为1/N
% 计算频率响应
[H,f] = freqz(b,1,[],fs);
% 计算幅频响应和相频响应
mag = 20*log10(abs(H)); % 幅频响应
phase = unwrap(angle(H))*(180/pi); % 相频响应,转化为角度制
% 绘制幅频响应曲线
figure;
subplot(211)
plot(f,mag,'LineWidth',1.1);
grid on
xlabel('频率 (Hz)');
ylabel('幅值 (dB)');
title(sprintf('N=%d的滑动平均滤波器幅频响应', N));
% 绘制相频响应曲线
subplot(212)
plot(f,phase,'LineWidth',1.1);
grid on
xlabel('频率 (Hz)');
ylabel('相位 (°)');
title(sprintf('N=%d的滑动平均滤波器相频响应', N));
end
%读取"Don_Giovanni.wav"并进行滑动平均滤波测试
[music0,Fs] = audioread('Don_Giovanni.wav');
player = audioplayer(music0,Fs);
play(player)
% stop(player) %停止播放
pause(abs(length(music0)/Fs))
disp('原始音频播放完毕')
music3 = movmean(music0,31);
player = audioplayer(music3,Fs); %N=3的效果
play(player)
pause(abs(length(music3)/Fs))
disp('N=3滑动平均滤波音频播放完毕')
music11 = movmean(music0,11);
player = audioplayer(music11,Fs); %N=11的效果
play(player)
pause(abs(length(music11)/Fs))
disp('N=11滑动平均滤波音频播放完毕')
music31 = movmean(music0,31);
player = audioplayer(music31,Fs); %N=31的效果
play(player)
pause(abs(length(music31)/Fs))
disp('N=31滑动平均滤波音频播放完毕')



clc,clear,close all;
%**************************************************************
% Filtering a signal by a filter with a transfer function
% B(z) b0+b1z^{-1}+b2z^{-2}+ ... +bMz^{-M}
% H(z)= ----- = -------------------------------------
% A(z) a0+a1z^{-1}+a2z^{-2}+ ... +aNz^{-N}
%**************************************************************
% General parameters
Fe=44100; %采样频率
Te=1/Fe; %周期
N=Fe; %序列长度
t=0:Te:(N-1)*Te;
f=-Fe/2:Fe/N:Fe/2-Fe/N;
% 1st order IIR filter
a=-0.5;
b0 = 1+a;
b1 = 0;
a0 = 1;
a1 = a;
B1=[b0 b1];
A1=[a0 a1];
% Filtering a step signal
echelon=ones(1,N);
rep_echelon=filter(B1,A1,echelon);
% Synthesis of a square signal
fc=250;
carre=1+square(2*pi*fc*t);
rep_carre=filter(B1,A1,carre);
Npts=floor(Fe/fc);
figure(1)
%阶跃信号滤波
stem(t(1:Npts),echelon(1:Npts),'k');
grid on; hold on;
stem(t(1:Npts),rep_echelon(1:Npts),'r');
xlabel('Time (s)');
ylabel('Amplitude');
legend('Input signal','Output signal')
title('1st order filter : step signal filtering')
%方波信号滤波
figure(2)
stem(t(1:Npts),carre(1:Npts),'k');
grid on; hold on;
stem(t(1:Npts),rep_carre(1:Npts),'r');
xlabel('Time (s)');
ylabel('Amplitude');
legend('Input signal','Output signal')
title('1st order filter : square signal filtering')


% Filtering a signal by a filter with a transfer function
% B(z) b0+b1z^{-1}+b2z^{-2}+ ... +bMz^{-M}
% H(z)= ----- = -------------------------------------
% A(z) a0+a1z^{-1}+a2z^{-2}+ ... +aNz^{-N}
%**************************************************************
% 1st order IIR filter 一阶IIR滤波器
分子分母系数:
b0 = 0.5
b1 = 0
a0 = 1
a1 = -0.5
系统函数:
H(z) = (1/2+0)/(1-1/2z^-1)
= 1/2/(1-1/2z^-1)
= 1/(2-z^-1)
= z/(2*z-1) = Y(z)/X(z) (传递函数:2y(n) - y(n-1) = x(n))
零极点分析
一阶零点:z = 0
一阶极点:z = 0.5
稳定性:系统稳定时的收敛域为z>0.5,收敛域包含单位圆,故该系统是稳定的。
IIR滤波器可以对信号进行低通、高通、带通或带阻滤波。
低通滤波器可以通过去除高频信号来平滑信号;
高通滤波器可以通过去除低频信号来突出高频信号;
带通和带阻滤波器可以选择性地保留或去除一定范围内的频率分量;
在本代码中,通过part2的波形图可以看出信号跳变部分变的更加平滑(变化更加缓慢),
使用的是一阶低通IIR滤波器(w=0(z=1)时,H(z)不为0.),其目的是平滑输入信号。
因此,滤波后的输出信号可能会比输入信号更平滑。
clc,clear,close all;
Fe=44100; %采样频率
Te=1/Fe; %周期
N=Fe; %序列长度
t=0:Te:(N-1)*Te;
f=-Fe/2:Fe/N:Fe/2-Fe/N;
% 2st order IIR filter 二阶IIR
a_1=-1.2;
a_2=.3;
b0 = 1+a_1+a_2;
b1 = 0;
b2 = 0;
a0 = 1;
a1 = a_1;
a2 = a_2;
B2=[b0 b1 b2];
A2=[a0 a1 a2];
% Filtering a step signal
echelon=ones(1,N);
rep_echelon=filter(B2,A2,echelon);
% Synthesis of a square signal
fc=250;
carre=1+square(2*pi*fc*t);
rep_carre=filter(B2,A2,carre);
% Sinusoidal signal
n = 1:N;
sin_signal=sin(2*pi*1000.*n./Fe);
rep_sin_signal=filter(B2,A2,sin_signal);
Npts=floor(Fe/fc);
figure(1)
%阶跃信号滤波
stem(t(1:Npts),echelon(1:Npts),'k');
grid on; hold on;
stem(t(1:Npts),rep_echelon(1:Npts),'r');
xlabel('Time (s)');
ylabel('Amplitude');
legend('Input signal','Output signal')
title('2st order filter : step signal filtering')
%方波信号滤波
figure(2)
stem(t(1:Npts),carre(1:Npts),'k');
grid on; hold on;
stem(t(1:Npts),rep_carre(1:Npts),'r');
xlabel('Time (s)');
ylabel('Amplitude');
legend('Input signal','Output signal')
title('2st order filter : square signal filtering')
%正弦信号滤波
figure(3)
stem(t(1:Npts),sin_signal(1:Npts),'k');
grid on; hold on;
stem(t(1:Npts),rep_sin_signal(1:Npts),'r');
xlabel('Time (s)');
ylabel('Amplitude');
legend('Input signal','Output signal')
title('2st order filter : Sinusoidal signal filtering')
% 音乐信号滤波
[music0,Fs] = audioread('Don_Giovanni.wav');
player = audioplayer(music0,Fs);
play(player)
% stop(player) %停止播放
pause(abs(length(music0)/Fs))
disp('原始音频播放完毕')
music1 = filter(B2,A2,music0);
player = audioplayer(music1,Fs);
play(player)
pause(abs(length(music1)/Fs))
disp('IIR滤波音频播放完毕')
figure(4)
subplot(2,1,1)
plot(music0)
grid on
title('原始音乐信号')
subplot(2,1,2)
plot(music1)
grid on
title('IIR滤波后的音乐信号')

