提取滤波

提取滤波（精选七篇）

提取滤波篇1

近十多年来学者们提出了许多水稻种植制度的遥感监测方法, 使用的遥感信息源来自光学传感器和星载雷达。由于星载雷达数据相干斑噪声严重, 且容易其他因素的影响, 遥感光学数据依旧是研究和应用的主要数据源。在大区域作物遥感监测中, MODIS数据得到广泛应用, 如Wardlow等使用250 m分辨率MODIS数据进行水稻制图;Sakamoto等[4]使用多年MODIS数据监测越南水稻年变化[5];Boschetti等利用MODIS数据获得NDVI时相曲线进行水稻物候期提取[6];徐岩岩等基于MODIS-EVI提取东北地区水稻物候期等[7]。

针对植物而言, 它的物候包括发芽、展叶、开花、果实成熟、叶变色、落叶等生命活动的季节现象[8]。物候信息不仅反映当时的天气条件, 还反映过去时间段内气象条件影响下的累积情况, 在预告农事活动, 作物引种, 农业气象预报, 农作物估产, 区域气候乃至古气候的研究中都有广泛的应用价值[9]。传统的物候期观测方法一般通过人工目测来获取相关信息[10], 需要大量人力物力, 且周期长, 覆盖面积较小;基于积温进行预测是另一种较为传统的方法[11], 但这种方法需要大量的先验知识, 如播种日期, 日温度资料等, 这二者都因其明显的弊端限制着物候期提取方法的发展和进步。而随着遥感的大力发展, 卫星遥感数据可获得能够反映植被生长和覆盖的植被指数[12], 这为通过遥感监测物候期信息提供了条件。

近年来, 通过遥感结合滤波方法来获取植被的物候期信息已取得了初步成果。其中, 辛景峰等利用以旬为单位的合成AVHRR NDVI时序数据, 采用条件时间内插平滑法对数据进行平滑后, 对黄淮海的冬麦区的物候期进行了估测[13]。吴文斌等利用NDVI数据, 采用Asymmetric Gaussian方法重构ND-VI曲线, 提取了全国过去20年耕地第一生长季起始期时空变化[14]。鹿琳琳等基于SPOT/VEGETA-TION时序数据, 使用TIMESAT软件, 使用Savitzky–Golay滤波、Asymmetric Gaussian、Double Logistic曲线3种方法重构NDVI时序曲线, 进行物候参数提取[15]。孙华生等[16]、李正国等[17]、肖江涛[18]结合遥感数据也做了许多研究。通过前人的研究和结论, 可见遥感数据在物候期提取方面的作用和效果都非常显著, 这为以后的研究方向和方法提供了支持和参考。

作者在通过对研究区域水稻面积进行提取的研究基础上 (大气科学学报《基于Fast ICA算法和MO-DIS数据的水稻面积提取》一文) , 提取水稻的EVI数据, 结合不同的滤波方法和新颖的移栽期提取方法, 对研究区域进一步进行了水稻物候期的提取, 检验了前人的方法并结合新的种植作物和种植制度得出新的结论和问题。

1 研究区域和数据

1.1 研究区域概况

江苏省位于长江下游地区, 全省各地平均气温在15℃左右, 由东北向西南逐渐增高, 全省年平均降水量约为1 100 mm。省内平原面积70 000 km2, 占全省面积的70%以上, 主要作物有水稻、玉米、油菜、蔬菜等, 有“鱼米之乡”之称。水稻种植主要集中在苏南平原、苏中江淮平原、苏北黄淮平原以及沿海一些大型农场, 多为一季粳稻。全省水稻种植制度各不相同且水稻种植条件差异较大, 有着不同的水稻品种以及不同的水稻生育期。一般5月播种, 6月移栽, 10~11月收获, 2010年总种植面积达3 300万亩 (1亩=0.000 666 7 km2, 总产量达1 807×104t。

1.2 研究数据

从UGGS EROS数据中心 (http://glovis.usgs.gov/) 下载了覆盖整个研究区的2010年46景8 d合成地表反射率数据 (MOD09A1) 产品数据。MOD09A1数据经过大气校正、气溶胶校正以及卷云处理, 空间分辨率为500 m, 包含7个波段, 波段中心分别为648 nm、858 nm、470 nm、555 nm、1 240 nm、1 640 nm和2 130 nm, 包括两个表征质量的波段QA (quality assurance) 。另外, 获取了2010年江苏省近10年农业气象观测站获取了各观测点水稻生育期资料如表1所示。

2 研究方法

2.1 植被指数计算

本文使用ENVI软件中的波段运算功能, 在获取水稻像元后, 对文中需要使用的增强型植被指数EVI以及土壤含水量指数LSWI进行计算, 按照如下公式进行[19]:

式 (2) 中, ρnir、ρred、ρswir分别表示近红外波段 (841~875) nm、红光波段 (620~670) nm、短波红外波段 (1 628~1 652) nm、蓝光波段 (459~479) nm对应的反射率。

2.2 HANTS分析法简介和参数设置

HANTS分析法是以傅里叶变换为基础的谐波分析法, 它将NDVI时间序列表示为不同相位、频率和幅度的正弦函数组合[20], 能很好的利用遥感图像时空性的特点, 将其空间上的分布规律和时间上的变化规律联系起来。

算法的具体实现过程如下:首先, 通过傅里叶变换得到非零频率的振幅和相位, 其中各频率的振幅和相位都是通过不停的迭代获得[21];然后将所有的点进行最小二次方拟合。通过观测资料与拟合曲线的比较, 将那些差值明显大于设定阈值的点作为云污染点, 通过把它们的权重赋为零而拒绝参与曲线的拟合。基于剩余点进行新的曲线拟合, 通过反复迭代上述过程实现图像的重构。HANTS分析法在进行植被指数时序处理时, 除了要针对不同数据进行数据参数设定外, 还需要设置滤波重构过程中多个参数, 在参数的设置中, 由于不同的研究具有不同的需求, 只能根据一定先验知识进行多次试验来确定。由于使用的软件只支持8字节与16字节的数据格式, 所以需要在前期对需进行滤波重构的数据进行预处理, 本文中采用的其他各参数值为如图1所示。

2.3 小波变换简介、小波分解和重构

小波变换的概念是由法国从事石油信号处理的工程师Morlet J在1974年首先提出的。传统的傅里叶分析中, 信号完全是在频域展开的, 而小波变换与傅里叶变换相比, 它是时间和频率的局域变换, 因此, 它能解决傅里叶变换许多不能解决的困难问题[22]。

小波变换中包括连续小波变换、离散小波变换、多分辨率分析, Mallat算法以及小波包技术。由于本文需要的只是针对MODIS-EVI数据的小波去噪分析, 对于小波变换的需求仅限于较为基础的层面。下面将详细介绍文中使用到的小波变换知识。

小波变换主要分为信号分解和信号重构两大部分, 在信号分解中, 一般先使用连续小波变换 (CWT) 离散成离散小波变换 (DWT) , 在此过程中, 主要信息并不无丢失, 而针对离散小波变换进行分析后, 舍去不需要的信息, 逆向重构小波即可得到除去了无关信息 (噪声) 的重构信号。连续小波变换公式表示为:

式 (3) 中, c是连续小波变换后产生的小波系数, 称为基小波或者母小波 (mother wavelet) , 它是缩放因子a (scale) 和平移因子b (position) 的函数。连续小波变换中的缩放因子a和平移因子b都是连续变化的实数, 处理数字信号时很不方便, 实际问题的数值计算中为了减少计算量, 往往采用离散形式, 即离散小波变换, DWT可以通过离散化CWT中的缩放因子a和平移因子b得到, a与b通常选2j (j>0的整数) 的倍数。

离散小波变换将小波分为低频部分和离散部分, 进行滤波重构时往往去除离散部分, 对低频部分保留并对其进行小波逆变换, 重构低频部分, 因为低频部分代表了时间序列的主要特征。但在实际应用中, 往往还需针对高频部分进行阈值量化, 保留部分高频系数, 来达到阈值去噪的目的。一般阈值去噪处理分为3种:默认阈值去噪处理、给定软阈值去噪处理和强制去噪处理。常用的阈值方案有Rigrsure, Sqtwolog, Heursure和Minimaxi 4种[23]。在本文的研究中, 主要基于Matlab软件, 对EVI的时相曲线利用Wavedec函数进行了离散小波变换进行三层小波分解, 经多次比较后, 选择了“Sym11”小波函数, 再利用给定阈值的阈值方案重构了EVI时相曲线, 用以去除噪声。根据江苏省水稻的生育期长度, 节选了处于水稻生长期内23个时相的EVI数据。

2.4 物候期提取方法

水稻的生育期一般包括:播种、出苗、移栽、分蘖、拔节、孕穗、抽穗、灌浆和成熟。EVI时序曲线上的变化实时反映出水稻的生长情况, 所以, 可以通过识别EVI时序曲线上的某些特征点来获取作物的物候情况 (图2) , 一些关键物候期的识别方法如下:

移栽期:水稻移栽前后会进行灌溉泡田, 这会导致EVI的降低, 而随着移栽后水稻根部的发育, 水稻开始生长, EVI开始增长。因此, 一般都会针对EVI时相曲线中, 水稻生长季前期出现的最小值 (即曲线的一阶导数为0, 且由负变正) 后的拐点 (曲线的二阶导数为0, 且由负变正) 作为对应的移栽期, 或者选择这两个点靠后的点作物移栽期。而实际情况中, 这种方法仍然无法有效的获取移栽日期。

本文在相似性指数提取含水稻像元的基础上, 利用自编Matlab小程序实现水稻移栽期的提取。通过Xiao等提到的水稻移栽期间LSWI和EVI的特有关系[24], 结合多年的水稻生长发育期数据, 提取出水稻可能移栽日期所有的LSWI-EVI (定义为T) 的波段数据, 江苏地区即为时相数21~23的T数据, 根据T的正负情况, 归纳为不同的移栽日期, 如在22、23时相出现正值, 23时相代表期间即为移栽日期, 将这样的正负情况设定为一种目标结构, 符合这样的结构即输出对应的目标日期, 从而实现移栽期的区域性提取。

抽穗期:抽穗期的开始代表着水稻由营养生长的停滞向生殖生长的发育开始转变。水稻在幼穗分化后进入生殖生长阶段, 并在抽穗前完成最后一片叶的生长, 达到EVI指数的最大值, 之后叶片逐渐衰老枯死。EVI时相曲线的极大值一般出现在抽穗期附近[25]。故可将水稻生长季中期EVI时相曲线达到最大值时 (曲线一阶导数为0) 的日期估计为水稻的抽穗期。

成熟期:在水稻进入成熟期前, 开始灌浆, 先达到乳熟期, 此时穗粒中淀粉不断积累, 干、鲜重持续增加;而之后达到蜡熟期时, 穗壳则开始变黄并且由于籽粒灌浆对同化物的竞争, 叶片开始衰老枯黄, 植被指数开始迅速下降, 水稻成熟期也就是EVI指数下降速率达到最大的时期。以此为依据, 在EVI时相曲线上表现为曲线的二阶导数为0, 并且从负变为正, 成熟期往往位于抽穗期的后40 d左右。

3 结果和分析

3.1 HANTS与小波分析结果对比

分别使用HANTS分析法和小波变换对研究区域的水稻EVI时相曲线进行了滤波重构。通过对比分析发现, 在针对云影响较小的区域, 二者都有较好的去噪效果如图3 (a) ;而针对云影响较大的区域, HANTS分析法相比小波变换的去噪, 小波滤波的效果更加接近真实情况, 如图3 (b) 。由于HANTS分析法中剔除的值是用临近的值来替换, 去噪后该片段的数据无法反应真实情况, 会出现向前偏移或者向后偏移, 偏移方向取决于抽穗期前后数据的权重。

3.2 物候期提取结果与验证

利用Matlab软件, 基于前文提取获得的水稻像元, 针对移栽期使用自编程序进行了批对比提取, 获取了江苏地区2010年水稻移栽期空间分布图, 如图4 (a) 。针对抽穗期和成熟期, 使用小波变换的去噪重构结果, 针对上述物候期提取方法进行提取, 获取了江苏地区2010年水稻抽穗期和成熟期空间分布图, 如图4 (b) 、图4 (c) 。

其中, 移栽期:江苏地区的水稻移栽期主要集中在6月, 其中盐城北部和沿海地区、启东地及泰州东部地区移栽较早, 盐城南部、泰州北部地区和连云港东海、赣榆等地次之, 而苏北黄淮平原地区和苏南平原地区移栽日期较迟;抽穗期:苏北和启东地区较早, 苏中江淮平原次之, 盐城部分地区和泰州北部较迟, 而苏南地区呈现零散分布, 总体较苏北地区晚;成熟期:江苏地区水稻成熟期的空间分布与抽穗期大体接近, 值得提出的是, 由于水稻品种的差异以及中东部2010年9月长期云覆盖的影响, 给提取结果带来了一定的误差, 连云港东海县部分水稻的呈现较迟的成熟期, 盐城东台和启东如皋部分地区水稻同样如此。

江苏地区水稻种植历史久远, 相同地区各个年份之间水稻种植日期差距很小, 故选取徐州、赣榆、淮安、兴化、镇江、宜兴、昆山、高淳8个地区多个站点的2007年水稻发育期数据与结果进行对比, 进一步验证提取水稻物候期的准确性, 对比结果如图5。

从图中可以看出, 基于MODIS数据提取的所有站点的模拟结果相比实测结果的误差几乎全部都在±16 d内, 且大部分地区的误差在±8 d以内, 这说明由Matlab编程实现的移栽日期提取方法准确度较好, 能够很好的表现江苏地区水稻真实的移栽日期;而由小波变换进行滤波重构后进行的抽穗期、成熟期提取方法同样真实的反应了研究区域的实际物候情况, 这给利用遥感数据提取作物物候期进一步提供了研究基础, 也为研究遥感与作物模型结合提供了一定的数据基础。

4 讨论和结论

4.1 小波去噪阈值处理探讨

由前文可知, 小波分析过程分为分解和重构两个过程, 针对信号分解后的高频和低频部分, 对高频部分的阈值判断, 往往分为默认阈值处理, 给定软阈值处理以及强制阈值处理。其中, 默认阈值处理一般基于一个函数生成默认阈值, 再利用一个函数进行去噪处理, 给定阈值中的阈值往往可通过经验公式或者先验知识获得, 而强制去噪即为直接舍去高频部分, 这能得到较为平滑的信号, 却也丢失了有用的信息。随机提取研究区域中一个水稻像元, 选取23个 (19~41) 时相的EVI数据, 利用“sym11”小波对其进行了三层分解并分别进行了上述三种阈值去噪方法对比结果如图6。

由图6 (a) 可知, 强制阈值去噪可以得到非常平滑的信号曲线, 但重构后的曲线无法反应原始信号的任何特殊信息;默认阈值去噪可以得到较为平滑的信号曲线同时也能反应一定的原始信号;给定阈值去噪的结果完全可以由给定的阈值来决定, 从图6 (b) 可以看出, 分解得到的第三层高频系数能够反应出原始信号的突变。重构时对其部分保留可以较好的还原原始信息。为了更好地制定阈值方案, 随机抽取了部分出现突变的像元, 集成不同突变时期的像元作为实验样本。试验发现, 通过对第一层信号的部分保留以及第二、三层信号的去除可以很好的重构原始信号, 对比结果如图7。图中3种不同信号都进行如下阈值处理:第一层高频信号给定0.1为阈值进行部分保留, 第二、三层高频信号去除, 结合低频信号进行重构。结果表明, 该阈值方法能够较好的对原始信号进行滤波去噪并保留有效信息。因此文中, 整个研究区域的小波去噪采用了该阈值方法。

4.2 结论

提取滤波篇2

关键词：齿轮,故障诊断,卡尔曼滤波,局部均值

0 引言

随着科学技术的发展, 旋转机械在社会生活和生产中所占的比例越来越大, 当旋转机械发生故障时会影响生产, 带来巨大经济损失[1]。根据最近的研究显示, 旋转机械中齿轮故障所占的比重很大, 不但影响生产, 带来巨大的经济损失, 严重时还会危及人身安全。

齿轮的失效形式主要是断齿、点蚀、胶合、磨损、疲劳裂纹等, 当齿轮发生严重故障时, 会在噪声上产生异常, 人们可以通过噪声的变化很容易判断出故障, 但是在破坏初期或环境噪声过大时, 对于是否产生故障以及故障位置并不能很好把握, 这时就需借助信号处理方法来进行有效判断。当齿轮发生故障时, 其振动通常会产生幅值和相位的变化, 即其故障信号具有调幅和调频作用。因此, 如何准确找出这些变化是对齿轮故障诊断的关键, 而窗口傅里叶变换、小波变换和Winger分布等信号处理方法就能被很好地用于故障诊断中。然而, 大多齿轮故障振动信号都是非线性、非平稳的时变信号, 上述方法的局限性也就体现出来了。因此适用于非平稳、非线性信号的处理方法掀起了新的研究热潮。如经验模态分解法[2] (Empirical Mode Decomposition, EMD) 自提出后取得了广泛运用, 它能自适应地将复杂的多分量信号分解成多个单分量信号IMF之和, 通过求出每个IMF分量的瞬时幅值和瞬时频率, 得到完整的时频分布图。基于以上优点, 不少学者将此方法运用到齿轮故障诊断中[3,4]。局部均值分解 (Local Mean Decomposition, LMD) 是由Jonathan S.Smith[5]提出来的一种适合于非线性、非稳定性信号方法, 该方法最先用于脑电图信号的研究, 取得了不错的成绩。同时, 研究者发现, LMD方法同样适合用于振动故障信号检测, 并通过实验证明在一些方面LMD较EMD存在优势。

尽管LMD方法在故障诊断上已不断完善, 并取得大量成果, 但实际采集的信号中的噪声信号也会进行LMD分解, 产生模态混淆, 造成故障特征难以提取, 严重时会影响齿轮故障诊断的准确性。考虑到噪声的影响, 一些学者也想出了相应的解决方法, 如提出适合振动信号的小波降噪方法[6,7], 该方法中阈值的选择是关键, 直接影响到结果的准确性。因此, 为提高诊断的准确性, 提出了一种将卡尔曼滤波和LMD相结合的故障诊断方法, 试验研究说明了该方法的有效性。

1 卡尔曼滤波

卡尔曼滤波是由R.E.Kalman提出的去除噪声还原真实数据的一种数据处理技术。与维纳滤波根据当前估计值x (n) 和过去全部观测值x (n-1) , x (n-2) , …来估计信号的当前值s (n) 不同, 卡尔曼滤波是根据前一个估计值s (n) 和最近一个观测值x (n) 来估计信号的当前值s赞 (n) , 它的解的形式是以均方误差最小为原则下的系统的传递函数, 用系统方程和观测方程对信号进行估计[8]。因此在通信与信号处理等方面得到了很大的利用[9]。

设采集的齿轮振动信号为y (t) , 真实信号为x (t) , t=1, 2, 3…, N。噪声信号为Q, 则滤波方程可以表示为[10]

2 LMD

任意振动信号经LMD方法分解都能得到若干个PF分量和一个余量u, 每个PF分量可由对应的包络信号和纯调频信号的乘积得到, 其中包络信号即瞬时幅值, 而由纯调频信号可直接求出其瞬时频率, 整合所有的瞬时幅值和瞬时频率便可以得到原信号的时频分布[11]。其基本分解如下:

1) 找到信号x (t) 的全部局部极值点ni、相邻极值点ni+1, 由式mi= (ni+ni+1) /2求出相邻两极值点的平均值mi, 由式ai=|ni-ni+1|/2求出包络估计值ai, 用滑动平均法分别连接所有的平均值点mi和包络估计值ai, 得到局部均值函数m11 (t) 和包络估计函数a11 (t) 。

2) 从信号x (t) 中分离出局部均值函数, 得到h11 (t) =x1 (t) -m11 (t) , 用h11 (t) 除以包络估计函数a11 (t) 进行解调, 如式s11 (t) =h11 (t) /a11 (t) , 反复迭代, 直到s1n (t) 为纯调频信号止, 即a1 (n+1) (t) =1。把迭代过程中产生的所有包络估计函数相乘便得到包络信号。

3) 第一个PF分量可由包络信号a1 (t) 和纯调频信号s1n (t) 相乘得到, 即PF1 (t) =a1 (t) s1n (t) 。包络信号a1 (t) 就是其瞬时幅值, 由纯调频信号s1n (t) 可求出其瞬时频率。

4) 将PF1 (t) 从信号x (t) 中分离出来, 得到新信号u1 (t) , 再将u1 (t) 作为原始数据, 重复以上步骤, 循环K次, 直到uk (t) 为一个单调函数为止。

则信号x (t) 完整的时频分布可由所有的瞬时幅值和瞬时频率组合得到。

3 基于卡尔曼滤波的局部齿轮故障特征提取的应用

齿轮故障诊断中难以避免噪声的影响, 从而降低了结果的准确性。因此结合卡尔曼滤波和LMD方法原理, 在对信号用LMD方法分解之前, 采用卡尔曼滤波对有故障加速度信号进行降噪处理, 获取更真实加速度信号;通过LMD方法分解降噪后的信号, 获取PF分量的瞬时幅值和瞬时频率, 进一步对PF分量的瞬时幅值进行分析, 从而提取出具有故障特征的频率。

试验采用减速箱、数据采集仪和加速度传感器等装置, 并在试验前人为设置了齿轮故障。设置采样频率为1k Hz, 在不含负载情况下, 测得转速为n为500 r/min, 计算得到转频f为8.3 Hz。图1为采集的原加速度信号和将原加速度信号经过卡尔曼滤波降噪处理后的加速度信号。

利用LMD方法对原加速度信号和经过卡尔曼滤波降噪后的原加速度信号进行分解, 得到三个分量和一个余量。原加速度信号经LMD分解得到第一个PF分量的瞬时幅值a1, 如图2所示。用卡尔曼滤波降噪后经LMD分解得到第一个PF分量的瞬时幅值a1′, 如图3所示。进一步对a1和a1′进行谱分析, 分别得到如图4、图5所示的谱图。

从图4的幅值谱中可以看出, 在转频f的10倍频处有较为明显的冲击, 且存在一些虚假信息, 但初步可以判断齿轮存在故障, 对比降噪后的幅值谱图5可以看出, 在5f、10f、18f处存在明显冲击, 更加清晰地表明齿轮存在局部故障, 充分说明降噪后有利于判断齿轮的故障诊断。

4 结论

振动信号难免包含不少的干扰噪声信号, 本文的创新点在于先降噪, 减少一部分干扰信息, 再运用适用于非平稳、非线性的信号处理方法, 即LMD方法, 根据信号本身进行自适应分解, 得到具有物理意义的PF分量和瞬时幅值, 由幅值谱判断故障特征。将卡尔曼滤波降噪和LMD方法结合, 减少了虚假频率成分, 使获得的瞬时频率更加真实, 从而故障诊断的结果更具有准确性。

参考文献

[1]HE Zhengjia, CHEN Jin, WANG Taiyong, et al.Theory and application of mechanical fault diagnosis[M].Beijing:Higher Education Press, 2010.

[2]程军圣, 张亢, 杨宇, 等.局部均值分解与经验模式分解的对比研究[J].振动与冲击, 2009, 28 (5) :13-16.

[3]于德介, 程军圣, 杨宇.Hilbert-Huang变换在齿轮故障诊断中的应用[J].机械工程学报, 2005, 41 (6) :102-107.

[4]张强.基于EMD的齿轮箱故障特征提取方法研究[D].大连:大连理工大学, 2013.

[5]SMITH J S.The local mean decomposition and its application to EEG perception data[J].Journal of the Royal Society Interface, 2005, 2 (5) :443-454.

[6]寿海飞.基于小波分析的齿轮故障诊断[D].杭州:浙江工业大学, 2007.

[7]李鹏, 陈永当, 刘广梅, 等.小波分析在机械故障诊断中的应用[J].机电一体化, 2013 (10) :85-89.

[8]彭丁聪.卡尔曼滤波的基本原理及应用[J].软件导刊, 2009, 8 (11) :32-34.

[9]邱凤云.Kalman滤波理论及其在通信与信号处理中的应用[D].济南:山东大学, 2008.

[10]刘志川, 唐力伟, 曹立军.基于卡尔曼滤波和数学形态学分形分析的齿轮箱故障诊断[J].机械传动, 2014 (7) :141-144.

一种基于滤波器组的信号提取方法篇3

当今的电磁环境越来越复杂,信号的个数及类型越来越多,侦察接收机不得不考虑多个信号同时到达的情形,数字信道化接收机应用而生。数字信道化接收机具有宽频率覆盖范围、较好的灵敏度和实时信号处理等性能。侦察技术的发展对雷达造成了极大的威胁,迫使雷达不断发展,超宽带雷达的出现对数字信道化接收机提出了挑战[1]。

信道化的基本原理是将输入的全带信号进行频带分割,又称为子频段或子信道,然后对各子信道分别处理。为了得到更高的频率分辨率,各子信道可分别再进行第2次分割、第3次分割,直到满足频率分辨率的要求为止。

由于信道化接收机的信道数要预先确定,信道数的确定要同时考虑分辨率与监视频带范围的影响。面对超宽带雷达的侦察接收,跨道问题成为不得不解决的问题。许多专家学者已经研究了很多信道融合的办法,其中有基于滤波器特征研究频域加窗的信道拼接方法、基于接收信道的数学模型研究频域消去的信道融合方法等。

为了避免超宽带雷达信号出现信道化中的跨道问题,本文提出了一种利用全景监视与滤波器组相结合的方法来完成将多个信号中的有效信号进行完整提取[2]。全景监视即将监视频带内的所有信号进行检测,可采取FFT算法持续检测各个频点是否有信号存在。全景监视的目的是发现信号的到来与获得粗略的频率及带宽信息,然后根据各个信号的频率及带宽信息进行滤波器的选择,原始信号经过不同的下变频及滤波处理,各个信号则可被完整提取出来。该方法不仅能处理同时到达的多个信号又能解决超宽带信号的跨道截断问题。

1 整体方法研究

在信道化原理中是将监视频带内的所有信号都进行信道划分,然后分别进行各信道的信号检测来判断信号的到来,而实际中感兴趣的信号仅仅出现在1个信道或少数几个信道,并且当信号带宽比信道带宽宽时就会出现跨道问题,为后续处理带来困难。本方法利用信道化中滤波器组的思想,通过下变频与滤波器相结合的结构将不同信道的信号进行分隔[3,4]。

整体方案实现流程如图1所示。图中全景监视测得各个信号的起始频率及信号带宽,然后根据测得的起始频率进行下变频,并根据带宽选择滤波器带宽,最后将经过下变频的信号进行滤波处理,从而将有用信号从多个信号中提取出来。当关注多个信号时,根据不同信号起始频率与带宽信息,进行不同的滤波处理,从而将各个信号完整提取出来。

2 全景监视

为了便于表述,定义fk为各个信号的起始频率,Bk为各个信号的带宽,即为所关注信号所处频带( k = 0,1,2…) 。

全景监视即监视整个关注频带,通常采用的方法是做FFT在频域进行全频带监视,当有信号出现时,经过FFT运算,通过峰值搜索即可将相应频点搜索出来,从而实现了全景监视。由于全景监视的目的只是获得信号的到来与大概频率,不需要特别精细的频率测量,因此可选特定的数据点数做FFT[5]。

全景监视流程图如图2所示。首先对输入数据不断地进行特定点数N截取,然后将截取的N点数据进行N点的FFT运算,再将FFT的输出进行峰值搜索,搜索出信号的起始频率与信号带宽,当多个信号同时到达时,搜索出不同的起始频率fk与信号带宽Bk。

3 滤波处理

输入信号经过全景监视测得各个信号的中心频率fk与带宽Bk,然后将下变频后的信号经过滤波器组即可将不同的信号分离开来,实现多个信号同时到达的接收处理,通过选择不同的滤波器通带宽度从而避免了宽带信号的跨道问题[6]。

滤波处理模块图如图3所示,下变频量fk与滤波器通带宽度Bk的不同组合能够实现动态的信道划分与处理。为了便于硬件实现,需将fk与Bk进行量化,由于利用的是fk与Bk的组合来确定动态滤波器的作用频带,因而fk与Bk并不必要很精确,这种特性使得该方法适应性更好[7,8,9]。

由图3可知,该滤波处理主要由fk与Bk的量化模块、下变频模块和低通滤波模块3个模块组成。

1 fk与Bk的量化模块: fk与Bk根据实际需要可以量化为不同分辨率的有限个值fs1、Bk1,该量化模块可按fk/ M1、Bk/ M2取整的方式实现,其中M1、M2根据量化分辨率而定。

2下变频模块: 根据量化后的fs1产生下变频所需的下变频信号,然后将输入信号与下变频信号进行相乘,从而实现对输入信号的下变频操作。由信号处理理论可知,下变频的过程即进行频谱向下搬移的过程,由傅里叶变换的性质可知,频谱的向下搬移即对应时域信号与exp( - j* 2* pi* fk1/ fs* t) 相乘,即din_dds = din* exp( - j* 2* pi* fk1/ fs* t) ,din为AD采样后信号,din_dds为经过下变频后信号。在FPGA中实现时可直接利用丰富的IP核资源,在FPGA开发软件ISE中有DDS模块的IP核,该核能方便快捷地实现下变频操作。该核的实际功能是根据输入的下变频量,选择产生的exp( - j*2* pi* fk1/ fs* t) ,当输入信号到来时进行相乘操作。

3低通滤波模块: 根据量化后的Bk1进行低通滤波器系数的选取,使得经过下变频后的输入信号带宽完全处在滤波器的通带内,从而完全检出所关心的信号。利用Matlab中firpm能较方便地产生不同带宽的128阶低通滤波器系数,所产生的低通滤波器仅通带宽度与截止频率不同。将产生的每个滤波器系数进行量化,便于FPGA中实现,量化位数可根据实际情况改变。将量化后的所有系数存储起来,按顺序排列在FPGA的rom模块中,根据接收到的受保护雷达的带宽Bk的量化值Bk1进行选择提取不同的系数从而应用不同通带宽度的低通滤波器。

经过上述3个模块的处理,完整、有效地从多个信号中提取出所关注的信号[10,11,12]。

4 仿真结果分析

为了验证方法的有效性,进行Matlab仿真验证。仿真输入1个单频信号和1个线性调频信号的混合信号,fs= 500 MHz,单频信号频率为44 MHz,线性调频信号其实频率为20 MHz,带宽为10 MHz,信噪比选为10 d B,其中线性调频幅度为6,单频为1。

2个信号经过全景监视测得起始频率点分别为9. 8 MHz和44. 9 MHz,带宽分别为9. 77 MHz和0 MHz。经过量化后,信号1: fs1= 10,Bk1= 10;信号2: fs2= 44,Bk2= 0。为了将信号1提取出来,首先进行20 MHz的下变频,根据信号然后选择通带为12 MHz的低通滤波器进行滤波提取。经过下变频之后的2个信号如图4所示,频率起始频率分别为0 MHz和34 MHz,线性调频带宽保持为10 MHz。

经过通带为12 MHz的低通滤波器滤波之后的结果如图5所示。由图5可知,此时信号2单频信号已经被滤除,仅剩下单一的线性调频信号,且带宽保持为10 MHz。

经过仿真实现可以看出,将同时到达的多个信号经过全景监视,测得各个信号的起始频率与带宽,再根据起始频率进行下变频,将下变频后数据经过不同通带宽度的滤波器进行滤波处理即可完整提取出各个信号。仿真结果验证了算法能够较好地处理多个信号同时到达且避免了超宽带信号的跨道问题。

5 结束语

信道化能较好地解决多个信号同时到达的问题,但在信号带宽较宽时易出现不便解决的跨道问题,影响信号的完整性。本文采用先进行全景监视,获得信号的粗略信息后进行下变频及滤波器组的滤波处理,从而实现将有用信号分别从同时到达的多个信号中完整提取出来,避免信号跨道问题。通过仿真验证了该方法的有效性。

提取滤波篇4

多目标跟踪旨在分析量测以获取个数时变的目标位置和速度等状态信息,量测包含真实目标量测和由于杂波、噪声引起的虚警量测。由于无法预先获知目标个数,传统的基于数据关联的多目标跟踪方法在跟踪多目标时适应性不强[1,2,3,4]。Mahler基于严格的随机有限集数学理论[5],提出概率假设密度[6](PHD)滤波方法,通过递推多目标PHD,估计目标个数和位置,避免了数据关联问题。但PHD方法中存在高维积分,难以得到解析解。高斯混合概率假设密度滤波(GM-PHD)可对PHD滤波具体实现[7]。多目标PHD由一组高斯序列和,即所谓混合高斯近似表示。通过预测、更新步骤递推高斯混合函数中每一高斯分量的权值、均值和协方差,间接递推多目标的后验PHD。根据高斯分量信息估计多目标状态,当目标过于密集,噪声、杂波等干扰较严重时,估计精度会降低。为解决此问题,提出一种改进的基于高斯混合概率假设密度滤波的多目标状态提取方法,以满足较为复杂场景下的多目标状态提取要求。

1 高斯混合概率假设密度滤波

1.1 随机有限集与PHD滤波

随机集即为一个随机变量的集合,是指取值为集合的随机元,是概率论中随机向量的推广。当随机集的取值是有限值时就成为随机有限集(RFS)。多目标跟踪系统中,目标个数是一随时间变化的离散随机变量,因而将多目标集表示为随机有限集时,集合中随机变量的个数本身就是一个随机变量,导致多目标的状态空间和量测空间的维数也会随着时间变化。多目标的状态空间和量测空间分别以Es、Eo表示,则k时刻多目标状态集Xk和量测集Zk为:

其中Mk、Nk分别表示k时刻的目标个数和量测个数,其中量测来自于真实目标和虚警。

以随机有限集Ξk和Σk分别对多目标的状态和量测建模。假设Xk-1是RFSΞk-1在k-1时刻的具体实现,则k时刻多目标的状态以随机集表示为:

其中,Sk(Xk-1)表示k时刻继续存活目标的RFS,Bk(Xk-1)代表由k-1时刻的目标状态Xk-1衍生出的目标的RFS,Γk代表在k时刻自然产生的目标的RFS。多目标的量测以随机有限集可表示为:

其中Ek(Xk)表示产生于Xk的观测值的RFS,Ck(Xk)代表杂波或误警生成的RFS。

根据多目标Bayes滤波理论,可得:

其中,pk|k-1(Xk|Z1:k-1)和pk|k(Xk|Z1:k)分别表示多目标联合先验概率密度函数和后验概率密度函数。

随机有限集的PHD D(x)与随机变量的期望类似,即PHD D(xk|Z1:k)即为多目标后验概率密度pk|k(Xk|Z1:k)的一阶矩。在多目标跟踪系统中,∫sD(x)dx是出现在S内的期望的目标数,而D(x)的极值给出目标的状态估计。以PHD代替后验概率密度以贝叶斯滤波方式递推,便得到PHD的预测和更新方程如式(5)、式(6)所示[6]:

PHD滤波仍然包含多重高维积分,得不到闭合形式的解析解。高斯混合PHD滤波是PHD滤波的一种简单近似,解决了难以具体实现的问题。

1.2 GM-PHD滤波

高斯混合PHD滤波采用一组权值wk(i)、均值mk(i)和协方差Pk(i)随时间推移不断更新的高斯分量的加权和来近似后验PHD[7]。即

其中Jk-1表示k-1时刻高斯分量的个数。N代表正态分布。假设目标的运动观测均为线性模型。其转移概率密度与观测似然函数为

其中Fk-1、Hk分别为线性状态转移矩阵和线性观测矩阵,Qk-1、Rk为过程噪声和量测噪声协方差。

将高斯混合形式的PHDDk(x)代入PHD预测、更新方程,可得到高斯混合形式的PHD递推形式。k时刻的预测、更新PHD以高斯混合可分别表示为:

其中:

PHD更新完成后保留权值较高的高斯分量,将权值大于某一阈值的高斯分量的均值作为多目标估计状态。多目标个数可由权值之和求得。

2 改进的高斯混合PHD滤波

2.1 权值修正

由GM-PHD更新式(10),更新后的高斯分量总个数为Jk|k-1+Jk|k-1·Mk,其中Jk|k-1表示预测高斯分量个数,Mk为k时刻量测个数。第一项表示漏检高斯分量个数,第二项表示由量测更新的高斯分量个数。由式(11)可由量测zn更新Jk|k-1个高斯分量的权值、均值和协方差矩阵。均值和协方差的更新方法参照文献[7]。对于权值更新,令:

代入式(11),则第m个高斯分量由第n个量测zn更新后的权值为:

同理,Jk|k-1个高斯分量可分别由Mk个量测更新权值,由此组成一Jk|k-1×Mk维权值更新矩阵W。行对应某一高斯分量,列代表某一量测,第m行n列矩阵元素即为由量测zn更新的第m个高斯分量的权值wm,nk,update。由式(13),权值更新带有归一化操作,使得:

对于量测zn会出现以下两种情形:

第一种情形对应第n个量测zn由某一真实目标生成的情况。当观测空间没有虚警杂波量测时等号成立。第二种情形对应第n个量测zn为虚警杂波的情况,根据其与各个高斯分量均值的似然关系,其对应更新权值极小,表明该量测与当前预测高斯分量相关性极低。

以上两种情形体现在更新矩阵W中即为第n列的权值和不大于1。其现实意义是一个量测至多只由一个目标生成(考虑虚警情形)。每时刻完成高斯分量权值更新后,每一列的权值和必定满足上述性质,这是由GM-PHD滤波的递推特性决定的。

当虚警量测较多或目标较为密集时,多个量测集中于某一高斯分量均值m(m)k|k-1位置附近,有:

体现在更新矩阵W中即为第m行的权值和大于1,表示第m个高斯分量(隐指一个目标)生成了多个量测。在现实中,考虑漏检的影响,一个目标至多生成一个量测。因此这种更新方式会错误地把多个量测分配给一个目标并以较高权值进行更新,导致状态估计结果错误或冗余。

为此对矩阵W的行权值修正,使其行权值和不大于1,即满足一个目标至多生成一个量测的实际条件。

行权值修正后,该行中最大的更新权值wm*,n*k,update得以保留,保持第n*个量测相对第m*个高斯分量的高似然度,该行其余权值按各自比重调整,减小原有权值使得该行的权值和为1。

行权值修正操作会改变原有列权值的“归一化”秩序。对于第一种情形,即第n个量测zn由某一真实目标生成的情况,W矩阵中第m*行n列的权值由于行权值修正操作而减小。为使第n列的权值和满足式(14)的条件,也即使该列符合真实量测必定源于某一目标的实际条件,需要对第n(n≠n*)列进行列权值修正:

第n(n≠n*)列列权值修正过程保留因第m*行行权值修正而减小的(m*,n)处的权值,使得该列其他权值按各自比重相应增加,提高第n个量测相对其他高斯分量(不包含第m*个分量)的似然度。因为由行权值修正过程已得知,仅第n*个量测最有可能源于第m*个高斯分量或目标。

对于第二种情形,由于量测zn本身为杂波虚警,由其更新的各个高斯分量的权值极小,可视为零。故对第m*行行权值修正时不会改变(m*,n)处零权值的性质,对应第n列的权值和仍然满足式(14),不需要进行列权值修正。故在进行列权值修正前,应先大致判别该量测是否为虚警杂波,再决定是否修正列权值。当权值更新矩阵W中第n列的权值和高于某一阈值wth时可施行列权值修正操作。文献[8,9]中对所有的量测均进行重归一化处理会会导致资源的浪费,增加运算时间。权值修正方法流程如图1所示。

进行行权值与列权值修正的目的在于使更新权值满足以下条件:一个量测至多源自一个目标,一个目标生成至多一个量测。从而提高状态估计的精度。

2.2 高斯分量选择性合并

每时刻更新后的高斯分量个数会爆炸式增长,因此Vo提出在高斯分量权值更新后引入剪枝操作[7],将低于阈值T的高斯分量删除,仅保留高权值高斯分量,用以控制高斯分量个数,减小运算量。

对于剪枝后的高斯分量采取合并的方式,将多个相互靠近即满足如下条件的高斯分量合而为一[10]:

但当目标较为密集时,按照原有的合并规则会导致代表相近不同目标的多个高斯分量合并,即原本代表多个不同目标的高斯分量被合并后的一个高斯分量所代替,则由合并后的高斯分量估计出的目标状态会丢失部分真实信息。

为防止代表不同相近目标的高斯分量合并为一,对高斯分量的合并条件进一步限制。对于剪枝后的高斯分量集合,将权值高于T2的高斯分量{wk(j),mk(j),Pk(j)}认定为能独立代表某一真实目标,可直接根据其均值估计目标状态,不参与合并操作,防止因与其他高斯分量过近而被误合并,导致该真实目标信息的丢失。即满足如下条件的高斯分量才能进行合并操作:

选择性合并后的高斯分量的均值、权值和协方差矩阵的计算与文献[7]一致。对合并后的高斯分量集合进行状态提取,将权值较高的高斯分量的均值视为估计状态,至此完成该时刻的状态估计。将合并后的混合高斯PHD传入下一时刻。

假设新生目标和衍生目标的概率密度均可用如下高斯混合的形式表示:

k-1时刻的PHD如式(7)所示,与式(21)、式(22)一并代入式(5),得到如式(9)所示的预测后高斯混合形式。对于更新步骤,在原有更新权值的基础上[7]采用本文提出的权值修正策略,对权值调整,使其符合现实情况。针对更新后的高斯分量,加强合并条件,防止错误合并,权值高于设定阈值(通常取0.5)的高斯分量的均值即为k时刻多目标的状态估计。

3 仿真结果及分析

设定目标在某一二维杂波区域[-1000,1000]m×[-1000,1000]m内运动,目标的运动服从高斯线性模型为:

目标的状态向量包含位置和速度信息。采样周期T=1 s,过程噪声vk~N(·;0,σv2I2)为高斯白噪声,σv=0.3 m/s2。为简化讨论,不考虑目标衍生的情况。新生目标的RFS服从Poisson分布,其强度为:

各个参数为m(1)γ=[250,250,0,0]T,m(2)γ=[-250,-250,0,0]Tm(3)γ=[-250,-500,0,0]T,Pγ=diag([100,100,25,25])目标的存活概率为ps=0.99。

系统的观测模型为:

其中,wk是不相关的高斯白噪声,wk~N(·;0,σw2I2),σw=10 m。目标的检测概率pd=0.98。杂波均匀分布在观测区域。每一采样时刻杂波的个数服从均值λk为10的Poisson分布,杂波均匀分布在观测区域内。T=10-5,U=4,T2=0.5为相关的修剪与合并参数。仿真中共有6个目标,各个目标的运动信息如表1所示。目标的真实轨迹和量测分别如图2和图3所示。

