parameter_jakes.asv

《移动衰落信道》Mobiel_Fading_Channels一书后面的相关仿真程序。

%--------------------------------------------------------------------
% parameter_Jakes.m -------------------------------------------------
%
% Program for the computation of the discrete Doppler frequencies,
% Doppler coefficients and Doppler phases by using the Jakes power
% spectral density.
%
% Used m-files: LPNM_opt_Jakes.m, fun_Jakes.m,
% grad_Jakes.m, acf_mue.m
%--------------------------------------------------------------------
% [f_i_n,c_i_n,theta_i_n]=parameter_Jakes(METHOD,N_i,sigma_0_2,f_max,PHASE,PLOT)
%--------------------------------------------------------------------
% Explanation of the input parameters:
%
% METHOD:
% |----------------------------------------------|------------------|
% | Methods for the computation of the discrete  |       Input      |
% | Doppler frequencies and Doppler coefficients |                  |
% |----------------------------------------------|------------------|
% |----------------------------------------------|------------------|
% | Method of equal distances (MED)              |      'ed_j'      |
% |----------------------------------------------|------------------|
% | Mean square error method (MSEM) | 'ms_j' |
% |----------------------------------------------|------------------|
% | Method of equal areas (MEA) | 'ea_j' |
% |----------------------------------------------|------------------|
% | Monte Carlo method (MCM) | 'mc_j' |
% |----------------------------------------------|------------------|
% | Lp-norm method (LPNM) | 'lp_j' |
% |----------------------------------------------|------------------|
% | Method of exact Doppler spread (MEDS) | 'es_j' |
% |----------------------------------------------|------------------|
% | Jakes method (JM) | 'jm_j' |
% |----------------------------------------------|------------------|
%
% N_i: number of harmonic functions
% sigma_0_2: average power of the real deterministic Gaussian process mu_i(t)
% f_max: maximum Doppler frequency
%
% PHASE:
% |----------------------------------------------|------------------|
% | Methods for the computation of the Doppler | Input |
% | phases | |
% |----------------------------------------------|------------------|
% |----------------------------------------------|------------------|
% | Random Doppler phases | 'rand' |
% |----------------------------------------------|------------------|
% | Permuted Doppler phases | 'perm' |
% |----------------------------------------------|------------------|
%
% PLOT: plot of the ACF and the PSD of mu_i(t), if PLOT==1
function [f_i_n,c_i_n,theta_i_n]=parameter_Jakes(METHOD,N_i,sigma_0_2,f_max,PHASE,PLOT)
if nargin<6, % 函数中参数输入个数
   error('Not enough input parameters')
end

sigma_0=sqrt(sigma_0_2);

% Method of equal distances (MED,等距法)
if     METHOD=='ed_j',
       n=(1:N_i)';
       f_i_n=f_max/(2*N_i)*(2*n-1);
       c_i_n=2*sigma_0/sqrt(pi)*(asin(n/N_i)-asin((n-1)/N_i)).^0.5;
       K=1;

% Mean square error method (MSEM,均方误差法)
elseif METHOD=='ms_j',
       n=(1:N_i)';
       f_i_n=f_max/(2*N_i)*(2*n-1);
       Tp=1/(2*f_max/N_i);
       t=linspace(0,Tp,5E3);
       Jo=besselj(0,2*pi*f_max*t);
       c_i_n=zeros(size(f_i_n));
       for k=1:length(f_i_n),
           c_i_n(k)=2*sigma_0*...
                    sqrt(1/Tp*( trapz( t,Jo.*...
                    cos(2*pi*f_i_n(k)*t )) ));
       end
       K=1;
       
% Method of equal areas (MEA,等面积法)
elseif METHOD=='ea_j'
       n=(1:N_i)';
       f_i_n=f_max*sin(pi*n/(2*N_i));
       c_i_n=sigma_0*sqrt(2/N_i)*ones(size(n));
       K=1;
       
% Monte Carlo method (MCM,蒙特卡洛法)
elseif METHOD=='mc_j'
       n=rand(N_i,1);
       f_i_n=f_max*sin(pi*n/2);
       c_i_n=sigma_0*sqrt(2/N_i)*ones(size(n));
       K=1;
       
% Lp-norm method (LPNM)
elseif METHOD=='lp_j',
       if exist('fminu')~=2
          disp([' =====> This method requires ',...
                'the Optimization Toolbox !!'])
          return
       else
          N=1E3;
          p=2; % Norm
          s_o=1;
          [f_i_n,c_i_n]=LPNM_opt_Jakes(N,f_max,sigma_0,p,N_i,s_o);
          K=1;
       end
       
% Method of exact Doppler spread (MEDS,精确多普勒扩展法)
elseif METHOD=='es_j',
       n=(1:N_i)';
       f_i_n=f_max*sin(pi/(2*N_i)*(n-1/2));
       c_i_n=sigma_0*sqrt(2/(N_i))*ones(size(f_i_n));
       K=1;
       
% Jakes method (JM,Jakes 法)
elseif METHOD=='jm_j',
       n=1:N_i-1;
       f_i_n=f_max*[[cos(pi*n/(2*(N_i-1/2))),1]',...
                    [cos(pi*n/(2*(N_i-1/2))),1]'];
       c_i_n=2*sigma_0/sqrt(N_i-1/2)*[[sin(pi*n/(N_i-1)),1/2]',...
                                     [cos(pi*n/(N_i-1)),1/2]'];
       K=1;
       theta_i_n=zeros(size(f_i_n));
       PHASE='none';
else
       error('Method is unknown')
end

% Computation of the Doppler phases():
if     PHASE=='rand',
       theta_i_n=rand(N_i,1)*2*pi;
       
elseif PHASE=='perm',
       n=(1:N_i)';
       Z=rand(size(n));
       [dummy,I]=sort(Z);
       theta_i_n=2*pi*n(I)/(N_i+1);
end;

if PLOT==1,
   if METHOD=='jm_j'
      subplot(2,3,1)
      stem([-f_i_n(N_i:-1:1,1);f_i_n(:,1)],...
            1/4*[c_i_n(N_i:-1:1,1);c_i_n(:,1)].^2)
      title('i=1')
      xlabel('f (Hz)')
      ylabel('PSD')
      subplot(2,3,2)
      stem([-f_i_n(N_i:-1:1,2);f_i_n(:,2)],...
      1/4*[c_i_n(N_i:-1:1,2);c_i_n(:,2)].^2)
      title('i=2')
      xlabel('f (Hz)')
      ylabel('PSD')
      tau_max=N_i/(K*f_max);
      tau=linspace(0,tau_max,500);
      r_mm=sigma_0^2*besselj(0,2*pi*f_max*tau);
      
      r_mm_tilde1=acf_mue(f_i_n(:,1),c_i_n(:,1),tau);
      subplot(2,3,4)
      plot(tau,r_mm,'r-',tau,r_mm_tilde1,'g--')
      title('i=1')
      xlabel('tau (s)')
      ylabel('ACF')
      r_mm_tilde2=acf_mue(f_i_n(:,2),c_i_n(:,2),tau);
      subplot(2,3,5)
      plot(tau,r_mm,'r-',tau,r_mm_tilde2,'g--')
      title('i=2')
      xlabel('tau (s)')
      ylabel('ACF')
      subplot(2,3,3)
      stem([-f_i_n(N_i:-1:1,1);f_i_n(:,1)],...
            1/4*[c_i_n(N_i:-1:1,1);c_i_n(:,1)].^2+...
            1/4*[c_i_n(N_i:-1:1,2);c_i_n(:,2)].^2)
      title('i=1,2')
      xlabel('f (Hz)')
      ylabel('PSD')
      subplot(2,3,6)
      plot(tau,2*r_mm,'r-',tau,r_mm_tilde1+r_mm_tilde2,'g--')
      title('i=1,2')
      xlabel('tau (s)')
      ylabel('ACF')
   else
      subplot(1,2,1)
      stem([-f_i_n(N_i:-1:1);f_i_n],...
            1/4*[c_i_n(N_i:-1:1);c_i_n].^2)
      xlabel('f/Hz')
      ylabel('LDS')
      tau_max=N_i/(K*f_max);
      tau=linspace(0,tau_max,500);
      r_mm=sigma_0^2*besselj(0,2*pi*f_max*tau);
      r_mm_tilde=acf_mue(f_i_n,c_i_n,tau);
      subplot(1,2,2)
      plot(tau,r_mm,'r-',tau,r_mm_tilde,'g--')
      xlabel('tau (s)')
      ylabel('ACF')
    end
end