*Letter* **0.5 V Fifth-Order Butterworth Low-Pass Filter Using Multiple-Input OTA for ECG Applications**

#### **Montree Kumngern 1, Nattharinee Aupithak 1, Fabian Khateb 2,3,\* and Tomasz Kulej <sup>4</sup>**


Received: 15 November 2020; Accepted: 14 December 2020; Published: 21 December 2020

**Abstract:** This paper presents a 0.5 V fifth-order Butterworth low-pass filter based on multiple-input operational transconductance amplifiers (OTA). The filter is designed for electrocardiogram (ECG) acquisition systems and operates in the subthreshold region with nano-watt power consumption. The used multiple-input technique simplifies the overall structure of the OTA and reduces the number of active elements needed to realize the filter. The filter was designed and simulated in the Cadence environment using a 0.18 μm Complementary Metal Oxide Semiconductor (CMOS) process from Taiwan Semiconductor Manufacturing Company (TSMC). Simulation results show that the filter has a bandwidth of 250 Hz, a power consumption of 34.65 nW, a dynamic range of 63.24 dB, attaining a figure-of-merit of 0.0191 pJ. The corner (process, voltage, temperature: PVT) and Monte Carlo (MC) analyses are included to prove the robustness of the filter.

**Keywords:** fifth-order low-pass filter; operational transconductance amplifier; multiple-input bulk-driven technique; subthreshold region; nanopower

#### **1. Introduction**

Continuous-time filters are widely used in biomedical systems devoted to applications in electroencephalographic (EEG), electromyographic (EMG), and electrocardiographic (ECG) systems. The biological signals processed in these systems typically occupy the frequency range of 0.05–250 Hz, with an amplitude of 15 μV–5 mV [1]. In more detail, the frequency/amplitude ranges for EEG, EMG, and ECG signals are 0.05–60 Hz/15−100 μV, 10−200 Hz/0.1−5 mV, and 0.05−250 Hz/100 μV−5 mV, respectively. Figure 1 shows a typical data acquisition system for ECG signal processing. The pre-amplifier stage amplifies a low-amplitude ECG signal, then the low-pass filter selects the frequency range and eliminates out-of-band noise. The filtered analog signal is converted into digital form by an analog-to-digital converter (ADC) and then it is further processed by a digital signal processing (DSP) block. This work focused on the design of a low-pass filter with the cutoff frequency of 250 Hz. The analog low-pass filters for ECG acquisition systems should be designed to meet specific requirements, such as high dynamic range, low-power consumption, and small chip area. There are many low-pass filters for ECG acquisition systems described in the literature [2–10]. The Butterworth approximation is usually used because it provides a better linear phase and flat response within each bandwidth. Considering the analog filters in [2–10], one can distinguish two main techniques that

have been used to realize the low-pass Butterworth filters: the cascade approach [2–6] and the ladder simulation approach [7–10]. The cascade structure can be obtained by cascading several biquad filters, which leads to a simple and easy-to-tune realization.

**Figure 1.** Electrocardiogram acquisition system.

The present work focused on the second approach, i.e., the ladder simulation of a prototype filter. In particular, we designed a fifth-order low-pass Butterworth filter based on the RLC prototype shown in Figure 2. As it is widely known, the high-order filters based on the RLC prototypes have lower pass-band sensitivity to the variation of passive elements, compared with that of the cascade approach.

**Figure 2.** Prototype of a fifth-order low-pass filter.

The fifth-order low-pass Butterworth filters derived from the LC ladder-type filter were reported in [7–10]. The fifth-order Butterworth low-pass filter using fully differential operational transconductance amplifiers (FD-OTAs) is shown in Figure 3a [7]. The floating inductors L2 and L4 are simulated using OTA-based gyrators. The resistors RS and RL are simulated using OTAs as well. It should be noted that the filter in [7] employs eleven FD-OTAs and consumes 453 nW of power. The number of active devices that are used to realize this fifth-order Butterworth filter can be reduced by using multiple-output fully differential OTA (MOFD-OTA) as shown in Figure 3b [8,9], or fully differential-difference transconductance (FDDA) (a multiple-input active device) as shown in Figure 3c [10]. The structures in [8,9] employ six MOFD-OTA while the structure in [10] employs five FDDAs and one OTA. The filter in [8] consumes 350 nW of power and offers a 49.9 dB dynamic range while the filter in [9] consumes 41 nW of power and offers a 61.2 dB dynamic range. The filter in [10] consumes 453 nW of power and offers a 50 dB dynamic range.

This paper proposes a fifth-order Butterworth low-pass filter based on multiple-input operational transconductance amplifiers. It is clearly shown that the number of active devices needed to realize the fifth-order low-pass filter can be reduced by using the multiple-input OTA and results in reducing the power consumption and the active chip area. A novel technique with a multiple-input gate-driven (MIGD) transistor is used to realize multiple-input OTA with an internal CMOS structure as simple as a conventional OTA, hence, no additional current branches or cascade connections of multiple OTAs is needed. Unlike the floating-gate technique, the multiple-input technique does not require any additional processing steps to eliminate the trapped charge effect on the isolated gate nor any auxiliary circuit. Another advantage is that the multiple-input gate-driven PorN-MOS transistors can be realized with any CMOS process. It is worth noting that the results presented in this work are based on pre-layout simulation and this work does not include the physical realization of the filter, nor the experimental testing in the context of ECG applications. However, the principle of multiple-input transistors, as multiple-input bulk-driven and multiple-input bulk-driven quasi-floating-gates, have been confirmed experimentally by Khateb et al. in previous

works [11–13]. The paper is organized as follows: Section 2 shows the principle of multiple-input gate-driven OTA and the filter design based on it, Section 3 the simulation results, and finally Section 4 the conclusion.

**Figure 3.** Fifth-order Butterworth low-pass filters, (**a**) FD-OTA-C filter [7], (**b**) MOFD-OTA-C filter [8,9], (**c**) FDDA-based filter [10].

#### **2. Fifth-Order Butterworth Low Pass Filter**

#### *2.1. Multiple-Input Gate-Driven OTA*

The active filter proposed in this work exploits multiple-input OTAs, which allows for simplifying its overall structure [14]. The multiple-input OTA is realized using a concept of a multiple-input MOS transistor. The symbol and CMOS realization of this element are shown in Figure 4a,b, respectively. As it is seen in Figure 4b, the multiple-input MOS can be seen as a connection of an "internal" MOS transistor and a voltage divider/analog summing circuit, composed of capacitances *CGi* (*i* = 1 ... N). The capacitors *CGi* are shunted by the large resistances *RLi*, which ensures proper biasing of the gate terminal of the internal MOS for DC. The large resistances can be realized using an anti-parallel connection of two minimum-size MOS transistors operating in a cutoff region, as shown in Figure 4b. The small-signal equivalent circuit of the resulting multiple-input MOS is shown in Figure 4c. Assuming 1/ω*CGi* -*RLi*, the gate potential *VG* is given by

