A fusion terminal and electric energy meter power line harmonic content algorithm based on phase difference method and prony algorithm

Abstract

The application discloses a kind of fusion terminal and electric energy meter electric power line harmonic content algorithm based on phase difference method and Prony algorithm, it is related to power quality monitoring technical field, including steps S1: signal synchronous acquisition and pre-processing, fusion terminal passes through control electric energy meter acquisition signal, eliminate direct current component and high-frequency noise;Step S2: phase difference method preliminary examination, adopt double threshold zero-crossing detection positioning fundamental wave period, calculate adjacent period phase difference by hilbert transform, deduce fundamental wave frequency;Step S3: Prony algorithm fine correction, construct Prony model, by Hankel matrix, SVD decomposition and characteristic equation solution;Step S4: adaptive fusion and output.The application adopts four-level progressive architecture of signal pre-processing, phase difference method preliminary examination, Prony algorithm fine correction and adaptive fusion output, and Prony algorithm is accurately corrected to preliminary examination result, makes up the deficiency of phase difference method in high-order harmonic and anti-interference, realizes the unity of speed and accuracy.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of power quality monitoring, in particular to a fusion terminal and electric energy meter power line harmonic content algorithm based on a phase difference method and a Prony algorithm. BACKGROUND

[0002] Harmonic pollution of power systems has become a key factor affecting power quality and safe and stable operation of power grids. With the wide access of nonlinear loads (such as frequency converters, electric arc furnaces, new energy power generation equipment, etc.), the harmonic components in the power grid are becoming more and more complex, which not only leads to increased line loss, equipment overheating, and protection misoperation, but also interferes with the normal work of precision instruments. Therefore, accurate and rapid harmonic detection of power lines is a necessary prerequisite for realizing fine management of smart grids and ensuring power supply reliability.

[0003] Traditional harmonic detection methods, such as Fast Fourier Transform (FFT), are widely used due to their high computational efficiency and simple implementation. However, the FFT method is limited by spectral leakage and the fence effect, and its detection accuracy significantly decreases in non-synchronous sampling or frequency fluctuation scenarios, especially for high-order harmonics or interharmonics. The phase difference method is based on time-domain signal zero-crossing points and phase information, and has the advantages of fast initial estimation of fundamental and main harmonic frequencies, and acceptable noise resistance. However, its accuracy is easily affected by noise and signal distortion, and it lacks the ability to detect high-order weak harmonics. The Prony algorithm is sensitive to noise, and in the solving process, it is easy to fall into local optimization or appear false components, and its convergence depends heavily on the initial parameters, resulting in the existing detection method being unable to achieve both accuracy and speed. SUMMARY

[0004] The present application aims to provide a fusion terminal and electric energy meter power line harmonic content algorithm based on a phase difference method and a Prony algorithm to solve the problems raised in the background.

[0005] To achieve the above-mentioned purpose, the present application provides the following technical solution: a fusion terminal and electric energy meter power line harmonic content algorithm based on a phase difference method and a Prony algorithm, comprising the following specific steps:

[0006] Step S1: signal synchronous acquisition and preprocessing, the fusion terminal acquires signals by controlling the electric energy meter, eliminates direct current components and high-frequency noise, and performs SNR verification to ensure that the signal-to-noise ratio is greater than or equal to 30 dB;

[0007] Step S2: initial detection by phase difference method, double-threshold zero-crossing detection is used to locate the fundamental period, the Hilbert transform is used to calculate the phase difference between adjacent periods, the fundamental frequency is derived, and the initial parameters of harmonics are estimated based on the integer multiple relationship of the fundamental;

[0008] Step S3: fine calibration by Prony algorithm, constrained by initial detection results, build Prony model, solve by Hankel matrix, SVD decomposition and characteristic equation, iterative optimization of harmonic parameters, for improving the detection accuracy of high-order harmonics;

[0009] Step S4: adaptive fusion and output, dynamically allocate weights according to signal-to-noise ratio and harmonic order, fuse the advantages of the two algorithms, and output the final parameters.

[0010] Preferably, the steps of step S1 are as follows:

[0011] S11: signal acquisition, the fusion terminal acquires 3-cycle voltage or current signals at a sampling frequency of 2560Hz by controlling the electric energy meter, and the original voltage or current discrete signals acquired by the electric energy meter are represented by ;

[0012] S12: remove the direct current component, use the sliding window average method to remove the direct current component, for eliminating the interference of direct current on phase calculation, the sliding window average method takes 64 points per cycle as the sliding window, calculates the direct current component in real time and removes it from the original signal, and the calculation method is as follows:

[0013]

[0014]

[0015] When k < 63, the window is supplemented with the previous period historical data, the direct current component calculation is effective for the first period, and is the original signal from which the direct current component is removed;

[0016] S13: anti-aliasing FIR filter, suppress high-frequency noise by an 8th-order Hann window FIR filter, the FIR filter cutoff frequency The filter has a passband ripple ≤0.1dB and a stopband attenuation ≥60dB, wherein it covers 50th harmonic 2500Hz, reserves 100Hz protection band, and the specific parameters are as follows:

[0017] The unit impulse response is:

[0018]

[0019] Filter implementation: based on FPGA parallel convolution operation, convolution operation is synchronized by 8 multipliers, and the filter formula is:

[0020]

[0021] Based on FPGA parallel convolution implementation, the processing delay is ≤4μs;

[0022] S14: Check, check the signal by SNR, verify the signal effective, get the determined de-noising and de-DC signal , provide clean input signal with SNR≥30dB for subsequent algorithm.

