
Publication History:
This article is based on
"Crain's Seismic Petrophysics" by E. R. (Ross) Crain, P.Eng., first
published in 2003, and updated annually until 2016.
This
webpage version is the copyrighted intellectual
property of the author.
Do not copy or distribute in any form without explicit
permission. 
Calculating Seismic
Wavelets
If a wavelet can be extracted by autocorrelation of a real seismic
trace, it should be used to make the synthetic. If this cannot
be done, wavelets are generated from equations which describe
the frequency content of the wavelet. There are two types  Ricker
wavelets are generated directly in the time domain and Klauder
wavelets are generated in the frequency domain. Both are called
zero phase wavelets and are a minimum delay. That is, the maximum
energy is at the beginning of the wavelet. Ricker wavelets are
shown below.
Standard computed wavelet shapes
There
are four different Ricker wavelets, called Ricker Type 1, Type
2, Type 3, and Type 4 respectively. The higher types have more
peaks and valleys, representing more "ringing", and
are less symmetrical.
Ricker wavelets are used when
the synthetic seismogram is to match dynamite or airgun data.
They are computed from the following:
Ricker
Type 1:
1: G = 3.4641 * CFREQ / 1000
2: FOR I = 1 TO WVLTLENGTH / SAMPRATE
3: X = (I * SAMPRATE  1) * G  3.75
4: Wvlt [I] = (X * X  1) / EXP (X * X / 2)
5: NEXT I
Ricker
Type 2:
1: G = 3.65 * CFREQ / 1000
2: FOR I = 1 TO WVLTLENGTH / SAMPRATE
3: X = (I * SAMPRATE  1) * G  3.85
4: Wvlt [I] = (X * X  X  1) / EXP (X * X / 2)
5: NEXT I
Ricker
Type 3:
1: G = 2.9 * CFREQ / 1000
2: FOR I = 1 TO WVLTLENGTH / SAMPRATE
3: X = (I * SAMPRATE  1) * G  3.85
4: Wvlt [I] = (X^3 + 3 * X) / EXP (X * X / 2)
5: NEXT I
Ricker
Type 4:
1: G = 2.35 * CFREQ / 1000
2: FOR I = 1 TO WVLTLENGTH / SAMPRATE
3: X = (I * SAMPRATE  1) * G  3.2
4: Wvlt [I] = 0.1*(X^6+X^5+15*X^410*X^345*X*X+15*X+15)/EXP(X*X/2)
5: NEXT I
When
comparing synthetic seismograms with different wavelets, it is
important that each wavelet have the same initial energy. The
energy can be found by squaring and summing each of the Wvlt[I]
terms. The sum is then compared to a standard energy of 1.0 and
a scale factor derived, which is multiplied against each Wvlt[I]
term to obtain normalized wavelet values.
The
Klauder wavelet is generated in the frequency domain, which gives
the amplitude and phase spectra of the wavelet. The autocorrelation
of the spectra produces the wavelet shape in the time domain.
This wavelet is used when the synthetic is to match Vibroseis
data.
Klauder
Wavelet:
STEP
1: Generate chirp signal amplitude and phase
1: CFREQ = 0.5 * (LFREQ + HFREQ)
2: KLFREQ = (HFREQ  LFREQ) / SWEEPLENGTH
3: K = SAMPRATE / 1000
4: FOR I = 1 TO SWEEPLENGTH / K
Calculate
phase
5: I[I] = 2 * PI * (LFREQ * I * K + KLFREQ * I * I * K * K / 2)
+ INITPHASE
Calculate
amplitude
6: R[I] = 1250 * COS(I[I])
Put
phase in all 4 quadrants
7: I[I] = TAN (I[I])
8: I[I] = ATN (I[I])
9: I[I] = I[I] * 180 / PI + 180 * (SIN (I[I]) < 0) + 90
10: NEXT I
STEP
2: Calculate autocorrelation to get wavelet in time domain. This
double nested loop takes the first 100 data points and cross correlates
them with the last 100 data points of the sweep. Since the number
of data points is small, there is little advantage to using the
FFT for the correlation.
1: I = SWEEPLENGTH / SAMPRATE  1
2: K = 1
3: FOR J = SWEEPLENGTH/SAMPRATE TO SWEEPLENGTH/SAMPRATE  100
BY 1
4: Auto[K] = 0
5: FOR L = 1 TO J
6: Auto[K] = R[L] * R[I  J + L] + Auto[K]
7: NEXT L
8: K = K + 1
9: NEXT J
To
see the wavelet shape, the correlation string is reversed and
paired with itself:
10: FOR I = 1 to 100
11: Wvlt[I] = Auto[100  I + 1]
12: NEXT I
13: FOR I = 101 TO 199
14: Wvlt[I] = Auto[I  99]
15: NEXT I
