CN110873895A - Variable grid micro-seismic reverse-time interference positioning method - Google Patents

Variable grid micro-seismic reverse-time interference positioning method Download PDF

Info

Publication number
CN110873895A
CN110873895A CN201811013189.2A CN201811013189A CN110873895A CN 110873895 A CN110873895 A CN 110873895A CN 201811013189 A CN201811013189 A CN 201811013189A CN 110873895 A CN110873895 A CN 110873895A
Authority
CN
China
Prior art keywords
grid
micro
seismic
time
reverse
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
CN201811013189.2A
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.)
China Petroleum and Chemical Corp
Sinopec Geophysical Research Institute
Original Assignee
China Petroleum and Chemical Corp
Sinopec Geophysical Research Institute
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 China Petroleum and Chemical Corp, Sinopec Geophysical Research Institute filed Critical China Petroleum and Chemical Corp
Priority to CN201811013189.2A priority Critical patent/CN110873895A/en
Publication of CN110873895A publication Critical patent/CN110873895A/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/307Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/60Analysis
    • G01V2210/63Seismic attributes, e.g. amplitude, polarity, instant phase

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Acoustics & Sound (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

The invention relates to a variable grid micro-seismic reverse-time interference positioning method, which comprises the following steps: establishing a grid speed model based on stratum speed, and dividing a global grid speed model grid and a local grid speed model grid; based on the grid velocity model, establishing a microseism wave field reverse time continuation equation by using a variable grid finite difference method, and inputting a microseism record serving as the reverse time continuation equation to reversely delay the microseism wave field; and performing interference imaging processing on the reverse time continuation micro seismic wave field based on interference imaging conditions, and picking up an interference imaging maximum value point so as to determine the position of a micro seismic source point. The method can still realize the accurate positioning of the micro-seismic source under the conditions of strong ground noise interference and slight disturbance of the velocity model, does not increase excessive calculation amount while improving the reverse continuation precision of the local wave field of the micro-seismic, and improves the calculation efficiency of the micro-seismic reverse time interference positioning in the future.

Description

Variable grid micro-seismic reverse-time interference positioning method
Technical Field
The invention relates to the field of micro-seismic source positioning research, in particular to a variable grid micro-seismic reverse-time interference positioning method.
Background
The source location of a microseism is a very important step in microseism monitoring. The positioning of the microseism event can describe the shape of the hydraulic fracturing fracture and further guide the implementation of the hydraulic fracturing construction process. The microseismic source location method is mainly derived from natural earthquake at first, and the classical Gieger location method is proposed by Gieger (1912), wherein an objective function is firstly constructed, and then the seismic source location is located by an inversion iteration method by utilizing the picked seismic signal first arrival. The joint epicenter determination method (JED) proposed by Douglas (1967) combines a plurality of seismic events to carry out positioning, and overcomes the influence of stratum velocity model errors on seismic positioning to a certain extent. Experts and scholars both at home and abroad (Dewey, 1972; King Ailanthus variegatus et al, 1993; Spence, 1980; Zhou Shi Yong et al, 1999) have developed improvements to the JED method and have evolved relative positioning algorithms. The relative positioning algorithm well solves the problem of seismic positioning errors caused by stratum model errors. However, the positioning method based on the first arrival of the seismic signals is susceptible to the signal-to-noise ratio of the seismic signals in the ground micro-seismic positioning because the micro-seismic magnitude is small and the signal-to-noise ratio is low.
In recent years, many experts and scholars have proposed microseismic waveform-based positioning methods based on this problem. Wangwave et al (2012) propose a ground microseism emission tomography (SET) positioning algorithm, which utilizes travel time information to shift ground microseism signals with low signal-to-noise ratio and perform multi-channel correlation processing, thereby reducing the interference of noise to ground monitoring microseism positioning. The method is a micro-seismic positioning method based on a ray theory, and due to the limitation of the ray theory, the method only utilizes travel time information of seismic signals, and the micro-seismic signals are not comprehensively utilized.
Based on the above problems, in recent years, many experts and scholars have developed the research of microseism positioning method based on wave theory, and the method is mainly from the popular reverse time migration theory at present. Gajewski et al (2005,2007) used for reference to the reverse time continuation method in the reverse time migration imaging theory to locate the micro-seismic source by recording the reverse time reverse continuation wave field by the micro-seismic. Saenger (2010) uses the acoustic wave equation to perform a reverse time localization of microseismic events. Sava (2008,2011) suppresses wave field crosstalk in a micro-seismic reverse time-delay wave field by utilizing an interference imaging condition, and can be well suitable for micro-seismic positioning under a sparse observation condition. Wanchenlong et al (2013) and li zhen et al (2014) conducted reverse time localization studies of microseismic events under combined well and surface observations. The method mainly applies the finite difference method to carry out the reverse continuation of the microseism wave field, and because the calculation efficiency of the finite difference method is mainly influenced by the size of the global grid division, the method inevitably reduces the calculation efficiency due to the fine division of the global grid when the precision of the reverse continuation of the wave field is ensured.
In summary, most of the existing microseism positioning method researches are based on a wave theory, and a finite difference method is mainly used for carrying out reverse continuation on a microseism wave field to realize positioning of a microseism seismic source. The method cannot simultaneously ensure the reverse continuation precision and the calculation efficiency of the microseism wave field.
Disclosure of Invention
In view of the above, the main objective of the present invention is to provide a variable-grid finite difference method-based microseism reverse-time interferometric positioning method, so as to solve the problem that the prior art cannot simultaneously ensure the reverse continuation precision and the calculation efficiency of a microseism wave field.
The purpose of the invention is realized by the following technical scheme:
a variable grid micro-seismic reverse-time interference positioning method comprises the following steps:
establishing a grid speed model based on stratum speed, and dividing a global grid speed model grid and a local grid speed model grid;
based on the grid velocity model, establishing a microseism wave field reverse time continuation equation by using a variable grid finite difference method, and inputting a microseism record serving as the reverse time continuation equation to reversely delay the microseism wave field;
and performing interference imaging processing on the reversely time-delayed microseism wave field based on interference imaging conditions, and picking up an interference imaging maximum value point so as to determine the position of the microseism seismic source point.
Preferably, a grid velocity model is established based on the stratum velocity, and a global grid velocity model grid and a local grid velocity model grid are divided, specifically comprising the following steps:
establishing a global grid region velocity model based on the stratum velocity structure, taking the x direction of the velocity structure as an example, setting the grid interval of the global grid region in the x direction as delta x, and completing the establishment of the global grid region velocity model in the x direction;
according to the requirement of a variable coefficient algorithm, the grid interval of the global grid region in the x direction is any odd number multiple b of the local grid interval delta x, the grid interval of the local grid region is set to delta x/b, and the establishment of a local grid region speed model in the x direction is completed according to a stratum speed structure;
establishing a global grid region velocity model based on the stratum velocity structure, taking a velocity structure z direction as an example, setting grid intervals in the z direction of a global grid region as delta z, and completing the establishment of the global grid region velocity model in the z direction;
according to the requirement of a variable coefficient algorithm, the grid interval of the global grid region in the z direction is any odd number b times of the local grid interval delta z, the grid interval of the local grid region is set to delta z/b, and the establishment of a local grid region speed model in the z direction is completed according to a stratum speed structure.
Preferably, the x-direction is a direction extending along a surface of the formation and the z-direction is a direction extending along a depth of the formation.
Preferably, based on the grid velocity model, a variable grid finite difference method is used for establishing a microseism wave field reverse time continuation equation, a microseism record is used as the input of the reverse time continuation equation, and the microseism wave field is reversely delayed, and the method specifically comprises the following steps:
taking the differential calculation of the wave field u as an example to illustrate the partial derivative of the wave field u with respect to x, first assuming that u (x) has a derivative of order 2N +1, then x is set as:
Figure BDA0001785562700000031
in the above formula, Δ x represents a discrete interval in the x direction, and m is x and x0The value range of m is 1,2, …, N, then u (x) at x the taylor expansion of 2N +1 order is:
Figure BDA0001785562700000032
under a varying grid, the 2N order differential formula of the first derivative at x of the finite difference method is represented as:
Figure BDA0001785562700000033
the formula is simplified to obtain:
Figure BDA0001785562700000034
undetermined coefficient a in the above formulamThe calculation is performed by the following matrix equation:
Figure BDA0001785562700000035
calculating the matrix equation to obtain undetermined coefficient am
When N is equal to 1, find a1=1,
When N is 2, find a1=1.125,a2=-0,04166667,
When N is present>2, undetermined coefficient amThe solution of (a) is given by:
Figure BDA0001785562700000041
in the case of the two-dimensional acoustic wave equation, the first order velocity stress equation of the two-dimensional acoustic wave equation is expressed in the form:
Figure BDA0001785562700000042
in the above formula, vxAnd vzThe velocities of the subsurface particles in the x-direction and z-direction, p the normal stress of the particles, c the medium velocity, ρ the medium density,
carrying out difference calculation on the formula to obtain a variable-grid high-order finite difference format with second-order difference precision in time and any-order difference precision in space, carrying out variable-grid finite difference on an acoustic wave equation, and setting stress, a mass point velocity transverse component and a mass point velocity longitudinal component on points (i, j, k) of space grid points as
Figure BDA0001785562700000043
Figure BDA0001785562700000044
The discrete differential representation of the particle stress, the particle velocity transverse component and the longitudinal component is expressed as:
Figure BDA0001785562700000045
Figure BDA0001785562700000046
Figure BDA0001785562700000047
in the above formula, anThe difference coefficient of the variable grid finite difference method is represented, wherein delta x is a space grid interval in the x direction, delta z is a space grid interval in the z direction, and delta t is a time step interval;
and taking the micro-seismic record as input, and performing reverse time continuation on the micro-seismic wave field by using the formula to obtain the micro-seismic wave field at each moment.
Preferably, the medium density ρ is 1g/cm3 in the acoustic wave equation.
Preferably, the method comprises the following steps of performing interference imaging processing on the reversely time-delayed microseism wave field based on interference imaging conditions, and picking up an interference imaging maximum value point to determine the position of a microseism seismic source point:
WDF filtering is carried out on each space grid point in the reverse time continuation wave field of the global and local grid areas near the seismic source excitation time according to the interference positioning imaging condition to obtain an interference imaging section, wherein the WDF filtering formula is shown as follows,
Figure BDA0001785562700000051
in the above equation, Y represents the spatial position of an imaging point, T represents the imaging time, W represents the interference imaging result after WDF filtering, V represents the reverse time-lapse wave field, T represents the size of the time window in the interference condition, Y represents the size of the spatial window in the interference condition, ThRepresenting any two symmetric wavefield time intervals, y, around the excitation timehRepresenting any two symmetrical space point space intervals near the imaging point;
searching a grid point with the maximum numerical value on the interference imaging profile to obtain the spatial position of the grid point with the maximum numerical value; and the position is used as the position of the micro seismic source to finish the positioning of the micro seismic source.
Preferably, the interference imaging profile used to determine the location of the microseismic source is an interference imaging profile of the local grid area.
Compared with the prior art, the invention has the following advantages or beneficial effects:
the method can still realize the accurate positioning of the micro-seismic source under the conditions of strong ground noise interference and slight disturbance of the velocity model, does not increase excessive calculation amount while improving the reverse continuation precision of the local wave field of the micro-seismic, and improves the calculation efficiency of the reverse time interference positioning of the micro-seismic in the future.
Additional advantages, objects, and features of the invention will be set forth in part in the description which follows and in part will become apparent to those having ordinary skill in the art upon examination of the following or may be learned from practice of the invention. The objectives and other advantages of the invention will be realized and attained by the structure particularly pointed out in the written description and claims hereof as well as the appended drawings.
Drawings
The accompanying drawings are included to provide a further understanding of the technology or prior art of the present application and are incorporated in and constitute a part of this specification. The drawings expressing the embodiments of the present application are used for explaining the technical solutions of the present application, and should not be construed as limiting the technical solutions of the present application.
FIG. 1 is a flowchart illustrating a first embodiment of the present invention;
FIG. 2(a) is a diagram of a deep anticline velocity model in accordance with a second embodiment of the present invention;
FIG. 2(b) is a diagram of a grid triple-encryption-position velocity model in the second embodiment of the present invention;
FIG. 3(a) is a non-interfering microseismic ground record diagram according to the second embodiment of the present invention;
FIG. 3(b) is a recording diagram of a non-interfering microseismic well in accordance with a second embodiment of the present invention;
FIG. 4(a) is a cross-sectional view of the positioning of the non-interfering microseismic global grid region in the second embodiment of the present invention;
FIG. 4(b) is a cross-sectional view of the location of the non-interfering microseismic local grid area in the second embodiment of the present invention;
FIG. 5(a) is a ground noise microseismic ground record diagram according to a second embodiment of the present invention;
FIG. 5(b) is a plot of a surface noisy microseismic well in accordance with a second embodiment of the present invention;
FIG. 6(a) is a section view of a ground noise-added microseismic global grid region reverse time extension in accordance with a second embodiment of the present invention;
FIG. 6(b) is a sectional view of a local grid region of a ground noise-added microseism according to a second embodiment of the present invention;
FIG. 7(a) is a sectional view of a ground noise-added microseismic recording reverse time interferometric global grid area location in accordance with a second embodiment of the present invention;
FIG. 7(b) is a sectional view of a ground noise-added microseismic recording positioning of a local grid area of reverse time interference in accordance with a second embodiment of the present invention;
FIG. 8(a) is a model diagram of the velocity of the back-tilt disturbance at the noisy deep portion in the second embodiment of the present invention;
FIG. 8(b) is a diagram of a disturbance velocity model at a triple-encryption position of a grid according to a second embodiment of the present invention;
fig. 9(a) is a cross-sectional view of the global grid area extension when the micro-seismic record after ground noise is subjected to micro-seismic reverse time extension under the condition of velocity model disturbance in the second embodiment of the present invention;
fig. 9(b) is a sectional view of a local grid area continuation when the micro-seismic record after ground noise is subjected to micro-seismic reverse time continuation under the condition of velocity model disturbance in the second embodiment of the present invention;
FIG. 10(a) is a cross-sectional view of the global grid region after applying interferometric imaging conditions to the microseismic reverse time-lapse wavefield in accordance with a second embodiment of the present invention;
FIG. 10(b) is a cross-sectional view of the local grid region after applying interference imaging conditions to the microseismic reverse time-lapse wavefield in accordance with a second embodiment of the present invention.
Detailed Description
The following detailed description of the embodiments of the present invention will be provided with reference to the accompanying drawings and examples, so that how to apply the technical means to solve the technical problems and achieve the corresponding technical effects can be fully understood and implemented. The embodiments and the features of the embodiments can be combined without conflict, and the technical solutions formed are all within the scope of the present invention.
In the following description, for purposes of explanation, numerous specific details are set forth in order to provide a thorough understanding of the embodiments of the invention. It will be apparent, however, to one skilled in the art that the present invention may be practiced without some of these specific details or with other methods described herein.
The first embodiment is as follows:
the invention provides a variable grid micro-seismic reverse-time interference positioning method, which comprises the following steps:
the method comprises the following steps: and establishing a grid speed model based on the stratum speed, and dividing a global grid speed model grid and a local grid speed model grid.
Further, taking the x direction of the speed structure as an example, the specific steps of dividing the global mesh speed model mesh are as follows: and setting the grid interval of the global grid region in the x direction as delta x, and finishing the establishment of the global grid region speed model in the x direction. The specific steps for dividing the local grid speed model grid are as follows: according to the requirement of a variable coefficient algorithm, the grid interval of the global grid region in the x direction is any odd number multiple b of the local grid interval delta x, the grid interval of the local grid region is set to delta x/b, and the establishment of a local grid region speed model in the x direction is completed according to a stratum speed structure.
Taking the speed structure z direction as an example, the specific steps of dividing the global grid speed model grid are as follows: and setting the grid interval in the z direction of the global grid area as delta z, and finishing the establishment of the speed model of the global grid area in the z direction. The specific steps for dividing the local grid speed model grid are as follows: according to the requirement of a variable coefficient algorithm, the grid interval of the global grid region in the z direction is any odd number b times of the local grid interval delta z, the grid interval of the local grid region is set to delta z/b, and the establishment of a local grid region speed model in the z direction is completed according to a stratum speed structure.
Where the x-direction is the direction extending along the surface of the formation and the z-direction is the direction extending along the depth of the formation.
Step two: and based on the grid velocity model established in the step S1, establishing a microseism wave field reverse time continuation equation by using a variable grid finite difference method, and inputting the microseism record serving as the reverse time continuation equation to reversely delay the microseism wave field. The method specifically comprises the following steps:
taking the differential calculation of the wave field u as an example to illustrate the partial derivative of the wave field u with respect to x, first assuming that u (x) has a derivative of order 2N +1, then x is set as:
Figure BDA0001785562700000071
in the above formula (1), Δ x represents a discrete interval in the x direction, and m is x and x0The value range of m is 1,2, …, N, then u (x) at x the taylor expansion of 2N +1 order is:
Figure BDA0001785562700000081
in equation (2), where m is 1,2, …, N, the 2N order differential equation of the first derivative at x of the finite difference method under the variable grid is expressed as:
Figure BDA0001785562700000082
the formula (3) is simplified to obtain:
Figure BDA0001785562700000083
undetermined coefficient a in expressions (3) and (4)mThe calculation is performed by the following matrix equation:
Figure BDA0001785562700000084
calculating the matrix equation (5) to obtain a undetermined coefficient am
When N is equal to 1, find a1=1,
When N is 2, find a1=1.125,a2=-0,04166667,
When N is present>2, undetermined coefficient amIs obtained by the following equation (6):
Figure BDA0001785562700000085
in the case of the two-dimensional acoustic wave equation, the first order velocity stress equation of the two-dimensional acoustic wave equation is expressed in the form:
Figure BDA0001785562700000086
in the above formula (7), vxAnd vzThe velocity of underground particles in the x direction and the z direction is shown, p is the normal stress of the particles, c is the medium velocity, rho is the medium density, and the medium density rho is 1g/cm in the acoustic wave equation3
Carrying out differential calculation on the formula (7) to obtain the second-order differential precision in time and the variable grid high-order with any order differential precision in spaceAnd in a limited difference cellular mode, variable-grid finite difference is carried out on the sound wave equation. Let stress, the transverse component of particle velocity, and the longitudinal component of particle velocity at points (i, j, k) of the spatial grid points be
Figure BDA0001785562700000091
The discrete differential representation of the particle stress, the particle velocity transverse component and the longitudinal component is expressed as:
Figure BDA0001785562700000092
Figure BDA0001785562700000093
Figure BDA0001785562700000094
in the above formulae (8), (9) and (10), anAnd the difference coefficient represents a variable grid finite difference method, wherein delta x is a space grid interval in an x direction, delta z is a space grid interval in a z direction, delta t is a time step interval, p is the normal stress of a particle point, c is a medium speed, and rho is a medium density.
And then, taking the micro-seismic record as input, and performing reverse time continuation on the micro-seismic wave field by using the above formulas (8), (9) and (10) to obtain the micro-seismic wave field at each moment.
Step three: and performing interference imaging processing on the reversely time-delayed microseism wave field based on interference imaging conditions, and picking up an interference imaging maximum value point so as to determine the position of the microseism seismic source point.
Specifically, WDF filtering is carried out on each spatial grid point in the reverse time continuation wave field of the global and local grid areas near the seismic source excitation time according to the interference positioning imaging condition to obtain an interference imaging section, the WDF filtering formula is shown as follows,
Figure BDA0001785562700000095
in the above formulaY denotes a spatial position of an imaging point, T denotes an imaging time, W denotes an interference imaging result after WDF filtering processing, V denotes a reverse time-lapse wave field, T denotes a size of a time window in an interference condition, Y denotes a size of a space window in the interference condition, T denotes a size of a time window in the interference condition, andhrepresenting any two symmetric wavefield time intervals, y, around the excitation timehRepresenting any two symmetric spatial point spatial separations in the vicinity of the imaging point.
Searching a grid point with the maximum numerical value on the interference imaging profile to obtain the spatial position of the grid point with the maximum numerical value; and the position is used as the position of the micro seismic source to finish the positioning of the micro seismic source.
In addition, since the wavefield delineation for the local grid regions is finer, the microseismic source locations are typically determined by interferometric imaging profiles of the local grid regions. Therefore, in the present embodiment, the interferometric imaging profile used to determine the location of the microseismic source is preferably an interferometric imaging profile of a local grid area.
Example two
The following operations are carried out in combination with the steps in the first embodiment:
firstly, model trial calculation is carried out under the condition of no interference, a deep anticline speed model shown in a part a of figure 2 is established, the transverse range of the model is 0-2000m, the longitudinal range is 0-2000m, the size of a speed model grid is 10m in the transverse direction, the longitudinal size of the speed model grid is 10m, 3 times of encryption is carried out on the area range grid of 1200m in the transverse direction and 1900m in the depth direction, the encryption area speed model is shown in a part b of figure 2, and the density is set to be 1g/cm 3. With coordinates (0,0) as a base point, 101 ground detectors are respectively arranged at intervals of 20m along the positive direction of the transverse coordinate axis, and 11-stage borehole detectors are arranged at intervals of 10m in a borehole (left side of the model) from 1600m to 1700 m. The source point location is loaded at coordinates (1004,1704) with a source firing time of 0 s. Setting the time sampling interval to be 0.1ms, and selecting a seismic source wavelet with the wavelet frequency of 30Hz to perform micro-seismic recording variable-staggered grid acoustic forward numerical simulation.
As shown in FIG. 3, the non-interfering microseismic recordings are recorded in the surface in section a of FIG. 3 and in the well in section b of FIG. 3 (0-1S). Since the interference imaging processing is performed under the condition of no interference, the noise interference and the wave field crosstalk are not obvious, so the interference imaging processing is not performed under the interference imaging condition.
As shown in fig. 4, the section a of fig. 4 is a global grid area positioning section, and the section b of fig. 4 is a local grid area positioning section. And (1003.3,1703.3) extracting positioning coordinates for the microseism reverse time positioning local grid area positioning section. It can be seen that the microseism reverse-time positioning result under the theoretical model is very accurate.
Next, model trial calculation is performed to simulate strong ground noise, and noise is added to the acoustic forward modeling numerical simulation ground micro-seismic record to make the signal-to-noise ratio of the micro-seismic record be-1, the micro-seismic record is shown in fig. 5, wherein part a of fig. 5 is the ground record, and part b of fig. 5 is the well record (0-1S).
The microseismic reverse time continuation profile at the seismic source excitation time (ground plus noise) is shown in fig. 6. It can be seen that the noise has little influence on the microseism reverse-time positioning result. Where part a of fig. 6 is a global mesh region continuation section, and part b of fig. 6 is a local mesh region continuation section.
An interference imaging condition (ground plus noise) was applied to the wavefield and the microseismic reverse time localization profile after applying the interference condition is shown in figure 7. Where part a of fig. 7 is a global mesh region localization section and part b of fig. 7 is a local mesh region localization section. As can be seen from FIG. 7, the background noise of the reverse time positioning section of the micro-earthquake is effectively suppressed after the interference imaging condition is applied, and the focusing performance at the seismic source is greatly improved. The local grid area location profile shown in part b of fig. 7 is extracted with location coordinates, and as a result (1003.3,1703.3), the location result is still accurate.
And finally, on the basis of the noise addition of the ground micro-seismic record, adding noise to the velocity model in a disturbance mode, wherein the disturbance range is 5%, and the velocity model after the noise addition and disturbance is shown in a figure 8. Part a of fig. 8 is a deep anticline disturbance velocity model, and part b of fig. 8 is a grid triple-encryption-point disturbance velocity model.
And (3) carrying out micro-seismic reverse time continuation processing on the micro-seismic records subjected to ground noise under the condition of velocity model disturbance, wherein the micro-seismic reverse time continuation section at the seismic source excitation time is shown in figure 9. Part a of fig. 9 is a global mesh region continuation section, and part b of fig. 9 is a local mesh region continuation section.
An interference imaging condition is applied to the micro-seismic reverse-time continuation wave field, and a micro-seismic reverse-time positioning section after the interference condition is applied is shown in fig. 10. Part a of fig. 10 is a global mesh region localization section, and part b of fig. 10 is a local mesh region localization section. The local grid area positioning profile shown in part b of fig. 10 is subjected to positioning coordinate extraction, and the result is (1003.3,1703.3), and the positioning result shows that in the case of slight velocity model disturbance, the positioning precision of the method is not greatly influenced and is still more accurate.
In conclusion, the method can still realize the accurate positioning of the micro-seismic source under the conditions of strong ground noise interference and slight disturbance of the velocity model, does not increase excessive calculation amount while improving the reverse continuation precision of the local wave field of the micro-seismic, and improves the calculation efficiency of the reverse time interference positioning of the micro-seismic in the future.
The above description is only a preferred embodiment of the present invention, but the scope of the present invention is not limited thereto, and any changes or substitutions that can be easily conceived by those skilled in the art within the technical scope of the present invention are included in the scope of the present invention. Therefore, the protection scope of the present invention shall be subject to the protection scope of the claims.

Claims (7)

1. A variable grid micro-seismic reverse-time interference positioning method comprises the following steps:
establishing a grid speed model based on stratum speed, and dividing a global grid speed model grid and a local grid speed model grid;
based on the grid velocity model, establishing a microseism wave field reverse time continuation equation by using a variable grid finite difference method, and inputting a microseism record serving as the reverse time continuation equation to reversely delay the microseism wave field;
and performing interference imaging processing on the reversely time-delayed microseism wave field based on interference imaging conditions, and picking up an interference imaging maximum value point so as to determine the position of the microseism seismic source point.
2. The variable-grid micro-seismic reverse-time interferometric positioning method of claim 1, characterized in that a grid velocity model is established based on formation velocity, and a global grid velocity model grid and a local grid velocity model grid are divided, specifically comprising the following steps:
establishing a global grid region velocity model based on the stratum velocity structure, taking the x direction of the velocity structure as an example, setting the grid interval of the global grid region in the x direction as delta x, and completing the establishment of the global grid region velocity model in the x direction;
according to the requirement of a variable coefficient algorithm, the grid interval of the global grid region in the x direction is any odd number multiple b of the local grid interval delta x, the grid interval of the local grid region is set to delta x/b, and the establishment of a local grid region speed model in the x direction is completed according to a stratum speed structure;
establishing a global grid region velocity model based on the stratum velocity structure, taking a velocity structure z direction as an example, setting grid intervals in the z direction of a global grid region as delta z, and completing the establishment of the global grid region velocity model in the z direction;
according to the requirement of a variable coefficient algorithm, the grid interval of the global grid region in the z direction is any odd number b times of the local grid interval delta z, the grid interval of the local grid region is set to delta z/b, and the establishment of a local grid region speed model in the z direction is completed according to a stratum speed structure.
3. The variable-grid microseismic reverse-time interferometric positioning method of claim 2 wherein the x-direction is a direction extending along the surface of the formation and the z-direction is a direction extending along the depth of the formation.
4. The variable-grid micro-seismic reverse-time interference positioning method according to claim 2 or 3, wherein a micro-seismic wave field reverse time continuation equation is established by using a variable-grid finite difference method based on the grid velocity model, a micro-seismic record is input as the reverse time continuation equation, and a micro-seismic wave field is reversely time-delayed, and the method specifically comprises the following steps:
taking the differential calculation of the wave field u as an example to illustrate the partial derivative of the wave field u with respect to x, first assuming that u (x) has a derivative of order 2N +1, then x is set as:
Figure FDA0001785562690000011
in the above formula, Δ x represents a discrete interval in the x direction, and m is x and x0The value range of m is 1,2, …, N, then u (x) at x the taylor expansion of 2N +1 order is:
Figure FDA0001785562690000021
under a varying grid, the 2N order differential formula of the first derivative at x of the finite difference method is represented as:
Figure FDA0001785562690000022
the formula is simplified to obtain:
Figure FDA0001785562690000023
undetermined coefficient a in the above formulamThe calculation is performed by the following matrix equation:
Figure FDA0001785562690000024
calculating the matrix equation to obtain undetermined coefficient am
When N is equal to 1, find a1=1,
When N is 2, find a1=1.125,a2=-0,04166667,
When N is present>2, undetermined coefficient amThe solution of (a) is given by:
Figure FDA0001785562690000025
in the case of the two-dimensional acoustic wave equation, the first order velocity stress equation of the two-dimensional acoustic wave equation is expressed in the form:
Figure FDA0001785562690000026
in the above formula, vxAnd vzThe velocities of the subsurface particles in the x-direction and z-direction, p the normal stress of the particles, c the medium velocity, ρ the medium density,
carrying out difference calculation on the formula to obtain a variable-grid high-order finite difference format with second-order difference precision in time and any-order difference precision in space, carrying out variable-grid finite difference on an acoustic wave equation, and setting stress, a mass point velocity transverse component and a mass point velocity longitudinal component on points (i, j, k) of space grid points as
Figure FDA0001785562690000031
Figure FDA0001785562690000032
The discrete differential representation of the particle stress, the particle velocity transverse component and the longitudinal component is expressed as:
Figure FDA0001785562690000033
Figure FDA0001785562690000034
Figure FDA0001785562690000035
in the above formula, anDifferential coefficient representing variable-mesh finite difference methodΔ x is a spatial grid interval in the x direction, Δ z is a spatial grid interval in the z direction, and Δ t is a time step interval;
and taking the micro-seismic record as input, and performing reverse time continuation on the micro-seismic wave field by using the formula to obtain the micro-seismic wave field at each moment.
5. The variable-grid micro-seismic reverse-time interferometric positioning method of claim 4, characterized in that:
the density rho of the medium is 1g/cm in an acoustic wave equation3
6. The variable-grid micro-seismic reverse-time interference positioning method according to claim 1, wherein the interference imaging processing is performed on the reversely delayed micro-seismic wave field based on the interference imaging condition, and an interference imaging maximum value point is picked up, so as to determine the position of the micro-seismic source point, and the method specifically comprises the following steps:
WDF filtering is carried out on each space grid point in the reverse time continuation wave field of the global and local grid areas near the seismic source excitation time according to the interference positioning imaging condition to obtain an interference imaging section, wherein the WDF filtering formula is shown as follows,
Figure FDA0001785562690000036
in the above equation, Y represents the spatial position of an imaging point, T represents the imaging time, W represents the interference imaging result after WDF filtering, V represents the reverse time-lapse wave field, T represents the size of the time window in the interference condition, Y represents the size of the spatial window in the interference condition, ThRepresenting any two symmetric wavefield time intervals, y, around the excitation timehRepresenting any two symmetrical space point space intervals near the imaging point;
searching a grid point with the maximum numerical value on the interference imaging profile to obtain the spatial position of the grid point with the maximum numerical value; and the position is used as the position of the micro seismic source to finish the positioning of the micro seismic source.
7. The variable-grid micro-seismic reverse-time interferometric positioning method of claim 6, characterized in that the interferometric imaging profile is a local grid area interferometric imaging profile.
CN201811013189.2A 2018-08-31 2018-08-31 Variable grid micro-seismic reverse-time interference positioning method Pending CN110873895A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811013189.2A CN110873895A (en) 2018-08-31 2018-08-31 Variable grid micro-seismic reverse-time interference positioning method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201811013189.2A CN110873895A (en) 2018-08-31 2018-08-31 Variable grid micro-seismic reverse-time interference positioning method

Publications (1)

Publication Number Publication Date
CN110873895A true CN110873895A (en) 2020-03-10

Family

ID=69715869

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201811013189.2A Pending CN110873895A (en) 2018-08-31 2018-08-31 Variable grid micro-seismic reverse-time interference positioning method

Country Status (1)

Country Link
CN (1) CN110873895A (en)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112114362A (en) * 2020-09-07 2020-12-22 中北大学 Method for reconstructing space-time field of underground shallow layer explosion
CN114089419A (en) * 2020-08-24 2022-02-25 中国石油化工集团有限公司 Optimized variable grid earthquake forward modeling method

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20150134308A1 (en) * 2012-09-14 2015-05-14 Institute Of Geology And Geophysics, Chinese Academy Of Sciences Method and device for acquiring optimization coefficient, and related method and device for simulating wave field
CN105652320B (en) * 2015-12-30 2018-05-04 中国石油天然气集团公司 Reverse-time migration imaging method and device

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20150134308A1 (en) * 2012-09-14 2015-05-14 Institute Of Geology And Geophysics, Chinese Academy Of Sciences Method and device for acquiring optimization coefficient, and related method and device for simulating wave field
CN105652320B (en) * 2015-12-30 2018-05-04 中国石油天然气集团公司 Reverse-time migration imaging method and device

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
DIRK GAJEWSKI AND EKKEHART TESSMER: "Reversemodelling for seismic event characterization", 《GEOPHYS. J. INT.》 *
周德山: "基于变网格逆时偏移的微地震定位", 《中国优秀硕士学位论文全文数据库》 *
朱生旺,等: "变网格有限差分弹性波方程数值模拟方法", 《石油地球物理勘探》 *
李振春,等: "井地联合观测多分量微地震逆时干涉定位", 《石油地球物理勘探》 *
王晨龙,等: "地面与井中观测条件下的微地震干涉逆时定位算法", 《地球物理学报》 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114089419A (en) * 2020-08-24 2022-02-25 中国石油化工集团有限公司 Optimized variable grid earthquake forward modeling method
CN114089419B (en) * 2020-08-24 2024-04-30 中国石油化工集团有限公司 Optimized variable grid earthquake forward modeling method
CN112114362A (en) * 2020-09-07 2020-12-22 中北大学 Method for reconstructing space-time field of underground shallow layer explosion

Similar Documents

Publication Publication Date Title
CN107505654B (en) Full waveform inversion method based on earthquake record integral
CN108345031B (en) Full waveform inversion method for elastic medium active source and passive source mixed mining seismic data
CN106353792B (en) Method suitable for positioning micro-seismic source of hydraulic fracturing
CN101334483B (en) Method for attenuating rayleigh wave scattered noise in earthquake data-handling
CN110133715B (en) Microseism seismic source positioning method based on first-arrival time difference and waveform superposition
CN111158049B (en) Seismic reverse time migration imaging method based on scattering integration method
CN111045077B (en) Full waveform inversion method of land seismic data
CN111290017B (en) Surface wave exploration method for jointly extracting Rayleigh wave frequency dispersion characteristics through seismic electric wave field
CN110579795B (en) Joint velocity inversion method based on passive source seismic waveform and reverse-time imaging thereof
CN105388520B (en) Seismic data prestack reverse time migration imaging method
CN108845351A (en) A kind of VSP seismic data converted wave full waveform inversion method
CN105549068A (en) 3D anisotropic micro seismic interference inverse positioning method and 3D anisotropic micro seismic interference inverse positioning system
CN109669212A (en) Seismic data processing technique, interval quality factors evaluation method and device
CN103758511A (en) Method and device for recognizing hidden reservoirs by underground reverse-time migration imaging
CN103926623A (en) Method for suppressing reverse time migration low frequency noise
CN104570116A (en) Geological marker bed-based time difference analyzing and correcting method
CN104330826A (en) A method for removing various noises under the condition of complex surface
CN105652319A (en) Estimation method of near-surface stratum Q value of complex mediums
CN111781639B (en) Shot-geophone reciprocal elastic wave full waveform inversion method for OBS multi-component data
CN110873895A (en) Variable grid micro-seismic reverse-time interference positioning method
CN113885079A (en) Elastic wave field decoupling-based high-precision multi-azimuth reverse-time seismic source imaging method
CN104730576A (en) Curvelet transform-based denoising method of seismic signals
CN106443793A (en) Space-time bivariant forward modeling method
US11754744B2 (en) Surface wave prospecting method for jointly extracting Rayleigh wave frequency dispersion characteristics by seismoelectric field
CN105182414A (en) Direct wave removing method based on wave equation forward modeling

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

Application publication date: 20200310

WD01 Invention patent application deemed withdrawn after publication