Transcript Dmitri Lokshtanov, Multiple removal with local plane waves

Multiple Removal with Local Plane Waves
Dmitri Lokshtanov
1-
Content
• Motivation
• WE multiple suppression operator
• Fast 2D/3D WE approach for simple sea-floor
• 2D/3D WE approach for irregular sea-floor
• Conclusions
2-
Motivation
• Seismic processing and imaging - main challenges:
− Velocity model building for sub-salt and sub-basalt imaging
− Removal of multiples from strong irregular boundaries
3-
Near offset section (no AGC)
Depth migration with water velocity
Input shot gathers (no AGC)
Multiple suppression
• For multiples from complex boundaries the methods based on periodicity
or kinematic discrimination usually don’t work or are not sufficient.
• In such cases the main demultiple tools are based on the Surface
Related Multiple Elimination (SRME) or Wave-Equation (WE) techniques.
8-
SRME (Berkhout, 1982; Verschuur, 1991) –
advantages and limitations
• Does not require any structural information. Predicts all free-surface
multiples
• As a rule becomes less efficient with increased level of interference of
multiples of different orders
• Requires the same dense sampling between sources as between receivers
• Noise in data and poor sampling significantly degrade the prediction quality
• Missing traces required by 3D SRME are reconstructed with least-square
Fourier or Radon interpolation; residual NMO correction; DMO/inverse
DMO; migration/demigration
9-
WE approach versus SRME
• SRME is the method of preference for data from areas with deep sea-floor,
especially when a thick package of strong reflectors is present below the sea-floor
• WE approach is especially efficient when the main free-surface multiples are just
‘pure’ water-layer multiples and peg-legs. Gives usually better results than SRME
when several orders of multiples are involved
• 3D WE approach has less sampling problems than 3D SRME and it gives a
flexiblility in methods for wavefield extrapolation depending on complexity of
structure
10 -
The operator Pg transforms the primary reflection event recorded at receiver 1 into the
multiple event recorded at receiver 2 (Wiggins, 1988; Berryhill & Kim, 1986).
Principles of WE approach
1. Remove pure multiples and receiver - side peg - legs : F1  (1  Pg ) D , where
D are input data; Pg is the receiver - side extrapolation operator.
2. Remove source - side peg - legs : F2  (1  Ps ) F1 , where Ps is the source - side
extrapolation operator.
3. Correct for the first - order water - layer multiple : F (1  Ps )(1  Pg ) D  Ps Dw ,
where Dw is the primary reflection from the water - bottom.
12 -
Adaptive subtraction of predicted multiples
F (1  Ps )(1  Pg ) D  Ps Dw ,
(1)
We apply the ‘scaled version’ of (1) trace by trace to the tau-p transformed CMP
or CS gathers. For each p-trace the operator has the form:
f (t )  d (t )  rg (t )  d g (t )  rs (t )  d s (t )  rsg (t )  d sg (t ),
( 2)
where d (t ), d g (t ), d s (t ), d sg (t ) are p - traces for the input data and the results
of extrapolation through the water-layer from the receiver-side, source-side (of
muted input data) and source-side after receiver-side respectively.
13 -
Wave-equation approach – main features
• All predicted multiples are split into 3 terms, where each term requires
the same amplitude correction
• All source-side and receiver-side multiples of all orders are suppressed
simultaneously in one consistent step
• The prediction and the adaptive subtraction of multiples are performed
in the same domain
• Fast version (WEREM) for a simple sea-floor. Slower version for
irregular sea-floor
14 -
Why in the tau-p domain
• Easier to apply antialiasing protection
• No problems with muting of direct arrival
• Easier to define ‘multiple’ zone of tau-p domain and mute it away
• Estimated reflection coefficients are explicitly angle dependent
15 -
2D WEREM – prediction of multiples - 1
The input Radon transformed CMP gathers D( p, y ) can be represented as follows:

D ( p, y ) 
R ( p, pd ) expi pd ydpd ,

2
(3)
where p g  p  p d 2 , p s  p  p d 2 . For ‘locally’ 1D water-bottom and arbitrary
2D structure below it the results of receiver- side extrapolation D g ( p, y ) is:
Dg ( p , y ) 

R ( p, pd ) expi ( pd y  2q g h)dpd 

2



D ( p, x )  expi  pd ( y  x )  2q g h dpd dx,

2
(4)
where h is the ‘local’ water-bottom depth, while qg is the vertical slowness.
2D WEREM – prediction of multiples - 2
st
For each pair x, y the stationary point pd corresponds to a simple relation:
x  y  h tg g where pg  c sin  g and pg  p  pdst 2 .
R2
R1


y

x

S
(5)
CMP

*
g
water bottom
Geometrical illustration of the result (5). x and y are CMP positions of the primary and
multiple events respectively. They are related by : x  y  0.5 ( R1  R2 )  h tg  g .
Velocity model used to generate synthetic FD data
Constant P sections (angle at the surface is about 3º)
Input
After Werem
After Remul
Constant P sections (angle at the surface is about 15º)
Input
After Werem
After Remul
Velocity model 2 used to generate synthetic FD data
Constant P sections (angle at the surface is about 3º)
 m. residual
Input
After Werem
T. Shetland
 T. Draupne
 T. Brent
Stack before multiple suppression
Stack after Werem
Constant P sections (angle at the surface is about 10º
Input
After Werem
Stack before multiple suppression (left) and after Werem multiple suppression (right).
The pink line shows the expected position of the first-order water-layer peg-leg from the
Top Cretaceous (black line). The multiple period is about 140 msec.
Constant P sections (angle at the surface is about 8º
Input
After Werem multiple suppression
Difference
raw stack
stack after WEREM
500 m
input
WEREM
4000 m
input
WEREM
Improving the results - local prediction /
subtraction of multiples
• Within the same prediction term, for the same CMP and the same p we
have events reflected at different positions along the water bottom
• Inconsistency between prediction and subtraction in case of rapid
variation of sea-floor reflectivity
• The problem is partly solved by applying adaptive subtraction in different
time windows
• Or by making prediction dependent on both p and offset (window)
3D WEREM – basic features
• 3D data can be represented as a sum of plane waves with different vertical
angles and azimuths from the source-side and receiver-side.
• Current quasi 3D marine acquisition does not allow full 4D decomposition
• Decomposition uniquely defines the direction of propagation from the receiverside and is an integral over crossline slownesses from the source-side
• The result of decomposition are used for exact prediction of multiples from the
receiver-side and approximate prediction from the source-side
• The approximation is that the crossline slowness from the source-side is the
same as from the receiver-side (the same azimuth for 1D structures). The
approximation allows us to mix data for flip flop shooting
Decomposition of point source data into plane waves (Weyl Integral):
1
2
 R
exp  i  
2
R
 c  2 
 

exp i  p x x  p y y  q z z 
iq z
  
1
2
 1
2
2
where q z   2  p x  p y  ,
c

Imq z  0.
dp x dp y ,
Reflection from a 3D Earth as a sum of plane waves
(positive direction of p xs , p xr are opposite to the direction of increasing xs )
d ( xs , ys  0, xr , yr ) 
4

2 4
   

~
R
 ( pxr , p yr , pxs , p ys ) S ( pxs , p ys ) exp i  pxs xs  pxr xr  p yr yr dp xs dp ys dp xr dp yr 
    
3

R( p xr , p yr , p xs ) exp i  p xs xs  p xr xr  p yr yr dp xs dp xr dp yr ,
3 
2 
where R( p xr , p yr , p xs ) 

~
S
(
p
,
p
)
R
( p xr , p yr , p xs , p ys ) dp ys .
xs
ys
2 
Radon transformed CMP gathers :
3
d ( x , h, y ) 
2 3
3

2 3
 
h
h

 R exp i  p xs ( x  2 )  p xr ( x  2 )  p yr y  dp xs dp xr dp yr 
 R expi  ph  p
d
x  p yr y dp xs dp xr dp yr ,
p xr  p xs
where p 
, p d  p xs  p xr .
2
2
D ( x, p , y ) 
2
2 
and R ( p, p d , p yr
 R expi  p x  p y dp dp ,
)   D exp i  p x  p y dx dy .
d
yr
d
d
yr
yr
3D prediction from the receiver-side :
2

Dr ( ~
x , p, ~
y) 
2 2
~
~
 R expi  p x  p y  2q z dp dp 
  D( p, x, y ) expi  p ( ~
x  x)  p ( ~
y  y )  2q z dp dp
d
d
yr
r
yr
d
yr
r
d
yr
 dx dy .
  p d ( ~x  x)  p yr ( ~y  y )  2q r z ,
1

2
 1
2 
where q r   2  p r  , p r  p xr2  p yr2
c

Stationary phase condition :

 0  (x  ~
x )  z tg r cos  r ,
p d

1
2
, p xr  p 
pd
2
.

 0  (~
y  y )  z tg r sin  r .
p yr
R2
R1


y

x

S
CMP

*
g
water bottom
Illustration of the stationary phase result: x and y are CMP positions of
the primary and multiple events respectively. They are related by:
x  y  0.5 ( R1  R2 )  h tg  g
37 -
Approximate 3D prediction from the source-side :
2
~
~
Ds ( x , p, y ) 
2 2

~
~


R
exp
i

p
x

p
y  2qs z dpd dp yr 
d
yr

 D( p, x, y) expi  p
d

(~
x  x)  p yr ( ~
y  y )  2qs z dpd dp yr dx dy .
  pd ( ~x  x)  p yr ( ~y  y )  2qs z ,
1

2
1
2
where qs   2  ps  , ps  p xs2  p ys2
c


1
2
, p xs  p 
pd
2
and p ys  p yr .
Stationary phase condition :

 0  (~
x  x)  z tg s cos  s ,
pd

0  (y  ~
y )  z tg s sin  s .
p yr
Input constant P section (small angles)
Predicted multiples – R-side (small angles)
Input constant P section (small angles)
Constant P section – after prediction / subtraction (small angles)
Difference (Input – 3D WEREM), small angles
Input constant P section (larger angles)
Constant P section – after 3D WEREM (larger angles)
Werem - conclusions
• Very efficient when the main assumptions are met: strongest multiples are
water-layer multiples and peg-legs and the sea-floor is simple
• Very fast - each predicted p trace is simply obtained as a sum of timedelayed input traces with the same p from the neighbour CMPs
46 -
WE for irregular sea-floor
• Kinematic prediction of multiples (extrapolation through the water layer)
takes into account coupling between incident and reflected / scattered
plane waves with different slownesses
• Both multiple reflections and diffractions are predicted
• The procedure starts from the Radon transformed CS gathers (no
interleaving is required)
• In 3D exact prediction from the receiver side; approximate prediction from
the source side
47 -
2D prediction of multiples from the receiver
side for irregular sea-floor
1. Extrapolate Radon transformed CS gather D ( p r , x s , ) down to the sea-floor:

W  x, z ( x), xs  
D( pr , xs ) expi  pr ( x  xs )  qr z ( x)dpr ,

2
2. Calculate the amplitude Dg ( psc , xs ) of the reflected/scattered plane wave with
slowness psc (Wenzel et al., 1990)
 dz p sc 
Dg ( p sc , x s )   W  x, z ( x), xs   1 

expi  p sc ( x  x s )  q sc z ( x)dx.

 dx q sc 
48 -
2D prediction of multiples from the source side for
irregular sea-floor
1.
Use FFT to decompose the input data D ( p r , x s ,  ) into contributions with
different propagation angles (wavenumbers) from the source side:
R ( p r , k ,  )   D ( p r , x s ,  ) exp  ikxs dx s ,
where the source-side wavenumber k s is defined as: k s  p r  k .
2.
Extrapolate the results of decomposition down to the sea-floor:
W  x, z ( x ), p r  
3.
49 -
1
R ( p r , k ) exp ik s x  k sz z ( x )dk .

2
Calculate reflected / scattered responses for each k and then use inverse FFT to define
Ds ( p r , x s , )
. All steps are performed in a double loop over pr and over .
Velocity model for FD modelling
Input P-section (zero angle)
Receiver-side prediction
Input P-section (zero angle)
Receiver-side prediction
Input P-section (20 degrees)
Receiver-side prediction
Source-side prediction
Input CS gather
After prediction + subtraction
Input CS gather
After prediction + subtraction
Raw stack with final velocity / mute libraries
Stack after WE + VF
Raw CMP (no AGC)
CMP after WE + VF (no AGC)
Input Constant P section
R-side prediction
Input Constant P section
R-side prediction
Input Constant P section
After adaptive subtraction
3D prediction from the receiver-side - 1
Decompose the recorded CS data into plane waves and then extrapolate each plane
wave down to the sea-floor:
2
W ( x, y , z ) 
2 2
2

2 2
 R ( p , p
x
y
) expi  p x x  p y y  q z z dpx dpy 
 D( p , y ) expip x expi  p
x
r
x
y

( y  yr )  q z z dpy dpx dyr
Innerint egralis calculat ed by t he st at ionaryphase approximation.For st at ionarypoint

  p y ( y  yr )  q z z  q~z z 2  ( y  yr )
62 -

1
2 2
1
1
2
1
1
2
2
2
~
, where q z   2  p x  p y  , q z   2  p x  .
c

c

2
3D prediction from the receiver-side - 2
Extrapolate the wavefield along the sea-floor W(x,y) up to the free-surface and
calculate Radon transformed CS gathers after prediction:
Dg ( p x , yr ) 


2
 W ( x, y) exp ip x C expi  p
x
y

( yr  y )  q z z ( x, y ) dpy dx dy
As above, innerintegralis calculated by the stationaryphase approximation.
63 -
Model for 3D ray tracing of primaries and sea-floor multiples
primary
peg-legs
Modelled CS gathers
Radon transformed CS gathers
Constant P sections (small angle) for line with crossline offset 250m.
Input
3D prediction
2D prediction
Constant P sections (larger angle) for line with crossline offset 250m.
Input
3D prediction
2D prediction
Constant P sections (small angle) for line with crossline offset 250m.
Input
quasi 3D prediction
2D prediction
Constant P sections (larger angle) for line with crossline offset 250m.
Input
quasi 3D prediction
2D prediction
WE for irregular sea-floor
• Both multiple reflections and diffractions are predicted
• Exact 3D prediction of pure water-layer multiples and peg-legs from the
receiver-side
• Quasi 3D prediction of peg-legs from the the source side
• 3-5 times slower than WEREM
71 -
Conclusions
• SRME is the method of preference for data from areas with deep seafloor, especially when a thick package of strong reflectors is present
below the sea-floor
• As a rule the method becomes less efficient when several orders of
multiples are involved
• For such data we use the wave-equation schemes, especially when the
main free-surface multiples are just water-layer multiples and peg-legs
• The 3D WE approach has fewer sampling problems than 3D SRME and it
allows us to use different WE extrapolation schemes for different
complexities of sea-floor and structure below it
72 -
Thank you
73 -