一种基于回波能量积分的气体超声波流量计信号处理方法和系统与流程

本发明涉及流量测量领域，是一种气体超声波流量计变送器，特别是一种基于回波能量积分的气体超声波流量计信号处理方法和系统。采用本发明提出的气体超声波流量计信号处理方法处理超声波回波信号，实现气体流量的测量。

背景技术：

气体超声波流量计由于其测量精度高、量程比大、无压损等优点，被广泛应用于天然气等大口径管道的气体流量测量领域。在测量方法上，因传播时间差法受温度变化引起的测量误差较小，是最常用、最有效的测量方法。采用传播时间差法测量流量时，首先需要驱动发射换能器产生超声波信号，然后再由接收换能器将得到的超声波信号转换为回波信号；通过测量发射和接收信号之间的时间差来计算出相应的顺流、逆流传播时间，获得管道横截面的面速度，进而计算气体流量。超声波信号的起始时刻是已知的，而回波信号到达时刻却是未知的，需要根据回波信号上某一稳定的特征点来确定。由于超声波信号在气体中传播时，其能量衰减严重，且随着气体流速的增加会产生相应的路径偏移，进一步削弱能量，使得回波信号幅值微弱，易受噪声干扰。而直接通过最大峰值点等简单方法确定的特征点，都会产生较大的波动，影响测量精度。

国内外学者对气体超声流量计信号处理方法进行了研究，主要有以下几种方法。

(1)基于固定阈值的方法

从回波信号幅值的角度出发，固定阈值的测量方法是通过设定一个固定的阈值，当超声波回波信号的幅值大于这个值时，就把这一个时刻作为超声波回波信号的到达时间，并以此求得超声波回波信号的传播时间。tamotsukobayashi等人(kobayashit.,shigetak.,fujiit.transit-timedifferencetypeultrasonicflowmeter[p].u.s.:20020062690a1,2002-05-30.)介绍了固定阈值方法的硬件实现，即利用比较器以及参考电压值来确定超声波回波信号的到达时间。r.bates等人(batesr.,etal.developmentofacustomon-lineultrasonicvapouranalyzer/flowmeterfortheatlasinnerdetector,withapplicationtogaseoustrackingandcherenkovdetectors[j].journalofinstrumentation,2013,8:c01002.)提到了上述固定阈值方法的软件实现，通过在超声波信号发出的同时启动时钟，并在回波信号幅值超过预设值之后停止时钟来计算超声波信号的传播时间。然而，噪声污染、漂移或者超声波换能器的使用年限等都会引起超声波回波信号幅值的变化，一旦超声波回波信号幅值波动较大，根据固定阈值法获得的回波信号接收时间就会出错。

(2)基于能量变化率的方法

从回波信号能量变化的角度出发，沈子文等人也提出了一种基于能量变化率的信号处理方法(沈子文,徐科军,方敏,等.基于能量变化率的气体超声波流量计信号处理方法[j].仪器仪表学报,2015,36(9):2138-2144.)，首先对采样的回波信号进行滤波等预处理，接着找到回波信号的各峰值点，并对其幅值进行平方，作为该点峰值点的能量(信号的能量与信号的幅值平方呈比例关系，fangm,xukj,zhuwj,etal.energytransfermodelanditsapplicationsofultrasonicgasflow-meterunderstaticanddynamicflowrates.[j].reviewofscientificinstruments,2016,87(1):015107.)，得到每个峰值点的能量变化率，并对其进行包络线拟合，为了减少计算量，保证信号处理的实时性，采用线性插值的方法求取超声回波信号能量的包络线。在得到各流量下的能量变化曲线后，根据能量变化曲线设定用于区分特征波的阈值。在实际测量时，根据预先设定的阈值和实时求取回波信号的能量变化曲线，判别对应的特征点。由于该方法需要在测量流量之前进行阈值参数选取，操作相对复杂。

(3)基于能量突变的方法

美国daniel公司采用检测能量突变的方法查找特征点的方法(williamfreund,winsorletton,jamesmc-clellan,baocangjia,anniwey,wenchang.methodandapparatusformeasuringthetimeofflightofasignal,uspatentno.5983730,nov.16,1999)。由于回波信号的能量在整个激励到传播再到接收的过程中，从弱到强，再从强到弱，因此回波能量变化率将先增大再减少，依据能量变化率的临界变化点，即可确定超声波传播时间。具体地说，首先求取回波信号各点的幅值平方；然后，利用滑动平均求取每一点的平均能量值，再绘制出前、后两点平均能量之比的变化曲线，即能量变化率曲线；最后，找到能量变化率的临界点，确定超声波传播时间。但是，该专利没有披露方法的关键技术，如滑动窗口的时间间隔选取，以及求取能量变化率时前、后两点之间的时间间隔等。

技术实现要素：

针对气体超声波流量计回波信号难以确定稳定的特征点，判断回波信号的到达时刻，而仪表又要求具有较高实时性的问题，本发明提出了基于回波能量积分的气体超声波流量计信号处理方法。因为信号的能量与信号的幅值平方呈比例关系，因此，为了简化计算，本发明使用回波信号幅值的平方表征回波信号的能量。在以dsp和fpga为双核心的硬件平台上(徐科军，方敏，汪伟，朱文姣，沈子文.一种基于fpga和dsp的气体超声波流量计，申请中国发明专利，公开号cn104697593a，申请日2015.06.10)，编写了软件程序，实时实现基于回波能量积分的信号处理方法。

基于回波能量积分的气体超声波流量计信号处理方法的计算步骤为：首先根据采样的超声波回波信号，通过比较得到回波能量信号的最大正峰值点；然后向前确定区间起始点，以大于回波信号归一化幅值的0.1倍为判断条件，找到积分起始点；从积分起始点(积分起始点也是一个峰值点)到最大正峰值点之间，每找到一个峰值点，将该点幅值平方后归一化，计算相邻峰值点间包围的梯形面积并累加；判断累加面积值和预设面积值的大小；当累加面积值超过预设面积值后停止积分累加，通过累加面积值减去预设面积值得到面积差值(小梯形面积)；已知面积差值(小梯形面积)和相邻峰值点的纵坐标(小梯形上下底边值)，求解小梯形上底边与x轴的交点(梯形的高)，作为特征点；确定超声波回波信号的顺流、逆流传播时间，并计算气体流量。

附图说明

图1是基于回波能量积分的信号处理方法流程图。

图2是0m³/h下的回波能量信号及其包络线和不同流量下的回波能量信号及其包络线图。

图3是确定积分区间起始点示意图。

图4是确定选取积分起始点示意图。

图5是基于回波能量积分方法的积分面积示意图。

图6是确定预设面积值原理图。

图7是基于回波能量积分的特征点计算方法。

图8是气体超声波流量计的硬件框图。

图9是气体超声流量计系统主控芯片dsp软件流程图。

具体实施方式

图1是基于回波能量积分的信号处理方法流程图。将回波信号进行滤波等处理后，便开始计算特征点。

(1)寻找回波信号最大正峰值点和最小负峰值点。在采集的回波信号数据中,通过逐点比较得到满足条件的两个峰值点。

(2)将两个峰值点平方后做比较。根据信号的能量与信号的幅值平方呈比例关系，将最大正峰值点和最小负峰值点分别平方后比较大小。若要完整地表示回波能量信号，则需要将采集的回波信号每一点的幅值都进行平方，这样处理会增加算法的运算量，影响程序的整体效率。本发明采用找到回波信号的上升段各个峰值点，对每个峰值点的幅值进行平方处理，既可以充分利用回波能量信号所携带的信息，又可以减少整个算法的计算量。并且，在寻找到一个峰值点后，便进行后续的面积积分计算，当达到预设条件(累加面积值大于预设面积值sa)后即可完成计算，进一步减少计算过程的重复循环。采用这样的流程减少了数据的运算量，有助于提高系统的实时性。

(3)得到回波能量信号最大正峰值点。为确定积分起始点做准备，回波能量信号最大正峰值点的幅值(纵坐标)用amax表示。

(4)确定积分起始点，幅值平方归一化。在已知最大正峰值点的情况下，首先找到积分的区间起始点，然后根据设定的要求(大于回波信号归一化幅值的0.1倍)，找到积分起始点，从积分起始点开始对面积进行积分。具体设定的要求在图3进行详细叙述。

(5)求积分面积。找到积分起始点后开始计算面积值，具体包括6个步骤：

求积分面积的第一步：寻找积分起始点后第一个峰值点。

求积分面积的第二步：对该峰值点幅值平方后归一化。该峰值点幅值用an表示，归一化后的幅值都用表示，该峰值点的回波能量信号幅值是其平方，用an²表示，n表示从积分起始点到最大正峰值点间的各峰值点下标。归一化公式为

求积分面积的第三步：计算相邻两个峰值点间梯形面积，并累加。如图5所示，相邻的两个峰值点都可以和横坐标轴围成一个小梯形，通过梯形面积累加得到回波信号的能量值。相邻两个峰值点(幅值归一化)坐标分别用和(xn+1,)表示，n和n+1表示从积分起始点到最大正峰值点间的各峰值点下标。小梯形面积sn计算公式为

将所有小梯形的面积进行累加，得到累加面积值。累加面积s的计算公式为

求积分面积的第四步：是否大于预设面积值sa？判断累加面积值与预设面积值的大小：

若累加面积值未达到预设面积值，则保存该峰值点平方后归一化的幅值(梯形下底)，向后寻找下一个峰值点，然后重复第二、三、四步骤，继续判断是否大于预设面积值sa；

若累加面积值达到了预设面积值，则开始进行下面一步。

求积分面积的第五步：计算面积差值s'。累加面积值减去预设面积值得到面积差值s'。

求积分面积的第6步：通过面积差值求特征点橫坐标。该步骤在图7中详细叙述。

(6)得到特征点的值。通过特征点的值则可以得到顺流、逆流传播时间，获得管道横截面的面速度，进而计算气体流量。以某一声道，顺流传播时间td为例，计算公式为：

td＝xt/f(4)式中，xt为计算得到的特征点横坐标值，f为回波信号的采样频率；

气体流量q的计算公式为：

式中，为声道上的平均速度，s为管道的横截面积，α表示一个与声道数量、声道分布及雷诺数等因素有关的修正系数，l表示声道长度，θ表示换能器与气体管道之间夹角，td、tu分别表示同一声道的超声波顺流、逆流传播时间。

图2是0m³/h下的回波能量信号及其包络线和不同流量下的回波能量信号及其包络线图。上半部分的图是0m³/h下的回波能量信号及其包络线，这是通过将50组0m³/h下的回波信号进行滤波、幅值平方、求平均和归一化等预处理，利用matlab中曲线拟合工具对回波能量信号的所有峰值点进行包络线拟合得到的。下半部分的图是不同流量下的回波能量信号及其包络线，将0m³/h、200m³/h、400m³/h、600m³/h、800m³/h、1000m³/h下的各50组回波信号进行同样的处理拟合得到的。通过不同流量回波能量信号包络线可以看出，当流量增大时，回波能量信号下降阶段会发生畸变，(有研究表明，当被测流体流速增大到一定范围时，回波信号会产生严重的畸变。sabatiniam.adigital-signal-processingtechniqueforultrasonicsignalmodelingandclassification[j].ieeetransactionsoninstrumentationmeasurement,2001,50(1):15-21.)，回波信号的幅值会经历从零开始快速增加，到达最大值后缓慢衰减至零的过程，即在释放能量的过程中，是缓慢且不平稳的释放。但是，在回波能量信号的上升阶段，由于受激励换能器发射的超声波直接作用，基本不存在反射声波等影响，能量集中且稳定，因此，可以看到，在图2中1、2、3、4、5、6所表示的各个流量下的上升阶段的波形轮廓基本不变，且随着流量的增大，上升阶段的波形轮廓发生平移，但各流量下波形轮廓基本保持平行。由于每次发射的超声波所携带的能量是相同的(激励信号的频率、幅值、相位都相同)，所以，回波能量信号上升阶段所携带的能量在同一流量下是稳定的，通过对上升段所包围的面积进行积分，可以得到稳定的特征点，进而求得被测气体的流量。

图3是确定积分区间起始点图。当被测气体的流速超过一定的范围时，即在大流量下，回波信号中所掺杂的噪声信号明显变大，可以看到，图3中归一化后的噪声信号幅值已经超过0.1，在图2中，已经说明同一流量下上升段的能量是稳定的，所以积分的区间结束点◇就可以确定为最大峰值点，那么，确定面积开始积分的区间起始点就十分关键。如果区间起始点选取在采样的第一个点，那么噪声信号就会被当作回波信号计算到面积中，不但增加计算量，将干扰引入计算，而且不能反映回波信号的真实信息，降低测量精度。而通过最大峰值点向前确定区间起始点，虽然最大峰值的位置(橫坐标)可能有所波动，但是可以有效去除噪声干扰，通过算法处理可以减小位置波动的影响。兼顾大流量和小流量上升段的峰值数，选取最大峰值点向前206到220点(每个回波信号采样2048个点，每个周期25个点)的位置，作为积分区间起始点，如图2中的o所示。积分起始点●的选取在图4中详细说明。

图4是确定选取积分起始点图。本发明运用穷举的方法,通过大量实验，在找到区间起始点后，从区间起始点到区间结束点(最大峰值点)之间，设定以回波信号幅值归一化的0.1为基准(平方后表示回波能量信号，则幅值归一化为0.01，根据图1所述的计算流程，先找到回波信号的峰值点，再平方即为回波能量信号的峰值点)，当归一化的峰值点幅值满足大于等于0.1时，就是所需的峰值点，这些峰值点作为后续面积计算效果最佳。如图4所示，在积分起始点●和区间结束点◇之间的峰值都是满足条件的峰值点。从图4中可以看出，在回波能量信号归一化幅值为0.01时，在区间起始点o前面与x轴所覆盖的面积很小，即信号前段的噪声信号干扰已经很小，可以忽略不计。

图5是基于回波能量积分方法的积分面积示意图。在确定面积积分的区间起始点、进而找到积分起始点以后，既可以进行具体的面积计算。图中的●所代表的就是我们已经找到满足设定条件的积分起始点和上升段各峰值点，包括图中所示的最大峰值点b、次大峰值点a、以及中间段的峰值点c。各峰值点相邻两点之间采用线性拟合的方法，构成一条线段，然后与这两点的纵坐标所代表的相邻的梯形底边(图中的实线)、横坐标轴共同组成一个小梯形，对每个小梯形的面积进行计算并累加，在保证精度的同时，简化计算过程，避免了对信号逐点计算所带来的计算量。对每个小梯形的面积计算如上述公式(2)所述。

图6是确定预设面积值原理图。找到所需要的积分起始点和上升段各峰值点后，我们需要设置一个面积的预设值，当累加面积值超过预设面积值时，便停止累加计算面积，开始求解所需的特征点的横坐标。那么，这个面积预设值同样十分关键。当我们选取的面积值较小时，对应的峰值点数就少，即所有覆盖的面积在中间段的峰值点c处，点3就是我们求解得到的特征点；当我们选取的面积值较大时，对应的峰值点数就多，即所有覆盖的面积在次大峰值点a处，点1就是我们求解得到的特征点。当回波能量信号波形发生波动时，回波能量信号波动造成峰值点幅值会变化，进而计算得到的累加面积值会发生变化，波动造成累加面积值增加或减少的部分，需要用累加面积值中最后一部分梯形的面积来补偿，即图中梯形s1或梯形s3的面积值(图中的两部分阴影面积)来补偿，因为面积预设值是一定的，所以需要补偿部分的面积值是相等的，即梯形s1的面积值等于梯形s3面积值，又因为梯形s3的上下底边明显小于梯形s1的上下底边，所以，梯形s3的高(中间段的峰值点c的横坐标与特征点3横坐标之间的距离)要比梯形s1的高(次大峰值的a的横坐标与特征点1横坐标之间的距离)大的多，即所求解得到的特征点1橫坐标的波动要明显小于特征点3横坐标的波动，那么，特征点1横坐标的稳定性要高于特征点3横坐标的稳定性。特征点的稳定性是保证气体流量测量精度的关键。所以，面积预设值应选取在次大峰值点a和最大峰值点b(图5中所示)之间，且不能超过最大峰值点b，兼顾大流量和小流量上升段的峰值数所带来的面积差异，选取归一化后幅值计算后的面积值，确定预设面积值sa为60到63之间时，效果最佳。

图7是基于回波能量积分的特征点计算方法。面积预设值选取在了次大峰值点a和最大峰值点b，当累加面积值大于预设面积值sa时，停止累加，已知面积差(梯形面积)和相邻峰值点的纵坐标(上下底边值)，就可以求解梯形的高度。次大峰值点a和最大峰值点b(如图5所示)的坐标分别为(x1,y1)、(x2,y2),其中，y1和y2都是归一化的回波能量信号的幅值，那么，这两点之间线性拟合的斜率k为

累加面积值与预设面积值sa之间的面积差值s'为

特征点1横坐标与最大峰值点b横坐标之间的距离为x,即梯形s'的高，那么,s'的面积值为

根据式(6)可以求出未知数x(只有一个解符合条件)：

则最终所求解的特征点1横坐标的值为(x2-x)。

图8是气体超声波流量计的硬件框图。主要包括一次仪表、dsp最小系统、fpga最小系统、驱动信号生成与放大电路、发射/接收信号通道切换电路、回波信号调理和采集电路、电源管理模块。

所述一次仪表由四个换能器和流量管组成，四个换能器分别固定安装在流量管上；四个换能器分别为换能器1、换能器2、换能器3和换能器4，组建双声道直射式结构；换能器为收发两用型，通过发射/接收信号通道切换电路控制，既可作为发射换能器用于发射激励信号，又可作为接收换能器用于接收回波信号。

所述fpga最小系统由fpga芯片、时钟、复位电路和flash(闪存)组成，其中，fpga芯片型号为ep2c8q208c8n。fpga内部包括dac(数模转换器)控制模块、rom(只读存储器)、分频模块、ram_2sport。

所述dsp最小系统由dsp芯片、时钟、外部看门狗、fram(铁电存储器)、人机接口、串口通信组成，其中，dsp芯片型号为tms320f28335。

所述电源管理模块为系统提供所需幅值的数字电源和模拟电源。

图9是气体超声流量计系统主控芯片dsp软件流程图。系统上电后，dsp进行初始化，完成对gpio(通用输入/输出)口的分配、定时器中断初始化、液晶显示初始化、流量计显示初始化等初始化任务。dsp向fpga发送激励起始信号，并切换发射、接收通道。检测定时器时间标志位是否为0，若为0，则直接进入dsp收到转存起始信号等待状态；若为1，则更新液晶显示，与上位机通信，并将定时器时间标志位置0。若dsp收到转存起始信号后，则开始读取fpga中已存的回波信号，计算积分面积，求得特征点横坐标，得到传播时间，求各声道的流量，进而输出最终流量q。