notebook.community

Edit and run 





In [1]:



    


%pylab inline



    


















    






Populating the interactive namespace from numpy and matplotlib

















In [2]:



    


from pkprocess import *



    











In [3]:



    


sd=read_su('data/marmousi.su') # Seismic Data



    











In [4]:



    


get_ns(sd),get_ntr(sd),get_dt(sd)



    


















    
Out[4]:








(723, 23040, 0.004)

















In [5]:



    


print_range(sd)



    


















    






keyword:: min(loc) ~ max(loc): [first - last]// # of unique keys

tracl:: 1(0) ~ 23040(23039): [1 - 23040]// 23040
tracr:: 1(0) ~ 23040(23039): [1 - 23040]// 23040
fldr:: 3000(0) ~ 8975(22944): [3000 - 8975]// 240
tracf:: 1(0) ~ 96(95): [1 - 96]// 96
cdp:: 1(0) ~ 574(23039): [1 - 574]// 574
cdpt:: 1(0) ~ 96(95): [1 - 96]// 96
trid:: 1(0) ~ 1(0): [1 - 1]// 1
duse:: 1(0) ~ 1(0): [1 - 1]// 1
offset:: -2575(0) ~ -200(95): [-2575 - -200]// 96
sx:: 3000(0) ~ 8975(22944): [3000 - 8975]// 240
gx:: 425(0) ~ 8775(23039): [425 - 8775]// 335
ns:: 723(0) ~ 723(0): [723 - 723]// 1
dt:: 4000(0) ~ 4000(0): [4000 - 4000]// 1

Header information





In [6]:



    


plot(ntr_per_shot(sd),'o')
xlabel('Shot gather numbers')
ylabel('Number of traces/shot gather')
grid()



    


















    
























In [7]:



    


plot(get_key(sd,'sx'))
xlabel('Trace number')
ylabel('Sources x-axis locations (m)')



    


















    
Out[7]:








Text(0,0.5,'Sources x-axis locations (m)')










    
























In [10]:



    


print(get_key(sd,'sx'))
print(get_key(sd,'gx'))



    


















    






[3000 3000 3000 ... 8975 8975 8975]
[ 425  450  475 ... 8725 8750 8775]

















In [11]:



    


plot(get_key(sd,'sdepth'))
xlabel('Trace number')
ylabel('Sources depth (m)')



    


















    
Out[11]:








Text(0,0.5,'Sources depth (m)')










    
























In [12]:



    


plot(get_key(sd,'gelev'))
xlabel('Trace number')
ylabel('Receivers depth (m)')



    


















    
Out[12]:








Text(0,0.5,'Receivers depth (m)')










    
























In [14]:



    


stacking_chart(sd)
savefig('stacking_chart.png')

Displaying seismic data





In [15]:



    


get_key_unique(sd,'fldr')



    


















    
Out[15]:








array([3000, 3025, 3050, 3075, 3100, 3125, 3150, 3175, 3200, 3225, 3250,
       3275, 3300, 3325, 3350, 3375, 3400, 3425, 3450, 3475, 3500, 3525,
       3550, 3575, 3600, 3625, 3650, 3675, 3700, 3725, 3750, 3775, 3800,
       3825, 3850, 3875, 3900, 3925, 3950, 3975, 4000, 4025, 4050, 4075,
       4100, 4125, 4150, 4175, 4200, 4225, 4250, 4275, 4300, 4325, 4350,
       4375, 4400, 4425, 4450, 4475, 4500, 4525, 4550, 4575, 4600, 4625,
       4650, 4675, 4700, 4725, 4750, 4775, 4800, 4825, 4850, 4875, 4900,
       4925, 4950, 4975, 5000, 5025, 5050, 5075, 5100, 5125, 5150, 5175,
       5200, 5225, 5250, 5275, 5300, 5325, 5350, 5375, 5400, 5425, 5450,
       5475, 5500, 5525, 5550, 5575, 5600, 5625, 5650, 5675, 5700, 5725,
       5750, 5775, 5800, 5825, 5850, 5875, 5900, 5925, 5950, 5975, 6000,
       6025, 6050, 6075, 6100, 6125, 6150, 6175, 6200, 6225, 6250, 6275,
       6300, 6325, 6350, 6375, 6400, 6425, 6450, 6475, 6500, 6525, 6550,
       6575, 6600, 6625, 6650, 6675, 6700, 6725, 6750, 6775, 6800, 6825,
       6850, 6875, 6900, 6925, 6950, 6975, 7000, 7025, 7050, 7075, 7100,
       7125, 7150, 7175, 7200, 7225, 7250, 7275, 7300, 7325, 7350, 7375,
       7400, 7425, 7450, 7475, 7500, 7525, 7550, 7575, 7600, 7625, 7650,
       7675, 7700, 7725, 7750, 7775, 7800, 7825, 7850, 7875, 7900, 7925,
       7950, 7975, 8000, 8025, 8050, 8075, 8100, 8125, 8150, 8175, 8200,
       8225, 8250, 8275, 8300, 8325, 8350, 8375, 8400, 8425, 8450, 8475,
       8500, 8525, 8550, 8575, 8600, 8625, 8650, 8675, 8700, 8725, 8750,
       8775, 8800, 8825, 8850, 8875, 8900, 8925, 8950, 8975], dtype=int32)

















In [16]:



    


shot5000=window(sd,'fldr',5000)
plot_wiggle(shot5000,perc=99)



    


















    






min=-562.14453125 max=593.078125










    
























In [17]:



    


plot_image(shot5000,perc=99)



    


















    






min=-562.14453125 max=593.078125










    
























In [18]:



    


plot_image(shot5000,perc=99,cmap='bwr')



    


















    






min=-562.14453125 max=593.078125










    
























In [19]:



    


plot_image(sd,figsize=[15,10],perc=99)



    


















    






min=-719.484375 max=729.921875










    
























In [20]:



    


get_key_unique(sd,'offset')



    


















    
Out[20]:








array([-2575, -2550, -2525, -2500, -2475, -2450, -2425, -2400, -2375,
       -2350, -2325, -2300, -2275, -2250, -2225, -2200, -2175, -2150,
       -2125, -2100, -2075, -2050, -2025, -2000, -1975, -1950, -1925,
       -1900, -1875, -1850, -1825, -1800, -1775, -1750, -1725, -1700,
       -1675, -1650, -1625, -1600, -1575, -1550, -1525, -1500, -1475,
       -1450, -1425, -1400, -1375, -1350, -1325, -1300, -1275, -1250,
       -1225, -1200, -1175, -1150, -1125, -1100, -1075, -1050, -1025,
       -1000,  -975,  -950,  -925,  -900,  -875,  -850,  -825,  -800,
        -775,  -750,  -725,  -700,  -675,  -650,  -625,  -600,  -575,
        -550,  -525,  -500,  -475,  -450,  -425,  -400,  -375,  -350,
        -325,  -300,  -275,  -250,  -225,  -200], dtype=int32)

Nearest-offset traces





In [21]:



    


near=window(sd,'offset',-200)



    











In [22]:



    


plot_wiggle(near,figsize=[10,10],perc=99)



    


















    






min=-709.8984375 max=728.31640625

Gain control





In [23]:



    


shot5000_tpow=gain(shot5000,tpow=1)
shot5000_epow=gain(shot5000,epow=1)
plot_comp((shot5000_tpow,shot5000_epow),plot='wiggle',perc=99)



    


















    






min=-871.3752890625 max=863.01325
min=-3087.258500663747 max=3069.501727809401










    
























In [24]:



    


sd_gained=gain(sd,tpow=1)
sd_gained.print_log()



    


















    






read_su: data/marmousi.su
gain: tpow=1 epow=0 agc=False agc_gate=0.5 norm=rms

















In [25]:



    


seis_env_dB(shot5000,shot5000_tpow)

FX, FK Spectrum





In [26]:



    


specfx(shot5000_tpow)



    


















    






dt=0.004, fmax=125.0










    
























In [27]:



    


specfk(shot5000_tpow)



    


















    






dt=0.004, fmax=125.0
dx=0.025, kmax=20.0

Autocorrelograms





In [28]:



    


sd3=window(sd_gained,'fldr',[5000, 5025, 5050])
max_lag=0.2
auto=auto_correlation_map(sd3,max_lag)
plot_wiggle(auto,figsize=[10,10],perc=99)
ylabel('Time lag (s)')



    


















    






min=-41984683.13342782 max=54597445.083550796










    
Out[28]:








Text(0,0.5,'Time lag (s)')










    
























In [29]:



    


plot_wiggle(sd3,figsize=[15,10],perc=99)



    


















    






min=-892.1047500000001 max=890.817890625

Spiking deconvolution





In [30]:



    


mu=0.1
sd_decon=spiking_decon(sd_gained,max_lag,mu)



    


















    






Calculating the filter
Applying the filter

















In [32]:



    


#%timeit -r 1 -n 1 spiking_decon(sd_gained,max_lag,mu)



    











In [33]:



    


sd3=window(sd_decon,'fldr',[5000, 5025, 5050])
plot_wiggle(sd3,figsize=[15,10],perc=99)



    


















    






min=-3.2932219421432025e-09 max=3.3076127049752527e-09

Automatic Gain Control





In [34]:



    


sd_agc=gain(sd_decon,agc=True,agc_gate=0.5,norm='rms')



    











In [35]:



    


plot_wiggle(window(sd_agc,'fldr',[5000,5025,5050]),figsize=[15,10],perc=99)



    


















    






min=-3.1169822500678457 max=3.179162786660583










    
























In [36]:



    


sd_agc.print_log()



    


















    






read_su: data/marmousi.su
gain: tpow=1 epow=0 agc=False agc_gate=0.5 norm=rms
spiking_decon: max_lag=0.2 mu=0.1
gain: tpow=0 epow=0 agc=True agc_gate=0.5 norm=rms

CDP Sorting





In [37]:



    


sd_sort=trace_sort(sd_agc,('+cdp','-offset'))



    











In [38]:



    


cmps=get_key(sd_sort,'cdp')
cmpu=get_key_unique(sd_sort,'cdp')
fold_num = [sum(icmp==cmps) for icmp in cmpu]

plot(cmpu,fold_num)
xlabel('CMP numbers')
ylabel('Fold')



    


















    
Out[38]:








Text(0,0.5,'Fold')










    
























In [39]:



    


sd_sort.print_log()



    


















    






read_su: data/marmousi.su
gain: tpow=1 epow=0 agc=False agc_gate=0.5 norm=rms
spiking_decon: max_lag=0.2 mu=0.1
gain: tpow=0 epow=0 agc=True agc_gate=0.5 norm=rms
sort: keys=+cdp,-offset

Output for the velocity analysis





In [40]:



    


sd_sort.write("marm_gain_decon_agc_sort.sd")

Velocity analysis

  1. Open and edit Par_velan.py file.
  2. Run Run_velan.py. ($ python Run_velan.py)
    • Left click: Add a point
    • Right click: Remove a point




In [41]:



    


from IPython.display import YouTubeVideo
display(YouTubeVideo("zWt5lbWEx3g",width=775,height=542))

After the velocity analysis





In [42]:



    


sd_nmo=read("marm_gain_decon_agc_sort_nmo.sd")



    











In [43]:



    


sd_nmo.print_log(nmo=True)



    


















    






