don't use global configuration

This commit is contained in:
Roman Zeyde
2014-11-09 12:14:54 +02:00
parent c8f5924c12
commit ceb826728a
15 changed files with 326 additions and 300 deletions

View File

@@ -19,7 +19,7 @@ log = logging.getLogger('__name__')
from amodem import config from amodem import config
from amodem import recv from amodem import recv
from amodem import send from amodem import send
from amodem import wave from amodem import audio
from amodem import calib from amodem import calib
null = open('/dev/null', 'wb') null = open('/dev/null', 'wb')
@@ -32,10 +32,11 @@ def FileType(mode, process=None):
fname = '-' fname = '-'
if fname is None: if fname is None:
assert process is not None
if 'r' in mode: if 'r' in mode:
return process(stdout=wave.sp.PIPE, stderr=null).stdout return process.launch(stdout=audio.sp.PIPE, stderr=null).stdout
if 'w' in mode: if 'w' in mode:
return process(stdin=wave.sp.PIPE, stderr=null).stdin return process.launch(stdin=audio.sp.PIPE, stderr=null).stdin
if fname == '-': if fname == '-':
if 'r' in mode: if 'r' in mode:
@@ -81,7 +82,7 @@ def main():
sender.set_defaults( sender.set_defaults(
main=run_send, main=run_send,
input_type=FileType('rb'), input_type=FileType('rb'),
output_type=FileType('wb', wave.play) output_type=FileType('wb', audio.play(Fs=config.Fs))
) )
# Demodulator # Demodulator
@@ -104,7 +105,7 @@ def main():
help='plot results using pylab module') help='plot results using pylab module')
receiver.set_defaults( receiver.set_defaults(
main=run_recv, main=run_recv,
input_type=FileType('rb', wave.record), input_type=FileType('rb', audio.record(Fs=config.Fs)),
output_type=FileType('wb') output_type=FileType('wb')
) )
@@ -127,6 +128,7 @@ def main():
if getattr(args, 'plot', False): if getattr(args, 'plot', False):
import pylab import pylab
args.plot = pylab args.plot = pylab
args.config = config
args.main(args) args.main(args)
@@ -148,18 +150,18 @@ def run_modem(args, func):
def run_send(args): def run_send(args):
if args.calibrate: if args.calibrate:
calib.send(verbose=args.verbose) calib.send(config=config, verbose=args.verbose)
elif args.wave: elif args.wave:
join_process(wave.play(fname=args.input)) join_process(audio.play(Fs=config.Fs).launch(fname=args.input))
else: else:
run_modem(args, send.main) run_modem(args, send.main)
def run_recv(args): def run_recv(args):
if args.calibrate: if args.calibrate:
calib.recv(verbose=args.verbose) calib.recv(config=config, verbose=args.verbose)
elif args.wave: elif args.wave:
join_process(wave.record(fname=args.output)) join_process(audio.record(Fs=config.Fs).launch(fname=args.output))
else: else:
run_modem(args, recv.main) run_modem(args, recv.main)

30
amodem/audio.py Normal file
View File

@@ -0,0 +1,30 @@
import subprocess as sp
import logging
import functools
log = logging.getLogger(__name__)
class ALSA(object):
def __init__(self, tool, Fs):
self.Fs = int(Fs) # sampling rate
self.bits_per_sample = 16
self.bytes_per_sample = self.bits_per_sample / 8.0
self.bytes_per_second = self.bytes_per_sample * self.Fs
# PCM signed little endian
self.audio_format = 'S{}_LE'.format(self.bits_per_sample)
self.audio_tool = tool
def launch(self, fname=None, **kwargs):
if fname is None:
fname = '-' # use stdin/stdout if filename not specified
args = [self.audio_tool, fname, '-q',
'-f', self.audio_format,
'-c', '1',
'-r', str(self.Fs)]
log.debug('Running: %r', ' '.join(args))
p = sp.Popen(args=args, **kwargs)
return p
# Use ALSA tools for audio playing/recording
play = functools.partial(ALSA, tool='aplay')
record = functools.partial(ALSA, tool='arecord')

View File

@@ -3,20 +3,21 @@ import logging
log = logging.getLogger(__name__) log = logging.getLogger(__name__)
from subprocess import PIPE
from . import common from . import common
from . import config from . import audio
from . import wave
CALIBRATION_SYMBOLS = int(1.0 * config.Fs)
ALLOWED_EXCEPTIONS = (IOError, KeyboardInterrupt) ALLOWED_EXCEPTIONS = (IOError, KeyboardInterrupt)
def send(wave_play=wave.play, verbose=False): def send(config, audio_play=audio.play, verbose=False):
t = np.arange(0, CALIBRATION_SYMBOLS) * config.Ts calibration_symbols = int(1.0 * config.Fs)
t = np.arange(0, calibration_symbols) * config.Ts
signal = [np.sin(2 * np.pi * f * t) for f in config.frequencies] signal = [np.sin(2 * np.pi * f * t) for f in config.frequencies]
signal = common.dumps(np.concatenate(signal)) signal = common.dumps(np.concatenate(signal))
p = wave_play(stdin=wave.sp.PIPE) player = audio_play(Fs=config.Fs)
p = player.launch(stdin=PIPE)
fd = p.stdin fd = p.stdin
try: try:
while True: while True:
@@ -27,16 +28,17 @@ def send(wave_play=wave.play, verbose=False):
p.kill() p.kill()
FRAME_LENGTH = 200 * config.Nsym def run_recorder(config, recorder):
FRAME_LENGTH = 200 * config.Nsym
process = recorder.launch(stdout=PIPE)
frame_size = int(recorder.bytes_per_sample * FRAME_LENGTH)
def recorder(process):
t = np.arange(0, FRAME_LENGTH) * config.Ts t = np.arange(0, FRAME_LENGTH) * config.Ts
scaling_factor = 0.5 * len(t) scaling_factor = 0.5 * len(t)
carriers = [np.exp(2j * np.pi * f * t) for f in config.frequencies] carriers = [np.exp(2j * np.pi * f * t) for f in config.frequencies]
carriers = np.array(carriers) / scaling_factor carriers = np.array(carriers) / scaling_factor
frame_size = int(wave.bytes_per_sample * FRAME_LENGTH)
fd = process.stdout fd = process.stdout
states = [True] states = [True]
@@ -80,13 +82,15 @@ def recorder(process):
fmt = '{freq:6.0f} Hz: {message:s}{extra:s}' fmt = '{freq:6.0f} Hz: {message:s}{extra:s}'
fields = ['peak', 'total', 'rms', 'coherency'] fields = ['peak', 'total', 'rms', 'coherency']
def recv(wave_record=wave.record, verbose=False, output=None): def recv(config, audio_record=audio.record, verbose=False):
extra = '' extra = ''
if verbose: if verbose:
extra = ''.join(', {0}={{{0}:.4f}}'.format(f) for f in fields) extra = ''.join(', {0}={{{0}:.4f}}'.format(f) for f in fields)
for r in recorder(wave_record(stdout=wave.sp.PIPE)):
msg = fmt.format(extra=extra.format(**r), **r) recorder = audio_record(Fs=config.Fs)
if not r['error']: for result in run_recorder(config=config, recorder=recorder):
msg = fmt.format(extra=extra.format(**result), **result)
if not result['error']:
log.info(msg) log.info(msg)
else: else:
log.error(msg) log.error(msg)

