[发明专利]一种基于联合字典的孔隙度反演方法

权利要求书

1.一种基于联合字典的孔隙度反演方法，其特征在于，包括以下步骤：

S1、数据预处理：将测井数据进行分块处理，去除地震数据中的噪声信号；

S2、将分块处理后的测井数据进行字典学习，得到阻抗字典和孔隙度字典；测井数据进行分块处理的具体实现方法为：将测井曲线分成N个相同长度的数据块，每个数据块的长度为L，得到训练样本集Y_Ip表示阻抗数据块集合，作为阻抗字典D_Ip学习的训练样本集；Y_por表示孔隙度数据块集合，作为孔隙度字典D_por学习的训练样本集；

具体实现方法为：将字典的学习过程描述为下述问题的优化过程：

其中，矩阵为期望学习得到的稀疏字典；矩阵A＝[α₁ α₂ … α_N]是稀疏表示矩阵，向量α_i为阻抗字典D_Ip和孔隙度字典D_por所共有的稀疏系数；λ是正则化参数；表示的是向量的F范数，即Frobenius范数，该式表示矩阵元素绝对值的平方和；表示的是向量的0范数，K表示的是预设的稀疏度最大值；

具体字典学习的过程为：

S21、设置正则化参数λ、最大迭代次数T，初始化字典D和稀疏系数矩阵A，当前迭代次数置为t＝1；

S22、进行稀疏编码：假设字典D固定，使用匹配正交追踪算法对训练样本集Y进行稀疏编码，得到稀疏系数矩阵A；

S23、进行字典更新：采用逐个更新策略对字典D中的原子进行更新：k＝1,2,3,…,K，对于第k个原子，进行如下步骤：

S231、为当前要更新的原子，记I_k＝{i|α_i(k)≠0,1≤i≤N}，α_i(k)为α_i中的第k个元素，则I_k表示所有样本中使用到第k列原子d_k的索引集；

S232、用表示去除字典中第k个原子后的字典，表示去除第k行系数以后的系数矩阵，然后计算新的字典与系数矩阵与原信号之间的残差矩阵E^k表示去除第k个原子后对所有样本造成的误差；

S233、根据I_k中的索引值选取E_k中相应的列向量，构成新的误差矩阵并对进行奇异值分解，得到

S234、选取矩阵U中的第一列，即最大特征值对应的特征向量，作为更新后的原子，并更新相应的系数；

S24、令t＝t+1，如果t≤T，则返回步骤S22；否则停止迭代，输出字典

S3、将地震数据进行基于字典的孔隙度反演，得到地震数据的阻抗和孔隙度；具体实现方法为：

S31、提取地震子波ω(t)；

S32、进行井震标定；

S33、建立初始模型：反演方法是建立在地震褶积模型的基础上的，反演模型为：

d(t)＝ω(t)*r(t) (2)

其中，d(t)是地震记录，ω(t)是地震子波，*是褶积运算操作符，r(t)是地层反射系数序列；

当地震波沿着地层的垂直方向入射时，地层反射系数与地层波阻抗具有如下关系：

令地震子波序列为ω(t)＝[ω₁ ω₂ … ω_m]^T，反射系数序列为r(t)＝[r₁ r₂ … r_n-1]^T，地震记录为d(t)＝[d₁ d₂ … d_n-1]^T，其中m为子波序列长度，将式(2)改写成如下的矩阵形式：

简写为：

d＝Wr (5)

其中W为子波褶积矩阵，由下式进行定义：

对于叠后波阻抗反演，根据Russell波阻抗近似公式，当波阻抗变化量相对比波阻抗的绝对数值来说非常小时，将式(3)改写为：

根据式(7)，将反射系数序列r表示为如下的矩阵形式：

将式(8)带入式(4)，得到：

令式(9)表示为：

d＝WDm (10)

令G＝WD，式(10)进一步简写成：

d＝Gm (11)

其中d为(n-1)×1维向量，代表单道地震叠后记录；G为(n-1)×n维矩阵，代表正演模型矩阵；m为n×1维向量，代表单道波阻抗系数的对数；

S34、将作为正则化约束项加入到反演目标函数中，得到最终目标求解函数：

其中，m_Ip表示阻抗，表示通过目标函数得到的阻抗最优解，R_i是一个由0和1形成的对阻抗进行取块的矩阵，表示矩阵元素绝对值的平方和；

S35、通过L曲线法选取正则化参数λ和μ；并设置稀疏度K；

S36、利用最优化相关算法最小化目标函数，求解公式(12)。