In [1]:
using Seismic, PyPlot # required deps
In [ ]:
####################################################################
### FLAT LAYER MODEL ###############################################
####################################################################
In [2]:
# Let's first record some seismic data:
include("../vel/flat.vel_correct.jl")
include("../mod/modeler.jl");
modeler("../dat/vel_correct","../dat/image_correct")
In [3]:
# now let's view an inline out of the "real" seismic data:
seismic, seismic_h = SeisRead("../dat/image_correct"); # load data
SeisPlot(seismic[:,:,50], cmap = "seismic") # plot data
Out[3]:
In [4]:
# Then let's assume we've picked a velocity model
# and converted it to interval velocities.
# We'll assume there are interpretive errors in
# the model.
include("../vel/flat.vel_incorrect.jl")
# An inline from our picked vel model looks like this:
vel_initial, vel_h = SeisRead("../dat/vel_incorrect"); # read data
SeisPlot(vel_initial[:,:,50], cmap = "viridis_r") # plot data
Out[4]:
In [5]:
# Obviously, this model is incorrect or we'd see three
# seismic events. Let's forward model this
# incorrect vel model and view the output modeled seismic:
modeler("../dat/vel_incorrect","../dat/image_incorrect") # model seismic from our bad vels
seis_incorrect, seismic_h = SeisRead("../dat/image_incorrect"); # load data
SeisPlot(seis_incorrect[:,:,50], cmap = "seismic") # plot data
Out[5]:
In [6]:
# Therefore, if we tried to interpret on our bad velocity
# model, we would never find our drilling target. Hence
# we need to improve our velocity model somehow. So let's
# use full waveform inversion.
# Let's input our poorly picked velocity model and the
# seismic data into an FWI algorithm and take a look
# at the updated velocity model to see if there was any
# improvement:
# NOTE: RUNNING FWI MAY TAKE 5 MINUTES OR MORE...
tic() # start timer
include("../fwi/fwi.jl") # prep deps
fwi("../dat/vel_incorrect", "../dat/image_correct",
"../dat/wav", "../dat/updated_vel", .03, 10, 0); # FWI algorithm
toc() # end timer + print
vel_updated, vel_h = SeisRead("../dat/updated_vel") # load data
SeisPlot(vel_updated[:,:,50], cmap = "viridis_r") # plot data
Out[6]:
In [7]:
# We can see here that the vel model has been improved considerably
# by full waveform inversion. FWI has corrected the velocity
# error in out model to the tune of about 3%.
# For comparison, let's view the insitu velocities.
# (This is the Earth model we used to "record" the "real"
# seismic, above):
vel_true, vel_h = SeisRead("../dat/vel_correct"); # read data
SeisPlot(vel_true[:,:,50], cmap = "viridis_r") # plot data
Out[7]:
In [ ]:
# As you can see, the FWI updated vel model is much closer
# to the true velocity model.
In [ ]:
####################################################################
### PYRAMID DOME MODEL #############################################
####################################################################
In [8]:
# Let's first record some seismic data:
include("../vel/pyramid.vel_correct.jl")
include("../mod/modeler.jl");
modeler("../dat/vel_correct","../dat/image_correct")
In [9]:
# First let's view an inline out of the "real" seismic data:
seismic, seismic_h = SeisRead("../dat/image_correct"); # load data
SeisPlot(seismic[:,:,50], cmap = "seismic") # plot data
Out[9]:
In [10]:
# Then let's assume we've picked a velocity model
# and converted it to interval velocities.
# We'll assume there are interpretive errors in
# the model.
include("../vel/pyramid.vel_incorrect.jl")
# An inline from our picked vel model looks like this:
vel_initial, vel_h = SeisRead("../dat/vel_incorrect"); # read data
SeisPlot(vel_initial[:,:,50], cmap = "viridis_r") # plot data
Out[10]:
In [11]:
# Obviously, this model is incorrect. Hypothetically the
# geophysicist who picked this section was afraid to pick
# the velocity inversion. Let's forward model this
# incorrect vel model and view the output modeled seismic:
modeler("../dat/vel_incorrect","../dat/image_incorrect") # model seismic from our bad vels
seis_incorrect, seismic_h = SeisRead("../dat/image_incorrect"); # load data
SeisPlot(seis_incorrect[:,:,50], cmap = "seismic") # plot data
Out[11]:
In [13]:
# Therefore, if we tried to interpret on our bad velocity
# model, we would never find our drilling target. Hence
# we need to improve our velocity model somehow. So let's
# use full waveform inversion.
# Let's input our poorly picked velocity model and the
# seismic data into an FWI algorithm and take a look
# at the updated velocity model to see if there was any
# improvement:
# NOTE: RUNNING FWI MAY TAKE 10 MINUTES OR MORE...
tic() # start timer
fwi("../dat/vel_incorrect", "../dat/image_correct",
"../dat/wav", "../dat/updated_vel", .03, 10, 0); # FWI algorithm
toc() # end timer + print
vel_updated, vel_h = SeisRead("../dat/updated_vel") # load data
SeisPlot(vel_updated[:,:,50], cmap = "viridis_r") # plot data
Out[13]:
In [12]:
# We can see here that the vel model has been improved considerably
# by full waveform inversion. Though the edges of the reservoir
# are smeared (because the wavelet envelope comparison algorithm
# isn't sophisticated enough to handle tuning), you are clearly
# able to delinate where the reservoir is.
# For comparison, let's view the insitu velocities.
# (This is the model the author used to create the "real" seismic)
vel_true, vel_h = SeisRead("../dat/vel_correct"); # read data
SeisPlot(vel_true[:,:,50], cmap = "viridis_r") # plot data
Out[12]:
In [7]:
# As you can see, the FWI updated vel model is much closer
# to the true velocity model.
In [ ]: