瞬时频率的初步认识

shuihai707

        说到瞬时频率，这是搞时频分析的同学不可绕过的话题。我们的教科书上往往讲的是对实信号x(t)，求出其Hilbert变换，构成解析信号，对解析信号的相位求导，便得到瞬时频率。这些表达让我们只会机械的使用Hilbert变换构造解析信号，而不知所以然。
        换一种思路想这个问题，调制信号是我们日常生活中经常使用的信号，它可以传递大量信息。对于一实信号x(t)，它能否表示成调制信号的形式呢，即x(t)=a(t)cos[b(t)]，并且a(t)和b(t)可以唯一的确定x(t)，a(t)能表示信号x(t)的瞬时幅值，b(t)能表示信号x(t)的瞬时相位，且瞬时相位求导能得到有意义的瞬时频率呢？答案是：可以的。我们接下来从解析信号入手，慢慢剥开这层迷雾。
        怎样得到这样的一对[a(t),b(t)]，解析信号法给我们提供了一个途径。H(x(t))=x(t)*(1/t乘以pai)=-isgn(f)X(f)（时域卷积等于频域相乘），则解析信号的频谱Z(f)=X(f)+i[-isgn(f)X(f)]=2X(f),对于f>0;Z(f)=0,对于f<0。这也给出一个复信号是否是解析信号的一个充分条件：如果一个复信号的频谱当f<0时，Z(f)=0，那么它就是一个解析信号。接下来的一个问题是：是否任意一个实信号都可以通过Hilbert变换构造成一个解析信号呢？显然不是，前提条件是实信号必须是一个单分量信号（这里的单分量信号可以简单描述成一个信号在任意一个时刻点只有一个频率值的信号，这不是一个很严格的概念，暂且这样理解）。我们构造解析信号的目的是找到这样一对[a(t),b(t)]，a(t)能表示信号x(t)的瞬时幅值，b(t)能表示信号x(t)的瞬时相位。如果实信号满足单分量的条件，通过Hilbert变换，我们可以得到瞬时幅值a(t)=根号下x²(t)+H²(x(t))，瞬时相位b(t)=arctanH(x(t))/x(t)，这样一对[a(t),b(t)]就可以唯一确定x(t)了，瞬时频率就可以通过瞬时相位求导得到了。事实真的是这样吗？不是，通过这个方法求出的瞬时频率有时会为负值，相关的例子很多。
        这是怎么回事呢？Bedrosian定理出场了，它指出通过解析信号法求得的瞬时频率如果有物理意义的话，瞬时幅值a(t)的频谱和cos[b(t)]的频谱不能重叠，并且a(t)的频谱要位于低频端，cos[b(t)]的频谱要位于高频端，一个著名的表达式就是H[a(t)cos(b(t))]=a(t)H[cos(b(t))]，原来是这样的，解析信号法不是万能的，构造完还需要通过傅里叶变换检验下构造的瞬时幅值和相位的频谱是否满足Bedrosian定理，如果满足了，才能求出真实的瞬时频率啊。事实真的会这样吗？继续看看。
        我们知道一个余弦函数cos(t)它的Hilbert变换是sin(t)，再回过头看看Bedrosian等式右边的式子，对于H[cos(b(t))]=sin(b(t))这个式子也会成立吗？答案是否定的，对于任意的b(t)是不会成立的，这就是著名的Nuttall定理（如果信号的相位是线性的，这个式子是成立的）。这个定理说明什么呢？说明我们通过Hilbert变换构造的解析信号的虚部与原信号的正交信号不相等，不相等怎么了？！不相等就要推翻Hilbert变换构造解析信号求取瞬时频率！！！你解析信号的虚部与原信号的正交信号不相同，或者说解析信号的虚部与原信号不正交，你能用b(t)=arctanH(x(t))/x(t)求瞬时相位吗？？说了一大顿原来解析信号法求瞬时频率不靠谱啊，从Bedrosian定理Nuttall定理的角度看确实是这样。
        但是，凡事都有个但是，我们能不能找到这样的一对[a(t),b(t)]，满足Bedrosian定理和Nuttall定理，并且瞬时频率有意义呢？答案是肯定的，有并且有比较多，这里就要涉及到单分量信号的概念了，我国澳门大学数学教授钱涛提出了单分量信号的数学上的概念，单分量信号是指满足H[a(t)cos(b(t))]=a(t)sin(b(t))【这个式子将Bedrosian定理和Nuttall定理两个式子合二为一了】,并且b(t)的导数大于等于0的一类信号。他提出了自适应傅里叶分解算法，就是将信号分解成这样一类单分量信号，有扎实的理论基础，相关的文献自己找找吧，里面涉及到复变函数相关知识，需要有一定的数学功底，我也说不清楚，只看个皮毛。
        最后再说说imf，imf只是调幅调频信号，采用Hilbert变换求瞬时频率肯定会受到Bedrosian定理和Nuttall定理的限制了，出现负频率也不足为奇了。所以黄提出了直接正交法、跨零点法解调瞬时频率，其余一些方法如能量算子法，论文都写烂了，并且效果很差（只对线性相位的效果好些），不知为什么有些论文说的那么好。还有一些求取瞬时频率的方法，我国王金良博士提出的直接插值法，还有alfa-counting瞬时频率法，都是最近1年提出的新方法，效果怎样大家可以试试。
        以上就是我对瞬时频率一些粗浅认识，欢迎大家一起讨论！

