Multi-position north seeking method of north seeker

Abstract

The invention discloses a multi-position north seeking method of a north seeker, which comprises the following steps of: calibrating an Allan variance curve of a gyroscope, and determining a characteristic time constant corresponding to each noise type; controlling the north seeker to sequentially rotate to a plurality of preset positions, identifying a dominant noise type in real time at each position through power spectral density analysis, and dynamically adjusting a sampling or signal processing strategy; meanwhile, acquiring internal multi-point temperature information to calculate a temperature gradient, and compensating zero offset of the gyroscope in real time by using a compensation model containing the temperature gradient; and finally, calculating the northbound angle based on the sampling data of all positions and the compensation result. According to the method, through adaptive noise suppression and multivariable temperature compensation, north-seeking precision and temperature adaptability in a complex environment are remarkably improved, and measurement efficiency and system stability are both considered by means of an optimization algorithm.

Description

Technical Field

[0001] This invention belongs to the field of inertial navigation technology, specifically a multi-position north-finding method for a north-finding instrument. Background Technology

[0002] North finders, as inertial measurement devices that determine the northward reference based on the Earth's rotational angular velocity vector, are widely used in surveying, military, aerospace, and other fields. Their core performance indicators include north-finding accuracy, measurement time, and environmental adaptability, which are primarily affected by the noise characteristics of the gyroscope sensor, temperature drift, and signal processing algorithms.

[0003] In existing technologies, the measurement modes of north-finding instruments are mainly divided into two categories: multi-position cyclic measurement and continuous rotation measurement. The multi-position cyclic measurement mode collects data by rotating the north-finding instrument to multiple fixed positions and uses position difference to suppress fixed zero bias. However, the traditional method uses fixed sampling time and uniform signal processing strategy, which cannot adapt to the changes in noise type at different positions (such as white noise, 1 / f noise, and vibration noise dominating alternately), resulting in insufficient noise suppression and limited north-finding accuracy. The continuous rotation measurement mode achieves rapid north-finding by rotating at a constant angular velocity. However, the angular velocity value lacks the correlation constraint with the gyroscope noise correlation time and temperature drift characteristic time, making it susceptible to the superposition of noise and temperature drift. In addition, the zero bias compensation is not accurate during the phase-locked amplification calculation, which also affects the north-finding accuracy.

[0004] Furthermore, gyroscope zero-bias drift is temperature-sensitive, and existing temperature compensation methods mostly consider only single-point temperature or static temperature, ignoring the influence of temperature gradients and rates of temperature change between the chip and the casing on zero bias, resulting in poor temperature adaptability. Meanwhile, traditional power spectral density analysis relies on software calculations, which lacks real-time performance and fails to meet the needs of dynamic noise identification; the data fusion process does not fully utilize the differences in noise characteristics at each measurement location, further limiting the improvement of north-finding accuracy.

[0005] Therefore, there is an urgent need for a north-finding method that can adapt to noisy environments, suppress temperature drift, and balance measurement accuracy and efficiency to address the shortcomings of existing technologies. Summary of the Invention

[0006] The purpose of this invention is to provide a multi-location north-finding method for a north-finding instrument, so as to solve the problems of insufficient noise suppression, large temperature drift, single measurement mode, and low calculation efficiency of existing north-finding instruments mentioned in the background art.

[0007] To solve the above-mentioned technical problems, the technical solution adopted by the present invention is as follows:

[0008] A multi-location north-finding method for a north-finding instrument includes the following steps:

[0009] Step S1: Calibration and initialization. The Allan variance curve of the gyroscope is calibrated, and the white noise coefficient Q, angular random walk coefficient N, and velocity random walk coefficient K are extracted. The noise correlation time τ is obtained by performing autocorrelation analysis on the gyroscope's static output data. c (By fitting the autocorrelation function R(τ)=R(0)·exp(-τ / τ) c The time constant for the autocorrelation function to decay to 1 / e is obtained by feature extraction from long-term sampling data of the gyroscope under varying temperature conditions, which is defined as the time constant for the autocorrelation function to decay to 1 / e. T ;

[0010] Step S2, multi-position cyclic measurement: control the north finder to rotate sequentially around the vertical axis to multiple preset measurement positions; at each measurement position, perform the following operations:

[0011] Step S201: Perform power spectral density analysis on the real-time acquired gyroscope output signal;

[0012] Step S202: Based on the power spectral density analysis results, identify the dominant noise type of the current gyroscope output signal in real time;

[0013] Step S203: Based on the identified dominant noise type, dynamically adjust the sampling time or signal processing strategy at the current measurement location;

[0014] Step S3, temperature compensation: obtain temperature information from at least two different locations inside the north finder, calculate the temperature gradient based on the temperature information, and use a compensation model that includes the temperature gradient as an independent variable to compensate for the zero-bias drift of the gyroscope.

