Angular momentum forecasting method and device for lunar orbit satellite and electronic equipment

By obtaining the satellite's angular momentum telemetry value and attitude quaternion, and combining the angular momentum increment model that fits the sun's direction angle with historical data, the accuracy of angular momentum prediction of satellites in lunar orbit is solved, and the accurate prediction of diagonal momentum increment is achieved, ensuring the orbital stability of the satellite.

CN119872928AActive Publication Date: 2025-04-25BEIJING AEROSPACE CONTROL CENT
View PDF 3 Cites 0 Cited by

Patent Information

Application Number
CN202510118583.6
Authority / Receiving Office
CN · China
Patent Type
Applications(China)
Current Assignee / Owner
Filing Date
2025-01-24
Publication Date
2025-04-25
Estimated Expiration
2045-01-24

AI Technical Summary

Technical Problem

In the angular momentum prediction of lunar orbit satellites, it is difficult for the prior art to accurately predict the changes in momentum wheel speed, resulting in orbit divergence effects and affecting the stable operation of the satellite.

Method used

By obtaining the satellite's angular momentum telemetry value and attitude quaternion, calculating the transformation matrix, determining the satellite's angular momentum vector under the inertia system and this system, using the angular momentum increment model that fits the angle included in the sun's direction with historical data, predicting the angular momentum increment in real time, and improving the accuracy of angular momentum forecasting.

Benefits of technology

Accurate prediction of the angular momentum increment of satellites in lunar orbit is achieved, the accuracy of angular momentum forecasting is improved, and the orbital stability of satellites and the reliability of execution of complex tasks is ensured.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN119872928A_ABST
    Figure CN119872928A_ABST
Patent Text Reader

Abstract

The invention discloses an angular momentum forecasting method and device for a lunar orbit satellite and electronic equipment. The method comprises the following steps: acquiring satellite historical angular momentum and attitude quaternion telemetering data, dividing the satellite historical angular momentum and attitude quaternion telemetering data into a plurality of periods by taking a rapid angular momentum change process as a node, and performing polynomial fitting on the angular momentum and attitude quaternion in each period to complement and smooth missing data caused by measurement and control constraints. And based on the complemented historical angular momentum and attitude quaternion, calculating the historical angular momentum increment of the satellite in the system within the specified step length. And calculating a historical included angle between the body axis with the maximum solar light pressure area of the satellite and the sun position vector in the system based on the satellite ephemeris. And fitting the historical included angle and the historical angular momentum increment corresponding to the historical included angle by using a sine and a model of a preset order to obtain a three-axis direction fitting coefficient of the system, and constructing an angular momentum increment fitting model. The model is adopted to predict the angular momentum forecast value of the satellite in the system, and accurate prediction of the angular momentum is achieved.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] The present application relates to the field of aerospace technology, and in particular to a method, device and electronic equipment for predicting angular momentum of a lunar orbit satellite. Background Art

[0002] The attitude control of satellites is mainly achieved by driving the momentum wheel to rotate through a motor, thereby exchanging angular momentum with the satellite. However, the rotation speed of the momentum wheel is not infinite. Once it reaches the critical value, it needs to be restored to the initial rotation speed through unloading. Although the speed increment generated by the unloading of the momentum wheel is relatively limited, in some delicate and sensitive three-body orbital environments in the Earth-Moon space, such a small change will also cause the orbit to produce a divergent effect that cannot be ignored. Therefore, it is particularly important to improve the accuracy of satellite angular momentum prediction. Summary of the invention

[0003] The embodiments of the present application provide a method, device and electronic equipment for predicting the angular momentum of a lunar orbit satellite, which are used to improve the accuracy of satellite angular momentum prediction.

[0004] In a first aspect, an embodiment of the present application provides a method for predicting angular momentum of a lunar orbit satellite, the method comprising:

[0005] Obtain the angular momentum telemetry values ​​of the satellite's three axes in this system and the satellite attitude quaternion at the start of the prediction. Based on the satellite attitude quaternion, calculate the transformation matrix of the satellite from the inertial system to this system. The satellite attitude quaternion is used to describe the satellite's attitude. The elements in the transformation matrix represent the projection components of the unit vector in the inertial system in this system.

[0006] Based on the angular momentum telemetry value and the reference angular momentum, the angular momentum vector of the satellite in this system is obtained. Based on the conversion matrix and the angular momentum vector in this system, the angular momentum vector of the satellite in the inertial system is determined through coordinate transformation. The reference angular momentum represents the indicated value of the angular momentum in the three axes when the angular momentum of the satellite is zero after the satellite is unloaded.

[0007] Input the solar direction angle into the angular momentum increment fitting model to obtain the angular momentum increment of the satellite in this system. The solar direction angle represents the angle between the axis with the largest solar light pressure area of ​​the satellite and the solar position vector in this system. The angular momentum increment fitting model is obtained by fitting the corresponding relationship between the historical solar direction angle and the historical angular momentum increment. The axis with the largest solar light pressure area represents the coordinate axis with the largest solar light pressure area of ​​the satellite in this system.

[0008] Based on the angular momentum increment, conversion matrix, angular momentum telemetry value and reference angular momentum under this system, the predicted angular momentum of the satellite's next step under this system is determined. If the conversion matrix and the reference angular momentum remain unchanged, the predicted angular momentum of the next step is positively correlated with the angular momentum increment under this system.

[0009] In the above method, by fitting the corresponding relationship between the historical solar direction angle and the historical angular momentum increment, the angular momentum increment of the satellite can be predicted more accurately, and the prediction accuracy of the angular momentum prediction value of the next step of the satellite under this system can be improved, so as to provide more accurate data support for the subsequent complex mission execution of the satellite. By obtaining the angular momentum telemetry value of the satellite under this system at the start of the prediction, the current angular momentum state of the satellite can be understood in real time. The reference angular momentum can ensure that when the angular momentum of the satellite is zero after the satellite is unloaded, the indicated values ​​of the three-axis directions of the angular momentum are accurate. By calculating the conversion matrix of the satellite from the inertial system to the present system, the unit vector in the inertial system can be projected into the present system to realize the coordinate system transformation. By fitting the corresponding relationship between the historical solar direction angle and the historical angular momentum increment, the future angular momentum increment of the satellite can be predicted more accurately. The angular momentum increment is calculated by the angular momentum increment fitting model and the real-time solar direction angle data, so that the angular momentum increment fitting model can be predicted in real time, and the accuracy of the angular momentum prediction value of the next step of the satellite under this system is improved.

[0010] Optionally, the above solar direction angle is obtained by:

[0011] By looking up the planetary ephemeris, calculate the solar position vector of the sun in the inertial frame at the start of the prediction;

[0012] Obtain satellite ephemeris data from the database to obtain the satellite position vector in the inertial system at the start time of prediction;

[0013] Based on the vector difference between the sun position vector and the satellite position vector, determine the relative position vector of the sun relative to the satellite in the inertial system;

[0014] Based on the transformation matrix and the relative position vector in the inertial system, the relative position vector of the sun relative to the satellite in this system is determined through coordinate transformation;

[0015] The inverse cosine function operation is performed on the ratio of the component of the relative position vector under this system on the axis with the largest solar light pressure area and the modulus of the relative position vector under this system to obtain the solar direction angle.

[0016] In the above method, high-precision solar position information can be obtained by searching the planetary ephemeris. The high precision of the satellite position vector is ensured by obtaining the satellite ephemeris data, thereby ensuring the accuracy of the relative position vector. By determining the precise solar direction angle through the relative position vector, the accuracy of the angular momentum prediction value of the next step of the satellite under this system can be improved.

[0017] Optionally, the angular momentum increment fitting model is fitted in the following way:

[0018] In each direction of the three axes of the satellite under this system, a sine sum model of preset different orders is used to fit each historical solar direction angle in the historical solar direction angle set and the historical angular momentum increment corresponding to the historical solar direction angle in the system, and the fitting coefficient of the satellite in each direction of the three axes is obtained. The historical solar direction angle set includes the historical angle between the axis with the largest solar light pressure area and the solar position vector of the satellite under this system at the time node corresponding to each step length in the historical time period;

[0019] The angular momentum increment fitting model is constructed based on the fitting coefficients.

[0020] In the above method, by using each historical solar direction angle in the historical time period and the historical angular momentum increment corresponding to the historical solar direction angle in this system for fitting, the model can better learn the law between the historical solar direction angle and the historical angular momentum increment, thereby enhancing the generalization ability of the model.

[0021] Optionally, the historical angular momentum increment corresponding to the above historical solar direction angle in this system is obtained by the following method:

[0022] Obtaining a historical telemetry data set and performing standardization on the historical telemetry data set, the historical telemetry data set includes multiple satellite times of the satellite in the system within a historical period, historical angular momentum telemetry values ​​corresponding to the satellite times, and historical satellite attitude quaternions corresponding to the satellite times;

[0023] For each moment in the standardized historical telemetry data set, the vector difference between the historical angular momentum telemetry value corresponding to each moment and the reference angular momentum is calculated to obtain the instantaneous angular momentum of the satellite in this system at each moment;

[0024] Based on the transformation matrix and instantaneous angular momentum of the satellite from the inertial system to the local system at each moment, the instantaneous angular momentum of the satellite in the inertial system is determined by coordinate transformation. The transformation matrix of the satellite from the inertial system to the local system at each moment is obtained by calculating the historical satellite attitude quaternion at each moment.

[0025] Calculate the vector difference between the instantaneous angular momentum of the next moment in the inertial system and the instantaneous angular momentum of each moment in the inertial system, and obtain the historical angular momentum increment of the satellite in the inertial system;

[0026] Based on the transformation matrix of the satellite from the inertial system to the current system at each moment and the historical angular momentum increment of the satellite in the inertial system, the historical angular momentum increment of the satellite in the current system is determined by coordinate transformation.

[0027] In the above method, by standardizing the historical telemetry data set, the data can be made more unified, which helps to improve the efficiency and accuracy of subsequent data processing. By calculating the vector difference between the historical angular momentum telemetry value and the reference angular momentum corresponding to each moment, the instantaneous angular momentum of the satellite in this system at each moment can be obtained, which helps to more accurately understand the angular momentum state of the satellite. Using the conversion matrix, flexible coordinate transformations can be performed between the inertial system and this system, and the angular momentum state of the satellite can be analyzed from multiple angles, improving the accuracy of the subsequent prediction of the angular momentum forecast value of the next step of the satellite in this system.

[0028] Optionally, the above historical telemetry data set includes a historical angular momentum telemetry data set and a historical attitude quaternion telemetry data set, and the historical telemetry data set is standardized, specifically including:

[0029] Delete the relevant data of the rapid change process of angular momentum in the historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set. The rapid change process of angular momentum at least includes the satellite's rapid attitude adjustment, momentum wheel unloading, and orbit control;

[0030] Taking the rapid change process of angular momentum as a node, the historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set are divided into multiple angular momentum periods, each angular momentum period includes a subset of historical angular momentum telemetry data and a subset of historical attitude telemetry data in different time periods;

[0031] The historical angular momentum telemetry data subset and the historical attitude telemetry data subset within each angular momentum cycle are fitted to complete and smooth the missing data points of the historical angular momentum telemetry data subset and the historical attitude telemetry data subset.

[0032] In the above method, by deleting the relevant data of the rapid change process of angular momentum, such as the data of the period of rapid satellite attitude adjustment, momentum wheel unloading, orbit control, etc., the interference of special states on data analysis can be avoided, thereby improving the stability of the data. By dividing the processed data into multiple angular momentum periods, it is helpful to identify and analyze the periodic change law of satellite angular momentum and improve the accuracy of subsequent data analysis. The missing data points can be supplemented and smoothed through fitting processing, making the historical telemetry data set more complete and continuous, and reducing the analysis errors caused by missing or abnormal data.

[0033] Optionally, the above fitting processing is performed on the historical angular momentum telemetry data subset and the historical attitude telemetry data subset in each angular momentum period to complete and smooth the missing data points of the historical angular momentum telemetry data subset and the historical attitude telemetry data subset, specifically including:

[0034] For each angular momentum period, a preset first-order polynomial is used to fit the components of each historical angular momentum telemetry value in each direction of the three axes of the satellite in this system to obtain a first fitting coefficient, and a first fitting formula is determined based on the first fitting coefficient;

[0035] Using a preset second-order polynomial to fit each component of each historical satellite attitude quaternion to obtain a second fitting coefficient, and determining a second fitting formula based on the second fitting coefficient;

[0036] Using a first fitting formula, based on known data points of the historical angular momentum telemetry data subset and a preset step size, the missing data points of the historical angular momentum telemetry data subset are determined, and the completed historical angular momentum telemetry data subset is smoothed;

[0037] The second fitting formula is used to determine the missing data points of the historical attitude quaternion telemetry data subset based on the known data points of the historical attitude quaternion telemetry data subset and the preset step size, and the completed historical attitude quaternion telemetry data subset is smoothed.

[0038] In the above method, the missing values ​​in the historical angular momentum telemetry data subset and the historical attitude telemetry data subset can be supplemented by polynomial fitting, making the historical angular momentum telemetry data subset and the historical attitude telemetry data subset more complete. At the same time, since the fitting process itself has the characteristics of smoothing data, it can reduce the noise and outliers in the data and improve the continuity and smoothness of the data.

[0039] Optionally, the above-mentioned determination of the predicted angular momentum value of the next step of the satellite under the present system based on the angular momentum increment, the conversion matrix, the angular momentum telemetry value and the reference angular momentum under the present system specifically includes:

[0040] The transpose of the conversion matrix is ​​multiplied by the transpose of the angular momentum increment in this system to obtain the angular momentum increment of the satellite in the inertial system;

[0041] Based on the angular momentum increment in the inertial system and the vector sum of the angular momentum vector in the inertial system, the angular momentum vector of the satellite's next step in the inertial system is determined;

[0042] The transformation matrix is ​​multiplied by the angular momentum vector of the next step length in the inertial system, and the angular momentum vector of the next step length of the satellite in the current system is determined by coordinate transformation;

[0043] Based on the vector sum of the angular momentum vector of the next step length in this system and the reference angular momentum, the predicted value of the angular momentum of the next step length of the satellite in this system is determined.

[0044] In the above method, the angular momentum increment in this system can be accurately converted to the inertial system by multiplying the transpose of the conversion matrix and the transpose of the angular momentum increment, ensuring the consistency of the angular momentum vector between different coordinate systems and improving the accuracy of subsequent data analysis. Since the angular momentum vector in the inertial system represents the absolute motion state of the satellite in space, and the angular momentum vector in this system represents the motion state of the satellite relative to its own coordinate system. By understanding the predicted angular momentum of the satellite at different times and in different coordinate systems, it is convenient to control the attitude and orbit of the satellite more accurately in the future and ensure the stable operation of the satellite.

[0045] In a second aspect, an embodiment of the present application provides an angular momentum prediction device for a lunar orbit satellite, the device comprising:

[0046] The transceiver module is used to obtain the angular momentum telemetry value of the satellite's three axes in this system at the start time of the prediction and the satellite attitude quaternion, which is used to describe the satellite attitude;

[0047] A processing module is used to calculate the transformation matrix of the satellite from the inertial system to the current system based on the satellite attitude quaternion, where the elements in the transformation matrix represent the projection components of the unit vector in the inertial system in the current system;

[0048] The processing module is further used to obtain the angular momentum vector of the satellite in the present system based on the angular momentum telemetry value and the reference angular momentum, and determine the angular momentum vector of the satellite in the inertial system through coordinate transformation based on the conversion matrix and the angular momentum vector in the present system, wherein the reference angular momentum represents the indicated values ​​of the angular momentum in the three axes when the angular momentum of the satellite is zero after the satellite is unloaded;

[0049] The processing module is also used to input the solar direction angle into the angular momentum increment fitting model to obtain the angular momentum increment of the satellite in this system. The solar direction angle represents the angle between the axis with the largest solar light pressure area of ​​the satellite and the solar position vector in this system. The angular momentum increment fitting model is obtained by fitting the corresponding relationship between the historical solar direction angle and the historical angular momentum increment. The axis with the largest solar light pressure area represents the coordinate axis with the largest solar light pressure area of ​​the satellite in this system.

[0050] The processing module is also used to determine the predicted angular momentum of the next step of the satellite under the present system based on the angular momentum increment, conversion matrix, angular momentum telemetry value and reference angular momentum under the present system. If the conversion matrix and the reference angular momentum remain unchanged, the predicted angular momentum of the next step is positively correlated with the angular momentum increment under the present system.

[0051] In a third aspect, an embodiment of the present application provides an electronic device, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein when the computer program is executed by the processor, the processor implements any one of the methods in the first aspect.

[0052] In a fourth aspect, an embodiment of the present application further provides a computer-readable storage medium, in which a computer program is stored. When the computer program is executed by a processor, any one of the methods in the first aspect is implemented.

[0053] In a fifth aspect, an embodiment of the present application further provides a computer program product, including a computer program, which is executed by a processor to implement any method in the first aspect described above.

[0054] The technical effects brought about by any one of the implementation methods in the second to fifth aspects can refer to the technical effects brought about by the corresponding implementation method in the first aspect, and will not be repeated here. BRIEF DESCRIPTION OF THE DRAWINGS

[0055] Figure 1 This is a flow chart of a method for predicting angular momentum of a lunar orbit satellite provided in an embodiment of the present application;

[0056] Figure 2 A schematic diagram of the relationship between the angular momentum increment of a satellite on the X-axis of the system and the angle to the sun provided in an embodiment of the present application;

[0057] Figure 3 A schematic diagram of the relationship between the angular momentum increment of a satellite on the Y-axis of the system and the angle to the sun provided in an embodiment of the present application;

[0058] Figure 4 A schematic diagram of the relationship between the angular momentum increment of a satellite on the Z axis of the system and the angle to the sun provided in an embodiment of the present application;

[0059] Figure 5 A schematic diagram of the structure of an angular momentum prediction device for a lunar orbit satellite provided in an embodiment of the present application;

[0060] Figure 6 A schematic diagram of the structure of a control device provided in an embodiment of the present application. DETAILED DESCRIPTION

[0061] Below, some terms in the embodiments of the present application are explained to facilitate understanding by those skilled in the art.

[0062] 1. "And / or" describes the relationship between related objects, indicating that there can be three relationships. For example, A and / or B can mean: A exists alone, A and B exist at the same time, and B exists alone. The character " / " generally indicates that the related objects are in an "or" relationship.

[0063] 2. The terms "first", "second", etc. in the specification and claims of this application and the above-mentioned drawings are used to distinguish similar objects, and are not necessarily used to describe a specific order or sequence. It should be understood that the terms used in this way can be interchanged where appropriate, so that the embodiments described herein can be implemented in an order other than that illustrated or described herein.

[0064] 3. Inertial frame, full name is inertial reference frame (Inertial Frame of Reference), is a crucial concept in physics. It refers to a reference frame in which an object will remain stationary or in uniform linear motion when not acted upon by external forces. Inertial frame is a reference frame in which Newton's laws of motion are valid, also known as inertial coordinate system, or inertial frame for short. In an inertial frame, the passage of time is uniform, and space is uniform and isotropic.

[0065] 4. This system, also called the body coordinate system or satellite body coordinate system, is a reference system established relative to the satellite itself. In this reference system, the origin of the coordinates is usually selected at the center of mass of the satellite or a specific point, and the coordinate axes are determined according to the structure, function and design requirements of the satellite. This system is mainly used to describe the attitude, orbit and motion state of the satellite itself.

