CN107702885B - 床面颗粒起悬时间的提取方法 - Google Patents

Abstract

本发明公开了一种床面颗粒起悬时间的提取方法，包括：将高速摄像机连拍的彩色图像做预处理，获得一系列仅凸出显示起悬颗粒轮廓的二值图；以二值图中的床面颗粒起始位置为原点，建立垂面直角坐标系，对二值图进行逐帧批处理，确定各帧二值图上的起悬颗粒在垂面上的位置坐标；根据起悬颗粒在垂面上的位置坐标计算颗粒起悬高度，并逐帧比较颗粒起悬高度与颗粒粒径的大小，直至找到颗粒起悬高度与颗粒粒径之差小于阈值的二值图，将相应二值图的拍摄时间与最后时刻图像的拍摄时间做差，获得床面颗粒起悬一个颗粒高的时间。该方法可以准确计算的床沙颗粒的起动时间，从而准确计算水流挟沙能力S^*，进而合理预测河流演变、解决工程泥沙问题。

Description

床面颗粒起悬时间的提取方法

技术领域

本发明涉及水力学及河流动力学领域，尤其涉及一种床面颗粒起悬时间的提取方法。

背景技术

在天然河流中运动着的泥沙，存在两种运动形式，一是随水流悬浮运动的悬移质，二是在床面滚动、滑动、跳跃的推移质。推移质和悬移质可以互相转换，当水流强度足够大，床面上的颗粒可被水流冲起，以悬移状态运动；若水流强度很小，含沙水流中的悬沙颗粒会在重力作用下沉积于河床表面。

天然河流是水、沙输运的通道，河流中泥沙的运移，会改变河道形态，引起河道的冲淤演变。为预测河道冲淤演变，需要计算水流的挟沙能力S^*，所谓挟沙能力，是指单位水体能够输运的泥沙量，其单位是kg/m³(千克/立方米)，在水沙一维数学模型中，它是一个断面平均值，它与断面平均含沙量S(单位水体中所含的泥沙量)的量纲相同，但含义不同，挟沙能力S^*是期望值，含沙量S则是实测值。当水流的实际含沙量S大于水流挟沙能力S^*时，水流中的部分泥沙就要淤积在河床；当水流的实际含沙量小于S^*时，河床就要发生冲刷。

对挟沙能力S^*的估算是河流泥沙运动力学的重要问题。对于粒径为D_l的均匀沙，基于泥沙运动统计理论的挟沙能力计算公式如下：

上式中，S^*(l)表示粒径为D_l的均匀沙的挟沙能力，l表示非均匀沙的粒径分组号。是床面上第l组颗粒由静止起动到脱离床面一个粒径高D_l的平均时间，目前它是上式中最难确定的参数，式中其他各量的物理意义和确定方法相对简单，可参考韩其为院士专著《水库淤积》；同时，韩其为院士也给出数学推导得到的参数t₀₄的方式，但是形式较为复杂，属于理论公式，尚待实测数据验证。

目前，主要有如下两个方案利用高速摄影技术研究颗粒运动的成果。

方案1：惠遇甲、胡春宏利用高速摄影技术研究了水流中颗粒跃移的运动学特征。

在他们的研究中，采用不同容重的颗粒材料和不同的水流强度，利用高速摄影技术系统测量了水流中颗粒的跃移过程。他们的工作是在1980年代末开展的，研究过程采用的高速摄影仪器是当时非常先进的Pentazet-35和GSJ型高速摄影机，摄影频率为60～300f/s。基于此设备记录的1000多条单颗粒泥沙在水槽中的运动轨迹，他们提炼了泥沙颗粒的常规运动形式，并统计了不同运动形式出现的概率与水流强度之间的关系。这是利用高速摄影技术研究颗粒运动的经典成果。在他们的研究中，并未涉及对床面颗粒起动时间t₀₄的推求。

方案2：苗蔚、陈启刚、李丹勋等人通过高速摄影技术研究了推移质泥沙的起动概率。

苗蔚等人对泥沙起动概率的研究，是通过在试验水槽正上方布设俯视高速摄影设备获取床沙影像的，所获得的床沙影像，是床面一个平面区域内的所有床沙分布图，为了增加图像的清晰度和对比度，避免背景照明不均匀、杂质和气泡等对识别结果产生影响，她们对原始图像进行了顶帽变换、对比度线性拉伸和高斯滤波处理，见图1，其中左侧为原始图，右侧为处理后的图。

