CN114966874B - Forward-modeling-based frequency domain electromagnetic method time sequence data synthesis method - Google Patents

Abstract

The present invention is applicable to the field of electromagnetic sounding technology, and provides a method for synthesizing time series data of frequency domain electromagnetic method based on forward modeling, comprising the following steps: Step 1: Perform frequency domain electromagnetic forward modeling to obtain field strength values of several frequencies; Step 2: Segment the time series of each frequency, and calculate the field value that changes with time; Step 3: Synthesize the segmented time series of each frequency, and splice the segmented time series; Step 4: Superimpose the spliced time series of each frequency to form a complete time series data and output it. The time series synthesis technology described in this method can be used to synthesize time series corresponding to the forward modeling results of various theoretical models. The results of the magnetotelluric method example show that the synthesized time series meets the basic characteristics of the frequency domain electromagnetic method, and the signal processing results are completely consistent with the model response.

Description

A method for synthesizing time series data of frequency domain electromagnetic method based on forward modeling

Technical Field

The invention belongs to the technical field of electromagnetic sounding and relates to a method for synthesizing time series data of a frequency domain electromagnetic method based on forward modeling.

Background Art

Frequency Domain Electromagnetic Method (FDMEM) is a general term for a series of geophysical methods that use frequency domain electromagnetic fields from natural and artificial sources to study the internal electrical structure of the earth, including magnetotelluric (MT), audio magnetotelluric (AMT), and controlled source audio magnetotelluric (CSAMT).

Its basic principle is: based on the principle that electromagnetic waves of different frequencies have different skin depths in conductive media, the electromagnetic field response sequence of natural sources (artificial sources) from high frequency to low frequency is measured on the surface, and the electrical structure of the earth from shallow to deep is obtained through relevant data processing. It has the advantages of not being shielded by high-resistance layers, strong resolution of high-conductivity layers, and strong lateral resolution, but also has the disadvantages of rapid weakening of vertical resolution with increasing depth. Its main steps are: data acquisition, data processing, and inversion interpretation.

The time series of electromagnetic fields collected by frequency domain electromagnetic method needs to be obtained by Fourier transform before the frequency domain field value can be processed and interpreted in the next step. The electromagnetic field time series is susceptible to interference, which affects the accuracy of subsequent data processing and interpretation. Existing data processing methods urgently need theoretical time series data as a standard for verification. At present, various denoising methods in frequency domain electromagnetic method are booming, but the effects of various processing methods are different for different data, and the processing results of different methods for the same data are also different. It is necessary to find a reasonable method to establish standard data to verify the processing effects of various data processing methods. In different frequency domain electromagnetic methods, the signals of some methods are noise to other methods. For example, line source signals are signals in controlled source electromagnetic method, but noise in natural source magnetotelluric method. This technology can synthesize signals from natural sources and artificial sources at the same time, so it can be used for the research of signal-to-noise separation technology.

Summary of the invention

The purpose of the embodiment of the present invention is to provide a method for synthesizing time series data of frequency domain electromagnetic method based on forward modeling, aiming to simulate the time series signal of frequency domain electromagnetic method, simulate the influence characteristics of various noises, and provide technical support for studying the distribution characteristics and propagation laws of electromagnetic method in time domain, frequency domain and space. For different data, the effects of various processing methods are not the same, and the processing results of different methods for the same data are also different. The present invention can reasonably establish standard data to test the processing effects of various data processing methods.

The embodiment of the present invention is implemented as follows: a method for synthesizing time series data of frequency domain electromagnetic method based on forward modeling, comprising the following steps:

Step 1: Perform frequency domain electromagnetic forward calculation to obtain field strength values of several frequencies;

Step 2: Segment the time series of each frequency and calculate the field value that changes with time;

Step 3: Synthesize the segmented time series of each frequency and splice the segmented time series together;

Step 4: Superimpose the spliced time series of each frequency to form a complete time series data and output it.

In a further technical solution, the specific steps of step 1 include:

Step 1.1: Design the theoretical model to be studied;

Step 1.2: For simple models, directly use analytical methods to solve them. For complex models that cannot be solved analytically, follow the steps below;

Step 1.2.1: Mesh the model, design measuring points, and record the number of measuring points as M;

Step 1.2.2: According to the frequency band and sampling rate of the study, design a frequency table with N frequency points, denoted as ω[j], j=1…N;

