Calculate TD shot gather from FD modelling results
Daniel Köhn
Kiel, 01.09.2017
Import necessary packages
In [ ]:
import numpy as np
import matplotlib as plt
import matplotlib.pyplot as plt
from matplotlib.pyplot import gca
import matplotlib as mpl
from pylab import rcParams
from matplotlib import rc
from scipy import interpolate
from scipy import signal
#%matplotlib notebook
Define parameters for data visualization
In [ ]:
clip = 5e-3 # data clip
nsrc = 1 # number of shots
ntr = 117 # number of receivers
shotno = 1 # extract FDFD data for shot shotno
# number of frequencies nf and frequency interval df
FC_low = 1e6
FC_high = 800e6
nf = 400
df = (FC_high - FC_low) / (nf-1)
nfnsrc = nf * nsrc
# TD parameters
TmaxTD = 150e-9 # maximum time TD
TmaxFD = 0.25 / df # maximum time FD
dt = TmaxFD / nf # time sampling FD
# maximum time sample of FD2TD corresponding to TmaxTD
nmaxFD = np.int(TmaxTD / dt)
# define time samples in time and frequency domain
tFD = np.arange(0.0, TmaxTD, dt)
offset = np.arange(0, ntr, 1)
Font properties
In [ ]:
FSize = 25
font = {'color': 'black',
'weight': 'bold',
'size': FSize}
mpl.rc('xtick', labelsize=FSize)
mpl.rc('ytick', labelsize=FSize)
rcParams['figure.figsize'] = 10, 10
Time Damping
In [ ]:
S = 0.0
Read 2 cross TE-mode FD data
In [ ]:
name_FDFD = "../seis/2_cross_TE_p_stage_1.bin"
f = open (name_FDFD)
data_type = np.dtype ('float32').newbyteorder ('<')
data_FDFD = np.fromfile(f, dtype=data_type)
data_FDFD = data_FDFD.reshape(nf*nsrc,2*ntr)
data_FDFD = np.transpose(data_FDFD)
data_FDFD = np.flipud(data_FDFD)
Extract FD data for shot shotno
In [ ]:
FD_real = data_FDFD[0:2*ntr-1:2,shotno-1:nfnsrc:nsrc]
FD_imag = data_FDFD[1:2*ntr:2,shotno-1:nfnsrc:nsrc]
data_FDFD = None
# taper FD data
window = signal.tukey(nf,1)
for i in range (ntr):
FD_real[i,:] *= window
FD_imag[i,:] *= window
In [ ]:
# assemble complex data
FDFD = FD_real + 1j*FD_imag
FD_real = None
FD_imag = None
FDFD = np.transpose(FDFD)
tmp1 = np.concatenate((FDFD, np.zeros((nf, ntr)), np.zeros((nf, ntr)), np.flipud(FDFD)))
#tmp1 = np.concatenate((np.zeros((1, ntr)),FDFD,np.flipud(FDFD)))
FDFD = None
FDFD = tmp1
In [ ]:
#plt.imshow(np.real(FDFD))
#plt.show()
IFFT of FD data
In [ ]:
# transformation from FD to TD
FD2TD = np.fft.ifft(FDFD,axis=0)
FD2TD = np.real(FD2TD)
tmp1 = None
Extract FD2TD up to TmaxTD and resample FD2TD from dtFD to dt
In [ ]:
# extract FD2TD data up to TmaxTD
tmp1 = FD2TD[1:nmaxFD,:]
FD2DTD = None
# clean memory
tmp = None
In [ ]:
#plt.imshow(FD2TD, vmin=-clip, vmax=clip)
#plt.show()
In [ ]:
#print(np.shape(tFD), np.shape(offset), np.shape(tmp1))
Damp Time Domain Data
In [ ]:
#tmp2 = np.exp(-S * tTD)
#tdamp = np.tile(tmp2,(ntr,1))
#tdamp = np.flipud(tdamp.T)
#traces = traces * tdamp
#tmp2 = None
#tdamp = None
Normalize data
In [ ]:
#max_TD = np.max(traces)
max_FD = np.max(tmp1)
tmp1 = tmp1 / max_FD
FD2TD = np.flipud(tmp1)
tmp1 = None
Shuffle 2-cross TE model data
In [ ]:
tmp1 = FD2TD[:,0:29]
tmp2 = FD2TD[:,30:57]
tmp3 = FD2TD[:,58:87]
tmp4 = FD2TD[:,88:116]
FD2TD = None
tmp = np.concatenate((tmp1.T, tmp4.T, np.flipud(tmp2.T), np.flipud(tmp3.T)))
FD2TD = tmp.T
tmp = None
tmp1 = None
tmp2 = None
tmp3 = None
tmp4 = None
#plt.imshow(FD2TD)
#plt.show()
Plot shot gather
In [ ]:
def do_plot(n, model, cm, an, title, clip):
extent = [np.min(offset), np.max(offset), 0.0, dt * nmaxFD * 1e9]
plt.rc('font',**{'family':'sans-serif','sans-serif':['Helvetica']})
plt.rc('text', usetex=True)
im = plt.imshow(model, cmap=cm, interpolation='nearest', extent=extent, vmin=-clip, vmax=clip, aspect=0.5)
a = gca()
a.set_xticklabels(a.get_xticks(), font)
a.set_yticklabels(a.get_yticks(), font)
plt.title(title, fontdict=font)
if n==1:
plt.ylabel('Time [ns]', fontdict=font)
plt.xlabel('Trace no.', fontdict=font)
plt.gca().invert_yaxis()
In [ ]:
plt.close('all')
plt.figure()
do_plot(1, FD2TD, 'gray', '(b)', r"2-cross TE-mode model (GERMAINE)", clip)
plt.tight_layout()
output_file = "test.pdf"
plt.savefig(output_file, bbox_inches='tight', format='pdf')
plt.show()
In [ ]: