jakes_model_1.m

using Jakes model to simulate rayleigh channel fading

% Simulation Of Jakes Model 
clear all; 
f_max = 10; 
M = 8;                                 % # of low frequency oscillators - 1 
N = 4*M+2; 
Ts=1e-03; 
sq = 2/sqrt(N); 
sigma = 1/sqrt(2); 
theta = 0;                              % Fixed Phase 
count = 0; 
t0=0.24; %(0.24, 0.45); 
    for t = 0:Ts:3            % Varying time 
               
        count = count + 1 
        g(count) = 0; 
        for n = 1 : M+1, 
            if n <= M 
               c_q(count,n) = 2*sigma*sin(pi*n/M);  % Gain associated with quadrature component 
               c_i(count,n) = 2*sigma*cos(pi*n/M);  % Gain associated with inphase component 
               f_i(count,n) = f_max*cos(2*pi*n/N);  % Discrete doppler frequencies of inphase component 
               f_q(count,n) = f_max*cos(2*pi*n/N);  % Discrete doppler frequencies of quadrature component 
           else 
               c_i(count,n) = sqrt(2)*cos(pi/4); 
               c_q(count,n) = sqrt(2)*sin(pi/4); 
               f_i(count,n) = f_max; 
               f_q(count,n) = f_max; 
           end;  % end if 
                
               g_i(count,n) =  c_i(count,n)*cos(2*pi*f_i(count,n)*(t-t0) + theta); % Inphase component for one oscillator 
              
               g_q(count,n) =  c_q(count,n)*cos(2*pi*f_q(count,n)*(t-t0) + theta); % Quadrature component for one oscillator 
                
       end;   %end n 
       tp(count) =  sq*sum(g_i(count,1:M+1));% Total Inphase component  
       tp1(count) = sq*sum(g_q(count,1:M+1));% Total quadrature component  
       tpp = tp + j*tp1;
    end;  % end count no n again 
   envelope=sqrt(tp.^2+tp1.^2); 
   rmsenv=sqrt(sum(envelope.^2)/count); 
   
   [auto_real,lag_i] = xcorr(tpp,'coeff') ; 
%    [auto_i,lag_i] = xcorr(tp,'coeff') ;   % Auto-correlation associated with inphase component 
%    [auto_q,lag_q] = xcorr(tp1,'coeff');  % Auto-correlation associated with quadrature component 
   
    len=length(lag_i); 
   [corrx2,lag2] = xcorr(tp,tp1,'coeff'); 
   
% Cross Correlation between inphase and quadrature components 
   
      aa=-(len-1)/2:1:(len-1)/2;%total duration for lag 
      bb=(len-2001)./2;%mid ... points for drawing figures 
      cc=bb+1:1:bb+2001;%for getting the mid-values  
      dd=-1000:1:1000; 
  %_______________ 
  tdd=dd*Ts; 
  z=2.*pi.*f_max*tdd; 
  sigma0=1; 
  T_bessel=sigma0.^2.*besselj(0,z); 
      
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
   figure; 
   plot(tdd,auto_real(cc),'-',tdd,real(T_bessel),'*');%in-phase 
   xlabel('t(Second)'); 
   ylabel('Auto-correlation'); 
   legend('I-Simulated','I-Theoretic'); 
    
%       figure; 
%    plot(tdd,auto_q(cc),'-',tdd,real(T_bessel),'*');%quadrature 
%    xlabel('t(Second)'); 
%    ylabel('Auto-correlation'); 
%    legend('Q-Simulated','Q-Theoretic'); 
   
   %figure; 
   %plot(tdd,corrx2(cc)); 
   %xlabel('t(Second)'); 
   %ylabel('cross-correlation'); 
    
    
figure 
co1=1:1000; 
 
semilogy(co1*Ts,envelope(1:1000)); 
 xlabel('t(Second)'); 
 ylabel('Rayleigh Coef.'); 
%%__________________________________________ 
length_r=length(envelope); 
 
%____________________________ 
pdf_env=zeros(1,501); 
 
count=0; 
temp=round(100.*envelope); 
for k=1:length_r 
    if temp(k)<=500 
        count=count+1; 
        pdf_env(1,temp(k)+1)=pdf_env(1,temp(k)+1)+1; 
    end 
end     
count 
pdf_env=pdf_env./count./0.01; 
 
sgma2=0.5; 
x=[0:0.01:5]; 
pdf_theory=(x./sgma2).*exp(-1.*x.^2./(2.*sgma2)); 
     
figure 
plot(x,pdf_env,'-',x,pdf_theory,'*'); 
legend('Simulated','Theoretic'); 
xlabel('r'); 
ylabel('PDF of r');