[Top][All Lists]

RE: [Discuss-gnuradio] Compare Phase block

 From: Geof Nieboer Subject: RE: [Discuss-gnuradio] Compare Phase block Date: Sun, 6 Apr 2008 19:56:53 +0900

My VOR receiver is almost there, but not quite...

My code is below.  I've demodulated the first 30 Hz signal from AM, and the
second from the FM subcarrier.  Both appear pretty solid on the O-scope
(except for a DC offset on the AM signal), and visually I can see a constant
difference in phase on the scope.

However, my attempts to get relative phase aren't working.  First I tried
exactly what was listed below, and I got either 0 or pi as a phase
difference.  Then I realized that the demodulators were returning a real
signal, and type converting it back to complex was leaving off the imaginary
portion.  So I used a Hilbert transform to re-create the imaginary portion,
and then fed those to the phase difference part of the graph.

Unfortunately, I'm not getting a constant return, it's basically 'spinning',
the returned value is constantly changing in a cyclical manner.
I tried another method to determine the relative phase, but it didn't do
anything different.

I figure it's either got to do with the DC offset, and I'm not sure the best
way to get rid of this because my signal is so low frequency, or I'm using
the Hilbert filter wrong.

Below is the source code, it's not cleaned up yet, but it's functional.  To
use you'd need to find a VOR station and get near it.  They really aren't
designed for ground use, so you may need to get close.
But like I said, the signals look clean on the scope, so I'm baffled.

Thoughts?

Geof

------------------------

> You can get the phase difference by multiplying one signal by the
> complex conjugate of the other, then take the arg (arctan) of that.
> There are blocks for each of those operations.
>
> sig0 = ...
> sig1 = ..
>
> conj = gr.conjugate_cc()
> mult = gr.multiply_cc()
> arg = gr.complex_to_arg()
>
> dst = ... # stream of float radians
>
> connect(sig0, conj, mult)
> connect(sig1, (mult, 1))
> connect(mult, arg)
> connect(arg, dst)

--------------------

#!/usr/bin/env python
#
# Copyright 2005,2006,2007,2008 Free Software Foundation, Inc.
#
# This file is part of GNU Radio
#
# GNU Radio is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 3, or (at your option)
# any later version.
#
# GNU Radio is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with GNU Radio; see the file COPYING.  If not, write to
# the Free Software Foundation, Inc., 51 Franklin Street,
# Boston, MA 02110-1301, USA.
#
# This module tunes to VOR receivers and determines the bearing from the
# to the receiver
#
# See http://www.doe.carleton.ca/courses/ELEC4504/nav_r7.pdf for a quick
# VOR theory discussion

from gnuradio import gr, gru, eng_notation, optfir
from gnuradio import audio
from gnuradio import usrp
from gnuradio import blks2
from gnuradio.eng_option import eng_option
from gnuradio.wxgui import slider, powermate, numbersink2
from gnuradio.wxgui import stdgui2, fftsink2, form, scopesink2
from optparse import OptionParser
from usrpm import usrp_dbid
import sys
import math
import wx

def pick_subdevice(u):
"""
The user didn't specify a subdevice on the command line.
Try for one of these, in order: BASIC_RX,TV_RX, BASIC_RX, whatever is on
side A.

@return a subdev_spec
"""
return usrp.pick_subdev(u, (usrp_dbid.BASIC_RX,
usrp_dbid.LF_RX,
usrp_dbid.TV_RX,
usrp_dbid.TV_RX_REV_2,
usrp_dbid.TV_RX_REV_3))

class vor_rx_block (stdgui2.std_top_block):
def __init__(self,frame,panel,vbox,argv):
stdgui2.std_top_block.__init__ (self,frame,panel,vbox,argv)

parser=OptionParser(option_class=eng_option)
default=None,
help="select USRP Rx side A or B (default=A)")
default=1008.0e3,
help="set frequency to FREQ", metavar="FREQ")
default=False,
help="use intermediate freq (compensates DC
problems in quadrature boards)" )
parser.add_option("-g", "--gain", type="eng_float", default=None,
help="set gain in dB (default is maximum)")
parser.add_option("-V", "--volume", type="eng_float", default=None,
help="set volume (default is midpoint)")
parser.add_option("-O", "--audio-output", type="string", default="",
help="pcm device name.  E.g., hw:0,0 or surround51
or /dev/dsp")

(options, args) = parser.parse_args()
if len(args) != 0:
parser.print_help()
sys.exit(1)

self.frame = frame
self.panel = panel
self.use_IF=options.use_if_freq
if self.use_IF:
self.IF_freq=64000.0
else:
self.IF_freq=0.0

self.vol = 0
self.state = "FREQ"
self.freq = 0

# build graph

#TODO: add an AGC after the channel filter and before the AM_demod

self.u = usrp.source_c()                    # usrp is data source

usrp_decim = 250
self.u.set_decim_rate(usrp_decim)
usrp_rate = adc_rate / usrp_decim           # 256 kS/s
chanfilt_decim = 8
demod_rate = usrp_rate / chanfilt_decim     # 32 kHz
audio_decimation = 128
audio_rate = demod_rate / audio_decimation  # 500 Hz
# This is the range of
data which holds the carrier AM 30 Hz signal
fm_decimation = 1
self.fm_rate = demod_rate / fm_decimation
fm_audio_decimation = 128
self.fm_audio_rate = self.fm_rate / fm_audio_decimation# 500 Hz

if options.rx_subdev_spec is None:
options.rx_subdev_spec = pick_subdevice(self.u)

self.u.set_mux(usrp.determine_rx_mux_value(self.u,
options.rx_subdev_spec))
self.subdev = usrp.selected_subdev(self.u, options.rx_subdev_spec)
print "Using RX d'board %s" % (self.subdev.side_and_name(),)

# First bring in all the data we need to use, which includes the 30
Hz signal (AM modulated),
# the 1020 Hz morse code identifier
# and the 9960 Hz subcarrier which has a 480Hz width, so we'll cut
off at 11000 Hz
chan_filt_coeffs = optfir.low_pass (1,           # gain
usrp_rate,   # sampling rate
11e3,        # passband cutoff
12e3,       # stopband cutoff
1.0,         # passband ripple
60)          # stopband
attenuation
if self.use_IF:
# Turn If to baseband and filter.
self.chan_filt = gr.freq_xlating_fir_filter_ccf (chanfilt_decim,
chan_filt_coeffs, self.IF_freq, usrp_rate)
else:
self.chan_filt = gr.fir_filter_ccf (chanfilt_decim,
chan_filt_coeffs)
self.am_demod = gr.complex_to_mag()

# Now we split off the 30 Hz signal after it's been demodulated
audio_filt_coeffs = gr.firdes.band_pass (.5,           # gain
demod_rate,   # sampling rate
20,        # low cutoff
50,        # high cutoff
150)       # stopband cutoff

self.audio_filt=gr.fir_filter_fff(audio_decimation,audio_filt_coeffs)

# convert it back to complex with a Hilbert filter
am_hilbert_coeffs = gr.firdes.hilbert (27)
am_hilbert_filt = gr.fir_filter_fff(1, am_hilbert_coeffs)

# that detected the main AM 30 Hz signal, now we need to decode the
second signal
# It's on a 9960Hz subcarrier, FM modulated

second_chan_coeffs = gr.firdes.low_pass (1.0,                # gain
demod_rate,       #
sampling rate
1000,               # low
pass cutoff freq
1000,                #
width of trans. band
gr.firdes.WIN_HANN) #
filter type

self.fm_filt = gr.freq_xlating_fir_filter_fcc(fm_decimation,
# decimation rate
second_chan_coeffs,
# taps
-9960,
# frequency translation amount
demod_rate)         #
input sample rate

# Once the signal is brought down to baseband, then we can
demodulate it

# convert it back to complex with a Hilbert filter
fm_hilbert_coeffs = gr.firdes.hilbert (27)
fm_hilbert_filt = gr.fir_filter_fff(1, fm_hilbert_coeffs)

# and once demodulated, then we can filter out everything but the
30Hz signal
fm_audio_filt_coeffs = optfir.low_pass (10,           # gain
demod_rate / fm_decimation,   #
sampling rate
50,        # passband cutoff
150,       # stopband cutoff
0.1,         # passband ripple
60)          # stopband
attenuation

self.fm_audio_filt=gr.fir_filter_fff(fm_audio_decimation,fm_audio_filt_coeff
s)

# Now there are 2 30 Hz signal, one out of fm_audio_filt, another
out of audio_filt.
# The directional signal is found by comparing their phase

self.fm_conv = gr.float_to_complex()
self.am_conv = gr.float_to_complex()
conj = gr.conjugate_cc()
self.mult = gr.multiply_cc()
self.arg = gr.complex_to_arg()
self.fm_arg = gr.complex_to_arg()
#self.sub = gr.sub_ff()

# For the moment, no sound is output.  However, we really should
also take the chan_filt output,
# and run it through a different filter/AM demod to get the 1020 Hz
ID signal, and send that to the sound card
# The code has been chopped for the moment, but we'll keep part of
it here for future use

self.volume_control = gr.multiply_const_ff(self.vol)

# sound card as the final sink
audio_sink = audio.sink (int (audio_rate),
options.audio_output,
False)  # ok_to_block

# now wire it all together
# obtain our 2 real signals
self.connect (self.u, self.chan_filt, self.am_demod,
self.audio_filt) #, (self.volume_control, 0), (audio_sink,0))
self.connect (self.am_demod, self.fm_filt, self.fm_rx,
self.fm_audio_filt) #, (audio_sink,1))

# convert the real signals back to complex w/ a hilbert filter
self.connect(self.audio_filt, (self.am_conv,0))
self.connect(self.audio_filt, am_hilbert_filt, (self.am_conv,1))
self.connect(self.fm_audio_filt, (self.fm_conv,0))
self.connect(self.fm_audio_filt, fm_hilbert_filt, (self.fm_conv,1))

# now run them through the phase comparator
self.connect(self.am_conv, conj, self.mult)
self.connect(self.fm_conv,  (self.mult, 1))
self.connect(self.mult, self.arg)

# alternative phase comparator (comment out the mult/conj
definitions above to use)
#self.connect(self.am_conv, self.arg, (self.sub,0))
#self.connect(self.fm_conv, self.fm_arg, (self.sub,1))

self._build_gui(vbox, usrp_rate, demod_rate, audio_rate)

if options.gain is None:
g = self.subdev.gain_range()
if True:
# if no gain was specified, use the maximum gain available
# (usefull for Basic_RX which is relatively deaf and the most
probable board to be used for AM)
# TODO: check db type to decide on default gain.
#options.gain = float(g)
options.gain = float(g+g)/2 #Changed with to 50% because
the TV Tuner board doesn't need the full gain and that's what I'm using.
else:
# if no gain was specified, use the mid-point in dB
options.gain = float(g+g)/2

if options.volume is None:
g = self.volume_range()
options.volume = float(g*3+g)/4

if abs(options.freq) < 1e3:
options.freq *= 1e3

# set initial values

self.set_gain(options.gain)
self.set_vol(options.volume)
if not(self.set_freq(options.freq)):
self._set_status_msg("Failed to set initial frequency")

def _set_status_msg(self, msg, which=0):
self.frame.GetStatusBar().SetStatusText(msg, which)

def _build_gui(self, vbox, usrp_rate, demod_rate, audio_rate):

def _form_set_freq(kv):
return self.set_freq(kv['freq'])

if 0:
self.src_fft = fftsink2.fft_sink_c(self.panel, title="Data from
USRP",
fft_size=512,
sample_rate=usrp_rate,
ref_scale=32768.0,
ref_level=0.0, y_divs=12)
self.connect (self.u, self.src_fft)
vbox.Add (self.src_fft.win, 4, wx.EXPAND)

if 0:
self.post_filt_fft = fftsink2.fft_sink_c(self.panel, title="Post
Channel filter",
fft_size=512,
sample_rate=demod_rate)
self.connect (self.chan_filt, self.post_filt_fft)
vbox.Add (self.post_filt_fft.win, 4, wx.EXPAND)

if 0:
pre_demod_fft = fftsink2.fft_sink_f(self.panel, title="Alt
Signal pre-FM Demod",
fft_size=512,
sample_rate=demod_rate,
y_per_div=10, ref_level=40)
self.connect (self.am_demod, pre_demod_fft)
vbox.Add (pre_demod_fft.win, 4, wx.EXPAND)
if 0:
post_demod_fft = fftsink2.fft_sink_f(self.panel, title="Alt
Signal",
fft_size=512,
sample_rate=self.fm_audio_rate,
y_per_div=10, ref_level=40)
self.connect (self.fm_audio_filt, post_demod_fft)
vbox.Add (post_demod_fft.win, 4, wx.EXPAND)
if 1:
oscope = scopesink2.scope_sink_c(self.panel,
sample_rate=audio_rate,
frame_decim=1,
v_scale=2,
t_scale=.05,
num_inputs=2)
#self.connect (self.mult, (oscope,1))
self.connect (self.fm_conv, (oscope,1))
self.connect (self.am_conv, (oscope,0))
vbox.Add (oscope.win, 4, wx.EXPAND)
if 0:
audio_fft = fftsink2.fft_sink_f(self.panel, title="Carrier AM
signal",
fft_size=512,
sample_rate=audio_rate,
y_per_div=10,
ref_level=40)
self.connect (self.audio_filt, audio_fft)
vbox.Add (audio_fft.win, 4, wx.EXPAND)
if 1:
num_sink = numbersink2.number_sink_f (self.panel,
avg_alpha=1.0e-5,average=False,
sample_rate=audio_rate,
factor=1.0,base_value=0,
minval=-7, maxval=7,
ref_level=0,
decimal_places=5,number_rate=15)
self.connect(self.arg, num_sink)
vbox.Add (num_sink.win, 1, wx.EXPAND)

# control area form at bottom
self.myform = myform = form.form()

hbox = wx.BoxSizer(wx.HORIZONTAL)
myform['freq'] = form.float_field(
parent=self.panel, sizer=hbox, label="Freq", weight=1,
callback=myform.check_input_and_call(_form_set_freq,
self._set_status_msg))

myform['freq_slider'] = \
form.quantized_slider_field(parent=self.panel, sizer=hbox,
weight=3,
range=(108.0e6, 118.0e6, 1.0e5),
callback=self.set_freq)

hbox = wx.BoxSizer(wx.HORIZONTAL)

myform['volume'] = \
form.quantized_slider_field(parent=self.panel, sizer=hbox,
label="Volume",
weight=3, range=self.volume_range(),
callback=self.set_vol)

myform['gain'] = \
form.quantized_slider_field(parent=self.panel, sizer=hbox,
label="Gain",
weight=3,
range=self.subdev.gain_range(),
callback=self.set_gain)

try:
self.knob = powermate.powermate(self.frame)
self.rot = 0
powermate.EVT_POWERMATE_ROTATE (self.frame, self.on_rotate)
powermate.EVT_POWERMATE_BUTTON (self.frame, self.on_button)
except:
print "FYI: No Powermate or Contour Knob found"

def on_rotate (self, event):
self.rot += event.delta
if (self.state == "FREQ"):
if self.rot >= 3:
self.set_freq(self.freq + .1e6)
self.rot -= 3
elif self.rot <=-3:
self.set_freq(self.freq - .1e6)
self.rot += 3
else:
step = self.volume_range()
if self.rot >= 3:
self.set_vol(self.vol + step)
self.rot -= 3
elif self.rot <=-3:
self.set_vol(self.vol - step)
self.rot += 3

def on_button (self, event):
if event.value == 0:        # button up
return
self.rot = 0
if self.state == "FREQ":
self.state = "VOL"
else:
self.state = "FREQ"
self.update_status_bar ()

def set_vol (self, vol):
g = self.volume_range()
self.vol = max(g, min(g, vol))
self.volume_control.set_k(10**(self.vol/10))
self.myform['volume'].set_value(self.vol)
self.update_status_bar ()

def set_freq(self, target_freq):
"""
Set the center frequency we're interested in.

@param target_freq: frequency in Hz
@rypte: bool

Tuning is a two step process.  First we ask the front-end to
tune as close to the desired frequency as it can.  Then we use
the result of that operation and our target_frequency to
determine the value for the digital down converter.
"""
r = usrp.tune(self.u, 0, self.subdev, target_freq  + self.IF_freq)
#TODO: check if db is inverting the spectrum or not to decide if we
should do + self.IF_freq  or - self.IF_freq

if r:
self.freq = target_freq
self.myform['freq'].set_value(target_freq)         # update
displayed value
self.myform['freq_slider'].set_value(target_freq)  # update
displayed value
self.update_status_bar()
self._set_status_msg("OK", 0)
return True

self._set_status_msg("Failed", 0)
return False

def set_gain(self, gain):
self.myform['gain'].set_value(gain)     # update displayed value
self.subdev.set_gain(gain)

def update_status_bar (self):
msg = "Volume:%r  Setting:%s" % (self.vol, self.state)
self._set_status_msg(msg, 1)
try:
self.src_fft.set_baseband_freq(self.freq)
except:
None

def volume_range(self):
return (-40.0, 0.0, 0.5)

if __name__ == '__main__':
app = stdgui2.stdapp (vor_rx_block, "USRP VOR RX")
app.MainLoop ()

reply via email to