# Multisensors Cooperative Detection Task Scheduling Algorithm Based on Hybrid Task Decomposition and MBPSO

Authors: Changyun Liu; Xiangke Guo; Zhihui Li; Yingying Wang; Gang Wei
Journal: Mathematical Problems in Engineering
 (2017)

---

## Abstract

A multisensor scheduling algorithm based on the hybrid task decomposition and modified binary particle swarm optimization (MBPSO) is proposed. Firstly, aiming at the complex relationship between sensor resources and tasks, a hybrid task decomposition method is presented, and the resource scheduling problem is decomposed into subtasks; then the sensor resource scheduling problem is changed into the match problem of sensors and subtasks. Secondly, the resource match optimization model based on the sensor resources and tasks is established, which considers several factors, such as the target priority, detecting benefit, handover times, and resource load. Finally, MBPSO algorithm is proposed to solve the match optimization model effectively, which is based on the improved updating means of particle’s velocity and position through the doubt factor and modified Sigmoid function. The experimental results show that the proposed algorithm is better in terms of convergence velocity, searching capability, solution accuracy, and efficiency.


---

## Body

## 1. Introduction

The early warning detection system of the ballistic missile defense system is composed of the space-based satellite, the ground radar, the sea-based radar, and so on. Due to the early warning equipment deployment and sensor detecting capability, multitactical ballistic missile (TBM) emerging randomly cannot be tracked by a single early warning equipment stably and continuously in the entire trajectory of TBM.The goal of the multisensors cooperative detection is to ascertain the detection sequence and the detection time about TBM under the constraint of limited sensor resources and visualization period and to achieve the optimized match of sensor resources and detecting tasks, so that the optimum detecting profit is obtained.Many researchers have done lots of researches on the resource scheduling of the early warning satellites in depth. In [1, 2], the parallel tabu genetic algorithm and improved particle swarm optimization algorithm are used, respectively, to resolve early warning satellites scheduling problems. A variable neighborhood heuristic algorithm for the dynamic task planning of the Low Earth Orbit warning system is studied in [3]. The task decomposition method for early warning satellites scheduling is studied in [4, 5]. The above researches on the early warning satellites task scheduling is focused on the method of numerical solution and task decomposition, while the sensor resources schedule in BMDS should be the combined schedule containing the early warning satellites and ground radars and so on. Aiming at antimissile resource scheduling, the multitarget optimization model including many factors is established in [6], but the task decomposition method is different from the practical problem and the difference between early warning satellites and radars is not reflected in its scheduling model.Sensor resources scheduling is a nonlinear multitarget combination optimization problem in essence. Because of the high velocity of TBM, it is vital to resolve the sensor scheduling problem in real time. Particle swarm optimization (PSO) [7–9] is easy and effective to use, so it has been widely applied in the sensor resources scheduling and the weapon-target allocation (WTA) [10]. Aiming at the optimization problems of the binary representations solution, BPSO is proposed by Kennedy and Eberhart in [11]. And many improved BPSO algorithms aiming at the shortcomings of BPSO are presented in literature [12–15]. Among them, BPSO algorithm is modified in generation mechanism, update mechanism, and selection mechanism for multiobjective resource allocation problem in [13]. And a memetic binary particle swarm optimization algorithm is introduced in [15] to improve the local and global search performance.Multisensors task scheduling is not only a nonlinear multitarget combination optimization problem, but also a real-time algorithm, so a multisensors task scheduling algorithm based on hybrid task decomposition and modified binary particle swarm optimization (MBPSO) is proposed. Aiming at the complex relationship between sensor resources and tasks, a hybrid task decomposition method is presented, so that the resource scheduling problem is decomposed into subtasks, and then the sensor resource scheduling problem is changed into the match problem of sensors and subtasks. Through the analysis of subtasks, the resource match optimization model based on the sensor resources and tasks is established, which contains several factors, such as the target priority, detecting benefit, handover times, and resource load. For the purpose of global optimization and real-time calculation, a match optimization model for MBPSO algorithm is proposed, which is based on the improved updating means of particle’s through the doubt factor and modified Sigmoid function, multisensors cooperative detecting tasks optimization, and real-time solving are realized.


## 2. Cooperative Detecting Task Decomposition Based on Hybrid Task Decomposition

Antimissile multisensors cooperative detecting task has the following characteristics:(1)
There are many different kinds of sensor resources in the ballistic missile defense system, so the visualization relationship between sensors and targets is complex when there are many TBM.
(2)
The flight range of the ballistic target is large, and it is easy to be out of the visual field of one sensor and enters into the visual field of another sensor.
(3)
The ballistic target has a low possibility to maneuver in the midcourse flight, so its motion trajectory can be predicted precisely.
(4)
The emergence time, position, number, and type of ballistic targets are difficult to predict, so the cooperative detecting tasks are random and dynamic.
(5)
The threat degree of multiple ballistic targets is different, so the priority of cooperative detection tasks is also different.It can be concluded from the above task characteristics that one single cooperative detection task cannot be involved independently in the scheduling, the cooperative detection task should be decomposed into subtasks, and these subtasks can be used in scheduling.


### 2.1. Task Decomposition’s Requirements

Through the task decomposition, the multisensors cooperative detecting tasks are decomposed into many basic and independent subtasks, which cannot be decomposed any longer. Some requirements to subtasks are as follows:(1)
Completeness: cooperative detecting tasks are decomposed into many subtasks, and the sum of subtasks is equal to the cooperative detecting tasks.
(2)
Independence: the subtasks cannot be replaced by each other.
(3)
Granularity: the length of subtasks should be moderate. If the length of the subtasks is too long, it will occupy a sensor resource for a long time, which does not meet the requirement for optimization. If the length of subtasks is too short, it will not meet the requirements of estimation and prediction, which improves the complexity of optimization. Assume that the longest subtask isD
max
⁡ and the shortest is D
min
⁡.
(4)
Interval: a certain interval is allowed among subtasks. If the interval is too long, it increases the trajectory prediction error and is more difficult for sensors to capture the target in the next subtask’ time. So the interval cannot exceed the upper limit value.


### 2.2. Method of Hybrid Task Decomposition

The basic methods [4–6] of task decomposition include “the longest observation time” decomposition method, “evenly split” decomposition method, and “start-stop time” decomposition method. The characteristics of these methods are mainly as follows:(1)
“The longest observation time” decomposition method: a sensor will be occupied for a long time, which not only decreases the detection efficiency of sensor resources but also increases the resource conflict.
(2)
“Evenly split” decomposition method: it is difficult to determine the length of subtasks. If the length of subtasks is too long, some tasks will not be completed by proper sensors, and if the length of subtasks is too short, the handover of sensors will be frequent.
(3)
“Start-stop time” decomposition method: a sensor detects the same target for a long time, which reduces the detection efficiency, and the number of subtasks would increase significantly with the increasing number of targets.Due to the characteristics of the above task decomposition methods, the subtasks in a single task decomposition method would conflict each other, which results in the local optimization of resource scheduling, so the single task decomposition method cannot meet the requirement of task decomposition. Then combining “evenly split” decomposition method with “start-stop time” decomposition method, a hybrid task decomposition method for multisensors detecting multitarget is designed. Its main idea is as follows: firstly, using the “start-stop time” decomposition method, the whole task is decomposed into many subtasks according to the visualization start time and end time of all sensors to all targets; secondly, using “evenly split” decomposition method, the subtasks are selected according to the constraint condition, if the length of subtasks is too long, it will continue to be split again, and if the length of subtasks is too short, it will be integrated. The specific diagram of hybrid task decomposition method is shown in Figure1.Figure 1
Process of hybrid task decomposition.In Figure1, S1 denotes sensor 1, S2 denotes sensor 2, S3 denotes sensor 2, and T1 to T13 denote subtask 1 to subtask 13, respectively. Green-1, orange-2, and blue-3 denote the interval of time; that is, target-1, target-2, and target-3 are in the detecting area of different sensor, respectively. For example, the green-1 in the area of S1 means that the interval of time of the target-1 in the sensor 1 detection area is from zero seconds to forty seconds, which denotes [
0
,
40
]. The hybrid task decomposition is the basis of the match optimization model.The detailed steps of hybrid task decomposition method are as follows.Step 1 (calculating). According to the visual relationship between all sensors and all ballistic targets, the visualization interval of sensors and targets is determined.Step 2 (sorting). The start time and end time in visualization interval of all sensors to all targets are marked, and the marked time is sorted in time from small to large, and then the two adjacent moments are synthesized into an interval; as a result, the whole tasks are decomposed into subtasks.Step 3 (selecting). According to the third constraint condition of task decomposition, the subtasks between D
m
a
x and D
min
⁡ are selected. If the length of subtasks is longer than D
max
⁡, go to Step 4.Step 4 (splitting). If the length of subtasks is longer than D
max
⁡, the subtasks would be split continually, so the length of subtasks is in D
min
⁡
,
D
max
⁡. It is supposed that the length of subtasks isM, the task is split into N subtasks whose length is D
z, the detailed method is as follows:(1)
N
=
M
/
D
min
⁡
if
M
mod
D
min
⁡
≠
0
M
/
D
min
⁡
if
M
mod
D
min
⁡
=
0
D
z
=
M
/
N
,where M
/
D
min
⁡ denotes the largest integer beyond M
/
D
min
⁡.Step 5. If the length of subtasks is less than D
min
⁡, the processing step is as follows: (
1
) if the subtask is in existence independently, namely, the subtask length of the previous, and next subtask exceeds D
min
⁡, the subtask would be deleted; (
2
) if the number of adjacent subtasks is more than two, it is supposed that the start time and the end time of ith subtask are t
start
i and t
end
i, respectively, and the time of the integrated subtask is t
start and t
end, respectively. So the integration method is as follows:(1)
initialize:t
start
=
t
start
1
,
i
=
2;
(2)
ift
end
i
-
t
start
∉
[
D
min
⁡
,
D
max
⁡
], go to (
3
); otherwise, go to (
4
);
(3)
i
=
+
+, go to (
2
);
(4)
k
=
i, the k adjacent subtasks are integrated, and the visual interval of the integrated subtask whose length is less than D
min
⁡ would be deleted, as it be shown in slash part of Figure 1.Step 6 (deleting). After subtasks are processed by Steps 4 and 5, the subtasks whose length cannot meet the third constraint condition is deleted.


## 2.1. Task Decomposition’s Requirements

Through the task decomposition, the multisensors cooperative detecting tasks are decomposed into many basic and independent subtasks, which cannot be decomposed any longer. Some requirements to subtasks are as follows:(1)
Completeness: cooperative detecting tasks are decomposed into many subtasks, and the sum of subtasks is equal to the cooperative detecting tasks.
(2)
Independence: the subtasks cannot be replaced by each other.
(3)
Granularity: the length of subtasks should be moderate. If the length of the subtasks is too long, it will occupy a sensor resource for a long time, which does not meet the requirement for optimization. If the length of subtasks is too short, it will not meet the requirements of estimation and prediction, which improves the complexity of optimization. Assume that the longest subtask isD
max
⁡ and the shortest is D
min
⁡.
(4)
Interval: a certain interval is allowed among subtasks. If the interval is too long, it increases the trajectory prediction error and is more difficult for sensors to capture the target in the next subtask’ time. So the interval cannot exceed the upper limit value.


## 2.2. Method of Hybrid Task Decomposition

The basic methods [4–6] of task decomposition include “the longest observation time” decomposition method, “evenly split” decomposition method, and “start-stop time” decomposition method. The characteristics of these methods are mainly as follows:(1)
“The longest observation time” decomposition method: a sensor will be occupied for a long time, which not only decreases the detection efficiency of sensor resources but also increases the resource conflict.
(2)
“Evenly split” decomposition method: it is difficult to determine the length of subtasks. If the length of subtasks is too long, some tasks will not be completed by proper sensors, and if the length of subtasks is too short, the handover of sensors will be frequent.
(3)
“Start-stop time” decomposition method: a sensor detects the same target for a long time, which reduces the detection efficiency, and the number of subtasks would increase significantly with the increasing number of targets.Due to the characteristics of the above task decomposition methods, the subtasks in a single task decomposition method would conflict each other, which results in the local optimization of resource scheduling, so the single task decomposition method cannot meet the requirement of task decomposition. Then combining “evenly split” decomposition method with “start-stop time” decomposition method, a hybrid task decomposition method for multisensors detecting multitarget is designed. Its main idea is as follows: firstly, using the “start-stop time” decomposition method, the whole task is decomposed into many subtasks according to the visualization start time and end time of all sensors to all targets; secondly, using “evenly split” decomposition method, the subtasks are selected according to the constraint condition, if the length of subtasks is too long, it will continue to be split again, and if the length of subtasks is too short, it will be integrated. The specific diagram of hybrid task decomposition method is shown in Figure1.Figure 1
Process of hybrid task decomposition.In Figure1, S1 denotes sensor 1, S2 denotes sensor 2, S3 denotes sensor 2, and T1 to T13 denote subtask 1 to subtask 13, respectively. Green-1, orange-2, and blue-3 denote the interval of time; that is, target-1, target-2, and target-3 are in the detecting area of different sensor, respectively. For example, the green-1 in the area of S1 means that the interval of time of the target-1 in the sensor 1 detection area is from zero seconds to forty seconds, which denotes [
0
,
40
]. The hybrid task decomposition is the basis of the match optimization model.The detailed steps of hybrid task decomposition method are as follows.Step 1 (calculating). According to the visual relationship between all sensors and all ballistic targets, the visualization interval of sensors and targets is determined.Step 2 (sorting). The start time and end time in visualization interval of all sensors to all targets are marked, and the marked time is sorted in time from small to large, and then the two adjacent moments are synthesized into an interval; as a result, the whole tasks are decomposed into subtasks.Step 3 (selecting). According to the third constraint condition of task decomposition, the subtasks between D
m
a
x and D
min
⁡ are selected. If the length of subtasks is longer than D
max
⁡, go to Step 4.Step 4 (splitting). If the length of subtasks is longer than D
max
⁡, the subtasks would be split continually, so the length of subtasks is in D
min
⁡
,
D
max
⁡. It is supposed that the length of subtasks isM, the task is split into N subtasks whose length is D
z, the detailed method is as follows:(1)
N
=
M
/
D
min
⁡
if
M
mod
D
min
⁡
≠
0
M
/
D
min
⁡
if
M
mod
D
min
⁡
=
0
D
z
=
M
/
N
,where M
/
D
min
⁡ denotes the largest integer beyond M
/
D
min
⁡.Step 5. If the length of subtasks is less than D
min
⁡, the processing step is as follows: (
1
) if the subtask is in existence independently, namely, the subtask length of the previous, and next subtask exceeds D
min
⁡, the subtask would be deleted; (
2
) if the number of adjacent subtasks is more than two, it is supposed that the start time and the end time of ith subtask are t
start
i and t
end
i, respectively, and the time of the integrated subtask is t
start and t
end, respectively. So the integration method is as follows:(1)
initialize:t
start
=
t
start
1
,
i
=
2;
(2)
ift
end
i
-
t
start
∉
[
D
min
⁡
,
D
max
⁡
], go to (
3
); otherwise, go to (
4
);
(3)
i
=
+
+, go to (
2
);
(4)
k
=
i, the k adjacent subtasks are integrated, and the visual interval of the integrated subtask whose length is less than D
min
⁡ would be deleted, as it be shown in slash part of Figure 1.Step 6 (deleting). After subtasks are processed by Steps 4 and 5, the subtasks whose length cannot meet the third constraint condition is deleted.


## 3. Match Optimization Model of Resources and Subtasks

After the task is decomposed, the implementation of multisensor cooperative detecting task is transformed into the implementation of its subtasks. In the interval of each subtask, there are two relationships between sensors and targets:(
1
) whether the target is in the visualization interval or not; (
2
) if the target is in the visualization interval of sensors, whether it is detected by the corresponding sensors or not. The multitargets subtasks optimization of multisensors resource is resolved through establishing the match optimization model of sensor resources and subtasks.


### 3.1. Description of the Problem

In order to facilitate the description of the problem, the following definitions are given.Definition 1 (subtask detection matrixD
k).
The subtask detection matrix denotes the matrix whose component denotes whether the target can be detected by sensor or not in the subtask interval of sensor, one denotes that the target can be detected by sensor, and zero denotes that the target cannot be detected by sensor.Definition 2 (subtask scheduling matrixE
k).
The subtask scheduling matrix denotes the matrix whose component denotes whether the detecting task of sensor is implemented or not in the subtask interval of sensor at the basis of subtask detection matrix, zero denotes that the sensor is to detect target, and zero denotes that the sensor is not to detect target.
It is supposed that, in [T
s,T
e], the number of ballistic targets ism, in addition, the number of LEO satellites is n
1, and the number of ground radars is n
2 in BMDS. After the task is decomposed,TNsubtasks are obtained in time sequence, and thekth subtask is in [ST
k,ET
k], so the subtask detection matrixD
k and subtask scheduling matrixE
k are as follows:(2)
D
k
=
d
i
j
k
n
1
+
n
2
m
,
E
k
=
e
i
j
k
n
1
+
n
2
m
.
In [ST
k,ET
k], if thejth ballistic target is in the visualization interval of the ith sensor, then d
i
j
k
=
1; otherwise d
i
j
k
=
0. In like manner, in [ST
k,ET
k], if the ith ballistic target is detected by thejth sensor, then e
i
j
k
=
1; otherwise e
i
j
k
=
0 and 0
≤
e
i
j
k
≤
d
i
j
k.
In the scheduling interval [T
s,T
e], the whole detection matrix D and the whole scheduling matrix E, respectively, are(3)
D
=
D
1
,
…
,
D
k
,
…
,
D
T
N
(4)
E
=
E
1
,
…
,
E
k
,
…
,
E
T
N
.
So the solution to the whole scheduling problem is to determine the value of each element in the whole scheduling matrixE.


### 3.2. Objective Function

The goal of multisensors cooperative detecting is to detect ballistic target stably and continuously. When the cooperative detecting task is scheduled, not only the property of targets but also the feasibility of detecting task to target is considered. So when the match optimization objective function between sensors and subtasks is established, such factors as the target priority, detecting benefit, handover times, and resource load are considered comprehensively.(1) Target Priority. P
r
j denotes the priority of the jth target and 0
≤
P
r
j
≤
1.(2) Detection Benefit of Sensor to Target [16]. Detection benefit includes two factors: one is the distance of target to sensor, the closer the distance between sensor and the target is, the more advantage the target is detected by sensor; the other is detection angle of sensor to target, there is the optimal detection angle, and the larger the difference between the detection angle and the best detection angle is, the poorer the detection efficiency is. The distance detection benefit A
i
j
k and angle detection benefit P
i
j
k of the ith sensor to thejth target in the interval of the kth subtask are defined, respectively, as follows:(5)
A
i
j
k
=
D
max
⁡
i
-
D
¯
i
j
k
D
max
⁡
i
(6)
P
i
j
k
=
max
⁡
θ
i
opt
-
θ
¯
i
j
k
-
θ
i
opt
-
θ
¯
i
j
k
max
⁡
θ
i
opt
-
θ
¯
i
j
k
,where D
m
a
x
i is the maximum detection distance of the ith sensor, D
-
i
j
k is the average detection distance of the ith sensor to thejth target in the S
T
k
,
E
T
k, θ
opt
i is the best detection angle of the ith sensor, and θ
¯
i
j
k is the average detection angle of the ith sensor to thejth target in the S
Time
k
,
E
Time
k.The normalized distance detection benefitA and normalized angle detection benefit P are as follows:(7)
A
=
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
A
i
j
k
·
e
i
j
k
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
d
i
j
k
P
=
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
P
i
j
k
·
e
i
j
k
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
d
i
j
k
,where 0
≤
A
≤
1, 0
≤
P
≤
1.(3) Handover Rate of Sensor. In the scheduling interval, the larger the number of total handover is, the more frequent the handover is, which would cause the loss of target more easily. The whole resources scheduling matrix (namely,E)’s corresponding element of adjacent subtask is calculated by “XOR,” and then the number of “1” is added up as the total number of handover. The handover rate in the scheduling interval is defined as(8)
B
=
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
-
1
e
i
j
k
⊕
e
i
j
k
+
1
n
1
+
n
2
·
m
·
T
N
-
1
,where 0
≤
B
≤
1.(4) The Total Load Rate of Resource [17]. The total load rate of sensor resource is defined as the ratio of total detection time and total visualization time of sensor to target in the Time
s
,
Time
e, which is marked asL.(9)
L
=
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
E
T
k
-
S
T
k
e
i
j
k
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
E
T
k
-
S
T
k
d
i
j
k
,where 0
≤
L
≤
1; the smaller L is, the lighter the load of sensor resource is.It is known from the above analysis that the matching optimization model of resources and subtasks contains such several factors asPr
j,A,P,B, and L
.Then the objective function is established through linear weighting of these factors:(10)
Z
=
α
1
P
r
j
+
α
2
A
+
α
3
P
-
α
4
B
-
α
5
L
,where α
1
,
α
2
,
α
3
,
α
4, and α
5 are the weight coefficient, which are adjusted adaptively, and α
1
+
α
2
+
α
3
+
α
4
+
α
5
=
1.


### 3.3. Constraint Condition Analysis

The factors of matching optimization objective function in (10) include the target, sensor resource, and subtask. In the course of optimization solution of the objective function, it must be subject to the corresponding constraint conditions. So analysis of constraint conditions is the prerequisite for the optimization solution of objective function.(1)
The tracking target capacity constraint of sensor.The tracking target capacity of sensors is(11)
∑
j
=
1
m
e
i
j
k
≤
R
i
,
∀
i
,
k
,where R
i denotes the tracking target capacity of the ith sensor.(2)
Binary detection constraint of LEO satellites.The TBM must be detected by two satellites at the same time; namely,(12)
∑
i
=
1
n
1
e
i
j
k
=
0,2
,
4
,
…
,
∀
j
,
k
.
(3)
The number constraint of sensor.(13)
r
j
s
1
≤
∑
i
=
1
n
1
e
i
j
k
≤
r
j
e
1
∩
r
j
s
2
≤
∑
i
=
n
1
n
1
+
n
2
e
i
j
k
≤
r
j
e
2
,
∀
j
,
k
,where r
j
s
1 and r
j
e
1 are the lower limit and upper limit of the satellites, respectively, and r
j
s
2 and r
j
e
2 are the lower limit and upper limit of the radars, respectively.(4)
Feasibility constraint.(14)
0
≤
e
i
j
k
≤
d
i
j
k
,
∀
i
,
j
,
k
.
(5)
Solution space constraint.(15)
d
i
j
k
=
0,1
e
i
j
k
=
0,1
.


## 3.1. Description of the Problem

In order to facilitate the description of the problem, the following definitions are given.Definition 1 (subtask detection matrixD
k).
The subtask detection matrix denotes the matrix whose component denotes whether the target can be detected by sensor or not in the subtask interval of sensor, one denotes that the target can be detected by sensor, and zero denotes that the target cannot be detected by sensor.Definition 2 (subtask scheduling matrixE
k).
The subtask scheduling matrix denotes the matrix whose component denotes whether the detecting task of sensor is implemented or not in the subtask interval of sensor at the basis of subtask detection matrix, zero denotes that the sensor is to detect target, and zero denotes that the sensor is not to detect target.
It is supposed that, in [T
s,T
e], the number of ballistic targets ism, in addition, the number of LEO satellites is n
1, and the number of ground radars is n
2 in BMDS. After the task is decomposed,TNsubtasks are obtained in time sequence, and thekth subtask is in [ST
k,ET
k], so the subtask detection matrixD
k and subtask scheduling matrixE
k are as follows:(2)
D
k
=
d
i
j
k
n
1
+
n
2
m
,
E
k
=
e
i
j
k
n
1
+
n
2
m
.
In [ST
k,ET
k], if thejth ballistic target is in the visualization interval of the ith sensor, then d
i
j
k
=
1; otherwise d
i
j
k
=
0. In like manner, in [ST
k,ET
k], if the ith ballistic target is detected by thejth sensor, then e
i
j
k
=
1; otherwise e
i
j
k
=
0 and 0
≤
e
i
j
k
≤
d
i
j
k.
In the scheduling interval [T
s,T
e], the whole detection matrix D and the whole scheduling matrix E, respectively, are(3)
D
=
D
1
,
…
,
D
k
,
…
,
D
T
N
(4)
E
=
E
1
,
…
,
E
k
,
…
,
E
T
N
.
So the solution to the whole scheduling problem is to determine the value of each element in the whole scheduling matrixE.


## 3.2. Objective Function

The goal of multisensors cooperative detecting is to detect ballistic target stably and continuously. When the cooperative detecting task is scheduled, not only the property of targets but also the feasibility of detecting task to target is considered. So when the match optimization objective function between sensors and subtasks is established, such factors as the target priority, detecting benefit, handover times, and resource load are considered comprehensively.(1) Target Priority. P
r
j denotes the priority of the jth target and 0
≤
P
r
j
≤
1.(2) Detection Benefit of Sensor to Target [16]. Detection benefit includes two factors: one is the distance of target to sensor, the closer the distance between sensor and the target is, the more advantage the target is detected by sensor; the other is detection angle of sensor to target, there is the optimal detection angle, and the larger the difference between the detection angle and the best detection angle is, the poorer the detection efficiency is. The distance detection benefit A
i
j
k and angle detection benefit P
i
j
k of the ith sensor to thejth target in the interval of the kth subtask are defined, respectively, as follows:(5)
A
i
j
k
=
D
max
⁡
i
-
D
¯
i
j
k
D
max
⁡
i
(6)
P
i
j
k
=
max
⁡
θ
i
opt
-
θ
¯
i
j
k
-
θ
i
opt
-
θ
¯
i
j
k
max
⁡
θ
i
opt
-
θ
¯
i
j
k
,where D
m
a
x
i is the maximum detection distance of the ith sensor, D
-
i
j
k is the average detection distance of the ith sensor to thejth target in the S
T
k
,
E
T
k, θ
opt
i is the best detection angle of the ith sensor, and θ
¯
i
j
k is the average detection angle of the ith sensor to thejth target in the S
Time
k
,
E
Time
k.The normalized distance detection benefitA and normalized angle detection benefit P are as follows:(7)
A
=
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
A
i
j
k
·
e
i
j
k
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
d
i
j
k
P
=
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
P
i
j
k
·
e
i
j
k
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
d
i
j
k
,where 0
≤
A
≤
1, 0
≤
P
≤
1.(3) Handover Rate of Sensor. In the scheduling interval, the larger the number of total handover is, the more frequent the handover is, which would cause the loss of target more easily. The whole resources scheduling matrix (namely,E)’s corresponding element of adjacent subtask is calculated by “XOR,” and then the number of “1” is added up as the total number of handover. The handover rate in the scheduling interval is defined as(8)
B
=
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
-
1
e
i
j
k
⊕
e
i
j
k
+
1
n
1
+
n
2
·
m
·
T
N
-
1
,where 0
≤
B
≤
1.(4) The Total Load Rate of Resource [17]. The total load rate of sensor resource is defined as the ratio of total detection time and total visualization time of sensor to target in the Time
s
,
Time
e, which is marked asL.(9)
L
=
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
E
T
k
-
S
T
k
e
i
j
k
∑
i
=
1
n
1
+
n
2
∑
j
=
1
m
∑
k
=
1
T
N
E
T
k
-
S
T
k
d
i
j
k
,where 0
≤
L
≤
1; the smaller L is, the lighter the load of sensor resource is.It is known from the above analysis that the matching optimization model of resources and subtasks contains such several factors asPr
j,A,P,B, and L
.Then the objective function is established through linear weighting of these factors:(10)
Z
=
α
1
P
r
j
+
α
2
A
+
α
3
P
-
α
4
B
-
α
5
L
,where α
1
,
α
2
,
α
3
,
α
4, and α
5 are the weight coefficient, which are adjusted adaptively, and α
1
+
α
2
+
α
3
+
α
4
+
α
5
=
1.


## 3.3. Constraint Condition Analysis

The factors of matching optimization objective function in (10) include the target, sensor resource, and subtask. In the course of optimization solution of the objective function, it must be subject to the corresponding constraint conditions. So analysis of constraint conditions is the prerequisite for the optimization solution of objective function.(1)
The tracking target capacity constraint of sensor.The tracking target capacity of sensors is(11)
∑
j
=
1
m
e
i
j
k
≤
R
i
,
∀
i
,
k
,where R
i denotes the tracking target capacity of the ith sensor.(2)
Binary detection constraint of LEO satellites.The TBM must be detected by two satellites at the same time; namely,(12)
∑
i
=
1
n
1
e
i
j
k
=
0,2
,
4
,
…
,
∀
j
,
k
.
(3)
The number constraint of sensor.(13)
r
j
s
1
≤
∑
i
=
1
n
1
e
i
j
k
≤
r
j
e
1
∩
r
j
s
2
≤
∑
i
=
n
1
n
1
+
n
2
e
i
j
k
≤
r
j
e
2
,
∀
j
,
k
,where r
j
s
1 and r
j
e
1 are the lower limit and upper limit of the satellites, respectively, and r
j
s
2 and r
j
e
2 are the lower limit and upper limit of the radars, respectively.(4)
Feasibility constraint.(14)
0
≤
e
i
j
k
≤
d
i
j
k
,
∀
i
,
j
,
k
.
(5)
Solution space constraint.(15)
d
i
j
k
=
0,1
e
i
j
k
=
0,1
.


## 4. Fast Computation for Resource Matching Model Based on MBPSO

In order to realize fast computation for the match optimization model, two more steps are involved:(1)
For the optimization process of objective function is affected by the above constraints, these five constraints should be processed further to eliminate the corresponding effect.
(2)
A modified binary particle swarm optimization (MBPSO) is devised by modifying the updates on particle’s velocity and position. On one hand, the modified update on particle’s velocity is used to solve the problem that particle is easily trapped into local optimum in BPSO; on the other hand, the modified updates on particle’s position are used to solve the problem that particles in BPSO cannot converge to the optimal solution.


### 4.1. Processing of Constraint Condition

#### 4.1.1. Binary Matrix Encoding

The mapping relationship between solution matrixE and particle’s position X is established according to resource matching model. A 0-1 binary matrix is built by binary encoding matrixE, whose dimension is C
×
T
N
·
m (where C
=
R
1
+
R
2
+
⋯
+
R
n
1
+
n
2). The encoding process is shown in Figure 2.Figure 2
Encoding process of particle.This encoding structure transforms the sensor with multitargets detecting capability into several sensors with the single-target detecting capability. The particle’s position isX
=
x
i
j
C
×
T
N
·
m(if x
i
j is equal to 1, this denotes that thejth target is detected by the ith sensor in the interval of subtask, while if x
i
j is equal to 0, this denotes that thejth target is not detected by sensor ith sensor). The particles created with this encoding method are very intuitive, which are the one-to-one mapping relationship with total scheduling matrixE and denote the scheduling results directly without decoding computation.The analysis shows that this binary matrix encoding will satisfy the constraint conditions (1) and (5) of resources optimization objective function.


#### 4.1.2. Further Processing of Constraint Condition

In order to eliminate the influence of constraint condition (2), (3), and (4) to objective function, these three constraints should be further processed, respectively.(1)
Processing the Constraint (2). The deployment location of satellites is considered, the two satellites whose location are close to each other are integrated as one sensor, and the target capacity of combined satellites is the minimum of the two satellites.
(2)
Processing the Constraint (4). The penalty function method is adopted, which is to check the position where the element is equal to 0 in total detection matrixD, and then the corresponding position’ value in the distance detection benefit matrix and angle detection benefit matrix is set zero, which can eliminate the infeasible solution.
(3)
Processing the Constraint (3). In the course of encoding and updating the particle, constraint (3) is inspected, if it is not conformed to constraint (3), then particle will be modified, and the modified steps are as follows.Step 1. The number of “1” in the column of particle is added up, which is defined asN.Step 2. If N
<
r
j
s
1
+
r
j
s
2, the value of r
j
s
1
+
r
j
s
2
-
N in this column is set to zero randomly, and then Step 4 is run; if N
≥
r
j
s
1
+
r
j
s
2, Step 3 is run.Step 3. If N
>
r
j
e
1
+
r
j
e
2, the value of N
-
r
j
e
1
+
r
j
e
2 in this column is set to zero randomly.Step 4. If all columns of particle are corrected, the step is jumped to Step 1; otherwise, the modified step is completed.


### 4.2. Modified BPSO Algorithm

#### 4.2.1. Modified Method for Particle Velocity Updating Means

The particle velocity updating means in basic BPSO is as follows:(16)
v
i
j
k
+
1
=
ω
v
i
j
k
+
c
1
r
1
p
i
j
k
-
x
i
j
k
+
c
2
r
2
p
g
j
k
-
x
i
j
k
,where k is thekth iteration, ω is the inertia weight, c
1 and c
2 are the cognition learning rate and social learning rate, respectively, r
1and r
2 are random number between 0 and 1, P
i
best
=
p
i
1
,
p
i
2
,
…
,
p
i
D is the optimal value of personal, and P
g
best
=
p
g
1
,
p
g
2
,
…
,
p
g
D is the swarm optimal value.The useful information is not made full use of in the iterative process of basic BPSO; namely, such factors as the particle’s velocity, personal optimal solution, and global optimal solution are contained. The whole swarm would fall into local optima, when global optimal solution is invalid at a certain moment or global optimal solution is trapped into local optima. The personal optimal solution and global optimal solution are doubted through introducing the doubt factor-
ω
k
r
3
p
i
j
k
+
p
g
j
k, which is to reduce the probability of the premature convergence when particle is fallen into local optima. And the modified velocity updating means is as follows:(17)
v
i
j
k
+
1
=
ω
k
v
i
j
k
+
c
1
r
1
p
i
j
k
-
x
i
j
k
+
c
2
r
2
p
g
j
k
-
x
i
j
k
-
ω
k
r
3
p
i
j
k
+
p
g
j
k
,where r
3 is random number between 0 and 1.The modified velocity updating means of the particle are shown in Figure3.Figure 3
Graphical representation of the particle updating.


#### 4.2.2. Modified Method of Particle Position Updating Means

