基于MATLAB图像处理的水泵特性曲线自动提取

丁 聪，邸南思

(上海勘测设计研究院有限公司，上海 200434)

1 概述

水泵、水轮机特性曲线是泵站、电站设计的重要基础资料，广泛用于原型机组参数计算、水力过渡过程分析[1-2]、经济效益评价[3]等。日常设计需要从期刊、书籍等资料中收集水泵、水轮机特性曲线。但资料多为纸质、图片格式，需进行人工读点、数值化。而人工数值化费时费力，且主观误差较大。为减少图像曲线数据提取工作量，提高曲线数据提取准确性，已有图像曲线数据提取的软件开发，如Getdata Dragh Digitizer、FindGraph、Digitize。但此类软件一般为手动选点或无法较好地适用于含纵横坐标轴网格线的图像。

Matlab因其强大的矩阵处理能力，在图像处理方面具有较大优势。MATLAB图像处理研究已经涉及到医学[4]、食品[5]、交通[6]等各行各业。并且已有研究利用Matlab进行图像曲线数据提取，但多数仍停留在手动选点、曲线拟合或仅适用含单条曲线的图像，或未考虑曲线走向[7-9]。为此提出利用已识别曲线像素，寻求周围曲线像素并结合曲线斜率不可突变作为判断条件，甄别取舍寻求到的周围像素点进行曲线自动追踪，并根据像素坐标与真实坐标一一对应关系进行曲线数据换算，以避免手工操作的繁琐和曲线追踪的偏离。

2 图像处理难点及思路

若曲线与背景绝对分开，可降低曲线追踪难度。收集的曲线图像多为彩色图像，因此，需先将图像灰度化、二值化处理。

若曲线很细，仅为曲线中心线，一可提高曲线数据提取准确性；二可方便曲线追踪。因此，需先对曲线细化处理。

收集的曲线图像可能有褪色、变质、偏转等问题，致使影响曲线提取的因素较多(噪声大)，需先进行图像预处理(如：去噪、旋转、裁剪)。

多数收集的曲线图像内含有坐标轴、网格线，且曲线、网格线交叉较多。较多的交叉点增加了曲线追踪难度，因此，需将图像中的坐标轴、网格线先行删除。

处理好的曲线图像便可利用遍历法及斜率预测法进行曲线的追踪识别，并可根据坐标、像素点位置对应关系进行曲线数据的换算导出。

水泵特性曲线一般比较复杂、且有多条，智能分条识别难度较大，可采用人为指定需识别曲线减小识别难度。

依据上述图像处理难点及解决方法，曲线处理可采用图1流程。

图1 曲线提取处理流程示意

3 图像处理过程

3.1 图像读入

图像处理前，需将图像读入Matlab工作空间。Matlab提供有图像文件的读入和显示函数imread、imshow，用来读取与显示如bmp、jpg、tiff、hdf等格式文件。

3.2 图像灰度化及灰度增强

图像灰度化就是由RGB三通道数据的彩色图像变为单通道数据图像。灰度图像上每个像素的颜色值称为灰度，范围为0～255。图像灰度化常用方法有3种：① 分量法，将彩色图像的R、G、B三分量的亮度之一作为灰度图像的灰度值；② 最大值法，将彩色图像的R、G、B三分量亮度最大值作为灰度图像的灰度值；③ 加权平均法，将彩色图像的R、G、B三分量以不同的权重进行加权平均。Matlab提供的灰度化函数为rgb2gray，使用的是加权平均法。由于人眼对绿色敏感最高，对蓝色敏感最低，故采用心理学灰度公式：

Gray= 0.114B+ 0.587G+ 0.299R

(1)

灰度增强是指按一定变换关系逐点改变原图中每一个像素的灰度值，以改善图像视觉效果。灰度增强常用的方法有线性灰度变化、非线性灰度变化、直方图均衡化等。直方图是图像的重要统计特征，是表示数字图像中每一灰度级与该灰度级出现的频数间的统计关系。直方图均衡化是通过灰度变化将一幅图像转换为另一幅均衡直方图，即在每个灰度级上都具有相同的像素点数的过程。采用直方图均衡化可使原始图像的灰度直方图从比较集中的某个灰度区间变成在全部灰度范围内的均匀分布。Matlab直方图均衡化函数为histeq。

