IS-95_仿真源程序

本 科 毕 业 设 计

IS-95前向链路通信系统的设计与仿真 (附录）

姓 名 李家蓬 学 院 信息与电气工程学院 专 业 电子信息工程专业 年 级 2008级 学 号 *********** 指导教师 王丽丽

2012年 5月 21日

IS-95系统仿真源程序

%main_IS95_forward.m

%此函数用于IS-95前向链路系统的仿真，包括扩 %频调制，匹配滤波，RAKE接收等相关通信模块。 %仿真环境: 加性高斯白噪声信道. %数据速率 = 9600 KBps %

clear all close all clc

disp('--------------start-------------------'); global Zi Zq Zs show R Gi Gq

clear j; show = 0;

SD = 0; % 选择软/硬判决接收

%-------------------主要的仿真参数设置------------------

BitRate = 9600;

ChipRate = 1228800; N = 184;

MFType = 1; % 匹配滤波器类型--升余弦 R = 5;

%+++++++++++++++++++Viterbi生成多项式++++++++++++++++++ G_Vit = [1 1 1 1 0 1 0 1 1; 1 0 1 1 1 0 0 0 1]; K = size(G_Vit, 2); L = size(G_Vit, 1);

%++++++++++++++++++++++++++++++++++++++++++++++++++++++

%++++++++++++++++++++++Walsh矩阵++++++++++++++++++++++++

WLen = ;

Walsh = reshape([1;0]*ones(1, WLen/2), WLen , 1); %Walsh = zeros(WLen ,1);

%++++++++++++++++++++++++++++++++++++++++++++++++++++++

%++++++++++++++++++扩频调制PN码的生成多项式++++++++++++++ %Gi = [ 1 0 1 0 0 0 1 1 1 0 1 0 0 0 0 1]'; %Gq = [ 1 0 0 1 1 1 0 0 0 1 1 1 1 0 0 1]';

Gi_ind = [15, 13, 9, 8, 7, 5, 0]';

Gq_ind = [15, 12, 11, 10, 6, 5, 4, 3, 0]';

Gi = zeros(16, 1);

Gi(16-Gi_ind) = ones(size(Gi_ind)); Zi = [zeros(length(Gi)-1, 1); 1]; % I路信道PN码生成器的初始状态

Gq = zeros(16, 1);

Gq(16-Gq_ind) = ones(size(Gq_ind)); Zq = [zeros(length(Gq)-1, 1); 1]; % Q路信道PN码生成器的初始状态

%++++++++++++++++++++++++++++++++++++++++++++++++++++++

%+++++++++++++++++++扰码生成多项式++++++++++++++++++++++

Gs_ind = [42, 35, 33, 31, 27, 26, 25, 22, 21, 19, 18, 17, 16, 10, 7, 6, 5, 3, 2, 1, 0]'; Gs = zeros(43, 1);

Gs(43-Gs_ind) = ones(size(Gs_ind)); Zs = [zeros(length(Gs)-1, 1); 1]; % 长序列生成器的初始状态

%++++++++++++++++++++++++++++++++++++++++++++++++++++++

%++++++++++++++++++++++AWGN信道++++++++++++++++++++++++ EbEc = 10*log10(ChipRate/BitRate); EbEcVit = 10*log10(L); EbNo = [-2 : 0.5 : 6.5]; %仿真信噪比范围(dB) % EbNo = [2 : 0.5 : 2.5];

%++++++++++++++++++++++++++++++++++++++++++++++++++++++

%------------------------------------------------------

%-------------------------主程序-------------------------

ErrorsB = []; ErrorsC = []; NN = []; if (SD == 1)

fprintf('\\n SOFT Decision Viterbi Decoder\\n\\n'); else

fprintf('\\n HARD Decision Viterbi Decoder\\n\\n'); end

for i=1:length(EbNo)

fprintf('\\nProcessing %1.1f (dB)', EbNo(i)); iter = 0; ErrB = 0; ErrC = 0;

while (ErrB <300) & (iter <150) drawnow;

%++++++++++++++++++++++发射机+++++++++++++++++++++++ TxData = (randn(N, 1)>0); % 速率为19.2Kcps

[TxChips, Scrambler] = PacketBuilder(TxData, G_Vit, Gs); % 速率为1.2288Mcps

[x PN MF] = Modulator(TxChips, MFType, Walsh);

%++++++++++++++++++++++++++++++++++++++++++++++++++++++

%++++++++++++++++++++++++信道+++++++++++++++++++++++++++

noise = 1/sqrt(2)*sqrt(R/2)*( randn(size(x)) + j*randn(size(x)))*10^(-(EbNo(i) EbEc)/20);

r = x+noise;

%++++++++++++++++++++++++++++++++++++++++++++++++++++++

%+++++++++++++++++++++++++接收机++++++++++++++++++++++++ RxSD = Demodulator(r, PN, MF, Walsh); %软判决，速率为19.2 Kcps RxHD = (RxSD>0); % 定义接收码片的硬判决 if (SD)

[RxData Metric]= ReceiverSD(RxSD, G_Vit, Scrambler); %软判决 else

[RxData Metric]= ReceiverHD(RxHD, G_Vit, Scrambler); %硬判决 end

%++++++++++++++++++++++++++++++++++++++++++++++++++++++

if(show)

subplot(311); plot(RxSD, '-o'); title('Soft Decisions');

subplot(312); plot(xor(TxChips, RxHD), '-o'); title('Chip Errors'); subplot(313); plot(xor(TxData, RxData), '-o');

title(['Data Bit Errors. Metric = ', num2str(Metric)]); pause; end

if(mod(iter, 50)==0) fprintf('.');

save TempResults ErrB ErrC N iter end

ErrB = ErrB + sum(xor(RxData, TxData)); ErrC = ErrC + sum(xor(RxHD, TxChips)); iter = iter+ 1; end

-

ErrorsB = [ErrorsB; ErrB]; ErrorsC = [ErrorsC; ErrC]; NN = [NN; N*iter]; save SimData * end

%+++++++++++++++++++++++++误码率计算++++++++++++++++++++++++ PerrB = ErrorsB./NN;

%PerrB1 = ErrorsB1./NN1; PerrC = ErrorsC./NN;

Pbpsk= 1/2*erfc(sqrt(10.^(EbNo/10)));

PcVit= 1/2*erfc(sqrt(10.^((EbNo-EbEcVit)/10))); Pc = 1/2*erfc(sqrt(10.^((EbNo-EbEc)/10)));

%+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

%%+++++++++++++++++++++++++性能仿真显示++++++++++++++++++++++ figure;

semilogy(EbNo(1:length(PerrB)), PerrB, 'b-*'); hold on;

% %semilogy(EbNo(1:length(PerrB1)), PerrB1, 'k-o'); hold on; % semilogy(EbNo(1:length(PerrC)), PerrC, 'b-o'); grid on; % semilogy(EbNo, Pbpsk, 'b-.^');

% %semilogy(EbNo, PcVit, 'k-.x'); ylabel('BER'); % semilogy(EbNo, Pc, 'b-.x'); xlabel('信噪比/dB'); ylabel('误码率'); grid on;

% legend('Pb of System (HD)', 'Pb of System (SD)', 'Pc before Viterbi of System', % ... 'Pb of BPSK with no Viterbi (theory)', 'Pc on Receiver (theory)'); % % %

% legend('Pb of System', 'Pc before Viterbi of System', ... %'Pb of BPSK with no Viterbi (theory)',

%'Pc before Viterbi (theory)', 'Pc on Receiver (theory)');

%+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

disp('--------------end-------------------');

%------------------------------------------------------

% ************************beginning of file***************************** %PacketBuilder.m

function [ChipsOut, Scrambler] = PacketBuilder(DataBits, G, Gs);

%

%此函数用于产生IS-95前向链路系统的发送数据包

%+++++++++++++++++++++++variables++++++++++++++++++++++++++++ % DataBits 发送数据（二进制形式） % G Viterbi编码生成多项式

% Gs 长序列生成多项式（扰码生成多项式） % ChipsOut 输入到调制器的码序列（二进制形式） % Scrambler 扰码

%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

global Zs

K = size(G, 2); L = size(G, 1);

N = *L*(length(DataBits)+K-1);% 码片数 (9.6 Kbps -> 1.288 Mbps)

chips = VitEnc(G, [DataBits; zeros(K-1,1)]); % Viterbi编码

% 交织编码

INTERL = reshape(chips, 24, 16); % IN：列, OUT：行 chips = reshape(INTERL', length(chips), 1); %速率=19.2 KBps

% 产生扰码

[LongSeq Zs] = PNGen(Gs, Zs, N); Scrambler = LongSeq(1::end);

ChipsOut = xor(chips, Scrambler);

%************************end of file***********************************

% ************************beginning of file***************************** %VitEnc.m

function y = VitEnc(G, x);

% 此函数根据生成多项式进行Viterbi编码 %

%+++++++++++++++++++++++variables++++++++++++++++++++++++++++ % G 生成多项式的矩阵

% x 输入数据（二进制形式） % y Viterbi编码输出序列

%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

K = size(G, 1); L = length(x);

yy = conv2(G, x'); yy = yy(:, 1:L);

y = reshape(yy,K*L, 1);

y = mod(y, 2);

% ************************end of file***********************************

% ************************beginning of file***************************** %PNGen.m

function [y, Z] = PNGen(G, Zin, N); %

% 此函数是根据生成多项式和输入状态产生长度为N的伪随机序列

%+++++++++++++++++++++++variables++++++++++++++++++++++++++++ % G 生成多项式

% Zin 移位寄存器初始化 % N PN序列长度

% y 生成的PN码序列 % Z 移位寄存器的输出状态

%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

L = length(G);

Z = Zin; % 移位寄存器的初始化

y = zeros(N, 1); for i=1:N

y(i) = Z(L);

Z = xor(G*Z(L), Z); Z = [Z(L); Z(1:L-1)]; end

%yy = filter(1, flipud(G), [1; zeros(N-1, 1)]); %yy = mod(yy, 2);

%************************end of file***********************************

% ************************beginning of file***************************** %Modulator.m

function [TxOut, PN, MF] = Modulator(chips, MFType, Walsh); %

%此函数用于实现IS-95前向链路系统的数据调制

%+++++++++++++++++++++++variables++++++++++++++++++++++++++++ % chips 发送的初始数据

% MFType 成型滤波器的类型选择 % Walsh walsh码

% TxOut 调制输出信号序列

% PN 用于扩频调制的PN码序列 % MF 匹配滤波器参数

%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

global Zi Zq show R Gi Gq

N = length(chips)*length(Walsh);

% 输入速率 = 19.2 KBps, 输出速率= 1.2288 Mcps tmp = sign(Walsh-1/2)*sign(chips'-1/2); chips = reshape(tmp, prod(size(tmp)), 1);

[PNi Zi] = PNGen(Gi, Zi, N); [PNq Zq] = PNGen(Gq, Zq, N);

PN = sign(PNi-1/2) + j*sign(PNq-1/2);

chips_out = chips.*PN;

chips = [chips_out, zeros(N, R-1)]; chips = reshape(chips.' , N*R, 1);

%成型滤波器 switch (MFType) case 1

%升余弦滤波器 L = 25;

L_2 = floor(L/2); n = [-L_2:L_2]; B = 0.7;

MF = sinc(n/R).*(cos(pi*B*n/R)./(1-(2*B*n/R).^2)); MF = MF/sqrt(sum(MF.^2)); case 2

%矩形滤波器 L = R;

L_2 = floor(L/2); MF = ones(L, 1);

MF = MF/sqrt(sum(MF.^2)); case 3

%汉明滤波器 L = R;

L_2 = floor(L/2); MF = hamming(L);

MF = MF/sqrt(sum(MF.^2)); end

MF = MF(:);

TxOut = sqrt(R)*conv(MF, chips)/sqrt(2); TxOut = TxOut(L_2+1: end - L_2);

if (show) figure;

subplot(211); plot(MF, '-o'); title('Matched Filter'); grid on; subplot(212); psd(TxOut, 1024, 1e3, 113); title('Spectrum'); end

% ************************end of file***********************************

% ************************beginning of file***************************** %Demodulator.m

function [SD] = Demodulator(RxIn, PN, MF, Walsh); %

% 此函数是实现基于RAKE接收机的IS-95前向信链路系统的数据包的解调

%+++++++++++++++++++++++variables++++++++++++++++++++++++++++ % RxIn 输入信号

% PN PN码序列（用于解扩） % MF 匹配滤波器参数 % Walsh 用于解调的walsh码

% SD RAKE接收机的软判决输出

%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

global R

N = length(RxIn)/R;

L = length(MF); L_2 = floor(L/2);

rr = conv(flipud(conj(MF)), RxIn);

rr = rr(L_2+1: end - L_2);

Rx = sign(real(rr(1:R:end))) + j*sign(imag(rr(1:R:end)));

Rx = reshape(Rx, , N/);

Walsh = ones(N/, 1)*sign(Walsh'-1/2); PN = reshape(PN, , N/)'; PN = PN.*Walsh;

% 输入速率 = 1.2288 Mpbs, 输出速率 = 19.2 KBps SD= PN*Rx; SD= real(diag(SD));

% ***********************end of file***********************************

% ************************beginning of file***************************** %ReceiverSD.m

function [DataOut, Metric] = ReceiverSD(SDchips, G, Scrambler); %

% 此函数用于实现基于Viterbi译码的发送数据的恢复

%+++++++++++++++++++++++variables++++++++++++++++++++++++++++ % SDchips 软判决RAKE接收机输入符号 % G Viterbi编码生成多项式矩阵 % Scrambler 扰码序列

% DataOut 接收数据（二进制形式） % Metric Viterbi译码最佳度量

%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

if (nargin == 1)

G = [1 1 1 1 0 1 0 1 1; 1 0 1 1 1 0 0 0 1]; end

% 速率=19.2 KBps

SDchips = SDchips.*sign(1/2-Scrambler);

INTERL = reshape(SDchips, 16, 24);

SDchips = reshape(INTERL', length(SDchips), 1); % 速率=19.2 KBps

[DataOut Metric] = SoftVitDec(G, SDchips, 1);

% ************************end of file**********************************

% ************************beginning of file***************************** %SoftVitDec.m

function [xx, BestMetric] = SoftVitDec(G, y, ZeroTail); %

% 此函数是实现软判决输入的Viterbi译码

%+++++++++++++++++++++++variables++++++++++++++++++++++++++++ % G 生成多项式的矩阵 % y 输入的待译码序列 % ZeroTail 判断是否包含‘0’尾 % xx Viterbi译码输出序列 % BestMetric 最后的最佳度量

%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

L = size(G, 1); % 输出码片数

K= size(G, 2); % 生成多项式的长度 N = 2^(K-1); % 状态数

T = length(y)/L; % 最大栅格深度

OutMtrx = zeros(N, 2*L); for s = 1:N

in0 = ones(L, 1)*[0, (dec2bin((s-1), (K-1))-'0')]; in1 = ones(L, 1)*[1, (dec2bin((s-1), (K-1))-'0')];

out0 = mod(sum((G.*in0)'), 2); out1 = mod(sum((G.*in1)'), 2);

OutMtrx(s, :) = [out0, out1]; end

OutMtrx = sign(OutMtrx-1/2);

PathMet = [100; zeros((N-1), 1)]; % 初始状态 = 100 PathMetTemp = PathMet(:,1);

Trellis = zeros(N, T); Trellis(:,1) = [0 : (N-1)]';

y = reshape(y, L, length(y)/L); for t = 1:T

yy = y(:, t); for s = 0:N/2-1

[B0 ind0] = max( PathMet(1+[2*s, 2*s+1]) + [OutMtrx(1+2*s, 0+[1:L]) * yy;

OutMtrx(1+(2*s+1), 0+[1:L])*yy] );

[B1 ind1] = max( PathMet(1+[2*s, 2*s+1]) + [OutMtrx(1+2*s, L+[1:L]) * yy; OutMtrx(1+(2*s+1), L+[1:L]) * yy] );

PathMetTemp(1+[s, s+N/2]) = [B0; B1];

Trellis(1+[s, s+N/2], t+1) = [2*s+(ind0-1); 2*s + (ind1-1)]; end

PathMet = PathMetTemp; end

xx = zeros(T, 1); if (ZeroTail) BestInd = 1; else

[Mycop, BestInd] = max(PathMet); end

BestMetric = PathMet(BestInd); xx(T) = floor((BestInd-1)/(N/2));

NextState = Trellis(BestInd, (T+1)); for t=T:-1:2

xx(t-1) = floor(NextState/(N/2));

NextState = Trellis( (NextState+1), t); end

if (ZeroTail)

xx = xx(1:end-K+1); end

% ************************end of file***********************************

% ************************beginning of file***************************** %ReceiverHD.m

function [DataOut, Metric] = ReceiverHD(HDchips, G, Scrambler); %

% 此函数用于实现基于Viterbi译码的硬判决接收机

%+++++++++++++++++++++++variables++++++++++++++++++++++++++++ % SDchips 硬判决RAKE接收机输入符号 % G Viterbi编码生成多项式矩阵 % Scrambler 扰码序列

% DataOut 接收数据（二进制形式）

% Metric Viterbi译码最佳度量

%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

if (nargin == 1)

G = [1 1 1 1 0 1 0 1 1; 1 0 1 1 1 0 0 0 1]; end

% 速率=19.2 KBps

HDchips = xor(HDchips, Scrambler);

INTERL = reshape(HDchips, 16, 24);

HDchips = reshape(INTERL', length(HDchips), 1);

[DataOut Metric] = VitDec(G, HDchips, 1);

%************************end of file***********************************

% ************************beginning of file***************************** %VitDec.m

function [xx, BestMetric] = VitDec(G, y, ZeroTail); %

% 此函数是实现硬判决输入的Viterbi译码

%+++++++++++++++++++++++variables++++++++++++++++++++++++++++ % G 生成多项式的矩阵 % y 输入的待译码序列 % ZeroTail 判断是否包含‘0’尾 % xx Viterbi译码输出序列 % BestMetric 最后的最佳度量

%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

L = size(G, 1); % 输出码片数 K= size(G, 2); % 生成多项式长度 N = 2^(K-1); % 状态数

T = length(y)/L; % 最大栅格深度

OutMtrx = zeros(N, 2*L); for s = 1:N

in0 = ones(L, 1)*[0, (dec2bin((s-1), (K-1))-'0')]; in1 = ones(L, 1)*[1, (dec2bin((s-1), (K-1))-'0')];

out0 = mod(sum((G.*in0)'), 2); out1 = mod(sum((G.*in1)'), 2);

OutMtrx(s, :) = [out0, out1]; end

PathMet = [0; 100*ones((N-1), 1)]; PathMetTemp = PathMet(:,1);

Trellis = zeros(N, T); Trellis(:,1) = [0 : (N-1)]';

y = reshape(y, L, length(y)/L); for t = 1:T

yy = y(:, t)'; for s = 0:N/2-1

[B0 ind0] = min( PathMet(1+[2*s, 2*s+1]) + [sum(abs(OutMtrx(1+2*s, 0+[1:L]) - yy).^2); sum(abs(OutMtrx(1+(2*s+1), 0+[1:L]) - yy).^2)] );

[B1 ind1] = min( PathMet(1+[2*s, 2*s+1]) + [sum(abs(OutMtrx(1+2*s, L+[1:L]) - yy).^2); sum(abs(OutMtrx(1+(2*s+1), L+[1:L]) - yy).^2)] );

PathMetTemp(1+[s, s+N/2]) = [B0; B1];

Trellis(1+[s, s+N/2], t+1) = [2*s+(ind0-1); 2*s + (ind1-1)]; end

PathMet = PathMetTemp; end

xx = zeros(T, 1); if (ZeroTail) BestInd = 1; else

[Mycop, BestInd] = min(PathMet); end

BestMetric = PathMet(BestInd); xx(T) = floor((BestInd-1)/(N/2));

NextState = Trellis(BestInd, (T+1)); for t=T:-1:2

xx(t-1) = floor(NextState/(N/2));

NextState = Trellis( (NextState+1), t); end

if (ZeroTail)

xx = xx(1:end-K+1); end

% ************************end of file***********************************