CN107044947B - 一种基于图像特征的pm2.5污染指数的识别方法 - Google Patents

一种基于图像特征的pm2.5污染指数的识别方法 Download PDF

Info

Publication number
CN107044947B
CN107044947B CN201710301867.4A CN201710301867A CN107044947B CN 107044947 B CN107044947 B CN 107044947B CN 201710301867 A CN201710301867 A CN 201710301867A CN 107044947 B CN107044947 B CN 107044947B
Authority
CN
China
Prior art keywords
image
value
pixel
weighted average
queue
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.)
Active
Application number
CN201710301867.4A
Other languages
English (en)
Other versions
CN107044947A (zh
Inventor
白鹤翔
李艳红
李德玉
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
BEIJING JIAHUA ZHILIAN TECHNOLOGY CO.,LTD.
Original Assignee
Shanxi University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Shanxi University filed Critical Shanxi University
Priority to CN201710301867.4A priority Critical patent/CN107044947B/zh
Publication of CN107044947A publication Critical patent/CN107044947A/zh
Application granted granted Critical
Publication of CN107044947B publication Critical patent/CN107044947B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume or surface-area of porous materials
    • G01N15/06—Investigating concentration of particle suspensions
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/08—Learning methods
    • G06N3/084—Backpropagation, e.g. using gradient descent
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V20/00—Scenes; Scene-specific elements
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume or surface-area of porous materials
    • G01N15/06—Investigating concentration of particle suspensions
    • G01N15/075—Investigating concentration of particle suspensions by optical means

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Chemical & Material Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Artificial Intelligence (AREA)
  • Evolutionary Computation (AREA)
  • Pathology (AREA)
  • Biochemistry (AREA)
  • Analytical Chemistry (AREA)
  • Biomedical Technology (AREA)
  • Biophysics (AREA)
  • Computational Linguistics (AREA)
  • Data Mining & Analysis (AREA)
  • Immunology (AREA)
  • Molecular Biology (AREA)
  • Computing Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Software Systems (AREA)
  • Dispersion Chemistry (AREA)
  • Multimedia (AREA)
  • Image Analysis (AREA)
  • Image Processing (AREA)

Abstract

本发明属于图像处理及气象学交叉领域，具体涉及一种基于图像特征的PM2.5污染指数的识别方法，主要解决了现有PM2.5的测定方法存在操作复杂和测定结果不精确的技术问题，本发明将标注了PM2.5实测值的历史图像统一缩放为XS×YS像素大小并转换为灰度图；对每个灰度图计算其图像特征；将图像特征作为输入使用合适的回归方法得到PM2.5的预测模型；将目标图像全部缩放为XS×YS像素大小并转换为灰度图；针对缩放后的目标图像计算其各种特征；将目标图像的图像特征作为已训练好预测模型的输入计算得到目标图像PM2.5值。本发明通过对一幅图像的特征进行分析来估计出图像所描绘场景的PM2.5值，从而可以应用到日常的出行、健康、天气预测等人类活动，并可用于大规模PM2.5污染指数的快速测量。

Description

一种基于图像特征的PM2.5污染指数的识别方法
技术领域
本发明属于图像处理及气象学交叉领域，具体涉及一种基于图像特征的PM2.5污染指数的识别方法。
背景技术
空气中漂浮着的颗粒物大小不同，PM2.5是其中较细小的那部分。测定PM2.5的浓度一般分为两步：(1)把PM2.5与较大的颗粒物分离，几乎所有的测定方法都需要这一步；(2)测定分离出来的PM2.5的重量。目前，各国环保部门广泛采用的PM2.5测定方法有三种：重量法、β射线吸收法和微量振荡天平法。
PM2.5手动监测方法又称标准称重法(滤膜称重法)属于手工监测方法，主要用于研究，或对自动分析方法的校准。这种方法较简便，只需要一个PM2.5切割头、一台泵和膜架及其滤膜，采集24小时样品后，取下滤膜称重即可，必要时可以平行采集3个样品，经恒温恒湿后再称重。
滤膜称重法的优点是经济成本低，容易实施。缺点：(1)气流长时间不断通过采样滤膜，滤膜上采集到的物质随着气流和温度的变化也会造成挥发性和半挥发性物质的损失；(2)一些极细小的颗粒穿过滤膜，造成结果偏低；(3)气态物质可能被滤膜吸附，造成结果偏高。
自动监测方法的测站点无人值守，数据通过网络直接获得，包括β射线吸收法和微量振荡天平法。β射线吸收法的基本原理是利用堆积在石英滤膜上的颗粒物对C14释放的β射线衰减量的变化检测大气颗粒物质量的变化。环境空气由采样泵经切割器吸入采样管，经过滤膜后排出。颗粒物沉淀在条状石英滤膜上，当β射线通过沉积着颗粒物的滤膜时，β射线强度发生衰减，通过对衰减量的测定计算出颗粒物的浓度。微量振荡天平(TEOM)法又称微量石英振荡天平法，该方法是在质量传感器内使用一个石英空心锥形管，在空心锥形管振荡端上安放可更换的滤膜，振荡频率取决于石英锥形管特性和它的质量。当采样气流通过滤膜，其中的颗粒物沉积在滤膜上，滤膜质量变化导致振荡频率变化，通过测量振荡频率的变化计算出沉积在滤膜上颗粒物的质量，再根据采样流量、采样现场环境温度和气压计算出该时段的颗粒物标态质量浓度。
振荡天平法的优点是定量关系明确。但其存在下面两个缺点：(1)样品加热后挥发性和半挥发性物质的损失，导致测定结果偏低；(2)需要加装膜动态测量系统(FilterDynamic Measurement System，)对偏低的结果进行校准。
β射线法是基于两个假设：仪器的石英采样滤膜条带均一；采集下来的PM2.5粒子物理特性均一，对β射线强度衰减率相同。但是现实条件下，以上两个假设通常很难成立，因此测定数据一般被认为也存在偏差。这种方法的缺点是在相对干净和干燥的地区故障率低，而在潮湿高温区域故障率很高。
发明内容
本发明的目的是解决现有PM2.5的测定方法存在操作复杂和测定结果不精确的技术问题，提供一种基于图像特征的PM2.5污染指数的识别方法。
为解决上述技术问题，本发明采用的技术方案为：
一种基于图像特征的PM2.5污染指数的识别方法，包括以下步骤：
步骤1.搜集标注了PM2.5实测值的历史图像，将搜集的历史图像全部统一缩放为XS×YS像素大小；如果图像是彩色的，那么就将其转换为灰度图；
步骤2.计算每个缩放并转换为灰度图后图像的图像特征；
步骤3.将每个图像的特征作为输入，使用回归方法得到PM2.5的预测模型；
步骤4.采集目标图像，按照步骤1中同样的方法将目标图像缩放为XS×YS像素大小，如果图像是彩色的，将其转换为灰度图；
步骤5.针对缩放后的目标图像，按照步骤2中的同样的方法，计算缩放并转换为灰度图后图像的图像特征；
步骤6.将目标图像的图像特征作为步骤3中已经训练好的预测模型的输入，计算得到目标图像的PM2.5值。
所述步骤2中在计算每个缩放并转换为灰度图后图像的图像特征之后,还需要对每个图像的图像特征以及PM2.5的值进行归一化；所述步骤5中针对缩放后的目标图像，在计算缩放并转换为灰度图后图像的图像特征之后，还需要对目标图像的图像特征按照步骤2中同样的归一化方法，进行归一化。在机器学习中，使用同训练时一样的方法对数据进行预处理，即使用同训练时一样的方法将图像转换为灰度值，计算图像特征并进行归一化，在使用先前训练好的预测模型进行预测时才能得到最佳的预测结果。
所述步骤1或步骤4中采用公式GRAY＝0.299×R+0.587×G+0.114×B将彩色图像转换为灰度图像，其中R、G、B分别表示彩色图像的红色、绿色和蓝色波段；
所述步骤2或步骤5中的图像特征为图像的平均像素值、局部极值点的个数NMax、局部方差的均值μLVAR、图像水平方向和垂直方向二阶差分加权平均值、图像水平方向和垂直方向的二阶梯度和的加权平均值、Qx、Qy和Qg折线差分的加权平均值以及图像中的斑点数；之所以选择这些特征，是因为在实验过程中发现这些特征能够很好的反映实际拍摄图像中PM2.5高低。
所述图像的平均像素值为所有像元的像素值加起来之后除以XS×YS；
所述图像局部极值点的个数NMax和局部方差的均值μLVAR的计算方法为：
(1)设定NMax和μLVAR的初值为0；
(2)对于图像中的每一个非边界像元找出邻近的8个像元，记8个像元的构成集合为PIXA；如果当前像元坐标为(i,j),i和j分别表示该像元和图像最左上角像元分别在垂直和水平方向同当前像元间隔的像元个数，那么其8个近邻的坐标分别为(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)；这里的非边界像元是指图像中所有能够找到上述8个近邻，并且任一近邻像元像素值不为空的像元；
(3)如果当前像元值大于周边的8个像元，那么NMax值加1；并且
这里，p表示当前像元，I(p)表示像元p的灰度值，表示PIXA中的8个像元和当前像元像素值的均值；
(4)当所有像元都处理完之后，μLVAR＝μLVAR/((XS-2)×(YS-2))；
所述图像水平方向和垂直方向二阶差分加权平均值、图像水平方向和垂直方向二阶梯度和的加权平均值的计算方法为：
(1)分别计算图像水平方向x和垂直方向y每个像元p的二阶差分,分别记作和并分别将两个方向的二阶差分按照从小到大的序列排序，排序过程中去掉重复值，并过滤掉小于水平方向给定阈值αx或垂直方向给定阈值αy的所有差分值，由此得到两个大小分别为nx和ny的水平方向队列Qx和垂直方向队列Qy，两个队列中第i的元素分别为和即为队列中排序为i的像元所对应的二阶差分值；
(2)计算图像每个像元p的二阶梯度和并将二阶梯度和值按照从小到大的序列排序，排序过程中去掉重复值，并过滤掉小于给定阈值αg的所有二阶梯度和值，由此得到大小为ng的二阶梯度和队列Qg，队列中第i的元素分别为Qg(i)＝grad(pi)，即为队列中排序为i的像元所对应的二阶梯度和值；
(3)分别对三个队列Qx、Qy和Qg计算对应折线总长度和以及总体加权变化度 和
(4)分别计算三个队列Qx、Qy和Qg的加权平均值WVLx＝Vx/Lx、WVLy＝Vy/Ly和WVLg＝Vg/Lg，这三个值就分别为水平方向二阶差分加权平均值、垂直方向二阶差分加权平均值以及二阶梯度和的加权平均值；
所述图像中的斑点数的计算方法为：
(1)给定的一系列有序尺度值σ1,σ2,…,σs,s≥3，并选定某个半径r，根据规范化拉普拉斯高斯算子生成对应各个尺度的卷积模板，其中每个卷积模板的中心坐标为(0，0)，所有和卷积模板中心水平或者垂直方向相差r个像元的所有网格都属于模板；模板中所有网格的坐标集合为{(tx,ty)|tx＝-r,-r+1,…,r-1,r∧ty＝-r,-r+1,…,r-1,r}，对应尺度σ∈{σ1,σ2,…,σs}的卷积模板中每个网格(tx,ty)上的值使用下面公式计算：
对每个尺度σ∈{σ1,σ2,…,σs}，使用对应的卷积模板对图像做卷积运算，得到一系列卷积图像I1,I2,…,Is；
(2)对于每个卷积图像Ii,i＝2,…,s-1中的每个像元p，假设其在图像中的坐标为(i,j)，其8个最近邻的像元坐标分别为(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)；如果Ii(p)是所有近邻中最小的，并且Ii(p)也比Ii-1和Ii+1中所有坐标为(i,j)、(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)的像元的值都小，那么该像元对应的就是一个斑点；如果Ii(p)是所有近邻中最大的，并且Ii(p)也比Ii-1和Ii+1中所有坐标为(i,j)、(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)的像元的值都大，那么该像元对应的也是一个斑点；在所有卷积图像Ii,i＝2,…,s-1上找到的斑点总数即为图像中的斑点数。
所述Qx、Qy和Qg折线差分的加权平均值的计算方法为：
(1)首先针对Qx、Qy和Qg的折线分别计算相邻两个元素之间的差值，并且将这一差值按照从小到大排序，分别形成一个队列，分别记为QQx、QQy和QQg，其第i个元素分别为QQx(i)、QQy(i)和QQg(i)。
(2)分别对三个队列QQx、QQy和QQg计算对应折线总长度和以及总体加权变化度 和
(3)分别计算三个队列QQx、QQy和QQg的加权平均值WVLQx＝VQx/LQx、WVLQy＝VQy/LQy和WVLQg＝VQg/LQg，这三个值就分别为Qx、Qy和Qg折线差分的加权平均值。
所述步骤2和步骤5中按照如下公式对图像的特征和PM2.5的值进行归一化：
其中F表示某个特征，F(I)表示缩放并转换为灰度图后的历史图像或者目标图像I在特征F上的取值，minF表示特征F在所有缩放并转换为灰度图后的历史图像上的最小值，maxF表示特征F在所有缩放并转换为灰度图后的历史图像上的最大值,Fnorm(I)表示I在特征F上归一化后的值，并且在对转换为灰度图后的目标图像进行归一化时若Fnorm(I)＞1时将Fnorm(I)置为1，若Fnorm(I)＜0时将Fnorm(I)置为0,Inew为待归一化的新值。这一方法实现起来比较简单高效，因此选择这一方法进行归一化。
所述步骤3中使用回归方法得到PM2.5的预测模型的方法为前馈神经网络、高斯过程或常用的最小二乘等方法。
所述采用前馈神经网络方法得到PM2.5的预测模型的具体操作步骤为：
(1)在前馈神经网络中设定好激活函数，这一步可以使用逻辑斯特函数或者双曲正切函数；
(2)设定神经网络的层次和每层的神经元数量；
(3)使用反向传播算法学习神经网络；
(4)将学习到的网络保存下来，需要保存的内容有使用了何种激活函数，神经网络的层数、每层的神经元数量以及每个神经元的权重；
所述步骤6中采用前馈神经网络方法得到PM2.5的预测模型时，需要读入保存的基于前馈神经网络的预测模型参数，包括激活函数，神经网络的层数、每层的神经元数量以及每个神经元的权重，然后按照前馈神经网络的前馈传播方法将归一化后的图像特征作为输入计算出Inew所对应的归一化后的PM2.5值PM2.5norm(Inew)，最后将其按照下面公式反归一化后作为新拍摄图像的PM2.5值PM2.5(Inew)：
PM2.5(Inew)＝PM2.5norm(Inew)×(maxF-minF)+minF。
本发明充分利用历史PM2.5的测量值和对应的实景照片，综合使用图像处理技术和数据挖掘技术，只需要使用目标区域景物照片就能够快速估计其PM2.5指数，克服了传统方法通常只能定点测量的缺点，提供了一种能够快速估计大范围区域PM2.5的技术方法。本发明采用以上技术方案，与现有技术相比的优点在于：
(1)不需要制作或者使用特别的仪器以及化学试剂，仅使用目标区域的照片就能够识别PM2.5指数；
(2)当只有目标区域照片时，无法重建或者重复观察目标区域的环境因素，本方法仍然可以识别其PM2.5指数；
(3)能够直接识别大面积区域的平均PM2.5值，而不是定点观测。
附图说明
图1为本发明的实现流程图；
图2-5为本发明实例中的四个输入图像。
具体实施方式
本实施例中的一种基于图像特征的PM2.5污染指数的识别方法，包括以下步骤：
步骤1.选择一个区域，使用摄像头拍摄不同时间的彩色照片并记录这一时间的PM2.5值以搜集标注了PM2.5实测值的历史图像，根据季节和时间将所有的照片区分为白天和晚上两大类，下面计算都以白天为例；将搜集的历史图像全部统一缩放为XS×YS像素大小，这里使用imagemagick软件的convert命令：“convert–resize 640x320”将搜集的历史图像全部统一缩放为640×320像素大小；这里使用只使用彩色照片的红色波段的方式将图像转换为灰度图，也可以采用公式GRAY＝0.299×R+0.587×G+0.114×B将图像转换为灰度图，其中R、G、B分别表示彩色图像的红色、绿色和蓝色波段；
步骤2.计算每个缩放并转换为灰度图后图像的平均像素值、局部极值点的个数NMax、局部方差的均值μLVAR、图像水平方向和垂直方向二阶差分加权平均值、图像水平方向和垂直方向的二阶梯度和的加权平均值、Qx、Qy和Qg折线差分的加权平均值以及图像中的斑点数；
所述图像的平均像素值为所有像元的像素值加起来之后除以XS×YS；
所述图像局部极值点的个数NMax和局部方差的均值μLVAR的计算方法为：
(1)设定NMax和μLVAR的初值为0；
(2)对于图像中的每一个非边界像元找出邻近的8个像元，记8个像元的构成集合为PIXA；如果当前像元坐标为(i,j),i和j分别表示该像元和图像最左上角像元分别在垂直和水平方向同当前像元间隔的像元个数，那么其8个近邻的坐标分别为(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)；这里的非边界像元是指图像中所有能够找到上述8个近邻，并且任一近邻像元像素值不为空的像元；
(3)如果当前像元值大于周边的8个像元，那么NMax值加1；并且
这里，p表示当前像元，I(p)表示像元p的灰度值，表示PIXA中的8个像元和当前像元像素值的均值；
(4)当所有像元都处理完之后，μLVAR＝μLVAR/((XS-2)×(YS-2))；
所述图像水平方向和垂直方向二阶差分加权平均值、图像水平方向和垂直方向二阶梯度和的加权平均值的计算方法为：
(1)分别计算图像水平方向x和垂直方向y每个像元p的二阶差分,分别记作和并分别将两个方向的二阶差分按照从小到大的序列排序，排序过程中去掉重复值，并过滤掉小于水平方向给定阈值αx或垂直方向给定阈值αy的所有差分值，由此得到两个大小分别为nx和ny的水平方向队列Qx和垂直方向队列Qy，两个队列中第i的元素分别为和即为队列中排序为i的像元所对应的二阶差分值；
(2)计算图像每个像元p的二阶梯度和并将二阶梯度和值按照从小到大的序列排序，排序过程中去掉重复值，并过滤掉小于给定阈值αg的所有二阶梯度和值，由此得到大小为ng的二阶梯度和队列Qg，队列中第i的元素分别为Qg(i)＝grad(pi)，即为队列中排序为i的像元所对应的二阶梯度和值；
(3)分别对三个队列Qx、Qy和Qg计算对应折线总长度和以及总体加权变化度 和
(4)分别计算三个队列Qx、Qy和Qg的加权平均值WVLx＝Vx/Lx、WVLy＝Vy/Ly和WVLg＝Vg/Lg，这三个值就分别为水平方向二阶差分加权平均值、垂直方向二阶差分加权平均值以及二阶梯度和的加权平均值；
所述Qx、Qy和Qg折线差分的加权平均值的计算方法为：
(1)首先针对Qx、Qy和Qg的折线分别计算相邻两个元素之间的差值，并且将这一差值按照从小到大排序，分别形成一个队列，分别记为QQx、QQy和QQg，其第i个元素分别为QQx(i)、QQy(i)和QQg(i)。
(2)分别对三个队列QQx、QQy和QQg计算对应折线总长度和以及总体加权变化度 和
(3)分别计算三个队列QQx、QQy和QQg的加权平均值WVLQx＝VQx/LQx、WVLQy＝VQy/LQy和WVLQg＝VQg/LQg，这三个值就分别为Qx、Qy和Qg折线差分的加权平均值。
所述图像中的斑点数的计算方法为：
(1)给定的一系列有序尺度值σ1,σ2,…,σs,s≥3，并选定某个半径r，根据规范化拉普拉斯高斯算子生成对应各个尺度的卷积模板，其中每个卷积模板的中心坐标为(0，0)，所有和卷积模板中心水平或者垂直方向相差r个像元的所有网格都属于模板；模板中所有网格的坐标集合为{(tx，ty)|tx＝-r，-r+1，…，r-1，r∧ty＝-r，-r+1，…，r-1，r}，对应尺度σ∈{σ1,σ2,…,σs}的卷积模板中每个网格(tx，ty)上的值使用下面公式计算：
对每个尺度σ∈{σ1,σ2,…,σs}，使用对应的卷积模板对图像做卷积运算，得到一系列卷积图像I1,I2,…,Is；
(2)对于每个卷积图像Ii,i＝2,…,s-1中的每个像元p，假设其在图像中的坐标为(i,j)，其8个最近邻的像元坐标分别为(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)；如果Ii(p)是所有近邻中最小的，并且Ii(p)也比Ii-1和Ii+1中所有坐标为(i,j)、(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)的像元的值都小，那么该像元对应的就是一个斑点；如果Ii(p)是所有近邻中最大的，并且Ii(p)也比Ii-1和Ii+1中所有坐标为(i,j)、(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)的像元的值都大，那么该像元对应的也是一个斑点；在所有卷积图像Ii,i＝2,…,s-1上找到的斑点总数即为图像中的斑点数；
本实例选择了计算局部极值点的个数NMax、局部方差的均值μLVAR、平均像素值，基于拉普拉斯算子的二阶梯度和加权平均值和Qg折线差分的加权平均值、基于Sobel算子的Qy折线差分的加权平均值以及基于高斯拉普拉斯算子的图像中的斑点数这些特征，计算出的图像部分特征如下表所示：
计算高斯拉普拉斯图像斑点数时使用了OpenCV库，其部分代码如下：
Mat image＝img；
Mat descriptors；
vector<KeyPoint>keypoints；
SimpleBlobDetector::Params params；
params.minThreshold＝10；
params.maxThreshold＝100；
params.thresholdStep＝10；
params.minArea＝10；
params.minConvexity＝0.3；
params.minInertiaRatio＝0.01；
params.maxArea＝8000；
params.maxConvexity＝10；
params.filterByColor＝false；
params.filterByCircularity＝false；
SimpleBlobDetector blobDetector(params)；
blobDetector.create("SimpleBlob")；
blobDetector.detect(image,keypoints)；
drawKeypoints(image,keypoints,image,Scalar(255,0,0))；
cout<<keypoints.size()<<""；
计算基于拉普拉斯算子的二阶梯度和加权平均值和Qg折线差分的加权平均值也使用了OpenCV库，其部分代码如下：
这里αg阈值设定为50。
并对每个图像的各图像特征以及PM2.5的值按照如下公式进行归一化：
其中F表示某个特征，F(I)表示缩放并转换为灰度图后的历史图像或者目标图像I在特征F上的取值，minF表示特征F在所有缩放并转换为灰度图后的历史图像上的最小值，maxF表示特征F在所有缩放并转换为灰度图后的历史图像上的最大值,Fnorm(I)表示I在特征F上归一化后的值，并且在对转换为灰度图后的目标图像进行归一化时若Fnorm(I)＞1时将Fnorm(I)置为1，若Fnorm(I)＜0时将Fnorm(I)置为0,Inew为待归一化的新值。以平均像素值为例，训练数据中这一特征的最大值为147.323，最小值为49.6862，因此，当这一值为113.21时，其归一化后的值为(113.21-49.6862)/(147.323-49.6862)＝0.181324。下表为各个特征和PM2.5归一化之后的结果：
步骤3.将每个图像的特征作为输入，使用前馈神经网络方法得到PM2.5的预测模型，其具体操作步骤为：
(1)在前馈神经网络中设定好激活函数，这一步可以使用逻辑斯特函数或者双曲正切函数，这里设置激活函数为逻辑斯特函数，也叫做sigmoid函数；
(2)设定神经网络的层次和每层的神经元数量，这里设置网络层数为3，其中输入层神经元个数为7,输出层神经元个数为1，中间层神经元个数为1000；
(3)使用反向传播算法学习神经网络，这里使用的参数最低训练速率为0.001，最高训练速率为0.01，误差区间为[0.001,0.02],迭代次数为10000；
神经网络采用了lwneuralnetplus 0.998版本的库，部分训练代码如下：
iomanager＝new iomanagelwnnfann()；
iomanager->info_from_file(argv[1],&npatterns,&ninput,&noutput)；
net＝new network(network::LOGISTIC,3,ninput,1000,noutput)；
net->set_momentum(0)；
net->set_learning_rate(0.001)；
net->jolt(0.03,0.22)；
train＝new trainer(net,"","")；
train->set_iomanager(iomanager)；
cout<<"Loading training data..."<<endl；
train->load_training(argv[1])；
train_go(argv,train,net)；
其中train_go部分代码如下所示：
(4)将学习到的网络保存下来，需要保存的内容有使用了何种激活函数，神经网络的层数、每层的神经元数量以及每个神经元的权重；
步骤4.采集目标图像，按照步骤1中同样的方法将目标图像缩放为XS×YS像素大小，如果图像是彩色的，将其转换为灰度图；
步骤5.针对缩放后的目标图像，按照步骤2中的同样的方法，计算缩放并转换为灰度图后图像的平均像素值、局部极值点的个数NMax、局部方差的均值μLVAR、图像水平方向和垂直方向二阶差分加权平均值、图像水平方向和垂直方向的二阶梯度和的加权平均值、Qx、Qy和Qg折线差分的加权平均值以及图像中的斑点数，并对目标图像的各图像特征按照步骤2中同样的归一化方法，进行归一化；
步骤6.将目标图像的图像特征作为步骤3中已经训练好的预测模型的输入，读入保存的基于前馈神经网络的预测模型参数，包括激活函数，神经网络的层数、每层的神经元数量以及每个神经元的权重，然后按照前馈神经网络的前馈传播方法将归一化后的图像特征作为输入计算出Inew所对应的归一化后的PM2.5值PM2.5norm(Inew)，最后将其按照下面公式反归一化后作为新拍摄图像的PM2.5值PM2.5(Inew)：
PM2.5(Inew)＝PM2.5norm(Inew)×(maxF-minF)+minF。
我们选用了下图中的四个图像作为新拍摄图像进行识别。这里需要先对下图中的各个图像的特征进行归一化，然后带入先前训练好的神经网络计算归一化后的PM2.5值，最后进行反归一化即可。识别效果如下表所示：
实测值 预测值
图像1 37 49
图像2 223 247
图像3 323 306
图像4 117 103
所述步骤2和步骤5中不论使用何种归一化方法，只要能够将特征值或者PM2.5值变为[0,1]区间或者(0,1)区间的数，就属于本专利权利要求范围；
所述步骤3中的预测模型的可选建立方法包括神经网络、高斯过程以及常用的最小二乘等方法等，不论使用何种方法，只要是在预测过程中使用了神经网络、高斯过程以及常用的最小二乘等方法中的一种或者多种特征都在本专利的权利要求范围；
所述步骤2和步骤5中，计算图像水平方向和垂直方向二阶差分加权平均值以及二阶梯度和的加权平均值可以使用多种方法，重复计算并获得多种特征，例如分别使用Sobel算子和拉普拉斯算子计算这三个特征；
所述步骤2和步骤5的归一化为可选步骤，如果不使用步骤2和步骤5，则需要再步骤6中预测出的PM2.5值无需进行反归一化操作。
本专利中的预测目标PM2.5的值也可以为能见度等可以使用实数进行度量的其他实际预测目标，不论使用何种预测目标，只要使用了图像水平方向和垂直方向二阶差分加权平均值、二阶梯度和的加权平均值以及Qx、Qy和Qg折线差分的加权平均值这几个特征就都属于本专利的权利要求范围。
以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的权利要求范围之内。

Claims (6)

1.一种基于图像特征的PM2.5污染指数的识别方法，其特征在于：包括以下步骤：
步骤1.搜集标注了PM2.5实测值的历史图像，将搜集的历史图像全部统一缩放为XS×YS像素大小；如果图像是彩色的，那么就将其转换为灰度图；
步骤2.计算每个缩放并转换为灰度图后图像的图像特征；
步骤3.将每个图像的特征作为输入，使用回归方法得到PM2.5的预测模型；
步骤4.采集目标图像，按照步骤1中同样的方法将目标图像缩放为XS×YS像素大小，如果图像是彩色的，将其转换为灰度图；
步骤5.针对缩放后的目标图像，按照步骤2中的同样的方法，计算缩放并转换为灰度图后图像的图像特征；
步骤6.将目标图像的图像特征作为步骤3中已经训练好的预测模型的输入，计算得到目标图像的PM2.5值；
所述步骤2中在计算每个缩放并转换为灰度图后图像的图像特征之后,还需要对每个图像的图像特征以及PM2.5的值进行归一化；所述步骤5中针对缩放后的目标图像，在计算缩放并转换为灰度图后图像的图像特征之后，还需要对目标图像的图像特征按照步骤2中同样的归一化方法，进行归一化；
所述步骤1或步骤4中采用公式GRAY＝0.299×R+0.587×G+0.114×B将彩色图像转换为灰度图像，其中R、G、B分别表示彩色图像的红色、绿色和蓝色波段；
所述步骤2或步骤5中的图像特征为图像的平均像素值、局部极值点的个数NMax、局部方差的均值μLVAR、图像水平方向和垂直方向二阶差分加权平均值、图像水平方向和垂直方向的二阶梯度和的加权平均值、Qx、Qy和Qg折线差分的加权平均值以及图像中的斑点数,所述Qx、Qy、Qg的定义如下述步骤a-b中所述；
所述图像的平均像素值为所有像元的像素值加起来之后除以XS×YS；
所述图像的局部极值点的个数NMax和局部方差的均值μLVAR的计算方法为：
(1)设定NMax和μLVAR的初值为0；
(2)对于图像中的每一个非边界像元找出邻近的8个像元，记8个像元的构成集合为PIXA；如果当前像元坐标为(i,j),i和j分别表示该像元和图像最左上角像元分别在垂直和水平方向同当前像元间隔的像元个数，那么其8个近邻的坐标分别为(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)；这里的非边界像元是指图像中所有能够找到上述8个近邻，并且任一近邻像元像素值不为空的像元；
(3)如果当前像元值大于周边的8个像元，那么NMax值加1；并且
这里，p表示当前像元，I(p)表示像元p的灰度值，表示PIXA中的8个像元和当前像元像素值的均值；
(4)当所有像元都处理完之后，μLVAR＝μLVAR/((XS-2)×(YS-2))；
所述图像水平方向和垂直方向二阶差分加权平均值、图像水平方向和垂直方向二阶梯度和的加权平均值的计算方法为：
a.分别计算图像水平方向x和垂直方向y每个像元p的二阶差分,分别记作和并分别将两个方向的二阶差分按照从小到大的序列排序，排序过程中去掉重复值，并过滤掉小于水平方向给定阈值αx或垂直方向给定阈值αy的所有差分值，由此得到两个大小分别为nx和ny的水平方向队列Qx和垂直方向队列Qy，两个队列中第i的元素分别为和即为队列中排序为i的像元所对应的二阶差分值；
b.计算图像每个像元p的二阶梯度和并将二阶梯度和值按照从小到大的序列排序，排序过程中去掉重复值，并过滤掉小于给定阈值αg的所有二阶梯度和值，由此得到大小为ng的二阶梯度和队列Qg，队列中第i的元素分别为Qg(i)＝grad(pi)，即为队列中排序为i的像元所对应的二阶梯度和值；
c.分别对三个队列Qx、Qy和Qg计算对应折线总长度和以及总体加权变化度 和
d.分别计算三个队列Qx、Qy和Qg的加权平均值WVLx＝Vx/Lx、WVLy＝Vy/Ly和WVLg＝Vg/Lg，这三个值就分别为水平方向二阶差分加权平均值、垂直方向二阶差分加权平均值以及二阶梯度和的加权平均值；
所述图像中的斑点数的计算方法为：
(1)给定的一系列有序尺度值σ1,σ2,…,σs,s≥3，并选定某个半径r，根据规范化拉普拉斯高斯算子生成对应各个尺度的卷积模板，其中每个卷积模板的中心坐标为(0，0)，所有和卷积模板中心水平或者垂直方向相差r个像元的所有网格都属于模板；模板中所有网格的坐标集合为{(tx,ty)|tx＝-r,-r+1,…,r-1,r∧ty＝-r,-r+1,…,r-1,r}，对应尺度σ∈{σ1,σ2,…,σs}的卷积模板中每个网格(tx,ty)上的值使用下面公式计算：
对每个尺度σ∈{σ1,σ2,…,σs}，使用对应的卷积模板对图像做卷积运算，得到一系列卷积图像I1,I2,…,Is；
(2)对于每个卷积图像Ii,i＝2,…,s-1中的每个像元p，假设其在图像中的坐标为(i,j)，其8个最近邻的像元坐标分别为(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)；如果Ii(p)是所有近邻中最小的，并且Ii(p)也比Ii-1和Ii+1中所有坐标为(i,j)、(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)的像元的值都小，那么该像元对应的就是一个斑点；如果Ii(p)是所有近邻中最大的，并且Ii(p)也比Ii-1和Ii+1中所有坐标为(i,j)、(i+1,j)、(i,j+1)、(i-1,j)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)和(i+1,j+1)的像元的值都大，那么该像元对应的也是一个斑点；在所有卷积图像Ii,i＝2,…,s-1上找到的斑点总数即为图像中的斑点数。
2.根据权利要求1所述的一种基于图像特征的PM2.5污染指数的识别方法，其特征在于：所述Qx、Qy和Qg折线差分的加权平均值的计算方法为：
(1)首先针对Qx、Qy和Qg的折线分别计算相邻两个元素之间的差值，并且将这一差值按照从小到大排序，分别形成一个队列，分别记为QQx、QQy和QQg，其第i个元素分别为QQx(i)、QQy(i)和QQg(i)；
(2)分别对三个队列QQx、QQy和QQg计算对应折线总长度和以及总体加权变化度 和
(3)分别计算三个队列QQx、QQy和QQg的加权平均值WVLQx＝VQx/LQx、WVLQy＝VQy/LQy和WVLQg＝VQg/LQg，这三个值就分别为Qx、Qy和Qg折线差分的加权平均值。
3.根据权利要求2所述的一种基于图像特征的PM2.5污染指数的识别方法，其特征在于：
所述步骤2和步骤5中按照如下公式对图像的特征和PM2.5的值进行归一化：
其中F表示某个特征，F(I)表示缩放并转换为灰度图后的历史图像或者目标图像I在特征F上的取值，minF表示特征F在所有缩放并转换为灰度图后的历史图像上的最小值，maxF表示特征F在所有缩放并转换为灰度图后的历史图像上的最大值,Fnorm(I)表示I在特征F上归一化后的值，并且在对转换为灰度图后的目标图像进行归一化时若Fnorm(I)＞1时将Fnorm(I)置为1，若Fnorm(I)＜0时将Fnorm(I)置为0,Inew为待归一化的新值。
4.根据权利要求3所述的一种基于图像特征的PM2.5污染指数的识别方法，其特征在于：所述步骤3中使用回归方法得到PM2.5的预测模型的方法为前馈神经网络、高斯过程或常用的最小二乘方法。
5.根据权利要求4所述的一种基于图像特征的PM2.5污染指数的识别方法，其特征在于：所述使用回归方法得到PM2.5的预测模型的方法为前馈神经网络，具体操作步骤为：
(1)在前馈神经网络中设定好激活函数，这一步可以使用逻辑斯特函数或者双曲正切函数；
(2)设定神经网络的层次和每层的神经元数量；
(3)使用反向传播算法学习神经网络；
(4)将学习到的网络保存下来，需要保存的内容有使用了何种激活函数，神经网络的层数、每层的神经元数量以及每个神经元的权重。
6.根据权利要求5所述的一种基于图像特征的PM2.5污染指数的识别方法，其特征在于：采用前馈神经网络方法得到PM2.5的预测模型时，需要读入保存的基于前馈神经网络的预测模型参数，包括激活函数，神经网络的层数、每层的神经元数量以及每个神经元的权重，然后按照前馈神经网络的前馈传播方法将归一化后的图像特征作为输入计算出Inew所对应的归一化后的PM2.5值PM2.5norm(Inew)，最后将其按照下面公式反归一化后作为新拍摄图像的PM2.5值PM2.5(Inew)：
PM2.5(Inew)＝PM2.5norm(Inew)×(maxF-minF)+minF。
CN201710301867.4A 2017-05-02 2017-05-02 一种基于图像特征的pm2.5污染指数的识别方法 Active CN107044947B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710301867.4A CN107044947B (zh) 2017-05-02 2017-05-02 一种基于图像特征的pm2.5污染指数的识别方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710301867.4A CN107044947B (zh) 2017-05-02 2017-05-02 一种基于图像特征的pm2.5污染指数的识别方法

Publications (2)

Publication Number Publication Date
CN107044947A CN107044947A (zh) 2017-08-15
CN107044947B true CN107044947B (zh) 2019-11-19

Family

ID=59546232

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710301867.4A Active CN107044947B (zh) 2017-05-02 2017-05-02 一种基于图像特征的pm2.5污染指数的识别方法

Country Status (1)

Country Link
CN (1) CN107044947B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109087277B (zh) * 2018-06-11 2021-02-26 北京工业大学 一种空气细颗粒物pm2.5测量方法
CN111488804B (zh) * 2020-03-19 2022-11-11 山西大学 基于深度学习的劳保用品佩戴情况检测和身份识别的方法
CN111912755B (zh) * 2020-08-07 2021-08-10 山东中煤工矿物资集团有限公司 一种矿用粉尘浓度传感器、传感器系统及方法

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6310346B1 (en) * 1997-05-30 2001-10-30 University Of Central Florida Wavelength-tunable coupled antenna uncooled infrared (IR) sensor
CN103903273A (zh) * 2014-04-17 2014-07-02 北京邮电大学 一种基于手机终端的pm2.5等级快速评价系统
CN103954542A (zh) * 2014-05-12 2014-07-30 中国计量学院 基于无参考图像清晰度评价的pm2.5浓度检测仪
CN104462778A (zh) * 2014-11-06 2015-03-25 华北电力大学 一种基于深度学习的pm2.5污染物测量方法
CN103593660B (zh) * 2013-11-27 2016-08-17 青岛大学 一种稳定特征图像下交叉梯度编码的掌纹识别方法
CN106295516A (zh) * 2016-07-25 2017-01-04 天津大学 基于图像的雾霾pm2.5值估计方法
CN106401359A (zh) * 2016-08-31 2017-02-15 余姚市泗门印刷厂 基于红外摄像的窗体操控平台

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6310346B1 (en) * 1997-05-30 2001-10-30 University Of Central Florida Wavelength-tunable coupled antenna uncooled infrared (IR) sensor
CN103593660B (zh) * 2013-11-27 2016-08-17 青岛大学 一种稳定特征图像下交叉梯度编码的掌纹识别方法
CN103903273A (zh) * 2014-04-17 2014-07-02 北京邮电大学 一种基于手机终端的pm2.5等级快速评价系统
CN103954542A (zh) * 2014-05-12 2014-07-30 中国计量学院 基于无参考图像清晰度评价的pm2.5浓度检测仪
CN104462778A (zh) * 2014-11-06 2015-03-25 华北电力大学 一种基于深度学习的pm2.5污染物测量方法
CN106295516A (zh) * 2016-07-25 2017-01-04 天津大学 基于图像的雾霾pm2.5值估计方法
CN106401359A (zh) * 2016-08-31 2017-02-15 余姚市泗门印刷厂 基于红外摄像的窗体操控平台

Also Published As

Publication number Publication date
CN107044947A (zh) 2017-08-15

Similar Documents

Publication Publication Date Title
EP3614308B1 (en) Joint deep learning for land cover and land use classification
CN107909109B (zh) 基于显著性和多尺度深度网络模型的sar图像分类方法
CN108921203B (zh) 一种指针式水表的检测与识别方法
CN108334847A (zh) 一种真实场景下的基于深度学习的人脸识别方法
CN108615046A (zh) 一种储粮害虫检测识别方法及装置
CN107044947B (zh) 一种基于图像特征的pm2.5污染指数的识别方法
CN107665492A (zh) 基于深度网络的结直肠全景数字病理图像组织分割方法
CN109064462A (zh) 一种基于深度学习的钢轨表面缺陷检测方法
CN108629369A (zh) 一种基于Trimmed SSD的尿沉渣有形成分自动识别方法
CN109509187A (zh) 一种针对大分辨率布匹图像中的小瑕疵的高效检验算法
CN104239902B (zh) 基于非局部相似性和稀疏编码的高光谱图像分类方法
CN111860459B (zh) 一种基于显微图像的禾本科植物叶片气孔指数测量方法
CN110517311A (zh) 基于叶片病斑面积的病虫害监测方法
CN113435282B (zh) 基于深度学习的无人机影像麦穗识别方法
CN111724355B (zh) 一种鲍鱼体型参数的图像测量方法
CN110569747A (zh) 一种利用图像金字塔与Faster-RCNN快速统计大田水稻稻穗数的方法
CN107978110A (zh) 基于图像匹配的电子围栏智能在位识别系统及识别方法
CN109815953A (zh) 一种基于车辆年检标的车辆识别匹配系统
CN109635634A (zh) 一种基于随机线性插值的行人再识别数据增强方法
CN111426637A (zh) 一种利用无人机获取城市黑臭河流水质参数的方法
KR101933856B1 (ko) 콘벌루션 신경망을 이용한 영상 처리 시스템 및 이를 이용한 영상 처리 방법
CN106650823A (zh) 一种基于概率极限学习机集成的泡沫镍表面缺陷分类方法
CN110298410A (zh) 基于深度学习的低对比度图像中弱目标检测方法及装置
CN113469097B (zh) 一种基于ssd网络的水面漂浮物多相机实时检测方法
CN112085725B (zh) 基于启发式迭代算法的残膜残留量检测方法及预警系统

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
GR01 Patent grant
GR01 Patent grant
TR01 Transfer of patent right

Effective date of registration: 20211222

Address after: 101100 4th floor, building 6, No.10 yard, Jiachuang Road, Tongzhou District, Beijing

Patentee after: BEIJING JIAHUA ZHILIAN TECHNOLOGY CO.,LTD.

Address before: 030006 No. 92, Hollywood Road, Taiyuan, Shanxi

Patentee before: SHANXI University

TR01 Transfer of patent right