% 本数据参数：Recording-10  1500r/min，外圈，径向加载1kN,轴向加载2kN  所有通道可用
% 经过计算 外圈故障特征频率：理论值（5.6Fr=5.6*1500/60=140）
clc
clear all;

load Recording_10.mat

ch1=Channel_1_Data;
ch2=Channel_2_Data;

Fs = 65536;            % Sampling frequency 
T = 1/Fs;             % Sampling period       
L = length(ch1)           % Length of signal
t = (0:L-1)*T;        % Time vector
x = (0:L-1);        

figure(1)
%subplot(2,1,1)
%plot(t,ch1)   % t轴单位是s  
plot(x,ch1)
xlabel('Time/s','fontsize',10.5)
ylabel('Amp','fontsize',10.5)

%axis([0 0.1 -30 30])
axis([0 4000 -50 50])
%subplot(2,1,2)
%plot(ch2)

axis([0 4000 -50 50])

%================FFT===================
                  


Y = fft(ch1);   %此时的Y为复数 a+bi 
%计算双侧频谱 P2。然后基于 P2 和偶数信号长度 L 计算单侧频谱 P1。


P2 = abs(Y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
%定义频域 f 并绘制单侧幅值频谱 P1。一般情况下，较长的信号会产生更好的频率近似值。

f = Fs*(0:(L/2))/L;
figure(2)
plot(f,P1) 

title('x-axis vibration spectrogram','fontsize',10.5)

xlabel('Freq','fontsize',10.5)
ylabel('Amp','fontsize',10.5)

% set(gca, 'LineWidth',1)  
% set(gca,'FontSize',10.5);
% set(gcf,'Position',[100 100 1300 280]);

% axis([0 200 0 2])


%================Envelope===================
%  'analytic' 效果不佳
%  40,'rms   效果不错 （推荐）
%  peak 运算较慢   50,'peak 效果较好 特征峰值较高
[up,lo]=envelope(ch1,50,'rms');  
 
figure(3)
subplot(2,1,1)
plot(t,ch1)
axis([0 0.1 -30 30])

hold on;
plot(t,up,'linewidth',1.5)
subplot(2,1,2)
plot(t,up)
axis([0 0.1 -30 30])

%================envelope FFT===================

Y = fft(up);   %此时的Y为复数 a+bi 
%计算双侧频谱 P2。然后基于 P2 和偶数信号长度 L 计算单侧频谱 P1。


P2 = abs(Y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
%定义频域 f 并绘制单侧幅值频谱 P1。一般情况下，较长的信号会产生更好的频率近似值。

f = Fs*(0:(L/2))/L;
figure(4)
plot(f,P1,'r') 

title('x-axis vibration spectrogram','fontsize',10.5)

xlabel('Freq','fontsize',10.5)
ylabel('Amp','fontsize',10.5)

% set(gca, 'LineWidth',1)  
% set(gca,'FontSize',10.5);
% set(gcf,'Position',[100 100 1300 280]);

axis([0 200 0 6])