Filtration of UAV Movement Parameters Based on the Received Signal Strength Measurement Sensor Networks in the Presence of Anomalous Measurements of Unknown Power at the Transmitter

Methods based on received signal strength measurements (RSS measurements) are used to determine the unmanned aerial vehicle (UAV) location using a wireless sensor network. The UAV transmitter power is usually unknown. In real conditions, it often becomes necessary to consider existence of anomalous measurement results. The reasons for the violation of the measurement process can be: the influence of interference, errors in the identification of signals during primary processing, failures of the equipment and similar. The optimum and quasi-optimal adaptive algorithms of UAV movement parameters filtration based on the RSS-measurement sensor networks in the presence of anomalous measurements at the unknown power of the transmitter are developed. These algorithms are obtained using Bayes’ theorems and the Markov property of a mixed process, including a vector of target movement parameters and a discrete component characterizing the type of measurement. Analysis of developed algorithm performance was carried out by Monte Carlo method on 2D plane. The quasi-optimal adaptive filtering algorithm detects the appearance of anomalous measurements with probabilities close to unity and allows one to eliminate their influence on the accuracy of UAV movement parameters estimation and also to estimate the UAV unknown transmitter power.


INTRODUCTION
In recent years, small unmanned aerial vehicles (UAVs) (both mini and micro) have become a popular surveillance tool in the field of defense and security. The military also increased interest in small UAVs, which can be used to solve tactical reconnaissance, electronic warfare or to deliver small bombs (Baker 2014). In many countries, continuous research and development work is being carried out in this direction. According to IC Insights, taking into account all potential areas of use, the global sales of devices in 2025 may exceed $ 10 billion (Business Insider Intelligence 2020; Hindle 2017).
On the other hand, the use of UAVs has led to new potential threats to national and public security. Such UAVs can carry explosives, biological or chemical weapons to carry out terrorist acts. Devices can also be used to transport smuggling, drugs, jamming GPS signals or Wi-Fi, which will lead to interruption in communication and data transfer (Wallace and Loffi 2015;Nasr et al. 2019).
The greatest vulnerability of UAVs is due to the presence of electromagnetic radiation. The standard radio frequency UAV bands are ISM 2.4 GHz and ISM 5.8 GHz, in which most commercial Wi-Fi, Bluetooth and internet of things (IoT) systems (i.e., ZigBee, Z-Wave, LoRa) operate. The signals in these bands are freely regulated using free access rules.
Recent advances in wireless sensor networks (WSN) are opening new possibilities in solving the problem of determining radio sources location. Wireless sensor network is a set of miniature and inexpensive devices equipped with various types of sensors, a small microcontroller and a receiver, which are connected via a Wireless Local Area Network (WLAN) and uses radio channels for data transmission. This task has a wide range of applications, such as: rescue operations, autonomous surveillance and monitoring of industrial processes and the environment (monitoring of the animal world), monitoring and control of moving objects, etc.
One of the important features of wireless sensor networks is the ability to track moving objects (targets), including UAVs (Chu and Han 2019;Zhang et al. 2018). To determine the UAV location using a wireless sensor network, methods based on measuring received signal strength (RSS measurements) are the most common. One of the localization methods used in such networks is fingerprinting. In this method, the correspondence of the current signal strength is searched for with available reference values in the calibration database (Zegeye et al. 2016). Alternative to the fingerprinting method are various mathematical algorithms based on the correspondence of the received RSS measurements to the signal propagation model. The accuracy of these methods largely depends on the accuracy of the signal propagation model (Azmi et al. 2018;Uluskan et al. 2017;Li and Huang 2016;Song and Yu 2008;Zou et al. 2016). However, today there is not even a single standard that defines the types of UAVs and their application. Therefore, the UAV transmitter power is, as a rule, unknown. This makes it difficult to use signal propagation model-based techniques to locate an unknown UAV.
When solving the problem of estimating UAV movement parameters (the coordinates of position, velocity, acceleration), it is usually assumed that the measurement results contain information about the process being evaluated with some errors. Of considerable practical interest are cases in which measurements may appear that either do not contain useful information at all or contain results with significantly larger errors. Such measurements are called abnormal (rough), as opposed to the normal measurement, formed during the usual mechanism of results formation (Kupriyanov and Sakharov 2007). The reasons for the violation of the measurement process can be: the influence of interference, errors in the identification of signals during primary processing, failures of the equipment etc. The appearance of anomalous measurements in the case blinds the observer temporarily, leading to a significant decrease in the accuracy of determining the UAV location and movement parameters (Sirota and Kirsanov 2006).
The article is devoted to the development of an adaptive estimation algorithm of UAV movement parameters based on RSS measurements sensor networks at unknown transmitter power and the presence of anomalous measurements.

PROBLEM DEFINITION
The UAV movement in the rectangular coordinate system can be described by using a stochastic dynamic system with random structure having the form (Tovkach and Zhuk 2019): (1) where is the state vector including parameters of UAV movement along coordinates of the rectangular coordinate system; , are known matrices: (2) (3) a is the root mean square (RMS) of random fluctuations of the rate of change of the UAV acceleration; T is the period of receipt of measurements from network sensors; ω(k) T = (ω x (k), ω y (k), ω z (k)) is the uncorrelated sequence of Gaussian vectors with zero expected value and identity correlation matrix Q ω = I; k -is a discrete time instant. Wireless sensor network has S sensors. For determination of UAV location on the 2D plane, the wireless sensor network has to consist of at least three sensors. The use of the RSS method must take into account the dependence of the received signal strength on the distance between the sensor and UAV. Wide range of applications use the model of signal forward propagation, which takes into account only its attenuation (Tovkach and Zhuk 2017a). In this case, the observation/measurement equation describing the measurement of UAV coordinates has the form: where are the measured UAV signal strength by s-th WSN sensors, ; a js (k)υ s (k) are the errors of UAV signal strength measurement by s-th WSN sensors; υ s (k) are the uncorrelated sequence of Gaussian values with zero expected value and dispersion σ 2 , which are independent.
Switching variables determine the type of measurement error. The first type of measurement error is normal the second is abnormal a 2s (k)=γ, s= 0 -,S, where γ determines how many times the RMS of normal measurement error σ increases. Switching variables are independent and take values with probabilities . Nonlinear function in Eq. 6 is described by the expression: where are UAV position coordinates; are coordinates of the s-th WSN sensors position, ; a is the signal attenuation coefficient (close to 2); p 0 is the signal strength at a radial distance r 0 to the sensor, which is an unknown parameter (Tomic et al. 2015). At unknown transmitter power, the parameter p 0 is unknown.

Extending of the model state vector
For the decision of the formulated task, the bayesian method of the adaptive filtering is applicable (Aoki 1971), whose main idea is to include unknown parameters in the state vector of the process. Following this method, the state vector of the process subject to filtering in the considered task can be extended to the form: It includes a vector and UAV transmitter power (Tovkach et al. 2018). In this case, the model of the dynamical system (1) is described by the expression: (9) in which the extended matrices F and G are defined as: The equation of the RSS measurements eq. 6 does not change and the state vector is replaced by the extended state vector . Equation 6 described the observation of a single sensor. Combine the switching variables into a vector Γ(k) . Also combine Eq. 6 for all sensors into a single equation: where is joint measurement vector; is joint vector-valued measurement function; v(Γ(k)) = [a j1 (k)υ 1 (k),...,a jS (k)υ S (k)] T is joint vector of measurement errors with a zero mean vector and covariance matrix σ 2 ,...,a 2 js (k)σ 2 ).
A posteriori probability density function of a mixed Markov process In Zhuk (1989), it was shown that an extended process involving a continuous component and discrete components is a mixed Markov process in discrete time. The a posteriori probability density function (p.d.f.) of a mixed Markov process is described by the expression , where . Based on the Markov property of the mixed process, Bayes' theorem and following the procedure given in Zhuk et al. (2019) p.d.f. is recurrently calculated by the equation: (11) where is the conditional p.d.f. defined by the eq. (5) [N(*,*) -the notation for the Gaussian probability density with the corresponding mathematical expectation and correlation matrix]; is one-step likelihood function; -conditional p.d.f. defined based on the Eq. 6 and 7; is the conditional p.d.f., executing role of a normalizing multiplier. Equation 11 describes an algorithm for nonlinear cofiltering of the state vector and discrete components . Using the probability multiplication theorem, Eq. 11 can be represented as: where is the conditional predicted p.d.f. of vector ; is the conditional a posteriori p.d.f. of the vector on condition Γ(k); is unconditional a posteriori p.d.f. of the vector , defined on k-th step; are predicted and a posteriori probabilities discrete components; is the conditional p.d.f., determined by the formula: is determined by the formula: The initial conditions for algorithm Eqs. 12 to 16 have the form . By means of the Eqs. 13, 14 and 16 filtering of continuous component and Eqs. 12 and 15, filtration of discrete components is executed. An inseparable connection of filtration equations for discrete and continuous components is a specific feature of this algorithm.
The optimal device realizing algorithm Eqs. 12 to 16, is multichannel with the number of channels 2 S . Each channel is coordinated with a certain combination of parameters of the switching variables Γ(k). In this case, the unconditional a posteriori p.d.f. of the vector is determined by averaging the conditional a posteriori p.d.f. using a posteriori probabilities , determined in each channel of the optimal device.

