Method for frequency domain conversion of time domain electromagnetic data of magnetic source at specified measuring point in specified work area
By determining the typical resistivity value and calculating the induced electromotive force vector, using the random resistivity vector and the fitness threshold optimization, the electromagnetic response of the magnetic source frequency domain is finally calculated, which solves the matrix singularity and initial value dependence problems in the time domain data conversion of magnetic source, and improves the accuracy and reliability of the conversion.
Patent Information
- Application Number
- CN202310589376.X
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2023-05-23
- Publication Date
- 2025-08-29
- Estimated Expiration
- 2043-05-23
AI Technical Summary
When the prior art performs frequency domain conversion of the magnetic source time domain electromagnetic data at designated measurement points in the specified work area, there are strong singularity of the solution matrix, strong dependence on the initial value, and easy to cause local small problems, resulting in insufficient conversion reliability and accuracy.
By determining the typical resistivity value of the specified work area, collecting time-domain electromagnetic data and calculating the induced electromotive force vector, determining the optimal resistivity vector using the random resistivity vector and the fitness threshold, and finally calculating the electromagnetic response of the magnetic source frequency domain, and optimizing using the smooth constraint least squares method and gradient algorithm.
The frequency domain conversion accuracy and reliability of magnetic source time domain electromagnetic data is improved, and the matrix singularity and initial value dependence problems exist in the prior art are solved, thereby avoiding the occurrence of local extremely small points.
Smart Images

