Sub-unit light pressure perturbation calculation method for antenna satellite

Abstract

The invention discloses a sub-unit light pressure perturbation calculation method for an antenna satellite, and relates to a light pressure perturbation calculation method for the antenna satellite. The invention provides a subunit light pressure perturbation calculation method for an antenna satellite in order to solve the problem of poor satellite orbit prediction and attitude control precision caused by inaccurate solar light pressure calculation of a large satellite antenna. The sub-unit light pressure perturbation calculation method for the antenna satellite specifically comprises the following steps: step 1, determining a sun exposure factor in an ECI coordinate system according to the position relation of the sun, the earth and the antenna satellite; 2, calculating the light pressure projection area of the satellite; step 3, calculating light pressure perturbative force of each part according to the sun exposure factor in the step 1 and the light pressure projection area in the step 2; and step 4, judging whether a set time length is reached or not, if not, repeatedly executing the step 1 to the step 3, and otherwise, ending. The method is applied to the field of antenna satellite light pressure perturbation acquisition.

Description

Technical Field

[0001] The invention relates to a method for calculating light pressure perturbation of an antenna satellite. Background Art

[0002] With the dramatic increase in global communication needs and the rapid development of space technology, building a high-performance communication infrastructure has become a top priority. In this context, large satellite antennas, as a core component of space infrastructure, have become increasingly important. These antennas not only support data and communication services around the earth, but are also key technologies for achieving global wide-area coverage, deep space exploration missions, and disaster emergency response. In addition, with the continuous expansion of earth observation, space science research, and commercial satellite services, the capabilities of large satellite antennas in data collection and other aspects have become the key to promoting progress in these fields. By improving the capacity and stability of space communication networks, large antennas enable satellites to process and transmit massive amounts of data, supporting complex scientific research missions and commercial operations. From climate monitoring to astrophysical research, from global positioning systems to broadband Internet services, the role of large satellite antennas is everywhere, and the innovation and application of their technology is directly related to the operating efficiency of all walks of life and the well-being of all mankind. The development of large satellite antenna technology is not only a reflection of the country's scientific and technological strength, but also an indispensable part of the global communication development trend. As the core of satellite antenna design and control, dynamic models play a vital role in improving mission success, reducing risks, and optimizing orbital operations.

[0003] For large space satellites equipped with large structures such as antennas, the accuracy of the dynamic model becomes an important factor in ensuring their structural stability and functionality. Due to the uniqueness and complexity of the structure of such satellites and their large surface area, solar pressure becomes the main factor of their perturbation. Traditional solar pressure calculations usually set the satellite surface mass ratio to a constant value. Although this simplifies the design and calculation in some aspects, the fixed surface mass ratio also means that the satellite's sensitivity to solar pressure is also fixed, sacrificing the accuracy of the dynamic model. Overly simplified models may not be able to accurately predict and respond to the dynamic conditions actually encountered, which may bring additional challenges to orbit maintenance and attitude control. Summary of the invention

[0004] The purpose of the present invention is to solve the problem of inaccurate calculation of solar radiation pressure of large satellite antennas, which leads to poor accuracy of satellite orbit prediction and attitude control, and propose a sub-unit solar radiation pressure perturbation calculation method for antenna satellites.

[0005] A method for calculating the perturbation of light pressure of a satellite antenna comprises the following steps:

[0006] Step 1: In the ECI coordinate system, determine the insolation factor based on the positional relationship between the sun, the earth, and the antenna satellite;

[0007] Step 2: Calculate the light pressure projection area of ​​the satellite;

[0008] Step 3, calculating the light pressure perturbation force on each part according to the exposure factor of step 1 and the light pressure projection area of ​​step 2;

[0009] Step 4: Determine whether the set time length has been reached. If not, repeat steps 1 to 3. Otherwise, end.

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

[0011] 1. The present invention aims to solve the problem that the large satellite antenna structure is complex and it is difficult to accurately calculate its light pressure perturbation. By decomposing the satellite into several key units and calculating the solar light pressure received by each unit separately, the overall light pressure perturbation can be predicted more accurately. This method allows the specific attributes of each unit (such as area, material, shape and orientation) to be considered in detail, thereby improving the accuracy and reliability of the overall model.

[0012] 2. The present invention is conducive to improving the accuracy of orbit prediction and attitude control of satellites and spacecraft. By calculating the influence of solar pressure with high precision, the accuracy of orbit prediction can be significantly improved, the attitude control capability can be enhanced, and the use of fuel can be optimized. This is especially critical for tasks that require high-precision pointing, such as earth observation, astronomical observation, and intersatellite communication. In addition, the use of solar pressure for orbit and attitude control can also reduce the consumption of propellant and extend the life of the satellite, thereby improving the economy and sustainability of the mission.

[0013] 3. The present invention is conducive to improving mission planning and execution. High-precision solar pressure calculations help to formulate more reasonable mission plans, especially in satellite formation flying and collaborative observation missions. For long-term missions, accurate solar pressure models help predict and respond to long-term orbital drift and attitude changes, and improve mission success rates. At the same time, accurate assessment and prediction of orbital changes caused by solar pressure can help enhance risk management and collision avoidance, improve the accuracy of collision risk assessment, and plan effective collision avoidance maneuvers to reduce collision risks in space environments.

[0014] 4. The present invention is beneficial to improving the quality of scientific experiments and data, especially for scientific satellites equipped with precision instruments, and can correct small disturbances caused by solar pressure to ensure the accuracy and reliability of data. BRIEF DESCRIPTION OF THE DRAWINGS

[0015] Figure 1 is a flow chart of the method of the present invention;

[0016] Figure 2is the antenna parabola fitting diagram, with the center point M of the antenna canopy as the origin, the speed direction of the antenna satellite as the x” axis, the direction pointing to the center of the earth as the z” axis, and the y” axis determined by the right-hand rule; points L, F, R, and B are four points symmetrically located on the boundary of the antenna canopy; θ represents the direction of the sun’s rays; Represents the curve on the antenna umbrella surface that is tangent to the sun's rays;

[0017] Figure 3 It is the projection diagram of the antenna parabola at different solar incidence angles;

[0018] Figure 4 The incident angle variation diagram of the antenna projection area is calculated for the macro model MM (Macro Model) and the high precision analytical light pressure perturbation model HPAM (High Precision Analytical Model) of the present invention;

[0019] Figure 5 This is a graph showing the changes in the projection area of ​​the satellite antenna on different dates;

[0020] Figure 6 At the vernal equinox, the antenna has a roll angle of 6° and is subjected to the light pressure perturbation diagram. Fx, Fy, and Fz are the three components of the solar light pressure perturbation force exerted on the antenna satellite in the orbital coordinate system, and Tx, Ty, and Tz are the three components of the solar light pressure perturbation moment exerted on the antenna satellite in the orbital coordinate system. DETAILED DESCRIPTION

[0021] Specific implementation method 1: Combination Figure 1 The present embodiment is described. The present embodiment is a method for calculating the optical pressure perturbation of a satellite antenna unit, which specifically includes the following steps:

[0022] The technical solution adopted by the present invention to solve the above technical problems is: a unit-by-unit light pressure perturbation calculation method that comprehensively considers the shape, structure, and posture of the satellite antenna, and the method specifically includes the following steps:

[0023] Step 1: In the ECI coordinate system, determine the insolation factor based on the positional relationship between the sun, the earth, and the antenna satellite;

[0024] Step 2: Calculate the light pressure projection area of ​​the satellite;

[0025] Step 3, calculating the light pressure perturbation force on each part according to the exposure factor of step 1 and the light pressure projection area of ​​step 2;

[0026] Step 4: Determine whether the set time length has been reached. If not, repeat steps 1 to 3. Otherwise, end.

[0027] Specific implementation method 2: This implementation method is different from the specific implementation method 1 in that: in step 1, in the ECI coordinate system, the insolation factor is determined according to the positional relationship between the sun, the earth, and the antenna satellite; the specific process is:

[0028] The ECI coordinate system is a geocentric inertial coordinate system. The origin of the ECI coordinate system is located at the center of mass of the earth. The z-axis is along the direction of the earth's rotation axis, pointing from the center of the earth to the North Pole. The x-axis is in the equatorial plane, pointing from the center of the earth to the vernal equinox. The y-axis, x-axis, and z-axis satisfy the right-hand rule.

[0029] Let the ECI coordinate of the antenna satellite be r sat =[r satx r saty r satz ] T , the ECI coordinate of the sun is r sun =[r sunx r suny r sunz ] T , based on the ECI coordinates of the antenna satellite and the ECI coordinates of the sun, the angle α between the antenna satellite and the sun relative to the center of the earth is obtained;

