Reconfigurable Satellite Constellation Design for Disaster Monitoring Using Physical Programming

Data collection by satellites during and after a natural disaster is of great significance. In this work, a reconfigurable satellite constellation is designed for disaster monitoring, and satellites in the constellation are made to fly directly overhead of the disaster site through orbital transfer. By analyzing the space geometry relations between satellite orbit and an arbitrary disaster site, a mathematical model for orbital transfer and overhead monitoring is established. Due to the unpredictability of disasters, target sites evenly spaced on the Earth are considered as all possible disaster scenarios, and the optimal reconfigurable constellation is designed with the intention to minimize total velocity increment, maximum and mean reconfiguration time, and standard deviation of reconfiguration times for all target sites. To deal with this multiobjective optimization, a physical programming method together with a genetic algorithm is employed. Numerical results are obtained through the optimization, and different observation modes of the reconfigurable constellation are analyzed by a specific case. Superiority of our design is demonstrated by comparing with the existing literature, and excellent observation performance of the reconfigurable constellation is demonstrated.


Introduction
A global or regional observation satellite constellation is used for some emergency or significant missions such as disaster monitoring, and usually, the observation pattern of a specific constellation is restricted by both on-board instruments and geometric configuration of satellites [1,2]. However, natural disasters such as earthquakes, forest fires, or extreme weathers occur randomly, and usually, data collection is achieved only when satellites happen to fly over the disaster location [3,4]. To get timely and adequate observation for the disaster area, orbital maneuvers of satellites are required. Therefore, configuration of a constellation needs to be adjusted due to the change of observation mission, and reconfiguration of the constellation is demanded.
The satellite constellation reconfiguration is defined as "a deliberate change of the relative arrangements of satellites in a constellation by addition or subtraction of satellites and orbital maneuvering in order to achieve desired changes in coverage or capacity" [5]. Usually, optimization of constellation reconfiguration includes final constellation configuration design and optimal orbital transfer, and the purposes are to minimize the total fuel consumption as well as the reconfiguration time [6,7]. Based on the satellite constellation reconfiguration, the concept of the reconfigurable constellation (ReCon) is proposed [8,9]. Operation of a ReCon comprises the following two modes: global observation mode (GOM) for normal operations and regional observation mode (ROM) for contingent responses, and due to the flexibility of the ReCon, timely and adequate observation for uncertain areas of interest can be achieved [10].
In recent years, plenty of research about constellation reconfiguration has been carried out. De Weck et al. [5] used the auction algorithm to optimize the reconfiguration of constellation, and the optimal orbital transfer strategy for all satellites was designed with the purpose of minimizing the total fuel consumption. Reconfiguration time is not considered as a design metric in their work. Hong et al. [11] introduced a concept and a baseline design of a dual-mode disaster monitoring constellation, and a normal mode as well as a disaster mode is considered. However, the constellation needs extra on-orbit propellant depots to provide plenty of fuel for the mode transformation. Chen et al. [12] investigated the observation for some particular targets on the ground when an emergent requirement arises in a short period and proposed a variable-size multiobjective differential evolution method to deal with the reconfiguration of constellation, where the main idea of the reconfiguration is to adjust the mean anomaly of satellites. The optimization is implemented with known initial constellations and target sites in their work. Fakoor et al. [13,14] used the hybrid Invasive Weed Optimization/Particle Swarm Optimization algorithm to deal with the constellation reconfiguration. Optimal assignment of the satellites and optimal orbital transfer were combined in one step to simplify the reconfiguration process. Yet, reconfiguration time is not considered as an optimization objective in their work. Hoskins et al. [3] aimed at forest fire monitoring and introduced the L-shaped method and stochastic programing approach to design the initial constellation and maneuver sequences for reconfiguration. The direct overhead monitoring for the disaster site was achieved after the reconfiguration, and results proved that their design performs better than the existing observation satellites. As for regional or specific target observation, repeating ground track orbit is a usual choice, which has been proved to have better partial coverage properties than the unrepeated one [15,16]. A ground track is considered as repeating or periodic if it repeats at a fixed time interval, and the interval is usually an integral multiple of one nodal day [17][18][19][20]. Hence, sites which are located on the ground track can be regularly observed. In addition, since the on-board sensors benefit from minimizing the distance between the satellite and the site to be monitored, an optimal way to achieve data collection for the disaster area is to make satellites fly directly over the target site [3]. Therefore, to achieve the timely and adequate observation for a natural disaster site, the repeating ground track orbit, of which the ground track passes the site directly, can be chosen as the target orbit for each satellite in the initial constellation. After the reconfiguration, each satellite in the final constellation is able to fly directly over the target site for data collection in a fixed period, and the period can be designed as one or more nodal days according to the requirement of a specific mission.
Previous investigations mentioned above considered plenty of approaches to optimize the reconfiguration process and also the use of repeating orbit for partial coverage. However, most of them were carried out with a known initial constellation, a known final constellation, or a set of known observation locations, and few of them considered the reconfiguration for monitoring the unpredictable disaster site as well as the design of a reconfigurable initial constellation. Besides, most of the previous researches took fuel consumption as the optimization target, and few of them considered the optimization of reconfiguration time or tradeoff between them. In addition, none of the current constellations with the task of observing the Earth can achieve active maneuvers to collect data for the natural disaster or other significant events [3]. Therefore, the design of a ReCon according to the specific constraint of fuel consumption, reconfiguration time, and unpredicted disaster sites is of great significance.
In this paper, we aim to design a ReCon, which can achieve timely and adequate observation for an arbitrary and unpredictable disaster site through orbital maneuvers. Satel-lites in the constellation flying directly over the disaster area for data collection are the purpose of the reconfiguration. Due to the unpredictability of disaster, target sites evenly spaced on the Earth are considered as all possible disaster scenarios. For the optimal ReCon, several optimization objectives are considered, which are (1) minimizing fuel consumption, which is equivalent to minimizing the total velocity increment during reconfiguration; (2) minimizing maximum and mean reconfiguration time to ensure timely observation for the target site; and (3) minimizing standard deviation of reconfiguration times to ensure similar reconfiguration ability for all target sites. To deal with this multiobjective optimization, a physical programming method together with a genetic algorithm is employed. The preferred design metrics in physical programming are considered according to the specific requirement or preference of the designer, which makes the design much more feasible, practical, and reasonable. Furthermore, GOM and ROM of the ReCon are analyzed to demonstrate the superior observation performance for normal global scan as well as some emergency disaster monitoring. As a whole, the main contributions of our work can be summarized as follows: firstly, a mathematical model for orbital transfer and overhead monitoring is established by analyzing the space geometry relations between satellite orbit and an arbitrary disaster site. Secondly, a ReCon for unpredictable disaster monitoring is designed, and several optimization objectives are considered to optimize the reconfiguration process as well as to ensure that the ReCon has similar reconfiguration ability to observe any possible disaster site. Thirdly, two observation modes are presented in the ReCon, and normal global scan as well as emergency disaster monitoring with excellent observation performance can be achieved, respectively.
The rest of this paper is organized as follows: The detailed problem formulation process is described in Section 2. An overview of the physical programming method is proposed in Section 3. Numerical results of a design example and comparison with existing literature are given in Section 4. Finally, conclusions are presented in Section 5.

