CN114966874A - A Synthesis Method of Frequency Domain Electromagnetic Method Time Series Data Based on Forward Modeling - Google Patents

Abstract

The invention is suitable for the technical field of electromagnetic sounding, and provides a forward-modeling-based method for synthesizing time series data of a frequency domain electromagnetic method, which comprises the following steps: step 1: performing electromagnetic forward calculation in a frequency domain to obtain field intensity values of a plurality of frequencies; step 2: segmenting the time sequence of each frequency, and calculating to obtain a field value which changes along with time; and step 3: synthesizing the segmented time sequences of each frequency, and splicing the segmented time sequences; and 4, step 4: and overlapping the spliced time series of each frequency to form complete time series data and outputting the complete time series data. By adopting the time sequence synthesis technology, time sequences corresponding to forward modeling results of various theoretical models can be synthesized. The result of the geoelectromagnetic method calculation shows that the time sequence obtained by synthesis meets the basic characteristics of the frequency domain electromagnetic method, and the processing result of the signal is completely consistent with the model response.

Description

A Synthesis Method of Frequency Domain Electromagnetic Method Time Series Data 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 based on a forward modeling frequency domain electromagnetic method.

Background technique

Frequency Domain Electromagnetic Method (Frequency Domain Electromagnetic Method) is a general term for a series of geophysical methods that use the frequency domain electromagnetic fields of natural and artificial sources to study the electrical structure of the earth's interior, including magnetotelluric (MT), audio geodetic methods. Electromagnetic bathymetry (Audio Magnetotelluric, AMT), Controlled Source Audio Magnetotelluric bathymetry (Controlled Source Audio Magnetotelluric, CSAMT) and so on.

The basic principle is: According to 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 earth is obtained through relevant data processing. Electrical structure from shallow to deep. It has the advantages of not being shielded by high-resistance layers, strong ability to distinguish high-conductivity layers, and strong lateral resolution. The main steps are: data acquisition, data processing, inversion interpretation.

The time series of the electromagnetic field collected by the frequency domain electromagnetic method needs to be obtained by Fourier transform to obtain the frequency domain field value before the next step of processing and interpretation. The electromagnetic field time series is easily affected by interference, which affects the accuracy of subsequent data processing and interpretation. Existing data processing methods urgently need theoretical time series data to be tested as a standard. At present, various denoising methods in the frequency domain electromagnetic method are booming, but 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. It is necessary to find a reasonable method to establish standard data to test The processing effect of various data processing methods. In different frequency domain electromagnetic methods, the signal of some methods is noise to other methods. For example, a line source signal is a signal in the controlled source electromagnetic method, but is noise in the natural source magnetotelluric method. This technology can simultaneously synthesize natural and artificial signals, so it can be used for the research of signal-to-noise separation technology.

SUMMARY OF THE INVENTION

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

The embodiments of the present invention are implemented in this way, a method for synthesizing time series data based on a forward modeling frequency domain electromagnetic method, comprising the following steps:

Step 1: Perform forward calculation of frequency domain electromagnetic method 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 splicing the segmented time series together;

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

Further technical solutions, the specific steps of the step 1 include:

Step 1.1: Design the theoretical model to be studied;

Step 1.2: For simple models, use the analytical method to solve them directly. 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 researched frequency band and sampling rate, design a frequency table with a total of N frequency points, denoted as ω[j], j=1...N;

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

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

使得每段时间序列中的电磁场的强度和极化方向都不同。In a further technical solution, the specific steps of step 2 include: for each frequency f, dividing the time series to be synthesized into multiple segments, and the length of each segment is not fixed, so as to simulate the variation of the field source with time. The present invention does not limit the specific length of each segment, which can be set according to specific research objectives. In order to ensure that the time series can completely describe the spectrum information, and also ensure that the spectrum information can be recovered during data processing, the length of each segment should be greater than the length of a cycle corresponding to the frequency (1/f), and if the length is too small, it is regarded as noise. In the best example of the present invention, the segment length is a random number conforming to a Gaussian distribution, the mean value is set to 8/f, and the variance is set to 4/f. The parameters of random length are also not fixed. In order to ensure the stable response of synthetic data, it is recommended that the average value 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... ^{C E and RH} [k,j], k=1... ^CH , and then calculate A ^E [k for each measurement point ,i,j]=A ^E [i,j]*R ^E [k,j] and A ^H [k,i,j]=A ^H [i,j]*R ^H , obtain the E pole of each segment and H-polarization field values, and calculate their amplitude |A| and phase

The intensity and polarization direction of the electromagnetic field in each time series are different.

按照式

计算每个分段时间序列，然后把这些分段时间序列收尾相连拼接起来，得到a^E(j，t)和a^H(j，t) 二者叠加得到a(j，t)，计算公式为a(j，t)＝a^E(j，t)+a^H(j，t)。为了减小时间序列拼接产生的信号阶跃，可以在拼接前对每段时间序列进行加窗处理，本发明不限定具体的窗函数，本发明的最佳实例中的使用的是汉宁窗。计算公式为a_w(j，t)＝a(j，t)*w(j，t)，其中w(j，t)为窗函数，a_w(j，t)为加窗后的分段时间序列。In a further technical solution, the specific steps of step 3 are: for each frequency, use the amplitude |A| and phase calculated in step 2

According to the formula

Calculate each segmented time series, and then splicing these segmented time series at the end to obtain a ^E (j, t) and a ^H (j, t) The superposition of both 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 time series splicing, each time series can be windowed before splicing. The present invention does not limit the specific window function, and the Hanning window is used in the best example of the present invention. The calculation formula is a _w (j, t)=a(j, t)*w(j, t), where w(j, t) is the window function, and a _w (j, t) is the windowed segment sequentially.

In a further technical solution, the steps 2 and 3 are a method of synthesizing an electromagnetic field in which the polarization direction changes with time by using segmented splicing. The present technology does not limit how to realize the change of the polarization direction.

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

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

Description of drawings

Fig. 1 is a kind of magnetotelluric time series signal and frequency spectrum based on a method for synthesizing time series data of frequency domain electromagnetic method provided by an embodiment of the present invention;

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

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

Fig. 4 is the grid division and measuring point design of a two-dimensional model in a forward modeling-based frequency-domain electromagnetic method time series data synthesis method provided by an embodiment of the present invention;

FIG. 5 is a sinusoidal signal of each frequency point and a superimposed signal thereof in a method for synthesizing time series data based on a forward modeling frequency domain electromagnetic method provided by an embodiment of the present invention;

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

Detailed ways

In order to make the objectives, technical solutions and advantages of the present invention clearer, the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention, but not to limit the present invention.

The specific implementation of the present invention will be described in detail below with reference to specific embodiments of magnetotelluric method time series synthesis.

Preconditions:

The synthetic time series data conforms to the basic characteristics of natural magnetotelluric signals. The impedance response obtained after the synthetic data is processed by the conventional data processing process is consistent with the theoretical response of the design model. The synthetic data should be similar to the existing measured data and can be used by mature data. Handle software processing and get a consistent response. Build a theoretical model for forward modeling, synthesize time series, and compare and verify the synthetic time series after data processing. The response of the synthetic data is consistent with the theoretical response.

working principle:

In the method for synthesizing time series data using the forward response of the theoretical model, the theoretical electrical structure model is designed, and the numerical simulation method (finite difference method, finite element method, etc.) is used to calculate the frequency domain at the measuring point. Electromagnetic field response.

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

和

为相位。利用欧拉恒等式e^iφ＝cos(φ)+isin(φ)，其中

为虚数单位，式(1)可以被写作：In formula (1), r is the position vector, A _e and A _h are the amplitudes, ω is the angular frequency,

and

for the phase. Using Euler's identity e ^iφ =cos(φ)+isin(φ), where

is an imaginary unit, equation (1) can be written as:

表示取实部。现定义复变量：in the formula

Represents the real part. Now define a complex variable:

Equation (2) can be written as:

The Maxwell equations of the time-domain electromagnetic method are:

In formula (5), μ is the magnetic permeability, ε is the permittivity, σ is the electrical conductivity, and J(r, t) is the current source, and its form is the same as the electromagnetic field in formula (3) and formula (1). Frequency domain electromagnetic method of natural sources, this item is 0. Putting Equation (4) into Equation (5), the Maxwell equation in the time domain of a single frequency point can be obtained as:

The solution result of formula (6) takes the real part, which is formula (4). According to the arithmetic of complex numbers:

In formula (6), divide both sides of each formula by e ^tωt , and finally we can get:

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

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

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

Step 1: Carry out the magnetotelluric forward calculation to obtain the field strength values of several frequencies. Design a suitable model to be studied (Fig. 3), the model of the present invention can be a one-dimensional, two-dimensional or three-dimensional model, mesh the model, design measuring points (Fig. 4), and the number of measuring points is M. According to the frequency band and sampling rate studied, a frequency table is designed, with a total of N frequency points, denoted as ω[j], j=1...N. Perform forward calculation for each frequency point of the model to obtain the electromagnetic field value at each measurement 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 present invention does not limit the specific forward modeling method, and adopts analytical method, or numerical simulation method such as finite difference method and finite element method, or uses existing open source software such as ModEM to calculate, or uses existing commercial software such as COMSOL to calculate, to obtain Electromagnetic field value at the measuring point. 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 two polarization modes of 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 of 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 variation of the field source with time. The present invention does not limit the specific length of each segment, which can be set according to specific research objectives. In order to ensure that the time series can completely describe the spectrum information, and also ensure that the spectrum information can be recovered during data processing, the length of each segment should be greater than the length of one cycle of the corresponding frequency (1/f), and if the length is too small, it is regarded as noise. In the best example of the present invention, the segment length is a random number conforming to a Gaussian distribution, the mean value is set to 8/f, and the variance is set to 4/f. The parameters of random length are also not fixed. In order to ensure the stable response of synthetic data, it is recommended that the average value 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 , and then calculate ^{A E [} k for each measurement point ,i,j]=A ^E [i,j]*R ^E [k,j] and A ^H [k,i,j]=A ^H [i,j]*R ^H , obtain the E pole of each segment and H-polarization field values, and calculate their amplitude |A| and phase

The intensity and polarization direction of the electromagnetic field in each time series are different.

按照式

计算每个分段时间序列，然后把这些分段时间序列收尾相连拼接起来，得到 a^E(j，t)和a^H(j，t)，二者叠加得到a(j，t)，计算公式为a(j，t)＝a^E(j，t)+a^H(j，t)。为了减小时间序列拼接产生的信号阶跃，可以在拼接前对每段时间序列进行加窗处理，本发明不限定具体的窗函数，本发明的最佳实例中的使用的是汉宁窗。计算公式为a_w(j，t)＝a(j，t)*w(j，t)，其中w(j，t)为窗函数,a_w(j，t)为加窗后的分段时间序列。Step 3: Synthesize the segmented time series of each frequency, and splicing the segmented time series together. For each frequency, use the amplitude |A| and phase calculated in step 2

According to the formula

Calculate each segmented time series, and then splicing these segmented time series together to obtain a ^E (j, t) and a ^H (j, t), which 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 time series splicing, each time series can be windowed before splicing. The present invention does not limit the specific window function, and the Hanning window is used in the best example of the present invention. The calculation formula is a _w (j, t)=a(j, t)*w(j, t), where w(j, t) is the window function, and a _w (j, t) is the windowed segment sequentially.

附图4展示的是多个频点的正弦信号叠加。Step 4: Superimpose the spliced time series of each frequency to form a complete time series data and output, the calculation formula is

FIG. 4 shows the superposition of sinusoidal signals of multiple frequency points.

If you want to synthesize the time series of the measured data format, you need to add the observation system response for each electric field and magnetic field before step 3. The observed system response means that after the actual input signal enters the instrument system, the recorded signal will be affected by the distance between the electrodes of the channel, the self-response of the magnetic probe, the self-response of the instrument, and the signal amplifier. and phase advance or lag characteristics. Therefore, before synthesizing the time series of the measured data format, the observation system of the measurement point should be defined, the response of the observation system should be obtained, and the system response should be added before calculating the electric field and magnetic field signals of each frequency point. Let the input signal be S ₀ (ω), the instrument recording signal is S(ω), and the system response is r(ω), then S(ω)=S ₀ (ω)r(ω), where r(ω) is Complex number, its modulus is the magnification of the amplitude of the signal after passing through the system, and the phase is the lag of the phase after the signal passes through the system.

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

The design of this model is shown in Figure 3. The model is that the uniform half space contains a rectangular abnormal body. The background resistivity of the uniform half space is 100Ω·m, the resistivity of the abnormal body is 10Ω·m, and the abnormal body is 10km wide and 5km thick. , the top is buried at a depth of 5km. Figure 4 shows the core area (-20km<X<20km, -40km<Y<0km) of the mesh division of the model. After the grid of the core area expands outward, the left and right boundaries of the actual model are respectively -200km and 200km. The boundary is -100km, and the minimum layer thickness of the model is 200m. The designed frequency table (50 frequency points in total, frequency range from 320Hz to 0.00001Hz) is used for forward calculation. The electromagnetic field value of each measuring point is calculated by the finite element method. For the time series segments of the two polarization directions of each frequency point, the segment length is a random number of Gaussian distribution, the mean value is set to 8/f, and the variance is set to 4/f, and the field strength value of each segment is calculated. Calculate the sine signal corresponding to each electric field and magnetic field in the two polarization directions of each frequency point of measuring point A. Windowing is performed on each time series, and the window function is Hanning window. Splicing the segmented time series, and then superimposing the time series of the two polarization directions to form the time series of each frequency point. The time domain signals calculated by each frequency point are superimposed together 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 the measuring point O. Accompanying drawing 6 is the time series (A) corresponding to the MTU-5A instrument of Phoenix Company synthesized at this point, the synthetic data adopts the result (B) of mature commercial software SSMT2000 processing, the contrast of processing result and theoretical response (C: apparent resistivity, D: Phase). Through comparison, it is found that the results are basically the same, indicating that the synthetic data of the present invention is correct.

The above descriptions are only preferred embodiments of the present invention and are not intended to limit the present invention. Any modifications, equivalent replacements and improvements made within the spirit and principles of the present invention shall be included in the protection of the present invention. within the range.

Claims (7)

1. a frequency domain electromagnetic method time series data synthesis method based on theoretical model forward calculation, is characterized in that, comprises the following steps: Step 1: Perform forward calculation of frequency domain electromagnetic method 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 splicing the segmented time series together; Step 4: Superimpose the spliced time series of each frequency to form complete time series data and output. 2. the frequency domain electromagnetic method time series data synthesis method based on theoretical model forward calculation according to claim 1, is characterized in that, the concrete steps of described step 1 comprise: Step 1.1: Design the theoretical model to be studied; Step 1.2: For simple models, use the analytical method to solve them directly. 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 researched frequency band and sampling rate, design a frequency table with a total of N frequency points, denoted as ω[j], j=1...N; Step 1.2.3: Perform forward calculation of frequency domain electromagnetic method for each frequency point of the model to obtain the electromagnetic field value at each measurement point. Each frequency point needs to be forward modeled twice, once for the E-polarized field source , the field source of primary H polarization, the electric field and magnetic field values obtained by forward modeling are A ^E [i,j] and A ^H [i,j], where superscript ^E and ^H represent E polarization and H pole in forward modeling , i=1...M, j=1...N, A∈{Ex, Ey, Hx, Hy, Hz}. 3. The frequency domain electromagnetic method time series data synthesis method based on theoretical model forward calculation according to claim 1 The frequency domain electromagnetic method time series data synthesis method based on theoretical model forward calculation is characterized in that, the frequency Domain electromagnetic methods include but are not limited to all frequency-domain electromagnetic detection methods, such as magnetotelluric method, audio-frequency magnetotelluric method, controllable-source audio-frequency magnetotelluric method, and frequency-domain airborne electromagnetic method. 4. the frequency domain electromagnetic method time series data synthesis method based on theoretical model forward calculation according to claim 2, is characterized in that, the forward method in described step 1.2 includes but is not limited to analytical method, numerical simulation method, Open source software computing and commercial software computing. 5. The frequency domain electromagnetic method time series data synthesis method based on theoretical model forward calculation according to claim 1, is characterized in that, the concrete steps of described step 2 comprise: for each frequency f, combine the time to be synthesized The sequence is divided into multiple segments, and the length of each segment is not fixed. It is used to simulate the change of the field source with time. In order to ensure that the time series can completely describe the spectrum information, and also ensure that the spectrum information can be recovered during data processing, the length of each segment should be greater than the frequency specified. The length of a corresponding period (1/f), if the length is too small, it is regarded as noise, and the parameters of random length are not fixed. The segments of E polarization and H polarization can be different, and the final time series is divided into C ^E segment and ^CH segment. For each time series, generate random complex numbers RE[k,j], k ⁼ 1... ^CE and ^RH [k,j], k=1... ^CH , j=1...N, and then for each measurement Point calculation A ^E [k,i,j]=A ^E [i,j]*R ^E [k,j] and A ^H [k,i,j]=A ^H [i,j]*R ^H , get The E-polarization field value and H-polarization field value of each segment make the intensity and polarization direction of the electromagnetic field in each time series different. 6. the frequency domain electromagnetic method time series data synthesis method based on theoretical model forward calculation according to claim 1, is characterized in that, the concrete step of described step 3 is: for each frequency, utilize the E obtained in step 2 Polarization field value and H polarization field value (both complex numbers), respectively calculate the respective amplitude |A| and phase

According to the formula

Calculate each segmented time series, and then splicing these segmented time series together to obtain a ^E (j, t) and a ^H (j, t), which 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 time series splicing, each time series can be windowed before splicing. 7. the frequency domain electromagnetic method time series data synthesis method based on theoretical model forward calculation according to claim 1, is characterized in that, the concrete step of described step 4 is: the splicing time series signal corresponding to all frequency points is superimposed To get the total time series, the formula is

Images