MATLAB滤波程序

% 这是时域积分的程序，原始信号为速度信号y1
% x_velocity  速度信号
% D 数据长度
% time_series 分析时间序列
% total_time  采样时间
%mean_velocity  速度均值
% spec_area 速度频谱
%freq_area 分析频率序列

%有效值 Xrms
close all
clear all
clc;
%=========================================================================%
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  信号产生   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%=========================================================================%
N=1024;
fs=250;
fc=40;
fb=80;
t=（0:1/fs:（N-1）/fs）;
y0=10*sin（2*pi*fc*t+200*pi/180）+20*cos（2*pi*fb*t+150*pi/180）+6*sin（2*pi*0.5*t+150*pi/180）+0.15;  %模拟信号作为加速度信号
%加白噪声，a=0;b=0.5;  %均值为a，方差为b^2
%y2=a+b*randn（1length（y1）） 高斯白噪声 title（‘N（0，0.25）的高斯白噪声‘）;
n=length（y0）;
y1=wgn（1n4）;                                                             %产生高斯白噪声
y=y0+y1;                                                                   %加入噪声之后的信号
x_accr=y;
D=length（x_accr）;                                                          %数据长度
time_series=（1/fs:1/fs:D/fs）;                                              %分析时间序列也可time_series=（0:D-1）/fs;                                                %分析时间序列的另一种表示
time_series=time_series‘;                                         
total_time=D/fs;                                                           %采样时间，作为横坐标的范围
N=2^nextpow2（D）;                                                           %确定分析频率范围的点数:取最接近D的较大2次幂

spec_y0=fft（y0N）;
spec_y=fft（yN）;
freq_series=fs*（0:N/2-1）/N;

figure（4）
subplot（311）;plot（ty1‘r‘）;
axis（[0 1.1 -5  5]）;
ylabel（‘y（um）‘）;xlabel（‘t（s）‘）;
title（‘强度为4dwb的高斯白噪声‘）;
grid on;
figure（1）
subplot（421）;plot（ty0）;
ylabel（‘y0（um）‘）;xlabel（‘t（s）‘）;
title（‘理想信号‘）;
grid on;
subplot（422）;plot（freq_seriesabs（spec_y0（1:N/2））*2/N）;
xlabel（‘频率（Hz）‘）;
ylabel（‘幅值（mm/g）‘）;
title（‘理想加速度幅频谱‘）;
grid on;
subplot（423）;plot（ty1）;
ylabel（‘y（um）‘）;xlabel（‘t（s）‘）;
title（‘原始信号（含噪）‘）;
grid on;
subplot（424）;plot（freq_seriesabs（spec_y（1:N/2））*2/N）;
xlabel（‘频率（Hz）‘）;
ylabel（‘幅值（mm/g）‘）;
title（‘原始加速度幅频谱‘）;
grid on;

%=========================================================================%
 %%%%%%%%%%%%%%%%%%%%%%%%%  加速度信号滤波处理   %%%%%%%%%%%%%%%%%%%%%%%%%%%
%=========================================================================%                                                         
mean_accr=mean（x_accr）;                                                    %去直流
disp（‘实际的直流量‘）;                         
disp（mean_accr）;
x_accr=x_accr-mean_accr;    
spec_accr=fft（x_accrN）;                                   
freq_series=fs*（0:N/2-1）/N;                                                %建立分析频率序列
% 输出加速度时域波形和频谱图
subplot（425）;
plot（time_seriesx_accr‘r‘）;
ax=max（abs（x_accr））;
ra=ax+1;
axis（[0 total_time -ra ra]）;
xlabel（‘采样时间（s）‘）;
ylabel（‘幅值（mm/g）‘）;
title（‘振动加速度，已去直流‘）;                    
grid on;
subplot（426）;
plot（freq_seriesabs（spec_accr（1:N/2））*2/N‘r‘）; 
xlabel（‘频率（Hz）‘）;
ylabel（‘幅值（mm/g）‘）;
title（‘去直流的加速度幅频谱‘）;
grid on;
%设计高通滤波器