源码 https://blog.csdn.net/qq_44394952/article/details/124490764?spm=1001.2014.3001.5502
1.实验目的
(1)学习并熟悉Bellohop模型的原理 (2)使用Bellohop模型产生水声信道,如通信距离、水深、浪高等,进行仿真 (3)实现水声信道调制技术仿真,采用相干检测的方法进行PSK、QAM调制解调,分析水声通信的特点
2.Bellohop模型
Bellohop模型是通过高斯波束跟踪方法,计算水平非均匀环境中的声场。高斯波束跟踪方法对于高频水平变化问题特别有吸引力,这是简正波、波数积分和抛物线模型不可替代的。高斯束射线跟踪法的基本思想是将高斯强度分布与每条声线联系起来,该声线为高斯声束的中心声线,这些声线能较平滑的过渡到声影区,也能较平滑的穿过焦散线,所提供的结果与全波动模型的结果更为一致。在频率范围为600Hz-30kHz时与实验数据及理论模型符合得很好。 Bellohop模型产生信道的过程如下: (1)BELLHOP在通过相关命令来读取可以描述发送端和接收端的环境和几何特性文件。mallab画图程序(plotssp,plotbty,plotbrc,etc.)用来显示每一个输入文件。输入文件有: *.env:描述水体环境,必选; *.ati:描述海面形状,可选; *.bty:描述海底形状,可选; *.trc:描述顶部反射系数,可选; *.brc:描述底部反射系数,可选。 (2)BELLHOP在整体环境文件中通过相关选项来产生不同的输出文件: *.ray:描述声线和本征声线; *.shd:描述声线的传播损失; *.arr:描述声线到达的时间-幅度序列等参数 系统的冲激响应仿真结果图: ![在这里插入图片描述](https://img-blog.csdnimg.cn/49a9c569e55f44a1a641d78bfba8fa6a.png)
3.系统组成及原理
(1)信源产生待传输的序列 信源发射端使用randsrc()产生一定长度的串行的二进制序列 (2)调制 1)16QAM调制 2/4电平转换由二进制变为四进制数字基带信号并进行格雷码映射,取值为-3,-1,1,3。格雷码的特点是邻近的两个信号,其标识仅有一位不同,这样当误符号率一定时,误比特率最小。 插值:在QAM符号之间插入一些0点,和后面的低通滤波联合起来的作用是模拟信号高频调制之后的频谱。分别在I路和Q路中,相邻的两个码字之间添加7个0。 低通:用于调整QAM符号的波形,使两路QAM符号分别经过平方根升余弦滤波器,形成等待调制到高频的信号。 两路QAM多进制信号分别和两路正交载波相乘。 2)4PSK调制 4PSK利用4种不同的相位表征数字信息,将二进制数字序列中每两个比特分为一组,00、01、10、11分别对应相位1/4π、3/4π、5/4π、7/4π。
(3)通过信道 将信号与Bellohop模型水声信道产生的单位冲激响应信号卷积,并通过awgn()函数加入高斯白噪声
(4)解调 1)16QAM解调 接入端接收到的QAM信号分别使用两路相干载波进行解调: 抽样:用于还原为调制时插值之前的QAM符号序列数目。经过插值,4QAM序列的点数是原本序列的4倍,16QAM序列的点数是原本序列的8倍,需要去除这些冗余点。 判决:还原QAM符号序列为二进制码元序列。通过根据最大后验概率准则,得到最小检测距离。 按照调制中的格雷码映射规则,经过“L-2电平变换”转换成还原为原来的I路和Q路码元信息。 对两路信号进行并/串变换,将两路码元数据合并后形成完整的比特流数据, I路信号作为还原后信号中序号为奇数的码元,Q路则作为还原后信号中序号为偶数的码元。 2)4PSK解调 接收端根据判断出的相位,确定发送端发送的数字序列。若收到的相位靠近1/4π,则认为发送端发送的是00;若收到的相位靠近3/4π,则认为发送端发送的是01;若收到的相位靠近5/4π,则认为发送端发送的是10;若收到的相位靠近7/4π,则认为发送端发送的是11。
(5) 计算误码率 使用biterr()函数比较解调后输出的比特流与发送端产生的比特流,计算误码率,绘制仿真误符号率与理论误符号率曲线,进行对比分析。 ![在这里插入图片描述](https://img-blog.csdnimg.cn/0dbd2325bf30465da6fc403dab8453fd.png)
4.模块源代码及各模块产生波形
(1)Bellohop模型水声信道产生单位冲激响应信号 代码:
filename = 'C:\CMST Software\AcTUP v2.2L\AcTUP\Output\Bounce&Bellhop
\Test20211109_00100.arr';
Minimum_range=100 %(接收水听器的水平方向上接收范围最小值,m)
Maximum_range=1000 %(接收水听器的水平方向上接收范围最大值,m)
[ amp1, delay, SrcAngle, RcvrAngle, NumTopBnc, NumBotBnc, narrmat, Pos ]... = read_arrivals_asc(filename);
%%单位冲激响应
[m,n]=size(amp1);
amp = abs(amp1); %取模
x = delay(m,:); %获取第50个接收机的时延和幅值
y = amp(m,:);
%%归一化冲激响应
Amp_Delay = [x;y];
Amp_Delay(:,all(Amp_Delay==0,1))=[]; %去掉0值
Amp_Delay=sortrows(Amp_Delay',1); %按照时延从小到大排序
normDelay = Amp_Delay(:,1)-Amp_Delay(1,1);%归一化时延
normAmp = Amp_Delay(:,2)/Amp_Delay(1,2);%归一化幅度
%距离时延
figure(3)
mum=1:m;
ReRange = Minimum_range+(Maximum_range-Minimum_range)/m*mum;
for i=1:min(narrmat)
plot(delay(:,i),ReRange,'o')
hold on
end
hold off
grid on
colorbar
xlabel('时延(sec)')
ylabel('Range(m)')
title(filename)
波形: (2)16QAM调制解调 1)调制 产生随机二进制序列
%定义待仿真序列的长度 N
global N
N=10000;
% 首先产生随机二进制序列
source=randsrc(1,N,[1,0]);
串并转换
I路:第1,3,5,7...元素
Q路:第2,4,6,8...元素
N=length(source);
a=1:2:N;
source1=source(a);%I路
source2=source(a+1);%Q路
2/4电平变换
a=1:2:N/2;
temp1=source1(a);
temp2=source1(a+1);
y11=temp1*2+temp2;
%I路 用十进制表示二进制
temp1=source2(a);
temp2=source2(a+1);
y22=temp1*2+temp2;
%Q路 用十进制表示二进制
a=1:N/4;
source1=(y11*2-1-4)*1.*cos(2*pi*a);
source2=(y22*2-1-4)*1.*cos(2*pi*a);
格雷码映射
%格雷码映射
source1(y11==0)=-3;
source1(y11==1)=-1;
source1(y11==3)=1;
source1(y11==2)=3;%I路
source2(y22==0)=-3;
source2(y22==1)=-1;
source2(y22==3)=1;
source2(y22==2)=3;%Q路
插值
%两路信号source1、source2进行插值,插值的比例为8
sig_insert1=zeros(1,8*length(source1));
a=1:8:length(sig_insert1);
sig_insert1(a)=source1;
sig_insert2=zeros(1,8*length(source2));
a=1:8:length(sig_insert2);
sig_insert2(a)=source2;
低通
[sig_rcos1,sig_rcos2]=rise_cos(sig_insert1,sig_insert2,0.25,2);
function [y1,y2]=rise_cos(x1,x2,fd,fs)
%生成平方根升余弦滤波器
[yf, tf]=rcosine(fd,fs, 'fir/sqrt');
%使用平方根升余弦滤波器对两路信号进行滤波
[yo1, to1]=
rcosflt(x1, fd,fs,'filter/Fs',yf);
[yo2, to2]=
rcosflt(x2, fd,fs,'filter/Fs',yf);
y1=yo1;
y2=yo2;
将基带信号调制到高频上
f=0.25; %输入信号的频率
hf=2.5; %载波的频率
%用于存储插值后的输入信号
yo1=zeros(1,length(sig_rcos1)*hf/f*10);
yo2=zeros(1,length(sig_rcos2)*hf/f*10);
n=1:length(yo1);
yo1(n)=sig_rcos1(floor((n-1)/(hf/f*10))+1);
yo2(n)=sig_rcos2(floor((n-1)/(hf/f*10))+1);
t=(1:length(yo1))/hf*f/10;
sig_modulate
sig_modulate=yo1.*cos(2*pi*hf*t)-yo2.*sin(2*pi*hf*t);
2)解调 判决
%对x1路信号进行判决
xx1(sig_pick1>=2)=3;
xx1((sig_pick1=0))=1;
xx1((sig_pick1>=-2)&(sig_pick1=-2)&(sig_pick2=0)
demodulated_X(k*2-1)=0;
demodulated_X(k*2)=0;
elseif(real(noised_X(k))=0)
demodulated_X(k*2-1)=0;
demodulated_X(k*2)=1;
elseif(real(noised_X(k)) |