CN104155030B - 一种考虑声线弯曲的声学ct温度场重建方法 - Google Patents

Abstract

本发明公开了一种考虑声线弯曲的声学CT温度场重建方法。将若干声波收发器布置在被测层面的周围，形成m条有效穿过被测层面的声波路径，并将该层面均匀地划分成n个像素。测量声波在各路径上的传播时间构成向量t。由预先算好的不考虑声线弯曲的重建矩阵l和向量t，重建出各像素的几何中心点处声速。依此，再用三角形前向展开法追踪出各路径的声线轨迹。用图形裁剪算法求出各声线轨迹在n个像素内的长度，形成考虑声线弯曲的重建矩阵l′。再由矩阵l′和向量t，为每个像素重建一声速。利用声速与温度的关系和插值算法，实现被测层面考虑声线弯曲地、准确快速温度场重建。该方法可显著提高非均匀温度场重建精度，且具有易编程实现、速度快等优点。

Description

一种考虑声线弯曲的声学CT温度场重建方法

技术领域

本发明属于声学CT温度场重建技术领域，通过考虑声波在非均匀温度场中会发生的声线弯曲现象，提高非均匀温度场声学CT重建精度。

背景技术

声学法测温是根据声波在介质中的传播速度间接地得到介质温度。气体介质中的声速c(m/s)与该气体介质的绝对温度T(K)、由气体组成决定的气体声音常数z间的关系可表示为对于烟道混合气而言该值为19.08，对于空气而言为20.05。

声学CT温度场重建技术需要在被测区域周围布置多个声波收发器。任一发射器发射的声波，可被所有接收器接收。这些收发器之间所形成的声波路径应尽可能均匀地覆盖被测区域。测量声波在各有效路径上的传播时间，便可在已知收发器位置的前提下，用适当的重建算法推算出被测区域的温度分布。

声学CT温度场重建技术具有非接触不干扰被测温场、测温范围广(零度～2000℃)、测量对象空间范围大(可达数十米)等优点。目前的声学CT温度场重建方法普遍将声波在发射/接收器间的传播轨迹(本征声线)按直线处理。但声波在非均匀温度场中传播时，其轨迹会因声波折射而产生弯曲。考虑声线弯曲现象，能有效地提高温度场的重建精度。

有学者采用费马原理与数学变分原理相结合的方法获取声线轨迹，实现了考虑声线弯曲的温度场重建。但费马原理结合数学变分原理的声线追踪法需要求解复杂的偏微分方程组，计算复杂、运算速度慢。

发明内容

1、发明目的

本发明提出了一种快速地、考虑声线弯曲的声学CT温度场重建方法，其目的在于提高非均匀温度场的声学CT重建精度。

2、技术方案：

本发明是通过以下技术方案来实现的：

一种考虑声线弯曲的声学CT温度场重建方法，其特征在于该方法由以下步骤构成：

步骤(1)：预先将若干个声波收发器布置在被测层面的周围，形成m条有效穿过被测层面的声波路径，并将被测层面均匀地划分成n个像素，n＜m。收发器的位置应使收发器间所形成的声波路径尽可能均匀地穿过被测区域。测量声波在m条有效路径上的传播时间并将它们组合成m维声波传播时间向量t。由预先计算好的不考虑声线弯曲的重建矩阵l和向量t，用式(1)求出每个像素的几何中心点处声速的倒数，进而求出各像素几何中心点处的声速。

其中s₁,...,s_n分别表示n个像素几何中心点处声速的倒数；t₁,...,t_m分别表示声波在m条路径上的传播时间；l_ki表示第k条路径在第i个像素内的长度，k∈{1,2,...,m},i∈{1,2,...n}。由于声波收发器位置已知，所以被测区域像素划分后，各声波路径视为直线即不考虑声线弯曲时的l_ki可求，不考虑声线弯曲的重建矩阵l可预先算好。

步骤(2)：用三角形前向展开法追踪出步骤(1)给出的重建场中各有效声波路径的声线轨迹。当已知入射点声速及声线出射角时，用三角形前向展开法可在给定等腰三角形底和高且各三角形顶点声速可获得的前提下，求出声线轨迹；下面只给出入射点声速、声线出射角、等腰三角形底和高及各三角形顶点声速的确定方法。

