地震反演基础知识2(代码演示)

数据集代码演示

1. SEG盐真实数据

# 绘制SEG盐层数据的地震图像
def pain_seg_seismic_data(para_seismic_data):
    '''
    Plotting seismic data images of SEG salt datasets
    :param para_seismic_data:   Seismic data (400 x 301) (numpy)
    :param is_colorbar:         Whether to add a color bar (1 means add, 0 is the default, means don't add)
    '''

    fig, ax = plt.subplots(figsize=(6.2, 8), dpi=120)

    im = ax.imshow(para_seismic_data, extent=[0, 300, 400, 0], cmap=plt.cm.seismic, vmin=-0.4, vmax=0.44)

    ax.set_xlabel('Position (km)', font21)
    ax.set_ylabel('Time (s)', font21)

    ax.set_xticks(np.linspace(0, 300, 5))
    ax.set_yticks(np.linspace(0, 400, 5))
    ax.set_xticklabels(labels=[0, 0.75, 1.5, 2.25, 3.0], size=21)
    ax.set_yticklabels(labels=[0.0, 0.50, 1.00, 1.50, 2.00], size=21)

    plt.rcParams['font.size'] = 14  # Set colorbar font size
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("top", size="3%", pad=0.32)
    plt.colorbar(im, ax=ax, cax=cax, orientation='horizontal')

    plt.subplots_adjust(bottom=0.08, top=0.98, left=0.11, right=0.99)

    plt.show()

在这里插入图片描述

2. SEG盐速度模型

def pain_seg_velocity_model(para_velocity_model):
    '''
    :param para_velocity_model: Velocity model (200 x 301) (numpy)
    :param min_velocity:        Upper limit of velocity in the velocity model
    :param max_velocity:        Lower limit of velocity in the velocity model
    :return:
    '''

    fig, ax = plt.subplots(figsize=(5.8, 4.3), dpi=150)
    im = ax.imshow(para_velocity_model, extent=[0, 3, 2, 0])

    ax.set_xlabel('Position (km)', font18)
    ax.set_ylabel('Depth (km)', font18)
    ax.tick_params(labelsize=14)

    plt.rcParams['font.size'] = 14  # Set colorbar font size
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("top", size="3%", pad=0.32)
    plt.colorbar(im, ax=ax, cax=cax, orientation='horizontal')
    plt.subplots_adjust(bottom=0.12, top=0.95, left=0.11, right=0.99)

    plt.show()

在这里插入图片描述

3. SEG盐模拟地震数据

def pain_seg_seismic_data(para_seismic_data):
    '''
    Plotting seismic data images of SEG salt datasets
    :param para_seismic_data:   Seismic data (400 x 301) (numpy)
    :param is_colorbar:         Whether to add a color bar (1 means add, 0 is the default, means don't add)
    '''

    fig, ax = plt.subplots(figsize=(6.2, 8), dpi=120)

    im = ax.imshow(para_seismic_data, extent=[0, 300, 400, 0], cmap=plt.cm.seismic, vmin=-0.4, vmax=0.44)

    ax.set_xlabel('Position (km)', font21)
    ax.set_ylabel('Time (s)', font21)

    ax.set_xticks(np.linspace(0, 300, 5))
    ax.set_yticks(np.linspace(0, 400, 5))
    ax.set_xticklabels(labels=[0, 0.75, 1.5, 2.25, 3.0], size=21)
    ax.set_yticklabels(labels=[0.0, 0.50, 1.00, 1.50, 2.00], size=21)

    plt.rcParams['font.size'] = 14  # Set colorbar font size
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("top", size="3%", pad=0.32)
    plt.colorbar(im, ax=ax, cax=cax, orientation='horizontal')

    plt.subplots_adjust(bottom=0.08, top=0.98, left=0.11, right=0.99)

    plt.show()

在这里插入图片描述

4. SEG盐模拟速度模型

def pain_seg_velocity_model(para_velocity_model):
    '''
    :param para_velocity_model: Velocity model (200 x 301) (numpy)
    :param min_velocity:        Upper limit of velocity in the velocity model
    :param max_velocity:        Lower limit of velocity in the velocity model
    :return:
    '''

    fig, ax = plt.subplots(figsize=(5.8, 4.3), dpi=150)
    im = ax.imshow(para_velocity_model, extent=[0, 3, 2, 0])

    ax.set_xlabel('Position (km)', font18)
    ax.set_ylabel('Depth (km)', font18)
    ax.tick_params(labelsize=14)

    plt.rcParams['font.size'] = 14  # Set colorbar font size
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("top", size="3%", pad=0.32)
    plt.colorbar(im, ax=ax, cax=cax, orientation='horizontal')
    plt.subplots_adjust(bottom=0.12, top=0.95, left=0.11, right=0.99)

    plt.show()

在这里插入图片描述

5. openfwi地震数据