[0015] Step S4, data calculation and fusion: Based on the sampling data of all measurement locations in step S2 and the compensation results in step S3, calculate and output the north angle.

[0016] According to the above technical solution, in step S2, the dominant noise type is identified in real time based on the power spectral density analysis results, specifically including:

[0017] Extract at least two characteristic parameters from the power spectral density (PSD). The characteristic parameters are selected from: PSD flatness F and the slope k of the logarithmic fit of the PSD in the 1-10Hz frequency band. PSD And the peak ratio R of the line spectrum;

[0018] Based on the comparison between the extracted feature parameter values ​​and the preset threshold, the dominant noise type is determined to be at least two of the following: white noise dominant type, 1 / f noise dominant type, or vibration noise dominant type.

[0019] According to the above technical solution, in step S2, the sampling time or signal processing strategy at the current measurement position is dynamically adjusted, specifically as follows:

[0020] When it is determined to be the white noise dominant type, the first sampling time t is adopted s1 for sampling, and its value is associated with the white noise coefficient of the gyroscope. The optimal relationship is:

[0021]

[0022] where k1 ∈ [0.5, 1.0];

[0023] When it is determined to be the 1 / f noise dominant type, wavelet threshold denoising processing is started. The db4 wavelet is used for 5-layer decomposition, and the soft threshold processing is performed on the detail coefficients of the 3rd to 5th layers using the Donoho-Johnstone universal threshold; the calculation formula for the universal threshold λ is:

[0024]

[0025] where σ is the noise standard deviation, and N data is the data length;

[0026] When it is determined to be the vibration noise dominant type, digital notch filtering for the vibration line spectrum frequency is started. The filter is a second-order IIR structure, and its transfer function is:

[0027]

[0028] where ω0 is the vibration line spectrum angular frequency, r is the damping coefficient and satisfies 0 < r < 1; the center frequency of the filter is configured as the vibration line spectrum frequency, and its bandwidth is 5% to 10% of the center frequency. The quality factor Q f (dimensionless, defined as ), and Δω is the 3dB bandwidth) takes values from 5 to 20.

[0029] According to the above technical solution, step S3 specifically includes:

[0030] Measure the surface temperature T of the gyroscope chip through the first temperature sensing unit chip ;

[0031] Measure the inner wall temperature T of the north-seeking instrument housing through the second temperature sensing unit shell ;

[0032] Calculate the temperature gradient ΔT = T chip - T shell ;

[0033] The expression of the compensation model includes the temperature gradient ΔT and the housing temperature T shell as independent variables.

[0034] According to the above technical solution, the expression of the compensation model is:

[0035]

[0036] Where, ΔB comp The zero-bias drift to be compensated. K is the rate of change of the outer shell temperature, and k1, k2, and k3 are compensation coefficients determined through calibration experiments.

[0037] According to the above technical solution, the multi-position cyclic measurement in step S2 can be switched to a continuous rotation measurement mode, specifically as follows:

[0038] Control the north finder to maintain a constant angular velocity ω r Continuous rotation around the vertical axis;

[0039] angular velocity ω r The value of satisfies the following constraints:

[0040] ω min <ω r <ω max

[0041] in, ω max =min(0.2f Nyquist ,ω limit );f c 1 / f is the noise corner frequency; f Nyquist ω is the Nyquist frequency. limit The upper limit of safe angular velocity is determined based on the gyroscope bandwidth or the system's mechanical resonant frequency;

[0042] The dynamically adjusted signal processing strategy is phase-locked amplification, which involves generating a reference signal based on the rotational position signal, extracting northward information from the gyroscope signal through coherent integration, and selecting the integration time T of the phase-locked amplification to cover 2-5 rotation cycles. n∈{2,3,4,5}, where n is the number of rotation cycles covered by the integration time, ensuring that northward information is fully extracted.

[0043] According to the above technical solution, the phase-locked amplification calculation specifically includes:

[0044] A sinusoidal reference signal sin(ω) of the same frequency is generated based on the rotational position signal. r t) and cosine reference signal cos(ω) r t);

[0045] The gyroscope output signal Ω(t) is multiplied by the two reference signals and then integrated to obtain the orthogonal integral components X and Y:

[0046]

[0047] Where T is the integration time;

[0048] The north angle α0 is calculated using the following formula:

[0049] α0 = arctan 2(Y-Y0, X-X0)

[0050] Where X0 and Y0 are the zero-bias compensation terms obtained from calibration.

[0051] According to the above technical solution, in step S2, the power spectral density analysis is accelerated by hardware. The field programmable gate array (FPGA) is used to perform fast Fourier transform (FFT) on the acquired gyroscope signal. The signal sampling adopts a sliding window method with overlap between windows.

[0052] According to the above technical solution, in step S4, data calculation and fusion includes: constructing opposing position differences based on azimuth data from multiple measurement locations; and calculating the north angle result {α} from M cyclic measurements. 01 α 02 ,...α 0M The variance of the i-th measurement is calculated by weighting and fusing the data based on the dispersion of each measurement. Calculated based on n independent sampling points within this loop:

[0053]

[0054] Where, α ij The north angle calculated for the j-th sampling point in the i-th iteration. Assign weights to the mean north-facing angle for this cycle:

[0055]

[0056] Final North Angle:

[0057]

[0058] According to the above technical solution, the elevation angle θ of the north-finding instrument is obtained in real time during the north-finding calculation process. p With roll angle θ r The local horizontal projection component of the Earth's rotational angular velocity is corrected based on the pitch and roll angles.

[0059] Compared with the prior art, the present invention has the following beneficial effects:

[0060] This invention introduces a noise suppression strategy based on power spectral density feature identification and targeted processing. It employs adaptive dynamic sampling, wavelet threshold denoising, and precise notch filtering to address different noise characteristics, significantly improving orientation accuracy in complex environments. By constructing a compensation model based on multiple temperature parameters, it effectively suppresses zero-bias drift over a wide temperature range. The invention designs a working mode compatible with multiple positions and continuous rotation, and sets reasonable angular velocity constraints, achieving a flexible adaptation and balance between high precision and high efficiency in various scenarios. Finally, by integrating algorithms such as opposing differential, weighted fusion, and attitude correction, it effectively improves the long-term stability and repeatability accuracy of the system. Attached Figure Description

[0061] Figure 1 This is a flowchart of the multi-location north-finding method of the north-finding instrument of the present invention. Detailed Implementation

[0062] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0063] Example 1

[0064] like Figure 1 As shown, a multi-location north-finding method for a north-finding instrument in this invention includes the following steps:

[0065] Step S1, Calibration and Initialization: Calibrate the Allan variance curve of the gyroscope, and extract the white noise coefficient Q, angular random walk coefficient N, and velocity random walk coefficient K; Allan variance calibration must be performed in a constant temperature environment (25℃±1℃), collecting at least 1 hour of static angular rate data, and the fitting formula is:

[0066]

[0067] Where σ(τ) is the Allan standard deviation (unit: ° / h), τ is the average time (unit: h), and Q, N, and K are obtained by least squares fitting.

[0068] By performing autocorrelation analysis on the static output data of the gyroscope, the noise correlation time τ is obtained. c (By fitting the autocorrelation function R(τ)=R(0)·exp(-τ / τ) c (This is obtained by defining the time constant for the autocorrelation function to decay to 1 / e).

[0069] By extracting features from long-term sampling data of gyroscope zero bias under variable temperature conditions (-40℃~60℃, heating / cooling rate 5℃ / min), the temperature drift characteristic time τ was obtained. T (Unit: s, defined as the time it takes for the zero bias drift to reach 63% of the steady-state value), complete the initialization preparation.

[0070] Step S2, Multi-position cyclic measurement: Control the north finder to rotate sequentially around the vertical axis to multiple preset measurement positions (preferably 4-8, with the included angle between adjacent positions). (Total number of locations); at each measurement location, perform the following operations:

[0071] Step S201: Power spectral density analysis: The acquired gyroscope signal is subjected to a Fast Fourier Transform (FFT) using a Field Programmable Gate Array (FPGA), with N FFT points. FFT =1024, sampling frequency f s Based on the gyroscope bandwidth setting (e.g., f) s =100Hz), the analysis frequency band covers 0.1Hzf s / 2; Signal sampling uses a sliding window method, with a window length of N. win =512-1024 sampling points, window overlap rate ( 75%, such as η = 60%), balancing spectral estimation accuracy (the higher the overlap rate, the smaller the spectral estimation variance) and data update rate.

[0072] Step S202, Dominant Noise Type Identification: Extract at least two characteristic parameters of the Power Spectral Density (PSD). The characteristic parameters are selected from: PSD flatness F (defined as the ratio of the maximum to minimum PSD value within the selected analysis frequency band) and the logarithmic fitting slope k of the PSD in the 1-10Hz frequency band. PSD 1. Peak Spectrum Ratio (PSD) (the ratio of the peak value of the PSD to the background PSD value). Furthermore, to enable those skilled in the art to more clearly understand and implement the noise adaptive identification and processing process in this invention, the following provides a more detailed explanation of the dominant noise type identification in step 2:

[0073] Detailed implementation of the noise adaptive identification and processing unit: The core of this unit is to perceive the noise characteristics of the gyroscope output signal in real time and dynamically optimize the sampling or processing strategy accordingly. Its process includes three stages: quantitative extraction of feature parameters, decision logic judgment, and adaptive strategy execution.

[0074] First, the real-time acquired gyroscope signals are windowed (e.g., using a Hanning window), and their power spectral density (PSD(f)) is calculated using a Fast Fourier Transform (FFT) accelerated by a Field-Programmable Gate Array (FPGA). Then, the following three key feature parameters are extracted:

[0075] PSD flatness F is used to assess the degree of fluctuation in the noise power spectrum within a selected frequency band, and is calculated using the following formula:

[0076]

[0077] Where, f∈[f low ,f high ], f low and f high Based on the effective bandwidth setting of the gyroscope, it is usually set to f. low =0.1Hz, f high = Sampling frequency / 2, ideal white noise F≈1, actual set threshold F th When F≤F th At that time, it is assumed that the spectrum is flat and that F > 1 for colored noise (such as 1 / f noise).

[0078] Low-frequency band logarithmic spectrum slope k PSD This parameter is used to identify Class 1 / f colored noise. A log-linear fit is performed on the PSD(f) within the 1Hz to 10Hz frequency band:

[0079] log 10 (PSD(f))=k PSD ·log 10 (f)+b

[0080] Using the least squares method for linear fitting, the slope of the fitted line is k. PSD k of white noise PSD k ≈0, 1 / f noise PSD ≈-1, k of random walk noise PSD ≈-2 is an effective characteristic for identifying low-frequency colored noise.

[0081] Detection and calculation of the peak value ratio (PSR) of the line spectrum: This parameter is used to detect periodic vibration disturbances. First, significant line spectrum peaks in the PSD(f) are detected. peak and its frequency f peak And estimate the background noise level PSD near that frequency. bg (e.g., f) peak Median PSD in the region excluding ±5Hz). Peak-to-peak ratio is defined as:

[0082]

[0083] By calculating the peak ratio, a threshold R is set. th (e.g., 10), when R > R th When the vibration noise is determined to be dominant, F > F th And k PSD <k PSD,th The noise is dominated by 1 / f noise, and the rest are dominated by white noise.

[0084] To address the problem that traditional methods cannot adapt to dynamic changes in noise type, this invention employs a rule-based multi-feature fusion decision logic. This logic uses the feature vector (F,k) calculated in real time. PSD The input is R, and the output is a dominant noise type indicator. The specific judgment rules and priorities are as follows (thresholds were determined through prior calibration experiments):

[0085] Input: Feature parameter triples (F, k) calculated in real time PSD ,R);

[0086] Preset threshold: The flatness threshold F is set through experimental calibration. th1 (e.g., 1.5), F th2 (e.g., 3.0), slope threshold k PSD,th (e.g., -0.5), the line spectrum peak ratio threshold R th (e.g., 10), and F th1 >F th2 .

[0087] Decision rules (executed in order of priority):

[0088] Step 1, Vibration and noise determination: If R > R th If the noise is immediately identified as dominant, this rule has the highest priority to avoid strong line spectrum noise interfering with the judgment of other features.

[0089] Step 2, 1 / f noise determination: If vibration noise determination is not triggered, and F>F is satisfied simultaneously. th2 and k PSD <k PSD,th It was determined to be dominated by 1 / f noise, and the dual conditions ensured the accuracy of low-frequency colored noise identification.

[0090] Step 3, White Noise Determination: If none of the above conditions are met, and F≤F th1 It was determined to be dominated by white noise, indicating that the noise power distribution was uniform across the entire frequency band.

[0091] Step 4, Mixed or Transitional State: If neither satisfies F th1 <F≤F th2 And k PSD If the value is close to 0, it is determined to be mixed noise, and is treated as 1 / f noise by default or fine time-frequency analysis is started.

[0092] Output: Dominant noise type flag, used to trigger subsequent dynamic adjustment strategies.

[0093] Dynamic adjustment strategy: Based on the identified dominant noise type, dynamically adjust the sampling time or signal processing strategy at the current measurement location.

[0094] When it is determined to be the white noise dominant type, the first sampling time t is adopted. s1 Sampling is carried out, and its value is associated with the white noise coefficient Q. The optimal relationship is:

[0095]

[0096] where k1 ∈ [0.5, 1.0]. By appropriately extending the sampling time, the low-frequency fluctuation characteristics of the partial smoothing 1 / f noise are reduced.

[0097] When it is determined to be the 1 / f noise dominant type, when it is determined to be the 1 / f noise dominant type, wavelet threshold denoising processing is started. The db4 wavelet is used for 5-layer decomposition, and the soft threshold processing is carried out on the detail coefficients of the 3rd - 5th layers by using the Donoho-Johnstone universal threshold; the calculation formula of the universal threshold λ is:

[0098]

[0099] where σ is the noise standard deviation, N data is the data length, and the soft threshold function is:

[0100] T λ T(x) = sign(x)·max(|x| - λ, 0)

[0101] When it is determined to be the vibration noise dominant type, digital notch filtering for the vibration line spectrum frequency is started. The filter is a second-order IIR structure, and its transfer function is:

[0102] [[ID=3P1]]

[0103] where ω0 is the vibration line spectrum angular frequency, ω0 = 2πf0, and f0 is the vibration line spectrum frequency; r is the damping coefficient and satisfies 0 < r < 1, usually r = 0.95; the bandwidth (Q f is the quality factor, Q f = 5 - 20) The center frequency of the filter is configured as the vibration line spectrum frequency, and its bandwidth is 5% to 10% of the center frequency.

[0104] To avoid system instability caused by frequent switching of strategies near the noise feature boundary, the following mechanism is introduced:

[0105] Hysteresis mechanism: The noise type determination result needs to remain consistent within N consecutive processing cycles (for example, N = 3) before triggering the strategy switch.

[0106] / / 这里原文中ID=48的内容似乎有误,推测应该是ID=48和ID=49的内容合并为一句,翻译如下 Minimum sampling time constraint: In any case, the actual sampling time shall not be lower than the minimum response time T min (for example, 0.1 s) to ensure the basic data update rate.

[0107] The multi-position cyclic measurement in step S2 can be switched to continuous rotation measurement mode, which is suitable for rapid north-finding scenarios. The specific operation is as follows:

[0108] Control the north finder to maintain a constant angular velocity ω r Continuous rotation about the vertical axis, angular velocity ω r The value of satisfies the following constraints:

[0109] ω min <ω r <ω max

[0110] in, ω max =min(0.2f Nyquist ,ω limit );f c 1 / f is the noise corner frequency; f Nyquist ω is the Nyquist frequency. limit This is a safety angular velocity upper limit determined based on the gyroscope bandwidth or the system's mechanical resonant frequency; this constraint ensures the rotation period Shorter than the characteristic time τ of temperature drift T (Suppresses temperature drift) and is longer than the noise correlation time τ c (Smoothed noise).

[0111] The gyroscope angular rate signal Ω(t) and rotation position signal are acquired synchronously during the rotation process;

[0112] The dynamic signal processing strategy is lock-in amplification, specifically including:

[0113] A sinusoidal reference signal sin(ω) of the same frequency is generated based on the rotational position signal. r t) and cosine reference signal cos(ω) r t);

[0114] The gyroscope output signal Ω(t) is multiplied by the two reference signals and then integrated. The integration time T satisfies:

[0115]

[0116] n is a positive integer, and n≥2, ensuring that the integration period covers an integer number of rotation periods, thus obtaining orthogonal integral components X and Y:

[0117]

[0118] The north angle α0 is calculated using the following formula:

[0119] α0 = arctan 2(Y-Y0, X-X0)

[0120] X0 and Y0 are zero-bias compensation terms obtained through calibration, used to eliminate the effects of fixed zero bias.

[0121] Step 3, Temperature Compensation: Obtain temperature information from at least two different locations inside the north-finding instrument, calculate the temperature gradient based on the temperature information, and use a compensation model to compensate for the gyroscope's zero-bias drift. Specifically, this includes:

[0122] The surface temperature T of the gyroscope chip is measured by the first temperature sensing unit. chip The first temperature sensing unit is preferably a miniature thermocouple mounted on the edge of the gyroscope chip die, with a measurement accuracy higher than ±0.1℃;

[0123] The temperature T of the inner wall of the north-finding instrument's outer casing is measured by the second temperature sensing unit. shell The second temperature sensing unit is preferably a thermistor fixed to the inner wall of the housing, with a response time of less than 10ms;

[0124] Temperature information is collected at a frequency of 10-100Hz, and the temperature gradient ΔT between the chip and the casing is calculated. grad =T chip -T shell ;

[0125] Calculate the rate of change of the outer casing temperature By continuously collecting T shell The data is obtained by first-order difference calculation, as shown in the following formula:

[0126]

[0127] Where Δt=1 / f T , where is the temperature sampling interval (unit: s).

[0128] The zero-bias drift ΔB to be compensated is calculated using a compensation model. comp The compensation model expression is:

[0129]

[0130] Among them, k1, k2, and k3 are compensation coefficients determined through calibration experiments. The specific calibration method is as follows: the north finder is subjected to a preset temperature cycle (temperature range -40℃ to 60℃, heating / cooling rate 5℃ / min), gyroscope zero bias data at different temperatures are collected, and the compensation coefficients are obtained by fitting using a recursive least squares algorithm.

[0131] Step 4, Data processing and fusion, and attitude correction:

[0132] Data processing and fusion: In the case of a multi-location loop mode, based on the sampling data from all measurement locations and the compensation result (ΔB) from step 3. compFirst, construct the opposing position difference to offset part of the fixed error. Specifically, for the azimuth angle α between the i-th position and the (i+m)-th position (m = M / 2, where M is the total number of positions). i α i+m Calculate the difference Δα i =α i -α i+m This offsets the fixed zero bias error.

[0133] In continuous rotation mode, the initial north angle is obtained directly based on the lock-in amplification calculation results. Then, the north angle results from multiple cyclic measurements are weighted and fused according to the dispersion (e.g., variance) of each measurement. Specifically, the north angle results from M cyclic measurements {α 01 ,α 02 ,...α 0M Based on the dispersion (variance) of each measurement Weight allocation: The variance of the i-th measurement Calculated based on n independent sampling points within this loop:

[0134]

[0135] Where α ij The north angle calculated for the j-th sampling point in the i-th iteration. This is the mean of the northward angles for this cycle. Weights are assigned as follows:

[0136]

[0137] The smaller the dispersion, the larger the weight. The weighted least squares algorithm is used to calculate the final north angle.

[0138]

[0139] If it is a continuous rotation mode, the preliminary north angle is obtained directly based on the phase-locked amplification solution result in step 2 and the compensation result in step 3.

[0140] Attitude correction: During the north-finding calculation, the pitch angle θ is obtained in real time through the accelerometer or tilt sensor built into the north-finding instrument. p With roll angle θ r The local horizontal projection component of the Earth's rotation angular velocity is corrected based on the pitch and roll angles. The correction formula is as follows:

[0141] Ω h =Ω e ·cos L·cosθ p ·cosθ r

[0142] Among them, Ω hΩ represents the corrected projection of the Earth's rotational angular velocity onto the local horizontal plane. e The Earth's rotational angular velocity constant (value 7.292115 × 10⁻⁶) -5 rad / s), L is the latitude of the north-finding instrument; combined with Ω h The north angle calculation result is corrected, and the final north angle is output with an accuracy better than ±0.01° (1σ).

[0143] Example 2

[0144] This embodiment is a further refinement of Embodiment 1. This embodiment provides a specific implementation method for a multi-location north-finding method.

[0145] In multi-location cyclic measurement mode, the following steps are included:

[0146] Calibration and initialization: The MEMS gyroscope was calibrated, and the Allan variance characteristic time t0 = 10s and the noise correlation time τ were measured. c =2s, characteristic time of temperature drift τ t =100s;

[0147] Multiple position settings: 6 preset measurement positions, with an angle of 60° between adjacent positions. The north finder rotates around the vertical axis to each position in sequence and then remains stationary.

[0148] Noise identification and dynamic adjustment:

[0149] At position 1, PSD analysis shows F = 1.2 ≤ F th1 =1.5, k PSD =-0.2∈[k low1 =-0.3,k high1 =0.1], which is determined to be dominated by white noise, and t is used s1 =8s sampling;

[0150] At position 3, PSD analysis shows F = 3.5 > F th2 =2.0, k PSD =-1.1∈[k low2 =-1.5,k high2 =-0.8], determined to be dominated by 1 / f noise, using t s2 =Sampling time is 30 seconds;

[0151] At position 5, PSD analysis shows R = 15 > R th =10, indicating vibration noise as the dominant factor, and activate the digital notch filter with a center frequency of 50Hz, bandwidth of 3Hz, and Q=10;

[0152] Temperature compensation: T measured by miniature thermocouple chip =45.2℃, T measured by the thermistor case=42.5℃, calculate ΔT grad =2.7℃, outer casing temperature change rate The compensation coefficients k1 = 0.002° / ℃, k2 = 0.005° / ℃, and k3 = 0.01°·s / ℃ are used to calculate ΔB. comp =0.002×42.5+0.005×2.7+0.001×0.3=0.1025°, to compensate for the zero bias of the gyroscope.

[0153] Data processing and fusion: Constructing opposing position differences, and using a weighted least squares algorithm to fuse data from six positions, combined with the pitch angle θ. p =0.5°, roll angle θ r =0.3°, latitude L=30°N, attitude correction is performed, and the final output north angle α0=89.763°, north finding accuracy ±0.008°(1σ).

[0154] In continuous rotation measurement mode, the following steps are included:

[0155] Calibration and initialization: Same as multi-position cyclic measurement mode, τ c =2s,τ t =100, take C1=1.5, C2=3.0, calculate the angular velocity constraint range:

[0156]

[0157] That is, 0.021 rad / s < ω r <1.57 rad / s, select ω r =0.5 rad / s = 0.5 (corresponding to a rotation period of 12.57 s);

[0158] Continuous rotation measurement: Control the north finder to rotate continuously around the vertical axis at 0.5 rad / s, and simultaneously acquire the gyroscope angular rate signal Ω(t) and rotation position signal;

[0159] Phase-locked amplification calculation: Generate sin(0.5t) and cos(0.5t) reference signals, take the integration time T = 25.14s (2 rotation cycles), calculate X = 0.0035rad, Y = 0.0012rad, calibrate the zero bias compensation terms X0 = 0.0002rad, Y0 = 0.0001rad, and preliminarily calculate the north angle arctan2(0.0012-0.0001,0.0035-0.0002) = 18.1°.

[0160] Temperature compensation: The temperature compensation process is the same as that of the multi-position cyclic measurement mode. After compensation, the north angle is corrected to 18.092°.

[0161] Attitude correction: based on pitch angle θp =0.2°, roll angle θ r The coordinates were corrected to 0.4° and latitude L=30°N, and the final output north angle was 18.089°, with a north-finding accuracy of ±0.015° (1σ) and a measurement time of 28s.

[0162] Compared with the prior art, the present invention has the following advantages:

[0163] Strong adaptive noise suppression capability: By identifying the dominant noise type through real-time power spectral density analysis and dynamically adjusting the sampling time or signal processing strategy, it can accurately suppress white noise, 1 / f noise and vibration noise, significantly improving the north-finding accuracy in complex noise environments.

[0164] Accurate temperature drift compensation: Considering the temperature gradient between the chip and the casing and the temperature change rate of the casing, a multivariable temperature compensation model is constructed to effectively reduce the impact of temperature changes on the gyroscope's zero bias, and the temperature adaptability range is widened to -40℃~60℃.

[0165] Dual-mode compatibility, balancing accuracy and efficiency: the multi-position cyclic mode is suitable for high-precision scenarios, while the continuous rotation mode is suitable for fast measurement scenarios. The angular velocity value is constrained by noise correlation time and temperature drift characteristics to balance the measurement performance of both modes.

[0166] High calculation efficiency: Power spectral density analysis is achieved through FPGA hardware acceleration, sliding window sampling improves the data update rate, and phase-locked amplification combined with coherent integration shortens the measurement time. The north-finding time in multi-position mode is ≤5min, and the north-finding time in continuous rotation mode is ≤30s.

[0167] Data fusion and attitude correction optimization: Random errors are reduced by using opposing position differential and weighted fusion algorithms, and attitude correction compensates for projection errors caused by tilt, further improving north-finding accuracy and stability.

[0168] It should be noted that, in this document, relational terms such as "first" and "second" are used only to distinguish one entity or operation from another, and do not necessarily require or imply any such actual relationship or order between these entities or operations. Furthermore, the terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such process, method, article, or apparatus.

[0169] Finally, it should be noted that the above descriptions are merely preferred embodiments of the present invention and are not intended to limit the present invention. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art can still modify the technical solutions described in the foregoing embodiments or make equivalent substitutions for some of the technical features. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A multi-location north-finding method for a north-finding instrument, characterized in that: Includes the following steps: Step S1, Calibration and Initialization: Calibrate the Allan variance curve of the gyroscope, and extract the white noise coefficient Q, angular random walk coefficient N, and velocity random walk coefficient K; obtain the noise correlation time τ by performing autocorrelation analysis on the static output data of the gyroscope. c By extracting features from long-term sampling data of gyroscope zero bias under varying temperature conditions, the temperature drift characteristic time τ was obtained. T ; Step S2, Multi-position Cyclic Measurement: Control the north finder to rotate sequentially around the vertical axis to multiple preset measurement positions; at each measurement position, perform the following operations: Step S201: Perform power spectral density analysis on the real-time acquired gyroscope output signal; Step S202: Based on the power spectral density analysis results, identify the dominant noise type of the current gyroscope output signal in real time; Step S203: Based on the identified dominant noise type, dynamically adjust the sampling time or signal processing strategy at the current measurement location; Step S3, Temperature Compensation: Obtain temperature information from at least two different locations inside the north finder, calculate the temperature gradient based on the temperature information, and use a compensation model that includes the temperature gradient as an independent variable to compensate for the zero-bias drift of the gyroscope. Step S4, Data Calculation and Fusion: Based on the sampling data of all measurement locations in Step S2 and the compensation results in Step S3, calculate and output the north angle.

2. The multi-location north-finding method for a north-finding instrument according to claim 1, characterized in that: In step S2, the dominant noise type is identified in real time based on the power spectral density analysis results, specifically including: Extract at least two characteristic parameters from the power spectral density (PSD). The characteristic parameters are selected from: PSD flatness F and the logarithmic fitting slope k of the PSD in the 1-10Hz frequency band. PSD And the peak ratio R of the line spectrum; Based on the comparison between the extracted feature parameter values ​​and the preset threshold, the dominant noise type is determined to be at least two of the following: white noise dominant type, 1 / f noise dominant type, or vibration noise dominant type.

3. The multi-location north-finding method for a north-finding instrument according to claim 2, characterized in that: In step S2, the sampling time or signal processing strategy at the current measurement location is dynamically adjusted, specifically as follows: When the noise type is determined to be white noise dominant, the first sampling time t is used. s1 The sampling is performed, and its value is correlated with the white noise coefficient of the gyroscope. The optimal relationship is as follows: Where k1∈[0.5,1.0]; When the noise type is determined to be 1 / f, wavelet thresholding denoising is initiated. A 5-level decomposition is performed using the db4 wavelet, and the detail coefficients of levels 3-5 are soft-thresholded using the Donoho-Johnstone universal threshold. The formula for calculating the universal threshold λ is as follows: Where σ is the noise standard deviation, N data For data length; When the vibration noise is determined to be the dominant type, a digital notch filter is activated for the vibration line spectrum frequency. The filter is a second-order IIR structure with the following transfer function: where ω0 is the angular frequency of the vibration line spectrum, r is the damping coefficient and satisfies 0 < r < 1; the center frequency of the filter is configured as the vibration line frequency, and its bandwidth is 5% to 10% of the center frequency, and the quality factor Q f (dimensionless, defined as ), and Δω is the 3dB bandwidth) takes values from 5 to 20.

4. The multi-location north-finding method for a north-finding instrument according to claim 1, characterized in that: Step S3 specifically includes: The surface temperature T of the gyroscope chip is measured by the first temperature sensing unit. chip ; The temperature T of the inner wall of the north-finding instrument's outer casing is measured by the second temperature sensing unit. shell ; Calculate the temperature gradient ΔT = T between the chip and the casing. chip -T shell ; The expression for the compensation model includes the temperature gradient ΔT and the shell temperature T. shell As an independent variable.

5. The multi-position north-finding method for a north-finding instrument according to claim 4, characterized in that: The expression for the compensation model is: Where, ΔB comp The zero-bias drift to be compensated. K is the rate of change of the outer shell temperature, and k1, k2, and k3 are compensation coefficients determined through calibration experiments.

6. The multi-location north-finding method for a north-finding instrument according to claim 1, characterized in that: The multi-position cyclic measurement in step S2 can be switched to continuous rotation measurement mode, specifically as follows: Control the north finder to maintain a constant angular velocity ω r Continuous rotation around the vertical axis; angular velocity ω r The value of satisfies the following constraints: oh min <ω r <ω max in, ω max =min(0.2f Nyquist ,ω limit );f c 1 / f is the noise corner frequency; f Nyquist ω is the Nyquist frequency. limit The upper limit of safe angular velocity is determined based on the gyroscope bandwidth or the system's mechanical resonant frequency; The dynamically adjusted signal processing strategy is phase-locked amplification, which involves generating a reference signal based on the rotational position signal, extracting northward information from the gyroscope signal through coherent integration, and selecting the integration time T of the phase-locked amplification to cover 2-5 rotation cycles. n represents the number of rotation cycles covered by the integration time, ensuring that northbound information is fully extracted.

7. The multi-position north-finding method for a north-finding instrument according to claim 6, characterized in that: The phase-locked amplification solution specifically includes: A sinusoidal reference signal sin(ω) of the same frequency is generated based on the rotational position signal. r t) and cosine reference signal cos(ω) r t); The gyroscope output signal Ω(t) is multiplied by the two reference signals and then integrated to obtain the orthogonal integral components X and Y: Where T is the integration time; The north angle α0 is calculated using the following formula: α0 = arctan 2(Y-Y0, X-X0) Where X0 and Y0 are the zero-bias compensation terms obtained from calibration.

8. The multi-position north-finding method for a north-finding instrument according to claim 1, characterized in that: In step S2, the power spectral density analysis is accelerated by hardware. A field-programmable gate array (FPGA) is used to perform a fast Fourier transform (FFT) on the acquired gyroscope signal. The signal sampling adopts a sliding window method with an overlap rate between windows.

9. The multi-location north-finding method for a north-finding instrument according to claim 1, characterized in that: In step S4, data processing and fusion include: constructing opposing position differences based on azimuth data from multiple measurement locations; and processing the north angle results {α} from M cyclic measurements. 01 α 02 ,...α 0M The variance of the i-th measurement is calculated by weighting and fusing the data based on the dispersion of each measurement. Calculated based on n independent sampling points within this loop: Where, α ij The north angle calculated for the j-th sampling point in the i-th iteration. Assign weights to the mean north-facing angle for this cycle: Where M is the total number of cycles for multi-location cyclic measurement, and n is the number of independent sampling points in a single cycle; Final North Angle:

10. The multi-position north-finding method for a north-finding instrument according to claim 1, characterized in that: During the north-finding calculation, the elevation angle θ of the north-finding instrument is acquired in real time. p With roll angle θ r The local horizontal projection component of the Earth's rotational angular velocity is corrected based on the pitch and roll angles.

Images