$$V\_G = \sum\_{i=1}^{N} \frac{C\_{Gi}}{C\_{\Sigma}} V\_{mi} \tag{1}$$

where *C*<sup>Σ</sup> is the sum of the capacitances *CGi* and the input capacitance of an internal MOS seen from its gate terminal *Cin*:

$$\mathbb{C}\_{\Sigma} = \mathbb{C}\_{\text{in}} + \sum\_{i=1}^{N} \mathbb{C}\_{\text{Gi}} \tag{2}$$

**Figure 4.** MIGD MOS transistor, (**a**) symbol, (**b**) realization, (**c**) small-signal model.

Since the AC signal at the gate of the internal MOS transistor is attenuated by the capacitive divider, the transconductance of the multiple-input device seen from its *i*-th input, and operating in the subthreshold region, can be expressed as:

$$\mathbf{g}\_{mi} = \frac{I\_D}{n\_p L I\_T} \cdot \frac{\mathbf{C}\_{Gi}}{\mathbf{C}\_{\Sigma}} \tag{3}$$

where *ID* is the DC drain current, *np*. is the subthreshold slope, and *UT*. is the thermal potential. As it is seen from (3), the transconductance seen from the *i*-th input is equal to the transconductance of the internal MOS, multiplied by the voltage gain of the capacitive voltage divider.