3.3 图像去噪

图像去噪是减少数字图像中噪声的过程。图像去噪的方法有均值滤波、自适应维纳滤波、中值滤波、形态学噪声滤波、小波变换等方法。小波变换在时域、频域同时具有良好的局部化性质和多分辨率，可有效地把信号和噪声区别开来。

小波分析进行图像去噪复原主要有3个步骤[10]：① 对图像信号进行小波分解。选择合适的小波和恰当的分解层次，然后对含有噪声的图形进行N层小波分解。② 对经过层次分解后的高频系数进行阈值量化。对分解的每一层，选择一个恰当的阈值，并对该层高频系数进行阈值量化处理。③ 利用二维小波重构图像信号。同样，根据小波分解后的第N层近似(低频系数)和经过闭值量化处理后的各层细节(高频系数)，来计算二维信号的小波重构。

Matlab小波去噪函数为ddencmp(阈值选取函数)、wdencmp(去噪函数)。

3.4 图像二值化

灰度阈值法是最古老的图像二值化方法[11]。选取适当阈值，像素灰度值超过阈值则分配到最大灰度，低于阈值分配到最小灰度。

阈值的选取方法可分为全局阈值法和局部阈值法。全局阈值法指对整幅图像选取单一的阈值进行二值化，在物体和背景的灰度差别明显时效果较好。曲线图像中曲线与背景的灰度差别较大，全局阈值法较为合适。典型的全局阈值法包括Otsu 法、全局迭代法、最大墒方法等。Otsu法主要思想是从图像灰度直方图中把灰度值的集合用阈值T分为两类，然后根据两个类别的平均值方差(类间方差)和各类的方差(类内方差)的比为最大来确定阈值T。

假定给定的图像有L个灰度级，设阈值为T，把具有T以下灰度值的像素和具有比T大的灰度值的像素分为两类，设为类1和类2，类1中的像素数量为W1(T)，灰度值的均值为M1(T)，方差为σ1(T)，同样类2中的像素量为W2(T)，灰度值的均值为M2(T)，方差为σ2(T)。全体像素的平均值为MT，则类内方差由下式计算：

(2)

类间方差由下式计算：

(3)

Matlab最大类间差法全局阈值求解函数为graythesh。

3.5 曲线细化

图像的细化算法较多，数学形态学细化算法提取的图像骨架较连续，断点少，对噪声不敏感，且提取的边缘较光滑[12-13]。细化过程根据击中或击不中变换定义，图像二值文件A除去被结构元素序列B击中的元素。

A⊗B=A-(A#B)=A∩(A#B)C

(4)

图像细化是对图像文件进行基于结构元素序列{B}的细化

A⊗{B}=((…((A⊗B1)⊗B2)…)⊗Bn

(5)

即连续使用{B}对A细化。

Matlab图像数学形态学函数为bwmorph，其中细化操作为thin(细化)和skel(骨架化)。

3.6 图像旋转、裁剪、坐标轴网格线删除

由于曲线图像中的坐标轴、网格线一般均为直线。可通过寻找计算坐标轴、网格线所在直线与水平、竖直存在的夹角进行图像旋转。

直线寻找可采用Hough变换，其原理是利用图像空间和Hough参数空间的点、线对偶性，把图像空间中的检测问题转换到参数空间，通过在参数空间里进行简单的累加统计，然后在Hough参数空间寻找累加器峰值的方法检测直线。

Hough变换的参数空间为极坐标空间，极坐标中的方程为：

ρ=xcosθ+ysinθ

(6)

其中ρ是直线离原法线距离，θ是法线对x轴的角度(如图2所示)。

图2 直角坐标、极坐标对应关系示意

可见直线的Hough变换在极坐标中的一点，而点的Hough变化是正弦函数。由此可知在图像空间中共线的点对应在参数空间中相交的线，参数空间中相交于同一点的所有曲线在图像空间里都有共线的点与之对应。利用相交点的累加统计便可检测直线。

依次利用Matlab图像函数hough(Hough变换)、houghpeaks(提取峰值点)、houghlines(提取峰值点对应的直线)提取图像文件中的直线。然后根据夹角采用imrotate函数进行图像旋转。

图像中若含有纵横网格线，对曲线追踪干扰较大。可先将网格线删除，具体方法为：首先利用Hough变化识别图像中的直线，然后对直线长度进行排序，最后根据给定的纵、横网格线条数删除识别出的长度较长的直线。

3.7 曲线区域特征点获得

图像文件中曲线数量可能为多条，若采用遍历算法寻找图像中曲线，其计算量较大，并且曲线易混淆。因此，提出在手动选择需要提取数据曲线的基础上进行曲线追踪。

MATALB函数getpts可获得用户在当前图像上使用鼠标点击点的坐标，利用此函数便可获得鼠标点击曲线处的像素点位置。为避免点击出现偏差，可以点击像素点位置为圆心逐渐向外寻找曲线像素点。

3.8 曲线追踪识别

获得区域特征点后可利用遍历法寻找曲线像素点，并确定为曲线追踪的起始点。然后利用遍历法在已识别曲线像素周围寻求下一曲线像素并结合曲线斜率不可突变作为判断条件，甄别取舍寻求到的周围像素点进行曲线追踪。曲线斜率可通过综合比较统计方向码、已得曲线拟合计算端点斜率获得。

方向码是边界跟踪中描述行走方向规定的概念，可规定如图3所示：东、东北、北、西北、西、西南、南、东南，分别使用数字0、1、2、3、4、5、6、7表示。统计已得曲线最近几个像素点的方向码，并计算其平均值，便可预测下一像素点的方向码。

图3 方向码示意

曲线拟合方法可利用Matlab多项式拟合函数polyfit对已得曲线像素点进行拟合，然后利用拟合曲线得到端点斜率或下一像素点位置。

3.9 数据换算导出

图像是由像素点构成，像素点位置与实际坐标值(x轴、y轴)存在一一对应关系。因此，将识别出曲线的像素点位置与坐标范围进行比例换算便可得到曲线坐标值。假设曲线图像像素尺寸(长度×宽度)为(A，B)；曲线的像素点位置，行、列值为(a、b)；曲线X轴坐标范围为(Xmin，Xmax)，Y轴坐标范围为(Ymin，Ymax)。则数据换算公式如下。

(7)

最后将换算得到的曲线数据导出。Matlab导出的数据文件可为mat(Matla数据)、Excel等形式。可利用uiputfile函数保存为mat文件，xlswrite函数保存为Excel文件。

4 水泵特性曲线提取实例

根据上述图像处理过程，利用MATLAB的GUI编辑界面，编写曲线数据提取软件包。利用此软件包进行水泵特性曲线提取。编制的软件界面如图4所示。

图4 曲线提取软件GUI界面示意

假设需处理的水泵特性曲线如图5所示，二值化后图像如图6所示，细化并删除坐标网格线的图像如图7所示。图像文件中有多条曲线，可依次选择需提取曲线，得到曲线数据绘制图形如图8所示，其中红色曲线是依据提取出的数据点重新绘制而成。

图5 水泵特性曲线原图像示意

图6 二值化处理后图像示意

图7 细化并删除坐标网格线后示意

图8 读取所有曲线后图像示意

比较原始图片和提取绘制的图像可以得出，提取绘制图像与原图像无偏差。

5 结语

1)利用已识别曲线像素寻求周围曲线像素并结合曲线斜率不可突变作为判断条件，甄别取舍寻求到的周围像素点可保证曲线追踪的连续性。

2)利用MATLAB软件进行图像曲线数据提取，能克服肉眼判断带来的误差，提取数据与原图像一致，大大提高了准确性及提取效率。

3)MATLAB语句简单，GUI软件包可实现软件的迅速调用及推广应用。