Q factor estimation for seismic data with localized phase space method

Report
Q factor estimation for seismic data with
localized phase space method
Jinghuai Gao
Institute of Wave and Information
Xi’an Jiaotong University, China
Outline
1. Why we need to estimate the Q ?
2. Methods and challenges
3. WEPIF principle
4. WEPIF methods and applications
5. Conclusions
2
Why we need to estimate the Q ?
Seismic waves propagating through the earth
suffer attenuation and dispersion due to the
viscosity of the media (Ricker, 1953; Futterman,
1962; Kneib and Shapiro, 1995). Seismic
attenuation is usually quantified by the quality
factor Q.
3
Measurements in both the field and the laboratory demonstrate that Q
is associated with factors such as physical properties of rocks, types of
fluids, and saturation of fluids (Winkler and Nur, 1982; Sheriff and
Geldart, 1995). Thus,

Q is a diagnostic tool for reservoir characterization and hydrocarbon
detection (Toksöz et al., 1979; Frisillo and Stewart, 1980).

Q is very important when interpreting the effects of AVO, improving
the resolution of seismic imaging, and advancing the study of material
properties.
Methods and difficulties

Approach 1:Time-domain methods: pulse
amplitude decay, pulse rising time, pulse
broadening.
Problem: amplitude information of seismic pulses is often
influenced by scattering, geometric spreading, and other
factors.
5
Methods and difficulties

Approach 2: Fourier-frequency domain
methods: logarithm spectral ratio (LSR), centroid
frequency shift (CFS), peak frequency shift
methods.
Problem:To properly select the window function
and window length (Harris, 1978) is difficulty.

Approach 3:Instantaneous frequency (IF) methods

Tonn (1991), Barnes (1991 and 1993), and Engelhard (1996) obtained
the relationship between the measured instantaneous spectra and
seismic attenuation.
Assuming that the source wavelet is an idealized band-pass wavelet,
Barnes (1991) derived the relation between Q and IF variations of
seismic waves, establishing a new approach for Q estimation.

Problem: These reletionships are qualitative or implicit and cannot
be used to estimate the Q.
In this work, we further develop Barnes’s work
(Barnes, 1991), and propose a method for
estimating Q based on our new development.
WEPIF-method principle
Definition: Instantaneous frequency at the wavelet envelope peak is
called WEPIF
80
source
60
40
Receiver 1
0
Receiver 1
20
0
t=0
80
-20
-40
z
-60
10
20
Receiver 2
Z
30
40
50
60
70
80

40
z
P(c, , z, Q)
60
0
20
Receiver 2
0
t 
-20
A slab
-40
-60
0
10
20
30
40
50
60
70
80
9
Propagator
Considering a one-way plane wave propagating in a
homogeneous anelastic medium with a frequencyindependent Q, the propagator is as follows:
 i z
z 
P(c, Q, z, )  exp 


c
(

)
2
c
(

)
Q


where , z denotes the travel-distance, c(ω) the phase
velocity, Q is the quality factor in the medium.
Source wavelet
We assume that the source wavelet and seismic wavelet can be
approximated in the frequency domain by
Uˆ (0, )  A
 
4

2
1/ 4
 (   )2

exp 

i


2 2


(2)
where, σ is the modulating frequency, A and φ are the amplitude
and phase factor, δ is the energy decay factor.
At Z  0 the seismic wavelet recorded by receiver 1 can be
approximated by
Uˆ (0, )  A
 
4

2
1/ 4
 (   )2

exp 

i

.
2

2


(2)
At Z  z the seismic wavelet recorded by receiver 2 is
Uˆ ( z,  )  P(c, Q, z,  )Uˆ (0,  )
2






 z   z 
 4 
 A    exp  

 i    

2
2Qc 
c  
 

1/4
(3)
In order to find the relation between Q and the EPIF, we
first discuss the EPIF of a seismic wavelet. For a constantphase wavelet propagated for time τ= z/c through a
homogeneous anelastic medium, its EPIF is exactly equal
to the average Fourier frequency weighted by its amplitude
spectrum.

f p ( ) 


0
f Uˆ ( , f )df

0
Uˆ ( , f )df
(4)

Using eq.(4 ), we get the EPIF of the source wavelet :
 2 2    
2
exp   2 
 
4 2

2


 

2
f p (0) 


2


0
 2 
 
exp   2  f 

2 
  
2
2
(5)

 df

and the EPIF of the wavelet at τ :
 2 2  
2
  
exp   2 


4 2

2

4 Q 





  

f p ( )  



2
 2 2 

 2 4 Q 
    
 0 exp   2  f  2  4 Q   df


2
(6)
We apply Taylor series expansion to Eq.6, and obtain the
approximate expression ,Eq.7, with its error estimation, based
on analyzing the characteristics of the seismic wavelets which
are extracted from zero-offset VSP data and poststack seismic
data, and investigating the ranges of the Q-factor in various
media.
15
Q
 2
2 

1

D

 2 
4  f p (0)  f p ( )  
1 2
2 2
D
exp

2

 

2
I1 4
I1 

2
(2 )

(7)

2
1
( x) 
2
 t2 
   exp   2  dt
x
(Jinghuai Gao, et al, 2011, Journal of Computational
Acoustic, Vol. 19,No.2)
IF extraction
Conventional method
s(t ) 
 s (t ), sc (t )  s(t )  is (t )
H
*
*
sc is the complex signal of analyzed signal
ds  (t )
ds (t )

s (t )
 s (t )
1
dt
dt
f (t ) 
2
e(t )   em
e(t )  s 2 (t )  s2 (t )
(8)
18
Nosy signal
soft threshold,
hard threshold
Extracting IF of a noisy signal in wavelet domain

Step 1: Choose a proper basic wavelet g(t) and
transform the signal s(t) into the wavelet domain. The
wavelet transform of s(t) is defined as (Mallat, 1999):
1 
 t b 
S  b, a    s(t ) g 
dt

a 
 a 

(9)
Step 2: Determine the range Ω of the energy
distribution of the useful contents of the signal s(t) in
the wavelet domain.

Step 3: Compute the analytic signal sc (t ) of the s(t) in
Ω,let sI(t) represent the imaginary parts of sc (t ) :
1
sc (t ) 
Cg

da
  S t, a  a
(10)
Step 4. Calculate the IF using eq.(8) by replacing s*(t)
with sI (t )  Im(sc (t )).
22
Testing Wavelet-Dependence
Infinite viscoelestic medium, zero-offset VSP record
The shape of four different
seismic wavelet
Receiver
Depth
A slab
The source and the first receiver
are located at the same localization
The first three wavelet have constant-phase,
and the last wavelet is Ricker wavelet.
Wavelet-Dependence
Trace number
The Q estimated by the four methods when the parameter  of the
wavelet is 0.65. The interval between adjacent two traces is 5m.
Wavelet-Dependence
Trace number
The Q estimated by the four methods
when the parameter of the wavelet is 0.5
Wavelet-Dependence
WEPIF
LSR
CFS
PFS
Trace number
The Q estimated by the four methods
when the parameter  of the wavelet is 0.35
26
Wavelet-Dependence testing
tracenumber
The Q estimated by the four methods
when the wavelet is Ricker
27
WEPIF methods and applications-VSP
Data
Figure 1. Forming the slabs for zero-offset VSP data.
28
Calculation procedure
a)
Split the medium into slabs with thicknesses equal to the
space between two successive VSP geophones;
b)
Approximate the direct wave at the top of the ith slab by
a constant-phase wavelet and determine its energy decay
factor δi and modulated frequency σi; then calculate ratio
ηi = σi /2 πδi
c)
Calculate the instantaneous amplitude (IA) of the direct
waves at the top and bottom of the ith slab, respectively,
pick their envelope-peak and determine the traveltime τi
in the ith slab;
d)
e)
f)
g)
Calculate the IF of the direct waves recorded at the slabtop (denoted by ft(t) ) and the slab-bottom (denoted by
fb(t) )
Pick the EPIFs of direct waves with fp(0)= ft(tt) and fp(τi )
= ft(tb), then compute the EPIF variation by fp(0)- fp(τi )
the ith slab.
Insertσi , δi , ηi and fp(0)- fp(τi ) into eq.(7) and compute
the Q-factor for the ith slab.
For all slabs, loop from step (b) to (f) until the Q-factor
of the last slab is obtained.
Examples of Synthetic VSP data
b)
0
density: 2.0 g/cc
c=1.0 km/s
Q=30
Depth (m)
source
layer 1
200
density: 2.5 g/cc
c=1.25 km/s
Q=60
layer 2
400
Depth (m)
a)
density: 3.0 g/cc
c=2.0 km/s
Q=100
layer 3
600
geophones
0
20
40
60
80
100
120
140
160
180
200
220
240
260
280
300
320
340
360
380
400
440
480
520
560
600
0
50
100
150
200
250
300
Time (ms)
Figure 3. 3-layer depth model (a) and synthetic seismograms (b).
Depth (m)
0
20
40
60
80
100
120
140
160
180
200
220
240
260
280
300
320
340
360
380
400
420
440
460
480
500
520
540
560
580
600
-50
T rue Q
LSR
CFS
CFS
WEP IF
LSR
WEPIF
CFS
LSR
WEPIF
-25
0
25
50
75
100
125
150
Q
Figure 4. Theoretical Q-factors (black line). Estimated Q-factors from the
synthetic VSP data using the LSR method (green line), the CFS method
(blue line), and the WEPIF method (red line).
Real zero-offset VSP data with six levels
c)
a)
3380
3380
3400
3400
3420
3420
3440
3440
3460
3460
b)
3380
3380
LSR
3420
3440
3460
T arget
reservoirs
3420
3440
3460
3480
3480
-10
800
Depth (m)
Depth (m)
Depth (m)
CFS
WEPIF
3400
3400
850
900
950 1000 1050 1100 1150 1200
Time (ms)
0
10
20
Q
30
40
3480
0
100
GR (API)
200
3480
170 230 290
350
AC (us/m)
Figure 5. (a) One-shot zero-offset VSP data with six levels, (b) Q-factors estimated
by the LSR (black solid line), CFS (green solid line), and WEPIF methods (red line:
IF calculated with Hilbert transform method). The GR and AC logs are given in (c).
The depth ranges of target reservoirs are indicated by blue crosses in (b)
WEPIF methods and applicationsPoststack Data
Figure 2. Forming the slabs for poststack reflection seismic data and
finding the EPIF: (a) CDP trace, (b) IA of (a), (c) IF of (a), and (d) EPIF.
34
Example of Real reflection seismic data
Known Excellent Clastic Reservoirs
Known
prediction
No gas
Known
Motivation
How do we Identify
Excellent Clastic Reservoirs
Moderate Clastic Reservoirs
Poor Clastic Reservoirs ?
37
WEPIF methods and applications-Prestack
Data
 Response of seismic waves for the gas-bearing
reservoir and gas-not bearing reservoir
 Criterion for distinguishing three classes reservoir
38
Wave field from the complex reservoir
Laterally inhomogeneous media
4800m
Source
Receiver
2500m
Vp=4600m/s Q=100
600m
100m
Vp=4500m/s Q=80
Vp=4200m/s Q=20
Vp=4500m/s Q=80
Vp=4600m/s Q=100
Total 280 shot records with 20m space interval
Partially stacked seismic dada
A
Q=20
B
A
A
B
B
Q=20
Q=20
Partial stacked profiles for small(top-left)、medium(top-right)、wide(bottom)angle gathers
QVA (Q versus angle)
attenuation
Small-angle
Mid-angle
Wide-angle
Well
apparent absorption versus
prestack angles
Excellent Clastic Reservoirs
Class Ⅰ
Small angle: strong absorption
Wide angle: weak absorption
Moderate Reservoirs
Class Ⅱ
The attenuation estimated for the Partially
stacked profile from the small, medium and
wide angle
Small angle: moderate absorption
Wide angle: weak absorption
Poor Reservoirs
Class Ⅲ
Small angle: none absorption
Wide angle: none absorption
Testing the Criterion
Zh10
10.6 104 m3 /Day
召33
Zh17
1.1104 m3 /Day
召17
E
小角度吸收
Small- angle
大角度吸收
The most advantageous
gas-bearing
Apparent attenuation profiles for S061292
Wide -angle
The apparent attenuation of the small-angle profiles near the two
wells, Well zhao10 and well zhao17, is strong. But in the wide-angle
profiles, apparent attenuation near well zhao 10 and zhao 17 is very
different. The former is weak , but the latter is still strong.
Our criterion state that the well zhao 10 has an advantage over the
zhao 17.
Gas well testing productivity shows that well zhao10 has good gasbearing ,the output is ten times as much as well zhao17.
Testing the Criterion
Apparent attenuation in small-angle is stronger than wide-angle . Our criterion
states that the well belong to ClassⅠ.
Well-drilling-results: Su46 drills across a thickness of 24m for sands of He8, well log reveals
ClassⅠ, Gas well testing productivity is
5.4 104 m3 /Day
.
SD46
Small angle Medium angle
Wide angle
Apparent attenuation profiles (S077003)
E
Apparent attenuation profiles(S076272 )
鄂20
N
Small angle
Wide angle
Apparent attenuation in small-angle is stronger than that in wide-angle, but the difference of
them is not larger. For this case, our criterion states the well is Class Ⅱ.
Well-drilling-results: Well E20 drills across a thickness of 18.1m for sands of He8,well log
reveals Class Ⅱ, Gas well testing productivity is 2.5 104 m3 /Day
. .
S60
SE
Small
angle
Medium
angle
Wide angle
Absorption attenuation profiles ( L051735 )
Well-drilling-results: Well Su60 is water well(25m3/d), not gas.
The apparent attenuation is more stronger in small-angle than that in the wide-angle. In this
case, our method does not give an identification.
Pre-drilling-prediction:The apparent attenuation profiles show that the absorption at small-angle
is strong ,but the apparent attenuation at wide-angle is weak. Based on our criterion, the well has
good Physical Properties and gas-bearing, and it belongs to Class Ⅰ.
Well-drilling-results: Su172 drills across a thickness of 23.6m for sands of He8,well log revealed
that it belong to Class Ⅰ. Gas well testing productivity is 4.14 104 m3 /Day
E
Su172
Small angle
Attenuation at small
angle in He8
Wide angle
Absorption at wide
angle in He8
Absorption attenuation profiles (L076763 )
Field result in angle domain –prediction before drilling
As the angle increases, the absorption is weak. It is predicted as Class Ⅰby post-stack.
But it is predicted as Class Ⅱby pre-stack data.
Drilling-result shows that the reservoir belong to Class Ⅱ.
S134
Small
angle
Attenuation at
small angle in H8
Medium
angle
Attenuation at medium
angle in H8
Wide
angle
Attenuation at
Wide angle in H8
Absorption attenuation profiles(L 051735 )
SE
Results
Workload for processing:
13 lines, total 831km
Workload for prediction: 55 wells
Testing:26 wells(Correct rate of 70%)
Predict before drill: 29 wells (Correct rate of 69%)
5 Conclusions



We propose the WEPIF method for estimating Q. A test using
synthetic VSP data shows that the WEPIF method is more stable, more
convenient for selecting parameters, and less sensitive to the
interference from interface reflections than the LSR and CFS methods.
The results, in which the WEPIF method is used for zero-offset VSP
data and post-stack seismic data, show that this method is valid and has
great potential for seismic attenuation estimation and gas reservoir
characterization.
The Q analyzing method for pre-stack multi-angle partially stack
profile can be used to identify the ClassⅠ, Class Ⅱand Class Ⅲ of
reservoir which can not be done with post-stack data.
50

similar documents