View File

@@ -28,11 +28,9 @@ def loads(data):
return x return x
def dumps(sym, n=1): def dumps(sym):
sym = sym.real * scaling sym = sym.real * scaling
sym = sym.astype('int16') return sym.astype('int16').tostring()
data = sym.tostring()
return data * n
def iterate(data, size, func=None, truncate=True, enumerate=False): def iterate(data, size, func=None, truncate=True, enumerate=False):

View File

@@ -4,7 +4,6 @@ import logging
log = logging.getLogger(__name__) log = logging.getLogger(__name__)
from . import config
from . import common from . import common
@@ -64,9 +63,9 @@ def estimate(x, y, order, lookahead=0):
class Demux(object): class Demux(object):
def __init__(self, sampler, freqs): def __init__(self, sampler, omegas, Nsym):
Nsym = config.Nsym self.Nsym = Nsym
self.filters = [exp_iwt(-f, Nsym) / (0.5*Nsym) for f in freqs] self.filters = [exp_iwt(-w, Nsym) / (0.5*self.Nsym) for w in omegas]
self.filters = np.array(self.filters) self.filters = np.array(self.filters)
self.sampler = sampler self.sampler = sampler
@@ -74,8 +73,8 @@ class Demux(object):
return self return self
def next(self): def next(self):
frame = self.sampler.take(size=config.Nsym) frame = self.sampler.take(size=self.Nsym)
if len(frame) == config.Nsym: if len(frame) == self.Nsym:
return np.dot(self.filters, frame) return np.dot(self.filters, frame)
else: else:
raise StopIteration raise StopIteration
@@ -83,18 +82,17 @@ class Demux(object):
__next__ = next __next__ = next
def exp_iwt(freq, n): def exp_iwt(omega, n):
iwt = 2j * np.pi * freq * np.arange(n) * config.Ts return np.exp(1j * omega * np.arange(n))
return np.exp(iwt)
def norm(x): def norm(x):
return np.sqrt(np.dot(x.conj(), x).real) return np.sqrt(np.dot(x.conj(), x).real)
def coherence(x, freq): def coherence(x, omega):
n = len(x) n = len(x)
Hc = exp_iwt(-freq, n) / np.sqrt(0.5*n) Hc = exp_iwt(-omega, n) / np.sqrt(0.5*n)
norm_x = norm(x) norm_x = norm(x)
if norm_x: if norm_x:
return np.dot(Hc, x) / norm_x return np.dot(Hc, x) / norm_x
@@ -151,9 +149,3 @@ class MODEM(object):
if error_handler: if error_handler:
error_handler(received=received, decoded=decoded) error_handler(received=received, decoded=decoded)
yield bits yield bits
def __repr__(self):
return '<{:.3f} kbps, {:d}-QAM, {:d} carriers>'.format(
config.modem_bps / 1e3, len(self.symbols), len(config.carriers))
__str__ = __repr__

View File