[0066] 5. Satellite ephemeris data, also known as two-line orbital data (TLE, Two-Line Orbital Element), is an expression invented by Celestrak in the United States to describe the position and speed of space flight objects (such as satellites, spacecraft, etc.), that is, the two-line orbital data system. It determines the time, coordinates, direction, speed and other parameters of the flight object based on the mathematical relationship between the six orbital parameters of Kepler's law.

[0067] 6. Satellite time usually refers to the time record on the satellite, which may be synchronized or calibrated based on a specific time standard (such as UTC time, GPS time, etc.). In the satellite system, satellite time is used to record key information such as the satellite's operating status and data collection time.

[0068] In order to make the purpose, technical solutions and advantages of the present application clearer, the present application will be further described in detail below in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present application.

[0069] The application scenarios described in the embodiments of the present application are intended to more clearly illustrate the technical solutions of the embodiments of the present application, and do not constitute a limitation on the technical solutions provided by the embodiments of the present application. It is known to those skilled in the art that with the emergence of new application scenarios, the technical solutions provided by the embodiments of the present application are also applicable to similar technical problems. In the description of the present application, unless otherwise specified, the meaning of "multiple" is two or more.

[0070] It can be understood that the present application is applicable to the case where the solar light pressure area of ​​a certain axis of the satellite in this system is much larger than the solar light pressure area of ​​the other two axis directions. Among them, the axis with the largest solar light pressure area represents the coordinate axis with the largest solar light pressure area in the satellite system, which can be the X-axis, Y-axis or Z-axis of the system. For example, the present application is applicable to the case where the solar light pressure area of ​​the lunar orbit satellite in the Z-axis direction is much larger than the solar light pressure area in the X-axis and Y-axis directions in this system. In this case, the axis with the largest solar light pressure area is the Z-axis.

[0071] The letters A, B, C, etc. in the formulas of this application, whether in regular or italic form, have the same meaning.

[0072] It should be noted that, in the satellite's body coordinate system, the axis with the largest satellite solar light pressure area can be any one of the X-axis, Y-axis or Z-axis in this system. The following takes the Z-axis as an example to illustrate the embodiment of the present application.

[0073] like Figure 1 As shown, a flow chart of a method for predicting angular momentum of a lunar orbit satellite provided in an embodiment of the present application may specifically include the following steps.

[0074] Step S101, obtaining the angular momentum telemetry values ​​of the three axes of the satellite in the current system and the satellite attitude quaternion at the prediction start time, and calculating the transformation matrix of the satellite from the inertial system to the current system based on the satellite attitude quaternion.

[0075] The satellite attitude quaternion is used to describe the satellite attitude. The elements in the transformation matrix represent the projection components of the unit vector in the inertial system in this system.

[0076] For example, assuming that the satellite attitude quaternion Q at the start of the prediction is [Q c ,Q1,Q2,Q3]. Among them, Q c Q1, Q2, and Q3 are the components of the satellite attitude quaternion. Based on the satellite attitude quaternion, the transformation matrix M from the inertial system to the local system is calculated. bi The values ​​of each element in satisfy the following formula:

[0077]

[0078] Among them, M bi (0,0) represents the element in the first row and first column of the transformation matrix. bi (0,1) represents the element in the first row and second column of the transformation matrix. bi (0,2) represents the element in the first row and third column of the transformation matrix. bi (1,0) represents the element in the second row and first column of the transformation matrix. bi (1,1) represents the element in the second row and second column of the transformation matrix. bi (1,2) represents the element in the second row and third column of the transformation matrix. bi (2,0) represents the element in the third row and first column of the transformation matrix. bi (2,1) represents the element in the third row and second column of the transformation matrix. bi (2,2) represents the element in the third row and third column of the transformation matrix.

[0079] Step S102: obtaining the angular momentum vector of the satellite in the present system based on the angular momentum telemetry value and the reference angular momentum, and determining the angular momentum vector of the satellite in the inertial system through coordinate transformation based on the conversion matrix and the angular momentum vector in the present system.

[0080] Among them, the reference angular momentum refers to the indicated values ​​of the angular momentum in the three-axis directions when the satellite angular momentum is zero after the satellite is unloaded.

[0081] Since satellites usually use a redundant design of four momentum wheels, the X, Y, and Z axis angular momentum telemetry is not zero when the satellite's overall angular momentum returns to zero. Suppose the X, Y, and Z axis angular momentum indication value when the satellite's angular momentum returns to zero, that is, the reference angular momentum is H inib =[h inibx h iniby h inibz ], the angular momentum telemetry value at the beginning of the prediction is H yb0 =[h ybx0 h yby0 h ybz0 ]. Then the actual angular momentum vector H of the satellite in this system is b0 Satisfies the following formula:

[0082] H b0 =H yb0 -H inib

[0083] For example, the actual angular momentum vector is (0, 0, 0). The angular momentum telemetry value displayed by the satellite should be (0, 0, 0). However, since satellites usually use a redundant design of four momentum wheels, the overall angular momentum of the satellite is not displayed as zero when it returns to zero, such as (1, 1, 1). Then the reference angular momentum is (1, 1, 1). For another example, assuming the reference angular momentum is (2, 2, 2), when the actual angular momentum vector of the satellite is (1, 1, 1), the angular momentum telemetry value displayed by the satellite is (3, 3, 3).

[0084] It is understandable that if the satellite does not adopt a redundant design, the angular momentum telemetry value displayed by the satellite is the same as the actual angular momentum vector of the satellite. At this time, the reference angular momentum is (0, 0, 0).

[0085] In an optional embodiment, after obtaining the angular momentum vector of the satellite in this system, the transpose of the transformation matrix from the inertial system to the present system can be multiplied by the angular momentum vector of the satellite in this system, and the angular momentum vector of the satellite in the inertial system can be determined by coordinate transformation.

[0086] For example, suppose the transformation matrix from the inertial system to the present system is M bi The angular momentum vector of the satellite in this system is H b0 , then the angular momentum vector H of the satellite in the inertial system is i0 Satisfies the following formula:

[0087]

[0088] Step S103: input the sun direction angle into the angular momentum increment fitting model to obtain the angular momentum increment of the satellite in this system.

[0089] Among them, the solar direction angle represents the angle between the axis with the largest solar light pressure area of ​​the satellite and the solar position vector in this system. The angular momentum increment fitting model is obtained by fitting the corresponding relationship between the historical solar direction angle and the historical angular momentum increment. The axis with the largest solar light pressure area represents the coordinate axis with the largest solar light pressure area of ​​the satellite in this system.

[0090] The following is an explanation of how to obtain the sun direction angle:

[0091] In an optional embodiment, the sun position vector of the sun in the inertial system at the start time of the prediction can be calculated by looking up the planetary ephemeris. The satellite ephemeris data is obtained from the database to obtain the satellite position vector of the satellite in the inertial system at the start time of the prediction. Based on the vector difference between the sun position vector and the satellite position vector, the relative position vector of the sun relative to the satellite in the inertial system is determined. Based on the conversion matrix and the relative position vector in the inertial system, the relative position vector of the sun relative to the satellite in this system is determined by coordinate transformation. The ratio of the component of the relative position vector in this system on the axis with the largest solar light pressure area to the modulus of the relative position vector in this system is subjected to an inverse cosine function operation to obtain the solar direction angle.

[0092] For example, by looking up the JPL planetary ephemeris, we can calculate the solar position vector r of the sun in the inertial frame at the start of the prediction. sun =[r sunx r suny r sunz ] Search the satellite's ephemeris data from the database to obtain the satellite position vector r in the inertial system at the start of the prediction sat =[r satx r saty r satz ].

[0093] Based on the vector difference between the sun position vector and the satellite position vector, determine the relative position vector Δr of the sun relative to the satellite in the inertial system i , satisfying the following formula:

[0094] Δr i =r sun -r sat

[0095] The transformation matrix M that transforms the satellite from the inertial system to this system bi Multiply it with the relative position vector in the inertial system and determine the relative position vector Δr of the sun relative to the satellite in this system through coordinate transformation b , satisfying the following formula:

[0096]

[0097] The inverse cosine function operation is performed on the ratio of the component of the relative position vector on the axis with the largest solar light pressure area under this system to the modulus of the relative position vector under this system to obtain the solar direction angle Satisfies the following formula:

[0098]

[0099] The following is an explanation of how to obtain the angular momentum increment fitting model:

[0100] In an optional embodiment, in each direction of the three axes of the satellite under this system, a sine sum model of preset different orders can be used to fit each historical solar direction angle in the historical solar direction angle set and the historical angular momentum increment corresponding to the historical solar direction angle in this system, and the fitting coefficient of the satellite in each direction of the three axes can be obtained. The historical solar direction angle set includes the historical angle between the axis with the largest solar light pressure area of ​​the satellite and the solar position vector under this system at the time node corresponding to each step in the historical time period. The angular momentum increment fitting model is constructed based on the fitting coefficient.

[0101] Among them, the historical angular momentum increment corresponding to the historical solar direction angle in this system is obtained by the following method:

[0102] Obtain a historical telemetry data set and perform standardization on the historical telemetry data set. The historical telemetry data set includes multiple satellite star times of the satellite in the current system during the historical period, historical angular momentum telemetry values ​​corresponding to the satellite star times, and historical satellite attitude quaternions corresponding to the satellite star times. For each moment in the standardized historical telemetry data set, calculate the vector difference between the historical angular momentum telemetry value corresponding to each moment and the reference angular momentum to obtain the instantaneous angular momentum of the satellite in the current system at each moment. Based on the transformation matrix and instantaneous angular momentum of the satellite from the inertial system to the current system at each moment, determine the instantaneous angular momentum of the satellite in the inertial system through coordinate transformation. The transformation matrix of the satellite from the inertial system to the current system at each moment is obtained by calculating the historical satellite attitude quaternions at each moment. Calculate the vector difference between the instantaneous angular momentum of the next moment in the inertial system and the instantaneous angular momentum of the satellite in the inertial system at each moment to obtain the historical angular momentum increment of the satellite in the inertial system. Based on the transformation matrix of the satellite from the inertial system to the current system at each moment and the historical angular momentum increment of the satellite in the inertial system, the historical angular momentum increment of the satellite in the current system is determined by coordinate transformation. The time difference between any two adjacent moments in the time series is a preset step size.

