大气PM2.5颗粒的图像分析研究

吴晶晶，王 柯，和 森

(1. 西安邮电大学 电子工程学院，陕西 西安 710061；2. 西安邮电大学 通信与信息工程学院，陕西 西安 710061)



大气PM2.5颗粒的图像分析研究

吴晶晶1，王 柯2，和 森1

(1. 西安邮电大学 电子工程学院，陕西 西安 710061；2. 西安邮电大学 通信与信息工程学院，陕西 西安 710061)

基于对空气中PM2.5颗粒群图像特征分析，通过图像处理的方法，将图像中的颗粒进行分割，再采用八连通判别算法识别、标记这些颗粒物，测量出颗粒物的数目及大小，并在Matlab 7.0软件平台上实现了标记颗粒物和测量颗粒物粒径的算法，这种处理与识别方法不受颗粒形状的影响，为PM2.5等微米级颗粒材料的图像分析提供了参考。

PM2.5；图像分割；颗粒标记；粒径测量

科学证实，PM2.5是霾的主要组成部分[1]，是指空气动力学直径≤2.5 μm的大气颗粒物质[2]，由于其富集有毒重金属、有机污染物、细菌和病毒，可通过呼吸道进入肺泡，进而经过血液循环到达其他器官，从而导致人体呼吸系统和其他功能系统的损害，因此已受到全世界普遍关注[3]。本文基于Matlab7.0平台研究了计算颗粒物数量与大小的算法，为控制PM2.5空气污染提供了参考。

1 PM2.5显微图像的预处理与分割

1.1 自适应阈值分割

图像分割可便于计算机自动对图像进行识别和处理。图像分割是指按照一定原则，将图像分成若干个互不相交的小区[4]。阈值法是一种简单有效的图像分割方法，此法用一个或几个阈值将图像的灰度级分为几个部分，将属于同一部分的像素视为相同的物体[5]。

本文采用自适应阀值法对颗粒图像进行分割。该方法将原始颗粒图像分割成若干个子图像，并采用最大类间方差法求取每一个子图像的阀值[6]，再用该阈值分割对应的子图像。采用自适应阈值法分割图像能较好地消除图像光照不均的影响。

建立一个平面、菱形结构元素，使用二值结构要素矩阵SE对原始图像数据矩阵im1执行腐蚀操作，将目标的边缘“毛刺”剔除。使用二值结构要素矩阵SE对腐蚀图像BW2执行膨胀操作，将目标的边缘或者是内部的坑填掉，利用原始图像im1和膨胀图像BW3，将图像转换成灰度图像I_gray。根据灰度图像，统计256个灰度级中每个灰度级的个数，再根据图像灰度大小求出归一化直方图。首先计算出总的平均灰度，然后分别计算前景与背景的平均灰度比例，最后计算前景与背景图像的方差，并将方差值存储。每次计算出方差后，与上一次计算的方差作比较，保留最大方差为 finalT。以最佳阈值finalT为界，将前景像素点赋值为1，背景像素点赋值为0。图2是自适应阈值分割后的二值图像，图1是分割前的二值图像，对照表明经自适应阈值法分割后，图像中重叠的颗粒得到很好的分离，降低了PM2.5颗粒尺寸分布的误差。

图1 分割前二值图像

图2 分割后二值图像

1.2 中值滤波

如图2所示，图像上存在诸多较小的黑点，这属于椒盐噪声，本文选用中值滤波对图像进行滤波处理。中值滤波是一种非线性滤波，可较好地避免线性滤波带来的图像边缘模糊问题，一般用来处理图像的椒盐噪声[7]。中值滤波是将领域内像素灰度的中值代替该像素的值,适用于处理脉冲噪声[8]。

图3是图2采用中值滤波处理后的结果。如图所示，图片中滤除了许多较小的黑点，中值滤波后的二值图像比图2更加清晰。

图3 滤波后图像

1.3 与运算

