\addcontentsline{toc}{chapter}{Dictionary of acronyms} \chapter*{Dictionary of acronyms} \begin{itemize} \item \emph{ARFCN} - Absolute Radio Frequency Channel Number - The channel number specifies the physical frequency channel used for transmission and reception of radio waves inside of an BTS covered area. \item \emph{BTS} - Base Transceiver Station - \item \emph{DC} - Direct Current \item \emph{GNSS} - Global Navigation Satellite System - A satellite navigation system that allows a specialized receive to determine its location on Earth. \item \emph{PCB} - Printed Circuit Board - The board where electronic components are soldered onto and wired through conductive tracks. \item \emph{RRLP} - Radio Resource Location Protocol - The employed protocol in GSM, UMTS and other wireless networks for providing and exchange of geolocation information. \item \emph{SMA} - SubMiniature version A - SMA is a connector used for interconnecting coaxial cables or PCB electronics that work in the frequency range between 0-18 GHz. \item \emph{TRX} - \item \emph{UART} - Universal Asynchronous Receiver Transmitter - A serial communication interface used by computers or other peripheral devices to communicate. \item \emph{UMTS} - Universal Mobile Telecommunications System - Third generation mobile network based on the GSM standards. \end{itemize} \addchap{Appendix} \numberwithin{equation}{subsection} \numberwithin{table}{subsection} \captionsetup[figure]{list=no} \captionsetup[table]{list=no} \section{Installation and configuration guide} In order to evaluate the localization system, it is required to install OpenBSC and to modify the proper source files and compile the system. The aim of this section is to describe that process in such detail that the presented material is sufficient to reproduce equivalent or similar results. The guide was successfully tested out on the following operating systems: Ubuntu 10.04 LTS 64 bit and Ubuntu 12.04 LTS 64 bit. A self-bootable test USB system is supplied with the thesis and it can be evaluated without executing the given steps. There is a marking difference between text given in light and dark grey background color, the first ought to be typed in into the terminal window or it may be an output produced by an application, whereas the later emphasizes a file modification case. \subsection{Installation of OpenBSC} In order to compile OpenBSC it is required to install the following precompiled packages\footnote{If more details are required for the installation process a guide can be found at \citep{openbscInstall}.}: \begin{itemize}\addtolength{\itemsep}{-0.8\baselineskip} \item libdbi0 \item libdbi0-dev \item libdbd-sqlite3 \item libortp-dev \item build-essential \item libtool \item autoconf \item automake \item git-core \item pkg-config \end{itemize} Before installing the required packages and libraries, to keep the installation process clean and free of modifying other files, the author will create a new directory. \begin{lstlisting}[backgroundcolor=\color{light-gray}] mkdir gsm_localization cd gsm_localization \end{lstlisting} By executing the following instructions the required libraries will be installed. \begin{lstlisting}[backgroundcolor=\color{light-gray}] sudo apt-get install libdbi0-dev libdbd-sqlite3 build-essential sudo apt-get install libtool autoconf automake git-core sudo apt-get install pkg-config libortp-dev \end{lstlisting} After the packages were installed, \textit{libosmocore} library must be downloaded, compiled and installed. By executing the following instructions: \begin{lstlisting}[backgroundcolor=\color{light-gray}][numbers = none] git clone git://git.osmocom.org/libosmocore.git cd libosmocore autoreconf -fi ./configure make sudo make install sudo ldconfig cd .. \end{lstlisting} In the next step \textit{libosmo-abis} will be installed. \begin{lstlisting}[backgroundcolor=\color{light-gray}][numbers = none] git clone git://git.osmocom.org/libosmo-abis.git cd libosmo-abis autoreconf -fi ./configure make sudo make install sudo ldconfig cd .. \end{lstlisting} After the previous steps have finished successfully, the author will proceed with downloading, compiling and installing OpenBSC. \begin{lstlisting}[backgroundcolor=\color{light-gray}][numbers = none] git clone git://git.osmocom.org/openbsc.git cd openbsc/openbsc autoreconf -i sudo export PKG_CONFIG_PATH=/usr/local/lib/pkgconfig ./configure make \end{lstlisting} At this point, OpenBSC should be successfully compiled. \newpage \subsection{Configuring nanoBTS for OpenBSC} To enable the nanoBTS and OpenBSC to be fully operational, the last configuration steps have to be made. It is necessary to inform the nanoBTS of the IP address of the server that is running OpenBSC since it must connect to OpenBSC. We need to find a free ARFCN channel where our system is expected to operate\footnote{A licence has to be obtained from the Federal Network Agency (German: \textit{Bundesnetzagentur}), otherwise it is ilegal and may be considered as a criminal act.}. To find the ID and the IP address of the nanoBTS it is required to start \textit{ipaccess-find}\footnote{The nanoBTS ought to be blinking in orange color before starting \textit{ipaccess-find}.}. \begin{lstlisting}[backgroundcolor=\color{light-gray}][numbers = none] cd ~/gsm_localization/openbsc/openbsc/src/ipaccess ./ipaccess-find \end{lstlisting} \textit{ipaccess-find} will produce an output similar to the one given: \begin{lstlisting}[backgroundcolor=\color{light-gray}][numbers = none] Trying to find ip.access BTS by broadcast UDP... MAC_Address='00:02:95:00:61:70' IP_Address='132.230.4.63' Unit_ID='1801/0/0' Location_1='' Location_2='BTS_NBT131G' Equipment_Version='165g029_73' Software_Version='168a352_v142b30d0' Unit_Name='nbts-00-02-95-00-61-70' Serial_Number='00110533' \end{lstlisting} In the next step, the nanoBTS is informed of the OpenBSC IP address by typing the following commands (the first IP address belongs to the server running OpenBSC and the second to the nanoBTS): \begin{lstlisting}[backgroundcolor=\color{light-gray}][numbers = none] cd ~/gsm_localization/openbsc/openbsc/src/ipaccess ./ipaccess-config -o 132.230.4.65 132.230.4.63 -r \end{lstlisting} It is required to create the directory where the configuration file will be located and to modify the configuration file. \begin{lstlisting}[backgroundcolor=\color{light-gray}][numbers = none] sudo mkdir /usr/local/lcr cd ~/gsm_localization/openbsc/openbsc/doc/ cd examples/osmo-nitb/nanobts sudo cp openbsc.cfg /usr/local/lcr sudo vim /usr/local/lcr/openbsc.cfg \end{lstlisting} A free ARFCN channel can be found using a spectrum analyzer and by setting the frequency range to the GSM frequency band. One has to slide through the frequencies shown on the X-axis, and by looking at the Y-axis with appropriate frequency resolution\footnote{The frequency resolution must be set to $f_{CB}=200 \,\mathrm{kHz}$ or higher values for faster movement in the frequency spectrum.}, where the received power is represented\footnote{ Dependent of the manufacturer and settings of the spectrum analyzer, it can show signal amplitude, magnitude and power.}. By patiently observing the Y-axis it can be easily seen on the X-axis which channels are taken by other GSM service providers and which are free. The chosen channel ought to be peak free. Once a free frequency channel has been found, it is necessary to instruct the nanoBTS to operate in that frequency range. The line, numbered 58, has to be modified with the correct free ARFCN channel,in this case 877. \begin{lstlisting} arfcn 877 \end{lstlisting} The ARFCN channel value can be calculated using the given formula in \eqref{eq:arfcn}, where $f_{start}$ is the starting frequency of the uplink bandwitdh for DCS1800, $f_{CB}$ is the channel bandwidth and \textit{Offset} is the offset\footnote{ A table with frequency channels can be found at the following URL: \url{https://gsm.ks.uni-freiburg.de/arfcn.php}}. \begin{equation} \label{eq:arfcn} \centering \begin{array}{l} \displaystyle f_{up}(\mathrm{ARFCN}) = f_{start}+f_{CB}\cdot(\mathrm{ARFCN}-\mathrm{Offset}) \\ \displaystyle \\ \displaystyle where \left\{ \begin{array}{rcl} f_{start} & = & 1710.2 \,\mathrm{MHz} \\ f_{CB} & = & 200 \,\mathrm{kHz} \\ \mathrm{Offset} & = & 512 \end{array}\right. \end{array} \end{equation} %Multiple aligned equation %\begin{equation} %\label{eq:15} %\centering %where \left\{ \begin{array}{rcl} % f_{start} & = & 1710.2 \,\mathrm{MHz} \\ % f_{CB} & = & 200 \,\mathrm{KHz} \\ % \mathrm{Offset} & = & 512 %\end{array}\right. %\end{equation} On line numbered 53, the last configuration file modification has to be made for the final configuration of the OpenBSC software. The Unit ID from the output above has to be set\footnote{Indentation has to match the one of the configuration file.}. \begin{lstlisting} ip.access unit_id 1801 0 \end{lstlisting} At this point the nanoBTS and OpenBSC configuration is done. \newpage \subsection{Installation and configuration of RRLP assistance software} \label{sec:appendSoft} To install the RRLP software that generates assistance data, several libraries are required to be installed, \textit{cURL}\footnote{It may happen that the given download URLs are incorrect and have changed in the meantime, but one can easily find the latest versions on \url{http://curl.haxx.se/} and \url{http://www.hyperrealm.com/libconfig/}}, \textit{libconfig} and \textit{SQLite}. \textit{cURL} was used for the purpose of safely downloading assistance data from the Navigation Center of the US Coast Guard and Trimble server. \textit{libconfig} library is used for reading in the configuration file, this way compiling of the software whenever one changes the settings was avoided. The \textit{SQLite} library was employed to access the database used by OpenBSC to store the respondence data from the mobile stations. \begin{lstlisting}[backgroundcolor=\color{light-gray}][numbers = none] cd ~/gsm_localization sudo apt-get install libsqlite3-dev wget http://curl.haxx.se/download/curl-7.25.0.tar.gz wget http://www.hyperrealm.com/libconfig/libconfig-1.4.8.tar.gz tar -xvzf curl-7.25.0.tar.gz tar -xvzf libconfig-1.4.8.tar.gz cd curl-7.25.0 make sudo make install cd .. cd libconfig-1.4.8/ ./configure make sudo make install \end{lstlisting} Once the libraries have been successfully installed, the user may proceed with the configuration and compiling the RRLP assistance software, which is the key software produced in this thesis. The configuration file can be found in the same directory as the RRLP modules under the name: ``gnssrrlp.cfg''. The sample configuration file is already preconfigured for the location of ``Angewandte Mathematik und Rechenzentrum'' building. Latitude and longitude of the BTS are expressed in decimal degrees and are bounded by \textpm90\textdegree and \textpm180\textdegree respectively. Positive latitudes are north of the equator, whereas negative are south of the equator. It is alike for longitude coordinates, positive longitudes are east of Prime Meridian and negative are west of the Prime Meridian. If the position in decimal degrees of the BTS is unknown, it is straightforward to derive them using the formula given in \eqref{eq:dd}, where $D$ are degrees, $M$ are minutes and $S$ are seconds\footnote{An online converter of the Federal Communication Commission can be used as well to convert from degrees, minutes and seconds to decimal degrees and vice versa \url{http://transition.fcc.gov/mb/audio/bickel/DDDMMSS-decimal.html}}. \begin{equation} \label{eq:dd} \centering DD = D + \frac{M}{60} + \frac{S}{3600} \end{equation} The altitude may be left as it is, set to 0, since it is not used in the current measurement technique\footnote{If the value is set to zero, it is important to set it to 0.0 because \textit{libconfig} would otherwise convert it to an integer however it is a floating point number.}. The boolean variables can take \textit{true} or \textit{false} values. \begin{lstlisting} // An example configuration file for the GNSS RRLP software. name = "Configuration for RRLP"; // Change the settings if required: settings = {config = ( { ephemeris_url = "ftp://ftp.trimble.com/pub/eph/CurRnxN.nav"; almanac_url = "http://www.navcen.uscg.gov/ ?pageName=currentAlmanac&format=yuma"; latitude_of_BTS = 48.003601; longitude_of_BTS = 7.848056; altitude_of_BTS = 0.0; uncertainty_of_lat_long = 7; uncertainty_of_alt = 7; confidence_level = 0; ephemeris_repair = false; use_reference_time = false; extra_seconds_to_add = 7; timezone_of_BTS = 1; time_to_refresh_ephem = 1; time_to_refresh_alm = 1 ; } );}; \end{lstlisting} The uncertainty of the latitude and longitude correctness can be described using equation \eqref{eq:unclatlong} \citep{3gppequations}. The uncertainty of $r$ is expressed in meters, it defines how accurate is the specified location of the BTS. In the configuration file, $K$ is set to 7, which corresponds to $r$ = 9.4872 m. \begin{equation} \label{eq:unclatlong} \centering \begin{array}{l} \displaystyle r=C((1+x)^{K}-1)\;\; \displaystyle where \left\{ \begin{array}{rcl} C & = & 10 \\ x & = & 0.1 \\ K & \in & [0,127] \cap \mathbb{N}_{0} \end{array}\right. \end{array} \end{equation} A set of uncertainties $r$ is given in table \ref{tab:unclatlong} for various integer values of $K$. \begin {table}[ht] \caption{Example uncertainties (latitude and longitude) for various integer values of $K$} \label{tab:unclatlong}\centering %\rowcolor{2}{light-gray}{} \scriptsize\fontfamily{iwona}\selectfont \begin{tabular}{llll} \toprule %$D$&&$P_u$&$\sigma_N$\\ \textbf{Value of $K$}&\textbf{Value of uncertainty $r$}&\textbf{Value of $K$}&\textbf{Value of uncertainty $r$}\\\toprule 0 & 0 m&20 & 57.3 m\\\midrule 1 & 1 m&60 & 3.0348 km %2 & 2.1 m&100 & 137.8 km\\\midrule %3 & 3.3 m&- & - \\\bottomrule \end {tabular} \end {table} Altitude uncertainty can be described using the same Binomial expansion method, as given in \eqref{eq:uncalt}, however with altered constant values \citep{3gppequations}. The altitude uncertainty ranges between 0 m and 990.5 m ($h\in[0,990.5]\, \mathrm{m}$). Although the same constant name $K$ is used, it describes the altitude uncertainty. A set of uncertainties $h$ is given in table \ref{tab:uncalt} for various integer values of $K$. \begin{equation} \label{eq:uncalt} \centering \begin{array}{l} \displaystyle h=C((1+x)^{K}-1) \;\; \displaystyle where \left\{ \begin{array}{rcl} C & = & 45 \\ x & = & 0.025 \\ K & \in & [0,127] \wedge \|K\| \end{array}\right. \end{array} \end{equation} \begin {table}[] \caption{Example uncertainties (altitude) for various integer values of $K$} \label{tab:uncalt}\centering %\rowcolor{2}{light-gray}{} \scriptsize\fontfamily{iwona}\selectfont \begin{tabular}{llll} \toprule %$D$&&$P_u$&$\sigma_N$\\ \textbf{Value of $K$}&\textbf{Value of uncertainty $h$}&\textbf{Value of $K$}&\textbf{Value of uncertainty $h$}\\\toprule 0 & 0 m&20 & 28.74 m\\\midrule 1 & 1.13 m&60 & 152.99 m %2 & 2.28 m&100 & 486.62 m\\\midrule %3 & 3.46 m&- & - \\\bottomrule \end {tabular} \end {table} Confidence level can take any integer value between 0 and 127. The confidence level defines the percentage of the confidence that the target entity, the GSM user one wants to locate, is within the geometric shape defined earlier. A value between 0 and 100; 127 may be interpreted as ``no information'' \citep{3gppequations}. Ephemeris repair is a variable of the boolean type. Ephemeris data may contain errors or miss some satellite information \citep{NASA-Ephem-Errors} \citep{Stanford-Ephem-Errors} and the ephemeris repair function, if set to true, will take data of the previous measurement report. This introduces an error as well. Reference time can be used to provide extra information for the A-GPS in the MS of target entity. This field is of boolean type, if set to true, reference time is included in the sent packets. Since the sent packets are not transmitted in real time but put on a stack and then sent to the MS, a time delay exists. The reference time being sent to the MS is Coordinated Universal Time (UTC). The GPS device receives UTC time from the satellites and adjusts the computer time. To set the correct time, time zone offset of the BTS ought to be set correctly. Finally, the refresh time of downloading new almanac and ephemeris data has to be set. The variable uses the hour unit, how often the data are being refreshed and downloaded. The almanac data are valid for up to 180 days \citep{GPS-Guide} but are updated usually every day\footnote{Almanac update times can be found here: \url{http://www.navcen.uscg.gov/?pageName=currentNanus&format=txt}} \citep{GPS-Pentagon}. \clearpage \section{Troubleshooting the BTS} While the work has been performed on OpenBSC, to open a data channel (SDCCH), the BTS was sometimes sent in erroneous states. These states are reported through a LED light on the BTS. Based on the color and flash type of the LED one can find out the state of the BTS. These states are given in table \ref{tbl:LEDStatus} with their appropriate meaning. They may help the developer to troubleshoot and find the bug. \begin {table}[ht] \caption{Indicator LED status on the nanoBTS. Table courtesy of \citep{installnanoBTS}.} \label{tbl:LEDStatus}\centering %\rowcolor{2}{light-gray}{} \scriptsize\fontfamily{iwona}\selectfont \begin{tabular}{llll} \toprule %$D$&&$P_u$&$\sigma_N$\\ \textbf{State}&\textbf{Color \& Pattern}&\textbf{When}&\textbf{Precedence}\\\toprule Self-test failure&Red - Steady &In boot or application code when a power&1 (High) \\ &&on self-test fails\\\midrule Unspecified failure&Red - Steady &On software fatal errors&2\\\midrule No ethernet&Orange - Slow flash &Ethernet disconnected&3\\\midrule Factory reset&Red - Fast blink &Dongle detected at start up and the&4\\ &&factory defaults have been applied\\\midrule Not configured&Alternating Red/&The unit has not been configured&5\\ &Green Fast flash\\\midrule Downloading code&Orange - Fast flash &Code download procedure is in progress&6\\\midrule Establishing XML&Orange - Slow blink &A management link has not yet been established&7\\ &&but is needed for the TRX to become operational.\\ &&Specifically: for a master a Primary OML or\\ &&Secondary OML is not yet established; for a\\ &&slave an IML to its master or a Secondary \\ &&OML is not yet established. \\\midrule Self-test &Orange - Steady &From power on until end of backhaul&8\\ &&power on self-test\\\midrule NWL-test &Green - Fast flash& OML established, NWL test in progress&9\\\midrule OCXO Calibration &Alternating Green/& The unit is in the fast calibrating state [SYNC]&10\\ &Orange - Slow blink\\\midrule Not transmitting &Green - Slow flash & The radio carrier is not being transmitted &11\\\midrule Operational &Green - Steady & Default condition if none of the above apply&12 (Low)\\\bottomrule \end {tabular} \end {table} \clearpage \section{GPS assistance data descriptions} Description of assistance data converted and sent inside the RRLP protocol. \begin {table}[ht!] \caption{Almanac message. Table courtesy of \citep{harper2010server-side}.} \label{tbl:almanacMessage}\centering %\rowcolor{2}{light-gray}{} \scriptsize\fontfamily{iwona}\selectfont \begin{tabular}{lll} \toprule %$D$&&$P_u$&$\sigma_N$\\ \textbf{Field (IE)}& \textbf{Description}\\\toprule SatelliteID&This is the satellite ID that is in the range of 0 to 63. PRN=SatelliteID + 1\\\midrule SV Health&Satellite health (e.q. 000 means the satellite is fully operational)\\\midrule $e$&``Eccentricity shows the amount of the orbit deviation from circular (orbit). It is the distance\\ &between the foci divided by the length of the semi-major axis'' \citep{ubxGPSDict}\\\midrule TOA&Time of applicability, reference time for orbit and clock parameters (seconds). ``The number of\\ &seconds in the orbit when the almanac data were generated'' \citep{ubxGPSDict}\\\midrule OI&Orbital inclination (radians). The angle to which the SV orbit meets the equator \citep{ubxGPSDict}\\\midrule RORA&Rate or right ascension (radians/second). ``Rate of change of the angle of right ascension as\\ &defined in the Right Ascension mnemonic'' \citep{ubxGPSDict}\\\midrule $A^{1/2}$& Square root of semi-major axis (meters$^{1/2}$). `` This is defined as the measurement from the center\\ &of the orbit to either the point of apogee or the point of perigee'' \citep{ubxGPSDict}\\\midrule $\Omega_0$& Right Ascension at Week (radians). Longitude of ascending node of orbit plane at weekly epoch\\\midrule $\omega$&Argument of perigee (semicircles). ``An angular measurement along the orbital path measured from\\ &the ascending node to the point of perigee, measured in the direction of the SV's motion'' \citep{ubxGPSDict}\\\midrule $M_0$&Mean anomaly (radians)\\\midrule $a_{f0}$&Satellite clock bias (seconds). Satellite clock error at reference time\\\midrule $a_{f1}$&Satellite clock drift (seconds per second). Satellite clock error rate\\\midrule Week&Week number since the last reset (i.e. since year 1980 modulo 1024 weeks) \\\bottomrule \end {tabular} \end {table} \begin{table}[hc] \scriptsize\fontfamily{iwona}\selectfont \begin{minipage}[b]{.49\textwidth} \centering \begin{tabular}{ll} \toprule %$D$&&$P_u$&$\sigma_N$\\ \textbf{Field (IE)} & \textbf{Description}\\\toprule $\alpha_{0}$&Coefficient 0 of vertical delay\\\midrule $\alpha_{1}$&Coefficient 1 of vertical delay\\\midrule $\alpha_{2}$&Coefficient 2 of vertical delay\\\midrule $\alpha_{3}$&Coefficient 3 of vertical delay\\\midrule $\beta_{0}$&Coefficient 0 of period of the model\\\midrule $\beta_{1}$&Coefficient 1 of period of the model\\\midrule $\beta_{2}$&Coefficient 2 of period of the model\\\midrule $\beta_{3}$&Coefficient 3 of period of the model \\\bottomrule \end {tabular} \caption{GPS Ionosphere Model.} \label{tbl:ionoModel} \end{minipage} \begin{minipage}[b]{.43\textwidth} \centering \begin{tabular}{ll} \toprule %$D$&&$P_u$&$\sigma_N$\\ \textbf{Field (IE)} & \textbf{Description}\\\toprule $A_{1}$&Drift coefficient of GPS time scale relative\\ &to UTC time scale\\\midrule $A_{0}$&Bias coefficient of GPS time scale relative\\ &to UTC time scale\\\midrule $t_{ot}$&Time data reference time of week\\\midrule $\Delta t_{LS}$&Current or past leap second count\\\midrule $WN_{0}$&Time data reference week number\\\midrule $WN_{LSF}$&Leap second reference week number\\\midrule $DN$&Leap second reference day number\\\midrule $\Delta t_{LSF}$&Current of future leap second count \\\bottomrule \end {tabular} \caption{GPS UTC Model.} \label{tbl:utcModel} \end{minipage} \end{table} \newpage \begin {table}[ht!] \caption{Navigation message (ephemeris). Table courtesy of \citep{harper2010server-side}.} \label{tbl:navMessage}\centering %\rowcolor{2}{light-gray}{} \scriptsize\fontfamily{iwona}\selectfont \begin{tabular}{llll} \toprule %$D$&&$P_u$&$\sigma_N$\\ \textbf{Field (IE)} & \textbf{Description}\\\toprule Satellite ID&This is the satellite ID that is in the range of 0 to 63. PRN=SatelliteID + 1\\\midrule Satellite status&This is an indicator of whether this is a new or existing satellite and whether\\ &the navigation model is new or the same.\\\midrule C/A or P on L2&Code(s) on L2 channel\\\midrule URA Index&User range accuracy\\\midrule SV Health&Satellite health\\\midrule IODC&Issue of data, clock\\\midrule L2 P Data flag& \\\midrule SF 1 Reserved& \\\midrule $T_{GD}$&Estimated group delay differential\\\midrule $t_{oc}$&Apparent clock correction\\\midrule $a_{f2}$&Apparent clock correction\\\midrule $a_{f1}$&Apparent clock correction\\\midrule $a_{f0}$&Apparent clock correction\\\midrule $C_{rs}$&Ampltitude of the sine harmonic correction term to the orbit radius (meters)\\\midrule $\Delta n$&Mean motion difference from computed value (semicircles/second)\\\midrule $M_{0}$&Mean anomaly at reference time (semicircles)\\\midrule $C_{uc}$&Ampltitude of the cosine harmonic correction term to the\\ &argument of latitude (radians)\\\midrule $e$&Eccentricity\\\midrule $C_{us}$&Amplitude of the sine harmonic correction term to the argument of latitude\\ &(radians)\\\midrule $A^{1/2}$&Square root of semi-major axis (meters)\\\midrule $t_{oe}$&Reference time ephemeris\\\midrule Fit Interval Flag&\\\midrule AODO&Age of data offset\\\midrule $C_{ic}$&Amplitude of the cosine harmonic correction term to the angle of inclination\\ &(radians)\\\midrule $\Omega_0$&Longitude of ascending node of orbit plane at weekly epoch (semicircles)\\\midrule $C_{is}$&Amplitude of the cosine harmonic correction term to the angle of inclination\\ &(radians)\\\midrule $i_{0}$&Inclination angle at reference time (semicircles)\\\midrule $C_{rc}$&Amplitude of the cosine harmonic correction term to the orbit radius (meters)\\\midrule $\omega$&Argument of perigee (semicircles)\\\midrule OMEGAdot&Rate of right ascension (semicircles/second)\\\midrule Idot&Rate of inclination angle (semicircles/second) \\\bottomrule \end {tabular} \end {table} \newpage \clearpage \section{GPS distance and position estimation} \label{sec:distanceAndPosition} In this section the focus is set on distance and position estimation inside of the GPS receiver. GPS system, as discussed earlier, takes advantage of the TOA ranging concept to determine user position. Time is measured how long it takes for a signal to arrive from a known location. \begin{figure}[ht!] \centering \includegraphics[scale=0.50]{img/Localization.pdf} \caption{Basic distance estimation principle for one satellite. Image courtesy of \citep{understandGPS}.} \label{img:SatLocalization} \end{figure} In figure \ref{img:SatLocalization}, an example concept can be seen, where $\vec{u}=(x_u,y_u,z_u)$ represents the unknown GPS user position vector with respect to Earth-Centered, Earth-Fixed\footnote{ECEF is a Cartesian coordinate system where the point $(0,0,0)$ is defined as the center of mass of the Earth \citep{earthCoordinates}.} (ECEF) coordinate system, $\vec{r}$ is the distance vector from the satellite to the user and $\vec{s}=(x_s,y_s,z_s)$ represents the GPS satellite position with respect to ECEF at a timepoint. Vector $\vec{s}$ is computed from ephemeris data broadcasted by the satellite. The distance vector $\vec{r}$, distance between the satellite and user, can be computed using equation \eqref{eq:r} and its magnitude is given in equation \eqref{eq:rMag}. \begin{equation} \label{eq:r} \vec{r}=\vec{s}-\vec{u} \end{equation} \begin{equation} \label{eq:rMag} r=\Vert s-u\Vert \end{equation} The geometric distance of $r$ is computed by measuring the signal propagation time, this is illustrated in figure \ref{img:TimingLoc} and it was discussed in section \ref{sec:CAdemod}. The PRN code generated on the GPS satellite at time $t_1$ arrives at the time $t_2$, the difference between these two time stamps, $\Delta t$, represents the propagation time. By multiplying the propagation time, $\Delta t$, with the speed of light, $c$, the geometric distance $r$ is computed, as given in equation \eqref{eq:rDist}. \begin{figure}[ht!] \centering \includegraphics[scale=0.50]{img/TimingLoc.pdf} \caption{Estimating the distance by phase shift $\Delta t =t_2 - t_1 =\tau$. Image courtesy of \citep{understandGPS}.} \label{img:TimingLoc} \end{figure} \begin{equation} \label{eq:rDist} r=c\Delta t \end{equation} Since the clocks are not synchronized, as described in sections \ref{sec:SigDemod} and \ref{sec:2dSearch}, clock error offsets have to be added to the geometric distance $r$. This new distance is called \textit{pseudorange}, $\rho$, because the range is determined using the difference of two nonsynchronized clocks (one on the GPS satellite and the other one on the receiver) that generate PRN codes\footnote{pseudo - Not genuine; sham; not perfect.}. Pseudorange is calculated as given in equation \eqref{eq:rho}, where $t_{u}$ is the advance of the receiver clock with respect to the system time\footnote{System time is the exact time on Earth and it is the most precise time known!} and $\delta t$ is the offset of the satellite clock from the system time \citep{understandGPS}. \begin{equation} \label{eq:rho} \rho=r + c(t_{u}-\delta t) \end{equation} Equation \eqref{eq:rMag} can be rewritten as \eqref{eq:rhoR} with respect to equation \eqref{eq:rho}. \begin{equation} \label{eq:rhoR} \rho - c(t_{u}-\delta t) = \Vert s-u\Vert \end{equation} Offset of the satellite clock from the system time, $\delta t$, is updated from Earth, as discussed in \ref{sec:SigDemod} and for that reason it can be removed for sake of simplicity, i.e. it is not an unknown term anymore, then the eqaution \eqref{eq:rhoR} can be rewritten as \eqref{eq:rhoNew}. \begin{equation} \label{eq:rhoNew} \rho - ct_{u} = \Vert s-u\Vert \end{equation} In order to estimate the user (GPS receiver) position, advance of the receiver clock with respect to the system time, $t_u$, has to be found, in other words equation \eqref{eq:rhoSats} has to be solved, where $i$ is the index of visible satellites at the moment of signal reception \citep{understandGPS}. \begin{equation} \label{eq:rhoSats} \rho_i= \Vert s_i-u\Vert + ct_u \end{equation} The estimated position of the user, $\vec{u}=(x_u,y_u,z_u)$, is a three dimensional vector and as stated above the clock offset, $t_u$, is unknown as well. This four dimensional space requires to have at least four pseudorange equations \eqref{eq:rhoSats} to find all the four unknown terms. As a result of this fact, at least four satellites have to be visible at the same time to estimate the position of the target user. Equation given in \eqref{eq:rhoSats} takes the form in \eqref{eq:rhoSatsNew} because the coordinate system is Cartesian and $\rho_i$ is nothing else but Euclidean distance where $i=1,2,...,n$ such that $n\geq4$ and $\vec{s_i}=(x_i,y_i,z_i)$ is the satellite position estimated from the ephemeris data. \begin{equation} \label{eq:rhoSatsNew} \rho_i= \sqrt{(x_i-x_u)^2+(y_i-y_u)^2+(z_i-z_u)^2} + ct_u \end{equation} Undoubtedly, the given equation in \eqref{eq:rhoSatsNew} is a nonlinear equation\footnote{Nonlinear equations, also known as polynomial equations, are equations that can not satisfy both of the linearity properties: additivity $f(x+y)=f(x)+f(y)$ and homogeneity $f(\alpha x) = \alpha f(x)$, $\alpha \in \mathbb{R}$ \citep{nonlinear}.}. It is not straightforward to find explicit solutions of nonlinear equations, it is more difficult than compared to linear equations. There are different techniques to solve sets of nonlinear equations \citep[Chapter 7]{understandGPS} but in this work the linearization method\footnote{Linear approximation is a technique where a function is approximated using a linear function.} shall be presented to find the unknown terms $(x_u,y_u,z_u,t_u)$, i.e. out of an approximate position and clock offset the true user position and the true clock offset shall be calculated. \begin{equation} \label{eq:rhoSatsNewFun} \rho_i= \sqrt{(x_i-x_u)^2+(y_i-y_u)^2+(z_i-z_u)^2} + ct_u = f(x_u,y_u,z_u,t_u) \end{equation} Let the equation \eqref{eq:rhoSatsNew} for pseudoranges, be rewritten as a function $f$ of four unknown terms $x_u$, $y_u$, $z_u$ and $t_u$, as given in \eqref{eq:rhoSatsNewFun} \citep[Chapter 2]{understandGPS}. Suppose that the approximation of the position and the clock offset are known, denoted as $\hat{x_u}$, $\hat{y_u}$, $\hat{z_u}$ and $\hat{t_u}$, then equation \eqref{eq:rhoSatsNewFun} can be rewritten as an approximate pseudorange \eqref{eq:rhoSatsNewFunApprox}. \begin{equation} \label{eq:rhoSatsNewFunApprox} \hat{\rho_i}= \sqrt{(x_i-\hat{x_u})^2+(y_i-\hat{y_u})^2+(z_i-\hat{z_u})^2} + c\hat{t_u} = f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u}) \end{equation} In other words, the unknown true position terms $x_u$, $y_u$, $z_u$ and the clock offset term $t_u$, of the GPS receiver, shall be expressed by the approximate values and an incremental component as shown in equation \eqref{eq:userCoordinates} \citep{understandGPS}. \begin{equation} \label{eq:userCoordinates} \begin{array}{l} x_u = \hat{x_u}+\Delta x_u \\ y_u = \hat{y_u}+\Delta y_u \\ z_u = \hat{z_u}+\Delta z_u \\ t_u = \hat{t_u}+\Delta t_u \end{array} \end{equation} By inserting the terms from \eqref{eq:userCoordinates} into equation \eqref{eq:rhoSatsNewFun}, a new equation is derived as in \eqref{eq:rhoSatsNewFunwithApprox}. \begin{equation} \label{eq:rhoSatsNewFunwithApprox} f(x_u,y_u,z_u,t_u) = f(\hat{x_u}+\Delta x_u, \hat{y_u}+\Delta y_u, \hat{z_u}+\Delta z_,\hat{t_u}+\Delta t_u) \end{equation} In the next step the pseudorange function shall be approximated using Taylor series\footnote{Taylor series ``is a representation of a function as an infinite sum of terms that are calculated from the values of the function's derivatives at a single point'' \citep[Chapter 11]{taylor}.} (linearization of the nonlinear equation). Taylor series for a function $f(x)$ is given in equation \eqref{eq:taylor}, where as $a$ approaches $x$ the estimation error shall be smaller and smaller, i.e. $f(x) = f(a)$ when $x=a$. The approximation error depends on Taylor polynomial degree (the amount of terms or taken derivatives of the function) and how far away the point $a$ is from $x$ \citep[Chapter 11.9]{taylor}. The basic idea of the principle can be seen in figure \ref{img:taylorSeries}. \begin{equation} \label{eq:taylor} f(x) = \sum_{n=0}^{\infty}\frac{f^{(n)}(a)}{n!}(x-a)^n = f(a) + \frac{f'(a)}{1!}(x-a)+\frac{f''(a)}{2!}(x-a)^2+... \end{equation} \begin{figure}[ht!] \centering \includegraphics[scale=0.50]{img/TaylorSeries.pdf} \caption{Taylor series approximation for a point $a=0.5$ where $n$ is the Taylor polynomial degree.} \label{img:taylorSeries} \end{figure} Due to the four unknown terms, Taylor series for multivariables have to be used. The general formula is given in equation \eqref{eq:Multitaylor}, where vector $\mathbf{x}\in\mathbb{R}^n$ denotes $n$ variables, $\nabla$ (nabla) is the Del\footnote{Del, $\nabla$, is the vector differential operator.} operator given in \eqref{eq:Del} and $\mathbf{a}$ is the linearization point of interest \citep{multiTaylor}. \begin{equation} \label{eq:Multitaylor} f(\mathbf{x}) \approx f(\mathbf{a}) + \nabla f |_{\mathbf{x=a}} \cdot (x-a) \end{equation} \begin{equation} \label{eq:Del} \nabla^{T} = \left[\frac{\partial}{\partial x_{1}}...\frac{\partial}{\partial x_{n}}\right] \end{equation} One can note that in equation \eqref{eq:Multitaylor} the Taylor series polynomial is of the first degree. This is because of one reason, it linearizes the approximation of the function $f(\mathbf{x})$ at point $\mathbf{a}$ and as a consequence it removes the nonlinearities \citep{understandGPS} \citep[Chapter 11.10]{taylor}, as seen in figure \ref{img:taylorSeries}, for $n=1$ the resulting function is linear. In the previously described step, one would calculate a hyperplane tangent to a point $a$ in a $n$-Dimensional space. By inserting equation \eqref{eq:rhoSatsNewFunwithApprox} in equation \eqref{eq:Multitaylor}, it yields equation \eqref{eq:MultitaylorFour} where $\mathbf{x}=(x_u,y_u,z_u,t_u)$ and $\mathbf{a}=(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u})$. \begin{equation} \label{eq:MultitaylorFour} \begin{array}{l} f(\hat{x_u}+\Delta x_u, \hat{y_u}+\Delta y_u, \hat{z_u}+\Delta z_,\hat{t_u}+\Delta t_u) \approx f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u}) \\[0.5em] + \dfrac{\partial f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u})}{\partial \hat{x_u}}\Delta x_u +\dfrac{\partial f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u})}{\partial \hat{y_u}}\Delta y_u \\ +\dfrac{\partial f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u})}{\partial \hat{z_u}}\Delta z_u +\dfrac{\partial f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u})}{\partial \hat{t_u}}\Delta t_u \end{array} \end{equation} The terms from equation \eqref{eq:MultitaylorFour} are solved individually in equations \eqref{eq:MultitaylorDeriv} where $\sqrt{(x_i-\hat{x_u})^2+(y_i-\hat{y_u})^2+(z_i-\hat{z_u})^2}$ has been subsituted with $\hat{r_i}$. \begin{equation} \label{eq:MultitaylorDeriv} \begin{array}{l} \dfrac{\partial f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u})}{\partial \hat{x_u}} = \dfrac{1}{2}\dfrac{-2(x_{i}-\hat{x_{u}})}{\sqrt{(x_i-\hat{x_u})^2+(y_i-\hat{y_u})^2+(z_i-\hat{z_u})^2}} =-\dfrac{x_i-\hat{x_u}}{\hat{r_i}}\\[0.9em] \dfrac{\partial f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u})}{\partial \hat{y_u}} = \dfrac{1}{2}\dfrac{-2(y_{i}-\hat{y_{u}})}{\sqrt{(x_i-\hat{x_u})^2+(y_i-\hat{y_u})^2+(z_i-\hat{z_u})^2}} =-\dfrac{y_i-\hat{y_u}}{\hat{r_i}}\\[0.9em] \dfrac{\partial f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u})}{\partial \hat{z_u}} = \dfrac{1}{2}\dfrac{-2(z_{i}-\hat{z_{u}})}{\sqrt{(x_i-\hat{x_u})^2+(y_i-\hat{y_u})^2+(z_i-\hat{z_u})^2}} =-\dfrac{z_i-\hat{z_u}}{\hat{r_i}}\\[0.9em] \dfrac{\partial f(\hat{x_u},\hat{y_u},\hat{z_u},\hat{t_u})}{\partial \hat{t_u}} = c \end{array} \end{equation} Then by substituting the equation terms from \eqref{eq:MultitaylorDeriv}, \eqref{eq:rhoSatsNewFun} and \eqref{eq:rhoSatsNewFunApprox} into \eqref{eq:MultitaylorFour}, the resulting equation is given in \eqref{eq:MultitaylorDerivAfter}. \begin{equation} \label{eq:MultitaylorDerivAfter} \begin{array}{l} \rho_i = \hat{\rho_i} -\dfrac{x_i-\hat{x_u}}{\hat{r_i}}\Delta x_u -\dfrac{y_i-\hat{y_u}}{\hat{r_i}}\Delta y_u -\dfrac{z_i-\hat{z_u}}{\hat{r_i}}\Delta z_u + c\Delta t_u \end{array} \end{equation} At this step, by solving equation \eqref{eq:MultitaylorFour}, the linearization of the nonlinear equations is completed. \begin{equation} \label{eq:MultitaylorDerivAfterRearange} \begin{array}{l} \hat{\rho_i} - \rho_i = \dfrac{x_i-\hat{x_u}}{\hat{r_i}}\Delta x_u +\dfrac{y_i-\hat{y_u}}{\hat{r_i}}\Delta y_u +\dfrac{z_i-\hat{z_u}}{\hat{r_i}}\Delta z_u - c\Delta t_u \end{array} \end{equation} \begin{equation} \label{eq:SubsTerms1} \Delta\rho = \hat{\rho_i} - \rho_i \\[0.7em] \end{equation} \begin{equation} \label{eq:SubsTerms2} \alpha_{xi} = \dfrac{x_i - \hat{x_u}}{\hat{r_i}} \hspace{1.5em} \alpha_{yi} = \dfrac{y_i - \hat{y_u}}{\hat{r_i}} \hspace{1.5em} \alpha_{zi} = \dfrac{z_i - \hat{z_u}}{\hat{r_i}} \end{equation} By rearanging the equation \eqref{eq:MultitaylorDerivAfter} one derives equation \eqref{eq:MultitaylorDerivAfterRearange}. And then by substituting the terms in \eqref{eq:SubsTerms1} and \eqref{eq:SubsTerms2} into \eqref{eq:MultitaylorDerivAfterRearange}, the equation resembles the one given in \eqref{eq:userPosition}. \begin{equation} \label{eq:userPosition} \Delta\rho_i = \alpha_{xi}\Delta x_u + \alpha_{yi}\Delta y_u + \alpha_{zi}\Delta z_u - c\Delta t_u \end{equation} There are four unknowns, $\Delta x_u$, $\Delta y_u$, $\Delta z_u$ and $\Delta t_u$, in equation \eqref{eq:userPosition}. By solving this set of linear equations, which shall result in finding $\Delta x_u$, $\Delta y_u$, $\Delta z_u$ and $\Delta t_u$, the GPS receiver position $(x_u, y_u, z_u)$ and clock offset $t_u$ is computed by replacing the same into equations in \eqref{eq:userCoordinates}. Equation \eqref{eq:userPosition} can be rewritten for four satellites in the matrix form as in \eqref{eq:userPositionMatrix}. \begin{equation} \label{eq:userPositionMatrix} \Delta\boldsymbol{\rho} = \boldsymbol{\alpha} \Delta \boldsymbol{x} \end{equation} \begin{equation} \Delta\boldsymbol{\rho}= \begin{bmatrix} \Delta \rho_1 \\ \Delta \rho_2 \\ \Delta \rho_3 \\ \Delta \rho_4 \end{bmatrix} \hspace{1.5em} \boldsymbol{\alpha}= \begin{bmatrix} \alpha_{x1} & \alpha_{y1} & \alpha_{z1} & 1 \\ \alpha_{x2} & \alpha_{y2} & \alpha_{z2} & 1 \\ \alpha_{x3} & \alpha_{y3} & \alpha_{z3} & 1 \\ \alpha_{x4} & \alpha_{y4} & \alpha_{z4} & 1 \end{bmatrix} \hspace{1.5em} \Delta \boldsymbol{x}= \begin{bmatrix} \Delta x_u \\ \Delta y_u \\ \Delta z_u \\ -\Delta ct_u \end{bmatrix} \end{equation} Finally, by multiplying both left sides\footnote{Matrix multiplication is not communitative, $\mathbf{AB\neq BA}$.} of the equation \eqref{eq:userPositionMatrix} with the inverse term of $\boldsymbol{\alpha}$, it yields the result of the unknown terms, as given in equation \eqref{eq:userPositionMatrixFinal}. \begin{equation} \label{eq:userPositionMatrixInverseMult} \boldsymbol{\alpha}^{-1}\Delta\boldsymbol{\rho} = \boldsymbol{\alpha}^{-1}\boldsymbol{\alpha} \Delta \boldsymbol{x} \end{equation} \begin{equation} \label{eq:userPositionMatrixFinal} \Delta \boldsymbol{x} = \boldsymbol{\alpha}^{-1} \Delta\boldsymbol{\rho} \end{equation} Linearization is repeated in a loop, where in the next round the approximate positions are set to the just derived position values, that is, $\hat{x_u}=x_u$, $\hat{y_u}=y_u$, $\hat{z_u}=z_u$ and $\hat{t_u}=t_u$. This process is repeated until the approximated positions converge to their final values. It is not necessarily required that the initial positions are very accurate and the results are usually obtained by 4-5 itterations \citep{pseudorangeError}. Risks exist that the solutions shall still be corrupted but there are different error avoiding mechanisms to solve these problems, like minimizing the error contribution using more than four satellite measurements \citep{pseudorangeError} \citep[Chapter 7]{understandGPS}. %\clearpage %\section{GPS Constants and equations} %\label{sec:gpsConsAndEq} %\begin{alignat}{4} % & A & = & \; (\sqrt{A})^2 \nonumber \\ % & n_{0} & = &\; \sqrt{\frac{\mu}{A^3}} \nonumber \\ % & t_{k} & = &\; t-t_{oe} \nonumber \\ % & n & = &\; n_{0} + \Delta n \nonumber \\ % & M_{k} & = &\; M_{0} + nt_{k} \nonumber \\ % & M_{k} & = &\; E_{k} - e\sin E_{k} \nonumber \\ % & v_{k} & = & \tan ^{-1} \left( \frac{\sin v_{k}}{\cos v_{k}} \right) = \tan ^{-1} \left( \frac{\frac{\sqrt{1-e^2} \sin E_{k}}{1-e \cos E_{k}}}{\frac{\cos E_{k}-e}{1-e\cos E_{k}}} \right) \nonumber \\ % & v_{k} & = & \tan ^{-1} \left( \frac{\sin v_{k}}{\cos v_{k}} \right) = \tan ^{-1} \left( \frac{\sqrt{1-e^2} \sin E_{k}/(1-e \cos E_{k})}{(\cos E_{k}-e)/(1-e\cos E_{k})} \right) = \tan ^{-1} \left( \frac{\sqrt{1-e^2} \sin E_{k}}{\cos E_{k} - e} \right) \nonumber \\ % & E_{k} & = & \cos ^{-1} \left( \frac{e+\cos v_{k}}{1+e \cos v_{k}} \right) \nonumber \\ % & \Phi_{k} & = &\; v_{k} + \omega \nonumber \\ % & \delta u_{k} & = &\; c_{us} \sin{2\Phi_{k}} + C_{us} \cos{2\Phi_{k}} \\ % & \delta r_{k} & = &\; c_{rc} \cos{2\Phi_{k}} + C_{rs} \sin{2\Phi_{k}} \nonumber \\ % & \delta i_{k} & = &\; c_{ic} \cos{2\Phi_{k}} + C_{is} \sin{2\Phi_{k}} \nonumber \\ % & u_{k} & = &\; \Phi_{k} + \delta u_{k} \nonumber \\ % & r_{k} & = &\; A(1-e\cos{E_{k}})+\delta r_{k} \nonumber \\ % & i_{k} & = &\; i_{0} + \delta i_{k} + (IDOT)t_{k} \nonumber \\ % & x_{k}^{'} & = &\; r_{k} \cos{u_{k}} \nonumber \\ % & y_{k}^{'} & = &\; r_{k} \sin{u_{k}} \nonumber \\ % & \Omega_{k} & = &\; \Omega_{0} + (\Omega - \Omega_{e})t_{k} - \Omega_{e}t_{oe} \nonumber \\ % & x & = &\; x_{k}^{'} \cos{\Omega_{k}}-y_{k}^{'}\cos{i_{k}}\sin{\Omega_{k}} \nonumber \\ % & y & = &\; x_{k}^{'} \sin{\Omega_{k}}-y_{k}^{'}\cos{i_{k}}\cos{\Omega_{k}} \nonumber \\ % & z & = &\; y_{k}^{'} \sin{i_{k}} \nonumber %\end{alignat} %\begin{equation} %\label{eq:paramconst1} % \begin{split} % \mu_{e} = 3.986004418\cdot 10^{14} \frac{m^3}{s^2} % \end{split} %\quad\Longleftarrow\quad % \begin{split} % \mbox{Geocentric gravitational constant} % \end{split} %\end{equation} %\begin{equation} %\label{eq:paramconst2} % \begin{split} % c= 2.99792458\cdot 10^{8} \frac{m}{s} % \end{split} %\quad\Longleftarrow\quad % \begin{split} % \mbox{speed of light} % \end{split} %\end{equation}