Indoor Positioning Using Angle of Departure Information

(1)

Department of Science and Technology

Institutionen för teknik och naturvetenskap

Linköping University

Linköpings universitet

g

n

i

p

ö

k

r

o

N

4

7

1

0

6 n

e

d

e

w

S

,

g

n

i

p

ö

k

r

o

N

4

7

1

0

6 -E

S

LiU-ITN-TEK-A--15/061--SE

Indoor Positioning Using Angle

of Departure Information

Erica Gunhardson

(2)

LiU-ITN-TEK-A--15/061--SE

Indoor Positioning Using Angle

of Departure Information

Examensarbete utfört i Elektroteknik

vid Tekniska högskolan vid

Linköpings universitet

Erica Gunhardson

Handledare Adriana Serban

Examinator Qin-Zhong Ye

(3)

Upphovsrätt

Detta dokument hålls tillgängligt på Internet – eller dess framtida ersättare –

under en längre tid från publiceringsdatum under förutsättning att inga

extra-ordinära omständigheter uppstår.

Tillgång till dokumentet innebär tillstånd för var och en att läsa, ladda ner,

skriva ut enstaka kopior för enskilt bruk och att använda det oförändrat för

ickekommersiell forskning och för undervisning. Överföring av upphovsrätten

vid en senare tidpunkt kan inte upphäva detta tillstånd. All annan användning av

dokumentet kräver upphovsmannens medgivande. För att garantera äktheten,

säkerheten och tillgängligheten finns det lösningar av teknisk och administrativ

art.

Upphovsmannens ideella rätt innefattar rätt att bli nämnd som upphovsman i

den omfattning som god sed kräver vid användning av dokumentet på ovan

beskrivna sätt samt skydd mot att dokumentet ändras eller presenteras i sådan

form eller i sådant sammanhang som är kränkande för upphovsmannens litterära

eller konstnärliga anseende eller egenart.

För ytterligare information om Linköping University Electronic Press se

förlagets hemsida

http://www.ep.liu.se/

Copyright

The publishers will keep this document online on the Internet - or its possible

replacement - for a considerable time from the date of publication barring

exceptional circumstances.

The online availability of the document implies a permanent permission for

anyone to read, to download, to print out single copies for your own use and to

use it unchanged for any non-commercial research and educational purpose.

Subsequent transfers of copyright cannot revoke this permission. All other uses

of the document are conditional on the consent of the copyright owner. The

publisher has taken technical and administrative measures to assure authenticity,

security and accessibility.

According to intellectual property law the author has the right to be

mentioned when his/her work is accessed as described above and to be protected

against infringement.

For additional information about the Linköping University Electronic Press

and its procedures for publication and for assurance of document integrity,

please refer to its WWW home page:

http://www.ep.liu.se/

(4)

Sammanfattning

I detta examensarbete undersöks möjligheten att kunna använda en positione-ringsmetod som inte enbart förlitar sig på den uppmätta signalstyrkan. Istället används en metod som bestämmer från vilken vinkel en signal uppkommer ifrån. Den här tekniken kallas för direction-finding. När informationen om signalens vinkel fastställts används den i ett positioningsfilter som uppskattar positionen. Två tillvägagångssätt har använts i den här rapporten, ett där enbart vinkeln an-vänds och ett där både signalstyrka och vinkel anan-vänds.

Simuleringar där direction-finding algoritmen tillsammans med positionerings-filtret har använts, visar goda resultat. Verklig data behövs för att vidare kunna analysera prestandan hos det framtagna positioneringssystemet.

(5)

(6)

Abstract

This master thesis investigates the possibility to use an indoor positioning ap-proach that does not completely rely on the strength of the received signal. In-stead, information from what angle the incoming signal originates is utilized. This technique is called direction-finding. When the angle information is deter-mined, it is used in a position filter. Two approaches for the estimation filter has been conducted, one which relies only on the angle information and one that re-lies both on the received signal strength and the angle information.

Simulations using the direction-finding algorithm together with the estimation filter generates promising results. Real data are required to further analyze the performance of the positioning system proposed in this thesis.

(7)

(8)

Acknowledgments

I would like to thank my supervisor at SenionLab, David Törnqvist, for his con-stant support in my work. I would like to thank Per Skoglar from SenionLab for great input along the course of the project.

I would also like to give a special thank you to my supervisor at Linköping Uni-versity, Adriana Serban, for her valuable insights, guidance and encouragement. Another special thank you to Johan Gustafsson who never stops believing in me and supports me no matter what.

Erica Gunhardson Linköping, June 2015

(9)

(10)

Notation

Abbreviations

Abbreviation Description aoa Angle-of-Arrival aod Angle-of-Departure

ble Bluetooth Low Energy df Direction Finding ekf Extended Kalman Filter

esprit Estimation Parameters via Rotational Invariance Tech-nique

gps Global Positioning System ips Indoor Positioning System kf Kalman Filter

mcu MCU

mems Micro Electro Mechanical Sensor mimo Multiple Input Multiple Output

miso Multiple Input Single Output music Multiple Signal Classification

od Output Data odr Output Data Rate

rf Radio Frequency

rssi Received Signal Strength Indicator simo Single Input Multiple Output

siso Single Input Multiple Output snr Signal to Noise Ratio

vhf Very High Frequency hf High Frequency

(13)

(14)

1

Introduction

This thesis covers the development of a direction-finding algorithm together with a position filter. The project was carried out as a master thesis on the Electronic Design program at Linköping University. The foundation of the thesis was pro-vided by SenionLab AB.

1.1 Background and Motivation

The interest in identifying ones geographical position is ancient and various tech-niques have been used over the years. Modern technology became a part of it when the first step towards the Global Positioning System (gps) was made in the 1950s [9]. While the positioning system for outdoor use has been widely re-searched and developed, the positioning system for indoor use has been lagging behind. There is no definite standard for indoor positioning today, but the most-frequently used technique is to use the received signal strength [19]. Together with a technique called fingerprinting, in which a radio map over the indoor en-vironment is created, the strength of the signal makes it possible to estimate ones location. A disadvantage with this technique is when it is used in larger areas where the changes in the received signal strength are too small to be recognized due to the weakness of the signal strength.

In the light of this, a localization approach called Direction-Finding (df) will be investigated. Instead of measuring the received strength of a signal, the direction from which the signal originated will be utilized. With this proposed technique, an estimation filter will be implemented to enable a localization estimate, using the result from the df algorithm.

(15)

2 1 Introduction

1.2 Objectives

The main objective of this thesis has been to investigate and develop a df algo-rithm as a complement with the received signal strength . The main potential benefit would be to eliminate the problem of finding ones location in larger ar-eas. As a complement to the result of the algorithm, a position estimate filter has been investigated and used.

1.3 Limitations

A pre-determined hardware design for a beacon device has acted as a design foundation for the df algorithm. To minimize the current consumption of the beacon, all intended calculations will be performed within the receiving device. Specifications for the antenna constellation for the receiver is pre-determined and the microcontroller within the receiver is pre-determined to be the Texas Instruments radio chip CC2650.

Figure 1.1:Block diagram of the system.

Figure 1.1 presents the antenna constellation design for the system. A re-search foundation for the functionality of such device will be presented as a fu-ture work suggestion. Hardware design of the receiving device lies outside the scope of this thesis.

All data are simulated, which limits the results to only be theoretical and are not validated with real data.

1.4 Approach

The possibility of determining a position based on the Angle of Departure (aod) information is investigated. The technique is called df and requires an array of antennas, either on the transmitting end, the receiver end or both. Generally, the antenna array is located at the receiving end of the system. A df algorithm will provide the angle of the signal but this information is not sufficient to determine the position of a target. To get an estimation of its position, two methods called bearings-range and bearing-only localization, which is based on the Extended Kalman Filter (ekf), is investigated.

(16)

1.5 Outline of the Thesis 3

Both algorithms are simulated in MATLAB, in separate files to get a better un-derstanding of each method individually.

1.5 Outline of the Thesis

Chapter 2 outlines the theory and background regarding different df methods. It also describes the differences in various antenna constellations.

Chapter 3 provides a deeper understanding of a specific df algorithm called Mul-tiple Signal Classification (music). This algorithm will serve as a foundation for a further development called M:1 music algorithm that suit the specifications of this thesis.

Chapter 4 presents simulation of the df algorithm.

Chapter 5 outlines the theory behind the bearings-range and bearing-only local-ization, it also includes the theory of the ekf.

Chapter 6 explains the modeled scenario of the bearing-range, bearing-only and ekfmethods, it also includes simulations of these methods.

Chapter 7 provides a research foundation for a receiver device called position-ing tag.

(17)

(18)

2

Direction-Finding Theory

This chapter includes the theory behind different df methods, together with the theory behind different antenna arrangements.

2.1 Direction-Finding Techniques

dfis a method to provide the functionality of Angle of Arrival (aoa) for a re-ceived signal, or the Angle of Depature (aod) of a transmitted signal. The method uses a single or multiple transmit antennas and a single or multiple receiving an-tennas to determine the azimuth angle of a distant emitter [23].

The df method has existed for as long as electromagnetic waves have been known. Heinrich Hertz discovered the directivity of antennas in 1888.

2.1.1 Watson-Watt Method

One of the first methods of df was the polarisation df, frequently used in the First World War. After the war, Sir Watson-Watt developed a non-mechanical df system using crossed loop antennas. These type of system consist of four spatially displaced monopole or vertical dipole antennas. The angle of the incoming signal is determined by the differences in amplitude of the received waveform. One additional antenna is usually placed in the centre to resolve uncertainties in the bearing, Figure 2.1.

(19)

6 2 Direction-Finding Theory

Figure 2.1:Typical Watson-Watt antenna arrangement.

The system compares the amplitudes received by each antenna, which are ob-tained as voltages out from the transmitting antenna. The antenna constellation creates two axes, one in the y-direction that gives the y-axis angle and one in the x-direction that gives the x-axis angle. By changing the spacing of the antennas, the system is capable to operate over a wide range of frequencies.

2.1.2 Doppler Method

In 1941, the Doppler principle of df systems appeared, consisting of a circular array of antennas. In 1950, airports were equipped with vhf/hf Doppler df systems for air traffic control. Doppler df systems use phase differences in the received antenna array. Early designs used rotating antenna arrays to obtain the Doppler shift, but more modern designs called Pseudo-Doppler systems use elec-trical methods to simulate the rotation. These systems use the changes in the ve-locity of a signal introduced by the rotating antenna array to induce the Doppler effect. Pseudo-Doppler df uses four equally spaced antennas positioned on the circumference on a circle. The antennas are being switched between in order to simulate the rotation.

Figure 2.2:Typical Doppler antenna arrangement.

If a signal would be incoming from the North, Figure 2.2, the switching makes the East-antenna to rotate away from the source and the West-antenna to rotate towards the source. There is no Doppler effect on the South-antenna. This shows that the signal is coming from the North.

(20)

2.2 Antenna Constellations 7

2.1.3 Interferometry

From the 1970s, the systems started to be digitalized and due to digital signal pro-cessing advances in the early 1980s, interferometry df systems were developed. Interferometry df systems uses phase differences to determine the angle of the transmitted signal. This method consists of an antenna array which is equally spaced in a straight line. Since the frequency of transmission is known, the phase differences can be calculated for each of the antennas. By knowing the phase differences and the distance between each antenna the angle can be calculated.

Figure 2.3:Typical interferometry antenna arrangement.

A typical interferometry antenna arrangement is shown in Figure 2.3, where the arrows represents the incoming signals and the dashed lines represents the advancing wavefront [22].

2.2 Antenna Constellations

In wireless communication, the constellation of the antenna systems are of im-portance. This is to know how to create a mathematical model of the system. The communication systems can be divided into four main types called single-input single-output (siso) system, single-single-input multiple-output (simo) system, multiple-input multiple-output (mimo) system and multiple-input single-output (miso) system. A siso system is the most commonly used consisting of one trans-mitting antenna and one receiving antenna, Figure 2.4.

Figure 2.4: sisosystem.

The simo system consist of one transmitting antenna and N receiving anten-nas, Figure 2.5.

In a simo system, due to the receiving antenna array, it is the aoa that is of importance when the phase difference will be determined.

(21)

8 2 Direction-Finding Theory

Figure 2.5: simosystem.

The miso system consist of M transmitting antennas and one receiving an-tenna, Figure 2.6.

Figure 2.6: misosystem.

When the phase difference of a miso system is to be determined, it is the aod that is of importance. The most complex system is the mimo, which consist of M transmitting antennas and N receiving antennas, Figure 2.7.

Figure 2.7: mimosystem.

In a mimo system, both the aoa and aod information can be used to deter-mine the phase difference of signals [6]. The transmitting and receiving end can be categorized depending on the number of rf channels. When each antenna got one rf channel, its called a multi-channel transmitter, and when there is only one rf channel independent of the number of antennas, its called a single-channel transmitter. Figure 2.8 illustrates a multi-channel transceiver.

(22)

2.2 Antenna Constellations 9

Figure 2.9 illustrates the single-channel transceiver.

Figure 2.9:An example of a single-channel transceiver.

This kind of transceiver use some form of switching among the antenna ele-ments or combining them to present the receiver with a single signal [24].

(23)

(24)

3

Direction-Finding Algorithm

A pre-determined hardware design for the beacon and the antenna constellation regarding the receiver are taken into consideration when an existing df algorithm will be investigated. Specifications of interest within this project are that the bea-con possess four transmitting antennas which are arranged in an uniform linear array, and the receiver possess one receiving antenna.

The interferometry df method described in Section 2.1.3 will serve as a design foundation for the algorithm due to the existents of an antenna array. Within the interferometry method, various techniques are available. These can be classified into conventional beamforming techniques, subspace-based techniques and max-imum likelihood techniques. The subspace-based technique was chosen after an evaluation regarding these methods, [31]. Within the subspace-based tech-nique, there are two main approaches which are called Multiple Signal Classifi-cation (music) and Estimation Parameters via Rotational Invariance Techniques (esprit). Here, the music approach was chosen, due to its proven stability and accurate results, [26].

First, a general case of the music algorithm will be presented, followed by a mod-ified development called M:1 music Algorithm, which is suited for this specific case.

(25)

12 3 Direction-Finding Algorithm

3.1 MUSIC Algorithm

The classic Multiple Signal Classification (music) algorithm is a high performance algorithm using a number of source signals and an array of receiving antennas. The algorithm achieves its high resolution from the evaluation of an input ance matrix derived from the input data model. The eigenvalues of this covari-ance matrix are determined and partitioned into two sets, the signal and the noise sub-spaces [24]. In general, the music algorithm acts as a simo system, thus the location of the antenna array is at the receiving end.

3.1.1 Uniform Linear Array

The following assumptions are being made for a mathematical model: The an-tenna array is located far from the signal sources such that far-field characteris-tics are applied. This means that the wavefront generated by each signal source arrives at all the array elements at an equal direction of propagation and the wave-front is planar. Each source is narrow band with the same center frequency, ω0, and the noise is assumed to be Gaussian white noise. The noise at every array ele-ment have a common standard deviation of σn. The number of signal sources are

D(d = 0, 1, .., D). The antenna array consists of M(m = 0, 1, .., M) elements and each element has the same characteristics. The antennas are aligned and equally spaced. [29]. An illustration of the scenario is displayed in Figure 3.1.

Figure 3.1:D signal sources and M receiving antennas.

Figure 3.2 illustrates the antenna array and one incoming signal. A signal gener-ated by signal source d are approaching the array with an angle θd. The signal

travel distance differs with a factor of ds for each antenna due to the antenna

(26)

3.1 MUSIC Algorithm 13

Figure 3.2:The receiving antenna array.

3.1.2 Single Source Model

To be able to evaluate the spectrum of a received signal, a general mathematical model for the input data is constructed, Equation (3.1). This function is based on the fact that the signal belongs to the narrow band spectrum.

x(t) =                   x0(t) x1(t) .. . x_M(t)                   =                   e−jζ0 e−jζ1 .. . e−jζM                   s(t) = a(θ)s(t) (3.1)

where the vector x(t) is called input data vector, a(θ) is called the steering vector and s(t) is the vector of incident signals derived from the narrow band signal model [27]. The a(θ) is a general form of the steering vector which is a function of the individual element response. The ζi is the phase shift factor that is the effect

of the propagation delay which occurs due to the placement of the antennas. To get a more realistic model of the input data, Equation (3.1) is rewritten into:

x(t) = a(θ)s(t) + n(t) (3.2) which can be expressed as

                  x0(t) x1(t) .. . x_M(t)                   =                   e−jζ0 e−jζ1 .. . e−jζM                                  s0(t) s1(t) . . . s_D(t)                +                   n0(t) n1(t) .. . n_M(t)                   (3.3)

(27)

where the vector n(t) represents the interference and noise components. As illus-trated in Figure 3.2, the incoming plane wave travels a longer distance to reach the mth_{antenna element, than it does to reach the first element. The difference}

in distance between two elements is ds = dasin θ, and the corresponding phase

shift differs between the elements with respect to the first element. Between the first and second element, the phase shift is ζ1= βds, and for the mthelement, the

phase shift is:

ζ_m= mβds = mβdasin θ (3.4)

where the β is a propagation factor equal to 2π/λ where λ is the wavelength. By using Equation (3.5), the a(θ) vector can be rewritten into:

a(θ) =                   1 e−jζ1 .. . e−jζM                   (3.5)

3.1.3 General Model

The input data model in Equation (3.2) can be extended into a more general case of multiple transmitting signals [27]:

x(t) = A(Θ)s(t) + n(t) (3.6) The A(Θ) in Equation (3.6) can be expanded into:

A(Θ) =                   1 1 . . . 1

e−jβdasin θ0 _e−jβdasin θ1 _{. . .} _e−jβdasin θD

..

. ... . .. ... e−jβMdasin θ0 _e−jβMdasin θ1 _{. . . e}−jβMdasin θD

                  (3.7)

To better understand how the A(Θ) matrix is composed, Figure 3.3 illustrates how the transmitting signals is received by the antenna array.

(28)

Figure 3.3:Illustration of the transmitting signals and the receiving array.

Each row of the A(Θ) matrix in Equation (3.7) is derived from how one receiv-ing antenna receives signals from each transmitter. Figure 3.4 illustrates how the first (m = 0) receiving antenna receives signals from each transmitter.

