CN112347823B - 沉积相边界识别方法及装置 - Google Patents

沉积相边界识别方法及装置 Download PDF

Info

Publication number
CN112347823B
CN112347823B CN201910733457.6A CN201910733457A CN112347823B CN 112347823 B CN112347823 B CN 112347823B CN 201910733457 A CN201910733457 A CN 201910733457A CN 112347823 B CN112347823 B CN 112347823B
Authority
CN
China
Prior art keywords
pixel
super
seismic attribute
result
attribute map
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
CN201910733457.6A
Other languages
English (en)
Other versions
CN112347823A (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.)
Petrochina Co Ltd
Original Assignee
Petrochina Co Ltd
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 Petrochina Co Ltd filed Critical Petrochina Co Ltd
Priority to CN201910733457.6A priority Critical patent/CN112347823B/zh
Priority to US16/988,019 priority patent/US11327194B2/en
Publication of CN112347823A publication Critical patent/CN112347823A/zh
Application granted granted Critical
Publication of CN112347823B publication Critical patent/CN112347823B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V20/00—Scenes; Scene-specific elements
    • G06V20/10—Terrestrial scenes
    • G06V20/13—Satellite images
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/30—Analysis
    • G01V1/307—Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/30—Analysis
    • G01V1/301—Analysis for determining seismic cross-sections or geostructures
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G06F18/232—Non-hierarchical techniques
    • G06F18/2321—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions
    • G06F18/23213—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions with fixed number of clusters, e.g. K-means clustering
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/25—Fusion techniques
    • G06F18/254—Fusion techniques of classification results, e.g. of results related to same input data
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/20—Image preprocessing
    • G06V10/26—Segmentation of patterns in the image field; Cutting or merging of image elements to establish the pattern region, e.g. clustering-based techniques; Detection of occlusion
    • G06V10/267—Segmentation of patterns in the image field; Cutting or merging of image elements to establish the pattern region, e.g. clustering-based techniques; Detection of occlusion by performing operations on regions, e.g. growing, shrinking or watersheds
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/40—Extraction of image or video features
    • G06V10/56—Extraction of image or video features relating to colour
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V30/00—Character recognition; Recognising digital ink; Document-oriented image-based pattern recognition
    • G06V30/10—Character recognition
    • G06V30/24—Character recognition characterised by the processing or recognition method
    • G06V30/248—Character recognition characterised by the processing or recognition method involving plural approaches, e.g. verification by template match; Resolving confusion among similar patterns, e.g. "O" versus "Q"
    • G06V30/2504—Coarse or fine approaches, e.g. resolution of ambiguities or multiscale approaches
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/34—Displaying seismic recordings or visualisation of seismic data or attributes
    • G01V1/345—Visualisation of seismic data or attributes, e.g. in 3D cubes
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/64—Geostructures, e.g. in 3D data cubes
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/64—Geostructures, e.g. in 3D data cubes
    • G01V2210/641—Continuity of geobodies
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/66—Subsurface modeling
    • G01V2210/661—Model from sedimentation process modeling, e.g. from first principles

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Remote Sensing (AREA)
  • Theoretical Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Data Mining & Analysis (AREA)
  • Multimedia (AREA)
  • Geophysics (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • Acoustics & Sound (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Astronomy & Astrophysics (AREA)
  • Artificial Intelligence (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Probability & Statistics with Applications (AREA)
  • Image Analysis (AREA)
  • Information Retrieval, Db Structures And Fs Structures Therefor (AREA)

Abstract

本发明公开了一种沉积相边界识别方法及装置，该方法包括获取地震属性图中沉积相的初步标记结果；利用最大类间方差及K‑means聚类获取K‑means分类结果；根据SLIC超像素分割获取超像素分类结果；以初步标记结果和K‑means分类结果为约束，对超像素分类结果进行区域生长融合，确定沉积相边界的识别结果。本发明中，初步标记结果仅需少量的人工标签，SLIC超像素分割利用少量的超像素代替大量的像素，同时以初步标记结果作为约束大幅减少人工标签，提高沉积相边界识别的效率。利用最大类间方差及K‑means聚类获取K‑means分类结果提高聚类的准确性，以初步标记结果为约束提高沉积相边界识别的准确率。

Description

沉积相边界识别方法及装置
技术领域
本发明涉及数字图像处理技术领域，尤其涉及一种沉积相边界识别方法及装置。
背景技术
本部分旨在为权利要求书中陈述的本发明实施例提供背景或上下文。此处的描述不因为包括在本部分中就承认是现有技术。
沉积相是当前沉积环境的岩性特征(岩石的颜色、物质成分、结构、构造、岩石类型与组合)、古生物特征(生物的种属和生态)和地球化学特征的总和。沉积相图的绘制是油气勘探开发的重要基础工作。
现今结合地震资料的沉积相图绘制主要是用户依靠人工手段，进行手工绘制。绘制人员凭借专业知识以及地质背景判断沉积相类别，然后基于地震属性图，结合钻井、测井资料，人工绘制出相带边界。由于沉积相展布非常复杂，不同的沉积相带常相互交错，基于地震属性图手工解释沉积相边界存在周期长、主观性强、难度较大等问题。因此，解释人员更希望应用图像分割算法来提取地震数据中潜在的沉积信息。
数字图像处理技术是通过计算机对图像进行去除噪声、增强、复原、分割、提取特征等处理的方法和技术。从20世纪50年代发展至今，已经有不少算法开始应用于地质图像的分析处理当中。如LeYu在地质遥感图的层位识别中使用了SVM机器学习算法，结合ASTER图像和DEM模型，航磁数据三种数据提高了识别的正确率。Abdul-Qadir等人使用了最大似然分类算法进行岩性分类。C.A.Perez等使用SVM算法提取纹理信息进行岩性分类。同时在Yathunanthan Vasuki论文中，讨论了一种交互式岩性边界检测算法(ILBD算法)，在这种算法中，要求用户根据自己的经验判断输入初始层位样本标签，先使用超像素算法对图像进行分割，最后再基于区域生长法从标记的区域即可得到分类区域。虽然将上述算法应用于航拍图片中的岩性分类可以取得很好的效果，但应用于基于地震属性图的沉积相带划分中，会存在诸多问题。因地震属性图相对噪音较大、各区域之间无明显边界、图像结构复杂，沉积相带识别和划分精度对标签的数量及准确度要求极高。超像素块划分不准确会极大影响最终的划分效果，从而会极大增加后期沉积边界修改的工作量。虽然有很多数字图像处理技术已经应用于地质图像的分析处理当中，但是这些方法一般着重于航拍图中岩性边界的判别，对于更加复杂的沉积相展布图并没有一种高效、准确的处理方法。
因此，现有的沉积相边界识别方法存在效率低、准确率低的问题。
发明内容
本发明实施例提供一种沉积相边界识别方法，用以提高沉积相边界识别的效率和准确率，该方法包括：
获取地震属性图中沉积相的初步标记结果；
利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果；
根据SLIC超像素分割获取地震属性图的超像素分类结果；
以初步标记结果和K-means分类结果为约束，对超像素分类结果进行区域生长融合，确定地震属性图中沉积相边界的识别结果。
本发明实施例还提供一种沉积相边界识别装置，用以提高沉积相边界识别的效率和准确率，该装置包括：
初步标记结果获取模块，用于获取地震属性图中沉积相的初步标记结果；
K-means分类结果获取模块，用于利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果；
超像素分类结果获取模块，用于根据SLIC超像素分割获取地震属性图的超像素分类结果；
识别结果确定模块，用于以初步标记结果和K-means分类结果为约束，对超像素分类结果进行区域生长融合，确定地震属性图中沉积相边界的识别结果。
本发明实施例还提供一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述计算机程序时实现上述沉积相边界识别方法。
本发明实施例还提供一种计算机可读存储介质，所述计算机可读存储介质存储有执行上述沉积相边界识别方法的计算机程序。
本发明实施例中，初步标记结果仅需地质专家添加少量的人工标签，SLIC超像素分割利用少量的超像素代替大量的像素表达地震属性图的特征，同时将初步标记结果和K-means分类结果作为约束大幅减少人工标签的数量，可以极大的提高沉积相边界识别的效率。另外，利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果，提高聚类的准确性和效率；同时以地质专家的初步标记结果为约束，可以抵消数据的不准确性，极大的提高沉积相边界识别的准确率。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。在附图中：
图1为本发明实施例提供的沉积相边界识别方法的实现流程图；
图2为本发明实施例提供的沉积相边界识别方法的另一实现流程图；
图3为本发明实施例提供的沉积相边界识别方法中步骤102的实现流程图；
图4为本发明实施例提供的沉积相边界识别方法中步骤302的实现流程图；
图5为本发明实施例提供的沉积相边界识别方法中步骤103的实现流程图；
图6为本发明实施例提供的沉积相边界识别方法中步骤104的实现流程图；
图7为本发明实施例提供的沉积相边界识别方法中步骤201的实现流程图；
图8为本发明实施例提供的沉积相边界识别方法中步骤701的实现流程图；
图9为本发明实施例提供的沉积相边界识别方法中步骤703的实现流程图；
图10为本发明实施例提供的沉积相边界识别装置的功能模块图；
图11为本发明实施例提供的沉积相边界识别装置的另一功能模块图；
图12为本发明实施例提供的沉积相边界识别装置中K-means分类结果获取模块1002的结构框图；
图13为本发明实施例提供的沉积相边界识别装置中K-means分类结果获取单元1202的结构框图；
图14为本发明实施例提供的沉积相边界识别装置中超像素分类结果获取模块1003的结构框图；
图15为本发明实施例提供的沉积相边界识别装置中识别结果确定模块1004的结构框图；
图16为本发明实施例提供的沉积相边界识别装置中最终识别结果确定模块1101的结构框图；
图17为本发明实施例提供的沉积相边界识别装置中全局能量函数确定单元1601的结构框图；
图18为本发明实施例提供的沉积相边界识别装置中迭代估计单元1603的结构框图；
图19为本发明实施例提供的沉积相边界识别方法中的初步标记结果示意图；
图20为本发明实施例提供的沉积相边界识别方法中的K-means分类结果示意图；
图21为本发明实施例提供的沉积相边界识别方法中的超像素分类结果示意图；
图22为本发明实施例提供的沉积相边界识别方法中沉积相边界的识别结果示意图；
图23(a)—图23(b)为本发明实施例提供的沉积相边界识别方法中沉积相边界的最终识别结果示意图；
图24为本发明实施例提供的K-means分类结果与超像素分类结果叠加效果示意图；
图25(a)—图25(b)分别为本发明实施例提供的实际地震属性图在不使用超像素分割需要输入的人工标记和在使用超像素分割需要输入的人工标记示意图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚明白，下面结合附图对本发明实施例做进一步详细说明。在此，本发明的示意性实施例及其说明用于解释本发明，但并不作为对本发明的限定。
图1示出了本发明实施例提供的沉积相边界识别方法的实现流程，为便于描述，仅示出了与本发明实施例相关的部分，详述如下：
如图1所示，沉积相边界识别方法，其包括：
步骤101，获取地震属性图中沉积相的初步标记结果；
步骤102，利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果；
步骤103，根据SLIC超像素分割获取地震属性图的超像素分类结果；
步骤104，以初步标记结果和K-means分类结果为约束，对超像素分类结果进行区域生长融合，确定地震属性图中沉积相边界的识别结果。
沉积相是沉积物的生成环境、生成条件和其特征的总和。此处首先要求地质专家根据实际情况和通过测井及钻井实际数据得到的经验认识对不同种类的沉积相进行人工标记，添加人工标签以区分不同种类的沉积相。例如，地质专家可以使用任意长度的不连续的线段大致地对不同种类的沉积相进行标记，可以采用不同的颜色进行区分，例如不同颜色的线代表不同种类的沉积相。鉴于后续采用的K-means聚类及SLIC超像素分割等处理可以大幅减少地质专家输入的人工标签的数量，此处地质专家仅需输入少量的人工标签即可。
最大类间方差，是由日本学者大津于1979年提出的一种适合于双峰情况的自动求取阈值的方法，又叫做大津法(简称Ostu)。它是按照图像的灰度特性，将图像分成背景和目标两部分，背景和目标之间的类间方差越大说明构成图像的两部分的差别越大，当部分目标错分为背景或者背景错分为目标时都会导致两部分的差别变小。因此，使类间方差最大的分割意味着错分概率最小。
K-means聚类(英文全称：k-means clustering algorithm)，是一种迭代求解的聚类分析算法，先随机选取多个对象作为初始的聚类中心，然后把每个对象分配给距离它最近的聚类中心，对聚类中心进行更新直至聚类中心不再发生变化，此时没有对象重新分配给不同的聚类，或者误差平方和局部最小。
在本发明实施例中，利用最大类间方差对地震属性图中的沉积相自动进行分类，并将分类结果作为K-means聚类的输入，即最大类间方差的分类结果作为K-means聚类的初始聚类中心进行K-means聚类，得到基于颜色的K-means分类结果，这样可以提高K-means聚类的准确性和效率。
SLIC(英文全称：simple linear iterative clustering)超像素分割，即简单线性叠加聚类，在本发明实施例中，其实质是对地震属性图的图像像素进行局部聚类的过程。SLIC超像素分割可以形成紧凑的、近似均匀的超像素。SLIC超像素分割利用了地震属性图中像素之间的颜色和空间位置特征的相似性将像素分组，用少量的超像素代替大量的像素来表达地震属性图的特征，因此，可以在很大程度上降低地震属性图处理的复杂度，极大的提高后续区域生长融合的效率和沉积相边界识别的效率。
区域生长融合，是指将组成的像素或者区域发展成更大区域的过程。从种子点的集合开始，从这些点的区域增长是通过将与每个种子点有着相似属性，例如强度、灰度级、纹理颜色等的相邻像素合并到此区域。在本发明实施例中，是指在以初步标记结果和K-means分类结果为约束的前提下，对超像素分类结果进行区域生长融合，最终确定地震属性图中沉积相边界的识别结果。
在本发明实施例中，初步标记结果仅需地质专家添加少量的人工标签，SLIC超像素分割利用少量的超像素代替大量的像素表达地震属性图的特征，同时将初步标记结果和K-means分类结果作为约束可以大幅减少人工标签的数量，极大的提高沉积相边界识别的效率。另外，利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果，提高聚类的准确性和效率；同时以地质专家的初步标记结果为约束，可以抵消数据的不准确性，极大的提高沉积相边界识别的准确率。
图2示出了本发明实施例提供的沉积相边界识别方法的另一实现流程，为便于描述，仅示出了与本发明实施例相关的部分，详述如下：
在本发明的一实施例中，如图2所示，沉积相边界识别方法，在上述图1所示方法步骤的基础上，还包括：
步骤201，利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相边界的最终识别结果。
申请人通过进一步的研究发现，由于地震属性图结构复杂，在经过上述方法步骤的区域生长融合后得到的沉积相边界的识别结果可能会出现区域不连续(即存在缺失区域)的问题。针对此问题，申请人创造性的提出基于先验形状的偏置场变分水平集算法，将上述方法步骤得到的地质专家的初步标记结果，作为先验形状带入到水平集运算过程中，进而得到地震属性图沉积相边界的最终识别结果。利用本发明实施例的上述步骤，可以得到准确、光滑的沉积相边界。
传统的偏置场变分水平集算法通过设定一个缓慢变化的参数(即偏置场)来表示图像中灰度不均匀的部分，其与灰度不均匀的原始图像和灰度均匀的理想图像之间满足特定的关系。本发明实施例为了弥补缺失区域，在传统偏置场变分水平集算法的全局能量函数中引入先验形状以提高分割精度。其中，先验形状需要地质专家手动添加，地质专家添加的先验形状代表了地质专家认为需要填补的缺失区域。
在本发明实施例中，利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相边界的最终识别结果，将先验形状带入到水平集运算过程中，可以弥补缺失区域，得到准确、光滑的沉积相边界。
图3示出了本发明实施例提供的沉积相边界识别方法中步骤102的实现流程，为便于描述，仅示出了与本发明实施例相关的部分，详述如下：
在本发明一的实施例中，如图3所示，步骤102，利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果，包括：
步骤301，利用最大类间方差获取地震属性图的最大类间方差分类结果；
步骤302，以最大类间方差分类结果中每类沉积相的RGB均值为初始簇类中心，利用K-means聚类获取地震属性图基于颜色的K-means分类结果。
在获取地震属性图基于颜色的K-means分类结果时，首先利用最大类间方差方法获取地震属性图的最大类间方差分类结果，进而将最大类间方差分类结果作为K-means聚类的输入，即以最大类间方差分类结果中每类沉积相的RGB均值为K-means聚类的初始簇类中心，对地震属性图进行基于颜色的K-means聚类，获得基于颜色的K-means分类结果。
在本发明实施例中，首先利用最大类间方差获取地震属性图的最大类间方差分类结果，进而以最大类间方差分类结果中每类沉积相的RGB均值为初始簇类中心，利用K-means聚类获取地震属性图基于颜色的K-means分类结果，不仅可以提高分类的准确性，还可以提高分类的效率。
在本发明的一实施例中，步骤301，利用最大类间方差获取地震属性图的最大类间方差分类结果，包括：
步骤：对地震属性图进行灰度化处理获得地震属性灰度图；
步骤：分别确定地震属性灰度图中灰度值不小于初始像素阈值的像素比例和像素平均灰度；及灰度值小于初始像素阈值的像素比例和像素平均灰度；
步骤：根据灰度值不小于初始像素阈值的像素比例和像素平均灰度、灰度值小于初始像素阈值的像素比例和像素平均灰度，及地震属性灰度图的总平均灰度，确定类间方差；
步骤：在预设灰度区间中确定使得地震属性灰度图的类间方差最大的灰度值，更新灰度区间及类间方差分类数；
步骤：若更新后的类间方差分类数小于预设类间方差分类数，重复上述步骤，继续在更新后的灰度区间确定使得更新后的灰度区间对应的类间方差最大的灰度值，继续更新灰度区间及类间方差分类数；
步骤：直至更新后的类间方差分类数等于预设类间方差分类数，获得最大类间方差分类结果。
在本发明实施例中，可以采用分量法、最大值发，平均值法及加权平均法对图像进行灰度化处理。灰度图像是由RGB三个分量相同的一种特殊的彩色图像，其一个像素点的变化范围为255种，黑白图像中点的颜色深度，范围一般从0到255，白色为255，黑色为0。所谓灰度值是指色彩的浓淡程度。对地震属性图进行灰度化处理后，得到地震属性灰度图。
其中，初始像素阈值为预先设定的像素阈值。假设初始像素阈值采用T表示，在首次进行最大类间方差分类时，可以预先设定初始像素阈值T为0。在确定类间方差时，具体可以通过如下公式确定：
g＝ω0(μ0-μ)2+ω1(μ1-μ)2；
其中，g表示类间方差，ω0表示地震属性灰度图中灰度值不小于初始像素阈值的像素比例，μ0表示地震属性灰度图中灰度值不小于初始像素阈值的像素平均灰度，ω1表示地震属性灰度图中灰度值小于初始像素阈值的像素比例，μ1表示地震属性灰度图中灰度值小于初始像素阈值的像素平均灰度，μ表示地震属性灰度图的总平均灰度。
其中，预设灰度区间[0,T]可以预先设定为[0,255]，类间方差分类数与像素阈值的个数有关，在更新灰度区间及类间方差分类时，至少存在以下两种情况：
(1)在初始像素阈值T为0或255的情况下，鉴于预设灰度区间为[0,255]，意味着没有对灰度值进行分类，初始像素阈值个数(或最大类间方差迭代次数)为0，此时类间方差分类数为0+1＝1。
在确定类间方差g表达式后，继而在灰度区间[0,255]中确定使得类间方差g取得最大值时的灰度值T0，即更新后的用于分类的像素阈值，此时像素阈值为T0，像素阈值的个数(或最大类间方差迭代次数)由0个更新为1个，灰度区间[0,255]更新为灰度区间[0,T0]和灰度区间[T0,255]两个灰度区间，更新后的类间方差分类数为1+1＝2。此时，判断更新后的类间方差分类数2与预先设定的预设类间方差分类数是否相等。若更新后的新的类间方差分类数2小于预设类间方差分类数，重复上述类间方差分类步骤，在像素比例较大的灰度区间继续进行类间方差分类步骤。假设灰度区间[0,T0]的像素比例大于灰度区间[T0,255]的像素比例，继续在灰度区间[0,T0]进行类间方差分类，确定使得灰度区间[0,T0]的类间方差最大的灰度值T1，此时灰度区间[0,T0]的像素阈值为T1，像素阈值的个数(或最大类间方差迭代次数)由1个更新为2个，灰度区间[0,T0]更新为灰度区间[0,T1]和灰度区间[T1,T0]，灰度区间的个数由2个更新为2+1＝3个。即类间方差分类数由2个更新为2+1＝3个。
继续判断更新后的类间方差分类数与预先设定的预设类间方差分类数是否相等。若更新后的类间方差分类数等于预先设定的预设类间方差分类数，则终止分类过程，获得最大类间方差分类结果。
若更新后的类间方差分类数小于预先设定的预设类间方差分类数，在灰度区间[0,T1]和灰度区间[T1,T0]中像素比例较大的灰度区间继续进行类间方差分类。假设灰度区间[0,T1]的像素比例小于灰度区间[T1,T0]的像素比例，在灰度区间[T1,T0]继续进行类间方差分类，确定使得灰度区间[T1,T0]的类间方差最大的像素阈值T2，以继续更新灰度区间及类间方差分类数。直至更新后的类间方差分类数等于预先设定的预设类间方差分类数时终止分类过程，获得最大类间方差分类结果，最大类间方差分类结果将地震属性图中的沉积相分为差别最大的多类沉积相。
因此，在本发明实施例中，假设像素阈值个数(或最大类间方差迭代次数)为X，对应的类间方差分类数为Y，则像素阈值个数(或最大类间方差迭代次数)X与对应的类间方差分类数Y之间满足：Y＝X+1。
本领域技术人员可以理解的是，可以根据实际情况和具体需求预先设定类间方差分类数。例如，预先设定该预设类间方差分类数为2，或者预先设定该预设类间方差分类数为4。本领域技术人员可以理解的是，还可以预先设定该预设类间方差分类数为除上述2和4以外的其他数值，例如，3或6或9等等，本发明实施例对此不做特别的限制。
(2)在初始像素阈值T为0或255的情况下，鉴于预设灰度区间为[0,255]，意味着没有对灰度值进行分类，初始像素阈值个数(或最大类间方差迭代次数)为0，此时类间方差分类数为20＝1。
在确定类间方差g表达式后，继而在灰度区间[0,255]中确定使得类间方差g取得最大值时的灰度值T0，即更新后的用于分类的像素阈值，此时像素阈值为T0，像素阈值的个数(或最大类间方差迭代次数)由0个更新为1个，灰度区间[0,255]更新为灰度区间[0,T0]和灰度区间[T0,255]两个灰度区间，更新后的类间方差分类数为21＝2。此时，判断更新后的类间方差分类数2与预先设定的预设类间方差分类数是否相等。若更新后的新的类间方差分类数2小于预设类间方差分类数，重复上述类间方差分类步骤，在灰度区间[0,T0]和灰度区间[T0,255]继续进行类间方差分类步骤。分别确定使得灰度区间[0,T0]对应的类间方差取得最大值时的像素阈值T1，及确定使得灰度区间[T0,255]对应的类间方差取得最大值时的像素阈值T2。此时像素阈值为T1和T2，像素阈值的个数由1个(T0)更新为2个(T1和T2)，灰度区间由两个，即灰度区间[0,T0]和[T0,255]，更新为22＝4个，即灰度区间[0,T1]和灰度区间[T1,T0]，及灰度区间[T0,T2]和灰度区间[T2,255]。
继续判断更新后的类间方差分类数与预先设定的预设类间方差分类数是否相等。若更新后的类间方差分类数等于预先设定的预设类间方差分类数，则终止分类过程，获得最大类间方差分类结果。
若更新后的类间方差分类数小于预先设定的预设类间方差分类数，在灰度区间[0,T1]和灰度区间[T1,T0]，及灰度区间[T0,T2]和灰度区间[T2,255]，继续进行类间方差分类，继续更新灰度区间及类间方差分类数，直至更新后的类间方差分类数等于预先设定的预设类间方差分类数时终止分类过程，获得最大类间方差分类结果，最大类间方差分类结果将地震属性图中的沉积相分为差别最大的多类沉积相。
因此，在本发明实施例中，假设像素阈值个数(或最大类间方差迭代次数)为X，对应的类间方差分类数为Y，则像素阈值个数(或最大类间方差迭代次数)X与对应的类间方差分类数Y之间满足：Y＝2X。
本领域技术人员可以理解的是，可以根据实际情况和具体需求预先设定类间方差分类数。例如，预先设定该预设类间方差分类数为2，或者预先设定该预设类间方差分类数为4。本领域技术人员可以理解的是，还可以预先设定该预设类间方差分类数为除上述2和4以外的其他数值，例如，8和16(即2X形式)等，本发明实施例对此不做特别的限制。
在本发明一实施例中，可以利用最大类间方差分类将地震属性图分割为差别最大的两类：河道沉积相和非河道沉积相。
在本发明实施例中，不断的更新类间方差，进而获得最大类间方差分类结果，不仅可以进一步提高分类结果的准确性，还可以进一步提高分类的效率。
图4示出了本发明实施例提供的沉积相边界识别方法中步骤302的实现流程，为便于描述，仅示出了与本发明实施例相关的部分，详述如下：
在本发明一的实施例中，如图4所示，步骤302，以最大类间方差分类结果中每类沉积相的RGB均值为初始簇类中心，利用K-means聚类获取地震属性图基于颜色的K-means分类结果，包括：
步骤401，分别确定最大类间方差分类结果中每类沉积相包含的所有像素点的RGB均值，并将其作为K-means聚类的聚类中心；
步骤402，确定地震属性图中每一个像素点到所有聚类中心的距离，将该像素与距离其最近的聚类中心合并，形成新的聚类中心，直至聚类中心不再发生变化，获得地震属性图基于颜色的K-means分类结果。
在进行K-means聚类时，将最大类间方差分类结果作为K-means聚类的输入，即最大类间方差分类结果中每类沉积相的RGB均值作为K-means聚类的初始聚类中心，进而针对地震属性图中的每一个像素，分别确定每一个像素距离所有初始聚类中心的距离，将该像素与距离其最近的初始聚类中心进行合并，形成新的聚类中心，进而对聚类中心进行更新。对所有像素都进行上述如此的操作，直至最终聚类中心不再发生变化，即可获得地震属性图基于颜色的K-means分类结果。
在本发明实施例中，利用基于最大类间方差的K-means聚类对地震属性图中的沉积相进行分类，可以进一步提高沉积相分类结果的准确性。
图5示出了本发明实施例提供的沉积相边界识别方法中步骤103的实现流程，为便于描述，仅示出了与本发明实施例相关的部分，详述如下：
在本发明一的实施例中，如图5所示，步骤103，根据SLIC超像素分割获取地震属性图的超像素分类结果，包括：
步骤501，将地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心；
步骤502，遍历地震属性图中每一像素，搜索以该每一像素为中心的第一预设邻域范围内的每个超像素块中心；
步骤503，分别确定每个超像素块中心与所述每一像素之间的距离，将所述每一像素和与其距离最小的超像素块中心聚类为超像素块，直至超像素块不再变化，获得超像素分类结果。
在对地震属性图进行SLIC超像素分割时，利用规则采样网格对地震属性图进行采样，将采样结果作为初始超像素块中心。其中，该规则采样网格为间隔为S个像素的采样网格，即超像素块分块步长为S，在本发明的一实施例中，分块步长S为25，即地震属性图每隔25个像素划分一个超像素。
其中，该第一预设邻域范围为预先设定的邻域范围，本领域技术人员可以理解的是，可以根据实际情况和具体需求设定该预设邻域范围。例如，在本发明的一实施例中，预先设定该预设邻域范围为2S×2S。本领域技术人员还可以预先设定该预设邻域范围为，除上述2S×2S以外的其他邻域范围，本发明实施例对此不做特别的限制。
在以采样结果作为初始超像素块中心后，对地震属性图中的每一像素进行遍历，搜索以该每一像素为中心的2S×2S邻域范围内的每个超像素块，分别确定该每一像素与该2S×2S邻域范围内的每个超像素块中心的距离，将该每一像素和与其距离最小的超像素块进行合并，聚类为新的超像素块，以对超像素块进行更新，对地震属性图中所有的像素都进行上述操作，直至超像素块不再发生变化，即可获得超像素分类结果。
在本发明的一实施例中，可以通过如下公式确定像素与地震属性图超像素块中心之间的距离：
其中，D表示像素与地震属性图超像素块中心之间的距离，dc表示像素与地震属性图超像素块中心颜色近似度，m表示权衡颜色相似度和空间相近度的重要性参数，ds表示像素与超像素块中心空间近似度，lj表示像素j的明度，li表示超像素块中心的明度，aj表示像素j颜色a通道的值，ai表示超像素块中心颜色a通道的值，bj表示像素j颜色b通道的值，bi表示超像素块中心颜色b通道的值，xj表示像素j的横坐标，xi表示超像素块中心的横坐标，yj表示像素j的纵坐标，yi表示超像素块中心的纵坐标。
m表示权衡颜色相似度和空间相近度的重要性参数，m的值越大，表明在超像素块距离度量中空间相近度所占的比重越大，产生的超像素块越接近长方形；m的取值越小，则表明在超像素块距离度量中颜色相似度所占的比重更大，此时颜色相似度更重要，即所得的超像素块更紧密的依附于图像原有边界。在本发明的一实施例中，m一般在[1,40]范围内取值。本领域技术人员可以理解的是，可以根据实际情况和具体需求预先设定m的取值，在本发明实施例中，例如，预先设定m取值为10，本领域技术人员可以理解的是，还可以预先设定m取值为除10以外的其他数值，本发明实施例对此不做特别的限制。
在本发明实施例中，不断的更新超像素块，进而获得超像素分类结果，利用像素之间颜色和空间位置特征的相似性，用少量的超像素块代替大量的像素，来表达地震属性图的特征，可以很大程度上降低图像处理后的复杂度，极大的提高沉积相识别的效率。
本发明的一实施例中，在上述图5所示方法步骤501之前，步骤103，根据SLIC超像素分割获取地震属性图的超像素分类结果，还包括：
步骤：将地震属性图由RGB颜色空间转化为LAB颜色空间；
相应的，步骤501，将地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心，包括：
步骤：将转化为LAB颜色空间的地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心。
鉴于LAB颜色空间更接近于人的视觉感应，在LAB颜色空间上对地震属性图进行SLIC超像素分割，可以提高SLIC超像素分块的效果。
在本发明实施例中，将地震属性图由RGB颜色空间转化为LAB颜色空间，对转化为LAB颜色空间的地震属性图进行SLIC超像素分割，可以提高SLIC超像素分块的效果。
本发明的一实施例中，在上述图5所示方法步骤501之后，步骤103，根据SLIC超像素分割获取地震属性图的超像素分类结果，还包括：
步骤：将初始超像素块中心调整至其第二预设邻域范围内梯度最低的位置，获得调整后的超像素块中心；
相应的，步骤502，遍历地震属性图中每一像素，搜索以该每一像素为中心的第一预设邻域范围内的每个超像素块中心，包括：
步骤：遍历地震属性图中每一像素，搜索以该每一像素为中心的第一预设邻域范围内的每个调整后的超像素块中心。
在本发明实施例中，第二预设邻域范围为预先设定的邻域范围，本领域技术人员可以理解的是，第二预设邻域范围可以根据实际情况和具体需求预先设定该第二预设邻域范围。在本发明一实施例中，可以预先设定该第二预设邻域范围为3×3的邻域范围，还可以预先设定该第二预设邻域范围为3×3的邻域范围以外的其他邻域范围，本发明实施例对此不做特别的限制。
在本发明实施例中，将初始超像素块中心调整至其第二预设邻域范围内梯度最低的位置，可以有效避免初始超像素块中心被定为在边界上或者噪音点上。
图6示出了本发明实施例提供的沉积相边界识别方法中步骤104的实现流程，为便于描述，仅示出了与本发明实施例相关的部分，详述如下：
在本发明一的实施例中，如图6所示，步骤104，以初步标记结果和K-means分类结果为约束，对超像素分类结果进行区域生长融合，确定地震属性图中沉积相边界的识别结果，包括：
步骤601，获取每个被标记的超像素块与所有相邻的超像素块构成的超像素块集；每个被标记的超像素块和每个未被标记的超像素块通过将地震属性图的初步标记结果映射至超像素分类结果得到；
步骤602，若每个被标记的超像素块与其对应的超像素块集满足预设条件，将每个被标记的超像素块与其对应的超像素块集融合为新的超像素块；
重复上述步骤601至步骤602，直至超像素块集不再发生变化，步骤603，获得地震属性图中沉积相边界的识别结果。
在本发明实施例中，通过将地震属性图的初步标记结果映射至超像素分类结果，即可得到被标记的超像素块和没有被标记的超像素块。进而对于被标记的超像素块来说，分别确定每个被标记的超像素块与所有与该被标记的超像素块相邻的超像素块，共同构成的超像素块集，在该被标记的超像素块与其对应的超像素块集满足预设条件的情况下，将该被标记的超像素块和与其对应的超像素块集进行融合，形成新的超像素块，以此对超像素块进行更新。重复上述更新的过程，直至超像素块不再发生变化，上述过程即对超像素分类结果进行区域生长融合，最终获得的即为地震属性图中沉积相边界的识别结果。
在本发明实施例中，通过生长区域融合不断的更新超像素块，通过少量的超像素块表达地震属性图的特征，可以极大的提高沉积相识别的效率。
在本发明的一实施例中，满足预设条件包括：
与每个被标记的超像素块对应的超像素块集所包括的超像素块均为未被标记的超像素块；且
每个被标记的超像素块的颜色类别与超像素块集中心点在K-means分类结果中的颜色类别相同；且
每个被标记的超像素块的中心点与其对应的超像素块集中心点之间的距离不大于预设距离阈值。
在本发明实施例中，在将被标记的超像素块与其对应的超像素块集进行融合时，首先应当满足被标记的超像素块对应的超像素块集所包括的超像素块均为未被标记的超像素块；其次应当满足被标记的超像素块的颜色类别，与超像素块集中心点在K-means分类结果中的颜色类别相同；最后，每个被标记的超像素块的中心点与其对应的超像素块集中心点之间的距离不大于预设距离阈值。其中，预设距离阈值为预先设定的距离阈值，本领域技术人员可以理解的是，可以根据实际情况和具体需求预先设定该预设距离阈值，本发明实施例对此不做特别的限制。只有在同时满足上述三个条件时，才将被标记的超像素块与其对应的超像素块集进行生长区域融合，更新为新的超像素块。
在本发明的一实施例中，可以通过如下公式确定被标记的超像素块的中心点与对应的超像素块集中心点之间的距离：
其中，d表示被标记的地震属性图超像素块的中心点与对应的地震属性图超像素块集中心点之间的距离，l1表示被标记的地震属性图超像素块的中心点的亮度值，l2表示对应的地震属性图超像素块集中心点的亮度值，a1表示被标记的地震属性图超像素块的中心点的颜色a通道值，a2表示对应的地震属性图超像素块集中心点的亮度值，b1表示被标记的地震属性图超像素块的中心点的颜色b通道值，b2表示对应的地震属性图超像素块集中心点的亮度值。
图7示出了本发明实施例提供的沉积相边界识别方法中步骤201的实现流程，为便于描述，仅示出了与本发明实施例相关的部分，详述如下：
在本发明一的实施例中，如图7所示，步骤201，利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相的最终识别结果，包括：
步骤701，根据沉积相边界识别结果的先验形状确定基于先验形状的全局能量函数；
步骤702，根据基于先验形状的全局能量函数确定变分水平集能量泛函；
步骤703，利用迭代法估计变分水平集能量泛函最小时的最终偏置场和最终水平集函数；
步骤704，根据最终偏置场和最终水平集函数确定地震属性图沉积相的最终识别结果。
为了弥补沉积相边界识别结果中的缺失区域，本发明在传统的偏置场变分水平集算法的全局能量泛函中引入先验形状以提高分割精度。其中，先验形状为地质专家手动添加的，代表了地质专家认为需要填补的缺失区域。假设先验形状用ψ表示，下文将详述如何采用符号距离函数表示先验形状ψ。此处不再赘言。因此，本发明实施例中的全局能量函数为基于先验形状的全部能量函数，此处基于先验形状的全部能量函数可以采用E表示。
假设变分水平集能量泛函采用F表示，在确定基于先验形状的全局能量函数E后，将该基于先验形状的全局能量函数E作为数据驱动引入到变分水平集能量泛函中，可以得到：
F＝E+vL(φ)+μR(φ)；
其中，F表示变分水平集能量泛函，E表示基于先验形状的全局能量函数，φ表示水平集函数，L(φ)表示地震属性图的水平集函数零水平面的轮廓长度，作用是用于平滑演化趋势轮廓，v表示地震属性图的水平集函数零水平面的轮廓长度的权重参数，R(φ)表示距离正则项，μ表示距离正则项的权重参数。
在确定变分水平集能量泛函F的表达式后，通过迭代法不断的更新演化趋势轮廓(假设演化趋势轮廓采用c0表示)，确定在变分水平集能量泛函F取得最小值时的偏置场b和水平集函数φ，作为最终偏置场和最终水平集函数，进而根据最终偏置场和最终水平集函数确定地震属性图沉积相的最终识别结果。
在本发明实施例中，首先确定基于先验形状的全局能量函数，进而确定变分水平集函数，然后利用迭代法估计最终偏置场和最终水平集函数，进而获得地震属性图沉积相的最终识别结果，本发明实施例在传统的偏置场变分水平集算法中引入先验形状，填补了地震属性图中的缺失部分，识别得到的沉积相边界更加光滑，因而可以进一步提高地震属性图沉积相的识别结果。
图8示出了本发明实施例提供的沉积相边界识别方法中步骤701的实现流程，为便于描述，仅示出了与本发明实施例相关的部分，详述如下：
在本发明一的实施例中，如图8所示，步骤701，根据沉积相边界识别结果的先验形状确定基于先验形状的全局能量函数，包括：
步骤801，将初步结果地震属性图划分为相互正交的第一区域和第二区域，基于相互正交的第一区域和第二区域确定初步结果地震属性图的初始全局能量函数；初步结果地震属性图包含沉积相边界识别结果；
步骤802，根据先验形状确定基于先验形状的约束能量函数；
步骤803，根据地震属性图的初始全局能量函数和基于先验形状的约束能量函数，确定基于先验形状的全局能量函数。
传统的偏置场变分水平集算法通过设定一个缓慢变化的参数b(即偏置场)来表示图像中灰度不均匀的部分，I表示观测图像(原始图像)，J表示灰度均匀的理想图像，三者关系如公式(5-5)所示：
I＝b·J；
在本发明实施例地震属性图划分中，可以用b表示地震属性图中灰度不均匀的组成部分，I表示灰度不均匀的原始地震属性图，J表示灰度均匀的理想地震属性图(即理想图像)，并且缓慢变化的参数，即偏置场b和地震属性图的理想图像J存在以下的假设：
(1)地震属性图上的偏置场b缓慢变化，这意味着在一定区域范围内，可以认为偏置场b是常数，便于后续计算；
(2)地震属性图的理想图像J可以被分为n个相互正交的区域，每个区域可以近似地使用常数c1,c2,…cn-1,cn来表示每个区域的灰度值。
根据以上(1)和(2)的假设和公式I＝b·J，可以得到关于地震属性图局部特性的有用信息。考虑地震属性图中某一个半径为ρ的圆形邻域Oy，它的圆心为y∈Ω，其中Ω为整个地震属性图，定义为其中x为圆形邻域Oy内任意一点，y为圆形邻域Oy的圆心。假设偏置场b在该圆形区域Oy的范围内是缓慢变化的，所以Oy区域中偏置场b(y)和圆心的偏置场b(x)相近似，即b(x)≈b(y),x∈Oy。则在每一个正交区域中，偏置场b(x)和理想图像J的灰度值J(x)的乘积就近似等于常数b(y)ci,ci＝c1,c2,…cn-1,cn，即b(x)J(x)≈b(y)ci,x∈Oy。
因此，在区域Oy∩Ωi,Ωi＝Ω1,Ω2,…Ωn-1,Ωn中，可以以mi＝b(x)ci,i＝1,2,…n-1,n为中心将该圆形邻域Oy有效地分为n类，这样就将地震属性图的灰度信息纳入到水平集算法中，并将地震属性图划分为n个部分，每个部分的局部能量定义为：
其中，K(y-x)表示一个非负的窗口函数，即圆形区域，代表了地震属性图中的一局部区域，并且在的区域满足K(y-x)＝0，其表达公式如下：
其中，u为圆形邻域Oy内任一点x和圆心y的距离，a为一归一化常数，其使得∫K(u)＝1，σ表示高斯函数的标准差，ρ是区域Oy的半径大小。
在本发明实施例中，只需将地震属性图划分为目标区域(河道沉积相)和背景区域(非河道沉积相)，所以设定上述分类数n为2，即将地震属性图Ω划分为相互正交的第一区域Ω1和第二区域Ω2，假设水平集函数用φ表示，那么地震属性图中相互正交的第一区域Ω1和第二区域Ω2可以分别定义为：
Ω1＝{x:φ(x)＞0}；
Ω2＝{x:φ(x)＜0}；
另外，相互正交的第一区域Ω1和第二区域Ω2还可以使用成员函数表示为：
M1(φ)＝H(φ)；
M2(φ)＝1-H(φ)；
其中，H(·)表示heaviside函数。
基于上述，基于相互正交的第一区域Ω1和第二区域Ω2，确定包含沉积相边界识别结果的初步结果地震属性图的初始全局能量函数：
其中，Eimage表示初始全局能量函数，b表示偏置场，ci表示区域i的灰度信息向量。
为了填补地震属性图中沉积相识别结果的缺失区域，本发明在传统偏置场变分水平集算法的能量泛函中引入先验形状以提高分割精度，其中先验形状需地质专家手动添加，代表了地质专家认为需要填补的缺失区域。首先采用符号距离函数代表先验形状ψ：
其中d(·)函数表示地震属性图中的坐标值为(x,y)的像素点和演化趋势轮廓c0的距离。
在确定先验形状ψ的表达式后，根据先验形状ψ确定基于先验形状的约束能量函数：
Eshape＝∫Ω(H(φ)-H(ψ))2dx；
在确定初始全局能量函数Eimage和基于先验形状的约束能量函数Eshape后，根据地震属性图的初始全局能量函数Eimage和基于先验形状的约束能量函数Eshape，通过下述公式确定基于先验形状的全局能量函数：
E＝Eimage+λEshape；
其中，λ表示基于先验形状的约束能量函数Eshape的权重参数，且λ≥0。在本发明实施例中，基于相互正交的第一区域和第二区域确定初始全局能量函数，再根据先验形状确定基于先验形状的约束能量函数，进而确定基于先验形状的全局能量函数，将先验形状引入能量函数中，可以使得地震属性图中沉积相相带边界的识别结果更加光滑，进一步提高沉积相相带边界识别结果的准确性。
图9示出了本发明实施例提供的沉积相边界识别方法中步骤703的实现流程，为便于描述，仅示出了与本发明实施例相关的部分，详述如下：
在本发明一的实施例中，如图9所示，步骤703，利用迭代法估计变分水平集能量泛函最小时的最终偏置场和最终水平集函数，包括：
步骤901，通过变分法得到变分水平集能量泛函的演化趋势轮廓；
步骤902，根据初始偏置场和演化趋势轮廓的灰度信息向量，利用标准梯度下降更新水平集函数；
步骤903，根据初始偏置场和更新后的水平集函数，确定变分水平集能量泛函最小时灰度信息向量的最优值，更新灰度信息向量；
步骤904，根据更新后的水平集函数和更新后的灰度信息向量，确定变分水平集能量泛函最小时偏置场的最优值，更新偏置场；
步骤905，利用更新后的偏置场和更新后的灰度信息向量重复上述迭代，直至满足预设迭代次数，将最后一次迭代的偏置场和水平集函数分别作为最终偏置场和最终水平集函数。
在根据上述步骤702确定变分水平集能量泛函的表达式后，可以通过变分法得到变分水平集能量泛函的演化方程。其中，变分法是一种使得泛函取得极大值或者极小值的方法。本发明实施例中，变分法使得变分水平集能量泛函取得极小值。其中，演化方程表示了在每一次迭代过程中演化曲线的变化趋势，演化曲线即为演化趋势轮廓c0，可以定义为任意初始形状。在本发明实施例中，将演化趋势轮廓c0的初始形状定义为初步结果地震属性图中的任意圆形区域，通常选取图像左上角某一半径为30像素点的圆形区域，该圆形区域的半径及位置选取对最终沉积相边界识别结果没有任何影响。
通过变分法得到变分水平集能量泛函的演化方程为：
其中，t表示偏置场迭代的时间，δ(·)表示Dirac函数，即 表示梯度算子，div表示散度算子，dp(·)函数表示v表示地震属性图的水平集函数零水平面的轮廓长度的权重参数，且v≥0；μ表示距离正则项的权重参数，且μ≥0；λ表示基于先验形状的约束能量函数Eshape的权重参数。在本发明的一实施例中，设定ν＝0.06，μ＝0.001，λ＝0.001。
ei表示地震属性图的灰度约束信息，其表达式如下：
ei(x)＝∫K(y-x)|I(x)-b(y)ci|2dy,(i＝1,2)；
在利用迭代法估计变分水平集能量泛函最小时的最终偏置场和最终水平集函数时，通过变分水平集能量泛函的演化趋势轮廓c0，得到演化趋势轮廓的灰度信息向量c，进而在给定的初始偏置场b和演化趋势轮廓c0的灰度信息向量c，利用标准梯度下降更新水平集函数φ。其中，初始偏置场b为预先设定的偏置场，本领域技术人员可以根据实际情况和具体需求预先设定，本发明实施例对此不做特别的限制。另外，梯度下降是迭代法中一种用于求解最小二乘问题的方法，标准梯度下降是梯度下降的一种。
在更新水平集函数φ后，根据初始偏置场b和更新后的水平集函数φ，确定在变分水平集能量泛函取得最小值时，灰度信息向量c的最优值，进而对灰度信息向量c进行更新，更新后的灰度信息向量c的表达式为：
其中，K表示地震属性图K-means聚类中聚类中心个数。进而，在确定更新后的灰度信息向量c后，根据更新后的水平集函数φ和更新后的灰度信息向量c，确定变分水平集能量泛函取得最小值时偏置场b的最优值，进而对偏置场b进行更新，更新后的偏置场b的表达式为：
其中，i表示某个聚类区域，μi(φ)表示某个聚类区域i的隶属函数。在确定更新后的偏置场b和更新后的灰度信息向量c后，利用更新后的偏置场b和更新后的灰度信息向量c，重复上述迭代步骤902至904，直至在满足预设迭代次数时，将最后一次迭代的偏置场和水平集函数分别作为最终偏置场和最终水平集函数。
其中，预设迭代次数为预先设定的迭代次数，本领域技术人员可以理解的是，预设迭代次数可以根据实际情况和具体需求预先设定，例如，预先设定迭代次数为100次，还可以预先设定迭代次数为100次以外的其他数值，例如50次或者200次等，本发明实施例对此不做特别的限制。
在本发明实施例中，利用迭代法不断的更新水平集函数、灰度信息向量及偏置场，直至在满足迭代次数时，确定迭代后的最终偏置场和最终水平集函数，在变分水平集算法中引入先验形状，并利用迭代法估计最终偏置场和最终水平集函数，可以进一步提高沉积相边界识别结果的准确性。
在本发明的一实施例中，在上述图7所示方法步骤的基础上，步骤201，利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相的最终识别结果，还包括：
步骤：在迭代的过程中，每隔N次迭代更新先验形状为当前演化趋势轮廓；当前演化趋势轮廓为该N次迭代后的演化趋势轮廓。
在迭代的过程中，为了使得人工标注的最初的先验形状附近的曲线演化得更加光滑，并且使得人工标注的人工连接区域的曲线形状可以很好的符合地震属性图中沉积相(例如河道沉积相)之间的形状演化趋势，需要每隔N次迭代过程，就更新先验形状为当前演化得到的趋势轮廓。
在本发明实施例中，每隔N次迭代就将先验形状更新为当前由演化所得到的趋势轮廓，可以使得人工标注的人工连接区域附近的曲线演化的更加光滑，并且使得人工标注的人工连接区域的曲线形状可以很好的符合目标区域之间的形状演化趋势。
在本发明的一实施例中，在上述方法步骤的基础上，步骤201，利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相的最终识别结果，还包括：
步骤：在迭代的过程中，每隔M次迭代更新先验形状为最初的先验形状；其中，N和M均为正整数，且N小于M。
申请人经过进一步的研究发现，每隔N次更新先验形状为当前演化趋势轮廓，虽然可以使得人工标注的人工连接区域附近的曲线演化的更加光滑，但也会使得目标区域的曲线演化变得不准确。因此，为了使得目标区域的曲线演化变得更加准确，需要每隔M次迭代更新先验形状为最初的先验形状，这样可以目标区域的曲线演化的更加符合地震属性图中沉积相的边界。其中，M和N满足N小于M的关系。
在本发明的一实施例中，预先设定N的取值为10，M的取值为15。本领域技术人员可以理解的是，还可以预先设定N和M的取值为除上述N为10、M为15的其他数值，例如，N的取值为8，M的取值为12，或者N的取值为12，M的取值为20等，本发明实施例对此不做特别的限制。
在本发明实施例中，在迭代的过程中，每隔M次迭代更新先验形状为最初的先验形状，可以使得目标区域的曲线演化变得更加准确。
本发明实施例中还提供了一种沉积相边界识别装置，如下面的实施例所述。由于这些装置解决问题的原理与沉积相边界识别方法相似，因此这些装置的实施可以参见方法的实施，重复之处不再赘述。
图10示出了本发明实施例提供的沉积相边界识别装置的功能模块，为便于说明，仅示出了与本发明实施例相关的部分，详述如下：
参考图10，所述沉积相边界识别装置所包含的各个模块用于执行图1对应实施例中的各个步骤，具体请参阅图1以及图1对应实施例中的相关描述，此处不再赘述。本发明实施例中，所述沉积相边界识别装置包括初步标记结果获取模块1001、K-means分类结果获取模块1002、超像素分类结果获取模块1003及识别结果确定模块1004。
初步标记结果获取模块1001，用于获取地震属性图中沉积相的初步标记结果。
K-means分类结果获取模块1002，用于利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果。
超像素分类结果获取模块1003，用于根据SLIC超像素分割获取地震属性图的超像素分类结果。
识别结果确定模块1004，用于以初步标记结果和K-means分类结果为约束，对超像素分类结果进行区域生长融合，确定地震属性图中沉积相边界的识别结果。
在本发明实施例中，初步标记结果获取模块1001得到的初步标记结果仅需地质专家添加少量的人工标签，超像素分类结果获取模块1003中的SLIC超像素分割利用少量的超像素代替大量的像素表达地震属性图的特征，同时将初步标记结果和K-means分类结果作为约束可以大幅减少人工标签的数量，极大的提高沉积相边界识别的效率。另外，K-means分类结果获取模块1002利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果，提高聚类的准确性和效率；同时以地质专家的初步标记结果为约束，可以抵消数据的不准确性，极大的提高沉积相边界识别的准确率。
图11示出了本发明实施例提供的沉积相边界识别装置的另一功能模块，为便于说明，仅示出了与本发明实施例相关的部分，详述如下：
参考图11，所述沉积相边界识别装置所包含的各个模块用于执行图2对应实施例中的各个步骤，具体请参阅图2以及图2对应实施例中的相关描述，此处不再赘述。本发明实施例中，在上述图10所示模块结构的基础上，所述沉积相边界识别装置，还包括最终识别结果确定模块1101。
最终识别结果确定模块1101，用于利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相边界的最终识别结果。
在本发明实施例中，最终识别结果确定模块1101利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相边界的最终识别结果，将先验形状带入到水平集运算过程中，可以弥补缺失区域，得到准确、光滑的沉积相边界。
图12示出了本发明实施例提供的沉积相边界识别装置中K-means分类结果获取模块1002的结构示意，为便于说明，仅示出了与本发明实施例相关的部分，详述如下：
参考图12，所述K-means分类结果获取模块1002所包含的各个单元用于执行图3对应实施例中的各个步骤，具体请参阅图3以及图3对应实施例中的相关描述，此处不再赘述。本发明实施例中，所述K-means分类结果获取模块1002包括最大类间方差分类结果获取单元1201和K-means分类结果获取单元1202。
最大类间方差分类结果获取单元1201，用于利用最大类间方差获取地震属性图的最大类间方差分类结果；
K-means分类结果获取单元1202，用于以最大类间方差分类结果中每类沉积相的RGB均值为初始簇类中心，利用K-means聚类获取地震属性图基于颜色的K-means分类结果。
在本发明实施例中，最大类间方差分类结果获取单元1201首先利用最大类间方差获取地震属性图的最大类间方差分类结果，进而K-means分类结果获取单元1202以最大类间方差分类结果中每类沉积相的RGB均值为初始簇类中心，利用K-means聚类获取地震属性图基于颜色的K-means分类结果，不仅可以提高分类的准确性，还可以提高分类的效率。
在本发明的一实施例中，最大类间方差分类结果获取单元1201包括：地震属性灰度图获取子单元、类间方差参数确定子单元、类间方差确定子单元、类间方差更新子单元及最大类间方差分类结果获取子单元。
地震属性灰度图获取子单元，用于对地震属性图进行灰度化处理获得地震属性灰度图；
类间方差参数确定子单元，用于分别确定地震属性灰度图中灰度值不小于初始像素阈值的像素比例和像素平均灰度；及灰度值小于初始像素阈值的像素比例和像素平均灰度；
类间方差确定子单元，用于根据灰度值不小于初始像素阈值的像素比例和像素平均灰度、灰度值小于初始像素阈值的像素比例和像素平均灰度，及地震属性灰度图的总平均灰度，确定类间方差；
类间方差更新子单元，用于在预设灰度区间中确定使得地震属性灰度图的类间方差最大的灰度值，更新灰度区间及类间方差分类数；
类间方差更新子单元，还用于若更新后的类间方差分类数小于预设类间方差分类数，重复上述步骤，继续在更新后的灰度区间确定使得更新后的灰度区间对应的类间方差最大的灰度值，继续更新灰度区间及类间方差分类数；
最大类间方差分类结果获取子单元，用于直至更新后的类间方差分类数等于预设类间方差分类数，获得最大类间方差分类结果。
在本发明实施例中，类间方差更新子单元不断的更新类间方差，进而获得最大类间方差分类结果，不仅可以进一步提高分类结果的准确性，还可以进一步提高分类的效率。
图13示出了本发明实施例提供的沉积相边界识别装置中K-means分类结果获取单元1202的结构示意，为便于说明，仅示出了与本发明实施例相关的部分，详述如下：
参考图13，所述K-means分类结果获取单元1202所包含的各个子单元用于执行图4对应实施例中的各个步骤，具体请参阅图4以及图4对应实施例中的相关描述，此处不再赘述。本发明实施例中，所述K-means分类结果获取单元1202包括：聚类中心确定子单元1301和聚类中心更新子单元1302。
聚类中心确定子单元1301，分别确定最大类间方差分类结果中每类沉积相包含的所有像素点的RGB均值，并将其作为K-means聚类的聚类中心。
聚类中心更新子单元1302，用于确定地震属性图中每一个像素点到所有聚类中心的距离，将该像素与距离其最近的聚类中心合并，形成新的聚类中心，直至聚类中心不再发生变化，获得地震属性图基于颜色的K-means分类结果。
在本发明实施例中，聚类中心更新子单元1302利用基于最大类间方差的K-means聚类更新聚类中心，获得地震属性图中的沉积相进行分类，可以进一步提高沉积相分类结果的准确性。
图14示出了本发明实施例提供的沉积相边界识别装置中超像素分类结果获取模块1003的结构示意，为便于说明，仅示出了与本发明实施例相关的部分，详述如下：
参考图14，所述超像素分类结果获取模块1003所包含的各个单元用于执行图5对应实施例中的各个步骤，具体请参阅图5以及图5对应实施例中的相关描述，此处不再赘述。本发明实施例中，所述超像素分类结果获取模块1003包括采样单元1401、遍历搜索单元1402及超像素块更新单元1403。
采样单元1401，用于将地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心。
遍历搜索单元1402，用于遍历地震属性图中每一像素，搜索以该每一像素为中心的第一预设邻域范围内的每个超像素块中心。
超像素块更新单元1403，用于分别确定每个超像素块中心与所述每一像素之间的距离，将所述每一像素和与其距离最小的超像素块中心聚类为超像素块，直至超像素块不再变化，获得超像素分类结果。
在本发明实施例中，超像素块更新单元1403不断的更新超像素块，进而获得超像素分类结果，利用像素之间颜色和空间位置特征的相似性，用少量的超像素块代替大量的像素，来表达地震属性图的特征，可以很大程度上降低图像处理后的复杂度，极大的提高沉积相识别的效率。
本发明的一实施例中，在上述图14所述结构示意的基础上，超像素分类结果获取模块1003还包括：
转化单元，用于将地震属性图由RGB颜色空间转化为LAB颜色空间；
相应的，采样单元1401，具体用于将转化为LAB颜色空间的地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心。
在本发明实施例中，转化单元将地震属性图由RGB颜色空间转化为LAB颜色空间，对转化为LAB颜色空间的地震属性图进行SLIC超像素分割，可以提高SLIC超像素分块的效果。
本发明的一实施例中，在上述所述结构示意的基础上，超像素分类结果获取模块1003还包括：
调整单元，用于将初始超像素块中心调整至其第二预设邻域范围内梯度最低的位置，获得调整后的超像素块中心。
相应的，遍历搜索单元1402，具体用于遍历地震属性图中每一像素，搜索以该每一像素为中心的第一预设邻域范围内的每个调整后的超像素块中心。
在本发明实施例中，调整单元将初始超像素块中心调整至其第二预设邻域范围内梯度最低的位置，可以有效避免初始超像素块中心被定为在边界上或者噪音点上。
图15示出了本发明实施例提供的沉积相边界识别装置中识别结果确定模块1004的结构示意，为便于说明，仅示出了与本发明实施例相关的部分，详述如下：
参考图15，所述识别结果确定模块1004所包含的各个单元用于执行图6对应实施例中的各个步骤，具体请参阅图6以及图6对应实施例中的相关描述，此处不再赘述。本发明实施例中，所述识别结果确定模块1004包括超像素块集获取单元1501、超像素块集更新单元1502及识别结果确定单元1503。
超像素块集获取单元1501，用于获取每个被标记的超像素块与所有相邻的超像素块构成的超像素块集；每个被标记的超像素块和每个未被标记的超像素块通过将地震属性图的初步标记结果映射至超像素分类结果得到。
超像素块集更新单元1502，用于若每个被标记的超像素块与其对应的超像素块集满足预设条件，将每个被标记的超像素块与其对应的超像素块集融合为新的超像素块；
识别结果确定单元1503，用于在超像素块集不再发生变化时，获得地震属性图中沉积相边界的识别结果。
在本发明实施例中，超像素块集更新单元1502通过生长区域融合不断的更新超像素块，通过少量的超像素块表达地震属性图的特征，可以极大的提高沉积相识别的效率。
在本发明的一实施例中，满足预设条件包括：
与每个被标记的超像素块对应的超像素块集所包括的超像素块均为未被标记的超像素块；且
每个被标记的超像素块的颜色类别与超像素块集中心点在K-means分类结果中的颜色类别相同；且
每个被标记的超像素块的中心点与其对应的超像素块集中心点之间的距离不大于预设距离阈值。
图16示出了本发明实施例提供的沉积相边界识别装置中最终识别结果确定模块1101的结构示意，为便于说明，仅示出了与本发明实施例相关的部分，详述如下：
参考图16，所述最终识别结果确定模块1101所包含的各个单元用于执行图7对应实施例中的各个步骤，具体请参阅图7以及图7对应实施例中的相关描述，此处不再赘述。本发明实施例中，所述最终识别结果确定模块1101包括全局能量函数确定单元1601、变分水平集能量泛函确定单元1602、迭代估计单元1603及最终识别结果确定单元1604。
全局能量函数确定单元1601，用于根据沉积相边界识别结果的先验形状确定基于先验形状的全局能量函数。
变分水平集能量泛函确定单元1602，用于根据基于先验形状的全局能量函数确定变分水平集能量泛函。
迭代估计单元1603，用于利用迭代法估计变分水平集能量泛函最小时的最终偏置场和最终水平集函数。
最终识别结果确定单元1604，用于根据最终偏置场和最终水平集函数确定地震属性图沉积相的最终识别结果。
在本发明实施例中，全局能量函数确定单元1601首先确定基于先验形状的全局能量函数，进而变分水平集能量泛函确定单元1602确定变分水平集函数，然后迭代估计单元1603利用迭代法估计最终偏置场和最终水平集函数，进而最终识别结果确定单元1604获得地震属性图沉积相的最终识别结果，本发明实施例在传统的偏置场变分水平集算法中引入先验形状，填补了地震属性图中的缺失部分，识别得到的沉积相边界更加光滑，因而可以进一步提高地震属性图沉积相的识别结果。
图17示出了本发明实施例提供的沉积相边界识别装置中全局能量函数确定单元1601的结构示意，为便于说明，仅示出了与本发明实施例相关的部分，详述如下：
参考图17，所述全局能量函数确定单元1601所包含的各个子单元用于执行图8对应实施例中的各个步骤，具体请参阅图8以及图8对应实施例中的相关描述，此处不再赘述。本发明实施例中，所述全局能量函数确定单元1601包括初始全局能量函数确定子单元1701、约束能量函数确定子单元1702及全局能量函数确定子单元1703。
初始全局能量函数确定子单元1701，用于将初步结果地震属性图划分为相互正交的第一区域和第二区域，基于相互正交的第一区域和第二区域确定初步结果地震属性图的初始全局能量函数；初步结果地震属性图包含沉积相边界识别结果。
约束能量函数确定子单元1702，用于根据先验形状确定基于先验形状的约束能量函数。
全局能量函数确定子单元1703，用于根据地震属性图的初始全局能量函数和基于先验形状的约束能量函数，确定基于先验形状的全局能量函数。
在本发明实施例中，初始全局能量函数确定子单元1701基于相互正交的第一区域和第二区域确定初始全局能量函数，约束能量函数确定子单元1702再根据先验形状确定基于先验形状的约束能量函数，进而全局能量函数确定子单元1703确定基于先验形状的全局能量函数，将先验形状引入能量函数中，可以使得地震属性图中沉积相相带边界的识别结果更加光滑，进一步提高沉积相相带边界识别结果的准确性。
图18示出了本发明实施例提供的沉积相边界识别装置中迭代估计单元1603的结构示意，为便于说明，仅示出了与本发明实施例相关的部分，详述如下：
参考图18，所述迭代估计单元1603所包含的各个子单元用于执行图9对应实施例中的各个步骤，具体请参阅图9以及图9对应实施例中的相关描述，此处不再赘述。本发明实施例中，所述迭代估计单元1603包括演化趋势轮廓迭代更新子单元1801、水平集函数迭代更新子单元1802、灰度信息向量迭代更新子单元1803、偏置场迭代更新子单元1804及确定子单元1805。
演化趋势轮廓迭代更新子单元1801，用于通过变分法得到变分水平集能量泛函的演化趋势轮廓。
水平集函数迭代更新子单元1802，用于根据初始偏置场和演化趋势轮廓的灰度信息向量，利用标准梯度下降更新水平集函数。
灰度信息向量迭代更新子单元1803，用于根据初始偏置场和更新后的水平集函数，确定变分水平集能量泛函最小时灰度信息向量的最优值，更新灰度信息向量。
偏置场迭代更新子单元1804，用于根据更新后的水平集函数和更新后的灰度信息向量，确定变分水平集能量泛函最小时偏置场的最优值，更新偏置场。
确定子单元1805，用于利用更新后的偏置场和更新后的灰度信息向量重复上述迭代，直至满足预设迭代次数，将最后一次迭代的偏置场和水平集函数分别作为最终偏置场和最终水平集函数。
在本发明实施例中，水平集函数迭代更新子单元1802利用迭代法不断的更新水平集函数、灰度信息向量迭代更新子单元1803不断的更新灰度信息向量及偏置场迭代更新子单元1804不断的更新偏置场，直至在满足迭代次数时，确定子单元1805确定迭代后的最终偏置场和最终水平集函数，在变分水平集算法中引入先验形状，并利用迭代法估计最终偏置场和最终水平集函数，可以进一步提高沉积相边界识别结果的准确性。
在本发明的一实施例中，在上述图16所示结构示意的基础上，最终识别结果还包括：
第一更新单元，用于在迭代的过程中，每隔N次迭代更新先验形状为当前演化趋势轮廓；当前演化趋势轮廓为该N次迭代后的演化趋势轮廓。
在本发明实施例中，第一更新单元每隔N次迭代更新先验形状为当前演化趋势轮廓，可以使得人工标注的人工连接区域附近的曲线演化的更加光滑，并且使得人工标注的人工连接区域的曲线形状可以很好的符合目标区域之间的形状演化趋势。
在本发明的一实施例中，在上述所示结构示意的基础上，最终识别结果还包括：
第二更新单元，用于在迭代的过程中，每隔M次迭代更新先验形状为最初的先验形状；其中，N和M均为正整数，且N小于M。
在本发明实施例中，在迭代的过程中，第二更新单元每隔M次迭代更新先验形状为最初的先验形状，可以使得目标区域的曲线演化变得更加准确。
为更加清楚的说明本发明提供的沉积相边界识别的原理和效果，以下实施例以实际地震属性图为例进行说明。
图19示出了本发明实施例提供的实际地震属性图的初步标记结果示意，为便于说明，仅示出了与本发明实施例现相关的部分，详述如下：
如图19所示，图中所示任意长度的不连续的线段为地质专家基于经验认识输入的多个人工标签，多个人工标签用于对不同种类的沉积相进行大致的区分和标记。例如，可以用不同颜色的线段代表不同类型的沉积相。
图20示出了本发明实施例提供的实际地震属性图的K-means分类结果示意，为便于说明，仅示出了与本发明实施例现相关的部分，详述如下：
如图20所示，图中所示为经过K-means聚类得到的K-means分类结果。在本发明实施例中，K-means分类结果包含两类差别较大的沉积相：河道沉积相(目标区域)和非河道沉积相(背景区域)。
图21示出了本发明实施例提供的实际地震属性图的超像素分类结果示意，为便于说明，仅示出了与本发明实施例现相关的部分，详述如下：
如图21所示，图中所示为经过SLIC超像素分割后得到的超像素分类结果。从图中可以看出，每个超像素块都包含大量的像素点，少量的超像素块可以替代大量的像素点，表达整个地震属性图的特征，这样可以极大的降低沉积相识别的复杂度，提高沉积相识别的效率。
图22示出了本发明实施例提供的实际地震属性图的沉积相边界的识别结果示意，为便于说明，仅示出了与本发明实施例现相关的部分，详述如下：
如图22所示，图中示出了实际地震属性图的沉积相边界的识别结果，从图中可以看出，地震属性图沉积相的初步识别结果中，已经能够较好的识别地震属性图中的河道沉积相和非河道沉积相，但仍然存在少量缺失的不连续区域，例如图中圆圈所示位置。
图23(a)示出了本发明实施例提供的实际地震属性图的沉积相边界的最终识别结果示意，为便于说明，仅示出了与本发明实施例现相关的部分，详述如下：
如图23(a)所示，在偏置场变分水平集算法中引入先验形状后，填补了地震属性图中河道沉积相的缺失、不连续部分，使得地震属性图最终识别结果中沉积相的边界更加连续和光滑。
图23(b)示出了本发明实施例提供的实际地震属性图的沉积相边界的最终识别结果示意，为便于说明，仅示出了与本发明实施例现相关的部分，详述如下：
如图23(b)所示，图示为将图23(a)所示的地震属性图最终识别结果中沉积相的边界，提取并输出的沉积相边界，从图中可以看出，最终输出的沉积相边界光滑、连续。
图24示出了本发明实施例提供的实际地震属性图的K-means分类结果与超像素分类结果叠加效果示意，为便于说明，仅示出了与本发明实施例现相关的部分，详述如下：
如图24所示，图示为将K-means分类结果映射在超像素分类结果上，得到的叠加效果示意。地震属性图像素点的数量级为10万，从图中可以看出，在经过K-means聚类和SLIC超像素分割后，将地震属性图缩减为约1000个左右的超像素块，可以极大的提高生长区域融合阶段的效率，进而提高沉积相边界的识别效率。
图25(a)示出了本发明实施例提供的实际地震属性图在不使用超像素分割需要输入的人工标记示意，图25(b)示出了本发明实施例提供的实际地震属性图在使用超像素分割需要输入的人工标记示意，为便于说明，仅示出了与本发明实施例现相关的部分，详述如下：
如图25(a)及图25(b)所示，通过对比可以发现，图25(b)中使用SLIC超像素分割后需要输入的人工标记的数量，相较于图25(a)不使用SLIC超像素分割后需要输入的人工标记的数量要大大减少。由此可见，对地震属性图使用SLIC超像素分割算法，可以大幅减少地质专家输入的人工标记的数量，提高沉积相识别的效率。
本发明实施例还提供一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述计算机程序时实现上述沉积相边界识别方法。
本发明实施例还提供一种计算机可读存储介质，所述计算机可读存储介质存储有执行上述沉积相边界识别方法的计算机程序。
综上所述，本发明实施例中，初步标记结果仅需地质专家添加少量的人工标签，SLIC超像素分割利用少量的超像素代替大量的像素表达地震属性图的特征，同时将初步标记结果和K-means分类结果作为约束大幅减少人工标签的数量，可以极大的提高沉积相边界识别的效率。另外，利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果，提高聚类的准确性和效率；同时以地质专家的初步标记结果为约束，可以抵消数据的不准确性，极大的提高沉积相边界识别的准确率。
本领域内的技术人员应明白，本发明的实施例可提供为方法、系统、或计算机程序产品。因此，本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (20)

1.一种沉积相边界识别方法，其特征在于，包括：
获取地震属性图中沉积相的初步标记结果；
利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果；
根据SLIC超像素分割获取地震属性图的超像素分类结果；
以初步标记结果和K-means分类结果为约束，对超像素分类结果进行区域生长融合，确定地震属性图中沉积相边界的识别结果；
其中，利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果，包括：
利用最大类间方差获取地震属性图的最大类间方差分类结果；
以最大类间方差分类结果中每类沉积相的RGB均值为初始簇类中心，利用K-means聚类获取地震属性图基于颜色的K-means分类结果；
其中，以初步标记结果和K-means分类结果为约束，对超像素分类结果进行区域生长融合，确定地震属性图中沉积相边界的识别结果，包括：
获取每个被标记的超像素块与所有相邻的超像素块构成的超像素块集；每个被标记的超像素块和每个未被标记的超像素块通过将地震属性图的初步标记结果映射至超像素分类结果得到；
若每个被标记的超像素块与其对应的超像素块集满足预设条件，将每个被标记的超像素块与其对应的超像素块集融合为新的超像素块；
重复上述步骤，直至超像素块集不再发生变化时，获得地震属性图中沉积相边界的识别结果；
满足预设条件包括：
与每个被标记的超像素块对应的超像素块集所包括的超像素块均为未被标记的超像素块；且
每个被标记的超像素块的颜色类别与超像素块集中心点在K-means分类结果中的颜色类别相同；且
每个被标记的超像素块的中心点与其对应的超像素块集中心点之间的距离不大于预设距离阈值。
2.如权利要求1所述的沉积相边界识别方法，其特征在于，还包括：
利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相边界的最终识别结果。
3.如权利要求1所述的沉积相边界识别方法，其特征在于，利用最大类间方差获取地震属性图的最大类间方差分类结果，包括：
对地震属性图进行灰度化处理获得地震属性灰度图；
分别确定地震属性灰度图中灰度值不小于初始像素阈值的像素比例和像素平均灰度；及灰度值小于初始像素阈值的像素比例和像素平均灰度；
根据灰度值不小于初始像素阈值的像素比例和像素平均灰度、灰度值小于初始像素阈值的像素比例和像素平均灰度，及地震属性灰度图的总平均灰度，确定类间方差；
在预设灰度区间中确定使得地震属性灰度图的类间方差最大的灰度值，更新灰度区间及类间方差分类数；
若更新后的类间方差分类数小于预设类间方差分类数，重复上述步骤，继续在更新后的灰度区间确定使得更新后的灰度区间对应的类间方差最大的灰度值，继续更新灰度区间及类间方差分类数；
直至更新后的类间方差分类数等于预设类间方差分类数，获得最大类间方差分类结果。
4.如权利要求1或3所述的沉积相边界识别方法，其特征在于，以最大类间方差分类结果中每类沉积相的RGB均值为初始簇类中心，利用K-means聚类获取地震属性图基于颜色的K-means分类结果，包括：
分别确定最大类间方差分类结果中每类沉积相包含的所有像素点的RGB均值，并将其作为K-means聚类的聚类中心；
确定地震属性图中每一个像素点到所有聚类中心的距离，将该像素与距离其最近的聚类中心合并，形成新的聚类中心，直至聚类中心不再发生变化，获得地震属性图基于颜色的K-means分类结果。
5.如权利要求1所述的沉积相边界识别方法，其特征在于，根据SLIC超像素分割获取地震属性图的超像素分类结果，包括：
将地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心；
遍历地震属性图中每一像素，搜索以该每一像素为中心的第一预设邻域范围内的每个超像素块中心；
分别确定每个超像素块中心与所述每一像素之间的距离，将所述每一像素和与其距离最小的超像素块中心聚类为超像素块，直至超像素块不再变化，获得超像素分类结果。
6.如权利要求5所述的沉积相边界识别方法，其特征在于，在将地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心之前，还包括：
将地震属性图由RGB颜色空间转化为LAB颜色空间；
相应的，将地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心，包括：
将转化为LAB颜色空间的地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心。
7.如权利要求5或6所述的沉积相边界识别方法，其特征在于，在将地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心之后，还包括：
将初始超像素块中心调整至其第二预设邻域范围内梯度最低的位置，获得调整后的超像素块中心；
相应的，遍历地震属性图中每一像素，搜索以该每一像素为中心的第一预设邻域范围内的每个超像素块中心，包括：
遍历地震属性图中每一像素，搜索以该每一像素为中心的第一预设邻域范围内的每个调整后的超像素块中心。
8.如权利要求2所述的沉积相边界识别方法，其特征在于，利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相的最终识别结果，包括：
根据沉积相边界识别结果的先验形状确定基于先验形状的全局能量函数；
根据基于先验形状的全局能量函数确定变分水平集能量泛函；
利用迭代法估计变分水平集能量泛函最小时的最终偏置场和最终水平集函数；
根据最终偏置场和最终水平集函数确定地震属性图沉积相的最终识别结果。
9.如权利要求8所述的沉积相边界识别方法，其特征在于，根据沉积相边界识别结果的先验形状确定基于先验形状的全局能量函数，包括：
将初步结果地震属性图划分为相互正交的第一区域和第二区域，基于相互正交的第一区域和第二区域确定初步结果地震属性图的初始全局能量函数；初步结果地震属性图包含沉积相边界识别结果；
根据先验形状确定基于先验形状的约束能量函数；
根据地震属性图的初始全局能量函数和基于先验形状的约束能量函数，确定基于先验形状的全局能量函数。
10.如权利要求8所述的沉积相边界识别方法，其特征在于，利用迭代法估计变分水平集能量泛函最小时的最终偏置场和最终水平集函数，包括：
通过变分法得到变分水平集能量泛函的演化趋势轮廓；
根据初始偏置场和演化趋势轮廓的灰度信息向量，利用标准梯度下降更新水平集函数；
根据初始偏置场和更新后的水平集函数，确定变分水平集能量泛函最小时灰度信息向量的最优值，更新灰度信息向量；
根据更新后的水平集函数和更新后的灰度信息向量，确定变分水平集能量泛函最小时偏置场的最优值，更新偏置场；
利用更新后的偏置场和更新后的灰度信息向量重复上述迭代，直至满足预设迭代次数，将最后一次迭代的偏置场和水平集函数分别作为最终偏置场和最终水平集函数。
11.如权利要求8所述的沉积相边界识别方法，其特征在于，还包括：
在迭代的过程中，每隔N次迭代更新先验形状为当前演化趋势轮廓；当前演化趋势轮廓为该N次迭代后的演化趋势轮廓。
12.如权利要求11所述的沉积相边界识别方法，其特征在于，还包括：
在迭代的过程中，每隔M次迭代更新先验形状为最初的先验形状；其中，N和M均为正整数，且N小于M。
13.一种沉积相边界识别装置，其特征在于，包括：
初步标记结果获取模块，用于获取地震属性图中沉积相的初步标记结果；
K-means分类结果获取模块，用于利用最大类间方差及K-means聚类获取地震属性图基于颜色的K-means分类结果；
超像素分类结果获取模块，用于根据SLIC超像素分割获取地震属性图的超像素分类结果；
识别结果确定模块，用于以初步标记结果和K-means分类结果为约束，对超像素分类结果进行区域生长融合，确定地震属性图中沉积相边界的识别结果；
其中，K-means分类结果获取模块包括：
最大类间方差分类结果获取单元，用于利用最大类间方差获取地震属性图的最大类间方差分类结果；
K-means分类结果获取单元，用于以最大类间方差分类结果中每类沉积相的RGB均值为初始簇类中心，利用K-means聚类获取地震属性图基于颜色的K-means分类结果；
其中，识别结果确定模块包括：
超像素块集获取单元，用于获取每个被标记的超像素块与所有相邻的超像素块构成的超像素块集；每个被标记的超像素块和每个未被标记的超像素块通过将地震属性图的初步标记结果映射至超像素分类结果得到；
超像素块集更新单元，用于若每个被标记的超像素块与其对应的超像素块集满足预设条件，将每个被标记的超像素块与其对应的超像素块集融合为新的超像素块；
识别结果确定单元，用于在超像素块集不再发生变化时，获得地震属性图中沉积相边界的识别结果；
满足预设条件包括：
与每个被标记的超像素块对应的超像素块集所包括的超像素块均为未被标记的超像素块；且
每个被标记的超像素块的颜色类别与超像素块集中心点在K-means分类结果中的颜色类别相同；且
每个被标记的超像素块的中心点与其对应的超像素块集中心点之间的距离不大于预设距离阈值。
14.如权利要求13所述的沉积相边界识别装置，其特征在于，还包括：
最终识别结果确定模块，用于利用基于先验形状的偏置场变分水平集算法对沉积相边界的识别结果进行追踪，确定地震属性图沉积相边界的最终识别结果。
15.如权利要求13所述的沉积相边界识别装置，其特征在于，超像素分类结果获取模块包括：
采样单元，用于将地震属性图在间隔S个像素的规则采样网格上的采样结果作为初始超像素块中心；
遍历搜索单元，用于遍历地震属性图中每一像素，搜索以该每一像素为中心的第一预设邻域范围内的每个超像素块中心；
超像素块更新单元，用于分别确定每个超像素块中心与所述每一像素之间的距离，将所述每一像素和与其距离最小的超像素块中心聚类为超像素块，直至超像素块不再变化，获得超像素分类结果。
16.如权利要求14所述的沉积相边界识别装置，其特征在于，最终识别结果确定模块包括：
全局能量函数确定单元，用于根据沉积相边界识别结果的先验形状确定基于先验形状的全局能量函数；
变分水平集能量泛函确定单元，用于根据基于先验形状的全局能量函数确定变分水平集能量泛函；
迭代估计单元，用于利用迭代法估计变分水平集能量泛函最小时的最终偏置场和最终水平集函数；
最终识别结果确定单元，用于根据最终偏置场和最终水平集函数确定地震属性图沉积相的最终识别结果。
17.如权利要求16所述的沉积相边界识别装置，其特征在于，全局能量函数确定单元包括：
初始全局能量函数确定子单元，用于将初步结果地震属性图划分为相互正交的第一区域和第二区域，基于相互正交的第一区域和第二区域确定初步结果地震属性图的初始全局能量函数；初步结果地震属性图包含沉积相边界识别结果；
约束能量函数确定子单元，用于根据先验形状确定基于先验形状的约束能量函数；
全局能量函数确定子单元，用于根据地震属性图的初始全局能量函数和基于先验形状的约束能量函数，确定基于先验形状的全局能量函数。
18.如权利要求16所述的沉积相边界识别装置，其特征在于，迭代估计单元包括：
演化趋势轮廓迭代更新单元，用于通过变分法得到变分水平集能量泛函的演化趋势轮廓；
水平集函数迭代更新单元，用于根据初始偏置场和演化趋势轮廓的灰度信息向量，利用标准梯度下降更新水平集函数；
灰度信息向量迭代更新单元，用于根据初始偏置场和更新后的水平集函数，确定变分水平集能量泛函最小时灰度信息向量的最优值，更新灰度信息向量；
偏置场迭代更新单元，用于根据更新后的水平集函数和更新后的灰度信息向量，确定变分水平集能量泛函最小时偏置场的最优值，更新偏置场；
确定单元，用于利用更新后的偏置场和更新后的灰度信息向量重复上述迭代，直至满足预设迭代次数，将最后一次迭代的偏置场和水平集函数分别作为最终偏置场和最终水平集函数。
19.一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述处理器执行所述计算机程序时实现权利要求1至12任一所述方法。
20.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储有计算机程序，所述计算机程序被处理器执行时实现权利要求1至12任一所述方法。
CN201910733457.6A 2019-08-09 2019-08-09 沉积相边界识别方法及装置 Active CN112347823B (zh)

Priority Applications (2)

Application Number Priority Date Filing Date Title
CN201910733457.6A CN112347823B (zh) 2019-08-09 2019-08-09 沉积相边界识别方法及装置
US16/988,019 US11327194B2 (en) 2019-08-09 2020-08-07 Method for identifying boundary of sedimentary facies, computer device and computer readable storage medium

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910733457.6A CN112347823B (zh) 2019-08-09 2019-08-09 沉积相边界识别方法及装置

Publications (2)

Publication Number Publication Date
CN112347823A CN112347823A (zh) 2021-02-09
CN112347823B true CN112347823B (zh) 2024-05-03

Family

ID=74367643

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910733457.6A Active CN112347823B (zh) 2019-08-09 2019-08-09 沉积相边界识别方法及装置

Country Status (2)

Country Link
US (1) US11327194B2 (zh)
CN (1) CN112347823B (zh)

Families Citing this family (19)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US11321915B2 (en) * 2020-06-12 2022-05-03 Chevron U.S.A. Inc. Conversion of rock classification maps for modeling
CN112883898A (zh) * 2021-03-11 2021-06-01 中国科学院空天信息创新研究院 一种基于sar影像的地物分类方法及装置
CN113139302B (zh) * 2021-05-20 2023-03-10 电子科技大学 基于区域生长的断溶体识别方法
CN113419278B (zh) * 2021-06-21 2022-03-08 大庆油田有限责任公司 一种基于状态空间模型与支持向量回归的井震联合多目标同时反演方法
CN113780259B (zh) * 2021-11-15 2022-03-15 中移(上海)信息通信科技有限公司 一种路面缺陷检测方法、装置、电子设备和可读存储介质
US20230220770A1 (en) * 2022-01-07 2023-07-13 Schlumberger Technology Corporation Systems and methods for measuring physical lithological features based on calibrated photographs of rock particles
CN115018828B (zh) * 2022-08-03 2022-10-25 深圳市尹泰明电子有限公司 一种电子元器件的缺陷检测方法
CN115578732B (zh) * 2022-11-21 2023-03-28 山东爱福地生物股份有限公司 一种肥料生产线的标签识别方法
CN115661669B (zh) * 2022-12-13 2023-03-21 山东省土地发展集团有限公司 一种基于视频监控的耕地违法占用监测方法及系统
CN116109662B (zh) * 2023-04-13 2023-06-23 中国科学院国家空间科学中心 一种红外图像的超像素分割方法
CN116452613B (zh) * 2023-06-14 2023-08-29 山东省国土空间生态修复中心(山东省地质灾害防治技术指导中心、山东省土地储备中心) 一种地质调查中裂缝轮廓提取方法
CN116823826B (zh) * 2023-08-29 2023-11-03 无锡康贝电子设备有限公司 一种数控机床刀具崩刃异常检测方法
CN116993947B (zh) * 2023-09-26 2023-12-12 光谷技术有限公司 一种三维场景可视化展示方法及系统
CN117689657B (zh) * 2024-02-02 2024-04-16 深圳市福源晖科技有限公司 一种柔性线路板贴片异常检测方法及系统
CN117690030B (zh) * 2024-02-02 2024-04-26 陕西仙喜辣木茯茶有限公司 基于图像处理的多面发花识别方法及系统
CN117853493B (zh) * 2024-03-08 2024-05-24 山东天意机械股份有限公司 一种基于图像处理的墙板生产质量检测方法
CN118298123B (zh) * 2024-06-04 2024-09-13 中国科学院地质与地球物理研究所 一种大规模数据耦合三维地质精细建模方法、系统
CN118334656B (zh) * 2024-06-12 2024-08-20 济南市计量检定测试院 基于机器视觉的计量检测数据识别方法
CN118395210B (zh) * 2024-07-01 2024-09-10 广州海洋地质调查局三亚南海地质研究所 一种微体古生物沉积环境识别方法、装置、设备及介质

Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102053270A (zh) * 2009-10-30 2011-05-11 中国石油化工股份有限公司 一种基于沉积地层单元的地震相分析方法
CN106056155A (zh) * 2016-05-30 2016-10-26 西安电子科技大学 基于边界信息融合的超像素分割方法
CN106545337A (zh) * 2016-11-25 2017-03-29 西南石油大学 一种基于支持向量机的沉积微相测井识别方法
CN107229917A (zh) * 2017-05-31 2017-10-03 北京师范大学 一种基于迭代聚类的多幅遥感影像共性显著目标检测方法
CN107545577A (zh) * 2017-08-23 2018-01-05 电子科技大学 基于神经网络的沉积相图像分割方法
CN108846381A (zh) * 2018-07-06 2018-11-20 西安电子科技大学 基于最大边缘度量学习的sar图像变化检测方法
CN109801295A (zh) * 2018-12-26 2019-05-24 西安电子科技大学 基于二次限制区域生长法的图像分割方法
WO2019104767A1 (zh) * 2017-11-28 2019-06-06 河海大学常州校区 基于深度卷积神经网络与视觉显著性的织物缺陷检测方法
WO2019110561A1 (en) * 2017-12-06 2019-06-13 Ventana Medical Systems, Inc. Method of storing and retrieving digital pathology analysis results

Family Cites Families (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB0226787D0 (en) * 2002-11-18 2002-12-24 Qinetiq Ltd Measurement of mitotic activity
JP2006266858A (ja) 2005-03-24 2006-10-05 Railway Technical Res Inst 硬質堆積岩地域の災害発生危険度評価方法
US20080285860A1 (en) * 2007-05-07 2008-11-20 The Penn State Research Foundation Studying aesthetics in photographic images using a computational approach
US8472718B2 (en) * 2011-04-27 2013-06-25 Sony Corporation Superpixel segmentation methods and systems
CN102831598B (zh) 2012-07-04 2014-10-01 西安电子科技大学 多分辨率NMF和Treelet融合的遥感图像变化检测方法
US9424337B2 (en) * 2013-07-09 2016-08-23 Sas Institute Inc. Number of clusters estimation
CN103488971B (zh) 2013-09-06 2017-02-08 电子科技大学 生物礁储层的几何形态识别方法
US9412174B2 (en) * 2013-10-16 2016-08-09 3M Innovative Properties Company Note recognition for overlapping physical notes
US9202178B2 (en) * 2014-03-11 2015-12-01 Sas Institute Inc. Computerized cluster analysis framework for decorrelated cluster identification in datasets
US20150302252A1 (en) * 2014-04-16 2015-10-22 Lucas A. Herrera Authentication method using multi-factor eye gaze
JP2017142196A (ja) 2016-02-12 2017-08-17 公益財団法人鉄道総合技術研究所 表層地盤増幅特性を考慮した、ｐ波初動立ち上がりによる震央距離推定方法
US11166670B1 (en) * 2019-04-30 2021-11-09 Verily Life Sciences Llc Managing meal excursions in blood glucose data

Patent Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102053270A (zh) * 2009-10-30 2011-05-11 中国石油化工股份有限公司 一种基于沉积地层单元的地震相分析方法
CN106056155A (zh) * 2016-05-30 2016-10-26 西安电子科技大学 基于边界信息融合的超像素分割方法
CN106545337A (zh) * 2016-11-25 2017-03-29 西南石油大学 一种基于支持向量机的沉积微相测井识别方法
CN107229917A (zh) * 2017-05-31 2017-10-03 北京师范大学 一种基于迭代聚类的多幅遥感影像共性显著目标检测方法
CN107545577A (zh) * 2017-08-23 2018-01-05 电子科技大学 基于神经网络的沉积相图像分割方法
WO2019104767A1 (zh) * 2017-11-28 2019-06-06 河海大学常州校区 基于深度卷积神经网络与视觉显著性的织物缺陷检测方法
WO2019110561A1 (en) * 2017-12-06 2019-06-13 Ventana Medical Systems, Inc. Method of storing and retrieving digital pathology analysis results
CN108846381A (zh) * 2018-07-06 2018-11-20 西安电子科技大学 基于最大边缘度量学习的sar图像变化检测方法
CN109801295A (zh) * 2018-12-26 2019-05-24 西安电子科技大学 基于二次限制区域生长法的图像分割方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
结合卷积神经网络和超像素聚类的细胞图像分割方法;杨金鑫;杨辉华;李灵巧;潘细朋;刘振丙;周洁茜;;计算机应用研究;20170614(第05期);全文 *

Also Published As

Publication number Publication date
US11327194B2 (en) 2022-05-10
CN112347823A (zh) 2021-02-09
US20210041588A1 (en) 2021-02-11

Similar Documents

Publication Publication Date Title
CN112347823B (zh) 沉积相边界识别方法及装置
Ibrahim et al. Image segmentation methods based on superpixel techniques: A survey
US9042648B2 (en) Salient object segmentation
CN109522908B (zh) 基于区域标签融合的图像显著性检测方法
Huang et al. Morphological building/shadow index for building extraction from high-resolution imagery over urban areas
CN109919944B (zh) 一种复杂场景建筑物变化检测的联合超像素图割优化方法
Derivaux et al. Supervised image segmentation using watershed transform, fuzzy classification and evolutionary computation
Yi et al. A scale-synthesis method for high spatial resolution remote sensing image segmentation
CN110119728A (zh) 基于多尺度融合语义分割网络的遥感图像云检测方法
CN108564085B (zh) 一种自动读取指针式仪表读数的方法
CN105096347B (zh) 图像处理装置和方法
CN105930815A (zh) 一种水下生物检测方法和系统
Hormese et al. Automated road extraction from high resolution satellite images
CN106846322B (zh) 基于曲线波滤波器和卷积结构学习的sar图像分割方法
CN106909902A (zh) 一种基于改进的层次化显著模型的遥感目标检测方法
CN106611423B (zh) 基于脊波滤波器和反卷积结构模型的sar图像分割方法
Schoening et al. Fully automated image segmentation for benthic resource assessment of poly-metallic nodules
CN107464249A (zh) 一种羊只无接触体尺测量方法
Aguilar et al. Optimizing multiresolution segmentation for extracting plastic greenhouses from WorldView-3 imagery
CN107067037B (zh) 一种使用llc准则定位图像前景的方法
CN111259808A (zh) 一种基于改进ssd算法的交通标识的检测识别方法
CN115690086A (zh) 一种基于对象的高分辨率遥感影像变化检测方法及系统
CN112037230B (zh) 一种基于超像素和超度量轮廓图的林区图像分割方法
CN107292268A (zh) 快速脊波反卷积结构学习模型的sar图像语义分割方法
Jenicka et al. A textural approach for land cover classification of remotely sensed image

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