CN114254253A - 一种超声波流量测量中互相关算法的轻量化实现方法 - Google Patents
一种超声波流量测量中互相关算法的轻量化实现方法 Download PDFInfo
- Publication number
- CN114254253A CN114254253A CN202111148649.4A CN202111148649A CN114254253A CN 114254253 A CN114254253 A CN 114254253A CN 202111148649 A CN202111148649 A CN 202111148649A CN 114254253 A CN114254253 A CN 114254253A
- Authority
- CN
- China
- Prior art keywords
- cross
- discrete
- signals
- frequency domain
- echo
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 45
- 238000005259 measurement Methods 0.000 title claims abstract description 36
- 238000005314 correlation function Methods 0.000 claims abstract description 91
- 238000011144 upstream manufacturing Methods 0.000 claims abstract description 91
- 238000001914 filtration Methods 0.000 claims abstract description 37
- 238000005070 sampling Methods 0.000 claims abstract description 25
- 239000012530 fluid Substances 0.000 claims abstract description 23
- 238000012545 processing Methods 0.000 claims abstract description 15
- 238000006073 displacement reaction Methods 0.000 claims description 13
- 230000001131 transforming effect Effects 0.000 claims description 3
- 238000004364 calculation method Methods 0.000 abstract description 16
- 230000001186 cumulative effect Effects 0.000 abstract description 2
- 238000013461 design Methods 0.000 description 7
- 238000005516 engineering process Methods 0.000 description 6
- 108010076504 Protein Sorting Signals Proteins 0.000 description 4
- 238000012360 testing method Methods 0.000 description 4
- 230000008569 process Effects 0.000 description 3
- 230000000644 propagated effect Effects 0.000 description 3
- 238000010586 diagram Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 238000009776 industrial production Methods 0.000 description 2
- 238000009434 installation Methods 0.000 description 2
- 238000003860 storage Methods 0.000 description 2
- 230000009466 transformation Effects 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 150000001875 compounds Chemical class 0.000 description 1
- 230000021615 conjugation Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 238000007405 data analysis Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 230000007123 defense Effects 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 239000000843 powder Substances 0.000 description 1
- 230000001902 propagating effect Effects 0.000 description 1
- 238000003908 quality control method Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 230000007306 turnover Effects 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
- G06F17/141—Discrete Fourier transforms
- G06F17/142—Fast Fourier transforms, e.g. using a Cooley-Tukey type algorithm
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01F—MEASURING VOLUME, VOLUME FLOW, MASS FLOW OR LIQUID LEVEL; METERING BY VOLUME
- G01F1/00—Measuring the volume flow or mass flow of fluid or fluent solid material wherein the fluid passes through a meter in a continuous flow
- G01F1/66—Measuring the volume flow or mass flow of fluid or fluent solid material wherein the fluid passes through a meter in a continuous flow by measuring frequency, phase shift or propagation time of electromagnetic or other waves, e.g. using ultrasonic flowmeters
- G01F1/667—Arrangements of transducers for ultrasonic flowmeters; Circuits for operating ultrasonic flowmeters
Landscapes
- Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Mathematical Optimization (AREA)
- Data Mining & Analysis (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Theoretical Computer Science (AREA)
- Pure & Applied Mathematics (AREA)
- Algebra (AREA)
- Discrete Mathematics (AREA)
- Electromagnetism (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Fluid Mechanics (AREA)
- Measuring Volume Flow (AREA)
- Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)
Abstract
本发明公开了一种超声波流量测量中互相关算法的轻量化实现方法，包括步骤：流体管道的上、下游两个位置点，布置一对超声波换能器，所述超声波换能器兼具发送和接收超声波的功能；分别对两组所述回波信号进行信号采样，得到离散的数字信号；将信号由时域转到频域，再进行带通滤波处理；对两组信号在频域内进行离散互相关运算；对离散互相关函数进行傅里叶逆变换；对离散互相关函数进行三次样条插值处理，拟合关键部分曲线，求取最大值点，并据此得到飞行时间差；从而计算出瞬时及累积流量值。本发明提高了测量抗干扰能力、大幅降低了计算量，实现了互相关算法的轻量化，并使得互相关算法计算结果能够用于实际测量，且提高了测量精度。
Description
技术领域
本发明属于流量检测技术领域，涉及互相关算法在超声波流量计上的应用，特别是基于时差法原理的超声波流量计。
背景技术
计量是工业生产过程中的眼睛，能够实时准确的得知生产过程中的各个参数变化情况，对工业生产的安全与质量的把控有着重要作用，流量计量是计量科学技术的组成部分之一，它与国民经济、国防建设、科学研究有着密切的关系。现如今芯片技术以及流体力学的迅速发展使得流量测量技术迈向新台阶，不再单单要求测量过程中的稳定，而是更加追求测量结果的准确性，因此，仪器的测量精度和稳定性备受关注。
如今，一些发达国家在超声波流量测量技术处于国际领先水平，与之相比，我国尤其是在测量精度及稳定性方面仍然还有着较大差距。所以创新流量测算的方法，提升仪器准确度和稳定性迫在眉睫。
基于时差法的超声波流量计关键在于超声波的下游飞行时间与上游飞行时间两者之差即飞行时间差的准确程度，只有飞行时间差精确才能根据它与被测流体流速、流量之间的关系计算出相应的流量值。
目前测量飞行时间差的方法以阈值比较法为主，通过比较上、下游回波信号前沿到达某一阈值电压的时间来确定首波到达的时刻，期间通过高精度的时间数字芯片来计算时间获得上、下游飞行时间，再据此计算出飞行时间差，此方法在实际的超声波流量计应用中，尤其是在接收信号时，由于超声换能器的谐振特性以及超声波在传播过程中的衰减和噪声干扰，使得难以准确判断接收信号何时到达阈值电压，而且还可能会出现波形信号串波现象，导致流量计的适用场合较为局限、重复性差以及测量精度低下等一系列问题。
发明内容
针对时差式的超声波流量计现有技术中存在的缺陷，本发明的目的在于提供了一种更为轻量化设计的互相关算法来计算飞行时间差，进而避免了现有技术阈值比较法测量飞行时间差出现的较大误差。本发明采用轻量化设计的数字相关法来计算飞行时间差，提升超声波流量计的精度和抗噪声干扰能力。
为了达到上述目的，本发明一种超声波流量测量中互相关算法的轻量化实现方法，包括如下步骤：
S1、分别位于直管型流体管道的上、下游两个位置点布置一对超声波换能器，所述超声波换能器兼具发送和接收超声波的功能；
其中，两个位置点，分时交替发出超声波信号，经管道内流体传播，到达相对应的位置点分别接收回波信号；所述回波信号为电信号。再具体的说，两个位置点，通过微控制器分时交替输出的脉冲信号来驱动超声波换能器进而激发超声波信号；分时交替发出超声波信号，经管道内流体传播，到达相对应的位置点分别接收回波信号；所述回波信号为电信号。
这对换能器，通过微控制器输出的脉冲驱动信号给上游超声波换能器，上游换能器就激发出超声波信号，下游换能器则就会接收到上游回波信号；之后，两者交换一下，微控制器的驱动信号给下游换能器，下游换能器激发超声波信号，上游换能器就会接收到下游回波信号，至此得到了上下游的回波信号。
S2，分别对两组所述回波信号进行信号采样，得到离散的数字信号；
S3，对所述数字信号进行快速傅里叶变换，将两组所述回波信号由时域转到频域；
S4，在频域内对步骤S3所得两组所述回波信号进行带通滤波处理，分别得到下游滤波离散频域信号Y′(n)和上游滤波离散频域信号X′(n)；
S5，对步骤S4所得两组所述回波信号在频域内进行离散互相关运算；
具体说，对下游滤波离散频域信号Y′(n)取共轭运算获得共轭输出
再与上游滤波离散频域信号X′(n)做乘积运算，以得到互相关函数频域离散输出序列；
S6，对离散互相关函数进行傅里叶逆变换，将步骤S5所得两组所述回波信号由频域反变换回到时域；
具体为：将互相关函数频域离散输出序列进行快速傅里叶逆变换得到互相关函数时域离散输出序列Ryx(τ)；
S7，对离散互相关函数进行三次样条插值处理，拟合关键部分曲线，求取最大值点，并据此得到飞行时间差；
具体为：采用三次样条插值的方法对互相关函数时域离散输出序列进行插值拟合，得到互相关函数时域连续输出多项式，再求取三次拟合多项式二阶导数为零时对应的位移时间，即互相关函数输出最大值点对应的位移时间，从而得到飞行时间差的准确估计；
S8，根据得到飞行时间差与流体流速的关系，并根据流量定义计算出瞬时及累积流量值。
优选方式下，步骤S2的方法为：使用高速模数转换器对接收到的上游回波信号x(t)和下游回波信号y(t)分别进行高速采样，得到上游原始离散回波时域信号x(n)和下游原始离散回波时域信号y(n)。
步骤S3的方法为：采用快速傅里叶变换分别将上、下游的原始离散回波时域信号变换得到上、下游离散回波频域信号X(n)和Y(n)。
步骤S4的方法为：在频域内对X(n)和Y(n)进行带通滤波，以滤除杂波信号，从而得到上游滤波离散频域信号X′(n)和下游滤波离散频域信号Y′(n)。
步骤S8的方法为：依据飞行时间差估计依次计算得到流体流速、瞬时流量和累积流量。具体，可参考文献<<基于时差法的单声道气体超声波流量计的研究>>此篇文献的方法(危鄂元.基于时差法的单声道气体超声波流量计的研究 [D].浙江大学,2014.)。计算飞行时间差与流体流速的关系：
并根据流量定义计算出瞬时流量值
及累积流量值
最优方式下，所述发出和接收超声波信号的超声波换能器为集成式设计，即同一位置点为一个既能发出又能接收超声波信号的超声波换能器。
本发明实现了互相关算法的应用，不同于传统的阈值比较法，互相关法计算结果的准确与否只与上、下游两组回波信号的相似程度有关系，而且对上、下游回波信号进行了高速采样、滤波等处理，不会受信号失真和畸变的影响，提高了测量抗干扰能力，可应用于大流速环境下以及电磁环境较为复杂的场合，大幅度提高流量计的稳定性。
本发明采用轻量化设计的互相关法计算飞行时间差，在时域进行互相关函数计算，实际运算量十分庞大，在超声波流量计的硬件平台上实现起来困难，为了实现轻量化设计，避免在时域进行卷积运算，将两组回波信号转换到频域内进行乘积运算来获得互相关函数频域输出，通过将时域的卷积运算转换到频域进行乘积处理，大幅度降低互相关函数的计算量，提高运算效率，为在嵌入式设备中实现该算法提供了基础。
由上述可知，得到的互相关函数实际是一个离散序列，这样很难辨别出互相关函数的最大值，真实的互相关函数最大值较大概率落在某两个离散点之间，所以要对互相关函数离散序列进行三次样条插值来将离散序列进行数据拟合为连续曲线，便可以确定互相关函数的最大值，并求取三次拟合多项式二阶导数为零时对应的位移时间，即互相关函数输出最大值点对应的位移时间，从而得到飞行时间差的准确估计，使得互相关算法计算结果能够用于实际测量，仪器精度得到大幅度提升。
本发明一种超声波流量测量中互相关算法的轻量化实现与应用方法，具体包括以下内容：微控制器输出脉冲电信号，驱动超声波换能器产生超声波信号，经流体传播，到达接收换能器转换为回波电信号，超声波信号从上游换能器到下游换能器的时间为上游飞行时间，从下游换能器到上游换能器的时间为下游飞行时间，二者之差为飞行时间差；使用高速模数转换器分别对上、下游接收到的回波信号进行采样，得到上、下游原始离散回波时域信号；采用快速傅里叶变换，得到上、下游离散回波频域信号，分别进行频域带通滤波，得到上、下游离散回波频域滤波信号，提高了测量抗干扰能力；对下游离散回波频域滤波信号取共轭，再与上游离散回波频域滤波信号进行乘积运算，得到互相关函数频域输出，进而采用快速傅里叶逆变换，得到互相关函数时域离散输出，大幅降低计算量，实现了互相关算法的轻量化；采用三次样条插值方法，对互相关函数时域离散输出进行拟合，得到互相关函数时域连续输出，再求取三次拟合多项式二阶导数为零时对应的位移时间，即互相关函数输出最大值点对应的位移时间，从而得到飞行时间差的准确估计，使得互相关算法计算结果能够用于实际测量，实现了互相关算法的应用，且提高了测量精度。本发明将高速采样得到的离散回波时域信号转换至频域信号，开展频域带通滤波，提高了测量抗干扰能力，且通过频域的共轭和乘积运算，获得互相关函数频域输出，大幅降低了计算量，实现了互相关算法的轻量化，并进一步采用三次样条插值方法，拟合得到了互相关函数连续输出，从而获得准确的飞行时间差估计，使得互相关算法计算结果能够用于实际测量，且提高了测量精度。
总之，本发明的有益效果为：
本发明使用轻量化设计的互相关算法来计算超声波流量计的飞行时间差，与传统的阈值比较法相比而言，不仅仅提高了飞行时间差的准确程度，而且提高了超声波流量计的抗杂波信号干扰能力，使得应用场合更为丰富。同时，大幅度降低相关算法的运算量，从而降低对微控制器等硬件性能的要求，三次样条插值函数引入使得互相关函数的最大值点更加准确，进而飞行时间差的估计也更为真实。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1是超声波流量计测量模型；
图2是上、下游接收到的回波信号图；
图3是上、下游飞行时间差的示意图；
图4是互相关函数曲线图；
图5是轻量化互相关算法实现的系统流程框图；
图6是离散互相关序列三次样条插值后的连续输出拟合曲线图；
图7是上、下游原始离散回波信号的时域图像；
图8是上、下游原始离散回波信号的频域图像；
图9是经过频域带通滤波处理后的上、下游离散回波信号的时域图像；
图10是经过频域带通滤波处理后的上、下游离散回波信号的频域图像；
图11是互相关函数的时域离散序列输出图像。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的其他实施例，都属于本发明保护的范围。
超声波流量计测量模型如图1所示，超声波换能器兼具信号发送和接收的功能，当发送声波信号时称之为发送换能器，当接收回波信号时称之为接收换能器，测量过程中，微控制器输出脉冲电信号交替驱动上、下游换能器产生超声波信号，经流体传播后，到达接收换能器转换为回波信号，超声波信号从上游换能器到下游换能器的时间为上游飞行时间，从下游换能器到上游换能器的时间为下游飞行时间，下游飞行时间减去上游飞行时间即为飞行时间差。
实际接收到的两组回波信号图像如图2所示，可以看出上游回波信号和下游回波信号两组信号的相似程度极高，它们都是相应的超声波换能器发出相同频率的超声波信号经管道内流体传播后接收得到的正弦信号，只不过上游信号是顺流传播，而下游信号是逆流传播，所以它们的最大区别在于时间轴上有一定的时间延迟，而该时间延迟为对应下游回波信号与上游回波信号之间的飞行时间差。
上、下游飞行时间差的示意图如图3所示，微控制器的驱动信号均为方波脉冲信号，接收换能器接收得到的上、下游回波信号x(t)和y(t)是正弦信号，可以发现x(t)和y(t)的波形曲线几乎可以完全一样，但且两组信号之间有一定的时间延迟Δt。
互相关函数可以很好的表征两组信号的相关程度，当互相关函数值达到最大值时，其所对应的横坐标即为两组信号的时间延迟，因此通过计算下游回波信号与上游回波信号的互相关函数就可以得到飞行时间差。
采用互相关函数Ryx(τ)描述了上游回波信号x(t)与下游回波信号y(t)在不同时间位移τ下的相似程度，上游回波信号为x(t)和下游回波信号为y(t)的互相关函数Ryx(τ)可定义为
其中，τ代表时间位移；
互相关函数曲线如图4所示，当Ryx(τ)的横坐标时间位移τ达到时间位移τ0时，且此时出现互相关函数输出曲线最大值，则说明此时上游回波信号x(t)和下游回波信号y(t)最为相似，则时间位移τ0就是两个回波信号之间的飞行时间差Δt。
利用了高速模数转换器对接收到的上游回波信号x(t)和下游回波信号y(t) 分别进行高速采样，得到便于数据处理的离散数字信号：上游原始离散回波时域信号x(n)和下游原始离散回波时域信号y(n)，所以测量过程中实际处理的是离散信号，则离散互相关函数Ryx(m)定义为：
其中，m为x(n)与y(n)两个离散信号序列之间的数据位移个数，M为采样点数，n为数据序号。
根据式2可知，在时域进行M点的Ryx(m)离散互相关函数运算，需要M2次乘法运算和M2-M次加法运算，运算量十分庞大，为此采用在频域内来实现互相关函数的计算。
根据离散卷积定理，离散信号序列y(n)与x(n)进行M点的卷积运算可表示为：
其中，x(m)代表离散信号序列x(n)，y(m)代表离散信号序列y(n)翻转后向右偏移m点所得的离散序列，y(-m)代表离散序列y(m)翻转后的离散序列。
离散互相关函数Ryx(m)与Rxy(m)之间的关系为：
Ryx(m)＝Rxy(-m) (5)
根据式4、5可知，离散互相关函数Ryx(m)与y(n)和x(n)的卷积之间的关系为：
离散互相关函数可通过离散卷积进行计算，根据式6以及时域的卷积对应频域的乘积、时域的翻转对应频域的共轭性质，将y(n)和x(n)的卷积运算转化为频域信号的乘积运算方可简化运算。轻量化设计的互相关算法实现的具体流程如图5所示，采用快速傅里叶变换分别将上、下游的原始离散回波时域信号变换得到上、下游离散回波频域信号X(n)和Y(n)；在频域内对X(n)和Y(n)进行带通滤波，以滤除杂波信号，从而得到上游滤波离散频域信号X′(n)和下游滤波离散频域信号Y′(n)；对下游滤波离散频域信号Y′(n)取共轭运算获得共轭输出
再与上游滤波离散频域信号X′(n)做乘积运算，以得到互相关函数频域离散输出序列；将互相关函数频域离散输出序列进行快速傅里叶逆变换得到互相关函数时域离散输出序列Ryx(τ)。
完成上述过程M点的互相关运算，共需要2次M点快速傅里叶变换运算、M 次复数乘法运算和1次M点的快速傅里叶逆变换运算。借助快速傅里叶变换函 数和快速傅里叶逆变换函数完成互相关运算，只需要
次复数乘法运 算和3Mlog2M次复数加法运算，对比在时域进行M点的Ryx(m)离散互相关函数 运算，一共需要M2次乘法运算和M2-M次加法运算，尤其是当M越大时，可以 明显发现轻量化设计的互相关算法计算量有大幅度降低，该算法更加轻量化， 提高运算效率，降低了硬件要求以及仪器成本。
上述计算得到的互相关函数实际上是一个离散互相关函数序列，由于互相关函数离散输出最大值点与互相关函数连续输出最大值点不一定能够碰巧重合，二者对应的延迟时间一般有着微小的偏差，如果直接将互相关函数离散序列的最大值点所对应的时间延迟作为飞行时间差Δt的测量结果，会产生一定程度的测量误差。
为了可以找到真实的互相关函数最大值，可以使用三次样条插值函数来拟合离散序列点获得互相关函数的连续输出。互相关函数连续输出拟合曲线如图6 所示，由于进行1024点的三次样条插值函数的运算量也是比较巨大，为了继续降低微控制器的运算量，只在关键部分进行函数插值运算，只需找出互相关函数离散输出最大值点C0与其左邻点C-1、右邻点C1，它们之间存在一定间距，只需在这三个关键点进行两次三次样条插值运算即可，根据互相关函数图像可知互相关函数连续输出最大值点必将会分布在这三点之内的区间。
根据三次样条插值函数定义及C-1、C0和C1三个插值节点信息计算出区间
和
上的唯一多项式S1(x)和S2(x)，这样互相关函数离散输出序列进行了关键点位置的曲线拟合，得到了互相关函数关键部分的连续输出曲线。
由图6可知，实际的连续的互相关函数曲线最大值点C0′不同于原离散序列 C0点，互相关函数Ryx(τ0)在C0′达到最大值，再求取三次拟合多项式二阶导数为零时对应的位移时间，即C0′点所对应的时间位移τ0就是上游连续回波信号与下游连续回波信号之间的飞行时间差估计值Δt，通过引入三次样条插值函数，飞行时间差的估计结果更加准确。最后，就可以根据得到飞行时间差与流体流速的关系，并根据流量定义计算出瞬时及累积流量值。
本发明具体实施以管段DN100为例进行实验，流体介质为纯水，流速v为 2.0m/s，超声波换能器采用直接对射式的Z型安装法，正如图1所示，一对超声波换能器以安装角度θ为45°相向安装，理论谐振点为1MHz，驱动信号为20 个方波脉冲信号，通过高速模数转换器对回波信号进行数据采样，采样频率fs为 7.2MHz，采样点数为1024点。按着上述要求进行实验测试并导出一组采样信号的测试数据，将采样数据借助MATLAB软件仿真并结合实际流量计硬件设备进行数据分析。
根据香农采样定理，为使采样后的离散信号可以完整地保留原始回波信号的信息，并可以无失真的恢复原始波形信号，需保证采样频率fs高于被测信号有效频率序列中的最高频率fnmax的2倍。所使用的采样频率为7.2MHz，而驱动脉冲信号的频率为1MHz，完全满足采样定理的要求。
通过高速模数转换器采样得到的上、下游原始离散回波信号对应的时域图如图6所示，为了提高算法的运算效率，首先利用快速傅里叶变换函数将上游原始离散回波信号x(n)和下游原始离散回波信号y(n)从时域转换到频域，上、下游原始离散回波信号对应的频域图如图8所示。
在图7中，可以看到上、下游原始离散回波时域信号都接收到了20个波形，回波信号接收较为完整，离散采样信号能够较为准确地还原出连续回波信号曲线。但是从图8中可以看到，上、下游原始离散回波频域信号在0.8MHzˉ0.9MHz 之间以及1.1MHzˉ1.2MHz之间有一些噪声干扰信号。为增强测量抗干扰能力，需要对采样得到的上游原始离散回波信号和下游原始离散回波信号在频域中进行了带通滤波，剔除杂波信号的干扰。
频域带通滤波就是将带通滤波范围之外的杂波信号全部归为零，只保留处于带通频率之内的信号，即将f1≤|fn|≤f2频段内的信号保持不变，却将除此之外的信号全部归零。其中，f1、f2分别为带通滤波器的下限截止频率和上限截止频率；采样频率记为fs，fn为被测信号有效频率序列。
设定带通滤波的上限截止频率设置为1.08MHz，下限截止频率为0.92MHz，将快速傅里叶变换输出的上游离散回波频域信号X(n)和下游离散回波频域信号 Y(n)分别与理想带通滤波器频谱函数做乘积运算，从而完成了频域滤波，得到的上、下游滤波后的离散回波信号对应的时域图如图9所示，上、下游滤波后的离散回波信号对应的频域图如图10所示，通过观察频域内的信号可发现干扰信号已几乎被剔除。
将频域输出的上游滤波离散频域信号X′(n)与下游滤波离散频域信号Y′(n) 的共轭输出
做乘积运算，就可以得到互相关函数的频域离散输出，再利用快速傅里叶逆变换函数将互相关函数的频域离散输出转换为互相关函数的时域离散输出。从而完成了上游离散回波信号x(n)和下游离散回波信号y(n)的互相关函数运算，得到的互相关函数的时域离散输出如图11所示。
在图11中，互相关函数离散输出所对应的横坐标时间位移τ等于由式2中上游离散回波信号x(n)与下游离散回波信号y(n)的离散时间位移m除以采样频率fs所得到。
当互相关函数达到最大值时，两组回波信号x(n)和y(n)具有最大相似性，所对应的时间位移
即为下游离散回波信号y(n)与上游离散回波信号x(n)之间的飞行时间差Δt，借助MATLAB软件可计算出互相关函数离散序列计算的飞行时间差Δt为277.778ns。三个插值节点C-1、C0、C1的坐标值分别为(0.138889， 2674)、(0.277778，3420)、(0.4166667，1719)，根据三次样条插值定义可分别计算出区间(τC-1,τC0)和(τC0,τC1)上的唯一的三次多项式(即拟合曲线)S1(x)和 S2(x)，计算结果如下：
S1(x)＝-742444x3+431888x2-74326x+6655.3
S2(x)＝1034600x3-1049000x2+337022x-31432 (7)
如图6所示是根据S1(x)和S2(x)来拟合互相关函数关键离散点为连续曲线，可以看出，新峰值C0′点为三次样条插值后的互相关函数连续曲线的最大值点，求取三次拟合多项式二阶导数为零时对应的位移时间τ0即为下游连续回波信号与上游连续回波信号之间的飞行时间差Δt，根据三次样条插值结果可以重新计算出互相关函数连续曲线所的最大值所对应的时间位移(即飞行时间差)为 252.689ns，对比三次样条插值前后飞行时间差结果对比可知，两者并不相等，确实存在一定的偏差，数值上精确了25.089ns。因此，得到了更加准确的飞行时间差估计，流量值也随之相应的更为准确。
根据实际的计量检定测试实验，并与标准表对比，记录的实验表最大相对误差的绝对值为0.77％，实验误差在1.0％以内，且实验结果的重复性误差在 0.25％ˉ0.5％之间，由此可见基于轻量化设计的互相关算法的超声波流量计测量精度较高和测量重复性也较为稳定。
综上，本发明一种超声波流量测量中互相关算法的轻量化实现与应用方法，具体包括以下步骤：
S1，微控制器驱动超声波换能器激发出超声波信号，该信号经由管道流体后被相应的接收换能器转为回波信号；
S2，对接收到的回波信号进行信号采样，得到便于存储和处理的离散数字信号；
S3，对采样后的数字信号进行快速傅里叶变换，进而上、下游回波信号由时域转到频域；
S4，在频域内对上、下游回波信号进行带通滤波处理；
S5，对上、下游两组回波信号在频域内进行离散互相关运算；
S6，对离散互相关函数进行傅里叶逆变换，上、下游回波信号由频域反变换回到时域；
S7，对离散互相关函数进行三次样条插值处理，拟合关键部分曲线，求取最大值点，并据此得到飞行时间差的估计；
S8，根据得到飞行时间差与流体流速的关系，并根据流量定义计算出瞬时及累积流量值。
步骤S1微控制器驱动超声波换能器激发出超声波信号，该信号经由管道流体后被相应的接收换能器转为回波信号具体指的是：时差式超声波流量计的微控制器输出一组一定频率的脉冲信号，分时交替驱动上、下游超声波换能器，产生超声波发送信号，该超声波信号经管道内流体传播，到达相对应的接收换能器，进而由接收换能器将转换为回波信号，该回波信号为电信号。
步骤S2对接收到的回波信号进行信号采样，得到便于存储和处理的离散数字信号具体指的是：使用高速模数转换器对接收到的上游回波信号x(t)和下游回波信号y(t)分别进行高速采样，得到上游原始离散回波时域信号x(n)和下游原始离散回波时域信号y(n)。
步骤S3对采样后的数字信号进行快速傅里叶变换，上、下游回波信号由时域转到频域具体指的是：采用快速傅里叶变换分别将上、下游的原始离散回波时域信号变换得到上、下游离散回波频域信号X(n)和Y(n)。
步骤S4在频域内对上、下游回波信号进行带通滤波处理具体指的是：
在频域内对X(n)和Y(n)进行带通滤波，以滤除杂波信号，从而得到上游滤波离散频域信号X′(n)和下游滤波离散频域信号Y′(n)。
步骤S5对上、下游两组回波信号在频域内进行离散互相关运算具体指的是：
对下游滤波离散频域信号Y′(n)取共轭运算获得共轭输出
再与上游滤波离散频域信号X′(n)做乘积运算，以得到互相关函数频域离散输出序列。
步骤S6对离散互相关函数进行傅里叶逆变换，上、下游回波信号由频域反变换回到时域具体指的是：将互相关函数频域离散输出序列进行快速傅里叶逆变换得到互相关函数时域离散输出序列Ryx(τ)。
步骤S7对离散互相关函数进行三次样条插值处理，拟合关键部分曲线，求取最大值点，并据此得到飞行时间差的估计具体指的是：采用三次样条插值的方法对互相关函数时域离散输出序列进行插值拟合，得到互相关函数时域连续输出多项式，再求取三次拟合多项式二阶导数为零时对应的位移时间，即互相关函数输出最大值点对应的位移时间，从而得到飞行时间差的准确估计。
步骤S8根据得到飞行时间差与流体流速的关系，并根据流量定义计算出瞬时及累积流量值具体指的是：最后依据飞行时间差估计依次计算得到流体流速、瞬时流量和累积流量
本发明未尽事宜为公知技术。
上述实施例只为说明本发明的技术构思及特点，其目的在于让熟悉此项技术的人士能够了解本发明的内容并据以实施，并不能以此限制本发明的保护范围。凡根据本发明精神实质所作的等效变化或修饰，都应涵盖在本发明的保护范围之内。
Claims (4)
1.一种超声波流量测量中互相关算法的轻量化实现方法，其特征在于，包括如下步骤：
S1、分别位于直管型流体管道的上、下游两个位置点布置一对超声波换能器，所述超声波换能器兼具发送和接收超声波的功能；
其中，两个位置点，分时交替发出超声波信号，经管道内流体传播，到达相对应的位置点分别接收回波信号；所述回波信号为电信号；
S2，分别对两组所述回波信号进行信号采样，得到离散的数字信号；
S3，对所述数字信号进行快速傅里叶变换，将两组所述回波信号由时域转到频域；
S4，在频域内对步骤S3所得两组所述回波信号进行带通滤波处理，分别得到下游滤波离散频域信号Y′(n)和上游滤波离散频域信号X′(n)；
S5，对步骤S4所得两组所述回波信号在频域内进行离散互相关运算；
具体说，对下游滤波离散频域信号Y′(n)取共轭运算获得共轭输出
再与上游滤波离散频域信号X′(n)做乘积运算，以得到互相关函数频域离散输出序列；
S6，对离散互相关函数进行傅里叶逆变换，将步骤S5所得两组所述回波信号由频域反变换回到时域；
具体为：将互相关函数频域离散输出序列进行快速傅里叶逆变换得到互相关函数时域离散输出序列Ryx(τ)；
S7，对离散互相关函数进行三次样条插值处理，拟合关键部分曲线，求取最大值点，并据此得到飞行时间差；
具体为：采用三次样条插值的方法对互相关函数时域离散输出序列进行插值拟合，得到互相关函数时域连续输出多项式，再求取三次拟合多项式二阶导数为零时对应的位移时间，即互相关函数输出最大值点对应的位移时间，从而得到飞行时间差的准确估计；
S8，根据得到飞行时间差与流体流速的关系，并根据流量定义计算出瞬时及累积流量值。
2.根据权利要求1所述超声波流量测量中互相关算法的轻量化实现方法，其特征在于，步骤S2的方法为：
使用高速模数转换器对接收到的上游回波信号x(t)和下游回波信号y(t)分别进行高速采样，得到上游原始离散回波时域信号x(n)和下游原始离散回波时域信号y(n)。
3.根据权利要求1所述超声波流量测量中互相关算法的轻量化实现方法，其特征在于，步骤S3的方法为：
采用快速傅里叶变换分别将上、下游的原始离散回波时域信号变换得到上、下游离散回波频域信号X(n)和Y(n)。
4.根据权利要求1所述超声波流量测量中互相关算法的轻量化实现方法，其特征在于，步骤S4的方法为：
在频域内对X(n)和Y(n)进行带通滤波，以滤除杂波信号，从而得到上游滤波离散频域信号X′(n)和下游滤波离散频域信号Y′(n)。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111148649.4A CN114254253B (zh) | 2021-09-29 | 2021-09-29 | 一种超声波流量测量中互相关算法的轻量化实现方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111148649.4A CN114254253B (zh) | 2021-09-29 | 2021-09-29 | 一种超声波流量测量中互相关算法的轻量化实现方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114254253A true CN114254253A (zh) | 2022-03-29 |
| CN114254253B CN114254253B (zh) | 2024-09-20 |
Family
ID=80790338
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202111148649.4A Active CN114254253B (zh) | 2021-09-29 | 2021-09-29 | 一种超声波流量测量中互相关算法的轻量化实现方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114254253B (zh) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114910132A (zh) * | 2022-04-12 | 2022-08-16 | 上海理工大学 | 一种基于超声回波信号轮廓的特征点识别方法 |
| CN118010116A (zh) * | 2024-04-09 | 2024-05-10 | 潍坊奥博仪表科技发展有限公司 | 一种抗干扰的超声波流量计量方法 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20170138773A1 (en) * | 2014-07-01 | 2017-05-18 | Pcme Limited | Method Of Measuring Time Of Flight Of An Ultrasound Pulse |
| CN108534855A (zh) * | 2018-05-11 | 2018-09-14 | 大连理工大学 | 一种新型超声波气体流量计及其测量方法 |
| CN110646042A (zh) * | 2019-10-16 | 2020-01-03 | 上海交通大学 | 一种用于低功耗超声流量计飞行时间差计算的互相关插值方法 |
| US20200141777A1 (en) * | 2017-06-27 | 2020-05-07 | Sagemcom Energy & Telecom Sas | Method for measuring a speed of a fluid |
| CN113375737A (zh) * | 2020-06-15 | 2021-09-10 | 郑州大学 | 时差式超声气体流量计的流速计量方法 |
-
2021
- 2021-09-29 CN CN202111148649.4A patent/CN114254253B/zh active Active
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20170138773A1 (en) * | 2014-07-01 | 2017-05-18 | Pcme Limited | Method Of Measuring Time Of Flight Of An Ultrasound Pulse |
| US20200141777A1 (en) * | 2017-06-27 | 2020-05-07 | Sagemcom Energy & Telecom Sas | Method for measuring a speed of a fluid |
| CN108534855A (zh) * | 2018-05-11 | 2018-09-14 | 大连理工大学 | 一种新型超声波气体流量计及其测量方法 |
| CN110646042A (zh) * | 2019-10-16 | 2020-01-03 | 上海交通大学 | 一种用于低功耗超声流量计飞行时间差计算的互相关插值方法 |
| CN113375737A (zh) * | 2020-06-15 | 2021-09-10 | 郑州大学 | 时差式超声气体流量计的流速计量方法 |
Non-Patent Citations (1)
| Title |
|---|
| 董雅琳: "管道传输气体超声波流量计的研制", 《中国硕士论文数据库》, 15 February 2019 (2019-02-15) * |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114910132A (zh) * | 2022-04-12 | 2022-08-16 | 上海理工大学 | 一种基于超声回波信号轮廓的特征点识别方法 |
| CN118010116A (zh) * | 2024-04-09 | 2024-05-10 | 潍坊奥博仪表科技发展有限公司 | 一种抗干扰的超声波流量计量方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN114254253B (zh) | 2024-09-20 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN114254253A (zh) | 一种超声波流量测量中互相关算法的轻量化实现方法 | |
| CN107290564B (zh) | 一种基于相位差的超声波流速测量方法 | |
| CN103868555B (zh) | 一种用于超声波流量计的循环时差检测方法 | |
| CN105890685B (zh) | 一种基于累积相位差的超声波流量测量装置 | |
| CN107860430B (zh) | 一种基于时差法的超声波气体流量计时间差测量方法 | |
| CN103353550A (zh) | 一种测量电力系统信号频率及谐波参数的方法 | |
| Sun et al. | A novel signal processing method based on cross-correlation and interpolation for ToF measurement | |
| CN110646042A (zh) | 一种用于低功耗超声流量计飞行时间差计算的互相关插值方法 | |
| CN105180853B (zh) | 一种电磁超声金属材料测厚方法 | |
| CN104748704A (zh) | 薄壁结构超声共振测厚频谱分析内插校正方法 | |
| CN114001804A (zh) | 一种基于时间差法超声波计量装置的校准方法及系统 | |
| CN113916719A (zh) | 一种流体密度与流速在线同步检测系统与检测方法 | |
| CN110987102A (zh) | 一种高抗干扰高精度气体超声流量计及其测量方法 | |
| CN113124948A (zh) | 一种基于fpga与互相关法的高精度时差测量方法 | |
| CN103744066A (zh) | 数字正交鉴相与匹配滤波的优化方法 | |
| CN106054154A (zh) | 一种基于icpf的机动目标的步进频率信号高分辨率成像方法 | |
| CN113917184B (zh) | 一种超声波测速方法及装置 | |
| CN112304376B (zh) | 基于数据融合的超声波流量计流量测量方法 | |
| CN105738482A (zh) | 一种fir多相内插滤波超声相控阵波束延时方法 | |
| CN212206222U (zh) | 一种高抗干扰高精度气体超声流量计 | |
| CN109059815B (zh) | 一种基于超声测量气液两相流液膜厚度的方法 | |
| CN115808213A (zh) | 基于tdc方案的改进型气体计量方法 | |
| CN111157066A (zh) | 基于第一包络重合度的气体超声波流量计渡越时间计算方法 | |
| CN103486987B (zh) | 一种利用曲线拟合提高超声测厚精度的方法 | |
| CN105759255A (zh) | 一种cic多相内插滤波超声相控阵波束延时方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |