A Fast Carrier Frequency Offset Position Detection Algorithm for Passive Backscatter Communication System
-
摘要: 近来无源反向散射通信技术作为绿色物联网的关键技术,引起了广泛的关注。在无源反向散射通信系统中,收发节点与反射节点的振荡器差异、相对运动以及环境变化,导致接收端和发送端存在载波频率偏移(CFO)。CFO对信号检测和系统性能有重要影响,而当前多数无源反向散射通信系统研究忽略了CFO。该文设计了一种适用于频移键控调制(FSK)的CFO快速检测方法,不需要导频就能快速有效检测出CFO是否存在并找出存在的位置。首先,根据信号在CFO存在与否的不同,采用去直流后取模方法对信号进行处理,而后根据处理后的信号的特点,基于累积和(CUSUM)设计快速检测算法对载波频偏出现的位置进行检测,并对理论分析结果进行仿真验证,仿真结果表明取模检测可以有效检测出CFO出现的位置。Abstract: Recently, passive backscatter communication technology has attracted extensive attention as one key technology for green Internet of Things. In one typical backscatter communication system, there may exist Carrier Frequency Offset (CFO) between the receiver and the transmitter due to the relative motion or the difference in oscillators or environmental changes. CFO has an important impact on signal detection and system performance, but most current studies ignore the CFO. In this paper, a fast CFO detection method suitable for Frequency Shift Keying (FSK) modulation is designed, which can quickly and effectively detect without pilots whether there exists CFO and find out the location where the CFO begins. Specifically, this paper designs one detector amplitude-taking method. Next, based on CUmulative SUM (CUSUM) algorithm a fast detection algorithm is designed to detect the location of CFO. Finally, simulation results are provided to corroborate the proposed studies. The simulation results show that the designed detector can effectively locate the position where CFO appears.
-
1. 引言
合成孔径雷达(Synthetic Aperture Radar , SAR)是最早提出并投入实用的成像雷达。合成孔径雷达可实现全天时、全天候的区域监测成像,且对植被等介质具有穿透能力,在灾害评估、环境监测和军事等领域得到了广泛的应用[1-3]。
由雷达分辨率理论和奈奎斯特采样定理可知,合成孔径雷达系统性能的提高通常伴随着采样数据量的显著增加,这给系统设计和实现带来困难。近年来,随着压缩感知[4]和相关理论的发展,稀疏信号处理在SAR成像中得到广泛应用,这表明当回波信号具有稀疏性或可压缩性时,便能够以远低于奈奎斯特的采样频率,用较少的观测数据高概率地恢复出原信号。现有的压缩感知SAR成像模型中,大致分为两大类,其一,单独考虑距离向和方位向的稀疏模型,然后与匹配滤波相结合进行成像。例如文献[5]仅在距离向实现了压缩感知技术,文献[6]仅在方位向实现了压缩感知技术,上述方法都可以实现SAR成像,但也存在一些问题。首先,仍然需要处理较大的数据量且没有充分利用空间的稀疏性;其次,成像过程中用到的匹配滤波方式导致很大的局限性。其二,在距离向和方位向同时构建压缩感知模型进行求解[7-9],其中文献[7]将接收到的雷达回波和2维后向散射系数矩阵叠加成非常大的列向量,然后进行方位向距离向压缩感知处理,虽然比仅对距离向和方位向处理更进一步减少了采样数据量,但求解过程需要对大的列向量进行处理,大幅增加大数据运算量与信号的重构时间。针对出现的问题,文献[8,9]将接收到的回波与后向散射系数保持为2维矩阵,而不是堆叠成向量,因此可以在降低采样率的同时很好地重建目标信号。不管是求解第1类成像模型还是第2类成像模型,都需要用到重构算法,最终目的都是解决l0范数最小化问题。传统的压缩感知重构算法包括凸优化算法、贪婪算法和统计优化算法3大类,例如文献[10]实现了基于随机孔径贝叶斯压缩感知理论的高分辨合成孔径雷达成像,这里用到了统计优化类算法,文献[11]将隶属于贪婪算法类的正交匹配追踪(Orthogonal Matching Pursuit, OMP)算法应用到合成孔径雷达成像中,文献[12]将求解凸优化问题的快速迭代收缩算法(Fast Iterative Shrinkage Thresholding Algorithm, FISTA)应用到压缩感知ISAR成像的过程中。另有一些学者将非凸函数族中的平滑l0算法应用到合成孔径雷达方位向稀疏模型求解过程中。但这些重构算法都存在重构精度不高、抗噪性能差等问题。
针对上述问题,本文在文献[8,9]构建的压缩感知SAR成像模型的基础上,将非凸非平滑优化算法引入到压缩感知SAR成像模型的求解问题中,提出一种基于迭代近端投影的2维欠采样合成孔径雷达成像重建算法。首先将压缩感知SAR成像模型转化为近端函数[13,14]稀疏信号模型,然后采用SCAD判罚函数[15]来进一步求解该模型,进而重构出图像。仿真与实测数据处理结果表明,本文方法具有很好的重构性能,具有较好的工程应用前景。
2. SAR回波信号2维稀疏表示模型
图1为SAR成像几何模型,假设雷达平台沿
y 轴以恒定的速度v 前行,若雷达发射线性调频(Linearly Frequency Modulated, LFM)信号,则理想点目标P(x,y) 的基带回波为s(tr,ta;r)=∬Dg(x,y)aa(ta)ar(tr−2R(ta;r)c)⋅exp[jπγ(tr−2R(ta;r)c)2]⋅exp[−j4πλR(ta;r)]dxdy (1) 其中,
D={(x,y)|xmin ,R({t_{\text{a}}};r) 为{t_{\text{a}}} 时刻雷达与目标间的距离,其表达式为R({t_{\rm{a}}};r){\text{ = }}\sqrt {{r^2} + {v^2} \cdot {{({t_{\rm{a}}} - {t_{{{\rm{a}}_0}}})}^2}} \approx r + \frac{{{v^2}}}{{2r}}{\left({t_{\rm{a}}} - {t_{{{\rm{a}}_0}}}\right)^2} (2) 其中,
{t_{\rm{r}}} 和{t_{\rm{a}}} 分别为距离向快时间和方位向慢时间;r 为点目标P(x,y) 到雷达航线的垂直距离;{\boldsymbol{g}}\left( {x,y} \right) 为目标的后向散射系数矩阵;a{}_{\rm{a}}( \cdot ) 为方位向窗函数;{a_{\rm{r}}}({t_{\rm{r}}}) = {\rm{rect}}\left(\dfrac{{{t_{\rm{r}}}}}{{{T_{\rm{p}}}}}\right) 为发射信号包络,{T_{\rm{p}}} 为脉冲持续时间;c 为光速;\gamma 为发射的LFM信号的线性调频率;\lambda 为波长。根据文献[8]和文献[9]知,理想点目标的回波经采样后可以表示为{\boldsymbol{S}} = {{\boldsymbol{\varPhi}} _{\rm{r}}}{\boldsymbol{G}}{{\boldsymbol{\varPhi}} _{\rm{a}}} (3) 其中,
{\boldsymbol{S}} \in {\mathbb{C}^{M \times N}} 表示均匀采样回波信号;{\boldsymbol{G}} \in {\mathbb{C}^{M \times N}} 表示观测场景目标点的后向散射矩阵;矩阵{{\boldsymbol{\varPhi}} _{\text{r}}} \in {\mathbb{C}^{M \times M}} 为距离向的卷积矩阵,其列是距离向参考函数的共轭与快时间域中采样数的卷积,其中距离向参考函数可以表示为{s_{\rm{r}}}({t_{\rm{r}}}) = {a_{\rm{r}}}\left( {{t_{\rm{r}}}} \right)\exp ({\rm{j}}\pi \gamma t_r^2) ;矩阵{{\boldsymbol{\varPhi}} _{\text{a}}} \in {\mathbb{C}^{N \times N}} 为方位向的卷积矩阵,其行是方位向参考函数的共轭与慢时间域中采样数的卷积,其中方位向参考函数可以表示为{s_{\rm{a}}}({t_{\rm{a}}}) = {a_{\rm{a}}}\left( {{t_{\rm{a}}}} \right)\exp \left( {{\rm{j}}\pi {k_{\rm{a}}}t_{\rm{a}}^2} \right) ,其中{k_{\rm{a}}} = \dfrac{{2{v^2}r}}{\lambda } 为多普勒调频率;\mathbb{C} 表示复数域;M 为距离向采样点数;N 为方位向采样点数。对飞机、桥梁等具有空域稀疏性的场景,首先,构建稀疏随机测量矩阵来减少2维采样量,随机测量矩阵构建方式为随机选取方位向和距离向部分采样点信号,如图2所示,其中灰色部分代表采样点信号。此时均匀采样回波信号
{\boldsymbol{S}} \in {\mathbb{C}^{M \times N}} 变为随机采样的回波信号{{\boldsymbol{S}}_{2{\rm{CS}}}} \in {\mathbb{C}^{P \times Q}} ,P\left( {P \le M} \right) 为距离向随机采样点数,Q\left( {Q \le N} \right) 为方位向随机采样点数。与此同时,距离向卷积矩阵{{\boldsymbol{\varPhi}} _{\rm{r}}} 变为{{\boldsymbol{\varPhi}} _{{\rm{RCS}}}} \in {\mathbb{C}^{P \times M}} ,方位向卷积矩阵{{\boldsymbol{\varPhi}} _{\rm{a}}} 变为{{\boldsymbol{\varPhi}} _{{\rm{ACS}}}} \in {\mathbb{C}^{N \times Q}} ,通过上述整理可得,2维随机降采样合成孔径雷达回波信号的欠定方程表示为{{\boldsymbol{S}}_{2{\rm{CS}}}}{\text{ = }}{{\boldsymbol{\varPhi}} _{{\rm{RCS}}}}{\boldsymbol{G}}{{\boldsymbol{\varPhi}} _{{\rm{ACS}}}} (4) 考虑到现实场景中加性噪声的影响,式(4)可以表示为
{{\boldsymbol{S}}_{2{\rm{ZCS}}}}{\text{ = }}{{\boldsymbol{\varPhi}} _{{\rm{RCS}}}}{\boldsymbol{G}}{{\boldsymbol{\varPhi}} _{{\rm{ACS}}}} + {{\boldsymbol{N}}_Z} (5) 其中,
{{\boldsymbol{S}}_{2{\rm{CS}}}} \in {\mathbb{C}^{P \times Q}},{{\boldsymbol{S}}_{2{\rm{ZCS}}}} \in {\mathbb{C}^{P \times Q}} 为随机采样的回波信号,{{\boldsymbol{\varPhi}} _{{\rm{RCS}}}} \in {\mathbb{C}^{P \times M}} 为距离维观测矩阵,{{\boldsymbol{\varPhi}} _{{\rm{ACS}}}} \in {\mathbb{C}^{N \times Q}} 为方位维观测矩阵,{\boldsymbol{G}} \in {\mathbb{C}^{M \times N}} 为待重构的目标图像,{{\boldsymbol{N}}_Z} \in {\mathbb{C}^{P \times Q}} 为复高斯白噪声信号。3. 基于迭代近端投影的合成孔径雷达成像方法
由上述章节描述可知,式(5)所表示的信号模型是2维卷积模型,它可以用两个连续的1维方程表示。首先,令
{\boldsymbol{X}} \in {\mathbb{C}^{M \times Q}} 作为方位维压缩数据。其表达式为{\boldsymbol{X}} = {\boldsymbol{G}}{{\boldsymbol{\varPhi}} _{{\rm{ACS}}}} (6) 则式(5)可以写为
{{\boldsymbol{S}}_{{\text{2ZCS}}}} = {{\boldsymbol{\varPhi}} _{{\rm{RCS}}}}{\boldsymbol{X}} + {{\boldsymbol{N}}_Z} (7) 式(6)表示方位向压缩感知模型,式(7)表示距离向压缩感知模型。
求解式(7)和式(6)的过程如下,首先,求解欠定方程式(7)得到距离向压缩后的信号,其次,将接收到的回波
{{\boldsymbol{S}}_{{\text{2ZCS}}}} 分解为列向量的形式,其表达式为{{\boldsymbol{S}}_{{\text{2ZCS}}}} = \left[ {{s_1},{s_2}, \cdots ,{s_Q}} \right] (8) 同时,
{\boldsymbol{X}} 也可以分解为列向量形式{\boldsymbol{X}} = \left[ {{{\boldsymbol{x}}_1},{{\boldsymbol{x}}_2}, \cdots ,{{\boldsymbol{x}}_Q}} \right] (9) 那么,式(7)可以表示为
{{\boldsymbol{s}}_q} = {{\boldsymbol{\varPhi}} _{{\rm{RCS}}}}{{\boldsymbol{x}}_q}{\kern 1pt} {\kern 1pt} + {{\boldsymbol{n}}_z},{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} q = 1,2, \cdots ,Q (10) 其中,
{{\boldsymbol{s}}_q} \in {\mathbb{C}^P} ,{{\boldsymbol{x}}_q} \in {\mathbb{C}^M} ,{{\boldsymbol{n}}_z} \in {\mathbb{C}^P} ,Q 为方位向接收到的脉冲数。由于式(10)是欠定方程,其对应的优化问题为式(11)的拉格朗日(Lagrangian)形式[16],\mathop {\min }\limits_{{x_q}} {\left\| {{x_q}} \right\|_0}, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\text{s}}{\text{.t}}{\text{.}}{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\left\| {{s_q} - {{\boldsymbol{\varPhi}} _{{\rm{RCS}}}}{x_q}} \right\|_2} \le \xi (11) 其中,
{\left\| \cdot \right\|_0} 是l0范数,它表示向量{x_q} 中非零元素个数;{\Vert ·\Vert }_{2} 是l2范数;\xi 是一个与噪声水平有关的参数。众所周知,求解式(11)最小化问题是一个非凸非平滑优化问题,传统基于压缩感知的SAR成像方法是在一定条件下将l0范数最小化的非凸非平滑优化问题转化为凸或平滑问题进行求解,这会使稀疏模型不可避免地存在一定误差,进而导致SAR成像效果不佳。而本文将能够解决非凸非平滑问题的迭代近端投影的方法应用到合成孔径雷达成像中,具体过程如下。由文献[13,14]知,为求得式(10)稀疏解,可将其转化为近端函数优化问题,其稀疏信号表示模型为\mathop {\min }\limits_{{x_q},z} F\left( z \right) + {\delta _\mathcal{A}}_{_\varepsilon }\left( {{x_q}} \right),{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\text{s}}{\text{.t}}{\text{.}}{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} z = {x_q}{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} (12) 其中,
F\left( \cdot \right) 是一个类似l0范数的非凸非平滑函数,z 是一个列向量的辅助变量,{\delta _{{\mathcal{A}_\varepsilon }}} 定义为{\mathcal{A}_\varepsilon } \triangleq \{ {x_q}:{\left\| {{s_q} - {{\boldsymbol{\varPhi}} _{{\rm{RCS}}}}{{\boldsymbol{x}}_q}} \right\|_2} \le \varepsilon \} 可行集的指示函数[13],即{\delta _{{\mathcal{A}_\varepsilon }}} = \left\{ {\begin{array}{*{20}{c}} {0,}&{{{\boldsymbol{x}}_q} \in {\mathcal{A}_\varepsilon }} \\ { + \infty ,}&{{{\boldsymbol{x}}_q} \notin {\mathcal{A}_\varepsilon }} \end{array}} \right. (13) 其中,
{\Vert ·\Vert }_{2} 表示l2范数,式(12)的稀疏信号表示模型中引入惩罚函数方法,即\mathop {\min }\limits_{{x_q},z} F\left( z \right) + {\delta _\mathcal{A}}_{_\varepsilon }\left( {{x_q}} \right){\kern 1pt} + \frac{1}{{2\alpha }}\left\| {{x_q} - z} \right\|{\kern 1pt} _2^2{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} (14) 其中,
\alpha > 0 为惩罚因子,通过交替最小化方法[17]解决式(14)所表示问题。表达式为\left. \begin{aligned} & {{{\boldsymbol{z}}_{k + 1}} = \arg {{\min }_z}\left\{ \alpha F\left( {\boldsymbol{z}} \right) + \frac{1}{2}\left\| {{\boldsymbol{z}} - {{\boldsymbol{x}}_{q,z}}} \right\|_2^2\right\} } \\ & {{{\boldsymbol{x}}_{q,z + 1}} = \arg {{\min }_{{x_q}}}\left\{ {\delta _{{\mathcal{A}_\varepsilon }}} + \frac{1}{2}\left\| {{{\boldsymbol{x}}_q} - {{\boldsymbol{z}}_{k + 1}}} \right\|_2^2\right\} } \end{aligned} \right\} (15) 根据文献[13]可知,式(15)可以进一步写为式(16)和式(17)
\left. \begin{aligned} \quad & {{{\boldsymbol{z}}_{k + 1}} = {\text{pro}}{{\text{x}}_{\alpha ,F}}({{\boldsymbol{x}}_{q,k}})} \\ & {{{\boldsymbol{x}}_{q,k + 1}} = {\text{pro}}{{\text{x}}_{{\delta _{{\mathcal{A}_\varepsilon }}}}}({{\boldsymbol{z}}_{k + 1}})} \end{aligned} \right\} (16) \begin{split} \qquad\;\;\quad {{\boldsymbol{x}}_{q,k + 1}} & = {\text{pro}}{{\text{x}}_{{\delta _{{\mathcal{A}_\varepsilon }}}}}\left( {{\text{pro}}{{\text{x}}_{\alpha ,F}}({{\boldsymbol{x}}_{q,k}})} \right) \\ & = {\mathcal{P}_{{\mathcal{A}_\varepsilon }}}({\text{pro}}{{\text{x}}_{\alpha ,F}}({{\boldsymbol{x}}_{q,k}})) \end{split} (17) 其中,
{\text{pro}}{{\text{x}}_{\alpha ,F}}({{\boldsymbol{x}}_{q,k}}) = \mathop {\arg \min }\limits_{{x_q}} \Bigr\{ \dfrac{1}{2}\left\| {{\boldsymbol{z}} - {{\boldsymbol{x}}_{q,k}}} \right\|_2^2 + \alpha F\left( {\boldsymbol{z}} \right)\Bigr\} 是F\left( z \right) 的近端算子,即求解{{\boldsymbol{x}}_{q,k}} 在满足约束条件下的最小的变量值z ,其核心思想就是将凸或非凸函数的优化问题转化为易求解的近端映射函数,从而实现近似解。{\mathcal{P}_{{\mathcal{A}_\varepsilon }}} 代表可行集{\mathcal{A}_\varepsilon } 的投影。此外,为了提高算法的性能和更好地解决非凸问题存在多个局部极小值的问题,本文利用外推步骤来改善算法性能,即\left. \begin{aligned} & {{{\bar {\boldsymbol{x}}}_{q,k}} = {{\boldsymbol{x}}_{q,k}} + \omega \left( {{{\boldsymbol{x}}_{q,k}} - {{\boldsymbol{x}}_{q,k - 1}}} \right)} \\ & {{{\bar {\boldsymbol{x}}}_{q,k + 1}} = {\mathcal{P}_{{\mathcal{A}_\varepsilon }}}\left( {{\text{pro}}{{\text{x}}_{\alpha ,F}}({{\bar {\boldsymbol{x}}}_{q,k}})} \right)} \end{aligned} \right\} (18) 其中,
\omega \ge 0 为权重常数。由文献[13]可知,非平滑函数的近端算子{\text{pro}}{{\text{x}}_{\alpha ,F}}({\bar {\boldsymbol{x}}_{q,k}}) 可由SCAD惩罚函数F\left( {{{\bar {\boldsymbol{x}}}_{q,k}}} \right) 产生相应的SCAD阈值函数T_{\lambda ,a}^{{\rm{scad}}} 来获得,该阈值函数可以进一步减小软、硬阈值收缩函数带来的偏差和促进解的稀疏性。因此,为进一步求得更优解,非平滑函数的近端算子{\text{pro}}{{\text{x}}_{\alpha ,F}}({\bar {\boldsymbol{x}}_{q,k}}) 可以表示为\begin{split} {\text{pro}}{{\text{x}}_{\alpha ,F}}({\bar x_{q,k}}) &= \mathop {\arg \min }\limits_{{{\bar x}_q}} \left\{ {\frac{1}{2}\left\| {{{\bar x}_q} - {{\bar x}_{q,k}}} \right\|_2^2 + \alpha F\left( {{{\bar x}_q}} \right)} \right\} \\ &= T_\lambda ^{{\rm{scad}}}\left( {{{\bar x}_q}} \right)\\[-10pt] \end{split} (19) 其中,惩罚函数
F\left( {{{\bar {\boldsymbol{x}}}_{q,k}}} \right) 的表达式为\begin{split} & F\left( {{{\bar {\boldsymbol{x}}}_{q,k}}} \right){\text{ = }}\\ & \left\{ \begin{aligned} & \lambda \left| {{{\bar {\boldsymbol{x}}}_{q,k}}} \right|,\qquad\qquad\qquad\qquad\qquad 0 \le \;\left| {{{\bar {\boldsymbol{x}}}_{q,k}}} \right| < \lambda \hfill \\ & - \frac{{\left( {{{\left| {{{\bar {\boldsymbol{x}}}_{q,k}}} \right|}^2} - 2a\lambda \left| {{{\bar {\boldsymbol{x}}}_{q,k}}} \right| + {\lambda ^2}} \right)}}{{2\left( {a - 1} \right)}},\; \lambda \le \;\left| {{{\bar {\boldsymbol{x}}}_{q,k}}} \right| < a\lambda \hfill \\ & \frac{{\left( {a + 1} \right){\lambda ^2}}}{2},\qquad\qquad\qquad\qquad\;\; \left| {{{\bar {\boldsymbol{x}}}_{q,k}}} \right| \ge a\lambda \end{aligned} \right. \end{split} (20) SCAD函数的近端映射为
\begin{split} & {T}_{\lambda ,a}^{{\rm{scad}}}\left({\overline{{\boldsymbol{x}}}}_{q}\right)=\\ & \left\{ \begin{aligned} & {\rm{sign}}\left({\overline{{\boldsymbol{x}}}}_{q}\right){\left(\left|{\overline{{\boldsymbol{x}}}}_{q}\right|-\lambda \right)}_+,\qquad \left|{\overline{{\boldsymbol{x}}}}_{q}\right|\le 2\lambda \\ & \frac{\left(a-1\right){\overline{{\boldsymbol{x}}}}_{q}-{\rm{sign}}\left({\overline{{\boldsymbol{x}}}}_{q}\right)a\lambda }{a-2},2\lambda < \left|{\overline{{\boldsymbol{x}}}}_{q}\right|\le a\lambda \\ & {\overline{{\boldsymbol{x}}}}_{q},\qquad \qquad \qquad \qquad \quad\;\; \left|{\overline{{\boldsymbol{x}}}}_{q}\right| > a\lambda \end{aligned}\right. \end{split} (21) 其中,
\lambda 为调整参数,a 为常量且a > 2 ,{\rm{sign}}\left( \cdot \right) 为符号函数,{\left(\beta \right)}_{\text{+}}\text{=}\mathrm{max}\left(\beta ,0\right) 。通过求解式(18)得到稀疏向量{\bar {\boldsymbol{x}}_q} ,至此,由{\bar{\boldsymbol{ x}}_q} 构成的{\boldsymbol{X}} 就是经过距离向压缩感知得到的1维距离像。然后,再对方位向压缩感知模型式(6)进行求解,为了更好地表述,将式(6)进行转置得{{\boldsymbol{X}}^{\rm{T}}} = {\boldsymbol{\varPhi}} _{{\rm{ACS}}}^{\rm{T}}{{\boldsymbol{G}}^{\rm{T}}} (22) 将
{{\boldsymbol{X}}^{\rm{T}}} 表示为列向量形式即{{\boldsymbol{X}}^{\rm{T}}} = \left[ {{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{{\boldsymbol{x}}} }_1},{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{{\boldsymbol{x}}} }_2} \cdot \cdot \cdot ,{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{{\boldsymbol{x}}} }_M}} \right] ,将{{\boldsymbol{G}}^{\rm{T}}} 也表示为列向量形式即{{\boldsymbol{G}}^{\rm{T}}} = \left[ {{{\boldsymbol{g}}_{_1}},{{\boldsymbol{g}}_2}, \cdots ,{{\boldsymbol{g}}_M}} \right] ,那么欠定方程式(22)可以进一步写为{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{{\boldsymbol{x}}} _m} = {\boldsymbol{\varPhi}} _{{\rm{ACS}}}^{\rm{T}}{{\boldsymbol{g}}_{m,{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} }}{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} m = 1,2, \cdots ,M (23) 其中,
{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{{\boldsymbol{x}}} _m} \in {\mathbb{C}^Q} ,{{\boldsymbol{g}}_m} \in {\mathbb{C}^N} ,M 为距离向采样点。对方位向压缩感知模型也使用本文提出的迭代近端投影的方法求解,求解得到稀疏向量{{\boldsymbol{g}}_m} ,进一步得到2维的SAR图像。4. 仿真结果与分析
为综合比较分析本文方法的性能,将正交匹配追踪(OMP)算法、平滑l0范数(Smoothed l0 norm, SL0)算法、贝叶斯压缩感知(Bayesian Compressive Sensing, BCS)算法与本文提出的迭代近端投影(Iterative Proximal Projection, IPP)算法作为对比;同时,为进一步定量分析成像结果,本文利用目标杂波比(Target to Clutter Ratio, TCR)[18]和图像熵(Image Entropy, IE)[18]来定量评价成像质量,用成像时间来评定算法的运算效率。同时,为验证算法在噪声情况下的性能,对回波加入高斯白噪声进行仿真实验。
4.1 不同采样率成像实验
为进一步验证算法的性能,对30个点目标仿真产生的数据在不同采样率和不同信噪比条件下进行仿真实验。主要仿真参数设置如表1所示。为方便分析,设置相同的距离采样率和方位采样率
\delta 。图3给出了目标散射点模型与不同采样率下各算法成像结果图,其中图3(a)为散射点模型,图3(b-e)表示采样率为原采样率1/2时4种算法成像结果。图3(f-i)表示采样率为原采样率1/4时4种算法成像结果。表 1 雷达仿真参数参数 数值 雷达信号载频 3 GHz 雷达信号带宽 150 MHz 采样频率 300 MHz 雷达距目标区域中心点 4200 m 由图3可以看出,在采样率为原采样率1/2时,几种算法都可以清晰地实现目标点SAR成像;当采样率为1/4时,OMP算法的成像效果变差,其他算法依然能够清晰成像。
为了进一步比较上述算法的成像质量,利用TCR和IE作为衡量指标。图4(a)、图4(b)分别给出了不同采样率条件下几种算法的TCR及IE对比曲线,从图4(a)可以看出,随着采样率的下降,SL0算法、BCS算法与本文算法的TCR相对保持稳定,而OMP算法在采样率较低时产生较大波动。同时,可以看出,本文算法的TCR一直大于其他算法,说明本文算法成像结果聚焦性能更好,虚假散射点更少。从图4(b)可以看出,本文算法的IE一直低于其他算法,进一步说明本文算法具有很好的成像性能。
为进一步验证算法的运算效率,表2给出了采样率同为原采样率的1/2时几种算法成像时间对比,从中进一步得出本文算法成像时间也相对较短。说明本文算法在保证成像质量的前提下能够实现快速成像。
表 2 采样率为原采样率1/2时各算法成像时间算法 时间(s) IPP 26.7381 OMP 18.7749 SL0 27.3091 BCS 178.2851 4.2 抗噪声性能实验
首先在不同信噪比条件下进行数据仿真实验,将距离向采样率和方位向采样率同时设置为原采样率的1/2,即
\delta = \dfrac{1}{2} ,并在回波信号处添加高斯白噪声使SNR分别为–10 dB和20 dB。成像结果如图5,其中图5(a)—图5(d)为SNR=20 dB时4种算法成像结果,图5(e)—图5(h)为SNR=–10 dB时4种算法成像结果。由图5可以看出,当信噪比为20 dB时,几种算法都可以清晰成像;当采样率为–10 dB时,OMP算法与BCS算法成像结果出现严重失真,SL0算法成像结果虽可以成像,但是出现较多虚假目标,只有本文算法受噪声影响较小,可以较为清晰地成像。为了进一步验证算法的鲁棒性,以TCR, IE作为衡量指标。图6(a)、图6(b)分别为不同信噪比条件下各算法的TCR与IE曲线,由图6(a)可以看出,随着信噪比的上升,各算法的TCR逐渐上升,但在较低信噪比时本文算法依然高于其他算法,说明本文算法有一定的抗噪声能力。同时,由图6(b)可以得出,本文算法的IE随着信噪比的上升逐渐下降且一直低于其他算法,进一步说明本文算法具有良好的鲁棒性。综合实验对比可知,本文算法相较于其他算法具有更好的欠采样成像能力。
5. 实测数据处理
为了充分验证本文算法对实际场景处理的有效性,本文将采用来自加拿大温哥华RADARSAT-1精细模式2(加拿大航天局版权所有)的实测数据进行实验。主要参数设置如表3所示。图7给出各算法成像结果,其中图7(a)为基于距离多普勒算法的温哥华地区英吉利海湾附近区域的实测场景图。本文将选取图7(a)中红色矩形区域中6个目标点(即6艘货船)进行实验,其原因是这6个目标点相较于海平面来说是强散射点,而又相较于英吉利海湾是稀疏的。图7(b)—图7(e)是4种算法在距离采样率为0.5,方位采样率为0.3时的成像结果。成像结果表明,当采样数据量减少后,本文算法仍能重建出目标,且相较于其他算法能够消除图像中的模糊,进一步验证了本文算法的有效性。
表 3 温哥华场景RADARSAT-1参数参数 数值 距离带宽 30.3 MHz 距离向采样频率 32.317 MHz 脉冲宽度 30.111 MHz 卫星轨道半径 7186029 m 雷达波长 0.05657 m 6. 结束语
针对RD等成像算法存在数据量大、采样率高及OMP、SL0和BCS等算法存在成像精度低及抗噪性能差等问题,本文提出一种迭代近端投影的2维欠采样SAR成像方法,首先对SAR回波信号进行分析,构建基于距离向与方位向的稀疏表示模型,在此基础上,利用近端函数优化模型来表示距离向与方位向稀疏模型,并采用SCAD函数获得近端算子来求解该模型。仿真和实际数据结果验证了该方法的有效性。结果表明,该文方法在低采样率以及低信噪比条件下成像结果优于OMP算法、SL0算法与BCS算法。
-
[1] 张晓茜, 徐勇军. 面向零功耗物联网的反向散射通信综述[J]. 通信学报, 2022, 43(11): 199–212. doi: 10.11959/j.issn.1000-436x.2022199ZHANG Xiaoxi and XU Yongjun. Survey on backscatter communication for zero-power IoT[J]. Journal on Communications, 2022, 43(11): 199–212. doi: 10.11959/j.issn.1000-436x.2022199 [2] XING Chengwen, JING Yindi, WANG Shuai, et al. New viewpoint and algorithms for water-filling solutions in wireless communications[J]. IEEE Transactions on Signal Processing, 2020, 68: 1618–1634. doi: 10.1109/TSP.2020.2973488 [3] CAO Yinwen, YU Song, SHEN Jing, et al. Frequency estimation for optical coherent MPSK system without removing modulated data phase[J]. IEEE Photonics Technology Letters, 2010, 22(10): 691–693. doi: 10.1109/LPT.2010.2044170 [4] LEVEN A, KANEDA N, KOC U V, et al. Frequency estimation in intradyne reception[J]. IEEE Photonics Technology Letters, 2007, 19(6): 366–368. doi: 10.1109/LPT.2007.891893 [5] FATADIN I and SAVORY S J. Compensation of frequency offset for 16-QAM optical coherent systems using QPSK partitioning[J]. IEEE Photonics Technology Letters, 2011, 23(17): 1246–1248. doi: 10.1109/LPT.2011.2158994 [6] HUANG Dezhao, CHENG T H, and YU Changyuan. Accurate two-stage frequency offset estimation for coherent optical systems[J]. IEEE Photonics Technology Letters, 2013, 25(2): 179–182. doi: 10.1109/LPT.2012.2232288 [7] YANG Tao, SHI Chen, CHEN Xue, et al. Hardware-efficient multi-format frequency offset estimation for M-QAM coherent optical receivers[J]. IEEE Photonics Technology Letters, 2018, 30(18): 1605–1608. doi: 10.1109/LPT.2018.2863739 [8] KIM J W, LEE Y S, JIN M Y, et al. Carrier frequency offset estimation for OFDM system with large oscillator phase noise[C]. 2021 International Conference on Information and Communication Technology Convergence (ICTC), Jeju Island, Korea, 2021: 368–370. [9] ZHENG Shuai, CHEN Jian, KUO Yonghong, et al. Statistical histogram-based blind frequency offset estimation for MPSK and MQAM signal[C]. 2022 International Conference on Machine Learning and Knowledge Engineering (MLKE), Guilin, China, 2022: 125–129. [10] JAYAPRAKASH A and REDDY G R. Covariance-fitting-based blind carrier frequency offset estimation method for OFDM systems[J]. IEEE Transactions on Vehicular Technology, 2016, 65(12): 10101–10105. doi: 10.1109/TVT.2016.2542181 [11] PAGE E S. A test for a change in a parameter occurring at an unknown point[J]. Biometrika, 1955, 42(3/4): 523–527. doi: 10.2307/2333401 [12] KESHAVARZ H, SCOTT C, and NGUYEN X. Optimal change point detection in Gaussian processes[J]. Journal of Statistical Planning and Inference, 2018, 193: 151–178. doi: 10.1016/j.jspi.2017.09.003 [13] KO S I M, CHONG TT L, and GHOSH P. Dirichlet process hidden Markov multiple change-point model[J]. Bayesian Analysis, 2015, 10(2): 275–296. doi: 10.1214/14-BA910 [14] KOKOSZKA P and LEIPUS R. Change-point in the mean of dependent observations[J]. Statistics & Probability Letters, 1998, 40(4): 385–393. doi: 10.1016/S0167-7152(98)00145-X [15] NA O, LEE Y, and LEE S. Monitoring parameter change in time series models[J]. Statistical Methods & Applications, 2011, 20(2): 171–199. doi: 10.1007/s10260-011-0162-3 [16] LING Jin, LI Xiaoqin, YANG Wenzhi, et al. The CUSUM statistic of change point under NA sequences[J]. Applied Mathematics-A Journal of Chinese Universities, 2021, 36(4): 512–520. doi: 10.1007/s11766-021-4015-z [17] YU Yuncai and CHEN Zhicheng. Strong convergence rates of multiple change-point estimator for ρ-mixing sequence[J]. Communications in Statistics - Theory and Methods, 2023, 52(13): 4605–4621. doi: 10.1080/03610926.2021.1998532 [18] KIM K, PARK J H, LEE M, et al. Unsupervised change point detection and trend prediction for financial time-series using a new CUSUM-based approach[J]. IEEE Access, 2022, 10: 34690–34705. doi: 10.1109/ACCESS.2022.3162399 期刊类型引用(0)
其他类型引用(1)
-