Demo entry 6741580

matlab

   

Submitted by tt on May 16, 2018 at 09:19
Language: Matlab. Code size: 4.6 kB.

fc = 1000;  %载波频率
fs = len*fc;%采样频率
A = 1;      %载波幅度
SNR = 10;   

num = 10e3; %码元个数
len = 20;   %码元长度

bt = [];    %基带信号
bt1 = [];  %恢复后的基带信号
n = 20;     %绘图周期个数
N = 2^18;   %FFT点数
ts=1/fs;    %采样周期
t0 = 10;     %总时长
t=ts:ts:t0; %时间序列
ct = A*sin(2*pi*fc*t); %载波

code = randi([0,1],1,num); %等概分布的基带信号
decode = zeros(1,num);     %判决后的基带信号

%生成基带信号
for i = 1:1:num   
    if (code(i) ==0)
        cp = zeros(1,len);
    else 
        cp = ones(1,len);
    end
    bt = horzcat(bt,cp);  %向量拼接
end
bt2 = 2*bt-1; %双极性基带信号
st = ct.*bt2; %BPSK调制后的信号
figure(1);
subplot(2,1,1);
plot(t(1:n*len),bt(1:n*len));       %基带信号
hold on 
grid
axis([0,n*len/fs,-0.1,1.1])
xlabel('t');
ylabel('b(t))'); 
title('基带信号');

subplot(2,1,2);
plot(t(1:n*len),st(1:n*len));       %BPSK调制信号
axis([0,n*len/fs,-1,1])
hold on 
grid
xlabel('t');
ylabel('s(t)'); 
title('BPSK调制信号');

%BPSK经过AWGN后相干解调
st_channel = add_noise(st,SNR); %加性高斯白噪声
st_channel = filterdesign(ceil(2*fc*N/fs),N,len*num,st_channel); %接收信号经过带通滤波器
dmst1 = st_channel.*ct; %乘以相干载波
dmst = filterdesign(ceil(fc*N/fs),N,len*num,dmst1); %经过低通滤波器
figure(2);
subplot(3,1,1);       %BPSK调制信号
plot(t(1:n*len),st(1:n*len));       
axis([0,n*len/fs,-1,1])
hold on 
grid
xlabel('t');
ylabel('s(t)'); 
title('BPSK调制信号');
subplot(3,1,2);       %BPSK调制信号经过AWGN信道
plot(t(1:n*len),st_channel(1:n*len));       
axis([0,n*len/fs,-1,1])
hold on 
grid
xlabel('t');
ylabel('s(t)'); 
title('通过AWGN信道后的BPSK信号(SNR=10dB)');
subplot(3,1,3);      %解调后的信号
plot(t(1:n*len),dmst(1:n*len));       
axis([0,n*len/fs,-1,1])
hold on 
grid
xlabel('t');
ylabel('s(t) '); 
title('BPSK解调信号');

%抽样判决
for k=1:num
    if(dmst(-len/2+len*k)>0) %取周期中点进行抽样判决
        decode(k) = 1;
    else 
        decode(k) = 0;
    end
end
for j=1:num                  %恢复基带信号成形
    if(decode(j) == 0)       
        cp2 = zeros(1,len) ;
    else 
        cp2 = ones(1,len);   
    end
    bt1 = horzcat(bt1,cp2);   
end
figure(3);
subplot(2,1,1);   
plot(t(1:n*len),bt(1:n*len));       
hold on 
grid
axis([0,n*len/fs,-0.1,1.1])
xlabel('t');
ylabel('b(t))'); 
title('基带信号');
subplot(2,1,2);
plot(t(1:n*len),bt1(1:n*len));      
hold on 
grid
axis([0,n*len/fs,-0.1,1.1])
xlabel('t');
ylabel('bt1(t))'); 
title('复原的基带信号');


%误码率与信噪比的关系曲线
figure(4);
spot = 60; %点数
range = 30; %信噪比范围
error_rate = zeros(1,spot);
SNR_array  = zeros(1,spot);
for m = 1:1:spot  %SNR步进
    SNR_array(m) = -30+range/spot*m;
    error_rate(m) = error_SNR(SNR_array(m),code,st); 
    hold on
end
plot(SNR_array,error_rate,'b',SNR_array,qfunc(2*sqrt(10.^(SNR_array/10))),'r');  %在误码率与SNR关系曲线上描点
legend('实际曲线','理论曲线');
hold on
set(gca,'YScale','log')           %纵坐标对数坐标
grid
xlabel('信噪比(dB)');
ylabel('误码率');
title('误码率与信噪比之间的的关系曲线');


%各种函数
%1  AWGN加性高斯白噪声
function [ Y,noise ] = add_noise( X,SNR ) 
noise = randn(size(X));
noise = noise - mean(noise);
signal_power = 1/length(X)*sum(X.^2);
noise_variance = signal_power/(10^(SNR/10));
noise = sqrt(noise_variance)/std(noise)*noise;
Y = X + noise;
end

%2误码率计算函数
function [ error_rate] = error_counter( code,decode,num )
error = abs(code-decode); %误码序列
error_num = 0;
for j = 1:1:num
    if error(j)~=0
        error_num = error_num+1; %对误码计数
    end
end
 error_rate = error_num/num; %计算误码率
end

%3  确定SNR下的误码率计算函数
function [error_rate] = error_SNR( SNR , code , st )

num = 10e3;%码元个数
len = 20; %码元长度
fc = 1000; %载波频率(Hz)
fs = len*fc;%采样频率(Hz)
A = 1;     %载波幅度
ts=1/fs;   %采样的时间步长
t0 = 10; %仿真时长(秒)
t=ts:ts:t0;%时间轴序列
ct = A*sin(2*pi*fc*t); %载波
N = 2^18;   %FFT分析点数
decode = zeros(1,num); %判决后的基带信号序列

st_channel = add_noise(st,SNR); %加性高斯白噪声
passband = ceil(2*fc*N/fs);     %数字滤波器的等效截止频率
st_channel = filterdesign(passband,N,len*num,st_channel); %接收信号先经过带通滤波器滤除带外噪声
dmst1 = st_channel.*ct; %乘以本地的相干载波
dmst = filterdesign(passband,N,len*num,dmst1); %低通滤波器

for k=1:num
    if(dmst(-len/2+len*k)>0) %取码元周期的中点进行抽样判决
        decode(k) = 1;
    else 
        decode(k) = 0;
    end
end
error_rate = error_counter(code,decode,num); %计算该信噪比下的误码率
end

%4  滤波器设计
function [ output ] = filterdesign(N_filter,N_FFT,data_length,input)
filter1 = ones(1,N_filter);
filter2 = zeros(1,N_FFT/2 - N_filter);
filter0 = [filter1,filter2,filter2,filter1];
input_spectrum = fft(input,N_FFT);%输入频谱
output_spectrum = input_spectrum.*filter0;%输出频谱
output_spectrum_ifft = real(ifft(output_spectrum,N_FFT));%ifft
output = real(output_spectrum_ifft(1:round(data_length)));
end

This snippet took 0.02 seconds to highlight.

Back to the Entry List or Home.

Delete this entry (admin only).