Matlab簡單的OFDM仿真,信道估計,有BER曲線

1、clear all;close all;fprin tf( 'n OFDM仿真 n n');% % 參數(shù)定義 % %IFFT_bin_length = 1024;carrier_count = 200;bits_per_symbol = 2;symbols_per_carrier = 50;% 子載波數(shù)200% 位數(shù) / 符號2% 符號數(shù) / 載波 50% 訓(xùn)練 符號數(shù) 10% 循環(huán)前綴長度 T/4 (作者注明) All-zero CP% 調(diào)制方式QDPSK% 多徑信道數(shù) 2、 3、 4(缺省)% 信道最大時延 7 ( 單位數(shù)據(jù) 符號)% 仿真條件收發(fā)之間嚴(yán)格同步%SNR=in

2、put('SNR='); % 輸入 信噪比參數(shù)SNR=3:14;淀 義信噪比范圍BER=zeros(1,length(SNR);baseband_out_length = carrier_count * symbols_per_carrier * bits_per_symbol;% 計算發(fā)送的二進(jìn)制序列長度carriers =(1:carrier_count) + (floor(IFFT_bin_length/4)-floor(carrier_count/2); % 坐標(biāo)： (1 to 200) + 156 , 157 - 356 conjugate_carriers=IFFT

3、_bin_length-carriers+2; %坐 標(biāo) ： 1024 -(157:356) + 2 = 1026 - (157:356) = (869:670)% 構(gòu)造共軛 時間 載波 矩陣 ，以便 應(yīng)用 所謂的 RCC，Reduced ComputationalComplexity 算法，即ifft之后結(jié)果為實數(shù)% Define the conjugate time-carrier matrix% 也可以用 flipdim 函數(shù) 構(gòu)造對稱共軛矩陣% %信號 發(fā)射% % %out = rand(1,baseband_out_length);%baseband_out1 = round(out

4、) ;%baseband_out2 = floor(out*2) ;%baseband_out3 = ceil(out*2)-1 ;%baseband_out4 = randint(1,baseband_out_length);% 四種生成發(fā)送的二進(jìn)制序列的方法，任取一種產(chǎn)生要發(fā)送的二進(jìn)制序列%if (baseband_out1 = baseband_out2 & baseband_out1 = baseband_out3 )% fprintf('Transmission Sequence Generated n n');% baseband_out = baseban

5、d_out1 ;%else% fprintf('Check Code! n n');%end% 驗證四種生成發(fā)送的二進(jìn)制序列的方法baseband_out=round( rand(1,baseband_out_length); convert_matrix reshape(baseband_out,bits_per_symbol,length(baseband_out)/bits_per_sy mbol);for k = 1length(baseband_out)/bits_per_symbol),modulo_baseband(k) = 0;for i = 1:bits_pe

6、r_symbolmodulo_baseband(k) = modulo_baseband(k) + convert_matrix(i,k)*2A(bits_per_symbol - i);endend% 每2個比特轉(zhuǎn)化為整數(shù) 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_

7、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');% 函數(shù)說明：% BI2DE Convert binary vectors to decimal numbers.% D = BI2DE(B) converts a bina

8、ry vector B to a decimal value D. When Bis% 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

9、_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 &q

10、uot;modulo_baseband = Test_convert_matrix1".% carrier_matrix = reshape(modulo_baseband,carrier_count,symbols_per_carrier)'% 生成時間載波矩陣% % QDPSK 調(diào)制 % %carrier_matrix = zeros(1,carrier_count); carrier_matrix; % 添加一個 差分調(diào)制的初始相位，為 0for i = 2symbols_per_carrier + 1)carrier_matrix(i, = rem(carrier_m

11、atrix(i, + carrier_matrix (i-1,2Abits_per_symbol) ; % 差分調(diào)制endcarrier_matrix = carrier_matrix*(2*pi)/(2Abits_per_symbol); % 產(chǎn)生差分相位X,Y=pol2cart(carrier_matrix,ones(size(carrier_matrix,1),size(carrier_matrix,2); %由極坐標(biāo)向復(fù)數(shù)坐標(biāo)轉(zhuǎn)化 第一參數(shù)為相位 第二參數(shù)為幅度%Carrier_matrix contains all the phase information and all the

12、amplitudes are the same 1'.complex_carrier_matrix = complex(X, Y) ;%添加訓(xùn)練序列'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 -11 j j 1 -1 .-j -j -1 1 j j 1 -1 -j -j -11 j j 1 -1 -j-j -1 1 j j 1-1 -j -j -1 1j j 1 -1 -j -j -1 1 j j 1 -1 -j -j -1 .1

13、j j 1 -1 -j -j -1 1 j j 1-1 -j -j -1 1jj 1 -1 -j -j-1 1 j j 1 -1-j -j -1 1 j j 1 -1 -j -j -11 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 -11 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 t

14、imes "1 j j 1" , 25 times "-1 -j -j -1", totally 200symbols 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

15、,complex_carrier_matrix) ; % 訓(xùn)練序列與數(shù)據(jù)合并% 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_modulati

16、on" 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

17、_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(I

18、FFT_modulation') ; %進(jìn)行 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 + 1windowed_time_wave_matrix( i, : ) = real(time_wave_matrix( i, : ) ; end% get the real part of t

19、he result of IFFT% 這一步可以省略，因為 IFFT 結(jié)果都是實數(shù)% 由此可以看出，只是取了 IFFT 之后載波上的點，并未進(jìn)行 CP 的復(fù)制和添加 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_

20、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 = of

21、dm_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) ;end copy2 = zeros(size(Tx_data) ) ;for i = 1 + d2: length( Tx_data) co

22、py2(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-paths Tx_s

23、ignal_power = var(Tx_data);for idx=1:length(SNR)%monte carlo 仿真模擬lin ear_SNR = 10A( 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(R

24、x_Data, IFFT_bin_length, symbols_per_carrier + 1) ;Rx_spectrum = fft(Rx_Data_matrix) ;% Suppose precise synchronazition between Tx and Rx Rx_carriers = Rx_spectrum( carriers, : )'Rx_training_symbols = Rx_carriers( (1: 4) , : ) ;Rx_carriers = Rx_carriers(5: 55), : ) ;4 +%信道估計Rx_training_symbols =

25、 Rx_training_symbols./ training_symbols;Rx_trai nin g_symbols_de no = Rx_trai nin g_symbols42;Rx_training_symbols_denoRx_training_symbols_deno(1,+Rx_training_symbols_deno(2,+Rx_training_s ymbols_deno(3,+Rx_training_symbols_deno(4, ;Rx_training_symbols_nume = Rx_training_symbols(1,:)+Rx_training_symb

26、ols(2, :)+ Rx_training_symbols(3, :)+Rx_training_symbols(4, : ) ;Rx_training_symbols_nume = conj(Rx_training_symbols_nume) ;% 取 4 個向量的導(dǎo)頻符號是為了進(jìn)行平均 優(yōu)化%都是針對 “行向量”即單個的OFDM符號 進(jìn)行操作%原理：尋求1/H，對FFT之后的數(shù)據(jù)進(jìn)行頻域補償 % 1/H = conj(H)/HA2 because 屮2 = H * conj(H) Rx_training_symbols Rx_training_symbols_nume./Rx_traini

27、ng_symbols_deno;Rx_training_symbolsRx_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_

28、4) ;Rx_training_symbols_8,Rx_training_symbols_16,Rx_training_symbols_32,Rx_training_symbols_48,cat(1,進(jìn)行頻域單抽頭均衡Rx_training_symbols_16= cat(1,Rx_training_symbols_8) ;Rx_training_symbols_32= cat(1,Rx_training_symbols_16) ;Rx_training_symbols_48=cat(1,Rx_training_symbols_16) ;Rx_training_symbols_50=cat(

29、1,Rx_training_symbols_2) ;Rx_training_symbols=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) ;Using%TestofIIrem"%Rx_phase1 = Rx_phase;%Rx_phase2 = Rx_phase;%Rx_phase1(p

30、hase_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('nThere is no need using rem in negative phase transition.n')%else% fprintf('nWeneed to use rem in

31、negative phase transition.n')%end%Rx_phase(phase_negative) = rem(Rx_phase(phase_negative) + 360, 360) ; % 把負(fù)的相位轉(zhuǎn)化為正的相位Rx_decoded_phase = diff(Rx_phase) ;% 這也是為什么要在前面加上初始相位的原因% “ Here a row vector of zeros is between training symbols and datasymbols! ”phase_negative = find(Rx_decoded_phase < 0

32、) ;Rx_decoded_phase(phase_negative)= rem(Rx_decoded_phase(phase_negative) + 360, 360) ; %再次把負(fù)的相位轉(zhuǎn)化為正的相位% % QDPSK 解調(diào) % % base_phase = 360 /2Abits_per_symbol;delta_phase = base_phase /2;Rx_decoded_symbols zeros(size(Rx_decoded_phase,1),size(Rx_decoded_phase,2) ; for i = 1: (2Abits_per_symbol - 1)center_phase = base_phase*i;plus_delta = center_phase + delta_phase; % Decision threshold 1minus_delta = center_phase - delta_phase; % Decision threshold 2 decoded = find(Rx_decoded_phase <= plus_delta)&(Rx_decoded_phase > minus_delta) ;Rx_decoded_symbols(decoded) = i;end% 僅僅對