一种太赫兹时域光谱稀疏成像方法

摘要

本发明公开了一种太赫兹时域光谱稀疏成像方法，利用移不变小波变换对太赫兹图像进行分解，另外为了增强太赫兹图像小波变换系数的稀疏性，即增强获得的小波系数中大系数值的同时抑制小系数值的影响，本发明对经过移不变小波变换后获得的小波系数利用指数函数进行非线性变换；本发明采用分裂迭代算法将构建的太赫兹时域光谱稀疏成像优化代价函数分裂为一系列子问题分别进行迭代优化，提高重建质量的同时能快速获得太赫兹时域光谱成像结果。

主权项

一种太赫兹时域光谱稀疏成像方法，其特征在于，包括以下步骤：步骤1：构建太赫兹时域光谱成像系统模型1a)假设太赫兹二维图像对象x有m×n个像素点，令N＝m×n，从这N个像素点中随机选取M个位置的频域数据进行太赫兹检测，则在每个检测位置，都将获得一个完整的太赫兹时域波形数据；1b)选取每个太赫兹时域波形的峰值作为对应检测像素点的像素值，则对应M个检测位置将获得一个M维的检测数据y；1c)由于太赫兹检测数据y是对成像对象的频域数据进行随机采样获得的，因此可构建太赫兹成像系统模型为y＝RFx   (1)其中y为M维太赫兹检测数据，x为二维成像对象，其m×n个元素排列在一个N维列向量里；F代表傅里叶变换矩阵，R为由随机采样构成的测量矩阵；步骤2：利用指数移不变小波对二维太赫兹图像进行稀疏变换2a)由于正交小波基是由基本小波函数经过平移与伸缩变换得到的函数集，平移是非一致取样，随着尺度增大，位移取样间隔随2的幂次增大，不能从多尺度角度很好的匹配信号的局部结构特征，因此在信号急剧变化部分可能会产生振荡现象和块效应；为了有效地消除这种人为的振荡现象，利用移不变小波变换对太赫兹图像进行分解，得到太赫兹图像的稀疏表示为r＝W(x)＝Ψ^‑1x   (2)其中W(x)表示对图像x进行移不变小波变换，Ψ表示移不变小波变换矩阵，r表示小波系数，其为稀疏信号；2b)为增强太赫兹图像小波变换系数的稀疏性，即增强获得的小波系数中大系数值的同时抑制小系数值的影响，对经过移不变小波变换后获得的小波系数利用指数函数进行非线性变换；首先将获得的小波系数范围进行归一化处理，之后利用指数函数进行运算，得到指数移不变小波变换系数为$r_e=W_e\left(x\right)=\frac{a^r-1}{a-1}=\frac{a^{{\Psi }^{-1}x}-1}{a-1}---\left(3\right)$其中W_e(x)表示对图像x进行指数移不变小波变换；a为大于1的常数，本发明中a的取值为10；步骤3：构建太赫兹时域光谱稀疏成像优化代价函数由于经过指数移不变小波变换后获得的小波系数具有稀疏性，这里采用1范数对小波系数进行约束，可构建太赫兹时域光谱稀疏成像代价函数为如下有约束最优化问题$\begin{array}{ccc}\underset{x}{min}\left|\right|W_{e}x|{|}_1& s.t.& y=RFx\end{array}---\left(4\right)$通过拉格朗日乘子法，将该有约束最优化问题转化为无约束最优化问题$\underset{x}{min}\left|\right|W_{e}x|{|}_1+\frac{\mu }{2}||RFx-y|{|}_2^2---\left(5\right)$步骤4：基于分裂迭代算法，对构建的太赫兹时域光谱稀疏成像优化代价函数进行求解，得到太赫兹时域光谱稀疏成像结果为提高重建质量的同时快速获得太赫兹时域光谱成像结果；本发明采用分裂迭代算法将无约束最优化问题(5)分裂为子问题分别进行迭代优化；具体的成像步骤为：4a)输入太赫兹检测数据y，测量矩阵R，傅里叶变换矩阵F和指数移不变小波变换函数W_e；令r_e＝W_ex，基于分裂迭代算法，对太赫兹时域光谱稀疏成像代价函数(5)进行优化，等价为如下一系列子问题分别求解$x^{i+1}=\underset{x}{min}\frac{\mu }{2}\left|\right|RFx-y|{|}_2^2+\frac{\lambda }{2}||r_e^i-W_{e}x-b_r^i|{|}_2^2---\left(6\right)$$r_e^{i+1}=\underset{r_e}{min}\left|\right|r_e|{|}_1+\frac{\lambda }{2}||r_e-W_e{x}^{i+1}-b_r^i|{|}_2^2---\left(7\right)$$b_r^{i+1}=b_r^i+(W_e{x}^{i+1}-r_e^{i+1})---\left(8\right)$4b)初始化信号估计x⁰＝F^‑1y，迭代次数i＝0；4c)第i次迭代时，求解最优化问题(6)对信号x的估计值进行更新；令式(6)的代价函数对x进行微分并等于0，可得到$({\mathrm{\mu F}}^T{R}^{T}RF+\lambda I)x^{i+1}={\mathrm{\mu F}}^T{R}^{T}y+{\mathrm{\lambda W}}_e^T(r_e-b_r)$则$x^{i+1}={({\mathrm{\mu F}}^T{R}^{T}RF+\lambda I)}^{-1}({\mathrm{\mu F}}^T{R}^{T}y+{{\mathrm{\lambda W}}_e}^T(r_e-b_r\left)\right)$4d)利用更新的估计值xⁱ⁺¹，求解最优化问题(7)对信号r_e的估计值进行更新；利用shrinkage操作获得式(7)的解为：$r_e^{i+1}=shrink(W_e{x}^{i+1}+b_r^i,1/\lambda )$其中$shrink(x,\lambda )=\frac{x}{\left|x\right|}\times max\left(\right|x|-\lambda ,0)$4e)将更新后的值xⁱ⁺¹和代入式(8)对b_r的取值进行更新，获得新的估计值4f)判断||xⁱ⁺¹‑xⁱ||₂是否小于预先设置的迭代终止门限δ，若大于，则令i＝i+1，返回步骤(4c)；否则，迭代终止，输出最终估计结果xⁱ⁺¹；4g)将获得的列向量xⁱ⁺¹重新排列为矩阵形式，即为最终的太赫兹时域光谱成像结果。