[0030] Based on the angle α between the antenna satellite and the sun relative to the center of the earth, the determination condition for the satellite's exposure factor k=0 in the earth's shadow area is obtained, and the expression is:

[0031]

[0032] Among them, r e is the radius of the Earth; the superscript T indicates transposition;

[0033] If the above determination condition (2) is not met, the exposure factor k=1.

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

[0035] Specific implementation method three: This implementation method is different from specific implementation methods one or two in that: the angle α between the antenna satellite and the sun relative to the center of the earth is obtained based on the ECI coordinates of the antenna satellite and the ECI coordinates of the sun; the specific process is:

[0036] Then the angle α between the antenna satellite and the sun relative to the center of the earth is:

[0037]

[0038] In the formula, |r sat ∣ represents the modulus of the satellite position vector in the ECI coordinate system; · represents the multiplication sign, r satx represents the x-axis coordinate of the antenna satellite in the ECI coordinate system, r satyrepresents the y-axis coordinate of the antenna satellite in the ECI coordinate system, r satz represents the z-axis coordinate of the antenna satellite in the ECI coordinate system, r sunx represents the x-axis coordinate of the sun in the ECI coordinate system, r suny represents the y-axis coordinate of the sun in the ECI coordinate system, r sunz Represents the z-axis coordinate of the sun in the ECI coordinate system.

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

[0040] Specific implementation method 4: This implementation method is different from any one of the specific implementation methods 1 to 3 in that: in step 2, the light pressure projection area of ​​the satellite is calculated; the specific process is:

[0041] Step 21: According to the antenna satellite structure, the antenna satellite is divided into a solar wing unit, a star unit and an antenna unit;

[0042] Step 22: Consider the solar wing unit as a flat plate structure with a length of m and a width of n, and calculate the light pressure projection area of ​​the solar wing;

[0043] Step 2: Consider the star unit as a rectangular parallelepiped structure with a length, width and height of l, w and h respectively, and calculate the light pressure projection area of ​​the star unit in the illumination direction;

[0044] Step 24: Consider the antenna unit as a parabolic structure. Due to the complexity of the parabolic structure, obtain the light pressure projection area of ​​the antenna unit in the illumination direction.

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

[0046] Specific implementation method 5: This implementation method is different from specific implementation methods 1 to 4 in that: in step 22, the solar wing unit is regarded as a flat plate structure with a length of m and a width of n, and the light pressure projection area of ​​the solar wing is calculated;

[0047] The specific process is:

[0048] Since the main structure of the large satellite antenna is composed of a star body, two symmetrical solar wings and an antenna, it is divided into solar wing units, star body units and antenna units for separate light pressure modeling, and the light pressure projection area in the direction of illumination is calculated.

[0049] Let the antenna satellite velocity vector be v, then the unit vector of the antenna satellite orbital coordinate system is as follows:

[0050]

[0051] In the formula, (x o ,y o ,z o) represents the three-axis unit vector of the antenna satellite orbit coordinate system in the ECI coordinate system; |v| represents the modulus of the antenna satellite operating velocity vector;

[0052] The origin of the antenna satellite orbital coordinate system x'y'z' is the center of mass of the antenna satellite, the x' axis is the velocity direction of the antenna satellite, the z' axis points to the center of the earth, and the y' axis is determined by the right-hand rule;

[0053] Based on the unit vector of the antenna satellite orbit coordinate system, the incident angle β of the sunlight is obtained;

[0054] Based on the incident angle β of sunlight, the light pressure projection area of ​​the solar wing is obtained, which is expressed as

[0055] S sol =m×n×|sinβ| (4).

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

[0057] Specific implementation method 6: This implementation method is different from any one of the specific implementation methods 1 to 5 in that: the incident angle β of the sunlight is obtained based on the unit vector of the antenna satellite orbit coordinate system, and the expression is:

[0058]

[0059] In the formula, · represents the dot product.

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

[0061] Specific implementation method 7: This implementation method is different from specific implementation methods 1 to 6 in that: in steps 2 and 3, the star unit is regarded as a rectangular parallelepiped structure with a length, width and height of l, w and h respectively, and the light pressure projection area of ​​the star unit in the illumination direction is calculated; the specific process is:

[0062] When the sun is at the vernal equinox, for any incident angle β of sunlight, calculate the light pressure projection area of ​​each sun-receiving surface of the rectangular structure; the expression is:

[0063] S WX =l×w×|sinβ|+w×h×|cosβ| (6)

[0064] In the formula, S WX It represents the light pressure projection area of ​​the star unit in the illumination direction.

[0065] The other steps and parameters are the same as those in Specific Implementations 1 to 6-1.

[0066] Specific implementation eight: This implementation differs from any one of specific implementations one to seven in that: in step 24, the antenna unit is regarded as a parabolic structure, and due to the complexity of the parabolic structure, the light pressure projection area of ​​the antenna unit in the illumination direction is obtained;

[0067] The specific process is:

[0068] 1) The antenna surface can be modeled as a parabola, which is a special form of a quadratic surface. The general equation of a quadratic surface is:

[0069] a0+a1x”+a2y”+a3z”+a4x”y”+a5y”z”+a6z”x”+a7x” 2 +a8y” 2 +a9z” 2 =0 (7)

[0070] In the formula, a0, a1, a2, a3, a4, a5, a6, a7, a8, and a9 represent parameters; x', y', and z' represent the coordinate values ​​of the three axes in the Cartesian coordinate system;

[0071] The Cartesian coordinate system x"y"z" is:

[0072] Taking the center of the antenna canopy as the origin, the speed direction of the antenna satellite is the x” axis, the direction pointing to the center of the earth is the z” axis, and the y” axis is determined by the right-hand rule;

[0073] The surface equation contains 10 unknown parameters;

[0074] Parameter a0 determines the position offset of the surface;

[0075] Parameters a1 to a3 are linear terms that affect the surface translation;

[0076] Parameters a4 to a6 represent the correlation between the two directions and affect the shear effect between the two directions;

[0077] Parameters a7 to a9 affect the curvature in three directions.

[0078] Solving the antenna surface equation requires more than ten parameters. Finite element software can be used to analyze the various parts of the antenna. Then fit the antenna parabola, such as Figure 2 As shown; the parabola equation of the antenna can be obtained;

[0079] Finite element software is used to fit the general equation of the quadratic surface (Formula 7), and the parabolic equation of the antenna is obtained as follows:

[0080] 0 = ax” 2 +ay” 2 -z" (8)

[0081] Where a represents a parameter; N is the total height of the antenna structure; z”≤N;

[0082] For the form 0 = ax" 2 +ay” 2 -z" antenna unit, assuming the direction of the sun's rays is Find the light pressure projection area of ​​the antenna unit in the direction of illumination, that is, find the area of ​​the parabola with the normal vector The projected area of ​​the plane y”cosθ+z”sinθ=0;

[0083] In the formula, θ represents the direction of sunlight;

[0084] The solar susceptibility force is mainly determined by the projected area of ​​the satellite exposed to the incident sunlight. The effective irradiated area is related to the solar incident angle θ. When θ = 0°, the sunlight directly illuminates the side of the antenna, resulting in the smallest projected area. When θ = 90°, the sunlight directly illuminates the antenna surface, resulting in the largest projected area. The incident angle changes continuously with the movement of the geosynchronous satellite and the sun, such as Figure 3 shown.

[0085] Parabola equation 0 = ax" 2 +ay” 2 The surface gradient of "-z" is

[0086] 2) When When , the light pressure projection area of ​​the antenna unit in the illumination direction is calculated based on the light pressure projection area formula (19);

[0087] 3) When θ is equivalent to π-θ, ​​substitute it into the light pressure projection area formula (19) to calculate the light pressure projection area of ​​the antenna unit in the illumination direction;

[0088] where θ is given by get;

[0089] 4) When θ is equivalent to θ-π, substitute it into the light pressure projection area formula (19) to calculate the light pressure projection area of ​​the antenna unit in the illumination direction;

[0090] where θ is given by get;

[0091] 5) When θ is equivalent to 2π-θ, substitute it into the light pressure projection area formula (19) to calculate the light pressure projection area of ​​the antenna unit in the illumination direction;

[0092] where θ is given by get.

[0093] The other steps and parameters are the same as those in Specific Implementation 1 to 7-1.

[0094] Specific implementation method 9: This implementation method is different from the specific implementation methods 1 to 8 in that: When , the light pressure projection area of ​​the antenna unit in the illumination direction is calculated based on the light pressure projection area formula (19); the specific process is:

[0095] 21) When When , since the parabola surface is continuous, the contour curve of the projected area is obtained by projecting the gradient and the normal vector perpendicularly.

