MIMO天线发射技术matlab程序

单入单出SISO系统：

1. clear all;
   
2. close all;
   
3. clc;
   
4. M = 1000;
   
5. Nt = [1];           % Nt = [2]; Nt = [4]; Nt = [8];
   
6. Nr = [1];
   
7. s={'b-hexagram'};
   
8. SNR = [0:2:50]; %dB
   
9. figure(1);
   
10. title('单入单出系统');
    
11. xlabel('SNR [dB]');
    
12. ylabel('CU');
    
13. grid on;
    
14. hold on;
    
15. for n = 1:length(Nt)
    
16. t = Nt(n);
    
17. r = Nr(n);
    
18. for m = 1:M
    
19. H = raylrnd(1,r,t);
    
20. for l = 1:length(SNR)
    
21. rou = 10^(SNR(l)/10);
    
22. % Find the capacity while CSIT unknown
    
23. CU(m,l) = log2(real(det(eye(r)+rou*H*H'/t)));
    
24. end
    
25. end
    
26. C_unknown(:,n) = mean(CU)'
    
27. for l = 1:length(SNR)
    
28. [cdf_u,co_u] = hist(CU(:,l),100);
    
29. cdf_u = cumsum(cdf_u);%数组累加
    
30. idx_ten_percent = find(abs(cdf_u-100)==min(abs(cdf_u-100)));
    
31. C_unknown_outage(l,n) = co_u(idx_ten_percent(1));
    
32. end
    
33. figure(1);
    
34. linetype = s(n);
    
35. plot(SNR,C_unknown(:,n),deal(linetype{:}));
    
36. legend('Mt=1');
    
37. end

复制代码

二发单收系统：

四发单收系统：

八发单收系统：

多发多收和单发单收系统比较：

1. clear all;
   
2. 
   
3. close all;
   
4. 
   
5. clc;
   
6. 
   
7. M = 1000;
   
8. 
   
9. Nt = [1,2,4,8];
   
10. 
    
11. Nr = [1,2,4,8];
    
12. 
    
13. s={'b-o','b-*','b-square','b-hexagram'};
    
14. 
    
15. SNR = [0:2:50]; %dB
    
16. 
    
17. figure(1);
    
18. 
    
19. title('多发多收和单发单收系统比较')
    
20. 
    
21. xlabel('SNR [dB]');
    
22. 
    
23. ylabel('CU');
    
24. 
    
25. grid on;
    
26. 
    
27. hold on;
    
28. 
    
29. for n = 1:length(Nt)
    
30. 
    
31. t = Nt(n);
    
32. 
    
33. r = Nr(n);
    
34. 
    
35. for m = 1:M
    
36. 
    
37. H = raylrnd(1,r,t);
    
38. 
    
39. for l = 1:length(SNR)
    
40. 
    
41. rou = 10^(SNR(l)/10);
    
42. 
    
43. % Find the capacity while CSIT unknown
    
44. 
    
45. CU(m,l) = log2(real(det(eye(r)+rou*H*H'/t)));
    
46. 
    
47. end
    
48. 
    
49. end
    
50. 
    
51. C_unknown(:,n) = mean(CU)';
    
52. 
    
53. for l = 1:length(SNR)
    
54. 
    
55. [cdf_u,co_u] = hist(CU(:,l),100);
    
56. 
    
57. cdf_u = cumsum(cdf_u);%数组各行累加
    
58. 
    
59. idx_ten_percent = find(abs(cdf_u-100)==min(abs(cdf_u-100)));
    
60. 
    
61. C_unknown_outage(l,n) = co_u(idx_ten_percent(1));
    
62. 
    
63. end
    
64. 
    
65. figure(1);
    
66. 
    
67. linetype = s(n);
    
68. 
    
69. plot(SNR,C_unknown(:,n),deal(linetype{:}));
    
70. 
    
71. legend('Mt=1','Mt=2','Mt=4','Mt=8');
    
72. 
    
73. end

复制代码

STBC编码原理：

x及其共轭的线性组合。一个编码码字共有P个时刻，并按行由N副天线同时发送，即在第一个时刻发送第一行，第二个时刻发送第二行，依此类推。在第t个时刻发送第t行，总共需P个时刻才可完成一个编码码字的发送。因此，矩阵的每一列符号实际是由同一副发送天线在不同时刻发送的。

Clear

N = 10^6; % number of bits or symbols

Eb_N0_dB = [0:25]; % multiple Eb/N0 values多个信噪比，信号传输的衰弱程度不一

nRx = 2;%接收天线的个数为2

for ii = 1:length(Eb_N0_dB)%信噪比的个数有多少，就做多少次循环

% Transmitter发射机

ip = rand(1,N)>0.5; % generating 0,1 with equal probability（rand产生一个1行N列的随机数）后面的">0.5"的意思是指当随机取得的数大于0.5时，那个数直接变为一，当小于0.5时为0

s = 2*ip-1; % BPSK 调制 0 -> -1; 1 -> 0（另外的一种调试方法，估计是仿真的时候用来对比的）

% Alamouti STBC

sCode = 1/sqrt(2)*kron(reshape(s,2,N/2),ones(1,2)) ;%kron（A，B）是求矩阵A和矩阵B的乘积

reshape（s，2，N/2）是将s矩阵改成2行，N/2列的矩阵；ones（1,2）生成一个1行2列的全1矩阵

% channel

h = 1/sqrt(2)*[randn(nRx,N) + 1i*randn(nRx,N)]; % Rayleigh channel (sqrt平方根)

n = 1/sqrt(2)*[randn(nRx,N) + 1i*randn(nRx,N)]; % white gaussian noise, 0dB variance

%

y = zeros(nRx,N);

yMod = zeros(nRx*2,N);

hMod = zeros(nRx*2,N);

for kk = 1:nRx %从1一直循环到nRx

hMod = kron(reshape(h(kk,:),2,N/2),ones(1,2)); % repeating the same channel for two symbols (重复相同信道的两个符号)

temp = hMod;

hMod(1,[2:2:end]) = conj(temp(2,[2:2:end]));

hMod(2,[2:2:end]) = -conj(temp(1,[2:2:end]));

% Channel and noise Noise addition

y(kk,:) = sum(hMod.*sCode,1) + 10^(-Eb_N0_dB(ii)/20)*n(kk,:);

% Receiver

yMod([2*kk-1:2*kk],:) = kron(reshape(y(kk,:),2,N/2),ones(1,2));

% forming the equalization matrix（形成均衡矩阵）

hEq([2*kk-1:2*kk],:) = hMod;

hEq(2*kk-1,[1:2:end]) = conj(hEq(2*kk-1,[1:2:end]));

hEq(2*kk, [2:2:end]) = conj(hEq(2*kk, [2:2:end]));

end

% equalization

hEqPower = sum(hEq.*conj(hEq),1);

yHat = sum(hEq.*yMod,1)./hEqPower; % [h1*y1 + h2y2*, h2*y1 -h1y2*, ... ]

yHat(2:2:end) = conj(yHat(2:2:end));

% receiver - hard decision decoding

ipHat = real(yHat)>0;

% counting the errors

nErr(ii) = size(find([ip- ipHat]),2);

end

simBer = nErr/N; % simulated ber

EbN0Lin = 10.^(Eb_N0_dB/10);

theoryBer_nRx1 = 0.5.*(1-1*(1+1./EbN0Lin).^(-0.5));

p = 1/2 - 1/2*(1+1./EbN0Lin).^(-1/2);

theoryBerMRC_nRx2 = p.^2.*(1+2*(1-p));

pAlamouti = 1/2 - 1/2*(1+2./EbN0Lin).^(-1/2);

theoryBerAlamouti_nTx2_nRx1 = pAlamouti.^2.*(1+2*(1-pAlamouti));

close all

figure

semilogy(Eb_N0_dB,theoryBer_nRx1,'bp-','LineWidth',2); %semilogy：绘制半对数坐标图形

hold on;

semilogy(Eb_N0_dB,theoryBerMRC_nRx2,'kd-','LineWidth',2);

semilogy(Eb_N0_dB,theoryBerAlamouti_nTx2_nRx1,'c+-','LineWidth',2);

semilogy(Eb_N0_dB,simBer,'mo-','LineWidth',2);

axis([0 25 10^-5 0.5]);

grid on;

legend('theory (nTx=1,nRx=1)', 'theory (nTx=1,nRx=2, MRC)', 'theory (nTx=2, nRx=1, Alamouti)', 'sim (nTx=2, nRx=2, Alamouti)');

xlabel('Eb/No, dB');

ylabel('Bit Error Rate');

title('BER for BPSK modulation with 2Tx, 2Rx Alamouti STBC (Rayleigh channel)');