CN113375737A - Flow velocity metering method of time difference type ultrasonic gas flowmeter - Google Patents

Flow velocity metering method of time difference type ultrasonic gas flowmeter Download PDF

Info

Publication number
CN113375737A
CN113375737A CN202010543939.8A CN202010543939A CN113375737A CN 113375737 A CN113375737 A CN 113375737A CN 202010543939 A CN202010543939 A CN 202010543939A CN 113375737 A CN113375737 A CN 113375737A
Authority
CN
China
Prior art keywords
upstream
downstream
echo signal
time
signal
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
CN202010543939.8A
Other languages
Chinese (zh)
Other versions
CN113375737B (en
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.)
Zhengzhou University
Original Assignee
Zhengzhou University
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 Zhengzhou University filed Critical Zhengzhou University
Priority to CN202010543939.8A priority Critical patent/CN113375737B/en
Publication of CN113375737A publication Critical patent/CN113375737A/en
Application granted granted Critical
Publication of CN113375737B publication Critical patent/CN113375737B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01FMEASURING VOLUME, VOLUME FLOW, MASS FLOW OR LIQUID LEVEL; METERING BY VOLUME
    • G01F1/00Measuring the volume flow or mass flow of fluid or fluent solid material wherein the fluid passes through a meter in a continuous flow
    • G01F1/66Measuring 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
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01FMEASURING VOLUME, VOLUME FLOW, MASS FLOW OR LIQUID LEVEL; METERING BY VOLUME
    • G01F1/00Measuring the volume flow or mass flow of fluid or fluent solid material wherein the fluid passes through a meter in a continuous flow
    • G01F1/66Measuring 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/667Arrangements of transducers for ultrasonic flowmeters; Circuits for operating ultrasonic flowmeters
    • G01F1/668Compensating or correcting for variations in velocity of sound

Landscapes

  • Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Fluid Mechanics (AREA)
  • General Physics & Mathematics (AREA)
  • Measuring Volume Flow (AREA)

Abstract

The embodiment of the disclosure discloses a flow rate metering method of a time difference type ultrasonic gas flowmeter. One embodiment of the method comprises: acquiring attitude information of an upstream ultrasonic transducer and a downstream ultrasonic transducer; acquiring the propagation distance of an ultrasonic signal from an upstream ultrasonic transducer to a downstream ultrasonic transducer; determining a time difference through a cross-power spectral density function and a coherence function, wherein the time difference is the difference between the transit time of the upstream ultrasonic signal and the transit time of the downstream ultrasonic signal; analyzing the upstream ultrasonic signal and the downstream ultrasonic signal based on Hilbert transform, and determining a variable threshold; determining a transit time of the upstream ultrasonic signal and a transit time of the downstream ultrasonic signal based on a variable threshold; and determining the flow velocity of the medium detected by the flowmeter according to the propagation distance, the attitude information, the time difference, the transition time of the upstream ultrasonic signal and the transition time of the downstream ultrasonic signal. This embodiment enables an accurate calculation of the flow rate of the medium.

Description

Flow velocity metering method of time difference type ultrasonic gas flowmeter
Technical Field
The embodiment of the disclosure relates to the field of flowmeters, in particular to a flow rate metering method of a time difference type ultrasonic gas flowmeter.
Background
Ultrasonic gas flow meters have been increasingly used in many fields of industry and life. For example, ultrasonic gas flow meters can measure the flow of natural gas, chlorine in water treatment, and multi-component gases.
The existing ultrasonic gas flowmeter is usually calculated by adopting a time difference method. As an example, dual threshold comparison methods are commonly used to determine flow rate measurements for time-difference ultrasonic flow meters, which focus on calculating the exact time for a signal to pass a threshold point. However, the above method is susceptible to noise interference. It is difficult to identify when the starting point of the echo signal is low in amplitude or two very close maximum points. Thereby affecting the accuracy of the measured value.
As another example, a method of determining the position of the echo signal by setting a transformation ratio threshold and finding out the feature points may also be used, and then a threshold proportionality coefficient may be calculated according to the peak value of the echo signal and the corresponding peak value of the echo signal. However, when the flow rate changes rapidly, the change speed of the echo energy integral is easily interfered by noise in the gas, so that the change of the characteristic point has instability, and the accuracy of the measured value is further influenced.
As another example, the transit time of the ultrasonic signal may be calculated by using a cross-correlation method in the time domain, which can solve the above two problems. However, the cross-correlation method requires interpolation to improve accuracy, searches for a suitable reference signal, and has a large calculation amount.
Accordingly, there is a need in the art for a new method for determining the flow rate of a transit time ultrasonic gas meter.
Disclosure of Invention
This summary is provided to introduce a selection of concepts in a simplified form that are further described below in the detailed description. This summary is not intended to identify key features or essential features of the claimed subject matter, nor is it intended to be used to limit the scope of the claimed subject matter.
Some embodiments of the present disclosure provide a flow rate metering method, apparatus, electronic device and computer readable medium for a time difference ultrasonic gas meter to solve the technical problems mentioned in the background section above.
Some embodiments of the present disclosure provide a flow rate metering method of a transit time ultrasonic gas flowmeter, the method comprising: acquiring attitude information of an upstream ultrasonic transducer and a downstream ultrasonic transducer; acquiring the propagation distance of an ultrasonic signal from an upstream ultrasonic transducer to a downstream ultrasonic transducer; determining a time difference through a cross-power spectral density function and a coherence function, wherein the time difference is the difference between the transit time of the upstream ultrasonic signal and the transit time of the downstream ultrasonic signal; analyzing the upstream ultrasonic signal and the downstream ultrasonic signal based on Hilbert transform to determine a variable threshold; determining a transit time of the upstream ultrasonic signal and a transit time of the downstream ultrasonic signal based on the variable threshold; and determining the flow velocity of the medium detected by the flow meter based on the propagation distance, the attitude information, the time difference, and the transit times of the upstream ultrasonic signal and the downstream ultrasonic signal.
In some embodiments, the flow rate is determined by the following equation:
Figure BDA0002540014280000021
Figure BDA0002540014280000022
wherein v represents a flow rate of the medium detected by the flow meter; l represents the propagation distance; theta represents the included angle between the upstream ultrasonic transducer and the pipe wall of the carried medium and the included angle between the downstream ultrasonic transducer and the pipe wall of the carried medium; Δ T represents a time difference; t isuRepresenting a transit time of the upstream ultrasonic signal; t isdRepresenting the transit time of the downstream ultrasonic signal.
In some embodiments, the determining the time difference by the cross-power spectral density function and the coherence function includes: acquiring an upstream echo signal and a downstream echo signal; determining the correlation degree of the upstream echo signal and the downstream echo signal in the time domain through a cross-correlation function; the cross-correlation function of the upstream echo signal and the downstream echo signal in the time domain is:
Figure BDA0002540014280000023
wherein t represents the selected time; τ represents a time delay of the upstream echo signal and the downstream echo signal; indicating a degree of correlation between the upstream echo signal and the downstream echo signal at a time of a time delay τ; x (t) represents the amplitude of the upstream echo signal or the downstream echo signal at time t; y (t + τ) represents the amplitude of the downstream echo signal or the upstream echo signal at time t + τ; t represents the period of sampling; for discrete time, the above correlation function translates into:
Figure BDA0002540014280000024
wherein m represents the selected discrete time; n represents a discrete time delay between the downstream echo signal and the upstream echo signal; n denotes the period of discrete sampling.
In some embodiments, the determining the time difference by the cross-power spectral density function and the coherence function includes: the cross-power spectral density function of the upstream echo signal and the downstream echo signal is:
Figure BDA0002540014280000031
wherein e is-j2πωτRepresenting a complex variable function; the cross-power spectral density function is expressed in complex form, with complex polar coordinates as follows:
Figure BDA0002540014280000032
wherein, Im [ P ]xy(ω)]An imaginary number value part representing a cross-power spectral density function; re [ P ]xy(ω)]Real number representing cross-power spectral density functionA value portion; j represents an imaginary number; the correlation between the upstream echo signal and the downstream echo signal in the frequency domain is determined by the following formula:
Figure BDA0002540014280000033
wherein, γ2 xy(ω) is a function of the power spectral density of the amplitude-squared coherence, which represents the correlation between the upstream echo signal and the downstream echo signal in the frequency domain, also called correlation coefficient, and the maximum value is 1; px(ω) represents the self-power spectral density function of x; py(ω) represents the self-power spectral density function of y; wherein, wmIs shown as
Figure BDA0002540014280000034
When the value is maximum, the value of w is obtained; omegamAlso indicates the signal frequency when the phase difference between the upstream echo signal and the downstream echo signal is thetaxy(wm) According to the relationship between the cross-power spectrum phase and the frequency, the phase difference can be obtained by a trigonometric function relationship;
Figure BDA0002540014280000035
the time difference is determined according to the following formula:
Figure BDA0002540014280000036
in some embodiments, the determining the transit time of the upstream ultrasonic signal and the transit time of the downstream ultrasonic signal based on the variable threshold includes: the above variable threshold is determined by the following formula: δ ═ 0.1max (a (t)); wherein δ represents a variable threshold; a (t) represents an envelope of the upstream echo signal or the downstream echo signal; wherein a (t) is determined by the following formula:
Figure BDA0002540014280000037
where x (t) represents the upstream or downstream signal: the above x (t) is determined by the following formula:
Figure BDA0002540014280000038
wherein the content of the first and second substances,
Figure BDA0002540014280000039
represents the convolution of x (t) and 1/π t, with the formula:
Figure BDA00025400142800000310
wherein, H [ x (t)]Denotes the Hilbert transform of x (t).
In some embodiments, the determining the transit time of the upstream ultrasonic signal and the transit time of the downstream ultrasonic signal based on the variable threshold includes: determining the maximum value of the envelope rising stage of the upstream echo signal or the downstream echo signal; the local maximum point is determined by the following formula: nn0O + δ; wherein nn0Representing a local maximum point Z0The ordinate of (a); δ represents a variable threshold; o represents the minimum difference between the variable threshold and the vertical coordinate of the local maximum; wherein said O is determined by the formula: o ═ min (| δ -nn)iL) wherein nniRepresents the ith local maximum ZiThe ordinate of (a); by aligning the local maximum point Z0Two adjacent sampling points S-1,S+1Performing parabolic interpolation processing on the three points to obtain a peak value; determining the peak value as an upstream feature point or a downstream feature point Sp(mp,np) (ii) a Determining a transit time of the upstream echo signal or the downstream echo signal according to the following formula: t isu=Tp+Tu’;Td=Tp+Td'; wherein, TPA fixed time obtained according to an experiment from the time when the upstream ultrasonic signal or the downstream ultrasonic signal is transmitted to the time when the upstream echo signal or the downstream echo signal is received; t isu' represents the time from the sampling start of the upstream echo signal to the upstream feature point; t isd' represents the time from the sampling start of the downstream echo signal to the downstream feature point; wherein, the above-mentioned Td' is determined by the following equation: t isd’=mp1/Fs;Wherein m isp1Indicating a downstream feature point Sp1The abscissa of (a); fs represents the sampling rate; wherein Tu' is determined by the following formula: tu ═ mp2/Fs, wherein mp2Represents the upstream feature point Sp2The abscissa of (a).
One of the above-described various embodiments of the present disclosure has the following advantageous effects: firstly, the time delay is directly calculated by adopting a cross-power spectrum method, the time difference of the upstream and downstream flight signals is determined, and compared with a related algorithm for respectively calculating the time of the upstream and downstream flight signals, the calculation amount is reduced. Meanwhile, the mutual power method is an unstable process, and nanosecond-level accurate time delay can be directly obtained. And further improves the accuracy of the flow velocity calculation.
In addition, the variable threshold is determined by adopting Hilbert transform, different thresholds are set according to waveforms of different signals, and the position of the echo signal is further accurately determined. Thereby improving the accuracy of the method.
Drawings
The above and other features, advantages and aspects of various embodiments of the present disclosure will become more apparent by referring to the following detailed description when taken in conjunction with the accompanying drawings. Throughout the drawings, the same or similar reference numbers refer to the same or similar elements. It should be understood that the drawings are schematic and that elements and features are not necessarily drawn to scale.
FIG. 1 is a schematic representation of a reflective structure of a transit time ultrasonic gas meter according to some embodiments of the present disclosure;
FIG. 2 is a flow chart of some embodiments of a flow rate metering method of a transit time ultrasonic gas meter according to some embodiments of the present disclosure;
FIG. 3 is a schematic illustration of upstream and downstream ultrasonic transmit signals and upstream and downstream echo signals of the present disclosure;
FIGS. 4 and 5 are cross-power spectral phase diagrams of an upstream echo signal and a downstream echo signal, respectively;
FIGS. 6 and 7 are schematic diagrams of the correlation degree of the upstream echo signal and the downstream echo signal, respectively;
FIG. 8 is a schematic diagram of a first echo signal or a second echo signal;
FIG. 9 is a schematic diagram of a functional block diagram of portions of the ultrasonic flow meter;
FIG. 10 is a time of flight T of an upstream echo signaluTransit time T of the downstream signaldA data table of the relationship with the time difference Δ T at 9 flow points;
FIG. 11 is a graphical illustration of the approximate linear relationship of time-of-flight difference to standard flow point;
FIG. 12 is a schematic illustration of an ultrasonic flow meter calibration test;
fig. 13 and 14 are graphs showing the results of zero drift experiments performed at high temperatures of 25.9 c and 70.0 c, respectively, for a filter with a moving average order of 30;
FIG. 15 is a schematic diagram of the error of K value verification for 9 dynamic flow points based on the above method;
FIG. 16 shows the repeatability error results for dynamic flow points;
FIG. 17 is a schematic diagram illustrating an indication error of a measurement of the present disclosure;
fig. 18 is a schematic illustration of the repeatability of the measurement of the present disclosure and the conventional measurement.
Embodiments of the present disclosure will be described in more detail below with reference to the accompanying drawings. While certain embodiments of the present disclosure are shown in the drawings, it is to be understood that the disclosure may be embodied in various forms and should not be construed as limited to the embodiments set forth herein. Rather, these embodiments are provided for a more thorough and complete understanding of the present disclosure. It should be understood that the drawings and embodiments of the disclosure are for illustration purposes only and are not intended to limit the scope of the disclosure.
It should be noted that, for convenience of description, only the portions related to the related invention are shown in the drawings. The embodiments and features of the embodiments in the present disclosure may be combined with each other without conflict.
It should be noted that the terms "first", "second", and the like in the present disclosure are only used for distinguishing different devices, modules or units, and are not used for limiting the order or interdependence relationship of the functions performed by the devices, modules or units.
It is noted that references to "a", "an", and "the" modifications in this disclosure are intended to be illustrative rather than limiting, and that those skilled in the art will recognize that "one or more" may be used unless the context clearly dictates otherwise.
The names of messages or information exchanged between devices in the embodiments of the present disclosure are for illustrative purposes only, and are not intended to limit the scope of the messages or information.
The present disclosure will be described in detail below with reference to the accompanying drawings in conjunction with embodiments.
Referring first to fig. 1, fig. 1 is a schematic view of a reflection type structure of a transit time ultrasonic gas flowmeter according to some embodiments of the present disclosure. As shown in fig. 1, an ultrasonic gas flow meter generally includes an upstream ultrasonic transducer 1 and a downstream ultrasonic transducer 2. The upstream ultrasonic transducer 1 and the downstream ultrasonic transducer 2 are symmetrically inserted into the pipeline 3. In the working state, the upstream ultrasonic transducer 1 sends out an upstream ultrasonic signal, and the upstream ultrasonic signal is reflected by the inner wall of the pipeline 3 and received by the downstream ultrasonic transducer 2. The time from the upstream ultrasonic transducer 1 to the downstream ultrasonic transducer 2 of the upstream ultrasonic signal can be called the transit time T of the downstream ultrasonic signald. Meanwhile, the downstream ultrasonic transducer 2 sends out a downstream ultrasonic signal, and the downstream ultrasonic signal is reflected by the inner wall of the pipeline 3 and received by the upstream ultrasonic transducer 1. The time from the downstream ultrasonic transducer 2 to the upstream ultrasonic transducer 1 of the downstream ultrasonic signal can be called the transit time T of the upstream ultrasonic signalu. The propagation distance L of the upstream ultrasonic signal can be determined by the diameter D of the pipeline and the placement posture of the upstream ultrasonic transducer or the downstream ultrasonic transducer. The placing posture can be an included angle theta between the upstream ultrasonic transducer or the downstream ultrasonic transducer and the pipe wall of the carried medium. Next, the flow rate of the medium in the pipe can be determined by the following formula:
Figure BDA0002540014280000061
where Δ T represents the upstream ultrasonic signal transit time TuTransit time T with downstream ultrasonic signaldThe time difference of (a).
Referring next to fig. 2, fig. 2 is a flow chart 200 of some embodiments of a flow rate metering method of a transit time ultrasonic gas meter according to some embodiments of the present disclosure. The method comprises the following steps:
step 201, obtaining the posture information of the upstream ultrasonic transducer and the downstream ultrasonic transducer.
In some embodiments, the attitude information may be an angle between the upstream ultrasonic transducer and the downstream ultrasonic transducer when the upstream ultrasonic transducer and the downstream ultrasonic transducer are inserted into the pipeline. Typically, the upstream ultrasonic transducer and the downstream ultrasonic transducer are symmetrically disposed in the conduit. The execution body can be a data processing module of the flowmeter. The staff can pass through the man-machine interface with the contained angle of upstream ultrasonic transducer and above-mentioned downstream ultrasonic transducer and pipeline, inputs this contained angle data into this data processing module.
Step 202, acquiring the propagation distance of the ultrasonic signal from the upstream ultrasonic transducer to the downstream ultrasonic transducer.
In some embodiments, the propagation distance may be calculated from the inner diameter of the pipe and the angle data obtained in step 202. Likewise, the inner diameter of the tubing may be used to input this data into the actuator body via a human interface. Specifically, the description will be given with reference to fig. 1. The propagation distance L may be twice the quotient of D and sin θ in the pipeline. Namely, it is
Figure BDA0002540014280000071
Step 203, determining a time difference through a cross-power spectral density function and a coherence function.
In some embodiments, the time difference represents the upstream ultrasonic signal transit time TuTransit time T with downstream ultrasonic signaldThe time difference of (a).
Specifically referring to fig. 3, fig. 3 is a schematic diagram of the upstream and downstream ultrasonic signals and the upstream and downstream echo signals of the present disclosure. The time difference may be determined by:
in a first step, an upstream echo signal 12 and a downstream echo signal 22 are acquired. Specifically, the upstream echo signal 12 and the downstream echo signal 22 may be acquired by a downstream ultrasonic transducer or an upstream ultrasonic transducer. The upstream echo signal 12 and the downstream echo signal 22 may be preprocessed by those skilled in the art, so that the upstream echo signal 12 or the downstream echo signal 22 is obtained by the downstream ultrasonic transducer or the receiving array element of the upstream ultrasonic transducer.
And secondly, determining the correlation degree of the upstream echo signal and the downstream echo signal in the time domain through a cross-correlation function. The cross-correlation degree is determined by a cross-correlation function, and the specific formula is as follows:
Figure BDA0002540014280000072
wherein t represents a selected time, wherein the selected time may be set by a technician or set by the default of the execution main body; τ represents the time delay of the upstream echo signal and the downstream echo signal, and is also the abscissa of the cross-correlation function, and the value range is between 0 and T, and the selection of τ may be selected by a technician or may be set by default by an executing entity; representing the degree of correlation of the upstream echo signal and the downstream echo signal at the time of time delay τ; x (t) represents the amplitude of the upstream echo signal or the downstream echo signal at time t; y (t + τ) represents the amplitude of the downstream echo signal or the upstream echo signal at time t + τ; t represents a sampling period, which may be set by a technician or set by the execution agent by default, for example; for discrete time, the above correlation function translates into:
Figure BDA0002540014280000073
wherein, m represents the selected discrete time, wherein the discrete time can be set by a technician or can be set by the default of the execution main body; n represents the discrete time delay of the downstream echo signal and the upstream echo signal, is also the abscissa of the discrete cross-correlation function, and has a value range of an integer from 0 to N, and the selection of N can be selected by a technician or can be set by executing a main body default; n represents a period of discrete sampling, which may be set by a technician or set by the execution agent as a default.
Next, the cross-power spectral density function of the upstream echo signal and the downstream echo signal is:
Figure BDA0002540014280000081
wherein e is-j2πωτRepresenting a complex variable function; the cross-power spectral density function is expressed in complex form, with complex polar coordinates as follows:
Figure BDA0002540014280000082
Figure BDA0002540014280000083
wherein, Pxy(ω) represents a cross-power spectral density function; im [ P ]xy(ω)]An imaginary number value part representing a cross-power spectral density function; re [ P ]xy(ω)]A real number value part representing a cross-power spectral density function; j represents an imaginary number. Referring to fig. 4 and 5, fig. 4 and 5 are graphs of the cross-spectral phase of the upstream echo signal and the downstream echo signal, respectively. The correlation between the upstream echo signal and the downstream echo signal in the frequency domain is determined by the following formula:
Figure BDA0002540014280000084
wherein, γ2 xy(ω) is an amplitude-squared coherence function, which represents the correlation between the upstream echo signal and the downstream echo signal in the frequency domain, and may also be referred to as a correlation coefficient, and the maximum value is 1; px(ω) represents the self-power spectral density function of x, where x represents the upstream echoA signal; py(ω) represents the self-power spectral density function of y, where y represents the downstream echo signal. Wherein, wmIs shown as
Figure BDA0002540014280000085
When the value is maximum, the value of w is obtained; omegamAlso representing the signal frequency. Referring to fig. 6 and 7, fig. 6 and 7 are schematic diagrams illustrating the correlation degree of the upstream echo signal and the downstream echo signal, respectively. The phase difference between the upstream echo signal and the downstream echo signal is thetaxy(wm) (ii) a According to the relation between the cross-power spectrum phase and the frequency, the phase difference can be obtained by a trigonometric function relation:
Figure BDA0002540014280000086
the time difference is determined according to the following formula:
Figure BDA0002540014280000087
in some embodiments, the determining the transit time of the upstream ultrasonic signal and the transit time of the downstream ultrasonic signal based on the variable threshold includes:
the above variable threshold is determined by the following formula:
δ ═ 0.1max (a (t)); wherein δ represents a variable threshold; a (t) represents an envelope of the upstream echo signal or the downstream echo signal; wherein a (t) is determined by the following formula:
Figure BDA0002540014280000088
where x (t) represents an upstream echo signal or a downstream echo signal: wherein, the above
Figure BDA0002540014280000089
Is determined by the following formula:
Figure BDA00025400142800000810
wherein the content of the first and second substances,
Figure BDA00025400142800000811
represents the convolution of x (t) and 1/π t, with the formula:
Figure BDA00025400142800000812
wherein, H [ x (t)]Denotes the Hilbert transform of x (t).
In some embodiments, the determining the transit time of the upstream ultrasonic signal and the transit time of the downstream ultrasonic signal based on the variable threshold includes:
and determining a maximum value having a positive slope and an amplitude greater than 0 at the rising stage of the envelope of the upstream echo signal or the downstream echo signal. Specifically, the description will be given with reference to fig. 8. Fig. 8 is a schematic diagram of the first echo signal or the second echo signal. As shown in fig. 8, the local maxima are the values of the locations characterized by the "triangles" in fig. 8; the local maximum point is determined by the following formula: nn0O + δ; wherein nn0Representing a local maximum point Z0The ordinate of (a); δ represents a variable threshold; o represents the minimum difference between the variable threshold and the local maximum; wherein the above O is determined by the following formula: o ═ min (| δ -nn)iI)); wherein nniRepresents the ith local maximum ZiThe ordinate of (a); by aligning the local maximum point Z0Two adjacent sampling points S-1,S+1Performing parabolic interpolation processing on the three points to obtain a peak value; determining the peak value as an upstream feature point or a downstream feature point Sp(mp,np)。
Next, determination of the transit time of the upstream echo signal or the downstream echo signal will be described with reference to fig. 8. Specifically, the transit time of the upstream echo signal or the downstream echo signal may be determined according to the following formula: t isu=Tp+Tu’;Td=Tp+Td'; wherein, TPAnd a fixed time obtained according to an experiment from the sending of the upstream ultrasonic signal or the downstream ultrasonic signal to the receiving of the upstream echo signal or the downstream echo signal, wherein the fixed time can be adjusted according to the size of the pipeline. T isu' represents the time from the sampling start of the upstream echo signal to the upstream feature point; t isd' represents the time from the sampling start of the downstream echo signal to the downstream feature point, wherein T isd' is determined by the following equation: t isd’=mp1(ii) Fs; wherein m isp1Indicating a downstream feature point Sp1The abscissa of (a); fs represents the sampling rate, where T is the aboveu' is determined by the following equation: t isu’=mp2/Fs, wherein mp2Represents the upstream feature point Sp2The abscissa of (a).
The flow velocity metering method of the time difference type ultrasonic gas flowmeter disclosed by the invention has the advantages that firstly, the time delay is directly calculated by adopting a cross-power spectrum method, the time difference of upstream and downstream flight signals is determined, and compared with a related algorithm for respectively calculating the time of the upstream and downstream flight signals, the calculated amount is reduced. Meanwhile, the mutual power method is an unstable process, and nanosecond-level accurate time delay can be directly obtained. And further improves the accuracy of the flow velocity calculation.
In addition, the variable threshold is determined by adopting Hilbert transform, different thresholds are set according to waveforms of different signals, and the position of the echo signal is further accurately determined. Thereby improving the accuracy of the method.
The following are calibration experiments performed.
First, a calibration apparatus is described. The calibration device is a set of sonic nozzle type calibration device, and the uncertainty of the calibration device is 0.5%. When calibrating the test meter, the test meter is mounted on the pipeline. The sonic nozzle type calibration device provides standard flow for a test instrument. The test instrument measures the pulse rate according to the pulse coefficient and converts the frequency value into the pulse number. The pulse equivalent K of the test meter was set to 500 pulses/m3. 500 pulses represent 1m3/h。
Specifically, the ultrasonic flowmeter will be described with reference to fig. 9. Fig. 9 is a schematic diagram of functional block diagrams of various parts of the ultrasonic flow meter. As shown in FIG. 9, the device comprises a flowmeter prototype, an ultrasonic sensor, a pulse emission driving module, an amplifier filtering module, a comparator module, and a statorThe device comprises a time module, a power module, a control data processing module and a storage display module. Wherein, the ultrasonic transducer can send an excitation signal and receive an echo signal, and the center frequency of the ultrasonic transducer is 200 KHz. The JK type transducer is arranged at a very small caliber (176.7 mm)2) Experimental data acquisition is carried out, and the feasibility of the algorithm is verified. And in medium size pipes (1963.5 mm)2) Experimental verification was performed above. It has very good portability.
In the formula
Figure BDA0002540014280000101
Figure BDA0002540014280000102
Is a constant number of times, and is,
Figure BDA0002540014280000103
the value of (d) determines the final average flow rate, and experiments show that the delta T time difference is decisive in the calculation.
Δ T is taken as a numerator, which is about 1 to 400 times the dynamic range of the denominator. The denominator is the product of the transit time of the upstream echo signal and the transit time of the downstream echo signal, and does not vary much from numerator over the full denominator. The sharp change in the time difference reflects a sharp change in the flow rate. The significance of the time difference is further verified and can be explained with reference to fig. 10. FIG. 10 is a time of flight T of an upstream echo signaluTransit time T of the downstream signaldAnd the time difference Δ T at 9 flow points. The cross-sectional area S is 176.7mm2Time, 100 sets of Td,TuAnd the average value of Δ T. At the same time, the effect of temperature on this conclusion is also discussed. When the temperature rises to 70 ℃, TuAnd TdWill simultaneously decrease a little, however, Td×TuThe contribution of the value of (a) to the overall flow calculation is also much less than at.
Next, the importance of the accuracy Δ T can also be seen by deriving from the relationship between the calculated data Δ T and the average flow velocity v, and the measurement of Δ T can be seen in fig. 10. Fig. 11 is a schematic diagram of the approximate linear relationship between the time difference parameter and the standard flow point, thereby verifying the significance of accurately calculating the time difference. There are 100 sets of time reference values at 9 standard flow points in the graph, which gives the concentration of time delays at each standard flow point.
Specifically, the experimental results were divided into two parts, one being static results and the other being dynamic indicating error and repeatability results at 9 standard flow points. Static experiments were performed at ambient (25.9 ℃) and elevated (70.0 ℃) temperatures. In the static condition, the rubber cover is used for carrying out spiral connection on the pipeline, so that the static condition is ensured.
At this time, the zero point of the ultrasonic gas flowmeter is detected at the static flow rate. According to the ultrasonic sound velocity, the ultrasonic sound velocity at the ambient temperature and the length of the ultrasonic propagation channel, the propagation time of the ultrasonic wave can be calculated. By this method and experiment, we can roughly estimate the propagation time T of the pulsepAnd is used to determine a fixed value between the ADC (Analog-to-Digital Converter) acquisition window and the pulse emission. Then, based on the cumulative flow over a period of time, 0m is obtained3Average measurement error of/h, and record K in software, complete zero calibration. The experimental measurement is performed by using a sonic nozzle type gas meter calibration device, which is specifically shown in fig. 12. FIG. 12 is a schematic illustration of an ultrasonic flow meter calibration test. The proposed ultrasonic flow meter system has been designed for use with a reflective structure flow meter. The algorithm has been tested in 2 pairs of 500kHz transducers and 5 pairs of 200kHz transducers.
Referring next to fig. 13 and 14, fig. 13 and 14 are schematic diagrams of near-zero measurements at high temperatures of 25.9 ℃ and 70.0 ℃, respectively, at a filter order of 30 for a moving average.
The maximum zero drift flow at 25.9 ℃ without filtering and after filtering is respectively: q is 0.01567m3H and Q0.00272 m3H is used as the reference value. At a high temperature of 70.0 ℃, the maximum zero drift flow rate of the unfiltered and the maximum zero drift flow rate of the filtered are respectively Q-0.01935 m3H and Q0.00295 m3H is used as the reference value. The raw measured flow rates measured in FIGS. 13 and 14 are broken line portions, dashed line portions of the graphThe inching flow, the curved portion represents the flow after filtering. And adds a visual comparison. As can be seen from fig. 13 and 14, the flow calculation model formed by the analysis model and the averaging method satisfies the zero point error requirement.
Before the start of the zero point test experiment, the instruments of the experimental prototype were placed in a high and low temperature laboratory box and the temperature was gradually increased to 70 ℃ over 2 hours. This temperature was maintained for 2 hours. During this time, the air temperature was stabilized at 70 ℃ with 1% fluctuation. Zero drift is achieved at high temperatures. At the zero flow point, the calculated flow needs to be smaller than the initial flow q according to the national standardsAt this time, no flow is defaulted. The pick-up flow of the G4 meter was 0.003m3H is used as the reference value. FIG. 15 shows the results of the test at a high temperature of 70 ℃. From the experimental results, it can be seen that the zero flow drift problem is satisfied even at a high temperature of 70 ℃.
For the dynamic result, an air dynamic experiment is carried out on the standard sonic nozzle type gas meter calibration device, and the measurement performance of the prototype gas flowmeter is tested. The center frequency and diameter of the cylindrical transducer were 200kHz and 18mm, respectively. The AD sampling rate is 1 MHz. The system is suitable for accurate gas metering at extremely small flow rates. The flow velocity detection range is 0.04-6.00 m3H is used as the reference value. The volume of the reference gas was measured by a sonic nozzle type calibration apparatus, and its uncertainty was 0.5% and the temperature was about 25 ℃. According to the regulations of the ultrasonic gas flowmeter verification regulation, three equivalent limit values (minimum gas power qmin, transition gas power and maximum gas power q) are determined in the measurement power rangemax) Q in this textt=0.1qmax. Thereby, the gas flow rate is controlled within the range of 0.04 to 6.00m3Per hour) at 0.6m3The/h dividing flow is divided into two groups, and the requirements on indicating value error and repeatability are different under the two flow ranges.
The final flow value needs to be verified with the value of K. The flow calculation formula is as follows: KvS where S is the cross-sectional area of the pipe and different flow points have different values of K. Fig. 15 is a schematic diagram of the indicating error of K value verification for 9 dynamic flow points based on the above method. Since the K value verification can correct the indicating value error correctly, the key effect of repeatability is illustrated. The computation of the transit time of the echo signal uses the Hilbert transform to set a variable threshold. The results of comparing the conventional measurement method with the measurement method of the present disclosure at 9 dynamic flow points are shown in fig. 16. Figure 16 shows the repeatability error results for dynamic flow points.
Fig. 15 and 16 show the recognition error and the repeatability under the condition of extremely small aperture for k value verification, respectively, and compare the results of the two principle calculations. With reference to fig. 17, fig. 17 shows the indication error of the measurement result of the present disclosure, and in the case of a small flow, the repeatability of the conventional measurement method is hard to guarantee and the proposed method meets the precision requirement. In modern ultrasonic flow measurement, accurate measurement of high-precision small flow is the key to realizing overall high-precision dynamic measurement results.
With reference to fig. 18, fig. 18 is a schematic illustration of the repeatability of the measurement results of the present disclosure and the conventional measurement results. The subject of this test is a very small diameter flowmeter with a cross-sectional diameter of 176.7mm2. It has high requirements for repeatability at small flow rates. In this case, the time difference needs to be measured more accurately. This method recalculates a time difference rather than subtracting absolute time of flight directly.
The flow velocity metering method of the time difference type ultrasonic gas flowmeter verifies the importance of the time difference in flow calculation by using the approximate linear relation between the time difference and the standard flow. The influence of the dynamic change of the time difference value on the flow calculation of the full time difference formula is analyzed. This makes accurate measurement of transit-time ultrasonic gas meters critical to time difference estimation. And establishing a separation model of the time difference and the time difference of the upstream echo signal and the downstream echo signal. The cross-power spectrum method is used for calculating the phase difference of the uplink signal and the downlink signal, the time difference is obtained through phase difference analysis, then the variable threshold is set through Hilbert transform, the absolute flight time is obtained, and the calculated amount is reduced. After the K value is detected, the indicating value error is reduced, and the repeatability is the most important index of the whole ultrasonic measurement system. The above experiment compares the calculated delay after subtraction by the variable threshold method with the repetitive results of the proposed method. The experimental results show the superiority of the algorithm and the superiority of the model.
The algorithm has higher precision, and the experimental result also proves that the algorithm is suitable for small pipelines capable of accurately detecting the flow. Under the condition of large flow experiment, the method has stable result and meets the repeatability requirement. The algorithm is transplanted to the pipeline with the cross section area of 1963mm2And good repeatability results are obtained.
The foregoing description is only exemplary of the preferred embodiments of the disclosure and is illustrative of the principles of the technology employed. It will be appreciated by those skilled in the art that the scope of the invention in the embodiments of the present disclosure is not limited to the specific combination of the above-mentioned features, but also encompasses other embodiments in which any combination of the above-mentioned features or their equivalents is made without departing from the inventive concept as defined above. For example, the above features and (but not limited to) technical features with similar functions disclosed in the embodiments of the present disclosure are mutually replaced to form the technical solution.

Claims (6)

1. A flow rate metering method of a time difference type ultrasonic gas flowmeter comprises the following steps:
acquiring attitude information of an upstream ultrasonic transducer and a downstream ultrasonic transducer;
acquiring the propagation distance of an ultrasonic signal from the upstream ultrasonic transducer to the downstream ultrasonic transducer;
determining a time difference by a cross-power spectral density function and a coherence function, wherein the time difference is a difference between a transit time of the upstream ultrasonic signal and a transit time of the downstream ultrasonic signal;
analyzing the upstream ultrasonic signal and the downstream ultrasonic signal based on Hilbert transform, and determining a variable threshold;
determining a transit time of the upstream ultrasonic signal and a transit time of the downstream ultrasonic signal based on the variable threshold;
and determining the flow velocity of the medium detected by the flowmeter according to the propagation distance, the attitude information, the time difference, the transition time of the upstream ultrasonic signal and the transition time of the downstream ultrasonic signal.
2. The method of claim 1, wherein the flow rate is determined by the formula:
Figure FDA0002540014270000011
wherein v represents a flow rate of the medium detected by the flow meter;
l represents the propagation distance;
theta represents the included angle between the upstream ultrasonic transducer and the downstream ultrasonic transducer and the pipe wall of the carried medium;
Δ T represents a time difference;
Turepresenting a transit time of the upstream ultrasonic signal;
Tdrepresenting the transit time of the downstream ultrasonic signal.
3. The method of claim 2, wherein the determining the time difference by a cross-power spectral density function and a coherence function comprises:
acquiring an upstream echo signal and a downstream echo signal;
determining the correlation degree of the upstream echo signal and the downstream echo signal in the time domain through a cross-correlation function;
the cross-correlation function of the upstream echo signal and the downstream echo signal in the time domain is:
Figure FDA0002540014270000021
wherein t represents the selected time;
τ represents a time delay of the upstream echo signal and the downstream echo signal;
Rxy(τ) represents the degree of correlation of the upstream echo signal and the downstream echo signal at a time delay τ;
x (t) represents the amplitude of the upstream echo signal or the downstream echo signal at time t;
y (t + τ) represents the amplitude of the downstream echo signal or the upstream echo signal at time t + τ;
t represents the period of sampling;
for discrete time, the correlation function translates into:
Figure FDA0002540014270000022
wherein m represents the selected discrete time;
n represents a discrete time delay of the downstream echo signal and the upstream echo signal;
n denotes the period of discrete sampling.
4. The method of claim 3, wherein the determining the time difference by a cross-power spectral density function and a coherence function comprises:
the cross-power spectral density function of the upstream echo signal and the downstream echo signal is:
Figure FDA0002540014270000023
wherein e is-j2πωτRepresenting a complex variable function;
the cross-power spectral density function is expressed in complex form, complex polar coordinate form as follows:
Figure FDA0002540014270000024
wherein, Im [ P ]xy(ω)]An imaginary number value part representing a cross-power spectral density function;
Re[Pxy(ω)]a real number value part representing a cross-power spectral density function;
j represents an imaginary number;
the correlation of the upstream echo signal and the downstream echo signal in the frequency domain is determined by the following formula:
Figure FDA0002540014270000031
wherein, γ2 xyAnd (ω) is the magnitude squared coherence as a function of the power spectral density. The correlation of the upstream echo signal and the downstream echo signal on the frequency domain is represented, namely the correlation coefficient, and the maximum value is 1;
Px(ω) represents the self-power spectral density function of x;
Py(ω) represents the self-power spectral density function of y;
wherein, ω ismIs shown as
Figure FDA0002540014270000032
When the value is maximum, the value of omega is obtained; omegamRepresents the signal frequency when the phase difference between the upstream echo signal and the downstream echo signal is thetaxym) According to the relation between the cross-power spectrum phase and the frequency, the phase difference can be obtained by a trigonometric function relation:
Figure FDA0002540014270000033
determining the time difference according to the following formula:
Figure FDA0002540014270000034
5. the method of claim 4, wherein said determining a time of flight of an upstream ultrasonic signal and a time of flight of a downstream ultrasonic signal based on the variable threshold comprises:
the variable threshold is determined by the following formula:
δ=0.1max(A(t));
wherein δ represents a variable threshold;
a (t) represents an envelope of the upstream echo signal or the downstream echo signal;
wherein A (t) is determined by the following formula:
Figure FDA0002540014270000035
x (t) represents the upstream or downstream signal:
wherein the Z (t) is determined by the following formula:
Figure FDA0002540014270000036
wherein the content of the first and second substances,
Figure FDA0002540014270000037
represents the convolution of x (t) and 1/t, and has the formula
Figure FDA0002540014270000038
Wherein H [ x (t) ] represents Hilbert transform of x (t).
6. The method of claim 5, wherein said determining a time of flight of an upstream ultrasonic signal and a time of flight of a downstream ultrasonic signal based on the variable threshold comprises:
determining the maximum value of the envelope rising stage of the upstream echo signal or the downstream echo signal;
the local maximum point is determined by the following formula:
nn0=O+δ;
wherein nn0Representing a local maximum point Z0The ordinate of (a);
δ represents a variable threshold;
o represents the minimum difference value between the variable threshold and the maximum value ordinate of the envelope rising stage;
wherein said O is determined by the formula:
O=min(|δ-nni|);
wherein nniRepresents the ith local maximum ZiThe ordinate of (a);
by aligning the local maximum point Z0Two adjacent sampling points S-1,S+1Performing parabolic interpolation processing on the three points to obtain a peak value;
determining the peak value as an upstream feature point or a downstream feature point Sp(mp,np);
Determining a transit time of the upstream echo signal or the downstream echo signal according to the following equation:
Tu=Tp+Tu’;
Td=Tp+Td’;
wherein, TPRepresenting a fixed time obtained according to an experiment from the sending of the upstream ultrasonic signal or the downstream ultrasonic signal to the receiving of the upstream echo signal or the downstream echo signal;
Tu' represents the time from the sampling start of the upstream echo signal to the upstream feature point;
Td' denotes the time from the start of sampling to the downstream feature point of the downstream echo signal,
wherein, T isd' is determined by the following equation:
Td’=mp1/Fs;
wherein m isp1Indicating a downstream feature point Sp1The abscissa of (a);
fs represents the sampling rate of the sample,
wherein, T isu' is determined by the following equation:
Tu’=mp2/Fs;
wherein m isp2Represents the upstream feature point Sp2The abscissa of (a).
CN202010543939.8A 2020-06-15 2020-06-15 Flow velocity metering method of time difference type ultrasonic gas flowmeter Active CN113375737B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010543939.8A CN113375737B (en) 2020-06-15 2020-06-15 Flow velocity metering method of time difference type ultrasonic gas flowmeter

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010543939.8A CN113375737B (en) 2020-06-15 2020-06-15 Flow velocity metering method of time difference type ultrasonic gas flowmeter

Publications (2)

Publication Number Publication Date
CN113375737A true CN113375737A (en) 2021-09-10
CN113375737B CN113375737B (en) 2024-01-12

Family

ID=77568924

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010543939.8A Active CN113375737B (en) 2020-06-15 2020-06-15 Flow velocity metering method of time difference type ultrasonic gas flowmeter

Country Status (1)

Country Link
CN (1) CN113375737B (en)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113916719A (en) * 2021-10-12 2022-01-11 北京航空航天大学 Fluid density and flow rate online synchronous detection system and detection method
CN114397475A (en) * 2022-03-25 2022-04-26 青岛鼎信通讯股份有限公司 Water flow velocity measuring method suitable for ultrasonic water meter

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20020143479A1 (en) * 2001-02-15 2002-10-03 Satoshi Fukuhara Ultrasonic flowmeter
CN103542901A (en) * 2012-07-09 2014-01-29 德克萨斯仪器股份有限公司 Flow meter
CN105891326A (en) * 2015-02-16 2016-08-24 森萨克申公司 Method for determining properties of medium and device for determining properties of medium
JP2017187310A (en) * 2016-04-01 2017-10-12 株式会社ソニック Ultrasonic flowmeter
US20170350842A1 (en) * 2016-06-03 2017-12-07 Mohr and Associates Method for Indentifying and Measuring Volume Fraction Constituents of a Fluid

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20020143479A1 (en) * 2001-02-15 2002-10-03 Satoshi Fukuhara Ultrasonic flowmeter
CN103542901A (en) * 2012-07-09 2014-01-29 德克萨斯仪器股份有限公司 Flow meter
CN105891326A (en) * 2015-02-16 2016-08-24 森萨克申公司 Method for determining properties of medium and device for determining properties of medium
JP2017187310A (en) * 2016-04-01 2017-10-12 株式会社ソニック Ultrasonic flowmeter
US20170350842A1 (en) * 2016-06-03 2017-12-07 Mohr and Associates Method for Indentifying and Measuring Volume Fraction Constituents of a Fluid

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
SIHAO SUN: "A Novel Signal Processing Method Based on Cross-correlation and Interpolation for ToF Measurement", 《 2019 IEEE 4TH INTERNATIONAL CONFERENCE ON SIGNAL AND IMAGE PROCESSING (ICSIP)》 *
王飞: "数字式时差法超声流量计的设计与实现", 《中国优秀博硕士学位论文全文数据库(硕士) 工程科技Ⅱ辑》 *
王飞: "数字式时差法超声流量计的设计与实现", 《中国优秀博硕士学位论文全文数据库(硕士) 工程科技Ⅱ辑》, no. 01, 15 January 2015 (2015-01-15), pages 47 - 49 *
王飞: "数字式时差法超声流量计的设计与实现", 《自动化仪表》 *
王飞: "数字式时差法超声流量计的设计与实现", 《自动化仪表》, vol. 35, no. 9, 30 September 2014 (2014-09-30) *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113916719A (en) * 2021-10-12 2022-01-11 北京航空航天大学 Fluid density and flow rate online synchronous detection system and detection method
CN114397475A (en) * 2022-03-25 2022-04-26 青岛鼎信通讯股份有限公司 Water flow velocity measuring method suitable for ultrasonic water meter

Also Published As

Publication number Publication date
CN113375737B (en) 2024-01-12

Similar Documents

Publication Publication Date Title
CN107003332B (en) Improved signal travel time flow meter
US9689726B2 (en) Flow meter
Zhu et al. Variable ratio threshold and zero-crossing detection based signal processing method for ultrasonic gas flow meter
CN101078640B (en) Ultrasonic wave air flow-meter and device for measuring waste gas flow of internal combustion engine and method for obtaining gas flow
RU2478919C2 (en) Calibration method under operating conditions of ultrasonic flow metres - counters of flow rate and volume of liquid single-phase media
CN107655533B (en) A kind of Ultrasonic Wave Flowmeter signal processing method and system based on backward energy integral
JP2002243514A (en) Ultrasonic flow meter
Tian et al. Energy peak fitting of echo based signal processing method for ultrasonic gas flow meter
CN113375737A (en) Flow velocity metering method of time difference type ultrasonic gas flowmeter
Zhu et al. Mathematical modeling of ultrasonic gas flow meter based on experimental data in three steps
CN111157065A (en) Acoustic time delay measuring method in ultrasonic signal transmission loop of gas ultrasonic flowmeter
US6816808B2 (en) Peak switch detector for transit time ultrasonic meters
CN112304376B (en) Ultrasonic flowmeter flow measuring method based on data fusion
CN109632025B (en) Data processing method for flow velocity of acoustic Doppler flowmeter
CN116558587A (en) Ultrasonic fluid flow measurement method and system based on flow prediction
Ma et al. Signal processing method based on connection fitting of echo peak point with a large slope for ultrasonic gas flow meter
Bucci et al. Numerical method for transit time measurement in ultrasonic sensor applications
JP4904099B2 (en) Pulse signal propagation time measurement device and ultrasonic flow measurement device
Li et al. Research on transit-time ultrasonic flowmeter with signal characteristic analysis
Li et al. A novel differential time-of-flight algorithm for high-precision ultrasonic gas flow measurement
Li et al. Study on transit-Time ultrasonic flow meter with waveform analysis
Zheng et al. A new characteristic peaks group judgement method for the accurate measurement of time‐of‐flight in the ultrasonic gas flowmeter
Ma Ultrasonic Gas Flowmeter Based on Dynamic Reference Waveform Cross-Correlation Algorithm
Kong et al. Time Delay Study of Ultrasonic Gas Flowmeters Based on VMD–Hilbert Spectrum and Cross-Correlation
Lucena et al. An innovative ultrasonic time of flight method based on extended Kalman filter for wind speed measurement

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