BPSK матлаб модулятор демодулятор Опции

[Acvarif]

Имеется код BPSK модема, передатчик и приемник.

CODE

clear all;
Br = 117.1875*10^3;
Fd = Br*8;

Fif = 20*10^6;
Fdif = 60*10^6;
Fdopl = 100;

%% === TRANSMITTER ===
M = pngen(5);
D = rand(1, 200*8) > 0.5;
CRC = zeros(1, 16);
for i=1:length(D),
CRC(16) = xor(CRC(16), D(i));
CRC = [CRC(2:16) CRC(1)];
end;
FILL = ones(1, 200)*0.5;
DATA = [dec2bin(length(D)/8, 11)=='1' D CRC];
packet = [FILL M M diff_enc([0 DATA]) FILL]*2-1;

% redescr
Ndkr = Fd/Br;
pkt = packet(ceil(0.01:1/Ndkr:length(packet)));

hch = firls(31, [0 .24 .3 .5]*2, [1 1 0 0]);
pktf = filter(hch, [1], pkt);

% upsample
t1 = upsample4(pktf, 31);
t2 = upsample4(t1, 31);
pktfu = upsample4(t2, 31)';

% -> IF
l = length(pktfu);
s = sin(2*pi*(Fif+Fdopl)*(0:1/Fdif:l/Fdif));
c = cos(2*pi*(Fif+Fdopl)*(0:1/Fdif:l/Fdif));
pkts = pktfu.*s(1:l);
pktc = pktfu.*c(1:l);

%% === TRANSMIT LINE ===
Ir = awgn(awgn(pkts, 0.05), 0.01);
Qr = awgn(awgn(pktc, 0.05), 0.01);

%% === RECEIVER ===
% -> 0
l = length(Ir);
s = sin(2*pi*Fif*(0:1/Fdif:l/Fdif));
c = cos(2*pi*Fif*(0:1/Fdif:l/Fdif));

Rreu = Ir.*s(1:l) + Qr.*c(1:l);
Rimu = Qr.*s(1:l) - Ir.*c(1:l);

tr00 = downsample4(Rreu);
tr01 = downsample4(tr00);
Rre_ = downsample4(tr01);
Rre = filter(hch, [1], Rre_);

tr10 = downsample4(Rimu);
tr11 = downsample4(tr10);
Rim_ = downsample4(tr11);
Rim = filter(hch, [1], Rim_);

%figure; hold on; plot(Rre); plot(Rim, 'r');

% detection
Ms = M(ceil(0.01:1/Ndkr:length(M)));
cfh = Ms(length(Ms):-1:1);
cre = filter(cfh, [1], Rre)/length(cfh);
cim = filter(cfh, [1], Rim)/length(cfh);
%cr = sqrt(cre.^2 + cim.^2);
cr = (7/8)*max(abs(cre), abs(cim)) + (1/8)*min(abs(cre), abs(cim));

%figure; hold on; plot(cre); plot(cim); plot(cr, 'r');

env(length(cr)) = 0;
for j=1:length(cr),
f = max(1, j-Ndkr*length(M)+2*Ndkr);
t = min(length(cr), max(2, j-Ndkr));
env(j) = sum(cr(f:t))/(t-f);
end;

%figure; hold on; plot(cr, 'r'); plot(env, 'g'); plot(env*5);

%
for i=length(M)*Ndkr+1:length(cr)-1,
c1 = cr(i) > env(i)*5;
% c2 = (cr(i) > cr(i-1)) && (cr(i) > cr(i+1));
c2 = (cr(i-1) > cr(i-2)) && (cr(i-1) > cr(i));
c3 = abs(cr(i)-cr(i-length(M)*Ndkr)) < 4*env(i);
c4 = cr(i-length(M)*Ndkr) > env(i-length(M)*Ndkr)*5;

if c1 && c2 && c3 && c4,
break;
end;
end;
i
re = Rre(i:length(Rre));
im = Rim(i:length(Rim));

%figure; hold on; plot(re); plot(im, 'r');
figure; hold on; plot(cr, 'r'); plot(env, 'g'); plot(env*5); plot((c1&c2&c3&c4)*1.3, 'k');

%
l = length(re);
lb = length(packet);
Yre(lb) = 0;
Yim(lb) = 0;
fre(lb) = 0;
fim(lb) = 0;
are(lb) = 0;
aim(lb) = 0;
mark(l) = 0;
bit(lb) = 0;

k = 1;
gamma = 0.1;
Yre(1) = re(1);
Yim(1) = im(1);

for i = 8:8:length(re)-10,
mark(i) = 1;
are(k) = sum(re(i-6:i-1));
aim(k) = sum(im(i-6:i-1));

% okre(k) = are(k)*Yim(k) - aim(k)*Yre(k);
okim(k) = aim(k)*Yim(k) + are(k)*Yre(k);

bit(k) = sign(okim(k));

fre(k) = are(k)*bit(k);
fim(k) = aim(k)*bit(k);

Yre(k+1) = (1-gamma)*Yre(k) + gamma*fre(k);
Yim(k+1) = (1-gamma)*Yim(k) + gamma*fim(k);

k = k + 1;
end;

bit_dec = diff_dec(sign(bit+1));
bits = bit_dec(2:length(bit_dec));

%figure; hold on; plot(okim); plot(okre, 'r');
%figure; hold on; plot(re); plot(Yre(ceil(0.01:1/Ndkr:length(bit))), 'r');
%figure; hold on; plot(im); plot(Yim(ceil(0.01:1/Ndkr:length(bit))), 'r');

figure; hold on; plot(im); plot(mark, 'k');

%
len = bin2dec(int2str(bits(1:11)));
data = bits(12:12+len*8-1);
rcrc = bits(12+len*8:12+len*8+15);
crc = zeros(1, 16);
for i=1:length(D),
crc(16) = xor(crc(16), data(i));
crc = [crc(2:16) crc(1)];
end;

disp('---');
disp(sprintf('Received packet length: %d bytes', len));
disp(sprintf('Data ok: %d', sum(abs(data-D))==0));
disp(sprintf('CRC ok: %d', sum(abs(crc-rcrc))==0));

При компиляции возникает ошибка:
Undefined function 'pngen' for input arguments of type 'double'.

Error in Modem (line 10)
M = pngen(5);
Я так понимаю, что функция pngen должна быть где-то описана. Или такая функция уже встроена в Матлаб последних версий.
Использую R2014a.
Подскажите пожалуйста в чем проблема.
Если кто сталкивался, подсобите пожалуйста простейшим рабочим кодом BPSK модема (передатчик + приемник).

[Grizzzly]

Нет, это не встроенная функция. Ее написал автор этого кода. Надо искать реализацию.

[Acvarif]

Нет, это не встроенная функция. Ее написал автор этого кода. Надо искать реализацию.

Спасибо. Понятно.

CODE