[0023] Preferably, the SNR checking step in step S14 is as follows:

[0024] S141: Fundamental preliminary fitting: assuming the fundamental frequency is 50Hz, the fitting expression is:

[0025]

[0026] Its Take 0.707 times of the maximum value, 0.707 times is the effective value estimation;

[0027] S142: Power calculation:

[0028]

[0029] S143: SNR calculation:

[0030]

[0031] S143: Abnormal processing, the abnormal processing mechanism process is:

[0032] If SNR<30dB, send resampling instruction;

[0033] If SNR<30dB for 3 times in a row, trigger device abnormal alarm, and the alarm information is uploaded to the master station through the 4G module, and the event record is stored into the database.

[0034] Preferably, the phase difference method initial detection in step S2 is based on the linear change characteristic of signal phase, locks the fundamental frequency and harmonic initial parameters, and provides initial constraints for Prony algorithm, wherein the detection step is as follows:

[0035] S21: Input signal pre-processing filtered signal ;

[0036] S22: Time domain zero-crossing detection and period calculation, the detection process is as follows:

[0037] Threshold setting for double threshold zero-crossing detection: positive zero-crossing threshold , negative zero-crossing threshold ;

[0038] The zero-crossing time is calculated by linear interpolation: when detecting , the accurate zero-crossing time calculation formula is:

[0039]

[0040] wherein, is the zero-crossing time;

[0041] Period calculation: take the adjacent two positive zero-crossing time , The fundamental period is ;

[0042] The fundamental period T1 is constrained for effectiveness: ; if it exceeds the range, expand the detection window to 5 periods to take the average value to ensure the reliability of period calculation;

[0043] S23: Hilbert transform and analytic signal construction, the Hilbert transform equation is:

[0044]

[0045] Wherein, the Hilbert transform is realized by FFT acceleration to reduce the computational complexity;

[0046] Analytic signal and instantaneous phase extraction:

[0047]

[0048] Wherein A(k) is the instantaneous amplitude, is the instantaneous phase;

[0049] The phase unwrapping formula is:

[0050]

[0051] Fundamental frequency calculation: take the phase difference of adjacent periods , then the fundamental frequency is:

[0052]

[0053] The fundamental frequency resolution is ≤0.01Hz;

[0054] S24: parameter estimation, based on the characteristics of power grid harmonics, the initial parameters of harmonics are estimated, and the integer frequency is erected, the erection process is as follows:

[0055] The frequency calculation of the initial frequency of harmonics is:

[0056] ;

[0057] Based on the statistical calculation of the massive measured data of State Grid, the initial amplitude of harmonics is calculated under general scenarios , wherein The initial amplitude of the harmonic is used for the estimation of the effective value of the fundamental wave.

[0058] S25: result output, output the preliminary inspection parameters to the Prony module.

[0059] Preferably, the step S3 comprises the following steps:

[0060] S31: Prony model construction, the model is established by taking the preliminary inspection parameters as constraints, and the filtered signal is input and the phase difference method preliminary inspection parameters , , a P=50 order Prony model is constructed, and the model formula is:

[0061] wherein p is the model order, is a complex amplitude, is a fitting residual;

[0062] S32: Hankel matrix and SVD decomposition, the Hankel matrix is constructed by selecting the first 2p points of the filtered signal, and the Hankel matrix is constructed. The matrix is as follows:

[0063]

[0064]

[0065]

[0066]

[0067] The singular value SVD decomposition is the characteristic equation solution of the Prony model. The characteristic equation solution is realized by calling the LAPACK library dgesvd function to realize matrix decomposition, separating the effective signal and the noise, and the decomposition formula is:

[0068]

[0069] wherein is a diagonal matrix, and the first p non-zero singular values are retained, and the small singular values corresponding to the noise are removed;

[0070] S33: characteristic equation solution and parameter iteration, used for extracting accurate harmonic parameters;

[0071] S34: parameter accuracy verification, the amplitude error of the 1-13th harmonic is less than or equal to 3%, and the error of the 21-50th harmonic is less than or equal to 12%.

[0072] Preferably, the step of extracting accurate harmonic parameters in the step S33 is:

[0073] S331: coefficient vector extraction, from Matrix last column extraction coefficient , build characteristic equation ;

[0074] S332: eigenvalue solution, through polynomial root algorithm, get eigenvalue :

[0075]

[0076] Where no attenuation, so ;

[0077] S333: parameter update, frequency, complex amplitude and amplitude and phase calculation, where the frequency The formula is:

[0078] Complex amplitude for building Vandermonde matrix , through the least square method to solve ;

[0079] Amplitude is , the phase is ;

[0080] S334: iterative optimization of harmonic parameters, the iteration termination condition is , or the number of iterations ≥ 30 times.

[0081] Preferably, the step S4 comprises the following steps:

[0082] S41: input phase difference method initial detection parameters And and Prony fine calibration parameters , and ;

[0083] S42: weight calculation, calculate the signal-to-noise ratio weight factor, SNR ≥ 20dB is effective, the signal-to-noise ratio weight factor calculation formula is:

[0084]

[0085] Where the higher the SNR, the stronger the reliability of the phase difference method, the greater the weight; SNR = 40dB = 0.67, SNR = 25dB ;

[0086] Calculate the harmonic number weight factor, the formula is:

[0087]

[0088] The larger n is, the weaker the harmonic signal is, the more obvious the noise resistance of Prony algorithm is, and the larger the weight is;

[0089] The maximum weight is the weight of the phase difference method and the weight of Prony algorithm, and the weight calculation formula of the phase difference method is:

