Activation array element rapid calculation method based on trajectory included angle solution

Abstract

The invention relates to the technical field of radio waves, in particular to an activated array element rapid calculation method based on trajectory included angle solution. Comprising the steps of performing slice preprocessing on a visible arc section of a future moving trajectory of a satellite; discrete time sampling is carried out according to the sliced visible arc segments; calculating azimuth pitch angle information of connection between the satellite and the ground measurement and control station at the sampling moment; converting the azimuth pitching information into three-dimensional rectangular coordinate information; according to the three-dimensional rectangular coordinate information, fitting the relationship between the satellite coordinates in the visible arc section of the slice and the time t; for any array element, calculating the minimum included angle of the rotation track of the array element in the visible arc section of the slice; and comparing the minimum included angle with a threshold value to finish judgment. According to the method, the precision equivalent to that of a discrete sampling method is realized with relatively low calculation complexity which is close to constant or linearly increased, and the complete activated array plane in the whole visible arc section is directly output, so that the calculation complexity is remarkably reduced.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of radio waves, in particular to a fast calculation method of active elements based on trajectory angle solution. BACKGROUND

[0002] In the future-oriented field of spaceflight TT&C (Tracking, Telemetry and Command), especially in the service of the increasingly large LEO (Low Earth Orbit) satellite constellation system, the full-sky phased array technology shows incomparable advantages. A single ground TT&C station equipped with a full-sky phased array can continuously and agilely track and control multiple satellites (or other spacecraft) flying in its visible range (limited by the curvature of the earth). This capability is a key cornerstone for large satellite constellations with high frequency and multi-target collaborative operation to improve TT&C efficiency, ensure operation safety, optimize ground station layout, and reduce system construction and operation cost.

[0003] One of the core performances of this technology is the dynamic adjustment of the "active array surface". The traditional method is to determine the active state based on the angle between the spacecraft and the ground TT&C station connecting line (line-of-sight vector) and the array normal vector. This method is conceptually clear and direct in calculation. It is feasible and efficient for calculating the active array surface required by a single satellite at a specific time (i.e., the subset of array elements required for instantaneous beam pointing). However, space missions often require the TT&C system to have the ability to continuously serve throughout the visible arc segment of the satellite. For the planning of the "active array surface" within the visible arc segment, all possible working array elements within the entire time period need to be covered. Currently, the discrete time sampling method is widely used in engineering practice: the visible arc segment of the satellite is divided into a series of discrete time points, and the corresponding active array surface is calculated independently at each sampling time point. Finally, the union of the obtained active array surfaces is taken as the final active array surface estimate covering the entire visible arc segment.

[0004] Obviously, although this discrete sampling method is logically simple, its calculation accuracy completely depends on the sampling density (time step). In order to capture the subtle changes of the active array elements during the rapid movement of the satellite and the change of the beam pointing (especially at perigee or high dynamic maneuver), a smaller sampling interval needs to be set in theory to ensure that the final union can cover all possible active array elements within the entire visible arc segment without omission. However, the reality is contradictory. High sampling accuracy will lead to a dramatic increase in the number of discrete points, followed by an explosive increase in the number of instantaneous active array surface calculations. For advanced phased array arrays composed of thousands or even tens of thousands of elements, there is a certain computational load in calculating the instantaneous active array surface. When multiplied by a large number of sampling points, the overall computational complexity will become unacceptable, which cannot meet the real-time planning or large-scale constellation automatic scheduling requirements. The theoretical calculation complexity is O(MN), where M represents the number of sampling points and N represents the number of phased array elements. SUMMARY

[0005] The application discloses a method for calculating active array elements based on trajectory included angle solving, and has the following specific method:

[0006] The visible arc segment of the future orbit of the satellite is sliced and preprocessed;

[0007] The sliced visible arc segment is discretely time-sampled;

[0008] The azimuth and elevation angle information of the satellite connected with the ground tracking station at the sampling time is calculated;

[0009] The azimuth and elevation information is converted into three-dimensional rectangular coordinate information;

[0010] The relationship between the satellite coordinates in the sliced visible arc segment and the time t is fitted according to the three-dimensional rectangular coordinate information;

[0011] For any array element, the minimum included angle of the rotating trajectory in the sliced visible arc segment is calculated;

[0012] The minimum included angle is compared with a threshold value, if the minimum included angle is smaller than the threshold value, the array element is in an active state, otherwise the array element is in an inactive state.

[0013] Further, the time slice length of the slicing preprocessing is determined, and the specific method is as follows:

[0014] The solving of the appropriate time slice length is modeled as a multi-objective optimization problem, and the specific formula is as follows:

[0015]

[0016] s.t.t th ∈T

[0017] In the formula, P s (e i ,e j ,t th ) represents the conflict calculation accuracy rate of the i th task and the j th task when the time slice length is t th ; T = [0:t s :T tot ] represents a series of candidate values of the time slice obtained with t s as the step length, T tot represents the duration of the visible arc segment; and C(t th ) represents the calculation complexity function corresponding to the conflict between the calculation tasks under the condition of the time slice with the length of t th .

[0018] Through multi-objective weighting technology, the following single-objective optimization problem is converted:

[0019]

[0020] s.t.t th ∈T

[0021] w1+w2=1

[0022] In the formula, w1, w2 are weight coefficients for balancing the two objectives;

[0023] By setting the weight ratio of the two objectives, the optimal solution of different quality, i.e. time slice length, is obtained by using reinforcement learning technology.

[0024] Further, when performing discrete time sampling according to the sliced visible arc segment, if the satellite flight trajectory is stable, the number of sampling points is reduced, and if the satellite flight is in a high dynamic maneuvering state, the number of sampling points is increased.

[0025] The sampling time is recorded as [0, t1, …, t max ].

[0026] Further, the azimuth and pitch information is converted into three-dimensional rectangular coordinate information, and the specific method is as follows:

[0027] The azimuth and pitch angle information of the satellite and the ground tracking station connecting line at each sampling time is calculated and saved in an angle matrix, recorded as

[0028] The corresponding azimuth and pitch information in the angle array is converted into three-dimensional rectangular coordinates, and the specific conversion formula is:

[0029]

[0030] Further, according to the three-dimensional rectangular coordinate information, the relationship between the satellite coordinates in the sliced visible arc segment and the time t is fitted, and the specific method is as follows:

[0031] Using t i ,1≤i≤N, Fitting the relationship between the x-axis coordinate of the satellite in the sliced visible arc segment and the time t, the corresponding function relationship is obtained, and the specific formula is as follows:

[0032] x=f x (t)

[0033] y=f y (t)

[0034] z=f z (t)

[0035] The cubic polynomial is used to f x (t), f y (t), f z(t) is performed, expressed as:

[0036] d(t) = [f x (t), f y (t), f z (t)].

[0038] Further, for any array element, the minimum included angle of its rotating track in the slice visible arc segment is calculated, and the specific method is as follows:

[0039] Suppose its unitized coordinates in three-dimensional space are The included angle of its rotating track in the visible arc segment can be expressed as

[0040] cos(θ) = e·d(t) T

[0041] Since the cosine function is a decreasing function between 0-180 degrees, the minimum included angle corresponds to the maximum value in the above formula, so solving the minimum included angle of the array element and the rotating track is equivalent to solving the following optimization problem:

[0042] min f(t) = -e·d(t) T

[0043] s.t.0 ≤ t ≤ t max

[0044] The objective function is a nonlinear function and is a one-dimensional function, and is solved by using an optimization algorithm, and the optimal solution is denoted as t * , and the minimum included angle calculation formula is:

[0045] θ * = arccos(-f(t * )).

[0047] Due to the adoption of the above technical solutions, the present application has the following beneficial effects:

[0048] 1. The present application realizes the accuracy comparable to the discrete sampling method with lower calculation complexity close to constant or linear growth, directly outputs the complete active array in the entire visible arc segment, and significantly reduces the calculation complexity.

[0049] 2. The present application reduces the invalid activation time of the array element and improves the array element utilization efficiency by optimizing the length of the slice preprocessing while ensuring the calculation complexity.

[0050] Additional advantages, objects, and features of the application will be apparent to those skilled in the art upon examination of the following specification. It is intended to BRIEF DESCRIPTION OF DRAWINGS

[0051] The drawings accompanying the specification are included as part of, and are to be taken in conjunction with, the specification. The drawings illustrate various

[0052] Figure 1 No sa4111 satellite visible arc segment 50 seconds active array surface calculation;

