CN117635874A - 基于图像多特征提取与融合算法的破片场三维快速反演系统 - Google Patents
基于图像多特征提取与融合算法的破片场三维快速反演系统 Download PDFInfo
- Publication number
- CN117635874A CN117635874A CN202311687287.5A CN202311687287A CN117635874A CN 117635874 A CN117635874 A CN 117635874A CN 202311687287 A CN202311687287 A CN 202311687287A CN 117635874 A CN117635874 A CN 117635874A
- Authority
- CN
- China
- Prior art keywords
- fragment
- perforation
- dimensional
- image
- target plate
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 239000012634 fragment Substances 0.000 title claims abstract description 254
- 238000000605 extraction Methods 0.000 title claims abstract description 63
- 238000004422 calculation algorithm Methods 0.000 title claims abstract description 44
- 230000004927 fusion Effects 0.000 title claims abstract description 25
- 238000004364 calculation method Methods 0.000 claims abstract description 50
- 238000000034 method Methods 0.000 claims abstract description 49
- 238000004880 explosion Methods 0.000 claims abstract description 47
- 238000005259 measurement Methods 0.000 claims abstract description 46
- 238000009826 distribution Methods 0.000 claims abstract description 18
- 238000004458 analytical method Methods 0.000 claims abstract description 9
- 239000011159 matrix material Substances 0.000 claims description 49
- 238000012937 correction Methods 0.000 claims description 36
- 230000009466 transformation Effects 0.000 claims description 30
- 238000003384 imaging method Methods 0.000 claims description 19
- 238000012545 processing Methods 0.000 claims description 19
- 238000001914 filtration Methods 0.000 claims description 16
- 230000011218 segmentation Effects 0.000 claims description 12
- 230000000877 morphologic effect Effects 0.000 claims description 11
- 238000006243 chemical reaction Methods 0.000 claims description 8
- 238000003708 edge detection Methods 0.000 claims description 8
- 229910052704 radon Inorganic materials 0.000 claims description 8
- SYUHGPGVQRZVTB-UHFFFAOYSA-N radon atom Chemical compound [Rn] SYUHGPGVQRZVTB-UHFFFAOYSA-N 0.000 claims description 8
- 238000005286 illumination Methods 0.000 claims description 7
- 238000005457 optimization Methods 0.000 claims description 7
- 238000001514 detection method Methods 0.000 claims description 6
- 238000007781 pre-processing Methods 0.000 claims description 5
- 230000009172 bursting Effects 0.000 claims description 3
- 230000008030 elimination Effects 0.000 claims description 3
- 238000003379 elimination reaction Methods 0.000 claims description 3
- 230000002708 enhancing effect Effects 0.000 claims description 3
- 230000001131 transforming effect Effects 0.000 claims description 3
- 238000012360 testing method Methods 0.000 abstract description 18
- 238000011156 evaluation Methods 0.000 abstract description 4
- 230000006870 function Effects 0.000 description 14
- 238000010586 diagram Methods 0.000 description 12
- 230000008569 process Effects 0.000 description 8
- 230000007797 corrosion Effects 0.000 description 7
- 238000005260 corrosion Methods 0.000 description 7
- 238000005516 engineering process Methods 0.000 description 6
- 230000001629 suppression Effects 0.000 description 5
- 101100446506 Mus musculus Fgf3 gene Proteins 0.000 description 4
- 230000000694 effects Effects 0.000 description 4
- 230000002147 killing effect Effects 0.000 description 4
- 238000013519 translation Methods 0.000 description 4
- 101000767160 Saccharomyces cerevisiae (strain ATCC 204508 / S288c) Intracellular protein transport protein USO1 Proteins 0.000 description 3
- 230000009286 beneficial effect Effects 0.000 description 3
- 239000002131 composite material Substances 0.000 description 3
- 238000005192 partition Methods 0.000 description 3
- 238000011160 research Methods 0.000 description 3
- 230000003068 static effect Effects 0.000 description 3
- 235000002566 Capsicum Nutrition 0.000 description 2
- 239000006002 Pepper Substances 0.000 description 2
- 241000722363 Piper Species 0.000 description 2
- 235000016761 Piper aduncum Nutrition 0.000 description 2
- 235000017804 Piper guineense Nutrition 0.000 description 2
- 235000008184 Piper nigrum Nutrition 0.000 description 2
- 230000006978 adaptation Effects 0.000 description 2
- 238000004891 communication Methods 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 239000006185 dispersion Substances 0.000 description 2
- 239000004744 fabric Substances 0.000 description 2
- 238000007499 fusion processing Methods 0.000 description 2
- 230000005484 gravity Effects 0.000 description 2
- 230000007246 mechanism Effects 0.000 description 2
- 230000003287 optical effect Effects 0.000 description 2
- 238000003672 processing method Methods 0.000 description 2
- 150000003839 salts Chemical class 0.000 description 2
- 238000003860 storage Methods 0.000 description 2
- 230000001360 synchronised effect Effects 0.000 description 2
- 238000012795 verification Methods 0.000 description 2
- 230000000007 visual effect Effects 0.000 description 2
- 241000270295 Serpentes Species 0.000 description 1
- 229910000831 Steel Inorganic materials 0.000 description 1
- 230000009471 action Effects 0.000 description 1
- 230000003044 adaptive effect Effects 0.000 description 1
- 230000003321 amplification Effects 0.000 description 1
- 238000005422 blasting Methods 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 150000001875 compounds Chemical class 0.000 description 1
- 230000006835 compression Effects 0.000 description 1
- 238000007906 compression Methods 0.000 description 1
- 238000004590 computer program Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000012942 design verification Methods 0.000 description 1
- 238000005474 detonation Methods 0.000 description 1
- 229910003460 diamond Inorganic materials 0.000 description 1
- 239000010432 diamond Substances 0.000 description 1
- 230000010339 dilation Effects 0.000 description 1
- 230000003628 erosive effect Effects 0.000 description 1
- 230000000670 limiting effect Effects 0.000 description 1
- 230000014759 maintenance of location Effects 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000013178 mathematical model Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000010606 normalization Methods 0.000 description 1
- 238000003199 nucleic acid amplification method Methods 0.000 description 1
- 230000035515 penetration Effects 0.000 description 1
- 230000002829 reductive effect Effects 0.000 description 1
- 238000007493 shaping process Methods 0.000 description 1
- 230000008054 signal transmission Effects 0.000 description 1
- 239000010959 steel Substances 0.000 description 1
- 238000010998 test method Methods 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
Landscapes
- Image Analysis (AREA)
Abstract
本发明适用于弹药测试技术领域,提供了基于图像多特征提取与融合算法的破片场三维快速反演系统,所述系统包括:破片穿孔位置和面积计算模块,用于获取战斗部爆炸后的等效靶板图像;破片穿孔位置三维反演模块,用于建立破片穿孔三维坐标系;破片飞散参数计算模块,用于根据所述破片穿孔三维坐标系并参照国标计算得到破片飞散角参数、破片飞散方向角参数和破片分布密度参数;系统测量模块,用于根据破片穿孔识别率、破片穿孔面积以及破片穿孔位置三个指标参数进行精度计算。本发明建立的模型与提出方法将使破片战斗部试验结果的分析、评估在效率、精度上得到显著提高。
Description
技术领域
本发明涉及弹药测试技术领域,具体是涉及基于图像多特征提取与融合算法的破片场三维快速反演系统。
背景技术
现代战争对战斗部的杀伤威力和毁伤效应提出了更高要求,世界各国都在积极研制与发展各种毁伤机理的新型战斗部技术,以提高弹药和导弹的高效毁伤能力。战斗部是导弹武器系统的重要组成部分和完成战斗使命的主要部件,根据形成毁伤的特点及毁伤机理,可分为杀伤战斗部、爆破战斗部、侵彻战斗部、聚能战斗部等。杀伤战斗部主要是破片毁伤,破片的穿孔数量、速度、飞散角、飞散方向角、密度分布以及破片穿孔面积等参数是判定杀伤战斗部作用威力的重要指标。而在杀伤战斗部设计研发中,破片的参数一般采用靶场静态爆轰试验来测定与考核,从而为战斗部的方案评估、设计校验及定型验收环节提供直接依据。
目前杀伤战斗部破片飞散特性的测量,主要利用人工数着靶点、划刻线等方式进行统计分析与计算,随着测试量的增大,人工方法测试参数误差大、效率低。本专利采用数字近景摄影测量技术及图像处理方法,解决爆炸破片飞散特性快速准确获取与破片场的三维反演问题,以期弥补人工方法的不足,从而对战斗部威力全面、准确、快速评估提供测试方法和技术。
发明内容
针对现有技术存在的不足,本发明的目的在于提供基于图像多特征提取与融合算法的破片场三维快速反演系统,以解决上述背景技术中存在的问题。
本发明是这样实现的,基于图像多特征提取与融合算法的破片场三维快速反演系统,所述系统包括:
破片穿孔位置和面积计算模块,用于获取战斗部爆炸后的等效靶板图像,并通过几何畸变校正、靶板区域提取以及图像处理得到破片穿孔中心定位参数与破片穿孔面积参数;
破片穿孔位置三维反演模块,用于根据所述破片穿孔中心定位参数建立破片穿孔二维坐标系并将其转化为破片穿孔三维坐标系,所述破片穿孔三维坐标系为破片穿孔相对于爆心建立;
破片飞散参数计算模块,用于根据所述破片穿孔三维坐标系并参照国标计算得到破片飞散角参数、破片飞散方向角参数和破片分布密度参数;
系统测量模块,用于根据破片穿孔识别率、破片穿孔面积以及破片穿孔位置三个指标参数进行精度计算。
作为本发明进一步的方案:所述破片穿孔位置和面积计算模块包括:
几何畸变校正单元,用于根据畸变校正矩阵H对所述等效靶板图像进行校正;
破片穿孔分割提取单元,用于根据破片穿孔与非穿孔部分灰度值的差别,并结合边缘和纹理特征对等效靶板图像进行提取、识别得到破片穿孔图像,所述破片穿孔图像表现形式为等效靶板图像上的破片穿孔位置具有标注;
破片位置与面积计算单元,用于根据破片穿孔图像对每个破片穿孔位置进行提取和统计,并基于灰度边缘轮廓与灰度加权计算得到破片穿孔中心定位参数与破片穿孔面积参数。
作为本发明进一步的方案:所述几何畸变校正单元包括:
靶板轮廓提取单元,用于根据Otsu算法按图像的灰度特性将背景和靶板区分,靶板变为黑色,分割之后的结果是一幅二值图,再基于形态学边界提取的方法得到为白色四边形的靶板边界线;
边界线解析单元,用于根据Radon变换检测直线的方法得到边界线方程,并根据边界线方程获取顶点坐标参数,所述顶点坐标参数为靶板边界线四个顶点的位置左边;
畸变校正坐标变换单元,用于根据所述畸变校正矩阵H将顶点坐标参数进行变换得到校正坐标参数。
作为本发明进一步的方案:所述破片穿孔分割提取单元包括:
图像预处理单元,用于根据自适应中值滤波方法对等效靶板图像进行去噪处理,用于保留等效靶板图像的边缘细节信息;
对比度增强单元,用于采用灰度变换的方法将等效靶板图像灰度等级扩大;
目标增强单元,用于根据Canny边缘检测的目标增强算法将等效靶板图像中的破片穿孔目标的边缘进行增强,所述增强后的破片穿孔目标的边缘线为零散状;
破片穿孔边缘优化处理单元,用于根据主动轮廓法对等效靶板图像进行优化处理;
特征提取单元,用于基于纹理特征对等效靶板图像的破片穿孔区域进行特征提取;
干扰清除单元,用于基于数学形态学对等效靶板图像进行干扰目标滤除;
目标识别单元,用于根据所述边缘特征和纹理特征对破片穿孔进行识别,分离出破片穿孔位置并进行标注得到破片穿孔图像。
作为本发明进一步的方案:所述破片穿孔位置三维反演模块包括:
坐标系变换单元,用于根据所述破片穿孔中心定位参数建立破片穿孔二维坐标系,并基于旋转矩阵R将穿孔二位坐标系并将其转化为破片穿孔三维坐标系;
三维反演单元,用于基于变换矩阵M将不同的等效靶板图像的破片穿孔二维坐标系反演到统一的已战斗部爆心为基准的破片穿孔三维坐标系。
作为本发明进一步的方案:所述破片飞散参数计算模块包括:
破片分布密度计算单元,用于基于破片穿孔三维坐标系并参照球面投影的方法计算得到破片分布密度参数;
联合计算单元,用于对等效靶板图像进行虚拟网格划分并进行计算得到破片飞散角参数和破片飞散方向角参数。
作为本发明进一步的方案:所述系统测量模块包括:
识别率核算单元,用于根据相机成像分辨率和光照角度两个方面对破片穿孔识别率进行核算,使得破片穿孔识别率大于98%;
面积核算单元,用于根据破片穿孔的像元分辨率对穿孔面积进行核算,使得破片穿孔面积大于4.2%;
位置核算单元,用于根据战斗部爆炸后破片穿孔和单点标志的二维位置和战斗部爆炸前单点标志的三维反演坐标三个指标计算得到破片穿孔位置测量精度。
作为本发明进一步的方案:所述位置核算单元包括:
二维位置精度计算单元,用于根据破片穿孔和单点标志的像元位置的提取精度、靶板的长、宽尺寸测量的精度确定战斗部爆炸后破片穿孔和单点标志的二维位置的定位误差;
三维反演精度计算单元,用于根据三维反演原理确定战斗部爆炸前单点标志三维坐标的位置测量精度。
与现有技术相比,本发明的有益效果是:
本发明提出爆炸试验测试过程中的靶板几何畸变校正方法,破片穿孔和凹坑区域精准分割和提取方法,提出了采用灰度加权重心法进行了破片穿孔二维坐标的定位,采用光速平差方法破片穿孔和凹坑位置实现了破片穿孔三维位置坐标的测量,结合爆炸场景反演结果,实现破片场的三维反演,同时开展了系统测量精度的分析研究,形成数据处理平台系统,自动输出破片飞散特性参数,并可以利用历史外场试验结果进行了试验验证,本发明建立的模型与提出方法将使破片战斗部试验结果的分析、评估在效率、精度上得到显著提高。
附图说明
图1为基于图像多特征提取与融合算法的破片场三维快速反演系统的结构示意图。
图2为基于图像多特征提取与融合算法的破片场三维快速反演系统中破片穿孔位置和面积计算模块的结构示意图。
图3为基于图像多特征提取与融合算法的破片场三维快速反演系统中几何畸变校正单元的结构示意图。
图4为基于图像多特征提取与融合算法的破片场三维快速反演系统中破片穿孔分割提取单元的结构示意图。
图5为基于图像多特征提取与融合算法的破片场三维快速反演系统中破片穿孔位置三维反演模块的结构示意图。
图6为基于图像多特征提取与融合算法的破片场三维快速反演系统中破片飞散参数计算模块的结构示意图。
图7为基于图像多特征提取与融合算法的破片场三维快速反演系统中系统测量模块的结构示意图。
图8为基于图像多特征提取与融合算法的破片场三维快速反演系统中位置核算单元的结构示意图。
图9为图像变形的水平倾斜的示意图。
图10为图像变形的垂直倾斜的示意图。
图11为图像变形的水平垂直倾斜的示意图。
图12为Radon变换定义的示意图。
图13为二值图像的结构元素的示意图。
图14为爆炸场景坐标系的示意图。
图15为不同坐标系间的坐标转换的示意图。
图16为爆心坐标系下某靶板上破片穿孔位置的三维反演图。
图17为等效靶板球面投影示意图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清晰,以下结合附图及具体实施例,对本发明进行进一步详细说明。应当理解,此处所描述的具体实施例仅仅用以解释本发明,并不用于限定本发明。
以下结合具体实施例对本发明的具体实现进行详细描述。
如图1所示,本发明实施例提供了基于图像多特征提取与融合算法的破片场三维快速反演系统,所述系统包括:
破片穿孔位置和面积计算模块100,用于获取战斗部爆炸后的等效靶板图像,并通过几何畸变校正、靶板区域提取以及图像处理得到破片穿孔中心定位参数与破片穿孔面积参数;
破片穿孔位置三维反演模块200,用于根据所述破片穿孔中心定位参数建立破片穿孔二维坐标系并将其转化为破片穿孔三维坐标系,所述破片穿孔三维坐标系为破片穿孔相对于爆心建立;
破片飞散参数计算模块300,用于根据所述破片穿孔三维坐标系并参照国标计算得到破片飞散角参数、破片飞散方向角参数和破片分布密度参数;
系统测量模块400,用于根据破片穿孔识别率、破片穿孔面积以及破片穿孔位置三个指标参数进行精度计算。
需要说明的是,目前杀伤战斗部破片飞散特性的测量,主要利用人工数着靶点、划刻线等方式进行统计分析与计算,随着测试量的增大,人工方法测试参数误差大、效率低。本专利采用数字近景摄影测量技术及图像处理方法,解决爆炸破片飞散特性快速准确获取与破片场的三维反演问题。
本发明实施例中,本发明提出爆炸试验测试过程中的靶板几何畸变校正方法,破片穿孔和凹坑区域精准分割和提取方法,提出了采用灰度加权重心法进行了破片穿孔二维坐标的定位,采用光速平差方法破片穿孔和凹坑位置实现了破片穿孔三维位置坐标的测量,结合爆炸场景反演结果,实现破片场的三维反演,同时开展了系统测量精度的分析研究,形成数据处理平台系统,自动输出破片飞散特性参数,并可以利用历史外场试验结果进行了试验验证,此专利建立的模型与提出方法将使破片战斗部试验结果的分析、评估在效率、精度上得到显著提高。
如图2所示,作为本发明一个优选的实施例,所述破片穿孔位置和面积计算模块100包括:
几何畸变校正单元101,用于根据畸变校正矩阵H对所述等效靶板图像进行校正;
破片穿孔分割提取单元102,用于根据破片穿孔与非穿孔部分灰度值的差别,并结合边缘和纹理特征对等效靶板图像进行提取、识别得到破片穿孔图像,所述破片穿孔图像表现形式为等效靶板图像上的破片穿孔位置具有标注;
破片位置与面积计算单元103,用于根据破片穿孔图像对每个破片穿孔位置进行提取和统计,并基于灰度边缘轮廓与灰度加权计算得到破片穿孔中心定位参数与破片穿孔面积参数。
本发明实施例中,视觉测量系统是利用计算机来代替人眼对图像进行理解和识别。一方面由于拍摄角度和拍摄距离的原因,特别是在拍摄长方形或正方形物体时,拍摄的图像存在一定的倾斜形变,竖线变成倾斜,正方形变成菱形;另一方面因摄像机镜头不对称或者摄像机的组装存在偏差,引起的拍摄图像的畸变。但是,计算机不能准确的从图像中识别出物体的形状信息,需要对图像进行畸变校正,再通过计算机识别出图像中的物体形状,就能更真实的反应物体的实际尺寸,针对上述已校正的靶板对其上的破片穿孔目标进行分割与提取。根据破片穿孔与非穿孔部分灰度值存在较大的差别,利用最大类间方差确定分割破片穿孔区域的阈值,结合边缘和纹理特征进行提取、识别。
基于连通区域对获取的每个破片穿孔区域进行提取和统计,并基于灰度边缘轮廓与灰度加权计算其中心坐标位置和面积。
破片穿孔位置定位:对于光亮强度不均匀的破片穿孔区域,采用破片穿孔灰度分布计算灰度权重质心坐标认为是破片穿孔位置定位点。设破片穿孔区域为m*n大小的图像f,灰度值超过阈值T的像素均参与重心的计算,破片穿孔坐标为:
其中,
xi表示像素第i行的坐标,yi表示像素第j列的坐标,fij表示第i行第j列的像素值。
破片穿孔面积计算:根据统计提取的每个破片穿孔区域的像素数量,及校正等效靶板图像与其空间几何图形的尺寸比例关系,进行面积计算。
根据计算的面积对破片穿孔噪声进行滤除,并对有效破片穿孔进行统计,经过以上处理与计算,破片穿孔坐标与面积计算的结果如表1所示,靶板上实际有9个破片穿孔,识别出9个破片穿孔,并且多次进行了重复大量的数据识别,破片穿孔识别率为100%;破片穿孔面积计算精度提高20%。
表1破片穿孔位置与面积示例统计表
如图3所示,以上一实施例为基础进行补充,所述几何畸变校正单元101包括:
靶板轮廓提取单元1011,用于根据Otsu算法按图像的灰度特性将背景和靶板区分,靶板变为黑色,分割之后的结果是一幅二值图,再基于形态学边界提取的方法得到为白色四边形的靶板边界线;
边界线解析单元1012,用于根据Radon变换检测直线的方法得到边界线方程,并根据边界线方程获取顶点坐标参数,所述顶点坐标参数为靶板边界线四个顶点的位置左边;
畸变校正坐标变换单元1013,用于根据所述畸变校正矩阵H将顶点坐标参数进行变换得到校正坐标参数。
在本实施例中,存在畸变的形式主要有拍摄角度引起的图像倾斜畸变和镜头引起的畸变两种;
①对于拍摄角度引起的图像倾斜畸变,摄影测量系统实际测量时,由于拍摄方位的问题出现不同程度的变形,当像平面与被拍摄物平面平行时,如果被拍摄的实物放置不正,与镜头有一定的角度则会产生变形,即倾斜形变。一般情况下,图像变形分为三种形式:水平倾斜、垂直倾斜、水平垂直倾斜,如图9、10、11所示。
②对于镜头引起的畸变,镜头畸变的是一般是通过线性小孔成像模型校准,但是由于引起镜头畸变的原因较多,也受相机内部参数的影响,变成非线性模型。
相机成像的原理实质上是通过坐标系之间的转换实现的。首先将物方空间中任一点都可以转换到相机成像的坐标系下,然后将其投射到像平面上,最后将像平面上的数据点变换到像素坐标系下。但是载透镜制造、组装的过程中会引入一定的畸变,致使拍摄的图像失真,镜头的畸变包含径向畸变和切向畸变两种。
③对于径向畸变,径向畸变是顺着椭圆形透镜半径方向产生的畸变,原因是在透镜中心地方的光线比其它地方的更加弯曲。当成像光轴中心处没有畸变,沿着透镜半径方向由内向外移动,畸变会越来越严重。畸变的数学模型是由镜头上某点分布位置处径向方向上的调节公式为:
x0=x(1+k1r2+k2r4) (1)
y0=y(1+k1r2+k2r4) (2)
式中:(x0,y0)为图像坐标系畸变点的位置,(x,y)为图像坐标系理想点的位置,k1、k2为径向形变系数,
对于切向畸变,切向畸变是由于透镜自身和图像平面不平行而产生的。建立的畸变模型如下:
x0=x+[2p1xy+p2(r2+2x2)] (3)
y0=y+[2p2xy+p1(r2+2x2)] (4)
式中:(x0,y0)为图像坐标系畸变点的位置,(x,y)为图像坐标系理想点的位置,p1、p2为切向形变系数,
关于畸变校正方法的说明
当拍摄的时候,因拍摄角度和拍摄距离引起的图像产生的梯形畸变,形状失真。通过形态学提取边缘轮廓和阈值分割的方法,识别出图像的任意形状的边缘轮廓,利用Radon变换算法检测直线的特性,解析出四条边缘直线方程,求出四个顶点位置坐标;利用已知的顶点坐标和长宽比计算出校正之后的顶点位置坐标,利用这两组顶点坐标计算畸变校正矩阵,然后利用畸变校正矩阵可以对原图进行校正。
包括具体实施步骤如下。
①靶板边缘轮廓的提取
靶板的颜色和背景差别相对较大,采用Otsu算法按图像的灰度特性将背景和靶板区分,靶板变为黑色,分割之后的结果是一幅二值图,采用形态学边界提取的方法,提取白色四边形的边界直线。
边界提取的算法公式如下:
β(A)=A-(AΘB) (5)
其中,β(A)是二值图A的边界二值图,Θ是腐蚀算子,B是一个适当的结构元素。
处理过程:建立一个像素值全为1的3×3矩阵,即结构元素B。用结构元素B对上面获得的二值图进行腐蚀,得到一幅新的二值图。假定原来的二值图为BW,新的二值图为BW1。两幅图相减,β(BW)=BW-BW1,得到原图的边界轮廓二值图。
②边界直线方程的解析
提取轮廓之后的二值图上,只有四变形的4条边界直线是白色的,其余都是黑色的。利用Radon变换检测直线的特点,将4条边界直线的方程解析出来。
Radon变换可在任意维空间定义,给出二维空间的定义式:
式中,D为整个图像xy平面;f(x,y)为在图像点(x,y)的灰度;ρ为坐标原点到直线的距离;θ为距离与x轴的夹角;δ为狄雷克莱函数。它使f(x,y)沿直线ρ=x cosθ+y sinθ进行积分,如图12所示。
经过Radon变换后,图像中的每条直线会在ρ-θ空间形成一个亮点,直线的检测转化成为在ρ-θ变化域中亮点的检测。上面得到的边缘轮廓图中有4条直线,分别对应变换图中4个亮点,坐标(ρ,θ)可从图中读出,每个亮点对应一条直线,由已知(ρ,θ),根据公式ρ=xcosθ+y sinθ可以将4条直线的数学方程解析出来。利用两条直线方程联立获得方程组,求出4个顶点的位置坐标。
关于畸变校正的说明
根据图像顶点间的模型关系和图像的长宽比,计算出经过校正后的顶点位置坐标。假设校正之前的4个顶点的坐标为(x1,y1)、(x2,y2)、(x3,y3)、(x4,y4)和畸变校正后的顶点坐标为(x1',y1')、(x'2,y'2)、(x'3,y'3)、(x'4,y'4),令畸变校正矩阵H为:
引入中间参数h,h=(h11 h12 h13 h21 h22 h23 h31 h32 h33)T。
利用投影关系和坐标转换公式,得到下面公式:
只要将校正之前的4个顶点坐标和校正之后的4个顶点坐标带入上面公式,求出h,那么畸变校正矩阵H就可以求出,利用畸变校正矩阵就可以对原图像进行校正。
如图4所示,作为本发明一个优选的实施例,所述破片穿孔分割提取单元102包括:
图像预处理单元1021,用于根据自适应中值滤波方法对等效靶板图像进行去噪处理,用于保留等效靶板图像的边缘细节信息;
对比度增强单元1022,用于采用灰度变换的方法将等效靶板图像灰度等级扩大;
目标增强单元1023,用于根据Canny边缘检测的目标增强算法将等效靶板图像中的破片穿孔目标的边缘进行增强,所述增强后的破片穿孔目标的边缘线为零散状;
破片穿孔边缘优化处理单元1024,用于根据主动轮廓法对等效靶板图像进行优化处理;
特征提取单元1025,用于基于纹理特征对等效靶板图像的破片穿孔区域进行特征提取;
干扰清除单元1026,用于基于数学形态学对等效靶板图像进行干扰目标滤除;
目标识别单元1027,用于根据所述边缘特征和纹理特征对破片穿孔进行识别,分离出破片穿孔位置并进行标注得到破片穿孔图像。
本发明实施例中,对于等效靶板图像预处理技术,对于爆炸后所采集的等效靶板图像,在图像传感器成像、信号传输等过程中,由于多种因素的影响,会产生黑白相间的亮暗点会污染图像,噪声来源主要是椒盐噪声,这些小区域或者小的噪声点将影响系统对破片在等效靶板上着靶点位置的检测,因此采用自适应中值滤波对图像进行去噪处理,在滤除椒盐噪声时很好地保留边缘细节信息。
自适应中值滤波的方法是先通过图像处理技术检测图像中的像素是否存在噪声,若没有噪声则不用滤波。若含有噪声,则检测滤波器窗口中的值是否含有噪声,若没有噪声则输出其值。若其值是噪声,说明噪声密度比较大,变大滤波器窗口再次滤波;若增大到最大尺寸时,还无法实现滤波,则噪声密度非常的大,则求取窗口外的平均灰度值,用此均值代替像素输出。若窗口以外没有像素点,则采用折中取值,极小灰度像素和极大灰度像素的平均灰度值代替当前像素输出。
具体操作步骤如下:
(1)设Wij为需要处理的、中心在(i,j)处的滤波窗口,W0为初始窗口,Wmax为最大窗口。(i,j)处的灰度值为fij,极小灰度值为fmin,极大灰度值为fmax,灰度中值为fmed。
(2)如果fmin<fij<fmax,输出fij,否则转2);
(3)如果fmin<fmed<fmax,输出fmax,否则转3);
(4)如果Wij+2≤Wmax,增大窗口:Wij=Wij+2,转1),否则转4);
(5)最大窗口中需要输出极小灰度值、极大灰度值和灰度均值favg(其中k为灰度极值的像素数目):
如果favg不存在,即Wij*Wij=k,则转e);
(6)输出favg=(fmin+fmax)/2。
分段线性灰度变换增强图像对比度
由于相机成像系统的性能以及外界光照等影响,会出现所拍摄的等效靶板图像对比度低的现象,如很高的灰度值把暗区的信号都掩盖了或是很低的灰度值模糊了亮区的信号,极大地影响图像视觉效果。灰度变换则能使等效靶板图像的灰度等级扩大,对比度增强,明显改善了图像视觉效果。
为了突出等效靶板的破片穿孔所在的灰度区间,相对抑制那些不感兴趣的等效靶板背景灰度区间,采用三段线性灰度变换,对应的数学表达式为4-9。其对灰度区间[a,b]进行了线性拉伸,而灰度区间[0,a]和[b,Mf]被压缩。通过调整得到合适的两个拐点,实现任意灰度区间的拉伸或压缩。
通过以上的图像预处理,可以得到抑制椒盐噪声,对比度较好的图像,提高破片穿孔目标检测的精度。
基于边缘检测的等效靶板目标增强算法
采集的等效靶板图像上需要检测破片穿孔目标,为了更好地分割破片穿孔所在的区域,采用边缘检测提取目标区域的算法。Canny边缘检测是拥有高精度和信噪比的最优化思想的算法,基于Canny边缘检测的目标增强算法步骤如下:
(1)设置滤波器。利用滤波器去除一部分噪声平滑图像。大多数使用的滤波器有高斯滤波、均值滤波、中值滤波等,噪声主要为高斯白噪声,采用高斯滤波器。
(2)边缘增强一阶差分算子。通过采用一阶差分算子来增强图像的边缘,常使用的一阶差分算子包含有Roberts算子、Sobel算子、Prewitt算子等,其中Prewitt算子表达式为:
Gx(i,j)=f(i-1,j+1)+f(i,j+1)+f(i+1,j+1)-f(i-1,j-1)-f(i,j-1)-f(i+1,j-1) (4-10)
Gy(i,j)=f(i+1,j-1)+f(i+1,j)+f(i+1,j+1)-f(i-1,j-1)-f(i-1,j)-f(i-1,j+1) (4-11)
式中Gx(i,j)和Gy(i,j)分别是点(i,j)在X和Y方向的梯度函数,f(i,j)为坐标(i,j)点处的灰度值,其梯度幅值为通过一阶差分算子计算后的是梯度幅值图。
(3)抑制极大值。抑制极大值是通过计算获得的灰度梯度幅值,其值是不是沿梯度方向中最大值,若是最大值,输出梯度幅值,否则为0。
(4)区间阀值检测。设T1为低阀值,T2为高阀值,抑制极大值后的梯度幅值再与低阀值进行比较,大的为强边缘,小的为弱边缘。先搜索强边缘,发现存在任意一个非零点时,以此点为开始点顺着其8邻域追踪边缘检测,当出现断点时,采用从该断点开始在弱边缘处的8领域内开始搜索检测,寻找非零点记录到强边缘图中,不断循环采用此方法跟踪连接,直到全部连接完成。
通过以上步骤的操作,可以将等效靶板图像中的破片穿孔目标的边缘进行增强,但边缘检测增强的是零散的边缘段,没有形成目标封闭的轮廓,有必要对零散的边缘段进一步进行处理,以有利于后续的破片穿孔目标的分割与提取。基于主动轮廓的破片穿孔图像边缘提取算法
对于边缘检测到的是破片穿孔零散的边缘段,需要形成闭合的轮廓,以利于内部区域参数的计算;同时由于破片穿孔边缘灰度对比度的不一致性,Canny算子利用局部特征获取的边缘段也不是全局最优的位置,从而影响破片穿孔位置计算的准确性,所以需要对提取的边缘进一步进行优化处理。
作为主动轮廓的候选点,给出每个破片穿孔的初始轮廓,采用主动轮廓自适应地搜索出每个破片穿孔的最优边缘轮廓。主动轮廓法的能量函数为:
式(12)中,n为控制点数。公式的第一项为Snake的内量,Eint(i)=T(i)。Eint1(i)+U(i)Eint2(i)。Eint1(i)和Eint2(i)分别对应于式(12)中的一阶和二阶连续性强制力。T(i)与U(i)分别为离散的一阶和二阶项的系数,在Kass的算法中,使用|Vi-Vi-1|作为Eint1(i)的离散近似会导致曲线收缩,会将控制点挤到曲率高的轮廓线部分,所以采用离散近似形式为:
/>
式(13)中:是各控制点之间的平均空间距离。平均距离附近的点能量值较小,促使轮廓线上的点符合均匀分布。二阶连续性强制力Eint2和一阶的形式一样,即:
Eint2=|Vi-1-2Vi+Vi+1|2 (14)
为了各项的影响达到平衡,对一阶和二阶连续性强制力进行归一化处理,其值在[0,1]之间。
由于破片穿孔图象中穿孔和靶板之间存在比较明显的灰度差异,穿孔与靶板在其边缘连接处存在灰度过渡区域,采用双门限法判断穿孔的边缘。
基于纹理特征的破片穿孔等区域特征提取
(1)基于Gabor小波的纹理特征提取
纹理特征的提取通常采用Gabor小波变换。纹理特征识别的原理如下:
一个2-DGabor函数为:
式(15)中:ω为调制频率,Gaussian函数是由ψ(x,y)经过复数正弦函数调制。小波函数ψ(x,y)经过旋转、伸缩及平移获取一些自相似的Gabor小波函数ψmn(x,y)=a-mψ(x,y),m、n分别为小波函数的尺度和方向,m=0,1,...,M-1;n=0,1...,N-1。M、N分别为小波函数的尺度数和方向数。式中x=a-m(x cosθ+y sinθ),y=a-m(-x sinθ+y cosθ),其中:a为伸缩因子,a>1,θ=nπ/N。假设给定图像I(x,y)的尺寸为P×Q,离散型Gabor小波变换为:
式(16)中:是自相似函数ψmn(s,t)的共轭复数函数,s,t是用于设置滤波器模板大小的参数。
在不同尺度和方向上计算变换后的幅度序列其中m=0,1,...,M-1;n=0,1...,N-1。在不同尺度和方向上这些系数幅度值的大小表征图像含有的能量。图像的纹理特征由系数幅度序列均值μmn和标准方差σmn表示。
(2)基于灰度共生矩阵的纹理特征提取
能量、惯量、熵和相关性为灰度共生矩阵特征参数,分别表示纹理一致性、纹理对比等。要想获得提取纹理的参数,图像中相距(Δx,Δy)中的两个灰度像素,在同一时刻出现的分布联合频率,需要用灰度共生矩阵表示。若某图像灰度为N级,伴随着N×N的共生矩阵,用函数M(Δx,Δy)(h,k)表示。其中,(h,k)表示两个灰度像素,mh,k表示两个灰度像素相距为(Δx,Δy)时出现的频次。
图像亮度分量定义成64个灰度级,构造4个方向的灰度共生矩阵,然后分别计算M(1,0)M(0,1)M(1,1)M(1,-1),然后通过计算共生矩阵获取4个纹理特征参数。
能量J可以判定图像灰度分布是否均匀,其表示为:
当mh,k数值分布在主对角线附近时,J值越大;反之,越小。惯量G为。
对于粗纹理,当mh,k的数值集中在主对角线附近,(h,k)灰度值较小时,G值也小;对于细纹理,则G值大。熵S为。
灰度共生矩阵中的各mh,k数值相差很小且比较分散时,S值大;反之,S值小。相关性C为
式(21)中:μxμyσxσy分别为mxmy的均值和标准差,表示是对矩阵M中每行元素求和,/>表示是对矩阵M中每列元素求和。
求解上述4个纹理特征参数后,联合其均值和标准差,令μG,σG,μJ,σJ,μS,σS,μC,σC为纹理特征参数向量中的分量,对8个分量进行内部归一化处理。
基于数学形态学的干扰目标滤除
性能是由结构元素矩阵里面的元素决定。二值图像的结构元素只有0与1,结构元素H如图13所示。
在膨胀放大的运算过程中,选用半径值为r的圆盘形结构元素,将圆盘形结构元素放在图像的任意前景结构上,反之称为腐蚀运算。基于数字近景图像开展形态学滤波,其二维坐标为点集形态运算和结构元素的表现方法,二值图像I(u,v)可以表述为点集形式。与之相对应的Q点集涵盖坐标p=(u,v)的全部像素点。
Q={p|I(p)=1} (22)
膨胀是直观性的形态学集合运算,定义如下:
I⊕H={(p+q)|p∈I且q∈H} (23)
集合I与H均为坐标点对向量膨胀运算产生的点集,腐蚀运算与之相反。
形态滤波复合运算为膨胀和腐蚀运算有机融合在一起的产物,通过二者的半对偶性完成形态滤波的效果,利用结构元素膨胀和腐蚀的先后顺序,融合成复合运算。
二值开运算就是先腐蚀后膨胀复合运算。
先提取腐蚀运算里面的结构比结构元素更小的前景结构,后经过膨胀其结构被平滑,恢复到原来的尺寸大小。
闭运算开运算相反,先膨胀后腐蚀的复合运算,用I·H表示。
前景结构里面中比结构元素小的缝隙和孔洞会被闭运算填补剔除。
基于融合特征的目标识别算法
根据破片穿孔轮廓区域内的边缘特征和纹理特征,对破片穿孔进行识别,分离出破片穿孔。
首先选取破片穿孔的位置为中心点,再选取其周边24个点构成一个大小5×5的窗口区域,然后对这个窗口区域内所有的点以加权的方式计算作为这个窗口的区域能量值,得到当前融合系数的重要性测度:
式(26)中,ω(l,t)为权值系数,(2L+1)代表窗口的长,(2T+1)代表窗口的宽,I=A/B代表原图像A或B。然后可以利用重要性测度进一步得到相应窗口区域内系数的归一化相关度/>如式(4-27)所示:
给定相似度阈值α=0.9,如果小于设定的相似度域值α,则第i个区域的融合处理如式(28)所示:
大于α,则第i个区域的融合处理如公式(4-29)所示:
式(29)中,ωA(x,y)、ωB(x,y)为相应的权值系数,且ωB(x,y)=-ωA(x,y)。
通过以上算法,得到破片穿孔图像,采用连通区域提取每个破片穿孔区域,并进行标注。
如图5所示,作为本发明一个优选的实施例,所述破片穿孔位置三维反演模块200包括:
坐标系变换单元201,用于根据所述破片穿孔中心定位参数建立破片穿孔二维坐标系,并基于旋转矩阵R将穿孔二位坐标系并将其转化为破片穿孔三维坐标系;
三维反演单元202,用于基于变换矩阵M将不同的等效靶板图像的破片穿孔二维坐标系反演到统一的已战斗部爆心为基准的破片穿孔三维坐标系。
本发明实施例中,对于获得的战斗部爆炸后单块靶板上破片穿孔的中心定位,通过光速平差法实现破片穿孔位置三维反演,并编制数据处理软件,实现破片飞散参数的快速计算。
选取同一块等效靶板上四个反演的单点标志数据,拟合平面边缘获取每块等效靶板。依据布设靶板的方位数据,以及已知的高程,以左右两侧等效靶板平面的中分面为YOZ平面,作为第一基准;以正前方为Y轴,基于调平后的高度向上偏移至已知高度,从正前方至爆心的距离已知,基于点线面原则建立理论爆心坐标系O′-X′Y′Z′,如图14所示。
假设原坐标系O-XYZ在爆心坐标系O′-X′Y′Z′中的位姿已知,原点O对原点O′的平移向量T=[Δx Δy Δz]T,X轴对X′轴的旋转角为ψ,Y轴对Y′轴的旋转角为θ,Z对Z′轴的旋转角为φ,那么原坐标系O-XYZ相对爆心坐标系O′-X′Y′Z′的旋转变换矩阵R为:
其中,
如果等效靶板形状反演数据在原坐标系O-XYZ中的坐标为(X1,Y1,Z1),那么转换到爆心坐标系O′-X′Y′Z′中的坐标(X1′,Y1′,Z1′)为:
从而将原来坐标系O-XYZ的坐标转换到爆心坐标系O′-X′Y′Z′。
破片穿孔位置的二维坐标的三维反演
①靶板坐标系到爆心坐标系的变换模型建立
等效靶板上的破片穿孔二维坐标是以当前每块等效靶板建立的坐标系,为o-xyz坐标系,其中z=0。三维反演的坐标系是以战斗部爆心建立的坐标系,为O′-X′Y′Z′。任意两个空间直角坐标系都可以通过旋转和平移转换到同一坐标系下。设战斗部爆心坐标系O′-X′Y′Z′先平移(X0,Y0,Z0)),再依次绕X轴、旋转后的Y轴和旋转后的Z轴旋转ω、κ后,与等效靶板坐标系o-xyz重合,如图15所示:
对于破片穿孔点P和p其实质是一个三维坐标变换,爆心坐标系和所在等效靶板坐标系破片穿孔点P和p之间的向量关系,由射影几何知识得。
设R是由旋转角ω、κ确定的旋转变换矩阵,又因为:
将(36)式代入(35)式,即得到爆心坐标系O′-X′Y′Z′下破片穿孔三维反演点P(X,Y,Z)与其在等效靶板坐标系o-xyz中对应的点p(x,y,z)有以下关系。
其中:旋转矩阵R是一个3×3的矩阵:
旋转矩阵R为正交矩阵,R中的9个结构元素是3个独立的旋转角ω、和κ的函数。本项目采用的是ω、/>κ转角顺序,经推导R中的各元素值如下:
式(39)中旋转矩阵R与平移向量t=[X0,Y0,Z0]-1共12个未知数,可组成3×4的变换矩阵。但由于破片穿孔在等效靶板坐标系的z=0,旋转矩阵的第3列没有作用,去掉,记旋转矩阵R的第i列为ri,得到3×3的矩阵M=[r1 r2 t]。对M矩阵进行归一化后,变换矩阵M=[r1r2 t]共8个未知数,采用4对单点标志的等效靶板坐标系坐标,以及相应的三维反演的战斗部爆心坐标系的坐标,建立如下变换矩阵计算模型:
通过模型(33),利用破片穿孔在爆心坐标系和等效靶板坐标系下的一对坐标(X,Y,Z)和(x,y)可确定关于M矩阵参数m1-m8的两个方程:
于是,利用单点标志在爆心坐标系和等效靶板坐标系下的四对坐标(X,Y,Z)和(x,y)即可求解自由度为8的变换矩阵M,从而可求得两个坐标系的旋转矩阵R和平移向量t。
②破片穿孔二维坐标的三维反演快速算法提出
在已知变换矩阵M的情况下,通过式(41)矩阵变换,可以将不同组等效靶板图像的破片穿孔的二维坐标反演到统一的战斗部爆心坐标系中的三维坐标。
在破片穿孔二维坐标的三维反演模块中,首先基于爆心坐标系下等效靶板边缘四个已知单点标志三维坐标和相应的等效靶板坐标系下二维坐标,计算破片穿孔二维坐标的三维反演变换矩阵M;然后将表1中破片穿孔中心二维xy坐标加入“Z=0”这一列,利用最小二乘原理进行求解,通过“分析-最佳拟合转化-点到点”,将等效靶板破片穿孔二维坐标转换至爆心坐标系下的三维坐标,如图16所示,三维反演的破片穿孔位置的三维坐标如表2所示。
表2破片穿孔位置三维反演示例
如图6所示,作为本发明一个优选的实施例,所述破片飞散参数计算模块300包括:
破片分布密度计算单元301,用于基于破片穿孔三维坐标系并参照球面投影的方法计算得到破片分布密度参数;
联合计算单元302,用于对等效靶板图像进行虚拟网格划分并进行计算得到破片飞散角参数和破片飞散方向角参数。
本发明实施例中,根据破片穿孔相对于爆心的三维坐标,按照中华人民共和国国家军用标准GJB5232.2-2004,战术导弹战斗部靶场试验方法第2部分:静爆试验破片飞散特性测试计算破片飞散角、破片飞散方向角和破片分布密度,计算破片场飞散的3个参数。
球面投影有效破片密度计算
根据战斗部放置方向,对数据坐标进行坐标系转换,给出破片总数N和分区角度结合等效靶板统计每一垂直分格内相邻两条分划线之间的有效破片穿孔数Ni,按照下式计算出该分格的破片密度νi,按照表3进行统计。
式中:ΔSi:为第i个铅垂分格的面积,单位为m2。
给出布靶半径R,等效靶板高度h,设等效靶板的内接球如图17所示;
图示中上面的圆环为等效靶板在球面上的投影,圆环半径为r,与平行的赤道面的夹角为圆环宽度对应的圆心角为/>圆环宽度为/>圆环周长为2πr,所以圆环上分区面积为:
其中,m为分区的数量。
分格序号 | 1 | 2 | …… | i | …… | n-1 | n |
对应方向角(°) | φ1 | φ2 | …… | φi | …… | φn-1 | φn |
有效破片数 | N1 | N2 | …… | Ni | …… | Nn-1 | Nn |
破片密度(枚/米2) | ν1 | ν2 | …… | νi | …… | νn-1 | νn |
表3破片场密度分布计算表
破片飞散角和飞散方向角的计算
对等效靶板进行虚拟网格划分,按照表一统计每一铅垂分格内的破片穿孔数Ni,并计算Ni与Na的比值,按照下式计算破片比率积分函数:
并将计算结果光滑连接成破片比率积分曲线。破片飞散角Ω为:
Ω=φ0.95-φ0.05 (45)
式中:φ0.95:F(φ)等于0.95所对应的角度,单位为度;φ0.05:F(φ)等于0.05所对应的角度,单位为度;
破片飞散方向角为:
ψ=90°-φ0.5 (46)
式中:φ0.5:F(φ)等于0.5所对应的角度,单位为度;
按照表3统计靶板的破片数据,表4为破片飞散角、飞散角方向数据处理表:
表4破片飞散角、飞散角方向数据处理表
则φ0.05、φ0.5、φ0.95分别按下式计算:
式中:
φ0--从靶板的正航向一端数起第一枚破片所在分格的正向边线所对应的φ角,单位为度;
k--φ0.05所在分格的序数;
m--φ0.5所在分格的序数;
r--φ0.95所在分格的序数;
Nk--φ0.05所在铅垂分格内的破片穿孔数;
Nm--φ0.5所在铅垂分格内的破片穿孔数;
Nr--φ0.95所在铅垂分格内的破片穿孔数;
Ni--第i个铅垂分格内的破片穿孔数;
Δ--每一分格的两条边线与爆心组成的两个面夹角的圆心角,单位为度。
按照上述计算公式得到破片分散角、破片分散方向角。
如图7所示,作为本发明一个优选的实施例,所述系统测量模块400包括:
识别率核算单元401,用于根据相机成像分辨率和光照角度两个方面对破片穿孔识别率进行核算,使得破片穿孔识别率大于98%;
面积核算单元402,用于根据破片穿孔的像元分辨率对穿孔面积进行核算,使得破片穿孔面积大于4.2%;
位置核算单元403,用于根据战斗部爆炸后破片穿孔和单点标志的二维位置和战斗部爆炸前单点标志的三维反演坐标三个指标计算得到破片穿孔位置测量精度。
如图8所示,以上一实施例为基础进行补充说明,所述位置核算单元403包括:
二维位置精度计算单元4031,用于根据破片穿孔和单点标志的像元位置的提取精度、靶板的长、宽尺寸测量的精度确定战斗部爆炸后破片穿孔和单点标志的二维位置的定位误差;
三维反演精度计算单元4032,用于根据三维反演原理确定战斗部爆炸前单点标志三维坐标的位置测量精度。
在本实施例中,测量系统衡量的指标参数有破片穿孔识别率、破片穿孔面积和破片穿孔位置3个参数,所述破片穿孔识别率也被称为捕获率,表明转化为图像后破片穿孔相比对实际保留程度。
关于破片穿孔识别率/捕获率如何获得
测量系统对破片穿孔的识别率/捕获率受靶板尺寸、相机成像距离、相机成像分辨率、光照角度、靶板材质、相机成像角度,以及破片穿孔图像处理算法等影响。其中相机成像分辨率以及光照角度起决定性作用,下面主要从这两个方面进行分析。
(1)相机成像分辨率
比如采用靶板的尺寸一般为2m*1.2m、6m*1.5m、2m*1m等规格,而相机的分辨率为6016×4016像元,为了将单块等效靶板全部成像在图像上,需要相机的成像视场略大于等效靶板的尺寸。由于等效靶板的长宽比例与相机的两个方向的分辨率不一致,所以按照等效靶板的长边计算相机视场。假设等效靶板为最大尺寸6m*1.5m,采集视场尺寸范围为6.5m*2m,那么相机在两个方向的物体空间分辨率取较低的分辨率为:Re1=6.5*1000/6016=1.08mm/像元。最小的破片穿孔目标直径为6个像素,所以其成像像元数量为:Count=pi*(6/2)^2/(1.08)^2=24像元
根据以上计算,在所选择相机条件下,对于最小的破片穿孔,从采集到破片穿孔的像元分辨率上,可以有效地采集到破片穿孔目标,为可靠地识别破片穿孔提供了可能性。
(2)光照角度计算
为了保证成像的质量,以便于后续的图像处理算法识别,将等效靶板底部放置在离地面1m~1.5m的高度,相对于太阳光为逆光拍摄;如果等效靶板不便移动,采用在等效靶板后面放置反光布进行逆光照明,使得破片穿孔的灰度具有高度的一致性,并辅以手动选择方式进行判断,使得破片穿孔识别率/捕获率达到98%以上。
关于破片穿孔面积如何测量计算
破片穿孔面积测量精度与破片穿孔的像元分辨率有关,识别的最小破片穿孔是6mm,其像元分辨率为24个,那么其面积测量的分辨率为:
当面积测量的绝对精度为±1.2m2,其对应的面积相对测量精度为:
如果破片穿孔直径增大,那么其成像的像元数就会增多,但面积测量的分辨率Resarea是不变的。当面积测量的分辨率Resarea不变的情况下,相对测量精度会随着破片穿孔直径的增大而减小,从而其相对测量精度就会提高。所以破片穿孔直径6mm以上时,面积测量的精度最低为4.2%。
关于破片穿孔位置如何测量
破片穿孔位置测量是以爆心为原点的坐标系的三维坐标,所以破片穿孔位置测量精度与破片穿孔的二维位置,二维位置的反演变换矩阵有关系。破片穿孔的二维位置的提取精度是由工业测量相机的分辨率和检测算法确定的。破片穿孔二维位置的反演变换矩阵是由单点标志的二维坐标,以及单点标志的三维反演坐标计算出来的。单点标志的三维反演坐标是通过单点的二维坐标,与相机的位姿计算出来的。通过以上分析可知,单个破片穿孔位置测量精度的主要影响因素包括战斗部爆炸后破片穿孔、单点标志的二维位置和战斗部爆炸前单点标志的三维反演坐标。
(1)战斗部爆炸后破片穿孔、单点标志的二维位置计算精度
对于破片穿孔的二维位置测量,是由破片穿孔和单点标志的像元位置(u,v)的提取精度,靶板的长、宽尺寸测量的精度共同决定的。如果等效靶板的尺寸L*H的测量精度为Δxbmm,成像的像元数量为uL*vH个,破片穿孔和单点标志的二维像素位置为(u,v)。那么由于像素位置(u,v)的提取精度引入的破片穿孔和单点标志二维位置定位误差Δxh1=1个像元,对应的实际物理尺寸也就是相机成像时的空间分辨率:
那么由于尺寸测量引入的破片穿孔和单点标志二维位置定位误差为:
如果等效靶板的尺寸6m*1.5m的测量精度为Δxb=1mm,成像的像元数量为uL*vH=4000*1000个,那么相机成像时的空间分辨率为Δxh1=1.5mm,破片穿孔和单点标志二维位置定位误差Δxh2的最大值为1mm。
(2)战斗部爆炸前单点标志的三维反演计算精度
对于爆炸前等效靶板中单点标志的三维反演坐标,根据三维反演原理,设相机的有效焦距为f,光轴与x轴的夹角为α1、α2,ω1、ω2为小于相机视场角的偏航方向的投影角。β1、β2为俯仰方向的投影角,单点标志的三维反演坐标P(x,y,z)为:
单点标志的三维反演坐标P(x,y,z)分别在偏航方向ω1、ω2反演的精度为:
/>
单点标志的三维反演坐标P(x,y,z)分别在俯仰方向β1、β2反演的精度,类似偏航方向ω1、ω2反演精度的计算,可以推导出。
设相机在偏航方向提取精度为δω1、δω2,在俯仰方向提取精度为δβ1、δβ2,则反演点P(x,y,z)的测量精度为:
通过式(35)爆心坐标系转换,以及式(40)二维坐标的三维反演变换,可知:破片穿孔三维反演位置坐标P(XP,YP,ZP)与战斗部爆炸后破片穿孔二维位置坐标p(xp,yp)、单点标志二维位置坐标b(xb,yb)和爆炸前单点标志三维反演坐标B(XB,YB,ZB)有关。由于方程的非线性,它们之间的关系无法写成解析式。为了分析破片穿孔三维反演位置的测量误差,它们之间的关系简写为:
P(XP,YP,ZP)=F[p(xp,yp),b(xb,yb),B(XB,YB,ZB)] (62)
结合式(61)得到的战斗部爆炸前单点标志三维反演的精度Δp,式(52)和(53)可得到战斗部爆炸后破片穿孔和单点标志的二维位置计算精度Δx1和Δx2,根据式(62)即可得到破片穿孔三维反演位置测量精度。
以上仅对本发明的较佳实施例进行了详细叙述,并不用以限制本发明,凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等,均应包含在本发明的保护范围之内。
应该理解的是,虽然本发明各实施例的流程图中的各个步骤按照箭头的指示依次显示,但是这些步骤并不是必然按照箭头指示的顺序依次执行。除非本文中有明确的说明,这些步骤的执行并没有严格的顺序限制,这些步骤可以以其它的顺序执行。而且,各实施例中的至少一部分步骤可以包括多个子步骤或者多个阶段,这些子步骤或者阶段并不必然是在同一时刻执行完成,而是可以在不同的时刻执行,这些子步骤或者阶段的执行顺序也不必然是依次进行,而是可以与其它步骤或者其它步骤的子步骤或者阶段的至少一部分轮流或者交替地执行。
本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程,是可以通过计算机程序来指令相关的硬件来完成,所述的程序可存储于一非易失性计算机可读取存储介质中,该程序在执行时,可包括如上述各方法的实施例的流程。其中,本申请所提供的各实施例中所使用的对存储器、存储、数据库或其它介质的任何引用,均可包括非易失性和/或易失性存储器。非易失性存储器可包括只读存储器(ROM)、可编程ROM(PROM)、电可编程ROM(EPROM)、电可擦除可编程ROM(EEPROM)或闪存。易失性存储器可包括随机存取存储器(RAM)或者外部高速缓冲存储器。作为说明而非局限,RAM以多种形式可得,诸如静态RAM(SRAM)、动态RAM(DRAM)、同步DRAM(SDRAM)、双数据率SDRAM(DDRSDRAM)、增强型SDRAM(ESDRAM)、同步链路(Synchlink)DRAM(SLDRAM)、存储器总线(Rambus)直接RAM(RDRAM)、直接存储器总线动态RAM(DRDRAM)、以及存储器总线动态RAM(RDRAM)等。
本领域技术人员在考虑说明书及实施例处的公开后,将容易想到本公开的其它实施方案。本申请旨在涵盖本公开的任何变型、用途或者适应性变化,这些变型、用途或者适应性变化遵循本公开的一般性原理并包括本公开未公开的本技术领域中的公知常识或惯用技术手段。说明书和实施例仅被视为示例性的,本公开的真正范围和精神由权利要求指出。
Claims (8)
1.基于图像多特征提取与融合算法的破片场三维快速反演系统,其特征在于,所述系统包括:
破片穿孔位置和面积计算模块,用于获取战斗部爆炸后的等效靶板图像,并通过几何畸变校正、靶板区域提取以及图像处理得到破片穿孔中心定位参数与破片穿孔面积参数;
破片穿孔位置三维反演模块,用于根据所述破片穿孔中心定位参数建立破片穿孔二维坐标系并将其转化为破片穿孔三维坐标系,所述破片穿孔三维坐标系为破片穿孔相对于爆心建立;
破片飞散参数计算模块,用于根据所述破片穿孔三维坐标系并参照国标计算得到破片飞散角参数、破片飞散方向角参数和破片分布密度参数;
系统测量模块,用于根据破片穿孔识别率、破片穿孔面积以及破片穿孔位置三个指标参数进行精度计算。
2.根据权利要求1所述的基于图像多特征提取与融合算法的破片场三维快速反演系统,其特征在于,所述破片穿孔位置和面积计算模块包括:
几何畸变校正单元,用于根据畸变校正矩阵H对所述等效靶板图像进行校正;
破片穿孔分割提取单元,用于根据破片穿孔与非穿孔部分灰度值的差别,并结合边缘和纹理特征对等效靶板图像进行提取、识别得到破片穿孔图像,所述破片穿孔图像表现形式为等效靶板图像上的破片穿孔位置具有标注;
破片位置与面积计算单元,用于根据破片穿孔图像对每个破片穿孔位置进行提取和统计,并基于灰度边缘轮廓与灰度加权计算得到破片穿孔中心定位参数与破片穿孔面积参数。
3.根据权利要求2所述的基于图像多特征提取与融合算法的破片场三维快速反演系统,其特征在于,所述几何畸变校正单元包括:
靶板轮廓提取单元,用于根据Otsu算法按图像的灰度特性将背景和靶板区分,靶板变为黑色,分割之后的结果是一幅二值图,再基于形态学边界提取的方法得到为白色四边形的靶板边界线;
边界线解析单元,用于根据Radon变换检测直线的方法得到边界线方程,并根据边界线方程获取顶点坐标参数,所述顶点坐标参数为靶板边界线四个顶点的位置左边;
畸变校正坐标变换单元,用于根据所述畸变校正矩阵H将顶点坐标参数进行变换得到校正坐标参数。
4.根据权利要求2所述的基于图像多特征提取与融合算法的破片场三维快速反演系统,其特征在于,所述破片穿孔分割提取单元包括:
图像预处理单元,用于根据自适应中值滤波方法对等效靶板图像进行去噪处理,用于保留等效靶板图像的边缘细节信息;
对比度增强单元,用于采用灰度变换的方法将等效靶板图像灰度等级扩大;
目标增强单元,用于根据Canny边缘检测的目标增强算法将等效靶板图像中的破片穿孔目标的边缘进行增强,所述增强后的破片穿孔目标的边缘线为零散状;
破片穿孔边缘优化处理单元,用于根据主动轮廓法对等效靶板图像进行优化处理;
特征提取单元,用于基于纹理特征对等效靶板图像的破片穿孔区域进行特征提取;
干扰清除单元,用于基于数学形态学对等效靶板图像进行干扰目标滤除;
目标识别单元,用于根据所述边缘特征和纹理特征对破片穿孔进行识别,分离出破片穿孔位置并进行标注得到破片穿孔图像。
5.根据权利要求1所述的基于图像多特征提取与融合算法的破片场三维快速反演系统,其特征在于,所述破片穿孔位置三维反演模块包括:
坐标系变换单元,用于根据所述破片穿孔中心定位参数建立破片穿孔二维坐标系,并基于旋转矩阵R将穿孔二位坐标系并将其转化为破片穿孔三维坐标系;
三维反演单元,用于基于变换矩阵M将不同的等效靶板图像的破片穿孔二维坐标系反演到统一的已战斗部爆心为基准的破片穿孔三维坐标系。
6.根据权利要求1所述的基于图像多特征提取与融合算法的破片场三维快速反演系统,其特征在于,所述破片飞散参数计算模块包括:
破片分布密度计算单元,用于基于破片穿孔三维坐标系并参照球面投影的方法计算得到破片分布密度参数;
联合计算单元,用于对等效靶板图像进行虚拟网格划分并进行计算得到破片飞散角参数和破片飞散方向角参数。
7.根据权利要求1所述的基于图像多特征提取与融合算法的破片场三维快速反演系统,其特征在于,所述系统测量模块包括:
识别率核算单元,用于根据相机成像分辨率和光照角度两个方面对破片穿孔识别率进行核算,使得破片穿孔识别率大于98%;
面积核算单元,用于根据破片穿孔的像元分辨率对穿孔面积进行核算,使得破片穿孔面积大于4.2%;
位置核算单元,用于根据战斗部爆炸后破片穿孔和单点标志的二维位置和战斗部爆炸前单点标志的三维反演坐标三个指标计算得到破片穿孔位置测量精度。
8.根据权利要求7所述的基于图像多特征提取与融合算法的破片场三维快速反演系统,其特征在于,所述位置核算单元包括:
二维位置精度计算单元,用于根据破片穿孔和单点标志的像元位置的提取精度、靶板的长、宽尺寸测量的精度确定战斗部爆炸后破片穿孔和单点标志的二维位置的定位误差;
三维反演精度计算单元,用于根据三维反演原理确定战斗部爆炸前单点标志三维坐标的位置测量精度。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311687287.5A CN117635874A (zh) | 2023-12-11 | 2023-12-11 | 基于图像多特征提取与融合算法的破片场三维快速反演系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311687287.5A CN117635874A (zh) | 2023-12-11 | 2023-12-11 | 基于图像多特征提取与融合算法的破片场三维快速反演系统 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN117635874A true CN117635874A (zh) | 2024-03-01 |
Family
ID=90028706
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202311687287.5A Pending CN117635874A (zh) | 2023-12-11 | 2023-12-11 | 基于图像多特征提取与融合算法的破片场三维快速反演系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN117635874A (zh) |
-
2023
- 2023-12-11 CN CN202311687287.5A patent/CN117635874A/zh active Pending
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN111260731B (zh) | 一种棋盘格亚像素级角点自适应检测的方法 | |
US9123112B2 (en) | Method for the pre-processing of a three-dimensional image of the surface of a tyre for use in the inspection of said surface | |
CN110033516B (zh) | 基于双目相机图像采集识别的针片状颗粒含量检测方法 | |
CN109859226B (zh) | 一种图形分割的棋盘格角点亚像素的检测方法 | |
CN112017223B (zh) | 基于改进SIFT-Delaunay的异源图像配准方法 | |
CN111369605B (zh) | 一种基于边缘特征的红外与可见光图像的配准方法和系统 | |
WO2010116885A1 (ja) | データ処理装置、画像照合方法、プログラムおよび画像照合システム | |
CN110222661B (zh) | 一种用于运动目标识别及跟踪的特征提取方法 | |
CN113256653B (zh) | 一种面向高层地物的异源高分遥感影像配准方法 | |
CN108647580B (zh) | 基于改进sift引导isar图像特征点提取匹配方法 | |
CN112435252B (zh) | 一种战斗部破片穿孔和凹坑检测方法 | |
JP5023238B2 (ja) | 薬莢ベースの自動領域セグメンテーションと薬莢比較用のベスト痕跡領域の選択の方法 | |
CN115880373A (zh) | 基于新型编码特征的立体视觉系统的标定板及标定方法 | |
CN116358449A (zh) | 一种基于双目面结构光的飞机铆钉凹凸量测量方法 | |
CN111950559A (zh) | 一种基于径向灰度的指针仪表自动读数方法 | |
Wang | Automatic extraction of building outline from high resolution aerial imagery | |
CN104268550A (zh) | 特征提取方法及装置 | |
CN107765257A (zh) | 一种基于反射强度辅助外部校准的激光探测与测量方法 | |
CN116883446B (zh) | 一种车载摄像头镜片碾磨程度实时监测系统 | |
CN112734816B (zh) | 基于CSS-Delaunay的异源图像配准方法 | |
CN117635874A (zh) | 基于图像多特征提取与融合算法的破片场三维快速反演系统 | |
CN116206139A (zh) | 一种基于局部自卷积的无人机图像升尺度匹配方法 | |
CN116052020A (zh) | 基于无人机的图像快速解译方法 | |
CN115330832A (zh) | 一种基于计算机视觉的输电铁塔全自由度位移监测系统和方法 | |
CN113822361A (zh) | 一种基于汉明距离的sar图像相似程度度量方法和系统 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination |