#### Abstract

Grid point technique is a classical method in computing satellite constellation coverage to the ground regions. Aiming at improving the low computational efficiency of the conventional method, a method using latitude stripe division is proposed, which has high efficiency, and we name it latitude stripe method. After dividing the target region into several latitude stripes, the coverage status of each latitude stripe is computed by means of the spherical geometry relationship in the first orbital period. The longitude coverage intervals in the remaining orbital periods are computed by sliding the coverage status in the first orbital period. Based on this method, the instantaneous and cumulative coverage in simulation time can be calculated more efficiently. As well, the relationship between the cumulative coverage and altitude can be computed fast by this method, which could be used in the optimized design of repeating sun-synchronous orbits. The comparison between the conventional grid point method and the latitude stripe method shows that the latitude stripe method has high efficiency and accuracy. Through various case studies, the optimization in repeating sun-synchronous orbits design is successfully represented.

#### 1. Introduction

With the development of civilian satellite technology, constellation composed of a number of satellites plays an important role in remote sensing [1], communication [2], navigation [3], and environmental monitoring [4]. Remote sensing satellites for earth imaging and telecommunication satellites for uplink/downlink data exchange are both involving the analysis of the satellite constellation coverage to ground. Constellation coverage performance to ground is one of the indicators to judge the merits of a satellite constellation design. Rapidity and accuracy of the coverage performance evaluation have an extremely important meaning in optimized design of satellite constellations.

Grid point method is proposed by Morrison [5]; he divided the ground region into some grid points. By judging whether the satellite can cover each grid point, the instantaneous and cumulative coverage was given. Wang et al. [6] investigated the Common Aero Vehicle (CAV) constellation based on the analysis of earth coverage of a single CAV according to the grid point method. Deng et al. [7] divided the space regions into some grid points and analyze the space coverage performance of infrared LEO constellation. However, all of them divide the target area into several grid points to analyze the coverage performance. When the target region is large and the high division accuracy is required, it costs a lot of time to calculate the coverage performance of satellite constellation.

In addition to the grid point method, there are other ways to solve the coverage problem. For the satellite constellation design problem with different kinds of coverage requirements, Walker [8] and Ballard [9] proposed the circular orbit Walker Constellation and Rosette Constellation respectively for the continuous global coverage. Palmerini [10] and Ulybyshev [11] designed the elliptical orbit constellations for single continuous coverage. Rider and Adams [12, 13] analyzed the coverage region of the ground track of the satellite by introducing the concept of “Streets-of-Coverage” and investigated the polar and inclined orbits satellite constellation design for the single and multiple coverage requirements. The methods mentioned above are only suitable for a certain given area. Song et al. [14] developed the concept of “longitude stripe,” then divided the target area into several longitude stripes, and proposed a method for the instantaneous and cumulative coverage capacity calculation. The accurate coverage rate can be obtained by this method, but it based on the precise division of the target region. Aiming at the design of repeating orbit, an analytic algorithm is addressed for global coverage. The proposed method is mapping every pass into its ascending and descending nodes on the specified latitude circle, and then it accumulates the projected width on the circle by the field of view of the satellite [15]. However, this method can only be used to judge complete coverage of the globe, and it can not derive an accurate coverage. Some of above methods have low calculation accuracy, while others can only be used to special cases, resulting in limited practical applications.

When designing the remote sensing constellations, the percentage area of coverage is one of the satellite performances to be considered [16]. Currently, the theoretical method of orbit calculation of the repeating orbits is quite mature [17, 18]. As well, the coverage percentage is considered in the optimization of repeating orbits [19]. However, all the researches mentioned above are achieved by the numerical calculation, which costs a lot of computation. As a design tool, it is preferable to a semianalytic method due to the ergodic performance on the orbit altitude.

This paper presents a new method to calculate the coverage between satellite constellation and the area target based on latitude stripe division, which is named latitude stripe method. First, the globe is divided by several narrow latitude stripes, and then the longitude intervals of each latitude stripe in the area target are obtained. By computing the coverage status of each latitude stripe, the coverage of satellite is acquired. Through analyzing the ground track of satellite, the relationship of coverage between the same latitude stripe and satellite orbital period is derived. Consider the rotation of the earth in earth-centered inertial coordinate (ECI) system. The ground track of satellite is slid, as well as satellite coverage longitude interval of the same latitude stripe. Therefore in coverage calculation, the coverage status of the latitude stripe is computed in first orbital period, then the latitude coverage intervals in the remaining orbital periods are computed by sliding the coverage status in the first orbital period, and, finally, the comprehensive statistics is made to get the coverage rate.

#### 2. Basic Concepts and Definitions

First, a series of terms used in this paper are defined. The single satellite is recorded as . Single satellite coverage of the ground is recorded as . The constellation is recorded as , the global area is recorded as , and the ground area target is recorded as . The longitude and latitude of a point on the surface of the earth are denoted as and .

*Definition 1 (latitude stripe). *Latitude stripe of which central latitude is refers to a region on the surface of the earth. Suppose stripe width is ; then it indicates a region of which latitude interval is and longitude interval is .

*Definition 2 (longitude interval). *The longitude range of a latitude stripe is and the longitude interval of a latitude stripe represents a subrange of the longitude range, which is recorded as . The latitude stripe with the longitude interval represents a region, which is surrounded by longitude range and latitude range .

*Definition 3 (regional longitude interval). *Regional longitude interval is the part of the latitude stripe which is in the region ; it is a longitude interval. The regional longitude interval is recorded as ; when the latitude stripe width is not emphasized, it is denoted as . That is, .

*Definition 4 (coverage longitude interval). *Coverage longitude interval is the intersection of latitude stripe and the coverage range of a single satellite to ground is a longitude interval. The coverage longitude interval is recorded as ; when the latitude stripe width is not emphasized, it is denoted as . That is, .

*Definition 5 (regional coverage longitude interval). *Regional coverage longitude interval is the intersection of the regional longitude interval and the coverage longitude interval, which is recorded as . , where denotes a longitude interval. This item represents the coverage ranges between the single satellite and latitude stripe in the target region.

For each satellite in the constellation, the regional coverage longitude intervals can be calculated, and the coverage longitude intervals of constellation are the unions of all satellite coverage longitude intervals:

The width of latitude stripe is and the longitude interval of latitude stripe is ; the area corresponding to the latitude stripe is

The regional coverage intervals corresponding to a latitude stripe may be a single interval or a range of intervals. So the coverage area of the constellation to a latitude stripe is

The total coverage area of the constellation to target region is

The coverage rate of the satellite constellation to target region at a certain time is where is the area of target region. In this paper, the target region can be an arbitrary region, the latitude zone region, or the global region.

#### 3. Basic Ideas of Latitude Stripe Method

Here are basic ideas of latitude stripe method: according to a given division granularity, the globe is divided into several small latitude stripes. By computing the intersections of the target region and latitude stripes, the longitude intervals of each latitude stripe in the target are obtained. Figure 1 is the global division by latitude stripe method.

Then, the simulation time is divided into a number of orbital periods. Firstly, the coverage status of latitude stripes in the first satellite orbital period is analyzed. According to the characteristics of the earth’s rotation, the latitude coverage intervals in the remaining orbital periods are computed by sliding the coverage status in the first orbital period. Finally, the coverage of each latitude stripe by satellites in the constellation is combined; comprehensive statistics is made to get the instantaneous coverage rate and accumulated coverage rate.

The calculation procedure of the instantaneous coverage of a single satellite on the target area at time is given as follows.

*Step 1. *According to the user’s demand, the division granularity of latitude stripes is set as and then the globe is divided into several latitude stripes. These stripes form a set of latitude stripes.

*Step 2. *With the intersection of the latitude stripes and the target , the latitude stripes in target are obtained; if the target area is the global region, this step can be ignored. Through this step, a set of latitude stripes which has intersection with target can be obtained.

*Step 3. *Choose a latitude stripe which is not calculated in the set of latitude stripes.

*Step 4. *Calculate the intersection of and the target ; then regional longitude interval is obtained.

*Step 5. *Calculate the coverage area of satellite at time .

*Step 6. *Get the intersection of latitude stripe and coverage area and the coverage longitude interval is obtained.

*Step 7. *Get the intersection of regional longitude interval and the coverage longitude interval and the regional coverage longitude interval of satellite is obtained.

*Step 8. *Turn to Step 3 when a latitude stripe which is not calculated in the set of latitude stripes exists.

*Step 9. *Find out the area corresponding to the set of regional coverage longitude intervals; thus the coverage rate of target region of the satellite at time is obtained.

This method is similar to conventional grid point method in some aspects, but the two methods are different in essence. The essential idea of grid point is a method based on probability sampling, but the latitude method is a method based on calculus infinitesimal. This method divides the region by several latitude stripes. When the latitude stripe width is narrow, the change in the longitude direction is small. At the same time, this method is very similar to the longitude stripe method [14]; both of them use stripes to divide the target region; while the longitude stripe method uses longitude stripes to divide the target area, the latitude stripe method uses latitude stripes. The latitude circles do not change with the rotation of the earth, so the coverage calculation of satellite constellation is simplified by the combination of the characteristics of earth rotation and satellite orbital period. According to this method, computational efficiency is improved.

#### 4. Coverage Computations for Satellite Constellation to Ground

##### 4.1. Coverage Computations in One Orbital Period

Since the target region is given, so a set of latitude stripes in this region can be obtained. The longitude intervals of these latitude stripes are . Take a latitude stripe of which central latitude is in the set of latitude stripe as an example. The geometric relationship between the satellite and the latitude strip is shown in Figure 2.

If the satellite does not cover the at the moment , the coverage rate of satellite to the latitude stripe is 0. The longitude and latitude of subsatellite point are recorded as . is a segmental arc in the latitude circle of which latitude is . Points and are the intersection of the satellite coverage circle and arc segment . Assume that the longitude of point is and the coverage radius of satellite is . By taking into account the spherical triangle where is the North pole, the equation can be obtained:

And then (6) can be simplified to

Make , the longitude of point is , and the longitude of point is :

According to (8), the longitude and latitude of points and can be calculated and then the coverage longitude intervals of latitude stripe can be obtained. Assume that the regional longitude intervals of latitude stripe are . Taking into account that the longitude interval is the subinterval of , so the regional coverage longitude intervals of satellite to ground can be expressed as

##### 4.2. Coverage Computations in Simulation Time

Consider the rotation of the earth in earth-centered inertial coordinate (ECI) system. The ground track of satellite is slid, as well as satellite coverage longitude interval of the same latitude stripe. After a satellite orbital period, the satellite can revisit the same latitude stripe; namely, the time difference of visiting same latitude stripe is a satellite orbit period. According to this theory, the coverage status of latitude stripes in the first satellite orbital period is computed, and then the latitude coverage intervals in the remaining orbital periods are computed by sliding the coverage status in the first orbital periods.

This method only needs to calculate the intersection of coverage longitude intervals and regional longitude intervals when computing the coverage in remaining orbital period, which saves the time of computing the coverage longitude intervals.

The inclined line in Figure 3 is the ground track of subsatellite points and the ellipse is the coverage of satellite to ground on earth at a certain time.

In Figure 3, and are the intersections of the satellite coverage and the latitude stripe in the first orbital period. and are the intersections of the satellite coverage and the latitude stripe in the second orbital period. and are the intersections of the satellite coverage and the latitude stripe in the third orbital period. is obtained by translating with an interval; is obtained by translating with two intervals.

Assume that the longitude of is , namely, one of intersections of satellite coverage, and the latitude stripe is . After th orbital periods, the intersection of latitude stripe iswhere is the nodal period. Taking into account that the longitude interval range of a latitude stripe is , is needed to be adjusted into . is mean angular velocity of the earth. is, respectively, the temporal variations of the Right Ascension of Ascending Node in perturbation model. is the difference of and ; it can be expressed as

#### 5. Numerical Simulations for Instantaneous Coverage

##### 5.1. Instantaneous Coverage Computation

This chapter verifies the instantaneous coverage of the target region by different methods. Take the Walker constellation [8] with parameter as an example. The satellite orbit altitude is 1300 km, the inclination of the orbit is , and every satellite of the constellation is equipped with a simple conic sensor whose half-angle is . The target region is defined as a latitude zone ranging from to . Simulation start time is 2016-1-1 00:00:00 (UTCG) and end time is 2016-1-3 00:00:00 (UTCG). The simulation step is set as 1 s.

The target region is divided by different division granularities which, respectively, are , , and . The division granularity of the latitude strip method is defined as the distance (unit: degree) between the two adjacent stripes. The division granularity of the grid point method is defined as the distance (unit: degree) between the two consecutive grid points when the target area is evenly divided by the grid points. In Figures 4, 5, and 6, the constellation coverage rate to target region is computed in comparisons by two methods.

The horizontal axis represents time (unit: second), and the vertical axis represents the constellation coverage rate. Red solid line is the coverage rate by latitude stripe method and black dotted line by grid point method.

Comparing the results shown in Figures 4, 5, and 6 and computing the difference of the coverage rates between the two methods, Figures 7, 8, and 9 can be obtained. The horizontal axis represents time (unit: second), and the vertical axis represents the difference of constellation coverage rate by two methods.

Comparing the computation results in Figures 4–9, there are some analyses as follows:(1)The division granularity is smaller, the results of two methods are more accurate, and the difference of constellation coverage rate is smaller.(2)When the division granularity is , the difference of coverage rate by two methods is in the order of and the accuracy of latitude strip method is demonstrated.

In order to compare the computational efficiency of the latitude stripe method with the conventional grid method, longitude stripe method is introduced which is proposed in the paper [14]. The calculation accuracy is defined as the reciprocal of the division granularity. Figure 10 is the computation time for target region of three methods in simulation time. Three methods are grid point, latitude stripe, and longitude stripe.

From Figure 10, when the division granularity is the same, the time spent by the latitude strip method is the least in three methods. When the division granularity is , it costs 26 s to compute the coverage rate by latitude stripe, which is far below the time spent by grid point, at the same time below the time spent by longitude method. It shows the efficiency of the latitude stripe method.

##### 5.2. Result Analysis

From Section 5.1, when the division granularity is small, the computation result difference of the grid point method and the longitude strip method is small as well. When computing the coverage rate by grid point method, the target region is divided along the longitude and latitude direction into several grids. However, when using latitude stripe method, the target region is divided along the longitude direction into several strips. The coverage longitude intervals in the latitude direction are accurate, which eliminates the error in the latitude direction when dividing the target region. When the calculation accuracy is large, the calculation results of latitude strip method and grid point method are in agreement.

When being in the same division granularity, the calculation time of the latitude strip method is shorter than that of the grid point method. The reason is that the latitude strip method makes full use of the characteristics of the earth’s rotation. The coverage of the target in the simulation time can be obtained by only calculating the coverage of the target in the first orbital period, which can improve calculation efficiency. Grid point method and longitude stripe method do not consider the characteristics of the earth’s rotation, so they spend more time to compute the coverage.

#### 6. Numerical Simulations for Cumulative Coverage

##### 6.1. Design of Repeating Orbit

The repeating orbits are also known as the repeating ground track orbits, which are defined as orbits with periodic repeating ground tracks. Their ground track will repeat after a whole numbers of revolutions in nodal days. These orbits have good appearances for earth’s coverage, which is good for remote sensing. In engineering practice, the repeating factor represents the number of orbits completed per day [17, 18]. determines the location and sequence of all ground traces, which is defined as

In (12), can be written as an integer number plus a fractional part , where is an integer number which is prime to and .

For remote sensing satellites, the global coverage to the earth is one of most important parameters to be required. Currently, repeating orbits are most designed by finding out the optimal orbit altitude, but it costs a lot of time.

Based on the latitude stripe method proposed in this paper, the relationship between the cumulative coverage rate and the orbit altitude can be quickly obtained, which makes it more convenient for designer to select the best orbit. According to this method, the set of different orbit altitudes can be obtained, which satisfy the global coverage to earth. By selecting different orbit values in the altitude set, the designers can save a lot of time of designing repeating orbit and make the decisions quickly.

Orbit data of HJ-1A is used in this chapter, HJ-1A is a satellite of satellite constellation in the environment monitoring and forecasting [4], the design goal is to cover the earth within four days. HJ-1A is equipped with a simple conic sensor whose half-angle is , and this satellite only takes photographs in the periods of the descending orbits. Consider the goal of HJ-1A to cover the earth within 4 days. So simulation start time is 2016-1-1 00:00:00 (UTCG) and end time is 2016-1-5 00:00:00 (UTCG). The simulation step is set as 1 s and the orbit altitude of the HJ-1A is chosen in the range from 600 to 700 km.

The relationship between the coverage percentage and the altitude for HJ-1A is demonstrated in Figure 11. With the different orbit altitudes, the coverage rate of the HJ-1A to the globe is different as well within four days. When the range of altitude is 639–650 km, the coverage rate to the globe for the satellite is 100%. So when the goal is the global coverage, the altitude of the satellite can be chosen in the range of altitude 639–650 km.

Actually, the orbit altitude of HJ-1A is designed as 649.053 km according to the article [4], which is in the range of 639–650 km. It shows the accuracy of the calculation results by latitude stripe method. When the division accuracy is , it costs 30 s to compute the coverage rate by latitude stripe.

##### 6.2. Design of Repeating Orbit for CFOSAT

CFOSAT (Chinese-French Oceanic Satellite orbit) is a French-Chinese mission embedding two wave and wind measuring instruments [15]. CFOSAT owns a narrow field of view of which the half-angle is . It is sun-synchronous polar orbit, and it imagines in both ascending and descending orbits. The design goal is to cover the earth with 11–14 days.

The relationship between the cumulative coverage percentage and the altitude of CFOSAT during several periods is demonstrated in Figure 12. The figure obtained by latitude method is very similar to the figure in the article [15]. It makes more convenient for designer to choose the best orbits. There are some analyses for Figure 12 as follows:(1)Within 11 days, the orbit altitude of the satellite about 502 km, 512 km, 520 km, and 540 km will have good coverage characteristics. When the orbit altitude is 502 km, the cumulative coverage rate is 96.68%.(2)Within 14 days, when the orbital altitude changes in 501–503 km, 520–521 km, 526-527 km, and 538–548 km, the global coverage rate of the satellite constellation is 100%. Therefore, the altitude of the satellite can be chosen in the ranges above for a global full coverage problem. Table 1 shows the parameters of some repeating sun-synchronous orbits when the orbital altitude range is 501–503 km. In Table 1, is the ground track repetition parameter. These orbits have different repeat cycles and ground tracks, but all of them can achieve the global full coverage within 14 days. According to paper [20], when the orbital altitude is approximately the same, one orbit can transfer to another orbit using a low technique. This technology guarantees the fulfillment of several objectives in the course of the same mission, which is useful for different task.(3)The influence of transmit errors of the coverage analysis: take an orbital altitude 557 km within 12 days as an example; the cumulative coverage rate is 72.80%. When the orbital altitude is increased to 567 km, the coverage rate reduced to 20.7%. Meanwhile, according to the relationship between the cumulative coverage rate and the orbital altitude, it can develop detail strategies of orbit control. In addition, the overall design could be combined with the working condition of payload at the different altitude and weigh the choice of them.

When the division granularity is and the simulation time is 14 days, it costs 50 s to compute the coverage rate by latitude stripe.

#### 7. Conclusions

This paper proposes a new method for the satellite constellation coverage to the ground. This method divides the target region by several narrow stripes. One contribution is that this method eliminates the errors in the latitude direction when dividing the target region, so the calculation result is more accurate. Another contribution is that it makes full uses of the characters of earth’s rotation, so calculation efficiency is improved. According to this method, instantaneous and cumulative coverage rate to the region can be quickly computed. When this method is applied to the design of repeating sun-synchronous orbit, the relationship between coverage percentage and orbit altitude can be quickly obtained. Compared with the traditional method, this method makes it more convenient for the designers to select the best orbit and it is possible to develop detailed strategies of orbit control as well.

#### Competing Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

#### Acknowledgments

This work is supported by Natural Science Foundation of China under Grants no. 41571403, 61472375, and 61103144. The project was supported by the Fundamental Research Founds for National University, China University of Geosciences (Wuhan), under Grant no. 1610491B21. The corresponding author also acknowledges the encouragement and care by Qin Luo (father) and E. Liu (mother).