图4和图5分别给出了GM-PHD滤波算法和本文提出的改进算法的多目标状态提取结果,实线代表目标的真实轨迹,“·”代表估计出的目标位置。从图中可以看出,改进算法的状态提取结果相对于原有算法的结果更为准确,估计出的目标位置要更加集中,与目标的真实轨迹贴合更加紧密。

为进一步分析状态提取的精度,采用OSPA准则[11]评估算法性能。OSPA评价指标引入了截断常数c,能够处理两个随机集之间目标个数不一致时的性能评价。在本文中截断常数c为100。图6给出了两种算法在20次MC仿真条件下的OSPA均值距离。距离越小表示估计结果越精确。由图可知,改进算法的OSPA曲线大体低于GM-PHD滤波状态提取算法,表明本文改进算法具有良好的估计精度。

当每时刻的杂波虚警个数的均值为100时,由真实目标生成的量测会被大量不相干噪声点湮没。即在不进行权值修正的情况下,会经常出现一个目标与多个量测相对应的情况,导致多目标状态估计的错误,为说明此问题,图7给出了两种算法在高浓度虚警条件下的平均OSPA距离。相比于低杂波浓度的情形(图6),两种算法的OSPA距离均有明显增大,这是由大量虚警量测导致的必然结果。但总体而言,改进算法的OSPA距离依然要低于GM-PHD滤波状态提取算法,说明了在目标、虚警密集等严苛条件下,改进算法的性能有所提升。

4 结语

本文在原有GM-PHD滤波方法的基础上,修正更新高斯分量的权值,并加强高斯分量合并规则,以减轻目标密集和虚警杂波对算法造成的影响。下一步可针对非线性非高斯模型进一步改进高斯混合PHD滤波状态提取算法,扩展其应用范围。

摘要：高斯混合概率假设密度滤波(GM-PHD)方法可有效解决线性高斯模型下的多目标跟踪问题,在估计目标个数的同时提取多目标状态。但当杂波浓度过高、目标过于密集时,GM-PHD的状态提取精度有所下降。针对GM-PHD滤波算法在复杂环境下性能下降的问题,提出一种改进的GM-PHD滤波多目标状态提取方法,通过修正高斯分量更新权值,强化合并规则,降低密集目标和杂波造成的干扰。仿真实验表明该方法能在不同杂波环境下提高多目标状态估计的准确度。

关键词：概率假设密度,高斯混合,多目标跟踪,状态提取

参考文献

[1]Samuel Blackman,Robert Popoli.Design and analysis of modern tracking systems[M].London:Artech House,1999:7-12.

[2]Poore A B.Multidimensional assignment formulation of data association problems arising from multitarget and multisensor tracking[J].Computational Optimization and Applications,1994,2(3):27-57.

[3]Barshalom,Rong Li,T Kirubarajan.Estimation with applications to tracking and navigation[M].USA:John Wiley&Sons,Inc,2001:33-39.

[4]林两魁.天基红外传感器对中段弹道目标群的跟踪与超分辨技术研究[D].长沙:国防科学技术大学,2013.

[5]Mahler R.Multitarget bayes filtering via first order multi-target moments[J].IEEE Trans.on Aerospace and Electronic Systems,2003,39(4):1152-1178.

[6]Mahler R.Statistical multisource-multitarget information fusion[M].Boston:Artech House,2007:711-715.

[7]Vo B N,Ma W K.The gaussian mixture probability hypothesis density filter[J].IEEE Transactions on Signal Processing,2006,54(11):4091-4104.

[8]Mahdi Yazdian.An improvement on GM-PHD filter for occluded target tracking[C].ICASSP,2013:1773-1776.

[9]Hongwei Ge.Improved gaussian mixture PHD for close multi-target tracking[C].Information Technology and Artificial Intelligence Conference,2014:311-315.

[10]Fuliang Yin,Liming Chen,Zhe Chen.A novel merging algorithm in gaussian mixture probability hypothesis density filter for close proximity targets tracking[J].Journal of Information and Computational Science,2013,8(12)2283-2299.

提取滤波篇5

齿轮传动具有传动力矩大、传动精度高、结构紧凑等优点, 是机械设备中必不可少的动力传动部件, 旋转机械约有10%的故障是由齿轮故障引发的, 因此, 齿轮故障特征参数的提取是旋转机械故障诊断的关键[1,2]。顺序形态变换是排序统计理论与形态学的结合, 引入了“循环统计学”的思想, 通过合适的组合可以改进形态学方法的缺陷。文献[3, 4]已将顺序形态滤波器用于旋转机械振动信号的降噪和转子轴心轨迹提纯, 均取得了较好的效果。对现场采集的含有大量噪声干扰的齿轮故障信号进行顺序形态滤波将有助于故障特征的提取。

奇异熵是对一维时间序列进行相空间重构和奇异值分解来计算熵源, 反映振动信号能量在奇异谱划分下的不确定性。信号成分越简单, 能量越集中于少数几个分量;反之, 信号成分越复杂, 能量就越分散。因此, 奇异熵可以作为振动信号非线性的一种度量[5,6]。

针对实测齿轮振动信号含有大量的噪声干扰而无法准确反映故障特征的问题, 本文采用顺序形态滤波与奇异熵相结合, 提出了一种新的齿轮故障特征提取方法。采用顺序形态滤波器对原始振动信号进行降噪预处理, 计算齿轮四种工况下信号的奇异熵, 并以此作为区分不同故障状态的特征, 实验结果验证了该方法的有效性。

1 顺序形态滤波原理

传统形态滤波处理过程是基于Minkowski (明可夫斯基) 运算的, 实质上就是一种极值运算, 不可避免地导致信号和噪声之间的误处理。顺序形态运算引入了“循环统计学”的思想, 通过合适的组合可以改进形态学方法的缺陷[7]。

设有集合A和B, 且0<μ (B) =k<+∞ (μ (·) 为测度) , 则集合A关于结构元素B的顺序形态变换A (p) B, (p=0, 1/ (k-1) , …, 1) 定义为:

式中, Bx={x-b|b∈A}表示结构元素B关于原点对称后沿向量x的平移, 变量p称为顺序形态变换的百分位。

式 (1) 的含义表示A (p) B是由其中至少含有A的[k- (k-1) p]个点的那些x组成的集合。

此时可定义复合顺序形态变换为:

设0≤p, q≤1, p, q=0, 1/ (k-1) , …, 1。令:

很明显, 传统的形态运算就是顺序形态变换的特例, 比如:

显然, 当对集合A使用不同的百分位值进行顺序形态变换时, 便可构造极值滤波、中值滤波、多尺度形态滤波等不同的传统形态变换。

为了同时去除信号中的正、负噪声干扰, 仿照Maragos定义的形态开-闭 (open-closing) 和闭-开 (close-opening) 滤波器进行组合, 可设计如下的顺序组合形态滤波器[8]:

式中, f (1/2, 1) B表示先采用结构元素B对采样信号f进行中值滤波, 以初步消除各种噪声的干扰, 然后进行膨胀运算。此时由于先进行的中值滤波, 可以较大修正膨胀作用的缺陷。同理, f (1/2, 0) B也可以修正腐蚀作用的缺陷[7]。

顺序形态滤波的效果与所采用的结构元素有着密切的关系, 相对而言, 结构元素越复杂, 滤波能力就越强, 但所耗费的时间就越长。综合考虑待分析信号特点与运算复杂度, 本文采用最简单的直线结构元素来处理含有噪声的振动信号, 结构元素的长度取3。

2 奇异熵

设离散信号x (i) , i=1, 2, …, N为观测的时间序列, N为采样点数, 基于相空间重构理论构造矩阵A为:

其中:1

矩阵A进行奇异值分解得到的奇异值为:

由此可定义观测时间序列的奇异谱为:

根据文献[6]可知, 信号能量可表示为:

3 齿轮故障特征提取实例分析

为了验证本文提出方法在齿轮故障特征提取中的有效性, 在齿轮试验台上分别对正常、齿面轻度磨损、齿面中度磨损和断齿四种工况下的齿轮进行了试验。被试齿轮转频为fr=23.6Hz, 啮合频率为fz=686Hz, 振动信号的采样频率为16384Hz。对齿轮四种工况分别采样, 各取10个样本, 先选取直线结构元素B={0, 0, 0}, 采用式 (3) 构造的顺序组合形态滤波器将原始信号进行降噪预处理以消除噪声的影响, 然后根据式 (8) 计算滤波后信号的奇异熵。以齿面轻度磨损为例, 图1给出了信号降噪前后的时域波形及其频谱。

对比降噪前后的图形可知, 信号经过顺序形态滤波降噪处理后, 原信号中含有的高频噪声得到了很好的抑制, 这对于信号进行后续的特征提取具有十分重要的意义。

图2给出了降噪后各种工况信号计算得到的奇异熵分布。从图中可以看出, 不同的故障状态对应的奇异熵具有明显区别, 而且奇异熵值分布比较平缓, 分类效果比较理想, 可用于齿轮故障诊断的特征提取与故障分类。

其中, “·”代表齿轮正常工况的奇异熵;“O”代表齿面轻度磨损的奇异熵;“+”代表齿面中度磨损的奇异熵;“x”代表断齿工况的奇异熵 (以下皆同) 。

为便于比较, 图3给出了各种工况信号降噪前的奇异熵。对比图2、3可知, 未经降噪处理的信号, 不同故障状态的奇异熵存在交叉, 而且奇异熵波动范围较大, 不能有效用于故障分类。这也充分说明了对原始信号进行降噪预处理的重要性。

4 结束语

本文将奇异值分解引入到齿轮的故障诊断中, 并定义了奇异熵用于对齿轮进行故障特征提取, 同时对传统形态滤波方法进行改进, 定义了顺序形态滤波方法, 然后将两者进行有效结合, 对实测齿轮故障信号进行了特征提取和故障分类。结果表明该方法具有较好的实用性, 为齿轮故障特征提取提供了一条新思路。

参考文献

[1]张业林, 程刚, 成钰龙, 等.基于小波包和SOM神经网络的齿轮故障诊断仿真研究[J].制造业自动, 2012, 34 (7) :82-84.

[2]张培林, 李兵, 徐超, 等.齿轮箱故障诊断的油液、振动信息融合方法[M].北京:机械工业出版社, 2011.

[3]Zhang W B, Wang H J, Teng R J, et al.Application of rank-order morphological filter in vibration signal denoising[C].Proceedings-2010 3rd International Congress on Image and Signal Processing, 2010, 8:4025-4027.

[4]Zhang W B, Su Y P, Zhou Y J, et al.Application of rankorder morphological filter in refinement of rotor center’s orbit[C].Proceedings-2010 4th International Congress on Image and Signal Processing, 2011, 4:2278-2280.

[5]王林鸿, 吴波, 杜润生, 等.用奇异谱和奇异熵研究数控工作台动态特性[J].振动、测试与诊断, 2012, 32 (1) :116-119.

[6]张小鹏, 范影乐, 杨勇.基于奇异谱熵的脑电意识任务识别方法的研究[J].计算机工程与科学, 2009, 31 (12) :117-120.

[7]欧阳森.基于顺序形态方法的电力信号处理方法[J].电工电能新技术, 2005, 24 (1) :45-48.

[8]张文斌.顺序形态滤波与样本熵在转子故障特征提取中的应用[J].制造业自动化, 2013, 35 (4) :79-81.

[9]赵学智, 叶邦彦, 陈统坚.基于小波-奇异值分解差分谱的弱故障特征提取方法[J].机械工程学报, 2012, 48 (7) :37-48.

提取滤波篇6

近年来利用先进的信息技术、电子控制技术和计算机处理技术开展的智能交通技术研究已成为我国大力扶持的一个极具应用前景和挑战性的重要研究方向。交通视频分析是智能交通的重要组成部分, 交通视频分析由三个步骤组成, 分别是局部特征提取和表示, 车辆检测和跟踪以及交通视频分析, 其中局部特征提取和表示是交通视频分析的基础和关键。

Gabor函数能够兼顾时域与频域的要求, 能充分描述图像的纹理信息, 并且具有其他提取方法所不具备的优越性, 因此Gabor变换很适用于应用在分析图像局部信息上。本文将研究基于2D Gabor滤波器的车辆局部特征提取方法。

1图像局部特征提取方法

为了满足车辆检测的实时性要求, 首先要针对车辆局部特征建立提取与表示的方法。边缘作为图像最基本的特征, 其包含了图像中用于识别的有用信息, 因此边缘检测是图像特征提取中的重要手段。边缘是图像分割所依赖的重要特征, 也是纹理分割的重要信息源和形状特征的基础;而图像的纹理形状特征的提取又常常依赖于图像分割。简单的边缘检测方法, 如soble算子, 是对原始图像按像素的某邻域构造边缘算子。这些特征提取方法一方面不够稳定, 容易受光照以及噪声等因素影响, 另一方面这些方法丢失的线条和纹理等其它局部特征对图像检测同样具有重要价值。近年来很多新的目标检测方法应用于人脸识别, 并取得了很好的效果, 例如DCT和LBP特征融合的人脸识别及Gabor小波和局部二值模式结合的人脸识别等, 车辆图像的特征提取以及车牌识别都可以借鉴和采用这些方法。

2基于2DGabor滤波器的车辆局部特征提取

2.1 Gabor滤波器简介

Gabor变换是英国物理学家Gabor提出来的, 具有最小的时频窗, 即Gabor函数能做到具有最精确的时间-频率的局部化。Gabor的小波特性表明图像的纹理信息可以用Gabor滤波来抽取, 这是Gabor滤波器的广泛应用之处。Gabor核函数由于去掉了直流分量, 对局部光照y影响不敏感[1]。Gabor滤波结果可以描述不同方向上灰度的分布信息, 并可以容忍图像有一定的平移、旋转、亮度不均、尺度变化等情况[2]。恰当地选择其参数, Gabor变换可以出色地进行图像分割、识别与理解。

一般的二维Gabor函数[3]可以表示为:

式 (1) 中,

(x′, y′) = (xcosφ+ysinφ, -xsinφ+ycosφ) 是直角坐标系旋转角度φ后新坐标下的坐标;U, V表示滤波器径向中心频率两个分量;g (x, y) 代表高斯函数:

从式 (2) 可以看出二维Gabor由二维高斯函数调制的复数正弦栅栏, 主轴相对于x轴旋转角度φ, σ为相对于y轴的方差。二维Gabor函数的空间频率响应为:

式 (3) 中, u, v表示空间频率变量,

(3) , H (u, v) , u旋转角度为径向中心频率

2.2图像预处理

首先把采集到的车辆图像归一化到256×256像素, 这样就可以使用相同的Gabor滤波器组参数进行滤波, 利用Gabor滤波器对车辆局部特征进行提取。灰度图像包含的信息更多, 所以在车辆局部特征提取中采用灰度图, 为减少非理想光照的影响, 凸显主要特征, 还要进行直方图平衡。

2.3 Gabor滤波器的选择

Gabor滤波器参数优化设计问题始终是难点, 总结前人的研究成果, 决定采用下面的方法。据前面的2.1分析可知, Gabor滤波器的位置由方向θ和径向中心频率F两个参数决定。在方向上取4个:0°、45°、90°和135°, 这四个方向能满足一般处理的需要。当然取的方向数越多处理的效果就会越好, 但这会使滤波器的个数大大增加, 其响应速度就会下降。

径向中心频率F的选择根据Jain[4]给出的方法:设一幅图像的长宽为Nc, Nc是2的幂次方, 则径向中心频率F可取如下数值:

据上面对方向和径向中心频率的分析可知, 选取了4个方向, 每个方向上有多个中心频率, 滤波器的数量还是很多, 计算量很大。这要求在使用时, 一方面根据图像的特征等选取最合适的方向和中心频率, 另一方面根据仿真处理的图像结果选取最合适的参数。

3实验结果

通过Matlab进行仿真实验。对输入的图像先进行归一化、灰度化和直方图处理, 分别送入到0°、45°、90°和135°等4个方向的2D Gabor滤波器中, 在整个实验仿真过程中要测试方向、径向中心频率以及滤波窗口三者对结果的影响。原图像如下图1。

首先设定高斯窗口为8×8, 频率 $F = (256 / 4) \sqrt{2} \approx 90$ , 测试方向对结果的影响。图2至图5是仿真的结果。

由图2至图5可以看出方向的变化并不会对最终的仿真结果造成过多影响, 四幅图像基本相同。出现这种结果的原因主要是

(1) Gabor滤波器允许图像有一定的平移、旋转, 图像的平移和旋转相对于滤波器变换方向是相同的。

(2) 对图像进行了归一化、灰度化和直方图处理, 这也会使图像最终的处理结果较好而减少滤波器方向变化造成的影响。当然这并不说明在所有图像的处理中都不需考虑滤波器的方向, 要根据不同图像及仿真结果具体处理。

下面看径向中心频率的变化会对仿真结果造成何种影响, 设定方向为45°, 高斯窗口为8×8, 径向中心频率取为F=45、23、11, 图6至图8为其仿真的结果。