[0103] For example, the historical telemetry data set may include a historical angular momentum telemetry data set and a historical attitude quaternion telemetry data set. After obtaining the historical telemetry data set, the historical telemetry data set may be standardized in the following manner:

[0104] Delete the relevant data of the rapid change process of angular momentum in the historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set. The rapid change process of angular momentum at least includes rapid attitude adjustment of the satellite, momentum wheel unloading, and orbit control. Taking the rapid change process of angular momentum as a node, the processed historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set are divided into multiple angular momentum cycles. Each angular momentum cycle includes a subset of historical angular momentum telemetry data and a subset of historical attitude telemetry data in different time periods. Fit the subset of historical angular momentum telemetry data and the subset of historical attitude telemetry data in each angular momentum cycle to complete and smooth the missing data points of the subset of historical angular momentum telemetry data and the subset of historical attitude telemetry data.

[0105] For example, a historical telemetry data set after the satellite enters the lunar mission orbit is obtained from the database. The historical telemetry data set may include a historical angular momentum telemetry data set S H0 And the historical attitude quaternion telemetry data set S Q0 .

[0106] Among them, S H0 Contains satellite time (i.e., the time stamp of satellite operation) and the historical angular momentum telemetry value corresponding to the satellite time in this system. Q0 Contains satellite time and historical satellite attitude quaternions.

[0107] Historical angular momentum telemetry data set S H0 And the historical attitude quaternion telemetry data set S Q0 Satisfy the following formulas respectively:

[0108]

[0109] Where i is an integer. i is the satellite time. bxi ,H byi ,H bzi They are the components of the historical angular momentum telemetry values ​​in the three directions of the satellite in the system, namely, X, Y, and Z. ci ,Q 1i ,Q 2i ,Q 3i is the component of historical attitude quaternion telemetry data.

[0110] In the historical angular momentum telemetry data set S H0 And the historical attitude quaternion telemetry data set S Q0 In the process of identifying the rapid change of angular momentum, the relevant data of the process of identifying the rapid change of angular momentum can be included in the process of identifying the rapid change of angular momentum. Ai (i=1,...,k), orbit control process T Ci(i=1,...,m), momentum wheel unloading T Wi (i=1,...,n). The relevant data of the angular momentum rapid change process are collected from the original angular momentum telemetry data set S H0 And the original attitude telemetry data collection S Q0 , and the processed historical angular momentum telemetry data set S is obtained. H and the historical attitude quaternion telemetry data set S Q .

[0111] The processed historical angular momentum telemetry data set S H and the historical attitude quaternion telemetry data set S Q Sort by time order to get the time node sequence T i (i=1,...,t,t=k+m+n). According to the time periods divided by different processes of rapid change of angular momentum, the processed historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set are divided into multiple angular momentum periods.

[0112] Each angular momentum period satisfies the following formula:

[0113]

[0114] Where t is an integer. In the time period from T1 to T2, the first angular momentum cycle includes the historical angular momentum telemetry data subset S Q1 And the historical attitude telemetry data subset S H1 During the time period from T2 to T3, the second angular momentum cycle includes the historical angular momentum telemetry data subset S Q2 And the historical attitude telemetry data subset S H2 The tth angular momentum cycle includes the historical angular momentum telemetry data subset S Qt And the second historical attitude telemetry data subset S Ht .

[0115] For example, suppose that the time series in the historical telemetry data set is from 0:00 to 24:00 on a certain day. In this event sequence, the rapid change process of angular momentum includes rapid satellite attitude adjustment, momentum wheel unloading, and orbit control. Assume that the satellite rapid attitude adjustment is performed from 5:00 to 5:02, the momentum wheel unloading is performed from 12:00 to 12:02, and the orbit control is performed from 18:00 to 18:05. Then the time series from 0:00 to 24:00 can be divided into 4 angular momentum cycles. That is, 0:00 to 5:00 is the time series corresponding to the first angular momentum cycle. 5:02 to 12:00 is the time series corresponding to the second angular momentum cycle. 12:02 to 18:00 is the time series corresponding to the third angular momentum cycle. 18:05 to 24:00 is the time series corresponding to the fourth angular momentum cycle.

[0116] Due to the influence of tracking arc segment limitation and telemetry minimum division value, the historical angular momentum telemetry data subset and historical attitude telemetry data subset within the angular momentum cycle will have discontinuity and data rebound problems. In order to obtain more accurate, continuous and smooth data, it is necessary to fit the historical angular momentum telemetry data subset and historical attitude telemetry data subset within each angular momentum cycle, and complete and smooth the missing data points of the historical angular momentum telemetry data subset and historical attitude telemetry data subset.

[0117] For example, the missing data points of the historical angular momentum telemetry data subset and the historical attitude telemetry data subset can be completed in the following manner:

[0118] For each angular momentum cycle, a preset first-order polynomial is used to fit the components of each historical angular momentum telemetry value in each direction of the three axes of the satellite under this system to obtain a first fitting coefficient, and a first fitting formula is determined based on the first fitting coefficient. A preset second-order polynomial is used to fit each component of each historical satellite attitude quaternion to obtain a second fitting coefficient, and a second fitting formula is determined based on the second fitting coefficient. The first fitting formula is used to determine the missing data points of the historical angular momentum telemetry data subset based on the known data points and preset step size of the historical angular momentum telemetry data subset, and the completed historical angular momentum telemetry data subset is smoothed. The second fitting formula is used to determine the missing data points of the historical attitude quaternion telemetry data subset based on the known data points and preset step size of the historical attitude quaternion telemetry data subset, and the completed historical attitude quaternion telemetry data subset is smoothed.

[0119] For example, a fifth-order polynomial is used to calculate the historical angular momentum telemetry value [H bx H by H bz ] Under this system, the components of the satellite's three-axis directions (x, y, z) are fitted, and the interpolation coefficients of the x, y, and z directions are [a1 b1 c1 d1 e1 f1], [a2 b2 c2 d2e2 f2] and [a3 b3 c3 d3 e3 f3].

[0120] Among them, the fifth-order polynomial satisfies the following formula:

[0121]

[0122] The historical satellite attitude quaternion [Q c Q1 Q2 Q3] each component is fitted, fitting Q c, the interpolation coefficients of the four elements Q1, Q2, Q3 are [a0 b0 c0 d0], [a1 b1 c1 d1], [a2 b2 c2 d2] and [a3 b3 c3d3].

[0123] Among them, the 4th order polynomial satisfies the following formula:

[0124]

[0125] After the interpolation coefficient fitting is completed, the first fitting formula can be used to determine the missing data points of the historical angular momentum telemetry data subset based on the known data points of the historical angular momentum telemetry data subset and the preset step size, and the completed historical angular momentum telemetry data subset is smoothed. The second fitting formula is used to determine the missing data points of the historical attitude quaternion telemetry data subset based on the known data points of the historical attitude quaternion telemetry data subset and the preset step size, and the completed historical attitude quaternion telemetry data subset is smoothed.

[0126] It is understandable that if the ground has attitude quaternion prediction data, if there is prediction data, the data can be directly used for subsequent calculations. If there is no prediction data, the telemetry attitude quaternion can be interpolated and fitted.

[0127] In an optional embodiment, after the historical telemetry data set is standardized, for each moment in the standardized historical telemetry data set, the vector difference between the historical angular momentum telemetry value corresponding to each moment and the reference angular momentum is calculated to obtain the instantaneous angular momentum of the satellite in the present system at each moment. The instantaneous angular momentum of the satellite in the inertial system is determined by multiplying the instantaneous angular momentum by the transpose of the transformation matrix of the satellite from the inertial system to the present system at each moment, and by coordinate transformation.

[0128] For example, suppose the historical angular momentum telemetry value corresponding to a certain historical moment is H ybq =[h ybxq h ybyq h ybzq ], the reference angular momentum is H inib =[h inibx h iniby h inibz ]. Where q is an integer. Then the instantaneous angular momentum H of the satellite in this system at this historical moment is bq , satisfying the following formula:

[0129] H bq =H ybq -H inib

[0130] The transpose of the transformation matrix that transforms the satellite from the inertial system to this system at that historical moment and the instantaneous angular momentum H bq Multiply them and determine the instantaneous angular momentum H of the satellite in the inertial system through coordinate transformation iq , satisfying the following formula:

[0131]

[0132] After obtaining the instantaneous angular momentum at each moment in the inertial system, the vector difference between the instantaneous angular momentum at the next moment in the inertial system and the instantaneous angular momentum at each moment in the inertial system can be calculated to obtain the historical angular momentum increment of the satellite in the inertial system.

[0133] For example, suppose that the instantaneous angular momentum in the inertial system at a certain moment in history is H i1 The instantaneous angular momentum in the inertial system at the next moment is H i2 Calculate H i2 With H i1 The vector difference is used to obtain the historical angular momentum increment ΔH of the satellite in the inertial system at that moment. i , satisfying the following formula:

[0134] ΔH i =H i2 -H i1

[0135] Based on the transformation matrix M of the satellite from the inertial system to the current system at each moment bi The historical angular momentum increment ΔH of the satellite in the inertial system i Multiply them and determine the historical angular momentum increment ΔH of the satellite in this system through coordinate transformation b , satisfying the following formula:

[0136] ΔH b =M bi ΔH i =[ΔH bx ΔH by ΔH bz ]

