A method and device for detecting interference of a magnetometer based on multi-feature fusion

Abstract

The application relates to the technical field of magnetometer interference detection, and discloses a magnetometer interference detection method and device based on multi-feature fusion, which comprises the following steps: collecting three-axis magnetometer data and constructing a sampling matrix to calculate static baseline statistics and a static threshold value by using initial static samples of the sampling matrix; sequentially calculating a magnetic field module value, a magnetic inclination angle, a magnetic azimuth angle and corresponding magnetometer feature instantaneous increments according to the sampling matrix; judging an anomaly by using a candidate indicator function, and weighting and fusing an anomaly proportion into a comprehensive detection score by using a cause-effect confirmation window, and then calculating a center estimate and a deviation, combining the static threshold value to construct a main threshold value and upper and lower hysteresis threshold values; and comparing the comprehensive detection score with the main threshold value and the upper and lower hysteresis threshold values and historical results to output a magnetic interference judgment result. The application can consider adaptive statistical discrimination ability and avoid significant time delay, can accurately perform magnetometer interference detection, and has the characteristics of high elimination precision and low method time delay.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of magnetometer interference detection, and particularly relates to a magnetometer interference detection method and device based on multi-feature fusion. BACKGROUND

[0002] As an important part of inertial / attitude solution system, magnetometers are widely used in unmanned aerial vehicle navigation, wearable motion sensing, vehicle heading estimation, robot positioning and mobile terminal attitude measurement. By measuring the geomagnetic field vector, it can provide heading angle and attitude reference in the absence of satellite signals or indoor environment. However, in the actual application process, the magnetometer is often affected by the complex magnetic environment, such as electromagnetic interference generated by surrounding motors, power lines, metal structures and adjacent electronic equipment. These interferences will cause the measured value of the magnetic field to change suddenly, distort or shift, resulting in attitude estimation drift, heading angle jump and even navigation system failure. Therefore, in order to effectively use the magnetic field data, the original data of the magnetometer needs to be further processed by an effective magnetometer interference detection method.

[0003] The existing magnetometer interference detection method mainly relies on the principle of multi-magnetometer differential discrimination, that is, the difference in magnetic field distribution is sensed by multiple magnetometers at the same position to identify the interference source. This method can suppress the influence of local magnetic distortion to a certain extent, but inevitably increases the hardware cost and system complexity. For example, the output difference of multiple magnetometers is used to calculate the target magnetic field strength to determine whether there is an interference magnetic field in the environment; the difference between the magnetic field information measured by multiple magnetometers is used for magnetic interference detection. Although the above method can improve the detection accuracy, additional hardware synchronization, spatial calibration and signal matching process are required in actual engineering application, resulting in complex system implementation and insufficient real-time performance.

[0004] In summary, in the engineering scene with high real-time response requirements, a multi-feature fusion adaptive detection method based on single magnetometer signal is still needed to accurately identify magnetic interference without additional hardware, by statistical features and dynamic threshold discrimination, so as to balance the detection accuracy and system real-time performance. SUMMARY

[0005] The purpose of the present application is to solve the problem of insufficient detection accuracy of the existing magnetometer interference detection method, and to provide a magnetometer interference detection method and device based on multi-feature fusion.

[0006] The present application provides a magnetometer interference detection method based on multi-feature fusion, comprising:

[0007] S1, collecting three-axis magnetometer data and constructing a sampling matrix M, calculating a static baseline statistic and a static threshold value using the initial static sample of the sampling matrix M ; wherein the static baseline statistics include a median and a median absolute deviation ;

[0008] S2, sequentially calculating magnetic field module values, magnetic inclination angles, magnetic azimuth angles and corresponding magnetometer feature instantaneous increments according to the sampling matrix M ;

[0009] S3, judging abnormalities through a candidate indicator function according to the magnetometer feature instantaneous increments, and statistically calculating an abnormality proportion through a causal confirmation window, and weighting and fusing the abnormality proportion into a comprehensive detection score ;

[0010] S4, calculating a central estimate of the comprehensive detection score according to the comprehensive detection score , a deviation , and constructing a main threshold value and upper and lower hysteresis threshold values in combination with a static threshold value ;

[0011] S5, outputting a magnetic interference judgment result through comparison of the comprehensive detection score with the main threshold value and the upper and lower hysteresis threshold values and historical results .

