cornaka_noise.m

MIMO系统

%相关衰落信道和噪声--Sfunction
function [sys,x0,str,ts]=cornaka_noise(t,x,u,flag,L,N,BitRate,m,Rou_gamma_Tx,Rou_gamma_Rx,omegap,sgma)
j=sqrt(-1);
%Rou_gamma_Tx: correlation coefficient matrix of transmitter
%Rou_gamma_Rx: correlation coefficient matrix of receiver
%m: fading parameter
%omegap:  mean fading power omegap=E(x^2)
%N: number of transmit branches
%L: number of receive branches
%sgma： noise 标准差 
%h=[ h(1,1) h(1,2) ... h(1,N)
%    h(2,1) h(2,2) ... h(2,N)
%    ... ... ... ...
%    h(L,1) ... ... h(L,N)]
%g--Gaussian gam--Gamma x--Nakagami(|h|)
tic;
Rou_gamma=kron(Rou_gamma_Tx,Rou_gamma_Rx);%MIMO correlation coefficient matrix ((N*L)*(N*L))
%correlation coefficient matrix of x_vector=[x(1,1) x(2,1) ...x(L,1) x(1,2)...x(L,2)......x(1,N) ... x(L,N)]'
C_gamma=Rou_gamma*(omegap^2/m);
%covariance matrix of Gamma
C_g=sqrt(C_gamma)/(2*sqrt(m));
%covariance matrix of Gaussian

L1=chol(C_g)';%cholesky factorization.
%lower triangular matrix L

if (2*m-floor(2*m))==0 %%2*m is integer
    g=L1*randn(L*N,2*m*BitRate);%2*m*BitRate组独立的高斯随机向量(L*N*(2*m*BitRate))（每一组向量的协方差矩阵为C_g）
  for i=1:L*N
       gam(i,:)=sum(reshape(g(i,:).^2,2*m,BitRate));%Gamma RVs 
  end

elseif (floor(2*m)~=1)
  p=floor(2*m);
  alpha=(2*p*m+sqrt(2*p*m*(p+1-2*m)))/(p*(p+1));
  beta=2*m-p*alpha;
  g=L1*randn(L*N,p*BitRate);%p*BitRate组独立的高斯随机向量(L*N*(p*BitRate))（每一组向量的协方差矩阵为C_g）
  for i=1:L
       gam1(i,:)=alpha*sum(reshape(g(i,:).^2,p,BitRate));%the first of Gamma RVs 
  end
 gam=gam1+beta*(L1*randn(L*N,BitRate)).^2;%Gamma RVs
else
    p=1;
    alpha=(2*p*m+sqrt(2*p*m*(p+1-2*m)))/(p*(p+1));
    beta=2*m-p*alpha;
    g=L1*randn(L*N,BitRate);%BitRate组独立的高斯随机向量((L*N)*BitRate)（每一组向量的协方差矩阵为C_g）
    gam=alpha*g.^2+beta*(L1*randn(L*N,BitRate)).^2;%Gamma RVs((L*N)*BitRate) 
end

delta1=sum((reshape(gam',L*BitRate,N))');%(1*(L*BitRate))
delta=2*sum((reshape(delta1,BitRate,L))');%(1*BitRate)
epsilon=normrnd(0,sqrt(delta)*sgma)+j*normrnd(0,sqrt(delta)*sgma);%(1*BitRate))(noise)

Nakanew=[delta';(real(epsilon))';(imag(epsilon))'];%((3*BitRate)*1)符合S函数的模式

switch flag,
    case 0
        [sys,x0,str,ts]=mdlInitializeSizes(Nakanew);%初始化模块
    case 3
        sys=mdlOutputs(t,x,u,Nakanew);%计算输出
    case {1,2,4,9}
        sys=[];%无用的标志
    otherwise
        error(['Unhandled flag =',num2str(flag)]);%错误处理
end
%模块结束
%%%%%%%初始化状态、采样周期、Size结构%%%%%%%%%
function [sys,x0,str,ts]=mdlInitializeSizes(Nakanew)
%调用simsizes创建sizes结构
sizes=simsizes;
%填写sizes结构
sizes.NumContStates=0;
sizes.NumDiscStates=0;
sizes.NumOutputs=length(Nakanew);
sizes.NumInputs=1;
sizes.DirFeedthrough=1;
sizes.NumSampleTimes=1;
sys=simsizes(sizes);
x0=[];%没有连续状态
str=[];%没有状态次序
ts=[0 0];%继承性的采样周期

%%%%%模块输出%%%%%%%%%%
function sys=mdlOutputs(t,x,u,Nakanew)
sys=Nakanew.*u;