Step 1.2.3: Perform frequency domain electromagnetic forward modeling calculations on each frequency point of the model to obtain the electromagnetic field value at each measuring point. Each frequency point needs to be forward modeled twice, once for the E-polarized field source and once for the H-polarized field source. The electric field and magnetic field values obtained by the forward modeling are A ^E [i, j] and A ^H [i, j], where the superscripts ^E and ^H represent the E polarization and H polarization in the forward modeling, i = 1…M, j = 1…N, A∈{Ex, Ey, Hx, Hy, Hz};

A further technical solution is that the specific forward modeling method is not limited in step 1.2, and an analytical method for a simple model, or a numerical simulation method such as a finite difference method and a finite element method, or calculations using existing open source software such as ModEM, or calculations using existing commercial software such as COMSOL, are used to obtain the electromagnetic field value at the measuring point.

A further technical solution, the specific steps of step 2 include: for each frequency f, the time series to be synthesized is divided into multiple segments, and the length of each segment is not fixed, which is used to simulate the change of the field source over time. The present invention does not limit the specific length of each segment, which can be set according to the specific research objectives. In order to ensure that the time series can fully describe the spectrum information and also ensure that the spectrum information can be restored during data processing, the length of each segment should be greater than the length of a cycle corresponding to the frequency (1/f). If the length is too small, it is regarded as noise. The segment length in the best example of the present invention is a random number that conforms to the Gaussian distribution, with a mean value of 8/f and a variance of 4/f. The parameters of the random length are also not fixed. In order to ensure the stability of the synthetic data response, it is recommended that the mean value is above 4/f. The segments of E polarization and H polarization can be different, and the final time series is divided into ^CE segments and ^CH . For each time series, generate random complex numbers ^RE [k,j], k=1… ^CE and ^RH [k,j], k=1… ^CH , then calculate ^AE [k,i,j]= ^AE [i,j]* ^RE [k,j] and ^ＡＨ [k,i,j]= ^ＡＨ [i,j]* ^ＲＨ for each measuring point, obtain the E polarization field value and H polarization field value of each segment, and calculate their amplitude |A| and phase This makes the intensity and polarization direction of the electromagnetic field in each time series different.

In a further technical solution, the specific steps of step 3 are: for each frequency, using the amplitude |A| and phase calculated in step 2 According to the formula Calculate each segmented time series, then connect and splice these segmented time series at the end to obtain ^aE (j, t) and ^aH (j, t). The two are superimposed to obtain a(j, t). The calculation formula is a(j, t) = ^aE (j, t) + ^aH (j, t). In order to reduce the signal step generated by the splicing of time series, each time series can be windowed before splicing. The present invention does not limit the specific window function. The best embodiment of the present invention uses a Hanning window. The calculation formula is _aw (j, t) = a(j, t) * w(j, t), where w(j, t) is the window function and _aw (j, t) is the segmented time series after windowing.

In a further technical solution, step 2 and step 3 are a method of synthesizing an electromagnetic field whose polarization direction changes with time by using segmented splicing. This technology does not limit how to achieve the change of polarization direction.

In a further technical solution, the specific steps of step 4 are: superimposing the spliced time series signals corresponding to each frequency point to obtain the total time series, and the calculation formula is:

The embodiment of the present invention provides a method for synthesizing frequency domain electromagnetic method time series data based on forward modeling, which uses the time series synthesis technology described in the method to synthesize the time series of various theoretical models. The results show that the synthesized time series meets the basic characteristics of the frequency domain electromagnetic method signal, and the signal processing result is completely consistent with the model response. Therefore, the time series synthesized by this method can be used as a data standard, which is of great significance to the research on frequency domain electromagnetic method data processing methods.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG1 is a magnetotelluric time series signal and spectrum of a method for synthesizing frequency domain electromagnetic method time series data based on forward modeling provided by an embodiment of the present invention;

FIG2 is a flow chart of a time series synthesis method for synthesizing time series data of a frequency domain electromagnetic method based on forward modeling provided by an embodiment of the present invention;

FIG3 is a two-dimensional theoretical model of a method for synthesizing time series data of a frequency domain electromagnetic method based on forward modeling provided by an embodiment of the present invention;

FIG4 is a diagram showing the mesh generation and measurement point design of a two-dimensional model in a method for synthesizing time series data of a frequency domain electromagnetic method based on forward modeling provided by an embodiment of the present invention;

FIG5 is a diagram showing sinusoidal signals at various frequency points and their superimposed signals in a method for synthesizing time series data of a frequency domain electromagnetic method based on forward modeling provided by an embodiment of the present invention;

FIG. 6 shows a time series and processing results synthesized by a method for synthesizing frequency domain electromagnetic method time series data based on forward modeling provided in an embodiment of the present invention.

DETAILED DESCRIPTION

In order to make the purpose, technical solution and advantages of the present invention more clearly understood, the present invention is further described in detail below in conjunction with the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention and are not intended to limit the present invention.

The specific implementation of the present invention is described in detail below in conjunction with a specific embodiment of magnetotelluric time series synthesis.

Prerequisites:

The synthetic time series data conforms to the basic characteristics of natural source magnetotelluric signals. The impedance response obtained after the synthetic data is processed by conventional data processing procedures is consistent with the theoretical response of the design model. The synthetic data should be similar to the existing measured data and can be processed by mature data processing software to obtain a consistent response. The theoretical model is constructed for forward modeling, and the time series is synthesized. After data processing, the synthetic time series is compared and verified. The response of the synthetic data is consistent with the theoretical response.

Working principle:

The method for synthesizing time series data using the forward response of a theoretical model described in the present invention designs a theoretical electrical structure model and uses numerical simulation methods (finite difference method, finite element method, etc.) to calculate the frequency domain electromagnetic field response at the measuring point.

Frequency domain electromagnetic method studies the stable relationship between any single frequency electromagnetic field signals. Single frequency signals can be expressed as sinusoidal signals with a certain amplitude and phase. The electric field and magnetic field signals of a certain frequency can be written as:

In formula (1), r is the position vector, _Ae and _Ah are the amplitudes, ω is the angular frequency, and is the phase. Using the Euler identity e ^iφ = cos(φ) + isin(φ), where is an imaginary unit, formula (1) can be written as:

In the formula Denotes taking the real part. Now define the complex variable:

Formula (2) can be written as:

The Maxwell equations for the time domain electromagnetic method are:

In formula (5), μ is the magnetic permeability, ε is the dielectric constant, σ is the conductivity, and J(r, t) is the current source. Its form is the same as the electromagnetic field in formula (3) and formula (1). For the frequency domain electromagnetic method of natural sources, this term is 0. Substituting formula (4) into formula (5), the time domain Maxwell equation of a single frequency point can be obtained as follows:

The real part of the solution of formula (6) is formula (4). According to the operation rules of complex numbers:

Dividing both sides of the equations in equation (6) by e ^tωt , we can finally get:

Formula (8) is the Maxwell equations in the frequency domain that are independent of time, and is also the set of equations that need to be solved for the forward modeling of the frequency domain electromagnetic method. The above derivation process shows that if E(r) and H(r) are calculated according to formula (8), they can be substituted into formula (3) and formula (1) to obtain the time-varying field strength value of the frequency ω/2π. The process of calculating H(r) and E(r) with a frequency of ω/2π by formula (8) is a frequency domain forward modeling process, which is already very mature. Forward modeling can obtain electromagnetic fields of a finite number (N) of frequencies. According to formula (1), the field strength signals of each frequency that vary with time can be obtained. By superimposing them together, the frequency domain electromagnetic field signals of multiple frequency points can be obtained:

As shown in FIG1 , the magnetotelluric measured time series signal (A) and the spectrum (B) of the signal, the spectrum (B) can be calculated by equation (9) to obtain the time series (A).

Referring to FIG. 2 , the time series data synthesis process of the present invention includes the following steps:

Step 1: Perform magnetotelluric forward modeling to obtain field strength values of several frequencies. Design a suitable model to be studied (Figure 3). The model of the present invention can be a one-dimensional, two-dimensional or three-dimensional model. The model is meshed and the measuring points are designed (Figure 4). The number of measuring points is recorded as M. According to the frequency band and sampling rate of the study, a frequency table is designed, with a total of N frequency points, recorded as ω[j], j=1…N. Forward modeling is performed on each frequency point of the model to obtain the electromagnetic field value at each measuring point. Each frequency point needs to be forward modeled twice, once for an E-polarized field source and once for an H-polarized field source. The present invention does not limit the specific forward modeling method, and adopts an analytical method, or a numerical simulation method such as a finite difference method and a finite element method, or uses existing open source software such as ModEM for calculation, or uses existing commercial software such as COMSOL for calculation to obtain the electromagnetic field value at the measuring point. The electric field and magnetic field values obtained by forward modeling are ^AE [i,j] and ^ＡH [i,j], where the superscripts ^E and ^H represent the two polarization modes, E polarization and H polarization, in the forward modeling, i=1…M, j=1…N, A∈{Ex, Ey, Hx, Hy, Hz}.

Step 2: Segment the time series for each frequency and calculate the field value that changes with time. For each frequency f, the time series to be synthesized is divided into multiple segments, and the length of each segment is not fixed, which is used to simulate the change of the field source over time. The present invention does not limit the specific length of each segment, which can be set according to the specific research objectives. In order to ensure that the time series can fully describe the spectrum information and also ensure that the spectrum information can be restored during data processing, the length of each segment should be greater than the length of one cycle of the corresponding frequency (1/f). If the length is too small, it is regarded as noise. The segment length in the best example of the present invention is a random number that conforms to the Gaussian distribution, with a mean of 8/f and a variance of 4/f. The parameters of the random length are also not fixed. In order to ensure the stability of the synthetic data response, it is recommended that the mean be above 4/f. The segments of E polarization and H polarization can be different, and the final time series is divided into ^CE segments and ^CH . For each time series, generate random complex numbers ^RE [k,j], k=1… ^CE and ^RH [k,j], k=1… ^CH , then calculate ^AE [k,i,j]= ^AE [i,j]* ^RE [k,j] and ^ＡＨ [k,i,j]= ^ＡＨ [i,j]* ^ＲＨ for each measuring point, obtain the E polarization field value and H polarization field value of each segment, and calculate their amplitude |A| and phase This makes the intensity and polarization direction of the electromagnetic field in each time series different.

Step 3: Synthesize the segmented time series of each frequency and concatenate the segmented time series. For each frequency, use the amplitude |A| and phase calculated in step 2 According to the formula Calculate each segmented time series, then connect and splice these segmented time series at the end to obtain a ^E (j, t) and a ^H (j, t), and superimpose the two to obtain a (j, t), and the calculation formula is a (j, t) = a ^E (j, t) + a ^H (j, t). In order to reduce the signal step generated by the splicing of time series, each time series can be windowed before splicing. The present invention does not limit the specific window function. The best embodiment of the present invention uses a Hanning window. The calculation formula is _aw (j, t) = a (j, t) * w (j, t), where w (j, t) is the window function, and _aw (j, t) is the segmented time series after windowing.

Step 4: Superimpose the spliced time series of each frequency to form a complete time series data and output it. The calculation formula is: FIG. 4 shows the superposition of sinusoidal signals at multiple frequencies.

If you want to synthesize a time series in the measured data format, you need to add the observation system response to each electric field and magnetic field before step 3. The observation system response refers to the fact that after the actual input signal enters the instrument system, the recorded signal will be affected by the electrode distance of the electric track, the self-response of the magnetic probe, the self-response of the instrument, and the signal amplifier. These influences are manifested in the frequency domain as the characteristics of amplitude scaling and phase advance or lag at different frequencies. Therefore, before synthesizing the time series in the measured data format, you need to define the observation system of the measuring point, obtain the response of the observation system, and add the system response before calculating the electric field and magnetic field signals at each frequency point. Let the input signal be S ₀ (ω), the instrument recorded signal be S(ω), and the system response be r(ω), then S(ω)＝S ₀ (ω)r(ω), where r(ω) is a complex number, its modulus is the amplitude magnification after the signal passes through the system, and its phase is the phase lag after the signal passes through the system.

Figure 6 (A) shows the time series synthesized at a certain measuring point. The data is saved in the time series format of Phoenix formula MTU-5A, and the results are compared using the commercial software SSMT2000. After the test of synthetic time series data processing, it is shown that the theoretical time series synthesized by this method can be processed to obtain accurate responses, which can be used as a research standard for magnetotelluric time series processing technology.

The model is designed as shown in Figure 3. The model is a rectangular anomaly in a uniform half space. The background resistivity of the uniform half space is 100Ω·m, the resistivity of the anomaly is 10Ω·m, the width of the anomaly is 10km, the thickness is 5km, and the top is buried at a depth of 5km. Figure 4 is the core area of the model grid division (-20km<X<20km, -40km<Y<0km). After the core area grid is expanded outward, the left and right boundaries of the actual model are -200km and 200km respectively, the lower boundary is -100km, and the minimum layer thickness of the model is 200m. The designed frequency table (a total of 50 frequency points, the frequency range is from 320Hz to 0.00001Hz) is used for forward calculation. The finite element method is used to calculate the electromagnetic field value of each measuring point. The time series of the two polarization directions of each frequency point is segmented, and the segment length is a random number of Gaussian distribution, the mean is set to 8/f, and the variance is set to 4/f. The field strength value of each segment is calculated. Calculate the sinusoidal signals corresponding to each electric field and magnetic field in each polarization direction of each frequency point at measuring point A. Perform windowing on each time series, and the window function is a Hanning window. Splice the segmented time series, and then superimpose the time series in the two polarization directions to form a time series for each frequency point. Superimpose the time domain signals calculated for each frequency point to form a time series containing all frequency point information. Figure 5 shows the sinusoidal curves (Index=1~15) and superimposed curves (Index=16) corresponding to the first 15 frequency points of Ex, Ey, Hx and Hy at measuring point O. Figure 6 shows the time series (A) corresponding to the Phoenix MTU-5A instrument synthesized at this point, the result of processing the synthesized data using the mature commercial software SSMT2000 (B), and the comparison of the processing results with the theoretical response (C: apparent resistivity, D: phase). After comparison, it was found that the results were basically consistent, indicating that the synthesized data of the present invention is correct.