[0053] Figure 2 No sa4111 satellite visible arc segment 100 seconds active array surface calculation;

[0054] Figure 3 No sa4111 satellite visible arc segment 200 seconds active array surface calculation;

[0055] Figure 4 No se108 satellite visible arc segment 200 seconds active array surface calculation;

[0056] Figure 5 The whole process is shown in the schematic diagram. DETAILED DESCRIPTION

[0057] The application will be further described below in conjunction with the drawings and examples.

[0058] A method for calculating active array elements based on trajectory angle solution, as shown in Figure 5 The specific steps are as follows:

[0059] S1, slice preprocessing is performed on the visible arc segment of the future orbit of the satellite.

[0060] In step S1, in order to reduce the idle time of the phased array elements in the satellite visible arc segment as much as possible, i.e., the difference between the reserved time and the actual used time, the time slice length needs to be set to a suitable value before the calculation of the active array elements. If the time slice is too long, the satellite task processing granularity is large, there are more redundant array elements, and the false alarm rate of the active array surface conflict detection between satellite tasks is higher. If the time slice is too small, although the active array surface conflict detection between satellite tasks is more accurate, the calculation complexity will increase rapidly. In order to solve the appropriate time slice length, the following optimization steps are adopted:

[0061] S11, modeling the solution of the appropriate time slice length as the following multi-objective optimization problem, and the multi-objective optimization model formula is as follows:

[0062]

[0063] s.t.t th ∈T

[0064] where P s (e i ,e j ,t th ) represents the conflict calculation accuracy rate of the ith task and the jth task at the time slice length t th . T = [0:t s :T tot ] represents a series of candidate values of time slices obtained in steps of t s , and T tot represents the visible arc segment duration. C(t th ) represents the calculation complexity function corresponding to the calculation conflict between tasks under the condition of a time slice of t th length.

[0065] S12, the above problem can be further converted into the following single-objective optimization problem through a multi-objective weighting technique. By setting the weight proportion of the two objectives, an optimal solution of different quality can be obtained by using reinforcement learning technology, and finally an ideal time slice length is selected. The formula of the single-objective optimization model is as follows:

[0066]

[0067] s.t.t th ∈T

[0068] w1+w2=1

[0069] In the formula, w1 and w2 are weight coefficients for balancing the two objectives.

[0070] S2, discrete time sampling is performed according to the sliced visible arc segment.

[0071] In step S2, according to the length of the time slice duration to be calculated, moderate discrete time sampling is performed in [0, t max ]. If the satellite flight trajectory is relatively stable, a small number of sampling points can be used, for example, 5. If the satellite performs high dynamic maneuvering, a larger number of sampling points is required, for example, 20. N sampling time instants [0, t1, …, t max ] are obtained. The start time and the end time of the visible arc segment are generally contained in the sampling points.

[0072] S3, the azimuth-elevation angle information of the satellite connected with the ground TT&C station at the sampling time instants is calculated and saved in an angle matrix,

[0073] S4, transform the azimuth-elevation information in the angle array obtained in step S3 into three-dimensional rectangular coordinate information, and the transformation formula is

[0074] S5, use t i ,1≤i≤N, The relationship between the satellite x-axis coordinate and time t in the fitting visible arc segment is obtained, and the corresponding function relationship is x=f x (t). Similarly, the function relationships of y-axis and z-axis coordinates with respect to time are y=f y (t) and z=f z (t). A cubic polynomial is used to fit f x (t), f y (t), and f z (t) with high precision. Therefore, the trajectory in the visible arc segment can be parameterized as d(t)=[f x (t), f y (t), f z (t)].

[0075] S6, for any array element, assuming its unitized coordinate in three-dimensional space is The angle between the array element and the rotating trajectory of the satellite in the visible arc segment can be expressed as:

[0076] cos(θ)=e·d(t) T

[0077] Since the cosine function is a decreasing function between 0 and 180 degrees, the minimum angle corresponds to the maximum value in the above formula. Therefore, solving the angle between the array element and the closest point of the rotating trajectory can be equivalent to solving the following optimization problem:

[0078] min f(t)=-e·d(t) T

[0079] s.t.0≤t≤t max

[0080] The objective function of the above problem is a nonlinear function, but since it is a one-dimensional function, it can be efficiently solved using existing optimization algorithms. In this embodiment, the fminbnb function of the Matlab optimization toolbox is used.

[0081] The optimal solution of the above optimization problem is denoted as t * , and the minimum angle can be calculated by the following formula, i.e., θ * =arccos(-f(t * )).

[0082] S7, if θ * ≤θ thIf the condition is met, the array element is in an active state; otherwise, the array element is in an inactive state.

[0083] It is important to note that although discrete-time sampling is used in this embodiment, the sampling point information is only used for satellite trajectory fitting within the time slice. Each array element actually participates in step S7 only once, which significantly reduces computational complexity compared to traditional activation array calculation methods based on discrete-time sampling.

[0084] Simulation experiment:

[0085] To verify the effectiveness of the fast calculation method proposed in this invention, the following content provides a specific example of phased array activation element calculation in a satellite telemetry and control scenario.

[0086] The STK software was used to generate GPS satellite orbit information and Starlink low-Earth orbit satellite orbit information. Based on the generated satellite information, a set of visible arc segments and a set of satellite missions were created. The set of visible arc segments includes: satellite number, visible arc segment number, start time, end time, and duration of the visible arc segment. The set of missions includes: mission number, satellite ID, mission execution duration, and mission execution status. The scheduling cycle starts at 2024-03-46:0:0 and ends at 2024-03-46:20:0. The tracking and control station is located at longitude 103°, latitude 31°. A full-space phased array is used. Figure 1 The array shown is a combination of spherical and cylindrical surfaces, with 26469 array elements. For the Starlink low-Earth orbit satellite mission set, the phased array activation angle is set to 20 degrees. To statistically analyze the differences between the proposed algorithm and traditional algorithms in calculating the activated array elements, a binary vector a of length 26469 is first defined. st a st [i] = 1 indicates that the i-th element is active. The error ratio is defined as follows:

[0087]

[0088] This embodiment takes the trajectory fitting of a 100-second visible arc segment of satellite number "se108" as an example. The specific steps are as follows:

[0089] Firstly, a time slice is obtained as the time processing range for solving the activated elements, and the satellite visible arc segment in this embodiment is 10 minutes, and the time slice size is determined to be 100 seconds. In the calculation process, the proportion used is w1=0.7, w2=0.3. The reinforcement learning method used is deep Q-learning, DQN algorithm. The time slice candidate set is T=[0:2:600], and the final calculated time slice length is 100 seconds. For a certain time slice, discrete time sampling is performed every 10, and the azimuth and elevation angle information is obtained and stored in the angle array. The contents of the corresponding angle array are:

[0090]

[0091]

[0092] Convert to rectangular coordinates and normalize to obtain the corresponding rectangular coordinate values, which are:

[0093]

[0094] The formula fitting of the x, y, and z axis coordinates can be obtained as

[0095] f x (t)=0.001t-0.845

[0096] f y (t)=0.001t+0.527

[0097] f z (t)=0.001t+0.098

[0098] For a certain element, the rectangular coordinates are e=[0 100 100], and the solution of the corresponding optimization problem can obtain the cosine value of the nearest point on the trajectory to the point as 0.6294, which is less than the cosine value corresponding to 20 degrees, so the point is determined as an inactive state.

[0099] Figures 1 to 3 The phased array activation array calculation of the low-orbit satellite numbered sa411 at 50 seconds, 100 seconds and 200 seconds is given, respectively, and the visible arc segment starts at 2024-3-46:9:20. From the figure, it can be seen that the algorithm can accurately calculate the phased array activation array set in different time periods. With the increase of the time period length, the error ratio performance of the calculated activated element set gradually decreases, but in the test range, the error ratio performance is still very good. The error ratio in 50 seconds, 100 seconds and 300 seconds is 0, 0 and 0.0112, respectively.

[0100] In order to verify the effectiveness of the algorithm, the application further simulates the phased array active element calculation of the satellite numbered se 108 within 200 seconds, and the visible arc segment starts at 2024-3-46:12:40. As shown in FIG. 8, within the 200 second time period, the proposed algorithm can still accurately calculate the active element set of the phased array within the time period, and the corresponding error ratio is 0.0179, which still remains at a very low level. Figure 4