DEVELOPMENT OF A QUASI-OPTIMAL ESTIMATION ALGORITHM
Large computational costs related to the need of integrating multidimensional p.d.f. and also the considerable number of channels in case of big S are a disadvantage of the optimal algorithm. A quasi-optimal adaptive filtering algorithm can be obtained by linearizing Eq. 6 and a Gaussian approximation of the a posteriori p.d.f.
. At the same time, predicted p.d.f. is also Gaussian and its expected value and covariance matrix are calculated by the formulas (Zhuk 1989): where , are expected value and correlation matrix of a posteriori p.d.f . Determination of a posteriori p.d.f.
(13) comes down to calculating its expected value and covariance matrix (Evlanov and Zhuk 1990) using notation : In synthesis of the quasi-optimal algorithm, a sequential method of processing incoming data is used ). The procedure is recurrent when measurements are received from the s-th sensor. In the first step, the expected value and correlation matrix of conditional a posteriori p.d.f. are calculated, provided that the value of the switching variable is known precisely. These equations coincide with the well-known equations of the extended Kalman filter: where -quasi-optimal filter channel gain; is the dispersion of a discrepancy determined by a formula: Estimated û s (k) and the correlation matrix of the estimation error P s (k) are the mathematical expectation and the correlation matrix of the Gaussian approximation of the a posteriori p.d.f.
. This allows to obtain a quasi-optimal device with a number of channels equal to two.
Initial conditions for the procedure Eqs. 17 to 24 at s = 0 have the form and, respectively, resultant values in case at s = S have the form . In case of implementation of the procedure Eqs. 19 to 24 on each l-th step, Gaussian approximation of a posteriori p.d.f. is executed and only the moments remain , .

ANALYSIS OF PERFORMANCE OF THE ALGORITHM
Analysis of algorithm performance Eqs. 17 to 24 was carried out by Monte Carlo method on 2D-plane. The wireless sensor network, represented in Fig. 1, consists of eight sensors located at the following:  For clarity of the algorithm, a test UAV motion trajectory (Fig. 1) was formed. It contains steps. The trajectory consists of five sections (Tovkach and Zhuk 2017b). Abnormal measurement errors occur in the following time instants: for the first sensor k = 14, 85, 87; for the second sensor k = 26, 27, 28; for the sixth sensor k = 55; for the seventh sensor k = 105. Error RMS of measurement is σ dBm, the rate of receipt of measurements from network sensors Т = 1 s, number of realizations is L = 100. The simulation was conducted with The probability of appearance of anomalous sensor measurements at any sensor is equal to and parameter γ=6. Using the Monte Carlo method means, , and RMS errors , of estimation UAV locations were calculated by the formulas: where , are estimations of the X,Y coordinates on the k-th step in the l-th test (trajectory); , are true coordinates; ε , ε are X,Y coordinates estimation errors at the k-th step in the l-th test (trajectory); L is the number of tests. Figure 2 shows RMS (curve 3) errors of estimation , of the coordinates X,Y, obtained when algorithm Eqs. 17 to 24 is executed. Also Fig. 2 shows mean (curve 1) , and RMS (curve 2) , errors of estimation, obtained by statistical simulation. Theoretical and actual RMS errors estimation are close, which indicates the correct operation of algorithm. Also, Fig. 2 shows RMS (curve 4) , errors of estimation when using the Kalman filter obtained by statistical modeling. Figure 3 shows the relationships of the probabilities of detecting anomalous errors of measurements for all sensors. The adaptive filter provides a high probability of recognition of abnormal measurements.   (Tovkach et al. 2018). After completing the adaptation process, occupying 30 steps, the developed filter allows to determine the location of the UAV with the accuracy characteristics close to an algorithm, in which the known power of the transmitter is used. As follows from Fig. 2, the appearance of anomalous measurement errors leads to a significant increase in the estimation errors of the Kalman filter. The application of the developed quasi-optimal adaptive algorithm of trajectory filtering allows one to eliminate their influence on the accuracy of estimating UAV movement parameters in comparison with the Kalman filtering algorithm, as shown in Fig.  2. When abnormal measurement errors appear, the RMS error of the estimation increases by 2 -5 times when using the Kalman filter.

CONCLUSIONS
To identify the unknown transmitter power, the method of expanding the state vector is applied, according to which the unknown parameter is included in the state vector of the filtered process.
Using the optimal adaptive filtering algorithm (Eqs. 12 to 16), the calculation of the joint a posteriori p.d.f. state vector and switching variables is performed. It is recurrent, due to the Markov property of the extended mixed process. The optimal filter is multichannel. Each channel is matched to a specific combination of switching variable parameters. The disadvantage of the optimal algorithm is the large computational costs associated with the need to integrate multidimensional p.d.f., as well as a significant number of filter channels.
The quasi-optimal adaptive filtering algorithm (Eqs. 17 to 24) was obtained by linearizing the measurement in Eq. 6 and Gaussian approximation of the a posteriori p.d.f. state vector. It performs sequential processing of incoming measurements. The number of channels of the quasi-optimal device is equal to two. The inclusion of unknown transmitter power in the state vector allows simultaneously with solving the main task of estimating UAV movement parameters, also to estimate its power, which can be useful for UAVs identification.
The quasi-optimal adaptive filtering algorithm detects the appearance of anomalous measurements with probabilities close to unity and allows one to eliminate their influence on the accuracy of UAV movement parameters estimation in comparison with the Kalman filtering algorithm. After completing the adaptation process of RMS (curve 3) of UAV, transmitter power error is 0.3 dBm.
A further area of research is the development of trajectory filtering algorithms for maneuvering UAVs with the integrated use of data from the sensor network, using the time difference of arrival (TDOA), angle of arrival (AOA) and RSS methods.

DATA AVAILABILITY STATEMENT
All data sets were generated or analyzed in the current study.