CN102288968A - 小孔径高分辨相控阵超声探测成像方法 - Google Patents

Abstract

本发明公开了一种小孔径高分辨相控阵超声探测成像方法，主要解决目前小孔径超声探测横向分辨率低的问题，其实现步骤为：根据设定的分辨率对探测区域进行空间离散化，并对离散化的探测区域进行扇形扫描；将得到的回波数据和扫描矩阵排列成回波列向量和操作实现矩阵；对回波列向量和操作实现矩阵进行随机抽取，得到观测向量和观测矩阵，并由观测向量和观测矩阵构造基于CS的高分辨超声探测成像模型；利用迭代加权最小二乘法求解出该模型重构出场景目标散射点系数；对该系数取模值并排列成一个二维矩阵即为小孔径高分辨超声探测成像。本发明相对传统的相控阵成像方法，成像的分辨率显著提高；适用于医用超声内窥镜、细小窄内壁材料探测。

Description

小孔径高分辨相控阵超声探测成像方法

技术领域

本发明属于超声探测领域，是一种小孔径高分辨相控阵超声探测成像方法，适用于医用超声内窥镜、细小窄内壁材料探测以及在安装探头空间有限的潜艇和鱼雷等场合，以获取高分辨的探测成像，从而为决策者提供有效的信息。

技术背景

众所周知，超声成像因无创伤，经济且安全已被广泛应用于医疗诊断、工业探伤、海洋水下探测等。从传统探测测量观点而言，超声探测中的辐射波的波长犹如一把尺子，它是度量单位，用于度量物体的空间尺度信息。传统的超声成像中认为横向分辨率取决于换能器阵列形成的孔径。换能器阵列孔径大则方位分辨率高。但从实际系统来分析，孔径的增大是有限的。大孔径不仅使得换能器尺寸变大，使得在很多场合，如在内壁狭小的管道探测成像难以应用，而且造成设备制造复杂成本增加。人们期望使用小孔径探头进行探测成像，但是由于小孔径探头不仅探测成像的横向分辨率低，而且还难以满足应用的要求。为了提高超声探测成像的横向分辨率，目前人们提出了如下成像方法：

1.频域法：该方法实际上是在频域内解决图像内插问题。1964年，Harris的“Diffraction and Resolving Power”文章奠定了超分辨率重建的数学基础，他首次将带限信号外推的方法应用于提高成像分辨率，但是这类方法的理论前提过于理想化，不能应用于小孔径超声探测成像。

2.空域法：利用图像在空域的先验特征为约束，如最小方差，全局和局部运动、空间可变模糊点扩散函数、非理想亚采样等建立观测模型，利用优化方法提高成像分辨率，如非均匀空域样本内插、迭代反投影、凸集投影、统计复原等方法，这类方法对于在换能器阵列孔径很小，表征对象信号发生混叠时，均无法取得横向高分辨率的探测成像。

上述两种高分辨率超声探测成像方法，均是针对于常规孔径或大孔径的。对于小孔径换能器阵列探测成像，由于在信号获取阶段，信息发生了混叠，上述两种的方法均不能从混叠信号中准确的重构原始信号。

发明内容

本发明的目的在于针对上述已有的技术问题，提出一种小孔径高分辨率相控阵超声探测成像方法，以从目标的混叠信号中，重构出高分辨率的影像。

实现本发明的技术方案是：将近来兴起的压缩感知CS理论应用于超声探测成像中，该压缩感知理论，是信号处理领域的重大成果，其核心思想是可以从低于奈奎斯特采样频率的有限采样信号中重构信号，而低于奈奎斯特采样频率的信号其频域会发生混叠，也就是说CS可以从混叠的采样数据中恢复原始信号，这就为解决小孔径相控阵超声探测成像横向分辨率低的难题带来了新机遇。根据这一理论，本发明首先根据设定的分辨率对探测区域进行空间离散化，然后用小孔径换能器阵列对离散化的探测区域进行扇形扫描，得到回波数据和扫描矩阵，将回波数据排列成一个列向量，同时将对应的扫描矩阵排列成一个操作实现矩阵，由此得到回波数据，操作实现矩阵与探测区域原始信息三者之间的关系式：y＝Ax，其中y是由回波数据排列成的回波列向量，A是操作实现矩阵，x为探测区域数据按行排列成的一个信号列向量，对回波数据和实现操作矩阵进行随机观测，得到观测向量b和观测矩阵Ψ；接着利用迭代加权最小二乘法求解出x，最后将x重新排列成一幅高分辨率图像，其具体实现步骤包括如下：

(1)根据超声波的波长λ，设定探测区域横向采样间隔为0.1λ和纵向采样间隔为0.2λ，并对探测区域进行离散化，得到离散化的二维探测区域；

(2)对二维离散化的探测区域进行相控阵超声扇形扫描，分别得到各个不同角度的探测目标的回波数据

同时保存各个探测角度的声束指向性函数矩阵

根据指向性函数与发射信号之间的关系计算得到各个角度的扫描矩阵

其中θ₁，θ₂，L，θ_P为相应的扫描角度；

(3)将各个角度的回波数据

排成一个回波列向量y，同时将对应的扫描矩阵

排成一个操作实现矩阵A，将离散化的二维探测区域按照行优先的原则排列成的一个信号列向量x，得到三者之间的关系为：y＝Ax；

(4)用随机矩阵Φ分别与列向量y和操作实现矩阵A进行相乘，实现对回波数据和扫描实现矩阵的随机抽取，得到观测向量b和观测矩阵Ψ，由观测向量b和观测矩阵Ψ构建基于CS的小孔径相控阵高分辨超声探测成像的数学模型为：

s.t b＝Ψx.

其中min表示最小化，

0≤p≤1，i表示信号列向量x中的第i个元素，x_i表示信号列向量x中第i个元素的值，b＝Φy，Ψ＝ΦA；

(5)按如下迭代加权最小二乘的方法对上述数学模型进行求解，得到信号列向量x：

(5a)初始化：

k＝0，e_k，e_min＝10^-8

其中x_k表示当前迭代得到的系数向量，k表示第k次迭代，

表示Ψ的伪逆，e_k是一个可调参数，初始设置在0～1之间，e_min表示e_k在迭代中的最小值；

(5b)由x_k构造系数向量的权值矩阵：

0≤p≤1，并更新观测矩阵Ψ_k：

Ψ_k+1＝ΨW_k+1

(5c)由Ψ_k+1更新系数向量x_k，得到x_k+1：

$x_{k+1}=W_{k+1}{\mathrm{\Ψ}}_{k+1}^{+}b$

(5d)判断如下条件是否成立：

${\left|\right|x_{k+1}-x_k\left|\right|}_2<\sqrt{{\mathrm{\&epsiv;}}_k}/100$

其中

i表示系数向量x_k+1和x_k中的第i个元素，x_k+1[i]表示系数向量x_k+1的第i个元素的值，x_k[i]表示系数向量x_k的第i个元素的值，N表示系数向量x_k+1的长度；

如果成立，判断条件ε_k＜ε_min是否满足，如果满足，则停止迭代，得到信号列向量x＝x_k+1；如果不满足则更新ε_k＝ε_k/10，返回步骤(5b)继续迭代；

如果不成立，返回步骤(5b)继续迭代；直到满足条件ε_k＜ε_min，停止迭代，返回得到的信号列向量x＝x_k+1：

(6)对上面得到的信号列向量x取模值，并排列成一个二维矩阵，即为得到的高分辨率探测战像。

本发明具有如下优点：

A.相比于传统的相控阵超声探测成像，本发明由于采用基于压缩感知的非凸函数最小化模型和迭代加权最小二乘的工程优化方法，所以能够得到高分辨率的小孔径相控阵超声探测成像。

B.本发明从导致信号混叠的原理入手，建立基于混叠信号的高分辨数学模型，提高了小孔径超声探测成像的分辨率。

附图说明

图1为本发明的总流程图；

图2为本发明对探测区域离散化的示意图；

图3为本发明对离散化探测区域进行扇形扫描示意图；

图4为本发明使用的换能器阵列与探测区域位置示意图；

图5为本发明将回波数据排列成列向量的示意图；

图6为本发明将扫描矩阵排列成操作实现矩阵的示意图；

图7为本发明将离散化探测区域排列成一个列向量的示意图；

图8为本发明中的迭代加权最小二乘算法子流程图；

图9为本发明将一个列向量排列成一个矩阵的示意图；

图10为本发明仿真使用的7个理想点目标的二维场景图；

图11为本发明小孔径相控阵超声探测高分辨率成像与传统成像对比图。

具体实施方式

参见图1，本发明的具体实现步骤如下：

步骤一，对探测区域进行离散化。

参见图2，本发明对探测区域进行离散化是根据超声波的波长λ，设定探测区域横向采样间隔为0.1λ和纵向采样间隔为0.2λ，并对探测区域进行等间隔采样，得到离散化的二维探测区域L×M，其中L为横向采样点数，M为纵向采样点数，L，M的取值与探测区域大小有关。

步骤二，对二维离散化的探测区域进行相控阵超声扇形扫描，得到回波数据

和扫描矩阵

参见图3，本发明的具体实现如下：

(2.1)根据图4所示的换能器阵列和探测区域的相对位置建立直角坐标系x-y，其中x表示横坐标，y表示纵坐标，

(2.2)在直角坐标系中，将换能器阵列固定在纵坐标为零的位置，即y＝0，且将阵列中心与探测区域中心对齐，换能器阵列的长度为Nd，则第l个阵元的横坐标x_l为：

$x_l=(l-\frac{N-1}{2})d$

其中N为换能器个数，d为相邻两换能器阵元之间的间隔，l＝0，1，L，N-1；

(2.3)当超声声束主瓣偏转角为θ_p时，计算位于探测区域角度θ_k，辐射距离r_k处的第k个目标散射点P_k(x_k，y_k)接收到的声场强度D(r_k，θ_k)：

$D(r_k,{\mathrm{\&theta;}}_k)=\frac{\exp (-j\frac{2\mathrm{\&pi;}}{\mathrm{\&lambda;}}{r}_k)}{r_k}p(t-\frac{r_k}{c})\exp (-j\frac{(N-1)}{\mathrm{\&lambda;}}\mathrm{\&pi;}d\sin {\mathrm{\&theta;}}_k)\underset{l=0}{\overset{N-1}{\mathrm{\&Sigma;}}}\left|a_l\right|\exp \left(\mathrm{jl}\frac{2\mathrm{\&pi;}}{\mathrm{\&lambda;}}d({\sin \mathrm{\&theta;}}_k-{\sin \mathrm{\&theta;}}_p)\right)$

其中λ为超声波的波长，c为声速，p(t)为激励信号的包络，该包络是缓变的，)是激励信号包络的延时，|a_l|表示第l个换能器阵元的复加权系数的模值，l表示第l个换能器阵元，-80°≤θ_p≤80°，-60°≤θ_k≤60°；

(2.4)当超声声束主瓣偏转角度为θ_p时，计算位于探测区域角度θ_k，辐射距离r_k处的第k个目标散射点P_k(x_k，y_k)反射的声场强度

其中，

表示换能器阵列的声束指向性函数；

(2.5)当超声声束主瓣偏转角为θ_p时，计算回波信号

$y_{{\mathrm{\&theta;}}_p}\left(t\right)=\underset{k=1}{\overset{\mathrm{LM}}{\mathrm{\&Sigma;}}}{a}_{{\mathrm{\&theta;}}_p}(t,k)\&CenterDot;x\left(k\right)=a_{{\mathrm{\&theta;}}_p}\left(t\right)x$

其中k表示第个k目标散射点，LM表示目标散射点的总个数，

为扫描矩阵函数，

表示扫描矩阵的行向量，

为待恢复的信号，γ_k表示第k个目标散射点的强度，x＝[x(1)，x(2)，L，x(LM)]^T表示待恢复的信号列向量，T表示转制；

(2.6)对回波信号

进行等间隔采样，采样间隔为T_s，采样点数为K，将回波信号写成矩阵与向量乘积的形式：

$y_{{\mathrm{\&theta;}}_p}=A_{{\mathrm{\&theta;}}_p}x$

其中

为离散回波列向量，t_k＝kT_s，k＝1，2，L，K，表示第k个采样时刻，为扫描矩阵，且满足：

$A_{{\mathrm{\&theta;}}_p}=\left[\begin{array}{cccc}{a}_{{\mathrm{\&theta;}}_p}(t_1,1)& a_{{\mathrm{\&theta;}}_p}(t_1,2)& L& a_{{\mathrm{\&theta;}}_p}(t_1,\mathrm{LM})\\ a_{{\mathrm{\&theta;}}_p}(t_2,1)& a_{{\mathrm{\&theta;}}_p}(t_2,2)& L& a_{{\mathrm{\&theta;}}_p}(t_2,\mathrm{LM})\\ M& M& O& M\\ a_{{\mathrm{\&theta;}}_p}(t_K,1)& a_{{\mathrm{\&theta;}}_p}(t_K,2)& L& a_{{\mathrm{\&theta;}}_p}(t_K,\mathrm{LM})\end{array}\right].$

