CN104181593B - Three-dimensional ray-free tracing contraflexure wave tomographic imaging method and device - Google Patents

Abstract

The invention relates to a three-dimensional ray-free tracing contraflexure wave tomographic imaging method and device. The three-dimensional ray-free tracing contraflexure wave tomographic imaging method comprises the steps that single-shot first arrivals are picked up, first arrivals at multiple azimuthal angles are selected according to an azimuthal angle interval among all the single-shot first arrivals, offset-distance progressive increase sequencing is carried out on the first arrivals at all the azimuthal angles, and therefore a first arrival curve within the range of 0-x is formed; low-order polynomial least square method global fitting is carried out on the first arrival curve within the range of 0-x with the offset distance as the horizontal coordinate and the first arrival time as the vertical coordinate, and multinomial coefficients corresponding to the azimuthal angle of each single-shot first arrival are obtained; an azimuthal angle slowness curve corresponding to each single-shot first arrival is obtained according to the corresponding multinomial coefficients; the velocity and depth of a contraflexure wave contraflexure point corresponding to the azimuthal angle selected by each single-shot first arrival which is picked up are obtained according to the corresponding azimuthal angle slowness curve of each single-shot first arrival; a three-dimensional velocity field is obtained through interpolation according to the velocity and depth of the contraflexure wave contraflexure point corresponding to the azimuthal angle selected by each single-shot first arrival which is picked up.

Description

Three-dimensional ray-tracing-free refracted wave tomography method and device

Technical Field

The invention relates to the technical field of seismic data processing, in particular to a three-dimensional ray-tracing-free refracted wave tomography method and device.

Background

In seismic data processing, the travel-time tomography method is an important velocity imaging method. The travel time tomography can be divided into two main steps of forward traveling and slowness inversion according to a calculation process. Forward travel includes travel time calculations, and typically also includes ray tracing. Travel time tomography can be simply divided into two categories depending on whether ray tracing is performed after forward travel time:

(1) conventional ray tracing travel time tomography methods. The method sets the speed distribution on discrete grid nodes, and the surface element speed between the nodes is constant or gradually changed. The ray path of the grid after the travel time is evolved can be obtained by using a gradient method and an interchange method, so that the error of the travel time and the first arrival of each receiving point is obtained, and then the velocity field is updated by iterative inversion.

The existing three-dimensional inflection wave tomography mainly adopts a traditional ray tracing travel time tomography method, the restraint speed is transversely unchanged at the moment, and the vertical direction gradually increases gradually, so that the condition of inflection wave propagation is met, and then three-dimensional ray tracing and inversion are carried out iteratively to update the slowness.

The method has the advantages that no assumption is made on the underground speed distribution, and the method can adapt to any complex or simple speed field. The method has the disadvantages that the complex velocity field has higher requirements on a ray path calculation method, and complex problems of travel time, ray caustic, discontinuous velocity gradient and the like need to be processed. In addition, ray paths are generally curves, so that the inversion presents strong nonlinearity, the solution is difficult, and multiple solutions exist. When the appointed speed is only gradually changed along the vertical direction, the speed field with the wave gradually changed along the vertical direction is transmitted in a form of a reverse wave and is reversely reflected to the ground surface, and the traditional reverse wave travel time tomography method is obtained.

(2) A non-ray tracing echo traveling time tomography method. This type of method also provides for the speed to be ramped only vertically. However, according to the propagation rule of the inflection wave, the depth and the speed of the inflection point of the inflection wave can be calculated by derivation formula analysis, so that the complex processes of forward traveling and ray tracing and inversion updating slowness are avoided.

The existing three-dimensional inflection wave tomography mainly adopts a ray-tracing-free inflection wave travel time tomography method, which is essentially a two-dimensional algorithm and cannot carry out three-dimensional velocity field imaging.

In summary, the conventional method is to perform three-dimensional decomposition (ospypov, 2000) on the first arrival time of three-dimensional seismic data and then apply a two-dimensional algorithm to the decomposed travel time for calculation. The decomposition algorithm is complex and requires multiple iterations.

Disclosure of Invention

In order to solve the problems in the prior art, the invention provides a three-dimensional ray-tracing-free inflection wave tomography method and a three-dimensional ray-tracing-free inflection wave tomography device, which extend a two-dimensional ray-tracing-free inflection wave tomography technology to three-dimensional seismic data processing and avoid a complex three-dimensional travel time iterative decomposition process.

In order to achieve the above object, the present invention provides a three-dimensional no-ray tracing folded wave tomography method, which comprises:

picking up single cannon first arrivals, selecting the first arrivals of at least two azimuth angles for each picked single cannon first arrival according to the azimuth angle interval, and sequencing each azimuth angle first arrival according to the incremental mode of offset distance to form a first arrival curve in the range of 0-x; wherein x is the maximum offset distance;

taking the offset as a horizontal coordinate and the first arrival time as a vertical coordinate, and performing low-order polynomial least square method integral fitting on the first arrival curve in the range of 0-x to obtain a polynomial coefficient corresponding to the azimuth angle of each single shot first arrival;

obtaining an azimuth slowness curve of the first arrival of each single shot correspondingly according to the polynomial coefficients;

obtaining the speed and the depth of a corresponding inflection point of the selected azimuth angle of each single shot first arrival according to the azimuth angle slowness curve of each single shot first arrival;

and interpolating into a three-dimensional velocity field according to the velocity and the depth of the corresponding inflection point of the inflection wave of each picked single shot first arrival selected azimuth angle.

Preferably, the speed obtaining method of the inflection point of the corresponding inflection wave of the selected azimuth angle of each picked single shot first arrival comprises the following steps:

randomly taking M points as inflection points of the inflection waves within the range of 0-x for the azimuth slowness curve of the first arrival of each picked single cannon;

calculating the reciprocal of the slowness curve value of the position corresponding to the point i to obtain the speed of the point i; wherein i is 1,2, …, M.

Preferably, the depth obtaining method of the inflection point corresponding to the selected azimuth angle of each picked single shot first arrival comprises:

the depth of the point i is obtained by adopting an integral formula; wherein the integral formula is:

$z_i=\frac{1}{\mathrm{\π}}{\∫}_0^{x_i}\cos h^{-1}\left(\frac{\frac{\mathrm{dt}}{\mathrm{dx}}}{{\left(\frac{\mathrm{dt}}{\mathrm{dx}}\right)}_i}\right)\mathrm{dx}$

wherein z is_iDepth, x, corresponding to inflection point i_iIs the offset corresponding to point i.

Preferably, the method for obtaining the azimuth slowness curve of each single shot first arrival comprises the following steps:

obtaining a polynomial according to the polynomial coefficient;

and (4) carrying out derivation on the polynomial to obtain an azimuth slowness curve of the first arrival of each single shot.

Preferably, the azimuthal angular interval is a comparable or unequal azimuthal angular interval.

In order to achieve the above object, the present invention further provides a three-dimensional raytrace-free refracted-wave tomography apparatus, including:

the system comprises a first arrival curve acquisition unit, a second arrival curve acquisition unit and a third arrival curve acquisition unit, wherein the first arrival curve acquisition unit is used for picking single cannon first arrivals, selecting first arrivals of at least two azimuth angles for each picked single cannon first arrival according to azimuth angle intervals, and sequencing each azimuth angle first arrival according to an offset distance increasing mode to form a first arrival curve in a range of 0-x; wherein x is the maximum offset distance;

the fitting unit is used for performing low-order polynomial least square method overall fitting on the first arrival curve in the range of 0-x by taking the offset distance as a horizontal coordinate and the first arrival time as a vertical coordinate to obtain a polynomial coefficient corresponding to the azimuth angle of each single shot first arrival;

the azimuth slowness curve acquisition unit is used for acquiring an azimuth slowness curve of the first arrival of each single shot correspondingly according to the polynomial coefficients;

the velocity and depth acquisition unit of the inflection point is used for acquiring the velocity and depth of the corresponding inflection point of the selected azimuth angle of each single shot first arrival according to the azimuth angle slowness curve of each single shot first arrival;

