Plot elastic FWI Results

Daniel Köhn Kiel, 10/06/2016

Import Libraries


In [33]:
import matplotlib.pyplot as plt
import matplotlib as mpl
import numpy as np
from matplotlib.colors import LightSource, Normalize
from matplotlib.pyplot import gca
from pylab import rcParams
from matplotlib import rc
from matplotlib.ticker import FormatStrFormatter
import pickle

Import Colormap


In [2]:
fp = open('cmap_cm.pkl', 'rb')
my_cmap_cm = pickle.load(fp)
fp.close()

FD grid dimensions


In [3]:
DH = 20.0;
NX = 500;
NY = 174;
vpmin=800.0

Define fonts


In [65]:
FSize = 20
font = {'color':  'black',
        'weight': 'normal',
        'size': FSize}
mpl.rc('xtick', labelsize=FSize) 
mpl.rc('ytick', labelsize=FSize) 
rcParams['figure.figsize'] = 12, 11

Read $V_p$, $V_s$ and density model


In [110]:
stage = 4;

In [121]:
#basename = "elastic_FWT/OBC_FWT/LBFGS/10_06_2016_vx_vy_1D_GRAD_FORM_1_2Hz/"
basename = "elastic_FWT/OBC_FWT/LBFGS/21_05_2015_vx_vy_1D/"

In [122]:
name_vp = basename + "modelTest_vp_stage_" + "%0.*f" %(0,np.fix(stage)) + ".bin"
f = open (name_vp)
data_type = np.dtype ('float32').newbyteorder ('<')
vp = np.fromfile (f, dtype=data_type)
vp = vp.reshape(NX,NY)
vp = np.transpose(vp)
vp = np.flipud(vp)

In [123]:
name_vs = basename + "modelTest_vs_stage_" + "%0.*f" %(0,np.fix(stage)) + ".bin"
f = open (name_vs)
data_type = np.dtype ('float32').newbyteorder ('<')
vs = np.fromfile (f, dtype=data_type)
vs = vs.reshape(NX,NY)
vs = np.transpose(vs)
vs = np.flipud(vs)

In [124]:
name_rho = basename + "modelTest_rho_stage_" + "%0.*f" %(0,np.fix(stage)) + ".bin"
f = open (name_rho)
data_type = np.dtype ('float32').newbyteorder ('<')
rho = np.fromfile (f, dtype=data_type)
rho = rho.reshape(NX,NY)
rho = np.transpose(rho)
rho = np.flipud(rho)

Define Axis


In [14]:
x = np.arange(0.0, DH*NX, DH)
y = np.arange(0.0, DH*NY, DH)
x = np.divide(x,1000.0);
y = np.divide(y,1000.0);

Define SubPlot


In [166]:
def do_plot(n, model, cm, an, title, vpmin, vpmax):
    
    ax=plt.subplot(3, 1, n)
    ax.set_xticks([1, 2, 3, 4, 5, 6, 7, 8, 9, 10])
    ax.set_yticks([0.5, 1, 1.5, 2, 2.5, 3, 3.5])
    
    #plt.rc('font',**{'family':'sans-serif','sans-serif':['Helvetica']})
    rc('font',**{'family':'sans-serif','sans-serif':['Helvetica']})
    ## for Palatino and other serif fonts use:
    #rc('font',**{'family':'serif','serif':['Palatino']})
    #plt.rc('text', usetex=True)
    rc('text', usetex=True)
    
    # plt.pcolor(x, y, vp, cmap=cm, vmin=vpmin)
    plt.imshow(model, cmap=cm, interpolation='none', extent=[0.0,NX*DH/1000.0,0.0,NY*DH/1000.0], vmin=vpmin, vmax=vpmax)
    a = gca()
    a.set_xticklabels(a.get_xticks(), font)
    a.set_yticklabels(a.get_yticks(), font)
    plt.axis('scaled')
    plt.ylabel('Depth [km]', fontdict=font)
    if n==3:
        plt.xlabel('Distance [km]', fontdict=font)
    plt.gca().invert_yaxis()
    cbar=plt.colorbar(aspect=8, pad=0.02)
    cbar.set_label(title, fontdict=font, labelpad=10)
    plt.text(0.1, 0.32,an,fontdict=font,color='white')
    plt.tight_layout()

Plot SubPlots


In [167]:
plt.close('all')
plt.figure()
do_plot(1, vp, my_cmap_cm, '(a)', r"$\rm{V_p [m/s]}$", 950, 4600)
do_plot(2, vs, my_cmap_cm, '(b)', r"$\rm{V_s [m/s]}$", 550, 2700)
do_plot(3, rho, my_cmap_cm, '(c)', r"$\rm{\rho [kg/m^3]}$", 1200, 2800)
#plt.savefig('test.png', format='png', dpi=100)
plt.savefig('test.pdf', format='pdf')
plt.show()

In [ ]: