|
9 v7 o+ t1 f( L* M3 E9 R
《测绘学报》 5 v; t, j4 \) ]" U
构建与学术的桥梁 拉近与权威的距离 3 b3 }0 |5 M5 U4 F! n: H3 P
海底底质分类反向散射强度三维概率密度法
8 |/ W! H5 {4 _5 K( W2 u) d; c6 E 金绍华1,2, 李家彪1, 吴自银1 , 边刚2, 崔杨2 5 x2 [! M( h6 `/ E" l2 Z
1. 国家海洋局海底科学重点实验室, 浙江 杭州 310012; 2. 海军大连舰艇学院军事海洋与测绘系, 辽宁 大连 116018 * o" ~9 r/ e) A1 X1 z! Z; I/ e
收稿日期：2017-11-10；修回日期：2018-06-15
, s$ d! ?+ b* `/ {* o 基金项目：国家自然科学基金（41876103；41576105；41374018；41830540） H/ h! c; z7 ^4 ^1 K* u
第一作者简介：金绍华(1978-), 男, 博士, 副教授, 研究方向为海洋大地测量及海洋测绘数据处理。E-mail:jsh_1978@163.com ( i9 [3 P) T2 S9 m( e: ~% u4 d* l- Z
通信作者：吴自银, E-mail: ziyinwu@163.com 1 E8 c/ s: a- K, `% e% k9 v
摘要：反向散射强度随入射角变化信息蕴含海底底质特性，是多波束海底底质分类的重要数据源。目前，多波束海底底质分类主要提取平均反向散射强度随入射角的变化信息或基于声呐图像，没有综合考虑两者的相关信息。本文综合考虑平均反向散射强度随入射角变化信息及其概率密度分布，绘制三维概率密度图，提出了基于反向散射强度三维概率密度图的海底底质分类方法。实例计算结果表明：该方法能直观表达多波束条带覆盖区内不同海底底质类型的数量，有效判断不同海底底质类型的边界，实现海底底质类型的快速识别。
# L( g' j- v9 h 关键词：多波束 反向散射强度 海底入射角 三维概率密度图 海底底质分类 5 w* D- g9 z- o5 R* G8 C' F
3D histogram of backscatter strength for seafloor substrates classification ) D$ p; x6 U+ R: @ V) l
JIN Shaohua1,2, LI Jiabiao1, WU Ziyin1 , BIAN Gang2, CUI Yang2 6 ]5 a, Z8 ~! _2 j: P# s
1. Key Laboratory of Submarine Geosciences of SOA, Hangzhou 310012, China;2. Department of Military Oceanography and Hydrography and Cartography, Dalian Naval Academy, Dalian 116018, China ' U: |) y2 R) X+ {
Foundation support: The National Natural Science Foundation of China (Nos. 41876103; 41576105; 41374018; 41830540)
0 ?4 p/ J' K; W/ n: O2 } First author: JIN Shaohua (1978-), male, PhD, associate professor, majors in and data processing of hydrographic surveying and charting. E-mail: jsh_1978@163.com
4 G% T! ]& G4 Y8 ?$ [ Corresponding author: WU Ziyin, E-mail: ziyinwu@163.com
@1 X0 q- W8 ^ Abstract: Backscatter strength angular response embodies the seafloor substrates classification, which are the important data sources for multibeam seafloor classification. At present, multibeam seafloor classification mainly extract mean backscatter strength angular response or sonar image without considering the relevant information of two aspects. This paper comprehensively analyzes backscatter strength angular response and probability distribution of backscatter strength, portraits 3D histogram, and proposes ways of seafloor classification based on 3D histogram of backscatter strength. Results show that the method can directly express different numbers of seafloor classification within multibeam swath, effectively judge the boundaries and realize the fast identification of different seafloor substrates.
9 V6 a/ S1 _) W' C$ x Key words: multibeam backscatter strength seafloor incidence angle 3D histogram seafloor substrates classification 1 f! v( G; n7 g, @. i
多波束测深系统实现了从“点”“线”测量到条带式、全覆盖、“面”测量的转变，在精度和效率方面都具有传统单波束测深系统无法比拟的优势[1-2]。作为当前海底地形测量的主流产品，多波束测深系统在获取水深数据的同时，还记录了反向散射强度数据。对于多波束测深数据的处理，经过近30多年的发展已经日趋成熟，形成了较为完善的理论体系和数据处理方法[3-5]。而对于反向散射强度的数据处理及其应用方面，尽管国内外学者进行了相关的研究并取得了一定的研究成果，但由于海底反向散射强度的随机性和复杂性，利用多波束反向散射强度进行海底底质分类仍然是国内外研究的热点问题[6-9]。
% O1 R! p g: j8 M 海底反向散射强度的研究主要集中在两个方面：一是从能量的观点出发探索海底平均反向散射强度的变化规律；二是采用随机过程的数学方法研究散射强度的统计特性，如分布函数和能量谱等。多波束海底底质分类也是基于这两方面展开的。在平均反向散射强度变化规律分析方面，文献[10]根据声散射在不同角度下的自然特性，将海底反向散射强度随入射角的变化划分为3个区域，并解释了每一个区域海底反向散射强度的主要贡献量。基于文献[10]的工作，文献[11-12]对实测的多波束海底反向散射强度数据进行了处理，提取平均反向散射强度随入射角变化曲线上的显著特征，以此估计海底的底质类型。针对文献[11]方法容差性差的问题，文献[13]提出了综合描述反向散射强度随入射角变化整体信息的参数拟合法。在散射强度统计特性分析方面，主要采用基本统计量、分位数和直方图、功率谱率、灰度共生矩阵和纹理分析等算法从声呐图像上提取特征参数以描述海底底质类型变化[14-18]。这两种方法各有优缺点，第1种方法突出了平均反向散射强度随入射角的变化信息，忽略了反向散射强度的统计特性；第2种方法反映了反向散射强度的统计特性，但没有顾及与入射角的相关性。如果能在一幅图上既能反映平均反向散射强度随入射角的变化信息，又能观察反向散射强度在不同入射角处的统计特性，将增加海底底质分类的特征信息，有望提高海底底质分类的精度和可靠性。 / ]6 j4 Y' n$ U) R( A4 }' j
基于此，本文探讨综合反映平均反向散射强度随入射角变化信息及反向散射强度统计信息的强度表示方法，构建反向散射强度随入射角变化的三维概率密度图。针对条带上存在不同海底底质类型时，平均反向散射强度和概率密度极大值会出现系统性偏差或较大波动的情况，提出基于概率密度极大值和次极大值比率的底质边界迭代逼近技术，并通过实例验证分类方法的有效性。 ; A* L5 I' ~+ u
另外，与最新研究成果相比[19-21]，本文提出的沿航迹方向的相同入射角的反向散射强度概率密度统计分析方法，通过相同入射角下声强点云的概率曲线分析，规避了海底地形改正模型误差等因素的影响，凸显了海底入射角与反向散射强度的响应关系。该方法更大的优势是可以较为精确地识别不同底质类型的边界，对于海底混合型底质的识别与分类具有重要价值。
0 W W$ A f- M 1 反向散射强度随入射角变化的三维概率密度表示 - P. c4 M* C% @$ ]3 l
一般而言，多波束采用“波束-强度式”记录反向散射强度数据[22]，也就是一个波束内记录多个反向散射强度数据。图 1为Simrad EM3000多波束测深系统实测的1 ping反向散射强度数据。随着波束入射角的增加，波束采样点数量增加。由于反向散射强度的随机特性，在相同入射角处，反向散射强度的幅度变化较大，这导致即使将多ping(10 ping)反向散射强度以点云的形式表示(图 2)，也很难看出其规律性信息。 
9 u% l1 L4 y3 u" r6 ^% U1 D2 `5 i 图 1 1 ping数据Fig. 1 1 ping data 
' Y3 {' D9 I7 M. i 图 2 反向散射强度随入射角变化点云数据Fig. 2 Point cloud data of backscatter strength with incidence angle . r+ M% o( q' j7 D$ C5 t
为了直观显示反向散射强度随入射角变化的趋势性信息和分布情况，构建三维概率密度图，如图 3(a)所示。横坐标表示海底入射角，纵坐标表示反向散射强度，在指定入射角处，反向散射强度出现的频率以图右侧的颜色条表示。由图 3(a)既可以直观地看出反向散射强度随入射角变化的趋势性信息，也可以通过不同入射角处反向强度的宽度分析分布特征。图 3(b)分别绘制了入射角为0、20°、40°和60°时不同反向散射强度出现的频率(横坐标为反向散射强度，纵坐标为出现的频率)。由图上可以看出，在指定入射角处，反向散射强度基本呈正态分布，且随着入射角的增加，方差减小。 
# {% D- s8 P! F" }0 R% n 图 3 反向散射强度概率密度Fig. 3 Probability density graph of backscatter strength 9 {( K W+ Z/ t
当多波束条带覆盖区域存在不同的海底底质类型时，反向散射强度随入射角变化的趋势性信息和分布都会发生变化，体现在三维概率密度图上，会出现两组趋势性信息，如图 4(a)所示。体现在图 4(b)上，在入射角为20°、40°和60°处，会出现多个峰值，以60°入射角为例，存在主峰fmax主和次峰fmax次。这说明反向散射强度概率密度图能直观反映海底底质类型的变化。  5 a$ {) C4 m; N1 x
图 4 不同底质类型反向散射强度的概率密度Fig. 4 Probability density graph of backscatter strength in different types of substrates
/ v2 N3 c- P+ {8 u; c& N9 H5 } 2 理论模型 & Y X( z6 }3 D7 n) T$ x$ F f
反向散射强度随入射角变化的三维概率密度图显示反向散射强度在给定入射角θ处的数量q，用坐标表示为(θ,bs,q)。有研究结果表明[23-25]：在给定入射角处，海底反向散射强度服从K分布，其概率密度函数为 
( H. Y0 Q" ^+ U& l: Q0 I+ n- b (1) ' P4 I; \: n9 I, g6 t
式中，Kv－1(·)为第2类v－1阶修正的贝塞尔函数；v为形状参数；λ为尺度参数。
5 j2 U& x: m) u+ e1 |3 w) O; D 在不同入射角处，反向散射强度的采样率不同，如图 1所示，边缘波束采样率较高(每个波束有上百个采样点)，而中央波束采样率较低(有的波束仅采样1个点)。为了在统一尺度上分析反向散射强度随入射角的变化规律，同时简化研究的分布模型，取每个波束的反向散射强度均值作为该波束的反向散射强度值。设每个波束的采样数为N，则均值为  % j; m. i( v$ o' i) @' X
(2)
2 ^' A' l% ?" q% J, f 根据中心极限定理，只要母体分布的生矩函数存在，且母体分布偏度不过大，则当N>30时，子样均值 
4 F9 L/ ~+ I1 G (3) ! \6 k: H B% Q( B9 i% L
式中，u和σ分别为反向散射强度的均值和方差。假设将多ping反向散射强度以一定间隔划分为M个部分，则在给定入射角θ处，其均值BSmean(θ)和概率密度极大值对应的反向散射强度BSmax(θ)用公式分别表示为  0 X+ u! t$ V: a2 J3 M
(4) 
8 D8 p3 G# n/ Z+ g (5)
1 T9 M/ ?+ n# c# w" S. a$ f. a 由以上分析知，如果多波束条带覆盖区海底底质类型相同， mean(θ)与BSmax(θ)应相等。但实际测量中，由于样本数量有限且受误差的影响，BSmax(θ)相对于BSmean(θ)有微小的波动，如图 3(a)所示。 - A. l& D3 q+ U3 @8 \* |* P7 B+ z
如果多波束条带覆盖区海底底质类型不同，将不再服从单一的正态分布，在给定入射角处，反向散射强度概率密度曲线会出现多个峰值，如图 4(b)所示。此时，BSmean(θ)和BSmax(θ)会有较大的差值，峰值比率(fmax主:fmax次)可近似表示条带覆盖区内不同海底底质类型所占的ping数。 * X) W+ D1 V& S
由以上分析知，可利用概率密度极大值和反向散射强度的差值Δ(θ)统计信息作为海底底质类型划分的依据 . j- `. T' s& m: C+ Y' R0 t
 (6)
9 m/ L' d$ }- B4 W
如果Δ(θ)随入射角变化呈随机波动且振幅较小，则认为ping覆盖区域为相同的底质类型，如果Δ(θ)呈系统性偏差或波动幅度较大，则认为ping覆盖区域底质类型不同，逐步调整ping范围，以划分海底底质类型边界。
+ h- g9 T, @+ b/ y& L 3 海底底质分类的方法流程 1 p8 K0 W6 F8 H# m+ k
海底底质分类流程图如图 5所示，具体计算步骤如下：  ' N. r: z7 ?1 J# c
图 5 海底底质分类流程Fig. 5 Flow chart of seafloor classification
$ f f7 Y, `! V) @. c2 ]" i* n @6 L) Q (1) 绘制多ping反向散射强度随入射角变化三维概率密度图，结合Δ(θ)，判断ping覆盖海底是否为同一底质类型。
7 L3 ~0 _! k0 \ (2) 如果为不同的底质类型，通过在相同入射角θ(一般取40°~60°的某一数值)处的反向散射强度频率峰值fmax主和fmax次的比率确定不同底质类型所占的ping数N1和N2 
* R% I @6 L7 ]0 P' D$ M- S (7) $ z: L) M5 m$ \; P& L/ m
式中，N为ping的总数，假设ping的取值为n1~n2，则N=n2－n1。 % C" q& Q# |+ p( k" |
(3) 以ping数较少的底质类型为基准，即以N2为基准，调整并缩减ping的显示范围，取值为(n1+N－2N2)~n2或n1~(n2－N+2N2)。
) }7 \9 O+ L2 K5 F (4) 以调整后的ping为基础，重复步骤(1)、(2)、(3)，直到fmax主:fmax次=1：1。 ' }, _4 W9 B# @8 z4 ?5 H( ^
(5) 设步骤(4)后，最终显示的ping的总数为N′，显示的ping号为n′1~n′2，则可近似认为底质类型的分界ping号为
% {6 g7 L0 K0 Z (6) 由于刚开始选取的ping数较多，频率峰值并不能完全反映两种底质类型所占的ping数，因此，以上确定的底质类型分界ping号精度相对较低。为了精确确定底质类型的分界线，需减小显示ping的数量，通常取分界ping号 
9 f1 ^6 ]/ f( W P 前后少数ping(一般取100 ping)的数据。重复步骤(1)-(5)，逐步逼近底质类型的分界线。 9 y8 {7 y& ?6 s, {0 ~
4 实例计算
; N- D4 E( O7 ` 图 4所示数据为EM3000多波束测深系统在青岛胶州湾海域获取的某条带1~3000 ping的数据，将其与海底底质直接采样点数据比较结果为：该覆盖区内海底底质类型为淤泥和粗砂，分界ping号为1820，其中淤泥占1820 ping，粗砂占1180 ping。 & g9 ~. [6 i. w
选择图 4(b)中60°区域的反向散射强度频率峰值(fmax主=0.107，fmax次=0.081)计算两种底质类型的ping数比，计算结果为N2=1292，根据步骤(3)，重新显示416~3000 ping的数据，其三维概率密度图如图 6所示。  8 {0 C* N" z- ^( e5 m
图 6 不同海底底质反向散射强度随入射角变化概率密度(415~3000 ping)Fig. 6 Probability density graph of backscatter strength with incidence angle in different types of substrates(415~3000 ping)
6 q3 u' v4 Y1 J) L% M& l 由图 6(b)可以看出，60°倾角处两个峰值频率之间仍存在差值。重复步骤(2)和(3)，依据图 6(b)中的频率峰值((fmax主=0.110，fmax次=0.093))，计算得N2=1184，经计算显示ping号为632~3000 ping，三维概率密度图如图 7所示。  2 [' O q* \2 ^, ~
图 7 不同海底底质反向散射强度随入射角变化概率密度(632~3000 ping)Fig. 7 Probability density graph of backscatter strength with incidence angle in different types of substrates(632~3000 ping)
7 T% U h" j8 S7 v9 e7 j9 ] 经过以上计算，可概略得到两种底质类型的分界ping号为1816。取1816 ping前后100 ping的数据，生成的概率密度图如图 8所示。依据图 8(b)中的频率峰值((fmax主=0.14，fmax次=0.12))，计算得N2=92，分界ping号为1824。与已知两种海底的分类ping号1820相比，本次分类结果只有4 ping的误差，相当于2 m(本次测量船速为5 m/s，采样率为10 ping/s，4 ping相当于航迹方向上的2 m)的误差。  b& `. W7 A# J4 j; y
图 8 不同海底底质反向散射强度随入射角变化概率密度(1716~1916 ping)Fig. 8 Probability density graph of backscatter strength with incidence angle in different types of substrates(1716~1916 ping)
& r% X/ G* u) L% n; d. h 为了进一步说明分类结果的有效性，取分界ping前后50 ping的数据生成概率密度图，如图 9所示。尽管显示的ping数(100 ping)，但仍可看出有两个概率密度极大值，证明存在两种底质类型，这也间接说明了该分类方法具有较高的分辨率和精度。由图 9也可以看出，条带左侧与右侧底质类型分界线可能略有不同，因此，也可以将左侧波束和右侧波束分开进行单独讨论。 
7 l- k; q2 r- X( Q0 I% ? 图 9 不同海底底质反向散射强度随入射角变化概率密度(1774~1874 ping)Fig. 9 Probability density graph of backscatter strength with incidence angle in different types of substrates(1774~1874 ping) 1 i! p# C- ~+ D1 c& ?5 C
图 10中的(a)-(e)分别为图 4、图 6-图 9显示数据Δ(θ)的变化情况。由图中数据可以看出，随着分界ping的不断优化，Δ(θ)呈现0均值的随机变化。  e$ e* @8 y3 ^, ?5 [9 e
图 10 Δ(θ)随入射角变化Fig. 10 Δ(θ) with incidence angle - ?9 Q" L+ m, K/ L& P! |1 `8 R* u% k9 q
为了进一步说明本文分类方法的有效性，选择已知海底底质类型的区域进行多波束底质分类结果验证。利用离散的海底采样点数据格网化后形成的海底底质类型范围如图 11(a)所示，分别为砂、淤泥质粉砂和泥粉砂质砂。0007和0011为EM3000多波束测量的两条测线，利用本文方法确定海底底质边界线后，分别用红色、绿色、天蓝色表示不同的海底底质类型，如图 11(b)所示。由图 11(b)可以看出，多波束底质分类结果与已知海底类型分布基本一致，这充分说明了本文分类方法的有效性。  ) v" ]5 X8 l$ T$ X
图 11 多波束海底底质分类结果Fig. 11 Multibeam seafloor substracts classification results 7 F# d. `" b7 O$ f
5 结论 & r- f; o& s. w+ b9 u
本文构建了多波束反向散射强度随入射角变化的三维概率密度图。该图像既能表达平均反向散射强度随入射角变化的趋势性信息，又能反映反向散射强度的分布特性，可直观表达多波束条带覆盖区内海底底质的种类，便于快速划分海底底质类型数量。对于条带覆盖区内存在多种底质类型的情况，提出了利用给定入射角处反向散射强度频率峰值比率确定不同底质类型ping数的方法，并采用迭代技术，逐步逼近底质类型边界。实例验证结果表明，该方法能有效判别航迹方向上海底底质类型的界线，精度达到10 ping以内。 / \1 V5 h1 \3 C2 [% l0 M/ E8 k F
【引文格式】金绍华, 李家彪, 吴自银, 等. 海底底质分类反向散射强度三维概率密度法. 测绘学报，2019，48(1)：124-131. DOI: 10.11947/j.AGCS.2019.20170631
+ P# b* Y) P+ ^3 z( G0 Y5 F) F  ' z. f% [6 g$ J: j* J$ n8 ?$ M
精
" ]! D7 }- a5 f. L' W5 v: i 彩 1 R& Q; U. V2 f. T8 v5 [7 a% N$ _
回
5 ]- X4 u! S, Q) _+ F 顾
S( {2 y' b: | V3 J6 ? n 测绘学科2019年第1期网刊发布 " P& E; n- p9 r. _0 Z
测绘学科18人上榜，第64批中国博士后基金面上项目出炉， ! Q6 @0 Y5 z$ ?2 B7 Z- C! }
174人！四川、陕西、黑龙江、海南四大省测绘局招聘地学人才！ , {' S Y/ |; g9 ]# h) a2 E5 b
重磅！中共自然资源部党组关于深化科技体制改革提升科技创新效能的实施意见
$ r6 F2 K: n4 w( X# y; M0 K9 q2 | Nature排行榜：中国地学院所TOP 50！
6 \4 i9 V+ v% h 李德仁院士：论遥感卫星商业化的道路—如何从事业走向商业 + K+ V5 a* I+ T! D5 N- n" U3 O
2018年最新自然指数公布，哪些高校表现优异？ 6 K3 J/ p& l h2 ^) J
2018年中科院地学类期刊分区名单（完整版）！
7 z* G3 ?4 \$ c5 G. W 重磅！关于2019年度国家自然科学基金项目申请与结题等有关事项的通告
$ M B. h- v2 J9 d+ b 权威 | 专业 | 学术 | 前沿 0 G2 {' Z# D3 W- P! T: p
微信投稿邮箱 | song_qi_fan@163.com ! G; H0 C6 ], o0 [9 O$ r8 z
微信公众号中搜索「测绘学报」，关注我们，长按上图二维码，关注学术前沿动态。
{+ |; O* R+ z+ l 欢迎加入《测绘学报》作者QQ群： 297834524
8 H+ X$ W9 o" L$ @ ` 进群请备注：姓名+单位+稿件编号
6 C8 S; U! C1 A% }7 |7 j1 `8 s' P, e8 \- @, d( t. s
0 k) k, V( \" A% f3 }8 j
! z# A: J8 o' F( p) C, ^) v q' _0 w" O' r
| 
