信号处理——Hilbert端点效应浅析

作者:桂。

时间:2017-03-05  19:29:12

链接:http://www.cnblogs.com/xingshansi/p/6506405.html

声明:转载请注明出处,谢谢。



前言


本文为Hilbert变换分析的补充,主要介绍Hilbert变换中的端点效应,内容拟分两部分展开:

  1)Gibbs现象介绍;

  2)端点效应分析;

内容为自己一家之言,其中不合理的地方,希望各位告诉我,反正改不改是我的事(●‘?‘●)

一、Gibbs现象

关于Gibbs的理论推导,可以参考郑君里的《信号与系统》上册(第二版,P97~101),里边有详细的理论分析。也可参考:维基百科,里边有一张图很形象:

给出对应的代码(用该程序,借助MATLAB便可以生成.gif图):

%% 代码介绍========================================================================
%% Title:
%  1) Fourier series approximation of square wave.
%  2) Demonstration of Gibbs phenomenon (verification of Fig. 3.9 of [1])
%% Author:
%  Ankit A. Bhurane ([email protected])
%% Expression:
% The Fourier series approximation of a square wave signal existing between
% -Tau/2 to Tau/2 and period of T0 will have the form:
%
% Original signal to be approximated:
%                         :
%               __________:__________  A
%              |          :          |
%              |          :          |
%              |          :          |
%    __________|          :          |__________
% -T0/2     -Tau/2        0        Tau/2       T0/2
%                         :
%                         :
%
% Its Fourier series approximation:
%
%                  Inf
%                  ___
%         A*Tau   \     / sin(pi*n*Tau/T0)                    % r(t) = -------   |   | ------------------ exp^(j*n*2*pi/T0*t) |
%          T0     /___  \  (pi*n*Tau/T0)                      /
%               n = -Inf
%
% The left term inside summation are the Fourier series coeffs (Cn). The
% right term is the Fourier series kernel.
% Tau: range of square wave, T: period of the square wave,
% t: time variable, n: number of retained coefficients.
%
%% Observations:
% 1) As number of retained coefficients tends to infinity, the approximated
% signal value at the discontinuity converge to half the sum of values on
% either side.
% 2) Ripples does not decrease with increasing coefficients with
% approximately 9% overshoot.
% 3) Energy in the error between original and approximated signal, reduces
% as the number of retained coefficients are increased.
%
%% References:
% [1] Oppenheim, Willsky, Nawab, "Signals and Systems", PHI, Second edition
% [2] Dean K. Frederick and A. Bruce Carlson, "Fourier series" section in
%     Linear systems in communication and control
%% Last Modified: Sept 24, 2013.
%% Copyright (c) 2013-2014 | Ankit A. Bhurane
%%

  

clc; clear all; close all;

% Specification
A = 1; % Peak-to-peak amplitude of square wave
Tau = 10; % Total range in which the square wave is defined (here -5 to 5)
T0 = 20; % Period (time of repeatation of square wave), here 10
C = 200; % Coefficients (sinusoids) to retain
N = 1001; % Number of points to consider
t = linspace(-(T0-Tau),(T0-Tau),N); % Time axis
X = zeros(1,N); X(t>=-Tau/2 & t<=Tau/2) = A; % Original signal
R = 0; % Initialize the approximated signal
k = -C:C; % Fourier coefficient number axis
f = zeros(1,2*C+1); % Fourier coefficient values

% Loop for plotting approximated signals for different retained coeffs.
for c = 0:C % Number of retained coefficients
    for n = -c:c % Summation range (See equation above in comments)

        % Sinc part of the Fourier coefficients calculated separately
        if n~=0
            Sinc = (sin(pi*n*Tau/T0)/((pi*n*Tau/T0))); % At n NOTEQUAL to 0
        else
            Sinc = 1; % At n EQUAL to 0
        end
        Cn = (A*Tau/T0)*Sinc; % Actual Fourier series coefficients
        f(k==n) = Cn; % Put the Fourier coefficients at respective places
        R = R + Cn*exp(1j*n*2*pi/T0.*t); % Sum all the coefficients
    end

    R = real(R); % So as to get rid of 0.000000000i (imaginary) factor
    Max = max(R); Min = min(R); M = max(abs(Max),abs(Min)); % Maximum error
    Overshoot = ((M-A)/A)*100; % Overshoot calculation
    E = sum((X-R).^2); % Error energy calculation

    % Plots:
    % Plot the Fourier coefficients
    subplot 211; stem(k,f,‘m‘,‘LineWidth‘,1); axis tight; grid on;
    xlabel(‘Fourier coefficient index‘);ylabel(‘Magnitude‘);
    title(‘Fourier coefficients‘);

    % Plot the approximated signal
    subplot 212; plot(t,X,t,R,‘m‘,‘LineWidth‘,1); axis tight; grid on;
    xlabel(‘Time (t)‘); ylabel(‘Amplitude‘);
    title([‘Approximation for N = ‘, num2str(c),...
    ‘. Overshoot = ‘,num2str(Overshoot),‘%‘,‘. Error energy: ‘,num2str(E)])

    pause(0.1); % Pause for a while
    R = 0; % Reset the approximation to calculate new one
end

  简而言之:

  • 对于跳变的点,傅里叶变换的分量只能是能量收敛,而不是一致(幅度)收敛;
  • 对于跳变的点,如门函数/方波,信号由低频到高频众多分量组合而成。

故对于连续时域信号,Fourier级数只能无限项逼近,不能完全一致。

二、端点效应

  A-理论解释

借用之前博文的一张图:

连续信号到离散信号需要进行采样,由于有无穷多项,采样率对于高频的部分(因为无穷多项)难以满足Nyquist采样定理,因此会出现失真,这也就是端点效应的理论解释。失真的根本原因在时域采样,频域采样步骤只不过影响频域分辨率,即所谓的栅栏效应,但不是造成失真的根本原因。

  B-现象分析

为了验证理论,此处采样一个正弦信号进行分析,按连续信号来看其Hilbert变换应该是余弦信号。

给出代码:

clc;clear all;close all;
fs = 2000;
f0 = 80;
t = 0:1/fs:.1;
sig = sin(2*pi*f0*t);
sig_ref = -cos(2*pi*f0*t);
sig_hilbert = hilbert(sig);
figure;
subplot 311
plot(t,sig,‘k‘,‘linewidth‘,2);hold on;
plot(t,sig_ref,‘r.-‘,‘linewidth‘,2);hold on;
plot(t,imag(sig_hilbert),‘b‘,‘linewidth‘,2);hold on;
legend(‘原信号‘,‘余弦信号‘,‘Hilbert变换信号‘,‘localization‘,‘best‘);
ylim([-2,2]);
grid on;
%频谱
f = linspace(0,fs,length(t));
subplot 312
plot(f,abs(fft(sig)),‘k‘);hold on;
plot(f,abs(fft(sig_hilbert)),‘r‘);hold on;
legend(‘原信号‘,‘Hilbert变换信号‘,‘localization‘,‘best‘);
grid on;

  对应的结果图:

时域

可以看到明显的端点效应(说是端点效应,如果间断点在中间位置,一样会有该效应,毕竟都是Nyquist定理下的误差嘛)。

频域

采样率取得够大了,理论上正弦信号Hilbert之后应该是一个单边的冲击。但频谱显示,峰值旁边的一个点,幅度达到160.6,这也印证了能量的泄露。

发现端点效应,与前文理论对应:Hilbert就是一个频域的变相器,本质也是fourier变换,是不是也因为跳变点(即信号中的Gibbs问题)?观察频谱,峰值旁边的点,幅值有160.5大小,我们计算原信号

2*sig(1)-sig(2)

  得到结果是:-0.4820,而对于没有的点,应该是默认是0的,因此端点算是一个跳变点了
下面再做一组实验,验证我的推断
思路——只要保证2sig(1)-sig(2)逼近0,理论上是不是就没有边界效应了?
此时:2*sig(1)-sig(2)=1.5873e-05,再来看看结果图:

 再来看看频谱,没错,峰值旁边的点压的很小(6.183e-13),符合预期。

此处代码:

clc;clear all;close all;
fs = 20000;
f0 = 80;
t = 0.01255:1/fs:.1;%t = 0.0252:1/fs:.1则没有跳变。
sig = sin(2*pi*f0*t);
2*sig(1)-sig(2)
sig_ref = -cos(2*pi*f0*t);
sig_hilbert = hilbert(sig);
figure;
subplot 211
plot(t,sig,‘k‘,‘linewidth‘,2);hold on;
plot(t,sig_ref,‘r.-‘,‘linewidth‘,2);hold on;
plot(t,imag(sig_hilbert),‘b‘,‘linewidth‘,2);hold on;
legend(‘原信号‘,‘余弦信号‘,‘Hilbert变换信号‘,‘localization‘,‘best‘);
ylim([-2,2]);
grid on;
%频谱
f = linspace(0,fs,length(t));
subplot 212
plot(f,abs(fft(sig)),‘k‘);hold on;
plot(f,abs(fft(sig_hilbert)),‘r‘);hold on;
legend(‘原信号‘,‘Hilbert变换信号‘,‘localization‘,‘best‘);
ylim([0 6000]);
grid on;

基于这个基本要点,很多方法对信号:预测、延拓、镜像等操作,以消除或者降低端点效应。为什么用2*sig(1)-sig(2)呢,这是因为实验的函数为sin,在0附近sinx~x ,严格意义上应该是2*f(1)-f(2)≈ 0 。f是函数表达式,而且实际中信号不是单一频点,仅仅通过两个点判断也是不够的。
只是消除端点方法不同,但端点现象的本质一致。
多说一句:Hilbert将双边谱压成单边,这让频域滤波等操作方便了不少。

参考:

Gibbs phenomen:https://en.wikipedia.org/wiki/Gibbs_phenomenon

时间: 2024-09-29 18:02:15

信号处理——Hilbert端点效应浅析的相关文章

Linux——浅析信号处理

信号及其处理 信号处理是Unix和LInux系统为了响应某些状况而产生的事件,通常内核产生信号,进程收到信号后采取相应的动作. 例如当我们想强制结束一个程序的时候,我们通常会给它发送一个信号,然后该进程会捕捉到信号,紧接着该进程执行一定操作后最终被终止掉.不仅仅如此,通常下面几种情况 ①键盘事件(ctrl+c.ctrl+\) ②访问非法内存 ③硬件故障(如算术运算执行除以0操作) ④ 环境切换 都会有信号的产生,而对这些产生的信号是需要让进程来处理的,进而信号也被作为进程间通信或修改行为的一种方

SylixOS中pthread_cancel函数浅析

1 知识简介 1.1 概述 取消一个线程要确保该线程能够释放其所持有的任何锁.分配的内存,使整个系统保持一致性.在很多复杂情况下要保证这种正确性是有一定困难的. 一种简单的线程取消:取消线程调用一个取消线程的函数,被取消线程死亡.在这种情况下,被取消线程所持有的的资源得不到释放.取消线程负责保证被取消者处于可安全取消状态,在一个要求可靠性高的系统中,这种保证非常困难或者无法实现.这种取消称为不受限制的异步取消. 还存在另外一种更安全的线程取消机制.一个线程可以以可靠的受控制的方式向进程的其他线程

Matlab信号处理工具箱函数

波形产生和绘图chirp 产生扫描频率余弦diric 产生Dirichlet函数或周期Sinc函数gauspuls 产生高斯调制正弦脉冲pulstran 产生脉冲串rectpuls 产生非周期矩形信号sawtooth 产生锯齿波或三角波sinc 产生sinc函数square 产生方波strips 产生条图tripuls 产生非周期三角波 滤波器分析和实现abs 绝对值(幅值)angle 相位角conv 卷积和多项式乘法conv2 二维卷积fftfilt 基于FFT重叠加法的数据滤波filter

浅析 Linux 中的时间编程和实现原理一—— Linux 应用层的时间编程【转】

本文转载自:http://www.cnblogs.com/qingchen1984/p/7007631.html 本篇文章主要介绍了"浅析 Linux 中的时间编程和实现原理一—— Linux 应用层的时间编程",主要涉及到浅析 Linux 中的时间编程和实现原理一—— Linux 应用层的时间编程方面的内容,对于浅析 Linux 中的时间编程和实现原理一—— Linux 应用层的时间编程感兴趣的同学可以参考一下. 简介: 本文试图完整地描述 Linux 系统中 C 语言编程中的时间问

boost.asio源码剖析(二) ---- 架构浅析

* 架构浅析 先来看一下asio的0层的组件图.                     (图1.0) io_object是I/O对象的集合,其中包含大家所熟悉的socket.deadline_timer等对象,主要功能是提供接口给用户使用. services服务是逻辑功能的实现者,其中包含提供定时功能的deadline_timer_service.提供socket相关功能的win_iocp_socket_service(windows平台)/reactive_socket_service(其他

linux进程状态浅析

linux中的进程状态: ◆运行状态(TASK_RUNNING)(R状态) 指正在被CPU运行或者就绪的状态.这样的进程被成为runnning进程.运行态的进程可以分为3种情况:内核运行态.用户运行态.就绪态. 只有在该状态的进程才可能在CPU上运行.而同一时刻可能有多个进程处于可执行状态,这些进程的task_struct结构(进程控制块)被放入对应CPU的可执行队列中(一个进程最多只能出现在一个CPU的可执行队列中).进程调度器的任务就是从各个CPU的可执行队列中分别选择一个进程在该CPU上运

Python之encode与decode浅析

 Python之encode与decode浅析 在 python 源代码文件中,如果你有用到非ASCII字符,则需要在文件头部进行字符编码的声明,声明如下: # code: UTF-8 因为python 只检查 #.coding 和编码字符串,为了美观等原因可以如下写法: #-*-coding:utf-8-*- 常见编码介绍: GB2312编码:适用于汉字处理.汉字通信等系统之间的信息交换. GBK编码:是汉字编码标准之一,是在 GB2312-80 标准基础上的内码扩展规范,使用了双字节编码.

浅析PHP的开源产品二次开发的基本要求

浅析PHP的开源产品二次开发的基本要求 第一, 基本要求:HTML(必须要非常熟悉),PHP(能看懂代码,能写一些小系统,如:留言板,小型CMS),Mysql(至少会一种数据库),Javascript(能看懂,能改现成的一些代码),Div+Css(能进行界面的调整,明白CSS是怎么使用的) 第二, 熟悉开源产品的使用,比如 Dedecms,你要知道怎么登录,怎么新建栏目,怎么添加文章,模板标签的使用方法,模型的概念和使用方法等等一些功能 第三, 要熟悉这个开源产品的数据库结构,还要理解里面核心文

word-break|overflow-wrap|word-wrap——CSS英文断句浅析

---恢复内容开始--- word-break|overflow-wrap|word-wrap--CSS英文断句浅析 一 问题引入 今天在再次学习 overflow 属性的时候,查看效果时,看到如下结果,内容在 div 中国换行了,可是两个 P 元素的内容并没有换行,搜索一番没有找到系统的答案,截图到群里请教大神,才知道是英文断句的问题,但是还是不太明白.之前没有遇到这种情况,为了彻底搞清楚,英文断句,又开始学习英文断句到底是怎么回事. 二 换行 每种语言里都有换行,就中文而言,我们最小语言单位