CN110703327A - Full-band imaging method - Google Patents

Full-band imaging method Download PDF

Info

Publication number
CN110703327A
CN110703327A CN201910970424.3A CN201910970424A CN110703327A CN 110703327 A CN110703327 A CN 110703327A CN 201910970424 A CN201910970424 A CN 201910970424A CN 110703327 A CN110703327 A CN 110703327A
Authority
CN
China
Prior art keywords
imaging
point
seismic
dip angle
dip
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.)
Pending
Application number
CN201910970424.3A
Other languages
Chinese (zh)
Inventor
吴吉忠
付晓飞
石颖
王维红
赵小青
吴吉厚
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Northeast Petroleum University
Original Assignee
Northeast Petroleum University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Northeast Petroleum University filed Critical Northeast Petroleum University
Priority to CN201910970424.3A priority Critical patent/CN110703327A/en
Publication of CN110703327A publication Critical patent/CN110703327A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28Processing seismic data, e.g. analysis, for interpretation, for correction
    • G01V1/30Analysis
    • G01V1/301Analysis for determining seismic cross-sections or geostructures
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28Processing seismic data, e.g. analysis, for interpretation, for correction
    • G01V1/36Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy

Abstract

The invention relates to a full-band imaging method, which is an imaging method capable of carrying out full-band frequency expansion processing on seismic data before migration through pulse deconvolution and improving the signal-to-noise ratio of a seismic imaging section through a self-adaptive aperture; secondly, generating an inclination angle gather corresponding to the sample CDP, and determining an inclination angle imaging area on the inclination angle gather; and thirdly, carrying out pre-stack time migration imaging by using the inclination angle imaging area to obtain a final imaging result. The invention can broaden effective seismic signals in a full frequency band, fully excavate available information in seismic data, determine an accurate imaging area by utilizing the dip angle gather, realize self-adaptive aperture imaging, and achieve the purposes of suppressing migration noise and improving the signal-to-noise ratio of a seismic stack section.

Description

Full-band imaging method
The technical field is as follows:
the invention relates to the technical field of reflection seismic data processing in seismic exploration, in particular to a full-band imaging method.
Background art:
in the conventional imaging processing process of the reflection seismic data, only the seismic data in a frequency band range with a higher signal-to-noise ratio is generally extracted and imaged, and low-frequency and high-frequency data which are in a frequency suppression band and have a relatively lower signal-to-noise ratio are not processed, so that effective signals of seismic original data are not fully mined and utilized, the effective bandwidth of the seismic imaging data is reduced, and the precision and the accuracy of subsequent seismic exploration work are influenced. Full frequency band imaging can bring the problem that the noise is violently enlargied at the in-process of fully excavating the effective signal of earthquake, generally adopts the mode of the post-stack data body filtering and denoising to improve the SNR of seismic section, but can lose, harm the effective signal of earthquake, has violated the original intention of full frequency processing, has weakened the effect of full frequency processing.
In summary, the conventional imaging method often has the disadvantage that data in a frequency suppression band and with a low signal-to-noise ratio cannot be effectively utilized, and the adopted post-stack filtering denoising can bring adverse effects on effective seismic signals aiming at the problem of severe noise amplification caused by full-frequency processing. In view of the above problems, no effective solution has been proposed.
The invention content is as follows:
the invention aims to provide a full-band imaging method, which is used for solving the problems that the conventional imaging method cannot effectively utilize data which are positioned in a frequency suppression band and have a low signal-to-noise ratio and the problem of severe noise amplification caused by full-band processing.
The technical scheme adopted by the invention for solving the technical problems is as follows: this full band imaging method:
step one, performing conventional noise suppression processing on seismic data before stacking;
carrying out full-band frequency extension processing on the seismic data subjected to noise suppression by using a single-channel pulse deconvolution technology, wherein a pulse function is used as expected output of pulse deconvolution;
step three, aiming at a plurality of common center point positions on the seismic imaging section, extracting common center point gathers, recording the gathers as sample CDP, carrying out conventional dynamic correction speed pickup on the gathers, and carrying out interpolation smoothing on the obtained result to be used as an initial migration speed field;
step four, performing offset impulse response calculation on each channel of input data by using an initial offset velocity field, classifying and partially overlapping the offset impulse response according to the offset distance of the input data to form a common reflection point offset distance gather for each common center point position, and determining a final offset velocity field according to whether the in-phase axis in the common reflection point gather is straight or not;
step five, calculating the seismic wave travel time, the migration amplitude and the dip angle of each imaging point of the pre-stack data at the sample CDP position by using the final migration velocity field, arranging according to the dip angle, and partially overlapping to obtain a dip angle gather of the sample CDP position;
the calculation of the travel time, the offset amplitude and the dip angle of the seismic wave is realized as follows: starting from a one-way wave equation and a phase-stabilizing point principle, defining the horizontal coordinate of an excitation point as xsHorizontal coordinate of receiving point is xgPhi is the azimuth angle, the coordinates of the imaging point are (x, T), and the offset velocity of the imaging point is vrmsAnd seismic wave travel time from the excitation point to the imaging point and then to the receiving point:
defining integers
Figure BDA0002231894970000022
Wherein delta t is the sampling interval of the seismic data, the amplitude values of four points i-1, i, i +1 and i +2 on the seismic data are used for quadratic curve interpolation, the seismic amplitude value delta at the travel time tau is solved, and the dimensionless imaging weight coefficient is calculated
Figure BDA0002231894970000023
The magnitude of the shift at the imaging point (x, T) is then
Figure BDA0002231894970000024
And (3) obtaining the inclination angle: the following variables, x1=xs-x,x2=xs+xg-2 x; calculating the distance variable a ═ x1cosφ;b=x1x2+2(vrmsT)2;c=x2cosφ;Then k is equal to (ad)2-bc)/(c2-d2) The tilt angle theta is arctan [ (x)1+kcosφ)/(vrmsT)];
Determining dip angle imaging areas aiming at the dip angle gather corresponding to the CDP positions of the sample, and determining the dip angle imaging areas corresponding to the rest CDP positions by utilizing the determined dip angle imaging areas in an interpolation continuation mode;
step seven, utilizing the final offset velocity field to calculate the dip angle of each imaging point on all CDP positions of the data before the stack, if the dip angle is positioned in the dip angle imaging area, continuously calculating the travel time and the amplitude of the imaging point to finish imaging, and if the dip angle is positioned outside the dip angle imaging area, repeating the operation of the step on the next imaging point;
and step eight, finishing the calculation of all CDP position imaging points to obtain a final seismic imaging section.
In the scheme, the step three is realized by the following steps of extracting a common-center-point gather from a plurality of common-center-point positions on the seismic imaging section: and selecting the position of the common center point at equal intervals on the conventional section according to the change condition of the structure, and distributing the offset distances at equal intervals from small to large when extracting the common center point gather.
In the scheme, the method for determining the dip angle imaging area for the dip angle gather corresponding to the CDP position of the sample in the step six comprises the following steps: a decision condition is introduced to the system that,
Figure BDA0002231894970000031
psi (t, x) is an imaging value corresponding to an angle x of an imaging point at a time depth t position in a dip gather,
Figure BDA0002231894970000032
Angle of imaging point at time depth t position of dip trace concentration from theta-eta1To theta + eta2The imaging superposition value in the range, theta is the true formation dip of the dip gather time depth t position, gamma is the threshold value, and is set as
Figure BDA0002231894970000033
AmaxWhen the condition is satisfied, the imaging result at the time depth t position tends to be stable and presents a phase-stable state [ theta-eta ] for the maximum absolute amplitude of the dip gather1,θ+η2]Is the dip imaging area to be determined at the time depth t position of the dip gather.
The invention has the following beneficial effects:
1. the invention carries out full-band frequency extension processing on the seismic data before migration, so that the data with low signal-to-noise ratio in the frequency suppression band of the seismic data is fully mined and effectively utilized.
2. The invention utilizes the dip angle channel set to determine the accurate imaging area, suppresses signal noise in the imaging process and solves the problem of severe noise amplification caused by full frequency processing.
Drawings
FIG. 1 is a conventional deconvolution stack cross-section.
FIG. 2 is a pulse deconvolution stack cross section.
FIG. 3 compares conventional deconvolution to pulse deconvolution spectra.
Fig. 4 is a full-frequency imaging method using a conventional post-stack denoising profile.
Fig. 5 is a cross-section of full-frequency imaging using the method of the present invention.
FIG. 6 shows a comparison of the cross-section spectrum obtained by the method of the present invention and the conventional post-stack denoising method used in full-frequency imaging.
Detailed Description
The invention is further illustrated below:
the full-band imaging method takes a certain seismic data of No. 1 block of the Jidong oilfield, and specifically comprises the following steps:
1. the pre-stack seismic data is subjected to conventional noise suppression processing, the data sampling interval is 2ms, the recording duration of the seismic signals is 5000ms, the CDP interval is 25m, the minimum offset is 100m, the maximum offset is 5400m, and the offset interval is 50 m.
(2) And carrying out full-band frequency broadening processing on the seismic data after the noise suppression by utilizing a single-channel pulse deconvolution technology, wherein a pulse function is taken as expected output of the pulse deconvolution.
(3) And (3) extracting a common-center offset gather aiming at a plurality of common-center positions on the section, recording the common-center offset gather as a sample CDP, performing conventional dynamic correction speed pickup on the gathers, and performing interpolation smoothing on the obtained result to be used as an initial offset speed field.
(4) And (3) performing offset impulse response calculation on each channel of data by using an initial offset velocity field, classifying and partially superposing impulse responses according to the offset distance of the data at the position of each common center point to form a common reflection point offset distance gather.
(5) And (3) calculating the seismic wave travel time, the migration amplitude and the dip angle of each imaging point of the pre-stack data at the sample CDP position by using the final migration velocity field, arranging according to the dip angle, and partially stacking to obtain the dip angle gather of the sample CDP position.
(6) And determining dip angle imaging areas aiming at the dip angle gather corresponding to the sample CDP, and determining the dip angle imaging areas corresponding to the rest CDP positions by utilizing the determined plurality of dip angle imaging areas in an interpolation continuation mode.
(7) And calculating the dip angle of each imaging point on all CDP positions of the data before the stacking by using the final offset velocity field, continuously calculating the travel time and amplitude of the imaging point if the dip angle is positioned in the dip angle imaging area, finishing imaging, and repeating the operation of the step for the next imaging point if the dip angle is positioned outside the dip angle imaging area.
(8) And finishing the calculation of all CDP position imaging points to obtain a final seismic imaging section.
Fig. 1 is a conventional deconvolution cross section, fig. 2 is a pulse deconvolution cross section obtained after the use, it can be seen that the signal-to-noise ratio of the deep layer in fig. 2 is significantly low, and high-frequency noise is developed violently, fig. 3 is a comparison of conventional deconvolution and pulse deconvolution spectra, the dotted line on the graph represents the conventional deconvolution spectrum, the solid line represents the pulse deconvolution spectrum, and it can be seen that the high-frequency part of the spectrum is significantly broadened after the pulse deconvolution processing. Fig. 4 is a cross section of full-frequency imaging after conventional post-stack denoising, fig. 5 is a cross section obtained by full-frequency imaging after the method of the present invention is adopted, it can be seen that interlayer information is richer in fig. 5, details are depicted more clearly, fig. 6 is a cross section frequency spectrum obtained by full-frequency imaging after conventional post-stack denoising and the method of the present invention, a dotted line on the graph represents a frequency spectrum of the conventional post-stack denoising method, a solid line represents a frequency spectrum processed by the method of the present invention, it can be seen that the frequency spectrum obtained by the method of the present invention is obviously superior to the conventional post-stack denoising method in terms of bandwidth, although the conventional post-stack denoising suppresses high-frequency noise, improves the signal-to-noise ratio, but damages the high-frequency effective signal, and the method of the present invention does not damage the high-frequency effective.

