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

воскресенье, 19 мая 2019 г.

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

Привет Хабр.

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

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

SoapySDR

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

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

Библиотека SoapySDR является кроссплатформенной, т.е. написанный под неё код будет работать под Windows, OSX, Linux, и даже на Raspberry Pi. Для Windows нужные библиотеки входят в состав пакета PothosSDR, для остальных платформ скомпилировать SoapySDR придется самостоятельно. Нужно скомпилировать две части — собственно библиотеку, и «плагин» для нужного приемника, в нашем случае это будет SoapyRTLSDR (под 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-данные с фотокамеры. Кому интересно подробнее, могут почитать например здесь. Для нас достаточно знать, что мы эти данные можем записать, а другие 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, https://www.n2yo.com/passes/?s=28654 и https://www.n2yo.com/passes/?s=33591 соответственно.

Чтобы не сидеть у компьютера, добавим в программу ожидание нужного времени. Это позволит также запускать программу на 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м звездной величины, предел человеческого глаза в идеальных условиях 6м. Т.е. при реально темном небе, далеко за городом, пролет спутника NOAA теоретически заметить можно, в городе конечно, шансов нет (как писали на Хабре, уже выросло поколение людей, никогда в жизни не видевших Млечный Путь).

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

Мы видим вполне различимый сигнал, хотя конечно со специальной антенной для приема NOAA качество было бы гораздо лучше. Формат сигнала называется APT (Automatic Picture Transmission), из него можно получить изображение земной поверхности, если кому интересно, можно отдельно рассмотреть его декодирование. Но есть разумеется, и готовые программы, декодировать такие сигналы можно с помощью 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 на примере создания нескольких виртуальных приемников на базе одного «железного».

Let's block ads! (Why?)

Парсер Хабра

...

воскресенье, 19 мая 2019 г.

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

SoapySDR

Прием NOAA

Комментариев нет:

Отправить комментарий

...

воскресенье, 19 мая 2019 г.

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

SoapySDR

Прием NOAA

Комментариев нет:

Отправить комментарий

воскресенье, 19 мая 2019 г.