WO2017024523A1 - 一种射线弹性参数的反演方法 - Google Patents
一种射线弹性参数的反演方法 Download PDFInfo
- Publication number
- WO2017024523A1 WO2017024523A1 PCT/CN2015/086647 CN2015086647W WO2017024523A1 WO 2017024523 A1 WO2017024523 A1 WO 2017024523A1 CN 2015086647 W CN2015086647 W CN 2015086647W WO 2017024523 A1 WO2017024523 A1 WO 2017024523A1
- Authority
- WO
- WIPO (PCT)
- Prior art keywords
- data
- seismic
- profile
- elastic impedance
- impedance
- Prior art date
Links
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V11/00—Prospecting or detecting by methods combining techniques covered by two or more of main groups G01V1/00 - G01V9/00
Definitions
- the invention relates to the field of oil and gas field exploration technology, in particular to reservoir prediction and oil and gas detection technology in oil and gas field development process, and specifically relates to a method and system for inverting ray elasticity parameters.
- the prestack elastic inversion is more based on the EI inversion of the angle-dependent elastic impedance inversion method proposed by Connolly in 1999 based on the Zoeppritz equation approximation.
- the method assumes that the longitudinal and shear wave velocity ratios of the whole seismic section are constant at different depths, which is not consistent with reality; it is assumed that the seismic wave propagates in the formation according to the normal incidence angle; as the incident angle of the longitudinal wave increases, the error of the EI inversion gradually increases.
- the amplitude of the elastic impedance EI varies greatly with the angle of incidence, and it is difficult to control by the normalization method.
- the change of the amplitude is easy to cause errors in lithology and fluid identification.
- Wang Yanghua and Ma Jinfeng respectively disclose a ray path elastic impedance inversion method.
- the elastic impedance is expressed as a function of the longitudinal wave impedance and the longitudinal and transverse wave velocity ratio.
- the ray elastic impedance does not need to assume that the longitudinal and transverse wave velocity ratio is constant and does not need to be incident.
- Angle normalization is unstable and the anti-noise ability is poor, thus limiting the practical application of the ray elastic impedance.
- Liu Lihui is equivalent to the binomial expansion of the ray-elastic impedance method in 2011, and expresses the ray-elastic impedance as a function of the longitudinal wave impedance and the shear wave impedance.
- Liu Lihui's approximate ray-elastic impedance method has improved stability and anti-noise ability compared with the ray-elastic impedance method at a small incident angle (less than 30°).
- the incident angle is greater than 30°
- the error with the ray elastic impedance is much larger than 10%, and the error becomes larger as the incident angle increases.
- this method cannot effectively distinguish reservoir physical properties.
- the above-mentioned ray-elastic impedance inversion has many advantages, but all of them have the problem of analytical instability.
- the present invention provides an inversion method and system for ray elasticity parameters, which is an accurate ray elastic parameter inversion scheme.
- the ray-elastic impedance inversion formula is further deduced by the ray-elastic impedance of multiple incident angles, which can accurately and quickly solve the shot.
- An object of the present invention is to provide a method for inverting a ray elasticity parameter, comprising: collecting seismic data and logging data of a prospecting area, the logging data including depth data, density data, longitudinal wave velocity data, and transverse wave Speed data; converting the logging data from the depth domain to the time domain according to the seismic data; determining the ray elastic impedance corresponding to the plurality of incident angles of the well point position according to the logging data in the converted time domain; The ray elastic impedance and the seismic data determine a ray elastic impedance profile corresponding to the plurality of incident angles; and determine a reservoir elastic parameter according to the ray elastic impedance profile corresponding to the plurality of incident angles.
- One of the objectives of the present invention is to provide an inversion system for radioelastic parameters, comprising: a data acquisition device for collecting seismic data and logging data of the exploration area, the logging data including depth data, density Data, longitudinal wave velocity data, shear wave velocity data; time domain conversion means for converting said logging data from a depth domain to a time domain based on said seismic data; and a ray elastic impedance determining means for converting
- the logging data in the time domain determines the ray elastic impedance corresponding to the plurality of incident angles of the well point position;
- the elastic impedance profile determining device is configured to determine the radiant elastic impedance profile corresponding to the plurality of incident angles according to the ray elastic impedance and the seismic data.
- the reservoir elastic parameter determining device is configured to determine the reservoir elastic parameter according to the ray elastic impedance profile corresponding to the plurality of incident angles.
- the invention has the beneficial effects of providing a method and system for inverting ray elastic parameters, belonging to reservoir prediction and oil and gas detection technology in oil and gas field exploration and development process, and is an accurate ray elastic parameter inversion scheme, utilizing A plurality of incident angle ray elastic impedances are used to further derive the ray elastic impedance inversion formula, and an accurate expression for solving the longitudinal and transverse wave velocity ratio, longitudinal wave impedance and shear wave impedance is obtained, which can accurately and quickly solve the elastic parameters without generating The error, and there is no multi-solution, and can be adapted to the seismic gather data of large incident angle to meet the requirements of oil and gas seismic exploration.
- FIG. 1 is a flowchart of a method for inverting a ray elasticity parameter according to an embodiment of the present invention
- step S102 in FIG. 1 is a specific flowchart of step S102 in FIG. 1;
- FIG. 3 is a specific flowchart of step S104 in FIG. 1;
- FIG. 5 is a specific flowchart of Embodiment 2 of step S304 in FIG. 3;
- FIG. 6 is a specific flowchart of step S105 in FIG. 1;
- FIG. 7 is a structural block diagram of an inversion system for a ray elasticity parameter according to an embodiment of the present invention.
- FIG. 8 is a block diagram showing a specific structure of a time domain conversion apparatus 200 in a system for inverting a ray elasticity parameter according to an embodiment of the present invention
- FIG. 9 is a block diagram showing a specific structure of an elastic impedance profile determining apparatus 400 in a system for inverting a radiation elastic parameter according to an embodiment of the present invention.
- FIG. 10 is a block diagram showing a specific structure of a first embodiment of a relative elastic impedance profile determining module 404 in a system for inverting a ray elastic parameter according to an embodiment of the present invention
- FIG. 11 is a block diagram showing a specific structure of a second embodiment of a relative elastic impedance profile determining module 404 in an inversion system for a ray elastic parameter according to an embodiment of the present invention
- FIG. 12 is a block diagram showing a specific structure of a reservoir elastic parameter determining apparatus 500 in a system for inverting a ray elasticity parameter according to an embodiment of the present invention
- Figure 13 is a schematic view showing a seismic section of the S087017 line in the specific embodiment
- Figure 14 is a schematic diagram of a logging curve of a su116 well in a specific embodiment
- 15 is a schematic diagram of layer calibration of a seismic section by using a synthetic record in a specific embodiment
- 16 is a schematic diagram of a ray-elastic impedance synthetic recording model as a function of incident angle in a specific embodiment
- 17 is a schematic diagram of a ray elastic impedance model with an incident angle of 5° in a specific embodiment
- 18 is a schematic diagram of a radiation elastic impedance model of an incident angle of 15° in a specific embodiment
- 19 is a schematic diagram of a ray elastic impedance model of an incident angle of 30° in a specific embodiment
- 20 is a schematic cross-sectional view of an angular gather of an incident angle of 5° in a specific embodiment
- 21 is a schematic cross-sectional view of an angle gather of an incident angle of 15° in a specific embodiment
- Figure 22 is a schematic cross-sectional view of an angular gather of an incident angle of 30° in a specific embodiment
- Figure 23 is a schematic cross-sectional view showing a radiation elastic impedance of an incident angle of 5° in a specific embodiment
- Figure 24 is a schematic cross-sectional view showing a radiation elastic impedance of an incident angle of 15° in a specific embodiment
- Figure 25 is a schematic cross-sectional view showing a radiation elastic impedance of an incident angle of 30° in a specific embodiment
- 26 is a schematic cross-sectional view showing a longitudinal-to-transverse wave velocity ratio obtained by inversion of a ray elastic impedance parameter in a specific embodiment
- FIG. 27 is a schematic cross-sectional view of a longitudinal wave impedance obtained by inversion of a ray elastic impedance parameter in a specific embodiment
- FIG. 28 is a schematic cross-sectional view of a transverse wave impedance obtained by inversion of a ray elastic impedance parameter in a specific embodiment.
- the invention belongs to reservoir prediction and oil and gas detection technology in oil and gas field exploration and development process, and is an accurate method for inversion of ray elastic parameters, which can be solved by precise expression in the inversion process from ray elastic impedance to elastic parameter. Out of the elastic parameters, there is no error, and there is no multiplicity.
- FIG. 1 is a specific flowchart of a method for inverting a ray elasticity parameter according to the present invention. As shown in FIG. 1, the method includes:
- S101 Collect seismic data and logging data of the exploration area, and the logging data includes depth data, density data, longitudinal wave velocity data, and shear wave velocity data.
- seismic data is collected in the exploration area, and longitudinal and transverse waves of the seismic are excited in the exploration well to obtain longitudinal wave velocity, shear wave velocity data, and depth data, and the density data of the formation is obtained in the exploration zone.
- FIG. 14 it is a schematic diagram of the logging curve of the su116 well in the specific embodiment, and the longitudinal wave velocity data, the shear wave velocity data, the density data, and the depth data are visible.
- FIG. 2 is a specific flowchart of step S102.
- S103 Determine a ray elastic impedance corresponding to the plurality of incident angles of the well point position according to the log data of the converted time domain.
- the incident angle may be from 2 to 6.
- the following is an example of three.
- the ray elastic reflection coefficient is Increasing the incident angle ⁇ increases the type of sandstone belonging to the third type. Therefore, the three incident angles can be selected to be 5°, 15°, 30°.
- REI (5), REI (15), and REI (30) are radiation elastic impedances corresponding to incident angles of 5°, 15°, and 30°.
- the ray elastic impedance corresponding to the three incident angles is calculated from the log data using the following formula:
- ⁇ is the incident angle
- ⁇ is the longitudinal wave velocity
- ⁇ is the shear wave velocity
- ⁇ is the density
- REI( ⁇ i ) is the ray elastic impedance of the incident angle ⁇
- k is the ratio between the density change rate and the shear wave velocity change rate.
- Table 1 is a table of elastic impedance data.
- FIG. 3 is a specific flowchart of step S104.
- FIG. 6 is a specific flowchart of step S105.
- FIG. 2 is a specific flowchart of step S102. As shown in FIG. 2, the step specifically includes:
- the seismic data is collected in the exploration area, and the seismic profile is obtained, as shown in FIG. 13 , which is a schematic diagram of the seismic profile of the S087017 line.
- a standard seismic wavelet is selected or a seismic wavelet is extracted from the well seismic profile.
- S203 Synthesize a seismic record according to the longitudinal wave velocity data, the density data, and the seismic wavelet.
- a 30hz Ray's wavelet is selected, and the seismic velocity is synthesized by combining the longitudinal wave velocity and the density velocity with the seismic wavelet.
- S204 Perform layer calibration on the seismic section by using the seismic record to obtain a layer calibration result and a time-depth relationship.
- the seismic section is calibrated by a synthetic seismic record to obtain a layer calibration result, a time-depth relationship, a target layer segment is determined, and a seismic reflection layer is picked up, as shown in FIG. 15 in the specific embodiment.
- a schematic diagram of the horizon calibration of the seismic section is performed by synthetic recording, and Tp7, Tp8 and Tc2 are the seismic reflection horizons picked up.
- S205 Interpret the seismic horizon according to the layer calibration result
- S206 convert the logging data from the depth domain to the time domain according to the time-depth relationship.
- the logging data is converted from the depth domain to the time domain by uniform sampling or non-uniform sampling.
- Table 2 shows the deep conversion and re-sampling data table of the logging curve. Depth domain log curve data, time domain log curve data obtained by the right side frame conversion.
- FIG. 3 is a specific flowchart of step S104. As shown in FIG. 3, the step specifically includes:
- S301 Establish a plurality of ray elastic impedance models according to the ray elastic impedance corresponding to the plurality of incident angles of the well point position and the seismic horizon.
- the three incident angles are Three ray elastic impedance models are established for the ray elastic impedance values and the seismic horizon.
- Figure 17 is a schematic diagram of a ray-elastic impedance model with an incident angle of 5°
- Figure 18 is a schematic diagram of a ray-elastic impedance model with an incident angle of 15°
- Figure 19 is a schematic diagram of a ray-elastic impedance model with an incident angle of 30°.
- the three ray elastic impedance models are established by using the inverse distance interpolation modeling algorithm and the Kriging geological modeling algorithm.
- three ray elastic impedance models are filtered by a digital low pass filter to obtain three ray elastic impedance low frequency models.
- the low pass frequency of the digital low pass filter is determined according to the range of low frequency missing in the spectrum of the seismic profile, such as 5-8 hz.
- FIG. 4 is a specific flowchart of Embodiment 1 of Step S304
- FIG. 5 is a specific flowchart of Embodiment 2 of Step S304.
- S305 Determine a ray elastic impedance profile according to the ray elastic impedance low frequency model and the relative elastic impedance profile.
- the beam-elastic impedance low-frequency model and the relative elastic impedance profile of the corresponding incident angle are added under the constraint of the well-elastic impedance of the well to obtain a three-angle beam-elastic impedance profile.
- Figure 23 shows the incident angle as Schematic diagram of a 5° beam-elastic impedance profile
- FIG. 24 is a schematic diagram of a beam-elastic impedance profile with an incident angle of 15°
- FIG. 25 is a schematic diagram of a beam-elastic impedance profile with an incident angle of 30°.
- FIG. 4 is a specific flowchart of Embodiment 1 of step S304. As shown in FIG. 4, in the first embodiment, the step specifically includes:
- S401 Determine pre-stack gather data according to the seismic data.
- the seismic data arranged by the track before the completion of the basic data processing superposition is the pre-stack gather data.
- S402 Determine a range of each angle superposition according to the pre-stack gather data.
- S403 Extract an angle gather superposition profile corresponding to the ray elastic impedance low frequency model from the pre-stack gather data according to the range of the angle superposition.
- the angular gathers of the three incident angles corresponding to the model are extracted from the seismic pre-stack gather data.
- 20 is a schematic diagram of a superposition profile of an angle gather of an incident angle of 5°
- FIG. 21 is a schematic diagram of an angle gather superimposed profile of an incident angle of 15°
- FIG. 22 is a schematic diagram of an angular gather superimposed profile of an incident angle of 30°.
- S404 compressing the seismic wavelets into the angle gather superposition profile to obtain a reflection coefficient profile.
- the seismic gather wavelet is compressed by a pulse deconvolution algorithm for the angle gather superposition profile of the three incident angles to obtain a reflection coefficient profile.
- S405 Obtain a relative elastic impedance profile according to the reflection coefficient profile.
- three relative elastic impedance profiles are obtained using a recursive inversion algorithm based on S403.
- FIG. 5 is a specific flowchart of Embodiment 2 of step S304. As shown in FIG. 5, in the second embodiment, the step specifically includes:
- S501 Determine pre-stack gather data according to the seismic data.
- the seismic data arranged by the track before the completion of the basic data processing superposition is the pre-stack gather data.
- S502 Determine a range of each angle superposition according to the pre-stack gather data.
- S503 Extract an angle gather superposition profile corresponding to the ray elastic impedance low frequency model from the pre-stack gather data according to the range of the angle superposition.
- three models corresponding to the model are extracted from the seismic pre-stack gather data.
- the angle of the incident angle gathers the profile.
- 20 is a schematic diagram of a superposition profile of an angle gather of an incident angle of 5°
- FIG. 21 is a schematic diagram of an angle gather superimposed profile of an incident angle of 15°
- FIG. 22 is a schematic diagram of an angular gather superimposed profile of an incident angle of 30°.
- S504 Perform time difference adjustment on the angle gather superposition profile, and eliminate time difference of reflection main phase axis of the main target interval on the three sections by time difference adjustment.
- S505 Perform amplitude matching on the angle gather superposition section after the time difference adjustment, and make the total energy of the main target interval on the three sections uniform by the amplitude matching method, thereby further eliminating the influence of spherical diffusion and the like. That is, the time-division and amplitude consistency correction is performed on the angle gather superposition profiles of the three incident angles.
- S506 compressing the seismic wavelet according to the amplitude-matched angular gather superposition profile to obtain a reflection coefficient profile.
- the seismic gather wavelet is compressed by a pulse deconvolution algorithm for the angle gather superposition profile of the three incident angles to obtain a reflection coefficient profile.
- S507 Obtain a relative elastic impedance profile according to the reflection coefficient profile.
- three relative elastic impedance profiles are obtained using a recursive inversion algorithm based on S505.
- FIG. 6 is a specific flowchart of step S105. As shown in FIG. 6, the step specifically includes:
- S601 Determine a longitudinal-to-transverse wave velocity ratio according to a ray elastic impedance profile corresponding to the plurality of incident angles.
- the step is performed by using the following formula:
- ⁇ 1 , ⁇ 2 , and ⁇ 3 are three incident angles
- ⁇ is the ratio of longitudinal to transverse wave velocities
- ⁇ 1 and ⁇ 2 are intermediate parameters
- k is a proportional coefficient between the rate of change of density and the rate of change of shear wave velocity
- REI ( ⁇ 1 ) and REI ( ⁇ 2 ) are the beam elastic impedance profiles corresponding to the incident angles ⁇ 1 and ⁇ 2 , respectively.
- S602 Determine a longitudinal wave impedance according to the aspect ratio.
- S603 Determine a transverse wave impedance according to the longitudinal wave-to-surface wave velocity ratio and the longitudinal wave impedance.
- S604 Determine a reservoir elastic parameter according to the aspect ratio, the longitudinal wave impedance, and the transverse longitudinal wave impedance, wherein the reservoir elastic parameter comprises a Poisson's ratio, a Young's modulus, and a shear modulus.
- 26 is a schematic cross-sectional view showing longitudinal and transverse wave velocity ratios obtained by inversion of a ray elastic impedance parameter in a specific embodiment
- FIG. 27 is a schematic diagram showing a longitudinal wave impedance profile obtained by inversion of a ray elastic impedance parameter in a specific embodiment
- FIG. 28 is a ray in a specific embodiment.
- Schematic diagram of the shear wave impedance profile obtained from the inversion of the elastic impedance parameters shows that the longitudinal wave impedance near the well 116 of the box 8 is low, and the longitudinal wave impedance near the T26 well is high; the shear wave impedance does not change much near the two wells; the longitudinal and transverse wave velocity ratio profile is The vicinity of the S116 well is low, the T26 well is near high value, and the longitudinal and transverse wave velocity ratio caused by the decrease of the longitudinal wave velocity after gas inclusion is low, which is consistent with the test gas result S116 well 41,000 square meters/day and T26 well no output.
- the inversion method of the ray elastic parameter provided by the present invention is an accurate ray elastic parameter inversion scheme, and the ray elastic impedance inversion formula is calculated by using the ray elastic impedance of multiple incident angles. Further derivation, the elastic parameters can be solved accurately and quickly, and there is no multi-solution, and it can be adapted to the seismic gather data of large incident angle to meet the requirements of oil and gas seismic exploration.
- FIG. 7 is a structural block diagram of a system for inverting a ray elasticity parameter according to an embodiment of the present invention. As shown in FIG. 7, the system includes:
- the data collection device 100 is configured to collect seismic data and logging data of the exploration area, and the logging data includes depth data, density data, longitudinal wave velocity data, and shear wave velocity data.
- seismic data is collected in the exploration area, and longitudinal and transverse waves of the seismic are excited in the exploration well to obtain longitudinal wave velocity, shear wave velocity data, and depth data, and the density data of the formation is obtained in the exploration zone.
- FIG. 14 it is a schematic diagram of the logging curve of the su116 well in the specific embodiment, and the longitudinal wave velocity data, the shear wave velocity data, the density data, and the depth data are visible.
- the time domain conversion device 200 is configured to convert the logging data from a depth domain to a time domain according to the seismic data.
- FIG. 8 is a block diagram showing a specific structure of the time domain conversion device 200.
- the ray elastic impedance determining device 300 is configured to determine a ray elastic impedance corresponding to the plurality of incident angles of the well point position according to the converted time domain logging data.
- the incident angle may be from 2 to 6.
- the following is an example of three.
- the ray elastic reflection coefficient is Increasing the incident angle ⁇ increases
- the three incident angles can be selected to be 5°, 15°, 30°.
- REI (5), REI (15), and REI (30) are radiation elastic impedances corresponding to incident angles of 5°, 15°, and 30°.
- the ray elastic impedance corresponding to the three incident angles is calculated from the log data using the following formula:
- ⁇ is the incident angle
- ⁇ is the longitudinal wave velocity
- ⁇ is the shear wave velocity
- ⁇ is the density
- REI( ⁇ i ) is the ray elastic impedance of the incident angle ⁇
- k is the ratio between the density change rate and the shear wave velocity change rate.
- Table 1 is a table of elastic impedance data.
- the elastic impedance profile determining device 400 is configured to determine a beam elastic impedance profile corresponding to the plurality of incident angles according to the beam elastic impedance and the seismic data.
- FIG. 9 is a block diagram showing a specific structure of the elastic impedance profile determining device 400.
- the reservoir elasticity parameter determining device 500 is configured to determine a reservoir elastic parameter according to a ray elastic impedance profile corresponding to the plurality of incident angles.
- FIG. 12 is a block diagram showing a specific structure of the reservoir elasticity parameter determining device 500.
- FIG. 8 is a block diagram showing a specific structure of a time domain conversion apparatus 200 in a system for inverting a ray elasticity parameter according to an embodiment of the present invention. As shown in FIG. 8, the time domain conversion apparatus 200 specifically includes:
- the seismic section determining module 201 is configured to determine a seismic section according to the seismic data.
- the seismic data is collected in the exploration area, and the seismic profile is obtained, as shown in FIG. 13 , which is a schematic diagram of the seismic profile of the S087017 line.
- the seismic wavelet selection module 202 is configured to select a seismic wavelet.
- a standard seismic wavelet is selected or a seismic wavelet is extracted from the well seismic profile.
- the seismic record synthesis module 203 is configured to synthesize the seismic record based on the longitudinal wave velocity data, the density data, and the seismic wavelet.
- a 30hz Ray's wavelet is selected, and the seismic velocity is synthesized by combining the longitudinal wave velocity and the density velocity with the seismic wavelet.
- the layer calibration module 204 is configured to perform layer calibration on the seismic section by using the seismic record to obtain a layer calibration result and a time-depth relationship.
- the seismic section is calibrated by a synthetic seismic record to obtain a layer calibration result, a time-depth relationship, a target layer segment is determined, and a seismic reflection layer is picked up, as shown in FIG. 15 in the specific embodiment.
- a schematic diagram of the horizon calibration of the seismic section is performed by synthetic recording, and Tp7, Tp8 and Tc2 are the seismic reflection horizons picked up.
- a seismic horizon determining module 205 configured to interpret the seismic horizon according to the layer calibration result
- a time domain conversion module 206 configured to convert the logging data from a depth domain to a time according to the time-depth relationship area.
- the logging data is converted from the depth domain to the time domain by uniform sampling or non-uniform sampling.
- Table 2 shows the deep conversion and re-sampling data table of the logging curve. Depth domain log curve data, time domain log curve data obtained by the right side frame conversion.
- FIG. 9 is a block diagram showing a specific structure of an elastic impedance profile determining apparatus 400 in a system for inverting a radioelastic parameter according to an embodiment of the present invention.
- the elastic impedance profile determining apparatus 400 specifically includes:
- the ray elastic impedance model establishing module 401 is configured to establish a plurality of ray elastic impedance models according to the ray elastic impedance corresponding to the plurality of incident angles of the well point position and the seismic horizon. In a specific embodiment, the constraint of the seismic reflection horizon Next, three ray elastic impedance models are established from the ray elastic impedance values of the three incident angles and the seismic horizon.
- Figure 17 is a schematic diagram of a ray-elastic impedance model with an incident angle of 5°
- Figure 18 is a schematic diagram of a ray-elastic impedance model with an incident angle of 15°
- Figure 19 is a schematic diagram of a ray-elastic impedance model with an incident angle of 30°.
- the three ray elastic impedance models are established by using the inverse distance interpolation modeling algorithm and the Kriging geological modeling algorithm.
- the low frequency missing range determining module 402 is configured to determine a range of low frequency missing in the spectrum of the seismic section;
- the filtering module 403 is configured to filter the ray elastic impedance model according to the range of the low frequency missing to obtain a ray elastic impedance low frequency model.
- three ray elastic impedance models are filtered by a digital low pass filter to obtain three ray elastic impedance low frequency models.
- the low pass frequency of the digital low pass filter is determined according to the range of low frequency missing in the spectrum of the seismic profile, such as 5-8 hz.
- the relative elastic impedance profile determining module 404 is configured to determine a relative elastic impedance profile corresponding to the plurality of incident angles according to the seismic data.
- the ray elastic impedance profile determining module 405 is configured to determine a ray elastic impedance profile according to the ray elastic impedance low frequency model and the relative elastic impedance profile.
- the beam-elastic impedance low-frequency model and the relative elastic impedance profile of the corresponding incident angle are added under the constraint of the well-elastic impedance of the well to obtain a three-angle beam-elastic impedance profile.
- 23 is a schematic diagram of a beam elastic impedance profile with an incident angle of 5°
- FIG. 24 is a schematic diagram of a beam elastic impedance profile with an incident angle of 15°
- FIG. 25 is a schematic diagram of a beam elastic impedance profile with an incident angle of 30°.
- the relative elastic impedance profile determining module 404 specifically includes:
- the pre-stack gather data determining unit 4041 is configured to determine the pre-stack gather data according to the seismic data.
- the seismic data arranged by the track before the completion of the basic data processing superposition is the pre-stack gather data.
- the range determining unit 4042 is configured to determine a range of each angle superposition according to the pre-stack gather data.
- the range of each angle superposition is selected according to the magnitude of the signal-to-noise ratio of the pre-stack gather data. For example, in a concrete In the example, when the angle gathers with the incident angle of 5° are extracted, the gathers of 3° to 7° are selected for superposition; when the angle gathers with the incident angle of 15° is extracted, the lanes of 13° to 17° are selected. The set is superimposed; when the angle gathers with an incident angle of 30° is extracted, the gathers of 28° to 32° are selected for superposition.
- the angle gather superposition section extracting unit 4043 is configured to extract an angle gather set superposition section corresponding to the ray elastic impedance low frequency model from the pre-stack gather data according to the range of the angle superposition.
- the angular gathers of the three incident angles corresponding to the model are extracted from the seismic pre-stack gather data.
- 20 is a schematic diagram of a superposition profile of an angle gather of an incident angle of 5°
- FIG. 21 is a schematic diagram of an angle gather superimposed profile of an incident angle of 15°
- FIG. 22 is a schematic diagram of an angular gather superimposed profile of an incident angle of 30°.
- the reflection coefficient profile determining unit 4044 is configured to compress the seismic wavelets into the angular gathers superimposed profile to obtain a reflection coefficient profile.
- the seismic gather wavelet is compressed by a pulse deconvolution algorithm for the angle gather superposition profile of the three incident angles to obtain a reflection coefficient profile.
- the relative elastic impedance profile determining unit 4045 is configured to obtain a relative elastic impedance profile according to the reflection coefficient profile.
- three relative elastic impedance profiles are obtained by using a recursive inversion algorithm based on the reflection coefficient profile determining unit 4043.
- FIG. 11 is a block diagram of a specific structure of the second embodiment of the relative elastic impedance profile determining module 404. As shown in FIG. 11, the second elastic impedance profile determining module 404 further includes:
- the time difference adjustment unit 4046 is configured to perform time difference adjustment on the angle gather set superimposed section, and eliminate the time difference of the main target layer reflection in the same phase axis on the three sections by the time difference adjustment.
- the amplitude matching unit 4047 is configured to perform amplitude matching on the angle gather superimposed profile after the time difference adjustment, and the total energy of the main target intervals on the three sections is matched by the amplitude matching method, thereby further eliminating the influence of factors such as spherical diffusion. That is, the time-division and amplitude consistency correction is performed on the angle gather superposition profiles of the three incident angles.
- FIG. 12 is a block diagram showing a specific structure of a reservoir elastic parameter determining apparatus 500 in an inversion system of a ray elastic parameter according to an embodiment of the present invention.
- the reservoir elastic parameter determining apparatus 500 specifically includes:
- the longitudinal and transverse wave velocity ratio determining module 501 is configured to determine the longitudinal and transverse wave velocity ratio according to the ray elastic impedance profile corresponding to the plurality of incident angles. In a specific embodiment, when the incident angle is three, the step is performed by the following formula:
- ⁇ 1 , ⁇ 2 , and ⁇ 3 are three incident angles
- ⁇ is the ratio of longitudinal to transverse wave velocities
- ⁇ 1 and ⁇ 2 are intermediate parameters
- k is a proportional coefficient between the rate of change of density and the rate of change of shear wave velocity
- REI ( ⁇ 1 ) and REI ( ⁇ 2 ) are the beam elastic impedance profiles corresponding to the incident angles ⁇ 1 and ⁇ 2 , respectively.
- the longitudinal wave impedance determining module 502 is configured to determine the longitudinal wave impedance according to the aspect ratio.
- the shear wave impedance determining module 503 is configured to determine the shear wave impedance according to the aspect ratio and the longitudinal wave impedance.
- the reservoir elasticity parameter determining module 504 is configured to determine a reservoir elastic parameter according to the aspect ratio, the longitudinal wave impedance, and the transverse longitudinal wave impedance, wherein the reservoir elastic parameter comprises a Poisson's ratio, a Young's modulus, and a shear Modulus.
- 26 is a schematic cross-sectional view showing longitudinal and transverse wave velocity ratios obtained by inversion of a ray elastic impedance parameter in a specific embodiment
- FIG. 27 is a schematic diagram showing a longitudinal wave impedance profile obtained by inversion of a ray elastic impedance parameter in a specific embodiment
- FIG. 28 is a ray in a specific embodiment.
- Schematic diagram of the shear wave impedance profile obtained from the inversion of the elastic impedance parameters shows that the longitudinal wave impedance near the well 116 of the box 8 is low, and the longitudinal wave impedance near the T26 well is high; the shear wave impedance does not change much near the two wells; the longitudinal and transverse wave velocity ratio profile is The vicinity of the S116 well is low, the T26 well is near high value, and the longitudinal and transverse wave velocity ratio caused by the decrease of the longitudinal wave velocity after gas inclusion is low, which is consistent with the test gas result S116 well 41,000 square meters/day and T26 well no output.
- the inversion system for the ray elasticity parameter provided by the present invention is an accurate ray elastic parameter inversion scheme, and the ray elastic impedance inversion formula is calculated by using the ray elastic impedance of multiple incident angles. Further derivation, the elastic parameters can be solved accurately and quickly, and there is no multi-solution, and it can be adapted to the seismic gather data of large incident angle to meet the requirements of oil and gas seismic exploration.
- Seismic data is collected in the exploration area, and the seismic section is processed, as shown in Fig. 13, which is a schematic diagram of the seismic section of the S087017 line.
- the longitudinal and transverse waves of the earthquake are excited in the exploration area, and the longitudinal wave velocity and the shear wave velocity data are obtained.
- the formation density is obtained in the exploration area, as shown in Fig. 14.
- It is a schematic diagram of the logging curve of the su116 well in the specific example, and the longitudinal wave velocity data is visible.
- Select 30hz Rake wavelet, synthesize seismic record with longitudinal wave velocity and formation density use the synthetic seismic record to stratify the seismic profile, determine the target interval, and pick up the seismic reflection horizon.
- Figure 15 shows the synthesis in the specific example.
- a schematic diagram of the horizon calibration of the seismic section is recorded, and Tp7, Tp8, and Tc2 are the seismic reflection horizons picked up.
- the radiation elastic reflection coefficient increases with the increase of the incident angle ⁇ and belongs to the third type of sandstone. Therefore, the three incident angles can be selected to be 5°, 15°, 30°.
- the time-depth relationship obtained by the layer calibration is used to convert the logging data from the depth domain to the time domain by uniform sampling or non-uniform sampling.
- the left frame in Table 1 is the depth domain logging curve data, and the time in the right frame conversion is obtained. Domain logging curve data.
- the ray elastic impedance corresponding to the three incident angles is calculated from the log data in the time domain. Under the constraint of the seismic reflection horizon, three ray elastic impedance models are respectively established from the ray elastic impedance values of the three incident angles.
- Fig. 17 is a schematic diagram of the ray elastic impedance model with an incident angle of 5°
- Fig. 18 is an incident angle of 15°.
- Fig. 19 is a schematic diagram of the ray elastic impedance model with an incident angle of 30°.
- the three ray elastic impedance models are established by using the inverse distance interpolation modeling algorithm and the Kriging geological modeling algorithm.
- the three ray elastic impedance models are filtered by a digital low-pass filter to obtain three ray elastic impedance low-frequency models, and the low-pass frequency of the digital low-pass filter is determined according to the range of the low-frequency missing in the spectrum of the seismic section;
- the low pass frequency of the digital low pass filter is selected from 5 to 8 hz.
- FIG. 20 is a schematic diagram of the angle gather superposition profile of the incident angle of 5°
- FIG. 21 is the angle gather superposition profile of the incident angle of 15°.
- Fig. 22 is a schematic diagram of the superposition of the angular gathers with an incident angle of 30°.
- the range of angle superposition is selected: when the angle gathers with an incident angle of 5° is extracted, the gathers of 3° to 7° are selected for superposition; the incident is extracted.
- the angle gathers with an angle of 15° is superimposed, the gathers of 13° to 17° are selected for superposition; when the angle gathers with an incident angle of 30° is extracted, the gathers of 28° to 32° are selected for superposition.
- the time difference and amplitude are corrected for the angle gather superposition profiles of the three incident angles.
- the time difference is used to eliminate the time difference of the main target interval reflection in the three sections, and the main target layer on the three sections is obtained by the amplitude matching method.
- the total energy of the segments is consistent, further eliminating the effects of spherical diffusion and other factors.
- the angled superposition profile of the three incident angles is compressed by the pulse deconvolution algorithm to obtain the reflection coefficient profile.
- the recursive inversion algorithm is used to obtain three relative elastic impedance profiles; Under the constraint of ray elastic impedance, the low-frequency model of the ray-elastic impedance and the relative elastic-impedance profile of the corresponding incident angle are added; the ray-elastic impedance profile of three angles is obtained.
- 23 is a schematic diagram of a beam elastic impedance profile with an incident angle of 5°
- FIG. 24 is a schematic diagram of a beam elastic impedance profile with an incident angle of 15°
- FIG. 25 is a schematic diagram of a beam elastic impedance profile with an incident angle of 30°.
- the longitudinal-to-transverse-wave velocity ratio, the longitudinal-wave impedance, and the transverse-wave impedance are calculated.
- the other elastic parameters of the reservoir are obtained from the longitudinal-to-transverse-wave velocity ratio, the longitudinal-wave impedance, and the transverse-longitudinal-wave impedance. Flexibility Parameter inversion.
- Other reservoir elastic parameters are Poisson's ratio, Young's modulus, and shear modulus.
- FIG. 26 is a schematic cross-sectional view showing longitudinal and transverse wave velocity ratios obtained by inversion of a ray elastic impedance parameter in a specific embodiment
- FIG. 27 is a schematic diagram showing a longitudinal wave impedance profile obtained by inversion of a ray elastic impedance parameter in a specific embodiment
- FIG. 28 is a ray in a specific embodiment.
- Schematic diagram of the shear wave impedance profile obtained from the inversion of the elastic impedance parameters shows that the longitudinal wave impedance near the well 116 of the box 8 is low, and the longitudinal wave impedance near the T26 well is high; the shear wave impedance does not change much near the two wells; the longitudinal and transverse wave velocity ratio profile is The vicinity of the S116 well is low, the T26 well is near high value, and the longitudinal and transverse wave velocity ratio caused by the decrease of the longitudinal wave velocity after gas inclusion is low, which is consistent with the test gas result S116 well 41,000 square meters/day and T26 well no output.
- the inversion method and system for ray elastic parameters proposed by the present invention belong to reservoir prediction and oil and gas detection technology in oil and gas field exploration and development process, and are an accurate ray elastic parameter inversion scheme, which utilizes A plurality of incident angle ray elastic impedances are used to further derive the ray elastic impedance inversion formula, and an accurate expression for solving the longitudinal and transverse wave velocity ratio, longitudinal wave impedance and shear wave impedance is obtained, which can accurately and quickly solve the elastic parameters without generating The error, and there is no multi-solution, and can be adapted to the seismic gather data of large incident angle to meet the requirements of oil and gas seismic exploration.
- the storage medium may be a magnetic disk, an optical disk, a read-only memory (ROM), or a random access memory (RAM).
Landscapes
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
一种射线弹性参数的反演方法,所述方法包括:采集探区的地震资料、测井资料(S101),所述的测井资料包括深度数据、密度数据、纵波速度数据、横波速度数据;根据所述的地震资料将所述的测井资料由深度域转换为时间域(S102);根据转换得到的时间域的测井资料确定井点位置的多个入射角对应的射线弹性阻抗(S103);根据所述的射线弹性阻抗以及地震资料确定多个入射角对应的射线弹性阻抗剖面(S104);根据多个入射角对应的射线弹性阻抗剖面确定储层弹性参数(S105)。该反演方法可以准确快速的求解出弹性参数,而且不存在多解性,同时可适合于大入射角的地震道集数据,满足油气地震勘探的要求。
Description
本发明关于油气田勘探技术领域,特别是关于油气田开发过程中的储层预测、油气检测技术,具体的讲是一种射线弹性参数的反演方法及系统。
目前,叠前弹性反演较多采用的是Connolly于1999年提出的基于Zoeppritz方程近似的、与角度有关的弹性阻抗反演方法EI反演。该方法假设整个地震剖面不同深度上纵、横波速度比为常数,不符合实际;假设地震波在地层中按照常入射角度传播;随着纵波入射角度的增大,EI反演的误差逐渐增大。
在实际应用中,弹性阻抗EI的幅值随入射角度变化很大,用归一化方法难以控制,幅值的变化易导致岩性和流体识别的错误。王仰华和马劲风分别公开了一种射线路径弹性阻抗反演方法,将弹性阻抗表示为纵波阻抗和纵横波速度比的函数关系式,射线弹性阻抗不需要假设纵横波速度比为常数,也无需进行入射角归一化处理。但是在求解弹性参数时,解析式不稳定,抗噪能力差,因此限制了射线弹性阻抗的实际应用。
刘力辉等于2011年对射线弹性阻抗方法进行二项式展开,将射线弹性阻抗表达成纵波阻抗和横波阻抗的函数。通过误差推导和试验分析,在小入射角度(小于30°)的情况下,与射线弹性阻抗方法相比,刘力辉的近似射线弹性阻抗方法在稳定性和抗噪能力方面有所改善。但是在入射角大于30°时,与射线弹性阻抗误差远大于10%,而且随着入射角的增大误差越来越大。而对于入射角大于30°的弹性反演,这种方法不能有效区分储层物性特征。上述射线弹性阻抗反演具有诸多优势,但都存在解析式不稳定的问题。
因此,如何提出一种新的射线弹性参数反演方案,其能够准确快速的求解出弹性参数,同时可适合于大入射角的地震道集数据、满足油气地震勘探的要求是本领域亟待解决的技术难题。
发明内容
为了克服现有技术中的射线弹性阻抗反演方案存在的解析式不稳定的问题,本发明提供了一种射线弹性参数的反演方法及系统,是一种精确的射线弹性参数反演方案,利用多个入射角的射线弹性阻抗对射线弹性阻抗反演公式进行了进一步推导,可以准确快速的求解出弹
性参数,而且不存在多解性,同时可适合于大入射角的地震道集数据,满足油气地震勘探的要求。
本发明的目的之一是,提供一种射线弹性参数的反演方法,包括:采集探区的地震资料、测井资料,所述的测井资料包括深度数据、密度数据、纵波速度数据、横波速度数据;根据所述的地震资料将所述的测井资料由深度域转换为时间域;根据转换得到的时间域的测井资料确定井点位置的多个入射角对应的射线弹性阻抗;根据所述的射线弹性阻抗、地震资料确定多个入射角对应的射线弹性阻抗剖面;根据多个入射角对应的射线弹性阻抗剖面确定储层弹性参数。
本发明的目的之一是,提供了一种射线弹性参数的反演系统,包括:资料采集装置,用于采集探区的地震资料、测井资料,所述的测井资料包括深度数据、密度数据、纵波速度数据、横波速度数据;时间域转换装置,用于根据所述的地震资料将所述的测井资料由深度域转换为时间域;射线弹性阻抗确定装置,用于根据转换得到的时间域的测井资料确定井点位置的多个入射角对应的射线弹性阻抗;弹性阻抗剖面确定装置,用于根据所述的射线弹性阻抗以及地震资料确定多个入射角对应的射线弹性阻抗剖面;储层弹性参数确定装置,用于根据多个入射角对应的射线弹性阻抗剖面确定储层弹性参数。
本发明的有益效果在于,提供了一种射线弹性参数的反演方法及系统,属于油气田勘探、开发过程中的储层预测和油气检测技术,是一种精确的射线弹性参数反演方案,利用多个入射角射线弹性阻抗,对射线弹性阻抗反演公式进行了进一步推导,得到了求解纵横波速度比、纵波阻抗、横波阻抗的精确表示式,可以准确快速的求解出弹性参数,不会产生误差,而且不存在多解性,同时可适合于大入射角的地震道集数据,满足油气地震勘探的要求。
为让本发明的上述和其他目的、特征和优点能更明显易懂,下文特举较佳实施例,并配合所附图式,作详细说明如下。
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1为本发明实施例提供的一种射线弹性参数的反演方法的流程图;
图2为图1中的步骤S102的具体流程图;
图3为图1中的步骤S104的具体流程图;
图4为图3中的步骤S304的实施方式一的具体流程图;
图5为图3中的步骤S304的实施方式二的具体流程图;
图6为图1中的步骤S105的具体流程图;
图7为本发明实施例提供的一种射线弹性参数的反演系统的结构框图;
图8为本发明实施例提供的一种射线弹性参数的反演系统中的时间域转换装置200的具体结构框图;
图9为本发明实施例提供的一种射线弹性参数的反演系统中的弹性阻抗剖面确定装置400的具体结构框图;
图10为本发明实施例提供的一种射线弹性参数的反演系统中的相对弹性阻抗剖面确定模块404的实施方式一的具体结构框图;
图11为本发明实施例提供的一种射线弹性参数的反演系统中的相对弹性阻抗剖面确定模块404的实施方式二的具体结构框图;
图12为本发明实施例提供的一种射线弹性参数的反演系统中的储层弹性参数确定装置500的具体结构框图;
图13为具体实施例中S087017测线的地震剖面示意图;
图14为具体实施例中su116井的测井曲线示意图;
图15为具体实施例中用合成记录对地震剖面进行层位标定示意图;
图16为具体实施例中随入射角变化的射线弹性阻抗合成记录模型示意图;
图17为具体实施例中入射角为5°的射线弹性阻抗模型示意图;
图18为具体实施例中入射角15°的射线弹性阻抗模型示意图;
图19为具体实施例中入射角30°的射线弹性阻抗模型示意图;
图20为具体实施例中入射角5°的角度道集叠加剖面示意图;
图21为具体实施例中入射角15°的角度道集叠加剖面示意图;
图22为具体实施例中入射角30°的角度道集叠加剖面示意图;
图23为具体实施例中入射角5°的射线弹性阻抗剖面示意图;
图24为具体实施例中入射角15°的射线弹性阻抗剖面示意图;
图25为具体实施例中入射角30°的射线弹性阻抗剖面示意图;
图26为具体实施例中射线弹性阻抗参数反演得到的纵横波速度比剖面示意图;
图27为具体实施例中射线弹性阻抗参数反演得到的纵波阻抗剖面示意图;
图28为具体实施例中射线弹性阻抗参数反演得到的横波阻抗剖面示意图。
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
本发明属于油气田勘探、开发过程中的储层预测和油气检测技术,是一种精确的射线弹性参数反演方法,在由射线弹性阻抗到弹性参数的反演过程中可以由精确的表达式求解出弹性参数,不会产生误差,而且不存在多解性。
图1为本发明提出的一种射线弹性参数的反演方法的具体流程图,由图1可知,所述的方法包括:
S101:采集探区的地震资料、测井资料,所述的测井资料包括深度数据、密度数据、纵波速度数据、横波速度数据。
在具体的实施例中,在探区内采集地震资料,在探区井中激发地震纵波和横波,得到纵波速度、横波速度数据、深度数据,在探区测井得到地层的密度数据。如图14所示,为具体实施例中su116井的测井曲线示意图,可见纵波速度数据、横波速度数据、密度数据、深度数据。
S102:根据所述的地震资料将所述的测井资料由深度域转换为时间域。图2为步骤S102的具体流程图。
S103:根据转换得到的时间域的测井资料确定井点位置的多个入射角对应的射线弹性阻抗。
在具体的实施方式中,所述的入射角可为2个至6个。下面以3个为例进行说明。在具体的实施方式中,确定出目的层段后(如图16中的虚线范围),确定射线弹性反射系数随入射角变化的类型,如图16所示的实施例中,射线弹性反射系数随入射角θ的增大而增大属于第三类砂岩类型。因此,三个入射角可以选择为5°、15°、30°。REI(5)、REI(15)、REI(30)为入射角为5°、15°、30°对应的射线弹性阻抗。
在时间域采用下式由测井数据计算三个入射角对应的射线弹性阻抗:
公式中,θ为入射角,α为纵波速度,β为横波速度,ρ为密度,REI(θi)为入射角为θ的射线弹性阻抗,k是密度变化率与横波速度变化率之间的比例系数。表1为弹性阻抗数据表。
表1
S104:根据所述的射线弹性阻抗以及地震资料确定多个入射角对应的射线弹性阻抗剖面。图3为步骤S104的具体流程图。
S105:根据多个入射角对应的射线弹性阻抗剖面确定储层弹性参数。图6为步骤S105的具体流程图。
图2为步骤S102的具体流程图,由图2可知,该步骤具体包括:
S201:根据所述的地震资料确定地震剖面。
在具体的实施方式中,在探区内采集地震资料,处理得到地震剖面,如图13所示,为S087017测线的地震剖面示意图。
S202:选择地震子波。
在具体的实施方式中,选择标准地震子波或者从井旁地震剖面中提取地震子波。
S203:根据所述的纵波速度数据、密度数据以及地震子波合成地震记录。
在具体的实施方式中,选择30hz的雷克子波,用纵波速度、密度速度结合地震子波合成地震记录。
S204:用所述的地震记录对所述的地震剖面进行层位标定,得到层位标定结果、时深关系。
在具体的实施方式中,用合成的地震记录对地震剖面进行层位标定,得到层位标定结果、时深关系,确定目的层段,并拾取地震反射层位,如图15为具体实施例中用合成记录对地震剖面进行层位标定示意图,Tp7、Tp8、Tc2即为拾取的地震反射层位。
S205:根据所述的层位标定结果解释得到地震层位;S206:根据所述的时深关系将所述的测井资料由深度域转换为时间域。
根据层位标定得到的时深关系,用均匀采样或者非均匀采样把测井资料由深度域转换为时间域,如表2为测井曲线深时转换及重新采样数据表,其中左侧框为深度域测井曲线数据,右侧框转换得到的时间域测井曲线数据。
表2
图3为步骤S104的具体流程图,由图3可知,该步骤具体包括:
S301:根据井点位置的多个入射角对应的射线弹性阻抗以及地震层位建立多个射线弹性阻抗模型,在具体的实施方式中,在地震反射层位的约束下,由三个入射角的射线弹性阻抗值以及地震层位分别建立三个射线弹性阻抗模型。图17为入射角5°的射线弹性阻抗模型示意图,图18为入射角15°的射线弹性阻抗模型示意图,图19为入射角30°的射线弹性阻抗模型示意图。所述的建立三个射线弹性阻抗模型采用反距离内插建模算法、克里金地质建模算法。
S302:确定所述地震剖面的频谱中低频缺失的范围;
S303:根据所述低频缺失的范围对所述射线弹性阻抗模型进行滤波,得到射线弹性阻抗低频模型。
在具体的实施方式中,用数字低通滤波器对三个射线弹性阻抗模型进行滤波,得到三个射线弹性阻抗低频模型。根据地震剖面的频谱中低频缺失的范围确定数字低通滤波器的低通频率,低通频率诸如5~8hz。
S304:根据所述的地震资料确定多个入射角对应的相对弹性阻抗剖面;图4为步骤S304的实施方式一的具体流程图,图5为步骤S304的实施方式二的具体流程图。
S305:根据所述射线弹性阻抗低频模型以及所述的相对弹性阻抗剖面确定射线弹性阻抗剖面。在具体的实施方式中,在已知井的射线弹性阻抗约束下,将射线弹性阻抗低频模型和对应入射角的相对弹性阻抗剖面相加,获得三个角度的射线弹性阻抗剖面。图23为入射角为
5°的射线弹性阻抗剖面示意图,图24为入射角为15°的射线弹性阻抗剖面示意图,图25为入射角30°的射线弹性阻抗剖面示意图。
图4为步骤S304的实施方式一的具体流程图,由图4可知,在实施方式一中,该步骤具体包括:
S401:根据所述的地震资料确定叠前道集资料。在完成基本资料处理叠加之前的按道排列的地震数据即为叠前道集资料。
S402:根据所述叠前道集资料确定每个角度叠加的范围,在具体的实施方式中,根据叠前道集资料的信噪比的大小选择每个角度叠加的范围。诸如,在具体的实施例中,提取入射角为5°的角度道集叠加时,选择3°~7°的道集进行叠加;提取入射角为15°的角度道集叠加时,选择13°~17°的道集进行叠加;提取入射角为30°的角度道集叠加时,选择28°~32°的道集进行叠加。
S403:根据所述角度叠加的范围从所述叠前道集资料中提取所述射线弹性阻抗低频模型对应的角度道集叠加剖面。在具体的实施方式中,从地震叠前道集数据提取模型对应的三个入射角的角度道集叠加剖面。图20为入射角5°的角度道集叠加剖面示意图,图21为入射角15°的角度道集叠加剖面示意图,图22为入射角30°的角度道集叠加剖面示意图。
S404:对所述的角度道集叠加剖面压缩地震子波,得到反射系数剖面。在具体的实施方式中,对三个入射角的角度道集叠加剖面采用脉冲反褶积算法压缩地震子波,获得反射系数剖面。
S405:根据所述的反射系数剖面得到相对弹性阻抗剖面。在具体的实施方式中,在S403的基础上采用递推反演算法得到三个相对弹性阻抗剖面。
图5为步骤S304的实施方式二的具体流程图,由图5可知,在实施方式二中,该步骤具体包括:
S501:根据所述的地震资料确定叠前道集资料。在完成基本资料处理叠加之前的按道排列的地震数据即为叠前道集资料。
S502:根据所述叠前道集资料确定每个角度叠加的范围,在具体的实施方式中,根据叠前道集资料的信噪比的大小选择每个角度叠加的范围。诸如,在具体的实施例中,提取入射角为5°的角度道集叠加时,选择3°~7°的道集进行叠加;提取入射角为15°的角度道集叠加时,选择13°~17°的道集进行叠加;提取入射角为30°的角度道集叠加时,选择28°~32°的道集进行叠加。
S503:根据所述角度叠加的范围从所述叠前道集资料中提取所述射线弹性阻抗低频模型对应的角度道集叠加剖面。在具体的实施方式中,从地震叠前道集数据提取模型对应的三个
入射角的角度道集叠加剖面。图20为入射角5°的角度道集叠加剖面示意图,图21为入射角15°的角度道集叠加剖面示意图,图22为入射角30°的角度道集叠加剖面示意图。
S504:对所述的角度道集叠加剖面进行时差调整,通过时差调整消除三个剖面上主要目的层段反射同相轴的时差。
S505:对时差调整后的角度道集叠加剖面进行振幅匹配,通过振幅匹配方法使三个剖面上主要目的层段的总能量一致,进一步消除球面扩散等因素的影响。即对三个入射角的角度道集叠加剖面进行时差、振幅的一致性校正。
S506:对振幅匹配后的角度道集叠加剖面压缩地震子波,获得反射系数剖面。在具体的实施方式中,对三个入射角的角度道集叠加剖面采用脉冲反褶积算法压缩地震子波,获得反射系数剖面。
S507:根据所述的反射系数剖面得到相对弹性阻抗剖面。在具体的实施方式中,在S505的基础上采用递推反演算法得到三个相对弹性阻抗剖面。
图6为步骤S105的具体流程图,由图6可知,该步骤具体包括:
S601:根据多个入射角对应的射线弹性阻抗剖面确定纵横波速度比,在具体的实施方式中,当入射角为3个时,该步骤通过如下公式进行:
γ=(γ1+γ2)/2
其中,θ1、θ2、θ3为三个入射角,γ为纵横波速度比,γ1、γ2为中间参数,k为密度变化率与横波速度变化率之间的比例系数,REI(θ1)、REI(θ2)分别为入射角θ1、θ2对应的射线弹性阻抗剖面。
S602:根据所述的纵横波速度比确定纵波阻抗。在具体的实施方式中,当入射角为3个时,将横纵波速度比γ=β/α带入下面任意一个公式,求解纵波阻抗αρ。
S603:根据所述的纵横波速度比、纵波阻抗确定横波阻抗。在具体的实施方式中,将横纵波速度比γ=β/α和纵波阻抗αρ带入下式求出横波阻抗βρ。
βρ=αρ*β/α
S604:根据所述的纵横波速度比、纵波阻抗、横纵波阻抗确定储层弹性参数,所述的储层弹性参数包括泊松比、杨氏模量、剪切模量。图26为具体实施例中射线弹性阻抗参数反演得到的纵横波速度比剖面示意图,图27为具体实施例中射线弹性阻抗参数反演得到的纵波阻抗剖面示意图,图28为具体实施例中射线弹性阻抗参数反演得到的横波阻抗剖面示意图,可以看出盒8目的层116井附近纵波阻抗低、T26井附近纵波阻抗高;横波阻抗在两口井附近变化不大;段纵横波速度比剖面在S116井附近为低值,T26井附近为高值,具有含气后纵波速度下降造成的纵横波速度比低特点,与试气结果S116井4.13万方/日、T26井无产量相吻合。
如上所述,即为本发明提供的一种射线弹性参数的反演方法,是一种精确的射线弹性参数反演方案,利用多个入射角的射线弹性阻抗对射线弹性阻抗反演公式进行了进一步推导,可以准确快速的求解出弹性参数,而且不存在多解性,同时可适合于大入射角的地震道集数据,满足油气地震勘探的要求。
图7为本发明实施例提供的一种射线弹性参数的反演系统的结构框图,由图7可知,所述的系统包括:
资料采集装置100,用于采集探区的地震资料、测井资料,所述的测井资料包括深度数据、密度数据、纵波速度数据、横波速度数据。
在具体的实施例中,在探区内采集地震资料,在探区井中激发地震纵波和横波,得到纵波速度、横波速度数据、深度数据,在探区测井得到地层的密度数据。如图14所示,为具体实施例中su116井的测井曲线示意图,可见纵波速度数据、横波速度数据、密度数据、深度数据。
时间域转换装置200,用于根据所述的地震资料将所述的测井资料由深度域转换为时间域。图8为时间域转换装置200的具体结构框图。
射线弹性阻抗确定装置300,用于根据转换得到的时间域的测井资料确定井点位置的多个入射角对应的射线弹性阻抗。
在具体的实施方式中,所述的入射角可为2个至6个。下面以3个为例进行说明。在具体的实施方式中,确定出目的层段后(如图16中的虚线范围),确定射线弹性反射系数随入射角变化的类型,如图16所示的实施例中,射线弹性反射系数随入射角θ的增大而增大属于
第三类砂岩类型。因此,三个入射角可以选择为5°、15°、30°。REI(5)、REI(15)、REI(30)为入射角为5°、15°、30°对应的射线弹性阻抗。
在时间域采用下式由测井数据计算三个入射角对应的射线弹性阻抗:
公式中,θ为入射角,α为纵波速度,β为横波速度,ρ为密度,REI(θi)为入射角为θ的射线弹性阻抗,k是密度变化率与横波速度变化率之间的比例系数。表1为弹性阻抗数据表。
弹性阻抗剖面确定装置400,用于根据所述的射线弹性阻抗以及地震资料确定多个入射角对应的射线弹性阻抗剖面。图9为弹性阻抗剖面确定装置400的具体结构框图。
储层弹性参数确定装置500,用于根据多个入射角对应的射线弹性阻抗剖面确定储层弹性参数。图12为储层弹性参数确定装置500的具体结构框图。
图8为本发明实施例提供的一种射线弹性参数的反演系统中的时间域转换装置200的具体结构框图,由图8可知,时间域转换装置200具体包括:
地震剖面确定模块201,用于根据所述的地震资料确定地震剖面。
在具体的实施方式中,在探区内采集地震资料,处理得到地震剖面,如图13所示,为S087017测线的地震剖面示意图。
地震子波选择模块202,用于选择地震子波。在具体的实施方式中,选择标准地震子波或者从井旁地震剖面中提取地震子波。
地震记录合成模块203,用于根据所述的纵波速度数据、密度数据以及地震子波合成地震记录。
在具体的实施方式中,选择30hz的雷克子波,用纵波速度、密度速度结合地震子波合成地震记录。
层位标定模块204,用于用所述的地震记录对所述的地震剖面进行层位标定,得到层位标定结果、时深关系。
在具体的实施方式中,用合成的地震记录对地震剖面进行层位标定,得到层位标定结果、时深关系,确定目的层段,并拾取地震反射层位,如图15为具体实施例中用合成记录对地震剖面进行层位标定示意图,Tp7、Tp8、Tc2即为拾取的地震反射层位。
地震层位确定模块205,用于根据所述的层位标定结果解释得到地震层位;
时间域转化模块206,用于根据所述的时深关系将所述的测井资料由深度域转换为时间
域。
根据层位标定得到的时深关系,用均匀采样或者非均匀采样把测井资料由深度域转换为时间域,如表2为测井曲线深时转换及重新采样数据表,其中左侧框为深度域测井曲线数据,右侧框转换得到的时间域测井曲线数据。
图9为本发明实施例提供的一种射线弹性参数的反演系统中的弹性阻抗剖面确定装置400的具体结构框图,由图9可知,弹性阻抗剖面确定装置400具体包括:
射线弹性阻抗模型建立模块401,用于根据井点位置的多个入射角对应的射线弹性阻抗以及地震层位建立多个射线弹性阻抗模型,在具体的实施方式中,在地震反射层位的约束下,由三个入射角的射线弹性阻抗值以及地震层位分别建立三个射线弹性阻抗模型。图17为入射角5°的射线弹性阻抗模型示意图,图18为入射角15°的射线弹性阻抗模型示意图,图19为入射角30°的射线弹性阻抗模型示意图。所述的建立三个射线弹性阻抗模型采用反距离内插建模算法、克里金地质建模算法。
低频缺失范围确定模块402,用于确定所述地震剖面的频谱中低频缺失的范围;
滤波模块403,用于根据所述低频缺失的范围对所述射线弹性阻抗模型进行滤波,得到射线弹性阻抗低频模型。
在具体的实施方式中,用数字低通滤波器对三个射线弹性阻抗模型进行滤波,得到三个射线弹性阻抗低频模型。根据地震剖面的频谱中低频缺失的范围确定数字低通滤波器的低通频率,低通频率诸如5~8hz。
相对弹性阻抗剖面确定模块404,用于根据所述的地震资料确定多个入射角对应的相对弹性阻抗剖面。
射线弹性阻抗剖面确定模块405,用于根据所述射线弹性阻抗低频模型以及所述的相对弹性阻抗剖面确定射线弹性阻抗剖面。在具体的实施方式中,在已知井的射线弹性阻抗约束下,将射线弹性阻抗低频模型和对应入射角的相对弹性阻抗剖面相加,获得三个角度的射线弹性阻抗剖面。图23为入射角为5°的射线弹性阻抗剖面示意图,图24为入射角为15°的射线弹性阻抗剖面示意图,图25为入射角30°的射线弹性阻抗剖面示意图。
图10为相对弹性阻抗剖面确定模块404的实施方式一的具体结构框图,由图10可知,在实施方式一中,相对弹性阻抗剖面确定模块404具体包括:
叠前道集资料确定单元4041,用于根据所述的地震资料确定叠前道集资料。在完成基本资料处理叠加之前的按道排列的地震数据即为叠前道集资料。
范围确定单元4042,用于根据所述叠前道集资料确定每个角度叠加的范围,在具体的实施方式中,根据叠前道集资料的信噪比的大小选择每个角度叠加的范围。诸如,在具体的实
施例中,提取入射角为5°的角度道集叠加时,选择3°~7°的道集进行叠加;提取入射角为15°的角度道集叠加时,选择13°~17°的道集进行叠加;提取入射角为30°的角度道集叠加时,选择28°~32°的道集进行叠加。
角度道集叠加剖面提取单元4043,用于根据所述角度叠加的范围从所述叠前道集资料中提取所述射线弹性阻抗低频模型对应的角度道集叠加剖面。在具体的实施方式中,从地震叠前道集数据提取模型对应的三个入射角的角度道集叠加剖面。图20为入射角5°的角度道集叠加剖面示意图,图21为入射角15°的角度道集叠加剖面示意图,图22为入射角30°的角度道集叠加剖面示意图。
反射系数剖面确定单元4044,用于对所述的角度道集叠加剖面压缩地震子波,得到反射系数剖面。在具体的实施方式中,对三个入射角的角度道集叠加剖面采用脉冲反褶积算法压缩地震子波,获得反射系数剖面。
相对弹性阻抗剖面确定单元4045,用于根据所述的反射系数剖面得到相对弹性阻抗剖面。在具体的实施方式中,在反射系数剖面确定单元4043的基础上采用递推反演算法得到三个相对弹性阻抗剖面。
图11为相对弹性阻抗剖面确定模块404的实施方式二的具体结构框图,由图11可知,在实施方式二中,相对弹性阻抗剖面确定模块404还包括:
时差调整单元4046,用于对所述的角度道集叠加剖面进行时差调整,通过时差调整消除三个剖面上主要目的层段反射同相轴的时差。
振幅匹配单元4047,用于对时差调整后的角度道集叠加剖面进行振幅匹配,通过振幅匹配方法使三个剖面上主要目的层段的总能量一致,进一步消除球面扩散等因素的影响。即对三个入射角的角度道集叠加剖面进行时差、振幅的一致性校正。
图12为本发明实施例提供的一种射线弹性参数的反演系统中的储层弹性参数确定装置500的具体结构框图,由图12可知,储层弹性参数确定装置500具体包括:
纵横波速度比确定模块501,用于根据多个入射角对应的射线弹性阻抗剖面确定纵横波速度比,在具体的实施方式中,当入射角为3个时,该步骤通过如下公式进行:
γ=(γ1+γ2)/2
其中,θ1、θ2、θ3为三个入射角,γ为纵横波速度比,γ1、γ2为中间参数,k为密度变
化率与横波速度变化率之间的比例系数,REI(θ1)、REI(θ2)分别为入射角θ1、θ2对应的射线弹性阻抗剖面。
纵波阻抗确定模块502,用于根据所述的纵横波速度比确定纵波阻抗。在具体的实施方式中,当入射角为3个时,将横纵波速度比γ=β/α带入下面任意一个公式,求解纵波阻抗αρ。
横波阻抗确定模块503,用于根据所述的纵横波速度比、纵波阻抗确定横波阻抗。在具体的实施方式中,将横纵波速度比γ=β/α和纵波阻抗αρ带入下式求出横波阻抗βρ。
βρ=αρ*β/α
储层弹性参数确定模块504,用于根据所述的纵横波速度比、纵波阻抗、横纵波阻抗确定储层弹性参数,所述的储层弹性参数包括泊松比、杨氏模量、剪切模量。图26为具体实施例中射线弹性阻抗参数反演得到的纵横波速度比剖面示意图,图27为具体实施例中射线弹性阻抗参数反演得到的纵波阻抗剖面示意图,图28为具体实施例中射线弹性阻抗参数反演得到的横波阻抗剖面示意图,可以看出盒8目的层116井附近纵波阻抗低、T26井附近纵波阻抗高;横波阻抗在两口井附近变化不大;段纵横波速度比剖面在S116井附近为低值,T26井附近为高值,具有含气后纵波速度下降造成的纵横波速度比低特点,与试气结果S116井4.13万方/日、T26井无产量相吻合。
如上所述,即为本发明提供的一种射线弹性参数的反演系统,是一种精确的射线弹性参数反演方案,利用多个入射角的射线弹性阻抗对射线弹性阻抗反演公式进行了进一步推导,可以准确快速的求解出弹性参数,而且不存在多解性,同时可适合于大入射角的地震道集数据,满足油气地震勘探的要求。
下面结合具体的实施例,详细介绍本发明的技术方案。
在探区采集地震资料,处理得到地震剖面,如图13所示,为S087017测线的地震剖面示意图。在探区井中激发地震纵波和横波,得到纵波速度、横波速度数据,在探区测井得到地层密度,如图14所示,为具体实施例中su116井的测井曲线示意图,可见纵波速度数据、横
波速度数据、密度数据、深度数据。选择30hz雷克子波,用纵波速度、地层密度合成地震记录,用合成地震记录对地震剖面进行层位标定,确定目的层段,并拾取地震反射层位,如图15为具体实施例中用合成记录对地震剖面进行层位标定示意图,Tp7、Tp8、Tc2即为拾取的地震反射层位。
确定目的层段(如图16中的虚线范围),射线弹性反射系数随入射角θ的增大而增大属于第三类砂岩类型,因此,三个入射角可以选择为5°、15°、30°。
由层位标定得到的时深关系,用均匀采样或者非均匀采样把测井数据由深度域转换为时间域,表1中左侧框为深度域测井曲线数据,右侧框转换得到的时间域测井曲线数据。在时间域由测井数据计算三个入射角对应的射线弹性阻抗。在地震反射层位的约束下,由三个入射角的射线弹性阻抗值分别建立三个射线弹性阻抗模型,图17为入射角5°的射线弹性阻抗模型示意图,图18为入射角15°的射线弹性阻抗模型示意图,图19为入射角30°的射线弹性阻抗模型示意图。所述的建立三个射线弹性阻抗模型采用反距离内插建模算法、克里金地质建模算法。
用数字低通滤波器对三个射线弹性阻抗模型进行滤波,得到三个射线弹性阻抗低频模型,根据地震剖面的频谱中低频缺失的范围确定数字低通滤波器的低通频率;所述的确定数字低通滤波器的低通频率选择5~8hz。
从地震叠前道集数据提取模型对应的三个入射角的角度道集叠加剖面,图20为入射角5°的角度道集叠加剖面示意图;图21为入射角15°的角度道集叠加剖面示意图;图22为入射角30°的角度道集叠加剖面示意图,角度叠加的范围选择:提取入射角为5°的角度道集叠加时,选择3°~7°的道集进行叠加;提取入射角为15°的角度道集叠加时,选择13°~17°的道集进行叠加;提取入射角为30°的角度道集叠加时,选择28°~32°的道集进行叠加。
对三个入射角的角度道集叠加剖面进行时差、振幅的一致性校正,通过时差调整消除三个剖面上主要目的层段反射同相轴的时差,通过振幅匹配方法使三个剖面上主要目的层段的总能量一致,进一步消除球面扩散等因素的影响。对三个入射角的角度道集叠加剖面采用脉冲反褶积算法压缩地震子波,获得反射系数剖面,在此基础上采用递推反演算法得到三个相对弹性阻抗剖面;在已知井的射线弹性阻抗约束下,将射线弹性阻抗低频模型和对应入射角的相对弹性阻抗剖面相加;获得三个角度的射线弹性阻抗剖面。图23为入射角为5°的射线弹性阻抗剖面示意图,图24为入射角为15°的射线弹性阻抗剖面示意图,图25为入射角30°的射线弹性阻抗剖面示意图。
在已知三个角度的射线弹性阻抗剖面的情况下,计算纵横波速度比、解纵波阻抗、横波阻抗,由纵横波速度比、纵波阻抗、横纵波阻抗求出其他储层弹性参数,完成射线弹性
参数反演。其他储层弹性参数为泊松比、杨氏模量、剪切模量。
图26为具体实施例中射线弹性阻抗参数反演得到的纵横波速度比剖面示意图,图27为具体实施例中射线弹性阻抗参数反演得到的纵波阻抗剖面示意图,图28为具体实施例中射线弹性阻抗参数反演得到的横波阻抗剖面示意图,可以看出盒8目的层116井附近纵波阻抗低、T26井附近纵波阻抗高;横波阻抗在两口井附近变化不大;段纵横波速度比剖面在S116井附近为低值,T26井附近为高值,具有含气后纵波速度下降造成的纵横波速度比低特点,与试气结果S116井4.13万方/日、T26井无产量相吻合。
综上所述,本发明提出的一种射线弹性参数的反演方法及系统,属于油气田勘探、开发过程中的储层预测和油气检测技术,是一种精确的射线弹性参数反演方案,利用多个入射角射线弹性阻抗,对射线弹性阻抗反演公式进行了进一步推导,得到了求解纵横波速度比、纵波阻抗、横波阻抗的精确表示式,可以准确快速的求解出弹性参数,不会产生误差,而且不存在多解性,同时可适合于大入射角的地震道集数据,满足油气地震勘探的要求。
本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程,可以通过计算机程序来指令相关的硬件来完成,所述的程序可存储于一般计算机可读取存储介质中,该程序在执行时,可包括如上述各方法的实施例的流程。其中,所述的存储介质可为磁碟、光盘、只读存储记忆体(Read-Only Memory,ROM)或随机存储记忆体(Random Access Memory,RAM)等。
本领域技术人员还可以了解到本发明实施例列出的各种功能是通过硬件还是软件来实现取决于特定的应用和整个系统的设计要求。本领域技术人员可以对于每种特定的应用,可以使用各种方法实现所述的功能,但这种实现不应被理解为超出本发明实施例保护的范围。
本发明中应用了具体实施例对本发明的原理及实施方式进行了阐述,以上实施例的说明只是用于帮助理解本发明的方法及其核心思想;同时,对于本领域的一般技术人员,依据本发明的思想,在具体实施方式及应用范围上均会有改变之处,综上所述,本说明书内容不应理解为对本发明的限制。
Claims (8)
- 一种射线弹性参数的反演方法,其特征是,所述的方法具体包括:采集探区的地震资料、测井资料,所述的测井资料包括深度数据、密度数据、纵波速度数据、横波速度数据;根据所述的地震资料将所述的测井资料由深度域转换为时间域;根据转换得到的时间域的测井资料确定井点位置的多个入射角对应的射线弹性阻抗;根据所述的射线弹性阻抗、地震资料确定多个入射角对应的射线弹性阻抗剖面;根据多个入射角对应的射线弹性阻抗剖面确定储层弹性参数。
- 根据权利要求1所述的方法,其特征是,所述的入射角为2个至6个。
- 根据权利要求2所述的方法,其特征是,根据所述的地震资料将所述的测井资料由深度域转换为时间域具体包括:根据所述的地震资料确定地震剖面;选择地震子波;根据所述的纵波速度数据、密度数据以及地震子波合成地震记录;用所述的地震记录对所述的地震剖面进行层位标定,得到层位标定结果、时深关系;根据所述的层位标定结果解释得到地震层位;根据所述的时深关系将所述的测井资料由深度域转换为时间域。
- 根据权利要求3所述的方法,其特征是,根据所述的射线弹性阻抗、地震资料确定多个入射角对应的射线弹性阻抗剖面具体包括:根据井点位置的多个入射角对应的射线弹性阻抗以及地震层位建立多个射线弹性阻抗模型;确定所述地震剖面的频谱中低频缺失的范围;根据所述低频缺失的范围对所述射线弹性阻抗模型进行滤波,得到射线弹性阻抗低频模型;根据所述的地震资料确定多个入射角对应的相对弹性阻抗剖面;根据所述射线弹性阻抗低频模型以及所述的相对弹性阻抗剖面确定射线弹性阻抗剖面。
- 根据权利要求4所述的方法,其特征是,根据所述的地震资料确定多个入射角对应的相对弹性阻抗剖面具体包括:根据所述的地震资料确定叠前道集资料;根据所述叠前道集资料确定每个角度叠加的范围;根据所述角度叠加的范围从所述叠前道集资料中提取所述射线弹性阻抗低频模型对应的角度道集叠加剖面;对所述的角度道集叠加剖面压缩地震子波,得到反射系数剖面;根据所述的反射系数剖面得到相对弹性阻抗剖面。
- 根据权利要求4所述的方法,其特征是,根据所述的地震资料确定多个入射角对应的相对弹性阻抗剖面具体包括:根据所述的地震资料确定叠前道集资料;根据所述叠前道集资料确定每个角度叠加的范围;根据所述角度叠加的范围从所述叠前道集资料中提取所述射线弹性阻抗低频模型对应的角度道集叠加剖面;对所述的角度道集叠加剖面进行时差调整;对时差调整后的角度道集叠加剖面进行振幅匹配;对振幅匹配后的角度道集叠加剖面压缩地震子波,获得反射系数剖面;根据所述的反射系数剖面得到相对弹性阻抗剖面。
- 根据权利要求5或6所述的方法,其特征是,根据多个入射角对应的射线弹性阻抗剖面确定储层弹性参数具体包括:根据多个入射角对应的射线弹性阻抗剖面确定纵横波速度比;根据所述的纵横波速度比确定纵波阻抗;根据所述的纵横波速度比、纵波阻抗确定横波阻抗;根据所述的纵横波速度比、纵波阻抗、横纵波阻抗确定储层弹性参数,所述的储层弹性参数包括泊松比、杨氏模量、剪切模量。
Priority Applications (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
PCT/CN2015/086647 WO2017024523A1 (zh) | 2015-08-11 | 2015-08-11 | 一种射线弹性参数的反演方法 |
PCT/CN2015/096205 WO2017024702A1 (zh) | 2015-08-11 | 2015-12-02 | 一种射线弹性参数的反演系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
PCT/CN2015/086647 WO2017024523A1 (zh) | 2015-08-11 | 2015-08-11 | 一种射线弹性参数的反演方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
WO2017024523A1 true WO2017024523A1 (zh) | 2017-02-16 |
Family
ID=57983964
Family Applications (2)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
PCT/CN2015/086647 WO2017024523A1 (zh) | 2015-08-11 | 2015-08-11 | 一种射线弹性参数的反演方法 |
PCT/CN2015/096205 WO2017024702A1 (zh) | 2015-08-11 | 2015-12-02 | 一种射线弹性参数的反演系统 |
Family Applications After (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
PCT/CN2015/096205 WO2017024702A1 (zh) | 2015-08-11 | 2015-12-02 | 一种射线弹性参数的反演系统 |
Country Status (1)
Country | Link |
---|---|
WO (2) | WO2017024523A1 (zh) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110542924A (zh) * | 2019-09-02 | 2019-12-06 | 成都理工大学 | 一种高精度的纵横波阻抗反演方法 |
Families Citing this family (25)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112147696B (zh) * | 2019-06-27 | 2023-05-23 | 中国石油化工股份有限公司 | 基于储集性能指数的页岩储集性能地震评价方法及系统 |
CN112180442A (zh) * | 2019-07-03 | 2021-01-05 | 中国石油天然气集团有限公司 | 一种基于crp道集的岩性反演方法及系统 |
CN112230279B (zh) * | 2019-07-15 | 2024-03-01 | 中国石油天然气集团有限公司 | 增强纵波地震数据品质的方法及装置 |
CN110568489B (zh) * | 2019-07-23 | 2024-06-25 | 中国石油化工股份有限公司 | 一种块状介质的宽频反演方法 |
CN112305602B (zh) * | 2019-08-01 | 2023-02-24 | 中国石油天然气股份有限公司 | 基于叠前多属性与古地貌融合技术的碳酸岩储层预测方法 |
CN112346116A (zh) * | 2019-08-09 | 2021-02-09 | 中国石油天然气集团有限公司 | 储层预测方法及装置 |
CN112394393B (zh) * | 2019-08-13 | 2023-03-21 | 中国石油化工股份有限公司 | 一种crp道集数据体重构的方法 |
CN112649862B (zh) * | 2019-10-12 | 2024-06-18 | 中国石油化工股份有限公司 | 基于地层结构信息分离的断溶体识别方法和装置 |
CN110850473B (zh) * | 2019-11-19 | 2021-06-15 | 怀化学院 | 一种基于稀疏变换学习的地震波阻抗反演方法 |
CN113156495B (zh) * | 2020-01-07 | 2024-06-25 | 中国石油天然气集团有限公司 | 网格层析反演反射点确定方法及装置 |
CN111474576B (zh) * | 2020-04-21 | 2022-10-21 | 中海石油(中国)有限公司 | 一种保持地层结构的叠前地震道集反演初始模型构建方法 |
CN113589378B (zh) * | 2020-04-30 | 2023-04-28 | 中国石油化工股份有限公司 | 基于三维地震资料的断层封堵性评价方法 |
CN111914609B (zh) * | 2020-05-08 | 2023-12-01 | 中国石油天然气集团有限公司 | 井震联合叠前地质统计学弹性参数反演方法及装置 |
CN113759419B (zh) * | 2020-06-04 | 2024-06-18 | 中国石油化工股份有限公司 | 一种储层预测方法、装置、存储介质及电子设备 |
CN113917529B (zh) * | 2020-07-07 | 2024-06-18 | 中国石油化工股份有限公司 | 波阻抗反演方法、装置、存储介质及电子设备 |
CN113933910B (zh) * | 2020-07-13 | 2024-04-09 | 中国石油化工股份有限公司 | 物性参数同时计算的方法、系统、存储介质与电子设备 |
CN112363222A (zh) * | 2020-10-28 | 2021-02-12 | 中国石油天然气集团有限公司 | 叠后自适应宽带约束波阻抗反演方法及装置 |
CN112180464B (zh) * | 2020-11-03 | 2024-05-24 | 中国石油化工股份有限公司 | 一种储层物性的识别方法 |
CN113671565B (zh) * | 2021-09-17 | 2022-11-04 | 中国海洋石油集团有限公司 | 一种针对巨厚储层气藏开发的地震多尺度储层预测方法 |
CN114035229B (zh) * | 2021-10-26 | 2023-11-10 | 西安石油大学 | 叠前地震数据小波阈值去噪最优小波基选取方法 |
CN115993649B (zh) * | 2023-02-21 | 2024-03-19 | 中国石油大学(华东) | 基于等效方位杨氏模量的裂缝参数预测方法及系统 |
CN116299706A (zh) * | 2023-03-23 | 2023-06-23 | 福瑞升(成都)科技有限公司 | 基于叠前拟横波反射率属性的混合阻抗砂体识别方法 |
CN116609830A (zh) * | 2023-04-11 | 2023-08-18 | 福瑞升(成都)科技有限公司 | 基于砂体avo响应特征的河道砂体识别方法 |
CN116699684B (zh) * | 2023-06-09 | 2024-08-02 | 中海石油(中国)有限公司深圳分公司 | 面向储层预测的地震数据拼接方法、装置、设备及介质 |
CN117687095B (zh) * | 2023-12-14 | 2024-08-27 | 中海石油(中国)有限公司上海分公司 | 一种基于方向分量数据提取与子波谱拓展的时变反褶积方法 |
Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6058073A (en) * | 1999-03-30 | 2000-05-02 | Atlantic Richfield Company | Elastic impedance estimation for inversion of far offset seismic sections |
CN102466816A (zh) * | 2010-11-04 | 2012-05-23 | 中国石油天然气集团公司 | 一种叠前地震数据地层弹性常数参数反演的方法 |
CN104184610A (zh) * | 2013-05-24 | 2014-12-03 | 株式会社日立制作所 | 信息管理装置和信息管理方法 |
CN104516021A (zh) * | 2013-09-27 | 2015-04-15 | 中国石油天然气集团公司 | 一种同时提高解析式稳定性和精度的射线弹性参数反演方法 |
CN104597490A (zh) * | 2015-01-28 | 2015-05-06 | 中国石油大学(北京) | 基于精确Zoeppritz方程的多波AVO储层弹性参数反演方法 |
CN104614763A (zh) * | 2015-01-19 | 2015-05-13 | 中国石油大学(北京) | 基于反射率法的多波avo储层弹性参数反演方法及系统 |
Family Cites Families (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101907729B (zh) * | 2010-06-11 | 2012-09-05 | 北京诺克斯达石油科技有限公司 | 射线弹性参数反演方法 |
US9213122B2 (en) * | 2010-09-24 | 2015-12-15 | Schlumberger Technology Corporation | Single well anisotropy inversion using velocity measurements |
CN103760598A (zh) * | 2013-12-29 | 2014-04-30 | 中国石油大学(华东) | 一种射线弹性参数反演方法 |
CN104181610B (zh) * | 2014-08-07 | 2017-02-08 | 中国石油天然气集团公司 | 一种射线路径弹性反演方法以及系统 |
-
2015
- 2015-08-11 WO PCT/CN2015/086647 patent/WO2017024523A1/zh active Application Filing
- 2015-12-02 WO PCT/CN2015/096205 patent/WO2017024702A1/zh active Application Filing
Patent Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6058073A (en) * | 1999-03-30 | 2000-05-02 | Atlantic Richfield Company | Elastic impedance estimation for inversion of far offset seismic sections |
CN102466816A (zh) * | 2010-11-04 | 2012-05-23 | 中国石油天然气集团公司 | 一种叠前地震数据地层弹性常数参数反演的方法 |
CN104184610A (zh) * | 2013-05-24 | 2014-12-03 | 株式会社日立制作所 | 信息管理装置和信息管理方法 |
CN104516021A (zh) * | 2013-09-27 | 2015-04-15 | 中国石油天然气集团公司 | 一种同时提高解析式稳定性和精度的射线弹性参数反演方法 |
CN104614763A (zh) * | 2015-01-19 | 2015-05-13 | 中国石油大学(北京) | 基于反射率法的多波avo储层弹性参数反演方法及系统 |
CN104597490A (zh) * | 2015-01-28 | 2015-05-06 | 中国石油大学(北京) | 基于精确Zoeppritz方程的多波AVO储层弹性参数反演方法 |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110542924A (zh) * | 2019-09-02 | 2019-12-06 | 成都理工大学 | 一种高精度的纵横波阻抗反演方法 |
Also Published As
Publication number | Publication date |
---|---|
WO2017024702A1 (zh) | 2017-02-16 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
WO2017024523A1 (zh) | 一种射线弹性参数的反演方法 | |
CN113759424B (zh) | 基于频谱分解和机器学习的岩溶储层充填分析方法和系统 | |
AU2010236999B2 (en) | Interferometric seismic data processing | |
CN106368691B (zh) | 基于岩石物理地震信息三维异常孔隙压力预测方法 | |
Tkalčić et al. | Multistep modelling of teleseismic receiver functions combined with constraints from seismic tomography: Crustal structure beneath southeast China | |
CN110954943B (zh) | 被动源地震频率谐振勘探方法 | |
CN106094027B (zh) | 一种垂直地震剖面vsp钻前压力预测方法和系统 | |
RU2550770C1 (ru) | Способ определения геометрических характеристик трещины гидроразрыва пласта | |
WO2016008105A1 (zh) | 一种基于柯西分布的叠后波阻抗反演方法 | |
US9952341B2 (en) | Systems and methods for aligning a monitor seismic survey with a baseline seismic survey | |
CN113031068B (zh) | 一种基于反射系数精确式的基追踪叠前地震反演方法 | |
CN111722284B (zh) | 一种基于道集数据建立速度深度模型的方法 | |
CN106597545B (zh) | 一种水平裂缝地震叠前反演方法和装置 | |
CN103954992A (zh) | 一种反褶积方法及装置 | |
CN104570116A (zh) | 基于地质标志层的时差分析校正方法 | |
CN109444959A (zh) | 全频高精度层速度场建立方法 | |
CN104316966A (zh) | 一种流体识别方法及系统 | |
CN104422955A (zh) | 一种利用旅行时变化量进行各向异性参数提取的方法 | |
Cichostępski et al. | Simultaneous inversion of shallow seismic data for imaging of sulfurized carbonates | |
Schaff | Placing an upper bound on preseismic velocity changes measured by ambient noise monitoring for the 2004 M w 6.0 Parkfield earthquake (California) | |
CN113552624B (zh) | 孔隙度预测方法及装置 | |
CN102914797A (zh) | 一种获得地层各向异性系数的方法及装置 | |
CN104181610B (zh) | 一种射线路径弹性反演方法以及系统 | |
RU2011148308A (ru) | Способ комплексной обработки геофизических данных и технологическая система "литоскан" для его осуществления | |
WO2023123971A1 (zh) | 基于vsp的深度域地震剖面层位标定方法及装置 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
121 | Ep: the epo has been informed by wipo that ep was designated in this application |
Ref document number: 15900721 Country of ref document: EP Kind code of ref document: A1 |
|
NENP | Non-entry into the national phase |
Ref country code: DE |
|
122 | Ep: pct application non-entry in european phase |
Ref document number: 15900721 Country of ref document: EP Kind code of ref document: A1 |