Nuclear magnetic resonance echo denoising method and device and electronic equipment

Abstract

The application provides a method and device for echo denoising of nuclear magnetic resonance and electronic equipment, and belongs to the technical field of nuclear magnetic resonance application. The method comprises the following steps: performing rotation correction on the real part echo in the original complex nuclear magnetic echo to obtain a corrected real part echo; determining a baseline correction echo; generating a pretreatment echo and performing variational mode decomposition on the pretreatment echo to obtain a plurality of mode components; screening effective mode components from the plurality of mode components based on the correlation, energy proportion and center frequency between each mode component and the pretreatment echo, and superimposing the screened effective mode components to obtain a reconstructed echo; determining a noise level based on a sampling point, determining a wavelet shrinkage threshold according to the noise level, and performing wavelet threshold shrinkage and reconstruction on the reconstructed echo based on the wavelet shrinkage threshold to obtain a target denoised echo. The above scheme can effectively improve the echo quality of the nuclear magnetic resonance echo.

Description

Technical Field

[0001] This invention relates to the field of nuclear magnetic resonance (NMR) application technology, specifically to NMR echo noise reduction methods, apparatus, and electronic equipment. Background Technology

[0002] In oil and gas exploration operations, nuclear magnetic resonance (NMR) logging is widely used to evaluate reservoir pore structure and fluid properties by analyzing the relaxation signals of hydrogen nuclei in formation fluids. For low signal-to-noise ratio NMR echoes obtained in complex environments such as ultra-deep formations, existing technologies often employ preprocessing methods such as smoothing filtering, wavelet denoising, or mode decomposition to suppress noise.

[0003] However, existing technologies struggle to simultaneously handle phase distortion, baseline shift, local spikes, and broadband random noise in echoes. Existing smoothing filters tend to weaken short-relaxation components; wavelet thresholding denoising is limited in effectiveness under strong noise and the threshold is fixed; mode decomposition suffers from difficulties in component selection, easily losing weak signals or causing mode aliasing; and existing processes often ignore the effects of phase and baseline errors, leading to echo distortion and producing spurious peaks and porosity deviations during inversion.

[0004] Therefore, there is an urgent need for a multi-stage conformal method that can sequentially eliminate phase distortion and baseline shift, remove local spikes, reliably distinguish between signal and noise modes, and reasonably determine the denoising threshold to improve the echo quality of nuclear magnetic resonance echoes. Summary of the Invention

[0005] To address the aforementioned problems, this application provides a method for denoising nuclear magnetic resonance (NMR) echoes to improve the echo quality of NMR signals. The technical solution is as follows: In a first aspect, the present invention provides a method for nuclear magnetic resonance echo denoising, comprising: The original complex NMR echo is obtained, and the real part echo in the original complex NMR echo is rotated and corrected to obtain the corrected real part echo. The rotation angle of the rotation correction is the negative angle of the principal phase angle. The principal phase angle is the arctangent angle between the accumulated imaginary part value and the accumulated real part value of the previous preset number of echo points of the original complex NMR echo. The baseline value is determined based on the sampling points at a preset ratio in the tail of the corrected real part echo, and the baseline value is subtracted from the corrected real part echo to obtain the baseline corrected echo. Local spikes are identified and replaced in the baseline-corrected echo to obtain a preprocessed echo. Variational mode decomposition is then performed on the preprocessed echo to obtain multiple modal components. Based on the correlation, energy ratio, and center frequency between each modal component and the preprocessed echo, effective modal components are selected from multiple modal components, and the selected effective modal components are superimposed to obtain the reconstructed echo; The noise level is determined based on the sampling points, the wavelet shrinkage threshold is determined based on the noise level, and the reconstructed echo is subjected to wavelet threshold shrinkage and reconstruction based on the wavelet shrinkage threshold to obtain the target denoised echo.

[0006] Combining the first aspect and the above implementation methods, in some possible implementation methods, the real part echo in the original complex NMR echo is rotated and corrected to obtain the corrected real part echo, including: Complex echo points are constructed based on the real and imaginary echoes in the original complex NMR echoes, and each complex echo point is rotated around the origin by the negative angle of the principal phase angle to obtain the rotated echo points; Extract the real component of each rotated echo point to obtain the corrected real echo, and extract the imaginary component of each rotated echo point to obtain the corrected imaginary echo.

[0007] Combining the first aspect and the above implementation methods, in some possible implementation methods, the baseline value is determined based on the sampling points of a preset proportion in the tail of the corrected real part echo, and the baseline value is subtracted from the corrected real part echo to obtain the baseline corrected echo, including: The tail median is determined based on the median of the sampling points, and the tail fluctuation scale is determined based on the median of the absolute values ​​of the differences between each sampling point and the tail median. Sampling points whose absolute value of the difference from the tail median is greater than a preset multiple of the tail fluctuation scale are removed, and the median of the remaining sampling points is determined as the baseline value. The baseline-corrected echo is obtained by subtracting the baseline value from each echo point of the corrected real part echo.

[0008] Combining the first aspect and the above implementation methods, in some possible implementation methods, local spike identification and replacement are performed on the baseline correction echo to obtain a preprocessed echo, including: Based on the median and dispersion of echo point values ​​within the neighborhood window of each echo point in the baseline-corrected echo, peak anomalies in the baseline-corrected echo are identified, and the values ​​of the identified peak anomalies are replaced with the median of the corresponding echo point values ​​within the neighborhood window to obtain the preprocessed echo.

[0009] Combining the first aspect and the above implementation methods, in some possible implementation methods, based on the correlation, energy ratio, and center frequency between each modal component and the preprocessed echo, effective modal components are selected from multiple modal components, and the selected effective modal components are superimposed to obtain the reconstructed echo, including: The correlation coefficient between each modal component and the preprocessed echo is obtained. The correlation coefficient is used to characterize the consistency between the attenuation trend of the modal component and the preprocessed echo. The ratio of the energy of each modal component to the total energy of all modal components is obtained to obtain the energy proportion of each modal component, and the center frequency of each modal component is determined based on the spectrum of each modal component. Modal components with correlation coefficients greater than or equal to the correlation threshold, energy proportions greater than or equal to the energy proportion threshold, and center frequencies less than the center frequency threshold are identified as effective modal components. If the number of effective modal components is lower than the preset number, modal components are added from the remaining modal components in descending order of energy percentage as effective modal components. The determined effective modal components are superimposed to obtain the reconstructed echo.

[0010] Combining the first aspect and the above implementation methods, in some possible implementation methods, the noise level is determined based on the sampling points, including: Peak removal is performed on the sampling points to obtain the removed sampling points; The first noise estimate is determined based on the median of the absolute value of the difference between the median of each sampled point after elimination and the median of the sampled points after elimination, and the second noise estimate is determined based on the median of the absolute value of the difference between the sampled points with a preset proportion in the tail of the corrected imaginary part echo and the median of the median. The maximum value between the first noise estimate and the weighted second noise estimate is determined as the noise level.

[0011] Combining the first aspect and the above implementation methods, in some possible implementation methods, a wavelet shrinkage threshold is determined based on the noise level, and wavelet threshold shrinkage and reconstruction are performed on the reconstructed echo based on the wavelet shrinkage threshold to obtain the target denoised echo, including: Wavelet decomposition is performed on the reconstructed echo to obtain low-frequency approximate components and high-frequency detail components at each layer. A general threshold is determined based on the noise level and the total number of echo points in the reconstructed echo, and a candidate threshold is determined based on the ratio of the square of the noise level to the fluctuation scale of the high-frequency detail component of the current layer. Determine the lower threshold, take the minimum value between the general threshold and the lower threshold, and multiply it by the scaling coefficient of the current layer to obtain the wavelet shrinkage threshold of the current layer; Based on the wavelet shrinkage threshold of each layer, the high-frequency detail components of each layer are subjected to soft threshold shrinkage, and the low-frequency approximate components are reconstructed with the shrunken high-frequency detail components of each layer to obtain the target denoised echo.

[0012] Combining the first aspect and the above implementation methods, in some possible implementation methods, after performing local spike identification and replacement on the baseline correction echo to obtain the preprocessed echo, the following steps are also included: If the mean of the early echo points of the preprocessed echo is negative, then the sign of the preprocessed echo is flipped. The robust signal-to-noise ratio (SNR) is determined by the ratio of the truncated mean of the absolute values ​​of the early echo points in the preprocessed echo after sign inversion to the median absolute deviation of the sampling points at a predetermined proportion in the tail of the corrected imaginary part echo. The robust SNR is used to characterize the signal-to-noise level of the original complex NMR echo.

[0013] Secondly, the present invention also provides a nuclear magnetic resonance echo denoising device, comprising: The data acquisition unit is used to acquire the original complex NMR echo, perform rotation correction on the real part echo in the original complex NMR echo to obtain the corrected real part echo. The rotation angle of the rotation correction is the negative angle of the principal phase angle. The principal phase angle is the arctangent angle between the accumulated imaginary part value and the accumulated real part value of the previous preset number of echo points of the original complex NMR echo. The data correction unit is used to determine the baseline value based on the sampling points at a preset ratio in the tail of the corrected real echo, and to subtract the baseline value from the corrected real echo to obtain the baseline corrected echo. The data decomposition unit is used to identify and replace local spikes in the baseline correction echo to obtain a preprocessed echo, and to perform variational mode decomposition on the preprocessed echo to obtain multiple mode components. The data reconstruction unit is used to select effective modal components from multiple modal components based on the correlation, energy ratio and center frequency between each modal component and the preprocessed echo, and to superimpose the selected effective modal components to obtain the reconstructed echo; The data generation unit is used to determine the noise level based on the sampling points, determine the wavelet shrinkage threshold based on the noise level, and perform wavelet threshold shrinkage and reconstruction on the reconstructed echo based on the wavelet shrinkage threshold to obtain the target denoised echo.

[0014] Thirdly, the present invention also provides an electronic device, including a data acquisition unit, a memory, and a processor, wherein, Data acquisition unit, used to acquire raw complex NMR echoes; Memory, used to store executable program code; The processor performs rotational correction on the real part of the original complex NMR echo to obtain a corrected real part echo; determines the baseline value based on a preset proportion of sampling points at the tail of the corrected real part echo, and subtracts the baseline value from the corrected real part echo to obtain a baseline corrected echo; performs local peak identification and replacement on the baseline corrected echo to obtain a preprocessed echo; performs variational mode decomposition on the preprocessed echo to obtain multiple modal components; selects effective modal components from the multiple modal components based on the correlation, energy ratio, and center frequency between each modal component and the preprocessed echo, and superimposes the selected effective modal components to obtain a reconstructed echo; determines the noise level based on the sampling points, determines the wavelet shrinkage threshold based on the noise level, and performs wavelet threshold shrinkage and reconstruction on the reconstructed echo based on the wavelet shrinkage threshold to obtain the target denoised echo.

[0015] Fourthly, the present invention also provides a computer-readable storage medium for storing a computer-readable program or instruction, which, when executed by a processor, can implement the steps of the nuclear magnetic resonance echo denoising method in any of the above implementations.

[0016] The beneficial effects of this invention are as follows: By rotating and correcting the original complex NMR echo using the negative angle of the principal phase angle, corrected real and imaginary echoes are obtained, effectively eliminating phase distortion. By selecting a preset proportion of sampling points at the tail of the corrected real echo to solve for the baseline value and subtracting it, a baseline-corrected echo is obtained, thus preventing amplitude misjudgment caused by baseline interference. Preprocessed echoes are obtained by identifying and replacing local spikes in the baseline-corrected echoes, eliminating abnormal spikes caused by occasional circuit pulses and avoiding spike contamination of subsequent decomposition results. Variational mode decomposition is performed on the preprocessed echoes to obtain multiple sets of modal components. Effective modes are selected and superimposed and reconstructed based on correlation, energy proportion, and center frequency, eliminating noise-dominated invalid modes and suppressing some noise interference. A wavelet shrinkage threshold is determined by estimating the noise level based on the sampling points, and wavelet threshold shrinkage processing is performed on the reconstructed echo based on the wavelet shrinkage threshold, thereby effectively generating high-quality target denoised echoes. Attached Figure Description

[0017] To more clearly illustrate the technical solutions in the embodiments of the present invention, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0018] Figure 1 A schematic flowchart of the nuclear magnetic resonance echo denoising method provided by the present invention; Figure 2 For the present invention Figure 1 A schematic diagram of an embodiment of S101; Figure 3 A comparison diagram of the distribution of robust signal-to-noise ratio after phase correction along the normalized depth position provided by the present invention; Figure 4 For the present invention Figure 1 A schematic diagram of an embodiment of S102; Figure 5 A comparison diagram of the distribution of the baseline of the corrected real part echo tail along the normalized depth position after phase correction provided by the present invention. Figure 6 For the present invention Figure 1 A schematic diagram of an embodiment of S104; Figure 7 For the present invention Figure 1 A schematic diagram of an embodiment of S104; Figure 8 This is a comparison image of VMD, wavelet, and VMD-wavelet fusion denoising echo provided by the present invention; Figure 9 Comparison of tail baseline correction and spike removal provided by the present invention; Figure 10 For the present invention Figure 1 A schematic diagram of an embodiment of S105; Figure 11 Comparison of NMR echo, T2 spectrum and porosity before and after VMD-wavelet fusion denoising provided by the present invention; Figure 12 A schematic diagram of the nuclear magnetic resonance echo noise reduction device provided by the present invention; Figure 13 This is a schematic diagram of the structure of the electronic device provided by the present invention. Detailed Implementation

[0019] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of the present invention, and not all of them. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without creative effort are within the scope of protection of the present invention.

[0020] In the description of the embodiments of the present invention, unless otherwise stated, "multiple" means two or more. "And / or" describes the relationship between related objects, indicating that there can be three relationships. For example, A and / or B can represent three situations: A exists alone, A and B exist simultaneously, and B exists alone.

[0021] The terms "first," "second," etc., used in the embodiments of this invention are for descriptive purposes only and should not be construed as indicating or implying their relative importance or implicitly specifying the number of technical features indicated. Therefore, a technical feature defined with "first" or "second" may explicitly or implicitly include at least one of that feature.

[0022] In this document, the term "embodiment" means that a particular feature, structure, or characteristic described in connection with an embodiment may be included in at least one embodiment of the invention. The appearance of this phrase in various places throughout the specification does not necessarily refer to the same embodiment, nor is it a separate or alternative embodiment mutually exclusive with other embodiments. It will be explicitly and implicitly understood by those skilled in the art that the embodiments described herein can be combined with other embodiments.

[0023] This application provides a method for nuclear magnetic resonance (NMR) echo denoising, wherein the execution subject of the method is an NMR echo denoising device or an electronic device equipped with an NMR echo denoising device. The following is a detailed description; it should be noted that the order of description of the following embodiments is not intended to limit the preferred order of the embodiments.

[0024] Please see Figure 1 , Figure 1 This is a schematic flowchart of the nuclear magnetic resonance echo denoising method provided in the embodiments of this application. Figure 1 As shown, the method in this application embodiment may include the following steps S101-S105: S101, acquire the original complex NMR echo, perform rotation correction on the real part echo in the original complex NMR echo to obtain the corrected real part echo.

[0025] In this embodiment, the original complex NMR echo refers to the echo signal, represented in complex form, acquired by a NMR logging instrument downhole. This signal consists of a series of echo points arranged in echo interval order, with each echo point containing both an in-phase component and a quadrature component. The real part of the echo refers to the sequence of in-phase components of each echo point in the original complex NMR echo, arranged in echo time order; that is, the projection sequence of the values ​​of each echo point onto the real axis of the complex plane. The imaginary part of the echo refers to the sequence of quadrature components of each echo point in the original complex NMR echo, arranged in echo time order; that is, the projection sequence of the values ​​of each echo point onto the imaginary axis of the complex plane.

[0026] Rotation correction refers to the process of rotating the complex numbers corresponding to each echo point of the original complex NMR echo around the origin by the same angle in the complex plane. This angle is the negative angle of the principal phase angle. The purpose of rotation correction is to rotate the effective relaxation signal in the original complex NMR echo as close to the real axis as possible. The corrected real part echo refers to the sequence of real components of each echo point in the rotation-corrected complex echo, and it carries the main effective NMR attenuation signal. The corrected imaginary part echo refers to the sequence of imaginary components of each echo point in the rotation-corrected complex echo, and it mainly reflects the quadrature residual components and noise.

[0027] In one feasible implementation, the echo time corresponding to the i-th echo point is taken as... The i-th echo point of the original complex NMR echo is represented by the value of the real part of the echo at that echo point. The real part and the imaginary part are the values ​​of the echo at that echo point. The imaginary part is a complex number. Ideally, the effective NMR attenuation signal should be mainly concentrated in the real axis direction, and the imaginary part should be close to zero mean noise. However, due to instrument phase error, downhole environmental interference, and low signal-to-noise ratio during actual acquisition, the effective NMR attenuation signal has a phase deviation relative to the real axis direction, so the effective NMR attenuation signal may not fall completely in the real axis direction. Therefore, rotation correction is required.

[0028] For example, the specific form of the formula for calculating the principal phase angle in step S101 above is as follows:

[0029] in, The principal phase angle is used to characterize the phase deviation of the effective NMR attenuation signal relative to the real axis direction. It is a two-parameter arctangent function, which can determine the quadrant of the principal phase angle based on the signs of the two parameters, and its value range is (-π, π]. The number of pre-set echo points used for principal phase angle estimation. The range is 20 to 80, preferably 50. The signal amplitude of the early echo points is relatively high and the signal-to-noise ratio is relatively good. Therefore, selecting the first N echo points as the phase estimation interval can improve the estimation accuracy of the principal phase angle. Let the imaginary part of the echo at the i-th echo point be taken as the value; Let be the real part echo value of the i-th echo point.

[0030] For example, the specific form of the calculation formula for rotation correction in step S101 above is as follows:

[0031] in, Let i be the complex echo point corresponding to the i-th echo point. The rotated echo point is obtained after rotation correction for the i-th echo point. It is a natural constant. The imaginary unit, This represents the negative angle of the principal phase angle rotated around the origin. After performing the above rotation on all echo points of the original complex NMR echo, the real component of each rotated echo point is extracted to obtain the corrected real part echo, and the imaginary component of each rotated echo point is extracted to obtain the corrected imaginary part echo.

[0032] For example, for an original complex NMR echo with 400 echo points, the first 50 echo points are selected as the phase estimation interval. If the sum of the imaginary echo values ​​of these 50 echo points is 120 and the sum of the real echo values ​​is 600, then the principal phase angle φ = arctan2(120, 600) ≈ 0.197 radians. After rotating the complex echo points corresponding to each echo point around the origin by -0.197 radians, the effective NMR attenuation signal is rotated to the real axis direction. After rotation correction, if the corrected real echo shows a clear attenuation trend, it indicates that the effective NMR attenuation signal has been basically rotated to the real axis direction; if the corrected imaginary echo still has strong oscillations or local spikes, it indicates that there are still orthogonal components remaining or noise interference.

[0033] S102, determine the baseline value based on the sampling points of the tail of the corrected real part echo with a preset ratio, and subtract the baseline value from the corrected real part echo to obtain the baseline corrected echo.

[0034] In this embodiment, the tail preset ratio refers to the proportion of echo points located at the end of the echo time and used for baseline estimation out of the total number of echo points in the corrected real part echo. The value of the tail preset ratio ranges from 0.20 to 0.35, preferably 0.25. The baseline value is a numerical value characterizing the non-zero offset level of the tail of the corrected real part echo. Theoretically, NMR echoes should gradually decay to zero as the echo time increases. The values ​​of the sampling points of the tail preset ratio should approach zero after deducting the influence of noise. If its values ​​deviate from zero overall, then the deviation is the baseline value. The baseline-corrected echo refers to the echo sequence obtained by subtracting the baseline value from the values ​​of each echo point of the corrected real part echo. The tail of the baseline-corrected echo is pulled back to near zero.

[0035] In one feasible implementation, the tail sampling points participating in baseline estimation are first determined. Assuming there are N echo points in the corrected real part echo, the number of tail sampling points, Ntail, is determined according to the following formula:

[0036] in, The number of tail sampling points; To correct the total number of echo points in the real part of the echo; Preset proportions for the tail section; To maximize the function, this formula ensures that when the total number of echo points is small, at least 20 tail sampling points participate in baseline estimation to guarantee the statistical stability of the baseline estimation. For example, when =400 and When the value is 0.25, =max(20, 100)=100, which means taking the last 100 echo points of the corrected real part echo as the tail sampling points.

[0037] The median of all tail sampling points is calculated as the tail median, and the median absolute deviation (MAD) is used to estimate the tail volatility. The median absolute deviation is the median of a dataset, calculated as the absolute values ​​of the differences between each data point and the median of that dataset. The median absolute deviation relies primarily on median statistics and is insensitive to isolated peaks. For example, the specific form of the formula for calculating the tail volatility is as follows:

[0038] in, is the tail fluctuation scale, used to characterize the dispersion of the tail sampling point values; 1.4826 is the scaling coefficient for converting the absolute deviation of the median into a scaling estimator consistent with the standard deviation under Gaussian noise conditions. To perform the median calculation; To correct the value of the real part echo at the i-th echo point; The median of the tail; It is the set of all tail sampling points.

[0039] Furthermore, to avoid potential spikes or outliers in the tail sampling points affecting the baseline estimation, sampling points whose absolute value differs from the tail median by a preset multiple κ are removed; that is, points satisfying the following criteria are retained. The tail sampling points are selected, where the value of κ ranges from 2.5 to 4.0, preferably 3.0. The median of the remaining tail sampling point values ​​after removal is determined as the baseline value, and this baseline value is subtracted from the values ​​of each echo point of the corrected real echo. That is, the value of the baseline corrected echo at the i-th echo point is the difference between the value of the corrected real echo at the i-th echo point and the baseline value.

[0040] For example, if the median of 100 tail sampling points is 35, the tail fluctuation scale is 20, and κ is 3.0, then tail sampling points whose values ​​fall outside the interval [35-3.0×20, 35+3.0×20] (i.e., [-25, 95]) are identified as spike anomalies and removed. If the median of the remaining 92 tail sampling points after removal is 32, then the baseline value is 32. Subtracting 32 from the values ​​of the 400 echo points of the corrected real part echo yields the baseline-corrected echo. After this step, the tail of the baseline-corrected echo is pulled back to near zero, which helps to weaken the false tail of the long transverse relaxation time caused by baseline drift, providing a more stable input for subsequent variational mode decomposition, wavelet thresholding denoising, and transverse relaxation time spectrum inversion.

[0041] S103 performs local spike identification and replacement on the baseline correction echo to obtain a preprocessed echo, and performs variational mode decomposition on the preprocessed echo to obtain multiple mode components.

[0042] In this embodiment, local spike identification and replacement refers to identifying local sudden jump points in the baseline-corrected echo that are short in duration, high in amplitude, and without obvious regularity, based on the local statistical characteristics of the echo point values ​​within the sliding neighborhood window, and replacing the values ​​of the identified jump points with the median of the corresponding echo point values ​​within the neighborhood window. The preprocessed echo refers to the echo sequence obtained after local spike identification and replacement; in the preprocessed echo, local spikes have been eliminated while the main attenuation envelope remains unchanged.

[0043] Variational Mode Decomposition (VMD) is a non-recursive, adaptive signal decomposition method. It constructs and solves a variational problem to decompose the input signal into several modal components, each with finite bandwidth and revolving around a central frequency. VMD adaptively determines the center frequency and bandwidth of each modal component. A modal component refers to a sub-signal output by VMD; each modal component represents the wave component within a specific frequency band in the preprocessed echo, and the superposition of all modal components equals the preprocessed echo.

[0044] In one feasible implementation, local spike identification and replacement are first performed on the baseline correction echo. Specifically, for the i-th echo point in the baseline correction echo, a neighborhood window with length Ω centered on that echo point is taken. Select 7 to 15 echo points, preferably 9, meaning the neighborhood window includes the i-th echo point and the four echo points before and after it; calculate the median of the values ​​of each echo point within the neighborhood window as the local median mi, and use the absolute deviation of the median to estimate the local robust standard deviation. ,Right now ,in Let y be the value of the j-th echo point within the neighborhood window; when the value of the i-th echo point is y i satisfy When this happens, the echo point is identified as a spike anomaly, and its value is replaced with the local median. Otherwise, the original value of the echo point remains unchanged, where The value ranges from 3.0 to 5.0, with 4.0 being preferred.

[0045] For example, if the values ​​of the nine echo points within the neighborhood window of a certain echo point are 10, 12, 11, 13, 300, 12, 11, 10, and 12 respectively, then the local median is... =12, Locally robust standard deviation Since |300-12|=288>4.0×1.48, this echo point is identified as a spike anomaly and replaced with 12, while the remaining echo points within the neighborhood window retain their original values. This processing only replaces local anomalous spikes and does not perform mean smoothing on the entire echo, thus avoiding damage to the main attenuation envelope of the NMR echo. After performing the above identification and replacement on all echo points of the baseline correction echo, the preprocessed echo is obtained.

[0046] Variational mode decomposition (VMD) is performed on the preprocessed echo. VMD represents the preprocessed echo as a superposition of K modal components, i.e., the preprocessed echo. satisfy:

[0047] in, To preprocess the echo values ​​at echo time t; K is the number of modal components, which can be selected between 3 and 10 according to the echo complexity, preferably between 5 and 7; (t) represents the value of the k-th mode component at echo time t. Variational mode decomposition iteratively solves for the center frequency and corresponding waveform of each mode component, minimizing the sum of the estimated bandwidths of each mode component and reconstructing the preprocessed echo from the sum of all mode components. After decomposition, K mode components are obtained.

[0048] It should be noted that traditional variational mode decomposition denoising often adopts a strategy of retaining low-order modes and discarding high-order modes. However, this strategy is not entirely suitable for low signal-to-noise ratio NMR echoes because the effective decay signal with short transverse relaxation time itself has a rapid change characteristic and may be distributed in higher frequency mode components. If high-frequency mode components are simply deleted, the short transverse relaxation time signal is easily weakened. If too many mode components are retained, random noise will remain. Therefore, this application further introduces an effective mode screening mechanism after variational mode decomposition, as detailed in step S104 below.

[0049] S104. Based on the correlation, energy ratio and center frequency between each modal component and the preprocessed echo, effective modal components are selected from multiple modal components, and the selected effective modal components are superimposed to obtain the reconstructed echo.

[0050] In this embodiment, correlation refers to the degree of consistency between each modal component and the preprocessed echo in terms of overall attenuation trend, and correlation is quantified by correlation coefficient. Energy ratio refers to the ratio of the energy of each modal component to the total energy of all modal components, used to characterize the contribution of the corresponding modal component to the reconstruction of the preprocessed echo. Center frequency refers to the frequency corresponding to the energy centroid of the spectrum of each modal component, used to characterize the frequency band position of the corresponding modal component. Effective modal component refers to the modal component containing effective NMR attenuation information selected from multiple modal components. Reconstructed echo refers to the echo sequence obtained by superimposing all effective modal components point by point. The reconstructed echo retains the main attenuation trend of the preprocessed echo and removes noise modal components.

[0051] In one feasible implementation, for the K modal components obtained by variational mode decomposition, three screening features are calculated for each modal component.

[0052] First, calculate the modal components. With preprocessed echo Correlation coefficient between :

[0053] in, Let be the correlation coefficient between the k-th modal component and the preprocessed echo, with a value range of [-1, 1]. The mean of the k-th modal component; The mean value of the preprocessed echo; This involves traversing all echo points. A higher correlation coefficient indicates that the modal component is more consistent with the overall attenuation trend of the preprocessed echo, and is more likely to contain a valid signal.

[0054] Second, calculate the energy percentage of each modal component. The energy of the k-th modal component is the sum of the squares of its values ​​at all echo points. The total energy of all modal components is the sum of the energies of each modal component. The energy percentage is... for:

[0055] in, For the first The energy percentage of each modal component; For the first Energy of each modal component; This represents the total energy of all modal components. Modal components with a larger energy proportion usually contribute more significantly to echo reconstruction.

[0056] Third, the center frequency is estimated using the spectra of the modal components. A Fourier transform is performed on each modal component to obtain its corresponding spectrum, and the center frequency is the weighted average frequency of the squares of the spectrum moduli.

[0057] in, For the first The center frequency of each modal component; For the first Modal components (t) The value of the spectrum obtained by Fourier transform at frequency f; f traverses all frequency points of the spectrum. Generally speaking, modal components with extremely high center frequencies and low correlation with preprocessed echoes are more likely to be random noise, while modal components with moderate center frequencies and high energy proportions may contain effective attenuation information with short transverse relaxation times.

[0058] Based on the three screening features mentioned above, a relevance threshold is set. Energy percentage threshold and center frequency threshold , will simultaneously satisfy ≥ , ≤ and ≥ The modal components are determined as effective modal components, that is, the effective mode set formed by the subscripts of the effective modal components is:

[0059] in, For effective mode set, For correlation threshold, For the center frequency threshold, The three thresholds are the energy percentage thresholds. They can be determined based on the actual work area calibration data. For example, the correlation threshold can be set to 0.3, the energy percentage threshold can be set to 0.05, and the center frequency threshold can be set to the lower limit frequency of the noise main frequency band.

[0060] It should be noted that the modal screening does not adopt a single low-frequency retention principle, but comprehensively considers three characteristics: correlation, energy proportion, and center frequency, in order to reduce the risk of accidentally deleting information with rapid decay due to short transverse relaxation time. When the number of effective modal components selected is too small, the unselected modal components are further sorted in descending order of energy proportion and added as effective modal components in turn until the number of effective modal components reaches the preset number, so as to ensure that the main decay information is not excessively weakened.

[0061] By superimposing all valid mode components in the valid mode set point by point, the reconstructed echo is obtained:

[0062] in, To reconstruct the echo values ​​at echo time t; S is the set of effective modes; (t) represents the value of the kth effective mode component at echo time t. This step retains the main mode components related to the preprocessed echo attenuation trend and removes high-frequency, low-correlation, and low-energy noise mode components, thereby achieving the first stage of random noise suppression.

[0063] For example, suppose variational mode decomposition yields six modal components with correlation coefficients of 0.92, 0.85, 0.61, 0.40, 0.12, and 0.05, and energy proportions of 0.45, 0.30, 0.12, 0.06, 0.04, and 0.03, respectively. The center frequencies correspond to frequency bands from low to high, and the center frequencies of the 5th and 6th modal components are higher than the center frequency threshold. Based on a correlation threshold of 0.3 and an energy proportion threshold of 0.05, the 1st to 4th modal components are determined as valid modal components, while the 5th and 6th modal components are discarded as noise modal components due to their excessively high center frequencies and low correlation. The reconstructed echo is obtained by superimposing the 1st to 4th valid modal components point by point.

[0064] S105, determine the noise level based on the sampling points, determine the wavelet shrinkage threshold based on the noise level, and perform wavelet threshold shrinkage and reconstruction on the reconstructed echo based on the wavelet shrinkage threshold to obtain the target denoised echo.

[0065] In this embodiment, the noise level refers to a robust estimate of the noise intensity in the original complex NMR echo, which forms the basis for subsequent wavelet shrinkage threshold calculation and signal-to-noise ratio calculation. The wavelet shrinkage threshold is a threshold value used in wavelet thresholding to determine whether wavelet coefficients in the high-frequency detail components should be preserved or compressed. Wavelet coefficients with absolute values ​​less than the wavelet shrinkage threshold are considered noise and compressed to zero. Wavelet threshold shrinkage and reconstruction involves performing wavelet decomposition on the reconstructed echo to obtain low-frequency approximate components and high-frequency detail components at each level. The wavelet coefficients of each high-frequency detail component are then shrunk according to the wavelet shrinkage threshold, and finally, the low-frequency approximate components and the shrunk high-frequency detail components are reconstructed using wavelet technology. The target denoised echo refers to the denoised echo obtained after wavelet thresholding and reconstruction. The target denoised echo is used as input for transverse relaxation time spectrum inversion to recover the transverse relaxation time (T2) spectrum and calculate interpretation parameters such as porosity. The transverse relaxation time refers to the time constant in nuclear magnetic resonance that characterizes the rate of decay of the transverse magnetization vector of the pore fluid. The transverse relaxation time spectrum refers to the signal amplitude distribution corresponding to different transverse relaxation times.

[0066] In one feasible implementation, the noise level is first determined. After phase correction, the effective NMR signal is mainly located in the corrected real part echo, while the corrected imaginary part echo should mainly consist of orthogonal residuals and noise. Therefore, in the presence of imaginary part data, the noise level is estimated by using a predetermined proportion of sampling points at the tail of the corrected imaginary part echo. Specifically, a predetermined proportion of sampling points at the tail of the corrected imaginary part echo is taken to form the imaginary part tail sampling set. First, Hampel peak removal is performed on the imaginary part tail sampling set (the specific process of Hampel peak removal is described in step S103 above, and will not be repeated here). Then, the robust standard deviation of the removed sampling points is calculated as the imaginary part noise estimate. :

[0067] in, This is the estimated value of the imaginary part noise; 1 represents the value of each sampling point in the imaginary tail segment sampling set after removing spikes; 1.4826 is the scaling coefficient for converting the absolute deviation of the median into a scaling estimator consistent with the standard deviation under Gaussian noise conditions.

[0068] Meanwhile, to avoid underestimation of noise due to abnormally low corrected imaginary echo or phase distortion, a robust standard deviation of sampling points at a preset proportion of the baseline corrected echo tail is calculated as an estimate of the real noise. :

[0069] in, This is the estimated value of the real part noise; The set of sampling points with a preset ratio for the tail of the baseline correction echo.

[0070] The final noise level is the maximum of three values: the imaginary noise estimate, the weighted real noise estimate, and a small quantity to prevent the noise level from being zero.

[0071] in, Noise level; is the weighting coefficient for the real part noise estimate. The value ranges from 0.3 to 0.8, preferably 0.5; To prevent small values ​​of zero noise, ε is set between 10^-12 and 10^-8, preferably 10^-12. This noise estimation method is more robust than the standard deviation because ultra-deep NMR echo tails often contain spikes and outliers. Directly using the standard deviation would be amplified by extreme values, while the median absolute deviation estimation mainly relies on the median and is not sensitive to isolated spikes, making it more suitable for echo data with low signal-to-noise ratios and many outliers.

[0072] The wavelet shrinkage threshold is determined based on the noise level. The reconstructed echo undergoes J-level wavelet decomposition, yielding the J-th level low-frequency approximation component and the 1st to J-th level high-frequency detail components. The low-frequency approximation component primarily reflects the overall attenuation trend of the reconstructed echo, while the high-frequency detail components contain random noise and local oscillations. The focus of wavelet thresholding is to shrink the high-frequency detail components at each level while preserving the low-frequency attenuation trend of the low-frequency approximation component as much as possible. The number of wavelet decomposition levels J can be determined based on the total number of echo points in the reconstructed echo; for example, 3 to 5 levels can be used. Wavelet basis functions such as sym8 or db4, which have good compact support and symmetry, can be selected.

[0073] For determining the wavelet shrinkage threshold, this application adopts a combination of noise estimation and scale-adaptive thresholding.

[0074] For example, the specific form of the formula for calculating the general threshold is as follows:

[0075] in, This is a general threshold; The noise intensity estimate can be directly obtained from the above noise levels. It can also be estimated from the median absolute deviation of the high-frequency detail coefficients of the highest frequency layer; The total number of echo points for reconstructing the echo.

[0076] Meanwhile, to avoid the universal threshold being too strong in local low signal-to-noise ratio echoes, this application constructs candidate thresholds by incorporating the Bayesian thresholding concept:

[0077] in, The candidate threshold; The square of the noise intensity estimate; The fluctuation scale of the high-frequency detail components of the current layer can be the standard deviation of the high-frequency detail coefficients of the current layer or a robust scale estimate based on the median absolute deviation. To prevent small quantities with a denominator of zero.

[0078] Finally, the wavelet shrinkage threshold for the j-th layer high-frequency detail components is:

[0079] in, The wavelet shrinkage threshold for the high-frequency detail components of the j-th layer; is the scaling factor for the j-th layer, used to fine-tune the threshold according to the differences in noise distribution at different decomposition levels; This is the lower limit control coefficient, used to prevent noise residue caused by excessively small candidate thresholds; min is the minimum value operation; max is the maximum value operation. This threshold form is achieved through... A fixed ratio is maintained to ensure that the threshold is greater than or equal to the general threshold to avoid noise residue caused by an excessively small threshold. Furthermore, the general threshold is used as the upper limit through a min operation to avoid loss of effective attenuation information caused by an excessively large threshold.

[0080] Soft thresholding is applied to the high-frequency detail components of each layer based on the wavelet shrinkage threshold for each layer. Soft thresholding refers to setting the absolute values ​​of wavelet coefficients that do not exceed the wavelet shrinkage threshold to zero, and compressing the absolute values ​​of wavelet coefficients that exceed the wavelet shrinkage threshold towards zero by the wavelet shrinkage threshold size. For example, the specific form of the calculation formula for soft thresholding is as follows:

[0081] in, For the high-frequency detail components after the j-th layer is shrunk; For the high-frequency detail components before the j-th layer shrinks; This is the sign function, used to preserve the sign of wavelet coefficients; Let be the wavelet shrinkage threshold for the j-th layer. Soft thresholding can continuously compress small-amplitude high-frequency disturbances, avoiding discontinuous oscillations caused by hard thresholding. After thresholding, the low-frequency approximation component of the j-th layer and the shrunken high-frequency detail components of each layer are reconstructed using wavelets to obtain the target denoised echo. The target denoised echo maintains the main attenuation trend while suppressing high-frequency random noise and local small-amplitude oscillations.

[0082] For example, suppose the total number of echo points in the reconstructed echo is N=400, and the noise level is... =15, then the general threshold =15×√(2ln400)≈15×3.46≈51.9; If the fluctuation scale σD of the high-frequency detail component in the second layer is 30 and the lower limit control coefficient c=0.3, then the candidate threshold is... =15² / 30=7.5, the lower threshold is max(7.5, 0.3×51.9)=max(7.5, 15.57)=15.57, the wavelet shrinkage threshold of the second layer is the scaling factor multiplied by min(51.9, 15.57)=15.57. If the scaling factor is 1, then the wavelet shrinkage threshold of the second layer is 15.57. Wavelet coefficients with an absolute value less than 15.57 in the high-frequency detail components of the second layer are set to zero, and wavelet coefficients with an absolute value not less than 15.57 are compressed to zero by 15.57.

[0083] In summary, this application performs rotational correction on the original complex NMR echo using the negative angle of the principal phase angle to obtain corrected real and imaginary echoes, effectively eliminating phase distortion. By selecting a preset proportion of sampling points at the tail of the corrected real echo to solve for the baseline value and subtracting it, a baseline-corrected echo is obtained, thus preventing amplitude misjudgment caused by baseline interference. Preprocessed echoes are obtained by identifying and replacing local spikes in the baseline-corrected echoes, eliminating abnormal spikes caused by occasional circuit pulses and avoiding spike contamination of subsequent decomposition results. Variational mode decomposition is performed on the preprocessed echoes to obtain multiple sets of modal components. Effective modes are selected and superimposed and reconstructed based on correlation, energy proportion, and center frequency, eliminating noise-dominated invalid modes and suppressing some noise interference. A wavelet shrinkage threshold is determined by estimating the noise level based on the sampling points, and wavelet threshold shrinkage processing is applied to the reconstructed echo based on the wavelet shrinkage threshold, thereby effectively generating high-quality target denoised echoes.

[0084] Please see Figure 2 , Figure 2 This is for the present invention Figure 1 A schematic diagram of an embodiment of S101. (See attached diagram.) Figure 2 As shown, the method in this application embodiment may include the following steps S201-S202: S201, construct complex echo points based on the real and imaginary echoes in the original complex nuclear magnetic resonance echoes, and rotate each complex echo point around the origin by the negative angle of the principal phase angle to obtain the rotated echo points.

[0085] In this embodiment, an echo point refers to the signal value corresponding to each sampling time obtained by the nuclear magnetic resonance logging instrument sampling the nuclear magnetic resonance signal according to the echo interval. A complex echo point is a complex number formed by taking the real part of the echo at the corresponding echo point as the real part and the imaginary part of the echo at the corresponding echo point as the imaginary part. A complex echo point corresponds to a point in the complex plane, where the x-coordinate is the value of the real part echo and the y-coordinate is the value of the imaginary part echo. The principal phase angle is the angle characterizing the phase deviation of the effective nuclear magnetic attenuation signal in the original complex nuclear magnetic echo relative to the real axis direction of the complex plane. The principal phase angle is the arctangent angle between the accumulated imaginary part value and the accumulated real part value of the previous preset number of echo points in the original complex nuclear magnetic echo. The specific form of the formula for calculating the principal phase angle is given in step S101 above and will not be repeated here. The rotated echo point refers to the complex number obtained by rotating the complex echo point around the origin in the complex plane by the negative angle of the principal phase angle. The effective nuclear magnetic attenuation signal component of the rotated echo point is rotated to the real axis direction.

[0086] In one feasible implementation, for the i-th echo point of the original complex NMR echo, the value of the real part of the echo at that echo point is... The values ​​of the real and imaginary parts of the echo at this echo point. As the imaginary part, construct complex echo points. = +j Where j is the imaginary unit. Then, the arctangent angle between the accumulated imaginary part and the accumulated real part of the original complex NMR echo points is calculated to obtain the principal phase angle. The preset number ranges from 20 to 80, preferably 50. Multiplying each complex echo point by e raised to the power of -jφ rotates each complex echo point around the origin in the complex plane by the negative angle of the principal phase angle, resulting in the rotated echo point. For example, if the principal phase angle φ = 0.197 radians, then each complex echo point is uniformly rotated clockwise by 0.197 radians around the origin, and the effective NMR attenuation signal component of the rotated echo point is rotated to the real axis direction.

[0087] S202, extract the real component of each rotated echo point to obtain the corrected real echo, and extract the imaginary component of each rotated echo point to obtain the corrected imaginary echo.

[0088] In the embodiments of this application, the real component refers to the projection value of the rotated echo point onto the real axis of the complex plane, that is, the real part value of the rotated echo point. The imaginary component refers to the projection value of the rotated echo point onto the imaginary axis of the complex plane, that is, the imaginary part value of the rotated echo point.

[0089] In one feasible implementation, for each rotated echo point obtained in step S201, the operations of extracting the real part and extracting the imaginary part are performed respectively. Specifically, the real components of all echo points are arranged in the order of echo time to form the corrected real part echo; the imaginary components of all echo points are arranged in the order of echo time to form the corrected imaginary part echo.

[0090] For example, if the rotated echo point is 580 + j25, then the real component of the echo point is 580, and the imaginary component is 25. After correction, the corrected real echo carries the main effective NMR attenuation signal and exhibits a clear attenuation trend. The corrected imaginary echo mainly reflects the orthogonal residual component and noise. The tail of the corrected imaginary echo will serve as an important basis for noise estimation in subsequent noise estimation.

[0091] Please refer to the following: Figure 3 , Figure 3 This is a comparison chart showing the distribution of robust signal-to-noise ratio along the normalized depth position after phase correction. (See figure) Figure 3 As shown, the horizontal axis represents the normalized depth position, ranging from 0 to 1, used to characterize the relative depth of each sampling point within the measurement well section; the vertical axis represents the robust signal-to-noise ratio (Robust SNR). For the specific form of the robust SNR calculation formula, please refer to the above implementation method of robust SNR, which will not be repeated here.

[0092] Figure (a) corresponds to the Tarim ultra-deep interval, Figure (b) corresponds to the Junggar mid interval, and Figure (c) corresponds to the Junggar bot interval. Each sub-figure contains three types of curves: the gray point-wise SNR curve, which represents the robust SNR at each depth sampling point without smoothing; the green solid line is the moving median curve (window n=31), which represents the local trend of robust SNR obtained by taking the median of the point-wise SNR with 31 sampling points as the moving window; and the two horizontal dashed lines are the reference lines for SNR equal to 1 (SNR=1) and SNR equal to 2 (SNR=2), respectively, used to define the low SNR interval.

[0093] Depend on Figure 3 As can be seen, the point-to-point signal-to-noise ratios (SNRs) of the three well sections all exhibit drastic random fluctuations and numerous isolated upward spikes, indicating that the SNR distribution along the depth in ultra-deep NMR echoes is extremely uneven, and conventional point-to-point estimation is easily affected by spike anomalies. In contrast, the moving median curve fluctuates gently and can stably characterize the overall SNR level of each well section: In Figure (a), the moving median of the ultra-deep Tarim well section is generally located near or below the reference line with an SNR of 2, indicating that the well section is in a low SNR state; In Figure (b), the moving median of the central Junggar well section is mostly below the reference line with an SNR of 1, only slightly rising after about 0.9 normalized depth; In Figure (c), the moving median of the bottom Junggar well section shows a significant rise in the normalized depth range of about 0.3 to 0.8 and exceeds the reference line with an SNR of 2, indicating that there is a reservoir section with a relatively high SNR in this range, while the rest of the positions are still at a low SNR level.

[0094] Figure 3 This indicates that ultra-deep NMR data are generally in a low signal-to-noise ratio (SNR) state and the SNR varies significantly along the depth. The robust SNR calculation method based on median statistics can effectively suppress the interference of isolated spikes on the SNR estimation, thereby truly reflecting the SNR level of each well section and providing a basis for subsequent adaptive estimation of noise level and quantitative evaluation of denoising effect.

[0095] In summary, this application constructs complex echo points from the echo points of the original complex NMR echo and rotates them around the origin by a negative angle of the principal phase angle, so that the effective NMR attenuation signal is concentrated in the real part of the rectified echo, while the imaginary part of the rectified echo mainly reflects the noise component. This provides the prerequisite for subsequent robust noise estimation and signal-to-noise ratio calculation based on the tail of the imaginary part of the rectified echo, and improves the ability to distinguish between signal and noise components under low signal-to-noise ratio conditions.

[0096] In one feasible implementation, when performing local spike identification and replacement on the baseline correction echo to obtain the preprocessed echo, the following is also specifically performed: In this embodiment, the neighborhood window refers to a local echo interval centered on the echo point to be judged in the baseline-corrected echo, extending forward and backward by the same number of echo points. The length of the neighborhood window is 7 to 15 echo points, preferably 9. The dispersion refers to the degree of deviation of each echo point value within the neighborhood window from the median value within that neighborhood window. The dispersion is quantified using a local robust standard deviation based on the absolute deviation of the median. A spike outlier is an echo point whose absolute value differs from the median value of the corresponding echo point within the neighborhood window by a preset multiple K times the dispersion. Spike outliers are characterized by short duration, high amplitude, and no obvious regularity, resulting in sudden, localized jumps.

[0097] In one feasible implementation, for the i-th echo point in the baseline correction echo, firstly, its neighborhood window is taken, which includes the i-th echo point and k echo points before and after it. When the neighborhood window length is 9, k is 4. Then, the median of the values ​​of each echo point within the neighborhood window is calculated as the local median, and the median of the absolute values ​​of the differences between the values ​​of each echo point within the neighborhood window and the local median is calculated. This median is then multiplied by a scaling factor of 1.4826 to obtain the local robust standard deviation characterizing the degree of dispersion. The value of the i-th echo point is compared with the local median. If the absolute value of the difference between the two is greater than a preset multiple K times the local robust standard deviation, the i-th echo point is determined to be a peak anomaly, and the value of the i-th echo point is replaced with the local median. Otherwise, the original value of the i-th echo point is kept unchanged. The preset multiple K ranges from 3.0 to 5.0, preferably 4.0.

[0098] The specific forms of the calculation formulas for the local median, local robust standard deviation, and peak outlier detection mentioned above are provided in step S103 above and will not be repeated here.

[0099] For example, if the value of an echo point in the baseline-corrected echo is 300, and the values ​​of the nine echo points in its neighborhood window are 10, 12, 11, 13, 300, 12, 11, 10, and 12 respectively, then the local median is 12, and the local robust standard deviation is approximately 1.48. Since |300-12|=288 is greater than 4.0×1.48, this echo point is identified as a spike anomaly and replaced with 12. Performing the above identification and replacement on all echo points in the baseline-corrected echo yields the preprocessed echo. This processing only replaces local anomalous spikes and does not perform mean smoothing on the entire echo, thus avoiding disruption of the main attenuation envelope of the NMR echo.

[0100] In summary, this application identifies and replaces spike anomalies based on the median and dispersion of echo points within a neighborhood window. This makes the discrimination of spike anomalies dependent on the median statistic, which is insensitive to isolated anomalies. This can eliminate local spikes without destroying the main attenuation envelope of the NMR echo, avoiding the amplified impact of spike anomalies on subsequent variational mode decomposition, wavelet shrinkage threshold estimation, and transverse relaxation time spectrum inversion, and reducing the risk of false peaks or spectral anomalies in the inversion results.

[0101] Please see Figure 4 , Figure 4 This is the present invention. Figure 1 A schematic diagram of an embodiment of S102. (See attached diagram.) Figure 4 As shown, the method in this application embodiment may include the following steps S401-S403: In this embodiment, the tail median refers to the median of all sampling points at a predetermined tail proportion in the corrected real echo. The tail median is used to initially characterize the offset level of the tail of the corrected real echo. The tail ripple scale is a robust scale estimate obtained by multiplying the median of the absolute values ​​of the differences between each sampling point and the tail median by a scaling factor of 1.4826. The tail ripple scale is used to characterize the dispersion of the tail sampling point values. The specific form of the formula for calculating the tail ripple scale is described in step S102 above and will not be repeated here.

[0102] In one feasible implementation, tail sampling points are first extracted from the end of the corrected real echo according to a preset tail ratio. The preset tail ratio ranges from 0.20 to 0.35, preferably 0.25. The number of tail sampling points is the product of the total number of echo points of the corrected real echo and the preset tail ratio. When the product is less than 20, it is taken as 20 to ensure the statistical sample size for estimation.

[0103] For example, if the real echo correction has 400 echo points and the tail preset ratio is 0.25, then the last 100 echo points are taken as sampling points. Then, the values ​​of the 100 sampling points are sorted from smallest to largest, and the value in the middle position after sorting is taken as the tail median, for example, the tail median is 35; then the absolute value of the difference between each sampling point and the tail median 35 is calculated to obtain 100 absolute deviation values. The median of these 100 absolute deviation values ​​is taken and multiplied by the scaling factor 1.4826 to obtain the tail fluctuation scale, for example, the tail fluctuation scale is 20.

[0104] S402, the sampling points whose absolute value of the difference from the tail median is greater than a preset multiple of the tail fluctuation scale are removed, and the median of the remaining sampling points is determined as the baseline value.

[0105] In the embodiments of this application, the preset multiple refers to the threshold coefficient used to determine whether the tail sampling point is a peak anomaly point. The value range of the preset multiple is 2.5 to 4.0, preferably 3.0.

[0106] In one feasible implementation, a retention interval is constructed with the tail median as the center and a preset multiple of the tail fluctuation scale as the half-width. That is, the lower limit of the retention interval is the tail median minus the preset multiple of the tail fluctuation scale, and the upper limit of the retention interval is the tail median plus the preset multiple of the tail fluctuation scale. Sampling points whose values ​​fall outside the retention interval are identified as peak outliers and removed from the tail sampling points. The median of the remaining sampling points within the retention interval is determined as the baseline value.

[0107] For example, using the values ​​from step S401 above, with a tail median of 35, a tail fluctuation scale of 20, and a preset multiple of 3.0, the retention interval is [-25, 95]. If 8 out of 100 sampling points fall outside this interval, these 8 sampling points are removed, and the median of the remaining 92 sampling points is 32. Therefore, the baseline value is determined to be 32. Since the baseline value is taken from the median of the remaining sampling points after removing the peak outliers, the baseline estimation is less susceptible to the influence of tail peak outliers.

[0108] S403, subtract the baseline value from each echo point of the corrected real echo to obtain the baseline corrected echo.

[0109] In one feasible implementation, a subtraction operation is performed on all echo points of the corrected real part echo one by one, that is, the value of each echo point is subtracted from the baseline value determined in step S402, and the results of the subtraction of each echo point are arranged in the order of echo time to form the baseline corrected echo.

[0110] For example, using the values ​​from step S402 above, the baseline value is 32. If the first echo point of the corrected real part echo is 612, then the first echo point of the baseline corrected echo is 612-32=580. If a certain echo point at the tail of the corrected real part echo is 30, then the corresponding echo point of the baseline corrected echo is 30-32=-2, that is, the tail echo is pulled back to near zero.

[0111] Please refer to the following: Figure 5 , Figure 5 A comparison diagram showing the distribution of the baseline of the corrected real part echo tail along the normalized depth position after phase correction, as provided by this invention. (See diagram for reference.) Figure 5As shown, the horizontal axis represents the normalized depth position, with a value range of 0 to 1, used to characterize the relative depth of each sampling point within the measurement well section; the vertical axis represents the tail baseline of the corrected real echo, which is the baseline value estimated by the preset proportion of the tail sampling points of the corrected real echo at each depth sampling point according to the baseline estimation method of steps S102 to S403 above.

[0112] Figure (a) corresponds to the Tarim ultra-deep interval, Figure (b) corresponds to the Junggar mid interval, and Figure (c) corresponds to the Junggar bot interval. Each sub-figure contains four types of information: the gray curve represents the point-wise tail baseline, indicating the unsmoothed tail baseline estimate at each depth sampling point; the blue solid line represents the moving median curve (window n=31), indicating the local trend of the baseline obtained by taking the median of the point-wise tail baseline using 31 sampling points as a moving window; the black dashed line represents the median baseline, indicating the median of all point-wise tail baselines within the entire interval; and the gray shaded band represents the baseline ±1σ band, indicating the overall dispersion of the tail baseline along the depth. The legend of each sub-figure indicates the median and standard deviation (σ) of the corresponding interval.

[0113] As shown in the figure, the tail baselines of all three well segments are not consistently zero, indicating that ultra-deep NMR echoes generally exhibit non-zero baseline shifts. Without baseline correction, this shift will be misinterpreted as a long transverse relaxation time component in subsequent transverse relaxation time spectrum inversion. Furthermore, the point-by-point tail baselines fluctuate dramatically around the moving median curve and exhibit numerous isolated spikes, suggesting that tail baseline estimation is susceptible to spike anomalies and random noise. A robust estimation method combining the median and median absolute deviation should be employed to suppress these effects.

[0114] Comparing the three well sections, we can see that: In Figure (a), the median baseline of the ultra-deep Tarim section is 2.77, and the standard deviation is 3.07. The baseline shift and dispersion are relatively small, and the moving median curve is generally flat, with only a slight rise around the normalized depth of about 0.85. In Figure (b), the median baseline of the central Junggar section is 2.97, and the standard deviation is 2.34. The baseline is the most stable overall. In Figure (c), the median baseline of the bottom Junggar section is 3.16, and the standard deviation is 6.36. The standard deviation is significantly larger than the previous two sections. Moreover, the moving median curve shows a significant rise in the normalized depth ranges of about 0.2 to 0.4 and 0.7 to 0.8, and the baseline shift is significantly enhanced. This indicates that the baseline drift is severe and changes drastically along the depth in the corresponding depth range of this section.

[0115] The figure shows that the tail baseline shift of ultra-deep NMR echoes is common and non-uniformly distributed along the depth. The baseline shift amplitude and dispersion vary significantly in different well sections. The robust baseline estimation method, which uses the median of the tail sampling points with a preset proportion and the absolute deviation of the median to remove peak anomalies, can adapt to well sections with different baseline shift degrees and provide a stable input near the tail regression zero value for subsequent denoising and lateral relaxation time spectrum inversion.

[0116] In summary, this application determines the retention interval by combining the tail median and the tail fluctuation scale and removes tail spike anomalies. Then, the median of the remaining sampling points is used as the baseline value and subtracted from the corrected real part echo as a whole. This makes the baseline estimation unaffected by tail spike anomalies, and the echo tail is pulled back to near zero. This weakens the false tail of long transverse relaxation time caused by baseline drift and avoids the tail offset being misinterpreted as a long transverse relaxation time component in the subsequent transverse relaxation time spectrum inversion, which would lead to spectral long tail and porosity anomalies.

[0117] Please see Figure 6 , Figure 6 This is the present invention. Figure 1 A schematic diagram of an embodiment of S104. (See attached diagram.) Figure 6 As shown, the method in this application embodiment may include the following steps S601-S604: In this embodiment of the application, the correlation coefficient is a statistical measure that characterizes the degree of linear correlation between the modal component and the preprocessed echo. The correlation coefficient ranges from [-1, 1]. The higher the correlation coefficient, the more consistent the overall attenuation trend of the modal component and the preprocessed echo, and the more likely the modal component is to contain effective NMR attenuation information. The specific form of the formula for calculating the correlation coefficient is described in step S104 above, and will not be repeated here.

[0118] In one feasible implementation, for each modal component obtained from variational mode decomposition, the mean of the modal component and the mean of the preprocessed echo are calculated. Then, the mean of the modal component is subtracted from the value of the modal component at each echo point to obtain a first deviation sequence. The mean of the preprocessed echo is subtracted from the value of the preprocessed echo at each echo point to obtain a second deviation sequence. The first deviation sequence and the second deviation sequence are multiplied point by point and summed to obtain a covariance term. Then, the arithmetic square root of the sum of squares of the values ​​of the first deviation sequence and the sum of squares of the values ​​of the second deviation sequence are calculated. The two arithmetic square roots are multiplied to obtain a normalization term. Finally, the covariance term is divided by the normalization term to obtain the correlation coefficient between the modal component and the preprocessed echo. The above calculation is repeated for all K modal components to obtain K correlation coefficients.

[0119] S602, obtain the ratio of the energy of each modal component to the total energy of all modal components, obtain the energy proportion of each modal component, and determine the center frequency of each modal component based on the spectrum of each modal component.

[0120] In this embodiment, energy refers to the sum of the squares of the values ​​of the modal components at all echo points, and energy is used to characterize the signal strength carried by the modal components. Total energy refers to the sum of the energies of all modal components. The spectrum refers to the function obtained after the modal component undergoes a Fourier transform, which characterizes the amplitude distribution of each frequency component in that modal component.

[0121] In one feasible implementation, for the k-th mode component, the energy Ek of the mode component is obtained by summing the squares of the values ​​of the mode component at each echo point; the total energy is obtained by summing the energies of all K mode components; and the energy of each mode component is divided by the total energy to obtain the energy proportion of each mode component. Simultaneously, a Fourier transform is performed on each mode component to obtain the corresponding spectrum. The center frequency of the mode component is obtained by using the squares of the moduli at each frequency point of the spectrum as weights. The specific form of the formula for calculating the center frequency is given in step S104 above and will not be repeated here.

[0122] For example, suppose the variational mode decomposition yields 6 modal components, with energies of 450, 300, 120, 60, 40, and 30 respectively. The total energy is 1000, and the energy percentages of each modal component are 0.45, 0.30, 0.12, 0.06, 0.04, and 0.03 respectively.

[0123] S603 defines modal components whose correlation coefficient is greater than or equal to the correlation threshold, whose energy percentage is greater than or equal to the energy percentage threshold, and whose center frequency is less than the center frequency threshold as effective modal components.

[0124] In this embodiment, the correlation threshold refers to the minimum correlation coefficient used to determine whether the correlation between the modal component and the preprocessed echo meets the requirements, and can be exemplarily set to 0.3. The energy percentage threshold refers to the minimum energy percentage used to determine whether the contribution of the modal component to the echo reconstruction meets the requirements, and can be exemplarily set to 0.05. The center frequency threshold refers to the upper frequency limit used to distinguish between the effective attenuation frequency band and the high-frequency noise frequency band, and the center frequency threshold can be set according to the lower limit frequency of the main noise frequency band.

[0125] In one feasible implementation, for each modal component, its correlation coefficient is compared with the correlation threshold, its energy percentage is compared with the energy percentage threshold, and its center frequency is compared with the center frequency threshold. The modal component that simultaneously satisfies the three conditions of correlation coefficient greater than or equal to the correlation threshold, energy percentage greater than or equal to the energy percentage threshold, and center frequency less than the center frequency threshold is determined as a valid modal component.

[0126] For example, using the values ​​from steps S601 and S602 above, the correlation coefficients of the six modal components are 0.92, 0.85, 0.61, 0.40, 0.12, and 0.05, respectively, and the energy proportions are 0.45, 0.30, 0.12, 0.06, 0.04, and 0.03, respectively. Furthermore, the center frequencies of the 5th and 6th modal components are higher than the center frequency threshold. With a correlation threshold of 0.3 and an energy proportion threshold of 0.05, the 1st to 4th modal components simultaneously meet all three conditions and are determined to be valid modal components, while the 5th and 6th modal components are identified as noise modal components. This screening method does not employ a single low-frequency retention principle but comprehensively considers three characteristics: correlation, energy proportion, and center frequency, which can reduce the risk of mistakenly deleting information with rapidly decaying short transverse relaxation times.

[0127] S604, if the number of effective modal components is lower than the preset number, modal components are added from the remaining modal components in descending order of energy percentage as effective modal components, and the determined effective modal components are superimposed to obtain the reconstructed echo.

[0128] In this embodiment, the preset quantity refers to the minimum number of effective modal components set to ensure that the main attenuation information is not excessively weakened. The preset quantity can be determined according to actual data, and can be set to 2 or 3 for example. The remaining modal components refer to the modal components that have not been selected as effective modal components among multiple modal components.

[0129] In one feasible implementation, the number of effective modal components selected in step S603 is first counted and compared with a preset number. If the number of effective modal components is greater than or equal to the preset number, the process proceeds directly to the superposition step. If the number of effective modal components is less than the preset number, the remaining modal components are sorted in descending order of energy percentage, and modal components are selected from the remaining sorted modal components from front to back to supplement the effective modal components until the number of effective modal components reaches the preset number.

[0130] For example, if the preset quantity is 3, and step S603 only selects 2 valid modal components, then the modal component with the largest energy proportion among the remaining modal components is added as a valid modal component. Finally, the values ​​of all valid modal components at each echo point are added point by point, that is, the value of the reconstructed echo at each echo point is the sum of the values ​​of each valid modal component at that echo point. After superposition, the reconstructed echo is obtained. The reconstructed echo retains the main modal components related to the attenuation trend of the preprocessed echo, and removes high-frequency, low-correlation, and low-energy noise modal components, thus achieving random noise suppression in the first stage.

