最小二乘偏移算法核心是根据反偏移模拟数据与观测数据的匹配程度来判定成像的准确性，并根据残差对成像结果进行修正.用于理论模型时，速度模型、数据正演模拟算法是已知的，并且震源子波都是已知的，这些都提高了数据匹配程度，但是在实际资料情况下是做不到的.对于实际资料来说，除了成像不精确引起的反偏移模拟数据与实际观测数据之间的误差之外，还有很多其他因素的影响，需要尽量消除，为基于数据残差进行成像更新提供相对纯净的环境.

为消除由于子波及其他因素带来的观测数据与反偏移模拟数据之间的振幅差异，对反偏移模拟数据进行振幅校正，使其振幅与观测数据处于同一量级，再求取数据残差.常规方法是采用归一化残差计算方法，即

这种方法是各自将数据调整到统一的数量级上，但还是不满足两套实用数据匹配要求.本文研究并采用能量一致性数据残差求取算法，即

采用最小平方匹配滤波方法来求取匹配因子，使得反偏移模拟数据与观测数据振幅能量数量级一致.通过莱文森(Levinson)快速递推算法求解式(9)所示方程求得匹配因子为(李振春和张军华，2004)

这一步骤也是本文所述方法的创新点所在.相较于归一化数据残差算法，本文能量一致性数据残差可以避免由于残差振幅奇异值对反演成像的影响，也更有利于实现数据匹配，求取残差并进行下一步更新.

目标泛函的梯度如式(12)，

具体实现时，目标泛函的梯度即基于数据残差的逆时偏移，即

① 对每一炮，计算检波点位置的波场残差；将波场残差反传，并与震源激发的正传波场内积，得到单炮的梯度；

② 叠加所有炮求取的梯度值，得到与模型空间的全局梯度.

共轭梯度法是最优化里边最常用的方法之一，它是求解无约束最优化问题中的一种非常重要的方法.本文利用共轭梯度法进行迭代时，模型更新的表达式可以写为

从上述式子当中可以看出共轭梯度法的表达形式比较简单，而且易于实现.仅需存储目标函数的两次梯度值，因此所需内存相对较小.为了使共轭梯度法全局收敛，需要求解出最优化步长α为

迭代最小二乘偏移实现流程如图 1所示.

图 1

Fig.1

图 1 最小二乘偏移流程图 Fig. 1 Flow chart of least square migration for real data

① 面向最小二乘偏移技术的观测数据预处理：最小二乘偏移是基于Born近似波场为基础，因此，直达波、多次波、折射波、面波等反偏移算法无法重构模拟的波场以及噪声都要消除掉.此预处理环节我们是结合商业软件相关流程来完成的.

② 基于观测数据进行初始逆时偏移成像.

③ 基于成像结果进行反偏移数据重构.

④ 能量一致性求取观测数据与重构模拟数据的残差.

⑤ 计算梯度和步长，并求取成像更新量，对成像结果进行更新.

⑥ 循环流程③、④、⑤.

⑦ 判断是否满足终止迭代条件，若满足退出保存成像结果，不满足回到第③步继续迭代.

用模型数据测试方法技术的可行性和有效性.图 2是模拟溶洞地质现象设计的模型，其中包含半径为10 m、20 m、30 m、40 m以及50 m不同尺度散射体，共有四个主要层，自上而下速度分别为2000 m/s、2200 m/s、2500 m/s以及2800 m/s，散射体速度为2000 m/s.

图 2

Fig.2

图 2 原始速度模型 Fig. 2 The velocity model

图 3a是最小二乘偏移初次RTM成像结果，图 3b是最小二乘偏移反演成像结果.图 4是迭代最小二乘偏移数据残差收敛图，说明此方法技术在模型测试过程中数据残差稳步收敛，残差逐渐变小，成像精度越来越高.从图 3展示的两种成像方法效果对比来看，最小二乘反演成像结果剖面无论是横向还是纵向的振幅均衡性提高，中深层构造成像更加清晰，整体照明更加均匀，而且具有更高的成像分辨率，散射体边界更加清晰.图 5是频谱对比图，说明了最小二乘偏移成像提高了主频，拓宽了频带，进一步呼应了提高了成像分辨率的应用效果.图 6a和图 6b分别是逆时偏移成像结果与最小二乘逆时偏移结果与真实反射系数单道振幅对比图，从中可以看出，最小二乘逆时偏移与真实反射系数振幅更加接近，尤其是中深层部分，说明最小二乘逆时偏移技术可以更好地对由于几何扩散或吸收衰减引起的中深层振幅能量衰减进行振幅和照明补偿，从而具有更高的保幅性.

图 3

Fig.3

图 3 RTM成像结果(a)和最小二乘偏移成像结果(b) Fig. 3 The comparison of (a)RTM imaging; (b)Least-squares migration imaging

图 4

Fig.4

图 4 最小二乘偏移数据残差收敛图 Fig. 4 Data residual curve versus the iteration numbers

图 5

Fig.5

图 5 RTM与LSRTM频谱对比 Fig. 5 The comparison of frequency spectrum of RTM and LSRTM imaging

图 6

Fig.6

图 6 RTM(a)与LSRTM(b)与真实反射系数单道振幅对比 Fig. 6 The comparison of reflection coefficient of RTM and LSRTM

为了更好地说明此方法技术对实际资料的适用性和有效性，结合某区二维实际资料进行试处理.实际地震数据为我国某区陆上资料.资料中面波、折射波、不规则强干扰发育.根据上文所述技术流程，需先对其进行预处理.在剔除异常振幅和直达波之后，选用保真度高的去噪方法，采用噪声自动识别与衰减、炮集自动统计道编辑、时空域相干噪声衰减等技术，压制原始资料中存在的上述干扰波，确保各种强能量干扰得到很好的压制，同时不损失有效波.

在迭代最小二乘偏移实现过程中，测试了不同方法求取数据残差的效果.图 7为第6次迭代时归一化方法和能量一致性方法求取残差的效果对比.其中，本文技术流程中所用的能量一致性残差可以更有效进行数据匹配，求取数据残差进行下一步成像更新.

图 7

Fig.7

图 7 不同方法求取残差对比 (a)归一化方法求取残差；(b)能量一致性求取残差. Fig. 7 The comparison of (a) Normalized error residual; (b) Energy consistency error residual

图 8为随着迭代次数增加某单炮数据残差示意图，图 9为残差收敛曲线.从中可以看到，随着迭代次数增加，数据残差越来越弱直至趋于稳定，说明随着逐次迭代，成像过程中逐步加入了更丰富的有效信息，实现有效信息利用最大化，从而有利于在常规偏移中不易利用的一些弱有效信号的利用，提高成像分辨率和保幅性.

图 8

Fig.8

图 8 最小二乘偏移不同迭代次数数据残差示意图 (a)第2次迭代；(b)第4次迭代；(c)第5次迭代；(d)第6次迭代；(e)第8次迭代. Fig. 8 The error residual versus the iteration numbers (a)The second iteration；(b)The fourth; (c)The fifth; (d)The sixth; (e)The eighth iteration.

图 9

Fig.9

图 9 最小二乘偏移不同迭代次数数据残差曲线 Fig. 9 Data residual curve versus the iteration numbers

图 10a和图 10b分别是逆时偏移成像与最小二乘逆时偏移成像结果.从两者效果对比看出，最小二乘偏移成像结果浅层成像更加清晰，分辨率更高，中深层部分如左下角椭圆框内小断块成像更加干脆清晰.图 11是逆时偏移成像和最小二乘偏移成像频谱分析，可以看出，最小二乘偏移成像提高了主频，拓宽了有效频带，提高了分辨率，展示了最小二乘偏移在小尺度构造如小断块等构造刻画方面的应用潜力.

图 10

Fig.10

图 10 RTM成像结果(a)和最小二乘逆时偏移成像结果(b) Fig. 10 The comparison of (a)RTM imaging; (b)Least-squares migration imaging

图 11

Fig.11

图 11 成像频谱分析 Fig. 11 The comparison of frequency spectrum of RTM and LSRTM imaging

4.1      研究和讨论了最小二乘逆时偏移成像方法原理及具体实现关键技术，包括有限差分解微分方程进行反偏移数据重构，基于数据残差逆时偏移求取梯度及共轭梯度法反演优化算法，实现了最小二乘逆时偏移技术.并研究了在实际处理时还需要注意的关键技术细节，包括面向最小二乘偏移技术的观测数据预处理及能量一致性数据残差求取.模型测试和实际资料试处理说明了方法技术流程的有效性和适应性.数值测试结果表明，与常规逆时偏移方法相比，该方法技术流程可以有效应用于实际资料处理，提高成像的分辨率与保幅性，从而更有利于后续岩性储层成像.

4.2     最小二乘偏移是基于反演理论的方法，算法核心是根据反偏移模拟数据与观测数据的匹配程度来判定成像的准确性，并根据残差对成像结果进行修正.在用于实际资料处理时，由于地下介质、传播过程及采集条件的复杂性，要真正实现完全模拟地震数据是做不到的，有多方面因素影响.下面是笔者分析的把最小二乘方法技术推向实用化接下来还需要改进优化的技术点：

(1)数据模拟算法因素：例如我们目前通常所用的是基于声波方程进行波场模拟，而实际地下介质可能是黏声、各向异性等特性，从而导致无法模拟出观测记录包含的全部波场，因此实际处理时要根据处理探区的实际地下介质情况，建立相应的数据模拟算法，才能模拟出观测记录中包含的全部有效波场；

(2)噪声因素：野外采集条件复杂多样，会产生很多干扰和噪声，这些噪声无法通过反偏移算法模拟出来，就会影响数据残差的准确求取，影响方法技术的收敛和稳定性，需要从观测数据中消除；

(3)子波因素：子波的估计影响反偏移模拟数据的振幅、相位，从而影响数据匹配精度.而实际震源子波尤其是在陆上是空变的，每炮都可能各不相同，强吸收衰减时，震源子波时空变化剧烈，难以准确估计，从而导致地震仪器记录的波场与反偏移模拟的波场的振幅相位不匹配.因此，子波估计会是最小二乘偏移的优化比较挑战性的一环.

(4)速度场因素：在实际资料处理时，地震速度误差是不可避免的.速度场的误差也会引起反偏移模拟记录的误差，因此，背景速度估计要尽量精确.

(5)反偏移模拟方法：目前常用的反偏移都是基于叠加成像剖面进行反偏移数据模拟，这就会导致中远偏移距模拟数据存在比较严重的由于大的速度误差引起的数据时差问题，研发基于道集的反偏移算法，是完善最小二乘偏移方法所需要攻关的方向之一；

(6)效率因素：当遇到大规模的尤其是三维实际资料时，计算量就是很大难题.引入GPU加速可以大幅地提高最小二乘偏移的计算效率，从而尽量满足实用需求.