The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included in the protection scope of the present invention.

Claims (5)

1. A method for synthesizing time series data of frequency domain electromagnetic method based on theoretical model forward calculation, characterized in that it comprises the following steps: Step 1: Perform frequency domain electromagnetic forward calculation to obtain field strength values of several frequencies; Step 2: Segment the time series of each frequency and calculate the field value that changes with time; Step 3: Synthesize the segmented time series of each frequency and splice the segmented time series together; Step 4: Superimpose the spliced time series of each frequency to form a complete time series data and output it. 2. The method for synthesizing frequency domain electromagnetic method time series data based on theoretical model forward calculation according to claim 1 is characterized in that the specific steps of step 1 include: Step 1.1: Design the theoretical model to be studied; Step 1.2: For simple models, directly use analytical methods to solve them. For complex models that cannot be solved analytically, follow the steps below; Step 1.2.1: Mesh the model, design measuring points, and record the number of measuring points as M; Step 1.2.2: According to the frequency band and sampling rate of the study, design a frequency table with N frequency points, denoted as ω[j], j=1…N; Step 1.2.3: Perform frequency domain electromagnetic forward modeling calculations at each frequency point of the model to obtain the electromagnetic field value at each measuring point. Each frequency point requires two forward models, one for the E-polarized field source and one for the H-polarized field source. The electric and magnetic field values obtained by the forward modeling are ^AE [i,j] and ^ＡＨ [i,j], where the superscripts ^E and ^H represent E polarization and H polarization in the forward modeling, i＝1…M, j＝1…N, A∈{Ex，Ey，Hx，Hy，Hz}. 3. The method for synthesizing frequency domain electromagnetic method time series data based on theoretical model forward calculation according to claim 2 is characterized in that the specific steps of step 2 include: for each frequency f, the time series to be synthesized is divided into multiple segments, the length of each segment is not fixed, and is used to simulate the change of the field source over time. In order to ensure that the time series can fully describe the spectrum information and also ensure that the spectrum information can be restored during data processing, the length of each segment should be greater than the length of a cycle corresponding to the frequency (1/f). If the length is too small, it is regarded as noise. The parameters of the random length are also not fixed. The E polarization and H polarization segments can be different. Finally, the time series is divided into ^CE segments and ^CH segments. For each time series, random complex numbers ^RE [k,j], k=1... ^CE and ^RH [k,j], k=1... ^CH are generated, and then ^AE [k,i,j]= ^AE [i,j]* ^RE [k,j] and ^ＡＨ [k,i,j]= ^ＡＨ [i,j]* ^ＲＨ are calculated for each measuring point. [k, j], the E polarization field value and H polarization field value of each segment are obtained, so that the intensity and polarization direction of the electromagnetic field in each time series are different. 4. The method for synthesizing frequency domain electromagnetic time series data based on theoretical model forward calculation according to claim 3 is characterized in that the specific steps of step 3 are: for each frequency, using the E polarization field value and the H polarization field value obtained in step 2, wherein the E polarization field value and the H polarization field value are both complex numbers, respectively calculate the respective amplitude |A| and phase According to the formula Calculate each segmented time series, and then connect the ends of these segmented time series to obtain a ^E (j, t) and a ^H (j, t). The two are superimposed to obtain a(j, t). The calculation formula is a(j, t) = a ^E (j, t) + a ^H (j, t). In order to reduce the signal step generated by time series splicing, each time series is windowed before splicing. 5. The method for synthesizing frequency domain electromagnetic method time series data based on theoretical model forward calculation according to claim 4 is characterized in that the specific steps of step 4 are: superimposing the spliced time series signals corresponding to all frequency points to obtain the total time series, and the calculation formula is: