CN105954695B - Synchronization-based homogeneous sensor mutation parameter identification method and device - Google Patents

Synchronization-based homogeneous sensor mutation parameter identification method and device Download PDF

Info

Publication number
CN105954695B
CN105954695B CN201610248401.8A CN201610248401A CN105954695B CN 105954695 B CN105954695 B CN 105954695B CN 201610248401 A CN201610248401 A CN 201610248401A CN 105954695 B CN105954695 B CN 105954695B
Authority
CN
China
Prior art keywords
signal
mutation
data
sensor
homogeneous sensor
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.)
Active
Application number
CN201610248401.8A
Other languages
Chinese (zh)
Other versions
CN105954695A (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.)
Yinchuan Power Supply Co Of State Grid Ningxia Electric Power Co
State Grid Corp of China SGCC
Xuji Group Co Ltd
Henan Xuji Instrument Co Ltd
Original Assignee
Yinchuan Power Supply Co Of State Grid Ningxia Electric Power Co
State Grid Corp of China SGCC
Xuji Group Co Ltd
Henan Xuji Instrument Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Yinchuan Power Supply Co Of State Grid Ningxia Electric Power Co, State Grid Corp of China SGCC, Xuji Group Co Ltd, Henan Xuji Instrument Co Ltd filed Critical Yinchuan Power Supply Co Of State Grid Ningxia Electric Power Co
Priority to CN201610248401.8A priority Critical patent/CN105954695B/en
Publication of CN105954695A publication Critical patent/CN105954695A/en
Application granted granted Critical
Publication of CN105954695B publication Critical patent/CN105954695B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R35/00—Testing or calibrating of apparatus covered by the other groups of this subclass

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Testing Or Calibration Of Command Recording Devices (AREA)

Abstract

The invention relates to a synchronous homogeneous sensor mutation parameter identification method and a synchronous homogeneous sensor mutation parameter identification device. And preprocessing the obtained data samples based on maximum likelihood estimation. Performing empirical mode decomposition on the signals of the sensor containing the sudden change based on the combined preprocessing to obtain IMF modal components of each order; and performing Hilbert-yellow transformation on the signal leading IMF component to obtain an instantaneous frequency graph and an instantaneous amplitude graph, and performing parameter detection and identification on abnormal signals of the sensor according to the catastrophe points, amplitude change and instantaneous frequency trend characteristic information in the HHT graph. The invention carries out undifferentiated fusion on different types of sensors to obtain equal-period data samples, reduces the difference of different sensors and random errors in the signal acquisition process, and improves the identification precision of the mutation parameters of the sensors.

Description

Synchronization-based homogeneous sensor mutation parameter identification method and device
Technical Field
The invention provides an effective method for identifying the detection and identification of sudden change parameters of a sensor after researching a fusion process of electrical data measured by the sensor and providing a method for engineering calculation and implementation by taking a homogeneous sensor signal commonly used by power equipment as an object in research.
Background
On one hand: the traditional power supply and utilization information collection is usually completed by a single sensor, even if a plurality of sensors are adopted, the sensors are used in a time sharing mode, and therefore the information of the power grid is reflected from a plurality of side surfaces in an isolated mode. As technology advances, these measurement data require a fusion process, i.e., the output of multiple sensors is used to infer a valid information.
On the other hand: the sensors in the power system are subject to different production environments and different technologies of manufacturers, and the physical data of the same object measured by instruments of different manufacturers or even instruments of different batches of the same manufacturer are different, particularly, the difference in amplitude is obvious, so that difficulties are caused in mass data comparison, automatic analysis and the like.
Various sensors in the power system are homogeneous sensors for the reasons described above. A homogeneous sensor refers to several sensors observing the same physical phenomenon, which may be different in time and location, but the characteristics of the detected or collected signals are the same. Due to individual differences of different sensors and random errors in the signal acquisition process, the mass data of the homogeneous sensor are difficult to detect and identify, the data quality is poor, and the identification precision of the mutation parameters of the sensor is low.
Disclosure of Invention
The invention aims to provide a synchronous homogeneous sensor mutation parameter identification method and a synchronous homogeneous sensor mutation parameter identification device, which are used for reducing the difference of different sensors and random errors in the signal acquisition process, improving the data quality and improving the identification precision of the sensor mutation parameters.
In order to achieve the above object, the scheme of the invention comprises:
a synchronous homogeneous sensor mutation parameter identification method comprises the following steps:
step A: the central processing unit takes T as a sampling period, the homogeneous sensor samples and quantizes a measured signal of the system at regular time, and sampling data under the same sampling frequency is obtained;
and B: intercepting the sampling data according to cycles to obtain m cycle sequences as follows: y is1(t),Y2(t),…Ym(t) N data points, Y, are contained in each periodic samplei＝[Xi(1),Xi(2),…Xi(N)]Wherein i is 1,2, …, m;
and C: constructing a periodic data sequence constructed according to periodic synchronization: each time point data in each period data is regarded as a multiple measurement result of the same object after synchronization, namely the time point data corresponding to each period is a sample, each time point value corresponding to the corresponding time in the multiple periods data forms a data sequence, and therefore the data sequence can be regarded as a one-dimensional data sequence measured by a single sensor multiple times, namely a period data sequence Y constructed according to period synchronization is constructed1′(t),Y2′(t),…Y′N(t), each set of one-dimensional data sequences includes m data points, i.e. Yi′＝[X1(1),X2(2),…Xm(N)]；
Step D: for data sample Yt' performing optimal combination preprocessing based on maximum likelihood estimation;
step E: performing EMD on the mutation-containing sensor signals based on combined preprocessing to obtain IMF components of each order;
step F: adopting HHT (Hilbert-Huang transform) to transform the dominant IMF component of the abnormal signal of the sensor to obtain an instantaneous frequency graph and an instantaneous amplitude graph;
step G: and detecting and identifying the sudden change parameters of the sensor according to the sudden change points, amplitude change and instantaneous frequency trend characteristic information in the HHT transformation diagram.
Further, in the step D, the data samples Y are processedtThe mode of performing the optimal combination preprocessing based on the maximum likelihood estimation is as follows:
s101: assuming that the environmental characteristic to be measured is X and the value of the sensor is Y at a given time, the measurement model of the sensor is: y ═ f (x) + V where V is a gaussian-distributed noise term; the data fusion is to obtain the measured value Y by N sensors1、Y2、…、YNObtaining the optimal estimation of the characteristic parameter X from the measured values according to a certain estimation criterion;
s102: finding a suitable criterion function that produces minimal loss when X is estimated to be X (y); taking the loss function as the uniform loss:
Figure BDA0000970259010000031
s103: on the basis of the loss function L, a function R of the corresponding estimated risk is defined:
Figure BDA0000970259010000032
wherein p (x), p (x | y) represent probability distributions;
s104: taking the minimum risk as an estimation criterion, i.e.
Figure BDA0000970259010000033
The optimum estimate in conformity with the formula (1) can be obtained as
Figure BDA0000970259010000034
S105: in a system with N sensors, the corresponding information fusion can be seen as the observation Y1、Y2、…、YNNext, the value X has an estimate of the maximum a posteriori of
Figure BDA0000970259010000035
S106: to find the maximum a posteriori estimate, p (X) is 1 for all possible parameters X, derived from the formula:
Figure BDA0000970259010000036
at this time, the maximum posterior estimation is simplified into maximum likelihood estimation, and the corresponding fusion calculation formula is as follows:
Figure BDA0000970259010000041
Figure BDA0000970259010000042
when the sensor is one-dimensional, and the coordinate transformation is not considered, equations (3) and (4) can be simplified as follows:
Figure BDA0000970259010000043
Figure BDA0000970259010000044
further, the step E of obtaining each order of IMF component by performing EMD decomposition on the mutation-containing sensor signal based on the combined preprocessing includes:
s201: interpolating all local maximum points and all local minimum points of the signal x [ t ] by a cubic spline function, and fitting an upper envelope line and a lower envelope line; x < t > is the sensor signal containing the mutation characteristic after combined pretreatment;
s202: calculating the average curve M of the upper and lower envelope lines1(t) (t is 1,2, 3 … m), the signal x [ t ] is sampled]And M1The difference of (t) is P1(t)：P1(t)＝x[t]-M1(t)；
S203: if P is1(t) satisfy at the same timeThe following two conditions of IMF are then the first IMF component, otherwise, it is repeated as the original signal S201 to S202, resulting in P11(t)：P11(t)＝P1(t)-M11(t), wherein: m11(t) is P1(t) an average curve of the upper and lower envelope lines;
the IMF component satisfies two conditions: in the whole time process, the zero crossing times are equal to the number of extreme points or have at most 1 difference; at any point on the signal, the average value of an upper envelope line defined by a local maximum value and a lower envelope line defined by a local minimum value point is 0, namely the signal is locally symmetrical about a time axis;
s204: repeating the above steps until P obtained by formula (1) at the k-th screening1k(t) two conditions for the IMF component are satisfied: p1k(t)＝P1(1-k)(t)-M1k(t) (1)
S205: in addition to the signal, the threshold value S can be obtained by the equation (2) in the actual calculationDJudging whether each screening result is an IMF component:
Figure BDA0000970259010000051
wherein: r is the number of sampling points of the sensor signal, threshold value SDUsually 0.2 to 0.3;
s206: let C1(t)＝P1k(t) then C1(t) is the first IMF component, which contains the original signal x [ t ]]The least intermediate-cycle IMF component; c is to be1(t) from x [ t ]]Separating out: r1(t)＝x[t]-C1(t)；
S207: r is to be1(t) repeating the above steps S201 to S205n times as a new value, a signal x [ t ] can be obtained]N IMF components: rn(t)＝Rn-1(t)-Cn(t)；
S208: when R isn(t) is a monotonic function from the signal x [ t ]]When other components can not be solved, the whole decomposition process is finished; obtaining:
Figure BDA0000970259010000052
further, the step F of obtaining the instantaneous frequency map and the instantaneous amplitude map by applying HHT transformation to the dominant IMF component of the sensor abnormal signal includes:
s301: performing Hilbert transform on all IMF components obtained after EMD decomposition; the Hilbert form of given C (t) is:
Figure BDA0000970259010000053
where λ is the integral variable.
S302: constructing an analysis signal Z (t):
Z(t)＝C(t)+iH(t)＝A(t)eiθ(t)；
s303: amplitude function:
Figure BDA0000970259010000054
s304: argument function:
Figure BDA0000970259010000055
s305: instantaneous frequency:
Figure BDA0000970259010000056
further, the process of detecting and identifying according to the mutation point, amplitude variation and instantaneous frequency trend characteristic information sensor mutation parameter in the HHT transformation graph in the step G is as follows:
s306: classifying the mutation parameters in the sensor: dip, swell, harmonic, pulse;
s307: according to different characteristic information such as mutation points, amplitude changes, instantaneous frequency trends and the like of different sensor abnormal signals in the HHT graph, the sensor is classified and identified when the sensor undergoes any mutation.
The invention also provides a synchronous homogeneous sensor mutation parameter identification device, which comprises:
a module A: the central processing unit takes T as a sampling period, the homogeneous sensor samples and quantizes a measured signal of the system at regular time, and sampling data under the same sampling frequency is obtained;
and a module B: intercepting the sampling data according to cycles to obtain m cycle sequences as follows: y is1(t),Y2(t),…Ym(t) N data points, Y, are contained in each periodic samplei＝[Xi(1),Xi(2),…Xi(N)]Wherein i is 1,2, …, m;
and a module C: constructing a periodic data sequence constructed according to periodic synchronization: each time point data in each period data is regarded as a multiple measurement result of the same object after synchronization, namely the time point data corresponding to each period is a sample, each time point value corresponding to the corresponding time in the multiple periods data forms a data sequence, and therefore the data sequence can be regarded as a one-dimensional data sequence measured by a single sensor multiple times, namely a period data sequence Y constructed according to period synchronization is constructed1′(t),Y2′(t),…Y′N(t), each set of one-dimensional data sequences includes m data points, i.e. Yi′＝[X1(1),X2(2),…Xm(N)]；
A module D: for data sample Yt' performing optimal combination preprocessing based on maximum likelihood estimation;
and a module E: performing EMD on the mutation-containing sensor signals based on combined preprocessing to obtain IMF components of each order;
and a module F: adopting HHT (Hilbert-Huang transform) to transform the dominant IMF component of the abnormal signal of the sensor to obtain an instantaneous frequency graph and an instantaneous amplitude graph;
a module G: and detecting and identifying the sudden change parameters of the sensor according to the sudden change points, amplitude change and instantaneous frequency trend characteristic information in the HHT transformation diagram.
Further, in the module D, the data samples Y are processedtThe mode of performing the optimal combination preprocessing based on the maximum likelihood estimation is as follows:
s101: assuming that the environmental characteristic to be measured is X and the value of the sensor is Y at a given time, the measurement model of the sensor is: y ═ f (x) + V where V is a gaussian-distributed noise term; the data fusion is formed by N sensorsThe device obtains a measured value Y1、Y2、…、YNObtaining the optimal estimation of the characteristic parameter X from the measured values according to a certain estimation criterion;
s102: finding a suitable criterion function that produces minimal loss when X is estimated to be X (y); taking the loss function as the uniform loss:
Figure BDA0000970259010000071
s103: on the basis of the loss function L, a function R of the corresponding estimated risk is defined:
Figure BDA0000970259010000072
wherein p (x), p (x | y) represent probability distributions;
s104: taking the minimum risk as an estimation criterion, i.e.
Figure BDA0000970259010000073
The optimum estimate in conformity with the formula (1) can be obtained as
Figure BDA0000970259010000074
S105: in a system with N sensors, the corresponding information fusion can be seen as the observation Y1、Y2、…、YNNext, the value X has an estimate of the maximum a posteriori of
Figure BDA0000970259010000075
S106: to find the maximum a posteriori estimate, p (X) is 1 for all possible parameters X, which can be derived from the formula:
Figure BDA0000970259010000081
at this time, the maximum posterior estimation is simplified into maximum likelihood estimation, and the corresponding fusion calculation formula is as follows:
Figure BDA0000970259010000082
Figure BDA0000970259010000083
when the sensor is one-dimensional, and the coordinate transformation is not considered, equations (3) and (4) can be simplified as follows:
Figure BDA0000970259010000084
Figure BDA0000970259010000085
further, in the module E, the process of obtaining each order of IMF component by performing EMD decomposition on the mutation-containing sensor signal based on the combined preprocessing includes:
s201: interpolating all local maximum points and all local minimum points of the signal x [ t ] by a cubic spline function, and fitting an upper envelope line and a lower envelope line; x < t > is the sensor signal containing the mutation characteristic after combined pretreatment;
s202: calculating the average curve M of the upper and lower envelope lines1(t) (t is 1,2, 3 … m), the signal x [ t ] is sampled]And M1The difference of (t) is P1(t)：P1(t)＝x[t]-M1(t)；
S203: if P is1(t) if both of the following IMF conditions are satisfied, it is the first IMF component, otherwise, it is repeated as the original signal S201 to S202, resulting in P11(t)：P11(t)＝P1(t)-M11(t), wherein: m11(t) is P1(t) an average curve of the upper and lower envelope lines;
the IMF component satisfies two conditions: in the whole time process, the zero crossing times are equal to the number of extreme points or have at most 1 difference; at any point on the signal, the average value of an upper envelope line defined by the local maximum value and a lower envelope line defined by the local minimum value point is 0, namely the signal is locally symmetrical about a time axis;
s204: repeating the above steps until P obtained by formula (1) at the k-th screening1k(t) two conditions for the IMF component are satisfied: p1k(t)＝P1(1-k)(t)-M1k(t) (1)
S205: in addition to the signal, the threshold value S can be obtained by the equation (2) in the actual calculationDJudging whether each screening result is an IMF component:
Figure BDA0000970259010000091
wherein: r is the number of sampling points of the sensor signal, threshold value SDUsually 0.2 to 0.3;
s206: let C1(t)＝P1k(t) then C1(t) is the first IMF component, which contains the original signal x [ t ]]The least intermediate-cycle IMF component; c is to be1(t) from x [ t ]]Separating out: r1(t)＝x[t]-C1(t)；
S207: r is to be1(t) repeating the above steps S201 to S205n times as a new value, a signal x [ t ] can be obtained]N IMF components: rn(t)＝Rn-1(t)-Cn(t)；
S208: when R isn(t) is a monotonic function from the signal x [ t ]]When other components can not be solved, the whole decomposition process is finished; obtaining:
Figure BDA0000970259010000092
further, the process that the module F obtains the instantaneous frequency map and the instantaneous amplitude map by applying HHT transformation to the dominant IMF component of the sensor abnormal signal is as follows:
s301: subjecting all IMF components obtained by EMD decomposition to Hilbert transform
(Hilbert) transform. The Hilbert form of given C (t) is:
Figure BDA0000970259010000093
where λ is the integral variable.
S302: constructing an analysis signal Z (t):
Z(t)＝C(t)+iH(t)＝A(t)eiθ(t)；
s303: amplitude function:
Figure BDA0000970259010000094
s304: argument function:
Figure BDA0000970259010000101
s305: instantaneous frequency:
Figure BDA0000970259010000102
further, the process of detecting and identifying in the module G according to the abrupt change point, the amplitude change, and the instantaneous frequency trend characteristic information sensor abrupt change parameter in the HHT transformation diagram is as follows:
s306: classifying the mutation parameters in the sensor: dip, swell, harmonic, pulse;
s307: according to different characteristic information such as mutation points, amplitude changes, instantaneous frequency trends and the like of different sensor abnormal signals in the HHT graph, the sensor is classified and identified when the sensor undergoes any mutation.
The invention discloses a synchronous homogeneous sensor mutation parameter identification method and a synchronous homogeneous sensor mutation parameter identification device, wherein a central processing unit takes T seconds as a sampling period, a homogeneous sensor (the same physical phenomenon is observed by the sensor) samples and quantizes a detected signal at regular time and obtains a data sequence X (N) related to time, sampled data is intercepted according to the period, and m period sequences are obtained as follows: y is1(t),Y2(t),…Ym(t) N data points, Y, are contained in each periodic samplei＝[Xi(1),Xi(2),…Xi(N)]Where i is 1,2, …, m. One point value is extracted from each time point of each period sequence to form a one-dimensional data sequence measured by a single sensor for multiple times, and a period-based one-dimensional data sequence is constructedThe periodic data sequence being constructed synchronously, i.e. Yi′(t),Y2′(t),…Y′N(t), each set of one-dimensional data sequences includes m data points, i.e. Yi′＝[X1(1),X2(2),…Xm(N)]. For the obtained data sample Yt' performing optimization combination preprocessing based on maximum likelihood estimation and least square estimation, and guiding periodic sampling data of the power system. Performing empirical mode decomposition on the signals of the sensor containing the sudden change based on the combined preprocessing to obtain IMF modal components of each order; and performing Hilbert-yellow transformation on the signal leading IMF component to obtain an instantaneous frequency graph and an instantaneous amplitude graph, and performing parameter detection and identification on abnormal signals of the sensor according to the catastrophe points, amplitude change and instantaneous frequency trend characteristic information in the HHT graph. The invention carries out undifferentiated fusion on different types of sensors to obtain equal-period data samples, reduces the difference of different sensors and random errors in the signal acquisition process, and improves the identification precision of the mutation parameters of the sensors.
Drawings
FIG. 1 is a flow chart of a synchronization-based homogeneous sensor mutation parameter identification method;
FIG. 2 is a graph of HHT instantaneous frequency and amplitude of the IMF1 component of the analog sensor dip signal of the present invention;
FIG. 3 is a graph of HHT instantaneous frequency and amplitude of the IMF1 component of the analog sensor transient signal of the present invention;
FIG. 4 is a graph of HHT instantaneous frequency and amplitude of the IMF1 component of the analog sensor pulse signal of the present invention;
FIG. 5 is a graph of the HHT instantaneous frequency and amplitude of the IMF1 component of the analog sensor harmonic signal of the present invention.
Detailed Description
The present invention will be described in further detail with reference to the accompanying drawings.
Method embodiment
As shown in fig. 1, a homogeneous sensor mutation parameter identification method based on synchronization includes the following steps:
step A: the central processing unit takes T (second) as a sampling period, the homogeneous sensor samples and quantifies a measured signal of the system at regular time, and the specific steps of obtaining data samples X (N) under the same sampling frequency are as follows:
the sampling period of the central processing unit is T (second), m sensors (which can be different and different in position, but have the same characteristics of detected or collected signals) contained in the system sample and quantize the detected signals in the system at regular time, and then data samples X under the same sampling frequency are obtainedi(N), wherein i ═ 1,2, …, m.
Further, the step C: constructing periodic data sequences constructed synchronously by period, i.e. Y1′(t),Y2′(t),…Y′N(t), each set of one-dimensional data sequences includes m data points, i.e. Yi′＝[X1(1),X2(2),…Xm(N)]The method comprises the following specific steps:
the method comprises the steps of intercepting collected power system data according to periods, regarding each point-in-time data in each period data as a multiple measurement result of the same synchronized object, namely, the point-in-time data corresponding to each period is a sample, and each point-in-time data corresponding to the time in the multiple periods forms a data sequence, so that the data sequences can be regarded as one-dimensional data sequences of multiple measurements of a single sensor, namely, a period data sequence synchronously constructed according to the periods is constructed, namely Y1′(t),Y2′(t),…Y′N(t), each set of one-dimensional data sequences includes m data points, i.e. Yi′＝[X1(1),X2(2),…Xm(N)]。
Further, the step D: for data sample Yt' the specific steps for preprocessing based on maximum likelihood estimation are:
s101: assuming that the environmental characteristic to be measured is X and the value of the sensor is Y at a given time, the measurement model of the sensor is: y ═ f (x) + V where V is a gaussian-distributed noise term. The data fusion is to obtain the measured value Y by N sensors1、Y2、…、YNAnd an optimal estimation of the characteristic parameter X is obtained from these measurements according to some estimation criterion.
S102: find the appropriate criterion function that produces the least loss when X is estimated to be X (y). Taking the loss function as the uniform loss:
Figure BDA0000970259010000121
s103: on the basis of the loss function L, a function R of the corresponding estimated risk is defined:
Figure BDA0000970259010000122
wherein p (x), p (x | y) represent probability distributions;
s104: taking the minimum risk as an estimation criterion, i.e.
Figure BDA0000970259010000123
An optimum estimate (maximum a posteriori estimate) can be obtained which conforms to equation (1) as
Figure BDA0000970259010000124
S105: in a system with N sensors, the corresponding information fusion can be seen as the observation Y1、Y2、…、YNNext, the value X has an estimate of the maximum a posteriori of
Figure BDA0000970259010000131
S106: in order to find the maximum a posteriori estimate, under certain conditions, we cannot determine the prior distribution of the characteristic parameter X, and we use the concept of "fuzzy prior", that is, p (X) -1 is used for all possible parameters X, and we can obtain by formula derivation:
Figure BDA0000970259010000132
at this time, the maximum posterior estimation is simplified into maximum likelihood estimation, and the corresponding fusion calculation formula is as follows:
Figure BDA0000970259010000133
Figure BDA0000970259010000134
when the sensor is one-dimensional, and the coordinate transformation is not considered, equations (3) and (4) can be simplified as follows:
Figure BDA0000970259010000135
Figure BDA0000970259010000136
further, the step E:
note that: and x [ t ] is the sensor signal containing the mutation characteristic after combined pretreatment.
S201: interpolating all local maximum points and all local minimum points of the signal x [ t ] by a cubic spline function, and fitting an upper envelope line and a lower envelope line;
s202: calculating the average curve M of the upper and lower envelope lines1(t) (t is 1,2, 3 … m), the signal x [ t ] is sampled]And M1The difference of (t) is P1(t)：P1(t)＝x[t]-M1(t)；
S203: if P is1(t) if both of the following IMF conditions are satisfied, it is the first IMF component, otherwise, it is repeated as the original signal S201 to S202, resulting in P11(t)：P11(t)＝P1(t)-M11(t), wherein: m11(t) is P1(t) an average curve of the upper and lower envelope lines;
the IMF component satisfies two conditions: in the whole time process, the zero crossing times are equal to the number of extreme points or have at most 1 difference; at any point on the signal, the average of the upper envelope defined by the local maxima and the lower envelope defined by the local minima is 0, i.e. the signal is locally symmetric about the time axis.
S204: repeating the above steps(ii) step selection until P obtained by the formula (1) at the k-th selection1k(t) two conditions for the IMF component are satisfied: p1k(t)＝P1(1-k)(t)-M1k(t) (1)
S205: in addition to the signal, the threshold value S can be obtained by the equation (2) in the actual calculationDJudging whether each screening result is an IMF component:
Figure BDA0000970259010000141
wherein: r is the number of sampling points of the power system signal, and the threshold value SDUsually 0.2 to 0.3;
s206: let C1(t)＝P1k(t) then C1(t) is the first IMF component, which contains the original signal x [ t ]]The IMF component with the shortest mid-period. C is to be1(t) from x [ t ]]Separating out: r1(t)＝x[t]-C1(t)；
S207: r is to be1(t) repeating the above steps S201 to S205n times as a new value, a signal x [ t ] can be obtained]N IMF components: rn(t)＝Rn-1(t)-Cn(t)；
S208: when R isn(t) is a monotonic function from the signal x [ t ]]When the other components can not be solved again, the whole decomposition process is finished. Obtaining:
Figure BDA0000970259010000142
further, the step F: the specific steps of obtaining an instantaneous frequency map and an instantaneous amplitude map by adopting HHT transformation on the dominant IMF component of the sensor abnormal signal are carried out according to two schemes:
s301: performing Hilbert transform on all IMF components obtained after EMD decomposition; the Hilbert form of given C (t) is:
Figure BDA0000970259010000143
where λ is the integral variable.
S302: constructing an analysis signal Z (t):
Z(t)＝C(t)+iH(t)＝A(t)eiθ(t)；
s303: amplitude function:
Figure BDA0000970259010000151
s304: argument function:
Figure BDA0000970259010000152
s305: instantaneous frequency:
Figure BDA0000970259010000153
further, the step G: the process of detecting and identifying the sudden change parameters of the sensor according to the sudden change points, the amplitude change and the instantaneous frequency trend characteristic information in the HHT transformation diagram comprises the following steps:
s306: classifying the mutation parameters in the sensor: dip, swell, harmonic, pulse.
S307: according to different characteristic information such as mutation points, amplitude changes, instantaneous frequency trends and the like of different sensor abnormal signals in the HHT graph, the sensor is classified and identified when the sensor undergoes any mutation.
The invention is introduced according to simulation tests: MATLAB was used to generate various sensor burst signals (dip, swell, harmonic, pulse) at a sampling frequency of 1kHz, 1000 points at a sampling frequency, and 50Hz at the fundamental frequency of the voltage, which are plotted in Table 1.
TABLE 1 Voltage disturbance types
Figure BDA0000970259010000154
Figure BDA0000970259010000161
According to table 1, different sensor mutation types are shown, and a simulation experiment is performed according to the present invention, and fig. 2 to 5 show the detection results of the algorithm of the present embodiment.
From the results fig. 2 and fig. 3 it can be seen that: the amplitude and the instantaneous frequency jump at 0.45s and 0.55s, the instantaneous frequency is slowly increased and then decreased, then is maintained for a period of stable time, then is slowly increased and then decreased, and the instantaneous frequency is maintained at about 50Hz at normal time; the amplitude is greatly increased and then is subjected to a stable process and then is greatly reduced, thereby indicating that the abnormal signal is suddenly reduced. And the sudden rising abnormal signal characteristic information is just opposite to the sudden falling abnormal signal. The abnormal signals of the sensors for sudden drop and sudden rise can be well distinguished by utilizing the different information.
From the results fig. 4 to 5, it can be seen that: the instantaneous frequency and the amplitude show spike information at a certain moment, and correspond to the pulse abnormal signal; instantaneous frequency over a certain period of time; the instantaneous frequency experiences violent jitter in a certain time period, and the amplitude of the instantaneous frequency firstly undergoes temporary rise in the same time period and then returns to a steady state after the amplitude is violently changed, and corresponds to a harmonic abnormal signal.
In summary, the following steps: aiming at different abnormal signals of the sensor, the algorithm can detect that different instantaneous frequency and amplitude characteristic information correspond to different abnormal signals; the abnormal signals of the sensor can be accurately classified, and the automatic analysis of the sudden change characteristics of the sensor is realized. The algorithm has good system identification and is convenient for hardware implementation.
The invention discloses a synchronous homogeneous sensor mutation parameter identification method, wherein a central processing unit takes T seconds as a sampling period, a homogeneous sensor (the same physical phenomenon is observed by the sensor) samples and quantizes a detected signal at regular time and obtains a data sequence X (N) related to time, and the sampled data is intercepted according to the period to obtain m period sequences as follows: y is1(t),Y2(t),…Ym(t) N data points, Y, are contained in each periodic samplei＝[Xi(1),Xi(2),…Xi(N)]Where i is 1,2, …, m. One point value is extracted from each time point of each period sequence to form a one-dimensional data sequence measured by a single sensor for multiple times, and a period data sequence constructed synchronously according to periods, namely Y1′(t),Y2′(t),…Y′N(t), each set of one-dimensional data sequences includes m data points, i.e. Yi′＝[X1(1),X2(2),…Xm(N)]. For the obtained data sample Yt' performing optimization combination preprocessing based on maximum likelihood estimation and least square estimation, and guiding periodic sampling data of the power system. Performing empirical mode decomposition on the signals of the sensor containing the sudden change based on the combined preprocessing to obtain IMF modal components of each order; and performing Hilbert-yellow transformation on the signal leading IMF component to obtain an instantaneous frequency graph and an instantaneous amplitude graph, and performing parameter detection and identification on abnormal signals of the sensor according to the catastrophe points, amplitude change and instantaneous frequency trend characteristic information in the HHT graph. The invention carries out undifferentiated fusion on different types of sensors to obtain equal-period data samples, reduces the difference of different sensors and random errors in the signal acquisition process, and improves the identification precision of the mutation parameters of the sensors.
Device embodiment
A synchronization-based homogeneous sensor mutation parameter identification device comprises:
a module A: the central processing unit takes T as a sampling period, the homogeneous sensor samples and quantizes a measured signal of the system at regular time, and sampling data under the same sampling frequency is obtained;
and a module B: intercepting the sampling data according to cycles to obtain m cycle sequences as follows: y is1(t),Y2(t),…Ym(t) N data points, Y, are contained in each periodic samplei＝[Xi(1),Xi(2),…Xi(N)]Wherein i is 1,2, …, m;
and a module C: constructing a periodic data sequence constructed according to periodic synchronization: each time point data in each period data is regarded as a multiple measurement result of the same object after synchronization, namely the time point data corresponding to each period is a sample, each time point value corresponding to the corresponding time in the multiple periods data forms a data sequence, and therefore the data sequence can be regarded as a one-dimensional data sequence measured by a single sensor multiple times, namely a period data sequence Y constructed according to period synchronization is constructed1′(t),Y2′(t),…Y′N(t), each set of one-dimensional data sequences includes m data points, i.e. Yi′＝[X1(1),X2(2),…Xm(N)]；
A module D: for data sample Yt' performing maximum likelihood estimation-based preprocessing;
and a module E: performing EMD on the mutation-containing sensor signals based on combined preprocessing to obtain IMF components of each order;
and a module F: adopting HHT (Hilbert-Huang transform) to transform the dominant IMF component of the abnormal signal of the sensor to obtain an instantaneous frequency graph and an instantaneous amplitude graph;
a module G: and detecting and identifying the sudden change parameters of the sensor according to the sudden change points, amplitude change and instantaneous frequency trend characteristic information in the HHT transformation diagram.
The apparatus in this embodiment is actually a software architecture for implementing the above method embodiment, and various modules in the apparatus are software functional modules, stored in the memory, and executed by the processor.
The present invention has been described in relation to particular embodiments thereof, but the invention is not limited to the described embodiments. In the thought given by the present invention, the technical means in the above embodiments are changed, replaced, modified in a manner that is easily imaginable to those skilled in the art, and the functions are basically the same as the corresponding technical means in the present invention, and the purpose of the invention is basically the same, so that the technical scheme formed by fine tuning the above embodiments still falls into the protection scope of the present invention.

Claims (6)

1. A synchronous homogeneous sensor mutation parameter identification method is characterized by comprising the following steps:
step A: the central processing unit takes T as a sampling period, the homogeneous sensor samples and quantizes a measured signal of the system at regular time, and sampling data under the same sampling frequency is obtained;
and B: intercepting the sampling data according to cycles to obtain m cycle sequences as follows: y is1(t),Y2(t),…Ym(t) N data points, Y, are contained in each periodic samplei＝[Xi(1),Xi(2),…Xi(N)]Wherein i is 1,2, …, m;
and C: constructing a periodic data sequence constructed according to periodic synchronization: the data of each time point in each period data is regarded as a multiple measurement result of the same object after synchronization, namely the data of the time point corresponding to each period is a sample, each point value of the corresponding time point in the data of multiple periods forms a data sequence, and the data is regarded as a one-dimensional data sequence measured by a single homogeneous sensor for multiple times, namely a period data sequence Y constructed according to period synchronization is constructed1′(t),Y2′(t),…YN' (t) each one-dimensional data sequence includes m data points, i.e. Yi′＝[X1(1),X2(2),…Xm(N)]；
Step D: for data sample Yt' performing maximum likelihood estimation-based preprocessing;
step E: performing EMD on the mutation-containing homogeneous sensor signal based on combined preprocessing to obtain IMF components of each order;
step F: adopting HHT (Hilbert-Huang transform) to obtain an instantaneous frequency graph and an instantaneous amplitude graph for the dominant IMF component of the abnormal signal of the homogeneous sensor;
step G: detecting and identifying the mutation parameters of the homogeneous sensor according to the mutation points, amplitude variation and instantaneous frequency trend characteristic information in the instantaneous frequency graph and the instantaneous amplitude graph;
in the step E, the EMD decomposition is adopted for the mutation-containing homogeneous sensor signals based on the combined pretreatment to obtain IMF components of each order:
s201: interpolating all local maximum points and all local minimum points of the signal x [ t ] by a cubic spline function, and fitting an upper envelope line and a lower envelope line; x < t > is the homogeneous sensor signal containing the mutation characteristic after combined pretreatment;
s202: calculating the average curve M of the upper and lower envelope lines1(t), if t is 1,2, 3 … m, the signal x [ t ] is sampled]And M1The difference of (t) is P1(t)：P1(t)＝x[t]-M1(t)；
S203: if P is1(t) if the following two conditions of IMF are satisfied simultaneously, it is the first IMF component, otherwise, it is repeated as original signal S201 to S202, resulting in P11(t)：P11(t)＝P1(t)-M11(t), wherein: m11(t) is P1(t) an average curve of the upper and lower envelope lines;
the IMF component satisfies two conditions: in the whole time process, the zero crossing times are equal to the number of extreme points or have at most 1 difference; at any point on the signal, the average value of an upper envelope line defined by the local maximum value and a lower envelope line defined by the local minimum value point is 0, namely the signal is locally symmetrical about a time axis;
s204: repeating the above steps S201 to S203 until P obtained by formula (1) at the k-th screening1k(t) two conditions for the IMF component are satisfied: p1k(t)＝P1(1-k)(t)-M1k(t) (1)
S205: calculation of threshold value S by formula (2)DJudging whether each screening result is an IMF component:
Figure FDA0002311775810000021
wherein: r is the number of sampling points of the homogeneous sensor signal, threshold value SDTaking 0.2 to 0.3;
s206: let C1(t)＝P1k(t) then C1(t) is the first IMF component, which contains the original signal x [ t ]]The least intermediate-cycle IMF component; c is to be1(t) from x [ t ]]Separating out: r1(t)＝x[t]-C1(t)；
S207: r is to be1(t) repeating the above steps S201 to S205n times as a new value, obtaining a signal x [ t ]]N IMF components: rn(t)＝Rn-1(t)-Cn(t)；
S208: when R isn(t) is a monotonic function and the signal x [ t ]]When other components can not be solved, the whole decomposition process is finished; obtaining:
Figure FDA0002311775810000022
2. the synchronous-based homogeneous sensor mutation parameter identification method according to claim 1, wherein the step F of obtaining the instantaneous frequency map and the instantaneous amplitude map by applying HHT transformation to the dominant IMF component of the homogeneous sensor abnormal signal comprises the following steps:
s301: performing Hilbert transform on all IMF components obtained after EMD decomposition; the Hilbert transform of all IMF components C (t) is of the form:
Figure FDA0002311775810000031
wherein λ is an integral variable;
s302: constructing an analysis signal Z (t):
Z(t)＝C(t)+iH(t)＝A(t)eiθ(t)；
s303: amplitude function:
Figure FDA0002311775810000032
s304: argument function:
Figure FDA0002311775810000033
s305: instantaneous frequency:
Figure FDA0002311775810000034
3. the homogeneous sensor mutation parameter identification method based on synchronization as claimed in any one of claims 1-2, wherein the step G of detecting and identifying the homogeneous sensor mutation parameter according to the mutation point, amplitude variation, and instantaneous frequency trend feature information in the instantaneous frequency map and the instantaneous amplitude map comprises:
s306: classifying the mutation parameters in the homogeneous sensor: dip, swell, harmonic, pulse;
s307: and classifying and identifying when the homogeneous sensor undergoes any sudden change according to the sudden change points, amplitude change and instantaneous frequency trend characteristic information of the abnormal signals of different homogeneous sensors in the instantaneous frequency diagram and the instantaneous amplitude diagram.
4. A synchronous homogeneous sensor mutation parameter identification device is characterized by comprising:
a module A: taking T as a sampling period, sampling and quantizing a system detected signal at regular time, and obtaining sampling data under the same sampling frequency;
and a module B: intercepting the sampling data according to cycles to obtain m cycle sequences as follows: y is1(t),Y2(t),…Ym(t) N data points, Y, are contained in each periodic samplei＝[Xi(1),Xi(2),…Xi(N)]Wherein i is 1,2, …, m;
and a module C: constructing a periodic data sequence constructed according to periodic synchronization: the data of each time point in each period data is regarded as a multiple measurement result of the same object after synchronization, namely the data of the time point corresponding to each period is a sample, each point value of the corresponding time point in the data of multiple periods forms a data sequence, and the data is regarded as a one-dimensional data sequence measured by a single homogeneous sensor for multiple times, namely a period data sequence Y constructed according to period synchronization is constructed1′(t),Y2′(t),…YN' (t) each one-dimensional data sequence includes m data points, i.e. Yi′＝[X1(1),X2(2),…Xm(N)]；
A module D: for data sample Yt' performing maximum likelihood estimation-based preprocessing;
and a module E: performing EMD on the mutation-containing homogeneous sensor signal based on combined preprocessing to obtain IMF components of each order;
and a module F: adopting HHT (Hilbert-Huang transform) to obtain an instantaneous frequency graph and an instantaneous amplitude graph for the dominant IMF component of the abnormal signal of the homogeneous sensor;
a module G: detecting and identifying the mutation parameters of the homogeneous sensor according to the mutation points, amplitude variation and instantaneous frequency trend characteristic information in the instantaneous frequency graph and the instantaneous amplitude graph;
the module E adopts EMD to decompose the mutation-containing homogeneous sensor signals based on combined preprocessing to obtain IMF components of each order:
s201: interpolating all local maximum points and all local minimum points of the signal x [ t ] by a cubic spline function, and fitting an upper envelope line and a lower envelope line; x < t > is the homogeneous sensor signal containing the mutation characteristic after combined pretreatment;
s202: calculating the average curve M of the upper and lower envelope lines1(t), if t is 1,2, 3 … m, the signal x [ t ] is sampled]And M1The difference of (t) is P1(t)：P1(t)＝x[t]-M1(t)；
S203: if P is1(t) if both of the following IMF conditions are satisfied, it is the first IMF component, otherwise, it is repeated as the original signal S201 to S202, resulting in P11(t)：P11(t)＝P1(t)-M11(t), wherein: m11(t) is P1(t) an average curve of the upper and lower envelope lines;
the IMF component satisfies two conditions: in the whole time process, the zero crossing times are equal to the number of extreme points or have at most 1 difference; at any point on the signal, the average value of an upper envelope line defined by the local maximum value and a lower envelope line defined by the local minimum value point is 0, namely the signal is locally symmetrical about a time axis;
s204: repeating the above steps S201 to S203 until P obtained by formula (1) at the k-th screening1k(t) two conditions for the IMF component are satisfied: p1k(t)＝P1(1-k)(t)-M1k(t) (1)
S205: in addition to the signal, the threshold value S can be obtained by the equation (2) in the actual calculationDJudging whether each screening result is an IMF component:
Figure FDA0002311775810000051
wherein: r is the number of sampling points of the homogeneous sensor signal, threshold value SDTaking 0.2 to 0.3;
s206: let C1(t)＝P1k(t) then C1(t) is the first IMF component, which contains the original signal x [ t ]]The least intermediate-cycle IMF component; c is to be1(t) from x [ t ]]Separating out: r1(t)＝x[t]-C1(t)；
S207: r is to be1(t) repeating the above steps S201 to S205n times as a new value, obtaining a signal x [ t ]]N IMF components: rn(t)＝Rn-1(t)-Cn(t)；
S208: when R isn(t) is a monotonic function and the signal x [ t ]]When other components can not be solved, the whole decomposition process is finished; obtaining:
Figure FDA0002311775810000052
5. the synchronous-based homogeneous sensor mutation parameter identification device according to claim 4, wherein the module F obtains the instantaneous frequency map and the instantaneous amplitude map by HHT transformation on the dominant IMF component of the homogeneous sensor abnormal signal by:
s301: performing Hilbert transform on all IMF components obtained after EMD decomposition; the hubert transform form of all IMF components c (t) is:
Figure FDA0002311775810000053
wherein λ is an integral variable;
s302: constructing an analysis signal Z (t):
Z(t)＝C(t)+iH(t)＝A(t)eiθ(t)；
s303: amplitude function:
Figure FDA0002311775810000054
s304: argument function:
Figure FDA0002311775810000055
S305：instantaneous frequency:
Figure FDA0002311775810000056
6. the homogeneous sensor mutation parameter identification apparatus based on synchronization as claimed in any one of claims 4 to 5, wherein the process of detecting and identifying the homogeneous sensor mutation parameter according to the mutation point, amplitude variation, and instantaneous frequency trend characteristic information of the instantaneous frequency map and the instantaneous amplitude map in the module G is as follows:
s306: classifying the mutation parameters in the homogeneous sensor: dip, swell, harmonic, pulse;
s307: and classifying and identifying when the homogeneous sensor undergoes any sudden change according to the sudden change points, amplitude change and instantaneous frequency trend characteristic information of the abnormal signals of different homogeneous sensors in the instantaneous frequency diagram and the instantaneous amplitude diagram.
CN201610248401.8A 2016-04-20 2016-04-20 Synchronization-based homogeneous sensor mutation parameter identification method and device Active CN105954695B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201610248401.8A CN105954695B (en) 2016-04-20 2016-04-20 Synchronization-based homogeneous sensor mutation parameter identification method and device

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201610248401.8A CN105954695B (en) 2016-04-20 2016-04-20 Synchronization-based homogeneous sensor mutation parameter identification method and device

Publications (2)

Publication Number Publication Date
CN105954695A CN105954695A (en) 2016-09-21
CN105954695B true CN105954695B (en) 2020-04-17

Family

ID=56917825

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201610248401.8A Active CN105954695B (en) 2016-04-20 2016-04-20 Synchronization-based homogeneous sensor mutation parameter identification method and device

Country Status (1)

Country Link
CN (1) CN105954695B (en)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107885940A (en) * 2017-11-10 2018-04-06 吉林大学 A kind of signal characteristic extracting methods for distributed optical fiber vibration sensing system
CN110823560A (en) * 2018-08-07 2020-02-21 上海华依科技集团股份有限公司 Data acquisition method for automatic transmission offline test system
CN111092875B (en) * 2019-12-12 2022-05-24 南京国电南自电网自动化有限公司 Transmission and transformation operation and inspection platform Internet of things edge information transmission compression method and system
CN111242366B (en) * 2020-01-08 2023-10-31 广东技术师范大学 EMD method and device capable of being used for processing signals in real time
CN115659160B (en) * 2022-12-28 2023-06-16 北京中航路通科技有限公司 Data quality measurement method for digital twin model optimization

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6862558B2 (en) * 2001-02-14 2005-03-01 The United States Of America As Represented By The Administrator Of The National Aeronautics And Space Administration Empirical mode decomposition for analyzing acoustical signals
CN101551433B (en) * 2009-05-05 2011-04-06 昆明理工大学 Distribution network feed out circuit fault circuit selection method by using HHT detection technology
KR101077541B1 (en) * 2010-08-20 2011-10-27 주식회사 후크앤타임 A sensor system based on time synchronization using interactive wireless communication
CN103576167B (en) * 2013-10-28 2015-09-23 中国科学院国家授时中心 Based on the cycle-slip detection and repair method of HHT and support vector machine
CN105426583B (en) * 2015-11-03 2019-11-26 国网江西省电力科学研究院 It is a kind of based on synchronous homogeneity sensor method for amalgamation processing
CN105510783A (en) * 2015-12-30 2016-04-20 华北电力大学 Switchgear partial discharge detecting system based on ultrasonic signal

Also Published As

Publication number Publication date
CN105954695A (en) 2016-09-21

Similar Documents

Publication Publication Date Title
CN105954695B (en) Synchronization-based homogeneous sensor mutation parameter identification method and device
CN106845010B (en) Low-frequency oscillation dominant mode identification method based on improved SVD noise reduction and Prony
Zygarlicki et al. A reduced Prony's method in power-quality analysis—parameters selection
Ji et al. Disturbance detection, location and classification in phase space
CN114153888A (en) Abnormal value detection method and device for time series data
CN107238756B (en) A kind of intelligence distribution transformer terminals impact load Identifying Methods of Harmonic Source
CN110648088A (en) Electric energy quality disturbance source judgment method based on bird swarm algorithm and SVM
CN112213687B (en) Gateway electric energy meter data anomaly detection method and system based on pseudo-anomaly point identification
CN110632546A (en) Electronic transformer credibility evaluation method and device based on whole-network-domain evidence set
CN106137184B (en) Electrocardiosignal QRS complex detection method based on wavelet transformation
CN110940933B (en) Comprehensive calculation method for measuring rising edge starting time of steep pulse
CN109410178A (en) A kind of workpiece crack detection method and system
CN115993511A (en) Partial discharge type high-precision detection and identification device, method and equipment
KR101807280B1 (en) Apparatus and Method for detecting accidents in Power Systems using PMU signal
CN111273126B (en) Power distribution network topology rapid sensing method
Ivanov et al. Determination of frequency in three-phase electric circuits with autocorrelated noise
CN117330906A (en) Equipment arc fault detection method, device, equipment and storage medium
CN108334822B (en) Kalman and modified wavelet transform filtering method based on electric vehicle charging nonlinear load characteristics
Held et al. Parameter optimized event detection for NILM using frequency invariant transformation of periodic signals (FIT-PS)
CN111075660B (en) Frequency domain analysis method, device and equipment for monitoring variables of wind turbine generator
CN114091593A (en) Network-level arc fault diagnosis method based on multi-scale feature fusion
CN114184838A (en) Power system harmonic detection method, system and medium based on SN mutual convolution window
CN107506824B (en) Method and device for detecting bad observation data of power distribution network
CN108614147A (en) Voltage fluctuation detection method and its Source of Gateway Meter
Lugnani et al. Real-time Coherency Identification using a Window-Size-Based Recursive Typicality Data Analysis

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant