Site menu Radioescuta com Python

Radioescuta com Python

Neste artigo, vamos descrever como funciona um programa de radioescuta automatizada, desenvolvido com a linguagem Python, que funciona com um dongle SDR comum.

Pode ser um brinquedo bem interessante, seja em uso como está (quem não quer ouvir conversas "proibidas" ou "secretas"?) quanto pelos aspectos de processamento digital de sinais que ele embute. Há potencial quase infinito de melhoria e ajustes, que você pode fazer.

Em primeiro lugar, uma demonstração do que nosso programinha pode fazer:

Demonstração do programa capturando alguns canais de intercomunicadores Talkabout das redondezas.

Na sessão de exemplo, o programa fica corujando diversas freqüências GMRS (utilizadas por comunicadores Talkabout e similares). O áudio é gravado apenas quando um sinal "bom" é detectado.

No repositório GitHub, há diversos scripts que utilizam o mesmo programa-base nfm.py para monitorar canais VHF de radioamador, UHF de polícia, comunicadores GMRS ("Talkabout"), freqüências de aviação e ATIS. Apesar do programa chamar-se "nfm", ele também decodifica AM.

Antes que alguém pergunte, radioescuta é legal no Brasil. Divulgar gravações entra numa zona cinza, com chance maior de dar problema se violar a privacidade dos interlocutores (identificando pessoas, lugares e/ou datas).

Banda estreita

As rádios FM comerciais ocupam uma banda de 200kHz cada, a fim de obter boa qualidade musical com tecnologia analógica. Nas faixas de intercomunicação (aviação, polícia, serviços públicos, radioamadorismo em fonia FM) dedicadas à transmissão de voz, basta que o áudio seja inteligível. Cada canal ocupa uma banda muito estreita, em torno de 10kHz.

Nesta seara, a modulação analógica mais comum é FM, embora aviação use AM. A tendência geral é migrar para modos digitais (DMR, P25) que usam menos banda e suportam criptografia. As agências reguladoras (ANATEL, FCC) insistem que novas licenças usem modos digitais para acomodar o maior número possível de usuários nas disputadas bandas sub-GHz.

Também é importante ressaltar que essas faixas de uso comercial são "canalizadas", ou seja, cada canal tem uma frequência e uma banda fixas. Em radioamadorismo, existem faixas não canalizadas, por exemplo Morse (CW) e fonia AM/SSB podem usar qualquer frequência da faixa destinada a esses modos. O ouvinte tem de procurar por um sinal, e fazer a sintonia fina para demodular perfeitamente. (Claro que os radioamadores tentam convencionar frequências e horários para facilitar essa procura.)

Para um serviço comercial, usar canais fixos é muito mais fácil e prático. Os radioamadores também usam mais fonia FM do que qualquer outro modo. Por motivos técnicos, modulação analógica AM e principalmente FM combinam melhor com canais fixos. SSB não combina, por exigir sintonia fina.

Conforme mencionamos no artigo sobre FM comercial, o RTL-SDR não é capaz de filtrar uma banda muito estreita. É preciso capturar uma banda maior e filtrar por software. Por outro lado, como há uma concentração de canais numa faixa relativamente estreita, transformamos o limão numa limonada: temos a possibilidade de escutar muitos canais ao mesmo tempo, com um único dongle.

Script de radioescuta

O script abaixo é o utilizado no vídeo de demonstração:

BW=250000
CHANNEL_BW=12500
FREQS="462562500 462587500 462612500 462637500 462662500 462687500 \
	462712500 462550000 462575000 462600000 462625000 462650000 \
	462675000 462700000 462725000"
STEP=2500

# Determine center automatically
s=0
n=0
for f in $FREQS; do
	s=$(($s + $f))
	n=$(($n + 1))
done
CENTR=$(($s / $n))
CENTR=$(($CENTR / $STEP))
CENTR=$(($CENTR * $STEP))

rtl_sdr -f $CENTR -g 30 -s $BW - | \
	./nfm.py $CENTR $BW $STEP $CHANNEL_BW $FREQS . -a $*

A variável BW é a largura de banda total que desejamos monitorar. Essa banda deve ser compatível com o suportado pelo RTL-SDR (200-300k, 900k-2800k) e deve ser 20% maior que a faixa que realmente nos interessa. Ela também deve ser múltipla de 25000, por questões internas do programa.

A variável CHANNEL_BW é a largura do canal. No caso do Talkabout analógico, é 12.5kHz.

A variável FREQS contém a lista de frequências, dos canais que desejamos monitorar. Quanto mais canais, mais CPU o programa vai usar, até chegar num ponto que ele não dá mais conta de trabalhar real-time. Sugiro testar com poucos canais e ir aumentando o número, sem deixar o uso de CPU chegar a 100% em todos os núcleos.

Note ainda que a distância entre a menor e a maior freqüência de FREQS deve "caber" dentro de BW. Essa distância deve ser uns 20% menor que BW, pois a filtragem do RTL-SDR não é perfeita nos extremos da banda capturada.

A variável STEP é um "denominador comum". Todas as demais variáveis devem ser múltiplos exatos de STEP, e STEP deve ser um divisor de 25000. Isto vai permitir uma série de otimizações. Em geral, STEP será escolhido em função de CHANNEL_BW.

A variável CENTR, que no script nós calculamos automaticamente, é a frequência central de captura para o RTL-SDR, em torno da qual a banda BW se distribui. (Teria sido mais elegante calcular esse valor dentro do Python, mas o cacoete de usar shell script para tudo é difícil de perder.)

Dois parâmetros passados ao rtl_sdr que podem exigir um pouco de experimentação: o ganho (-g) e o desvio de frequência (-p). Para captura automatizada, eu prefiro usar ganho fixo; falaremos disto novamente, no tópico de squelch. Quanto ao desvio, ele será negligenciável se você usa um dongle RTL-SDR genuíno. Se houver desvio, é 1ppm ou 2ppm no máximo. Dongles xing-ling podem ter desvios de 50, 70ppm e variando o tempo todo conforme esquentam.

Agora, vamos ver como funciona o programa nfm.py, onde a mágica realmente acontece.

Blocos funcionais do nfm.py

A demodulação acontece da seguinte forma:

Utilizamos extensivamente a biblioteca NumPy, que nos entrega uma performance muito boa sem precisar apelar para Cython.

A primeira fase não tem muita novidade. As amostras I/Q do RTL-SDR são convertidas para uma lista de números complexos.

# Convert to complex numbers
iqdata = numpy.frombuffer(data, dtype=numpy.uint8)
iqdata = iqdata - 127.5
iqdata = iqdata / 128.0
iqdata = iqdata.view(complex)

# Forward I/Q samples to all channels
for k, d in demodulators.items():
	d.ingest(iqdata)

Usamos threads nos demoduladores de canal pois isto permite aproveitar melhor as CPUs multi-núcleos que todo computador possui hoje em dia. O famigerado "GIL" não é um problema pois o trabalho pesado acontece nas entranhas no NumPy.

Demodulação complexa

O princípio básico da demodulação complexa é o mesmo de qualquer demodulação: mover a banda de interesse para que ela fique em torno do zero, tornando-se assim banda-base, tal qual o sinal (áudio) original.

# Get a cosine table in correct phase
carrier = self.carrier_table[self.if_phase:self.if_phase \
				+ len(iqsamples)]
# Advance phase
self.if_phase = (self.if_phase + len(iqsamples)) % self.if_period
# Demodulate
ifsamples = iqsamples * carrier

A multiplicação de iqsamples por carrier não é entre dois números, e sim entre duas listas, otimizada pelo NumPy.

Porém, em vez de multiplicar por uma portadora composta de números reais, usamos uma sequência complexa. Além disso, precalculamos as amostras da portadora no início da execução e reutilizamos a mesma lista continuamente:

self.if_freq = freq - CENTER
self.carrier_table = \
	[ cmath.exp((0-1j) * tau * t * (self.if_freq / INPUT_RATE)) \
		for t in range(0, INGEST_SIZE * 2) ]

A fórmula com cmath.exp() intimida, mas tudo o que ela faz é criar duas portadoras ao mesmo tempo, defasadas de 90º, uma alojada na parte real das amostras, outra alojada na parte imaginária.

Para que seja possível usar uma portadora precalculada, o período dela deve caber num número inteiro de amostras, e essa é uma razão pela qual todas as bandas e frequências têm de ser múltiplas de STEP.

O sinal I/Q vindo do RTL-SDR também é complexo, assim como essa portadora. A demodulação complexa de um sinal complexo tem uma grande vantagem: ela apenas desloca o espectro, sem criar cópias adicionais do sinal:

Figura 1: De cima para baixo: 1) Espectro de um sinal complexo, que pode ser diferente em frequências positivas e negativas. 2) Espectro de uma portadora complexa. 3) Espectro do sinal modulado, que é uma translação perfeita do original.

Na modulação AM, sempre tínhamos que nos preocupar com as cópias adicionais, eliminando-as com filtros ou mediante técnicas mais elaboradas.

Por que essa diferença? Existem dois motivos, com um mesmo fundamento.

O primeiro é que a portadora puramente real pode ser modelada como a soma de duas portadoras complexas, "girando" em direções opostas, uma com frequência positiva e outra com frequência negativa. O espectro de uma portadora real é simétrico em torno do zero.

Um sinal complexo modulado por essa portadora real produzirá duas cópias, uma delas no espectro negativo:

Figura 2: De cima para baixo: 1) Espectro de um sinal complexo. 2) Espectro de uma portadora real. 3) Espectro do sinal modulado, em duplicata mas não-simétrico em torno de f=0.

O segundo motivo é que todo sinal puramente real possui espectro simétrico (*). Ou seja, o espectro nas frequências negativas é um espelho perfeito do espectro positivo. De um ponto de vista espectral, tal sinal já possui duas cópias, uma delas "oculta".

Quando um sinal real é modulado por uma portadora complexa, as duas cópias são transladadas juntas:

Figura 3: De cima para baixo: 1) Espectro de um sinal real, com um `fantasma` do lado negativo que o torna simétrico. 2) Espectro de uma portadora complexa. 3) Espectro do sinal modulado; o fantasma materializou-se no lado positivo.

Podemos imaginar então que modular um sinal real com uma portadora real produz nada menos que quatro cópias do sinal:

Figura 4: De cima para baixo: 1) Espectro de um sinal real. 2) Espectro de uma portadora real. 3) Espectro do sinal modulado, simétrico e em duplicata.

Pelo aspectro simétrico, temos que o resultado da modulação é puramente real. Isso era esperado, já que sinal e portadora eram ambos reais. (Lembrando que apenas sinais reais podem materializar-se no mundo real, seja na forma de ondas de rádio ou ondas sonoras. Uma tecnologia de modulação precisa sempre cuidar que o resultado final seja realizável.)

Outra característica "mágica" dos sinais complexos é conter informação diferente nas bandas positiva e negativa do espectro. Portadoras positivas ou negativas podem ser usadas para deslocar o espectro para a direita ou para a esquerda, respectivamente.

Se o sinal possui diversos canais, o primeiro passo para isolar o canal que nos interessa é deslocá-lo para a banda-base.

Figura 5: De cima para baixo: 1) Espectro complexo de um sinal captado pelo SDR, contendo diversos canais, com o canal de interesse em azul. 2) Portadora complexa negativo. 3) Espectro do sinal demodulado, com o canal de interesse sobre f=0 (banda-base).

(*) Um sinal real possui espectro simétrico, ou "par" no tocante à amplitude. Já no tocante à fase, ele apresenta espectro anti-simétrico ou "impar".

Filtro passa-baixas

Neste ponto o canal de interesse está na banda-base, porém os demais canais ainda estão presentes no sinal. Removemos os canais indesejados com um filtro passa-baixas.

self.if_filter = filters.low_pass(INPUT_RATE, IF_BANDWIDTH / 2, 48)
...
ifsamples = self.if_filter.feed(ifsamples)

Utilizamos nossa própria implementação de filtro. É um filtro FIR tradicional, que faz uso da operação de convolução do NumPy para obter performance aceitável.

O fato das amostras serem complexas não altera o princípio de funcionamento do filtro, que remove frequências (positivas ou negativas) cujo valor absoluto ultrapasse a frequência de corte.

Redução da taxa de amostragem

Neste ponto, podemos reduzir a taxa de amostragem a um patamar intermediário, a fim de otimizar o processamento do único canal que sobrou.

self.if_decimator = filters.decimator(INPUT_RATE // IF_RATE)
...
ifsamples = self.if_decimator.feed(ifsamples)

Existem duas formas básicas de converter a taxa de amostragem de um sinal digital: decimação para reduzir, e enchimento para aumentar. Uma e outra só fazem conversão inteira (2:1, 3:1, etc.). Para uma conversão racional (3:2, 4:3, etc.) usa-se uma combinação das duas.

Para aumentar a taxa, digamos de 100 para 300 amostras/s, inserimos duas amostras zero depois de cada amostra original. Isto vai adicionar ruído da alta frequência, que é removido simplesmente passando um filtro passa-baixas de 50. (Segundo o teorema de Nyquist, uma taxa de amostragem 100 só comporta frequências até 50.)

Para diminuir a taxa, digamos de 300 para 100, conservamos 1 a cada 3 amostras, descartando as outras 2. Isto transladaria sinais de alta frequência para baixo, pelo fenômeno de "aliasing". A fim de evitar isso, passamos o filtro passa-baixas de 50 antes da decimação.

Já tínhamos feito uma filtragem para eliminar os canais que não nos interessavam, então a decimação pode ser executada logo a seguir.

A taxa intermediária IF_RATE é de 25000, e esta taxa não é configurável na versão atual do programa. Portanto, a taxa de amostragem vinda do RTL-SDR (varíavel BW no script) deve ser múltipla de 25000.

Extração de áudio

O processo de extração de áudio FM é bem descrito neste artigo. Consiste em transformar variações de fase do sinal em amostras de áudio, sem prestar atenção especial à amplitude.

No caso de sinal AM, a detecção é ainda mais simples: a amplitude (ou seja, o valor absoluto das amostras complexas) é interpretado como áudio.

if am:
	output_raw = numpy.absolute(ifsamples)[1:]
	output_raw *= 0.25 / strength
else:
	angles = numpy.angle(ifsamples)
	rotations = numpy.ediff1d(angles)
	rotations = (rotations + numpy.pi) % (2 * numpy.pi) - numpy.pi
	output_raw = numpy.multiply(rotations, DEVIATION_X_SIGNAL)

O único detalhe é que o volume de áudio é escalonado proporcionalmente à portadora adicionada, que no sinal AM demodulado aparece como um grande DC Bias. (Vide o artigo sobre modulação AM para saber porque existe essa portadora adicionada.)

A saída deste estágio é um sinal de áudio puramente real, não é mais complexo. Como dissemos antes, apenas um sinal real pode "materializar-se", no caso em ondas sonoras.

Squelch

"Squelch" é aquela função dos rádios intercomunicadores que corta o ruído branco quando não há recepção, ou quando ela é muito fraca. Neste nosso programa, implementamos algo parecido para gravar áudio apenas quando há recepção real.

if use_autocorrelation:
	vote = self.vote_by_autocorrelation(output_raw)
else:
	vote = self.vote_by_dbfs(dbfs)

A maioria dos rádios implementa squelch com base na potência do sinal recebido, e nós também implementamos essa modalidade (ausência do parâmetro -a). Isto é possível quando a relação sinal/ruído é muito alta, que é o caso dos modos de banda estreita.

strength = numpy.sum(numpy.absolute(ifsamples)) / len(ifsamples)
dbfs = 20 * math.log10(strength)
...
def vote_by_dbfs(self, dbfs):
	vote = -1
	if dbfs > (self.dbfs_off + THRESHOLD_SNR):
		vote = +1
	if not self.is_recording and vote < 0:
		# Use sample to find background noise level
		if dbfs < self.dbfs_off:
			self.dbfs_off = 0.05 * dbfs \
				+ 0.95 * self.dbfs_off
		else:
			self.dbfs_off = 0.005 * dbfs \
				+ 0.995 * self.dbfs_off
	return vote

Como pode ser visto acima, usamos uma média móvel do nível de sinal, e tentamos identificar o nível de sinal correspondente ao ruído de fundo automaticamente.

É importante que o rtl_sdr opere com ganho fixo (parâmetro -g com valor acima de 0). Do contrário, na ausência de recepção, o RTL-SDR amplificaria o ruído de fundo, que pareceria ter a mesma potência de uma recepção boa.

Essa modalidade de squelch funciona particularmente bem com sinais AM, onde a portadora adicionada é indício inconfundível de presença.

A outra modalidade de squelch, que funciona muito bem com sinais FM, é a autocorrelação.

def vote_by_autocorrelation(self, output_raw):
	ac_r = numpy.abs(numpy.correlate(output_raw, output_raw, \
			'same'))
	ac_metric = numpy.sum(ac_r) / numpy.max(ac_r) \
		/ len(output_raw)

	if self.ac_avg is None:
		self.ac_avg = ac_metric
	else:
		self.ac_avg = 0.2 * ac_metric + 0.8 * self.ac_avg

	vote = -1
	if ac_metric > THRESHOLD_AC:
		vote = +1
	return vote

Também consideramos a média móvel da autocorrelação, não o valor instantâneo. Esta técnica merece uma explicação mais alongada.

Correlação é uma medida de semelhança entre dois sinais. Uma forma de calcular a correlação é multiplicar um sinal pelo outro, amostra por amostra, e sumarizar os produtos. Quanto maior o valor absoluto obtido, mais parecidos são os sinais.

Por exemplo, se multiplicarmos cos(x) por sin(x), amostra por amostra, a soma dos produtos será zero, pois são sinais 100% ortogonais, podem inclusive coexistir sem interferência.

Por outro lado, se multiplicarmos cos(x) por uma versão deslocada de sin(x), o resultado será diferente de zero, pois haverá alguma correlação. O caso extremo é multiplicar cos(x) por cos(x), em que todos os produtos são positivos, redundando numa soma total bem grande.

Para calcular de fato a correlação entre dois sinais, devemos repetir a multiplicação por todos os deslocamentos de um sinal — ou seja, calculamos a convolução entre os sinais. Depois, estudamos a área total da convolução para determinar a semelhança.

(Num sinal de comprimento infinito, a convolução também demoraria um tempo infinito para calcular. Mas em DSP tratamos sinais bloco a bloco. A convolução de dois blocos de 400 amostras implicaria, grosso modo, em "apenas" 16 mil multiplicações.)

Autocorrelação é calcular a correlação de um sinal com ele mesmo.

A ideia é a seguinte: se o sinal possui alguma ordem interna, ele apresenta grande autocorrelação. Um sinal que contenha voz ou dados apresenta uma regularidade ao longo do tempo, uma periodicidade, ele se parece com versões defasadas de si mesmo.

Por outro lado, se o sinal é apenas ruído branco, cada bloco é completamente diferente de qualquer outro. O sinal não se parece com nenhuma versão deslocada de si mesmo. A autocorrelação tende a zero, e entendemos isso como ausência de recepção.

Resta a pergunta: por que autocorrelação não funciona tão bem com AM? Provavelmente porque AM é muito sensível ao QRM, o ruído artificial (motores, lâmpadas LED, fontes, etc.). O QRM apresenta uma "regularidade" que ludibria o algoritmo.

Já a modulação FM é notoriamente resistente ao QRM; uma estação fora de sintonia ouve apenas QRN, o ruído natural (térmico e aleatório). Ainda que o controle automático de ganho amplifique esse ruído, ele ainda é detectado como QRN pela autocorrelação.

Segunda redução da taxa, remoção de DC Bias

A taxa de amostragem é novamente reduzida de 25000 para 12500. Como dito antes, o sinal deve ser filtrado antes da decimação. Neste estágio usamos um filtro de 4000Hz, pois o áudio de intercomunicador é muito limitado em qualidade, nunca ultrapassando 3400Hz.

output_raw = self.audio_filter.feed(output_raw)
output = self.audio_decimator.feed(output_raw)
output = self.dc_filter.feed(output)
output = numpy.clip(output, -0.999, +0.999)

Também removemos o DC Bias usando um filtro passa-altas. Ele aparece como um valor constante adicionado ao áudio, principalmente em sinais AM.

Gravação de áudio WAV

Se o algoritmo de squelch determinou que estamos recebendo sinal válido, as amostras são finalmente gravadas num arquivo WAV. Cada canal vai para um arquivo diferente.

if self.insert_timestamp:
	output = numpy.concatenate((self.gen_timestamp(output), \
					output))
	self.insert_timestamp = False
	
bits = struct.pack('<%dh' % len(output), *output)
if not self.wav:
	self.create_wav()
self.wav.writeframes(bits)

Um arquivo WAV conterá inúmeros fragmentos de recepção de um canal. Para que os fragmentos sejam separáveis a posteriori, é adicionado um preâmbulo com algumas pseudo-amostras de áudio, que também embutem um timestamp. Por serem poucas, são inaudíveis.