[0137] In an optional embodiment, after obtaining each historical solar direction angle and the historical angular momentum increment corresponding to the historical solar direction angle in this system, each historical solar direction angle in the historical solar direction angle set and the historical angular momentum increment corresponding to the historical solar direction angle in this system can be fitted using preset sine sum models of different orders to obtain the fitting coefficients of the satellite in each direction of the three axes, and a fitting model for the angular momentum increment is constructed based on the fitting coefficients.

[0138] For example, suppose the angular momentum increment of the satellite in this system is ΔH b =[ΔH bx ΔH byΔH bz ]. For the component ΔH of the angular momentum increment on the X axis bx , we can use the fifth-order sine sum model to calculate ΔH bx Perform fitting and obtain a xi ,b xi ,c xi (i=1,...,5) coefficient values. Where i is a positive integer.

[0139] The fifth-order sine sum model satisfies the following formula:

[0140]

[0141] For the angular momentum increment on the Y axis, ΔH by , the seventh-order sine sum model can be used to calculate ΔH by Perform fitting and obtain a yi ,b yi ,c yi (i=1,...,7) coefficient values.

[0142] The seventh-order sine sum model satisfies the following formula:

[0143]

[0144] For the component of angular momentum increment in the Z axis, ΔH bz , the third-order sine sum model can be used to calculate ΔH bz Perform fitting and obtain a zi ,b zi ,c zi (i=1,...,3) coefficient values.

[0145] The third-order sine sum model satisfies the following formula:

[0146]

[0147] After obtaining the fitting coefficients in each direction of the three axes, the above-mentioned angular momentum increment fitting model is constructed.

[0148] In an optional embodiment, it is assumed that the preset step length is Δt. The angle of the sun direction at the current moment is Enter the sun angle To the angular momentum increment fitting model, the angular momentum increment ΔH of the satellite in this system is obtained b0 =[ΔH bx0 ΔH by0 ΔH bz0 ], the components of the angular momentum increment satisfy the following formula:

[0149]

[0150] Step S104: determining the predicted angular momentum value of the next step of the satellite in the present system based on the angular momentum increment, the conversion matrix, the angular momentum telemetry value and the reference angular momentum in the present system.

[0151] Among them, if the conversion matrix and the reference angular momentum remain unchanged, the predicted value of the angular momentum of the next step length is positively correlated with the angular momentum increment under this system.

[0152] In an optional embodiment, the transpose of the conversion matrix and the transpose of the angular momentum increment in this system can be multiplied to obtain the angular momentum increment of the satellite in the inertial system. Based on the angular momentum increment in the inertial system and the vector sum of the angular momentum vector in the inertial system, the angular momentum vector of the next step of the satellite in the inertial system is determined. The conversion matrix is ​​multiplied by the angular momentum vector of the next step of the satellite in the inertial system, and the angular momentum vector of the next step of the satellite in this system is determined by coordinate transformation. Based on the vector sum of the angular momentum vector of the next step of the satellite in this system and the reference angular momentum, the predicted angular momentum value of the next step of the satellite in this system is determined.

[0153] For example, suppose the angular momentum increment of the satellite in this system is ΔH b0 =[ΔH bx0 ΔH by0 ΔH bz0 ]. The transformation matrix is ​​M bi Then the angular momentum increment of the satellite in the inertial system is ΔH i0 . ΔH i0 Satisfies the following formula:

[0154]

[0155] Assume that the angular momentum vector in the inertial system is H i0 . Based on the angular momentum increment ΔH in the inertial system i0 And the angular momentum vector H in the inertial frame i0 The vector sum of the satellite's next step in the inertial system determines the angular momentum vector H i1 .H i1 Satisfies the following formula;

[0156] H i1 =H i0 +ΔH i0

[0157] The transformation matrix M bi The angular momentum vector H of the next step length in the inertial frame i1 Perform product operation and determine the angular momentum vector H of the satellite's next step in this system through coordinate transformation b1 .H b1 Satisfies the following formula:

[0158] H b1 =M bi ×H i1

[0159] Assume the reference angular momentum is H inib =[h inibx h iniby h inibz ]. Based on the next step length angular momentum vector H in this system b1 With the reference angular momentum H inib =[h inibx h iniby h inibz ] to determine the predicted angular momentum H of the satellite's next step in this system. yb1 .H yb1 Satisfies the following formula:

[0160] H yb1 =H b1 +H inib

[0161] After calculating the predicted angular momentum of the satellite in the future, it can be used to predict the saturation time of the satellite's momentum wheel, plan its unloading time, and avoid unexpected unloading of the momentum wheel. Because the orbital stability of some Earth-Moon systems is poor, unexpected unloading of the momentum wheel may cause a large deviation in the orbit, resulting in the failure of some on-orbit tests with high orbital accuracy requirements.

[0162] In another possible scenario, the angular momentum prediction method can be used for satellite angular momentum management. For example, a Halo orbit satellite at the Earth-Moon L2 point can maintain its orbital configuration by unloading the momentum wheel. Based on the obtained angular momentum increment fitting model, the satellite's demand for angular momentum change is considered, the satellite's on-orbit attitude is planned, and the satellite's angular momentum is controlled to increase. Unloading the orbital configuration at the appropriate time can greatly reduce the orbit maintenance frequency of the Halo orbit satellite.

[0163] In another possible scenario, satellite fault troubleshooting can also be performed based on the angular momentum prediction value. For example, by comparing the angular momentum prediction value with the actual angular momentum telemetry value, it can be determined whether there is a potential fault. If there is a large deviation between the satellite's angular momentum prediction value and the actual angular momentum telemetry value, the satellite's attitude control system or related sensors are faulty.

[0164] It is understandable that when the axis with the largest satellite solar pressure area is not the Z axis in this system, the coordinate axis related items in the formula also need to be replaced accordingly.

[0165] If the X-axis is the axis with the largest solar radiation pressure area, the variables involving the Z-axis in the original formula should be replaced by X, and the Y-axis remains unchanged.

[0166] If the Y axis is the axis with the largest solar radiation pressure area, the variables involving the Z axis in the original formula should be replaced by Y, and the X axis remains unchanged.

[0167] The angular momentum prediction method when the axis with the largest satellite solar light pressure area under this system is the X-axis or the Y-axis is consistent with the angular momentum prediction method when the axis with the largest satellite solar light pressure area under the above-mentioned system is the Z-axis. It only requires the coordinate axes in the formula to be replaced accordingly, which will not be repeated here.

[0168] Below Figure 1 The following examples are used to illustrate:

[0169] Taking the Chang'e-4 relay satellite as an example, this satellite has been operating in the Halo orbit of the Earth-Moon L2 point for a long time. The Z-axis of the satellite in this system is a large umbrella-shaped parabolic antenna, and the Y-axis is a solar sail. The light pressure area of ​​the satellite in the Z-axis direction is much larger than that in the X-axis and Y-axis directions, which is suitable for this application.

[0170] The historical angular momentum telemetry data set is obtained by searching the database for the satellite's angular momentum telemetry data from January to December 2022, and the historical attitude quaternion telemetry data set is obtained through the historical attitude prediction file. Due to the huge amount of data, the telemetry data of the first hour of the first two arcs of the historical angular momentum telemetry data set are displayed at 5-minute intervals as shown in Table 1. The unit of angular momentum is Nms, the same below.

[0171] Table 1: Angular momentum raw telemetry data

[0172]

[0173]

[0174] In the entire data set, we search for the rapid attitude adjustment process TAi (i=1,...,k), the track control process TCi (i=1,...,m), the momentum wheel unloading TWi (i=1,...,n), etc. Taking January 2022 as an example, the results are shown in Table 2:

[0175] Table 2: Events of rapid change in angular momentum (January 2022)

[0176] Beijing Time Key Events 2022-01-03T09:00 Momentum wheel unloading 2022-01-05T19:00 Track maintenance 2022-01-08T18:00 Momentum wheel unloading 2022-01-11T21:00 Momentum wheel unloading 2022-01-15T14:30 Momentum wheel unloading 2022-01-19T08:00 Momentum wheel unloading 2022-01-21T22:30 Track maintenance 2022-01-23T23:00 Momentum wheel unloading 2022-01-26T09:00 Momentum wheel unloading 2022-01-27T17:10 Momentum wheel unloading 2022-01-30T09:20 Quick posture adjustment

[0177] Delete the relevant data of the rapid change process of angular momentum in the historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set. Taking the rapid change process of angular momentum as the node, the historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set are divided into multiple angular momentum periods. According to the rapid change events of angular momentum in Table 2, the angular momentum data in January 2022 can be divided into 10 angular momentum periods as shown in Table 3.

[0178] Table 3: Events of rapid change in angular momentum (January 2022)

[0179]

[0180]

[0181] For each angular momentum cycle, the historical angular momentum telemetry data subset and the historical attitude telemetry data subset are respectively obtained for the angular momentum data of the satellite in the three axes x, y, and z directions under this system [H bx H by H bz ]Use the fifth-order polynomial to fit and get the fitting coefficients corresponding to the three-axis x, y, and z directions of the satellite. The fitting coefficients are as follows:

[0182]

[0183] The fitting coefficients corresponding to the three-axis x, y, and z directions of the satellite are used to construct the fitting formula. The fitting formula is used to determine the missing data points of the historical angular momentum telemetry data subset based on the known data points and the preset step size of the historical angular momentum telemetry data subset, and the completed historical angular momentum telemetry data subset is smoothed. Taking the third angular momentum period as an example, the fitted data is shown in Table 4.

[0184] Table 4: Angular momentum data after interpolation fitting

[0185]

[0186]

[0187] Since the ground has satellite attitude quaternion prediction data, there is no need to perform interpolation fitting processing on the telemetry attitude quaternion.

[0188] In each angular momentum cycle, the angle between the Z axis of the satellite system and the solar vector and the corresponding increment of angular momentum in the system are calculated. When the angular momentum of the satellite is zero, the angular momentum values ​​of the X, Y, and Z axes are (1.5, 1.5, -1.5). With a calculation step of 10 minutes and the fourth angular momentum cycle as an example, the angle between the Z axis and the solar vector and the increment of angular momentum of the X, Y, and Z axes of the system are calculated as shown in Table 5. The unit of angle is degree, the same below.

[0189] Table 5: Angle between Z axis and sun vector and angular momentum increment

[0190]

[0191]

[0192] The angle between the Z axis and the sun vector is As the independent variable, the sine sum model is used to fit the three-direction angular momentum increment ΔH bx , ΔH by , ΔH bz , the fitting coefficient values ​​are as follows, only the Y-axis coefficient a is shown yi ,b yi ,c yi (i=1,...,7), the rest will not be repeated here.

[0193]

[0194] The relationship between the satellite's angular momentum increment on the X-axis of this system and the angle with the sun is as follows: Figure 2 shown.

[0195] The relationship between the satellite's angular momentum increment on the Y axis of this system and the angle with the sun is as follows: Figure 3 shown.

[0196] The relationship between the satellite's angular momentum increment on the Z axis of the system and the angle with the sun is as follows: Figure 4 shown.

[0197] Based on the angular momentum increment, conversion matrix, angular momentum telemetry value and reference angular momentum in this system, the angular momentum forecast value of the satellite's next step in this system is determined. Let the forecast time period be the start time T s =2023-02-01T10:00:00.000, end time T e =2023-02-01T23:00:00.000, start time T s The satellite's three-axis angular momentum telemetry value H under this system yb0 =[1.51.5-1.5]. The satellite attitude during the forecast period is exemplified by the following two common attitude conditions.

[0198] Example 1:

[0199] During the forecast period, the satellite attitude is oriented toward the sun with the positive Z axis deviating 37 degrees from the positive X axis, and the earth is constrained in the XZ plane. The satellite angular momentum forecast results are shown in Table 6.

[0200] Table 6: Angular momentum prediction results for example 1

[0201]

[0202] Example 2:

[0203] During the forecast period, the satellite attitude is oriented toward the moon with the positive Z axis, and the sun is constrained in the XZ plane. The satellite angular momentum forecast results are shown in Table 7.

[0204] Table 7: Example two-angular momentum prediction results

[0205]

[0206]

[0207] The above method is used to predict the angular momentum of the satellite in the lunar relay tracking attitude mode and the biased solar inertial flight attitude mode to verify the correctness of the method.

[0208] Verification 1: Angular momentum prediction in lunar relay tracking attitude mode

[0209] The satellite adjusted its attitude and unloaded on November 23, 2024, and then entered the lunar relay tracking mode. After unloading, the satellite telemetry angular momentum data was collected at 18:00:00.457 on November 23, 2024, and the angular momentum telemetry value was [1.4351.435 -1.435]Nms. Using this set of data as the initial value, the subsequent angular momentum was predicted and compared with the telemetry results. The results are shown in Tables 8, 9, and 10.

[0210] Table 8 0.69-day forecast deviation (2024-11-24T10:27:52)

[0211] type X-axis angular momentum Y-axis angular momentum Z-axis angular momentum Telemetry results 1.435 0.665 -1.505 Forecast results 1.469 0.558 -1.520 Forecast bias 0.034 -0.107 -0.015

[0212] Table 9 1.96-day forecast deviation (2024-11-25T16:59:53)

[0213] type X-axis angular momentum Y-axis angular momentum Z-axis angular momentum Telemetry results 1.680 -1.120 -1.680 Forecast results 1.734 -1.276 -1.662 Forecast bias 0.054 -0.156 0.018

[0214] Table 10 2.96-day forecast deviation (2024-11-26T17:00:00)

[0215] type X-axis angular momentum Y-axis angular momentum Z-axis angular momentum Telemetry results 2.345 -2.310 -1.820 Forecast results 2.368 -2.089 -1.776 Forecast bias 0.023 0.221 0.044

[0216] It can be seen from the above table that the angular momentum deviation is mainly in the Y-axis direction. The forecast deviation for 0.7 days is about -0.11Nms, the forecast deviation for 2.0 days is about -0.16Nms, and the forecast deviation for 3.0 days is about 0.22Nms.

[0217] Verification 2: Angular momentum prediction in the day-biased inertial flight mode

[0218] The satellite was unloaded on December 17, 2024, and then entered the inertial flight attitude with a sun bias. The +Z axis was offset from the sun by 32 degrees, and the corresponding attitude quaternion was [0.159148512 -0.188430480.82269712 -0.51218664]. After unloading, the satellite telemetry angular momentum data was collected at 11:40:03.510 on December 17, 2024, and the angular momentum telemetry value was [1.470 1.505 -1.575]Nms. This set of data was used as the initial value to predict the subsequent angular momentum and compare it with the telemetry results. The results are shown in Tables 11, 12, and 13.

[0219] Table 11 0.94-day forecast deviation (2024-12-18T10:09:55)

[0220] type X-axis angular momentum Y-axis angular momentum Z-axis angular momentum Telemetry results 1.365 0.245 -1.610 Forecast results 1.357 0.231 -1.635 Forecast bias -0.008 -0.014 -0.025

[0221] Table 12 1.98-day forecast deviation (2024-12-19T11:10:04)

[0222] type X-axis angular momentum Y-axis angular momentum Z-axis angular momentum Telemetry results 1.260 -1.085 -1.680 Forecast results 1.240 -1.038 -1.695 Forecast bias -0.020 0.047 -0.015

[0223] Table 13 3.12-day forecast deviation (2024-12-20T14:31:08)

[0224]

[0225]

[0226] As can be seen from the table above, in the inertial flight attitude with a bias toward the sun, the forecast error is relatively small overall. The forecast error for 0.94 days is about 0.03Nms, the forecast error for 2.0 days is about 0.05Nms, and the forecast error for 3.1 days is about 0.10Nms.

[0227] Figure 5 The following is a schematic diagram of the structure of an angular momentum prediction device for a lunar orbit satellite provided in an embodiment of the present application. Figure 5 As shown, the device includes: a transceiver module 501 and a processing module 502.

[0228] The transceiver module 501 is used to obtain the angular momentum telemetry value of the satellite's three axes in this system and the satellite attitude quaternion at the prediction start time, and the satellite attitude quaternion is used to describe the satellite attitude;

[0229] The processing module 502 is used to calculate the transformation matrix of the satellite from the inertial system to the current system based on the satellite attitude quaternion, where the elements in the transformation matrix represent the projection components of the unit vector in the inertial system in the current system;

[0230] The processing module 502 is further used to obtain the angular momentum vector of the satellite in the present system based on the angular momentum telemetry value and the reference angular momentum, and determine the angular momentum vector of the satellite in the inertial system through coordinate transformation based on the conversion matrix and the angular momentum vector in the present system, wherein the reference angular momentum represents the indicated values ​​of the angular momentum in the three-axis directions when the angular momentum of the satellite is zero after the satellite is unloaded;

[0231] The processing module 502 is further used to input the solar direction angle into the angular momentum increment fitting model to obtain the angular momentum increment of the satellite in the present system. The solar direction angle represents the angle between the axis with the largest solar light pressure area of ​​the satellite and the solar vector in the present system. The angular momentum increment fitting model is obtained by fitting the corresponding relationship between the historical solar direction angle and the historical angular momentum increment. The axis with the largest solar light pressure area represents the coordinate axis with the largest solar light pressure area of ​​the satellite in the present system.

[0232] The processing module 502 is also used to determine the predicted angular momentum of the next step of the satellite in the present system based on the angular momentum increment, the conversion matrix, the angular momentum telemetry value and the reference angular momentum in the present system. If the conversion matrix and the reference angular momentum remain unchanged, the predicted angular momentum of the next step is positively correlated with the angular momentum increment in the present system.

[0233] Optionally, the sun direction angle is obtained in the following manner, and the processing module 502 is further used for:

[0234] By looking up the JPL planetary ephemeris, calculate the solar position vector of the Sun in the inertial frame at the start of the prediction;

[0235] Obtain satellite ephemeris data from the database to obtain the satellite position vector in the inertial system at the start time of prediction;

[0236] Based on the vector difference between the sun position vector and the satellite position vector, determine the relative position vector of the sun relative to the satellite in the inertial system;

[0237] Based on the transformation matrix and the relative position vector in the inertial system, the relative position vector of the sun relative to the satellite in this system is determined through coordinate transformation;

[0238] The inverse cosine function operation is performed on the ratio of the component of the relative position vector under this system on the axis with the largest solar light pressure area and the modulus of the relative position vector under this system to obtain the solar direction angle.

[0239] Optionally, the angular momentum increment fitting model is fitted in the following manner, and the processing module 502 is further used for:

[0240] In each direction of the three axes of the satellite under this system, a sine sum model of preset different orders is used to fit each historical solar direction angle in the historical solar direction angle set with the historical angular momentum increment corresponding to the historical solar direction angle in the system, and the fitting coefficient of the satellite in each direction of the three axes is obtained. The historical solar direction angle set includes the historical angle between the axis with the largest solar light pressure area of ​​the satellite and the solar position vector under this system at the time node corresponding to each step in the historical time period;

[0241] The angular momentum increment fitting model is constructed based on the fitting coefficients.

[0242] Optionally, the historical angular momentum increment corresponding to the historical solar direction angle in the present system is obtained by the following method:

[0243] The transceiver module 501 is also used to: obtain a historical telemetry data set;

[0244] The processing module 502 is also used to: perform standardization processing on a historical telemetry data set, where the historical telemetry data set includes multiple satellite times of the satellite in the system within a historical period, historical angular momentum telemetry values ​​corresponding to the satellite times, and historical satellite attitude quaternions corresponding to the satellite times;

[0245] For each moment in the standardized historical telemetry data set, the vector difference between the historical angular momentum telemetry value corresponding to each moment and the reference angular momentum is calculated to obtain the instantaneous angular momentum of the satellite in this system at each moment;

[0246] Based on the transformation matrix and instantaneous angular momentum of the satellite from the inertial system to the local system at each moment, the instantaneous angular momentum of the satellite in the inertial system is determined by coordinate transformation. The transformation matrix of the satellite from the inertial system to the local system at each moment is obtained by calculating the historical satellite attitude quaternion at each moment.

[0247] Calculate the vector difference between the instantaneous angular momentum of the next moment in the inertial system and the instantaneous angular momentum of each moment in the inertial system, and obtain the historical angular momentum increment of the satellite in the inertial system;

[0248] Based on the transformation matrix of the satellite from the inertial system to the current system at each moment and the historical angular momentum increment of the satellite in the inertial system, the historical angular momentum increment of the satellite in the current system is determined by coordinate transformation.

[0249] Optionally, the historical telemetry data set includes a historical angular momentum telemetry data set and a historical attitude quaternion telemetry data set, and the historical telemetry data set is standardized. The processing module 502 is specifically used for:

[0250] Delete the relevant data of the rapid change process of angular momentum in the historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set. The rapid change process of angular momentum at least includes the satellite's rapid attitude adjustment, momentum wheel unloading, and orbit control;

[0251] Taking the rapid change process of angular momentum as a node, the historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set are divided into multiple angular momentum periods, each angular momentum period includes a subset of historical angular momentum telemetry data and a subset of historical attitude telemetry data in different time periods;

[0252] The historical angular momentum telemetry data subset and the historical attitude telemetry data subset within each angular momentum cycle are fitted to complete and smooth the missing data points of the historical angular momentum telemetry data subset and the historical attitude telemetry data subset.

[0253] Optionally, fitting processing is performed on the subset of historical angular momentum telemetry data and the subset of historical attitude telemetry data in each angular momentum period to complete and smooth the missing data points of the subset of historical angular momentum telemetry data and the subset of historical attitude telemetry data. The processing module 502 is specifically used to:

[0254] For each angular momentum period, a preset first-order polynomial is used to fit the components of each historical angular momentum telemetry value in each direction of the three axes of the satellite in this system to obtain a first fitting coefficient, and a first fitting formula is determined based on the first fitting coefficient;

[0255] Using a preset second-order polynomial to fit each component of each historical satellite attitude quaternion to obtain a second fitting coefficient, and determining a second fitting formula based on the second fitting coefficient;

[0256] Using a first fitting formula, based on known data points of the historical angular momentum telemetry data subset and a preset step size, the missing data points of the historical angular momentum telemetry data subset are determined, and the completed historical angular momentum telemetry data subset is smoothed;

[0257] The second fitting formula is used to determine the missing data points of the historical attitude quaternion telemetry data subset based on the known data points of the historical attitude quaternion telemetry data subset and the preset step size, and the completed historical attitude quaternion telemetry data subset is smoothed.

[0258] Optionally, based on the angular momentum increment, the conversion matrix, the angular momentum telemetry value and the reference angular momentum in the present system, the angular momentum prediction value of the satellite for the next step in the present system is determined, and the processing module 502 is specifically used to:

[0259] The transpose of the conversion matrix is ​​multiplied by the transpose of the angular momentum increment in this system to obtain the angular momentum increment of the satellite in the inertial system;

[0260] Based on the angular momentum increment in the inertial system and the vector sum of the angular momentum vector in the inertial system, the angular momentum vector of the satellite's next step in the inertial system is determined;

[0261] The angular momentum vector of the next step length in the inertial system is multiplied by the conversion matrix, and the angular momentum vector of the next step length in the satellite in the current system is determined by coordinate transformation;

[0262] Based on the vector sum of the angular momentum vector of the next step length in this system and the reference angular momentum, the predicted angular momentum value of the next step length of the satellite in this system is determined.

[0263] Based on the same technical concept, an electronic device is also provided in an embodiment of the present application, and the electronic device can realize the function of the angular momentum prediction device of the aforementioned lunar orbit satellite.

[0264] Figure 6 A schematic diagram of the structure of an electronic device provided in an embodiment of the present application.

[0265] At least one processor 601, and a memory 602 connected to the at least one processor 601. The specific connection medium between the processor 601 and the memory 602 is not limited in the embodiment of the present application. Figure 6 In the example, the processor 601 and the memory 602 are connected via a bus 600. The bus 600 is Figure 6 The connection between other components is shown by bold lines, and is not intended to be limiting. The bus 600 can be divided into an address bus, a data bus, a control bus, etc. For ease of representation, Figure 6 Only one thick line is used in the figure, but it does not mean that there is only one bus or one type of bus. Alternatively, the processor 601 can also be called a controller, and there is no limitation on the name.

[0266] In the embodiment of the present application, the memory 602 stores instructions that can be executed by at least one processor 601. The at least one processor 601 can execute the angular momentum prediction method of a lunar orbit satellite discussed above by executing the instructions stored in the memory 602. The processor 601 can implement Figure 6 The functions of each module in the device shown.

[0267] Among them, the processor 601 is the control center of the device, and can use various interfaces and lines to connect the various parts of the entire control device. By running or executing instructions stored in the memory 602 and calling the data stored in the memory 602, the various functions of the device and processing data, the device can be monitored as a whole.

[0268] In one possible design, the processor 601 may include one or more processing units, and the processor 601 may integrate an application processor and a modem processor, wherein the application processor mainly processes an operating system, a driver interface, and application programs, and the modem processor mainly processes wireless communications. It is understandable that the modem processor may not be integrated into the processor 601. In some embodiments, the processor 601 and the memory 602 may be implemented on the same chip, and in some embodiments, they may also be implemented separately on separate chips.

[0269] The processor 601 may be a general-purpose processor, such as a central processing unit (CPU), a digital signal processor, an application-specific integrated circuit, a field programmable gate array or other programmable logic device, a discrete gate or transistor logic device, or a discrete hardware component, and may implement or execute the methods, steps, and logic block diagrams disclosed in the embodiments of the present application. A general-purpose processor may be a microprocessor or any conventional processor, etc. The steps of the angular momentum prediction method for a lunar orbit satellite disclosed in the embodiments of the present application may be directly embodied as being executed by a hardware processor, or may be executed by a combination of hardware and software modules in the processor.

[0270] The memory 602 is a non-volatile computer-readable storage medium that can be used to store non-volatile software programs, non-volatile computer executable programs and modules. The memory 602 may include at least one type of storage medium, such as a flash memory, a hard disk, a multimedia card, a card-type memory, a random access memory (Random Access Memory, RAM), a static random access memory (Static Random Access Memory, SRAM), a programmable read-only memory (Programmable Read Only Memory, PROM), a read-only memory (Read Only Memory, ROM), an electrically erasable programmable read-only memory (Electrically Erasable Programmable Read-Only Memory, EEPROM), a magnetic memory, a disk, an optical disk, etc. The memory 602 is any other medium that can be used to carry or store a desired program code in the form of an instruction or data structure and can be accessed by a computer, but is not limited thereto. The memory 602 in the embodiment of the present application can also be a circuit or any other device that can realize a storage function, for storing program instructions and / or data.

[0271] By programming the processor 601, the code corresponding to the method for predicting the angular momentum of a lunar orbit satellite described in the above embodiment can be fixed into the chip, so that the chip can execute the method when running. Figure 1The embodiment shown is a method for predicting the angular momentum of a lunar orbit satellite. How to design and program the processor 601 is a technique known to those skilled in the art and will not be described in detail here.

[0272] It should be noted here that the above-mentioned electronic device provided in the embodiment of the present application can implement all the method steps implemented in the above-mentioned method embodiment, and can achieve the same technical effect. The parts and beneficial effects of this embodiment that are the same as the method embodiment will not be described in detail here.

[0273] An embodiment of the present application also provides a computer-readable storage medium, which stores computer-executable instructions. The computer-executable instructions are used to enable a computer to execute a method for predicting angular momentum of a lunar orbit satellite in the above embodiment.

[0274] Those skilled in the art will appreciate that the embodiments of the present application may be provided as methods, systems, or computer program products. Therefore, the present application may adopt the form of a complete hardware embodiment, a complete software embodiment, or an embodiment in combination with software and hardware. Moreover, the present application may adopt the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) that contain computer-usable program code.

[0275] The present application is described with reference to the flowcharts and / or block diagrams of the methods, devices (systems), and computer program products according to the present application. It should be understood that each process and / or block in the flowchart and / or block diagram, as well as the combination of the processes and / or blocks in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to generate a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowchart and / or block diagram. Figure 1 Flow or flows and / or boxes Figure 1 A device that performs the functions specified in a box or multiple boxes.

[0276] These computer program instructions may also be stored in a computer-readable memory capable of directing a computer or other programmable data processing device to operate in a specific manner, so that the instructions stored in the computer-readable memory produce an article of manufacture comprising an instruction device, which implements the process Figure 1 Flow or flows and / or boxes Figure 1 The functions specified in one or more boxes.

[0277] These computer program instructions can also be loaded onto a computer or other programmable data processing device so that a series of operating steps are executed on the computer or other programmable device to produce a computer-implemented process, thereby providing instructions for implementing the process. Figure 1 Flow or flows and / or boxes Figure 1 A step that specifies a function in one or more boxes.

[0278] Obviously, those skilled in the art can make various changes and modifications to the present application without departing from the spirit and scope of the present application. Thus, if these modifications and variations of the present application fall within the scope of the claims of the present application and their equivalents, the present application is also intended to include these modifications and variations.

Claims

