Software Defined Radio — как это работает? Часть 3

Во второй части ^[1] были рассмотрены практические аспекты использования SDR. В этой части мы разберемся, как принять данные метеоспутника NOAA с помощью Python и недорогого (30$) приемника RTL-SDR. Рассмотренный код будет работать везде — на Windows, OSX, Linux и даже на Raspberry Pi.

Кому интересно, продолжение под катом.

SoapySDR

Производителей различных SDR-устройств довольно много, и поддерживать каждое по отдельности было бы весьма неудобно, да и дорого в плане покупки «железа» для тестирования. В принципе, для унифицированного доступа существует две библиотеки, которые по сути стали стандартом. Первая — это уже довольно старый интерфейс ExtIO DLL ^[2], которому наверно не менее 10 лет, вторая — более современная библиотека SoapySDR ^[3], которую мы и рассмотрим.

SoapySDR — это набор кроссплатформенных библиотек, написанных на C++, предоставляющих унифицированный доступ к SDR-устройствам, как к приемникам, так и трансиверам. Если производитель делает такой интерфейс, то его устройство «автоматом» будет работать с достаточно большим числом популярных программ (GQRX, GNU Radio, CubicSDR и пр). Практически все адекватные производители, кроме некоторых, (пользуясь случаем, передаю привет компании ЕЕ) имеют поддержку SoapySDR, список поддерживаемых устройств можно посмотреть на странице проекта ^[3]. Как можно видеть, он довольно большой, и включает HackRF, USRP, SDRPlay, LimeSDR, RTL-SDR, Red Pitaya и многие другие.

Библиотека SoapySDR является кроссплатформенной, т.е. написанный под неё код будет работать под Windows, OSX, Linux, и даже на Raspberry Pi. Для Windows нужные библиотеки входят в состав пакета PothosSDR ^[4], для остальных платформ скомпилировать SoapySDR придется самостоятельно. Нужно скомпилировать две части — собственно библиотеку ^[5], и «плагин» для нужного приемника, в нашем случае это будет SoapyRTLSDR ^[6] (под Windows библиотеку тоже можно собрать из исходников, для этого понадобится Visual Studio, Cmake, и SWIG). Теперь все готово, и можно писать код.

Импортируем библиотеку и получаем список приемников:

from __future__ import print_function
import SoapySDR

# Enumerate devices
print("SDR devices:")
for d in SoapySDR.Device.enumerate(''):
    print(d)
print()

Подключаем приемник, запускаем код и видим список устройств, среди которых есть наш rtlsdr.

Остальные устройства это звуковые карты, как мы помним, исторически первые SDR работали именно через линейный вход ПК, и библиотека их тоже поддерживает. Получаем информацию об устройстве — число доступных каналов, частотный диапазон и пр:

soapy_device = "rtlsdr"
device = SoapySDR.Device(dict(driver = soapy_device))

channels = list(range(device.getNumChannels(SoapySDR.SOAPY_SDR_RX)))
print("Channels:", channels)

ch = channels[0]

sample_rates = device.listSampleRates(SoapySDR.SOAPY_SDR_RX, ch)
print("Sample rates:n", sample_rates)

bandwidths = list(map(lambda r: int(r.maximum()), device.getBandwidthRange(SoapySDR.SOAPY_SDR_RX, ch)))
print("Bandwidths:n", bandwidths)

print("Gain controls:")
for gain in device.listGains(SoapySDR.SOAPY_SDR_RX, ch):
    print("  %s: %s" % (gain, device.getGainRange(SoapySDR.SOAPY_SDR_RX, ch, gain)))

frequencies = device.listFrequencies(SoapySDR.SOAPY_SDR_RX, ch)
print("Frequencies names:", frequencies)

frequency_name = frequencies[0]
print("Frequency channel name:", frequency_name)

print("Frequency range:", device.getFrequencyRange(SoapySDR.SOAPY_SDR_RX, ch, frequency_name)[0])

Запускаем программу и видим информацию о приемнике:

Мы видим что приемник имеет один входной канал с названием «RF», возможные частоты дискретизации [250000.0, 1024000.0, 1536000.0, 1792000.0, 1920000.0, 2048000.0, 2160000.0, 2560000.0, 2880000.0, 3200000.0] и частотный диапазон 24МГц-1.7ГГц.

Лайфхак — те же данные можно получить и из командной строки, набрав команду SoapySDRUtil --probe=«driver=rtlsdr».

Зная это, мы можем записать поток данных в WAV. Как говорилось в предыдущей части, данные с SDR представлены потоком сигналов, называемых I и Q, представляющих собой отсчеты с АЦП, грубо их можно представить как RAW-данные с фотокамеры. Кому интересно подробнее, могут почитать например здесь ^[7]. Для нас достаточно знать, что мы эти данные можем записать, а другие SDR-программы потом могут с ними работать.

Сама по себе запись довольно проста — функция readStream заполняет буфер если есть данные, если данных еще нет, то вернется -1. Ниже показан код записи 10 отсчетов (несущественные части кода опущены).

device.setFrequency(SoapySDR.SOAPY_SDR_RX, channel, "RF", frequency)
device.setGain(SoapySDR.SOAPY_SDR_RX, channel, "TUNER", gain)
device.setGainMode(SoapySDR.SOAPY_SDR_RX, channel, False)
device.setSampleRate(SoapySDR.SOAPY_SDR_RX, channel, sample_rate)

# Number of blocks to save
block, max_blocks = 0, 10

block_size = device.getStreamMTU(stream)
print("Block size:", block_size)

buffer_format = np.int8
buffer_size = 2*block_size  # I+Q
buffer = np.empty(buffer_size, buffer_format)

while True:
    d_info = device.readStream(stream, [buffer], buffer_size)
    if d_info.ret > 0:
        wav.write(buffer[0:2*d_info.ret])
        print("Bytes saved:", 2*d_info.ret)
        block += 1
        if block > max_blocks:
            break

Результат на скриншоте:

Как можно видеть, мы получаем от устройства блоки данных, размер одного блока составляет 131072 байта, что при частоте дискретизации 250000 дает нам длительность примерно полсекунды. В общем, кто работал ранее со звуковыми картами в Windows, найдет много общего.

Для теста записываем файл и проверяем, что все нормально — его можно воспроизвести в SDR#. Тут есть еще одна хитрость — чтобы SDR# корректно показал частоты станций, имя файла нужно записать в формате, совместимом с HDSDR, вида «HDSDR_20190518_115500Z_101000kHz_RF.wav» (как нетрудно догадаться, в начале идет дата и время в GMT, затем частота в килогерцах). Это нетрудно записать на Python:

frequency = 101000000
file_name = "HDSDR_%s_%dkHz_RF.wav" % (datetime.datetime.utcnow().strftime("%Y%m%d_%H%M%SZ"), frequency/1000)

Для начала проверяем на FM-диапазоне. Все нормально, станции видно, музыка играется, RDS работает.

Можно приступать к записи NOAA.

Прием NOAA

Итак, у нас есть приемник и есть программа записи. Нам будут интересны метеоспутники NOAA 15, NOAA 18 и NOAA 19, передающие изображения поверхности Земли на частотах 137.620, 137.9125 и 137.100МГц. Основная сложность здесь в том, что нужно «поймать» момент, когда спутник пролетает над нами. Узнать время пролета можно онлайн по ссылкам https://www.n2yo.com/passes/?s=25338 ^[8], https://www.n2yo.com/passes/?s=28654 ^[9] и https://www.n2yo.com/passes/?s=33591 ^[10] соответственно.

Чтобы не сидеть у компьютера, добавим в программу ожидание нужного времени. Это позволит также запускать программу на Raspberry Pi, без дисплея и клавиатуры.

import datetime

def wait_for_start(dt):
    # Wait for the start
    while True:
        now = datetime.datetime.now()
        diff = int((dt - now).total_seconds())
        print("{:02d}:{:02d}:{:02d}: Recording will be started after {}m {:02d}s...".format(now.hour, now.minute, now.second, int(diff / 60), diff % 60))
        time.sleep(5)
        if diff <= 1:
            break

wait_for_start(datetime.datetime(2019, 5, 18, 21, 49, 0))

Кстати, чтобы запустить скрипт на Raspberry Pi и оставить его работать после закрытия консоли, нужно ввести команду «nohup python recorder.py &».

Все готово, запускаем скрипт и можем заниматься другими делами, запись длится примерно 20 минут. Параллельно может возникнуть вопрос — можно ли увидеть пролет спутника невооруженным глазом? Согласно таблице, его максимальная яркость составляет порядка 5.5м звездной величины ^[11], предел человеческого глаза в идеальных условиях 6м. Т.е. при реально темном небе, далеко за городом, пролет спутника NOAA теоретически заметить можно, в городе конечно, шансов нет (как писали на Хабре ^[12], уже выросло поколение людей, никогда в жизни не видевших Млечный Путь ^[13]).

Результатом работы скрипта является записанный wav-файл, его спектр показан на скриншоте.

Мы видим вполне различимый сигнал, хотя конечно со специальной антенной ^[14] для приема NOAA качество было бы гораздо лучше. Формат сигнала называется APT (Automatic Picture Transmission ^[15]), из него можно получить изображение земной поверхности, если кому интересно, можно отдельно рассмотреть его декодирование. Но есть разумеется, и готовые программы, декодировать такие сигналы можно с помощью WxToImg или MultiPSK.

Интересно видеть допплеровский сдвиг на спектре, который происходит оттого, что спутник пролетает мимо нас. Наверное, несложно посчитать и его скорость, зная время и сдвиг частоты в герцах.

Разумеется, программу можно использовать не только для записи NOAA, любую ширину полосы и частоту можно указать в настройках. Для тех кто захочет поэкспериментировать с SoapySDR самостоятельно, код программы целиком размещен под спойлером.

from __future__ import print_function
import SoapySDR
import numpy as np
import struct
import sys
import time
import datetime


def wait_for_start(dt):
    # Wait for the start
    while True:
        now = datetime.datetime.now()
        diff = int((dt - now).total_seconds())
        print("{:02d}:{:02d}:{:02d}: Recording will be started after {}m {:02d}s...".format(now.hour, now.minute, now.second, int(diff / 60), diff % 60))
        time.sleep(5)
        if diff <= 1:
            break


def sdr_enumerate():
    # Enumerate SDR devices
    print("SDR devices:")
    for d in SoapySDR.Device.enumerate(''):
        print(d)
    print()


def sdr_init():
    soapy_device = "rtlsdr"
    device = SoapySDR.Device({"driver": soapy_device})

    channels = list(range(device.getNumChannels(SoapySDR.SOAPY_SDR_RX)))
    print("Channels:", channels)

    ch = channels[0]

    sample_rates = device.listSampleRates(SoapySDR.SOAPY_SDR_RX, ch)
    print("Sample rates:n", sample_rates)

    print("Gain controls:")
    for gain in device.listGains(SoapySDR.SOAPY_SDR_RX, ch):
        print("  %s: %s" % (gain, device.getGainRange(SoapySDR.SOAPY_SDR_RX, ch, gain)))

    frequencies = device.listFrequencies(SoapySDR.SOAPY_SDR_RX, ch)
    print("Frequencies names:", frequencies)

    frequency_name = frequencies[0]
    print("Frequency channel name:", frequency_name)

    print("Frequency range:", device.getFrequencyRange(SoapySDR.SOAPY_SDR_RX, ch, frequency_name)[0])

    print()
    return device


def sdr_record(device, frequency, sample_rate, gain, blocks_count):
    print("Frequency:", frequency)
    print("Sample rate:", sample_rate)
    print("Gain:", gain)

    channel = 0  # Always for RTL-SDR
    device.setFrequency(SoapySDR.SOAPY_SDR_RX, channel, "RF", frequency)
    device.setGain(SoapySDR.SOAPY_SDR_RX, channel, "TUNER", gain)
    device.setGainMode(SoapySDR.SOAPY_SDR_RX, channel, False)
    device.setSampleRate(SoapySDR.SOAPY_SDR_RX, channel, sample_rate)

    data_format = SoapySDR.SOAPY_SDR_CS8 # if 'rtlsdr' in soapy_device or 'hackrf' in soapy_device else SoapySDR.SOAPY_SDR_CS16
    stream = device.setupStream(SoapySDR.SOAPY_SDR_RX, data_format, [channel], {})
    device.activateStream(stream)

    block_size = device.getStreamMTU(stream)
    print("Block size:", block_size)
    print("Data format:", data_format)
    print()

    # IQ: 2 digits ver variable
    buffer_format = np.int8
    buffer_size = 2 * block_size # I + Q samples
    buffer = np.empty(buffer_size, buffer_format)

    # Number of blocks to save
    block, max_blocks = 0, blocks_count

    # Save to file
    file_name = "HDSDR_%s_%dkHz_RF.wav" % (datetime.datetime.utcnow().strftime("%Y%m%d_%H%M%SZ"), frequency/1000)
    print("Saving file:", file_name)
    with open(file_name, "wb") as wav:
        # Wav data info
        bits_per_sample = 16
        channels_num, samples_num = 2, int(max_blocks * block_size)
        subchunk_size = 16  # always 16 for PCM
        subchunk2_size = int(samples_num * channels_num * bits_per_sample / 8)
        block_alignment = int(channels_num * bits_per_sample / 8)

        # Write RIFF header
        wav.write('RIFF'.encode('utf-8'))
        wav.write(struct.pack('<i', 4 + (8 + subchunk_size) + (8 + subchunk2_size)))  # Size of the overall file
        wav.write('WAVE'.encode('utf-8'))
        # Write fmt subchunk
        wav.write('fmt '.encode('utf-8'))  # chunk type
        wav.write(struct.pack('<i', subchunk_size))  # subchunk data size (16 for PCM)
        wav.write(struct.pack('<h', 1))  # compression type 1 - PCM
        wav.write(struct.pack('<h', channels_num))  # channels
        wav.write(struct.pack('<i', int(sample_rate)))  # sample rate
        wav.write(struct.pack('<i', int(sample_rate * bits_per_sample * channels_num/ 8)))  # byte rate
        wav.write(struct.pack('<h', block_alignment))  # block alignment
        wav.write(struct.pack('<h', bits_per_sample))  # sample depth
        # Write data subchunk
        wav.write('data'.encode('utf-8'))
        wav.write(struct.pack('<i', subchunk2_size))
        while True:
            d_info = device.readStream(stream, [buffer], buffer_size)
            if d_info.ret > 0:
                data = buffer[0:2*d_info.ret]
                fileData = data
                if data_format == SoapySDR.SOAPY_SDR_CS8:
                   fileData = data.astype('int16')
                wav.write(fileData)
                print("Block %d saved: %d bytes" % (block, 2*d_info.ret))
                block += 1
                if block > max_blocks:
                    break

    device.deactivateStream(stream)
    device.closeStream(stream)

if __name__ == "__main__":
    print("App started")

    # Forecast for active NOAA satellites
    # NOAA 15: 137.620, https://www.n2yo.com/passes/?s=25338
    # NOAA 18: 137.9125, https://www.n2yo.com/passes/?s=28654
    # NOAA 19: 137.100, https://www.n2yo.com/passes/?s=33591
    # Wait for the start: 18-May 21:49 21:49:
    wait_for_start(datetime.datetime(2019, 5, 18, 21, 49, 0))

    device = sdr_init()

    t_start = time.time()

    sdr_record(device, frequency=137912500, sample_rate=250000, gain=35, blocks_count=2100)

    print("Recording complete, time = %ds" % int(time.time() - t_start))
    print()

Плюс SoapySDR в том, что эта же программа с минимальными изменениями будет работать и с другими приемниками, например с SDRPlay или HackRF. Ну и про кроссплатформенность тоже уже упоминалось.

Если у читателей еще остался интерес к теме радиоприема, можно рассмотреть пример использования SDR с GNU Radio на примере создания нескольких виртуальных приемников на базе одного «железного».

Автор: DmitrySpb79

Источник ^[16]