Claims (3)

1. A full band imaging method, comprising the steps of:
step one, performing conventional noise suppression processing on seismic data before stacking;
carrying out full-band frequency extension processing on the seismic data subjected to noise suppression by using a single-channel pulse deconvolution technology, wherein a pulse function is used as expected output of pulse deconvolution;
step three, aiming at a plurality of common center point positions on the seismic imaging section, extracting common center point gathers, recording the gathers as sample CDP, carrying out conventional dynamic correction speed pickup on the gathers, and carrying out interpolation smoothing on the obtained result to be used as an initial migration speed field;
step four, performing offset impulse response calculation on each channel of input data by using an initial offset velocity field, classifying and partially overlapping the offset impulse response according to the offset distance of the input data to form a common reflection point offset distance gather for each common center point position, and determining a final offset velocity field according to whether the in-phase axis in the common reflection point gather is straight or not;
step five, calculating the seismic wave travel time, the migration amplitude and the dip angle of each imaging point of the pre-stack data at the sample CDP position by using the final migration velocity field, arranging according to the dip angle, and partially overlapping to obtain a dip angle gather of the sample CDP position;
the calculation of the travel time, the offset amplitude and the dip angle of the seismic wave is realized as follows: starting from a one-way wave equation and a phase-stabilizing point principle, defining the horizontal coordinate of an excitation point as xsHorizontal coordinate of receiving point is xgPhi is the azimuth angle, the coordinates of the imaging point are (x, T), and the offset velocity of the imaging point is vrmsAnd seismic wave travel time from the excitation point to the imaging point and then to the receiving point:
Figure FDA0002231894960000011
defining integers
Figure FDA0002231894960000012
Wherein delta t is the sampling interval of the seismic data, the amplitude values of four points i-1, i, i +1 and i +2 on the seismic data are used for quadratic curve interpolation, the seismic amplitude value delta at the travel time tau is solved, and the dimensionless imaging weight coefficient is calculated
Figure FDA0002231894960000013
The magnitude of the shift at the imaging point (x, T) is then
Figure FDA0002231894960000014
And (3) obtaining the inclination angle: the following variables, x1=xs-x,x2=xs+xg-2 x; calculating the distance variable a ═ x1cosφ;b=x1x2+2(vrmsT)2;c=x2cosφ;
Figure FDA0002231894960000021
Then k is equal to (ad)2-bc)/(c2-d2) The inclination angle theta is arctan[(x1+kcosφ)/(vrmsT)];
Determining dip angle imaging areas aiming at the dip angle gather corresponding to the CDP positions of the sample, and determining the dip angle imaging areas corresponding to the rest CDP positions by utilizing the determined dip angle imaging areas in an interpolation continuation mode;
step seven, utilizing the final offset velocity field to calculate the dip angle of each imaging point on all CDP positions of the data before the stack, if the dip angle is positioned in the dip angle imaging area, continuously calculating the travel time and the amplitude of the imaging point to finish imaging, and if the dip angle is positioned outside the dip angle imaging area, repeating the operation of the step on the next imaging point;
and step eight, finishing the calculation of all CDP position imaging points to obtain a final seismic imaging section.
2. The full band imaging method of claim 1, wherein: the step three is realized by that a plurality of common center point positions on the seismic imaging section are extracted, and the common center point gather is extracted: and selecting the position of the common center point at equal intervals on the conventional section according to the change condition of the structure, and distributing the offset distances at equal intervals from small to large when extracting the common center point gather.
3. The full band imaging method of claim 1, wherein: the sixth step is a method for determining the dip angle imaging area for the dip angle gather corresponding to the sample CDP position:
a decision condition is introduced to the system that,
Figure FDA0002231894960000022
psi (t, x) is an imaging value corresponding to an angle of an imaging point x at a time depth t position in the dip gather,angle of imaging point at time depth t position of dip trace concentration from theta-eta1To theta + eta2Imaging stack value in the range, θ is tiltThe true formation dip angle at the time depth t position of the angle trace concentration, gamma is a threshold value and is set as
Figure FDA0002231894960000024
AmaxWhen the condition is satisfied, the imaging result at the time depth t position tends to be stable and presents a phase-stable state [ theta-eta ] for the maximum absolute amplitude of the dip gather1,θ+η2]Is the dip imaging area to be determined at the time depth t position of the dip gather.
CN201910970424.3A 2019-10-13 2019-10-13 Full-band imaging method Pending CN110703327A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910970424.3A CN110703327A (en) 2019-10-13 2019-10-13 Full-band imaging method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910970424.3A CN110703327A (en) 2019-10-13 2019-10-13 Full-band imaging method

Publications (1)

Publication Number Publication Date
CN110703327A true CN110703327A (en) 2020-01-17

Family

ID=69200220

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910970424.3A Pending CN110703327A (en) 2019-10-13 2019-10-13 Full-band imaging method

Country Status (1)

Country Link
CN (1) CN110703327A (en)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111538077A (en) * 2020-05-07 2020-08-14 中国石油天然气集团有限公司 Pre-stack depth migration method and device based on dip angle constraint
CN115857007A (en) * 2022-11-15 2023-03-28 东北石油大学 Full-signal imaging method and device

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20090257308A1 (en) * 2008-04-11 2009-10-15 Dimitri Bevc Migration velocity analysis methods
CN101957455A (en) * 2010-09-20 2011-01-26 中国海洋石油总公司 Method of three-dimensional preserved-amplitude pre-stack time migration
CN104297789A (en) * 2014-10-23 2015-01-21 中国科学院地质与地球物理研究所 Three-dimensional dip angle domain stationary phase pre-stack time migration method and system
CN106896408A (en) * 2017-03-23 2017-06-27 中国石油天然气股份有限公司 A kind of angle domain prestack time migration method
CN109254324A (en) * 2018-10-19 2019-01-22 中国石油天然气股份有限公司 Full range protects width seismic data processing technique and device

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20090257308A1 (en) * 2008-04-11 2009-10-15 Dimitri Bevc Migration velocity analysis methods
CN101957455A (en) * 2010-09-20 2011-01-26 中国海洋石油总公司 Method of three-dimensional preserved-amplitude pre-stack time migration
CN104297789A (en) * 2014-10-23 2015-01-21 中国科学院地质与地球物理研究所 Three-dimensional dip angle domain stationary phase pre-stack time migration method and system
CN106896408A (en) * 2017-03-23 2017-06-27 中国石油天然气股份有限公司 A kind of angle domain prestack time migration method
CN109254324A (en) * 2018-10-19 2019-01-22 中国石油天然气股份有限公司 Full range protects width seismic data processing technique and device

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
吴吉忠 等: "倾角域自适应孔径叠前时间偏移", 《石油地球物理勘探》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111538077A (en) * 2020-05-07 2020-08-14 中国石油天然气集团有限公司 Pre-stack depth migration method and device based on dip angle constraint
CN115857007A (en) * 2022-11-15 2023-03-28 东北石油大学 Full-signal imaging method and device

Similar Documents

Publication Publication Date Title
CN108196305B (en) Mountain land static correction method
WO2006054181A1 (en) Method for processing at least two sets of seismic data
AU2014280832B2 (en) Seismic data spectrum restoring and broadening
CN110703327A (en) Full-band imaging method
CN111103621A (en) Analysis method for superposition of active source common imaging points and multiple surface waves
CN109765615A (en) A kind of inversion method for stratigraphic quality factor and device
CN102073064A (en) Method for improving velocity spectrum resolution by using phase information
CN110389377B (en) Microseism offset imaging positioning method based on waveform cross-correlation coefficient multiplication
CN109884693B (en) Self-adaptive trend velocity spectrum solving method and system
CN105093318B (en) A kind of adaptive wave equation wave field extrapolation static correcting method
CN111239814B (en) Shallow profile data mechanical interference suppression method based on same-phase axis frequency division tracking smoothing
CN110244383B (en) Geological lithology comprehensive model establishing method based on near-surface data
CN104155688A (en) High precision weighted stack method
CN111766626A (en) Algorithm for realizing three-dimensional space fitting of seismic data
CN112213784B (en) One-time processing fast static correction method for complex surface seismic data
CN113671566B (en) Method for calculating crack parameters based on depth domain seismic data
CN110879413A (en) Ray parameter domain converted wave static correction method and system
CN112946742B (en) Method for picking up accurate superposition velocity spectrum
CN110568491B (en) Quality factor Q estimation method
CN109884701B (en) Geologic body scattering angle guiding depth imaging method
CN110554435A (en) method for constructing quality factor body by using micro-logging data
CN102890288B (en) Interval velocity inversion method for earthquake waves
CN110531409B (en) Method for extracting azimuth information of marine narrow-azimuth seismic data
CN110737019B (en) Processing method and system for improving velocity spectrum precision
CN114594517B (en) CRS common reflection surface element superposition imaging method

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
WD01 Invention patent application deemed withdrawn after publication
WD01 Invention patent application deemed withdrawn after publication

Application publication date: 20200117