and the interpolation unit is used for interpolating a three-dimensional velocity field according to the velocity and the depth of the inflection point of the corresponding inflection wave of the selected azimuth angle of the first arrival of each picked single shot.

Preferably, the velocity and depth acquiring unit of the inflection point of the inflection wave comprises: the device comprises an inflection wave inflection point acquisition module and a first calculation module; wherein,

the inflection wave inflection point acquisition module is used for arbitrarily taking M points as inflection points of the inflection wave for the azimuth angle slowness curve of the first arrival of each picked single cannon within the range of 0-x;

the first calculation module is used for calculating the reciprocal of the slowness curve value of the position corresponding to the point i to obtain the speed of the point i; wherein i is 1,2, …, M.

Preferably, the speed and depth of the inflection point of the inflection wave further comprises a second calculation module; wherein,

the second calculation module is used for obtaining the depth of the point i by adopting an integral formula; wherein the integral formula is:

$z_i=\frac{1}{\mathrm{\π}}{\∫}_0^{x_i}\cos h^{-1}\left(\frac{\frac{\mathrm{dt}}{\mathrm{dx}}}{{\left(\frac{\mathrm{dt}}{\mathrm{dx}}\right)}_i}\right)\mathrm{dx}$

wherein z is_iDepth, x, corresponding to inflection point i_iIs the offset corresponding to point i.

Preferably, the azimuth slowness curve acquiring unit includes: the device comprises a polynomial acquisition module and a derivation module; wherein,

the polynomial acquisition module is used for acquiring a polynomial according to the polynomial coefficients;

and the derivation module is used for deriving the polynomial to obtain an azimuth slowness curve of the first arrival of each single cannon.

Preferably, the azimuth angle intervals adopted by the first arrival curve acquiring unit are equal or unequal azimuth angle intervals.

The technical scheme has the following beneficial effects:

the traditional two-dimensional tomography method without ray tracing during the traveling of the refracted wave is applied in different azimuth angles along a plurality of azimuth angles of the three-dimensional seismic data, so that a three-dimensional calculation result is obtained. The technical scheme avoids a complex iterative decomposition process during three-dimensional travel, so the method is easy to realize, and the calculation speed is high when the number of the selected azimuth angles is less, such as less than 40.

Drawings

In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings used in the description of the embodiments or the prior art will be briefly described below, it is obvious that the drawings in the following description are only some embodiments of the present invention, and for those skilled in the art, other drawings can be obtained according to the drawings without creative efforts.

FIG. 1 is a flow chart of a three-dimensional non-ray tracing refracted wave tomography method according to the present invention;

FIG. 2 is a block diagram of a three-dimensional non-ray tracing refracted wave tomography device according to the present invention;

FIG. 3 is a flowchart of a three-dimensional raytrace-free echo tomography method according to the present embodiment;

fig. 4 is a tomographic effect chart of the present embodiment.

Detailed Description

The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.

The working principle of the technical scheme is as follows: the two-dimensional single-shot first arrival is extracted to obtain two-dimensional lines of a plurality of azimuth angles, a two-dimensional ray-free tracking echo travel time tomography algorithm is used for calculation, the effect of three-dimensional calculation is achieved according to the combination of the plurality of azimuth angles, and complex iterative three-dimensional travel time decomposition calculation is avoided.

Fig. 1 is a flow chart of a three-dimensional non-ray tracing refracted wave tomography method according to the present invention. The method comprises the following steps:

step 101) picking single shot first arrivals, selecting the first arrivals of at least two azimuth angles for each picked single shot first arrival according to the azimuth angle interval, and sequencing the first arrivals of each azimuth angle according to the incremental offset mode to form a first arrival curve in the range of 0-x; wherein x is the maximum offset distance;

in step 101), the single shot first arrival picking is obtained by processing software in an automatic or interactive picking mode, and the obtained first arrival data must contain the azimuth angle information of the current path. And extracting the first arrivals of a plurality of azimuth angles from the picked single-shot first arrival data according to equal or unequal azimuth angle intervals, and gradually sorting the first arrivals of each azimuth angle according to the offset distance to form a first arrival curve in the range of 0-x.

The azimuth of the single shot first arrival can be one or a plurality of azimuths. Only one azimuth is equal to the effect of two-dimensional tomographic calculations. When a plurality of azimuth angles are specified, the required angles can be flexibly selected at equal or unequal azimuth angle intervals.

Step 102), taking the offset distance as a horizontal coordinate and the first arrival time as a vertical coordinate, and performing low-order polynomial least square method integral fitting on the first arrival curve in the range of 0-x to obtain a polynomial coefficient corresponding to the azimuth angle of each single shot first arrival;

in step 102), the first-motion curve is integrally fitted to a first-motion curve based on the selected single azimuth angle: and (3) performing low-order polynomial least square method integral fitting on the first arrivals in the range of 0-x by taking the offset as a horizontal coordinate and the first arrival time as a vertical coordinate. During fitting, a polynomial which can be quickly integrated by weighted summation of polynomial roots is selected, such as Legendre polynomial and Chebyshev polynomial. The polynomial is of lower order, below 5 th order, usually 2 nd order. And obtaining a polynomial coefficient corresponding to the first-motion data after the fitting is finished.

Step 103) obtaining an azimuth slowness curve of each single shot first arrival correspondingly according to the polynomial coefficients;

in step 103), the coefficients of the polynomial in step 102) are obtained as a polynomial, and the slowness curve is obtained by taking the derivative of the polynomial. The calculation method comprises the following steps: and calculating the value of dt/dx of the derivative p in the range of 0 to x by using the derivative formula of the polynomial calculated in the step 102), thereby obtaining the slowness curve. Wherein t is the first arrival time, x is the offset, and d is the derivation operation symbol. Since the slowness curve calculation directly uses the polynomial derivative formula, the coordinates used for calculation may be arranged regularly or irregularly.

Step 104) obtaining the speed and the depth of a corresponding inflection point of the selected azimuth angle of each picked single shot first arrival according to the azimuth angle slowness curve of each single shot first arrival;

in step 104), the speed calculation methodThe method comprises the following steps: randomly taking M points from the slowness curve of the specified azimuth angle obtained in the step 103) within the range of 0-x, and obtaining the reciprocal 1/(dt/dx) of the slowness curve value at the position of the point i_i(i ═ 1,2, …, M), the velocity of the inflection point of this point was obtained.

In step 104), the depth of the point i adopts the following integral formula:

$z_i=\frac{1}{\mathrm{\π}}{\∫}_0^{x_i}\cos h^{-1}\left(\frac{\frac{\mathrm{dt}}{\mathrm{dx}}}{{\left(\frac{\mathrm{dt}}{\mathrm{dx}}\right)}_i}\right)\mathrm{dx}$

wherein z is_iDepth, x, corresponding to inflection point i_iIs an offset. The depth is calculated to be 0-x_iWithin the range, the integral kernel function is calculated point by point, and then, for 0 to x_iIntegration gives z_i. Since the polynomial in step 103) is selected to be a polynomial capable of performing fast integration by polynomial root weighted summation, the integration is performed by gaussian integration, and for M integral calculation points, an accurate solution of the polynomial of order 2M-1 can be obtained. Since the first-arrival curve fitting in step 102) uses a low-order polynomial, such as 2 orders, M can be very small, such as M ═ 2, and an accurate solution of a 3-order polynomial can be obtained. Therefore, when the depth of the inflection point is calculated, the method has the advantages of high integration precision, small integration calculation amount and high calculation efficiency.

In the step 104), the selected M point coordinates may be at equal intervals or at unequal intervals. The M dots may cover the entire offset range or only a portion of the offset range.

And 105) interpolating into a three-dimensional velocity field according to the velocity and the depth of the corresponding inflection point of the inflection wave of each picked single shot first arrival selected azimuth angle.

Fig. 2 is a block diagram of a three-dimensional non-ray tracing refracted wave tomography apparatus according to the present invention. The device includes:

a first arrival curve obtaining unit 201, configured to pick up single shot first arrivals, select first arrivals at least two azimuth angles for each picked single shot first arrival according to azimuth angle intervals, and sort each azimuth angle first arrival in an increasing offset manner to form a first arrival curve in a range from 0 to x; wherein x is the maximum offset distance;

the fitting unit 202 is configured to perform low-order polynomial least square method overall fitting on the first arrival curve in the range from 0 to x by using the offset as a horizontal coordinate and the first arrival time as a vertical coordinate, so as to obtain a polynomial coefficient corresponding to the azimuth angle of each single shot first arrival;

an azimuth slowness curve obtaining unit 203, configured to obtain an azimuth slowness curve of the first arrival of each single shot correspondingly according to the polynomial coefficient;

a velocity and depth obtaining unit 204 for obtaining the velocity and depth of the inflection point of the echo corresponding to the selected azimuth angle of each single shot first arrival according to the azimuth angle slowness curve of each single shot first arrival;

and the interpolation unit 205 is used for interpolating a three-dimensional velocity field according to the velocity and the depth of the inflection point of the corresponding inflection wave of each picked single shot first arrival.

Preferably, the velocity and depth acquiring unit 204 for the inflection point of the inflection wave comprises: the device comprises an inflection wave inflection point acquisition module and a first calculation module; wherein,

the inflection wave inflection point acquisition module is used for arbitrarily taking M points as inflection points of the inflection wave for the azimuth angle slowness curve of the first arrival of each picked single cannon within the range of 0-x;

the first calculation module is used for calculating the reciprocal of the slowness curve value of the position corresponding to the point i to obtain the speed of the point i; wherein i is 1,2, …, M.

Preferably, the speed and depth obtaining unit 204 of the inflection point of the inflection wave further includes a second calculating module; wherein,

the second calculation module is used for obtaining the depth of the point i by adopting an integral formula; wherein the integral formula is:

$z_i=\frac{1}{\mathrm{\π}}{\∫}_0^{x_i}\cos h^{-1}\left(\frac{\frac{\mathrm{dt}}{\mathrm{dx}}}{{\left(\frac{\mathrm{dt}}{\mathrm{dx}}\right)}_i}\right)\mathrm{dx}$

wherein z is_iDepth, x, corresponding to inflection point i_iIs the offset corresponding to point i.

Preferably, the azimuth slowness curve acquiring unit 203 includes: the device comprises a polynomial acquisition module and a derivation module; wherein,

the polynomial acquisition module is used for acquiring a polynomial according to the polynomial coefficients;

and the derivation module is used for deriving the polynomial to obtain an azimuth slowness curve of the first arrival of each single cannon.

Preferably, the azimuth angle interval adopted by the first arrival curve acquiring unit 201 is a comparable or unequal azimuth angle interval.

For a better understanding of the invention, embodiments thereof are described in further detail below with reference to the accompanying drawings and examples of application.

Fig. 3 is a flowchart of the three-dimensional raytracing-free echo tomography method of the present embodiment. The specific implementation mode is as follows:

step 1) picking up single shot first arrivals and sorting according to azimuth angles: first arrival data are picked up, a grid speed model is established according to the work area range and the earth surface elevations of the shot and inspection points, and the initialization speed is 0. Sorting the picked single shot first arrivals according to azimuth angles, and sequencing the sorted first arrivals in an increasing mode according to offset distances to form a first arrival curve within the range of 0-3000 m, wherein 3000 m is the selected maximum offset distance.

Step 2) carrying out integral fitting on the first-motion curve of the specified azimuth angle: in a single gun, for the first arrival of the designated azimuth angle, a 2-order Legendre polynomial is adopted, the offset distance is used as a horizontal coordinate, the first arrival time is used as a vertical coordinate, the least square method is carried out to carry out overall fitting, and the polynomial coefficient corresponding to the azimuth angle is solved.

Step 3) calculating a slowness curve by using the fitted polynomial curve for the specified azimuth angle: calculating a slowness curve p ═ dt/dx within the range of 0-3000 m by using the polynomial coefficient obtained in the step 2) and according to a 2-order Legendre polynomial derivative formula and by taking the offset distance in the step 2) as an abscissa.

