clear all;
close all;
clc;
%%
Ldata = 2^10; %data部分总符号长度 800ms
LUW = 2^8; %UW部分长度 水声信道最大时延扩展不超过50ms
NumFrame = 100; %帧数一般情况下仿10000帧,时变信道下仿5000帧
SNR = [0:0.5:20]; %信噪比dB
index = 1; %进制数转换
% index: modulation index
% 1---bpsk
% 2---qpsk
% 4---16qam
% 6---64qam
% else is error
%载波频率fc = 10kHz
%带宽是B = 2kHz
%采样速率fs = 10kHz
%采样间隔ts = 0.1ms
%速度v = 10km/h
%最大多普勒频移fd = fc*v/1500 = 18.7Hz
%假设船只不动,水流运动造成的fd = 1.5Hz
%一个生成多径信道的函数,信道的每一径的衰落都独立的服从Rayleigh分布chan=rayleighchan(Ts,fd,tau,pdb)
%Ts:采样时间,如果考虑基带信号,这个和接收机要处理的数据速率是一样的,要考虑过采样的影响
%fd:就是Doppler频偏,以Hz为单位,与速率的换算关系为v×fc/c,fc是载频
%tau:输入的信道参数,一个向量,包含了各径的延时,以s为单位
%pdb:输入的信道参数,一个向量,包含了各径的功率(当然是均值啦,实际产生的能量都是以此为均值的随机量),以dB为单位
fs = 10*10^3; %采样速率10KHz
ts = 1/fs; %抽样间隔0.1ms
fc = 5*10^3; %载波频率5kHz
fd = 0.2; %多普勒频移20Hz
tau = [0 0.13 0.13 0.5 0.5 1.2 1.2 2.1 2.1 3.3 3.3 4.8 4.8 6.5 6.5 8.5 8.5 10.8 10.8 13.3]*10^-4; %20径延时
pdb = [0 0 -0.967 -0.967 -0.967 -0.967 -1.933 -1.933 -1.933 -1.933 -2.900 -2.900 -2.900 -2.900 -3.86 -3.86....
-3.86 -3.86 -4.84 -4.84]; %20径衰落dB
% tau = [0,0.2]*10^-4;
% pdb = [0,-0.38];
%%
%信道模型
% pedAchannel = [1 10^(-9.7/20) 10^(-22.8/20)];
% pedAchannel = pedAchannel/sqrt(sum(pedAchannel.^2));
% vehAchannel = [1 0 10^(-1/20) 0 10^(-9/20) 10^(-10/20) 0 0 0 10^(-15/20) 0 0 0 10^(-20/20)];
% vehAchannel = vehAchannel/sqrt(sum(vehAchannel.^2));
% idenChannel = 1;
RayleighChan = rayleighchan(ts,fd,tau,pdb);
%选择信道
%InputSamplePeriod:直接设置即可,无线信道传输的基带信号的符号周期。
%PathDelays:每条径的延时,单位是seconds。
%MaxDopplerShift:最大多普勒频移:对于rayleighchan函数的使用方法研究
%AvgPathGaindB:每条路径增益的平均功率,单位是分贝。
%1) ChannelType:这个就没什么好说了,只是说明使用函数生成的信道是瑞利的还是莱斯的。
%2) PathGains:信道各路径的真实增益,这个增益不是通过AvgPathGaindB换算过来的,AvgPathGaindB是平均增益,所以生成的各路径增益的功率值是围绕这个值生成的。
%3) ChannelFilterDelay:这个值没有什么好说的,原因是这个值没有什么用,它一般是和信道的PathDelays的第一个值是对应的。
%4) NumSamplesProcessed:这个值是说信道处理了多少个数据,如果输入信号矢量的长度是N,则这个值就是N。}
%5) StorePathGains:将chan2.StorePathGains设置为1时,各路径的增益将会被记录下来,但是不能使用Plot(chan2)画图。
%在这里额外说一点,信道的StorePathGains设为1时,查看chan2.PathGains可以发现这个值是一个N×M的矩阵,在这里N是信道处理的数据长(也就是输入信号的长度),M是多径数。这里的chan2.PathGains每列是线性变化的,原因是在处理每个数据时信道由于受到多普勒谱的影响是时变的,但信道的状态是连续的。
% Channel = pedAchannel;
% Channel = vehAchannel;
% Channel = idenChannel;
Channel = RayleighChan;
Channel.NormalizePathGains = false;
%6)将Channel.ResetBeforeFiltering在使用filter函数前设置为0时,可以看到信道各路径的斜率会是一个值,也就是无论处理多少数据信道的实际状态都是连续的,无论其他的条件(比如噪声和输入信号)发生怎样的变化,
%信道的状态是连续的,但不是说信道的增益值是恒定的数值,只是信道增益的实际值线性增长的,信道的状态时连续的。(每次使用filter函数一次,信道增益的实际数值都是会发生变化的。只要信道的状态不变
%就相当于信道没有发生质的变化,在仿真时我们可以认为使用的是同一个信道。如果设置为1,则是说明每次生成的信道增益的状态是不连续的,也就是信道发生了变化,我认为这在仿真时可以控制信道变化的快慢。
Channel.ResetBeforeFiltering = false;
Channel.StoreHistory = 1; %将Channel.StoreHistory设置为1时,可以把信道的信息记录下来,并可以通过Plot(Channel.StoreHistory)画出当前的信道的时域IR,频域相应等各种图。
Channel.NormalizePathGains = 1; %将Channel.NormalizePathGains设置为1时,每条路径增益绝对值的平方和(即功率和)的平均值为1。这个可以通过多次使用rayleighchan生成信道,求取信道各路径增益绝对值的平方和即可看到其平均值为1。
%所选择信道的FFT变换
% H_channel0 = fft(Channel,Ldata);
H_channel0 = fft(Channel.PathGains./sqrt(sum((abs(Channel.PathGains)).^2)),Ldata);
%%
%生成UW序列(CHU序列)
UW = zeros(1,LUW);
u = LUW;
for k = 0:u-1
Q1(k+1) = pi*k^2./u;
end
I = cos(Q1);
Q = sin(Q1);
UW = I+i*Q;
% 生成UW序列(FRANK序列)
% UW=zeros(1,LUW);
% u = LUW;
% for(p = 0:sqrt(u)-1)
% for(q = 0:sqrt(u)-1)
% Q(p+q*sqrt(u)+1) = 2*pi*p*q/sqrt(u);
% end
% end
% I = cos(Q);
% Q = sin(Q);
% uw=I;
% UW = Q+i*I;
%%
%不同信噪比下发送NumFrame帧数据统计相应的误码率
for n = 1:length(SNR) %对每一snr计算一次
tic %SC-FDE计时开始
%误码数量
Errcount1_ZF = 0;
Errcount1_MMSE = 0;
Errcount0 = 0;
for k = 1:NumFrame
%随机生成data
Inputsymbols = randint(1,(Ldata-LUW)*index);
%调制
Data = modulation(Inputsymbols,index);
% if SNR(n)==20 && k==NumFrame
% scatterplot(Data);
% title('调制星座图');
% grid on;
% end;
%帧合成
% Txsymbols = [UW,Data]
Txsymbols = [UW,Data,UW];
% Txsymbols = [UW,UW,Data,UW,UW];
%过时变多径瑞利信道模拟
% Rxsymbols = filter(Channel,Txsymbols); %多径瑞利信道
% Rxsymbols = filter(Channel.PathGains,1,Txsymbols);
% Channel.PathGains
% H_channel0 = fft(Channel.PathGains(LUW+Ldata,:),Ldata);
% if SNR(n)==20 && k==NumFrame
% plot(Channel);
% plot(abs(H_channel0))
% end
%
Channel0 = Channel.PathGains./sqrt(sum((abs(Channel.PathGains)).^2));
Rxsymbols = filter(Channel0,1,Txsymbols); %非时变水声信道
% Rxsymbols = filter(Channel,1,Txsymbols); %普通信道
%噪声功率
Pnoise = 10^(-SNR(n)/10);
Rxsymbols1 = awgn(Rxsymbols,SNR(n),'measured');%过瑞利信道后添加加性高斯白噪声
Rxsymbols0 = awgn(Txsymbols,SNR(n),'measured');%不过瑞利信道添加加性高斯白噪声
%去UW 这里是如何实现的????
%RxData = Rxsymbols(LUW+1:LUW+Ldata);
RxData1 = Rxsymbols1(LUW+1:LUW+Ldata);
RxData0 = Rxsymbols0(LUW+1:Ldata);
% RxData1 = Rxsymbols1(LUW*2+1:LUW*2+Ldata);
% RxData0 = Rxsymbols0(LUW*2+1:LUW*2+Ldata);
if SNR(n)==20 && k==NumFrame
scatterplot(RxData1);
title('均衡前星座图');
grid on;
end;
%信道估计
% RxUW1 = Rxsymbols(1:LUW); %接收前缀
% % RxUW11 = Rxsymbols1(1:LUW); %接收前缀
% % RxUW12 = Rxsymbols1(LUW+1:LUW*2); %接收前缀
%
% % RxUW2 = Rxsymbols(Ldata+LUW+1:Ldata+2*LUW); %接收后缀
% % RxUW21 = Rxsymbols1(Ldata+LUW*2+1:Ldata+3*LUW); %接收后缀
% % RxUW22 = Rxsymbols1(Ldata+LUW*3+1:Ldata+4*LUW); %接收后缀
%
% HRxUW1 = fft(RxUW1); %接收前缀FFT
% % HRxUW11 = fft(RxUW11); %接收前缀FFT
% % HRxUW12 = fft(RxUW12); %接收前缀FFT
%
% % HRxUW2 = fft(RxUW2); %接收后缀FFT
% % HRxUW21 = fft(RxUW21); %接收后缀FFT
% % HRxUW22 = fft(RxUW22); %接收后缀FFT
%
% HTxUW = fft(UW); %发送UW的FFT
%
% HRxUW = HRxUW1; %接收UW的FFT
% % HRxUW = (HRxUW1+HRxUW2)/2; %接收UW的FFT
% % HRxUW = (HRxUW11+HRxUW12+HRxUW21+HRxUW22)/4; %接收UW的FFT
%
% H_esti = HRxUW./HTxUW; %初步估计信道FFT
% h_esti = ifft(H_esti);
% h_estimate = [h_esti,zeros(1,Ldata-LUW)];
% H_estimate = fft(h_estimate); %插值法得到的
评论7