多尺度模拟（精选7篇）
多尺度模拟 篇1
0 引言
随着地震勘探面临的勘探任务越来越复杂, 地震勘探成效在很大程度上, 取决于对基于符合实际介质模型的方法理论研究。本文就详细介绍了有限差分法对交错网格地震波场进行的模拟仿真研究。一般的有限差分地震模拟方法是规则的交错网格。但是当地质体强纵横向变速介质时计算效率却很低。为此, 人们研究了基于可变网格和不规则网格的地震波数值模拟方法。它对地质模型的离散化更为合理, 同时, 对于保持模型计算的灵活性也非常重要。从空间采样的角度考虑, 最有效的提高模拟精度同时又降低计算机内存需求的方法, 就是在模型的不同区域采用不同的网格步长, 即可变网格 (也称为不规则网格) [1]。
1 正演模拟方法原理
1.1 均匀各向同性完全弹性介质的弹性波动方程
由本质上讲, 地震波就是在地下岩层中传播的机械弹性波, 但反映弹性波传播的基本规律是由弹性波的波动方程[2]。由研究岩石的不同物理参数得到均匀各项同性完全弹性介质的波动方程, 即
式中:t为时间分量;u, v, ω为介质位移在x, y, z三个方向上的分量;σxx、σyy、σzz、τxy、τyz、τzx为应力分量;gx、gy、gz为体力密度分量;ρ为介质密度。
1.2 各向同性介质中的一阶弹性波方程
弹性波理论是弹性体受力和应变的关系, 根据各向同性介质表示的的本构方程和表示的柯西方程, 可以推导出一阶速度—应力弹性波方程。二维二分量各向同性介质中的一阶应力—速度弹性波方程 (假定体力为零) [4], 即
式中λ、μ表示拉梅弹性常数, ρ表示弹性体密度, σxx、σzz表示正应力, σxz表示切应力, vx、vz分别表示x方向、z方向速度分量。
2 网格有限差分法
2.1 有限差分法基本原理
有限差分法 (Finite Differential Method) 是基于差分原理的一种数值计算方法, 基本思想是:将波场域离散为许多小网格, 应用差分原理, 用差商来近似代替微商, 将求解的连续函数问题转换为求解网格节点上的离散的差分方程组的问题[5]。
2.2 高阶方程化为低阶方程的一般方法
解决高阶方程化低阶方程的问题, 其基本思路是将高阶差分分裂成低阶方程 (二阶) 。在每一延拓步长中交替使用低阶方程进行延拓, 即可达到高阶差分延拓的目的。
当色散方程Muir展开式的项数取得比较多时, 一般式为
并可 (在zj, zj+1间隔内) 分裂成下列方程组
频率———波数域方程组
时间———空间域方程组
在p1 (x, τ, t) 间隔内交替使用上式中的各二阶方程组, 即可完成高阶方程的延拓、偏移。
2.3 有限差分方程的应用过程三步骤
首先, 通过对研究体网格分解, 将连续函数空间分成许多不重合的子空间, 子空间采用正方形网格, 减少后面计算的复杂。其次, 由泰勒公式展开, 用差商代替微商, 将连续变化变量离散化处理, 得到有限差分方程组。最后, 编写程序, 求解有限差分方程组, 显示输出结果。
3 各向同性介质弹性波正演模拟结果
通过模拟地震波在均匀介质中的传播过程来分析可变网格差分正演模拟的精确度。
均匀介质模型的计算区域为均匀介质情况下弹性波脉冲, 纵波速度Vp=2000m/s, 横波速度Vs=800m/s, 震源在平面中心点。采样间隔Δt=2ms。
图1各向同性介质中网格:301×301, 空间步长Δx=Δz=5m时二阶差分模拟波场快照 (左图:水平分量;右图:垂直分量)
图2各向同性介质中网格:301×301, 空间步长Δx=Δz=10m时二阶差分模拟波场快照 (左图:水平分量;右图:垂直分量)
图3各向同性介质中网格:301×301, 空间步长Δx=Δz=20m;时二阶差分模拟波场快照 (左图:水平分量;右图:垂直分量)
通过模拟和仿真实验, 得到对应于典型地质模型的波场快照, 并依据这两种图像进行了对比, 分析了其相似和区别, 得到了与具体地质模型相符合的结论。
4 结论
依据上述研究, 本文通过对地震波动方程的有限差分化推导和fortran编程, 基于可变网格的技术, 对典型地质体模型进行了模拟实验, 得到了相对的波场快照, 并依据图像进行了对比, 分析了其相似和区别。证明了可变网格多尺度有限差分模拟方法是一种实用性好、精度高、计算效率高的地震波场正演模拟方法。
高阶有限差分法在对波动方程正演模拟时, 具有计算速度快、频散小的特征, 在对有限差分正演过程中必须注意震源、地质体、边界条件等因素的影响。
摘要：针对典型各项同性地质模型, 地震波波动方程和有限差分法对交错网格地震波场进行模拟仿真研究。为了得到一阶速度应力弹性波方程, 首先对地震波波动方程原始形式进行推演, 然后用时空二元高阶交错差分网格的方法对其进行离散化, 最后对地质模型模拟仿真实验, 并获取了对应的波场快照。最后由仿真结果表明, 算法计算精度准确, 计算效率高, 而且可变网格的变换对地质模型的研究起到了很重要的作用。
关键词：有限差分,地震模型,波动方程,正演模型
参考文献
[1]Smith D N, Ferguson J F.Constrained inversion of seis—mic refraction data using the controlled random search[J].Geophysics, 2000, 65 (5) :1622-1630.
[2]杜世通.地震波动力学[M].山东东营:石油大学出版社, 1999:238.
[3]董良国, 马在田, 曹景忠.一阶弹性波方程交错网格高阶差分解法[J].地球物理学报, 2000, 43 (3) :411-419.
[4]刘军迎, 雍学善, 高建虎, 杨午阳.多波多分量地震波场数值模拟及分析[J].石油物探, 2007, 46 (5) :451-456.
[5]左莹 (导师:李庆春) .基于高阶交错网格的有限差分地震波场数值模拟[D].长安大学, 2009-05-25.
多尺度模拟 篇2
1 颗粒接触理论
1. 1 黏结颗粒模型( BPM)
BPM模型的基本思想是将一群装配好的球形颗粒集合体,在各球的黏结点处通过添加平行键,形成可破碎的矿石模型如图1 所示[11]。颗粒间的力-位移行为由法向、切向的刚度kn和ks及颗粒间摩擦系数 μ 描述;平行键的力-位移行为由单位面积上的法向、切向的刚度、抗拉、抗压强度及键的半径倍数进行描述。其中平行键半径为两颗粒半径的 λmin(R(A)、R(B))。两个颗粒间的黏结键是一种弹性行为关系,键的总力Fi和总力矩Mi,可分解为切向和法向,即
式(1)中,分别为法向、切向的力、力矩; ni、ti为接触平面的单位矢量。
当平行键建立时,总力Fi和总力矩Mi初始值赋值均为零,随后的相对位移和旋转增量(ΔUn,ΔUs,Δθn,Δθs) 在每时步产生相应弹性力和力矩增量
式(2)中,A、I、J分别为平行键的初始面积、惯性矩和极矩,计算公式见文献[11]。
作用在平行键上的最大拉应力和剪应力为:
当最大拉应力超过抗拉强度或是最--大剪应力大于剪切强度时,平行键发生破裂失效。
1. 2 Hertz Mindlin无滑动接触模型
接触模型就是描述颗粒碰撞时接触的产生和发生作用过程。其中Hertz-Mindlin无滑动接触模型是在Mindlin取得的研究成果的基础上建立的[12]。当两颗粒球发生碰撞重叠时,法向接触力和切向接触力的弹簧阻尼简化系统如图2 所示。模型中设半径为RA颗粒球A和半径为RB的颗粒球B发生弹性接触碰撞,其法向重叠量 ψ 计算公式
式(4)中,r1,r2是两颗粒球心位置矢量。
颗粒碰撞后重叠面为圆形,其接触半径c采用下式计算
因此可求出颗粒碰撞重叠的法向力Fn:
接触模型中的法向阻尼力Fnd采用下式计算:
式(7)中,m*为等效质量。计算公式如下式所示
假设两颗粒碰撞前速度分别为vA、vb,则碰撞时法向单位矢量n为
式(7)中法向相对速度vn_rel、系数 α 和法向刚度Sn可由下式求得:
式(11)中,e为弹性恢复系数。
模型中颗粒间的切向力由式(13)求出:
式(13)中,Sτ、φ 分别为两颗粒切向刚度和切向重叠量,Sτ由式(14)计算而得:
两颗粒间的切向阻尼力
式(15)中,vτ_rel为切向相对速度。
合力矩为:
式(16)中,ωi为接触点处单位角速度矢量,μr为滚动摩擦系数,Ri为质点到接触点间距离。
式(5)、式(6)、式(14)中,R*为颗粒等效半径E*为等效弹性模量,G*为等效剪切模量,采用计算公式见文献[12—14]。
2 模型构建
离散元软件EDEM提供了生成多种粒度分布(如正态分布、对数正态分布、随机分布)的颗粒,现结合文献[15]选用双正态分布生成的具有双峰分布特性的颗粒集合体能较好地达到较高的配位数与较低的孔隙率。由于岩石内部结构的复杂性,同种矿物内部甚至不同矿物内部的粘结强度不同,基于Potyondy & Cundall提出的颗粒粘结模型[16]结合EDEM二次开发添加随机的粘结参数值来体现矿石内部同种矿物及不同矿物之间内聚力的差异性,本模型粘结参数如表1 所示。
随机选取30 mm矿石几何体建立颗粒模型如图3(a)所示,并建立矿石多尺度内聚颗粒模型如图3( b) 所示。图3( b) 即为构建的符合矿石特性同时符合冲击破碎特点的多尺度内聚颗粒模型。其中四种不同颜色的球,分别代表四种相同分布结构、不同类型颗粒集合体;图中颜色、大小、长都不同的连接线表示不同粘结键。
本征参数是指材料本身的特性参数,即泊松比、剪切模量和密度,这类参数一般可从文献和手册中查取。实验矿石为石灰石,冲击摆锤为碳钢材料,其特性参数如表2 所示,接触参数如表3 所示。
3 数值模拟及分析
3. 1 多尺度内聚颗粒模型破碎过程分析
由于冲击破碎过程时间短暂,采用实验方法分析十分困难,而数值模拟,能够较好地反映矿料破碎的每个时步破碎特征。根据构建的矿石颗粒模型分析冲击角30°、35°、40°、45°下矿料的破碎结果。通过EDEM后处理提取破碎后每个时步图片,其中有代表性的图片如图4 所示。
图4 为30 mm矿料颗粒模型在冲击角为40°下冲击破碎过程中的部分截图,其中图4( a) 至图4(g)颜色表示模型受挤压力的变化情况,挤压力随蓝、绿、红依次递增。图4(a) 中摆锤与矿料模型未接触,颗粒颜色为蓝色;图4(b)中,摆锤与矿料模型接触,内部颗粒挤压力剧增,并逐渐传递至整个矿料模型;图4(c)至图4(d)可以看出矿料在冲击破碎过程中具有一定的弹性变形阶段,内部黏结键挤压变形并开始断裂,同时由两图中颜色可以看出矿石模型中间部位受力最大。图4(e)至图4(g)为矿料塑性变形至断裂的阶段,此阶段矿料模型逐渐变形破碎,并伴随着内部黏结键快速断裂。图4(h)为破碎后的结果,这与实验中矿料破碎的结果较为接近。
3. 2 破碎黏结键及几何体破碎力分析
3. 2. 1 黏结键分析
为了对比不同冲击破碎能下模型内部黏结键的断裂数目变化情况,将黏结键刚初始断裂时作为断裂起始点,绘出不同冲击角下未断裂黏结键数目随时间的变化情况(图5)。
图5 给出了不同冲击角下未断裂黏结键数量与时间的曲线图,由图中可以看出,多尺度内聚颗粒模型的初始形成时黏结键数目为2 270,0. 001 s之前黏结键已开始发生断裂,且冲击角越大,黏结键剩余数越少,这符合冲击能越大矿料越碎的现象;0. 001s之后随着冲击角越大,黏结键断裂数量越大,这说明破碎时外界能量的大小与裂纹的扩展程度有关;0. 005 s之后模型内部剩余黏结键数量逐渐趋于稳定,且冲击角越大稳定的趋势越滞后;对比四条曲线斜率可知,在黏结键数目趋于稳定前,冲击破碎过程中黏结键随着时间变化具有一定的线性特征,且冲击角越大黏结键断裂越迅速,由此可得外界能量越大矿料破碎时内部裂纹扩展速率越大,并且黏结键数量和冲击能有某种线性关系。
根据剩余黏结键数与冲击破碎能数据,得到黏结键数目与冲击破碎能之间的曲线关系如图6所示。
图6 为冲击破碎能与剩余黏结键数目之间的拟合关系曲线,拟合度R2= 0. 993。由图中可看出冲击破碎能与破碎后剩余黏结键数目之间的数学关系可视为线性的。
3. 2. 2 破碎力分析
通过EDEM后处理器提取不同冲击角下双摆锤冲击过程中力和时间的数据,以摆锤与颗粒模型相接触时为起点,绘制不同冲击角下摆锤破碎力随时间的变化曲线如图7 所示。
图7 给出摆锤几何体破碎力随时间的变化曲线图,由图可看出0. 001 s开始,摆锤破碎力随着时间迅速增大,至0. 005 s左右达到最大值后急剧下降,0. 006 s时下降速率减小; 对比不同冲击角可以看出,在开始阶段冲击角越大,破碎力增加的速率也越大,0. 005 s之后冲击角对破碎力下降的速率影响不大。为分析最大破碎力与冲击破碎能关系,根据提取数据得出最大冲击力与冲击角度之间的关系曲线如图8 所示。
图8 给出了最大破碎力与冲击破碎能之间的拟合关系曲线,拟合度R2= 0. 994。由图可以看出,冲击破碎能与摆锤的最大破碎力具有较好的线性特性。
3. 3 破碎能耗及粒度分布分析
当前在EDEM仿真中,粒度分布主要在于获取碎后仍黏结以及未黏结的颗粒集合体在不同粒径筛网上的含量。为了近似给出其破碎后的粒度分布,通过PROE构建2、5、8、10、12、14 mm的物理筛网对其筛分,获取破碎后的粒度分布。筛分后通过EDEM中manual selection确定不同筛网上的每个球颗粒质量,得出各筛网上所占质量分数;根据所得质量分数得到35°、40°、45°冲击角下,矿料破碎后粒度分布概率密度曲线图如图9 所示,其中30°碎为4 块及部分碎屑。
图9 为不同冲击角下矿料颗粒模型破碎后的粒度分布曲线。图中可看出,6 mm以下的颗粒冲击角45° 的质量分数较高,而中间粒度颗粒质量分数较小,由此可知构建的矿料颗粒模型外部颗粒黏结不如内部紧密,在外力作用下易碎,同时由于破碎后仍黏结的颗粒在筛分时由于黏结不够紧密,易发生再次碎裂现象;曲线具有一定的单峰偏正态特性与真实岩石破碎后的粒度分布特征比较接近。
3. 4 冲击破碎能与冲击粒度分析
为研究实际与模拟下矿料冲击破碎后粒度分布与冲击破碎能之间关系,选取粒度为30 mm石灰石在不同冲击角下进行双摆锤冲击破碎试验。根据试验及模拟结果求得实验和模拟破碎后矿料的平均粒径并与冲击能建立曲线关系,结果如图10 所示。
图10 为30 mm矿料在冲击角为30° ~ 45°下的实验与仿真破碎结果对比图,仿真与实验所得曲线具有相同趋势。当冲击角为30°时,矿料与岩石颗粒模型都破碎为4 块,破碎后粒度接近;由图中可以看出,当冲击角大于35°时仿真得到的平均粒度略小于实验结果,仿真破碎结果与实验结果误差在10% 左右,这是由于筛分时矿石模型发生了再次断裂的缘故。
4 结论
利用EAlink将EDEM与ADAMS实行双向耦合,模拟了双摆锤冲击破碎过程,进一步验证了该颗粒模型的合理性,并得到以下结论。
(1)通过建立多尺度内聚模型数值模拟矿石冲击破碎过程中的弹性阶段、塑性变形以及断裂破碎阶段。
(2)通过提取EDEM后处理中黏结键断裂数目,建立黏结键数目与冲击能及破碎时间的变化关系曲线,为进一步研究岩石破坏过程中断裂裂纹与能量间的关系提供参考。
(3)通过冲击破碎能与破碎后的平均粒度之间的仿真和实验曲线趋势的一致性,表明建立的多尺度内聚模型能较好的模拟矿料的冲击破碎过程。
摘要：基于离散元方法采用Hertz-Mindlin接触模型和黏结颗粒模型,构建石灰岩的多尺度内聚颗粒模型。利用EDEM与ADAMS耦合将该模型应用于双摆锤冲击破碎过程,分析冲击破碎能对石灰岩冲击破碎的影响。结果表明:冲击破碎能与黏结键的数目有一定的线性关系,且冲击破碎能越大黏结键断裂越迅速,矿料破碎得越彻底;冲击角度越大,破碎力增加的速率也越大。最后该模型给出与试验相吻合的颗粒粒度分布曲线,说明该模型能够较好地模拟矿料的冲击破碎过程。
多尺度模拟 篇3
旋转机械在工业领域中的使用十分广泛, 振动故障是旋转机械各类故障中出现频次较高, 后果较为严重的类型之一。在旋转机械振动信号的各类参数图形中, 存在着大量反映设备运行状态的信息, 如二维幅频或相频特性曲线、小波图、趋势图、三维谱图、三维阶比图等。通过分析图形信息实现对旋转机械的状态监测是目前该领域的研究热点[1,2,3,4]。
在旋转机械故障诊断中, 振动信号参数图形的有用信息没有得到很好的利用, 这主要是由于施工环境较为复杂, 噪声干扰严重, 参数图形的边缘特征提取困难所致。数学形态学的主要研究对象是图像的形态特征, 此类特征可通过某种结构元素的形态与图像相应形态的对比方式来确定, 以此完成对图像的分析、滤波、识别、边缘检测、分割和重建等处理过程[5]。使用数学形态学相关方法分析处理旋转机械振动信号参数图形时, 可通过调节结构元素尺度来剔除环境噪声等干扰信息, 有效提取参数图形的边缘特征, 为进一步的故障诊断扫清障碍。
为此, 本文在文献[1]实验的基础上, 根据数学形态学Top-Hat变换和Bottom-Hat变换理论, 对旋转机械振动信号参数图形进行多尺度滤波增强处理;通过多结构元边缘检测方法对滤波处理后的旋转机械参数图形进行边缘检测。
1 数学形态学的基本原理
数学形态学的基本原理是通过一整套的变换来描述图像的基本特征和结构。数学形态学最基本的2种变换是腐蚀和膨胀, 其他变换都是由这2种变换的组合来定义的[6]。
1.1 形态学腐蚀运算、膨胀运算
设A为待处理的灰度图像, B为结构元素, 则结构元素B关于图像A的腐蚀与膨胀运算定义为
其中, DA和DB分别是A和B的定义域, 位移参数则必须包含在灰度图像A的定义域内。
腐蚀运算可以消除图形中的小成分, 可从内部对图形进行滤波;膨胀运算可以填充图形边缘处小的凹陷部分以及图形中比结构元素小的孔洞, 可从外部对图形进行滤波[7]。
1.2 形态学开运算、闭运算
形态学开运算、闭运算分别定义为
开运算在纤细处分离物体和平滑较大物体边界, 具有消除散点、毛刺和小桥等细小物体的作用;闭运算连接两个邻近的区域和平滑边界, 具有填充物体内细小孔洞的作用[8]。
1.3 形态学变换
形态学Top-Hat变换是对灰度图像做减去其开运算结果处理, 该变换可以提取亮度较高的背景中的较暗区域;形态学Bottom-Hat变换是对灰度图像的闭运算结果做减去原始图像处理, 该变换可提取亮度较低的背景中的较亮区域。形态学变换可用来提取目标图像中尺度小于结构元素的峰值和谷值[9]。
Top-Hat变换定义为
Bottom-Hat变换定义为
2 多尺度滤波增强处理
形态学腐蚀、膨胀、开、闭4种运算中的1种或2种串联或并联的组合就是形态学滤波运算。多尺度形态学滤波增强处理是通过不同尺寸的结构元素多次对图像进行滤波的, 其中多尺度开闭滤波在消除噪声、保持图像细节和提高信噪比等方面优于多尺度腐蚀膨胀滤波, 从而在一定程度上优化了灰度图像的有用信息, 令后续边缘检测结果更加真实可靠, 因此在形态学滤波中应用较多。
多尺度结构元素定义为
其中, B为十字形3×3结构元素, n为滤波尺度, 式 (7) 含义即为大尺度结构元素由小尺度结构元素多次膨胀得到。
为了得到足够平滑的图像, 本文采用最大尺度的结构元素Bn对图像进行多尺度开闭滤波增强处理, 其表达式为
其中, 权值ω对最后的滤波增强结果有较大影响, 一般取为0.5, 本文根据滤波增强处理结果的优劣, 取0.3。
图像经过多尺度开闭滤波增强处理后得到足够平滑的低频图像, 为获得更全面的有用信息, 还需提取图像的高频细节信息。在多尺度滤波增强处理方法中, 由于噪声在经小尺度结构元素处理的图像中出现几率较大, 并且随着尺度的增加其影响逐渐消失[10], 故本文选用带有修正系数的Top-Hat变换 (FT (i) ) 和Bottom-Hat变换 (FB (i) ) 来提取图像的高频细节信息。为减小噪声对图像的影响, 修正系数设定为公比为0.5的等比数列, 此过程完成了不同尺度间小尺度图像特征的平滑处理, 具体的表达式如下:
由多尺度开闭滤波增强处理的图像最终由三部分组成:第一部分是图像经最大尺度结构元素开闭滤波增强以后生成的低频平滑图像, 该部分包含图像中的大尺度图像信息;第二部分是提取比该滤波增强尺度还小的亮点图像高频特征;第三部分是提取比该滤波增强尺度还小的暗点图像高频特征。至此, 一幅灰度图像经多尺度滤波增强处理后生成的图像为[11]
3 多结构元边缘检测算子
在图像边缘检测处理中存在着多种梯度, 若在某一像素点处梯度值大, 则表示在该像素点处图像的灰度值变化迅速, 从而认定该点可能是图像的边缘点。数学形态学边缘检测方法主要是利用形态学梯度来完成图像的边缘检测。若将数学形态学的腐蚀、膨胀、开、闭等基本运算用于图像处理, 可构造出合适的形态学梯度算子 (经典边缘检测算子) 用于图像的边缘检测[12]。
腐蚀型边缘检测算子:
膨胀型边缘检测算子:
膨胀腐蚀型边缘检测算子:
上述3种形态学边缘检测算子是一种非线性的差分算子, 这些算子容易实现, 在实际中有一定的应用。但是, 这些算子对噪声都很敏感, 不能在保持较高检测精度的同时又不损失抗噪性能。由于旋转机械振动信号中普遍存在噪声, 虽然已经过多尺度滤波增强处理, 但仍有少量残留, 而且噪声信号和参数图形的边缘又均为频域中的高频分量, 因此, 为了更好地提取旋转机械振动信号参数图形的边缘特征, 应选择抗噪性能优于经典边缘检测算子的方法对参数图形进行边缘检测。根据腐蚀、膨胀、开、闭4种运算抑制噪声的相关特性, 本文对式 (12) ~式 (14) 做如下改进。
抗噪腐蚀型边缘检测算子:
抗噪膨胀型边缘检测算子:
抗噪膨胀腐蚀型边缘检测算子:
数学形态学边缘检测方法不仅与所使用的边缘检测算子有关, 还与结构元素自身特点密切相关, 如大小、方向、形状等。在边缘检测过程中, 不同结构元素对图像不同边缘细节信息的敏感性各不相同, 一种结构元素只能提取图像的一种边缘信息, 这不利于保持图像边缘的有用信息。因此, 应尽量选用具有不同特征的结构元素对图像进行边缘检测, 让每个结构元素都发挥作用, 提取出具有其自身特征的边缘信息, 这样可以充分保持图像的各种边缘信息, 达到既能检测出图像的各种边缘纹理, 又能抑制噪声的目的[13]。本文利用抗噪膨胀腐蚀型边缘检测算子 (式 (17) ) 构造多结构元边缘检测算子, 其表达式如下:
其中, B1、B2、B3为结构元素, 尺寸固定不变 (3×3正方形) , B1, B2可取为同一种结构元素, 也可取为不同的结构元素。
4 多尺度多结构元边缘检测仿真
为验证多尺度多结构元边缘检测方法的正确性与有效性, 本文选取结构元素B1=[1 2 1;2 62;1 2 1], B2=[0 1 0;1 1 1;0 1 0], B3=[1 0 1;0 1 0;1 0 1], 对含有5%椒盐噪声的Lenna灰度图像进行多尺度多结构元边缘检测, 其中多尺度滤波增强处理使用结构元素B1作为初始结构元素, 滤波尺度n取4, 多结构元边缘检测算子使用结构元素B1、B2、B3进行检测。图1a为原始灰度图像, 图1b为边缘检测结果。从图中可以看出:多尺度多结构元边缘检测方法滤除了Lenna图像中的椒盐噪声, 检测出的图像边缘轮廓清晰、纹理明确, 信噪比有所提高。该方法边缘检测效果优于经典边缘检测算子边缘检测效果, 更适用于含有噪声污染图像的边缘检测。
5 旋转机械参数图形边缘检测实例
5.1 实验
旋转机械故障模拟实验在600MW超临界汽轮发电机组轴系试验台上完成, 分别进行了转子正常、转子不对中和轴承松动故障的实验。试验台主要包括5个部分, 即发电机组轴系、润滑系统、动力系统、供气系统和信号采集分析系统。其中发电机组轴系由9个轴承5跨组成;润滑系统用独立的油路系统对各个轴承供油, 每个轴承座均安装BENTLY3000 XL8 mm电涡流传感器, 输出为7.87V/mm;动力装置采用55k W变频电机经过FRENIC变频器输出转速和功率, 并采用HG0G-C2型变速箱, 试验台详细结构布置如图2所示。在实验过程中, 采样时间为0.64s, 采样频率为转速的32倍, 实验时转子最高工作转速为3200r/min, 采集的信号经A/D卡传送到计算机, 为后续的数据分析做准备[7]。
实验中对转子正常、转子不对中及轴承松动故障, 每种采集40个启停机样本, 共计120个。首先将每个原始振动信号的采集样本进行处理, 生成各自的振动三维谱图, 如图3所示。
5.2 多尺度多结构元边缘检测
根据三维谱图倍频特征明显的特点, 将频率作为横轴, 转速作为纵轴, 像素点灰度值作为该转速下、该频率下幅值的大小, 将其转化为二维灰度图形, 结果如图4所示。灰度图中明显的竖线为倍频线, 与三维谱图中的倍频线相对应。
为了有效地提取旋转机械振动信号参数图形的边缘特征, 本文对图4各种状态下的参数图形进行量化、直方图均衡化等预处理, 选取结构元素B4=[1 3 1;3 5 3;1 3 1], B5=[0 1 0;1 1 1;0 10], B6=[1 0 1;0 1 0;1 0 1], 应用上述多尺度多结构元边缘检测方法对其进行边缘检测。其中多尺度滤波增强处理使用结构元素B4作为初始结构元素, 滤波尺度n取4;多结构元边缘检测算子使用结构元素B4、B5、B6进行检测, 最终的多尺度多结构元边缘检测结果如图5所示。从图5可以看出:旋转机械振动信号参数图形经多尺度多结构元边缘检测处理后, 噪点大幅降低, 环境污染噪声基本被滤除干净, 有用信息得到保持的同时信噪比大幅提高, 边缘鲜明, 轮廓清晰, 充分保持了图形的细节特征。至此已说明多尺度多结构元边缘检测方法能够有效地提取旋转机械振动信号参数图形的边缘特征, 具有较强的抗噪声干扰能力, 适合在环境比较复杂、噪声污染较为严重的情况下对旋转机械实施状态监测。
6 结论
(1) 依据数学形态学多尺度图形处理方法, 结合Top-Hat变换和Bottom-Hat变换处理方法, 选取合适的结构元素, 在对旋转机械振动信号参数图形进行有效滤波的同时, 可以保持图形的高频细节特征, 增强参数图形的有用信息, 提高参数图形的信噪比。
(2) 运用多结构元边缘检测算子检测旋转机械振动信号参数图形的边缘, 能够有效剔除多尺度滤波增强处理过程残留的噪点信息, 提取的参数图形边缘特征质量较高。
(3) 在实际应用中, 结合旋转机械振动信号参数图形及其噪声的特点, 多尺度多结构元边缘检测方法可以较好地解决边缘检测精度与抗噪声性能的协调问题, 为基于振动三维图形的旋转机械故障诊断奠定基础。
多尺度模拟 篇4
标准Kalman滤波是用于线性系统的最小均方意义上的最优状态估计。Kalman滤波对非平稳信号具有较强的估计能力,但标准最优Kalman滤波的缺点和局限性是要求知道系统的数学模型和噪声统计特性。所以,在单一尺度上,如果由于多种不确定因素干扰,导致使用不精确的模型和噪声统计,设计出的Kalman滤波器会导致滤波器性能变坏[1]。
多分辨率的数据进行多级别、多层次的处理能够获得更有价值的信息,而这种信息是单一尺度上所无法获得的。通过多尺度分解,在不同的尺度空间上分别刻画目标的特征属性再加以融合,可以有效地降低不确定因素的干扰,从而提高滤波器的滤波性能[2,3,4]。小波变换的多尺度特点非常适合多尺度信号的处理,特别是依靠小波变换的重构算法对低分辨率补充一定信息可以获得高分辨率近似补充。文献[5]中提出基于小波变换,寻找基于全局信息的一步预测值和估计误差协方差矩阵,通过综合不同尺度的信息来降低Kalman滤波对数学模型和噪声统计特性的依赖,这种方法比起单一尺度上的Kalman滤波效果是明显的,但计算量较大,由于算法侧重于横向数据更新,全局最优预测值和误差矩阵对滤波效果的影响变得较大,容易受到干扰。本文就此提出一种基于多尺度Kalman的融合滤波方法。充分利用初始估计序列来有效地降低向量和矩阵维数,减少运算量,把标准Kalman滤波只在单一尺度和时间轴上,对状态估计值和误差协方差进行数据更新,改进为基于小波变换在尺度轴和时间轴上,通过不同尺度的观测数据,纵向和横向双向数据更新,算法侧重于纵向数据更新,充分利用多尺度数据进行多层次融合处理,从而获得更好的滤波效果。
本文主要思想可以概括为:首先,在最细尺度N上通过标准Kalman滤波得到初始估计序列,然后通过小波变换将状态估计值和状态估计误差协方差的分块序列分解到最粗尺度,并根据相应尺度上的观测数据进行数据更新,然后通过小波重构,实现不同尺度上的数据更新,当纵向回归到原始尺度即最细尺度时,再跳至下一个数据块进行,直至滤波结束。在实际中,也可以利用此方法解决多传感器多分辨率数据的融合,更好地利用不同分辨率数据的互补信息,达到更佳的融合效果。
2 多尺度动态系统的描述及Kalman最优滤波
通常情况,动态系统状态方程在尺度N上的描述为
其中:k为离散时间变量,x(N,k)是状态向量,Φ(N,k)是系统状态转移矩阵,w(N,k)是系统噪声,且服从均值为零,方差为Q(N,k)的正态分布,初始值x(N0,)为一随机变量,其均值和方差是E[x(N,0)]=x0和E{[x(N,O)-x0][x(N,O)-x0]T}=P0。
动态系统量测方程在尺度N上的描述为
其中:C(N)系统观测阵,z(N,k)是对x(N,k)的观测值,观测噪声v(N,k)服从均值为零,方差为R(N,k)的正态分布,初始值x(N0,)与观测噪声v(N,k)以及系统噪声w(N,k)之间是互不相关的。
在尺度N上标准Kalman最优滤波的基本方程,如式(3)~式(7)所示:
3 基于小波变换的多尺度动态系统描述
我们由小波理论可知[6],通过一个脉冲响应为h(l)的低通滤波器可以从尺度i上获得粗尺度i-1上的平滑信号xL(i-,1k),脉冲响应为g(l)的高通滤波器可获得细节信号xH(i-,1k):
通过小波变换将状态方程从尺度i分解到粗尺度i-1,我们可以发现不同尺度之间的参数矩阵之间的内在联系,为了计算方便,一般我们取状态方程中的G(N,k)为单位阵,则i-1尺度上系统状态方程为
其中上标i表示结果由尺度i上小波分解得到的,比较式(10)和(11)可得式(12)~式(13):
由式(13)可推导出:
同样,用小波变换将量测方程从尺度i分解到尺度i-1,则i-1尺度上系统量测方程为
比较式(15)和式(16)可得式(17)~式(18):
从式(12)∼式(13)和式(17)∼式(18)可以看出,通过尺度i上的系统状态和量测方程的参数矩阵,可以很容易地获得尺度i-1上的相应的参数矩阵,研究标准Kalman最优滤波的基本方程可知,在尺度i-1上进行标准Kalman滤波估计还需要和Pk+1/k(i-1),我们将在多尺度Kalman滤波算法中给出计算方法。
4 基于小波变换的多尺度Kalman滤波
首先,假设一长度为T的信号,我们在单一尺度N上用标准Kalman滤波根据观测值进行估计。由于多种不确定因素干扰,滤波效果并不会很理想,在尺度N上,经Kalman滤波后,我们将得到估计序列和相应的估计误差协方差阵序列Pk/k(N),序列长度均为T,我们把得到的这一系列估计序列定义为初始估计序列,将这些序列分割成长度为M的数据块,其中T和M均设定为2的整数次幂,方便后续的小波变换。这样分割后的序列可表示为Xm(N)和Pm(N),则:
对于数据块Xm(N)中的每一个数据元素均来自于初始估计序列的估计值,由于我们要在尺度轴上对此估计值进一步修正,所以我们可以将此值看作一步预测状态估计值,为了保持与标准Kalman最优滤波方程中的符号的一致性,我们用符号来代替数据块Xm(N)的表示,此时是由初始状态估计值组成的数据块。
由式(8)~式(9)可知:经过从尺度N到尺度i的小波塔式分解,信号可分解为尺度i上的平滑信号和相应的各尺度l(i≤l≤N-)1的细节信号,这样我们把初始状态估计序列通过小波变换向尺度i分解,可得到尺度i上的平滑信号和相应各尺度l上的细节信号,如式(21)~式(22)所示:
对于初始状态估计序列,分解后不同尺度上的序列集合为
接下来,我们要获得相应尺度上的预测误差协方差阵。由估计误差协方差阵的定义可知:
其中:为尺度N上的一步预测误差,由于我们是对长度为M的数据块操作,所以此时估计误差协方差为一矩阵块。根据初始估计误差协方差阵序列,我们可以获得数据块Xm(N)的估计误差协方差阵为,如果用表示第i行第j列的元素,则:
同样,为了保持与标准Kalman最优滤波方程中的符号的一致性,我们用符号Pk/k-1(N)来代替数据块表示。这样对Pk/k-1(N)做小波变换实际上是一个二维离散小波变换[6],通过小波变换我们可以得到尺度i上的Pk/k-1(i)。
若在单一尺度上测量,我们可以通过下采样Zk(N)获得尺度i上的测量序列Zk(i),若是多传感器多分辨率系统,则尺度i上的测量值Zk(i)可以通过不同传感器的观测获得。在尺度i上,我们利用Kalman滤波方程对数据块进行数据更新。由式(17)、(18)、(23)和(26),我们可以得到尺度i上的Ck(i),Rk(i),和Pk/k-1(i),并通过式(27)~式(31)在尺度i上Kalman滤波:
数据更新结束后,通过小波重构我们可以得到尺度i+1上的以及Pk/k-1(i+1),在通过尺度i+1上观测值Zk(i+1)用同样的方法对序列进行数据更新,以此类推,直至尺度N,我们便得到了尺度N上的。同时在小波重构时,还可根据实际的研究情况,结合一定的小波域值去噪方法来去除部分噪声影响[7,8]。
5 仿真算例
实验一
系统在尺度4上的状态模型和尺度j(j=4,3,2,1)上的量测模型,如式(32)~式(33)所示:
其中:、w(,4k)~N[,0Q(,4k)]、v(j,k)~N[,0R(j,k)];式中Q(,4k)=diag1(1,),测量矩阵C(j,k)是单位矩阵,测量误差协方差矩阵R(j,k)为
尺度4为最细尺度,初始值x0=[10,]0T和P0=diag(4,4),选用d B4小波分解和重构。图1(a)信号在尺度4上的量测值(T=512,含两个状态分量)。图1(b)为最细尺度(i=4)上Kalman滤波结果。图1(c)为融合4层信息后的多尺度Kalman滤波结果(M=32),可以很直观地看出,多尺度Kalman滤波的有效性。
从表1中可以看出,融合的尺度越多,测量值的改进越明显,但当融合尺度数目过多时,测量值的改进程度明显减少且运算量随之增大,这是由于测量值提供的有用信息变少了,在粗尺度上的数据更新也就相应变少。从表2中可以看出,不同数据块大小对滤波效果也有明显影响,当数据块选取过大时状态估计误差协方差矩阵较大,数据之间的影响变大,导致滤波性能降低。所以,我们可以确定实际需求,确定适当融合尺度数目和数据块大小有效地提高滤波性能和效率。
实验二
继续采用式(32)∼式(33)的系统状态模型和量测模型,尺度4为最细尺度。假定此系统为多传感器多分辨率系统,我们用三个分辨率不同的传感器来观测系统信号,并进行数值仿真,三个传感器对信号观测的结果如图2(a)∼图2(c)所示,其信噪比(SNR)分别为13.4d B、9.2d B、6.7d B,图中为信号的第一个状态分量,初始值x0=[10,0]T和P0=diag(4,4)。这样式(27)的尺度i上的测量值Zk(i)可以通过不同传感器的观测获得,我们选用d B4小波分解和重构,传感器1、2、3的量测分别对应尺度2、3、4上的观测,将三个尺度上的数据融合后的结果图如图2(d)所示,其信噪比为19.3d B,可以看出基于多尺度Kalman的多传感器数据融合滤波算法具有良好的滤波效果,可用于多分辨率多传感器数据融合。
6 结论
本文通过深入分析基于小波变换的动态系统模型,提出一种基于多尺度Kalman的数据融合滤波的方法,利用小波的多尺度特点,把初始估计序列多尺度分解,并进行分层Kalman滤波估计,通过小波重构进行估计融合。该算法将小波多尺度分解和Kalman滤波结合起来,同时还有效地利用了多分辨率的数据信息。通过实验可以看出,该算法对实际中含较强噪声的动态系统状态估计效果较好,同时也能用于多分辨率多传感器数据融合。
摘要：本文通过分析基于小波变换的动态系统模型,提出一种基于小波多尺度的Kalman数据滤波方法,本文利用小波的多尺度特点,把初始估计序列多尺度分解,并在不同尺度层上进行Kalman滤波估计,再利用小波重构来融合各层的估计信息,把标准Kalman滤波只在单一尺度和时间轴上对状态估计值和误差协方差进行数据更新,改进为基于小波变换的尺度轴和时间轴上的双向数据更新,该算法将小波多尺度分解去噪和Kalman滤波相结合,对实际中含较强噪声的动态系统的状态估计效果较好。算法也可用于多分辨率多传感器数据融合。
关键词：多尺度,Kalman滤波,小波变换,数据融合
参考文献
[1]邓自立,许燕.基于Kalman滤波的白噪声估值理论[J].自动化学报,2003,29(1):21-23.DENG Zi-li,XU Yan.White Noise Estimation Theory Based on Kalman Filtering[J].Acta Automatica Sinica,2003,29(1):21-23.
[2]HONG L,CHEN G R,CHUI C K.A Filter-bank-based Kalman Filtering Technique for Wavelet Estimation and Decomposition of Random Signal[J].IEEE Trans.On Circuits and systems II:Analog and Digital Signal Processing,1998,45(2):237-241.
[3]WEN X B,TIAN Z,LIN W.Robust Kalman Filter and Smoothing Recursive Estimator for Multiscale Autoregressive Process[C]//ICSP’04.7th International Conf.Beijing:IEEE,2004,1(1):364-367.
[4]雷明,韩崇昭,元向辉,等.非线性系统的小波分频的扩展Kalman滤波[J].光电工程,2006,33(6):68-72.LEI Ming,HAN Chong-zhao,YUAN Xiang-hui,et al.Wavelet Decomposing based Subbandwise Extended Kalman Filtering for the Nonlinear System[J].Opto-Electronic Engineering,2006,33(6):68-72.
[5]文成林.多传感器单模型动态系统多尺度数据融合[J].电子学报,2001,29(3):341-345.WEN Cheng-lin.Multiscale Data Fusion for Multi-Sensor Single Model Dynamic Systems[J].Acta Electronica Sinica,2001,29(3):341-345.
[6]MALLAT S.A Wavelet Tour of Signal Processing(Second Edition)[M].San Diego:Academic Press,1999.
[7]彭玉华.一种改进的小波变换阈值去噪方法[J].通信学报,2004,25(8):119-123.PENG Yu-hua.An Improved Thresholding Method in Wavelet Transform Domain for Denosing[J].Journal on Communications,2004,25(8):119-123.
基于多尺度特征学习的阴影检测 篇5
阴影是许多图像的组成部分,有时候反映了物体形状、大小等重要信息,然而图像中的阴影会对图像的后续处理造成不利的影响。例如:在智能视频监控中,阴影会严重影响目标的跟踪。阴影检测率会直接影响阴影消除的效果,因此阴影检测也是近些年的研究热点之一。
目前阴影检测方法大体可以分为两类:基于模型[1,2,3]的方法和基于特征[4,5,6]的方法。基于模型的方法利用一定的先验信息建立光照模型来检测阴影;基于特征的方法根据阴影的色度、亮度等特征来识别和检测阴影。基于模型的方法通常需要利用多幅图像的先验信息进行建模,该类方法有一定的局限性,只适合于特定的场景。基于特征的方法不需要场景的先验知识,仅仅通过图像的特征检测阴影,时间复杂度较小。当阴影区被覆盖时,虽然亮度变低却仍然保留着一定的色度信息,依据这一特点,在阴影检测时往往将图像转换为亮度与色度相互独立的颜色空间[7,8]。基于色度的阴影检测方法计算量小,但是对噪声和光照比较敏感。近几年来,基于统计属性[9,10,11,12]的阴影检测方法受到了更高的关注。Zhu等[9]基于图像的可变特征和不可变特征、Lalonde等[10]基于图像的色彩特征和纹理属性、Jiang等[11]基于图像的色彩分割与照明估计、Guo等[12]基于成对区域,这类方法都采用条件随机场标识阴影区域。该类方法对不同的场景与不同光照条件具有较强的适应性,代表了当前的研究水平与未来的发展方向。但是该类方法都需要花费大量精力小心的设计这些特征。Zhu等[9]提取可变特征与不可变特征,这些特征用来训练一个分类器,该方法消耗资源大训练时间长。Lalonde等[10]使用了48维的特征向量,虽然增加特征维数能够达到很好的效果,但是阴影检测的时间大大的增长了。Jiang等[11]对Lalonde等[10]方法进行了改进,取得了更好的效果,却花费了更长的时间。Guo等[12]首先采用均值漂移的方法对图像进行分割,把提取的特征作为支持向量机(SVM)的输入建立分类模型。该方法对局部细节噪声较为敏感并且当场景比较复杂时大大增加了SVM的测试时间。
针对上述方法的缺点,本文提出一种新的多尺度特征学习的阴影检测方法。与上述基于统计属性的方法相比,本文方法是一种在场景中,通过特征学习的方法学习大部分相关特征,这些特征用来进行阴影检测。卷积神经网络CNNs(Convolutional Deep Neural Networks)[13]是人工神经网络的一种,已成为当前图像处理领域的研究热点。它的权值共享网络结构使之更类似于生物神经网络,降低了网络模型的复杂度,减少了权值的数量。它模拟人脑进人脑的机制来解释、处理图像等。因此,本文采用卷积神经网络进行特征学习。一个像素是否属于阴影区像素,不仅与像素自身的亮度、色度等特征有关,还与其周围相邻区域与不相邻区域有紧密的关系,多尺度(图像的拉普拉斯金字塔)的特征学习很好地解决了这个问题。条件随机场CRF(Conditional Random Fields)[14]是一种基于统计的模型。2001年,由John等人首次提出,可在给定全局条件的前提下,计算全局最优输出标记的条件概率。它是一种判别式概率模型,是随机场的一种,常用于标注或分析序列资料。卷积神经网络特征学习之后的后验分布反馈给条件随机场,条件随机场对图像进行标注达到阴影检测的目的。
1 自学习阴影检测框架
给定一副图像,本文在像素级别下检测和定位阴影。对于一个像素来说,不能简单地根据与周围相邻像素之间的关系来判断它是否属于阴影区像素。一个像素与其周围较远处不相邻的像素也有较为密切的关系。人眼能够根据图像中阴影区与其周围非阴影区的不同快速定位阴影。本文对输入的图像经过拉普拉斯金字塔变换,分别以确定的聚类中心为中心进行窗口提取达到考虑像素周围较远处区域的目的。合成训练样本后在卷积神经网络中进行训练,用小的测试集进行测试,最后产生的后验分布反馈给条件随机场,不同标签集合的交集形成阴影检测的结果。本文方法流程如图1所示。
1.1 确定聚类中心
超像素是指许多相似的像素点组合在一起,作为一个整体来处理,这个整体就称之为超像素。近年来,超像素已经被越来越多地应用于图像预处理过程中。超像素利用像素与像素之间特征的相似程度对像素进行分组,从而获取图像的冗余信息,在很大程度上降低了后续图像处理任务的复杂程度。SLIC算法[15]是一种思想简单、效率很高并且运行速度快的算法,通常只需迭代10次就完成聚类。该方法将彩色图像转换为CIELAB颜色空间和XY坐标下的五维特征向量,然后对五维特征向量构造度量标准,对图像像素进行局部聚类的过程。假设图像有N个像素,预分割为K个超像素,则每个超像素的平均像素数目为N/K,聚类中心间的近似距离为。该聚类方法是对k-means聚类算法的一种改进,为提高k-means算法的运算速度,在进行聚类时是以聚类中心的2S×2S正方形区域内搜索相似的像素,而不是在整张图像搜索。本文进行聚类时在聚类中心周围半径为S的区域内进行搜索,如图2所示。
在聚类中心点数目相同时,该方法能够减少重复计算的时间,因为圆形具有严格的几何不变性[16],确定的聚类中心为:
其中,lk、ak与bk表示坐标为(xk,yk)的像素点在CIELAB颜色空间下各通道的值,k∈[1,K]。
1.2 窗口提取与特征学习
输入图像经过拉普拉斯金字塔变换,以每个聚类中心为中心进行窗口提取,则对于坐标为(xk,yk)的像素点,此操作可以描述为:
其中,Ωk为所提取的图像块Ωk∈Ω,W(·)为窗口操作,t为提取窗口的大小,由于所选图像集阴影区域的不平衡性,采用合成少数类过取样算法[17]合成训练样本:
其中,Λ(·)为合成样本操作,Ο为样本集,样本集在卷积神经网络中进行特征学习。本文卷积神经网络结构如图3所示。
对于卷积层来说,上一层的特征图被一个可学习的卷积核进行卷积,然后通过一个激活函数,就可以得到输出特征图。每一个输出特征图可能是组合卷积多个输入的特征图:
其中,xjl和xil-1分别为当前层l和上一层l-1的特征图(下同),k是卷积核,f(·)是非线性激活函数,b为偏置,Mj代表输入特征图的一个选择。对于降采样层来说,有多少个输入特征图,就有多少个输出特征图,只是每个输出特征图都变小了:
其中,down(·)表示降采样函数,kjl为权值,神经元的非线性作用函数为f(x)=(1+e-x)-1。神经元的输出层:
对于式(6)j∈[shadow,n-shadow],卷积神经网络以一个图像块作为输入,经过训练之后得到一个后验分布。在本文中,池操作的步幅等于池邻域的均值,通过降采样层执行池操作有助于学习不变性的特征描述。全连接层工作原理是有一层隐藏层的多层传感器,这个隐藏层后边跟着一个逻辑回归的输出层,逻辑回归的输出层提供类的分布。图像转化为二值图像,二值输出变量的后验分布:
图像块训练之前先进行预处理,采用卷积神经网络可以把概率分布的模型PCNNs(Yi|Οk)定义为:
其中,θ(·)是预处理,F(·)是有5层隐藏层的卷积神经网络。预处理之后,通过在线学习(随机梯度下降),卷积神经网络进行有监督的训练。在训练过程中,梯度通过反向传播的方式计算,交叉熵损失函数被最小化[17]。采用交叉验证的方法确定训练参数,训练样本在训练之前被打乱,这是因为卷积神经网络对未知样本学习得更快。卷积神经网络的初值是从0均值的高斯分布中随机取样来初始化的。
卷积神经网络在训练过程的每个阶段用一个小的验证集来评估正在训练的网络,一旦验证集在执行时不能达到κ步,训练过程就停止(本文κ=5)。在验证集上表现好的网络被用来在图像库上进行测试,最初的学习率是试探性的选择能使误差收敛的最大学习率。
1.3 条件随机场
给定输入图像,阴影检测的任务实际上就是把整幅图像的像素分为了两类,也是对每一个像素进行标签的过程,一个像素要么被标记为阴影像素要么被标记为非阴影像素。这种二分类问题从有监督的特征学习中建立概率估计并提供给条件随机场。条件随机场定义在网格拓扑结构图上,图的节点对应于图像的像素,如式(9),条件随机场是对像素标签最常用的方法之一。由于训练空间(标签图像)的大小增加了计算似然函数梯度的难度,因此条件随机场的参数不能被简单的手工标签阴影区的最大似然估计来得到,采用一种最大化利润学习方法[18]来学习条件随机场的参数。把阴影检测定义为条件概率分布:
其中,ω是该模型的权值,Z(ω)是归一化函数,ν表示节点(比如单个像素),ε表示边(比如相邻像素),X表示整幅图像的像素。对吉布斯能量函数而言,上述分布形式为:
由式(7)与式(8)可知,存在两个隐函数:一元隐函数和二元隐函数。一元函数从卷积神经网络的概率估计得到:
二元函数由类转移隐函数和空间转移隐函数结合得到:
当Yi≠Yj时,否则如下:
其中,α和β通过在每个数据集上交叉验证得到,对随机变量Y∈ΓN,通过采用最大后验估计Y'来确定像素的标签,由于归一化函数Z(ω)不依赖与Y',这种估计转化为能量最小化问题:
根据Y'求标签集合SJ最后求交集:
2 结果与讨论
2.1 实验方案
为了验证提出方法的有效性和高效性,本文实验方案如下:
(1)图像库:本文实验图像库来源于文献[9]图像库(大部分阴影图像背景复杂,多为硬阴影,阴影面积小)和文献[12]图像库(大大部分阴影图像背景简单,多为软阴影,阴影面积大)。
(2)实验参数:本文聚类数目K取150,拉普拉斯金字塔J取5。实验条件为Intel(R)Core(TM)i3-2370@2.40 GHz处理器,2 GB内存,32位操作系统,Visual Studio2012开发平台编程实现。
K的值小于图像大小与样本大小的商,这样能够保证在窗口提取时能够提取每个超像素的大部分像素,从而能够得到更可靠的样本。而当K等于150时,对于图像库中的所有图像,窗口提取都能够满足条件。对于一些特定的实时领域,K的值需要简单的计算与测试确定。实验的精度并不随J的增大而增大,当J=5时,阴影的检出率最高,而当J>5时,阴影的检出率反而下降,因此J取5。
(3)实验方法:宏观上交叉训练,交叉测试;微观上选取9幅场景复杂度不同的图像进行测试。
(4)计算方法:
对于一副有阴影标签的图像来说,假设A表示阴影区像素的集合,B表示非阴影区像素的集合,A'、B'分别表示经过实验得到的阴影区与非阴影区像素的集合。阴影检测效果采用阴影检出率SD、阴影误检率SE来评价,TP计算方法如下(其中N(·)表示计算像素的数目):
2.2 实验结果
(1)宏观测试
本文采用在其中一个图像库上训练然后在另一图像库进行测试的方法,然后比较本文在不同图像集上所能达到的最高精度(阴影检出率SD)。在文献[9]图像库上进行训练,在文献[12]图像库测试的精度达到81.9%;而在文献[12]图像库上训练,在文献[9]图像库上测试的精度为79.8%。这是因为文献[9]图像库的图像场景比较复杂能够训练得到一个更好的网络。本文方法在文献[9]图像库与文献[12]图像库上测试达到的最高精度分别为90.71%与93.37%。由阴影检测计算公式可知:在相同条件下误检率相同时,阴影区的面积越大阴影的检出率相对较高,而文献[12]图像库的阴影区面积较大,因此本文方法在文献[12]图像库上取得较高的精度。
(2)微观测试
为了验证本文方法的场景适应性,选取了9幅场景复杂度不同的图像作为本位阴影检测结果的实例,如图4所示。其中第一排、第三排和第五排为源图像,第二排、第四排和第六排为阴影检测效果图。
定性上,从阴影检测的效果来看,本文方法能够很好地检测软阴影、硬阴影、自阴影以及场景非常复杂的阴影。(a)、(b)、(d)与(f)场景较为简单且无自阴影,但是(a)与(b)阴影多为软阴影,(d)与(f)阴影接近硬阴影;(c)、(d)、(g)、(h)与(i)都包含自阴影,其中(c)与(d)场景较为简单,(g)、(h)与(i)场景较为复杂,(g)与(h)为遥感图像,(i)为复杂的室外自然场景。实验表明,本文提出方法在各种不同场景下的图像的阴影检测效果较好,能够准确检测出由于遮挡产生的自阴影。
定量上,本文对选取的图像进行测试并将本文方法的阴影检出率与误检率和传统的方法进行比较,如表1、表2所示。
从表1与表2可以看出,本文方法阴影检出率较传统的方法高,阴影误检率较传统方法低。而SD与SE并没有完全的相关性,这与SD与SE的计算方法有关。传统方法阴影的检出率随场景的变化变化较大且检出率较低,而本文方法相对较为稳定且能够达到更高的检出率。
在训练资源的消耗上,文献[9]所采用的方法需要更大的内存(9 GB)和更大的时间消耗(大约10小时训练125幅图像)。本文的多尺度特征学习方法只需要2 GB内存,训练时间也更短(大约4小时能训练两个图像库全部的图像)。在实验中当用网络i的权值初始化网络i+1时(1<=i,i+1<=5),本文方法的训练时间能够缩短到2.5小时左右。文献[12]采用的支持向量机方法只适合小数据,卷积神经网络模拟人的可视域可以用多个卷积核映射出多个特征图,更适合图像像素的分类。当图像扩展到三维,或图像场景非常复杂,卷积神经网络的优势就体现出来了。就目前来说,图像越来越复杂,卷积神经网络的适应性越强。
3 结语
本文主要阐述了一种新的阴影检测方法,利用卷积神经网络进行特征学习,结合条件随机场对图像进行标签。该方法场景适应性强,在阴影面积较大、软阴影较多的图像中也取得了很好的效果,而且训练的时间也比现有的基于统计属性的方法大大缩短,阴影检出率也较现有传统的算法高。
多尺度模拟 篇6
Marching Cubes算法(简称MC算法)是由William E.Lorensen和Harvey E.Cline在1987年的Computer Graphics中第一次提出来的,它是面显示算法中的经典算法,也被称为“等值面提取”(Isosurface Extraction)算法,它根据阈值在一个3维数据场中抽取等值面实现3维重建。自从MC算法提出以来,通过不断改进,已广泛地应用在可视化的各个领域,特别是在医学图像三维重建、立体地图测绘等领域。MC算法提出以后,一直就被认为是迄今为止应用最广泛的面重建算法之一。
2. MC算法简介
在医学图像三维重建中,一张CT图像为一个断层切面,两张CT图像组成上下两层切面,每张断层切面取相邻的4个点,共8个点组成了一个立方体(Cube),如图1所示。Marching Cubes算法先单独计算一个立方体与等值面相交的情况,得到若干个三角面,然后再遍历所有的立方体,得到所有的相交三角面,最后将这些三角面拟合为等值面输出。
在原始数据(比如CT、MRI三维图像)中,它们的三维空间数据通常是离散的,沿着立方体边的方向数据场呈连续线性变化,如果立方体的一条边的两个顶点分别大于、小于等值面的值,则在该边上有且仅有一个交点。如果边上的两个顶点都大于或都小于等值面的值,则这条边和等值面没有交点。它的基本算法是逐个处理数据场中的立方体,找出等值面与这些立方体的交点,进而找出包含在这些体元中的等值面。
Lorensen在他们的论文中第一次提出MC算法时,具体分析了Cube与等值面相交的各种情况:因为一个立方体有8个顶点,每个顶点有2种状态,即表面内和表面外,因而Cube与表面的关系有2的8次方,即256种情况,通过枚举这256种情况,可得到一个状态表以得到等值面与Cube各边相交的信息以及Cube内表面的拓扑结构。因为Cube是立方体,利用2种不同的对称性,即表面与Cube上每个顶点与表面的关系同时逆转,不会影响Cube内三角化表面的拓扑结构,也就是说256种情况是两两互补的,所以只需要考虑0-4种顶点在表面内的情况。另外,根据立方体3个方向旋转上的对称性,即Cube向3个方向任意旋转90度的倍数,Cube内表面的三角形拓扑结构不变。这样经过分析,最终将256种情况减少为15种。这15种情况如图2所示:
Marching Cubes具体的算法步骤:
(1)根据对称关系构建一个256种相交关系的索引表。该表指明等值面与体素的哪条边相交。
(2)取相邻两层图片中相邻的8个象素,构成一个体素并把这8个象素编号。
(3)比较每个象素与阈值确定该象素是l还是0。
(4)把这8个象素构成的01串组成一个8位的索引值。
(5)用索引值在上边的索引表里查找对应关系,并求出与立方体每条边的交点。
(6)交点构成三角形面片。
(7)遍历三维图像的所有体素,重复执行(2)到(6)。
(8)将得到的所有的三角面片输出。
3. MC算法的缺点
Marching Cube算法的优点是明显的,但它也存几个缺点:
(1)MC算法局限于通过阈值分割来提取等值面,而医学图像由于成像机理不同,所适用的分割方法也不同,这就造成了MC算法适用的局限性;
(2)MC算法中关于等值面与立方体交点的计算是采用插值法的,计算量比较大;
(3)MC算法存在二义性。如果在体元的一个边界面上,标号为0和标号为1的角点分别位于对角线的两端,则该表面上存在四个等值点,那么就会有两种可能的等值线的连接方式,这就是二义性;
(4)MC算法的输出是三角形网格,由于得到的数据分布十分密集,由此建立的三维网格数量庞大。
4. MC算法在不同尺度上的改进
根据上面的分析,MC算法由于计算量非常大,所以显示三维面时会适成延时。本文根据MC算法可以选择不同的像素点组成立方体进行求交点提出了一种多尺度的MC算法。其原理如下:
设每相隔N个点取一个顶点,上下两层共8个顶点组成一个大间距的立方体,从这个大间距立方体开始计算。下图3为同一个方向上的两个顶点a1,a2。求这个两个顶点组成的边与等值面的交点K。
根据MC算法求交点的公式,K的坐标为:
K=(a2-a1)(阈值-v1)/(v2-v1)(1)
遍历完所有大间距的立方体后,就得到了一个大尺度的交点集,由这个交点集所组成的三角面片会明显减少,计算量下降,使得三维图像绘制速度加快。
计算完上述大尺度的立方体后,缩小顶点之间的间距,将立方体的边长变为原大尺度的一半。下图4为缩小一半间距所同一方向上的三个顶点。
则以a1,a3组成边的两个顶点,a3,a2组成另一组边的两个顶点。则这两条边和等值面的交点求法如下:
K1=(a3-a1)(阈值-v1)/(v 3-v1)(2)
K2=(a2-a3)(阈值-v3)/(v 2-v3)(3)
观察(2)式和(3)式有
K=(a2-a1)(阈值-v1)/(v2-v1)
K1=(a3-a1)(阈值-v1)/(v3-v1)
其中a3-a1=(a2-a1)/2
所以有:
K(v2-v1)=(a2-a1)(阈值-v1)
K1(v3-v1)=(a2-a1)(阈值-v1)/2
即
K1(v3-v1)=K(v2-v1)/2
K1=K(v2-v1)/2(v3-v1)
同样K2=K(v2-v1)/2(v2-v3)
根所上述理论,计算一次小尺度的交点K1,K2,都可以利用到上一次大尺度的计算结果K,而v2,v1,v3都是已知值,所以上述计算过程将会非常之简单。这样就从理论上找到了一个可以提高MC算法计算速度的方法。以下表1是本文根据这个理论所做的实验,实验的数据是来自田捷的中科院自动化研究所医学影像研究室提供的一个老鼠的CT图,该数据的大小是512层,每层是512×512。每个像素点是32bit的值。
由表1可以看出,采用不同尺度的方向后,生成的三解形面片少了,速度也快了。对比图5和图6可以看出,用两个尺度相差数倍组成立方体的MC算法计算出的三维图像的效果差别是不明显的。
5. 结论
本文根据MC算法计算量大,难以实时绘图的缺点,提出了一个多尺度的MC算法,经过实验对比,新的算法比原算法速度更快,生成三维面的效率更高,而且生成三维面的质量并没有明显的下降。
参考文献
[1]徐毅,李晓梅.对体可视化Marching Cube算法的改进[J].计算机工程,1999,25(11):52—54.
[2]谢小棉,李树祥.基于MC的医学三维等值面的平滑与归并[J].中国图像图形学报,2001,6(8):806—809.
[3]William E.Lorensen,Harvey E.Cline,Marching Cubes:A High Resolution3D Surface Construction Algorithm,1987
[4]Doi A,Koide A.An Efficient Method of Triangulating Equi-Valued Surfaces hy Using Tetrahedral Cells,1988.
多尺度模拟 篇7
本文立足于提高基于HOG的人体检测方法的性能,并增强存在遮挡情况下的检测率,提出了一种基于多尺度方向直方图的特征提取方法提取图像的粗特征和精特征。在此基础上使用级联SVM模型和Adaboost算法进行检测,形成一种多尺度多分类器的人体检测模型,并对其效能进行测试。
人体检测作为计算机视觉研究应用领域的重要研究方向之一,有着巨大的市场应用潜力和研究意义,可广泛运用于人机交互、自动控制、视频检索和智能视频监控等方面[1]。这也进一步推动了学者对视频帧和图像中的人体目标检测,但迄今为止人体检测仍然存在着效率低下且易受遮挡、视角、光照等因素影响,这些方面的问题成为人体检测领域的难点亟待解决[2]。
在目前的研究成果中,主要有三类人体检测方法,分别是模板匹配、基于Haar特征的方法和基于梯度方向直方图HOG(Histograms of oriented gradients)的方法[3,4,5]。其中模板匹配是提取图像并与已有人体形状模板进行匹配检测是否是人体,较为典型的有头肩模型用于匹配检测[6],但人体形态的多样性导致这种方法效果较差,已较少使用;基于Haar特征的方法使用Haar特征和分类器材识别,由于人脸特征能与Haar特征较好契合从而在人脸检测中得到较好应用,但因背景颜色、衣服色泽等因素影响,在人体检测中的效果不如预期[7];基于HOG的方法经描述局部图像的方向分布情况和梯度获取该区域的HOG作为检测特征并结合训练的分类器进行识别,该方法精度较高因此得到了较多的重视和研究[8]。但基于HOG的也存在训练速度较慢、特征窗口提取单一等缺陷,因此许多研究对其进行改进,如多尺度HOG方法[9]、局部二值模式结合HOG方法[10]、加权块[11]等方法来增强局部特征的表达能力从而提高基于HOG方法的效率。但这类方法都是在固定尺度上计算,没有将人体的多尺度轮廓纳入特征提取中,而多尺度特征表达的提出通过使用Gabor小波、灰度协方差、组合块等方法提高了基于HOG方法的性能,进一步增强了基于HOG方法的应用价值。
1 本文人体检测方法
本文在已有研究的基础上,综合考虑复杂背景下相互遮挡且不易检测的难点,以背景区分、多视角检测、遮挡去除等方法用于人体目标检测。训练特征块以多尺度方法从训练样本中提取和选择,包含粗特征和精特征,其中粗特征使用扩展多尺度方向特征提取并使用粗分类器选择,精特征使用WTA hash编码的扩展多尺度方向直方图提取并使用贪婪算法选择,对特征集使用多层级联的Adaboost算法和SVM分类器进行多视角分类器进行分类训练。将训练好的分类器运用于待检测图像,并结合复杂背景处理、多视角样本生成和遮挡处理得到最后的检测结果,算法流程如图1所示。
2 特征提取
本文分别使用扩展多尺度方向特征提取和WTA hash编码的多尺度梯度方向直方图特征提取方法获得图像的粗特征和精特征。
2.1 HOG原理
HOG是按照图像的外观、形状等特征的边缘方向和梯度来获取图像特征,这种特征并没有从图像整体去提取特征,而是将图像分割成许多个名为cell的小细胞单元,再计算所有cell中的梯度方向直方图。同时为了提高算法效率将多个cell融合成一个名为block的大单元并将其梯度进行归一化,计算梯度方向向量[12]。
在HOG特征计算过程中,首先使用一维离散微分模板从垂直和水平两个方向将180°分成N个方向计算内部各个cell的梯度,研究得出,当N为9时梯度精度最大[13]。然后取各cell单元中所有像素点的梯度幅值进行加权求得其比重,依照比重为直方图通道投票加权。由于图像背景的存在和影响,各cell的梯度变化幅度较大,分类器难以适应其变化,因此再进行归一化处理,使用L2范数为因子,计算如式(1)所示。
式中e为常数,V表示归一化前的梯度向量值。这样可求得图像的HOG特征。
2.2 扩展多尺度方向特征提取
多尺度方向MSO(Multi-scale orientation)特征是从水平和垂直两个方向综合描述直立或者姿态变化小的人体轮廓信息,但却不能有效识别姿态变化明显的人体轮廓信息[14]。为此设计一种扩展多尺度方向特征EMSO(Extended MSO)提取方法,以扩展的Haar-like特征为基础,在MSO的基础上再提取两个特征块,使人体轮廓信息特征提取的姿态变化适应性更强,其原理如图2所示。
EMSO的方法是将图像分解为n个单元块,再将各单元块划分为4个子块并计算各子块颜色积分值,常规MSO已有垂直梯度Dd和左旋45°梯度Ddl,而水平梯度Dl和右旋45°梯度Ddr需计算得出,Dl计算方法为:
Ddr计算方法为:
式中I(·)是该点的像素值。计算4个子块的颜色积分值后,可得EMSO单元块方向He的值为:
在此基础上可得量化块方向Fi的值为:
即将180°以内方向的连续值进行8值离散化处理,将所有Fi进行连接后可得EMSO特征,该特征可用方块表征任何矩形,其方向计算可在方向上进行,可降低EMSO特征维数,且不易受光照、形变影响,因此可提高运算效率。
2.3 扩展梯度方向直方图提取
设存在图像,其尺寸为(W,H),对应的多尺度图像金字塔度数为K=3,其中各度在k=4个相互独立级的量化将使用的尺度因子。在这种条件下,各级图像将被分为8×8的cell,在2×2的重叠block上提取HOG特征,各级HOG级联可得多尺度梯度方向直方图MHOG(Multi-scale Histograms of oriented gradients)特征[15],对应的维数为4644,维数值较大,计算复杂度较高,因此使用WTA hash编码予以改进,即WM-HOG。
图像上的X和Y两点相似性在某种程度上和对对应的特征维数排序类似,其pairwise-order量度计算方法为:
式中,xi、xj分别是X和Y对应的第i特征维数,T(·)是阈值函数,其公式为:
式(6)可变换为:
式中有:
WTA hash编码[16]如图3所示,图3(a)中L(X,0)={1,2,3,5}对应式(10)的低于0下标元素的元素下标集合,类似的L(Y,0)={3,4,5},代入式(9)可得R0(X,Y)=|{3,5}|=2,由于维数为6,因此K=4,θ={1,4,2,5,0,3},图3(a)和(b)的X不相关,编码结果分别为1和2,即不相等。而图3中的(c)、(d)是对(a)做线性运算而得,因此对应的输出结果也为1。通过实例可知,WTA hash编码可抗干扰,且具备线性稳定性,能将MHOG转换到对特征维数不敏感的空间中得到稀疏特征,避免原特征空间中的冗余度大的缺陷,提高运算速度和运算效率。
WTA hash编码的步骤如图4所示,实际上是将向量X中的值随机排列,通过滑动窗口法扫描,将最大值保留,其他值设置为0,同时利用随机数组θ消除向量中元素的相关性,这就是本文所用的WMHOG特征提取方法。
图5所示为图2的例图经MHOG和WMHOG两种特征提取方法后的效果对比,分析发现,WMHOG特征比MHOG特征更为稀疏,降低了约17%。
3 特征选择
对EMSO提取的粗特征块的选择使用Adaboost算法[17],对WMHOG提取的精特征块使用贪婪算法[18]。
使用Adaboost算法选择粗特征先需构建弱分类器,由于人体的轮廓普遍存在一定的角度范围,对应的方向特征fn在人体轮廓的表征上具有较好的性能。因此对弱分类器的需求是在特征选择时将处于一定角度范围内的特征作为正例样本,这是选择Adaboost算法的主要原因。弱分类器的构建方法为:
式中,C(n)是粗特征经粗分类器分类之后的结果,fn是通过式(1)所得的EMSO方向特征。Tmin(n)、Tmax(n)为弱分类器的最小和最大角度阈值,阈值可使用贪婪算法获取。
WMHOG特征块选择使用贪婪算法,在其第t次求解过程中,由于局部最优解的优势,新加入的特征ft+1可进一步增大分类精度,整个求解过程是从待选特征集A中选出局部最优解的特征值加入到已选集S中,即:
式中,I(u;C)=∑Cp(C,f)log(p(C|f))为f的条件熵,p(C|f)为C样例正确的概率。这样经t次求解后,A、S的更新如下:
即从集合A中去除ft+1加入到集合S中,持续这一过程直到分类精度不增长或增长幅度小于一定阈值时为止。
4 目标检测
4.1 多视角检测
分类器训练使用级联Adaboost和线性SVM,共分四层,第一层全部视角样本共同训练单一模型,第二、三层分别使用3个、5个视角,前三层使用Adaboost,第四层使用SVM,其框架如图6所示。在人体目标检测时,对图像窗口进行逐一分类以判定图像中是否存在目标人体,并对其位置进行确定,当4层分类器都通过时判定存在目标人体,反之判定不存在目标人体。图6的分类器框架经过实验验证确定前三层只用3个级联就能实现较高的查全率和较快的检测速率,对第四层的SVM分类器,在超过128个精特征块可使SVM分类精度超过97%。
但是这种多视角检测中存在一定的问题,即人体的正面(0°)和侧面(90°)的样本易获取,但是获得中间角度(30°、45°、60°)的样本较困难,对此采取特征块位置线性平移的方法合成中间角度样本。由于人体的对称性完成0°到90°的样本获取即可,设特征块n的视角样本度数为D,对应的位置为XD(n),则可得0°,90°内的样本线性合成为:
式中α(n)为合成系数,由于不进行精确的姿态识别,因此可假定人体为圆柱体,可得α(n)为合成系数α(n)为:
类似地,块特征也可由式(16)获取:
式中FD(n)表示特征块n在D°的特征值,这样可实现通过正面样本及侧面样本获取中间各视角的样本从而完成分类器训练。
4.2 复杂背景处理
图像人体识别过程中复杂背景处理主要解决人体轮廓与背景图层灰度相似情况下检测目标出现空洞的问题[15]。因此在复杂背景存在的情况下,定义函数f(a,b):
式(17)中a(x,y)是当前帧在像素点(x,y)的像素值,b(x,y)是对应点的背景像素值,这样根据式(17)可知:当(x,y)是背景点时,f(a,b)=0,且f(a,b)分散分布;当(x,y)为人体目标区域内点时,f(a,b)≥0且f(a,b)集中分布。计算f(a,b)集中分布的范围即为目标人体区域范围。
4.3 遮挡处理
解决了复杂背景问题,还需解决其他目标或其他人体对人体目标遮挡的问题,在遮挡存在的情况下会使人体目标局部轮廓缺失导致分类错误[17],在此使用特征重装方法。先设是被遮挡的特征之一,可由式(18)合成:
式中f0m(n)是训练样本集中样本m的第n特征块,M是正例样本个数,wm是样本m的权值向量,wm的计算方法为:
式(19)是限定与遮挡人体相似度高的样本在样本重装时权值较大的基础上得出,fu(n)是需分类目标的未被遮挡特征,fum(n)是第m个样本的第n特征块且该特征块未被遮挡,最后对其进行归一化从而完成被遮挡样本块样本的重装再进行分类。
5 实验及分析
为对本文的人体检测算法进行验证,使用法国国家信息与自动化研究所的公共测试集INRIA行人数据库进行实验。该数据库是目前最常用的静态行人检测数据库,提供了原始的图片和相应的标注文件。训练时采用了2000正例样本和1000反例样本,反例样本主要是建筑、树干等和人体轮廓相似的圆柱形物体,其中正面和侧面视角的样本由人工划分形成,中间视角(30°、45°、60°)的样本能过多视角样本合成算法生成。最后使用测试集INRIA数据库进行实验,部分人体检测结果如图7所示。可以看出本文算法可以检测出大部分人体,图7(a)-(g)包含了复杂背景、与类人形圆柱体较为贴近、非正面或侧面多视角人体、多姿态人体等多种情况,本文算法都能成功检测,这是由于使用了复杂背景处理和特征重装成功检测出了中间复杂背景中有遮挡的人像。但在图7(h)和图7(i)中漏检了后方复杂背景下的人体目标,分析发现是WTA hash编码导致部分有用信息被删除从而无法检测出复杂背景下较小像素的遮挡人体目标。
算法性能定量分析比较使用单位图片误检/图片FPPI(False Positive per Image)条件下的漏检率,分别与常规的HOG特征提取检测[9]、HOG和局部二值模式HOG+LBP(HOG+Local Binary Pattern)相结合的方法[10]及多尺度HOG(MHOG)特征方法[11]相比较,其结果如图8所示。从图8可知,在相同FPPI的条件下本文的多尺度多视角算法的漏检率低于HOG、HOG+LBP、MHOG等方法,但在FPPI接近1的时候,本文算法的漏检率与MHOG逐渐接近,经分析这是由于WTA hash编码后部分有用信息被滤除导致漏检率增加。
上述各种算法在复杂背景下的检测效率如表1所示。由表1可知,本文的多尺度多视检测算法的检测率比HOG等算法高,检测时间比HOG等算法少得多,这主要是因为使用WTA hash编码后滤除图像特征中的冗余信息使图像特征变得较为稀疏,从而使后续算法较为迅速,降低了检测时间。
6 结语