Addition of Signal Analyzer tool
"signal_analyzer.py" takes in a raw complex64 file as input and it
extracts period, duty-cycle, and amplitude.
Change-Id: I319c8b90e7cd0ffec4966a6664091d96c23efb59
diff --git a/bandpass_recorder/README b/bandpass_recorder/README
index f796401..44c3072 100644
--- a/bandpass_recorder/README
+++ b/bandpass_recorder/README
@@ -1,5 +1,5 @@
Interference Recorder Script 1.0
-=====================================================================
+===============================================================================
This script takes in a filename as in input and opens a gui which allows
the user to adjust the frequency, bandwidth, gain, and record time. Upon
@@ -9,13 +9,13 @@
Usage:
./bandpass_recorder.py filename
----------------------------------------------------------------------
+-------------------------------------------------------------------------------
Requirements:
- USRP Software Defined Radio
- GNURadio and UHD Driver
- "bandpass_recorder_file.py" and "bandpass_recorder_gui.py"
----------------------------------------------------------------------
+-------------------------------------------------------------------------------
Changing radio behaviour:
The recorder behaviour can be changed by modifying "bandpass_recorder_file.py"
or "bandpass_recorder_gui.py" in 2 ways:
@@ -27,3 +27,15 @@
this is used by the script to specify the output filename. The script
will output the data stream to recorder.out and it will then be renamed
to the specified filename.
+
+===============================================================================
+Signal Analysis Tool
+
+"signal_analyzer.py" can be used to analyze recordings in order to extract
+parameters from the data files created by the recorder script.
+
+Use "./signal_analyzer.py -h" to view options
+
+Note that large window/display values can crash your system as they are very
+memory intensive. Window values should be less than or equal to 1e8 and
+display values should be less than or equal to 1e7.
diff --git a/bandpass_recorder/bandpass_recorder.py b/bandpass_recorder/bandpass_recorder.py
index fb5ae30..053df1d 100755
--- a/bandpass_recorder/bandpass_recorder.py
+++ b/bandpass_recorder/bandpass_recorder.py
@@ -1,5 +1,11 @@
#!/usr/bin/python
-import ctypes
+"""Bandpass Recorder Tool.
+
+This script uses the URSP and records signals from the air. The user can
+specify the center-frequency, gain value, width of the recording and the
+duration of the recording.
+"""
+
from distutils.version import StrictVersion
import os
import sys
@@ -8,67 +14,88 @@
from gnuradio import gr
from PyQt4 import Qt
-# Import modules from gnuradio-companion generated files.
from bandpass_recorder_file import bandpass_recorder_file
from bandpass_recorder_gui import bandpass_recorder_gui
+import options
+optspec = """
+bandpass_recorder.py [options...] filename
+--
+v,visual Displays a GUI for viewing the region to be recorded.
+f,frequency= Specify the center frequency to record from in MHz [2412]
+g,gain= Specify the decibel gain of the recording. [50.0]
+w,width= Specify the width of the interval of the recording in KHz. [20000]
+d,duration= Specify the length of the recording in seconds. [5]
+"""
+
+
+def record(gain, freq, width, duration, filename):
+ """Record a specified signal and store it into a file.
+
+ Args:
+ gain: Amplification of the signal. (dB)
+ freq: Center Frequency of the recording (MHz).
+ width: Width of the bandpass filter. (KHz)
+ duration: Length of the recording. (s)
+ filename: Path to the file that the recording is stored into.
+ """
+ print 'test'
+ file_block = bandpass_recorder_file()
+ print 'recording... ({0} seconds)\n'.format(duration)
+ file_block.set_gain(gain)
+ file_block.set_freq(freq * 1e6)
+ file_block.set_width(width * 1e3)
+
+ file_block.start()
+ time.sleep(duration)
+
+ file_block.stop()
+ file_block.wait()
+
+ # Rename output file
+ try:
+ os.rename('recorder.out', filename)
+ except OSError as e:
+ print '\nOSError [{0}]: {1}'.format(e.errno, e.strerror)
+ print 'Recording stored in \'recorder.out\'\n'
+
+
+def main():
+ o = options.Options(optspec)
+ opt, flags, extra = o.parse(sys.argv[1:])
+
+ if len(extra) != 1:
+ o.fatal('Filename missing or extra arguments present')
+ sys.exit(1)
+
+ recording_filename = extra[0]
+
+ if StrictVersion(Qt.qVersion()) >= StrictVersion('4.5.0'):
+ prefs = gr.prefs().get_string('qtgui', 'style', 'raster')
+ Qt.QApplication.setGraphicsSystem(prefs)
+ qapp = Qt.QApplication(sys.argv)
+
+ if opt.visual:
+ gui_block = bandpass_recorder_gui()
+ gui_block.start()
+ gui_block.show()
+
+ # Define callback behaviour for the gui_block: starts the file_block and
+ # begins recording once the gui_block is closed.
+ def gui_callback():
+ gain = gui_block.get_gain()
+ freq = gui_block.get_freq() / 1000000.0
+ width = gui_block.get_width()
+ record_time = gui_block.get_record_time() / 1000.0
+
+ gui_block.stop()
+ gui_block.wait()
+ record(gain, freq, width, record_time, recording_filename)
+ qapp.connect(qapp, Qt.SIGNAL('aboutToQuit()'), gui_callback)
+ qapp.exec_()
+ else:
+ record(float(opt.gain), float(opt.frequency), float(opt.width),
+ float(opt.duration), recording_filename)
if __name__ == '__main__':
- # Checks if libX11 works properly.
- if sys.platform.startswith('linux'):
- try:
- x11 = ctypes.cdll.LoadLibrary('libX11.so')
- x11.XInitThreads()
- except:
- print 'Warning: failed to XInitThreads()'
-
- # Check and parse arguments.
- if len(sys.argv) != 2:
- print 'usage: ./bandpass_recorder_final.py filename\n'
- sys.exit(1)
- if StrictVersion(Qt.qVersion()) >= StrictVersion('4.5.0'):
- Qt.QApplication.setGraphicsSystem(gr.prefs().get_string('qtgui', 'style', 'raster'))
- qapp = Qt.QApplication(sys.argv)
-
- # Start the GUI
- gui_block = bandpass_recorder_gui()
- gui_block.start()
- gui_block.show()
-
- # Define callback behaviour for the gui_block: starts the file_block and begins
- # recording once the gui_block is closed.
- def begin_recording():
- # Retrieve parameters
- gain = gui_block.get_gain()
- freq = gui_block.get_freq()
- width = gui_block.get_width()
- record_time = gui_block.get_record_time() / 1000.0
-
- # Stop the gui
- gui_block.stop()
- gui_block.wait()
-
- # Start the file recorder
- file_block = bandpass_recorder_file()
- print 'recording... ({0} seconds)\n'.format(record_time)
-
- # Set parameters
- file_block.set_gain(gain)
- file_block.set_freq(freq)
- file_block.set_width(width)
-
- file_block.start()
- # Run the recording for the specified time.
- time.sleep(record_time)
-
- file_block.stop()
- file_block.wait()
-
- # Rename output file
- try:
- os.rename('recorder.out', sys.argv[1])
- except OSError as e:
- print '\nOSError [{0}]: {1}'.format(e.errno, e.strerror)
- print 'Recording stored in \'recorder.out\'\n'
- qapp.connect(qapp, Qt.SIGNAL('aboutToQuit()'), begin_recording)
- qapp.exec_()
+ main()
diff --git a/bandpass_recorder/bandpass_recorder_gui.py b/bandpass_recorder/bandpass_recorder_gui.py
index 29ffd63..6e2a57b 100755
--- a/bandpass_recorder/bandpass_recorder_gui.py
+++ b/bandpass_recorder/bandpass_recorder_gui.py
@@ -106,16 +106,16 @@
self.qtgui_sink_x_0.set_update_time(1.0/10)
self._qtgui_sink_x_0_win = sip.wrapinstance(self.qtgui_sink_x_0.pyqwidget(), Qt.QWidget)
self.top_layout.addWidget(self._qtgui_sink_x_0_win)
-
+
self.qtgui_sink_x_0.enable_rf_freq(False)
-
-
-
+
+
+
##################################################
# Connections
##################################################
- self.connect((self.uhd_usrp_source_0, 0), (self.qtgui_sink_x_0, 0))
+ self.connect((self.uhd_usrp_source_0, 0), (self.qtgui_sink_x_0, 0))
def closeEvent(self, event):
self.settings = Qt.QSettings("GNU Radio", "bandpass_recorder_gui")
diff --git a/bandpass_recorder/signal_analyzer.py b/bandpass_recorder/signal_analyzer.py
new file mode 100755
index 0000000..a3a4ce5
--- /dev/null
+++ b/bandpass_recorder/signal_analyzer.py
@@ -0,0 +1,130 @@
+#!/usr/bin/python
+"""Signal Analysis Program.
+
+This program takes in a recording of an interference pattern and extracts the
+frequency, duty-cycle, and square of the signal so that the signal can be
+regenerated using parameters.
+"""
+import sys
+
+import matplotlib.pyplot as plt
+import numpy as np
+from scipy.signal import decimate
+from scipy.signal import hilbert
+from scipy.signal import symiirorder1
+
+import options
+
+SAMPLE_RATE = 30e6
+DOWNSAMPLE_FACTOR = 1000
+MAX_DUTY_PERIOD_USEC = 10000.0
+MIN_DUTY_PERIOD = SAMPLE_RATE / DOWNSAMPLE_FACTOR / MAX_DUTY_PERIOD_USEC
+optspec = """
+signal_analyzer [options...] filename
+--
+w,window= Specify the number of samples to be analyzed. [100000000]
+d,display Display the raw data.
+o,output= Specify the output filename for the profile.
+"""
+
+
+def get_active_rms(raw_data, active_indices):
+ """Calculates the rms of the signal when its on."""
+ square_data = np.square(np.absolute(raw_data))
+ raw_active_indices = np.repeat(active_indices, DOWNSAMPLE_FACTOR)
+ active_values = square_data[raw_active_indices]
+
+ # TODO(sunnyd): There's a lot of values here that are not part of the active
+ # signal, need to find a way to parse them out in the future
+ return np.sqrt(np.average(active_values))
+
+
+def get_period(switch_indices):
+ """Calculates the frequency of the signal peaks."""
+ period_samples = []
+ for i in range(0, len(switch_indices) - 1):
+ period_samples.append(switch_indices[i + 1] - switch_indices[i])
+ return np.median(period_samples) * DOWNSAMPLE_FACTOR / SAMPLE_RATE
+
+
+def get_switch_indices(active_indices):
+ # Locate the indices where the signal is on
+ off_switch_indices = []
+ on_switch_indices = []
+ for i in range(len(active_indices) - 1):
+ if active_indices[i + 1] - active_indices[i] > MIN_DUTY_PERIOD:
+ off_switch_indices.append(active_indices[i])
+ on_switch_indices.append(active_indices[i+1])
+ return (off_switch_indices, on_switch_indices)
+
+
+def save_profile(raw_data, off_switch_indices, on_switch_indices, filename):
+ """Extract the fourier transform of the signal."""
+ # TODO(sunnyd): Find an average of many peaks rather than choosing an
+ # arbitrary one.
+ cur_slice = raw_data[off_switch_indices[3] * DOWNSAMPLE_FACTOR :
+ on_switch_indices[3] * DOWNSAMPLE_FACTOR]
+ fourier = np.fft.fft(cur_slice, int(SAMPLE_RATE))
+ fourier_down = decimate(fourier, DOWNSAMPLE_FACTOR, 2, ftype='fir')
+
+ # Filter out some of the noise
+ fourier_mean = np.average(np.absolute(fourier_down))
+
+ # TODO(sunnyd): 1.5 is a bit arbitrary, need a good metric for this. STD Dev
+ # might work but could filter out a lot of useful signal if the entire band
+ # is high
+ fourier_down[np.absolute(fourier_down) < fourier_mean * 1.5] = 0
+ np.save(filename, np.array(fourier_down))
+
+
+def get_duty_cycle(active_bitmap, switch_indices):
+ """Calculate the Duty cycle of the signal."""
+ duty_cycles = []
+ for i in range(len(switch_indices) - 1):
+ duty_cycles.append(np.average(active_bitmap[switch_indices[i] :
+ switch_indices[i + 1]]))
+ return np.median(duty_cycles)
+
+
+def main():
+ o = options.Options(optspec)
+ opt, _, extra = o.parse(sys.argv[1:])
+
+ window_size = int(opt.window)
+ if len(extra) != 1:
+ o.fatal('Missing filename or extra arguments')
+ if not opt.output:
+ o.fatal('Missing output option required')
+
+ f = open(extra[0], 'rb')
+ raw_data = np.fromfile(f, dtype=np.complex64, count=window_size)
+ downsampled_data = decimate(raw_data, DOWNSAMPLE_FACTOR, 2, ftype='fir')
+
+ hilbert_data = hilbert(np.absolute(downsampled_data))
+ envelope = symiirorder1(hilbert_data, 1, 0.5)
+ abs_envelope = abs(envelope)
+
+ if opt.display:
+ display_size = min(opt.window, 10e6)
+ if opt.window > 10e6:
+ sys.stderr.write('Warning: window size too large to display. Using '
+ 'default value of 10e6\n')
+ plt.plot(downsampled_data[0: display_size].real, color='red', hold=True)
+ plt.plot(abs_envelope)
+ plt.show()
+
+ # Calculate useful values for processing.
+ mean = np.average(abs_envelope)
+ active_bitmap = abs_envelope > mean
+ active_indices = active_bitmap.nonzero()[0]
+
+ off_switch_indices, on_switch_indices = get_switch_indices(active_indices)
+ duty = get_duty_cycle(active_bitmap, off_switch_indices)
+ period = get_period(off_switch_indices)
+ rms = get_active_rms(raw_data, active_bitmap)
+ save_profile(raw_data, off_switch_indices, on_switch_indices, opt.output)
+ print 'Duty-Cycle: %f, Period: %fs, Active RMS: %f' % (duty, period, rms)
+
+
+if __name__ == '__main__':
+ main()
diff --git a/options.py b/options.py
new file mode 100644
index 0000000..7fb5bcf
--- /dev/null
+++ b/options.py
@@ -0,0 +1,274 @@
+# Copyright 2010-2012 Avery Pennarun and options.py contributors.
+# All rights reserved.
+#
+# (This license applies to this file but not necessarily the other files in
+# this package.)
+#
+# Redistribution and use in source and binary forms, with or without
+# modification, are permitted provided that the following conditions are
+# met:
+#
+# 1. Redistributions of source code must retain the above copyright
+# notice, this list of conditions and the following disclaimer.
+#
+# 2. Redistributions in binary form must reproduce the above copyright
+# notice, this list of conditions and the following disclaimer in
+# the documentation and/or other materials provided with the
+# distribution.
+#
+# THIS SOFTWARE IS PROVIDED BY AVERY PENNARUN ``AS IS'' AND ANY
+# EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
+# PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL <COPYRIGHT HOLDER> OR
+# CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
+# EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+# PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+# PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+# LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+# NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+# SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+#
+"""Command-line options parser.
+With the help of an options spec string, easily parse command-line options.
+
+An options spec is made up of two parts, separated by a line with two dashes.
+The first part is the synopsis of the command and the second one specifies
+options, one per line.
+
+Each non-empty line in the synopsis gives a set of options that can be used
+together.
+
+Option flags must be at the begining of the line and multiple flags are
+separated by commas. Usually, options have a short, one character flag, and a
+longer one, but the short one can be omitted.
+
+Long option flags are used as the option's key for the OptDict produced when
+parsing options.
+
+When the flag definition is ended with an equal sign, the option takes one
+string as an argument. Otherwise, the option does not take an argument and
+corresponds to a boolean flag that is true when the option is given on the
+command line.
+
+The option's description is found at the right of its flags definition, after
+one or more spaces. The description ends at the end of the line. If the
+description contains text enclosed in square brackets, the enclosed text will
+be used as the option's default value.
+
+Options can be put in different groups. Options in the same group must be on
+consecutive lines. Groups are formed by inserting a line that begins with a
+space. The text on that line will be output after an empty line.
+"""
+import sys, os, textwrap, getopt, re, struct
+
+
+def _invert(v, invert):
+ if invert:
+ return not v
+ return v
+
+
+def _remove_negative_kv(k, v):
+ if k.startswith('no-') or k.startswith('no_'):
+ return k[3:], not v
+ return k,v
+
+
+class OptDict(object):
+ """Dictionary that exposes keys as attributes.
+
+ Keys can be set or accessed with a "no-" or "no_" prefix to negate the
+ value.
+ """
+ def __init__(self, aliases):
+ self._opts = {}
+ self._aliases = aliases
+
+ def _unalias(self, k):
+ k, reinvert = _remove_negative_kv(k, False)
+ k, invert = self._aliases[k]
+ return k, invert ^ reinvert
+
+ def __setitem__(self, k, v):
+ k, invert = self._unalias(k)
+ self._opts[k] = _invert(v, invert)
+
+ def __getitem__(self, k):
+ k, invert = self._unalias(k)
+ return _invert(self._opts[k], invert)
+
+ def __getattr__(self, k):
+ return self[k]
+
+
+def _default_onabort(msg):
+ sys.exit(97)
+
+
+def _intify(v):
+ try:
+ vv = int(v or '')
+ if str(vv) == v:
+ return vv
+ except ValueError:
+ pass
+ return v
+
+
+def _atoi(v):
+ try:
+ return int(v or 0)
+ except ValueError:
+ return 0
+
+
+def _tty_width():
+ s = struct.pack("HHHH", 0, 0, 0, 0)
+ try:
+ import fcntl, termios
+ s = fcntl.ioctl(sys.stderr.fileno(), termios.TIOCGWINSZ, s)
+ except (IOError, ImportError):
+ return _atoi(os.environ.get('WIDTH')) or 70
+ (ysize,xsize,ypix,xpix) = struct.unpack('HHHH', s)
+ return xsize or 70
+
+
+class Options:
+ """Option parser.
+ When constructed, a string called an option spec must be given. It
+ specifies the synopsis and option flags and their description. For more
+ information about option specs, see the docstring at the top of this file.
+
+ Two optional arguments specify an alternative parsing function and an
+ alternative behaviour on abort (after having output the usage string).
+
+ By default, the parser function is getopt.gnu_getopt, and the abort
+ behaviour is to exit the program.
+ """
+ def __init__(self, optspec, optfunc=getopt.gnu_getopt,
+ onabort=_default_onabort):
+ self.optspec = optspec
+ self._onabort = onabort
+ self.optfunc = optfunc
+ self._aliases = {}
+ self._shortopts = 'h?'
+ self._longopts = ['help', 'usage']
+ self._hasparms = {}
+ self._defaults = {}
+ self._usagestr = self._gen_usage() # this also parses the optspec
+
+ def _gen_usage(self):
+ out = []
+ lines = self.optspec.strip().split('\n')
+ lines.reverse()
+ first_syn = True
+ while lines:
+ l = lines.pop()
+ if l == '--': break
+ out.append('%s: %s\n' % (first_syn and 'usage' or ' or', l))
+ first_syn = False
+ out.append('\n')
+ last_was_option = False
+ while lines:
+ l = lines.pop()
+ if l.startswith(' '):
+ out.append('%s%s\n' % (last_was_option and '\n' or '',
+ l.lstrip()))
+ last_was_option = False
+ elif l:
+ (flags,extra) = (l + ' ').split(' ', 1)
+ extra = extra.strip()
+ if flags.endswith('='):
+ flags = flags[:-1]
+ has_parm = 1
+ else:
+ has_parm = 0
+ g = re.search(r'\[([^\]]*)\]$', extra)
+ if g:
+ defval = _intify(g.group(1))
+ else:
+ defval = None
+ flagl = flags.split(',')
+ flagl_nice = []
+ flag_main, invert_main = _remove_negative_kv(flagl[0], False)
+ self._defaults[flag_main] = _invert(defval, invert_main)
+ for _f in flagl:
+ f,invert = _remove_negative_kv(_f, 0)
+ self._aliases[f] = (flag_main, invert_main ^ invert)
+ self._hasparms[f] = has_parm
+ if f == '#':
+ self._shortopts += '0123456789'
+ flagl_nice.append('-#')
+ elif len(f) == 1:
+ self._shortopts += f + (has_parm and ':' or '')
+ flagl_nice.append('-' + f)
+ else:
+ f_nice = re.sub(r'\W', '_', f)
+ self._aliases[f_nice] = (flag_main,
+ invert_main ^ invert)
+ self._longopts.append(f + (has_parm and '=' or ''))
+ self._longopts.append('no-' + f)
+ flagl_nice.append('--' + _f)
+ flags_nice = ', '.join(flagl_nice)
+ if has_parm:
+ flags_nice += ' ...'
+ prefix = ' %-20s ' % flags_nice
+ argtext = '\n'.join(textwrap.wrap(extra, width=_tty_width(),
+ initial_indent=prefix,
+ subsequent_indent=' '*28))
+ out.append(argtext + '\n')
+ last_was_option = True
+ else:
+ out.append('\n')
+ last_was_option = False
+ return ''.join(out).rstrip() + '\n'
+
+ def usage(self, msg=""):
+ """Print usage string to stderr and abort."""
+ sys.stderr.write(self._usagestr)
+ if msg:
+ sys.stderr.write(msg)
+ e = self._onabort and self._onabort(msg) or None
+ if e:
+ raise e
+
+ def fatal(self, msg):
+ """Print an error message to stderr and abort with usage string."""
+ msg = '\nerror: %s\n' % msg
+ return self.usage(msg)
+
+ def parse(self, args):
+ """Parse a list of arguments and return (options, flags, extra).
+
+ In the returned tuple, "options" is an OptDict with known options,
+ "flags" is a list of option flags that were used on the command-line,
+ and "extra" is a list of positional arguments.
+ """
+ try:
+ (flags,extra) = self.optfunc(args, self._shortopts, self._longopts)
+ except getopt.GetoptError, e:
+ self.fatal(e)
+
+ opt = OptDict(aliases=self._aliases)
+
+ for k,v in self._defaults.iteritems():
+ opt[k] = v
+
+ for (k,v) in flags:
+ k = k.lstrip('-')
+ if k in ('h', '?', 'help', 'usage'):
+ self.usage()
+ if (self._aliases.get('#') and
+ k in ('0','1','2','3','4','5','6','7','8','9')):
+ v = int(k) # guaranteed to be exactly one digit
+ k, invert = self._aliases['#']
+ opt['#'] = v
+ else:
+ k, invert = opt._unalias(k)
+ if not self._hasparms[k]:
+ assert(v == '')
+ v = (opt._opts.get(k) or 0) + 1
+ else:
+ v = _intify(v)
+ opt[k] = _invert(v, invert)
+ return (opt,flags,extra)
diff --git a/signal_generator/options.py b/signal_generator/options.py
new file mode 120000
index 0000000..3508154
--- /dev/null
+++ b/signal_generator/options.py
@@ -0,0 +1 @@
+../options.py
\ No newline at end of file
