MATLAB编程并绘制函数或者信号包络线

大家好，欢迎来到IT知识分享网。

MATLAB绘制包络线可以通过hilbert或者envelope实现。Hilbert变换是一个很有用的变换，用它来做包络分析更是一种有效的数据处理方法。

实例1

程序

clc;
clear all;
close all;
t = 0:0.01:6*pi;
y = sin(t).*sin(9.*t);
h = hilbert(y); % Hilbert变换
y1 = 3*cos(t).*sin(9.*t);
h1 = hilbert(y1); % Hilbert变换
figure;
subplot(2,1,1);
plot(t,y);
hold on
plot(t,abs(h),'r--','linewidth',1.5)
hold on
plot(t,-1*abs(h),'r--','linewidth',1.5)
grid on
subplot(2,1,2);
plot(t,y1);
hold on
plot(t,abs(h1),'g-.','linewidth',1.5)
hold on
plot(t,-1*abs(h1),'g-.','linewidth',1.5)
grid on

运行结果

实例2

程序

% Hilbert变换测试
clc
clear all
close all

ts = 0.001;
fs = 1/ts;
N = 200;
f = 50;
k = 0:N-1;
t = k*ts;

% 信号变换
% 结论：sin信号Hilbert变换后为cos信号
y = sin(2*pi*f*t);
yh = hilbert(y);    % matlab函数得到信号是合成的复信号
yi = imag(yh);      % 虚部为书上定义的Hilbert变换

figure
subplot(211)
plot(t, y)
title('原始sin信号')
subplot(212)
plot(t, yi)
title('Hilbert变换信号')

% 检验两次Hilbert变换的结果（理论上为原信号的负值）
% 结论：两次Hilbert变换的结果为原信号的负值
yih = hilbert(yi);
yii = imag(yih);
max(y + yii)

% 信号与其Hilbert变换的正交性
% 结论：Hilbert变换后的信号与原信号正交
sum(y.*yi)

% 谱分析
% 结论：Hilbert变换后合成的复信号的谱没有大于奈氏频率的频谱，即其谱为单边的
NFFT = 2^nextpow2(N);
f = fs*linspace(0,1,NFFT);
Y = fft(y, NFFT)/N;
YH = fft(yh, NFFT)/N;

figure
subplot(211)
plot(f,abs(Y))
title('原信号的双边谱')
xlabel('频率f (Hz)')
ylabel('|Y(f)|')
subplot(212)
plot(f,abs(YH))
title('信号Hilbert变换后组成的复信号的双边谱')
xlabel('频率f (Hz)')
ylabel('|YH(f)|')

运行结果

实例3

程序

% 包络分析(高中心频率的窄带信号分析)
% 基于：两个信号乘积的Hilbert变换取决于高频信号的Hilbert变换
clc
clear all
close all

ts = 0.001;
fs = 1/ts;
N = 200;
k = 0:N-1;
t = k*ts;

% 原始信号
f1 = 10;
f2 = 70;
% a = cos(2*pi*f1*t);       % 包络1
a = 2 + cos(2*pi*f1*t);     % 包络2
% a = 1./(1+t.^2*50);       % 包络3
m = sin(2*pi*f2*t);         % 调制信号
y = a.*m;  % 信号调制

figure
subplot(241)
plot(t, a)
title('包络')
subplot(242)
plot(t, m)
title('调制信号')
subplot(243)
plot(t, y)
title('调制结果')

% 包络分析
% 结论：Hilbert变换可以有效提取包络、高频调制信号的频率等
yh = hilbert(y);
aabs = abs(yh);                 % 包络的绝对值
aangle = unwrap(angle(yh));     % 包络的相位
af = diff(aangle)/2/pi;         % 包络的瞬时频率，差分代替微分计算

% NFFT = 2^nextpow2(N);
NFFT = 2^nextpow2(1024*4);      % 改善栅栏效应
f = fs*linspace(0,1,NFFT);

YH = fft(yh, NFFT)/N;           % Hilbert变换复信号的频谱
A = fft(aabs, NFFT)/N;          % 包络的频谱

subplot(245)
plot(t, aabs, t, a, '.')
title('包络的绝对值')
legend('包络分析结果', '真实包络')
subplot(246)
plot(t, aangle)
title('调制信号的相位')
subplot(247)
plot(t(1:end-1), af*fs)
title('调制信号的瞬时频率')
subplot(244)
plot(f,abs(YH))
title('原始信号的Hilbert谱')
xlabel('频率f (Hz)')
ylabel('|YH(f)|')
subplot(248)
plot(f,abs(A))
title('包络的频谱')
xlabel('频率f (Hz)')
ylabel('|A(f)|')

运行结果

envelope绘制包络线

程序

clc;
clear all;
close all;
t = 0:0.01:6*pi;
y = sin(t).*sin(9.*t);
figure;
plot(t,y);
hold on
% 包络线绘制
[up,lo] = envelope(y);
plot(t,up,t,lo,'linewidth',1.5)
legend('signal','上包络线','下包络线')
hold off

运行结果

本文部分内容来源于长笛人倚楼Gloria的博客，仅供参考学习，如内容、图片有任何版权问题，请联系处理，24小时内删除。

作 者 | 郭志龙
编 辑 | 郭志龙
校 对 | 郭志龙