资源介绍
平方包络信号的仿真 帮助大家学者仿真平方包络信号
% squaredd1
% 调制函数呈衰减趋势
% 滤波法
% 每阶啮合谐波加入两阶调制谐波
clear all
% 冲击形状按照文章中的表达式
N=5120;
t=1:N;
fs=1000;
t=t/fs;
u=0.925;
deta=0.05;
d0=exp(-(log(t)-u).^2./(deta.^2))./(t.*deta.*sqrt(2*pi));
a11=0.12.*exp(-0.15.*t).*cos(2*pi.*t);
a12=0.08.*exp(-0.15.*t).*cos(2*pi.*t);
a1=a11+a12;
b11=0.12.*exp(-0.15.*t).*cos(2*pi.*t);
b12=0.08.*exp(-0.15.*t).*cos(2*pi.*t);
b1=b11+b12;
d=1.45.*d0./max(abs(d0));
a1=a1+d/121;
b1=b1-d/121;
a21=0.12.*exp(-0.15.*t).*cos(4*pi.*t);
a22=0.08.*exp(-0.15.*t).*cos(4*pi.*t);
a2=a21+a22;
b21=0.12.*exp(-0.15.*t).*cos(4*pi.*t);
b22=0.08.*exp(-0.15.*t).*cos(4*pi.*t);
b2=b21+b22;
a2=a2+d/121;
b2=b2-d/121;
a31=0.12.*exp(-0.15.*t).*cos(6*pi.*t);
a32=0.08.*exp(-0.15.*t).*cos(6*pi.*t);
a3=a31+a32;
b31=0.12.*exp(-0.15.*t).*cos(6*pi.*t);
b32=0.08.*exp(-0.15.*t).*cos(6*pi.*t);
b3=b31+b32;
a41=0.12.*exp(-0.15.*t).*cos(8*pi.*t);
a42=0.08.*exp(-0.15.*t).*cos(8*pi.*t);
a4=a41+a42;
b41=0.12.*exp(-0.15.*t).*cos(8*pi.*t);
b42=0.08.*exp(-0.15.*t).*cos(8*pi.*t);
b4=b41+b42;
a51=0.12.*exp(-0.15.*t).*cos(10*pi.*t);
a52=0.08.*exp(-0.15.*t).*cos(10*pi.*t);
a5=a51+a52;
b51=0.12.*exp(-0.15.*t).*cos(10*pi.*t);
b52=0.08.*exp(-0.15.*t).*cos(10*pi.*t);
b5=b51+b52;
r1=10*(a1.*cos(54*pi*t)-b1.*sin(54*pi*t));
r2=8*(a2.*cos(108*pi*t)-b2.*sin(108*pi*t));
r3=6*(a3.*cos(162*pi*t)-b3.*sin(162*pi*t));
r4=5*(a4.*cos(216*pi*t)-b4.*sin(216*pi*t));
r5=2*(a5.*cos(270*pi*t)-b5.*sin(270*pi*t));
yr=r1+r2+r3+r4+r5+d.*cos(242*pi*t);
Y=fft(yr);
fil=[zeros(1,610) Y(611:624) zeros(1,4496)];%一共5120个点
fy1=ifft(fil);
n1=10;wn1=[110 135]/(fs/2); % 带通滤波
[b,a]=butter(n1,wn1);
fy1=filter(b,a,yr);
y=fy1.*fy1;
Y=fft(y);
fil=[Y(1:50) zeros(1,5070)];%一共5120个点
y=abs(ifft(fil));
t=t*360*1000/5120;
plot(t,y);
xlabel('轴转角/度');
ylabel('Am/mm');
title('平方包络信号的时域波形');