# Judgement Theorems and an Approach for Solving the Constellation-to-Ground Coverage Problem

Authors: Zhiming Song; Xiangyun Hu; Maocai Wang; Guangming Dai
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

The satellite constellation-to-ground coverage problem is a basic and important problem in satellite applications. A group of judgement theorems is given, and a novel approach based on these judgement theorems for judging whether a constellation can offer complete single or multiple coverage of a ground region is proposed. From the point view of mathematics, the constellation-to-ground coverage problem can be regarded as a problem entailing the intersection of spherical regions. Four judgement theorems that can translate the coverage problem into a judgement about the state of a group of ground points are proposed, thus allowing the problem to be efficiently solved. Single- and multiple-coverage problems are simulated, and the results show that this approach is correct and effective.


---

## Body

## 1. Introduction

Currently, satellite constellations technology are widely used in many practical applications, including communication, navigation, remote sensing, and science missions [1, 2], and these applications provide very important support to many fields, such as the observation of the land, sea and air [3, 4], disaster monitoring [5, 6], and resource exploration and development [7]. To meet the requirements of these applications, some typical problems encountered when designing satellite constellations have been studied by scholars. These problems include constellation design [8, 9], constellation performance evaluation [10], satellite-ground communication [11], and satellite scheduling [12, 13]. In these typical problems, the service of the constellation to the ground region, which depends on the visibility between the constellation satellite and the ground region, is a basic and important problem. This problem is called the constellation-to-ground coverage problem.Solutions to the constellation-to-ground coverage problem aim to calculate the coverage region and other information about the ground regions that are covered by satellite constellations during a time period or at a point in time. This problem is essentially a problem involving Boolean operations on 2D graphs [14, 15], which is commonly encountered in the fields of computer graphics and geology. However, compared with the traditional problems involving the Boolean operations on 2D graphs, the constellation-to-ground coverage problem has some differences because of the following three aspects: (1) the background manifold is a spherical or ellipsoidal surface in three dimensions rather than a 2D Euclidean space, (2) the shape of ground regions and the view field of the sensor can be arbitrary, and (3) the temporal characteristics during long simulation periods frequently arise in calculating the coverage region of sensor. These three aspects make this problem more difficult than traditional Boolean operation problems.The solution of the coverage problem depends on the configuration of the constellation and the coverage problem itself. Different constellation configurations, including the star constellation [16], the Walker-Delta constellation [17, 18], the rosette constellation [19], and the flower constellation [20], may have some particular configuration characteristics, such as spatial symmetry [21, 22], common-track [20], regressive orbits characteristic [23], and polar orbits [24]. Proper usage of these characteristics can greatly increase computational efficiency. From the aspect of the coverage problem, owing to the different requirements of the time and space factors for different practical applications, constellation-to-ground coverage problems are different. For navigation and communication applications, the timeliness requirement is critical; however, for some remote sensing applications such as resource exploration, the timeliness requirement is minimal. Thus, the concepts of continuous coverage problems [8, 25] and discontinuous coverage problems [24, 26] have been proposed by many scholars, and thus, different constellation-to-ground coverage problems have different characteristics and different strategies for achieving a solution.Until now, to solve constellation-to-ground coverage problems, there have been three typical methods, including the grid point method [27], the triangle mesh subdivision method [17], and the street-of-coverage method [28]. The triangle mesh subdivision method is a method that is favored by many scholars [18, 19]. This method can quickly judge whether the constellation can cover the globe at a given time. However, more detailed coverage information, such as the coverage rate, cannot be obtained. The street-of-coverage method is a method for solving the continuous coverage problem and can obtain the approximate analytic relationship between coverage and the constellation configuration so as to guide constellation design. The disadvantage of this approach is that the result is an approximate value and the method ignores phase information.Additionally, some mathematical strategies have been applied in constellation-to-ground coverage problems, such as probability statistics [29], iterative methods [30], the differential quadrature method [31], the method based on the transformation group [21, 22], the polyhedron surrounding method [32], the method based on auxiliary space [26, 33], and the numerical fitting method [34]. These mathematical strategies can be used in some of the three aforementioned methods, which can accelerate computing efficiency. However, most of these mathematical strategies have some restrictions in that they can only work in some certain conditions and can only solve some certain problems.In this paper, we focus on a subproblem of the constellation-to-ground coverage problem, the complete coverage problem, which judges whether a constellation can cover a ground region completely using one or multiple satellites. If the constellation cannot cover a ground region completely, we do not care about the coverage rate or other coverage properties.In this paper, we will discuss this problem from the aspect of mathematics. From the point view of mathematics, the instantaneous coverage region of the ground region from a satellite can be seen as a spherical region. Thus, the constellation-to-ground coverage problem can be translated into a problem on the intersections of spherical regions. According to the basic characteristics of the coverage problem, a group of definitions are given, and a group of theorems as well as proofs are put forward. Then, according to the judgement theorems, highly accurate, efficient, and novel approaches for solving this problem are proposed.The paper is organized as follows. Section2 gives a brief review of the traditional methods. Section 3 introduces the basic concepts and appointments of this paper. Section 4 proposes a group of judgement theorems for this coverage problem. Section 5 discusses the detailed procedures of the approaches for solving this coverage problem. In Section 6, a group of experiments, including single-coverage problems and multiple-coverage problems, are conducted. Section 7 concludes the paper.


## 2. Background

The traditional methods for solving the complete coverage problem include the grid point method and the triangle mesh subdivision method. A brief introduction to these two methods will be given.


### 2.1. The Process of the Grid Point Method

The grid point method is one of the most commonly used approaches for solving the constellation-to-ground coverage problem. The main idea of this approach is discretization, which is the division of a continuous object into a group of discrete objects and then calculating the characteristics of these discrete objects to express the characteristics of the coverage problem.For a coverage problem, two types of objects need to be discretized: the first is space, and the other is time. The discretization of space means discretizing the ground region into a group of grid points according to a regional discretization strategy, and the discretization of time means discretizing a time range into a set of time points.For the grid point method, the process is as shown in Figure1, and the basic steps are as follows:Figure 1
The basic process of traditional GPA.Step 1.
For ground regionR, a minimum longitude and latitude box that can completely contain R is found (the longitude and latitude box corresponds to the blue rectangle in the second panel of Figure 1).Step 2.
Given a space partition precisionε and the chosen sampling strategy, a group of feature points are chosen and are called grid characteristic points, which are denoted by {Nk}, where k=1,2,…,nR. {Nk} corresponds to the red and greed dots in the second panel of Figure 1.Step 3.
For every grid characteristic point belonging to{Nk}, whether this point belongs to the ground region R is judged. Then, all grid characteristic points that belong to R are gathered, and these points are written as {N′k}. {N′k} corresponds to the green dots in the second panel of Figure 1, while the red dots represent the grid characteristic points that do not belong to R.Step 4.
For every grid characteristic point belonging to{N′k}, whether this point can be covered by the constellation C is judged. Then, all grid characteristic points that belong to R are gathered, and these points are written as {N~′k}. {N~′k} corresponds to the black dots in the third panel of Figure 1, while the red circles stand for the coverage region of constellation C.Step 5.
The number of points in{N~′k} and {N′k} are counted, and the ratio of these two sets is regarded as the coverage rate.This process of the grid point method does not address the temporal characterization. The temporal characterization is hidden in Step5, the process of judging whether a grid characteristic point can be covered by the constellation. The basic approach for solving the temporal characterization also involves discretization. For the grid point method, the basic steps of discretizing the temporal elements are as follows:Step 1.
Given a time partition precision corresponding to time intervalΔt, a group of time points t0,t1,…,tn can be obtained from the time interval [tS,tE], where tS=t0, tE=tn, and tk-tk-1=Δt for any k=1,2,…,n.Step 2.
For any time pointtk, the coverage state of every grid characteristic point belonging to {N′k} is checked for whether each point can be covered by the constellation C at the current time.Step 3.
The coverage characteristics of all time points are analyzed, and statistics are calculated for these time points. Then, the temporal characteristics of the constellation-to-ground region can be obtained.The greatest advantage of the grid point method is that it has a wide response domain, which means that it can solve nearly all coverage problems. However, when high computational accuracy is required, this approach is highly temporally and spatially complexity. The temporal and spatial complexity grows exponentially as the required calculation accuracy increases. This feature is one of the biggest drawbacks of the grid point method and greatly restricts the practical application of this approach.Additionally, the grid point method is good at computing the value of the coverage rate, but when solving complete coverage problems, the grid point method is not efficient as it must compute a great amount of information that will not be used when solving complete coverage problems.


### 2.2. The Spherical Triangle Mesh Subdivision Method

The main idea of the spherical triangle mesh subdivision method is to divide the globe into a set of spherical Delaunay triangles, where the vertexes of the set of spherical Delaunay triangles are the points beneath the satellites in the constellation, and then checking whether the constellation can completely cover a spherical Delaunay triangle given the radius of the spherical circumcircle of the spherical Delaunay triangle.Spherical Delaunay triangulation is an important pretreatment technique in both graphical and numerical analyses. By following certain rules, the globe can be divided into adjacent and mutually disjoint triangles that have the following characteristics:(i)
The minimum angle of the spherical triangles is maximized. This characteristic can make the subdivision result more uniform and regular.(ii)
The triangles have empty interiors. The characteristic of empty interiors is that the interior of every spherical triangle circumcircle does not contain other vertexes.The subdivision result is generally unique. ForN vertexes, the globe can be divided into 2N-4 spherical triangles, where each has 3N-6 edges.The basic process of the spherical triangle mesh subdivision method is as follows:Step 1.
Calculate the instantaneous location of the subsatellite point of all the satellites in the constellation.Step 2.
Take the instantaneous subsatellite position of all the satellites in the constellation as the vertex set and then construct a spherical Delaunay triangle subdivision of the globe.Step 3.
For any Delaunay triangle, calculate the radius of the spherical circumcircle and obtain the maximum value, which is denoted byDmax.Step 4.
Denote the instantaneous coverage radius on the earth of the satellite byrS. If Dmax>rS, this indicates that the constellation cannot fully cover the global. If Dmax<rS, then the constellation can completely cover the globe.The spherical triangle mesh subdivision method is a method for solving complete coverage problems. However, the target of complete coverage problems must be the globe. However, the method does not work well for other types of targets, such as zonal targets, polar cap targets, and local area targets. Additionally, this method can only solve problems in whichrS of every satellite is equal.


## 2.1. The Process of the Grid Point Method