The particle position updating means in basic BPSO are as follows:(18)
x
i
j
k
+
1
=
1
,
a
<
S
v
i
j
k
+
1
0
,
otherwise,where a is random number between 0 and 1 and S
v
i
j
k
+
1 denotes a Sigmoid fuzzy function (where S
v
i
j
k
+
1
=
1
/
1
+
e
x
p
-
v
i
j
k
+
1). v
i
j
k
+
1
∈
v
min
⁡
,
v
max
⁡ denotes the probability of x
i
j
k
+
1 whose value is 1 or 0.Sigmoid function is modified referring methods in [15], and the contrast diagram of the basic Sigmoid function and modified Sigmoid function is shown in Figure 4.Figure 4
The contrast diagram of the basic Sigmoid function and modified Sigmoid function.Analysis of Figure4 shows that, in the Sigmoid function of BPSO, the particle velocity is larger and positive, which means that probability of 1 for the particle position is large; in the same way, the particle velocity is less and negative, which means that probability of 0 for the particle position is large, and the particle velocity is 0, which means that the probability of 0 and 1 for particle position is 0.5. While there is no difference between positive velocity and negative velocity in the modified Sigmoid function, a large absolute value of velocity means that the current particle’s position is not suitable and a great movement is required to make the particle reach the optimal position. Also, a small absolute value of the velocity indicates that the current particle’s position is close to the optimal position.So the modified Sigmoid function is used to update particle’s position in MBPSO algorithm; the modified velocity updating means is(19)
S
v
i
j
k
+
1
=
E
+
1
-
E
×
t
a
n
h
v
i
j
k
+
1
(20)
x
i
j
k
+
1
=
x
-
i
j
k
,
a
<
S
v
i
j
k
+
1
x
i
j
k
,
other,where E
=
erf
⁡
N
F
/
T
′
=
2
/
π
∫
0
N
F
/
T
′
e
-
t
2
d
t is a small value. If the optimal position of particles is not changed in several sequential iterations,Eis added to accelerate the convergence rate and to escape from a local optima.NFis the number of generations that the optimal position is not changed in the sequential iterations, T
′ is a time constant, and x
-
i
j
k is the complement of x
i
j
k.


### 4.3. Algorithm Process

MBPSO algorithm process is as follows.Step 1. Set parameter. Set ω, c
1, c
2, c
3, the largest number of iterations Iter
max
⁡,NF, T
′, population size, particle’s dimension, and so on.Step 2. Initialize the particle swarm, the current personal optimal fitness value, and global optimal fitness value. The position and velocity of initial swarm are generated randomly and k
=
1.Step 3. Judge whether the particle meets constraint (3) or not; if constraint (3) is met, the processing step is gone to Step 4; otherwise, the particle is modified with the modified method.Step 4. The fitness value of each particle is calculated by (10), then P
ibest and P
gbest are updated.Step 5. The velocity of particle is updated by (17), and the position is updated by (19) and (20); then it is judged that whether the particle meets constraint (3) or not, if it meets, the process step is gone to Step 6; otherwise, the particle is modified with the modified method.Step 6. If k
<
Iter
max
⁡, then k
=
k
+
1, and the process step is gone to Step 4; otherwise, the process step is gone to Step 7.Step 7. Output P
g
b
e
s
t. The optimal scheduling matrix scheme is obtained; then the algorithm is ended.


## 4.1. Processing of Constraint Condition

### 4.1.1. Binary Matrix Encoding

The mapping relationship between solution matrixE and particle’s position X is established according to resource matching model. A 0-1 binary matrix is built by binary encoding matrixE, whose dimension is C
×
T
N
·
m (where C
=
R
1
+
R
2
+
⋯
+
R
n
1
+
n
2). The encoding process is shown in Figure 2.Figure 2
Encoding process of particle.This encoding structure transforms the sensor with multitargets detecting capability into several sensors with the single-target detecting capability. The particle’s position isX
=
x
i
j
C
×
T
N
·
m(if x
i
j is equal to 1, this denotes that thejth target is detected by the ith sensor in the interval of subtask, while if x
i
j is equal to 0, this denotes that thejth target is not detected by sensor ith sensor). The particles created with this encoding method are very intuitive, which are the one-to-one mapping relationship with total scheduling matrixE and denote the scheduling results directly without decoding computation.The analysis shows that this binary matrix encoding will satisfy the constraint conditions (1) and (5) of resources optimization objective function.


### 4.1.2. Further Processing of Constraint Condition

In order to eliminate the influence of constraint condition (2), (3), and (4) to objective function, these three constraints should be further processed, respectively.(1)
Processing the Constraint (2). The deployment location of satellites is considered, the two satellites whose location are close to each other are integrated as one sensor, and the target capacity of combined satellites is the minimum of the two satellites.
(2)
Processing the Constraint (4). The penalty function method is adopted, which is to check the position where the element is equal to 0 in total detection matrixD, and then the corresponding position’ value in the distance detection benefit matrix and angle detection benefit matrix is set zero, which can eliminate the infeasible solution.
(3)
Processing the Constraint (3). In the course of encoding and updating the particle, constraint (3) is inspected, if it is not conformed to constraint (3), then particle will be modified, and the modified steps are as follows.Step 1. The number of “1” in the column of particle is added up, which is defined asN.Step 2. If N
<
r
j
s
1
+
r
j
s
2, the value of r
j
s
1
+
r
j
s
2
-
N in this column is set to zero randomly, and then Step 4 is run; if N
≥
r
j
s
1
+
r
j
s
2, Step 3 is run.Step 3. If N
>
r
j
e
1
+
r
j
e
2, the value of N
-
r
j
e
1
+
r
j
e
2 in this column is set to zero randomly.Step 4. If all columns of particle are corrected, the step is jumped to Step 1; otherwise, the modified step is completed.


## 4.1.1. Binary Matrix Encoding

The mapping relationship between solution matrixE and particle’s position X is established according to resource matching model. A 0-1 binary matrix is built by binary encoding matrixE, whose dimension is C
×
T
N
·
m (where C
=
R
1
+
R
2
+
⋯
+
R
n
1
+
n
2). The encoding process is shown in Figure 2.Figure 2
Encoding process of particle.This encoding structure transforms the sensor with multitargets detecting capability into several sensors with the single-target detecting capability. The particle’s position isX
=
x
i
j
C
×
T
N
·
m(if x
i
j is equal to 1, this denotes that thejth target is detected by the ith sensor in the interval of subtask, while if x
i
j is equal to 0, this denotes that thejth target is not detected by sensor ith sensor). The particles created with this encoding method are very intuitive, which are the one-to-one mapping relationship with total scheduling matrixE and denote the scheduling results directly without decoding computation.The analysis shows that this binary matrix encoding will satisfy the constraint conditions (1) and (5) of resources optimization objective function.


## 4.1.2. Further Processing of Constraint Condition

In order to eliminate the influence of constraint condition (2), (3), and (4) to objective function, these three constraints should be further processed, respectively.(1)
Processing the Constraint (2). The deployment location of satellites is considered, the two satellites whose location are close to each other are integrated as one sensor, and the target capacity of combined satellites is the minimum of the two satellites.
(2)
Processing the Constraint (4). The penalty function method is adopted, which is to check the position where the element is equal to 0 in total detection matrixD, and then the corresponding position’ value in the distance detection benefit matrix and angle detection benefit matrix is set zero, which can eliminate the infeasible solution.
(3)
Processing the Constraint (3). In the course of encoding and updating the particle, constraint (3) is inspected, if it is not conformed to constraint (3), then particle will be modified, and the modified steps are as follows.Step 1. The number of “1” in the column of particle is added up, which is defined asN.Step 2. If N
<
r
j
s
1
+
r
j
s
2, the value of r
j
s
1
+
r
j
s
2
-
N in this column is set to zero randomly, and then Step 4 is run; if N
≥
r
j
s
1
+
r
j
s
2, Step 3 is run.Step 3. If N
>
r
j
e
1
+
r
j
e
2, the value of N
-
r
j
e
1
+
r
j
e
2 in this column is set to zero randomly.Step 4. If all columns of particle are corrected, the step is jumped to Step 1; otherwise, the modified step is completed.


## 4.2. Modified BPSO Algorithm

### 4.2.1. Modified Method for Particle Velocity Updating Means

The particle velocity updating means in basic BPSO is as follows:(16)
v
i
j
k
+
1
=
ω
v
i
j
k
+
c
1
r
1
p
i
j
k
-
x
i
j
k
+
c
2
r
2
p
g
j
k
-
x
i
j
k
,where k is thekth iteration, ω is the inertia weight, c
1 and c
2 are the cognition learning rate and social learning rate, respectively, r
1and r
2 are random number between 0 and 1, P
i
best
=
p
i
1
,
p
i
2
,
…
,
p
i
D is the optimal value of personal, and P
g
best
=
p
g
1
,
p
g
2
,
…
,
p
g
D is the swarm optimal value.The useful information is not made full use of in the iterative process of basic BPSO; namely, such factors as the particle’s velocity, personal optimal solution, and global optimal solution are contained. The whole swarm would fall into local optima, when global optimal solution is invalid at a certain moment or global optimal solution is trapped into local optima. The personal optimal solution and global optimal solution are doubted through introducing the doubt factor-
ω
k
r
3
p
i
j
k
+
p
g
j
k, which is to reduce the probability of the premature convergence when particle is fallen into local optima. And the modified velocity updating means is as follows:(17)
v
i
j
k
+
1
=
ω
k
v
i
j
k
+
c
1
r
1
p
i
j
k
-
x
i
j
k
+
c
2
r
2
p
g
j
k
-
x
i
j
k
-
ω
k
r
3
p
i
j
k
+
p
g
j
k
,where r
3 is random number between 0 and 1.The modified velocity updating means of the particle are shown in Figure3.Figure 3
Graphical representation of the particle updating.


### 4.2.2. Modified Method of Particle Position Updating Means

The particle position updating means in basic BPSO are as follows:(18)
x
i
j
k
+
1
=
1
,
a
<
S
v
i
j
k
+
1
0
,
otherwise,where a is random number between 0 and 1 and S
v
i
j
k
+
1 denotes a Sigmoid fuzzy function (where S
v
i
j
k
+
1
=
1
/
1
+
e
x
p
-
v
i
j
k
+
1). v
i
j
k
+
1
∈
v
min
⁡
,
v
max
⁡ denotes the probability of x
i
j
k
+
1 whose value is 1 or 0.Sigmoid function is modified referring methods in [15], and the contrast diagram of the basic Sigmoid function and modified Sigmoid function is shown in Figure 4.Figure 4
The contrast diagram of the basic Sigmoid function and modified Sigmoid function.Analysis of Figure4 shows that, in the Sigmoid function of BPSO, the particle velocity is larger and positive, which means that probability of 1 for the particle position is large; in the same way, the particle velocity is less and negative, which means that probability of 0 for the particle position is large, and the particle velocity is 0, which means that the probability of 0 and 1 for particle position is 0.5. While there is no difference between positive velocity and negative velocity in the modified Sigmoid function, a large absolute value of velocity means that the current particle’s position is not suitable and a great movement is required to make the particle reach the optimal position. Also, a small absolute value of the velocity indicates that the current particle’s position is close to the optimal position.So the modified Sigmoid function is used to update particle’s position in MBPSO algorithm; the modified velocity updating means is(19)
S
v
i
j
k
+
1
=
E
+
1
-
E
×
t
a
n
h
v
i
j
k
+
1
(20)
x
i
j
k
+
1
=
x
-
i
j
k
,
a
<
S
v
i
j
k
+
1
x
i
j
k
,
other,where E
=
erf
⁡
N
F
/
T
′
=
2
/
π
∫
0
N
F
/
T
′
e
-
t
2
d
t is a small value. If the optimal position of particles is not changed in several sequential iterations,Eis added to accelerate the convergence rate and to escape from a local optima.NFis the number of generations that the optimal position is not changed in the sequential iterations, T
′ is a time constant, and x
-
i
j
k is the complement of x
i
j
k.


## 4.2.1. Modified Method for Particle Velocity Updating Means

The particle velocity updating means in basic BPSO is as follows:(16)
v
i
j
k
+
1
=
ω
v
i
j
k
+
c
1
r
1
p
i
j
k
-
x
i
j
k
+
c
2
r
2
p
g
j
k
-
x
i
j
k
,where k is thekth iteration, ω is the inertia weight, c
1 and c
2 are the cognition learning rate and social learning rate, respectively, r
1and r
2 are random number between 0 and 1, P
i
best
=
p
i
1
,
p
i
2
,
…
,
p
i
D is the optimal value of personal, and P
g
best
=
p
g
1
,
p
g
2
,
…
,
p
g
D is the swarm optimal value.The useful information is not made full use of in the iterative process of basic BPSO; namely, such factors as the particle’s velocity, personal optimal solution, and global optimal solution are contained. The whole swarm would fall into local optima, when global optimal solution is invalid at a certain moment or global optimal solution is trapped into local optima. The personal optimal solution and global optimal solution are doubted through introducing the doubt factor-
ω
k
r
3
p
i
j
k
+
p
g
j
k, which is to reduce the probability of the premature convergence when particle is fallen into local optima. And the modified velocity updating means is as follows:(17)
v
i
j
k
+
1
=
ω
k
v
i
j
k
+
c
1
r
1
p
i
j
k
-
x
i
j
k
+
c
2
r
2
p
g
j
k
-
x
i
j
k
-
ω
k
r
3
p
i
j
k
+
p
g
j
k
,where r
3 is random number between 0 and 1.The modified velocity updating means of the particle are shown in Figure3.Figure 3
Graphical representation of the particle updating.


## 4.2.2. Modified Method of Particle Position Updating Means

The particle position updating means in basic BPSO are as follows:(18)
x
i
j
k
+
1
=
1
,
a
<
S
v
i
j
k
+
1
0
,
otherwise,where a is random number between 0 and 1 and S
v
i
j
k
+
1 denotes a Sigmoid fuzzy function (where S
v
i
j
k
+
1
=
1
/
1
+
e
x
p
-
v
i
j
k
+
1). v
i
j
k
+
1
∈
v
min
⁡
,
v
max
⁡ denotes the probability of x
i
j
k
+
1 whose value is 1 or 0.Sigmoid function is modified referring methods in [15], and the contrast diagram of the basic Sigmoid function and modified Sigmoid function is shown in Figure 4.Figure 4
The contrast diagram of the basic Sigmoid function and modified Sigmoid function.Analysis of Figure4 shows that, in the Sigmoid function of BPSO, the particle velocity is larger and positive, which means that probability of 1 for the particle position is large; in the same way, the particle velocity is less and negative, which means that probability of 0 for the particle position is large, and the particle velocity is 0, which means that the probability of 0 and 1 for particle position is 0.5. While there is no difference between positive velocity and negative velocity in the modified Sigmoid function, a large absolute value of velocity means that the current particle’s position is not suitable and a great movement is required to make the particle reach the optimal position. Also, a small absolute value of the velocity indicates that the current particle’s position is close to the optimal position.So the modified Sigmoid function is used to update particle’s position in MBPSO algorithm; the modified velocity updating means is(19)
S
v
i
j
k
+
1
=
E
+
1
-
E
×
t
a
n
h
v
i
j
k
+
1
(20)
x
i
j
k
+
1
=
x
-
i
j
k
,
a
<
S
v
i
j
k
+
1
x
i
j
k
,
other,where E
=
erf
⁡
N
F
/
T
′
=
2
/
π
∫
0
N
F
/
T
′
e
-
t
2
d
t is a small value. If the optimal position of particles is not changed in several sequential iterations,Eis added to accelerate the convergence rate and to escape from a local optima.NFis the number of generations that the optimal position is not changed in the sequential iterations, T
′ is a time constant, and x
-
i
j
k is the complement of x
i
j
k.


## 4.3. Algorithm Process

MBPSO algorithm process is as follows.Step 1. Set parameter. Set ω, c
1, c
2, c
3, the largest number of iterations Iter
max
⁡,NF, T
′, population size, particle’s dimension, and so on.Step 2. Initialize the particle swarm, the current personal optimal fitness value, and global optimal fitness value. The position and velocity of initial swarm are generated randomly and k
=
1.Step 3. Judge whether the particle meets constraint (3) or not; if constraint (3) is met, the processing step is gone to Step 4; otherwise, the particle is modified with the modified method.Step 4. The fitness value of each particle is calculated by (10), then P
ibest and P
gbest are updated.Step 5. The velocity of particle is updated by (17), and the position is updated by (19) and (20); then it is judged that whether the particle meets constraint (3) or not, if it meets, the process step is gone to Step 6; otherwise, the particle is modified with the modified method.Step 6. If k
<
Iter
max
⁡, then k
=
k
+
1, and the process step is gone to Step 4; otherwise, the process step is gone to Step 7.Step 7. Output P
g
b
e
s
t. The optimal scheduling matrix scheme is obtained; then the algorithm is ended.


## 5. Simulation Analysis

### 5.1. Simulation Conditions

