WO2014063635A1 - 数据采样方法与系统及其在参数辨识中的应用方法与系统 - Google Patents

数据采样方法与系统及其在参数辨识中的应用方法与系统 Download PDF

Info

Publication number
WO2014063635A1
WO2014063635A1 PCT/CN2013/085813 CN2013085813W WO2014063635A1 WO 2014063635 A1 WO2014063635 A1 WO 2014063635A1 CN 2013085813 W CN2013085813 W CN 2013085813W WO 2014063635 A1 WO2014063635 A1 WO 2014063635A1
Authority
WO
WIPO (PCT)
Prior art keywords
sampling
value
sequence
steady state
parameter identification
Prior art date
Application number
PCT/CN2013/085813
Other languages
English (en)
French (fr)
Inventor
郝玉山
Original Assignee
Hao Yushan
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 Hao Yushan filed Critical Hao Yushan
Priority to JP2015537134A priority Critical patent/JP6069809B2/ja
Priority to US14/437,972 priority patent/US9438449B2/en
Priority to CA2889334A priority patent/CA2889334C/en
Priority to EP13849426.5A priority patent/EP2913930A4/en
Publication of WO2014063635A1 publication Critical patent/WO2014063635A1/zh

Links

Classifications

    • H—ELECTRICITY
    • H04—ELECTRIC COMMUNICATION TECHNIQUE
    • H04L—TRANSMISSION OF DIGITAL INFORMATION, e.g. TELEGRAPHIC COMMUNICATION
    • H04L25/00—Baseband systems
    • H04L25/02—Details ; arrangements for supplying electrical power along data transmission lines
    • H04L25/0262—Arrangements for detecting the data rate of an incoming signal
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F7/00—Methods or arrangements for processing data by operating upon the order or content of the data handled
    • G06F7/38—Methods or arrangements for performing computations using exclusively denominational number representation, e.g. using binary, ternary, decimal representation
    • G06F7/48—Methods or arrangements for performing computations using exclusively denominational number representation, e.g. using binary, ternary, decimal representation using non-contact-making devices, e.g. tube, solid state device; using unspecified devices
    • G06F7/52—Multiplying; Dividing
    • G06F7/535—Dividing only
    • H—ELECTRICITY
    • H03—ELECTRONIC CIRCUITRY
    • H03M—CODING; DECODING; CODE CONVERSION IN GENERAL
    • H03M1/00—Analogue/digital conversion; Digital/analogue conversion
    • H03M1/06—Continuously compensating for, or preventing, undesired influence of physical parameters
    • H03M1/0617—Continuously compensating for, or preventing, undesired influence of physical parameters characterised by the use of methods or means not specific to a particular type of detrimental influence
    • H03M1/0626—Continuously compensating for, or preventing, undesired influence of physical parameters characterised by the use of methods or means not specific to a particular type of detrimental influence by filtering
    • H—ELECTRICITY
    • H04—ELECTRIC COMMUNICATION TECHNIQUE
    • H04L—TRANSMISSION OF DIGITAL INFORMATION, e.g. TELEGRAPHIC COMMUNICATION
    • H04L25/00—Baseband systems
    • H04L25/02—Details ; arrangements for supplying electrical power along data transmission lines
    • H04L25/06—Dc level restoring means; Bias distortion correction ; Decision circuits providing symbol by symbol detection
    • H04L25/068—Dc level restoring means; Bias distortion correction ; Decision circuits providing symbol by symbol detection by sampling faster than the nominal bit rate

Definitions

  • the present invention relates to computer application techniques, particularly sampling, parameter identification, and digital control techniques.
  • sampling The process of converting analog quantities into digital quantities is called sampling.
  • the object of the present invention is to solve the above problems. After being digitally sampled by a stable analog controller and converted to a digital controller, it is still stable and has the same characteristics, and accurate parameters can be identified based on the sampled data.
  • the inventor discovered through repeated trials and summaries and theoretical derivation analysis that the cause of the above problems was not caused by the identification algorithm or the digital control method, but by the sampling data; the sampling data problem is not sampling. The result is too slow, but the sampling is too fast; the sampling frequency is not as fast as possible as generally recognized by those skilled in the art, but there is an upper limit of the sampling frequency, and the sampling frequency exceeding the upper limit is caused by a large error in the sampling data.
  • the invention provides a data sampling method, which comprises the following steps: 51. Perform analog low-pass filtering on the physical quantity of the electrical signal y to obtain the filtered electrical signal;
  • the sampling system has a domain error ⁇ ⁇ ⁇ , and the ⁇ is a truncation error of the sampling system; the truncation error ⁇ of the sampling system is related to the sum, which is half of the resolution of the analog-to-digital converter, and the value is The truncation error of the number of bits is represented in the computer.
  • the maximum error ⁇ of the S domain is a maximum error that can be accepted by the S domain determined by the specific application, for example, ⁇ determined by the acceptable error of the zero pole of the application system S domain, or by the application system time domain differential equation
  • the acceptable maximum error 0 corresponds to the determined ⁇
  • the sampling frequency tends to be high in practical applications, and the upper limit of the sampling frequency is not satisfied. Therefore, the present invention also proposes a second data sampling method, which includes the following step:
  • the electrical signal y of the physical quantity is sampled to satisfy the sampling frequency of the Nyquist theorem, and the sampling sequence Y1 is obtained;
  • the cutoff frequency of the low pass filtering is, the / c ⁇ 0.5x, to prevent mixing when resampling.
  • the error of the S domain can be reduced to a predetermined range ⁇ .
  • f sl ⁇ f s ⁇ ⁇ of the first method there may be a contradiction in / s3 ⁇ 4 ⁇ / s , which results in the inability to select the sampling frequency according to f sl ⁇ f s ⁇ ⁇ of the first method.
  • the second data sampling method described above that is, low-pass filtering and heavy The sampling method can solve this contradiction. Since in all the poles ⁇ p n ⁇ and zero ⁇ p ra ⁇ , the main pole is the main zero point P.
  • m plays a leading role, the higher the frequency, the lower the pole or zero effect, so when the contradiction of / s3 ⁇ 4 ⁇ f sl occurs, the high-frequency pole and zero point are filtered by the low-pass filter, and the sampling frequency is reduced by re-sampling, not only The contradiction can be solved, and the main pole P rn ⁇ P main zero point P can also be guaranteed.
  • the accuracy of m where the main pole or the main zero is the pole or zero near the origin on the S domain.
  • the invention also proposes a data sampling system comprising:
  • a sampling unit for sampling the physical quantity of the electrical signal y to satisfy the sampling frequency of the Nyquist theorem to obtain a sampling sequence Y1;
  • a low-pass filter for digitally low-pass filtering the sample sequence Y1 to obtain a filtered sample sequence ⁇ ' J Y2;
  • a resampling unit for resampling the filtered sample sequence Y2 to obtain a sample sequence Y;
  • ⁇ is a domain of the sampling system Error, which is the maximum error of the S domain;
  • the cutoff frequency of the low pass filter is , y ⁇ o. 5x, to prevent mixing when resampling.
  • Digital sampling is used for measurement or for identification parameter or digital control.
  • Parameters include static parameters and dynamic parameters, static parameters are the coefficients of system equations, and dynamic parameters are system differential equations.
  • static parameters refer to the admittance matrix of the power flow equation, that is, the resistance, reactance, and susceptance of components such as transformers, wires, and generators.
  • the dynamic parameters refer to the excitation time constant, the time constant of the speed regulation, and the time constant of the rotor.
  • a parameter identification method using the above two data sampling methods characterized in that it further comprises steps S5, comprising performing a dynamic parameter identification step S51 and/or a static parameter identification step S52 according to the sampling sequence Y.
  • the step S51 of the dynamic parameter identification includes, according to the sampling sequence Y, the step S511 of identifying the order and parameters of the ARMAX model by using the identification method in the adaptive control, or the step S512 of identifying the parameter by the reference model method, or Kalman Step S513 of the filter identification parameter.
  • the static parameter identification step S52 includes a step S521 of first calculating a steady state value for the sampling sequence Y, and a parameter estimation method, a correlation coefficient method, a linear regression method, a linearizable linear regression method according to the steady state value, Or step S522 of identifying a static parameter by a least square method; and the calculating the steady state value step S521 includes:
  • the average value is calculated according to the sampling sequence Y as a steady state value
  • the forced component is calculated as the steady state value according to the sampling sequence Y.
  • the step S5211 determines that the current sample value is in a steady state process or a transient process according to the t distribution or its binization, or the filter output result; the filter is a Kalman filter or an ⁇ filter. For example, performing Kalman filtering or ⁇ filtering on the sampling sequence to obtain components of the state variable; determining whether each component of the state variable exceeds a corresponding set value; if one component or several components exceeds the set value, It is judged that the current sampling value is in a transient process; otherwise, if it is not exceeded, it is judged that the current sampling value is in a steady state process.
  • the Kalman filter is at least 1st order, or higher order, for example, 2nd order, 3rd order, etc., and the invention is not limited.
  • the step S5213 calculates the forcing component according to the following steps:
  • step S52135 is as follows: ⁇ . + ⁇ ⁇ + ⁇ 2 + ⁇ 3 ⁇ 3 + ⁇ + ⁇ resort ⁇ "Calculate the forcing component, where ⁇ 3 ⁇ 4 , ⁇ , ⁇ ⁇ , ⁇ are in accordance with the derivative formula and the physical quantity y and its derivatives The constant determined by the value.
  • the step S52133 determines whether the current value is in a steady state process or a transient process according to the t distribution or its binization, or the filter output result; the filter is a Kalman filter or an ⁇ filter. For example, Kalman filtering or ⁇ filtering is performed on the sampling sequence ⁇ ⁇ to obtain components of the state variable; and it is determined whether each component of the state variable exceeds a corresponding set value; if one component or several components exceeds the set value, Then judge the current value is in the transient process; otherwise, if it is not exceeded, it is judged that the current value is in the steady state process.
  • Kalman filtering or ⁇ filtering is performed on the sampling sequence ⁇ ⁇ to obtain components of the state variable; and it is determined whether each component of the state variable exceeds a corresponding set value; if one component or several components exceeds the set value, Then judge the current value is in the transient process; otherwise, if it is not exceeded, it is judged that the current value is in the steady state process.
  • the Kalman filter is at least 1st order, or higher order, for example, 2nd order, 3rd order, etc., and the invention is not limited.
  • a parameter identification system comprising the above data sampling system, further comprising: a parameter identification unit, configured to perform dynamic parameter identification and/or static parameter identification according to the sampling sequence.
  • the parameter identification unit is a dynamic parameter identification unit, configured to identify the order and parameters of the ARMAX model, or the reference model identification parameter according to the sampling sequence ⁇ , using the identification method in the adaptive control, Or the Kalman filter to identify the parameters.
  • the parameter identification unit is a static parameter identification unit, including a steady state value calculation unit, configured to calculate a steady state value according to the sampling sequence, and a static parameter estimation unit, configured to use the steady state value according to the steady state value.
  • the static parameters are identified by parameter estimation, correlation coefficient method, linear regression method, linear linear regression method, or least squares method.
  • the steady state value calculation unit includes:
  • the determining unit determines, according to the sampling sequence ⁇ , that the current sampling value is in a steady state process or a transient process according to the t distribution or its binization, or the filter output result;
  • the filter is a Kalman filter or an ⁇ filter;
  • the steady state unit ⁇ is calculated for calculating the forced component as the steady state value according to the sampling sequence ⁇ when determining that the current sample value is in the transient process.
  • the steady state unit B is calculated, including:
  • a determining unit configured to receive a data sequence "[ ⁇ , to determine whether the current value is in a transient process or a steady state process
  • the data sampling method of the present invention reduces the truncation error in the digital sampling and digital system to the S domain or the time domain by making the sampling frequency satisfy the upper limit of the sampling frequency or reducing the frequency by resampling to meet the upper limit of the sampling frequency.
  • the error ensures that the error of the S domain is within the acceptable error range, which makes the digital control stable and meets the needs of specific applications.
  • the parameter identification method and system of the present invention can accurately identify the dynamic parameters and the static parameters by using the data sampling method and system as described above, especially after calculating the steady state value and then identifying the static parameters, so that the finally obtained static parameters The parameters are more accurate.
  • FIG. 1 is a schematic diagram of a first data sampling method according to an embodiment of the present invention.
  • FIG. 4 is a schematic diagram of a second data sampling method according to an embodiment of the present invention.
  • FIG. 5 is a schematic structural diagram of a data sampling system according to an embodiment of the present invention.
  • Figure 6 is a diagram showing the zero-order point on the S domain of an application system. detailed description
  • Figure 1 is a schematic diagram of the first data sampling method.
  • the physical quantity of the electrical signal y is first filtered by the analog low-pass filter 1 and then sampled by the analog-to-digital converter 2 to output a sampling sequence Y.
  • the sampling process is controlled by the sampling control signal.
  • the frequency/ s of the sampling control signal satisfies: >
  • ⁇ is the Z-domain error of the sampling system
  • the dish is the maximum error of the S-domain.
  • the mathematical description is the differential equation of the time domain; the differential equation of the time domain is transformed by the Laplace transform to the transfer function G(s) on the S domain; it has been proved in the cybernetics that G (s) has a strict correspondence with the differential equations in the time domain, without any distortion, so the error ⁇ of the time domain closely corresponds to the error S of the S domain, so there is an acceptable maximum error ⁇ at the time, which is acceptable.
  • the maximum error of the S domain is the differential equation of the time domain.
  • the signal is sampled in the time domain analog quantity, and the analog quantity is sampled and converted to obtain the sampled data.
  • the time domain is transformed into the ⁇ domain, and the differential equation in the time domain is transformed into the ⁇ domain. Function 3 ⁇ 42).
  • a minimum grid is generated, which is customarily called resolution, and half of the resolution is the truncation error ⁇ ; and the number is also represented in the computer, and the truncation error ⁇ is also generated.
  • There is an error ⁇ and when the hardware and software of the sampling system used are determined, the system's domain error ⁇ is also determined.
  • the number of bits in the ADC of the analog-to-digital converter determines the resolution of the analog-to-digital conversion, resulting in a truncation error of half the resolution.
  • the 8-bit ADC has a Sadc of 0.004 and a 10-bit of 0.001.
  • the value is in the computer. represented by a binary number, digit binary number also generates truncation error, 8 16-bit binary number is 2x10- 5, single-precision floating-point number is 10-6 8, 8, double precision floating point
  • the Z-domain error ⁇ can be analyzed by the truncation error ⁇ .
  • the sampling domain has a ⁇ ⁇ ⁇ ⁇
  • the ⁇ is the truncation error of the sampling system
  • the truncation error ⁇ of the sampling system In relation to ⁇ , the tube can be taken as the smaller of ⁇ ⁇ bcpu , or ⁇ .
  • T s is the sampling interval.
  • the S-domain grid will become larger and larger.
  • Half of the grid is the quantization error, so the error in the S domain will increase.
  • the same error analysis method as above can be used to derive the relationship between the time domain error ⁇ and the ⁇ domain error ⁇ , and the sampling frequency / s . Since the frequency domain, that is, the S domain, is used in practical applications, It will be described in detail, but it should also fall within the scope of protection of the present invention.
  • the sampling frequency / s tends to be high, and can not meet the requirement of less than the upper limit of the sampling frequency / 53 ⁇ 4 . If not processed, the S domain and time domain errors will be large, and the application requirements cannot be met. . Therefore, the present invention further proposes a second data sampling method and system, so that the sampled data can satisfy the upper limit of the sampling frequency, thereby improving the accuracy of the sampled data.
  • an embodiment of a second data sampling method includes the following steps:
  • the electrical signal y of the physical quantity is sampled to satisfy the sampling frequency of the Nyquist theorem, and the sampling sequence Y1 is obtained;
  • the Z-domain error ⁇ ⁇ ⁇ of the sampling system is a truncation error of the sampling system; the truncation error ⁇ of the sampling system is related to ⁇ , and the single ground can take ⁇ ⁇ , h cpu , or The smaller of ⁇ and ⁇ .
  • FIG. 5 is a schematic diagram of an embodiment of a data sampling system, including:
  • the sampling unit 3 is configured to sample the physical quantity of the electrical signal y to satisfy the sampling frequency of the Nyquist theorem, and obtain the sampling sequence Y1;
  • a low pass filter 4 configured to perform low pass filtering on the sampling sequence Y1 to obtain a filtered sampling sequence Y2;
  • the resampling unit 5 is configured to resample the filtered sample sequence Y2 to obtain the sample sequence Y, and then output to the output unit 6;
  • ⁇ is the area error of the sampling system, and the dish is the maximum error of the S domain;
  • Output unit 6 output sample sequence ⁇ .
  • the cutoff frequency of the low pass filter 4 is, / c ⁇ 0.5x, to prevent mixing when resampling.
  • the Z-domain error ⁇ ⁇ ⁇ of the sampling system is a truncation error of the sampling system; the truncation error ⁇ of the sampling system is related to ⁇ , and the single ground can take ⁇ ⁇ , h cpu , or The smaller of ⁇ and ⁇ .
  • the error of the S domain can be reduced to a predetermined range ⁇ .
  • the contradiction of / s3 ⁇ 4 ⁇ may occur in practical applications, resulting in failure to follow the first method.
  • f sl ⁇ f s ⁇ ⁇ Select the sampling frequency, which can be solved by the second data sampling method described above, that is, the method of low-pass filtering and re-sampling. Since in all the poles ⁇ p n ⁇ and zero ⁇ p ra ⁇ , the main pole is the main zero point P. m plays a leading role, the higher the frequency, the lower the pole or zero effect, so when the contradiction occurs, the high-frequency pole and zero point are filtered by the low-pass filter, and the sampling frequency is reduced by resampling, which can not only solve the problem.
  • the above sampling methods and systems have many applications.
  • the first application is parameter identification.
  • both dynamic and static parameters have been rarely successful. Now the sampling data is accurate. The error is within a given range, which is very beneficial for parameter identification.
  • Parameters include static parameters and dynamic parameters, static parameters are the coefficients of the system equation, and dynamic parameters are the coefficients of the system differential equation.
  • static parameters refer to the power flow equation.
  • Admittance matrix that is, resistance, reactance and susceptance of components such as transformers, wires, generators, etc.
  • Dynamic parameters refer to excitation time constant, speed regulation time constant, rotor time constant, etc.
  • a parameter identification method using the above two data sampling methods further comprising the step S5, comprising performing a dynamic parameter identification step S51 and/or a static parameter identification step S52 according to the sampling sequence ⁇ .
  • the step S51 of the dynamic parameter identification includes, according to the sampling sequence Y, the step S511 of identifying the order and parameters of the ARMAX model by using the identification method in the adaptive control, or the step S512 of identifying the parameter by the reference model method, or Kalman Step S513 of the filter identification parameter.
  • the static parameter identification step S52 includes a step S521 of first calculating a steady state value for the sampling sequence Y, and a parameter estimation method, a correlation coefficient method, a linear regression method, a linearizable linear regression method according to the steady state value, Or step S522 of identifying a static parameter by a least squares method;
  • Value step S521 includes:
  • the average value is calculated according to the sampling sequence Y as a steady state value
  • the forced component is calculated as the steady state value according to the sampling sequence Y.
  • the step S5211 determines whether the current sample value is in a steady state process or a transient process according to the t distribution or its binization, or the filter output result; the filter is a Kalman filter or an ⁇ filter.
  • the step S52135 calculates the forced component according to: ⁇ + ⁇ + ⁇ + ⁇ +...+ ⁇ , where ⁇ 3 ⁇ 4, ⁇ , ⁇ ⁇ , ⁇ are according to the derivative formula and the physical quantity y and The constant determined by the initial value of each derivative.
  • a parameter identification system comprising the above data sampling system, further comprising: a parameter identification unit, which is decomposed into a dynamic parameter identification unit 51 and a static parameter identification unit 52 according to application requirements, as shown in FIG.
  • the dynamic parameter identification unit 51 identifies the order and parameters of the ARMAX model, or the reference model estimation parameter, or the Kalman filter estimation parameter according to the sampling sequence Y, using the identification method in the adaptive control.
  • the static parameter identification unit 52 includes, including a steady state value calculation unit 521, Calculating a steady state value according to the sampling sequence Y; a static parameter estimating unit 522 for using a parameter estimation method, or a correlation coefficient method, or a linear regression method, or a linearizable linear regression method, or a minimum two according to a steady state value Multiply estimates static parameters.
  • FIG. 8 is a schematic diagram of a steady-state value calculation unit, the steady-state value calculation unit 521, comprising: a determination unit 5211, according to the sampling sequence ⁇ , determining whether the current sample value is in a steady state process or a transient process according to the t-distribution or its cylinderization;
  • Calculating a steady state unit 5212 configured to calculate an average value as a steady state value according to the sampling sequence Y when determining that the current sampling value is in a steady state process
  • the steady state unit 5213 is calculated to calculate the forced component as the steady state value according to the sampling sequence Y when it is determined that the current sample value is in the transient process.
  • the steady state unit 5213 is calculated, including:
  • n + l to determine n, where ⁇ is a constant close to 0;
  • the determining unit 52133 is configured to receive the data sequence "[ ⁇ , and determine that the current value is in a steady state process or a transient process according to the t distribution or its cylinderization;
  • the determining unit 52133 the t distribution criterion:
  • t(k) is the t distribution with a degree of freedom of k.
  • the criterion can be further reduced to:
  • ⁇ is a given constant
  • Xe is the nominal value of the physical quantity X.
  • the ⁇ is between 0.1% and 10%.
  • the average and standard deviation 3 ⁇ 4 calculations can be performed as follows:
  • the filtering algorithm can also be used to judge whether the current value x k is in a steady state process or a transient process, specifically: performing Kalman filtering on ⁇ ⁇ to obtain each component of the state variable; and determining whether each component of the state variable exceeds the corresponding setting. Value; if there is one component or several components exceeding the set value, it is judged that the current value is in the transient process; otherwise, if it is not exceeded, it is judged that the current value is in the steady state process.
  • the determining unit 5211 according to the sampling sequence Y, determining the current sampling value in the steady state process or the transient process according to the t distribution or the sizing thereof, and the steps of calculating the average value and the variance are similar to the determining unit 52133, I won't go into details here.
  • the data sampling method and system of the invention reduces the truncation error in the digital sampling and digital system to the S domain and the time domain by making the sampling frequency satisfy the upper limit of the sampling frequency or reducing the frequency by resampling to meet the upper limit of the sampling frequency.
  • the error makes the final S-domain error within the acceptable maximum error, meeting the needs of the specific application.
  • the parameter identification method and system of the invention adopts the data sampling method and system as described above, and can accurately identify the dynamic parameters and the static parameters, especially after calculating the steady state value, and then identifying the static parameters, so that the finally obtained static parameters More accurate. Since telemetry collects local digital measurement data by communication, which is consistent with the principle and nature of data obtained by in-situ digital measurement, the data sampling method and system of the present invention are equally applicable to telemetry methods and systems, and should be equally protected.

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Physics & Mathematics (AREA)
  • Power Engineering (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Signal Processing (AREA)
  • Mathematical Analysis (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computing Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Feedback Control In General (AREA)
  • Emergency Protection Circuit Devices (AREA)
  • Testing And Monitoring For Control Systems (AREA)
  • Analogue/Digital Conversion (AREA)

Abstract

本发明涉及一种数据采样方法：对物理量按照采样频率fs≤fsh进行采样，fsh=ζmax/ε为采样频率上限；以及一种数据采样方法及系统：对物理量以满足奈奎斯特定理的采样频率进行采样，对得到的采样序列先进行低通滤波，再进行重抽样，重抽样的频率fs≤ζmax/ε；其中ε为采样系统的Ζ域误差，ζmax为S域的最大误差。本发明还涉及一种先采用上述数据采样方法及系统得到采样数据，再利用采样数据进行动态和/或静态参数辨识的参数辨识方法及系统。本发明的数据采样方法及系统，以及参数辨识方法及系统，解决了采样数据误差较大、数字控制失稳、和参数辨识失败的技术难题。

Description

说 明 书 数据采样方法与系统及其在参数辨识中的应用方法与系统
技术领域
本发明涉及计算机应用技术， 尤其是采样、 参数辨识和数字控制技术。
背景技术
将模拟量转换为数字量的过程， 称之为采样， 早在 1924年 Nyquist就给 出了采样频率的下限 = 2· ， 即奈奎斯特采样定理； 但是在模拟控制器下稳 定的系统， 经过采样后成为数字控制系统， 尽管遵循了采样定理， 反而系统不 再稳定； 利用采样数据进行参数辨识也鲜有成功。
面对不成功的情况， 按照采样定理， 通常认为采样频率越快越好， 技术人 员往往会选择加快采样频率， 但是， 实际情况证明， 采样频率再怎么快， 也不 能解决问题。 于是， 数字控制和参数辨识需要进行不断修改、 调试、 再修改、 再调试的不断反复摸索。
发明内容
本发明的目的就是要解决上述问题，由稳定的模拟控制器进行数字采样后 转换为数字控制器后仍然稳定并特性相当，并且依据采样数据能够辨识出准确 的参数。
发明人经过长达三十年的悉心研究，通过反复试验总结、再到理论推导分 析， 发现上述问题的原因不是辨识算法、数字控制方法造成的， 而是采样数据 造成的； 采样数据问题不是采样太慢造成的， 而是采样太快造成的； 采样频率 并不是像本领域技术人员通常认为的那样越快越好，而是存在一个采样频率上 限 ， 采样频率超过上限是造成采样数据误差较大、 数字控制失稳、 和参数 辨识失败的根源所在。
本发明提出一种数据采样方法， 包括以下步骤： 51、 对物理量的电信号 y进行模拟低通滤波， 得到滤波后的电信号；
52、对滤波后的电信号按照采样频率/进行采样， 所述采样频率 /满足奈 奎斯特定理， 即 > ^；
53、 输出采样序列 Y;
其中， 所述步骤 S2 的采样频率 /s还满足 fs≤ fsh=U； 其中 ε为采样系 统的 Ζ域误差， ζ S域的最大误差。
其中， 所述采样系统的 Ζ域误差 ε≥δ , 所述 δ为采样系统的截断误差； 所 述采样系统的截断误差 δ与 和 有关， 为模数转换器的分辨率的一 半， 为数值在计算机中表示位数的截断误差。
进一步的， 所述采样系统的 Ζ域误差 ε = · δ , 所述 为大于 1的保险系 数， 其取值由采样时的随机干扰水平、数值在采样系统中的运算误差或其组合 决定。
其中， 所述 S域的最大误差 ζ , 是由具体应用决定的 S域可以接受的最 大误差， 例如， 由应用系统 S域零极点的可接受误差决定的 ζ , 或者由应用 系统时域微分方程可以接受的的最大误差0 对应确定的 ζ 为了缩小体积、 降低成本， 实际应用中采样频率往往较高， 不满足采样频 率的上限 ， 因此本发明还提出了第二种数据采样方法， 其包括以下步骤：
51、对物理量的电信号 y以满足奈奎斯特定理的采样频率进行采样，得到 采样序列 Y1 ;
52、 对采样序列 Y1进行数字低通滤波， 得到滤波后的采样序列 Y2;
53、 对滤波后的采样序列 Y2进行重抽样， 得到采样序列 Y; 所述重抽样 的频率为 ， 所述/„≤ =ζ^ /ε ; 其中 ε为采样系统的 Ζ域误差， ζ匪为 S域 的最大误差；
54、 输出采样序列 Y
进一步的， 所述低通滤波的截止频率为 ， 所述 /c≤0.5x ， 以防止重抽 样时产生混频。
经过重抽样后， 由于重抽样的频率 frs能够满足小于采样频率的上限 fsh , 所以能够将 S域的误差降低到预先设定的范围 ζ 内。 实际应用有中可能出现 /s¾ < /s,的矛盾， 导致无法按照第一种方法的 fsl < fs≤ ^选取采样频率， 通过上述第二种数据采样方法， 即经过低通滤波和 重抽样的方法能够解决该矛盾。 由于在所有的极点 {pn}和零点 {pra}中， 主极点 主零点 P。m起主导作用， 频率越高的极点或零点作用越小， 所以， 当发 生/ s¾ < fsl的矛盾时， 通过所述低通滤波滤除高频极点和零点、 重抽样降低采样 频率， 不仅能够解决所述矛盾， 还能够保证主极点 Prn^P主零点 P。m的精度， 其中主极点或主零点是 S域上靠近原点的极点或零点。 本发明还提出了一种数据采样系统， 其包括：
采样单元，用于对物理量的电信号 y以满足奈奎斯特定理的采样频率进行 采样， 得到采样序列 Y1 ;
低通滤波器， 用于对采样序列 Y1进行数字低通滤波， 得到滤波后的采样 序歹' J Y2;
重抽样单元，用于对滤波后的采样序列 Y2进行重抽样，得到采样序列 Y; 所述重抽样的频率为 frs , 所述 frs≤ fsh= I ； 其中 ε为采样系统的 Ζ域误差， 为 S域的最大误差；
进一步的， 所述低通滤波器的截止频率为 ， 所述 y≤o.5x ， 以防止重 抽样时产生混频。
本发明通过重抽样降低采样频率， 重抽样的频率/„≤/5¾= ζ /ε , 从而大 大降低了因数字采样和数字系统中的截断误差给时域、 S域造成更大的误差， 使 S域和时域的误差满足应用要求。 数字采样要么用于测量，要么用于辨识参数或数字控制。 参数包括静态参 数和动态参数，静态参数即系统方程的系数，动态参数即系统微分方程的系数。 对于电力系统， 静态参数是指潮流方程的导纳矩阵， 即变压器、 导线、 发电机 等元件的电阻、 电抗和电纳， 动态参数是指励磁时间常数、 调速时间常数、 转 子时间常数等。 为了更好地理解所述采样方法的实用性， 本发明给出所述采样 方法及系统在参数辨识方面的应用。
一种采用上述两种数据采样方法的参数辨识方法，其特征在于还包括步骤 S5 , 包括根据采样序列 Y进行动态参数辨识步骤 S51和 /或静态参数辨识步骤 S52。
其中， 所述动态参数辨识的步骤 S51包括， 根据采样序列 Y, 采用自适应 控制中的辨识方法辨识 ARMAX模型的阶数和参数的步骤 S511、 或参考模型 法辨识参数的步骤 S512、 或卡尔曼滤波器辨识参数的步骤 S513。
其中， 所述静态参数辨识步骤 S52包括， 对采样序列 Y先计算稳态值的 步骤 S521 , 再根据稳态值采用参数估计法、 相关系数法、 线性回归法、 可线 性化的线性回归法、 或最小二乘法辨识静态参数的步骤 S522; 所述计算稳态 值步骤 S521包括：
55211、 依据采样序列 Y, 判断当前采样值处于稳态过程还是暂态过程；
55212、 当判定当前采样值处于稳态过程， 根据采样序列 Y计算平均值作 为稳态值；
55213、 当判定当前采样值处于暂态过程， 根据采样序列 Y计算强迫分量 作为稳态值。
其中， 所述步骤 S5211 , 依据 t分布或其筒化、 或滤波器输出结果判断当 前采样值处于稳态过程或暂态过程； 所述滤波器为卡尔曼滤波器或 αβγ 滤波 器。 例如， 对采样序列 Υ进行卡尔曼滤波或 αβγ滤波， 得到状态变量的各分 量； 再判断状态变量的各分量是否超出对应的设定值； 若有一个分量或几个分 量超出设定值， 则判断当前采样值处于暂态过程； 否则， 均未超出， 则判断当 前采样值处于稳态过程。 所述卡尔曼滤波器至少为 1 阶， 或更高阶， 例如 2 阶、 3阶等， 本发明不做限定。
其中， 所述步骤 S5213 , 按照以下步骤计算强迫分量：
S52131 , 根据采样序列 Υ计算 y 的 n+1 阶导数， η=0,1,...， 直到满足 1
d y
n+l ≤ε , 从而确定 n, 其中 ε为接近于 0的常数;
dt
552132, 做变量替换; c = ^Z , 由采样序列 Y计算出数据序列 { }；
dtn
552133 ,对于数据序列 { } ,判断当前数值 处于暂态过程还是稳态过程；
552134, 当判定当前数值 处于稳态过程， 根据 { }计算平均值 作为稳 态值 , 进行步骤 S52135; 当判定当前数值 处于暂态过程， 令 k=l ; S52135、 对 x的稳态值 做积分反变换， 计算出 y的强迫分量 。
其中， 所述步骤 S52135 , 按照 = α。 + αι · + · 2 + α3 · 3 +〜 + α„· "计算强迫 分量 , 其中， 丄 ·¾ , ^ , α α^ ,α^是按照导数公式和物理量 y及其各阶 导数的初值确定的常数。
其中， 所述步骤 S52133 ,依据 t分布或其筒化、或滤波器输出结果判断当 前数值 处于稳态过程或暂态过程； 所述滤波器为卡尔曼滤波器或 αβγ 滤波 器。 例如， 对采样序列 { }进行卡尔曼滤波或 αβγ滤波， 得到状态变量的各分 量； 再判断状态变量的各分量是否超出对应的设定值； 若有一个分量或几个分 量超出设定值， 则判断当前数值 处于暂态过程； 否则， 均未超出， 则判断当 前数值 处于稳态过程。 所述卡尔曼滤波器至少为 1 阶， 或更高阶， 例如 2 阶、 3阶等， 本发明不做限定。 一种包括上述数据采样系统的参数辨识系统， 其还包括： 参数辨识单元， 用于根据采样序列 Υ进行动态参数辨识和 /或静态参数辨识。
其中， 需要辨识动态参数时， 所述参数辨识单元为动态参数辨识单元， 用 于根据采样序列 Υ, 采用自适应控制中的辨识方法辨识 ARMAX模型的阶数 和参数、 或参考模型法辨识参数、 或卡尔曼滤波器辨识参数。
其中，需要进行静态参数辨识时，所述参数辨识单元为静态参数辨识单元， 包括稳态值计算单元， 用于根据采样序列 Υ计算稳态值； 静态参数估计单元， 用于根据稳态值， 采用参数估计法、 相关系数法、 线性回归法、 可线性化的线 性回归法、 或最小二乘法辨识静态参数。
所述稳态值计算单元， 包括：
判断单元， 根据采样序列 Υ, 依据 t分布或其筒化、 或滤波器输出结果判 断当前采样值处于稳态过程或暂态过程； 所述滤波器为卡尔曼滤波器或 αβγ 滤波器；
计算稳态单元 Α, 用于当判定当前采样值处于稳态过程时，根据采样序列 Υ计算平均值作为稳态值；
计算稳态单元 Β , 用于当判定当前采样值处于暂态过程时， 根据采样序列 Υ计算强迫分量作为稳态值。 其中， 计算稳态单元 B, 包括：
求导单元， 用于根据采样序列 Y计算 y的 n+1阶导数， η=0,1,...， 直到满 1
足 d y
n+l ≤ε , 从而确定 n, 其中 ε为接近于 0的常数;
dt
替换单元， 用于做变量替换 c = 2, 由采样序列 Y计算出数据序列 { }； dtn
判断单元， 用于接收数据序列" [ }, 判断当前数值 处于暂态过程还是稳 态过程；
平均单元， 用于当判定当前数值 处于稳态过程， 根据 { }计算平均值 作为稳态值 , 进入还原单元； 当判定当前数值 处于暂态过程， 令 k=l, 到 出口；
还原单元，按照公式 ^+^ + ^^ + ^^+…+ ^^"计算强迫分量 , 其 中， 《„=丄 ·¾, , 是按照导数公式和物理量 y及其各阶导数的初值 确定的常数。
本发明的数据采样方法， 通过令采样频率满足采样频率上限，或者通过重 抽样降低频率以满足采样频率上限，来降低因数字采样和数字系统中的截断误 差给 S域或时域造成更大的误差， 保障 S域的误差在可以接受的误差 的范 围内， 使数字控制稳定， 满足了具体应用的需求。
本发明的参数辨识方法及系统，由于采用了如上所述的数据采样方法及系 统， 能够准确辨识出动态参数和静态参数，特别是计算出稳态值后再辨识静态 参数， 使最后得到的静态参数更为准确。
附图说明
图 1为本发明实施例第一种数据采样方法的示意图。
图 2为 /s=10Hz时， S域和 Z域的栅格受采样频率的影响示意图。
图 3为 /s=30Hz时， S域和 Z域的栅格受采样频率的影响示意图。
图 4为本发明实施例第二种数据采样方法的示意图。
图 5为本发明实施例数据采样系统的结构示意图。
图 6为假定某个应用系统的 S域上零级点示意图。 具体实施方式
下面详细描述本发明的实施例， 所述实施例的示例在附图中示出， 其中自 始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元 件。 下面通过参考附图描述的实施例是示例性的， 仅用于解释本发明， 而不能 解释为对本发明的限制。
图 1为第一种数据采样方法的示意图，将物理量的电信号 y先经过模拟低 通滤波器 1滤波，再经过模数转换器 2采样后输出采样序列 Y, 采样过程受采 样控制信号的控制， 采样控制信号的频率/ s除了满足 > ^外， 还满足：
f - f ζ
其中， ε为采样系统的 Z域误差， ζ皿为 S域的最大误差。
对于模拟量随时间变化的动态过程， 用数学描述是时域的微分方程； 时域 的微分方程经过拉普拉斯变换到 S域上的传递函数 G(s); 控制论中已经证明， G(s)与时域的微分方程有严密的对应关系， 没有任何失真， 所以时域的误差 σ 严密对应着 S域的误差 ζ , 所以当时域存在可以接受的的最大误差 σ , 对应 有可以接受的 S域的最大误差 皿。
对时域模拟量进行信号采样， 将模拟量经过模数转换采样后获得采样数 据， 数字自动控制理论中， 认为将时域变换到 Ζ域， 时域上的微分方程变换 为 Ζ域上的传递函数 ¾2)。2域中因为数字计算机的位数有限会产生最小栅格， 习惯上称为分辨率， 分辨率的一半为截断误差 δ； 并且数字在计算机中表示时 也会产生截断误差 δ , 于是， Ζ域存在误差 ε , 并且当采用的采样系统的软硬 件确定后， 系统的 Ζ域误差 ε也会随之确定。
例如， 模数转换器 ADC的位数决定模数转换的分辨率， 产生的截断误差 为分辨率的一半, 8位 ADC的 S adc为 0.004、 10位的 为 0.001 ;再例如， 数值在计算机中用二进制数表示， 二进制数的位数也产生截断误差 , 16 位二进制数的8 为 2x10— 5、单精度浮点数的8 为 10— 6、双精度浮点数的8 经过误差分析， 由截断误差 δ可以分析出 Z域误差 ε , 具体来说， 所述采 样系统的 Ζ域误 ε≥δ , 所述 δ为采样系统的截断误差； 所述采样系统的截断误 差 δ与 和 δ 有关， 筒单地， 可以取 δα bcpu , 或 和 δ 中的较小者。
进一步的， 所述采样系统的 Ζ域误差 ε= ·δ , 所述 为大于 1的保险系 数，取决于采样时的随机干扰水平、 以及数值在采样系统中的运算误差或其组 合。
经过发明人悉心研究， 发现系统的 Ζ域误差 ε和 S域的误差 ζ存在着一定 关系， 具体来说：
S域到 Ζ域的变换为 z = e s7 rs为采样周期； 反之 Z域到 S域的变换为 s = _· ln(z) i己 s = JC+ 和 z = , 有：
Figure imgf000010_0001
x = fs-Hr), = ·Θ (2)
Figure imgf000010_0002
由于 S域的靠近虚轴的左半平面映射到 Z域接近半径为 1的圆环上， r-1, ··· dx^ fs · dr (3) dy = fs-dQ (4) 于是， S域 x方向的误差 ζχ和 y方向的误差 分别为：
ζ Λ·ε (5) ζ,=Λ·ε (6) 由式 (5)和式 (6), 得， S域的误差 ζ为：
ζ-Λ·ε (7)
Ts为采样间隔。
Figure imgf000010_0003
所以， 为了保证 S域的误差在可以接受的范围 皿内， 采样频率/应当满 足:
fs≤fsh=U (8) 即采样频率的上限/ 5¾=ζ
参考图 2所示， 为 /s=10Hz时， S域和 Z域的栅格受采样频率的影响示意 图。 当 /s=10Hz时， 将 S域的 (-10 ~ 0, -31.4 ~ 31.4)映射到 Z域的半径 (0.37 ~
1)圆环上。
参考图 3所示， 为/ =30Hz时， S域和 Z域的栅格受采样频率的影响示意 图。 当 /s=30Hz时， 将 S域的 (-10 ~ 0, -94.2 ~ 94.2)映射到 Z域的半径 (0.72 ~ 1)圆环上。
对比图 2和图 3 , 可见 /s由 10Hz加快 3倍到 30Hz, S域扩大了 6.67倍， 而 Z域缩小了 1.92倍， 也就是说， 随着采样频率 的增加 Z域的一小块面积 将代表 S域越来越大的面积； 图 3中的 Z域栅格明显比图 2的 Z域栅格小、 图 3中的 S域栅格明显比图 2的 S域栅格大，可以推断， 随着采样频率 /s的加 快， Z域的栅格会越来越小， 反而， 对应的 S域栅格会越来越大。 当采样系统 的软硬件一定后， Z域误差 ε和栅格也随之确定， 这时随着采样频率/ 5的加快，
S域栅格会变的越来越大。 栅格的一半即量化误差， 所以 S域的误差会加大。
因此，发明人不同于通常观念中认为的采样频率越快越好， 而是创造性的 提出了采样频率有上限， 为 =ζ /ε。 同样道理， 也可以通过与上述一样的 误差分析办法， 推导出时域误差 σ与 Ζ域误差 ε、 以及采样频率 /s的关系， 由 于实际应用中多采用频域， 即 S域， 所以这里不再详述， 但也应属于本发明的 保护范围内。
在实际应用中， 为了缩小体积、 降低成本， 采样频率/ s往往较高， 不能够 满足小于采样频率上限/ 5¾的要求， 如果不处理将导致 S域和时域误差很大， 不能满足应用要求。所以本发明进一步提出了第二种数据采样方法及系统，使 采样后的数据能够满足采样频率上限， 从而提高采样数据的准确性。
如图 4所示为第二种数据采样方法的实施例， 其包括以下步骤：
51、对物理量的电信号 y以满足奈奎斯特定理的采样频率进行采样，得到 采样序列 Y1 ;
52、 对采样序列 Y1进行低通滤波， 得到滤波后的采样序列 Y2;
53、 对滤波后的采样序列 Y2进行重抽样， 得到采样序列 Y; 所述重抽样 的频率为 ， 所述/„≤ =ζ^ /ε ; 其中 ε为采样系统的 Ζ域误差， ζ匪为 S域 的最大误差；
54、 输出采样序列 Y。 进一步的， 所述低通滤波的截止频率为 ， 所述 y ≤o.5x ， 以防止重抽 样时产生混频。
其中， 所述采样系统的 Z域误 ε≥δ , 所述 δ为采样系统的截断误差； 所述 采样系统的截断误差 δ与 和 δ 有关， 筒单地， 可以取 δα 、 h cpu , 或 和 δ 中的较小者。
进一步的， 所述采样系统的 Ζ域误差 ε = · δ , 所述 为大于 1的保险系 数，取决于采样时的随机干扰水平、 以及数值在采样系统中的运算误差或其组 合。 如图 5所示为一种数据采样系统的实施例示意图， 其包括：
采样单元 3 , 用于对物理量的电信号 y以满足奈奎斯特定理的采样频率进 行采样， 得到采样序列 Y1 ;
低通滤波器 4, 用于对采样序列 Y1进行低通滤波， 得到滤波后的采样序 列 Y2;
重抽样单元 5 , 用于对滤波后的采样序列 Y2进行重抽样， 得到采样序列 Y后， 到输出单元 6; 所述重抽样的频率为 ， 所述/„≤ =ζΜ∞/ε ; 其中 ε为 采样系统的 Ζ域误差， ζ皿为 S域的最大误差；
输出单元 6, 输出采样序列 Υ。
进一步的， 所述低通滤波器 4的截止频率为 ， 所述 /c≤0.5x ， 以防止 重抽样时产生混频。
其中， 所述采样系统的 Z域误 ε≥δ , 所述 δ为采样系统的截断误差； 所 述采样系统的截断误差 δ与 和 δ 有关，筒单地，可以取 δα 、 h cpu ,或 和 δ 中的较小者。
进一步的， 所述采样系统的 Ζ域误差 ε = · δ , 所述 为大于 1的保险系 数，取决于采样时的随机干扰水平、 以及数值在采样系统中的运算误差或其组 合。
经过重抽样后， 由于重抽样的频率 frs能够满足小于采样频率的上限 fsh , 所以能够将 S域的误差降低到预先设定的范围 ζ 内。
此外， 实际应用中可能出现/ s¾ < 的矛盾， 导致无法按照第一种方法的 fsl < fs≤ ^选取采样频率， 通过上述第二种数据采样方法， 即经过低通滤波和 重抽样的方法能够解决该矛盾。 由于在所有的极点 {pn}和零点 {pra}中， 主极点 主零点 P。m起主导作用， 频率越高的极点或零点作用越小， 所以， 当发 生 < 的矛盾时， 通过所述低通滤波滤除高频极点和零点、 重抽样降低采样 频率， 不仅能够解决所述矛盾， 还能够保证主极点 Prn^P主零点 P。m的精度， 其中主极点或主零点是 S域上靠近原点的极点或零点。如图 6所示， X为极点、 0为零点， 区域 Q1中有主极点 Prn^主零点 P。m, 区域 Q2中为高频极点和高 频零点， Q1 中的零极点比 Q2 中的零极点重要， 一旦发生/ s¾ < /s,的矛盾， 应 当采用图 5中的低通滤波器滤除图 6中的 Q2中靠上， ω较大和靠左， X负值 的绝对值较大的高频极点。
本发明通过重抽样降低采样频率， 重抽样的频率/„≤/5¾= ζ /ε , 从而大 大降低了因数字采样和数字系统中的截断误差给 S域或时域造成更大的误差， 使最后的 S域和时域的误差满足应用要求。 上述采样方法和系统有多种应用， 首个应用就是参数辨识， 以往无论是动 态参数还是静态参数， 都鲜有成功， 现在由于采样数据准确了， 误差在给定范 围内， 非常有利于参数辨识。 参数包括静态参数和动态参数， 静态参数即系统 方程的系数， 动态参数即系统微分方程的系数。 对于电力系统， 静态参数是指 潮流方程的导纳矩阵， 即变压器、 导线、 发电机等元件的电阻、 电抗和电纳， 动态参数是指励磁时间常数、 调速时间常数、 转子时间常数等。
一种采用上述两种数据采样方法的参数辨识方法， 还包括步骤 S5 , 包括 根据采样序列 Υ进行动态参数辨识步骤 S51和 /或静态参数辨识步骤 S52。
其中， 所述动态参数辨识的步骤 S51包括， 根据采样序列 Y, 采用自适应 控制中的辨识方法辨识 ARMAX模型的阶数和参数的步骤 S511、 或参考模型 法辨识参数的步骤 S512、 或卡尔曼滤波器辨识参数的步骤 S513。
其中， 所述静态参数辨识步骤 S52包括， 对采样序列 Y先计算稳态值的 步骤 S521 , 再根据稳态值采用参数估计法、 相关系数法、 线性回归法、 可线 性化的线性回归法、 或最小二乘法辨识静态参数的步骤 S522; 所述计算稳态 值步骤 S521包括：
55211、 依据采样序列 Y判断当前采样值处于稳态过程还是暂态过程；
55212、 当判定当前采样值处于稳态过程， 根据采样序列 Y计算平均值作 为稳态值；
55213、 当判定当前采样值处于暂态过程， 根据采样序列 Y计算强迫分量 作为稳态值。
其中， 所述步骤 S5211, 依据 t分布或其筒化、 或滤波器输出结果判断当 前采样值处于稳态过程或暂态过程； 所述滤波器为卡尔曼滤波器或 αβγ 滤波
≤ε , 从而确定 n, 其中 ε为接近于 0的常数;
Figure imgf000014_0001
552132, 做变量替换; c = ^， 由采样序列 Y计算出数据序列 { }；
dtn
552133,对于数据序列 { },判断当前数值 处于暂态过程还是稳态过程； S52134, 当判定当前数值 处于稳态过程， 根据 { }计算平均值 作为稳 态值 , 进行步骤 S52135; 当判定当前数值 处于暂态过程， 令 k=l;
S52135, 对 X的稳态值 做积分反变换， 计算出 y的强迫分量 。
其中， 所述步骤 S52135, 按照 ：^ + ^ + ^^ + ^^+…+ ^^计算强迫 分量 , 其中， 丄 ·¾, ^,α α^,α^是按照导数公式和物理量 y及其各阶 导数的初值确定的常数。 一种包括上述数据采样系统的参数辨识系统， 其还包括： 参数辨识单元， 根据应用要求， 分解为动态参数辨识单元 51和静态参数辨识单元 52, 参见图 7所示。
图 7中， 所述动态参数辨识单元 51, 根据采样序列 Y, 采用自适应控制 中的辨识方法辨识 ARMAX模型的阶数和参数、 或参考模型法估计参数、 或 卡尔曼滤波器估计参数。
图 7中， 所述静态参数辨识单元 52包括， 包括稳态值计算单元 521, 用 于根据采样序列 Y计算稳态值； 静态参数估计单元 522, 用于根据稳态值， 采 用参数估计法、 或相关系数法、 或线性回归法、 或可线性化的线性回归法、 或 最小二乘法估计静态参数。
图 8为稳态值计算单元示意图， 所述稳态值计算单元 521, 包括： 判断单元 5211, 根据采样序列 Υ,依据 t分布或其筒化判断当前采样值处 于稳态过程或暂态过程；
计算稳态单元 5212, 用于当判定当前采样值处于稳态过程时， 根据采样 序列 Y计算平均值作为稳态值；
计算稳态单元 5213, 用于当判定当前采样值处于暂态过程时， 根据采样 序列 Y计算强迫分量作为稳态值。
其中， 计算稳态单元 5213, 包括：
求导单元 52131, 用于根据采样序列 Y计算 y的 n+1阶导数， η=0,1,...，
7/2 + 1
直到满足 d y
≤ε ,
n+l 从而确定 n, 其中 ε为接近于 0的常数;
dt
替换单元 52132, 用于做变量替换 c = 2, 由采样序列 Y计算出数据序 dtn
列 { };
判断单元 52133, 用于接收数据序列" [ }, 依据 t分布或其筒化判断当前 数值 处于稳态过程或暂态过程；
平均单元 52134, 用于当判定当前数值 处于稳态过程， 根据 { }计算平 均值 作为稳态值 ,进入还原单元 52135; 当判定当前数值 处于暂态过程， 令 k=l, 到出口；
还原单元 52135, 按照 二 ^？^^？ + ^？^^ + ^^^+…+ ^?^"计算强迫分量 , 其中， =丄.¾, ,Α, ,···,^^是按照导数公式和物理量 y及其各阶导数的初 值确定的常数。
其中， 判断单元 52133, 所述 t分布判据：
服从 分布
F
0, 如果^ ^不服从 分布
F=l时， 处于稳态过程； F=0时， 处于暂态过程； 其中， 为平均值、 ¾ 为标准方差、 t(k)是自由度为 k的 t分布。
所述判据可以筒化为：
如果 | - |≤G-4
Figure imgf000016_0001
如果 | — |>G- 其中， 为平均值、 ¾为标准方差、 G为给定常数。所述 G位于 2.5~15之间。
所述判据还可以进一步筒化为：
如果 | _ |≤δ ·
Figure imgf000016_0002
如果 \xk -χΑ>δ■ xe
其中， 为平均值， δ为给定常数， Xe为物理量 X的额定值。所述 δ位于 0.1% ~ 10%之间。
所述平均值 和标准方差 ¾计算可以按照以下方法进行：
当 F=l变为 F = 0时，即当前数值 由稳态过程进入暂态过程时，令 k=l;
Figure imgf000016_0003
还可以采用滤波算法判断当前数值 xk处于稳态过程或暂态过程， 具体为： 对 { }进行卡尔曼滤波， 得到状态变量的各分量； 再判断状态变量的各分量是 否超出对应的设定值； 若有一个分量或几个分量超出设定值， 则判断当前数值 处于暂态过程； 否则， 均未超出， 则判断当前数值 处于稳态过程。
其中， 所述判断单元 5211, 根据采样序列 Y, 依据 t分布或其筒化判断当 前采样值处于稳态过程或暂态过程的步骤， 以及计算平均值和方差的步骤， 与 判断单元 52133类似， 这里不再赘述。
本发明数据采样方法及系统，通过令采样频率满足采样频率上限，或者通 过重抽样降低频率以满足采样频率上限，来降低因数字采样和数字系统中的截 断误差给 S域和时域造成更大的误差， 使最后的 S域的误差在可以接受的最 大误差 的范围内， 满足了具体应用的需求。
本发明参数辨识方法及系统， 由于采用了如上所述的数据采样方法及系 统， 能够准确辨识出动态参数和静态参数，特别是计算出稳态值后再辨识静态 参数， 使最后得到的静态参数更为准确。 由于遥测通过通信收集就地的数字测量数据 ,与就地数字测量获得数据的 原理和性质一致， 所以， 本发明的数据采样方法及系统同样适用于遥测方法及 系统， 理应受到同样的保护。
以上依据图式所示的实施例详细说明了本发明的构造、 特征及作用效果， 以上所述仅为本发明的较佳实施例，但本发明不以图面所示限定实施范围， 凡 是依照本发明的构想所作的改变，或修改为等同变化的等效实施例，仍未超出 说明书与图示所涵盖的精神时， 均应在本发明的保护范围内。

