将上一篇中的程序改成函数文件：

[ a f y] = afyimf( x )，其中x是IMF单分量信号输入，a、f、y分别是单分量瞬时振幅、瞬时频率、对应复信号虚部输出。现在看看imfMBrd各个分量的瞬时频率情况。

imf=imfMBrd;

sr=size(imf,1)-1;
sc=size(imf,2);

%A=zeros(sr,sc);
F=zeros(sr,sc);
%Y=zeros(sr,sc);

for k=1:sr
    x=imf(k,:);
    [ a f y] = afyimf( x );
    %A(i,:)=a;
    F(k,:)=f;
    %Y(i,:)=y;
end

zs=zeros(1,sc);

fg1=figure(1);
set(fg1,'position',[50,1,700,1050])

li=fix(sr/2);
lj=1;

for i=1:li
for j=1:lj
n=(i-1)*lj+j;
s(n)=subplot(li,lj,n);

set(s(n),'position',[0.06+(j-1)/lj,(li-0.94*(0.999^i)*i)/li,0.92/lj,0.65/li]);

plot(F(n,:))
hold on
plot(zs,'k--')
xlim([0 sc+100])
title(['imfMBrd__tf',int2str(n)])
end
end

    运行 得

    图32-1

fg2=figure(2);
set(fg2,'position',[50,1,700,1050])

li=ceil(sr/2);%
lj=1;

for i=1:li
for j=1:lj
n=(i-1)*lj+j;
s(n)=subplot(li,lj,n);

set(s(n),'position',[0.06+(j-1)/lj,(li-0.94*(0.999^i)*i)/li,0.92/lj,0.65/li]);

plot(F(n+fix(sr/2),:))
hold on
plot(zs,'k--')
xlim([0 sc+10])
title(['imfMBrd__tf',int2str(n+fix(sr/2))])
end
end

    运行 得：

    图32-2

    从上二图可以看到，最大的频率取值范围是-0.5～0.5Hz，这也是我的程序算法所决定的频率取值域。上面13个IMF分量，只有前两个分量出现了负频，其余都是正频。

    将前两个分量出现负频的时段放大看，如下：

fg3=figure(3);
set(fg3,'position',[50,250,1050,700])

li=2;%
lj=1;

for i=1:li
for j=1:lj
n=(i-1)*lj+j;
s(n)=subplot(li,lj,n);

set(s(n),'position',[0.03+(j-1)/lj,(li-0.94*(0.999^i)*i)/li,0.95/lj,0.82/li]);

plot(F(n,:))
hold on
plot(zs,'k--')
xlim([12000 12200])
title(['imfMBrd__tf',int2str(n)])
end
end
    
     运行 得

    图32-3

    可以看到，负频部分与正频部分的连接，是连续的，再没有以前我一直抱怨的“跳变”。可以认为，这种负频，是信号本身的属性，而不是人为地制造出来的。负频部分，其绝对值绝大部分小于0.25Hz，但有少量大于0.25Hz，这都取决于信号本身的情况。

    换一组IMF看看。将上述最前面的imf=imfMBrd改成imf=imf_TW486012(体温序列IMF，具见博文《13、EMD分解－－生理信号的Hilbert-Huang变换（一）》)，修改下绘图句柄，即得：

    图32-4

    图32-5

    图32-6
 

    可以看到，情况与imf=imfMBrd是一样的。特别说明一点：上面两组IMF，在前两个分量的瞬时频率中，很容易看到，距数据右端点将近20000点的地方，频率变化范围有一个突变。经我查实，在那个时间点，我的数据测量工作中发生最重大的一件事，是采样频率，由一个时辰测量一次，改成了两个时辰测量一次，中间点数据由样条插值补上。这除了说明，高频信号必须用高频采样，也说明，实际信号当中的所谓“高频”，其实是“宽频”的代名词。IMF序号靠前的分量，与其说是“高频分量”，还不如说是“宽频分量”。

    这一篇就到这里，下一篇续谈。

（本文首发于：http://blog.sina.com.cn/s/blog_6ad0d3de01016gw5.html

首发时间：2012-06-13 13:45:56）