In order to validate efficiency of the resources match model and MBPSO algorithm proposed in this paper, simulation scenario is supposed as follows:(1)
Command and control system received early warning information of 4 TBMs at 18:00:30. Relevant information is shown in Table1.
(2)
Two LEO satellites and three ground radars are deployed in early warning system. LEO satellite’s parameters are as follows: orbit altitude is 1600 km, detection range is as far as 1400 km, field size is 20° × 10°, the best detection angle is 90°, target capacity is 2. Ground radar’s parameters are as follows: two works onP-band, one works onX-band, and the detection range of three radars are, respectively, 3000 km, 4000 km, and 2000 km, the area of azimuth-angle is ±60°, the area of pitch-angle is 1~90°, the best detection angle is 45°, 60°, and 80°, respectively, and target capacity of the radar is 3, 4,, and 8, respectively.Table 1
Relevant parameters of TBM.
TBM
Launch time
Flight time (s)
Maximum height (km)
Priority
TBM1
18:00:00
930
551
2
TBM2
18:00:00
1032
615
1
TBM3
18:00:00
879
529
4
TBM4
18:00:00
911
543
3Task scheduling time is from 18:06:30 to 18:08:30, in which sensors will detect targets and provide subsequent intercepting with significant instructions of TBM. The visualization interval of sensor to target is shown in Table2.Table 2
Visualization interval of sensor to TBM.
TBM1
TBM2
TBM3
TBM4
Sat1
[
0,85
]
[
0,102
]
[
0,68
]
[
0,89
]
Sat2
[
0,120
]
[
0,120
]
[
0,120
]
[
0,120
]
PR1
[
0,85
]
[
0,89
]
[
0,72
]
[
0,83
]
PR2
[
40,120
]
[
46,120
]
[
24,120
]
[
29,120
]
XR1
[
97,120
]
[
110,120
]
[
72,120
]
[
83,120
]


### 5.2. Simulation Process

AssumeD
max
⁡ is 20 s, D
min
⁡ is 10 s, andMacLong is 15 s. Nine subtasks are obtained by the task decomposition method proposed in this paper, whose corresponding interval, respectively, is [
0
,
12
], [
12
,
24
], [
29
,
40
], [
46
,
57
], [
57
,
68
], [
72
,
83
], [
83
,
97
], [
97
,
110
], and [
110
,
120
]. Meanwhile, subtasks’ detection matrix is as follows:(21)
D
1
=
1
1
1
1
1
1
1
1
1
1
1
1
0
0
0
0
0
0
0
0
D
2
=
1
1
1
1
1
1
1
1
1
1
1
1
0
0
0
0
0
0
0
0
D
3
=
1
1
1
1
1
1
1
1
1
1
1
1
0
0
1
1
0
0
0
0
D
4
=
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
0
0
0
0
D
5
=
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
0
0
0
0
D
6
=
1
1
0
1
1
1
1
1
1
1
0
1
1
1
1
1
0
0
1
0
D
7
=
0
1
0
0
1
1
1
1
0
0
0
0
1
1
1
1
0
0
1
1
D
8
=
0
0
0
0
1
1
1
1
0
0
0
0
1
1
1
1
1
0
1
1
D
9
=
0
0
0
0
1
1
1
1
0
0
0
0
1
1
1
1
1
1
1
1
.Assume thatIter
max
⁡ = 50,NF= 5, c
1 = c
2 = c
3 = 2.05, T
′ = 20, swarm size is 30, weight coefficient α
1
=
0.3, α
2
=
0.25, α
3
=
0.2, α
4
=
0.15, α
5
=
0.1, r
j
s
1
+
r
j
s
2
=
1, and r
j
e
1
+
r
j
e
2
=
4. According to resources match model in Section 2, the total scheduling matrix is solved by MBPSO algorithm. So subtasks’ scheduling matrix is as follows:(22)
D
1
=
1
1
0
0
1
1
0
0
1
1
0
1
0
0
0
0
0
0
0
0
D
2
=
1
1
0
0
1
1
0
0
1
1
0
1
0
0
0
0
0
0
0
0
D
3
=
1
1
0
0
1
1
0
0
1
1
0
1
0
0
1
1
0
0
0
0
D
4
=
1
1
0
0
1
1
0
0
1
1
0
1
1
1
1
1
0
0
0
0
D
5
=
1
1
0
0
1
1
0
0
1
1
0
1
1
1
1
1
0
0
0
0
D
6
=
1
1
0
0
1
1
0
0
1
1
0
1
1
1
1
1
0
0
1
0
D
7
=
0
1
0
0
0
1
0
0
0
0
0
0
1
1
1
1
0
0
1
1
D
8
=
0
0
0
0
0
0
0
0
0
0
0
0
1
1
1
1
1
0
1
1
D
9
=
0
0
0
0
0
0
0
0
0
0
0
0
1
1
1
1
1
1
1
1
.


### 5.3. Simulation Conclusion

#### 5.3.1. Optimizing Scheduling Scheme

The total scheduling matrixD is transformed, and then the optimal match scheme is obtained, as shown in Figure 5.Figure 5
Optimal scheduling scheme of multisensor for four-TBM detection.Analysis of Figure5 shows the following:(1)
The optimal scheduling scheme meets the constraint of total detection matrix.
(2)
It is guaranteed that there is at least one sensor to detect each TBM, in addition to TBM3 which is not detected in period[
0
,
29
]. Under the constraint of binary satellite detection, the TBM is detected by sensor 1 and sensor 2 for the meantime.
(3)
More sensors are needed to detect TBM2, 4 sensors at the most time, while less sensors are needed to detect TBM3, 2 at the most time. That is, TBM2’s priority is the highest, and TBM3’s priority is the lowest.
(4)
Besides intervals [24
,
29] and [40
,
46], detection of each sensor to TBM is not interrupted, which denotes that detection to TBM is continuous and has fewer handover times.


#### 5.3.2. Algorithm Performance Contrast

The algorithm is verified in MATLAB R2013a on Pentium 3.2 GHz platform. The above match optimization model is solved by BPSO [11], NBPSO [18], MBPSO-2, respectively, and the proposed algorithm, and the performance of each algorithm is compared through calculating 100 times. MBPSO-2 is the algorithm which is obtained by improving the particle’s velocity updating means in BPSO.The optimum fitness value contrast of four kinds of algorithm is shown in Figure6.Figure 6
Contrast of optimal fitness value curves.The performance contrast of four kinds of algorithm is shown in Table3.Table 3
Performance contrast of algorithms.
Algorithm
Average convergence iteration
Average computation time/s
BPSO
16
2.5382
NBPSO
32
3.0759
MBPSO-2
22
2.5968
MBPSO
18
2.6125It can be analyzed from Figure6 and Table 3 that the optimal fitness value of the MBPSO algorithm is better than that of the other three algorithms, and the iteration number of MBPSO is the least, which is 18, while the iteration number of MBPSO-2 and NBPSO is 22 and 32, respectively, which suggest that the real time of the MBPSO algorithm is the best among the BPSO, NBPSO, MBPSO-2, and MBPSO.On one hand, the Sigmoid function in NBPSO is modified, and the optimal fitness value of NBPSO is better than BPSO, but NBPSO still easily traps into local optima due to its velocity updating means. On the other hand, the velocity updating means in MBPSO-2 are modified, and the optimal fitness value of MBPSO-2 is greatly improved, but the particle cannot be converged to the best solution because of the problem of basic Sigmoid function. With both the velocity and position updating means of particle modified in MBPSO, which makes it easier for particle to jump out of local optima, the convergence speed and search capability of MBPSO are better than those of the other three algorithms.


## 5.1. Simulation Conditions

In order to validate efficiency of the resources match model and MBPSO algorithm proposed in this paper, simulation scenario is supposed as follows:(1)
Command and control system received early warning information of 4 TBMs at 18:00:30. Relevant information is shown in Table1.
(2)
Two LEO satellites and three ground radars are deployed in early warning system. LEO satellite’s parameters are as follows: orbit altitude is 1600 km, detection range is as far as 1400 km, field size is 20° × 10°, the best detection angle is 90°, target capacity is 2. Ground radar’s parameters are as follows: two works onP-band, one works onX-band, and the detection range of three radars are, respectively, 3000 km, 4000 km, and 2000 km, the area of azimuth-angle is ±60°, the area of pitch-angle is 1~90°, the best detection angle is 45°, 60°, and 80°, respectively, and target capacity of the radar is 3, 4,, and 8, respectively.Table 1
Relevant parameters of TBM.
TBM
Launch time
Flight time (s)
Maximum height (km)
Priority
TBM1
18:00:00
930
551
2
TBM2
18:00:00
1032
615
1
TBM3
18:00:00
879
529
4
TBM4
18:00:00
911
543
3Task scheduling time is from 18:06:30 to 18:08:30, in which sensors will detect targets and provide subsequent intercepting with significant instructions of TBM. The visualization interval of sensor to target is shown in Table2.Table 2
Visualization interval of sensor to TBM.
TBM1
TBM2
TBM3
TBM4
Sat1
[
0,85
]
[
0,102
]
[
0,68
]
[
0,89
]
Sat2
[
0,120
]
[
0,120
]
[
0,120
]
[
0,120
]
PR1
[
0,85
]
[
0,89
]
[
0,72
]
[
0,83
]
PR2
[
40,120
]
[
46,120
]
[
24,120
]
[
29,120
]
XR1
[
97,120
]
[
110,120
]
[
72,120
]
[
83,120
]


## 5.2. Simulation Process