Figure 3.4:One receiver receives signals from different transmitters. Due to the considerable distance between each transmitter, the aoa of each signal in Figure 3.4 is unique. Each column of the A(Θ) matrix in Equation (3.7) is derived from how each receiving antenna receives a signal from the same trans-mitter, illustrated in Figure 3.5.

(29)

Figure 3.5:One transmitter sending signals for the receiver array to receive.

The signal traveling distance is considerably larger than the antenna spacing, d_a, which results in the aoa being assumed to be equal. The A(Θ) matrix in Equation (3.7) is usually denoted as:

A(Θ) =ha(θ0) . . . a(θD)

i

(3.8)

For the vector x, the covariance matrix R_xxcan be expressed as:

Rx,x= E[xxH] = E[(As+n)(As+n)H] = AE[ssH]AH+ E[nnH] (3.9)

where E[ssH_{] = R}

ss is the signal correlation matrix and E[nnH] = σn2Ithe noise

correlation matrix, rewriting Equation (3.9) to:

Rx,x= ARssAH+ σn2I (3.10)

The vector space of R_x,xspans one set composed by noise and one set composed by both noise and signal. It is the latter that is interesting for further calculations. In practical applications, the sample covariance ˆR_xxis usually used [29]:

ˆR_x,x= 1 N N X i=1 x(ti)x(ti)H (3.11)

where ˆR_x,xis the maximum likelihood estimation of R_x,xand N is the number of samples obtained at some time instances.

Equation (3.10) has positive eigenvalues that occur D times that correspond to signals and small (close to zero) eigenvalues that occur NN = M − D times,

(30)

belong to either of the two orthogonal subspaces, the signal subspace or the noise subspace. The steering vectors of A(Θ) lie in the signal subspace and the noise eigenvalues, V_N lie in the noise subspace. Due to orthogonality of the two sub-spaces, it is possible to search though all possible array steering vectors to find those which are orthogonal to the noise eigenvectors. By doing this its possible to determine the angle of arrival, θ.

To get an estimation with this method, the number of transmitters must be greater than the number of receivers in the antenna array;

(31)

3.2 M:1 MUSIC Algorithm

A modified version of the music algorithm, suited for the specifications of this thesis, is presented in this chapter.

3.2.1 Uniform Linear Array

The considered transmitter contains the antenna array instead of the receiving end of the system, illustrated in Figure 3.6. As described in Section 2.2, it is the aod that will be determined in this case. The antenna array consists of four antennas, which are being switched between. When the switch is connected to an antenna, it sends out a signal during a known period of time, illustrated in Figure 3.7. The spacing between the transmitting antennas is daand it is assumed

that the transmitter is at a distance from the receiver so that the aod from each transmitter antenna is constant during one switching period.

Figure 3.6:Antenna constellation for the system.

The signals sent from each transmitting antenna are illustrated in Figure 3.7 and the four signals form one switching period. The appearance and properties of the signals are purely for illustrative purpose.

The receiver receives a combined signal of the four signals, and due to the antenna spacing daon the transmitter, the signals from respectively transmitting

antenna gets an unique phase shift. As illustrated in Figure 3.8. The received signal with respect to the phase can be described as:

x(t) = e−jζ(t)_{s(t) + n(t)} _(3.13) where

ζ(t) = ζ0u0T(t) + ζ1u2T_T (t) + ζ2u3T_2T(t) + ζ3u_3T4T(t) (3.14) and ut1

t0(t) is a rectangle function defined as:

ut1 t0(t) =        1, if t0< t ≤ t1 0, otherwise

(32)

3.2 M:1 MUSIC Algorithm 19

Figure 3.7:Signals sent from each TX-antenna.

Figure 3.8:Illustrative signal received by the RX-antenna.

In the original music algorithm D represents the number of source signals. In the M:1 music algorithm, it represents one switching period which means that D = 1 for this system. Instead of looking at the number of receiving antennas, M here represents the number of transmitting antennas.

3.2.2 Model

The general case of the input data model is equal to the expression in Equa-tion (3.2), except from the notaEqua-tion of a(θ) that is denoted a(ζ) due to the constant DoA. Here, s(t) is a scalar since it is one signal and thus D = 1.

(33)

which can be expressed as                x(t) x(t + T ) x(t + 2T ) x(t + 3T )                =                 e−jζ0 e−jζ1 e−jζ2 e−jζ3                 s(t) +                n(t) n(t + T ) n(t + 2T ) n(t + 3T )                , 0 ≤ t ≤ T (3.16)

where the noise, n(t) can be rewritten into:

n(t) =                n(t) n(t + T ) n(t + 2T ) n(t + 3T )                ∆ =                n0(t) n1(t) n2(t) n3(t)                (3.17)

due to the fact that ni _{is independent noise. This result in:}

               x(t) x(t + T ) x(t + 2T ) x(t + 3T )                =                 e−jζ0 e−jζ1 e−jζ2 e−jζ3                 s(t) +                n0(t) n1(t) n2(t) n3(t)                , 0 ≤ t ≤ T (3.18)

Between the first and second element the phase shift is ζ1 = βdasin θ and for the

mthelement the phase shift is:

ζ_m= mβd_asin θ (3.19) where the β is the propagation factor, Equation (3.5).

There are five steps to implement the M:1 music algorithm. The process is described in Algorithm 1.

(34)

3.2 M:1 MUSIC Algorithm 21

Algorithm 1M:1 music algorithm

1. Collect N (n = 0, 1, .., N ) input samples to form the array x in Equa-tion (3.15) and estimate the input covariance matrix ˆR_x,x using Equa-tion (3.11).

2. Determine the eigenvectors of ˆR_x,xusing eigenvalue decomposition: ˆRx,xV= VΛ (3.20)

where Λ = diag[λ0, λ1, .., λM], λ0 ≤ λ1 ≤, .., ≤ λM are the eigenvalues of

ˆRx,xand V = [v0, v1, .., vM] are the corresponding eigenvectors

[27].

3. Evaluate the music function, P(θ):

P(θ) = 1 aH_(ζ)V NVHNa(ζ) ,−π 2 ≤θ ≤ π 2 (3.21) where V_N = [v_D+1, v_D+2, .., v_M]. The product V_NVH_N represents the pro-jection matrix on the subspace. Orthogonality between a(θ) and VN will

minimize the denominator.

4. Find the D maximum peaks of P(θ). This occurs when the denominator is at its minimum. [20] [25]

(35)

(36)

4

Direction-Finding Simulations

This chapter presents the simulation results for the 4:1 music algorithm. To fully understand the impact of different system design aspects, simulations with varying values of parameters have been performed. The evaluated parameters is the number of samples taken for each phase, N , antenna spacing for the beacon device, d_a, and the number of transmitting antennas, M.

4.1 Algorithm Simulations

Calculations for this system are intended to be performed in the receiver, as de-scribed in Section 1.3. The microcontroller within the receiver is pre-determined to be the Texas Instrument CC2650 radio chip. The incoming signal frequency is 2.45 GHz, and the microcontroller has a sample rate of 200 kbit/s [4], which indicates a signal downconvertion to a maximum of 100 kHz [18]. It is this fre-quency that are used in the df algorithm. The simualtions are being made with the software tool MATLAB.

4.2 Simulation with Pre-Determined Parameters

Simulations of the 4:1 music algorithm is presented in Figure 4.1, the parame-ters for these are set to match the pre-determined values. The graph to the left shows a signal being recognized by the algorithm at -70° and the rightmost graph presents a signal being recognized at 10°.

(37)

24 4 Direction-Finding Simulations

Figure 4.1:Simulations with pre-determined parameters.

The number of samples per phase is set to 100, the antenna spacing is λ/4, number of transmitting antennas is 4 and the snr noise level is 20 dB. The noise is ideal white Gaussian noise. The music function does not estimate the signal power associated with each recognized angle. Instead, the peaks of P(θ) corre-spond to the true angles of departure.

4.3 Number of Samples

The simulations for varying number of samples per phase are illustrated in Figure 4.2. Other conditions remain unchanged. The graph to the left represents a aod of -70° and the right represents a aod of 10°.

(38)

4.4 Antenna Spacing 25

For the solid line in the figure, the number of samples are 50 which provides the highest noise floor but requires the least amount of processing data. The dashed line represents 100 samples and the dash-dotted line represents 300 sam-ples. The number of samples is a software design aspect that can be varied to obtain a suitable result.

When the number of samples are varied, the highest number of performed sam-ples provide the best angle detection. But a high number of samsam-ples result in a high amount of processing data. Figure 4.2 proves that an angle detection is pos-sible with a considerably low number of samples. Which also result in a lower amount of computation resources.

4.4 Antenna Spacing

The simulations for varying antenna spacing are illustrated in Figure 4.3. Other conditions remain unchanged compared to the basic simulation. The graph to the left represents a aod of -70° and the right represents a aod of 10°.

Figure 4.3:Varying antenna spacing.

The solid line represents an antenna spacing of λ/4, the dashed line repre-sents an antenna spacing of λ/2 and the dash-dotted line reprerepre-sents a spacing of λ.

Figure 4.3 shows that the antenna spacing is an important hardware design as-pect. When the antenna spacing are half the wavelength or less, the beam width becomes narrow. And when the antenna spacing is more than half the wave-length, a false peak will emerge. This phenomenon is created due to the angle of departure condition of −π

2 ≤ θ ≤ π2. For each aod angle θ there is a corre-sponding phase shift ζ. In order to determine θ uniquely from ζ, a one-to-one-correspondence is desired between them. This results in a phase shift condition of −π ≤ ζ ≤ π and thus

(39)

26 4 Direction-Finding Simulations da,max = ζ_maxλ 2π sin θ_max (4.1) which yields d_a ≤ λ₂ (4.2) So when da exceeds this value there will be an ambiguity in the detection.

4.5 Number of Transmitting Antennas

The simulations for varying number of transmitting antennas are illustrated in Figure 4.4. Other conditions remain unchanged compared to the basic simulation. The graph to the left represents a aod of -70° and the right represents a aod of 10°.

Figure 4.4:Varying number of transmitting antennas.

The solid line represents 2 transmitting antennas, the dashed line represents 4 antennas, which conforms with the actual hardware design. The dash-dotted line represents 8 antennas.

Figure 4.4 shows that the aod estimation beam width becomes narrower for a higher amount of antennas. By increasing the number of antennas, the signal gets easier to distinguish but the processing data are also increasing for every added antenna. This would also be a practical aspect, though every antenna requires a certain amount of spacing on the beacon device.

(40)

4.6 Noise Level 27

4.6 Noise Level

The simulations for varying noise level are illustrated in Figure 4.5. Other condi-tions remain unchanged compared to the basic simulation. The graph to the left represents a aod of -70° and the right represents a aod of 10°.

Figure 4.5:Varying noise level.

The solid line in the figure represents a snr of 5 dB, the dashed line represents a snr of 20 dB and the dash-dotted line represents a snr of 35 dB. The snr is defined as:

snr= Psignal

P_noise (4.3)

where Psignal is the power of the signal and Pnoiseis the power of the background

noise. The noise of the system affects the peak definition, and a high snr value proves to have less of an impact on the system.

(41)

(42)

5

Localization Background

There are many different techniques to wirelessly locate objects with unknown positions. This section describes two methods, bearing-range and bearings-only. An estimation algorithm is presented in which these techniques will be used in Chapter 6.

5.1 Bearing-Range Localization

When both bearing and range measurements are available for position estimation of a target, it is called bearing-range localization. The bearing measurement are the aod data calculated in the df algorithm and the range data could be obtained by using the signal strength of a received signal. To calculate the received signal strength, Friis transmission formula is utilized:

P_r = PtGtGrλ 2

4πR2 (5.1)

where Ptis the power at which the signal was transmitted, λ is the wavelength

and R is the range between the transmitter and receiver. The Gt and Gr are the

gains of the transmitter and receiver respectively [7]. With bearing-range mea-surements, one beacon is sufficient to estimate a position of a target.

5.2 Bearing-Only Localization

Bearing-only localization is a technique used to determine the location of a target by using bearings measurements. The bearings are calculated using:

θ = tan−1 Ry Rx

!

(5.2)

(43)

30 5 Localization Background

where Ryis the y-term range from the beacon to the target, and Rxis the x-term

range from the beacon to the target [5].

When the line of bearings are intersecting, it is possible to calculate the position of the target. With bearing-only measurement, one beacon would be insufficient to determine a position. Here, at least two beacons are required to get the inter-secting point.