@@ -2,7 +2,6 @@ import numpy as np
from numpy.linalg import lstsq from numpy.linalg import lstsq
from amodem import dsp from amodem import dsp
from amodem import config
from amodem import sampling from amodem import sampling
import itertools import itertools
@@ -10,76 +9,76 @@ import random
_constellation = [1, 1j, -1, -1j] _constellation = [1, 1j, -1, -1j]
class Equalizer(object):
def train_symbols(length, seed=0, Nfreq=config.Nfreq): def __init__(self, config):
r = random.Random(seed) self.carriers = config.carriers
choose = lambda: [r.choice(_constellation) for j in range(Nfreq)] self.omegas = 2 * np.pi * np.array(config.frequencies) / config.Fs
return np.array([choose() for i in range(length)]) self.Nfreq = config.Nfreq
self.Nsym = config.Nsym
def train_symbols(self, length, seed=0):
r = random.Random(seed)
choose = lambda: [r.choice(_constellation) for j in range(self.Nfreq)]
return np.array([choose() for i in range(length)])
def modulator(symbols): def modulator(self, symbols):
carriers = config.carriers gain = 1.0 / len(self.carriers)
gain = 1.0 / len(carriers) result = []
result = [] for s in symbols:
for s in symbols: result.append(np.dot(s, self.carriers))
result.append(np.dot(s, carriers)) result = np.concatenate(result).real * gain
result = np.concatenate(result).real * gain assert np.max(np.abs(result)) <= 1
assert np.max(np.abs(result)) <= 1 return result
return result
def demodulator(self, signal, size):
signal = itertools.chain(signal, itertools.repeat(0))
symbols = dsp.Demux(sampler=sampling.Sampler(signal),
omegas=self.omegas, Nsym=self.Nsym)
return np.array(list(itertools.islice(symbols, size)))
def demodulator(signal, size): def equalize_symbols(self, signal, symbols, order, lookahead=0):
signal = itertools.chain(signal, itertools.repeat(0)) assert symbols.shape[1] == self.Nfreq
symbols = dsp.Demux(sampling.Sampler(signal), config.frequencies) length = symbols.shape[0]
return np.array(list(itertools.islice(symbols, size)))
matched = np.array(self.carriers) / (0.5*self.Nsym)
matched = matched[:, ::-1].transpose().conj()
signal = np.concatenate([signal, np.zeros(lookahead)])
y = dsp.lfilter(x=signal, b=matched, a=[1])
def equalize_symbols(signal, symbols, order, lookahead=0): A = []
Nsym = config.Nsym b = []
Nfreq = config.Nfreq
carriers = config.carriers
assert symbols.shape[1] == Nfreq for j in range(self.Nfreq):
length = symbols.shape[0] for i in range(length):
offset = (i+1)*self.Nsym
row = y[offset-order:offset+lookahead, j]
A.append(row)
b.append(symbols[i, j])
matched = np.array(carriers) / (0.5*Nsym) A = np.array(A)
matched = matched[:, ::-1].transpose().conj() b = np.array(b)
signal = np.concatenate([signal, np.zeros(lookahead)]) h, residuals, rank, sv = lstsq(A, b)
y = dsp.lfilter(x=signal, b=matched, a=[1]) h = h[::-1].real
A = [] return h
b = []
for j in range(Nfreq): def equalize_signal(self, signal, expected, order, lookahead=0):
for i in range(length): signal = np.concatenate([np.zeros(order-1), signal, np.zeros(lookahead)])
offset = (i+1)*Nsym length = len(expected)
row = y[offset-order:offset+lookahead, j]
A.append(row)
b.append(symbols[i, j])
A = np.array(A) A = []
b = np.array(b) b = []
h, residuals, rank, sv = lstsq(A, b)
h = h[::-1].real
return h for i in range(length - order):
offset = order + i
row = signal[offset-order:offset+lookahead]
A.append(np.array(row, ndmin=2))
b.append(expected[i])
A = np.concatenate(A, axis=0)
def equalize_signal(signal, expected, order, lookahead=0): b = np.array(b)
signal = np.concatenate([np.zeros(order-1), signal, np.zeros(lookahead)]) h, residuals, rank, sv = lstsq(A, b)
length = len(expected) h = h[::-1].real
return h
A = []
b = []
for i in range(length - order):
offset = order + i
row = signal[offset-order:offset+lookahead]
A.append(np.array(row, ndmin=2))
b.append(expected[i])
A = np.concatenate(A, axis=0)
b = np.array(b)
h, residuals, rank, sv = lstsq(A, b)
h = h[::-1].real
return h

View File

@@ -12,87 +12,97 @@ from . import dsp
from . import sampling from . import sampling
from . import train from . import train
from . import common from . import common
from . import config
from . import framing from . import framing
from . import equalizer from . import equalizer
modem = dsp.MODEM(config.symbols) class Detector(object):
# Plots' size (WIDTH x HEIGHT) COHERENCE_THRESHOLD = 0.99
HEIGHT = np.floor(np.sqrt(config.Nfreq))
WIDTH = np.ceil(config.Nfreq / float(HEIGHT))
COHERENCE_THRESHOLD = 0.99 CARRIER_DURATION = sum(train.prefix)
CARRIER_THRESHOLD = int(0.99 * CARRIER_DURATION)
SEARCH_WINDOW = 10 # symbols
CARRIER_DURATION = sum(train.prefix) def __init__(self, config):
CARRIER_THRESHOLD = int(0.99 * CARRIER_DURATION) self.freq = config.Fc
SEARCH_WINDOW = 10 # symbols self.omega = 2 * np.pi * self.freq / config.Fs
self.Nsym = config.Nsym
self.Tsym = config.Tsym
self.maxlen = config.baud # 1 second of symbols
def run(self, samples):
counter = 0
bufs = collections.deque([], maxlen=self.maxlen)
for offset, buf in common.iterate(samples, self.Nsym, enumerate=True):
bufs.append(buf)
def report_carrier(bufs, begin): coeff = dsp.coherence(buf, self.omega)
x = np.concatenate(tuple(bufs)[-CARRIER_THRESHOLD:-1]) if abs(coeff) > self.COHERENCE_THRESHOLD:
Hc = dsp.exp_iwt(-config.Fc, len(x)) counter += 1
Zc = np.dot(Hc, x) / (0.5*len(x)) else:
amp = abs(Zc) counter = 0
log.info('Carrier detected at ~%.1f ms @ %.1f kHz:'
' coherence=%.3f%%, amplitude=%.3f',
begin * config.Tsym * 1e3 / config.Nsym, config.Fc / 1e3,
np.abs(dsp.coherence(x, config.Fc)) * 100, amp)
return amp
if counter == self.CARRIER_THRESHOLD:
def detect(samples, freq): length = (self.CARRIER_THRESHOLD - 1) * self.Nsym
counter = 0 begin = offset - length
bufs = collections.deque([], maxlen=config.baud) # 1 second of symbols amplitude = self.report_carrier(bufs, begin=begin)
for offset, buf in common.iterate(samples, config.Nsym, enumerate=True): break
bufs.append(buf)
coeff = dsp.coherence(buf, config.Fc)
if abs(coeff) > COHERENCE_THRESHOLD:
counter += 1
else: else:
counter = 0 raise ValueError('No carrier detected')
if counter == CARRIER_THRESHOLD: log.debug('Buffered %d ms of audio', len(bufs))
length = (CARRIER_THRESHOLD - 1) * config.Nsym
begin = offset - length
amplitude = report_carrier(bufs, begin=begin)
break
else:
raise ValueError('No carrier detected')
log.debug('Buffered %d ms of audio', len(bufs)) bufs = list(bufs)[-self.CARRIER_THRESHOLD-self.SEARCH_WINDOW:]
trailing = list(itertools.islice(samples, self.SEARCH_WINDOW*self.Nsym))
bufs.append(np.array(trailing))
bufs = list(bufs)[-CARRIER_THRESHOLD-SEARCH_WINDOW:] buf = np.concatenate(bufs)
trailing = list(itertools.islice(samples, SEARCH_WINDOW*config.Nsym)) offset = self.find_start(buf, self.CARRIER_DURATION*self.Nsym)
bufs.append(np.array(trailing)) log.debug('Carrier starts at %.3f ms',
offset * self.Tsym * 1e3 / self.Nsym)
buf = np.concatenate(bufs) return itertools.chain(buf[offset:], samples), amplitude
offset = find_start(buf, CARRIER_DURATION*config.Nsym)
log.debug('Carrier starts at %.3f ms',
offset * config.Tsym * 1e3 / config.Nsym)
return itertools.chain(buf[offset:], samples), amplitude
def find_start(buf, length): def report_carrier(self, bufs, begin):
N = len(buf) x = np.concatenate(tuple(bufs)[-self.CARRIER_THRESHOLD:-1])
carrier = dsp.exp_iwt(config.Fc, N) Hc = dsp.exp_iwt(-self.omega, len(x))
z = np.cumsum(buf * carrier) Zc = np.dot(Hc, x) / (0.5*len(x))
z = np.concatenate([[0], z]) amp = abs(Zc)
correlations = np.abs(z[length:] - z[:-length]) log.info('Carrier detected at ~%.1f ms @ %.1f kHz:'
return np.argmax(correlations) ' coherence=%.3f%%, amplitude=%.3f',
begin * self.Tsym * 1e3 / self.Nsym, self.freq / 1e3,
np.abs(dsp.coherence(x, self.omega)) * 100, amp)
return amp
def find_start(self, buf, length):
N = len(buf)
carrier = dsp.exp_iwt(self.omega, N)
z = np.cumsum(buf * carrier)
z = np.concatenate([[0], z])
correlations = np.abs(z[length:] - z[:-length])
return np.argmax(correlations)
class Receiver(object): class Receiver(object):
def __init__(self, plt=None): def __init__(self, config, plt=None):
self.stats = {} self.stats = {}
self.plt = plt or common.Dummy() self.plt = plt or common.Dummy()
self.modem = dsp.MODEM(config.symbols)
self.frequencies = np.array(config.frequencies)
self.omegas = 2 * np.pi * self.frequencies / config.Fs
self.Nsym = config.Nsym
self.Tsym = config.Tsym
self.iters_per_report = config.baud # report once per second
self.modem_bitrate = config.modem_bps
self.equalizer = equalizer.Equalizer(config)
self.carrier_index = config.carrier_index
def _prefix(self, sampler, freq, gain=1.0, skip=5): def _prefix(self, symbols, gain=1.0, skip=5):
symbols = dsp.Demux(sampler, [freq]) S = common.take(symbols, len(train.prefix))
S = common.take(symbols, len(train.prefix)).squeeze() * gain S = S[:, self.carrier_index] * gain
sliced = np.round(np.abs(S)) sliced = np.round(np.abs(S))
self.plt.figure() self.plt.figure()
self.plt.subplot(121) self.plt.subplot(121)
@@ -116,7 +126,7 @@ class Receiver(object):
self.plt.plot(indices, phase, ':') self.plt.plot(indices, phase, ':')
self.plt.plot(indices, a * indices + b) self.plt.plot(indices, a * indices + b)
freq_err = a / (config.Tsym * config.Fc) freq_err = a / (self.Tsym * self.frequencies[self.carrier_index])
last_phase = a * indices[-1] + b last_phase = a * indices[-1] + b
log.debug('Current phase on carrier: %.3f', last_phase) log.debug('Current phase on carrier: %.3f', last_phase)
@@ -125,16 +135,16 @@ class Receiver(object):
return freq_err return freq_err
def _train(self, sampler, order, lookahead): def _train(self, sampler, order, lookahead):
gain = config.Nfreq Nfreq = len(self.frequencies)
train_symbols = equalizer.train_symbols(train.equalizer_length) train_symbols = self.equalizer.train_symbols(train.equalizer_length)
train_signal = equalizer.modulator(train_symbols) * gain train_signal = self.equalizer.modulator(train_symbols) * Nfreq
prefix = postfix = train.silence_length * config.Nsym prefix = postfix = train.silence_length * self.Nsym
signal_length = train.equalizer_length * config.Nsym + prefix + postfix signal_length = train.equalizer_length * self.Nsym + prefix + postfix
signal = sampler.take(signal_length + lookahead) signal = sampler.take(signal_length + lookahead)
coeffs = equalizer.equalize_signal( coeffs = self.equalizer.equalize_signal(
signal=signal[prefix:-postfix], signal=signal[prefix:-postfix],
expected=train_signal, expected=train_signal,
order=order, lookahead=lookahead order=order, lookahead=lookahead
@@ -147,7 +157,7 @@ class Receiver(object):
equalized = list(equalization_filter(signal)) equalized = list(equalization_filter(signal))
equalized = equalized[prefix+lookahead:-postfix+lookahead] equalized = equalized[prefix+lookahead:-postfix+lookahead]
symbols = equalizer.demodulator(equalized, train.equalizer_length) symbols = self.equalizer.demodulator(equalized, train.equalizer_length)
sliced = np.array(symbols).round() sliced = np.array(symbols).round()
errors = np.array(sliced - train_symbols, dtype=np.bool) errors = np.array(sliced - train_symbols, dtype=np.bool)
error_rate = errors.sum() / errors.size error_rate = errors.sum() / errors.size
@@ -160,17 +170,15 @@ class Receiver(object):
SNRs = 20.0 * np.log10(signal_rms / noise_rms) SNRs = 20.0 * np.log10(signal_rms / noise_rms)
self.plt.figure() self.plt.figure()
for i, freq, snr in zip(range(config.Nfreq), config.frequencies, SNRs): for (i, freq), snr in zip(enumerate(self.frequencies), SNRs):
log.debug('%5.1f kHz: SNR = %5.2f dB', freq / 1e3, snr) log.debug('%5.1f kHz: SNR = %5.2f dB', freq / 1e3, snr)
self.plt.subplot(HEIGHT, WIDTH, i+1)
self._constellation(symbols[:, i], train_symbols[:, i], self._constellation(symbols[:, i], train_symbols[:, i],
'$F_c = {} Hz$'.format(freq)) '$F_c = {} Hz$'.format(freq), index=i)
assert error_rate == 0, error_rate assert error_rate == 0, error_rate
return equalization_filter return equalization_filter
def _demodulate(self, sampler, freqs): def _demodulate(self, sampler, symbols):
streams = [] streams = []
symbol_list = [] symbol_list = []
errors = {} errors = {}
@@ -178,52 +186,51 @@ class Receiver(object):
def error_handler(received, decoded, freq): def error_handler(received, decoded, freq):
errors.setdefault(freq, []).append(received / decoded) errors.setdefault(freq, []).append(received / decoded)
symbols = dsp.Demux(sampler, freqs) generators = common.split(symbols, n=len(self.omegas))
generators = common.split(symbols, n=len(freqs)) for freq, S in zip(self.frequencies, generators):
for freq, S in zip(freqs, generators):
equalized = [] equalized = []
S = common.icapture(S, result=equalized) S = common.icapture(S, result=equalized)
symbol_list.append(equalized) symbol_list.append(equalized)
freq_handler = functools.partial(error_handler, freq=freq) freq_handler = functools.partial(error_handler, freq=freq)
bits = modem.decode(S, freq_handler) # list of bit tuples bits = self.modem.decode(S, freq_handler) # list of bit tuples
streams.append(bits) # stream per frequency streams.append(bits) # stream per frequency
self.stats['symbol_list'] = symbol_list self.stats['symbol_list'] = symbol_list
self.stats['rx_bits'] = 0 self.stats['rx_bits'] = 0
self.stats['rx_start'] = time.time() self.stats['rx_start'] = time.time()
log.info('Starting demodulation: %s', modem) log.info('Starting demodulation: %s', self.modem)
for i, block in enumerate(common.izip(streams)): # block per frequency for i, block in enumerate(common.izip(streams), 1):
for bits in block: for bits in block:
self.stats['rx_bits'] = self.stats['rx_bits'] + len(bits) self.stats['rx_bits'] = self.stats['rx_bits'] + len(bits)
yield bits yield bits
if i > 0 and i % config.baud == 0: if i % self.iters_per_report == 0:
err = np.array([e for v in errors.values() for e in v]) err = np.array([e for v in errors.values() for e in v])
err = np.mean(np.angle(err))/(2*np.pi) if len(err) else 0 err = np.mean(np.angle(err))/(2*np.pi) if len(err) else 0
errors.clear() errors.clear()
duration = time.time() - self.stats['rx_start'] duration = time.time() - self.stats['rx_start']
sampler.freq -= 0.01 * err / config.Fc sampler.freq -= 0.01 * err * self.Tsym
sampler.offset -= err sampler.offset -= err
log.debug( log.debug(
'Got %8.1f kB, realtime: %6.2f%%, drift: %+5.2f ppm', 'Got %8.1f kB, realtime: %6.2f%%, drift: %+5.2f ppm',
self.stats['rx_bits'] / 8e3, self.stats['rx_bits'] / 8e3,
duration * 100.0 / (i*config.Tsym), duration * 100.0 / (i*self.Tsym),
(1.0 - sampler.freq) * 1e6 (1.0 - sampler.freq) * 1e6
) )
def start(self, signal, freqs, gain=1.0): def start(self, signal, gain=1.0):
sampler = sampling.Sampler(signal, sampling.Interpolator()) sampler = sampling.Sampler(signal, sampling.Interpolator())
symbols = dsp.Demux(sampler=sampler, omegas=self.omegas, Nsym=self.Nsym)
freq_err = self._prefix(sampler, freq=freqs[0], gain=gain) freq_err = self._prefix(symbols, gain=gain)
sampler.freq -= freq_err sampler.freq -= freq_err
filt = self._train(sampler, order=11, lookahead=6) filt = self._train(sampler, order=11, lookahead=6)
sampler.equalizer = lambda x: list(filt(x)) sampler.equalizer = lambda x: list(filt(x))
bitstream = self._demodulate(sampler, freqs) bitstream = self._demodulate(sampler, symbols)
self.bitstream = itertools.chain.from_iterable(bitstream) self.bitstream = itertools.chain.from_iterable(bitstream)
def run(self, output): def run(self, output):
@@ -237,7 +244,7 @@ class Receiver(object):
def report(self): def report(self):
if self.stats: if self.stats:
duration = time.time() - self.stats['rx_start'] duration = time.time() - self.stats['rx_start']
audio_time = self.stats['rx_bits'] / float(config.modem_bps) audio_time = self.stats['rx_bits'] / float(self.modem_bitrate)
log.debug('Demodulated %.3f kB @ %.3f seconds (%.1f%% realtime)', log.debug('Demodulated %.3f kB @ %.3f seconds (%.1f%% realtime)',
self.stats['rx_bits'] / 8e3, duration, self.stats['rx_bits'] / 8e3, duration,
100 * duration / audio_time if audio_time else 0) 100 * duration / audio_time if audio_time else 0)
@@ -247,13 +254,18 @@ class Receiver(object):
self.plt.figure() self.plt.figure()
symbol_list = np.array(self.stats['symbol_list']) symbol_list = np.array(self.stats['symbol_list'])
for i, freq in enumerate(config.frequencies): for i, freq in enumerate(self.frequencies):
self.plt.subplot(HEIGHT, WIDTH, i+1) self._constellation(symbol_list[i], self.modem.symbols,
self._constellation(symbol_list[i], config.symbols, '$F_c = {} Hz$'.format(freq), index=i)
'$F_c = {} Hz$'.format(freq))
self.plt.show() self.plt.show()
def _constellation(self, y, symbols, title): def _constellation(self, y, symbols, title, index=None):
if index is not None:
Nfreq = len(self.frequencies)
height = np.floor(np.sqrt(Nfreq))
width = np.ceil(Nfreq / float(height))
self.plt.subplot(height, width, index + 1)
theta = np.linspace(0, 2*np.pi, 1000) theta = np.linspace(0, 2*np.pi, 1000)
y = np.array(y) y = np.array(y)
self.plt.plot(y.real, y.imag, '.') self.plt.plot(y.real, y.imag, '.')
@@ -267,6 +279,7 @@ class Receiver(object):
def main(args): def main(args):
config = args.config
reader = stream.Reader(args.input, data_type=common.loads) reader = stream.Reader(args.input, data_type=common.loads)
signal = itertools.chain.from_iterable(reader) signal = itertools.chain.from_iterable(reader)
@@ -275,12 +288,13 @@ def main(args):
reader.check = common.check_saturation reader.check = common.check_saturation
receiver = Receiver(plt=args.plot) detector = Detector(config=config)
receiver = Receiver(config=config, plt=args.plot)
success = False success = False
try: try:
log.info('Waiting for carrier tone: %.1f kHz', config.Fc / 1e3) log.info('Waiting for carrier tone: %.1f kHz', config.Fc / 1e3)
signal, amplitude = detect(signal, config.Fc) signal, amplitude = detector.run(signal)
receiver.start(signal, config.frequencies, gain=1.0/amplitude) receiver.start(signal, gain=1.0/amplitude)
receiver.run(args.output) receiver.run(args.output)
success = True success = True
except Exception: except Exception:

View File

@@ -5,78 +5,72 @@ import itertools
log = logging.getLogger(__name__) log = logging.getLogger(__name__)
from . import train from . import train
from . import wave
from . import common from . import common
from . import config
from . import stream from . import stream
from . import framing from . import framing
from . import equalizer from . import equalizer
from . import dsp from . import dsp
modem = dsp.MODEM(config.symbols) class Sender(object):
def __init__(self, fd, config):
class Writer(object):
def __init__(self, fd):
self.offset = 0 self.offset = 0
self.fd = fd self.fd = fd
self.modem = dsp.MODEM(config.symbols)
self.carriers = config.carriers / config.Nfreq
self.pilot = config.carriers[config.carrier_index]
self.silence = np.zeros(train.silence_length * config.Nsym)
self.iters_per_report = config.baud # report once per second
self.padding = [0] * config.bits_per_baud
self.equalizer = equalizer.Equalizer(config)
def write(self, sym, n=1): def write(self, sym):
sym = np.array(sym) sym = np.array(sym)
data = common.dumps(sym, n) data = common.dumps(sym)
self.fd.write(data) self.fd.write(data)
self.offset += len(data) self.offset += len(sym)
def start(self): def start(self):
carrier = config.carriers[config.carrier_index]
for value in train.prefix: for value in train.prefix:
self.write(carrier * value) self.write(self.pilot * value)
silence = np.zeros(train.silence_length * config.Nsym) symbols = self.equalizer.train_symbols(train.equalizer_length)
symbols = equalizer.train_symbols(train.equalizer_length) signal = self.equalizer.modulator(symbols)
signal = equalizer.modulator(symbols) self.write(self.silence)
self.write(silence)
self.write(signal) self.write(signal)
self.write(silence) self.write(self.silence)
def modulate(self, bits): def modulate(self, bits):
padding = [0] * config.bits_per_baud bits = itertools.chain(bits, self.padding)
bits = itertools.chain(bits, padding) Nfreq = len(self.carriers)
symbols_iter = modem.encode(bits) symbols_iter = common.iterate(self.modem.encode(bits), size=Nfreq)
carriers = config.carriers / config.Nfreq for i, symbols in enumerate(symbols_iter, 1):
for i, symbols in common.iterate(symbols_iter, self.write(np.dot(symbols, self.carriers))
size=config.Nfreq, enumerate=True): if i % self.iters_per_report == 0:
symbols = np.array(list(symbols)) total_bits = i * Nfreq * self.modem.bits_per_symbol
self.write(np.dot(symbols, carriers)) log.debug('Sent %8.1f kB', total_bits / 8e3)
data_duration = (i / config.Nfreq + 1) * config.Tsym
if data_duration % 1 == 0:
bits_size = data_duration * config.modem_bps
log.debug('Sent %8.1f kB', bits_size / 8e3)
def main(args): def main(args):
writer = Writer(args.output) sender = Sender(args.output, config=args.config)
Fs = args.config.Fs
# pre-padding audio with silence # pre-padding audio with silence
writer.write(np.zeros(int(config.Fs * args.silence_start))) sender.write(np.zeros(int(Fs * args.silence_start)))
writer.start() sender.start()
training_size = writer.offset training_duration = sender.offset
training_duration = training_size / wave.bytes_per_second log.info('Sending %.3f seconds of training audio', training_duration / Fs)
log.info('Sending %.3f seconds of training audio', training_duration)
reader = stream.Reader(args.input, bufsize=(64 << 10), eof=True) reader = stream.Reader(args.input, bufsize=(64 << 10), eof=True)
data = itertools.chain.from_iterable(reader) data = itertools.chain.from_iterable(reader)
bits = framing.encode(data) bits = framing.encode(data)
log.info('Starting modulation: %s', modem) log.info('Starting modulation: %s', sender.modem)
writer.modulate(bits=bits) sender.modulate(bits=bits)
data_size = writer.offset - training_size data_duration = sender.offset - training_duration
log.info('Sent %.3f kB @ %.3f seconds', log.info('Sent %.3f kB @ %.3f seconds',
reader.total / 1e3, data_size / wave.bytes_per_second) reader.total / 1e3, data_duration / Fs)
# post-padding audio with silence # post-padding audio with silence
writer.write(np.zeros(int(config.Fs * args.silence_stop))) sender.write(np.zeros(int(Fs * args.silence_stop)))

View File

@@ -1,27 +0,0 @@
import subprocess as sp
import logging
import functools
log = logging.getLogger(__name__)
from . import config
Fs = int(config.Fs) # sampling rate
bits_per_sample = 16
bytes_per_sample = bits_per_sample / 8.0
bytes_per_second = bytes_per_sample * Fs
audio_format = 'S{}_LE'.format(bits_per_sample) # PCM signed little endian
def launch(tool, fname=None, **kwargs):
fname = fname or '-'
args = [tool, fname, '-q', '-f', audio_format, '-c', '1', '-r', str(Fs)]
log.debug('Running: %r', args)
p = sp.Popen(args=args, **kwargs)
return p
# Use ALSA tools for audio playing/recording
play = functools.partial(launch, tool='aplay')
record = functools.partial(launch, tool='arecord')

View File

@@ -1,27 +1,29 @@
from amodem import wave from amodem import audio
import subprocess as sp import subprocess as sp
import signal import signal
def test_launch(): def test_launch():
p = wave.launch(tool='true', fname='fname') p = audio.ALSA(tool='true', Fs=32000).launch(fname='fname')
assert p.wait() == 0 assert p.wait() == 0
def test_exit(): def test_exit():
p = wave.launch(tool='python', fname='-', stdin=sp.PIPE) p = audio.ALSA(tool='python', Fs=32000).launch(fname='-', stdin=sp.PIPE)
s = b'import sys; sys.exit(42)' s = b'import sys; sys.exit(42)'
p.stdin.write(s) p.stdin.write(s)
p.stdin.close() p.stdin.close()
assert p.wait() == 42 assert p.wait() == 42
def test_io(): def test_io():
p = wave.launch(tool='python', fname='-', stdin=sp.PIPE, stdout=sp.PIPE) p = audio.ALSA(tool='python', Fs=32000)
p = p.launch(fname='-', stdin=sp.PIPE, stdout=sp.PIPE)
s = b'Hello World!' s = b'Hello World!'
p.stdin.write(b'print("' + s + b'")\n') p.stdin.write(b'print("' + s + b'")\n')
p.stdin.close() p.stdin.close()
assert p.stdout.read(len(s)) == s assert p.stdout.read(len(s)) == s
def test_kill(): def test_kill():
p = wave.launch(tool='python', fname='-', stdin=sp.PIPE, stdout=sp.PIPE) p = audio.ALSA(tool='python', Fs=32000)
p = p.launch(fname='-', stdin=sp.PIPE, stdout=sp.PIPE)
p.kill() p.kill()
assert p.wait() == -signal.SIGKILL assert p.wait() == -signal.SIGKILL

View File

@@ -1,4 +1,5 @@
from amodem import calib from amodem import calib
from amodem import config
from io import BytesIO from io import BytesIO
@@ -8,10 +9,13 @@ class ProcessMock(object):
self.buf = BytesIO() self.buf = BytesIO()
self.stdin = self self.stdin = self
self.stdout = self self.stdout = self
self.bytes_per_sample = 2
def __call__(self, *args, **kwargs): def launch(self, *args, **kwargs):
return self return self
__call__ = launch
def kill(self): def kill(self):
pass pass
@@ -26,9 +30,9 @@ class ProcessMock(object):
def test_success(): def test_success():
p = ProcessMock() p = ProcessMock()
calib.send(p) calib.send(config, p)
p.buf.seek(0) p.buf.seek(0)
calib.recv(p) calib.recv(config, p)
def test_errors(): def test_errors():
@@ -37,11 +41,11 @@ def test_errors():
def _write(data): def _write(data):
raise IOError() raise IOError()
p.write = _write p.write = _write
calib.send(p) calib.send(config, p)
assert p.buf.tell() == 0 assert p.buf.tell() == 0
def _read(data): def _read(data):
raise KeyboardInterrupt() raise KeyboardInterrupt()
p.read = _read p.read = _read
calib.recv(p, verbose=True) calib.recv(config, p, verbose=True)
assert p.buf.tell() == 0 assert p.buf.tell() == 0

View File

@@ -59,11 +59,12 @@ def test_estimate():
def test_demux(): def test_demux():
freqs = [1e3, 2e3] freqs = np.array([1e3, 2e3])
carriers = [dsp.exp_iwt(f, config.Nsym) for f in freqs] omegas = 2 * np.pi * freqs / config.Fs
carriers = [dsp.exp_iwt(2*np.pi*f/config.Fs, config.Nsym) for f in freqs]
syms = [3, 2j] syms = [3, 2j]
sig = np.dot(syms, carriers) sig = np.dot(syms, carriers)
res = dsp.Demux(sampling.Sampler(sig.real), freqs) res = dsp.Demux(sampling.Sampler(sig.real), omegas, config.Nsym)
res = np.array(list(res)) res = np.array(list(res))
assert np.max(np.abs(res - syms)) < 1e-12 assert np.max(np.abs(res - syms)) < 1e-12

View File

@@ -13,8 +13,9 @@ def assert_approx(x, y, e=1e-12):
def test_training(): def test_training():
L = 1000 L = 1000
t1 = equalizer.train_symbols(L) e = equalizer.Equalizer(config)
t2 = equalizer.train_symbols(L) t1 = e.train_symbols(L)
t2 = e.train_symbols(L)
assert (t1 == t2).all() assert (t1 == t2).all()
@@ -35,10 +36,11 @@ def test_commutation():
def test_modem(): def test_modem():
L = 1000 L = 1000
sent = equalizer.train_symbols(L) e = equalizer.Equalizer(config)
sent = e.train_symbols(L)
gain = config.Nfreq gain = config.Nfreq
x = equalizer.modulator(sent) * gain x = e.modulator(sent) * gain
received = equalizer.demodulator(x, L) received = e.demodulator(x, L)
assert_approx(sent, received) assert_approx(sent, received)
@@ -46,23 +48,24 @@ def test_symbols():
length = 100 length = 100
gain = float(config.Nfreq) gain = float(config.Nfreq)
symbols = equalizer.train_symbols(length=length) e = equalizer.Equalizer(config)
x = equalizer.modulator(symbols) * gain symbols = e.train_symbols(length=length)
assert_approx(equalizer.demodulator(x, size=length), symbols) x = e.modulator(symbols) * gain
assert_approx(e.demodulator(x, size=length), symbols)
den = np.array([1, -0.6, 0.1]) den = np.array([1, -0.6, 0.1])
num = np.array([0.5]) num = np.array([0.5])
y = dsp.lfilter(x=x, b=num, a=den) y = dsp.lfilter(x=x, b=num, a=den)
lookahead = 2 lookahead = 2
h = equalizer.equalize_symbols( h = e.equalize_symbols(
signal=y, symbols=symbols, order=len(den), lookahead=lookahead signal=y, symbols=symbols, order=len(den), lookahead=lookahead
) )
assert norm(h[:lookahead]) < 1e-12 assert norm(h[:lookahead]) < 1e-12
assert_approx(h[lookahead:], den / num) assert_approx(h[lookahead:], den / num)
y = dsp.lfilter(x=y, b=h[lookahead:], a=[1]) y = dsp.lfilter(x=y, b=h[lookahead:], a=[1])
z = equalizer.demodulator(y, size=length) z = e.demodulator(y, size=length)
assert_approx(z, symbols) assert_approx(z, symbols)
@@ -72,9 +75,10 @@ def test_signal():
den = np.array([1, -0.6, 0.1]) den = np.array([1, -0.6, 0.1])
num = np.array([0.5]) num = np.array([0.5])
y = dsp.lfilter(x=x, b=num, a=den) y = dsp.lfilter(x=x, b=num, a=den)
e = equalizer.Equalizer(config)
lookahead = 2 lookahead = 2
h = equalizer.equalize_signal( h = e.equalize_signal(
signal=y, expected=x, order=len(den), lookahead=lookahead) signal=y, expected=x, order=len(den), lookahead=lookahead)
assert norm(h[:lookahead]) < 1e-12 assert norm(h[:lookahead]) < 1e-12

View File

@@ -1,6 +1,7 @@
import numpy as np import numpy as np
from amodem import config from amodem import config
from amodem import dsp
from amodem import recv from amodem import recv
from amodem import train from amodem import train
from amodem import sampling from amodem import sampling
@@ -10,29 +11,34 @@ def test_detect():
P = sum(train.prefix) P = sum(train.prefix)
t = np.arange(P * config.Nsym) * config.Ts t = np.arange(P * config.Nsym) * config.Ts
x = np.cos(2 * np.pi * config.Fc * t) x = np.cos(2 * np.pi * config.Fc * t)
samples, amp = recv.detect(x, config.Fc)
detector = recv.Detector(config)
samples, amp = detector.run(x)
assert abs(1 - amp) < 1e-12 assert abs(1 - amp) < 1e-12
x = np.cos(2 * np.pi * (2*config.Fc) * t) x = np.cos(2 * np.pi * (2*config.Fc) * t)
try: try:
recv.detect(x, config.Fc) detector.run(x)
assert False assert False
except ValueError: except ValueError:
pass pass
def test_prefix(): def test_prefix():
symbol = np.cos(2 * np.pi * config.Fc * np.arange(config.Nsym) * config.Ts) omega = 2 * np.pi * config.Fc / config.Fs
symbol = np.cos(omega * np.arange(config.Nsym))
signal = np.concatenate([c * symbol for c in train.prefix]) signal = np.concatenate([c * symbol for c in train.prefix])
sampler = sampling.Sampler(signal) def symbols_stream(signal):
r = recv.Receiver() sampler = sampling.Sampler(signal)
freq_err = r._prefix(sampler, freq=config.Fc) return dsp.Demux(sampler=sampler, omegas=[omega], Nsym=config.Nsym)
r = recv.Receiver(config)
freq_err = r._prefix(symbols_stream(signal))
assert abs(freq_err) < 1e-16 assert abs(freq_err) < 1e-16
try: try:
silence = 0 * signal silence = 0 * signal
r._prefix(sampling.Sampler(silence), freq=config.Fc) r._prefix(symbols_stream(silence))
assert False assert False
except ValueError: except ValueError:
pass pass
@@ -40,6 +46,7 @@ def test_prefix():
def test_find_start(): def test_find_start():
sym = np.cos(2 * np.pi * config.Fc * np.arange(config.Nsym) * config.Ts) sym = np.cos(2 * np.pi * config.Fc * np.arange(config.Nsym) * config.Ts)
detector = recv.Detector(config)
length = 200 length = 200
prefix = postfix = np.tile(0 * sym, 50) prefix = postfix = np.tile(0 * sym, 50)
@@ -48,6 +55,6 @@ def test_find_start():
prefix = [0] * offset prefix = [0] * offset
bufs = [prefix, prefix, carrier, postfix] bufs = [prefix, prefix, carrier, postfix]
buf = np.concatenate(bufs) buf = np.concatenate(bufs)
start = recv.find_start(buf, length*config.Nsym) start = detector.find_start(buf, length*config.Nsym)
expected = offset + len(prefix) expected = offset + len(prefix)
assert expected == start assert expected == start

View File

@@ -8,6 +8,7 @@ from amodem import recv
from amodem import common from amodem import common
from amodem import dsp from amodem import dsp
from amodem import sampling from amodem import sampling
from amodem import config
import logging import logging
logging.basicConfig(level=logging.DEBUG, logging.basicConfig(level=logging.DEBUG,
@@ -27,7 +28,7 @@ class Args(object):
def run(size, chan=None, df=0, success=True): def run(size, chan=None, df=0, success=True):
tx_data = os.urandom(size) tx_data = os.urandom(size)
tx_audio = BytesIO() tx_audio = BytesIO()
send.main(Args(silence_start=1, silence_stop=1, send.main(Args(config=config, silence_start=1, silence_stop=1,
input=BytesIO(tx_data), output=tx_audio)) input=BytesIO(tx_data), output=tx_audio))
data = tx_audio.getvalue() data = tx_audio.getvalue()
@@ -43,7 +44,8 @@ def run(size, chan=None, df=0, success=True):
rx_audio = BytesIO(data) rx_audio = BytesIO(data)
rx_data = BytesIO() rx_data = BytesIO()
result = recv.main(Args(skip=0, input=rx_audio, output=rx_data)) result = recv.main(Args(config=config,
skip=0, input=rx_audio, output=rx_data))
rx_data = rx_data.getvalue() rx_data = rx_data.getvalue()
assert result == success assert result == success
@@ -61,7 +63,7 @@ def test_small(small_size):
def test_error(): def test_error():
skip = 1 * send.config.Fs # remove trailing silence skip = 32000 # remove trailing silence
run(1024, chan=lambda x: x[:-skip], success=False) run(1024, chan=lambda x: x[:-skip], success=False)