The grid point method is one of the most commonly used approaches for solving the constellation-to-ground coverage problem. The main idea of this approach is discretization, which is the division of a continuous object into a group of discrete objects and then calculating the characteristics of these discrete objects to express the characteristics of the coverage problem.For a coverage problem, two types of objects need to be discretized: the first is space, and the other is time. The discretization of space means discretizing the ground region into a group of grid points according to a regional discretization strategy, and the discretization of time means discretizing a time range into a set of time points.For the grid point method, the process is as shown in Figure1, and the basic steps are as follows:Figure 1
The basic process of traditional GPA.Step 1.
For ground regionR, a minimum longitude and latitude box that can completely contain R is found (the longitude and latitude box corresponds to the blue rectangle in the second panel of Figure 1).Step 2.
Given a space partition precisionε and the chosen sampling strategy, a group of feature points are chosen and are called grid characteristic points, which are denoted by {Nk}, where k=1,2,…,nR. {Nk} corresponds to the red and greed dots in the second panel of Figure 1.Step 3.
For every grid characteristic point belonging to{Nk}, whether this point belongs to the ground region R is judged. Then, all grid characteristic points that belong to R are gathered, and these points are written as {N′k}. {N′k} corresponds to the green dots in the second panel of Figure 1, while the red dots represent the grid characteristic points that do not belong to R.Step 4.
For every grid characteristic point belonging to{N′k}, whether this point can be covered by the constellation C is judged. Then, all grid characteristic points that belong to R are gathered, and these points are written as {N~′k}. {N~′k} corresponds to the black dots in the third panel of Figure 1, while the red circles stand for the coverage region of constellation C.Step 5.
The number of points in{N~′k} and {N′k} are counted, and the ratio of these two sets is regarded as the coverage rate.This process of the grid point method does not address the temporal characterization. The temporal characterization is hidden in Step5, the process of judging whether a grid characteristic point can be covered by the constellation. The basic approach for solving the temporal characterization also involves discretization. For the grid point method, the basic steps of discretizing the temporal elements are as follows:Step 1.
Given a time partition precision corresponding to time intervalΔt, a group of time points t0,t1,…,tn can be obtained from the time interval [tS,tE], where tS=t0, tE=tn, and tk-tk-1=Δt for any k=1,2,…,n.Step 2.
For any time pointtk, the coverage state of every grid characteristic point belonging to {N′k} is checked for whether each point can be covered by the constellation C at the current time.Step 3.
The coverage characteristics of all time points are analyzed, and statistics are calculated for these time points. Then, the temporal characteristics of the constellation-to-ground region can be obtained.The greatest advantage of the grid point method is that it has a wide response domain, which means that it can solve nearly all coverage problems. However, when high computational accuracy is required, this approach is highly temporally and spatially complexity. The temporal and spatial complexity grows exponentially as the required calculation accuracy increases. This feature is one of the biggest drawbacks of the grid point method and greatly restricts the practical application of this approach.Additionally, the grid point method is good at computing the value of the coverage rate, but when solving complete coverage problems, the grid point method is not efficient as it must compute a great amount of information that will not be used when solving complete coverage problems.


## 2.2. The Spherical Triangle Mesh Subdivision Method

The main idea of the spherical triangle mesh subdivision method is to divide the globe into a set of spherical Delaunay triangles, where the vertexes of the set of spherical Delaunay triangles are the points beneath the satellites in the constellation, and then checking whether the constellation can completely cover a spherical Delaunay triangle given the radius of the spherical circumcircle of the spherical Delaunay triangle.Spherical Delaunay triangulation is an important pretreatment technique in both graphical and numerical analyses. By following certain rules, the globe can be divided into adjacent and mutually disjoint triangles that have the following characteristics:(i)
The minimum angle of the spherical triangles is maximized. This characteristic can make the subdivision result more uniform and regular.(ii)
The triangles have empty interiors. The characteristic of empty interiors is that the interior of every spherical triangle circumcircle does not contain other vertexes.The subdivision result is generally unique. ForN vertexes, the globe can be divided into 2N-4 spherical triangles, where each has 3N-6 edges.The basic process of the spherical triangle mesh subdivision method is as follows:Step 1.
Calculate the instantaneous location of the subsatellite point of all the satellites in the constellation.Step 2.
Take the instantaneous subsatellite position of all the satellites in the constellation as the vertex set and then construct a spherical Delaunay triangle subdivision of the globe.Step 3.
For any Delaunay triangle, calculate the radius of the spherical circumcircle and obtain the maximum value, which is denoted byDmax.Step 4.
Denote the instantaneous coverage radius on the earth of the satellite byrS. If Dmax>rS, this indicates that the constellation cannot fully cover the global. If Dmax<rS, then the constellation can completely cover the globe.The spherical triangle mesh subdivision method is a method for solving complete coverage problems. However, the target of complete coverage problems must be the globe. However, the method does not work well for other types of targets, such as zonal targets, polar cap targets, and local area targets. Additionally, this method can only solve problems in whichrS of every satellite is equal.


## 3. Basic Appointment and Definition