def pain_openfwi_seismic_data(para_seismic_data):
    '''
    Plotting seismic data images of openfwi dataset
    :param para_seismic_data:   Seismic data (1000 x 70) (numpy)
    '''
    data = cv2.resize(para_seismic_data, dsize=(400, 301), interpolation=cv2.INTER_CUBIC)

    fig, ax = plt.subplots(figsize=(6.1, 8), dpi=120)
    im = ax.imshow(data, extent=[0, 0.7, 1.0, 0], cmap=plt.cm.seismic, vmin=-18, vmax=19)

    ax.set_xlabel('Position (km)', font21)
    ax.set_ylabel('Time (s)', font21)

    ax.set_xticks(np.linspace(0, 0.7, 5))
    ax.set_yticks(np.linspace(0, 1.0, 5))
    ax.set_xticklabels(labels=[0, 0.17, 0.35, 0.52, 0.7], size=21)
    ax.set_yticklabels(labels=[0, 0.25, 0.5, 0.75, 1.0], size=21)

    plt.rcParams['font.size'] = 14  # Set colorbar font size
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("top", size="3%", pad=0.3)
    plt.colorbar(im, ax=ax, cax=cax, orientation='horizontal')
    plt.subplots_adjust(bottom=0.08, top=0.98, left=0.11, right=0.99)

    plt.show()

在这里插入图片描述

6. openfwi速度模型

def pain_openfwi_velocity_model(para_velocity_model):
    '''
    Plotting seismic data images of openfwi dataset
    :param para_velocity_model: Velocity model (70 x 70) (numpy)
    :param min_velocity:        Upper limit of velocity in the velocity model
    :param max_velocity:        Lower limit of velocity in the velocity model
    :param is_colorbar:         Whether to add a color bar (1 means add, 0 is the default, means don't add)
    :return:
    '''

    fig, ax = plt.subplots(figsize=(5.8, 6), dpi=150)

    im = ax.imshow(para_velocity_model, extent=[0, 0.7, 0.7, 0])

    ax.set_xlabel('Position (km)', font18)
    ax.set_ylabel('Depth (km)', font18)
    ax.set_xticks(np.linspace(0, 0.7, 8))
    ax.set_yticks(np.linspace(0, 0.7, 8))
    ax.set_xticklabels(labels=[0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7], size=18)
    ax.set_yticklabels(labels=[0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7], size=18)

    plt.rcParams['font.size'] = 14  # Set colorbar font size
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("top", size="3%", pad=0.35)
    plt.colorbar(im, ax=ax, cax=cax, orientation='horizontal')
    plt.subplots_adjust(bottom=0.10, top=0.95, left=0.13, right=0.95)

    plt.show()

在这里插入图片描述

  • 1
    点赞
  • 4
    收藏
    觉得还不错? 一键收藏
  • 5
    评论
地震波阻抗反演是地球物理勘探中的一种重要手段,其原理是通过测量地震波在不同介质中的传播速度和衰减程度,推断出地下介质的物理参数。下面是一个简单的地震波阻抗反演的MATLAB代码: ```matlab % 地震波阻抗反演代码 clear;clc; % 假设地下介质为两层,分别为地表到100m深度的砂岩和100m到200m深度的泥岩 vp=[1500 2000]; % 砂岩和泥岩的纵波速度,单位 m/s vs=[500 650]; % 砂岩和泥岩的横波速度,单位 m/s rho=[2000 2300]; % 砂岩和泥岩的密度,单位 kg/m3 h=[100 100]; % 砂岩和泥岩的厚度,单位 m f=[5:5:100]; % 反演频率范围,单位 Hz w=2*pi*f; % 角频率,单位 rad/s z=zeros(length(f),1); % 存储阻抗 for i=1:length(f) M=[1 1 0 0;-rho(1)*w(i) rho(1)*w(i) -w(i)^2 0;0 0 1 1;rho(1)*w(i)/vs(1) -rho(1)*w(i)/vs(1) 0 -w(i)^2]; N=[1 1 0 0;-rho(2)*w(i) rho(2)*w(i) -w(i)^2 0;0 0 1 1;rho(2)*w(i)/vs(2) -rho(2)*w(i)/vs(2) 0 -w(i)^2]; M=M*diag([1i*h(1) 1i*h(1) 1i*h(1) 1i*h(1)]); N=N*diag([1i*h(2) 1i*h(2) 1i*h(2) 1i*h(2)]); Z=M\N*[vp(1) vs(1) vp(2) vs(2)]'; z(i)=Z(3); end % 绘制阻抗随频率变化的曲线 figure; plot(f,real(z),'r',f,imag(z),'b'); % 实部为红色,虚部为蓝色 xlabel('Frequency (Hz)');ylabel('Impedance (kg/m^2s)'); legend('Real Part','Imaginary Part'); grid on; ``` 代码中假设地下介质为两层,分别为地表到100m深度的砂岩和100m到200m深度的泥岩,通过输入砂岩和泥岩的纵波速度、横波速度、密度和厚度,以及反演频率范围,计算得到阻抗随频率变化的曲线。
评论 5
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值