基于三边结构导向滤波的三维地震数据图像降噪方法

摘要

本发明公开了一种基于三边结构导向滤波的三维地震数据图像降噪方法，基于传统结构导向滤波的框架，改善结构张量的构建，并结合保边效果较好的三边滤波，大幅提升了原有的高斯核各向异性扩散滤波的性能，使得去噪效果好，结构信息损失少，在地震图像的保边降噪中发挥了极大功效，能使地震图像在去噪的同时留存更多结构信息，便于地震图像的解释及后续处理。

主权项

一种基于三边结构导向滤波的三维地震数据图像降噪方法，其特征在于，该方法为：1)设置对三维地震图像数据的三边结构导向滤波参数：计算结构张量的滤波参数σ₁，扩散滤波的时间步长timeSetp，三边结构导向滤波的迭代次数iterTimes；2)导入三维地震图像数据dataIn，并计算三维地震图像中每一点处的梯度，得到三个梯度分量体g_x、g_y和g_z：三维地震图像点(m,n,k)处的梯度计算式如下：$\left\{\begin{array}{c}{g}_x(m,n,k)=\frac{dataIn(m+1,n,k)-dataIn(m-1,n,k)}{2};\\ g_y(m,n,k)=\frac{dataIn(m,n+1,k)-dataIn(m,n-1,k)}{2};\\ g_z(m,n,k)=\frac{dataIn(m,n,k+1)-dataIn(m,n,k-1)}{2};\end{array}\right.$其中，dataIn(m,n,k)表示三维地震图像中，点(m,n,k)的像素值；g_x(m,n,k)表示点(m,n,k)处梯度的第一分量体，即三维坐标系中的x方向分量，g_y(m,n,k)表示点(m,n,k)处梯度的第二分量体，即三维坐标系中的y方向分量，g_z(m,n,k)表示(m,n,k)处梯度的第三分量体，即三维坐标系中的z方向分量；dataIn(m+1,n,k)、dataIn(m‑1,n,k)、dataIn(m,n+1,k)、dataIn(m,n‑1,k)、dataIn(m,n,k+1)、dataIn(m,n,k‑1)为点(m,n,k)的邻近点对应的像素值；3)利用上述梯度分量体g_x、g_y和g_z，计算三维地震图像数据dataIn的梯度张量GT：$GT=g·g^T=\left[\begin{array}{ccc}{g}_x^2& g_x{g}_y& g_x{g}_z\\ g_y{g}_x& g_y^2& g_y{g}_z\\ g_z{g}_x& g_z{g}_y& g_z^2\end{array}\right]=\left[\begin{array}{ccc}{\mathrm{GT}}_{11}& {\mathrm{GT}}_{12}& {\mathrm{GT}}_{13}\\ {\mathrm{GT}}_{21}& {\mathrm{GT}}_{22}& {\mathrm{GT}}_{23}\\ {\mathrm{GT}}_{31}& {\mathrm{GT}}_{32}& {\mathrm{GT}}_{33}\end{array}\right],$g＝[g_x,g_y,g_z]^T,其中，GT₁₁,GT₁₂,GT₁₃,GT₂₂,GT₂₃,GT₃₃分别为梯度张量GT的六个分量体；4)对上述计算梯度张量GT的六个分量体GT₁₁,GT₁₂,GT₁₃,GT₂₂,GT₂₃,GT₃₃分别进行三边滤波，得到三维地震图像数据dataIn的结构张量ST；5)根据上述结构张量ST构建三维地震图像的扩散张量DT，并计算三维地震图像每一点处的流量；其中：扩散张量DT的表达式为：$\left\{\begin{array}{c}DT=ϵ·(v_2·v_2^T+v_3·v_3^T)\\ ϵ=\frac{Tr\left(J_{0}J\right)}{Tr\left(J_0\right)Tr\left(J\right)}\\ [v_1,v_2,v_3]=eigV\left(ST\right)\end{array}\right.$三维地震图像点(m,n,t)处的流量flux(m,n,t)的三个分量体flux_x，flux_y，flux_z表达式为：$flux=ϵ·DT·\left[\begin{array}{c}{g}_x\\ g_y\\ g_z\end{array}\right]$其中，ε取值范围为0～1；Tr(·)表示求矩阵的迹；J₀为初始三维地震图像结构张量，即ST₀；J为当前迭代的三维地震图像的结构张量，即ST；eigV(·)表示特征向量求解，且特征向量v₁,v₂,v₃对应的特征值满足λ₁>＝λ₂>＝λ₃；flux＝[flux_x，flux_y，flux_z]^T；6)计算三维地震图像点处流量的散度，获得每一点处在扩散过程中单位时间上的增加量I_Δ；其中所述散度即流量各分量在对应方向上的差分的和；7)根据步骤1)中设定的时间步长timeSetp，计算出当前迭代的三边结构导向滤波的结果，即计算出输出的三维地震图像数据dataOut：dataOut＝dataIn+timeSetp·I_Δ；8)判断是否达到步骤1)中设定的三边结构导向滤波的迭代次数iterTimes，若达到，则输出；否则，将dataOut作为待导入的三维地震图像数据，返回步骤2)。