高级检索

ISSN1001-3806CN51-1125/TN 网站地图

留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

基于FPGA的光纤光斑中心定位算法研究

张睛 吴友宇

引用本文:
Citation:

基于FPGA的光纤光斑中心定位算法研究

    作者简介: 张睛(1991-), 女, 硕士研究生, 研究方向为嵌入式图像处理.
    通讯作者: 吴友宇, wuyouyu1@whut.edu.cn
  • 中图分类号: TN911.73;TN47

Locating algorithm of optical fiber spot center based on FPGA

    Corresponding author: WU Youyu, wuyouyu1@whut.edu.cn
  • CLC number: TN911.73;TN47

  • 摘要: 为了解决传统数字图像处理算法中数据运算量大、复杂度高、耗时长的问题,提出一种基于可编程门阵列(FPGA)光纤光斑中心定位的方法。采用数字信号处理系统,利用开发工具(DSP builder),设计了光斑图像预处理算法和边缘检测算法,用最小二乘法拟合光斑边界,采用流水线设计,增强了数据处理的并行能力,提高了处理速度。在Cyclone V平台上进行理论分析和实验验证,取得了光斑图像边界、中心坐标数据。结果表明,在保证对光斑中心定位的绝对误差小于0.1pixel的条件下,使用FPGA比计算机运算速度能提高21倍以上。该研究能够在FPGA平台上快速准确定位光斑中心。
  • Figure 1.  Spot image before and after filtering

    Figure 2.  Spot image before and after edge detection with Sobel operator

    Figure 3.  Structure diagram of decomposition module

    Table 1.  Deviation between FPGA and PC

    coordinate of center absolute deviation
    visual studio (150.06907, 166.097) (0.06907, 0.097)
    FPGA (150.06881, 166.048) (0.06881, 0.048)
    下载: 导出CSV

    Table 2.  Processing time

    FPGA visual studio proportion
    decomposition/μs 36.38 620 1/17.04
    solving/μs 37.06 995 1/26.85
    total/μs 73.44 1615 1/21.99
    下载: 导出CSV
  • [1]

    FAN Ch L, LIU Ch, DING G, et al.A method of circle curve fitting based on the cumulative error of the radius error[C]//International Conference on Computational Intelligence and Security. New York, USA: IEEE, 2015: 211-214.
    [2]

    WANG Q Q, LIU J, PENG Zh, et al.Measurement system for laser divergence angle based on LabView[J]. Chinese Journal of Lasers, 2012, 39(11):1108005(in Chinese). doi: 10.3788/CJL
    [3]

    ZHANG L. System of measuring laser spot[D]. Xi'an: Xidian University, 2010: 19-45(in Chinese).
    [4]

    WANG L L, HU Zh W, JI H X. Laser spot center location algorithm based on Gaussian fitting[J]. Journal of Applied Optics, 2012, 33(5):985-990(in Chinese). 
    [5]

    LIU H L, HOU W, FAN Y L, et al. An improved algorithm of laser spot center location[J]. Computer Measurement & Control, 2014, 22(1):139-141(in Chinese). 
    [6]

    LAN Zh L, YANG X F. Practical improvement of laser spot center location algorithm[J]. Computer Engineering, 2008, 34(6):7-9(in Chinese). 
    [7]

    NING Sh N, ZHU M, SUN H H, et al. Realization of improved sobel adaptive edge detection algorithm based on FPGA[J]. Chinese Journal of Liquid Crystals and Displays, 2014, 29(3):395-402(in Chinese). doi: 10.3788/YJYXS
    [8]

    SHI Y, CHENG X. Laser spot center detection based on the geometric feature[C]//International Symposium on Information Science & Engineering. New York, USA: IEEE, 2010: 322-325.
    [9]

    LIU J S, YUAN S C. Subpixel location algorithm for circular targets center based on Zernike moments and curvature[J]. Computer Engineering & Applications, 2010, 46(29):153-156. 
    [10]

    ZHONG X J. Research on FPGA calculation method of typical matrix decomposition[D]. Harbin: Harbin Institute of Technology, 2012: 14-66(in Chinese).
    [11]

    TAI Y G, LO C T D, PSARRIS K. Applying out-of-core QR decomposition algorithms on FPGA-based systems[C]//International Conference on Field Programmable Logic and Applications. New York, USA: IEEE, 2007: 86-91.
    [12]

    HAN B, YANG Z, ZHENG Y R. FPGA implementation of QR decomposition for MIMO-OFDM using four CORDIC cores[C]//IEEE International Conference on Communications. New York, USA: IEEE, 2013: 4556-4560.
    [13]

    WU G M, DOU Y, PETERSON G D. Blocking LU decomposition for FPGAs[C]//IEEE International Symposium on Field-Programmable Custom Computing Machines. New York, USA: IEEE, 2010: 109-112.
    [14]

    TANG B P, JIANG Y H, ZHANG X Ch. Feature extraction method of rolling bearing fault based on singular value decomposition-morphology filter and empirical mode decomposition[J]. Chinese Journal of Mechanical Engineering, 2010, 46(5):37-42(in Chinese). doi: 10.3901/JME.2010.05.037
    [15]

    WU G M, DOU Y, WANG M. A fine-grained parallel algorithm for the Cholesky decomposition[J]. Computer Engineering & Science, 2010, 32(9):102-106(in Chinese). 
    [16]

    GUO L, TANG Y H, ZHOU J, et al. Research on parallel architecture for LDLT decomposition -processor[J]. Computer Engineering, 2011, 37(21):241-243(in Chinese).
  • [1] 蔡旭明李笑刘玉县何春华林俊杰 . 基于灰度直方图的激光光斑中心定位算法. 激光技术, 2023, 47(2): 273-279. doi: 10.7510/jgjs.issn.1001-3806.2023.02.018
    [2] 王芳荣赵丁选廖宗建张宇 . 激光光斑中心空间定位方法的研究. 激光技术, 2005, 29(1): 87-89.
    [3] 王博刘晓东李君豪陈泰宇刘容麟 . 基于邻域特征的网点激光打孔定位算法研究. 激光技术, 2019, 43(5): 591-596. doi: 10.7510/jgjs.issn.1001-3806.2019.05.001
    [4] 张海庄姚梅雷萍李鹏曾庆平 . 远场激光光斑图像处理方法研究. 激光技术, 2013, 37(4): 460-463. doi: 10.7510/jgjs.issn.1001-3806.2013.04.010
    [5] 鲍鸿曾海涛白玉磊胡忠向志聪周延周申作春 . 基于概率密度最小二乘拟合的叶片后缘轮廓. 激光技术, 2016, 40(4): 555-559. doi: 10.7510/jgjs.issn.1001-3806.2016.04.021
    [6] 孙立环赵霄洋高凌妤李杏华 . 基于亚像素定位技术的激光光斑中心位置测量. 激光技术, 2017, 41(4): 511-514. doi: 10.7510/jgjs.issn.1001-3806.2017.04.011
    [7] 王杰飞刘洁瑜赵晗沈强 . 一种改进的激光光斑中心亚像素定位方法. 激光技术, 2015, 39(4): 476-479. doi: 10.7510/jgjs.issn.1001-3806.2015.04.010
    [8] 王国军黄亚新赵启林张冬冬 . 基于自适应区域的光斑中心鲁棒性研究. 激光技术, 2020, 44(5): 616-622. doi: 10.7510/jgjs.issn.1001-3806.2020.05.015
    [9] 张健李白燕 . 基于图论最小割集算法的图像分割研究. 激光技术, 2014, 38(6): 863-866. doi: 10.7510/jgjs.issn.1001-3806.2014.06.030
    [10] 周中亮周冰何永强王斌 . 成像型激光探测系统中光斑精确定位方法研究. 激光技术, 2008, 32(3): 248-251.
    [11] 邓晓鹏文伟 . 基于干涉的二值图像逻辑运算加密技术. 激光技术, 2010, 34(3): 401-404. doi: 10.3969/j.issn.1001-3806.2010.03.033
    [12] 李世阳禹延光叶会英付广春 . 一种半导体激光自混合效应模型参数的测量方法. 激光技术, 2005, 29(5): 519-521.
    [13] 谢珊珊王哲强黄河陈宝宝汪培李劲松 . 随机抽样一致性算法在激光光谱中的应用研究. 激光技术, 2017, 41(1): 133-137. doi: 10.7510/jgjs.issn.1001-3806.2017.01.027
    [14] 李玉瑶张婉怡刘喆李美萱付秀华S-on-1测量方式下薄膜激光损伤的累积效应. 激光技术, 2018, 42(1): 39-42. doi: 10.7510/jgjs.issn.1001-3806.2018.01.008
    [15] 杨友良刘爱旭马翠红连畅 . 基于红外CCD的钢水红外测温模型分析. 激光技术, 2018, 42(4): 562-566. doi: 10.7510/jgjs.issn.1001-3806.2018.04.024
    [16] 张春光张玉钧韩道文刘文清陈臻懿 . 机动车颗粒物的激光雷达监测. 激光技术, 2009, 33(2): 130-133.
    [17] 李文龙戈海龙任远成巍 . 图像处理技术在激光熔池温度检测的应用. 激光技术, 2018, 42(5): 599-604. doi: 10.7510/jgjs.issn.1001-3806.2018.05.004
    [18] 汤敏王惠南 . 激光扫描共聚焦显微镜图像的计算机处理. 激光技术, 2007, 31(5): 558-560.
    [19] 张羽鹏王开福 . LabVIEW和MATLAB在电子散斑干涉图像处理中的应用. 激光技术, 2009, 33(6): 582-585,589. doi: 10.3969/j.issn.1001-3806.2009.06.007
    [20] 冯煦张瑞瑛周萍李松 . 大功率半导体线激光图像处理方法研究. 激光技术, 2010, 34(5): 624-627. doi: 10.3969/j.issn.1001-3806.2010.O5.013
  • 加载中
图(3) / 表(2)
计量
  • 文章访问数:  4820
  • HTML全文浏览量:  3724
  • PDF下载量:  331
  • 被引次数: 0
出版历程
  • 收稿日期:  2016-10-21
  • 录用日期:  2016-12-08
  • 刊出日期:  2017-09-25

基于FPGA的光纤光斑中心定位算法研究

    通讯作者: 吴友宇, wuyouyu1@whut.edu.cn
    作者简介: 张睛(1991-), 女, 硕士研究生, 研究方向为嵌入式图像处理
  • 武汉理工大学 信息工程学院, 武汉 430070

摘要: 为了解决传统数字图像处理算法中数据运算量大、复杂度高、耗时长的问题,提出一种基于可编程门阵列(FPGA)光纤光斑中心定位的方法。采用数字信号处理系统,利用开发工具(DSP builder),设计了光斑图像预处理算法和边缘检测算法,用最小二乘法拟合光斑边界,采用流水线设计,增强了数据处理的并行能力,提高了处理速度。在Cyclone V平台上进行理论分析和实验验证,取得了光斑图像边界、中心坐标数据。结果表明,在保证对光斑中心定位的绝对误差小于0.1pixel的条件下,使用FPGA比计算机运算速度能提高21倍以上。该研究能够在FPGA平台上快速准确定位光斑中心。

English Abstract

    • 随着科学技术的进步,光纤被广泛应用在各个领域,光纤光斑中心位置检测是光电测量手段中的关键技术,其准确性、精度和抗干扰性直接影响测量结果。目前光斑中心定位算法主要有Hough变换算法、重心法和拟合法。鉴于Hough变换算法需要逐点进行投票、记录,花费时间较长;重心法则要求光斑图像光强分布均匀,且其受光斑形状的影响较大,会产生误差;而拟合法时间复杂度较小、运算精度高[1],但容易受噪声影响。

      随着集成电路的发展,可编程门阵列(field-programmable gate array, FPGA)器件有着以下特点:优秀的可编程能力、并行执行能力、硬件资源丰富、低功耗。因此FPGA被广泛应用于图像处理算法的领域。现有测量激光光束质量分析的系统都已经相当成熟,但都是通过计算机对激光光斑图像进行处理来计算光斑的参量。例如,参考文献[2]中基于LabVIEW平台开发了一套快速检测激光发散角的测量系统。而参考文献[3]中将FPGA中采集到的激光光斑数据经USB线与上位机通信,并对数据二次处理,得到光斑的能量分布图,并解析出光斑的中心。参考文献[4]中对Hough变换法、重心法以及高斯拟合算法定位光斑中心相比较,得到高斯拟合方法效果最好,误差小于0.1pixel。参考文献[5]中对传统最小二乘法拟合圆以及改进的最小二乘法拟合圆进行比较,改进的圆拟合方法处理得到的中心坐标误差在0.1pixel之内。因此,为保证误差小于0.1pixel并能够提高计算效率,研究基于FPGA平台处理光斑数据的系统是非常必要的。针对这一问题,相对于传统软件进行光斑图像处理的算法,本文中提出一种基于FPGA的光纤光斑中心定位的算法,提高了运算速率,并且增加了图像预处理部分,可高效地过滤图像中的噪声,采用Sobel算子提取光斑轮廓,获得光斑边界像素坐标,再利用代数距离的最小二乘逼近光斑的轮廓,从而得到光斑的中心位置坐标,为求取光斑质量参量奠定基础。

    • 目前,常见的中心定位方法有Hough变换、重心法和最小二乘拟合[6]。鉴于重心法和Hough变换法具有复杂度高和精度低的缺陷,因此,本文中选用最小二乘法拟合光斑边界。最小二乘法拟合算法对噪声敏感,所以,问题的根本在于提高该算法的抗干扰能力,以提高光斑中心定位的精度。

    • 中值滤波是一种非线性平滑技术,通常用于处理椒盐噪声。其实现过程为:首先对图像中的某个采样窗口中的所有数据进行排序,然后用排序后的中值作为该区域的中心位置的像素值。传统的中值滤波算法运算复杂度高,故本文中采用快速中值滤波算法:在3×3的滑动窗口中,先比较每列值中的最大值、最小值、中值;再取出3个最大值中的最小值A0,3个中值的中值A1,3个最小值中的最大值A2;最后取出A0, A1, A2的中值即为最终结果。与线性滤波相比,中值滤波可以保护图像的边缘信息并且有效滤除噪声。经滤波以后的图像如图 1所示。

      Figure 1.  Spot image before and after filtering

    • 图像经中值滤波后,需对图像进行阈值分割,即二值化处理——选取适当的阈值[7],用0或255将灰度图像分成黑与白两种颜色,本文中选取的阈值为图像中灰度最大值的一半,而后采用边缘检测算法提取光斑的轮廓。Sobel算子与Prewitt算子相比,对像素位置的影响做了加权,降低了噪声的影响,用其进行目标光斑的边缘检测。换言之,同时测量水平方向和竖直方向的灰度差,将得到的数据叠加。经Sobel检测后,会得到所需的边缘像素坐标,再对数据进行处理。经过Sobel算子边缘检测后的图像如图 2所示。

      Figure 2.  Spot image before and after edge detection with Sobel operator

    • 对图像数据进行预处理后,得到平滑的光斑边缘,进而采用曲线拟合对光斑的中心精确定位。圆可以看成是一个特殊的椭圆,本文中研究一般椭圆的光斑中心定位[8-9]。由图像预处理得到的所有的边缘数据对圆进行最小二乘拟合,利用求解出的方程系数可直接计算出光斑的中心位置。

    • 在2维坐标系中,平面椭圆方程的隐式方程表示为F(a, x)=0,其中$a = (A, B, C, D, E, F) $为椭圆方程的参量,即椭圆隐式方程为:

      $ f\left( {x, y} \right) = A{x^2} + Bxy + C{y^2} + Dx + Ey + F $

      (1)

      式中,$ f\left( {x, y} \right)$为椭圆方程,其中xy分别为椭圆上点的横纵坐标值; $ A, B, C, D, E, F$为椭圆方程参量。

      比较常见的椭圆拟合方法有基于代数距离的最小二乘法。由最小二乘法可知,椭圆的拟合可以近似看成是求光斑边缘上的i个数据点$\left( {{x_i}, {y_i}} \right) $到曲线$ f\left( {x, y} \right) = 0$的残差平方和的最小值,即求目标函数$f\left( {B, C, D, E, F} \right) = \sum\limits_{i = 1}^m {} {\left( {{x_i}^2 + B{x_i}{y_i} + C{y_i}^2 + D{x_i} + E{y_i} + F} \right)^2} $的最小值。再由极值原理,欲使$f\left( {B, C, D, E, F} \right) $值为最小,必有$ \frac{{\partial f}}{{\partial B}} = \frac{{\partial f}}{{\partial C}} = \frac{{\partial f}}{{\partial D}} = \frac{{\partial f}}{{\partial E}} = \frac{{\partial f}}{{\partial F}} = 0$。由此可得一个线性方程组,将其转换成求线性矩阵方程Ax=b的解。这里A为5×5维实矩阵,x为待求5维列向量,bn维实向量。

      $ \begin{array}{l} \left[ {\begin{array}{*{20}{c}} {\sum\limits_{i = 1}^m {} {x_i}^2{y_i}^2}&{\sum\limits_{i = 1}^m {} {x_i}{y_i}^3}&{\sum\limits_{i = 1}^m {} {x_i}^2{y_i}}&{\sum\limits_{i = 1}^m {} {x_i}{y_i}^2}&{\sum\limits_{i = 1}^m {} {x_i}{y_i}}\\ {\sum\limits_{i = 1}^m {} {x_i}{y_i}^3}&{\sum\limits_{i = 1}^m {} {y_i}^4}&{\sum\limits_{i = 1}^m {} {x_i}{y_i}^2}&{\sum\limits_{i = 1}^m {} {y_i}^3}&{\sum\limits_{i = 1}^m {} {y_i}^2}\\ {\sum\limits_{i = 1}^m {} {x_i}^2{y_i}}&{\sum\limits_{i = 1}^m {} {x_i}{y_i}^2}&{\sum\limits_{i = 1}^m {} {x_i}^2}&{\sum\limits_{i = 1}^m {} {x_i}{y_i}}&{\sum\limits_{i = 1}^m {} {x_i}}\\ {\sum\limits_{i = 1}^m {} {x_i}{y_i}^2}&{\sum\limits_{i = 1}^m {} {y_i}^3}&{\sum\limits_{i = 1}^m {} {x_i}{y_i}}&{\sum\limits_{i = 1}^m {} {y_i}^2}&{\sum\limits_{i = 1}^m {} {y_i}}\\ {\sum\limits_{i = 1}^m {} {x_i}{y_i}}&{\sum\limits_{i = 1}^m {} {y_i}^2}&{\sum\limits_{i = 1}^m {} {x_i}}&{\sum\limits_{i = 1}^m {} {y_i}}&{\sum\limits_{i = 1}^m {} 1} \end{array}} \right]\left[ \begin{array}{l} B\\ C\\ D\\ E\\ F \end{array} \right] = \\ \left[ \begin{array}{l} - \sum\limits_{i = 1}^m {} {x_i}^3{y_i}\\ - \sum\limits_{i = 1}^m {} {x_i}^2{y_i}^2\\ - \sum\limits_{i = 1}^m {} {x_i}^3\\ - \sum\limits_{i = 1}^m {} {x_i}^2{y_i}\\ - \sum\limits_{i = 1}^m {} {x_i}^2 \end{array} \right] \end{array} $

      (2)

      式中,$ \left( {{x_i}, {y_i}} \right)$为光斑边界上点的坐标,$\left( {{x_i}, {y_i}} \right) \subset E $表示光斑图像中边缘上所有的点的集合。通过求解该线性方程组,便可得到椭圆方程参量B, C, D, EF,则中心坐标$\left( {{x_{\rm{c}}}, {y_{\rm{c}}}} \right) $为:

      $ \left\{ \begin{array}{l} {x_{\rm{c}}} = \frac{{BE - 2CD}}{{4C - {B^2}}}\\ {y_{\rm{c}}} = \frac{{BD - 2E}}{{4C - {B^2}}} \end{array} \right. $

      (3)
    • 在数字信号处理(digital signal processing,DSP)系统设计开发工具中以高斯消元法为基础对线性方程组进行求解,鉴于DSP builder的除法模块只能计算出被除数和除数相除以后的商和余数,故而采用移位寄存器将被除数移位倍增,再进行相除,经验证,其精度远远小于单精度浮点数,计算的结果与理论比较误差大于10%,因此选用Quartus Ⅱ的单精度浮点运算知识产权(intellectual property,IP)核,运用IP核设计运算处理单元。

      一般地,线性最小二乘问题的求解主要有迭代法和矩阵分解法两类。其中矩阵分解法可以最大限度实现并行计算,从而更适用于FPGA实现[10]。矩阵分解的方法主要有QR分解[11-12]LU分解[13]、奇异值分解[14]以及Cholesky分解[15]QR分解的计算过程复杂;而奇异值分解与QR分解的计算难度相当;大规模矩阵一般采用LU分解来进行分块分解,但是LU分解过程中需要大量的乘法和除法运算,会造成一定的延时,并且消耗过多的硬件资源;LDLT分解适用于对称的正定矩阵,分解过程中避免了开方运算,其计算的复杂度仅为LU分解的一半,故而本文中采取LDLT分解算法来实现矩阵的分解。

    • LDLT分解适用于对称的正定矩阵,对任意一个对称正定矩阵A,将其分解成A=LDLT的形式,其中L为对角线为1的下三角矩阵,D为对角阵,LTL的转置,即:

      $ \begin{array}{l} \mathit{\boldsymbol{A}} = \mathit{\boldsymbol{LD}}{\mathit{\boldsymbol{L}}^{\rm{T}}} = \left[ {\begin{array}{*{20}{c}} 1&{}&{}&{}&{}\\ {{l_{21}}}&1&{}&{}&{}\\ {{l_{31}}}&{{l_{32}}}&1&{}&{}\\ \ldots&\ldots&\ldots&\ldots &{}\\ {{l_{n1}}}&{{l_{n2}}}&{{l_{n3}}}& \ldots &1 \end{array}} \right] \cdot \\ \left[ {\begin{array}{*{20}{c}} {{d_1}}&{}&{}&{}&{}\\ {}&{{d_2}}&{}&{}&{}\\ {}&{}&{{d_3}}&{}&{}\\ {}&{}&{}& \ldots &{}\\ {}&{}&{}&{}&{{d_n}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} 1&{{l_{21}}}&{{l_{31}}}& \ldots &{{l_{n1}}}\\ {}&1&{{l_{32}}}& \ldots &{{l_{n2}}}\\ {}&{}&1& \ldots&\ldots \\ {}&{}&{}&1& \ldots \\ {}&{}&{}&{}&1 \end{array}} \right] \end{array} $

      (4)

      矩阵DL中的元素djlij根据下式计算:

      $ \left\{ \begin{array}{l} {d_j} = {a_{jj}} - \sum\limits_{k = 1}^{j - 1} {} {l_{jk}}^2{d_k}\\ {l_{ij}} = \frac{{{a_{ij}} - \sum\limits_{k = 1}^{j - 1} {} {l_{ik}}{d_k}{l_{jk}}}}{{{d_j}}} \end{array} \right. $

      (5)

      式中,$ j = 1, 2, \ldots , n;i = j + 1, j + 2, \ldots , n$。aijajj分别为矩阵A中的第ij列和第jj列所对应的元素,djdk为对角阵D中对角线上所对应的的元素(k=j-1),lijlikljk分别为矩阵L中第ij列、第ik列和第jk列所对应的元素。

      分解后即将求解的线性方程组Ax=b转换为求LDLTx=b的解,等同于求解3个三角线性方程组。

      $ \left\{ \begin{array}{l} \mathit{\boldsymbol{Lz}} = \mathit{\boldsymbol{b}}\\ \mathit{\boldsymbol{Dr}} = \mathit{\boldsymbol{z}}\\ {\mathit{\boldsymbol{L}}^{\rm{T}}}\mathit{\boldsymbol{x}} = \mathit{\boldsymbol{r}} \end{array} \right. $

      (6)

      式中, zr为引入的n维待求列向量。根据以上步骤即可求得方程的解向量x

    • 为了实现在FPGA硬件平台上高效的矩阵分解,必须合理采用流水线技术,设计并行计算[16]。算法并行化的根本是数据之间的依赖关系,下面将分析LDLT分解算法的数据之间的关系。

      根据(5)式,每一个内层循环用于计算矩阵L的一列元素和矩阵D对应的一个元素,以第j个内层循环中计算的具体数据依赖关系举例说明。将(5)式展开如下所示:

      $ \left\{ \begin{array}{l} {d_j} = {a_{jj}} - \left[ {{l_{i1}}^2{d_1} + {l_{i2}}^2{d_2} + \ldots + {l_{i(j - 1)}}^2{d_{(j - 1)}}} \right]\\ {l_{ij}} = [{a_{ij}} - ({l_{i1}}{d_1}{l_{j1}} + {l_{i2}}{d_2}{l_{j2}} + \ldots + \\ {l_{i(j - 1)}}{d_{(j - 1)}}{l_{j(j - 1)}})]/{d_j} \end{array} \right. $

      (7)

      由(7)式可知,矩阵D中元素dj的计算依赖于矩阵A中的元素ajj、矩阵L的第j行元素以及矩阵D的前(j-1)个元素。L矩阵的元素lij的计算依赖于矩阵A中的元素aij、矩阵D中的元素dj以及矩阵L的第i行元素。而lijl(i+1)j之间的计算没有依赖关系,且djlij的计算都要计算中间结果ljkdk

      通过分析数据依赖关系可知,LDLT分解算法的外层数据循环计算需要依存已经获得的数据,内层数据循环计算则没有依赖关系,直接并行计算。

    • 通过对LDLT分解算法及相关数据分析,可设计如图 1所示的分解结构框图。分解模块主要分为计算对角阵D元素的单元PED、下三角阵L元素的单元PEL和数据选择器,PED与PEL共享PED的计算结果ljkdk,所以可以将此结果缓存于先入先出缓存器(first in, first out, FIFO)中的直接用于lij的计算,从而避免重复运算,节省资源。PED的输入ljk通过四选一数据选择器,按照一定的顺序将4组存储矩阵L中元素的随机存取存储器(random access memory, RAM)中的数据读出,用于计算dj。矩阵A的数据通过数据分配器从外部传输到PEL与PED。其中,矩阵中aijdj, lij各自关联,但是PED与PEL计算存在延迟,故将数据存于FIFO_a0~FIFO_a4中。当前lij的计算结果与当前dj值关联很大,且lij的计算还会用到计算dj时的中间值ljkdk,所以PED要在PEL之前执行。在理想的情况下,为了实现运算速度最快,PEL的数量越多,效率越高。输入矩阵A为5×5的矩阵,在每一个内层循环中,需要计算的矩阵D的对角线元素dj仅为1个,采用一个PED单元即可;而计算的矩阵L的元素lij的个数最多为4个,且互不依赖,故采用4个PEL并行计算。分解模块具体结构如图 3所示。其中RAMD为存储矩阵D中元素的随机存取存储器。

      Figure 3.  Structure diagram of decomposition module

      图 3中,PED和PEL模块均采取流水化设计方案,运算单元主要包括乘法器(multiplication)、除法器(division)、加法器(addition)和减法器(subtraction),均使用Quartus Ⅱ提供的IEEE 754规范的浮点数运算IP核实现。

      在设计PED和PEL模块时,鉴于乘法器在计算时存在固有延时,所以将输入的数据ljk写入FIFO1中实现缓存,此时FIFO1就作为第2个乘法器的输入端,与第1个乘法器的计算结果进行第二次乘法运算。由于ljkdk为PED和PEL计算的共有部分,为减少乘法器的使用,故将计算中间结果ljkdk存于FIFO_ld中。累加器是由加法器和锁存器组成,每隔加法器计算的延时时间就将FIFO2中数据读出一个,累加结束后,将结果传递给减法器输入端,并使能减法器。最终的计算结果存入RAM中。内层循环结束后,循环计算5次,即可获得对角阵D和下三角阵L矩阵的元素djlij

    • 分解结束后,获得下三角矩阵L和对角阵D,根据(LDLT)x=b,可以等价成3个简单的方程组Lz=bDr=zLTx=r,由此可得求解公式:

      $ \left\{ \begin{array}{l} {z_m} = {b_m} - \sum\limits_{p = 1}^{m - 1} {} {l_{mp}}{z_p}, \left( {m = 1, 2, \ldots , q} \right)\\ {r_m} = \frac{{{z_m}}}{{{d_m}}}, \left( {m = 1, 2, \ldots , q} \right)\\ {x_m} = {r_m} - \sum\limits_{u = m + 1}^q {} {l_{um}}{x_u}, \left( {m = q, q - 1, \ldots , 1} \right) \end{array} \right. $

      (8)

      式中,lmplum分别表示矩阵L中第mp列和第um列所对应的元素;dm表示矩阵D对角线上对应的元素;zm, rm, xm, bm分别代表列向量z, r, x, b的第m个元素;zpxu分别代表列向量z的第p个元素和列向量x的第u个元素;q为待求向量的维数(q=5)。

      求解模块要用到PEL模块中计算的数据lij,故每算出一个lij的值,求解模块中就要有一个计算单元与之相对应。求解模块的结构框图与PED和PEL的框图基本相同。

    • 本文中基于Quartus Ⅱ 14.0的开发环境,将上述设计在Altera公司的DE1-SOC开发板上进行了验证。鉴于作者无法准确地得到光斑图像真实的中心,因此通过MATLAB程序生成一幅人工光斑图像,并附着一些噪声及干扰,该光斑图像的理论中心坐标为(150,166),分别在FPGA平台和计算机Visual Studio(VS)环境下利用本文中的算法对光斑进行中心定位,根据(3)式计算的最终结果对比如表 1所示。其中,FPGA的工作时钟设定为50MHz,并采取IEEE 754规范下的单精度浮点制IP核来实现设计中的各浮点运算。

      Table 1.  Deviation between FPGA and PC

      coordinate of center absolute deviation
      visual studio (150.06907, 166.097) (0.06907, 0.097)
      FPGA (150.06881, 166.048) (0.06881, 0.048)

      表 1可知,本方法对圆的中心定位的精度小于0.1pixel,但FPGA平台与计算机VS环境下实现相同的算法计算的结果有差异,出现此结果的原因是,在FPGA平台使用单精度浮点IP核进行计算时,首先要将数据转换为单精度浮点数(32bit),但对于较大的数据,尾数部分会有部分省略,进而对结果产生一定的误差。不过两个平台计算的中心位置的误差都小于0.1pixel,所以该方法保证了中心定位算法的精度。

    • 保证中心定位算法精度的前提下,为测试方法的效率,将实现LDLT解最小二乘问题的方法在FPGA平台和计算机VS下进行了计算效率比较。具体结果如表 2所示。

      Table 2.  Processing time

      FPGA visual studio proportion
      decomposition/μs 36.38 620 1/17.04
      solving/μs 37.06 995 1/26.85
      total/μs 73.44 1615 1/21.99

      表 2可知,在FPGA平台下实现解5维最小二乘问题时,本文中方法在保证精度的前提下,可达到相对于计算机21倍以上的速度。

    • 光斑中心位置检测是光纤光束测量系统中的一个关键步骤。针对现有的基于计算机平台实现光斑参量测量实时性差,数据吞吐量不高的问题,本文中提出一种基于FPGA光纤光斑中心定位的方法:先对光斑图像进行中值滤波、二值化和边缘提取等预处理,提高算法的抗干扰能力,而后在FPGA平台上将获得的光斑边界数据进行分析处理。实验结果表明,在单精度浮点制条件下,本文中方法对光斑中心定位的绝对误差小于0.1pixel,在保证精度的同时可以实现21倍以上的速度提升。

参考文献 (16)

目录

    /

    返回文章
    返回