The lower input transconductance *gmi* entails a lower intrinsic voltage gain of the multiple-input MOS, as well as an increased input-referred noise. Both parameters are degraded by the factor of *C*Σ/*CGi*. However, it is worth noting that the linear range for such a device is also increased by the factor of *C*Σ/*CGi*, therefore, its dynamic range (DR) remains the same as that of the internal MOS.

The multiple-input MOS transistors were used to design a multiple-input OTA. The symbol and CMOS realization of the circuit are shown in Figures 5 and 6, respectively. The multiple-input MOS transistors M1 and M2 were used to create a multiple-input differential pair, biased by the self-cascode current sources M7,7c and M8,8c. The drain currents of the input differential pair are transferred to the outputs (Io<sup>+</sup> and Io-) through the current mirrors composed of the self-cascode transistors M3/3c-M4/4c and M5,5c-M6,6c. The current mirrors are loaded with the self-cascode current sources M10,10c and M9,9c. Note that the tail node that supplies the differential pair in Figure 6 is drawn with two branches for esthetic reasons. The application of self-cascode connections in this design allows for an increase in the output resistance of the OTA, which entails increasing the DC voltage gain of this circuit. The transistors M9c-M11c form a simple common-mode feedback circuit (CMFB) circuit, which forces the output common-mode level to be equal to the reference potential VCM. All the transistors operate in a subthreshold triode region. If the common-mode level is increasing/decreasing, the channel resistances of M10C1,c2 are increasing/decreasing as well, thus lowering the currents flowing through M10 and M9, and consequently, decreasing/increasing the common-mode level to the desired value. The transistors M9c and M10c are divided into two devices, which makes the circuit insensitive to the output differential signals of the OTA, at least for small amplitudes of the signal. For larger amplitudes of the output

signals, one can observe nonlinear components of the drain currents ID9 and ID10, caused by the differential output voltage of the OTA. However, this nonlinear effect is not apparent at the differential output of OTA, since variation of ID9 and ID10 are identical. This effect, however, causes variation of the output common-mode level. Figure 7 illustrates the large signal transfer characteristics and the common-mode level variation for unloaded OTA in Figure 6 controlled with differential signals. Note, moderate nonlinear effects are caused by the nonlinear output conductance of the OTA rather than that of the CMFB. Variations of the common-mode output voltage are maintained at an acceptable level.

**Figure 5.** Symbol of a multiple-input operational transconductance amplifier (OTA).

**Figure 6.** CMOS implementation for an MIGD OTA.

**Figure 7.** Output differential voltage and common-mode level versus input differential voltage for unloaded OTA in Figure 6.

One can say that the applied CMFB has a simple structure and does not consume additional power from supply rails. On the other hand, it slightly limits the maximum output voltage swing due to nonzero voltage drops across transistors M9c–M11c and variations of the output common-mode level caused by differential signals. However, the negative effects can be maintained at an acceptable level.

Assuming 1/ω*CGi* -*RLi*, the differential output current of the OTA can be expressed as:

$$I\_{o+} - I\_{o-} = I\_B \tanh\left(\sum\_{i=1}^{N} \frac{V\_{+\text{ini}} - V\_{-\text{ini}}}{n\_p l I\_T} \cdot \frac{\mathbb{C}\_{Gi}}{\mathbb{C}\_{\Sigma}}\right) \tag{4}$$

where *IB* is the biasing current (it was assumed that *ID*<sup>7</sup> = *ID*<sup>8</sup> = *ID*11). From (4), the small-signal transconductance from *i*-th input is given by:

$$\mathcal{g}\_{mi} = \frac{I\_D}{n\_p \mathcal{U}\_T} \cdot \frac{\mathbb{C}\_{Gi}}{\mathbb{C}\_{\Sigma}} \tag{5}$$

The DC voltage gain of the OTA from the *i*-th input can be expressed as:

$$A\_{vd} = \mathcal{g}\_{mi} r\_{out} \tag{6}$$

where *rout* is the output resistance of the OTA, given by:

$$\mathcal{F}\_{\text{out}} \triangleq \mathcal{G}\_{\text{5}} \mathfrak{a}\_{\text{4}} \delta \mathcal{F}\_{\text{ds}4} \mathfrak{a}\_{\text{4}} \mathcal{F}\_{\text{ds}4} \mathfrak{a}\_{\text{,}\text{6c}} \left\| \mathcal{G}\_{\text{m}9} \mathfrak{a}\_{\text{,}10} \mathcal{F}\_{\text{ds}9} \mathfrak{a}\_{\text{,}10} \left( r\_{\text{ds}9,10 \text{c}} / 2 \right) \right\|\tag{7}$$

Thanks to the self cascode connections, the voltage gain of the OTA can be at an acceptable level, despite the lower transconductance of the input differential pair.

From (4), the third order harmonic distortion of the OTA for a sinusoidal signal applied to one pair of input terminals, while the other pairs are shorted to ground the AC signals, can be expressed as:

$$HD\_3 = \frac{1}{48} \left( \frac{V\_{+i} - V\_{-i}}{n\_P Ul\_T} \cdot \frac{C\_{Gi}}{C\_\Sigma} \right)^2 \tag{8}$$

Thus, as it is seen from (8), the input linear range is increased by the factor of *C*Σ/*CGi*, i.e., the voltage attenuation factor introduced by the input capacitive divider.

The input referred noise of the OTA, including both thermal and flicker noise components, can be expressed as:

ı

$$\overline{v\_{n\text{th}}^2} = 2 \left(\frac{\text{l}I\_T}{I\_B}\right)^2 \left(\frac{\text{C}\_\Sigma}{\text{C}\_{Gi}}\right)^2 \left[\overline{i\_{1,2}^2} + 2\overline{i\_{3-6}^2} + 2\overline{i\_{9,10}^2}\right] \tag{9}$$

where:

$$\overline{I\_{1,2}^2} = 2qI\_B + \frac{K\_{Fp} \left(\frac{I\_B}{UL}\right)^2}{f \mathbb{C}\_{OX} (WL)\_{1,2}} \tag{10}$$

$$\overline{\sigma\_{3-6}^2} = 4kT \mathcal{g}\_{d\text{-}6} - \left(1 + \frac{2}{3} \frac{\mathcal{g}\_{d\text{-}6} - \text{c}}{\mathcal{g}\_{m3-6}}\right) \frac{\left(\mathcal{g}\_{m3-6}r\_{d\text{-}6}\right)^2}{\left(1 + \mathcal{g}\_{m3-6}r\_{d\text{-}6}\right)^2} \\ \quad + \frac{1}{4} \cdot \frac{K\_{\text{Fn}}\left(\frac{l\_R}{l\Gamma\_{\text{F}}}\right)^2}{f \cdot \mathbb{C}\_{\text{OX}}\left[\left(\mathcal{WL}\_{\text{/3-}6\text{cf}}\right)\right]} \tag{11}$$

$$\overline{\sigma\_{10-14}^2} = 4kT g\_{d\Phi 9-10c} \Big( 1 + \frac{2}{3} \frac{g\_{d\Phi 9-10c}}{g\_{m\Phi -10}} \Big) \frac{\left(g\_{m\Phi -10r} r\_{d\Phi 9-10c}\right)^2}{\left(1 + g\_{m\Phi -10r} r\_{d\Phi 9-10c}\right)^2} \\ \quad + \frac{1}{4} \cdot \frac{K\_{Fp} \left(\frac{I\_B}{II\_T}\right)^2}{f \mathbb{C}\_{\rm OX} \left[\left(WL\right)\_{9, 10cff}\right]} \Big) \tag{12}$$

where *gds*<sup>9</sup>−10*<sup>c</sup>* = *gds*<sup>9</sup>−10*c*1//*gds*<sup>9</sup>−10*c*1, WLieff = (WLi·WLic)/(WLi + WLic), i = 3 ... 10, WL9,10c = WL9,10c1 + WL9-10c2, KFn and KFp are the flicker noise constants for n- and p-channel transistors, respectively, and COX is the oxide capacitance per unit area.

As it is easy to note from (9), the input referred noise is increased by the factor of *C*Σ/*CGi*, as compared with the input noise of a single-input OTA biased with the same current. However, if the multiple input OTA is realized with N identical OTAs, each biased with the current of *IB*/N, then the input transconductance from each input and the input referred noise would be the same as that for the proposed realization (see the Appendix A). Since the linear range in the proposed design is increased *C*Σ/*CGi* times, then the DR of the proposed solution is also increased in the same proportion. The improved DR can be considered as the most important advantage of the proposed approach. Note that a similar capacitive attenuation approach that increase the dynamic range of OTAs has been presented before [15].

#### *2.2. Proposed Filter*

The proposed fifth-order Butterworth low-pass filter is shown in Figure 8a. It was developed from the LC-ladder filter based on the OTA-C topology. Its signal flow graph is shown in Figure 8b, where τ<sup>1</sup> = *C*1/*gm*1, τ<sup>2</sup> = *C*2/*gm*2, τ<sup>3</sup> = *C*3/*gm*3, τ<sup>4</sup> = *C*4/*gm*4, and τ<sup>5</sup> = *C*5/*gm*5. The filter comprises five MIGD OTAs and five capacitors. The number of active devices is reduced from 6 to 5, as compared with [8–10], which allows for the reduction of the active area and power.

**Figure 8.** (**a**) Proposed fifth-order Butterworth low-pass filter, (**b**) signal flow graph.

Considering OTA0, OTA1 in Figure 3b and OTA0, FDDA1 in Figure 3c, it can be noted that these devices are used to realize a floating resistor [9]. In this work these components together with the capacitor C1 create a lossy integrator as shown in Figure 9a [8], Figure 9b [10]. The ideal transfer function of these circuits can be expressed as:

$$\frac{V\_{op1} - V\_{on1}}{V\_{ip} - V\_{in}} = \frac{g\_{m0} / g\_{m1}}{(sC\_1 / g\_{m1}) + 1} \tag{13}$$

It is evident that the circuits work as lossy integrators, where the voltage gain can be controlled by *gmo*. Usually, all transconductances are set to be equal for easy tuning. Figure 9c shows the lossy integrator based on the three-input OTA that is proposed in this paper. The ideal transfer function of the circuit in Figure 9c can be expressed as:

$$\frac{V\_{op1} - V\_{on1}}{V\_{ip} - V\_{in}} = \frac{1}{\left(\text{sC}\_1/\text{g}\_{m1}\right) + 1} \tag{14}$$

Thus, the circuit works as a lossy integrator with unity gain. Assuming that *gmo* = *gm*1, Equations (13) and (14) will be identical. Thus, it can be concluded that the OTA0 in Figure 3b,c can be removed by using multiple-input OTA. This application can only be realized using multiple-input OTA and it is not possible by using conventional OTA. It should be noted that only the parts mentioned above in Figure 9a of [8], Figure 9b of [10] are modified, the other parts (OTA2-5 or FDDA2-5) are not changed and the feedback connection is still similar to the filters in [8,10].

**Figure 9.** Lossy integrator, (**a**) circuit in [8], (**b**) circuit in [10], (**c**) proposed circuit.

#### **3. Results and Discussion**

The circuit was designed in the Cadence environment using a TSMC 0.18 μm CMOS process with a metal-insulator-metal (MIM) capacitor. The OTA with bias current IB = 3.3 nA consumes 8.25 nW under a 0.5 V supply voltage. The isolation between OTA inputs is assured by the large value resistance of the MOS transistor operating in a cutoff region. The input currents are well below 100 pA for input range rail-to rail.

The RLC filter in Figure 2 was designed for the cut-off frequency of 250 Hz. The prototype element values were chosen as follows: RS = RL = 1 Ω, C1 = C5 = 393.4 μF, C3 = 1.27 mF, and L2 = L4 = 1.03 mH. For the OTA-C filter C1 = C5 = 5.43 pF, C2 = C4 = 14.2 pF, C3 = 17.57 pF, and the bias current for each OTA was IB = 3.3 nA. Note that the bias current circuit serves to bias all OTAs hence the maximum power consumption of the filter is 34.65 nW. Figure 10 shows the frequency responses of the RLC and the proposed filter. The gain magnitude at low frequency was −6 dB and −6.4 dB and the cut-off frequency (*fc*) was 250.2 Hz and 250.4 for the RLC and OTA filters, respectively. Both curves are in good agreement up to −70 dB. Figure 11 shows the frequency response of the filter with different bias currents ranging from 0.1 nA to 3.3 nA while the *fc* was in the range of 17.11 Hz to 250.4 Hz. The tuning capability and the linear relation between *fc* and IB are demonstrated in Figure 12. The transient response of the filter for the input sine wave of Vinpp = 100 mV and 10-Hz frequency are illustrated in in Figure 13. The total harmonic distortion (THD) was 1%.

To check the influence of the process, voltage, and temperature (PVT) variations on the filter performance, the corner analysis was performed. TheMOS transistor corners (ss, sf, fs, ff), MIM capacitor corners (ss, ff), voltage supply corners (490 mV, 510 mV), and temperature corners (0 ◦C, 60 ◦C) were used. The variation of the gain was in the range of −7.2 dB to −6.13 dB while the variation of the cut-off frequency was in the range of 100.6 Hz to 326.7 Hz, as shown in Figure 14. Note that the temperature corner has the most effect of the variation of the frequency response since the circuit operates in a subthreshold region. However, since the circuit is proposed for biomedical applications it is expected

that the temperature variation will be less than the chosen temperature corners. Although the variation of the cut-off frequency is large, the needed value can be simply re-adjusted by the bias current. Note that the amplitudes of the bumps at low bias currents in Figure 11 and at higher frequencies in Figure 14 do not exceed 1.6 dB and do not affect stability of the circuit in a significant manner.

The Monte Carlo analysis with 200 runs was performed for the filter gain and cut-off frequency as shown in Figures 15 and 16, respectively. The mean value of the gain was −6.23 dB with standard deviation of 0.14 dB, while the mean value of the cut-off frequency was 251.7 Hz with standard deviation of 4.9 Hz. Figure 17 shows the output referred noise density of the filter. The integrated in-band noise between 0.1 Hz to 250 Hz shows that the output referred noise is 77 μVrms. Figure 18 shows the performance of the proposed filter in processing the ECG signal where (a) depicts the ECG signal with a distortion signal (5 mV/500 Hz) that was applied at the input of the filter and (b) depicts the filtered output signal.

**Figure 10.** The frequency response of the RLC and the proposed filter.

**Figure 11.** The frequency response of the proposed filter with different bias currents.

**Figure 12.** The cut-off frequency versus the bias current.

**Figure 13.** The transient response of the filter for input sine wave with Vinpp = 100 mV and 10 Hz.

**Figure 14.** The frequency response of the proposed filter under process, voltage and temperature (PVT) corners.

**Figure 15.** Monte Carlo simulation of the voltage gain.

**Figure 16.** Monte Carlo simulation of the cut-off frequency.

**Figure 17.** The output referred noise density of the proposed filter.

**Figure 18.** Transient response of the filter for ECG signal: (**a**) input; (**b**) output.

The summary and comparison between the proposed filter and some previous works are shown in Table 1. Only the fifth-order Butterworth low-pass filters simulated by the LC-ladder type filter and suitable for ECG signal acquisition [7–10] have been selected for comparison. From Table 1, it is clear that the proposed filter has a lower number of active devices, power consumption, and figure-of-merit (FOM). Finally, the FOM versus VDD of fifth-order low-pass filters are shown in Figure 19. Compared with the works in [7,8,10], the proposed filter offers clearly better FOM. The FOM is even slightly lower than the one in [9] with half the value of VDD. It is worth noting that the estimated chip area of 2-inputs and 3-inputs OTA based on the MIGD technique is increased by approximately 5% and 8%, respectively, compared to that of a single-input conventional OTA with the same transistor dimensions. This confirms the advantage of this technique of saving chip area. Note, a similar conclusion of this advantage based on experimental results is stated in [11]. The small chip area of the proposed filter is evident in Table 1 compared with that of [10] that used off-chip capacitors for filter realization.


**Table 1.** Performance comparison between the proposed filter and other fifth-order low-pass filters for ECG signal acquisition.

**Figure 19.** Figure-of-merit (FOM) against VDD of the fifth-order low-pass filters.

#### **4. Conclusions**

In this paper, a fifth-order Butterworth low-pass filter using multiple-input OTA was proposed. The design proves that the number of OTAs for realizing the fifth-order low-pass filter architecture can be reduced using multiple-input OTAs. This entails the reduction of both the power consumption and the active area. Comparison with other designs in the literature shows that the proposed structure is the most beneficial, regarding the number of active devices and power consumption. The proposed filter was simulated with a 0.18 μm CMOS process and supplied with 0.5 V, which entailed operation in a subthreshold region. Simulation results including PVT corner and Monte Carlo (MC) analyses confirmed the robustness of the design.

**Author Contributions:** Conceptualization, F.K. and M.K.; methodology, M.K. and T.K.; software, F.K.; validation, F.K., N.A., and M.K.; formal analysis, M.K. and T.K.; investigation, F.K., M.K., and T.K.; writing—original draft preparation, M.K. and F.K.; writing—review and editing, M.K., F.K., and T.K. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work was supported by King Mongkut's Institute of Technology Ladkrabang under grant KREF026201. For the research, the infrastructure of the SIX Center was used.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Appendix A**

Let us compare noise properties of a fully-differential OTA in Figure A1a, biased with current *IB*, and a multiple-input OTA composed of *n* identical OTAs of the same structure, but biased with currents of *IB*/*n* (Figure A1b). For simplicity, let us consider only the thermal noise.

The mean-square value of the output noise current of the OTA in Figure A1a, operating in a weak-inversion region can be expressed as:

$$
\overline{\mathbf{u}\_n^2} = 2qI\_B\mathbf{A} \tag{A1}
$$

where *q* is the electron charge and *A* is a constant depending on the particular structure of the OTA. Consequently, the input referred noise is given by:

$$\overline{v\_n^2} = \frac{2qI\_BA}{\mathcal{g}\_m^2} \tag{A2}$$

The output noise current of each OTA in Figure A1b is:

$$
\overline{n\_n^2} = 2q \frac{I\_B}{N} A \tag{A3}
$$

However, the total output noise current, equal to the sum of N output currents, is the same as that for the reference OTA in Figure A1a. If the total output noise is referred to one input, we obtain:

$$
\overline{\upsilon\_m^2} = \frac{\overline{\upsilon\_n^2}}{\left(\frac{\mathcal{S}\_m}{\mathcal{N}}\right)^2} = \frac{2qI\_B A}{\mathcal{S}\_m^2} \mathcal{N}^2 \tag{A4}
$$

Thus, the rms value of the input noise is given by:

$$\sqrt{\overline{v\_{\rm un}^2}} = N \sqrt{\frac{2qI\_B A}{\mathcal{g}\_{\rm un}^2}}\tag{A5}$$

If a noiseless passive voltage divider, with N inputs, and a voltage gain of 1/N from each input is added at the input of the OTA in Figure A1a, then the rms value of the i-th input referred noise voltage is given by the same equation, namely, the i-th input-referred noise is the same as that for the OTA in Figure A1b.

**Figure A1.** Single-input fully-differential OTA (**a**) and multiple-input fully-differential OTA (**b**).

If we define the dynamic range as the ratio of the maximum input rms voltage, limited by an assumed level of nonlinear distortion (*Vinmax*) to the *i*-th input referred noise, then for the multiple-input OTA in Figure A1b we have:

$$DR\_{mi} = \frac{V\_{inmax}}{\sqrt{v\_m^2}}\tag{A6}$$

