第三章随机过程
本节首先介绍利用matlab现有的库函数根据实际需要直接产生均分分布和高斯分布随机变量的方法，然后重点讲解蒙特卡罗算法。
一、均匀分布的随机数
利用MATLAB库函数rand产生。
rand函数产生（0，1）内均匀分布的随机数，使用方法如下：
1）x=rand(m)；产生一个m×m的矩阵，所含元素取值均为在（0，1）内均匀分布的随机数。
2）x=rand(m，n)；产生一个m×n的矩阵，所含元素取值均为在（0，1）内均匀分布的随机数。
3）x=rand；产生一个随机数。
举例：1、产生一个5×5服从均匀分布的随机矩阵，所含元素取值均为在（0，1）内均匀分布的随机数。
x=rand(5)
2、产生一个5×3服从均匀分布的随机矩阵，所含元素取值均为在（0，1）内均匀分布的随机数。
x=rand(5,3)
二、高斯分布的随机数
randn函数产生均值为0，方差为1的高斯分布的随机数，使用方法如下：
1）x=randn(m)；产生一个m×m的矩阵，所含元素都是均值
为0，方差为1的高斯分布的随机数。
2）x=randn(m，n)；产生一个m×n的矩阵，所含元素都是均值为0，方差为1的高斯分布的随机数。
3）x=randn；产生一个均值为0，方差为1的高斯分布的随机数。
举例：1、产生一个5×5的矩阵，所含元素都是均值为0，方差为1的高斯分布的随机数。
x=randn(5)
2、产生一个5×3的矩阵，所含元素都是均值为0，方差为1的高斯分布的随机数。
x=randn(5,3)
3、产生一个5×3的矩阵，所含元素都是均值为0，方差为4的高斯分布的随机数。
x=2×randn(5,3)
三、蒙特卡罗仿真
1、蒙特卡罗算法
蒙特卡罗估计是指通过随机实验估计系统参数值的过程。
蒙特卡罗算法的基本思想：由概率论可知，随机实验中实验的结果是无法预测的，只能用统计的方法来描述。
故需进行大量的随机实验，如果实验次数为N，以
N表示事件A发
A
生的次数。
若将A发生的概率近似为相对频率，定义为
N N。
A
这样，在相对频率的意义下，事件A发生的概率可以通过重
复无限多次随机实验来求得，即：
()lim
A N N P A N
→∞
=
在二进制数字通信系统中，若N 是发送端发送的总码元数，A N 是差错发生的次数，则总误码率可通过蒙特卡罗算法计算。
2、举例
本节用蒙特拉罗仿真研究一个简单的二进制双极性数字基带通信系统的误比特率。
数字基带信号传输系统模型如图1所示：
图1 数字基带信号传输系统模型
假设该通信系统满足以下条件：
① 信源输出的数据符号是相互独立和等概的双极性基带信号
② 发送端没有发送滤波器，接收端没有接收滤波器，满足无码间串扰条件
③ 信道是加性高斯白噪声信道，即只考虑噪声对误比特率的影响。
⑴理论分析：
由通信理论可知，对于二进制双极性数字基带通信系统，当1，0出现概率相同，即
()()1012P P ==
时，最佳判决门限： 0d V *=
误码率：1
122e P erfc erfc ⎛⎫==
，利用2/2
()t x
Q x e
dt ∞
-=
⎰
，2
()t x
erfc x e dt ∞
-=
⎰
，可得
1()2Q x erfc =，
故可用Q
函数表示误码率e P Q =。
此时，在接收端，抽样判决器输入信噪比：22
n r E σ=
抽样判决器输入信号为：111000r E n r E n =+⎧⎨=-+⎩
，发，发
E ±为判决器输入有用信号电压，1n ，0n 为信道输入的均值为0，方差为2
n σ高斯噪声。
依据上述分析，可得通信系统的蒙特卡罗仿真模型如图2所示。
图2 通信系统的蒙特卡罗仿真模型
⑵仿真流程：
①规定信号电压E=1
②将信噪比从dB数转化为信号与噪声的功率比
③计算噪声方差
④用均匀分布的随机数产生二进制数字信号，若随机数大于或等于0.5，则产生1，用高电平表示；否则产生0，用低电平表示
⑤将所产生的数字信号送入信道，叠加高斯白噪声（均值为0，方差由③产生）
⑥在接收端，对收到的信号按最佳判决门限进行判决
⑦比较原始数字信号和判决后的数字信号
⑧计算不一样的码元的个数，得到误比特率
⑨设定不同的信噪比，重复②～⑧，得到不同信噪比下的误码率，画出曲线，并和这些信噪比下的理论误码率相比较。
⑶源代码：
function [p]=smldPe54(snr_in_dB)%计算误码率
%信噪比与误码率的互换
E=1;
SNR=10^snr_in_dB/10;
sgma=sqrt(1/SNR);
%二进制序列的产生
N=10000;
for i=1:N
temp=rand;
if(temp<0.5)
dsource(i)=0;
else
dsource(i)=1;
end
end;
%计算误码率
numoferr=0;
for i=1:N
if(dsource(i)==0)
r=-E+gngauss(sgma); else
r=E+gngauss(sgma); end
if(r<0)
decis=0;
else
decis=1;
end
if(decis~=dsource(i))
numoferr=numoferr+1;
end
end
p=numoferr/N;
高斯随机数发生器
function [gsrv1,gsrv2]=gngauss(m,sgma) if nargin==0
m=0;
sgma=1;
elseif nargin==1
sgma=m;
m=0;
end
u=rand;
z=sgma*(sqrt(2*log(1/(1-u))));
gsrv1=m+z*cos(2*pi*u);
gsrv2=m+z*sin(2*pi*u);
理论误码率计算
function [y]=Qfunct(x)
y=(1/2)*erfc(x/sqrt(2));
主程序：
echo on;
SNRindB1=0:1:10;
SNRindB2=0:0.1:10;
%计算实际误码率
for i=1:length(SNRindB1)
smld_err_prb(i)=smldPe54(SNRindB1(i)); end
%计算理论误码率
for i=1:length(SNRindB2)
SNR=exp(SNRindB2(i)*log(10)/10);
theo_err_prb(i)= Qfunct(sqrt(SNR));
end
semilogy(SNRindB1,smld_err_prb,'r*');
hold
semilogy(SNRindB2,theo_err_prb);
仿真结果：
012345678910
10
10
10
10
10
解法：2 clear all close all clc
EbN0dB=1:0.5:10 N0=10.^(-EbN0dB/10); sigma=sqrt(N0/2); %理论计算的误码率 Pb=0.5*erfc(sqrt(1./N0)); %仿真误码率
numberror=zeros(1,length(EbN0dB)) for n=1:length(EbN0dB)
a=sign(rand(1,100000));%产生等概信源+1，-1
rk=a+sigma(n)*randn(1,100000);%离散等效接收模型dec_a=sign(rk);%判决
ber(n)=sum(abs(a-dec_a)/2)/length(a);
end
semilogy(EbN0dB,Pb);
hold;
semilogy(EbN0dB,ber,'rd-');
legend('理论值','仿真结果');
xlabel('Eb/N0(dB)');
Ylabel('Pb');。