【轴承故障检测】【借助倒谱预白化技术在变速条件下诊断轴承故障的应用】带通滤波后的倒谱预白化的平方包络谱用于轴承故障检测（Matlab代码实现）

 💥💥💞💞欢迎来到本博客❤️❤️💥💥

🏆博主优势：🌞🌞🌞博客内容尽量做到思维缜密，逻辑清晰，为了方便读者。

⛳️座右铭：行百里者，半于九十。

📋📋📋本文目录如下：🎁🎁🎁

目录

💥1 概述

📚2 运行结果

🎉3 参考文献

🌈4 Matlab代码及文章

---

💥1 概述

文献来源：

这个函数可以用于常速滚动轴承的基于振动的故障诊断。这是一个三步骤的过程：
(i) 倒谱预白化：减少其他周期性来源（如齿轮）的贡献。
(ii) 带通滤波：提高信噪比，特别是在系统共振附近执行。
(iii) 平方包络谱：允许检测（伪）循环平稳贡献，其特征是在特定循环频率处具有大的分量。
摘要
滚动轴承的诊断涉及信号增强和分析的各种技术的组合。最常见的方法是首先进行阶次跟踪和同步平均，能够从信号中去除与轴谐波同步的不需要的分量，然后进行包络分析以获得平方包络谱。这个指标已经得到了深入研究，并获得了基于统计的准则，以便识别受损轴承。统计阈值只有在信号中的所有确定性分量已被去除的情况下才有效。不幸的是，在各种工业应用中，由异质振动源所特征化，同步平均的第一步并不足以完全消除确定性分量，需要在包络分析之前进行额外的预白化步骤。过去提出了不同的技术来实现这个目的：最广泛应用的是线性预测滤波器和谱峭度方法。最近，提出了一种基于倒谱分析的新技术用于预白化：所谓的倒谱预白化。由于其低计算需求和简单性，它似乎是在自动损伤识别算法中执行中间预白化步骤的一个不错的选择。在本文中，将对在全尺寸工业轴承测试台上测得的数据进行新技术的有效性测试，该测试台能够重现严酷的运行条件。将与传统的预白化技术进行基准比较，作为验证倒谱预白化潜力的最后一步。

📚2 运行结果

部分代码：

function [xSES,alpha,th] = SES(x,fs,bpf,plotFlag,p,cpswFlag)
%% Estimation of the Squared Envelope Spectrum
%     this function can be used for detecting bearing faults under constant
%     working speed
%
% INPUTS
%    x = input signal
%    fs = sampling frequency
%    bpf = band-pass filter frequencies, use a vector as [f lower, f higher]
%       put and empty vector if band-pass filtering is not needed
%        bearing fault detection can be improved if performed in a frequency band
%        wher the SNR is high (typically about a system resonance)
%    plotFlag = display the SES, 0 -> no (default), 1 -> yes 
%    p = threshold significance level, default p = .999 (99.9%)
%    cpswFlag = cesptrum pre-whitening, 0 -> no (default), 1 -> yes 
%        bearing fault detection is affected by periodic contribution due to
%        external sources such as gears. This effect can be reduced by whitening
%        the signal before SES
%
% OUTPUTS
%    SES = squared envelope spectrum
%    alpha = cyclic frequencies
%    th = threshold
%
 

if nargin < 4
  plotFlag = 0;
  p = .999;
  cpswFlag = 0;
end
if nargin < 5
  p = .999;
  cpswFlag = 0;
end
if nargin < 6
  cpswFlag = 0;
end

L = length(x);
k = (0:L-1);

% cepstrum pre-whitening
if cpswFlag == 1;
  x = real(ifft(fft(x) ./ abs(fft(x))));
end

% band-pass filtering and ses estimation
if isempty(bpf)
  l = 1
  h = floor(L/2)+1;
  wfilt = zeros(size(x)); wfilt(l:h) = 1;
  xf = ifft(2 .* fft(x) .* wfilt); % filtered analytic signal  
else
  l = floor(bpf(1)*L/fs); % lower freq. index
  h = floor(bpf(2)*L/fs); % higher freq. index
  wfilt = zeros(size(x)); wfilt(l:h) = 1;
  xf = ifft(2 .* fft(x) .* wfilt); % filtered analytic signal
end
ENV = abs(xf).^2; % squared envelope
xSES = abs(1/L .* fft( ENV )) .^ 2; % squared envelope spectrum

% threshold
S0 = (h - l - k) ./ (2 * (h - l)^2 ) .* (mean(abs(xf).^2)).^2;
th = chi2inv(p,2) .* S0;

% keep only meaningful cyclic frequencies
alpha = k .* fs ./ L; % cyclic frequencies vector
alpha = alpha(1:h - l);
xSES = xSES(1:h - l); xSES(1) = 0; % put to zero the DC-term of SES in order to 
th = th(1:h - l);                  % improve its visualization 
if plotFlag == 1
% display results
  tt = k ./ fs; % time vector
  figure
  subplot(211)
  plot(tt,ENV,'k')
  title('squared envelope')
  xlabel('time (s)')
  box off
  subplot(212)
  plot(alpha,xSES,'k')
  title('squared envelope spectrum')
  hold on, plot(alpha,th,'r')
  legend('SES',[num2str(p .* 100) '% threhsold'  ])
  xlabel('cyclic frequency (Hz)')
  box off
end

🎉3 参考文献

文章中一些内容引自网络，会注明出处或引用为参考文献，难免有未尽之处，如有不妥，请随时联系删除。

🌈4 Matlab代码及文章