基于形态完备性的非均质岩心三维结构重建算法的制作方法

本发明涉及利用计算机数学建模方法从一张二维(2d)训练图像随机模拟非均质三维(3d)岩心结构技术，属于图像处理技术领域。

背景技术：

由于在岩心形成过程中的构造作用、沉积作用和成岩作用等综合因素使得岩心结构存在非均质性。这种非均质性决定着储层质量的好坏，影响油田开发的效果。在不同储层岩心中，这种非均质表现程度不同。如常规油气开采的砂岩岩心中，这种非均质性较小；而在非常规油气开采的页岩岩心中，这种非均质特性会非常明显。不同的学者对储层非均质性提出不同的划分方案，一般地将储层非均质性分为宏观非均质性(层间非均质性、平面非均质性和层内非均质性)和微观非均质性(孔隙非均质性、颗粒非均质性和填隙物非均质性)。在数字岩心三维重建中，非均质性是指岩心微观孔隙非均质性。

可以通过物理成像和数学建模两种方式获取岩心三维结构。本发明是针对借助计算机数学建模方法，从单幅二维训练图像重建出非均质岩心三维结构。

目前，在计算机数学建模过程中通常利用训练图像同时控制三维结构三个正交切面，这样操作的目的使重建结构三个正交切面同时具有与训练图像相似的形态。一般来讲当训练图像具有较好的均质性时，采用这种重建方式能够得到较为理想的三维重建结构。另一方面，非均质岩心图像往往会具有孔隙多尺寸、分布不均匀等特征，使得岩心三维结构各个正交切面形态差异较大。即说明对于非均质岩心图像利用常规的重建方式无法得到理想的三维重建结果。

对于非均质岩心的三维重建研究,目前学术界处于探索阶段。tahmasebi在2012年提出的ccsim算法，虽然相关文献描述可以重建具有非平稳特征的岩心结构，但是却没能够给出重建出这种非平稳特征三维结构可信服的具体理论，因此该算法实际上是无法稳定而有效重建非平稳三维结构。tahmasebi在2016年继而提出利用改进ccsim算法从单幅图像重建具有多尺寸特征页岩三维结构。在考虑到满足气体渗流特性条件下的页岩分辨率大小，实际上能够重建的图像仅仅是岩心一个尺寸很小的视域。而该尺寸大小下的视域实际上是无法有效代表岩心的渗流特性。同时该尺寸大小下的二维视域图像本身无法有效代表三维结构的整体特征，即该视域图像并不适合用于三维重建。总之非均质岩心的三维重建仍处于初始发展阶段，存在众多需要继续深入研究的内容。

本发明受国家自然科学基金项目“岩石微观非均质结构三维图像重建及分辨率提升技术研究(61372174)”资助，以解决上述现有技术中所存在的问题，为非均质岩心三维重建提供一种新的算法。

技术实现要素：

本发明从基于单幅训练图像重建三维结构方式出发，首先给出了正交弦的定义，然后通过分析训练图像形态、三维结构各正交切面形态和正交弦之间的关系，从同时考虑三个正交切面和两个正交切面两种方式分别给出了训练图像形态完备性测评方法，最后依据该测评分析方法给非均质岩心三维重建算法，在算法设计中又给出特征区域提取和特征区域修复的具体措施。本发明具体原理内容包括下述四个方面：

1、同时考虑三个正交切面时训练图像形态完备性分析方法

为了使得重建三维结构各正交切面都具有训练图像的统计特征和形态特征，在基于单幅训练图像的三维重建过程中，利用同一训练图像对目标三维结构的三个正交切面同时进行约束，如附图1所示。

为了说明在附图1所示的重建过程中各正交切面之间的满足关系，附图2-6给出了线、弦与正交弦的定义。附图2是一幅二维岩心图像，附图3是从附图2中所提取的一个局部形态。将附图4中d方向上开始与终止于同一相，且所有像素点均属于该相的有向线段称为弦，用c^d表示。将弦的任意部分称为线，用l^d表示，即满足(1)式。附图4中c^y是y方向的一条弦，而l^x是x方向的一条线。将c^d上所包含像素点的个数称之为c^d的弦长，用l(c^d)表示，相应地将l^d上所包含像素点的个数称之为l^d的线长，用l(l^d)表示。

将附图5中相正交于点的弦c^x和c^y所构成的结构c^x⊥c^y称为正交弦，标记为(2)式。相应地将附图6中相正交于点的线l^x和l^y所构成的结构l^x⊥l^y称为正交线，标记为(3)式。

为了显示三维岩心中不同正交切面正交弦长之间的关系，用附图7(a)所示的一个尺寸大小为128×128×128micro-ct真实岩心结构进行说明。附图7(b)为附图7(a)中所标注的局部区域，尺寸大小为31×31×31。附图7(c)为过附图7(b)所示三维结构中点的正交切面图。附图7(d)为附图7(c)中所对应的sxy、szx和syz三个相正交切面。

通过分析附图7发现，正交弦在三维重建中具有的关键作用：第一，正交弦能够有效表述三维重建中训练图像具体形态。重建结构各正交切面中具有与训练图像相似的形态是由训练图像所对应的正交弦所构成。第二，由于每个像素点对应一个正交弦，因此当一个局部形态看作是由一系列像素点构成时，该局部形态也可以看作是由一系列正交弦所组合而成。因此可以通过每个点来测量训练图像的具体形态。

继续深入研究附图7发现，由点在正交切面sxy、szx和syz分别对应投影点和所确定的正交弦满足(4)式：

以为例，其表示由szx切面内点所确定正交线的x方向第i条线的线长等于由sxy切面内点所确定正交线的x方向第i条线的线长，同时要求每个切面x方向的两条线能够组合成一条弦。

根据附图1所示的重建方式和正交弦所代表的物理含义，以sxy切面中点为例，若在szx和syz分别存在满足(4)式的点和说明和可以构成三维空间中的一点并且点在三个正交切面对应的形态均属于训练图像，如附图8所示。简言之，该训练图像所提供用于重建由点所确定的形态是完备的。反之，若点无法在szx和syz找到满足(4)式的点和时，则说明点无法构成点当在sxy切面所对应的形态为点所确定形态(正交弦)时，在szx和syz切面所对应的形态(正交弦)仍同时属于该训练图像。简言之，该训练图像所提供用于重建由点所确定的形态是不完备的。

2、同时考虑两个正交切面时训练图像形态完备性分析方法

在实际中对于某些二维图像，虽然不存在三个正交切面都具有与之相似形态的三维结构，但是存在两个正交切面都具有与之相似形态的三维结构。即说明该二维图像在重建两个正交切面都具有与之相似形态三维结构时，形态是完备的。以三维结构szx和sxy正交切面为例，此时在三维重建中需要满足(5)式所示的完备性。具体含义如附图8所示。

本发明引入完备率参量rc用来表示某相满足完备性要求点数ns与总点数nall之间的比值，即

3、特征区域提取

特征区域是该区域内绝大数形态只能够出现在三维结构一个切面上。主要从下述几个方面来确定该特征区域。

(1)在控制两个正交切面重建方式中形态完备性判断结果。这里从三个判决结果中选用完备性最低的切面作为参考图像，以附图9所示的非均质二维岩心图像为例，此时选择的参考图像如附图10所示。

(2)在控制三个正交切面重建方式中形态完备性判断结果。这里选用完备性最高的那个切面作为参考图像，如附图11所示。

(3)判决区域的尺寸大小ix×iy。判决区域尺寸大小是指搜索不满足完备性要求点的范围大小。很明显此区域大小直接决定于训练图像两相平均颗粒尺寸大小。以岩心二维图像为例，此判决区域直接取决于岩石相和孔隙的平均尺寸，即两相的平均弦长大小。而岩石相和孔隙相平均弦长可以用(7)和(8)表示。而最终判决区域的尺寸ix×iy为两相平均尺寸之和，即可以用(9)式表示。

(4)控制两个正交切面方式中不满足完备性要求的点数与该区域大小的比值r2，该参数表示判决区域只可能出现在三维结构一个切面的概率大小，本发明中设置该参数的门限值为表示当某判决区域不满足完备性要求的比率大于该门限值时，就将该区域认定为只可能出现在三维结构的一个切面上。

(5)控制三个正交切面方式中满足完备性要求的点数与该区域大小的比值r3，该参数表示判决区域可能同时出现在三维结构三个正交切面的概率大小。本发明中设置该参数的门限值为表示当某判决区域满足完备性要求的比率大于该门限值时，就认定该区域被排除在只可能出现在三维结构的一个切面上范围外。最终提取的特征区域如附图12所示。

4、特征区域修复

特征区域修复是指利用原始二维图像中除特征区域以外的区域修复特征区域，如附图13所示。本发明中具体修复方法采用块匹配的算法，该算法主要包括这么几个关键环节。

(1)修复模板尺寸大小。该修复模板决定了每次修复模式和整体修复区域的形态。为了有效均衡相邻修复模块之间的连续性和变化性，本发明采用已知区域的孔隙相平均弦长和岩石相平均弦长的和作为该修复模板尺寸大小。即sx和sy通过式(10)计算。其中(i＝x,y；j＝p,g)表示i方向j相的平均弦长。

(2)匹配区域的尺寸大小与在本发明中采用修复模板尺寸的1/4大小，即满足(11)式。需要说明这里匹配区域并不是由构成的一个矩形区域。由于待填补区域没有固定的形态，使得匹配区域也没有固定形态。实际操作中该区域是指以某个待填补区域左上角第一个点为基准，将修复模板向上回退大小和向左回退大小之后所包括已知区域大小。

(3)最佳匹配区域。最佳匹配区域是指利用修复模板在原始训练图像中进行全图逐点扫描已知区域，当某扫描区域中匹配区域与待修复区域中所对应匹配区域相似度最高时，将该扫描区域称为最佳匹配区域。当获得最佳匹配块后，就可以将最佳匹配区域填补到待修复区域中对应区域中。填补时只需要对待修复的区域进行填补就可以，也就是说已知的匹配区域保持不变。最终修复后的图像为附图14。

附图说明

图1基于二维图像的三维重建方式；

图2二维岩心图像；

图3从图2提取的一个局部形态；

图4二维图像中的弦与线；

图5二维图像中的正交弦；

图6二维图像中的正交线；

图7三维结构中点与正交弦关系；

图8三维重建过程中各正交切面正交弦关系；

图9本发明实施例中二维训练图像；

图10是本发明实施例中控制两个切面形态完备性；

图11是本发明实施例中控制三个切面形态完备性；

图12是本发明实施例中提取的特征区域；

图13是本发明实施例中图像修复原理；

图14是本发明实施例中修复后的图像；

图15是本发明实施例中重建的三维结构；

图16是本发明实施例中参照的三维结构；

图17是本发明实施例中重建结构和参照结构x方向两点相关函数对比结果；

图18是本发明实施例中重建结构和参照结构y方向两点相关函数对比结果；

图19是本发明实施例中重建结构和参照结构x方向线性路径函数对比结果；

图20是本发明实施例中重建结构和参照结构y方向线性路径函数对比结果；

具体实施方式

下面结合一个具体案例对本发明具体实施过程作进一步详细说明，但所述实施案例仅对本发明实现方法作一个详细说明，而不应理解为是对本发明保护内容的任何限制。本发明具体实施过程包括下述几个步骤：

第一步，根据(4)式计算训练图像在控制三个正交切面情况下的完备性数据。

第二步，根据(5)式计算训练图像在控制两个正交切面情况下的完备性数据。

第三步，结合前两步分析结果，确定对应给定训练图像在目标三维结构各个正交切面形态重现情况，从而确定所对应的重建方式。以附图9为例，孔隙相在三个正交切面形态完备性都比较高，可以采用同时控制三个正交切面的方式。但是存在一个特征区域在控制两种重建方式中都非常低，即说明该特征区域只能够出现在三维结构的一个正交切面上。因此先确定该特征区域，然后将该特征区域采用图像修复技术进行修复，再采用同时控制三个正交切面的重建方式以修复后的图像为训练图像，重建三维结构，最后再在所重建三维结构的一个正交切面上形成该特征区域。在此例中是sxy切面。

第四步，提取图像的特征区域，即由考虑两个正交切面完备性结果中选取完备性最低的结果，如附图10，和考虑三个正交切面完备性结果中选取完备性最高的结果，如附图11，确定出具有特征区域的结果图，如附图12。

第五步，特征区域修复，即利用附图13中已知区域对对所标注的特征区域进行填充，在具体填充过程中采用前面所陈述的块匹配方式。所修复后的结果如附图14所示。

第六步，以修复后的图像，附图14，作为三维重建中的训练图像，采用同时控制三个正交切面的重建方式，重建一个三维结构。

第七步，根据前面描述，可以知道训练图像中的特征区域形态只能够出现在三维结构的一个正交切面上，因此再在重建三维结构的一个切面上构成该特征区域形态，最终的重建结构如附图15所示，附图16是所参照的三维结构。

为了显示本发明效果，在该实施方案中利用随机图像重建中常用的两点相关函数和线性路径函数对重建的结果进行量化比较。比较结果如附图17-20所示。由图可见，重建图像与原始图像的统计函数曲线具有较好的吻合。从而说明采用本发明能够有效地实现非均质岩心三维重建工作。

上述实施例只是本发明的优选实施案例，并不是对本发明所述技术方案的限制，只要是不经过创造性劳动即可在上述实施案例的基础上实现的技术方案，均应视为落入本发明内容的保护范围内。