基于Open CV的瓦斯爆炸火焰形态演化与传播速度关系的研究
2023-09-11姚礼琳王学强席国军徐景德周振兴
姚礼琳,王学强,席国军,徐景德,周振兴
(1.华北科技学院,北京 东燕郊 065201;2.陕西煤业股份有限公司,陕西 西安 710000; 3.陕西彬长胡家河矿业有限公司,陕西 咸阳 712000)
0 引言
甲烷作为常见的可燃爆炸性气体,在工业生产与生活领域广泛应用,也是矿井瓦斯爆炸的根源之一。瓦斯爆炸带来的巨大损失,严重威胁我国煤矿安全生产。基于气体爆燃的“两波三区”理论,瓦斯爆炸爆燃波为冲击波波阵面传播提供能量。受激波与障碍物影响,流场紊流度急剧增大,火焰波面形态出现不规则形状。火焰形态与速度时空分布表征了燃烧波及爆炸波传播演化规律,因此研究燃烧波形态与速度变化有助于掌握火焰传播过程,掌握瓦斯爆炸机理及规律。
利用高速相机拍摄爆炸图像[1]是研究爆炸流场的重要手段,高速相机配合纹影技术[2-4]能获得清晰的流场结构演变图像。但由于外界因素影响与干扰[5],难以获得清晰直观的图像,进而无法深刻认识流场中火焰的特征结构。传统火焰速度测量方法多为沿爆炸主体壁面布置传感器测点[6,7],通过火焰抵达传感器时刻计算火焰传播速度。但受传感器动态响应时间限制,管道壁面无法密集布置火焰传感器,因此无法捕捉到精细的火焰传播特性。爆炸图像包含火焰速度信息,多位学者对爆炸图像速度信息的提取进行了研究。Park[8]自编程序检测方形空间内爆炸火焰波面位置变化。聂百胜[9]采用图像相关系数法,获得了小尺度激波管内火焰运动速度。余明高[10]使用Matlab对小尺度有机玻璃管中纹影图像进行处理,求得了火焰面轴向运动速度。
但针对中尺度激波管内爆炸火焰图像处理的相关研究较少,因此本文使用Python 3.6 自编程序调用Open CV 3.4.1对中尺度激波管内火焰图像进行图像处理。Open CV是开源的计算机视觉库,能直接使用代码调用不同的边缘检测算子,检测复杂条件下火焰锋面边缘。进而计算出火焰传播速度,通过边缘检测图像分析爆炸火焰形态演化与火焰传播速度变化规律。
1 原始图像采集
1.1 实验系统
实验系统由管道系统、激光纹影系统、点火系统组成。管道系统主体选用6段截面尺寸为200mm×200mm 的方形激波管,如图1所示;总长度为17.05m,管路末端使用盲板密封。K9玻璃观察窗直径为250mm,设置于管道主体9.05m 处。激光纹影系统使用绿激光模块(激光强度5W,532nm)作为纹影光源;使用双凹球反射镜形成Z 字形光路;高速相机选用加拿大Mega Speed 公司相机,采样频率为5000Hz,满足爆炸纹影图像采集要求。点火系统为自行研制,同步装置输入TTL 电平产生高压火花,引燃激波管内预混气体。
图1 实验装置
本实验共计 3 种工况。工况1 全管道充入9.5%体积分数的甲烷-空气预混气(当量比为1)点火能量为255mJ。工况 2 与工况3 以工况1为基础,分别在观察窗处(9.05m)设置1个与2个阻塞率为16.74%的圆柱形障碍物(直径40.34mm、高166mm),如图2所示。
图2 工况示意图
1.2 原始纹影图像分析
图3 为部分原始纹影图像。流场的密度分布与纹影图像中明暗变化存在一定的映射关系[5],如图3所示,密闭管道中激波震荡致使未燃气体密度分布不均,纹影图像中产生散乱分布的黑斑状结构,对火焰形态的精细分析造成干扰。高速相机曝光值通常不能在爆炸实验过程中很好地进行动态调节,纹影图像部分区域亮度过高导致火焰形态显示不完整,如图3(c)所示。爆炸火焰在管道三维空间内运动,高速相机获取的图像是垂直于视窗不同位置火焰波面叠加后的二维图像。因此纹影图像中火焰锋面沿轴向最前端一点的速度大小,能在一定程度上表征火焰轴向速度的快慢。为精确分析火焰锋面形态演变、计算火焰轴向运动速度,需要进行图像处理检测火焰锋面边缘位置信息。
图3 各工况下部分原始纹影图像
2 图像处理
2.1 图像边缘检测原理
二维图像中,物体边缘位于像素灰度值剧烈变化的位置[11]。Open CV能将二维灰度图像转化为数字矩阵,矩阵元素的列数与行数对应灰度图像中像素的横、纵位置坐标,矩阵中元素的数值大小表征对应位置像素灰度值的强度。因此,对图像的数字矩阵做微分运算能检测图像中的边缘信息。梯度幅度较大的位置为图像中物体边缘的位置,梯度幅度的输出值为边缘点的强度数值。
图4 为不同边缘检测算子与相应算子的边缘检测图像。基于梯度幅度与卷积运算,学者们提出了许多边缘检测算子。如图4(b)(c)(d)(e)所示,Open CV内置了Canny、Sobel、Shcarr、Laplacian边缘检测算子。其中,Canny算子先计算像素点灰度值沿x,y方向的一阶导数,将其组合成4个梯度。随后设置梯度值阈值区间,保留区间内的梯度值作为图像中物体边缘的输出值。Sobel算子计算x,y方向像素灰度值的梯度数值,随后叠加两个方向梯度数值得出图像中物体边缘的输出值。Shcarr算子在Sobel算子基础上提升了算法精度,而 Laplacian使用二阶梯度运算来检测物体边缘。
图4 不同边缘检测算子与相应算子的边缘检测图像
由图4(b’)可得,Canny算子过滤了部分无关像素,对比原始纹影图像中的火焰锋面边缘,检测精度最高。如图4(c’)所示,Sobel算子边缘检测图像中,火焰边缘对比原始纹影图像明显外扩,边缘检测结果存在明显误差。Scharr算子保留了大量噪点,如图4(d’)所示,无法对其深入分析。如图4(e’)所示,Laplacian算子为二阶算子,处理后的边缘检测图像暗淡。因此,本文使用Canny算子,进行火焰锋面形态、位置的检测。
以Canny边缘检测图像为基础,继续通过二值化运算提取火焰锋面边缘位置。二值化运算依据预先设定的全局阈值,将图像中低于较小阈值的无关像素灰度值设置为‘0’,阈值区间内的像素灰度值设置为‘255’。多次进行二值化运算,并与原始火焰图像对比确定全局阈值范围。
2.2 图像处理步骤
(1) 火焰边缘检测
如图5(a)所示,首先对纹影图像执行灰度化与ROI操作,截取209×286大小的灰度图像,使用Canny算子进行边缘检测。对比二值化运算结果与原始纹影图像,确定二值化全局阈值的最佳取值,如图5(b)所示;随后使用Photoshop去除二值化运算后的图像噪点,如图5(c)所示。
图5 图像处理过程
(2) 火焰位移计算
使用Python将二值化图像转化为像素矩阵,并对所有列元素求和,将求和结果不为0 的元素数值设置为‘255’。如图5(d)为处理后像素矩阵折线图,横轴为横向像素位置坐标,纵轴为列像素灰度值加和,首个不为零像素点的位置坐标即为火焰锋面像素位置。图像横向像素为‘286’,纹影视窗直径为250mm,因此每个像素实际长度为0.8741mm/pixel。两相邻图像的行矩阵相减得到火焰锋面像素位移数值,乘以像素实际长度换算为火焰锋面位移。
(3) 火焰轴向速度计算
对火焰位移运用式(1)计算平均火焰速度V。
V=ΔL/Δt
(1)
式中,ΔL为火焰锋面位移;Δt为相邻两纹影图像时间间隔。
3 火焰形态演化与速度变化计算
3.1 火焰形态(不加障碍物)演化分析
图6、7为工况 1 原始纹影图像与Canny算子边缘检测图像,图中激波均为向右移动,点火端位于右侧。对比原始纹影图像,Canny算子能降低纹影图像中黑斑结构的影响,清晰的显示火焰锋面的边缘形态。如图6(a)所示,698ms盲板反射的一道激波与火焰前锋同时抵达观察视窗。反射激波与未燃气体冲击火焰锋面,致使火焰锋面变得光滑,如图6(b’)所示。同时降低火焰轴向传播速度,诱导火焰锋面开始发生形变失稳,如图6(c’)所示。
图6 原始纹影图像与Canny算子边缘检测图像(工况1)
如图7(d’)所示,激波诱导火焰锋面发生RM(Richtmyer-Meshkov)不稳定;球状褶皱内部密度场分布不均,表明反射激波穿过火焰锋面致使管内流体密度场与压力场斜交[12],产生涡流并快速湍流化,加速火焰向前运动,如图7(e’)所示。
图7 原始纹影图像与Canny算子边缘检测图像(工况1)
3.2 置障条件下火焰形态演化分析
图8为工况 2火焰纹影图像与Canny算子边缘检测图像,视野中央黑色矩形为圆柱形障碍物。如图8(b’)所示,587.2ms 火焰锋面越过障碍物时发生绕射,速度梯度明显降低;障碍物后方产生孔洞状火焰面褶皱,其规模与数量的不断增加表明障碍物附近流场湍流程度不断提升,如图8(b’)、(c’)所示。
图8 原始纹影图像与Canny算子边缘检测图像(工况2)
图9、10为工况 3火焰纹影图像与Canny算子边缘检测图像,Canny算子削弱了纹影图像中亮度过高对火焰边缘显示的影响,如图9(a’),(b’)所示。707.8ms 此时一道反射激波经过障碍物运动至火焰锋面前方,如图9(a’)所示;由图9(b’),(c’)可知,反射激波致使火焰向来流方向加速回退。火焰完全退出纹影视窗视野后,在首个障碍物表面产生黑色流体涡团,如图9(c)所示。
图9 原始纹影图像与Canny算子边缘检测图像(工况3)
相较火焰锋面首次越过第一个障碍物,如图9(a’)所示,激波冲击后,火焰锋面快速扩张并向前喷射传播,如图10(e’),(f’)所示,表明反射激波驱动未燃气体平动能增强障碍物附近伴随流动场[13]湍流强度,加速火焰锋面向前运动。
图10 原始纹影图像与Canny算子边缘检测图像(工况3)
3.3 火焰速度变化计算
图11为工况1、2、3 纹影视窗内火焰速度图,图中横坐标为相对时间,横轴 0 时刻为火焰进入视窗边缘时刻。在Canny算子边缘检测图像基础上进行数字化处理并结合式(1),对各个工况下火焰速度进行计算,求得不同障碍物数量下火焰轴向传播速度曲线,如图11所示。工况1中,激波与伴流阻滞作用使火焰速度降低了93.3%,火焰运动近乎停滞。随后RM(Richtmyer-Meshkov)不稳定使火焰速度“突跃式”升高,峰值速度提升了13.3%。工况2中,因障碍物遮挡,2~3.6ms,火焰速度为障碍物直径与火焰通过时间的比值,如图11 左侧黄色区域(a)所示。火焰绕行障碍物消耗了部分能量,随后障碍物间湍流使火焰加速向前运动,峰值速度提升了14.2%。工况3中,反射激波冲击致使火焰锋面加速回退,视窗视野中仅可观察到流体涡团。如图11 右侧黄色区域(b)所示,相对时间7.2~11ms火焰回退后重新抵达观察视窗边缘。流体涡团处发育出新的火焰锋面向前喷射传播,峰值速度较初始速度提升了25.9%。通过分析可知,计算所得速度变化规律与爆炸火焰演化规律一致,证明该方法能准确获取爆炸火焰轴向运动速度。
图11 工况1、2、3 纹影视窗内火焰速度图
4 结论
(1) 基于Open CV的边缘检测方法对复杂条件下爆炸火焰边缘的显示效果清晰,火焰边缘位置检测结果准确。
(2) 边缘检测算子对火焰边缘检测效果有重要影响,四种算子中,Canny算子对火焰边缘位置检测最精确。
(3) 计算所得火焰速度时空分布规律与爆炸图像吻合,证明该方法能计算复杂条件下火焰锋面轴向速度时空分布。
(4) 障碍物对火焰运动具有阻滞与加速作用;随着障碍物数量的增加,障碍物附近流场湍流程度提升,显著提升火焰锋面轴向传播速度。
(5) 结合井下爆炸事故现场毁伤程度,障碍物间爆炸流场演化可为事故调查提供判定依据。