ANALYTICAL EXPRESSION FOR UNCERTAINTY PROPAGATION OF AERIAL COOPERATIVE NAVIGATION

. In this paper, the propagation of uncertainty in a cooperative navigation algorithm (CNA) for a group of flying robots (FRs) is investigated. Each FR is equipped with an inertial measurement unit (IMU) and range-bearing sensors to measure the relative distance and bearing angles between the agents. In this regard, an extended Kalman filter (EKF) is implemented to estimate the position and rotation angles of all the agents. For further studies, a relaxed analytical performance index through a closed-form solution is derived. Moreover, the effects of the sensors noise covariance and the number of FRs on the growth rate of the position error covariance is investigated. Analytically, it is shown that the covariance of position error in the vehicles equipped with the IMU is proportional to the cube of time. However, the growth rate of the navigation error is, considerably more rapid compared to a mobile robot group. Furthermore, the covariance of position error is independent of the path and noise resulting from the relative position measurements. Further, it merely depends on both the size of the group and noise characteristics of the accelerometers. Lastly, the analytical results are validated through comprehensive Guidance, Navigation, and Control (GNC) in-the-loop simulations.


Introduction
In recent years, autonomous robots have achieved significant impacts on business, industry and people's social and private affairs as well. Needless to say, FRs have made a substantial contribution to field the field of autonomous robotics which have led to significant improvements. In this regard, the IOD, one of the elements of smart cities, is a newly emerging branch of aerial robotics. The IOD, a layered network control architecture, is designed primarily for managing aerial robots' (AR) access to controlled airspace (Chakraborty et al., 2016). ARs can also be regarded as an aerial base station to collect data service from sensors to the consumers, deploying in a smart city (Chakraborty et al., 2019).
By the advent of a new generation of technologies in robotics, a new idea called "Cooperative Robotics" (CR) has recently emerged to achieve complex missions (Enayatollahi & Atashgah, 2018;Kemsaram et al., 2017). This concept was rapidly embraced and developed in the field of AR as well. In this regard, the authors of (Lee, 2020) have proposed a CNA using the technology of IODs. A group of FRs is more robust than a single one and by the employment of a team of inexpensive and miniature FRs, complex missions can be performed (Ryan et al., 2004). One of the most important subsystems of autonomous robots is the navigation subsystem providing the position, velocity, and attitude (PVA) of a vehicle. Moreover, this subsystem feeds guidance and control (GC) systems to guidee the vehicle towards the desired state/trajectory. Consequently, the role of this subsystem, specifically in ARs, is is so vital that its failure can halt the mission.
Most of the current navigation methods are dedicated to the optimal fusion of the measurements of proprioceptive and exteroceptive sensors. The proprioceptive sensors monitor the motion of the vehicle, while the exteroceptive sensors detect the surrounding environment and its signals (Roumeliotis & Bekey, 2000). One of the most ubiquitous proprioceptive sensors in FRs is IMU which conventionally involves three accelerometers and three gyros. By integrating the outputs of the IMU sensors, the three-dimensional (3D) position and attitude of the vehicle are calculated, that's why the INS has an incremental error over time. To resolve this problem, sensors and other technologies with non-dead-reckoning nature are fused with the inertial sensors. One of the most widely used technologies in almost all navigation systems is the GNSS. However, GNSS signals may be unavailable due to a variety of reasons including: the state of the ionosphere, towers, and buildings in urban areas, mountain ranges, or jamming conditions. There are many tools available for estimation and identification of the state variables and parameters of a vehicle by fusion of the measurements and motion data. The family of Kalman filters (KFs) are refered to as the most pervasive real-time estimators in the field of navigation (Simon, 2006).
As mentioned before, sometimes exteroceptive data like GNSS signals are not available, as a result, the mission may fail. To enhance the reliability and robustness in the field of robotics, the concept of cooperation and collaboration was also developed for the navigation subsystems. For the first time, Kurazume et al. (1994) proposed a new localization method called "cooperative positioning". In this model, the agents share their position and measurement data with the other ones to minimize the overall position error by integrating all the data. Later, Roumeliotis and Bekey (2000) presented a new method called "Collective Localization" based on which the localization problem in a group of mobile robots (MR) is solved in a distributed configuration and without any external source like GNSS or map data. Each robot collects data from its motion sensor (proprioceptive sensor) and relative sensor (exteroceptive sensor) then shares this information with the rest of the team. Linear observability analysis indicates that the system is not observable unless at least one of the robots has access to the absolute position (i.e. GNSS). The authors in (Rekleitis et al., 2002) investigated the effects of different observations like relative range, bearing, and orientation, besides the number of agents on the performance of the cooperative localization. Moreover, Rekleitis and  and Roumeliotis and Rekleitis (2003), presented an analytical solution of the upper bound of the increase rate of the uncertainty, as a function of the number of the agents and the odometer error. Furthermore, in Rekleitis and Roumeliotis (2003) all robots are restricted to move along one dimension while in Roumeliotis and Rekleitis (2003) the robots can perform planar motions. For the first time,  presented a nonlinear observability analysis method (Lie algebra) for two MRs, based on polar coordinates. Moreover, it is shown that a relative bearing measurement provides the most observability between the robots. It is worth noting that the method presented in Roumeliotis and Bekey (2000), was developed with the most general relative measurements between two robots. The method used by , is practical for real application; due to the existence of sensors to measure the data.
Succinctly put, the first results indicate that the relative bearing can improve localization accuracy better than other relative measurements. Later, Mourikis et al. (2006) presented a closed-form relation for calculating steady-state covariance of the position error when only one of the robots has access to a piece of absolute posi-tion information. The paper presents a closed-form solution for calculating the growth rate of the covariance of the position error in the absence of absolute position data. For the first time, Sharma and Taylor (2008) used the concept of the CN for a group of miniature aerial vehicles (MAV). Each agent receives all IMU and relative position measurements from adjacent agents and integrates them by using an EKF. In this method, each MAV estimates the PVA of itself and the nearby MAVs. Likewise, by a nonlinear observability analysis, it is shown that, in this method, the biases of the IMU can be observable. In 2009, Sharma presented a distributed CNA for a MAV group using a camera mounted on each agent, which measures only the relative bearing angle (Sharma & Taylor, 2009). The simulation results were presented for two modes: the first one with two known landmarks, and the other without any landmark. Afterward, in Bahr et al. (2009) the authors presented a technique for a consistent CN of multiple autonomous underwater vehicles (AUV's). The paper addresses the problem of overconfidence in the localization of robots while sharing navigation information for the CN. Each AUV accurately traces the source of the measurements and prevents them from being used more than once. Multiple estimates are consistently combined and provide conservative covariance estimates. Gao and Chitre (2010) presented a CNA for two AUVs with different navigational capabilities using range-only measurements; where AUV with higher positioning accuracy plays the role of a beacon for the one with lower accuracy. Subsequently, Sharma et al. (2012) presented a graph-based method to analyze the nonlinear observability of the bearing-only CN. they used graph properties between two nodes to explain the observability of the system. In Sharma et al. (2013), the authors developed a test-bed of three ground robots equipped with the wheel encoders and omnidirectional cameras, to implement the bearing-only cooperative localization. Simulation and experimental results are presented to validate the observability conditions (Sharma et al., 2012) for the complete observability of the bearing-only cooperative localization problem. Furthermore, Wanasinghe et al. (2016) used the split covariance intersection algorithm for decentralized CN. Similarly, Goel et al. (2016) developed an algorithm to estimate the 3D position of a group of UAVs using low-cost sensors. The paper investigated the performance of the algorithm under available and unavailable GNSS conditions in both centralized and decentralized configurations.
At the same time, Chakraborty et al. (2016) presented a centralized cooperative localization for a group of fixedwing UAVs, based on relative measurements between the agents and some landmarks with a known position, in a GNSS denied condition. However, Goel et al. (2017) used a UWB module as a relative sensor between UAVs, in which the estimated state vector includes some IMU parameters, such as gyro and accelerometer biases to improve accuracy. In Sullivan et al. (2018) Sullivan investigated how a CL algorithm performs under different conditions. This paper investigated a CL under varying sensor specifications (position accuracy, yaw accuracy, sample rate), communication rates, and the number of robots for both homogeneous and heterogeneous multirobot systems. In this paper it is founded that yaw accuracy has a substantial effect on CL performance, also a fast communication rate can be detrimental, and heterogeneous systems are greater candidates for cooperative localisation than homogeneous systems. Later, Chakraborty et al. (2019) demonstrated that each FR group can localize itself using only two known landmarks and only relative bearing observations. He used two bearing angles in 3D space for the CNA. The conditions required to maintain observability and perform consistent estimations, using Lie algebra were also presented.
To the best of the authors' knowledge, the analysis of the behavior of error covariance variations in the cooperative navigation for a group of FRs still demands more elaborations; though, the reviews admit this statement. One of the contributions of this work is that in this research a closed-form solution is derived as a function of time to predict the covariance of the position error due to the CN which is applied in a group of FRs. Furthermore, the validation process is conducted via numerical GNC in-the-loop simulations. It is considered that all of the FRs comprise a guidance and control system (GCS) to follow a predefined path, planned via arbitrarily distributed waypoints. It is worthwhile to mention that, the GCS includes a waypoint navigation system (WNS) for path following and a cascade control system as the controller, which performs a convenient apparatus for the applicability, validity, and performance evaluation purposes. Besides, the simulations exploit 6DoF equations of motion along with the GC and CNA systems for each of the FRs. The related toolsets have been accompanied meticulously under the Matlab environment, for the supplementary software-inthe-loop (SIL), hardware-in-the-loop (HIL), and practical filed tests. This paper is organized as follows. Section 1 describes the algorithm of the CN in a group of FRs in which the formulations of the problem are rendered. In Section 2 a closed-form solution is derived to describe the behavior of the covariance of the position errors. Subsequently, the simulation results are presented in Section 3. Lastly, the conclusions and future work are presented in the last section.

Cooperative navigation
The group of FRs here consists of 'N' quadrotors which fly inside a GNSS denied environment and employs an EKF to perform the CN algorithm. Each FR is equipped with an IMU to perform dead-reckoning process (propagation of the state variables) and a relative sensors to measure the inter-agent relative range and bearing angles. In Figure 1, the notations θ and ψ signify pitch and heading angles, while the symbols α and β stand for the relative elevation and relative azimuthangles, respectively. Here, we extend the cooperative localization method introduced in (Roumeliotis & Bekey, 2000), for a group of FRs. In this method, all agents of the group are considered as a unified system. The agents exchange their information including measurements and the estimated position data. Accordingly, the exchange of each exteroceptive measurement leads to an improvement in the estimation of the overall position error. The exteroceptive sensors measure the relative range and bearing angles between the agents.
As stated, an EKF estimator is used to performe the method. As usual, the EKF is divided into the prediction/ propagation and update cycles (Bayat & Atashgah, 2017). In the prediction phase, the knowledge concerning the state variables of the system is propagated to the next step; based on the evolution of the state variable equations, the measured control inputs, and the statistical behavior of the system noise. In the update phase, the relative pose measurements are processed to update the propagated state variables during the prior phase.

Inertial navigation state-space equations
As mentioned before, all agents are quadrotor type. Since the operating conditions are at low speed and range, the non rotating NED reference frame is utilized to resolve the state variables. The state vector of the i th agent consists of 3D position, velocity and attitude vectors (PVA), which is defined as below: where the superscript 'T' indicates the transpose mode of the state vector and subscript 'i' denotes the number of FRs in the group. The notations r i v i and ψ i respectively are the PVA vectors of the i th agent. Subsequently, the navigation equations for a single flying robot in a nonerotating NED frame can be written as: where e indicates the NED frame and b stands for the body frame. Also, b i a , e g , and b i ω respectively stand for the measured specific force vector (the linear acceleration except for gravity) in the body frame, the gravity acceleration vector in the Earth frame, and the measured rotational velocity vector in the body frame. C are determined as follows (Titterton & Weston, 2004): in which e b C indicates the body to Earth frame transformation matrix, and c and s notations signify the mathematical cosine and sine functions. As mentioned earlier, b i a and b i ω are the measured linear specific force and rotational velocity vectors, which are corrupted by a zeromean white Gaussian noise ( a ε and ω ε ) with known covariances: ,

Cooperative EKF prediction
This section describes the prediction state equations established by purturbing the state and inputs variablesin Eqs.
where '^' signifies the estimated parameters. Consequently, Eq. (9) in discrete-time form is given by: represents the i th error state vector, t δ is the time step, i F and i G are the system dynamic and input noise matrices, and i U  indicates the system input due to IMU measurement noise. In the following, i F and i G are given by: Therefore, the system noise covariance is given by: 3 3  2  3 3   3 3  3 3  3 3  2  3 3  3 3  2  3 3 3 3 Subsequently, the error propagation equation for the augmented system including the error dynamic model of all agents can be written as: (16) Consequently, the error covariance can be predicted through the following equation: where ( ) k P t + and ( ) 1 k P t − + are respectively the priorly corrected and the currently updated covariance matrices of the whole system can be defined according to the following relation:

Cooperative EKF update
Assuming that each flying robot is equipped with a UWB module or a camera to measure the relative range and bearing angles; the measurement equation between the th i and th j agents can be written as below: where ij ρ , β and α indicate the relative range, azimuth, and elevation angles between the agents, respectively. Figure 1 depicts the relevant parameters in defining the azimuth and elevation angles between the j th and i th agents of the group. These measurements can be utilized to update the state variables and error covariance matrices; the linearized equation of observation used in the update step is as follows: while ij n  is the measurement zero-mean white noise with the covariance of Also ij H is formed as follows: in which all derivatives are calculated by the following relations:

Closed-form solution on uncertainty propagation
In this section, the performance of the CN algorithm in a group of FR is analyzed and then a closed form solution is derived. A: Prediction Step Uncertainty Propagation.

Flying robot case
As stated before, this work employs EKF as an estimator in the CN method. To evaluate EKF performance, changes in the error covariance matrix over time can be a proper index. As mentioned before, the EKF comprises two primary steps: Prediction and Update, where each section is analyzed separately. For simplicity and coherence with other works, it is assumed that all agents perform a planar motion. Consequently, the state vector of each agent can be redefined as: in the following, the Riccati equation is used to evaluate the error covariance in each section (Simon, 2006): where i A denotes the system dynamic matrix of the th i agent in planar motion, which is defined as: where e y a and e x a are the agent's acceleration in y-and x-direction at the flat earth NED reference frame. By substituting Eq. (31) into Eq. (30) the differential equation of error covariance can be written in a state space form as below: (32) Hence, the closed-form solution of the above equations in constant acceleration mode is obtained as follows: where all 0 P −− are initial values of the covariance matrix. The third relation in Eq. (33) demonstrates that the heading error covariance grows linearly with time and slope 2 ω σ . In the following, through the numerical simulations, it is shown that, the effect of the acceleration terms on the first two relations of Eq. (33) is negligible compared to other terms. To illustrate this claim, the percentage of error covariance difference in both accelerated and non-accelerated trajectory is defined as a criterion, which can be calculated through the following relation: where a P and Na P are respectively the position error covariance in an accelerated and non-accelerated trajectory. In the following, Figure 2 depicts the numerical simulation results of a quadrotor flying in an accelerated trajectory. In this simulation, the quadrotor flies to pass through a waypoint at 40 meters away in the x-direction. acceleration, position, and position error covariance due to the dead reckoning are noticeable in Figure 2. Subsequently, Figure 3 exhibits the time-history of the performance criteria as stated in Eq. (34). As stated earlier, Figure 3 reveals that the effect of body acceleration on the position error covariance is negligible.
Besides, Table 1 presents the range of standard deviation of IMU sensors noises in different grades including "control", "tactical" and "navigation"; which have been used in different simulaitions to evaluate the position error covariance. It is concluded that the magnitude of the third term, in the first two relations of Eq. (33), dominates the rest. It should be noted that in subsequent simulations, the IMU sensor with a "control" grade is utilized. Assuming that the initial state errors are zero, the relations of Eq. (33) can be approximated as follows: Eq. (35) shows that the position error covariance in the X-and Y-varies in the same way over time; while xx P and yy P are proportional with 3 t and P ψψ is proportional with t . Furthermore, the relations of Eq. (35) are always positive and are related to 2 a σ and 2 ω σ respectively. Therefore, the position error covariance does not depend on the path. This claim is illustrated in Figure 4. The figure depicts the real and dead-reckoning-based 2D position of a quadrotor that follows a predefined path via several waypoints. Similarly, the ellipse of the position error covariance is depicted over the flight time. Based on the presented results, it is concluded that the growth of the covariance ellipse of position error in the vehicles equipped with the IMU does not depend on the path. However, the position error covariance in the ground mobile robots equipped with the encoder, depends on the path and varies with the path direction.

Comparison to the mobile robot
Similar to the previous approach, the position error covariance of a mobile robot, equipped with an encoder, can be solved as follows: sin 3 where X V and y V are components of robot velocity vector in the reference coordinate system, Φ indicates the estimated heading angle, and lastly ω σ and V σ stand for the standard deviations of the rotational and translational velocities of the robot, respectively. Assuming that the initial error state vector is zero, Eq. (36) is rewritten as follows: Figure 5 shows the 2D trajectory of three mobile robots moving randomly. As depicted in figure 5 the position errors grow due to the dead reckoning. In Figure 5, the ellipse of the position error covariance is depicted alongside the trajectory, as well. As expected from the relations of Eq. (37), position error covariance of a ground mobile robot equipped with an encoder depends on the robot trajectory.

Update step uncertainty propagation
Like the prediction/propagation step, the Riccati equation is used to analyze the error covariance in the update step: It is assumed that each agent only applies the data that it receives directly from other agents to update its estimation. Hence, by performing several algebraic operations in Eq. (38), based on the definitions made so far, the following equation is obtained: where ii P is the covariance sub-matrix of the th i agent, ij P is the correlation covariance sub-matrix between the th i agent and th j agent. Hence, the following relation is obtained via Eqs. (26)-(28): Substituting Eq. (41) with Eq. (40), the following relation is derived: The last term is the update term which includes the other agents covariance matrices. In other words, because of each agent's contribution to the position error of the other agents, its error covariance matrix is reduced by the correlation covariance matrix. Assuming that the absolute position of the th j agent is known and it plays the role of a landmark in the group, thereby the correlation covariance matrix becomes zero and Eq. (42) is rewritten as below:

Closed-form solution to update uncertainty
To derive a closed-form solution for position error covariance, the following assumptions are considered: -The group consists of two agents -The agents only move in only one dimension -The cooperative navigation is only based on relative range measurements -Each agent is only equipped with an accelerometer -Initial position error and initial error covariance matrices are zero -The group comprises of homogeneous agents According to the above-mentioned assumptions, the equation of motion of each agent is regarded as follows: Subsequently, Eq. (44) can be rewritten in the statespace form, for each agent, as below: As before, by linearizing Eq. (45), the error state equations for each agent can be obtained as follows: where i A stands for the system sub-matrix, a ε denotes the system noise due to the acceleration measurements, and i G indicates the system noise input sub-matrix. Afterward, the system noise covariance of each agent is obtained as: Hence, the error state equation for a two-agent system can be written as follows: The relative range between agent 1 and agent 2 measured by agent 1 can be modeled as: where D is the cross distance between the two agents and is always constant according to the assumptions mentioned before. This measurement model is used to update the position estimate and the error covariance matrix estimations. By linearizing Eq. (49), the error measurement model is obtained as below: L is defined as a simplification parameter related to agent 1. ν is the measurement noise and involves a white Gaussian probability density function (pdf). The covariance of ν is given as below: where ρ σ is the standard deviation of the relative range measurement noise. Based on the measurements of the agent 1 and using Eq. (38), the covariance matrix of the system, is obtained as: The equations for the covariance sub-matrices in Eq. (52) are separately given as below: where xx P and vv P are the covariance matrices of the position and velocity errors, respectively. For better comprehension and more abstraction, the following parameters, 11 L P = , 12 M P = , and 22 N P = are defined. Therefore, based on these definitions, Eq. (54) is re-shaped as follows: Consequently, according to the previous analysis, the growth of the position error covariance in the vehicles equipped with an IMU, does not depend on the path. Hence, it can be assumed that the error covariance of both agents, yields approximately the same relations as a function of time; which causes the operator L to be equal with the N one. Based on this statement, Eq. (55) can be rewritten as below: Furthermore, in Eq. (56), there are common terms between the two sets of the equations ( L  and M  ). By eliminating these terms, the following relationships are obtained: ( ) In short, based on the assumptions given in this section and according to Eq. (60), the covariance of the position error, resulting from the implementation of the CNA in a group of FRs, involves the following features: it is proportional to the cube of time ( 3 t ), independent of the covariance of measurement noise, and dependent on the covariance of the system noise (accelerometer noise). Also, the covariance of the agents equals the correlation covariance ( 11 12 P P = ), while it is inversely proportional to the number of agents , where N is the number of agents, subscript CN denotes "cooperative navigation", and subscript DR stands for "deadreckoning".

Simulation experiments
In this section, to evaluate the performance of the proposed method, a series of simulations are performed in the Matlab environment. The results of the 6DoF dynamic GNC in-the-loop simulations confirm the theoretical analysis, deliberated in the previous sections. In the conducted simulations, a tri-agent group of FRs is considered to fly in a 20*20*15 (m) space. Employing the GNC system, each of the FRs follows a pre-planned path through arbitrarily chosen waypoints. As mentioned before, in this work, the agents are quadrotors and the simulation environment encompasses a 6DOF dynamic model of them. Each FR is equipped with an IMU in the grade of "control" class and white Gaussian noise (WGN) is used as a measurement noise model of the sensors. Accordingly, the standard deviation of the acceleration measurements noise ( a σ ) is 2 0.012m / s and that of the rotational velocity ( ω σ ) equals 100 deg / hr. Likewise, the exteroceptive sensor measures the range and the two bearing angles, degraded by a WGN with 5 cm ρ σ = and 4 deg α β σ = σ = , while a 1-Hz update rate is considered in the relative measurements between the agents. In the following, Figures 6 and 7 present the actual trajectory (solid line), the dead-reckoning based trajectory (dashed-dotted line), and the CN (dashed line) for the FRs. In Figure 6, the agents have crossed a curved trajectory ,while in Figure 7, they have crossed a hexagon one. In 6 and 7 figures, the left sub-figure depicts the 2D trajectory and the right one illustrates the 3D path at different altitudes (1.5, 2.5, 3.5 (m)).
Since all measurements are relative and there isn't any absolute position data, the CN merely improves the relative geometry between the agents. As expected, if the initial absolute position data is available accurately, the CN can enhance overall position error. In the following, the analytical outcomes obtained in the previous section, are validated by the 6DoF dynamic GNC in-the-loop simulations. In these numerical experiments, the agents employ only the relative range measurements among each other. Observations indicate that the position error covariance is inversely proportional to the number of agents.
Afterward, Figure 8 depicts the position error covariance vs. time for different number of agents, obtained from the EKF covariance matrix, during the GNC inthe-loop simulations. As expected, the position error covariance grows with the cube of time. Correspondingly, as the number of FRs increases, the growth rate of the covariance decreases. To validate the theoretical solution obtained in Eq. (60), we substitute 100 t = (s) in the following relationship: where N is the number of FRs. In the following, Table 2 presents the results pertaining to the numerical experiments. It is worth mentioning that, the values determined by Eq. (61) in Table 2, are the same as those of the curves at 100 t = (s), attainable in Figure 7. Compared with the case of the MR group (Mourikis et al., 2006), as depicted in Figure 8, the covariance of position error grows more rapidly with time. Also despite the analytical results in an MR group , in an FR one equipped with the IMU, the position error covariance does not depend on the noise characteristics of the relative position measurements.
More on the prior observations, the position error covariance in an FR group is independent of the trajectory and grows uniformly among the agents. Actually, in an MR group, it depends on the trajectory of each agent and grows According according to the path. Therefore the non-uniformity property of the position error covariance among the agents of an MR group helps to minimize the navigation errors and their covariance in a CN algorithm. Contrarily, in an FR group, the lack of this property causes navigation errors to diverge much faster. These issues can be resolved based on the incorporation of the dynamic model of the vehicle inside the navigation system (Bayat & Atashgah, 2017;Koifman & Bar-Itzhack, 1999). In the proposed method, all measurements including relative measurements can be used to estimate the model parameters. By this method, an accurate estimate of the model of at least one agent of the group will be sufficient to increase the performance indices of the CN algorithm.

Conclusions
In this paper, a theoretical analysis for the position error covariance of a group of flying robots, equipped with an IMU and relative range and bearing angle sensors. An analytical formula, as the main contribution of the work, was derived to express the covariance of position error as a function of time and the noise characteristics of the IMU. The validations indicate that the derived relation (Eq. (61)) can be used to analyze and also perform trade-  Figure 6. 2D/3D trajectory plots; (a) the 2D, and (b) the 3D curved trajectories; "*" signifies the initial position offs in the preliminary design phase of the CNA for a group of robots, equipped with inertial sensors. However, the simulation results exhibit that the covariance of position error grows rapidly over time and consequently the CNA can't be used in this case for long-time operations without absolute position information of at least one of the agents. Moreover, regarding the group size, smaller groups seem to be more effective in terms of decreasing the growth rate of the covariance and hence have a greater cost-benefit; compared with the larger groups. An interesting result here is that the covariance of position error is approximately independent of the noise characteristics of the relative position measurements, while it depends on the size of the group and the noise covariance of the accelerometers. In GNSS denied condition, to promote the robustness and stability of the CN algorithm in a group of vehicles, equipped with IMU, a model-based estimation method is proposed as future work. In this method, some parameters of the agents' dynamic model, which aren't accurately available, can be estimated properly. Hence a model-based INS can be more accurate than the unaided INS and helps the CN algorithm to achieve higher performance indices.