一种基于摄影测量的小水槽波浪要素的提取方法

：
1.本发明涉及一种摄影测量的技术领域，尤其涉及一种基于摄影测量的小水槽波浪要素提取方法。


背景技术：

2.海浪在海洋活动中是非常常见的一种表现形式。在大多数情况下，海况都受到风浪和涌浪的共同影响，但风浪和涌浪的外在表现形式和内部机制差别很大。
3.在不同沿海海域和受到不同气象环境作用下，风浪和涌浪在海浪谱中通常以双峰谱或者多峰谱的表现形式居多。海况由于受到风浪和涌浪的共同影响，将导致海浪参数中的波向、波高和波周期等参数改变了原来的物理意义，因此风涌的研究就变得尤为迫切。
4.风浪在风直接作用下产生的水面波动是一个极度复杂的过程，要严格加以定量分析是十分困难的。风浪的生成和成长的机制是海浪研究中最基本的问题。
5.涌浪是风浪离开风吹的区域后所形成的波浪。涌浪与风浪相比，具有较规则的外形，排列比较整齐，波峰线较长，波面较平滑，比较接近于正弦波的形状，而涌浪可以看作是由许多振幅不等、频率不等、传播方向不同并具有随机初相位的正弦波的分量叠加而成。
6.随着海洋科学及海洋技术的发展，海洋开发的热度提升，越来越多的海洋工程及水动力实验得到了重视，许多实验室模拟真实海洋的波浪海流及潮汐的实验系统也在逐步完善，其中小型的波浪水槽相对于大型造波水槽有着占地面积小、造价低、移动方便等优点，越来越多的小型开敞式波浪水槽走进实验室，它们可以通过机械造波板的往复运动模拟出不同形态的波浪，并提供了模拟真实海洋进行实验探究的可能。
7.一种波浪视频测量装置及方法的测量原理是通过视频图像采集模块同步拍摄水面预置的浮体运动序列图像，利用图像识别处理方法提取图像坐标系浮体中心，基于双目视觉摄影测量原理计算浮体世界坐标系中三维坐标，进而按照《海洋观测规范》中海浪观测相关要求反演波浪参数。
8.其中波高﹑波周期波浪参数反演是通过获取一段时间内浮体世界坐标系三维坐标，绘制随时间变化的模拟曲线，利用所述波浪参数显示模块实时显示数据接收处理模块计算得到的波高、波周期的波浪参数。
9.海面测量系统、海面测量方法以及存储介质中的海面测量方法是通过设置于不同地点处的多个摄像机同时拍摄包括同一海面区域的图像，从通过多个摄像机同时拍摄到的各个图像中提取波浪，针对提取出的各波浪，从通过多个摄像机同时拍摄到的各个图像中找出相同的波浪并进行对应，通过对被对应起来的图像中的各波浪进行三维分析，来计算海面的高度和波浪的高度。
10.现今实验室波浪水槽在实验中应用广泛，能够在一定程度上模拟海面的环境，波浪水槽的造波器能够产生不同形态的波浪，此外，还配备有产生风的仪器。
11.对于小型开敞式波浪水槽所制造出的波浪，根据实验所需，研究者有时需要得到所造波浪的波浪、波长、周期和频率等波浪要素数据进行分析研究。在具体实施过程中，由
于波浪各要素会随着造波器频率的变化而变化，获取实时的波浪形态并不好采集，这些仪器所造出的波浪可以做到“定性”，却无法达到“定量”的要求，且波浪要素随着水槽的规模、水深、水体的性质、造波器振幅等实验条件的改变而改变。而在实验中，往往需要获取准确的波浪要素信息才能为实验提供支持。


技术实现要素：

12.本发明针对现有技术的不足和缺陷，提供了一种基于摄影测量的小水槽波浪要素的提取方法，以便于实现对风浪和涌浪的观察、数据获取与分析。
13.为实现上述目的，本发明提供了如下技术方案：一种基于摄影测量的小水槽波浪要素的提取方法，其特征在于，包括水槽、造波器、吹风设备、消波板和摄像机；造波器制造波浪，吹风设备制造风浪，通过水槽的一侧的造波器造波，造波器下方是一个楔形块体，造波器以设定频率上下移动楔形块体，楔形块体的斜面与水面接触，对水面起力的作用，做上下运动，形成波浪；在水槽远离造波器的另一侧放置消波板，消波板削弱波浪遇水槽壁后产生的反射波所造成的误差，在水中加入颜色指示剂，以便水体与空气界面的灰度产生差异。
14.包括如下步骤：
15.100：视频录制与图像提取；
16.200：rgb数据提取；
17.300：灰度处理；
18.400：欧拉法提取波形、波高和周期；
19.500：拉格朗日法提取波速和波长；
20.600：海浪频谱提取并获取谱峰周期、平均周期及平均波高。
21.在一个实施例中，在步骤100中还包括，在所述水槽造出的波浪通过架设的所述摄像机完成采集，拍摄视频时从水槽的一侧横向水平拍摄，尽量避免左右倾斜所造成的误差。
22.在一个实施例中，在步骤200中还包括，通过读取视频，保存单帧图像，获得一定时间间隔内的波浪现场图像。为了提取其中的信息，将读取到的单帧图像转化为rgb的三分量灰度图。
23.在一个实施例中，在步骤300中还包括，利用灰度图像中波浪的水体与空气界面的灰度差异，对图像灰度进行边缘处理，使得图像颜色呈现两极化，突出边界，获得清晰的波线。判断水体与空气的灰度阈值，选取临界值为40，将颜色较深的水体，即rgb值小于40的像素点，均令为黑色为0；将颜色较浅的空气，即rgb值大于40的像素的，均令为白色为255。将颜色进行两极化处理之后，即可以得到边界分明的波浪图像。对灰度处理后的图像进行分析，从图片的左上角开始，从上到下从左到右读取每个像素点的灰度值，找到第一个黑色像素点的坐标，得到第一列数据的水体与空气分界点，依次类推得到其他列的第一个黑色像素点坐标，将各个第一黑色像素点连接得到水体与空气的分界线，即可得到单帧图像的准确波线。
24.在一个实施例中，在步骤400中还包括：选定图像中的一列像素点，寻找水体与空气的界面点，得到间隔δn张图像的时间间隔δt的离散点，δt＝fps*δn，通过等比例转换，选取每张单帧图像的中心轴线来观察波面位置，记录相应的像素坐标，即可绘制中心轴线处的波形随时间变化图。fps是frames per second的缩写，为帧速率。
25.在拍摄视频时，预先测得了静止水平面高度，当波面从静止水平面以下上升至静止水平面以上时，波面与静止水平面的交点即为上跨零点，反之，当波面从静止水平面以上下降至静止水平面以下时，波面与静止水平面的交点为下跨零点。利用上跨零点或下跨零点计算波周期：相邻两个上跨零点所对应的时间分别为t
01
、t
02
，即可得到上跨零点周期t0＝t
02-t
01
。
26.利用波峰和波谷计算波周期：相邻两个波峰或波谷所对应的时间分别为t
11
、t
12
，即可得到周期t1＝t
12-t
11
。波高是相邻的波峰和波谷间的垂直距离，获取波峰、波谷的像素坐标(x1,y1)、(x2,y2)，通过作差的方式得到某一波高h1＝|y1–
y2|，将得到的所有波按照波高从大到小排列，取最大的前1/3个波的波高做平均即可得到有效波高h
1/3
，其对应的周期做平均即为有效周期t
1/3
，取最大的前1/10个波的波高做平均即可得到1/10大波波高，其对应的周期平均为1/10大波周期。对所有的波高hi进行算数平均即可得到平均波高hm，对所有的周期ti进行算数平均即可得到平均周期tm。
27.在一个实施例中，在步骤500中还包括：根据在步骤400中计算得到的平均周期tm，在单个周期时间内选取两帧不相邻的图像，时间分别为t
21
、t
22
，其波峰对应横轴坐标分别为x
p1
、x
p
2，使用步骤400中的方法计算得到t
21
、t
22
所在的波周期t2，通过计算不同帧波线同一波峰或波谷的横坐标差值与时间差的比值，提取波速c、波长l：
[0028][0029]
l＝c
×
t2[0030]
再结合步骤400中的方法计算得到t
21
、t
22
所在波的波高h2，计算得到波陡：
[0031][0032]
自第一帧图片起，通过分别取间隔为1/2平均周期t
m1
的两帧图片计算波速与波长，即使用第1帧与第1 fps*t
m1
/2帧、第2帧与第2 fps*t
m1
/2帧，以此类推，计算得到相应波速ci与波长li后取平均值，得到平均波速cm与平均波长lm：
[0033][0034][0035]
fps-帧速率。
[0036]
在一个实施例中，在步骤600中还包括：海浪频谱提取并获取谱峰周期、平均周期及平均波高，海浪可以看作由多个正弦波叠加而成，每个波具有不同的圆频率与能量。以圆频率为自变量，能量为因变量可以绘制海浪频谱，从而展示海浪分解出的组成波的能量关于频率的信息，得到海浪频谱并获取谱峰周期、平均周期及平均波高。
[0037]
本发明的主要有益效果在于：通过实验室波浪水槽和运用造波器和吹风设备，在不同条件下模拟波浪，通过方便有效的录制视频的手段。
[0038]
利用计算机读取物理模型水槽视频中的波浪曲线及其有效数据，使得实验数据采集更规范准确。
[0039]
对风浪涌浪的物理模拟进行可视化处理，直接获取波浪要素数据。以便于实现对风浪、涌浪的观察、数据获取与分析，及对波浪生成、成长、衰减破碎的过程，进行深入探讨研究。
附图说明：
[0040]
图1揭示了本发明一实施例中，一种基于摄影测量的小水槽波浪要素的提取方法的流程图；
[0041]
图2揭示了本发明一实施例中，一种基于摄影测量的小水槽波浪要素的提取方法的整体结构布置示意图；其中，图2a为主视图，图2b为图2a的侧视图；
[0042]
图3揭示了本发明一实施例中，rgb值绘制频率统计图；
[0043]
图4揭示了本发明一实施例中，颜色两极化处理过程的波浪图像；其中图4a为原始图像，图4b为灰度处理后的图像；
[0044]
图5揭示了本发明一实施例中，单帧图像的准确波线图；
[0045]
图6揭示了本发明一实施例中，定点观测波形图；
[0046]
图7揭示了本发明一实施例中，计算平均上跨零点周期的示意图；
[0047]
图8揭示了本发明一实施例中，不同帧波线的同一波峰或波谷提取波速和波长的示意图；
[0048]
图9揭示了本发明一实施例中，风浪组成波的能量关于频率的信息，横坐标为频率，纵坐标为能量的风浪频谱图。
[0049]
在图示中：
[0050]
1-造波器，
[0051]
2-水槽，
[0052]
3-消波板，
[0053]
4-相机。
具体实施方式：
[0054]
以下由特定的具体实施例说明本发明的实施方式，本领域技术人员可由本说明书所揭示的内容轻易地了解本发明的其他优点及功效。虽然本发明的描述将结合较佳实施例一起介绍，但这并不代表此发明的特征仅限于该实施方式。恰恰相反，结合实施方式作发明介绍的目的是为了覆盖基于本发明的权利要求而有可能延伸出的其它选择或改造。为了提供对本发明的深度了解，以下描述中将包含许多具体的细节。本发明也可以不使用这些细节实施。此外，为了避免混乱或模糊本发明的重点，有些具体细节将在描述中被省略。需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。
[0055]
在本实施例的描述中，需要说明的是，术语“上”、“下”、“内”、“底”等指示的方位或位置关系为基于附图所示的方位或位置关系，或者是该发明产品使用时惯常摆放的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限制。
[0056]
术语“第一”、“第二”等仅用于区分描述，而不能理解为指示或暗示相对重要性。
[0057]
在本实施例的描述中，还需要说明的是，除非另有明确的规定和限定，术语“设
置”、“相连”、“连接”应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以具体情况理解上述术语在本实施例中的具体含义。
[0058]
参考图1并结合图2-图9，在一实施例中的一种基于摄影测量的小水槽波浪要素的提取方法，包括如下步骤；
[0059]
100：视频录制。
[0060]
使用150cm*10cm*20cm的波浪水槽、bmud30-c2造波器、2000w吹风设备、65cm*10cm的消波板和gopro7运动摄像机等器材。使用造波器制造波浪，使用吹风设备制造风浪，消波板防止水槽壁的反射波对波的影响，使用运动摄像机进行视频的录制。
[0061]
首先向水槽加水，水位在水槽高度的1/3到2/3之间，向水槽内加入颜色指示剂，以便水体与空气界面的灰度差异更加明显，本实验中以蓝色为例，该颜色不易在水槽内壁留下痕迹，也可尝试其他颜色达到相同效果。
[0062]
然后架设拍摄设备，该设备需与波面在同一水平面，且摄像机镜头与水槽纵截面平行，从波浪水槽一侧横向水平拍摄，尽量避免左右倾斜所造成的误差，拍摄得到静止状态下的水平面高度。保持摄像机位置不变，为了后期计算得到真实的波浪参数，在视频开头将3cm*3cm的方格纸粘贴至水槽外壁，获取带有比例尺信息的视频。在视频中提取方格纸所占像素数量num，设置rule为3cm，得到视频中单个像素长度所代表的实际长度pixel＝rule/num，从而用于对由视频获得的波浪参数进行比例转换，得到真实数值。然后，在不挪动摄像机的前提下将比例尺移开，继续拍摄水槽波浪视频。
[0063]
通过水槽一侧的造波器造波，造波器下方是一个楔形的块体，造波器以设定频率上下移动块体，楔形块体斜面与水面接触，由此对水面起到力的作用，从而使水做既定波运动，形成波浪。在水槽的另一侧放置一块消波板，通过消波板削弱波浪遇水槽壁后产生的反射波所造成的误差。
[0064]
采集视频资料时，读取视频图像的尺寸大小(长*宽)、时长t、帧速率frame rate等参数，并获取视频的总帧数，继而利用视频采集图像，将所需帧数间隔的图像采集出来，该采集图像的原理就是从视频文件中的所有帧读取数据，显示每一帧，然后以视频的帧速率播放所有帧，并保存每一帧。设置视频分辨率为1080p，帧速率为29.98，设置一个循环来将视频的每一帧保存成jpg格式的图片，图片按顺序编号从第一帧到最后一帧，每次读取一帧可防止内存不足，每一张图片对应视频中相应的一帧。
[0065]
200：rgb数据提取。
[0066]
通过读取图像，获得r、g、b三个图层的信息，用n*m*3的矩阵表示出来。
[0067]
矩阵中n*m表示图像大小，矩阵内的数字表示的是各图层的rgb数值大小。数值的大小范围为0-255，指的是图像的灰度，其数值越大，灰度越大，即当三个图层均为255时为白色。为了提取波浪信息，将从步骤100中读取到的单帧图像转化为rgb三分量灰度图。
[0068]
300：灰度处理(颜色极化)。
[0069]
利用灰度图像中波浪的水体与空气界面的灰度差异，对图像灰度进行边缘处理，使得图像颜色呈现两极化，突出边界，获得清晰的波线。初步判断水体与空气的灰度阈值，选取临界值为40，将颜色较深的水体，即rgb值小于40的像素点，均令为黑色为0；将颜色较
浅的空气，即rgb值大于40的像素的，均令为白色为255。将颜色进行两极化处理之后，即可以得到边界分明的波浪图像。对灰度处理后的图像进行分析，从图片的左上角开始，从上到下从左到右读取每个像素点的灰度值，找到第一个黑色像素点的坐标，得到第一列数据的水体与空气分界点，依次类推得到其他列的第一个黑色像素点坐标，将各个第一黑色像素点连接得到水体与空气的分界线，即可得到单帧图像的准确波线。
[0070]
400：欧拉法提取波形、波高、周期。
[0071]
欧拉法是指不关注于单个流体质点，而是以某一流场为对象，观察该流场内的流体运动状态与特性。因为拍摄范围内很难获取同一时间点上完整的一个或多个波形，因此使用欧拉法对视频图像进行分析，观测波面位移随时间变化。选定图像中的一列像素点，寻找水体与空气的界面点，得到间隔δn张图像的时间间隔δt的离散点，δt＝fps*δn，fps为帧速率(frames per second)。通过等比例转换，即可获得定点观测波形随时间变化图。例如，选取每张单帧图像的中心轴线来观察波面位置，记录相应的像素坐标，即可绘制中心轴线处的波形随时间变化图。在步骤100中预先测得了静止水平面高度，当波面从静止水平面以下上升至静止水平面以上时，波面与静止水平面的交点即为上跨零点，反之，当波面从静止水平面以上下降至静止水平面以下时，波面与静止水平面的交点为下跨零点。利用上跨零点或下跨零点计算波周期：相邻两个上跨零点所对应的时间分别为t
01
、t
02
，即可得到上跨零点周期t0＝t
02-t
01
。
[0072]
对所有的上跨零点周期t0进行算数平均即可得到平均上跨零点周期t
01
：
[0073][0074]
利用波峰和波谷计算波周期：相邻两个波峰或波谷所对应的时间分别为t
11
、t
12
，即可得到周期t1＝t
12-t
11
。波高是相邻的波峰和波谷间的垂直距离，获取波峰、波谷的像素坐标(x1,y1)、(x2,y2)，通过作差的方式得到某一波高h1＝|y1–
y2|，将得到的所有波按照波高从大到小排列，取最大的前1/3个波的波高做平均即可得到有效波高h
1/3
：
[0075][0076]
其中n1为总个数。
[0077]
有效周期t
1/3
：
[0078][0079]
取最大的前1/10个波的波高做平均即可得到1/10大波波高，其对应的周期平均为1/10大波周期。对所有的波高hi进行算数平均即可得到平均波高hm：
[0080]
[0081]
对所有的周期ti进行算数平均即可得到平均周期tm：
[0082][0083]
500：拉格朗日法提取波速、波长。
[0084]
拉格朗日法指的是以某一流体质点为观测对象，观察水质点的运动状态以及运动中物理量随时间的变化。为获得波速与波长等波浪要素，使用拉格朗日法观察图像，即专注于波形上某一点不变，观察波形随时间传播在图像上的位置变化情况。
[0085]
本发明所关注的是单帧图像上的波峰/波谷点。根据在步骤400中计算得到的平均周期tm，在单个周期时间内选取两帧不相邻的图像，其对应时间分别为t
21
、t
22
，其波峰对应横轴坐标分别为x
p1
、x
p2
，使用步骤400中的方法计算得到t
21
、t
22
所在的波周期t2，通过计算不同帧波线同一波峰或波谷的横坐标差值与时间差的比值，提取波速c、波长l：
[0086][0087]
l＝c
×
t2[0088]
再结合步骤400中的方法计算得到t
21
、t
22
所在波的波高h2，计算得到波陡：
[0089][0090]
本发明自第一帧图片起，通过分别取间隔为1/2平均周期t
m1
的两帧图片计算波速与波长，即使用第1帧与第1 fps*t
m1
/2帧、第2帧与第2 fps*t
m1
/2帧，以此类推，计算得到相应波速ci与波长li后取平均值，得到平均波速cm与平均波长lm：
[0091][0092][0093]
fps-帧速率。
[0094]
600：海浪频谱提取并获取谱峰周期、平均周期及平均波高。
[0095]
海浪可以看作由多个正弦波叠加而成，每个波具有不同的圆频率与能量。以圆频率为自变量，能量为因变量可以绘制海浪频谱，从而展示海浪分解出的组成波的能量关于频率的信息。
[0096]
根据步骤300获取的定点波面随时间变化的离散数据记录，利用中轴线上的波面位移随时间的变化y(t)，进行离散傅里叶变换y1(f)：
[0097][0098]
其中，为n1次单位根之一，n1为离散波面位移数据个数，i是虚数单位；f＝1:n2，为了能够更好的进行傅里叶变换，数据n2选择为距离数据个数n1最近的二
的幂次的数。由于傅里叶级数对应时域幅值，其中已经包含了1/n2项，而傅里叶变换中没有该系数，因此，进行傅里叶变换后，需除以n2/2才能与时域对上，从而得到y(f)：
[0099][0100]
利用摄像机采样时的频率即fps帧速率，计算得到傅里叶变换后的频率fv：
[0101][0102]
其中，v＝0:n2/2。
[0103]
傅里叶变换后的振幅为a＝|yk|，只选取小于采样频率的部分，k＝1:n2/2 1。
[0104]
能量s：
[0105][0106]
谱的零阶矩m0：
[0107][0108]
圆频率ω0：
[0109]
ω0＝2πfv[0110]
谱的二阶矩m2：
[0111][0112]
利用零阶矩求得平均波高
[0113][0114]
利用零阶矩和二阶矩求得平均周期
[0115][0116]
由于本实验拍摄设备位于水槽中间，因此只能获得衰减部分的风浪频谱，若将相机向左移，即可获得完整的风浪频谱，从而获取谱峰周期。获得海浪最大谱值所对应的频率ω＝2πf
p
，并通过求倒数换算出谱峰周期t
p
＝2π/ω＝1/f
p
，f
p
是能量s最大值时的f。
[0117]
本发明介绍的方法对通过实验室小型开敞式波浪水槽造出来的波浪进行了真实数据的采集，计算了有效波高、有效波长、上跨零点周期、有效周期、谱峰周期、波速、波陡和波谱等波浪要素，其数值大小及特征与大水槽造出的波浪及外海波浪非常接近或相似，此方法在实验室中基本试验成功。
[0118]
对于本领域技术人员而言，显然本发明不限于上述示范性实施例的细节，而且在不背离本发明的精神或基本特征的情况下，能够以其他的具体形式实现本发明。因此，无论
从哪一点来看，均应将实施例看作是示范性的，而且是非限制性的，本发明的范围由所附权利要求而不是上述说明限定，因此旨在将落在权利要求的等同要件的含义和范围内的所有变化囊括在本发明内。也不应将权利要求中的任何附图标记视为限制所涉及的权利要求。