图像的逻辑运算在一定程度上也可以用来优化图像效果[9]。如图3所示，中值滤波后图像中颗粒物边缘扩大，尤其是图片边界的颗粒物。颗粒物边缘扩大，白色增多(二值图像“1”增多)。利用自适应图像阈值分割后二值图像与中值滤波后的二值图像进行“与”运算将优化图像效果[10]。如图4所示是优化后的图像效果。

图4 与运算优化后图像

1.4 Canny边缘检测

步骤1 去除噪声。将原始图像数据与高斯平滑模板进行卷积处理，对比原始图像后发现，经过卷积的图像变得轻微模糊。这样，经过高斯平滑模板处理后的图像将不再受到单独一个像素噪声的影响[11]。

步骤2 寻找图像中的亮度梯度。针对图像边缘可能会指向不同方向的这种情况，Canny 算法使用 4 个 mask对水平、垂直以及对角线方向的边缘分别进行检测[12]。将原始图像与4个 mask分别进行卷积处理并存储。对于每个点都标识在这个点上的最大值以及生成的边缘方向[12]。这样就从原始图像中生成了每个点亮度梯度图以及亮度梯度的方向[12]。

步骤3 在图像中跟踪边缘。在一定程度上，边缘的亮度梯度有可能是较高的亮度梯度，但不确定多高的亮度梯度可作为一个门限来判定是否为边缘，所以 Canny算法使用了滞后阈值。滞后阈值需要高低两个阈值。假定图像中的重要边缘由连续曲线组成，在避免将噪声像素当成边缘的前提下跟踪给定曲线的模糊部分[13]。于是可从其中一个较大的阈值开始，使用前面导出的方向信息，标识出比较确信的真实边缘，从这些真正的边缘开始在图像中跟踪整个边缘[14]。再使用一个较小的阈值跟踪曲线的模糊部分直到回到起点[14]。当完成此过程后，便得到一个二值图像，其中每个点都表示一个边缘点。

图5 边缘检测后的图像

1.5 孔洞填充

利用自动阀值选择法对图像进行Canny算子检测，提取边界，但处理后有些PM2.5颗粒内部存在孔隙，可利用imfill函数对这些孔洞进行填充，避免在计算颗粒数目过程中产生误差。

如图3所示，分割后的图像黑色的背景上有众多白色颗粒，而有些白色颗粒内有黑色的孔洞，本文采用孔洞填充来填补二值图像中的孔洞区域。图5所示是图3经过孔洞填补后的结果。

1.6 开操作

采用开操作分开图像中粘连在一起的颗粒，首先将粘连颗粒“腐蚀”，可断开颗粒间狭窄的连接部分并消除细的突出物，这样不仅可使颗粒轮廓光滑，且也适当缩小颗粒尺寸；然后对经过“腐蚀”操作后的颗粒再进行“膨胀”操作，恢复颗粒的大小到原来的尺寸。但是开操作具有一定的局限性，只能分开连接长度较小的颗粒群，对于分离连接长度较大的颗粒群以及直径小于或等于结构元素大小的颗粒都不适用于此操作[15]。

本文针对采用开操作无法分离的个别颗粒，起始阶段使用了自适应阀值分割法来处理。图6显示了对图5进行开操作的结果。

图6 开操作后图像

2 颗粒数目

图6中相互连接的白色象素的集合表示一个颗粒，颗粒的标记过程即给颗粒标上记号。八连通判别算法是指在由“0”象素(颗粒信息)和“1”象素(背景信息)组成的点阵图像的二值图像中，把颗粒信息从相互邻接的8个区域的“0”值象素集合中提取出来。

本程序采用八连通判别算法，通过对图像每个联通的白色区域进行标记操作，从而求得颗粒数目[15]。如衅7所示，设有4个不联通的颗粒，0表示背景，1,2,3,4表示颗粒[15]。在算法中，定义count记录颗粒数目，定义数组arr(count)表示第count个颗粒内部包含的象素个数。如第1个颗粒，arr(1)=5；第2个颗粒，arr(2)=5。下面以标记颗粒1为例，说明计算颗粒数目和面积的具体算法。

(1)初始化count=0；从左到右、上到下逐个扫描图像；若没有遇到颗粒，则继续扫描。若遇到颗粒，如(1,1)，则依次判断像素点的左方(1,0)，左上方(0,0)，上方(0,1)，右上方(0,2)，右方(1,2)，右下方(2,2)，下方(2,1)，左下方(2,0)是否为颗粒，若均不为颗粒，则count值加1；

(2)若左方点为颗粒，则当前点跟左方点均为白点，白点入队列，并标记为已访问。如：当前点(1,2)，其左方点(1,1)为颗粒，所以当前点(1,2)与左方点(1,1)相同的值，即count相同，arr(count)的值加1；

(3)若左方点不为颗粒，则判断左上方点。如当前点(2,4)，其左方点(2,3)不为颗粒，而左上方点(1,3)为颗粒，所以当前点(2,4)与左上方点(1,3)相同的值，即count相同，但arr(count)的值加1；

(4)若左上方点不为颗粒，则判断上方点。如当前点(5,6)，其左方点(5,5)和左上方点(4,5)均不为颗粒，而上方点(4,6)为颗粒，所以当前点(5,6)与正上方点相同点值，即count相同，但arr(count)的值加1；

(5)同理，若当前点的左方点、左上方点、正上方点均不为颗粒，则以同样的方法一次判断右上方点，若右上方点也不为颗粒，则判断右方点；

(6)由于当前点(1,7)的左方点(1,6)，左上方点(0,6)，正上方点(0,7)，右上方点(0,8)均不为颗粒，则当前的值count加1；

在具体实现中，首先从左到右，从上到下，依次检测每个象素，如果发现某象素值为1，则依次检测该点的左方点、左上方点、正上方点、右上方点，根据前面介绍的连通性判断，记录颗粒数目，最后count的值就是统计的颗粒的总个数，arr(count)就是每个颗粒的象素总数。对图6计算出的颗粒总个数是302，其中PM2.5颗粒总数是247。

0列1列2列3列4列5列6列7列8列0行0000330001行0110330402行0110300443行0010000444行0000004445行222200444

图7 八连通判别法示意图

3 颗粒大小

本文研究颗粒大小时，先处理原始图像，获得二值图像，二值图象中，相互连接的白色象素的集合成为一个颗粒。颗粒的边界内部有众多像素，通过这些像素来计算颗粒的大小。本文先逐行逐列地遍历寻找首个前景(白点)p，标记为1，并调用连通域函数,将该点传递给形参，返回值为每个domain中包含的总象素。调用函数中首个白点入队列，以该点为中心，该点左方的点为出发点，顺时针方向遍历寻找白点，并将寻找到的白点入队列，对队列中的点继续进行查找。同时将白点周围领域已访问的点标上标记，避免重复访问，降低效率。当该点领域的点全部被访问到时,该点出队列，对应的指针后移,直到指针指向空为止，将该域的象素总数count赋值给返回值c,c为该域即该颗粒的大小。

在实际测量中，往往需要对颗粒进行实际尺寸的测量。利用了Nikon公司提供的标称值为0.01 mm的测微尺，使用与图1～图6相同的放大倍数[15]，计算得出1 μm对应3.015个象素。因此，通过转换，可得颗粒的实际面积，把颗粒近似看成圆形，可计算出其的等效直径。从而在测出的所有颗粒中，利用PM2.5颗粒直径的限制，可以统计出PM2.5颗粒的数目，再统计出PM2.5颗粒面积和粒径大小，在实验结果中将显示其运行的PM2.5颗粒尺寸分布结果。

3 实验结果

在Matlab 7.0软件平台上运行程序，在PM2.5颗粒显微图像预处理后得到图1～图6所示图像，充分将重叠颗粒进行分离，并测得颗粒物的总个数是302颗，每个颗粒物的大小如图8所示(单位：象素)。而又因为PM2.5颗粒的直径是≤2.5 μm的，利用直径大小，可从302粒颗粒中筛选出247个PM2.5颗粒。在Matlab7.0软件平台上运行程序，为更客观地看到颗粒数目及大小，本文截取5张运行结果的图片。如图8～图12所示，显示了247颗PM2.5颗粒的粒径大小及颗粒面积大小。

图8 颗粒总数及大小

图9 PM2.5颗粒数目

图10 PM2.5颗粒粒径

图11 PM2.5颗粒大小

图12 PM2.5颗粒大小

4 结束语

对上述颗粒物图像的处理过程是在Matlab 7.0平台上进行的，共处理了302个颗粒,每个颗粒的大小存储在数组arr(count)中，其中247个PM2.5颗粒，并将每个PM2.5颗粒的大小存储在数组arry (l)中，将每个颗粒的粒径大小存储在arrys(l)。由此可见，本文颗粒处理、识别、颗粒数目计数和颗粒面积计算等均是有效的。本文的程序适用于分析与处理微米级颗粒物的尺寸分布，为控制PM2.5颗粒污染提供了理论基础。

[1] 谭译,陈默.城市雾霾治理政策措施探讨[J].环球人文地理,2014(2):88-90.

[2] 叶秋霞.颗粒物PM10和PM2.5对体育锻炼和身体健康的影响[J].河南农业,2014(8):57-58.

[3] 杨洪斌,邹旭东,汪宏宇,等.大气环境中PM2.5的研究进展与展望[J].气象与环境学报,2012,28(3):77-82.

[4] 盛仲飙.基于对比度增强的分水岭分割算法[J].计算机技术与发展,2012(12):90-92,96.

[5] 刘超,蔡文华,陆玲.图像阈值法分割综述[J].电脑知识与技术,2015(1):140-142,145.

[6] 黄丽丽,杨帆,王东强,等.基于改进型最大类间方差法的瞳孔定位方法[J].计算机工程与应用,2013(23):137-140.

[7] 杨明,李晶.一种新的椒盐噪声去除方法[J]. 科技视界,2014(25):111.

[8] 王梅,黄华,应大力.基于隶属度函数的自适应加权中值滤波[J].电脑知识与技术,2012(31):7565-7567.

[9] 樊雅平,黄生学.一种新的基于像素交迭的图像插值方法[J].广西师范学院学报:自然科学版,2003,20(4):72-76.

[10] 蒋爱如.基于Matlab的激光图像背景噪声的处理[J].大众科技,2012(11):28-29.

[11] 关榆君,杨帆.图像处理在管道机器人焊缝自动定位中的研究[J].工业控制计算机,2012,25(12):66-67.

[12] 李武.多景近景摄影测量的特征匹配[J]. 北京测绘,2013(5):46-49.

[13] 安然,王辉麟,史天运.基于智能视频分析的铁路弓网状态监测技术研究[C].合肥:第八届中国智能交通年会,2013.

[14] 施诚诚.基于视频流的弹点捕捉和跟踪[D].上海:东华大学,2014.

[15] 李仪芳,刘景琳.大气PM2.5显微图像粒径大小的测量[J].佛山科学技术学院学报:自然科学版,2008,26(2):68-70.

Image Analysis of Atmospheric PM2.5 Particles

WU Jingjing1，WANG Ke2，HE Sen1

(1. School of Electronic Engineering, Xi’an University of Posts and Telecommunications, Xi’an 710061, China; 2. School of Communication and Information Engineering, Xi’an University of Posts and Telecommunications,Xi’an 710061, China)

Based on the analysis of PM2.5 particle group image feature in the air,through pictures processing.Use eight communication discrimination algorithm to identify,so as to label particles,measure the number and size of the particles and achieve marking particles and measuring particle diameter of the algorithm in the Matlab7.0 software platform,such kind of processing and recognition is not affected by particle shape,and it provides reference for the image analysis of PM2.5 micron-sized particles.

PM2.5; image segmentation; particle marker; particle size distribution measurement

10.16180/j.cnki.issn1007-7820.2016.12.045

2016- 03- 01

吴晶晶(1990-)，女，硕士研究生。研究方向：图像处理。

TP391.41；TP312

A

1007-7820(2016)12-162-05