平方包络信号仿真

平方包络信号的仿真 帮助大家学者仿真平方包络信号 % 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('平方包络信号的时域波形');