AbstractPrecise navigation, as a method for guiding vehicles from one point to another, is an important subject these days especially in navigation of aircraft. Global navigation satellite systems (GNSSs) are capable tools for such a purpose. Any intentional or unintentional interference in satellite signals may cause risks of deadly accidents. Therefore, it is tremendously important to control airports or harbors and locate any existing radio frequency interference device. This localization can be done based on measuring time of arrival (TOA), angle of arrival (AOA), or time difference of arrival (TDOA) of signals from the device to sensors or receivers at some basepoints. In this article, a method is proposed based on these arrivals for optimizing the configuration created by these basepoints from a large grid of points covering a control area. Furthermore, a simulation test was performed to verify the theory, and after that a control network was designed and optimized for the international Landvetter Airport of Sweden. Our simulation studies show that when the AOA is used, our optimization is more robust with respect to the control grid resolution. In addition, optimization based on the TDOA improves the coverage over the control area with a significant reduction of error of control points, but because of the special geometric shape of the Landvetter Airport, such an optimization was not successful.IntroductionToday, almost all vehicles at the surface of the Earth are navigated by global navigation satellite systems (GNSSs). Electromagnetic signals, carrying satellite coordinates, are transmitted from satellites constantly toward all GNSS receivers mounted on vehicles. However, navigation of aircrafts needs more attention because any radio frequency interference in their navigation signal may lead to catastrophes for crew, passengers, and even people on the ground. Taking off and landing are the two most risky times of a flight and are close to an airport. Weather conditions, connection with the airport traffic control tower, and navigation signal are important factors for a successful landing or takeoff. Any intentional or unintentional interference in aircraft navigation and connection with a control tower might lead to serious problems and risk human life. Therefore, some sensors or receivers with the possibility of providing information about interference are needed. The geometric configuration of the basepoints, at which these sensors and receivers are placed, plays an important role in successful localization of an interference device. In this article, a method is developed and applied for optimization of such configuration around an area in such a way that any interference device can be localized with a higher precision.Jamming and spoofing are two well-known types of signal interference. Jamming means transmitting a radio frequency signal into the same band as, or a band near to, the satellite navigation band of interest to prevent the signal to reach to the sensors or receivers, and spoofing is the transmission of a fake GNSS signal (Dempster 2016). There are studies showing that a simple and relatively cheap GNSS spoofer can be used to overtake, for example, a ship navigation without being detected; see Humphreys et al. (2008) and Divis (2013). Because the power level of GNSS signals is low, such signals are susceptible to interference; therefore, a relatively weak interference signal can jam a receiver (Dempster 2016). There are real examples that this interference affected operational infrastructures; see, for example, Balaei et al. (2007), Clynch et al. (2003), Grant et al. (2009), Hambling (2011), Motella et al. (2008), and Pullen et al. (2012). Specifically, we can point to unintentional cases such as a faulty TV amplifier, which jammed the global positioning system global positioning system (GPS) operation at a harbor in Monterey, California, for 37 days (Clynch et al. 2003). A small jammer that was used in a delivery van disrupted the ground-based augmentation system (GBAS) system aiding aircraft approaches at Newark Airport while driving on a nearby highway in 2009 (Hambling 2011; Pullen et al. 2012; Warburton and Tedeschi 2011). The Central Radio Management Office of South Korea reported several disruptions from 2010 to 2012 due to GPS jammers (Seo and Kim 2013). In Australia, Balaei et al. (2007) detected some interference and in Italy Motella et al. (2008) found some from TV signals in the GNSS band, disrupting GPS. Recognition of an interference signal among all scattered signals is a complicated process and requires skills in signal processing, which is outside the scope of this paper.Location of an interference device can be determined from basepoints with sensors or receivers, which are able to detect and analyze an interference signal. The time of arrival (TOA), angle of arrival (AOA), and time difference of arrival (TDOA) are known measurements being used for localization. Drake and Dogancay (2004) performed the localization problem on prolate spheroidal coordinates and stated that the mathematical equations of the TDOA will be greatly simplified by them. The equations will be linear corresponding to the hyperbolic asymptotes of the TDOA. Ananthasubramanian and Madlhow (2008) investigated AOA measurements and developed a sequential algorithm; they concluded that the localization error is proportional to the AOA error variance and coverage area, and reducible by increasing the number of estimates. Thompson et al. (2009) studied the configuration of the sensors’ positions for localization of interfering devices, presented a method using differences of received signal strength measurements, and concluded that they can be alternatives to the TOA, AOA, and TDOA. Thompson (2013) investigated the interference device detection and localization by analyzing the dilution of precision (DOP), from received signal strength, and concluded that TDOA is superior to the received signal strength measurements. Numerous literature studies exist regarding localization of GNSS interference using the mentioned quantities; however, optimization of configuration of the basepoints is a novel idea, which is presented in this paper.In geodetic network optimization, one of the purposes is to determine the optimal configuration of networks for maximizing the precision and reliability because the configuration has significant influence on the quality of the network; for more details, see Koch (1982, 1985), Xu (1989), Kuang (1996), Eshagh and Kiamehr (2007), and Eshagh and Alizadeh-Khameneh (2014). Stochastic and evolutionary methods such as genetic algorithm (GA) (Batilović et al. 2021), particle swarm optimization (PSO) (Yetkin et al. 2009; Singh et al. 2016), generalized particle swarm optimization (GPSO) (Batilović et al. 2022), and simulated annealing (SA) (Berné and Baselga 2004) are applicable for geodetic network optimization as well.The configuration, or the geometry formed by the basepoints for localization, is also a type of geodetic network. However, in the geodetic network optimization, the configuration is optimized by varying the unknown control points in such a way that the desired precision for these points is achieved, while in optimizing a wireless network for localization of interference device the configuration of the known basepoints is optimized by varying the basepoints to reach the desired precision for the control points. So far, no optimization has been used to determine optimal configuration of wireless security networks around airports, harbors, or important infrastructures.In this article, a new approach is proposed to determine the optimal configuration of the basepoints, at which the sensors or receivers are placed. The TOA, AOA, and TDOA of signals to these sensors and receivers are considered as observables and a criterion matrix is selected for required precision of the points over the control area. The basepoints move until the estimated variance-covariance (VC) matrices of the control points are fitted to the criterion matrix subjected to some constraints. A grid of points is designed over the control area and each one of its cells is the probable location of an interfering device. The location of basepoints is determined in such a way that the location errors of the grid points are minimized based on the type of observables and fit to the predefined criterion matrix for the grid points.Control Network and ObservablesPrecision of a network depends on both the quality of measurements and the network configuration. In this section, a two-dimensional control network for localization is defined and after that the observables TOA, AOA, and TDOA as well as their mathematical models are presented.Control NetworkThe control network is a grid of points over an area in addition to some basepoints from which an interfering device is localized by some measurements. The grid points, which are called hereafter control points, can be regarded as the probable position of this device. Fig. 1 illustrates a schematic control network; the basepoints and control points are respectively shown by triangles and small circles. Configuration means the geometric form, created by the basepoints; e.g., in Fig. 1, the three basepoints form a triangle as the configuration. Considering more basepoints leads to a more complicated configuration and a harder optimization process, but it is possible. The resolution of the grid, or the distance between the control points, is defined as the resolution of design. The spacing between the points along the x- and y-axes may not be same, it depends on the designer’s opinion.The main goal of considering a grid of points to cover an area is that the positions of the control points can be determined locally from the grid, e.g., the point at the left-down corner of the grid can be regarded as a local coordinate system origin and the first column as the y-axis and the lowest row as the x-axis. Because the grid resolution is known, coordinates of all control points can be simply determined. The important issue in optimization is the precision of these points when their coordinates are determined from the basepoints using some observables.ObservablesOur goal is to determine the horizontal position of an interfering device from some basepoints based on two-dimensional observables of TOA, AOA, and TDOA. If the velocity of the interfering signal is known, by measuring the TOA, its distances from the basepoints can be determined. Generally, localization with the TOA requires precise time synchronization of transmitter and receivers so that distances can be computed from the measured TOAs. However, in a two-dimensional (2D) localization by using at least three sensors or receivers, the transmission time can be considered as an extra unknown in the system of equations and approximated simultaneously with the coordinates of the transmitter. Today’s GNSS environmental monitoring systems (GEMSs) are probably able to measure this signal arrival. An advantage of using the TDOA is that the transmitter needs neither synchronizing with the receivers and sensors (see Gustafsson 2018, p. 78) nor sending the transmission time. The TDOAs between the basepoints are estimated by cross-correlation processes among the received signals (e.g., Lindstrom et al. 2007). In addition, there are new GEMSs consisting of several low-cost sensors to monitor GNSS system performance in a specific area (Trinkle et al. 2012). The AOA can also be determined at each basepoint using antenna arrays; see Trinkle et al. (2012) for mathematical modeling and estimation of AOA from these arrays (Trinkle et al. 2012; Huang et al. 2022).The mathematical formula of a distance between the interfering device at the point i and the jth basepoint is (1a) lij=(xj−xi)2+(yj−yi)2where xj and yj=x- and y-coordinates of the basepoint; and xi and yi = coordinates for the interfering device.The AOA is, in fact, the direction from which the interfering signal enters a sensor. This angle has the following mathematical expression in terms of the interfering device and the basepoint coordinates (e.g., Trinkle et al. 2012; Gustafsson 2018): (1b) θij=tan−1xj−xiyj−yiThe range difference (dijk) from an interfering device at point i and two basepoints at j and k has the following mathematical formula (Trinkle et al. 2012; Gustafsson 2018): (1c) dijk=(xj−xi)2+(yj−yi)2−(xk−xi)2+(yk−yi)2where xi and yi = coordinates of the interfering device; and xj, yj and xk, yk = pair coordinates of the jth and kth basepoints.Our goal is to find an optimal configuration for the basepoints in such a way that the error of the position of the interfering device is minimized. No measurement, speed of signal, and errors in the TOA, AOA, and TDOA are needed for such a design problem. When the location of the basepoints is optimally determined, stations equipped by sensors or receivers can be established for measuring the TOA, AOA, and TDOA. However, the risk of lack of connection between the basepoints and grid points always exists, which means the sensors and receivers at these points should be able to receive interference signal; otherwise, localization is not possible. Therefore, enough of such basepoints should be considered for covering the control area.Solution of the Coordinates and from the BasepointsFirst, assume that the basepoints have known coordinates and we want to see how large the errors of the control points are over the area based on our observables. However, because the observables have nonlinear mathematical models, they need to be linearized. Generally, a linearized model is of Gauss-Markov type for a control point over the area (see Koch 2010) (2a) Ax=L−ε,E{ε}=0E{εεT}=CL=σ02Qwhere A = coefficient matrix, whose elements are partial derivatives of observables with respect to coordinates of the control point; x = vector of the coordinate updates to their initial values; L = vector of difference between actual and computed observations from the initial coordinates; ε = vector of random errors with E{ε}=0, where E{} stands for the statistical expectation; CL = VC matrix of the observations; Q = cofactor matrix; and finally σ02 = a priori variance of unit weight.In Eq. (2a), the only parameter that carries the geometrical properties of the basepoints and the control point is the matrix A. This matrix contains the partial derivatives of observables with respect to the coordinates of the unknown control point, which are, in fact, sensitivities of the observables with respect to the coordinates. For example, in the case of 2D localization of one interference device from m basepoints, the structure of matrix A for TOA, AOA, or TDOA is (2b) A=[a1wb1wa2wb2w⋮⋮amwbmw],where  w=TOA, AOA, or TDOA and i=1,2,…,nwhere the elements of the matrix A are (2c) aiTOA=∂lij∂xi=−xj−xilijandbiTOA=∂lijyi=−yj−yilij j=1,2,…,m(2d) aiAOA=∂θij∂xi=−yj−yilij2andbiAOA=∂θij∂yi=−(xj−xi)(yj−yi)lij2 j=1,2,…,m(2e) aiTDOA=∂djik∂xi=−xj−xilij+xk−xilikandbiTDOA=∂dijk∂yi=−yj−yilij+yk−yilik j,k=1,2,…,mThe least-squares solution of Eq. (2a) is (e.g., Cooper 1987) (2f) x^=(ATQ−1A)−1ATQ−1Land the VC matrix of the estimated coordinates (2g) No observation is needed to compute the VC matrix of the estimated coordinates and the configuration can be optimized using Eq. (2g) and some a priori values of the measurement quality.Eq. (2a) describes a system of equations connecting the observables to two-dimensional coordinates of an unknown control point. Obviously, at least two observables from two different basepoints are needed to have a unique solution for the system. In the case of having more basepoints, the least-squares method is applicable.Suppose that three basepoints and three observables are used, and then the system in Eq. (2a) will have three equations and two unknowns for one control point. Our goal is to optimize the configuration of these basepoints so they cover the control area. Consequently, many such points should be considered over this area, which are probable locations of an interfering device.Optimal Configuration of the Basepoints Based on the Desired Errors of Control PointsAs shown, the VC matrices of control points are determined from the observables and location of the basepoints. The coordinates of these control points are known as the size of the grid and its cells are designed by us. Then the basepoints are considered as unknowns but an initial configuration is needed, which can be derived from the local coordinate system defined by the grid. Therefore, initial VC matrices for all control points are computed and after that the configuration is optimized by varying the coordinates of the basepoints.A criterion matrix is needed for the optimization process, which should be a priori known from the expected variances of the control points. In a two-dimensional network, each VC matrix has four elements, two variances as diagonal elements and two equal covariances as off-diagonal ones. Normally, a diagonal matrix with the desired variances for the coordinates of the unknown points on its diagonal elements is considered as the criterion. The initial VC matrices should be fitted to this criterion in a least-squares sense by varying the coordinates of the basepoints. In other words, by optimization the locations of the basepoints are determined in such a way that the VC matrix of the control points is fitted to the criterion.To present this idea mathematically, the VC matrix of one point [Eq. (2g)] is expanded for the basepoints by the Taylor series (e.g., Xu 1989; Koch 1982, 1985; Kuang 1996) (3a) Cx^=σ02(ATQ−1A)−1=Cx0+∑j=1m[∂Cx0∂xj∂Cx0∂yj][ΔxjΔyj]where A = initial coefficient matrix derived from the approximate positions of the basepoints; Δxj and Δyj = coordinate updates to the basepoints’ positions; and m = number of basepoints. It will not be difficult to show that (see also Kuang 1996; Eshagh and Alizadeh-Khameneh 2014) (3b) {∂Cx^0∂xj∂Cx^0∂yj}={∂∂xj(ATQ−1A)−1∂∂yj(ATQ−1A)−1}=−(ATQ−1A)−1{∂(ATQ−1A)∂xj∂(ATQ−1A)∂yj}(ATQ−1A)−1(3c) {∂(ATQ−1A)∂xj∂(ATQ−1A)∂yj}={∂AT∂xj∂AT∂yj}Q−1A+ATQ−1{∂A∂xj∂A∂yj}where the structure of (∂A)/(∂xj) or (∂A)/(∂yj), which we show that by (∂A)/[∂xj(yj)] for the TOA and AOA measurement in a 2D case of m basepoints and one interference device, will be (3d) ∂A∂x1(y1)=[∂a1w∂x1(y1)∂b1w∂x1(y1)0000⋮⋮],∂A∂x2(y2)=[00∂a2w∂x2(y2)∂b2w∂x2(y2)00⋮⋮],and∂A∂xm(ym)=[0000⋮⋮∂amw∂xm(ym)∂bmw∂xm(ym)]However, when TDOA is applied this structure for the first three basepoints is (3e) ∂A∂x1(y1)=[∂ai12TDOA∂x1(y1)∂ai12TDOA∂x1(y1)00⋮⋮∂ai1kTDOA∂x1(y1)∂ai1kTDOA∂x1(y1)00⋮⋮],∂A∂x2(y2)=[00∂ai21TDOA∂x2(y2)∂ai21TDOA∂x2(y2)00⋮⋮∂ai2kTDOA∂x2(y2)∂ai2kTDOA∂x2(y2)00⋮⋮],and∂A∂x3(y3)=[00⋮⋮∂ai31TDOA∂xk(yk)∂ai31TDOA∂xk(yk)00⋮⋮∂ai3kTDOA∂x3(y3)∂ai3kTDOA∂x3(y3)⋮⋮]The elements of these matrices are(3f) ∂aiTOA∂xj=∂2lij∂xj∂xi=−(yj−yi)2lij3and∂biTOA∂xj=∂2lij∂xj∂yi=(yj−yi)(xj−xi)lij3(3g) ∂aiAOA∂xj=∂2θij∂xj∂xi=2(yj−yi)(xj−xi)lij4and∂biAOA∂xj=∂2θij∂xj∂yi=1lij2−2(xj−xi)2lij4(3h) ∂aiTDOA∂xj=∂2dijk∂xj∂xi=−(yj−yi)2lij3and∂biTDOA∂xj=∂2lij∂xj∂yi=(yj−yi)(xj−xi)lij3(3i) ∂aiTOA∂yj=∂2lij∂yj∂xi=(yj−yi)(xj−xi)lij3and∂aiTOA∂yj=∂2lij∂yj∂yi=−(xj−xi)2lij3(3j) ∂aiAOA∂yj=∂2θij∂yj∂xi=1lij2−2(xj−xi)2lij4and∂bjAOA∂yi=∂2θij∂yi∂yj=−2(xj−xi)(yj−yi)lij4(3k) ∂aiTDOA∂yj=∂2dijk∂yj∂xi=(yj−yi)(xj−xi)lij3and∂biTDOA∂yj=∂2lij∂yj∂yi=−(xj−xi)2lij3Cx^ on the left-hand side of Eq. (3a) is the criterion matrix, and σ02(ATQ−1A)−1, which is the initial VC matrix, should be fitted to Cx^. Let us rewrite Eq. (3a) in the following form: (3l) Cx^−Cx0=∑j=1m[∂Cx0∂xj∂Cx0∂yj][ΔxjΔyj]Cx^ is a 2×2 matrix having the variances of the x- and y-coordinates on its diagonal elements and covariance between them on the off-diagonal. Eq. (3l) means that the updates Δxj and Δyj to the initial coordinates xj and yj are estimated in such a way that the updated Cx0 is fitted to the desired Cx^. Because in this study three basepoints are considered and each has two coordinates (x and y), six unknown parameters exist in the system of equations. Our VC matrices are 2×2, meaning they contain four elements and our system has four equations for each control points. However, there are many points already designed over the control area, and each one adds four equations to our systems, while the number of the unknown parameters remains constant. In this case, a huge system of equations is created and the coordinates of the basepoints are estimated in such a way that the best fit to all criterion matrices for the control points are achieved.Let us present this system in the following form: (4a) where ε′ = vector of residuals; ΔL = vector of differences between the elements of the criterion and initial VC matrices; B = coefficients matrix containing the partial derivatives of the VC matrix with respect to the basepoints’ coordinates; and Δx = basepoints’ coordinate updates. The mathematical descriptions of them are (4b) ΔL=[vec(Cx^1)−vec(Cx10)vec(Cx^2)−vec(Cx20)⋯vec(Cx^n)−vec(Cxn0)]T(4c) B=[vec(∂Cx10∂x1)vec(∂Cx10∂y1)vec(∂Cx10∂x2)vec(∂Cx10∂y2)⋯vec(∂Cx10∂xm)vec(∂Cx10∂ym)vec(∂Cx20∂x1)vec(∂Cx20∂y1)vec(∂Cx20∂x2)vec(∂Cx20∂y2)⋯vec(∂Cx20∂xm)vec(∂Cx20∂ym)⋮⋮⋮⋮⋱⋮⋮vec(∂Cxn0∂x1)vec(∂Cxn0∂y1)vec(∂Cxn0∂x2)vec(∂Cxn0∂y2)⋯vec(∂Cxn0∂xm)vec(∂Cxn0∂ym)]n×2m(4d) Δx=[Δx1Δy1Δx2Δy2⋯ΔxmΔym]Twhere operator vec inserts the columns of the VC matrices below each other and converts the 2×2 matrices to 4×1 vectors; ()T stands for transposition operator of matrix algebra; and n = number of control points.ConstraintsAt first glance, one expects to solve the system of Eq. (4a) by the least-squares method. However, unrealistic results may be obtained, for example, some basepoints may be located far from the area or they may be located close to the center and create a weak geometrical configuration. Hence, some constraints are needed for solving this system of equations, like: •The basepoints should be around the area.•Their average x- and y-coordinates should be at the center of the area.•They should be well distributed around the area.The use all these constraints depends on the study area and how much freedom our basepoints must move during the optimization process. These constraints act like inner constraints (e.g., Cooper 1987; Tan 2005) that are applied in geodetic network adjustment for overcoming rank deficiency of the established system of equations due to datum defect. However, our system does not have this deficiency because the grid points already have known coordinates. These constraints are used to limit the movements of the basepoints.Limiting Search Area for the BasepointsTo keep the basepoints outside the control area, some inequality constraints are applied for limiting the search domain of the coordinate variations; such inequality constraints are (compare Kuang 1996; Eshagh and Alizadeh-Khameneh 2014) (5a) (5b) where wjL and wjU = lower and upper bounds of the inequality that contain limiting x-coordinates of the jth basepoint, respectively; and vjL and vjU = corresponding limits for the y-coordinates. To write these constraints in terms of the coordinate updates being estimated from the optimization process, Eqs. (5a) and (5b) are written in the following forms (compare Eshagh and Alizadeh-Khameneh 2014): (5c) wjL−xj≤Δxj≤wjU−xj(5d) vjL−yj≤Δyj≤vjU−yjAccording to Eqs. (5c) and (5d), the coordinate updates in such a way that the coordinates remain in the specified interval [Eqs. (5a) and (5b)] outside the area.These inequality constraints can be written in the following vector form: (5e) where (5f) Lb=[w1L−x1v1L−y1w2L−x2v2L−y2⋯wmL−xmvmL−ym]T(5g) Ub=[w1U−x1v1U−y1w2U−x2v2U−y2⋯wmU−xmvmU−ym]TAverage Coordinates of the Basepoints at the Center of AreaTo keep the average of the basepoints’ coordinates at the center of control area, the following equations can be used (compare Cooper 1987; Kuang 1996): (6a) (6b) where xG and yG = coordinates of the center of the control area, which can be determined from the control points; and m=3 is the number of the basepoints.Because the nonlinear mathematical models are linearized in our optimization process, Δxj and Δyj are estimated and not the coordinates directly. Therefore, Eqs. (6a) and (6b) need to be rewritten in terms of these coordinate updates (6c) 1m∑j=1mΔxj=xG−1m∑j=1mxj=[1010⋯10]1×2mx=xG−1m∑j=1mxj(6d) 1m∑j=1mΔyj=yG−1m∑j=1myj=[0101⋯01]1×2mx=yG−1m∑j=1myjor (6e) where (6f) D1=1m[1010⋯100101⋯01]andd1=[xGyG]−1m[∑j=1mxj0∑j=1myj0]Distribution of the Basepoints around the AreaBearing angles from the center of the area to the basepoints can be used for making a constraint for distributing the basepoints around the area. By assuming that the summation of all these bearing angles is zero, we can write (compare Kuang 1996) (7a) ∑j=1mtan−1xj−xGyj−yG=0Again, Eq. (7a) should be written in terms of the coordinate updates; therefore, its linearized form will be applied (7b) ∑j=1m[yj−yGljG2−xj−xGljG2][ΔxjΔxj]=0−∑j=1mtan−1xj−xGyj−yGwhere ljG = distance from the center of control area and the basepoint j. Eq. (7b) can be written in the following matrix form (compare Kuang 1996): (7c) where (7d) [y1−yGl1G2−x1−xGl1G2y2−yGl2G2−x2−xGl2G2⋯ym−yGlmG2−xm−xGlmG2](7e) d2=−∑j=1mtan−1xj−xGyj−yGOptimization ModelThe system of equations in Eq. (4a) should be solved for the coordinate updates in a least-squares sense but subjected to the mentioned constraints. Such an optimization model is (compare Koch 1982, 1985; Kuang 1996) (8) min(12ΔxTBTBΔx−BTΔL)subject  toD1Δx=d1D2Δx=d2Lb≤Δx≤UbConsidering all these constraints, Eq. (8) might not be possible in practice. As shown in the synthetic test, all basepoints have the freedom to move, but their movements can be limited by the three constraints so they can move around the area using the search area constraint and around the center of the area with a good distribution. However, in the real study for the international Landvetter Airport, based on the presence of forests around the area, keeping the basepoints around the area is not practically possible and only the limiting search area around each basepoint is suitable to keep them in area.Today, there are different software programs for solving the optimization problem in Eq. (8). The theory of solving this problem is known (e.g., Bazaraa and Shetty 1976; Grafarend and Sanso 1985). However, the optimization problem in Eq. (8) is not global because of the involvement of different constraints. The unconstrained optimization may diverge, or the basepoints may move outside the area or toward themselves or become colinear. Applying at least the search area constraint is a necessity for obtaining a convergent solution. The resolution of the grid of control points plays a role in the rate of convergence and computational time. In addition, a threshold should be a prior defined to stop the iterative optimization; in this study, the norm of the coordinates updates is computed and when it becomes smaller than 10 cm during the iterations, the optimization is stopped.Synthetic Numerical TestTo test the presented theory, a rectangular area of 4×5  km is considered with resolutions of 10, 20, 40, and 80 m. Three basepoints, M, N, and O, are chosen manually around the area; see their initial coordinates in Table 1 and their locations in Fig. 2(a) by small circles. Fig. 2(a) shows the pattern of the square root of the trace of the VC matrix of the control points, i.e., horizontal dilution of precision (HDOP), based on the TOA observables and a resolution of 40×40  m. Because the HDOP patterns for different resolutions were similar, only one is presented. Generally, the HDOP based of the TOA has no units, but we considered a priori precision of 1 m for the TOA distances to give a unit to the computed HDOP. Such a precision seems to be large for geodetic networks, but it is acceptable in localization of the interference device by wireless networks (e.g., Trinkle et al. 2012). As observed, the largest error reaches to about 2.2 m, close the farthest point O at the upper-right corner of the area. A diagonal criterion matrix with equal diagonal elements to 2 m was considered for fitting the VC matrices of the control points to it. The basepoints were kept 400 m, as an example, outside the area from the margins by using the inequality constraint in Eq. (5e). Fig. 2(b) shows the optimized location of the basepoints on the map of the HDOP. A significant reduction of HDOP is seen with the maximum reducing to 1.6 m. In addition, the configuration forms almost an equilateral triangle, which is more stable than the initial design forming an isosceles triangle. However, the weak points having larger HDOPs are seen around the basepoints and elongated toward the center of the area. This could be expected because two of the three distances from each control point and to the three basepoints are considerably longer than the other one, forming a weaker geometry than the points located at the center of the area.Table 1. Optimal coordinates for the basepoints M, N, and O based on resolution of design and type of observableTable 1. Optimal coordinates for the basepoints M, N, and O based on resolution of design and type of observableObservableBasepointx (m)y (m)x (m)y (m)x (m)y (m)x (m)y (m)TOAM0.03,605.60.03,605.60.03,950.10.03,575.9N1,600.1−4001,600−4002,000.1−4001,600−400O4,4004,294.54,4004,294.54,0003,9504,4004,264.1AOAM−400.03,950.1−400.03,950.1−400.03,950.1−400.03,920.1N2,000.1−400.02,000.1−400.02,000.12,000.12,000.1−400.0O4,400.03,950.04,400.03,950.04,400.03,950.04,400.03,920.0TDOAM−355.73,711.5−363.93,718.7−380.73,733.2−400.03,720.1N1,955.70.01,964.00.01,980.70.02,000.10.0O4,400.03,788.54,400.03,781.44,400.03,766.84,400.03,720.1A similar HDOP pattern is seen when the AOA are considered. The same criterion matrix was considered for the control points and a priori variance factor of 60 arc-seconds as the accuracy the measured AOAs. Fig. 2(c) shows that the HDOP before optimization reaches to about 3.5 m, and after optimization, Fig. 2(d), it significantly reduced to about 1.8 m with a symmetric configuration.Fig. 2(e) is a similar HDOP for the TDOA, with a priori variance factor of 1 m, based on the same resolution before optimization. It shows that the HDOP reaches to about 3.5 m and the weak points are seen at the corners. Generally, a better coverage over the area can be seen using the TDOA. The HDOP after optimization is shown in Fig. 2(f). As observed, a symmetric geometry is created by the optimization process and the weak points are around the corners of the area, but with a proper coverage over the inner area meanwhile reducing the HDOP to 2.5 m.The HDOP derived from the TDOA gives a good coverage over the area, and the weak points have large HDOPs and are not close to the basepoints. However, the mathematical models are more complicated because of their differential nature.To test sensitivity of our optimization method based on the observables, the process was repeated with different resolutions. Table 1 gives the optimal horizontal positions of the three basepoints after optimization, while the initial positions of these points were M (0.5, 2,500.5), N (2,000.5, 0.5), and O (4,000.5, 5,000.5) prior to optimization.According to Table 1, when the TOA is used, the y-coordinate of M remains the same for resolutions 10×10  m and 20×20  m, but the x-coordinate remains constant for all. There are variations in the x-coordinate of N but it is almost constant for these resolutions as is the x-coordinate of O. Generally, with resolutions less than 20×20  m, the optimized coordinates of the points are almost the same as the observable.Similar results are obtained when the AOA is used. A resolution of 30×30  m is good enough for optimization, but for TDOA no convergence is seen by increasing the resolution. Nevertheless, the differences between the coordinates obtained from optimization using the two successive resolutions are in order of the resolutions. Therefore, the TDOA is more sensitive to the selected resolution and no specific resolution can be determined for optimization based on it.Generally, applying a dense grid for the control points is suggested, but it takes a huge computational burden and long time to perform the optimization process. One important issue that needs to be highlighted is that our study area is rectangular 4×5  km, close to a square.Real Case Study: Landvetter AirportSo far, our simulation study showed that the presented theory works with some constraints. Designing a security network is not as simple as what was done in the simulation study unless the area is flat with no signal hindrance. Here, an international airport is considered and positions of three basepoints are optimized in the designed control network with the assumption of the presence of the interference device inside the airport. Not all our proposed constraints can be applied for this study area. First, there are forests around, which are not suitable places for our sensors. Therefore, spreading the basepoint outside and around the airport is not realistic. Second, search areas around the basepoints are limited because there is not much freedom for them to move.The international Landvetter Airport of Sweden in Gothenburg was selected, and a local planar coordinate system was defined for the airport. The Gaussian radius of curvature at the system origin was computed from its latitude at the surface of the WGS84 reference ellipsoid with a=6,378,137  m, and e2=0.0068. A grid with a resolution of 40×40  m was considered over the airport; was rotated based on azimuth of the most left takeoff way (extracted from Google Maps), and later their coordinates were transformed to the geodetic coordinates where φi and λi = geodetic coordinates of the ith point of the grid; φO and λO = geodetic coordinates of the local system origin; xi and yi = local Cartesian coordinates of the point; and R = Gaussian radius of the curvature and the local system’s origin (compare Vanicek and Krakiwsky 1986; Jekeli 2012) R=MN,with  N=a(1−e2sin2φO)1/2  and  M=a(1−e2)(1−e2sin2φO)3/2where N and M = well-known radius of the prime vertical and curvature of the local meridian at the point with the latitude φO, respectively.Landvetter Airport’s runway is elongated in southwest–northeast direction with an approximate bearing angle of 26°. The area is almost rectangular, with length of 5 km and width of 2 km. Based on the geometry of the airport, no symmetric configuration can be considered for the basepoints. Therefore, from the satellite photo and Google Earth, we selected three suitable places; Fig. 3 shows thee selected areas by the squares for optimal locations of the basepoint. In fact, an initial position was considered for each basepoint, which varied in a square with the size of 400 m.Fig. 4(a) shows the maps of HDOP for the TOA over Landvetter Airport and the small circles represent the basepoints before optimization. The points having larger HDOPs were extended along the runway of the airport. The HDOP reaches to about 6 m over the area. Fig. 4(b) illustrates the optimized positions for the basepoints, and it shows that the optimization has reduced the HDOP to 4 m with slight change in the initial positions; see Table 2. Fig. 4(c) is a similar map to Fig. 4(a) but for the AOA with the same initial positions for the basepoints. After optimization, the map in Fig. 4(d) is obtained, and small changes in the positions of the basepoints lead to reduction of the HDOP to less than 4 m. Because of the special narrow rectangular form of the airport, the initial design based on the TDOA was not successful; the reason could be that the hyperbolic curves, which are created from each two points, do not make a good geometry. The HDOP at some points became extremely large and unrealistic, and even if the optimization was performed the constraints were violated and the basepoints moved outside the area. For this reason, no result is presented for design using the TDOA.Table 2. Initial and optimized coordinates of basepoints before and after optimizationTable 2. Initial and optimized coordinates of basepoints before and after optimizationObservableBasepointx (m)y (m)x (m)y (m)TOAM600.54,000.5400.53,800.5N600.5200.5533.4229.26O1,299.53,199.51,499.52,999.5AOAM600.54,000.5400.53,800.5N600.5200.5400.50.5O1,299.53,199.51,499.52,999.5Table 2 lists the initial and optimized positions of the basepoints M, N, and O. Optimization based on the TOA and AOA leads to the same positions for the basepoints M and O, but N is more influenced by the type of observables. Generally, M and O are closer to each other, but N is far from both. These points form a triangle having a small angle at N. It is known that an equilateral or equiangular triangle is a strong geometric form. When a triangle has a different form, it has weaker geometry, such as when M, N, and O form a triangle having such a small angle at N. This is the reason N moves during the optimization process, because it has direct influence on the strength of the geometry of the basepoints.ConclusionsTOA, AOA, or TDOA of signals into sensors or receivers can be used for localization of a signal interference device. This article showed that the optimal geometrical configuration formed by basepoints, at which these signal arrivals are measured, strongly depends on the type of arrival and the ability of sensors or receivers to measure these arrivals. A design using the AOA showed a better stability and robustness with respect to the resolution of the designed grid compared to the rest of its rivals, and optimization based on the TDOA improves coverage over the control area. Generally, our optimization method needs some constraints for stabilizing the solution; otherwise, the basepoints may form improper configuration during the iterative optimization process, e.g., they may move toward each other or become collinear. However, in our real case study, the search domains for the basepoints were limited and subjecting the optimization process to extra constraints, as done in the simulation test, led to infeasible solutions. In other words, the choice of the constraints is highly dependent on the shape of the area.Data Availability StatementThe data that support the findings of this study are available from the author upon reasonable request.AcknowledgmentsThe author is thankful to both reviewers and their fruitful comments enhancing the quality of the paper.References Ananthasubramanian, B., and U. Madlhow. 2008. “Cooperative localization using angle of arrival measurements in non-line-of-sight environments.” In Proc., Melt’08 September 19, New York: Association for Computing Machinery Digital Library. Balaei, A. T., B. Motella, and A. G. Dempster. 2007. “GPS interference detected in Sydney Australia.” In Proc., Int. Symp. GPS/GNSS IGNSS2007. Sydney, Australia: IGNSS Society. Batilović, M., Z. Sušić, Ž. Kanović, M. Z. Marković, D. Vasić, and V. Bulatović. 2021. “Increasing efficiency of the robust deformation analysis methods using genetic algorithm and generalised particle swarm optimization.” Surv. Rev. 53 (378): 193–205. Bazaraa, M. S., and C. M. Shetty. 1976. Vol. 122 of Foundations of optimization. Lecture notes in economics and mathematical systems book series. New York: Springer. Clynch, J. R., A. A. Parker, R. W. Adler, and W. R. Vincent. 2003. “System challenge—The hunt for RFI—Unjamming a coast harbor.” GPS World 14 (1): 16–22. Cooper, M. A. R. 1987. Control surveys in civil engineering. 1st ed. New York: Nicols Publication. Divis, D. A. 2013. “GPS spoofing experiment knocks ship off course.” InsideGNSS 31 (Jul): 1–3. Drake, S. P., and K. Dogancay. 2004. “Geolocation by time difference of arrival using hyperbolic asymptotes.” In Proc., ICASSP 2004. New York: IEEE. Eshagh, M., and M. A. Alizadeh-Khameneh. 2014. “The effect of constraints on bi-objective optimisation of geodetic networks.” Acta Geod. Geophys. 50 (4): 449–459. Eshagh, M., and R. Kiamehr. 2007. “A strategy for optimum designing of the geodetic networks from the cost, reliability and precision views.” Acta Geod. Geophys. Hung. 42 (3): 297–308. Grafarend, E. W., and F. Sanso. 1985. Optimization and design of geodetic networks. Berlin: Springer. Gustafsson, F. 2018. Statistical sensor fusion. 3rd ed. Lund, Sweden: Studentlitteratur AB. Hambling, D. 2011. “GPS chaos: How a $30 box can jam your life.” New Sci. 2803 (Mar): 44–47. Huang, L., Z. Lu, Z. Xiao, C. Ren, J. Song, and B. Li. 2022. “Suppression of jammer multipath in GNSS antenna array receiver.” Remote Sens. 14 (2): 350. Humphreys, T. E., B. M. Ledvina, M. L. Psiaki, B. W. O’Hanlon, and P. M. Kintner. 2008. “Assessing the spoofing threat: Development of a portable GPS civilian spoofer.” In Proc., 21st Int. Technical Meeting Satellite Division of the Institute of Navigation (ION GNSS 2008), 2314–2325. Manassas, VA: Institute of Navigation. Jekeli, C. 2012. Geometric reference systems in geodesy, lecture notes. Columbus, OH: Ohio State Univ. Koch, K. R. 1982. “Optimization of the configuration of geodetic networks.” Dtsch. Geodaetische Kommission 3 (258): 82–89. Koch, K. R. 1985. “First order design: Optimization of the configuration of a network by introducing small position changes.” In Optimization and design of geodetic networks, edited by E. W. Grafarend and F. Sanso, 56–73. Berlin: Springer. Koch, K. R. 2010. Parameter estimation and hypothesis testing in linear models. 2nd ed. Berlin: Springer. Kuang, S. 1996. Geodetic network analysis and optimal design: Concepts and applications. Chelsea, MI: Ann Arbor Press. Lindstrom, J., D. M. Akos, O. Isoz, and M. Junered. 2007. “GNSS interference detection and localization using a network of low cost front-end modules.” In Proc., 20th Int. Technical Meeting of the Satellite Division of The Institute of Navigation (ION GNSS 2007), 1165–1172. Manassas, VA: Institute of Navigation. Motella, B., M. Pini, and F. Dovis. 2008. “Investigation on the effect of strong out-of-band signals on global navigation satellite systems receivers.” GPS Solut. 12 (2): 77–86. Pullen, S., G. Gao, C. Tedeschi, and J. Warburton. 2012. “The impact of uninformed RF interference on GBAS and potential mitigations.” In Proc., 2019 Int. Technical Meeting of the Institute of Navigation, 780–789. Manassas, VA: Institute of Navigation. Seo, J., and M. Kim. 2013. “eLoran in Korea—Current status and future plans.” In Proc., European Navigation Conf. (ENC-GNSS). Bonn, Germany: German Institute of Navigation. Thompson, R. J. R. 2013. “Detection and localisation of radio frequency interference to GNSS reference stations.” Ph.D. thesis, School of Electrical Engineering and Telecommunications, Univ. of New South Wales. Thompson, R. J. R., A. Tabatabaei Balaei, and A. Dempster. 2009. “Dilution of precision for GNSS interference localisation systems.” In Proc., European Navigation Conf. (ENC GNSS2009). Naples, Italy: Parthenope Univ. Trinkle, M., E. Cetin, R. J. R. Thompson, and A. G. Dempster. 2012. “Interference localisation within the GNSS environmental monitoring system (GEMS)—Initial field test results.” In Proc., 25th Int. Technical Meeting of the Satellite Division of The Institute of Navigation (ION GNSS 2012), 2930–2939. Manassas, VA: Institute of Navigation. Vaníček, P., and E. J. Krakiwsky. 1986. Geodesy: The concepts. Amsterdam, Netherlands: Elsevier. Warburton, J., and C. Tedeschi. 2011. “GPS privacy jammers and RFI at Newark: Navigation team AJP-652 results.” In Proc., 12th Int. GBAS Working Group Meetings (I-GWG-12). Washington, DC: Federal Aviation Administration. Yetkin, M., C. Inal, and C. O. Yigit. 2009. “Use of the particle swarm optimization algorithm for second order design of levelling networks.” J. Appl. Geod. 3 (3): 171–178.

Source link

Leave a Reply

Your email address will not be published. Required fields are marked *