留言板

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

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

连续波太赫兹成像系统的单幅图像超分辨重建

王欢 郎利影 庞亚军 张雷 郑伟 席思星

王欢, 郎利影, 庞亚军, 张雷, 郑伟, 席思星. 连续波太赫兹成像系统的单幅图像超分辨重建[J]. 红外与激光工程, 2023, 52(1): 20220292. doi: 10.3788/IRLA20220292
引用本文: 王欢, 郎利影, 庞亚军, 张雷, 郑伟, 席思星. 连续波太赫兹成像系统的单幅图像超分辨重建[J]. 红外与激光工程, 2023, 52(1): 20220292. doi: 10.3788/IRLA20220292
Wang Huan, Lang Liying, Pang Yajun, Zhang Lei, Zheng Wei, Xi Sixing. Single-image super-resolution reconstruction for continuous-wave terahertz imaging systems[J]. Infrared and Laser Engineering, 2023, 52(1): 20220292. doi: 10.3788/IRLA20220292
Citation: Wang Huan, Lang Liying, Pang Yajun, Zhang Lei, Zheng Wei, Xi Sixing. Single-image super-resolution reconstruction for continuous-wave terahertz imaging systems[J]. Infrared and Laser Engineering, 2023, 52(1): 20220292. doi: 10.3788/IRLA20220292

连续波太赫兹成像系统的单幅图像超分辨重建

doi: 10.3788/IRLA20220292
基金项目: 国家自然科学基金(11904073);河北省重点研发项目(20371802D);河北省自然科学基金( F2019402351);河北省教育厅青年拔尖人才项目(BJ2020028)
详细信息
    作者简介:

    王欢,女,硕士生,主要研究方向为太赫兹图像处理

    郎利影,女,教授,博士,主要研究方向为人工智能与模式识别及光信息应用

  • 中图分类号: TN219; TP399

Single-image super-resolution reconstruction for continuous-wave terahertz imaging systems

  • 摘要: 针对现有的太赫兹成像系统所需硬件设备复杂且昂贵的问题,设计了基于单幅图像超分辨重建的连续波太赫兹成像系统,降低设备复杂度和硬件成本。通过对该成像系统生成的太赫兹图像进行双维度预处理,降低图像处理的占用内存,提高后续处理速度。引入限制对比度自适应直方图均衡方法对太赫兹图像进行分区域对比度提升,有效解决太赫兹图像对比度低的问题。利用稀疏表示和字典学习实现太赫兹图像的超分辨重建,提出了反余割拟牛顿平滑零范数的算法解决零范数优化问题,提高了重建精度。通过对该成像系统采集的单幅太赫兹图像进行超分辨重建,在边缘强度上提高了3.232,在平均梯度对比中提高了0.300,验证了对单幅太赫兹图像超分辨重建的有效性与优越性。
  • 图  1  成像系统光路图

    Figure  1.  Optical path schematic of imaging system

    图  2  成像系统实物图

    Figure  2.  Physical diagram of imaging system

    图  3  太赫兹扫描图像

    Figure  3.  Scanned raw terahertz image

    图  4  三种函数的函数值分布($ \sigma = 0.1 $

    Figure  4.  Distribution of the three functions($ \sigma = 0.1 $

    图  5  样品实物图。(a)样品大小示意图;(b)样品位置示意图

    Figure  5.  Physical image of the sample. (a) Sample size diagram; (b) Sample location diagram

    图  6  重建结果。(a) lanczos2插值算法;(b) NEDI算法;(c) NSL0算法;(d) ICNSL0算法;(e) lanczos2细节图;(f) NEDI细节图;(g) NSL0细节图;(h) ICNSL0细节图

    Figure  6.  Reconstructions images. (a) lanczos2 interpolation algorithm; (b) NEDI algorithm; (c) NSL0 algorithm; (d) ICNSL0 algorithm; (e) Detail of lanczos2 interpolation algorithm; (f) Detail of NEDI algorithm; (g) Detail of NSL0 algorithm; (h) Detail of ICNSL0 algorithm

    表  1  ICNSL0算法的步骤

    Table  1.   Steps of ICNSL0 algorithm

    输入:图像Y,过完备字典D,递减数列$ \sigma [{\sigma _1},{\sigma _2},...,{\sigma _n}] $,步长d$ {r_0} = 0 $
    输出:稀疏表示系数$ \alpha $
    1. 初始化计算公式$\alpha = \alpha {}^{\rm{T} }{(\alpha \alpha {}^{\rm{T} })^{ - 1} }{{Y} }$
    2. 最速下降法取函数一阶导数的相反数,设为阻尼牛顿法初始值
    3. 阻尼牛顿法更新优化方向:$ \alpha = \alpha + d $
    4. 梯度投影更新$ \alpha $值:$\alpha = \alpha - \alpha {}^{\rm{T} }{(\alpha \alpha {}^{\rm{T} })^{ - 1} }(D\alpha - {{Y} })$
    5. 计算$r = {{Y} } - D\alpha$, 若$ r - {r_0} < \varepsilon $或超出最大迭代次数,输出结果$ \alpha $
    否则跳至步骤3
    下载: 导出CSV

    表  2  不同算法重建结果的NIQE、PSNR、SSIM值、边缘强度和平均梯度

    Table  2.   NIQE, PSNR, SSIM, edge intensity and average gradient value of the reconstruction results of different algorithms

    Measureslanczos2NEDINSL0ICNSL0
    NIQE8.6486.3527.4996.481
    PSNR27.71314.43922.872 421.433
    SSIM0.897 10.2660.3040.864 2
    ES16.70220.74119.59622.828
    AG1.5251.9611.7962.096
    下载: 导出CSV
  • [1] Ouchi T, Kajiki K, Koizumi T, et al. Terahertz imaging system for medical applications and related high efficiency Terahertz devices [J]. J Infrared Millim Terahertz Waves, 2014, 35(1): 118-130. doi:  10.1007/s10762-013-0004-5
    [2] Chen Z, Ma X Y, Zhang B, et al. A survey on Terahertz communications [J]. China Commun, 2019, 16(2): 1-35.
    [3] Ge H Y, Lv M, Lu X J, et al. Applications of THz spectral imaging in the detection of agricultural products [J]. Photonics, 2021, 8(11): 518. doi:  10.3390/photonics8110518
    [4] Liu Xinyuan, Zeng Haomin, Tian Xin, et al. Transmission simulation and safety analysis of terahertz radiation in skin [J]. Optics and Precision Engineering, 2021, 29(5): 999-1007. (in Chinese) doi:  10.37188/OPE.20212905.0999
    [5] Wang Yuye, Jiang Bozhou, Xu Degang, et al. Continuous terahertz wave biological tissue imaging technology based on focal plane array [J]. Acta Optica Sinica, 2021, 41(7): 0711001. (in Chinese) doi:  10.3788/AOS202141.0711001
    [6] Li Yashang, Zhao Guozhong, Wei Qingyun, et al. Study on the performance of terahertz passive imaging system [J]. Infrared and Laser Engineering, 2020, 49(4): 0404005. (in Chinese) doi:  10.3788/IRLA202049.0404005
    [7] Li Y D, Hu W D, Zhang X, et al. Adaptive terahertz image super-resolution with adjustable convolutional neural network [J]. Opt Express, 2020, 28(15): 22200-22217. doi:  10.1364/OE.394943
    [8] Zhang Saiwen, Deng Yaqi, Wang Chong, et al. Research on super-resolution fluorescence microscopy imaging based on multiple measurement vector compressed sensing [J]. Infrared and Laser Engineering, 2021, 50(11): 20210484. (in Chinese) doi:  10.3788/IRLA20210484
    [9] Zhao Haoguang, Qu Hanshi, Wang Xin, et al. Super-resolution reconstruction of micro-scanning images [J]. Optics and Precision Engineering, 2021, 29(10): 2456-2464. (in Chinese) doi:  10.37188/OPE.20212910.2456
    [10] Liu Mingxin, Zhang Xin, Wang Lingjie, et al. Optimization of matching coded aperture with detector based on compressed sensing spectral imaging technology [J]. Chinese Optics, 2020, 13(2): 290-301. (in Chinese) doi:  10.3788/CO.20201302.0290
    [11] Li W, Li B, Li P F. Image super-resolution via sparse representation and local texture constrain[C]//12th IEEE Conference on Industrial Electronics and Applications (ICIEA), 2017: 1044-1049.
    [12] Wu Y Q, Pan S H, Bi S J, et al. Kurdyka-Lojasiewicz property of zero-norm composite functions [J]. J Optim Theory Appl, 2020, 188(1): 94-112.
    [13] Xu M X, Yang Y, Sun Q S, et al. Image super-resolution reconstruction based on adaptive sparse representation [J]. Concurr Comput -Pract Exp, 2018, 30(24): e4968. doi:  10.1002/cpe.4968
    [14] Wang M J, Feng Z D, Wu J G. Pixel super-resolution lensless in-line holographic microscope with hologram segmentation [J]. Chin Opt Lett, 2019, 17(11): 110901.
    [15] Aamir M, Rehmanp Z, Pu Y F, et al. Image enhancement in varying light conditions based on wavelet transform[C]//16th IEEE International Computer Conference on Wavelet Active Media Technology and Information Processing (ICCWAMTIP), 2019: 317-322.
    [16] Chen Tianze, Ge Baozhen, Luo Qijun. Pose estimation for free binocular cameras based on reprojection error optimization [J]. Chinese Optics, 2021, 14(6): 1400-1409. (in Chinese) doi:  10.37188/CO.2021-0105
    [17] Wu L Y, Zhang X G, Deng J F. Making a “Completely blind” image quality analyzer [J]. IEEE Signal Process Lett, 2013, 20(3): 209-212. doi:  10.1109/LSP.2012.2227726
    [18] Li Q H, Fang Y M, Lin W S, et al. Gradient-weighted structural similarity for image quality assessments[C]//IEEE International Symposium on Circuits and Systems (ISCAS), Lisbon, IEEE, 2015: 2165-2168.
    [19] Shi Z F, Zhang J P, Cao Q J, et al. Full-reference image quality assessment based on image segmentation with edge feature [J]. Signal Process, 2018, 145: 99-105. doi:  10.1016/j.sigpro.2017.11.015
  • [1] 苏英蔚, 田震.  基于相位梯度光栅介电超表面的高效太赫兹波异常反射器 . 红外与激光工程, 2023, 52(2): 20220304-1-20220304-7. doi: 10.3788/IRLA20220304
    [2] 武军安, 郭锐, 刘荣忠, 柯尊贵, 赵旭.  边缘区域约束的导向滤波深度像超分辨率重建算法 . 红外与激光工程, 2021, 50(1): 20200081-1-20200081-11. doi: 10.3788/IRLA20200081
    [3] 佘荣斌, 祝永乐, 刘文权, 鲁远甫, 李光元.  太赫兹单像素计算成像原理及其应用(特邀) . 红外与激光工程, 2021, 50(12): 20210717-1-20210717-19. doi: 10.3788/IRLA20210717
    [4] 贺敬文, 董涛, 张岩.  太赫兹波前调制超表面器件研究进展 . 红外与激光工程, 2020, 49(9): 20201033-1-20201033-11. doi: 10.3788/IRLA20201033
    [5] 刘朝阳, 刘力源, 吴南健.  采用CMOS太赫兹波探测器的成像系统 . 红外与激光工程, 2017, 46(1): 125001-0125001(6). doi: 10.3788/IRLA201746.0125001
    [6] 吴俊政, 严卫东, 倪维平, 张晗.  圆周阵列太赫兹干涉成像仿真 . 红外与激光工程, 2017, 46(8): 825002-0825002(8). doi: 10.3788/IRLA201746.0825002
    [7] 张智, 林栩凌, 张建兵, 何红艳.  太赫兹成像质量提升方法 . 红外与激光工程, 2017, 46(11): 1126002-1126002(5). doi: 10.3788/IRLA201746.1126002
    [8] 李晗, 余晨.  太赫兹波对肾癌组织的光谱检测 . 红外与激光工程, 2016, 45(5): 525001-0525001(5). doi: 10.3788/IRLA201645.0525001
    [9] 刘伟伟, 赵佳宇, 张逸竹, 王志, 储蔚, 曾斌, 程亚.  飞秒激光成丝过程中的太赫兹波超光速传输现象研究 . 红外与激光工程, 2016, 45(4): 402001-0402001(7). doi: 10.3788/IRLA201645.0402001
    [10] 张学迁, 张慧芳, 田震, 谷建强, 欧阳春梅, 路鑫超, 韩家广, 张伟力.  利用介质超材料控制太赫兹波的振幅和相位 . 红外与激光工程, 2016, 45(4): 425004-0425004(6). doi: 10.3788/IRLA201645.0425004
    [11] 王花, 孙晓红, 王真, 齐永乐, 王毅乐.  太赫兹波超材料吸波体的特性分析 . 红外与激光工程, 2016, 45(12): 1225003-1225003(5). doi: 10.3788/IRLA201645.1225003
    [12] 王启超, 汪家春, 赵大鹏, 林志丹, 苗雷.  太赫兹波对烟幕的透射能力研究 . 红外与激光工程, 2015, 44(12): 3696-3700.
    [13] 许文忠, 钟凯, 梅嘉林, 徐德刚, 王与烨, 姚建铨.  太赫兹波在沙尘中衰减特性 . 红外与激光工程, 2015, 44(2): 523-527.
    [14] 杨晶, 赵佳宇, 郭兰军, 刘伟伟.  超快激光成丝产生太赫兹波的研究 . 红外与激光工程, 2015, 44(3): 996-1007.
    [15] 李琦, 杨永发, 胡佳琦.  一种用于太赫兹共焦扫描图像复原的复合算法 . 红外与激光工程, 2015, 44(1): 321-326.
    [16] 江舸, 成彬彬, 杨陈, 蔡英武, 张健.  提高0.14 THz 雷达成像质量的新方法 . 红外与激光工程, 2014, 43(9): 2912-2918.
    [17] 郭澜涛, 牧凯军, 邓朝, 张振伟, 张存林.  太赫兹波谱与成像技术 . 红外与激光工程, 2013, 42(1): 51-56.
    [18] 李泽, 王民钢, 刘小华, 赵跃进, 张存林.  基于压缩传感的太赫兹成像 . 红外与激光工程, 2013, 42(6): 1523-1527.
    [19] 李乾坤, 李德华, 周薇, 马建军, 鞠智鹏, 屈操.  单缝双环结构超材料太赫兹波调制器 . 红外与激光工程, 2013, 42(6): 1553-1556.
    [20] 罗俊, 公金辉, 张新宇, 季安, 谢长生, 张天序.  基于超材料的连续太赫兹波透射特性研究 . 红外与激光工程, 2013, 42(7): 1743-1747.
  • 加载中
图(6) / 表(2)
计量
  • 文章访问数:  284
  • HTML全文浏览量:  105
  • PDF下载量:  100
  • 被引次数: 0
出版历程
  • 收稿日期:  2022-04-10
  • 修回日期:  2022-05-20
  • 刊出日期:  2023-01-18

连续波太赫兹成像系统的单幅图像超分辨重建

doi: 10.3788/IRLA20220292
    作者简介:

    王欢,女,硕士生,主要研究方向为太赫兹图像处理

    郎利影,女,教授,博士,主要研究方向为人工智能与模式识别及光信息应用

基金项目:  国家自然科学基金(11904073);河北省重点研发项目(20371802D);河北省自然科学基金( F2019402351);河北省教育厅青年拔尖人才项目(BJ2020028)
  • 中图分类号: TN219; TP399

摘要: 针对现有的太赫兹成像系统所需硬件设备复杂且昂贵的问题,设计了基于单幅图像超分辨重建的连续波太赫兹成像系统,降低设备复杂度和硬件成本。通过对该成像系统生成的太赫兹图像进行双维度预处理,降低图像处理的占用内存,提高后续处理速度。引入限制对比度自适应直方图均衡方法对太赫兹图像进行分区域对比度提升,有效解决太赫兹图像对比度低的问题。利用稀疏表示和字典学习实现太赫兹图像的超分辨重建,提出了反余割拟牛顿平滑零范数的算法解决零范数优化问题,提高了重建精度。通过对该成像系统采集的单幅太赫兹图像进行超分辨重建,在边缘强度上提高了3.232,在平均梯度对比中提高了0.300,验证了对单幅太赫兹图像超分辨重建的有效性与优越性。

English Abstract

    • 太赫兹指的是频率在0.1~10 THz(波长30~3000 μm)之间的电磁波辐射,位于电磁波谱中微波和红外波之间。它的位置使其存在宽带性、光子能量低、吸水性和瞬时性等优良特性。太赫兹成像的相关技术在医学成像研究[1]、安全检查[2]和环境监测[3]和生物研究[4]等方面都存在着广泛应用。受光学衍射极限和成像环境、成像核心元件发展现状的影响,太赫兹图像存在亮度、对比度和分辨率低,纹理特征不明显的问题[5]。为解决上述问题,对成像系统的核心硬件直接升级是最直接有效的方式,但该方法技术难度高且实际应用中会面临成本约束的问题[6]。近年来,研究人员将注意力集中在开发新的算法来提高太赫兹图像的质量[7]

      目前,提高图像分辨率的常用方法是超分辨重建算法[8],利用同一场景的低分辨率(Low resolution,LR)图像得到对应的高分辨率(High resolution,HR)图像,提高现有图像空间分辨率[9]。传统单幅图像的超分辨重建算法中,基于插值的重建算法会忽略图像的高频边缘信息。基于模型的算法需要同一场景的多帧图像进行重建,而太赫兹图像的原始数据无法满足这一条件,不适用于太赫兹图像的超分辨重建。基于学习的重建算法中,利用稀疏表示和字典学习的原理实现超分辨重建,对原始数据的要求较低,同时在重建中不会忽略高频信息。这类方法[10]通过求解稀疏表示系数结合过完备字典实现超分辨重建。针对这一过程,Li[11]提出了一种基于稀疏编码改进的超分辨重建算法。为了解决稀疏表示系数求解过程中的耗时问题,Wu[12]将L1参数约束变为L2约束,以模型精度为代价简化求解过程同时缩短时间。Xu[13]添加了重建的先验性约束,利用自适应结构相似性提高重建精度。

      太赫兹图像质量受成像系统和光学衍射极限限制,图像分辨率受限,存在低对比度、边缘信息少,信息量不足的问题,无法满足后续研究需要。因此超分辨重建算法被提出并用于突破硬件限制提高图像分辨率。实际成像系统获取的太赫兹图像往往是单幅图像,其数据量不足以进行超分辨重建。因此双维度预处理构造高分辨率图像方法被提出,可以突破原始数据量的限制。由于太赫兹图像所含信息少,使用联合字典学习方法,丰富原始数据的信息量,用于图像的超分辨重建。此外,针对稀疏表示和字典学习算法的零范数优化问题,提出新的拟合函数提高图像重建精度,优化重建后图像的纹理特征,丰富太赫兹图像的边缘信息。

    • 文中搭建的连续波太赫兹快速扫描成像系统以美国Tera Sense公司的TeraFast-256-HS硬件设备为基础,该设备由线阵太赫兹相机与连续波太赫兹源组成。太赫兹源向空气中发射出电磁波后,不同的环境中太赫兹波传播时衰减率不同,不同样品对太赫兹波的透射率或反射率也不同,使用探测器获取携带了样品信息的太赫兹波,利用上位机对该电磁波进行信息处理从而获取该样品的太赫兹图像。

      文中搭建的太赫兹成像系统中使用的线阵探测器大小为256 pixel×1 pixel,像素尺寸为1.5 mm × 3 mm,成像速率为64 fps。由于线阵探测器接收像素大小有限,因此需要通过样品的移动实现整体扫描。传统方法中利用传送带带动样品沿某一方向位移实现太赫兹波的扫描[14]。为了降低系统的复杂性从而降低设备成本,仅使用单个的太赫兹源、线阵探测器结合可调控的旋转台实现连续波太赫兹线性扫描成像,通过预先对不同厚度的标准样品进行测试模拟,利用不同厚度样品生成的强度值图像拟合求得衰减率和透射率。综上所述,成像系统的光路图如图1所示。

      成像系统中的太赫兹源为0.1 THz的频率辐射,功率为20 μW,成像动态范围为24 dB。太赫兹源和线性探测器分别放置于待测样品两侧,并保持在同一水平高度。待测样品放置在旋转台上,旋转台通过驱动控制匀速旋转。成像系统的实物图如图2所示。

      图  1  成像系统光路图

      Figure 1.  Optical path schematic of imaging system

      图  2  成像系统实物图

      Figure 2.  Physical diagram of imaging system

      太赫兹源发射出的电磁波从样品左侧(图2中标注的样品B面)开始进行扫描,可调控的旋转台以6 (°)/s速度匀速旋转一周,与线性探测器间距离为5 cm。线性探测器获取到太赫兹图像的扫描视频,将其连接到太赫兹TERA Fast®图像处理软件。图像视频经处理拼接得到原始太赫兹扫描图像。扫描图像为从侧视视角观察的样品旋转一周后得到的结果图像。

    • 待测样品如图2被放置在可调控旋转台上,太赫兹波对样品持续扫描,固定发射太赫兹波的角度,从三角架B面进行扫描,样品台匀速旋转一周,得到完整的太赫兹侧面扫描图像如图3所示。该图像的行像素值大小由线阵探测器的单次探测能力决定,列像素值大小等于成像速率与成像时间乘积。

      图  3  太赫兹扫描图像

      Figure 3.  Scanned raw terahertz image

      采用双维度行列分离的预处理方法,将原始太赫兹图像(256 pixel× 256 pixel)中每一行视为独立的图像,对其进行行维度的超分辨重建(放大倍数为2)后得到重建结果图像(512 pixel× 256 pixel)。同时, 对原始图像进行列维度重建得到重建图像(256 pixel× 512 pixel)。但是,如果不进行该预处理,在原始图像为256 pixel×256 pixel的情况下直接进行超分辨率重建,处理过程中所占用内存将超过30 G,而内存占用率过大时会导致重建速度缓慢,从而影响重建的效率。经过上述分离操作优化后,内存占用能够减少至650~700 M。因此,双维度分离的预处理方法对太赫兹扫描图像的超分辨重建非常重要。

      太赫兹图像存在对比度低和边缘模糊的缺点。对太赫兹图像进行超分辨重建前,对太赫兹图像进行对比度有限的自适应直方图均衡化[15],分区域进行对比度提升,抑制图像中无用信息,增强样品与背景的对比度。而成像系统利用透射样品的光谱强度成像,样品的内部损耗和不同材料边缘使得太赫兹波受散射效应影响了强度值的分布,因此图像为该样品的伪彩色图像,而处理彩色图像时会增加通道信息,造成大量的时间浪费。在对图像进行处理时将其转换成灰度值图像,降低处理过程的内存占比,同时灰度图像中的明暗分布能够确定并识别出样品内部的损失、不规则缺陷和边缘信息。

    • 基于稀疏表示和字典学习的重建算法是基于过完备字典的约束,利用较少的信息对原始图像进行线性表示。该类算法包括字典学习和图像重建两个过程。字典学习利用学习到的图像特征对字典约束进行训练,使其满足稀疏先验约束,即HR图像与LR图像存在相同的稀疏表达系数。图像重建要求LR图像与HR图像生成模型一致,传统算法中将HR图像下采样和模糊化后获得LR图像,满足这一约束条件。

      实际成像系统获取的太赫兹图像本身分辨率低无法作为HR图像直接进行超分辨重建,因此无法直接满足上述两个约束条件。文中利用NEDI插值算法对LR图像插值构建类HR图像,使其满足生成模型一致的约束。为了防止LR图像与类HR图像像素矩阵相似度过高,导致字典学习过程中图像特征单一使得模型泛化能力差,将LR图像像素矩阵进行随机值浮动,进行差异化处理,由此使得原始太赫兹图像满足重建约束。为了满足稀疏先验约束,将LR图像与HR图像进行联合字典训练,用于提高图像稀疏表示系数的相似性。基于稀疏表示和字典学习的图像重建原理可以表示为:

      $$ {{Y}} = D\alpha = \sum\limits_{i = 1}^k {{d_i}{\alpha _i}} $$ (1)

      式中:$Y \in {R^n}$表示原始太赫兹图像;$D = \left[ {{d_1},{d_2},\cdots,{d_i}} \right] \in {R^{{{n}} \times {{m}}}}$表示过完备字典;$\alpha = {\left[ {{\alpha _1},{\alpha _2},\cdots,{\alpha _{\text{i}}}} \right]^{\rm{T}}} \in {R^m}$表示稀疏表示系数。每个向量$ {d_i} $代表过完备字典中第i个字典原子,字典中的大部分元素为零。由于稀疏性限制的存在,稀疏表示系数的大部分系数也为零,图像Y在过完备字典D的约束下可以利用稀疏表示系数$ \alpha $表示。

      在进行联合字典训练时,需要大量的原始数据用于获取丰富的图像特征,使其可以完整地表示出待重建图像的所有特征。同时要求数据有一定的相关性或重叠性。然而太赫兹图像所含信息少、原始数据量不足,因此将图像中的每一行和每一列作为独立数据参与过完备字典的建立,利用这些信息构建过完备字典,进行联合字典训练的表达式为:

      $$ \mathop {\min }\limits_{\alpha \in {R_n}} {\left\| {{{Y}} - D\alpha } \right\|_0} + \lambda \left( {\frac{1}{L} + \frac{1}{H}} \right){\left\| \alpha \right\|_0} $$ (2)

      式中:LH分别表示低分辨率图像块和高分辨率图像块的维度,该参数可以作为两个字典表示误差的权值;$ \lambda $为正则化参数,对训练字典的结果进行规范化。该式为联合字典训练模型,求解得到对应的字典D后,经过拆分得到LR图像字典和HR图像字典。获取了过完备字典后再求解对应的稀疏表示系数,其计算公式为:

      $$ \mathop {\min }\limits_{\alpha \in {{\text{R}}_{\text{n}}}} {\left\| \alpha \right\|_0},{\text{s}}.{\text{t}}.{{Y}} = D\alpha $$ (3)

      式中:$ {\left\| \alpha \right\|_0} $表示非零元素个数,用于衡量矩阵的稀疏程度。因此求解稀疏表示系数的问题转化为零范数优化,然而直接求解零范数属于NP-Hard问题。通常利用连续函数对其进行近似拟合求解最优值来解决该问题。因此文中提出反余割函数近似拟合零范数,函数表达式为:

      $$ {f_\sigma }({x_i}) = \frac{2}{\pi }{\rm{{acsc}}} \left( {\frac{{x_i^{^{\frac{1}{4}}}}}{{{{\mu }}{\sigma ^2}}}} \right) $$ (4)

      式中:$ \sigma > 0 $${\mu}$取大于1的常数,$ {x_i} $表示$ x $的第i个分量。该拟合函数的函数族表示为:

      $$ {\left\| x \right\|_0} = \mathop {\lim }\limits_{\sigma \to 0} {F_\sigma }\left( x \right) = \mathop {\lim }\limits_{\sigma \to 0} \sum\limits_1^N {{f_\sigma }({x_i})} $$ (5)

      由上式可以看出,当$ \sigma $趋近于0时,该函数族$ {F_\sigma }\left( x \right) $代表$ {x_i} $中不为0的元素个数,因此零范数优化的问题转化为:$ \sigma \to 0 $$ F\left( x \right) $的最小化问题。从数学意义的角度分析,文中提出的反余割函数相比SL0算法中使用的高斯函数和NSL0算法中使用的双曲正切函数更具“陡峭性”,三种函数在[−1,1]区间(其中$ \sigma = 0.1 $)的仿真对比结果如图4所示。因此反余割函数对零范数的拟合度更高,更逼近于零范数,对零范数的估计更精确,因而寻找到的最优解也更为精准。

      图  4  三种函数的函数值分布($ \sigma = 0.1 $

      Figure 4.  Distribution of the three functions($ \sigma = 0.1 $

      函数进行迭代优化时,最速下降法在初始值的选择上存在很大优势,但迭代过程会产生锯齿效应。而拟牛顿法收敛速度快,但其结果受初始值选择影响较大。因此选择利用最速下降法选择迭代初始值,利用拟牛顿法更新迭代步长。其中拟牛顿法的步长求解公式表示为:

      $$ d = - {\nabla ^2}{F_\sigma }{\left( x \right)^{ - 1}}\nabla {F_\sigma }\left( x \right) $$ (6)

      其中$ \nabla {F_\sigma }\left( x \right) $$ {\nabla ^2}{F_\sigma }\left( x \right) $的求解结果如公式(7)和(8)所示,可表示为:

      $$ \begin{gathered} \nabla {F_\sigma }\left( x \right) = {\left[ {\frac{{\partial {f_\sigma }({x_1})}}{{\partial {x_1}}},...,\frac{{\partial {f_\sigma }({x_N})}}{{\partial {x_N}}}} \right]^{\text{T}}}=\\ {\left[ {\dfrac{{\dfrac{1}{2}{\mu}{\sigma ^2}x_{^1}^{ - \frac{3}{4}}}}{{\pi \left( {{\mu}{}^2{\sigma ^4} + x_1^{\frac{1}{2}}} \right)}},...,\frac{{\dfrac{1}{2}{\mu}{\sigma ^2}x_N^{ - \frac{3}{4}}}}{{\pi \left( {{\mu}{}^2{\sigma ^4} + x_N^{\frac{1}{2}}} \right)}}} \right]^{\text{T}}} \\ \end{gathered} $$ (7)
      $$ {\nabla ^2}{F_\sigma }\left( x \right) = \left( {\begin{array}{*{20}{c}} {\dfrac{{{\partial ^2}{f_\sigma }({x_1})}}{{\partial x_1^2}}}&0&{...}&0 \\ 0&{\dfrac{{{\partial ^2}{f_\sigma }({x_2})}}{{\partial x_2^2}}}&{...}&0 \\ {...}&{...}&{...}&{...} \\ 0&0&{...}&{\dfrac{{{\partial ^2}{f_\sigma }({x_{{N}}})}}{{\partial x_{{N}}^2}}} \end{array}} \right) $$ (8)

      其中矩阵对角线的元素表达式为:

      $$ \frac{{{\partial ^2}{f_\sigma }({x_{{k}}})}}{{\partial x_{{k}}^2}} = \frac{{2{{\mu}^3}{\sigma ^6}x_{^{{k}}}^{ - 2} - \frac{5}{2}{\mu}{\sigma ^2}x_{^1}^{ - \frac{3}{2}}}}{{4\pi {{\left( {x_{{N}}^{\frac{1}{2}} - {\mu}{}^2{\sigma ^4}} \right)}^{\frac{3}{2}}}}} $$ (9)

      由于$ {\nabla ^2}{F_\sigma }\left( x \right) $(Hesse矩阵)只有满足非奇异且正定时,才能满足牛顿方向为下降方向。为了使得Hesse矩阵的特征根均为正数,即公式(9)的值恒大于零,需要构造一个新矩阵对其进行修正,修正结果表示为:

      $$ {G_k}{\text{ = }}\frac{{{{\mu}^3}{\sigma ^6}x_{^k}^{ - 2}}}{{2\pi {{\left( {x_N^{\frac{1}{2}} - {\mu}{}^2{\sigma ^4}} \right)}^{\frac{3}{2}}}}} $$ (10)

      结合公式(6)对修正牛顿方向的定义,计算得出修改后牛顿方向如公式(11)所示,迭代步长取该方向的最小值。

      $$ d = {\left[ {\frac{{x_1^{\frac{3}{2}} - {x_1}{\mu ^2}{\sigma ^4}}}{{{{\mu}^2}{\sigma ^4}}},...,\frac{{x_N^{\frac{3}{2}} - {x_N}{{\mu}^2}{\sigma ^4}}}{{{{\mu}^2}{\sigma ^4}}}} \right]^{\rm{T}}} $$ (11)

      综上所述,基于稀疏表示和字典学习的太赫兹图像超分辨重建ICNSL0算法的步骤可描述为如表1所示。

      表 1  ICNSL0算法的步骤

      Table 1.  Steps of ICNSL0 algorithm

      输入:图像Y,过完备字典D,递减数列$ \sigma [{\sigma _1},{\sigma _2},...,{\sigma _n}] $,步长d,$ {r_0} = 0 $
      输出:稀疏表示系数$ \alpha $
      1. 初始化计算公式$\alpha = \alpha {}^{\rm{T} }{(\alpha \alpha {}^{\rm{T} })^{ - 1} }{{Y} }$
      2. 最速下降法取函数一阶导数的相反数,设为阻尼牛顿法初始值
      3. 阻尼牛顿法更新优化方向:$ \alpha = \alpha + d $
      4. 梯度投影更新$ \alpha $值:$\alpha = \alpha - \alpha {}^{\rm{T} }{(\alpha \alpha {}^{\rm{T} })^{ - 1} }(D\alpha - {{Y} })$
      5. 计算$r = {{Y} } - D\alpha$, 若$ r - {r_0} < \varepsilon $或超出最大迭代次数,输出结果$ \alpha $;
      否则跳至步骤3

      其中递减序列$ \sigma [{\sigma _1},{\sigma _2},...,{\sigma _n}] $的取值设置为${\sigma _n}{\text{ = 0}}{\text{.8}} × {\sigma _{n{\text{ - 1}}}}$,初始参数$ {\sigma _1}{\text{ = 1}} $,最小值不低于0.0001。利用ICNSL0算法得到的稀疏表示系数结合HR图像的过完备字典进行超分辨重建。对单幅原始太赫兹图像进行行列双维度的超分辨重建得到HR1图像和HR2图像。将双维度重建图像进行混合插值进行融合重建:边缘部分使用NEDI插值算法,图像平滑区域使用双立方插值算法降低融合图像的时间成本。图像融合后得到大小为512 pixel×7680 pixel图像。最后对这一图像进行滤波反投影处理[16],得到大小为512 pixel×512 pixel的样品太赫兹图像(三角架俯视图)。

    • 为了说明文中算法的重建效果,使用连续波太赫兹快速线阵扫描成像系统对样品进行扫描成像,该样品的形状为边长4 cm的等腰三角形如图5(a)所示,其材质为铜质三角架,在旋转台上的相对位置如图5(b)所示。

      图  5  样品实物图。(a)样品大小示意图;(b)样品位置示意图

      Figure 5.  Physical image of the sample. (a) Sample size diagram; (b) Sample location diagram

      太赫兹源对该样品进行扫描后得到单幅低分辨率图像,对该图像使用文中提出的算法,并与传统的lanczos2插值算法、基于边缘指导的NEDI算法和基于稀疏表示的NSL0算法的重建结果进行比较。各个算法重建图像的对比结果如图6(a)~(d)所示,其大小为512 pixel×512 pixel。

      图  6  重建结果。(a) lanczos2插值算法;(b) NEDI算法;(c) NSL0算法;(d) ICNSL0算法;(e) lanczos2细节图;(f) NEDI细节图;(g) NSL0细节图;(h) ICNSL0细节图

      Figure 6.  Reconstructions images. (a) lanczos2 interpolation algorithm; (b) NEDI algorithm; (c) NSL0 algorithm; (d) ICNSL0 algorithm; (e) Detail of lanczos2 interpolation algorithm; (f) Detail of NEDI algorithm; (g) Detail of NSL0 algorithm; (h) Detail of ICNSL0 algorithm

      通过成像结果对比可以得出,lanczos2算法和NEDI算法重建图像均属于基于插值的重建算法,重建图像对比度较低,样品与成像背景边缘模糊、不易区分。基于稀疏表示的NSL0算法重建结果中高频信息(图像边缘部分)比较模糊。而ICNSL0算法重建图像的对比度有所提高,样品边缘部分的清晰度有所提升。截取部分重建图像(图中蓝色方框部分)进行对比如图6(e)~(h)所示,该图为同比例放大的三脚架左边缘部分。通过细节对比可以看出lanczos算法、NEDI算法和NSL0算法重建图像中三角架的边缘较为模糊,而ICNSL0算法重建图像对比度提高后使得重建后的三角架厚度(图中白色部分)的边缘细节更加清晰。

      由不同算法重建后的图像设立统一衡量标准,利用主观和客观双重评价方法对重建图像的质量比较。主观评价中使用自然图像质量评估(Natural Image Quality Evaluator,NIQE)[17]作为衡量标准。客观衡量标准分为两类:第一类以统计学为基础的峰值信噪比(Peak Signal-to-Noise Ratio, PSNR)和结构相似度(Structural Similarity, SSIM)[18]。利用双三次线性插值算法构建“参考图像”。第二类以边缘信息为主的边缘强度(Edge strength, ES)和平均梯度(Average gradient, AG)[19]。对不同算法的重建结果进行评价结果如表2所示,表中有效数字为三位。

      表 2  不同算法重建结果的NIQE、PSNR、SSIM值、边缘强度和平均梯度

      Table 2.  NIQE, PSNR, SSIM, edge intensity and average gradient value of the reconstruction results of different algorithms

      Measureslanczos2NEDINSL0ICNSL0
      NIQE8.6486.3527.4996.481
      PSNR27.71314.43922.872 421.433
      SSIM0.897 10.2660.3040.864 2
      ES16.70220.74119.59622.828
      AG1.5251.9611.7962.096

      对于NIQE的主观评价指标,该指标越低表示人眼对图像主观感受越好,从表中数据可以看出ICNSL0算法相比lanczos2算法降低了2.167。这表明ICNSL0算法在太赫兹图像重建的结果中具有更好的主观感知。

      基于统计学的客观评价指标PSNR和SSIM对比结果显示,lanczos2算法在这两项评价指标中表现出较好的特性。这两种评价指标的原理为计算“参考图像”与待测图像之间的区别,数值越高表示两个图像之间的相似度越高。而参考图像是利用双三次线性插值算法构造而来,因此在这两项指标的对比中ICNSL0算法相比其他算法不存在明显优势。对于以边缘信息为主的边缘强度和平均梯度的两种衡量标准的对比下,基于插值的lanczos2算法和NEDI算法降低了重建图像的锯齿效应,但削弱了图像的高频信息。ICNSL0算法与NSL0算法相比在边缘强度上提高了3.232,在平均梯度上提高了0.300,在一定程度上提高了重建精度。ICNSL0算法针对不同区域的图像进行混合重建融合时避免了削弱边缘部分的表示,因此可以在一定程度上提高图像边缘强度与平均梯度。

      综合上述三类评价指标的结果表明,ICNSL0算法与其他算法相比,对单幅太赫兹图像的超分辨重建存在一定的有效性与优越性。

    • 以搭建精简化、低成本的太赫兹成像系统为目的,文中提出利用单个连续波太赫兹源和线性探测器结合可调控旋转台实现太赫兹波的快速线阵扫描,利用该旋转台匀速移动对样品进行扫描获得单幅低分辨率太赫兹图像。由于原始图像的像素值较大,提出了双维度分离的预处理方法,有效地减少了图像处理中的占用内存。引入CLAHE算法对原始图像进行分区域对比度提升从而解决太赫兹图像对比度低的问题。利用稀疏表示和字典学习的方法对单幅太赫兹图像进行超分辨重建,解决原始图像分辨率较低的问题。文中提出ICNSL0算法解决稀疏表示系数的求解问题,使用反余割函数拟合零范数进行优化求解。与传统的lanczos2插值算法、以边缘为指导的NEDI算法和基于稀疏表示的NSL0算法相比,ICNSL0算法重建的太赫兹图像边缘细节更加清晰,样品与背景旋转台的对比度有所提升。从主观和客观的双重衡量标准对图像重建结果进行评价,ICNSL0算法重建的图像更符合人类主观视觉的观察习惯,其边缘信息更加细腻突出,更适合于实际成像系统采集的单幅太赫兹图像的超分辨率重建。

参考文献 (19)

目录

    /

    返回文章
    返回