[0096]

[0097] In the formula, Represents the parabola equation 0 = ax" 2 +ay” 2 - surface gradient of z”, z”≤N; represents a normal vector; represents a dot product;

[0098] The tangent plane is

[0099]

[0100] Substituting (10) into the parabola equation 0 = ax" 2 +ay” 2 -z", get the curve

[0101]

[0102] To find the curve The projection of (x”, y”, z”) is replaced by Substituting into equation (11), eliminating k, we get the curve The equation of the projection contour curve is

[0103]

[0104] Since the parabola is not a complete closed surface, consider the non-tangent case;

[0105] 22) The parabola boundary Γ is

[0106]

[0107] 23) When When , the contour curve of the projection of the antenna unit is the projection of the boundary Γ on the projection plane y”cosθ+z”sinθ=0;

[0108] When θ = 0, the curve Go to Curve

[0109] Find the critical angle θ0; the specific process is:

[0110] when and the parabola equation 0 = ax" 2 +ay” 2 -z” has no intersection point, the projected contour curve is the projection of the boundary Γ on y”cosθ+z”sinθ=0, and the critical value Critical Angle

[0111] when When , the light ray has no tangent point with the parabola, then the projected contour curve is the projection of the boundary Γ on y”cosθ+z”sinθ=0;

[0112] When 0≤θ≤θ0, the projection is complicated and the projection is a curve With curve Projection on y”cosθ+z”sinθ=0;

[0113] 24) Find the coordinates of point O and point P; the specific process is:

[0114] Simultaneous equations

[0115]

[0116] Get the coordinates of point O and point P coordinates

[0117] 25) Based on the coordinates of point O and point P in 24), calculate The projected area S1 of the curve formed with the straight line OP is expressed as:

[0118]

[0119] In the formula, t1 represents the intermediate variable; · represents the dot product;

[0120] 26) Based on the coordinates of point O and point P in 24), for The projected area S2 of the curve formed with the straight line OP is expressed as:

[0121]

[0122] curve With curve The projected area of ​​the enclosed curve at ycosθ+zsinθ=0 is

[0123] S=S1+S2 (18)

[0124] 27) Based on 23), 25) and 26) and the light pressure projection area formula (19), calculate the light pressure projection area of ​​the antenna unit in the illumination direction; the expression is

[0125]

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

[0127] Specific implementation method 10: This implementation method is different from any one of specific implementation methods 1 to 9 in that: in step 3, the light pressure perturbation force received by each part is calculated according to the exposure factor of step 1 and the light pressure projection area of ​​step 2; the specific process is:

[0128] If k = 0, the pressure force of sunlight is 0;

[0129] If k = 1, the solar pressure force is not 0, and the solar pressure force F is calculated. SPP , the calculation formula is:

[0130]

[0131] Among them, C R is the satellite surface reflection coefficient, which characterizes the satellite surface material's ability to absorb and reflect sunlight; ρ SR is the solar light pressure constant, generally 4.56×10 -6 N / m 2 ; k is the exposure factor; during the satellite operation, when the satellite is in the shadow area on the back of the earth, it will not be exposed to sunlight. In this case, k = 0; when the satellite is not in the shadow area, k = 1; AU is 1 astronomical unit, 1AU = 1.495978707×10 11 ;|r sun | is the modulus of the sun's position vector in the ECI coordinate system; S Ref is the satellite's light pressure projection area; e D is the unit vector pointing from the satellite to the sun.

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

[0133] Example:

[0134] According to the first embodiment, the algorithm of the present invention is experimented in a simulation environment. Using the projection area calculation method of the macroscopic model and the high-precision analytical light pressure perturbation dynamic model, the relationship between the projection area and the incident angle is obtained through simulation analysis, such as Figure 4 As shown in Figure 2, the two methods show similar trends, with the maximum error reaching 1000m at 30°. 2The present invention is more scientific and rigorous, and more in line with the actual situation of light pressure perturbation. On March 21 and April 21, the light pressure projection area of ​​the antenna satellite was simulated for 24 hours using a macroscopic model and a high-precision analytical light pressure perturbation model. Figure 5 On March 21, 2020, the maximum and minimum projection areas of the antenna surface calculated by the two models were equal, with a maximum difference of 1314.7m. 2 The difference between the entire antenna satellite and the antenna surface is very small, mainly due to the influence of the satellite body and the solar wing, indicating that their shading area is negligible. After one month of in-orbit operation, the maximum difference between the two models is 1314.8m 2 , which is not much different from March 21. Although the maximum projected area on the antenna surface is equal on both dates, the difference in the minimum value is about 694.2m 2 This shows that as the position of the sun changes along the ecliptic plane, the photoperturbation dynamics of the antenna satellite becomes more complex, increasing the error in the macroscopic model. In contrast, the high-precision analytical photoperturbation dynamics model improves the calculation accuracy.

[0135] Figure 6 This is the perturbation force and torque diagram of the light pressure obtained from the on-orbit simulation of the vernal equinox. The satellite has a roll angle of 6°. At the beginning of the simulation, the earth, the satellite and the sun are in a straight line, and the Z axis of the satellite always points to the center of the earth. From 6 to 18 hours, the angle between the Z axis and the solar radiation pressure is less than 90°, the perturbation force is positive, and it reaches a peak at 12 hours. From 18 hours to 0 hours and then to 6 hours, the angle is greater than 90°, so the perturbation force is negative, reaching a peak at 0 / 24 hours. From 0 to 12 hours, the satellite moves from aphelion to perihelion, and the solar radiation pressure behaves as resistance along the velocity direction (X axis), which is a negative value. At 6 hours, the solar radiation pressure is the smallest because the incident angle of the antenna surface is about 90°. From 18 to 0 hours and then to 6 hours, the satellite moves from perihelion to aphelion, and the solar radiation pressure along the velocity direction (X axis) is a positive thrust. Since the flip angle of the satellite along the Z axis is 6°, the solar radiation pressure will not pass through the center of mass, thus generating a Z-axis torque. The solar radiation pressure curve shows an equal-length shadow effect around 0 hours.

[0136] The present invention may also have many other embodiments. Without departing from the spirit and essence of the present invention, those skilled in the art may make various corresponding changes and modifications based on the present invention, but these corresponding changes and modifications should all fall within the scope of protection of the claims attached to the present invention.

Claims

1. A method for calculating the optical pressure perturbation of a satellite antenna, characterized in that: The method specifically comprises the following steps: Step 1: In the ECI coordinate system, determine the insolation factor based on the positional relationship between the sun, the earth, and the antenna satellite; Step 2: Calculate the light pressure projection area of ​​the satellite; Step 3, calculating the light pressure perturbation force on each part according to the exposure factor of step 1 and the light pressure projection area of ​​step 2; Step 4: Determine whether the set time length has been reached. If not, repeat steps 1 to 3. Otherwise, end.

2. The method for calculating the optical pressure perturbation of a satellite antenna according to claim 1, characterized in that: In the step 1, in the ECI coordinate system, the insolation factor is determined according to the positional relationship among the sun, the earth, and the antenna satellite; the specific process is: The ECI coordinate system is a geocentric inertial coordinate system. The origin of the ECI coordinate system is located at the center of mass of the earth. The z-axis is along the direction of the earth's rotation axis, pointing from the center of the earth to the North Pole. The x-axis is in the equatorial plane, pointing from the center of the earth to the vernal equinox. The y-axis, x-axis, and z-axis satisfy the right-hand rule. Let the ECI coordinate of the antenna satellite be r sat =[r satx r saty r satz ] T , the ECI coordinate of the sun is r sun =[r sunx r suny r sunz ] T , based on the ECI coordinates of the antenna satellite and the ECI coordinates of the sun, the angle α between the antenna satellite and the sun relative to the center of the earth is obtained; Based on the angle α between the antenna satellite and the sun relative to the center of the earth, the determination condition for the satellite's exposure factor k=0 in the earth's shadow area is obtained, and the expression is: Among them, r e is the radius of the Earth; the superscript T indicates transposition; If the above determination condition (2) is not met, the exposure factor k=1.

3. The method for calculating the optical pressure perturbation of a satellite antenna according to claim 2, characterized in that: The angle α between the antenna satellite and the sun relative to the center of the earth is obtained based on the ECI coordinates of the antenna satellite and the ECI coordinates of the sun. The specific process is: Then the angle α between the antenna satellite and the sun relative to the center of the earth is: In the formula, |r sat ∣ represents the modulus of the satellite position vector in the ECI coordinate system; ·Represents the multiplication sign, r satx represents the x-axis coordinate of the antenna satellite in the ECI coordinate system, r saty represents the y-axis coordinate of the antenna satellite in the ECI coordinate system, r satz represents the z-axis coordinate of the antenna satellite in the ECI coordinate system, r sunx represents the x-axis coordinate of the sun in the ECI coordinate system, r suny represents the y-axis coordinate of the sun in the ECI coordinate system, r sunz Represents the z-axis coordinate of the sun in the ECI coordinate system.

4. The method for calculating the optical pressure perturbation of a satellite antenna according to claim 3, characterized in that: In the step 2, the light pressure projection area of ​​the satellite is calculated; the specific process is: Step 21: According to the antenna satellite structure, the antenna satellite is divided into a solar wing unit, a star unit and an antenna unit; Step 22: Consider the solar wing unit as a flat plate structure with a length of m and a width of n, and calculate the light pressure projection area of ​​the solar wing; Step 2: Consider the star unit as a rectangular parallelepiped structure with a length, width and height of l, w and h respectively, and calculate the light pressure projection area of ​​the star unit in the illumination direction; Step 24: Consider the antenna unit as a parabolic structure and obtain the light pressure projection area of ​​the antenna unit in the illumination direction.

5. The method for calculating the optical pressure perturbation of a satellite antenna according to claim 4, characterized in that: In step 22, the solar wing unit is regarded as a flat plate structure with a length of m and a width of n, and the light pressure projection area of ​​the solar wing is calculated; the specific process is: Let the antenna satellite velocity vector be v, then the unit vector of the antenna satellite orbital coordinate system is as follows: In the formula, (x o ,y o ,z o ) represents the three-axis unit vector of the antenna satellite orbit coordinate system in the ECI coordinate system; |v| represents the modulus of the antenna satellite operating velocity vector; The origin of the antenna satellite orbital coordinate system x'y'z' is the center of mass of the antenna satellite, the x' axis is the velocity direction of the antenna satellite, the z' axis points to the center of the earth, and the y' axis is determined by the right-hand rule; Based on the unit vector of the antenna satellite orbit coordinate system, the incident angle β of the sunlight is obtained; Based on the incident angle β of sunlight, the light pressure projection area of ​​the solar wing is obtained, which is expressed as S sol =m×n×|sinβ| (4)。 6. The method for calculating the optical pressure perturbation of antenna satellites according to claim 5, characterized in that: The unit vector based on the antenna satellite orbit coordinate system is used to obtain the incident angle β of sunlight, which is expressed as: In the formula, · represents the dot product.

7. The method for calculating the optical pressure perturbation of a satellite antenna according to claim 6, characterized in that: In the steps 2 and 3, the star unit is regarded as a rectangular parallelepiped structure with a length, width and height of l, w and h respectively, and the light pressure projection area of ​​the star unit in the illumination direction is calculated; the specific process is: When the sun is at the vernal equinox, for any incident angle β of sunlight, calculate the light pressure projection area of ​​each sun-receiving surface of the rectangular structure; the expression is: S WX =l×w×|sinβ|+w×h×|cosβ| (6) In the formula, S WX It represents the light pressure projection area of ​​the star unit in the illumination direction.

8. The method for calculating the optical pressure perturbation of antenna satellites according to claim 7, characterized in that: In the step 24, the antenna unit is regarded as a parabolic structure to obtain the light pressure projection area of ​​the antenna unit in the illumination direction; The specific process is: 1) The antenna surface is modeled as a parabola, which is a special form of a quadratic surface. The general equation of a quadratic surface is: <h2 style=";text-align:left;direction:ltr">a0+a1x”+a2y”+a3z”+a4x”y”+a5y”z”+a6z”x”+a7x”<h2 style=";text-align:left;direction:ltr"> 2 <h2 style=";text-align:left;direction:ltr"> +a8y”<h2 style=";text-align:left;direction:ltr"> 2 <h2 style=";text-align:left;direction:ltr"> +a9z”<h2 style=";text-align:left;direction:ltr"> 2 <h2 style=";text-align:left;direction:ltr"> =0 (7) In the formula, a0, a1, a2, a3, a4, a5, a6, a7, a8, and a9 represent parameters; x', y', and z' represent the coordinate values ​​of the three axes in the Cartesian coordinate system; The Cartesian coordinate system x"y"z" is: Taking the center of the antenna canopy as the origin, the speed direction of the antenna satellite is the x-axis, the direction pointing to the center of the earth is the z-axis, and the y-axis is determined by the right-hand rule; Finite element software is used to fit the general equation of the quadratic surface, and the parabolic equation of the antenna is obtained as follows: 0=ax” 2 +month” 2 -z” (8) Where a represents a parameter; N is the total height of the antenna structure; z”≤N; For the form 0 = ax" 2 +ay” 2 -z" antenna unit, assuming the direction of the sun's rays is Find the light pressure projection area of ​​the antenna unit in the direction of illumination, that is, find the area of ​​the parabola with the normal vector The projected area of ​​the plane y”cosθ+z”sinθ=0; In the formula, θ represents the direction of sunlight; Parabola equation 0 = ax" 2 +ay” 2 The surface gradient of "-z" is 2) When When , the light pressure projection area of ​​the antenna unit in the illumination direction is calculated based on the light pressure projection area formula; 3) When θ is equivalent to π-θ, ​​substitute it into the light pressure projection area formula to calculate the light pressure projection area of ​​the antenna unit in the illumination direction; where θ is given by get; 4) When θ is equivalent to θ-π, substitute it into the light pressure projection area formula to calculate the light pressure projection area of ​​the antenna unit in the illumination direction; where θ is given by get; 5) When θ is equivalent to 2π-θ, substitute it into the light pressure projection area formula to calculate the light pressure projection area of ​​the antenna unit in the illumination direction; where θ is given by get.

9. The method for calculating the optical pressure perturbation of antenna satellites according to claim 8, characterized in that: In the above 2) When , the light pressure projection area of ​​the antenna unit in the illumination direction is calculated based on the light pressure projection area formula; the specific process is: 21) When When , since the parabola surface is continuous, the contour curve of the projected area is obtained by projecting the gradient and the normal vector perpendicularly. In the formula, Represents the parabola equation 0 = ax" 2 +ay” 2 - surface gradient of z”, z”≤N; represents a normal vector; represents a dot product; The tangent plane is Substituting (10) into the parabola equation 0 = ax" 2 +ay” 2 -z", get the curve Replace (x”, y”, z”) with Substituting into equation (11), eliminating k, we get the curve The equation of the projection contour curve is 22) The parabola boundary Γ is 23) When When , the contour curve of the projection of the antenna unit is the projection of the boundary Γ on the projection plane y”cosθ+z”sinθ=0; When θ = 0, the curve Go to Curve Find the critical angle θ0; the specific process is: when and the parabola equation 0 = ax" 2 +ay” 2 -z” has no intersection point, the projected contour curve is the projection of the boundary Γ on y”cosθ+z”sinθ=0, and the critical value Critical Angle when When , the light ray has no tangent point with the parabola, then the projected contour curve is the projection of the boundary Γ on y”cosθ+z”sinθ=0; When 0≤θ≤θ0, the projection is a curve With curve Projection on y”cosθ+z”sinθ=0; 24) Find the coordinates of point O and point P; the specific process is: Simultaneous equations Get the coordinates of point O and point P coordinates 25) Based on the coordinates of point O and point P in 24), calculate The projected area S1 of the curve formed with the straight line OP is expressed as: In the formula, t1 represents the intermediate variable; · represents the dot product; 26) Based on the coordinates of point O and point P in 24), for The projected area S2 of the curve formed with the straight line OP is expressed as: curve With curve The projected area of ​​the enclosed curve at ycosθ+zsinθ=0 is S=S1+S2 (18) 27) Based on 23), 25) and 26) and the light pressure projection area formula (19), calculate the light pressure projection area of ​​the antenna unit in the illumination direction; the expression is 10. The method for calculating the optical pressure perturbation of antenna satellites according to claim 9, characterized in that: In step 3, the light pressure perturbation force received by each part is calculated according to the exposure factor of step 1 and the light pressure projection area of ​​step 2. The specific process is: If k = 0, the pressure force of sunlight is 0; If k = 1, the solar pressure force is not 0, and the solar pressure force F is calculated. SPP , the calculation formula is: Among them, C R is the satellite surface reflection coefficient; ρ SR is the solar light pressure constant; k is the insolation factor; AU is 1 astronomical unit, 1AU=1.495978707×10 11 ;|r sun | is the modulus of the sun's position vector in the ECI coordinate system; S Ref is the satellite's light pressure projection area; e D is the unit vector pointing from the satellite to the sun.

Images