SURFACE WAVE PROSPECTING METHOD FOR JOINTLY EXTRACTING RAYLEIGH WAVE FREQUENCY DISPERSION CHARACTERISTICS BY SEISMOELECTRIC FIELD

Abstract

A surface wave prospecting method for jointly extracting Rayleigh wave frequency dispersion characteristics in a seismoelectric field. A surface wave prospecting method includes following steps of: acquiring jointly acquired data, where the jointly acquired data includes seismic wave data and electric field data; carrying out jointly imaging processing on jointly acquired data to obtain a superposed frequency dispersion spectrum; carrying out extraction processing on superposed frequency dispersion spectrum to obtain a frequency dispersion curve, outperforming inversion processing on frequency dispersion curve to obtain a stratum structure profile. As seismic wave data and electric field data are adopted to carry out combined imaging processing to obtain superposed frequency dispersion spectrum, multi-mode frequency dispersion curve is extracted, multiplicity of solutions of inversion is greatly reduced during inversion, precision and stability of surface wave prospecting are greatly improved.

Claims

1-10. (canceled)

11. A surface wave prospecting method for jointly extracting Rayleigh wave frequency dispersion characteristics by a seismoelectric field, the method comprising: acquiring a plurality of jointly collected data, wherein the jointly collected data comprises a seismic wave data and an electric field data; performing a jointly imaging processing onto the jointly collected data to obtain a superposed dispersion spectrum; performing an extraction processing onto the superposed dispersion spectrum to obtain a dispersion curve; and performing an inversion processing onto the dispersion curve to obtain a stratum structure profile.

12. The method according to claim 11, wherein the seismic wave data comprises a seismic wave component; the electric field data composes an electric field component; the performing of jointly imaging processing onto the jointly collected data to obtain the superposed dispersion spectrum, comprises: obtaining a dispersion spectrum of the seismic wave component and a dispersion spectrum of the electric field component according to the seismic wave component and the electric field component, respectively; and superposing the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component, to obtain a superposed dispersion spectrum.

13. The method according to claim 12, wherein the obtaining of the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component according to the seismic wave component and the electric field component, respectively, comprises: extracting respectively an imaginary part of the seismic wave component and an imaginary part of the electric field component, before performing a frequency scanning to obtain the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component.

14. The method according to any one of claim 12 or 13, wherein the seismic wave component comprises: a seismic wave radial component and/or a seismic wave vertical component; and the electric field component comprises: an electric field radial component and/or an electric field vertical component.

15. The method according to claim 14, wherein the superposing of the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component, to obtain the superposed dispersion spectrum, comprises: superposing the dispersion spectrum of the seismic wave vertical component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum; or superposing the dispersion spectrum of the seismic wave radial component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum; or superposing the dispersion spectrum of the seismic wave vertical component, the dispersion spectrum of the seismic wave radial component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum.

16. The method according to claim 11, wherein the performing of the extraction processing onto the superposed dispersion spectrum to obtain the dispersion curve, comprises: extracting the superposed dispersion spectrum by an energy maximum to obtain the dispersion curve; and the dispersion curve is a multi-mode dispersion curve comprising the basic-order mode and the high-order mode.

17. The method according to claim 11, wherein performing of the inversion processing onto the dispersion curve to obtain the stratum structure profile, comprises: establishing an initial inversion stratum model, and using a plurality of inversion methods to fit the dispersion curve to obtain the stratum structure profile; and the inversion methods comprise at least one of: a genetic algorithm, a swarm algorithm, and a quasi-Newton algorithm.

18. A surface wave prospecting device for jointly extracting Rayleigh wave frequency dispersion characteristics by the seismoelectric field, comprising: a collection device, applied to acquiring the plurality of jointly collected data, wherein the plurality of jointly collected data comprises: a seismic wave data and an electric field data; a jointly imaging device, applied to performing a jointly imaging processing onto the jointly collected data to obtain the superposed dispersion spectrum; and an inversion device, applied to performing the extraction processing onto the superposed dispersion spectrum to obtain the dispersion curve, and performing the inversion processing onto the dispersion curve to obtain the stratum structure profile.

19. A terminal device, comprising: a processor and a storage coupled to the processor, the storage has a surface wave prospecting program for extracting the Rayleigh wave dispersion characteristics jointly by the seismoelectric field stored, and when the program is executed by the processor, a plurality of steps of the surface wave prospecting method for jointly extracting the Rayleigh wave frequency dispersion characteristics by the seismoelectric field according to any one of claims 1-7 are implemented.

20. A storage medium, wherein the storage medium has a surface wave prospecting program for jointly extracting the Rayleigh wave dispersion characteristics by the seismoelectric field stored, and when program is executed by the processor, a plurality of steps of the surface wave prospecting method for jointly extracting the Rayleigh wave frequency dispersion characteristics by the seismoelectric field according to any one of claims 1-7 are implemented.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0041] FIG. 1 illustrates a shot set record of a linear observation system in the present invention receiving the seismic wave components u.sub.r, u.sub.z, and the electric field components E.sub.r, E.sub.z.

