CN101916439A - 基于希尔伯特-黄变换的空间碎片高光谱序列探测方法 - Google Patents

基于希尔伯特-黄变换的空间碎片高光谱序列探测方法 Download PDF

Info

Publication number
CN101916439A
CN101916439A CN 201010236974 CN201010236974A CN101916439A CN 101916439 A CN101916439 A CN 101916439A CN 201010236974 CN201010236974 CN 201010236974 CN 201010236974 A CN201010236974 A CN 201010236974A CN 101916439 A CN101916439 A CN 101916439A
Authority
CN
China
Prior art keywords
curve
spectrum
imf component
amplitude
hyperspectral
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN 201010236974
Other languages
English (en)
Other versions
CN101916439B (zh
Inventor
沈毅
张淼
王强
王艳
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Harbin Institute of Technology
Original Assignee
Harbin Institute of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Harbin Institute of Technology filed Critical Harbin Institute of Technology
Priority to CN2010102369741A priority Critical patent/CN101916439B/zh
Publication of CN101916439A publication Critical patent/CN101916439A/zh
Application granted granted Critical
Publication of CN101916439B publication Critical patent/CN101916439B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Image Analysis (AREA)

Abstract

基于希尔伯特-黄变换的空间碎片高光谱序列探测方法，属于图像处理领域，本发明为解决采用传统的危险空间碎片探测方法时，必须依赖有关空间碎片图像的样本信息，算法适应性差的问题，本发明方法包括：一、连续T次采样疑似目标中心位置的高光谱曲线，并将获取的T段高光谱曲线处理合成为待处理曲线；二、一维经验模态分解，获取两个本征模态函数分量；三、对二阶IMF分量进行希尔伯特变换，得到其幅值和瞬时频率；四、保留幅值高于均值一半的部分，保留瞬时频率高于均值的部分，并分割形成特征波段集合；五、循环搜索特征波段集合，判断目标是否旋转，如旋转，则确定疑似目标为空间碎片，如不旋转，则确定疑似目标不是空间碎片。

Description

基于希尔伯特-黄变换的空间碎片高光谱序列探测方法
技术领域
本发明涉及基于希尔伯特-黄变换的空间碎片高光谱序列探测方法，属于图像处理领域。
背景技术
航天技术的日益发展有利地促进了世界各国的经济建设和社会进步，但与此同时，由于人类对太空开发的狂热，致使大量的人造物体被发射到太空之中。随着航天器发射数量的增加，失效的航天器，以及脱落、碰撞、爆炸等形成的空间碎片（space debris）的数量也在逐年增加。大、中型空间碎片主要包括废弃的卫星、火箭及其脱落物等；小型空间碎片主要包括卫星表面脱落物、发动机喷出的碎屑、碰撞或爆炸产生的碎片等等。
探测是掌握空间碎片分布情况的前提和基础，也是空间碎片观测研究中至关重要的工作。质量很大的大型碎片主要在地面进行观测，也较为容易识别，测量方法主要是可见光观测、雷达观测等。数量众多的中型空间碎片（尺寸在1mm-10cm之间），由于尺度较小，地基望远镜和雷达一般无法观测或测定其轨道，而航天器表面采样分析方法对其撞击后再进行分析又可能为时已晚，因此他们常被称为危险空间碎片。目前探测危险空间碎片的研究方法集中在天基可见光遥感探测，无论是探测方法还是识别方法都处于发展的初级阶段，且可见光相机较高光谱成像仪在视域信息提供量上存在较大差距。在传统的危险空间碎片（中型空间碎片）探测手段中，主要依靠天基可见光相机的成像序列并利用空间碎片的运动特性进行探测，必须依赖有关空间碎片图像的样本信息，对疑似目标只能获取其位置及灰度信息，天基可见光相机无法利用空间碎片围绕它的主轴进行自旋运动，算法适应性差。
发明内容
本发明目的是为了解决采用传统的危险空间碎片探测方法时，天基可见光相机无法利用空间碎片围绕它的主轴进行自旋运动，必须依赖有关空间碎片图像的样本信息，算法适应性差的问题，提供了基于希尔伯特-黄变换的空间碎片高光谱序列探测方法。
本发明包括如下步骤：
步骤一、连续T次采样疑似目标中心位置的高光谱曲线，并将获取的T段高光谱曲线处理合成为待处理曲线；
步骤二、对步骤一获取的待处理曲线数据进行一维经验模态分解，获取两个本征模态函数分量和残差，两个本征模态函数分量为一阶IMF分量和二阶IMF分量；
步骤三、对所述二阶IMF分量进行希尔伯特变换，得到二阶IMF分量的幅值和瞬时频率；
步骤四、保留同时满足幅值高于均值一半和瞬时频率高于均值两个条件的二阶IMF分量，并分割形成特征波段集合；
步骤五、循环搜索特征波段集合，判断目标是否旋转，进而确定疑似目标是否为空间碎片，如旋转，则确定疑似目标为空间碎片，如不旋转，则确定疑似目标不是空间碎片。
本发明的优点：
1）本发明利用高光谱成像仪进行天基采样，能够对目标因旋转所产生的光谱微弱变化进行检测，克服了天基可见光相机无法利用空间碎片围绕它的主轴进行自旋运动这一事实来进行探测，使得空间碎片的判定能够基于对象的特殊运动信息，而无需任何有关空间碎片图像的样本信息，因而算法适应性强。
2）本发明所提出的探测方法，在确定空间碎片的同时，不但可以通过高光谱图像信息定位对象的位置，而且可以估计对象的自旋周期，可为空间碎片的天基探测搜集更多的目标信息。
附图说明
图1为本发明方法流程图；
图2为一维经验模态分解流程图；
图3为实施例中观测像素的20段原始采样数据合成后的高光谱曲线；
图4为合成高光谱曲线的一阶本征模态函数（IMF1）曲线；
图5为合成高光谱曲线的二阶本征模态函数（IMF2）曲线；
图6为合成高光谱曲线的残差函数（RES）曲线；
图7为IMF2经希尔伯特变换后的幅值曲线；
图8为IMF2经希尔伯特变换后的瞬时频率曲线；
图9为本发明方法筛选后的特征波段定位图。
具体实施方式
具体实施方式一：下面结合图1和图2说明本实施方式，
1998年NASA的黄锷（Norden E Huang）提出了一种依据数据本身的时间尺度特征来对信号进行分解的新算法——希尔伯特-黄变换（Hilbert-Huang Transform，HHT）。其中经验模态分解法（Empirical Mode Decomposition，EMD）是HHT算法的核心步骤，利用信号内部时间尺度的变化做能量与频率的解析，将信号展开成若干个本征模态函数（Intrinsic Mode Function，IMF），再利用希尔伯特变换（Hilbert Transform，HT）获得IMF的瞬时频率及振幅，最后求得时间-频率-振幅的三维谱分布，即Hilbert谱。
其中，IMF必须满足下列条件：
1）在整个函数中，极值点的数目与穿越零点的数目相等或者相差1；
2）在任何时刻，由局部极值包络线所定义的包络线局部均值为零。
依托这两个条件构建起来的EMD及HHT被认为是对以傅立叶变换为基础的线性及稳态谱分析的重大突破，是近年来强有力地求解非线性、非平稳信号的自适应方法。
在传统的危险空间碎片（中型空间碎片）探测手段中，主要依靠天基可见光相机的成像序列并利用空间碎片的运动特性进行探测，对疑似目标只能获取其位置及灰度信息。而遥感探测器若改用高光谱成像仪，则可以增加捕获对象的光谱信息，进而可利用空间碎片通常会围绕它的主轴进行简单的自旋运动这一事实来进行探测。由于目标像素的光谱序列属于非线性、非平稳信号，因此采用HHT进行自旋特征提取进而检测出空间碎片是非常适合的。
本实施方式方法通过以下技术方案实现的：选定观察对象中心位置处目标像素的                                                
Figure 553860DEST_PATH_IMAGE001
段连续高光谱数据序列，合成1段数据后进行一维经验模态分解（Empirical Mode Decomposition，EMD），对得到的二阶本征模态函数（Intrinsic Mode Function，IMF）进行希尔伯特变换（Hilbert Transform，HT）获得其瞬时频率及振幅，进而综合判断特征位置并分割为
Figure 471000DEST_PATH_IMAGE001
段，各段反映的是每个采样时刻得到的高光谱数据的高频特征信息所在波段，通过循环搜索判断目标是否旋转从而检测出空间碎片。具体方法包括以下步骤：
步骤一、连续T次采样疑似目标中心位置的高光谱曲线，并将获取的T段高光谱曲线处理合成为待处理曲线；
步骤二、对步骤一获取的待处理曲线数据进行一维经验模态分解，获取两个本征模态函数分量和残差，两个本征模态函数分量为一阶IMF分量和二阶IMF分量；
步骤三、对所述二阶IMF分量进行希尔伯特变换，得到二阶IMF分量的幅值和瞬时频率；
步骤四、保留同时满足幅值高于均值一半和瞬时频率高于均值两个条件的二阶IMF分量，并分割形成特征波段集合；
步骤五、循环搜索特征波段集合，判断目标是否旋转，进而确定疑似目标是否为空间碎片，如旋转，则确定疑似目标为空间碎片，如不旋转，则确定疑似目标不是空间碎片。
高光谱成像仪每次采样会得到观测区域的若干幅连续光谱图像，选择其中疑似目标的中心位置像素，可以得到一条光谱曲线。连续
Figure 562322DEST_PATH_IMAGE001
次观测采样，则可得到该位置像素的
Figure 19848DEST_PATH_IMAGE001
段光谱曲线。如果疑似目标为空间碎片，大多数情况下都会发生自旋，表现在同一位置像素的光谱曲线上就是会有周期性的变化。下面在步骤一中详细说明。
步骤一中的连续T次采样疑似目标中心位置的高光谱曲线，并将获取的T段高光谱曲线处理合成为待处理曲线的过程为：
步骤11、采用高光谱成像仪采集观测区域的多幅连续高光谱图像，选择其中疑似目标的中心位置像素，获取一条疑似目标中心位置的高光谱曲线，所述疑似目标中心位置的高光谱曲线的长度为m，
步骤12、按照步骤11的方法连续T次采样疑似目标中心位置的高光谱曲线，获取具有连续数据的T段高光谱曲线，
步骤13、对每段高光谱曲线进行左右翻转，并连接在原始高光谱曲线之后，形成T段长度为2m的新的高光谱曲线，
步骤14、将所述T段长度为2m的新的高光谱曲线依次拼接，合成为待处理曲线。待处理曲线的长度为
Figure 834220DEST_PATH_IMAGE002
。
步骤二中的获取两个本征模态函数分量和残差的过程为：
设定输入的待处理曲线为
Figure 469732DEST_PATH_IMAGE003
，
Figure 799082DEST_PATH_IMAGE004
，
步骤21、IMF分解过程初始化：，且满足关系式
Figure 978446DEST_PATH_IMAGE006
成立，其中
Figure 34126DEST_PATH_IMAGE007
为第
Figure 601505DEST_PATH_IMAGE008
次分解后剩余的残差函数；
步骤22、筛选过程初始化，
Figure 338517DEST_PATH_IMAGE009
，且满足关系式
Figure 924219DEST_PATH_IMAGE010
成立，其中
Figure 354063DEST_PATH_IMAGE011
为第
Figure 907273DEST_PATH_IMAGE012
次本征模态函数分解中经过第
Figure 244714DEST_PATH_IMAGE013
次筛选后的剩余函数；
步骤23、根据筛选程序获取输入的待处理曲线
Figure 888185DEST_PATH_IMAGE003
经过第次本征模态函数分解的剩余的残差函数中经过第
Figure 827639DEST_PATH_IMAGE015
次筛选后的剩余函数
Figure 906453DEST_PATH_IMAGE016
；
步骤24、采用标准偏差准则判断输入的待处理曲线
Figure 716015DEST_PATH_IMAGE003
经过第
Figure 549979DEST_PATH_IMAGE014
次本征模态函数分解的剩余的残差函数中经过第
Figure 500618DEST_PATH_IMAGE015
次筛选后的剩余函数
Figure 196172DEST_PATH_IMAGE017
是否满足本征模态函数的条件，即
Figure 610973DEST_PATH_IMAGE018
是否小于阈值
Figure 553521DEST_PATH_IMAGE019
，
Figure 303040DEST_PATH_IMAGE020
；
判断结果为是，执行步骤25，判断结果为否，则
Figure 785974DEST_PATH_IMAGE021
，然后执行步骤23，
步骤25、提取一个本征模态函数分量
Figure 258544DEST_PATH_IMAGE022
；
步骤26、获取输入的待处理曲线
Figure 185043DEST_PATH_IMAGE003
经过第
Figure 172590DEST_PATH_IMAGE014
次本征模态函数分解的剩余的残差函数
Figure 396898DEST_PATH_IMAGE023
；
步骤27、判断是否满足下述关系式成立：
Figure 769980DEST_PATH_IMAGE024
；
判断结果为否，则
Figure 116647DEST_PATH_IMAGE025
，然后执行步骤22，判断结果为是，完成提取过程，获得两个本征模态函数分量：一阶IMF分量、二阶IMF分量
Figure 529174DEST_PATH_IMAGE026
；和
Figure 370222DEST_PATH_IMAGE027
个经过第2次本征模态函数分解的剩余的残差RES：
Figure 614122DEST_PATH_IMAGE028
。即
Figure 803795DEST_PATH_IMAGE029
。
步骤三所述的获取二阶IMF分量的幅值和瞬时频率的过程为：
步骤31、对二阶IMF分量
Figure 15202DEST_PATH_IMAGE030
进行离散褶积，获得其希尔伯特变换；
Figure 945298DEST_PATH_IMAGE032
，
步骤32、获取二阶IMF分量的幅值：
Figure 118921DEST_PATH_IMAGE030
的解析信号为
Figure 568357DEST_PATH_IMAGE033
，将所述解析信号的包络振幅
Figure 938159DEST_PATH_IMAGE034
作为二阶IMF分量的幅值，解析信号的包络振幅
Figure 140339DEST_PATH_IMAGE034
按如下公式计算：
Figure 999710DEST_PATH_IMAGE035
，
步骤33、获取二阶IMF分量的瞬时频率
Figure 608546DEST_PATH_IMAGE036
：先求取所述解析信号的相角
Figure 595088DEST_PATH_IMAGE037
：
Figure 668086DEST_PATH_IMAGE038
，
然后，根据下式获取二阶IMF分量的瞬时频率
Figure 636042DEST_PATH_IMAGE036
：
Figure 43759DEST_PATH_IMAGE039
。
合成高光谱曲线记录了同一目标位置像素的连续采样光谱序列，如果目标发生自旋，则该序列会存在循环的光谱。但是同样环境下采集的光谱曲线是具有高度相关性的，只在一部分高频信息所在波段存在差异，因此自适应地提取高频信息所在波段是非常有意义的，可以用这些波段的集合作为各采样时刻的特征位置，下面在步骤四中详细说明。
步骤四所述的特征波段集合的获取过程为：
步骤41、计算二阶IMF分量的幅值平均值和瞬时频率的平均值
Figure 682867DEST_PATH_IMAGE041
：
Figure 634774DEST_PATH_IMAGE042
，其中，
Figure 530235DEST_PATH_IMAGE044
，
步骤42、舍弃二阶IMF分量中幅值小于幅值平均值
Figure 827093DEST_PATH_IMAGE040
一半的波段，同时舍弃瞬时频率小于瞬时频率平均值
Figure 136851DEST_PATH_IMAGE041
，获取特征时刻集合
Figure 269893DEST_PATH_IMAGE045
：
Figure 136348DEST_PATH_IMAGE046
，
步骤43、在范围内依次截取长度为的集合，获得
Figure 976129DEST_PATH_IMAGE001
段采样数据各自对应的特征波段集合
Figure 519105DEST_PATH_IMAGE047
。
如果对象为空间碎片，在自旋一周之后，相同位置的观察像素会具有一定相关性，即便观察像素在碎片上的位置稍有偏移也是如此，只是相关性有所减弱。下面在步骤五进行详细说明。
步骤五所述的循环搜索特征波段集合，判断目标是否旋转，进而确定疑似目标是否为空间碎片的过程为：
步骤51、计算特征波段集合
Figure 123131DEST_PATH_IMAGE048
两两之间特征集合的相异波段数目
Figure 980228DEST_PATH_IMAGE049
，
步骤52、计算所述相异波段数目的均值
Figure 736832DEST_PATH_IMAGE050
；
步骤53、计算等间隔采样数据之间的相异波段数目采用下式计算均值
Figure 201442DEST_PATH_IMAGE051
：
Figure 246759DEST_PATH_IMAGE052
，
由于采样数据量为
Figure 704285DEST_PATH_IMAGE001
，循环周期的有效数值经验证应在
Figure 588100DEST_PATH_IMAGE053
到
Figure 676142DEST_PATH_IMAGE054
之间。
步骤54、获取中最小值所对应的间隔
Figure 17441DEST_PATH_IMAGE056
：
，
在到
Figure 807915DEST_PATH_IMAGE054
之间选取多个k值，计算出每个k值对应的均值
Figure 544926DEST_PATH_IMAGE051
，并将其中最小均值找出来，这个最小均值对应的间隔
Figure 130629DEST_PATH_IMAGE056
是用于衡量相关性的一个周期参数。
步骤55、判断
Figure 560473DEST_PATH_IMAGE058
是否成立，
判断结果为是，断定目标旋转导致观测像素在间隔
Figure 349568DEST_PATH_IMAGE059
的周期具有较高的光谱相关性，使得特征集合相异波段数目减少，进而判断疑似目标为空间碎片；判断结果为否，判断疑似目标不是空间碎片。
具体实施方式二、下面结图1至图9，结合某空间碎片模型的高光谱旋转成像序列给出一个具体实施例：
执行步骤一：选定观察对象中心位置处目标像素的多个连续高光谱数据序列进行数据合成，形成待处理曲线。
对连续采集的20个高光谱数据，调用目标中心位置像素的高光谱曲线（共计20段，每段长度m为220），接着对每段高光谱曲线进行左右翻转，并连接在原始光谱曲线之后，形成新的高光谱曲线（共计20段，每段长度为2m=440），最后将20段新的光谱曲线依次连接得到一条合成曲线为待处理曲线（如图3所示，长度
Figure 952588DEST_PATH_IMAGE002
=8800）作为后续步骤的一维输入光谱数据——待处理曲线。
执行步骤二：对合成的一维输入光谱数据进行一维经验模态分解，分解至一阶、二阶本征模态函数以及残差函数即停止分解，具体描述如下：
设定输入的待处理曲线为，
Figure 242810DEST_PATH_IMAGE060
，
步骤21、IMF分解过程初始化：
Figure 34048DEST_PATH_IMAGE005
，且满足关系式
Figure 112863DEST_PATH_IMAGE006
成立，其中
Figure 423890DEST_PATH_IMAGE007
为第
Figure 992274DEST_PATH_IMAGE008
次分解后剩余的残差函数；
步骤22、筛选过程初始化，
Figure 208492DEST_PATH_IMAGE009
，且满足关系式
Figure 402582DEST_PATH_IMAGE010
成立，其中
Figure 817383DEST_PATH_IMAGE011
为第
Figure 759931DEST_PATH_IMAGE012
次本征模态函数分解中经过第
Figure 10915DEST_PATH_IMAGE013
次筛选后的剩余函数；
步骤23、根据筛选程序获取输入的待处理曲线
Figure 493849DEST_PATH_IMAGE003
经过第次本征模态函数分解的剩余的残差函数中经过第
Figure 125873DEST_PATH_IMAGE015
次筛选后的剩余函数；
步骤24、采用标准偏差准则判断输入的待处理曲线
Figure 337729DEST_PATH_IMAGE003
经过第
Figure 477854DEST_PATH_IMAGE014
次本征模态函数分解的剩余的残差函数中经过第
Figure 824522DEST_PATH_IMAGE015
次筛选后的剩余函数
Figure 237048DEST_PATH_IMAGE017
是否满足本征模态函数的条件，即
Figure 311053DEST_PATH_IMAGE018
是否小于阈值
Figure 554952DEST_PATH_IMAGE019
，
Figure 10204DEST_PATH_IMAGE061
；
判断结果为是，执行步骤25，判断结果为否，则
Figure 723076DEST_PATH_IMAGE021
，然后执行步骤23，
步骤25、提取一个本征模态函数分量
Figure 351504DEST_PATH_IMAGE022
；
步骤26、获取输入的待处理曲线
Figure 387593DEST_PATH_IMAGE003
经过第次本征模态函数分解的剩余的残差函数；
步骤27、判断是否满足下述关系式成立：
Figure 878989DEST_PATH_IMAGE024
；
判断结果为否，则
Figure 848213DEST_PATH_IMAGE025
，然后执行步骤22，判断结果为是，完成提取过程，获得两个本征模态函数分量：一阶IMF分量（见图4）、二阶IMF分量（见图5）：
Figure 707585DEST_PATH_IMAGE026
；和
Figure 316421DEST_PATH_IMAGE027
个经过第2次本征模态函数分解的剩余的残差RES（见图6）：
Figure 801497DEST_PATH_IMAGE028
。即
Figure 874496DEST_PATH_IMAGE029
。
执行步骤三：对IMF2进行希尔伯特变换，得到二阶IMF分量的幅值和瞬时频率。
首先，按如下算式对二阶IMF分量：
Figure 655501DEST_PATH_IMAGE062
进行离散褶积得到其希尔伯特变换
Figure 751633DEST_PATH_IMAGE031
：
Figure 259975DEST_PATH_IMAGE063
；
接着，按如下算式计算的解析信号
Figure 841184DEST_PATH_IMAGE065
的包络振幅
Figure 486929DEST_PATH_IMAGE066
：
Figure 549694DEST_PATH_IMAGE067
；
然后，按如下算式计算解析信号的相角
Figure 597284DEST_PATH_IMAGE068
和瞬时频率
Figure 907043DEST_PATH_IMAGE069
：
Figure 289351DEST_PATH_IMAGE070
，
Figure 405075DEST_PATH_IMAGE071
；
该步骤完毕，得到二阶IMF分量的幅值和瞬时频率
Figure 538564DEST_PATH_IMAGE072
，分别如图7和图8所示。
执行步骤四：综合判断高光谱数据的高频信息所在波段，得到特征波段集合。
对于二阶IMF分量的瞬时频率中小于其平均值（
Figure 893322DEST_PATH_IMAGE073
）的波段予以排除以保留具有特征意义的高频信息所在波段；同时由于二阶IMF分量的幅值和瞬时频率是依波段对应的，而过小的幅值其相应波段以噪声成分居多，因此对于幅值小于其平均值（
Figure 750420DEST_PATH_IMAGE074
）二分之一的波段也要予以排除，最后得到特征时刻集合
Figure 756291DEST_PATH_IMAGE045
（见图9）：
Figure 470169DEST_PATH_IMAGE075
，
然后对
Figure 515485DEST_PATH_IMAGE076
依次截取长度为440的集合，即得到20段采样数据各自对应的特征波段集合。
执行步骤五：循环搜索判断目标是否旋转，进而确定是否为空间碎片。
首先，计算
Figure 334854DEST_PATH_IMAGE077
两两之间特征集合的相异波段数目；
接着，计算相异波段数目的均值
Figure 1513DEST_PATH_IMAGE079
；
由于采样数据量为
Figure 997151DEST_PATH_IMAGE080
，循环周期的有效数值经验证应在到之间，对等间隔采样数据之间的相异波段数目采用下式计算均值：
Figure 289089DEST_PATH_IMAGE083
；
找出
Figure 291680DEST_PATH_IMAGE084
（见表1）中最小值所对应的间隔
Figure 126650DEST_PATH_IMAGE085
；
由于
Figure 290915DEST_PATH_IMAGE086
大于阈值2，因此可判断目标旋转导致观测像素在间隔10的周期具有较高的光谱相关性，使得特征集合相异波段数目减少，进而判断目标为空间碎片。
表1  等间隔采样数据之间的相异波段数目均值
间隔k 相异波段数目均值Mk
2 46.23
3 47.86
4 50.05
5 36.17
6 48.17
7 51.95
8 48.56
9 48.69
10 16.30
11 52.33
12 55.13
13 49.43
14 56.83
15 44.20
16 39.50
17 57.33
18 52.00

Claims (7)

1.基于希尔伯特-黄变换的空间碎片高光谱序列探测方法，其特征在于，它包括如下步骤：
步骤一、连续T次采样疑似目标中心位置的高光谱曲线，并将获取的T段高光谱曲线处理合成为待处理曲线；
步骤二、对步骤一获取的待处理曲线数据进行一维经验模态分解，获取两个本征模态函数分量和残差，两个本征模态函数分量为一阶IMF分量和二阶IMF分量；
步骤三、对所述二阶IMF分量进行希尔伯特变换，得到二阶IMF分量的幅值和瞬时频率；
步骤四、保留同时满足幅值高于均值一半和瞬时频率高于均值两个条件的二阶IMF分量，并分割形成特征波段集合；
步骤五、循环搜索特征波段集合，判断目标是否旋转，进而确定疑似目标是否为空间碎片，如旋转，则确定疑似目标为空间碎片，如不旋转，则确定疑似目标不是空间碎片。
2.根据权利要求1所述的基于希尔伯特-黄变换的空间碎片高光谱序列探测方法，其特征在于，步骤一中的连续T次采样疑似目标中心位置的高光谱曲线，并将获取的T段高光谱曲线处理合成为待处理曲线的过程为：
步骤11、采用高光谱成像仪采集观测区域的多幅连续高光谱图像，选择其中疑似目标的中心位置像素，获取一条疑似目标中心位置的高光谱曲线，所述疑似目标中心位置的高光谱曲线的长度为m，
步骤12、按照步骤11的方法连续T次采样疑似目标中心位置的高光谱曲线，获取具有连续数据的T段高光谱曲线，
步骤13、对每段高光谱曲线进行左右翻转，并连接在原始高光谱曲线之后，形成T段长度为2m的新的高光谱曲线，
步骤14、将所述T段长度为2m的新的高光谱曲线依次拼接，合成为待处理曲线。
3.根据权利要求1所述的基于希尔伯特-黄变换的空间碎片高光谱序列探测方法，其特征在于，步骤二中的获取两个本征模态函数分量和残差的过程为：
设定输入的待处理曲线为                                                
Figure DEST_PATH_IMAGE001
，
Figure 2010102369741100001DEST_PATH_IMAGE002
，
步骤21、IMF分解过程初始化：
Figure 2010102369741100001DEST_PATH_IMAGE003
，且满足关系式
Figure 2010102369741100001DEST_PATH_IMAGE004
成立，其中
Figure 2010102369741100001DEST_PATH_IMAGE005
为第
Figure 2010102369741100001DEST_PATH_IMAGE006
次分解后剩余的残差函数；
步骤22、筛选过程初始化，
Figure 2010102369741100001DEST_PATH_IMAGE007
，且满足关系式
Figure 2010102369741100001DEST_PATH_IMAGE008
成立，其中
Figure 697460DEST_PATH_IMAGE009
为第
Figure DEST_PATH_IMAGE010
次本征模态函数分解中经过第
Figure 2010102369741100001DEST_PATH_IMAGE011
次筛选后的剩余函数；
步骤23、根据筛选程序获取输入的待处理曲线
Figure 2010102369741100001DEST_PATH_IMAGE012
经过第
Figure DEST_PATH_IMAGE013
次本征模态函数分解的剩余的残差函数中经过第
Figure DEST_PATH_IMAGE014
次筛选后的剩余函数
Figure DEST_PATH_IMAGE015
；
步骤24、采用标准偏差准则判断输入的待处理曲线
Figure 175453DEST_PATH_IMAGE012
经过第
Figure 17507DEST_PATH_IMAGE013
次本征模态函数分解的剩余的残差函数中经过第
Figure 412716DEST_PATH_IMAGE014
次筛选后的剩余函数是否满足本征模态函数的条件，即是否小于阈值
Figure DEST_PATH_IMAGE018
，
Figure DEST_PATH_IMAGE019
；
判断结果为是，执行步骤25，判断结果为否，则
Figure DEST_PATH_IMAGE020
，然后执行步骤23，
步骤25、提取一个本征模态函数分量
Figure DEST_PATH_IMAGE021
；
步骤26、获取输入的待处理曲线经过第
Figure 111868DEST_PATH_IMAGE013
次本征模态函数分解的剩余的残差函数
Figure DEST_PATH_IMAGE022
；
步骤27、判断是否满足下述关系式成立：
Figure DEST_PATH_IMAGE023
；
判断结果为否，则
Figure DEST_PATH_IMAGE024
，然后执行步骤22，判断结果为是，完成提取过程，获得两个本征模态函数分量：一阶IMF分量、二阶IMF分量
Figure DEST_PATH_IMAGE025
；和
Figure DEST_PATH_IMAGE026
个经过第2次本征模态函数分解的剩余的残差RES：
Figure DEST_PATH_IMAGE027
。
4.根据权利要求3所述的基于希尔伯特-黄变换的空间碎片高光谱序列探测方法，其特征在于，步骤24的中
Figure 628169DEST_PATH_IMAGE018
=0.25。
5.根据权利要求1所述的基于希尔伯特-黄变换的空间碎片高光谱序列探测方法，其特征在于，步骤三所述的获取二阶IMF分量的幅值和瞬时频率的过程为：
步骤31、对二阶IMF分量
Figure DEST_PATH_IMAGE028
进行离散褶积，获得其希尔伯特变换
Figure DEST_PATH_IMAGE029
；
Figure DEST_PATH_IMAGE030
，
步骤32、获取二阶IMF分量的幅值：
Figure 263287DEST_PATH_IMAGE028
的解析信号为
Figure DEST_PATH_IMAGE031
，将所述解析信号的包络振幅
Figure DEST_PATH_IMAGE032
作为二阶IMF分量的幅值，解析信号的包络振幅
Figure 807532DEST_PATH_IMAGE032
按如下公式计算：
Figure DEST_PATH_IMAGE033
，
步骤33、获取二阶IMF分量的瞬时频率
Figure DEST_PATH_IMAGE034
：先求取所述解析信号的相角
Figure DEST_PATH_IMAGE035
：
Figure DEST_PATH_IMAGE036
，
然后，根据下式获取二阶IMF分量的瞬时频率
Figure 987846DEST_PATH_IMAGE034
：
Figure DEST_PATH_IMAGE037
。
6.根据权利要求1所述的基于希尔伯特-黄变换的空间碎片高光谱序列探测方法，其特征在于，步骤四所述的特征波段集合的获取过程为：
步骤41、计算二阶IMF分量的幅值平均值
Figure DEST_PATH_IMAGE038
和瞬时频率的平均值
Figure DEST_PATH_IMAGE039
：
Figure DEST_PATH_IMAGE040
，其中
Figure DEST_PATH_IMAGE041
，
Figure DEST_PATH_IMAGE042
，
步骤42、舍弃二阶IMF分量中幅值小于幅值平均值
Figure 991443DEST_PATH_IMAGE038
一半的波段，同时舍弃瞬时频率小于瞬时频率平均值
Figure 790772DEST_PATH_IMAGE039
，获取特征时刻集合
Figure DEST_PATH_IMAGE043
：
Figure DEST_PATH_IMAGE044
，
步骤43、在范围内依次截取长度为的集合，获得段采样数据各自对应的特征波段集合
Figure DEST_PATH_IMAGE046
。
7.根据权利要求1所述的基于希尔伯特-黄变换的空间碎片高光谱序列探测方法，其特征在于，步骤五所述的循环搜索特征波段集合，判断目标是否旋转，进而确定疑似目标是否为空间碎片的过程为：
步骤51、计算特征波段集合
Figure DEST_PATH_IMAGE047
两两之间特征集合的相异波段数目
Figure DEST_PATH_IMAGE048
，
步骤52、计算所述相异波段数目的均值
Figure DEST_PATH_IMAGE049
；
步骤53、计算等间隔采样数据之间的相异波段数目采用下式计算均值
Figure DEST_PATH_IMAGE050
：
Figure DEST_PATH_IMAGE051
，
步骤54、获取
Figure DEST_PATH_IMAGE052
中最小值所对应的间隔
Figure DEST_PATH_IMAGE053
：
Figure DEST_PATH_IMAGE054
，
步骤55、判断
Figure DEST_PATH_IMAGE055
是否成立，
判断结果为是，断定目标旋转导致观测像素在间隔的周期具有较高的光谱相关性，使得特征集合相异波段数目减少，进而判断疑似目标为空间碎片；判断结果为否，判断疑似目标不是空间碎片。
CN2010102369741A 2010-07-27 2010-07-27 基于希尔伯特-黄变换的空间碎片高光谱序列探测方法 Expired - Fee Related CN101916439B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN2010102369741A CN101916439B (zh) 2010-07-27 2010-07-27 基于希尔伯特-黄变换的空间碎片高光谱序列探测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN2010102369741A CN101916439B (zh) 2010-07-27 2010-07-27 基于希尔伯特-黄变换的空间碎片高光谱序列探测方法

Publications (2)

Publication Number Publication Date
CN101916439A true CN101916439A (zh) 2010-12-15
CN101916439B CN101916439B (zh) 2012-08-08

Family

ID=43323943

Family Applications (1)

Application Number Title Priority Date Filing Date
CN2010102369741A Expired - Fee Related CN101916439B (zh) 2010-07-27 2010-07-27 基于希尔伯特-黄变换的空间碎片高光谱序列探测方法

Country Status (1)

Country Link
CN (1) CN101916439B (zh)

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105571412A (zh) * 2015-12-11 2016-05-11 中国人民解放军63850部队 基于Hilbert变换的弹丸进动周期提取方法
CN106373149A (zh) * 2015-07-23 2017-02-01 中国人民解放军海军大连舰艇学院 基于辅助图像分量和Hilbert变换的图像分解方法
WO2019068650A1 (fr) 2017-10-02 2019-04-11 Ecole Polytechnique Radar laser pour la detection des debris spatiaux
CN111198385A (zh) * 2019-12-26 2020-05-26 北京旷视机器人技术有限公司 障碍物检测方法、装置、计算机设备和存储介质
CN114298106A (zh) * 2021-12-29 2022-04-08 中国铁道科学研究院集团有限公司铁道建筑研究所 路基碾压中特征波识别方法、碾压状态判别方法及其应用
CN116992293A (zh) * 2023-09-26 2023-11-03 北京豪迈生物工程股份有限公司 一种用于化学发光仪器的智能化数据处理方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2003228707A (ja) * 2001-10-19 2003-08-15 Stmicroelectronics Sa 同心情報の符号化
CN101685435A (zh) * 2008-09-26 2010-03-31 财团法人工业技术研究院 用于图像纹理分析的多维度经验模态分析方法
CN101763336A (zh) * 2010-01-18 2010-06-30 哈尔滨工业大学 一种具有激励噪声添加参数选取功能的集合经验模态分解方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2003228707A (ja) * 2001-10-19 2003-08-15 Stmicroelectronics Sa 同心情報の符号化
CN101685435A (zh) * 2008-09-26 2010-03-31 财团法人工业技术研究院 用于图像纹理分析的多维度经验模态分析方法
CN101763336A (zh) * 2010-01-18 2010-06-30 哈尔滨工业大学 一种具有激励噪声添加参数选取功能的集合经验模态分解方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
《红外与毫米波学报》 20081031 陈志刚 等 高光谱图像光谱域噪声去除的经验模态分解方法 第27卷, 第5期 2 *
《遥感学报》 20070131 王坚 等 基于经验模态分解的高分辨率影像融合 第11卷, 第1期 2 *

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106373149A (zh) * 2015-07-23 2017-02-01 中国人民解放军海军大连舰艇学院 基于辅助图像分量和Hilbert变换的图像分解方法
CN106373149B (zh) * 2015-07-23 2019-03-22 中国人民解放军海军大连舰艇学院 基于辅助图像分量和Hilbert变换的图像分解方法
CN105571412A (zh) * 2015-12-11 2016-05-11 中国人民解放军63850部队 基于Hilbert变换的弹丸进动周期提取方法
WO2019068650A1 (fr) 2017-10-02 2019-04-11 Ecole Polytechnique Radar laser pour la detection des debris spatiaux
CN111198385A (zh) * 2019-12-26 2020-05-26 北京旷视机器人技术有限公司 障碍物检测方法、装置、计算机设备和存储介质
CN114298106A (zh) * 2021-12-29 2022-04-08 中国铁道科学研究院集团有限公司铁道建筑研究所 路基碾压中特征波识别方法、碾压状态判别方法及其应用
CN114298106B (zh) * 2021-12-29 2022-07-22 中国铁道科学研究院集团有限公司铁道建筑研究所 路基碾压中特征波识别方法
CN116992293A (zh) * 2023-09-26 2023-11-03 北京豪迈生物工程股份有限公司 一种用于化学发光仪器的智能化数据处理方法
CN116992293B (zh) * 2023-09-26 2023-12-08 北京豪迈生物工程股份有限公司 一种用于化学发光仪器的智能化数据处理方法

Also Published As

Publication number Publication date
CN101916439B (zh) 2012-08-08

Similar Documents

Publication Publication Date Title
CN108921885B (zh) 一种综合三类数据源联合反演森林地上生物量的方法
Jin et al. Stem–leaf segmentation and phenotypic trait extraction of individual maize using terrestrial LiDAR data
Falkowski et al. The influence of conifer forest canopy cover on the accuracy of two individual tree measurement algorithms using lidar data
CN102959354B (zh) 用于利用LiDAR数据来分析树冠层的方法和装置
CN109031344B (zh) 一种全波形激光雷达和高光谱数据联合反演森林结构参数的方法
CN101916439A (zh) 基于希尔伯特-黄变换的空间碎片高光谱序列探测方法
Hauglin et al. Geo-referencing forest field plots by co-registration of terrestrial and airborne laser scanning data
CN101807301A (zh) 一种基于高阶统计量的高光谱图像目标检测方法
CN106778640B (zh) 一种三维可视化环境下地表植被覆盖模型的生成方法
Sun et al. Remote estimation of grafted apple tree trunk diameter in modern orchard with RGB and point cloud based on SOLOv2
CN103090946A (zh) 果树单树产量测量的方法和系统
CN1996044B (zh) 基于高空间分辨率遥感影像的林冠空间统计学定量估计方法
CN108898070A (zh) 一种基于无人机平台的高光谱遥感提取薇甘菊装置及方法
Yang et al. Forest canopy height mapping over China using GLAS and MODIS data
CN107220628A (zh) 红外干扰源检测的方法
Su et al. Tree skeleton extraction from laser scanned points
CN103926203A (zh) 一种针对地物光谱不确定性的光谱角度制图方法
Lantini et al. A deep learning approach for tree root detection using GPR spectrogram imagery
Wei et al. Spiral band model for locating tropical cyclone centers
Yao et al. Corn area extraction by the integration of MODIS-EVI time series data and China’s environment satellite (HJ-1) data
Celik et al. Identification of corn and cotton fields using multi-temporal Spot6 NDVI data
Zhang et al. Extracting planning areas of paddy rice in Southern China by using EOS/MODIS data
Liu et al. Fish-pond change detection based on short term time series of RADARSAT images and object-oriented method
CN106780532A (zh) 一种面向df图的互调频率方程快速检测方法
CN114813587B (zh) 一种花椒树冠层氮浓度无损诊断方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20120808

CF01 Termination of patent right due to non-payment of annual fee