read_su: data/marmousi.su
gain: tpow=1 epow=0 agc=False agc_gate=0.5 norm=rms
spiking_decon: max_lag=0.2 mu=0.1
gain: tpow=0 epow=0 agc=True agc_gate=0.5 norm=rms
sort: keys=+cdp,-offset
write: marm_gain_decon_agc_sort.sd
read: marm_gain_decon_agc_sort.sd
velan_nmo: vmin=1500.0 dv=100 nv=41 cmp_start=25 cmp_end=525 cmp_step=50 max_stretch=10.0 R=4 L=8
write: marm_gain_decon_agc_sort_nmo.sd
read: marm_gain_decon_agc_sort_nmo.sd
nmo picks
{25: (array([1722.3655914 , 1868.60215054, 2117.20430108, 2512.04301075,
       3257.84946237]), array([0.62836364, 1.34649351, 1.78036364, 2.24041558, 2.69298701])), 75: (array([1722.3655914 , 1839.35483871, 2146.4516129 , 2570.53763441,
       2994.62365591, 3798.92473118]), array([0.53485714, 1.22680519, 1.63075325, 2.23293506, 2.61444156,
       2.70420779])), 125: (array([1736.98924731, 1853.97849462, 2234.19354839, 2409.67741935,
       2541.29032258, 2716.77419355, 3067.74193548, 3696.55913978,
       4515.48387097]), array([0.50493506, 1.05475325, 1.51106494, 1.75792208, 1.9598961 ,
       2.14690909, 2.42742857, 2.58077922, 2.81641558])), 175: (array([1663.87096774, 2058.70967742, 2117.20430108, 2541.29032258,
       3009.24731183, 4983.44086022]), array([0.25433766, 0.76675325, 1.17444156, 1.59709091, 2.07584416,
       2.43490909])), 225: (array([2000.21505376, 2877.6344086 , 3067.74193548, 3608.8172043 ]), array([0.864     , 1.68685714, 2.07584416, 2.16187013])), 275: (array([1824.7311828 , 2292.68817204, 2438.92473118, 2687.52688172,
       2804.51612903, 4603.22580645]), array([0.49371429, 1.16322078, 1.70181818, 2.09828571, 2.3601039 ,
       2.72664935])), 325: (array([1751.61290323, 1941.72043011, 2336.55913978, 2497.41935484,
       2994.62365591, 4720.21505376]), array([0.16831169, 0.79667532, 1.27542857, 1.98607792, 2.44987013,
       2.84633766])), 375: (array([1722.3655914 , 1853.97849462, 2380.43010753, 2424.30107527,
       2877.6344086 , 5275.91397849]), array([0.26555844, 0.62462338, 1.57090909, 1.9598961 , 2.46109091,
       2.83511688])), 425: (array([1693.11827957, 2117.20430108, 2541.29032258, 2658.27956989,
       3594.19354839]), array([0.61714286, 1.31657143, 1.56716883, 1.90005195, 2.10202597])), 475: (array([1707.74193548, 2087.95698925, 2468.17204301, 2775.2688172 ,
       3462.58064516, 4837.20430108]), array([0.5198961 , 1.26420779, 1.66067532, 1.85890909, 2.13568831,
       2.73038961])), 525: (array([1751.61290323, 2248.8172043 , 2877.6344086 , 4939.56989247]), array([0.44883117, 1.28664935, 1.70555844, 2.54711688]))}

















In [44]:



    


plot_image(sd_nmo,figsize=[15,10],perc=99)



    


















    






min=-2.7117873033959325 max=2.7374865904961885

Stacking





In [45]:



    


sd_stk=stack(sd_nmo)



    











In [46]:



    


plot_image(sd_stk,figsize=[15,10],perc=99,key='cdp')



    


















    






min=-52.22541936884288 max=52.232315642273306










    
























In [47]:



    


print_range(sd_stk)



    


















    






keyword:: min(loc) ~ max(loc): [first - last]// # of unique keys

cdp:: 1(0) ~ 574(573): [1 - 574]// 574
ns:: 723(0) ~ 723(0): [723 - 723]// 1
dt:: 4000(0) ~ 4000(0): [4000 - 4000]// 1
shortpad:: 1(0) ~ 48(94): [1 - 1]// 48

















In [48]:



    


get_key(sd_stk,'shortpad') # shortpad keyword contains the fold number



    


















    
Out[48]:








array([ 1,  1,  2,  2,  3,  3,  4,  4,  5,  5,  6,  6,  7,  7,  8,  8,  9,
        9, 10, 10, 11, 11, 12, 12, 13, 13, 14, 14, 15, 15, 16, 16, 17, 17,
       18, 18, 19, 19, 20, 20, 21, 21, 22, 22, 23, 23, 24, 24, 25, 25, 26,
       26, 27, 27, 28, 28, 29, 29, 30, 30, 31, 31, 32, 32, 33, 33, 34, 34,
       35, 35, 36, 36, 37, 37, 38, 38, 39, 39, 40, 40, 41, 41, 42, 42, 43,
       43, 44, 44, 45, 45, 46, 46, 47, 47, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48, 48,
       48, 48, 48, 48, 47, 47, 46, 46, 45, 45, 44, 44, 43, 43, 42, 42, 41,
       41, 40, 40, 39, 39, 38, 38, 37, 37, 36, 36, 35, 35, 34, 34, 33, 33,
       32, 32, 31, 31, 30, 30, 29, 29, 28, 28, 27, 27, 26, 26, 25, 25, 24,
       24, 23, 23, 22, 22, 21, 21, 20, 20, 19, 19, 18, 18, 17, 17, 16, 16,
       15, 15, 14, 14, 13, 13, 12, 12, 11, 11, 10, 10,  9,  9,  8,  8,  7,
        7,  6,  6,  5,  5,  4,  4,  3,  3,  2,  2,  1,  1], dtype=int16)

















In [49]:



    


plot_image(window(sd_stk,'shortpad',48),figsize=[15,10],perc=99,key='cdp')



    


















    






min=-52.11353700118288 max=51.574066660876326

Stolt migration





In [50]:



    


v=3000.
dx=12.5
sd_mig=stolt_mig(sd_stk,v,dx)



    











In [51]:



    


plot_image(sd_mig,figsize=[15,10],perc=99,key='cdp')



    


















    






min=-39.275908898720076 max=39.69829313571435










    
























In [52]:



    


plot_image(window(sd_mig,'shortpad',48),figsize=[15,10],perc=99,key='cdp')



    


















    






min=-41.265648180822005 max=41.91458695484362










    
























In [53]:



    


plot_comp((sd_stk,sd_mig),figsize=[12,10],perc=99,key='cdp')



    


















    






min=-52.22541936884288 max=52.232315642273306
min=-39.275908898720076 max=39.69829313571435










    
























In [54]:



    


sd_mig.print_log()



    


















    






read_su: data/marmousi.su
gain: tpow=1 epow=0 agc=False agc_gate=0.5 norm=rms
spiking_decon: max_lag=0.2 mu=0.1
gain: tpow=0 epow=0 agc=True agc_gate=0.5 norm=rms
sort: keys=+cdp,-offset
write: marm_gain_decon_agc_sort.sd
read: marm_gain_decon_agc_sort.sd
velan_nmo: vmin=1500.0 dv=100 nv=41 cmp_start=25 cmp_end=525 cmp_step=50 max_stretch=10.0 R=4 L=8
write: marm_gain_decon_agc_sort_nmo.sd
read: marm_gain_decon_agc_sort_nmo.sd
stack
stold_mig: v=3000.0 dx=12.5

















In [51]:



    


#sd_mig.write('mig.sd')
#sd_mig.write_su('mig.su')

Content source: pkgpl/IPythonProcessing

Similar notebooks:

notebook.community | gallery | about