通信系统的计算机模拟第十三讲.ppt

1,通信系统的计算机模拟,第十三讲,2,通信系统的蒙特卡罗仿真,相移键控（Phase-shift Keying, PSK）数字通信系统，尽管非常简单，基本构建模块，差分QPSK系统的仿真，其中考虑了相位和符号同步误差的影响。 讨论半解析（SA）方法，即把蒙特卡罗仿真和解析的方法结合在一起。 蒙特卡罗方法几乎不需要任何的数学分析，只要系统方框图中每个功能模块的信号处理算法是已知的就能应用 蒙特卡罗仿真是一个很通用的工具，缺点是仿真运行时间长，半解析方法需要更高水平的分析，降低了仿真运行时间。 运行一次蒙特卡罗仿真得到的是单个误比特率估计值，而半解析仿真得到的是以为Eb/N0函数的完整的误比特率曲线，半解析仿真不是一个普遍适用的方法，仅适用于一类很有限的系统，半解析仿真消耗的机时是微不足道的，因此在能运用半解析仿真时，都应优先考虑这种方法。,3,两个蒙特卡罗实例,将蒙特卡罗方法用于估计数字通信系统的误比特率时，是让N个采样符号通过系统的仿真模型，并计算产生差错的个数来实现的 假设符号导致了Ne个差错，那么误比特率的估计值为 是一个随机变量，要获得误比特率的准确估计，估计器必须是无偏，并具有小的方差，小的方差要求有大的N，而这又会导致较长的计算时间。,4,例10-1（PSK）,假设为BPSK调制，信号星座中的两个信号点都位于直接（同相）信道上（回顾一下例9-3），根据这一假设，在仿真中我们可以不考虑正交信道。 假设调制器输出端的滤波器是三阶巴特沃思滤波器，其带宽等于比特率（BW=rb），该滤波器会产生码间干扰（ISI）。 仿真的目的是确定由滤波器带来的ISI所增加的误比特率。 这里采用了块级联（block-serial）的方法，迭代处理由1000个符号组成的块，直到处理完所有N个符号。 这样做主要是为了利用MATLAB内置的函数filter，它可以实现时域卷积。它作为内置函数具有非常高的运行效率，可大大地减少仿真时间。注意必须保证滤波器的输出在块与块之间是连续的，这是通过使用filter提供的初始条件参数来实现的。,5,Solution,第一个问题是要确定delay的值 最妙的方法是，把调制器的输入和接收机的输出进行互相关。 在进行半解析仿真时我们会采用这种方法。 为了说明正确选择时延值的重要性，具体来说，我们会选定一个值，用不同的delay值对系统进行仿真，并观察结果,6,Matlab Program,EbNodB = 6; % Eb/No (dB) value z = 10.(EbNodB/10); % convert to linear scale delay = 0:8; % delay vector BER = zeros(1,length(delay); % initialize BER vector Errors = zeros(1,length(delay); % initialize Errors vector BER_T = q(sqrt(2*z)*ones(1,length(delay); % theoretical BER vector N = round(100./BER_T); % 100 errors for ideal (zero ISI) system FilterSwitch = 1; % set filter switch (in=1 or out=0) for k=1:length(delay) BER(k),Errors(k) = c10_MCBPSKrun(N(k),z,delay(k),FilterSwitch) end semilogy(delay,BER,o,delay,BER_T,-); grid; xlabel(Delay); ylabel(Bit Error Rate);,7,Result,AWGN环境下的理想系统（没有ISI）在Eb/N0=6dB时的性能,8,function BER,Errors=MCBPSKrun(N,EbNo,delay,FilterSwitch) SamplesPerSymbol = 10; % samples per symbol BlockSize = 1000; % block size NoiseSigma = sqrt(SamplesPerSymbol/(2*EbNo); % scale noise level DetectedSymbols = zeros(1,BlockSize); % initialize vector NumberOfBlocks = floor(N/BlockSize); % number of blocks processed BTx,ATx = butter(5,2/SamplesPerSymbol); % compute filter parameters TxOutput,TxFilterState = filter(BTx,ATx,0); % initialize state vector BRx = ones(1,SamplesPerSymbol); ARx=1; % matched filter parameters Errors = 0; % initialize error counter % % Simulation loop begine here. % for Block=1:NumberOfBlocks,9,% % Generate transmitted symbols. % SymbolSamples,TxSymbols = random_binary(BlockSize,SamplesPerSymbol); % % Transmitter filter if desired. % if FilterSwitch=0 TxOutput = SymbolSamples; else TxOutput,TxFilterState = filter(BTx,ATx,SymbolSamples,TxFilterState); end % % Generate channel noise. % NoiseSamples = NoiseSigma*randn(size(TxOutput); % % Add signal and noise. %,10,RxInput = TxOutput + NoiseSamples; % % Pass Received signal through matched filter. % IntegratorOutput = filter(BRx,ARx,RxInput); % % Sample matched filter output every SamplesPerSymbol samples, % compare to transmitted bit, and count errors. % for k=1:BlockSize, m = k*SamplesPerSymbol+delay; if (m length(IntegratorOutput) DetectedSymbols(k) = (1-sign(IntegratorOutput(m)/2; if (DetectedSymbols(k) = TxSymbols(k) Errors = Errors + 1; end end end end BER = Errors/(BlockSize*NumberOfBlocks); % calculate BER,11,function x, bits = random_binary(nbits,nsamples) % This function genrates a random binary waveform of length nbits % sampled at a rate of nsamples/bit. x = zeros(1,nbits*nsamples); bits = round(rand(1,nbits); for m=1:nbits for n=1:nsamples index = (m-1)*nsamples + n; x(1,index) = (-1)bits(m); end end % End of function file.,12,讨论,正确的时延值极有可能在5和6个采样周期之间。要更精确地确定delay的估计值，可以用更高的采样频率（更小的采样周期）再次运行仿真程序。 既然已经知道了合适的delay值，现在就可以运行仿真程序，并确定以Eb/N0为函数的Pe值。,13,EbNodB = 0:8; % vector of Eb/No (dB) values z = 10.(EbNodB/10); % convert to linear scale delay = 5; % enter delay value (samples) BER = zeros(1,length(z); % initialize BER vector Errors = zeros(1,length(z); % initialize Errors vector BER_T = q(sqrt(2*z); % theoretical (AWGN) BER vector N = round(20./BER_T); % 20 errors for ideal (zero ISI) system FilterSwitch = 1; % Tx filter out (0) or in (1) for k=1:length(z) N(k) = max(1000,N(k); % ensure at least one block processed BER(k),Errors(k) = c10_MCBPSKrun(N(k),z(k),delay,FilterSwitch) end semilogy(EbNodB,BER,o,EbNodB,BER_T) xlabel(E_b/N_0 - dB); ylabel(Bit Error Rate); grid legend(System Under Study,AWGN Reference,0),14,结果,15,讨论,16,例10-2,（QPSK）在前面关于BPSK调制的例子中，为了使仿真代码简洁也为了易于分析，做了一些简化的假设。 而对此例的OPSK系统，我们对一些新的误差源进行了建模，并且包含了一些新的参数。以便比较容易地将仿真结果和由实际通信系统获得的结果关联起来。 例如，描述发射机和接收机之间传输损耗的信道衰减包括在仿真参数中 符号速率的实际（未经缩放的）数值和采样频率也包括在仿真中，该系统的方框图如图10-4所示,17,例10.2,18,载波和符号同步的仿真,对于相干射频（RF）系统，接收机必须具备提供载波和符号同步的功能 但噪声和信道失真使接收机不可能完美地实现载波和符号同步 不正确的载波同步将会导致发送信号相对于接收信号产生相位误差或相位旋转。 把相位误差仿真成随机过程，符号同步误差会导致积分-清除，检测器在不正确的时间区间上处理接收信号。 全例的仿真也使用户能够检验这种影响所造成的误差。,19,相位模糊问题,QPSK系统存在相位模糊的问题 由于信道造成了未知的信号时延，所以接收机不可能确定发送信号的绝对相位 假设信道所造成的时间时延为100.75个RF载波周期，接收机将会把45到-45 其他还会犯类似的错误。 不把信息位的编码包含在绝对相位中，而是包含在符号间的相位差中。 相位从45到135增加了90，接收机所检测到信号为-135到45，仍然保持相位增加90。 差分-附录B,20,Delay 估计,Eb = 23; No = -50; % Eb (dBm) and No (dBm/Hz) ChannelAttenuation = 70; % channel attenuation in dB N = 1000; delay = -0.1:0.1:0.5; EbNo = 10.(Eb-ChannelAttenuation)-No)/10); BER_MC = zeros(size(delay); for k=1:length(delay) BER_MC(k) = c10_MCQPSKrun(N,Eb,. No,ChannelAttenuation,delay(k),0,0,0); disp(Simulation ,. num2str(k*100/length(delay),% Complete); end BER_T = 0.5*erfc(sqrt(EbNo)*ones(size(delay); % Theoretical BER semilogy(delay,BER_MC,o,delay,2*BER_T,-) %Plot BER vs Delay xlabel(Delay (symbols); ylabel(Bit Error Rate); legend(MC BER Estimate,Theoretical BER);grid;,21,结果,因为没有使用信道滤波器，最优时延是零符号周期，周期来测量时延，而不是像前面的例子那样用采样周期。,22,误比特率对静态同步相位误差的灵敏度 当相位误差以10度为增量从0度变化到90度,PhaseError = 0:10:90; % Phase Error at Receiver Eb = 24; No = -50; % Eb (dBm) and No (dBm/Hz) ChannelAttenuation = 70; % dB EbNo = 10.(Eb-ChannelAttenuation-No)/10); BER_T = 0.5*erfc(sqrt(EbNo)*ones(size(PhaseError); N = round(100./BER_T); BER_MC = zeros(size(PhaseError); for k=1:length(PhaseError) BER_MC(k) = c10_MCQPSKrun(N(k),Eb,No,ChannelAttenuation,0,0,. PhaseError(k),0); disp(Simulation ,num2str(k*100/length(PhaseError),% Complete); end semilogy(PhaseError,BER_MC,o,PhaseError,2*BER_T,-) xlabel(Phase Error (Degrees); ylabel(Bit Error Rate); legend(MC BER Estimate,Theoretical BER); grid;,23,结果,误比特率在相位误差为45的时候达到最大，然后在相位误差为0或90时下降到最优值（零同步相位误差时的值），这是差分编码所带来的结果,24,BER:最优相移和信道时延已知,Eb = 22:0.5:26; No = -50; % Eb (dBm) and No (dBm/Hz) ChannelAttenuation = 70; % Channel attenuation in dB EbNodB = (Eb-ChannelAttenuation)-No; % Eb/No in dB EbNo = 10.(EbNodB./10); % Eb/No in linear units BER_T = 0.5*erfc(sqrt(EbNo); % BER (theoretical) N = round(100./BER_T); % Symbols to transmit BER_MC = zeros(size(Eb); % Initialize BER vector for k=1:length(Eb) % Main Loop BER_MC(k) = c10_MCQPSKrun(N(k),Eb(k),No,ChannelAttenuation,0,0,0,0); disp(Simulation ,num2str(k*100/length(Eb),% Complete); end semilogy(EbNodB,BER_MC,o,EbNodB,2*BER_T,-) xlabel(Eb/No (dB); ylabel(Bit Error Rate); legend(MC BER Estimate,Theoretical BER); grid;,25,仿真结果,26,相位抖动对系统误比特率的影响 用白噪声来对相位误差过程进行建模,PhaseBias = 0; PhaseJitter = 0:2:30; Eb = 24; No = -50; % Eb (dBm) and No (dBm/Hz) ChannelAttenuation = 70; % dB EbNo = 10.(Eb-ChannelAttenuation-No)/10); BER_T = 0.5*erfc(sqrt(EbNo)*ones(size(PhaseJitter); N = round(100./BER_T); BER_MC = zeros(size(PhaseJitter); for k=1:length(PhaseJitter) BER_MC(k) = c10_MCQPSKrun(N(k),Eb,No,ChannelAttenuation,0,0,. PhaseBias,PhaseJitter(k); disp(Simulation ,num2str(k*100/length(PhaseJitter),% Complete); end semilogy(PhaseJitter,BER_MC,o,PhaseJitter,2*BER_T,-) xlabel(Phase Error Std. Dev. (Degrees); ylabel(Bit Error Rate); legend(MC BER Estimate,Theoretical BER); grid;,27,讨论,当相位抖动的标准差增加时，误比特率也会增加 在许多系统的仿真中，用白噪声来建模相位抖动是不恰当的 如果是这样，可以设计一个有限冲激响应滤波器（FIR）来实现相位抖动过程所需要的功率谱密度（PSD）。,28,结果,29,误比特率对符号同步的灵敏度,SymJitter = 0:0.02:0.2; Eb = 24; No = -50; % Eb (dBm) and No (dBm/Hz) ChannelAttenuation = 70; % channel attenuation in dB EbNo = 10.(Eb-ChannelAttenuation-No)/10); BER_T = 0.5*erfc(sqrt(EbNo)*ones(size(SymJitter); N=round(100./BER_T); BER_MC = zeros(size(SymJitter); for k=1:length(SymJitter) BER_MC(k) = c10_MCQPSKrun(N(k),Eb,No,ChannelAttenuation,0,SymJitter(k),0,0); disp(Simulation ,num2str(k*100/length(SymJitter),% Complete); end semilogy(SymJitter,BER_MC,o,SymJitter,2*BER_T,-) xlabel(Symbol Timing Error Std. Dev. (Symbols); ylabel(Bit Error Rate); legend(MC BER Estimate,Theoretical BER); grid;,30,建模精度,如果要对符号抖动过程的记忆效应进行精确建模FIR 符号是互相关的vxcorr，以便正确地确定误比特率。要采用互相关方法来计算合适的时延值。,31,结果,32,10.2 半解析方法,蒙特卡罗方法是完全通用的，除了实现子系统所要求的外，不需要用到任何分析知识。 运行时间？ 半解析方法就是其中最有效的方法之一。,33,半解析的条件,34,10.2.1 基本考虑,35,差错概率,36,仿真方法,37,仿真V1,38,仿真方差,39,差错概率,40,非线性影响,41,42,10.2.2 等效噪声源,在应用半解析方法时，我们采用了等效噪声源的思想 我们已经看到判决统计是三个分量的函数，也就是说 其中是 Sk 由于信号产生的 Dk 来自于由系统因素引起的失真（如ISI） 而 Nk 是由噪声引起的。 Sk和Dk 的影响可由蒙特卡罗仿真确定，而如前所见，由Nk表示的噪声的影响则用解析的方法处理。 如果执行无噪声的仿真，则所得的充分统计量只是Sk和Dk的函数。该统计量加上一个方差由式（10-11）给出的随机变量Nk，可得 可将随机变量Nk视作来自于如图10-13所示等价噪声源ne(t)的一个样本。,43,cont,该等价噪声源包含了反映在积分-清除检测器的积分器输出中的 由热噪声、干扰和其他信道损伤所产生的总影响。 如果信道噪声是白噪声，则可用冲激响应，或等价地用式（10-11）所定义的传递函数，来把信道噪声变换到积分器的输出中。,44,10.2.3 PSK系统误比特率的半解析估计方法,半解析仿真来确定BPSK系统中的误比特率 我们采用一种很容易就可护展到QPSK系统的方法进行处理,而接收信号落在区域,考虑如图10-14所示的信号星座，发送信号的点用 表示，相应的判决区域用 表示。如果发送的是 而接收信号落在区域 中，则接收机作出了正确的判决；不然，就发生了差错。,45,差错概率,46,cont,47,48,BPSK 附录C,49,10.2.4 QPSK系统误比特率的半解析估计方法,由于QPSK的信号星座图有四个点而不是两个 并且信号空间是二维而不是一维 因此QPSK半解析估计器和PSK估计器的不同之处在于前者必须为正交信道增加一维。,发送的信号点记作S1，S2，S3，S4，判决区域记作D1，D2，D3，D4，。 如果发送的是Si，而接收到的信号落在区域Di中，则接收机作出了正确的判决，否则就发生了差错。,50,假设发送的是 ，接收的无噪声信号记作 。由于存在的符号间干扰和失真， 。 因为仿真考虑了符号间干扰的影响而没有考虑噪声的影响，所以半解析仿真所确定的是 。的同相和正交分量分别用 表示 当考虑噪声时，加入 nx, ny 。 在发送 S1 的条件下，如果（ ），则作出了正确的判决； 否则，则发生了差错。 由于开发的是半解析仿真估计器，噪声的影响用解析方法作了处理，所以没有出现在图10-16中,51,nx,ny及分布,假设同相和正交噪声分量不相关并且是联合高斯分布的,52,上界,可以看出式（10-26）中的四个积分中有两个为1,53,总的上界,符号k,54,例10-4（QPSK）,QPSK系统半解析仿真的MATLAB执行代码在附录D中。 仿真的目的是考察发送滤波器引起的ISI所产生的影响。 滤波器的带宽设为符号率（比特率的一半） 因为信号星座图是对称的，如同前面所讨论的一样，所有接收到的信号都旋转到第一象限。