CN1971308A - Inverse operator approximation method of 3D seismics offset imaging - Google Patents
Inverse operator approximation method of 3D seismics offset imaging Download PDFInfo
- Publication number
- CN1971308A CN1971308A CN 200610047210 CN200610047210A CN1971308A CN 1971308 A CN1971308 A CN 1971308A CN 200610047210 CN200610047210 CN 200610047210 CN 200610047210 A CN200610047210 A CN 200610047210A CN 1971308 A CN1971308 A CN 1971308A
- Authority
- CN
- China
- Prior art keywords
- matrix
- imaging
- seismic
- accompanying
- extrapolation
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Images
Landscapes
- Geophysics And Detection Of Objects (AREA)
Abstract
The invention specially relates to a method that lifting, imaging and calculating the seismic reflection data with an easy, fast algorism, and then embed into seismic data processing software platform to form the final result of seismic data processing and displayed as profile chart of underground medium structure, so it is provides a reliable method for array of exploration and development plan, determination of oil-gas well site or other mineral resource survey. The invention adopts the break-in polynomials of partial differential operator matrix power series to approximate the inverse operator directly. The approximate solution in second order or forth order can be obtained at same time with coming at the equation coefficient matrix. The direction distortion can not be generated, because splitting calculation is not adopted in direction. The aim of the invention is resolving the problems of current big amount of calculation, low accuracy, inescapable chaotic phenomenon, out of overcoming of direction distortion, model experiment and large amount of actual data treatment achievement represents that the invention has resolved the problems very well.
Description
One, technical field: the invention belongs to the problem that the seismic migration imaging aspect is studied in the exploration geophysics field, seismic exploration data be by the field instrumentation equipment records to magnetic medium and in the seismic data processing software platform data volume also can be shown as seismic section image with specific software.The present invention proposes a kind of simple algorithm the earthquake reflected wave data is carried out migration imaging processing calculating, for understanding underground medium structure instructs the deployment of exploration scheme, and then determine that exactly the actual position of exploration well location or other subterranean resource provides a kind of practicable reliable method.
Two, background technology: utilizing the earthquake reflected wave data information that stratal configuration is carried out to similarly is one of core technology of seismic prospecting.This technology that is called migration imaging (migration imaging), the spitting image of the CT imaging in the medical science, but its complexity but is to be difficult to analogy.The difficulty of seismic migration imaging is: the data amount is big especially; Signal source is uncertain; The dielectric structure complexity; Interference noise is powerful in the seismic wave propagation, and the unfavorable factor that diffraction, shielding, absorption etc. can't be grasped all produces bad influence to final imaging precision.Meanwhile, with regard to cost, the industrial input of a bite prospecting borehole generally reaches tens million of round Renminbi.Therefore the industrial migration imaging achievement fast and accurately that always requires is as disposing the exploration scheme, determining the foundation of well location.Therefore the three-dimensional seismic offset imaging development of technology is long-expected, from the theory to the commercial Application, the three-dimensional seismic offset imaging technology is subjected to common concern for a long time, has attracted large quantities of outstanding scientists and engineering technical personnel to be engaged in the research of this respect in the geophysics field.
Over nearly two, 30 years, by means of rapid development of computer technology, handle on a large computer system or the network group of planes and set up the software platform that much differs from one another being specifically designed to seismic data, this is a requisite basic normal signal handling implement in the energy industry exploration.But any platform all presses for and improves speed and the precision that migration imaging is handled, and it is gratifying also not having a kind of platform at present, and large quantities of technical forces is being devoted to research and develop new method.Current, what be subjected to generally paying attention to is three-dimensional implicit schemes finite-difference migration method, and the advantage of this method is: adapt to the horizontal change of underground speed, the precision height; Shortcoming is that calculated amount is too big.
1, data volume problem: for the poststack data information of 10km * 10km, if line-spacing and CDP spacing are 25m, sampling interval is 4.00ms, data on each time point reach 400 * 400=160000 so, the equation exponent number that needs to find the solution on this time point is 160000 * 160000, is equation 3000-6000 time that 6 seconds data need be found the solution above-mentioned exponent number for record length.
2, counting yield problem:, calculate numerous and diverse known method for solving usually and can't be competent at because data volume is big.Existing theory, general only research once or the algorithm of finding the solution several times.If, no longer be the problem of efficient and cost, whether the horizon problem of commercial Application has been arranged for thousands of so huge equations of exponent number of asking for loan of the 3D seismic data of small size very inferiorly.
3, direction distortion problem: since thousands of inferior in the wave field extrapolation process, the data of having calculated used, and if can not control distortion in handling well, error transfer accumulation so is enough to make that end result has no practical significance.
So for a long time, people expectation be computing technique fast and accurately, the problem of computing technique aspect has solved, it is very simple that the problem of imaging aspect also just becomes.To invert is relatively more commonly used in division in the existing known computing technique, roughly has alternating direction implicit (ADI) method to be also referred to as along X-Y direction division (Brown, 1983 etc.), to decompose (Rickett andClaerbout, 1998) based on the dimensionality reduction LU on spiral border.The latter is a kind of nonlinear method, can't avoid chaos phenomenon in the calculating, and boundary treatment also has problem.The former because circular symmetry is destroyed, can't overcomes the direction distortion end result is distorted in extrapolation process.Legend and application examples in the back will provide and specify.
Three, summary of the invention:
1, goal of the invention: the invention provides a kind of inverse operator approach method of three-dimensional seismic offset imaging, its fundamental purpose is to solve in the existing seismic migration imaging field that existing calculated amount is big, low precision, can't avoid chaos phenomenon and can't overcome problem such as direction distortion.
2, technical scheme: the present invention is achieved through the following technical solutions:
A kind of inverse operator approach method of three-dimensional seismic offset imaging is characterized in that: this method mainly follows these steps to carry out:
A, at first carry out matrix power series calculating with block polynomial calculating:
The second order difference operator is expressed as a ultra-large symmetric positive definite matrix of coefficients in the matrix equation of wave field extrapolation, this matrix in fact is the matrix of coefficients of implicit schemes extrapolation equation: Au=p, and the extrapolation equation can also be expressed as form
(I-Q)u＝(I+Q)p，
When the spectral radius of matrix Q satisfies ρ (Q)＜1, (I-Q) arranged
-1=I+Q+Q
2+ ..., this is a matrix power series, blocks to be item number matrix polynomial seldom when reality is used;
B, the spectral radius of carrying out matrix are then estimated:
Because Q=α
xT
x+ α
yT
y, T
x, T
yRepresent respectively along the second order difference operator matrix of x and y direction, get Δ x=Δ y, α=α is promptly arranged
x=α
y, matrix Q is got row and norm
As long as
Convergence is just guaranteed;
C, the countermeasure during real data handled:
Because acoustic velocity v, the track pitch Δ y on stratum, CDP separation delta x, time sampling interval Δ t and extrapolation step delta l choose, convergence, counting yield to this method have certain influence, the processing stage gathering and give, do trace interpolation, or simplify the distributed model of stratum acoustic velocity, can make counting yield obtain to a certain degree raising;
In d, the software module in being embedded into the seismic data processing software platform, use above-mentioned inverse operator approximation technique, and select, can obtain the higher image of precision fast by suitable parameter.
3, advantage and effect: the present invention adopts a kind of brand-new method for solving for solving the three-dimensional seismic offset imaging precision problem of long-term puzzlement exploration geophysics industry member.Novel part of the present invention be it abandoned fully this area continued to use decades, be in the division of research in exploring now and approach the conventional thought of inverting.Adopt the matrix power series of partial differential operator and block polynomial expression and directly approach inverse operator.Minimum calculated amount with the explicit scheme method has realized absolute stable implicit schemes calculating.It is a kind of implicit method fully, so absolute stability.But its calculated amount but is surprisingly little, to such an extent as to just obtained the approximate of second order or quadravalence when obtaining matrix of coefficients, this completely new approach drops to the calculated amount of implicit method minimum.Before this because or be that calculated amount is too huge as for cannot bearing, or be to produce distortion in the wave field extrapolation, make the method for finite difference three-D migration rest on the state of theoretical research for a long time.The circular symmetry that the present invention keeps wave field to propagate in data handling procedure fully, final imaging precision is reliably guaranteed, because counting yield improves more than hundred times than classic method, make Three dimensional finite difference method implicit expression migration and imaging techniques become the industrial technology of real practicality again.Can confirm in the practical application that the computing method that adopted have tangible technique effect and remarkable economic efficiency in Liaohe Oil Field calculating is implemented.
Four, description of drawings: accompanying drawing 1 is the classification chart of 3 D earthquake deflection method;
The various correction charts that accompanying drawing 2 is done for the distortion at wave field extrapolation operator circular symmetry;
Accompanying drawing 3 is the impulse response figure of CG skew;
Accompanying drawing 4 is a wave field extrapolation circular symmetry design sketch;
Accompanying drawing 5 is another wave field extrapolation circular symmetry design sketch;
Accompanying drawing 6 is the present invention and alternating direction implicit (ADI) method wave field extrapolation circular symmetry effect comparison diagram;
Accompanying drawing 7 is the preceding stacked profile map of Q171 line skew;
Accompanying drawing 8 is Q171 line offset effect figure;
Accompanying drawing 9 is another Q171 line offset effect figure;
Accompanying drawing 10 is Q171 line offset effect figure of the present invention;
Accompanying drawing 11 is three kinds of Q171 line offset effect comparison diagrams;
Accompanying drawing 12 is the preceding stacked profile map of Q201 line skew;
Accompanying drawing 13 is Q201 line offset effect figure;
Accompanying drawing 14 is another Q201 line offset effect figure;
Accompanying drawing 15 is Q201 line offset effect figure of the present invention;
Accompanying drawing 16 is three kinds of Q201 line offset effect comparison diagrams;
Accompanying drawing 17 is GR software 0471 a line offset effect comparison diagram;
Accompanying drawing 18 is CG software 0471 a line offset effect comparison diagram.
Five, embodiment:
The present invention mainly is by a kind of brand-new computing technique fast and accurately, the earthquake reflected wave data information is carried out migration imaging to be calculated, the data processing product that obtains is at the seismic data processing software platform, can be shown as the migrated seismic section image, this is the core technology that is relied in hydrocarbon resources exploration and the exploitation.Accurate migration imaging achievement can improve the precision of exploration, and its effect in mineral resources exploration exploitation industry is vital.Existing three-D migration formation method has a variety of, sees shown in the accompanying drawing 1.But the precision of migration imaging is unsatisfactory, mainly be because exist calculated amount excessive, can't avoid chaos phenomenon, the problem of aspects such as boundary treatment in calculating.The normal at present computing method that adopt have alternating direction implicit (ADI) method and spiral border Nonlinear Dimension Reduction LU factorization, chaos phenomenon always appears in spiral border Nonlinear Dimension Reduction LU factorization in calculating, the image that accompanying drawing 2,3,4,5 has just adopted alternating direction implicit (ADI) method to form.9 points, 13 points, 17 points, and sexangle McCllalen conversion are adopted in the various corrections that accompanying drawing 2 is done for the distortion at wave field extrapolation operator circular symmetry at 25 in the data.Accompanying drawing 3 is the impulse response of foreign software CG software skew, and as can be seen from the figure Pian Yi noise is serious.Accompanying drawing 4 is the impulse response dropping cut slice of foreign software CG skew, and as can be seen from the figure circular symmetry obviously distorts, and migration noise is serious.Accompanying drawing 5 is the impulse response dropping cut slice of software GR skew, and as can be seen from the figure circular symmetry obviously distorts, and migration noise is serious.Below by embodiment the present invention is done nearly step explanation:
Embodiment: for the poststack data information of 10km * 10km, adopting line-spacing and CDP spacing is 25m, and sampling interval is 4.00ms, uses this method, and the technical essential that relates generally to has:
1, the matrix power series theory of computation with block polynomial calculating: the second order difference operator is expressed as a ultra-large symmetric positive definite matrix of coefficients in the matrix equation of wave field extrapolation, this matrix and right-hand vector are calculated many common parts, in traditional division approach method, these are actually redundant computation, in this method it is used, simplifies or omit a large amount of calculating.In fact implicit schemes extrapolation equation Au=p also has following form:
(I-Q)u＝(I+Q)p
When the spectral radius of matrix Q satisfies ρ (Q)＜1, (I-Q) arranged
-1=I+Q+Q
2+ ...,
What this was that in fact matrix power series represent is the difference operator power series.When reality is used, block and be item number matrix polynomial seldom.So just inversion operation has been become multiplying, inverting has just been become drill, meaning wherein is very great.Realize absolute stable implicit method with explicit calculating, avoided the problem that must solve in many calculating.
2, the spectral radius of matrix is estimated: because Q=α
xT
x+ α
yT
y, T
x, T
yThe second order difference operator matrix of representing x and y direction is respectively got Δ x=Δ y, and α=α is promptly arranged
x=α
y, matrix Q is got row and norm
As long as
Convergence is just guaranteed fully.
3, the consideration during real data is handled: in the practice of Data Processing, the acoustic velocity v on stratum, track pitch Δ y, CDP separation delta x, time sampling interval Δ t and extrapolation step delta l choose, and to the convergence of this method, counting yield has certain influence.This with common seismic migration in choose the extrapolation step-length experience be consistent.Certainly we also have a whole set of departure to improve the scheme of speed of convergence.The appropriate processing of being done the processing stage of gathering and give such as trace interpolation, or is suitably simplified the distributed model of stratum acoustic velocity, also will make counting yield obtain to a certain degree raising.
In d, the software module in being embedded into the seismic data processing software platform, use above-mentioned inverse operator approximation technique, and select, can obtain the precision offset achievement fast and be shown as the image of underground medium structure really by suitable parameters.
Image shown in the accompanying drawing 6 is a model data, and the figure middle and upper part is the section of alternating direction implicit (ADI) method shift pulse level of response, can find out clearly that from figure its circular symmetry obviously distorts; The figure middle and lower part is that inverse operator approaches the section of (IOA) method shift pulse level of response, can clearly find out its strict circular symmetry that keeps from figure, does not have migration noise.Accompanying drawing 7 is the preceding stacked section of Q171 line skew.Accompanying drawing 8 is CG software Q171 line offset effect figure, and as can be seen from the figure relatively the position is in circled, and imaging is unclear.Accompanying drawing 9 is GR software Q171 line offset effect figure, as can be seen from the figure compares the position in circled, image blur.Accompanying drawing 10 is the Q171 line offset effect figure of inverse operator approximatioss of the present invention, and as can be seen from the figure the circle position is more clear than other method imaging.Accompanying drawing 11 is three kinds of Q171 line offset effect comparison diagrams, the top is GR software imaging effect figure, the centre is inverse operator approximatioss imaging effect figure, bottom is foreign software CG imaging effect figure, by the intermediate image lines is clear more as can be seen, imaging effect is significantly better than the imaging effect of upper and lower two images.Accompanying drawing 12 is the preceding stacked section of Q201 line skew, and accompanying drawing 13 is GR software Q201 line offset effect figure, and the imaging of circle position is unintelligible in the image.Accompanying drawing 14 is a CG software migration result, and the circled imaging is unclear in the image.Accompanying drawing 15 is the migration result of inverse operator approximatioss of the present invention (IOA), circle position signal to noise ratio (S/N ratio) height as can be seen among the figure, and imaging is clear.Accompanying drawing 16 is three kinds of Q201 line offset effect comparison diagrams, the left side is the image that software GR forms, the centre is the image that the inverse operator approximatioss forms, and the right is the image that foreign software CG forms, and as can be seen from the figure intermediate image effect is significantly better than the image on left and right both sides.
Enforcement by the technical program, the migration imaging achievement precision height that obtains calculates very fast, can reflect underground medium structure truely and accurately, be exploratory development scheme cloth heat, accurately determining the oil well well location or reconnoitring other mineral resources provides reliable foundation.
Claims (1)
1, a kind of inverse operator approach method of three-dimensional seismic offset imaging is characterized in that: this method mainly follows these steps to carry out:
A, at first carry out matrix power series calculating with block polynomial calculating:
The second order difference operator is expressed as a ultra-large symmetric positive definite matrix of coefficients in the matrix equation of wave field extrapolation, this matrix in fact is implicit schemes extrapolation equation: Au=p, and it has following form
(I-Q)u＝(I+Q)p，
When the spectral radius of matrix Q satisfies ρ (Q)＜1, (I-Q) arranged
-1=I+Q+Q
2+ ..., this is a matrix power series, blocks to be item number matrix polynomial seldom when reality is used;
B, the spectral radius of carrying out matrix are then estimated:
Because Q=α
xT
x+ α
yT
y, T
x, T
yRepresent respectively along the second order difference operator matrix of x and y direction, get Δ x=Δ y, α=α is promptly arranged
x=α
y, matrix Q is got row and norm
, as long as
, convergence is just reliably guaranteed;
C, the countermeasure during real data handled:
Because acoustic velocity v, the track pitch Δ y on stratum, CDP separation delta x, time sampling interval Δ t and extrapolation step delta t choose, convergence, counting yield to this method have certain influence, need the processing stage gathering and give, do trace interpolation, or simplify the distributed model of stratum acoustic velocity, can make counting yield obtain to a certain degree raising;
In d, the software module in being embedded into the seismic data processing software platform, use above-mentioned inverse operator approximation technique, and, can obtain the higher image of precision fast by selected suitable parameters.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CNB2006100472101A CN100547434C (en) | 2006-07-17 | 2006-07-17 | A kind of inverse operator approach method of three-dimensional seismic offset imaging |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CNB2006100472101A CN100547434C (en) | 2006-07-17 | 2006-07-17 | A kind of inverse operator approach method of three-dimensional seismic offset imaging |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN1971308A true CN1971308A (en) | 2007-05-30 |
| CN100547434C CN100547434C (en) | 2009-10-07 |
Family
ID=38112220
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CNB2006100472101A Expired - Fee Related CN100547434C (en) | 2006-07-17 | 2006-07-17 | A kind of inverse operator approach method of three-dimensional seismic offset imaging |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN100547434C (en) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101086533B (en) * | 2007-07-06 | 2010-04-14 | 福州华虹智能科技开发有限公司 | Shallow seismic exploration instrument soft system structure device based on embedded Linux |
| CN103454681A (en) * | 2013-08-28 | 2013-12-18 | 中国石油集团川庆钻探工程有限公司地球物理勘探公司 | Method and equipment for evaluating imaging effect of three-dimensional earthquake observing system |
Family Cites Families (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| FR2691812B1 (en) * | 1992-05-29 | 1994-09-02 | Geophysique Cie Gle | Signal processing method for geophysical prospecting using an operator to extrapolate an improved wave field. |
| US5392255A (en) * | 1992-10-15 | 1995-02-21 | Western Atlas International | Wavelet transform method for downward continuation in seismic data migration |
-
2006
- 2006-07-17 CN CNB2006100472101A patent/CN100547434C/en not_active Expired - Fee Related
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101086533B (en) * | 2007-07-06 | 2010-04-14 | 福州华虹智能科技开发有限公司 | Shallow seismic exploration instrument soft system structure device based on embedded Linux |
| CN103454681A (en) * | 2013-08-28 | 2013-12-18 | 中国石油集团川庆钻探工程有限公司地球物理勘探公司 | Method and equipment for evaluating imaging effect of three-dimensional earthquake observing system |
| CN103454681B (en) * | 2013-08-28 | 2016-03-16 | 中国石油集团川庆钻探工程有限公司地球物理勘探公司 | Evaluate the method and apparatus of 3 D seismic observation system imaging effect |
Also Published As
| Publication number | Publication date |
|---|---|
| CN100547434C (en) | 2009-10-07 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN102395902B (en) | Seismic imaging systems and methods employing a fast target-oriented illumination calculation | |
| CN102169189B (en) | Method for eliminating deep water interlamination multiple | |
| CN104316965B (en) | Prediction method and system for fissure azimuth and intensity | |
| CN108181652A (en) | A kind of subsea node seismic data uplink and downlink wave field numerical method | |
| CN103926619B (en) | Reverse time migration method of three-dimensional VSP data | |
| CN102269820A (en) | Three-dimensional seismic pre-stack reverse-time migration imaging method based on GPU (graphics processing unit) staggered grid with small memory capacity | |
| CN103424777B (en) | A kind of method that improves seismic imaging resolution ratio | |
| CN102466816A (en) | Inversion method for stratum elasticity constant parameter of pre-stack seismic data | |
| CN101910871A (en) | Spectral shaping inversion and migration of seismic data | |
| CN103389513A (en) | Method for improving seismic data resolution by using acoustic logging data constraint inversion | |
| CN106353797A (en) | High-precision earthquake forward modeling method | |
| CN102841378B (en) | Method for predicting reservoir stratum by seismic inversion data | |
| CN102636809B (en) | Method for generating spreading angle domain common image point gathers | |
| CN101487898A (en) | Method for identifying oil, gas and water by longitudinal wave seismic exploration post-stack data | |
| CN103733089A (en) | System and method for subsurface characterization including uncertainty estimation | |
| CN104570066A (en) | Construction Method of Seismic Inversion Low-Frequency Model | |
| CN101290355A (en) | Method for determining carbonate rock hydrocarbon reservoir structure | |
| CN100547434C (en) | A kind of inverse operator approach method of three-dimensional seismic offset imaging | |
| Angus et al. | Seismic waveforms and velocity model heterogeneity: towards a full-waveform microseismic location algorithm | |
| Ikutama et al. | Source modeling for predicting ground motions and permanent displacements very close to the fault trace | |
| CN102466818B (en) | A kind of well earthquake data that utilizes is to the method for anisotropic medium imaging | |
| Taleghani et al. | Post-treatment assessment of induced fracture networks | |
| CN104597496A (en) | Three-dimensional space homing method for speed in two-dimensional seismic data | |
| Chen et al. | Adjoint-state traveltime tomography for azimuthally anisotropic media in spherical coordinates | |
| CN102841374A (en) | Pseudo three-dimensional fast microseism forward modeling method based on scanning surface forward modeling |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| C14 | Grant of patent or utility model | ||
| GR01 | Patent grant | ||
| C17 | Cessation of patent right | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20091007 Termination date: 20110717 |