Перейти к содержанию
    

Синхронизация Park

Добрый день!

 

Решил я поразбираться с синхрой по преамбуле(тренировочным символам) реализовал метод Парка, не знаю правильно или нет замоделил, подскажите если что не так.

Хочу померить MSE оценок по данной преамбуле, как это можно правильно сделать?

 

clear all; clc; close all;

SNR=-5:10;
N_FFT=1024;
NGI=1/4;
N_used=N_FFT/2;
N_GI=round(NGI*N_FFT);
Nofdm=N_FFT+N_GI;
Nsym = 500;

nSTOs=50;
CFO=0.4;

%% Формирование тренировочных символов
B=2;
reS=(-1).^round(rand(1,round(N_used/B)));
imS=1i.*(-1).^round(rand(1,round(N_used/B)));
s1=reS+imS;
s2=fliplr(s1);
tx_signal_noGI=[s1 s2 conj(s1) conj(s2)];
tx_signal = zeros(1,Nsym*N_FFT);
tx_signal_GI=[tx_signal_noGI(1,N_FFT-N_GI+1:N_FFT) tx_signal_noGI];
tx_signal=repmat(tx_signal_GI,1,Nsym);

%% Добавление смещения частоты и времени
if nSTOs>=0, y_STO=[tx_signal(1,nSTOs+1:end) zeros(1,nSTOs)]; 
else
   y_STO=[zeros(1,-nSTOs) tx_signal(1,1:end+nSTOs)];
end
nn=0:length(y_STO)-1; y_CFO_STO = y_STO.*exp(1i*2*pi*CFO*nn/N_FFT);

%%
buff = [];
for k = 1:length(SNR)
recvd_signal=awgn(y_CFO_STO, SNR(k), 'measured', 'dB');
   for n = 1:Nsym
   %%
   Len_all=N_FFT; 
   recvd_signal_zeropad=[zeros(1,Len_all) recvd_signal(1,(n-1)*Nofdm+1:n*Nofdm) zeros(1,Len_all) zeros(1,Len_all)];
       for d=1:length(recvd_signal_zeropad)-Len_all-1
           P_Park(d)=sum(recvd_signal_zeropad(d+(1:N_FFT/2)).*fliplr(recvd_signal_zeropad(d+N_FFT/2+(1:N_FFT/2)))); %
           E_Park(d)=sum(abs(recvd_signal_zeropad(d+(1:N_FFT))).^2);
           Q_Park(d)=sum(conj(recvd_signal_zeropad(d+(1:N_FFT/2))).*fliplr(recvd_signal_zeropad(d+N_FFT/2+(1:N_FFT/2))));
       end
   Park_metric=(abs(P_Park).^2)./(E_Park).^2;
   buff = [buff Park_metric];
   plot(Park_metric/max(Park_metric))

   [~, maxipos1]=max(Park_metric-N_GI);
   t_est_Park(1,n)=Nofdm - maxipos1;
   FreqOffset(1,n) = angle(Q_Park(maxipos1-N_GI+1))/pi;
   end
mse_CFO(k) = sqrt(sum(abs(FreqOffset - CFO).^2)/Nsym)
mse_STO(k) = sqrt(sum(abs(t_est_Park - nSTOs).^2)/Nsym)
end

Изменено пользователем Ivan55

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

Если быть точнее то я щас делаю вот, так

 

clear all; clc; close all;

flTimeShift = 1; % 1 - с уходом времени; 0 - без ухода времени
flChann = 0; % 1 - с каналом; 0 - без канала

%choose SNR level in dB
SNR=-5:10;
%choose FFT length
N_FFT=1024;
%choose guard interval length as a percentage of N_FFT
NGI=1/4;
%choose number of used carriers
N_used=N_FFT/2;
N_GI=round(NGI*N_FFT);
Nofdm=N_FFT+N_GI;
Nsym = 500;
%frequency offset
nSTOs=50;
CFO=0.4;

Fs = 11.2e6;         % Частота дискретизации сигнала
Tofdm = Nofdm/Fs;   % Длительность OFDM сигнала в сек
ppm = 50*10^-6;     % Стабильность тактового гененратора
Nppm = Fs*ppm;      % Уход в отсчетах в сек

%% Канал 1
% PathDelays = [0 50e-9 110e-9 170e-9 290e-9 310e-9]; % задержки лучей
% U = [0 -3 -10 -18 -26 -32];
% Fd= [4.6 23 139]; % Доплеровское смещение в Гц
%% Канал 2
% PathDelays = [0 110e-9 190e-9 410e-9]; % задержки лучей
% U = [0 -9.7 -19.2 -22.8];
% Fd= [4.6 23 139]; % Доплеровское смещение в Гц
%% Канал 3
PathDelays = [0 310e-9 710e-9 1090e-9 1730e-9 2510e-9]; % задержки лучей
U = [0 -1.0 -9.0 -10.0 -15.0 -20.0];
Fd= [4.6 23 139]; % Доплеровское смещение в Гц

%% preamble Minn
B=2;
reS=(-1).^round(rand(1,round(N_used/B)));
imS=1i.*(-1).^round(rand(1,round(N_used/B)));
s1=reS+imS;
s2=fliplr(s1);
tx_signal_noGI=[s1 s2 conj(s1) conj(s2)];
tx_signal = zeros(1,Nsym*N_FFT);
tx_signal_GI=[tx_signal_noGI(1,N_FFT-N_GI+1:N_FFT) tx_signal_noGI];
tx_signal=repmat(tx_signal_GI,1,Nsym);

%% Модель канала
OutChan = complex(zeros(Nsym*N_FFT,1),zeros(Nsym*N_FFT,1));
rand('seed',0);
if flChann == 1
   tic
   %T = [0 PathDelays(OptChan)];
   %U = [0 0];
   T = PathDelays;
   chan = rayleighchan(1/Fs,Fd(3)/2,T,U);
   chan.DopplerSpectrum=doppler.ajakes; %flat ajakes
   chan.ResetBeforeFiltering=0;
   OutChan=filter(chan,tx_signal);
   t = toc;
   fprintf('Сигнал искажен каналом...\n Время выполнения %4.2f секунд\n', t);
else
   OutChan = tx_signal;
end

%% Добавление смещение частоты и времени
tic;
delay(1,1) = Nppm*Tofdm;
L = Nofdm;
for n = 2:Nsym
   delay(1,n) = delay(1,n-1) + Nppm*Tofdm; 
   int_delay = fix(delay);
   fract_delay = delay - int_delay;
end

MaxDelay    = fix(max(int_delay));

OutChan(1,end:end+MaxDelay) = 0;
buff_delay = complex(zeros(1, 2*L), zeros(1, 2*L));
if flTimeShift == 1
   for n = 1:Nsym
       buff_delay = circshift(buff_delay, [0 -(L+int_delay(1,n))]);
       buff_delay(1,L-int_delay(1,n)+1:end) = OutChan(1,(n-1)*L+1:n*L+int_delay(1,n));

       re_buff_delay = sig_delay(real(buff_delay),1/Fs,fract_delay(1,n)/Fs);
       im_buff_delay = sig_delay(imag(buff_delay),1/Fs,fract_delay(1,n)/Fs);

       buff_delay = re_buff_delay + 1i*im_buff_delay;

       y_STO(1,(n-1)*L+1:n*L) = buff_delay(1,L+1:2*L);
   end
else
   y_STO = OutChan;
end

if nSTOs>=0, y_STO=[y_STO(1,nSTOs+1:end) zeros(1,nSTOs)]; 
else
   y_STO=[zeros(1,-nSTOs) y_STO(1,1:end+nSTOs)];
end

nn=0:length(y_STO)-1; y_CFO_STO = y_STO.*exp(1i*2*pi*CFO*nn/N_FFT);
t = toc;
fprintf('Осуществлен сдвиг времени и частоты сигнала...\n Время выполнения %4.2f секунд\n', t);

%%
buff = [];
for k = 1:length(SNR)
recvd_signal=awgn(y_CFO_STO, SNR(k), 'measured', 'dB');

for n = 1:Nsym
%%
Len_all=N_FFT; %
recvd_signal_zeropad=[zeros(1,Len_all) recvd_signal(1,(n-1)*Nofdm+1:n*Nofdm) zeros(1,Len_all) zeros(1,Len_all)];
for d=1:length(recvd_signal_zeropad)-Len_all-1
   P_Park(d)=sum(recvd_signal_zeropad(d+(1:N_FFT/2)).*fliplr(recvd_signal_zeropad(d+N_FFT/2+(1:N_FFT/2)))); %fliplr
   E_Park(d)=sum(abs(recvd_signal_zeropad(d+(1:N_FFT))).^2);
   Q_Park(d)=sum(conj(recvd_signal_zeropad(d+(1:N_FFT/2))).*fliplr(recvd_signal_zeropad(d+N_FFT/2+(1:N_FFT/2))));
end
Park_metric=(abs(P_Park).^2)./(E_Park).^2;
buff = [buff Park_metric];
plot(Park_metric/max(Park_metric))

[~, maxipos1]=max(Park_metric-N_GI);
t_est_Park(1,n)=Nofdm - maxipos1;
FreqOffset(1,n) = angle(Q_Park(maxipos1-N_GI))/pi;
end

mse_CFO(k) = sqrt(sum(abs(FreqOffset - CFO).^2)/Nsym)
if flTimeShift == 1
   mse_STO(k) = sqrt(sum(abs(diff(t_est_Park - nSTOs) - diff(int_delay(1,:))).^2)/Nsym) %
else
   mse_STO(k) = sqrt(sum(abs(t_est_Park - nSTOs).^2)/Nsym)
end

end

 

function sig_d = sig_delay(sig,T,dt)
% Задержка сигнала по времени
% sig - отсчеты сигнала
%T - шаг дискретизации
%dt - задержка по времени
%sig_d - задержанный сигнал
if dt==0
   sig_d=sig;
   return
end

N=numel(sig);
sig_d(1:N)=0;
a=fix(dt/T);
b=dt/T-a;

x=-b;
y0=0;y1=0;y2=0;y3=0;

for m=1:N
   if m-a-2>0
       y0=sig(m-a-2);
   end
   if m-a-1>0
       y1=sig(m-a-1);
   end
   if m-a>0
       y2=sig(m-a);
   end
   if (m-a+1>0)&&(m-a+1<=N)
       y3=sig(m-a+1);
   end
   a3=(y3-y0)/6+(y1-y2)/2;
   a1=(y3-y1)/2-a3;
   a2=y3-y2-a1-a3;
   a0=y2;

   sig_d(m)=((a3*x+a2)*x+a1)*x+a0;   
end
end

 

но чет у меня какие-то сомнения в адекватности модели

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

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

Параметр fd у объекта rayleighchan задает не сдвиг по частоте, а доплеровское рассеяние, моделируя одновременный приход множества лучей с разной доплеровской частотой, объединяя их в один кластер. По факту генерится комплексный гауссовский случайный процесс с СПМ, заданной параметром DopplerSpectrum и умножается на входной сигнал. Это я на всякий случай написал, может вам это и нужно было. В принципе дальше по коду частотный сдвиг моделируется отдельно.

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

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

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

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

Вот в статье Парка10.1109_LCOMM.2003.812181.pdf у него тоже этот пик тока с другой стороны, видимо потому что суммирование в другом порядке делал

отсюда я и предположил что префикс нужен

А в жизни преамбулы разве без защитного передают?

Изменено пользователем Ivan55

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

В вашей модели у меня не получилось. Один хороший пик. В статье нарисовали два, но дальше на эту тему не наспространялись, а жаль. Сам я этот алгоритм не использовал. Все, что сложнее автокоррелятора не очень люблю, очень сложная схема получается. Может на DSP нормально можно реализовать.

ЗИ в преамбуле нужен, если она же будет использоваться для оценивания канала. В остальных случаях он скорее мешает.

 

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

Присоединяйтесь к обсуждению

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

Гость
Ответить в этой теме...

×   Вставлено с форматированием.   Вставить как обычный текст

  Разрешено использовать не более 75 эмодзи.

×   Ваша ссылка была автоматически встроена.   Отображать как обычную ссылку

×   Ваш предыдущий контент был восстановлен.   Очистить редактор

×   Вы не можете вставлять изображения напрямую. Загружайте или вставляйте изображения по ссылке.

×
×
  • Создать...