Figure CN116626771B_ABST
Abstract
Description
Technical Field
[0001] The embodiments of the present application relate to the field of geophysical electromagnetic data processing, and specifically to a method for performing frequency domain conversion on time-domain electromagnetic data of a magnetic source at a specified measuring point in a specified work area. Background Art
[0002] Electromagnetic methods are geophysical exploration techniques that generate electromagnetic fields to detect specific points within a designated work area. These include magnetic electromagnetic methods, which generate electromagnetic fields using a transmitting coil. These methods offer advantages such as ease of use and minimal labor intensity, making them commonly used for ground, aerial, and marine surveys.
[0003] Magnetic source time-domain electromagnetic data refers to magnetic source electromagnetic data obtained with time as the independent variable. Magnetic source time-domain electromagnetic data can be converted to magnetic source frequency-domain electromagnetic data, that is, magnetic source electromagnetic data obtained with frequency as the independent variable. Currently, when measuring at a designated point in a designated work area, magnetic source time-domain electromagnetic data is typically obtained. However, the processing of magnetic source time-domain electromagnetic data is one-dimensional, which is not as mature as the three-dimensional forward and inversion processing used for processing magnetic source frequency-domain electromagnetic data. Therefore, the obtained magnetic source time-domain electromagnetic data can be converted to the frequency domain, and the resulting magnetic source frequency-domain electromagnetic data can be processed using frequency-domain three-dimensional forward and inversion. Summary of the Invention
[0004] In view of the above problems, the present application is proposed to provide a method for performing frequency domain conversion on time-domain electromagnetic data of a magnetic source at a specified measuring point in a specified work area.
[0005] An embodiment of the present application provides a method for frequency-domain conversion of time-domain electromagnetic data of a magnetic source at a specified measuring point in a specified work area, comprising the following steps: S1: determining a typical resistivity value of the specified work area based on the specified work area; S2: collecting time-domain electromagnetic data of the specified measuring point in the specified work area, and determining a time trace vector and an induced electromotive force vector; S3: determining a random resistivity vector of the specified work area based on the typical resistivity value; S4: determining a calculated induced electromotive force vector based on the random resistivity value vector and the time trace vector; S5: determining fitness based on the calculated induced electromotive force vector and the induced electromotive force vector; S6: determining an optimal resistivity vector based on the fitness and a fitness threshold; and S7: calculating a frequency-domain electromagnetic response of the magnetic source based on the optimal resistivity vector.
[0006] The method for performing frequency domain conversion on the time domain electromagnetic data of the magnetic source at the specified measuring point in the specified work area in the embodiment of the present application can improve the accuracy and reliability of the frequency domain conversion. BRIEF DESCRIPTION OF THE DRAWINGS
[0007] Figure 1This is a flow chart of a method for performing frequency domain conversion on time-domain electromagnetic data of a magnetic source at a specified measuring point in a specified work area according to an embodiment of the present application;
[0008] Figure 2 is a flow chart of updating the random resistivity vector p according to an embodiment of the present application;
[0009] Figure 3 FIG. 4 is a flow chart of updating the random resistivity vector p according to yet another embodiment of the present application.
[0010] It should also be noted that the drawings are only for the purpose of illustrating the preferred embodiments, not the application itself. The drawings do not illustrate every aspect of the described embodiments and do not limit the scope of the application. DETAILED DESCRIPTION
[0011] To make the purpose, technical solutions, and advantages of this application more clear, the technical solutions of this application will be clearly and completely described below in conjunction with the accompanying drawings of the embodiments of this application. Obviously, the described embodiment is only one embodiment of this application, not all embodiments. Based on the described embodiments of this application, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of this application.
[0012] It should be noted that, unless otherwise defined, the technical terms or scientific terms used in this application should have the usual meanings understood by people with ordinary skills in the field to which this application belongs. If the full text involves descriptions such as "first" and "second", the "first" and "second" descriptions are only used to distinguish similar objects, and cannot be understood as indicating or implying their relative importance, order of precedence, or implicitly indicating the number of technical features indicated. It should be understood that the data described by "first" and "second" can be interchangeable under appropriate circumstances. If "and / or" appears in the full text, its meaning includes three parallel solutions. Taking "A and / or B" as an example, it includes Solution A, Solution B, or solutions that meet both A and B.
[0013] The inventors of the present application have discovered that, when the time-domain electromagnetic data of the magnetic source measured at a designated measuring point in a designated work area are converted to the frequency domain, there are methods based on discrete Fourier transform and using the smooth constrained least squares method for solution, which directly converts the time-domain data into the frequency domain data. However, the above methods have the problem of strong singularity of the solution matrix; there are conventional optimization methods, which have the problem of strong dependence on initial values; and there are gradient-based algorithms, which are prone to local minima and reduce the reliability of the time-frequency data conversion.
[0014] To this end, the embodiment of the present application provides a more reliable and accurate method for converting the time-domain electromagnetic data of the magnetic source at a specified measuring point in a specified work area into the frequency domain, such as Figure 1 As shown, it shows a flowchart of a method for performing frequency domain conversion on time-domain electromagnetic data of a magnetic source at a specified measuring point in a specified work area according to an embodiment of the present application, which includes the following steps: S1: determining a typical resistivity value of the specified work area based on the specified work area; S2: collecting time-domain electromagnetic data of the specified measuring point in the specified work area, and determining a time trace vector and an induced electromotive force vector; S3: determining a random resistivity vector of the specified work area based on the typical resistivity value; S4: determining a calculated induced electromotive force vector based on the random resistivity value vector and the time trace vector; S5: determining fitness based on the calculated induced electromotive force vector and the induced electromotive force vector; S6: determining an optimal resistivity vector based on the fitness and the fitness threshold; S7: calculating the frequency domain electromagnetic response of the magnetic source based on the optimal resistivity vector.
[0015] In some embodiments, in step S1, a typical resistivity value for a designated work area can be obtained by querying local resistivity rock property data or geological data. In some embodiments, a DC resistivity meter or a high-density resistivity meter can also be used to measure the resistivity value in the designated work area, and the measurement result is the typical resistivity value for the designated work area.
[0016] In some embodiments, in step S2, the time channel vector is a vector of a series of moments, such as T = [0.00001, 0.0001, 0.001, 0.01] s, and the specific moment values in the time channel vector are the moment values set in the acquisition instrument used.
[0017] In some embodiments, in step S2, when determining the induced electromotive force vector, the coordinates corresponding to when the induced electromotive force vector was collected can also be obtained to determine the position when the induced electromotive force was collected.
[0018] In some embodiments, the length of the random resistivity vector may be 10, and the number of random resistivity vectors may be twice its length.
[0019] In some embodiments, in step S3, the random resistivity vector is determined using the following expression:
[0020] p i =p min +rand×(p max –p min ),
[0021] Among them, p i is the resistivity value of serial number i, p min is the minimum resistivity value, p maxis the maximum resistivity value, rand is a random number between 0 and 1, and p min The value is 0.01p0, p0 is the typical resistivity value, p max The value is 100p0, if p i Less than p min , then p i =p min , if p i Greater than p max , then p i =p max .
[0022] That is, when the typical resistivity value is 100 ohm-meter, the minimum resistivity value is 1 ohm-meter, and the maximum resistivity value is 10,000 ohm-meter. If the resistivity value of sequence number i is less than the minimum resistivity value, the minimum resistivity value is assigned to the resistivity of sequence number i; if the resistivity value of sequence number i is greater than the minimum resistivity value, the maximum resistivity value is assigned to the resistivity of sequence number i.
[0023] In some embodiments, in step S4, the induced electromotive force vector is calculated by the following expression:
[0024] in,
[0025] J1 is the first-order Bessel function of the first kind, t represents the time channel, d mod To calculate the induced electromotive force vector, I is the emission current, a is the radius of the magnetic source, ω is the circular frequency, m is the Hankel transform wave number, r TE is the resistivity reflection coefficient, and the random resistivity vector p is used to determine r TE .
[0026] In some embodiments, the circular frequency is specified by a cosine transform. The cosine transform process and Hankel transform wavenumber are referenced in the following literature: New digital linear filters for Hankel J0 and J1 transforms (D. Guptasarma, B. Singh, Geophysical prospecting, 1997, 45: 745-762). and Optimal digital filters for Sine and Cosine transforms (Zhao Yun-wei, Zhu Zi-qiang, Lu Guang-yin, Han Bo, Journal of Applied Geophysics 2018). The process of calculating the resistivity reflection coefficient is referenced in the following literature: Inversion of helicopter electromagnetic data to a magnetic conductive layered earth (Haoping Huang, Douglas C. Fraser, Geophysics, 2003, 68(4): 1211-1223).
[0027] In some embodiments, in step S5, the fitness is determined by the following expression:
[0028] f=1 / [1+(d obs –d mod ) T (d obs –d mod )+β(Dp) T (Dp)], where f represents fitness, d obs is the induced electromotive force vector, d mod To calculate the induced electromotive force vector, β is a predetermined value, D is a one-dimensional difference matrix of l × l, whose diagonal elements are 1 and the elements on the right side of the diagonal are -1, and l is the length of the resistivity vector;
[0029] n f constitute the fitness vector, and n is the number of resistivity vectors.
[0030] In some embodiments, β can be adjusted. In some embodiments, the value of β can be 1.
[0031] In some embodiments, the method further includes the step of updating the random resistivity vector p, determining an updated fitness value based on the updated random resistivity vector p, and comparing the updated fitness value with a fitness threshold value, wherein the random resistivity vector p before the update and the updated random resistivity vector p satisfy the following expression:
[0032] p new =p+α×(p–p k ), where p new represents the resistivity value in the random resistivity vector after update, α is a random number between -1 and 1, p represents the resistivity value in the random resistivity vector before update, p k represents the value of the resistivity before the kth update in p;
[0033] If the updated fitness is better than the current fitness, the random resistivity vector p is updated, and the updated random resistivity vector p is adopted in step S7.
[0034] Figure 2 FIG. 4 shows a flow chart of updating the random resistivity vector p according to an embodiment of the present application. Figure 2 As shown, the following steps are included: S61: updating the random resistivity vector p; S62: determining the updated fitness according to the updated random resistivity vector p; S63: comparing the updated fitness with the current fitness, if the updated fitness is better than the current fitness, executing step S641; if the updated fitness is not better than the current fitness, executing step S642; S641: updating the random resistivity vector p; S642: not updating the random resistivity vector p; S651: adopting the updated random resistivity vector p in step S7; S652: adopting the current random resistivity vector p in step S7.
[0035] In some embodiments, all current random resistivity vectors are updated.
[0036] In some embodiments, the random resistivity vector p is updated multiple times until the difference between the fitness after the update and the fitness before the update is less than a predetermined value, and the updating is stopped, thereby determining the optimal random resistivity vector p.
[0037] Figure 3 FIG. 4 shows a flow chart of updating the random resistivity vector p according to another embodiment of the present application, as shown in FIG. Figure 3As shown, the method includes the following steps: S61: updating the random resistivity vector p; S62: determining the updated fitness according to the updated random resistivity vector p; S63: comparing the updated fitness with the current fitness, if the updated fitness is better than the current fitness, executing step S641; if the updated fitness is not better than the current fitness, executing step S642; S641: updating the random resistivity vector p; S642: not updating the random resistivity vector p; S65: comparing the difference between the current fitness and the previous fitness with a predetermined value, if the difference between the current fitness and the previous fitness is greater than the predetermined value, executing step S61 again; if the difference between the current fitness and the previous fitness is less than the predetermined value, executing step S66; S66: adopting the current random resistivity vector p in step S7.
[0038] In some embodiments, the predetermined value is 0.0001.
[0039] In some embodiments, the following steps are also included: generating a random number between 0 and 1, and calculating the probability of updating the random resistivity vector p by the following expression: [0.9*f / max(f)]+0.1; if the random number is less than this probability, then updating is performed.
[0040] In some embodiments, if the resistivity value in the updated random resistivity vector is less than the minimum resistivity value, the minimum resistivity value is assigned to the corresponding resistivity in the updated random resistivity vector; if the resistivity value in the updated random resistivity vector is greater than the maximum resistivity value, the maximum resistivity value is assigned to the corresponding resistivity in the updated random resistivity vector.
[0041] In some embodiments, if the updated fitness is not better than the current fitness, the random resistivity vector is not updated, and the number of retentions is increased by 1. If the updated fitness is better than the current fitness, the random resistivity vector is updated, and the number of retentions is set to 0. If the number of retentions of the random resistivity vector exceeds the maximum number of retentions, the random resistivity vector is discarded, and step S3 is executed to regenerate a random resistivity vector. Steps S4-S6 are executed using the newly generated random resistivity vector until the optimal resistivity vector is determined.
[0042] In some embodiments, the initial retention number is 0, and the maximum retention number may be 5 times the length of the resistivity vector.
[0043] In some embodiments, if the value of any resistivity in the regenerated new random resistivity vector is less than the minimum resistivity value, the minimum resistivity value is assigned to the corresponding resistivity in the regenerated new random resistivity vector; if the value of any resistivity in the regenerated new random resistivity vector is greater than the maximum resistivity value, the maximum resistivity value is assigned to the corresponding resistivity in the regenerated new random resistivity vector.
[0044] In some embodiments, in step S7, the frequency domain electromagnetic response of the magnetic source satisfies the following expression:
[0045] where d result Represents the frequency domain electromagnetic response of the magnetic source.
[0046] In the above expression, r TE is the resistivity reflection coefficient determined using the optimal resistivity vector, m is the Hankel transform wave number, and a is the radius of the magnetic source.
[0047] In some embodiments, when the fitness does not reach the fitness threshold, the calculation parameters in the above embodiments may be adjusted and the random resistivity vector may be regenerated.
[0048] In some embodiments, when the fitness does not reach the fitness threshold and the resistivity value in the updated random resistivity vector is the minimum resistivity value or the maximum resistivity value, the minimum resistivity value and the maximum resistivity value can be adjusted, and step S3 is executed to regenerate the random resistivity vector. Steps S4-S6 are executed using the newly generated random resistivity vector until the optimal resistivity vector is determined. In some embodiments, if any resistivity value in the newly generated random resistivity vector is less than the minimum resistivity value, the minimum resistivity value is assigned to the corresponding resistivity in the newly generated random resistivity vector; if any resistivity value in the newly generated random resistivity vector is greater than the maximum resistivity value, the maximum resistivity value is assigned to the corresponding resistivity in the newly generated random resistivity vector.
[0049] In some embodiments, when adjusting the minimum resistivity value and the maximum resistivity value, the minimum resistivity value may be reduced by 10 times, and the maximum resistivity value may be increased by 10 times.
[0050] In some embodiments, when the fitness does not reach the fitness threshold and the resistivity value in the updated random resistivity vector is not the minimum resistivity value or the maximum resistivity value, the predetermined value may be adjusted, and step S3 may be executed to regenerate the random resistivity vector. Steps S4-S6 may be executed using the newly generated random resistivity vector until the optimal resistivity vector is determined. In some embodiments, if any resistivity value in the newly generated random resistivity vector is less than the minimum resistivity value, the minimum resistivity value is assigned to the corresponding resistivity in the newly generated random resistivity vector; if any resistivity value in the newly generated random resistivity vector is greater than the maximum resistivity value, the maximum resistivity value is assigned to the corresponding resistivity in the newly generated random resistivity vector.
[0051] In some embodiments, when adjusting the predetermined value, the predetermined value may be reduced by a factor of 10. Furthermore, since excessively low predetermined values may cause parameter oscillation in the random resistivity vector, the predetermined value may be adjusted multiple times, and the maximum predetermined value that enables the fitness to reach the fitness threshold may be selected.
[0052] In some embodiments, when the fitness has not reached the fitness threshold and has been reduced by a predetermined value multiple times, the length of the resistivity vector may be increased, and step S3 may be executed to regenerate the random resistivity vector. Steps S4-S6 may be executed using the newly generated random resistivity vector until the optimal resistivity vector is determined. In some embodiments, if any resistivity value in the newly generated random resistivity vector is less than a minimum resistivity value, the minimum resistivity value is assigned to the corresponding resistivity in the newly generated random resistivity vector; if any resistivity value in the newly generated random resistivity vector is greater than a maximum resistivity value, the maximum resistivity value is assigned to the corresponding resistivity in the newly generated random resistivity vector.
[0053] In some embodiments, since increasing the length of the resistivity vector will slow down the calculation, the length of the resistivity vector may be adjusted multiple times, and the minimum resistivity vector length that enables the fitness to reach the fitness threshold is selected.
[0054] In some embodiments, if the fitness still does not reach the fitness threshold after all the above adjustments are made, the fitness threshold is set too low, and steps S4-S6 are performed multiple times. After the predetermined number of executions is reached, the random resistivity vector with the best fitness is determined as the optimal resistivity vector, and step S7 is performed.
[0055] Those skilled in the art will appreciate that the embodiments described above are exemplary and that they may be improved upon. The structures described in the various embodiments may be freely combined without causing any conflict in structure or principle.
[0056] Although some embodiments according to the overall technical concept of the present disclosure have been shown and described, those skilled in the art will understand that changes may be made to these embodiments without departing from the principles and spirit of the overall technical concept of the present disclosure, and the scope of the present disclosure is defined by the claims and their equivalents.
Claims
1. A method for converting time-domain electromagnetic data of a magnetic source at a specified measuring point in a specified work area into a frequency-domain, wherein: The following steps are involved: S1: Based on the designated work area, determining a typical resistivity value of the designated work area; S2: collecting time domain electromagnetic data of a designated measuring point in the designated work area, and determining a time trace vector and an induced electromotive force vector; S3: determining a random resistivity vector of the designated work area according to the typical resistivity value; S4: determining and calculating an induced electromotive force vector according to the random resistivity value vector and the time track vector; S5: determining fitness based on the calculated induced electromotive force vector and the induced electromotive force vector; S6: determining an optimal resistivity vector according to the fitness and the fitness threshold; S7: Calculate the frequency domain electromagnetic response of the magnetic source based on the optimal resistivity vector.
2. The method according to claim 1, wherein In step S3, the random resistivity vector is determined using the following expression: p i =p min +rand×(p max -p min ), Among them, p i is the resistivity value of serial number i, p min is the minimum resistivity value, p max is the maximum resistivity value, rand is a random number between 0 and 1, And p min The value is 0.01p0, p0 is the typical resistivity value, p max The value is 100p0, if p i Less than p min , then p i =p min , if p i Greater than p max , then p i =p max .
3. The method according to claim 2, wherein: In step S4, the calculated induced electromotive force vector is determined by the following expression: in, J1 is the first-order Bessel function of the first kind, t represents the time channel, d mod To calculate the induced electromotive force vector, I is the emission current, a is the radius of the magnetic source, ω is the circular frequency, m is the Hankel transform wave number, r TE is the resistivity reflection coefficient, and the random resistivity vector p is used to determine r TE .
4. The method according to claim 3, wherein: In step S5, the fitness is determined by the following expression: f=1 / [1+(d obs -d mod ) T (d obs -d mod )+β(Dp) T (Dp)], where f represents fitness, d obs is the induced electromotive force vector, d mod To calculate the induced electromotive force vector, β is a predetermined value, D is a one-dimensional difference matrix of l × l, whose diagonal elements are 1 and the elements on the right side of the diagonal are -1, and l is the length of the resistivity vector; n f constitute the fitness vector, and n is the number of resistivity vectors.
5. The method according to claim 4, wherein The method includes the steps of updating the random resistivity vector p, determining an updated fitness according to the updated random resistivity vector p, and comparing the updated fitness with a fitness threshold, wherein the random resistivity vector p before the update and the updated random resistivity vector p satisfy the following expression: p new =p+α×(pp k ), where p new represents the resistivity value in the random resistivity vector after update, α is a random number between -1 and 1, p represents the resistivity value in the random resistivity vector before update, p k represents the value of the resistivity before the kth update in p; If the updated fitness is better than the current fitness, the resistivity vector p is updated, and the updated random resistivity vector p is adopted in step S7.
6. The method according to claim 5, wherein: The random resistivity vector p is updated multiple times until the difference between the fitness after the update and the fitness before the update is less than a predetermined value, and the updating is stopped, thereby determining the optimal resistivity vector p.
7. The method according to claim 1, wherein In step S7, the frequency domain electromagnetic response of the magnetic source satisfies the following expression: where d result represents the electromagnetic response of the magnetic source in the frequency domain, J1 is the first-order Bessel function of the first kind, r TE is the resistivity reflection coefficient determined using the optimal resistivity vector, m is the Hankel transform wave number, and a is the magnetic source radius.
8. The method according to claim 6, wherein: The predetermined value is 0.0001.
9. The method according to claim 4, wherein: The β can be adjusted.
10. The method according to claim 6, wherein: The following steps are also included: Generate a random number between 0 and 1, The probability of updating the random resistivity vector p is calculated by the following expression: [0.9*f / max(f)]+0.1; If the random number is less than the probability, an update is performed.
Citation Information
Patent Citations
Electric prospecting method and device
CN102426393A
Ground-air transient electromagnetic data three-dimensional frequency domain interpretation method based on one-dimensional inversion
CN115047530A