自适应迭代算法

自适应迭代算法（共10篇）

自适应迭代算法 篇1

摘要：医学图像分割是医学图像处理分析领域中的研究重点和难点问题,文中对医学CT图像的三维分割方法进行了深入研究,提出了一种医学CT图像的三维分割框架——三维自适应迭代分割算法(SO3DAISA)。试验结果表明,本文的分割方法在很大程度上减少了人工干预、执行效率高、图像分割效果好,并且具有很好的实用性。

关键词：医学图像,自适应迭代分割,图像分割,实用性

计算机断层扫描(CT)数据及核磁共振图像(MRI)数据的可视化如今已经作为医学图像处理中研究的热点问题之一。实现医学图像三维重建的首要前提就是对图像数据进行正确、合理的分割,从中提取出感兴趣器官、组织或病变体的三维重建,达到辅助治疗与手术规划的目的[1]。

CT、MRI等图像数据,又称三维医学图像,对其进行分割是在三维空间进行。在医学图像的获取过程中,由于影像设备中各电子器件的随机扰动和受周围环境的影响,使图像多少含有噪声和失真,这就给医学图像准确、快速分割带来了一定的难度,影响了组织的分割与提取[2]。由于传统分割方法的阈值设置不合理时,会产生过分割和欠分割问题,从而导致算法的鲁棒性不高。因此如何设置最佳阈值,如何检测和避免过分割和欠分割是三维医学图像在实际应用中遇到的难题[3]。

文中针对三维图像分割难题提出了一种基于3D区域生长的三维迭代自适应分割算法(SO3DAISA),该算法主要包括4部分:种子选择、生长准则的设计、自适应迭代分割算法设计。

1 种子的选择

区域生长从种子开始,因此种子选取优劣对最后的分割结果有很大影响。最理想的种子选择方法是全自动方法,由计算机自动从输入图像中提取最佳种子,但医学图像结构复杂、数据量大,实现比较困难,因此目前主要还是采用人工交互方法[4]。

为尽可能减少和方便人工交互,文中设计了一个简单的种子设置软件。通过该软件可以方便地从 3D 图像中抽取任何一张切片Sz(x,y), 通过鼠标操作,可以用线条、矩形或圆形等形状在切片图像中标识出种子候选区域,如图1所示。基于这些候选区域,提出了一种基于概率统计的种子选择方法。

(1)从3D图像中抽取多张切片(通常采取每隔 3、4 张抽一张的方法),利用灰度直方图的分布计算概率得出种子的候选区域,并从切片图像中标识出种子的候选区域。假设共有N个候选区域Ri,1≤ i≤N,通常 N>6;

(2)对每个区域R中的像素进行灰度值统计,由式(1)和式(2)计算其均值和方差,分别记为$\overline{X},\sigma {}_i^2$,其中1≤i≤N。

$\overline{X}{}_i=\frac{1}{n}{\displaystyle \sum _{k=1}^{n}X}{}_k$ (1)

$\sigma {}_i^2=\frac{1}{n-1}{\displaystyle \sum _{k=1}^n(}X{}_k-\overline{X}){}^2$ (2)

其中,n 为区域Ri,1≤i≤N中的像素数,Xk,1≤ K≤n代表区域中第K个元素的灰度值;

(3)对每一个候选区域Ri,1≤i≤N,依据正态分布“3σ原则”,视$\left|X{}_i-\overline{X}{}_i\right|＞3\sigma$的样本点为不合理种子点,予以剔除;

(4)对剔除不合理候选点后的区域重新计算均值和方差;

(5)舍弃具有最大和最小方差的区域,计算所有合理种子点的均值$\overline{X}$和方差σ2;

(6)将所有种子点的二维坐标映射为三维坐标,并加入种子队列中;

(7)得出所需要的三维种子侯选队列。

2 生长准则的设计

种子选择后,要根据种子的特征进行生长。通常,合理的生长准则既要考虑目标对象的总体特征,也要考虑其局部特征。文中设计了一种生长准则,它综合了对象的全局和局部信息,具有较好的性能[5]。

全局特征$F{}_G=\exp \left(-\frac{({\rm I}{}_V-\overline{X}{}_R){}^2}{2\sigma {}_R^2}\right)(3)$

其中,IV表示体素V的灰度值,R表示目标区域,$\overline{X}{}_R$表示区域R中体素的平均灰度值,σ${}_R^2$表示区域R中体素灰度值方差。当新的体素添加到区域R中,应更新$\overline{X}{}_R$和σ${}_R^2$,因此,FG准则反映了当前体素V与区域R的相似性大小,是一种全局的相似性评价准则。

局部特征$F{}_L=\{\begin{array}{l}1\text{i}\text{f}{\rm N}{}_R＞\lambda {\rm N}\\ \lambda {\rm N}{}_R/{\rm N}\text{o}\text{t}\text{h}\text{e}\text{r}\text{w}\text{i}\text{s}\text{e}\end{array}(4)$

其中N和NR分别表示当前体素邻域大小和邻域中已划分到区域中体素的个数。λ为一调节参数,该准则反映了当前体素V与其近邻体素的相似性大小,是一种局部相似性评价准则。

为综合考虑全局和局部的特征信息,采用加权因子W将FG和FL的线性组合作为最终的生长准则F,F的取值越大,表明当前体素属于目标区域的可能性越大。

F=WFG+(1-W)FL (5)

加权因子W的取值是制约生长准则F可靠性的一个关键的因素。从式(5)可以看出,FG越大,表明当前体素与目标区域的相似性越大,考虑到区域生长的邻接性,此时当前体素归并到目标区域的可能性应该较大,因此FG占主导地位,其权值W应较大;相反,FG越小,表明当前体素与目标区域在灰度值上差别越大,此时应考虑其局部特征FL,所以(1-W)应该越大,即W应取较小值。因此可以将W看作是的函数,在文中,采用式(6)函数来表示两者的关系。

$W=1-\exp \left(-\frac{F{}_G^2}{2\sigma {}^2}\right)0\le F{}_G\le 1$ (6)

其中σ2为控制参数。由式(5)和式(6)推导得出分割阈值

$F=\left[1-\exp \left(-\frac{F{}_G^2}{2\sigma {}^2}\right)\right]F{}_G+\exp \left(-\frac{F{}_G^2}{2\sigma {}^2}\right)F{}_L$ (7)

3 三维自适应迭代分割算法的实现

根据式(7)可以计算出图像的分割阈值,算法主要有3个控制参数:Δ、η和δ。其中,Δ代表生长准则阈值TF的减小步长;δ是用户设置参数,用于控制分割的精度,该参数通常为一常数,在迭代分割过程中不变。

假设Δ和TE在第k次迭代时的取值为Δ(k)和TE(k),则第k+1次迭代Δ和TE的相应值Δ(k +1)和TE(k+ 1)分别由式(8)和式(9)确定。

TE(k+1)=ηN-VTF(k)+β (9)

当有新的体素归并到区域R时,其均值$\overline{X}{}_R$和方差σ${}_R^2$由式(10)和式(11)计算

$\overline{X}{}_{R\text{n}\text{e}\text{w}}=(\overline{X}{}_{R\text{o}\text{l}\text{d}}\times {\rm N}{}_R+V{}_{\text{o}\text{x}\text{e}\text{l}}V{}_{\text{a}\text{l}})/({\rm N}{}_R+1)$ (10)

$\sigma {}_{R\text{n}\text{e}\text{w}}^2=((\sigma {}_{R\text{o}\text{l}\text{d}}^2+\overline{X}{}_{R\text{o}\text{l}\text{d}})\times {\rm N}{}_R+V{}_{\text{o}\text{x}\text{e}\text{l}}V{}_{\text{a}\text{l}{}^2})/({\rm N}{}_R+1)-\overline{X}{}_{R\text{n}\text{e}\text{w}}^2(11)$

算法重复迭代执行,直到(TF(k)-TF(k+1))<σ时算法终止。

参数$\overline{X}{}_R$和σ${}_R^2$有两种更新方式:第一种方式为逐元素更新方式,在这种方式下,每当有一个新体素归并到目标区域时,就对$\overline{X}{}_R$和σ${}_R^2$更新一次;第二种方式是成批更新方式,在一次迭代下,将所有种子的邻域判决完毕后再对$\overline{X}{}_R$和σ${}_R^2$进行更新。

现以一张256×256的CT图为分割对象,令计数k=0,计算$\overline{X}(0)$和方差σ2(0),初始化FG,FL;设置TF的初始减小步长Δ(0)以及η,初始化掩模AuxVMask,将种子置1,其它体素置0;执行第一次区域生长。

当迭代一次后,执行 k=k+1;更新$\overline{X}(k),\sigma {}^2(k),F{}_G(k),F{}_L(k)$和TF(k);如果(VTF(k)-VTF(k-1))≥TE或越界,重新恢复到前一次迭代状态,令Δ=Δ/2;如果Δ<δ,程序执行完毕,退出。

分割结果对比,如图2所示。

在该算法中,每次迭代总是在上一次迭代的基础上向外扩展,并采用生长过度判断准则防止生长过度,具有较好的鲁棒性。通过逐次迭代得到针对不同三维图象的局部分割阈值,最大程度提取出医学图象中的细节信息。根据上图实验对比,自适应迭代分割比别的一般算法具有更好的分割效果。

4 结束语

利用三维自适应迭代分割算法,很好地解决了医学图像分割中的过分割与欠分割问题,并对其直方图进行多次迭代算法,最大程度找出图像中具有的细节信息。文中设计的人机交替分割系统,很好地解决了数据运算大、计算复杂等问题。具有很好的操作性和实用性,为下一步医生的诊断提供了很好的依据。同时,由于参数的设置依赖于图像的先验性,因此在进行三维自适应迭代分割过程中,必须反复交替对参数的设计进行调节,也得到合适的分割阈值。

参考文献

[1]向日华,王润生.一种基于高斯混合模型的距离图像分割算法[J].软件学报,2003,14(7):1250-1257.

[2]刘伟强,陈鸿,夏德深.基于马尔可夫随机场的快速图像分割[J].中国图像图形学报,2001,6(A版)(3):228-233.

[3]汪俊,周来水,安鲁陵,等.基于网格模型的一种新的区域分割算法[J].中国机械工程,2005,16(9):796-801.

[4]Kun Chang Yu,Erilc L Ritman,William E Higgins.System for the Analysis and Visyalization of Large3D An-atomical Tress[J].IEEE Transactions on Medical Ima-ging,1990,19(8):384-395.

[5]Babaguchi N,Yamada K,Kise K,et al.Connectionist Mod-el Binarization[C].Proceding10th ICPR,1990:51-56.

自适应迭代算法 篇2

在噪声消除系统中,若观测噪声为有色噪声,则基于最小二乘准则(LS),提出了两段RLS-RELS算法,这是一种改进的.递推增广最小二乘法,该自适应算法能显著减小噪声的影响,提高信号质量.并在此基础上提出了计算噪声方差的估值方法.计算机数值仿真例子和信噪比的计算比较证明了算法的正确性和有效性.

作 者：许鹏 窦寅丰 XU Peng DOU Yin-feng 作者单位：许鹏,XU Peng(暨南大学珠海学院计算机系,珠海,519070)

窦寅丰,DOU Yin-feng(黑龙江大学电子工程,哈尔滨,150080)

自适应迭代算法 篇3

关键词：农业图像； 滤波增强；非下采样轮廓波变换；离散小波变换；阈值函数模型

中图分类号：TP391.41；S126 文献标志码： A文章编号：1002-1302（2015）09-0448-02

农业图像由于受到复杂的成像环境、成像设备自身缺陷等因素的影响，导致图像中含有一定程度的噪声并且对比度不高，影响了图像的使用。对于该类图像的处理，主要思路有：（1）去噪。采用去噪算法在空间域或者变换域滤除图像中的噪声，如中值滤波[1]、非局部均值滤波[2]、非下采样轮廓波变换（NSCT）[3]、离散小波变换（DWT）[4]等，该类算法在滤除噪声的同时，较少地顾及改善图像的对比度，提高滤波后图像的视觉效果。（2）增强。采用如直方图均衡化[5]、模糊增强[6]等方法改善图像的对比度，突出图像中感兴趣的目标信息，但该类方法很难滤除图像中的噪声，且容易在增强图像的同时放大噪声。农业图像作为一种细节信息较为丰富的图像，对其进行处理，应当有效兼顾滤波和增强这2个环节。因此，结合非下采样轮廓波变换（NSCT）和离散小波变换（DWT），通过对农业噪声图像进行NSCT分解与重构，实现背景图像和细节图像的有效分离，在此基础上，对2者分别进行模糊增强和小波去噪，将处理后的2类图像进行叠加，获得去噪后的农业图像。

1算法基本原理

相对于DWT来说，NSCT在对农业噪声图像分解过程中，能更为有效地捕捉到图像中的线、面等信息，能够实现对图像更为稀疏的表达，但具有相当大的数据冗余度。DWT对于刻画图像中点奇信息是十分有效的，能从水平、垂直、对角等3个方向来刻画图像中的高频信息，尽管其对图像的稀疏性表达能力略逊于NSCT，但具有较小的数据冗余度。因此，构建了NSCT-DWT图像分析框架，通过对农业图像进行NSCT变换和逆变换，获得背景图像和细节图像；然后对细节图像采用DWT进行进一步分析。

1.1背景图像自适应模糊增强

为了更为有效地突显农业图像中的细节信息（如植物的根茎叶边缘等信息），对Pal-King算法[6]中的模糊隶属度映射函数进行改进，改进后的函数模型定义为：

式中：f′（x，y）为公式（1）中像素点灰度值f（x，y）模糊增强后的结果；T-1（）为模糊逆变换计算函数。

1.2细节图像自适应小波去噪

关于小波阈值函数模型，学者们提出一些改进型的模型[7-10]，该类函数去噪模型尽管取得了一些效果，但不足之处在于：（1）自适应不强。函数模型中加入了调节参数，该类参数的取值只能经过反复试验获得，并且对于一类图像甚至是一幅图像来说，如此取值可以满足去噪要求，但对于其余类型图像，则需要重新确定参数值。（2）耗时较大。相当一部分去噪模型中融入了诸如指数、对数之类的非线性函数，在逐步提高函数整体去噪效果的同时，去噪过程的耗时也急剧增大。因此，本研究提出了如下函数模型：

0|w|≤t1。（12）

1.3算法实现步骤

步骤1：对农业噪声图像进行NSCT多尺度分解，得到1幅低频图像和1幅带通图像；低频图像主要包含了农业图像中背景信息，带通图像则包含了农业图像中图像轮廓边缘等细节信息以及绝大多数的噪声信息。

步骤2：对步骤1中获得的带通图像采用一种非下采样多方向滤波器进行多方向分解，获得代表各个方向上的农业像高频信息的方向带通图像序列。

步骤3：对步骤1中的低频图像继续执行步骤1和步骤2，完成对图像更为细致的分解。

步骤4：对上述低频图像和带通图像分别进行NSCT系数重构，获得空间域中的背景图像和细节图像。

步骤5：分别对背景图像和细节图像进行自适应增强和滤波处理，得到处理后的背景图像和细节图像。

步骤6：对经过步骤5处理后的背景图像和细节图像进行叠加，得到最终处理后的农业图像。

2试验测试

借助MATLAB语言对本算法进行编程，采用1幅水果图像对本算法去噪能力进行测验，引入文献[8]和文献[10]中的小波阈值去噪算法与本算法进行对比（分别记为算法1、算法2、算法3），试验结果如图1所示。采用峰值信噪比（peak signal to noise ratio，PSNR）[11]、平均结构相似度算子（means structural similarity，MSSIM）[12]以及算法耗时对上述试验结果进行客观评价，结果见表1。

表1试验结果的定量评价

噪声方差PSNR（dB）算法1算法2算法3MSSIM算法1算法2算法3算法耗时（s）算法1算法2算法3525.35125.60826.1880.7010.7760.8343.4746.1924.8231022.30924.88525.5520.6550.7350.7683.6196.7035.0021520.02823.61225.1730.6020.6840.7363.9276.8915.202

就图像清晰度而言，图1-d稍优于图1-c这，因为算法2中通过将所有的图像小波分解系数分成3类，对于每一类的小波分解系数设置不同的计算方式进行处理，相对于算法1而言，能够更为精细地去除残留在小波系数中的噪声信息，但该算法的计算复杂度明显高于算法1。而图1-d的清晰程度明显优于图1-b和图1-c，主要体现在：（1）图1-d背景中的噪声点基本不存在，而图1-b和图1-c的背景中仍然存在不同程度的噪声，该类噪声的存在影响对图中枇杷果实信息的准确识别；（2）图1-d中枇杷果实轮廓基本从噪声信息中恢复出来，并且图像的对比度比图1-a有了明显提高。

当图像中的高斯白噪声方差为5时，3种算法的PSNR和MSSIM指标值相差不大。当噪声方差达到15时，本算法的PSNR值高于其他2种算法，分别比算法2高1.561 dB比算法1高5.145 dB；本算法的MSSIM值高于其他2种算法，分别比算法2高0.052、比算法1高0.134。对于含有3种不同方差高斯白噪声的农业图像去噪耗时来说，本研究算法相对于文献[12]中的小波阈值去噪算法而言，仍有一定的优势。

3结束语

结合NSCT和DWT，提出了1种农业图像滤波增强算法。通过对含有不同方差高斯白噪声的农业图像进行试验，主观分析结果，PSNR、MSSIM以及算法执行时间等指标定性、定量证明了本研究算法的有效性。

参考文献：

[1]周仕友. 一种结合小波变换的豆类图像增强算法[J]. 江苏农业科学，2014，42（7）：433-435.

[2]陈明举. 基于统计特性的非局部均值去噪算法[J]. 液晶与显示，2014，29（3）：450-454.

[3]韩伟，刘强. 一种 NSCT 域改进阈值函数的杂草图像去噪方法[J]. 江苏农业科学，2013，41（11）：151-153.

[4]Luo F，Wu S J，Jiao L C，et al. Asic design of adaptive threshold denoise DWT chip[J]. Journal of Electronics，2002，19（1）：1-7.

[5]赵吉文，魏正翠，汪洋，等. 基于灰度带比例的优质西瓜子识别算法研究与实现[J]. 农业工程学报，2011，27（4）：340-344.

[6]Pal S K，King R A. One edge detection of X-ray images using fuzzy sets[J]. IEEE Trans，1983，5（1）：69-77.

[7]李雷，魏连鑫. 自适应小波阈值函数的指纹图像去噪[J]. 上海理工大学学报，2014，36（2）：154-157，162.

[8]黄玉昌，侯德文. 基于改进小波阈值函数的指纹图像去噪[J]. 计算机工程与应用，2012，50（6）：179-181.

[9]王蓓，张根耀，李智，等. 基于新阈值函数的小波阈值去噪算法[J]. 计算机应用，2014，34（5）：1499-1502.

[10]王艺龙，杨守志. 基于连续阈值函数的小波去噪方法[J]. 汕头大学学报：自然科学版，2013，28（4）：66-74.

[11]杨福增，田艳娜，杨亮亮，等. 基于杂交小波变换的农产品图像去噪算法[J]. 农业工程学报，2011，27（3）：172-178.

自适应迭代算法 篇4

关键词：卫星信道,因子图,迭代均衡,自适应均衡,最小均方误差

1 引言

随着卫星通信的发展,高阶调制和高码率越来越多地被应用在系统设计中,而卫星转发器中滤波器的线性失真特性(包括群时延特性和幅频特性)会对高阶调制的信号性能产生严重恶化。因此必须在接收端采用均衡技术,以降低失真影响,提高卫星通信系统的吞吐率。

均衡器当中的线性均衡器因其简单而被广泛应用,均衡器阶数(抽头个数,L)可以与信道多径数(J)相当,运算复杂度与系统所采用的调制方式(M)无关,仅为O(L),但是线性均衡器并不能完全消除码间干扰,只能降低码间干扰的影响。迭代均衡器主要包括将译码结构与均衡技术相结合的Turbo均衡器、基于因子图的置信度传播(BP,Belief Propagation)均衡器和基于MMSE原则的线性均衡器。Turbo均衡算法是基于状态转移模型的一种最大后验概率(MAP,Maximum A Posteriori)检测法;而基于因子图模型的BP迭代均衡算法,因线性失真信道形成的因子图中存在大量长为4的环,会恶化均衡器的性能,是一种次优检测算法。假设高阶调制中,星座点个数为Ω=2M,符号个数为Ns,则Turbo均衡器与BP均衡器的运算复杂度为O(Ns×ΩL),与均衡器的阶数呈指数关系。虽然BP均衡器与Turbo均衡器的复杂度相同,但是BP均衡器采用的和积算法(SPA,Sum-Product Algorithm)和最小和算法(MSA,Min-Sum Algorithm)的运算量要小很多,在实际应用中具有很大的优势。本文将线性均衡器与BP迭代均衡器结合,得到混合式迭代均衡算法(HIE,Hybrid Iterative Equalization)。利用复杂度较低的线性均衡器消除大量幅度较小的码间干扰,用性能较好的迭代均衡器消除少量幅度较大的码间干扰。

2系统模型

本节介绍卫星通信线性失真系统模型和其对应的因子图表示。在卫星通信中,二进制信息比特序列B经过LDPC编码器后变为分组码LDPC码C,本文采用Wi MAX标准中的LDPC码进行仿真说明。码字C经过星座映射变为符号序列S,每个符号用M比特表示,共有Ω=2M种取值。

符号序列S经过卫星信道中的滤波器和行波管放大器后,进入加性高斯白噪声信道(AWGN),其中滤波器特性包含发送端的成型滤波特性、转发器的滤波特性等。接收机收到的信号yn为:

下面根据典型卫星信道的线性失真,研究其冲击响应特性。假设群时延的最大时延为Td,则线性群时延(A)和抛物线群时延(B)特性对应的滤波器频域函数的相位值为下式。

通过实际测试可知,线性失真导致的多径中,除了主径外,其它多径(副径)的幅度都较小(<0.2),但是副径的数目较多(>10)。图1描述了系统的因子图表示,其中变量节点包括接收到的观察节点(Y)、符号节点(S)、以及组成符号的比特节点(B);函数节点包括均衡函数节点(Eq,Equalization)、星座映射节点(M,Mapping)和校验函数节点(C,Check)。

3迭代接收算法

本节将在介绍函数节点的基础上,根据和积算法引入迭代均衡、迭代调制和迭代均衡。

在图1的因子图中,均衡节点(Eq)对应的方程为:

上式中p(Y|S)表示接收序列对于发送的可能序列的似然函数。根据和积算法:

(1)符号节点s传播到函数节点Eq的外信息:

其中~{Sk}表示序列S中除第k个符号确定已知外,其它符号遍历其所有Ω=8状态的序列集合,集合的序列数量为ΩNs-1。

(3)映射函数节点M传递给符号节点s的外信息为:

(4)符号节点s传递给函数节点M的外信息为:

上式中N(sk)表示与符号节点sk相连接的均衡节点集合。

(5)变量节点b传递给映射函数节点M的外信息为:

(6)映射函数节点M传递给比特节点b的外信息为:

最后LDPC译码器根据校验关系更新比特节点的后验概率信息,至此完成整个接收过程。

4 仿真结果

假设卫星信道B为抛物线群时延,其中最大延时Td=Rs。图2显示了在该信道条件下不同均衡方式的误码性能,其中映射方式都为8PSK,线性均衡器的阶数为11,迭代均衡器的阶数为4。图2中LDPC码的码长为960,码率为5/6。可以看出,在误码率1E-5时,本文提出的HIE算法性能比LMS算法优0.5d B,与AWGN信道下的信噪比仅差0.2d B,有效地提高了均衡器的性能。而迭代均衡器因其阶数有限,均衡性能较差。

5 结论

本文首先分析了卫星信道的模型和线性失真冲击响应特性,即主径幅度大,其它多径数量多且幅度小。然后引入了基于因子图的迭代接收算法,包括迭代均衡算法和迭代调制算法。本文的创新点在于基于卫星信道线性失真特性,提出横向线性均衡器与迭代均衡器相结合的混合式均衡算法。该算法有效地利用了线性均衡器的抽头个数不受限但简单的特点和迭代均衡器抽头数受限但准确的特点,得到一种结构简单且性能良好的均衡器。仿真结果表明在典型群时延信道和8PSK映射方式下,HIE要优于最小均方误差(LMS)线性均衡器0.5d B。为了进一步提高混合均衡器的性能,将在之后工作中联合优化线性均衡器与迭代均衡器的系数。

参考文献

[1]Proakis J G,Salehi M.Digital Communications(Fifth Edition).Americas,New York:Mc Graw-Hill Companies,Inc,2008.691-696.

[2]Kschischang F R,Frey B J,Loeliger HA.Factor Graphs and the Sum-Product Algorithm.IEEE Transactions on Information Theory,2001,47(02):498-519.

[3]Wymeersch H.Iterative Receiver Design.New York:Cambridge university press,2007.

自适应迭代算法 篇5

针对GPS/DR组合导航Kalman滤波的异常扰动影响问题,引入了自适应滤波算法.给出了由预测残差确定自适应因子的过程.利用实测数据进行验证,结果表明无论是单因子自适应滤波还是多因子自适应滤波都能够很好地控制状态异常对滤波估值的影响,滤波精度均优于标准Kalman滤波导航解;而且因为多因子自适应滤波避免损失可靠的状态参数信息,较单因子自适应滤波,精度又有明显提高.

作 者：吴富梅 杨元喜 田育民 WU Fu-mei YANG Yuan-xi TIAN Yu-min 作者单位：吴富梅,WU Fu-mei(信息工程大学,测绘学院,河南,郑州,450052;西安测绘研究所,陕西,西安,710054)

杨元喜,田育民,YANG Yuan-xi,TIAN Yu-min(西安测绘研究所,陕西,西安,710054)

自适应迭代算法 篇6

目前, CT已经成为一种主要的医学放射暴露源。最新的流行病学证据表明, 高剂量头颅CT的辐射会增加患白血病和颅脑肿瘤的风险[1,2,3], 所以寻找一种降低CT剂量的新方法是目前迫切需要解决的问题。CT最常用的重建算法是滤过反投影 (filetered back projection, FBP) 重建算法, 该算法的特点是:当降低放射剂量时, 图像信噪比会随之下降。而自适应统计迭代重建 (adaptive statistical iterative reconstruction, ASIR) 算法是最新的一种技术, 与FBP算法相比, ASIR最大的优势是能在常规FBP算法表现为低信噪比时明显提高图像质量, 也就是说ASIR能在低剂量扫描时提高图像信噪比, 因此ASIR为临床低剂量扫描开辟了一种新途径。

对于ASIR算法, 前人在胸部[4,5,6]、腹部[7,8,9]、心脏[10,11]及大血管[12]方面低剂量扫描进行了大量的研究并报道, 而在颅脑方面的应用报道较少。本研究旨在发现ASIR在颅脑扫描中降低辐射剂量的证据。

1 资料与方法

1.1 研究对象

受试者来自于2011年3月至2013年5月解放军总医院门诊就诊患者, 实验经解放军总医院医学伦理委员会批准。受试者均签署知情同意书。

纳入标准:CT检查示脑内无明确病变, 能顺利并配合完成检查, 签署知情同意书。

排除标准:CT检查结果证实有以下疾病: (1) 脑血管病, 如脑梗塞或脑出血病史等; (2) 脑肿瘤, 如脑原发肿瘤、转移瘤; (3) 脑外肿瘤, 如脑膜瘤等; (4) 脑外伤; (5) 异物; (6) 其他颅内疾病。

1.2 研究方法

1.2.1 数据采集

采用型号为Discovery CT 750 HD宝石能谱CT, 对32例受检者进行头部扫描, 年龄 (74±16.48) 岁。扫描条件120 k V、220 m A, 球管旋转1 s。轴扫, 间隔:5 mm, 层厚:5 mm, 共24张图像。

1.2.2 数据测量与处理

扫描结束后在双侧侧脑室最大层面重建图像, FOV 25 cm、层厚5 mm、间隔20 mm, 共5层, 分别将ASIR值设定为0%、20%、40%、60%、80%、100%重建图像, 窗宽为100 HU、窗位35 HU。在一侧侧脑室内选取圆形感兴趣区 (region of interest, ROI) (面积为9.9 mm2, 周长为7.8 mm, 像素点为30) , 测量不同的ASIR值时侧脑室脑脊液的CT值与标准差 (如图1所示) 。由2名高年资放射诊断科医生对数据进行测量, 每一个数据测量2次后取平均值, 采用双盲的方法记录数据。

1.2.3 统计学分析

应用统计学软件SPSS 13对数据进行统计学分析, 首先检验病灶的CT值和标准差等数据是否符合正态分布, 之后对数据进行t检验。

2 结果

当ASIR值为0%、20%、40%、60%、80%、100%时, 侧脑室脑脊液的CT值SD分别为 (4.706±1.192) 、 (3.466±0.872) 、 (2.8±0.761) 、 (2.284±0.618) 、 (1.784±0.547) 、 (1.425±0.507) , 各组间具有明显的统计学差异, P<0.01, 随ASIR百分比的增加, 图像噪声水平明显降低。

3 讨论

医用辐射已成为人类接受X线照射最大的来源。特别是近年来CT技术的发展, 病灶的检出率大大增加, CT应用的范围也越来越广泛, 导致CT检查次数在不断增加, 从而引起辐射剂量的增加。CT产生的电离辐射具有潜在致癌的风险, 会增加患白血病和颅脑肿瘤等恶性肿瘤发生的几率[1,2], 同时因CT扫描速度比较快, 所以头部CT扫描临床应用也越来越多, 进而带来的是放射剂量的危害也越来越严重。因此, 当前降低CT辐射剂量、提高有效剂量是目前迫切需要解决的问题。

随着社会发展进步和人们生活水平的提高, 人们越来越意识到放射防护的重要性。人们之前用了很多方法尝试降低头部扫描剂量[13,14,15,16], 如降低曝光量、遮挡腺体和眼睛等重要部位、自动球管电流调节、倾斜扫描架避免扫描眼眶等, 这些方法能在一定程度上降低CT扫描剂量。此外, 学者们还尝试在重建技术上降低辐射剂量, 其中, 临床上采用FBP重建算法时是通过降低曝光量而达到降低剂量的目的。但这种方法的局限性在于:该方法会导致图像噪声的增加和图像对比噪声比的下降。这是因为:FBP重建方法虽然速度较快, 但是该方法要求每次投影测量数据都是精确定量和完全的, 而FBP重建方法对噪声和伪影都很敏感, 当辐射剂量降低或投影数据采集不足时, 重建出的图像质量就会很差。因此, 为了保证图像质量, 采用FBP重建方法就不能大幅度地降低CT扫描的辐射剂量。

ASIR重建方法是一种较新的重建方法。该方法使用从FBP算法获取的图像信息作为图像重建的原始结构单元, 通过建立系统噪声模型并且利用迭代方法加以抑制, 从而得到更加清晰的图像。该方法是通过将FBP数据和ASIR数据按不同比例加权融合而达到不同程度的降噪效果。与FBP重建算法相比, AISR算法重建图像时图像噪声更低[17]。本研究就是依据这个原理进行研究的, 研究结果证实了AISR重建算法能有效地降低图像噪声。结果显示, 随着ASIR数据按不同比例加权的增加, 图像对比越来越平滑, 噪声越来越小。也就是说ASIR百分比的增加, 图像噪声水平明显降低。提高ASIR权重值能够降低颅脑CT图像噪声, 因此扫描时采用较高ASIR值重建图像具有降低射线剂量的作用。本研究结果为临床采用ASIR方法降低颅脑扫描剂量提供了理论依据。

本研究中, 随着ASIR算法比例加权的增加, FBP算法比例加权降低, 图像噪声降低。从这个角度上也说明, ASIR算法比FBP算法在降低图像噪声上更优越, 是一种更好的降噪手段, 对于降低辐射剂量更有效。在临床使用过程中, ASIR算法重建需要的时间很短 (据文献报道不会超过5 min[18]) , 不会影响临床的效率, 说明该算法有很好的临床可行性。

本研究采用图像中脑脊液作为测量对象, 这是因为有研究对ASIR算法和FBP算法2组脑组织噪声差异进行对比发现, ASIR算法组脑脊液图像噪声明显下降, 而脑白质、颅骨和空气的噪声对比2组没有差异[19], 原因是脑脊液密度接近于水的密度, 且CT值相对比较恒定, 不容易受其他因素的干扰。这里唯一需要考虑的是脑脊液流动伪影造成的影响, 而双侧侧脑室内脑脊液相对流速较慢, 就有效地避免了这个问题。颅骨和空气CT值变化较大, 同时与其他组织相邻会有界面征, 会给CT值测量带来麻烦。脑白质虽然没有颅骨和空气对于CT值影响那么大, 但脑白质本身在不同部位和纤维束走行都会影响测量结果。综上所述, 选择侧脑室脑脊液来测量CT值结果更科学。

自适应迭代算法 篇7

在汽车工业领域, 许多电线端子上需要安装防水密封的圆柱形橡胶栓, 其对橡胶栓的外形质量要求较高, 圆柱面上不能有明显的缺陷。传统的橡胶栓的检测主要依靠人工分拣, 费时费力, 而且受人为因素的影响, 检测效果也不理想。因此, 利用数字图像处理技术对橡胶栓进行在线检测以保证橡胶栓的质量, 具有重要的意义。

橡胶栓是通过注塑机成形的, 常见的缺陷主要是圆周上的飞边和缺料的凹坑等, 而CCD一次拍摄只能取到一半的圆柱面, 就可能会漏过橡胶栓背面的缺陷。因此在在线检测中, 课题组巧妙地采用了一组螺杆传送机构[1], 该机构不仅推动橡胶栓往前运动通过CCD, 同时使橡胶栓随螺杆的旋转时本身也回转, 这样在通过CCD下方时隔一定时间连续拍摄两次或三次, 基本上可实现整个柱面的检测。

由于橡胶栓在螺杆上是自由的, 因此在CCD拍摄时其位置角度是随机的, 而在图像识别中, 采用的方法往往是将被检测图像与标准模板图像进行比对来识别缺陷, 因此在利用图像处理技术对橡胶栓的表面缺陷进行识别时, 首先必须对不同位置的橡胶栓图像进行快速自适应旋转, 使图像中橡胶栓的中心轴与标准模板中橡胶栓的中心轴一致。

图像的自适应旋转是数字图像处理中的一个重点和难点, 目前, 这方面的文献报道较少。图像自适应旋转采用的方法主要有惯性主轴法[2]和圆形编码匹配[3]的方法。前者先检测出图像的边缘, 再利用边缘信息计算图像的质心和惯性主轴角度, 然后绕质心旋转, 所以该方法受限于边缘检测的精度, 并且容易受光照等外界条件的影响。后者将图像在圆周方向等分, 并建立方向编码来解决图像旋转匹配的问题, 计算量很大, 速度较慢, 难以满足在线检测的图像处理实时性要求。

本研究提出一种基于迭代匹配的快速图像自适应旋转算法, 将最小距离匹配的方法引入图像旋转, 然后通过迭代计算, 精确快速地实现橡胶栓图像的自适应旋转。

1 最小距离匹配

1.1最小距离匹配的原理

最小距离匹配又称模板匹配[4], 是一种重要的图像识别方法。最小距离匹配主要用于判断两个物体的相似度, 越相似则距离越小, 反之越大。

在利用最小距离匹配前, 需要预先将n个标准对象的特征存储为n个标准模板, 然后将具体的待检测样品与已知的标准模板进行比较, 计算它们之间的距离, 若待测样品与第k (k∈[0, n]) 个模板之间的距离最短, 则判定待测样品与第k个标准对象是同一类的物体。

两个对象间的距离可以用如下的方法计算:

例如, 假设一个标准对象A可以抽象为一个m维的特征向量:

XA=[xA1xA2 … xAm]T (1)

同理, 任意的待检测对象B也可以抽象为一个m维的特征向量:

XB=[xB1xB2 … xBm]T (2)

则, 待检测对象与标准对象间的距离可表示为:

$d{}_{AB}=d(X{}_A,X{}_B)=\left[{\displaystyle \sum _{i=1}^m(}x{}_{Ai}-x{}_{Bi}){}^2\right]{}^{\frac{1}{2}}$ (3)

由于上式涉及到乘方、开方等复杂运算, 可将距离公式简化为:

$D{}_{AB}=D(X{}_A,X{}_B)={\displaystyle \sum _{i=1}^m\left|x{}_{Ai}-x{}_{Bi}\right|}$ (4)

1.2特征向量的选取

利用最小距离匹配需要将标准对象抽象为一个特征向量, 然后存储为模板。对实际的待检测样品, 也需要用同样的方法将其抽象为一个特征向量, 然后才可以计算它与模板间的距离。因此, 特征向量的设计具有重大的意义。

目前, 对特征向量的选取没有固定和统一的方法;相反, 随着研究的深入, 选取方法越来越多, 如利用图像的几何矩构成特征向量[5]、利用图像的灰度信息构成特征向量[6]等。本研究根据橡胶栓的实际特点和检测要求, 以橡胶栓在图像中所占比重为特征值来构成图像的特征向量。如图1所示, 将整幅图像等分为M×N个小块, 每块的像素数为X×Y。对图像中的任一点像素, 令:

则图像中第m行第n列小块的面积比重为:

$\frac{S{}_{m,n}}{S}=\frac{{\displaystyle \sum _{i=m\cdot X}^{(m+1)\cdot X-1}{\displaystyle \sum _{j=n\cdot Y}^{(n+1)\cdot Y-1}}}g(i,j)}{X\cdot Y}$ (6)

对任意的m、n, 令k=M×n+m, 于是, 面积比重亦可表示为:

$x{}_k=\frac{S{}_k}{S}=\frac{S{}_{{\rm M}\cdot n+m}}{S}=\frac{S{}_{m,n}}{S}$ (7)

其中, m∈[0, M-1], m∈[0, N-1], k∈[0, MN-1]。所以, 该对象的特征向量就可以表示为:

X=[x1x2 … xMN-1]T (8)

上述方法得到了对象的k维特征向量, 在检测过程中, 需要根据实际的需要, 选择适当维数的特征向量。特征向量的维数越多, 那么匹配的精度越高, 但计算的时间复杂度更大;反之, 维数越少, 则精度降低, 但计算速度更快。

1.3最小距离匹配在自适应旋转中的应用

根据最小距离匹配的原理, 在图像旋转中, 假设被存储为模板的标准对象的旋转角度为0°, 则当图像的旋转角度越接近0°, 则它与模板间的距离越小;反之, 而当旋转角度增大时, 它与模板间的距离也随之增大[7,8]。

在橡胶栓检测过程中, 假设图像中橡胶栓中心轴水平时橡胶栓的旋转角度为0°, 并将该图像存储为模板, 那么当被检测的橡胶栓图像的旋转角度从-90°～+90°时, 它与模板间的距离变化关系如图2所示。

从图2可见, 在图像的旋转角度从-90°～+90°的范围内, 用最小距离匹配的方法来实现橡胶栓图像的自适应旋转是可行的。

对旋转角度不在-90°～+90°范围内的橡胶栓图像, 可根据橡胶栓的具体特点, 在实际检测过程中, 先利用数字图像处理技术将橡胶栓图像水平镜像, 然后再利用最小距离匹配的方法来实现图像的自适应旋转。

2 迭代匹配过程

2.1原始迭代方法

根据最小距离匹配原则, 在-90°～+90°的范围内, 待检测对象旋转角度偏离0°越大, 它与模板间的距离越大, 所以, 在实际的检测过程中, 可根据检测精度要求, 指定一个正数C, 当图像与模板间的距离小于指定正数C时, 即可认为图像已经旋转到位。对当前的待检测图像, 可以先计算距离D0, 然后分别计算顺时针和逆时针转动θ角 (θ>0, 下文顺时针用正数表示, 逆时针用负数表示) 后的偏离距离D1和D2, 则D0、D1、D2的4种位置关系如图3所示。

根据D0、D1、D2的位置关系就可以进一步确定之后图像旋转的方向和角度, 直到图像与模板的距离小于指定正数C。具体的迭代步骤如下:

(1) 判断D0是否小于指定正数C, 若是, 则跳到第 (8) 步;

(2) 计算当前位置的距离D0和当前位置正反转动θ角 (θ>0) 后的位置1、位置2的距离D1、D2;

(3) 若D2>D0>D1, 则继续;若D2<D0<D1, 则θ=-θ, 继续;若D0最小, 则跳到第 (6) 步;

(4) 判断D0是否小于指定正数C, 若是, 则跳到第 (8) 步;否则, 转动θ角, 计算当前位置距离D0′;

(5) 若D0′<D0, 则令D0=D0′, 跳到第 (4) 步;否则令θ=-θ/2, 跳到第 (4) 步;

(6) 若D2>D1, θ=θ/2;否则, θ=-θ/2;

(7) 转动θ角, 计算当前位置距离, 并赋值给D0, 若D0小于指定正数C, 跳到第 (8) 步;否则令θ=-θ/2, 跳到第 (7) 步;

(8) 结束。

2.2简化的迭代方法

原始的迭代算法通过不断地调整旋转方向和旋转角度实现自适应旋转, 可确保旋转的精度, 但如果迭代次数非常多, 势必影响迭代运算的速度。所以, 本研究下面提出一种简化的迭代算法, 可有效提高运算速度。

(1) 计算当前的距离D0;

(2) 将图像旋转θ角 (θ>0) ;

(3) 计算新位置的距离D1;

(4) 若D1<D0, 沿原方向不断转动θ角;若D1>D0, 则沿反方向不断转动θ角;

(5) 不断旋转θ角, 直到新位置的距离D1又变大;

(6) 反方向转动θ/2角。

(7) 结束。

简化的迭代算法的实现过程如图4所示。仅经过7步的迭代运算即可将图像旋转到位。

简化的迭代算法比原始的迭代算法的迭代步骤少, 因此, 它具有更快的运算速度。但由图4可见, 简化的迭代算法最后得到的图像与模板的距离并不一定是最小的, 它与最小距离可能有一点偏差, 因此为保证旋转的精度, 在计算过程中需要合理地选取θ值, 保证最后的图像与模板的距离小于指定正数C。

3 实验与结果

为验证基于最小距离匹配的橡胶栓图像自适应旋转算法的有效性, 在实际检测条件下, 本研究分别采用传统的惯性主轴法和原始迭代法、简化迭代法对具有不同旋转角度的橡胶栓图像进行旋转匹配。

实验中相关参数如下:

(1) CCD摄像头的分辨率:640×480;

(2) CPU:Intel Pentium (R) 4, 主频2.80 GHz;

(3) 内存:512 MB;

(4) 光源:LED环形灯;

(5) 模板:旋转0°的标准图像, 100维特征向量。

在实验中, 笔者分别对旋转角度为0°、10°、30°、45°、60°和80°的待匹配的图像进行对比实验, 如图5所示 ( (a) ～ (e) 为不同旋转角度的待检测图像, (f) 为自适应旋转后旋转角度接近0°的图像) 。

若实验中取正数θ=0.5, 初始角度θ=5°, 则两种迭代法匹配所需要的时间、迭代次数及最后距离如表1所示。

由表1可知, 相对于传统的惯性主轴法, 本研究中的原始迭代法能有效地完成橡胶栓图像的自适应旋转, 且效率没有明显降低;而简化迭代法不仅能有效地完成了橡胶栓图像的自适应旋转, 而且计算速度也有显著的提高。但是, 相对于原始迭代法, 简化迭代法的自适应旋转匹配较快, 但精度稍低 (最后距离大) 。因此, 在实际工业检测中, 可根据检测的精度要求选择适当的迭代匹配方法。

4 结束语

本研究将模式识别中的最小距离匹配引入图像的自适应旋转, 并利用迭代法实现匹配过程, 在保证匹配精度的同时, 也具有较快的运算速度, 因而在橡胶栓的检测中得到了很好的旋转匹配效果, 为后续橡胶栓图像表面缺陷识别提供了基础。另外, 笔者提出了用原始的迭代方法和简化的迭代方法分别实现自适应旋转匹配, 两者的侧重点有所不同, 可以根据实际的要求选择合适的算法。

参考文献

[1]魏燕定, 周豪, 陈波扬.在线橡胶栓外形质量检查机:中国, 101349550[P].2009-01-21.

[2]李金泉, 徐正华, 胡晓飞, 等.染色体质心、惯性主轴及旋转算法的研究[J].信号处理, 2001, 17 (3) :274-277.

[3]FARHAN U, SHUN-ICHI K.Using orientation codes for ro-tation-invariant template matching[J].Pattern Recogni-tion, 2004, 37 (2) :201-209.

[4]MCQUEEN M P C.A generalization of template matchingfor recognition of real objects[J].Pattern Recognition, 1981, 13 (2) :139-145.

[5]刘明, 叶正麟, 陈作平.基于二维特征向量的快速分形编码方法[J].计算机工程与应用, 2007, 43 (8) :82-84.

[6]陈玲, 沈红标, 李咸伟, 等.改进的图像纹理检索方法在矿石识别中的应用[J].中国图象图形学报, 2006, 11 (11) :1700-1703.

[7]RAFAEL C G, RICHARD E W.Digital Image Processing[M].2nd ed.London:Prentice Hall, 2007.

自适应迭代算法 篇8

迭代学习控制不需要依赖动态系统的精确数学模型就可以使系统实际输出完全跟踪期望输出。少量的先验知识就可以使迭代学习控制得以实现,并且重复运行次数越多控制精度越高,因此在很多具有重复运行特性的被控对象中得到了应用[1,2,3]。PID迭代学习控制以其良好的可靠性、鲁棒性已成为迭代学习控制领域中最为常用的控制策略。为了使PID学习律达到更好的控制效果,常将其与其他一些智能控制方法相结合,以期对非线性、不确定的复杂控制系统达到更好的控制效果。模糊控制对解决不确定性高的复杂非线性系统表现出诸多优越之处,因此很多学者将模糊控制和迭代学习控制相结合,以求达到较为理想的控制效果。张航等[4]采用自整定模糊控制器作为迭代学习律,并对单臂机械手进行了仿真研究,取得了较好的控制效果;张丽萍等[5]运用并行分配补偿方法确定T-S模型的迭代学习控制器结构,并给出了误差收敛条件;蒋思中等[6]利用作为模糊控制器的输入,生成PD增益矩阵的调整因子,以达到提高收敛速度的目的。

本文提出一种模糊自适应PID迭代学习控制算法,利用模糊控制器对PID学习律的参数进行实时整定,可以提高系统的控制精度,使控制系统具有较好的稳定性和较强的鲁棒性,同时可以提高系统的收敛速度。

1 问题的提出

考虑如下形式的线性系统:

其中t∈[0,T],x(t)∈Rm为系统的状态,u(t)∈Rn为系统的输入,y(t)∈Rr为系统的输出,A、B、C为适当维数的矩阵。

迭代学习控制是对给定控制系统在时间区间t∈[0,T]内,利用被控对象的期望输出yd(t),通过迭代学习算法寻找优化的控制信号uk(t),使得系统响应yk(t)相对于y0(t)有所改善,寻找uk(t)的过程就是系统学习的过程,且使得当k趋于∞时,满足yk(t)与yd(t)的偏差为零,即系统实现实际输出完全跟踪期望输出[7]。

在第次运行时,系统的动态方程可表示为:

闭环PID型学习律如下:

式中Γ、Φ、Ψ为学习增益矩阵,ek+1(t)为实际输出与期望输出的偏差:

可见PID学习律中的参数是一成不变的,目前所需要的就是寻找到一种方法来对式中Γ、Φ、Ψ三个参数进行实时校正,实现系统的动态学习过程,以提高系统的控制精度和收敛速度。

2 模糊自适应PID迭代学习控制算法

2.1 控制系统基本结构

模糊自适应PID迭代学习控制算法结构如图1所示,控制系统采用闭环学习控制算法,利用传统的经验PID参数作为参考,通过模糊整定单元对经验PID参数进行实时校正,生成精确度更高的模糊PID学习律,提高系统收敛速度,加强抗干扰能力,最终实现系统的完全跟踪,即。

2.2 模糊自适应PID学习律

模糊控制器的输入选用系统输出误差e和输出误差的变化量ec,而模糊控制器的输出则为PID学习律的增益整定量,采用PID学习律时,Kp用于提高系统响应速度,调节系统控制精度,过大则导致系统不稳定,过小则导致系统响应速度缓慢;Ki用于消除系统稳态误差,过大则会导致积分饱和,超调量较大,过小则无法起到相应作用;Kd有助于提高系统的动态性能,模糊控制器生成的整定参数应以此为原则而确定,其控制规则形式如下:

IF ek(t)is E and dek(t)is EC THEN Kpk(t)is Kp。

IF ek(t)is E and dek(t)is EC THEN Kik(t)is Ki。

IF ek(t)is E and dek(t)is EC THEN Kdk(t)is Kd。

定义模糊控制器的输入、输出量的模糊子集均为NB,NM,NS,ZE,PS,PM,PB,即负大,负中,负小,零,正小,正中,正大,其模糊论域均为[-E,E],同时可以建立如表1所示控制规则表:

通过模糊规则对其进行判定,并利用重心法对其进行解模糊处理,以Kp为例,可得模糊控制器输出的的PID参数分别为:

式中μKp为输出量的隶属度函数,Kpq为输出的模糊化变量,Kp为解模糊后的输出量。

由此可以得到模糊PID迭代学习律为:

其中Γ、Φ、Ψ为PID学习律的固定参数。

3 仿真研究

利用仿真对算法的有效性进行验证,选取如下形式LTI系统作为被控对象:

期望跟踪轨迹选用:。各模糊控制器的输入输出论域均选择在[-3,3]范围内,经验PID学习律参数为:Γ=30,Φ=2,Ψ=3。学习周期为T=6s。各变量均采用如图2形式的隶属度函数,其中ZE选用三角形隶属度函数,其余均选择梯形隶属度函数。

控制系统仿真框图如图3所示:

通过仿真可以得到如图4所示的第一次迭代学习下的跟踪轨迹,图5则为系统跟踪误差。

从仿真结果可以看出,系统具有较好的跟踪性能,并随学习次数的增加系统误差不断趋近于零,使得控制系统具有良好的动态性能,达到高精度趋近期望轨迹的控制目的,由此可见模糊自适应PID迭代学习控制系统对于轨迹跟踪具有较好的控制效果。

4 结语

现实问题中常常会遇到需要高精度轨迹跟踪的控制问题,针对这一问题本文提出一种模糊自适应PID迭代学习控制算法,利用模糊推理整定PID学习律的参数,以达到精确控制和快速收敛的目的,并通过仿真验证了方法的有效性,对同类控制对象的控制系统设计具有一定的参考价值。

参考文献

[1]王仲鸿,王强,张东霞,姜齐荣.电力系统暂态稳定问题和迭代学习控制的研究[J].电力系统自动化,1999,23(8):6～10.

[2]徐敏,林辉.基于迭代学习控制理论的励磁控制器设计[J].电力自动化设备,2006,26(3):69～72.

[3]王岩,付永领.模糊滑模迭代学习控制算法在液压系统中应用[J].北京航空航天大学学报,2007,33(1):86～89.

[4]张航,罗大庸,黄浩江,罗熊.机器人模糊迭代学习控制及其仿真研究[J].自动化技术与应用,2002,(2):3～5.

[5]张丽萍,杨富文.基于T-S模型的迭代学习控制算法及其在机器人点位控制中的仿真研究[J].系统仿真学报,2005,17(1):166～169.

[6]蒋思中,朱芳来,王心开,王改云.模糊增益PD迭代学习算法及其应用[J].电光与控制,2009,16(8):72～74.

自适应迭代算法 篇9

自组织竞争网络是一种无导师学习的人工神经网络, 通过竞争机制自动地寻找样本中的内在规律与本质属性, 自组织自适应的改变网络参数。其基本思路是网络竞争层各神经元竞争对输入模式的响应机会, 最后的只有一个神经元成为竞争的胜利者。并对那些与获胜神经元有关的各连接权一同朝着更有利于它竞争的方向调整。获胜的神经元就表示对输入模式的分类。除了竞争方式外, 还可以通过抑制方式、侧抑制方式获胜[1]。

基本的自组织竞争网络由输入层和竞争层组成。可假设竞争层有m个神经元, 输入层有n个神经元, 网络拓扑结构如下图1-1所示。

所对应的网络权值为:

自组织竞争网络具有较好的动态性能, 能够自我学习, 并根据所记忆的学习方式对输入模式进行最接近的分类, 即以竞争层获胜神经元表示分类结果。

2 基于迭代算法的改进神经网络

2.1 算法提出

本文将描述输入模式的变量定义为特征变量。根据1, 自组织竞争网络能够通过自我学习进行模式分类, 在特征变量充足且所有输入模式类别未知的情况下可以得到唯一且较为准确的分类结果。经MATLAB仿真试验, 在部分模式类别已知的情况下训练自组织竞争网络, 可能出现分类结果无法复现的情况 (即两次或多次分类结果不一致) , 故该情况下会得到失准的分类结果。

为弥补这一缺陷, 本文将引入迭代算法对部分模式类别已知情况下的分类结果进行校准。

2.2 算法描述

定义未知类别的模式为未知模式, 已知类别的模式为已知模式。改进算法的基本思路是为每一个未知模式赋一个随机的初始类别, 并与已知模式混合进行第一次分类。若分类结果与初始分类结果不一致, 则将已知模式的类别还原后重新进行一次分类;不断重复分类过程, 直至所有模式的类别不在发生改变。此时得到最终的模式分类结果。

定义Si, k为第i号模式在第k次迭代中的类别。以将n个输入模式分成t类为例, 则Si, k满足下列约束条件:

根据上述描述可知, 得到最终分类结果的充要条件是迭代收敛。迭代收敛的条件方程为:

引入迭代算法, 可使自组织竞争网络自我校准, 满足收敛条件方程, 得到唯一而可靠的分类结果。

2.3 仿真结果

验证迭代自组织神经网络算法的分类功能。

以我国南方某城市为例, 共3个城区, 分别为L区、O区、W区, 为城市中的小学划定所属区域可便于教育管理与辖区入学。已知3所小学的分区情况, 根据经纬度坐标为其余9所学校划定所属区域。

分别采用自组织竞争网络与迭代自组织神经网络进行分区, 得到MATLAB仿真结果如下。

(1) 已知分区:

L区:西山小学O区:塘下小学W区:新川小学

(2) 自组织神经网络分类结果:

L区:西山小学, 广场小学, 东方小学, 实验小学;

O区:塘下小学, 曲溪小学, 源口村小学, 南汇小学, 茶山小学;

W区:新川小学, 光明小学, 桥头小学。

仿真时间:4.5秒

(3) 迭代自组织神经网络分类结果:

L区:西山小学, 广场小学, 实验小学, 南汇小学;

O区:塘下小学, 曲溪小学, 源口村小学, 东方小学;

W区:茶山小学, 新川小学, 光明小学, 桥头小学。

仿真时间:25分30秒

结果显示, 迭代自组织神经网络给出了更为准确的分区结果, 但耗时较长。说明引入迭代算法在提高模式分类准确度的同时, 将很大程度的降低程序运行速度。

3 迭代自组织神经网络算法应用分析

根据2, 迭代自组织神经网络算法适用于已知部分输入模式类别情况下的模式分类。故可利用该算法识别未知模式的类别, 进而应用于电气工程、环境工程等领域。下面将以电气工程为例, 分析迭代自组织神经网络算法在变压器故障检测中应用的可能性。

电力变压器是电网中的主要元件之一, 其安全运行是电网可靠运作的必要条件。变压器油溶解气体的色谱分析是变压器故障诊断的一种重要方法[2]。为提高故障诊断的准确度, 可尝试引入迭代自组织神经网络算法。

以油中5中特征气体 (H2、CH4、C2H2、C2H4、C2H6) 的体积含量作为特征变量。根据长期运行经验, 可得到某型号变压器在正常运行及故障运行时的检测数据, 作为已知模式。同时, 以实时监测数据作为未知模式。采用迭代自组织神经网络法可进行检测数据的模式分类, 得到实时监测数据所属类别, 进而得到变压器故障检测结果。

同理, 迭代自组织神经网络算法也可试用于水质检测、人脸识别等方面。

4 总结

迭代自组织神经网络算法可以在分类过程中进行自我校准, 适用于已知部分输入模式类别情况下的模式分类。改进算法有效提高了模式分类的准确度, 具有应用于电气工程、环境工程等领域的可能性。但本文仅局限于算法的提出和算法应用思路的分析, 具体工程应用及算法速度问题的解决有待于后继研究来实现。

摘要：本文在自组织竞争网络的基础上, 结合迭代算法, 提出了一种改进的神经网络算法, 即迭代自组织神经网络算法。改进算法可有效提高模式分类的准确度, 但会很大程度降低程序的运行速度。本文最后, 将简要分析该改进算法在电力、环境工程等领域推广与应用的可能性。

关键词：自组织竞争网络,迭代算法,模式分类

参考文献

[1]汪海波, 张海臣, 段雪丽.基于MATLAB的自组织竞争神经网络聚类研究[J].邢台技术学院学报, 2005 (1) :46-47.

自适应迭代算法 篇10

目前针对装载机的轨迹跟踪控制问题已经取得一定程度的进展。柳波,等基于模糊自整定PID控制实现了装载机自动铲掘时的轨迹控制问题[1]。龚捷,等在考虑作业阻力补偿的情况下,基于计算力矩法设计装载机铲掘作业的轨迹跟踪控制器,实现铲斗轨迹跟踪控制[2]。

自适应迭代学习控制整合自适应控制和迭代学习控制的优点,基于PD反馈控制器增益自适应调节,迭代学习控制器修正控制力矩,该算法具有快速的渐进收敛性能[3]。拉格朗日法与牛顿-欧拉法是机构动力学建模的两种常见方法,封闭形式的拉格朗日方程是控制器设计和综合的有效工具,而递推形式的牛顿-欧拉方程则是实时控制、求解逆动力学问题的有效方法[4]。采用牛顿-欧拉法推导了装载机工作装置的二自由度动力学模型,在matlab环境下设计自适应迭代学习控制器并进行稳定性证明,对装载机工作装置的轨迹跟踪问题进行仿真控制,以期达到高精度追踪。

1 基于牛顿-欧拉法的动力学方程

将装载机工作装置机器人化处理,视为平面二连杆操作臂机构,建立如图1所示的D-H坐标系。设定动臂与车架铰接点为O0,铲斗与动臂铰接点为O1,铲斗斗尖位置为O2,分别以O0、O1、O2为原点建立坐标系O0x0y0z0、O1x1y1z1、O2x2y2z2。

图1中,θ1、θ2分别为动臂和铲斗相对于车架和动臂的转角;l1=lO0O1、l2=lO1O2分别为动臂和铲斗长度;m1、m2分别为铲斗和动臂的质量;假定动臂和铲斗的质心处于连杆O0O1、O1O2的中心位置。

1.1 牛顿-欧拉动力学方程

根据旋转关节的运动推导关节力矩的过程可细分为外推法和内推法两个步骤。外推法采用牛顿欧拉方程,按连杆编号从小到大的顺序逐次计算连杆的速度和加速度;内推法则按连杆编号从大到小的顺序计算连杆间的力、力矩[5]。运动学方程如下。

外推法求解过程:i=0,1,…,n。

内推法求解过程:i=n,n-1,…,1。

式中:iωi、${}^i\omega {}_{\dot{i}}$、ivi、${}^{i}b{}_{\dot{i}}$为连杆i的角速度、角加速度、线速度、线加速度相对于坐标系{i}的表达式;i+1iR和${}_{i+1}^i$R为坐标系{i}和{i+1}之间的姿态变换矩阵,且i+1iR${}_{i+1}^i$R=I;i$\widehat{{\rm Z}}$i为Zi轴方向上的单位矢量;iPi+1为节点i到节点i+1的矢量;iPCi为节点i到连杆i质心位置处的矢量;CiIi为连杆i质心处的惯性张量;mi为连杆i的质量;为连杆i质心上的惯性力和力矩相对于坐标系{i}的表达式;为连杆i-1作用在连杆i上的力和力矩相对于坐标系{i}的表达式;τi为关节i上的旋转力矩。

1.2 装载机工作装置的动力学计算

以图1中装载机工作装置的二连杆机器人操作臂为研究对象,推导其牛顿-欧拉动力学方程。假设装载机铲斗斗尖处不受力的作用f3=0,n3=0;车架固定不旋转有$\omega {}_0=0,\stackrel{\cdot }{{\displaystyle \omega {}_0}}=0$;计及重力因素[6]的影响,在式(1)中令$v{}_{\stackrel{\cdot }{{\displaystyle i+1}}}{}_0=G-Gy{}_0‚g=9.8\text{m}/\text{s}{}^{-2}$。相邻坐标系之间的姿态变换矩阵为(其中ci+1=cos(θi+1)),(si+1=sin(θi+1))。

$\begin{array}{l}{}_{i+1}^{i}R=\left[\begin{array}{ccc}c{}_{i+1}& -s{}_{i+1}& 0\\ s{}_{i+1}& c{}_{i+1}& 0\\ 0& 0& 1\end{array}\right],\\ {}_i^{i+1}R=\left[\begin{array}{ccc}c{}_{i+1}& s{}_{i+1}& 0\\ -s{}_{i+1}& c{}_{i+1}& 0\\ 0& 0& 1\end{array}\right],i=0,1(3)\end{array}$

连杆在坐标系{i}中的惯性张量为:

根据平行移轴定理[7]得到以刚体质心{C}为原点的坐标系与{i}坐标系间的转换式

$\{\begin{array}{l}{}^i{\rm I}{}_{xx}={}^{C{}_i}{\rm I}{}_{xx}+m(y{}_c{}^2+z{}_c{}^2)\\ {}^i{\rm I}{}_{yy}={}^{C{}_i}{\rm I}{}_{yy}+m(x{}_c{}^2+z{}_c{}^2),{}^i{\rm I}{}_{zz}=\\ {}^{C{}_i}{\rm I}{}_{zz}+m(x{}_c{}^2+y{}_c{}^2)\\ {}^i{\rm I}{}_{xy}={}^{C{}_i}{\rm I}{}_{xy}-mx{}_{c}y{}_c,{}^i{\rm I}{}_{xz}=\\ {}^{C{}_i}{\rm I}{}_{xz}-mx{}_{c}z{}_c,{}^i{\rm I}{}_{yz}={}^{C{}_i}{\rm I}{}_{yz}-my{}_{c}z{}_c\end{array}(6)$

式(6)中pc=[xcyczc]T表示连杆质心在坐标系{i}中的位置矢量。

将式(3)—式(6)代入式(1)和式(2)计算装载机工作装置的动力学方程如下(式子太长未列出)。

$\begin{array}{l}{}^{C{}_1}{\rm I}{}_1=\left[\begin{array}{ccc}0& 0& 0\\ 0& \frac{m{}_{1}l{}_1{}^2}{12}& 0\\ 0& 0& \frac{m{}_{1}l{}_1{}^2}{12}\end{array}\right],\\ {}^{C{}_2}{\rm I}{}_2=\left[\begin{array}{ccc}0& 0& 0\\ 0& \frac{m{}_{2}l{}_2{}^2}{12}& 0\\ 0& 0& \frac{m{}_{2}l{}_2{}^2}{12}\end{array}\right](7)\end{array}$

设$\overrightarrow{i}$、$\overrightarrow{j}$、$\overrightarrow{k}$分别为坐标系O0-x0y0z0中的单位向量

$\{\begin{array}{l}{}^1\omega {}_1=\dot{\theta }{}_1\overrightarrow{k}\\ {}^1\dot{\omega }{}_1=\ddot{\theta }{}_1\overrightarrow{k}\\ {}^1\dot{v}{}_1=gs{}_1\overrightarrow{i}+gc{}_1\overrightarrow{j}\\ {}^1\dot{v}{}_{C{}_1}=\left(-\frac{l{}_1}{2}\ddot{\theta }{}_1+gs{}_1\right)\overrightarrow{i}+\left(\frac{l{}_1}{2}\ddot{\theta }{}_1+gc{}_1\right)\overrightarrow{j}\\ {}^2\omega {}_2=(\dot{\theta }{}_1+\dot{\theta }{}_2)\overrightarrow{k}\\ {}^1\dot{\omega }{}_2=(\ddot{\theta }{}_1+\ddot{\theta }{}_2)\overrightarrow{k}\\ {}^2\dot{v}{}_2=(l{}_1\ddot{\theta }{}_{1}s{}_2-l{}_1\dot{\theta }{}_1{}^{2}c{}_2+gs{}_{12})\overrightarrow{i}+(l{}_1\ddot{\theta }{}_{1}c{}_2+\\ l{}_1\dot{\theta }{}_1^{2}s{}_2+gc{}_{12})\overrightarrow{j}\\ {}^2\dot{v}C{}_2=(l{}_1\ddot{\theta }{}_{1}s{}_2-l{}_1\dot{\theta }{}_1^{2}c{}_2+gs{}_{12}-\frac{l{}_2}{2}(\dot{\theta }{}_1+\dot{\theta }{}_2){}^2)\overrightarrow{i}+\\ (l{}_1\ddot{\theta }{}_{1}c{}_2+l{}_1\dot{\theta }{}_1^{2}s{}_2+gc{}_{12}+\frac{l{}_2}{2}(\ddot{\theta }{}_1+\ddot{\theta }{}_2))\overrightarrow{j}\\ {}^1{\rm N}{}_1=\frac{m{}_{1}l{}_1^2}{12}\ddot{\theta }{}_1\overrightarrow{k}\\ {}^{1}F{}_1=\left(-\frac{m{}_{1}l{}_1}{2}\ddot{\theta }{}_1+m{}_{1}gs{}_1\right)\overrightarrow{i}+\left(\frac{m{}_{1}l{}_1}{2}\ddot{\theta }{}_1+m{}_{1}gc{}_1\right)\overrightarrow{j}\\ {}^2{\rm N}{}_2=\frac{m{}_{2}l{}_2^2}{12}(\ddot{\theta }{}_1+\ddot{\theta }{}_2)\overrightarrow{k}\\ {}^{2}F{}_2=m{}_2(l{}_1\ddot{\theta }{}_{1}s{}_2-l{}_1\dot{\theta }{}_1^{2}c{}_2+gs{}_{12}-\frac{l{}_2}{2}(\dot{\theta }{}_1+\dot{\theta }{}_2){}^2)\overrightarrow{i}+\\ m{}_2(l{}_1\ddot{\theta }{}_{1}c{}_2+l{}_1\dot{\theta }{}_1^{2}s{}_2+gc{}_{12}+\frac{l{}_2}{2}(\ddot{\theta }{}_1+\ddot{\theta }{}_2))\overrightarrow{j}\end{array}(8)$

$\{\begin{array}{l}\tau {}_1=\left(\frac{m{}_{1}l{}_1{}^2}{3}+\frac{m{}_{2}l{}_2{}^2}{3}+m{}_{2}l{}_1{}^2+m{}_{2}l{}_{1}l{}_{2}c{}_2\right)\ddot{{\displaystyle \theta {}_1}}+\\ \left(\frac{m{}_{2}l{}_2{}^2}{3}+\frac{m{}_{2}l{}_{1}l{}_{2}c{}_2}{2}\right)\theta {}_2{}^{\ddot{2}}-\frac{m{}_{2}l{}_{1}l{}_{2}s{}_2}{2}\stackrel{\cdot }{{\displaystyle \theta {}_2}}{}^2-\\ m{}_{2}l{}_{1}l{}_{2}s{}_2\stackrel{\cdot }{{\displaystyle \theta {}_1}}\stackrel{\cdot }{{\displaystyle \theta {}_2}}+\frac{m{}_{2}gl{}_2}{2}c{}_{12}+\left(\frac{m{}_1}{2}+m{}_2\right)gl{}_{1}c{}_1\\ \tau {}_2=\left(\frac{m{}_{2}l{}_2{}^2}{3}+\frac{m{}_{2}l{}_{1}l{}_{2}c{}_2}{2}\right)\ddot{{\displaystyle \theta {}_1}}+\frac{m{}_{2}l{}_2{}^2}{3}\ddot{{\displaystyle \theta {}_2}}+\\ \frac{m{}_{2}l{}_{2}s{}_2}{2}\stackrel{\cdot }{{\displaystyle \theta {}_1}}{}^2+\frac{m{}_{2}l{}_{2}gc{}_{12}}{2}\end{array}(9)$

2 工作装置轨迹的自适应迭代学习控制

2.1 问题描述

考虑实际中存在不确定性、摩擦及外部干扰的影响,式(9)描述的装载机工作装置动力学方程可写为

${\rm M}(q{}_i(t))\ddot{{\displaystyle q{}_i}}(t)+C(q{}_i(t),\ddot{{\displaystyle q{}_i}}(t))\ddot{{\displaystyle q{}_i}}(t)+G(q{}_i(t))=\tau {}_i(t)+d{}_i(10)$

式(10)中:qi(t)=θi(t),i=1,2,t∈[0,T]为时间变量,i∈Z+为迭代次数,$q{}_i\in \mathbf{R}{}^n,\stackrel{\cdot }{{\displaystyle q{}_i}}\in \mathbf{R}{}^n,\ddot{{\displaystyle q{}_i}}\in \mathbf{R}{}^n$表示第i次迭代的关节角度、角速度、角加速度量。M(qi)∈Rn×n为惯性矩阵,$C(q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}})\stackrel{\cdot }{{\displaystyle q{}_i}}\in \mathbf{R}{}^n$表示科氏力和离心力量,G(qi)∈Rn为重力项,τi(t)∈Rn为作用在关节位置处的控制力或力矩,di∈Rn为未建模动态和扰动。

2.2 收敛性分析

控制目的是∀t∈[0,T]及∀i∈Z+,设计有界的自适应迭代学习控制器τi(t)来保证$q{}_i\in \mathbf{R}{}^n,\stackrel{\cdot }{{\displaystyle q{}_i}}\in \mathbf{R}{}^n,\stackrel{\cdot \cdot }{{\displaystyle q{}_i}}\in \mathbf{R}{}^n$有界,且当i→∞时,$q{}_i(t),\stackrel{\cdot }{{\displaystyle q{}_i}}(t)$收敛到理想的参考位置$q{}_d(t),\stackrel{\cdot }{{\displaystyle q{}_d}}(t)$。为实现控制目的,给出以下基本假设:

(A1) qd(t)为可实现的关节参考轨迹角度。

(A2) 对于∀t∈[0,T]及$\forall i\in \mathbf{{\rm Z}}{}_{+}‚q{}_i(t),\stackrel{\cdot }{{\displaystyle q{}_i}}(t),\stackrel{\cdot \cdot }{{\displaystyle q{}_i}}(t)‚d{}_i(t)$有界。

(A3) ∀i∈Z+,初始条件满足:

$\stackrel{\cdot }{{\displaystyle q{}_d}}(0)-\stackrel{\cdot }{{\displaystyle q{}_i}}(0)=q{}_d(0)-q{}_i(0)=0$。

且满足如下4个基本特征:

(B1) M(qi)∈Rn×n为有界正定对称矩阵;

$(\text{B}{}_2)\dot{{\rm M}}(q{}_i)-2C(q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}})\in \mathbf{R}{}^{n\times n}$为对称矩阵,且$\omega {}^{\text{Τ}}(\dot{{\rm M}}(q{}_i)-2C(q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}}))\omega =0,\forall \omega \in \mathbf{R}{}^n$;