[0012] Further, the S1 specifically includes:

[0013] S101, collecting three-axis magnetometer data in real time, and constructing the three-axis magnetometer data into a sampling matrix in time sequence, denoted as , is a three-axis measurement matrix, is an original measurement value at time i in channel c, and N is a sampling point number;

[0014] S102, calculating static baseline statistics including a median and a median absolute deviation using initial static samples of the sampling matrix M, and the calculation formula is:

[0015] ;

[0016] wherein, is a median of initial static interval data; is a sample for static baseline estimation, and j is a time index within a window; is an initial window size; is a median absolute deviation of initial static interval data;

[0017] ​​​​S103、According to the median and the median absolute deviation Calculate the static threshold , the formula is:

[0018] ;

[0019] Wherein, is the static threshold, is the median of the initial static interval data, is the median absolute deviation of the initial static interval data, is the static threshold coefficient.

[0020] Further, the S2 specifically comprises:

[0021] S201, according to the sampling matrix M, the magnetic field modulus value is calculated in turn , the formula is:

[0022] ;

[0023] Wherein, is the three-axis magnetometer modulus value at time k , when i=k,

[0024] ;

[0025] is the magnetometer X axis output at time k , is the magnetometer Y axis output at time , is the magnetometer Z axis output at time k ; k S202, according to the magnetic field modulus value at time

[0026] , the magnetic inclination k at time , the magnetic course angle k is calculated, the formula is: ;

[0027] S203, according to the magnetic field modulus value , the magnetic inclination

[0028] , the magnetic course angle , the corresponding magnetometer characteristic instantaneous increment is calculated, the formula is: ;

[0029] Wherein,

[0030] ​​For at any time k The magnitude increment of the triaxial magnetometer For at any time k The magnetic tilt increment, For at any time k The magnetic flight angle increment, For a moment k-1 The magnetic field magnitude, For a moment k-1 Magnetic tilt angle, For a moment k-1 Magnetic flight angle.

[0031] Furthermore, S3 specifically includes:

[0032] S301. Anomalies are determined based on the instantaneous increment of the magnetometer characteristic using a candidate indication function; wherein the expression of the candidate indication function is:

[0033] ;

[0034] in, Here, j is the indicator function; j represents the historical time within the causal confirmation window, indicating the sub-time within the causal confirmation window at time k. For at any time j The magnitude increment of the triaxial magnetometer For at any time j The magnetic tilt increment, For at any time j The magnetic flight angle increment, The multiple factor for modulus difference discrimination. Let the standard deviation of the modulus over the causal local window be denoted as . The threshold for judging magnetic inclination angle is... The threshold for determining the magnetic flight angle;

[0035] S302. Based on the causal confirmation window, the abnormal occurrence ratio of the instantaneous increments of the three magnetometer characteristics is statistically analyzed, and the formula is as follows:

[0036] ;

[0037] in, , , For a moment k The proportion of abnormal occurrences of instantaneous increments in the three types of magnetometer characteristics. L For the length of the causal confirmation window, For indicator functions;

[0038] S303. The abnormal occurrence ratios of the instantaneous increments of the three magnetometer features are linearly fused according to preset weights to obtain a comprehensive detection score for magnetometer anomaly detection. Its formula is:

[0039] ;

[0040] wherein, is the integrated detection score at time k , , , is a coefficient for linearly weighting the three ratios to obtain the integrated score, , , is the proportion of the three magnetometer feature instantaneous incremental anomalies appearing at time k .

[0041] Further, the S4 specifically comprises:

[0042] S401, according to the integrated detection score calculating the center estimate of the integrated detection score , the calculation formula is:

[0043] ;

[0044] wherein, is the center trend estimate obtained by exponentially weighted average of the integrated detection score at time k , is the center trend estimate obtained by exponentially weighted average of the integrated detection score at time k-1 , is a weight coefficient, is the integrated detection score at time k ;

[0045] S402, according to the center estimate of the integrated detection score calculating the deviation , the calculation formula is:

[0046] ;

[0047] wherein, is the exponentially weighted average of the absolute deviation at time k, is the exponentially weighted average of the absolute deviation at time k-1, is a weight coefficient, is the integrated detection score at time k , is the center trend estimate obtained by exponentially weighted average of the integrated detection score at time k , is an approximate value for comparison with the standard deviation under normal distribution;

[0048] S403、According to the center estimation of the comprehensive detection score With bias , combined with static threshold The main threshold is constructed and the static lower limit protection is applied, and the formula is:

[0049] ;

[0050] Wherein, The main threshold at time k , The main threshold scale term is multiplied by a factor;

[0051] S404、According to the main threshold The upper and lower hysteresis threshold values are derived, and the formula is:

[0052] ;

[0053] Wherein, The entry threshold at time k , The exit threshold at time k , The upper boundary coefficient, The lower boundary coefficient.

[0054] Further, in the S5, the expression of magnetometer anomaly determination is:

[0055] ;

[0056] Wherein, The magnetic interference determination result at time k , The magnetic interference determination result at time k-1 , The comprehensive detection score at time k .

[0057] The application also provides a magnetometer interference detection device based on multi-feature fusion, based on the magnetometer interference detection method based on multi-feature fusion as described above, the device comprises:

[0058] The acquisition module is used for acquiring three-axis magnetometer data and constructing a sampling matrix M, and calculating static baseline statistics and static threshold Using the initial static sample of the sampling matrix M; wherein the static baseline statistics include median And the absolute deviation of median ;

[0059] The calculation module is used for calculating the magnetic field module value , magnetic inclination , magnetic course , and corresponding magnetometer feature instantaneous increment ;

[0060] Fusion module for determining abnormalities by candidate indicator function according to the magnetometer feature instantaneous increment, and statistically determining the proportion of abnormalities by the causal confirmation window, and weighting the proportion of abnormalities into a comprehensive detection score ;

[0061] Construction module for constructing a main threshold value and upper and lower hysteresis threshold values according to the comprehensive detection score , and combining a static threshold value ; , and combining a static threshold value ;

[0062] Output module for outputting a magnetic interference judgment result by comparing the comprehensive detection score , and the main threshold value and upper and lower hysteresis threshold values with historical results .

[0063] The present application also provides a computer device comprising a memory and a processor, wherein the memory stores a computer program, and the processor implements the steps of the above method when executing the computer program.

[0064] The present application also provides a computer readable storage medium having a computer program stored thereon, wherein the computer program is executed by a processor to implement the steps of the above method.

[0065] The present application has the following beneficial effects:

[0066] The present application solves the problem of insufficient detection precision of existing magnetometer interference detection. First, three-axis magnetometer data acquisition and static baseline estimation are performed, then three-axis magnetometer instantaneous feature extraction is performed, then magnetometer interference detection score calculation is performed, and finally robust adaptive threshold generation and magnetometer interference judgment are performed. The present application has high reliability, strong universality, high precision, and good practicability. The present application can accurately identify magnetic interference by statistical features and dynamic threshold discrimination without additional hardware, thereby balancing detection precision and system real-time performance. BRIEF DESCRIPTION OF DRAWINGS

[0067] Figure 1 Fig. 1 is a flowchart of a magnetometer interference detection method based on multi-feature fusion according to the present application.

[0068] Figure 2 Fig. 2 is a structural schematic diagram of a magnetometer interference detection device based on multi-feature fusion according to the present application.

[0069] Figure 3The internal structure schematic diagram of the computer equipment of an embodiment of the present application.

[0070] The implementation, functional features and advantages of the present application will be further described with reference to the embodiments and the accompanying drawings. DETAILED DESCRIPTION

[0071] It should be understood that the specific embodiments described herein merely serve to explain the present application and do not limit the present application.

[0072] As Figure 1 shown, the present application also provides a magnetometer interference detection method based on multi-feature fusion, comprising:

[0073] S1, multi-channel three-axis magnetometer data acquisition and static baseline estimation. The three-axis magnetometer data is collected and constructed into a sampling matrix M, and the static baseline statistics and static threshold are calculated by using the initial static samples of the sampling matrix M; wherein the static baseline statistics include the median and the median absolute deviation ; specifically comprising:

[0074] S101, real-time acquisition of three-axis magnetometer data, and construction of the three-axis magnetometer data into a sampling matrix in time sequence, denoted as , is a three-axis measurement matrix, is the original measurement value at time i in channel c, and N is the number of sampling points;

[0075] S102, calculation of static baseline statistics by using the initial static samples of the sampling matrix M, including the median and the median absolute deviation , and the calculation formula is:

[0076] ;

[0077] wherein, is the median of the initial static interval data; is the sample for static baseline estimation, and j is the time index in the window; is the initial window size; is the median absolute deviation of the initial static interval data;

[0078] S103, calculation of the static threshold according to the median and the median absolute deviation , and the calculation formula is:

[0079] ;

[0080] wherein, This is a static threshold. This is the median of the initial static interval data. The absolute deviation of the median of the initial static interval data. This is the static threshold coefficient. Simultaneously, algorithm parameters are initialized, including the local sliding window half-width, the base Z-score threshold, and smoothing parameters.

[0081] S2. Instantaneous feature extraction of a triaxial magnetometer.

[0082] The magnetic field magnitude is calculated sequentially based on the sampling matrix M. Magnetic tilt angle Magnetic flight angle and the corresponding instantaneous increment of the magnetometer characteristic Specifically, it includes:

[0083] S201. Calculate the magnetic field magnitude sequentially according to the sampling matrix M. The calculation formula is:

[0084] ;

[0085] in, For at any time k The modulus of the triaxial magnetometer, when i=k

[0086] ;

[0087] For at any time k The magnetometer's X-axis output, For at any time k The magnetometer's Y-axis output, For at any time k The Z-axis output of the magnetometer;

[0088] S202, according to time k The magnetic field magnitude Calculation time k magnetic inclination Magnetic flight angle The formulas for calculating the magnetic tilt angle and the magnetic flight angle are as follows:

[0089] ;

[0090] S203, Based on the magnetic field modulus Magnetic tilt angle Magnetic flight angle The corresponding instantaneous increment of the magnetometer characteristic is calculated using the following formula:

[0091] ;

[0092] wherein, is the three-axis magnetometer module value increment at time k , is the magnetic inclination angle increment at time k , is the magnetic heading angle increment at time k , is the magnetic field module value at time k-1 , is the magnetic inclination angle at time k-1 , is the magnetic heading angle at time k-1 .

[0093] S3, magnetometer interference detection score calculation.

[0094] According to the magnetometer feature instantaneous increment, the abnormality is judged by a candidate indication function, and the abnormality proportion is counted by a causal confirmation window, and the abnormality proportion is weighted and fused into a comprehensive detection score ; Specifically comprising:

[0095] S301, according to the magnetometer feature instantaneous increment, the abnormality is judged by a candidate indication function; wherein the expression of the candidate indication function is:

[0096] ;

[0097] wherein, is the indication function; j is the historical time in the causal confirmation window, indicating the sub-time in the causal confirmation window at time k; is the three-axis magnetometer module value increment at time j , is the magnetic inclination angle increment at time j , is the magnetic heading angle increment at time j , is the multiple factor of the module value difference discrimination, is the standard deviation of the module value in the causal local window, is the discrimination threshold of the magnetic inclination angle, is the discrimination threshold of the magnetic heading angle;

[0098] S302, based on the causal confirmation window, the abnormality occurrence proportion of the three kinds of magnetometer feature instantaneous increments is counted, and the formula is:

[0099] ;

[0100] wherein, , , is the abnormality occurrence proportion of the three kinds of magnetometer feature instantaneous increments at time k ,L is the length of the causal confirmation window, is an indicator function;

[0101] S303, linearly fuse the abnormal occurrence proportions of the three kinds of magnetometer feature instantaneous increments according to preset weights to obtain a comprehensive detection score of magnetometer anomaly detection , and the formula is:

[0102] ;

[0103] wherein, is the comprehensive detection score at time k , , , is a coefficient for linearly weighting the three proportions to obtain the comprehensive score, , , is the proportion of abnormal occurrence of the instantaneous increment of the three kinds of magnetometer features at time k .

[0104] S4, robust adaptive threshold generation and magnetometer interference determination.

[0105] According to the comprehensive detection score , a central estimate of the comprehensive detection score is calculated , and a bias is calculated, and a main threshold and upper and lower hysteresis thresholds are constructed in combination with a static threshold ; specifically comprising:

[0106] S401, according to the comprehensive detection score , a central estimate of the comprehensive detection score is calculated , and the calculation formula is:

[0107] ;

[0108] wherein, is a central trend estimate obtained by exponentially weighted averaging the comprehensive detection score at time k , is a central trend estimate obtained by exponentially weighted averaging the comprehensive detection score at time k-1 , is a weight coefficient, is the comprehensive detection score at time k .

[0109] S402, according to the central estimate of the comprehensive detection score , a bias is calculated, and the calculation formula is:

[0110] ;

[0111] wherein, is the exponentially weighted average of absolute deviation at time k, is the exponentially weighted average of absolute deviation at time k-1, is the weight coefficient, is the comprehensive detection score at time k is the center tendency estimation by exponentially weighted average of the comprehensive detection score at time k is the approximation value for comparison with the standard deviation under normal distribution.

[0112] S403, according to the center estimation of the comprehensive detection score and the deviation , a main threshold is constructed in combination with a static threshold and a static lower limit protection, and the formula is:

[0113] ;

[0114] wherein, is the main threshold at time k , and is the multiple factor of the main threshold mesoscale term.

[0115] S404, according to the main threshold an up-and-down hysteresis threshold is derived, and the formula is:

[0116] ;

[0117] wherein, is the entry threshold at time k , and is the exit threshold at time k , and is the upper boundary coefficient, is the lower boundary coefficient, is the approximation value for comparison with the standard deviation under normal distribution.

[0118] S5, by comparing the comprehensive detection score with the main threshold and the up-and-down hysteresis threshold and the historical results, a magnetic interference judgment result is output ; wherein, the expression of the magnetometer anomaly judgment is: ;

[0119] wherein, is the magnetic interference judgment result at time k , and is the magnetic interference judgment result at time k-1 .​​ a comprehensive detection score at time k t, an entry threshold at time k t, an exit threshold at time k t.

[0120] Subsequently, the running phase, only the calculation is performed according to steps S1-S5 to perform the magnetometer anomaly detection in real time.

[0121] The application is explained below with specific examples.

[0122] In order to verify the correctness of the application, actual measurement experiments are carried out, and the performance parameters of the magnetometer in the experiment are listed in Table 1.

[0123] Table 1 Performance parameters of the magnetometer in the actual measurement experiment

[0124]

[0125] Three groups of actual measurement experiments are carried out by collecting magnetometer interference detection data from the magnetometer corresponding to Table 1, and the magnetometer interference detection is carried out by using the method provided by the application and the single magnetometer module value-based magnetic interference detection algorithm. The comparison of the success rate of magnetometer interference detection is shown in Table 2.

[0126] Table 2 Comparison table of magnetometer interference detection success rate

[0127]

[0128] As can be seen from Table 2, the success rate of magnetometer interference detection by the method of the application reaches 93.6%, and is improved by about 20% compared with the single module value detection method, indicating that the method provided by the application can effectively perform magnetometer interference detection. Therefore, the effectiveness and correctness of the method provided by the application are verified.

[0129] As shown in Figure 2 , the application further provides a magnetometer interference detection device based on multi-feature fusion, based on the magnetometer interference detection method based on multi-feature fusion as described above, the device comprises:

[0130] The acquisition module 1 is used for acquiring three-axis magnetometer data and constructing a sampling matrix M, and calculating a static baseline statistic and a static threshold using the initial static sample of the sampling matrix M; wherein the static baseline statistic includes a median and a median absolute deviation ;

[0131] The calculation module 2 is used for sequentially calculating a magnetic field module value , a magnetic inclination angle , magnetic heading and its corresponding magnetometer feature instantaneous increment

[0132] fusion module 3, for judging the anomaly by candidate indicator function according to the magnetometer feature instantaneous increment, and for calculating the anomaly proportion by cause-effect confirmation window, and for weighting and fusing the anomaly proportion into comprehensive detection score

[0133] construction module 4, for constructing the main threshold value and the upper and lower hysteresis threshold value according to the comprehensive detection score central estimate of the comprehensive detection score and deviation , and combining with static threshold

[0134] output module 5, for outputting the magnetic interference judgment result by comparing the comprehensive detection score with the main threshold value and the upper and lower hysteresis threshold value and historical results .

[0135] The above-mentioned modules are used to correspondingly execute the steps in the above-mentioned magnetic interference detection method based on multi-feature fusion, and the specific implementation manners are described in the above-mentioned method embodiments, which will not be described here.

