文章目录
- 前言
- 一、Bellhop 简介
- 二、Bellhop 结构
- 1、输入文件
- 2、输出文件
- 三、Bellhop 环境文件
- 1、OPTIONS1
- 2、OPTIONS2
- 3、OPTIONS3
- 4、OPTIONS4
- 5、其他参数
- 四、BELLHOP（Matlab_GUI）实例
- 1、bellhop 工具箱自取
- 2、解压缩工具箱
- 3、设置路径
- 4、添加并包含子文件夹
- 5、运行 BELLHOP
- 五、声速剖面与声线轨迹
- 1、环境文件
- 2、第一步：绘制声速剖面
- 3、第二步：计算声场并检查输入文件
- 4、第三步：绘制声线轨迹
- 六、绘制本征声线
- 1、环境文件
- 2、第一步：计算声场并检查输入文件
- 3、第二步：绘制本征声线
- 总结
前言
本文介绍了 Bellhop 工具箱的具体使用，通过介绍官方 demo 实例，对绘制声速剖面、声线轨迹和本征声线进行了讲解，便于新手快速入门。
一、Bellhop 简介
Bellhop，是一种水声工具箱，可以进行水下声学仿真，通过 env 文件设置环境参数，可得到与之对应的声线信息、多途等仿真；
目前 BELLHOP 可以用 Fortran、MATLAB 和 Python 实现，并可以用于（MAC，Windows 和 Linux）的多种平台。我们在 BELLHOP（MATLAB GUI）中输入海洋环境参数：信道几何结构、声速剖面、海底地形、界面反射损失等得出多径的数目 N，入射角、幅度和延时，进而可以给出系统的冲激响应，从而可以在后续的信道建模的修正及信道估计领域应用提供有效的参照。
二、Bellhop 结构
1、输入文件
BELLHOP 的整体结构图如下图所示：
输入文件有以下几个：
- *.env：描述水体环境，必选。
- *.ati：描述海面形状，可选。
- *.bty：描述海底形状，可选。
- *.trc：描述顶部反射系数，可选。
- *.brc：描述底部反射系数，可选。
- *.ssp：描述声速剖面图，可选。
2、输出文件
输出文件有以下几个：
- *.ray：描述声线和本征声线；
- *.shd：描述声线的传播损失；
- *.arr：描述声线到达的时间-幅度序列等参数
BELLHOP 的结构可以用下图直观的表示出来。
三、Bellhop 环境文件
*.env 文件，这是 Bellhop 调用的环境参数文件，其中包括了设置的海洋环境、计算选项等。
*.env 文件主要有 TITLE 和 OPTIONS1 至 OPTIONS5，各个部分相互独立并完成对水体环境不同方面的描述
1、OPTIONS1
OPTION1 是单引号中的五个字母组成的，下面我们分别介绍其含义。
<1>、OPTIONS1(1)
注： OPTIONS1（1）所指为 OPTIONS1 的第一个字母
主要描述了 BELLHOP 为了计算声速以及沿着声线相关的其他参数而采用的一些插值方法。
下面是为声速剖面进行插值，所包括的字母如下：
| 选项 | 说明 |
|---|---|
| S | 三次样条插值（推荐）注：样条拟合产生外表更光滑的声线轨迹图 |
| C | C 型线性插值 |
| N | N2 线性插值 |
| A | 分析插值，需要 ssp 子程序的调整和模型的进一步重编译 |
| Q | 声场的二次逼近，需要另外编写 *.ssp 文件 |
<2>、OPTIONS1(2)
主要描述海水表面的类型
| 选项 | 说明 |
|---|---|
| V | 表面真空 |
| R | 表面完全刚性的硬物 |
| A | 声学半空间 |
| F | 从一个*.brc 文件读入反射系数 |
如果不选’A’， BELLHOP 只用到前两个参数（每行要以“/” 结尾， 其余参数用默认的），后面的参数是其他模型使用的。
其中在描述为“A”的声学半空间中，SURFACE-LINE 格式如下:
| z-surface | cp-surface | cs-surface | density-surface | AP-surface | AS-surface |
|---|---|---|---|---|---|
| 深度 | 纵波声速 | 横波声速 | 表面密度 | 纵波吸收系数 | 横波吸收系数 |
<3>、OPTIONS1(3)
描述底部的衰减所用的单位
| 选项 | 说明 |
|---|---|
| F | (dB/m)kHz，F 指与频率相关 |
| L | 衰减单位对应于参数损失 损失参量（又称损失切线） |
| M | dB/m，M 指每米 |
| N | Nepers/m |
| Q | 品质因子 |
| W | dB/λ， W 指波长 |
<4>、OPTIONS1(4)
可选参数， 如果描述声音的 Thope Volume 衰减， 要设置为“T”
<5>、OPTIONS1(5)
可选参数， 如果不具体说明，就认为表面是平的；如果具体说明，则设置为’*'。在后一种情况下，表面坐标需要用*.ati 文件来描述， 如可以描述为高斯波浪等所需的海面形状。其结构如下：
参数“ 插值类型” 是一个字符，等于‘L’（ 用于对表面进行线性插值） 或‘C’（ 用于对表面进行曲线插值）；表面的距离单位为 km，表面的深度单位为 m。
在 Bellhop 中不使用参数 nmesh 和 sigmas（ σ），参数 z(nssp)用于检测声速剖面的最后一行。 z( )和 cp( )的数值分别对应于以 m 为单位的深度和以 m/s 为单位的 p 波声速。
2、OPTIONS2
由但括号内的两个字母组成，用来描述水柱下面的环境状况。
<1>、OPTIONS2(1)
描述水柱下面的介质类型
| 选项 | 说明 |
|---|---|
| V | 水柱下为真空 |
| R | 刚性物质 |
| A | 声学半空间， BOTTOM-LINE 格式如下 |
| F | 从一个*.brc 文件读入反射系数 |
其中在描述为“A”的声学半空间中，BOTTOM-LINE 格式如下:
| z-surface | cp-bottom | cs-bottom | density-bottom | AP-bottom | AS-bottom |
|---|---|---|---|---|---|
| 深度 | 纵波声速 | 横波声速 | 表面密度 | 纵波吸收系数 | 横波吸收系数 |
<2>、OPTIONS2(2)
描述底部形状。可以空着（对应于平坦底面），也可以设置为 '’ 。在后一种情况下，底部坐标用.bty 文件来描述，结构如下：
同样， “ 插值类型” 可设置为‘L’或‘C’；底端的距离单位为 km，底端的深度单位为 m。
3、OPTIONS3
用五个字母来描述输出选项，很重要！！！
<1>、OPTIONS3(1)
设定将写入输出文件的信息类型。
| 选项 | 说明 |
|---|---|
| A | 写入振幅和传播时间 |
| E | 写入本征射线坐标 |
| R | 写入射线坐标 |
| C | 写入相干声压 |
| I | 写入非相干声压 |
| S | 写入半相干声压 |
<2>、OPTIONS3(2)
描述计算声压的近似方法。
| 选项 | 说明 |
|---|---|
| G | 采用几何波束（默认值） |
| C | 采用笛卡尔波束 |
| R | 采用中心射线波束 |
| B | 采用高斯波束 |
<3>、OPTIONS3(3)
设定包含波束移位效果等的选项。
| 选项 | 说明 |
|---|---|
| ‘ ’ | 不包括波束移位效果（默认值） |
| S | 包括波束移位效果 |
| * | 使用一个波束类型文件 *.sbp， 格式与 *.ati 和 *.bty 一样 |
<4>、OPTIONS3(4)
描述声源的类型。
| 选项 | 说明 |
|---|---|
| R | 圆柱坐标系中的点声源(default) |
| X | 笛卡尔坐标系中的线声源 |
<5>、OPTIONS3(5)
设置阵列类型。
| 选项 | 说明 |
|---|---|
| R | 直线网络的接收器 |
| I | 不规则网络的接收器 |
整数 nbeams 表示出射角的数目， theta(1) 和 theta(nbeams) 是以度（ °） 为单位设定的第一个和最后一个出射角， 以指向底部的 出 射 角 为 正 值 ， 指 向 表 面 的 出 射 角 为 是 负 值 。 参 数 ray-step(m)、 zmax(m) 和 rmax (km) 定义射线和动态方程积分中的射线步距 ds 和“声线追踪区域（box） ” 的边界范围， “ 追踪区域” 外就停止射线的追踪。
4、OPTIONS4
如果 OPTIONS3 只有一个字母的话，就没有必要在 *.env 文件继续写内容了.否则就要多加两行用来包含波束特性的额外信息，下面为字母说明：
<1>、OPTIONS4(1)
描述波束的类型。
| 选项 | 说明 |
|---|---|
| C | Cerveny 类型 |
| F | 空间填充 |
| M | 最小宽度 |
| W | WKB 波束 |
<2>、OPTIONS4(2)
描述波束曲率的类型。
| 选项 | 说明 |
|---|---|
| D | 曲率翻倍 |
| S | 标准曲率 |
| Z | 零曲率 |
参数 epmult 和 rloop 应该为正实数， isingl、 nImage 和 ibwin 应该为整数。 整数 nImage 可以取 1、 2 或 3。 Component 是单个字符，只有在采用中心射线坐标（OPTIONS3(2)=’R’） 计算声压时才使用；它可以空着（ 将声压写进输出文件） 、 等于’H’（ 将声压的水平分量写入输出文件） 或等于’V’（ 将声压的垂直分量写入输出文件）。
5、其他参数
| 选项 | 说明 |
|---|---|
| NSD | 声源在垂直方向上个数(<51)；超过 50 个声源，声线图将会变得杂乱无章； |
| SD() | 写入本征射线坐标 |
| NRD | 接收水听器的垂直方向上个数 |
| RD() | 接收水听器的深度 |
| NR | 接收水听器的水平方向上个数 |
| R() | 接收水听器的水平接收范围 |
| STEP | 声线跟踪的步长 |
| ZBOX | 接收水听器的最大深度 |
| RBOX | 接收水听器与声源的最大水平距离 |
四、BELLHOP（Matlab_GUI）实例
1、bellhop 工具箱自取
链接：bellhop 工具箱
2、解压缩工具箱
将压缩包解压至 D:\Matlab2019a\toolbox\matlab 目录下
3、设置路径
在 matlab 里面添加路径
4、添加并包含子文件夹
在添加路径的过程中，选择添加并包含子文件夹
5、运行 BELLHOP
运行 D:\Matlab2019a\toolbox\matlab\atWin10_2020_11_4\tests\Munk 目录下的 runtests.m 文件
% Munk profile test cases
% mbp
global units
units = 'km';
%%
figure
plotssp( 'MunkB_ray' ) % 绘制声速剖面
bellhop( 'MunkB_ray' ) % 计算声场并检查输入文件
figure
plotray( 'MunkB_ray' ) % 绘制声线轨迹
bellhop( 'MunkB_eigenray' ) % 计算声场并检查输入文件
figure
plotray( 'MunkB_eigenray' ) % 绘制本征声线
执行结果如下：
下图从左至右依次为 声速剖面，声线轨迹，本征声线。
五、声速剖面与声线轨迹
作为第一个例子， 我们来考虑一个深海场景， Munk 声速剖面。通常我们应该以画出“声速剖面”并计算“声线轨迹”作为开始。输入文件（也称作环境文件）是一个简单的文本文件，可以用任何标准文本编辑器来创建，但扩展名必须是“.env”
此处，考虑 D:\Matlab2019a\toolbox\matlab\atWin10_2020_11_4\tests\Munk\MunkB_ray.env：
1、环境文件
MunkB_ray.env
'Munk profile' ! TITLE
50.0 ! FREQ (Hz)
1 ! NMEDIA
'PVF' ! SSPOPT (Analytic or C-linear interpolation)
51 0.0 5000.0 ! DEPTH of bottom (m)
0.0 1548.52 /
200.0 1530.29 /
250.0 1526.69 /
400.0 1517.78 /
600.0 1509.49 /
800.0 1504.30 /
1000.0 1501.38 /
1200.0 1500.14 /
1400.0 1500.12 /
1600.0 1501.02 /
1800.0 1502.57 /
2000.0 1504.62 /
2200.0 1507.02 /
2400.0 1509.69 /
2600.0 1512.55 /
2800.0 1515.56 /
3000.0 1518.67 /
3200.0 1521.85 /
3400.0 1525.10 /
3600.0 1528.38 /
3800.0 1531.70 /
4000.0 1535.04 /
4200.0 1538.39 /
4400.0 1541.76 /
4600.0 1545.14 /
4800.0 1548.52 /
5000.0 1551.91 /
'A' 0.0
5000.0 1600.00 0.0 1.0 /
2 ! NSD
1000.0 4000 / ! SD(1:NSD) (m)
51 ! NRD
0.0 5000.0 / ! RD(1:NRD) (m)
1001 ! NR
0.0 100.0 / ! R(1:NR ) (km)
'R' ! 'R/C/I/S'
41 ! NBeams
-20.0 20.0 / ! ALPHA1,2 (degrees)
0.0 5500.0 101.0 ! STEP (m), ZBOX (m), RBOX (km)
- 第 2 行：声源频率对基本的声线轨迹并不十分重要。声线与频率无关，不过，频率对声线步长大小有影响，因为程序假定在更高的频率需要画出更精确的声线轨迹。
- 第 3 行：在 BELLHOP 中， NMedia 总是设定为 1。包含此参数是为了与声学工具箱中其他模型相兼容，那些模型能够处理多层介质问题。
- 第 4 行：接下来的顶端选项被设定为“PVF”， 表明使用样条拟合来对声速剖面进行插值； 海面以上模拟成真空；所有衰减值都以dB/mkHz 为单位。选择样条拟合是因为知道本例中的声速剖面变化平滑，在这种情况下，样条拟合将生成更加平滑的声线轨迹图。
- 第 5 行：唯一重要的参数是海底深度（5000 米），它标识了需要读取的声速剖面的最后一行。 BELHOP 不使用（该行的） 前两个参数。
- 第 6 行~第 32 行：接下来，我们看到顺序排列以“深度-声速” 来定义的海洋声速剖面。“声速剖面最后一行必须以上面指定的海底深度值为开始。为了保证与声学工具箱中其他模型相兼容，我们通常用“/”来结束每一行。其他模型还需要衰减、剪切速度和密度作为附加参数， “/”表示停止读取该行并使用默认值。
- 第 33 行~第 34 行：这两行说明底部边界，字母选项“A”表示底部模拟为声弹性半空间。第二行表示该半空间声速为 1600m/s，密度为 1.0（这并不是实际海底底质的情况）。
- 第 35 行~第 40 行：描述了声源深度、接收器深度和接收器的距离。
- 深度单位总是为米，距离单位总是为公里。 对于程序首次运行，我们只生成声线轨迹，因此与接收器位置无关；但是，必须提供接收器位置。注意，此处设定了 51 个接收器深度。通常，为了展示声场，用户只想把接收器深度简单设为均匀分布。为了避免强迫用户必须键入所有这些深度数值， 还可以选择只简单地输入第一个和最后一个数值，并用“/”结束该行。 程序检测到预置的端点，然后通过插值产生一组完整的接收器深度值。 声源和接收器必须位于波导内部。
- 第 41 行：接下来是运行类型（RunType）。对于求解声线轨迹的程序运行，选项设置为“R”。。
- 第 42 行~第 43 行：设定声线扇面， 以声线数（第 42 行）和限制角（第 43 行）（单位： 度(°)） 给出。角度倾斜依循惯例，即：水平发出的声线为零度，向底部发出的声线为正角度。
- 第 44 行：设定声线追踪的步长（单位：米），以及声线追踪的深度与距离范围（box）， 范围外的声线不予追踪。
- 通常， 步长设为 0， BELLHOP 将自动选择水深的 1/10 作为步长。 不管选择什么步长， BELLHOP 均会随着声线追踪动态调整步长，以确保每根声线精确着陆在所有给定声速的深度上。对于给定的声速剖面采样，可以通过设定小于默认值的步长来获得更精确的声线轨迹。
2、第一步：绘制声速剖面
现在已经创建好了输入文件，可以使用 Matlab 程序 plotssp.m，开始来绘制声速剖面曲线。 其运行的 Matlab 命令的语法是：
plotssp('MunkB_ray')
此处 ‘MunkB_ray’ 是 BELLHOP 输入文件的名称。运行后产生的 matlab 命令行窗口输出和绘图窗口输出如下表所示：
>> plotssp('MunkB_ray')
__________________________________________________________________
Munk profile
Frequency = 50 Hz
Number of media = 1
PCHIP approximation to SSP
Attenuation units: dB/mkHz
VACUUM
z (m) alphaR (m/s) betaR rho (g/cm^3) alphaI betaI
( Number of pts = 51 Roughness = 0.00 Depth = 5000.00 )
0.00 1548.52 0.00 1.00 0.0000 0.0000
200.00 1530.29 0.00 1.00 0.0000 0.0000
250.00 1526.69 0.00 1.00 0.0000 0.0000
400.00 1517.78 0.00 1.00 0.0000 0.0000
600.00 1509.49 0.00 1.00 0.0000 0.0000
800.00 1504.30 0.00 1.00 0.0000 0.0000
1000.00 1501.38 0.00 1.00 0.0000 0.0000
1200.00 1500.14 0.00 1.00 0.0000 0.0000
1400.00 1500.12 0.00 1.00 0.0000 0.0000
1600.00 1501.02 0.00 1.00 0.0000 0.0000
1800.00 1502.57 0.00 1.00 0.0000 0.0000
2000.00 1504.62 0.00 1.00 0.0000 0.0000
2200.00 1507.02 0.00 1.00 0.0000 0.0000
2400.00 1509.69 0.00 1.00 0.0000 0.0000
2600.00 1512.55 0.00 1.00 0.0000 0.0000
2800.00 1515.56 0.00 1.00 0.0000 0.0000
3000.00 1518.67 0.00 1.00 0.0000 0.0000
3200.00 1521.85 0.00 1.00 0.0000 0.0000
3400.00 1525.10 0.00 1.00 0.0000 0.0000
3600.00 1528.38 0.00 1.00 0.0000 0.0000
3800.00 1531.70 0.00 1.00 0.0000 0.0000
4000.00 1535.04 0.00 1.00 0.0000 0.0000
4200.00 1538.39 0.00 1.00 0.0000 0.0000
4400.00 1541.76 0.00 1.00 0.0000 0.0000
4600.00 1545.14 0.00 1.00 0.0000 0.0000
4800.00 1548.52 0.00 1.00 0.0000 0.0000
5000.00 1551.91 0.00 1.00 0.0000 0.0000
Number of pts = 51
ACOUSTO-ELASTIC half-space
5000.00 1600.00 0.00 1.00 0.0000 0.0000
_______________________
Number of source depths, NSz = 2
Source depths, Sz (m)
1000.00
4000.00
_______________________
Number of receiver depths, NRz = 51
Receiver depths, Rz (m)
0.00 ... 5000.00
_______________________
Number of receiver ranges, NRr = 1001
Receiver ranges, Rr (km)
0.00 ... 100.00
_______________________
Ray trace run
Geometric hat beams
Number of beams = 41
Beam take-off angles (degrees)
-20.000000 20.000000
Step length, deltas = 0 m
Maximum ray depth, zBox = 5500 m
Maximum ray range, rBox = 101 km
Default step length, deltas = 500 m
No beam shift in effect
Point source (cylindrical coordinates)
>> plotssp('MunkB_ray')
__________________________________________________________________
Munk profile
Frequency = 50 Hz
Number of media = 1
PCHIP approximation to SSP
Attenuation units: dB/mkHz
VACUUM
z (m) alphaR (m/s) betaR rho (g/cm^3) alphaI betaI
( Number of pts = 51 Roughness = 0.00 Depth = 5000.00 )
0.00 1548.52 0.00 1.00 0.0000 0.0000
200.00 1530.29 0.00 1.00 0.0000 0.0000
250.00 1526.69 0.00 1.00 0.0000 0.0000
400.00 1517.78 0.00 1.00 0.0000 0.0000
600.00 1509.49 0.00 1.00 0.0000 0.0000
800.00 1504.30 0.00 1.00 0.0000 0.0000
1000.00 1501.38 0.00 1.00 0.0000 0.0000
1200.00 1500.14 0.00 1.00 0.0000 0.0000
1400.00 1500.12 0.00 1.00 0.0000 0.0000
1600.00 1501.02 0.00 1.00 0.0000 0.0000
1800.00 1502.57 0.00 1.00 0.0000 0.0000
2000.00 1504.62 0.00 1.00 0.0000 0.0000
2200.00 1507.02 0.00 1.00 0.0000 0.0000
2400.00 1509.69 0.00 1.00 0.0000 0.0000
2600.00 1512.55 0.00 1.00 0.0000 0.0000
2800.00 1515.56 0.00 1.00 0.0000 0.0000
3000.00 1518.67 0.00 1.00 0.0000 0.0000
3200.00 1521.85 0.00 1.00 0.0000 0.0000
3400.00 1525.10 0.00 1.00 0.0000 0.0000
3600.00 1528.38 0.00 1.00 0.0000 0.0000
3800.00 1531.70 0.00 1.00 0.0000 0.0000
4000.00 1535.04 0.00 1.00 0.0000 0.0000
4200.00 1538.39 0.00 1.00 0.0000 0.0000
4400.00 1541.76 0.00 1.00 0.0000 0.0000
4600.00 1545.14 0.00 1.00 0.0000 0.0000
4800.00 1548.52 0.00 1.00 0.0000 0.0000
5000.00 1551.91 0.00 1.00 0.0000 0.0000
Number of pts = 51
ACOUSTO-ELASTIC half-space
5000.00 1600.00 0.00 1.00 0.0000 0.0000
_______________________
Number of source depths, NSz = 2
Source depths, Sz (m)
1000.00
4000.00
_______________________
Number of receiver depths, NRz = 51
Receiver depths, Rz (m)
0.00 ... 5000.00
_______________________
Number of receiver ranges, NRr = 1001
Receiver ranges, Rr (km)
0.00 ... 100.00
_______________________
Ray trace run
Geometric hat beams
Number of beams = 41
Beam take-off angles (degrees)
-20.000000 20.000000
Step length, deltas = 0 m
Maximum ray depth, zBox = 5500 m
Maximum ray range, rBox = 101 km
Default step length, deltas = 500 m
No beam shift in effect
Point source (cylindrical coordinates)
3、第二步：计算声场并检查输入文件
在实践中，建议首先尝试以 BELLHOP 运行输入文件。BELLHOP 会生成一个打印文件，通过检查打印文件，任何格式错误我们通常都能清楚地发现。
运行 BELLHOP 的 MATLAB 命令是：
bellhop( 'MunkB_ray' )
此处 “MunkB_ray.env” 是输入文件的名称，BELLHOP 调用的计算文件为：D:\Matlab2019a\toolbox\matlab\atWin10_2020_11_4\at\bin\bellhop.exe。
假设成功完成，BELHOP 生成一个名为 “MunkB_ray.prt” 的打印文件和一个名为 “MunkB_ray.ray” 的声线文件。仔细检查打印文件，以查验问题是按照预想场景设定， BELLHOP 也完成运行。后者可以通过检查打印文件中没有错误消息来验证，打印文件的最后一行显示的是 CPU耗时。
回传输入数据的打印文件 MunkB_ray.prt
BELLHOP/BELLHOP3D
BELLHOP- Munk profile
frequency = 50.00 Hz
Dummy parameter NMedia = 1
PCHIP approximation to SSP
Attenuation units: dB/mkHz
VACUUM
Depth = 5000.00 m
Sound speed profile:
z (m) alphaR (m/s) betaR rho (g/cm^3) alphaI betaI
0.00 1548.52 0.00 1.00 0.0000 0.0000
200.00 1530.29 0.00 1.00 0.0000 0.0000
250.00 1526.69 0.00 1.00 0.0000 0.0000
400.00 1517.78 0.00 1.00 0.0000 0.0000
600.00 1509.49 0.00 1.00 0.0000 0.0000
800.00 1504.30 0.00 1.00 0.0000 0.0000
1000.00 1501.38 0.00 1.00 0.0000 0.0000
1200.00 1500.14 0.00 1.00 0.0000 0.0000
1400.00 1500.12 0.00 1.00 0.0000 0.0000
1600.00 1501.02 0.00 1.00 0.0000 0.0000
1800.00 1502.57 0.00 1.00 0.0000 0.0000
2000.00 1504.62 0.00 1.00 0.0000 0.0000
2200.00 1507.02 0.00 1.00 0.0000 0.0000
2400.00 1509.69 0.00 1.00 0.0000 0.0000
2600.00 1512.55 0.00 1.00 0.0000 0.0000
2800.00 1515.56 0.00 1.00 0.0000 0.0000
3000.00 1518.67 0.00 1.00 0.0000 0.0000
3200.00 1521.85 0.00 1.00 0.0000 0.0000
3400.00 1525.10 0.00 1.00 0.0000 0.0000
3600.00 1528.38 0.00 1.00 0.0000 0.0000
3800.00 1531.70 0.00 1.00 0.0000 0.0000
4000.00 1535.04 0.00 1.00 0.0000 0.0000
4200.00 1538.39 0.00 1.00 0.0000 0.0000
4400.00 1541.76 0.00 1.00 0.0000 0.0000
4600.00 1545.14 0.00 1.00 0.0000 0.0000
4800.00 1548.52 0.00 1.00 0.0000 0.0000
5000.00 1551.91 0.00 1.00 0.0000 0.0000
( RMS roughness = 0.00 )
ACOUSTO-ELASTIC half-space
5000.00 1600.00 0.00 1.00 0.0000 0.0000
__________________________________________________________________________
Number of Source depths, Sz = 2
Source depths, Sz (m)
1000.00 4000.00
__________________________________________________________________________
Number of Receiver depths, Rz = 51
Receiver depths, Rz (m)
0.00000 100.000 200.000 300.000 400.000
500.000 600.000 700.000 800.000 900.000
...
5000.00
__________________________________________________________________________
Number of Receiver ranges, Rr = 1001
Receiver ranges, Rr (km)
0.00000 0.100000 0.200000 0.300000 0.400000
0.500000 0.600000 0.700000 0.800000 0.900000
...
100.000
Ray trace run
Geometric hat beams in Cartesian coordinates
Point source (cylindrical coordinates)
Rectilinear receiver grid: Receivers at Rr( : ) x Rz( : )
__________________________________________________________________________
Number of beams in elevation = 41
Beam take-off angles (degrees)
-20.0000 -19.0000 -18.0000 -17.0000 -16.0000
-15.0000 -14.0000 -13.0000 -12.0000 -11.0000
...
20.0000
__________________________________________________________________________
Step length, deltas = 0.000 m
Maximum ray depth, Box%z = 5500. m
Maximum ray range, Box%r = 101.0 km
No beam shift in effect
Step length, deltas = 500.0 m (automatically selected)
Tracing beam 1 -20.00
Tracing beam 2 -19.00
Tracing beam 3 -18.00
Tracing beam 4 -17.00
Tracing beam 5 -16.00
Tracing beam 6 -15.00
Tracing beam 7 -14.00
Tracing beam 8 -13.00
Tracing beam 9 -12.00
Tracing beam 10 -11.00
Tracing beam 11 -10.00
Tracing beam 12 -9.00
Tracing beam 13 -8.00
Tracing beam 14 -7.00
Tracing beam 15 -6.00
Tracing beam 16 -5.00
Tracing beam 17 -4.00
Tracing beam 18 -3.00
Tracing beam 19 -2.00
Tracing beam 20 -1.00
Tracing beam 21 0.00
Tracing beam 22 1.00
Tracing beam 23 2.00
Tracing beam 24 3.00
Tracing beam 25 4.00
Tracing beam 26 5.00
Tracing beam 27 6.00
Tracing beam 28 7.00
Tracing beam 29 8.00
Tracing beam 30 9.00
Tracing beam 31 10.00
Tracing beam 32 11.00
Tracing beam 33 12.00
Tracing beam 34 13.00
Tracing beam 35 14.00
Tracing beam 36 15.00
Tracing beam 37 16.00
Tracing beam 38 17.00
Tracing beam 39 18.00
Tracing beam 40 19.00
Tracing beam 41 20.00
Tracing beam 1 -20.00
Tracing beam 2 -19.00
Tracing beam 3 -18.00
Tracing beam 4 -17.00
Tracing beam 5 -16.00
Tracing beam 6 -15.00
Tracing beam 7 -14.00
Tracing beam 8 -13.00
Tracing beam 9 -12.00
Tracing beam 10 -11.00
Tracing beam 11 -10.00
Tracing beam 12 -9.00
Tracing beam 13 -8.00
Tracing beam 14 -7.00
Tracing beam 15 -6.00
Tracing beam 16 -5.00
Tracing beam 17 -4.00
Tracing beam 18 -3.00
Tracing beam 19 -2.00
Tracing beam 20 -1.00
Tracing beam 21 0.00
Tracing beam 22 1.00
Tracing beam 23 2.00
Tracing beam 24 3.00
Tracing beam 25 4.00
Tracing beam 26 5.00
Tracing beam 27 6.00
Tracing beam 28 7.00
Tracing beam 29 8.00
Tracing beam 30 9.00
Tracing beam 31 10.00
Tracing beam 32 11.00
Tracing beam 33 12.00
Tracing beam 34 13.00
Tracing beam 35 14.00
Tracing beam 36 15.00
Tracing beam 37 16.00
Tracing beam 38 17.00
Tracing beam 39 18.00
Tracing beam 40 19.00
Tracing beam 41 20.00
CPU Time = 0.781E-01s
4、第三步：绘制声线轨迹
下一步是应用 MATLAB 命令绘制声线：
plotray( 'MunkB_ray' )
请注意，距离轴单位为米。如果更喜欢千米，那么就简单地设置 MATLAB 全局变量：global units; units = 'km';
依据声线是触碰单侧边界还是触碰两侧边界， 在绘图中使用不同颜色来表达。表面反弹和底端反弹数量均写入声线文件中，因此方便人们简单地修改 plotray 文件中的着色代码，以便以最好的方式来描述声线传播的物理机理。
六、绘制本征声线
BELLHOP 也能生成本征声线绘图，用以展示仅仅只连接声源到接收器的声线。只需将“运行类型（RunType）” 更改为“E”，即可完成此任务。
对于实际的本征声线计算，我们应采用默认波束，它具有由相邻声线形成的声线管所定义的波束宽度，我们称其为几何波束。当采用默认波束类型时，所写入的声线将仅仅只是穿过接收器位置的声线。
通常需要采用更精细的扇面。 比如，在前面的例子中，如果我们设定只计算 41 根声线，那么经过远距离传播后，这些声线就很发散了。当要保存穿过接收器的声线时，这些声线可能基本上都错过了接收器的位置。因此，本例中，我们将声线数设置为 5001。声线数设置得越多，本征声线就计算得越准确。当然，运行时间也相应地递增。
我们通常只用一个声源和一个接收器来执行本征声线计算。否则，得到的声线图会很凌乱。 输入文件“MunkB_eigenray.env” 经过这些改变后列入下表。
1、环境文件
MunkB_eigenray.env
'Munk profile' ! TITLE
50.0 ! FREQ (Hz)
1 ! NMEDIA
'CVF' ! SSPOPT (Analytic or C-linear interpolation)
51 0.0 5000.0 ! DEPTH of bottom (m)
0.0 1548.52 /
200.0 1530.29 /
250.0 1526.69 /
400.0 1517.78 /
600.0 1509.49 /
800.0 1504.30 /
1000.0 1501.38 /
1200.0 1500.14 /
1400.0 1500.12 /
1600.0 1501.02 /
1800.0 1502.57 /
2000.0 1504.62 /
2200.0 1507.02 /
2400.0 1509.69 /
2600.0 1512.55 /
2800.0 1515.56 /
3000.0 1518.67 /
3200.0 1521.85 /
3400.0 1525.10 /
3600.0 1528.38 /
3800.0 1531.70 /
4000.0 1535.04 /
4200.0 1538.39 /
4400.0 1541.76 /
4600.0 1545.14 /
4800.0 1548.52 /
5000.0 1551.91 /
'A' 0.0
5000.0 1600.00 0.0 1.0 /
1 ! NSD
1000.0 / ! SD(1:NSD) (m)
1 ! NRD
800.0 / ! RD(1:NRD) (m)
1 ! NR
100.0 / ! R(1:NR ) (km)
'E' ! 'R/C/I/S'
5001 ! NBeams
-25.0 25.0 / ! ALPHA1,2 (degrees)
0.0 5500.0 101.0 ! STEP (m), ZBOX (m), RBOX (km)
2、第一步：计算声场并检查输入文件
运行 BELLHOP 的 MATLAB 命令是：
bellhop( 'MunkB_eigenray' )
此处 “MunkB_eigenray.env” 是输入文件的名称。
3、第二步：绘制本征声线
本征声线用 plotray 命令绘制：
figure
plotray( 'MunkB_eigenray' )
总结
以上就是全部内容，本文仅仅简单介绍了 bellhop 的简单使用，对绘制声速剖面、声线轨迹和本征声线进行了讲解，后面会继续深入学习研究并分享。
我的qq：2442391036，欢迎交流！
