OFDM-MATLAB仿真程序.doc

word完整版)OFDM MATLAB仿真程序
OFDM。m： OFDM Simulator (outer function)
clear all;

A = ［1 1/exp(1） 1/exp(2)］; ％ power delay profile
N = 64； % number of symbols in a single OFDM symbol
GI = 16； ％ guard interval
Mt = 1; % number of Tx antennas
Mr = 1； ％ number of Rx antennas
sig2 = 1e—3; % noise variance
M = 8； % max constellation bit number
Mgap = 10.^(1:（1.7/10）：2。7)； % gap
Btot = 100*Mt； % total ＃ bits per OFDM symbol
TransmitIter = 50; ％ # iterations of symbol transmissions for each channel instance
ChannelIter = 100； ％ # iterations of independent identically distributed channel instances
GapIter = length(Mgap）；

load ENC2.mat
load ENC4.mat
load ENC16.mat
load ENC64.mat
load ENC256.mat

TotEbNo = ［]；
Errors =［］；
EbNo = [];
for lGap = 1:GapIter
lGap
gap = Mgap（lGap)；
totalErrors = 0;
for lChan = 1：ChannelIter

％ create channel
[H h_f]=create_channel（Mt, Mr, A， N+GI）;

% decompose each subchannel in the frequency domain
[U S V］ = svd_decompose_channel（Mt， Mr, h_f， N）;

% bitloading
［bits_alloc，energy_alloc］ = BitLoad（S,Btot，Mt＊N,gap，sig2,M）;
%energy_alloc=energy_alloc/(mean(energy_alloc））；
％energy_alloc=ones（1,128);

for lTrans = 1：TransmitIter

％ bits to transmit
x = （randn（1，Btot）>0);

% modulate
x_mod = modulate(x,bits_alloc，energy_alloc， s2,s4，s16，s64，s256）；

％ precode modulated signal
x_pre = precode（Mt， x_mod, V， N);

% ifft, with cyclic prefix for each antenna
ofdm_symbol =[];
for i=1:Mt
ofdmsymbol = ［ofdm_symbol； ifft_cp_tx_blk（x_pre(i：Mt：Mt＊(N-1)+i）,N，GI)］;
end
ofdm_symbol2 = reshape(ofdm_symbol,Mt*(N+GI),1)；

％ channel
y = transpose（channel（sig2, Mt， Mr， ofdm_symbol2， H, N+GI)）；

％ fft
rec_symbol =[］;
for i=1:Mt
rec_symbol = ［rec_symbol； fft_cp_rx_blk（y（i:Mt:Mt*（N+GI—1)+i)，N,GI）］；
end

rec_symbol2 = reshape(rec_symbol，1，Mt*N）；

％ shape received signal
shaped_vals = shape（rec_symbol2, Mr, U, N);

% demodulate
y_demod = demodulate(shaped_vals, bits_alloc， energy_alloc， S, s2,s4,s16,s64,s256， c2,c4，c16,c64，c256）;

% comparison
totalErrors = totalErrors + sum(xor(y_demod，x)）；
end
EbNo = ［EbNo sum(energy_alloc）/Btot/sig2]；
end

Errors = [Errors totalErrors/Btot/ChannelIter/TransmitIter］
TotEbNo = [TotEbNo mean（EbNo)]
EbNo = ［];
end

semilogx（TotEbNo， Errors)；
xlabel(’Eb/No')；
ylabel('BER’）；
title（’SISO link, adaptive rate and power')

save SISO_adaptive2。mat Errors EbNo

create_channel。m: Generates a Rayleigh fading frequency—selective channel， parametrized by the antenna configuration， the OFDM configuration, and the power-delay profile.

function [H, H_f］=create_channel（Mt， Mr， A， N）;
% function ［H, H_f]=create_channel（Mt, Mr， A, N)；
％
% A - vector containing the power-delay profile (real values）
% Mt - number of Tx antennas
% Mr — number of Rx antennas
％ N - number of vector symbols to be sent in a single OFDM symbol Tx
％ ie: N MIMO transmissions in one OFDM symbol

% This is for Rayleigh frequency—selective fading, which assumes complex
％ Gaussian matrix elements with in-phase and quadrature components independent。
％ Assume iid matrix channel elements, and further， independent channel taps

% define the channel taps
H_int = 1/sqrt（2)*（randn（Mr*length(A),Mt) + j*randn(Mr*length（A)，Mt)）；

H_int2=［]；
for i = 1:length（A)
H_int2 = [H_int2;sqrt（A（i））*H_int(（i-1）＊Mr+1:i*Mr，:）］;
end

％h_f = fft(H_int2',64）；

%%H = H_int2';

H_int2 = ［H_int2;zeros（(N-length（A））＊Mr,Mt）］；

H_f = zeros（Mr,Mt＊（N—16)）；
for i = 1:Mt
for j = 1:Mr
h_f = fft(H_int2(j:Mr：（N—16—1)＊Mr+j,i）);
for k = 1：（N—16)
H_f（j，i+（k—1)＊Mt） = h_f（k);
end
end
end

H=[H_int2］；
for i = 1：N—1
H=［H,［zeros(Mr＊i,Mt）;H_int2(1:(N-i)*Mr，：）]]；
end
svd_decompose_channel。m: Since full channel knowledge is assumed, transmission is across parallel singular value modes。 This function decomposes the channel into these modes.

function ［U, S， V] = svd_decompose_channel（Mt, Mr, h_f, N);
％ [U S V］ = svd_decompose_channel（Mt, Mr， h_f, N);
%
％ Function decomposes the channel at each subcarrier into its SVD components
%
％ Mt — ＃ Tx antennas
% Mr - ＃ Rx antennas
％ h_f — MIMO impulse response - Mr rows, Mt＊L columns， where L is the number of
％ channel taps
% N — # subcarriers

U = [］;
S = [］;
V = ［］；
for i = 1:N
[Utmp Stmp Vtmp] = svd(h_f(：，(i—1）＊Mt+1:i*Mt)）；
U=［U Utmp］；
V=[V Vtmp］;
S=[S Stmp］;
end

S = sum（S，1)；
BitLoad.m： Apply the bit-loading algorithm to achieve the desired bit and energy allocation for the current channel instance.

function [bits_alloc，energy_alloc] = BitLoad（subchan_gains,total_bits，num_subc,gap,noise，M)

％ Bit Loading Algorithm
％ —-—----————-----————-
％
% Inputs :
% subchan_gains : SubCarrier Gains
％ total_bits ： Total Number of bits
% num_subc : Number of Subcarriers
％ gap ： Gap of the system
% noise ： Noise Power
％ M : Max Constellation Size
% Outputs:
％ bits_alloc ： Bits allocation for each subchannel
％ power_alloc : Total Power allocation
％ -----—-—-——————-----—--—-—-——-————---——----—-—---———-————--—--—

％ Compute SNR's for each channel
SNR = ComputeSNR(subchan_gains，noise,gap）;

% This function just initializes the system with a particular bit
% allocation and energy allocation using Chow's Algorithm。 This is
％ further efficientize using Campello’s Algorithm

［bits_alloc, energy_alloc］ = chow_algo(SNR,num_subc,M）;

% Form the Energy Increment Table based on the present channel
% gains for all the subchannels in order to be used by Campello
% Algorithm

energytable = EnergyTableInit（SNR，M)；

％ Efficientize the algorithm using the Campello's algorithm
［bits_alloc,energy_alloc] = campello_algo（bits_alloc，energy_alloc,energytable，total_bits，num_subc，M)；

ComputeSNR.m: Given the subcarrier gains， this simple function generates the SNR values of each channel （each singular value on each tone is a separate channel）。

function SNR = ComputeSNR（subcar_gains,noise，gap）

SNR = abs（(subcar_gains.^2）。/（noise＊gap)）;

chow_algo.m: Apply Chow's algorithm to generate a particular bit and energy allocation.

% Chow's Algorithm
％ -——--—-—-——-—-—-
％ This is based on the paper by Chow et al titled
％
% A Practical Discrete Multitone Transceiver Loading Algorithm
％ for Data Transmission over Spectrally Shaped Channels。IEEE Trans
% on Communications. Vol。 43, No 2/3/4， pp。 773—775, Feb/Mar/Apr 1995

function ［bits_alloc， energy_alloc] = chow_algo(SNR,num_subc，M)

for i = 1：num_subc
％ Assuming each of the subchannels has a flat fading， we get initial estimate
％ of the bits for each subchannel
tempbits = log2(1 + abs（SNR（i))）； % bits per two dimension.
roundtempbits = round(tempbits）； ％ round the bits

if (roundtempbits 〉 8) ％ Limit them between 2 and 15
roundtempbits = 8；
end

if (mod(roundtempbits，2)== 1 ＆ roundtempbits ～= 1)
roundtempbits = roundtempbits —1；
end

if roundtempbits 〉 0 ％ Calculate the Energy required for the subchannel
energy_alloc(i） = （2^roundtempbits—1）/SNR(i） ;
else
energy_alloc（i) = 0；
end

bits_alloc(i） = roundtempbits； ％ Update the BitSubChan
end

％ end of function

EnergyTableInit。m： Given the SNR values， form a table of energy increments for each channel。
function energytable = EnergyTableInit（SNR,M);

％ Inputs:
% subcar_gains ： Subcarrier Gains
% M ： max Constellation Size
％ Gap : Gap of the system
％ Noise ： Noise Power
% Outputs：
% energytable ： Energytable
％
％ Based on the Subcarrier Gains， we calculate the energy
％ increment required by each subcarrier for transmitting
％ 1,2 ,4 ,6，8 bits。
％ Energy = 2^（i-1)/subcar_gains；
％ —----——-—-——-—-—-——-——-—-———-------—-——--—---—-—-——-—-

%subcar_gains = （subcar_gains。^2）/（Gap*Noise)；

energytable = abs(（1./SNR）’＊（2.^(［1：M+1］-1）)）；

％ Increase the energy value for constellation size of more than M to
％ a very high value so that it is not assigned.

energytable(：,M+1） = Inf*ones（size(energytable（:，M+1)）);

for i = 3：2:M
energytable(：，i） = (energytable(：，i） +energytable（:,i+1））/2；
energytable(：,i+1) = energytable（：,i);
end

％energytable = ［ones(1,size（energytable，1))’ energytable］；
campello_algo。m: Apply Campello’s algorithm to converge to the optimal bit and energy allocation for the given channel conditions.
％ campello_algo.m
% ---——---—-———-
% This function is used by Campello's algorithm to allocate bits and energy for
% each subchannel optimally。

function ［bits_alloc， energy_alloc] = campello_algo(bits_alloc,energy_alloc,energytable，total_bitsnum_subc,M）

bt = sum(bits_alloc);

％ We can't transmit more than M*（Number of subchannel） bits
if total_bits > M*num_subc
total_bits = M＊num_subc;
end

while (bt ~= total_bits）
if (bt > total_bits）

max_val = 0；
max_ind = ceil（rand（1)*num_subc)；

for i = 1:num_subc

if bits_alloc（i) ~= 0
temp = energytable(i，bits_alloc（i）） ；
else
temp = 0;
end

if （temp > max_val）
max_val = temp;
max_ind = i；
end

end

if （bits_alloc(max_ind) > 0)
bits_alloc（max_ind） = bits_alloc(max_ind） —1;
energy_alloc(max_ind） = energy_alloc（max_ind） - max_val；
bt = bt-1；
end

else

min_val = Inf;
min_ind = ceil（rand（1）＊num_subc)；
for i = 1：num_subc
if bits_alloc（i) ~=0 & bits_alloc（i) 〈9
temp = energytable(i，bits_alloc(i） + 1）；
else
temp = Inf；
end
if (temp < min_val)
min_val = temp;
min_ind = i;
end
end

if （bits_alloc(min_ind) 〈 8)
bits_alloc(min_ind) = bits_alloc(min_ind) +1;
if (min_val==inf)
min_val = energytable（min_ind，bits_alloc（min_ind）)；
end

energy_alloc（min_ind） = energy_alloc(min_ind) +min_val;
bt = bt+1;
end

end

end

for i = 1：length（bits_alloc)
if (mod（bits_alloc(i）,2） == 1 ＆ bits_alloc(i) ～=1)
[bits_alloc,energy_alloc］ = ResolvetheLastBit(bits_alloc,energy_alloc,i,energytable，num_subc）；
end
end

ResolvetheLastBit.m： An optimal bit—loading of the last bit requires a unique optimization.
function [bits_alloc， energy_alloc] = ResolvetheLastBit（bits_alloc,energy_alloc,index，energytable，num_subc）

max_val = 0;

for i = 1：num_subc
if （i ~= index ＆ bits_alloc（i) == 1）
if bits_alloc(i） ~= 0
temp = energytable（i，bits_alloc(i)） ；
end
if （temp > max_val)
max_val = temp；
max_ind = i;
end
end
end

min_val = Inf;
for i = 1:num_subc
if （i~= index & bits_alloc（i） == 1）
if bits_alloc（i） ～=0
temp = energytable（i，bits_alloc（i) + 1）；
end
if （temp 〈 min_val）
min_val = temp；
min_ind = i;
end
end
end

if (min_val 〈 max_val）
bits_alloc（min_ind) = bits_alloc(min_ind) + 1；
bits_alloc（index) = bits_alloc（index） - 1；
energy_alloc(index） = energy_alloc(index) — min_val；
else
bits_alloc(max_ind) = bits_alloc(max_ind） - 1；
bits_alloc（index) = bits_alloc(index) + 1;
energy_alloc(index） = energy_alloc（index） + max_val；
end