CN111192239B - 遥感影像变化区域检测方法、装置、存储介质及电子设备 - Google Patents

遥感影像变化区域检测方法、装置、存储介质及电子设备 Download PDF

Info

Publication number
CN111192239B
CN111192239B CN201911306081.7A CN201911306081A CN111192239B CN 111192239 B CN111192239 B CN 111192239B CN 201911306081 A CN201911306081 A CN 201911306081A CN 111192239 B CN111192239 B CN 111192239B
Authority
CN
China
Prior art keywords
pixel
vector
points
remote sensing
change area
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
CN201911306081.7A
Other languages
English (en)
Other versions
CN111192239A (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.)
Stargis Tianjin Technology Development Co ltd
University of Science and Technology of China USTC
Original Assignee
Stargis Tianjin Technology Development Co ltd
University of Science and Technology of China USTC
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 Stargis Tianjin Technology Development Co ltd, University of Science and Technology of China USTC filed Critical Stargis Tianjin Technology Development Co ltd
Priority to CN201911306081.7A priority Critical patent/CN111192239B/zh
Publication of CN111192239A publication Critical patent/CN111192239A/zh
Application granted granted Critical
Publication of CN111192239B publication Critical patent/CN111192239B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/0002Inspection of images, e.g. flaw detection
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00Pattern recognition
    • G06F18/20Analysing
    • G06F18/24Classification techniques
    • G06F18/241Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06VIMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V20/00Scenes; Scene-specific elements
    • G06V20/10Terrestrial scenes
    • G06V20/13Satellite images
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10032Satellite or aerial image; Remote sensing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30181Earth observation

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Data Mining & Analysis (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • General Engineering & Computer Science (AREA)
  • Artificial Intelligence (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Astronomy & Astrophysics (AREA)
  • Remote Sensing (AREA)
  • Multimedia (AREA)
  • Quality & Reliability (AREA)
  • Image Processing (AREA)

Abstract

本发明提供了一种遥感影像变化区域检测方法、装置、存储介质及电子设备,所述方法包括:对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,所述像素向量用于表示当前像素点与其指定邻域范围内其他像素点的相关性;查找各个分类类别中包含像素向量最少的目标分类类别;将所述目标分类类别中各个像素向量对应的像素点作为目标变化区域,对所述目标变化区域进行图像优化处理。本发明能够高效、准确地从两张遥感影像的差异图中检测用户实际感兴趣的变化区域。

Description

遥感影像变化区域检测方法、装置、存储介质及电子设备
技术领域
本发明涉及遥感影像处理技术领域,尤其涉及一种遥感影像变化区域检测方法、装置、存储介质及电子设备。
背景技术
遥感影像的对比变化是指在不同时间点,相同地区的遥感图片中显示出来的改变,这些改变包含建筑物的增加、植被的覆盖、土地的开发利用等。变化检测就是通过处理比较图片,找出随着时间发生变化的位置,标识出发生变化的区域所对应的像素点。
变化检测(Change Detection,CD)在利用现有遥感影像数据方面,在准确理解地表变化方面有着非常重要的作用。由于遥感影像具有容易获得且更新周期短的特点,从遥感影像中准确、及时地检测出变化区域,对于城市管理和规划、城市化评估以及灾后重建等多项工作都具有重要意义。将检测到的各个区域中新建的建筑物,结合其所在的区域,可以有效地检测出违规建筑物,便于城市规范化管理。分析遥感影像中的变化密集程度便于了解各个区域发展速度以及建筑物改变的密集区域,这有助于估计地区的发展进程并能评估地区的发展程度。通过定位发生剧烈变化的区域,能够判断受到地震等自然灾害的地区是否按照规划方式重新建设,监督灾后重建的进程。
为了得到图像的变化情况及位置,测绘等相关行业以人工的方式用ARCGIS等工具通过对比两张遥感影像,并标注出发生变化的区域。由于建筑物分布密集,需要细致地对比每个区域中存在的变化,同时由于遥感影像本身具有对比度低、城区建筑物分布密集、需要检测的区域范围大等特点,导致现有遥感影像的变化区域检测虽然消耗了大量的时间与精力,但是由于差异图片中许多像素点是并不关心的些许变化,因此,往往消耗大量的时间与精力标记出的变化区域并不是实际感兴趣的变化,造成时间和精力的不必要浪费,且影响检测效果。
发明内容
本发明提出了一种遥感影像变化区域检测方法、装置、存储介质及电子设备,解决了现有技术中难以高效、准确地从两张遥感影像的差异图中检测用户实际感兴趣的变化区域的问题。
本发明的一个方面,提供了一种遥感影像变化区域检测方法,所述方法包括:
对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,所述像素向量用于表示当前像素点与其指定邻域范围内其他像素点的相关性;
查找各个分类类别中包含像素向量最少的目标分类类别;
将所述目标分类类别中各个像素向量对应的像素点作为目标变化区域,对所述目标变化区域进行图像优化处理。
可选地,所述对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,包括:
配置分类类别的类别数;
根据所述分类类别的类别数将所述差异图中各个像素点的像素向量按照预设划分标准划分为不同分类类别。
可选地,所述根据分类类别的类别数将所述差异图中各个像素点的像素向量按照预设划分标准划分为不同分类类别,包括:
计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的距离;
若计算出的距离最小值大于任意两个类别点的最大距离,则当前像素点的像素向量被分为一个新的类,并将对应的两个类别合并;
否则,当前像素点的像素向量被分到与其距离最小的像素向量所属的分类类别。
可选地,所述计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的距离,包括:
计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的L1距离、L2距离或者直方图距离。
可选地,所述对所述目标变化区域进行图像优化处理,包括:
采用预设大小的扫描单元依次对所述目标变化区域内的像素点进行扫描;
利用第一预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于或等于预设的第一阈值,则去除当前处于所述扫描单元中心位置的像素点,所述第一预设矩阵具体如下:
Figure BDA0002323110060000031
利用第二预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于预设的第二阈值,则根据当前处于所述扫描单元中心位置的像素点对所述扫描单元中的所有像素点进行像素信息重置,所述第二预设矩阵具体如下:
Figure BDA0002323110060000032
可选地,所述对所述目标变化区域进行图像优化处理,还包括:
对像素信息重置后的目标变化区域进行矢量化处理。
本发明的另一个方面,提供了一种遥感影像变化区域检测装置,所述装置包括:
配置单元,用于对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,所述像素向量用于表示当前像素点与其指定邻域范围内其他像素点的相关性;
检测单元,用于查找各个分类类别中包含像素向量最少的目标分类类别;
处理单元,用于将所述目标分类类别中各个像素向量对应的像素点作为目标变化区域,对所述目标变化区域进行图像优化处理。
可选地,所述配置单元,包括:
配置子单元,用于配置分类类别的类别数;
分类子单元,用于根据所述分类类别的类别数将所述差异图中各个像素点的像素向量按照预设划分标准划分为不同分类类别。
可选地,所述分类子单元,具体用于计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的距离;若计算出的距离最小值大于任意两个类别点的最大距离,则当前像素点的像素向量被分为一个新的类,并将对应的两个类别合并;否则,当前像素点的像素向量被分到与其距离最小的像素向量所属的分类类别。
可选地,所述处理单元,包括:
扫描子单元,用于采用预设大小的扫描单元依次对所述目标变化区域内的像素点进行扫描;
第一计算子单元,用于利用第一预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于或等于预设的第一阈值,则去除当前处于所述扫描单元中心位置的像素点,所述第一预设矩阵具体如下:
Figure BDA0002323110060000041
第二计算子单元,用于利用第二预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于预设的第二阈值,则根据当前处于所述扫描单元中心位置的像素点对所述扫描单元中的所有像素点进行像素信息重置,所述第二预设矩阵具体如下:
Figure BDA0002323110060000051
可选地,所述处理单元,还包括:
矢量化子单元,用于对像素信息重置后的目标变化区域进行矢量化处理。
此外,本发明还提供了一种计算机可读存储介质,其上存储有计算机程序,该程序被处理器执行时实现如上所述方法的步骤。
此外,本发明还提供了一种电子设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,所述处理器执行所述程序时实现如上所述方法的步骤。
本发明实施例提供的遥感影像变化区域检测方法、装置、存储介质及电子设备,通过对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,进而确定出用户实际感兴趣目标变化区域,实现高效、准确地从两张遥感影像的差异图中检测用户实际感兴趣的变化区域。
上述说明仅是本发明技术方案的概述,为了能够更清楚了解本发明的技术手段,而可依照说明书的内容予以实施,并且为了让本发明的上述和其它目的、特征和优点能够更明显易懂,以下特举本发明的具体实施方式。
附图说明
通过阅读下文优选实施方式的详细描述,各种其他的优点和益处对于本领域普通技术人员将变得清楚明了。附图仅用于示出优选实施方式的目的,而并不认为是对本发明的限制。而且在整个附图中,用相同的参考符号表示相同的部件。在附图中:
图1为本发明实施例的一种遥感影像变化区域检测方法的流程示意图;
图2为本发明实施例的一种遥感影像变化区域检测方法中步骤S11的实现流程示意图;
图3为本发明实施例的一种遥感影像变化区域检测装置的结构示意图;
图4为本发明实施例的一种遥感影像变化区域检测装置中配置单元的内部结构示意图。
具体实施方式
下面将参照附图更详细地描述本公开的示例性实施例。虽然附图中显示了本公开的示例性实施例,然而应当理解,可以以各种形式实现本公开而不应被这里阐述的实施例所限制。相反,提供这些实施例是为了能够更透彻地理解本公开,并且能够将本公开的范围完整的传达给本领域的技术人员。
本技术领域技术人员可以理解,除非另外定义,这里使用的所有术语(包括技术术语和科学术语),具有与本发明所属领域中的普通技术人员的一般理解相同的意义。还应该理解的是,诸如通用字典中定义的那些术语,应该被理解为具有与现有技术的上下文中的意义一致的意义,并且除非被特定定义,否则不会用理想化或过于正式的含义来解释。
图1示意性示出了本发明一个实施例的遥感影像变化区域检测方法的流程图。参照图1所示,本发明实施例提出的遥感影像变化区域检测方法具体包括步骤S11~S13,如下所示:
S11、对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,所述像素向量用于表示当前像素点与其指定邻域范围内其他像素点的相关性。
遥感影像的差异图指的是由两张三通道影像图片得到的一张单通道图片。得到差异图的方法有多种,包括直接进行差值计算、调整图片通道权重计算以及结合像素点通道间差值计算。
本发明实施例中,初步得到的差异图包含了几乎所有微小的变化区域,无法用来直接生成最终的变化图像,因此在实现了遥感影像的差异图初步获取之后,需要计算所述差异图中各个像素点与其指定邻域范围内其他像素点的相关信息,以便于后续变化检测中能够结合各个像素点的邻域信息来对像素点进行变化区域分离检测,进而提高基于遥感影像差异图进行变化区域检测的效率和泛化能力。
S12、查找各个分类类别中包含像素向量最少的目标分类类别。
S13、将所述目标分类类别中各个像素向量对应的像素点作为目标变化区域,对所述目标变化区域进行图像优化处理。
本发明实施例提供的遥感影像变化区域检测方法,通过对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,进而确定出用户实际感兴趣目标变化区域,实现高效、准确地从两张遥感影像的差异图中检测用户实际感兴趣的变化区域。
下面通过一个具体实施例对影像的差异图获取的具体实现方式进行详细说明。本发明实施例中,首先需要得到两张在不同时间拍摄地点相同的遥感影像,本实施例中,以图片的大小均为1200*1200*3为例进行说明,允许两张图片拍摄条件不完全相同,即可能存在色差、光照角度及建筑物倾斜角度的差异。影像的差异图指的是由两张三通道图片得到的一张单通道图片。得到差异图的方法有多种,包括直接进行差值计算、调整图片通道权重计算以及结合像素点通道间差值计算,可根据所需检测差异种类从以下方法中选择:
方法一:直接进行差值计算。计算量最小的方案是将每张图片的三个通道值取平均值。采用这种方法可以最小化计算量,并且无差别地检测到图像中所有的变化位置,但是考虑到实际图像中由于光照、云层等各种原因,导致两张图片在某些没有变化的位置出现数值上的差异,会给差异图带来较大的噪声干扰。
具体实现方法包括:分别提取两幅遥感影像的R、G、B三个通道的特征值;根据两幅遥感影像中各个像素点对应的R、G、B三个通道的特征值的均值,分别生成每一遥感影像对应的灰度图片;分别计算两幅遥感影像对应的灰度图片中位于相同位置的像素点的差值,根据两幅遥感影像中各个像素点的差值生成两幅遥感影像的差异图。
本实施例中,分别提取两幅遥感影像的R、G、B三个通道的特征值,得到彩色图片不同通道的数值。对于每一张彩色图片,分别分离图片的三个通道(R、G、B)。
将每张图的三个通道的值,取对应每个像素点的均值,并注意取得到均值的整数值。采用如下公式:
Sk=αR+βG+γB
α=0.33,β=0.33,γ=0.33
计算两张灰度图片的差异,将对应位置的两个像素点直接相减,并且取结果的绝对值,得到需要的差异图片。具体的,两张图片的S1与S2,计算差异图为S1与S2的差的绝对值:
S=|S1-S2|
方法二:调整图片通道权重计算。图像的通道权重指的是在形成灰度图时,三个通道的数值在最终结果中占的比例。两张不同图片通常具有不同的光照强度和色差,这种影响范围是覆盖整个图片的,为了解决这个问题,所采用的方法是计算颜色通道的对应直方图之间的距离,确定新的通道计算时的权重,这样可以在很大程度上减小色差、光照带来的干扰。
具体实现方法包括:分别提取两幅遥感影像的R、G、B三个通道的特征值;对每一遥感影像的R、G、B三个通道的特征值分别进行直方图统计;计算每一遥感影像对应的三个直方图的直方图距离;根据每一遥感影像对应的三个直方图的直方图距离的大小关系,确定R、G、B三个通道的特征值在生成对应的灰度图片时所占的权重值,并根据确定的权重值生成当前遥感影像对应的灰度图片;分别计算两幅遥感影像对应的灰度图片中位于相同位置的像素点的差值,根据两幅遥感影像中各个像素点的差值生成两幅遥感影像的差异图。
本实施例中,对每张图片的各个通道的像素值进行统计,像素点的取值范围为0~255的整数值,则统计不同取值的像素点个数,每个梯度直方图一共有256个柱,计算三个通道的共3个梯度直方图。
具体可根据以下公式计算三个对应直方图之间的距离,直方图距离的计算公式如下:
Figure BDA0002323110060000091
Figure BDA0002323110060000092
其中,f[j]为一个直方图中第j个柱状图的值。
距离最小的通道为G,将其权重设为0.2,距离最大的通道为R,将它的权重设为0.5,减小色差、光照带来的干扰。
根据以上权值计算,得到两张图片的S1与S2,计算差异图为S1与S2的差的绝对值。
方法三:结合像素点通道间差值计算。考虑到需要检测变化的目标的不同,可以调整对应通道差值所占的比例。在需要检测的物体具有一定的颜色特征的时候,根据这个特征调整权值。
具体实现方法包括:获取待检测的物体的目标颜色特征;分别提取两幅遥感影像的R、G、B三个通道的特征值;计算每一遥感影像中各个像素点的与所述目标颜色特征对应的颜色通道与其他通道间特征值的差值;对每一遥感影像中各个像素点的与所述目标颜色特征对应的颜色通道的特征值以及该颜色通道与其他通道间特征值的差值在生成灰度图片时所占的权重值进行调整,并根据调整后的权重值生成当前遥感影像对应的灰度图片;分别计算两幅遥感影像对应的灰度图片中位于相同位置的像素点的差值,根据两幅遥感影像中各个像素点的差值生成两幅遥感影像的差异图。
本实施例中,需要检测的是房屋变化,在图像中,通常在旷地上变化的房屋具有红色或者蓝色的屋顶,而绿色植被、地表的变化并不被需要。
计算每张图片的三个通道将相互的差值,差值小于0的值设为0,分别计算R、G、B三个通道间共6个差值。
确定通道差值的权重。现在对建筑的变化最感兴趣,那么对于水域与植被的变化就不需要关心。植被与水域都呈现较为明显的绿色,那么像素点中绿色通道(G)与其他两个蓝色通道(B)与红色通道(R)之间的差值大的点就不需要过于关注,减小其权重,并增加R通道与其他通道差值的权重,可选地,设置α为0.4,β为0.4,γ为0.2。
生成差异图。根据以上确定的权重,根据以下公式计算每张图的新的灰度图对于红色分量最感兴趣,结果可以写成:
S=α(R-B)+β(R-G)+γR
α+β+γ=1
本发明实施例中,在实现了遥感影像的差异图初步获取,初步得到的差异图包含了几乎所有微小的变化区域,无法用来直接生成最终的变化图像,因此需要基于邻域信息对像素点进行调整处理,具体包括如下附图中未示出的步骤:
S101、为待处理的差异图中的每一像素点划分对应的邻域范围。
S102、根据每一邻域范围中各个像素点的像素信息生成当前邻域范围对应像素点的像素向量。
S103、根据所述差异图中每一像素点的像素向量生成目标矩阵。
S104、对所述目标矩阵进行稀疏化处理,采用稀疏化处理后的矩阵体现所述差异图中各个像素点与其邻域范围内其他像素点的相关信息。
其中,对所述目标矩阵进行稀疏化处理,具体包括以下附图中未示出的步骤:
S1041、按照预设规则从所述差异图中选取用于变化检测的目标像素点。具体的,为了避免造成专图专用结果,仅从一张图中选取一部分点进行计算,选取的方法有许多,主要保证的原则是选取点均匀分布、不重复、数量适中,本发明采用均匀随机分布的取样方式实现目标像素点选取,具体选取方式包括但不限于以下两种:在所述差异图的同一像素行中,每次相隔相同数量像素列选取一个目标像素点,且在同一像素列中,每次相隔相同数量的像素行选取一个目标像素点;以及根据预设划分标准将所述差异图中划分为多个像素单元,选取各个像素单元中位于同一分布位置的像素点作为目标像素点。
S1042、根据各个目标像素点对应的像素向量生成源矩阵。
S1043、计算所述源矩阵的协方差矩阵的变换矩阵。本实施例中,为了实现目标矩阵的线性变换,本发明通过计算源矩阵对应的变换矩阵实现对所述目标矩阵进行变换处理,其中,变换矩阵的计算方法包括:对所述源矩阵中的每一行进行零均值化处理;计算零均值化处理后的矩阵的协方差矩阵,并计算所述协方差矩阵的特征值和特征向量;根据特征值的大小顺序将所述特征向量从上到下依次排列,得到所述变换矩阵。
S1044、根据所述变换矩阵对所述目标矩阵进行变换处理。
进一步地,在所述根据所述变换矩阵对所述目标矩阵进行变换处理之后,所述方法还包括:对变换处理后的矩阵进行稀疏化处理。具体实现方式为:查找变换处理后的矩阵中存在的小于预设阈值的像素信息;将查找到的像素信息设置为0。
下面通过一个具体实施例对计算差异图中各个像素点与其邻域范围内其他像素点的相关信息的具体实现方法进行解释说明。
本实施例中,单个像素点往往只能反映其所在位置的变化程度,但是对于图像来说,感兴趣的变化必然是具有连续性的,所以像素点周围点的信息需要纳入考虑范围。本实施例将一个像素点相邻范围内的像素点加入,选取的邻域范围为5*5大小的区域,并将其展平,以一个像素点的信息从一个值变成一个向量1*25为大小的向量,这个向量代表了这个像素点及其周围点对这个点的贡献值。
为了避免造成专图专用结果,仅从一张图中选取一部分点进行计算,选取的方法有许多,主要保证的原则是选取点均匀分布、不重复、数量适中。
确定选取方案。本实施例中采用均匀随机分布的取样方式,选取的邻域范围为5*5的区域,图片大小为1200*1200,在没有重复的情况下,一张图片具有240*240个5*5的区域,得到了240*240个像素点,每个像素点和其上、下、左、右相邻的像素点的距离均为5,像素点的分布覆盖了整个图片,且分布均匀,满足选取点均匀分布、不重复、数量适中的原则。
将选取的像素点及其邻域所组成的向量,形成一个向量集,可以看作是一个矩阵,矩阵的每一列表示一个像素点的向量。操作步骤是将之前得到的240*240个像素点的25维矩阵组成一个25行240*240列的矩阵X。
减少像素点向量的维数。向量之间通常都具有线性关系,为了减少向量的维数,找到这个矩阵的特征向量。具体的操作步骤如下:
均值化。将240*240个像素点的25维向量组成一个25行240*240列的矩阵X,对X的每一行进行零均值化,即减去这一行的均值。
计算协方差矩阵。对之前进行过均值化的矩阵X,按如下公式求出协方差矩阵C:
Figure BDA0002323110060000121
求出变换矩阵。求出C的特征值与特征向量,将特征向量按特征值大小从上到下排成新的矩阵,得到的就是需要找的变换矩阵K,矩阵的维度为25*25。
将变换应用到所有像素点。前述步骤所得到的K是从一张图片的部分点得到的,认为这个矩阵K包括了适用于所有点的为了减少数据维度而作的线性变换。
包含所有像素点。对差异图的每一个像素点,取其5*5的邻域范围得到一个向量,即一共有1200*1200个向量,对于边缘点,将其不存在的邻域范围补充为0。
组成新的矩阵。25为一个向量的维数,1200*1200为整个图像中像素点的总数,将所有的向量组成一个行数为25的1200*1200列新矩阵X_n。
应用变换。将变换K应用于X_n,根据如下公式:
X_n=K·X_n
得到的结果是一个25行1200*1200列的矩阵,但是这个矩阵应该比原来的X_n更加稀疏,或者说有更多的值接近于0。
稀疏化。为了达到稀疏的目的,设定一个阈值5,让结果中小于5的值都设为0,在线性变化中引入非线性变换。
通过以上步骤,将一张由1200*1200个像素点组成的差异图变成了一个25行1200*1200列的稀疏矩阵。
本发明实施例中,所述对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,如图2所示,具体实现流程如下:
S111、配置分类类别的类别数;
S112、根据分类类别的类别数将所述差异图中各个像素点的像素向量按照预设划分标准划分为不同分类类别。其中,根据所述分类类别的类别数将所述差异图中各个像素点的像素向量按照预设划分标准划分为不同分类类别,具体包括:计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的距离,所述距离包括但不限于L1距离、L2距离或者直方图距离;若计算出的距离最小值大于任意两个类别点的最大距离,则当前像素点的像素向量被分为一个新的类,并将对应的两个类别合并;否则,当前像素点的像素向量被分到与其距离最小的像素向量所属的分类类别。
本实施例中,由于差异图片中许多像素点是并不关心的些许变化,往往感兴趣的变化是大范围且稳定的,因此检测出来的像素点可以看作是两大类组成:具有聚合性与连续性的真正的变化点,以及离散的非变化点。
对向量进行分类。通过之前的步骤,已经得到每个像素点及其周围领域的相关信息,需要根据这些信息对所有的像素点进行分类。理论上,所有成熟可靠的分类方法都可以应用在这个过程。由于这个分类问题没有具体的标签,即无法知道哪一个像素点实际上应该属于哪一个分类,因此只能将无监督的分类方法应用到过程之中,具体采用的步骤是:
确定类别数目。事先确定需要分成的类的数量,因为不具有标签,分类数决定了最后每个类别中元素的数量,选择的类别的数量不会超过5个,最后确定选择类别数量目为3。
进行聚类分类。从第一个点的向量开始,每一次加入一个新的向量进行计算,计算新加入的向量与其他向量之间的距离,距离可以选择采用L1距离、L2距离或者直方图距离等。本实施例中距离选择卡方距离,对于向量值,卡方距离可以写成:
Figure BDA0002323110060000141
其中a,b为两个向量,N为向量维数,值为25,如果距离最小的值大于某两类点的最大距离,那么新的向量被分为一个新的类,并将对应的两类合并,否则,新加入的向量被分到与它最近的向量所属的那一类。
选择结果分类的类别。选择的正确的类别作为最终结果往往决定了最后结果的好坏程度,经验结果表明,真正需要关注的变化的点在分类之后数量最少的类别中。因此,本实施例中利用预设程序实现选择结果分类,即选择结果中3个大类中向量数量最少的分类集合,并认为这个集合中所有向量对应的像素点就是检测到的变化点。
本发明实施例中,得到的目标变化区域依然可能包括一些噪声的干扰,产生的原因是由于大面积的土地结构改变,表现在图像上是小的、分散的、不规则的点状。这些确实是变化的一部分,但是在应用中不需要关注,所以需要对所述目标变化区域进行图像优化处理,以将其去除。
具体实现流程如下:采用预设大小的扫描单元依次对所述目标变化区域内的像素点进行扫描;利用第一预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于或等于预设的第一阈值,则去除当前处于所述扫描单元中心位置的像素点;利用第二预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于预设的第二阈值,则根据当前处于所述扫描单元中心位置的像素点对所述扫描单元中的所有像素点进行像素信息重置。
下面通过一个具体实施例实现对结果图像进行优化处理的解释说明。
首先,利用卷积操作扫描整个图片,决定每一个像素点是否应该被保留,以去除小的噪音点。
确定卷积核。定义一个长度不大的正方形矩阵作为第一预设矩阵,该矩阵的取值为0或1。
具体如下:
Figure BDA0002323110060000151
卷积操作去除噪音点。对得到的图像的每一块5*5区域进行卷积操作,设定阈值为9,只有卷积结果大于9,那么对应正方形矩阵的中心点被认为是需要保留的,否则去除这个点。
其次,利用卷积操作扫描整个图片,决定每一个像素点周围空间是否应该被填充,以放大感兴趣的区域。
确定卷积核。定义一个长度不大的正方形矩阵作为第二预设矩阵,该矩阵的取值为0或1。
具体如下:
Figure BDA0002323110060000152
卷积操作填充。对得到的图像的每一块5*5区域进行卷积操作,设定阈值为5,只有卷积结果小于5,那么对应正方形矩阵的中心点被认为是需要填充的,将该矩形对应的所有点都填充,以放大感兴趣的区域。
最后,将区域进行矢量化。在经过去除与放大的步骤之后,得到的结果中每一个区域依然很可能是不规则的,为了实际使用中观看的方便,找到一块连续区域的边界点,并根据这些边界点将区域重组成一块相对规则的区域。
本发明实施例提出了根据所得原图像的色差、光照条件、拍摄角度以及所需要检测变化的物体类型,选择较好的差异图的方法,实现了根据图片条件选择尽量合适的参数。
本发明实施例提出了结合像素点本身及其邻域信息对像素点进行分类处理的流程,实现了提取像素点与其领域相关性程度的变换过程。
本发明实施例实现了根据像素点的高维信息,在无标签的条件下,对高维数据进行分类的过程。
本发明实施例提出了针对带有噪音的结果图像优化处理的后续流程,实现了利用去除和填充操作消除噪音并放大有用区域,使得实用性更高。
对于方法实施例,为了简单描述,故将其都表述为一系列的动作组合,但是本邻域技术人员应该知悉,本发明实施例并不受所描述的动作顺序的限制,因为依据本发明实施例,某些步骤可以采用其他顺序或者同时进行。其次,本邻域技术人员也应该知悉,说明书中所描述的实施例均属于优选实施例,所涉及的动作并不一定是本发明实施例所必须的。
图3示意性示出了本发明一个实施例的遥感影像变化区域检测装置的结构示意图。参照图3,本发明实施例的遥感影像变化区域检测装置具体包括配置单元201、检测单元202以及处理单元203,其中:
配置单元201,用于对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,所述像素向量用于表示当前像素点与其指定邻域范围内其他像素点的相关性;
检测单元202,用于查找各个分类类别中包含像素向量最少的目标分类类别;
处理单元203,用于将所述目标分类类别中各个像素向量对应的像素点作为目标变化区域,对所述目标变化区域进行图像优化处理。
本发明实施例,如图4所示,所述配置单元201,包括配置子单元2011和分类子单元2012,其中:
配置子单元2011,用于配置分类类别的类别数;
分类子单元2012,用于根据所述分类类别的类别数将所述差异图中各个像素点的像素向量按照预设划分标准划分为不同分类类别。
进一步地,所述分类子单元2012,具体用于计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的距离;若计算出的距离最小值大于任意两个类别点的最大距离,则当前像素点的像素向量被分为一个新的类,并将对应的两个类别合并;否则,当前像素点的像素向量被分到与其距离最小的像素向量所属的分类类别。
本实施例中,所述处理单元203,具体包括附图中未示出的扫描子单元、第一计算子单元和第二计算子单元,其中:
扫描子单元,用于采用预设大小的扫描单元依次对所述目标变化区域内的像素点进行扫描;
第一计算子单元,用于利用第一预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于或等于预设的第一阈值,则去除当前处于所述扫描单元中心位置的像素点,所述第一预设矩阵具体如下:
Figure BDA0002323110060000171
第二计算子单元,用于利用第二预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于预设的第二阈值,则根据当前处于所述扫描单元中心位置的像素点对所述扫描单元中的所有像素点进行像素信息重置,所述第二预设矩阵具体如下:
Figure BDA0002323110060000181
进一步地,所述处理单元203,还包括附图中未示出的矢量化子单元,所述矢量化子单元,用于对像素信息重置后的目标变化区域进行矢量化处理。
本发明实施例中,所述装置还包括附图中未示出的提取单元,所述提取单元,用于提取待比对的同一地区不同时相的两幅遥感影像的差异图。
本发明实施例中,所述装置还包括附图中未示出的计算单元,所述计算单元,用于计算所述差异图中各个像素点与其指定邻域范围内其他像素点的相关信息,所述相关信息通过像素点的像素向量进行表示。
对于装置实施例而言,由于其与方法实施例基本相似,所以描述的比较简单,相关之处参见方法实施例的部分说明即可。
以上所描述的装置实施例仅仅是示意性的,其中所述作为分离部件说明的单元可以是或者也可以不是物理上分开的,作为单元显示的部件可以是或者也可以不是物理单元,即可以位于一个地方,或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本实施例方案的目的。本领域普通技术人员在不付出创造性的劳动的情况下,即可以理解并实施。
本发明实施例提供的遥感影像变化区域检测方法、装置,通过对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,进而确定出用户实际感兴趣目标变化区域,实现高效、准确地从两张遥感影像的差异图中检测用户实际感兴趣的变化区域。
此外,本发明实施例还提供了一种计算机可读存储介质,其上存储有计算机程序,该程序被处理器执行时实现如上所述方法的步骤。
本实施例中,所述遥感影像变化区域检测装置集成的模块/单元如果以软件功能单元的形式实现并作为独立的产品销售或使用时,可以存储在一个计算机可读取存储介质中。基于这样的理解,本发明实现上述实施例方法中的全部或部分流程,也可以通过计算机程序来指令相关的硬件来完成,所述的计算机程序可存储于一计算机可读存储介质中,该计算机程序在被处理器执行时,可实现上述各个方法实施例的步骤。其中,所述计算机程序包括计算机程序代码,所述计算机程序代码可以为源代码形式、对象代码形式、可执行文件或某些中间形式等。所述计算机可读介质可以包括:能够携带所述计算机程序代码的任何实体或装置、记录介质、U盘、移动硬盘、磁碟、光盘、计算机存储器、只读存储器(ROM,Read-Only Memory)、随机存取存储器(RAM,Random Access Memory)、电载波信号、电信信号以及软件分发介质等。需要说明的是,所述计算机可读介质包含的内容可以根据司法管辖区内立法和专利实践的要求进行适当的增减,例如在某些司法管辖区,根据立法和专利实践,计算机可读介质不包括电载波信号和电信信号。
本发明实施例提供的电子设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,所述处理器执行所述计算机程序时实现上述各个遥感影像变化区域检测方法实施例中的步骤,例如图1所示的S11~S13。或者,所述处理器执行所述计算机程序时实现上述各遥感影像变化区域检测装置实施例中各模块/单元的功能,例如图3所示的配置单元201、检测单元202以及处理单元203。
示例性的,所述计算机程序可以被分割成一个或多个模块/单元,所述一个或者多个模块/单元被存储在所述存储器中,并由所述处理器执行,以完成本发明。所述一个或多个模块/单元可以是能够完成特定功能的一系列计算机程序指令段,该指令段用于描述所述计算机程序在所述遥感影像变化区域检测装置中的执行过程。例如,所述计算机程序可以被分割成配置单元201、检测单元202以及处理单元203。
所述电子设备可以是移动型计算机、笔记本、掌上电脑及手机等设备。所述电子设备可包括,但不仅限于,处理器、存储器。本领域技术人员可以理解,本实施例中的电子设备可以包括更多或更少的部件,或者组合某些部件,或者不同的部件,例如所述电子设备还可以包括输入输出设备、网络接入设备、总线等。
所述处理器可以是中央处理单元(Central Processing Unit,CPU),还可以是其他通用处理器、数字信号处理器(Digital Signal Processor,DSP)、专用集成电路(Application Specific Integrated Circuit,ASIC)、现成可编程门阵列(Field-Programmable Gate Array,FPGA)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件等。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等,所述处理器是所述电子设备的控制中心,利用各种接口和线路连接整个电子设备的各个部分。
所述存储器可用于存储所述计算机程序和/或模块,所述处理器通过运行或执行存储在所述存储器内的计算机程序和/或模块,以及调用存储在存储器内的数据,实现所述电子设备的各种功能。所述存储器可主要包括存储程序区和存储数据区,其中,存储程序区可存储操作系统、至少一个功能所需的应用程序(比如声音播放功能、图像播放功能等)等;存储数据区可存储根据手机的使用所创建的数据(比如音频数据、电话本等)等。此外,存储器可以包括高速随机存取存储器,还可以包括非易失性存储器,例如硬盘、内存、插接式硬盘,智能存储卡(Smart Media Card,SMC),安全数字(Secure Digital,SD)卡,闪存卡(Flash Card)、至少一个磁盘存储器件、闪存器件、或其他易失性固态存储器件。
本领域的技术人员能够理解,尽管在此的一些实施例包括其它实施例中所包括的某些特征而不是其它特征,但是不同实施例的特征的组合意味着处于本发明的范围之内并且形成不同的实施例。例如,在下面的权利要求书中,所要求保护的实施例的任意之一都可以以任意的组合方式来使用。
最后应说明的是:以上实施例仅用以说明本发明的技术方案,而非对其限制;尽管参照前述实施例对本发明进行了详细的说明,本领域的普通技术人员应当理解:其依然可以对前述各实施例所记载的技术方案进行修改,或者对其中部分技术特征进行等同替换;而这些修改或者替换,并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (11)

1.一种遥感影像变化区域检测方法,其特征在于,所述方法包括:
对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,所述像素向量用于表示当前像素点与其指定邻域范围内其他像素点的相关性;
查找各个分类类别中包含像素向量最少的目标分类类别;
将所述目标分类类别中各个像素向量对应的像素点作为目标变化区域,对所述目标变化区域进行图像优化处理;
所述对所述目标变化区域进行图像优化处理,包括:
采用预设大小的扫描单元依次对所述目标变化区域内的像素点进行扫描;
利用第一预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于或等于预设的第一阈值,则去除当前处于所述扫描单元中心位置的像素点,所述第一预设矩阵具体如下:
Figure FDA0004132574790000011
利用第二预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于预设的第二阈值,则根据当前处于所述扫描单元中心位置的像素点对所述扫描单元中的所有像素点进行像素信息重置,所述第二预设矩阵具体如下:
Figure FDA0004132574790000012
2.根据权利要求1所述的方法,其特征在于,所述对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,包括:
配置分类类别的类别数;
根据所述分类类别的类别数将所述差异图中各个像素点的像素向量按照预设划分标准划分为不同分类类别。
3.根据权利要求2所述的方法,其特征在于,所述根据分类类别的类别数将所述差异图中各个像素点的像素向量按照预设划分标准划分为不同分类类别,包括:
计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的距离;
若计算出的距离最小值大于任意两个类别点的最大距离,则当前像素点的像素向量被分为一个新的类,并将对应的两个类别合并;
否则,当前像素点的像素向量被分到与其距离最小的像素向量所属的分类类别。
4.根据权利要求3所述的方法,其特征在于,所述计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的距离,包括:
计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的L1距离、L2距离或者直方图距离。
5.根据权利要求1所述的方法,其特征在于,所述对所述目标变化区域进行图像优化处理,还包括:
对像素信息重置后的目标变化区域进行矢量化处理。
6.一种遥感影像变化区域检测装置,其特征在于,所述装置包括:
配置单元,用于对同一地区不同时相的两幅遥感影像的差异图中各个像素点的像素向量进行分类,所述像素向量用于表示当前像素点与其指定邻域范围内其他像素点的相关性;
检测单元,用于查找各个分类类别中包含像素向量最少的目标分类类别;
处理单元,用于将所述目标分类类别中各个像素向量对应的像素点作为目标变化区域,对所述目标变化区域进行图像优化处理;
所述处理单元,包括:
扫描子单元,用于采用预设大小的扫描单元依次对所述目标变化区域内的像素点进行扫描;
第一计算子单元,用于利用第一预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于或等于预设的第一阈值,则去除当前处于所述扫描单元中心位置的像素点,所述第一预设矩阵具体如下:
Figure FDA0004132574790000031
第二计算子单元,用于利用第二预设矩阵对扫描单元内的像素点进行卷积操作,若卷积结果小于预设的第二阈值,则根据当前处于所述扫描单元中心位置的像素点对所述扫描单元中的所有像素点进行像素信息重置,所述第二预设矩阵具体如下:
Figure FDA0004132574790000032
7.根据权利要求6所述的装置,其特征在于,所述配置单元,包括:
配置子单元,用于配置分类类别的类别数;
分类子单元,用于根据所述分类类别的类别数将所述差异图中各个像素点的像素向量按照预设划分标准划分为不同分类类别。
8.根据权利要求7所述的装置,其特征在于,所述分类子单元,具体用于计算当前像素点的像素向量与所述差异图中其他像素点的像素向量之间的距离;若计算出的距离最小值大于任意两个类别点的最大距离,则当前像素点的像素向量被分为一个新的类,并将对应的两个类别合并;否则,当前像素点的像素向量被分到与其距离最小的像素向量所属的分类类别。
9.根据权利要求6所述的装置,其特征在于,所述处理单元,还包括:
矢量化子单元,用于对像素信息重置后的目标变化区域进行矢量化处理。
10.一种计算机可读存储介质,其上存储有计算机程序,其特征在于,该程序被处理器执行时实现如权利要求1-5任一项所述方法的步骤。
11.一种电子设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,其特征在于,所述处理器执行所述程序时实现如权利要求1-5任一项所述方法的步骤。
CN201911306081.7A 2019-12-18 2019-12-18 遥感影像变化区域检测方法、装置、存储介质及电子设备 Active CN111192239B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911306081.7A CN111192239B (zh) 2019-12-18 2019-12-18 遥感影像变化区域检测方法、装置、存储介质及电子设备

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911306081.7A CN111192239B (zh) 2019-12-18 2019-12-18 遥感影像变化区域检测方法、装置、存储介质及电子设备

Publications (2)

Publication Number Publication Date
CN111192239A CN111192239A (zh) 2020-05-22
CN111192239B true CN111192239B (zh) 2023-04-25

Family

ID=70711008

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911306081.7A Active CN111192239B (zh) 2019-12-18 2019-12-18 遥感影像变化区域检测方法、装置、存储介质及电子设备

Country Status (1)

Country Link
CN (1) CN111192239B (zh)

Families Citing this family (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111784606B (zh) * 2020-06-30 2023-11-24 长江大学 一种遥感图像分类后处理方法、存储介质及系统
CN112634138A (zh) * 2020-12-31 2021-04-09 紫光展锐(重庆)科技有限公司 图像处理方法、装置、设备及存储介质、芯片、模组设备
CN112926532B (zh) * 2021-04-01 2024-05-10 深圳前海微众银行股份有限公司 信息处理方法、装置、设备、存储介质及计算机程序产品
CN112990112B (zh) * 2021-04-20 2021-07-27 湖南大学 边缘引导的循环卷积神经网络建筑物变化检测方法及系统
CN113591933B (zh) * 2021-07-07 2024-04-09 中国人民解放军海军航空大学 一种基于相关性度量的遥感图像变化检测方法及系统
CN115164845B (zh) * 2022-07-04 2024-04-05 新疆维吾尔自治区环境保护科学研究院 一种生态恢复遥感自动探测装置及方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US9239384B1 (en) * 2014-10-21 2016-01-19 Sandia Corporation Terrain detection and classification using single polarization SAR
CN106203521A (zh) * 2016-07-15 2016-12-07 西安电子科技大学 基于差异图自步学习的sar图像变化检测方法
CN107423771A (zh) * 2017-08-04 2017-12-01 河海大学 一种两时相遥感图像变化检测方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US9239384B1 (en) * 2014-10-21 2016-01-19 Sandia Corporation Terrain detection and classification using single polarization SAR
CN106203521A (zh) * 2016-07-15 2016-12-07 西安电子科技大学 基于差异图自步学习的sar图像变化检测方法
CN107423771A (zh) * 2017-08-04 2017-12-01 河海大学 一种两时相遥感图像变化检测方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Yu-long Qiao et.al.."Target Recognition in SAR Images via Graph Wavelet Transfrom and 2DPCA".《ICIGP 2019》.2019,全文. *
冯文卿等."高分辨率遥感影像的随机森林变化检测方法".《测绘学报》.2017,全文. *

Also Published As

Publication number Publication date
CN111192239A (zh) 2020-05-22

Similar Documents

Publication Publication Date Title
CN111179230B (zh) 遥感影像对比变化检测方法、装置、存储介质及电子设备
CN111192239B (zh) 遥感影像变化区域检测方法、装置、存储介质及电子设备
CN109886997B (zh) 基于目标检测的识别框确定方法、装置及终端设备
US10803554B2 (en) Image processing method and device
CN110660066B (zh) 网络的训练方法、图像处理方法、网络、终端设备及介质
CN110705583B (zh) 细胞检测模型训练方法、装置、计算机设备及存储介质
US9865063B2 (en) Method and system for image feature extraction
US9454836B2 (en) Object display device and object display method
CN112102164B (zh) 一种图像处理方法、装置、终端及存储介质
KR101618996B1 (ko) 호모그래피를 추정하기 위한 샘플링 방법 및 영상 처리 장치
CN109977952B (zh) 基于局部最大值的候选目标检测方法
CN111199188B (zh) 遥感影像差异图的像素处理方法、装置、存储介质及设备
CN112966548A (zh) 一种大豆地块识别方法及系统
CN117011304B (zh) 缺陷检测方法、装置、计算机设备及计算机可读存储介质
CN111709876B (zh) 一种图像拼接的方法及装置、设备、存储介质
CN112651953A (zh) 图片相似度计算方法、装置、计算机设备及存储介质
CN111199567B (zh) 车道线绘制方法、装置及终端设备
CN111127358A (zh) 图像处理方法、装置及存储介质
CN108932703B (zh) 图片处理方法、图片处理装置及终端设备
CN107770487B (zh) 一种特征提取并优选的方法、系统及终端设备
Touya et al. Comparing image-based methods for assessing visual clutter in generalized maps
CN113269752A (zh) 一种图像检测方法、装置终端设备及存储介质
CN111539975B (zh) 一种运动目标的检测方法、装置、设备及存储介质
CN111062341B (zh) 视频图像区域的分类方法、装置、设备及存储介质
CN116977671A (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