[0131] In summary, this application selects effective modal components by comprehensively considering three characteristics: correlation coefficient, energy ratio, and center frequency. When the number of effective modal components is insufficient, they are supplemented in descending order of energy ratio. This ensures that the mid-to-high frequency modal components of the effective decay signal with short transverse relaxation time are not mistakenly deleted due to a single low-frequency retention strategy. At the same time, high-frequency, low-correlation, and low-energy noise modal components are removed, thereby achieving the first stage of random noise suppression while retaining the main decay information.

[0132] Please refer to the following: Figure 7 , Figure 7 This is a comparison image of VMD, wavelet, and VMD-wavelet fusion denoising echoes provided by this invention. Figure 7 As shown in the figure, the gray curve represents the original noisy signal, the blue curve represents the VMD denoising result, the green curve represents the wavelet denoising result, and the red curve represents the VMD-wavelet fusion denoising result. It can be seen from the figure that the original signal is severely affected by noise. The VMD method (blue) effectively removes high-frequency noise while preserving the overall trend of the signal, but significant amplitude fluctuations still exist at certain depths. The wavelet method (green) effectively suppresses spike noise but is prone to over-smoothing, leading to loss of signal details. The VMD-wavelet fusion method (red) combines the advantages of both methods, effectively eliminating high-frequency noise and spike interference while better preserving the amplitude characteristics and trends of the echo signal. Its denoised curve is the smoothest and closest to the true signal shape. The advantages of the fusion method are particularly evident in deep, high-noise environments, verifying the superiority of the VMD-wavelet fusion denoising method in nuclear magnetic resonance logging data processing.

[0133] Please see Figure 8 , Figure 8 This is the present invention. Figure 1 A schematic diagram of an embodiment of S104. (See attached diagram.) Figure 7 As shown, the method in this application embodiment may include the following steps S801-S803: S801 performs spike removal on the sampling points to obtain the removed sampling points.

[0134] In this embodiment, spike removal refers to the process of identifying and removing spike anomalies in the sampling points that deviate abnormally from the overall level based on the local statistical characteristics of the sampling points. Spike removal can adopt the local spike identification method in step S103 above, that is, it is based on the median and dispersion of the values ​​within the neighborhood window. The sampling points after spike removal refer to the sampling points retained after spike removal and used for noise estimation.

[0135] In one feasible implementation, taking a preset proportion of sampling points in the tail of the corrected imaginary part echo as an example, neighborhood window statistical discrimination is performed on each of these sampling points: for each sampling point, the median of the values ​​of the sampling points in its neighborhood window and the local robust standard deviation based on the absolute deviation of the median are calculated. If the absolute value of the difference between the value of the sampling point and the median in the neighborhood window is greater than a preset multiple of the local robust standard deviation, the sampling point is determined to be a peak outlier and removed; otherwise, the sampling point is retained. After all sampling points are discriminated, the retained sampling points are the sampling points after removal.

[0136] S802, a first noise estimate is determined based on the median of the absolute value of the difference between the median of each sampled point after elimination and the median of the sampled points after elimination, and a second noise estimate is determined based on the median of the absolute value of the difference between the sampled points of the tail of the corrected imaginary part echo at a preset proportion and the median of the median.

[0137] In this embodiment, the first noise estimate is a noise estimate determined based on the robust dispersion of sampling points after removing a predetermined proportion of the tail in the corrected imaginary echo. Since the corrected imaginary echo after phase correction should primarily consist of orthogonal residuals and noise, the first noise estimate directly reflects the noise intensity. The second noise estimate is a noise estimate determined based on the robust dispersion of sampling points in the corrected real echo or baseline corrected echo, based on a predetermined proportion of the tail. The second noise estimate is used to prevent noise from being underestimated when the corrected imaginary echo is abnormally low or there is phase distortion.

[0138] In one feasible implementation, for the first noise estimate, the median of the values ​​of each removed sampling point is calculated, and then the absolute value of the difference between each removed sampling point and the median is calculated. The median of this absolute value is then multiplied by a scaling factor of 1.4826 to obtain the first noise estimate. For the second noise estimate, sampling points at a predetermined proportion of the tail portion of the corrected real echo are taken, and the same robust dispersion calculation method is used, i.e., the median of the absolute value of the difference between this portion of sampling points and their median is calculated and multiplied by a scaling factor of 1.4826 to obtain the second noise estimate. The specific forms of the calculation formulas for the first and second noise estimates are given in step S105 above and will not be repeated here.

[0139] S803, the maximum value between the first noise estimate and the weighted second noise estimate is determined as the noise level.

[0140] In the embodiments of this application, weighting refers to the process of multiplying the second noise estimate by a weighting coefficient less than 1. The weighting coefficient ranges from 0.3 to 0.8, preferably 0.5. The purpose of weighting is to consider the possibility that a small amount of effective signal may remain in the tail of the real echo, which may cause the second noise estimate to be too high, and to reduce the second noise estimate.

[0141] In one feasible implementation, the second noise estimate is multiplied by a weighting factor to obtain a weighted second noise estimate. Then, the first noise estimate is compared with the weighted second noise estimate, and the maximum value of the two is determined as the noise level.

[0142] Furthermore, the noise level can be compared with a small amount (e.g., 10^-12) to prevent the noise level from being zero, and the maximum value can be taken to avoid the subsequent wavelet shrinkage threshold calculation failing due to a zero noise level in extreme cases. For example, using the values ​​from step S802 above, the first noise estimate is 12.0, the weighting coefficient is 0.5, then the weighted second noise estimate is 0.5 × 20.0 = 10.0, and the noise level is max(12.0, 10.0) = 12.0. This noise determination method mainly relies on the robust estimation of the imaginary part echo tail segment and uses the weighted real part tail segment estimate as a lower limit protection, which can avoid the noise underestimation caused by abnormally low imaginary part or phase distortion when relying solely on the imaginary part estimation.

[0143] Please refer to the following: Figure 9 , Figure 9 A comparison diagram of tail baseline correction and spike removal provided by this invention. (See attached image.) Figure 8As shown, the baseline correction and spike removal processes and effects of nuclear magnetic resonance echo signals at three different depths (8474.964 m, 8205.216 m, and 8984.933 m) are compared. From left to right, the five columns are: (a) the original phase-corrected real echo signal, showing obvious baseline drift and spike noise; (b) a local magnification of the echo tail (approximately 2400-2600 echoes), where the blue curve represents the signal before correction, and the baseline parameters b (30.0, 27.7, and 21.2 respectively) are estimated through fitting. The green dashed line represents the estimated baseline, and the overall signal shifts down to near zero after correction; (c) the baseline-corrected signal. The echo signal has eliminated the attenuation trend, but the spike noise is more prominent in the green curve; (d) Comparison before and after spike removal: the gray curve is the signal before processing, and the purple curve is the signal after processing. The number of spikes detected is marked (125, 86, and 85 respectively), which shows that the pulse interference has been effectively removed; (e) Residual distribution after spike removal: the red scatter points basically fluctuate randomly around the zero line, indicating that the processed signal no longer contains systematic bias. The residual has an approximate white noise characteristic, which verifies the effectiveness of the tail baseline correction and spike removal method.

[0144] In summary, this application determines the first noise estimate by removing spikes from the sampling points of the corrected imaginary part echo tail segment, and then determines the noise level by taking the maximum value of the weighted second noise estimate. This ensures that the noise estimate is not affected by isolated spikes, nor is it underestimated due to abnormally low corrected imaginary part echoes or phase distortions. This provides a reliable noise benchmark for the subsequent adaptive calculation of wavelet contraction threshold and robust signal-to-noise ratio calculation.

[0145] Please see Figure 10 , Figure 10 This is the present invention. Figure 1 A schematic diagram of an embodiment of S105. (See attached diagram.) Figure 9 As shown, the method in this application embodiment may include the following steps S1001-S1004: S1001 performs wavelet decomposition on the reconstructed echo to obtain low-frequency approximate components and high-frequency detail components of each layer.

[0146] In this embodiment, wavelet decomposition refers to the process of using wavelet basis functions to expand a signal at multiple scales, decomposing the signal into a low-frequency component representing the overall trend and a high-frequency component representing local changes at different scales. The low-frequency approximation component refers to the J-th level low-frequency component obtained by wavelet decomposition, which mainly reflects the overall attenuation trend of the reconstructed echo. The high-frequency detail components refer to the high-frequency components of each level obtained by wavelet decomposition, and these components include random noise and local oscillations.

[0147] In one feasible implementation, a wavelet basis function with good tight support and symmetry is selected, for example, the sym8 wavelet basis or the db4 wavelet basis can be selected, and the number of decomposition layers J is determined, which can be determined according to the total number of echo points of the reconstructed echo.

[0148] For example, 3 to 5 layers can be selected; for instance, when reconstructing 400 echo points, 4 layers can be used for decomposition. Then, wavelet decomposition is performed layer by layer on the reconstructed echo: the first layer decomposition yields the first layer low-frequency approximation component and the first layer high-frequency detail component; the first layer low-frequency approximation component is then further decomposed to obtain the second layer low-frequency approximation component and the second layer high-frequency detail component, and so on, until the Jth layer, finally obtaining the Jth layer low-frequency approximation component and the high-frequency detail components from the 1st to the Jth layers. The key point of wavelet thresholding is to shrink the high-frequency detail components of each layer while trying to maintain the low-frequency attenuation trend of the low-frequency approximation components.

[0149] S1002, a general threshold is determined based on the noise level and the total number of echo points of the reconstructed echo, and a candidate threshold is determined based on the ratio of the square of the noise level to the fluctuation scale of the high-frequency detail component of the current layer.

[0150] In this embodiment, the general threshold refers to a wavelet threshold that, asymptotically, can suppress Gaussian white noise, determined by the noise level and the total number of echo points. The specific form of the general threshold calculation formula is given in step S105 above and will not be repeated here. The candidate threshold is an adaptive threshold determined by combining the Bayesian thresholding concept and the ratio of the square of the noise level to the fluctuation scale of the high-frequency detail components in the current layer. The specific form of the candidate threshold calculation formula is given in step S105 above and will not be repeated here. The fluctuation scale refers to the dispersion of the high-frequency detail components in the current layer. The fluctuation scale can be taken as the standard deviation of the high-frequency detail coefficients in the current layer or a robust scale estimate based on the absolute deviation of the median.

[0151] In one feasible implementation, a general threshold is obtained by multiplying the noise level by the square root of the natural logarithm of the total number of echo points; at the same time, for the high-frequency detail components of the current layer, the fluctuation scale of the high-frequency detail coefficients of the current layer is calculated, and the square of the noise level is divided by the sum of the fluctuation scale and a small amount to prevent the denominator from being zero to obtain a candidate threshold.

[0152] For example, assuming a noise level of 15 and a total of 400 echo points in the reconstructed echo, the general threshold is 15 × √(2ln400) ≈ 51.9; if the fluctuation scale of the high-frequency detail component in the second layer is 30, then the candidate threshold is 15² / 30 = 7.5. The general threshold depends only on the noise level and data length, and may be too strong under local low signal-to-noise ratio conditions. The candidate threshold can adaptively decrease according to the actual fluctuation level of the high-frequency detail component in the current layer. The combination of the two can balance noise suppression and effective signal preservation.

[0153] S1003, determine the lower threshold, take the minimum value between the general threshold and the lower threshold and multiply it by the scale adjustment coefficient of the current layer to obtain the wavelet shrinkage threshold of the current layer.

[0154] In this embodiment, the lower limit coefficient refers to a proportional coefficient used to control the lower limit of the wavelet shrinkage threshold. The lower limit coefficient is used to prevent noise residue caused by excessively small candidate thresholds. The lower limit threshold is the threshold determined by the maximum value among the candidate thresholds, calculated by multiplying the lower limit coefficient and the general threshold. The scale adjustment coefficient is a proportional coefficient used to fine-tune the wavelet shrinkage threshold based on the noise distribution differences at different decomposition levels. The higher the frequency band of the high-frequency detail components in the current layer, the more the scale adjustment coefficient can be adjusted. For example, the scale adjustment coefficient can be a value between 0.8 and 1.2.

[0155] In one feasible implementation, the product of the lower limit coefficient and the general threshold is first calculated. This product is then compared with the candidate threshold, and the maximum value is taken to obtain the lower limit threshold. This step ensures that the wavelet shrinkage threshold is greater than or equal to a fixed proportion of the general threshold, avoiding noise residue caused by an excessively small threshold. Next, the general threshold is compared with the lower limit threshold, and the minimum value is taken. That is, the general threshold is used as the upper limit of the wavelet shrinkage threshold to avoid loss of effective attenuation information due to an excessively large threshold. Finally, the result of taking the minimum value is multiplied by the scaling coefficient of the current layer to obtain the wavelet shrinkage threshold of the current layer. The specific form of the final wavelet shrinkage threshold calculation formula is given in step S105 above and will not be repeated here.

[0156] For example, using the values ​​from step S902 above, the general threshold is 51.9, the candidate threshold is 7.5, and the lower limit coefficient is 0.3. The product of the lower limit coefficient and the general threshold is 15.57, and the lower limit threshold is max(7.5, 15.57) = 15.57. Then, the minimum value between the general threshold and the lower limit threshold is min(51.9, 15.57) = 15.57. If the scaling coefficient of the current layer is 1.0, then the wavelet shrinkage threshold of the current layer is 15.57.

[0157] S1004, based on the wavelet shrinkage threshold of each layer, performs soft threshold shrinkage on the high-frequency detail components of each layer, and reconstructs the low-frequency approximate components with the shrunken high-frequency detail components of each layer to obtain the target denoised echo.

[0158] In this embodiment, soft thresholding refers to a nonlinear process that sets the portion of the wavelet coefficients of high-frequency detail components whose absolute value does not exceed the wavelet shrinking threshold to zero, and compresses the portion whose absolute value exceeds the wavelet shrinking threshold toward zero. Soft thresholding can continuously compress small-amplitude high-frequency disturbances and avoid discontinuous oscillations caused by hard thresholding. The specific form of the calculation formula for soft thresholding is described in step S105 above, and will not be repeated here.

[0159] In one feasible implementation, for each high-frequency detail component, the wavelet coefficients of that layer are extracted one by one, and the absolute value of each wavelet coefficient is compared with the wavelet shrinkage threshold of that layer: if the absolute value of the wavelet coefficient does not exceed the wavelet shrinkage threshold, the wavelet coefficient is set to zero; if the absolute value of the wavelet coefficient exceeds the wavelet shrinkage threshold, the sign of the wavelet coefficient is retained, and its absolute value is subtracted from the wavelet shrinkage threshold as the value after shrinkage.

[0160] For example, if the wavelet shrinkage threshold of the second layer is 15.57, then the wavelet coefficients with a value of 10 in the high-frequency detail components of the second layer are set to zero, the wavelet coefficients with a value of 40 are shrunk to 40-15.57=24.43, and the wavelet coefficients with a value of -40 are shrunk to -24.43. After the high-frequency detail components of all layers are shrunk, the low-frequency approximate component of the Jth layer and the shrunk high-frequency detail components of each layer are reconstructed layer by layer according to the inverse process of wavelet decomposition, thus obtaining the target denoised echo. The target denoised echo maintains the main attenuation trend and suppresses high-frequency random noise and local small-amplitude oscillations. It can be directly input into the transverse relaxation time spectrum inversion process to recover the transverse relaxation time spectrum and calculate interpretation parameters such as porosity.

[0161] In summary, this application determines the wavelet shrinkage threshold for each layer by combining a general threshold with candidate thresholds and setting a lower limit threshold and a scaling coefficient. This allows the wavelet shrinkage threshold to adaptively adjust with the noise level and the actual fluctuation level of the high-frequency detail components in each layer. Furthermore, by using soft threshold shrinkage, the high-frequency detail components are continuously compressed. This avoids noise residue caused by an excessively small threshold, while also preventing the loss of effective attenuation information due to an excessively large threshold. Thus, the main attenuation envelope of the NMR echo is maintained while suppressing residual high-frequency disturbances.

[0162] In one feasible implementation, after performing local spike identification and replacement on the baseline correction echo to obtain the preprocessed echo, the following is further performed: If the mean of the early echo points of the preprocessed echo is negative, then the sign of the preprocessed echo is flipped. The robust signal-to-noise ratio (SNR) is determined by the ratio of the truncated mean of the absolute values ​​of the early echo points in the preprocessed echo after sign inversion to the median absolute deviation of the sampling points at a predetermined proportion in the tail of the corrected imaginary part echo. The robust SNR is used to characterize the signal-to-noise level of the original complex NMR echo.

[0163] In this embodiment, early echo points refer to the first few echo points in the preprocessed echo that are located at the beginning of the echo time segment, have high signal amplitude, and good signal-to-noise ratio. The number of early echo points can be selected between 20 and 80, preferably 50. Sign flipping refers to multiplying the values ​​of each echo point in the preprocessed echo by -1. Sign flipping does not change the signal-to-noise characteristics of the echo, but only changes the overall polarity of the echo. Truncated mean refers to the mean obtained by sorting a set of data by size, removing the data at both ends with a preset proportion, and then averaging the remaining data. The truncated mean is not sensitive to isolated extreme values. Median absolute deviation refers to the median of the data, which is the absolute value of the difference between each data point in a set of data and the median of that set of data. Robust signal-to-noise ratio is the signal-to-noise ratio determined by the ratio of the truncated mean of the absolute values ​​of the early echo points of the preprocessed echo after sign flipping to the robust noise estimate of the sampling points in the tail segment of the corrected imaginary part echo. The numerator of the robust signal-to-noise ratio represents the effective NMR signal intensity, and the denominator represents the noise intensity. Signal-to-noise level refers to the relative strength of the effective NMR signal in the original complex NMR echo compared to the noise level.

[0164] In one feasible implementation, the mean value of the early echo points of the preprocessed echo is first calculated. If the mean value is negative, it indicates that the principal polarity of the effective NMR signal after phase correction is negative. In this case, the values ​​of all echo points of the preprocessed echo are multiplied by -1 and the sign is flipped to keep the early effective NMR signal positive, so as to satisfy the physical characteristic of the overall non-negative attenuation of the spin echo train (Carr-Purcell-Meiboom-Gill, CPMG) echo. Here, the spin echo train refers to a series of transversely magnetized attenuated echoes excited and collected by the pulse sequence in NMR logging. If the mean value is not negative, the preprocessed echo is kept unchanged.

[0165] To calculate the robust signal-to-noise ratio, for example, the specific form of the formula for calculating the robust signal-to-noise ratio is as follows:

[0166] in, For robust signal-to-noise ratio; For example, to calculate the truncated mean, you can remove 10% of the data from both ends after sorting before averaging. The sequence of early echo points formed by the preprocessed echoes after symbol inversion, and the truncated mean of their absolute values ​​are used to characterize the effective NMR signal intensity. To correct the sequence of sampling points with a preset ratio in the tail of the imaginary part echo; 1.4826 is the absolute deviation of the median; 1.4826 is the scaling coefficient for converting the absolute deviation of the median into a scaling estimate consistent with the standard deviation under Gaussian noise conditions, and the denominator as a whole is used to characterize the noise intensity.

[0167] For example, assuming the early echo points are the first 50 echo points, after sorting their absolute values ​​and removing the last 5 data points from each end, the average of the remaining 45 data points is 180. The median absolute deviation of the 100 sampling points in the corrected imaginary part echo tail is 8.1. Therefore, the robust signal-to-noise ratio (SNR) is 180 / (1.4826 × 8.1) = 180 / 12.0 = 15. Since the numerator of the robust SNR uses the truncated mean and the denominator uses the median absolute deviation of the corrected imaginary part echo tail, both are insensitive to isolated spikes. Therefore, the robust SNR is not abnormally inflated by isolated spikes like the conventional SNR, and can more realistically reflect the SNR level of ultra-deep low SNR NMR data. The calculated robust SNR can be used to evaluate the degree of improvement in echo quality before and after denoising. Combined with the transverse relaxation time spectrum inversion results, the denoising effect can be comprehensively evaluated from the echo domain, spectral domain, and interpretation parameter domain.

[0168] Please refer to the following: Figure 11 , Figure 11 This is a comparison of NMR echo, T2 spectrum, and porosity before and after VMD-wavelet fusion denoising, as provided in an embodiment of this application. Figure 11 It can be seen that there are obvious random fluctuations and local anomalous interferences in the original real and imaginary echoes. After phase rotation, baseline correction, and peak replacement, the phase residue, non-zero baseline shift, and local sudden anomalies in the preprocessed echoes are weakened. After further processing with the VMD-wavelet fusion method, the random oscillation amplitude of the target denoised echo is significantly reduced, while the main attenuation envelope is well preserved. Compared with the Schlumberger T2 spectrum, the peak shape of the target denoised T2 spectrum is more continuous and stable, the local pseudo-peaks and spectral tail anomalies caused by noise are suppressed, and the distribution position and relative amplitude of the main relaxation components do not change significantly. The porosity comparison results show that the variation trend and numerical range of the NMR porosity curve after denoising are more consistent with the core porosity measurement points. This indicates that the VMD-wavelet fusion denoising method provided in this application can retain effective NMR attenuation information while suppressing echo noise, thereby improving the stability and reliability of T2 spectrum inversion and porosity calculation results.

[0169] DEPTH represents the downhole depth, GR curve represents the natural gamma curve, ECHO_AMP_R represents the real component of the original complex NMR echo, ECHO_AMP_X represents the imaginary component of the original complex NMR echo, ECHO_PREPROCESS represents the preprocessed echo obtained after phase rotation correction, tail baseline correction, and local spike identification and replacement, ECHO_DENOISED represents the target denoised echo obtained after variational mode decomposition, effective mode screening, and wavelet threshold shrinkage, T2_DIST represents the T2 spectrum of the denoised inversion given by Schlumberger, T2_VMD+WAVELET represents the T2 spectrum of the VMD-wavelet fusion denoised echo obtained by inversion based on the target denoised echo, and POR_CORE represents the core porosity, POR_BRD represents the target denoised NMR porosity, and POR_DIST represents the porosity obtained by Schlumberger.

[0170] In summary, this application reverses the sign of the preprocessed echo when the mean of the early echo points is negative, so that the calculation of the robust signal-to-noise ratio satisfies the physical characteristic of the overall non-negative attenuation of the NMR echo. Furthermore, it uses the truncated mean of the absolute values ​​of the early echo points to characterize the effective NMR signal intensity and the median absolute deviation of the corrected imaginary part of the echo tail segment to characterize the noise intensity. This ensures that neither the numerator nor the denominator of the robust signal-to-noise ratio is affected by isolated spikes, thereby avoiding the problem of the conventional signal-to-noise ratio being abnormally raised by isolated spikes under ultra-deep low signal-to-noise ratio conditions. This provides a reliable basis for the quantitative evaluation of the denoising effect.

[0171] The following will combine Figure 12 This paper provides a detailed description of the nuclear magnetic resonance echo denoising device provided in the embodiments of this application. It should be noted that... Figure 12 The nuclear magnetic resonance echo noise reduction device in the present application is used to perform the following functions. Figures 1-8 The methods shown in the embodiments are for illustrative purposes only, illustrating the parts relevant to the embodiments of this application. For specific technical details not disclosed, please refer to this application. Figures 1-8 In the embodiment shown, the nuclear magnetic resonance echo denoising device 10 may include a data acquisition unit 11, a data correction unit 12, a data decomposition unit 13, a data reconstruction unit 14, and a data generation unit 15, as detailed below: Data acquisition unit 11 is used to acquire the original complex NMR echo, perform rotation correction on the real part echo in the original complex NMR echo to obtain the corrected real part echo, and the rotation angle of the rotation correction is the negative angle of the principal phase angle. The principal phase angle is the arctangent angle between the accumulated imaginary part value and the accumulated real part value of the previous preset number of echo points of the original complex NMR echo. Data correction unit 12 is used to determine the baseline value based on the sampling points at a preset ratio in the tail of the corrected real echo, and to subtract the baseline value from the corrected real echo to obtain the baseline corrected echo; The data decomposition unit 13 is used to identify and replace local spikes in the baseline correction echo to obtain a preprocessed echo, and to perform variational mode decomposition on the preprocessed echo to obtain multiple mode components. The data reconstruction unit 14 is used to select effective modal components from multiple modal components based on the correlation, energy ratio and center frequency between each modal component and the preprocessed echo, and to superimpose the selected effective modal components to obtain the reconstructed echo. The data generation unit 15 is used to determine the noise level based on the sampling points, determine the wavelet shrinkage threshold based on the noise level, and perform wavelet threshold shrinkage and reconstruction on the reconstructed echo based on the wavelet shrinkage threshold to obtain the target denoised echo.

[0172] The nuclear magnetic resonance echo denoising device 10 provided in the above embodiments can realize the technical solutions described in the above nuclear magnetic resonance echo denoising method embodiments. The specific implementation principles of each module or unit can be found in the corresponding content in the above nuclear magnetic resonance echo denoising method embodiments, and will not be repeated here.

[0173] like Figure 13 As shown, the present invention also provides an electronic device 20. This electronic device 10 includes a processor 22, a memory 22, a display 23, and a data acquisition unit 24. Figure 13 Only some of the components of the electronic device 20 are shown, but it should be understood that it is not required to implement all of the components shown, and more or fewer components may be implemented instead.

[0174] In some embodiments, processor 22 may be a central processing unit (CPU), microprocessor or other data processing chip, used to run program code stored in memory 22 or process data, such as the nuclear magnetic resonance echo denoising method of the present invention.

