ofdm.m

基于训练序列的OFDM的信道估计

echo off;                                     % 关闭回显 
clear all;                                    % 从内存中清除变量和函数
close all;                                    % 关闭所有图形
fprintf( '\n OFDM仿真\n \n') ;                % 设置显示格式
% --------------------------------------------- %
%                   参数定义                    %
% --------------------------------------------- %
IFFT_bin_length = 1024;                                % 发送端的IFFT变换长度, 接收端的FFT变换长度，R代表接受，T代表发送
carrier_count = 200;                                   % 子载波数
bits_per_symbol = 2;                                   % 位数/符号
symbols_per_carrier = 50;                              % 符号数/载波
cp_length = input('cp length = ');                     % 输入循环前缀长度   
d4 = input('d4 = ');                                   % 输入最大多径时延扩展
a4 = input('a4 = ');                                   % 输入最大多径时延扩展的系数
SNR = input('SNR = ');                                 % 输入信道信噪比(dB)                        
% --------------------------------------------- %
%                初始参数设置完毕               %
% --------------------------------------------- %
baseband_out_length = carrier_count * symbols_per_carrier * bits_per_symbol;        % 计算发送的二进制序列长度：基带传送长度=载波数× 符号数/载波×位数/符号
carriers = (1:carrier_count) + (floor(IFFT_bin_length/4) - floor(carrier_count/2));  % 载波（坐标）：(1:200) + 256 - 100 = 157:356
conjugate_carriers = IFFT_bin_length - carriers + 2;                                 % 载波变换（坐标）：1024 - (157:356) + 2 = 1026 - (157:356) = (869:670) 
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% 构造共轭时间－载波矩阵，以便应用所谓的RCC，Reduced Computational Complexity算法，即ifft之后结果为实数 
% 也可以用flipdim函数构造对称共轭矩阵
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% --------------------------------------------- %
%                    发送信号                   %
% --------------------------------------------- %
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% Generate a random binary output signal:
%   - a row of uniform random numbers (between 0 and 1), rounded to 0 or 1
%   - this will be the baseband signal which is to be transmitted.
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
baseband_out = round(rand(1,baseband_out_length));
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% round：朝最近的整数取整，rand:产生均匀分布的随机数矩阵（1×baseband_out_length阶）               
% Convert to 'modulo N' integers where N = 2^bits_per_symbol
%   - this defines how many states each symbol can represent
%   - first, make a matrix with each column representing consecutive bits 
%     from the input stream and the number of bits in a column equal to the
%     number of bits per symbol
%   - then, for each column, multiply each row value by the power of 2 that 
%     it represents and add all the rows
%   - for example:  input 0 1 1 0 0 0 1 1 1 0
%                   bits_per_symbol = 2
%                   convert_matrix = 0 1 0 1 1
%                                    1 0 0 1 0
%                   modulo_baseband = 1 2 0 3 2
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
convert_matrix = reshape(baseband_out,bits_per_symbol,length(baseband_out)/bits_per_symbol) ;
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% RESHAPE Change size. 把baseband_out变为M×N阶的矩阵
%    RESHAPE(X,M,N) returns the M-by-N matrix whose elements
%    are taken columnwise from X.  An error results if X does
%    not have M*N elements.
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
for k = 1:(length(baseband_out)/bits_per_symbol)
  modulo_baseband(k) = 0 ;
for i = 1:bits_per_symbol
  modulo_baseband(k) =  modulo_baseband(k) + convert_matrix(i,k)*2^(bits_per_symbol - i) ;
end
end
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% Serial to Parallel Conversion       串并转换
%   - convert the serial modulo N stream into a matrix where each column 
%     represents a carrier and each row represents a symbol
%   - for example:
%      serial input stream = a b c d e f g h i j k l m n o p
%      parallel carrier distribution =
%                           C1/s1=a  C2/s1=b  C3/s1=c C4/s1=d
%                           C1/s2=e  C2/s2=f  C3/s2=g C4/s2=h
%                           C1/s3=i  C2/s3=j  C3/s3=k C4/s3=l
%                              .        .        .       .
%                              .        .        .       .
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
carrier_matrix = reshape(modulo_baseband, carrier_count, symbols_per_carrier)';          % 生成时间－载波矩阵（carrier_count×symbols_per_carrier）
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% Apply differential coding to each carrier string
%   - append an arbitrary start symbol (let it be 0, that works for all 
%     values of bits_per_symbol) (note that this is done using a vertical
%     concatenation [x;y] of a row of zeros with the carrier matrix, sweet!)
%   - perform modulo N addition between symbol(n) and symbol(n-1) to get the 
%     coded value of symbol(n)
%   - for example:
%                   bits_per_symbol = 2 (modulo 4)
%                   symbol stream =  3 2 1 0 2 3
%                   start symbol = 0
%
%                   coded symbols = 0 + 3 = 3
%                                   3 + 2 = 11 = 1
%                                   1 + 1 = 2
%                                   2 + 0 = 2
%                                   2 + 2 = 10 = 0
%                                   0 + 3 = 3
%
%                   coded stream =  0 3 1 2 2 0 3
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% --------------------------------------------- %
%                   QDPSK调制                   %
% --------------------------------------------- %
carrier_matrix = [zeros(1,carrier_count);carrier_matrix];                                    % 添加一个差分调制的初始相位，为0
for i = 2:(symbols_per_carrier + 1)
    carrier_matrix(i,:) = rem(carrier_matrix(i,:)+carrier_matrix(i-1,:),2^bits_per_symbol);  % 差分调制（rem除后取余）
end
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% Convert the differential coding into a phase
%   - each phase represents a different state of the symbol
%   - for example:
%                 bits_per_symbol = 2 (modulo 4)
%                 symbols = 0 3 2 1
%                 phases =
%                           0 * 2pi/4 = 0 (0 degrees)
%                           3 * 2pi/4 = 3pi/2 (270 degrees)
%                           2 * 2pi/4 = pi (180 degrees)
%                           1 * 2pi/4 = pi/2 (90 degrees)
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
carrier_matrix = carrier_matrix * ((2*pi)/(2^bits_per_symbol));                              % 产生差分相位
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% Convert the phase to a complex number
%   - each symbol is given a magnitude of 1 to go along with its phase 
%     (via the ones(r,c) function)
%   - it is then converted from polar to cartesian (complex) form
%   - the result is 2 matrices, X with the real values and Y with the imaginary
%   - each X column has all the real values for a carrier, and each Y column 
%     has the imaginary values for a carrier
%   - a single complex matrix is then generated taking X for real and 
%     Y for imaginary
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
[X,Y] = pol2cart(carrier_matrix, ones(size(carrier_matrix,1),size(carrier_matrix,2)));       % 由极坐标向复数坐标转化 第一参数为相位 第二参数为幅度
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% Carrier_matrix contains all the phase information and all the amplitudes are the same,‘1’     
% 函数说明：极或柱坐标变为直角坐标
%  POL2CART Transform polar to Cartesian coordinates.
%    [X,Y] = POL2CART(TH,R) transforms corresponding elements of data
%    stored in polar coordinates (angle TH, radius R) to Cartesian
%    coordinates X,Y.  The arrays TH and R must the same size (or
%    either can be scalar).  TH must be in radians.
% [X,Y,Z] = POL2CART(TH,R,Z) transforms corresponding elements of
%    data stored in cylindrical coordinates (angle TH, radius R, height Z) 
%    to Cartesian coordinates X,Y,Z. The arrays TH, R, and Z must be
%    the same size (or any of them can be scalar).  TH must be in radians.
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
complex_carrier_matrix = complex(X,Y);                     
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% 函数说明：
%  COMPLEX Construct complex result from real and imaginary parts.
%    C = COMPLEX(A,B) returns the complex result A + Bi, where A and B are
%    identically sized real N-D arrays, matrices, or scalars of the same data type.
%
% Assign each carrier to its IFFT bin
%   - each row of complex_carrier_matrix represents one symbol period, with
%     a symbol for each carrier
%   - a matrix is generated to represent all IFFT bins (columns) and all 
%     symbols (rows)
%   - the phase modulation for each carrier is then assigned to the 
%     appropriate bin
%   - the conjugate of the phase modulation is then assigned to the 
%     appropriate bin
%      - the phase modulation bins and their conjugates are symmetric about 
%        the Nyquist frequency in the IFFT bins
%      - since the first bin is DC, the Nyquist Frequency is located 
%        at (number of bins/2) + 1
%      - symmetric conjugates are generated so that when the signal is 
%        transformed to the time domain, the time signal will be real-valued
%   - example
%      - 1024 IFFT bins
%      - bin 513 is the center (symmetry point)
%      - bin 1 is DC
%      - bin 514 is the complex conjugate of bin 512
%      - bin 515 is the complex conjugate of bin 511
%      - ....
%      - bin 1024 is the complex conjugate of bin 2 (if all bins 
%        were used as carriers)
%      - So, bins 2-512 map to bins 1024-514
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－ %
%                               添加训练序列                                  %
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－ % 
training_symbols = [ 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 ...
-j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 ...
1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 ...
-1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j ...
-1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 ]; 
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% 25 times "1 j j 1" 
% 25 times "-1 -j -j -1"
% totally 200 symbols as a row
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
training_symbols = cat(1, training_symbols, training_symbols) ;  
training_symbols = cat(1, training_symbols, training_symbols) ;   
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% Production of 4 rows of training_symbols
% 函数说明：串接成高维数组
% CAT Concatenate arrays.
%     CAT(DIM,A,B) concatenates the arrays A and B along the dimension DIM.  
%     CAT(2,A,B) is the same as [A,B].
%     CAT(1,A,B) is the same as [A;B].
%     B = CAT(DIM,A1,A2,A3,A4,...) concatenates the input
%     arrays A1, A2, etc. along the dimension DIM.
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
complex_carrier_matrix = cat(1, training_symbols, complex_carrier_matrix);
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% 训练序列与数据合并，串接成高维数组
% block-type pilot symbols
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
IFFT_modulation = zeros(4 + symbols_per_carrier + 1, IFFT_bin_length);
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
% Here a row vector of zeros is between training symbols and data symbols!!! 
% 4 training symbols and 1 zero symbol
% every OFDM symbol takes a row of "IFFT_modulation" 
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
IFFT_modulation(:,carriers) = complex_carrier_matrix;
IFFT_modulation(:,conjugate_carriers) = conj(complex_carrier_matrix);
%－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－             
% Find the indices of zeros
% Description