Claims

权 利 要 求 书
1. 一种数据采样方法， 包括以下步骤：
51、 对物理量的电信号 y进行模拟低通滤波， 得到滤波后的电信号；
52、对滤波后的电信号按照采样频率/进行采样， 所述采样频率 /满足奈 奎斯特定理， 即 > ^；
53、 输出采样序列 Y;
其特征在于： 所述步骤 S2 的采样频率 /s满足 fs≤ fsh=U； 其中 ε为采 样系统的 Ζ域误差， ζ皿为 S域的最大误差。
2. 一种数据采样方法， 包括以下步骤：
51、对物理量的电信号 y以满足奈奎斯特定理的采样频率进行采样，得到 采样序列 Y1 ;
其特征在于： 还包括以下步骤：
52、 对采样序列 Y1进行数字低通滤波， 得到滤波后的采样序列 Y2;
53、 对滤波后的采样序列 Y2进行重抽样， 得到采样序列 Y; 所述重抽样 的频率为 ， 所述/„≤ =ζ^ /ε ; 其中 ε为采样系统的 Ζ域误差， ζ匪为 S域 的最大误差；
54、 输出采样序列 Y。
3. 如权利要求 2所述的数据采样方法， 其特征在于： 所述数字低通滤波 的截止频率为 , 所述 ≤ 0.5 X 。
4. 如权利要求 1或 2所述的数据采样方法， 其特征在于： 所述采样系统 的 Ζ域误差 ε≥δ , 所述 δ为采样系统的截断误差。
5. 如权利要求 4所述的数据采样方法， 其特征在于： 所述采样系统的截 断误差 δ与 5^和8 有关， 其中 为模数转换器分辨率的一半， δ 为数值 在计算机中表示位数的截断误差。
6. 如权利要求 4所述的数据采样方法， 其特征在于： 所述采样系统的 Ζ 域误差 ε = ^ · δ , 所述 为大于 1的保险系数， 其取值由采样时的随机干扰水 平、 数值在采样系统中的运算误差或其组合决定。
7. 如权利要求 1或 2所述的数据采样方法， 其特征在于： 所述 S域的最 大误差 由具体应用决定， 包括： 由应用系统 S域零极点的可接受误差决定 的 , 或者由应用系统时域微分方程可接受的的最大误差 σ 对应确定的 ζ
8. 一种数据采样系统， 其特征在于包括：
采样单元，用于对物理量的电信号 y以满足奈奎斯特定理的采样频率进行 采样， 得到采样序列 Y1 ; 其特征在于还包括：
低通滤波器， 用于对采样序列 Y1进行数字低通滤波， 得到滤波后的采样 序歹' J Y2;
重抽样单元，用于对滤波后的采样序列 Y2进行重抽样，得到采样序列 Y; 所述重抽样的频率为 frs , 所述 frs≤ fsh=D ； 其中 ε为采样系统的 Ζ域误差， 为 S域的最大误差。
9. 如权利要求 8所述的数据采样系统， 其特征在于： 所述低数字通滤波 器的截止频率为 ， 所述 /≤0·5χ 。
10. 一种采用如权利要求 1或 2所述的数据采样方法的参数辨识方法， 其 特征在于还包括步骤 S5 , 包括根据采样序列 Υ进行动态参数辨识的步骤 S51 和 /或静态参数辨识的步骤 S52
11. 如权利要求 10所述的参数辨识方法， 其特征在于： 所述动态参数辨 识的步骤 S51 包括， 根据采样序列 Y, 采用自适应控制中的辨识方法辨识 ARMAX模型的阶数和参数的步骤 S511、或参考模型法辨识参数的步骤 S512 或卡尔曼滤波器辨识参数的步骤 S513
12. 如权利要求 10所述的参数辨识方法， 其特征在于： 所述静态参数辨 识步骤 S52包括， 对采样序列 Y先计算稳态值的步骤 S521 , 再根据稳态值采 用参数估计法、 相关系数法、 线性回归法、 可线性化的线性回归法、 或最小二 乘法辨识静态参数的步骤 S522; 所述计算稳态值的步骤 S521包括：
55211、 依据采样序列 Y判断当前采样值处于稳态过程还是暂态过程；
55212、 当判定当前采样值处于稳态过程， 根据采样序列 Y计算平均值作 为稳态值；
55213、 当判定当前采样值处于暂态过程， 根据采样序列 Y计算强迫分量 作为稳态值。
13. 如权利要求 12所述的参数辨识方法， 其特征在于： 所述步骤 S5211, 依据 t分布或其筒化、或滤波器输出结果判断当前采样值处于稳态过程或暂态 过程； 所述滤波器为卡尔曼滤波器或 αβγ滤波器。
14. 如权利要求 12所述的参数辨识方法， 其特征在于： 所述步骤 S5213, 按照以下步骤计算强迫分量：
S52131, 根据采样序列 Υ计算 y 的 n+1 阶导数， η=0,1,...， 直到满足
≤ε , 从而确定 n, 其中 ε为接近于 0的常数;
Figure imgf000020_0001
552132, 做变量替换; c = ^Z, 由采样序列 Y计算出数据序列 { }；
dtn
552133,对于数据序列 { },判断当前数值 处于暂态过程还是稳态过程； S52134, 当判定当前数值 处于稳态过程， 根据 { }计算平均值 作为稳 态值 , 进行步骤 S52135; 当判定当前数值 处于暂态过程， 令 k=l;
S52135, 对 X的稳态值 做积分反变换， 计算出 y的强迫分量 。
15. 如权利要求 14所述的参数辨识方法，其特征在于： 所述步骤 S52135, 按照公式 =fl。+i + . 2 + i 3+. + fl„. "计算强迫分量 , 其中， an =丄.¾ , α0,α1,α2,···,απ_1是按照导数公式和物理量 y及其各阶导数的初值确定的常数。
16. 一种包括如权利要求 8所述的数据采样系统的参数辨识系统， 其特征 在于： 还包括参数辨识单元， 用于根据采样序列 Y进行动态参数辨识和 /或静 态参数辨识。
17. 如权利要求 16所述的参数辨识系统， 其特征在于： 所述参数辨识单 元为动态参数辨识单元， 用于才艮据采样序列 Y, 采用自适应控制中的辨识方法 辨识 ARM AX模型的阶数和参数、 或参考模型法辨识参数、 或卡尔曼滤波器 辨识参数。
18. 如权利要求 16所述的参数辨识系统， 其特征在于： 所述参数辨识单 元为静态参数辨识单元， 包括：
稳态值计算单元， 用于根据采样序列 Y计算稳态值；
静态参数估计单元， 用于根据稳态值， 采用参数估计法、 相关系数法、 线 性回归法、 可线性化的线性回归法、 或最小二乘法估计静态参数。
19. 如权利要求 18所述的参数辨识系统， 其特征在于： 所述稳态值计算 单元包括：
判断单元， 根据采样序列 Y, 依据 t分布或其筒化、 或滤波器输出结果判 断当前采样值处于稳态过程或暂态过程； 所述滤波器为卡尔曼滤波器或 αβγ 滤波器；
计算稳态单元 Α, 用于当判定当前采样值处于稳态过程时，根据采样序列 Υ计算平均值作为稳态值；
计算稳态单元 Β, 用于当判定当前采样值处于暂态过程时， 根据采样序列 Υ计算强迫分量作为稳态值。
20. 如权利要求 19所述的参数辨识系统， 其特征在于： 计算稳态单元 Β 包括：
求导单元， 用于根据采样序列 Υ计算 y的 n+1阶导数， η=0,1,...， 直到满 1
足 d y
n+l ≤ε , 从而确定 n, 其中 ε为接近于 0的常数;
dt
替换单元， 用于做变量替换 c = 2, 由采样序列 Y计算出数据序列 { }； dtn
判断单元， 用于接收数据序列" [ }, 判断当前数值 处于暂态过程还是稳 态过程；
平均单元， 用于当判定当前数值 处于稳态过程， 根据 { }计算平均值 作为稳态值 , 进入还原单元； 当判定当前数值 处于暂态过程， 令 k=l, 到 出口；
还原单元，按照公式 =0。+0 + ^ + 03^+〜+ 0„^"计算强迫分量 , 其 中， 《„=丄 ·¾, , 是按照导数公式和物理量 y及其各阶导数的初值 确定的常数。
PCT/CN2013/085813 2012-10-23 2013-10-23 数据采样方法与系统及其在参数辨识中的应用方法与系统 WO2014063635A1 (zh)

Priority Applications (4)

Application Number Priority Date Filing Date Title
JP2015537134A JP6069809B2 (ja) 2012-10-23 2013-10-23 データサンプリング方法とシステム、及びパラメータ識別におけるその応用方法とシステム
US14/437,972 US9438449B2 (en) 2012-10-23 2013-10-23 Data sampling method and system, and application method and system thereof in parameter identification
CA2889334A CA2889334C (en) 2012-10-23 2013-10-23 Data sampling method and system, and application method and system thereof in parameter identification
EP13849426.5A EP2913930A4 (en) 2012-10-23 2013-10-23 Data sampling method and system, and application method and system thereof in parameter identification

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201210408534.9A CN102946253B (zh) 2012-10-23 2012-10-23 数据采样方法与系统及其在参数辨识中的应用方法与系统
CN201210408534.9 2012-10-23

Publications (1)

Publication Number Publication Date
WO2014063635A1 true WO2014063635A1 (zh) 2014-05-01

Family

ID=47729167

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2013/085813 WO2014063635A1 (zh) 2012-10-23 2013-10-23 数据采样方法与系统及其在参数辨识中的应用方法与系统

Country Status (6)

Country Link
US (1) US9438449B2 (zh)
EP (1) EP2913930A4 (zh)
JP (1) JP6069809B2 (zh)
CN (1) CN102946253B (zh)
CA (1) CA2889334C (zh)
WO (1) WO2014063635A1 (zh)

Families Citing this family (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102946253B (zh) * 2012-10-23 2016-06-08 保定市三川电气有限责任公司 数据采样方法与系统及其在参数辨识中的应用方法与系统
CN102928014B (zh) * 2012-10-23 2015-05-13 保定市三川电气有限责任公司 电力系统数字测量或遥测处理的方法及装置
CN103441495B (zh) * 2013-08-28 2015-07-29 三川电力设备股份有限公司 电力系统元件参数和功率修正系数的辨识方法及系统
CN103439607B (zh) * 2013-08-28 2016-11-23 三川电力设备股份有限公司 由故障录波辨识元件全参数的方法和系统及故障定位方法
CN105787219B (zh) * 2016-04-21 2018-08-24 北京航空航天大学 一种利用临近频点采样建立传导干扰耦合通道多元线性回归模型的方法
CN106992661A (zh) * 2017-04-17 2017-07-28 广西大学 用于pwm数字控制的过采样即时信号处理方法
CN107748497B (zh) * 2017-09-28 2019-10-15 中国科学院长春光学精密机械与物理研究所 运动控制系统的模型辨识和参数设计方法及系统
CN108036864A (zh) * 2017-11-06 2018-05-15 武汉航空仪表有限责任公司 一种尾桨温度传感器信号的fir滤波方法
CN108197381B (zh) * 2017-12-29 2019-09-27 河海大学 基于寻优空间形态分析的参数辨识方法
CN111293686A (zh) * 2020-02-29 2020-06-16 上海电力大学 基于armax系统辨识的电力系统惯量实时评估方法
CN114167113B (zh) * 2021-12-31 2024-05-17 上海市计量测试技术研究院 一种精确确定积分式数字多用表带宽的采样方法及系统

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101614554A (zh) * 2009-07-08 2009-12-30 保定市三川电气有限责任公司 连续物理量测量装置及方法
CN102331535A (zh) * 2011-06-09 2012-01-25 郝玉山 交流电物理量测量和数据采集装置和方法
CN102377436A (zh) * 2010-08-16 2012-03-14 Nxp股份有限公司 低功率高动态范围西格玛-德尔塔调制器
CN102393214A (zh) * 2011-06-09 2012-03-28 郝玉山 连续物理量数据采集方法和装置
CN102946253A (zh) * 2012-10-23 2013-02-27 保定市三川电气有限责任公司 数据采样方法与系统及其在参数辨识中的应用方法与系统

Family Cites Families (21)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5015934A (en) * 1989-09-25 1991-05-14 Honeywell Inc. Apparatus and method for minimizing limit cycle using complementary filtering techniques
US5224170A (en) * 1991-04-15 1993-06-29 Hewlett-Packard Company Time domain compensation for transducer mismatch
US5297166A (en) * 1992-07-02 1994-03-22 National Semiconductor Corporation Method and apparatus for decision feedback equalization with reduced convergence time
JP3134602B2 (ja) * 1993-06-17 2001-02-13 ヤマハ株式会社 波形記録装置
JP3223491B2 (ja) * 1993-07-12 2001-10-29 ソニー株式会社 再生装置
US6005377A (en) * 1997-09-17 1999-12-21 Lucent Technologies Inc. Programmable digital controller for switch mode power conversion and power supply employing the same
JP3942790B2 (ja) * 2000-02-24 2007-07-11 アンリツ株式会社 信号分析装置
US6751043B2 (en) * 2000-12-15 2004-06-15 Texas Instruments Incorporated Digital servo control system for a hard disc drive using a voice coil motor in voltage mode
US6820108B2 (en) * 2001-09-07 2004-11-16 Motorola, Inc. Method and apparatus to perform division in hardware
ATE445709T1 (de) * 2002-10-16 2009-10-15 Scripps Research Inst Glycoproteinsynthese
EP1447952B1 (en) * 2002-12-09 2011-06-22 Rohde & Schwarz GmbH & Co. KG Method and device for analysing an OFDM signal
EP1580882B1 (en) * 2004-03-19 2007-01-10 Harman Becker Automotive Systems GmbH Audio enhancement system and method
JP4534609B2 (ja) * 2004-06-07 2010-09-01 ヤマハ株式会社 ディジタルフィルタ設計システム及び方法
JP2006217120A (ja) * 2005-02-02 2006-08-17 Yokogawa Electric Corp データ収集装置
EP1884031A4 (en) * 2005-05-25 2014-01-22 Broadcom Corp PLL WITH PHASE CONTROL AND RESYNCHRONIZATION
US7436911B2 (en) * 2005-10-11 2008-10-14 L-3 Communications Integrated Systems L.P. Nyquist folded bandpass sampling receivers with narrow band filters for UWB pulses and related methods
CN100439926C (zh) * 2006-03-31 2008-12-03 北京万工科技有限公司 用于洛高夫斯基(Rogowski)线圈的积分器及其实现方法
JP4258545B2 (ja) * 2006-11-22 2009-04-30 トヨタ自動車株式会社 デジタルローパスフィルタ
JP5092580B2 (ja) * 2007-06-26 2012-12-05 ソニー株式会社 デジタル信号処理装置、デジタル信号処理方法及びデジタル信号処理プログラム
CN100504309C (zh) * 2007-09-30 2009-06-24 南京大学 基于快速傅立叶变换的布里渊光时域反射测量方法
TWI434494B (zh) * 2009-10-12 2014-04-11 Richtek Technology Corp 具適應性電壓位置控制之電源轉換器控制電路及其控制方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101614554A (zh) * 2009-07-08 2009-12-30 保定市三川电气有限责任公司 连续物理量测量装置及方法
CN102377436A (zh) * 2010-08-16 2012-03-14 Nxp股份有限公司 低功率高动态范围西格玛-德尔塔调制器
CN102331535A (zh) * 2011-06-09 2012-01-25 郝玉山 交流电物理量测量和数据采集装置和方法
CN102393214A (zh) * 2011-06-09 2012-03-28 郝玉山 连续物理量数据采集方法和装置
CN102946253A (zh) * 2012-10-23 2013-02-27 保定市三川电气有限责任公司 数据采样方法与系统及其在参数辨识中的应用方法与系统

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
See also references of EP2913930A4 *

Also Published As

Publication number Publication date
US9438449B2 (en) 2016-09-06
CN102946253B (zh) 2016-06-08
EP2913930A1 (en) 2015-09-02
JP6069809B2 (ja) 2017-02-01
CA2889334C (en) 2017-07-18
CN102946253A (zh) 2013-02-27
CA2889334A1 (en) 2014-05-01
JP2016500958A (ja) 2016-01-14
EP2913930A4 (en) 2017-02-08
US20150280943A1 (en) 2015-10-01

Similar Documents

Publication Publication Date Title
WO2014063635A1 (zh) 数据采样方法与系统及其在参数辨识中的应用方法与系统
Wang et al. Multi-stage holomorphic embedding method for calculating the power-voltage curve
US20170030958A1 (en) Transformer parameter estimation using terminal measurements
CN110968833A (zh) 一种用于模拟量校准的校准函数关系获取方法及装置
CN105115573B (zh) 一种洪水流量预报的校正方法和装置
WO2019165743A1 (zh) 确定对风角度偏差及修正对风角度的方法、装置和系统
CN113364408B (zh) 一种光伏跟踪支架阵列运行精度的测量方法和装置
US20200110115A1 (en) Systems and methods for operating generators based on generator steady state stability limits
CN105552885B (zh) 一种提高配电网状态估计可观测性的方法及其系统
CN114460526B (zh) 基于随动补偿的变电站电流互感器误差预测方法及系统
WO2014063634A1 (zh) 电力系统数字测量或遥测处理的方法及装置
CN110470941B (zh) 交流电的漏电流检测方法、装置、设备及存储介质
CN116995655A (zh) 一种光伏短期预测曲线调整方法及系统
CN114069667B (zh) 一种储能群功率分配方法、系统、处理设备及存储介质
RU2683999C1 (ru) Способ цифровой коррекции эффекта насыщения магнитопровода трансформатора тока
JP5739309B2 (ja) ディジタル保護制御装置
WO2018050741A1 (en) Improvements in or relating to the measurement of current within a conductor
JP2018152935A (ja) 電力系統状態推定装置
CN110648021B (zh) 一种两级电力负荷预测结果协调方法、装置及设备
CN117169800B (zh) 一种基于知识引导的电流互感器在线监测方法及装置
CN111141254B (zh) 建筑物变形监测边点沉降缺失信息的获取方法和装置
CN117764271A (zh) 基于扩展卡尔曼粒子滤波的电网动态状态估计方法与系统
TWI550416B (zh) 利用高維度子空間產生取樣之方法
Liu Check for Application of Hidden Markov Model on the Prediction of Hepatitis B Incidences Qiong Liu and Jianhua Yang School of Economics and Management, University of Science and Technology
CN117634443A (zh) 关口表数据矫正方法、装置、设备和存储介质

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 13849426

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2015537134

Country of ref document: JP

Kind code of ref document: A

ENP Entry into the national phase

Ref document number: 2889334

Country of ref document: CA

NENP Non-entry into the national phase

Ref country code: DE

WWE Wipo information: entry into national phase

Ref document number: 14437972

Country of ref document: US

WWE Wipo information: entry into national phase

Ref document number: 2013849426

Country of ref document: EP