[0042] FIG. 2 illustrates a dispersion spectrum obtained by the shot set record of the seismic wave components u.sub.r, u.sub.z, and the electric field components E.sub.r, E.sub.z in the present invention.

[0043] FIG. 3a illustrates a superposed dispersion spectrum D.sub.Seismic obtained by the dispersion spectrum of the seismic wave components u.sub.r, u.sub.z in the present invention.

[0044] FIG. 3b illustrates a superposed dispersion spectrum D.sub.Seismoelectric obtained by the dispersion spectrum of the radial electric field E.sub.r and the dispersion spectrum of the seismic wave components u.sub.r, u.sub.z in the present invention.

[0045] FIG. 4a illustrates a superposed dispersion spectrum obtained by a superposition of u.sub.z, E.sub.r in the present invention.

[0046] FIG. 4b illustrates a superposed dispersion spectrum obtained by a superposition of u.sub.r, u.sub.z, E.sub.r in the present invention.

[0047] FIG. 5 illustrates a schematic diagram on the dispersion curve in the present invention.

[0048] FIG. 6 illustrates a formation velocity structure obtained by an inversion to the multi-mode dispersion curve in the present invention.

[0049] FIG. 7 illustrates a velocity structure profile of a formation transverse wave drawn by an interpolation to the velocity structure inversed by a plurality of surface wave in a survey line.

[0050] FIG. 8 illustrates a schematic diagram on a jointly collection of the seismic and electrical collection of the present invention.

[0051] FIG. 9 illustrates a flow chart on the surface wave prospecting method for jointly extracting the Rayleigh wave frequency dispersion characteristics by a seismoelectric field.

DETAILED DESCRIPTION OF THE EMBODIMENTS

[0052] In order to make the purpose, technical solution and the advantages of the present invention clearer and more explicit, further detailed descriptions of the present invention are stated herein, referencing to the attached drawings and some embodiments of the present invention. It should be understood that the detailed embodiments of the invention described herein are used to explain the present invention only, instead of limiting the present invention.

[0053] Referencing from FIG. 1 to FIG. 9 together, the present invention provides a plurality of preferred embodiments on a surface wave prospecting method for jointly extracting a plurality of Rayleigh wave frequency dispersion characteristics by a seismoelectric field.

[0054] A seismoelectric effect was first discovered in 1939 by a former Soviet Union scientist Ivanov, without applying a voltage, an electric field caused by a seismic wave was measured, and it was pointed out that this may be related to an electrokinetic effect of an electric double layer. Since then, Frenkel established a wave theory in a porous medium that a solid phase motion is coupled with a fluid motion between a plurality of pores, which was further improved by Biot, before forming a well-known Biot theoretical model of elastic wave in porous medium. In the 1990s, following a development of electronic science and technology and an improvement in a weak signal extraction and a signal processing method, a study on the seismoelectric effect has made a great progress. In 1994, based on the work of Frenkel and Biot, Pride proposed a group of equations on a macroscopic control for seismoelectric coupling. The group of equations has coupled a group of Biot pores elastic equations with a group of Maxwell electromagnetic equations together through a plurality of electrokinetic coupling coefficients, which describes a coupling relationship between a wave field and an electromagnetic field in a saturated porous medium, which has now become a basis of a seismoelectric theory research.

[0055] In a wave field simulation study of a layered porous medium model, it is found that when a seismic wave enters a porous medium interface exceeding a critical angle, a kind of electromagnetic (EM) wave, which is called evanescent EM wave, is generated by a seismoelectric conversion and is featured by rapid amplitude decay when leaving away from the interface. A plurality of studies have found that an evanescent electromagnetic wave generated by a Rayleigh wave excited by a vertical single-force point source has a same dispersion characteristic as the Rayleigh wave, which has provided a theoretical basis for extracting a dispersion information of a surface wave by receiving an evanescent electromagnetic wave signal generated by the surface wave. Studies have shown that a dispersion spectrum calculated by the evanescent electromagnetic wave contains a high-order mode information in rich and having a high-quality, while a dispersion in the high-order modes cannot be provided by a dispersion spectrum calculated directly by using a seismic wave (a Rayleigh wave). In addition, for a basic-order mode, a frequency range of the evanescent electromagnetic wave having a significant energy in the dispersion spectrum is not totally as same as that in a frequency range of the seismic wave, while forming a complementary. Therefore, the evanescent electromagnetic wave (EM) generated by a seismoelectric conversion has a great potential to improve a quality of a Rayleigh wave dispersion imaging, which has a great significance for further improving an accuracy of a Rayleigh surface wave prospecting.

[0056] Shown as FIG. 9, a surface wave prospecting method for jointly extracting the Rayleigh wave frequency dispersion characteristics by a seismoelectric field disclosed in the present invention, comprises a plurality of following steps:

[0057] S100. acquiring a plurality of jointly collected data; wherein the jointly collected data comprises a seismic wave data and an electric field data.

[0058] Specifically, the jointly collected data may further comprises a magnetic field data. The seismic wave data comprises a seismic wave component; the electric field data composes an electric field component. More specifically, the seismic wave component comprises a seismic wave radial component u.sub.r and/or a seismic wave vertical component u.sub.z; the electric field component comprises a electric field radial component E.sub.r, and/or an electric field vertical component E.sub.z.