[0175] In some embodiments, processor 22 may be a single server or a group of servers. The server group may be centralized or distributed. In some embodiments, processor 22 may be local or remote. In some embodiments, processor 22 may be implemented on a cloud platform. In one embodiment, the cloud platform may include a private cloud, public cloud, hybrid cloud, community cloud, distributed cloud, intranet, multi-cloud, etc., or any combination thereof.

[0176] In some embodiments, memory 22 may be an internal storage unit of electronic device 20, such as a hard disk or memory of electronic device 20. In other embodiments, memory 22 may also be an external storage device of electronic device 20, such as a plug-in hard disk, smart media card (SMC), secure digital (SD) card, flash card, etc. equipped on electronic device 20.

[0177] Furthermore, the memory 22 may include both internal storage units of the electronic device 20 and external storage devices. The memory 22 is used to store application software and various types of data installed on the electronic device 20.

[0178] Display 23 may be an LED display, a liquid crystal display, a touch-screen liquid crystal display, or an OLED (Organic Light-Emitting Diode) touchscreen in some embodiments. Display 23 is used to display information from electronic device 20 and to display a visual user interface. Data acquisition unit 24 may be a nuclear magnetic resonance logging tool in some embodiments for acquiring raw complex nuclear magnetic resonance echoes. Components 22-24 of electronic device 20 communicate with each other via a system bus.

[0179] In one embodiment, when the processor 22 executes the nuclear magnetic resonance echo denoising program in the memory 22, the following steps can be implemented: The real part echo in the original complex NMR echo is rotated and corrected to obtain the corrected real part echo. The rotation angle of the rotation correction is the negative angle of the principal phase angle. The principal phase angle is the arctangent angle between the accumulated imaginary part value and the accumulated real part value of the previous preset number of echo points in the original complex NMR echo. The baseline value is determined based on the sampling points at a preset ratio in the tail of the corrected real part echo, and the baseline value is subtracted from the corrected real part echo to obtain the baseline corrected echo. Local spikes are identified and replaced in the baseline-corrected echo to obtain a preprocessed echo. Variational mode decomposition is then performed on the preprocessed echo to obtain multiple modal components. Based on the correlation, energy ratio, and center frequency between each modal component and the preprocessed echo, effective modal components are selected from multiple modal components, and the selected effective modal components are superimposed to obtain the reconstructed echo; The noise level is determined based on the sampling points, the wavelet shrinkage threshold is determined based on the noise level, and the reconstructed echo is subjected to wavelet threshold shrinkage and reconstruction based on the wavelet shrinkage threshold to obtain the target denoised echo.

[0180] It should be understood that when the processor 22 executes the nuclear magnetic resonance echo denoising program in the memory 22, in addition to the functions mentioned above, it can also perform other functions, as can be found in the description of the corresponding method embodiments above.

[0181] Furthermore, this embodiment of the invention does not specifically limit the type of electronic device 20 mentioned. Electronic device 20 can be a tablet computer, a personal digital assistant (PDA), a wearable device, a laptop computer, or other portable electronic devices. Exemplary embodiments of portable electronic devices include, but are not limited to, portable electronic devices running iOS, Android, Microsoft, or other operating systems. The aforementioned portable electronic device can also be other portable electronic devices, such as a laptop computer with a touch-sensitive surface (e.g., a touch panel). It should also be understood that in some other embodiments of the invention, electronic device 20 may not be a portable electronic device, but rather a desktop computer with a touch-sensitive surface (e.g., a touch panel).

[0182] Accordingly, this application also provides a computer-readable storage medium for storing computer-readable programs or instructions. When the programs or instructions are executed by a processor, they can implement the steps or functions of the nuclear magnetic resonance echo denoising methods provided in the above-described method embodiments.

[0183] Those skilled in the art will understand that all or part of the processes of the methods described in the above embodiments can be implemented by a computer program instructing related hardware (such as a processor, controller, etc.), and the computer program can be stored in a computer-readable storage medium. The computer-readable storage medium may be a disk, optical disk, read-only memory, or random access memory, etc.

[0184] The above provides a detailed description of the nuclear magnetic resonance echo denoising method, apparatus, and electronic equipment provided by the present invention. Specific examples have been used to illustrate the principles and implementation methods of the present invention. The description of the above embodiments is only for the purpose of helping to understand the method and core ideas of the present invention. At the same time, for those skilled in the art, there will be changes in the specific implementation methods and application scope based on the ideas of the present invention. Therefore, the content of this specification should not be construed as a limitation of the present invention.

Claims

1. A method for nuclear magnetic resonance echo noise reduction, characterized in that, include: The original complex NMR echo is obtained, and the real part echo in the original complex NMR echo is rotated and corrected to obtain the corrected real part echo. The rotation angle of the rotation correction is the negative angle of the principal phase angle. The principal phase angle is the arctangent angle between the accumulated imaginary part value and the accumulated real part value of the previous preset number of echo points of the original complex NMR echo. The baseline value is determined based on the sampling points at a preset ratio in the tail of the real part echo of the correction, and the baseline value is subtracted from the real part echo of the correction to obtain the baseline correction echo. The baseline-corrected echo is subjected to local spike identification and replacement to obtain a preprocessed echo. The preprocessed echo is then subjected to variational mode decomposition to obtain multiple mode components. Based on the correlation, energy ratio, and center frequency between each modal component and the preprocessed echo, effective modal components are selected from the multiple modal components, and the selected effective modal components are superimposed to obtain the reconstructed echo; The noise level is determined based on the sampling points, the wavelet shrinkage threshold is determined based on the noise level, and the reconstructed echo is subjected to wavelet threshold shrinkage and reconstruction based on the wavelet shrinkage threshold to obtain the target denoised echo.

2. The method according to claim 1, characterized in that, The step of rotating and correcting the real part echo in the original complex NMR echo to obtain the corrected real part echo includes: Complex echo points are constructed based on the real and imaginary echoes in the original complex nuclear magnetic resonance echoes, and each complex echo point is rotated around the origin by the negative angle of the principal phase angle to obtain the rotated echo points; Extract the real component of each of the rotated echo points to obtain the corrected real component echo; After performing rotational correction on the real part echo in the original complex NMR echo to obtain the corrected real part echo, the method further includes: Extract the imaginary component of each of the rotated echo points to obtain the corrected imaginary echo.

3. The method according to claim 1, characterized in that, The step of determining the baseline value based on sampling points at a preset proportion in the tail of the corrected real echo, and subtracting the baseline value from the corrected real echo to obtain the baseline corrected echo, includes: Based on the median of the sampling points, the tail median is determined, and the tail fluctuation scale is determined based on the median of the absolute values ​​of the differences between each sampling point and the tail median. Sampling points whose absolute value of the difference from the median of the tail is greater than a preset multiple of the tail fluctuation scale are removed, and the median of the remaining sampling points is determined as the baseline value. The baseline-corrected echo is obtained by subtracting the baseline value from each echo point of the corrected real part echo.

4. The method according to claim 1, characterized in that, The step of identifying and replacing local spikes in the baseline-corrected echo to obtain a preprocessed echo includes: Based on the median and dispersion of the echo point values ​​within the neighborhood window of each echo point in the baseline-corrected echo, peak anomalies in the baseline-corrected echo are identified, and the values ​​of the identified peak anomalies are replaced with the median of the echo point values ​​within the corresponding neighborhood window to obtain the preprocessed echo.

5. The method according to claim 1, characterized in that, The process involves selecting effective modal components from the multiple modal components based on the correlation, energy ratio, and center frequency between each modal component and the preprocessed echo, and then superimposing the selected effective modal components to obtain the reconstructed echo, including: Obtain the correlation coefficient between each modal component and the preprocessed echo, the correlation coefficient being used to characterize the consistency between the attenuation trends of the modal components and the preprocessed echo; The energy ratio of each modal component to the total energy of all modal components is obtained to obtain the energy ratio of each modal component, and the center frequency of each modal component is determined based on the spectrum of each modal component. Modal components whose correlation coefficient is greater than or equal to the correlation threshold, whose energy percentage is greater than or equal to the energy percentage threshold, and whose center frequency is less than the center frequency threshold are determined as effective modal components. If the number of effective modal components is less than the preset number, modal components are added from the remaining modal components in descending order of energy percentage as effective modal components. The determined effective modal components are superimposed to obtain the reconstructed echo.

6. The method according to claim 1, characterized in that, Determining the noise level based on the sampling points includes: Peak removal is performed on the sampling points to obtain the removed sampling points; A first noise estimate is determined based on the median of the absolute value of the difference between the median of each of the removed sampling points and the median of the removed sampling points, and a second noise estimate is determined based on the median of the absolute value of the difference between the median of the sampling points of the preset proportion in the tail of the corrected imaginary part echo and the median of the median. The maximum value between the first noise estimate and the weighted second noise estimate is determined as the noise level.

7. The method according to claim 1, characterized in that, The step of determining a wavelet shrinkage threshold based on the noise level, and performing wavelet threshold shrinkage and reconstruction on the reconstructed echo based on the wavelet shrinkage threshold to obtain the target denoised echo, includes: The reconstructed echo is decomposed by wavelet to obtain low-frequency approximate components and high-frequency detail components of each layer; A general threshold is determined based on the noise level and the total number of echo points of the reconstructed echo, and a candidate threshold is determined based on the ratio of the square of the noise level to the fluctuation scale of the high-frequency detail component of the current layer. Determine the lower threshold, take the minimum value between the general threshold and the lower threshold, and multiply it by the scaling coefficient of the current layer to obtain the wavelet shrinkage threshold of the current layer; Based on the wavelet shrinkage threshold of each layer, the high-frequency detail components of each layer are subjected to soft threshold shrinkage, and the low-frequency approximation components are reconstructed with the shrunken high-frequency detail components of each layer to obtain the target denoised echo.

8. The method according to claim 1, characterized in that, After performing local spike identification and replacement on the baseline-corrected echo to obtain the preprocessed echo, the method further includes: If the mean value of the early echo points of the preprocessed echo is negative, then the sign of the preprocessed echo is flipped. The robust signal-to-noise ratio (SNR) is determined by the ratio of the truncated mean of the absolute values ​​of the early echo points of the preprocessed echo after sign flipping to the median absolute deviation of the sampling points at the predetermined proportion in the corrected imaginary part echo. The robust SNR is used to characterize the signal-to-noise level of the original complex NMR echo.

9. A nuclear magnetic resonance echo noise reduction device, characterized in that, The device includes: The data acquisition unit is used to acquire the original complex NMR echo, perform rotation correction on the real part echo in the original complex NMR echo to obtain the corrected real part echo, wherein the rotation angle of the rotation correction is the negative angle of the principal phase angle, and the principal phase angle is the arctangent angle between the accumulated imaginary part value and the accumulated real part value of the previous preset number of echo points of the original complex NMR echo. The data correction unit is used to determine the baseline value based on the sampling points at a preset ratio in the tail of the corrected real part echo, and to subtract the baseline value from the corrected real part echo to obtain the baseline corrected echo; The data decomposition unit is used to identify and replace local spikes in the baseline-corrected echo to obtain a preprocessed echo, and to perform variational mode decomposition on the preprocessed echo to obtain multiple mode components. The data reconstruction unit is used to select effective modal components from the multiple modal components based on the correlation, energy ratio and center frequency between each modal component and the preprocessed echo, and to superimpose the selected effective modal components to obtain the reconstructed echo; The data generation unit is used to determine the noise level based on the sampling points, determine the wavelet shrinkage threshold according to the noise level, and perform wavelet threshold shrinkage and reconstruction on the reconstructed echo based on the wavelet shrinkage threshold to obtain the target denoised echo.

10. An electronic device, characterized in that, The electronic device includes a data acquisition unit, a memory, and a processor, wherein, The data acquisition device is used to acquire raw complex nuclear magnetic resonance echoes; The memory is used to store executable program code; The processor is configured to perform rotational correction on the real part echo in the original complex NMR echo to obtain a corrected real part echo; determine a baseline value based on a preset proportion of sampling points at the tail of the corrected real part echo, and subtract the baseline value from the corrected real part echo to obtain a baseline corrected echo; perform local peak identification and replacement on the baseline corrected echo to obtain a preprocessed echo; perform variational mode decomposition on the preprocessed echo to obtain multiple modal components; select effective modal components from the multiple modal components based on the correlation, energy ratio, and center frequency between each modal component and the preprocessed echo, and superimpose the selected effective modal components to obtain a reconstructed echo; determine the noise level based on the sampling points, determine a wavelet shrinkage threshold based on the noise level, and perform wavelet threshold shrinkage and reconstruction on the reconstructed echo based on the wavelet shrinkage threshold to obtain a target denoised echo.

Images