A method and device for detecting faults in a direct current current transformer with varying magnetic flux density

Abstract

This invention proposes a method, device, and equipment for detecting faults in DC current transformers based on changes in magnetic flux density. It obtains a family of magnetic field recovery curves through gradual pulse excitation and extracts characteristic spectra reflecting the fatigue level of the iron core. Time-scale decomposition technology is used to separate fast and slow recovery components, calculate the scale coupling factor, and establish a dynamic health benchmark. A bidirectional approximation test is employed to locate the hysteresis critical point, and fault precursors are identified from the hysteresis width parameter. Intersection points are extracted through curve family cross-analysis, an angle field is constructed, and the angle sequence is calculated to quantify the fault development speed. An adaptive excitation scheme based on the development speed is designed to identify fault resonance modes from chaotic transition characteristics. Phase space trajectory deviation analysis is used to generate a multi-resolution spatial fault distribution map. A comprehensive fault index is calculated by integrating spatial distribution and temporal evolution information to achieve four-level graded diagnosis, forming a complete fault detection technology solution for current transformers from early warning to precise location.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of power equipment fault detection, in particular to a DC current transformer fault detection method based on magnetic flux density variation, a device and equipment. BACKGROUND

[0002] As a key measuring device in DC power transmission system, the operation reliability of DC current transformer directly affects the safety and stability of power grid. With the increase of operation time, the magnetic properties of the transformer core material will degrade, the magnetic hysteresis characteristics will change, etc. The traditional periodic maintenance and simple electrical test cannot timely detect these progressive faults, and the performance of the device will be severely degraded or even fail, which will bring safety hazards to the operation of power grid.

[0003] The existing transformer detection technology mainly focuses on the measurement of static electrical parameters such as insulation resistance, DC resistance, and transformation ratio error. These parameters can only reflect the instantaneous state of the device, and cannot reveal the dynamic change rule of the core magnetic properties. The current fault judgment method relies too much on the absolute value comparison of a single index, lacks analysis of the internal relationship between multiple related parameters, and is easy to miss complex fault modes. In addition, the traditional detection needs to be carried out after the device is shut down, and cannot realize continuous monitoring and trend analysis. More importantly, the existing method can only give a general conclusion of "qualified" or "unqualified", which cannot accurately evaluate the severity of the fault, nor can it determine the specific location of the fault in the device, which brings difficulties to the decision of maintenance personnel.

[0004] Therefore, there is an urgent need for a method to solve at least one of the above problems. SUMMARY

[0005] The present application provides a DC current transformer fault detection method based on magnetic flux density variation, device and equipment, which realizes early detection and accurate evaluation of faults by analyzing the dynamic recovery characteristics of the magnetic field, establishes a complete detection chain from the microscopic magnetization process to the macroscopic fault diagnosis, generates a comprehensive diagnostic report containing fault type, severity, spatial location and development trend, and realizes accurate state evaluation of DC current transformer.

[0006] The first aspect of the present application provides a DC current transformer fault detection method based on magnetic flux density variation, comprising the following steps:

[0007] A gradual pulse excitation sequence is applied to the DC current transformer to collect magnetic flux density recovery process data, and a time sequence recovery curve family is generated according to the magnetic flux density recovery process data, and a magnetic field fatigue feature spectrum is extracted from the time sequence recovery curve family;

[0008] The magnetic field fatigue characteristic spectrum is used for time scale decomposition to obtain a fast recovery component and a slow recovery component, ratio analysis is performed on the fast recovery component and the slow recovery component to generate a scale coupling factor, and a magnetic field health state benchmark is constructed through the scale coupling factor;

[0009] Bidirectional approximation testing is performed according to the magnetic field health state benchmark to obtain an ascending critical point and a descending critical point, a magnetic hysteresis width parameter is determined based on the ascending critical point and the descending critical point, and a fault precursor feature is extracted from the magnetic hysteresis width parameter;

[0010] Curve family intersection analysis is performed on the fault precursor feature to obtain a set of intersection positions, an intersection angle sequence is extracted through the set of intersection positions, and a fault development speed index is generated according to the intersection angle sequence;

[0011] A parameter coupling excitation scheme is designed according to the fault development speed index, the parameter coupling excitation scheme is used to excite the transformer to obtain a coupling response spectrum, and a fault characteristic resonance mode is extracted from the coupling response spectrum;

[0012] The fault characteristic resonance mode and the scale coupling factor are analyzed to obtain an adaptive window parameter, a multi-resolution deviation degree is obtained through the adaptive window parameter, and a spatial fault distribution map is generated according to the multi-resolution deviation degree;

[0013] Fusion analysis is performed on the spatial fault distribution map and the fault precursor feature to obtain a comprehensive fault index, a hierarchical diagnosis result is generated according to the comprehensive fault index, and direct current current transformer fault detection is completed.

[0014] The second aspect of the present application proposes a direct current current transformer fault detection device with varying magnetic flux density, comprising:

[0015] A feature extraction module is configured to apply a gradual pulse excitation sequence to the direct current current transformer to collect magnetic flux density recovery process data, generate a time sequence recovery curve family according to the magnetic flux density recovery process data, and extract a magnetic field fatigue characteristic spectrum from the time sequence recovery curve family;

[0016] A scale analysis module is configured to use the magnetic field fatigue characteristic spectrum for time scale decomposition to obtain a fast recovery component and a slow recovery component, perform ratio analysis on the fast recovery component and the slow recovery component to generate a scale coupling factor, and construct a magnetic field health state benchmark through the scale coupling factor;

[0017] A bidirectional testing module is configured to perform bidirectional approximation testing according to the magnetic field health state benchmark to obtain an ascending critical point and a descending critical point, determine a magnetic hysteresis width parameter based on the ascending critical point and the descending critical point, and extract a fault precursor feature from the magnetic hysteresis width parameter;

[0018] The cross analysis module is used for obtaining a set of intersection positions through curve family cross analysis for the fault precursor feature, extracting a cross angle sequence through the set of intersection positions, and generating a fault development speed index according to the cross angle sequence.

[0019] The coupling excitation module is used for designing a parameter coupling excitation scheme according to the fault development speed index, exciting the mutual inductor to obtain a coupling response spectrum through the parameter coupling excitation scheme, and extracting a fault feature resonance mode from the coupling response spectrum.

[0020] The adaptive analysis module is used for analyzing the fault feature resonance mode and the scale coupling factor to obtain an adaptive window parameter, obtaining a multi-resolution deviation degree through the adaptive window parameter, and generating a spatial fault distribution map according to the multi-resolution deviation degree.

[0021] The comprehensive diagnosis module is used for performing fusion analysis on the spatial fault distribution map and the fault precursor feature to obtain a comprehensive fault index, generating a hierarchical diagnosis result according to the comprehensive fault index, and completing the fault detection of the DC current mutual inductor.

[0022] The third aspect of the present application provides a computer device, comprising a memory, a processor and a computer program stored in the memory and executable on the processor, wherein the processor implements the steps of the magnetic flux density change DC current mutual inductor fault detection method disclosed in the first aspect when executing the program.

[0023] The beneficial effects of this invention are reflected in the following aspects: First, by using gradual pulse excitation and time-scale decomposition techniques, the changing characteristics of rapid and slow magnetization processes can be captured separately. The scale coupling factor quantifies the interaction strength of processes at different time scales, and deviations from the healthy baseline can detect early anomalies. Bidirectional approximation testing combined with micro-perturbation analysis accurately locates the nonlinear transition point of the magnetic field response, and the hysteresis width growth rate, as a sensitive fault precursor feature, can issue early warnings in the early stages of a fault. Second, cross-analysis of fault precursor feature curves generates an angle sequence, achieving continuous characterization of fault evolution; the angle change rate directly reflects the speed of fault development. The adaptive excitation scheme dynamically adjusts test parameters according to the fault development speed, and temperature-current composite excitation fully excites potential defects. Chaotic transition characteristics reflect changes in the dynamic characteristics of the magnetization process, and resonance frequency shift and quality factor changes provide criteria for fault type identification. Finally, phase space trajectory deviation analysis quantifies the degree of local anomaly in the magnetic field, and the spatial fault distribution map generated by multi-resolution processing intuitively displays the fault location and impact range. A comprehensive fault index integrates spatial and temporal multi-dimensional information, achieving a quantitative assessment of the fault severity. The tiered diagnostic results directly correspond to differentiated maintenance measures ranging from routine inspections to immediate shutdown, transforming complex analysis results into clear maintenance plans.

[0024] It should be understood that the above general description and the following detailed description are exemplary and explanatory only, and do not limit this application. Attached Figure Description

[0025] The accompanying drawings illustrate specific examples of the technical solutions described in this invention and, together with the detailed embodiments, form part of the specification, serving to explain the technical solutions, principles, and effects of this invention.

[0026] Unless otherwise specified or defined, the same reference numerals in different figures represent the same or similar technical features, and different reference numerals may be used to represent the same or similar technical features.

[0027] Figure 1 This is a flowchart illustrating a method for detecting faults in a DC current transformer based on changes in magnetic flux density, according to the present invention.

[0028] Figure 2 This is a structural block diagram of a DC current transformer fault detection device based on magnetic flux density variation according to the present invention.

[0029] Figure 3 This is a schematic diagram of the structure of a computer device according to the present invention. Detailed Implementation

[0030] In the following description, specific details such as particular system architectures and techniques are set forth for illustrative purposes and not for limitation, in order to provide a thorough understanding of the embodiments of this application. However, those skilled in the art will understand that this application may also be implemented in other embodiments without these specific details. In other instances, detailed descriptions of well-known systems, apparatuses, circuits, and methods have been omitted so as not to obscure the description of this application with unnecessary detail.

[0031] It should be understood that, when used in this application specification and the appended claims, the term "comprising" indicates the presence of the described features, integrals, steps, operations, elements and / or components, but does not exclude the presence or addition of one or more other features, integrals, steps, operations, elements, components and / or a collection thereof.

[0032] It should also be understood that the term “and / or” as used in this application specification and the appended claims means any combination of one or more of the associated listed items and all possible combinations, and includes such combinations.

[0033] As used in this application specification and the appended claims, the term "if" may be interpreted, depending on the context, as "when," "once," "in response to determination," or "in response to detection." Similarly, the phrase "if determined" or "if detected [the described condition or event]" may be interpreted, depending on the context, as meaning "once determined," "in response to determination," "once detected [the described condition or event]," or "in response to detection [the described condition or event]."

[0034] Furthermore, in the description of this application and the appended claims, the terms "first," "second," "third," etc., are used only to distinguish descriptions and should not be construed as indicating or implying relative importance.

[0035] References to "one embodiment" or "some embodiments" as described in this specification mean that one or more embodiments of this application include a specific feature, structure, or characteristic described in connection with that embodiment. Therefore, the phrases "in one embodiment," "in some embodiments," "in other embodiments," "in still other embodiments," etc., appearing in different parts of this specification do not necessarily refer to the same embodiment, but rather mean "one or more, but not all, embodiments," unless otherwise specifically emphasized. The terms "comprising," "including," "having," and variations thereof mean "including but not limited to," unless otherwise specifically emphasized.

[0036] The technical solutions of the embodiments of this application will be described below.

[0037] like Figure 1As shown, this embodiment of the invention provides a method for detecting faults in a DC current transformer based on changes in magnetic flux density, including the following steps S110-S170:

[0038] Step S110: Apply a gradual pulse excitation sequence to the DC current transformer to collect magnetic flux density recovery process data, generate a family of time-series recovery curves based on the magnetic flux density recovery process data, and extract the magnetic field fatigue characteristic spectrum from the family of time-series recovery curves.

[0039] Specifically, a gradually increasing pulse excitation sequence was applied to the excitation coil of the DC current transformer. The excitation sequence was designed as a square wave pulse train with linearly increasing amplitude. Pulse parameters were set as follows: initial amplitude 0.1A, final amplitude 5A, step increment 0.1A, pulse width 100ms, pulse interval 500ms, forming 50 increasing pulses. Magnetic flux density measurement employed a Hall effect sensor array, with eight measuring points arranged at different locations on the iron core, forming a spatially distributed measurement. The sensor sampling rate was set to 10kHz to ensure the capture of rapidly changing magnetic field transients. After each pulse was applied, the complete recovery process of the magnetic flux density from peak value to steady state was recorded, with a recovery time window of 400ms. Data acquisition used a 24-bit high-precision ADC with a range of ±2T and a resolution of 0.1mT. Temperature compensation was achieved through real-time monitoring, with each Hall sensor equipped with a temperature probe. Precise control of the excitation current was achieved through a constant current source, with current stability better than 0.1%. A three-layer shielding structure was used for magnetic shielding to reduce environmental magnetic field interference, ultimately obtaining complete magnetic flux density recovery process data.

[0040] The collected magnetic flux density recovery process data were grouped according to excitation amplitude to generate a family of time-series recovery curves. Each recovery curve represents the time evolution of magnetic flux density under a specific excitation amplitude, denoted as B(t,I), where B is the magnetic flux density, t is time, and I is the excitation current amplitude. The curve family contains 50 independent curves, covering an excitation range from 0.1A to 5A. Data preprocessing employed a 5-point moving average filter to remove high-frequency noise while preserving the dynamic characteristics of the recovery process. The steady-state magnetic flux density B_ss(I) was determined by arithmetically averaging the data points at the end of each recovery curve, serving as the steady-state reference value for the excitation current I. Normalization standardized the peak value of each curve to 1, facilitating comparison of recovery characteristics under different excitation levels. Time axis alignment was performed with the pulse end time as the zero point to ensure a consistent time reference for all curves. The characteristic time constant of the recovery process was extracted through exponential fitting, with typical values ​​varying with excitation intensity within the range of 10-50ms. The 3D visualization of the curve family used a waterfall plot to visually demonstrate the evolution of magnetic flux density with time and excitation intensity. Abnormal curves are detected by calculating the similarity between adjacent curves; curves with a similarity lower than 0.9 are marked and removed. Curve interpolation uses a cubic spline method to convert discrete sampling points into continuous functions.

[0041] The magnetic field fatigue characteristic spectrum is extracted from the family of time-series recovery curves to reflect the magnetization fatigue degree of the transformer core. Feature extraction first calculates the integral area A(I) = ∫[B(t,I)-B_ss(I)]dt for each recovery curve, where A(I) is the integral area corresponding to the excitation current I, B(t,I) is the time-series recovery curve, and B_ss(I) is the measured steady-state magnetic flux density under that excitation current, with an integration interval of 0 to 400 ms. The recovery rate characteristic is extracted by calculating the half-recovery time T_50, defined as the time required for the magnetic flux density to decrease from its peak to 50%. The nonlinear exponent is obtained by fitting the initial segment of the recovery curve, using the constrained fitting model B(t) = B_0 × exp(-t / τ) + B_ss(I), where B_0 is the initial magnetic flux density amplitude, τ is the time constant, and B_ss(I) is the previously determined measured steady-state value as the fitting constraint. Frequency domain analysis of fatigue characteristics employed short-time Fourier transform with a window length of 20 ms and an overlap rate of 50%, extracting the spectral evolution within the 0-500 Hz range. Hysteresis was quantified by comparing the recovery characteristics of the rising and falling edges, using the hysteresis exponent H = |τ_up - τ_down| / τ_avg, where τ_up is the time constant of the excitation rising edge, τ_down is the time constant of the excitation falling edge, and τ_avg is the average of the two. Higher-order features included details such as the curvature change of the recovery curve, inflection point position, and oscillation decay rate. Through multi-dimensional feature extraction and analysis, the magnetic field fatigue characteristic spectrum characterizing the magnetic field fatigue state was finally obtained.

[0042] Step S120: Use the magnetic field fatigue characteristic spectrum to perform time-scale decomposition to obtain the fast recovery component and the slow recovery component, perform ratio analysis on the fast recovery component and the slow recovery component to generate a scale coupling factor, and construct a magnetic field health status benchmark through the scale coupling factor.

[0043] Specifically, the recovery process at different time scales is separated using the magnetic field fatigue characteristic spectrum. First, the characteristic parameters A(I), T_50(I), and τ(I) are rearranged into time series A(t), T_50(t), and τ(t) according to the excitation current sequence, where time t corresponds to the time sequence of pulse application. The time scale decomposition employs the Empirical Mode Decomposition (EMD) method, decomposing each characteristic time series into multiple intrinsic mode functions (IMFs). The decomposition process involves finding local extrema, constructing upper and lower envelopes, and iteratively extracting the mean until the IMF condition is met. Taking the time constant sequence τ(t) as the main decomposition object, the first three IMF components correspond to rapidly changing components, reflecting the rapid changes between adjacent pulses; subsequent IMF components correspond to slowly changing components, reflecting the slow evolution of the cumulative effect. The fast recovery component R_fast = Σ(IMF_i), i = 1, 2, 3, where IMF_i is the i-th intrinsic mode function, reflecting the rapid domain flipping process. The slow recovery component R_slow = Σ(IMF_j), j = 4, 5, ..., n, where n is the total number of IMFs, reflecting the slow rearrangement process of magnetic domains. The decomposition residue represents the trend component, which is related to the permanent magnetization of the iron core. The time scale is defined by calculating the average period of each IMF. The energy proportion of the components varies between 0.3 and 0.7, reflecting the relative contribution of the fast and slow processes. EMD decomposition achieves effective separation of the magnetic field recovery process on the time scale, obtaining the fast recovery component R_fast and the slow recovery component R_slow.

[0044] In some embodiments, the step of performing ratio analysis on the fast recovery component and the slow recovery component to generate a scale coupling factor includes: performing a logarithmic transformation on the fast recovery component to obtain nonlinear decay characteristics; projecting the slow recovery component onto phase space to extract trajectory density; generating an implicit coupling mode by cross-mapping the nonlinear decay characteristics and the trajectory density; and deriving the scale coupling factor using the implicit coupling mode.

[0045] A logarithmic transformation is performed on the fast recovery component to obtain nonlinear decay characteristics. The reference amplitude R_ref is taken as the root mean square value of the fast recovery component R_fast to ensure consistency of normalization. The logarithmic transformation formula is L_fast = ln(1 + |R_fast| / R_ref), where L_fast is the logarithmically transformed signal and R_fast is the fast recovery component. The transformed signal exhibits approximately linear decay characteristics, facilitating the extraction of the decay rate. Piecewise linear fitting divides the decay process into an initial fast segment (0-5ms) and a subsequent slow segment (5-20ms), calculating the slopes k_1 and k_2 respectively. The nonlinear exponent is defined as N_idx = |k_1 / k_2-1|, where k_1 is the slope of the initial segment and k_2 is the slope of the subsequent segment, reflecting the degree of nonlinearity of the decay process. The decay time constant is extracted through exponential fitting, and the applicability of the model is evaluated by the fitting residuals. The abrupt change point is determined by the zero point of the second derivative; the distribution pattern of the abrupt change point is related to the step characteristics of the magnetic domain wall motion. Spectral analysis of the decay energy revealed the dominant decay mode: low-frequency components correspond to the slow adjustment of large magnetic domains, while high-frequency components reflect the rapid response of small magnetic domains. The temperature dependence of the nonlinear decay characteristics was obtained through multi-temperature point testing, with a temperature coefficient of approximately -0.5% / ℃. The eigenvector [N_idx,k_1,k_2,τ_f] fully describes the nonlinear decay characteristics of the fast component.

[0046] The slow recovery component is projected onto the phase space to extract the trajectory density. Phase space reconstruction employs the delayed coordinate method, constructing an m-dimensional state vector X(t) = [R_slow(t), R_slow(t-τ), ..., R_slow(t-(m-1)τ)]. The embedding dimension m = 5 is determined using the spurious nearest neighbor method; as m increases to 5, the proportion of spurious nearest neighbors drops below 1%. The delay time τ = 15 ms is determined using the mutual information method, corresponding to the delay at which the mutual information first reaches its minimum value. The phase space trajectory forms a specific geometric structure in 5-dimensional space, exhibiting a quasi-periodic torus under normal conditions and strange attractors under abnormal conditions. The trajectory density is calculated using kernel density estimation, and the bandwidth selection is optimized using cross-validation. The peak positions of the density distribution correspond to the attractor centers, and the number of peaks reflects multi-stable characteristics; normal current transformers typically have 1-2 stable attractors. The fractal dimension of the trajectory is calculated using the correlation dimension, with typical values ​​between 2.3 and 2.8, reflecting the system complexity. The Lyapunov exponent assesses the divergence velocity of the trajectory; a maximum Lyapunov exponent λ_max > 0 indicates chaotic behavior. The residence time distribution of the trajectory in phase space follows a power law, with an exponent of approximately -1.5. The curl calculation of the density gradient reveals the vortex structure of the trajectory, and the maximum trajectory density ρ_max and its spatial distribution constitute the dynamic characteristics of the slow recovery component.

[0047] Implicit coupling patterns are generated through cross-mapping of nonlinear decay features and trajectory density. The cross-mapping constructs a two-dimensional feature space (N_idx, ρ_max), where N_idx is the nonlinear exponent and ρ_max is the maximum trajectory density. Radial basis function interpolation is used to ensure the continuity and smoothness of the feature space. Cluster analysis divides the mapping space into four coupling pattern regions: strongly coupled (high N_idx, high ρ_max), weakly coupled (low N_idx, low ρ_max), fast-dominant (high N_idx, low ρ_max), and slow-dominant (low N_idx, high ρ_max). Pattern boundaries are obtained through support vector machine training, achieving a classification accuracy of 92%. Pattern transition paths are obtained by connecting temporally adjacent mapping points; the transition speed reflects the stability of the system state. Markov chain analysis calculates the transition probability matrix between patterns, and the steady-state distribution reveals the long-term proportion of each pattern. Outliers are defined as those deviating from all cluster centers by more than 2.5 standard deviations, with an occurrence frequency of less than 5%. The calculation of mode entropy quantifies the degree of order in the coupling state; an increase in entropy indicates that the system tends towards disorder. Mapping analysis identifies four implicit coupling modes and their dynamic evolution patterns, revealing the complex coupling relationship between fast and slow components.

[0048] The scale coupling factor is derived using implicit coupling modes. The coupling factor is calculated as S_c = Σp_i × s_i × exp(-d_i / d_0), where p_i is the probability of the i-th mode, s_i is the mode strength coefficient (1.5 for strong coupling, 0.5 for weak coupling, 0.8 for fast dominance, and 0.7 for slow dominance), d_i is the distance to the mode center, and d_0 is the normalized distance (twice the average distance). The probability p_i is calculated using time proportions, considering the dynamic evolution of the mode. The exponential decay of the distance weight ensures that states closer to the mode center contribute more. The time-varying coupling factor S_c(t) is calculated using a sliding window with a window length of 100 ms and a step size of 10 ms, achieving dynamic tracking of the coupling state. Spectral analysis of the coupling factor reveals a dominant oscillation of 0.5-5 Hz, corresponding to the characteristic frequency of the magnetization process. Statistical distribution analysis shows that S_c approximately follows a log-normal distribution, with the logarithmic mean and variance reflecting the average level and stability of the coupling, respectively. The autocorrelation function of the coupling factor exhibits exponential decay, with a correlation time of approximately 200 ms. An improved 3σ criterion is used for outlier detection, taking into account the asymmetry of the distribution. The gradient rate of change of the coupling factor serves as an early warning indicator; rapid changes often foreshadow impending failures. Quantitative analysis of the implicit coupling mode ultimately generates the dynamic-scale coupling factor S_c(t).

[0049] Construct a magnetic field health state benchmark through the scale coupling factor. The health state benchmark is defined as H_ref = S_c0 × [1 - α × (T - T0)], where S_c0 is the coupling factor in the initial state, taken from the S_c value measured at the initial stage of the transformer operation, α is the aging coefficient, with a typical value of 0.001 / day, T is the operation time, and T0 is the initial moment. This benchmark curve describes the natural degradation process of the iron core material over time through a linear aging model, and the aging coefficient α is determined based on the historical degradation data of the same type of transformers. The benchmark envelope is obtained by statistically analyzing a large number of S_c values in the normal operation state. The upper limit H_upper = H_ref + 2σ_H, and the lower limit H_lower = H_ref - 2σ_H, where σ_H is the standard deviation of the health state, reflecting the reasonable fluctuation range of the coupling factor under normal operation conditions. The health margin M_h = (S_c - H_lower) / (H_upper - H_lower) quantifies the relative distance from the current coupling factor S_c to the abnormal lower limit, where S_c is the scale coupling factor measured in real time, and the value range of M_h is 0 - 1. The state space is divided into three regions according to the health margin M_h: the normal region (M_h > 0.7) indicates that the magnetic field state of the transformer is good and it can operate safely; the warning region (0.3 < M_h ≤ 0.7) indicates a decline in magnetic field performance and requires an increased monitoring frequency; the abnormal region (M_h ≤ 0.3) indicates severe magnetic field deterioration and it is recommended to arrange maintenance immediately. The initial calibration of S_c0, the degradation modeling of α, the statistical determination of σ_H, and the real-time calculation of M_h together constitute the complete magnetic field health state benchmark H_ref.

[0050] Step S130, perform a two-way approximation test according to the magnetic field health state benchmark to obtain the rising critical point and the falling critical point, determine the hysteresis width parameter based on the rising critical point and the falling critical point, and extract the fault precursor characteristics from the hysteresis width parameter.

[0051] Specifically, a bidirectional approximation test is performed based on the magnetic field health state benchmark H_ref to find the critical point where the magnetic field state undergoes a sudden change. The bidirectional approximation test is designed with two test sequences: an increasing sequence starting from 0.1 times the rated current and an increasing step of 0.05 times until a nonlinear change occurs in the magnetic field response; the decreasing sequence starts from the current operating point and decreases step by 0.05 times until the magnetic field recovers a linear response. A step excitation lasting 200ms is applied to each test point, and the dynamic response of the magnetic flux density is recorded and the corresponding coupling factor S_c_test is calculated. The response evaluation index is defined as R_eval = |S_c_test - H_ref| / H_ref, where S_c_test is the simplified coupling factor of the test point, calculated by the ratio of the fast and slow components of the magnetic flux density recovery process: S_c_test = A_fast / A_slow, where A_fast is the integral area of ​​the first 20ms of the recovery process, A_slow is the integral area of ​​the subsequent time, and H_ref is the health state benchmark. When R_eval exceeds 0.1, it is determined to be a deviation from the normal state. The rising critical point I_up is defined as the excitation current value where R_eval > 0.1 first appears in the increasing sequence, and the falling critical point I_down is defined as the excitation current value where R_eval recovers to below 0.1 in the decreasing sequence. Precise location of the critical points employs a bisection method for refined searching, halving the step size and repeating the test until the critical point location accuracy reaches 0.01 times the rated current. Temperature compensation during the test is achieved by real-time monitoring of the core temperature, and all measurements are normalized to a reference temperature of 25℃. The convergence of the bidirectional approximation is judged by comparing the difference in critical points between two adjacent tests; convergence is considered when the difference is less than 1%. The test sequence design considers the influence of magnetization history, and demagnetization is performed before each change in excitation direction. The bidirectional approximation test accurately locates the rising critical point I_up and the falling critical point I_down.

[0052] In some embodiments, determining the hysteresis width parameter based on the rising critical point and the falling critical point includes: triggering a micro-perturbation test at the rising critical point to record the response amplitude; repeating the perturbation test at the falling critical point to obtain the decay amplitude; performing a time-series comparison between the response amplitude and the decay amplitude to generate an amplitude deviation sequence; and generating the hysteresis width parameter based on the extreme value distribution of the amplitude deviation sequence.

[0053] The response amplitude was recorded by triggering a micro-perturbation test at the rising critical point. The micro-perturbation was designed as a small-amplitude sinusoidal signal superimposed on the critical point current I_up, with a perturbation amplitude ΔI = 0.01 × I_up and a frequency range of 1-100Hz. The perturbation response was measured using a high-precision magnetic flux sensor, with a sampling rate set to more than 100 times the perturbation frequency. The response amplitude A_up(f) was defined as the peak-to-peak value of the magnetic flux density response, forming a frequency response curve as it varied with the perturbation frequency. The phase delay was obtained by calculating the phase difference between the excitation and the response, reflecting the dynamic characteristics of the system. Nonlinear indicators were extracted through harmonic analysis to quantify the degree of distortion in the response. Transient response characteristics were obtained by applying a step perturbation, with rise time and overshoot reflecting the dynamic performance of the system. The time-frequency analysis of the response used a short-time Fourier transform with a window length of 10 perturbation periods and an overlap rate of 75%. Energy transfer characteristics quantified the energy transfer efficiency at different frequencies. The probability density distribution of the response amplitude was obtained through long-term statistics, with skewness and kurtosis reflecting the non-Gaussian nature of the distribution. The sensitivity at the critical point was assessed by calculating the derivative of the amplitude with respect to the current; high sensitivity predicts potential instability. Micro-perturbation testing yielded the complete response amplitude characteristic set A_up(f) at the rising critical point.

[0054] The decay amplitude was obtained through repeated perturbation tests at the falling critical point. The micro-perturbation tests at the falling critical point I_down used the same test parameters as those at the rising critical point to ensure comparability of results. The decay amplitude A_down(f) exhibited different frequency characteristics than A_up(f), with more pronounced decay in the low-frequency range. The phase characteristics showed greater hysteresis, reflecting the slow response of the magnetic domains at low excitation levels. The nonlinearity index was generally smaller than at the rising critical point, indicating better linearity at the falling critical point. The transient response at the falling critical point exhibited overdamped characteristics, with no significant overshoot and a longer settling time. Time-frequency analysis revealed mode transitions during the decay process, gradually converging from initial multi-frequency oscillations to a single dominant frequency. The energy transfer function decreased significantly in the low-frequency range, reflecting enhanced energy dissipation. The decay time constant was extracted through exponential fitting, describing the decay rate of the amplitude. The spatial distribution of the decay amplitude was obtained through multi-point measurements, revealing the non-uniformity of magnetic field decay. Noise level analysis showed a reduced signal-to-noise ratio at the falling critical point, requiring more refined signal processing. The complete set of decay amplitude features A_down(f) contrasts with the rising critical point.

[0055] An amplitude deviation sequence is generated by comparing the response amplitude with the decay amplitude over time. The amplitude deviation is defined as ΔA(f) = A_up(f) - A_down(f), where A_up and A_down are the response amplitudes at the rising and falling critical points, respectively. Time-series comparison is achieved by synchronizing two sets of test data, with a time alignment accuracy of 0.1 ms. Moving average processing of the deviation sequence eliminates random fluctuations, with a window length of 5 data points. The relative deviation sequence provides a normalized comparison benchmark, facilitating comparisons under different test conditions. Spectral analysis of the deviation sequence reveals the characteristic frequency at which the deviation reaches its maximum value. Time-varying deviations are obtained through continuous testing, demonstrating the evolution of the deviation over time. The autocorrelation function of the deviation sequence exhibits a periodic structure, with the period related to the magnetization process. Cross-correlation analysis reveals the correlation between the rising and falling processes. The statistical moments of the deviation sequence, including mean, variance, skewness, and kurtosis, comprehensively characterize the distribution properties. Trend analysis extracts the long-term trend of the deviation through linear regression. The amplitude deviation sequence ΔA(f) and its statistical characteristics constitute a quantitative description of hysteresis asymmetry.

[0056] The hysteresis width parameter is generated based on the extreme value distribution of the amplitude deviation sequence. Extreme value extraction is achieved by finding the local maximum and minimum values ​​of the deviation sequence ΔA(f). The hysteresis width parameter is defined as W_h = max(ΔA(f)) - min(ΔA(f)), reflecting the maximum deviation range across the entire frequency band. The relative hysteresis width W_rel = W_h / I_rated × 100%, where I_rated is the rated current, achieving a normalized representation. The frequency-dependent hysteresis width W_h(f) is obtained through piecewise fitting, dividing the frequency range into low-frequency, mid-frequency, and high-frequency segments. The hysteresis asymmetry is obtained by calculating the ratio of positive to negative deviations, reflecting the symmetry of the hysteresis loop. The hysteresis area is estimated by integrating the deviation sequence, characterizing energy loss. The dynamic hysteresis width tracks the change of W_h over time, identifying degradation trends. The temperature correction of the hysteresis width considers the influence of temperature on magnetic properties. The comprehensive hysteresis width index integrates multiple parameters to provide an overall assessment. The set of hysteresis width parameters {W_h,W_rel,W_h(f)} generated by extreme value distribution analysis fully characterizes the hysteresis properties.

[0057] Fault precursor features are extracted from the hysteresis width parameter. The time series W_h(t) of the hysteresis width is obtained through periodic, repeated bidirectional approximation tests. The test cycle is determined based on the importance of the equipment, typically once a month. The primary indicator of the fault precursor features is the hysteresis width growth rate r_W = (W_h(t) - W_h(t - Δt)) / Δt, where W_h(t) is the current hysteresis width, W_h(t - Δt) is the hysteresis width at the previous moment, and Δt is the measurement time interval. When r_W exceeds three times the standard deviation of the normal fluctuation range, it is marked as abnormal growth. Changes in hysteresis asymmetry reflect unidirectional degradation of magnetic properties; significant changes indicate potential faults. Frequency characteristic shifts are achieved by monitoring changes in characteristic frequencies; a decrease in frequency indicates a decrease in core permeability. The relative increase in energy loss density serves as a criterion for energy efficiency degradation. A decrease in the hysteresis loop shape factor indicates magnetic performance degradation. Asymmetric development of harmonic characteristics indicates deterioration of nonlinear characteristics. Time series forecasting is constructed based on historical data to predict the evolution trend of hysteresis width. The fault risk index integrates multiple precursor features to provide a unified risk assessment. The multidimensional feature vector [r_W,W_rel,W_h(f)] constitutes a complete set of fault precursor features, enabling early warning of transformer faults.

[0058] Step S140: Perform curve family cross analysis on the fault precursor features to obtain the set of intersection points, extract the cross angle sequence from the set of intersection points, and generate a fault development speed index based on the cross angle sequence.

[0059] Specifically, a curve family cross-analysis is performed on the fault precursor feature set [r_W, W_rel, W_h(f)] to find the intersection points of different feature curves. The curve family is constructed by plotting the evolution of each precursor feature over time as a continuous curve, with the time span being the measurement data from the most recent 30 days. The W_h(t) time series is obtained based on the aforementioned periodic bidirectional approximation test. r_W(t) is calculated by the difference between W_h values ​​at adjacent times, and W_rel(t) is obtained by normalizing W_h(t) relative to the rated current. The hysteresis width growth rate curve r_W(t) is smoothed using a 5-point method to eliminate the influence of measurement noise. The relative hysteresis width curve W_rel(t) is normalized to the 0-100% range for easy comparison with other features. The frequency-dependent hysteresis width W_h(f,t) forms a three-dimensional surface, which is then used to generate a two-dimensional curve family through contour projection. The cross-analysis overlays all feature curves in the same coordinate system, and coordinate axis normalization ensures the comparability of features with different dimensions. Intersection detection is achieved by solving F_i(t) - F_j(t) = 0, where F_i and F_j are any two characteristic curves. The numerical solution uses the Newton-Raphson iterative method with a convergence accuracy of 10^-6. Intersection filtering removes false intersections and tangent points at the endpoints, retaining only the true crossing intersections. The time coordinate t_cross and the eigenvalue F_cross of the intersection record the time and state of the intersection. Intersection density analysis counts the number of intersections per unit time; increased density indicates rapid changes in the system state. The complete set of intersection locations P_cross = {(t_1,F_1),(t_2,F_2),...,(t_n,F_n)} contains the spatiotemporal information of all valid intersections.

[0060] In some embodiments, extracting the cross angle sequence from the set of intersection locations includes: applying a rotation transformation to the set of intersection locations to construct an angle field; searching for gradient abrupt change lines in the angle field to identify angle transition paths; integrating along the angle transition path to obtain cumulative angle changes; and arranging the cumulative angle changes in chronological order to form a cross angle sequence.

[0061] An angle field is constructed by applying a rotation transformation to the set of intersection points. This transformation converts the intersection points from a time-eigenvalue coordinate system (t, F) to a polar coordinate system (r, φ), where r is the polar radius and φ is the polar angle. The origin is chosen at the center of the data time window to ensure the symmetry of the transformation. The angle field A(x, y) is constructed using the polar angle φ as the field value through spatial interpolation of the φ values ​​at discrete intersection points, employing a radial basis function. The angle field is smoothed using kernel density estimation, with a Gaussian kernel and bandwidth optimized through cross-validation. Weight coefficients are determined based on the rate of change of eigenvalues ​​at the intersection points. The rate of change is calculated using the slope of a linear fit of the preceding and following five data points; a larger rate of change results in a higher weight. The spatial resolution of the angle field is set to 0.1 × 0.1, and the coverage area is adaptively adjusted according to the intersection point distribution to ensure all intersection points are within the field. The field strength gradient is calculated using the finite difference method, indicating the direction and intensity of the angle field's change. Contour lines are plotted using the marching squares algorithm, with the contour line intervals automatically determined based on the field value range. The curl calculation of the field identifies the vortex structure. The vortex center corresponds to the key point of state transition, and the vortex intensity reflects the severity of the transition. Zero-gradient boundary conditions are used for boundary treatment to avoid edge effects. The construction of the angle field A(x,y) realizes the mapping from discrete intersection points to a continuous field, providing a continuous search space for subsequent path search.

[0062] Gradient abrupt change lines are searched within the angle field to identify angle transition paths. A gradient abrupt change line is defined as the trajectory where the gradient magnitude changes abruptly, reflecting regions of rapid change in the angle field. Abrupt change detection is achieved by calculating the second derivative of the gradient, using the Sobel operator to improve edge detection accuracy. Abrupt change points are marked when the absolute value of the second derivative exceeds an adaptive threshold. The threshold is determined based on the statistical characteristics of the gradient distribution, taking 1.5 times the 75th percentile. Consecutive abrupt change points are connected by 8-neighborhoods to form abrupt change line. The connection process considers the consistency of the gradient direction; points with a directional deviation exceeding 45 degrees are not connected. Breakpoint handling employs morphological closing operations, first dilating and then eroding, connecting breakpoints with a distance of less than 3 pixels. The direction of the abrupt change line is determined by the normal vector of the local gradient, ensuring it is perpendicular to the direction of the largest gradient change. Path search starts from the position of maximum field strength, using a greedy algorithm to move along the direction of the most drastic gradient change, selecting the direction that maximizes the increase in path cost at each step. Search termination conditions include reaching the field boundary, path self-intersection, or gradient change falling below the threshold. Multiple paths are obtained by starting searches at the top 10 peak positions of the field strength. The curvature of the path was calculated using a three-point circular arc fitting, with high curvature points (curvature radius less than 5 pixels) marked as locations of rapid state transitions. Path smoothing was performed using B-spline fitting, with the number of control points being 1 / 3 of the original number. The identified set of angular transition paths describes the main change channels in the angular field.

[0063] The cumulative angle change is obtained by integrating along the angle transition path. Path integration calculates the sum of angle field values ​​along the path, which physically represents the total angle change along the path. Numerical integration employs an adaptive Simpson's law, adjusting the integration step size according to local curvature, decreasing the step size in regions of high curvature. The integration step size varies between 0.01 and 0.1, with the error controlled within 0.1%. Piecewise integration divides each path into several segments based on curvature changes, with curvature changes within each segment not exceeding 20%. The accumulation process records the cumulative value at each integration point, forming a cumulative curve of angle versus path arc length. Local extrema on the path are detected by changes in the sign of the first derivative; maxima correspond to regions of rapid angle increase, while minima correspond to regions of slowing angle decrease. The integration directionality is determined by comparing the results of forward (from start to end) and reverse (from end to start) integration, selecting the direction with the larger absolute value of the cumulative amount as the principal direction. Path weights are determined based on path length and average field strength, with longer paths and higher field strength paths assigned higher weights. Multipath integral results are combined using a weighted average, with weight normalization ensuring a sum of 1. The physical dimension of the integral value is angle, converted to degrees using radians. Outlier integral values ​​(exceeding 360 degrees) are handled using modular arithmetic. The cumulative angle change for each path and its distribution along the path constitute a spatial description of the angle evolution.

[0064] The cumulative angle changes are arranged chronologically to form a cross-angle sequence. Time mapping is achieved through inverse coordinate transformation of each point on the path, mapping spatial coordinates (x, y) back to the original time coordinate t. The inverse transformation uses nearest neighbor interpolation to find the closest time marker in the original intersection set. The cumulative angles of multiple paths are combined at the same time by weighted averaging, with the weights determined by the field strength of the path at that time. Missing values ​​in the time series are filled using PCHIP (piecewise cubic Hermite interpolation) to ensure monotonicity and continuity. Resampling unifies irregular time intervals to 1 hour, and anti-aliasing filtering is used to prevent high-frequency distortion. The rate of angle change is calculated using a Savitzky-Golay differential filter, with smoothing and differentiation performed simultaneously. The trend component of the sequence is extracted using the STL (Seasonal Trend Decomposition) method, separating the long-term trend, periodic components, and random components. Outlier detection uses the Local Outlier Factor (LOF) algorithm, considering the local density of the time series. Phase expansion is achieved through cumulative phase difference, eliminating 2π jumps and ensuring the continuity of the angle sequence. Timestamps are converted to relative times, with the first valid data point as zero, in hours. Quality labels are assigned a confidence level to each data point, with low-quality data being downweighted in subsequent analyses. The final cross-angle sequence Θ(t) is arranged in ascending order of time, containing the complete time history of angle evolution.

[0065] The fault development speed index is generated based on the cross-angle sequence. The basic definition of fault development speed is V_f = dΘ / dt, where Θ is the cross-angle sequence and t is time, in degrees per hour. Instantaneous speed is calculated using a five-point center difference method to improve the accuracy and stability of numerical differentiation. Acceleration a_f = dV_f / dt reflects the acceleration or deceleration trend of fault development, with a positive value indicating accelerated fault deterioration. Statistical characteristics of the speed include average speed V_avg, maximum speed V_max, and speed standard deviation σ_v, comprehensively characterizing the speed distribution. Acceleration period identification is achieved by detecting continuous periods where a_f > 0, with the shortest duration of a single acceleration period set at 3 hours. For example, after 10 years of operation, the cross-angle sequence of a 500kV instrument transformer shows that V_avg increased from the initial 0.5 degrees / hour to 2.8 degrees / hour, and V_max reached 5.2 degrees / hour, exhibiting three distinct acceleration periods, each lasting 12-24 hours, indicating rapid deterioration of magnetic performance. Velocity spectrum analysis revealed periodic components, mainly concentrated in the daily and weekly frequency bands. Velocity mutation detection identified abnormal acceleration moments where |ΔV_f|>3σ_v, and the increased mutation frequency is a key indicator of impending fault. The fault development speed index comprehensively considers multiple velocity characteristics: V_index = w_1×V_avg + w_2×V_max + w_3×σ_v + w_4×N_acc, where the weighting coefficients w_1 = 0.3, w_2 = 0.3, w_3 = 0.2, w_4 = 0.2, and N_acc represents the number of acceleration periods. This index compresses multi-dimensional velocity information into a single value; a V_index > 3.0 triggers a level one warning, and a V_index > 5.0 triggers emergency maintenance. Through velocity analysis of cross-angle sequences and multi-feature fusion, the fault development speed index V_index, which quantifies the dynamic characteristics of fault evolution, was obtained.

[0066] Step S150: Design a parameter coupling excitation scheme according to the fault development speed index, use the parameter coupling excitation scheme to excite the current transformer to obtain the coupling response spectrum, and extract the fault characteristic resonance mode from the coupling response spectrum.

[0067] Specifically, a parameter-coupled excitation scheme is designed based on the fault development speed index V_index to achieve targeted excitation testing of the current transformer. The basic principle of the excitation scheme is to adjust the excitation intensity and frequency range according to the value of V_index. The larger the V_index, the faster the fault develops, requiring a wider excitation spectrum to capture potential resonance characteristics. The excitation parameter mapping relationship is as follows: when V_index < 3, standard excitation is used, with a current amplitude of 0.5-1.0 times the rated value and a frequency range of 10-100Hz; when 3 ≤ V_index < 5, enhanced excitation is used, with a current amplitude of 0.8-1.5 times the rated value and a frequency range of 5-200Hz; when V_index ≥ 5, strong excitation is used, with a current amplitude of 1.0-2.0 times the rated value and a frequency range of 1-500Hz. Temperature excitation is achieved through a controllable heating device, and the temperature variation range is determined according to V_index: ΔT = 10 + 5 × V_index (°C), ensuring that temperature-related fault characteristics are fully excited within a safe range. The excitation timing design involves first applying a steady-state excitation to establish a baseline, then superimposing a sweep excitation to find the resonance point, and finally applying a pulse excitation to test the transient response. The sweep rate is inversely proportional to V_index, r_sweep = 10 / V_index (Hz / s). When the fault develops rapidly, the sweep rate is reduced to improve resolution. Coupling modes include synchronous coupling, anti-phase coupling, and quadrature coupling. The excitation duration T_exc = 60 + 20 × V_index (seconds) ensures the capture of complete response characteristics. The parameter-coupled excitation scheme P_scheme = {I_range, f_range, ΔT, r_sweep, T_exc} forms the complete test scheme.

[0068] In some embodiments, the step of using the parameter coupling excitation scheme to excite the current transformer and obtain the coupling response spectrum includes: generating a temperature-current composite excitation signal based on the parameter coupling excitation scheme; driving the current transformer to generate a nonlinear response based on the temperature-current composite excitation signal; capturing chaotic transition characteristics based on the nonlinear response; and performing spectral analysis based on the chaotic transition characteristics to obtain the coupling response spectrum.

[0069] A temperature-current composite excitation signal is generated based on a parameter-coupled excitation scheme. The composite excitation signal design comprehensively utilizes various parameters in P_scheme. Current excitation employs a multi-frequency superposition method, with the DC bias set to the lower limit of I_range, and the AC amplitude uniformly distributed within I_range. The frequency is determined based on f_range and the sweep rate r_sweep, selected at logarithmic intervals within f_range, and the frequency switching rate follows the r_sweep setting. The temperature excitation T(t) = T_0 + ΔT × g(t), where T_0 is the reference temperature of 25℃, ΔT directly adopts the temperature variation range from the excitation scheme, and g(t) is the normalized modulation function. The modulation function is designed as a sawtooth waveform, with the rising segment lasting 0.7 × T_exc and the falling segment lasting 0.3 × T_exc, where T_exc is the excitation duration specified in the scheme. The sweep process is completed within T_exc time, ensuring coverage of the entire frequency range. The composite excitation signal S_comp(t) = I(t) × [1 + α × T(t)] achieves temperature-current modulation, with a coupling coefficient α = 0.02 / ℃. Timing control ensures that all excitation combinations are completed within the time T_exc. Power limiting is set according to the upper limit of I_range to prevent over-excitation. The signal generator is automatically configured according to the P_scheme parameter to generate the temperature-current composite excitation signal S_comp(t).

[0070] The current transformer generates a nonlinear response based on a composite excitation signal. The temperature-current composite excitation signal S_comp(t) is applied to the primary winding of the current transformer through a power amplifier, the bandwidth of which is configured according to the frequency range of S_comp(t). The current transformer's response to S_comp(t) exhibits significant nonlinear characteristics, mainly reflected in the non-proportional relationship between the output magnetic flux density and S_comp(t). When the temperature component in S_comp(t) increases, the core permeability decreases, leading to a reduction in magnetic flux density under the same current excitation, forming a family of temperature-modulated magnetization curves. During the frequency sweep phase of S_comp(t), different frequency components excite different response modes; low-frequency components mainly cause an overall shift in the hysteresis loop, while high-frequency components excite eddy current effects and local magnetic domain resonance. The nonlinear response R_nl(t) is measured synchronously by the secondary-side induced voltage and a flux sensor, with the sampling rate set according to the highest frequency component of S_comp(t). Temperature-current coupling in R_nl(t) manifests as a modulation envelope, the envelope frequency of which is equal to the temperature change frequency. The pulse component of S_comp(t) excites transient oscillations in the iron core, and the oscillation characteristics contain information about the dynamic magnetic properties. The saturation effect becomes significant when S_comp(t) reaches the upper limit of I_range, and R_nl(t) exhibits flat-top distortion. The cross-modulation products appear as combined frequency components in the spectrum of R_nl(t). The nonlinear response R_nl(t) generated under the drive of S_comp(t) fully records the complex dynamic behavior of the current transformer.

[0071] Chaotic transition characteristics are captured based on nonlinear response. The phase space of the nonlinear response R_nl(t) is reconstructed, and an attractor is constructed using the delayed coordinate method. The delay time is determined by the autocorrelation function of R_nl(t). R_nl(t) exhibits an evolution from periodic to chaotic under different excitation intensities: it is a simple limiting cycle under weak excitation, exhibits period-doubling bifurcation as excitation increases, and enters a chaotic state with further excitation. The Lyapunov exponent is calculated by tracking the divergence rate of adjacent trajectories of R_nl(t), with a positive value indicating a chaotic state. During the chaotic transition, the power spectrum of R_nl(t) evolves from a discrete line spectrum to a continuous spectrum. Bifurcation points are determined by monitoring the periodicity of R_nl(t) using the Poincaré section technique. The correlation dimension of R_nl(t) increases from an integer value to a fractional value during the chaotic transition, typically ranging from 2.3 to 2.8. Intermittent chaos in R_nl(t) manifests as alternation between periodic and chaotic windows. The self-similarity of R_nl(t) is analyzed using a recursive graph, and the disruption of the diagonal structure corresponds to the emergence of chaos. The captured chaotic transition feature C_trans includes parameters such as the Lyapunov exponent, bifurcation point, correlation dimension, and intermittency ratio. Changes in these chaotic features directly reflect the nonlinearity of the magnetization process of the transformer core: under normal conditions, it exhibits a stable periodic response with a negative Lyapunov exponent; when the core material experiences fatigue or local defects, the disorder of magnetic domain motion increases, gradually entering a chaotic state, and the Lyapunov exponent turns positive; the increase in correlation dimension indicates increased complexity, corresponding to the deterioration of the internal microstructure of the core; the early appearance of bifurcation points predicts a decrease in the stability of the magnetic material, providing a sensitive indicator for early fault diagnosis.

[0072] Spectral analysis is performed based on chaotic transition characteristics to obtain the coupled response spectrum. The spectral analysis strategy is adaptively adjusted according to the parameters in C_trans: wavelet transform is used to capture broadband features in chaotic states, and FFT is used to extract discrete frequencies in periodic states. The analysis window length is adjusted according to the intermittency ratio in C_trans, shortening the window length when the intermittency ratio is high. Time-frequency analysis is performed on data near bifurcation points, with the resolution set according to the correlation dimension. When approaching chaos, the number of FFT spectral lines is increased to 16384 points to improve frequency resolution. The wavelet basis function is selected according to the chaotic transition speed: Daubechies wavelet is used for fast transitions, and Morlet wavelet is used for slow transitions. Spectral entropy calculation uses different frequency bands for different dynamic states. The bifurcation point information in C_trans determines the frequency range for key analysis. Cross-spectral analysis is performed separately in chaotic and periodic windows. The spectral estimation method is selected according to the Lyapunov exponent: Welch method is used for chaotic states, and periodogram method is used for periodic states. The coupled response spectrum S_coupled(f) is obtained based on the spectral analysis guided by C_trans.

[0073] In some embodiments, extracting fault feature resonance modes from the coupled response spectrum includes: decomposing the coupled response spectrum into a main frequency component and side frequency components; performing phase demodulation on the main frequency component to obtain a phase drift sequence; performing correlation matching between the side frequency components and the phase drift sequence; and identifying fault feature resonance modes through the correlation matching results.

[0074] The coupled response spectrum is decomposed into a main frequency component and sideband components. Main frequency identification is achieved by finding spectral peaks in S_coupled(f). The peak detection algorithm uses local maximum search, requiring the peak value to be more than twice the neighborhood average. The main frequency component is defined as the spectral composition near the center frequency, with its bandwidth determined by the -3dB point. Sideband components are frequency components symmetrically distributed on both sides of the main frequency, with intervals that are integer multiples of the modulation frequency. Spectral decomposition is implemented using a bandpass filter bank, with the center frequency set at both the main frequency and each sideband position. The filter design ensures good frequency selectivity. The main frequency component is extracted using the main frequency filter, preserving the spectrum near the center frequency. Sideband components are extracted separately using the sideband filter bank, with each sideband processed individually. Energy distribution calculations show a significant increase in the sideband energy proportion during faults. Spectral line splitting occurs during severe faults, where the originally single main frequency splits into multiple sub-peaks. The decomposed sets of main frequency and sideband components retain the original amplitude and phase information.

[0075] Phase drift sequence is obtained by phase demodulating the dominant frequency component. First, a Hilbert transform is performed on the dominant frequency component to obtain the analytic signal, and then the instantaneous phase is calculated. A phase expansion algorithm detects the phase difference between adjacent sampling points to eliminate phase jumps. The expanded phase contains two parts: a linear growth trend and a nonlinear drift. The linear trend is obtained by least squares fitting. The phase drift Δφ(t) = φ(t) - φ_linear(t) reflects the nonlinear change characteristics of the phase, where φ(t) is the instantaneous phase and φ_linear(t) is the linear trend. The drift rate is calculated using central difference to improve numerical stability. A reasonable threshold is set for phase abrupt change detection to mark abnormal moments. The RMS value of the phase drift is used to evaluate overall stability. Allan variance analysis is used to analyze the phase stability at different time scales. The modulation depth is obtained by calculating the peak-to-peak value of the phase drift. Filtering removes high-frequency noise while retaining effective phase change information. The final complete phase drift sequence is recorded at fixed sampling intervals.

[0076] The sideband components are correlated and matched with the phase drift sequence. Correlation analysis first extracts the amplitude envelope of each sideband component and obtains the envelope line using Hilbert transform. Cross-correlation calculation uses a normalization method to search for correlation peaks to determine the time delay relationship. The position of the correlation peak determines the time delay τ between the sideband and phase drift; a positive value indicates that the sideband lags behind the phase change. The synchronization index SI calculates the maximum correlation coefficient; a high value indicates a strong correlation. Frequency consistency testing compares the matching degree between the sideband spacing and the phase modulation frequency. Amplitude-phase transfer characteristics are obtained through frequency domain analysis to describe the coupling mechanism. Nonlinear correlation uses mutual information quantization to capture higher-order statistical correlations. Sliding window analysis generates time-varying correlation curves to track coupling dynamics. Granger causality test determines the influence direction D; D=1 indicates that the phase causes the sideband, and D=-1 indicates that the sideband causes the phase. Matching strength is graded according to the correlation coefficient, distinguishing between strong, medium, and weak matches. The matching results of each sideband form an correlation matrix M_corr, containing three key parameters: correlation coefficient SI, time delay τ, and causal direction D.

[0077] Fault characteristic resonance modes are identified through correlation matching results. Resonance mode identification first analyzes the dominant frequency offset, calculating the difference between the measured resonance frequency and the normal value Δf_res = f_res - f_res0, where f_res is the measured frequency and f_res0 is the normal frequency. The quality factor is calculated using Q = f_res / Δf_3dB, where Δf_3dB is the -3dB bandwidth; a decrease in Q indicates increased damping. Mode classification is determined based on the frequency range and the time delay τ in the correlation matrix M_corr: a fast response with τ < 10ms corresponds to electromagnetic resonance, 10ms < τ < 100ms corresponds to mechanical resonance, and τ > 100ms corresponds to thermal resonance. Coupling strength is directly evaluated using the correlation coefficient SI in M_corr; SI > 0.8 indicates a strongly coupled fault. Fault source localization utilizes the causal direction D; D = 1 indicates the fault originates from the main system influencing sidebands, and D = -1 indicates that abnormal sidebands react on the main system. Fault propagation speed is calculated using the time delay τ and the distance between nodes, v = L / τ. For example, a shift in the mechanical resonant frequency from 45Hz to 52Hz was detected, with the Q value decreasing from 15 to 8. The correlation matrix showed SI = 0.85 (strong coupling), τ = 25ms (mechanical response), and D = 1 (main frequency drive) at the 104Hz sideband. This was diagnosed as a mechanical fault caused by winding loosening developing into electromagnetic coupling. Through feature extraction and correlation parameter application, the fault characteristic resonance mode R_fault = {f_res, Q, coupling strength SI} was obtained.

[0078] Step S160: Analyze the fault feature resonance mode and scale coupling factor to obtain adaptive window parameters, obtain multi-resolution deviation through the adaptive window parameters, and generate a spatial fault distribution map based on the multi-resolution deviation.

[0079] Specifically, an adaptive window parameter is obtained by comprehensively analyzing the fault characteristic resonance mode R_fault={f_res,Q,SI} and the scale coupling factor S_c generated by S120. The window parameter design needs to consider both resonance characteristics and time-scale characteristics. The window length is related to the quality factor Q and the resonance frequency f_res, ensuring that it includes at least enough resonance periods to accurately capture oscillation characteristics. For example, when the mechanical resonance frequency of a 500kV transformer is detected to be 45Hz and the quality factor Q=12, the analysis window needs to cover about 12 oscillation periods to accurately identify the resonance characteristics. For damped oscillations with lower Q values, the window length can be shortened accordingly. The multi-resolution layer design is logarithmically divided from low frequency to high frequency. The bottom layer analyzes low-frequency thermal effects and mechanical displacement, the middle layer analyzes power frequency and its harmonics, and the top layer analyzes high-frequency partial discharge and switching transients. The time resolution of each layer is determined according to the center frequency, maintaining an inverse relationship between time and frequency resolution. The low-frequency layer provides high time resolution to track slow changes, and the high-frequency layer provides high frequency resolution to identify harmonic components. Boundary handling employs mirror extension to avoid edge effects. When the analysis window approaches the data boundary, the data length is extended through mirror reflection. Taking a 110kV distribution transformer with winding deformation as an example, its resonant frequency shifts to 120Hz and the coupling strength significantly increases. In this case, the window parameters are automatically adjusted to a shorter window length and more analysis layers to accurately capture the complex vibration modes caused by the deformation. The adaptive window parameter set W_params includes the window length, the number of multi-resolution layers, and the time resolution parameters for each layer.

[0080] In some embodiments, obtaining multi-resolution deviation through the adaptive window parameters includes: constructing a delay coordinate system based on the adaptive window parameters; projecting the magnetic field signal onto a multi-dimensional phase space based on the delay coordinate system; extracting deviation features based on trajectory evolution in the multi-dimensional phase space; and obtaining multi-resolution deviation based on the spatial distribution of the deviation features.

[0081] A delay coordinate system is constructed based on adaptive window parameters. The magnetic field signal uses the original magnetic flux density recovery process data acquired in S110 as the time series input for constructing the delay coordinate system. The construction of the delay coordinate system utilizes the time resolution parameter in W_params to determine the delay time of each layer. The delay time is determined according to the window length and sampling rate to ensure that the main dynamics within the window are captured. For example, for a 220kV transformer, when the detection window is set to 200ms, the delay time interval is set according to the main oscillation frequency of the transformer to capture the complete periodic characteristics of magnetic flux density changes. The embedding dimension is determined according to the ratio of the window length and the resonant period to ensure that the phase space is fully expanded. When the transformer winding becomes loose, the vibration frequency is usually in the range of 100-300Hz, requiring a higher embedding dimension to capture rapidly changing dynamic characteristics, while slow changes such as core saturation are handled by a lower embedding dimension. Multi-scale delay is achieved by doubling the layers one by one, corresponding to a halving of the frequency analysis range. The first layer captures high-frequency transient characteristics such as switching impact response, the second layer analyzes mid-frequency mechanical vibration, and the third layer extracts low-frequency thermally induced changes. The delay vector contains the values ​​of the original signal at different historical moments, forming a high-dimensional state representation. This method uses multiple "snapshots" to reconstruct the complete dynamic picture, with each delay component representing a "memory" at a specific moment. Delay coordinates within a window maintain consistency, avoiding discontinuities across windows. Boundary processing employs cyclic delay, utilizing the quasi-periodic characteristics of the magnetization process. As the analysis window moves to the data boundary, cyclic extension is performed using the periodicity of the magnetic field signal. The orthogonality of the coordinate system is enhanced through the Gram-Schmi dt process, eliminating correlations between different delay components. Normalization eliminates dimensional effects, ensuring that each component has equal weight in phase space analysis. The constructed delay coordinate system provides a framework for phase space projection, transforming the signal characteristics in the time domain into the geometry of the spatial domain.

[0082] For example, projecting the magnetic field signal onto a multidimensional phase space based on the delay coordinate system includes: determining the embedding dimension and delay time based on the delay coordinate system; constructing a high-dimensional state vector based on the embedding dimension and the delay time; generating a phase space trajectory by performing a nonlinear mapping based on the high-dimensional state vector; and generating a multidimensional phase space projection based on the phase space trajectory.

[0083] The embedding dimension and delay time are determined based on the delayed coordinate system. The embedding dimension *m* is determined using the false nearest neighbor method, starting with *m*=2 and incrementally increasing, calculating the distance change of nearest neighbors in phase space for each dimension. False neighbors are points that are close in the low-dimensional projection but far apart in the high-dimensional projection; these points cause trajectory self-intersection and loss of dynamic information. The rate of change of nearest neighbor distance is calculated after increasing the dimension; a rate of change exceeding 10 times is considered a false neighbor. The proportion of false neighbors decreases with increasing dimension; the minimum *m* value when the proportion drops below 1% is taken as the optimal embedding dimension, typically between 3 and 7. The delay time *τ* is determined using the mutual information method, calculating the mutual information *I(τ)* between the original signal *s(t)* and the delayed signal *s(t+τ)*. Mutual information decreases with increasing *τ* starting from its maximum value; the *τ* corresponding to the first local minimum is selected, maintaining necessary correlation while retaining independent information. Verification uses the phase space expansion index, calculating the average neighborhood radius of the trajectory to ensure sufficient expansion. Multi-scale embedding determines independent (m_i, τ_i) parameter pairs for each resolution layer, with lower frequency layers requiring greater delay. The embedding window length must be greater than (m-1)×τ to ensure the integrity of the state vector construction.

[0084] A high-dimensional state vector is constructed based on embedding dimension and time delay. The state vector extracts m delayed components from the time series s(t), forming a vector representation X(t) = [s(t), s(t-τ), s(t-2τ), ..., s(t-(m-1)τ)]. Each component represents the system's state at different historical moments, reflecting the overall dynamic evolution process. Vector normalization first calculates the mean and standard deviation of each component, then normalizes it to zero mean and unit variance, eliminating amplitude differences between different delays. Integrity checks ensure all components have valid data, especially at sequence boundaries. Boundary processing uses symmetric extension, mirroring the signal about the endpoints to generate additional data. Temporal continuity is achieved by checking the Euclidean distance between adjacent state vectors; abrupt changes in distance indicate data anomalies. Vector orthogonality uses QR decomposition, transforming related delayed components into orthogonal bases to enhance independence. Sparse representation is achieved through principal component analysis; extracting the first 3-5 principal components typically retains over 95% of the information. Quality assessment calculates the error between the original signal and the signal reconstructed from the first component of the state vector. The constructed high-dimensional state vector sequence accurately characterizes the phase space structure of the system.

[0085] Phase space trajectories are generated based on nonlinear mapping of high-dimensional state vectors. The nonlinear mapping employs the Locally Linear Embedding (LLE) algorithm to preserve the local geometric structure of the data. First, k nearest neighbors are found for each state vector, where k is typically 2-3 times the embedding dimension. The reconstruction weight matrix W is calculated to minimize the error of each point represented by a linear combination of its nearest neighbors, with the weights satisfying normalization constraints. Coordinates Y are found in the low-dimensional space to minimize the reconstruction error under the same weights, preserving neighborhood relationships. Trajectory generation connects the mapped points in chronological order, with adjacent points connected by line segments to form continuous trajectories. Trajectory smoothing uses cubic spline interpolation, with the control points being the original mapped points, ensuring the continuity of the first and second derivatives. Geometric feature calculations include tangent vectors (reflecting velocity), curvature (reflecting the degree of bending), and torsion (reflecting three-dimensional distortion). Trajectory segments are classified by local dynamic features: periodic segments exhibit closed loops, and chaotic segments exhibit complex entanglement. Abnormal trajectory segments are identified by Fréchet distance with a normal trajectory library; distances exceeding a threshold are marked as abnormal. The generated phase space trajectory fully preserves the nonlinear characteristics of the system.

[0086] Multi-dimensional phase space projections are generated based on phase space trajectories. The projection design showcases different dynamic aspects of the trajectory. Position projection directly uses trajectory coordinates to display the occupancy of the state space. Velocity projection calculates the tangent vector through numerical differentiation; the vector length represents the velocity magnitude, and the direction represents the motion trend. Acceleration projection calculates the second derivative, reflecting the forces and energy conversions of the system. The position-velocity phase plane is a classic projection, with position on the horizontal axis and velocity on the vertical axis; conservative systems form closed curves, while dissipative systems form spirals. The projection resolution adapts to local trajectory density, refining the mesh in dense areas to capture details and coarsening it in sparse areas to save computation. Color coding uses a time gradient, with cool colors for early trajectories and warm colors for later ones, intuitively showing the evolution process. Other physical quantities such as energy and temperature can also be encoded to increase the information dimension. Statistical features are extracted, including the centroid (average state), covariance matrix (state dispersion), and higher-order moments (distribution shape). Multi-scale projection is achieved by adjusting the observation window; short windows display transients, and long windows display trends. Cross-projection simultaneously displays multiple physical quantity trajectories, analyzing coupling and phase relationships. The quantitative analysis of the projection calculates geometric parameters such as the area occupied, perimeter, and fractal dimension.

[0087] Deviation features are extracted based on trajectory evolution in multidimensional phase space. First, a normal trajectory reference library is established, containing trajectories under various healthy operating states such as light load, full load, and start-stop. The comparison between the current trajectory and the reference library uses the Hausdorff distance, defined as the maximum of the minimum distances between all point pairs, reflecting the overall degree of deviation. Local deviations are analyzed through a sliding window, calculating the average distance between a fixed-length trajectory segment and the most similar reference segment. The deviation direction is determined through principal component analysis, calculating the principal direction of the deviation vector; the first principal component typically corresponds to the main fault mode. Deviation velocity is the time derivative of the deviation degree, and deviation acceleration is the second derivative; rapid growth indicates fault deterioration. Time patterns are identified through template matching, with predefined templates for sudden (step), gradual (ramp), periodic (sine), and intermittent (pulse) patterns. Statistical features are calculated within the sliding window, including deviation mean (average level), standard deviation (stability), skewness (asymmetry), and kurtosis (extreme frequency). Geometric features focus on trajectory deformation; a sudden increase in curvature indicates a sharp turn, and changes in torsion indicate spatial distortion. Dynamic characteristics are assessed for stability changes using invariants such as the Lyapunov exponent and correlation dimension. The eigenvector F_dev integrates all deviation information, providing a quantitative basis for fault diagnosis.

[0088] Multi-resolution deviation is derived from the spatial distribution of deviation features. Spatial distribution analysis establishes a mapping between deviation features and physical locations, with each measurement point corresponding to a specific location on the transformer, such as the winding layer or core area. The transformer meshing considers actual geometry; cylindrical coordinates are used for the windings to accommodate the cylindrical shape, while Cartesian coordinates are used for the core to accommodate the square structure. Due to the limited number of sensors, spatial interpolation is needed to estimate the deviation features of unmeasured points. Kriging interpolation considers spatial correlation, assuming that nearby points are similar, and that the correlation decays exponentially with distance. The variogram is fitted from measured data to determine the correlation length and nugget value. Interpolation weights are obtained by solving the Kriging equations to ensure unbiased optimal estimation. Multi-resolution uses a mesh pyramid, with a 100×100 fine mesh at the bottom, merging layer by layer into a single mesh at the top. For each layer's deviation calculation, the feature vector F_dev is normalized and its modulus is taken, with reference values ​​derived from historical normal data. Inter-layer propagation uses a constraint operator for upward aggregation and an extension operator for downward refinement. Temporal evolution is obtained through repeated spatial analysis, with the sampling interval matched to the fault development rate. Spatial smoothing employs anisotropic diffusion, maintaining sharpness along structural boundaries and smoothing noise in the vertical direction. Anomaly identification is based on local statistics, with anomalies defined as exceeding the mean plus three standard deviations. The completed D_multi(x,y,scale,t) provides multi-scale spatiotemporal deviation information.

[0089] A spatial fault distribution map is generated based on multi-resolution deviation. The spatial deviation distribution is extracted from D_multi and mapped to a physical coordinate system consistent with the actual equipment. A perceptually linear color scheme is used: deviations of 0-0.3 are represented by blue-green (safe), 0.3-0.7 by yellow-orange (warning), and 0.7-1.0 by red-purple (dangerous). Contour lines are generated using a contour tracing algorithm; primary contour lines of 0.5 and 0.8 are thick and labeled with values, while secondary contour lines are drawn as thinner lines every 0.1. Fault regions are defined as connected regions with deviations exceeding 0.5, and connected component analysis is performed using the 8-adjacency criterion. Key features are extracted for each fault region: geometric center for location, area to quantify the impact range, maximum deviation indicating severity, boundary perimeter reflecting shape complexity, and principal axis direction indicating the expansion trend. The gradient field is calculated using finite difference; the gradient vector indicates the direction of fastest deviation growth, with larger gradients where contour lines are dense. The propagation path is traced along the gradient direction, pointing from low-deviation regions to the high-deviation center. The time stamps include the first occurrence time, peak time, and current state, tracking the fault lifecycle. Statistical information summarizes the current number of faulty areas, the percentage of total area, the most severe location and value, and an area distribution histogram. Through color mapping, contour plotting, fault area identification, and gradient analysis, the spatial fault distribution map G_fault(x,y,t) is finally obtained.

[0090] Step S170: The comprehensive fault index is obtained by integrating the spatial fault distribution map and the fault precursor characteristics. The graded diagnosis result is generated based on the comprehensive fault index to complete the fault detection of the DC current transformer.

[0091] Specifically, the spatial fault distribution map G_fault(x,y,t) and the fault precursor features [r_W,W_rel,W_h(f)] are fused and analyzed from multiple dimensions to generate a unified comprehensive fault index. Spatial feature extraction obtains three quantitative indicators from G_fault(x,y,t): the fault area ratio is obtained by statistically analyzing the percentage of grid points with a deviation greater than 0.5 in the total monitored area, reflecting the spatial diffusion degree of the fault; the maximum deviation is directly read from the peak value of G_fault, representing the local severity of the fault; the number of fault areas is determined through connected component analysis, with multiple independent areas suggesting the complexity of the fault's causes. Temporal feature processing standardizes the three precursor features of S130: the hysteresis width growth rate r_W is compared with historical data, with daily growth rates exceeding 1% given high weight; the relative hysteresis width W_rel is already in percentage form, directly reflecting the current degree of degradation; and the frequency-related hysteresis width W_h(f) is calculated using the average value of the main frequency points. The fusion strategy employs a hierarchical weighting approach. First, spatial features are weighted internally: the fault area ratio is weighted at 0.3, the maximum deviation at 0.4, and the number of regions at 0.3, yielding a spatial composite value. Temporal features are then weighted internally: the growth rate r_W is weighted at 0.4 to reflect dynamic trends, W_rel at 0.3 to reflect static levels, and W_h(f) at 0.3 to characterize structural state, yielding a temporal composite value. Finally, the comprehensive fault index is obtained by weighted summation of the spatial composite value (weight 0.4) and the temporal composite value (weight 0.6), with a numerical range of 0-10 for easy intuitive understanding and subsequent classification.

[0092] Based on the comprehensive fault index, a four-level diagnostic standard system corresponding to the actual operating status of the equipment is established. Level 1 Normal State: Comprehensive fault index less than 2.0, indicating that all performance indicators of the instrument transformer are within the design limits, the magnetic field distribution is uniform, the timing characteristics are stable, and the equipment can operate safely and reliably. It is recommended to maintain the existing maintenance cycle and monitoring frequency; no special attention is required. Level 2 Mild Abnormality: Comprehensive fault index in the range of 2.0-4.0, indicating early signs of equipment degradation, possibly manifested as slight local magnetic field distortion or a slow performance degradation trend. Although it does not affect normal operation, it requires attention. It is recommended to shorten the routine inspection cycle by half, add trend monitoring points, and establish a degradation tracking file to provide data support for subsequent maintenance decisions. Level 3 Moderate Fault: A comprehensive fault index between 4.0 and 7.0 indicates obvious fault characteristics in the equipment, with significant spatial anomalies or accelerated deterioration in timing. Equipment performance declines but basic functions can still be maintained. A detailed maintenance plan is required, including maintenance time windows, necessary spare parts, and technical personnel allocation. Before maintenance, the monitoring frequency should be increased to once daily, and automatic alarms should be set for parameters exceeding limits to ensure maintenance is completed before complete equipment failure. Level 4 Severe Fault: A comprehensive fault index greater than or equal to 7.0 indicates the equipment is in a dangerous operating state. The fault has developed to the point of affecting safe operation. Large-scale spatial anomalies or abrupt temporal changes may occur. Continued operation carries the risk of equipment damage, protection malfunction, or system failure. Immediate measures should be taken, including load transfer, equipment isolation, and emergency shutdown. A professional maintenance team should be organized for comprehensive inspection and repair. If necessary, the equipment should be replaced directly. The DC current transformer fault detection is then completed.

[0093] To implement the DC current transformer fault detection method based on magnetic flux density changes corresponding to the above method embodiments, in order to achieve the corresponding functions and technical effects. See Figure 2 , Figure 2 This diagram illustrates a structural block diagram of a DC current transformer fault detection device 200 based on a magnetic flux density variation according to an embodiment of this application. For ease of explanation, only the parts relevant to this embodiment are shown. The DC current transformer fault detection device 200 based on a magnetic flux density variation according to an embodiment of this application includes:

[0094] Feature extraction module 201 is used to apply a gradual pulse excitation sequence to a DC current transformer to collect magnetic flux density recovery process data, generate a family of time-series recovery curves based on the magnetic flux density recovery process data, and extract magnetic field fatigue feature spectrum from the family of time-series recovery curves.

[0095] The scale analysis module 202 is used to perform time-scale decomposition using the magnetic field fatigue characteristic spectrum to obtain the fast recovery component and the slow recovery component, perform ratio analysis on the fast recovery component and the slow recovery component to generate a scale coupling factor, and construct a magnetic field health status benchmark through the scale coupling factor.

[0096] The bidirectional testing module 203 is used to perform a bidirectional approximation test based on the magnetic field health status benchmark to obtain the rising critical point and the falling critical point, determine the hysteresis width parameter based on the rising critical point and the falling critical point, and extract fault precursor features from the hysteresis width parameter.

[0097] Cross-analysis module 204 is used to perform curve family cross-analysis on the fault precursor features to obtain a set of intersection points, extract a cross-angle sequence from the set of intersection points, and generate a fault development speed index based on the cross-angle sequence.

[0098] The coupling excitation module 205 is used to design a parameter coupling excitation scheme according to the fault development speed index, use the parameter coupling excitation scheme to excite the current transformer to obtain the coupling response spectrum, and extract the fault characteristic resonance mode from the coupling response spectrum.

[0099] Adaptive analysis module 206 is used to analyze the fault feature resonance mode and the scale coupling factor to obtain adaptive window parameters, obtain multi-resolution deviation through the adaptive window parameters, and generate a spatial fault distribution map based on the multi-resolution deviation.

[0100] The comprehensive diagnostic module 207 is used to perform fusion analysis on the spatial fault distribution map and the fault precursor features to obtain a comprehensive fault index, generate a graded diagnostic result based on the comprehensive fault index, and complete the fault detection of the DC current transformer.

[0101] The aforementioned DC current transformer fault detection device 200 based on magnetic flux density variation can implement the DC current transformer fault detection method based on magnetic flux density variation described in the above method embodiments. The options in the above method embodiments are also applicable to this embodiment and will not be detailed here. The remaining contents of this application embodiment can be referred to the contents of the above method embodiments, and will not be repeated in this embodiment.

[0102] like Figure 3 As shown, the third embodiment of the present invention also provides a computer device, including a memory 301, a processor 302, and a computer program stored in the memory 301 and executable on the processor 302. The feature is that when the processor 302 executes the program, it implements the steps of the DC current transformer fault detection method based on magnetic flux density variation described in the first embodiment of the present invention.

[0103] The purpose of the above embodiments is to reproduce and derive the technical solution of the present invention by way of example, and to fully describe the technical solution, purpose and effect of the present invention. The purpose is to enable the public to have a more thorough and comprehensive understanding of the disclosure of the present invention, and not to limit the scope of protection of the present invention.

[0104] The above embodiments are not an exhaustive list based on the present invention, and there may be many other embodiments not listed. Any substitutions and improvements made without departing from the concept of the present invention are within the protection scope of the present invention.

Claims

1. A method for detecting faults in a DC current transformer based on changes in magnetic flux density, characterized in that, include: A gradual pulse excitation sequence is applied to a DC current transformer to collect magnetic flux density recovery process data. A family of time-series recovery curves is generated based on the magnetic flux density recovery process data. The magnetic field fatigue characteristic spectrum is extracted from the family of time-series recovery curves. The process involves using the magnetic field fatigue characteristic spectrum to perform time-scale decomposition to obtain a fast recovery component and a slow recovery component, and then performing ratio analysis on the fast and slow recovery components to generate a scale coupling factor. This includes: performing a logarithmic transformation on the fast recovery component to obtain nonlinear decay characteristics; projecting the slow recovery component onto phase space to extract trajectory density; generating an implicit coupling mode through cross-mapping between the nonlinear decay characteristics and the trajectory density; deriving the scale coupling factor using the implicit coupling mode; and constructing a magnetic field health state benchmark using the scale coupling factor. Based on the magnetic field health status benchmark, a bidirectional approximation test is performed to obtain the rising critical point and the falling critical point. The hysteresis width parameter is determined based on the rising critical point and the falling critical point, including: triggering a micro-perturbation test at the rising critical point and recording the response amplitude; repeating the perturbation test at the falling critical point to obtain the decayed amplitude; comparing the response amplitude and the decayed amplitude in a time sequence to generate an amplitude deviation sequence; generating the hysteresis width parameter based on the extreme value distribution of the amplitude deviation sequence; and extracting fault precursor features from the hysteresis width parameter. Curve family cross analysis is performed on the aforementioned fault precursor features to obtain a set of intersection points. Cross angle sequences are extracted from the set of intersection points, and a fault development speed index is generated based on the cross angle sequences. Design a parameter coupling excitation scheme according to the fault development speed index, use the parameter coupling excitation scheme to excite the current transformer to obtain the coupling response spectrum, and extract the fault characteristic resonance mode from the coupling response spectrum; The fault characteristic resonance mode and the scale coupling factor are analyzed to obtain adaptive window parameters. Multi-resolution deviation is obtained through the adaptive window parameters, and a spatial fault distribution map is generated based on the multi-resolution deviation. By combining the spatial fault distribution map and the fault precursor features, a comprehensive fault index is obtained. Based on the comprehensive fault index, a graded diagnosis result is generated to complete the fault detection of the DC current transformer.

2. The method according to claim 1, characterized in that, The step of extracting the intersection angle sequence from the set of intersection positions includes: A rotational transformation is applied to the set of intersection points to construct an angle field; Searching for gradient abrupt change lines in the angle field to identify angle transition paths; The cumulative angle change is obtained by integrating along the angle transition path. The cumulative angle changes are arranged in chronological order to form a cross-angle sequence.

3. The method according to claim 1, characterized in that, The step of using the parameter coupling excitation scheme to excite the current transformer and obtain the coupling response spectrum includes: A temperature-current composite excitation signal is generated based on the aforementioned parameter coupling excitation scheme; The current transformer is driven to generate a nonlinear response based on the temperature-current composite excitation signal. Based on the aforementioned nonlinear response, chaotic transition characteristics are captured; Based on the chaotic transition characteristics, a spectral analysis is performed to obtain the coupled response spectrum.

4. The method according to claim 1, characterized in that, The extraction of fault feature resonance modes from the coupled response spectrum includes: The coupled response spectrum is decomposed into main frequency components and side frequency components; Phase demodulation is performed on the main frequency component to obtain the phase drift sequence; The sideband components are associated and matched with the phase drift sequence; The resonance patterns of fault features are identified through the correlation matching results.

5. The method according to claim 1, characterized in that, The process of obtaining multi-resolution deviation through the adaptive window parameters includes: Construct a delayed coordinate system based on the adaptive window parameters; The magnetic field signal is projected into a multidimensional phase space based on the aforementioned delay coordinate system; Deviation features are extracted based on trajectory evolution in the multidimensional phase space; Multi-resolution deviation is obtained based on the spatial distribution of the deviation features.

6. The method according to claim 5, characterized in that, The projection of the magnetic field signal onto the multidimensional phase space based on the delayed coordinate system includes: The embedding dimension and delay time are determined based on the aforementioned delay coordinate system; A high-dimensional state vector is constructed based on the embedding dimension and the delay time; A phase space trajectory is generated by nonlinear mapping based on the high-dimensional state vector; A multidimensional phase space projection is generated based on the phase space trajectory.

7. A fault detection device for a DC current transformer based on changes in magnetic flux density, characterized in that, include: The feature extraction module is used to apply a gradual pulse excitation sequence to a DC current transformer to collect magnetic flux density recovery process data, generate a family of time-series recovery curves based on the magnetic flux density recovery process data, and extract the magnetic field fatigue feature spectrum from the family of time-series recovery curves. The scale analysis module is used to perform time-scale decomposition of the magnetic field fatigue characteristic spectrum to obtain fast recovery components and slow recovery components, and to perform ratio analysis on the fast recovery components and the slow recovery components to generate a scale coupling factor. This includes: performing a logarithmic transformation on the fast recovery components to obtain nonlinear decay characteristics; projecting the slow recovery components onto phase space to extract trajectory density; generating implicit coupling modes through cross-mapping of the nonlinear decay characteristics and the trajectory density; deriving the scale coupling factor using the implicit coupling mode; and constructing a magnetic field health state benchmark using the scale coupling factor. A bidirectional testing module is used to perform bidirectional approximation tests based on the magnetic field health status benchmark to obtain rising and falling critical points, and to determine hysteresis width parameters based on the rising and falling critical points. This includes: triggering a micro-perturbation test at the rising critical point and recording the response amplitude; repeating the perturbation test at the falling critical point to obtain the decay amplitude; performing a time-series comparison between the response amplitude and the decay amplitude to generate an amplitude deviation sequence; generating hysteresis width parameters based on the extreme value distribution of the amplitude deviation sequence; and extracting fault precursor features from the hysteresis width parameters. The cross-analysis module is used to perform curve family cross-analysis on the fault precursor features to obtain a set of intersection points, extract a cross-angle sequence from the set of intersection points, and generate a fault development speed index based on the cross-angle sequence. The coupling excitation module is used to design a parameter coupling excitation scheme according to the fault development speed index, use the parameter coupling excitation scheme to excite the current transformer to obtain the coupling response spectrum, and extract the fault characteristic resonance mode from the coupling response spectrum. An adaptive analysis module is used to analyze the fault feature resonance mode and the scale coupling factor to obtain adaptive window parameters, obtain multi-resolution deviation through the adaptive window parameters, and generate a spatial fault distribution map based on the multi-resolution deviation. The comprehensive diagnostic module is used to perform fusion analysis on the spatial fault distribution map and the fault precursor features to obtain a comprehensive fault index, generate a graded diagnostic result based on the comprehensive fault index, and complete the fault detection of the DC current transformer.

8. A computer device, characterized in that, It includes a memory and a processor; the memory is used to store a computer program; the processor is used to execute the computer program and, when executing the computer program, implement the method as described in any one of claims 1 to 6.

Images