Dynamic magnetic disturbance compensation method and system for magnetic characteristic test of aircraft

Abstract

The invention relates to the field of magnetic characteristic testing, in particular to a dynamic magnetic disturbance compensation method for magnetic characteristic testing of an aircraft, and the method comprises the steps: (1) selecting a target site as a testing region; (2) arranging an environment monitoring device at the periphery of the target area, and performing magnetic source positioning based on a measurement result of the environment monitoring device; and (3) performing dynamic compensation based on the positioned magnetic source position and magnetic source intensity and the position of the target area. According to the method, a static magnetic source positioning analysis algorithm is adopted, a closed-loop structure of'magnetic gradient tensor + extended Kalman filter (EKF) 'is adopted, traditional particle swarm optimization and active excitation links are canceled, the problem of initial value sensitivity is avoided, external magnetic field interference is not introduced, millisecond-level static magnetic source positioning is achieved, and the positioning accuracy is improved. And the requirements of the magnetic characteristic test on rapid, accurate and disturbance-free compensation are met.

Description

Technical Field

[0001] This invention relates to a magnetic interference compensation system for magnetic performance testing, belonging to the field of precision magnetic field control and signal processing technology. Background Technology

[0002] During the navigation of a spacecraft, magnetic fields are frequently used for purposes such as confirming direction and searching for targets.

[0003] Existing magnetic property testing systems mostly use static background subtraction or uniform feedback compensation methods, which cannot accurately locate the position and movement of the magnetic source, nor can they distinguish different types of interference such as geomagnetic diurnal variation, mechanical disturbance and power grid harmonics, resulting in compensation lag, overcompensation or undercompensation, and poor test repeatability. Summary of the Invention

[0004] To address this issue, this invention proposes a dynamic magnetic interference compensation method and system for aircraft, used to maintain the stability of the aircraft's magnetic properties. The magnetic interference compensation system includes a magnetic field compensation coil winding and a compensation algorithm integrated into the control module. The focus of this invention is on the improvement and innovation of environmental magnetic interference compensation, ensuring the stability of the magnetic field in the aircraft testing system.

[0005] The present invention aims to solve the problem that existing uniform magnetic field generating systems are easily affected by external disturbances in complex environments, and proposes a dynamic magnetic interference compensation method and system.

[0006] On one hand, the present invention provides a dynamic magnetic interference compensation method for testing the magnetic characteristics of a spacecraft, the method comprising:

[0007] Step (1) Select the target site as the test area; Step (2) Set up an environmental monitoring device around the target area and locate the magnetic source based on the measurement results of the environmental monitoring device;

[0008] Step (2) includes: (2.1) Constructing a six-dimensional state vector, which includes the three-dimensional spatial position and three-axis magnetic moments of the magnetic source; (2.2) Establishing the state transition system and observation model of the magnetic source, and constructing the observation function of the magnetic source. ,in To observe the noise, X k It is a six-dimensional state vector, that is ,in, Let be the spatial coordinates of the magnetic source at time k. (2.3) Linearize the observation model to obtain the linearized observation matrix of the magnetic source: ;

[0009] (2.4) Use Kalman gain calculation to perform error covariance update iteration and state vector update of observations;

[0010] (2.5) Estimate and iteratively update the observation noise covariance and process noise covariance; (2.6) Repeat steps (2.1)-(2.5) above to output the optimal state estimate of the magnetic source. This state estimate includes location. and magnetic moment Step (3): Based on the location and intensity of the magnetic source after positioning and the location of the target area, determine the magnetic field generated by the magnetic source in the target area; Step (4): Based on the calculated magnetic field generated by the magnetic source at the target location, generate a reverse magnetic field at the target location.

[0011] Furthermore, the method includes: performing magnetic field measurement using a cross-shaped planar measurement array, which includes four or 4N vector magnetic sensors; the tensor matrix expression of the cross-shaped planar measurement array at the center point o is as follows:

[0012] , among which, (B) 1x B 1y B 1z ,) represents the triaxial magnetic field values ​​measured by the vector magnetic sensor on the first cross branch of the cross-shaped planar measurement array, (B) 2x B 2y B 2z ,) represents the triaxial magnetic field values ​​measured by the vector magnetic sensor on the second cross branch of the cross-shaped planar measurement array, (B) 3x B 3y B 3z ,) represents the triaxial magnetic field values ​​measured by the vector magnetic sensor on the third cross branch of the cross-shaped planar measurement array, (B) 4x B 4y B 4z ,) represents the three-axis magnetic field values ​​measured by the vector magnetic sensor on the fourth cross branch of the cross-shaped measurement array.

[0013] Furthermore, the method includes: collecting several sets of initial observation data to statistically analyze and determine the initial covariance. This step includes: Step (2.21), collecting N sets of initial magnetic field data Z1, Z2, ... Z N And estimate the initial state using the initial observation data. ;

[0014] Step (2.22): Calculate the variance of the estimated initial state:

[0015] ,

[0016] ,

[0017] ,

[0018] ,

[0019] ,

[0020] ;x i ,y i ,z i Let m represent the coordinates of the magnetic source measured in the i-th group. xi m yi , m zi Representing the magnetic moment of the magnetic source, step (2.23) combines the above variances into a 6×6 diagonal matrix:

[0021] .

[0022] Furthermore, the method also includes: collecting initial observation data of the magnetic field to determine the initial covariance of the magnetic field as the initial value of the prior covariance; and for subsequent time steps, using the posterior covariance of the previous time step. process noise covariance Q k Perform the prior error covariance at the next time step. predict, Q k The process noise covariance is obtained by collecting magnetic field information from the magnetic source and calculating the variance of the estimated magnetic moment.

[0023] Furthermore, the method also includes: utilizing the obtained observation matrix H k Substitute into the formula The Kalman gain was calculated. , To observe the noise covariance, the observation residuals and Kalman gain are used. To correct the prior state, we obtain the posterior state:

[0024] ,in It is the difference between the "actual observed value" and the "observation value predicted based on the prior state";

[0025] Updated posterior error covariance: , where I is the identity matrix.

[0026] Further, step (4) includes: based on the formula To calculate the compensation magnetic field strength at the target location, The position vector of the magnetic source relative to the center of the compensation coil. To extend the magnetic source moment of the Kalman filter output, is the vacuum permeability.

[0027] Furthermore, the observation function is: ,

[0028] Where, m x m y m z , , and z are the magnetic moments of the magnetic source in the x, y, and z directions, respectively, r is the distance from the target point to the magnetic source, and x, y, and z are the coordinates of the target point.

[0029] The method further includes: classifying magnetic field interference and compensating for different types of interference. This step includes classifying the interfering magnetic field into: near-field interference generated by moving objects, far-field interference generated by moving objects, geomagnetic field interference, and power grid interference. For near-field interference, step (2) of the method is used to locate the magnetic source and determine the location of the magnetic source based on the formula. To calculate the compensation magnetic field strength at the target location, The position vector of the magnetic source relative to the center of the compensation coil. To extend the magnetic source moment of the Kalman filter output, The magnetic permeability is given by the vacuum magnetic field. For far-field interference, compensation is made using the reverse magnetic field strength based on the magnetic field measured by the magnetic sensor located in the target area. For grid interference, the magnetic field data of the peak area of ​​the grid interference is brought into step (2) of the method described in claim 1 according to the measured grid cycle to locate the magnetic source of the grid interference. Based on the located magnetic source position and magnetic source strength of the grid interference, the magnetic source is located using the formula... The peak compensation magnetic field strength at the target location is calculated, and the peak compensation magnetic field strength with the same period is taken as the maximum value for reverse magnetic field compensation.

[0030] The process of determining the magnetic interference category includes: arranging N triaxial magnetic sensors at equal intervals in three dimensions at different locations within the target area; each magnetic sensor measures the magnetic field at its location to obtain the corresponding magnetic sensor signal.

[0031] (1)

[0032] The obtained magnetic sensor signal is preprocessed; the variance of the magnetic field within the measurement time period T is measured.

[0033] Measure the kurtosis of the magnetic field within a time period T;

[0034] Spectral features of the magnetic field signal were extracted using Fast Fourier Transform;

[0035] Through the magnetic field vector B(x, y, z, t)=(B x By B z Spatial sampling is used to calculate the magnetic field gradient tensor;

[0036] Different magnetic field disturbances are classified using the obtained time-domain features, frequency-domain features, and spatial sampling results.

[0037] The method also includes distinguishing between geomagnetic interference and white noise. If it is geomagnetic interference, the geomagnetic field change in the test area is considered to be the same as that in the environmental monitoring area, and reverse compensation is performed directly by the magnitude and direction of the magnetic field obtained by magnetic sensing. If it is white noise, no compensation is performed, and the conventional Wiener filtering method is used for noise reduction.

[0038] This invention also provides a dynamic magnetic interference compensation system for aircraft. The system includes a magnetic field generating coil winding, a magnetic field measuring device, an environmental monitoring device, a magnetic field compensation coil winding, a control device, a load-bearing module, and a power supply device. The magnetic field generating coil winding is used to stably generate a uniform magnetic field at the microtesla level within a local space, and the magnetic field strength and direction are adjustable. The magnetic field measuring device consists of a multi-channel triaxial magnetic sensor array for real-time measurement of the magnetic field in the test area. The environmental monitoring device uses a triaxial sensor array positioned at a distance from the coil winding to monitor multi-source information such as geomagnetic changes and environmental interference magnetic fields, and identifies the type of interference in real time through a data processing unit.

[0039] Magnetic field compensation coil winding: It consists of coils in the X, Y, and Z directions, which can generate a magnetic field inside the coil that is equal in magnitude and opposite in direction to the external disturbance, so as to cancel out external geomagnetic fluctuations and environmental magnetic interference;

[0040] Control device: controls the entire magnetic measurement system;

[0041] Power supply unit: Provides power to the entire measurement system;

[0042] Load-bearing device: the equipment to be tested, so as to simulate interference compensation and test the results. The control device uses the method described above to determine the compensation magnetic field strength, and determines the output magnetic field of the magnetic field compensation coil winding based on the compensation magnetic field strength.

[0043] Beneficial effects

[0044] The magnetic field compensation method of this invention can better classify and locate magnetic sources, and then perform corresponding compensations by accurately identifying different similar interferences. The various compensations are then superimposed to achieve better interference compensation. This invention employs a static magnetic source location analysis algorithm with a closed-loop structure of "magnetic gradient tensor + extended Kalman filter (EKF)," eliminating the traditional particle swarm optimization and active excitation stages. This avoids initial value sensitivity issues and does not introduce external magnetic field interference, achieving millisecond-level static magnetic source location, meeting the requirements of magnetic property testing for fast, accurate, and disturbance-free compensation. The dynamic magnetic source location method introduces a time dimension based on the static magnetic source location method, continuously sampling magnetic signals at different times. The analysis and processing methods at each time point are the same as those in the static magnetic source location method. Preferably, this invention also adds a magnetic interference classification function. Through the magnetic field classification of this invention, interference magnetic fields can be dynamically compensated. Attached Figure Description

[0045] Figure 1 This is a diagram showing the composition of the electromagnetic field control and testing system for the aircraft of the present invention;

[0046] Figure 2 This is a flowchart of the magnetic interference compensation method of the present invention;

[0047] Figure 3 This is a schematic diagram showing the positional relationship of a simple cross-shaped measurement array.

[0048] Figure 4 Flowchart for Extended Kalman Filter (EKF) localization;

[0049] Figure 5 This is a vehicle movement trajectory diagram, which shows the change in distance of the vehicle relative to the test array;

[0050] Figure 6 A comparison chart of magnetic field fluctuations;

[0051] Figure 7 The figure shows the effect of magnetic moment estimation. It can be seen from the figure that the real-time estimation of the equivalent magnetic moment of the vehicle by EKF basically follows the true value, which verifies the effectiveness of the algorithm of the present invention.

[0052] Figure 8 The results of the magnetic field after compensation based on the estimated magnetic field components during the test show that the overall magnetic field after compensation is at a low level and has almost no fluctuations. Detailed Implementation

[0053] The present invention will be further described in detail below with reference to the embodiments and accompanying drawings, but the embodiments of the present invention are not limited thereto.

[0054] I. Overall Equipment Architecture

[0055] To address the problem that existing uniform magnetic field generating systems are susceptible to external disturbances in complex environments, this invention proposes a dynamic magnetic interference compensation method and a system for implementing the method.

[0056] like Figure 1 As shown, the dynamic magnetic interference compensation system for a spacecraft according to an embodiment of the present invention includes a magnetic field generating coil winding, a magnetic field measuring device, an environmental monitoring device, a magnetic field compensation coil winding, a control device, a load-bearing module, and a power supply device. The test site can be selected near the spacecraft or other areas to be tested.

[0057] The functions of each module in the magnetic measurement and compensation device are as follows:

[0058] (1) Magnetic field generating coil winding: used to stably generate a uniform magnetic field at the micro-Tesla level in a local space, and the magnetic field strength and direction are adjustable.

[0059] (2) Magnetic field measuring device: It consists of a multi-channel triaxial magnetic sensor array and is used to measure the magnetic field in the test area (inside the coil winding) in real time.

[0060] (3) Environmental monitoring device: A three-axis sensor array is set at a distance from the coil winding to monitor multi-source information such as geomagnetic changes and environmental interference magnetic fields, and the type of interference is identified in real time through the data processing unit. The environmental monitoring data serves as the input for the subsequent magnetic field compensation system.

[0061] (4) Magnetic field compensation coil winding: It consists of coils in the X, Y and Z directions, which can generate a magnetic field inside the coil that is equal in magnitude and opposite in direction to the external disturbance, thereby canceling the external geomagnetic fluctuations and environmental magnetic interference.

[0062] (5) Control device: controls the entire magnetic measurement system.

[0063] (6) Power supply: to supply power to the entire measurement system.

[0064] (7) Load-bearing device: to support the aircraft waiting for testing equipment in order to conduct interference compensation simulation and result testing.

[0065] II. The Interference Compensation Process

[0066] Figure 2 The diagram shows a flowchart of magnetic interference compensation. The required compensation current is calculated using the magnetic compensation method of this invention to drive the magnetic field compensation coil winding for dynamic magnetic compensation. First, a cross-shaped magnetic gradient tensor measurement four-probe array is deployed around the test area (see...). Figure 3Compared to other hexahedral and tetrahedral measurement arrays, the cross-shaped planar measurement array uses the fewest triaxial sensors and has the simplest structure. This reduces costs without compromising positioning accuracy, resulting in cost reduction and efficiency improvement. Four magnetic sensors are arranged in the horizontal and vertical directions of the cross-shaped measurement array, and the magnetic field gradient tensor is approximated through differential operations.

[0067] Secondly, by passively receiving magnetic anomaly signals and combining magnetic gradient tensor measurement with an extended Kalman filter algorithm, environmental interference is effectively suppressed, achieving high-precision positioning of the magnetic source. Furthermore, utilizing the five independent components of the magnetic gradient tensor, only four triaxial magnetic sensors forming a cross-shaped planar array are needed to complete tensor estimation within milliseconds, significantly reducing hardware complexity and power consumption.

[0068] The specific algorithm for magnetic gradient tensor localization is as follows:

[0069] (1)

[0070] The symbols in this formula have the following meanings: G is the magnetic gradient tensor; , , These are the magnetic field components of the magnetic field vector B in the x, y, and z directions in space. , , Each represents a component. Rates of change along the x, y, and z directions, respectively; , , Represents weight Rates of change along the x, y, and z directions, respectively; , , Represents weight Rates of change along the x, y, and z directions, respectively.

[0071] Since the magnetic field is a divergence-free field, its divergence at any point in space is 0, and it is an irrotational field in the absence of spatial current density, therefore we can obtain...

[0072] Therefore, we have: (2)

[0073] As shown above, the magnetic gradient tensor matrix actually contains only 5 independent elements. During the measurement process, only these 5 elements are needed to obtain the complete magnetic gradient tensor information.

[0074] Right now (3)

[0075] For cross-shaped measurement arrays, as follows: Figure 3 As shown, a schematic diagram of the magnetic gradient tensor system with a cross-shaped structure is drawn.

[0076] The system in the figure consists of four vector magnetic sensors, labeled 1, 2, 3, and 4 in counterclockwise order. Sensors 1 and 3 are located on the x-axis and pass through the center point o. The baseline distance between the center points of the two three-component magnetic sensors is d. Sensors 2 and 4 are located on the y-axis and also pass through the center point o. The baseline distance between the center points of the two vector magnetic sensors is d.

[0077] Assume the component output of the magnetic sensor The tensor matrix expression of the cross-shaped structure at the center point o, representing the magnitude of the magnetic field of sensor i in the j direction, is as follows:

[0078] (4)

[0079] Based on the above theory and array model, the magnetic tensor is measured:

[0080] Assuming the target object of the magnetic source is a magnetic dipole, the magnetic target object to be measured can be described by six quantities, namely the three-dimensional magnetic moment vector of the magnetic dipole. The position vector from measurement point P to the magnetic dipole Then the distance from the magnetic target point The magnetic field at that location is:

[0081]

[0082] (5)

[0083] In the above formula, , , The magnetic field strength of the magnetic target is represented by the components along the x, y, and z axes. Let be the magnetic permeability of the medium. The magnetic moment mode and direction cosine can be used to represent the magnetic moment components. Let the declination of the geomagnetic field at the magnetic dipole be A, the tilt angle be I, and the magnetic moment mode be... Therefore, the formula for calculating the direction cosine is: , (6)

[0084] The corresponding formula for calculating the magnetic moment component is:

[0085] , (7)

[0086] The theoretical values ​​of the five independent quantities of the magnetic gradient tensor matrix are expressed as follows:

[0087]

[0088]

[0089]

[0090]

[0091] (8)

[0092] In the above equation, the vacuum permeability , .

[0093] Next, let the distance from the field source in the magnetic dipole field be... The magnitude of the magnetic field strength at point A is Then the distance from the field source The magnetic field strength at point B can be expressed as:

[0094] , (9)

[0095] The difference in magnetic field strength between points A and B is:

[0096] , (10)

[0097] At the same time, the difference in magnetic field strength between the two points can also be expressed as

[0098] (11)

[0099] Therefore there is ,

[0100] Simplifying the above equation yields the position vector of the magnetic target relative to the measuring device:

[0101] .

[0102] The above is the calculation process of the single-point inversion positioning algorithm. The single-point positioning algorithm only needs to obtain the magnetic gradient tensor value of a point in the magnetic field to calculate the location information of the target. Its algorithm implementation is simple and the measurement is convenient. Therefore, the positioning of this system is based on the single-point positioning algorithm.

[0103] In this invention, the extended Kalman filter algorithm is implemented as follows: First, the three-dimensional spatial position of the magnetic source and its three-axis magnetic moments are combined to construct a six-dimensional state vector. During a single sampling period, it can be assumed that the magnetic source remains stationary; therefore, the state remains unchanged during the state prediction stage. Subsequently, twelve magnetic field data are acquired using four three-axis magnetic sensors deployed on a cross-shaped planar array. Differential operations are used to obtain the three-component magnetic field at the array center and five independent gradient tensor components, thus forming the observation vector. A nonlinear observation function is established based on the magnetic dipole model and linearized to the first order at the current state estimation point to obtain the observation matrix.

[0104] In the filtering process, the prior covariance (calculated from sensor magnetic field data) is first used to predict the covariance from the previous time step. Then, combining the observation matrix, process noise, and observation noise, the Kalman gain is calculated to weight and correct the residuals between the measured and predicted values, thereby updating the state vector and the posterior covariance. To improve robustness, this invention introduces the Sage-Husa adaptive estimation method to adjust the statistical characteristics of process noise and observation noise in real time, effectively offsetting model errors caused by environmental factors such as temperature drift and sudden changes in power grid harmonics. The above "prediction-update" loop is automatically executed after each magnetic field data acquisition, and can output the optimal estimate of the magnetic source position and magnetic moment within a millisecond time scale, providing a reliable benchmark for subsequent interference mapping and compensation current calculation.

[0105] The specific implementation steps of the algorithm are as follows:

[0106] I. Definition of State and System Model

[0107] 1. Six-dimensional state vector

[0108] Combine the three-dimensional spatial position of the magnetic source with the three-axis magnetic moment into a state vector:

[0109] ,in Let be the spatial coordinates of the magnetic source at time k. For the three-axis magnetic moment components of the magnetic source.

[0110] 2. System state transition

[0111] Since it is assumed that the magnetic source is stationary during a single sampling period, the state transition function is "state remains unchanged", that is: ,in It is the prior state estimate at time k (based only on (Time information) yes Posterior state estimate at time step (updated optimal state).

[0112] II. EKF Extended Kalman Filter Prediction Stage

[0113] 1. Prior state prediction

[0114] As per the assumption of "magnetic source at rest" mentioned above, the posterior state from the previous time step is directly inherited:

[0115]

[0116] 2. Initial covariance

[0117] Before the algorithm starts, a small amount of initial observation data (such as the initial magnetic field measurements of a cross array sensor) can be collected to statistically analyze and determine the initial covariance.

[0118] [1] Initial observations were collected: N sets (e.g., N=10) of initial magnetic field data Z1, Z2, ... Z N And estimate the initial state using the initial observation data. ;

[0119] [2] Calculate the variance of the estimated initial state (position + magnetic moment):

[0120] ,

[0121] ,

[0122] ,

[0123] ,

[0124] ,

[0125] (12)

[0126] [3] Combine the above variances into a 6×6 diagonal matrix:

[0127] , as the initial value of the prior covariance; (13)

[0128] [4] Note: It is sufficient to ensure that the values ​​are within a reasonable engineering range (absolute precision is not required). The covariance will be dynamically adjusted using the Sage-Husa adaptive algorithm based on actual observation data. Therefore, small deviations in the initial values ​​will not affect the long-term convergence of the algorithm.

[0129] 3. Prediction of prior error covariance

[0130] Prior error covariance Describes the "statistical characteristics of the deviation between the prior state and the actual state", which is derived from the posterior error covariance at the previous time step. process noise covariance Q k(Calculation reflecting the uncertainty of state transition) The process noise covariance Q at this point k This corresponds to the uncertainty in the prediction of the magnetic source state in this system, such as sudden changes in vehicle speed in the external environment, transient disturbances in the external environment, etc.

[0131] For the process noise covariance Q k The initial value needs to be obtained based on initial experimental statistics, by simulating interference sources (such as moving vehicles or energized coils) in the experimental environment, and recording their state fluctuations.

[0132] Position components: Collect estimated position values ​​of the magnetic source in the x / y / z directions (e.g., 1000 sets), and calculate the variance of all position estimates. For example, if a slight movement of the vehicle causes a position fluctuation of <0.01 m, then the position variance is set to... ;

[0133] Magnetic moment components: the magnetic source is collected at Magnetic moment estimates in different directions (e.g., 1000 sets, calculated from magnetic field data) are used, and the variance is calculated for all magnetic moment estimates. For example, power grid fluctuations may cause magnetic moment fluctuations <0.1. Then the variance of the magnetic moment is set as

[0134]

[0135] Therefore, we obtain the initial Q. k =diag(10 -4 , 10 -4 , 10 -4 (0.01, 0.01, 0.01). This can also be set based on industry experience and experimental experience; the adaptive algorithm will automatically correct the deviation subsequently.

[0136] That is, given an initial process noise, usually a very small value, it will be automatically corrected later.

[0137] III. Observation Model and Linearization

[0138] 1. Magnetic dipole observation function

[0139] Treating the magnetic source as a magnetic dipole, for any magnetic dipole, the magnetic dipole at a point in space... The resulting magnetic field strength (taking the x-direction as an example) is:

[0140] ,in It is the distance from the magnetic source to the observation point. , is the vacuum permeability; similarly The expression is: , .

[0141] Define the observation function as Given the state vector: ,have ,in To account for observation noise (sensor measurement errors, etc.), it must comply with... Observation noise covariance R k The measurement errors corresponding to the magnetic sensor array, including temperature drift and electronic noise, are used to quantify the statistical characteristics of noise during the measurement process. These can be estimated using long-term measurement data of the sensor array under stable conditions. k It is used to correct sensor measurement results and reduce the measurement error of the sensor itself. It reduces the impact of sensor error by adding observation noise.

[0142] make Based on the derivation of the magnetic field formula for a magnetic dipole, we have:

[0143] ;

[0144] By taking the spatial partial derivative of the magnetic field component (with the partial derivative applied to the location of the magnetic source), we obtain:

[0145]

[0146] Integrating the 8 components yields the following observation function:

[0147] .

[0148] Furthermore, the measured value vector Z at the center k Includes: 3 magnetic field components (Magnetic field generated by the magnetic source, selected from magnetic field data measured by a specific sensor, such as magnetic sensor 1); 5 independent gradient tensor components Therefore, an eight-dimensional measurement vector is formed:

[0149] 2. Jacobian matrix (observation matrix) )

[0150] EKF requires the above observation functions In the prior state Performing a first-order Taylor expansion, the Jacobian matrix becomes the linearized observation matrix: .

[0151] For example, to calculate gradient components For example, the Jacobian element is:

[0152] Similarly, it can be calculated China The partial derivatives are used to obtain an 8×6 Jacobian matrix (8-dimensional observation vector and 6-dimensional state vector).

[0153] IV. EKF Update Phase (Posterior Estimation)

[0154] 1. Kalman gain calculation

[0155] Using the obtained observation matrix H k Substitute into the formula The Kalman gain was calculated. This is used to balance the "reliability of prior estimation" and the "reliability of observation," and corresponds to the weighting coefficient matrix of "predicted state" and "sensor measured value" applied by the control module when calculating the compensation current. When the sensor measurement environment has high noise, K... k Automatically reduce the weight of observation data in compensation calculations to avoid noise amplification; when the environment is relatively stable, K k The weights of the observation data are increased to quickly track the position of the magnetic source and changes in the magnetic moment. Through dynamic adjustment of this gain matrix, the compensation current is ensured to respond promptly to sudden disturbances while remaining stable during long-term operation.

[0156] 2. Status Update

[0157] The prior state is corrected using the observation residuals (the difference between measured and predicted observations) to obtain the posterior state:

[0158] ,in It is the difference between the "actual observed value" and the "observation value predicted based on the prior state" (observation residual).

[0159] 3. Error Covariance Update

[0160] Updated posterior error covariance: , where I is the identity matrix.

[0161] V. Sage-Husa Adaptive Noise Estimation (Robustness Improvement)

[0162] To address non-modeling errors such as temperature drift, the Sage-Husa method is introduced to adjust the process noise covariance in real time. (Reflecting the uncertainty of state transitions) and observation noise covariance The two parameters mentioned above are inherent core parameters of the Kalman filter itself.

[0163] 1. Definition of Cumulative Factor and Residual

[0164] Time averaging factor The α value is used to balance the weights of new and old data. It is a balance parameter between real-time data and historical noise statistics in the control module: a large α value indicates greater emphasis on historical noise (stability); a small α value indicates greater emphasis on new data (fast response).

[0165] Cumulative factor (initial =0), which is a weighting parameter used internally by the control module to store and update the noise statistical characteristics. Its function is to gradually accumulate the statistical characteristics of the magnetic field measurement in each sampling period, and to correct the zero-point drift and system bias of the magnetic sensor over a long period of time.

[0166] It determines the degree to which "system historical noise information" is retained in subsequent noise estimation: the larger the cumulative factor, the stronger the system's memory of long-term stable noise, which is suitable for compensating for slowly changing background noise; the smaller the cumulative factor, the more sensitive the system is to newly sampled data, which is suitable for quickly responding to transient interference.

[0167] Observation residuals Its essence is "measured value - predicted observation value", reflecting the difference between the sensor's measured magnetic field data and the model's predicted value, and reflecting the deviation between the prediction model and the actual situation. The optimal state vector is obtained when the observation residual meets the threshold condition.

[0168] State residual Its essence is "the state vector of the next moment - the state vector of the previous moment", which corresponds to the difference between the predicted and estimated values ​​of the interfering magnetic source at the next moment.

[0169] 2. Real-time estimation of noise covariance

[0170] Observation noise covariance : (14)

[0171] Process noise covariance (15), of which, This is the cumulative factor from the previous time step. This is the cumulative factor for the next time step.

[0172] 3. Iterative update of noise covariance

[0173] The estimated , Replace the original fixed (original An initial value can be given (and then corrected using the algorithm above), and it participates in the EKF operation at the next time step to achieve adaptive adjustment of noise statistical characteristics.

[0174] VI. Loop Execution and Time Characteristics

[0175] The above-described process of "prediction → observation linearization → Kalman gain calculation → state and covariance update → Sage-Husa noise adaptation → prediction" is executed automatically after each magnetic field data acquisition (sampling period is in the millisecond range, such as 1ms), and finally outputs the optimal state estimate of the magnetic source. (Location +Magnetic moment This provides a reliable benchmark for subsequent mapping of the interference magnetic field and calculation of the compensation current.

[0176] The compensation current can be calculated in real time based on the magnetic source positioning results.

[0177] Let the magnetic source moment of the extended Kalman filter output be... The position vector of the magnetic source relative to the center of the compensation coil is The reverse magnetic field that the compensation coil needs to generate at the center of the test area is:

[0178] (15)

[0179] For a radius of R c The number of turns is N c A circular compensating coil, the magnetic field and current I on its central axis c The relationship is

[0180] (16)

[0181] The compensation current can be obtained by combining equations (17) and (18).

[0182] (17)

[0183] The controller refreshes the DAC output in real time according to formula (19) and drives the compensation coil to generate magnetic fields with equal amplitude and opposite direction.

[0184] Experimental Test

[0185] During the magnetic characteristic testing of the aircraft's magnetic field compensation, the test site was selected in an open experimental area far away from strong electromagnetic interference sources. However, even so, it was difficult to avoid magnetic field interference. Therefore, an anti-interference support platform was laid on the ground to support the aircraft under test and its supporting test equipment. During the test, a dynamic magnetic interference compensation system was activated.

[0186] After the test begins, the system first establishes a uniform reference magnetic field with an intensity of 0.5 μT using the magnetic field generation module. This magnetic field serves as the reference background field for the test environment. The control unit monitors the field in real time to ensure that the reference magnetic field fluctuates within a stable range of ±0.1 nT, thus providing reliable initial conditions for the measurement of the vehicle's magnetic properties. At this point, the vehicle is positioned at the center of the test area, and its magnetic properties (including the magnitude and direction of the magnetic moment, magnetic anomaly distribution, etc.) are continuously collected by the magnetic field sensor array and uploaded to the data processing unit, enabling precise measurement of the vehicle's magnetic field leakage characteristics.

[0187] During the test, a vehicle slowly drove from north to south near the test site. Its engine current and metal casing constantly changed relative positions during the movement, introducing time-varying local magnetic interference signals into the test area. The environmental monitoring device detected abnormal fluctuations in the magnetic field background, with peak amplitudes reaching ±15 nT and exhibiting significant time-varying characteristics, indicating that the interference was caused by mechanical activity. To further pinpoint the source of the interference, this invention utilizes a cross-shaped magnetic gradient tensor sensor array to collect instantaneous magnetic fields, calculates the magnetic field gradient tensor, and then applies the Extended Kalman Filter (EKF) algorithm. Using the absolute value of the magnetic field from the magnetic sensor and the magnetic field gradient data as input, it dynamically estimates the equivalent magnetic source and position of the vehicle within a millisecond timescale. The EKF algorithm outputs the equivalent magnetic moment magnitude and spatial coordinates of the vehicle in each sampling period. The results show that the vehicle gradually approached from approximately 150 m from the array center to 100 m and then gradually moved away. Throughout the movement, the magnetic source intensity fluctuated slightly within the range of 99~102 A·m². Figure 5 As shown.

[0188] Based on the EKF output (according to the algorithm mentioned above), the environmental interference compensation module calculates the required reverse current for the compensation coil to be 50~70 mA. The controller refreshes the output current within 5 ms, driving the compensation coil to generate an equivalent reverse magnetic field in the opposite direction within the test area. After compensation, the magnetic field fluctuation in the test area rapidly decreases to ±0.5 nT, restoring a stable level suitable for testing the magnetic properties of targets (such as aircraft).

[0189] Comparative experiments show that, without compensation, the passage of a vehicle can cause magnetic field drift of up to ±15 nT, and the repeatability error of the magnetic anomaly test exceeds 12%. After adopting the EKF positioning and compensation scheme proposed in this invention, the magnetic field fluctuation is controlled within ±1 nT, and the test accuracy is improved by an order of magnitude.

[0190] This embodiment provides an effective method for compensating for dynamic magnetic interference.

[0191] Example 2

[0192] In another preferred implementation, when multiple interferences exist, the present invention adopts a method of first classifying and then classifying and compensating. This embodiment is a further enhancement based on embodiment 1. In addition to the structure and method of embodiment 1, the following method for distinguishing magnetic interference types is added. The measurement sensors used below are additional sensors.

[0193] Specifically, the steps in this embodiment include: first, classifying the magnetic field interference, and then compensating for different types of interference.

[0194] The classification process includes:

[0195] Step (1) Arrange N triaxial magnetic sensors at equal intervals in three dimensions at different locations in the target area (for example, in a cubic matrix, arranged symmetrically with one vertex as the center). Each magnetic sensor measures the magnetic field at its location to obtain the corresponding magnetic sensor signal.

[0196] (1)

[0197] Step (2): Preprocess the obtained magnetic sensor signal, including DC drift removal, bandpass filtering, and normalization.

[0198] (2)

[0199] in Let i be the time average value of sensor i. This indicates a bandpass filter that filters out high-frequency noise and low-frequency drift, while retaining a frequency range that covers geomagnetic fluctuations and power frequency interference (e.g., 0.01~200 Hz); Step (3): Calculate the time-domain characteristic quantity of the three-component magnetic field values ​​measured by each sensor using the following formula: the average value of the magnetic field within the measurement time period T: , (3)

[0200] The variance of the magnetic field within the measurement time period T: (4)

[0201] Kurokiness of the magnetic field within the measurement time period T: (5) Among them To measure the average value of the magnetic field within a time period T, To measure the variance of the magnetic field within time period T, To measure the kurtosis of the magnetic field within a time period T, where i is the sensor number, i=1,…,N, The preprocessed magnetic field signal is represented by α, which represents the x, y, and z axes, t is the time variable, and t0 is the start time of the measurement period. Step (4): Extract the spectral features of the magnetic field signal using Fast Fourier Transform.

[0202] Extract the main frequency amplitude of the magnetic field signal and phase n represents the index of the time-domain sampling point, from 0 to M-1; M represents the total number of sampling points; f represents the sampling interval. k Represents the k-th discrete frequency; Step (5): Calculate the magnetic field gradient tensor by spatial sampling of the magnetic field vector B(x, y, z, t) = (Bx, By, Bz).

[0203] (7)

[0204] And further calculate its norm. ;

[0205] Step (6) Use the obtained time-domain features, frequency-domain features and spatial sampling results to classify different magnetic field interferences.

[0206] Interference magnetic fields are classified and identified as follows: Determine if both the variance and kurtosis of the interference magnetic field are greater than 3. If both are greater than 3, and ||G|| is greater than 0.15, it can be identified as a near-field disturbance. If neither variance nor kurtosis reaches 3 simultaneously, then the magnitude of ||G|| is used for identification. If ||G|| is greater than 0.6, the interference magnetic field is identified as a near-field disturbance. If the magnetic field signal does not meet the above requirements, the interference magnetic field is identified as a far-field disturbance. Here, near-field and far-field disturbances mainly refer to sudden or mobile interference sources. Extract the main frequency amplitudes of the interference magnetic field signal. and phase This allows us to obtain the periodic characteristics of power / grid interference. When the difference between the main frequency and 50 or 60 Hz is less than the threshold, it is identified as a power / grid interference signal.

[0207] Calculate the autocorrelation length and variance stability of the interfering magnetic field signal. If the autocorrelation length is much greater than 10s, the spectral slope is significantly less than zero, and the variance stability is not zero, then the interfering magnetic field is determined to be geomagnetic interference. If the autocorrelation length is less than 1s and the spectral slope is significantly less than zero, then the interfering magnetic field is determined to be geomagnetic interference. If the variance stability is approximately 0, then the interference signal is determined to be white noise.

[0208] For far-field interference and geomagnetic field interference, a reverse magnetic field of equal intensity is used for compensation. The applicant found that white noise has an approximately flat distribution in the frequency domain power spectrum, i.e., S(f)≈C0, where C0 is a constant and the spectral slope is... This reflects that this type of signal has a uniform energy distribution across all frequency components and lacks a specific dominant frequency component. Geomagnetic disturbances, on the other hand, are concentrated in the low-frequency range and typically manifest as… , where α 1~2 demonstrates the signal's characteristic of concentrated energy in the low-frequency range and rapid attenuation with increasing frequency, while the spectral slope... The value is significantly less than zero, indicating strong time correlation and non-stationary characteristics. Furthermore, from a time-domain statistical perspective, the autocorrelation function R(τ) of white noise approximates the Dirac function form: , where σ 2 The noise variance indicates that the time samples are completely uncorrelated. When τ≠0, R(τ) decays rapidly to zero, and the autocorrelation length Lc is extremely small, usually much smaller than 1 second, failing to form a smooth temporal change; while geomagnetic disturbances have a longer correlation, and their autocorrelation function usually exhibits exponential or power-law decay. Its Lc can reach tens of seconds, significantly greater than that of white noise, demonstrating the continuity and slow variation of geomagnetic disturbances. For a discrete signal x[n], the sample autocorrelation function is defined as: Where N is the number of sampling points, k is the lag order, and R xx [0] represents the energy of the signal, and the autocorrelation length Lc refers to the autocorrelation function decaying to zero or less than a certain threshold (such as R). xx [k]<0.05R xx [0]) corresponds to the lag step size. Furthermore, statistically, white noise approximately follows a Gaussian distribution within a short time window, with its mean stable near zero and variance remaining essentially constant; while geomagnetic disturbances exhibit non-stationary characteristics, with their mean slowly drifting over time and variance σ 2 Significant changes exist on long-term timescales, which can be represented by the relative rate of change. Characterization is performed. For a signal x[n], the variance is defined as: , The variance here is different from the variance used in time-domain characteristic analysis. (Used to roughly distinguish between steady-state disturbances and sudden disturbances), specifically for finely distinguishing between white noise and geomagnetic disturbances.

[0209] Near-field and far-field perturbations are distinguished by the gradient norm ||G||, and further analyzed by the spectral slope κ, autocorrelation length Lc, and variance stability Δσ. 2 The combined criteria can effectively distinguish between white noise and geomagnetic disturbances, ensuring the accuracy and robustness of the judgment.

[0210] In summary, the focus of this embodiment is on the dynamic magnetic compensation system, and a three-level progressive technical solution of "interference classification - magnetic source localization - dynamic compensation" is proposed.

[0211] First, classify the types of magnetic field interference to determine whether they are: near-field interference, far-field interference, power grid interference, geomagnetic interference, or white noise.

[0212] If the interference type is near-field interference (generally caused by moving objects or sudden events), magnetic source localization and inversion are performed; for far-field interference, the measured magnetic field is directly used for reverse compensation. If it is geomagnetic interference, the geomagnetic field change in the test area can be considered the same as that in the environmental monitoring area, and reverse compensation can be directly performed using the magnetic field measured by magnetic sensing data. For power grid magnetic field interference, the magnetic field data of the peak area of ​​the power grid interference is brought into step (2) of the method described in claim 1 to locate the magnetic source of the power grid interference. Based on the located magnetic source position and magnetic source intensity of the power grid interference, the formula is used to perform the inversion. The peak compensation magnetic field strength at the target location is calculated. The peak compensation magnetic field strength with the same period is taken as the maximum value, and the reverse magnetic field compensation is performed. For the geomagnetic interference magnetic field, the reverse magnetic field value of the geomagnetic interference magnetic field is directly used for compensation.

[0213] In this embodiment, the accuracy of noise category identification is significantly improved compared to existing technologies. For example, in a mixed environment simulating 50 nT of geomagnetic disturbance and 10 dB of white noise signal-to-noise ratio, the traditional single-spectrum criterion method achieves a classification accuracy of only about 72% and a false positive rate exceeding 28%. However, the joint discrimination mechanism of this invention (spatial gradient norm + spectral features + autocorrelation length) improves the classification accuracy to 92% under the same conditions, while reducing the false positive rate to less than 8%. Furthermore, in a scenario simulating the coexistence of near-field magnetic metal disturbance (5 cm distance, intensity approximately 200 nT) and power line interference (50 Hz fundamental wave and 150 Hz subharmonics), traditional methods often confuse metal disturbance with low-frequency harmonics, resulting in an accuracy of less than 65%. However, this invention, by introducing a dual-threshold confirmation mechanism, can simultaneously distinguish between two types of interference sources, achieving a discrimination accuracy of over 95%. These comparative results demonstrate that the present invention not only achieves higher accuracy under single interference but also maintains stable discrimination performance in complex environments with multiple interference sources, exhibiting outstanding engineering application value. Furthermore, by classifying more accurately, compensation for complex magnetic fields can also be effectively improved.

[0214] Although the principles of the present invention have been described in detail above with reference to preferred embodiments, those skilled in the art should understand that the above embodiments are merely illustrative explanations of the implementation of the present invention and are not intended to limit the scope of the present invention. The details in the embodiments do not constitute a limitation on the scope of the present invention. Any obvious changes, such as equivalent transformations or simple substitutions, based on the technical solutions of the present invention without departing from the spirit and scope of the present invention fall within the protection scope of the present invention.

Claims

1. A dynamic magnetic interference compensation method for testing the magnetic characteristics of aircraft, characterized in that, The method includes: Step (1) Select the target site as the test area; Step (2) Set up an environmental monitoring device around the target area and locate the magnetic source based on the measurement results of the environmental monitoring device; Step (2) includes: (2.1) Constructing a six-dimensional state vector, which includes the three-dimensional spatial position and three-axis magnetic moments of the magnetic source; (2.2) Establishing the state transition system and observation model of the magnetic source, and constructing the observation function of the magnetic source. ,in To observe the noise, X k It is a six-dimensional state vector, that is ,in, Let be the spatial coordinates of the magnetic source at time k. (2.3) Linearize the observation model to obtain the linearized observation matrix of the magnetic source: ; (2.4) Use Kalman gain calculation to perform error covariance update iteration and state vector update of observations; (2.5) Estimate and iteratively update the observation noise covariance and process noise covariance; (2.6) Repeat steps (2.1)-(2.5) above to output the optimal state estimate of the magnetic source. This state estimate includes location. and magnetic moment Step (3): Based on the location and intensity of the magnetic source after positioning and the location of the target area, determine the magnetic field generated by the magnetic source in the target area; Step (4): Based on the calculated magnetic field generated by the magnetic source at the target location, generate a reverse magnetic field at the target location.

2. The dynamic magnetic interference compensation method for aircraft according to claim 1, characterized in that, The method includes: Magnetic field measurements are performed using a cross-shaped planar measurement array, which includes four or 4N vector magnetic sensors. The tensor matrix expression of the cross-shaped planar measurement array at the center point o is as follows: , among which, (B) 1x B 1y B 1z ,) represents the triaxial magnetic field values ​​measured by the vector magnetic sensor on the first cross branch of the cross-shaped planar measurement array, (B) 2x B 2y B 2z ,) represents the triaxial magnetic field values ​​measured by the vector magnetic sensor on the second cross branch of the cross-shaped planar measurement array, (B) 3x B 3y B 3z ,) represents the triaxial magnetic field values ​​measured by the vector magnetic sensor on the third cross branch of the cross-shaped planar measurement array, (B) 4x B 4y B 4z ,) represents the three-axis magnetic field values ​​measured by the vector magnetic sensor on the fourth cross branch of the cross-shaped measurement array.

3. The dynamic magnetic interference compensation method for aircraft according to claim 1, characterized in that, The method includes: collecting several sets of initial observation data to statistically analyze and determine the initial covariance. This step includes: Step (2.21), collecting N sets of initial magnetic field data Z1, Z2, ... Z N And estimate the initial state using the initial observation data. ; Step (2.22): Calculate the variance of the estimated initial state: , , , , , ； x i ,y i ,z i Let m represent the coordinates of the magnetic source measured in the i-th group. xi m yi , m zi Representing the magnetic moment of the magnetic source, step (2.23) combines the above variances into a 6×6 diagonal matrix: 。 4. The dynamic magnetic interference compensation method for aircraft according to claim 3, characterized in that, The method further includes: collecting initial observation data of the magnetic field to determine the initial covariance of the magnetic field, which serves as the initial value of the prior covariance; and for subsequent time steps, using the posterior covariance of the previous time step. process noise covariance Q k Perform the prior error covariance at the next time step. predict, Q k The process noise covariance is obtained by collecting magnetic field information from the magnetic source and calculating the variance of the estimated magnetic moment.

5. The dynamic magnetic interference compensation method for aircraft according to claim 4, characterized in that, The method further includes: utilizing the obtained observation matrix H k Substitute into the formula The Kalman gain was calculated. , To observe the noise covariance, the observation residuals and Kalman gain are used. To correct the prior state, we obtain the posterior state: ,in It is the difference between the "actual observed value" and the "observation value predicted based on the prior state"; Updated posterior error covariance: , where I is the identity matrix.

6. The dynamic magnetic interference compensation method for aircraft according to claim 5, characterized in that, Step (4) includes: based on the formula To calculate the compensation magnetic field strength at the target location, The position vector of the magnetic source relative to the center of the compensation coil. To extend the magnetic source moment of the Kalman filter output, is the vacuum permeability.

7. The dynamic magnetic interference compensation method for aircraft according to claim 6, characterized in that, The observation function is: , Where, m x m y m z , , and z are the magnetic moments of the magnetic source in the x, y, and z directions, respectively, r is the distance from the target point to the magnetic source, and x, y, and z are the coordinates of the target point.

8. The dynamic magnetic interference compensation method for aircraft according to claim 1, characterized in that, The method further includes: classifying the magnetic field interference and compensating for different types of interference. This step includes classifying the interfering magnetic field into: near-field interference, far-field interference, and power grid interference. For near-field interference, step (2) of the method in claim 1 is used to locate the magnetic source and determine the location of the magnetic source based on the formula. To calculate the compensation magnetic field strength at the target location, The position vector of the magnetic source relative to the center of the compensation coil. To extend the magnetic source moment of the Kalman filter output, The magnetic permeability is given by the vacuum magnetic field. For far-field interference, compensation is made using the reverse magnetic field strength based on the magnetic field measured by the magnetic sensor located in the target area. For grid interference, the magnetic field data of the peak area of ​​the grid interference is brought into step (2) of the method described in claim 1 according to the measured grid cycle to locate the magnetic source of the grid interference. Based on the located magnetic source position and magnetic source strength of the grid interference, the magnetic source is located using the formula... The peak compensation magnetic field strength at the target location is calculated, and the peak compensation magnetic field strength with the same period is taken as the maximum value for reverse magnetic field compensation.

9. The dynamic magnetic interference compensation method for aircraft according to claim 1, characterized in that, The process of determining the magnetic interference category includes: arranging N triaxial magnetic sensors at equal intervals in three dimensions at different locations within the target area; each magnetic sensor measures the magnetic field at its location to obtain the corresponding magnetic sensor signal. （1） The obtained magnetic sensor signal is preprocessed; the variance of the magnetic field within the measurement time period T is measured. Measure the kurtosis of the magnetic field within a time period T; Spectral features of the magnetic field signal were extracted using Fast Fourier Transform; The magnetic field gradient tensor is calculated by spatial sampling of the magnetic field vector B(x, y, z, t) = (Bx, By, Bz). Different types of magnetic field interference are classified using the obtained time-domain features, frequency-domain features, and spatial sampling results. The method also includes distinguishing between geomagnetic interference and white noise. If it is geomagnetic interference, the geomagnetic field change in the test area is considered to be the same as that in the environmental monitoring area, and reverse compensation is performed directly by the magnitude and direction of the magnetic field obtained by magnetic sensing. If it is white noise, no compensation is performed, and the conventional Wiener filtering method is used for noise reduction.

10. A dynamic magnetic interference compensation system for an aircraft, characterized in that, The system includes a magnetic field generating coil winding, a magnetic field measuring device, an environmental monitoring device, a magnetic field compensation coil winding, a control device, a load-bearing module, and a power supply device. The magnetic field generating coil winding is used to stably generate a uniform magnetic field at the microtesla level within a local space; the magnetic field strength and direction are adjustable. The magnetic field measuring device consists of a multi-channel triaxial magnetic sensor array used to measure the magnetic field of the test area in real time. The environmental monitoring device uses a triaxial sensor array positioned at a distance from the coil winding to monitor multi-source information such as geomagnetic changes and environmental interference magnetic fields, and identifies the type of interference in real time through a data processing unit. Magnetic field compensation coil winding: It consists of coils in the X, Y, and Z directions, which can generate a magnetic field inside the coil that is equal in magnitude and opposite in direction to the external disturbance, so as to cancel out external geomagnetic fluctuations and environmental magnetic interference; Control device: controls the entire magnetic measurement system; Power supply unit: Provides power to the entire measurement system; Load-bearing device: the device to be tested, for simulating interference compensation and testing the results. The control device uses the method described in any one of claims 1-9 to determine the compensation magnetic field strength, and determines the output magnetic field of the magnetic field compensation coil winding based on the compensation magnetic field strength.

Images