[0090]

[0091] The weight calculation formula of Prony algorithm is:

[0092]

[0093] S43: Weighted fusion calculation, through the connection with the weight calculation, the output of the optimal parameter is realized, and the calculation method is:

[0094]

[0095]

[0096]

[0097] Wherein, The initial value is 0, which is gradually corrected through iteration, and the phase accuracy after fusion is ≤±0.5°

[0098] S44: Post-processing parameter output, harmonic content rate The calculation method is:

[0099]

[0100] Wherein, the special treatment is carried out on When , It is recorded as 0.00% and marked as not detected.

[0101] Preferably, the total harmonic distortion rate in S44 is calculated as:

[0102]

[0103] Wherein, the results of the total harmonic distortion rate and the harmonic content rate are kept to 2 decimal places;

[0104] The total harmonic distortion rate and the harmonic content rate data are uploaded according to the specified format through the 4G module at a regular interval of 1 minute, and the total harmonic distortion rate When the limit is exceeded, it is uploaded in real time for data security.

[0105] A fusion terminal for the above-mentioned algorithm includes a main controller, a memory, a 2500V isolated RS485 interface, an ME909s-8214G module, and a BeiDou / GPS dual-mode timing module, with an operating temperature range of -40℃ to 85℃.

[0106] The smart energy meter, used in conjunction with the aforementioned converged terminal, employs an STM32F407 microprocessor, an AD7606-16 ADC chip, a 1000:1 voltage divider network, and a Rogowski coil current sensor. The energy meter integrates DS18B20 temperature compensation and supports 2560Hz synchronous sampling.

[0107] The technical effects and advantages of this invention are as follows:

[0108] (1) The present invention adopts a four-level progressive architecture of signal preprocessing, initial detection by phase difference method, fine correction by Prony algorithm and adaptive fusion output. The phase difference method provides the Prony algorithm with an initial constraint close to the true value, which effectively avoids the problem that the Prony algorithm gets stuck in local optima and converges slowly due to large initial value deviation. The Prony algorithm accurately corrects the initial detection result, which makes up for the shortcomings of the phase difference method in high harmonics and anti-interference. The two work together to achieve the unity of speed and accuracy and solve the problem of insufficient detection capability of traditional single algorithm.

[0109] (2) This invention precisely controls sampling synchronization through protocol and adopts preprocessing steps such as sliding window averaging to remove DC, anti-aliasing FIR filtering and SNR verification to provide high-quality, high signal-to-noise ratio input signals for subsequent analysis, ensuring detection accuracy from the source. The phase difference method initial detection obtains high-precision fundamental frequency through double threshold zero-crossing detection and Hilbert transform, and sets the initial amplitude of harmonics based on the integer multiple frequency assumption and measured data statistics, controlling the initial parameter deviation of the Prony algorithm within 5%, improving the convergence probability of the Prony algorithm, reducing the number of iterations, and significantly shortening the calculation time.

[0110] (3) This invention innovatively introduces a dynamic weight allocation mechanism based on signal-to-noise ratio and harmonic order. In high signal-to-noise ratio and low harmonic scenarios, the phase difference method is given higher weight to ensure speed. In low signal-to-noise ratio and high harmonic scenarios, the Prony algorithm is given higher weight to ensure accuracy. This mechanism enables the algorithm to adapt to different load characteristics and output the optimal harmonic parameters. The accuracy of the harmonic parameters after fusion processing is improved. Attached Figure Description

[0111] Figure 1 This is a logic diagram of the present invention.

[0112] Figure 2 This is a detailed flowchart of the preprocessing process of the present invention.

[0113] Figure 3 Flow chart for phase difference method detection of the present application.

[0114] Figure 4 Flow chart for time domain zero-crossing detection and period calculation of the present application.

[0115] Figure 5 Flow chart for Hilbert transform and analytic signal construction of the present application.

[0116] Figure 6 Flow chart for frequency calculation and parameter estimation of the present application.

[0117] Figure 7 Flow chart for Prony algorithm correction of the present application.

[0118] Figure 8 Flow chart for Prony algorithm fine correction step of the present application.

[0119] Figure 9 Flow chart for output and processing step of the present application.

[0120] Figure 10 Flow chart for adaptive fusion of the present application. DETAILED DESCRIPTION

[0121] The technical solutions in the embodiments of the present application will be clearly and completely described below with reference to the drawings in the embodiments of the present application. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative labor fall within the scope of protection of the present application.

[0122] The present application provides a fusion terminal and electric energy meter power line harmonic content algorithm based on phase difference method and Prony algorithm as shown in Figures 1-10 The algorithm adopts a four-level progressive architecture of "signal preprocessing → phase difference method preliminary detection → Prony algorithm fine correction → adaptive fusion output", and each link is closely connected through mathematical models, including the following specific steps:

[0123] Step S1: signal synchronous acquisition and preprocessing, the fusion terminal acquires signals by controlling the electric energy meter, eliminates the direct current component and high frequency noise, and passes the SNR check to ensure that the signal-to-noise ratio is greater than or equal to 30dB;

[0124] Step S2: phase difference method preliminary detection, double-threshold zero-crossing detection is used to locate the fundamental period, the phase difference of adjacent periods is calculated through Hilbert transform, the fundamental frequency is derived, and the initial parameters of harmonics are estimated based on the integer multiple relationship of the fundamental frequency. The phase difference method is used for preliminary detection to provide initial constraints for Prony algorithm and avoid falling into local optimum.

[0125] Step S3: fine calibration by Prony algorithm, with the initial detection result as the constraint, build Prony model, solve by Hankel matrix, SVD decomposition and characteristic equation, iterative optimization of harmonic parameters, for the improvement of high harmonic detection accuracy;

[0126] Step S4: adaptive fusion and output, dynamically allocate weights according to signal-to-noise ratio and harmonic order, fuse the advantages of the two algorithms, and output the final parameters.

[0127] The specific steps of step S1 signal acquisition and preprocessing are as follows:

[0128] S11: signal acquisition, the fusion terminal controls the electric energy meter to collect 3-cycle voltage or current signals at a sampling frequency of 2560Hz through DL / T645-2007 protocol, and the original voltage or current discrete signals collected by the electric energy meter are sent to the fusion terminal through the RS485 bus. The signal acquisition process is as follows:

[0129] Sampling control: the fusion terminal sends a sampling instruction to the electric energy meter, the instruction contains a UTC timestamp, and is accurate to 1ms, ensuring that the sampling synchronization deviation of the whole network is less than or equal to 1ms;

[0130] Sampling parameters: sampling frequency , collect 3 cycles continuously, sample N=64 points per cycle, sampling interval , the sampling frequency standard is according to DL / T1475-2015 clause 5.2.2;

[0131] The expression of the original signal of the collected voltage or current signal is as follows:

[0132]

[0133] Wherein is the fundamental amplitude / frequency / phase; is the n-th harmonic parameter; is the DC component; is the Gaussian white noise; =0,1,...,191.

[0134] S12: DC component removal, use sliding window average method to remove DC component, for eliminating the interference of DC component on phase calculation, since the DC component will cause the baseline deviation of phase difference calculation, and make the fundamental cycle positioning error increase by 3-5 times, which must be removed first;

[0135] The sliding window average method takes 64 points per cycle as the sliding window, calculates the DC component in real time and removes it from the original signal, and the calculation method is as follows:

[0136] ​​

[0137]

[0138] Wherein k < 63, the window fills the previous period history data, determines the first period direct current component calculation effective, and The processing effect of the sliding window average method for the original signal without the direct current component is that the error of the fundamental wave phase detection is reduced from ± 5° to within ± 0.5°, which is used to meet the accuracy requirement of subsequent phase difference calculation.

[0139] S13: Anti-aliasing FIR filtering, anti-aliasing FIR filtering is connected with the direct current removal, 8-order Hanning window FIR filter is used to suppress high-frequency noise and avoid spectrum aliasing, the cut-off frequency of the FIR filter is The filter has a passband ripple of ≤0.1dB and a stopband attenuation of ≥60dB, wherein the 50th harmonic 2500Hz is covered, and a protection band of 100Hz is reserved, and the specific parameters are as follows:

[0140] The unit impulse response is:

[0141]

[0142] Filtering implementation: based on FPGA parallel convolution operation, the convolution operation is synchronized by 8 multipliers, and the filtering formula is:

[0143]

[0144] The filtered signal is realized based on FPGA parallel convolution, and the performance index is a passband ripple of ≤0.1dB (0-2500Hz) and a stopband attenuation of ≥60dB (3000Hz and above), a processing delay of ≤4μs, and the anti-aliasing characteristic test of the State Grid Electric Power Research Institute.

[0145] S14: verification, SNR verification is connected with the filtering in S13, the effectiveness of the verification signal is improved, invalid calculation is avoided, the signal is verified by SNR, and after the verification signal is effective, the determined denoising and direct current removal signal is obtained, which provides a clean input signal with a signal-to-noise ratio of ≥30dB for subsequent algorithms.

[0146] The signal-to-noise ratio is calculated by the fundamental fitting residual, and if the threshold is not met, re-sampling is triggered to ensure the quality of the input signal, wherein the SNR verification steps are as follows:

[0147] S141: preliminary fitting of the fundamental wave: assuming that the fundamental wave frequency is 50Hz, the fitting expression is:

[0148]

[0149] The is taken 0.707 times of the maximum value, 0.707 times is the effective value estimation.

[0150] S142: Signal power calculation:

[0151]

[0152] S143: SNR calculation:

[0153]

[0154] S143: Abnormal processing, the abnormal processing mechanism flow is:

[0155] If SNR < 30dB, send resampling instruction, the sampling instruction is resampled at most 3 times;

[0156] If SNR < 30dB for 3 times in succession, trigger device abnormal alarm, and alarm information is uploaded to the master station through the 4G module, and the event record is stored into the database;

[0157] Through the denoising and removing direct current of the signal , it is convenient to meet the quality requirements of the subsequent algorithm to the input signal, so that its SNR can be greater than or equal to 30dB.

[0158] The phase difference method in step S2 is based on the linear variation characteristics of the signal phase, locks the fundamental frequency and the initial parameters of the harmonic, provides the initial constraint for the Prony algorithm, provides the initial parameters close to the true value for the Prony algorithm, and shortens the convergence time, wherein the phase difference method preliminary detection step is as follows:

[0159] S21: Input signal pre-processing filtered signal ;

[0160] S22: Time domain zero crossing detection and period calculation, in order to solve the misjudgment problem of single threshold under low signal amplitude, the detection process is as follows:

[0161] Threshold setting for double threshold zero crossing detection: positive zero crossing threshold , negative zero crossing threshold ;

[0162] The zero crossing time is calculated by linear interpolation: when is detected, the accurate zero crossing time calculation formula is:

[0163] ;

[0164] Period calculation: take the adjacent two positive zero crossing times , , and the fundamental period is ;

