前言

盐丘轮廓的绘制在地震偏移领域用途十分广泛，主要目的是为了分析偏移结果的好坏。当偏移结果和盐丘的轮廓基本吻合时，偏移质量就很高。最开始绘制轮廓的时候我采用的PS法，通过磁性套索把整个盐丘模型轮廓勾勒出来，然后再保存下来，与偏移结果进行对比。这里介绍一种直接通过 Madagascar 脚本的形式来绘制盐丘的结果，非常方便，可复用性很强，还可以用于其他的方面。此方法也借鉴于Prof. Sava 在 Madagascar book 文件夹下贡献的一个关于 Sigsbee 模型的脚本。

操作

下面通过脚本的形式来说明主要的操作步骤：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38

# 提取盐丘外形
Flow(  'salt','vraw','add add=-4.280 | clip clip=0.300 ')
# vraw为盐丘速度模型,其速度值单位为km/s

#　计算水平梯度和垂直梯度
Flow(  'edgez','salt','         igrad square=y')
Flow(  'edgex','salt','transp | igrad square=y | transp')
Flow(  'edge','edgez edgex',
       'math z=${SOURCES[0]} x=${SOURCES[1]} output="z+x"')

# 对于给定的阈值进行mask，大于阈值部分为1，小于则为0
Flow(  'mask','edge','mask min=0.04')

# 计算盐丘外壳的x与z坐标
for x in ('x1','x2'):
    if(x=='x1'): o='zs'
    if(x=='x2'): o='xs'

    Flow(o,'edge mask',
         '''
         math output=%s |
         put n1=1 n2=%d |
         headerwindow mask=${SOURCES[1]} |
         window
         ''' % (x,par['nz']*par['nx']))

# 合并x和z坐标
Flow('ss',['xs','zs'],
     '''
     cat axis=2 space=n
     ${SOURCES[0]} ${SOURCES[1]} | transp |
     window j2=2

     ''', stdin=0)

fdmod.param(par)
# 绘制外壳图像
Plot(  'ss',fdmod.ssplot('plotfat=2 plotcol=6 symbol=.',par))

plot-salt