入射点声速：用布置在各声波收发器处的接触式温度传感器，获得声波收发器处的温度；利用声速与温度的关系求出各声波收发器处的声速，即入射点声速。

声线出射角：以声波收发器A、B间的声线轨迹追踪为例，简述如下。①以声线起点A到终点B的直线方向所对应的角度θ₁为声线出射角初值，用三角形前向展开法获得声线与终点所在直线(与起点到终点的连线相垂直的直线)的交点B1。若B1与B的间距达到误差限要求，则θ₁即为所求本征声线出射角；否则调整出射角为θ₂使其对应声线在终点所在直线上的交点B2位于B的另一侧(即B2、B1分列于B的两侧)，如图1所示。②将区间[min(θ₁,θ₂),max(θ₁,θ₂)]p等分，得到p个出射角估值min[θ₁,θ₂]+j×|θ₂-θ₁|/p，j＝1,2,...,p-1，其中p为正整数且其取值保证等分间隔≤π/100。③用三角形前向展开法求出上述p个出射角估值所对应的声线在终点所在直线上的p个交点，找出其中与终点B最邻近的两交点Bn1和Bn2，则所求本征声线出射角为θ_si＝θ_n1+(θ_n2-θ_n1)×(B-Bn1)/(Bn2-Bn1)。其中θ_n1,θ_n2分别是对应于交点Bn1和Bn2的声线出射角。

等腰三角形底和高：追踪两收发器间的本征声线轨迹时，展开的等腰三角形的高设置为两收发器间距的3％，且三角形的底边与高相等。

各三角形顶点声速：将被测区域均匀地细化为20×20＝400个像素，以步骤(1)获得的n个像素中心点的声速及各收发器处声速为已知，用插值的方法获得400个细化像素几何中心点处的声速。再根据三角形展开过程中三角形顶点所在细化像素4个顶点的声速，用双线性插值算法求出三角形顶点声速。

步骤(3)：用图形裁剪算法求出步骤(2)所获得的各有效声波路径的声线轨迹在n个像素内的声线长度，形成考虑声线弯曲的重建矩阵l′如式(2)所示。再由矩阵l′和声波传播时间向量t，用式(3)为每个像素的几何中心点求出一声速的倒数。然后，利用声速与温度的关系为每个像素的几何中心点求出一考虑声线弯曲的温度。再用插值运算的方法，获得整个被测层面的温度分布，实现被测层面考虑声线弯曲地、准确快速声学CT温度场重建。

s′＝(l′^Tl′)^-1l′^Tt， (3)

其中l′_ki表示由步骤(2)获得的第k条声波路径的声线轨迹在第i个像素内的长度，k∈{1,2,...,m},i∈{1,2,...n}，s′＝(s′₁,...,s′_n)^T，s′₁,...,s′_n分别表示n个像素几何中心点处声速的倒数。

3、优点及效果：

声学CT温度场重建技术可用于在大型工业炉内的温度场重建与检测中。并被尝试用于仓储粮食温度分布的非接触监测等新的领域。声波在非均匀温度场中传播时会发生折射，即声波不再沿直线传播。而目前的声学CT温度场重建方法普遍将声波在发射、接收器间的传播轨迹---本征声线，按直线处理，这必然会给测量带来误差。考虑声线弯曲现象，能有效地提高温度场的重建精度。目前虽有学者采用费马原理与数学变分原理相结合的方法获取声线轨迹，实现了考虑声线弯曲的温度场重建。但费马原理结合数学变分原理的声线追踪法需要求解复杂的偏微分方程组，计算复杂、运算速度慢。

本发明提出一种用三角形前向展开法追踪声线的、考虑声线弯曲的温度场重建方法。该方法可显著提高非均匀温度场重建精度，且具有易编程实现、速度快等优点。

无论是工业炉内的烟气温度场，还是仓储粮食中的温度分布，通常存在较大的温度(声速)梯度。粮食可视为多孔介质，声波可通过粮食颗粒间的空隙而传播。储粮中温度(声速)梯度的存在会导致储粮中宏观声线的弯曲，宏观声线的追踪可等效为气体温度场中的声线追踪。因此本发明提出的重建方法对提高工业炉烟气温度和储粮温度声学CT监测精度都具有重要意义，能更好地适应实际温度场重建需求，涵含良好的经济效益和社会效益。

附图说明

图1为本征声线出射角搜索示意图。

图2为将8个声波收发器布置在一个正方形被测层面周围，形成24条有效声波路径，并将被测层面划分成16个像素的示意图。图中给出的是不考虑声线弯曲的有效声波路径。

图3为将16个声波收发器布置在一个正方形被测层面周围，形成96条有效声波路径，并将被测层面划分成64个像素的示意图。图中给出的是不考虑声线弯曲的有效声波路径。

图4给出了采用图2所示声波收发器布局在某一偏置单峰温度场中追踪到的弯曲的声线轨迹。

图5给出了采用图3所示声波收发器布局在某一对称双峰温度场中追踪到的弯曲的声线轨迹。

图6给出前向展开的三角形和其中的声线轨迹示意图。

图7给出三角形内声波路径和声波传播时间计算示意图。

图8给出了追踪到的由S1至S3、S4、S5、S6、S7和S8的本征声线以及该温度场的等温线描述。

具体实施方式：

下面结合附图对本发明做进一步的说明：

本发明提出了一种考虑声线弯曲的声学CT温度场重建方法。所采用的温度场重建系统系统由声波收发器、机架、电缆、驱动及信号调理系统、数据采集系统和计算机组成。安装在机架上的声波收发器通过电缆连接到驱动与信号调理系统，信号调理系统输出经电缆连接到数据采集系统。任一发射器发射的声波，可被所有接收器接收。这些收发器之间所形成的声波路径应尽可能均匀地覆盖被测区域。测量出声波在各有效路径上的传播时间，便可在已知收发器位置的前提下，用本发明所述方法，重建出被测区域的温度分布。

本发明方法具体步骤如下：

一种考虑声线弯曲的声学CT温度场重建方法，其特征在于该方法由以下步骤构成：

步骤(1)：预先将若干个声波收发器布置在被测层面的周围，形成m条有效穿过被测层面的声波路径，并将被测层面均匀地划分成n个像素，n＜m。收发器的位置应使收发器间所形成的声波路径尽可能均匀地穿过被测区域。例如将8个声波收发器布置在一个正方形被测层面周围，形成24条(m＝24)有效声波路径，并将被测层面划分成16个像素(n＝16)，如图2所示。再如将16个声波收发器布置在一个正方形被测层面周围，形成96条(m＝96)有效声波路径，并将被测层面划分成64个像素(n＝64)，如图3所示。图2、图3中给出的是不考虑声线弯曲现象的有效声波路径，即直线路径。在一个检测周期内顺序启闭各声波发射器使其轮流发射声波，并保证每次最多一个发射器发射声波。任一发射器发射声波时，所有接收器均接收此声波并转换为电信号。这些电信号经信号调理器和数据采集卡后，同步(采集卡的输入通道数不小于收发器数目)或分组同步(采集卡的输入通道数大于收发器数目)进入计算机。采用某种时延估计算法，例如结合小波抑噪的互相关时延估计法，就可由这些数据，测量出声波在m条有效路径上的传播时间并将它们组合成m维声波传播时间向量t。由预先计算好的不考虑声线弯曲的重建矩阵l和向量t，用式(1)求出每个像素的几何中心点处声速的倒数，进而求出各像素几何中心点处的声速。

其中s₁,...,s_n分别表示n个像素几何中心点处声速的倒数；t₁,...,t_m分别表示声波在m条路径上的传播时间；l_ki表示第k条路径在第i个像素内的长度，k∈{1,2,...,m},i∈{1,2,...n}。由于声波收发器位置已知，所以被测区域像素划分后，各声波路径视为直线即不考虑声线弯曲时的l_ki可求，不考虑声线弯曲的重建矩阵l可预先组合好。

步骤(2)：用三角形前向展开法追踪出步骤(1)给出的重建场中各有效声波路径的声线轨迹。当已知入射点声速及声线出射角时，用三角形前向展开法可在给定等腰三角形底和高且各三角形顶点声速可获得的前提下，求出声线轨迹，具体计算过程详述如下。

设声线轨迹入射点坐标为(x_in,y_in)，描述声线出射角的入射方向矢量为(A_x,A_y)，以(x_in,y_in)为顶点作一底边长为d、高为h、方向与入射方向相同的等腰三角形。该三角形要足够小使其内部的声速分布可做线性近似。则声线在此三角形内的轨迹为始点为(x_in,y_in)、终点为(x_out,y_out)、半径为r的圆弧。如图6所示。

三角形内任意点的声速可表示为C(x,y)＝ax+by+c，由三角形三个顶点的坐标及三个顶点处的声速值C₁、C₂、C₃可以确定系数a、b、c。由Snell定律可确定该圆弧的轨迹方程为：

(x-x_o)²+(y-y_o)²＝r² (2)

其中圆半径r＝d·C₁/(C₂-C₃)，而圆心(x_o,y_o)为直线y-y_in＝-A_x/A_y·(x-x_in)与直线ax+by+c＝0的交点。圆弧曲线与三角形底边交点到垂足(x_h,y_h)的距离为：

$l=r-\sqrt{r^2-h^2}---\left(3\right)$

圆弧曲线与三角形底边的交点，即声线在该三角形区域的声线终点坐标为：

(xo_ut,yo_ut)＝(x_h+A_y·l,y_h-A_x·l)(4)

声线终点方向矢量为：

$\{\begin{array}{c}{B}_x=-h\·A_y/r+(1-l/r)\·A_x\\ B_y=h\·A_x/r+(1-l/r)\·A_y\end{array}---\left(5\right)$

将声线终点(x_out,y_out)和方向矢量(B_x,B_y)作为下次三角形展开的起点和起始方向，使声线追踪能够反复执行。声波传播的路径长度L和传播时间t可由各三角形内的路径长度L_i和传播时间t_i求和得到，即L＝∑L_i；t＝∑t_i。当三角形三个顶点处的声速C₁≠C₂≠C₃时，三角形内的声线轨迹是圆弧曲线如图7所示，三角形内的路径长度和传播时间分别为：

L_i＝r·(θ₂-θ₁) (6)

$t_i=\frac{1}{2\·\sqrt{a^2+b^2}}\left|log(\frac{1+{\mathrm{cos\θ}}_1}{1-{\mathrm{cos\θ}}_1}\·\frac{1-{\mathrm{cos\θ}}_2}{1+{\mathrm{cos\θ}}_2})\right|---\left(7\right)$

当C₂＝C₃时，声线沿着三角形顶角平分线前进，三角形内的路径长度和传播时间分别为

L_i＝h (8)

下面给出入射点声速、声线出射角、等腰三角形底和高及各三角形顶点声速的确定方法。

入射点声速：声波在气体介质中的声速c(m/s)与该气体介质的绝对温度T(K)、由气体组成决定的气体声音常数z间的关系可表示为对于烟道混合气而言该值为19.08，对于空气而言为20.05；用布置在各声波收发器处的接触式温度传感器，获得声波收发器处的温度；利用声速与温度的关系求出各声波收发器处的声速，即入射点声速。

声线出射角：以声波收发器A、B间的声线轨迹追踪为例，简述如下：

①以声线起点A到终点B的直线方向所对应的角度θ₁为声线出射角初值，用三角形前向展开法获得声线与终点所在直线(与起点到终点的连线相垂直的直线)的交点B1。若B1与B的间距达到误差限要求，则θ₁即为所求本征声线出射角；否则调整出射角为θ₂使其对应声线在终点所在直线上的交点B2位于B的另一侧(即B2、B1分列于B的两侧)，如图1所示。一般该误差限可设置为距离最近的两收发器间距的0.03％。

②将区间[min(θ₁,θ₂),max(θ₁,θ₂)]p等分，得到p个出射角估值min[θ₁,θ₂]+j×|θ₂-θ₁|/p，j＝1,2,...,p-1，其中p为正整数且其取值保证等分间隔≤π/100。

③用三角形前向展开法求出上述p个出射角估值所对应的声线在终点所在直线上的p个交点，找出其中与终点B最邻近的两交点Bn1和Bn2，则所求本征声线出射角为θ_si＝θ_n1+(θ_n2-θ_n1)×(B-Bn1)/(Bn2-Bn1)。其中θ_n1,θ_n2分别是对应于交点Bn1和Bn2的声线出射角。

等腰三角形底和高：追踪两收发器间的本征声线轨迹时，展开的等腰三角形的高设置为两收发器间距的3％，且三角形的底边与高相等。

各三角形顶点声速：将被测区域均匀地细化为20×20＝400个像素，以步骤(1)获得的n个像素中心点的声速及各收发器处声速为已知，用插值的方法，例如双调和样条插值法，获得400个细化像素几何中心点处的声速。再根据三角形展开过程中三角形顶点所在细化像素4个顶点的声速，用双线性插值算法求出三角形顶点声速。

步骤(3)：用图形裁剪算法求出步骤(2)所获得的各有效声波路径的声线轨迹在n个像素内的声线长度，形成考虑声线弯曲的重建矩阵l′如式(2)所示。再由矩阵l′和声波传播时间向量t，用式(3)为每个像素的几何中心点求出一声速的倒数。然后，利用声速与温度的关系，为每个像素的几何中心点求出一考虑声线弯曲的温度。再用插值运算的方法，例如三次样条算法，获得整个被测层面的温度分布，实现被测层面考虑声线弯曲地、准确快速声学CT温度场重建。

s′＝(l′^Tl′)^-1l′^Tt， (3)

其中l′_ki表示由步骤(2)获得的第k条声波路径的声线轨迹在第i个像素内的长度，k∈{1,2,...,m},i∈{1,2,...n}，s′＝(s′₁,...,s′_n)^T，s′₁,...,s′_n分别表示n个像素几何中心点处声速的倒数。

图6给出前向展开的三角形和其中的声线轨迹示意图。

图7给出三角形内声波路径和声波传播时间计算示意图。

实施例1：图4给出了采用图2所示声波收发器布局在某一单峰温度场中追踪到的声线轨迹。表1给出了由声波收发器S1至收发器S3、S4、S5、S6、S7和S8的6条本征声线出射角计算结果。图8给出了追踪到的由S1至S3、S4、S5、S6、S7和S8的本征声线以及该温度场的等温线描述。未考虑声线弯曲时重建温度场的最大相对误差为7.23％，均方根误差为3.4％；考虑声线弯曲后则分别为4.93％和2.27％。图8给出了追踪到的由S1至S3、S4、S5、S6、S7和S8的本征声线以及该温度场的等温线描述。

表1本征声线出射角计算结果

实施例2：图5给出了采用图3所示声波收发器布局在某一双峰温度场中追踪到的声线轨迹。未考虑声线弯曲时重建温度场的最大相对误差为14.43％，均方根误差为6.46％；考虑声线弯曲后则分别为10.98％和5.60％。

结论：本发明提出的这种考虑声线弯曲的声学CT温度场重建方法，可显著提高非均匀温度场重建精度，且具有易编程实现、速度快等优点。无论是工业炉内的烟气温度场，还是仓储粮食中的温度分布，通常存在较大的温度(声速)梯度。因此本发明提出的重建方法对提高工业炉烟气温度和储粮温度声学CT监测精度都具有重要意义，能更好地适应实际温度场重建需求，涵含良好的经济效益和社会效益。

Claims (2)

1.一种考虑声线弯曲的声学CT温度场重建方法，其特征在于：该方法由以下步骤构成：

步骤(1)：预先将若干个声波收发器布置在被测层面的周围，形成m条有效穿过被测层面的声波路径，并将被测层面均匀地划分成n个像素，n＜m；测量声波在m条有效路径上的传播时间并将它们组合成m维声波传播时间向量t；由预先计算好的不考虑声线弯曲的重建矩阵l和向量t，用式(1)求出每个像素的几何中心点处声速的倒数，进而求出各像素几何中心点处的声速及温度；

s＝(l^Tl)^-1l^Tt

其中s₁,...,s_n分别表示n个像素几何中心点处声速的倒数；t₁,...,t_m分别表示声波在m条路径上的传播时间；l_ki表示第k条路径在第i个像素内的长度，k∈{1,2,...,m},i∈{1,2,...n}；由于声波收发器位置已知，被测区域像素划分后，各声波路径视为直线即不考虑声线弯曲时的l_ki可求，不考虑声线弯曲的重建矩阵l能够预先算好；

步骤(2)：用三角形前向展开法追踪出步骤(1)给出的重建场中各有效声波路径的声线轨迹；当已知入射点声速及声线出射角时，用三角形前向展开法可在给定等腰三角形底和高且各三角形顶点声速可获得的前提下，求出声线轨迹；下面只给出入射点声速、声线出射角、等腰三角形底和高及各三角形顶点声速的确定方法；

入射点声速：用布置在各声波收发器处的接触式温度传感器，获得声波收发器处的温度；利用声速与温度的关系求出各声波收发器处的声速，即入射点声速；

声线出射角：以声波收发器A、B间的声线轨迹追踪为例，简述如下：

①以声线起点A到终点B的直线方向所对应的角度θ₁为声线出射角初值，用三角形前向展开法获得声线与终点所在直线的交点B1；若B1与B的间距达到误差限要求，则θ₁即为所求本征声线出射角；否则调整出射角为θ₂使其对应声线在终点所在直线上的交点B2位于B的另一侧，即B2、B1分列于B的两侧；

②将区间[min(θ₁,θ₂),max(θ₁,θ₂)]p等分，得到p个出射角估值min[θ₁,θ₂]+j×|θ₂-θ₁|/p，j＝1,2,...,p-1，其中p为正整数且其取值保证等分间隔≤π/100；

③用三角形前向展开法求出上述p个出射角估值所对应的声线在终点所在直线上的p个交点，找出其中与终点B最邻近的两交点Bn1和Bn2，则所求本征声线出射角为θ_si＝θ_n1+(θ_n2-θ_n1)×(B-Bn1)/(Bn2-Bn1)；其中θ_n1,θ_n2分别是对应于交点Bn1和Bn2的声线出射角；

等腰三角形底和高：追踪两收发器间的本征声线轨迹时，展开的等腰三角形的高设置为两收发器间距的3％，且三角形的底边与高相等；

各三角形顶点声速：将被测区域均匀地细化为20×20＝400个像素，以步骤(1)获得的n个像素中心点的声速及各收发器处声速为已知，用插值的方法获得400个细化像素几何中心点处的声速；再根据三角形展开过程中三角形顶点所在细化像素4个顶点的声速，用双线性插值算法求出三角形顶点声速；

步骤(3)：用图形裁剪算法求出步骤(2)所获得的各有效声波路径的声线轨迹在n个像素内的声线长度，形成考虑声线弯曲的重建矩阵l′如式(2)所示；再由矩阵l′和声波传播时间向量t，用式(3)为每个像素的几何中心点求出一声速的倒数；然后，利用声速与温度的关系为每个像素的几何中心点求出一考虑声线弯曲的温度；再用插值运算的方法，获得整个被测层面的温度分布，实现被测层面考虑声线弯曲地、准确快速声学CT温度场重建；

s′＝(l′^Tl′)^-1l′^Tt， (3)

其中l′_ki表示由步骤(2)获得的第k条声波路径的声线轨迹在第i个像素内的长度，k∈{1,2,...,m},i∈{1,2,...n}，s′＝(s′₁,...,s′_n)^T，s′₁,...,s′_n分别表示n个像素几何中心点处声速的倒数。

2.根据权利要求1所述的考虑声线弯曲的声学CT温度场重建方法，其特征在于：步骤(1)中收发器的位置应使收发器间所形成的声波路径尽可能均匀地穿过被测区域。