Research of Visual H∞ Control of Inverted Pendulum With Time-varying Computational Time and Computational Error
-
摘要: 针对基于视觉传感的倒立摆实时控制系统中，通过每帧图像计算小车位移和摆杆偏角产生的时变计算时间和计算误差直接影响控制系统性能甚至导致系统失稳问题，不同于目前不考虑计算时间和计算误差或将计算时间视为定时滞进行研究倒立摆视觉实时控制方法，本文首先创新设计了新型的倒立摆视觉伺服控制实验平台，然后设计了一种基于事件触发机制的工业图像采集策略，提出了考虑小车和摆杆特征的小车位移和摆杆偏角计算方法，并统计分析了图像处理计算时间和计算误差特性；进一步建立了融合计算时间和计算误差的闭环控制系统模型，理论证明系统的稳定性并建立了反映计算时间与系统稳定性能之间关系，给出了系统对计算误差的H∞扰动抑制水平γ和控制器设计方法.最后，仿真和实时控制实验验证了所提方法可行且有效.Abstract: In visual-feedback-based real-time control of an inverted pendulum on a cart, some new challenges are as follows. On the one hand to get the values of cart position and pendulum angle from images needs more computation time. On the other hand, the values obtained have bigger errors than those in the traditional sensor approach. These will affect the system performance and even lead the system to instability. In this paper, we present an analytical study and experiments on:1) a novel platform; 2) event-triggered sampling and signal processing of the values from images, and statistical analysis of time-varying computational time and error; 3) modelling of a visual-based feedback closed-loop system; 4) conditions for system stability in terms of the statistics as described before; and 5) H∞ norm based controller design. At the end simulation and experiment are conducted to validate the design method.
-
Key words:
- Inverted pendulum /
- visual servoing /
- image processing /
- computational time /
- computational error
-
星载合成孔径雷达 (Synthetic aperture radar, SAR) 具有运行周期固定、地面覆盖范围广、分辨率高等特点, 因此被应用于地球表面各种资源的探测.星载SAR图像被广泛应用于制图、土地监测、灾害监测、地质、水文、海洋、冰川等许多方面[1]. SAR图像中建筑物轮廓的抽取是一个难点, 其分割精度直接影响着后续解译的质量.
目前, SAR图像的分割方法主要有聚类方法[2-4]、阈值方法[5-7]、基于马尔科夫随机场的方法[8-12]和水平集方法[13-23]等.马秀丽等[2]结合分水岭和谱聚类算法实现了对SAR图像的分割.邓晓政等[3]提出了一种基于非负矩阵分解的谱聚类SAR图像分割方法.徐海霞等[4]针对谱聚类方法应用SAR图像分割时Laplace矩阵的特征值和特征向量难以计算的问题, 结合SAR图像的多尺度统计信息, 提出了基于谱聚类和混合模型的SAR图像分割方法.薛景浩等[5]在Rayleigh分布假设下, 提出一种最小误差的阈值化分割算法.聚类方法和阈值方法难以有效地利用SAR图像的边缘信息和空间信息, 当目标外观有歧义时, 这类算法的鲁棒性有待提高.余航等[8]针对SAR图像的特点, 提出了一种基于上下文分析的无监督分层迭代算法, 该算法结合了聚类算法和区域增长算法的优点, 提高了SAR图像分割的准确率. Deng等[9]采用MRF方法实现了SAR图像的非监督分割.傅兴玉等[11]在传统的MRF方法基础上, 提出了一种基于像素Gabor纹理相似度的邻域势函数模型, 改进了高分辨率SAR图像建筑物的分割效果.然而, MRF方法计算量巨大, 得到的分割结果受到相干斑噪声影响严重.
与上述分割方法相比, 水平集方法更能够适应图像的拓扑关系变化, 具有较强的实用性. Chan-Vese (CV) 模型[13]以分片光滑的区域描述图像的目标和背景, 对光学图像取得了较好的分割结果.为有效分割复杂背景下的感兴趣区域, Bresson等[14]基于待分割目标轮廓的先验知识, 在模型中增加了先验形状约束. Chen等[15]通过深度Boltzman机提取先验形状的分层结构, 在变分框架下通过耦合优化Boltzman机参数和变形曲线实现感兴趣目标的分割.田昊等[16]将多个先验形状竞争模型引入水平集方法中, 在标记函数的指导下, 利用先验形状能量来约束曲线的演化, 在对图像进行分割的同时完成建筑物的检测和提取.然而, 由于SAR发射的是纯相干波, 相干波经过与地物的相干作用, 特别是地物的后向散射作用, 目标的回波信号产生了衰减, 在图像中表现为相干斑噪声, 为了提高分割精度, 需要在模型中考虑噪声对分割结果的影响.涂松等[17]综合分析了基于活动轮廓模型的SAR图像分割方法. Mao等[18]使用复小波变换获取SAR图像的边缘信息, 然后使用CV模型抽取图像中的水边线. Zhang等[19]基于支持向量机获取SAR图像目标的初始轮廓, 再使用几何活动轮廓模型 (Geometric active contour, GAC) 实现SAR图像的分割. Ben等[20]对SAR图像基于Gamma分布建模, 通过水平集方法分割出SAR图像中的目标. Shuai等[21]在此基础上提出了一种静态全局最小能量函数分割模型, 改善了分割结果. Sui等[22]针对CV模型的缺点, 将CV模型中的灰度均值项改为Gamma概率分布项, 有效地降低了相干斑噪声对分割结果的影响.
上述文献提出的分割方法基于SAR图像的低层特征, 当目标周围背景复杂或在目标被部分遮挡的情况下会出现一定的分割误差.苏娟等[24]提出了一种基于SAR和可见光图像融合的SAR图像建筑物分割算法, 该方法充分利用可见光图像解译性好的特点, 提高了分割精度.本文结合高分辨率光学遥感影像, 提出了一种基于能量最小化的星载SAR图像建筑物分割模型:基于Gamma分布拟合SAR图像数据以构造条件概率能量项推动演化曲线向目标边缘演化, 定义演化曲线长度能量项以保证曲线在变形过程中的平滑; 根据先验知识构造先验形状能量项, 演化曲线在该能量项的推动下不断变形, 有效克服了目标被遮挡或复杂背景对分割结果造成的影响.实验结果表明, 该方法能够有效地提高复杂背景下SAR图像中目标轮廓的分割精度.
1. 星载SAR图像分割方法
1.1 分割模型算法框架
本文提出的分割算法流程如图 1所示.首先, 基于Gamma分布拟合SAR图像数据, 定义几何活动轮廓模型的曲线演化速度, 基于水平集方法对SAR图像分割得到粗分割结果; 再利用高分辨率光学遥感图像提取建筑物目标的先验形状, 以粗分割结果作为初始的轮廓线, 在先验形状约束下优化能量泛函直至收敛.在每一步迭代过程中, 需要将当前的演化曲线与先验形状进行配准.
1.2 基于能量最小化的SAR图像分割模型
由于SAR图像相干斑噪声的存在, 传统的混合Gaussian模型难以描述SAR图像的灰度分布.本文采用Gamma分布拟合图像数据. Gamma分布的概率密度函数 (Probability density function, PDF) 为[25]
$ \begin{align} \label{eq1} P_{\Gamma}^L(I\vert R) = \frac{L^L}{\mu _R(L - 1)!}\left( {\frac{I}{\mu _R }} \right)^{L - 1}{\rm e}^{ - \frac{LI}{\mu _R }} \end{align} $
(1) 其中, $I$ 为SAR图像的灰度, $R$ 为待分割SAR图像的目标和背景区域, $L$ 为SAR图像的等效视数, 本文取等效视数 $L = 1$ , $\mu _R $ 为区域的平均灰度值.
根据SAR图像数据的概率分布构造条件概率能量项
$ \begin{align} \label{eq2} E_{\rm con\_prob} = \sum\limits_{i = 1}^2 {\int_{R_i } { - \log P_{\Gamma}^L(I\vert R_i ){\rm d}\Omega } } \end{align} $
(2) 式中, $R_1 $ 和 $R_2 $ 分别为SAR图像的目标及背景区域, $\Omega \in {\bf R}^2$ 为图像的定义域.条件能量项的大小由目标和背景区域中SAR数据的概率分布给出, 该能量项推动曲线朝着目标的边界演化.
为使曲线在变形过程中保持平滑定义长度能量项
$ \begin{align} \label{eq3} E_{\rm length} = \int_\gamma {\rm d}s \end{align} $
(3) 式中, $\gamma $ 为当前演化曲线, 长度能量项的大小为演化曲线的长度.
由于SAR图像的信噪比低, 加之目标背景复杂, 仅依据条件概率能量项和长度能量项难以驱动演化曲线获取清晰、完整的轮廓.为了提高分割精度, 本文从高分辨率光学遥感图像中根据建筑物目标的先验形状构造先验能量项约束曲线的演化, 使得抽取的建筑物轮廓与先验形状相似.将先验形状以水平集方式表示, 设为 $\phi _s $ , 定义先验形状能量项
$ \begin{align} \label{eq4} E_{\rm prior} = \int_\Omega {(H(\phi(x, y)) - H(\tilde {\phi }_s(x_s, y_s )))^2} {\rm d}\Omega \end{align} $
(4) 式中, $H$ 为Heaviside函数, 当 $\phi \ge 0$ 时, $H(\phi ) = 1$ , 当 $\phi < 0$ 时, $H(\phi ) = 0$ . $\phi(x, y)$ 为当前演化曲线的水平集函数, $\tilde {\phi }_s(x_s, y_s )$ 为经仿射变换后的先验形状的水平集函数, $\tilde {\phi }_s(x_s, y_s )$ 定义为
$ \begin{align} \label{eq5} \tilde {\phi }_s(x_s, y_s ) = \phi _s(T^{ - 1}(x, y)) \end{align} $
(5) 式中, $T^{ - 1}$ 为演化曲线到先验形状之间的仿射变换, 考虑形状变化的几种可能情形:平移、缩放、旋转, 本文定义仿射变换参数 $P = \left[{t_x t_y s_x s_y \theta } \right]^{\rm T}$ , 其中, $t_x $ 为水平位移, $t_y $ 为垂直位移, $s_x $ 为水平尺度变换, $s_y $ 为垂直尺度变换, $\theta $ 为旋转角度. $(x, y)$ 和 $(x_s, y_s )$ 满足下面关系:
$ \begin{align} \left[{\begin{array}{l} x_s \\ y_s \\ \end{array}} \right] = &\ T^{ - 1}\left[{\begin{array}{l} x \\ y \\ \end{array}} \right] =\notag\\[2mm] & \left[{\begin{array}{l} {\kern 1pt} \dfrac{1}{s_x }[(x-x_g-t_x )\cos (\theta ) +\\ \qquad\quad(y-y_g - t_y )\sin(\theta )] + x_g \\[1mm] \dfrac{1}{s_y }[-(x-x_g-t_x )\sin(\theta ) +\\ \qquad\quad(y -y_g - t_y )\cos (\theta )] + y_g \\ \end{array}} \right] \end{align} $
(6) 式中, $(x_g, y_g )$ 为先验形状的中心点坐标.
为获取仿射变换参数, 使用梯度下降流极小化 $E_{\rm prior}$ 得到下式:
$ \begin{align} \label{eq7}\left\{ \begin{gathered} \frac{{{\rm d}{t_x}}}{{{\rm d}t}} = \iint_\Omega {(H({{\tilde \phi }_s}) - H(\phi )) \times 2\delta({{\tilde \phi }_s})} \times \\ \qquad\quad \left(\frac{{{{\tilde \phi }_{sx}}}}{{{s_x}}}\cos (\theta ) - \frac{{{{\tilde \phi }_{sy}}}}{{{s_y}}}{\rm sin}(\theta )\right){\rm d}x{\rm d}y \\[2mm] \frac{{{\rm d}{t_y}}}{{{\rm d}t}} = \iint_\Omega {(H({{\tilde \phi }_s}) - H(\phi )) \times 2\delta({{\tilde \phi }_s})} \times \\ \qquad\quad \left(\frac{{{{\tilde \phi }_{sx}}}}{{{s_x}}}\sin(\theta ) + \frac{{{{\tilde \phi }_{sy}}}}{{{s_y}}}\cos(\theta )\right){\rm d}x{\rm d}y \\[2mm] \frac{{{\rm d}{s_x}}}{{{\rm d}t}} = \iint_\Omega {(H({{\tilde \phi }_s}) - H(\phi ))} \times \\ \qquad\quad 2\delta({{\tilde \phi }_s}) \times \frac{{{{\tilde \phi }_{sx}}}}{{{s_x}}} \times {x_s}{\rm d}x{\rm d}y \\[2mm] \frac{{{\rm d}{s_y}}}{{{\rm d}t}} = \iint_\Omega {(H({{\tilde \phi }_s}) - H(\phi ))} \times \\ \qquad\quad 2\delta({{\tilde \phi }_s}) \times \frac{{{{\tilde \phi }_{sy}}}}{{{s_y}}} \times {y_s}{\rm d}x{\rm d}y \\[2mm] \frac{{{\rm d}\theta }}{{{\rm d}t}} = \iint_\Omega {(H({{\tilde \phi }_s}) - H(\phi )) \times 2\delta({{\tilde \phi }_s})} \times \\ \qquad\quad \left(\frac{{{{\tilde \phi }_{sy}}{s_x}}}{{{s_y}}} \times {x_s} - \frac{{{{\tilde \phi }_{sx}}{s_y}}}{{{s_x}}}{y_s}\right){\rm d}x{\rm d}y \end{gathered} \right. \end{align} $
(7) 式中, ${\tilde \phi _{sx}} =\dfrac {{\partial {{\tilde \phi }_s}}}{{\partial x}}$ , ${\tilde \phi _{sy}} =\dfrac{{\partial {{\tilde \phi }_s}}}{{\partial y}}$ .
演化曲线在条件概率能量项、长度能量项和先验形状能量项的作用下逐渐向目标轮廓演化, 模型的能量泛函为
$ \begin{align} \label{eq8} E = {\mu _1 E_{\rm con\_pro} + \mu _2 E_{\rm length} + \mu _3 E_{\rm prior} } \end{align} $
(8) 式中, $\mu _1 $ , $\mu _2 $ 和 $\mu _3 $ 分别为条件概率能量项、长度能量项和先验能量项权重系数. $\mu _1 $ 越大, 演化曲线在条件概率能量的驱动下向目标边界演化的速度越快; 增大 $\mu _2 $ , 可保证曲线演化平滑; 当SAR图像的信噪比较低或建筑物背景复杂时, 可增大 $\mu _3 $ , 演化曲线在先验能量的作用下将向与先验形状经仿射变换后相似的轮廓逼近.
计算式 (8) 的一阶变分, 基于梯度下降流得到曲线的演化方程
$ \begin{align} \frac{\partial \phi }{\partial t} =&\ \mu _1V\left\| {\nabla \phi } \right\| + \mu _2 \kappa \left\| {\nabla \phi } \right\| + \notag\\[1mm] &\ 2\mu _3(H(\phi ) - H(\tilde {\phi }_s ))\delta(\phi ) \end{align} $
(9) 式中, $V = -( {\log \hat {\mu }_{R_1 } + {I}/{\hat {\mu }_{R_1 } } - \log \hat {\mu }_{R_2 } - {I}/{\hat {\mu }_{R_2 } }} )$ , $\hat {\mu }_{R_1 } $ 和 $\hat {\mu }_{R_2 } $ 分别为当前演化曲线所划分的内外区域的灰度均值, $\delta $ 为Dirac函数, $\delta(\phi ) = {H}'(\phi )$ , $\kappa $ 为演化曲线的曲率.
为提高数值计算的效率, 避免迭代过程中重新初始化水平集函数, 本文使用有符号距离约束[26]
$ \begin{align} \label{eq10} E_{\rm Disreg} = \int_\Omega {\frac{1}{2}}(\left\| {\nabla \phi } \right\| - 1)^2{\rm d}\Omega \end{align} $
(10) 则式 (9) 变为
$ \begin{align} \frac{\partial \phi }{\partial t} =&\ \mu _1 V\left\| {\nabla \phi } \right\| + \mu _2 \kappa \left\| {\nabla \phi } \right\| +\notag\\[1mm] &\ 2\mu _3(H(\phi ) - H(\tilde {\phi }_s ))\delta(\phi ) +\\[1mm] &\ \mu _4 \left[\Delta \phi-{\rm div}\left(\frac{\nabla \phi }{\left\| {\nabla \phi } \right\|}\right)\right] \end{align} $
(11) 在SAR图像分割过程中, 为提高算法的执行效率, 首先令 $\mu _3 = 0$ , 利用条件概率能量和长度能量驱动演化曲线进行初步分割, 得到待分割目标的初步轮廓, 之后加入先验形状能量项进行分割, 从而得到最终的分割结果.
1.3 基于能量最小化的SAR建筑物分割算法步骤
基于能量最小化的SAR建筑物分割算法主要步骤如下:
步骤1.基于高分辨率光学图像获取感兴趣目标的轮廓, 作为目标的先验形状 $\phi _s $ .
步骤2.初始化:演化曲线的水平集函数为 $\phi ^0$ , 算法第一阶段最大迭代次数为 $Max\_Iter1$ , $n = 0$ , 设置模型各项的权重系数 $\mu _1$ , $\mu _2 $ , $\mu _3 $ 和 $\mu _4 $ .
步骤3.曲线进行第一阶段演化, 令式 (11) 中 $\mu _3= 0$ , 迭代:
while $n \le Max\_Iter1$ do
$ \begin{align*} &{\phi ^{n + 1}} = {\phi ^n} + \Delta t \times \\[1mm] &\qquad \Bigg\{ {\mu _1}V\left\| {\nabla {\phi ^n}} \right\| + {\mu _2}\kappa \left\| {\nabla {\phi ^n}} \right\| + \\[1mm] &\qquad {\mu _4}\left[\Delta {\phi ^n}-{\rm div}\left(\frac{{\nabla {\phi ^n}}}{{\left\| {\nabla {\phi ^n}} \right\|}}\right)\right]\Bigg\} \\[1mm] & n = n + 1 \end{align*} $
零水平曲线若停止演化, 转步骤4, 否则继续循环;
end
将输出结果作为粗分割结果 $\phi '$ .
步骤4.取 $\phi ^0=\phi '$ , 令算法第二阶段最大迭代次数为 $ Max\_Iter2$ , $n = 0$ , 仿射变换最大迭代次数为 $M$ , $m = 0$ ;
while $n \le Max\_Iter2$ do
初始化仿射变换的参数: $t_x = 0$ , $t_y = 0$ , $s_x$ $=$ $1$ , $s_y= 1$ , $\theta = 0$ ;
while $m \le M$ do
$ \begin{align*} &t_x^{m + 1} = t_x^m + \Delta t \times\\ &\qquad \Bigg\{ (H(\tilde \phi _s^m) - H({\phi ^n})) \times2\delta(\tilde \phi _s^m) \times\\ &\qquad \left(\frac{{\tilde \phi _{sx}^m}}{{s_x^m}}\cos({\theta ^m}) - \frac{{\tilde \phi _{sy}^m}}{{s_y^m}}\sin({\theta ^m})\right)\Bigg\}\\[1mm] &t_y^{m + 1} = t_y^m + \Delta t \times \\ &\qquad \Bigg\{(H(\tilde \phi _s^m) - H({\phi ^n})) \times 2\delta(\tilde \phi _s^m) \times \\ &\qquad \left(\frac{{\tilde \phi _{sx}^m}}{{s_x^m}}\sin({\theta ^m}) + \frac{{\tilde \phi _{sy}^m}}{{s_y^m}}\cos({\theta ^m})\right)\Bigg\}\\[1mm] &s_x^{m + 1} =s_x^m + \Delta t \times\\ &\qquad \Bigg\{ \left(H(\tilde \phi _s^m) - H({\phi ^n})\right) \times \\ &\qquad 2\delta(\tilde \phi _s^m) \times \frac{{\tilde \phi _{sx}^m}}{{s_x^m}} \cdot x_s^m\Bigg\}\\[1mm] & s_y^{m + 1} = s_y^m + \Delta t \times\\ &\qquad \Bigg\{ \left(H(\tilde \phi _s^m) - H({\phi ^n})\right) \times \\ &\qquad 2\delta(\tilde \phi _s^m) \times \frac{{\tilde \phi _{sy}^m}}{{s_y^m}} \times y_s^m\Bigg\}\\[1mm] &{\theta ^{m + 1}} = {\theta ^m} + \Delta t \times\\ &\qquad \Bigg\{ \left(H(\tilde \phi_s^m) - H({\phi ^n})\right) \cdot 2\delta(\tilde \phi _s^m) \times \\ &\qquad \left(\frac{{\tilde \phi _{sy}^ms_x^m}}{{s_y^m}} \times x_s^m - \frac{{\tilde \phi _{sx}^ms_y^m}}{{s_x^m}} \times y_s^m\right)\Bigg\}\\[1mm] & m = m + 1 \end{align*} $
end
$ \begin{align*} &{\phi ^{n + 1}} = {\phi ^n} + \Delta t \times \Bigg\{ {\mu _1}V\left\| {\nabla {\phi ^n}} \right\| + \\ &\qquad 2{\mu _3}\left(H({\phi ^n}) - H({\tilde \phi _s})\right)\delta({\phi ^n}) + {\mu _2}\kappa \left\| {\nabla {\phi ^n}} \right\| +\\ &\qquad {\mu _4}\left[\Delta {\phi ^n}-{\rm div}(\frac{{\nabla {\phi ^n}}}{{\left\| {\nabla {\phi ^n}} \right\|}})\right]\Bigg\} \end{align*} $
$n = n + 1$
零水平曲线若停止演化, 转步骤5, 否则继续循环;
end.
步骤5.输出最终分割结果 $\phi $
1.4 数值计算
离散化式 (11) 得到
$ \begin{align} \phi _{ij}^{n + 1} =&\ \phi _{ij}^n + \Delta t \times \bigg\{ {\mu _1}(\max(V_{ij}^n, 0){\nabla ^ + } + \\[1mm] &\ \min(V_{ij}^n, 0){\nabla ^ - }) + {\mu _2}\kappa _{ij}^n({(D_{ij}^{0x})^2} + \notag\\[1mm] & \ {(D_{ij}^{0y})^2})^{\frac{1}{2}} +2{\mu _3}\left( {H(\phi _{ij}^n) - H(\tilde \phi _{{s_{ij}}}^n)} \right)\times\notag\\[1mm] & \ \delta(\phi _{ij}^n) +{\mu _4}\left( {\Delta \phi _{ij}^n - \kappa _{ij}^n} \right)\bigg\} \end{align} $
(12) 式中, $\Delta t$ 为时间步长, $n$ 为迭代次数, $i$ , $j$ 为图像坐标, $\nabla ^ + $ 和 $\nabla ^ - $ 定义如下:
$ \begin{align} {\nabla ^ + } =&(\max {(D_{ij}^{ - x}, 0)^2} + \min {(D_{ij}^{ + x}, 0)^2} + \\[1mm] &\ \max {(D_{ij}^{ - y}, 0)^2} + \min {(D_{ij}^{ + y}, 0)^2}{)^{\frac{1}{2}}} \notag\\[2mm] {\nabla ^ - } =&(\max {(D_{ij}^{ + x}, 0)^2} + \min {(D_{ij}^{ - x}, 0)^2} +\\[1mm] & \ \max {(D_{ij}^{ + y}, 0)^2} + \min {(D_{ij}^{ - y}, 0)^2}{)^{\frac{1}{2}}} \end{align} $
(13) 其中, $D_{ij}^{ + x} $ , $D_{ij}^{ - x} $ , $D_{ij}^{0x} $ 和 $D_{ij}^{ + y} $ , $D_{ij}^{ - y} $ , $D_{ij}^{0y} $ 分别为 $x$ 和 $y$ 的前向差分、后向差分和中心差分格式.
曲率 $\kappa $ 采用中心差分格式离散化:
$ \begin{align} \label{eq14} \kappa _{ij}^k = \frac{D_{ij}^{xx}(D_{ij}^{0y} )^2 + 2D_{ij}^{0x} D_{ij}^{0y} D_{ij}^{xy} - D_{ij}^{yy}(D_{ij}^{0x} )^2}{\left( {(D_{ij}^{0x} )^2 +(D_{ij}^{0y} )^2} \right)^{\frac{3}{2}}} \end{align} $
(14) 水平集函数 $\phi $ 的拉普拉斯算子 $\Delta \phi _{ij}^n $ 的差分格式为
$ \begin{align} \label{eq15} \Delta \phi _{ij}^n = \phi _{i - 1, j}^n + \phi _{i + 1, j}^n + \phi _{ij - 1}^n + \phi _{ij + 1}^n - 4\phi _{ij}^n \end{align} $
(15) 在数值计算中, Heaviside函数和Dirac函数采用下面的形式进行计算:
$ \begin{align} &{H_\varepsilon }(z) =\begin{cases} 1, &z > \varepsilon \\ 0, & z < - \varepsilon \\ \dfrac{1}{2}\left(1 + \dfrac{2}{\pi }\arctan \dfrac{z}{\varepsilon }\right), &\left| z \right| \leq \varepsilon \end{cases} \notag\\[2mm] & {\delta _\varepsilon }(z) = \begin{cases} 0, &\left| { z} \right| > \varepsilon \\[1mm] \dfrac{1}{{2\varepsilon }}\left(1 + \cos \dfrac{{\pi z}}{\varepsilon }\right), & \left| { z} \right| \leq \varepsilon \end{cases} \end{align} $
(16) 2. 实验结果及分析
为了验证本文方法的有效性, 分别在带斑点噪声的合成图像和真实的SAR图像上进行分割实验, 并与CV方法[13]、Ben方法[20](算法实现的代码通过作者主页提供的链接下载, 参数设置通过实验进行了优化) 和Sui方法[22]所得到的分割结果进行比较.采用Dice (Dice similarity coefficient, DSC) 系数[27]作为指标对分割结果进行定量评价, DSC定义如下:
$ \begin{align} \label{eq17} DSC = \frac{2\times \vert A \cap B\vert }{\vert A\vert + \vert B\vert }\times 100 \% \end{align} $
(17) 式中, $A$ 为自动分割结果, $B$ 为手工分割结果. $DSC$ 越大, 分割效果越好, 当 $DSC$ 为100%时, $A$ 和 $B$ 完全匹配, 达到最佳分割效果.
2.1 带斑点噪声的合成图像分割
图 2(a)为使用Matlab的imnoise函数生成的带斑点噪声的合成图像, 斑点噪声的均值为1, 方差为0.04, 图像大小为100像素 $\times$ 100像素, 正方形目标右下角被遮挡; 图 2(b)为目标的先验形状.
图 3(a)中蓝色的初始曲线在条件概率能量和长度能量的共同作用下逐渐向内收缩, 红色曲线为迭代250次之后演化曲线的位置.当条件概率能量和长度能量达到最小时, 分割曲线停止了演化. 图 3(a)中的绿色曲线为模型没有加入先验形状约束的分割结果.为获取目标的完整轮廓, 在模型中引入由先验形状构造的先验形状能量项, 权重参数取 $\mu _1$ $=$ $1$ , $\mu _2 = 1$ , $\mu _3 = 1.5$ , $\mu _4 = 0.1$ , 时间步长 $\Delta t$ $=$ $0.5$ , 由于受到先验能量的作用, 图 3(e)中右下角处的曲线逐渐向外扩张并向先验形状逼近, 绿色曲线为迭代45次的结果, 迭代190次后, 模型最终获得了目标的较为完整的轮廓, 如图 3(e)中的红色曲线所示. 图 3(f)为手工分割的结果, 与本文模型相比, 图 3(b)的CV方法、图 3(c)的Ben方法和图 3(d)的Sui方法的演化曲线向目标的边缘演化, 可以有效地将目标和背景分割出来, 但是由于出现了遮挡, 这些模型均未能获取完整的目标轮廓.
为了研究图像噪声对本文方法的影响, 选取斑点噪声方差分别为0.1, 0.2, 0.3和0.4的图像进行分割实验, 结果如图 4所示.
从图 4可以看出, 本文方法在不同方差斑点噪声的影响下仍然能够分割出较为完整目标轮廓.但随着图像噪声方差的逐渐增大, 条件概率能量项受到较大的影响, 演化曲线难以完整分割出目标右下角的缺失部分, 得到的目标轮廓精度降低.
2.2 真实SAR图像建筑物分割
本文选取三幅真实SAR图像进行分割实验, 实验数据来自TerraSAR-X于2008年10月采集的上海交通大学区域的单视斜距图像 (Single look slant range complex, SSC) 及IKONOS 2007年3月采集的上海交通大学区域图像.如图 5(a)所示, SAR图像大小为106像素 $\times$ 145像素, 从对应的分辨率为1 m的IKONOS全色波段图像 (图 5(b)) 中抽取" H "状楼的轮廓作为先验形状 (图 5(c)).由于在分割过程中需要将变形曲线与先验形状配准, 因此不需要预先配准SAR和IKONOS图像.
图 5 " H "状楼的SAR图像和全色波段图像及其先验形状Fig. 5 SAR image of the building with " H " contour, its corresponding panchromatic image and the prior shape图 6为不同方法的分割结果.由于SAR图像存在大量强散射斑块, 在同质区域内表现为部分像素灰度值较高, 部分像素灰度值较低; 建筑物内部纹理复杂且轮廓不完整. CV方法和Ben方法只能分割出图像中的一些高亮线条, 无法得到完整的轮廓, 分割结果中存在较多孤立像素. Gamma分布能够较好地反映SAR图像数据的统计特性, Sui模型结合了Gamma分布和CV模型, 因此Sui模型的分割结果优于前两种方法.图 6(e)为手工标注的分割结果, 本文模型引入了建筑物的先验形状信息, 模型参数设置为 $\mu _1 = 1$ , $\mu _2 = 1$ , $\mu _3 = 2$ , $\mu _4 = 1$ , 时间步长 $\Delta t = 0.1$ , 最终分割结果由SAR图像建筑物的灰度分布特性及先验形状确定, 与Sui模型相比, 抽取的轮廓更为完整, 分割精度更高.
图 6 不同方法对" H "字楼的分割结果Fig. 6 Segmentation results of the building with " H " contour by different methods第二组实验SAR图像如图 7(a)所示, 图像大小为136像素 $\times$ 73像素, 从对应的IKONOS全色波段图像 (图 7(b)) 中抽取目标建筑物的轮廓作为先验形状 (图 7(c)).
图 7 " L "状楼的SAR图像和全色波段图像及其先验形状Fig. 7 SAR image of the building with " L " contour, its corresponding panchromatic image and the prior shape图 8为不同方法的分割结果.由于SAR图像相干斑噪声严重, 目标建筑物在SAR图像中无明显轮廓, CV方法和Ben方法只能分割出图像中的小块的分块区域. Sui方法优于前两种方法, 可以得到较为完整的形状, 但未能有效地将目标与背景分割开来, 且存在着明显的误分割现象.本文模型的参数设置为 $\mu _1 = 1$ , $\mu _2 = 1$ , $\mu _3 = 2$ , $\mu _4 = 0.5$ , 时间步长 $\Delta t$ $=$ $0.1$ , 在先验形状的约束下分割得到的目标建筑物轮廓更加精确完整.
图 8 不同方法对" L "状楼的分割结果Fig. 8 Segmentation results of the building with " L " contour by different methods为进一步验证本文方法在复杂背景下的分割效果, 选择了体育馆作为待分割的目标如图 9(a)所示.图像的大小为221像素 $\times$ 201像素; 从对应的IKONOS全色波段图像 (图 9(b)) 中抽取建筑物的轮廓作为先验形状 (图 9(c)).
图 9 体育馆的SAR图像和全色波段图像及其先验形状Fig. 9 The SAR image, its corresponding panchromatic image and the prior shape of the gymnasium图 10为不同方法的分割结果.由于SAR图像中建筑物右上方边缘断裂, CV方法和Ben方法得到的轮廓在该处出现凹陷现象.加之建筑物内部纹理结构复杂, Ben方法得到的轮廓内部存在较多的孔洞. Sui方法虽然优于前两种方法, 但由于该建筑物周围背景较为复杂, Sui方法未能有效地将目标与背景分割开来, 存在大量误分割现象.本文模型的参数设置为 $\mu _1 = 1$ , $\mu _2 = 1$ , $\mu _3 = 2$ , $\mu _4 = 1$ , 时间步长 $\Delta t = 0.1$ , 在先验形状的约束下分割得到的建筑物边界更为清晰平滑.
分割实验的精度指标DSC及运行时间 (运行环境: Intel Core I7, 4 GB RAM; Matlab 2012 a, 操作系统为Windows 7) 如表 1所示.从表中可以看出, 在SAR仿真图像和真实图像分割实验中, 本文方法的分割精度接近或超过90%, 优于其他方法.引入待分割目标的先验信息, 显著提高了分割精度.在算法的时间效率上, 本文方法优于CV方法, 但低于Ben和Sui方法.综合考虑分割精度和运行时间, 本文的方法具有更强的分割性能.
表 1 实验结果比较Table 1 Comparison of segmentation results实验名称 CV方法 Ben方法 Sui方法 本文方法 精度 (%) 耗时 (s) 精度 (%) 耗时 (s) 精度 (%) 耗时 (s) 精度 (%) 耗时 (s) 合成图像 97.4 15.4 97.4 21.8 97.3 16.9 99.4 113.7 " H "状楼 28.7 133.9 57.3 82.9 84.0 52.9 89.6 116.1 " L "状楼 13.4 162.3 21.7 20.8 62.9 13.7 91.5 67.1 体育馆 82.0 326.2 76.1 88.0 86.3 101.9 94.1 206.1 平均值 55.4 159.5 63.1 53.3 82.6 46.4 93.7 125.8 3. 结论
根据星载SAR图像的特点, 本文提出了一种基于能量最小化的星载SAR图像建筑物分割方法, 能够有效地提取建筑物目标的完整轮廓.通过星载SAR图像的灰度特征和GAC模型, 获取初步的分割结果, 然后从高分辨率光学图像提取感兴趣建筑物的先验形状以构造先验能量, 约束变形曲线直至收敛到目标的边界.模型中的参数通过实验设置, 如何自适应地选取权重系数以及根据SAR图像的特点采用更合适的分布模型以进一步提高建筑物目标的分割精度是下一步研究的方向.
-
图 1 倒立摆视觉伺服控制实验平台结构
Fig. 1 Scheme of inverted pendulum visual servoing control platform
图 2 倒立摆视觉伺服控制系统控制框图
Fig. 2 Control block diagram of inverted pendulum visual servoing control system
图 16 在控制律${K_1}$下的小车位移和摆杆偏角(仿真)
Fig. 16 Cart position and pendulum angle in controller ${K_1}$ (simulation)
图 17 在控制律${K_2}$下的小车位移和摆杆偏角(仿真)
Fig. 17 Cart position and pendulum angle in controller ${K_2}$ (simulation)
图 18 在控制律${K_1}$下的小车位移和摆杆偏角(实时控制)
Fig. 18 Cart position and pendulum angle in controller ${K_1}$ (real-time control)
图 19 不同图像处理计算时间上界$\bar d$的小车位移
Fig. 19 Cart position in different upper bound of image processing computational time $\bar d$
图 20 不同图像处理计算时间上界$\bar d$的摆杆偏角
Fig. 20 Pendulum angle in different upper bound of image processing computational time $\bar d$
表 1 相机内部参数
Table 1 Intrinsic parameters of camera
参数 取值 $f$ 909.5 pixels $d$ 0.9263 m ${\alpha _u}$ $1.6241\times 106 $ pixels/m ${\alpha _v}$ $1.6241\times106 $ pixels/m $s$ 0 ${u_0}$ 334.85 pixels ${v_0}$ 267.46 pixels 
下载: 导出CSV
-
[1] 王瑶为, 邢科新, 马剑, 张文安.直线一级倒立摆的自抗扰控制方法及实现.控制工程, 2017, 24(4):002 http://d.old.wanfangdata.com.cn/Periodical/jczdh201704002Wang Yao-Wei, Xing Ke-Xin, Ma Jian, Zhang Wen-An. Implementation and design of active disturbance rejection control for the linear inverted pendulum. Control Engineering of China, 2017, 24(4):002 http://d.old.wanfangdata.com.cn/Periodical/jczdh201704002 [2] Ronquillo-Lomeli G, Ríos-Moreno G J, Gómez-Espinosa A, Morales-Hernández L A, Perea M T. Nonlinear identification of inverted pendulum system using Volterra polynomials. Mechanics Based Design of Structures and Machines, 2016, 44(1-2):5-15 doi: 10.1080/15397734.2015.1028551 [3] Li Z J, Zhang Y N. Robust adaptive motion/force control for wheeled inverted pendulums. Automatica, 2010, 46(8):1346-1353 doi: 10.1016/j.automatica.2010.05.015 [4] Ravichandran M T, Mahindrakar A D. Robust stabilization of a class of underactuated mechanical systems using time scaling and Lyapunov redesign. IEEE Transactions on Industrial Electronics, 2011, 58(9):4299-4313 doi: 10.1109/TIE.2010.2102318 [5] 王忠杰, 谢璐璐.信息物理融合系统研究综述.自动化学报, 2011, 37 (10):1157-1166 http://www.aas.net.cn/CN/abstract/abstract17604.shtmlWang Zhong-Jie, Xie Lu-Lu. Cyber-physical systems:a survey. Acta Automatica Sinica, 2016, 37(10):1157-1166 http://www.aas.net.cn/CN/abstract/abstract17604.shtml [6] Bradley J M, Atkins E M. Toward continuous state-space regulation of coupled cyber-physical systems. Proceedings of the IEEE, 2012, 100(1):60-74 doi: 10.1109/JPROC.2011.2161239 [7] Ge X H, Yang F W, Han Q L. Distributed networked control systems:a brief overview. Information Sciences, 2017, 380:117-131 doi: 10.1016/j.ins.2015.07.047 [8] 张永立, 程会锋, 李洪兴.三级倒立摆的自动摆起与稳定控制.控制理论与应用, 2011, 28(1):37-45 http://d.old.wanfangdata.com.cn/Periodical/kzllyyy201101005Zhang Yong-Li, Cheng Hui-Feng, Li Hong-Xing. The swing-up and stabilization of the triple inverted pendulum. Control Theory and Applications, 2011, 28(1):37-45 http://d.old.wanfangdata.com.cn/Periodical/kzllyyy201101005 [9] 李雪冰, 马莉, 丁世宏.一类新的二阶滑模控制方法及其在倒立摆控制中的应用.自动化学报, 2015, 41(1):193-202 http://www.aas.net.cn/CN/abstract/abstract18598.shtmlLi Xue-Bing, Ma Li, Ding Shi-Hong. A new second-order sliding mode control and its application to inverted pendulum. Acta Automatica Sinica, 2015, 41(1):193-202 http://www.aas.net.cn/CN/abstract/abstract18598.shtml [10] 武玉强, 朱成龙.车轨长度受限的并行双摆能量控制.控制理论与应用, 2015, 32(9):1254-1260 http://d.old.wanfangdata.com.cn/Periodical/kzllyyy201509015Wu Yu-Qiang, Zhu Cheng-Long. Energy control for parallel-type double inverted pendulums with restricted cart rail length. Control Theory and Applications, 2015, 32(9):1254-1260 http://d.old.wanfangdata.com.cn/Periodical/kzllyyy201509015 [11] Muralidharan V, Mahindrakar A D. Position stabilization and waypoint tracking control of mobile inverted pendulum robot. IEEE Transactions on Control Systems Technology, 2014, 22(6):2360-2367 doi: 10.1109/TCST.2014.2300171 [12] Huang S H, Pan Y C. Automated visual inspection in the semiconductor industry:a survey. Computers in Industry, 2015, 66:1-10 doi: 10.1016/j.compind.2014.10.006 [13] 王耀南, 陈铁健, 贺振东, 吴成中.智能制造装备视觉检测控制方法综述.控制理论与应用, 2015, 32(3):273-286 http://d.old.wanfangdata.com.cn/Periodical/kzllyyy201503001Wang Yao-Nan, Chen Tie-Jian, He Zhen-Dong, Wu Cheng-Zhong. Review on the machine vision measurement and control technology for intelligent manufacturing equipment. Control Theory and Applications, 2015, 32(3):273-286 http://d.old.wanfangdata.com.cn/Periodical/kzllyyy201503001 [14] Riggio M, Sandak J, Franke S. Application of imaging techniques for detection of defects, damage and decay in timber structures on-site. Construction and Building Materials, 2015, 101:1241-1252 doi: 10.1016/j.conbuildmat.2015.06.065 [15] Gao J, Proctor A A, Shi Y, Bradley C. Hierarchical model predictive image-based visual servoing of underwater vehicles with adaptive neural network dynamic control. IEEE Transactions on Cybernetics, 2016, 46(10):2323-2334 doi: 10.1109/TCYB.2015.2475376 [16] Ye W Q, Li Z J, Yang C G, Sun J J, Su C Y, Lu R Q. Vision-based human tracking control of a wheeled inverted pendulum robot. IEEE Transactions on Cybernetics, 2016, 46(11):2423-2434 doi: 10.1109/TCYB.2015.2478154 [17] Wang H P, Vasseur C, Koncar V, Chamroo A, Christov N. Modelling and trajectory tracking control of a 2-DOF vision based inverted pendulum. Journal of Control Engineering and Applied Informatics, 2010, 12(3):59-66 [18] Magana M E, Holzapfel F. Fuzzy-logic control of an inverted pendulum with vision feedback. IEEE Transactions on Education, 1998, 41(2):165-170 doi: 10.1109/13.669727 [19] Tu Y W, Ho M T. Design and implementation of robust visual servoing control of an inverted pendulum with an FPGA-based image co-processor. Mechatronics, 2011, 21(7):1170-1182 doi: 10.1016/j.mechatronics.2011.07.011 [20] Kizir S, Ocak H, Bingul Z, Oysu C. Time delay compensated vision based stabilization control of an inverted pendulum. International Journal of Innovative Computing, Information and Control, 2012, 8(12):8133-8145 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=cb3b1b8e6f3b7b2de73fe7983a59202d [21] Benitez-Morales A, Santos O, Romero H, Ramos-Velasco L E. Suboptimal robust linear visual servoing for a delayed underactuated system. Optimal Control Applications and Methods, 2013, 34(6):696-711 doi: 10.1002/oca.v34.6 [22] Van Hamme D, Goeman W, Veelaert P, Philips W. Robust monocular visual odometry for road vehicles using uncertain perspective projection. EURASIP Journal on Image and Video Processing, 2015, 2015:10 doi: 10.1186/s13640-015-0065-6 [23] Canny J. A computational approach to edge detection. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1986, 8(6):679-698 http://d.old.wanfangdata.com.cn/OAPaper/oai_doaj-articles_3d39d0b11988c5f90bf44b10d764f020 [24] Liu M Q, Zhang X G, Zhang Y L, Lyu S. Calibration algorithm of mobile robot vision camera. International Journal of Precision Engineering and Manufacturing, 2016, 17(1):51-57 doi: 10.1007/s12541-016-0007-y [25] Xu Z Z, Shin B S, Klette R. Accurate and robust line segment extraction using minimum entropy with Hough transform. IEEE Transactions on Image Processing, 2015, 24(3):813-822 doi: 10.1109/TIP.2014.2387020 [26] Xu S Y, Lam J, Zou Y. New results on delay-dependent robust H∞ control for systems with time-varying delays. Automatica, 2006, 42(2):343-348 doi: 10.1016/j.automatica.2005.09.013 [27] Gao H J, Wu J L, Shi P. Robust sampled-data H∞ control with stochastic sampling. Automatica, 2009, 45(7):1729-1736 doi: 10.1016/j.automatica.2009.03.004 [28] Yan H C, Qian F F, Zhang H, Yang F W, Guo G. H∞ fault detection for networked mechanical Spring-Mass systems with incomplete information. IEEE Transactions on Industrial Electronics, 2016, 63(9):5622-5631 doi: 10.1109/TIE.2016.2559454 [29] Jiang X S, Tian X M, Zhang T L, Zhang W H. Quadratic stabilizability and H∞ control of linear discrete-time stochastic uncertain systems. Asian Journal of Control, 2017, 19(1):35-46 doi: 10.1002/asjc.v19.1 [30] Han Q L. Absolute stability of time-delay systems with sector-bounded nonlinearity. Automatica, 2005, 41(12):2171-2176 doi: 10.1016/j.automatica.2005.08.005 [31] Peng C, Tian Y C. Delay-dependent robust stability criteria for uncertain systems with interval time-varying delay. Journal of Computational and Applied Mathematics, 2008, 214(2):480-494 doi: 10.1016/j.cam.2007.03.009 期刊类型引用(2)
1. 樊新航，徐建军，王怿瑾，刘俊生，董春. 回旋式倒立摆的H∞控制. 测控技术. 2024(05): 93-101 . 
百度学术2. 韩治国，陈能祥，许锦，李伟. 基于BP神经网络的倒立摆智能控制研究. 实验技术与管理. 2021(06): 101-106 . 百度学术其他类型引用(6)
-
下载:
计量
- 文章访问数: 2416
- HTML全文浏览量: 266
- PDF下载量: 580
- 被引次数: 8
