无线通信系统信道仿真源码

function signal=bin(N)

signal=rand(1,N)>=0.5;

end

//BpskModu.m

function Signal=BpskModu(signal)%信号波形调制

t=linspace(0,0.000001,32);%把0.000001分成32等分

fc=1000000000;

carrier=sin(2*pi*fc*t);

signal=signal*2-1;

Signal=[];

for i=1:length(signal)

Signal=[Signal,signal(i)*carrier];

end

figure(2);

plot(Signal)

title('调制后的波形')

end

//demobpsk.m

function bitstream=demobpsk(gussianhannel,N)

load pass1

load pass2

%gussianhannel=filter(bandpass,1,gussianhannel);

%figure

%plot(gussianhannel);%滤除噪声后的波形

%title('滤除噪声后的波形')

t=linspace(0,0.000001,32);

fc=1000000;

carrier=cos(2*pi*fc*t);

result=[];

for i=1:N

result=[result,carrier.*gussianhannel(((i-1)*32+1):i*32)]; end

figure(4);

subplot(2,1,1);

plot(result);

title('与相干波相乘后的波形')

result=[result];

result1=filter(lowpass,1,result);

subplot(2,1,2);

plot(result1);

title('通过低通滤波器滤波后的波形')

IN=fix(length(lowpass)/2);

bitstream=[];

i=IN+LL/2;

while (i<=length(result1))

bitstream=[bitstream,result1(i)>0];

i=i+LL;

end

axis([0,300,-50,50]);

end

//fuzhi.m

function uu=fuzhi %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Jakes fading simulator

% designed by Kevin Jungclear; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% Model the received complex low-pass envelope

N1=10; N2=100; % Number of the sinusoids D=1; % "D" denotes variance fmax=60; % Maxium Doppler frequency

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%Compute correlation of U,U1 and U2.

K=0;

rUU=zeros(1,0.2/0.001);

rU1U1=zeros(1,0.2/0.001);

rU2U2=zeros(1,0.2/0.001);

rU1U2=zeros(1,0.2/0.001);

for T=0:0.001:0.2,

K=K+1;

rUU(K)=2*(D^2)*besselj(0,2*pi*fmax*T);

rU1U1(K)=(D^2)*besselj(0,2*pi*fmax*T);

rU2U2(K)=rU1U1(K);

rU1U2(K)=rU1U1(K);

end

% N=10 ; compute discrete Doppler frequencies fi,n

f=zeros(1,N1);

for n=1:(N1-1),

f(n)=fmax*cos(n*pi/(2*N1-1));

f(N1)=fmax;

% compute the gains ci,n

% compute the gains c1,n

c1=zeros(1,N1);

for n=1:(N1-1),

c1(n)=(2*sqrt(D))*sin(pi*n/(N1-1))/sqrt(N1-0.5);

end

c1(N1)=sqrt(D)/sqrt(N1-0.5);

% compute the gains c2,n

c2=zeros(1,N1);

for n=1:(N1-1),

c2(n)=(2*sqrt(D))*cos(pi*n/(N1-1))/sqrt(N1-0.5);

end

c2(N1)=sqrt(D)/sqrt(N1-0.5);

% compute the complex envelpe

M=2000; % M denotes the time interval area

u1=zeros(1,M); % u1 denotes μ1(t)

u2=zeros(1,M); % u2 denotes μ2(t)

u=zeros(1,M); % u denotes μ(t)

k=0;

% compute μ1(t)

for t=0:0.01:20-0.01,

w1=cos(2*pi*f*t);

ut1=c1*w1';

k=k+1;

u1(k)=ut1;

end

% compute μ2(t)

k=0;

for t=0:0.01:20-0.01,

w2=cos(2*pi*f*t);

ut2=c2*w2';

k=k+1;

u2(k)=ut2;

end

% compute μ(t)

for n=1:M

u(n)=u1(n)+i*u2(n);

for n=3:202

uu=-abs(u);

end

end

end

//Gain.m

function[raleigh]=Gain(Fd,N)

%% Producing a raleigh distribution with given Doppler frequency Fd;

%% N>10000,small N will distroy the gaussion of WGN;

Fs=256000; %% channel sample rate;

sgma=1;

tmax=40/(2*pi*Fd); %% 20dB time point;

r=ceil(64/(2*tmax)); %% symbol rate of the gaussian data;

rn=ceil(r*N/Fs); %% the numble of gaussion noise;

L=Fs/r; % interpolation factor;

gnr=WGN(1,rn,0); % white gaussion noise;

gni=WGN(1,rn,0); % white gaussion noise;

gn=gnr+j*gni;

Jakes=Jfilter(Fd); % Jakes' filter with Doppler frequency Fd;

J=conv(gn,Jakes); % Doppler filter;

J=abs(J);

J1=J(33:(rn+32)); % deleting the head and the tail;

JN=sqrt(var(J1));

J1=J1/JN;

m=1:rn; %% interpolation with factor L;

n=1:1/L:rn;

raleigh=spline(m,J1,n); %% raleigh data with Fd and channel sample rate; size(raleigh)

//gussianhannel.m