[0101] Finally, it should be noted that: the above examples are only used to illustrate the technical solutions of the present application, but not to limit it, although the present application has been described in detail with reference to the above examples, those skilled in the art should understand that: the specific embodiments of the present application can still be modified or replaced by the equivalent, without departing from the spirit and scope of the present application, any modification or equivalent replacement, which should be covered within the protection scope of the claims of the present application.​

Claims

1. A method for calculating an active element based on the trajectory angle solution, characterized in that, The specific method is as follows: The visible arc segment of the future orbit of the satellite is sliced and pretreated; The sliced visible arc segment is sampled at discrete time; The azimuth-elevation angle information of the satellite connected with the ground tracking station at the sampling time is calculated; The azimuth-elevation information is converted into three-dimensional rectangular coordinate information; The relationship between the satellite coordinates in the sliced visible arc segment and the time t is fitted according to the three-dimensional rectangular coordinate information; For any array element, the minimum included angle of the rotating track in the sliced visible arc segment is calculated; The minimum included angle is compared with the threshold value, if the minimum included angle is less than the threshold value, the array element is in the active state, otherwise the array element is in the inactive state.

2. The method of claim 1, wherein the method is based on a trajectory angle. The time slice length of the slicing pretreatment is determined, and the specific method is as follows: The solving of the appropriate time slice length is modeled as a multi-objective optimization problem, and the specific formula is: s.t.t th ∈T In the formula, P s (e i , e j , t th ) represents the conflict calculation accuracy rate of the ith task and the jth task when the time slice length is t th ; T = [0:t s :T tot ] represents a series of candidate values of time slices obtained in steps of t s ; T tot represents the visible arc segment duration; C(t th ) represents the calculation complexity function corresponding to the calculation task conflict under the condition of a time slice with a length of t th . Through multi-objective weighting technology, it is converted into the following single-objective optimization problem: s.t.t th ∈T w1+w2=1 In the formula, w1, w2 are weight coefficients for balancing the two objectives; By setting the weight proportion of the two objectives, the optimal solution of different quality, i.e. the time slice length, is obtained by using reinforcement learning technology.

3. The method of claim 1, wherein the method is based on a trajectory angle. When the sliced visible arc segment is sampled at discrete time, if the satellite flight trajectory is stable, the number of sampling points is reduced, and if the satellite flight is in a high dynamic state, the number of sampling points is increased; The sampling time instants are recorded as [0, t1,..., t max ].

4. The method of claim 1, wherein the method is based on a trajectory angle. The azimuth-elevation information is converted into three-dimensional rectangular coordinate information, and the specific method is as follows: The azimuth-elevation angle information of the satellite and the ground TT&C station connection line at each sampling time is calculated and saved in an angle matrix, recorded as The corresponding azimuth-elevation information in the angle array is converted into three-dimensional rectangular coordinates, and the specific conversion formula is:

5. The method of claim 4, wherein the method is based on the trajectory angle. The relationship between the satellite coordinates in the sliced visible arc segment and the time t is fitted according to the three-dimensional rectangular coordinate information, and the specific method is as follows: Use t i ,1≤i≤N, The relationship between the satellite x-axis coordinate and time t in the fitting slice visible arc segment is obtained, and the corresponding function relationship is as follows: x = f x (t) y = f y (t) z = f z (t) A cubic polynomial is used to fit f x (t), f y (t), f z (t) with high accuracy, expressed as: d(t) = [f x (t),f y (t),f z (t)].

6. The method of claim 5, wherein the angle between the trajectory and the activation element is calculated by: ###0001### where θ is the angle between the trajectory and the activation element, and φ is the angle between the trajectory and the x-axis. For any array element, the minimum included angle of the rotating track in the sliced visible arc segment is calculated, and the specific method is as follows: Assume its unitized coordinates in three-dimensional space are The angle between it and the satellite's rotating track in the visible arc segment can be expressed as cos(θ) = e - d(t) T Since the cosine function is a decreasing function between 0-180 degrees, the minimum included angle corresponds to the maximum value in the formula, so solving the minimum included angle of the array element and the rotating track is equivalent to solving the following optimization problem: min f(t) = -e-d(t) T s.t. 0≤t≤t max The objective function is a nonlinear function and only a monomial function, and is solved by an optimization algorithm, and the optimal solution is denoted as t * The minimum angle calculation formula is: θ * = arccos(-f(t * )).

Images