[0059] In a practical engineering application, a seismic source used in a transient Rayleigh wave prospecting is generally generated by a hammer hit, a hammer drop and more seismic sources, which can be approximated as a single point force source in a vertical direction, and generates theoretically only a longitudinal wave (P wave) and a vertical component transverse wave. (SV wave) as well as an electromagnetic field of a transverse magnetic mode. The P wave and the SV wave after excitation are coherent before generating a Rayleigh wave and an evanescent electromagnetic wave correspondingly. Thus in an actual prospecting, a radial component and a vertical component of the Rayleigh wave (in a cylindrical coordinate system), as well as a radial component and a vertical component of the electric field, together with a transverse component of a magnetic field. Since the component of the magnetic field is weak, in the present embodiment, only the seismic wave components and the electric field components are considered.

[0060] Since a signal of the electric field component is weaker than a signal of the seismic wave, the collection device adopts a different noise reduction system and a different front amplification system. After a seismic source is excited by a hammer hit, two different sensors of cosmic and electric will conduct a high-frequency digital sampling at a same time.

[0061] S200. performing a jointly imaging processing onto the jointly collected data to obtain a superposed dispersion spectrum.

[0062] Specifically, the step S200 comprises:

[0063] S210. obtaining a dispersion spectrum of the seismic wave component and a dispersion spectrum of the electric field component according to the seismic wave component and the electric field component, respectively.

[0064] S210 comprises a plurality of steps:

[0065] S211. extracting respectively an imaginary part of the seismic wave component and an imaginary part of the electric field component, before performing a frequency scanning to obtain the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component.

[0066] In a horizontal layered formation model, a response wave field of a pulse seismic source time function, such as G(r,θ,z,ω) can be expressed as:

[00001]$\begin{array}{cc}G\left(r,\theta ,z,\omega \right)=\frac{1}{2\pi }\underset{m=-\infty }{\overset{+\infty }{.Math.}}{\int }_0^{+\infty }\left\{G_{T,m}\left(z,k^{\prime },\omega \right)T_k^m\left(r,\theta \right)+G_{S,m}\left(z,k^{\prime },\omega \right)S_k^m\left(r,\theta \right)+G_{R,m}\left(z,k^{\prime },\omega \right)R_k^m\left(r,\theta \right)\right\}{k}^{\prime }d{k}^{\prime },& \left(1\right)\end{array}$

[0067] Wherein, k′ is all possible horizontal wave numbers, T.sub.k.sup.m(r,θ) Sk (r,θ) R.sub.k.sup.m (r,θ) are a set of complete orthogonal vector basis functions, a specific expression is as follows:

