802.11n信道模型matlab代码.doc

function R, Q, sigma_deg = correlation(M, spacing, d_norm, . cluster_number, amplitude_cluster, . PAS_type, phi_deg, AS_deg, . delta_phi_deg, type)% R, Q, sigma_deg = correlation(M, spacing, d_norm,% cluster_number, amplitude_cluster, PAS_type,% phi_deg, AS_deg, delta_phi_deg, type)% Derives the correlation matrix R from the description of the% environment, namely the antenna element spacing and the way the% waves impinge (number of clusters, PAS type, AS and mean angle of% incidence). The generated correlation coefficients are either% (complex) field correlation coefficients (type = 0) or (real% positive) power correlation coefficients (type = 1).% Inputs% * Variable M, number of antenna elements of the ULA% * Variable spacing, spacing of the antenna elements of the ULA% * Vector d_norm, containing the relative spacings of the M% elements of the ULA with respect to the first one% * Variable cluster_number, number of impinging clusters% * Vector amplitude_cluster, containing the amplitude of% the cluster_number impinging clusters% * Variable PAS_type, defines the nature of the PAS, isotropic% uniform (1), isotropic Gaussian (2), isotropic Laplacian (3)% or directive Laplacian (6) according to 3GPP-3GPP2 SCM AHG% radiation pattern% * Vector phi_deg, containing the AoAs of the number_cluster% clusters% * Vector AS_deg, containing the ASs of the number_cluster% clusters% * Vector delta_phi_deg, containing the constraints limits% of the truncated PAS, if applicable (Gaussian and Laplacian% case)% * Variable type, defines whether correlation properties% should be computed in field (0) or in power (1)% Outputs% * 2-D Hermitian matrix R of size M x M, whose elements are% the correlation coefficients of the corresponding elements% of the ULA impinged by the described PAS %Vector Q of cluster_number elements, giving the normalisation% coefficients to be applied to the clusters in order to have% the PAS fulfilling the definition of a pdf% * Vector sigma_deg of cluster_number of elements, giving% the standard deviations of the clusters, derived from% their description. There is not necessarily equality between% the given AS and the standard deviation of the modelling PAS.% Revision history% April 2003 - Addition of case 6 (directive Laplacian)% February 2002 - Creation% STANDARD DISCLAIMER% CSys is furnishing this item as is. CSys does not provideany% warranty of the item whatsoever, whether express, implied, or% statutory, including, but not limited to, any warranty of% merchantability or fitness for a particular purpose or any% warranty that the contents of the item will be error-free.% In no respect shall CSys incur any liability for any damages,% including, but limited to, direct, indirect, special, or% consequential damages arising out of, resulting from, or any way% connected to the use of the item, whether or not based upon% warranty, contract, tort, or otherwise; whether or not injury was% sustained by persons or property or otherwise; and whether or not% loss was sustained from, or arose out of, the results of, %the item, or any services that may be provided by CSys.% (c) Laurent Schumacher, AAU-TKN/IES/KOM/CPK/CSys - February 2002% Normalisationswitch(PAS_type) case 1 % Uniform distribution Q = normalisation_uniform(cluster_number, . amplitude_cluster, . AS_deg); sigma_deg = -1; % meaningless in uniform PAS case 2 % Gaussian distribution Q, sigma_deg = . normalisation_gaussian(cluster_number, . amplitude_cluster, . AS_deg, . delta_phi_deg); case 3, 6 % Laplacian distribution Q, sigma_deg = . normalisation_laplacian(cluster_number, . amplitude_cluster, . AS_deg, . delta_phi_deg); otherwise disp(PAS type non supported. Exiting.); return %breakend if (M1) % Complex correlation computation switch(PAS_type) case 1 % Isotropic Uniform distribution Rxx = bessel(0,2.*pi.*d_norm); Rxy = zeros(size(Rxx); for k=1:cluster_number Rxx = Rxx + Q(k).*Rxx_uniform(d_norm, . phi_deg(k), . AS_deg(k); Rxy = Rxy + Q(k).*Rxy_uniform(d_norm, . phi_deg(k), . AS_deg(k); end; case 2 % Isotropic Gaussian distribution Rxx = bessel(0,2.*pi.*d_norm); Rxy = zeros(size(Rxx); for k=1:cluster_number Rxx = Rxx + Q(k).*Rxx_gaussian(d_norm, . phi_deg(k), . sigma_deg(k), . delta_phi_deg(k); Rxy = Rxy + Q(k).*Rxy_gaussian(d_norm, . phi_deg(k), . sigma_deg(k), . delta_phi_deg(k); end; case 3 % Isotropic Laplacian distribution Rxx = bessel(0,2.*pi.*d_norm); Rxy = zeros(size(Rxx); for k=1:cluster_number Rxx = Rxx + Q(k).*Rxx_laplacian(d_norm, . phi_deg(k), . sigma_deg(k), . delta_phi_deg(k); Rxy = Rxy + Q(k).*Rxy_laplacian(d_norm, . phi_deg(k), . sigma_deg(k), . delta_phi_deg(k); end; case 6 % Directive Laplacian distribution % Hard-coded radiation pattern characteristics tmp = rho_lpl_SCM_AHG_approx(d_norm, . phi_deg, . Q, . sigma_deg, . delta_phi_deg, . 12, 20, 70); Rxx = real(tmp); Rxy = imag(tmp); end; % Correlation coefficient computation if (type = 0) tmp = Rxx + i.* Rxy; else tmp = Rxx.2 + Rxy.2; end; R = toeplitz(tmp);else R = 1;endfunction E=erfcomp(z); % ERFCOMP (komplexes) Fehlerintegral % % E=ERFCOMP(Z) % % Die (komplexe) Matrix Z ist das Eingabeargument % fr diese Prozedur. Zurckgeliefert wird die % komplexe Error-Funktion. % % Entnommen aus Abramowitz, % Abschnitt Error Functions and Fresnel Integrals, % Seite 299 zx,sx=size(z); % Dimension der Eingabematrix x=real(z); y=imag(z); x=x(:); y=y(:); i=sqrt(-1); % E=erf(x); KOa=+0.3275911; % Hier erfolgt die Berechnung der KOA=+0.254829592; % Error-Funktion durch eine rationale KOB=-0.284496736; % Approximation fr reelle Argumente. KOC=+1.421413741; KOD=-1.453152027; KOE=+1.061405429; w=1 ./(1 + KOa .*abs(x); P=KOE .*w + KOD; P=P .*w + KOC; P=P .*w + KOB; P=P .*w + KOA; P=P .*w; Q=exp(-(x .*x); E=(1 - Q .*P) .*sign(x); S=zeros(size(E); abserr=100; n=1; while abserr1E-12, fn=2 .*x - 2 .*x .*cosh(n .*y) .*cos(2 .*x .*y) + n .*sinh(n .*y) .*sin(2 .*x .*y); gn=2 .*x .*cosh(n .*y) .*sin(2 .*x .*y) + n .*sinh(n .*y) .*cos(2 .*x .*y); fgn=(fn + i .*gn) .*exp(-0.25 .*(n .2) ./( (n .2) + (4 .*(x .2); Sakt=2 .*fgn .*exp(-(x .2) ./pi; abserr=max(abs(Sakt); n=n+1; S=S + Sakt; end ki=find(x=0 & y=0); if isempty(ki), E(ki)=zeros(size(ki); end ki=find(x=0 & y=0); if isempty(ki), E(ki)=i .*y(ki) ./pi; end ki=find(x=0 & y=0); if isempty(ki), xk=x(ki); yk=y(ki); E(ki)=erf(xk) + (exp(-(xk .2) .*( (1 - cos(2 .*xk .*yk) + i .*sin(2 .*xk .*yk) ./(2 .*pi .*xk); end E=E+S; E=reshape(E,zx,sx);geometry2correlation% User interface to the definition of the correlation matrices% of the MIMO radio channel model. For both Node B and UE,% the user is prompted for% * the number of elements% * the spacing between elements (Uniform Linear Array is assumed)% * the number of impinging clusters of waves% * their Power Azimuth Spectrum (PAS) type (Uniform, Gaussian or Laplacian)% * the mean angle of incidence (in degrees)% * the Azimuth Spread (AS)% STANDARD DISCLAIMER% CSys is furnishing this item as is. CSys does not provide any% warranty of the item whatsoever, whether express, implied, or% statutory, including, but not limited to, any warranty of% merchantability or fitness for a particular purpose or any% warranty that the contents of the item will be error-free.% In no respect shall CSys incur any liability for any damages,% including, but limited to, direct, indirect, special, or% consequential damages arising out of, resulting from, or any way% connected to the use of the item, whether or not based upon% warranty, contract, tort, or otherwise; whether or not injury was% sustained by persons or property or otherwise; and whether or not% loss was sustained from, or arose out of, the results of, the% item, or any services that may be provided by CSys.% (c) Laurent Schumacher, AAU/TKN/I8/KOM/CPK/CSys - September 2001 clear;clf;display = 1; % Interactive dialog% % Downlink or uplink string = sprintf(n);disp(string);Downlink = -1;while (Downlink = 0) & (Downlink =1) Downlink = input(Downlink (Yes = 1, no = 0)? );end; % Field or envelope/power string = sprintf(n);disp(string);Field = -1;while (Field = 0) & (Field =1) Field = input(Type of correlation coefficient (Field = 0, Envelope = 1)? );end; % Configuration of Node B string = sprintf(n* Node B *n);M, spacing_Node_B, d_norm_Node_B, cluster_number_Node_B, . amplitude_cluster_Node_B, PAS_type_Node_B, phi_deg_Node_B, . AS_deg_Node_B, delta_phi_deg_Node_B = dialog(string); % Configuration of UE string = sprintf(n* UE *n);N, spacing_UE, d_norm_UE, cluster_number_UE, . amplitude_cluster_UE, PAS_type_UE, phi_deg_UE, AS_deg_UE, . delta_phi_deg_UE = dialog(string); % Computation of the correlation matrices% % Node B R_Node_B, Q_Node_B, sigma_deg_Node_B = correlation(M, spacing_Node_B, . d_norm_Node_B, . cluster_number_Node_B, . amplitude_cluster_Node_B, . PAS_type_Node_B, . phi_deg_Node_B, . AS_deg_Node_B, . delta_phi_deg_Node_B, Field); if (display & (M=1) figure(1); clf; switch(PAS_type_Node_B) case 1 plot_uniform(cluster_number_Node_B, Q_Node_B, phi_deg_Node_B, . AS_deg_Node_B, 1, b-); case 2 plot_gaussian(cluster_number_Node_B, Q_Node_B, phi_deg_Node_B, . sigma_deg_Node_B, delta_phi_deg_Node_B, 1, b-); case 3 plot_laplacian(cluster_number_Node_B, Q_Node_B, phi_deg_Node_B, . sigma_deg_Node_B, delta_phi_deg_Node_B, 1, b-); otherwise break; end;end; % UE R_UE, Q_UE, sigma_deg_UE = correlation(N, spacing_UE, d_norm_UE, . cluster_number_UE, . amplitude_cluster_UE, . PAS_type_UE, phi_deg_UE, . AS_deg_UE, delta_phi_deg_UE, Field);if (display & (N=1) figure(1); switch(PAS_type_UE) case 1 plot_uniform(cluster_number_UE, Q_UE, phi_deg_UE, AS_deg_UE, 2, b-); case 2 plot_gaussian(cluster_number_UE, Q_UE, phi_deg_UE, sigma_deg_UE, . delta_phi_deg_UE, 2, b-); case 3 plot_laplacian(cluster_number_UE, Q_UE, phi_deg_UE, sigma_deg_UE, . delta_phi_deg_UE, 2, b-); otherwise break end;end; % Information exploitation% if Downlink nTx = M; nRx = N; R = kron(R_Node_B, R_UE)else nTx = N; nRX = M; R = kron(R_UE, R_Node_B)end;% SelfGeometry2TGnCorrelation% Computation of the correlation matrices of the MIMO radio% channel model described in IEEE 802 document 11-03/940r0% for a limited set of hard-coded scenarios.% STANDARD DISCLAIMER% The Computer Science Institute of the University of Namur (hereafter% FUNDP-INFO) is furnishing this item as is. FUNDP-INFO does not% provide any warranty of the item whatsoever, whether express,% implied, or statutory, including, but not limited to, any warranty% of merchantability or fitness for a particular purpose or any% warranty that the contents of the item will be error-free.% In no respect shall FUNDP-INFO incur any liability for any damages,% including, but not limited to, direct, indirect, special, or% consequential damages arising out of, resulting from, or any way% connected to the use of the item, whether or not based upon% warranty, contract, tort, or otherwise; whether or not injury was% sustained by persons or property or otherwise; and whether or not% loss was sustained from, or arose out of, the results of, the% item, or any services that may be provided by FUNDP-INFO.% (c) Laurent Schumacher, FUNDP-INFO - November 2003 clear all;close all;delete TGnCorrelationMatrices.mat; TGnScenario = A,B,C,D,E,F;NumberOfTxAntennas = 1,2,4;SpacingTx = .5,1,4;NumberOfRxAntennas = 1,2,4;SpacingRx = .5,1,4; qq = 1;for (kk=1:size(TGnScenario,2) for (ll=1:size(NumberOfTxAntennas,2) for (mm=1:size(NumberOfRxAntennas,2) for (nn=1:size(SpacingTx,2) for (oo=1:size(SpacingRx,2) Scenario = TGnScenario(kk),_,num2str(NumberOfTxAntennas(ll),. x,num2str(NumberOfRxAntennas(mm),_; if (SpacingTx(nn) 1) Scenario = Scenario,i,num2str(inv(SpacingTx(nn),_; else Scenario=Scenario,num2str(SpacingTx(nn),_; end; if (SpacingRx(oo) 1) Scenario = Scenario,i,num2str(inv(SpacingRx(oo); else Scenario=Scenario,num2str(SpacingRx(oo); end; time = clock; disp(num2str(time(4),:,num2str(time(5),.,num2str(floor(time(6),. - Iteration ,num2str(qq),/,num2str(size(TGnScenario,2)*. size(NumberOfTxAntennas,2)*size(NumberOfRxAntennas,2)*. size(SpacingTx,2)*size(SpacingRx,2), - ,Scenario); % % Description of scenarios % switch TGnScenario(kk) case A % Flat fading with 0 ns rms delay spread (one tap at 0 ns delay model) PDP_dB = 0; % Average power dB 0; % Relative delay (ns) % Power roll-off coefficients Power_per_angle_dB = 0; % Tx AoD_Tx_deg = 45; AS_Tx_deg = 40; Type_Tx = 3; % Rx AoA_Rx_deg = 45; AS_Rx_deg = 40; Type_Rx = 3; case B % Typical residential environment, LOS conditions, 15 ns rms delay spread, and 10 dB Ricean K-factor at first delay PDP_dB = 0 -5.4287 -2.5162 -5.8905 -9.1603 -12.5105 -15.6126 -18.7147 -21.8168; % Average power dB 0 10e-9 20e-9 30e-9 40e-9 50e-9 60e-9 70e-9 80e-9 ; % Relative delay (ns) % Power roll-off coefficients Power_per_angle_dB = 0 -5.4287 -10.8574 -16.2860 -21.7147 -Inf -Inf -Inf -Inf; -Inf -Inf -3.2042 -6.3063 -9.4084 -12.5105 -15.6126 -18.7147 -21.8168; % Tx AoD_Tx_deg = 225.1084.*ones(1,5) -Inf.*ones(1,4); -Inf.*ones(1,2) 106.5545.*ones(1,7); AS_Tx_deg = 14.4490.*ones(1,5) -Inf.*ones(1,4); -Inf.*ones(1,2) 25.4311.*ones(1,7); Type_Tx = 3.*ones(1, size(AoD_Tx_deg, 2); % Rx AoA_Rx_deg = 4.3943.*ones(1,5) -Inf.*ones(1,4); -Inf.*ones(1,2