{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Lab 4: Frequency Analysis -- Calibration of the SDR frequency using GSM signals
\n",
"\n",
"\n",
"\n",
"\n",
"### Michael Lustig \n",
"\n",
"#### Note: The instruction in this lab are now less detailed and require some thought and common sense with respect to choosing parameters and implementations. You should probably make sure you plot various temporary results (even though we may not request them) in order to guarentee your code is correct -- like the results of filter design or filtering....\n",
"\n",
"* The lab was inspired by the Kalibrate project by Joshua Lackey\n",
"\n",
"We hope you enjoy it!\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# Import functions and libraries\n",
"from __future__ import division\n",
"import numpy as np, matplotlib.pyplot as plt, time, IPython.display\n",
"from numpy import *\n",
"from numpy.fft import *\n",
"from matplotlib.pyplot import *\n",
"from rtlsdr import RtlSdr\n",
"from scipy import signal \n",
"import pyaudio\n",
"import time\n",
"%matplotlib inline"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# function to compute average power spectrum\n",
"def avgPS( x, N=256, fs=1):\n",
" M = floor(len(x)/N)\n",
" x_ = reshape(x[:M*N],(M,N)) * np.hamming(N)[None,:]\n",
" X = np.fft.fftshift(np.fft.fft(x_,axis=1),axes=1)\n",
" return r_[-N/2.0:N/2.0]/N*fs, mean(abs(X**2),axis=0)\n",
"\n",
"\n",
"# Plot an image of the spectrogram y, with the axis labeled with time tl,\n",
"# and frequency fl\n",
"#\n",
"# t_range -- time axis label, nt samples\n",
"# f_range -- frequency axis label, nf samples\n",
"# y -- spectrogram, nf by nt array\n",
"# dbf -- Dynamic range of the spect\n",
"\n",
"def sg_plot( t_range, f_range, y, dbf = 60, fig = None) :\n",
" eps = 10.0**(-dbf/20.0) # minimum signal\n",
" \n",
" # find maximum\n",
" y_max = abs(y).max()\n",
" \n",
" # compute 20*log magnitude, scaled to the max\n",
" y_log = 20.0 * np.log10( (abs( y ) / y_max)*(1-eps) + eps )\n",
" \n",
" # rescale image intensity to 256\n",
" img = 256*(y_log + dbf)/dbf - 1\n",
" \n",
" \n",
" fig=figure(figsize=(16,6));\n",
" \n",
" plt.imshow( np.flipud( 64.0*(y_log + dbf)/dbf ), extent= t_range + f_range ,cmap=plt.cm.gray, aspect='auto');\n",
" plt.xlabel('Time, s');\n",
" plt.ylabel('Frequency, Hz');\n",
" plt.tight_layout();\n",
" \n",
"def myspectrogram_hann_ovlp(x, m, fs, fc, dbf = 60, fig = None):\n",
" # Plot the spectrogram of x.\n",
" # First take the original signal x and split it into blocks of length m\n",
" \n",
" # check if x is real\n",
" isreal_bool = isreal(x).all()\n",
" \n",
" # pad x up to a multiple of m \n",
" lx = len(x);\n",
" nt = (lx + m - 1) // m\n",
" x = append(x,zeros(-lx+nt*m))\n",
" x = x.reshape((m//2,nt*2), order='F')\n",
" x = concatenate((x,x),axis=0)\n",
" x = x.reshape((m*nt*2,1),order='F')\n",
" x = x[r_[m//2:len(x),ones(m//2)*(len(x)-1)].astype(int)].reshape((m,nt*2),order='F')\n",
" \n",
" \n",
" xmw = x * hanning(m)[:,None];\n",
" \n",
" # frequency index\n",
" t_range = [0.0, lx / fs]\n",
" \n",
" if isreal_bool:\n",
" f_range = [ fc, fs / 2.0 + fc]\n",
" xmf = fft(xmw,len(xmw),axis=0)\n",
" fig = sg_plot(t_range, f_range, xmf[0:m//2,:], dbf, fig);\n",
" else:\n",
" f_range = [-fs / 2.0 + fc, fs / 2.0 + fc]\n",
" xmf = abs(fftshift( fft( xmw ,len(xmw),axis=0), axes=0 ))\n",
" fig = sg_plot(t_range, f_range, xmf, dbf, fig);\n",
" \n",
" return fig"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## RTL-SDR frequency Drift\n",
"\n",
"The rtl-sdr has has a 28.8MHz crystal oscillator, which sets the reference frequency used for demodulation on the tuner. The crystal is known to have poor quallity and therefore its frequency can drift considerably with temperature. The drift is measured in parts per millions (PPM) and can go up to $\\pm50$ PPM. The result of such a drift is that when you tune to a certain frequency, you will actually be tunning to an offset frequency. This offset is proportional to the center frequency, so the absolute frequency offset gets worse for higher frequencies. \n",
"\n",
"This sort of oscillator drifts are very common in other devices, such as your cellphone. The wireless industry has come up with many different techniques in which the clock, or the absolute frequency can be corrected for. In this lab we will use the GSM cellphone network to correct for the frequency drift of the rtl-sdr. GSM base stations are required to have an accurate clock within 0.05 ppm, so they are an excellent source of \"accurate\" signals. This is very much the first step that a GSM based phone would do in order to connect to the network. \n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Before we start talking about how to use the GSM network to calibrate our rtl-sdr, let's learn a bit about the system\n",
"\n",
"\n",
"## GSM : Global System for Mobile Communications\n",
"\n",
"This section is heavily based on \"GSM for dummies\" (http://www.pennula.de/datenarchiv/gsm-for-dummies.pdf)\n",
"\n",
"GSM is a popular digital cellular network. It is based on time-division-multiple-access (TDMA), which means that several phones use the same channel by time interleaving. GSM operates in several different carrier frequencies. In the US and Canada and number of other countries, the frequency bands are 850MHz and 1900MHz. In this lab we will use the 850MHz frequency band, since for the majority of you have the 820T tuner based rtl-sdrs, and the 1900MHz band is outside of its range (28-1700MHz). \n",
"\n",
"### Frequencies\n",
"For the GSM-850 the band is separated into Uplink frequencies 824-849MHz and downlink frequencies 869-894MHz. The base station transmits on the downlink frequencies and mobile phones transmit on the uplink frequencies. GSM operates in duplex mode in which the phones and base stations receive and transmit at the same time. This is the reason that the uplink and the downlink are separated by 50MHz!. A channel number, also known as Absolute Radio Frequency Channel Number (ARFCN) is assigned a pair of uplink and downlink frequencies. There are 124 channels numbered 128-251 in the GSM-850 system. Each channel is allocated a bandwidth of 200KHz. For the downlink frequencies $f_N = 869.2+0.2*(ARFCN-128)$ MHz. \n",
"\n",
"### Modulation\n",
"GSM uses Gaussian Minimum Shift Keying (GMSK in short) as a modulation method. We will cover this in part II of this lab. But, in short.... it is a binary digital Frequency Modulated scheme in which the phase between each bit period is continuouse. Bits are encoded as different frequencies separated by one-half the bit-rate. The modulation rate in GSM is $1625/6 \\approx 270.833$ kb/s. For more information read (http://en.wikipedia.org/wiki/Minimum-shift_keying) "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Time-Division-Multiple-Acess\n",
"\n",
"GSM uses time-division to share a frequency channel among users. Each frequency is divided into blocks of time that are known as time-slots. There are 8 time-slots that are numberes TS0 - TS7. Each time-slot lasts about 576.9 $\\mu s$. Given the bit-rate above, a total of 156.25 bits can be transmitted in each slot. Each slot allocates 8.25 of it \"bits time\" as guard-time, split between the beginning and the end of each time-slot. This time is necessary to prevent overlapping between the different time slots. In addition, 3 bits on both sides of the time-slot do not contain any data and are there to allow for ramping the amplifiers up and down. \n",
"\n",
"Data transmitted within each time-slot is called a burst. There are several type of burts, with different purposes. We are interested in the so called \"frequency correction\" burst, which is a burst of a pure frequency tone at 1/4th the bitrate of GSM or (1625000 / 6) / 4 = 67708.3 Hz. \n",
"\n",
"\n",
"