Problem Formulation
Problem formulation is presented in this section. To simplify the formulation process, as well as to make the model easier to establish and understand, some reasonable and necessary assumptions ought to be made. Firstly, maneuvers for orbital transfer are impulsive, and this is reasonably accurate during the initial mission planning [19]. Secondly, the specific sensors are not modeled, since all sensors (e.g., visual spectrum cameras, infrared cameras, and radar sounders) benefit from minimizing the range between the satellite and the area being monitored [3]. Mathematic models for orbital transfer and overhead monitoring, as well as the description of the optimization problem for the optimal ReCon design, are presented as follows. International Journal of Aerospace Engineering modes for the ReCon, which are GOM for normal operations and ROM for contingent responses [10]. When an unpredictable disaster occurs, the observation model of the constellation can be converted to ROM through reconfiguration, and satellites in the constellation are able to fly overhead the disaster site for monitoring through orbital maneuvers. In this work, constellation with GOM is defined as initial constellation, while constellation with ROM is the target constellation, and initial constellation is designed to achieve optimal reconfiguration for monitoring any possible disaster. Features of the initial and target constellation, as well as the orbital transfer, can be described as follows.
(1) Initial Constellation. For the initial constellation with GOM, the number of satellites is K, which is also the number of orbit planes. Eccentricity and argument of perigee are both set to be 0 for each orbit. Radius r 1 and inclination i are set to be the same for each orbit. Right ascension of the ascending node (RAAN) and argument of latitude (AOL) of orbit kðk = 1, 2, ⋯, KÞ are Ω k and u k , respectively. It should be mentioned that for all satellites, Ω k are chosen to be evenly distributed along the equator to ensure the uniformity of the observation, and for the first satellite, Ω 1 = 0; thus, we have Ω k = 2πðk − 1Þ/K. u k is the angle between the ascending node and the satellite, and 0 ≤ u k < 2π. As a result, the design variables of the initial constellation are K, i, r 1 , and u k , where k = 1, 2, ⋯, K. J 2 perturbation of the Earth is considered in this work. For the initial orbit, assume that change rates of the RAAN, the argument of perigee, and mean anomaly are _ Ω 1 , _ ω 1 , and _ M 1 , respectively. n 1 denotes the mean angular velocity of the satellite. Then, the nodal period T p1 of the initial orbit can be defined as Change rates _ Ω 1 , _ ω 1 , and _ M 1 can be expressed as where ξ = 3R 2 e J 2 /4r 2 1 , R e denotes the average radius of the Earth, and the value of J 2 is 1:082626 × 10 −3 [21].
Also, for the nodal period of initial orbit, the corresponding angular velocity n 1p , which is different from n 1 , is defined as (2) Target Constellation. As for the target constellation with ROM, the number of orbit planes, the RAAN, and the inclination stay the same as the initial constellation, for only coplanar orbital transfer is considered. The type of each target orbit is considered to be the repeating ground track circular orbit, and radius r 2 of the orbit is determined by the number of revolutions and number of nodal days in the repetition period. Also, circular orbit means constant distance between the satellite and the Earth, which can ensure the uniformity of observation.
A ground track is considered as repeating or periodic if it repeats at a fixed time interval allowing a cyclic observation of the Earth [1]. Let T p2 be the nodal period of the target orbit and T d be the nodal day. N p and N d are positive integers, which are prime to one another. For the repeating ground track orbit, it should be guaranteed that the period T p2 precisely matches the period T d , and this can be described as a mathematical condition, which is the following expression [20]: where T r is a fixed time interval, which is also the period of repetition. For the target repeating ground track orbit, assume that change rates of the RAAN, the argument of perigee, and mean anomaly are _ Ω 2 , _ ω 2 , and _ M 2 , respectively. n 2 denotes the mean angular velocity of the satellite, ω e denotes the angular velocity of the Earth, and μ denotes the geocentric gravitational constant. Then, the period T d and T p2 can be defined as Also, radius r 2 can be calculated by [1] r 3: For the nodal period of orbit, the corresponding angular velocity n p2 , which is different from n 2 , is defined as In this work, the repetition period of one nodal day is considered to ensure the high frequency of observation, which 3 International Journal of Aerospace Engineering also means that N d = 1. Also, it is convenient to adjust our model for other kinds of repetition periods.
It should be mentioned that n p1 and n p2 will be used as the angular velocity of the initial and target orbit in the following analysis.
(3) Orbital Transfer. In this work, Hohmann transfer, as the optimal two-impulse coplanar circular orbital transfer, is introduced. The total velocity increment Δv of Hohmann transfer is determined by the radius of the initial circular orbit r 1 and radius of the target circular orbit r 2 . Let n T = r 2 /r 1 , and μ denotes the geocentric gravitational constant; then, Δv can be expressed by The transfer time ΔT h is a half of the period of the elliptical transfer orbit, which is defined as where a T is the semimajor axis of the transfer orbit and a T = ðr 1 + r 2 Þ/2.

Mathematical
Modeling. Due to the unpredictability of disaster, target sites around the Earth's surface in step of 1 degree are set as possible scenarios, and all satellites in the constellation are considered to monitor these target sites through orbital maneuvers. To illustrate the model mathematically, one arbitrary target site on the ground and an orbit with one satellite are considered as the basic scenario. The main task of this subsection is to calculate the specific transition time of each satellite monitoring each target site in the basic scenario and then to form the whole constellation reconfiguration process. It should be noted that transition time refers to the time from the mission beginning to the satellite finishing the transfer, while reconfiguration time represents the maximum transition time of all satellites in the constellation. The ascending pass of a satellite is first considered for the monitoring mission. The basic scenario is presented in Figure 1. In the figure, O E XYZ is the J2000 inertial system. The intersection of the X-axis and the equator is O, and the intersection of the Z-axis and the Earth is P. Ground track of the orbit shown in the figure is defined on the nonrotational Earth in the inertial system, which is also the projection of the orbit. It is a circle and would rotate along the Z-axis at a regular angular velocity due to the perturbation. T is the target site.
Suppose that the longitude and latitude of the target site T are λ and φ, respectively, and the current time is assumed to be t 0 , which is also the mission start time. The corresponding right ascension and declination of the target site at t 0 moment are λ 0 and φ 0 , respectively, where 0 ≤ λ 0 < 2π and −π/2 ≤ φ 0 ≤ π/2. As for the satellite, inclination of the orbit is i, RAAN at t 0 moment is Ωðt 0 Þ, and eccentricity as well as argument of perigee are both set to be 0.
In the inertial system, target site T is rotating along the Z -axis at a regular angular velocity, which is equal to the Earth rotating angular velocity. Let i 0 = min ði, π/2 − iÞ be the latitude limit of the ground track; then, only target sites which satisfy −i 0 ≤ φ 0 ≤ i 0 can intersect the ground track and get the chance to be directly monitored by the satellite flying overhead, and hence, only this kind of target sites are considered in the following analysis.
After a specific time, target site T will first intersect the ground track at point T 1 due to the rotation of the Earth. In this paper, this specific time is defined as the "minimum distance 1" between the target site and the ground track, which is marked as ΔT l1 min , and is determined by TT 1 _ , where the number "1" represents the ascending pass (the same below). Only after at least ΔT l1 min then the satellite could possibly get the chance to fly directly over the target site T during the ascending pass. Obviously, the target site will intersect the ground track every single nodal day after the first intersection. Therefore, if the satellite cannot finish the orbital transfer before the target site intersects the ground track, then it will get another one nodal day time to wait and to perform the orbital transfer, and so on. Here, the number of nodal days is marked as D 1 , and it is defined as "waiting days" in this paper.
For point T 1 on the ground track, the declination is φ 0 and the right ascension is assumed to be λ 1 at t 0 moment. Make an orthodrome crossing point P and T and intersecting with the equator at point M. Also, make an orthodrome crossing point P and T 1 and intersecting with the equator at point M 1 . Let u be the geocentric angle of NT 1 _ . Then, in the spherical right triangle NT 1 M 1 , u and λ 1 can be expressed as Target Ground track Ascending node  International Journal of Aerospace Engineering where i ≠ π/2. For the case that i = π/2, it is obvious that u = φ 0 and λ 1 = Ωðt 0 Þ. Therefore, the defined "minimum distance 1" ΔT l1 min can be calculated by where the function ½X 2π ensures that X falls into the interval of ½0, 2πÞ. ω e is the angular velocity of the Earth, and _ Ω is the change rate of the RAAN.
Here, we take ΔT l1 as the time limit for the satellite to perform the orbital transfer. According to the analysis above, ΔT l1 can be expressed as where D 1 = 0, 1, 2, ⋯, is the number of waiting days. T d is the time interval of a nodal day.
As for the descending pass, similar analysis can be carried out. As shown in Figure 2, the "minimum distance 2" Δ T l2 min can be determined by TT 2 _ , where the number "2" represent the descending pass.
Let P 2 be the midpoint of T 1 T 2 _ . Make an orthodrome crossing point P and P 2 and intersecting with the ground track and equator at point P 1 and P 3 , respectively. Also, make an orthodrome crossing point P and T 2 and intersecting with the equator at point M 2 . Declination φ 0 and right ascension λ 1 of point T 1 are already obtained. For point T 2 on the ground track, the declination is φ 0 and the right ascension is assumed to be λ 2 at t 0 moment. Thus, for point P 1 on the ground track, the declination is i 0 and the right ascension is ðλ 1 + λ 2 Þ/2. For point P 2 , the declination is φ 0 and the right ascension is also ðλ 1 + λ 2 Þ/2.
Then, in the spherical right triangle NP 3 P 1 , it is easy to obtain that where i ≠ π/2, and for the case that i = π/2, it is obvious that Therefore, the defined "minimum distance 2" ΔT l2 min can be calculated by Similarly, we take ΔT l2 as the time limit for the satellite to perform the orbital transfer. According to the analysis above, ΔT l2 can be expressed as where D 2 = 0, 1, 2, ⋯, is the number of waiting days.
Once ΔT l1 and ΔT l2 are determined, the next task is to obtain the transition time of the satellite. In this paper, we define the required position of the satellite in the target orbit and the position of the satellite in the initial orbit as "target position" P T and "initial position" P I , respectively. Obviously, these two positions are changing over time, and the main idea of the transition is to wait for a suitable moment to achieve the Hohmann transfer.
As illustrated in Figure 3, r 1 denotes the radius of the initial circular orbit and r 2 denotes the radius of the target circular orbit. Let n p1 and n p2 be the angular velocity of the initial and target orbits, respectively, which are calculated by equations (3) and (8). To calculate the transition time, the direct overhead observation during the ascending pass of a satellite is analyzed below, and for the descending pass, analysis can be carried out in the same way.
At t 0 moment, AOL of the initial position is assumed to be u 1 , which is a given parameter. After ΔT l1 , the target position is right above the ground point T 1 presented in Figure 1, and the corresponding AOL is u. Therefore, at t 0 moment, the AOL u 2 of the target position can be expressed by Then, the phase difference Δu 0 between the initial and target position at t 0 moment can be obtained, which is As shown in Figure 3, if Hohmann transfer is to be performed, the phase difference Δu T of the initial and target position need to be restricted as where ΔT h is the transfer time defined by equation (10). Also, the target position is required to be P H , as illustrated in Figure 3. The difference between Δu 0 and Δu T is gradually changed due to the difference between n p1 and n p2 . Define Δu c as the phase difference that needs to be reduced, then For the case that r 1 > r 2 , the target position is "chasing" the initial position, until the target position is Δu T behind the initial position, and then, the first maneuver can be implemented. Define the "chasing" time as ΔT c , which is also the time for waiting before the orbital transfer, then ΔT c can be expressed by Therefore, the total transition time ΔT 1 can be expressed as For the case that r 1 < r 2 , the analysis can be carried out in a similar way. Also, for the descending pass observation, similar analysis can be done, and the total transition time ΔT 2 can be obtained. So far, time limits ΔT l1 and ΔT l2 and total transition times ΔT 1 and ΔT 2 are obtained, and the problem turns to be judging whether the satellite can achieve the orbital transfer within the time limit of ΔT l1 or ΔT l2 . It should be mentioned that in order to ensure the quality of observation, the satellite needs a time interval to prepare for the monitoring mission after the last impulse, and the interval in this work is set to be a half of the target orbit period 2π/n p2 [3].
In conclusion, for the basic scenario that one satellite overhead monitors one arbitrary target site during ascending pass, the calculation process can be summarized as the following seven steps.
Similarly, analysis for descending pass can be carried out, and the mathematical model for one satellite monitoring one target can be concluded as Figure 4.
As a result, when considering K satellites in the constellation maneuvering for monitoring S target sites around the Earth's surface, the same analysis given in Figure 4 can be carried out individually, and S × K sets of Δv, ΔT, and D can be obtained. Finally, the total reconfiguration model is established.

Optimization Problem.
Constellation design is generally a complicated optimization problem with plenty of design variables; especially when there are several design objectives, the computational process can be very time-consuming [22] [23]. In this paper, the optimal design is simplified by restricting some design variables. Reasons and details for the simplification are analyzed as follows.
Firstly, radius r 1 or altitude h 1 and inclination i of each initial circular orbit are the same. The RAAN Ω of each orbit is evenly distributed along the equator, and the inclination i is determined by the highest latitude of the region that needs to be monitored.
Secondly, radius r 2 or altitude h 2 of the target orbit is determined by the number of revolutions in the period of repetition, and it is a given reference for the design of the altitude h 1 of the initial orbit.
Thirdly, the repetition period of the target orbit is one nodal day, and there is only one satellite in each orbit plane. The number of satellites K, which is also the number of orbit planes, is determined by the specific observation frequency requirements and the economic conditions of the designer.
Therefore, once the number of satellites K, the inclination i, and the altitude h 2 of the target orbit are determined by the designer, design variables would be reduced to the altitude h 1 of the initial orbit and the AOL u k each initial orbit, where k = 1, 2, ⋯, K. For the first satellite, it is reasonable to set u 1 = 0, since only phase differences between each satellite matter.
Based on the analysis above, our design can be concluded as designing the optimal altitude of the initial orbit and phase differences between each orbit plane, and the purpose is to minimize the fuel consumption, maximum reconfiguration time, mean transition time of all satellites for all targets, the Before formulating the optimization problem, symbols of design variables and design metrics need to be defined. Assume that there are S (S ≥ 1) target sites on the Earth's surface and K (K ≥ 2) satellites in the constellation. Design variables are the altitude h 1 of the initial orbit and the phase differences Δu k between satellite k and the first satellite, where 0 ≤ Δu k ≤ 2π for k = 2, ⋯, K. For an arbitrary target site s (s = 1, 2, ⋯, S) and an arbitrary satellite k (k = 1, 2, ⋯, K), the impulse Δv, transition time ΔT s,k , and waiting days D s,k can be obtained based on the mathematical model described before. It should be mentioned that for a constellation, the reconfiguration time is defined as the maximum transition time of all satellites. Therefore, reconfiguration time ΔT s of the constellation for target site s can be defined as Also, the maximum waiting days of the constellation for the target site s can be expressed as In conclusion, the optimization problem for ReCon design can be formulated as find : min f i x ð Þ, i = 1, 2, 3, 4, 5 where design metrics f i are defined as where Δ T denotes the mean value of reconfiguration times for all target sites, which can be expressed by Calculate time limit T l2 0 , 0 , t 0 , i, (t 0 ), u 1 r 1 , r 2 , D 1 -D 2 = 0 x denotes the vector of design variables, and it is expressed as x LB and x UB are the lower and upper bound, and they are determined by the designer. For the altitude h 1 , h min and h max are defined as the lower and upper bound, while for the phase differences Δu k (k = 2, ⋯, K), 0 and 2π are the lower and upper bound. Thus, x LB and x UB can be expressed as where numbers of "0" and "2π" are both K − 1.

Physical Programming Method
A physical programming method is an efficient way to deal with multiobjective problems, of which the main idea is to convert the problem into a single-objective problem by using preference functions that capture the designer's preferences [24,25]. Tradeoffs between all design metrics are considered during problem formulation, and Pareto design which is in the preferred region of the Pareto surface is generated. Consequently, a compromise solution is obtained.

Preference Function Description.
Preference functions are defined with the intention to reflect the preferences of the designer, and each preference function corresponds to a design metric. Assume that there are n sc design metrics, and each of them is marked as g i ; also, value of the corresponding preference function is marked as H i ðg i Þ, where i = 1, 2, ⋯, n sc . The less the value of preference function H i ðg i Þ, the more satisfied the design metric g i . It should be noticed that g i is directly determined by the vector of design variables x, and hence, g i can also be expressed as function g i ðxÞ.
Generally, preferences in physical programming are divided into the following four main classes: Class 1, smaller is better; Class 2, larger is better; Class 3, value is better; and Class 4, range is better. Each class is divided into soft and hard types as well, such as Class 1S and Class 1H, where the suffix letters S and H represent soft type and hard type, respectively. For the soft one, value of the preference function is gradually changed in the feasible region, yet for the hard one, value of the preference function in the feasible region stays the same, which is set to be 0. As a whole, preference functions can be concluded to be eight basic types, as listed in Figure 5.
Furthermore, for soft types, design metrics can be divided into several regions. Here, we take Class 1S as an example, and the rest of soft types can be defined in the similar pattern. Regions for Class 1S can be listed as follows: Region 1: highly desirable range (g i ≤ g i1 ).
Region 6: unacceptable range (g i > g i5 ), where g i1 to g i5 are the region end points, which are given by the designer according to the specific preferences. These six regions are given in Figure 6.

Preference Function Formulation.
Preference functions are formulated based on the six interval divisions of the design metrics. By designing the value and the slope of preference function at each region end point, and fitting the value of preference function in each region, the specific preference function can be obtained.
Class 1S is illustrated as an example. Let k be the region number, and k = 1, 2, ⋯, 5. As for region 6, the preference function is a linear extension of the end point of region 5. Value and slope of the preference function at the k region end point g ik are H ik and s ik , respectively.
To calculate the value H ik , the modification ΔH ik of the preference function in the k region is introduced. For the first region, the value and the modification of the preference function are usually set to be 0.1, and the rests can be calculated by the following equations: where β is a constant determined by the designer and β > 5.
To calculate the slope s ik , the average slope s ik of the preference function in k region is introduced, which is defined as Then, the slope of the first region can be expressed as where α is a constant determined by the designer and α > 0. Slopes for the adjacent region can be expressed as Therefore, the value of preference function in each region can be obtained, which is illustrated as follows.
When k = 1, When k = 2, 3, 4, 5, International Journal of Aerospace Engineering and functions T 0 to T 3 are expressed as follows: As a result, the preference function of Class 1S is formulated. Other types of the preference functions can be formulated in a similar way.
Finally, the aggregate objective function is formed, and the multiobjective function is consequently converted to a single-objective function, which can be expressed as where x is the vector of design variables.
In this paper, n sc = 5, and g i is the design metric f i proposed in equation (27). For the classification of design metrics, velocity increment, maximum reconfiguration time, mean transition time, and the standard deviation of the reconfiguration time for all targets belong to Class 1S, while maximum waiting days belong to Class 1H.

Numerical Results
In this section, simulation based on the proposed model is carried out. Initial conditions and preferences of design metrics in physical programming are set up by referring to the literature [3], and also, results are compared with that obtained in the literature. Furthermore, the GOM and ROM of the ReCon are analyzed and compared through a specific case, and the superior observation properties of the ReCon are demonstrated.

Simulation Setup.
Initial conditions for the simulation are presented as follows. The inclination is selected as 90°f or the global monitoring, and the number of satellites K is set to 5, which are the same as the BLASTOFF constellation presented in literature [3]. Target sites around the Earth's surface in step of 1 degree are considered, ignoring the north and south poles, where satellites will pass by directly without any maneuvering. That is to say, 64440 (179 × 360) target sites in total are considered as possible disaster sites. For the final repeating ground track orbit, parameters introduced in Section 2 are listed in Table 1.
As for design variables and design metrics, specific ranges are determined according to literature [3]. In the referenced literature, minimizing the change in the orbital period is introduced as a means to minimizing the fuel consumption. The minimum period change is 12.125 minutes for 5 satellites, and the corresponding velocity increment for each satellite is about 62.961 m/s. Considering that the constellation in our design needs to go back to GOM after disaster monitoring, two Hohmann transfers are required for each satellite. Thus, if half of 62.961 m/s is used for the transfer, then the initial orbital altitude is about 605.8 km. Therefore, it is reasonable to set boundaries of the initial altitude as h min = 580 km and h max = 620 km, and the corresponding velocity increments are about 17.5 m/s and 39.2 m/s. Also, the RAAN and the AOL of  Figure 5: Preference function classification for physical programming.

Highly desirable Desirable Tolerable
Undesirable Highly undesirable Unacceptable H i Figure 6: Preference function regions for Class 1S. 9 International Journal of Aerospace Engineering the first satellite, which is also the anchor satellite, are both set to 0, for only the relative positions of satellites matter.
Based on the boundaries of the initial altitude, the differences of angular velocity of the initial and target orbits are determined, and the transfer ability can be expressed as equation (22). Therefore, it is reasonable to set the boundaries of the maximum reconfiguration time as 100 and 200 hours, and the maximum waiting days are set to 6 days. Also, preferred value of the average reconfiguration time is set to 48 hours according to the BLASTOFF constellation in the literature [3]. As for the standard deviation of the reconfiguration time, 0 is a preferred value to ensure that the constellation has the similar ability to monitor any target site through orbital maneuver. As a whole, feasible regions of design metrics can be designed, and details of the region end points for Class 1S are provided in Table 2.
A genetic algorithm is employed for the optimization. Parameters used for the genetic algorithm are listed in Table 3.