function channelsignal=gussianhannel(transmittedSignal,sig)

%sig为加人的噪声db数

channelsignal=awgn(transmittedSignal,sig);

subplot(2,1,2);

plot(channelsignal);

title('加上高斯白噪声后的信号波形')

end

//Jfilter.m

function[hnw]=Jfilter(Fd)

% generate the time domain response according the Jakes Spectrum Model;

% fd is the maximum Doppler frequency;

NU=0.25;

Xt=20/(2*pi*Fd);

C=1.457; %% normolizing factor; t=-Xt:(2*Xt/63):Xt;

x=2*pi*Fd*abs(t);

ht=C*Fd*BESSELJ(NU,x)./sqrt(sqrt(x)); %% Bessel function;

n=0:63; %% 64 taps Jakes' filter; hn=ht;

hnw=hn.*hamming(64)'; %% hamming windowed; //onemain.m

clear all;

clc;

N=20;%产生20个序列

sig=100;

signal=bin(N);%返回二进制的字符串

figure(1);

stem(signal)

title('显示原序列')

%显示原序列

axis([0,10,-1,2]);

modemsignal=BpskModu(signal);

h=xitong;

channelsignal1=conv(modemsignal,h);

figure(3);

subplot(2,1,1);

plot(channelsignal1);

title('加路径损耗和阴影后的信号波形');

channelsignal2=gussianhannel(channelsignal1,sig);

bitstream=demobpsk(channelsignal2,N);

figure(5);

stem(bitstream)

title('解调后的序列')

axis([0,10,-1,2]);

%显示恢复后的序列

SNRindB3=0:0.2:10;

for i=1:length(SNRindB3),

SNR=exp(SNRindB3(i)*log(10)/10); % signal to noise ratio

snr=sqrt(SNR)

theo_err_prb(i)=0.5*erfc(snr); % theoretical bit error rat end;

figure(6);

semilogy(SNRindB3,theo_err_prb);

title('误码性能');

xlabel('SNR');

ylabel('BER');

//sunhaoshadow.m

function [pathlossv,pathlossvv]=sunhaoshadow

a=0.8;

do=200;

k=20;

distance=1000;

figure(7);

for distancex=do:10:distance;

pathlossv =10*log(k)- 10*a*log(distancex/do);

plot(distancex,pathlossv,'b*');

hold on;

R=NORMRND(0,2);

pathlossvv = pathlossv+R;

plot(distancex,pathlossvv,'r*');

hold on;

end;

legend('Pass Loss','Shadow');

title('大尺度路径损耗(对数正态阴影)与T-R距离的关系');

xlabel('T-R距离(m)');

ylabel('路径损耗(dB)')

%

%%

%% Plot the Raleigh channel taps;

fd = 200;

%%

%

% $$e^{\pi i} + 1 = 0$$

%

N = 200000; % N >200000,or it will distory the gaussian of the WGN;

h = TDL(fd,N);

h_dB = 20*log10(abs(h)+eps);

y = h_dB(64001:128000,:);

x=[1:length(y)]/256;

h2=zeros(64000,4);

h3=1./h2;

h3([1600:1600:64000],:)=y([1600:1600:64000],:);

%plot(x,h3(:,[1]),'*',x,h3(:,[2]),'v',x,h3(:,[3]),'o',x,h3(:,[4]),'d',x,y(:,[1]),x,y(:,[2]),x,y(:,[3]),x,y(:,[4])); %xlabel('Time(ms)');ylabel('Channel Amp(dB)');

%title('Raleigh Channel');

%legend('Tap 1;','Tap 2','Tap 3','Tap 4');

%grid on;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%

%pause;

figure(8);

%subplot(2,2,1);

plot(x,y(:,[1]));

grid on; xlabel('Time(ms)'); ylabel('Channel Amp(dB)');

title('Raleigh Channel'); legend('Tap 1');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%

//TDL.m

% genorate the impulse respose for the time variant Channel;

% 4 path

function[h]=TDL(Fd,N) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% sgma is the Tap Weight of each Delay channle %% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Fd = Fd;

sgma1=1;

sgma2=0.3679;

sgma3=0.1353;

sgma4=0.0498;

NormSgma=sqrt(sgma1^2+sgma2^2+sgma3^2+sgma4^2);

sgma1 = sgma1/NormSgma;

sgma2 = sgma2/NormSgma;

sgma3 = sgma3/NormSgma;

sgma4 = sgma4/NormSgma;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%

g=Gain(Fd,N)'; %% produce a raleigh distribution array

size(g)

h(:,[1])=sgma1*g; %% get the first tap Gain;

g=Gain(Fd,N)'; %% reproduce a raleigh distribution array

h(:,[2])=sgma2*g; %% get the 2nd tap Gain;

g=Gain(Fd,N)'; %% reproduce a raleigh distribution array

h(:,[3])=sgma3*g; %% get the 3thd tap Gain;

g=Gain(Fd,N)'; %% reproduce a raleigh distribution array h(:,[4])=sgma4*g; %% get the last tap Gain;

//xitong.m

function h=xitong

[L,S]=sunhaoshadow;

uu=fuzhi;

h=(L+S)*uu;

end