1. A method for predicting the angular momentum of a lunar orbit satellite, characterized in that: The method comprises: Obtain the angular momentum telemetry values ​​of the three axes of the satellite and the satellite attitude quaternion of the satellite in the current system at the prediction start time, and calculate the transformation matrix of the satellite from the inertial system to the current system based on the satellite attitude quaternion. The satellite attitude quaternion is used to describe the satellite attitude, and the elements in the transformation matrix represent the projection components of the unit vector in the inertial system in the current system; The angular momentum vector of the satellite in the present system is obtained based on the angular momentum telemetry value and the reference angular momentum. The angular momentum vector of the satellite in the inertial system is determined by coordinate transformation based on the conversion matrix and the angular momentum vector in the present system. The reference angular momentum represents the indicated values ​​of the angular momentum in the three-axis directions when the angular momentum of the satellite is zero after the satellite is unloaded. Inputting the solar direction angle into the angular momentum increment fitting model, obtaining the angular momentum increment of the satellite in the present system, wherein the solar direction angle represents the angle between the axis with the largest solar light pressure area of ​​the satellite and the solar position vector in the present system, and the angular momentum increment fitting model is obtained by fitting the corresponding relationship between the historical solar direction angle and the historical angular momentum increment, and the axis with the largest solar light pressure area represents the coordinate axis with the largest solar light pressure area of ​​the satellite in the present system; Based on the angular momentum increment under the present system, the conversion matrix, the angular momentum telemetry value and the reference angular momentum, a predicted angular momentum value of the next step of the satellite under the present system is determined. If the conversion matrix and the reference angular momentum are unchanged, the predicted angular momentum value of the next step is positively correlated with the angular momentum increment under the present system.

2. The method according to claim 1, characterized in that The sun direction angle is obtained in the following way: By looking up the planetary ephemeris, the solar position vector of the sun in the inertial system at the prediction start time is calculated; Acquire satellite ephemeris data from a database to obtain a satellite position vector of the satellite in an inertial system at the start time of the prediction; Determining a relative position vector of the sun relative to the satellite in the inertial system based on a vector difference between a sun position vector and a satellite position vector; Based on the conversion matrix and the relative position vector in the inertial system, determining the relative position vector of the sun relative to the satellite in the system through coordinate transformation; The solar direction angle is obtained by performing an inverse cosine function operation on the ratio between the component of the relative position vector under the present system on the axis with the largest solar light pressure area and the modulus of the relative position vector under the present system.

3. The method according to claim 1, characterized in that The angular momentum increment fitting model is obtained by fitting in the following way: In each direction of the three axes of the satellite under the present system, a sine sum model of preset different orders is used to fit each historical solar direction angle in the historical solar direction angle set and the historical angular momentum increment corresponding to the historical solar direction angle in the present system, so as to obtain the fitting coefficient of the satellite in each direction of the three axes, wherein the historical solar direction angle set includes the historical angles of the axis with the largest solar light pressure area of ​​the satellite under the present system and the solar position vector at the time node corresponding to each step in the historical time period; The angular momentum increment fitting model is constructed based on the fitting coefficients.

4. The method according to claim 3, characterized in that The historical solar direction angle is obtained in the historical angular momentum increment corresponding to the present system by the following method: Acquire a historical telemetry data set and perform standardization processing on the historical telemetry data set, wherein the historical telemetry data set includes a plurality of satellite times of the satellite in the system during the historical time period, historical angular momentum telemetry values ​​corresponding to the satellite times, and historical satellite attitude quaternions corresponding to the satellite times; For each moment in the standardized historical telemetry data set, calculate the vector difference between the historical angular momentum telemetry value corresponding to each moment and the reference angular momentum, and obtain the instantaneous angular momentum of the satellite in the system at each moment; Based on the transformation matrix of the satellite from the inertial system to the local system at each moment and the instantaneous angular momentum, the instantaneous angular momentum of the satellite in the inertial system is determined by coordinate transformation, wherein the transformation matrix of the satellite from the inertial system to the local system at each moment is obtained by calculating the historical satellite attitude quaternion at each moment; Calculating the vector difference between the instantaneous angular momentum at the next moment of each moment in the inertial system and the instantaneous angular momentum at each moment in the inertial system, to obtain the historical angular momentum increment of the satellite in the inertial system; Based on the transformation matrix of the satellite from the inertial system to the local system at each moment and the historical angular momentum increment of the satellite in the inertial system, the historical angular momentum increment of the satellite in the local system is determined through coordinate transformation.

5. The method according to claim 4, characterized in that The historical telemetry data set includes a historical angular momentum telemetry data set and a historical attitude quaternion telemetry data set, and the standardization processing of the historical telemetry data set specifically includes: Deleting the historical angular momentum telemetry data set and the relevant data of the rapid change process of angular momentum in the historical attitude quaternion telemetry data set, wherein the rapid change process of angular momentum at least includes rapid attitude adjustment of the satellite, momentum wheel unloading, and orbit control; Taking the angular momentum rapid change process as a node, dividing the historical angular momentum telemetry data set and the historical attitude quaternion telemetry data set into a plurality of angular momentum cycles, each angular momentum cycle including a subset of historical angular momentum telemetry data and a subset of historical attitude telemetry data in different time periods; The historical angular momentum telemetry data subset and the historical attitude telemetry data subset within each angular momentum cycle are fitted to complete and smooth the missing data points of the historical angular momentum telemetry data subset and the historical attitude telemetry data subset.

6. The method according to claim 5, characterized in that The fitting process is performed on the subset of historical angular momentum telemetry data and the subset of historical attitude telemetry data in each angular momentum cycle, and the missing data points of the subset of historical angular momentum telemetry data and the subset of historical attitude telemetry data are supplemented and smoothed, specifically including: For each angular momentum period, a preset first-order polynomial is used to fit the components of each historical angular momentum telemetry value in each direction of the three axes of the satellite in this system to obtain a first fitting coefficient, and a first fitting formula is determined based on the first fitting coefficient; Using a preset second-order polynomial to fit each component of each historical satellite attitude quaternion to obtain a second fitting coefficient, and determining a second fitting formula based on the second fitting coefficient; Using the first fitting formula, based on the known data points of the historical angular momentum telemetry data subset and a preset step size, determining the missing data points of the historical angular momentum telemetry data subset, and smoothing the completed historical angular momentum telemetry data subset; The second fitting formula is used to determine the missing data points of the historical posture quaternion telemetry data subset based on the known data points of the historical posture quaternion telemetry data subset and the preset step size, and the completed historical posture quaternion telemetry data subset is smoothed.

7. The method according to claim 1, characterized in that Determining the predicted angular momentum value of the next step of the satellite in the present system based on the angular momentum increment in the present system, the conversion matrix, the angular momentum telemetry value and the reference angular momentum specifically includes: Performing a product operation on the transpose of the conversion matrix and the transpose of the angular momentum increment in the present system to obtain the angular momentum increment of the satellite in the inertial system; Determine an angular momentum vector of a next step of the satellite in the inertial system based on the angular momentum increment in the inertial system and the vector sum of the angular momentum vector in the inertial system; The transformation matrix is ​​multiplied by the angular momentum vector of the next step length in the inertial system, and the angular momentum vector of the next step length of the satellite in the current system is determined by coordinate transformation; Based on the angular momentum vector of the next step in the present system and the vector sum of the reference angular momentum, a predicted angular momentum value of the next step of the satellite in the present system is determined.

8. An angular momentum prediction device for a lunar orbit satellite, characterized in that: The device comprises: A transceiver module is used to obtain the angular momentum telemetry value of the satellite's three axes in this system and the satellite attitude quaternion at the prediction start time, and the satellite attitude quaternion is used to describe the satellite attitude; A processing module, used for calculating a transformation matrix of the satellite from the inertial system to the local system based on the satellite attitude quaternion, wherein the elements in the transformation matrix represent the projection components of the unit vector in the inertial system in the local system; The processing module is further used to obtain the angular momentum vector of the satellite in the present system based on the angular momentum telemetry value and the reference angular momentum, and determine the angular momentum vector of the satellite in the inertial system through coordinate transformation based on the conversion matrix and the angular momentum vector in the present system, wherein the reference angular momentum represents the indicated values ​​of the angular momentum in the three-axis directions when the angular momentum of the satellite is zero after the satellite is unloaded; The processing module is further used to input the solar direction angle into the angular momentum increment fitting model to obtain the angular momentum increment of the satellite in the present system, wherein the solar direction angle represents the angle between the axis with the largest solar light pressure area of ​​the satellite and the solar position vector in the present system, and the angular momentum increment fitting model is obtained by fitting the corresponding relationship between the historical solar direction angle and the historical angular momentum increment, and the axis with the largest solar light pressure area represents the coordinate axis with the largest solar light pressure area of ​​the satellite in the present system; The processing module is further used to determine the predicted angular momentum of the next step of the satellite in the present system based on the angular momentum increment in the present system, the conversion matrix, the angular momentum telemetry value and the reference angular momentum. If the conversion matrix and the reference angular momentum remain unchanged, the predicted angular momentum of the next step is positively correlated with the angular momentum increment in the present system.

9. An electronic device comprising a memory, a processor and a computer program stored in the memory and executable on the processor, characterized in that: When the processor executes the computer program, the method according to any one of claims 1 to 7 is implemented.

10. A computer-readable storage medium, characterized in that: The computer-readable storage medium stores a computer program, and the computer program is used to cause the computer to execute the method according to any one of claims 1 to 7.

11. A computer program product, characterized in that When the computer program product is called by a computer, the computer is caused to execute the method according to any one of claims 1 to 7.

Citation Information

Patent Citations

  • Satellite on-orbit angular momentum management and control method with attitude switching function

    CN112572835A

  • Three-axis angular momentum unloading method for geosynchronous orbit satellite

    CN115626304A

  • Tracking Neighboring Quasi-Satellite Orbits Around Mars's Moon Phobos

    US20210206517A1