$(\text{B}{}_3)G(q{}_i)+C(q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}})\stackrel{\cdot }{{\displaystyle q{}_d}}(t)=\psi (q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}})\xi {}^{\text{Τ}}(t),\psi (q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}})\in \mathbf{R}{}^{n\times (m-1)}$为已知矩阵,ξ(t)∈Rm-1为未知向量;

$(\text{B}{}_4)\forall q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}}‚t\in [0,{\rm T}]‚i\in \mathbf{{\rm Z}}{}_{+}‚||C(q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}})||\le k{}_c||\stackrel{\cdot }{{\displaystyle q{}_i}}||‚||G(q{}_i)||\le k{}_g,k{}_c$、kg均为正实数。

定理 定义关节位置跟踪误差和速度误差分别为:$e{}_i(t)=q{}_d(t)-q{}_i(t),\stackrel{\cdot }{{\displaystyle e{}_i}}(t)=\dot{q}{}_d(t)-\dot{q}{}_i(t)$,采用控制律

$\{\begin{array}{l}\tau {}_i(t)={\rm K}{}_{{\rm P}}e{}_i(t)+{\rm K}{}_D\stackrel{\cdot }{{\displaystyle e{}_i}}(t)+\phi (q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}},\stackrel{\cdot }{{\displaystyle e{}_i}}(t))\stackrel{⌢}{\theta {}_i}(t)(11\text{a})\\ \stackrel{⌢}{\theta {}_i}(t)=\stackrel{⌢}{\theta }{}_{i-1}(t)+\Gamma \phi {}^{\text{Τ}}(q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}},\stackrel{\cdot }{{\displaystyle e{}_i}}(t))\stackrel{\cdot }{{\displaystyle e{}_i}}(t)(11\text{b})\end{array}$

式(11)中:$\stackrel{⌢}{\theta }{}_{i-1}(0)=0‚\phi (q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}},\stackrel{\cdot }{{\displaystyle e{}_i}})\in \mathbf{R}{}^{n\times n}$,且$\phi (q{}_i,\stackrel{\cdot }{{\displaystyle q{}_i}},\stackrel{\cdot }{{\displaystyle e{}_i}})=[\psi (q{}_i,\dot{q}{}_i)\mathrm{sgn}(\stackrel{\cdot }{{\displaystyle e{}_i}})]$,矩阵KP、KD、Γ均为对称正定矩阵,则ei(t)、$\stackrel{\cdot }{{\displaystyle e{}_i}}(t)$有界,且

$\underset{i\to \infty }{{\displaystyle \lim }}e{}_i(t)=\underset{i\to \infty }{{\displaystyle \lim }}\stackrel{\cdot }{{\displaystyle e{}_i}}(t)=0,\forall t\in [0,{\rm T}]$。

证明 第i次迭代时,构造Lyapunov函数

$W{}_i(t)=V{}_i(t)+\frac{1}{2}{\displaystyle {\int }_0^t\stackrel{⌢}{\theta {}_i}}(t)\Gamma {}^{-1}\stackrel{⌢}{\theta {}_i}(t)\text{d}\tau (12)$