5.3 Extended Kalman Filter

To estimate the location of the moving target, an estimation filter called the Ex-tended Kalman Filter (ekf) will be used. The ekf is a non-linear derivation of the linear Kalman Filter (kf) which optimizes the state estimate by minimizing the state covariance [30]. A continuously estimation of the state is conducted by us-ing a non-linear motion and measurement model of the system. To linearize the system, Taylor expansion is performed around the current state estimate. State space models of the motion and measurements are given by:

x_t= f_t−1(x_t−1) + u_t (5.3)

y_t= h_t(x_t) + v_t (5.4) where x_t is the motion state which contains position and velocity information. This model calculates a prediction of how the state change over time when the previous state, x_t−1, is taken into consideration. Equation (5.4) is the measure-ment model. The u_t-vector represents the process noise, which is assumed to be a zero mean Gaussian white noise with covariance Q. Further, the measurement noise, v_t, is assumed to be zero mean Gaussian white noise with covariance R. The filter involves two stages, prediction and measurement update. Equations for the prediction stage are:

ˆx_t|t−1= f_t−1(ˆx_t−1|t−1) (5.5)

P_t|t−1= F_t−1P_t−1|t−1FT_t−1+ Q_t−1 (5.6) where ˆx_t|t−1is the predicted state estimate and P_t|t−1is the corresponding covari-ance. Equations for the measurement stage of the filter are given by:

P_t|t_{= (I − K}_tH_t)P_t|t−1 (5.9) where ˆx_t|t is the state estimation corrected by the measurements and P_t|t is the updated covariance. By utilize the difference in the received and predicted mea-surements, and the modeled uncertainty, the state prediction from the prediction

(44)

5.3 Extended Kalman Filter 31

stage can be corrected. Estimated states are updated to get an estimation closer to the true state.

Further, the F_tand H_tare defined as:

Ht= ∂ht(xt ) ∂x_t _x t= ˆxt|t−1 (5.10) F_t= ∂ft(xt) ∂x_t _x t= ˆxt|t (5.11) where H_t and F_t are defined as the jacobians of the motion and measurement model respectively, Equation (5.3) - (5.4) [21]. The ekf is summarized in Algo-rithm 2.

Algorithm 2EKF

1. Initialize with ˆx_0|0= x0and ˆP0|0= P0. 2. Predicition update: ˆx_t|t−1= f_t−1(ˆx_t−1|t−1) P_t|t−1= F_t−1P_t−1|t−1FT_t−1+ Q_t−1 3. Measurement update: K_t= P_t|t−1HT_t(H_tP_t|t−1HT_t + R_t)−1 ˆx_t|t= ˆx_t|t−1+ K_t(y_t_{− h}_t(ˆx_t|t−1)) P_t|t_{= (I − K}_tH_t)P_t|t−1

(45)

(46)

6

Localization Simulations

To be able to test the ekf, data for both the target trajectory and the measure-ment has to be generated. Environmeasure-ment settings for the beacons will be included. Two simulation models are taken into consideration, one with measurements gen-erated from two beacons and one with measurements gengen-erated from three bea-cons. These models is based on the results given by the df algorithm presented in Chapter 3.2 and the theory introduced in Chapter 5.

6.1 Setup Models

By using the motion model, an actual trajectory is generated, which is repre-sented as the dashed line in Figure 6.1. Initialized data for the actual trajectory is set to a constant velocity in both x and y-direction. In x-direction,

v_x,t = 1 m/s, and in y-direction, vy,t = 0 m/s, which result in a horizontal path

with respect to the x-axis. The general motion model are given as:

x_t=                p_x,t p_y,t vx,t v_y,t                + u_t (6.1)

where px,t and py,t represents the position of the target, and vx,t and vy,t

repre-sents the velocity of the target. The u_t-vector represents the process noise, which is assumed to be a zero mean Gaussian white noise with a covariance of 0.5 me-ters. The initial motion model are given as:

(47)

34 6 Localization Simulations x0 =                0 0 1 0                + u₀ (6.2)

The first model, Model 1, is simulated using two beacons, and the second model, Model 2, is simulated using three beacons, Figure 6.1. The generated trajectory is equal in both Model 1 and Model 2.

(48)

6.1 Setup Models 35

Figure 6.1:Model 1 and Model 2.

To generate the measurement data, a measurement vector is utilized which consists of positioning data with added noise. The positions of the beacons are determined and are represented as the wedges in Figure 6.1. These positions are known, denoted as (bx1, by1), (bx2, by2) and (bx3, by3) for beacon 1 (b1), beacon 2 (b2) and beacon 3 (b3) respectively. The arrows in Figure 6.1 represent the direction reference for the aod:s. This reference is set to match the aod derived from the df algorithm. It is assumed that the beacons are placed aligned on a wall.

(49)

36 6 Localization Simulations

6.2 Bearing-Range

The setup model, measurement vectors and simulations of the bearing-range method is presented in this section.

6.2.1 Bearing-Range Model

When bearing-range is used together with the ekf, the measurement parameters for the system are the aod values and the distance between each beacon and the target, Figure 6.2. The positions of the target are denoted as p = (px, py).

Figure 6.2: Measurement parameters for Model 1 and Model 2, using the bearing-range method.

In the ekf, the measurement vector, y_t, will consist of four parameters for Model 1, denoted y_BR1,t. This is due to the fact that two beacons are used in this case, and each beacon results in two statements, one aod value and one distance value. Equation (6.3) is the specified measurement vector for this model.

y_BR1,t= h_BR1(xt) + vt=                d_1,t θ1,t d_2,t θ_2,t                +                v_d vθ v_d v_θ                (6.3)

where h_BR1(xt) is given by:

hBR1(xt) =                     ||p − b1|| tan−1py−by1 px−bx1 ||p − b2|| tan−1py−by2 px−bx2                     (6.4) which yields

(50)

6.2 Bearing-Range 37 H_BR1,t =∂hBR1,t(xt) ∂x_t =                         px−bx1 √ (px−bx1)2+(py−by1)2 py−by1 √ (px−bx1)2+(py−by1)2 0 0 by−py1 (bx1−px)2+(by1−py)2 px−bx1 (bx1−px)2+(by1−py)2 0 0 px−bx2 √_(p x−bx2)2+(py−by2)2 py−by2 √_(p x−bx2)2+(py−by2)2 0 0 by−py2 (bx2−px)2+(by2−py)2 px−bx2 (bx2−px)2+(by2−py)2 0 0                         (6.5)

The measurement vector will consist of six parameters when using Model 2, de-noted y_BR2,t. Equation (6.9) is the specified measurement vector for this model.

y_BR2,t= h_BR2(xt) + vt=                           d_1,t θ1,t d_2,t θ2,t d_3,t θ_3,t                           +                           v_d vθ v_d vθ v_d v_θ                           (6.6)

where h_BR2(x_t) is given by:

h_BR2(xt) =                                  ||p − b1|| tan−1py−by1 px−bx1 ||p − b2|| tan−1py−by2 px−bx2 ||p − b3|| tan−1py−by3 px−bx3                                  (6.7) which yields H_BR2,t =∂hBR2,t(xt) ∂xt =                                        px−bx1 √_(p x−bx1)2+(py−by1)2 py−by1 √_(p x−bx1)2+(py−by1)2 0 0 by−py1 (bx1−px)2+(by1−py)2 px−bx1 (bx1−px)2+(by1−py)2 0 0 px−bx2 √ (px−bx2)2+(py−by2)2 py−by2 √ (px−bx2)2+(py−by2)2 0 0 by−py2 (bx2−px)2+(by2−py)2 px−bx2 (bx2−px)2+(by2−py)2 0 0 px−bx3 √_(p x−bx3)2+(py−by3)2 py−by3 √_(p x−bx3)2+(py−by3)2 0 0 by−py3 (bx3−px)2+(by3−py)2 px−bx3 (bx2−px)2+(by3−py)2 0 0                                        (6.8)

The measurement noise, vd and vθ, has a standard deviation of 0.5 meters

re-spectively 0.05 radians. These vectors represent the measured target positions in polar coordinates with respect to each beacon.

(51)

6.2.2 Bearing-Range Simulations

Each model is simulated with a varying beacon spacing. This is to distinguish differences in the localization estimation. Values of the beacon spacing are 4 and 12 meters. In the bearing-range simulation, the noise used for the measurements are normally distributed pseudorandom numbers. The noise corresponding to the range has a standard deviation of 0.50 meters and the noise corresponding to the aod has a standard deviation of 2.86° (0.05 radians). Figure 6.3 presents the simulations performed for Model 1 with a beacon spacing of 4 meters.

Figure 6.3: Model 1 with a distance of 4 meters between each beacon, using bearing-range measurements.

The dashed lines represent the actual trajectory, in the figure to the left, the solid line represent the ekf estimation. To the right, the dots represents the mea-surements.

Figure 6.4:Angle and distance values for Model 1 with a beacon distance of 4 meters, using bearing-range measurements.

(52)

6.2 Bearing-Range 39

Figure 6.4 presents the values for the angles and the distances for each bea-con during the tracking with the beabea-con placement equal the one in Figure 6.3. The measurements are represented by the solid lines and the actual values are represented by the dashed lines. The thick lines correspond to the first beacon, b1 in Figure 6.2, and the narrow lines correspond to the second beacon, b2 in Figure 6.2. The measurements i Figure 6.4 are symmetric between beacon 1 and beacon 2. This is due to the placement of the beacons. Figure 6.5 presents the simulations performed for Model 1 with a beacon spacing of 12 meters.

Figure 6.5:Model 1 with a distance of 12 meters between each beacon, using bearing-range measurements.

Figure 6.6 presents the values for the angles and the distances for each beacon during the tracking with the beacon placement equal the one in Figure 6.5.

(53)

The measurements in Figure6.6 are more stable than the measurements in Figure 6.4, this is because the estimation filter has a longer period of time to stabilize here.

Figure 6.7 presents the simulations performed for Model 2 with a beacon spacing of 4 meters.

Figure 6.7: Model 2 with a distance of 4 meters between each beacon, using bearing-range measurements.

(54)

6.2 Bearing-Range 41

Figure 6.8 presents the values for the angles and the distances for each beacon during the tracking with the beacon placement equal the one in Figure 6.7. Mea-surements and true values for the third beacon, b3in Figure 6.2, are represented by the lines in magenta. The thick and the magenta lines are similar to the graph in Figure 6.6.

Figure 6.9 presents the simulations performed for Model 2 with a beacon spacing of 12 meters. Figure 6.3, 6.5, 6.7 and Figure 6.9 shows that the estimation filter has an initial uncertainty but find the right path quickly. The graphs in these figures displays a small instability in between the beacons.

Figure 6.9:Model 2 with a distance of 12 meters between each beacon, using bearing-range measurements.

Figure 6.10: Angle and distance values for Model 2 with a beacon distance of 12 meters, using bearing-range measurements.

Figure 6.10 presents the values for the angles and the distances for each bea-con during the tracking with the beabea-con placement equal the one in Figure 6.9.

(55)

Measurements and true values for the third beacon, b3 in Figure 6.2, are repre-sented by the lines in magenta.

6.3 Bearing-Only

The setup model, measurement vectors and simulations of the bearing-only method is presented in this section.

6.3.1 Bearing-Only Model

When the bearing-only method is used together with the ekf, the measurement parameters for the system are the aod values, Figure 6.11.

Figure 6.11: Measurement parameters for Model 1 and Model 2, using the bearing-only method.

In the ekf, the measurement vector will consist of two parameters for Model 1, denoted y_BO1,t. One aod value for each beacon. Equation (6.9) is the specified measurement vector for this model.

y_BO1,t = h_BO1(xt) + vt =

"θ1,t

θ_2,t #

+ v_θ (6.9)

where h_BO1(xt) is given by:

h_BO1(x_t) =           tan−1py−by1 px−bx1 tan−1py−by2 px−bx2           (6.10) which yields

(56)

6.3 Bearing-Only 43 H_BO1,t =∂hBO1,t(xt) ∂x_t =           px−bx1 √ (px−bx1)2+(py−by1)2 py−by1 √ (px−bx1)2+(py−by1)2 0 0 px−bx2 √_(p x−bx2)2+(py−by2)2 py−by2 √_(p x−bx2)2+(py−by2)2 0 0           (6.11)

The measurement vector will consist of six parameters when using Model 2, de-noted y_BO2,t. Equation (6.12) is the specified measurement vector for this model.

y_BO2,t = h_BO2(xt) + vt=           θ_1,t θ2,t θ_3,t           + v_θ (6.12)

where h_BO2(xt) is given by:

h_BO2(xt) =                  tan−1py−by1 px−bx1 tan−1py−by2 px−bx2 tan−1py−by3 px−bx3                  (6.13) which yields H_BO2,t =∂hBO2,t(xt) ∂xt =                   px−bx1 √ (px−bx1)2+(py−by1)2 py−by1 √ (px−bx1)2+(py−by1)2 0 0 px−bx2 √_(p x−bx2)2+(py−by2)2 py−by2 √_(p x−bx2)2+(py−by2)2 0 0 px−bx3 √ (px−bx3)2+(py−by3)2 py−by3 √ (px−bx3)2+(py−by3)2 0 0                   (6.14)

The measurement noise, v_θ, has a standard deviation of 0.05 radians. These vec-tors represent the measured positions of the target in polar coordinates with re-spect to each beacon. When the polar coordinates are converted to cartesian, the following equation is used:

di = ||p − bi|| (6.15)

(57)

6.3.2 Bearing-Only Simulations

Each model is simulated with a varying beacon spacing. This is to distinguish differences in the localization estimation. Values of the beacon spacing are 4 and 12 meters. In the bearing-only simulation, the noise used for the measurements are normally distributed pseudorandom numbers. The noise corresponding to the aod has a standard deviation of 2.86° (0.05 radians). Figure 6.12 presents the simulations performed for Model 1 with a beacon spacing of 4 meters.

Figure 6.12:Model 1 with a distance of 4 meters between each beacon, using bearing-only measurements.

Figure 6.13: Angle and distance values for Model 1 with a beacon distance of 4 meters, using bearing-only measurements.

(58)

6.3 Bearing-Only 45

Figure 6.13 presents the values for the angles and the distances for each bea-con during the tracking with the beabea-con placement equal the one in Figure 6.12. The measurements are represented by the solid lines and the actual values are represented by the dashed lines. The thick lines correspond to the first beacon, b1 in Figure 6.2, and the narrow lines correspond to the second beacon, b2 in Figure 6.2. This graph is similar to the leftmost graph in Figure 6.4, just with a bigger uncertainty.

Figure 6.14 presents the simulations performed for Model 1 with a beacon spac-ing of 12 meters.

Figure 6.14: Model 1 with a distance of 12 meters between each beacon, using bearing-only measurements.

(59)

Figure 6.7 presents the simulations performed for Model 2 with a beacon spac-ing of 4 meters.

Figure 6.16:Model 2 with a distance of 4 meters between each beacon, using bearing-only measurements.

The dashed lines represent the actual trajectory, in the figure to the left, the solid line represent the ekf estimation. To the right, the dots represents the mea-surements. Figure 6.17 presents the values for the angles and the distances for each beacon during the tracking with the beacon placement equal the one in Fig-ure 6.16.

Measurements and true values for the third beacon, b3in Figure 6.2, are rep-resented by the lines in magenta.

(60)

6.3 Bearing-Only 47

spacing of 12 meters. The simulations for the bearing-only estimation presents a bigger uncertainty in the initial phase than in the bearing-range estimation. In Figure 6.12, 6.14, 6.16 and Figure 6.18, the uncertainty in between the beacons gets more difficult to distinguish due to this.

Figure 6.18: Model 2 with a distance of 12 meters between each beacon, using bearing-only measurements.

(61)

6.4 Monte Carlo Simulations

Monte Carlo diagrams of the estimation error are presented, with a beacon spac-ing of 4 meters as well as 12 meters. The graphs represents the mean error, µ, and the standard deviation, σ, over time. In every timestep of the position estimation, 1000 samples has been collected and then been averaged for the mean value and a standard deviation has been calculated.

Notation and corresponding line representation are consistent throughout this chapter. The mean error for bearing-range estimations are denoted as µBR and

corresponds to the solid thin line. The mean error for bearing-only estimations are denoted as µBO and corresponds to the solid thick line. Standard deviations

for the bearing-range estimations are denoted as σBR and corresponds to to the

dashed thin line. The standard deviation for the bearing-only estimations are denoted σBOand corresponds to the dashed thick line.

6.4.1 Model 1

In Model 1, two beacons are utilized. Figure 6.20 presents the mean estimation error and the standard deviation with a beacon spacing of 4 meters. The graph to the left corresponds to the x-component and the rightmost graph corresponds to the y-component.

Figure 6.20: Monte Carlo diagram of the estimation error for a beacon dis-tance of 4 meters, Model 1.

Figure 6.21 displays the mean estimation error and the standard deviation with a beacon spacing of 12 meters.

(62)

6.4 Monte Carlo Simulations 49

6.4.2 Model 2

Three beacons are utilized in Model 2, which means that these graphs are more comprehensive and makes it easier to distinguish estimation error patterns. Fig-ure 6.22 presents the mean estimation error and the standard deviation with a beacon spacing of 4 meters. The graph to the left corresponds to the x-component and the rightmost graph corresponds to the y-component.

Figure 6.21 displays the mean estimation error and the standard deviation with a beacon spacing of 12 meters.

(63)

Simulations for the estimation error indicates that the estimations get more accurate with a smaller beacon spacing, Figure 6.20 and Figure 6.22, in which a spacing of 4 meters is utilized. Evaluation of Figure 6.21 and Figure 6.23 displays a bigger uncertainty in between the beacons, where a beacon spacing of 12 me-ters is used. The error pattern in Figure 6.23 shows that the error for the x-term decreases when the target is close to a beacon, and that the error for the y-term increases when the target is close to a beacon. This can be explained by the us-age of a distance error, that does not take the x-and y-term into consideration separately. Figure 6.24 shows the distance error when the target is close and far, respectively.

Figure 6.24:Error evaluation. From the figure, the following can be stated:

(64)

6.4 Monte Carlo Simulations 51

e_x,close< e_{x,f ar} (6.16) ey,close> ey,f ar (6.17)

This causes a big uncertainty over time when looking at the Monte Carlo di-agram i Figure 6.23, but these errors are canceling each other out. The mean values of the estimated errors are the most stable with three beacons and when smaller beacon distances are utilized.

(65)