A satellite formation system fusion positioning optimization method based on genetic algorithm

Abstract

A satellite formation system fusion positioning optimization method based on genetic algorithm belongs to the field of satellite formation navigation and positioning technology. The present invention solves the problems of low accuracy of existing space-based positioning methods of satellite formation systems and limitations of inter-satellite observations. The method of the present invention is as follows: Step 1. Kalman filtering is used to determine the integer ambiguity based on inter-satellite observation data, and the integer ambiguity is substituted into the baseline observation equation to obtain the processed baseline observation value; Step 2. The initial position and clock error of the user star are determined by extended Kalman filtering according to the pseudo-range observation data; the pseudo-range observation value is obtained based on the pseudo-range observation equation and the clock error; Step 3. According to the navigation star position, the baseline observation value obtained in step 1, the pseudo-range observation value obtained in step 2 and the initial position of the user star, multi-observation information fusion and intelligent optimization are performed through genetic algorithm to obtain the user star position information. The method of the present invention can be applied to the fusion positioning and optimization of satellite formation systems.

Description

Technical Field

[0001] The invention belongs to the technical field of satellite formation navigation and positioning, and in particular relates to a satellite formation system fusion positioning optimization method based on a genetic algorithm. Background Art

[0002] The satellite formation system is a distributed multi-body formation system with a reference satellite as the center and multiple user satellites maintaining a specific relative configuration. The system achieves functional reorganization of multiple satellites through inter-satellite collaboration. It can not only complete most traditional single-satellite tasks, but also realize some advanced aerospace missions, such as space long-baseline interferometry, large-range stereo imaging, intermittent positioning and navigation, etc., and its performance is superior to that of traditional single-satellite systems.

[0003] In the satellite formation flight mission, formation capture, configuration maintenance, formation reconstruction and other technologies all need to be based on the position and velocity information of the satellite in space. Therefore, high-precision satellite navigation positioning is one of the important guarantees for the design of the satellite formation system. Satellite navigation positioning is to observe several parameters related to the satellite movement at several times or in several time intervals, and then calculate the position, velocity or orbital elements of the satellite according to the orbital dynamics principle and specific algorithms. The traditional satellite positioning method is to obtain a set of distance observation values ​​by observing a single satellite at least four ground observation stations at the same time, and then estimate the position and velocity through estimation methods such as least squares and Kalman filtering. Finally, the satellite's orbital element data is obtained through a certain coordinate transformation. The observation value of this method itself has low accuracy, usually at the m or dm level, which cannot meet the needs of high-precision positioning tasks. The phase difference between the time when the signal is sent and the time when it is received is measured by carrier phase ranging technology to perform high-precision inter-satellite ranging. The ranging accuracy can usually reach cm or even mm level, but there are usually integer ambiguity parameters in the ranging process, and a more complex algorithm is required to estimate the integer ambiguity.

[0004] In summary, the space-based positioning of satellite formations can obtain the absolute position of user satellites in space, but the accuracy is low, and although intersatellite relative observation has high measurement accuracy, it is difficult to obtain the absolute position of satellites. Therefore, there is an urgent need for a satellite formation positioning and optimization method that solves the problems of low accuracy of existing space-based positioning and limitations of intersatellite observation. Summary of the invention

[0005] The purpose of the present invention is to solve the problems of low accuracy of existing space-based positioning methods for satellite formation systems and limitations of inter-satellite observations, and to propose a satellite formation system fusion positioning optimization method based on genetic algorithms.

[0006] The technical solution adopted by the present invention to solve the above technical problems is:

[0007] A satellite formation system fusion positioning optimization method based on genetic algorithm, the method specifically comprises the following steps:

[0008] Step 1: Use pseudorange observation data and integer ambiguity to construct the state variables of Kalman filtering, use pseudorange observation data and baseline observation data to construct the observation variables of Kalman filtering, then determine the integer ambiguity through Kalman filtering, and substitute the integer ambiguity into the baseline observation equation to obtain the processed baseline observation value;

[0009] Step 2: The state variables of the extended Kalman filter are constructed using the user satellite position and its change rate, the clock difference between the user satellite and the navigation satellite and its change rate, and the observation variables of the extended Kalman filter are constructed using the pseudo-range observation data. Subsequently, the initial position and clock difference of the user satellite are determined by the extended Kalman filter.

[0010] Then, based on the pseudorange observation equation and clock error, the processed pseudorange observation value is obtained;

[0011] Step 3: According to the navigation satellite position, the baseline measurement value obtained in step 1, the pseudo-range observation value obtained in step 2 and the initial position of the user satellite, the user satellite position information is obtained by fusing multi-observation information through a genetic algorithm;

[0012] The specific process of step 3 is as follows:

[0013] Step 301: Given the optimization rate λ, the maximum number of iterations n, and the maximum number of baseline optimizations n L , the population size m, the crossover rate α, and the mutation rate β, and initialize the elements in the baseline optimization list to all baselines;

[0014] Input the navigation satellite position, the baseline observation data obtained in step 1, the user satellite initial position obtained in step 2, and the pseudo-range observation data into the genetic algorithm, and initialize the user satellite initial position to the optimal position;

[0015] Step 302: for any baseline in the baseline optimization list, calculate the Euclidean distance between the two user satellites connected to the baseline, subtract the observed value of the baseline from the calculated Euclidean distance, and use the difference as the baseline error;

[0016] Similarly, after processing each baseline in the baseline optimization list, the baseline corresponding to the largest baseline error is selected;

[0017] The selected baseline is recorded as the target baseline flagL, and the selected baseline is removed from the baseline optimization list. The number of times the baseline is optimized flaggn is initialized to 0, and the maximum baseline error is taken as the maximum error, recorded as L flagL; The two user satellites connected to the selected baseline are taken as target satellites, and the numbers of the obtained target satellites are recorded as flag1 and flag2;

[0018] Step 303: Initialize the population according to the coding rules;

[0019] Use the cost function to calculate the cost of each individual in the population, take the individual with the smallest cost as the optimal individual, take the cost of the optimal individual as the optimal cost, and record the initial optimal cost before optimization as error_ori;

[0020] Step 304, set the number of iterations to 1;

[0021] Step 305: determine whether the number of iterations is greater than the maximum number of iterations n;

[0022] If the number of iterations is not greater than n, execute step 306;

[0023] If the number of iterations is greater than n, the iteration process ends and the current optimal position is output;

[0024] Step 306, increment the number of times the current baseline is optimized, flagn, by 1;

[0025] Step 307, traverse the individuals in the population, and set the individual number a=1;

[0026] Step 308: determine whether the individual number is greater than m;

[0027] If the individual number is not greater than m, execute step 309;

[0028] If the individual number is greater than m, execute step 311;

[0029] Step 309, calculating the cost of the ath individual, and decoding the genetic code of the ath individual to obtain an adjustment amount;

[0030] The positions of the target satellites flag1 and flag2 are adjusted by using the adjustment amount to obtain the adjusted positions of the target satellites flag1 and flag2, the Euclidean distance between the target satellite flag1 after the position adjustment and each navigation star is calculated, and the pseudo-range calculation values ​​corresponding to the target satellite flag1 are obtained respectively, the Euclidean distance between the target satellite flag2 after the position adjustment and each navigation star is calculated, and the pseudo-range calculation values ​​corresponding to the target satellite flag2 are obtained respectively;

[0031] For any pseudorange calculation value of the target satellite flag1, the difference between the pseudorange calculation value and the corresponding pseudorange observation value is calculated, and the obtained difference is used as a pseudorange error value of the target satellite flag1. Similarly, each pseudorange error value of the target satellite flag1 and each pseudorange error value of the target satellite flag2 are obtained;

[0032] Step 310: judging whether the ath individual meets the individual optimal condition according to the cost of the ath individual and the pseudorange error value calculated in step 309;

[0033] If the ath individual satisfies the individual optimal condition, the ath individual is updated to the optimal individual, the cost of the ath individual is updated to the optimal cost, and the individual number a is incremented by 1, and then the process returns to step 308;

[0034] If the ath individual does not meet the individual optimal condition, the individual number a is directly incremented by 1, and then the process returns to step 308;

[0035] Step 311, determine whether the cost of the optimal individual is less than error_ori*λ;

[0036] If the cost of the optimal individual is less than error_ori*λ, the iteration ends and the final position is output;