Appointment 1 (ground region).
The ground region is appointed as a connected spherical domain and also as a closed domain.Denote a spherical domain byR, and denote the boundary of R by ∂R. Because the spherical domain is appointed as a closed domain, then ∂R⊂R.Appointment 2 (satellite coverage region).
The satellite coverage region is appointed as a connected spherical domain as well as an open domain.LetS be a satellite. Denote the coverage region of S by ΛS, and denote the boundary of ΛS by ∂ΛS. According to Appointment 2, ΛS is an open domain; then ∂ΛS⊄ΛS, and ΛS∩∂ΛS=⌀. The union set of ΛS and ∂ΛS is called the closure of S and is denoted by Λ¯S.In actual applications, we may also encounter the situation in which the ground region or satellite coverage region is not a connected domain but a nonconnected domain. Then, this nonconnected ground region or satellite coverage region can be equivalently seen as a group of connected ground regions or satellite coverage regions, and this treatment has no influence on the results of the analysis. Thus, in this paper, for the ground region or satellite coverage region, only the connected spherical domains are considered.WhenΛS is a single connected domain, ∂ΛS consists of a connected curve. However, when ΛS is a multiple connected domain, ∂ΛS consists of more than one single connected curve. Let ∂λSk be the kth boundary curve of ∂ΛS; then ∂ΛS=⋃k=1n∂λSk, where n is the number of connected boundary curves.Because the position of satellite changes with time,ΛS is a time-varying region and has hidden time parameters t, which can be written as ΛS(t) when the time parameters are emphasized.A constellation is a set of satellites that has a relatively stable spatial configuration. Denote a constellation byC. Then, similarly, the coverage region of the constellation is written as ΛC, and the boundary of the constellation coverage region is written as ∂ΛC. Then, the following equation holds:(1)ΛC=⋃S∈CΛS.However,∂ΛC≠⋃S∈C∂ΛS because some points of ∂ΛS may be covered by other satellites belonging to C, which means that these points are no longer boundary points. Then, ∂ΛC⊂⋃S∈C∂ΛS.Generally,ΛC is not a single connected domain but a multiple connected domain. A connected domain of ΛC is called the coverage branch of ΛC. Let m be the count of the connected domain of ΛC, and denote the kth coverage branch of ΛC by ΛCk; then,(2)ΛC=⋃k=1mΛCk.Then,m is less than or equal to the number of satellites in constellation C.Appointment 3 (cover).
A ground point denoted byP is said to be covered by satellite S when P is an inner point of ΛS, which can be expressed by P∈ΛS according to Appointment 2. Additionally, a ground point P is said to be covered by constellation C when P∈ΛC.According to Appointments 2 and 3, satelliteS cannot cover a point located on ∂ΛS.Generally, in actual applications, these three concepts do not need to be defined so rigorously. However, because we will solve this problem from the point view of mathematics, the rigorous definition of the concepts can make the forms of the expression and process of the proof simpler.Then, some definitions that will be used later in this paper are given.Definition 4 (connected satellites).
Lett be a time point. For two satellites S1 and S2, if ΛS1(t) and ΛS2(t) belong to the same branch of ΛC(t), then S1 is said to be connected to S2 at time t, which is denoted by S1~S2, where t is a hidden parameter of this relationship. S1 and S2 are called connected satellites.For set{ΛS∣S∈C}, we can see that the connectivity of satellites can be regarded as an equivalence relationship. The equivalence class of satellite S using the equivalence relationship is written as [S]. For any satellite S′ belonging to [S], ΛS′ belongs to the same connected branch of ΛC.Definition 5 (connected coverage region).
The connected coverage region of satelliteS in constellation C is the union of the coverage region of satellites belonging to [S] and is denoted by US. Similarly, the boundary of US is written as ∂US.Then, evidently,∂US⊂∂ΛC, and ∂ΛC=⋃S∈C∂US.Additionally, ifΛS∩ΛS′≠⌀, then US=US′.Definition 6 (complete coverage).
ConstellationC is said to completely cover the ground region R at time t when ∀P∈R, ∃S∈C such that P∈ΛS.Strictly speaking, complete coverage is called single-satellite complete coverage. There also exists another type of complete coverage, which is called multiple-satellite complete coverage.Definition 7 (multiple-satellite complete coverage).
ConstellationC is said to possess k-multiple complete coverage of the ground region R at time t when ∀P∈R, ∃{S1,S2,…,Sk}⊂C such that P∈ΛSj, where j=1,2,…,k.


## 4. Coverage Problem Judgement Theorems

Theorem 8.
The necessary and sufficient conditions for constellationC to completely cover ground region R are as follows: (1)
∃S∈C such that ΛS∩R≠⌀.(2)
R⊂US.Proof.
For the sufficiency,R⊂US⊂ΛC; then C can completely cover R.
Because constellationC can completely cover ground region R, then condition (1) is clearly satisfied. For condition (2), supposing R⊄US, then R∩∂US≠⌀. Then, ∃P∈R such that P∈∂US⊂∂ΛC⊄ΛC, and then P cannot be covered by constellation C. This condition contradicts the condition that constellation C completely covers ground region R, so the assumption is invalid. Thus, R⊂US.
Then the theorem is proved.Then, according to Theorem8, if constellation C can completely cover region R, then constellation C has one and only one coverage brunch that can contain R. Additionally, if constellation C can completely cover region R, then, for any coverage branch ΛCk, either R⊂ΛCk or R∩ΛCk=⌀ is held.Theorem8 reveals some characteristics of the complete coverage problem, and this theorem can be used to prune the branches.Theorem 9.
The necessary and sufficient conditions for constellationC to completely cover ground region R are as follows: (1)
∃S∈C such that ΛS∩R≠⌀.(2)
R∩∂US=⌀.Proof.
For the sufficiency, if there exists a satelliteS such that ΛS∩R≠⌀, then, supposing R⊄US, because ΛS⊂US, then R∩US≠⌀. Then, the relationship between R and US may have two statuses, which are US⊂R and R⊄US, but R∩US≠⌀. If US⊂R, because US is an open domain and R is a closed domain, then ∂US⊂R, and then R∩∂US≠⌀. If R⊄US but R∩US≠⌀, then R∩∂US≠⌀. Then, the statuses both contradict condition (2). So, this assumption is invalid. Then, R⊂US. According to Theorem 8, constellation C can completely cover ground region R.
For the necessity, when constellationC can completely cover ground region R, then condition (1) is obviously satisfied. Then, according to Theorem 8, R⊂US. According to Appointment 2, it satisfies that ∂US∩US=⌀; then R∩∂US=⌀.Generally,R and US are both 2-dimensional spherical domains, and the coverage problem of C to R is an intersecting problem of two two-dimensional spherical domains. Then, Theorem 9 translates this coverage problem into an intersecting problem of a two-dimensional domain and a one-dimensional domain, which is much easier to solve than the original problem.However, the problem of intersecting a two-dimensional domain and a one-dimensional domain is also somewhat complex. We then introduce the boundary feature point and give another theorem, thus making the problem easier.Before we discuss the theorem, a definition that will be used in the theorem should be defined.Definition 10 (boundary feature points).
The boundary feature points of satelliteS are a group of special points located on ∂ΛS and are denoted by K(∂ΛS). The boundary feature point on the kth boundary curve, which is denoted by K(∂λSk), is divided into two conditions:
Condition  1. For ∂λSk, if ∂λSk∩∂R≠⌀ or ∂λSk∩∂ΛC∖S∩R≠⌀, then we select ∂λSk∩∂R∪∂λSk∩∂ΛC∖S∩R as the boundary feature point.
Condition 2. For ∂λSk, if ∂λSk⊂R and ∂λSk∩∂ΛC∖S=⌀, we select an arbitrary point on ∂λSk as the boundary feature point.
Then,K(∂ΛS)=⋃k=1nK(∂λSk).
Similarly,(3)K∂ΛC=⋃S∈CK∂ΛSK∂ΛS=⋃S′∈SK∂ΛS′.
Per the definition of a boundary feature point, we can see thatK(∂ΛS) consists of finite points when we do not consider the critical condition that ∂λSk and ∂R overlap with each other.Theorem 11.
The necessary and sufficient condition for constellationC to completely cover ground region R is as follows: (1)
∃S∈C such that ΛS∩R≠⌀.(2)
∀P∈K(∂Λ[S]), ∃S′∈[S] such that P∈ΛS′.Proof.
For the necessity, the first condition is obvious. For the second condition, per Definition10, K(∂Λ[S])⊂R; then, because constellation C can completely cover ground region R, the second condition holds.
For the sufficiency, from Figure2, we can see that the boundary of the uncovered region consists of ∂ΛC and ∂R and that the vertexes of the uncovered region are made up of the intersection of the coverage boundary of two satellites as well as the intersection of the coverage boundary of a satellite and the boundary of the region. Then, we can see that the vertexes of the uncovered region approximate the boundary feature points of the constellation. Intuitively, supposing that constellation C does not completely cover ground region R, if the uncovered region has vertexes, then the vertexes must be the boundary feature points, and if the uncovered region has no vertex, then boundary of the region is just Condition 2 of Definition 10, and an arbitrary point on the boundary is selected as the boundary feature point. Then, the boundary feature points are not covered by another satellite, which contradicts Condition 2. Thus, the assumption is invalid. Then, constellation C can completely cover ground region R.Figure 2
Two typical conditions indicating that a constellation does not completely cover region.Then, by Theorem11, the complete coverage problem can be translated into checking the coverage state of a finite number of points, which greatly reduces the computational difficulty. More importantly, this approach has no calculation error.Additionally, this theorem can not only be used for solving the single-coverage problem but also can solve the multiple-coverage problem. Thus, a new theorem that is used to deal with the multiple-coverage problem will be proposed.Theorem 12.
The necessary and sufficient condition for constellationC to possess k-multiple complete coverage of ground region R is as follows: (1)
∃S∈C such that ΛS∩R≠⌀.(2)
∀P∈K(∂Λ[S]), ∃{S1,S2,…,Sk}⊂[S] such that P∈ΛSj, where j=1,2,…,k.The proof of this theorem is similar to that of Theorem11. Thus, we do not cover this process here.


## 5. Approach for Judging the Coverage State of the Constellation-to-Ground Region

### 5.1. Approach for Judging the Coverage State at a Certain Time

In the actual computing process, generally, we use Theorems11 and 12 to solve the practical problem because the computational efficiency with this approach is much higher than when using Theorems 8 and 9.The basic steps for judging whether constellationC can singly or multiply cover ground region R completely at a certain time t are as follows:Step 1.
Check whether there exists a satelliteS such that ΛS∩R≠⌀. If S does not exist, then it can be asserted that the constellation C cannot completely cover ground region R at time t.Step 2.
Calculate the boundary feature points of constellationC to ground region R at time t, which is denoted by K(∂ΛC).Step 3.
Judge whether all points belonging toK(∂ΛC) can be covered by a satellite (or k satellites in k-multiple-coverage problem). If the result is yes, then C can completely cover R at t; else, C cannot completely cover R at t.


### 5.2. Approach for Judging the Coverage State for a Period of Time

Theorems11 and 12 can only judge the coverage state at a certain time. However, in actual engineering applications, time factors are often considered, including continuous coverage [8, 25].Definition 13 (continuous coverage).
Let[tS,tE] be a period of time. If for any t∈[tS,tE], the constellation C can completely cover ground region R at t, then we say that the constellation C can continuously cover ground region R at [tS,tE].Because the relevant judgement theorems about time factors have not been included in this paper, when dealing with the time factors, the discretization method is adopted. The basic steps for solving continuous coverage problems are as follows:Step 1.
Given a time partition precision corresponding to time intervalΔt, [tS,tE] can be divided into a group of time points t0,t1,…,tn where tS=t0, tE=tn, and tk-tk-1=Δt for any k=1,2,…,n.Step 2.
For any time pointtk, if constellation C cannot singly (or multiply) cover a ground region R completely at time tk, then C cannot continuously cover R at [tS,tE]. If C can completely cover R completely at all of the time points, then C can continuously cover R at [tS,tE].However, different from the approach for judging the coverage state at a certain time, the approach for evaluating the continuous coverage problem may have calculation errors whenΔt is not small enough. Thus, in order to control calculation errors, Δt is chosen to be a sufficiently small number.


### 5.3. Calculation of the Critical Beam Angle

Theorems11 and 12 can only obtain a result of “yes” or “no,” which is a qualitative result. However, in actual applications, a quantitative result is more useful. So, in this subsection, we will use Theorems 11 and 12 to obtain a typical relevant quantitative result by calculating the critical beam angle of the sensor.Definition 14 (the critical beam angle).
Consider an isomorphic constellationC, which means that all satellites in the constellation have the same semimajor axis, eccentricity, and orbit inclination, and all satellites carry same type sensor. Then, for ground region R, the critical sensor beam angle is the minimum beam angle of sensor with which the constellation C can completely cover the ground region R.For constellationC, if the half-beam angle of sensor equipped on the satellite that belongs to C is α, then the coverage region of C at t is written as ΛC(α,t). Then, clearly, if α<β, then ΛC(α,t)⊂ΛC(β,t); that is, the results are monotonous relative to the half-beam angle of sensor. So, the critical beam angle for C to R can be obtained by the bisection method.Denote the radius of the Earth byRE, and denote the semimajor axis of every satellite by a; then, given an elevation constraint of the sensor, the maximum half-beam angle denoted by αmax is arcsin⁡RE/a, and the minimum half-beam angle denoted by αmin is 0.Step 1.
For constellationC, calculate αmax and αmin.Step 2.
IfR⊄ΛC(αmax,t), then the algorithm terminates and the critical half-beam angle does not exist.Step 3.
Letαmid be equal to (αmax+αmin)/2. Then, check the coverage state of C with half-beam angle αmid to R at t; then, if R⊂ΛC(αmid,t), then αmax=αmid; else, αmin=αmid.Step 4.
LetΔα be the allowable error, which is given according to the requirements of the actual problem. If αmax-αmin<Δα, then the algorithm terminates, and αmid is the ultimate critical half-beam angle; else, jump to Step 3.The critical half-beam angle ofC to R at t is denoted by η(C,R,t).Then, for a time period[tS,tE], the critical half-beam angle of C to R, which is denoted by η(C,R,[tS,tE]), is defined as the maximum value of the critical half-beam angle of C to R at every moment in [tS,tE]. Then,(4)ηC,R,tS,tE=maxt∈tS,tE⁡ηC,R,t.


## 5.1. Approach for Judging the Coverage State at a Certain Time

In the actual computing process, generally, we use Theorems11 and 12 to solve the practical problem because the computational efficiency with this approach is much higher than when using Theorems 8 and 9.The basic steps for judging whether constellationC can singly or multiply cover ground region R completely at a certain time t are as follows:Step 1.
Check whether there exists a satelliteS such that ΛS∩R≠⌀. If S does not exist, then it can be asserted that the constellation C cannot completely cover ground region R at time t.Step 2.
Calculate the boundary feature points of constellationC to ground region R at time t, which is denoted by K(∂ΛC).Step 3.
Judge whether all points belonging toK(∂ΛC) can be covered by a satellite (or k satellites in k-multiple-coverage problem). If the result is yes, then C can completely cover R at t; else, C cannot completely cover R at t.


## 5.2. Approach for Judging the Coverage State for a Period of Time

Theorems11 and 12 can only judge the coverage state at a certain time. However, in actual engineering applications, time factors are often considered, including continuous coverage [8, 25].Definition 13 (continuous coverage).
Let[tS,tE] be a period of time. If for any t∈[tS,tE], the constellation C can completely cover ground region R at t, then we say that the constellation C can continuously cover ground region R at [tS,tE].Because the relevant judgement theorems about time factors have not been included in this paper, when dealing with the time factors, the discretization method is adopted. The basic steps for solving continuous coverage problems are as follows:Step 1.
Given a time partition precision corresponding to time intervalΔt, [tS,tE] can be divided into a group of time points t0,t1,…,tn where tS=t0, tE=tn, and tk-tk-1=Δt for any k=1,2,…,n.Step 2.
For any time pointtk, if constellation C cannot singly (or multiply) cover a ground region R completely at time tk, then C cannot continuously cover R at [tS,tE]. If C can completely cover R completely at all of the time points, then C can continuously cover R at [tS,tE].However, different from the approach for judging the coverage state at a certain time, the approach for evaluating the continuous coverage problem may have calculation errors whenΔt is not small enough. Thus, in order to control calculation errors, Δt is chosen to be a sufficiently small number.


## 5.3. Calculation of the Critical Beam Angle

Theorems11 and 12 can only obtain a result of “yes” or “no,” which is a qualitative result. However, in actual applications, a quantitative result is more useful. So, in this subsection, we will use Theorems 11 and 12 to obtain a typical relevant quantitative result by calculating the critical beam angle of the sensor.Definition 14 (the critical beam angle).
Consider an isomorphic constellationC, which means that all satellites in the constellation have the same semimajor axis, eccentricity, and orbit inclination, and all satellites carry same type sensor. Then, for ground region R, the critical sensor beam angle is the minimum beam angle of sensor with which the constellation C can completely cover the ground region R.For constellationC, if the half-beam angle of sensor equipped on the satellite that belongs to C is α, then the coverage region of C at t is written as ΛC(α,t). Then, clearly, if α<β, then ΛC(α,t)⊂ΛC(β,t); that is, the results are monotonous relative to the half-beam angle of sensor. So, the critical beam angle for C to R can be obtained by the bisection method.Denote the radius of the Earth byRE, and denote the semimajor axis of every satellite by a; then, given an elevation constraint of the sensor, the maximum half-beam angle denoted by αmax is arcsin⁡RE/a, and the minimum half-beam angle denoted by αmin is 0.Step 1.
For constellationC, calculate αmax and αmin.Step 2.
IfR⊄ΛC(αmax,t), then the algorithm terminates and the critical half-beam angle does not exist.Step 3.
Letαmid be equal to (αmax+αmin)/2. Then, check the coverage state of C with half-beam angle αmid to R at t; then, if R⊂ΛC(αmid,t), then αmax=αmid; else, αmin=αmid.Step 4.
LetΔα be the allowable error, which is given according to the requirements of the actual problem. If αmax-αmin<Δα, then the algorithm terminates, and αmid is the ultimate critical half-beam angle; else, jump to Step 3.The critical half-beam angle ofC to R at t is denoted by η(C,R,t).Then, for a time period[tS,tE], the critical half-beam angle of C to R, which is denoted by η(C,R,[tS,tE]), is defined as the maximum value of the critical half-beam angle of C to R at every moment in [tS,tE]. Then,(4)ηC,R,tS,tE=maxt∈tS,tE⁡ηC,R,t.


## 6. Experiment

To verify the correctness and effectiveness of the theorem discussed in this paper, a group of experiments are designed.For a certain constellation, the time for judging whether the constellation can cover a ground point is nearly a constant, so the computing time is associated with the number of ground points that need to be judged. So, to measure the effectiveness of the algorithm, we establish a performance index, which is the number of boundary feature points.


### 6.1. Single-Coverage Problem

Consider a Walker-Delta constellation with a configuration of 45/5/1. All satellites in the constellation have circular orbits and have the same altitude, which is 1300 km. The ground region is a latitude zonal target, the upper and lower bounds of which are 50° and 0°, respectively. Consider a simulation time period that lasts for three minutes with the start time 2017-01-01 0:0:0 and the end time 2017-01-01 0:3:0.The variation curve of the average number of boundary feature points and the orbital inclination during the simulation time period are presented in Figure3. Then, the variation curve in Figure 3 represents single-peak characteristics. When the orbital inclination is y = 90°, the number of boundary feature points reaches the minimum. The maximum number of boundary feature points is lower than 500; that is to say, if we want to judge whether the constellation can completely cover the ground region, no more than 500 points need to be checked. Thus, the efficiency for judging whether the constellation can completely cover the ground region is very high.Figure 3
Variation of the average number of boundary feature points with orbital inclination.The variation curve of the critical sensor beam angle with orbital inclination is presented in Figure4. Then, from Figure 4, we can see that the curve is nearly symmetric along the axis of y = 90°. The curve has three distinct troughs and four distinct peaks. The global minimum value is achieved when the orbital inclination equals 37°, at which point the corresponding critical beam angle is 54.6°.Figure 4
Variation of the critical beam angle with orbital inclination.To verify the correctness of the algorithm, we calculated this coverage problem by means of the grid point approach. We considered an orbital inclination for all satellites of 37°, and then we calculated the relationship between the continuous coverage rate during the simulation time period and the sensor beam angle of all the satellites.Figure5 illustrates the variation of coverage rate over orbital inclination. The x-axis represents the orbital inclination and the y-axis gives the corresponding coverage rate in the graph. It denotes that the constellation can continuously cover the ground region, while the orbital inclination is over 54.6°. The result shows that this approach is correct for solving the single-coverage problem.Figure 5
Variation of coverage rate with sensor beam angle.


### 6.2. Multiple-Coverage Problem

Consider a Walker-Delta constellation with a configuration of 18/3/0. The altitude of the satellites belonging to the constellation is 25000 km. Consider the global quadruple continuous coverage problem with the simulation time period 2017-01-01 0:0:0–0:53:0, which is the PU of this Walker constellation.The variation of the critical beam angle with orbital inclination is presented in Figure6. For this multiple-coverage problem, the variation of critical beam angle with orbital inclination is very similar to that in single-coverage problems.Figure 6
Variation of the critical beam angle with orbital inclination.The average number of boundary feature points during the simulation time period with the change in orbital inclination is as in Figure7. The curve is very regular and symmetric. The average number of boundary feature points is between 270 and 290.Figure 7
Variation of the number of boundary feature points with orbital inclination.Figures6 and 7 reflect that the best orbital inclination is approximately 54.2°, where the critical sensor beam angle is 11.3°.Additionally, to verify the correctness of the algorithm for solving multiple-coverage problems, the grid point approach is also used. We considered that the orbital inclination of all satellites is 54.2°; then, we calculated the global quadruple continuous coverage problem during the simulation time period.From Figure8, we can see that when the sensor beam angle is greater than 11.3°, the constellation offers quadruple continuous coverage of the globe. This result shows that this approach is also correct for solving the multiple-coverage problem.Figure 8
Variation of the global quadruple continuous coverage rate with sensor beam angle.


## 6.1. Single-Coverage Problem

Consider a Walker-Delta constellation with a configuration of 45/5/1. All satellites in the constellation have circular orbits and have the same altitude, which is 1300 km. The ground region is a latitude zonal target, the upper and lower bounds of which are 50° and 0°, respectively. Consider a simulation time period that lasts for three minutes with the start time 2017-01-01 0:0:0 and the end time 2017-01-01 0:3:0.The variation curve of the average number of boundary feature points and the orbital inclination during the simulation time period are presented in Figure3. Then, the variation curve in Figure 3 represents single-peak characteristics. When the orbital inclination is y = 90°, the number of boundary feature points reaches the minimum. The maximum number of boundary feature points is lower than 500; that is to say, if we want to judge whether the constellation can completely cover the ground region, no more than 500 points need to be checked. Thus, the efficiency for judging whether the constellation can completely cover the ground region is very high.Figure 3
Variation of the average number of boundary feature points with orbital inclination.The variation curve of the critical sensor beam angle with orbital inclination is presented in Figure4. Then, from Figure 4, we can see that the curve is nearly symmetric along the axis of y = 90°. The curve has three distinct troughs and four distinct peaks. The global minimum value is achieved when the orbital inclination equals 37°, at which point the corresponding critical beam angle is 54.6°.Figure 4
Variation of the critical beam angle with orbital inclination.To verify the correctness of the algorithm, we calculated this coverage problem by means of the grid point approach. We considered an orbital inclination for all satellites of 37°, and then we calculated the relationship between the continuous coverage rate during the simulation time period and the sensor beam angle of all the satellites.Figure5 illustrates the variation of coverage rate over orbital inclination. The x-axis represents the orbital inclination and the y-axis gives the corresponding coverage rate in the graph. It denotes that the constellation can continuously cover the ground region, while the orbital inclination is over 54.6°. The result shows that this approach is correct for solving the single-coverage problem.Figure 5
Variation of coverage rate with sensor beam angle.


## 6.2. Multiple-Coverage Problem

Consider a Walker-Delta constellation with a configuration of 18/3/0. The altitude of the satellites belonging to the constellation is 25000 km. Consider the global quadruple continuous coverage problem with the simulation time period 2017-01-01 0:0:0–0:53:0, which is the PU of this Walker constellation.The variation of the critical beam angle with orbital inclination is presented in Figure6. For this multiple-coverage problem, the variation of critical beam angle with orbital inclination is very similar to that in single-coverage problems.Figure 6
Variation of the critical beam angle with orbital inclination.The average number of boundary feature points during the simulation time period with the change in orbital inclination is as in Figure7. The curve is very regular and symmetric. The average number of boundary feature points is between 270 and 290.Figure 7
Variation of the number of boundary feature points with orbital inclination.Figures6 and 7 reflect that the best orbital inclination is approximately 54.2°, where the critical sensor beam angle is 11.3°.Additionally, to verify the correctness of the algorithm for solving multiple-coverage problems, the grid point approach is also used. We considered that the orbital inclination of all satellites is 54.2°; then, we calculated the global quadruple continuous coverage problem during the simulation time period.From Figure8, we can see that when the sensor beam angle is greater than 11.3°, the constellation offers quadruple continuous coverage of the globe. This result shows that this approach is also correct for solving the multiple-coverage problem.Figure 8
Variation of the global quadruple continuous coverage rate with sensor beam angle.


## 7. Conclusion

In this paper, we discussed the problem of judging whether a satellite constellation can offer complete single or multiple coverage for a ground region. To solve this problem, a group of judgement theorems and definitions are proposed.The numerical experiments indicated that the proposed approach can solve the problem correctly and efficiently. Compared with traditional approaches to solving this problem, such as the grid point approach and the spherical triangulation method, the process of this algorithm is simple, and the amount of computation is much less, which means that this approach has a high computational efficiency. Moreover, this approach to solving coverage problems at a given time has no computational error, so computational accuracy does not need to be considered.However, one drawback is that this approach can only judge whether the constellation can offer complete single or multiple coverage for the ground region. If the constellation cannot cover the ground region completely, the coverage rate or other coverage characteristics cannot be obtained from this approach.Additionally, the group of judgement theorems can only judge the coverage state at a certain time. However, by means of discretization, the coverage state for a period of time can be calculated. Future studies include proposing judgement theorems for solving coverage problems for a period of time without using the discretization method. Additionally, the accumulated coverage problem as well as the discontinuous coverage problem [24, 26] should be addressed.



---

*Source: 5097324-2018-02-22.xml*