快速回转波近地表速度建模方法
脉冲多普勒高度计地面回波的建模与仿真
ｌ７ ３
本 研 究 的地 面 是 一 个山 日标 。对 于 运 动 面 目标 ， 般 进 一 行 离 散 化 处 理 ， 其 看 成 是 很 多 运 动 点 ｝标 的集 合 。 仿 真 时 将 ＝ ｊ 采 用 网格 法 ［ 、 ６ 将 波 束 照 射 的 区 域 在 时 域 上 以 距 离 环 ， １５、］ 频 域 上 以等 多 谱 勒 带来 划 分 ， 图 ２所 示 。 单 元 格 所 在 的距 如
合理性 。 关键词 ： 达高度计 ； 回波 ； 雷 地 网格 法 ； 普 勒 频 率 ； 多 建模 与仿 真
中 图分 类 号 ：Ｎ ５ Ｔ ９８ 文 献 标 识 码 ： Ａ 文章 编号 ：０６—００ （０ １ ０ — １６— ３ １０ ７ ７ ２ １ ） ３ ０ ３ ０
雷达高度计 是测 定飞行 器与地 面的高 度设 备。对雷达 高度计挂机测试 能 真实 接收 地 面回波 ， 也存 在许 多 的弊 但 端： 飞机不能提供高过载 、 不能很好拟合 火箭 、 弹等的飞行 导
抗 干扰能力 。文献 ［ ］ ３ 中详尽 的描述了调制脉冲的特 点并给
１
是一个用来做转换的中介 坐标 系 ， 它与 ０ ｙ 坐标 的区别 就 一ｘｚ
是 ｏ原 点 是 弹 头 ， 即 ＝ｚ 。０ 一ｘｙｚ 标 系 如 图 所 示 —ｈ ｖｖ ｖ坐
出了多普勒频偏的表达式 ： 。
４ ）雷 达方 向 图模 型 雷 达 天 线 的 方 向 图 Ｇ是 方 向 角 的 函 数 ， Ｇ（ 。方 向 即 ）
功率
Ｐｒ ： —ｔ ｔ Ｐ Ｇ Ｇ ｒ２ ｈ
—
—
图的曲线一股 由一 个主峰 （ 为主瓣 ） 称 和许 多其他 峰 （ 为 称
副瓣 ） 成 。 考 虑 到 引 信 主 要 工 作 目的 是 测 高 ， 采 用 主 波 组 故
使用回折波走时数据的起伏地表速度建模方法
! $%$"年#月
中国海洋大学学报
&'()*+),-.*/*,'-0 10)2'(3)45 */,6)0-
!"!#"#%?=!%7% 389:;$$%$"
使用回折波走时数据的起伏地表速度建模方法"
何万利杨华臣张建中""
!中国海洋大学海洋地球科学学院$海底科学与探测技术教育部重点实验室$山东 青岛 $??"%%"
! "!基金项目#山东省自然科学基金项目!J($%"#]+%%""资助 3N99QT:8FHO:E80K:NTKV3AD8BA8/QNBFK:DQBQS3EKBFQBL&TQXDBA8!J($%"#]+%%""
收稿日期#$%$%<""<$7%修订日期#$%$"<%"<$$ 作者简介#何万利!"##!<"$男$硕士生''<UKDV#6MV#!="""?>;AQU !""!通讯作者#'<UKDV#ZEKBL@Z"QNA;8FN;AB
!!近地表 低 降 速 带 速 度 和 厚 度 的 横 向 变 化$会 使 来 自其下部的 地 震 波 产 生 时 移$导 致 地 震 资 料 成 像 畸 变 及解释错误'静校正是消除低降速带对地震成像影响 的常规技术$其 效 果 极 大 地 依 赖 于 近 地 表 速 度 模 型 的 精度)"<>*'在目前的地震资料处理中$主要用走时层析 方法建立近地表速度模型)=*'走时层析通过反复迭代 更新初始速 度 模 型$使 模 型 计 算 走 时 与 实 际 观 测 走 时 的残差达到最小$从而建立近地表速度模型)!<Y*'但这 种方法的结 果 严 重 依 赖 于 初 始 速 度 模 型$且 需 要 进 行 射线追踪运 算$计 算 效 率 低' 而 基 于 速 度 随 深 度 增 加 假设的回折波走时反演方法不需要初始速度模型和射 线追踪$因而具有较高的效率$得到了广泛的研究和应 用)#<"$*'
一种三维地面场景SAR回波仿真的快速实现方法
一种三维地面场景SAR回波仿真的快速实现方法景国彬;张云骥;孙光才;邢孟道;保铮【摘要】Due to huge calculation of SAR echo simulation of a three-dimensional (3D) ground scene,a fast simulation method isproposed.First,the data of a 3D ground scene are fractionally interpolated based on the Fractional Brownian Motion (FBM) model.Second,the data interpolated are divided into lots of small facets whose precision meets the simulation requirement,the backscattering coefficients of the small facets are calculated.Third,the shaded area is judged by the comparison method of overlook angles.Three parallel levels are analyzed and the kernel function is designed.Finally,the computing method of reduction adding and external thread are also utilized under the framework of GPU in order to get high efficiency.With the methods mentioned above,the imaging result of experimental DEM data verifies the validity and superiority of the proposed method.%针对三维地面场景合成孔径雷达回波仿真中存在计算量巨大的问题,提出了一种基于图形处理单元的合成孔径雷达回波快速仿真方法.首先,采用分形布朗运动模型对数字高程模型数据进行分形插值处理;再对插值后的数据进行小面元剖分并计算了剖分后小面元的后向散射系数;接着利用改进的下视角比较法对三维场景的阴影遮挡进行快速判断;最后详细分析了回波仿真过程中的3个并行层次,设计了核函数,并利用线程外推和归约相加的计算方法,实现了图形处理单元编程架构下三维场景合成孔径雷达回波的快速仿真.利用这种方法对实测数字高程模型数据进行了回波仿真和成像处理,验证了该方法的正确性和高效性.【期刊名称】《西安电子科技大学学报（自然科学版）》【年(卷),期】2017(044)003【总页数】7页(P1-7)【关键词】三维地面;合成孔径雷达;分形布朗运动模型;阴影遮挡;图形处理单元【作者】景国彬;张云骥;孙光才;邢孟道;保铮【作者单位】西安电子科技大学雷达信号处理国家重点实验室,陕西西安710071;西安电子科技大学雷达信号处理国家重点实验室,陕西西安710071;西安电子科技大学雷达信号处理国家重点实验室,陕西西安710071;西安电子科技大学雷达信号处理国家重点实验室,陕西西安710071;西安电子科技大学雷达信号处理国家重点实验室,陕西西安710071【正文语种】中文【中图分类】TN957在对新型合成孔径雷达(Synthetic Aperture Radar,SAR)设计方案和成像算法验证和评估时,需要满足特定参数的合成孔径雷达回波信号,而这无法直接通过雷达录取回波得到,因此需要进行合成孔径雷达回波仿真工作.合成孔径雷达回波仿真分为逆向法和正向法,其中,逆向法是基于二维合成孔径雷达图像进行回波录取,其简化了后向散射系数的计算,但丢失了场景的高程信息,最终获取的信息量只是基准合成孔径雷达图像的子集,成像效果有限;正向法是基于原始地形结构的三维数据进行回波录取,可以实现场景的重构,在复杂地形探测、景象匹配制导中有重大作用[1].相比逆向法,正向法虽然复杂,但却拥有逆向法无法实现的优势[2],因此研究三维地面场景的正向法回波快速仿真是非常必要的.针对正向法仿三维场景合成孔径雷达回波,文献[3]提出了一种三维雷达后向散射模型;文献[4]研究了理想轨迹下的具有地面起伏特征的三维场景合成孔径雷达回波仿真,对计算量巨大等问题未进行更进一步的优化处理;文献[5]研究了三维森林场景合成孔径雷达回波数据的模拟方法,提出了等效散射模拟方法,但对其后期计算需进行优化处理,不能运用到大场景的快速处理.正向法仿三维场景合成孔径雷达回波的研究当前还存在以下技术问题:数字高程模型(Digital Elevation Model,DEM)数据如何精细化处理,如何计算三维场景的后向散射系数,阴影区域范围的有效判断以及三维地面大场景回波仿真如何快速实现.针对以上问题,笔者给出了一种基于图形处理单元(Graphics Processing Unit,GPU)架构的合成孔径雷达回波快速仿真方法,详细讨论了基于分形布朗运动模型(Fractional Brownian Motion,FBM)的插值方法、小面元后向散射系数计算和三维遮挡判断等内容,并利用图形处理单元对回波仿真算法进行了并行优化实现.1.1 数字高程模型数据精细化处理三维地面场景的数字高程模型数据的地平面分辨单元虽然较小,但仍然可以认为在一个地面分辨单元中包含大量的散射体.为了获得高质量的合成孔径雷达回波,首先对数字高程模型数据进行精细化处理.此外,根据三维地面场景在很大范围尺度上具有统计自相似性(即随机分形特性),满足分形布朗运动模型的统计特性[6],笔者将基于分形布朗运动模型对原始数字高程模型数据进行分形插值,得到精细化的数字高程模型数据,具体原理如下.设地理坐标系中二维地面上的任意一点坐标为(x,y),z(x,y)为该点所对应的数字高程,则根据分形布朗模型[6]得其中,P(t)是随机变量t的累积概率密度函数,服从N(0,σ2)的高斯分布,平坦区域的值σ较小,陡峭区域的σ值较大;参数H为自相似参数,表示地形的自相似程度,可以用来描述地面的粗糙程度,满足0＜H＜1,平滑区域的值较大,粗糙区域的值较小;,表示样本的采样.设D为分形维数,则在二维平面空间中D和参数H的关系为根据分形布朗的随身统计特性[4],可得出式(1)等价为其中,E(·)表示数学期望,常数C表示随机变量的均值.对上式两边同时取对数,可得利用线性回归模型,求出直线的斜率k和Y轴的截距ln C,继而得到分形参数H=k,σ=(2π)1/2C/2.在确定了分形参数后继续进行分形插值,其处理过程为:首先创建一个两维均为原始数字高程模型数据维数两倍的零矩阵,将原始数字高程模型数据放置在新矩阵的行列坐标均为偶数的位置上,接着通过内插公式(5)得到新矩阵中行列坐标均为奇数的位置上的数据;之后通过内插公式(6)得到新矩阵中行列坐标有一个为奇数的位置上的数据;重复以上插值操作,直至得到符合合成孔径雷达回波仿真要求的精细化的数字高程模型数据为止.当新矩阵中行列坐标均为奇数时,其高程值对应的插值公式为当新矩阵中行列坐标只有一个为奇数时,其高程值对应的插值公式为其中,G代表服从标准正态分布的随机变量,代表插值后的矩阵的样本间隔.通过分形插值处理得到精细化数字高程模型数据,既保持了地面的随机起伏特性,又满足了回波仿真的精度要求.1.2 小面元后向散射系数的计算获得精细采样间隔后的三维场景数字高程模型数据,接着需要计算其被电磁波照射后的后向散射系数.通常三维场景可以看成是由大量平面小面元组成的[7],三维场景的电磁散射特性正是这些小面元的后向散射场相干叠加的结果.常规的小面元剖分是三角形面元剖分法,其计算量较大.因此在保证仿真精度的情况下,为了降低计算量,笔者提出了一种四边形面元剖分法,并对剖分后的四边形小面元进行后向散射系数的计算,其中,四边形小面元可看成两个三角形小面元的组合,其剖分原理与三角形小面元剖分类似[7].此外,四边形小面元拥有偶数边的性质有助于图形处理单元并行运算,实现了后向散射系数的快速计算.如图1(a)所示,坐标系O-iXYZ是三维空间中任意一个小面元α的局部坐标系,Oi为小面元中心,X轴、Y轴、Z轴的方向与地理坐标系相同,S点为雷达载机位置.设四边形小面元由M×N个散射点组成,并设在(i,j)处的散射点的高度为Hi,j,则小面元所在的平面方程为其中,p表示小面元在X轴方向的斜率,q表示小面元在Y轴方向的斜率,z表示小面元中心的高度,Δx是小面元沿X轴方向的采样间隔,Δy是小面元沿Y轴方向的采样间隔.用最小二乘法求解参数p、q、z,拟合小面元的最佳方程,即其中,hi,j表示(i,j)处散射点的真实高度.求出参数p、q、z后即可确定小面元的平面方程,则小面元所在平面的单位法向矢量n为由雷达相位中心发出射向此小面元中心的电磁波的入射矢量ki为则雷达发射电磁波在小面元上的局部入射角θi为设小面元在地面上的投影面积为Sd,则小面元的实际面积Si为假定电磁波垂直照射单位面积的该植被地面的后向散射系数为σi(0),则该小面元的后向散射系数为因此,可得任意一个四边形小面元的后向散射系数.遍历完扫描场景中的所有小面元,就会得到整个场景的后向散射特性分布.1.3 阴影区域的判断雷达电磁波在空间中沿着直线传播,传播过程中会遇到山峰等地形结构的遮挡.然而对地形结构遮挡区域的判断是一个相当复杂的过程,它将导致回波仿真的计算量急剧增加,尤其是出现多重遮挡时,遮挡区域的判断变得更为复杂.目前,解决阴影区域判断问题的较好算法是光线追踪算法和下视角比较算法.光线追踪算法原理是由光源发出的每一条光线遍历所有目标以此来判断照亮区域与阴影区域.然而,当三维地面剖分的小面元数量达到上百万量级时,其算法运算效率较低;下视角比较法原理通过对波束进行细划分,对波束照射下不同目标的下视角进行比较,如不满足判断要求,则认为遮挡,其原理简单且计算量小.为了实现三维场景合成孔径雷达回波快速仿真,下视角比较算法无疑是一种合适的选择.然而,常规的下视角比较算法是建立在雷达发射的电磁波为平面波假设条件下的.在电磁波为球面波假设下,常规的下视角比较算法将不能适用,因此笔者提出了改进的下视角比较法.在电磁波为球面波假设下,雷达波束在三维空间的照射示意图如图1(b)所示,理论上雷达波束在地面上投影区域为椭圆形.为了后期图形处理单元并行处理的方便,采用如图1(c)所示的扇环区域来确定雷达波束边界,并对扇环区域的雷达波束进行两维划分.设雷达波束中心下视角为β,中心斜距为Rs,方位向波束角宽度为θBW,距离向宽度为θRW;小面元地面投影的方位向采样间隔为ρx,距离向采样间隔为ρg,则雷达波束方位向划分数量为其中,表示向上取整.雷达波束的距离向划分数量为划分后每个波束网格中可近似认为只有一束平行电波穿过[6].接着进行阴影区域的判断,假设当前点为Pm,计算得到当前点的下视角为βm,找到在此方位角度间隔内所有地距小于Pm点的最大下视角为βmax.若βm≤βmax,则当前点被遮挡;否则,当前点没有被遮挡.为了方便描述,引入遮挡函数:其中,φ为当前点相对波束中心线的方位角度,βm为当前点相对波束中心线的下视角.通过上一小节计算出当前小面元的后向散射系数为σi,则在考虑了阴影遮挡后,小面元最终实际的后向散射系数变为σiL(φ,βm).下面给出三维场景阴影区域遍历判断的详细过程:步骤1 利用上述雷达波束划分法对某一方位时刻的雷达波束进行两维划分.步骤2 选定雷达波束范围内的某一方位角度φ,计算地距最近点的下视角β,并令βmax=β,该点的遮挡函数值为1,m=1,指向下一个点.步骤3 计算当前点Pm的下视角βm,比较βmax与βm,如图1(d)所示.若βm≤βmax,则令遮挡函数值L(φ,βm)为0;否则,令遮挡函数值L(φ,βm)为1;最后令βmax=max(βm,βmax),m值累加1,指向下一点.步骤4 循环执行步骤3,直至m=Nr.选定下一方位角度φ.步骤5 循环执行步骤2～步骤4,直至遍历完所有的方位角度后结束,完成了三维场景阴影区域的判断.2.1 并行性分析针对三维地面场景合成孔径雷达回波仿真,常规仿真算法运算效率仍然很慢.为了提高运算效率,笔者结合图形处理单元的编程架构对回波仿真算法进行并行设计优化,以满足工程应用中的实时信号处理需求.图形处理单元的编程架构是以大量线程来实现高吞吐量的数据并行运算的,线程间的独立性越强,可并行效果越明显[9-11].笔者基于常规同心圆算法进行三维场景合成孔径雷达回波仿真,其处理不同方位采样时刻的雷达回波是相互独立的,处理每个点目标的回波也是相互独立的,因此非常适合用图形处理单元进行并行计算,可以极大地提高运算效率,满足实时处理的需求.通过分析笔者提出的算法,发现其存在3个层次的并行.2.1.1 不同方位时刻生成合成孔径雷达回波的并行性首先雷达波束照射三维场景,然后接收三维场景反射的雷达回波.不同时刻间雷达接收回波是完全相互独立的,即每个方位时刻都会彼此独立地产生一个雷达回波,可以采用图形处理单元中的流思想来加速不同方位时刻合成孔径雷达回波生成的过程.2.1.2 某方位时刻计算不同小面元后向散射系数的并行性三维地面场景可以剖分成上百万个小面元,计算每个小面元后向散射系数的过程是相互独立的,存在高度的可并行性,大大地节省了计算时间.2.1.3 某方位时刻不同雷达波束划分区域内遮挡区域判断的并行性当把雷达波束划分角度后,在每个角度划分区域内采用改进的下视角比较法判断阴影区域.不同角度区域内的阴影判断是相互独立的,存在较高程度的并行性.判断阴影后,得到了所有小面元实际的后向散射系数,将其累加得到场景的回波,这个过程可以采用“归约相加”的思想[8]实现并行处理.2.2 并行化方案通过以上并行性分析,可知第2、3层次的并行度要远远高于第1层次的并行度.因此,在图形处理单元处理核函数设计中,针对并行度较高的第2、3层次来设计3个核函数:第1个核函数负责计算小面元的后向散射系数;第2个核函数负责进行阴影区域的判断;第3个核函数则负责将所有小面元的回波“归约相加”为一个雷达回波.针对第1层次,采用了图形处理单元中的流处理,加速了第1层次并行处理.在三维地面场景的合成孔径雷达回波仿真过程中,为了进一步提高仿真效率,笔者以两个图形处理单元进行设计,并行化处理流程图如图2所示.首先在中央处理器(Central Processing Unit,CPU)中设定雷达参数并分配物理内存,读入数字高程模型数据并进行分形插值处理,计算得到三维场景中所有小面元的平面信息,并通过PCI-E总线传到图形处理单元的显存空间;接着在图形处理单元中用多线程并行扫描场景中所有的小面元,并计算其斜距以及后向散射系数(其中,图中编号为0的图形处理单元处理方位时刻是偶数的数据,编号为1的图形处理单元负责方位时刻是奇数的数据);然后,在图形处理单元中对雷达波束进行角度划分,完成不同波束角度内的阴影遮挡判断,并结合小面元的后向散射系数和遮挡判断结果,得到小面元的实际回波;采用“线程外推”[8]的方法实现部分小面元回波的同心圆累加,接着分别在两个图形处理单元中利用“归约相加”思想[8]实现所有的小面元回波的同心圆累加;在所有方位时刻累加完成后,转到距离频域乘上中央处理器产生的调频项,再通过快速傅里叶逆变换(Inverse Fast Fourier Transform,IFFT)转到距离时域,即得到雷达回波数据;最后,结合图形处理单元流技术计算所有方位时刻的雷达回波,并将两个图形处理单元设备端的回波数据传回中央处理器,并由中央处理器对处理结果进行保存和输出,从而得到三维地面场景的合成孔径雷达回波数据.为了验证笔者提出方法的有效性及并行设计的实效性,利用实测数字高程模型数据进行仿真分析.仿真选用机载雷达C波段参数和Tesla K20C显卡的仿真平台.实测数字高程模型数据分形插值后的三维图如图3所示,其水平维数是2 048×2 048,采样间隔是1 m×1 m,场景最大高度为357.4 m,最低高度为5.3 m,场景斜距是20 km.为了对比阴影判断的有效性,预设了4条雷达载机飞行轨迹,不同轨迹的基本信息如表1所示.表中的H代表载机飞行高度,β代表波束中心线下视角,θ代表飞行轨迹与X轴的夹角.为了说明合成孔径雷达回波仿真方法的高效性,对沿轨迹1飞行时照射真实数字高程模型数据的合成孔径雷达回波仿真进行了计时,此时合成孔径雷达录取回波方位向上有1 608个采样点,距离向上有2 048个采样点,每生成一个脉冲的回波数据,需扫描2 048×2 048(即4 194 304)个小面元并进行遮挡区域的判断.在双TeslaK20C显卡上,经过10次实验,生成这组回波数据的平均测量时间为131.01 s.因此,可以看出笔者提出的方法的快速处理性能强.为了验证笔者提出的方法产生回波的有效性,对得到的回波数据进行常规的成像处理,得到如图4所示的成像结果.其中,图4(a)～(d)分别是载机在飞行轨迹1～轨迹4下的合成孔径雷达回波的成像结果,4幅图像均做了地距校正.对比图4(a)与图4(b)、图4(c)与图4(d)发现,在雷达高度、波束中心的下视角相同的情况下,由于飞行方向不同,导致阴影效果不同,所以成像结果中有着细节差异;对比图4(a)与图4(c)、图4 (b)与图4(d)发现,由于飞行方向相同,两幅图的纹理基本相同,不过由于载机飞行高度不一致,成像结果的阴影效果有所差异,下视角为60°时的阴影效果更加明显.通过以上仿真实验,充分说明了笔者提出的方法的有效性.笔者主要讨论了三维地面场景合成孔径雷达回波仿真方法及其并行优化设计问题.首先,基于分形布朗运动模型对三维地面场景数据进行了分形插值处理,对插值后的数据进行四边形小面元剖分,并计算小面元的后向散射系数;接着利用改进的下视角比较法对三维场景的阴影遮挡进行遍历判断;最后结合同心圆算法和图形处理单元的并行设计实现了三维地面场景合成孔径雷达回波的快速仿真.笔者的并行设计方案优势明显,可快速地实现大范围的三维地面场景的回波仿真,很好地解决运算效率的难题.笔者提出的方法在复杂地形探测、景象匹配制导中具有重大应用前景.【相关文献】[1]邢孟道,保铮,李真芳,等.雷达成像算法进展[M].北京:电子工业出版社,2014:8-25.[2]van der LAAN W J,JALBA A C,ROERDINK J B T M.Accelerating Wavelet Lifting on Graphics Hardware Using CUDA[J].IEEE Transactions on Parallel and Distributed Systems,2011,22(1):132-146.[3]SUN G,RANSON K J.A Three-dimensional Radar Backscatter Model of Forest Canopies[J].IEEE Transactions on Geoscience and Remote Sensing,1995,33(2):372-382.[4]任三孩.三维场景SAR回波模拟技术研究[D].长沙:国防科学技术大学,2007.[5]孙晗伟,胡程,曾涛.一种三维森林场景极化SAR数据的快速模拟方法[J].电子与信息学报,2012,34(6): 1297-1304. SUN Hanwei,HU Cheng,ZENG Tao.A Fast Method of Polarimetric SAR Data Simulation for Three-dimension Forest Stand[J].Jounal of Electronics&Information Technology,2012,34(6):1297-1304.[6]KRIEGER G,MOREIRA A,FIEDLER H,et al.TanDEM-X:a Satellite Formation for High-resolution SAR Interferometry[J].IEEE Transactions on Geoscience and Remote Sensing,2007,45(11):3317-3341.[7]DUDGEON J E,GOPALAKRISHNAN R.Fractal-based Modeling of 3D Terrain Surfaces[C]//Proceedings of the 1996 IEEE SoutheastconConference.Piscataway:IEEE,1996:246-252.[8]LLOYD D B,BOYD C,GOVINDARAJU N.Fast Computation of General Fourier Transforms on GPUs[C]// Proceedings of the 2008 IEEE International Conference on Multimedia and Expo.Piscataway:IEEE,2008:5-8.[9]JIN X,KO S B.GPU-based Parallel Implementation of SAR Imaging[C]//Proceedings ofthe 2012 International Symposium on Electronic System Design.Washington:IEEE Computer Society,2012:125-129.[10]HWANG K,SUNG W.Load Balanced Resampling for Real-time Particle Filtering on Graphics Processing Units[J]. IEEE Transactions on Signal Processing,2013,61(2):411-419.[11]武勇,王俊,张培川,等.CUDA架构下外辐射源雷达杂波抑制并行算法[J].西安电子科技大学学报,2015,42(1): 104-111. WU Yong,WANG Jun,ZHANG Peichuan,et al.Parallel Clutter Suppression Algorithm for Passive Radar in CUDA [J].Journal of XidianUniversity,2015,42(1):104-111.。
初至反演近地表速度的建模方法分析与应用
ｄ ｏ ｉ ： １ ０ ． ３ ９ ６ ９ ￣ ． ｉ ｓ ｓ ｎ ． １ ６ ７ ４ － １ ８ ０ ３ ． ２ ０ １ ３ ． ０ ５ ． １ ２
文 章编 号 ： １ ６ ７ ４ 一 ｌ ８ ０ ３ （ ２ ０ １ ３ ） ０ ５ — ０ ０ ４ ６ — ０ ４
ｗｈ ｉ ｌ ｅ ｔ ｈ ｅ ｐ ｒ ｉ ｍａ ｒ ｙ ｗａ ｖ ｅ ｔ ｏ ｍｏ ｇ ｒ ａ ｐ ｈ ｉ ｃ ｉ ｎ ｖ ｅ ｒ ｓ ｉ ｏ ｎ ｍｅ ｔ ｈ ｏ ｄ ｃ ａ ｎ ｅ ｓ ｔ ｉ ｍａ ｔ ｅ ｓ ｔ ａ ｂ ｌ ｅ ａ ｎ ｄ ｒ ｅ ｌ ｉ ａ ｂ ｌ ｅ ｓ ｕ ｆ ｒ ａ ｃ ｅ ｖ ｅ ｌ ｏ ｃ ｉ ｔ ｙ ｔ ｅ ｘ ｔ ｕ ｒ ｅ ｂ ｅ ｔ ｔ ｅ ｒ ．Ａ ｐ ｒ ｏ ｓ ｐ ｅ ｃ ｔ ｉ ｎ ｇ ｌ ｉ ｎ ｅ ｏ ｆ Ｔ Ｂ ｇ ａ ｓ ｉ ｆ ｅ ｌ ｄ ｉ ｎ Ｉ ｒ ａ ｎ ｌ ｉ ｅ ｓ ｔ ｈ ｒ ｏ ｕ ｇ ｈ ａ ｍｏ ｕ ｎ ｔ ａ ｉ ｎ ｐ ｅ ａ ｋ ， ｗｉ ｔ ｈ ｅ ｌ ｅ ｖ ａ ｔ ｉ ｏ ｎ ｄ ｉ ｆ ｆ ｅ ｒ ｅ ｎ ｃ ｅ ｕ ｐ ｔ ｏ ｌ Ｏ ０ ０ ｍ， ａ ｎ ｄ ｓ ｔ ｒ ｏ ｎ ｇ ｔ ｒ ａ ｎ ｓ ｖ ｅ ｒ ｓ ｅ ｉ ｎ ｓ ｔ ａ ｂ ｉ ｌ ｉ ｔ ｙ ｏ ｆ ｌ ｏ ｗ ａ ｎ ｄ ｄ ｅ ｃ ｅ ｌ ｅ ｒ ａ ｔ ｉ ｎ ｇ ｖ ｅ ｌ ｏ ｃ ｉ ｔ ｙ ｚ ｏ ｎ ｅ ｓ ．Ｕｎ ｄ ｅ ｒ ｗ ｉ ｔ ｈ ｏ ｕ ｔ ａ ｐ ｒ ｉ ｏ ｉ ｒ ｇ ｅ ｏ ｌ ｏ ｇ ｉ ｃ ａ ｌ ｉ ｎ ｆ ｏ ｒ ｍａ ｔ ｉ ｏ ｎ ｓ ｉ ｔ ｕ ａ ｔ ｉ ｏ ｎ ，ｉ ｎ ｔ ｅ ｐｒ ｒ ｅ ｔ ａ ｔ ｉ ｏ ｎ ａ ｍｂ ｉ ｇ ｕ ｉ ｔ ｙ ｏ ｆ ｌ ｏ ｗ ａ ｎ ｄ ｄ ｅ ｃ ｅ ｌ ｅ ｒ ａ ｔ ｉ ｎ ｇ ｖ ｅ ｌ ｏ ｃ ｉ ｔ ｙ ｚ ｏ ｎ ｅ ｓ ｂ ｏ ｔ ｔ ｏ ｍ ｂ ｏ ｕ ｎ ｄ ａ ｒ ｙ ｒ ｆ ｏ ｍ ｒ ｅ ｆ ｒ ａ ｃ ｔ ｉ ｏ ｎ ｍｅ ｔ ｈ ｏ ｄ ｉ ｓ ｈ ａ ｒ ｄ ｔ ｏ ｅ ｌ ｉ ｍｉ ｎ ａ ｔ ｅ ；ｗｈ ｉ ｌ ｅ ｔ ｈ ｅ ｖ ｅ ｌ ｏ ｃ ｉ ｔ ｙ ｔ ｒ ａ ｎ ｓ ｖ ｅ ｒ ｓ ｅ ｖ ａ ｉ ｒ ａ ｔ ｉ ｏ ｎ
复杂地表地震勘探中的高精度近地表速度建模方法
复杂地表地震勘探中的高精度近地表速度建模方法作者：李稳来源：《中国科技博览》2016年第23期[摘 ;要]在勘察内蒙阿拉善某高放废物地质处置预选场址目标花岗岩体完整性与稳定性的地震勘探过程中，利用基于MSFM（Multistencils Fast Marching Methods）射线追踪算法的地震初至波走时反演方法进行了近地表高精度速度建模研究。
获得了准确、精细的近地表速度模型成果。
该成果不但可被用于解决对于陆上复杂地表地震勘探来说至关重要的层析静校正问题，而且可作为地震解释工作中重要的参考资料。
[关键词]高放废物地质处置; 复杂地表; 近地表速度建模; 射线追踪; 速度反演中图分类号：P5 文献标识码：A 文章编号：1009-914X（2016）23-0013-02根据以往的探测经验可知，近地表速度不准确会严重影响地下地质结构的地震成像效果，对于地表条件复杂的工区更是如此。
在勘察内蒙阿拉善某高放废物地质处置预选场址目标花岗岩体的过程中，为了解决对于陆上工区地震勘探至关重要的层析静校正问题，利用基于MSFM （Multistencils Fast Marching Methods，多模板快速步进）射线追踪算法的地震初至波走时反演技术，建立了沿地震测线准确、精细的近地表速度模型。
该成果不但被用于有效解决了后续的层析静校正问题，而且是最终地震地质解释工作中精确、可靠的参考资料。
1 方法原理复杂地表问题通常指包含剧烈起伏的崎岖地表和强速度变化的近地表速度结构。
在选择适用于复杂地表初至波走时反演的射线追踪算法时，通常需要考虑计算精度、计算效率、以及算法对于模型复杂度的适应能力三个方面[1]。
1.1 MSFM走时计算公式MSFM（Multistencils Fast Marching Methods）多模板快速步进算法是由Hassouna等提出的对于经典FMM方法[2]的改进算法[3]。
其原理是在标准FMM算法的基础上通过坐标旋转产生多个FMM计算模板，使对角方向邻点也能够参与计算。
直达波反射波联合层析深度域速度建模方法
直达波反射波联合层析深度域速度建模方法直达波反射波联合层析深度域速度建模方法是地球物理探测领域中的一种重要技术，其主要作用是通过分析地下介质的物理特征，实现地下结构的精细化解析，为资源勘探和生产提供重要依据。
下面将分步骤阐述直达波反射波联合层析深度域速度建模方法的具体流程。
1. 数据预处理首先，需要对数据进行预处理，以确保数据质量。
数据预处理包括对数据进行采样、滤波、去除噪声等处理，可通过MATLAB等工具进行实现。
2. 确定速度模型在进行层析速度建模之前，需要确定速度模型。
可以通过VTI模型或者HTI模型来建立速度模型。
其中，VTI模型是指在纵向方向上具有各向同性，而在水平方向上具有各向异性的模型。
HTI模型则是在VTI模型的基础上，增加了沿着某一方向的轴对称性。
在实际应用中，需要根据地质结构的实际情况，选择合适的速度模型。
3. 建立初始速度模型在确定速度模型之后，需要建立初始速度模型。
通常可以采用叠前时距叠加技术来建立初始速度模型。
4. 反射波速度层析建模在建立初始速度模型之后，将反射波数据输入层析程序中，采用反射波速度层析建模方法。
该方法的主要思想是将射线从注入点投射至接收器，形成数据矩阵。
然后，借助叠前时距叠加技术，对数据进行处理，得到反射波速度层析图像。
5. 直达波速度层析建模反射波速度层析建模得到的速度模型并不完整，因此需要借助直达波数据进行补充。
直达波速度层析建模的方法与反射波速度层析建模的方法类似，只是在数据的处理方式上有所不同。
6. 联合反射波和直达波速度层析建模通过反射波速度层析和直达波速度层析建模，得到了两种不同的速度模型。
为了得到更加精确的速度模型，需要将两个速度模型进行联合。
联合的方法通常采用全波形反演或者模型耦合技术。
7. 模型评价联合反射波和直达波速度层析建模得到的速度模型需要进行验证和评价。
通常采用类背景波形剖面对模型进行测试，并进行统计和分析，从而得到最终的速度模型。
直达波反射波联合层析深度域速度建模方法
直达波反射波联合层析深度域速度建模方法近年来，随着人类在地震研究方面的研究取得了长足的进步，地球物理学家可以获得地壳结构的深度域速度建模的更为准确的估计。
这主要是通过利用波场全息联合层析技术来实现的，该技术允许将多重数据（如波动破裂速度，地壳密度和波振幅转换率）合并起来，使其可以更准确地估计地壳的深度域速度建模。
然而，在此层析过程中，反射波和直达波将会重叠，而这些重叠可能会影响联合层析技术得出的结果。
因此，有必要研究如何准确地分离反射波和直达波，以便更准确地估计深度域的速度建模。
针对这一问题，我们提出了一种新的“直达波反射波联合层析深度域速度建模方法”，它可以将反射波和直达波的重叠效果进行精确的消除。
这一方法的基本思想是，首先，针对已知的反射波数据，使用层析法确定反射波的地质参数；其次，通过深度域建模方法求解该地质参数对应的潜在反射波面；最后，再利用给定的反射波和直达波数据，完成反射波和直达波重叠效果消除。
为了验证这一方法的有效性，我们编写了相应的计算程序，并进行了一系列的测试。
实验结果表明，此方法能够有效地消除反射波和直达波的重叠效果，因此可以更准确地估计深度域的速度建模。
综上所述，我们提出了一种新的“直达波反射波联合层析深度域速度建模方法”，该方法能够更准确地估计深度域的速度建模，从而为科学家印证地壳层析结构提供了准确的依据。
虽然这项研究初期仍然存在一些缺陷，但它仍是地震研究领域的重要成果。
未来，可以使用更复杂的方法或更精确地分离反射波和直达波，他们可以更好地研究地质结构，从而更准确地估计地壳层析结构。
本文提出了一种新的“直达波反射波联合层析深度域速度建模方法”，它能够准确地消除反射波和直达波重叠效果，从而更准确地估算深度域速度建模。
本文探讨了这一新技术的原理，提出了相应的具体算法，并利用给定的数据进行了测试，以证明该方法的可行性。
本文的研究可以为地震研究领域的实践提供指导，并为更准确地估算地壳层析结构提供科学依据。
相控阵主动源近地表速度建模方法
相控阵主动源近地表速度建模方法
相控阵主动源呢，就像是一个超级厉害的小助手，在探测近地表速度方面有它独特的本事。
那这个速度建模为啥重要呢？你想啊，如果不知道近地表的速度情况，就好像你在一个陌生的地方走路，却不知道哪里是快行道，哪里是慢行道一样。
对于相控阵主动源近地表速度建模，有一种方法是通过发射特殊的信号。
这些信号就像一个个小信使，跑到地下各个地方，然后再把信息带回来。
根据信号传播的时间和距离，就能大概算出速度啦。
比如说，信号跑出去再回来花了多少时间，我们知道它跑的路程，速度不就出来了嘛。
还有一种方法是利用多组数据的对比。
就像我们看一群小朋友赛跑，有的跑得快，有的跑得慢，通过对比不同小朋友的速度情况，我们就能对整个赛场的速度分布有个概念。
在相控阵主动源这里，就是对比不同方向、不同位置的信号反馈数据，然后构建出近地表速度模型。
不过呢，这个过程也不是一帆风顺的。
有时候会受到外界干扰，就像你在听别人说话的时候，突然有个大喇叭在旁边响，干扰你听清楚内容。
在速度建模的时候，像地下的一些复杂结构或者其他的电磁干扰之类的，都会影响数据的准确性。
但是科学家们可有办法啦。
他们会对数据进行筛选和修正。
就像我们挑水果一样，把坏的挑出去，留下好的。
把不准确的数据找出来，然后用一些巧妙的算法进行修正，这样就能得到更靠谱的近地表速度模型啦。
多信息融合的近地表速度建模技术及应用_赵玲芝
２０１７年２月第５２卷 第１期 ＊河北省涿州市华阳东路东方地球物理公司科技园研究院处理中心，０７２７５０。
Ｅｍａｉｌ：８４３４２２５０８＠ｑｑ．ｃｏｍ本文于２０１６年５月３日收到，最终修改稿于同年１０月２８日收到。
·处理技术·文章编号：１０００－７２１０（２０１７）０１－００３４－０８多信息融合的近地表速度建模技术及应用赵玲芝＊①② 谷跃民② 张建中①（①中国海洋大学海洋地球科学学院，山东青岛２６６１００；②东方地球物理公司研究院，河北涿州０７２７５０）赵玲芝，谷跃民，张建中．多信息融合的近地表速度建模技术及应用．石油地球物理勘探，２０１７，５２（１）：３４－４１．摘要 为了提高复杂区近地表速度模型的构建精度，提出并实现了基于多信息融合的近地表速度模型构建技术。
采用协同克里金的插值方法，根据初至反演的近地表速度与近地表调查速度的相关性及变化趋势，确定对应的变差函数和用于插值的权系数，对速度变量进行协方差的空间估计，使横向的速度变化趋势既与大炮初至层析反演的结果较吻合，数值上又与近地表调查的结果基本一致。
该方法可以提高复杂近地表速度模型的精度。
在中国西部复杂区地震资料的应用取得了较好的效果，进一步验证了该方法的有效性。
关键词 大炮初至 近地表调查 近地表模型 协同克里金插值中图分类号：Ｐ６３１ 文献标识码：Ａ ｄｏｉ：１０．１３８１０／ｊ．ｃｎｋｉ．ｉｓｓｎ．１０００－７２１０．２０１７．０１．００６１ 引言近地表速度模型的精度直接影响静校正及叠前深度偏移成像的效果，如何得到高精度的近地表速度模型一直以来是陆上地震勘探的难题之一。
为了获得高精度的近地表速度模型，许多工区开展了微测井和小折射等近地表调查，大炮初至层析反演技术也得到了广泛的应用［１］。
对复杂区近地表速度建模的结果研究表明：近地表调查的速度基本代表了该位置近地表的垂向速度，精度较高，但平面上呈散点分布，横向密度低，分布不均匀，不易刻画速度的空间变化［２］；而大炮初至层析反演技术所反演的速度不只受炮点位置速度的影响，还受炮检距大小、射线追踪的空间分布特征及速度的空间变化等影响，反演的速度模型可以作为真实模型的“等效模型”，能反映近地表速度的变化趋势，但数值上与微测井、小折射的速度存在一定的差异。
landmark中depthteam速度场的建立及变速成图
地球物理勘探中，地下介质的速度是个非常关键的因素。
准确的速度求取一直是地震勘探的核心问题，速度求取的准确与否直接影响着地震勘探的各个环节及最终成果。
只有得到准确的速度，才能准确的确定地层的深度位置，以及根据速度研究岩石和空隙液体的性质。
1.1，国内外研究现状在过去，国内外一直致力于依据叠加速度或偏移速度借助于Dix公式将其转换成平均速度或层速度建立速度场。
但Dix公式是在地下水平层状均匀介质，射线垂直入射的前提下建立起来的速度关系式。
均匀介质或水平层状介质是水平叠加和速度分析的前提条件，而各向异性较强地区的速度分析偏离了这个条件。
因此，当地层倾角较大，速度横向变化大时，用此方法建立的速度场精度很低。
同时，对变速成图的研究主要集中在速度场建模和叠前偏移成像领域。
变速成图的常规作法主要有量版法、叠偏剖面法、人工T0图空校法和分区空校法等。
改进型的方法主要有等效直射线法、曲射线法、射线追踪法、叠偏剖面成图法和倾斜地表空校法等。
速度研究中，地球物理学家曾经提出了用射线追踪相干反演法通过叠前CMP道集反演层速度建立速度场。
对于二维来讲这种方法比较经济适用，同时从很大程度上解决了二维工区中层速度及深度的三维空间归位问题。
对于三维，与Dix转换相比，反演出的层速度精度更高。
同时，近年来也提出了在多井约束下反演层速度，建立地下速度场，进行变速构造成图。
地质统计法综合考虑二维平面的速度趋势，同时在井点处又严格遵循井上的硬数据，即将地震和钻井的速度有机地结合在一起，从而提高模型精度，降低勘探成本。
它将大量的迭加速度谱与少量的井速度资料相结合，根据协克里金地质统计综合分析技术，推导出地下速度的空间分布规律。
目前速度反演方法，国内外广为采用的是基于斯奈尔定理的二维和三维空间的射线追踪逐层反演方法。
其追踪方法以试射法和弯曲法为主，以试射法最为流行。
但是这种方法受界面影响比较大，存在着射线不能达到接收点的可能性，同时对逆断层区内的反射点不能有效追踪。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
快速回转波近地表速度建模方法郭振波;孙鹏远;钱忠平;李培明;唐博文;熊定钰【摘要】Conventional first-arrival travel time tomography based on iterative inversion is computationally intensive and timeconsuming,especially in a large-scale data processing.In order to improve the efficiency of near-surface modeling,a fast near-surface model building with turning wave is developed by assuming that the near-surface velocity is transversely constant within the maximum offset range and varies linearly with the depth.In order to improve the accuracy of inversion,a multi-datum correction method is adopted to reduce the influence of surface undulation on the inversion results.In order to enhance the stability of inversion,a local weighted ray parameter estimation is developed to increase the robustness with respect to noise.Synthetic and real data tests verify the validity and efficiency of the proposed method.%常规基于迭代反演的初至波层析近地表速度建模方法在大规模数据处理时计算量较大、耗时多.为了提高近地表建模的效率,发展了一种回转波快速近地表建模方法.该方法假设在最大炮检距范围内近地表速度横向不变但随深度线性变化.为了提高反演精度,在算法内部采用多基准面校正方法降低地表起伏对反演结果的影响;为了增强反演的稳定性,应用基于局部加权的稳定射线参数估计方法降低初至拾取误差对反演结果的影响.理论模型数据与实际数据测试验证了方法的有效性及高效性.【期刊名称】《石油地球物理勘探》【年(卷),期】2019(054)002【总页数】8页(P261-267,前插1)【关键词】近地表;层析;速度建模;反演;静校正;回转波【作者】郭振波;孙鹏远;钱忠平;李培明;唐博文;熊定钰【作者单位】东方地球物理公司物探技术研究中心,河北涿州 072751;东方地球物理公司物探技术研究中心,河北涿州 072751;东方地球物理公司物探技术研究中心,河北涿州 072751;东方地球物理公司物探技术研究中心,河北涿州 072751;东方地球物理公司物探技术研究中心,河北涿州 072751;东方地球物理公司物探技术研究中心,河北涿州 072751【正文语种】中文【中图分类】P6310 引言地震勘探由于激发和接收均在地表，因此需要消除地表起伏及近地表速度横向变化对反射信号的影响[1]，其主要处理方式有两种：①通过静校正技术将炮点、检波点校正到统一基准面上[2-3]，在消除近地表对反射波影响的同时也可使一些基于水平地表假设的方法(如水平叠加等)可行；②直接进行基于起伏地表的处理，在算法内部考虑近地表的影响(如起伏地表条件下的速度建模、偏移成像等)[4-5]。
不管采用哪种方式，都需要进行近地表建模，因此近地表建模方法的研究对油气地震勘探具有重要意义。
近地表建模主要利用地震记录中的初至波，可分为基于波动理论和射线理论的层析，也可分为波形层析和旅行时层析[6]。
由于观测误差、噪声等因素的影响，很难直接进行波形匹配，因此目前以旅行时层析为主[7]。
由于计算能力的限制，旅行时层析以基于射线理论的方法为主，而基于波动理论的旅行时层析还未得到大规模应用[8]。
常规的基于射线理论的旅行时层析方法(常规射线层析)[9]是基于射线追踪方程或程函方程、通过迭代求解正演旅行时与观测旅行时的最优拟合，可获取较高精度的近地表模型。
但由于需要大量的正演及反演计算，计算量较大，特别是对于目前的高密度、宽方位地震数据，限制了方法的大规模应用。
除了常规射线层析反演，还有一类方法是通过对速度分布做进一步假设构建快速近地表建模方法，该类方法在天然地震研究中广泛应用，通常用来构建一维参考模型[10]。
Diebold等[11]采用截断时间求和方法利用折射波与反射波旅行时构建一维速度模型；Rühl等[12]利用最大深度方法构建近地表模型并将其应用于静校正，实际数据测试验证了由该方法计算的长波长静校正量的准确性；Nowack[13]、Novotn等[14]、胡自多等[15]、徐涛等[16]利用 Herglotz-Wiechert方法构建近地表模型； Shi等[17]在线性速度分布假设下采用层剥离的方法进行快速近地表建模。
本文通过假设速度随深度线性变化发展了一种快速的回转波近地表建模方法。
通过在算法内部采用多基准面校正降低地表起伏对反演结果的影响，提高反演精度；应用基于局部加权的稳定射线参数估计方法降低初至拾取误差对反演结果的影响，增强反演的稳定性。
理论模型数据与实际数据测试验证了方法的有效性及高效性。
1 方法原理1.1 基本原理假设在最大炮检距范围内近地表速度横向不变但随深度线性变化，此时近地表的速度变化趋势可表述为v(z)=v0+gz(1)式中： v0为地表速度； z为深度； g为速度梯度。
在该介质中传播的地震波称为回转波。
给定射线入射角度θ、起始坐标(x0，z0)可导出其射线路径、炮检距、初至时间的理论表达式(图1)。
具体来说，其射线路径为一圆弧，圆心位于点[18](2)半径为(3)式中p=sinθ/v为射线参数。
相应地可求得对应的炮检距H及初至时间t(4)对于速度反演，已知炮检距、初至时间，求取地下不同位置处的速度参数。
假设速度水平横向不变，地震波的传播射线参数保持不变，可据此求得回转点处的速度为图1 地震波在速度随深度线性变换介质中的传播路径(5)由式(2)和式(3)可求得回转点的深度为(6)射线参数p可由初至时间估算；地表速度v0由近炮检距数据计算或由先验信息给定。
速度随深度变化的梯度通过求解目标函数(7)获得。
式中W为加权因子，调整炮检距项(第一项)与时间项(第二项)的权重。
当炮检距单位为m、时间单位为s时，加权因子W设为1000，以均衡由于不同测量单位固有的加权效应。
利用式(5)～式(7)可获取一个炮检距—初至时间对所对应回转点的速度及其深度。
通过在由近及远的炮检距范围内进行初至反演，可获取由浅到深的速度参数；通过对不同空间位置点进行相同的处理，可获取整个三维空间的速度参数。
给定离散的空间网格，将相应的速度参数投影到对应的网格点上，然后进行对应的速度内插、平滑等处理，可得到规则网格点上最终的速度参数。
1.2 多基准面初至时间校正由于上述方法基于水平地表假设，起伏地表情况下需要将其校正到水平基准面上。
如图2所示，本文通过两方面的策略保持校正精度：①将数据抽至CMP道集后再进行后续的处理，保证炮检点的空间局部性；②根据地表起伏的情况，针对不同的炮检距范围利用不同的基准面进行高程校正，尽量缩小校正时间以减少高程校正带来的误差。
一般情况下，近炮检距范围采用CMP点处的高程作为反演的基准面，远炮检距范围内采用炮检点的平均高程作为反演基准面。
图2 多基准面校正示意图1.3 基于局部加权的稳定射线参数估计射线参数又称为视慢度，可通过估算炮检距—初至时间曲线不同点的斜率获取。
射线参数估计的精度及稳定性直接影响最终的反演结果。
本文提出了一种基于局部加权的稳定射线参数估计方法，可有效降低初至拾取误差、异常值的影响，实现相对稳定的射线参数估计。
对于一个CMP道集，利用一个固定宽度的窗口沿炮检距—初至时间进行滑动，对平滑窗内的数据进行直线拟合[19]，所拟合直线的斜率即为该平滑窗中心位置对应的射线参数，如图3a所示。
利用这种直线拟合的方式可消除小的拾取误差对射线参数估计的影响，但是对于大的拾取误差容错能力较弱。
借鉴Cleveland[20]局部加权平滑的相关思想，本文在直线拟合的过程中引入局部加权减少大的拾取误差对射线参数估计的影响。
基于局部加权的直线拟合，可以表述为求解目标函数[20](8)式中： N为平滑窗内的点数； tk、Hk为第k个点处的初至时间与炮检距； t0、p分别为所拟合直线的截距时间与斜率(即射线参数)； wk为第k个点处的加权系数。
稳定的参数估计方法主要分为两步：①将所有点处的加权系数wk 设为1，首先进行常规线性拟合；②计算所有点与拟合直线的距离，根据给定的可容许的最大时间差Δtmax，重新设置其加权系数为(9)利用新的加权系数再进行一次线性拟合，消除异常初至对射线参数估计的影响，如图3b所示。
图3 射线参数估计示意图(a)射线参数整体估计方法示意图； (b)基于局部加权的拟合方法(线2)与常规拟合方法(线1)的对比2 理论数据测试为了说明本文方法的有效性，选用Amoco静校正基准测试模型[21](图4)进行理论数据的测试。
该模型包含了大部分常见的近地表地质构造，如高速层出露(区域A)、局部高速、低速异常体(区域B)、浅层低速层(区域C)、近地表复杂构造(区域D)及极浅层低速体(区域E)等，可在一定程度上说明该方法对不同近地表地质构造的适应性。
图4 Amoco静校正基准测试模型原有模型横向范围为50km，纵向深6km，由于本文只进行近地表模型的建立，因此只截取了浅部3.2km的部分模型进行测试。
模型纵、横向采样间隔均为5m，横向共10001个采样点，纵向共647个采样点，最大速度为5760 m/s，最小速度约为800m/s。
观测数据共1998炮，炮点以25m的间隔均匀分布于地表以下10m处，第一个炮点位于横向位置10m处。
按照陆上采集方式布设观测系统，最大炮检距为7.5km，检波点间距为10m。
为了避免初至拾取巨大的工作量以及拾取误差对分析的影响，利用有限差分求解程函方程计算得到各个检波点处的初至时间作为反演的输入。
本文方法反演的速度模型如图5所示。
对比图4与图5可知，本文方法除了在复杂构造区域(图4中D标识区域)反演效果较差外，其余部分反演结果与真实模型形态均相似。