All examples for ACOMP
# Complete job example for Q-compensation of prestack data
*call debug on # enable debugging messages from the code
## This info is printed on the screen when the job is executed:
*info
Demo of QCPost
A-compensation with Panimaging data. The attenuation-compensation
is performed by tool 'acomp' (see "*call acomp" below). All options
available from 'acomp' are listed. To try the various options,
you need to uncomment them and comment those which are not needed.
The procedure below consists of three steps:
1) Loading the stacking velocity model from ASCII file Plain/vel.txt,
transforming it in a series of traces and storing in a buffer
(see '*call store').
2) Deriving several other models (interval velocity, 1/Q)
from this stacking velocity. The stacking velocity is interpreted
as RMS velocity.
3) While loading the seismic data, performing the required processing,
and plotting in a PostScript file.
###############################################################################
*call trcinfo cdp={cdp} chan={chan} # this is just creation of a useful
# trace label that can be handy in error output
###############################################################################
## Step 1: Load the stacking-velocity model and
## store it as seismic traces in a memory buffer.
###############################################################################
#### input stacking velocity from Vista database file:
*call rdascii 10 400 Plain/vel.txt
vista vel2d cdp cdp_x cdp_y
## store copies of RMS velocities in memory buffer
*call store model_traces
*end keep # the traces are blocked here from passing futher
# in the job
############################################in#################################
## Step 2: Form the required models
## based on stored stacking-velocity traces
###############################################################################
## tell the program to treat the stored traces as a RMS velocity model:
*call modeltr vrms_model
source 1 cdp model_traces
## Derive interval velocity model from RMS velocity:
*call modeltr vint_model
source 1 cdp model_traces
transf vrms->vint
## Derive the reference velocity model from interval velocity.
## The reference velocity is the velocity Vr used in modeling dispersion.
*call modeltr vref_model
source 1 cdp model_traces
transf vrms->vint # derive interval velocity from RMS velocity
# and use it as reference
## Derive a model for 1/Q by scaling interval velocity:
*call modeltr q_model
source 1 cdp model_traces
#transf vrms->vint
#transf vint->1/Q 0.10 1.5 1.0
transf const 0.010
## Construct a model for eta (frequency dependence of 1/Q)
## as a constant = 0.0 for this demo:
*call modeltr eta_model
source 1 cdp model_traces
transf const 0.0 # set eta = 0
## Construct a model for gamma (local geometric spreading)
## as a constant (for this demo):
*call modeltr gamma_model
source 1 cdp model_traces
transf const 0.02 # set fixed gamma = 0.02
## Construct a model for scattering (RMS reflection amplitude)
## as a constant for this demo:
*call modeltr refl_model
source 1 cdp model_traces
transf const 0.05 # set fixed reflectivity = 0.05
###############################################################################
## Step 3: Start over with trace data input and processing
###############################################################################
## load data traces from one or several SEGY file(s):
*call load Plain/brute_stack_PanImaging.sgy
## resample to 1 ms
*call modify 1.0
## Trace header list for quality control and information:
*call table *
cdp offset chan sou_x sou_y selev rec_x rec_y relev
## Attenuation compensation:
*call acomp stacked
# uncomment this line to load all parameters from file 'params.acomp' instead of
# entering them futher in the job:
#loadpar params.acomp
# file name in which all parameters below are saved
savepar params.acomp
# Optional files in which the Q and phase-velocity curves are output for info.
# These files can be displayed in Matlab, Octave, gnuplot, GMT, input in Excel, etc:
print qn(f) 200 300 DB/qn_f # at t0 = 200ms, 300 points sampled across the frequency band
print v(f) 200 300 DB/v_f # note that file names are without extension ".txt"
print vn(f) 200 300 DB/vn_f
print wp(t) 200 200 4 100 DB/wp_t # time-domain propagation responses for [-100,100] ms
# in model at t0 = 200 ms,
# modeled at at times = 0, 1000, 2000, 3000, and 4000 ms
print w(t) 200 200 4 100 DB/w_t # time-domain propagation responses for [-100,100] ms
# Specify trace headers used to find trace locations in models:
layout cdp espnum # headers used
# Specify trace headers used to find trace locations in models,
# header containing offsets, statics, etc.:
geom offset sou_x sou_y selev sdepth rec_x rec_y relev
# Selection of a range of CDPs and time range for processing:
seltrc cdp 200 220
seltime 200 1500
# Specify the models to use:
m-vrms 2.4 # RMS velocity model
m-vint vint_model # Interval velocity model
m-vref vref_model # Reference velocity model for dispersion
m-q q_model # Model of 1/Q
m-rrms refl_model # model for RMS reflection amplitude (used for modeling scattering)
# Source wavelet selections (the first one uncommented is used):
#wavelet berlage zero 30 60 0.5 # 30-ms Berlage wavelet
#wavelet gabor zero 30 90 # 30-ms Gabor wavelet with 90-Hz modulaion
#wavelet gauss zero 30 # 30-ms Gaussian wavelet
#wavelet ormsby zero 10 20 40 60 # Ormsby wavelet
wavelet ricker zero 30 # zero-phase 30-Hz Ricker wavelet
#wavelet klauder zero 10 60 8 # Klauder wavelet
#wavelet saw zero 30 # 30-ms sawtooth wavelet
#wavelet spike zero 30 # spike (single trace sample) wavelet.
# for spike, we still have to declare
# wavelet breadth (= 30 ms here) in order
# to determine the multi-tapering window
# Reference time by which the attenuation and geometric spreading
# are normalized:
reftime 100 # 100 ms
# Geometric-spreading options
# Uncomment one or several of these options.
# All geometric-spreading corrections are superimposed
#g-vt # Newman's geometric spreading
#g-t 0.2 # Additional weak time-dependent geometric spreading
#g-texp 2.0 # time-exponential spreading in dB/sec
#g-r 0.2 # power-law spreading based on ray lengths
#g-gamma gamma-model # geometric correction using local values of gamma
# attenuation modeling options. You can use multiple entries of each type:
#a-q ampl qonly 0.03 50 2.4 # using Q model,
# amplitude-only application
a-q full constq 0.003 50 2.4 # Constant-Q (Kjartansson's) model
#a-q full kolsky q_model 50 vref_model # using Kolsky's model
#a-q full azimi q_model 50 vref_model # using Azimi's model
#a-q full powlaw 0.1 q_model 50 vref_model # using power-law Q (Muller's) model
#a-q full dv(f) q_model 50 vref_model # using Q models
#10 -0.02 # for this type of attenuation, velocity dispersion is given as dv(f) here
#50 0 # dispersion=0 at reference frequency = 50 Hz
#200 0.05
#a-visc # using solid viscosity models
#a-poro # using poroelastic model (not implemeted yet)
#a-scatt refl_model # using scatering
# deconvolution parameters:
band 10 20 80 160 # frequency band for deconvolution
#decon zero 0.0 wiener 0.001 nosrc # Wiener deconvolution
#decon minph 0.0 wlevel 0.001 decsrc # Water-level deconvolution
#decon model iter 0.001 0.3 decsrc # Iterative time-domain deconvolution
# When using iterative deconvolution, we usually need to convolve with Gaussian wavelet.
# Note that this option can also replace the selection of the frequency band (list 'band')
#convol gauss zero 30 # 30-ms Gaussian wavelet
# The following are several options for output.
# Uncomment the required option (the first one uncommented is used).
#output data # output original data
#output tapers # output the tapers instead of data
#output tapdata # output the data records tapered in processing windows
output wavelet align # output the wavelet instead of data
#output dwavelet align # output the wavelet instead of data
#output decon # output deconvolved records
## plot by using program psxy from GMT (https://www.soest.hawaii.edu/gmt/)
###############################################################################
## select a range of CDPs for plotting (this is not really required)
*call edit pass
range cdp 200 220
## text substitutions:
*define
psxydir /tmp/psxy/ # temporary directory for GMT inputs
psfile PS/Plain_stack.ps # output PostScript file name with the trace plot
psdisp PS/disprel.ps # output PostScript file with dispersion relations
pswt PS/wt.ps # output PostScript file with attenuating waveforms
## plot:
*call plot
gain 0.8
setamp peak 0 1200 # trace scaling for plotting
psxy {psxydir}/traces
offset cdp # horizontal axis for the plot (CDP, i.e. CMP)
trange 0 2500
# this tool runs a shell to display PostScript
*call unix eof
#sh plot.sh {psxydir} test # a couple of options for trace plotting
#sh plot_200ms.sh {psxydir} test &
sh plot_disprel.sh DB/vn_f DB/qn_f {psdisp} & # plot dispersion relations
sh plot_wt.sh DB/wp_t DB/w_t {pswt} & # plot attenuating spikes and waveforms
set-env GMTREGION_wt {gmt.region()}