%**************************************************************************
%  Functionfile:  init_fading.m
%  Description:
%  使用Jake模型仿真瑞利衰落信道
%***************************************************************************

% ***************************************************************** 
% change history:                                      
% date                  author                 Descriptionr           
% 2008/10/20           dengjuan                  Create        
% ***************************************************************** 
%*******************Pseudo Code*******************************************
%input
%FadingSamplingTime,表示表示衰落采樣時間。
%NumberOfTxAntennas，表示發射天線數。
%NumberOfRxAntennas，表示接收天線數。
%FadingOversamplingFactor，表示采樣因子。
%NumberOfPaths，表示多徑數。
%FadingNumberOfIterations，表示求衰落系數的迭代次數。
%fd，表示多普勒頻率。

%output
%Fadingmatrix,Jake模型的衰落矩陣
%****************************************************************
function Fadingmatrix=init_fading(FadingSamplingTime, ...
					FadingOversamplingFactor, ...
					NumberOfTxAntennas, ...
					NumberOfRxAntennas, ...
					NumberOfPaths,FadingNumberOfIterations,fd)
t=0:FadingSamplingTime:(FadingNumberOfIterations-1)*FadingSamplingTime;   
s=35;%s=N0

for l=1:NumberOfPaths
    gama=2*pi*(l-1)/(s+1);
for i=1:NumberOfTxAntennas*NumberOfRxAntennas
        hc(i,:)=zeros(1,length(t));
        hs(i,:)=zeros(1,length(t));
        for k=1:s
            fai(k)=2*pi*rand;
            beita(k)=pi*k/(s+1);
            sita(k)=beita(k)+gama;
            womiga(k)=2*pi*fd*cos(pi*k/(2*s+1));
            hc(i,:)=hc(i,:)+sqrt(2)/sqrt(4*s+2)*(2*cos( womiga(k)*t+fai(k))*cos(beita(k)/2)*cos(sita(k))-2*sin( womiga(k)*t+fai(k))*cos(beita(k)/2)*sin(sita(k)));
            hs(i,:)=hs(i,:)+sqrt(2)/sqrt(4*s+2)*(2*cos( womiga(k)*t+fai(k))*sin(beita(k)/2)*cos(sita(k))-2*sin( womiga(k)*t+fai(k))*sin(beita(k)/2)*sin(sita(k)));
        end
            h(i,:)=hc(i,:)+j*hs(i,:);
    end
    hh(l,:,:)=h;
end
Fadingmatrix=zeros(NumberOfTxAntennas*NumberOfRxAntennas*NumberOfPaths,FadingNumberOfIterations);
for ii=1:NumberOfPaths
    for jj=1:NumberOfTxAntennas*NumberOfRxAntennas
        for kk=1:FadingNumberOfIterations
            Fadingmatrix((ii-1)*NumberOfTxAntennas*NumberOfRxAntennas+jj,kk)=hh(ii,jj,kk);
        end
    end
end
% Normalisation to 1
Fadingmatrix= Fadingmatrix./sqrt(mean(abs(Fadingmatrix).^2,2)*...
				      ones(1,size(Fadingmatrix,2)));