Seismic scattered wave field forward modeling method and device, medium, electronic equipment and product

Abstract

The invention relates to a seismic scattered wave field forward modeling method and device, a medium, electronic equipment and a product, and the method comprises the steps: dividing a speed model into a plurality of thin layers from the top of the model to the bottom of the model in a depth direction, acquiring downlink scattered wave data and primary uplink scattered wave data corresponding to each thin layer in a first propagation direction from the top of the model to the bottom of the model along a second propagation direction by using a high-order generalized screen continuation operator based on the Muller model, primary uplink scattered wave data reserved in the first direction is used as an initial value, primary uplink scattered wave data corresponding to each thin layer in the second propagation direction is obtained through the high-order generalized screen continuation operator, and the second propagation direction is the direction from the bottom of the model to the top of the model; and determining a forward modeling seismic record of the detector position according to the uplink scattered wave data corresponding to each thin layer in the second propagation direction. According to the method, the accuracy of forward modeling of the seismic scattered wave field in the viscoacoustic medium can be improved.

Description

Technical Field

[0001] The present disclosure relates to the field of geophysical exploration technology, and in particular to a seismic scattered wave field forward modeling method, device, medium, electronic equipment and product. Background Art

[0002] Seismic wave field forward modeling is an important tool for understanding and analyzing the propagation characteristics of seismic waves in underground media, and is also a prerequisite for studying seismic exploration and development technology. Generally, the process of describing seismic scattered waves based on scattering theory is mainly aimed at solving the Lippmann-Schwinger equation. In the prior art, the thin plate approximation and screen approximation methods can reveal the propagation law of scattered waves in viscoacoustic media (Futterman model). However, in the case of the presence of lateral strong scattering media in the viscoacoustic medium, the error of this method is relatively large. Summary of the invention

[0003] The purpose of the present disclosure is to provide a seismic scattered wave field forward modeling method, device, medium, electronic equipment and product to solve the above technical problems.

[0004] In order to achieve the above object, according to a first aspect of an embodiment of the present disclosure, a seismic scattered wave field forward modeling method is provided, the method comprising:

[0005] The velocity model is divided into multiple thin layers along the depth direction from the top to the bottom of the model;

[0006] Along the first propagation direction, a first scattered wave field calculation is performed using a high-order generalized screen extension operator based on the Müller model to obtain downlink scattered wave data and one uplink scattered wave data corresponding to each thin layer in the first propagation direction; the first propagation direction is the direction from the top of the model to the bottom of the model;

[0007] Along a second propagation direction, a second scattered wave field calculation is performed using the high-order generalized screen extension operator based on the Müller model to obtain upgoing scattered wave data corresponding to each thin layer in the second propagation direction; the second propagation direction is the direction from the bottom of the model to the top of the model, and the second propagation direction is opposite to the first propagation direction;

[0008] The forward seismic record of the detector position is determined according to the upgoing scattered wave data corresponding to each thin layer in the second propagation direction.

[0009] Optionally, the number of the thin layers is N, where N is an integer greater than 1, and the first scattering wave field calculation is performed along the first propagation direction using a high-order generalized screen extension operator based on a Müller model, including:

[0010] Along the first propagation direction, the downlink scattered wave data of the previous thin layer is used as the boundary condition of the next adjacent thin layer, and the operation of calculating the scattered wave field by using the high-order generalized screen extension operator based on the Müller model is iteratively performed to obtain the downlink scattered wave data and the uplink scattered wave data corresponding to each thin layer in the N thin layers;

[0011] The iteration is stopped when the downlink scattered wave data and one uplink scattered wave data corresponding to the Nth thin layer are obtained; wherein the Nth thin layer is a thin layer located at the bottom of the model.

[0012] Optionally, the method comprises:

[0013] The source parameters are used as the boundary conditions of the first thin layer among N thin layers, and the scattered wave field is calculated using a high-order generalized screen continuation operator based on the Müller model to obtain the downlink scattered wave data and the first uplink scattered wave data corresponding to the second thin layer; wherein the first thin layer is the thin layer located at the top of the model among the N thin layers.

[0014] Optionally, the method of using the source parameter as the boundary condition of the first thin layer among the N thin layers and using a high-order generalized screen extension operator based on the Müller model to calculate the scattered wave field includes:

[0015] A two-dimensional coordinate system is established along the tunnel; wherein the origin of the two-dimensional coordinate system is the center of the tunnel surface on the velocity model boundary, the x-axis direction is the horizontal direction, and the z-axis direction is the depth direction;

[0016] Acquire the earthquake source parameters; the earthquake source parameters include one or more of earthquake source location parameters, main frequency, sampling interval, sampling time length, number of detectors and detector location;

[0017] The source parameters are used as boundary conditions of the first thin layer, and the scattered wave field is calculated using a high-order generalized screen extension operator based on the Müller model.

[0018] Optionally, the number of thin layers is N, where N is an integer greater than 1, and the second scattered wave field calculation is performed along the second propagation direction using the high-order generalized screen continuation operator based on the Müller model to obtain primary upgoing scattered wave data corresponding to each thin layer in the second propagation direction, including:

[0019] Along the second propagation direction, the first upward scattered wave of the N-1th thin plate in the first propagation direction is used as the initial value condition in the second propagation direction, and the operation of performing scattered wave field calculation using the high-order generalized screen extension operator based on the Müller model is iteratively performed, wherein the first upward scattered wave field value retained at the same position in the first propagation direction is added to the calculation result of each scattered wave field calculation to obtain the first upward scattered wave data corresponding to each thin layer in the N thin layers;

[0020] The iteration is stopped when the up-going scattered wave data corresponding to the Mth thin layer is obtained; wherein the Mth thin layer is the thin layer located at the top of the model; and the values ​​of M and N are the same.

[0021] Optionally, along the second propagation direction, taking the first upward scattered wave of the N-1th thin plate in the first propagation direction as the initial value condition in the second propagation direction, and iteratively performing the operation of calculating the scattered wave field using the high-order generalized screen extension operator based on the Müller model, comprises:

[0022] Taking the N-1th thin layer in the first propagation direction as the first thin layer in the second propagation direction;

[0023] Acquire the primary upgoing scattered wave data corresponding to the N-1th thin layer in the first propagation direction;

[0024] The first upward scattered wave data corresponding to the N-1th thin layer is used as the boundary condition of the first thin layer in the second propagation direction, and the scattered wave field is calculated using a high-order generalized screen extension operator based on the Müller model. The first upward scattered wave data retained by the N-2th thin layer in the first propagation direction is added to the calculation result to obtain the first upward scattered wave data corresponding to the second thin layer in the second propagation direction.

[0025] Optionally, the Müller model calculation formula used in the scattered wave field calculation is:

[0026]

[0027] Where ω is the frequency, v 0 (ω) is the wave velocity at frequency ω, Q 0 is the background quality factor, ω r is the reference frequency, γ is the frequency-dependent attenuation index, γ=1 / πQ 0 ;

[0028] The calculation formula for the downlink scattered wave or the uplink scattered wave in the scattered wave field calculation is:

[0029]

[0030] Among them, p(x,z j+1 ω) represents the down-going scattered wave corresponding to the next thin layer or an up-going scattered wave, p(x,z j ;ω) represents the downward scattered wave or the upward scattered wave corresponding to the previous thin layer, FT[] is the Fourier transform, FT -1 [] is the inverse Fourier transform, Δz is the extension step, A is the coefficient of the high-order generalized screen extension operator, which is used to correct the propagation characteristics of the wave field, ω is the circular frequency, is the vertical wave number; k x is the transverse wave number;

[0031]

[0032] Where n is the order, s 0 (x,z j ;ω) is the background velocity field, s(x,z j ; ω) is the actual velocity field, represents the background velocity field, which is used to represent the velocity distribution in the medium without disturbance; v(x,z j ; ω) represents the actual velocity field, which contains disturbances or changes outside the background velocity field;

[0033] The calculation formula of an upward scattered wave in the scattered wave field calculation is:

[0034] b(x,z j+1 ;ω)=R*p(x,z j+1 ;ω)

[0035] Among them, b(x,z j+1 ; ω) represents the once upward scattered wave, and R is the reflection coefficient.

[0036] Optionally, determining the forward seismic record according to the first upward scattered wave data corresponding to each thin layer in the first propagation direction and the first upward scattered wave data corresponding to each thin layer in the second propagation direction includes:

[0037] Acquire total first-order upward scattered wave data in the first propagation direction according to the downward scattered wave data corresponding to each thin layer in the first propagation direction;

[0038] Acquire total upward scattered wave data in the second propagation direction according to the first upward scattered wave data corresponding to each thin layer in the first propagation direction and the first upward scattered wave data corresponding to each thin layer in the second propagation direction;

[0039] A forward seismic record is determined based on the total upgoing scattered wave data in the second propagation direction.

[0040] According to a second aspect of an embodiment of the present disclosure, a seismic scattered wave field forward modeling device is provided, the device comprising:

[0041] A model partitioning module is used to divide the velocity model into multiple thin layers from the top to the bottom of the model along the depth direction;

[0042] A first calculation module is used to perform a first scattered wave field calculation along a first propagation direction using a high-order generalized screen extension operator based on a Müller model to obtain downlink scattered wave data and one uplink scattered wave data corresponding to each thin layer in the first propagation direction; the first propagation direction is a direction from the top of the model to the bottom of the model;

[0043] A second calculation module is used to perform a second scattered wave field calculation along a second propagation direction using the high-order generalized screen extension operator based on the Müller model to obtain the first upward scattered wave data corresponding to each thin layer in the second propagation direction; the second propagation direction is the direction from the bottom of the model to the top of the model, and the second propagation direction is opposite to the first propagation direction;

[0044] A determination module is used to determine the forward seismic record of the detector position according to the first upgoing scattered wave data corresponding to each thin layer in the second propagation direction.

[0045] According to a third aspect of an embodiment of the present disclosure, a non-temporary computer-readable storage medium is provided, on which a computer program is stored, and when the program is executed by a processor, the steps of any one of the methods described in the first aspect of the embodiment of the present disclosure are implemented.

[0046] According to a fourth aspect of an embodiment of the present disclosure, there is provided an electronic device, including:

[0047] a memory having a computer program stored thereon;

[0048] A processor is used to execute the computer program in the memory to implement the steps of any one of the methods described in the first aspect of the embodiments of the present disclosure.

[0049] According to a fifth aspect of an embodiment of the present disclosure, a computer program product is provided, including a computer program, which, when executed by a processor, implements the steps of any one of the methods described in the first aspect of the embodiment of the present disclosure.

[0050] In the above technical solution, the velocity model is divided into a plurality of thin layers along the depth direction from the top to the bottom of the model, and a first scattering wave field calculation is performed along a first propagation direction using a high-order generalized screen extension operator based on a Müller model to obtain downlink scattering wave data and first uplink scattering wave data corresponding to each thin layer in the first propagation direction; the first propagation direction is the direction from the top to the bottom of the model, and along a second propagation direction, the first uplink wave retained in the first direction is used as an initial value, and a second scattering wave field calculation is performed using the high-order generalized screen extension operator based on the Müller model to obtain uplink scattering wave data corresponding to each thin layer in the second propagation direction; the second propagation direction is the direction from the bottom to the top of the model, and the second propagation direction is opposite to the first propagation direction, and the forward seismic record of the detector position is determined according to the first uplink scattering wave data corresponding to each thin layer in the second propagation direction. Through the above technical solution, the high-order generalized screen extension operator is combined with the Müller model, which can improve the ability of the operator to adapt to the lateral strong scattering medium, and thus improve the accuracy of seismic wave field forward simulation when there is a lateral strong scattering medium in the medium.

[0051] Other features and advantages of the present disclosure will be described in detail in the following detailed description. BRIEF DESCRIPTION OF THE DRAWINGS

[0052] The accompanying drawings are used to provide a further understanding of the present disclosure and constitute a part of the specification. Together with the following specific embodiments, they are used to explain the present disclosure but do not constitute a limitation of the present disclosure. In the accompanying drawings:

[0053] Figure 1 is a flow chart of a seismic scattered wave field forward modeling method according to an exemplary embodiment;

[0054] Figure 2 is a schematic diagram of velocity model division according to an exemplary embodiment;

[0055] Figure 3 is a schematic diagram of a velocity model according to an exemplary embodiment;

[0056] Figure 4 is a block diagram of a seismic scattered wave field forward modeling device according to an exemplary embodiment;

[0057] Figure 5 is a block diagram of an electronic device 500 according to an exemplary embodiment. DETAILED DESCRIPTION

[0058] The specific implementation of the present disclosure is described in detail below in conjunction with the accompanying drawings. It should be understood that the specific implementation described herein is only used to illustrate and explain the present disclosure, and is not used to limit the present disclosure.