Since in the subthreshold region the linear range of a differential pair does not depend on the biasing current, then for the OTA in Figure A1a, with an additional passive voltage divider, the linear range will be extended N times, and the DR will be:

$$DR\_{pl} = N \cdot \frac{V\_{immax}}{\sqrt{\upsilon\_m^2}} = N \cdot DR\_{mi} \tag{A7}$$

Hence, the dynamic range of the OTA with a passive, noiseless voltage divider at the input is N times as large as that for the OTA composed of N identical OTAs, biased with N-times lower current. Similar proof could be concluded for flicker noise, however, in such a case not only the total biasing current, but also the total areas of transistor channels should be equal for the two compared circuits, i.e., the transistor channel areas of each OTA in Figure A1b should be N times smaller than for that of the OTA in Figure A1a.

#### **References**


**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

## *Article* **Design of a High Precision Ultrasonic Gas Flowmeter**

#### **Jianfeng Chen, Kai Zhang \*, Leiyang Wang and Mingyue Yang**

College of Metrology & Measurement Engineering, China Jiliang University, Hangzhou 310018, China; p1802085205@cjlu.edu.cn (J.C.); p1702085249@cjlu.edu.cn (L.W.); s1902080445@cjlu.edu.cn (M.Y.)

**\*** Correspondence: zkzb3026@cjlu.edu.cn

Received: 23 June 2020; Accepted: 19 August 2020; Published: 26 August 2020

**Abstract:** Aiming at the problems of substantial pressure loss, small range ratio and contact measurement in traditional gas flowmeters, this paper designs a new type of data-filtering ultrasonic gas flowmeter. The flowmeter is composed of hardware circuits such as STM32F407 (ARM Cortex 32-bit microcontroller) main control chip and high-precision timing chip TDC-GP22 (time to digital converter). The software uses a new data-filtering algorithm combining Kalman filtering algorithm and arithmetic average algorithm to improve the measurement accuracy of ultrasonic gas flowmeter. Through experimental comparison, we find that the filtering algorithm effectively reduces the measurement error of the system. Within the flow range of 0.025–4 m3/h, the maximum relative error of the system measurement is 2.7404%, which meets the national standard for the measurement error of the 1.5-level instruments. Moreover, it reduces the zero-drift to about one half of the original, which significantly improves the stability of the system. The gas flowmeter has the characteristics of high accuracy, good stability, low power consumption, and the overall performance is significantly improved.

**Keywords:** ultrasonic gas flowmeter; the principle of time-difference method; data filtering; low-power measurement

#### **1. Introduction**

Nowadays, flow detection technology has been widely used in various fields of industrial production such as petroleum, chemical industry, energy, etc. [1,2]. As a flowmeter for detecting fluid, flowmeters can be divided into gas flowmeters and liquid flowmeters from the measurement medium. Traditional turbine flowmeters, vortex flowmeters, and orifice flowmeters play an important role in the measurement of liquid flow [3]. However, due to the particularity of the gas medium flow and the complexity of the problems encountered in the signal transmission in the gas, the severe issues in gas flow measurement still plague people [4]. In the field of gas measurement, most domestic gas meters currently used in the country are membrane gas meters. The metering principle of the traditional membrane gas meter is to introduce the gas into a metering chamber with a constant volume and discharge it when it is full. Through a specific transmission mechanism, the number of charges and exhaust cycles converts into a mass, which reflects in the counter of the gas meter on [5]. The measurement technology of membrane gas meter is mature, reliable measurement and stable quality. However, its interior is contact measurement. Due to long-term wear, the measurement accuracy, sensitivity and stability of the measuring element is reduced, and there are problems such as substantial pressure loss, small range ratio and difficult meter reading. It has been challenging to meet the needs of modern society and people's daily life [6].

With the development of electronic information and the Internet of things technology, ultrasonic flowmeters have become a new type of intelligent instrument. Ultrasonic flowmeters achieve remote data transmission in combination with the Internet of things technology, solving the problem of field meter reading [7]. Compared with the traditional membrane gas meter, the ultrasonic

gas flowmeter has the outstanding advantages of noncontact type, small pressure loss, wide range ratio, high accuracy, etc. [8].

In recent years, a large number of researchers have made significant achievements in the study of fluid pipeline structure and ultrasonic transducers [9,10]. In 2013, Zhao Xuesong conducted a flow field simulation on the pipeline of a mono ultrasonic flowmeter. The flow field characteristics of different shapes and sizes in the pipeline are studied, and the best tapered pipeline model is designed [11]. In 2015, Li Yuxi studied the influence of different installation methods of ultrasonic transducers on the measurement accuracy of flow meters and determined the optimal installation method of ultrasonic flow meters [12]. Moreover, the data filtering of the ultrasonic gas flowmeter mainly adopts the median filtering algorithm and wavelet algorithm to improve the measurement performance [13]. In 2018, Yao Ping used algorithms such as wavelet analysis and recursive average filtering to study how to improve the measurement accuracy of gas ultrasonic flowmeters in complex flow fields [14]. These research results have promoted the widespread application of ultrasonic gas flowmeters [15].

However, most of the existing ultrasonic gas flowmeters meet the measurement accuracy requirements in the high zone (From the transition flow Qt (including Qt) to the maximum flow Qmax: Qt ≤ Qi ≤ Qmax) and the accuracy requirements of small flow in low zone (From minimum flow Qmin to transition flow Qt: Qmin ≤ Qi < Qt) cannot be guaranteed [16]. To further improve the overall performance of the ultrasonic gas flowmeter, this paper designs a new data-filtering algorithm combining Kalman filtering and arithmetic average. The small flow in the low zone has especially high measurement accuracy, which significantly improves its measurement stability and can better meet the needs of practical applications.

#### **2. Method of Time-Di**ff**erence Measurement**

#### *Principle of Time Di*ff*erence Measurement*

In this paper, the ultrasonic gas flowmeter adopts the principle of time-difference method. The time-difference method is to indirectly obtain the average flow velocity of the fluid medium by measuring the transmission time interval of the ultrasonic signal in the fluid medium in the forward and reverse directions [17].

Figure 1 shows the measurement principle of the time-difference method [18]. Transducer Pup is a forward-flow transducer, and Pdn is a reverse-flow transducer. The fluid flows from left to right in Figure 1 as a positive flow direction. The installation angle of the transducers and the pipe is θ, the pipe diameter is *D*. The linear distance between the two transducers is *L*. The transmission speed of ultrasonic waves in a gas medium is *c*. The forward velocity of the gas is *v*. The upstream time transmitted by the Pup transducer and received by the Pdn transducer is [19]:

$$t\_{\text{hyp}} = \frac{L}{\text{c} + \upsilon \cos \theta} \tag{1}$$

The downstream time transmitted by the Pdn transducer and received by the Pup transducer is:

$$t\_{\text{down}} = \frac{L}{c - v \cos \theta} \tag{2}$$