Step 4) calculating the corresponding speed and depth of the inflection point of the inflection wave of the azimuth angle: firstly, calculating the slowness obtained in the step 3) point by point within the range of 0-3000Reciprocal 1/(dt/dx) of curve p_i(i is 1,2, …, M), obtaining the velocity of the inflection point corresponding to the ith point; secondly, the 3-point Gauss-Legendre integral is adopted to calculate the inflection point depth of the point i. Using the depth calculation formula:

$z_i=\frac{1}{\mathrm{\π}}{\∫}_0^{x_i}\cos h^{-1}\left(\frac{\frac{\mathrm{dt}}{\mathrm{dx}}}{{\left(\frac{\mathrm{dt}}{\mathrm{dx}}\right)}_i}\right)\mathrm{dx}$

first, point-by-point derivative values are taken for points 0 to i and divided by a constant (dt/dx)_iCalculate cosh^-1Function values; secondly, looking up a mathematic manual to obtain a 3-point Gauss-Legendre polynomial root and weight coefficient between 0 and x_iRange weighted summation to obtain precise integral, and calculating depth z of inflection point_i。

The Gaussian-Legendre integral can be used for integrating equidistant nodes and non-equidistant nodes, so that the offset distance from the beginning to the corresponding point is taken as a calculation coordinate point, and the method can adapt to the situations of regular arrangement and irregular arrangement of detection points. When calculating, if only partial offset distance needs to be calculated, the integral range is limited.

And 5) repeating the steps 2) to 4) for each azimuth angle of the single shot to obtain the corresponding inflection point depth and the corresponding speed: at this time, 36 azimuths are selected, equally spaced by 10 degrees.

And 6) repeating the steps 2) to 5) for all the single guns to obtain all the inflection point depths and corresponding speeds.

And 7) interpolating the speed points into a three-dimensional speed field.

As shown in fig. 4, this embodiment is a tomographic effect chart. For the mountain data shown in fig. 4, the three-dimensional velocity field obtained by the invention gradually changes from the earth surface to the deep layer, is uniformly distributed in all directions in space, changes naturally, has no obvious abnormal body, is relatively uniform near the earth surface, and has no non-uniform calculation result and abnormal velocity distribution caused by the violent change of the earth surface in the conventional ray chromatography. In addition, when the result of the graph 4 is calculated, only 1 node and 1 thread on a cluster are used, the time is consumed for 3 minutes, while the conventional ray chromatography 8 nodes on the same cluster are parallel, each node adopts 10 threads, and the time is consumed for 30 minutes to obtain a similar result. The efficiency advantages of the present invention are therefore very significant.

The above-mentioned embodiments are intended to illustrate the objects, technical solutions and advantages of the present invention in further detail, and it should be understood that the above-mentioned embodiments are merely exemplary embodiments of the present invention, and are not intended to limit the scope of the present invention, and any modifications, equivalent substitutions, improvements and the like made within the spirit and principle of the present invention should be included in the scope of the present invention.

Claims (8)

1. A three-dimensional ray-tracing-free refracted wave tomography method is characterized by comprising the following steps:

picking up single cannon first arrivals, selecting the first arrivals of a plurality of azimuth angles for each picked single cannon first arrival according to azimuth angle intervals, and sequencing each azimuth angle first arrival according to an offset distance increasing mode to form a first arrival curve in a range of 0-x; wherein x is the maximum offset distance;

taking the offset as a horizontal coordinate and the first arrival time as a vertical coordinate, and performing low-order polynomial least square method integral fitting on the first arrival curve in the range of 0-x to obtain a polynomial coefficient corresponding to the azimuth angle of each single shot first arrival;

obtaining an azimuth slowness curve of the first arrival of each single shot correspondingly according to the polynomial coefficients;

obtaining the speed and the depth of a corresponding inflection point of the selected azimuth angle of each single shot first arrival according to the azimuth angle slowness curve of each single shot first arrival;

and interpolating into a three-dimensional velocity field according to the velocity and the depth of the corresponding inflection point of the inflection wave of each picked single shot first arrival selected azimuth angle.

2. The method of claim 1, wherein the velocity of the corresponding inflection point of each single shot first arrival picked is obtained by:

randomly taking M points as inflection points of the inflection waves within the range of 0-x for the azimuth slowness curve of the first arrival of each picked single cannon;

calculating the reciprocal of the slowness curve value of the position corresponding to the inflection point i to obtain the speed of the inflection point i; wherein i is 1,2, …, M.

3. The method of claim 2, wherein the depth of the inflection point corresponding to the selected azimuth angle for each single shot first arrival picked is obtained by:

obtaining the depth of the inflection point i by adopting an integral formula; wherein the integral formula is:

$z_i=\frac{1}{\mathrm{\π}}{\∫}_0^{x_i}{\cosh }^{-1}\left(\frac{\frac{dt}{dx}}{{\left(\frac{dt}{dx}\right)}_i}\right)dx$

wherein z is_iDepth, x, corresponding to inflection point i_iThe offset distance corresponding to the inflection point i; t denotes the first arrival time and d is the derivative operator.

4. The method of claim 1, wherein the azimuthal slowness curve for each single shot first arrival is obtained by:

obtaining a polynomial according to the polynomial coefficient;

and (4) carrying out derivation on the polynomial to obtain an azimuth slowness curve of the first arrival of each single shot.

5. A three-dimensional raytrace-free refracted-wave tomography apparatus, the apparatus comprising:

the system comprises a first arrival curve acquisition unit, a second arrival curve acquisition unit and a third arrival curve acquisition unit, wherein the first arrival curve acquisition unit is used for picking single cannon first arrivals, selecting first arrivals of at least two azimuth angles for each picked single cannon first arrival according to azimuth angle intervals, and sequencing each azimuth angle first arrival according to an offset distance increasing mode to form a first arrival curve in a range of 0-x; wherein x is the maximum offset distance;

the fitting unit is used for performing low-order polynomial least square method overall fitting on the first arrival curve in the range of 0-x by taking the offset distance as a horizontal coordinate and the first arrival time as a vertical coordinate to obtain a polynomial coefficient corresponding to the azimuth angle of each single shot first arrival;

the azimuth slowness curve acquisition unit is used for acquiring an azimuth slowness curve of the first arrival of each single shot correspondingly according to the polynomial coefficients;

the velocity and depth acquisition unit of the inflection point is used for acquiring the velocity and depth of the corresponding inflection point of the selected azimuth angle of each single shot first arrival according to the azimuth angle slowness curve of each single shot first arrival;

and the interpolation unit is used for interpolating a three-dimensional velocity field according to the velocity and the depth of the inflection point of the corresponding inflection wave of the selected azimuth angle of the first arrival of each picked single shot.

6. The apparatus of claim 5, wherein the velocity and depth of the inflection point acquisition unit comprises: the device comprises an inflection wave inflection point acquisition module and a first calculation module; wherein,

the inflection wave inflection point acquisition module is used for arbitrarily taking M points as inflection points of the inflection wave for the azimuth angle slowness curve of the first arrival of each picked single cannon within the range of 0-x;

the first calculation module is used for calculating the reciprocal of the slow curve value of the position corresponding to the inflection point i to obtain the speed of the inflection point i; wherein i is 1,2, …, M.

7. The apparatus of claim 6, wherein the velocity and depth of the inflection point acquisition unit further comprises a second calculation module; wherein,

the second calculation module is used for obtaining the depth of the inflection point i by adopting an integral formula; wherein the integral formula is:

$z_i=\frac{1}{\mathrm{\π}}{\∫}_0^{x_i}{\cosh }^{-1}\left(\frac{\frac{dt}{dx}}{{\left(\frac{dt}{dx}\right)}_i}\right)dx$

wherein z is_iDepth, x, corresponding to inflection point i_iT represents the first arrival time as the offset corresponding to the inflection point i.

8. The apparatus of claim 5, wherein the azimuth slowness curve acquisition unit comprises: the device comprises a polynomial acquisition module and a derivation module; wherein,

the polynomial acquisition module is used for acquiring a polynomial according to the polynomial coefficients;

and the derivation module is used for deriving the polynomial to obtain an azimuth slowness curve of the first arrival of each single cannon.