[0136] As shown in Figure 3 , the present application further provides a computer device, which can be a server, and the internal structure thereof can be as shown in Figure 3 . The computer device comprises a processor, a memory, a network interface and a database connected through a system bus. The processor of the computer device is used to provide calculation and control capabilities. The memory of the computer device comprises a non-volatile storage medium and an internal memory. The non-volatile storage medium stores an operating system, a computer program and a database. The internal memory provides an environment for the operation of the operating system and the computer program in the non-volatile storage medium. The database of the computer device is used to store all data required by the process of the magnetic interference detection method based on multi-feature fusion. The network interface of the computer device is used to communicate with external terminals through network connection. The computer program is executed by the processor to implement the magnetic interference detection method based on multi-feature fusion.

[0137] Those skilled in the art can understand that Figure 3 the structure shown in is only a block diagram of part of the structure related to the scheme of the present application, and does not constitute a limitation on the computer device to which the scheme of the present application is applied.

[0138] ​​​An embodiment of the present application further provides a computer readable storage medium, which stores a computer program. The computer program is executed by a processor to implement any one of the above-mentioned multi-feature fusion based magnetometer interference detection methods.

[0139] Those skilled in the art can understand that all or part of the processes in the above-mentioned embodiment methods can be completed by computer program instruction related hardware. The computer program can be stored in a non-volatile computer readable storage medium, and the program can include the processes of the above-mentioned method embodiments when executed. Any reference to memory, storage, database or other medium provided by the present application and used in the embodiments can include non-volatile and / or volatile memory. The non-volatile memory can include read-only memory (ROM), programmable ROM (PROM), electrically programmable ROM (EPROM), electrically erasable programmable ROM (EEPROM) or flash memory (Flash Memory). The volatile memory can include random access memory (RAM), such as dynamic RAM used as main storage or static RAM commonly used in cache memory. As an illustration but not limitation, RAM has various forms, such as static RAM (SRAM), dynamic RAM (DRAM), synchronous DRAM (SDRAM), double data rate SDRAM (DDR SDRAM), and Rambus DRAM (RDRAM), etc.

[0140] It should be noted that in this document, the terms "comprising", "including", or any other variant thereof are intended to cover non-exclusive inclusions, so that processes, devices, articles or methods including a series of elements not only include those elements, but also include other elements not explicitly listed, or include elements inherent to such processes, devices, articles or methods. Without more limitations, the element defined by the statement "comprising a" does not exclude the presence of other identical elements in the process, device, article or method including the element.

[0141] The above description is only preferred embodiments of the present application, and does not limit the patent scope of the present application. Any equivalent structure or equivalent process transformation, or direct or indirect application in other related technical fields, as described in the specification and drawings of the present application, are also included in the patent protection scope of the present application.

Claims

1. A method for detecting interference of a magnetometer based on multi-feature fusion, characterized in that, Comprise: S1, collect triaxial magnetometer data and construct as a sampling matrix M, calculate static baseline statistics and static threshold value using initial static samples of the sampling matrix M ; wherein the static baseline statistics include median and median absolute deviation ; S2, sequentially calculate the magnetic field module value according to the sampling matrix M , magnetic inclination , magnetic heading and its corresponding magnetometer characteristic instantaneous increment ; S3, judging the abnormality by a candidate indicator function according to the magnetic meter feature instantaneous increment, and statistically calculating the abnormality proportion by a cause-effect confirmation window, and weighting and fusing the abnormality proportion into a comprehensive detection score ; specifically comprising: S301, according to the magnetometer characteristic instantaneous increment through candidate indicator function judge abnormal;Wherein, the expression of candidate indicator function is: ； wherein, is an indicator function; j is a historical time within the causal confirmation window, representing a sub-time within the causal confirmation window at time k; is the increment of the three-axis magnetometer module value at time j is the increment of the magnetic inclination angle at time j is the increment of the magnetic heading angle at time j is a multiple factor of the module value difference discrimination, is the standard deviation of the module value on the causal local window, is the discrimination threshold of the magnetic inclination angle, is the discrimination threshold of the magnetic heading angle;​​​ S302, based on the causal confirmation window statistics three kinds of magnetometer characteristic instantaneous increment abnormal appearance ratio, its formula is: ； wherein , , is the time k the proportion of abnormal occurrences of instantaneous increments of the three magnetometer features, L is the length of the causal confirmation window, is an indicator function; S303, linearly fuse the abnormal occurrence proportions of the three magnetometer feature transient increments according to preset weights to obtain a comprehensive detection score of magnetometer anomaly detection The formula is: ； wherein, is the integrated detection score at time k , , , is the coefficient for linearly weighting the three ratios to obtain the integrated score, , , is the integrated detection score at time k the proportion of the three magnetometer characteristic instantaneous incremental anomalies S4, the center estimate of the comprehensive detection score calculating the center estimate of the comprehensive detection score with a bias and combining with a static threshold constructing a main threshold and upper and lower hysteresis thresholds; wherein constructing the main threshold and the upper and lower hysteresis thresholds comprises: constructing the main threshold and the upper and lower hysteresis thresholds according to the center estimate of the comprehensive detection score with a bias combining with a static threshold constructing the main threshold and applying a static lower limit protection, and the formula is: ； wherein is a main threshold value at time k is a main threshold mesoscale term, is a multiple factor of the main threshold mesoscale term; According to the main threshold value The derived upper and lower hysteresis thresholds are given by the formulae: ； wherein is an entry threshold for time of day k is an exit threshold for time of day is an entry threshold for day of week k is an exit threshold for day of week is an upper bound coefficient is a lower bound coefficient S5, outputting a magnetic interference determination result by synthesizing the detection scores comparing with the main threshold value and the upper and lower hysteresis threshold values and the history result, outputting a magnetic interference determination result . 2.The multi-feature fusion based magnetometer interference detection method according to claim 1, characterized in that, The S1 specifically includes: S101, real-time acquisition three-axis magnetometer data, and three-axis magnetometer data is constructed as sampling matrix according to time sequence, and is recorded as ， is a three-axis measurement matrix, is the raw measurement value at time i in channel c, and N is the number of sampling points. S102、Calculate static baseline statistics, including median, using initial static samples of the sampling matrix M and the median absolute deviation The formula is: ； wherein, is the median of the initial static interval data; is the sample for static baseline estimation, j is the time index within the window; is the initial window size; is the absolute deviation of the median of the initial static interval data; S103、According to the median and the absolute deviation of the median Calculate the static threshold The formula is: ； wherein, is a static threshold, is a median of the initial static interval data, is an absolute deviation of the median of the initial static interval data, is a static threshold coefficient.

3. The multi-feature fusion based magnetometer interference detection method of claim 2, wherein, The S2 specifically includes: S201、According to the sampling matrix M, the magnetic field modulus value is calculated in sequence The calculation formula is: ； wherein, is the three-axis magnetometer module value at time k ; ; is the magnetometer X-axis output at time k ; is the magnetometer Y-axis output at time k ; is the magnetometer Z-axis output at time k ; S202, according to time k The magnetic field magnitude Calculation time k magnetic inclination Magnetic flight angle The calculation formula is: ； S203、According to the time k of the magnetic field modulus , magnetic inclination , magnetic course , calculate its corresponding magnetometer characteristic instantaneous increment, the calculation formula is: ； in, For at any time k The magnitude increment of the triaxial magnetometer For at any time k The magnetic tilt increment, For at any time k The magnetic flight angle increment, For a moment k-1 The magnetic field magnitude, For a moment k-1 Magnetic tilt angle, For a moment k-1 Magnetic flight angle.

4. The multi-feature fusion based magnetometer interference detection method of claim 3, wherein, The S4 specifically includes: S401、According to the comprehensive detection score Central estimate of the comprehensive detection score , the calculation formula is: ； wherein, is a center tendency estimate obtained by exponentially weighted averaging of the composite detection scores at time k is a center tendency estimate obtained by exponentially weighted averaging of the composite detection scores at time is a center tendency estimate obtained by exponentially weighted averaging of the composite detection scores at time k-1 is a center tendency estimate obtained by exponentially weighted averaging of the composite detection scores at time is a weight coefficient, is a composite detection score at time k is a composite detection score at time S402、estimating a center according to the comprehensive detection score Computing bias , the formula is: ； in, The exponentially weighted average of the absolute deviations at time k. The exponentially weighted average of the absolute deviations at time k-1. These are the weighting coefficients. For a moment k The comprehensive test score, For time k The central tendency estimate is obtained by using an exponentially weighted average of the comprehensive test scores. This is an approximation used for comparison with the standard deviation under a normal distribution.

5. The multi-feature fusion based magnetometer interference detection method of claim 4, wherein, In the S5, the expression of magnetometer abnormality determination is: ； wherein is the magnetic interference determination result at time k is the magnetic interference determination result at time is the magnetic interference determination result at time k-1 is the magnetic interference determination result at time is the comprehensive detection score at time k is the comprehensive detection score at time 6. A multi-feature fusion based magnetometer interference detection device, characterized in that, The multi-feature fusion-based magnetometer interference detection method according to any one of claims 1-5, the device comprises: a collection module for collecting tri-axial magnetometer data and constructing a sample matrix M, calculating a static baseline statistic and a static threshold using initial static samples of the sample matrix M ; wherein the static baseline statistic includes a median and a median absolute deviation ; a calculating module for sequentially calculating the magnetic field modulus , the magnetic inclination , the magnetic heading and their corresponding magnetometer characteristic instantaneous increments ; a fusion module for judging abnormality by candidate indicator functions according to the magnetometer feature instantaneous increment, and for statistically calculating abnormality proportion by a causal confirmation window, and for weighting and fusing the abnormality proportion into a comprehensive detection score ; specifically comprising: S301, according to the magnetometer characteristic instantaneous increment through candidate indicator function judge abnormal;Wherein, the expression of candidate indicator function is: ； wherein, is an indicator function; j is a historical time within the causal confirmation window, representing a sub-time within the causal confirmation window at time k; is an indicator function; j is a historical time within the causal confirmation window, representing a sub-time within the causal confirmation window at time k; j is a three-axis magnetometer module value increment at time is a three-axis magnetometer module value increment at time j is a magnetic inclination angle increment at time is a magnetic inclination angle increment at time j is a magnetic inclination angle increment at time is a module value difference discrimination multiplier, is a module value standard deviation over the causal local window, is a magnetic inclination angle discrimination threshold value, is a magnetic inclination angle discrimination threshold value; S302, based on the causal confirmation window statistics three kinds of magnetometer characteristic instantaneous increment abnormal appearance ratio, its formula is: ； wherein , , is the time k the proportion of abnormal occurrences of instantaneous increments of the three magnetometer features, L is the length of the causal confirmation window, is the indicator function; The abnormal occurrence proportion of the instantaneous increment of the three magnetometer features is linearly fused according to preset weights to obtain a comprehensive detection score of magnetometer anomaly detection The formula is: ； wherein, is the integrated detection score at time k , , , is the coefficient for linearly weighting the three ratios to obtain the integrated score, , , is the integrated detection score at time k the proportion of the three magnetometer characteristic instantaneous incremental anomalies a construction module for constructing a main threshold value and upper and lower hysteresis threshold values in accordance with the center estimate of the comprehensive detection score a center estimate of the comprehensive detection score with a bias and in combination with a static threshold value a construction module for constructing a main threshold value and upper and lower hysteresis threshold values in accordance with the center estimate of the comprehensive detection score with a bias in combination with a static threshold value a construction module for constructing a main threshold value and upper and lower hysteresis threshold values in accordance with the center estimate of the comprehensive detection score ； wherein is a main threshold value at time k is a main threshold value at time is a multiple factor for the main threshold mesoscale term; According to the main threshold value The derived upper and lower hysteresis thresholds are given by the formulae: ； wherein is an entry threshold for time k , is an exit threshold for time k , is an upper bound coefficient, is a lower bound coefficient; an output module for outputting a magnetic interference determination result based on the comparison of the main threshold value and the upper and lower hysteresis threshold values and the history result an output module for outputting a magnetic interference determination result based on the comparison of the main threshold value and the upper and lower hysteresis threshold values and the history result . 7.A computer device, comprising a memory and a processor, wherein the memory stores a computer program, and the computer device is configured to perform the method according to any one of claims 1-6 when the computer program is executed by the processor. The processor executes the computer program and realizes the steps of the method in any one of claims 1 to 5.

8. A computer-readable storage medium having stored thereon a computer program, characterized in that, The computer program is executed by the processor and realizes the steps of the method in any one of claims 1 to 5.

Images