Comprehensive pulsar time determination method and device based on time-frequency combination and dynamic weighted fusion

Abstract

The invention discloses a comprehensive pulsar time determination method and device based on time-frequency combination and dynamic weighted fusion, and relates to the field of comprehensive pulsar time calculation. The method comprises the following steps: acquiring a timing residual sequence of each pulsar; the timing residual sequence is determined based on timing residuals collected at set sampling intervals; the timing residual error is the difference value between the actual arrival time and the preset arrival time of the pulsar pulse; carrying out adaptive wavelet packet decomposition and empirical mode decomposition on the timing residual sequence of each pulsar, and carrying out time-frequency combined filtering processing to obtain a processing signal; and performing fusion processing according to the processed signal by adopting a dynamic weighted fusion algorithm to obtain comprehensive pulsar time. According to the invention, dynamic fusion is realized, the noise suppression precision is improved, and the anti-interference capability is improved.

Description

Technical Field

[0001] The present application relates to the field of integrated pulsar time calculation, and in particular to a method and device for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion. Background Art

[0002] The existing methods for calculating integrated pulsar time are mainly based on two technical routes: atomic clock-assisted pulsar time algorithm and pure pulsar time algorithm.

[0003] Atomic clock-assisted pulsar timing algorithms rely on the short-term stability of high-precision atomic clocks (such as hydrogen maser clocks or optical clocks) to calibrate short-term fluctuations in pulsar timing residuals using atomic clock signals. The typical process involves using atomic clocks to generate a short-term time reference and then correcting for the long-term drift of the atomic clocks by combining them with pulsar timing models (such as parameters such as rotation period and period derivatives). For example, the literature proposes fusing atomic clock and pulsar data through Kalman filtering to suppress the influence of pulsar timing noise. However, this approach has the following drawbacks and deficiencies:

[0004] (1) Long-term error accumulation: The long-term drift of atomic clocks (such as frequency aging effects) needs to be corrected by pulsar data, but the pulsar timing residuals themselves contain low-frequency noise (such as red noise). The coupling of the two may cause the correction deviation to accumulate over time.

[0005] (2) Strong external dependence: The system needs to regularly synchronize the atomic clock with pulsar data. Once the atomic clock fails or an abnormality occurs (such as equipment failure in deep space exploration), the continuity of the time base will be threatened.

[0006] (3) Low calibration efficiency: The sampling rates of atomic clocks and pulsar data differ significantly (e.g., atomic clocks are at the nanosecond level, while pulsar observations are at the hour / day level). Real-time fusion requires complex interpolation processing, which introduces additional calculation errors.

[0007] The pure pulsar time algorithm abandons the atomic clock and directly generates a synthetic pulsar by fusing the timing residuals of multiple pulsars. However, its shortcomings and deficiencies are as follows:

[0008] (1) Insufficient suppression of non-stationary noise: The traditional wavelet decomposition has a fixed basis function, which makes it difficult to adapt to the non-stationary noise in the pulsar timing residual (such as time-varying color noise caused by interstellar medium scintillation), resulting in residual high-frequency noise or over-smoothing of the effective signal. For example, interstellar scintillation noise exhibits a power-law spectrum in the low-frequency band and white noise characteristics in the high-frequency band. The fixed wavelet basis cannot dynamically track such characteristics.

[0009] (2) Rigid multi-star fusion strategy: Existing weighted averaging methods often use fixed weights (e.g., static allocation based on the long-term stability of pulsars), without considering short-term disturbance factors in a single observation (e.g., transient disturbances in the interstellar medium and fluctuations in the signal-to-noise ratio of the observation system). For example, if a pulsar is disturbed by the solar wind during a specific period, causing its signal-to-noise ratio to drop sharply, the fixed weight will still give it a high weight, which will significantly reduce the reliability of the fusion result.

[0010] (3) Information separation in the time and frequency domains: Existing methods typically model only in a single dimension, either the time or frequency domain. For example, wavelet decomposition focuses on frequency domain localization but ignores time domain correlation; Bayesian methods rely on time domain statistical properties but do not utilize prior knowledge of frequency domain noise distribution. This separation results in incomplete noise feature extraction, which is particularly limited in error suppression when dealing with mixed color and white noise.

[0011] (4) Lack of autonomy: During long-term operation, system drift (such as instrument thermal noise accumulation) must be corrected by external references (such as regular atomic clock calibration), making it impossible to achieve closed-loop self-consistency. For example, when a deep space probe is far away from the Earth, the transmission delay of the external calibration signal may be as long as several hours, resulting in a loss of real-time performance.

[0012] Whether it is atomic clock-assisted or pure pulsar time algorithm, the core issues can be summarized as follows:

[0013] (1) Mismatch between noise model and real environment: Existing methods assume that the noise is stationary or obeys a fixed statistical distribution, which is significantly different from the non-stationary and time-varying noise characteristics in actual pulsar observations (such as interstellar medium disturbances and instrument nonlinear response).

[0014] (2) Inefficient use of multi-source information: The joint time-frequency characteristics of pulsar timing residuals (such as the time-domain aggregation of high-frequency noise and the frequency-domain correlation of low-frequency drift) are not fully exploited, nor is the real-time linkage between the dynamic changes of observation conditions and weight distribution achieved.

[0015] (3) Lack of long-term autonomous capabilities: The system relies on external intervention or preset model parameters, and is difficult to adapt to the needs of long-term isolated operation in scenarios such as deep space exploration.

[0016] Existing methods have fundamental defects in noise suppression, fusion strategies and system autonomy. Therefore, it is crucial to improve the accuracy of noise suppression and enhance anti-interference capabilities. Summary of the Invention

[0017] The purpose of this application is to provide a method and device for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion, which can realize dynamic fusion, improve the accuracy of noise suppression and enhance anti-interference capability.

[0018] To achieve the above objectives, this application provides the following solutions:

[0019] In a first aspect, the present application provides a method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion, comprising:

[0020] Obtaining a timing residual sequence for each single pulsar; the timing residual sequence is determined based on the timing residuals collected at a set sampling interval; the timing residual is the difference between the actual arrival time of the pulsar pulse and the preset arrival time;

[0021] The timing residual sequence of each single pulsar is subjected to adaptive wavelet packet decomposition and empirical mode decomposition, and then subjected to time-frequency joint filtering to obtain the processed signal;

[0022] A dynamic weighted fusion algorithm is used to perform fusion processing on the processed signals to obtain a comprehensive pulsar time.

[0023] In one embodiment, the method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion further includes:

[0024] Based on a time-frequency domain joint control loop, feedback calibration is performed according to the integrated pulsar time; the time-frequency domain joint control loop includes: a time domain loop for drift correction and a frequency domain loop for noise feature tracking.

[0025] In one embodiment, the timing residual sequence of each single pulsar is subjected to adaptive wavelet packet decomposition and empirical mode decomposition, and then subjected to time-frequency joint filtering to obtain a processed signal, specifically including:

[0026] Determine a wavelet packet basis function; the wavelet packet basis function is determined based on signal spectrum characteristics; the objective function corresponding to the wavelet packet basis function is determined based on a weighted sum of spectrum energy entropy and reconstruction error;

[0027] Performing wavelet packet decomposition on the timing residual sequence of each single pulsar based on the wavelet packet basis function to obtain multiple time-frequency subband signals;

[0028] Perform empirical mode decomposition and noise component identification on each time-frequency subband signal to obtain the noise intrinsic mode function;

[0029] Frequency domain threshold filtering and time domain correlation correction processing are performed based on the noise intrinsic mode function to obtain a processed signal.

[0030] In one embodiment, the objective function corresponding to the wavelet packet basis function is expressed as:

[0031]

[0032] Among them, J(ψ) is the objective function; ψ is the candidate wavelet basis; λ is the balance factor; x is the original time series; is the signal reconstructed based on wavelet packet decomposition; E k is the energy proportion of the kth sub-band.

[0033] In one embodiment, performing empirical mode decomposition and noise component identification on each time-frequency subband signal to obtain a noise intrinsic mode function specifically includes:

[0034] Perform empirical mode decomposition on each time-frequency subband signal to obtain the set of intrinsic mode functions;

[0035] For each eigenmode function in the eigenmode function set, the noise component is identified by frequency domain power spectrum density and time domain kurtosis to obtain the noise eigenmode function.

[0036] In one embodiment, a dynamic weighted fusion algorithm is used to perform fusion processing on the processed signals to obtain a synthetic pulsar, specifically including:

[0037] The reinforcement learning method is used to optimize the fusion weights to obtain the optimized weights;

[0038] Based on the optimized weights, fusion processing is performed on the processed signals to obtain the integrated pulsar time.

[0039] In one embodiment, a reinforcement learning method is used to optimize the fusion weights to obtain optimized weights, specifically including:

[0040] Determine the variance corresponding to the timing residual based on the sliding window and determine the stability index;

[0041] Determine the linear weighted score based on the signal-to-noise ratio and the interstellar medium disturbance index;

[0042] A dynamic weight calculation model is determined according to the stability index and the linear weighted score; the mathematical expression corresponding to the dynamic weight calculation model is:

[0043] ω i (t) = αS i (t)+βC i (t);

[0044] Among them, ω i (t) is the function corresponding to the dynamic weight calculation model; α is the first fusion weight; S i (t) is the stability index; β is the second fusion weight; C i (t) is the linear weighted score;

[0045] The reinforcement learning method is used to optimize the fusion weight based on the state vector and action space by minimizing the reward function to obtain the optimized weight.

[0046] In one embodiment, ARIMA modeling is performed on the integrated pulsar time within the drift-corrected time domain loop, and drift correction is performed based on the drift correction amount;

[0047] After ARIMA modeling of the integrated pulsar time, the corresponding mathematical expression is:

[0048]

[0049] The mathematical expression corresponding to the drift correction amount is:

[0050]

[0051] Among them, φ i is the autoregressive coefficient; B is the backshift operator; B i is the i-th power of the backshift operator B; EPT is the integrated pulsar time; θ j is the moving average coefficient; B j is the jth power of the backshift operator (acting on the time series function to delay it by j units); ε(t) is white noise, that is, noise with zero mean, constant variance and no correlation in the sequence; i and j are both serial numbers; p is the order of the autoregressive part (AR order); q is the order of the moving average part (MA order); ΔT drift is the drift correction amount (i.e. the time drift amount that needs to be corrected); t is the current time; ε(tj) is the white noise value at time tj.

[0052] In one embodiment, within the frequency domain loop of noise feature tracking, the power spectral density of the integrated pulsar is determined, and the number of decomposition levels of the wavelet packet decomposition is adjusted based on a set parameter adjustment rule to enhance low-frequency noise suppression; the mathematical expression corresponding to the set parameter adjustment rule is:

[0053]

[0054] Where flow is the low-frequency start frequency; fmid is the intermediate-frequency cutoff frequency; S(f) is the power spectral density of the integrated pulsar timing residual; f is the frequency; Γ is the low-frequency energy accumulation threshold, and the integral value corresponds to the total noise energy. This frequency band corresponds to the red noise (low-frequency drift noise) in pulsar timing.

[0055] In a second aspect, the present application provides a comprehensive pulsar time determination device based on time-frequency combination and dynamic weighted fusion, comprising:

[0056] A data acquisition module is configured to acquire a timing residual sequence for each individual pulsar; the timing residual sequence is determined based on the timing residuals collected at a set sampling interval; the timing residual is the difference between the actual arrival time of the pulsar pulse and the preset arrival time;

[0057] The processing module is used to perform adaptive wavelet packet decomposition and empirical mode decomposition on the timing residual sequence of each single pulsar, and perform time-frequency joint filtering to obtain a processed signal;

[0058] The fusion processing module is used to adopt a dynamic weighted fusion algorithm to perform fusion processing according to the processed signals to obtain a comprehensive pulsar time.

[0059] According to the specific embodiments provided in this application, this application discloses the following technical effects:

[0060] The present application provides a method and apparatus for determining integrated pulsar time based on time-frequency joint and dynamic weighted fusion. The method first obtains the timing residual sequence of each individual pulsar. The timing residual sequence is determined based on the timing residuals collected at a set sampling interval. The timing residual is the difference between the actual arrival time and the preset arrival time of the pulsar pulse. The timing residual sequence of each individual pulsar is then subjected to adaptive wavelet packet decomposition and empirical mode decomposition, and time-frequency joint filtering is performed to obtain a processed signal. Finally, a dynamic weighted fusion algorithm is used to perform fusion processing based on the processed signal to obtain integrated pulsar time. The present application solves the problems of insufficient noise suppression and rigid weights in existing pulsar time algorithms through time-frequency joint filtering denoising and dynamic weighted fusion, thereby improving the accuracy of noise suppression and enhancing anti-interference capabilities. BRIEF DESCRIPTION OF THE DRAWINGS

[0061] In order to more clearly illustrate the embodiments of the present application or the technical solutions in the prior art, the following briefly introduces the drawings required for use in the embodiments. Obviously, the drawings described below are only some embodiments of the present application. For ordinary technicians in this field, other drawings can be obtained based on these drawings without creative work.

[0062] Figure 1 This is a flow chart of the integrated pulsar time determination method based on time-frequency joint and dynamic weighted fusion;

[0063] Figure 2 This is a structural diagram of the integrated pulsar time determination device based on time-frequency combination and dynamic weighted fusion. DETAILED DESCRIPTION

[0064] The following will be combined with the drawings in the embodiments of this application to clearly and completely describe the technical solutions in the embodiments of this application. Obviously, the embodiments described are only part of the embodiments of this application, not all of the embodiments. Based on the embodiments in this application, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of this application.

[0065] In order to make the above-mentioned purposes, features and advantages of the present application more obvious and easy to understand, the present application is further described in detail below with reference to the accompanying drawings and specific implementation methods.

[0066] In an exemplary embodiment, Figure 1 As shown, a comprehensive pulsar time determination method based on time-frequency combination and dynamic weighted fusion is provided, including:

[0067] Step 100: Obtain a timing residual sequence for each individual pulsar. The timing residual sequence is determined based on the timing residuals collected at a set sampling interval; the timing residual is the difference between the actual arrival time of the pulsar pulse and the preset arrival time.

[0068] Step 200: Perform adaptive wavelet packet decomposition and empirical mode decomposition on the timing residual sequence of each single pulsar, and perform time-frequency joint filtering to obtain a processed signal.

[0069] The timing residual sequence of each single pulsar is subjected to adaptive wavelet packet decomposition and empirical mode decomposition, and then subjected to time-frequency joint filtering to obtain the processed signal, specifically including:

[0070] Step 201: Determine the wavelet packet basis function. The wavelet packet basis function is determined based on the signal spectrum characteristics; the objective function corresponding to the wavelet packet basis function is determined based on the weighted sum of spectrum energy entropy and reconstruction error.

[0071] The expression of the objective function corresponding to the wavelet packet basis function is:

[0072]

[0073] Among them, J(ψ) is the objective function; ψ is the candidate wavelet basis; λ is the balance factor; x is the original time series; is the signal reconstructed based on wavelet packet decomposition; E k is the energy proportion of the kth sub-band.

[0074] Step 202: Perform wavelet packet decomposition on the timing residual sequence of each single pulsar based on the wavelet packet basis function to obtain multiple time-frequency sub-band signals.

[0075] Step 203: Perform empirical mode decomposition and noise component identification on each time-frequency subband signal to obtain a noise intrinsic mode function.

[0076] Among them, empirical mode decomposition and noise component identification are performed on each time-frequency subband signal to obtain the noise intrinsic mode function, which specifically includes:

[0077] Perform empirical mode decomposition on each time-frequency subband signal to obtain a set of intrinsic mode functions.

[0078] For each eigenmode function in the eigenmode function set, the noise component is identified by frequency domain power spectrum density and time domain kurtosis to obtain the noise eigenmode function.

[0079] Step 204: Perform frequency domain threshold filtering and time domain correlation correction processing based on the noise intrinsic mode function to obtain a processed signal.

[0080] Step 300: Using a dynamic weighted fusion algorithm, perform fusion processing based on the processed signals to obtain the integrated pulsar time.

[0081] Among them, the dynamic weighted fusion algorithm is used to perform fusion processing based on the processed signals to obtain the synthetic pulsar, specifically including:

[0082] The reinforcement learning method is used to optimize the fusion weights to obtain the optimized weights; based on the optimized weights, fusion processing is performed according to the processed signals to obtain the integrated pulsar time.

[0083] In one embodiment, a reinforcement learning method is used to optimize the fusion weights to obtain the optimized weights, specifically including:

[0084] Based on the sliding window, the variance corresponding to the timing residual is determined, and the stability index is determined. According to the signal-to-noise ratio and the interstellar medium disturbance index, the linear weighted score is determined. According to the stability index and the linear weighted score, the dynamic weight calculation model is determined. The mathematical expression corresponding to the dynamic weight calculation model is:

[0085] ω i (t) = αS i (t)+βC i (t).

[0086] Among them, ω i (t) is the function corresponding to the dynamic weight calculation model; α is the first fusion weight; S i (t) is the stability index; β is the second fusion weight; C i (t) is the linear weighted score;

[0087] The reinforcement learning method is used to optimize the fusion weight based on the state vector and action space by minimizing the reward function to obtain the optimized weight.

[0088] As an optional implementation, the method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion further includes:

[0089] Based on the time-frequency domain joint control loop, feedback calibration is performed according to the integrated pulsar time; the time-frequency domain joint control loop includes: a time domain loop for drift correction and a frequency domain loop for noise feature tracking.

[0090] Specifically, in the drift-corrected time domain loop, ARIMA modeling is performed on the integrated pulsar time, and drift correction is performed based on the drift correction amount.

[0091] After ARIMA modeling of the integrated pulsar time, the corresponding mathematical expression is:

[0092]

[0093] The mathematical expression corresponding to the drift correction is:

[0094]

[0095] Among them, φ i is the autoregressive coefficient; B is the backshift operator; B i is the i-th power of the backshift operator B; EPT is the integrated pulsar time; θ j is the moving average coefficient; B j is the jth power of the backshift operator (acting on the time series function to delay it by j units); ε(t) is white noise, that is, noise with zero mean, constant variance and no correlation in the sequence; i and j are both serial numbers; p is the order of the autoregressive part (AR order); q is the order of the moving average part (MA order); ΔT drift is the drift correction amount (i.e. the time drift amount that needs to be corrected); t is the current time; ε(tj) is the white noise value at time tj.

[0096] In the frequency domain loop of noise feature tracking, the power spectrum density of the integrated pulsar is determined, and the number of decomposition layers of the wavelet packet decomposition is adjusted based on the set parameter adjustment rules to enhance low-frequency noise suppression. The mathematical expression corresponding to the set parameter adjustment rules is:

[0097]

[0098] Among them, f low is the low-frequency starting frequency; f mid is the intermediate frequency cutoff frequency; S(f) is the power spectral density of the residual of the integrated pulsar timing; f is the frequency; Γ is the low-frequency energy accumulation threshold, and the integral value corresponds to the total noise energy. This frequency band corresponds to the red noise (low-frequency drift noise) in pulsar timing.

[0099] The core of this application is to build a joint time-frequency analysis framework and a dynamic weighted fusion algorithm to improve the accuracy and autonomy of integrated pulsar timing through multi-dimensional information fusion and a closed-loop feedback mechanism. The specific implementation steps are as follows:

[0100] 1. Joint time-frequency modeling and denoising of timing residuals.

[0101] Objective: To solve the problem of non-stationary noise suppression and avoid the limitations of traditional single time domain or frequency domain methods.

[0102] Step 1: Adaptive Wavelet Packet Decomposition (AWPD).

[0103] Input: Single pulsar timing residual sequence x(t), with sampling interval Δt.

[0104] Decomposition process:

[0105] (1) Adaptive basis selection: Dynamically select the optimal wavelet packet basis function based on the signal spectrum characteristics. The objective function is defined as minimizing the weighted sum of the reconstruction error and the band energy entropy:

[0106]

[0107] Among them, ψ is the candidate wavelet basis (such as db4, sym6), E k is the energy proportion of the kth subband, and λ is the balance factor (the default value is 0.7).

[0108] (2) Multi-layer decomposition: Perform 3-layer wavelet packet decomposition on x(t) to obtain 8 time-frequency subband signals

[0109] Step 2: Empirical Mode Decomposition (EMD) assists in noise extraction.

[0110] (1) EMD decomposition: for each subband signal s i (t) Perform EMD to obtain the set of intrinsic mode functions (IMF)

[0111] (2) Noise component identification: The noise-dominated IMF is determined by combining the frequency domain power spectral density (PSD) and the time domain kurtosis:

[0112]

[0113] Where fmax is the high frequency threshold (such as 1Hz), θ f NoiseFlag(c i,j ) is the noise intrinsic mode function; fhigh is the high-frequency noise judgment threshold; PSD(c i,j ) is the eigenmode function c i,j The corresponding power spectrum density; θf is a constant.

[0114] Step 3: Joint time-frequency filtering.

[0115] (1) Frequency domain threshold filtering: Soft threshold processing is performed on the noise IMF:

[0116]

[0117] Among them, the threshold τ(f) is set adaptively according to the sub-band frequency: σ j is the noise standard deviation of the j-th IMF. is the IMF frequency domain coefficient after filtering; sign(c i,j (f)) is the phase sign function; c i,j (f) is the original IMF frequency domain coefficient, and τ(f) is the frequency adaptive threshold.

[0118] Time domain correlation correction: Use the time domain correlation of adjacent sub-band signals to smooth the filtered signal:

[0119]

[0120] Weight ω k Determined by the cross-correlation coefficient between sub-bands. is the smoothed output signal of the i-th sub-band; is the output signal after smoothing the k-th sub-band. Compared with traditional wavelet decomposition, the AWPD-EMD joint algorithm can dynamically adapt to the noise band, and the suppression effect of interstellar scintillation noise will be more significant.

[0121] 2. Dynamic weighted fusion algorithm.

[0122] Objective: To solve the suboptimal problem of fixed-weight fusion and achieve dynamic adjustment of multi-satellite weights according to observation conditions.

[0123] Step 1: Short-term stability index S i (t) Calculation, sliding window Allan variance: in window T w Calculate the Allan variance of the pulsar timing residuals within (e.g. 10 hours):

[0124]

[0125] Take τ = T w Allan variance corresponding to / 2 The reciprocal of is used as the stability index. is the Allan variance; τ is the smoothing time scale, which is the time interval used to calculate the Allan variance; N is the total number of data points in the window; x k+1is the k+1th timing residual value; x k is the kth timing residual value.

[0126] Step 2: Observation condition score, i.e. linear weighted score C i (t)Calculation.

[0127] (1) Scoring model: Constructing a linear weighted scoring model based on real-time observation parameters:

[0128]

[0129] Among them, SNR i (t) is the current signal-to-noise ratio (normalized to [0,1]), D i (t) is the interstellar medium disturbance index (calculated from Faraday rotation measurements), γ1, γ2, and D0 are all constants. γ1 = 0.6, γ2 = 0.4, D0 = 10 rad / m2.

[0130] Step 3: Online optimization of weights (reinforcement learning).

[0131] State definition: state vector s t =[S1(t),C1(t),...,S n (t),C n (t)].

[0132] Among them S i (t) is a short-term stability indicator based on the sliding window Allan variance.

[0133] Action space: Adjust the increments of α and β by Δα∈{-0.1, 0, 0.1}, Δβ=1-Δα.

[0134] Reward function R t :

[0135] R t =-(Var(EPT)+η||s t -s t-1 || 2 ).

[0136] Where Var is the variance, EPT is the integrated pulsar, and η is the weight change penalty coefficient (default is 0.1).

[0137] Dynamic weights improve the stability of fusion results under sudden interference.

[0138] 3. Joint feedback calibration in time and frequency domains.

[0139] Goal: Achieve long-term autonomy and eliminate dependence on external calibration.

[0140] Time domain loop: long-term drift correction.

[0141] Trend fitting: ARIMA modeling of the integrated pulsar time:

[0142]

[0143] Drift correction:

[0144]

[0145] Frequency domain loops: noise signature tracking.

[0146] Power Spectral Density (PSD) feedback: Calculates the PSD of the integrated pulsar and detects abnormal frequency bands (such as 1 / f noise areas).

[0147] Parameter adjustment rules, if:

[0148]

[0149] The number of AWPD decomposition layers is increased to enhance low-frequency noise suppression.

[0150] Closed-loop calibration allows the system to operate without external intervention during a 180-day simulation, with long-term drift error controlled within 0.1ns.

[0151] This application utilizes three core technologies: joint time-frequency modeling (AWPD-EMD), dynamic weight optimization (driven by reinforcement learning), and autonomous feedback calibration (ARIMA+PSD tracking). This achieves efficient suppression of non-stationary noise, real-time optimization of multi-satellite fusion strategies, and closed-loop control of long-term stability. Compared to existing technologies, this system significantly improves key metrics such as noise suppression accuracy, system autonomy, and anti-interference capability. It is particularly suitable for time base generation scenarios in extreme environments such as deep space navigation.

[0152] Specifically, (1) time-frequency joint noise suppression technology.

[0153] A joint denoising framework of adaptive wavelet packet decomposition (AWPD) and empirical mode decomposition (EMD) is proposed to solve the problem of non-stationary noise suppression through the following methods:

[0154] 1) Dynamic wavelet basis selection: Adaptively matching the pulsar noise band characteristics based on the joint optimization objective of spectral energy entropy and reconstruction error.

[0155] 2) Time-frequency domain judgment conditions: For the IMF components after EMD decomposition, the noise components are determined by combining the frequency domain power ratio and the time domain kurtosis, thus avoiding the risk of misjudgment by traditional single criteria.

[0156] (2) Dynamic weighted fusion mechanism.

[0157] A dynamic weight model based on short-term stability and real-time observation conditions is constructed, and reinforcement learning is introduced to optimize weight parameters to achieve real-time adaptation of multi-satellite fusion strategies:

[0158] Weight factor formulation and online parameter optimization: Reinforcement learning algorithms, such as Q-learning, are used to adjust α and β in real time to minimize the variance of the fusion results and the weight mutation (reward function).

[0159] (3) Autonomous closed-loop feedback calibration system.

[0160] Design a joint time-frequency domain control loop to achieve long-term stability autonomous maintenance.

[0161] Long-term drift correction in the time domain: Predict system drift based on the ARIMA model and generate real-time corrections.

[0162] Frequency domain noise tracking: Through power spectral density (PSD) monitoring, the number of AWPD decomposition layers and EMD threshold parameters are dynamically adjusted to suppress low-frequency 1 / f noise.

[0163] Specifically, the pulsar timing residual processing method based on AWPD-EMD joint denoising includes the following steps:

[0164] (1) Adaptive wavelet packet decomposition (AWPD) is performed on the original timing residual signal, and the band energy entropy is optimized by the objective function Select a wavelet basis.

[0165] (2) Perform empirical mode decomposition (EMD) on the decomposed sub-band signals and extract the intrinsic mode function (IMF).

[0166] (3) The noise-dominated IMF component is determined jointly based on the frequency domain power ratio and time domain kurtosis.

[0167] (4) Apply frequency domain soft threshold filtering and time domain correlation smoothing to the noise IMF to reconstruct the denoised timing residual signal.

[0168] Dynamic weight fusion method and reinforcement learning optimization:

[0169] Dynamic weight calculation model ω i (t) = αS i (t)+βC i (t), where S i (t) is a short-term stability indicator based on the sliding window Allan variance, Score the real-time observation conditions.

[0170] The online parameter optimization method based on Q-learning includes the state definition state vector s t=[S1(t),C1(t),...,S n (t),C n (t)], action space (Δα, Δβ), and reward function R t =-(Var(EPT)+η||s t -s t-1 || 2 ).

[0171] Joint time-frequency feedback calibration:

[0172] The time domain calibration module uses the ARIMA model to predict long-term drift:

[0173]

[0174] And correct the system time base.

[0175] Frequency domain calibration module: detects noise characteristics through power spectral density analysis and dynamically adjusts the number of AWPD decomposition layers and EMD threshold parameters.

[0176] Feedback control logic: When , triggers the AWPD parameter adjustment instruction.

[0177] The method mentioned in this application can be applied in deep space navigation, gravitational wave detection, and autonomous time reference generation systems.

[0178] This application addresses the problems of insufficient noise suppression, rigid weights, and strong external dependencies in existing pulsar timing algorithms through three technological innovations: joint time-frequency denoising, dynamic weight fusion, and autonomous feedback calibration. This application covers methods, systems (devices), and application scenarios, forming a multi-layered defense system that provides core technical support for the generation of high-precision time bases in fields such as deep space exploration.

[0179] In an exemplary embodiment, Figure 2 As shown, a comprehensive pulsar time determination device based on time-frequency combination and dynamic weighted fusion is provided, comprising:

[0180] The data acquisition module is used to obtain the timing residual sequence of each individual pulsar. The timing residual sequence is determined based on the timing residuals collected at a set sampling interval; the timing residual is the difference between the actual arrival time of the pulsar pulse and the preset arrival time.

[0181] The processing module is used to perform adaptive wavelet packet decomposition and empirical mode decomposition on the timing residual sequence of each single pulsar, and perform time-frequency joint filtering to obtain a processed signal.

[0182] The fusion processing module is used to adopt a dynamic weighted fusion algorithm to perform fusion processing based on the processed signals to obtain the integrated pulsar time.

[0183] The technical features of the above embodiments can be combined arbitrarily. To make the description concise, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, they should be considered to be within the scope of this specification.

[0184] This document uses specific examples to illustrate the principles and implementation methods of this application. The description of the above examples is only intended to help understand the method and core concept of this application. At the same time, for those skilled in the art, based on the concept of this application, there may be changes in the specific implementation methods and application scope. In summary, the content of this specification should not be understood as limiting this application.

Claims

1. A method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion, characterized in that: The method for determining the integrated pulsar time based on time-frequency combination and dynamic weighted fusion includes: Obtaining a timing residual sequence for each single pulsar; the timing residual sequence is determined based on the timing residuals collected at a set sampling interval; the timing residual is the difference between the actual arrival time of the pulsar pulse and the preset arrival time; The timing residual sequence of each single pulsar is subjected to adaptive wavelet packet decomposition and empirical mode decomposition, and then subjected to time-frequency joint filtering to obtain the processed signal; A dynamic weighted fusion algorithm is used to perform fusion processing on the processed signals to obtain a comprehensive pulsar time.

2. The method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion according to claim 1, characterized in that: The method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion further includes: Based on a time-frequency domain joint control loop, feedback calibration is performed according to the integrated pulsar time; the time-frequency domain joint control loop includes: a time domain loop for drift correction and a frequency domain loop for noise feature tracking.

3. The method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion according to claim 1, characterized in that: The timing residual sequence of each single pulsar is subjected to adaptive wavelet packet decomposition and empirical mode decomposition, and then subjected to time-frequency joint filtering to obtain the processed signal, which specifically includes: Determine a wavelet packet basis function; the wavelet packet basis function is determined based on signal spectrum characteristics; the objective function corresponding to the wavelet packet basis function is determined based on a weighted sum of spectrum energy entropy and reconstruction error; Performing wavelet packet decomposition on the timing residual sequence of each single pulsar based on the wavelet packet basis function to obtain multiple time-frequency subband signals; Perform empirical mode decomposition and noise component identification on each time-frequency subband signal to obtain the noise intrinsic mode function; Frequency domain threshold filtering and time domain correlation correction processing are performed based on the noise intrinsic mode function to obtain a processed signal.

4. The method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion according to claim 3, characterized in that: The expression of the objective function corresponding to the wavelet packet basis function is: Among them, J(ψ) is the objective function; ψ is the candidate wavelet basis; λ is the balance factor; x is the original time series; is the signal reconstructed based on wavelet packet decomposition; E k is the energy proportion of the kth sub-band.

5. The method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion according to claim 3, characterized in that: Perform empirical mode decomposition and noise component identification on each time-frequency subband signal to obtain the noise intrinsic mode function, which specifically includes: Perform empirical mode decomposition on each time-frequency subband signal to obtain the set of intrinsic mode functions; For each eigenmode function in the eigenmode function set, the noise component is identified by frequency domain power spectrum density and time domain kurtosis to obtain the noise eigenmode function.

6. The method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion according to claim 1, characterized in that: The method adopts a dynamic weighted fusion algorithm to perform fusion processing on the processed signals to obtain a synthetic pulsar, specifically including: The reinforcement learning method is used to optimize the fusion weights to obtain the optimized weights; Based on the optimized weights, fusion processing is performed on the processed signals to obtain the integrated pulsar time.

7. The method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion according to claim 6, characterized in that: The reinforcement learning method is used to optimize the fusion weights to obtain the optimized weights, including: Determine the variance corresponding to the timing residual based on the sliding window and determine the stability index; Determine the linear weighted score based on the signal-to-noise ratio and the interstellar medium disturbance index; A dynamic weight calculation model is determined according to the stability index and the linear weighted score; the mathematical expression corresponding to the dynamic weight calculation model is: ω i (t)=αS i (t)+βC i (t); Among them, ω i (t) is the function corresponding to the dynamic weight calculation model; α is the first fusion weight; S i (t) is the stability index; β is the second fusion weight; C i (t) is the linear weighted score; The reinforcement learning method is used to optimize the fusion weight based on the state vector and action space by minimizing the reward function to obtain the optimized weight.

8. The method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion according to claim 2, characterized in that: In the drift-corrected time domain loop, ARIMA modeling is performed on the integrated pulsar time, and drift correction is performed based on the drift correction amount. After ARIMA modeling of the integrated pulsar time, the corresponding mathematical expression is: The mathematical expression corresponding to the drift correction amount is: Among them, φ i is the autoregressive coefficient; B is the backshift operator; B i is the i-th power of the backshift operator B; EPT is the integrated pulsar time; θ j is the moving average coefficient; B j is the jth power of the backshift operator; ε(t) is white noise; i and j are both serial numbers; p is the order of the autoregressive part; q is the order of the moving average part; ΔT drift is the drift correction; t is the current time; ε(tj) is the white noise value at time tj.

9. The method for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion according to claim 2, characterized in that: In the frequency domain loop of noise feature tracking, the power spectrum density of the integrated pulsar is determined, and the number of decomposition layers of the wavelet packet decomposition is adjusted based on the set parameter adjustment rules to enhance low-frequency noise suppression. The mathematical expression corresponding to the set parameter adjustment rules is: Where flow is the low-frequency start frequency; fmid is the intermediate-frequency cutoff frequency; S(f) is the power spectral density of the integrated pulsar residual; f is the frequency; and Γ is the low-frequency energy accumulation threshold.

10. A device for determining integrated pulsar time based on time-frequency combination and dynamic weighted fusion, characterized in that: The integrated pulsar time determination device based on time-frequency combination and dynamic weighted fusion includes: A data acquisition module is configured to acquire a timing residual sequence for each individual pulsar; the timing residual sequence is determined based on the timing residuals collected at a set sampling interval; the timing residual is the difference between the actual arrival time of the pulsar pulse and the preset arrival time; The processing module is used to perform adaptive wavelet packet decomposition and empirical mode decomposition on the timing residual sequence of each single pulsar, and perform time-frequency joint filtering to obtain a processed signal; The fusion processing module is used to adopt a dynamic weighted fusion algorithm to perform fusion processing according to the processed signals to obtain a comprehensive pulsar time.

Images