[0165] Effectiveness constraints are applied to the fundamental period T1: The validity constraint corresponds to the power grid frequency fluctuation range of 50Hz±2.5Hz; if it exceeds the range, the detection window is expanded to take the average of 5 cycles to ensure the reliability of the cycle calculation.

[0166] The fundamental frequency is linked to the period positioning; the fundamental frequency calculation derives the accurate fundamental frequency through phase difference. ;

[0167] S23: Hilbert Transform and Analytic Signal Construction, the Hilbert transform equation is:

[0168]

[0169] The Hilbert transform is accelerated using the FFT to reduce computational complexity.

[0170] Analytical signal and instantaneous phase extraction:

[0171]

[0172] Where A(k) is the instantaneous amplitude, Instantaneous phase (range of values) );

[0173] The signal phase expansion formula is:

[0174]

[0175] Signal phase is to be eliminated A jump is used to obtain a continuous phase curve;

[0176] Fundamental frequency calculation: Take the phase difference between adjacent periods (64 points apart). Then the fundamental frequency is:

[0177]

[0178] Fundamental frequency Resolution ≤0.01Hz, fundamental frequency The precision is retained to 4 decimal places, which meets the requirements of Clause 6.2 of DL / T1475-2015.

[0179] S24: Parameter estimation. Based on the derivation of power grid harmonic characteristics, the initial harmonic parameters are estimated to achieve integer multiple frequency setting. Integer multiple frequency setting conforms to the harmonic generation law of the power system. The harmonics generated by nonlinear loads are mostly integer multiples of the fundamental frequency. The process of integer multiple frequency setting is as follows:

[0180] The initial frequency of the harmonic is calculated as follows:

[0181] ;

[0182] Based on the calculation of the initial amplitude of the harmonic based on the massive measured data statistics of the State Grid, the initial amplitude of the harmonic is set , wherein , the initial amplitude of the harmonic is used for the estimation of the effective value of the fundamental wave;

[0183] The Prony algorithm is sensitive to the initial value, and if the initial value deviates from the true value by more than 20%, the convergence probability is reduced to 60%. The initial deviation is controlled within 5% in this step, the convergence probability is improved to 100%, the iteration number is reduced from 50 times to within 30 times, and the calculation time is reduced by 40%;

[0184] Special scene adaptation: in the industrial load containing the frequency converter, the initial amplitudes of the characteristic harmonics of 6, 12 and 18 times are adjusted to , so as to improve the convergence speed of the Prony algorithm for the characteristic harmonics;

[0185] S25: result output, outputting the initial detection parameters to the Prony module, outputting the fundamental frequency , the initial amplitude of the harmonic and the frequency , the resolution of the fundamental frequency provides an initial constraint for the Prony algorithm, so that the iteration number is reduced from 50 times to within 30 times.

[0186] The accurate correction of the Prony algorithm in step S3 includes the following steps:

[0187] S31: Prony model construction, the model is established by taking the initial detection parameters as constraints, the filtered signal and the initial detection parameters of the phase difference method , are input, a P=50 order Prony model is constructed, the accurate harmonic amplitude , frequency and phase after the accurate correction of the Prony algorithm, the model adapts to the steady-state harmonic characteristics of the power grid, the characteristics are very small, and the attenuation factor can be ignored , and the formula of the established model is:

[0188] , wherein p is the model order, the model order p can be dynamically adjusted, the basic value of the model order p is 50, covering the harmonics of 2-50 times; when , the order is increased to p=52, which is used to cope with the 3rd harmonic frequency division problem;

[0189] is a complex amplitude, which contains amplitude and phase information;

[0190] And To fit the residual, by , to ensure the fitting accuracy.

[0191] S32: Hankel matrix and SVD decomposition, Hankel matrix and Prony model for the interface, easy to convert nonlinear problems into linear solution;

[0192] Hankel matrix construction selected filter signal before 2p points, the construction of p order Hankel matrix, Hankel matrix construction for enhancing the anti-noise of signal, the matrix as follows:

[0193]

[0194]

[0195]

[0196]

[0197] In which singular value (SVD) decomposition, namely Prony model Hankel matrix characteristic equation solution, the characteristic equation solution by calling LAPACK library dgesvd function to realize matrix decomposition, separation of effective signal and noise, the decomposition formula is:

[0198]

[0199] In which is a diagonal matrix, and the singular values are arranged in descending order, and the first p non-zero singular values are retained, and the small singular values corresponding to the noise are removed.

[0200] S33: characteristic equation solution and parameter iteration, used to extract accurate harmonic parameters, and the SVD in S32 step, to extract accurate harmonic parameter step is:

[0201] S331: coefficient vector extraction, from the last column of the matrix to extract the coefficient , to build the characteristic equation: ;

[0202] S332: eigenvalue solution, by polynomial root algorithm, the eigenvalue is:

[0203]

[0204] In which the eigenvalue has no attenuation, so ;

[0205] S333: parameter update, to calculate the frequency, complex amplitude and amplitude and phase, in which the frequency The calculation formula is:

[0206] The complex amplitude is used to construct a Vandermonde matrix , and is solved by least square method ;

[0207] The amplitude and phase are respectively , ;

[0208] S334: Iterative optimization of harmonic parameters, and the iteration termination condition is or the iteration number is greater than or equal to 30 times, wherein the condition is used for frequency convergence, or the iteration number is greater than or equal to 30 times for avoiding infinite iteration;

[0209] S34: Parameter accuracy verification, the amplitude error of 1-13th harmonic is less than or equal to 3%, and the error of 21-50th harmonic is less than or equal to 12%, and the accuracy of harmonic parameters is better than the error upper limit of "1-13th ≤5%, 21-50th ≤15%" specified in DL / T1475-2015.

[0210] The Prony algorithm is used to accurately correct and compensate for the weak anti-interference ability of the phase difference method, and to improve the detection accuracy of high-order harmonics.

[0211] The adaptive fusion and output in step S4 include the following steps:

[0212] S41: Input the initial detection parameters of the phase difference method and and the fine correction parameters of Prony , , ;

[0213] S42: Weight calculation, the weight calculation is based on the dynamic allocation of weights based on signal characteristics, and is used to calculate the influence of signal-to-noise ratio and harmonic order on detection accuracy;

[0214] The signal-to-noise ratio weight factor is calculated, which is effective when SNR≥20dB, and the calculation formula of the signal-to-noise ratio weight factor is:

[0215]

[0216] The higher the SNR, the stronger the reliability of the phase difference method, and the greater the weight; when SNR=40dB =0.67, and when SNR=25dB ;

[0217] The harmonic order weight factor is calculated, and the formula is:

[0218]

[0219] The greater n is, the weaker the harmonic signal is, the more obvious the Prony algorithm is against noise, and the greater the weight is; when n=3 =0.94, when n=30 =0.55;

[0220] The maximum weight is the weight of the phase difference method and the weight of the Prony algorithm, and the weight calculation formula of the phase difference method is:

[0221]

[0222] The weight calculation formula of the Prony algorithm is:

[0223]

[0224] wherein, is the phase difference method weight of the harmonic signal, is the Prony algorithm weight of the harmonic signal.

[0225] Phase difference method weight and Prony algorithm weight embodiment one:

[0226] Residential load (SNR=40dB, n=3): , (residential load mainly uses the phase difference method, and the speed is considered);

[0227] Industrial load (SNR=25dB, n=30): =0.17*0.55≈0.09, =0.91 (industrial load mainly uses Prony, to ensure accuracy).

[0228] S43: Weighted fusion calculation, through the connection with the weight calculation, the output of the optimal parameter is realized, and the calculation mode is:

[0229]

[0230]

[0231]

[0232] wherein, is initially set to 0, and is gradually corrected through iteration, and the phase accuracy after fusion is ≤±0.5°

[0233] S44: Post-processing parameter output, through the weighted fusion calculation, used for the calculation of the harmonic containing filter, and the harmonic containing rate The calculation mode is:

[0234]

[0235] Wherein the is specially processed when , is recorded as 0.00% and marked as not detected, in line with Q / GDW11645-2016 Article 6.3.2;

[0236] The total harmonic distortion rate is calculated as follows:

[0237]

[0238] Wherein the total harmonic distortion rate and the harmonic content rate are rounded to two decimal places.

[0239] The total harmonic distortion rate and the harmonic content rate and other data are uploaded through the 4G module at a regular interval of 1 minute in the prescribed format, the format conforms to DL / T645-2007, the data identifiers are 0x0701 (voltage harmonic) and 0x0702 (current harmonic), each harmonic point is 3 bytes, 3 bytes are 2 bytes of effective value + 1 byte of times, the total harmonic distortion rate is uploaded in real time when it exceeds the limit, DL / T698.45 protocol + AES-128 encryption is adopted to ensure data security and the safety of harmonic signals;

[0240] Data storage: the original data is stored for 90 days through the memory, the fusion processing result is stored for 365 days, and SHA-1 check code is added to each record to ensure integrity.

[0241] A fusion terminal is used to implement the above algorithm, which comprises a main controller TIAM3352, a 1GB DDR3+8G BeMMC memory, a 2500V isolation RS485 interface, a ME909s-8214G module and a Beidou / GPS dual-mode timing module, and the working temperature range is -40℃~85℃.

[0242] An intelligent electric energy meter is used in cooperation with the fusion terminal, which adopts an STM32F407 microprocessor, an AD7606-16 ADC chip, a 1000:1 voltage division network and a Rogowski coil current sensor, the electric energy meter is integrated with a DS18B20 temperature compensation, and supports 2560Hz synchronous sampling.

[0243] In the second embodiment, the hardware platform configuration of the fusion terminal and the intelligent electric energy meter is as follows:

[0244] The fusion terminal is the core processing unit of the algorithm. The specific configuration is shown in Table 1 as follows:

[0245]

[0246] The intelligent electric energy meter is the signal acquisition unit. The specific configuration is shown in Table 2 as follows:

[0247]

[0248] Example three, equipment debugging and standard process as follows:

[0249] a, equipment calibration, FLUKE6100A harmonic standard source for 0.05 level, Kato WT3000 power analyzer for 0.1 level;

[0250] b, calibration point setting as follows:

[0251] Fundamental frequency: 49.5Hz, 50.0Hz, 50.5Hz;

[0252] Harmonic content: 2 times 5%, 3 times 3%, 5 times 2%, 7 times 1%;

[0253] Temperature condition: 25℃ (normal temperature), -25℃, +70℃ (extreme working condition);

[0254] c, calibration steps:

[0255] Equipment preheating: standard source preheating 1 hour, preheating 30 minutes of the equipment to be tested;

[0256] Data acquisition: set the standard source output point by point, collect 10 groups of data after 5s stabilization;

[0257] Error calculation: is the measured value, is the standard source output value);

[0258] Coefficient update: if Threshold, update the temperature segmented calibration coefficient table (every 5℃ a segment), realize full temperature range compensation through linear interpolation;

[0259] d, on-site debugging;

[0260] Communication debugging: use DL / T645-2007 debugging software to send sampling instructions, verify the response time of electric energy meter ≤10ms, data transmission packet loss rate ≤0.1%.

[0261] Algorithm debugging: test Prony algorithm iteration times (target ≤30 times) in industrial field (including frequency converter, electric arc furnace), ensure single cycle processing time ≤200ms.

[0262] Abnormal processing debugging: simulate signal interruption (disconnect voltage / current connection), verify that the device triggers alarm within 1s, and the breakpoint continuation function is normal (cache the latest 1000 records of data).

[0263] e, compliance and performance verification report

[0264] The standard compliance list is shown in Table 3 as follows:

[0265]

[0266] The performance test results are shown in Table 4 as follows:

[0267]

[0268] Through the initial detection by the phase difference method, the fine calibration by the Prony algorithm, and the adaptive fusion and output, the problems of low precision, poor real-time performance and anti-interference performance of traditional harmonic detection algorithms are solved, and the detection precision and detection efficiency are improved.

[0269] Finally, it should be noted that: the above is only the preferred embodiment of the present application, and is not intended to limit the present application, although the present application has been described in detail with reference to the foregoing embodiments, for those skilled in the art, the technical solutions recorded in the foregoing embodiments can still be modified, or some technical features therein can be replaced by equivalents, any modification, equivalent replacement, improvement, etc. made within the spirit and principles of the present application shall be included within the scope of protection of the present application.

Claims

1. A fusion terminal and electric energy meter power line harmonic content algorithm based on the phase difference method and the Prony algorithm, characterized in that, The method comprises the following specific steps: Step S1: signal synchronization acquisition and pretreatment, the fusion terminal acquires the signal through the control of the electric energy meter, eliminates the direct current component and high frequency noise, and passes the SNR check to ensure that the signal-to-noise ratio is greater than or equal to 30 dB; Step S2: initial detection by phase difference method, double threshold zero-crossing detection is used to locate the fundamental wave period, the phase difference of adjacent periods is calculated by Hilbert transform, the fundamental wave frequency is derived, and the initial parameters of harmonics are estimated based on the integer multiple relationship of the fundamental wave; Step S3: Prony algorithm fine calibration, the initial detection result is used as a constraint to construct a Prony model, the Hankel matrix, SVD decomposition and characteristic equation solving are used to iteratively optimize the harmonic parameters, and the high-order harmonic detection precision is improved; Step S4: adaptive fusion and output, the weights are dynamically allocated according to the signal-to-noise ratio and the number of harmonics, the advantages of the two algorithms are fused, and the final parameters are output.

2. The phase difference method and Prony algorithm fusion terminal and electric energy meter power line harmonic content algorithm according to claim 1, characterized in that, The step S1 is as follows: S11: signal acquisition, the fusion terminal acquires 3-cycle voltage or current signals by controlling the electric energy meter at a sampling frequency of 2560 Hz ; the original voltage or current discrete signals acquired by the electric energy meter are represented by . S12: DC component removal, sliding window average method is used to remove DC component, used to eliminate DC The interference of phase calculation, the sliding window average method takes 64 points per cycle as the sliding window, and the DC component is calculated in real time and removed from the original signal. The calculation method is as follows: When k < 63, the window is supplemented with the previous period history data, the first period direct current component calculation is effective, and The original signal with the direct current component removed; S13: Anti-aliasing FIR filter, suppress high frequency noise through 8 order Hann window FIR filter, FIR filter cutoff frequency Filtering of the filter, passband ripple ≤0.1 dB, stopband attenuation ≥ 60 dB, covering 50th harmonic 2500 Hz, reserving 100 Hz guard band, specific parameters as follows: unit impulse response is: Filtering: based on FPGA parallel convolution operation, convolution operation is synchronized through 8 multipliers, and the filtering formula is: Based on FPGA parallel convolution, the processing delay is less than or equal to 4us; S14: Check, through the SNR signal check, verify the signal effective, get the determined de-noising DC signal To provide clean input signal for subsequent algorithm with signal-to-noise ratio ≥ 30dB.

3. The phase difference method and Prony algorithm fusion terminal and electric energy meter power line harmonic content algorithm according to claim 2, characterized in that, The SNR check step in the step S14 is as follows: S141: initial fitting of the fundamental wave: assuming that the fundamental wave frequency is 50 Hz, the fitting expression is: its take a maximum of 0.707 times, 0.707 times being the effective value estimate, ; S142: power calculation: S143: signal-to-noise ratio calculation: S143: abnormal processing, the abnormal processing mechanism process is: If SNR<30dB, send the resampling instruction; If SNR<30dB is not satisfied for 3 times in succession, trigger the device abnormal alarm, the alarm information is uploaded to the master station through the 4G module, and the event record is stored into the database.

4. The phase difference method and Prony algorithm fusion terminal and electric energy meter power line harmonic content algorithm according to claim 3, characterized in that, The initial detection by phase difference method in the step S2 is based on the linear variation characteristic of the signal phase, the fundamental wave frequency and the initial parameters of harmonics are locked, and the initial constraint is provided for the Prony algorithm, wherein the detection steps are as follows: S21: input signal pre-processing filtered signal ; S22: time domain zero-crossing detection and period calculation, the detection process is as follows: Threshold setting for double threshold zero-crossing detection: positive zero-crossing threshold , negative zero-crossing threshold ; The zero-crossing time calculation is performed by linear interpolation: when detecting the accurate zero-crossing time calculation formula is obtained as follows: wherein is the zero-crossing time; Period calculation: take the adjacent two positive zero-crossing time , , fundamental period is ; The fundamental period T1 is subjected to an effectiveness constraint: If the range is exceeded, the detection window is expanded to 5 periods and the average is taken to ensure the reliability of the period calculation. S23: Hilbert transform and analytic signal construction, the Hilbert transform equation is: Wherein, the Hilbert transform is realized by FFT acceleration to reduce the calculation complexity; Analytic signal and instantaneous phase extraction: where A(k) is the instantaneous amplitude, is the instantaneous phase; The phase unfolding formula is: Fundamental frequency calculation: take the phase difference of adjacent cycles The fundamental frequency is then: The fundamental wave frequency resolution is less than or equal to 0.01 Hz; S24: parameter estimation, based on the deduction of the power grid harmonic characteristics, the initial parameters of harmonics are estimated, the integer multiple frequency is erected, and the erection process is as follows: The frequency calculation of the initial frequency of harmonics is: ； The initial amplitude of the harmonic is calculated based on a large amount of measured data of the State Grid, and the initial amplitude of the harmonic is set in a general scene wherein the initial amplitude of the harmonic is used for the estimation of the effective value of the fundamental wave; S25: result output, the initial detection parameters are output to the Prony module.

5. The algorithm for calculating the harmonic content of the power line of the fusion terminal and the electric energy meter based on the phase difference method and the Prony algorithm according to claim 1, characterized in that, The step S3 comprises the following steps: S31: Prony model construction, establish a model by constraining the initial detection parameters, input the filtered signal and phase difference method initial detection parameters , , construct a P=50 order Prony model, the formula of the established model is: where p is the model order, is the complex amplitude, is the fitting residual, sampling period ; S32: Hankel matrix and SVD decomposition, the Hankel matrix is constructed by selecting the first 2p points of the filtered signal, and the matrix is as follows: Wherein, the singular value SVD decomposition is the characteristic equation solving of the Prony model, the characteristic equation solving is realized by calling the dgesvd function of the LAPACK library to realize the matrix decomposition, separate the effective signal and the noise, and the decomposition formula is: wherein is a diagonal matrix, retaining the first p nonzero singular values and rejecting the small singular values corresponding to noise; S33: characteristic equation solving and parameter iteration, used for extracting accurate harmonic parameters; S34: parameter accuracy check, the amplitude error of the 1-13th harmonic is less than or equal to 3%, and the error of the 21-50th harmonic is less than or equal to 12%.

6. The phase difference method and Prony algorithm fusion terminal and electric energy meter power line harmonic content algorithm according to claim 5, characterized in that, The step S33 of extracting the accurate harmonic parameter includes the following steps: S331: Perform coefficient vector extraction from the last column of the matrix to construct the characteristic equation ; S332: Eigenvalue solving, get eigenvalue by polynomial root finding algorithm is: where there is no attenuation, so ; S333: Parameter update, frequency, complex amplitude and amplitude and phase calculation, wherein the frequency The calculation formula is: , sampling period ; Complex amplitudes are used to construct a Vandermonde matrix , which is solved by least squares ; amplitude is , phase is ; S334: iteratively optimize the harmonic parameters, with the iteration termination condition being or the iteration number is ≥ 30.

7. The algorithm for calculating the harmonic content of the power line of the fusion terminal and the electric energy meter based on the phase difference method and the Prony algorithm according to claim 1, characterized in that, The step S4 includes the following steps: S41: input initial detection parameters of phase difference method and and Prony fine calibration parameters , and ; S42: weight calculation, calculate the signal-to-noise ratio weight factor, SNR≥20dB is effective, the signal-to-noise ratio weight factor calculation formula is: Wherein the higher the SNR, the stronger the reliability of the phase difference method, and the greater the weight; when SNR=40dB =0.67, when SNR=25dB ; Calculate the harmonic number weight factor, the formula is: The larger n is, the weaker the harmonic signal is, the more obvious the Prony algorithm noise resistance advantage is, and the greater the weight is; The maximum weight is the phase difference method weight and the Prony algorithm weight, the phase difference method weight calculation formula is: The Prony algorithm weight calculation formula is: S43: weighted fusion calculation, through the connection with the weight calculation, the output of the optimal parameter is realized, and the calculation mode is: wherein, initially set to 0, is gradually corrected by iteration, and the accuracy of the fused phase is ≤ ± 0.5° S44: Post-processing parameter output, harmonic content The calculation is: wherein the is treated specially when , is reported as 0.00% and labeled not detected.

8. The algorithm for determining the harmonic content of the power line of the electricity meter based on the fusion of the phase difference method and the Prony algorithm according to claim 7, characterized in that, The total harmonic distortion in the S44 is also included The calculation method is: The total harmonic distortion rate and the harmonic content rate are kept to 2 decimal places; The total harmonic distortion and harmonic content data will be uploaded through the 4G module at a regular interval of 1 minute in the prescribed format. The total harmonic distortion Real-time upload beyond the limit, for data security.

9. A fusion terminal for implementing the algorithm of any of claims 1-8, characterized by It includes a main controller, a memory, a 2500V isolation RS485 interface, an ME909s-8214G module and a Beidou / GPS dual-mode timing module, and the working temperature range is-40℃~85℃.

10. A smart meter for use with the converged terminal of claim 9, characterized in that, STM32F407 microprocessor, AD7606-16 ADC chip, 1000:1 voltage division network and Rogowski coil current sensor are adopted, the electric energy meter integrates DS18B20 temperature compensation, and supports 2560Hz synchronous sampling.

Images