Optimization
Results. Based on the introduced physical programming method and genetic algorithm, optimization is carried out, and results are given as follows.
Change of the best and mean fitness value during the optimization are presented in Figure 7, and these two values tend to be convergent within the given tolerance after 83 generations.
Values of design variables are listed in Table 4. From the results presented in Table 4, we can find that each phase difference is nearly an integral of 72. Actually, due to the uniformly distributed target sites and the uniformly distributed orbits in this problem, it is easy to form an intuitive conclusion that satellites should be distributed uniformly, and the optimization results exactly verified this conclusion.
Furthermore, if we consider the walker-δ constellation, of which configuration parameters can be expressed as N/P/F, i, h [22], and set N = P = 5, F = 3, i = 90°, and h = h 1 , then we can find that results presented in Table 4 are quite close to the proposed walker-δ constellation. Actually, it is easy to verify that the proposed walker-δ constellation, which can also be marked as configuration 5/5/3 [26], is the optimal result for our design by our program. Besides, uniform performance of the walker-δ constellation for the whole global can also consequently be verified in our work.
Values of design metrics are given in Table 5.
From the results presented in Table 5, we can find that the final velocity increment for each satellite is 27.172 m/s, and considering that the constellation is required to turn back to the GOM, the velocity increment for each satellite during the whole process is 54.344 m/s, which is about 86.3% of that of the BLASTOFF. It should be mentioned that the velocity increment is optimized by the designer's preference, and hence, the value can be smaller at the cost of the increment of the reconfiguration time.
Also, the maximum reconfiguration time is about 150 hours, while the mean transition time is only 50 hours. That is to say, for any target site on the ground, the whole constellation can finish the reconfiguration in 6 days, and the average transition time for each satellite is about 2 days. Therefore, satellites in the constellation can achieve the orbital transfer gradually within 6 days for any target on the global, and the target site can be overhead monitored during the whole reconfiguration process. Besides, the standard deviation of reconfiguration times for all target sites is minimized to about 22.6 hours, which is much less than the maximum reconfiguration time, and similar reconfiguration ability for all target sites can be achieved.
In addition, a repeating ground track orbit is chosen to be the target orbit of all satellites in our design, which means the disaster site can always be directly overhead monitored by the whole constellation after the reconfiguration. Once the observation mission ends, the constellation can convert back to the GOM in an inverse process.

Case Study.
In this subsection, a practical case is studied to demonstrate the superiority and effectiveness of the designed ReCon. A monitoring mission for a forest fire in Liangshan, Sichuan, China, is introduced. The GOM and ROM of the Recon, as well as the BLASTOFF proposed in literature [3], are all considered for the monitoring mission, and specific observation performances of these constellations are compared.
The mission starting time is 2020.3.30.10:00:00 UTC, and the geographic latitude and longitude of the target site are 28.531°N and 101.241°E. For the ReCon, RAAN and AOL of the anchor satellite at the starting time are both set to 0, and orbital altitude and phase differences of the rest satellites are provided in Table 4. For the BLASTOFF, specific parameters at the starting time can be found in Table 2 of literature [3]. As for the camera on the satellite, specific parameters are set according to the GF-1 [27]. High-resolution modes of the camera are considered for the overhead fire monitoring, and the width of the satellite swath on the ground is 60 km.
Observation sequences of different constellations for the target site are presented in Figure 8. The X-axis represents the time from the mission beginning, and the unit is hour. Each triangle along the X-axis denotes an observation, and the x-value of the vertex of a triangle represents the monitoring moment, since the observation duration is quite short. The figure of GOM shows the observation sequences of the initial constellation, while the figure of ROM presents the observation sequences of the ReCon during and after the reconfiguration. The total number of observations varies with the time from the mission beginning for each constellation, and the specific change trend can be concluded based on the observation sequences, as illustrated in Figure 9.
From Figures 8 and 9, we can find that the total number of observations of the ROM and that of the BLASTOFF are much more than that of the GOM, and hence, observation for the disaster area after reconfiguration is more adequate and frequent. Moreover, the ROM can achieve 4 extra obser-vations within 2 days compared to the BLASTOFF, and thus, observation of the ROM is much more timely, although the total number of observations of the ROM is 2 less than that of the BLASTOFF. Besides, observation of the ROM becomes regular (5 times per day as the BLASTOFF) after about 100 hours, and this regular observation would remain until the mission ends without extra orbital transfer of any satellite in the constellation.
In the foregoing case, each satellite in the ROM monitors the target site during the ascending pass, and the orbit of each satellite is a repeating orbit with the same repeating period. Therefore, all satellites in the ROM are flying along the same repeating ground track, and this constellation is also the socalled Flower constellation [17], which is marked as ROM1 in this work. Actually, each satellite is able to monitor the target site during the descending pass as illustrated in our mathematical mode. If the initial condition of the constellation changes, the final ROM might be different from the ROM1. Here, we set the AOL of the anchor satellite at the starting time as 180°, and the rest conditions stay the same. Then, we can get a new observation mode, which is marked as ROM2. In the ROM1, all satellites observe the target site during ascending pass, and the AOL of each satellite stays the same with that of the anchor one, while in the ROM2, there is one satellite observing the target site during the descending pass, and its AOL is 56.740°behind the others'. Specific configurations of the GOM, ROM1, and ROM2 are presented in Figure 10.
Based on the presented configurations, ground tracks and average observation numbers per day of the GOM, the ROM1, and the ROM2 can be analyzed, as shown in Figure 11, where the red point represents the target site. It should be mentioned that average observation numbers are calculated by the method presented by Razoumny [28,29]. From the figure, we can find that in the GOM, satellites can evenly scan the entire region within a latitude band without bias, while in the ROM, satellites can scan the region near a fixed ground track frequently. Specifically, in the GOM, satellites can observe the target site about 0.25 times a day, which means that the target site could only get one observation in about 4 days, while all areas in the latitude band can be monitored. In the ROM, satellites can observe the target site 5 times a day, and particularly in the ROM1, only the region around the target site can get 5 observations in one day, and most of the areas cannot be monitored at all. As a whole, the ReCon can scan the Earth uniformly in GOM as a normal operation, and when disaster occurs, it can turn to ROM, and satellites can scan the disaster area adequately and frequently.     To further analyze different properties of the GOM and the ROM, the low-resolution mode of the GF-1 camera is considered, of which the swath width is 800 km. Similarly, average observation numbers per day of the GOM, the ROM1, and the ROM2 are calculated, and results are given in Figure 12. With a wider swath of each satellite, the GOM can achieve about 3 times of observations per day for most of areas globally, while the ROM can monitor the target site 10 times per day. For the GOM, observation numbers per day for a specific latitude is still almost the same, while for both ROM1 and ROM2, areas around the target are the main observation object. Therefore, different modes of the ReCon show excellent observation performances for both global observation and disaster monitoring.
In addition, to verify the accuracy of the proposed mathematical model, distances between the target site and the satellite ground tracks in the ROM1 and the ROM2 are presented in Figure 13, which are also the enlarged views of the ground tracks shown in Figures 11(b) and 11(c). It should be noted that the orbital elements and ground tracks of each satellite in the ROM1 and the ROM2 are calculated through numerical simulation, while in the numerical simulation, the transition times and the moments to implement the orbital transfer are given by the proposed mathematical model. From Figure 13, we can find that the minimum arc distance between the target site and the satellite ground tracks is less than 0.01 degrees, which is about a few hundred meters and does not affect the overhead observation. Therefore, the analytical calculations of the transition times and the moments to implement the orbital transfer are accurate enough for the initial observation mission planning.

Conclusion
In this paper, a reconfigurable constellation, which can achieve timely and adequate observation for an unpredictable disaster site through orbital maneuvers, is successfully designed. By analyzing the space geometry relations between an arbitrary satellite orbit and disaster site, a mathematical model for orbital transfer and overhead monitoring is established. Aiming at minimizing both total velocity increment and reconfiguration time, optimization of the initial constellation and the reconfiguration process are combined and carried out, and a physical programming method along with a genetic algorithm is employed for the multiobjective optimization problem. Preferences in physical programming are set by referring to literature [3], and results are compared to demonstrate the superiority of our design. Also, the GOM and ROM of the ReCon as well as the BLASTOFF are analyzed and compared, and excellent observation performances of the ReCon are demonstrated.
As for the results of the optimization, the optimal reconfigurable constellation is extremely close to the walker-δ constellation. The reconfiguration velocity increment for each satellite is 54.344 m/s, which is about 86.3% of that in literature [3]. The maximum reconfiguration time is about 150 hours, while the mean transition time is only 50 hours, which means satellites can finish transition gradually, and the target site can be overhead monitored during the whole reconfiguration process. Also, the standard deviation of reconfiguration times for all target sites is minimized to less than one day to achieve similar reconfiguration ability for all target sites. A specific disaster case is analyzed, and the observation

13
International Journal of Aerospace Engineering preferred optimal results. Future work would be to monitor 2 or more disaster sites, as well as to monitor the moving disaster site through constellation reconfiguration. This research can provide basic idea and method for the preliminary design of reconfigurable constellation for disaster monitoring.

Data Availability
The data used to support the findings of this study are included within the article.

Conflicts of Interest
The authors declare that there is no conflict of interest regarding the publication of this paper.