式(12)中,设定$\stackrel{⌢}{\theta }{}_i(t)=\theta (t)-\stackrel{⌢}{\theta }{}_i(t)‚\theta (t)=[\xi {}^{\text{Τ}}(t)\beta ]{}^{\text{Τ}},\widehat{\theta }{}_i(t)=[\widehat{\xi }{}_i{}^{\text{Τ}}(t)\widehat{\beta }{}_i(t)]{}^{\text{Τ}}$为θ(t)的估计值,由假设(A1)、和基本特征(B1)知,

$||{\rm M}(q{}_i)\stackrel{\cdot \cdot }{{\displaystyle q}}{}_d-d{}_i||\le \beta ‚\beta >0$

$\begin{array}{l}V{}_i(e{}_{\dot{i}}(t),e{}_i(t))=\frac{1}{2}e{}_i{}^{\dot{\text{Τ}}}(t){\rm M}(q{}_i)e{}_{\dot{i}}(t)+\\ \frac{1}{2}e{}_i(t){}^{\text{Τ}}{\rm K}{}_{p}e{}_i(t)(13)\end{array}$

由于$\overline{\theta {}_i}=-\theta +\stackrel{⌢}{\theta }{}_i-\stackrel{⌢}{\theta }{}_{i+1}+\theta =-\stackrel{⌢}{\theta }{}_i+\stackrel{⌢}{\theta }{}_{i+1}$,

即$\stackrel{⌢}{\theta }{}_{i+1}=\overline{\theta }{}_i+\stackrel{⌢}{\theta }{}_i$,则

$\begin{array}{l}\tilde{\theta }{}_i(t)\Gamma {}^{-1}\tilde{\theta }{}_i(t)-\tilde{\theta }{}_{i-1}(t)\Gamma {}^{-1}\tilde{\theta }{}_{i-1}(t)=\\ -2\overline{\theta }{}_i{}^{\text{Τ}}(t)\Gamma {}^{-1}\tilde{\theta }{}_i(t)-\overline{\theta }{}_i{}^{\text{Τ}}(t)\Gamma {}^{-1}\overline{\theta }{}_i(t)(14)\end{array}$

2.2.1 Wi(t)的非递增性证明

$\text{Δ}W{}_i=W{}_i-W{}_{i-1}=V{}_i-V{}_{i-1}-\frac{1}{2}{\displaystyle {\int }_0^t(}\overline{\theta }{}_i{}^{\text{Τ}}\Gamma {}^{-1}\overline{\theta }{}_i+$

$2\overline{\theta }{}_{i-1}{}^{\text{Τ}}(t)\Gamma {}^{-1}\overline{\theta }{}_{i-1}(t))\text{d}\tau$ (15)

式(15)中$\overline{\theta }{}_i=\stackrel{⌢}{\theta }{}_i-\stackrel{⌢}{\theta }{}_{i-1}$。由于${\displaystyle {\int }_0^t\dot{V}}(t)\text{d}\tau =V{}_i(t)-V{}_i(0)$,即$V{}_i(t)=V{}_i(0)+{\displaystyle {\int }_0^t\dot{V}}(t)\text{d}\tau$,

又由于

$\{\begin{array}{l}V{}_{\dot{i}}(t)=e{}_i{}^{\dot{\text{Τ}}}(t){\rm M}(q{}_i)e{}_{\stackrel{\cdot \cdot }{{\displaystyle i}}}(t)+\frac{1}{2}e{}_i{}^{\dot{\text{Τ}}}(t)\times \\ \dot{{\rm M}}(q{}_i)e{}_{\dot{i}}(t)+e{}_i{}^{\dot{\text{Τ}}}(t){\rm K}{}_{p}e{}_i\\ V{}_i(e{}_{\dot{i}}(t),e{}_i(t))=\\ V{}_i(e{}_{\dot{i}}(0),e{}_i(0))+\\ {\displaystyle {\int }_0^t[}e{}_i{}^{\dot{\text{Τ}}}{\rm M}(q{}_i)e{}_{\stackrel{\cdot \cdot }{{\displaystyle i}}}+\frac{1}{2}e{}_i{}^{\dot{\text{Τ}}}\dot{{\rm M}}(q{}_i)e{}_{\dot{i}}+e{}_i{}^{\dot{\text{Τ}}}{\rm K}{}_{p}e{}_{\dot{i}}]\text{d}\tau \end{array}(16)$

由式(10)和特性(B2)和特性(B3)得

$\{\begin{array}{l}e{}_i{}^{\dot{\text{Τ}}}{\rm M}e{}_{\stackrel{\cdot \cdot }{{\displaystyle i}}}=e{}_i{}^{\dot{\text{Τ}}}{\rm M}(q{}_{\stackrel{\cdot \cdot }{{\displaystyle d}}}-q{}_{\stackrel{\cdot \cdot }{{\displaystyle i}}})=e{}_i{}^{\dot{\text{Τ}}}{\rm M}q{}_{\stackrel{\cdot \cdot }{{\displaystyle d}}}-\\ e{}_i{}^{\dot{\text{Τ}}}(-Cq{}_{\dot{i}}-G+\tau {}_i+d{}_i)\\ \frac{1}{2}e{}_i{}^{\dot{\text{Τ}}}\dot{{\rm M}}(q{}_i)e{}_{\dot{i}}=e{}_i{}^{\dot{\text{Τ}}}Ce{}_{\dot{i}}=e{}_i{}^{\dot{\text{Τ}}}C(q{}_{\dot{d}}-q{}_{\dot{i}})=\\ e{}_i{}^{\dot{\text{Τ}}}Cq{}_{\dot{d}}-e{}_i{}^{\dot{\text{Τ}}}Cq{}_{\dot{i}}\end{array}(17)$

由已知得

$\begin{array}{l}e{}_i{}^{\dot{\text{Τ}}}({\rm M}(q{}_i)q{}_{\stackrel{\cdot \cdot }{{\displaystyle d}}}-d{}_i)\le \\ \beta ||e{}_i{}^{\dot{\text{Τ}}}||=e{}_i{}^{\dot{\text{Τ}}}\beta \mathrm{sgn}(e{}_{\dot{i}})\psi (q{}_i,q{}_{\dot{i}})\xi {}^{\text{Τ}}+\beta \mathrm{sgn}(e{}_{\dot{i}})=\\ [\psi (q{}_i,q{}_{\dot{i}})\mathrm{sgn}(e{}_{\dot{i}})][\xi {}^{\text{Τ}}\beta ]{}^{\text{Τ}}=\text{ϕ}(q{}_i,q{}_{\dot{i}},e{}_{\dot{i}})\theta (18)\end{array}$

则

$\begin{array}{l}V{}_i(e{}_{\dot{i}}(t),e{}_i(t))=\\ V{}_i(e{}_{\dot{i}}(0),e{}_i(0))+{\displaystyle {\int }_0^{t}e}{}_i{}^{\dot{\text{Τ}}}({\rm M}(q{}_i)q{}_{\stackrel{\cdot \cdot }{{\displaystyle d}}}-d{}_i+\\ C(q{}_i,q{}_{\dot{i}})q{}_{\dot{d}}+G(q{}_i)+{\rm K}{}_{{\rm P}}e{}_i-\tau {}_i)\text{d}\tau \le V{}_i(e{}_{\dot{i}}(0),e{}_i(0))+{\displaystyle {\int }_0^{t}e}{}_i{}^{\dot{\text{Τ}}}(\psi (q{}_i,q{}_{\dot{i}})\xi {}^{\text{Τ}}+{\rm K}{}_{{\rm P}}e{}_i+\beta \mathrm{sgn}(e{}_{\dot{i}})-\\ \tau {}_i)\text{d}\tau \le V{}_i(e{}_{\dot{i}}(0),e{}_i(0))+\\ {\displaystyle {\int }_0^{t}e}{}_i{}^{\dot{\text{Τ}}}(\text{ϕ}(q{}_i,q{}_{\dot{i}},e{}_i)\theta +{\rm K}{}_{{\rm P}}e{}_i-\tau {}_i)\text{d}\tau (19)\end{array}$

将控制律式(11a)代入式(19),得

$\begin{array}{l}V{}_i(e{}_{\dot{i}}(t),e{}_i(t))\le V{}_i(e{}_{\dot{i}}(0),e{}_i(0))+\\ {\displaystyle {\int }_0^{t}e}{}_i{}^{\dot{\text{Τ}}}(\text{ϕ}(q{}_i,q{}_{\dot{i}},e{}_i)\tilde{\theta }{}_i-{\rm K}{}_{D}e{}_{\dot{i}})\text{d}\tau (20)\end{array}$

据式(11b)得:$\overline{\theta }{}_i{}^{\text{Τ}}(t)=(\Gamma \text{ϕ}{}^{\text{Τ}}e{}_{\dot{i}}){}^{\text{Τ}}=e{}_i{}^{\dot{\text{Τ}}}\text{ϕ}\Gamma$,有

由假设(A2)得$V{}_i(\stackrel{\cdot }{{\displaystyle e{}_i}}(0),e{}_i(0))=0$,将式(21)代入式(15),有

$\begin{array}{l}\text{Δ}W{}_i=-V{}_{i-1}+V{}_i-\frac{1}{2}{\displaystyle {\int }_0^t(}\overline{\theta }{}_i{}^{\text{Τ}}\Gamma {}^{-1}\overline{\theta }{}_i+2\overline{\theta }{}_i{}^{\text{Τ}}\Gamma {}^{-1}\tilde{\theta }{}_i)\text{d}\tau \le \\ -V{}_{i-1}+{\displaystyle {\int }_0^{t}e}{}_i{}^{\dot{\text{Τ}}}(\text{ϕ}\tilde{\theta }{}_i-{\rm K}{}_{D}e{}_{\dot{i}})\text{d}\tau -\frac{1}{2}{\displaystyle {\int }_0^t(}e{}_i{}^{\dot{\text{Τ}}}\text{ϕ}\Gamma \text{ϕ}{}^{\text{Τ}}e{}_{\dot{i}}+2e{}_i{}^{\dot{\text{Τ}}}\text{ϕ}\tilde{\theta }{}_i)\text{d}\tau \le -V{}_{i-1}-\frac{1}{2}{\displaystyle {\int }_0^{t}e}{}_i{}^{\dot{\text{Τ}}}(\text{ϕ}\Gamma \text{ϕ}{}^{\text{Τ}}+\\ 2{\rm K}{}_D)e{}_{\dot{i}}\text{d}\tau \le 0(22)\end{array}$

因Vi-1、Γ、KD均为正定阵,ΔWi≤0,因此Wi为非递增数列,得结论如下:若W0有界,那么Wi必定有界。

2.2.2 W0(t)的连续有界性证明

由式(12)和式(20)得

$\begin{array}{l}W{}_{\dot{0}}(t)\le e{}_0{}^{\dot{\text{Τ}}}(\text{ϕ}(q{}_0,q{}_{\dot{0}},q{}_{\stackrel{\cdot \cdot }{{\displaystyle 0}}})\tilde{\theta }{}_0-{\rm K}{}_{D}e{}_{\dot{0}})+\\ \frac{1}{2}\tilde{\theta }{}_0{}^{\text{Τ}}\Gamma {}^{-1}\tilde{\theta }{}_0(23)\end{array}$

由于$\stackrel{⌢}{\theta }{}_{-1}(t)=0$,根据控制律式(11b)得

$\stackrel{⌢}{\theta }{}_0(t)=\Gamma \phi {}^{\text{Τ}}(q{}_0,\stackrel{\cdot }{{\displaystyle q{}_0}},\stackrel{\cdot }{{\displaystyle e{}_0}})\stackrel{\cdot }{{\displaystyle e{}_0}}(t)$,有

$\begin{array}{l}W{}_{\dot{0}}(t)\le -e{}_0{}^{\dot{\text{Τ}}}{\rm K}{}_{D}e{}_{\dot{0}}+\left(\widehat{\theta }{}_0{}^{\text{Τ}}+\frac{1}{2}\tilde{\theta }{}_0{}^{\text{Τ}}\right)\Gamma {}^{-1}\tilde{\theta }{}_0=\\ -e{}_0{}^{\dot{\text{Τ}}}{\rm K}{}_{D}e{}_{\dot{0}}-\frac{1}{2}\tilde{\theta }{}_0{}^{\text{Τ}}\Gamma {}^{-1}\tilde{\theta }{}_0+\theta {}^{\text{Τ}}\Gamma {}^{-1}\tilde{\theta }{}_0(24)\end{array}$

因为$\theta {}^{\text{Τ}}\Gamma {}^{-1}\stackrel{⌢}{\theta }{}_0\le {\rm K}||\Gamma {}^{-1}\stackrel{⌢}{\theta }{}_0||{}^2+\frac{1}{4{\rm K}}||\theta ||{}^2,{\rm K}>0$,可得

$\stackrel{\cdot }{{\displaystyle W{}_0}}(t)\le -\beta {}_1||\stackrel{\cdot }{{\displaystyle e{}_0}}||{}^2-\beta {}_2||\stackrel{⌢}{\theta }{}_0||{}^2+\frac{1}{4{\rm K}}||\theta ||{}^2(25)$

$\begin{array}{l}\text{式}(25)\text{中}:\beta {}_1=\lambda {}_{\min }({\rm K}{}_D),\\ \beta {}_2=\frac{1}{2}\lambda {}_{\min }(\Gamma {}^{-1})-{\rm K}\lambda {}_{\max }^2(\Gamma {}^{-1}),{\rm K}\le \frac{\lambda {}_{\min }(\Gamma {}^{-1})}{2\lambda {}_{\max }^2(\Gamma {}^{-1})}‚\text{因}\text{此}\text{有}\end{array}$

$\stackrel{\cdot }{{\displaystyle W{}_0}}(t)\le \frac{1}{4{\rm K}}||\theta ||{}^2‚\forall t\in [0,{\rm T}](26)$

因为θ(t)连续有界,推出W0(t)亦连续有界。

2.2.3 Wi(t)的连续有界性证明

Wi(t)可表示为$W{}_i(t)=W{}_0+{\displaystyle \sum _{j=1}^i\text{Δ}}W{}_j$。

则由式(15)可得

$\begin{array}{l}W{}_i\le W{}_0-{\displaystyle \sum _{j=1}^{i}V}{}_{j-1}\le W{}_0-{\displaystyle \sum _{j=1}^{i}V}{}_{j-1}\le W{}_0-\\ \frac{1}{2}{\displaystyle \sum _{j=1}^{i}e}{}_{j-1}{}^{\text{Τ}}{\rm K}{}_{{\rm P}}e{}_{j-1}-\frac{1}{2}{\displaystyle \sum _{j=1}^{i}e}{}_{j-1}{}^{\dot{\text{Τ}}}{\rm K}{}_{{\rm P}}e{}_{\stackrel{\cdot }{{\displaystyle j-1}}}(27)\end{array}$

$\begin{array}{l}\text{那}\text{么}({\displaystyle \sum _{j=1}^{i}e}{}_{j-1}{}^{\text{Τ}}{\rm K}{}_{{\rm P}}e{}_{j-1}+{\displaystyle \sum _{j=1}^{i}e}{}_{j-1}{}^{\dot{\text{Τ}}}{\rm K}{}_{{\rm P}}e{}_{\stackrel{\cdot }{{\displaystyle j-1}}})\le \\ 2(W{}_0-W{}_i)\le 2W{}_0,\text{即}W{}_i(t)\text{有}\text{界}。\text{可}\text{得}\text{结}\text{论}:\end{array}$

$\underset{i\to \infty }{{\displaystyle \lim }}e{}_i(t)=\underset{i\to \infty }{{\displaystyle \lim }}\stackrel{\cdot }{{\displaystyle e{}_i}}(t)=0,\forall t\in [0,{\rm T}]$,定理得证。

3 仿真研究

对于式(10)所示的装载机工作装置的动力学方程及式(11)所设计的控制器,仿真研究如下。设置动臂、铲斗质量分别为m1=1 793.080 7 kg、m2=2 706.447 8 kg;动臂、铲斗长度l1=2.668 2 m、l2=1.262 1 m。对式(10)描述的系统,可写为

$\left[\begin{array}{cc}m{}_{11}& m{}_{12}\\ m{}_{21}& m{}_{22}\end{array}\right]\left[\begin{array}{l}\stackrel{\cdot \cdot }{{\displaystyle \theta {}_1}}\\ \stackrel{\cdot \cdot }{{\displaystyle \theta {}_2}}\end{array}\right]+\left[\begin{array}{cc}c{}_{11}& c{}_{12}\\ c{}_{21}& c{}_{22}\end{array}\right]\left[\begin{array}{l}\stackrel{\cdot }{{\displaystyle \theta {}_1}}\\ \stackrel{\cdot }{{\displaystyle \theta {}_2}}\end{array}\right]+\left[\begin{array}{c}G{}_1\\ G{}_2\end{array}\right]=\left[\begin{array}{c}\tau {}_1\\ \tau {}_2\end{array}\right](28)$