shuihai707

A NOVEL SIGNAL DECOMPOSITION APPROACH — ADAPTIVE FOURIER DECOMPOSITION

ESMD_for_Nonlinear_and_Non-Stationary_Signal_Processing
网上都有

shuihai707

我们知道一个余弦函数cos(t)它的Hilbert变换是sin(t)，再回过头看看Bedrosian等式右边的式子，对于H[cos(b(t))]=sin(b(t))这个式子也会成立吗？答案是否定的，对于任意的b(t)是不会成立的，这就是著名的Nuttall定理（如果信号的相位是线性的，这个式子是成立的）。这个定理说明什么呢？说明我们通过Hilbert变换构造的解析信号的虚部与原信号的正交信号不相等，不相等怎么了？！不相等就要推翻Hilbert变换构造解析信号求取瞬时频率！！！你解析信号的虚部与原信号的正交信号不相同，或者说解析信号的虚部与原信号不正交，你能用b(t)=arctanH(x(t))/x(t)求瞬时相位吗？？说了一大顿原来解析信号法求瞬时频率不靠谱啊，从Bedrosian定理Nuttall定理的角度看确实是这样。
关于上面这段，有些说的不严密。我做如下修改：
我们知道一个余弦函数cos(t)它的Hilbert变换是sin(t)，再回过头看看Bedrosian等式右边的式子，对于H[cos(b(t))]=sin(b(t))这个式子也会成立吗？答案是否定的，对于任意的b(t)是不会成立的，这就是著名的Nuttall定理（如果信号的相位是线性的，这个式子是成立的）。这个定理说明什么呢？先引入正交信号的概念，正交信号是一个复信号（这篇文章《Quadrature Signals-Complex But Not Complicated》对正交信号说明的很详细），我们通过解析信号法将信号x(t)表示成x(t)=a(t)cos(b(t))的形式，它的正交信号可以表示成x_q(t)=a(t)[cos(b(t))+isin(b(t))]，如果对于任意的b(t),Nuttall定理不成立，就意味着x(t)的解析信号和正交信号不相等，即a(t){cos(b(t))+iH[cos(b(t))]}不等于a(t)[cos(b(t))+isin(b(t))]。既然这两个表达式不相等，那么就存在能量误差，这个误差Nuttall也给出了度量。这样通过瞬时相位b(t)=arctanH(x(t))/x(t)求瞬时频率就会存在偏差，所以从Bedrosian定理和Nuttall定理的角度看解析信号法求得的瞬时频率是存在问题的。

yghit08

shuihai707 发表于 2013-11-20 16:28

                               
登录/注册后可看大图

我们知道一个余弦函数cos(t)它的Hilbert变换是sin(t)，再回过头看看Bedrosian等式右边的式子，对于H[cos(b( ...

不错，从理论上给出解释。我在看Hilbert变换的时候，里面提到这样的要求，比如信号x=a(t)*y(t),如果a(t)是比y(t)变化慢的函数的话，那么H（x）=a(t)*H(y).H(y)表示y的Hilbert变换。
求取解析信号方便将幅值和相位角分开，进而很容易求解瞬时频率。
目前应该有在不需要求取分析信号得到瞬时频率的研究工作在进行。

shuihai707

yghit08 发表于 2013-11-29 15:23

                               
登录/注册后可看大图

不错，从理论上给出解释。我在看Hilbert变换的时候，里面提到这样的要求，比如信号x=a(t)*y(t),如果a(t)是 ...

版主的那个等式就是Bedrosian定理的意思。信号变化快慢就是对应频率的高低，这样频谱就可以分离。直接插值法和alfa-counting瞬时频率法都不需要进行Hilbert变换求取信号瞬时频率，但都属于一种经验算法，不像Hilbert变换这样有理论基础，并且能和傅里叶变换相兼容。

lqlcandy

大师讲得很好，有一个问题，就是我用了台湾国立中央大学的matlab代码ifndq去求得瞬时频率，但是还会出现个别负值出现，问题是哪里呢？如果我的采样时间是1分钟，是不是由于我的采样频率太低造成的？希望大师指导，分享一下求得瞬时频率的代码。

shuihai707

lqlcandy 发表于 2015-1-23 00:07
大师讲得很好，有一个问题，就是我用了台湾国立中央大学的matlab代码ifndq去求得瞬时频率，但是还会出现个 ...

采样频率最好是信号最高频率4倍以上。计算imf出现瞬时频率负值现象，黄等也讨论过此问题,主要原因帖子里也说过，两个定理造成的。解决方法：黄提出的经验AM-FM分解加上直接正交法。

lqlcandy

shuihai707 发表于 2015-1-26 09:31
采样频率最好是信号最高频率4倍以上。计算imf出现瞬时频率负值现象，黄等也讨论过此问题,主要原因帖子里 ...

大师，我并没有直接求，我也知道那两个定理的局限性，所以我求得时候是先用他的代码hilbertnormaize 提取了AM，然后将已经正交化的data在用ifndq 程序求得瞬时频率，还是出现个别负值。说白了就是用的NHHT

shuihai707

lqlcandy 发表于 2015-1-26 13:14
大师，我并没有直接求，我也知道那两个定理的局限性，所以我求得时候是先用他的代码hilbertnormaize 提取 ...

好像ifndq这个函数本质就是用Hilbert变换，具体函数情况我不太清楚，看看ifndq函数是如何求瞬时频率的。

lqlcandy

shuihai707 发表于 2015-1-27 17:12
好像ifndq这个函数本质就是用Hilbert变换，具体函数情况我不太清楚，看看ifndq函数是如何求瞬时频率的。

我认真看了一下那个ifndq函数。
%good original value,direct calculate instantaneous frequency  
         omgcos(i)=abs(nvimf(i+1)-nvimf(i-1))*0.5/sqrt(1-nvimf(i)*nvimf(i)
这个就是所说的Quadrature就是直接正交化，应该来说是最接近真实瞬时频率的。
但是无论是我这个Quadrature还是使用NHHT方法，对于我的数据都有负的瞬时频率，我想方法都是不完善的，也可能是我的数据避免不了负值出现。但是不管怎么样，谢谢大师。
我有另外一个问题，就是他求得IMF后，求得希尔伯特边际谱，但是这个2维， 如果推广到任意维上，需要一个联合密度函数pdf，这样积分，我才能最小二乘拟合出来有关于柯西主值的函数，不知道大师指导怎么算这个pdf吗？或者有没有代码可以算任意维度。

奋斗奋斗

shuihai707 发表于 2013-11-20 16:28
我们知道一个余弦函数cos(t)它的Hilbert变换是sin(t)，再回过头看看Bedrosian等式右边的式子，对于H[cos(b( ...

mxlzhenzhu

楼主老师太帅了，我好喜欢，不不不，已经爱上你了。。。。

mxlzhenzhu

shuihai707 发表于 2013-11-20 16:28
我们知道一个余弦函数cos(t)它的Hilbert变换是sin(t)，再回过头看看Bedrosian等式右边的式子，对于H[cos(b( ...

赞同不宜直接上来就H变换求瞬时频率，但是文献里说DQ算法是最精确的，请问老师，现在如何看这个问题？