[0059] It should be noted that all actions of acquiring signals, information or data in the present disclosure are performed with authorization given by the owner of the corresponding device.

[0060] Among the research methods of seismic wave propagation characteristics, the thin plate approximation and screen approximation are methods of alternating calculations in the frequency-wavenumber domain and the frequency-space domain. Phase shift is completed in the frequency-wavenumber domain, and time shift correction is completed in the frequency-space domain. The two are added to obtain the total downlink scattering wave field, and the backscattering field is calculated once. The total downlink scattering field is used as the initial condition for the next thin plate, and so on, until the bottom of the model. The first reflection field at the bottom of the model is used as the initial value of the upgoing wave and extended. The previously saved backscattering field needs to be added during the calculation process until the top of the model is received by the detector. However, the applicant found that the prior art cannot handle the situation of strong lateral scattering media during implementation. When there is strong disturbance locally, the compensation accuracy of the operator is insufficient. The reasons include: (1) The local Born approximation can only adapt to weak scattering conditions. When the disturbance in the thin plate is small, the conventional extension operator can handle it normally. However, when the lateral change is strong or the frequency is high (not meeting the Born approximation condition), the correction term in the spatial domain cannot correctly compensate for the error; (2) The background medium Green function is approximated by the uniform medium Green function. In order to improve the calculation efficiency, the thin plate approximation and the screen approximation both use the local uniform medium Green function value to approximate the renormalized background Green function in the process of calculating the disturbance. However, the condition of the approximation requires that the disturbance in the thin plate is small. Therefore, the errors of these two approximation methods are relatively large in strong scattering media. In order to solve the problems in the prior art, the present invention provides a seismic scattering wave field forward modeling method, device, medium, electronic equipment and product, which are described below in conjunction with specific embodiments.

[0061] Figure 1 is a flow chart of a seismic scattering wave field forward modeling method according to an exemplary embodiment. Figure 1 As shown, the method includes the following steps.

[0062] In step S11 , the velocity model is divided into a plurality of thin layers along the depth direction from the top to the bottom of the model.

[0063] In the forward modeling of seismic scattering wave fields, the velocity model refers to a three-dimensional model that describes the velocity distribution of the underground medium. The velocity model is the basis for the simulation of seismic wave propagation and can define how seismic waves propagate in the underground medium, including the refraction, reflection and diffraction of the waves. The velocity model can include information about the underground rock type, the state of the pore fluid, and other factors that may affect the velocity of seismic waves. When performing forward modeling, the accuracy of the velocity model will affect the credibility of the simulation results. Therefore, it is necessary to construct a velocity model with higher accuracy during seismic exploration and research.

[0064] In step S12, along the first propagation direction, a first scattered wave field calculation is performed using a high-order generalized screen extension operator based on the Müller model to obtain downlink scattered wave data and one uplink scattered wave data corresponding to each thin layer in the first propagation direction; the first propagation direction is the direction from the top of the model to the bottom of the model.

[0065] In step S13, along the second propagation direction, the second scattered wave field calculation is performed using the high-order generalized screen extension operator based on the Müller model to obtain the first upward scattered wave data corresponding to each thin layer in the second propagation direction; the second propagation direction is the direction from the bottom of the model to the top of the model, and the second propagation direction is opposite to the first propagation direction.

[0066] Among them, the high-order generalized screen extension operator based on the Müller model is a mathematical tool used in seismic wave field simulation. This extension operator can more accurately describe the propagation characteristics of waves in complex geological structures and can improve the limitations of existing traditional methods in dealing with wavefront diffraction and scattering phenomena.

[0067] In step S14, the forward seismic record of the detector position is determined according to the primary upgoing scattered wave data corresponding to each thin layer in the second propagation direction.

[0068] For example, in a possible implementation, when performing seismic wave field forward modeling, a known velocity model is selected, which describes the velocity distribution of the underground medium. For the convenience of subsequent calculations, the velocity model is evenly divided into multiple thin layers from the top to the bottom along the depth direction, and each thin layer represents a homogeneous medium within a certain depth range. Such a division method can ensure the accuracy of the calculation and improve the calculation efficiency. Along the first propagation direction from the top to the bottom of the model, a high-order generalized screen extension operator based on the Müller model is used to perform the first scattered wave field calculation, wherein the first propagation direction is from the surface (top of the model) to the bottom of the ground (bottom of the model) along the depth direction. Downward scattered wave data and one upward scattered wave data corresponding to each thin layer can be obtained through the first scattered wave field calculation. Then, along the second propagation direction from the bottom to the top of the model, the high-order generalized screen extension operator based on the Müller model is used again to perform calculations, thereby obtaining the upward scattered wave data corresponding to each thin layer. Among them, the second propagation direction is from the ground (bottom of the model) upward back to the surface (top of the model). By combining the data obtained from the two calculations, namely the downlink scattered wave data and the uplink scattered wave data in the first propagation direction, and the uplink scattered wave data in the second propagation direction, the forward seismic record can be determined.

[0069] It can be understood that the total downlink scattered wave data in the first propagation direction is obtained according to the downlink scattered wave data corresponding to each thin layer in the first propagation direction; the total uplink scattered wave data in the second propagation direction is obtained according to the uplink scattered wave data corresponding to each thin layer in the first propagation direction and the uplink scattered wave data corresponding to each thin layer in the second propagation direction; the forward seismic record is determined according to the total uplink scattered wave data in the first propagation direction and the total uplink scattered wave data in the second propagation direction. The forward seismic record can provide a behavior pattern of seismic wave propagation in a simulated velocity model, which can be used for subsequent data analysis and geological structure interpretation.

[0070] In the above technical solution, the velocity model is divided into a plurality of thin layers along the depth direction from the top to the bottom of the model, and a first scattering wave field calculation is performed along a first propagation direction using a high-order generalized screen extension operator based on a Müller model to obtain downlink scattering wave data and first uplink scattering wave data corresponding to each thin layer in the first propagation direction; the first propagation direction is a direction from the top to the bottom of the model, and a second scattering wave field calculation is performed along a second propagation direction using the high-order generalized screen extension operator based on a Müller model to obtain uplink scattering wave data corresponding to each thin layer in the second propagation direction; the second propagation direction is a direction from the bottom to the top of the model, and the second propagation direction is opposite to the first propagation direction, and a forward seismic record is determined based on the downlink scattering wave data and first uplink scattering wave data corresponding to each thin layer in the first propagation direction, and the uplink scattering wave data corresponding to each thin layer in the second propagation direction. Through the above technical solution, the high-order generalized screen extension operator is combined with the Müller model, which can improve the ability of the operator to adapt to the lateral strong scattering medium, and thus improve the accuracy of seismic wave field forward simulation when there is a lateral strong scattering medium in the medium.