[0037] Otherwise, continue to execute step 312;

[0038] Step 312: determine whether the optimal individual meets the single optimization condition;

[0039] If it is satisfied, the genetic code of the optimal individual is decoded, and the position of the target satellites flag1 and flag2 is adjusted using the adjustment amount obtained by decoding to obtain the adjusted optimal position; the errors of each baseline in the baseline optimization list are recalculated, and the target baseline flagL is updated according to the recalculated baseline errors, and the baseline corresponding to the maximum baseline error in the recalculated baseline errors is removed from the baseline optimization list, and the maximum error L is updated. flagL , target satellite flag1 and flag2, and set the current baseline optimization times flagn to 0; then execute step 313;

[0040] If not satisfied, directly execute step 313;

[0041] Step 313: Determine whether the current baseline optimization times flagn is greater than n L ,

[0042] If the current baseline optimization times is greater than n L, then the target satellite connected to the current target baseline is no longer optimized, the first element is selected from the baseline optimization list, the first element selected is used as the target baseline flagL, and it is removed from the baseline optimization list, and then the target satellite flag1 and flag2 are updated, and the current baseline optimization times flagn is set to 0; it is determined whether the baseline optimization list is empty, if it is empty, the elements of the baseline optimization list are reinitialized to all baselines, if it is not empty, no processing is required; then step 314 is executed;

[0043] If not greater than n L , then execute step 314;

[0044] Step 314: cross the population to obtain a crossover population, and then perform a mutation operation on the crossover population to obtain a mutated population;

[0045] Step 315: Initialize a pseudo-random number for each individual in the mutated population. If the pseudo-random number of an individual is less than half of the crossover probability, reinitialize the individual, otherwise do not perform any operation.

[0046] After reinitialization, a new population is obtained;

[0047] Step 316: Initialize a pseudo-random number for each individual in the new population. If the pseudo-random number of an individual is less than 0.05, initialize the individual as a local better individual near the best individual. If the pseudo-random number of an individual is greater than or equal to 0.05, do not operate on the individual.

[0048] After the operation, an updated population is obtained, the number of iterations is incremented by 1, and the process returns to step 305 .

[0049] Furthermore, the specific process of determining the integer ambiguity is:

[0050] Construct the state variable X of the Kalman filter:

[0051] X=[ρ1 ρ2 … ρ L N1 N2 … N L ] T

[0052] Where ρ1 represents the first pseudorange observation value, ρ2 represents the second pseudorange observation value, and ρ L represents the Lth pseudorange observation value, N1 represents the first integer ambiguity, N2 represents the second integer ambiguity, and N L represents the Lth integer ambiguity, L = M(M-1) / 2, M is the number of user satellites;

[0053] Construct the observed variable Y of Kalman filtering:

[0054]

[0055] In the formula, represents the first baseline observation, represents the second baseline observation, represents the Lth baseline observation;

[0056] The Kalman filter formula is used to solve the state variable X to obtain the integer ambiguity N1, N2, ..., N L :

[0057]

[0058] P(k+1|k)=FP(k)F T +ΓQΓ T

[0059] K(k+1)=P(k+1|k)H T [HP(k+1|k)H T +R] -1

[0060]

[0061] P(k+1)=[IK(k+1)H]P(k+1|k)

[0062] In the formula, is the state variable estimate of the baseline observation model at time k, F is the state transfer matrix of the baseline observation model, and the superscript T represents the transpose of the matrix. is the one-step prediction of the state variables of the baseline observation model at time k, Γ is the input noise driving matrix of the baseline observation model, P(k) is the state quantity covariance estimate of the baseline observation model at time k, Q is the system model noise matrix, H is the observation matrix of the baseline observation model, P(k+1|k) is the one-step prediction of the state quantity covariance of the baseline observation model at time k, R is the baseline observation model noise matrix, and the superscript -1 represents the inverse of the matrix. is the state variable estimate of the baseline observation model at time k+1, K(k+1) is the gain of the baseline observation model at time k+1, Y(k+1) is the observed variable of the baseline observation model at time k+1, P(k+1) is the state variable covariance estimate of the baseline observation model at time k+1, and I is the unit matrix.

[0063] Furthermore, the integer ambiguity is substituted into the baseline observation equation to obtain the processed baseline observation value, which is specifically:

[0064]

[0065] Where, d ijrepresents the processed baseline observation value between the i-th user star and the j-th user star, where i = 1, 2, ..., M, j = 1, 2, ..., M, x j ,y j ,z j is the coordinate value of user star j in the three-dimensional rectangular coordinate system, x i ,y i ,z i is the coordinate value of user star i in the three-dimensional rectangular coordinate system, c is the speed of light constant, t i is the clock difference between user star i and navigation star, t j is the clock difference between user star j and navigation star, N is the integer ambiguity corresponding to the baseline observation value between the i-th user star and the j-th user star, λ is the wavelength of the observation signal, ε ij is the baseline observation noise between user star i and user star j.

[0066] Furthermore, the method of determining the initial position and clock error of the user satellite is:

[0067] The state vector X′ of the extended Kalman filter is as follows:

[0068]

[0069] In the formula, x, y, and z represent the three-dimensional coordinates of the user satellite. represents the three-dimensional velocity of the user satellite, t represents the clock difference between the user satellite and the navigation satellite, Represents the rate of change of the clock difference between the user satellite and the navigation satellite;

[0070] The three-dimensional coordinates of the user satellite, i.e. the initial position of the user satellite, are obtained by solving the extended Kalman filter;

[0071]

[0072] P′(k+1|k)=Φ(k)P′(k)Φ T (k)+Γ′Q′Γ′ T

[0073] K′(k+1)=P′(k+1|k)H′ T (k+1)[H′(k+1|k)P′(k+1|k)H′ T (k+1|k)+R′] -1

[0074]

[0075] P′(k+1)=[IK′(k+1)H′(k+1|k)]P′(k+1|k)

[0076] in, is the final position estimation result at time k, f[·] is the linear system model, is the one-step prediction of the state quantity of the pseudorange observation model at time k, Φ(k) is the state transfer matrix at time k, the superscript T represents the transpose of the matrix, P′(k) is the state quantity covariance estimate of the pseudorange observation model at time k, Q′ is the input noise covariance matrix of the pseudorange observation model, Γ′ is the input noise driving matrix of the pseudorange observation model, P′(k+1|k) is the one-step prediction of the state quantity covariance of the pseudorange observation model at time k, H′(k+1) is the Jacobi matrix of the pseudorange observation equation at time k, R′ is the covariance matrix of the observation noise of the pseudorange observation model, the superscript -1 represents the inverse of the matrix, h[·] is the pseudorange observation equation, is the state variable estimate of the pseudorange observation model at time k+1, K′(k+1) is the gain of the pseudorange observation model at time k+1, P′(k+1) is the state variable covariance estimate of the pseudorange observation model at time k+1, I is the unit matrix, and Y′(k+1) is the observation variable of the pseudorange observation model at time k+1.

[0077] Furthermore, the Jacobi matrix of the pseudorange observation equation is:

[0078]

[0079] in, x l ,y l 、z l represents the position of the lth navigation star, r l Represents (x l ,y l ,z l ) and (x, y, z), q represents the total number of pseudo-range observations between the user satellite with three-dimensional coordinates (x, y, z) and the navigation satellite, that is, the total number of navigation satellites.

[0080] Furthermore, the specific process of obtaining the processed pseudorange observation value based on the pseudorange observation equation and the clock error is as follows:

[0081]

[0082] In the formula, represents the pseudo-range observation value between the i′th navigation satellite and the jth user satellite, x i′ ,y i′ ,z i′ is the three-dimensional coordinate value of the i′th navigation star, x j ,y j ,z j is the three-dimensional coordinate value of the jth user star, t j is the clock difference between the user star and the navigation star, is the distance observation noise between the i′th navigation star and the jth user star.

[0083] Furthermore, the cost f cost (X0,ΔX best_old ) is calculated as:

[0084]

[0085] Where X0 represents the user satellite position before adjustment, ΔX best_old represents the adjustment amount obtained after individual decoding, represents the baseline measurement value between the i-th and j-th user satellites obtained in step 1, d i,j (X0,ΔX best_old ) indicates the use of ΔX best_old The Euclidean distance between the i-th user star and the j-th user star obtained after adjusting X0;

[0086] After updating the position using the individual adjustment amount, the difference between the baseline calculation value and the baseline measurement value after Kalman filtering between any two user satellites is calculated, and then the absolute values ​​of all the calculated differences are taken respectively, and the absolute values ​​are summed to obtain the cost function value.

[0087] Furthermore, the individual optimal condition is:

[0088]

[0089] Where ΔX best_new is the adjustment amount corresponding to the new optimal individual, f cost (X0′,ΔX best_new ) is the cost of the new individual obtained after adjusting X0′, X0′ is the position of the original optimal individual, f cost (X0,ΔX best_old ) is the cost of the original optimal individual, is the pseudorange observation value obtained in step 2, The adjustment amount ΔX best_new After updating the position X0′ of user star j, the Euclidean distance between the position of user star j and the position of navigation star i′ is adjusted; The adjustment amount ΔX best_old After updating the position X0 of user star j, the Euclidean distance between the position of user star j and the position of navigation star i′ is adjusted; Δρ is the pseudorange observation error value, and q is the number of navigation stars.

[0090] Furthermore, the single optimization condition is:

[0091] f cost (X k,ΔX k+1 )<f cost (X k ,0)-0.01*(0.5) k′ L flagL

[0092] Among them, X k is the user star position at the end of the kth iteration, ΔX k+1 is the adjustment amount obtained by optimal individual decoding at the k+1th iteration, f cost (X k ,ΔX k+1 ) is to use ΔX k+1 X k The cost of making adjustments, f cost (X k , 0) is to use 0 to X k The cost of the adjustment, k′ represents the number of iterations divided by 20 and rounded down.

[0093] Furthermore, the crossover operation on the population is specifically performed as follows:

[0094] Initialize a pseudo-random number for each individual in the population, compare the pseudo-random number with the crossover probability, if the pseudo-random number of an individual is less than the crossover probability, then put the individual into the crossover population, if the pseudo-random number of an individual is greater than or equal to the crossover probability, then put the individual into the non-crossover population;

[0095] Perform crossover on the individuals in the crossover population. Assume that the number of individuals in the crossover population is m. c , if m c is an even number, then select the ath individual and the mth individual c -a individuals participate in crossover in pairs, if m c If is an odd number, the individuals in the middle of the population container do not participate in the crossover;

[0096] When crossing, use random points to cross. First, initialize n p Pseudo-random numbers, the obtained n p After multiplying the pseudo-random numbers by the individual code length, the obtained n p The product is rounded down, and the rounded n p The values ​​are used as crossover points, and then the individuals participating in the crossover in pairs are crossed according to the crossover points to obtain the crossover population after the crossover.

[0097] Merge the crossover population and the non-crossover population to obtain the crossover population;

[0098] The mutation of the population after the crossover is specifically as follows:

[0099] For any individual in the population after crossover, a pseudo-random number is randomly initialized for each binary code in the encoding of the individual. If the pseudo-random number of a binary code is less than the mutation rate, the binary code is inverted. If the pseudo-random number of a binary code is greater than or equal to the mutation rate, no mutation operation is performed. After processing each individual in the population after crossover, the mutated population is finally obtained.

[0100] For any individual in the mutated population, a pseudo-random number is randomly initialized for the individual. If the pseudo-random number is less than 0.05, the individual is initialized as a local better individual near the best individual. Otherwise, if the pseudo-random number is greater than or equal to 0.05, the individual does not need to be processed. Similarly, each individual in the mutated population is processed separately.

[0101] Initialize the individual as a local better individual near the optimal individual, which is specifically:

[0102] Generate an adjustment amount with the first n_low bit value of 0, and use the adjustment amount to adjust the optimal individual to obtain a local better individual near the optimal individual.

[0103] The beneficial effects of the present invention are:

[0104] (1) The present invention combines extended Kalman filtering with genetic algorithms to fuse and optimize multi-source information of multi-observation data, thus realizing a high-precision fusion positioning optimization method for a satellite formation system. The method of the present invention has strong optimization and search capabilities, and introduces the probabilistic concept in nature. The search process is more in line with objective laws. Compared with traditional positioning algorithms, the method of the present invention can achieve higher positioning accuracy while obtaining the absolute position of the user satellite. It also solves the problems of low accuracy of existing space-based positioning methods and limitations of inter-satellite observations, and is more suitable for scenarios where actual satellite-borne computers perform high-precision position estimation;

[0105] (2) The present invention fuses and optimizes multiple observation data such as pseudorange and baseline, and has high-precision positioning and online intelligent optimization functions. At the same time, the algorithms are relatively separated and can improve the positioning accuracy of the satellite formation system while ensuring the stability of the algorithm through intelligent search.

[0106] (3) The present invention introduces baseline observation values ​​into the space-based pseudo-range positioning results of the satellite formation and performs fusion positioning and optimization based on genetic algorithms. Compared with traditional positioning methods, the positioning accuracy of the satellite formation is significantly improved.

[0107] (4) The fusion positioning optimization algorithm of the present invention has high positioning accuracy and strong scalability. It can adapt to positioning tasks of different satellite formation sizes. In addition, the positioning accuracy is significantly improved compared with traditional positioning methods under undersampling conditions, and its engineering application value is higher.

[0108] The experimental results show that the variance of the X-axis positioning error is reduced by 84.5% under full sampling, the variance of the Y-axis positioning error is reduced by 85.7%, and the variance of the Z-axis positioning error is reduced by 83.5%. Under undersampling, the variance of the X-axis positioning error is reduced by 59.6%, the variance of the Y-axis positioning error is reduced by 66.3%, and the variance of the Z-axis positioning error is reduced by 65.5%. BRIEF DESCRIPTION OF THE DRAWINGS

[0109] Figure 1 It is a schematic diagram of the method of the present invention;

[0110] Figure 2 is a flow chart of the method of the present invention;

[0111] Figure 3a It is a schematic diagram of the X-axis positioning error under undersampling working conditions;

[0112] Figure 3b It is a schematic diagram of the Y-axis positioning error under undersampling working conditions;

[0113] Figure 3c This is a schematic diagram of the Z-axis positioning error under undersampling conditions;

[0114] In the figure, the dotted line is the error of the traditional Kalman filter positioning result, and the solid line is the error of the fusion positioning optimization result based on the genetic algorithm;

[0115] Figure 4a This is a schematic diagram of the X-axis positioning error under full sampling working conditions;

[0116] Figure 4b This is a schematic diagram of the Y-axis positioning error under full sampling working conditions;

[0117] Figure 4c This is a schematic diagram of the Z-axis positioning error under full sampling working conditions;

[0118] In the figure, the dotted line is the error of the traditional Kalman filter positioning result, and the solid line is the error of the fusion positioning optimization result based on the genetic algorithm. DETAILED DESCRIPTION

[0119] Specific implementation method 1. Combination Figure 1 and Figure 2 The present embodiment describes a satellite formation system fusion positioning optimization method based on a genetic algorithm, and the method specifically includes the following steps:

[0120] Step 1: Use pseudorange observation data and integer ambiguity to construct the state variables of Kalman filtering, use pseudorange observation data and baseline observation data to construct the observation variables of Kalman filtering, then determine the integer ambiguity through Kalman filtering, and substitute the integer ambiguity into the baseline observation equation to obtain the processed baseline observation value;

[0121] Step 2: construct the state variables of the extended Kalman filter using the user star position and its (referring to the user star position) change rate, the user star and navigation star clock difference and its (referring to the user star and navigation star clock difference) change rate, construct the observation variables of the extended Kalman filter using the pseudorange observation data, and then determine the user star initial position and clock difference through the extended Kalman filter;

[0122] Then, based on the pseudorange observation equation and clock error, the processed pseudorange observation value is obtained;

[0123] Step 3: Based on the known position of the navigation satellite, the baseline measurement value obtained in step 1, the pseudo-range observation value obtained in step 2, and the initial position of the user satellite, the user satellite position information is obtained by performing multi-observation information (multi-pseudo-range observation values, multi-baseline observation values) fusion and intelligent optimization through a genetic algorithm;

[0124] The specific process of step 3 is as follows:

[0125] Step 301: Given the optimization rate λ, the maximum number of iterations n, and the maximum number of baseline optimizations n L , the population size m, the crossover rate α, and the mutation rate β, and initialize the elements in the baseline optimization list to all baselines;

[0126] Input the navigation satellite position, the baseline observation data obtained in step 1, the user satellite initial position obtained in step 2, and the pseudo-range observation data into the genetic algorithm, and initialize the user satellite initial position to the optimal position;

[0127] Step 302: for any baseline in the baseline optimization list, calculate the Euclidean distance between the two user satellites connected to the baseline, subtract the observed value of the baseline from the calculated Euclidean distance, and use the difference as the baseline error;

[0128] Similarly, after processing each baseline in the baseline optimization list, the baseline corresponding to the largest baseline error is selected;

[0129] The selected baseline is recorded as the target baseline flagL, and the selected baseline is removed from the baseline optimization list. The number of times the baseline is optimized flaggn is initialized to 0, and the maximum baseline error is taken as the maximum error, recorded as L flagL; The two user satellites connected to the selected baseline are taken as target satellites, and the numbers of the obtained target satellites are recorded as flag1 and flag2;

[0130] Step 303: Initialize the population according to the coding rules;

[0131] Use the cost function to calculate the cost of each individual in the population, take the individual with the smallest cost as the optimal individual, take the cost of the optimal individual as the optimal cost, and record the initial optimal cost before optimization (that is, the optimal cost at this time) as error_ori;

[0132] Step 304, set the number of iterations to 1;

[0133] Step 305: determine whether the number of iterations is greater than the maximum number of iterations n;

[0134] If the number of iterations is not greater than n, execute step 306;

[0135] If the number of iterations is greater than n, the iteration process ends and the current optimal position is output;

[0136] Step 306, increment the number of times the current baseline is optimized, flagn, by 1;

[0137] Step 307, traverse the individuals in the population (after numbering the individuals, traverse from small to large numbers), and set the individual number a=1;

[0138] Step 308: determine whether the individual number is greater than m;

[0139] If the individual number is not greater than m, execute step 309;

[0140] If the individual number is greater than m, execute step 311;

[0141] Step 309, calculating the cost of the ath individual, and decoding the genetic code of the ath individual to obtain an adjustment amount;

[0142] Use the adjustment amount to adjust the positions of the target satellites flag1 and flag2 (the starting point of the adjusted position is the user satellite position after the last population optimization), obtain the adjusted positions of the target satellites flag1 and flag2, calculate the Euclidean distance between the target satellite flag1 after the position adjustment and each navigation star, and obtain the pseudo-range calculation values ​​corresponding to the target satellite flag1, calculate the Euclidean distance between the target satellite flag2 after the position adjustment and each navigation star, and obtain the pseudo-range calculation values ​​corresponding to the target satellite flag2;

[0143] For any pseudorange calculation value of the target satellite flag1, the difference between the pseudorange calculation value and the corresponding pseudorange observation value (referring to the pseudorange observation value between the user satellite and the navigation satellite corresponding to the pseudorange calculation value) is calculated, and the obtained difference is used as a pseudorange error value of the target satellite flag1. Similarly, each pseudorange error value of the target satellite flag1 and each pseudorange error value of the target satellite flag2 are obtained;

[0144] Step 310: judging whether the ath individual meets the individual optimal condition according to the cost of the ath individual and the pseudorange error value calculated in step 309;

[0145] If the ath individual satisfies the individual optimal condition, the ath individual is updated to the optimal individual, the cost of the ath individual is updated to the optimal cost, and the individual number a is incremented by 1, and then the process returns to step 308;

[0146] If the ath individual does not meet the individual optimal condition, the individual number a is directly incremented by 1, and then the process returns to step 308;

[0147] Step 311, determine whether the cost of the optimal individual is less than error_ori*λ;

[0148] If the cost of the optimal individual is less than error_ori*λ, the iteration ends and the final position is output;

[0149] Otherwise, continue to execute step 312;

[0150] Step 312: determine whether the optimal individual meets the single optimization condition;

[0151] If it is satisfied, the genetic code of the optimal individual is decoded, and the position of the target satellites flag1 and flag2 is adjusted using the adjustment amount obtained by decoding (the starting point of the adjusted position is the user satellite position after the last population optimization), and the adjusted optimal position is obtained; the error of each baseline in the baseline optimization list is recalculated, and the target baseline flagL is updated according to the recalculated baseline error, and the baseline corresponding to the maximum baseline error in the recalculated baseline error is removed from the baseline optimization list, and the maximum error L is updated. flagL , target satellite flag1 and flag2, and set the current baseline optimization times flagn to 0; then execute step 313;

[0152] If not satisfied, directly execute step 313;

[0153] Step 313: Determine whether the current baseline optimization times flagn is greater than n L ,

[0154] If the current baseline optimization times is greater than n L, then the target satellite connected to the current target baseline will no longer be optimized, the first element will be selected from the baseline optimization list, the first element will be used as the target baseline flagL, and it will be removed from the baseline optimization list, and then the target satellite flag1 and flag2 will be updated, and the current baseline optimization times flagn will be set to 0;

[0155] Determine whether the baseline optimization list is empty. If it is empty, reinitialize the elements of the baseline optimization list to all baselines. If it is not empty, no processing is required. Then execute step 314.

[0156] If not greater than n L , then execute step 314;

[0157] Step 314: cross the population to obtain a crossover population, and then perform a mutation operation on the crossover population to obtain a mutated population;

[0158] Step 315: Initialize a pseudo-random number for each individual in the mutated population. If the pseudo-random number of an individual is less than half of the crossover probability, reinitialize the individual, otherwise do not perform any operation.

[0159] After reinitialization, a new population is obtained;

[0160] Step 316: Initialize a pseudo-random number for each individual in the new population. If the pseudo-random number of an individual is less than 0.05, initialize the individual as a local better individual near the best individual. If the pseudo-random number of an individual is greater than or equal to 0.05, do not operate on the individual.

[0161] After the operation, an updated population is obtained, the number of iterations is incremented by 1, and the process returns to step 305 .

[0162] The present invention fuses and optimizes the baseline observation data and pseudo-range observation data through multi-source information, realizes a satellite formation system positioning algorithm with high-precision positioning and online intelligent optimization functions, and can complete high-precision estimation and optimization of the kinematic information of any satellite in the formation, which has high engineering application significance. The implementation process of the method of the present invention is shown in Table 1:

[0163] Table 1

[0164]

[0165]

[0166] Specific implementation method 2: This implementation method is different from the specific implementation method 1 in that the specific process of determining the integer ambiguity is:

[0167] Construct the state variable X of the Kalman filter:

[0168] X=[ρ1ρ2…ρ L N1 N2…N L ] T

[0169] Where ρ1 represents the first pseudorange observation value, ρ2 represents the second pseudorange observation value, and ρ L represents the Lth pseudorange observation value, N1 represents the first integer ambiguity, N2 represents the second integer ambiguity, and N L represents the Lth integer ambiguity, L = M(M-1) / 2, M is the number of user satellites;

[0170] Select L pseudo-range observations from all pseudo-range observations to participate in the calculation, and all the selected pseudo-range observations are pseudo-range observations between user satellites;

[0171] Construct the observed variable Y of Kalman filtering:

[0172]

[0173] In the formula, represents the first baseline observation, represents the second baseline observation, represents the Lth baseline observation;

[0174] The Kalman filter formula is used to solve the state variable X to obtain the integer ambiguity N1, N2, ..., N L :

[0175]

[0176] P(k+1|k)=FP(k)F T +ΓQΓ T

[0177] K(k+1)=P(k+1|k)H T [HP(k+1|k)H T +R] -1

[0178]

[0179] P(k+1)=[IK(k+1)H]P(k+1|k)

[0180] In the formula, is the state variable estimate of the baseline observation model at time k, F is the state transfer matrix of the baseline observation model, and the superscript T represents the transpose of the matrix. is the one-step prediction of the state variables of the baseline observation model at time k, Γ is the input noise driving matrix of the baseline observation model, P(k) is the state quantity covariance estimate of the baseline observation model at time k, Q is the system model (constant speed model) noise matrix, H is the observation matrix of the baseline observation model, P(k+1|k) is the one-step prediction of the state quantity covariance of the baseline observation model at time k, R is the baseline observation model noise matrix, and the superscript -1 represents the inverse of the matrix. is the state variable estimate of the baseline observation model at time k+1, K(k+1) is the gain of the baseline observation model at time k+1, Y(k+1) is the observed variable of the baseline observation model at time k+1, P(k+1) is the state variable covariance estimate of the baseline observation model at time k+1, and I is the unit matrix.

[0181] State transfer matrix F:

[0182]

[0183] In the formula, Represents the number of combinations;

[0184] Observation matrix:

[0185]

[0186] System model noise matrix Q:

[0187]

[0188] Observation model noise matrix R:

[0189]

[0190] Initial covariance matrix P:

[0191]

[0192] The other steps and parameters are the same as those in the first embodiment.

[0193] Specific implementation method three: This implementation method is different from specific implementation method one or two in that the integer ambiguity is substituted into the baseline observation equation to obtain the processed baseline observation value, which is specifically:

[0194]

[0195] Where, d ij represents the processed baseline observation value between the i-th user star and the j-th user star, where i = 1, 2, ..., M, j = 1, 2, ..., M, x j ,y j ,z jis the coordinate value of user star j in the three-dimensional rectangular coordinate system, x i ,y i ,z i is the coordinate value of user star i in the three-dimensional rectangular coordinate system, c is the speed of light constant, c = 299792458m / s, t i is the clock difference between user star i and navigation star, t j is the clock difference between user star j and navigation star, N is the integer ambiguity corresponding to the baseline observation value between the i-th user star and the j-th user star, λ is the wavelength of the observation signal, ε ij is the baseline observation noise between user star i and user star j, usually Gaussian white noise.

[0196] The other steps and parameters are the same as those in the first or second embodiment.

[0197] Specific implementation method 4: This implementation method is different from one of the specific implementation methods 1 to 3 in that the method for determining the initial position and clock error of the user satellite is:

[0198] The state vector X′ of the extended Kalman filter is as follows:

[0199]

[0200] In the formula, x, y, and z represent the three-dimensional coordinates of the user satellite. represents the three-dimensional velocity of the user satellite, t represents the clock difference between the user satellite and the navigation satellite, Represents the rate of change of the clock difference between the user satellite and the navigation satellite;

[0201] The three-dimensional coordinates of the user satellite, i.e. the initial position of the user satellite, are obtained by solving the extended Kalman filter (EKF);

[0202]

[0203] P′(k+1|k)=Φ(k)P′(k)Φ T (k)+Γ′Q′Γ′ T

[0204] K′(k+1)=P′(k+1|k)H′ T (k+1)[H′(k+1|k)P′(k+1|k)H′ T (k+1|k)+R′] -1

[0205]

[0206] P′(k+1)=[IK′(k+1)H′(k+1|k)]P′(k+1|k)