步骤三，排列回波数据和扫描矩阵，构造回波数据的矩阵形式。

(3.1)将各个角度的回波数据

按照扫描角度的顺序依次从上到下排列成一个回波列向量y，如图5所示；

(3.2)将各个角度的扫描矩阵

按照扫描角度的顺序依次从上到下排列成操作实现矩阵A，如图6所示；

(3.3)将二维离散化的探测区域按照行优先的原则，依次从上到下排列成一个信号列向量x，如图7所示；

(3.4)回波列向量y由操作实现矩阵A与信号列向量x的乘积得到，其关系式为：y＝Ax。

步骤四，构建小孔径相控阵高分辨超声探测成像的数学模型。

(4.1)用随机矩阵分别与回波列向量y和操作实现矩阵A进行相乘，得到观测向量b＝Φy和观测矩阵Ψ＝ΦA，其中

表示实数域，Q表示随机矩阵Φ的行数，其取值大小与回波列向量y的长度有关，b∈￡^Q×1，Ψ∈￡^Q×LM，LM表示信号列向量x的长度；

(4.2)由观测向量b和观测矩阵Ψ构建基于压缩感知CS的小孔径相控阵高分辨超声探测成像的数学模型为：

s.t b＝Ψx。

步骤五，用加权迭代最小二乘的方法对小孔径相控阵高分辨超声探测成像的数学模型进行求解，得到信号列向量x。

参见图8，本步骤的具体实现如下：

(5a)初始化：

k＝0，e_k，e_min＝10^-8

其中x_k表示当前迭代得到的系数向量，k表示第k次迭代，

表示Ψ的伪逆，e_k是一个可调参数，初始设置在0～1之间，e_min表示e_k在迭代中的最小值；

(5b)由x_k构造系数向量的权值矩阵：

0≤p≤1，并更新观测矩阵Ψ_k：

Ψ_k+1＝ΨW_k+1

(5c)由Ψ_k+1更新系数向量x_k，得到x_k+1：

$x_{k+1}=W_{k+1}{\mathrm{\&Psi;}}_{k+1}^{+}b$

(5d)判断如下条件是否成立：

${\left|\right|x_{k+1}-x_k\left|\right|}_2<\sqrt{{\mathrm{\&epsiv;}}_k}/100$

其中i表示系数向量x_k+1和x_k中的第i个元素，x_k+1[i]表示系数向量x_k+1的第i个元素的值，x_k[i]表示系数向量x_k的第i个元素的值，N表示系数向量x_k+1的长度；

如果成立，则判断是否满足条件ε_k＜ε_min，如果满足，则停止迭代，得到信号列向量x＝x_k+1；如果不满足则更新ε_k＝ε_k/10，返回步骤(5b)继续迭代；

如果不成立，返回步骤(5b)继续迭代；直到满足条件ε_k＜ε_min，停止迭代，返回得到的信号列向量x＝x_k+1。

步骤六，将得到的信号列向量x取模值，排列成高分辨率探测成像。

参见图9，信号列向量x的长度为LM，从信号列向量x的最上端开始，依次从上到下，每次取M个数据排成一行，然后将各行按照从上到下的顺序依次排列成一个二维矩阵，矩阵的大小为L×M，该矩阵即为得到的高分辨率探测成像。

本发明的效果通过以下仿真进一步说明：

1.仿真条件：

运行系统为Intel(R)Core(TM)2Duo CPU E45002.20GHz 2.20GHz的32位Windows操作系统，仿真程序采用MATLAB程序设计语言实现。仿真一个有7个换能器阵元的线性相控阵阵列，相邻阵元间隔为150um，激励信号为单频包络信号，信号的中心频率5MHz，包络形式为高斯窗函数，持续时间为1.2μs；

对于传统的相控阵超声探测成像，扇形扫描角度从负30°到30°，每隔1°扫描一次，回波数据的采样频率为5MHz；对于本发明的方法，扇形扫描角度从负20°到20°，每隔4°扫描一次，回波数据的采样频率为1MHz；p取1，超声波在介质中传播速度设定为1500m/s；

2.仿真内容：

参照图10，构造有7个理想点目标的二维探测场景，理想点目标分别位于探测区域中不同位置；构造的回波信号带宽为1.67MHz；针对构造的场景分别用传统的相控阵扇形扫描的方法和本发明的方法进行探测成像，仿真结果如图11，其中图11(a)为传统相控阵超声探测成像结果图，图11(b)为本发明压缩感知CS方法探测结果图。

从图11可见，在换能器阵列孔径小时，利用传统相控阵超声探测方法成像的结果完全模糊，分辨率极低；而本发明方法将所有的理想点目标位置都精确重构出来了，且没有模糊；本发明方法的分辨率比传统的方法有了明显的提升。

Claims (2)

1.一种小孔径高分辨相控阵超声探测成像方法，包括如下步骤：

(1)根据超声波的波长λ，设定探测区域横向采样间隔为0.1λ和纵向采样间隔为0.2λ，并对探测区域进行离散化，得到离散化的二维探测区域；

(2)对二维离散化的探测区域进行相控阵超声扇形扫描，分别得到各个不同角度的探测目标的回波数据

同时保存各个探测角度的声束指向性函数

根据指向性函数与发射信号之间的关系计算得到各个角度的扫描矩阵

其中θ₁，θ₂，L，θ_P为相应的扫描角度；

(3)将各个角度的回波数据排成一个回波列向量y，同时将对应的扫描矩阵

排成一个操作实现矩阵A，将离散化的二维探测区域按照行优先的原则排列成的一个信号列向量x，得到三者之间的关系为：y＝Ax；

(4)用随机矩阵Φ分别与列向量y和操作实现矩阵A进行相乘，实现对回波数据和扫描实现矩阵的随机抽取，得到观测向量b和观测矩阵Ψ，由观测向量b和观测矩阵Ψ构建基于CS的小孔径相控阵高分辨超声探测成像的数学模型为：

s.t b＝Ψx.

其中min表示最小化，

0≤p≤1，i表示信号列向量x中的第i个元素，x_i表示信号列向量x中第i个元素的值，b＝Φy，Ψ＝ΦA；

(5)按如下迭代加权最小二乘的方法对上述数学模型进行求解，得到信号列向量x：

(5a)初始化：

k＝0，e_k，e_min＝10^-8

其中x_k表示当前迭代得到的系数向量，k表示第k次迭代，

表示Ψ的伪逆，e_k是一个可调参数，初始设置在0～1之间，e_min表示e_k在迭代中的最小值；

(5b)由x_k构造系数向量的权值矩阵：

0≤p≤1，并更新观测矩阵Ψ_k：

Ψ_k+1＝ΨW_k+1

(5c)由Ψ_k+1更新系数向量x_k，得到x_k+1：

$x_{k+1}=W_{k+1}{\mathrm{\&Psi;}}_{k+1}^{+}b$

(5d)判断如下条件是否成立：

${\left|\right|x_{k+1}-x_k\left|\right|}_2<\sqrt{{\mathrm{\&epsiv;}}_k}/100$

其中

i表示系数向量x_k+1和x_k中的第i个元素，x_k+1[i]表示系数向量x_k+1的第i个元素的值，x_k[i]表示系数向量x_k的第i个元素的值，N表示系数向量x_k+1的长度；

如果成立，判断条件ε_k＜ε_min是否满足，如果满足，则停止迭代，得到列向量x＝x_k+1；如果不满足则更新ε_k＝ε_k/10，返回步骤(5b)继续迭代；

如果不成立，返回步骤(5b)继续迭代；直到满足条件ε_k＜ε_min，停止迭代，返回得到的信号列向量x＝x_k+1

(6)对上面得到的信号列向量x取模值，并排列成一个二维矩阵，即为得到的高分辨率探测成像。

2.根据权利要求书1所述的小孔径高分辨相控阵超声探测成像方法，其特征在于，步骤2所述的根据指向性函数与发射信号之间的关系，计算各个角度的扫描矩阵

是通过如下公式计算：

$A_{{\mathrm{\&theta;}}_q}=D_{{\mathrm{\&theta;}}_q}\&CenterDot;p(t-{\mathrm{\&tau;}}_k)$

其中

表示扫描角度为θ_q时的声束指向性函数，p(t-τ_k)表示发射信号的包络延时，τ_k表示探测区域中第k个点目标的回波延时，θ_q＝θ₁，θ₂，L，θ_P为扫描角度，P为扫描角度的个数。

Images