CN112484798A - 一种基于改进卡尔曼滤波器的时差超声波流量测量方法 - Google Patents

一种基于改进卡尔曼滤波器的时差超声波流量测量方法 Download PDF

Info

Publication number
CN112484798A
CN112484798A CN202011364228.0A CN202011364228A CN112484798A CN 112484798 A CN112484798 A CN 112484798A CN 202011364228 A CN202011364228 A CN 202011364228A CN 112484798 A CN112484798 A CN 112484798A
Authority
CN
China
Prior art keywords
time
time difference
ultrasonic
filtering
flow
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
Application number
CN202011364228.0A
Other languages
English (en)
Other versions
CN112484798B (zh
Inventor
不公告发明人
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Jiangxi Zhongke Intelligent Aquaculture Research Co Ltd
Original Assignee
Jiangxi Zhongke Intelligent Aquaculture Research Co Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Jiangxi Zhongke Intelligent Aquaculture Research Co Ltd filed Critical Jiangxi Zhongke Intelligent Aquaculture Research Co Ltd
Priority to CN202011364228.0A priority Critical patent/CN112484798B/zh
Publication of CN112484798A publication Critical patent/CN112484798A/zh
Application granted granted Critical
Publication of CN112484798B publication Critical patent/CN112484798B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • 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
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01P—MEASURING LINEAR OR ANGULAR SPEED, ACCELERATION, DECELERATION, OR SHOCK; INDICATING PRESENCE, ABSENCE, OR DIRECTION, OF MOVEMENT
    • G01P5/00—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft
    • G01P5/24—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft by measuring the direct influence of the streaming fluid on the properties of a detecting acoustical wave
    • G01P5/245—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft by measuring the direct influence of the streaming fluid on the properties of a detecting acoustical wave by measuring transit time of acoustical waves
    • 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/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • 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/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Theoretical Computer Science (AREA)
  • Computational Mathematics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Probability & Statistics with Applications (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Operations Research (AREA)
  • Acoustics & Sound (AREA)
  • Multimedia (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Electromagnetism (AREA)
  • Fluid Mechanics (AREA)
  • Computing Systems (AREA)
  • Measuring Volume Flow (AREA)

Abstract

本发明公开了一种基于改进卡尔曼滤波器的时差超声波流量测量方法，其步骤包括，步骤（1）时差信号测量：步骤（2）时差信号的滤波：步骤（3）累计流量的计算。本发明提出了测量噪声协方差和过程噪声协方差的更新方法，并由此提出了一种改进的卡尔曼滤波器。然后，针对时差法超声波流量计中时差信号的滤波，提出了基于改进卡尔曼滤波器的时差信号滤波方法，及给出了基于改进卡尔曼滤波器的时差法超声波流量测量方案。通过本发明提供的技术方案，明显地改善了时差法超声波流量计中时差信号的滤波效果，从而有效地提高了时差法超声波流量计的测量精度。

Description

一种基于改进卡尔曼滤波器的时差超声波流量测量方法
技术领域
本发明属于水表超声波流量计量领域，具体涉及一种基于改进卡尔曼滤波器的时差超声波流量测量方法，适用于各种管道中水流体的流量测量。
背景技术
随着智慧水务和物联网技术的发展，超声波水表应运而生。与传统的机械水表相比，超声波水表具有非接触式、压损小、精度高、功耗低、使用寿命长等突出优点，必将成为下一代智能水表的重要选择。超声波水表所采用的超声波测量技术有多种方法，常用的测量方法有时差法、相关法、噪声法、多普勒法和波束偏移法等。
时差法超声波测量的原理是：超声波信号在流体顺流与逆流传播时的传播速度不同，进而产生传播时间差值。由于时间差值的大小与流速大小有关，因此，通过测量时间差值的大小，就可以测量流体的流速。
在超声波水表系统中，由于测量过程会引入噪声，使得测量结果与实际流速值产生一定偏差。偏差值的大小与流体中是否有气泡或颗粒物、超声波信号质量、电路板PCB布局、温度补偿超声波声速值误差、外部电磁干扰等有关。为了减少外部噪声对超声波水表的干扰，提高水表的精度，需要对测得的时差数据进行滤波。目前使用较多的是中位值平均滤波方法，其功能是消除大幅度的脉冲干扰，算数平均滤波是将多个数据取平均值，可以使数据波动变缓，从而使数据稳定可靠。但是，中位值滤波方法会产生较大的误差。因此就有人提出了卡尔曼滤波方法，卡尔曼滤波后的数据波动要比中值均值滤波要小得多，而且卡尔曼滤波的快速性和快速响应性远大于中值均值滤波。
但是如果时差信号快速变化，则卡尔曼滤波并不能获得良好的滤波效果，因此又有人提出了卡尔曼滤波算法与算术平均算法结合的数据滤波算法，但是算术平均与卡尔曼滤波之间的切换条件不好确定(因为涉及到时差数据的采样率Fs、求最大和最小值时用到的采样数量、切换门限值δ等多个阴虚)；此外，算术平均滤波对于快速变化的时差值，滤波效果也不太理想；而且，卡尔曼滤波器对于快速变化的时差值，会跟不上时差的变化而导致较大的误差，导致流量的累积计算可能会存在较大的误差。
发明内容
本发明要解决的技术问题是提供一种基于改进卡尔曼滤波器的时差超声波流量测量方法，通过更新卡尔曼滤波器中的过程噪声协方差矩阵和测量噪声协方差矩阵，极大地改善了超声波流量计中时差信号的滤波效果，特别是时差信号快速变化的情况下的滤波效果的改善更为明显，从而有效地提高了时差法超声波流量计的测量精度。
本发明通过以下技术方案解决上述技术问题，
一种基于改进卡尔曼滤波器的时差超声波流量测量方法，其步骤包括，
步骤(1)、时差信号测量：
沿水流方向前后设置第一超声波换能器P1和第二超声波换能器P2；规定换能器与反射片的距离为s，超声波在水中的速度为c，反射片与水平面夹角为45°，两反射片的中心距离为L，管道直径为D，水的正向流速即线速度vl为v；
第一超声波换能器P1发出超声波脉冲，在第二超声波换能器P2上接收超声波脉冲，并计算接收到超声波时刻和发送超声波时刻的时间差t12，则有：
Figure BDA0002804955810000021
通过第二超声波换能器P2发出超声波脉冲，在第一超声波换能器P1上接收超声波脉冲，并计算接收到超声波时刻和发送超声波时刻的时间差t21，则有：
Figure BDA0002804955810000022
最后，获得时间差Δt＝t21-t12，根据式(1)和式(2)，可求得顺逆流时间差Δt：
Figure BDA0002804955810000023
步骤(2)、时差信号的滤波：
滤波处理是为了去掉原始数据中的随机误差，在超声波流量测量装置中，由于测量过程会引入外部噪声，使得测量结果与实际流速值产生一定偏差值，偏差值的大小与流体中是否有气泡或颗粒物、超声波信号质量、电路板PCB布局、温度补偿超声波声速值误差、外部电磁干扰等有关，为了减少外部噪声对超声波水表的干扰，提高水表的精度，需要对测得的时差数据进行滤波；用xk表示时刻k时时差Δt(kT]的状态值，其中T为采样周期，具体如下：
xk＝Δt(kT]；
在卡尔曼滤波器中，系统状态方程如下：
xk＝Φk-1xk-1+wk-1， (4)；
Figure BDA0002804955810000031
式中，xk为k时刻的系统状态，xk-1为k-1时刻的系统状态，Φk-1为k-1时刻的状态转移矩阵，wk-1为k-1时刻的系统噪声，Qk为k时刻的系统噪声；
用zk表示时刻k时时差Δt(kT]的测量值，则时差Δt的测量方程如下：
zk＝Hkxk+vk， (6)；
Figure BDA0002804955810000032
式中，Hk表示k时刻的测量系统的参数，vk为k时刻的测量噪声，Rk为k时刻的测量噪声的协方差矩阵；
在系统中定义Φk-1＝1、Hk＝1，那么，时差数据的改进卡尔曼滤波方法如下：
(2.1)、初始条件计算：
Figure BDA0002804955810000033
式中，E[···]表示数学期望运算，
Figure BDA0002804955810000034
表示0时刻的状态值，P0(+)表示0时刻的误差协方差；
(2.2)、利用k-1时刻的误差协方差Pk-1(+)、k-1时刻的状态转移矩阵Φk-1和k-1时刻的系统噪声Qk-1，计算k时刻的状态预测
Figure BDA0002804955810000035
和误差协方差预测Pk(-)，具体如下：
Figure BDA0002804955810000036
Figure BDA0002804955810000041
(2.3)、利用k时刻的误差协方差预测Pk(-)、k时刻的测量系统的参数Hk和k时刻的测量噪声的协方差矩阵Rk，计算k时刻的卡尔曼增益
Figure BDA0002804955810000042
Figure BDA0002804955810000043
(2.4)、利用式(10)和式(11)计算获得的k时刻的
Figure BDA0002804955810000044
和Pk(-)，计算误差协方差更新矩阵如下：
Figure BDA0002804955810000045
(2.5)、利用k时刻的卡尔曼增益
Figure BDA0002804955810000046
k时刻的状态预测
Figure BDA0002804955810000047
给定的初始估计值
Figure BDA0002804955810000048
和输入数据zk，递归计算k时刻的状态值
Figure BDA0002804955810000049
如下：
Figure BDA00028049558100000410
(2.6)、k时刻的测量噪声协方差Rk和k-1时刻的过程噪声协方差Qk-1更新如下：
Figure BDA00028049558100000411
Figure BDA00028049558100000412
Figure BDA00028049558100000413
Figure BDA00028049558100000414
Ratio＝R1/Q0 (18)；
Rk＝Rk(+)/Ratio (19)；
Qk-1＝Qk-1(+)/Ratio (20)；
(2.7)、将步骤(2.6)更新后的k时刻的测量噪声协方差Rk和k-1时刻的过程噪声协方差Qk-1代入步骤(2.2)-(2.6)中，递归计算k时刻的状态值
Figure BDA00028049558100000415
(2.8)、按步骤(2.1)-(2.7)对式(6)中测得的时差信号zk进行滤波，通过滤波，得到滤波后的时差数据：
Figure BDA0002804955810000051
步骤(3)、累计流量的计算：
在计算瞬时流速和流量时，用到的是线平均流速，由流体力学可知，面平均流速根据流体的不同状态，有着不同的修正系数；
设修正系数为K，计算流体的累计流量Q为：
Figure BDA0002804955810000052
式中，vl为流体线速度，根据式(3)得到，vl值为：
Figure BDA0002804955810000053
为了获得更好的技术效果，步骤(3)前还有步骤(2.9)超声波声速的温度补偿：
由于超声波的声速随着温度的变化而变化，因此，为了准确计算式(23)中流体线速度vl值，需要对超声波的声速进行温度补偿；
(2.9.1)通过测量声速与温度的关系，得到基于曲线拟合得到声速与温度的关系曲线；
(2.9.2)制作超声波速度与温度的关系表格，并利用表格进行基于查表法的温度补偿，利用查表法得到的声速c，结合计算获得的滤波后的时差信号
Figure BDA0002804955810000054
计算该温度下流体线速度校正值vl'：
Figure BDA0002804955810000055
在上述流量计算的方法中，其中步骤(2)中时差信号的滤波算法，利用了本发明提出的改进卡尔曼滤波方法，改善了时差信号的滤波效果，从而有效地提高了时差法超声波流量计的测量精度。
本发明提出一种改进的卡尔曼滤波器，并提出一种基于改进卡尔曼滤波的时差法超声波流量测量方案。在该方案中，针对，提出了测量噪声协方差和过程噪声协方差的更新方法，并由此提出了一种改进的卡尔曼滤波器。然后，针对时差法超声波流量计中时差信号的滤波，提出了基于改进卡尔曼滤波器的时差信号滤波方法，及给出了基于改进卡尔曼滤波器的时差法超声波流量测量方案。通过本发明提供的技术方案，明显地改善了时差法超声波流量计中时差信号的滤波效果，从而有效地提高了时差法超声波流量计的测量精度。
附图说明
图1为本发明实施例超声波流量测量结构示意图；
图2为本发明实施例水流体下0～40℃声速-温度拟合曲线；
图3为本发明实施例超声波流量测量流程图；
图4为本发明实施例和对比例的平稳时差数据时的滤波结果；
图5为本发明实施例和对比例的时差数据线性降低时的滤波结果；
图6为本发明实施例和对比例的时差数据存在快速变化时的滤波结果；
图7为本发明对比例超声波流量测量流程图；
图8为本发明对比例卡尔曼滤波算法与算术平均算法结合的数据滤波算法流程图。
具体实施方式
下面结合附图和实施例，进一步详细阐述本发明的内容。
实施例1
一种基于改进卡尔曼滤波器的时差超声波流量测量方法，其步骤包括，
步骤(1)、时差信号测量：
如图1所示，沿水流方向前后设置第一超声波换能器P1和第二超声波换能器P2；规定换能器与反射片的距离为s，超声波在水中的速度为c，反射片与水平面夹角为45°，两反射片的中心距离为L，管道直径为D，水的正向流速即线速度vl为v；
第一超声波换能器P1发出超声波脉冲，在第二超声波换能器P2上接收超声波脉冲，并计算接收到超声波时刻和发送超声波时刻的时间差t12，则有：
Figure BDA0002804955810000061
通过第二超声波换能器P2发出超声波脉冲，在第一超声波换能器P1上接收超声波脉冲，并计算接收到超声波时刻和发送超声波时刻的时间差t21，则有：
Figure BDA0002804955810000071
最后，获得时间差Δt＝t21-t12，根据式(1)和式(2)，可求得顺逆流时间差Δt：
Figure BDA0002804955810000072
步骤(2)、时差信号的滤波：
滤波处理是为了去掉原始数据中的随机误差，在超声波流量测量装置中，由于测量过程会引入外部噪声，使得测量结果与实际流速值产生一定偏差值，偏差值的大小与流体中是否有气泡或颗粒物、超声波信号质量、电路板PCB布局、温度补偿超声波声速值误差、外部电磁干扰等有关，为了减少外部噪声对超声波水表的干扰，提高水表的精度，需要对测得的时差数据进行滤波；
用xk表示k时刻的时差Δt(kT]的状态值，其中T为采样周期，具体表示如下：
xk＝Δt(kT]；
在卡尔曼滤波器中，系统状态方程如下：
xk＝Φk-1xk-1+wk-1 (4)；
Figure BDA0002804955810000073
式中，xk为k时刻系统的状态，xk-1为k-1时刻的系统状态，Φk-1为k-1时刻的状态转移矩阵，wk-1为k-1时刻的系统噪声，Qk为k时刻的系统噪声；
用zk表示时刻k时时差Δt(kT]的测量值，则时差Δt的测量方程如下：
zk＝Hkxk+vk (6)；
Figure BDA0002804955810000074
式中，Hk表示k时刻的测量系统的参数，vk为k时刻的测量噪声，Rk为k时刻的测量噪声的协方差矩阵；
在系统中定义Φk-1＝1、Hk＝1，那么，时差数据的改进卡尔曼滤波方法如下：
(2.1)、初始条件计算：
Figure BDA0002804955810000081
式中，E[···]表示数学期望运算，
Figure BDA0002804955810000082
表示0时刻的状态值，P0(+)表示0时刻的误差协方差；
(2.2)、利用k-1时刻的误差协方差Pk-1(+)、k-1时刻的状态转移矩阵Φk-1和k-1时刻的系统噪声Qk-1，计算k时刻的状态预测
Figure BDA0002804955810000083
和误差协方差预测Pk(-)，具体如下：
Figure BDA0002804955810000084
Figure BDA0002804955810000085
(2.3)、利用k时刻的误差协方差预测Pk(-)、k时刻的测量系统的参数Hk和k时刻的测量噪声的协方差矩阵Rk，计算k时刻的卡尔曼增益
Figure BDA0002804955810000086
Figure BDA0002804955810000087
(2.4)、利用k时刻的卡尔曼增益
Figure BDA0002804955810000088
和误差协方差预测Pk(-)，计算k时刻的误差协方差Pk(+)，更新矩阵如下：
Figure BDA0002804955810000089
(2.5)、利用k时刻的卡尔曼增益
Figure BDA00028049558100000810
k时刻的状态预测
Figure BDA00028049558100000811
给定的初始估计值
Figure BDA00028049558100000812
和输入数据zk，递归计算k时刻的状态值
Figure BDA00028049558100000813
如下：
Figure BDA00028049558100000814
(2.6)、k时刻的测量噪声协方差Rk和k-1时刻的过程噪声协方差Qk-1更新如下：
Figure BDA00028049558100000815
Figure BDA00028049558100000816
Figure BDA0002804955810000091
Figure BDA0002804955810000092
Ratio＝R1/Q0 (18)；
Rk＝Rk(+)/Ratio (19)；
Qk-1＝Qk-1(+)/Ratio (20)；
(2.7)、将步骤(2.6)更新后的k时刻的测量噪声协方差Rk和k-1时刻的过程噪声协方差Qk-1代入步骤(2.2)-(2.6)中，递归计算k时刻的状态值
Figure BDA0002804955810000093
(2.8)、按步骤(2.2)-(2.7)对式(6)中测得的时差信号zk进行滤波，通过滤波，得到滤波后的时差数据：
Figure BDA0002804955810000094
(2.9)、超声波声速的温度补偿：
由于超声波的声速随着温度的变化而变化，因此，为了准确计算步骤(4)中式(23)中线速度vl值，需要对超声波的声速进行温度补偿；
(2.9.1)通过测量声速与温度的关系，得到基于曲线拟合得到声速与温度的关系曲线；
在容器中，保持一定间隔设置两个超声波换能器，其中一个发送超声波，另一个接收超声波，获得发送与接收时间差，通过间隔距离与时间差获得超声波的声速，之后改变容器中流体的温度，获得超声波的声速与温度的关系曲线，并制作超声波速度与温度的关系表格；
(2.9.2)利用超声波速度与温度的关系表格进行基于查表法的温度补偿；
利用查表法得到的声速c，结合计算获得的滤波后的时差信号
Figure BDA0002804955810000095
计算该温度下流体线速度校正值vl'：
Figure BDA0002804955810000096
步骤(3)、累计流量的计算：
在计算瞬时流速和流量时，用到的是线平均流速，由流体力学可知，面平均流速根据流体的不同状态，有着不同的修正系数；
假设修正系数为K，计算流体的累计流量Q为：
Figure BDA0002804955810000101
式中，vl为流体线速度，根据式(3)得到，vl值为：
Figure BDA0002804955810000102
结合步骤(2.9.2)获得的该温度下流体线速度校正值vl'，
流体的累计流量Q修正为：
Figure BDA0002804955810000103
在上述流量计算的方法中，其中步骤(2)中时差信号的滤波算法，利用了本发明提出的改进卡尔曼滤波方法，改善了时差信号的滤波效果，从而有效地提高了时差法超声波流量计的测量精度，流程图见3。
(一)、用仿真来验证基于改进卡尔曼滤波器的时差信号滤波效果
图4给出了基于改进卡尔曼滤波器的平稳时差数据(真实时差＝2ns)时的滤波结果；图5给出了时差数据线性降低(其中包含一段抖动的数据)时的滤波结果，其中时差从2ns降至0ns；图6给出了时差数据存在快速变化时的滤波结果。
相比现有方案中，本发明在时差信号平稳时，基于改进卡尔曼滤波器的时差数据滤波效果稍差于基于传统卡尔曼滤波器的时差信号滤波方法；但是，对于快速变化的时差信号，本发明滤波效果有了显著的提高，从而能有效提升时差法超声波流量计的测量精度。
对比例
时差法超声波测量的原理是：超声波信号在流体顺流与逆流传播时的传播速度不同，进而产生传播时间差值。由于时间差值的大小与流速大小有关，因此，通过测量时间差值的大小，就可以测量流体的流速，流程图如图7所示；超声波流量测量结构示意图如图1所示。
图1中，沿水流方向前后设置第一超声波换能器P1和第二超声波换能器P2；规定换能器与反射片的距离为s，超声波在水中的速度为c，反射片与水平面夹角为45°，两反射片的中心距离为L，管道直径为D，水的正向流速即线速度vl为v；
第一超声波换能器P1发出超声波脉冲，在第二超声波换能器P2上接收超声波脉冲，并计算接收到超声波时刻和发送超声波时刻的时间差t12，则有：
Figure BDA0002804955810000111
通过第二超声波换能器P2发出超声波脉冲，在第一超声波换能器P1上接收超声波脉冲，并计算接收到超声波时刻和发送超声波时刻的时间差t21，则有：
Figure BDA0002804955810000112
最后，获得时间差Δt＝t21-t12，根据式(1)和式(2)，可求得顺逆流时间差Δt：
Figure BDA0002804955810000113
由此，可以得到线速度vl：
Figure BDA0002804955810000114
在计算瞬时流速和流量时，用到的是线平均流速。由流体力学可知，面平均流速根据流体的不同状态，有着不同的修正系数。假设修正系数为K，则流体的累计流量Q为：
Figure BDA0002804955810000115
由式(23)可知，当顺逆流时间差Δt越大时，水流线速度越大，于是一定时间内的累计流量越大；根据式(22)可以发现，在修正系数K、管道直径D、反射片距离L确定的情况下，一定时间内的累计流量Q只与时间差值Δt和超声波传播速度c有关，因此，时间差值Δt和超声波传播速度c的测量成为超声波流量计的关键。
数据滤波处理是为了去掉原始数据中的随机误差，以提高数据质量的一种方法。在超声波水表系统中，由于测量过程会引入噪声，使得测量结果与实际流速值产生一定偏差。偏差值的大小与流体中是否有气泡或颗粒物、超声波信号质量、电路板PCB布局、温度补偿超声波声速值误差、外部电磁干扰等有关。为了减少外部噪声对超声波水表的干扰，提高水表的精度，需要对测得的时差数据进行滤波。
现有技术中有的采用了中位值平均滤波方法。中位值平均滤波是中位值滤波和算数平均滤波结合起来的滤波方法。中位值滤波最重要的功能是消除大幅度的脉冲干扰，算数平均滤波是将多个数据取平均值，可以使数据波动变缓，从而使数据稳定可靠。但是，中位值滤波方法会产生较大的误差。
现有技术中还有的提出了卡尔曼滤波方法。而且，结果表明，中位值滤波和卡尔曼滤波这两种滤波方法都能有效防止脉冲干扰，使得数据更加平滑，但是卡尔曼滤波后的数据波动要比中值均值滤波要小得多，而且卡尔曼滤波的快速性和快速响应性远大于中值均值滤波。然而，如果时差信号快速变化，则卡尔曼滤波并不能获得良好的滤波效果。目前已有人提出了卡尔曼滤波算法与算术平均算法结合的数据滤波算法并投入了使用，对测得的时差信号Δt＝{[Δt(iT],i＝0,1,2,…}进行滤波，如图8所示。
一种基于卡尔曼滤波算法与算术平均算法结合的时差超声波流量测量方法，其步骤包括，
步骤(1)、时差信号测量：
当管道中有水流时，首先激发超声波换能器，使得超声波信号在管道内沿着水流方向顺流和逆流传播，超声波在流体中传播时，超声波传播的速度是超声波在流体中的波速与流体流速的叠加。时差法超声波测量的原理是：超声波信号在流体顺流与逆流传播时的传播速度不同，进而产生传播时间差值。由于时间差值的大小与流速大小有关，因此，通过测量时间差值的大小，就可以测量流体的流速，如图1所示；
图1中，假设向右为正向，换能器P1为顺流换能器，P2为逆流换能器；换能器发射出的超声波与流体流向垂直；换能器与反射片的距离为s，超声波在水中的速度为c，反射片与水平面夹角为45°，两反射片的中心距离为L，管道直径为D，水的正向流速即线速度vl为v，那么，由P1发射、P2接收超声波的顺流时间为：
Figure BDA0002804955810000131
由P2发射、P1接收超声波的逆流时间为：
Figure BDA0002804955810000132
根据式(1)和式(2)，可求得顺逆流时间差Δt：
Figure BDA0002804955810000133
在该步骤中，利用TDC－GP22高精度计时芯片，测量超声波换能器的时差信号Δt＝{[Δt(iT],i＝0,1,2,…}；
步骤(2)、时差信号的滤波：
该滤波方法具体如下：
(2.1)快速测量八组流量数据并储存；
(2.2)判断其最大值与最小值的差值；
当超声波水表检测到有流速增加大于或等于一定阈值δ时，此时认为水的流量处于快速变化阶段，提高超声波水表采样频率Fs，同时进入算术平均滤波算法阶段；
当检测到流速增加值小于一定阈值范围δ内，此时认为水流量处于稳定状态，则减小采样频率Fs，并进入卡尔曼滤波算法阶段；
(2.2.1)初始条件计算：
Figure BDA0002804955810000134
(2.2.2)利用k-1时刻的误差协方差Pk-1(+)、k-1时刻的状态转移矩阵Φk-1和k-1时刻的系统噪声Qk-1，计算k时刻的状态预测
Figure BDA0002804955810000135
和误差协方差预测Pk(-)，具体如下：
Figure BDA0002804955810000141
Figure BDA0002804955810000142
(2.2.3)利用k时刻的误差协方差预测Pk(-)、k时刻的测量系统的参数Hk和k时刻的测量噪声的协方差矩阵Rk，计算k时刻的卡尔曼增益
Figure BDA0002804955810000143
Figure BDA0002804955810000144
(2.2.4)利用k时刻的卡尔曼增益
Figure BDA0002804955810000145
和误差协方差预测Pk(-)，计算k时刻的误差协方差Pk(+)，更新矩阵如下：
Figure BDA0002804955810000146
(2.2.5)利用k时刻的卡尔曼增益
Figure BDA0002804955810000147
k时刻的状态预测
Figure BDA0002804955810000148
给定的初始估计值
Figure BDA0002804955810000149
和输入数据zk，递归计算k时刻的状态值
Figure BDA00028049558100001410
如下：
Figure BDA00028049558100001411
(2.3)按步骤(2.2.1)-(2.2.5)对式(6)中测得的时差信号zk进行滤波，得到滤波后的时差数据：
Figure BDA00028049558100001412
步骤(3)、累计流量的计算：
在计算瞬时流速和流量时，用到的是线平均流速，由流体力学可知，面平均流速根据流体的不同状态，有着不同的修正系数；
假设修正系数为K，计算流体的累计流量Q为：
Figure BDA00028049558100001413
在对比例中，步骤(2)中时差信号的滤波算法存在一些不足之处，首先，算术平均与卡尔曼滤波之间的切换条件不好确定(因为涉及到时差数据的采样率Fs、求最大和最小值时用到的采样数量、切换门限值δ等多个阴虚)；其次，算术平均滤波对于快速变化的时差值，滤波效果也不太理想，而且，卡尔曼滤波器对于快速变化的时差值，会跟不上时差的变化而导致较大的误差。因此，流量的累积计算可能会存在较大的误差。
为了说明现有超声波流量计中常规卡尔曼滤波器的不足之处，使用卡尔曼滤波器对变化的时差数据进行滤波。图4为平稳时差数据(真实时差＝2ns)时的滤波结果。图5给出了时差数据线性降低(其中包含一段抖动的数据)时的滤波结果，其中时差从4ns降至1ns。图6给出了时差数据存在快速变化时的滤波结果。由图5和图6可知，在时差存在变化时，常规卡尔曼滤波器跟不上时差的变化，滤波效果很差。
与对比例相比，本发明提出一种改进的卡尔曼滤波方法，并基于此给出超声波水表的流量计算方案。该滤波方法通过更新测量噪声和过程噪声的协方差矩阵，有效的解决了常规卡尔曼滤波方法，从而能快速跟上时差的变化，具有很好的滤波效果。而且，克服了现有卡尔曼滤波算法与算术平均算法结合的数据滤波算法中需要在两种滤波方法中切换的不足之处。这样，在超声波水流量计中使用该该改进的滤波方法，能有效地提高超声波流量计的精度。

Claims (2)

1.一种基于改进卡尔曼滤波器的时差超声波流量测量方法，其步骤包括，
步骤(1)、时差信号测量：
沿水流方向前后设置第一超声波换能器P1和第二超声波换能器P2；规定换能器与反射片的距离为s，超声波在水中的速度为c，反射片与水平面夹角为45°，两反射片的中心距离为L，管道直径为D，水的正向流速即线速度vl为v；
第一超声波换能器P1发出超声波脉冲，在第二超声波换能器P2上接收超声波脉冲，并计算接收到超声波时刻和发送超声波时刻的时间差t12，则有：
Figure FDA0002804955800000011
通过第二超声波换能器P2发出超声波脉冲，在第一超声波换能器P1上接收超声波脉冲，并计算接收到超声波时刻和发送超声波时刻的时间差t21，则有：
Figure FDA0002804955800000012
最后，获得时间差Δt＝t21-t12，根据式(1)和式(2)，可求得顺逆流时间差Δt：
Figure FDA0002804955800000013
步骤(2)、时差信号的滤波：
滤波处理是为了去掉原始数据中的随机误差，在超声波流量测量装置中，由于测量过程会引入外部噪声，使得测量结果与实际流速值产生一定偏差值，偏差值的大小与流体中是否有气泡或颗粒物、超声波信号质量、电路板PCB布局、温度补偿超声波声速值误差、外部电磁干扰等有关，为了减少外部噪声对超声波水表的干扰，提高水表的精度，需要对测得的时差数据进行滤波；
用xk表示k时刻的时差Δt(kT]的状态值，其中T为采样周期，具体如下：
xk＝Δt(kT]；
在卡尔曼滤波器中，系统状态方程如下：
xk＝Φk-1xk-1+wk-1， (4)；
Figure FDA0002804955800000021
式中，xk为k时刻的系统状态，xk-1为k-1时刻的系统状态，Φk-1为k-1时刻的状态转移矩阵，wk-1为k-1时刻的系统噪声，Qk为k时刻的系统噪声；
用zk表示时刻k时时差Δt(kT]的测量值，则时差Δt的测量方程如下：
zk＝Hkxk+vk， (6)；
Figure FDA0002804955800000022
式中，Hk表示k时刻的测量系统的参数，vk为k时刻的测量噪声，Rk为k时刻的测量噪声的协方差矩阵；
在系统中定义Φk-1＝1、Hk＝1，那么，时差数据的改进卡尔曼滤波方法如下；
(2.1)、初始条件计算：
Figure FDA0002804955800000023
式中，E[···]表示数学期望运算，
Figure FDA0002804955800000024
表示0时刻的状态值，P0(+)表示0时刻的误差协方差；
(2.2)、利用k-1时刻的误差协方差Pk-1(+)、k-1时刻的状态转移矩阵Φk-1和k-1时刻的系统噪声Qk-1，计算k时刻的状态预测
Figure FDA0002804955800000025
和误差协方差预测Pk(-)，具体如下：
Figure FDA0002804955800000026
Figure FDA0002804955800000027
(2.3)、利用k时刻的误差协方差预测Pk(-)、k时刻的测量系统的参数Hk和k时刻的测量噪声的协方差矩阵Rk，计算k时刻的卡尔曼增益
Figure FDA0002804955800000028
Figure FDA0002804955800000029
(2.4)、利用k时刻的卡尔曼增益
Figure FDA00028049558000000210
和误差协方差预测Pk(-)，计算k时刻的误差协方差Pk(+)，更新矩阵如下：
Figure FDA0002804955800000031
(2.5)、利用k时刻的卡尔曼增益
Figure FDA0002804955800000032
k时刻的状态预测
Figure FDA0002804955800000033
给定的初始估计值
Figure FDA0002804955800000034
和输入数据zk，递归计算k时刻的状态值
Figure FDA0002804955800000035
如下：
Figure FDA0002804955800000036
(2.6)、k时刻的测量噪声协方差Rk和k-1时刻的过程噪声协方差Qk-1更新如下：
Figure FDA0002804955800000037
Figure FDA0002804955800000038
Figure FDA0002804955800000039
Figure FDA00028049558000000310
Ratio＝R1/Q0 (18)；
Rk＝Rk(+)/Ratio (19)；
Qk-1＝Qk-1(+)/Ratio (20)；
(2.7)、将步骤(2.6)更新后的k时刻的测量噪声协方差Rk和k-1时刻的过程噪声协方差Qk-1代入步骤(2.2)-(2.6)中，递归计算k时刻的状态值
Figure FDA00028049558000000311
(2.8)、按步骤(2.2)-(2.7)对式(6)中测得的时差信号zk进行滤波，通过滤波，得到滤波后的时差数据：
Figure FDA00028049558000000312
步骤(3)、累计流量的计算：
在计算瞬时流速和流量时，用到的是线平均流速，由流体力学可知，面平均流速根据流体的不同状态，有着不同的修正系数；
设修正系数为K，计算流体的累计流量Q为：
Figure FDA0002804955800000041
式中，vl为流体线速度，根据式(3)得到，vl值为：
Figure FDA0002804955800000042
2.如权利要求1所述测量方法，其特征在于，步骤(3)前还有步骤(2.9)超声波声速的温度补偿：
(2.9.1)通过测量声速与温度的关系，得到基于曲线拟合得到声速与温度的关系曲线；
(2.9.2)制作超声波速度与温度的关系表格，并利用表格进行基于查表法的温度补偿，利用查表法得到的声速c，结合计算获得的滤波后的时差信号
Figure FDA0002804955800000043
计算该温度下流体线速度校正值vl'：
Figure FDA0002804955800000044
CN202011364228.0A 2020-11-27 2020-11-27 一种基于改进卡尔曼滤波器的时差超声波流量测量方法 Active CN112484798B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011364228.0A CN112484798B (zh) 2020-11-27 2020-11-27 一种基于改进卡尔曼滤波器的时差超声波流量测量方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011364228.0A CN112484798B (zh) 2020-11-27 2020-11-27 一种基于改进卡尔曼滤波器的时差超声波流量测量方法

Publications (2)

Publication Number Publication Date
CN112484798A true CN112484798A (zh) 2021-03-12
CN112484798B CN112484798B (zh) 2024-03-19

Family

ID=74936678

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011364228.0A Active CN112484798B (zh) 2020-11-27 2020-11-27 一种基于改进卡尔曼滤波器的时差超声波流量测量方法

Country Status (1)

Country Link
CN (1) CN112484798B (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113155210A (zh) * 2021-03-30 2021-07-23 合肥工业大学 一种基于改进卡尔曼滤波的电磁流量测量系统信号处理方法
CN113340363A (zh) * 2021-06-09 2021-09-03 河北大学 一种高精度宽量程的超声流量测量装置及测量方法
CN113865659A (zh) * 2021-09-17 2021-12-31 湖南威铭能源科技有限公司 一种高精度户用超声波水表抗干扰滤波方法及设备
CN114184724A (zh) * 2022-02-15 2022-03-15 华谱科仪(北京)科技有限公司 一种色谱仪载气流量补偿方法、装置及其存储介质
CN116642562A (zh) * 2023-07-27 2023-08-25 黑龙江惠达科技股份有限公司 一种植保无人机药液质量测量系统、方法和无人机

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105181997A (zh) * 2015-08-20 2015-12-23 天津市众中科技发展有限公司 非接触式超声波水流流速计及非接触式流速检测方法
WO2017004887A1 (zh) * 2015-07-03 2017-01-12 中国矿业大学 一种时差式超声波流量测量方法及装置
CN110346004A (zh) * 2019-08-16 2019-10-18 杭州山科智能科技股份有限公司 一种双声道超声时差法的流量测量数据融合方法
CN110609973A (zh) * 2019-08-27 2019-12-24 广东艾科技术股份有限公司 一种用于流量测量的卡尔曼滤波方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2017004887A1 (zh) * 2015-07-03 2017-01-12 中国矿业大学 一种时差式超声波流量测量方法及装置
CN105181997A (zh) * 2015-08-20 2015-12-23 天津市众中科技发展有限公司 非接触式超声波水流流速计及非接触式流速检测方法
CN110346004A (zh) * 2019-08-16 2019-10-18 杭州山科智能科技股份有限公司 一种双声道超声时差法的流量测量数据融合方法
CN110609973A (zh) * 2019-08-27 2019-12-24 广东艾科技术股份有限公司 一种用于流量测量的卡尔曼滤波方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
李芳;李艾华;吴朝军;: "超声波流量测量时差改进卡尔曼滤波", 计算机工程与设计, no. 18, 28 September 2009 (2009-09-28) *
邵欣;檀盼龙;高杰;王崇超;杨彬;: "超声波流量计的测量噪声滤波方法研究", 信息通信, no. 11, 15 November 2016 (2016-11-15) *

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113155210A (zh) * 2021-03-30 2021-07-23 合肥工业大学 一种基于改进卡尔曼滤波的电磁流量测量系统信号处理方法
CN113340363A (zh) * 2021-06-09 2021-09-03 河北大学 一种高精度宽量程的超声流量测量装置及测量方法
CN113865659A (zh) * 2021-09-17 2021-12-31 湖南威铭能源科技有限公司 一种高精度户用超声波水表抗干扰滤波方法及设备
CN113865659B (zh) * 2021-09-17 2024-03-22 湖南威铭能源科技有限公司 一种高精度户用超声波水表抗干扰滤波方法及设备
CN114184724A (zh) * 2022-02-15 2022-03-15 华谱科仪(北京)科技有限公司 一种色谱仪载气流量补偿方法、装置及其存储介质
CN114184724B (zh) * 2022-02-15 2022-05-17 华谱科仪(北京)科技有限公司 一种色谱仪载气流量补偿方法、装置及其存储介质
CN116642562A (zh) * 2023-07-27 2023-08-25 黑龙江惠达科技股份有限公司 一种植保无人机药液质量测量系统、方法和无人机
CN116642562B (zh) * 2023-07-27 2023-10-20 黑龙江惠达科技股份有限公司 一种植保无人机药液质量测量系统、方法和无人机

Also Published As

Publication number Publication date
CN112484798B (zh) 2024-03-19

Similar Documents

Publication Publication Date Title
CN112484798B (zh) 一种基于改进卡尔曼滤波器的时差超声波流量测量方法
AU2015294885B2 (en) Improved signal travel time flow meter
Zhou et al. Multipath ultrasonic gas flow-meter based on multiple reference waves
CN107727176B (zh) 一种超声波气体流量计及其流量测量方法
CN106643939A (zh) 用于超声波流量计计算超声波传播时间的方法
CN107860430B (zh) 一种基于时差法的超声波气体流量计时间差测量方法
CN107655533B (zh) 一种基于回波能量积分的气体超声波流量计信号处理方法和系统
CN115824331B (zh) 一种适用于超声水表的低功耗测量方法
CN111238587A (zh) 一种超声波换能器传输时间校准方法
JP2007187506A (ja) 超音波流量計
CN105115552B (zh) 基于侧音测相的超声波流量测量方法及装置
CN214583449U (zh) 一种高精度宽量程的超声流量测量装置
CN107576964B (zh) 线性变频信号的回波时间测量方法
US11137276B1 (en) All digital travel time flow meter using time reversed acoustics
CN112304376B (zh) 基于数据融合的超声波流量计流量测量方法
CN112785052B (zh) 基于粒子滤波算法的风速风向预测方法
CN111756353A (zh) 一种基于非线性融合滤波的液位仪噪声优化方法
Ma et al. Signal processing method based on connection fitting of echo peak point with a large slope for ultrasonic gas flow meter
JP2003279396A (ja) 超音波流量計
CN117629382B (zh) 一种水下声速和地形距离测量方法
CN116338240B (zh) 一种基于抛物线拟合的超声波液体流速测量方法及测量装置
EP4067833B1 (en) All digital travel time flow meter using time reversed acoustics
CN113191307B (zh) 一种基于贝叶斯估计的浅海波导简正波阶数判定方法
CN112763023B (zh) 基于优化数据模型的雷达物位计高精度测量输出处理方法
Catak et al. Self-Calibration of Ultrasonic Water Flow Meter

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
CB02 Change of applicant information

Address after: 330095 Gaohang building, 899 Huoju 5th Road, high tech Development Zone, Nanchang City, Jiangxi Province

Applicant after: Zhongke Shuiyan (Jiangxi) Technology Co.,Ltd.

Address before: 330095 Gaohang building, 899 Huoju 5th Road, hi tech Development Zone, Nanchang City, Jiangxi Province

Applicant before: Jiangxi Zhongke intelligent aquaculture research Co.,Ltd.

CB02 Change of applicant information
GR01 Patent grant
GR01 Patent grant