[0071] Optionally, in a possible implementation, the velocity model is divided into a plurality of thin layers, the number of the thin layers is N, and N is an integer greater than 1. The calculation of the first scattered wave field along the first propagation direction by using a high-order generalized screen extension operator based on the Müller model in the above step S12 may include:

[0072] A1. Along the first propagation direction, the downward scattered wave data of the previous thin layer is used as the boundary condition of the next adjacent thin layer, and the operation of calculating the scattered wave field using the high-order generalized screen extension operator based on the Müller model is iteratively performed to obtain the downward scattered wave data and the upward scattered wave data corresponding to each thin layer in the N thin layers.

[0073] Optionally, the source parameters are used as the boundary conditions of the first thin layer among the N thin layers, and the scattering wave field calculation is performed using the high-order generalized screen extension operator based on the Müller model to obtain the downlink scattering wave data and the first uplink scattering wave data corresponding to the second thin layer. The first thin layer is the thin layer located at the top of the model among the N thin layers.

[0074] It can be understood that based on Huygens' principle, when the wave propagates forward, when it propagates to a new position, it continues to propagate forward as a new source. Therefore, taking the downward scattered wave as the initial value of the next thin layer is equivalent to loading the "source term" on each thin layer. The difference is that the data at this time is different from the source term. The boundary condition of the first thin layer among the N thin layers is the source parameter.

[0075] A2. Stop the iteration when the downlink scattered wave data and one uplink scattered wave data corresponding to the Nth thin layer are obtained; wherein the Nth thin layer is the thin layer located at the bottom of the model.

[0076] Optionally, the step S13 described above, performing the second scattered wave field calculation along the second propagation direction using the high-order generalized screen extension operator based on the Müller model to obtain the upgoing scattered wave data corresponding to each thin layer in the second propagation direction, may include:

[0077] B1. Along the second propagation direction, the first upward scattered wave of the N-1th thin plate in the first propagation direction is used as the initial value condition of the second propagation direction, and the operation of calculating the scattered wave field using the high-order generalized screen extension operator based on the Müller model is iteratively performed, wherein the first upward scattered wave field value retained at the same position in the first propagation direction is added to the calculation result of each scattered wave field calculation to obtain the first upward scattered wave data corresponding to each thin layer in the N thin layers.

[0078] Optionally, in a possible implementation, the iterative step may include:

[0079] The N-1th thin layer in the first propagation direction is taken as the first thin layer in the second propagation direction, and the first upward scattered wave data corresponding to the N-1th thin layer in the first propagation direction is obtained. The first upward scattered wave data corresponding to the N-1th thin layer is used as the boundary condition of the first thin layer in the second propagation direction. The scattered wave field is calculated using a high-order generalized screen extension operator based on the Müller model, and the first upward scattered wave data retained from the N-2th thin layer in the first propagation direction are added to the calculation result to obtain the first upward scattered wave data corresponding to the second thin layer in the second propagation direction, and then the downward scattered wave data of the previous thin layer is used as the boundary condition of the next adjacent thin layer for iteration.

[0080] B2. Stop the iteration when the uplink scattered wave data corresponding to the Mth thin layer is obtained; wherein the Mth thin layer is the thin layer located at the top of the model; and the values ​​of M and N are the same.

[0081] Optionally, the step of using the source parameters as the boundary condition of the first thin layer among the N thin layers and using a high-order generalized screen extension operator based on the Müller model to calculate the scattered wave field may include:

[0082] a) Establish a two-dimensional coordinate system along the tunnel; wherein the origin of the two-dimensional coordinate system is the center of the tunnel surface on the velocity model boundary, the x-axis direction is the horizontal direction, and the z-axis direction is the depth direction.

[0083] For example, Figure 2is a schematic diagram of velocity model division according to an exemplary embodiment. Figure 2 As shown in the figure, in one possible division method, the velocity model is divided into i thin layers along the z-axis direction. After the scattered wave field is calculated using the high-order generalized screen extension operator based on the Müller model, all the downlink scattered wave data p in the first propagation direction can be obtained. fi , and all the uplink scattered wave data p in the second propagation direction bi .

[0084] b) Acquiring the source parameters, wherein the source parameters may include one or more of the source location parameters, main frequency, sampling interval, sampling time length, number of detectors and detector location.

[0085] c) The source parameters are used as the boundary conditions of the first thin layer, and the scattering wave field is calculated using a high-order generalized screen extension operator based on the Müller model.

[0086] For example, the Müller model calculation formula used in the scattered wave field calculation is:

[0087]

[0088] Where ω is the frequency, v 0 (ω) is the wave velocity at frequency ω, Q 0 is the background quality factor, ω r is the reference frequency, γ is the frequency-dependent attenuation index, γ=1 / πQ 0 ;

[0089] The calculation formula for the downlink scattered wave or uplink scattered wave in the scattered wave field calculation is:

[0090]

[0091] Among them, p(x,z j+1 ; ω) represents the down-going scattered wave or up-going scattered wave corresponding to the next thin layer, p(x,z j ;ω) represents the downward or upward scattered wave corresponding to the previous thin layer, FT[] is the Fourier transform, FT -1 [] is the inverse Fourier transform, Δz is the extension step, A is the coefficient of the high-order generalized screen extension operator, which is used to correct the propagation characteristics of the wave field, ω is the circular frequency, is the vertical wave number; k x is the transverse wave number;

[0092]

[0093] Where n is the order, s 0 (x,zj ;ω) is the background velocity field, s(x,z j ; ω) is the actual velocity field, represents the background velocity field, which is used to represent the velocity distribution in the medium without disturbance; v(x,z j ; ω) represents the actual velocity field, which contains disturbances or changes outside the background velocity field;

[0094] The calculation formula for an upward scattered wave in the scattered wave field calculation is:

[0095] b(x,z j+1 ;ω)=R*p(x,z j+1 ;ω) (Formula 4)

[0096] Among them, b(x,z j+1 ; ω) represents the uplink scattered wave, and R is the reflection coefficient.

[0097] For example, the seismic scattered wave field forward modeling method provided by the present disclosure may include the following steps:

[0098] S1. Obtain the velocity model v(x,z) and divide it into thin layers along the depth direction, with a total of j layers, where the number of lateral parameters is nx, the number of longitudinal depths is nz, and the grid spacings are dx and dz respectively.

[0099] For example, in a possible implementation, the number of lateral parameters is 401, the number of longitudinal depth parameters is 301, and the grid spacings are 5 m and 5 m respectively.

[0100] S2. Calculate the corresponding quality factor Q parameter through the velocity model. The calculation formula is as follows:

[0101] Q=14×v 2.2 (x,z)

[0102] The unit of v(x,z) is km / s.

[0103] S3. Establish a two-dimensional coordinate system along the tunnel, with the x direction as the horizontal direction and the z direction as the depth direction. The two-dimensional coordinate system is established with the center of the tunnel surface on the model boundary as the origin. The source position parameters are (sx, sz), the main frequency is fHz, the sampling interval is set to 1 ms, the sampling time length is ts, there are N detectors in total, the position of the first detector is (gx, gz), and the track spacing is dm.

[0104] For example, in one possible implementation, the source location parameters are (1.0km, 0.0), the main frequency is 25Hz, the sampling interval is set to 1ms, the sampling time length is 1500s, there are 200 detectors in total, the first detector position is (0.5km, 0.0), and the track spacing is 5m. Figure 3 is a schematic diagram of a velocity model according to an exemplary embodiment. Figure 3 As shown, 1 represents the source position and 2 represents the position of the detector.

[0105] S4. Perform edge expansion processing on the model parameters, use the Hanning window function as the absorbing boundary processing, expand the edge to px grids (for example, it can be 100 grids), and the expression of the window function is as follows:

[0106]

[0107] Where r = (x, z), r 0 =(sx+px,sz+px).

[0108] S5. Using the source parameters as the boundary conditions of the first thin layer, the high-order generalized screen extension operator based on the Müller model is used to calculate the downlink scattered wave p(x,z j+1 ;ω) and the first up-scattered wave b(x,z j+1 ;ω), and the first up-scattered wave b(x,z j+1 ;ω) for subsequent calculation. Among them, the downlink scattered wave p(x,z j+1 ;ω) and the first up-scattered wave b(x,z j+1 ; ω), are calculated by the above formula (2) and formula (4), respectively, and the expression of the source term can be:

[0109]

[0110] This formula transforms the Ricker wavelet into the frequency domain through Fourier transform, where f is the main frequency of the Ricker wavelet, ω is the circular frequency, and exp() is the e exponent.

[0111] S6. Use the value of the downward scattered wave as the boundary condition of the next thin layer, iterate step S5, and obtain the downward scattered wave data and the upward scattered wave data of each thin layer respectively until the bottom of the model.

[0112] S7. Taking the first up-scattering wave field data of the j-th layer as the boundary condition of the first thin layer in the propagation direction toward the surface, iteratively calculate the up-scattering wave field data until the surface, thereby obtaining all the forward seismic records.

[0113] It can be understood that the total upward scattered wave field data in the propagation direction toward the surface is the sum of the primary upward scattered wave data corresponding to each thin layer in the first propagation direction and the upward scattered wave data corresponding to each thin layer in the second propagation direction.

[0114] In addition, it should be noted that in the calculation process of the scattered wave field, the velocity in each thin layer needs to be decomposed into background velocity and disturbance velocity, that is, v(x,z)=v 0 (x,z)+v per (x,z), where v 0 (x,z) is the average velocity of each thin layer, Accordingly, the quality factor is divided into background quality factor and disturbance quality factor parameters, namely

[0115] This processing is a requirement of multi-domain seismic migration technology itself. Multi-domain seismic migration technology includes split-step Fourier method, phase screen method, extended local Born approximation method, extended local Rytov approximation method and generalized screen migration method. The basic idea of ​​this type of algorithm is to decompose the velocity field into the constant velocity background and variable velocity disturbance within the layer. The background field is equivalent to solving the constant velocity acoustic wave equation, which can be achieved through phase shift. For variable velocity disturbance, this inhomogeneity can be considered to be equivalent to a scattering source (secondary source). The incident wave field acts on these scattering sources, thereby generating a scattered wave field. These algorithms are implemented by first converting the problem to the frequency domain, and then using the fast Fourier transform (FFT) to convert and calculate between the space domain and the wave number domain to realize multi-domain (space domain and wave number domain) operations. These methods generally include three items in terms of operators: free propagation, phase screen and high angle compensation. Free propagation is carried out in the wavenumber domain through a uniform medium with a certain reference velocity, which can vary with depth; the other two terms generally reflect the interaction between the wave and the inhomogeneous part of the medium, the phase screen term is carried out in the space domain, and high-angle propagation compensation can be carried out in the wavenumber domain or by finite difference in the frequency-space domain. The main advantage of multi-domain migration technology is that it does not produce discrete grid dispersion phenomenon, and it has high computational efficiency and strong adaptability, and can perform high-angle and lateral strong variable speed imaging. The multi-domain migration operator can be regarded as an ultra-wide-angle one-way wave propagation operator, ignoring the reverberation of round-trip propagation between inhomogeneities, but can correctly handle wave phenomena related to multiple downlink scattering, such as wave focusing, diffraction, refraction and interference.

[0116] In the above technical solution, the velocity model is divided into a plurality of thin layers along the depth direction from the top to the bottom of the model, and a first scattering wave field calculation is performed along a first propagation direction using a high-order generalized screen extension operator based on a Müller model to obtain downlink scattering wave data and first uplink scattering wave data corresponding to each thin layer in the first propagation direction; the first propagation direction is a direction from the top to the bottom of the model, and a second scattering wave field calculation is performed along a second propagation direction using the high-order generalized screen extension operator based on a Müller model to obtain uplink scattering wave data corresponding to each thin layer in the second propagation direction; the second propagation direction is a direction from the bottom to the top of the model, and the second propagation direction is opposite to the first propagation direction, and a forward seismic record is determined based on the downlink scattering wave data and first uplink scattering wave data corresponding to each thin layer in the first propagation direction, and the uplink scattering wave data corresponding to each thin layer in the second propagation direction. Through the above technical solution, the high-order generalized screen extension operator is combined with the Müller model, which can improve the ability of the operator to adapt to the lateral strong scattering medium, and thus improve the accuracy of seismic wave field forward simulation when there is a lateral strong scattering medium in the medium.

[0117] Figure 4 is a block diagram of a seismic scattering wave field forward modeling device according to an exemplary embodiment. Figure 4 As shown, the device 400 includes: a model division module 410, a first calculation module 420, a second calculation module 430 and a determination module 440.

[0118] The model division module 410 is used to divide the velocity model into a plurality of thin layers along the depth direction from the top of the model to the bottom of the model;

[0119] The first calculation module 420 is used to perform a first scattered wave field calculation along a first propagation direction using a high-order generalized screen extension operator based on a Müller model to obtain downlink scattered wave data and one uplink scattered wave data corresponding to each thin layer in the first propagation direction; the first propagation direction is the direction from the top of the model to the bottom of the model;

[0120] The second calculation module 430 is used to perform a second scattered wave field calculation along a second propagation direction using the high-order generalized screen extension operator based on the Müller model to obtain the primary upgoing scattered wave data corresponding to each thin layer in the second propagation direction; the second propagation direction is the direction from the bottom of the model to the top of the model, and the second propagation direction is opposite to the first propagation direction;

[0121] The determination module 440 is used for determining the forward seismic record of the detector position based on the uplink scattered wave data corresponding to each thin layer in the second propagation direction.

[0122] Optionally, the number of the thin layers is N, where N is an integer greater than 1, and the first calculation module 420 is further configured to:

[0123] Along the first propagation direction, the downlink scattered wave data of the previous thin layer is used as the boundary condition of the next adjacent thin layer, and the operation of calculating the scattered wave field by using the high-order generalized screen extension operator based on the Müller model is iteratively performed to obtain the downlink scattered wave data and the uplink scattered wave data corresponding to each thin layer in the N thin layers;

[0124] The iteration is stopped when the downlink scattered wave data and one uplink scattered wave data corresponding to the Nth thin layer are obtained; wherein the Nth thin layer is the thin layer located at the bottom of the model.

[0125] Optionally, the first calculation module 420 is further configured to:

[0126] The source parameters are used as the boundary conditions of the first thin layer among N thin layers, and the scattered wave field is calculated using the high-order generalized screen continuation operator based on the Müller model to obtain the downgoing scattered wave data and the first upgoing scattered wave data corresponding to the second thin layer; wherein the first thin layer is the thin layer located at the top of the model among the N thin layers.

[0127] Optionally, the first calculation module 420 includes a coordinate system establishment submodule, an acquisition submodule and a calculation submodule;

[0128] The coordinate system establishment submodule is used to establish a two-dimensional coordinate system along the tunnel; wherein the origin of the two-dimensional coordinate system is the center of the tunnel surface on the velocity model boundary, the x-axis direction is the horizontal direction, and the z-axis direction is the depth direction;

[0129] The acquisition submodule is used to acquire the source parameters; the source parameters include one or more of the source location parameters, main frequency, sampling interval, sampling time length, number of detectors and detector location;

[0130] The calculation submodule is used to use the source parameter as the boundary condition of the first thin layer and perform scattered wave field calculation using a high-order generalized screen extension operator based on the Müller model.

[0131] Optionally, the number of the thin layers is N, where N is an integer greater than 1, and the second calculation module 430 is further used for:

[0132] Along the second propagation direction, the first upward scattered wave of the N-1th thin plate in the first propagation direction is used as the initial value condition of the second propagation direction, and the first upward scattered wave field value of the same position retained in the first propagation direction is added to the calculation result each time, and the operation of calculating the scattered wave field using the high-order generalized screen extension operator based on the Müller model is iteratively performed to obtain the first upward scattered wave data corresponding to each thin layer in the N thin layers;

[0133] The iteration is stopped when the up-going scattered wave data corresponding to the Mth thin layer is obtained; wherein the Mth thin layer is the thin layer located at the top of the model; and the values ​​of M and N are the same.

[0134] Optionally, the second calculation module 430 is further used for:

[0135] Using the N-1th thin layer in the first propagation direction as the first thin layer in the second propagation direction;

[0136] Acquire the primary upgoing scattered wave data corresponding to the N-1th thin layer in the first propagation direction;

[0137] The first upward scattered wave data corresponding to the N-1th thin layer is used as the boundary condition of the first thin layer in the second propagation direction, and the scattered wave field is calculated using a high-order generalized screen continuation operator based on the Müller model. The first upward scattered wave data retained by the N-2th thin layer in the first propagation direction is added to the calculation result to obtain the upward scattered wave data corresponding to the second thin layer in the second propagation direction.

[0138] Optionally, the Müller model calculation formula used in the scattered wave field calculation is:

[0139]

[0140] Where ω is the frequency, v 0 (ω) is the wave velocity at frequency ω, Q 0 is the background quality factor, ω r is the reference frequency, γ is the frequency-dependent attenuation index, γ=1 / πQ 0 ;

[0141] The calculation formula for the downlink scattered wave or the uplink scattered wave in the scattered wave field calculation is:

[0142]

[0143] Among them, p(x,z j+1 ω) represents the down-going scattered wave corresponding to the next thin layer or an up-going scattered wave, p(x,z j; ω) represents the downward scattered wave or the upward scattered wave corresponding to the previous thin layer, FT[] is the Fourier transform, FT -1 [] is the inverse Fourier transform, Δz is the extension step, A is the coefficient of the high-order generalized screen extension operator, which is used to correct the propagation characteristics of the wave field, ω is the circular frequency, is the vertical wave number; k x is the transverse wave number;

[0144]

[0145] Where n is the order, s 0 (x,z j ;ω) is the background velocity field, s(x,z j ; ω) is the actual velocity field, v 0 (x,z j ; ω) represents the background velocity field, which is used to represent the velocity distribution in the medium without disturbance; v(x,z j ; ω) represents the actual velocity field, which contains disturbances or changes outside the background velocity field;

[0146] The calculation formula for an upward scattered wave in the scattered wave field calculation is:

[0147] b(x,z j+1 ;ω)=R*p(x,z j+1 ;ω)

[0148] Among them, b(x,z j+1 ; ω) represents the uplink scattered wave, and R is the reflection coefficient.

[0149] Optionally, the determining module 440 is further configured to:

[0150] Acquire total downward scattered wave data in the first propagation direction according to downward scattered wave data corresponding to each thin layer in the first propagation direction;