由图6至图8的仿真结果可以看出, 径向中心频率的改变会对图像的处理结果造成很大的影响。在这幅图的处理中, 随着所选径向中心频率的下降, 仿真结果随之变差。寻找一种快速地确定径向中心频率的方法是以后需要深入研究的内容。

最后要确定滤波窗口的改变会对仿真结果造成何种影响。根据以上得到的结果, 设定方向为45°, 径向中心频率为F=90, 设定了两种极限的情况, 滤波窗口取为256×256和1×1, 仿真结果如图9、图10所示。

对比图9、图10可以发现, 滤波窗口的选择会对仿真结果造成很大影响。滤波窗口取的太大, 则不能提取图像的局部特征信息, 窗口取的太小, 会出现很多孤岛。经多次改变滤波窗口进行仿真发现, 针对本文中的此类图像滤波窗口选取为8×8至12×12效果较好。

窗口的选择不仅和图像本身的纹理特征有关, 也和图像所选取图像的像素有关, 寻找一种滤波窗口自适应的方法是以后研究中进一步讨论的问题。

应用上述仿真结果得到的结论, 对采集到的一系列车辆图像做了处理, 得到了较好好的结果。

4结论

Gabor函数能做到具有最精确的时间-频率的局部化, 很好地解决时频的“测不准”问题, 并且对非理想环境的适应性较强, 非常适用于图像的分割和提取。

本文通过构建2D Gabor滤波器来提取车辆图像中的局部特征, 仿真结果证明了该方法的有效性。

从实验的整个过程来看, 在以后的工作中应从以下两个方面做更深层次的研究。一是滤波窗口的选择会对最终结果有很大的影响, 如何根据不同图像选取相适应的滤波窗口要进一步研究。另外, 对于滤波器参数的选择也值得进一步研究, 本研究中主要以实验结果来确定参数, 如何根据图像快速准确地选择滤波器参数也是要进一步研究的问题。

参考文献

[1]黎涛, 罗代升, 吴炜, 等.Gabor变换的参数设计及其在车牌字符识别中的应用.中国测试技术, 2006;32 (1) ;127—129

[2]宋刚, 徐光佑.基于特征的人脸识别.北京:清华大学, 2001

[3] Leduc J P.Spatio-temporal wavelet transform for digital signal analy-sis.Signal Processing, 1997;60 (1) :23—41

提取滤波篇7

有源滤波器是治理电力系统谐波污染的有效手段之一,其基本原理是向电网中注入一个与负载谐波电流大小相等、方向相反的补偿电流,从而消除负载谐波电流对电网的污染[1]。有源滤波器能够有效地补偿电网中谐波的必要条件之一就是谐波电流的准确提取。迄今为止,已有多种谐波检测方法被提出[2,3,4,5,6]。如:基于瞬时无功功率理论的谐波检测方法、基于FFT(Fast Fourier Transformation)的谐波电流检测方法、基于小波变换的谐波检测方法、基于ip-iq变换的谐波检测方法、基于同步检测法的谐波检测方法[6]以及基于人工神经元网络的谐波检测方法等。本文采用的基于旋转坐标系下进行投影变换的方法,通过对电网电压和电流矢量在dq坐标系下进行投影变换,实现谐波和无功电流的提取,提高了谐波检测精度。该方法在电网电压畸变情况下仍能够很好地检测和补偿电网中的谐波电流和无功电流。

2 基于旋转坐标系下进行投影变换的谐波电流检测方法

2.1 理论分析

如图1,abc/dq代表旋转坐标变换,LPF为低通滤波器,PLL为锁相环。输入为待检测三相电流ai、bi、ci和电源电压ua、bu、cu,输出三相电流中除了基波正序有功分量之外的全部电流分量i*ca、i*cb、i*cc,有学者称其为有害电流。该电流即为有源滤波器的补偿电流指令。谐波提取过程如下:电流信号经d-q正变换、低通滤波器和d-q反变换,得到基波正序电压和正序电流分量,补偿电流参考指令为负载电流与基波正序有功电流之差。

在三相三线系统中,电压(ua,bu,cu)和电流(ai,bi,ci)量变换到d-q坐标系下,有

其中θ=ωt+θ0,θ0为电网同步参考正弦初相角,ω为电网基波角频率。上标“+”表示正序分量,“-”表示负序分量表示;下标n表示n次谐波分量。式(1)和(2)经低通滤波器滤波提取直流分量分别为

对式(3)和(4)进行d-q反变换,可得基波正序电压分量和基波正序分量:

其中

分别为基波正序电压矢量和基波正序电流矢ρ

量;ipd和ipq分别为基波有功电流分量ipρ在d轴和q轴上投影。对ipd和ipq进行d-q反变换,可求得三项基波正序有功电流分量ifa、ifb、ifc,如式(7)。

补偿电流参考指令为负载电流与基波正序有功电流之差:

要得到基波无功电流,只需将电流矢量iρ在与电压矢量uρ垂直方向上投影即可得到。

要提取任意次(如k次)谐波电流,坐标变换中基波角频率ω用kω(或-kω)来代替,则d-q坐标系将以k倍基波角频率作正向(或逆向)旋转。

2.2 仿真模型

采用MATLAB软件建立有源滤波系统数字仿真模型,见图2。模型主要包括:电源(Source)、非线性负载(Nonlinear Load)、谐波提取模块(Harmonic Extraction)、PWM信号产生模块(PWM Generator)。电路参数为:电源相电压220Vrms/50Hz,直流侧电容5000μF,直流侧参考电压800V,逆变器输出电感0.5m H,控制周期100μS,可控整流桥直流侧电阻1.5Ω,电感20m H。

2.3 仿真结果

图3为补偿前电源电流波形及其对应频谱分析,电流波形接近方波,总畸变率约为26.7%。图4为补偿后电源电流波形及其频谱分析,电流波形对称且接近正弦,总畸变率约为7%,谐波含量大大降低。

非理想电源电压波形如图5所示:电源中均含有10%的5次谐波和5%的7次谐波,形成非理想电压源,其它电路参数不变。负载电流波形及其频谱如图6所示,电流总畸变率约为23%,补偿后电源电流波形及其频谱分析如图7所示,电流总畸变率约为6.1%。

3 实验系统与实验结果

开发了容量为2k VA的并联型有源电力滤波器实验装置,主要参数如下:电源电压150V/50Hz,直流侧电压400V,逆变器输出电感1.25m H,控制周期100μS,谐波源为三相全桥二极管整流桥。实验结果如下图。注:以下波形图(图9～图12)横坐标均为10ms/格;其中频谱波形横坐标均为200Hz/格。

图8与图9由上而下波形依次是不对称三相输入电压(幅值不等)、原边电流、原边电流频谱。显然补偿后的原边电流接近正弦,图9与图8相比,5次、7次、11次、13次及更高次谐波已消除。实验结果表明,电网电压存在不对称情况下,采用该电流检测方法能很好地补偿谐波。

图10与图11由上而下波形依次是对称三相交流输入电压、原边电流、原边电流频谱。图11是只补偿了11次谐波的电流波形,从频谱可以看到,5次谐波分量非常接近零,与图11相比,补偿后5次谐波明显减小。说明该方法可以很好的补偿单次谐波。

4 结束语

本文介绍了有源滤波器的谐波电流检测方法。理论分析、仿真、试验结果均表明,基于该检测方法的有源电力滤波器在电网电压对称情况下能够准确地补偿非线性负载产生的谐波电流,具有良好的滤波性能,并且可以补偿任意次谐波分量。

参考文献

[1]顾建军,刘汉奎,公茂忠,徐殿国.基于电源电压瞬时相位检测的并联型有源滤波器研究[J].电气传动,2003,33(2):45-47.

[2]L.L.LAI,W.L.CHAN.Real-time frequency and har-monic evaluation using artificial neural networks[J].IEEETransactions on power delivery.1999,14(1):52-59.

[3]叶忠明,董伯藩,钱照明.谐波电流的提取方法比较[J].电力系统自动化.1997,21(12):21-24

[4]杨桦,任震.基于小波变换检测谐波的新方法[J].电力系统自动化.1997,21,(10):39-42.

[5]郝瑞祥.基于全数字控制的有源电力滤波器研究[D].博士论文,河北工业大学,2004.

[6]游小杰,李永东.Victor Valouch.郝瑞祥.并联型有源电力滤波器在非理想电源电压下的控制[J].中国电机工程学报,2004,(2):43-47.

本文来自古文书网(www.gwbook.cn)，转载请保留网址和出处

边缘滤波02-02

中频滤波02-02

滤波方法02-02

中国人的婚恋观02-02

直流系统滤波02-02

我的园子作文500字02-02

县全面学习贯彻落实干部监督四项制度措施总结02-02