[0207] in, is the final position estimation result at time k (for the first moment that needs to be positioned, Just take it as 0 vector. When you need to calculate the value at time k+1, based on the value at time k That is, f[·] is a linear system model, is the one-step prediction of the state quantity of the pseudorange observation model at time k, Φ(k) is the state transfer matrix at time k, the superscript T represents the transpose of the matrix, P′(k) is the state quantity covariance estimate of the pseudorange observation model at time k, Q′ is the input noise covariance matrix of the pseudorange observation model, Γ′ is the input noise driving matrix of the pseudorange observation model, P′(k+1|k) is the one-step prediction of the state quantity covariance of the pseudorange observation model at time k, H′(k+1) is the Jacobi matrix of the pseudorange observation equation at time k, R′ is the covariance matrix of the observation noise of the pseudorange observation model, the superscript -1 represents the inverse of the matrix, h[·] is the pseudorange observation equation, is the state variable estimate of the pseudorange observation model at time k+1, K′(k+1) is the gain of the pseudorange observation model at time k+1, P′(k+1) is the state variable covariance estimate of the pseudorange observation model at time k+1, I is the unit matrix, and Y′(k+1) is the observation variable of the pseudorange observation model at time k+1.

[0208] State transition matrix:

[0209]

[0210] In the formula, dt represents the time interval between two adjacent observation epochs, and q represents the number of observation equations at the same time;

[0211] Input noise matrix:

[0212]

[0213] In the formula, Indicates the three-dimensional acceleration of the user satellite, in m / s 2 , is Gaussian white noise.

[0214] The input noise driving matrix Γ′ of the pseudorange observation model is:

[0215]

[0216] The covariance matrix Q′ of the input noise matrix W is:

[0217]

[0218] The observation noise covariance matrix R′ of the pseudorange observation model is a q×q diagonal matrix:

[0219]

[0220] Initial covariance matrix P0:

[0221]

[0222] The other steps and parameters are the same as those in Specific Embodiments 1 to 3.

[0223] Specific implementation method 5: This implementation method is different from any one of specific implementation methods 1 to 4 in that the Jacobi matrix of the pseudorange observation equation is:

[0224]

[0225] in, x l ,y l 、z l represents the position of the lth navigation star, r l Represents (x l ,y l ,z l ) and (x, y, z), q represents the total number of pseudo-range observations between the user satellite with three-dimensional coordinates (x, y, z) and the navigation satellite, that is, the total number of navigation satellites.

[0226] The number of pseudo-range observation values ​​involved in the calculation of this implementation is q, and the q pseudo-range observation values ​​involved in the calculation of this implementation are all pseudo-range observation values ​​between the navigation satellite and the user satellite.

[0227] The other steps and parameters are the same as those in Specific Embodiments 1 to 4.

[0228] Specific implementation method 6: This implementation method is different from the specific implementation methods 1 to 5 in that the processed pseudo-range observation value is obtained based on the pseudo-range observation equation and the clock difference, and the specific process is:

[0229]

[0230] In the formula, represents the pseudo-range observation value between the i′th navigation satellite and the jth user satellite, x i′ ,y i′ ,z i′ is the three-dimensional coordinate value of the i′th navigation star, x j ,y j ,z j is the three-dimensional coordinate value of the jth user star, t j is the clock difference between the user star and the navigation star, is the distance observation noise between the i′th navigation star and the jth user star.

[0231] The other steps and parameters are the same as those in Specific Implementation Methods 1 to 5.

[0232] Specific implementation method 7: This implementation method is different from the specific implementation methods 1 to 6 in that the cost f cost (X0,ΔX best_old ) is calculated as:

[0233]

[0234] Where X0 represents the user satellite position before adjustment (the first iteration is the initial position, otherwise it is the intermediate position variable before the iterative adjustment is completed during the iteration process), ΔX best_old represents the adjustment amount obtained after individual decoding, represents the baseline measurement value between the i-th and j-th user satellites obtained in step 1, d i,j (X0,ΔX best_old ) indicates the use of ΔX best_old The Euclidean distance between the i-th user star and the j-th user star obtained after adjusting X0;

[0235] After updating the position using the individual adjustment amount, the difference between the baseline calculation value and the baseline measurement value after Kalman filtering between any two user satellites is calculated, and then the absolute values ​​of all the calculated differences are taken respectively, and the absolute values ​​are summed to obtain the cost function value.

[0236] The other steps and parameters are the same as those in Specific Embodiments 1 to 6.

[0237] Specific implementation eight: This implementation differs from any one of specific implementations one to seven in that the individual optimal conditions are:

[0238]

[0239] Where ΔX best_new is the adjustment amount corresponding to the new optimal individual, f cost (X′0,ΔX best_new ) is the cost of the new individual obtained after adjusting X′0, X′0 is the position of the original optimal individual (using ΔX best_old After adjusting X0), f cost (X0,ΔX best_old ) is the cost of the original optimal individual, is the pseudorange observation value obtained in step 2, The adjustment amount ΔX best_new After updating the position X′0 of user star j, the Euclidean distance between the position of user star j and the position of navigation star i′ is adjusted; The adjustment amount ΔX best_old After updating the position X0 of user star j, the Euclidean distance between the position of user star j and the position of navigation star i′ is adjusted; Δρ is the pseudorange observation error value (a rough observation error of the pseudorange observation equipment, and it is approximately assumed that the observation error of each pseudorange is the same), and q is the number of navigation stars.

[0240] When making the optimal individual judgment, it is necessary to satisfy the requirement that the cost of the individual is less than the cost of the optimal individual, and that after the user satellite position is adjusted using the individual adjustment amount, the error between all pseudorange calculation values ​​and observation values ​​of the target satellite is reduced or less than the pseudorange observation error.

[0241] The other steps and parameters are the same as those in Specific Embodiments 1 to 7.

[0242] Specific implementation method 9: This implementation method is different from any one of specific implementation methods 1 to 8 in that the single optimization condition is:

[0243] f cost (X k ,ΔX k+1 )<f cost (X k ,0)-0.01*(0.5) k′ L flagL

[0244] Among them, X k is the user star position at the end of the kth iteration, ΔX k+1 is the adjustment amount obtained by optimal individual decoding at the k+1th iteration, f cost (X k ,ΔX k+1 ) is to use ΔX k+1 X k The cost of making adjustments, f cost (X k , 0) is to use 0 to X k The cost of adjustment, k′, represents the number of iterations divided by 20 and rounded down, and increases with the number of iterations.

[0245] Through a single optimization condition, it is judged whether the adjustment of the target satellite has achieved the expected effect, so as to decide whether to replace the target satellite. At the same time, the k′ value increases with the increase of the number of iterations, which not only avoids the problem of falling into the local optimal solution due to excessive adjustment in one time during the iteration process, but also avoids the phenomenon of the optimization process being stuck when the optimization effect is not obvious when the baseline errors are small in the later stage of the iteration.

[0246] The other steps and parameters are the same as those in Specific Embodiments 1 to 8.

[0247] Specific implementation method 10: This implementation method is different from the specific implementation methods 1 to 9 in that the crossover operation is performed on the population, which is specifically:

[0248] Initialize a pseudo-random number for each individual in the population, compare the pseudo-random number with the crossover probability, if the pseudo-random number of an individual is less than the crossover probability, then put the individual into the crossover population, if the pseudo-random number of an individual is greater than or equal to the crossover probability, then put the individual into the non-crossover population;

[0249] Perform crossover on the individuals in the crossover population. Assume that the number of individuals in the crossover population is m. c , if m c is an even number, then select the ath individual and the mth individual c -a individuals participate in crossover in pairs, if m c If is an odd number, the individuals in the middle of the population container do not participate in the crossover;

[0250] When crossing, use random points to cross. First, initialize n p Pseudo-random numbers, the obtained n p After multiplying the pseudo-random numbers by the individual code length, the obtained n p The product is rounded down, and the rounded n p The values ​​are used as crossover points, and then the individuals participating in the crossover in pairs are crossed according to the crossover points to obtain the crossover population after the crossover.

[0251] Merge the crossover population and the non-crossover population to obtain the crossover population;

[0252] The mutation of the population after the crossover is specifically as follows:

[0253] For any individual in the population after the crossover, a pseudo-random number is randomly initialized for each binary code in the encoding of the individual. If the pseudo-random number of a binary code is less than the mutation rate, the binary code is inverted (that is, if the binary code is 1 before the mutation, it is changed to 0, otherwise, if the binary code is 0 before the mutation, it is changed to 1). If the pseudo-random number of a binary code is greater than or equal to the mutation rate, no mutation operation is performed. After processing each individual in the population after the crossover, the mutated population is finally obtained;

[0254] For any individual in the mutated population, a pseudo-random number is randomly initialized for the individual. If the pseudo-random number is less than 0.05, the individual is initialized as a local better individual near the best individual. Otherwise, if the pseudo-random number is greater than or equal to 0.05, the individual does not need to be processed. Similarly, each individual in the mutated population is processed separately.

[0255] Initialize the individual as a local better individual near the optimal individual, which is specifically:

[0256] Generate an adjustment amount whose first n_low digits are 0 (that is, some high-order digits of the adjustment amount are 0, and the low-order digits are assigned values, and the number of high-order digits n_low with a value of 0 can be taken as half of the total number of digits of the entire adjustment amount), and use the adjustment amount to adjust the optimal individual to obtain a locally superior individual near the optimal individual.

[0257] The other steps and parameters are the same as those in Specific Embodiments 1 to 9.

[0258] The encoding method of the present invention is:

[0259] The optimization variables are encoded using binary coding. The three-dimensional coordinates of two satellites are optimized each time, with a total of six variables to be optimized. The adjustment amount of the optimization variable is recorded as ΔX = [Δx flag1 Δy flag1 Δz flag1 Δx flag2 Δy flag2 Δz flag2 ], where the subscript represents the satellite number to be adjusted. The adjustment amount of each variable to be optimized is encoded using an 8-bit binary code, where the first 7 bits of the 8-bit binary code are digital bits and the last bit is a sign bit, that is, the encoding length of the adjustment amount corresponding to each individual in the population is 48 bits;

[0260] The decoding method of the present invention is:

[0261] Divide the 48-bit binary code corresponding to the adjustment amount into 6 8-bit binary codes, convert the binary code into a decimal number according to the rules of the digit and sign bits, and then divide the obtained decimal number by 100 to obtain the adjustment amount with an accuracy of 0.01m and a range of [-1.27, 1.27];

[0262] For each variable to be optimized, the relationship between the encoded binary number and the decoded decimal number is:

[0263]

[0264] Among them, value D Indicates the decimal number obtained after decoding, value B Represents the encoded binary number, value B [i] represents the i-th digit of a binary number.

[0265] Experimental Section

[0266] Taking the satellite formation system under space-based navigation as an example, the navigation satellite is an observation satellite with known ephemeris in the global navigation and positioning system, and the user satellite is a satellite to be positioned in the formation. The basic principle of the satellite formation system fusion positioning optimization method based on genetic algorithm designed by the present invention is as follows: Figure 1 As shown in the figure, the baseline observation data is first preprocessed by Kalman filtering to obtain the baseline observation value, and then the pseudo-range observation data is subjected to extended Kalman filtering to obtain the initial position and pseudo-range observation value. Finally, the processed baseline and pseudo-range observation value are used as the reference value, and the initial position of step 2 is used as the initial value for genetic algorithm fusion positioning and optimization. The genetic algorithm optimization process is as follows: Figure 2 As shown, a high-precision position estimate is output when the baseline error is less than a given value.

[0267] For the high-precision fusion positioning optimization method of the satellite formation system based on genetic algorithm designed by the present invention, the simulation is completed through the STK / MATLAB multi-platform simulation environment. A simulation scene is established in STK and a satellite object is inserted. First, the relevant data in the simulation environment is read into the MATLAB environment through the relevant interface for algorithm processing, and then the calculation results are input into the STK scene for verification. During the simulation, 4 satellites are selected from the scene as space-based navigation stars with known position information, and 8 satellites are selected as user satellites to be positioned with unknown position information, that is, the number of baseline observation values ​​is 28 during full sampling, and the number of baseline observation values ​​is less than 28 during undersampling. For data settings, the pseudorange observation value noise is set to 2, and the baseline observation value noise is set to 0.005.

[0268] Satellite No. 4 is selected from the 8 user satellites for simulation result analysis. The simulation results under undersampling are as follows: Figure 3a , Figure 3b and Figure 3c As shown, the simulation results under full sampling working condition are as follows Figure 4a , Figure 4b and Figure 4c As shown. From the simulation results, it can be seen that the high-precision fusion positioning optimization method for satellite formation system based on genetic algorithm designed by the present invention has significantly improved the three-dimensional positioning accuracy compared with the traditional positioning algorithm, and the error covariance is also significantly reduced. The variance of the X-axis positioning error is reduced by 84.5% under full sampling, the variance of the Y-axis positioning error is reduced by 85.7%, and the variance of the Z-axis positioning error is reduced by 83.5%. Under undersampling, the variance of the X-axis positioning error is reduced by 59.6%, the variance of the Y-axis positioning error is reduced by 66.3%, and the variance of the Z-axis positioning error is reduced by 65.5%, indicating that the high-precision fusion positioning optimization method for satellite formation system based on genetic algorithm of the present invention has achieved the expected effect.

[0269] The above calculation examples of the present invention are only used to explain the calculation model and calculation process of the present invention in detail, and are not intended to limit the implementation methods of the present invention. For ordinary technicians in the relevant field, other different forms of changes or modifications can be made based on the above description. It is impossible to list all the implementation methods here. All obvious changes or modifications derived from the technical solution of the present invention are still within the scope of protection of the present invention.

Claims

1. A satellite formation system fusion positioning optimization method based on genetic algorithm, characterized in that: The method specifically comprises the following steps: Step 1: Use pseudorange observation data and integer ambiguity to construct the state variables of Kalman filtering, use pseudorange observation data and baseline observation data to construct the observation variables of Kalman filtering, then determine the integer ambiguity through Kalman filtering, and substitute the integer ambiguity into the baseline observation equation to obtain the processed baseline observation value; Step 2: The state variables of the extended Kalman filter are constructed using the user satellite position and its change rate, the clock difference between the user satellite and the navigation satellite and its change rate, and the observation variables of the extended Kalman filter are constructed using the pseudo-range observation data. Subsequently, the initial position and clock difference of the user satellite are determined by the extended Kalman filter. Then, based on the pseudorange observation equation and clock error, the processed pseudorange observation value is obtained; Step 3: According to the navigation satellite position, the baseline measurement value obtained in step 1, the pseudo-range observation value obtained in step 2 and the initial position of the user satellite, the user satellite position information is obtained by fusing multi-observation information through a genetic algorithm; The specific process of step 3 is as follows: Step 301: Given the optimization rate λ, the maximum number of iterations n, and the maximum number of baseline optimizations n L , the population size m, the crossover rate α, and the mutation rate β, and initialize the elements in the baseline optimization list to all baselines; Input the navigation satellite position, the baseline observation data obtained in step 1, the user satellite initial position obtained in step 2, and the pseudo-range observation data into the genetic algorithm, and initialize the user satellite initial position to the optimal position; Step 302: for any baseline in the baseline optimization list, calculate the Euclidean distance between the two user satellites connected to the baseline, subtract the observed value of the baseline from the calculated Euclidean distance, and use the difference as the baseline error; Similarly, after processing each baseline in the baseline optimization list, the baseline corresponding to the largest baseline error is selected; The selected baseline is recorded as the target baseline flagL, and the selected baseline is removed from the baseline optimization list. The number of times the baseline is optimized flaggn is initialized to 0, and the maximum baseline error is taken as the maximum error, recorded as L flagL ; The two user satellites connected to the selected baseline are taken as target satellites, and the numbers of the obtained target satellites are recorded as flag1 and flag2; Step 303: Initialize the population according to the coding rules; Use the cost function to calculate the cost of each individual in the population, take the individual with the smallest cost as the optimal individual, take the cost of the optimal individual as the optimal cost, and record the initial optimal cost before optimization as error_ori; Step 304, set the number of iterations to 1; Step 305: determine whether the number of iterations is greater than the maximum number of iterations n; If the number of iterations is not greater than n, execute step 306; If the number of iterations is greater than n, the iteration process ends and the current optimal position is output; Step 306, increment the number of times the current baseline is optimized, flagn, by 1; Step 307, traverse the individuals in the population, and set the individual number a=1; Step 308: determine whether the individual number is greater than m; If the individual number is not greater than m, execute step 309; If the individual number is greater than m, execute step 311; Step 309, calculating the cost of the ath individual, and decoding the genetic code of the ath individual to obtain an adjustment amount; The positions of the target satellites flag1 and flag2 are adjusted by using the adjustment amount to obtain the adjusted positions of the target satellites flag1 and flag2, the Euclidean distance between the target satellite flag1 after the position adjustment and each navigation star is calculated, and the pseudo-range calculation values ​​corresponding to the target satellite flag1 are obtained respectively, the Euclidean distance between the target satellite flag2 after the position adjustment and each navigation star is calculated, and the pseudo-range calculation values ​​corresponding to the target satellite flag2 are obtained respectively; For any pseudorange calculation value of the target satellite flag1, the difference between the pseudorange calculation value and the corresponding pseudorange observation value is calculated, and the obtained difference is used as a pseudorange error value of the target satellite flag1. Similarly, each pseudorange error value of the target satellite flag1 and each pseudorange error value of the target satellite flag2 are obtained; Step 310: judging whether the ath individual meets the individual optimal condition according to the cost of the ath individual and the pseudorange error value calculated in step 309; If the ath individual satisfies the individual optimal condition, the ath individual is updated to the optimal individual, the cost of the ath individual is updated to the optimal cost, and the individual number a is incremented by 1, and then the process returns to step 308; If the ath individual does not meet the individual optimal condition, the individual number a is directly incremented by 1, and then the process returns to step 308; Step 311, determine whether the cost of the optimal individual is less than error_ori*λ; If the cost of the optimal individual is less than error_ori*λ, the iteration ends and the final position is output; Otherwise, continue to execute step 312; Step 312: determine whether the optimal individual meets the single optimization condition; If it is satisfied, the genetic code of the optimal individual is decoded, and the position of the target satellites flag1 and flag2 is adjusted using the adjustment amount obtained by decoding to obtain the adjusted optimal position; the errors of each baseline in the baseline optimization list are recalculated, and the target baseline flagL is updated according to the recalculated baseline errors, and the baseline corresponding to the maximum baseline error in the recalculated baseline errors is removed from the baseline optimization list, and the maximum error L is updated. flagL , target satellite flag1 and flag2, and set the current baseline optimization times flagn to 0; then execute step 313; If not satisfied, directly execute step 313; Step 313: Determine whether the current baseline optimization times flagn is greater than n L , If the current baseline optimization times is greater than n L , the target satellite connected to the current target baseline will no longer be optimized, the first element will be selected from the baseline optimization list, and the first element selected will be used as the target baseline flag L , and remove it from the baseline optimization list, then update the target satellite flag1 and flag2, and set the current baseline optimization times flagn to 0; determine whether the baseline optimization list is empty, if it is empty, reinitialize the elements of the baseline optimization list to all baselines, if it is not empty, no processing is required; then execute step 314; If not greater than n L , then execute step 314; Step 314: cross the population to obtain a crossover population, and then perform a mutation operation on the crossover population to obtain a mutated population; Step 315: Initialize a pseudo-random number for each individual in the mutated population. If the pseudo-random number of an individual is less than half of the crossover probability, reinitialize the individual, otherwise do not perform any operation. After reinitialization, a new population is obtained; Step 316: Initialize a pseudo-random number for each individual in the new population. If the pseudo-random number of an individual is less than 0.05, initialize the individual as a local better individual near the best individual. If the pseudo-random number of an individual is greater than or equal to 0.05, do not operate on the individual. After the operation, an updated population is obtained, the number of iterations is incremented by 1, and the process returns to step 305 .

2. The satellite formation system fusion positioning optimization method based on genetic algorithm according to claim 1 is characterized in that: The specific process of determining the integer ambiguity by Kalman filtering is: Construct the state variable X of the Kalman filter: X=[ρ1ρ2…ρ L N1 N2…N L ] T Where ρ1 represents the first pseudorange observation value, ρ2 represents the second pseudorange observation value, and ρ L represents the Lth pseudorange observation value, N1 represents the first integer ambiguity, N2 represents the second integer ambiguity, and N L represents the Lth integer ambiguity, L = M(M-1) / 2, M is the number of user satellites; Construct the observed variable Y of Kalman filtering: In the formula, represents the first baseline observation, represents the second baseline observation, represents the Lth baseline observation; The Kalman filter formula is used to solve the state variable X to obtain the integer ambiguity N1, N2, ..., N L : P(k+1|k)=FP(k)F T +ΓQΓ T K(k+1)=P(k+1|k)H T [HP(k+1|k)H T +R] -1 P(k+1)=[IK(k+1)H]P(k+1|k) In the formula, is the state variable estimate of the baseline observation model at time k, F is the state transfer matrix of the baseline observation model, and the superscript T represents the transpose of the matrix. is the one-step prediction of the state variables of the baseline observation model at time k, Γ is the input noise driving matrix of the baseline observation model, P(k) is the state quantity covariance estimate of the baseline observation model at time k, Q is the system model noise matrix, H is the observation matrix of the baseline observation model, P(k+1|k) is the one-step prediction of the state quantity covariance of the baseline observation model at time k, R is the baseline observation model noise matrix, and the superscript -1 represents the inverse of the matrix. is the state variable estimate of the baseline observation model at time k+1, K(k+1) is the gain of the baseline observation model at time k+1, Y(k+1) is the observed variable of the baseline observation model at time k+1, P(k+1) is the state variable covariance estimate of the baseline observation model at time k+1, and I is the unit matrix.

3. The satellite formation system fusion positioning optimization method based on genetic algorithm according to claim 2 is characterized in that: Substituting the integer ambiguity into the baseline observation equation to obtain the processed baseline observation value is specifically: Where, d ij represents the processed baseline observation value between the i-th user star and the j-th user star, where i = 1, 2, ..., M, j = 1, 2, ..., M, x j ,y j ,z j is the coordinate value of user star j in the three-dimensional rectangular coordinate system, x i ,y i ,z i is the coordinate value of user star i in the three-dimensional rectangular coordinate system, c is the speed of light constant, t i is the clock difference between user star i and navigation star, t j is the clock difference between user star j and navigation star, N is the integer ambiguity corresponding to the baseline observation value between the i-th user star and the j-th user star, λ is the wavelength of the observation signal, ε ij is the baseline observation noise between user star i and user star j.

4. The method for optimizing satellite formation system fusion positioning based on genetic algorithm according to claim 3 is characterized in that: The method for determining the initial position and clock error of the user satellite is: The state vector X′ of the extended Kalman filter is as follows: In the formula, x, y, and z represent the three-dimensional coordinates of the user satellite. represents the three-dimensional velocity of the user satellite, t represents the clock difference between the user satellite and the navigation satellite, Represents the rate of change of the clock difference between the user satellite and the navigation satellite; The three-dimensional coordinates of the user satellite, i.e. the initial position of the user satellite, are obtained by solving the extended Kalman filter; P′(k+1|k)=Φ(k)P′(k)Φ T (k)+Γ′Q′Γ′ T K′(k+1)=P′(k+1|k)H′ T (k+1)[H′(k+1|k)P′(k+1|k)H′ T (k+1|k)+R′] -1 P′(k+1)=[IK′(k+1)H′(k+1|k)]P′(k+1|k) in, is the final position estimation result at time k, f[·] is the linear system model, is the one-step prediction of the state quantity of the pseudorange observation model at time k, Φ(k) is the state transfer matrix at time k, the superscript T represents the transpose of the matrix, P′(k) is the state quantity covariance estimate of the pseudorange observation model at time k, Q′ is the input noise covariance matrix of the pseudorange observation model, Γ′ is the input noise driving matrix of the pseudorange observation model, P′(k+1|k) is the one-step prediction of the state quantity covariance of the pseudorange observation model at time k, H′(k+1) is the Jacobi matrix of the pseudorange observation equation at time k, R′ is the covariance matrix of the observation noise of the pseudorange observation model, the superscript -1 represents the inverse of the matrix, h[·] is the pseudorange observation equation, is the state variable estimate of the pseudorange observation model at time k+1, K′(k+1) is the gain of the pseudorange observation model at time k+1, P′(k+1) is the state variable covariance estimate of the pseudorange observation model at time k+1, I is the unit matrix, and Y′(k+1) is the observation variable of the pseudorange observation model at time k+1.

5. The method for optimizing satellite formation system fusion positioning based on genetic algorithm according to claim 4, characterized in that: The Jacobi matrix of the pseudorange observation equation is: in, x l ,y l 、z l represents the position of the lth navigation star, r l Represents (x l ,y l ,z l ) and (x, y, z), q represents the total number of pseudo-range observations between the user satellite with three-dimensional coordinates (x, y, z) and the navigation satellite, that is, the total number of navigation satellites.

6. The satellite formation system fusion positioning optimization method based on genetic algorithm according to claim 5 is characterized in that: The specific process of obtaining the processed pseudorange observation value based on the pseudorange observation equation and the clock difference is as follows: In the formula, represents the pseudo-range observation value between the i′th navigation satellite and the jth user satellite, x i′ ,y i′ ,z i′ is the three-dimensional coordinate value of the i′th navigation star, x j ,y j ,z j is the three-dimensional coordinate value of the jth user star, t j is the clock difference between the user star and the navigation star, is the distance observation noise between the i′th navigation star and the jth user star.

7. The method for optimizing satellite formation system fusion positioning based on genetic algorithm according to claim 1, characterized in that: The cost f cost (X0,ΔX best_old ) is calculated as: Where X0 represents the user satellite position before adjustment, ΔX best_old represents the adjustment amount obtained after individual decoding, represents the baseline measurement value between the i-th and j-th user satellites obtained in step 1, d i,j (X0,ΔX best_old ) indicates the use of ΔX best_old The Euclidean distance between the i-th user star and the j-th user star obtained after adjusting X0; After updating the position using the individual adjustment amount, the difference between the baseline calculation value and the baseline measurement value after Kalman filtering between any two user satellites is calculated, and then the absolute values ​​of all the calculated differences are taken respectively, and the absolute values ​​are summed to obtain the cost function value.

8. The method for optimizing satellite formation system fusion positioning based on genetic algorithm according to claim 7, characterized in that: The individual optimal condition is: f cost (X0′,ΔX best_new )<f cost (X0,ΔX best_old ) j=flag1,flag2 i′=1,2,...,q Where ΔX best_new is the adjustment amount corresponding to the new optimal individual, f cost (X0′,ΔX best_new ) is the cost of the new individual obtained after adjusting X0′, X0′ is the position of the original optimal individual, f cost (X0,ΔX best_old ) is the cost of the original optimal individual, is the pseudorange observation value obtained in step 2, The adjustment amount ΔX best_new After updating the position X0′ of user star j, the Euclidean distance between the position of user star j and the position of navigation star i′ is adjusted; The adjustment amount ΔX best_old After updating the position X0 of user star j, the Euclidean distance between the position of user star j and the position of navigation star i′ is adjusted; Δρ is the pseudorange observation error value, and q is the number of navigation stars.

9. The method for optimizing satellite formation system fusion positioning based on genetic algorithm according to claim 8, characterized in that: The single optimization condition is: f cost (X k ,ΔX k+1 )<f cost (X k ,0)-0.01*(0.5) k′ L flagL Among them, X k is the user star position at the end of the kth iteration, ΔX k+1 is the adjustment amount obtained by optimal individual decoding at the k+1th iteration, f cost (X k ,ΔX k+1 ) is to use ΔX k+1 X k The cost of making adjustments, f cost (X k , 0) is to use 0 to X k The cost of the adjustment, k′ represents the number of iterations divided by 20 and rounded down.

10. The method for optimizing satellite formation system fusion positioning based on genetic algorithm according to claim 9, characterized in that: The crossover operation on the population is specifically as follows: Initialize a pseudo-random number for each individual in the population, compare the pseudo-random number with the crossover probability, if the pseudo-random number of an individual is less than the crossover probability, then put the individual into the crossover population, if the pseudo-random number of an individual is greater than or equal to the crossover probability, then put the individual into the non-crossover population; Perform crossover on the individuals in the crossover population. Assume that the number of individuals in the crossover population is m. c , if m c is an even number, then select the ath individual and the mth individual c -a individuals participate in crossover in pairs, if m c If is an odd number, the individuals in the middle of the population container do not participate in the crossover; When crossing, use random points to cross. First, initialize n p Pseudo-random numbers, the obtained n p After multiplying the pseudo-random numbers by the individual code length, the obtained n p The product is rounded down, and the rounded n p The values ​​are used as crossover points, and then the individuals participating in the crossover in pairs are crossed according to the crossover points to obtain the crossover population after the crossover. Merge the crossover population and the non-crossover population to obtain the crossover population; The mutation of the population after the crossover is specifically as follows: For any individual in the population after crossover, a pseudo-random number is randomly initialized for each binary code in the encoding of the individual. If the pseudo-random number of a binary code is less than the mutation rate, the binary code is inverted. If the pseudo-random number of a binary code is greater than or equal to the mutation rate, no mutation operation is performed. After processing each individual in the population after crossover, the mutated population is finally obtained. For any individual in the mutated population, a pseudo-random number is randomly initialized for the individual. If the pseudo-random number is less than 0.05, the individual is initialized as a local better individual near the best individual. Otherwise, if the pseudo-random number is greater than or equal to 0.05, the individual does not need to be processed. Similarly, each individual in the mutated population is processed separately. Initialize the individual as a local better individual near the optimal individual, which is specifically: Generate an adjustment amount with the first n_low bit value of 0, and use the adjustment amount to adjust the optimal individual to obtain a local better individual near the optimal individual.

Images