AssumeD
max
⁡ is 20 s, D
min
⁡ is 10 s, andMacLong is 15 s. Nine subtasks are obtained by the task decomposition method proposed in this paper, whose corresponding interval, respectively, is [
0
,
12
], [
12
,
24
], [
29
,
40
], [
46
,
57
], [
57
,
68
], [
72
,
83
], [
83
,
97
], [
97
,
110
], and [
110
,
120
]. Meanwhile, subtasks’ detection matrix is as follows:(21)
D
1
=
1
1
1
1
1
1
1
1
1
1
1
1
0
0
0
0
0
0
0
0
D
2
=
1
1
1
1
1
1
1
1
1
1
1
1
0
0
0
0
0
0
0
0
D
3
=
1
1
1
1
1
1
1
1
1
1
1
1
0
0
1
1
0
0
0
0
D
4
=
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
0
0
0
0
D
5
=
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
0
0
0
0
D
6
=
1
1
0
1
1
1
1
1
1
1
0
1
1
1
1
1
0
0
1
0
D
7
=
0
1
0
0
1
1
1
1
0
0
0
0
1
1
1
1
0
0
1
1
D
8
=
0
0
0
0
1
1
1
1
0
0
0
0
1
1
1
1
1
0
1
1
D
9
=
0
0
0
0
1
1
1
1
0
0
0
0
1
1
1
1
1
1
1
1
.Assume thatIter
max
⁡ = 50,NF= 5, c
1 = c
2 = c
3 = 2.05, T
′ = 20, swarm size is 30, weight coefficient α
1
=
0.3, α
2
=
0.25, α
3
=
0.2, α
4
=
0.15, α
5
=
0.1, r
j
s
1
+
r
j
s
2
=
1, and r
j
e
1
+
r
j
e
2
=
4. According to resources match model in Section 2, the total scheduling matrix is solved by MBPSO algorithm. So subtasks’ scheduling matrix is as follows:(22)
D
1
=
1
1
0
0
1
1
0
0
1
1
0
1
0
0
0
0
0
0
0
0
D
2
=
1
1
0
0
1
1
0
0
1
1
0
1
0
0
0
0
0
0
0
0
D
3
=
1
1
0
0
1
1
0
0
1
1
0
1
0
0
1
1
0
0
0
0
D
4
=
1
1
0
0
1
1
0
0
1
1
0
1
1
1
1
1
0
0
0
0
D
5
=
1
1
0
0
1
1
0
0
1
1
0
1
1
1
1
1
0
0
0
0
D
6
=
1
1
0
0
1
1
0
0
1
1
0
1
1
1
1
1
0
0
1
0
D
7
=
0
1
0
0
0
1
0
0
0
0
0
0
1
1
1
1
0
0
1
1
D
8
=
0
0
0
0
0
0
0
0
0
0
0
0
1
1
1
1
1
0
1
1
D
9
=
0
0
0
0
0
0
0
0
0
0
0
0
1
1
1
1
1
1
1
1
.


## 5.3. Simulation Conclusion

### 5.3.1. Optimizing Scheduling Scheme

The total scheduling matrixD is transformed, and then the optimal match scheme is obtained, as shown in Figure 5.Figure 5
Optimal scheduling scheme of multisensor for four-TBM detection.Analysis of Figure5 shows the following:(1)
The optimal scheduling scheme meets the constraint of total detection matrix.
(2)
It is guaranteed that there is at least one sensor to detect each TBM, in addition to TBM3 which is not detected in period[
0
,
29
]. Under the constraint of binary satellite detection, the TBM is detected by sensor 1 and sensor 2 for the meantime.
(3)
More sensors are needed to detect TBM2, 4 sensors at the most time, while less sensors are needed to detect TBM3, 2 at the most time. That is, TBM2’s priority is the highest, and TBM3’s priority is the lowest.
(4)
Besides intervals [24
,
29] and [40
,
46], detection of each sensor to TBM is not interrupted, which denotes that detection to TBM is continuous and has fewer handover times.


### 5.3.2. Algorithm Performance Contrast

The algorithm is verified in MATLAB R2013a on Pentium 3.2 GHz platform. The above match optimization model is solved by BPSO [11], NBPSO [18], MBPSO-2, respectively, and the proposed algorithm, and the performance of each algorithm is compared through calculating 100 times. MBPSO-2 is the algorithm which is obtained by improving the particle’s velocity updating means in BPSO.The optimum fitness value contrast of four kinds of algorithm is shown in Figure6.Figure 6
Contrast of optimal fitness value curves.The performance contrast of four kinds of algorithm is shown in Table3.Table 3
Performance contrast of algorithms.
Algorithm
Average convergence iteration
Average computation time/s
BPSO
16
2.5382
NBPSO
32
3.0759
MBPSO-2
22
2.5968
MBPSO
18
2.6125It can be analyzed from Figure6 and Table 3 that the optimal fitness value of the MBPSO algorithm is better than that of the other three algorithms, and the iteration number of MBPSO is the least, which is 18, while the iteration number of MBPSO-2 and NBPSO is 22 and 32, respectively, which suggest that the real time of the MBPSO algorithm is the best among the BPSO, NBPSO, MBPSO-2, and MBPSO.On one hand, the Sigmoid function in NBPSO is modified, and the optimal fitness value of NBPSO is better than BPSO, but NBPSO still easily traps into local optima due to its velocity updating means. On the other hand, the velocity updating means in MBPSO-2 are modified, and the optimal fitness value of MBPSO-2 is greatly improved, but the particle cannot be converged to the best solution because of the problem of basic Sigmoid function. With both the velocity and position updating means of particle modified in MBPSO, which makes it easier for particle to jump out of local optima, the convergence speed and search capability of MBPSO are better than those of the other three algorithms.


## 5.3.1. Optimizing Scheduling Scheme

The total scheduling matrixD is transformed, and then the optimal match scheme is obtained, as shown in Figure 5.Figure 5
Optimal scheduling scheme of multisensor for four-TBM detection.Analysis of Figure5 shows the following:(1)
The optimal scheduling scheme meets the constraint of total detection matrix.
(2)
It is guaranteed that there is at least one sensor to detect each TBM, in addition to TBM3 which is not detected in period[
0
,
29
]. Under the constraint of binary satellite detection, the TBM is detected by sensor 1 and sensor 2 for the meantime.
(3)
More sensors are needed to detect TBM2, 4 sensors at the most time, while less sensors are needed to detect TBM3, 2 at the most time. That is, TBM2’s priority is the highest, and TBM3’s priority is the lowest.
(4)
Besides intervals [24
,
29] and [40
,
46], detection of each sensor to TBM is not interrupted, which denotes that detection to TBM is continuous and has fewer handover times.


## 5.3.2. Algorithm Performance Contrast

The algorithm is verified in MATLAB R2013a on Pentium 3.2 GHz platform. The above match optimization model is solved by BPSO [11], NBPSO [18], MBPSO-2, respectively, and the proposed algorithm, and the performance of each algorithm is compared through calculating 100 times. MBPSO-2 is the algorithm which is obtained by improving the particle’s velocity updating means in BPSO.The optimum fitness value contrast of four kinds of algorithm is shown in Figure6.Figure 6
Contrast of optimal fitness value curves.The performance contrast of four kinds of algorithm is shown in Table3.Table 3
Performance contrast of algorithms.
Algorithm
Average convergence iteration
Average computation time/s
BPSO
16
2.5382
NBPSO
32
3.0759
MBPSO-2
22
2.5968
MBPSO
18
2.6125It can be analyzed from Figure6 and Table 3 that the optimal fitness value of the MBPSO algorithm is better than that of the other three algorithms, and the iteration number of MBPSO is the least, which is 18, while the iteration number of MBPSO-2 and NBPSO is 22 and 32, respectively, which suggest that the real time of the MBPSO algorithm is the best among the BPSO, NBPSO, MBPSO-2, and MBPSO.On one hand, the Sigmoid function in NBPSO is modified, and the optimal fitness value of NBPSO is better than BPSO, but NBPSO still easily traps into local optima due to its velocity updating means. On the other hand, the velocity updating means in MBPSO-2 are modified, and the optimal fitness value of MBPSO-2 is greatly improved, but the particle cannot be converged to the best solution because of the problem of basic Sigmoid function. With both the velocity and position updating means of particle modified in MBPSO, which makes it easier for particle to jump out of local optima, the convergence speed and search capability of MBPSO are better than those of the other three algorithms.


## 6. Conclusions

The resources scheduling problem for the multisensor cooperative detecting the ballistic target is researched in this paper. The multisensor cooperative detecting task is decomposed as many subtasks by hybrid task decomposition method. Then the resource match optimization model based on the resources and tasks is established. Finally, MBPSO algorithm is proposed to solve the match optimization model efficiently. Simulation results show that the proposed algorithm is better in terms of convergence velocity, searching capability, solution accuracy, and efficiency.



---

*Source: 3453589-2017-02-12.xml*