lab025_wfm_audio_100_1mhz
This commit is contained in:
BIN
data/processed/lab025/lab025_audio_spectrum.png
Normal file
BIN
data/processed/lab025/lab025_audio_spectrum.png
Normal file
Binary file not shown.
|
After Width: | Height: | Size: 136 KiB |
BIN
data/processed/lab025/lab025_audio_waveform.png
Normal file
BIN
data/processed/lab025/lab025_audio_waveform.png
Normal file
Binary file not shown.
|
After Width: | Height: | Size: 237 KiB |
29
data/processed/lab025/lab025_report.txt
Normal file
29
data/processed/lab025/lab025_report.txt
Normal file
@@ -0,0 +1,29 @@
|
|||||||
|
Lab025. Приём и программная WFM-демодуляция
|
||||||
|
|
||||||
|
URI Pluto+: ip:192.168.2.1
|
||||||
|
Частота станции: 100100000 Гц
|
||||||
|
RX LO: 100350000 Гц
|
||||||
|
Цифровое смещение: +250000 Гц
|
||||||
|
Частота дискретизации RX: 2400000 Гц
|
||||||
|
Полоса RX: 1500000 Гц
|
||||||
|
Количество рабочих буферов: 128
|
||||||
|
Количество IQ-сэмплов: 8388608
|
||||||
|
Длительность IQ-записи: 3.495253 с
|
||||||
|
Частота канального сигнала после децимации: 240000 Гц
|
||||||
|
Длительность итогового WAV: 3.445250 с
|
||||||
|
Частота WAV: 48000 Гц
|
||||||
|
Тип PCM: signed PCM16, mono
|
||||||
|
RMS итогового аудио: 0.389148
|
||||||
|
Пиковая амплитуда: 0.974792
|
||||||
|
Отсчёты на границе PCM16: 0.000000 %
|
||||||
|
|
||||||
|
Выходные файлы:
|
||||||
|
WAV: data\processed\lab025\lab025_wfm_audio_100_1mhz.wav
|
||||||
|
Радиоспектр: data\processed\lab025\lab025_rf_spectrum.png
|
||||||
|
Звуковая волна: data\processed\lab025\lab025_audio_waveform.png
|
||||||
|
Спектр звука: data\processed\lab025\lab025_audio_spectrum.png
|
||||||
|
Отчёт: data\processed\lab025\lab025_report.txt
|
||||||
|
|
||||||
|
Использован только приёмный канал RX1.
|
||||||
|
Передатчики TX1 и TX2 не использовались.
|
||||||
|
Необработанные IQ-сэмплы на диск не сохранялись.
|
||||||
BIN
data/processed/lab025/lab025_rf_spectrum.png
Normal file
BIN
data/processed/lab025/lab025_rf_spectrum.png
Normal file
Binary file not shown.
|
After Width: | Height: | Size: 138 KiB |
BIN
data/processed/lab025/lab025_wfm_audio_100_1mhz.wav
Normal file
BIN
data/processed/lab025/lab025_wfm_audio_100_1mhz.wav
Normal file
Binary file not shown.
875
tests/lab025_pluto_wfm_demodulation.py
Normal file
875
tests/lab025_pluto_wfm_demodulation.py
Normal file
@@ -0,0 +1,875 @@
|
|||||||
|
"""
|
||||||
|
Lab025. Приём и программная WFM-демодуляция FM-радиостанции.
|
||||||
|
|
||||||
|
Схема подключения:
|
||||||
|
антенна 40–860 МГц -> RX1
|
||||||
|
|
||||||
|
Используется только приёмный канал RX1. Передатчики TX1 и TX2
|
||||||
|
не используются. Необработанные IQ-сэмплы на диск не сохраняются.
|
||||||
|
"""
|
||||||
|
|
||||||
|
from pathlib import Path
|
||||||
|
|
||||||
|
import adi
|
||||||
|
import matplotlib
|
||||||
|
import numpy as np
|
||||||
|
from scipy import signal
|
||||||
|
from scipy.io import wavfile
|
||||||
|
|
||||||
|
|
||||||
|
# Для лабораторной графики сохраняется в файлы без блокирующих окон.
|
||||||
|
matplotlib.use("Agg")
|
||||||
|
import matplotlib.pyplot as plt
|
||||||
|
|
||||||
|
|
||||||
|
# ---------------------------------------------------------------------
|
||||||
|
# Параметры приёмника Pluto+
|
||||||
|
# ---------------------------------------------------------------------
|
||||||
|
|
||||||
|
PLUTO_URI = "ip:192.168.2.1"
|
||||||
|
|
||||||
|
STATION_FREQUENCY_HZ = 100_100_000
|
||||||
|
LO_OFFSET_HZ = 250_000
|
||||||
|
|
||||||
|
RX_LO_FREQUENCY_HZ = (
|
||||||
|
STATION_FREQUENCY_HZ + LO_OFFSET_HZ
|
||||||
|
)
|
||||||
|
|
||||||
|
SAMPLE_RATE_HZ = 2_400_000
|
||||||
|
RX_BANDWIDTH_HZ = 1_500_000
|
||||||
|
|
||||||
|
RX_BUFFER_SIZE = 65_536
|
||||||
|
DISCARD_BUFFER_COUNT = 3
|
||||||
|
CAPTURE_BUFFER_COUNT = 128
|
||||||
|
|
||||||
|
|
||||||
|
# ---------------------------------------------------------------------
|
||||||
|
# Параметры обработки сигнала
|
||||||
|
# ---------------------------------------------------------------------
|
||||||
|
|
||||||
|
CHANNEL_DECIMATION = 10
|
||||||
|
CHANNEL_SAMPLE_RATE_HZ = 240_000
|
||||||
|
|
||||||
|
AUDIO_CUTOFF_HZ = 15_000
|
||||||
|
AUDIO_SAMPLE_RATE_HZ = 48_000
|
||||||
|
|
||||||
|
DEEMPHASIS_TIME_CONSTANT_SECONDS = 50e-6
|
||||||
|
|
||||||
|
AUDIO_TRANSIENT_DURATION_SECONDS = 0.05
|
||||||
|
AUDIO_REFERENCE_PERCENTILE = 99.5
|
||||||
|
AUDIO_TARGET_LEVEL = 0.85
|
||||||
|
MINIMUM_AUDIO_REFERENCE_PEAK = 1e-12
|
||||||
|
|
||||||
|
RF_WELCH_MAXIMUM_SAMPLE_COUNT = 1_048_576
|
||||||
|
RF_WELCH_SEGMENT_LENGTH = 8_192
|
||||||
|
RF_WELCH_OVERLAP_LENGTH = 4_096
|
||||||
|
|
||||||
|
AUDIO_WELCH_SEGMENT_LENGTH = 8_192
|
||||||
|
|
||||||
|
|
||||||
|
# ---------------------------------------------------------------------
|
||||||
|
# Выходные файлы
|
||||||
|
# ---------------------------------------------------------------------
|
||||||
|
|
||||||
|
OUTPUT_DIRECTORY = Path("data/processed/lab025")
|
||||||
|
|
||||||
|
WAV_FILE_PATH = (
|
||||||
|
OUTPUT_DIRECTORY / "lab025_wfm_audio_100_1mhz.wav"
|
||||||
|
)
|
||||||
|
RF_SPECTRUM_FILE_PATH = (
|
||||||
|
OUTPUT_DIRECTORY / "lab025_rf_spectrum.png"
|
||||||
|
)
|
||||||
|
AUDIO_WAVEFORM_FILE_PATH = (
|
||||||
|
OUTPUT_DIRECTORY / "lab025_audio_waveform.png"
|
||||||
|
)
|
||||||
|
AUDIO_SPECTRUM_FILE_PATH = (
|
||||||
|
OUTPUT_DIRECTORY / "lab025_audio_spectrum.png"
|
||||||
|
)
|
||||||
|
REPORT_FILE_PATH = OUTPUT_DIRECTORY / "lab025_report.txt"
|
||||||
|
|
||||||
|
|
||||||
|
def configure_receiver() -> adi.Pluto:
|
||||||
|
"""
|
||||||
|
Подключается к Pluto+ и настраивает только приёмный канал RX1.
|
||||||
|
|
||||||
|
Частота гетеродина смещена на 250 кГц выше частоты станции.
|
||||||
|
Благодаря этому полезный сигнал не совпадает с аппаратным DC-пиком
|
||||||
|
в центре комплексной полосы приёмника.
|
||||||
|
"""
|
||||||
|
|
||||||
|
print("Подключение к Pluto+...")
|
||||||
|
|
||||||
|
sdr = adi.Pluto(uri=PLUTO_URI)
|
||||||
|
|
||||||
|
# Используем только первый приёмный канал RX1.
|
||||||
|
sdr.rx_enabled_channels = [0]
|
||||||
|
|
||||||
|
sdr.sample_rate = SAMPLE_RATE_HZ
|
||||||
|
sdr.rx_lo = RX_LO_FREQUENCY_HZ
|
||||||
|
sdr.rx_rf_bandwidth = RX_BANDWIDTH_HZ
|
||||||
|
|
||||||
|
# Медленная АРУ подходит для приёма вещательной FM-станции.
|
||||||
|
sdr.gain_control_mode_chan0 = "slow_attack"
|
||||||
|
sdr.rx_buffer_size = RX_BUFFER_SIZE
|
||||||
|
|
||||||
|
return sdr
|
||||||
|
|
||||||
|
|
||||||
|
def receive_samples(sdr: adi.Pluto) -> np.ndarray:
|
||||||
|
"""
|
||||||
|
Получает последовательность комплексных IQ-буферов.
|
||||||
|
|
||||||
|
Первые буферы отбрасываются, поскольку после настройки приёмника
|
||||||
|
в них могут присутствовать переходные процессы АРУ и фильтров.
|
||||||
|
Рабочие буферы объединяются только в оперативной памяти.
|
||||||
|
"""
|
||||||
|
|
||||||
|
print()
|
||||||
|
print(
|
||||||
|
"Отбрасывание переходных буферов: "
|
||||||
|
f"{DISCARD_BUFFER_COUNT}..."
|
||||||
|
)
|
||||||
|
|
||||||
|
for _ in range(DISCARD_BUFFER_COUNT):
|
||||||
|
_ = sdr.rx()
|
||||||
|
|
||||||
|
print("Получение рабочих IQ-буферов...")
|
||||||
|
|
||||||
|
received_buffers: list[np.ndarray] = []
|
||||||
|
|
||||||
|
for buffer_number in range(1, CAPTURE_BUFFER_COUNT + 1):
|
||||||
|
received_buffer = np.asarray(
|
||||||
|
sdr.rx(),
|
||||||
|
dtype=np.complex64,
|
||||||
|
)
|
||||||
|
|
||||||
|
if received_buffer.size == 0:
|
||||||
|
raise RuntimeError(
|
||||||
|
"Pluto+ вернул пустой рабочий IQ-буфер."
|
||||||
|
)
|
||||||
|
|
||||||
|
received_buffers.append(received_buffer)
|
||||||
|
|
||||||
|
if (
|
||||||
|
buffer_number % 16 == 0
|
||||||
|
or buffer_number == CAPTURE_BUFFER_COUNT
|
||||||
|
):
|
||||||
|
print(
|
||||||
|
f" Принято {buffer_number:3d}/"
|
||||||
|
f"{CAPTURE_BUFFER_COUNT} буферов"
|
||||||
|
)
|
||||||
|
|
||||||
|
if not received_buffers:
|
||||||
|
raise RuntimeError("Pluto+ не вернул IQ-сэмплы.")
|
||||||
|
|
||||||
|
combined_samples = np.concatenate(received_buffers)
|
||||||
|
|
||||||
|
if combined_samples.size == 0:
|
||||||
|
raise RuntimeError("После объединения получен пустой IQ-массив.")
|
||||||
|
|
||||||
|
if not np.all(np.isfinite(combined_samples)):
|
||||||
|
raise RuntimeError("В IQ-сэмплах обнаружены NaN или Inf.")
|
||||||
|
|
||||||
|
return np.asarray(combined_samples, dtype=np.complex64)
|
||||||
|
|
||||||
|
|
||||||
|
def shift_station_to_baseband(
|
||||||
|
samples: np.ndarray,
|
||||||
|
sample_rate_hz: float,
|
||||||
|
frequency_shift_hz: float,
|
||||||
|
) -> np.ndarray:
|
||||||
|
"""
|
||||||
|
Переносит выбранную станцию в центр цифровой полосы.
|
||||||
|
|
||||||
|
При положительном смещении комплексный генератор переносит сигнал,
|
||||||
|
расположенный на отрицательной относительной частоте, к 0 Гц.
|
||||||
|
"""
|
||||||
|
|
||||||
|
if samples.size == 0:
|
||||||
|
raise RuntimeError(
|
||||||
|
"Невозможно выполнить цифровой перенос пустого IQ-массива."
|
||||||
|
)
|
||||||
|
|
||||||
|
sample_indices = np.arange(
|
||||||
|
len(samples),
|
||||||
|
dtype=np.float64,
|
||||||
|
)
|
||||||
|
|
||||||
|
digital_oscillator = np.exp(
|
||||||
|
1j
|
||||||
|
* 2.0
|
||||||
|
* np.pi
|
||||||
|
* frequency_shift_hz
|
||||||
|
* sample_indices
|
||||||
|
/ sample_rate_hz
|
||||||
|
).astype(np.complex64)
|
||||||
|
|
||||||
|
centered_samples = samples * digital_oscillator
|
||||||
|
|
||||||
|
if not np.all(np.isfinite(centered_samples)):
|
||||||
|
raise RuntimeError(
|
||||||
|
"После цифрового переноса обнаружены NaN или Inf."
|
||||||
|
)
|
||||||
|
|
||||||
|
return np.asarray(centered_samples, dtype=np.complex64)
|
||||||
|
|
||||||
|
|
||||||
|
def extract_wfm_channel(
|
||||||
|
centered_samples: np.ndarray,
|
||||||
|
) -> np.ndarray:
|
||||||
|
"""
|
||||||
|
Фильтрует WFM-канал и понижает частоту до 240 кГц.
|
||||||
|
|
||||||
|
Полигармоническая передискретизация одновременно выполняет
|
||||||
|
низкочастотную фильтрацию и децимацию в десять раз.
|
||||||
|
"""
|
||||||
|
|
||||||
|
calculated_sample_rate_hz = (
|
||||||
|
SAMPLE_RATE_HZ / CHANNEL_DECIMATION
|
||||||
|
)
|
||||||
|
|
||||||
|
if not np.isclose(
|
||||||
|
calculated_sample_rate_hz,
|
||||||
|
CHANNEL_SAMPLE_RATE_HZ,
|
||||||
|
):
|
||||||
|
raise RuntimeError(
|
||||||
|
"Частота канального сигнала после децимации "
|
||||||
|
"не равна 240 кГц."
|
||||||
|
)
|
||||||
|
|
||||||
|
channel_samples = signal.resample_poly(
|
||||||
|
centered_samples,
|
||||||
|
up=1,
|
||||||
|
down=CHANNEL_DECIMATION,
|
||||||
|
window=("kaiser", 8.0),
|
||||||
|
)
|
||||||
|
|
||||||
|
if channel_samples.size == 0:
|
||||||
|
raise RuntimeError("После выделения WFM-канала массив пуст.")
|
||||||
|
|
||||||
|
if not np.all(np.isfinite(channel_samples)):
|
||||||
|
raise RuntimeError(
|
||||||
|
"После выделения WFM-канала обнаружены NaN или Inf."
|
||||||
|
)
|
||||||
|
|
||||||
|
return np.asarray(channel_samples, dtype=np.complex64)
|
||||||
|
|
||||||
|
|
||||||
|
def demodulate_fm(
|
||||||
|
channel_samples: np.ndarray,
|
||||||
|
) -> np.ndarray:
|
||||||
|
"""
|
||||||
|
Выполняет частотную демодуляцию по фазовой разности отсчётов.
|
||||||
|
|
||||||
|
Результат представляет монофонический композитный FM-сигнал
|
||||||
|
до звуковой фильтрации и коррекции предыскажений.
|
||||||
|
"""
|
||||||
|
|
||||||
|
if channel_samples.size < 2:
|
||||||
|
raise RuntimeError(
|
||||||
|
"Недостаточно канальных отсчётов для FM-демодуляции."
|
||||||
|
)
|
||||||
|
|
||||||
|
phase_difference = np.angle(
|
||||||
|
channel_samples[1:]
|
||||||
|
* np.conj(channel_samples[:-1])
|
||||||
|
).astype(np.float64)
|
||||||
|
|
||||||
|
demodulated = phase_difference - np.mean(phase_difference)
|
||||||
|
|
||||||
|
if demodulated.size == 0:
|
||||||
|
raise RuntimeError("FM-дискриминатор вернул пустой массив.")
|
||||||
|
|
||||||
|
if not np.all(np.isfinite(demodulated)):
|
||||||
|
raise RuntimeError(
|
||||||
|
"После FM-демодуляции обнаружены NaN или Inf."
|
||||||
|
)
|
||||||
|
|
||||||
|
return demodulated
|
||||||
|
|
||||||
|
|
||||||
|
def lowpass_audio(
|
||||||
|
demodulated_samples: np.ndarray,
|
||||||
|
sample_rate_hz: float,
|
||||||
|
) -> np.ndarray:
|
||||||
|
"""
|
||||||
|
Выделяет монофонический звук L+R в полосе до 15 кГц.
|
||||||
|
|
||||||
|
Фильтр подавляет стереопилот 19 кГц, стереоразностную часть,
|
||||||
|
RDS и внеполосный высокочастотный шум.
|
||||||
|
"""
|
||||||
|
|
||||||
|
audio_sos = signal.butter(
|
||||||
|
6,
|
||||||
|
AUDIO_CUTOFF_HZ,
|
||||||
|
btype="lowpass",
|
||||||
|
fs=sample_rate_hz,
|
||||||
|
output="sos",
|
||||||
|
)
|
||||||
|
|
||||||
|
filtered_audio = signal.sosfiltfilt(
|
||||||
|
audio_sos,
|
||||||
|
demodulated_samples,
|
||||||
|
)
|
||||||
|
|
||||||
|
if filtered_audio.size == 0:
|
||||||
|
raise RuntimeError("Звуковой фильтр вернул пустой массив.")
|
||||||
|
|
||||||
|
if not np.all(np.isfinite(filtered_audio)):
|
||||||
|
raise RuntimeError(
|
||||||
|
"После звукового фильтра обнаружены NaN или Inf."
|
||||||
|
)
|
||||||
|
|
||||||
|
return np.asarray(filtered_audio, dtype=np.float64)
|
||||||
|
|
||||||
|
|
||||||
|
def apply_deemphasis(
|
||||||
|
audio_samples: np.ndarray,
|
||||||
|
sample_rate_hz: float,
|
||||||
|
time_constant_seconds: float,
|
||||||
|
) -> np.ndarray:
|
||||||
|
"""
|
||||||
|
Выполняет европейскую FM-коррекцию предыскажений 50 мкс.
|
||||||
|
|
||||||
|
Используется устойчивый однополюсный рекурсивный фильтр,
|
||||||
|
реализованный функцией scipy.signal.lfilter.
|
||||||
|
"""
|
||||||
|
|
||||||
|
alpha = np.exp(
|
||||||
|
-1.0
|
||||||
|
/ (
|
||||||
|
sample_rate_hz
|
||||||
|
* time_constant_seconds
|
||||||
|
)
|
||||||
|
)
|
||||||
|
|
||||||
|
deemphasized = signal.lfilter(
|
||||||
|
[1.0 - alpha],
|
||||||
|
[1.0, -alpha],
|
||||||
|
audio_samples,
|
||||||
|
)
|
||||||
|
|
||||||
|
deemphasized = deemphasized - np.mean(deemphasized)
|
||||||
|
|
||||||
|
if deemphasized.size == 0:
|
||||||
|
raise RuntimeError("De-emphasis вернул пустой массив.")
|
||||||
|
|
||||||
|
if not np.all(np.isfinite(deemphasized)):
|
||||||
|
raise RuntimeError(
|
||||||
|
"После de-emphasis обнаружены NaN или Inf."
|
||||||
|
)
|
||||||
|
|
||||||
|
return np.asarray(deemphasized, dtype=np.float64)
|
||||||
|
|
||||||
|
|
||||||
|
def resample_audio_to_output_rate(
|
||||||
|
audio_samples: np.ndarray,
|
||||||
|
) -> np.ndarray:
|
||||||
|
"""
|
||||||
|
Преобразует звуковой сигнал с 240 кГц в выходные 48 кГц.
|
||||||
|
"""
|
||||||
|
|
||||||
|
output_decimation = 5
|
||||||
|
calculated_output_rate_hz = (
|
||||||
|
CHANNEL_SAMPLE_RATE_HZ / output_decimation
|
||||||
|
)
|
||||||
|
|
||||||
|
if not np.isclose(
|
||||||
|
calculated_output_rate_hz,
|
||||||
|
AUDIO_SAMPLE_RATE_HZ,
|
||||||
|
):
|
||||||
|
raise RuntimeError(
|
||||||
|
"Рассчитанная частота WAV не равна 48 кГц."
|
||||||
|
)
|
||||||
|
|
||||||
|
output_audio = signal.resample_poly(
|
||||||
|
audio_samples,
|
||||||
|
up=1,
|
||||||
|
down=output_decimation,
|
||||||
|
)
|
||||||
|
|
||||||
|
if output_audio.size == 0:
|
||||||
|
raise RuntimeError(
|
||||||
|
"После преобразования в 48 кГц получен пустой массив."
|
||||||
|
)
|
||||||
|
|
||||||
|
if not np.all(np.isfinite(output_audio)):
|
||||||
|
raise RuntimeError(
|
||||||
|
"После преобразования звука обнаружены NaN или Inf."
|
||||||
|
)
|
||||||
|
|
||||||
|
return np.asarray(output_audio, dtype=np.float64)
|
||||||
|
|
||||||
|
|
||||||
|
def convert_audio_to_pcm16(
|
||||||
|
audio_samples: np.ndarray,
|
||||||
|
) -> np.ndarray:
|
||||||
|
"""
|
||||||
|
Удаляет переходный участок, нормализует звук и создаёт PCM16.
|
||||||
|
|
||||||
|
Опорный уровень определяется по 99,5-му процентилю модуля,
|
||||||
|
поэтому единичный выброс не делает весь WAV слишком тихим.
|
||||||
|
"""
|
||||||
|
|
||||||
|
if audio_samples.size == 0:
|
||||||
|
raise RuntimeError("Невозможно нормализовать пустой звук.")
|
||||||
|
|
||||||
|
centered_audio = audio_samples - np.mean(audio_samples)
|
||||||
|
|
||||||
|
transient_sample_count = int(
|
||||||
|
round(
|
||||||
|
AUDIO_TRANSIENT_DURATION_SECONDS
|
||||||
|
* AUDIO_SAMPLE_RATE_HZ
|
||||||
|
)
|
||||||
|
)
|
||||||
|
|
||||||
|
if centered_audio.size <= transient_sample_count:
|
||||||
|
raise RuntimeError(
|
||||||
|
"Звуковой массив короче переходного участка 0,05 с."
|
||||||
|
)
|
||||||
|
|
||||||
|
centered_audio = centered_audio[transient_sample_count:]
|
||||||
|
|
||||||
|
reference_peak = float(
|
||||||
|
np.percentile(
|
||||||
|
np.abs(centered_audio),
|
||||||
|
AUDIO_REFERENCE_PERCENTILE,
|
||||||
|
)
|
||||||
|
)
|
||||||
|
|
||||||
|
if (
|
||||||
|
not np.isfinite(reference_peak)
|
||||||
|
or reference_peak <= MINIMUM_AUDIO_REFERENCE_PEAK
|
||||||
|
):
|
||||||
|
raise RuntimeError(
|
||||||
|
"Сигнал слишком мал для безопасной нормализации PCM16."
|
||||||
|
)
|
||||||
|
|
||||||
|
normalized_audio = (
|
||||||
|
centered_audio
|
||||||
|
/ reference_peak
|
||||||
|
* AUDIO_TARGET_LEVEL
|
||||||
|
)
|
||||||
|
|
||||||
|
normalized_audio = np.clip(
|
||||||
|
normalized_audio,
|
||||||
|
-1.0,
|
||||||
|
1.0,
|
||||||
|
)
|
||||||
|
|
||||||
|
pcm16_audio = np.round(
|
||||||
|
normalized_audio * np.iinfo(np.int16).max
|
||||||
|
).astype(np.int16)
|
||||||
|
|
||||||
|
if pcm16_audio.size == 0:
|
||||||
|
raise RuntimeError("После преобразования PCM16 массив пуст.")
|
||||||
|
|
||||||
|
return pcm16_audio
|
||||||
|
|
||||||
|
|
||||||
|
def calculate_rf_spectrum(
|
||||||
|
samples: np.ndarray,
|
||||||
|
) -> tuple[np.ndarray, np.ndarray]:
|
||||||
|
"""
|
||||||
|
Рассчитывает ограниченный спектр исходного IQ методом Уэлча.
|
||||||
|
|
||||||
|
Для графика используется не более 1 048 576 отсчётов, поэтому
|
||||||
|
полная восьмимиллионная запись не передаётся в одну большую FFT.
|
||||||
|
"""
|
||||||
|
|
||||||
|
analysis_sample_count = min(
|
||||||
|
len(samples),
|
||||||
|
RF_WELCH_MAXIMUM_SAMPLE_COUNT,
|
||||||
|
)
|
||||||
|
|
||||||
|
analysis_samples = samples[:analysis_sample_count]
|
||||||
|
|
||||||
|
relative_frequencies_hz, power_density = signal.welch(
|
||||||
|
analysis_samples,
|
||||||
|
fs=SAMPLE_RATE_HZ,
|
||||||
|
window="hann",
|
||||||
|
nperseg=RF_WELCH_SEGMENT_LENGTH,
|
||||||
|
noverlap=RF_WELCH_OVERLAP_LENGTH,
|
||||||
|
detrend=False,
|
||||||
|
return_onesided=False,
|
||||||
|
scaling="density",
|
||||||
|
)
|
||||||
|
|
||||||
|
relative_frequencies_hz = np.fft.fftshift(
|
||||||
|
relative_frequencies_hz
|
||||||
|
)
|
||||||
|
power_density = np.fft.fftshift(power_density)
|
||||||
|
|
||||||
|
minimum_positive_value = np.finfo(np.float64).tiny
|
||||||
|
power_db = 10.0 * np.log10(
|
||||||
|
power_density + minimum_positive_value
|
||||||
|
)
|
||||||
|
power_db -= np.max(power_db)
|
||||||
|
|
||||||
|
return relative_frequencies_hz, power_db
|
||||||
|
|
||||||
|
|
||||||
|
def create_rf_spectrum_graph(samples: np.ndarray) -> None:
|
||||||
|
"""
|
||||||
|
Сохраняет спектр до цифрового переноса станции в центр.
|
||||||
|
"""
|
||||||
|
|
||||||
|
frequencies_hz, power_db = calculate_rf_spectrum(samples)
|
||||||
|
|
||||||
|
figure, axes = plt.subplots(figsize=(13, 7))
|
||||||
|
|
||||||
|
axes.plot(
|
||||||
|
frequencies_hz / 1e3,
|
||||||
|
power_db,
|
||||||
|
linewidth=0.8,
|
||||||
|
)
|
||||||
|
|
||||||
|
expected_station_offset_hz = -LO_OFFSET_HZ
|
||||||
|
|
||||||
|
axes.axvline(
|
||||||
|
expected_station_offset_hz / 1e3,
|
||||||
|
color="tab:red",
|
||||||
|
linestyle="--",
|
||||||
|
linewidth=1.4,
|
||||||
|
label=(
|
||||||
|
"Ожидаемая станция 100,100 МГц "
|
||||||
|
"(-250 кГц)"
|
||||||
|
),
|
||||||
|
)
|
||||||
|
|
||||||
|
axes.set_title(
|
||||||
|
"Lab025. Спектр принятого IQ до цифрового переноса\n"
|
||||||
|
"RX LO = 100,350 МГц; станция = 100,100 МГц"
|
||||||
|
)
|
||||||
|
axes.set_xlabel("Относительная частота относительно RX LO, кГц")
|
||||||
|
axes.set_ylabel("Относительная мощность, дБ")
|
||||||
|
axes.set_ylim(-90, 5)
|
||||||
|
axes.grid(True)
|
||||||
|
axes.legend()
|
||||||
|
|
||||||
|
figure.tight_layout()
|
||||||
|
figure.savefig(RF_SPECTRUM_FILE_PATH, dpi=160)
|
||||||
|
plt.close(figure)
|
||||||
|
|
||||||
|
|
||||||
|
def create_audio_waveform_graph(pcm16_audio: np.ndarray) -> None:
|
||||||
|
"""
|
||||||
|
Сохраняет первые 0,1 секунды нормированного звука.
|
||||||
|
"""
|
||||||
|
|
||||||
|
displayed_sample_count = min(
|
||||||
|
len(pcm16_audio),
|
||||||
|
int(round(0.1 * AUDIO_SAMPLE_RATE_HZ)),
|
||||||
|
)
|
||||||
|
|
||||||
|
displayed_audio = (
|
||||||
|
pcm16_audio[:displayed_sample_count].astype(np.float64)
|
||||||
|
/ np.iinfo(np.int16).max
|
||||||
|
)
|
||||||
|
time_seconds = (
|
||||||
|
np.arange(displayed_sample_count, dtype=np.float64)
|
||||||
|
/ AUDIO_SAMPLE_RATE_HZ
|
||||||
|
)
|
||||||
|
|
||||||
|
figure, axes = plt.subplots(figsize=(12, 5))
|
||||||
|
axes.plot(time_seconds, displayed_audio, linewidth=0.8)
|
||||||
|
axes.set_title("Lab025. Первые 0,1 секунды демодулированного звука")
|
||||||
|
axes.set_xlabel("Время, с")
|
||||||
|
axes.set_ylabel("Нормированная амплитуда")
|
||||||
|
axes.set_ylim(-1.05, 1.05)
|
||||||
|
axes.grid(True)
|
||||||
|
|
||||||
|
figure.tight_layout()
|
||||||
|
figure.savefig(AUDIO_WAVEFORM_FILE_PATH, dpi=160)
|
||||||
|
plt.close(figure)
|
||||||
|
|
||||||
|
|
||||||
|
def create_audio_spectrum_graph(pcm16_audio: np.ndarray) -> None:
|
||||||
|
"""
|
||||||
|
Рассчитывает методом Уэлча и сохраняет односторонний спектр WAV.
|
||||||
|
"""
|
||||||
|
|
||||||
|
normalized_audio = (
|
||||||
|
pcm16_audio.astype(np.float64)
|
||||||
|
/ np.iinfo(np.int16).max
|
||||||
|
)
|
||||||
|
|
||||||
|
segment_length = min(
|
||||||
|
AUDIO_WELCH_SEGMENT_LENGTH,
|
||||||
|
len(normalized_audio),
|
||||||
|
)
|
||||||
|
|
||||||
|
frequencies_hz, power_density = signal.welch(
|
||||||
|
normalized_audio,
|
||||||
|
fs=AUDIO_SAMPLE_RATE_HZ,
|
||||||
|
window="hann",
|
||||||
|
nperseg=segment_length,
|
||||||
|
noverlap=segment_length // 2,
|
||||||
|
detrend="constant",
|
||||||
|
return_onesided=True,
|
||||||
|
scaling="density",
|
||||||
|
)
|
||||||
|
|
||||||
|
minimum_positive_value = np.finfo(np.float64).tiny
|
||||||
|
power_db = 10.0 * np.log10(
|
||||||
|
power_density + minimum_positive_value
|
||||||
|
)
|
||||||
|
power_db -= np.max(power_db)
|
||||||
|
|
||||||
|
figure, axes = plt.subplots(figsize=(12, 6))
|
||||||
|
axes.plot(frequencies_hz / 1e3, power_db, linewidth=0.9)
|
||||||
|
axes.set_xlim(0, 20)
|
||||||
|
axes.set_ylim(-100, 5)
|
||||||
|
axes.set_title("Lab025. Спектр демодулированного монофонического звука")
|
||||||
|
axes.set_xlabel("Частота, кГц")
|
||||||
|
axes.set_ylabel("Относительная спектральная плотность, дБ")
|
||||||
|
axes.grid(True)
|
||||||
|
|
||||||
|
figure.tight_layout()
|
||||||
|
figure.savefig(AUDIO_SPECTRUM_FILE_PATH, dpi=160)
|
||||||
|
plt.close(figure)
|
||||||
|
|
||||||
|
|
||||||
|
def calculate_audio_statistics(
|
||||||
|
pcm16_audio: np.ndarray,
|
||||||
|
) -> tuple[float, float, float]:
|
||||||
|
"""
|
||||||
|
Возвращает нормированные RMS, пик и процент предельных отсчётов.
|
||||||
|
"""
|
||||||
|
|
||||||
|
normalized_audio = (
|
||||||
|
pcm16_audio.astype(np.float64)
|
||||||
|
/ np.iinfo(np.int16).max
|
||||||
|
)
|
||||||
|
|
||||||
|
rms_value = float(
|
||||||
|
np.sqrt(np.mean(normalized_audio ** 2))
|
||||||
|
)
|
||||||
|
peak_value = float(np.max(np.abs(normalized_audio)))
|
||||||
|
|
||||||
|
clipped_sample_count = int(
|
||||||
|
np.count_nonzero(
|
||||||
|
(pcm16_audio == np.iinfo(np.int16).min)
|
||||||
|
| (pcm16_audio == np.iinfo(np.int16).max)
|
||||||
|
)
|
||||||
|
)
|
||||||
|
clipping_percentage = (
|
||||||
|
100.0
|
||||||
|
* clipped_sample_count
|
||||||
|
/ len(pcm16_audio)
|
||||||
|
)
|
||||||
|
|
||||||
|
return rms_value, peak_value, clipping_percentage
|
||||||
|
|
||||||
|
|
||||||
|
def save_report(
|
||||||
|
iq_sample_count: int,
|
||||||
|
iq_duration_seconds: float,
|
||||||
|
pcm16_audio: np.ndarray,
|
||||||
|
audio_rms: float,
|
||||||
|
audio_peak: float,
|
||||||
|
clipping_percentage: float,
|
||||||
|
) -> None:
|
||||||
|
"""
|
||||||
|
Создаёт текстовый отчёт с параметрами приёма и WAV-файла.
|
||||||
|
"""
|
||||||
|
|
||||||
|
wav_duration_seconds = (
|
||||||
|
len(pcm16_audio) / AUDIO_SAMPLE_RATE_HZ
|
||||||
|
)
|
||||||
|
|
||||||
|
report_lines = [
|
||||||
|
"Lab025. Приём и программная WFM-демодуляция",
|
||||||
|
"",
|
||||||
|
f"URI Pluto+: {PLUTO_URI}",
|
||||||
|
f"Частота станции: {STATION_FREQUENCY_HZ} Гц",
|
||||||
|
f"RX LO: {RX_LO_FREQUENCY_HZ} Гц",
|
||||||
|
f"Цифровое смещение: +{LO_OFFSET_HZ} Гц",
|
||||||
|
f"Частота дискретизации RX: {SAMPLE_RATE_HZ} Гц",
|
||||||
|
f"Полоса RX: {RX_BANDWIDTH_HZ} Гц",
|
||||||
|
f"Количество рабочих буферов: {CAPTURE_BUFFER_COUNT}",
|
||||||
|
f"Количество IQ-сэмплов: {iq_sample_count}",
|
||||||
|
f"Длительность IQ-записи: {iq_duration_seconds:.6f} с",
|
||||||
|
(
|
||||||
|
"Частота канального сигнала после децимации: "
|
||||||
|
f"{CHANNEL_SAMPLE_RATE_HZ} Гц"
|
||||||
|
),
|
||||||
|
f"Длительность итогового WAV: {wav_duration_seconds:.6f} с",
|
||||||
|
f"Частота WAV: {AUDIO_SAMPLE_RATE_HZ} Гц",
|
||||||
|
"Тип PCM: signed PCM16, mono",
|
||||||
|
f"RMS итогового аудио: {audio_rms:.6f}",
|
||||||
|
f"Пиковая амплитуда: {audio_peak:.6f}",
|
||||||
|
f"Отсчёты на границе PCM16: {clipping_percentage:.6f} %",
|
||||||
|
"",
|
||||||
|
"Выходные файлы:",
|
||||||
|
f"WAV: {WAV_FILE_PATH}",
|
||||||
|
f"Радиоспектр: {RF_SPECTRUM_FILE_PATH}",
|
||||||
|
f"Звуковая волна: {AUDIO_WAVEFORM_FILE_PATH}",
|
||||||
|
f"Спектр звука: {AUDIO_SPECTRUM_FILE_PATH}",
|
||||||
|
f"Отчёт: {REPORT_FILE_PATH}",
|
||||||
|
"",
|
||||||
|
"Использован только приёмный канал RX1.",
|
||||||
|
"Передатчики TX1 и TX2 не использовались.",
|
||||||
|
"Необработанные IQ-сэмплы на диск не сохранялись.",
|
||||||
|
]
|
||||||
|
|
||||||
|
REPORT_FILE_PATH.write_text(
|
||||||
|
"\n".join(report_lines),
|
||||||
|
encoding="utf-8",
|
||||||
|
)
|
||||||
|
|
||||||
|
|
||||||
|
def main() -> None:
|
||||||
|
"""
|
||||||
|
Выполняет полный цикл приёма, WFM-демодуляции и сохранения WAV.
|
||||||
|
"""
|
||||||
|
|
||||||
|
OUTPUT_DIRECTORY.mkdir(parents=True, exist_ok=True)
|
||||||
|
|
||||||
|
sdr = None
|
||||||
|
|
||||||
|
try:
|
||||||
|
sdr = configure_receiver()
|
||||||
|
|
||||||
|
print()
|
||||||
|
print("Параметры RX:")
|
||||||
|
print(f" URI: {PLUTO_URI}")
|
||||||
|
print(
|
||||||
|
f" Станция: "
|
||||||
|
f"{STATION_FREQUENCY_HZ / 1e6:.3f} МГц"
|
||||||
|
)
|
||||||
|
print(
|
||||||
|
f" RX LO: "
|
||||||
|
f"{sdr.rx_lo / 1e6:.3f} МГц"
|
||||||
|
)
|
||||||
|
print(
|
||||||
|
f" Частота дискретизации: "
|
||||||
|
f"{sdr.sample_rate / 1e6:.3f} Мвыб/с"
|
||||||
|
)
|
||||||
|
print(
|
||||||
|
f" Полоса RX: "
|
||||||
|
f"{sdr.rx_rf_bandwidth / 1e6:.3f} МГц"
|
||||||
|
)
|
||||||
|
print(
|
||||||
|
f" Режим усиления: "
|
||||||
|
f"{sdr.gain_control_mode_chan0}"
|
||||||
|
)
|
||||||
|
print(f" Размер буфера: {sdr.rx_buffer_size}")
|
||||||
|
print(" Активный канал: RX1")
|
||||||
|
|
||||||
|
samples = receive_samples(sdr)
|
||||||
|
|
||||||
|
iq_sample_count = len(samples)
|
||||||
|
iq_duration_seconds = iq_sample_count / SAMPLE_RATE_HZ
|
||||||
|
|
||||||
|
print()
|
||||||
|
print(f"Количество IQ-сэмплов: {iq_sample_count}")
|
||||||
|
print(f"Длительность записи: {iq_duration_seconds:.3f} с")
|
||||||
|
|
||||||
|
# Удаляем остаточную комплексную постоянную составляющую.
|
||||||
|
samples = samples - np.mean(samples)
|
||||||
|
|
||||||
|
print()
|
||||||
|
print("Цифровой перенос станции в центр полосы...")
|
||||||
|
centered_samples = shift_station_to_baseband(
|
||||||
|
samples=samples,
|
||||||
|
sample_rate_hz=SAMPLE_RATE_HZ,
|
||||||
|
frequency_shift_hz=LO_OFFSET_HZ,
|
||||||
|
)
|
||||||
|
|
||||||
|
print("Фильтрация WFM-канала и децимация до 240 кГц...")
|
||||||
|
channel_samples = extract_wfm_channel(centered_samples)
|
||||||
|
|
||||||
|
print("FM-демодуляция...")
|
||||||
|
demodulated_samples = demodulate_fm(channel_samples)
|
||||||
|
|
||||||
|
print("Звуковой low-pass фильтр 0–15 кГц...")
|
||||||
|
filtered_audio = lowpass_audio(
|
||||||
|
demodulated_samples=demodulated_samples,
|
||||||
|
sample_rate_hz=CHANNEL_SAMPLE_RATE_HZ,
|
||||||
|
)
|
||||||
|
|
||||||
|
print("De-emphasis 50 мкс...")
|
||||||
|
deemphasized_audio = apply_deemphasis(
|
||||||
|
audio_samples=filtered_audio,
|
||||||
|
sample_rate_hz=CHANNEL_SAMPLE_RATE_HZ,
|
||||||
|
time_constant_seconds=(
|
||||||
|
DEEMPHASIS_TIME_CONSTANT_SECONDS
|
||||||
|
),
|
||||||
|
)
|
||||||
|
|
||||||
|
print("Преобразование звука в 48 кГц...")
|
||||||
|
output_audio = resample_audio_to_output_rate(
|
||||||
|
deemphasized_audio
|
||||||
|
)
|
||||||
|
pcm16_audio = convert_audio_to_pcm16(output_audio)
|
||||||
|
|
||||||
|
audio_rms, audio_peak, clipping_percentage = (
|
||||||
|
calculate_audio_statistics(pcm16_audio)
|
||||||
|
)
|
||||||
|
|
||||||
|
print("Сохранение WAV...")
|
||||||
|
wavfile.write(
|
||||||
|
WAV_FILE_PATH,
|
||||||
|
AUDIO_SAMPLE_RATE_HZ,
|
||||||
|
pcm16_audio,
|
||||||
|
)
|
||||||
|
|
||||||
|
if not WAV_FILE_PATH.is_file():
|
||||||
|
raise RuntimeError("Выходной WAV-файл не создан.")
|
||||||
|
|
||||||
|
print("Сохранение графиков...")
|
||||||
|
create_rf_spectrum_graph(samples)
|
||||||
|
create_audio_waveform_graph(pcm16_audio)
|
||||||
|
create_audio_spectrum_graph(pcm16_audio)
|
||||||
|
|
||||||
|
print("Создание текстового отчёта...")
|
||||||
|
save_report(
|
||||||
|
iq_sample_count=iq_sample_count,
|
||||||
|
iq_duration_seconds=iq_duration_seconds,
|
||||||
|
pcm16_audio=pcm16_audio,
|
||||||
|
audio_rms=audio_rms,
|
||||||
|
audio_peak=audio_peak,
|
||||||
|
clipping_percentage=clipping_percentage,
|
||||||
|
)
|
||||||
|
|
||||||
|
print()
|
||||||
|
print("Созданы файлы:")
|
||||||
|
print(f" WAV: {WAV_FILE_PATH}")
|
||||||
|
print(f" Радиоспектр: {RF_SPECTRUM_FILE_PATH}")
|
||||||
|
print(f" Звуковая волна: {AUDIO_WAVEFORM_FILE_PATH}")
|
||||||
|
print(f" Спектр звука: {AUDIO_SPECTRUM_FILE_PATH}")
|
||||||
|
print(f" Отчёт: {REPORT_FILE_PATH}")
|
||||||
|
print()
|
||||||
|
print(f"RMS аудио: {audio_rms:.6f}")
|
||||||
|
print(f"Пиковая амплитуда: {audio_peak:.6f}")
|
||||||
|
print(
|
||||||
|
f"Отсчёты на границе PCM: "
|
||||||
|
f"{clipping_percentage:.6f} %"
|
||||||
|
)
|
||||||
|
print()
|
||||||
|
print("Lab025 выполнена успешно.")
|
||||||
|
print()
|
||||||
|
print("Для прослушивания откройте:")
|
||||||
|
print(WAV_FILE_PATH)
|
||||||
|
|
||||||
|
finally:
|
||||||
|
if sdr is not None:
|
||||||
|
destroy_buffer = getattr(
|
||||||
|
sdr,
|
||||||
|
"rx_destroy_buffer",
|
||||||
|
None,
|
||||||
|
)
|
||||||
|
|
||||||
|
if callable(destroy_buffer):
|
||||||
|
try:
|
||||||
|
destroy_buffer()
|
||||||
|
except Exception as cleanup_error:
|
||||||
|
print(
|
||||||
|
"Предупреждение: не удалось освободить "
|
||||||
|
f"RX-буфер: {cleanup_error}"
|
||||||
|
)
|
||||||
|
|
||||||
|
|
||||||
|
if __name__ == "__main__":
|
||||||
|
main()
|
||||||
Reference in New Issue
Block a user