clear
% несущаЯ
Fc = 2400;
% частота выборки бит из вектора Х
Fd = 2400;
% частота выборки выходного модулированного сигнала
Fs = 19200;
% количество сэмплов на 1 бит
Sb = Fs/Fd;
% количество 0 сэмплов перед сигналом
ZeroS = 10;%round(Sb/0.8);
% количество 0 сэмплов перед сигналом
NoiseLev = 0.2;
% ограничение тракта
Amax = 5;
% Количество наборов фаз на символ
M = 2;
% МодулирующаЯ последовательность
x = repmat([0 1 0 0 0 0 0 1 1 1 1 1 0 1 0 1 0 1 0 0 1 0 0 0],1,1); x=x(;
% Modulate, keeping track of time.
[y,t] = dmod(x,Fc,Fd,[Fs pi/2],'psk',M);
% Всего выборок в выходом сигнале
samples_total = length(y);

if ZeroS > 0
y(ZeroS + 1: samples_total) = y (1 : samples_total - ZeroS);
y(1 : ZeroS) = 0;
end

noise = NoiseLev*randn(size(y));

% yn сигнал на входе демодулЯтора
yn = y + noise;

% Смоделировать выход за динамический диапазон по входу
for sn = 1 : samples_total
if yn(sn) > Amax
yn(sn) = Amax;
elseif yn(sn) < -Amax
yn(sn) = -Amax;
end
end

% ДемодулЯциЯ
% Вектор номеров выборок
samples_vector = 1 : samples_total;
%ГенерациЯ векторов sin и cos локальной несущей длЯ коррелЯции c Y
dFs = 2*pi*Fc/Fs;
% sin и cos на длине одного периода
sinFc = sin(0:dFs:2*pi-dFs)';
cosFc = cos(0:dFs:2*pi-dFs)';
samples_sin = length(sinFc);
% sin и cos на длине вектора входного сигнала
sinY = repmat(sinFc, samples_total/samples_sin, 1);
cosY = repmat(cosFc, samples_total/samples_sin, 1);

% ДвигаЯсь по выборкам скользЯщим окном, равным длине битового интервала, вычислЯем коррелЯцию с несущей
demod_y(1:samples_total) = 0;
demod_I(1:samples_total) = yn.*cosY;
demod_Q(1:samples_total) = yn.*sinY;
demod_y = demod_I + demod_Q*j;

% Усреднение I и Q на длине битового интервала
demod_avr_I(1:samples_total) = 0;
demod_avr_Q(1:samples_total) = 0;
avr_interval = samples_sin;

for sn = avr_interval : samples_total
demod_avr_I(sn) = sum(demod_I( sn - avr_interval + 1: sn ))/avr_interval*2;
demod_avr_Q(sn) = sum(demod_Q( sn - avr_interval + 1: sn ))/avr_interval*2;
end

% Перемножение усредненных I и Q с задержанной на 1 бит версией
demod_dly_I(1:samples_total) = 0;
demod_dly_Q(1:samples_total) = 0;
demod_dly_I(1+Sb:samples_total) = demod_avr_I(1:samples_total - Sb).* demod_avr_I(1+Sb:samples_total);
demod_dly_Q(1+Sb:samples_total) = demod_avr_Q(1:samples_total - Sb).* demod_avr_Q(1+Sb:samples_total);
demod_dly = demod_dly_I + demod_dly_Q;

% ФильтрациЯ демодулированного сигнала
[demod_b,demod_a]=butter(2,(Fd/Fs));
demod_filt = filter(demod_b,demod_a,sign(demod_dly));


subplot(5,1,1), plot (samples_vector, y, 'b', samples_vector, y, 'r'), grid on;
set (gca, 'XLimMode', 'manual',...
'XLim', [1,samples_total]',...
'XTickMode', 'manual',...
'XTick', [Sb:Sb:samples_total]');
subplot(5,1,2), plot (samples_vector, yn), grid on;
set (gca, 'XLimMode', 'manual',...
'XLim', [1,samples_total]',...
'XTickMode', 'manual',...
'XTick', [Sb:Sb:samples_total]');
subplot(5,1,3), plot (samples_vector, demod_I, samples_vector, demod_Q, samples_vector, abs(demod_y)), grid on;
set (gca, 'XLimMode', 'manual',...
'XLim', [1,samples_total]',...
'XTickMode', 'manual',...
'XTick', [Sb:Sb:samples_total]');
subplot(5,1,4), plot (samples_vector, demod_dly_I, samples_vector, demod_dly_Q, samples_vector, demod_dly), grid on;
set (gca, 'XLimMode', 'manual',...
'XLim', [1,samples_total]',...
'XTickMode', 'manual',...
'XTick', [Sb:Sb:samples_total]');
subplot(5,1,5), plot (samples_vector, demod_dly, samples_vector, sign(demod_dly), samples_vector, sign(demod_filt)), grid on;
set (gca, 'XLimMode', 'manual',...
'XLim', [1,samples_total]',...
'XTickMode', 'manual',...
'XTick', [Sb:Sb:samples_total]',...
'YLimMode', 'manual',...
'YLim', [-1.5, 1.5]',...
'YTickMode', 'manual',...
'YTick', [-1.5:0.5:1.5]');

Очевидно в этом коде функция dmod тоже принадлежит автору?
https://www.gaussianwaves.com/digital-modul...s-from-scratch/
Я так понимаю это за деньги. Можно-ли найти подобную литературу не за деньги?

Хотя dmod вроде существует http://matlab.exponenta.ru/communication/book2/6/dmod.php
Но начиная с какой версии Матлаб?
Даже код кое-какой имеется http://read.pudn.com/downloads40/sourcecod...MM/DMOD.M__.htm

Сообщение отредактировал Acvarif - Mar 19 2018, 13:01

[ataradov]

Ух ты. Я написал этот код много лет назад и все эти функции есть у меня на домашнем компе.

Я так предполагаю это из РТИ? Или это код дальше разошелся?

pngen() генерирует М-последовательность длинны 2^n, где n - это параметр функции. Можно просто найти любую М-последовательность длинны 32 и подставить напрямую.

Сообщение отредактировал ataradov - Mar 19 2018, 17:08

[ataradov]

Вот исходник этой pngen():

CODE

function x=pngen(p);
% x=pnegn(p) - Generate pseudo-random sequence
% p = either a scaler (2..12) to generate a PN sequence
% using a predefined polynomial of length p, or
% a binary generating polynomial.
% x = sequence, of length 2^length(p) - 1.

% Predefined seqence polynomials:
pl{2} = [1 1];
pl{3} = [1 1 0];
pl{4} = [1 0 0 1];
pl{5} = [1 0 1 1 1];
pl{6} = [1 1 1 0 0 1];
pl{7} = [1 1 0 0 1 0 1];
pl{8} = [1 0 1 0 1 1 1 1];
pl{9} = [1 0 1 0 1 1 0 1 1];
pl{10} = [1 0 1 0 1 0 1 0 1 1];
pl{11} = [1 0 1 0 1 0 1 0 1 0 1];
pl{12} = [1 1 0 1 0 0 0 0 1 1 0 1];
if length(p) == 1
p = pl{p};
end

L=length(p);
N=2^L-1;

s=ones(1,length(p));
x=zeros(1,N);
for i=1:N
snew = mod(sum(and(p,s)),2);
x(i) = snew;
s=[s(2:end),snew];
end

Если нужна помощь с остальным - обращайтесь, помогу.

[Acvarif]

Если нужна помощь с остальным - обращайтесь, помогу.

Большое спасибо. Помощь конечно-же нужна.
Для начала немного осмыслю как все работает...

[ataradov]

Какая-то из версий этой модели была на Verilog-е реализована, но исходников этого мне сейчас уже точно не найти.

Ну и не стоит эту модель воспринимать как полную истину. Это моя первая модель модема созданная по следам чтения Скляра и экспериментов.

[Acvarif]

Какая-то из версий этой модели была на Verilog-е реализована, но исходников этого мне сейчас уже точно не найти.

Ну и не стоит эту модель воспринимать как полную истину. Это моя первая модель модема созданная по следам чтения Скляра и экспериментов.

Да, я понимаю.
У меня пока совсем никакого представления. Только в общих чертах. Как работает ФАПЧ мне известно. На практике приходилось это делать, но только в бинарном варианте.
Но сейчас это все реализовывается в дискретном варианте на микроконтроллерах http://www.ti.com/lit/an/slaa681a/slaa681a.pdf, ПЛИС. Для начала пытаюсь хоть как-то понять это дело через Матлаб. Конечно Verilog (VHDL) код очень бы не помешал. Поскольку больше работаю с ПЛИС.

В коде появилась еще ошибка Undefined function 'diff_enc' for input arguments of type 'double'.
Похоже должна быть определена еще и функция diff_enc.

[ataradov]

Похоже должна быть определена еще и функция diff_enc.

CODE

function enc=diff_enc(data, start);
    if ~(exist('start'))
        start = 0;

    de(1)=xor(data(1), start);
    for i=2:length(data), de(i)=xor(de(i-1), data(i)); end;
    enc = de;
end

Этот модем скорее всего перебор по сравнению с той программной реализацией на МК.

Сообщение отредактировал ataradov - Mar 20 2018, 07:18

[Acvarif]

Этот модем скорее всего перебор по сравнению с той программной реализацией на МК.

Спасибо. Дело продвинулось.
Далее Undefined function 'upsample4' for input arguments of type 'double'.
Похоже дальше есть еще неопределенные функции...

Перебор это в смысле..?

Сообщение отредактировал Acvarif - Mar 20 2018, 07:29

[petrov]

Для начала пытаюсь хоть как-то понять это дело через Матлаб. Конечно Verilog (VHDL) код очень бы не помешал. Поскольку больше работаю с ПЛИС.

Лучше сразу смотреть в строну симулинка, облегчает понимание и постепенно модель можно конретизировать до полноценного рабочего HDL кода не прибегая к внешним HDL симуляторам.

[ataradov]

Похоже дальше есть еще неопределенные функции...

Вот все функции какие у меня есть.

Перебор это в смысле..?

Я возможно не совсем понял статью, но они там похоже уже ожидают прямоугольный сигнал и только восстанавливают битовый клок.

Этот модем еще борется с несовпадением частоты дискретизации на приемнике и передатчике (или эффектом Допплера). Плюс тут много всяких заморочек с конкретным железом. Для модема на 100 Кбод не нужна дискретизация на 60 МГц.

Сообщение отредактировал ataradov - Mar 20 2018, 07:37

Прикрепленные файлы

 functions.zip ( 8.04 килобайт ) Кол-во скачиваний: 0

 

[Acvarif]

Лучше сразу смотреть в строну симулинка, облегчает понимание и постепенно модель можно конретизировать до полноценного рабочего HDL кода не прибегая к внешним HDL симуляторам.

Да. Согласен. Но в симуленке я сильно плаваю.

Вот все функции какие у меня есть.

Архив не скачивается. Выскакивает на пустую страницу.

Я возможно не совсем понял статью, но они там похоже уже ожидают прямоугольный сигнал и только восстанавливают битовый клок.

Насколько я понял там задействован внутренний АЦП микроконтроллера. Далее цифровые библиотеки FIR и т. п.

Этот модем еще борется с несовпадением частоты дискретизации на приемнике и передатчике (или эффектом Допплера). Плюс тут много всяких заморочек с конкретным железом. Для модема на 100 Кбод не нужна дискретизация на 60 МГц.

Доплер, это мне какраз и нужно.

[ataradov]

Да. Согласен. Но в симуленке я сильно плаваю.

Может быть проще найти время и разобраться. Будет гораздо проще в итоге.

Архив не скачивается. Выскакивает на пустую страницу.

Добавил снова.

Насколько я понял там задействован внутренний АЦП микроконтроллера. Далее цифровые библиотеки FIR и т. п.

Одного АЦП в общем случае не достаточно нужно 2 квадратуры цифровать.

Прикрепленные файлы

 functions.zip ( 8.04 килобайт ) Кол-во скачиваний: 12

 

[Acvarif]

Может быть проще найти время и разобраться. Будет гораздо проще в итоге.

Согласен. Придется разбираться.
Спасибо за функции. Главный код заработал. Попытаюсь разобраться.

Одного АЦП в общем случае не достаточно нужно 2 квадратуры цифровать.

Микроконтроллер имеет на борту многоканальный АЦП. Ну а если даже использовать один канал то при 4_х выборках на период можно использовать соседние выборки как синус и косинус.

[petrov]

Да. Согласен. Но в симуленке я сильно плаваю.

Получается как в анекдоте:
Мужик что-то ищет под фонарем. Тут к нему под ходит милиционер и
спрашивает: "Что вы тут делаете?" Мужик отвечает: "Ключи от квартиры
ищу". "А где потерял?". "В парке". "А зачем здесь ищешь?".
"А здесь светлее ".

Искать надо там где надо. Все эти обрывки чужих исходников - тупик. Симулинк - отличнейшее окружение, чтобы разобраться с сутью, да ещё в конце и HDL код получить можно, готовый для прошивки в FPGA.

[Acvarif]

Получается как в анекдоте:
Мужик что-то ищет под фонарем. Тут к нему под ходит милиционер и
спрашивает: "Что вы тут делаете?" Мужик отвечает: "Ключи от квартиры
ищу". "А где потерял?". "В парке". "А зачем здесь ищешь?".
"А здесь светлее ".

Искать надо там где надо. Все эти обрывки чужих исходников - тупик. Симулинк - отличнейшее окружение, чтобы разобраться с сутью, да ещё в конце и HDL код получить можно, готовый для прошивки в FPGA.

Принимается.
Кстати Ваш симулинк код под названием bpsk_costas_2008_08_25.mdl надеюсь в этом поможет.

Сообщение отредактировал Acvarif - Mar 20 2018, 08:15

[petrov]

Кстати Ваш симулинк код под названием bpsk_costas_2008_08_25.mdl надеюсь в этом поможет.

Возможно, лишь как отправная точка, воспринимайте критически, ищите пробелы, недостатки и ошибки.

[Acvarif]

Возможно, лишь как отправная точка, воспринимайте критически, ищите пробелы, недостатки и ошибки.

Посмотрел модельку. Сложно для начинающего. Не понимаю, что такое скрамблер.
Google советует туториал по цифровой модуляции http://www.cengage.com/resource_uploads/do...2511_315335.pdf
Только не понятно где эти все примеры находятся. В самом Матлаб, или нужно искать в сети?

[ataradov]

Посмотрел модельку. Сложно для начинающего. Не понимаю, что такое скрамблер.

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

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

[Acvarif]

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

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

Спасибо. Понятно. Тогда по ходу в демодуляторе должен быть дескрамблер?
Хотя на этом думаю мне пока не стоит зацикливаться.
В простом случае принцип фазовой манипуляции ничего особенного их себя не представляет. Это все очень похоже, например, на принцип формирования кода Баркера, с чем плотно приходилось и приходится иметь дело в цифровой технике.
Что касается демодулятора с ФАПЧ (или модное название петля Костаса) это все тоже очень похоже на цифровую ФАПЧ с использованием синтезатора частоты, например, по Уолшу в обратной связи, что тоже приходилось делать, и работало оно на ура.
Думаю не стоит мне сейчас влезать в дебри цифровой связи. Это все придет по ходу.
Пока в Матлабе ищу простейший способ запустить простейший модулятор-демодулятор (кодированием или симулинк) чтобы понять как все перенести в ПЛИС. Имеется ввиду в цифровом виде (АЦП, выборки и пр.). В бинарном виде ФАПЧ ничего особенного из себя не представляет. Разве, что только пока не понимаю как подбираются промежуточная по отношению к несущей, как пересчитывать коэффициенты для управления синтезатором уолша. В свое время это все делалось на варикапах и подбиралось экспериментально, потом по мере появления новых технологий все переносилось в ПЛИС, но соблюдая все наработки.

Сообщение отредактировал Acvarif - Mar 22 2018, 07:22

[ataradov]

Спасибо. Понятно. Тогда по ходу в демодуляторе должен быть дескрамблер?

Да, конечно.

В простом случае принцип фазовой манипуляции ничего особенного их себя не представляет. Это все очень похоже, например, на принцип формирования кода Баркера, с чем плотно приходилось и приходится иметь дело в цифровой технике.

Это каша полная. Причем тут Баркер?

Что касается демодулятора с ФАПЧ (или модное название петля Костаса)

Петля Костаса - это не модное название ФАПЧ, это очень специфическая схема с ФАПЧ. ФАПЧ бывает и с другими схемами, ничего общего с Костасом не имеющие.

Думаю не стоит мне сейчас влезать в дебри цифровой связи. Это все придет по ходу.

Как раз таки наоборот, с пониманием цифровой связи почти автоматически придет понимание как это сделать в матлабе.

[Acvarif]

Это каша полная. Причем тут Баркер?

Баркер тоже фазовая манипуляция. Применялась для формирования сложного сигнала.
BPSK это тоже фазовая манипуляция.

Петля Костаса - это не модное название ФАПЧ, это очень специфическая схема с ФАПЧ. ФАПЧ бывает и с другими схемами, ничего общего с Костасом не имеющие.

Не согласен. Ничего специфического там нет. В бинарном виде эта самая ФАПЧ делается с применением синтезатора частоты Уолша как гетеродина.
Используем, например, в гидроакустических лагах.

Как раз таки наоборот, с пониманием цифровой связи почти автоматически придет понимание как это сделать в матлабе.

Согласен. Только просто чтение, по опыту, бесполезно. Пытаюсь все делать по ходу. Может и не правильно, но так уж сложилось...
Конечно хотелось-бы выяснить как все-же правильно, подобрать частоту гетеродина, например при несущей 100 кГ, если один бит при BPSK модуляции планируется передавать четырьмя периодами несущей.
Как рассчитывается шаг изменения частоты гетеродина в петле ФАПЧ.
Должны-же быть основополагающие формулы для таких расчетов...

Сообщение отредактировал Acvarif - Mar 22 2018, 08:10

[petrov]

Только просто чтение, по опыту, бесполезно. Пытаюсь все делать по ходу.

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

Конечно хотелось-бы выяснить как все-же правильно, подобрать частоту гетеродина, например при несущей 100 кГ, если один бит при BPSK модуляции планируется передавать четырьмя периодами несущей.
Как рассчитывается шаг изменения частоты гетеродина в петле ФАПЧ.
Должны-же быть основополагающие формулы для таких расчетов...

Демодуляцию делают на нулевой частоте, частота комплексного гетеродина минус 100 кГц.

[Acvarif]

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

Пробую. Трудно дается. Неспецифическая тема.

Демодуляцию делают на нулевой частоте, частота комплексного гетеродина минус 100 кГц.

Не совсем врубаюсь.
Комплексная это понятно.
Допустим несущая 100 кГц. Какая должна быть средняя частота управляемого синтезатора частоты DDS в петле ФАПЧ?
Тоже 100кГ?
Если взять одно квадратурное плечо-умножитель то на один вход будет подаваться модулированная несущая 100 кГц, а на другой, например sin 100 кГц. В этом случае какраз и получится нулевая чатота. Так?
Если так почему именно на нулевой, а не со сдвигом, например 100 и 102?
Мы когда-то ФАПЧ делали с небольшой разностной, например 7 и 5 кГц.
Или я чего-то не так понял?

[petrov]

Пробую. Трудно дается. Неспецифическая тема.

Не совсем врубаюсь.
Комплексная это понятно.
Допустим несущая 100 кГц. Какая должна быть средняя частота управляемого синтезатора частоты DDS в петле ФАПЧ?
Тоже 100кГ?
Если взять одно квадратурное плечо-умножитель то на один вход будет подаваться модулированная несущая 100 кГц, а на другой, например sin 100 кГц. В этом случае какраз и получится нулевая чатота. Так?
Если так почему именно на нулевой, а не со сдвигом, например 100 и 102?
Мы когда-то ФАПЧ делали с небольшой разностной, например 7 и 5 кГц.
Или я чего-то не так понял?

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

[Acvarif]

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

Понял. Спасибо.
Тогда если в цифровом виде (не бинарном).
Несущая 100 кгц. Частота дискретизации 400 кГц. Исходная частота NCO в петле ФАПЧ 100 кГц.
Количество периодов несущей на символ (в смысле на бит) 4.
От этого можно оттолкнуться?

Сообщение отредактировал Acvarif - Mar 22 2018, 11:53

[petrov]

Понял. Спасибо.
Тогда если в цифровом виде (не бинарном).
Несущая 100 кгц. Частота дискретизации 400 кГц. Исходная частота NCO в петле ФАПЧ 100 кГц.
Количество периодов несущей на символ (в смысле на бит) 4.
От этого можно оттолкнуться?

Нет никакой петли ФАПЧ 100 кГц, нет никакого количества периодов несущей на символ, просто переносим нужную частоту на нулевую.

[Acvarif]

Нет никакой петли ФАПЧ 100 кГц, нет никакого количества периодов несущей на символ, просто переносим нужную частоту на нулевую.

Понятно. Тогда в простейшем виде и малыми шагами.
Модулятор
Выход элемента Switch - BPSK сигнал
Далее, чтобы его демодулировать необходимо выделить из него комплексную огибающую. Для этого его нужно умножить на комплексную экспоненту (перенос спектра модулированного сигнала в область 0) и затем отфильтровать с помощью ФНЧ
Выглядит это так
Если, допустим, в модуляторе несущая равна Fn, то в демодуляторе чему должна быть равна омега0, тоже Fn?
Вот на этом моменте я почему-то застрял...

[petrov]

[Acvarif]

Показан процесс демодуляции - повторное умножение на комплексную экспоненту со знаком минус.
Грубо говоря это избавление от несущей. Так?
Значит если несущая 100 кгц то в демодуляторе омега0 тоже должна быть 100 кГц. Так?
Тоесть если в квадратурном гетеродине омега0 сделать 110 кГц то спектр раздвоится на 10 и 210 кГц.
Произойдет перенос модулированного сигнала на другую частоту. А нам, грубо говоря, это н. не нужно.
Нужно всего лишь избавиться от несущей. Так?

[petrov]

При умножении на комплексную экспоненту ничего не раздваивается, просто спектр сдвигается. Если мы берём реальную часть комплексного сигнала, то получаем симметричный спектр относительно нуля.

[Acvarif]

При умножении на комплексную экспоненту ничего не раздваивается, просто спектр сдвигается. Если мы берём реальную часть комплексного сигнала, то получаем симметричный спектр относительно нуля.

Да. Погорячился. Сдвигается.
Но всеравно остается вопрос почему в квадратурном демодуляторе(гетеродине) на входы умножителей обычно подается частота равная несущей?
Для того чтобы сразу от нее избавиться?
Почему не перенести спектр на другую частоту Например несущая 100 кГц Гетеродин 98 кГц. Так у нас обычно делалось. Дальше уже работа идет с 2 кГц и там уже и выделяется комплексная огибающая таким же способом перемножения комплексных експонент.

[petrov]

Но всеравно остается вопрос почему в квадратурном демодуляторе(гетеродине) на входы умножителей обычно подается частота равная несущей?

Вы видели созвездие BPSK, QPSK, QAM? По которым мы принимаем решения о передаваемых символах. Чтобы такую картинку получить нужно осуществить фазовую и символьную синхронизации. Если будет сдвиг по частоте, то созвездие будет вращаться с угловой частотой соответствующей сдвигу и решения принять не получится.

[Acvarif]

Вы видели созвездие BPSK, QPSK, QAM? По которым мы принимаем решения о передаваемых символах. Чтобы такую картинку получить нужно осуществить фазовую и символьную синхронизации. Если будет сдвиг по частоте, то созвездие будет вращаться с угловой частотой соответствующей сдвигу и решения принять не получится.

Спасибо. Кажется становится понятно. Нужно еще осмыслить...

[Acvarif]

Для большего понимания для себя происходящего в BPSK сворганил простую модель.

CODE

%% BPSK модель
% Количество бит данных
N = 4;
% Данные
data = [1 0 1 0];
%data = randi(1,N);
pData = data*2 - 1;
% Несущая частота (100 кГц)
f= 100000;
% частота выборок
fs = f*100;
% Период выборки
Ts = 1/fs;
% Период несущей
T = 1/f;
% Количество периодов несущей на бит-символ
M = 4;
% Длина пакета данных
n = M*length(data);
% Текущее время
t = 0:Ts:n*T;
% Несущая в пакете
car = sin(2*pi*f*t);
% График данных
subplot(2, 1, 1);
stem(pData);
title('Данные');
% График несущей
subplot(2, 1, 2);
plot(car);
title('Несущая');

%% Преобразование данных в прямоугольные импульсы
tp = 0:Ts:T*M;
exdata = [];
for (i = 1:length(data))
for(j = 1:length(tp) - 1)
exdata = [exdata pData(i)];
end
end
exdata = [exdata 0];
figure;
plot(exdata, 'r-', 'LineWidth', 4);
hold on;
grid on;
plot(car, 'g-');
hold on

%% Модуляция
% Перенос exdata на несущую
mSig = exdata.*car;
plot(mSig, 'h-', 'LineWidth', 2);
hold off;
title('Прям. импульсы + несущая + промод.несущая');

%% AWGN (аддитивный белый гауссов шум) канал передачи данных
figure;
% отношение сигнал/шум в децибелах
SNR = 50;
% мощность сигнала в децибелах
SIGPOWER = 10;
% смесь мод.сигнала с шумом
rx = awgn( mSig, SNR, SIGPOWER);
plot(mSig, 'r-', 'LineWidth', 3);
hold on;
plot(rx, 'g-', 'LineWidth', 1);
grid on;
%% Демодуляция
% умножение мод.зашумленого сигн.на несущую
dem = rx.*car;
figure;
plot(dem);
grid on;
%% Декодирование
k = 1;
rcv = [];
for(i = 1:length(data))
sm = 0;
for(j = 1:length(tp) - 1)
sm = sm + dem(k);
k = k + 1;
end
if(sm > 0)
rcv = [rcv 1];
else
rcv = [rcv 0];
end
end
clc;
% Переданные данные
data
% Принятые данные
rcv

В данном случае это идеальный, полностью синхронизированный, вариант.
Я так понимаю в реальности все далеко не так. Демодулятор может иметь свою частоту дискретизации отличную от модулятора и она не синхронна с частотой дискретизации в модуляторе. Но это не самое важное. Важно что несущая подаваемая на умножитель в демодуляторе тоже не синхронна с несущей пришедшего на него сигнала, не говоря уже о доплере если источник и приемник двигаются друг относительно друга.
В этой части вопрос. Какая основная задача петли Костаса? Подстройка частоты подаваемой на умножитель, вернее на умножители демодулятора? Или задача петли Костаса это фазовая и символьная синхронизация в демодуляторе? Если это так то не врубаюсь как подстройкой частоты NCO можно добиться еще и фазовой синхронизации...

Сообщение отредактировал Acvarif - Mar 26 2018, 05:03

[ataradov]

Костас только убирает несовпадение частот и допплер (они не различимы на приемнике). Костас ничего не делает для символьной синхронизации. Для этого в данные вставляют синхропоследовательности для начальной синхронизации и дальнейшая синхронизация осуществляется по смене фаз дынных. И чтобы этих смен было достаточно используют скремблеры.

Сообщение отредактировал ataradov - Mar 26 2018, 05:10

[Acvarif]

Костас только убирает несовпадение частот и допплер (они не различимы на приемнике). Костас ничего не делает для символьной синхронизации. Для этого в данные вставляют синхропоследовательности для начальной синхронизации и дальнейшая синхронизация осуществляется по смене фаз дынных. И чтобы этих смен было достаточно используют скремблеры.

Спасибо. Понятно. Стало более понятно значение нескольких байт преамбулы обычно используемой в пакете данных.
Попытаюсь встроить петлю Костаса в код
На демодулятор поступает rx (смесь модулированного сигнала с аддит. гаусовым шумом) которую далее необходимо умножать на
sin и cos поступающие с выхода VCO
Очевидно что в демодуляторе необходимо будет организовать
1. Управляемый генератор sin cos (VCO).
2. Два одинаковых ФНЧ
3. Фазовый дискриминатор
4. Еще фильтр перед управляемым генератором
Какова роль фильтра (loopfilter) после фазового дискриминатора?

[petrov]

Попытаюсь встроить петлю Костаса в код

В коде сложнее что-либо разглядеть, чем в симулинк модели.

2. Два одинаковых ФНЧ

Согласованный фильтр лучше поставить до петли, иначе задержка фильтра будет входить в петлю ФАПЧ, с соответствующим увеличением времени настройки, для обеспечения устойчивасти.

4. Еще фильтр перед управляемым генератором
Какова роль фильтра (loopfilter) после фазового дискриминатора?

Пропорционально интегрирующее звено ФАПЧ.

Ещё что-то надо делать с амплитудой сигнала, которую надо учитывать при рассчёте коэффициентов передачи ФАПЧ.

[Acvarif]

Согласованный фильтр лучше поставить до петли, иначе задержка фильтра будет входить в петлю ФАПЧ, с соответствующим увеличением времени настройки, для обеспечения устойчивасти.

До петли это как? На схеме разве не так? Фильтры стоят сразу после I Q умножителей.
Разве можно как-то иначе поставить фильтры и так чтобы не учитывалась их задерка в петле ФАПЧ?
Хотя Ваша схема ФАПЧ не похожа на привычную петлю Костаса
Если не сложно, поясните пожалуйста в нескольких фразах как она работает.

Ещё что-то надо делать с амплитудой сигнала, которую надо учитывать при рассчёте коэффициентов передачи ФАПЧ.

С этим я вообще пока не знаю что делать...
В свое время когда петля ФАПЧ делалась в бинарном виде этот вопрос вообще не стоял.

Сообщение отредактировал Acvarif - Mar 26 2018, 11:20

[petrov]

До петли это как? На схеме разве не так? Фильтры стоят сразу после I Q умножителей.

Уже выше обсуждали, что комплексный сигнал у вас будет на нулевой частоте, сначала фильтр согласованный, а потом уже ФАПЧ - умножение на управляемую комплексную экспоненту, поворачиваем вращающееся созвездие в обратную сторону, чтобы оно стояло неподвижно.

Хотя Ваша схема ФАПЧ не похожа на привычную петлю Костаса
Если не сложно, поясните пожалуйста в нескольких фразах как она работает.

Это и не петля Костаса, в оригинальной статье петля работает без символьной синхронизации, здесь же она работает после символьной синхронизации в домене одного отсчёта на символ. Согласованный фильтр здесь как раз включён в петлю.

С этим я вообще пока не знаю что делать...
В свое время когда петля ФАПЧ делалась в бинарном виде этот вопрос вообще не стоял.

Моделировать, моделировать. У вас же в общем случае короткие пакеты будут приходить минимального и максимального уровня друг за другом.

[Acvarif]

Уже выше обсуждали, что комплексный сигнал у вас будет на нулевой частоте, сначала фильтр согласованный, а потом уже ФАПЧ - умножение на управляемую комплексную экспоненту, поворачиваем вращающееся созвездие в обратную сторону, чтобы оно стояло неподвижно.

Запутался.
Вот петля Костаса
Входной сигнал m(t) sin(ωt) Далее он сразу умножается в одном плече на sin(ωt) в другом на cos(ωt) Как результат - на выходе верхнего плеча компонента I, нижнего Q. Произошел снос на нулевую частоту. Что понимается под согласованным фильтром? У нас, например, под согласованным фильтром понимается коррелятор при работе со сложным сигналом.

Сообщение отредактировал Acvarif - Mar 26 2018, 13:10

[petrov]

Запутался.
Вот петля Костаса
Входной сигнал m(t) sin(ωt) Далее он сразу умножается в одном плече на sin(ωt) в другом на cos(ωt) Как результат - на выходе верхнего плеча компонента I, нижнего Q. Произошел снос на нулевую частоту. Что понимается под согласованным фильтром? У нас, например, под согласованным фильтром понимается коррелятор при работе со сложным сигналом.

Статья Костаса древняя, эту петлю в аналоге делали, не нужно её тупо переносить в цифру.
Разобрались, что мы можем перенести сигнал на нулевую частоту.
Далее фильтруем его согласованным фильтром, согласованным с формой видео импульсов, которые мы передаём.
Коррелятор фильтром не является, он не инвариантен ко времени прихода.
Далее петлёй Костаса на нулевой частоте разворачиваем созвездие BPSK, чтобы оно стояло.
Filter1,3 не нужны, 2 - ПИ звено. VCO - NCO, комплексная экспонента управляемая, умножается на комплексный сигнал с выхода согласованного фильтра. Детектор ошибки - произведение реальной и мнимой частей после умножителя.

[Acvarif]

Спасибо что Вы мне подсказываете. Общее понимание почти сложилось. Но непонятые детали не дают покоя.

Статья Костаса древняя, эту петлю в аналоге делали, не нужно её тупо переносить в цифру.

Да. Картинка немного не та. Вот та картинка
Но суть вроде не меняется.

Разобрались, что мы можем перенести сигнал на нулевую частоту.

Железно.
Перенос на нулевую частоту это умножения входного сигнала на комплексную экспоненту.

Далее фильтруем его согласованным фильтром, согласованным с формой видео импульсов, которые мы передаём.

В цифровом варианте это будет FIR (LPF) Это понимается под согласованным фильтром?

Filter1,3 не нужны

Если по цифровой схеме Костаса то не нужны. Там стоят LPF.

VCO - NCO, комплексная экспонента управляемая, умножается на комплексный сигнал с выхода согласованного фильтра.

Тут опять тупик. Получается что согласованный фильтр должен стоять до умножителей.
Может в качестве согласованного фильтра имеется ввиду это?

Сообщение отредактировал Acvarif - Mar 26 2018, 19:11

[petrov]

Чем отличается согласованный фильтр от ФНЧ?

Тем, что максимально возможный сигнал/шум даёт на выходе. Читайте про приподнятый косинус и корень из приподнятого косинуса.

А как быть с компонентами на удвоенной несущей после переноса на нулевую?
Вроде ФНЧ их устраняет.

Да устраняет.

Тут опять тупик. Получается что
1. Согласованный фильтр должен стоять до умножителей.
2. ФНЧ не нужны.

Да, всё уже отфильтровано, осталось только устранить вращение созвездия.

Поймите, что от переноса частоты в ноль и даунсемплинга нужно абстрагироваться, с этим вопросом можно разобраться отдельно, всё равно никакого Костаса у вас не будет в итоге, потому что ФАПЧ медленно настраивается. Далее вся обработка осуществляется на нулевой частоте, в комплексном виде, с малым количеством отсчётов на символьный интервал.

[Acvarif]

Далее вся обработка осуществляется на нулевой частоте, в комплексном виде, с малым количеством отсчётов на символьный интервал.

Да. Где-то я это понимаю. В частности про малое количество отсчетов.
Но тем не менее для моделирования должна сложиться четкая картинка компонентов и их соединения. Схема Костаса эту картинку дает.

всё равно никакого Костаса у вас не будет в итоге, потому что ФАПЧ медленно настраивается.

Как-же тогда быть? Никакой другой схемы кроме Костаса пока не встречал...

Сообщение отредактировал Acvarif - Mar 27 2018, 04:44

[petrov]

Как-же тогда быть? Никакой другой схемы кроме Костаса пока не встречал...

Можно почитать для начала RF Architectures and Digital Signal Processing Aspects of Digital Wireless Transceivers - Nezami.

[Acvarif]

Если посчитать за какое время сработает петля Костаса, например, при доплере 150 Гц (в воде это движение ~2 м/сек) (несущая 100 кГц) то получится, что для полной синхронизации потребуется время равное 150 символам (битам). Это если шаг изменения частоты NCO будет 1 Гц.
150 бит это ~18 байт преамбулы. Если длительность одного символа 4 периода несущей то это 40 мкс. Всего 150*40 = 6000 мкс. Тоесть полная подстройка произойдет за 6 мс. Это вполне нормально. Но в реальности наверняка шаг подстройки будет меньше чем 1 Гц, например 0.5 Гц. Ну тогда полная синхронизация наступит через 12 мс. Тоже вроде неплохо. Далее пойдет синхросимвол и данные.
Для случая с несущей в 100 кГц вроде неплохо.

Можно почитать для начала RF Architectures and Digital Signal Processing Aspects of Digital Wireless Transceivers - Nezami.

Уже читаю. Глава 5.
Думаю для начала всеравно нужно построить петлю Костаса. Приступать сразу к FeedForward наверное неправильно. Ведь так?

[ataradov]

VCO получает сигал ошибки в качестве контрольного, зачем перестраивать фиксированным шагом? Вам и предлагают ПИ звено поставить чтобы шаг сам нашелся какой нужно для оптимальной настройки.

[Acvarif]

VCO получает сигал ошибки в качестве контрольного, зачем перестраивать фиксированным шагом? Вам и предлагают ПИ звено поставить чтобы шаг сам нашелся какой нужно для оптимальной настройки.

Об этом постоянно думаю. В смысле о самонастраиваемом шаге.
Что такое ПИ звено?

Сообщение отредактировал Acvarif - Mar 27 2018, 07:48

[ataradov]

Что такое ПИ звено?

https://ru.wikipedia.org/wiki/%D0%9F%D0%98%...%82%D0%BE%D1%80

ПИ это тоже самое,только без "Д".

[petrov]

Уже читаю. Глава 5.
Думаю для начала всеравно нужно построить петлю Костаса.

Это ещё месяц назад надо было сделать, в симулинке изучить все возможные воздействия и сочетания параметров, непонятно чего тянете, стройте уже. Такими темпами и 5 лет не хватит, на то, что на самом деле надо сделать для вашей задачи.

[Acvarif]

Это ещё месяц назад надо было сделать, в симулинке изучить все возможные воздействия и сочетания параметров, непонятно чего тянете, стройте уже. Такими темпами и 5 лет не хватит, на то, что на самом деле надо сделать для вашей задачи.

Ok. Попробую. В простейшем виде.
Сформировал BPSK сигнал

Я так понимаю сигнал необходимо оцифровать для того, чтобы потом умножать на сигнал формируемый Discrete-Time VCO.
Сделал примерно так

Ход верный или ерунда?
На FPGA я делал-бы именно так...

[petrov]

Ok. Попробую. В простейшем виде.
Сформировал BPSK сигнал

Я так понимаю сигнал необходимо оцифровать для того, чтобы потом умножать на сигнал формируемый Discrete-Time VCO.
Сделал примерно так

Ход верный или ерунда?
На FPGA я делал-бы именно так...

Не нужно моделировать сигнал на несущей, этим вы только замедляете модель, увеличивая количество отсчётов, необходимое для представления сигнала на несущей. Ни о какой оцифровке на данном этапе думать не нужно, и про FPGA тоже. Работайте с комплексными сигналами. В общем примеры я выкладывал, смотрите как там сделано.

[Acvarif]

Не нужно моделировать сигнал на несущей, этим вы только замедляете модель, увеличивая количество отсчётов, необходимое для представления сигнала на несущей. Ни о какой оцифровке на данном этапе думать не нужно, и про FPGA тоже. Работайте с комплексными сигналами. В общем примеры я выкладывал, смотрите как там сделано.

Рад-бы. Но мой уровень симулинк не позволяет сходу это сделать.
Ваш пример для меня сложен. А примеры имеющиеся в сети на мой взгляд бестолковые. Их даже в качестве рыбы нельзя использовать. Оторванные от реального железа.
Поскольку в симулинк я вообще не умею работать то он вызывает некоторое отторжение.
Петлю ФАПЧ попробую сделать в классическом матлабе.

CODE

%% BPSK модель
% Количество бит данных
N = 4;
% Данные
data = [1 0 1 0];
%data = randi(1,N);
pData = data*2 - 1;
% Несущая частота (100 кГц)
fn= 100000;
% частота выборок
fns = fn*10;
% Период выборки
Tns = 1/fns;
% Период несущей
Tn = 1/fn;
% Количество периодов несущей на бит-символ
M = 4;
% Длина пакета данных
n = M*length(data);
% Текущее время
tn = 0:Tns:n*Tn;
% Несущая в пакете
car = sin(2*pi*fn*tn);

%% Преобразование данных в прямоугольные импульсы
tpn = 0:Tns:Tn*M;
exdata = [];
for (i = 1:length(data))
for(j = 1:length(tpn) - 1)
exdata = [exdata pData(i)];
end
end
exdata = [exdata 0];

%% Модуляция
% Перенос exdata на несущую
mSig = exdata.*car;

Теория говорит, что mSig по сути комплексный сигнал. Но в демодуляторе нужно использовать только реальную его часть, которую затем в петле ФАПЧ нужно умножать на комплексную экспоненту.
Какой выглядит в коде реальная часть mSig?

[ataradov]

Нужно взять исходный бинарный сигнал, умножить его на синус и косинус частоты Допплера без всяких несущих. Все, получили I и Q. Это сигнал идущий на схему вращающую созвездие. Костас или не Костас - это как захочется.

Чтобы визуально увидеть созвездие нужно построить график I(Q). Для нулевого Допплера и в отсутствии шума и прочих искажений это будут 2 точки. При добавлении Допплера это будет окружность. Задача ФАПЧ вернуть эти точки на исходное место и следить, чтобы они там и остались.

Сообщение отредактировал ataradov - Mar 28 2018, 06:50

[Acvarif]

Нужно взять исходный бинарный сигнал, умножить его на синус и косинус частоты Допплера без всяких несущих. Все, получили I и Q. Это сигнал идущий на схему вращающую созвездие. Костас или не Костас - это как захочется.

Чтобы визуально увидеть созвездие нужно построить график I(Q). Для нулевого Допплера и в отсутствии шума и прочих искажений это будут 2 точки. При добавлении Допплера это будет окружность. Задача ФАПЧ вернуть эти точки на исходное место и следить, чтобы они там и остались.

Понял. Спасибо.

[Acvarif]

Не получается.
Пытаюсь сформировать синфазную и квадратурную составляющую модулирующего сигнала.

CODE

%% BPSK модель
% Количество бит данных
N = 4;
% Данные
data = [1 0 1 0];
%data = randi(1,N);
pData = data*2 - 1;
% Несущая частота (100 кГц)
fn = 100000;
% Доплер
fd = 100;
% частота выборок
fns = fn*8;
% Период выборки
Tns = 1/fns;
% Период несущей
Tn = 1/fn;
% Количество периодов несущей на бит-символ
M = 4;
% Длина пакета данных
n = M*length(data);
% Текущее время
tn = 0:Tns:n*Tn;
% Несущая в пакете
car_sin = sin(2*pi*fn*tn);
car_cos = cos(2*pi*fn*tn);

%% Преобразование данных в прямоугольные импульсы
tpn_i = 0:Tns:Tn*M;
exdata_i = [];
for (i = 1:length(data))
for(j = 1:length(tpn_i) - 1)
exdata_i = [exdata_i pData(i)];
end
end
exdata_i = [exdata_i 0];

tpn_q = 0:Tns:Tn*M;
exdata_q = [];
for (i = 1:length(data))
for(j = 1:length(tpn_q) - 1)
% exdata_q = [exdata_q pData(i)];
if (length(tpn_q) == 17)
exdata_q = [exdata_q 1];
else
exdata_q = [exdata_q -1];
end
end

end
exdata_q = [exdata_q 0];

%% Модуляция
% Перенос exdata на несущую
mSig_i = exdata_i.*car_cos;
mSig_q = exdata_q.*car_sin;

%%

figure;
plot(exdata_i, 'r-', 'LineWidth', 4);
hold on;
grid on;
plot(exdata_q, 'b-', 'LineWidth', 4);
hold on;
plot(car_sin, 'g-');
hold on
hold off;

Получается ерунда.
У Вас похоже это работает

CODE

clear all;
Br = 117.1875*10^3;
Fd = Br*8;

Fif = 20*10^6;
Fdif = 60*10^6;
Fdopl = 100;

%% === TRANSMITTER ===
M = pngen(5);
%M = 0.5;
D = rand(1, 200*8) > 0.5;
%D = [1 0 1 0];
CRC = zeros(1, 16);
for i=1:length(D),
CRC(16) = xor(CRC(16), D(i));
CRC = [CRC(2:16) CRC(1)];
end;
FILL = ones(1, 200)*0.5;
DATA = [dec2bin(length(D)/8, 11)=='1' D CRC];
packet = [FILL M M diff_enc([0 DATA]) FILL]*2-1;

% redescr
Ndkr = Fd/Br;
pkt = packet(ceil(0.01:1/Ndkr:length(packet)));

hch = firls(31, [0 .24 .3 .5]*2, [1 1 0 0]);
pktf = filter(hch, [1], pkt);

% upsample
t1 = upsample4(pktf, 31);
t2 = upsample4(t1, 31);
pktfu = upsample4(t2, 31)';

% -> IF
l = length(pktfu);
s = sin(2*pi*(Fif+Fdopl)*(0:1/Fdif:l/Fdif));
c = cos(2*pi*(Fif+Fdopl)*(0:1/Fdif:l/Fdif));
pkts = pktfu.*s(1:l);
pktc = pktfu.*c(1:l);

%% === TRANSMIT LINE ===
Ir = awgn(awgn(pkts, 0.05), 0.01);
Qr = awgn(awgn(pktc, 0.05), 0.01);

%% === GRAPHICS ===

figure;
plot(pkts, 'g-');
grid on;
hold on
plot(pktc, 'r-');
hold on

figure;
plot(Ir, 'g-');
grid on;
hold on
plot(Qr, 'r-');
hold on
plot(pktfu, 'r-');
hold on

figure;
plot(pkt, 'g-');
grid on;
hold on

clc;
% Переданные данные
%pkt

Но там слишком большое количество данных. Сложно рассмотреть подробности.

[ataradov]

Для отладки синхронизации не нужна несущая, ее можно добавить в самом конце.

Вот как будут выглядеть I и Q компоненты идеального сигнала после всех переносов частоты.

CODE

clear all;
Br = 400e3;
Fd = Br*8;

Fdopl = 100;

D = (rand(1, 10000) > 0.5) * 2 - 1;
% redescr
Ndkr = Fd/Br;
pkt = D(ceil(0.01:1/Ndkr:length(D)));
l = length(pkt);

s = sin(2*pi*Fdopl*(0:1/Fd:l/Fd));
c = cos(2*pi*Fdopl*(0:1/Fd:l/Fd));

I = s(1:l) .* pkt;
Q = c(1:l) .* pkt;

figure; hold on; plot(I); plot(Q, 'r');
figure; hold on; plot(I, Q, '*');

Начинайте с этого кода и добавляйте схемы, которые восстановят этот сигнал.

Если запустить с Fdopl = 0, то на второй картинке будет 2 звездочки. Это идеальный конечный результат. Если запустить Fdopl = 100, то будет окружность, образованная вращением этих звездочек.

Задача ФАПЧ вернуть звездочки на место и поддерживать их там.

Как только заработает с идеальным сигналом, то нужно будет обрезать ему полосу и добавить шум. И отлаживать это. И только после этого можно смотреть на несущие, если сильно хочется.

Сообщение отредактировал ataradov - Mar 30 2018, 05:31

[Acvarif]

Для отладки синхронизации не нужна несущая, ее можно добавить в самом конце.
Вот как будут выглядеть I и Q компоненты идеального сигнала после всех переносов частоты.

Спасибо. Попробую разобраться.

[Acvarif]

Сворганил BPSK приемо-передатчик с петлей Костаса на базе VCO.
Основной код петли заимствовал на просторах сети.

CODE

%% Петля Костаса на базе VCO

%% Модуляция
clear;

% Битовая картика профиля pacman
image = [ 1, 1, 1, 0, 0, 1, 1, 1;
1, 1, 0, 0, 0, 0, 1, 1;
1, 0, 0, 0, 1, 0, 0, 1;
0, 0, 0, 0, 0, 0, 1, 1;
0, 0, 0, 0, 0, 1, 1, 1;
1, 0, 0, 0, 0, 0, 0, 1;
1, 1, 0, 0, 0, 0, 1, 1;
1, 1, 1, 0, 0, 1, 1, 1 ];

% Начальный символ приема
sync_symbol = [ 1 ];
% Битовый поток с начальным символом
bit_stream = [ sync_symbol ];

% Генерация битового потока из матрицы профиля pacman
for j = 1:8
bit_stream = [bit_stream image(j,];
end

%bit_stream = (rand(1, 6) > 0.5);
%bit_stream = [1 1 0 1 0 1 0 1 0 0 1 0 1 0 1 0];
bit_stream = [1 1 0 1 0];
N = length(bit_stream);
fprintf('bit_stream %d\n', bit_stream);

% Частота семплирования 288 кГц
fs = 288000;
% Несущая частота относительно частоты выборок 36 кГц
f = 0.125;
% Частота доплера 300 Гц
fdop = f/1200;
% Символьная частота 100 Гц
fsim = 100;
% Количество выборок на период несущей 8
ns = 1/f;
% Количество выборок необх. для одного символа (360*8 периодов нес. f)
% Символьная скорость 100 Гц
tsymbol = fs*f*ns/fsim;
fprintf('tsymbol %d\n', tsymbol);
% Счетчик выборок в 0
t = 0;
% Количество символов переданных с нулевой фазой
nDummySymbols = 4;
% Общее количество выборок
%ttotal = N*tsymbol + nDummySymbols * tsymbol + tsymbol*3;
ttotal = N*tsymbol + nDummySymbols * tsymbol;
fprintf('ttotal %d\n', ttotal);
% Создать пустой сигнал
BPSK_signal = zeros(ttotal,1);

% Символы с нулевой фазой
for i = 1:nDummySymbols
for j = 1:tsymbol
BPSK_signal(t+1) = cos( 2*pi*(f+fdop)*t );
BPSK_signal_nomod(t+1) = cos(2*pi*(f+fdop)*t);
t = t + 1;
end
end

fprintf('t %d\n', t);

% BPSK сигнал
for i = 1:length(bit_stream)
phase = pi * bit_stream(i);
for j = 1:tsymbol
BPSK_signal(t+1) = cos(2*pi*(f+fdop)*t + phase);
BPSK_signal_mod(t+1) = cos(2*pi*(f+fdop)*t + phase);
BPSK_signal_nomod(t+1) = 0;
t = t + 1;
end
end

fprintf('length BPSK_signal %d\n', length(BPSK_signal));
fprintf('t %d\n', t);

%% AWGN (аддитивный белый гауссов шум) канал передачи данных
% отношение сигнал/шум в децибелах
SNR = 20;
% мощность сигнала в децибелах
SIGPOWER = 10;
% смесь мод.сигнала с шумом
BPSK_signal = awgn( BPSK_signal, SNR, SIGPOWER);

%% Демодуляция

% полученная картинка в 0
rec_image = zeros(8,8);

% Координаты картинки
x = 1;
y = 1;

% начальный шаг в 0
t = 0;

% Все матрицы в 0
carrier_inph = zeros(ttotal,1);
mixer_inph = zeros(ttotal,1);
demod_thr_inph = zeros(ttotal,1);

% Фильтр демодулятора
b = fir1(100, 0.1);
% Память состояния демодулятора
stateLowpassInph = zeros(numel(-1,1);

% переменные для PLL осциллятора
% cosine output
c = 1;
% delayed cosine by one timestep
c_delay = 0;
% sine output
s = 0;
% delayed sine by one timestep
s_delay = 0;

% Данные из VCO для целей отладки
sine = zeros(ttotal,1);
cosine =zeros(ttotal,1);
vco = zeros(ttotal,1);

% PLL loop фильтр
bPLL = fir1(40, 0.1);
% and the state memory
zfPLLa = zeros(numel(bPLL)-1,1);
zfPLLc = zeros(numel(bPLL)-1,1);
zfPLLs = zeros(numel(bPLL)-1,1);

% Обратный отсчет. При нулевом значении отбирается символ.
% Начальный обратный отсчет учитывает задержку
% fir lowpass и пропускает начальный символ
sampleMoment = tsymbol/2+tsymbol;

% Это время необходимое PLL для стабилизации на несущей
nPLLsettle = (nDummySymbols/2) * tsymbol;

% Игнорировать данные, пока не получен символ "1"
ignoreData = 1;

% VCO шаг
VCOgain = 0.005;

% VCO, которое управляет частотой. при v = 0 он точно равен f.
v = 0;

% Основной цикл
for step = 1:ttotal

% Получить значение BPSK_signal
sample = BPSK_signal(step);

% PLL - "voltage" управляемый осциллятор с +/- v входом
% w0 частота VCO
w0 = 2*pi*(f+v*VCOgain);
% Cохранить предыдущее состояние
c_delay = c;
s_delay = s;
% Подсчет нового состояния
c = c_delay * cos(w0) - s_delay * sin(w0);
s = s_delay * cos(w0) + c_delay * sin(w0);
% Cохраним все в удобных векторах для графики
sine(step) = s;
cosine(step) = c;
vco(step) = v;
% конец VCO

% Сохраняем sine/cosine. Обратите внимание, что инфазная волна
% теперь является синусоидальной волной, потому что PLL с мультипликацией
% в качестве фазового детектора имеет сдвиг фазы на 90 градусов между
% входом и выходом
carrier_inph(step) = c;

% Это решает, находимся ли мы в режиме синхронизации PLL или режима замораживания
% Примечание: это должно быть заменено блокировкой PLL в реальном приеме,
% потому что мы не знаем, когда передача фактически
% начинается. Здесь мы просто предполагаем, что она начинается с самого начала

% PLL in action: изменение частоты через v
% Фазовый детектор
pc = sample*c;
ps = sample*-s;

% filtering out the 2f term
% Фильтрация pc
[vc,zfPLLc] = filter(bPLL,1,pc,zfPLLc);

% filtering out the 2f term
% Фильтрация ps
[vs,zfPLLs] = filter(bPLL,1,ps,zfPLLs);

[v,zfPLLa] = filter(bPLL,1,vc*vs,zfPLLa);

% we sneak here a -1 in because this is ambigous
% шаг на -1
mixer_inph(step) = -sample*carrier_inph(step);

[lp_demod_inph(step),stateLowpassInph] = filter(b,1,mixer_inph(step),stateLowpassInph);

demod_thr_inph(step) = lp_demod_inph(step);

if (ignoreData)
% Ждем стартового символа 1
if ( lp_demod_inph(step) > 0.25 )
% got it!
% Стартовый символ прошел
ignoreData = 0;
end
else
% Ожидание середины символа для выборки демодулятора
sampleMoment = sampleMoment - 1;
if (sampleMoment == 0)
sampleMoment = tsymbol;
if ( demod_thr_inph(step) > 0.2 )
lsb = 1;
else
lsb = 0;
end

fprintf('lsb %d\n', lsb);

% debug samples: add small spikes
% отладка: добавка небольших пиков
lp_demod_inph(step) = lp_demod_inph(step) + 0.1;

% Сохранить результат в картинку
if (((x<9) && (y<9)) && (ignoreData == 0))
rec_image(y,x) = lsb;
x = x + 1;
if (x > 8)
x = 1;
y = y + 1;
end % if (x > 8)
end % x>9 && y>9
end % sampleMoment == 0
end % ignore data
% Следующий шаг в выборке
t = t + 1;
end

%% Графика

% Немодулированный сигнал
figure;
plot(BPSK_signal_nomod, 'b-', 'LineWidth', 2);
hold on;
grid on;
title('Немодулированный сигнал');

% Модулированный сигнал
figure;
plot(BPSK_signal_mod, 'b-', 'LineWidth', 2);
hold on;
grid on;
title('Модулированный сигнал');

% Полный BPSK_signal
figure;
plot(BPSK_signal, 'b-', 'LineWidth', 2);
hold on;
grid on;
title('BPSK сигнал');

% Демодулированный сигнал
figure;
subplot(2,1,1);
plot(lp_demod_inph,'LineWidth',2);
xlabel('Samples');
ylabel('Amplitude');
title('demod / lowpass filtered');
grid on;

% Картинка
subplot(2,1,2);
imshow(rec_image);

% Напряжение VCO
figure;
subplot(1,1,1);
plot(vco);

Работает при соотнощении сигнал-шум 20 на 10 децибел, доплере 30 Гц при несущей в 36 кГц и символьной скорости 100 Гц
Напрягает то, что при соотношении сигнал шум 1 приемник перестает работать. Все разваливается, даже при нулевом доплере.
Все так и должно быть или нет?
Получается что для нормальной работы BPSK приемника необходим хороший приемный усилитель с узкой полосой (определяется удвоенной символьной скоростью) и коэффициентом усиления достаточным, чтобы при минимально возможном сигнале на его входе он мог обеспечить необходимый уровень сигнала на входе BPSK приемника.
Но с другой стороны при максимально возможном уровне сигнала на входе приемного усилителя он не должен входить в ограничение выхода.
Не представляю как решается такая задачка.

Сообщение отредактировал Acvarif - Apr 13 2018, 06:29

[Acvarif]

Ниже код BPSK модулятор-демодулятор с петлей Костаса на базе VCO для акустического модема.
Просьба покритиковать. На базе данного кода предполагается писать код для ПЛИС либо микроконтроллера, только с использованием NCO.
Характеристики:
- Символьная скорость 100 бит в сек
- Несущая частота 36 000 Гц
- 8 выборок на период несущей
- Устойчиво работает при соотношении сигнал-шум 1
- Максимальный Доплер 60 Гц (в воде ~ 2.5 м/сек)
Есть сомнения в необходимости входного полосового фильтра поскольку наверняка некоторую селективность будет иметь аналоговый приемный тракт.

CODE

%% Петля Костаса на базе VCO

%% Модуляция
clear;

% Битовая картика профиля pacman
% image = [ 1, 1, 1, 0, 0, 1, 1, 1;
% 1, 1, 0, 0, 0, 0, 1, 1;
% 1, 0, 0, 0, 1, 0, 0, 1;
% 0, 0, 0, 0, 0, 0, 1, 1;
% 0, 0, 0, 0, 0, 1, 1, 1;
% 1, 0, 0, 0, 0, 0, 0, 1;
% 1, 1, 0, 0, 0, 0, 1, 1;
% 1, 1, 1, 0, 0, 1, 1, 1 ];

image = [ 1, 0, 1, 0, 1, 0, 1, 0;
0, 1, 0, 1, 0, 1, 0, 1;
1, 0, 1, 0, 1, 0, 1, 0;
0, 1, 0, 1, 0, 1, 0, 1;
1, 1, 1, 1, 1, 1, 1, 1;
0, 0, 0, 0, 0, 0, 0, 0;
1, 1, 1, 1, 1, 1, 1, 1;
0, 0, 0, 0, 0, 0, 0, 0 ];

% Начальный символ приема
sync_symbol = [ 1 ];
% Битовый поток с начальным символом
bit_stream = [ sync_symbol ];

% Генерация битового потока из матрицы профиля pacman
for j = 1:8
bit_stream = [bit_stream image(j,];
end

%bit_stream = (rand(1, 20) > 0.5);
%bit_stream = [1 1 0 1 0 1 0 1 0 0 1 0 1 0 1 0 1 1 0 1 0 1 0 1 0];
%bit_stream = [1 1 0 1 0 1 0 1 0];
N = length(bit_stream);
fprintf('bit_stream %d\n', bit_stream);

% Частота семплирования 288 кГц (8 выборок на период)
fs = 288000;
% Несущая частота относительно частоты выборок 36 кГц
f = 0.125;
% Частота доплера 60 Гц (~2.5 м.сек)
fdop = f/600;

% % Частота семплирования 144 кГц (4 выборки на период)
% fs = 144000;
% % Несущая частота относительно частоты выборок 36 кГц
% f = 0.25;
% % Частота доплера 60 Гц (~2.5 м/сек)
% fdop = -f/600;

% Частота среза НЧ фильтра в демодуляторе
fc = 100;
% Частота Найквиста
fn = fs/2;
% Порядок НЧ фильтра демодулятора
lowpfiltord = 50;
% Порядок loop фильтра демодулятора
loopfiltord = 50;

% Амплитуда
A = 1.0;
% Символьная скорость 100 Гц
fsim = 100;
% Количество выборок на период несущей
ns = 1/f;
% Количество выборок необх. для одного символа (360*8 периодов нес. f)
% при символьной скорости fsim
tsymbol = fs*f*ns/fsim;
fprintf('tsymbol %d\n', tsymbol);
% Счетчик выборок в 0
t = 0;
% Количество символов переданных с нулевой фазой (пилот)
nDummySymbols = 4;
% Общее количество выборок
%ttotal = N*tsymbol + nDummySymbols * tsymbol + tsymbol*3;
ttotal = N*tsymbol + nDummySymbols * tsymbol;
fprintf('ttotal %d\n', ttotal);
% Создать пустой BPSK сигнал
BPSK_signal_t = zeros(ttotal,1);

% Пилот сигнал
for i = 1:nDummySymbols
for j = 1:tsymbol
BPSK_signal_t(t+1) = A*sin( 2*pi*(f+fdop)*t );
BPSK_signal_nomod(t+1) = A*sin(2*pi*(f+fdop)*t);
t = t + 1;
end
end

fprintf('t %d\n', t);

% BPSK сигнал
for i = 1:length(bit_stream)
phase = pi * bit_stream(i);
for j = 1:tsymbol
BPSK_signal_t(t+1) = A*sin(2*pi*(f+fdop)*t + phase);
BPSK_signal_mod(t+1) = A*sin(2*pi*(f+fdop)*t + phase);
BPSK_signal_nomod(t+1) = 0;
t = t + 1;
end
end

fprintf('length BPSK_signal %d\n', length(BPSK_signal_t));
fprintf('t %d\n', t);

%% AWGN (аддитивный белый гауссов шум) + полосовая фильтрация

% отношение сигнал/шум в децибелах
SNR = 9;
% мощность сигнала в децибелах
SIGPOWER = 12;
% смесь мод.сигнала с шумом
BPSK_signal_awgn = awgn( BPSK_signal_t, SNR, SIGPOWER );

% Порядок полосового фильтра
bandpfiltord = 50;
% нижняя частота полосового фильтра 35КГц
W1=35e3/(0.5*fs);
% верхняя частота полосового фильтра 37КГц
W2=37e3/(0.5*fs);
% расчет окна Ханна
k=hann(bandpfiltord+1);
% расчет коэффициентов входного полосового фильтра
bp=fir1(bandpfiltord, [W1 W2], k);
% Память состояния фильтра
stateBandpassBpsk = zeros(numel(bp)-1,1);
% Полосовая фильтрация BPSK сигнала + AWGN шум
[BPSK_signal_awgn_filtr,stateBandpassBpsk] = filter(bp,1,BPSK_signal_awgn,stateBandpassBpsk);
% BPSK_signal_awgn_filtr = filter(bp, 1, BPSK_signal_awgn);

% BPSK_signal
%BPSK_signal = BPSK_signal_awgn;
BPSK_signal = BPSK_signal_awgn_filtr;


%% Демодуляция

% полученная картинка в 0
rec_image = zeros(8,8);

% Координаты картинки
x = 1;
y = 1;

% начальный шаг в 0
t = 0;

% Все матрицы в 0
carrier_inph = zeros(ttotal,1);
mixer_inph = zeros(ttotal,1);
demod_thr_inph = zeros(ttotal,1);

%расчет окна Кайзера
kb=kaiser(lowpfiltord+1);
% Коэффициенты НЧ фильтра демодулятора порядок lowpfiltord
b = fir1(lowpfiltord, fc/fn, kb);
% Память состояния для НЧ фильтра демодул.
stateLowpassInph = zeros(numel(-1,1);

% переменные для PLL осциллятора
% cosine output
c = 1;
% delayed cosine by one timestep
c_delay = 0;
% sine output
s = 0;
% delayed sine by one timestep
s_delay = 0;

% Данные из VCO для целей отладки
sine = zeros(ttotal,1);
cosine =zeros(ttotal,1);
vco = zeros(ttotal,1);

%расчет окна Ханна
kp=hann(loopfiltord+1);
% Коэффициенты PLL loop фильтра
bPLL = fir1(loopfiltord, fc/fn, kp);
% Обнулить память
zfPLLa = zeros(numel(bPLL)-1,1);
zfPLLc = zeros(numel(bPLL)-1,1);
zfPLLs = zeros(numel(bPLL)-1,1);

% Обратный отсчет. При нулевом значении отбирается символ.
% Начальный обратный отсчет учитывает задержку
% fir lowpass и пропускает начальный символ
sampleMoment = tsymbol/2+tsymbol;

% Это время необходимое PLL для стабилизации на несущей
nPLLsettle = (nDummySymbols/2) * tsymbol;

% Игнорировать данные, пока не получен символ "1" (синхробит)
ignoreData = 1;

% VCO шаг
VCOgain = 0.005;

% VCO, которое управляет частотой. при v = 0 он точно равен f.
v = 0;

% Основной цикл
for step = 1:ttotal

% Полосовая фильтрация BPSK сигнала + AWGN шум
% [BPSK_signal(step),stateBandpassBpsk] = filter(bp,1,BPSK_signal_awgn(step),stateBandpassBpsk);

% Получить значение BPSK_signal в sample
sample = BPSK_signal(step);

% PLL - управляемый осциллятор с +/- v входом
% w0 частота VCO
w0 = 2*pi*(f+v*VCOgain);
% Cохранить предыдущее состояние
c_delay = c;
s_delay = s;
% Подсчет нового состояния
c = c_delay * cos(w0) - s_delay * sin(w0);
s = s_delay * cos(w0) + c_delay * sin(w0);
% Cохраним все в удобных векторах для графики
sine(step) = s;
cosine(step) = c;
vco(step) = v;
% конец VCO

% Сохраняем sine/cosine. Обратите внимание, что инфазная волна
% теперь является синусоидальной волной, потому что PLL с мультипликацией
% в качестве фазового детектора имеет сдвиг фазы на 90 градусов между
% входом и выходом
carrier_inph(step) = c;

% Это решает, находимся ли мы в режиме синхронизации PLL или режиме замораживания
% Примечание: это должно быть заменено блокировкой PLL в реальном приеме,
% потому что мы не знаем, когда передача фактически
% начинается. Здесь мы просто предполагаем, что она начинается с самого начала

% PLL изменение частоты через v
% Фазовый детектор
pc = sample*c;
ps = sample*-s;

% Фильтрация pc
[vc,zfPLLc] = filter(bPLL,1,pc,zfPLLc);

% Фильтрация ps
[vs,zfPLLs] = filter(bPLL,1,ps,zfPLLs);

[v,zfPLLa] = filter(bPLL,1,vc*vs,zfPLLa);

% шаг на -1
mixer_inph(step) = sample*carrier_inph(step);

[lp_demod_inph(step),stateLowpassInph] = filter(b,1,mixer_inph(step),stateLowpassInph);

demod_thr_inph(step) = lp_demod_inph(step);

if (ignoreData)
% Ждем стартового символа 1
if ( lp_demod_inph(step) > 0.25 )
% Стартовый символ прошел
ignoreData = 0;
end
else
% Ожидание середины символа для выборки демодулятора
sampleMoment = sampleMoment - 1;
if (sampleMoment == 0)
sampleMoment = tsymbol;
if ( demod_thr_inph(step) > 0.1 )
lsb = 1;
else
lsb = 0;
end

fprintf('lsb %d\n', lsb);

% отладка: добавка небольших пиков
lp_demod_inph(step) = lp_demod_inph(step) + 0.1;

% Сохранить результат в картинку 8:8
if (((x<9) && (y<9)) && (ignoreData == 0))
rec_image(y,x) = lsb;
x = x + 1;
if (x > 8)
x = 1;
y = y + 1;
end % if (x > 8)
end % x>9 && y>9
end % sampleMoment == 0
end % ignoreData
% Следующий шаг в выборке
t = t + 1;
end

%% Графика

% % Немодулированный сигнал
% figure;
% subplot(2,1,1);
% plot(BPSK_signal_nomod, 'b-', 'LineWidth', 2);
% hold on;
% grid on;
% title('Немодулированный сигнал (пилот)');
%
% % Модулированный сигнал
% %figure;
% subplot(2,1,2);
% plot(BPSK_signal_mod, 'r-', 'LineWidth', 2);
% hold on;
% grid on;
% title('Модулированный сигнал (данные)');

% Полный BPSK_signal
figure;
subplot(2,1,1);
plot(BPSK_signal_awgn, 'b-', 'LineWidth', 2);
xlabel('Выборки');
hold on;
grid on;
title('Полный BPSK сигнал + шум AWGN');

subplot(2,1,2);
plot(BPSK_signal, 'r-', 'LineWidth', 2);
xlabel('Выборки');
hold on;
grid on;
title('Полный BPSK сигнал + шум AWGN + Фильтр');

% Демодулированный сигнал
figure;
subplot(2,1,1);
plot(lp_demod_inph,'LineWidth',2);
xlabel('Выборки');
ylabel('Амплитуда');
title('Демодулированный сигнал / НЧ фильтрация');
grid on;

% Картинка
subplot(2,1,2);
imshow(rec_image);
title('Битовая картинка 8:8 квадратиков 1-белый, 0-черный кв. ' );

% Напряжение VCO
% figure
% plot(vco);
% title('Напряжение VCO');
%
% % Фильтр демодулятора
% figure
% % АЧХ НЧ фильтра
% [h,f]=freqz(b,1);
% % Параметры являются соответственно частотой и амплитудой
% plot(f*fs/(2*pi),20*log10(abs(h)))
% xlabel('частота/Hz');ylabel('амплитуда/dB');title('АЧХ НЧ фильтра демодулятора');
%
% % Фильтр полосовой на 36 кГц
% figure
% [h,fp]=freqz(bp,1);%amplitude-frequency characteristic graph
% plot(fp*fs/(2*pi),20*log10(abs(h)))%parameters are respectively frequency and amplitude
% xlabel('частота/Hz');ylabel('амплитуда/dB');title('АЧХ полосового фильтра');

Сообщение отредактировал Acvarif - Apr 16 2018, 07:57