[Из песочницы] Разработка многоканального SDR
Расскажу о своем опыте разработки цифрового многоканального широкополосного приемника.
Много лет работаю в области создания средств захвата и обработки сигналов от бортовых и береговых навигационных РЛС. Примерно года два назад выпустил последнюю, на сегодняшний день, версию нашей платы RVAQ (Radar Video AcQuisition) и задумался, чем в жизни заняться дальше. Хотелось чего-то нового и неизведанного. Выбор пал на неохваченную мной до сих пор область — цифровой радиоприем с легким заходом в СВЧ-область.
Это первая глава, посвященная начальной постановке задачи.
С чего начать, если никогда профессионально не занимался цифровым радиоприемом, если не считать приемник «Мишутка», собранный в детстве? Конечно, с освежения памяти чтения Полякова и модели в MATLAB. Исходной идеей являлось создание многоканального приемника диапазона 156–162МГц для мониторинга и записи всех активных переговоров в морском УКВ-диапазоне. Перечислю желаемые свойства такого приемника:
1. Полоса приема не менее 6 МГц (162–156=6)
2. Чувствительность не хуже -110дБм, а то засмеют
3. Большой динамический диапазон, так как когда слушаешь на берегу судно в море за 30 миль, рядом обязательно кто-то заорет своими 25-ю ваттами. Для приличных приемников уровень блокирования должен быть не меньше 70дБ. Забегая немного вперед, скажу, что получилось реализовать более 90дБ по блокированию. Одним словом, RTL-SDR решительно не соответствовал ожиданиям. Прикидку начал, как ни странно, с выбора АЦП. Так как если в природе нет соответствующих (хотя бы в теории) АЦП, то и браться не стоит. Такой АЦП нашелся.
Теперь нужно выбрать архитектуру приемника. Обзор актуальных решений, изучение элементной базы и интуиция позволили остановиться на приемнике прямого преобразования. Также решено было перенести в ноль частоты интересующий участок спектра при помощи квадратурного демодулятора и работать в первой зоне Найквиста, чтобы максимально утилизировать все качества выбранного АЦП.
clear all;
k = 1.381e-23; % Joule/K - Boltzmann's Constant
T0 = 290; % K - temperature
% Encoding Windows-1251
% 1. Функциональная схема
% 1010 - Схема защиты входных каскадов
% 1011 - Преселектор
% 1012 - Управляемый МШУ (вкл/выкл.)
% 1013 - Фильтр основной селекции
% 1014 - Управляемый усилитель петли АРУ
% 1015 - Квадратурный смеситель
% 1021, 1022 - Фильтры низкой частоты.
% 1031, 1032 - Низкочастотные усилители
% 1041, 1042 - Фильтры низкой частоты АЦП
% 1051, 1052 - АЦП
% 2. Системная диаграмма радиотракта для анализа мощности шумов и чувствительности
Rrf_inp_ohm = 50; % Ohm - входной импеданс.
% 2011 - Преселектор
BWrf_ekv_prf_hz = 20.0e6; % Hz - Эквивалентная полоса пропускания преселектора
Lrf_max_prf_db = 1.0; % dB - Максимальный уровень затухания преселектора в рабочей полосе РПУ (insertion loss)
% 2012 - МШУ
% (для примера взяты параметры Agilent MGA-71543 с учетом реализационных потерь)
Grf_lna_db = 16.0; % dB - Коэффициент усиления МШУ
NFrf_lna_db = 1.0; % dB - Коэффициент шума МШУ
% 2013 - Фильтр основной селекции
BWrf_ekv_fms_hz = 6.0e6; % Hz - Полоса пропускания фильтра основной селекции
Lrf_max_fms_db = 4.0; % dB - Максимальный уровень затухания фильтра основной селекции в полосе пропускания (insertion loss)
% 2014 - Усилитель высокой частоты
% (для примера взяты параметры Agilent MGA-71543 с учетом реализационных потерь)
Grf_amp_db = 16.0; % dB - Коэффициент усиления УВЧ
NFrf_amp_db = 1.0; % dB - Коэффициент шума УВЧ
% 2015 - Квадратурный смеситель
% (для примера взяты параметры Analog Devices ADL5387)
Grf_mix_db = 4.5; % dB - Коэффициент усиления квадратурного смесителя
NFrf_mix_db = 15.0; % dB - Коэффициент шума квадратурного смесителя
IP1dBrf_mix_dbw = 13.0 - 30.0; % dBW - Input P1dB (IP1dB)
% 3. Системная диаграмма цифрового приемника для анализа чувствительности и максимально допустимого уровня блокирующего сигнала
% 3011 - Радиотракт
% Коэффициент шума радиотракта (db)
NFrf_sys_db = pow2db( ( db2pow( Lrf_max_prf_db ) ) + ...
( db2pow( NFrf_lna_db ) + 1 ) / ( db2pow( -Lrf_max_prf_db )) + ...
( db2pow( Lrf_max_fms_db ) + 1 ) / ( db2pow( -Lrf_max_prf_db ) * db2pow( Grf_lna_db )) + ...
( db2pow( NFrf_amp_db ) + 1 ) / ( db2pow( -Lrf_max_prf_db ) * db2pow( Grf_lna_db ) * db2pow( -Lrf_max_fms_db )) + ...
( db2pow( NFrf_mix_db ) + 1 ) / ( db2pow( -Lrf_max_prf_db ) * db2pow( Grf_lna_db ) * db2pow( -Lrf_max_fms_db ) * db2pow( Grf_amp_db )) ...
);
% Коэффициент усиления радиотракта (dB)
Grf_sys_db = ( Grf_lna_db + Grf_amp_db + Grf_mix_db ) - ( Lrf_max_prf_db + Lrf_max_fms_db );
% Полоса сигнала на выходе радиотракта (Hz)
BWrf_sys_hz = BWrf_ekv_fms_hz;
% 3012 - Baseband LPF
Lbb_lpf_db = 9; % dB - Максимальный уровень затухания Baseband LPF в полосе пропускания (insertion loss)
% 3013 - Baseband УНЧ (LTC6400-14)
Gbb_opa_db = 0; % dB - Коэффициент усиления УНЧ
NFbb_opa_db = 0; % dB - Коэффициент шума УНЧ
% Эквивалентный коэффициент усиления УНЧ, состоящего из фильтра и ОУ (dB)
Gbb_lfa_db = Gbb_opa_db - Lbb_lpf_db;
% Коэффициент шума системы от входа до выхода УНЧ (dB) - практически не поменялся.
NFbb_sys_db = pow2db( ( db2pow( NFrf_sys_db ) ) + ...
( db2pow( NFbb_opa_db ) + 1 ) / ( db2pow( Grf_sys_db )) ...
);
% Коэффициент усиления системы от входа до выхода УНЧ (dB)
Gbb_sys_db = Grf_sys_db + Gbb_lfa_db;
% Полоса сигнала в основной полосе (Hz) - не поменялась
BWbb_sys_hz = BWrf_sys_hz;
% Мощность шума на выходе УНЧ в основной полосе (dBW)
PNbb_out_dbw = pow2db( k * T0 * BWrf_sys_hz ) + NFbb_sys_db + Gbb_sys_db;
% 3014 - АЦП
% (для примера взяты параметры 1 канала Linear Technology LTC2271)
FSadc_hz = 20.0e6; % Hz - Sampling rate
SNadc_fs_db = 84; % dB - SNR
NBadc_fs_bits = 16; % bits - Full scale bits
Vadc_fs_v = 2; % V - Full scale voltage
Radc_inp_ohm = 1000; % Ohm - Input ADC resistance
% Мощность синусоидального сигнала при размахе на полную шкалу для пары АЦП (dBW)
PFSadc_inp_dbw = pow2db( 2.0 * (( Vadc_fs_v * 0.5 * sqrt( 0.5 )) ^ 2 ) / Radc_inp_ohm );
% Мощность собственного шума пары АЦП (dBW)
PNadc_snr_dbw = PFSadc_inp_dbw - SNadc_fs_db;
% Мощность шума квантования пары АЦП (dBW)
PNadc_qan_dbw = PFSadc_inp_dbw - ( NBadc_fs_bits * mag2db( 2 ) + mag2db( sqrt( 6 ) / 2 )); % adding correction factor for sinusoidal signal
% 3015 - Эквивалентная схема обработки сигнала с ЧМ
SNfm_min_db = 12.0; % dB - минимальное отношение сигнал-шум для приема голосовых передач в режиме ЧМ
BWfm_max_hz = 25.0e3; % Hz - максимальная ширина полосы канала, предназначенного для передачи речи
BWfm_min_hz = 6.25e3; % Hz - минимальная ширина полосы канала, предназначенного для передачи речи
% Чувствительность радиотракта (dBW) при идеальных АЦП и цифровой обработке канала передачи (минимальная мощность в полосе канала передачи)
Pfm_min_dbw = pow2db( k * T0 * BWfm_max_hz ) + NFbb_sys_db + SNfm_min_db;
% Мощность сигнала на выходе унч при синусоидальном входном сигнале, по уровню равным уровню чувствительности (dBW)
Pfm_min_bb_sys_dbw = Pfm_min_dbw + Gbb_sys_db;
% Мощность сигнала на входе смесителя. Такой сигнал после прохождения всех последующих каскадов
% на выход АЦП даст размах полной шкалы (dBW)
PFSmix_inp_dbw = PFSadc_inp_dbw - Gbb_lfa_db - Grf_mix_db;
% Запас линейного участка передаточной характеристики смесителя (по входу)
deltaPmix_inp_lin = IP1dBrf_mix_dbw - PFSmix_inp_dbw;
% Разрешение БПФ равно минимальной полосе канала
Nfft = 2 * (( FSadc_hz / 2 ) / BWfm_min_hz );
Nsmp = Nfft;
tmp_fft_buf = zeros( 1, Nfft );
tmp_acc_buf = zeros( 1, Nfft );
tmp_smp_buf = zeros( 1, Nsmp );
max_acc = 30;
for acc = 1:max_acc
% Принимаемый сигнал 1
% на входе приемного тракта - сигнал на уровне чувствительности приемного тракта
% на входе АЦП - сигнал прошедший все каскады преобразования и усиления
PS1 = db2pow( Pfm_min_bb_sys_dbw );
WS1 = 25.0e3;
FS1 = 1.0e6;
Fstart = FS1;
Fstop = Fstart + WS1 - BWfm_min_hz;
Pstep = PS1 / ( WS1 / BWfm_min_hz );
PS1_smp_buf = zeros( 1, Nsmp );
for f = Fstart:BWfm_min_hz:Fstop
phi_acc = 2.0 * pi * rand( 1 ); % random phase
phi_stp = pi * ( f / ( FSadc_hz / 2 ));
for k = 1:Nsmp
PS1_smp_buf( k ) = PS1_smp_buf( k ) + sqrt( Pstep ) * exp( j * phi_acc );
phi_acc = phi_acc + phi_stp;
if( phi_acc > ( 2.0 * pi ))
phi_acc = phi_acc - 2.0 * pi;
else
if( phi_acc < ( -2.0 * pi ))
phi_acc = phi_acc + 2.0 * pi;
end
end
end
end
% Принимаемый сигнал 2
% на входе приемного тракта - сигнал по уровню близок блокировке приемного тракта
% на входе АЦП - сигнал прошедший все каскады преобразования и усиления
PS2 = db2pow( PFSadc_inp_dbw - 1.0 ); % -1 dB back off
WS2 = 25.0e3;
FS2 = -2.0e6;
Fstart = FS2;
Fstop = Fstart + WS2 - BWfm_min_hz;
Pstep = PS2;
PS2_smp_buf = zeros( 1, Nsmp );
for f = Fstart:BWfm_min_hz:Fstop
phi_acc = 2.0 * pi * rand( 1 ); % random phase
phi_stp = pi * ( f / ( FSadc_hz / 2 ));
for k = 1:Nsmp
PS2_smp_buf( k ) = PS2_smp_buf( k ) + sqrt( Pstep ) * exp( j * phi_acc );
phi_acc = phi_acc + phi_stp;
if( phi_acc > ( 2.0 * pi ))
phi_acc = phi_acc - 2.0 * pi;
else
if( phi_acc < ( -2.0 * pi ))
phi_acc = phi_acc + 2.0 * pi;
end
end
end
end
% Шум на выходе унч основной полосы
PN1 = db2pow( PNbb_out_dbw );
WN1 = BWbb_sys_hz;
Pfull_bw = PN1 * ( FSadc_hz / WN1 );
PN1_smp_buf = sqrt( 0.5 * Pfull_bw ) * complex( randn( 1, Nsmp ), randn( 1, Nsmp ));
tmp_fft_buf = fftshift( fft( PN1_smp_buf ));
tmp_msk_buf = zeros( 1, Nfft );
tmp_msk_buf((( Nfft / 2 ) - (( WN1 / FSadc_hz ) * ( Nfft / 2 )) + 1 ) : (( Nfft / 2 ) + (( WN1 / FSadc_hz ) * ( Nfft / 2 )))) = ...
ones( 1, (( WN1 / FSadc_hz ) * Nfft ));
tmp_fft_buf = tmp_fft_buf .* tmp_msk_buf;
PN1_smp_buf = ifft( fftshift( tmp_fft_buf ));
% Шум аналоговой части АЦП
PN2 = db2pow( PNadc_snr_dbw ) - db2pow( PNadc_qan_dbw );
%PN2 = db2pow( PNadc_snr_dbw );
Pfull_bw = PN2;
PN2_smp_buf = sqrt( 0.5 * Pfull_bw ) * complex( randn( 1, Nsmp ), randn( 1, Nsmp ));
% Ограничитель и квантователь АЦП
QAN_smp_buf = PS1_smp_buf + PS2_smp_buf + PN1_smp_buf + PN2_smp_buf;
QAN_delta = Vadc_fs_v / ( 2 ^ NBadc_fs_bits );
QAN_smp_buf = round( QAN_smp_buf ./ QAN_delta ) .* QAN_delta;
QAN_smp_buf_re = real( QAN_smp_buf );
QAN_smp_buf_re( find( QAN_smp_buf_re > ( Vadc_fs_v / 2.0 ))) = Vadc_fs_v / 2.0;
QAN_smp_buf_re( find( QAN_smp_buf_re < ( -Vadc_fs_v / 2.0 ))) = -Vadc_fs_v / 2.0;
QAN_smp_buf_im = imag( QAN_smp_buf );
QAN_smp_buf_im( find( QAN_smp_buf_im > ( Vadc_fs_v / 2.0 ))) = Vadc_fs_v / 2.0;
QAN_smp_buf_im( find( QAN_smp_buf_im < ( -Vadc_fs_v / 2.0 ))) = -Vadc_fs_v / 2.0;
QAN_smp_buf = complex( QAN_smp_buf_re, QAN_smp_buf_im );
% Спектр после ADC
tmp_smp_buf = QAN_smp_buf;
%tmp_smp_buf = PS1_smp_buf + PS2_smp_buf + PN1_smp_buf + PN2_smp_buf;
tmp_fft_buf = fft( tmp_smp_buf ) / Nfft;
tmp_acc_buf = tmp_acc_buf + ( tmp_fft_buf .* conj( tmp_fft_buf ));
end
tmp_acc_buf = tmp_acc_buf ./ max_acc;
f = linspace(( -FSadc_hz / 2 ) + BWfm_min_hz, FSadc_hz / 2, Nfft );
plot( f, pow2db( fftshift( tmp_acc_buf )));
xlim( [( -FSadc_hz / 2 ), ( FSadc_hz / 2 )] );
ylim( [-150.0, -20.0] );
title('Power Spectrum')
xlabel('Frequency (Hz)')
ylabel('P(f) dBW')
drawnow;
Что ж, модель показала жизнеспособность идеи — чувствительность -115дБм, блокровка под 90дБ.
Далее, в ПЛИС, с помощью блока нормализации сигналов квадратур, мы уберем постоянную составляющую, справимся с зеркальным каналом и подадим сигнал на вход DDC. После сноса интересующей частоты в ноль, сигнал попадет на цепочку цифровых CIC — и FIR-фильтров, формирующих канальную полосу. Разумеется, если мы хотим одновременно принимать не один канал, у нас должна быть кучка DDC и фильтров.
В следующей статье, если публике будет интересно, расскажу о дальнейших шагах в моделировании, и оценке аппаратных ресурсов ПЛИС.
