path: root/doc/thesis/chapters/implementation.tex

% vim: set ts=2 sw=2 noet spell:

\chapter{Implementation} \label{chp:implementation}

\section{Overview}

To illustrate the effects of multipath fading on a wireless transmission link, two machines with an USRP B210 each (hardware) were set up to run a software defined radio (SDR) using the GNU Radio (GR) signal processing toolkit. The latter can also be used as standalone (without hardware) to simulate a channel model. Therefore, simulations of different scenarios affected by multipath fading were elaborated and constructed first using a static FIR model\footnote{As discussed in sections \ref{sec:discrete-time-model} and \ref{sec:fractional-delay}}, and then with a dynamic statistical model\footnote{See section \ref{sec:statistical-model}}. To present the results a graphical interface was made using the Dear IMGUI framework.

The rest of the chapter is structured as follows. First the tools used in this project are introduced. Then implementations of the transmitter and receiver chains is explained. Subsequently simulations and measurements of fading channels with their respective empirically computed bit error rate (BER) are presented. Finally, the current state of the implementation and issues are discussed.

\section{Software Stack}

\subsection{GNU Radio}

For both the signal processing and the simulations the GNU Radio (GR) toolkit was chosen, as it already had integrations with the USRP hardware drivers. GR is an open-source free software framework that can be used to build signal processing chains and SDRs. GR is composed of two parts: a C\texttt{++} library with Python bindings to write signal processing code, and GNU Radio Companion (GRC). GRC is a graphical user interface made to more easily construct signal processing chains. Signal processing algorithms are graphically shown as ``blocks'' that can be chained together with arrows, essentially drawing a diagram called ``flow graph''. An example of a flow graph can be seen in \figref{fig:sync-lock-flowgraph}.

Internally GR works by keeping multiple memory buffers of samples, that are passed as pointers to the signal processing algorithms' ``work functions''. When the signal processing is complete, the output buffer of one block is given to the next block as input according to how they were connected in the flow graph. The structure of a block is shown in the Python listing \ref{lst:gr-block-py}. To improve performance GR creates a thread for each work function to parallelize the workload of the concurrently running signal processing blocks. For more details see the GNU Radio Wiki and User Manual in \cite{GRWiki}.

{\newcommand{\placeholder}[1]{\textit{\(\langle\)\,\textrm{#1}\,\(\rangle\)}}
\begin{lstlisting}[
	language = python, escapechar = {`}, float,
	caption = {A minimal GNU Radio block in Python.},
	captionpos = b, label = {lst:gr-block-py}
]
class myblock(gr.sync_block):
	# there are also other types of blocks such as interpolators 
	# (more outputs that inputs), decimators (more inputs than 
	# outputs) sync blocks have a 1-to-1 input to output ratio
	def __init__(self, `\placeholder{parameters}`):
		gr.sync_block.__init__(self, name="My Block",
			# this block as one input port and one output port
			# with samples that are 64 bit complex numbers
			in_sig=[np.complex64], out_sig=[np.complex64]
		)

	def work(self, inputs, outputs):
		# signal processing goes here
		# inputs and outputs are k x n arrays, where each
		# of the k rows is a port that contains n samples
		return `\placeholder{N of inputs that were processed}`
\end{lstlisting}}

\subsection{Dear PyGUI}\label{sec:GUI}

To construct a graphical interface for a demonstration platform the Dear IMGUI (immediate mode graphical user interface) library was chosen, mainly for its ease of use, wide range of technical capabilites and high refresh rate. Dear PyGUI (DPG) are the Python bindings for the Dear IMGUI library.

The DPG front-end communicates with the GR flow graphs using the IP/UDP protocol. This decision to separate the project into two parts that communicate over the IP network was made because it is not easy to extend the graphical interface of GRC without interfering with the sophisticated multi-threaded architecture of GR. Furthermore, this allows to have multiple correctly configured flow graphs on disk and to choose which one to run and display on the graphical interface, instead of having a single flow graph whose parameters need to be changed each time. As a side effect, theoretically this setup allows to have one computer running the graphical interface, and another remote machine running just the flow graph. Though the latency caused by the UDP/IP could be substantial.

\section{Hardware}

\begin{table}[h]
	\centering
	\begin{tabular}{ll}
		\toprule
		Dimensions               & \(9.7 \times 15.5 \times 1.5\)\,cm \\
		Ports                    & 2\,TX, 2\,RX, Half or Full Duplex   \\
		RF frequencies           & \SI{70}{\mega\hertz} to \SI{6}{\giga\hertz} \\
		Bandwidth                & \SI{200}{\kilo\hertz} --- \SI{56}{\mega\hertz} \\
		External reference input & \SI{10}{\mega\hertz}                         \\
		\bottomrule
	\end{tabular}
	\caption{USRP B210 specifications \cite{EttusUSRPB210}. \label{tab:usrp-specs}}
\end{table}

As already mentioned, for the receiver and transmitter in the SDR setup two USRP B210 devices from Ettus Research were used. Some technical specifications of the device are shown in \tabref{tab:usrp-specs}. GR provides off the shelf blocks that interface with the official API provided from Ettus Research.

\section{Transmitter chain}

\subsection{Data frame} \label{sec:data-frame}

\begin{figure}
\centering
\input{figures/tikz/packet-frame}
	\caption{
		Structure of framed data packets used in the implementation.
		\label{fig:dataframe}
	}
\end{figure}

To compute the empirical bit error rate (BER) of the setup, the data has to be framed by the transmitter and the bitstream needs to be synchronized on the receiver side. The structure of a data packet used in the implementation is shown in \figref{fig:dataframe}. A frame begins with a user specified \(k\)-byte preamble, that in the current implementation serves as synchronization pattern. Another use case for the preamble sequence could be to introduce channel estimation pilot symbols (this will be discussed more in section \ref{sec:psam}). Following the preamble are 4 bytes encoded using a (31, 26) Hamming code (plus 1 padding bit), that contain metadata about the packet, namely payload ID and payload length. Because the payload length in bytes is encoded in 21 bits, the theoretical maximum payload size is 2 MiB, which together with 32 possible unique IDs gives a maximum data transfer with unique frame headers of 64 MiB. These constraints are a result of decisions made to keep the implementation simple. % TODO: explain why its simpler this way


\subsection{Modulation}

GR provides a constellation modulator block, that already implements several standard constellations (QPSK and 16-ary QAM being of interest for us). The block also already integrates a root raised cosine filter, whose phase bandwidth (roll-off factor) can be given as parameter; in all flow graphs the roll off factor is \(\alpha = 0.35\).

% TODO: Warum alpha 0.35
%
% Because we had no restrictions on bandwidth (except for the physical ones, which are unreachable). Though, we are in the 2.4 GHz spectrum which is pretty crowded. For a sanity check: we are using a very short symbol time of around T = 1 / 1 MHz = 1 us, then the bandwidth with \alpha = 0.35 is B = (1 + \alpha) / (2 * T) = 675 kHz. Pretty small, when compared for ex to WiFi 802.11g channels that are like 20 MHz wide.

\section{Receiver chain}

\subsection{Envelope detector}

What is here referred to as envelope detector has the purpose of synchronizing the symbols and equalizing the input signal amplitude. This is accomplished in GRC using two blocks: a polyphase clock sync and a constant modulus adaptive (CMA) filter equalizer. The input signal for the envelope detector has 4 samples per symbol, while the output has only one sample per symbol. The choice of the CMA equalizer later turned out to be a mistake in the QAM flow graphs, as it only works for envelopes with a constant amplitude. In the latest version least mean square decision-directed (LMS-DD) equalizers have been used.

% \paragraph{Polyphase Clock Sync}
% %TODO : nochmals anschauen ob dieese erklärung verständlich ist und richtig interpretiert wurde.
% With the the polyphase clock sync the symbols can be synchronized by preforming a time synchronization with the help of multiple filterbanks. For that the derivative of the filtered signal should be minimized which turns to a better SNR. 
% This works with the help of two filterbanks, one of them contains the filters of the signal adapted to the pulse shaping with several phases. The other contains its derivative. So in the time domain it has a sinc shape, for the output Signal the sinc peak should be on a sample, with the fact that sinc(0) = 1 and sinc(0)' = 0 an error signal can be generated which tells how far away from the peak it is. This error Signal should be zero this is possible with the help of a loop second order whish constants the number of the filterbank and the rate. This rate is generated because of the clock difference between the transmitter and receiver to synchronized the receiver the filter goes through the phases. For the output one sample per symbol is enough.

% \paragraph{Equalizer}
% \skelpar[2]{CMA equalizer.}

\subsection{Frame synchronization and phase correction} \label{sec:phasecorr}

Once the envelope's clock is synchronized in the processing chain the data stream has one sample per symbol. At this point it is necessary to find where each data frame starts or end in order to correctly decode their payloads. For such purpose a special sequence called \emph{access code} is put in front of each frame. Access codes are sequences of samples that are carefully constructed to have an autocorrelation with a high peak at zero, and that rapidly decreases for increasing shifts. In other words, the autocorrelation of an access code high only when the sequence is perfectly aligned. Thus, by cross correlating an envelope signal \(r(t)\), that periodically contains an access code \(a(t)\) with the access code itself, and looking for peaks in the result, it is possible to determine where each frame begin. Furthermore, by analyzing the values of the peaks it is possible to extract information about the phase and frequency offsets.

To understand how correlation peaks allow for fine phase correction, recall that the cross correlation (denoted here by \(\star\)) of two complex valued signals is
\begin{equation}
	R_{ra}
	= (r \star a)(t)
	= \int_\mathbb{R} r(\tau) a^*(\tau - t) \,d\tau
	= r(t) * a^*(-t),
\end{equation}
which is equivalent to a convolution, with the left term being time-reversed complex conjugated \cite{Gallager}. This last property is especially useful because it makes possible to implement cross correlation using FIR filters. Some interesting properties of the cross correlation are that correlation with itself (autocorrelation) at \(t = 0\) is
\begin{equation}
	R_{aa} = (a \star a)(0)
	= \int_\mathbb{R} a(\tau) a^*(\tau - 0) \,d\tau
	= \int_\mathbb{R} |a(\tau)|^2 \,d\tau \in \mathbb{R},
\end{equation}
which is a real number. And more importantly the correlation with an out of phase copy \(a'(t) = a(t) e^{j\varphi}\) at \(t = 0\) is
\begin{equation} \label{eqn:xc-oop-copy}
	% R_{a'a} =
	(a' \star a)(0) 
	= \int_\mathbb{R} a(\tau)e^{j\varphi}  a^*(\tau) \,d\tau
	= R_{aa} e^{j\varphi}.
\end{equation}
The relevant observation to make in \eqref{eqn:xc-oop-copy} is that since \(R_{aa}\) is a real number, the phase of the cross correlation at \(t = 0\) is the phase of \(a'(t)\). This fact can be exploited to implement fine phase correction for the received envelope in relatively few steps as follows:
\begin{enumerate}
	\item Compute the cross correlation \(R_{ra}\) of the envelope \(r(t)\) with the access code \(a(t)\),
	\item Find the maximum value of \(\hat{R}_{ra} = \max R_{ra}(t)\) (correlation peak),
	\item Extract the phase offset \(\varphi = \arg \hat{R}_{ra}\),
	\item Remove the phase offset in the envelope by multiplying it with the complex conjugate of the offset, that is \(\hat{r}(t) = r(t) e^{-j\varphi}\).
\end{enumerate}

\begin{figure}
	\centering
	\includegraphics[width = .95\linewidth]{figures/screenshots/sync_lock}
	\caption{
		Part of the GNU Radio flow graph for the QPSK modulated link (with hardware). The shown blocks are used to equalize and lock the envelope.
		\label{fig:sync-lock-flowgraph}
	}
\end{figure}

\subsubsection{Implementing fine phase and frequency correction} \label{sec:implement-phasecorr}

To implement in GR what was discussed above the two blocks shown in \figref{fig:sync-lock-flowgraph} were used: a correlator estimator block, and a custom block. The former essentially implements the first 3 of the steps discussed at the end of the previous subsection. The correlator estimator block is given a sequence of samples, and when the cross correlation between them and the input stream is higher than a certain threshold (90\% of the amplitude of a perfect autocorrelation), it produces a ``tag'' in the output stream, that contains the phase estimate.

Tags are GR's way of working with metadata that is attached to a sample. Internally tags are just polymorphic data structures containing a number indicating the absolute offset (in samples), and a pair of arbitrary values called ``key'' and ``value''. Tags are passed on from one block to the next like sample streams (unless the block specifies to do otherwise).

Thus, the tagged stream is processed with a custom block, of which a simplified version of its work function shown in listing \ref{lst:phasecorr-work}. The custom block also implements fine frequency correction (shown in listing \ref{lst:phasecorr-blockphase}) by linearly interpolating the phase estimates between each pair of tags (called chunk). This can be rather trivially be formulated for a chunk of \(N\) samples as the
\begin{subequations}
	\begin{align}
		k\text{-th chunk digital frequency} \quad  & \Omega_k = (\varphi_{k+1} - \varphi_k) / N, \text{ and the }\\
		k\text{-th chunk phase estimate} \quad & \Phi(n) = \varphi_k - \omega_k n/N.
	\end{align}
\end{subequations}

\subsubsection{Performance of the implementation}\label{sec:preforming-implementation}

The phase and frequency correction block was implemented with the design goal of being able to correct a maximum frequency offset of \(\hat{\epsilon} = 0.1\%\) under ideal conditions, which is sufficient to take into account small Doppler shifts at walking speed (\(v = \SI{2}{\meter\per\second}\)) with carrier at \(f_c = 2.4\) GHz. The USRP B210 devices have an internal clock frequency accuracy of \(\epsilon = 1\text{ ppm} = 10^{-6}\), which results in a total frequency offset of
\begin{equation}\label{eq:doppler}
	\Delta f = f_c \left( \frac{v}{c_0} + \epsilon \right)
	= \SI{2.4}{\giga\hertz} \left(
		\frac{\SI{2}{\meter\per\second}}{\SI{3e8}{\meter\per\second}} + 10^{-6} 
	\right) = \SI{2416}{\hertz}.
\end{equation}
Because the frequency estimate is linearly interpolated, the phase error may not exceed \(\pi\) (half rotation) during one data frame (chunk). These constraints imply that for frames of \(N'\) symbols of duration \(T\) each, using \(\kappa\) samples per symbol the relation
\begin{equation}\label{Doppler-shift}
	2\pi\Delta f N' T \kappa \leq \pi
	\implies T = 1/f_s \leq \frac{1}{2\Delta f N' \kappa},
	\iff N' \leq \frac{1}{2\Delta f T \kappa},
\end{equation}
must hold. By further setting \(\kappa = 4\) and \(N' = 32\) we obtain a minimum sampling frequency of approximately \(\SI{618.5}{\kilo\hertz}\), or conversely by letting \(f_s = \SI{1}{\mega\hertz}\) we have a maximum frame length of \(N' = 51\) symbols. In other words, roughly every 50 symbols the system must send an access code sequence. This result is rather unfortunate as it means that more processing power than expected is required.

\begin{figure}
	\centering
	\input{figures/tikz/phasecorr-blockprocessing-diagram}
	\caption{
		Graphical representation of the input samples for the work function of the fine phase and frequency correction block (shown in listing \ref{lst:phasecorr-work}). Roughly every \(N\) samples there is a tag containing the information of the phase error (computed using the cross correlation peak). The white `chunks' of samples can be corrected using their respective left and right tag values. The samples in the red chunk need phase information from the previous block processing. The samples in the blue chunk need a phase information from the future, which is not attainable. Thus for the blue chunk the frequency estimate of the previous chunk is used.
		\label{fig:phasecorr-chunks}
	}
\end{figure}


\begin{lstlisting}[
	texcl = true, language = python, escapechar = {`},
	float, captionpos = b, label = {lst:phasecorr-work},
	caption = {
		Simplified work function of fine phase correction block that corrects only samples `in the middle'. The version that is actually used does handle edge cases that have been removed here for readability. See also \figref{fig:phasecorr-chunks} for a graphical representation of the inputs and listing \ref{lst:phasecorr-blockphase} for the definition of the \texttt{block\_phase} function.
	},
]
def work(self, inputs, outputs):
	# alias for inputs of the first port
	inp = inputs[0]
	# read phase tags from stream
	is_phase = lambda tag: pmt.to_python(tag.key) == "phase_est"
	tags = filter(is_phase, self.get_tags_in_window(0, 0, len(inp)))
	# create a list of pairs \(((\varphi_0,\varphi_1), (\varphi_1, \varphi_2), \ldots, (\varphi_{k-1}, \varphi_k)))\)
	pairs = zip(tags, tags[1:])
	# compute phase correction between each pair of tags
	chunks = [self.block_phase(start, end) for (start, end) in pairs]
	# flatten array to get \(\Phi(n)\) and compute the correction
	phases = np.concatenate(chunks)
	correction = np.exp(-1j * phases) 
	# write to the first output port
	left = tags[0].offset - self.nitems_written(0)
	right = tags[-1].offset - self.nitems_written(0)
	outputs[0][left:right] = inp * correction
	# return how many samples were processed
	return len(outputs[0])
\end{lstlisting}

\begin{lstlisting}[
	texcl = true, language = python, escapechar = {`},
	float, captionpos = b, label = {lst:phasecorr-blockphase},
	caption = {
		Block phase function referenced in listing \ref{lst:phasecorr-work}.
	},
]
def block_phase(self, start, end):
	# compute number of samples between tags
	nsamples = end.offset - start.offset
	# unpack pmt values into start and end phase
	sphase = pmt.to_python(start.value)
	ephase = pmt.to_python(end.value)
	# compute frequency offset between start and end
	phasediff = (ephase - sphase)
	freq = phasediff / nsamples
	if phsediff > np.pi:
		phasediff -= np.pi
	elif phasediff < -np.pi:
		phasediff += np.pi
	# compute chunk values
	return sphase * np.ones(nsamples) + freq * np.arange(0, nsamples)
\end{lstlisting}

\subsection{GUI implementation}
\begin{figure}
	\centering
	\includegraphics[frame, width = \linewidth]{figures/screenshots/gui_screenshot}
	\caption{Screenshot of the graphical interface of receiver built using the DearPyGUI library.
	\label{fig:GUI}}
\end{figure}

The GUI is implemented with the Dear PyGUI library as introduced in section \ref{sec:GUI}. In \figref{fig:GUI} the graphical interface of it is shown. 
Four different constellation plots are illustrated: the channel itself, after the polyphase clock synchronization, after the equalizer and the locked constellation.
The GUI shows the BER of the constellation and a time domain plot. The interface also contains a block diagram.


\section{Channel simulations}

In order to study the effects of multipath fading, a series of simulations have been made under different conditions. To simulate a channel affected by multipath fading two blocks from the GR library, and a third custom block were used. The channel model can simulate AWGN, a frequency offset and either a Rayleigh (NLOS) oder Rice (LOS) fading.

\subsection{Fading with discrete time model} \label{sec:discrete-time-model-fir}

 To implement and illustrate a static fading effect corresponding to section \ref{sec:discrete-time-model}, a custom block was created. The work function of the block is shown in listing \ref{lst:fractional-delay-fir}. This block is based on a FIR filter. It can be configured to simulate a with direct path (LOS) or without (NLOS).

 In this block it is possible to simulate any number of paths with different attenuations. A special case is show in the first and second colum of \figref{fig:qpsk-simulations-static}, where the delay is a multiple of the symbol time. Another example with more taps is shown in the third colum.  These simulation values are however not realistic, because the static model is too simple. It does not account for other effects, such as changes in the environment. The block additionally implemented the method described in section \ref{sec:fractional-delay} to allow non-integer delays. Here the sinc pulse was limited using a simple boxcar function, the filter could be improved by using for example a Hann window.


\begin{lstlisting}[
	texcl = true, language = python, escapechar = {`},
	float, captionpos = b, label = {lst:fractional-delay-fir},
	caption = {
		Implementation of a static fractional delay FIR filter.
	},
	]
	def work(self, input_items, output_items):
		inp = input_items[0]
		oup = output_items[0]
		# find the length of the highest order filter
		max_order = 2 * np.floor(np.max(self.delays)) + 1
		max_samples = np.arange(0, max_order +1)
		max_len = len(max_samples)
		# total impulse response (of all taps)
		tot_h = np.zeros(int(max_len))
		# compute for each tap
		for (a, d) in zip(self.amplitudes, self.delays):
			# compute fir coefficients
			order = 2 * np.floor(d) + 1
			samples = np.arange(0, order + 1)
			# compute impulse response
			h = a * np.sinc(samples - d)
			# correct length
			h = np.concatenate([h, np.zeros(max_len - len(h))])
			# add to other filters
			tot_h += h
		# add a LOS path if necessary
		tot_h[0] = self.los
		# compute output
		y = np.convolve(inp, tot_h)
		# add values from previous block processing
		y += np.concatenate([self.temp, np.zeros(len(y) - len(self.temp))])
		# write to output
		oup[:] = y[:len(inp)]
		# save rest for next block processing
		self.temp = y[len(inp):]
		return len(oup)
	
\end{lstlisting}

\subsection{Fading with statistical model}

In order to represent the effect of multipath fading not only statically, a second statical model fro GR was used. The block named Frequency Selective Fading Model is implemented using the algorithm from \cite{Alimohammad2009}, with the help of the sum-of sinusoid principle (SOS). The algorithm in this block is implemented with the aim that only a small number of sinusoids are needed to simulate each ray. For the simulations shown, the recommended value of 8 has been chosen.

In the block it is possible to choose between Rayleigh or Rician for the statistical modeling. When the Rician model is chosen, a realistic value for the \(K\) factor (between zero and ten \cite{Mathworks}) needs to be given. As mentioned in section \ref{sec:statistical-model}, if \(K=0\) the distribution is the same as in the Rayleigh model. The delays of the power profile are specified in samples. Those vectors depend on the environment. The magnitudes of the pulses are given in their linear values. In practice the average path gain of a fading path is in the range from \SI{-20}{\decibel} to \SI{0}{\decibel}. To add movement, some Doppler shift can be introduced, but this frequency offset needs to be normalized with the sampling rate.


%When nothing else is mentioned, the number of FIR-filter taps used is eight.


\subsubsection{Issues}

It is difficult to check, whether the noise generated from the statistical model is correct. Especially when the Doppler effect is included. Then the simulation is difficult to recreate, when the amplitude and phase parameter are not in a special state, in which the amplitude and the phase shift can be seen exactly. 
To have a possibility to verify the plot, mainly whether the movement of the signal could be correct, a Matlab model with the same values as in the GR simulation was used. Using this method, the model turned out to be correct.

\subsubsection{Example with realistic values}

\begin{table}[hb]
	\centering
	\begin{tabular}{rr rr}
		\toprule
		\bfseries Excess tap delay & \bfseries Relative power &
		\bfseries Excess tap delay & \bfseries Relative power \\
		\midrule
		\SI{   0}{\nano\second} & \(\SI{-1.0}{\decibel} \approx 0.7943\) &
		\SI{1.6}{\micro\second} & \(\SI{-3.0}{\decibel} \approx 0.5011\) \\
		\SI{  50}{\nano\second} & \(\SI{-1.0}{\decibel} \approx 0.7943\) &
		\SI{2.3}{\micro\second} & \(\SI{-5.0}{\decibel} \approx 0.3162\) \\
		\SI{ 120}{\nano\second} & \(\SI{-1.0}{\decibel} \approx 0.7943\) &
		\SI{5.0}{\micro\second} & \(\SI{-7.0}{\decibel} \approx 0.1995\) \\
		\SI{ 200}{\nano\second} & \(\SI{ 0.0}{\decibel} = 1.0000\)       \\
		\SI{ 230}{\nano\second} & \(\SI{ 0.0}{\decibel} = 1.0000\)       \\
		\SI{ 500}{\nano\second} & \(\SI{ 0.0}{\decibel} = 1.0000\)       \\
		\bottomrule
	\end{tabular}
	\caption{Extended Typical Urban model (ETU) ETSI Standard PDP values for multipath fading propagation conditions \cite{ETSI}. \label{tab:etsi-tap-values}}
\end{table}

In order to obtain a realistic simulation, the values for multipath fading propagation conditions of an Extended Typical Urban (ETU) model, from the ETSI (European Telecommunication Standards Institute) were used \cite{ETSI}. The values shown in \tabref{tab:etsi-tap-values}. For these  values the maximum Doppler frequency possibilities are predefined. In the following examples in \figref{fig:qpsk-simulations-dynamic} either \(\SI{5}{\hertz}\) or \(\SI{70}{\hertz}\) were used, opposed to the are values calculated in \eqref{eq:doppler} for a walking speed of \(\SI{2}{\meter\per\second}\), where the Doppler frequency is \(\SI{16}{\hertz}\). Those predefined values correspond to a speed of
\begin{align}
	v_1 &= \frac{\Delta f}{f_c}\cdot c_0 = \frac{\SI{5}{\hertz}}{\SI{2.4}{\giga\hertz}}\cdot \SI{3e8}{\meter\per\second}= \SI{0.625}{\meter\per\second}, \text{ and} \\
	v_2 &= \frac{\Delta f}{f_c}\cdot c_0 = \frac{\SI{70}{\hertz}}{\SI{2.4}{\giga\hertz}}\cdot \SI{3e8}{\meter\per\second}= \SI{8.75}{\meter\per\second}.
\end{align}
The numbers of taps used in this case are the number of given values.

\subsection{Measurements}

%\begin{figure}
%	\centering
%	\includegraphics[frame, width = \linewidth]{figures/screenshots/Hardware_indoor.png}
%	\caption{
%		Plot from the GNU Radio sink for an indoor environment demonstrator.
%		\label{fig:GR-Hardware-indoor}
%	}
%\end{figure}

\begin{figure}
	\centering
	\hfill
	\begin{subfigure}[t]{.47\linewidth}
		\includegraphics[frame, width = \linewidth]{figures/picture/PC210002.JPG}
		\caption{
			Reviser SDR in the outdoor environment set up.
			\label{fig:sdr1}
		}
	\end{subfigure}
	\hfill
	\begin{subfigure}[t]{.47\linewidth}
		\centering
		\includegraphics[frame, width = \linewidth]{figures/picture/PC210011.JPG}
		\caption{
		 Transmitter SDR in the outdoor environment measurement set up.
			\label{fig:sdr2}
		}
	\end{subfigure}
	\hfill
	\caption{
		Outdoor measurement set up.
		\label{fig:mesurement-set-up-outside}
	}
\end{figure}


To demonstrate the fading effect, two SDRs were used. For that multiple different measurements were made.
For example one in an indoor environment, the Lab. Another in an outdoor environment: The set up is shown in \figref{fig:mesurement-set-up-outside}.
The results of those measurements are illustrated in \figref{fig:hardware-mesurement}. Because of the USB cables the distance between the two SDRs was only about \si{2} to \SI{3}{\meter}.
The signals were all sent with a gain value of 0.4. The phase and amplitude changes could be seen well. Especially when the transmitter or the receiver were moved, the change of them got faster.

The BER, which will be described in detail in the next section, on average was 2.37 for the outdoor environment and about 2.67 for indoors. It makes sense that the fading effect occurs more in an indoor environment, because there were more surfaces for reflections. 



\subsection{Empirical BER} \label{sec:ber}

To find out how accurate the simulations are compared with a simulation of the fading effect and measurements, the bit error rate of the system is calculated. This is done with the help of a user specified \(2\)-byte test payload followed by 2 randomly generated bytes. As seen in the simplified listing \ref{lst:ber-work}, every bit is compared with the test vector after the demodulation stage. Because of the fact that the test vector has some random bits at the end, the bit error rate has always an average value of 32, even if the test vector perfectly matches. To only focus on the true BER of the signal, this value is subtracted. The random byte at the ends were introduced to help the feedback loops in the synchronization and equalization stages, which would otherwise would drift and lose their operating point. The test vector used is \texttt{[0x12, 0x48]}, because these binary numbers are easy to recognize in a continuous bit stream. For generating the bit error rate a bit stream with a specific length is compared with the test vector. To make it simpler and to avoid mistakes, additionally an average BER over the last 200 test vectors is computed.


\begin{lstlisting}[
	texcl = true, language = python, escapechar = {`},
	float, captionpos = b, label = {lst:ber-work},
	caption = {
		Custom block to compute the empirical BER.
	},
	]
	def work(self, input_items, output_items):
		# input is a list of blocks of k-bytes
		inp = input_items[0]
		# for each block
		for i in inp:
			i = np.array(i, dtype=np.uint8)
			# XOR to compute the difference
			v = np.array(self.vgl, dtype=np.uint8) ^ i
			# compute how many bits differ
			ber = sum(np.unpackbits(v))
			# save BER value
			self.ber_samples.appendleft(ber)
		# compute statistics and send to GUI
		ber_max, ber_min, ber_avg = self.ber_stats()
		self.send(self.encode([ber_max, ber_min, ber_avg]))
		return len(inp)
\end{lstlisting}

\section{Issues in the current implementation}

\subsection{Non modulated access codes} \label{sec:access-code-issue}

Currently, as described in section \ref{sec:data-frame}, the access codes are put as bytes in front of the frame in the \(k\)-byte preamble. For this to work, the access code bytes must still have a good autocorrelation function after being modulated into symbols using the chosen modulation scheme. This works well with QPSK, because the constellation is quite simple and the length of the sequence is only halved after the modulation (since QPKS has 2 bits per symbol). Thus, in the QPSK flow graph the longest known Barker sequence  \texttt{0x1f35} (13 bits, left padded with zeros) is sufficient (\(\approx 7\) symbols).

With QAM however, the complexity of the constellation and the higher number of bits per symbol makes it increasingly difficult to find binary sequences that retain a good autocorrelation after being modulated. A better solution for example would be to use a \emph{constant amplitude zero autocorrelation waveform} (CAZAC) of length \(N\), which is computed with
\begin{equation}
	u_k = \exp\left(j\frac{M\pi K}{N}\right) \text{ where }
	K = \begin{cases}
		k^2 & \text{when } N \text{ is even} \\
		k(k+1) & \text{when } N \text{ is odd},
	\end{cases}
\end{equation}
and \(M\) is relatively prime to \(N\) \cite{Chu1972}. CAZAC waveforms are ideal because they have a Dirac delta as autocorrelation \cite{Chu1972}, i.e. \(R_{uu}(\tau) = \delta(\tau)\). Unfortunately, since these complex values are not on any constellation point they break some assumptions of the polyphase clock sync and the LMD DD equalizer (but not CMA, which unfortunately cannot be used for QAM). Thus, to use CAZAC waveforms, the transmitter needs to put them in front of the modulated symbols (for example using a correctly parametrized Stream MUX block in GR), and the receiver would need to synchronize with the sequence before the clock recovery or equalization. The latter is especially problematic because then it is no longer possible to identify the peak by comparing the autocorrelation value to a fixed threshold as done in section \ref{sec:implement-phasecorr}. Because of the aforementioned reasons, in its current state the QAM flow graph is unable to lock and decode any signals.


%TODO: Show figurs with QAM 


\subsection{Single threaded GUI application} \label{sec:gui-issue-single-threaded}

The current GUI prototype built with DearPyGUI has some issues, the most critical begin that it is a single-threaded program. The interprocess communications (with GR's flow graphs) should be on a separate thread from the graphics, what is currently not the case. The problem is not noticeable as long as the flow graphs in the background keep sending data, but as soon as the UDP/IP data stream stops, the timeout of the socket interface causes the interface to run at less than 20 frames per second.

\subsection{Clock synchronization issues}

Unfortunately the two SDR need an external clock generator. For that a Rubidium Frequency standard device (Model FS725) is used with the clock frequency of \SI{10}{\mega\hertz}. Two of them are used to make them more movable and independent. Those clock generators  where needed, because the synchronization does not work as planed in \ref{sec:preforming-implementation}. 
%TODO: Right squenz?
Without those only the amplitudes could be seen in the plots. 


% TODO: Plots from the Hardware

\section{Produced constellation plots}
In the next few pages the constellation plots from the simulation results and the hardware measurements are shown.
% TODO anayl
%TODO achsenbeschrieftung 


\newgeometry{
	top = 25mm, bottom = 25mm,
	inner = 15mm, outer = 15mm,
}
\begin{figure}
	\centering
	\input{figures/tikz/qpsk-simulations-static}
	\caption{
		Simulations of a static fading channel models with different tap values. The samples were generated using the custom block discussed in section \ref{sec:discrete-time-model-fir}. For the 1 tap model the fading tap was \(0.2\delta(n - 0.25)\), and for the 4 tap model uses \(0.2 \delta(n - 0.25) + 0.08 \delta(n - 3.25) + 0.5 \delta(n - 4) + 0.4 \delta(n - 6.3)\). In both cases the delays are given in samples.
		% delay = [0.25, 3.25, 4, 6.3]
		% ampl =  [0.2, 0.08, 0.5, 0.4]
		\label{fig:qpsk-simulations-static}
	}
\end{figure}
\newpage
\begin{figure}
	\centering
	\input{figures/tikz/qpsk-simulations-dynamic}
	\caption{
		Simulations with QPSK modulation and a dynamic fading channel model that uses PDP values of the Extended Typical Urban model (ETU) from the ETSI standard normative Annex B.2 in \cite{ETSI}. The PDP values from the standard are reported in \tabref{tab:etsi-tap-values}. The color gradient represents progression in time: yellow samples are more recent than the blue samples.
		\label{fig:qpsk-simulations-dynamic}
	}
\end{figure}
\newpage
\begin{figure}
	\centering
	\input{figures/tikz/qam-simulations-dynamic}
	\caption{
		Simulation using the same channel parameters as in \figref{fig:qpsk-simulations-dynamic}, but with a 16-ary QAM modulation scheme. Unfortunately, because of the issue discussed in section \ref{sec:access-code-issue} the receiver cannot lock and decode the envelope.
	}
\end{figure}
\begin{figure}
	\centering
	\input{figures/tikz/hardware}
	\caption{
		Constellation plots from the different measurements with the two SDRs.
		\label{fig:hardware-mesurement}
	}
\end{figure}
\newpage
\restoregeometry