为了在从连续两张图片中获取泥沙起动的信息，他们采用了Keshavarzy,A等人提出的方法，将两帧图像相减，根据图像灰度变化判读泥沙是否发生运动，这种处理方式的理论依据在于，颗粒运动导致初始位置和当前位置的灰度均发生改变。在泥沙起动试验中，由于运动泥沙暴露度大、距光源近，其对应的灰度值一般比静止泥沙大，当泥沙发生运动后，其初始位置的灰度变小而当前位置的灰度增大，当利用第二帧图片减去第一帧图片时，初始位置的灰度为负，当前位置的灰度为正。分析本文实际拍摄的图片时发现，有部分运动颗粒比静止颗粒暗，即当第二帧减第一帧时，初始位置的灰度为正，当前位置的灰度为负，因此统计全部的负值区域或者正值区域将得到运动颗粒总数，两种方法得到的结果误差率为5％。通过这种处理方式，可以得到不同时刻两帧图像的灰度差分图(见图2)，图2中成对出现的正负峰值清晰标识了运动泥沙的初始位置和当前位置。

显然，上述两个方案关注的重点是颗粒的运动形式和床面推移质的运动过程，均未涉及床沙颗粒的起动时间t₀₄。

发明内容

本发明的目的是提供一种床面颗粒起悬时间的提取方法，可以准确计算床沙颗粒的起动时间。

本发明的目的是通过以下技术方案实现的：

一种床面颗粒起悬时间的提取方法，包括：

将高速摄像机连拍的彩色图像做预处理，获得一系列仅凸出显示起悬床面颗粒轮廓的二值图；

以二值图中的起悬床面颗粒起始位置为原点，建立垂面直角坐标系，对二值图进行逐帧批处理，确定各帧二值图上的起悬床面颗粒在垂面上的位置坐标；

根据起悬床面颗粒在垂面上的位置坐标计算起悬床面颗粒起悬高度，并逐帧比较起悬床面颗粒起悬高度与起悬床面颗粒粒径的大小，直至找到起悬床面颗粒起悬高度与起悬床面颗粒粒径之差小于预先规定的精度阈值的二值图，将相应二值图的拍摄时间与颗粒静止于床面最后时刻图像的拍摄时间做差，获得起悬床面颗粒起悬时间。

由上述本发明提供的技术方案可以看出，可以基于高速摄影连拍图像准确的计算床沙起悬时间，计算过程均可编程实现，通过程序高效运行，可实现对高速摄影海量连拍图片的高效处理。同时，根据准确计算的床沙颗粒的起动时间，准确计算水流挟沙能力S^*，从而合理预测河流演变、解决工程泥沙问题。

附图说明

为了更清楚地说明本发明实施例的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域的普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他附图。

图1为本发明背景技术提供的高速摄影俯视床沙分布图；

图2为本发明背景技术提供的连续两帧图像的灰度差示意图；

图3为本发明实施例提供的一种床面颗粒起悬时间的提取方法的流程图；

图4为本发明实施例提供的灰度图二值化的示意图；

图5为本发明实施例提供的床面颗粒起动过程在平面上的不同位置的示意图；

图6为本发明实施例提供的逐行扫描获取颗粒中心点的C++代码；

图7为本发明实施例提供的获取颗粒脱离床面一个粒径高所用时间t₀₄的C++代码；

图8为本发明实施例提供的叠加连拍二值图后获取的颗粒垂向运动轨迹图；

图9为本发明实施例提供的二值图叠加的C++代码；

图10为本发明实施例提供的床面颗粒垂向运动连拍二值图；

图11为本发明实施例提供的叠加连拍二值图后获取的颗粒垂向运动轨迹图。

具体实施方式

下面结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明的保护范围。

某粒径床面泥沙颗粒的平均起悬时间t₀₄，是泥沙运动统计理论水流挟沙能力S^*的计算公式的重要参数，计算床面泥沙颗粒由静止起悬t₀₄，进而确定床面泥沙颗粒由静止起悬的平均时间t₀₄，是准确计算水流挟沙能力S^*的必要条件，准确计算水流挟沙能力S^*，是合理预测河流演变、解决工程泥沙问题的关键所在。为此，本发明实施例提供一种床面颗粒起悬时间的提取方法，如图3所示，其主要包括如下步骤：

步骤1、将高速摄像机连拍的彩色图像做预处理，获得一系列仅凸出显示起悬床面颗粒轮廓的二值图。

通过本步骤的处理，可消除连拍图像中的不均匀背景和噪声信号，具体如下：

1)由于高速摄像机获取的彩色图片，存储的是RGB信息，需要将高速摄像机连拍的彩色图像批量转换为灰度图，转换公式如下：

Y＝α₁R+α₂G+α₃B；

其中，Y为转换后的灰度图的灰度值，R、G、B分别对应为彩色图像中红色值、绿色值、蓝色值；α₁、α₂、α₃分别对应为R、G、B的权重，且α₁+α₂+α₃＝1；示例性的，α₁、α₂、α₃可分别设为0.299、0.587、0.114。

通过编制计算机程序进行批处理，可以对高速摄影的连拍图像由彩色图转为灰度图进行逐帧处理。

2)对灰度图做逐帧自适应阈值二值化批处理，获得一系列仅凸出显示起悬床面颗粒轮廓的二值图。

图像二值化，是将有256个亮度等级的灰度图像转换为仅包含0和255两个灰度值图像的操作，通过这一处理可减小图像数据量，获得仅凸出显示起悬床面颗粒轮廓的二值图。

图像二值化，一般是在给定一个阈值后，逐像素点比较各点灰度值，大于或等于阈值的像素被判定为属于特定物体，设置其灰度值为255，否则这些像素点被排除在物体区域以外，设置其灰度值为0。如图4所示，为灰度图二值化的示例，左侧为灰度图，右侧为对应的二值图。

对于阈值的选取，有全局阈值和自适应阈值两种方式，全局阈值适用于背景均匀的情况，自适应阈值分割则适用于图像背景不均匀的情况。由于颗粒起动过程水流流态会影响背景光线，多会形成不均匀背景，因此本发明采用自适应阈值(AdaptiveThreshold)算法进行连拍灰度图的二值化分割。

同样的，上述二值化过程也可以通过编程进行批处理。

步骤2、以二值图中的起悬床面颗粒起始位置为原点，建立垂面直角坐标系，对二值图进行逐帧批处理，确定各帧二值图上的起悬床面颗粒在垂面上的位置坐标。

本发明实施例中，对二值图进行逐帧批处理，主要是为了获取每帧二值图的颗粒中心点，其过程如下：对二值图进行高斯滤波，经高斯滤波后床面颗粒中心像素点的灰度值将变成局部灰度最大值点；对滤波后的二值图进行逐点扫描检查，寻找局部灰度最大值点，从而确定起悬床面颗粒中心点。之后，就可以根据起悬床面颗粒中心点计算在垂面上的位置坐标。

如图5所示，为起悬床面颗粒起动过程在平面上的不同位置的示例。

假设每幅二值图均由m行×n列个像素点组成，设每个像素点的灰度值为G_i,j，下标中的i，j值分别是像素点所在行号和列号，则进行逐行扫描获取颗粒中心点的C++代码如图6所示。

第i帧二值图的起悬床面颗粒中心点坐标记为(x_i,y_i)，图6所示代码中有两个参数需要人工测取设置。一是实际长度与图中像素个数的转换关系，即1毫米等于二值图中多少个像素直径之和，设此值为Coef，则Coef值需要人工测取。二是颗粒粒径值D，这也需要人工测取。

步骤13、根据起悬床面颗粒在垂面上的位置坐标计算起悬床面颗粒起悬高度，并逐帧比较起悬床面颗粒起悬高度与起悬床面颗粒粒径的大小，直至找到起悬床面颗粒起悬高度与起悬床面颗粒粒径之差小于预先规定的精度阈值eps的二值图，将相应二值图的拍摄时间与颗粒静止于床面最后时刻图像的拍摄时间做差，获得起悬床面颗粒起悬时间。

本步骤具体过程如下：

1)计算起悬床面颗粒起悬高度，公式为：

Y_i＝y_i-y₀；

其中，Y_i为第i帧二值图的起悬床面颗粒起悬高度(脱离床面的高度)，y_i第i帧二值图中起悬床面颗粒的y轴坐标，y₀为起悬床面颗粒静止于床面时，在二值图中起始位置处的y轴坐标。

2)比较Y_i与预先测量的起悬床面颗粒粒径D的大小，直至找到|Y_i-D|＜eps的二值图，其中，eps为预先规定的精度阈值；

3)将相应二值图拍摄时间T₄与颗粒在床面上静止的最后时刻图像的拍摄时间T₀做差，即得起悬床面颗粒起悬一个粒径高的起悬时间：t₀₄＝T₄-T₀。

假设连拍图像共有m＝100帧，且第1帧恰为颗粒在床面上静止的最后时刻所拍摄，则获取颗粒脱离床面一个粒径D高的时间t₀₄的批处理的代码图7所示。

另外，在执行本发明实施例步骤1之后，可以获得一系列的二值图，而叠加一系列的二值图，可以获取床面颗粒起动垂向轨迹。具体的过程如下：

1)以循环语句逐帧加载并叠加二值图，将每幅二值图中处于同一(x,y)位置的像素点的灰度值(G_x,y)_k相加，表示为：

其中，m为二值图的总数，为叠加后的灰度值；

2)再将每个像素点灰度叠加值大于0的点的灰度值均设置为255，表示为：

从而得到叠加后的床面颗粒垂向运动轨迹图，如图8所示，为叠加二值图后获取的颗粒垂向运动轨迹图。

二值图叠加的C++代码如9所示。

另一方面，在发明人开展的国家国际科技合作专项“复杂条件下的河流演变模拟关键技术联合研究”项目和中国水利水电科学研究院科研专项“非均匀沙运动理论前沿研究”中，采用高速摄影技术对床面泥沙起动过程进行了高速摄影测验，对床面泥沙起动过程高速连拍图按本发明所述方法进行处理后，分别得到了如图10所示的连拍二值图和如图11所示的叠加的运动轨迹。

以本发明方法提取如图10和图11示的颗粒起动过程的起动时间，得到颗粒的起动时间t₀₄介于1e-2～1e-4秒之间。

通过以上的实施方式的描述，本领域的技术人员可以清楚地了解到上述实施例可以通过软件实现，也可以借助软件加必要的通用硬件平台的方式来实现。基于这样的理解，上述实施例的技术方案可以以软件产品的形式体现出来，该软件产品可以存储在一个非易失性存储介质(可以是CD-ROM，U盘，移动硬盘等)中，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行本发明各个实施例所述的方法。

以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明披露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求书的保护范围为准。

Claims (6)

1.一种起悬床面颗粒起悬时间的提取方法，其特征在于，包括：

将高速摄像机连拍的彩色图像做预处理，获得一系列仅凸出显示起悬床面颗粒轮廓的二值图；

以二值图中的起悬床面颗粒起始位置为原点，建立垂面直角坐标系，对二值图进行逐帧批处理，确定各帧二值图上的起悬床面颗粒在垂面上的位置坐标；

根据起悬床面颗粒在垂面上的位置坐标计算起悬床面颗粒起悬高度，并逐帧比较起悬床面颗粒起悬高度与起悬床面颗粒粒径的大小，直至找到起悬床面颗粒起悬高度与起悬床面颗粒粒径之差小于预先规定的精度阈值的二值图，将相应二值图的拍摄时间与颗粒静止于床面最后时刻图像的拍摄时间做差，获得起悬床面颗粒起悬时间。

2.根据权利要求1所述的一种起悬床面颗粒起悬时间的提取方法，其特征在于，所述将高速摄像机连拍的彩色图像做预处理，获得一系列仅凸出显示起悬床面颗粒轮廓的二值图包括：

将高速摄像机连拍的彩色图像批量转换为灰度图，转换公式如下：

Y＝α₁R+α₂G+α₃B；

其中，Y为转换后的灰度图的灰度值，R、G、B分别对应为彩色图像中红色值、绿色值、蓝色值；α₁、α₂、α₃分别对应为R、G、B的权重，且α₁+α₂+α₃＝1；

再对灰度图做逐帧自适应阈值二值化批处理，获得一系列仅凸出显示起悬床面颗粒轮廓的二值图。

3.根据权利要求1所述的一种起悬床面颗粒起悬时间的提取方法，其特征在于，所述确定各帧二值图上的起悬床面颗粒在垂面上的位置坐标包括：

对二值图进行高斯滤波，经高斯滤波后起悬床面颗粒中心像素点的灰度值将变成局部灰度最大值点；

对滤波后的二值图进行逐点扫描检查，寻找局部灰度最大值点，从而确定起悬床面颗粒中心点；

根据起悬床面颗粒中心点计算在垂面上的位置坐标。

4.根据权利要求1所述的一种起悬床面颗粒起悬时间的提取方法，其特征在于，所述获得起悬床面颗粒起悬时间包括：

计算起悬床面颗粒起悬高度，公式为：

Y_i＝y_i-y₀；

其中，Y_i为第i帧二值图的起悬床面颗粒起悬高度，y_i第i帧二值图中起悬床面颗粒的y轴坐标，y₀为起悬床面颗粒静止于床面时，在二值图中起始位置处的y轴坐标；

比较Y_i与预先测量的起悬床面颗粒粒径D的大小，直至找到|Y_i-D|＜eps的二值图，其中，eps为预先规定的精度阈值；

将相应二值图拍摄时间T₄与颗粒在床面上静止的最后时刻图像的拍摄时间T₀做差，即得起悬床面颗粒起悬一个粒径高的起悬时间：t₀₄＝T₄-T₀。

5.根据权利要求1-4任一项所述的一种起悬床面颗粒起悬时间的提取方法，其特征在于，该方法还包括：叠加一系列的二值图，获取起悬床面颗粒起动垂向轨迹。

6.根据权利要求5所述的一种起悬床面颗粒起悬时间的提取方法，其特征在于，所述叠加一系列的二值图，获取起悬床面颗粒起动垂向轨迹包括：

以循环语句逐帧加载并叠加二值图，将每幅二值图中处于同一(x,y)位置的像素点的灰度值(G_x,y)_k相加，表示为：

其中，m为二值图的总数，为叠加后的灰度值；

再将每个像素点灰度叠加值大于0的点的灰度值均设置为255，表示为：

从而得到叠加后的起悬床面颗粒垂向运动轨迹图。