[0151] Acquire total upward scattered wave data in the second propagation direction according to the primary upward scattered wave data corresponding to each thin layer in the first propagation direction and the upward scattered wave data corresponding to each thin layer in the second propagation direction;

[0152] A forward seismic record is determined based on the total upgoing scattered wave data in the second propagation direction.

[0153] Regarding the device in the above embodiment, the specific manner in which each module performs operations has been described in detail in the embodiment of the method, and will not be elaborated here.

[0154] Figure 5FIG. 5 is a block diagram of an electronic device 500 according to an exemplary embodiment. Figure 5 As shown, the electronic device 500 may include: a processor 501 and a memory 502. The electronic device 500 may also include one or more of a multimedia component 503, an input / output (I / O) interface 504, and a communication component 505.

[0155] The processor 501 is used to control the overall operation of the electronic device 500 to complete all or part of the steps in the above-mentioned seismic scattered wave field forward modeling method. The memory 502 is used to store various types of data to support the operation of the electronic device 500, and these data may include, for example, instructions for any application or method used to operate on the electronic device 500, and application-related data, such as contact data, sent and received messages, pictures, audio, video, etc. The memory 502 can be implemented by any type of volatile or non-volatile storage device or a combination thereof, such as static random access memory (Static Random Access Memory, referred to as SRAM), electrically erasable programmable read-only memory (Electrically Erasable Programmable Read-Only Memory, referred to as EEPROM), erasable programmable read-only memory (Erasable Programmable Read-Only Memory, referred to as EPROM), programmable read-only memory (Programmable Read-Only Memory, referred to as PROM), read-only memory (Read-Only Memory, referred to as ROM), magnetic memory, flash memory, magnetic disk or optical disk. The multimedia component 503 may include a screen and an audio component. The screen may be, for example, a touch screen, and the audio component is used to output and / or input audio signals. For example, the audio component may include a microphone for receiving external audio signals. The received audio signal may be further stored in the memory 502 or sent through the communication component 505. The audio component also includes at least one speaker for outputting audio signals. The I / O interface 504 provides an interface between the processor 501 and other interface modules, and the above-mentioned other interface modules may be keyboards, mice, buttons, etc. These buttons may be virtual buttons or physical buttons. The communication component 505 is used for wired or wireless communication between the electronic device 500 and other devices. Wireless communication, such as Wi-Fi, Bluetooth, Near Field Communication (NFC), 2G, 3G, 4G, NB-IOT, eMTC, or other 5G, etc., or a combination of one or more of them, is not limited here. Therefore, the corresponding communication component 505 may include: Wi-Fi module, Bluetooth module, NFC module, etc.

[0156] In an exemplary embodiment, the electronic device 500 can be implemented by one or more application specific integrated circuits (ASIC), digital signal processors (DSP), digital signal processing devices (DSPD), programmable logic devices (PLD), field programmable gate arrays (FPGA), controllers, microcontrollers, microprocessors or other electronic components to execute the above-mentioned seismic scattered wave field forward modeling method.

[0157] In another exemplary embodiment, a computer-readable storage medium including program instructions is also provided, and when the program instructions are executed by a processor, the steps of the above-mentioned seismic scattered wave field forward modeling method are implemented. For example, the computer-readable storage medium may be the above-mentioned memory 702 including program instructions, and the above-mentioned program instructions may be executed by the processor 701 of the electronic device 700 to complete the above-mentioned seismic scattered wave field forward modeling method.

[0158] In another exemplary embodiment, a computer program product is also provided. The computer program product includes a computer program executable by a programmable device, and the computer program has a code portion for executing the above-mentioned seismic scattered wavefield forward modeling method when executed by the programmable device.

[0159] The preferred embodiments of the present disclosure are described in detail above in conjunction with the accompanying drawings; however, the present disclosure is not limited to the specific details in the above embodiments. Within the technical concept of the present disclosure, a variety of simple modifications can be made to the technical solution of the present disclosure, and these simple modifications all fall within the protection scope of the present disclosure.

[0160] It should also be noted that the various specific technical features described in the above specific embodiments can be combined in any suitable manner without contradiction. In order to avoid unnecessary repetition, the present disclosure will not further describe various possible combinations.

[0161] In addition, various embodiments of the present disclosure may be arbitrarily combined, and as long as they do not violate the concept of the present disclosure, they should also be regarded as the contents disclosed by the present disclosure.

Claims

1. A seismic scattering wave field forward modeling method, characterized in that: The method comprises: The velocity model is divided into multiple thin layers along the depth direction from the top to the bottom of the model; Along the first propagation direction, a first scattered wave field calculation is performed using a high-order generalized screen extension operator based on the Müller model to obtain downlink scattered wave data and one uplink scattered wave data corresponding to each thin layer in the first propagation direction; the first propagation direction is the direction from the top of the model to the bottom of the model; Along a second propagation direction, a second scattered wave field calculation is performed using the high-order generalized screen extension operator based on the Müller model to obtain first-order upgoing scattered wave data corresponding to each thin layer in the second propagation direction; the second propagation direction is the direction from the bottom of the model to the top of the model, and the second propagation direction is opposite to the first propagation direction; The forward seismic record of the detector position is determined according to the first upward scattered wave data corresponding to each thin layer in the second propagation direction.

2. The seismic scattered wave field forward modeling method according to claim 1, characterized in that: The number of the thin layers is N, where N is an integer greater than 1. The first scattering wave field calculation is performed along the first propagation direction using a high-order generalized screen extension operator based on the Müller model, including: Along the first propagation direction, the downlink scattered wave data of the previous thin layer is used as the boundary condition of the next adjacent thin layer, and the operation of calculating the scattered wave field by using the high-order generalized screen extension operator based on the Müller model is iteratively performed to obtain the downlink scattered wave data and the uplink scattered wave data corresponding to each thin layer in the N thin layers; The iteration is stopped when the downlink scattered wave data and one uplink scattered wave data corresponding to the Nth thin layer are obtained; wherein the Nth thin layer is a thin layer located at the bottom of the model.

3. The seismic scattered wave field forward modeling method according to claim 2, characterized in that: The method comprises: The source parameters are used as the boundary conditions of the first thin layer among N thin layers, and the scattered wave field is calculated using a high-order generalized screen continuation operator based on the Müller model to obtain the downlink scattered wave data and the first uplink scattered wave data corresponding to the second thin layer; wherein the first thin layer is the thin layer located at the top of the model among the N thin layers.

4. The seismic scattered wave field forward modeling method according to claim 3, characterized in that: The method uses the source parameters as the boundary condition of the first thin layer among the N thin layers and uses a high-order generalized screen extension operator based on the Müller model to calculate the scattered wave field, including: A two-dimensional coordinate system is established along the tunnel; wherein the origin of the two-dimensional coordinate system is the center of the tunnel surface on the velocity model boundary, the x-axis direction is the horizontal direction, and the z-axis direction is the depth direction; Acquire the earthquake source parameters; the earthquake source parameters include one or more of earthquake source location parameters, main frequency, sampling interval, sampling time length, number of detectors and detector location; The source parameters are used as boundary conditions of the first thin layer, and the scattered wave field is calculated using a high-order generalized screen extension operator based on the Müller model.

5. The seismic scattered wave field forward modeling method according to claim 1 or 2, characterized in that: The number of the thin layers is N, where N is an integer greater than 1, and the second scattering wave field calculation is performed along the second propagation direction using the high-order generalized screen extension operator based on the Müller model to obtain the first upward scattered wave data corresponding to each thin layer in the second propagation direction, including: Along the second propagation direction, the first upward scattered wave of the N-1th thin plate in the first propagation direction is used as the initial value condition in the second propagation direction, and the operation of performing scattered wave field calculation using the high-order generalized screen extension operator based on the Müller model is iteratively performed, wherein the first upward scattered wave field value retained at the same position in the first propagation direction is added to the calculation result of each scattered wave field calculation to obtain the first upward scattered wave data corresponding to each thin layer in the N thin layers; The iteration is stopped when the up-going scattered wave data corresponding to the Mth thin layer is obtained; wherein the Mth thin layer is the thin layer located at the top of the model; and the values ​​of M and N are the same.

6. The seismic scattered wave field forward modeling method according to claim 5, characterized in that: The operation of calculating the scattered wave field by using the high-order generalized screen extension operator based on the Müüller model is iteratively performed along the second propagation direction, using the first upward scattered wave of the N-1th thin plate in the first propagation direction as the initial value condition in the second propagation direction, and comprising: Taking the N-1th thin layer in the first propagation direction as the first thin layer in the second propagation direction; Acquire the primary upgoing scattered wave data corresponding to the N-1th thin layer in the first propagation direction; The first upward scattered wave data corresponding to the N-1th thin layer is used as the boundary condition of the first thin layer in the second propagation direction, and the scattered wave field is calculated using a high-order generalized screen extension operator based on the Müüller model. The first upward scattered wave data retained by the N-2th thin layer in the first propagation direction is added to the calculation result to obtain the first upward scattered wave data corresponding to the second thin layer in the second propagation direction.

7. The seismic scattered wave field forward modeling method according to any one of claims 1 to 6, characterized in that: The calculation formula of the Müüller model used in the scattered wave field calculation is: Where ω is the frequency, v0(ω) is the wave velocity at frequency ω, Q0 is the background quality factor, ω r is the reference frequency, γ is the frequency-dependent attenuation index, γ=1 / πQ0; The calculation formula for the downlink scattered wave or the uplink scattered wave in the scattered wave field calculation is: Among them, p(x, z j+1 ω) represents the down-going scattered wave corresponding to the next thin layer or the up-going scattered wave, p(x, z j ;ω) represents the downward scattered wave or the upward scattered wave corresponding to the previous thin layer, FT[] is the Fourier transform, FT -1 [] is the inverse Fourier transform, Δz is the extension step, A is the coefficient of the high-order generalized screen extension operator, which is used to correct the propagation characteristics of the wave field, ω is the circular frequency, is the vertical wave number; k x is the transverse wave number; Where n is the order, s0(x, z j ; ω) is the background velocity field, s(x, z j ; ω) is the actual velocity field, s(x,z j ;ω)=1 / v 2 (x,z j ;ω);v0(x,z j ; ω) represents the background velocity field, which is used to represent the velocity distribution in the medium without disturbance; v(x, z j ; ω) represents the actual velocity field, which contains disturbances or changes outside the background velocity field; The calculation formula for an upward scattered wave in the scattered wave field calculation is: b(x,z j+1 ;ω)=R*p(x,z j+1 ;ω) Among them, b(x, z j+1 ; ω) represents the once upward scattered wave, and R is the reflection coefficient.

8. The seismic scattered wave field forward modeling method according to any one of claims 1 to 6, characterized in that: The method of determining the forward seismic record according to the first upward scattered wave data corresponding to each thin layer in the first propagation direction and the first upward scattered wave data corresponding to each thin layer in the second propagation direction comprises: Acquire uplink scattered wave data in the first propagation direction according to the downlink scattered wave data corresponding to each thin layer in the first propagation direction; Acquire total upward scattered wave data in the second propagation direction according to the first upward scattered wave data corresponding to each thin layer in the first propagation direction and the first upward scattered wave data corresponding to each thin layer in the second propagation direction; A forward seismic record is determined based on the total upgoing scattered wave data in the second propagation direction.

9. A seismic scattering wave field forward modeling device, characterized in that: The device comprises: A model partitioning module is used to divide the velocity model into multiple thin layers from the top to the bottom of the model along the depth direction; A first calculation module is used to perform a first scattered wave field calculation along a first propagation direction using a high-order generalized screen extension operator based on a Müller model to obtain downlink scattered wave data and one uplink scattered wave data corresponding to each thin layer in the first propagation direction; the first propagation direction is a direction from the top of the model to the bottom of the model; A second calculation module is used to perform a second scattered wave field calculation along a second propagation direction using the high-order generalized screen extension operator based on the Müller model to obtain upgoing scattered wave data corresponding to each thin layer in the second propagation direction; the second propagation direction is the direction from the bottom of the model to the top of the model, and the second propagation direction is opposite to the first propagation direction; A determination module is used to determine the forward seismic record of the detector position according to the first upgoing scattered wave data corresponding to each thin layer in the second propagation direction.

10. A non-transitory computer-readable storage medium having a computer program stored thereon, characterized in that: When the program is executed by a processor, the steps of the method described in any one of claims 1 to 8 are implemented.

11. An electronic device, characterized in that: include: a memory having a computer program stored thereon; A processor, configured to execute the computer program in the memory to implement the steps of the method according to any one of claims 1 to 8.

12. A computer program product, comprising a computer program, characterized in that When the computer program is executed by a processor, the steps of the method according to any one of claims 1 to 8 are implemented.

Images