From Equations (1) and (2), through the downstream time transmitted *tdown* minus the upstream time transmitted *tup*, the forward and reverse flow time difference Δ*t* of the ultrasonic signal transmission can be obtained.

$$
\Delta t = t\_{\text{donu}} - t\_{\text{up}} = \frac{2Lv\cos\theta}{c^2 - v^2\cos^2\theta} \tag{3}
$$

Because the sound velocity *c* of the ultrasonic wave is much greater than the flow velocity *v* of the gas, so *c*<sup>2</sup> much larger than *v*<sup>2</sup> cos2 θ, the Equation (3) can be approximated as:

$$v = \frac{\Delta t c^2}{2L \cos \theta} \tag{4}$$

The linear velocity *v* of the gas can be obtained by Equation (4). When measuring the instantaneous flow rate of gas, it is necessary to use the average surface velocity to calculate. According to the relevant knowledge of fluid mechanics, the average surface velocity has different correction coefficients according to different states of the fluid. Let the correction coefficient be *K*, then the instantaneous flow *Q* of the gas is:

$$Q = K \times S \times v = \frac{K\pi D^2}{4}v\tag{5}$$

where *S* is the cross-sectional area of the pipeline; *K* is the instrument factor, which is related to the Reynolds number of the fluid state.

It can be seen from Equations (4) and (5), under the condition of ensuring the accuracy of other fixed quantities, accurate measurement of the time difference Δ*t* of the forward and reverse flow is the key to ensure the system's metering accuracy [20].

**Figure 1.** Schematic diagram of the principle of time-difference method.

#### **3. Hardware Design and Signal-Processing Process**

#### *3.1. Overall Hardware Design*

Figure 2 shows the system hardware block diagram of the ultrasonic gas flowmeter. The system is mainly composed of a control chip, a timing chip, an excitation signal amplifying circuit, a receiving conditioning circuit, an LCD (liquid crystal display) module and a power module.

**Figure 2.** System hardware block diagram.

This system uses a 32-bit single-chip STM32F407 as the calculation and control module, which meets the requirements of high operation accuracy and low-power-consumption measurement. This system selects the TDC-GP22 high-precision timing chip as the time measurement module. Through the SPI (serial peripheral interface) communication mode, the single-chip microcomputer controls the pulse generator inside the TDC-GP22 chip to generate an excitation signal and amplified by the amplifier circuit, which acts on the transmitting transducer to send out the transmit waveform and start the measuring time. Moreover, transmitting in the pipeline for some time, the transmit waveform is received by the receiving transducer. After being amplified and filtered, the ultrasonic echo signal is sent to the threshold comparison circuit, and the zero-crossing generates a signal to stop timing. Finally, it sends to the TDC-GP22 timing chip to complete the time measurement. The TDC-GP22 sends an interrupt signal to the microcontroller after the measurement is completed. The microcontroller reads the time-difference measurement result and brings into the flow calculation formula to display and store the flow data.

#### *3.2. TDC-GP22 Circuit Generates an Excitation Signal*

The transit time of the ultrasonic gas flowmeter is at the nanosecond level; the TDC-GP22 can measure at 90-picosecond resolution in the single-precision mode and at 45-picosecond resolution in double precision mode. All meet our requirements for time-difference measurement accuracy. To reduce power consumption, we choose a single-precision measurement mode. The TDC-GP22 high-precision timing chip used in this paper integrates a pulse-generating circuit, which greatly simplifies the configuration of peripheral circuits. The internal logic gate delay is used to measure the time interval to ensure the accuracy of the measurement. Figure 3 shows the excitation signal generated by the TDC-GP22 pulse-generating circuit.

**Figure 3.** High-precision timing chip (TDC-GP22) circuit generates an excitation signal.

*3.3. The Amplified Excitation Signal by the Amplifier Circuit*

Figure 3 shows that the amplitude of the excitation signal is 3.3 V. In actual design—because the ultrasonic signal generated by low-amplitude excitation is attenuated severely in the gas medium—an excitation source with an amplitude of more than 15 V is generally required to excite the gas ultrasonic transducer. Therefore, we used the TPS61085 (Booster converter) boost chip to amplify the original signal. Figure 4 shows that the amplifier circuit boosts the 3.3 V voltage of the initial excitation signal to 18 V, which acts on the transmitting transducer to send out the transmit waveform.

**Figure 4.** Amplified excitation signal by the amplifier circuit.

#### *3.4. The Receiving Circuit Processes the Echo Signal*

After propagating through the fluid medium, the transmit waveform signal is received by the receiving transducer. Figure 5a shows the original echo signal generated by the receiving transducer has an amplitude of several tens of millivolts, and it contains some noise interference signals. To obtain a stable and pure echo signal, the received original signal needs to be conditioned [21]. Therefore, we use the low-noise and high-precision op amp OPA320 (precision operational amplifier) to amplify the weak echo signal containing noise and then pass the active band-pass filter circuit to filter the amplified high and low-frequency interference signals. Finally, a stable and pure echo signal is obtained. Figure 5b shows the echo signal obtained after amplification and filtering by the receiving circuit.

**Figure 5.** (**a**) Signal processed before the receiving circuit; (**b**) signal obtained by the receiving circuit.

#### *3.5. Threshold and Zero-Crossing Comparison Circuit Generates Stop Timing Signal*

To avoid some interference signals will also trigger the TDC-GP22 chip timing. We first use a threshold comparison circuit to set a threshold voltage. When the amplitude of the echo signal is higher than the threshold voltage, the enable pin of the zero-crossing comparison circuit is triggered and generates a stop-timing signal and sends it to the stop pin of TDC-GP22 to complete timing [22]. In this way, the interference before the useful signal is shielded by a combination of threshold and a zero-crossing comparison circuit. Figure 6 shows the stop-timing signal generated by the threshold and zero-crossing comparison circuit.

**Figure 6.** Threshold and zero-crossing comparison circuit generate the stop-timing signal.

#### **4. Software Design**

#### *4.1. Software System Design*

Figure 7 shows the flow chart of ultrasonic gas flowmeter measurement. When the measurement starts, the single-chip microcomputer needs to be initialized, including initializing the corresponding peripheral configuration and system clock and then configuring the registers of the TDC-GP22. The single-chip microcomputer sends an instruction to start measurement to the TDC-GP22 through the SPI communication method and simultaneously enters the time measurement program. After completing the upstream and downstream time measurement, it enters the time-difference processing module and then runs a filtering algorithm to process the time-measurement data. Finally, the single-chip microcomputer calculates the gas flow and send the data to the LCD. At the same time, the flow rate data are stored in the main control chip's internal flash to prevent data loss after power failure. After completing measurement, the single-chip microcomputer enters the low power STANDBY mode and waits for the timer interrupt to wake up for the next measuring.

The ultrasonic gas flowmeter in this study is battery-powered, so its system must use low-power measurement. The power consumption of the STM32F407 microcontroller in STANDBY mode is less than 2.5 μA and the static power consumption of the TDC-GP22 timing chip is the only 2.2 μA. To reduce the measurement power consumption of the system, the sampling frequency of the system is related to the change of the gas flow rate. When the gas flow rate is considered to be stable, the sampling period sets to 1 s and the rest of the time is in a dormant state. When the gas flow rate changes, the system will increase the sampling frequency to reflect the amount of gas change.

**Figure 7.** System software flow chart.

#### *4.2. Data Filtering*

The purpose of data filtering is to eliminate random errors in the original data. A reasonable data-filtering algorithm must be based on a thorough study of the original data, taking into account data volatility, data distribution characteristics and sudden changes [23,24]. In the measurement of ultrasonic gas flowmeters, measurement deviations are caused by disturbing factors such as random noise and piping structure design errors. For the arithmetic average filtering algorithm, the smoother the filtering effect, the higher the lag of the algorithm and the poor ability to suppress random errors that occur in the system. For the application of the Kalman filter algorithm, it can achieve a strong ability to suppress the arbitrary interference appearing in the order and the filter effect will be worse when the response to sudden changes is realized. Therefore, in the case of frequent changes of gas flow rate, how to make the ultrasonic gas flowmeter respond quickly to the amount of gas change and achieve accurate measurement is another test of the filtering algorithm.

To balance the smoothing effect of filtering and the timeliness of data processing. In the design of the algorithm, it is necessary to discriminate whether the measurement deviations are caused by the interference of external random factors or the design error of the internal structure of the measuring pipe, then adopts different algorithms for different situations. Therefore, this study proposes a new data-filtering algorithm combining a Kalman filtering algorithm and an arithmetic averaging algorithm.

#### 4.2.1. Kalman Filtering Algorithm

Kalman filtering is a time–domain filtering method that is suitable for recursive solving. It can process the data obtained at each sampling instant immediately and based on the state estimates before that moment. The recursive equation gives the new state at any time estimate. Inflow measurement systems, the entire process can be represented by the following discrete models [25]:

$$x(k) = Ax(k-1) + Bu(k-1) + w(k-1)\tag{6}$$

$$z(k) = \mathbb{C}x(k) + v(k) \tag{7}$$

where the *x*(*k*) is state variable and the gain matrix *A* is the state of *k* − 1 is linearly mapped to the state of *k* at the current time, matrix *B* represents the gain of the optional control input *u*(*k* − 1) and the random signal *w*(*k* − 1) is the process excitation noise, because the state of the process has no control inputs, *u* = 0, the matrix *C* represents the gain of the state variable *x*(*k*) against the measured variable *z*(*k*) and *v*(*k*) represents the observed noise [26].

The Kalman filtering algorithm is as follows:

1. Predict the current state:

$$
\hat{\mathfrak{x}}(k \mid k-1) = A(k-1) \cdot \hat{\mathfrak{x}}(k-1 \mid k-1) \tag{8}
$$

2. The covariance of prior estimation errors:

$$P(k|k-1) = A(k-1)P(k-1|k-1)A^T(k-1) + Q(k-1) \tag{9}$$

3. Gain calculation:

$$K(k) = P(k|k-1)\mathbb{C}^T(k)\left[\mathbb{C}^T(k)P(k|k-1)\mathbb{C}(k) + R(k)\right] \tag{10}$$

4. Status estimate update:

$$\hat{\mathfrak{x}}(k|k) = \hat{\mathfrak{x}}(k|k-1) + \mathcal{K}(k)[z(k) - \mathcal{C}(k)\hat{\mathfrak{x}}(k|k-1)] \tag{11}$$

#### 5. The covariance of posterior estimation error:

$$P(k|k) = [1 - K(k)C(k)]P(k|k-1) \tag{12}$$

where *x*ˆ(*k*|*k* − 1) is the prior estimate, *x*ˆ(*k*|*k*) is the posterior estimate and *P*(*k*|*k* − 1) is the covariance of the previous estimate error, *P*(*k*|*k*) is the covariance of the following estimation error, *K*(*k*) is the filter gain, *R*(*k*) is the observation noise covariance matrix, *Q*(*k* − 1) is the process excitation noise covariance matrix [27].

#### 4.2.2. Filter Algorithm Combining Kalman and Arithmetic Average

Figure 8 shows a data-filtering algorithm combining the Kalman filtering algorithm and the arithmetic averaging algorithm. First, the system quickly measures eight sets of time-difference data and determine the difference value between the maximum and minimum values. When the ultrasonic gas flowmeter detects that the difference value higher than or equal to a particular threshold δ, the data have a wide range at this time, the gas flow rate is considered to be in a rapidly changing stage due to external interference. Hence, the system increases the sampling frequency Fs to track the change of the gas flow rate. At this time, the system selects a Kalman filter algorithm to process the time-difference data. When the difference value is less than a particular threshold range δ, the range of data changes is small. The gas flow rate is considered to be relatively stable at this time. The sampling frequency Fs is reduced, and we select the arithmetic averaging algorithm to solve the problem that the internal structure design error of the actual pipeline affects the gas flow rate.

**Figure 8.** Data-filtering algorithm.

The data-filtering algorithm can effectively improve the measurement accuracy of the system. It can also significantly reduce the zero-drift of the system. Especially in the measurement of small flow, the stability and measurement accuracy of the system was improved.

#### **5. Results and Discussion**

#### *5.1. Experimental Device*

To test the effect of the comparison data-filtering algorithm, we use a gas flow measurement system for experiments. The schematic diagram of the experimental measurement system is shown in Figure 9, which is mainly composed of the standard flowmeter, the ultrasonic gas flowmeter, the pressure measurement sensor, the temperature measurement sensor, a flow regulating valve, a fan and an outlet valve.

**Figure 9.** Schematic diagram of the experimental measurement system.

In the experimental measurement, the ultrasonic gas flowmeter was installed on the pipe section of the verified flowmeter and the fan started to blow in air, and the air outlet valve was closed. After the pressure between the standard flowmeter and the verified flowmeter was consistent and the measurement system was stable. We opened the air outlet valve [28]. In the experiment, we adjusted the flow rate by adjusting the opening of the flow regulating valve. After running for some time, we read the flow data of the standard flowmeter and the verified flowmeter.

To verify the measurement performance of ultrasonic gas flow, it is necessary to test the instantaneous flow at different flow points. The measuring range designed is 0.025–4 m3/h in this study. According to "JJG-1030-2007 Verification regulation of ultrasonic flowmeters", the detection flow points are Qmin, Qt, 0.4 Qmax, Qmax, for each flow point test is no less than 3 times; the experimental results are shown in Table 1.


**Table 1.** Flow experiment data at different flow points.

#### *5.2. Zero-Drift Discussion*

Zero-drift means that the flow rate is zero under the test conditions: the pipe section is filled with air, and the inlet and outlet of the pipe are sealed. Therefore, the time-difference data of the zero-drift is continuously collected when the flow rate is zero.

Theoretically, when the flow rate of the fluid is zero, the transmission time of the ultrasonic signal in the forward and reverse directions is equal. Actually, it is disturbed by factors such as pipe structure, the performance difference of upstream and downstream transducers, which causes the time difference to not be zero, so the calculated flow rate is not zero, and it will cause errors in the actual flow measurement. Especially for the measurement of small flow, because the time difference of the small flow is close to the time difference of the zero-drift, the error of the measured flow result is greater. Therefore, the zero-drift phenomenon will affect stability and accuracy of the ultrasonic flowmeter measurement.

To compare the effect before and after the data-filtering algorithm, we sampled the original zero-drift time-difference data and the zero-drift time-difference data after data-filtering algorithm. We labeled them as Sample 1, Sample 2, respectively. As shown in Figure 10, the black line represents the original zero-drift data. It can be seen from figure that the zero-drift range of the system remains within ±10 ns. The thick red line represents the zero-drift data obtained after data filtering. The zero-drift data range of the system is within ±5 ns. By comparison, it can be seen that the scope of the zero-drift data becomes significantly smaller.

**Figure 10.** Comparison of zero-drift before and after data filtering.

To describe the degree of dispersion of the zero-drift data after the data filtering, the definition of sample standard deviation is introduced:

$$S = \sqrt{\frac{1}{n-1} \sum\_{i=1}^{n} \left(\mathbf{X}\_i - \overline{\mathbf{X}}\right)^2} \tag{13}$$

where *S* is the samples standard deviation, *Xi* is the time-difference value of the zero-drift, *X* is the average value of samples and *n* is the number of samples.

Table 2 analyses the difference value between the maximum and minimum value and standard deviation of the two groups of samples. The difference value indicates the range of the sample data [29]. The smaller the scale of the sample data, the more stable the data; and the smaller the standard deviation, the lower the degree of sample dispersion [30]. It can be seen from Table 1 that the difference value and standard deviation of Sample 2 are the smallest, and the difference value of Sample 1 is the largest. The results in Table 1 show that the data of Sample 2 are the most stable and have the least dispersion.



From the above analysis, it can be concluded that the data-filtering algorithm can effectively reduce the zero-drift of the system. The smaller the change range of the zero-point drift data, the more stable the measurement result of the ultrasonic flowmeter, and the less affected by the interference of the actual measurement environment.

#### *5.3. Experimental Data*

The range of the ultrasonic gas flowmeter is 0.025–4 m3/h, the transition flow Qt is 0.4 m3/h, the range ratio R-value is 160, the accuracy level is 1.5 and the measuring pipe diameter is nominally diameter 20-mm in diameter. To verify and compare the effectiveness of data-filtering algorithms, we select two flow points distributed in the low zone (0.025 m3/h <sup>≤</sup> Qi < 0.4 m3/h) and the high zone

(0.4 m3/<sup>h</sup> <sup>≤</sup> Qi <sup>≤</sup> 4 m3/h) and test multiple groups time-difference data. Three sets of data were selected for each flow point, as shown in Tables 3 and 4, respectively. The time-difference data of a particular flow point is chosen randomly and the effect of the data-filtering algorithm is tested. The comparison of the results before and after the data filtering is shown in Figure 11.


**Table 3.** Time-difference signal samples without data filtering.


**Table 4.** Time-difference signal samples with data filtering.

**Figure 11.** Comparison of time-difference signals before and after data filtering.

This study uses the time-difference method to measure the flow rate of an ultrasonic gas flow meter. Combining Equations (4) and (5), it can be seen that the ultrasonic velocity *c* in a static fluid is taken as a constant. The installation angle θ of the transducer and the pipe and the linear distance *L* between the two transducers are directly measurable physical quantities, so the accurate measurement of the time difference between the forward and reverse flow can improve the measurement accuracy of the system. To verify the effect of the filtering algorithm and eliminate the influence of other physical constants on the verification results, this study compares the relative error data of the three sets of time difference before and after filtering and improves the measurement accuracy of the system by accurately measuring the time difference between the forward and reverse flow.

Figure 11 is a comparison of the effect of time-difference signals before and after data-filtering algorithm at a particular flow point. The black line represents the signal before data filtering, and the thick red line represents the signal after data filtering. By comparison, we can find that

the time-difference signal after data filtering is smoother than the time-difference signal before data filtering. It shows that the filtered signal is less affected by the random interference.

#### *5.4. Data Discussion*

Tables 2 and 3 show three sets of time-difference data before and after data filtering. The formula for calculating relative errors in table is:

$$E = \frac{\Delta T}{T} \times 100\% \tag{14}$$

where *E* is the actual relative error, generally given as a percentage, Δ*T* is the absolute error, that is the difference value between the measured value and the real value, and *T* is the real value. In the actual calculation, the average value is used instead of the real value.

It can be seen from Table 3 that when the flow point is in the high zone, the relative error of time-difference data is small; when the flow point is in the low zone, the relative error of time-difference data is large, and it is difficult to achieve accurate measurement. It can be seen from Table 4 that when the flow point is in the high zone, the relative error is small. Especially the relative error of the low zone flow is significantly reduced, which proves that the relative error of time-difference data after data filtering is smaller, and the measurement accuracy is higher.

To better compare the time-difference data in Tables 3 and 4, the sample averages of the time-difference data of the low zone and high zone flow and the difference Δ*E* between the maximum and minimum of the relative error were calculated. The results are shown in Table 5.


**Table 5.** Performance comparison with and without data filtering.

By comparing the data in Table 5, it can be found that the relative error Δ*E* of a particular flow point in the low zone has decreased from 17.7258% before data filtering to 1.7101% after data filtering. The relative error Δ*E* of a particular flow point in the high zone decreased from 1.5303% before data filtering to 0.4260% after data filtering. The comparison shows that the data using the filtering algorithm is more stable and more accurate.

#### **6. Conclusions**

This paper designs a high-precision ultrasonic gas flow meter and analyses the signal-processing process of the ultrasonic gas flow meter. A data-filtering algorithm combining Kalman filtering and arithmetic averaging is proposed to improve the measurement accuracy and stability of the system.

The zero-drift of the ultrasonic flowmeter is an essential manifestation of the measurement performance of the entire system. The smaller the zero-drift is, the higher the measurement stability of the system is. Therefore, by comparing the effects of the data-filtering algorithm, it is evident that the zero-drift scope of the system is reduced from within ±10 ns to within ±5 ns. The zero-drift is about 1/2 of the original zero-drift of the system, which effectively improves the stability of the system measurement.

To verify the effectiveness of the data-filtering algorithm, it can be seen from the experimental comparison results that the relative error of the time-difference data after the data-filtering algorithm is reduced, especially the relative error of the small flow is reduced more obviously. Within the flow range of 0.025–4 m3/h, the maximum relative error of the system measurement is 2.7404%, which meets the national standard for the measurement error of the 1.5-level instruments. It shows that the measurement accuracy of the system is higher after data filtering, and its comprehensive measurement performance was significantly improved.

**Author Contributions:** J.C. did the research, wrote and edited the paper. L.W. provided article writing ideas and data collection. M.Y. completed signal processing. K.Z. provided the experiment conditions and gave many important suggestions on article. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by National Natural Science Foundation of China, Grant Number 11472260. **Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