[00002]$\begin{array}{cc}\{\begin{array}{c}{T}_k^m\left(r,\theta \right)=\left(1/k^{\prime }\right)\nabla \times \left[e_z{Y}_{k^{\prime }}^m\left(r,\theta \right)\right],\\ S_{k^{\prime }}^m\left(r,\theta \right)=\left(1/k^{\prime }\right)\nabla Y_{k^{\prime }}^m\left(r,\theta \right),\\ R_{k^{\prime }}^m\left(r,\theta \right)=-e_z{Y}_{k^{\prime }}^m\left(r,\theta \right),\end{array}& \left(2\right)\end{array}$

[0068] Wherein Y.sub.k′.sup.m(r,θ)=J.sub.m(k′r)e.sup.imθ, i is an imaginary unit, m=0,±1, ±2, . . . , k′∈(0,+∞). J.sub.m(k′r) represents an m-order Bessel function of first kind. G.sub.T,m, G.sub.S,m and G.sub.R,m in formula 1 are a plurality of expending coefficients according to a basis function, represented specifically as:

[00003]$\begin{array}{cc}\begin{array}{cc}\{\begin{array}{cc}\begin{array}{c}{G}_{T,m}\left(z,k^{\prime },\omega \right)={\int }_0^{2\pi }{\int }_0^{+\infty }G\left(r,\theta ,z,\omega \right){▯\left[T_{k^{\prime }}^m\left(r,\theta \right)\right]}^{*}\mathrm{rdrd}\theta ,\\ G_{S,m}\left(z,k^{\prime },\omega \right)={\int }_0^{2\pi }{\int }_0^{+\infty }G\left(r,\theta ,z,\omega \right){▯\left[S_{k^{\prime }}^m\left(r,\theta \right)\right]}^{*}\mathrm{rdrd}\theta ,\\ G_{R,m}\left(z,k^{\prime },\omega \right)={\int }_0^{2\pi }{\int }_0^{+\infty }G\left(r,\theta ,z,\omega \right){▯\left[R_{k^{\prime }}^m\left(r,\theta \right)\right]}^{*}\mathrm{rdrd}\theta ,\end{array}& \end{array}& \end{array}& \left(3\right)\end{array}$

[0069] Wherein a symbol * represents a complex conjugate. In a cylindrical coordinate system, G(r,θ,z,ω) may also be expressed as:

G r,θ,z,ω)=G,(r,θ,z,ω)e.sub.r+G.sub.θ(r,θ,z,ω)e.sub.θ+G.sub.z(r,θ,z,ω)e.sub.z.  (4)

[0070] Combining (1), (2), (3), (4), and considering a derivative property of 0-order Bessel function of first kind: J.sub.0′(k′r)=−J.sub.1(k′r), we will get:

[00004]$\begin{array}{cc}\{\begin{array}{c}{G}_r\left(r,\omega ,z\right)=\frac{1}{2\pi }{\int }_0^{+\infty }\left[-G_{S,0}\left(k^{\prime },\omega ,z\right)J_1\left(k^{\prime }r\right)\right]k^{\prime }{\mathrm{dk}}^{\prime },\\ G_{\theta }\left(r,\omega ,z\right)=\frac{1}{2\pi }{\int }_0^{+\infty }\left[G_{T,0}\left(k^{\prime },\omega ,z\right)J_1\left(k^{\prime }r\right)\right]k^{\prime }{\mathrm{dk}}^{\prime },\\ G_z\left(r,\omega ,z\right)=\frac{1}{2\pi }{\int }_0^{+\infty }\left[-G_{R,0}\left(k^{\prime },\omega ,z\right)J_0\left(k^{\prime }r\right)\right]k^{\prime }{\mathrm{dk}}^{\prime }.\end{array}& \left(5\right)\end{array}$

[0071] According to an orthogonal property of the Bessel function:

[00005]${\int }_0^{+\infty }{J}_m\left(k^{\prime }r\right)J_m\left(\mathrm{kr}\right)\mathrm{rdr}=\frac{\delta \left(k^{\prime }-k\right)}{k^{\prime }}$

(wherein k is the horizontal wave number observed), from formula 5, one can get:

[00006]$\begin{array}{cc}{\int }_0^{+\infty }{G}_r\left(r,\omega ,z\right)J_1\left(\mathrm{kr}\right)\mathrm{rdr}=\frac{1}{2\pi }{\int }_0^{+\infty }\left[-G_{S,0}\left(k^{\prime },\omega ,z\right){\int }_0^{+\infty }{J}_1\left(k^{\prime }r\right)J_1\left(\mathrm{kr}\right)\mathrm{rdr}\right]k^{\prime }{\mathrm{dk}}^{\prime }=-\frac{1}{2\pi }{G}_{S,0}\left(k,\omega ,z\right),& \left(6\right)\end{array}$$\begin{array}{cc}{\int }_0^{+\infty }{G}_{\theta }\left(r,\omega ,z\right)J_1\left(\mathrm{kr}\right)\mathrm{rdr}=\frac{1}{2\pi }{\int }_0^{+\infty }\left[-G_{T,0}\left(k^{\prime },\omega ,z\right){\int }_0^{+\infty }{J}_1\left(k^{\prime }r\right)J_1\left(\mathrm{kr}\right)\mathrm{rdr}\right]k^{\prime }{\mathrm{dk}}^{\prime }=\frac{1}{2\pi }{G}_{T,0}\left(k,\omega ,z\right),& \left(7\right)\end{array}$$\begin{array}{cc}{\int }_0^{+\infty }{G}_z\left(r,\omega ,z\right)J_0\left(\mathrm{kr}\right)\mathrm{rdr}=\frac{1}{2\pi }{\int }_0^{+\infty }\left[-G_{R,0}\left(k^{\prime },\omega ,z\right){\int }_0^{+\infty }{J}_0\left(k^{\prime }r\right)J_0\left(\mathrm{kr}\right)\mathrm{rdr}\right]k^{\prime }{\mathrm{dk}}^{\prime }=-\frac{1}{2\pi }{G}_{R,0}\left(k,\omega ,z\right).& \left(8\right)\end{array}$

[0072] Defining a kernel function as follows:

[00007]$\begin{array}{cc}\{\begin{array}{c}{G}_S\left(k,\omega ,z\right)=-\frac{1}{2\pi }{G}_{S,0}\left(k,\omega ,z\right),\\ G_T\left(k,\omega ,z\right)=\frac{1}{2\pi }{G}_{T,0}\left(k,\omega ,z\right),\\ G_R\left(k,\omega ,z\right)=-\frac{1}{2\pi }{G}_{R,0}\left(k,\omega ,z\right).\end{array}& \left(9\right)\end{array}$

[0073] Taking the formula 9 into the formulas (6)-(8), and getting:

[00008]$\begin{array}{cc}\{\begin{array}{c}{\int }_0^{+\infty }{G}_r\left(r,\omega ,z\right)J_1\left(kr\right)rdr=G_S\left(k,\omega ,z\right),\\ {\int }_0^{+\infty }{G}_{\theta }\left(r,\omega ,z\right)J_1\left(kr\right)rdr=G_T\left(k,\omega ,z\right),\\ {\int }_0^{+\infty }{G}_z\left(r,\omega ,z\right)J_0\left(kr\right)rdr=G_R\left(k,\omega ,z\right).\end{array}& \left(10\right)\end{array}$

[0074] Since there is a following relationship between a response wave field G(r,θ,z,ω) of a pulse source time function and an actual source time function F(ω) as well as an actual wave field U(r,θ,z,ω):

U(r,θ,z,ω)=F(ω)G(r,θ,z,ω)  (11)

[0075] Therefore, the radial component u.sub.r and the vertical component u.sub.z of the seismic wave (in a cylindrical coordinate system), as well as the radial component E.sub.r and the vertical component E.sub.z of the electric field, according to the formulas of (10) and (11), may be expressed as:

[00009]$\begin{array}{cc}\{\begin{array}{c}{\int }_0^{+\infty }\left[u_r\left(r,\omega ,z\right)/F\left(\omega \right)\right]J_1\left(\mathrm{kr}\right)\mathrm{rdr}=u_S\left(k,\omega ,z\right),\\ {\int }_0^{+\infty }\left[u_z\left(r,\omega ,z\right)/F\left(\omega \right)\right]J_0\left(\mathrm{kr}\right)\mathrm{rdr}=u_R\left(k,\omega ,z\right),\\ {\int }_0^{+\infty }\left[E_r\left(r,\omega ,z\right)/F\left(\omega \right)\right]J_1\left(\mathrm{kr}\right)\mathrm{rdr}=E_S\left(k,\omega ,z\right),\\ {\int }_0^{+\infty }\left[E_z\left(r,\omega ,z\right)/F\left(\omega \right)\right]J_0\left(\mathrm{kr}\right)\mathrm{rdr}=E_R\left(k,\omega ,z\right),\end{array}& \left(12\right)\end{array}$

[0076] Wherein u.sub.S(k,ω,z), u.sub.R (k,ω,z), E.sub.S(k,ω,z), E.sub.R (k,ω,z) are the core functions of the wave and the electric field; F(ω) represents the actual source time function; u.sub.r(r,ω,z), u.sub.z (r,ω,z), E.sub.r(r,ω,z), E.sub.z (r,ω,z) are a plurality of Fourier transforms of a plurality of time domain wavefield components u.sub.r (r,t,z), u.sub.z (r,t,z), E.sub.r(r,t,z), E.sub.z (r,t,z); in a real observation, r and z represents individually a radial distance between a receiver and the source as well as a depth of the receiver.

[0077] Therefore, in a surface wave prospecting, an integral in (12) can be approximated as a summation, and an imaginary part thereof is taken for a frequency scanning before obtaining the dispersion spectrum. A specific formula is as follows:

[00010]$\begin{array}{cc}\{\begin{array}{c}\mathrm{Im}\left[\underset{j=1}{\overset{N}{.Math.}}{\int }_{r_{j-1}}^{r_j}\left[u_r\left(r,\omega ,z\right)/F\left(\omega \right)\right]J_1\left(kr\right)rdr\right]={\tilde{u}}_r\left(k,\omega ,z\right),\\ \mathrm{Im}\left[\underset{j=1}{\overset{N}{.Math.}}{\int }_{r_{j-1}}^{r_j}\left[u_z\left(r,\omega ,z\right)/F\left(\omega \right)\right]J_0\left(kr\right)rdr\right]={\tilde{u}}_z\left(k,\omega ,z\right),\\ \mathrm{Im}\left[\underset{j=1}{\overset{N}{.Math.}}{\int }_{r_{j-1}}^{r_j}\left[E_r\left(r,\omega ,z\right)/F\left(\omega \right)\right]J_1\left(kr\right)rdr\right]={\tilde{E}}_r\left(k,\omega ,z\right),\\ \mathrm{Im}\left[\underset{j=1}{\overset{N}{.Math.}}{\int }_{r_{j-1}}^{r_j}\left[E_z\left(r,\omega ,z\right)/F\left(\omega \right)\right]J_0\left(kr\right)rdr\right]={\tilde{E}}_z\left(k,\omega ,z\right),\end{array}& \left(15\right)\end{array}$

[0078] Wherein ũ.sub.r (k,ω,z), ũ.sub.z (k,ω,z), {tilde over (E)}.sub.r (k,ω,z) and {tilde over (E)}.sub.z (k,ω,z) represent the dispersion spectrums obtained by an approximation; Im[.Math.] represents calculating an imaginary part of a complex variable; k represents the horizontal wave number of an observation, ω is the angular frequency, and N represents a total number of recording tracks; r.sub.j represents the radial distance from a j-th geophone to the source, since a surface wave prospecting is based on a horizontal strata, all sensors and the seismic source, do not require a plurality of the sensors to be arranged in a linear equidistance or in a plurality of other regular shapes, or maintain a certain offset from the seismic source, therefore the geophones can be arranged arbitrarily according to an actual situation. That is, in the present embodiment, through a jointly acquisition and imaging of the seismic wave and the electric field data (an existing surface wave prospecting uses only the seismic wave components for a detection), a collection system may be arbitrarily arranged according to a detection area, requiring neither a linear or other regular shapes nor an offset. Therefore, they can be arranged in any shapes according to an environment of a detection site, which can reflect an underground velocity structure of the detection area more comprehensively, while comparing to a traditional surface wave observation system evenly spaced in a straight line and requiring a certain minimum offset, it has a higher applicability.

[0079] S220, superposing the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component, to obtain the superposed dispersion spectrum.

[0080] There is often an underground water level in a shallow surface medium of the Earth, and a medium above the water level, especially a medium close to the surface, can be approximated as an air-filled porous medium. Therefore, a four-layer pore half-space model is established according to an actual formation condition, and a plurality of parameters are shown in Table 1.

TABLE-US-00001 TABLE 1 Parameters of a four-layer half-space porous medium model with a high-speed interlayer First layer Second layer Third layer Fourth layer Property medium medium medium medium Type of pore fluid air water air air Porosity φ (%) 5 25 0.1 0.1 Darcy permeability κ.sub.0 (meters   1 × 10.sup.−14 1 × 10.sup.−12   1 × 10.sup.−20   1 × 10.sup.−20 square) Solid volume modulus K.sub.s (GPa) 4.4 7.0 4.46 6.83 Fluid volume modulus K.sub.f (GPa) 1.43 × 10.sup.−4 2.2 1.43 × 10.sup.−4 1.43 × 10.sup.−4 Solid skeleton volume modulus 3.928 2.99 4.458 6.829 K.sub.fr (GPa) Solid skeleton shear modulus G 0.0577 0.2266 0.1125 0.6984 (GPa) Fluid viscosity coefficient η (Pa .Math. s)  1.8 × 10.sup.−5 1 × 10.sup.−3.sub.   1.8 × 10.sup.−5  1.8 × 10.sup.−5 Solid density ρ.sub.s (Kg/m.sup.3) 1874 2133 1802 1942 Fluid density ρ.sub.f (Kg/m.sup.3) 1.21 1000 1.21 1.21 Fluid salinity C.sub.0 (Mole/liter) 0 0.1 0 0 Absolute temperature T (Kalvin) 298 298 298 298 Relative permittivity of fluid κ.sub.f 1 80 1 1 Relative permittivity of solid κ.sub.s 4 4 4 4 Curvature α.sub.∞ 10 3 50 50 Fluid conductivity σ.sub.f (Siemens/m) 0 0.9275 0 0 Volume conductivity σ 0.01 0.0773 0.01 0.01 (Siemens/m) Fast longitudinal wave velocity V.sub.p 1500 1700 1600 2000 (m/s) Shear wave velocity V.sub.s (m/s) 180 350 250 600 Stratum depth z (m) 0~10 10~20 20~40 40~

[0081] A second layer of the model is arranged as a porous medium saturated with water, while a plurality of other layers are arranged as a porous medium saturated with air, and a phreatic water surface locates at a depth of 10 meters underground. The seismic source is a vertical single force point source excited by a surface (a hammer hit, a drop hammer seismic source), and a source time function adopts a Ricker wavelet with a dominant frequency of 10 Hz and a delay time of 0.5 seconds; in order to facilitate a data collection, shown as FIG. 8, an observation system adopts a linear arrangement, having an inter-channel space of 2 meters, and a total of 91 channels. The shot set records of the seismic wave component received and the shot set records of the electric field component are shown in FIG. 1. According to the shot set records of each seismic wave component and the shot set records of each electric field component, a recorded dispersion spectrum is shown as FIG. 2, wherein a dotted line represents a theoretical dispersion curve of the Rayleigh wave calculated by the formation parameters.

[0082] S220 comprises steps:

[0083] S221. superposing the dispersion spectrum of the seismic wave vertical component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum; or

[0084] S222. superposing the dispersion spectrum of the seismic wave radial component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum; or

[0085] S223. superposing the dispersion spectrum of the seismic wave vertical component, the dispersion spectrum of the seismic wave radial component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum.

[0086] From FIG. 1, it can be seen that the Rayleigh wave dominates the shot set records, and in the surface wave area, the high mode of the electric field components E.sub.r, E.sub.z are more developed than that of the seismic wave components u.sub.r, u.sub.z. The black dotted line in the dispersion spectrum of FIG. 2 is a theoretical dispersion curve. It can be seen that a position of a maximum energy in the dispersion spectrum coincides with that of the theoretical dispersion curve, thus a Rayleigh wave dispersion curve may be extracted through the maximum energy value. From FIG. 2 it can be seen that, in the dispersion spectrum (|{tilde over (E)}.sub.r|,|{tilde over (E)}.sub.z|) of the electric field component E.sub.r, E.sub.z, an imaging quality of the Rayleigh wave in the high-order mode is significantly higher than the dispersion spectrum (|ũ.sub.r|,|ũ.sub.z|) of the seismic component u.sub.r, u.sub.z; an imaging frequency band of the basic-order mode of the dispersion spectrum of u.sub.r, u.sub.z is significantly wider than that of the dispersion spectrum of E.sub.r, E.sub.z. Therefore, performing an inversion on the multi-mode dispersion curve jointly extracted by the seismic wave components and the electric field components will greatly reduce a versatility of the inversion, thereby greatly improving an accuracy and a stability of the inversion.

[0087] From FIG. 3a, it can be seen that the high-order mode imaging in the dispersion spectrum D.sub.Seismic obtained by superposing the seismic wave components u.sub.r, u.sub.z only, is still not as much as the high-order mode obtained by the electric field components. Considering that the signal of the vertical component E.sub.z of the electric field is relatively weak and a data acquisition in an actual detection is hard, thus in a practical application, the dispersion spectrum of the radial component E.sub.r of the electric field and D.sub.Seismic are superposed to obtain a new dispersion spectrum D.sub.Seismoelectric (shown as FIG. 3b), it can be seen that, in a superposed dispersion spectrum, the imaging qualities of both the base-order mode and the high-order mode are improved.

[0088] An acquisition of the seismic wave radial component u.sub.r requires a transform of two horizontal components, in an actual detection, a three-component geophone is required for an acquisition, while an acquisition of the vertical component u.sub.z of the seismic wave requires only a vertical geophone, which is not only economical but efficient, thus the dispersion spectrums by a superposition of u.sub.z, E.sub.r is compared with a superposed dispersion spectrum of three components of u.sub.r, u.sub.z, E.sub.r. FIG. 4a illustrates the dispersion spectrum obtained by superposing u.sub.z, E.sub.r, FIG. 4b illustrates the dispersion spectrum obtained by superposing u.sub.r, u.sub.z, E.sub.r, it can be seen that there is no significant difference in between. Therefore, it is also possible to perform a superposition of u.sub.z, E.sub.r only, to avoid an acquisition and a superposition of multi-seismic wave, and the electric field component. Thus it is possible to reduce an acquisition cost effectively, on a premise of ensuring a quality of the dispersion spectrum imaging.

[0089] S300. performing an extraction processing onto the superposed dispersion spectrum to obtain a dispersion curve, and performing an inversion processing onto the dispersion curve to obtain a stratum structure profile.

[0090] Specifically, the S300 comprises a plurality of steps:

[0091] S310. extracting the superposed dispersion spectrum by an energy maximum to obtain the dispersion curve.

[0092] Specifically, the dispersion curve is a multi-mode dispersion curve comprising the basic-order mode and the high-order mode. By setting a region, then by a method of using a computer program to automatically search for a maximum point in the region, a multi-mode dispersion curve comprising the basic-order mode and the high-order mode is extracted from the dispersion spectrum (see FIG. 5) for a later inversion.

[0093] S320. establishing an initial inversion stratum model, and using a plurality of inversion methods to fit the dispersion curve to obtain the stratum structure profile.

[0094] Specifically, the inversion methods comprise at least one of: a genetic algorithm, a swarm algorithm, and a quasi-Newton algorithm. Establishing an initial inversion stratigraphic model based on a basic data, then a plurality of inversion methods including the genetic algorithm, the swarm algorithm (a particle swarm, a bee colony, and more.), the quasi-Newton algorithm and more, are applied to simulate the multi-mode dispersion curve extracted from FIG. 5, therefore the stratigraphic structure is inversed, shown as FIG. 6.

[0095] A velocity structure of multiple stratum positions will be obtained by moving the observation system, and the velocity structures of a plurality of different positions are interpolated, before being able to draw a stratum section according to a plurality of spatial positions, shown as FIG. 7.

[0096] Based on any one of the embodiments on the surface wave prospecting method for jointly extracting the Rayleigh wave frequency dispersion characteristics by the seismoelectric field, the present invention further provides a plurality of preferred embodiments on a surface wave prospecting device for jointly extracting the Rayleigh wave frequency dispersion characteristics by the seismoelectric field:

[0097] The surface wave prospecting device for jointly extracting the Rayleigh wave frequency dispersion characteristics by the seismoelectric field as described in an embodiment of the present invention comprises: a collection device, applied to acquiring the plurality of jointly collected data; wherein the plurality of jointly collected data comprises: a seismic wave data and an electric field data;

[0098] a jointly imaging device, applied to performing a jointly imaging processing onto the jointly collected data to obtain the superposed dispersion spectrum;

[0099] an inversion device, applied to performing the extraction processing onto the superposed dispersion spectrum to obtain the dispersion curve, and performing the inversion processing onto the dispersion curve to obtain the stratum structure profile.

[0100] The seismic wave data comprises: a seismic wave component; the electric field data composes: an electric field component.

[0101] The jointly imaging device is specifically applied to obtaining a dispersion spectrum of the seismic wave component and a dispersion spectrum of the electric field component according to the seismic wave component and the electric field component, respectively; and superposing the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component, before obtaining a superposed dispersion spectrum.

[0102] The jointly imaging device is further applied to extracting respectively an imaginary part of the seismic wave component and an imaginary part of the electric field component, before performing a frequency scanning to obtain the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component.

[0103] The seismic wave component comprises: a seismic wave radial component and/or a seismic wave vertical component; the electric field component comprises: an electric field radial component and/or an electric field vertical component.

[0104] The jointly imaging device is further applied to superposing the dispersion spectrum of the seismic wave vertical component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum; or

[0105] superposing the dispersion spectrum of the seismic wave radial component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum; or

[0106] superposing the dispersion spectrum of the seismic wave vertical component, the dispersion spectrum of the seismic wave radial component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum.

[0107] The inversion device is specifically applied to extracting the superposed dispersion spectrum by an energy maximum to obtain the dispersion curve; the dispersion curve is a multi-mode dispersion curve comprising the basic-order mode and the high-order mode.

[0108] The inversion device is further applied to establishing an initial inversion stratum model, and using a plurality of inversion methods to fit the dispersion curve to obtain the stratum structure profile; the inversion methods comprise at least one of: a genetic algorithm, a swarm algorithm, and a quasi-Newton algorithm.

[0109] Based on any one of the embodiments on the surface wave prospecting method for jointly extracting the Rayleigh wave frequency dispersion characteristics by the seismoelectric field, the present invention further provides a plurality of preferred embodiments on a terminal device:

[0110] the terminal device described in an embodiment of the present invention, comprises: a processor, and a storage connected to the processor.

[0111] the storage has a surface wave prospecting program for jointly extracting the Rayleigh wave dispersion characteristics by the seismoelectric field stored, and a plurality of following steps are implemented when the surface wave prospecting program for extracting the Rayleigh wave dispersion characteristics jointly by the seismoelectric field is executed by the processor:

[0112] acquiring the plurality of jointly collected data; wherein the plurality of jointly collected data comprises: a seismic wave data and an electric field data;

[0113] performing the jointly imaging processing onto the jointly collected data to obtain the superposed dispersion spectrum;

[0114] performing the extraction processing onto the superposed dispersion spectrum to obtain the dispersion curve, and performing the inversion processing onto the dispersion curve to obtain the stratum structure profile.

[0115] The seismic wave data comprises the seismic wave component; the electric field data composes the electric field component.

[0116] When the surface wave prospecting program for extracting the Rayleigh wave dispersion characteristics jointly by the seismoelectric field is executed by the processor, a plurality of following steps are further implemented:

[0117] obtaining a dispersion spectrum of the seismic wave component and a dispersion spectrum of the electric field component according to the seismic wave component and the electric field component, respectively;

[0118] superposing the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component, to obtain a superposed dispersion spectrum.

[0119] When the surface wave prospecting program for extracting the Rayleigh wave dispersion characteristics jointly by the seismoelectric field is executed by the processor, a plurality of following steps are further implemented:

[0120] extracting respectively an imaginary part of the seismic wave component and an imaginary part of the electric field component, before performing a frequency scanning to obtain the dispersion spectrum of the seismic wave component and the dispersion spectrum of the electric field component.

[0121] The seismic wave component comprises: a seismic wave radial component and/or a seismic wave vertical component; the electric field component comprises: an electric field radial component and/or an electric field vertical component.

[0122] When the surface wave prospecting program for extracting the Rayleigh wave dispersion characteristics jointly by the seismoelectric field is executed by the processor, a plurality of following steps are further implemented:

[0123] superposing the dispersion spectrum of the seismic wave vertical component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum; or

[0124] superposing the dispersion spectrum of the seismic wave radial component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum; or

[0125] superposing the dispersion spectrum of the seismic wave vertical component, the dispersion spectrum of the seismic wave radial component and the dispersion spectrum of the electric field radial component to obtain the superposed dispersion spectrum.

[0126] When the surface wave prospecting program for extracting the Rayleigh wave dispersion characteristics jointly by the seismoelectric field is executed by the processor, a plurality of following steps are further implemented:

[0127] extracting the superposed dispersion spectrum by an energy maximum to obtain the dispersion curve.

[0128] The dispersion curve is a multi-mode dispersion curve comprising the basic-order mode and the high-order mode.

[0129] When the surface wave prospecting program for extracting the Rayleigh wave dispersion characteristics jointly by the seismoelectric field is executed by the processor, a plurality of following steps are further implemented:

[0130] establishing an initial inversion stratum model, and using a plurality of inversion methods to fit the dispersion curve to obtain the stratum structure profile.

[0131] The inversion methods comprise at least one of: a genetic algorithm, a swarm algorithm, and a quasi-Newton algorithm.

[0132] Based on any one of the embodiments on the surface wave prospecting method for jointly extracting the Rayleigh wave frequency dispersion characteristics by the seismoelectric field, the present invention further provides a plurality of preferred embodiments on a storage medium:

[0133] The storage medium described in an embodiment of the present invention, wherein has the surface wave prospecting program for jointly extracting the Rayleigh wave dispersion characteristics by the seismoelectric field stored thereon, and a plurality of following steps are implemented when the surface wave prospecting program for extracting the Rayleigh wave dispersion characteristics jointly by the seismoelectric field is executed by the processor:

[0134] acquiring the plurality of jointly collected data; wherein the plurality of jointly collected data comprises: a seismic wave data and an electric field data;

[0135] performing the jointly imaging processing onto the jointly collected data to obtain the superposed dispersion spectrum;

[0136] performing the extraction processing onto the superposed dispersion spectrum to obtain the dispersion curve, and performing the inversion processing onto the dispersion curve to obtain the stratum structure profile.

[0137] It should be understood that, the application of the present invention is not limited to the above examples listed. Ordinary technical personnel in this field can improve or change the applications according to the above descriptions, all of these improvements and transforms should belong to the scope of protection in the appended claims of the present invention.