Matlab简单的OFDM仿真,信道估计,有BER曲线.doc

clear all;close all;fprintf( n OFDM仿真n n) ;% - % 参数定义 % - %IFFT_bin_length = 1024;carrier_count = 200;bits_per_symbol = 2;symbols_per_carrier = 50;% 子载波数 200% 位数/ 符号 2% 符号数/ 载波 50% 训练符号数 10% 循环前缀长度 T/4（作者注明）All-zero CP% 调制方式 QDPSK% 多径信道数 2、3、4（缺省）% 信道最大时延 7 (单位数据符号)% 仿真条件 收发之间严格同步 %SNR=input(SNR=); % 输入信噪比参数SNR=3:14;%定义信噪比范围BER=zeros(1,length(SNR);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 to 200) + 156 ,157 - 356conjugate_carriers=IFFT_bin_length-carriers+2;% 坐标 ：1024 - (157:356) + 2 = 1026 - (157:356) = (869:670) % 构造共轭时间载波矩阵，以便应用所谓的RCC，Reduced Computational Complexity算法，即ifft之后结果为实数 % Define the conjugate time-carrier matrix% 也可以用flipdim函数构造对称共轭矩阵% - % 信号发射 % - %out = rand(1,baseband_out_length);%baseband_out1 = round(out) ;%baseband_out2 = floor(out*2) ;%baseband_out3 = ceil(out*2)-1 ;%baseband_out4 = randint(1,baseband_out_length);% 四种生成发送的二进制序列的方法，任取一种产生要发送的二进制序列%if (baseband_out1 = baseband_out2 & baseband_out1 = baseband_out3 )% fprintf(Transmission Sequence Generated n n);% baseband_out = baseband_out1 ;%else % fprintf(Check Code! n n);%end% 验证四种生成发送的二进制序列的方法baseband_out=round( rand(1,baseband_out_length);convert_matrix = reshape(baseband_out,bits_per_symbol,length(baseband_out)/bits_per_symbol);for k = 1length(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% 每2个比特转化为整数 0至3% 采用left-msb方式%-%Test by lavabin%A built-in function of directly change binary bits into decimal numbers%-%convert_matrix1 = zeros(length(baseband_out)/bits_per_symbol,bits_per_symbol);%convert_matrix1 = convert_matrix ;%Test_convert_matrix1 = bi2de(convert_matrix1,bits_per_symbol,left-msb);%Test_convert_matrix2 = bi2de(convert_matrix1,bits_per_symbol,right-msb);% 函数说明：% BI2DE Convert binary vectors to decimal numbers.% D = BI2DE(B) converts a binary vector B to a decimal value D. When B is% a matrix, the conversion is performed row-wise and the output D is a% column vector of decimal values. The default orientation of thebinary% input is Right-MSB; the first element in B represents the least significant bit.%if (modulo_baseband = Test_convert_matrix1) % fprintf(modulo_baseband = Test_convert_matrix1 nnn);%else if (modulo_baseband = Test_convert_matrix2) % fprintf(modulo_baseband = Test_convert_matrix2 nnn);% else% fprintf(modulo_baseband = any Test_convert_matrix nnn); % end%end% we get the result modulo_baseband = Test_convert_matrix1.%-carrier_matrix = reshape(modulo_baseband,carrier_count,symbols_per_carrier);% 生成时间载波矩阵% - % QDPSK调制 % - %carrier_matrix = zeros(1,carrier_count); carrier_matrix;% 添加一个差分调制的初始相位，为0for i = 2symbols_per_carrier + 1) carrier_matrix(i, = rem(carrier_matrix(i, + carrier_matrix (i-1, 2bits_per_symbol) ;% 差分调制 endcarrier_matrix = carrier_matrix*(2*pi)/(2bits_per_symbol) ;% 产生差分相位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 same1. complex_carrier_matrix = complex(X, Y) ;% 添加训练序列 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 rowtraining_symbols = cat(1, training_symbols, training_symbols) ;training_symbols = cat(1, training_symbols, training_symbols) ; % Production of 4 rows of training_symbolscomplex_carrier_matrix = cat(1, training_symbols, complex_carrier_matrix) ; % 训练序列与数据合并 % block-type pilot symbolsIFFT_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) ;%-% Test by lavabin- Find the indices of zeros %index_of_zeros = zeros(symbols_per_carrier,IFFT_bin_length - 2*carrier_count);%IFFT_modulation1 = zeros(4 + symbols_per_carrier + 1,IFFT_bin_length);%IFFT_modulation2 = zeros(4 + symbols_per_carrier + 1,IFFT_bin_length);%IFFT_modulation1(6:symbols_per_carrier+5, = IFFT_modulation(6:symbols_per_carrier+5,=0 ;%for i = 1:symbols_per_carrier%index_of_zeros(i, = find(IFFT_modulation1(i+5,=1);%end%-time_wave_matrix = ifft(IFFT_modulation) ; % 进行IFFT操作time_wave_matrix = time_wave_matrix;% If X is a matrix, ifft returns the inverse Fourier transform of each column of the matrix.for i = 1: 4 + symbols_per_carrier + 1 windowed_time_wave_matrix( i, : ) = real(time_wave_matrix( i, : ) ;end% get the real part of the result of IFFT% 这一步可以省略，因为IFFT结果都是实数% 由此可以看出，只是取了IFFT之后载波上的点，并未进行CP的复制和添加endofdm_modulation = reshape(windowed_time_wave_matrix,1, IFFT_bin_length*(4 + symbols_per_carrier + 1) ) ;% P2S operation%-% Test by lavabin% Another way of matrix transition%ofdm_modulation_tmp = windowed_time_wave_matrix.;%ofdm_modulation_test = ofdm_modulation_tmp(;%if (ofdm_modulation_test = ofdm_modulation)% fprintf(ofdm_modulation_test = ofdm_modulation nnn);%else%fprintf(ofdm_modulation_test = ofdm_modulation nnn);%end % We get the result ofdm_modulation_test = ofdm_modulation .%-Tx_data=ofdm_modulation;% - % 信道模拟 % - %d1= 4; a1 = 0.2; d2 = 5; a2 = 0.3; d3 = 6; a3 = 0.4; d4 = 7; a4 = 0.5;%信道模拟copy1 = zeros(size(Tx_data) ;for i = 1 + d1: length(Tx_data)copy1(i) = a1*Tx_data( i - d1) ;endcopy2 = zeros(size(Tx_data) ) ;for i = 1 + d2: length( Tx_data)copy2(i) = a2*Tx_data( i - d2) ;endcopy3 = zeros(size(Tx_data) ) ;for i = 1 + d3: length(Tx_data)copy3(i) = a3*Tx_data ( i - d3) ;endcopy4 = zeros(size(Tx_data) ) ;for i = 1 + d4: length( Tx_data)copy4(i) = a4*Tx_data(i - d4) ;endTx_data = Tx_data + copy1 + copy2 + copy3 + copy4; % 4 multi-pathsTx_signal_power = var(Tx_data);for idx=1:length(SNR)%monte carlo 仿真模拟 linear_SNR = 10( SNR(idx) /10) ;noise_sigma = Tx_signal_power / linear_SNR;noise_scale_factor = sqrt(noise_sigma) ;noise = randn(1, length(Tx_data) )*noise_scale_factor;Rx_Data = Tx_data + noise;% - % 信号接收 % - % Rx_Data_matrix = reshape(Rx_Data, IFFT_bin_length, 4 + symbols_per_carrier + 1) ;Rx_spectrum = fft(Rx_Data_matrix) ; %Suppose precise synchronazition between Tx and RxRx_carriers = Rx_spectrum( carriers, : );Rx_training_symbols = Rx_carriers( (1: 4) , : ) ;Rx_carriers = Rx_carriers(5: 55), : ) ;% - % 信道估计 % - % Rx_training_symbols = Rx_training_symbols./ training_symbols;Rx_training_symbols_deno = Rx_training_symbols.2;Rx_training_symbols_deno = Rx_training_symbols_deno(1,+Rx_training_symbols_deno(2,+Rx_training_symbols_deno(3,+Rx_training_symbols_deno(4, ;Rx_training_symbols_nume = Rx_training_symbols(1, : ) +Rx_training_symbols(2, : ) + Rx_training_symbols(3, : ) +Rx_training_symbols(4, : ) ;Rx_training_symbols_nume = conj(Rx_training_symbols_nume) ;% 取4个向量的导频符号是为了进行平均优化% 都是针对 “行向量”即单个的OFDM符号 进行操作% 原理：寻求1/H，对FFT之后的数据进行频域补偿% 1/H = conj(H)/H2 because H2 = H * conj(H) Rx_training_symbols = Rx_training_symbols_nume./Rx_training_symbols_deno;Rx_training_symbols = Rx_training_symbols_nume./Rx_training_symbols_deno;Rx_training_symbols_2 = cat(1, Rx_training_symbols,Rx_training_symbols) ;Rx_training_symbols_4 = cat(1, Rx_training_symbols_2,Rx_training_symbols_2) ;Rx_training_symbols_8 = cat(1, Rx_training_symbols_4,Rx_training_symbols_4) ;Rx_training_symbols_16 = cat(1, Rx_training_symbols_8, Rx_training_symbols_8) ;Rx_training_symbols_32 = cat(1, Rx_training_symbols_16, Rx_training_symbols_16) ;Rx_training_symbols_48 = cat(1, Rx_training_symbols_32, Rx_training_symbols_16) ;Rx_training_symbols_50 = cat(1, Rx_training_symbols_48, Rx_training_symbols_2) ;Rx_training_symbols = cat(1, Rx_training_symbols_50,Rx_training_symbols) ;Rx_carriers = Rx_training_symbols.*Rx_carriers; % 进行频域单抽头均衡 Rx_phase = angle(Rx_carriers)*(180/pi) ;phase_negative = find(Rx_phase 0) ;%-Test of Using rem-%Rx_phase1 = Rx_phase; %Rx_phase2 = Rx_phase;%Rx_phase1(phase_negative) = rem(Rx_phase1(phase_negative) + 360, 360) ;%Rx_phase2(phase_negative) = Rx_phase2(phase_negative) + 360 ;%if Rx_phase2(phase_negative) = Rx_phase1(phase_negative)%fprintf(n There is no need using rem in negative phase transition.n)%else% fprintf(n We need to use rem in negative phase transition.n) %end%-Rx_phase(phase_negative) = rem(Rx_phase(phase_negative) + 360, 360) ;% 把负的相位转化为正的相位Rx_decoded_phase = diff(Rx_phase) ;%这也是为什么要在前面加上初始相位的原因 % “Here a row vector of zeros is between training symbols and data symbols!”phase_negative = find(Rx_decoded_phase 0) ;Rx_decoded_phase(phase_negative)= rem(Rx_decoded_phase(phase_negative) + 360, 360) ;% 再次把负的相位转化为正的相位% - % QDPSK解调 % - % base_phase = 360 /2bits_per_symbol;delta_phase = base_phase /2;Rx_decoded_symbols = zeros(size(Rx_decoded_phase,1),size(Rx_decoded_phase,2) ;for i = 1: (2bits_per_symbol - 1)center_phase = base_phase*i;plus_delta = center_phase + delta_phase;% Decision threshold 1minus_delta = center_phase - delta_phase; % Decision threshold 2decoded =