Robust electromagnetic inverse scattering inversion method based on maximum correlation entropy

Abstract

The invention discloses a robust electromagnetic inverse scattering inversion method based on maximum correlation entropy, and belongs to the technical field of electromagnetic inverse scattering and signal processing. According to the method, a maximum correlation entropy criterion in an information theory is introduced into a contrast source inversion framework, and a novel cost function which is composed of a data fidelity item and a state fidelity item and is based on MCC is constructed; equivalently converting the minimization problem of the non-convex cost function into a series of iterative weighted least square sub-problems by adopting a positive semi-definite optimization technology; in the iteration process, the influence of outliers in the data is dynamically suppressed through a self-adaptive weighting mechanism depending on a previous iteration residual error, and an alternating optimization strategy is adopted to efficiently solve a contrast source and a contrast function. According to the method, the strong outlier in the data can be self-adaptively suppressed, the reconstruction precision and robustness are far better than those of a traditional method in a non-Gaussian noise environment, and meanwhile, the method also shows superior or equivalent performance in a noiseless or Gaussian noise environment, and has extremely high practical value.

Description

Technical Field

[0001] This invention relates to the field of electromagnetic imaging technology, and in particular to a signal processing and image reconstruction method for electromagnetic inverse scattering problems (EISP). Background Technology

[0002] Electromagnetic inverse scattering (EIS) technology aims to deduce the spatial distribution of electromagnetic parameters (such as dielectric constant and conductivity) within a target by measuring the electromagnetic field data scattered by the target. It has broad application prospects in fields such as medical imaging, nondestructive testing, geophysical exploration, and radar detection. However, EIS is a typical nonlinear and severely ill-posed problem, making its inversion process extremely sensitive to measurement noise and data errors, posing significant challenges to the stability and accuracy of the reconstructed results.

[0003] To address the aforementioned nonlinear problems, the Contrast Source Inversion (CSI) method was proposed. This method introduces a "contrast source" as an intermediate variable, decomposing the original nonlinear inversion problem into a series of alternating iterative linear subproblems. This effectively alleviates the local minima problem in strong scattering scenarios. Compared to traditional linearization methods such as the Born approximation or Rytov approximation, it exhibits significant advantages in convergence and reconstruction performance, and has become a benchmark algorithm in the field of electromagnetic inverse scattering.

[0004] However, traditional CSI methods and their derivative algorithms (such as multiplicative regularized CSI and MR-CSI) still face significant technical bottlenecks in practical applications. The objective function of these methods is typically constructed based on the minimum mean squared error (MSE) criterion, which aims to minimize the L2 norm between the model's predicted data and the actual measured data. From an information theory perspective, the MSE criterion only utilizes the second-order statistic (variance) of the error signal, and its optimality strictly depends on the assumption that the noise follows a Gaussian distribution.

[0005] In real-world measurement environments, scattered field data is often contaminated by various complex non-Gaussian noises, such as multiplicative speckle noise generated by coherent imaging systems, heavy-tailed distribution noise caused by multipath propagation or sensor malfunctions, and strong outliers like impulse noise. When such non-Gaussian noise exists, especially outliers with large amplitudes, the squared term of the MSE criterion over-amplifies the penalty weight for outlier errors, causing the optimization process to be dominated by a few "bad points." This leads to a sharp decline in the performance of the inversion algorithm, resulting in numerous artifacts in the reconstructed image, or even complete failure to converge, severely impacting the robustness and reliability of the imaging.

[0006] Although existing technologies attempt to improve the noise resistance of CSI methods by introducing techniques such as total variational regularization, sparse constraints, or Bayesian frameworks, most of these improvements remain within the data fidelity term framework of MSE and fail to fundamentally solve the sensitivity of the MSE criterion to non-Gaussian noise, especially impulse noise.

[0007] Therefore, there is an urgent need for a new electromagnetic inverse scattering inversion method that can still achieve stable and high-precision imaging in complex non-Gaussian noise environments, especially in the presence of impulse noise and strong outliers, in order to break through the limitations of the traditional MSE framework. Summary of the Invention

[0008] The purpose of this invention is to overcome the above-mentioned defects of the prior art and provide a signal processing method for electromagnetic inverse scattering inversion. It aims to solve the technical problems of poor robustness, low reconstruction accuracy, and susceptibility to abnormal data interference leading to result distortion when the traditional contrast source inversion (CSI) method based on the minimum mean square error (MSE) criterion is used to process scattering field data containing non-Gaussian noise, especially strong outliers such as impulse noise and heavy tail noise.

[0009] To achieve the above-mentioned objectives, this invention introduces the Maximum Correlation Criterion (MCC) from information theory into the CSI framework, replacing the traditional Minimum Mean Square Error (MSE) criterion. By constructing a new cost function and adopting an iterative optimization strategy, robust and high-precision reconstruction of the electromagnetic parameters of the target medium is achieved.

[0010] The technical solution provided by this invention is a robust electromagnetic inverse scattering inversion method based on maximum correlation entropy, comprising the following steps:

[0011] Step 1: Establish a discretized operator model for the electromagnetic scattering problem;

[0012] The target region D to be inverted is discretized into N grid cells; the state equations and data equations describing the scattering phenomenon are discretized to obtain the following matrix-form operator model:

[0013] Step 2: Construct a cost function based on the maximum correlation entropy criterion;

[0014] The electromagnetic inverse scattering inversion problem is reformulated as an optimization problem with a cost function. J MCC It consists of a data fidelity term and a state fidelity term, both of which are measured based on the maximum relevance entropy criterion. The cost function... for:

[0015] (4)

[0016] Where i is the index of the incident wave, N i This represents the total number of incident antennas;

[0017] For the current comparison source Calculated scattering field Compared with the actual measured scattering field The difference between them, i.e. ;

[0018] For the current comparison source With current contrast The difference between the contrast source calculated from the incident field and the incident field, i.e. ;

[0019] and These are the Gaussian kernel function bandwidth parameters for the data and state terms, respectively, used to control sensitivity to errors. and To balance the regularization parameters for data items and state items;

[0020] Step 3: Establish the equivalent transformation and iterative weighted least squares framework of the cost function;

[0021] Directly minimize the non-convex, non-quadratic cost function in step 2. The computational cost is enormous. This invention employs a half-quadratic (HQ) optimization technique, introducing auxiliary variables to transform it into an iterative weighted least squares (IRLS) problem. In the (n+1)th iteration, the optimization objective is transformed into minimizing a weighted quadratic cost function J. (n+1) :

[0022] (5)

[0023] in, , Represents the weighted norm. and These represent the corresponding diagonal weight matrices;

[0024] This weighting mechanism can automatically identify outliers in the residuals: when the residuals are small (credible data), the weights approach 1; when the residuals are large (outliers), the weights rapidly decay to 0, thereby suppressing the influence of outliers in subsequent optimizations.

[0025] Step 4: Iterative optimization solution using alternating steps;

[0026] In the (n+1)th iteration, an alternating optimization strategy is adopted, fixing one variable and solving for the optimal solution of the other variable:

[0027] The contrast function is fixed as the result of the previous iteration. At this point, the cost function It is about A quadratic function; by letting it about gradient Solve for zero The gradient is represented as:

[0028] (8)

[0029] The gradients of the data item and the state item are respectively:

[0030] (9)

[0031] (10)

[0032] Wherein, the superscript H indicates the conjugate transpose, ¯ indicates conjugate, and ⊙ indicates the Hadamard product. , These represent the Gaussian kernel function bandwidth parameters for the data item and the state item in the nth iteration, respectively. Represents the identity matrix. , Let these represent the residuals of the data item and the state item in the nth iteration, respectively.

[0033] The weighted conjugate gradient method is used for iterative solution to obtain the updated comparison source. ;

[0034] Update contrast function The comparison source is fixed as the newly updated one. At this point, minimize the cost function The problem relates only to the state terms and is decomposed into solutions that are independent for each discrete element r; for each element r, The contrast function representing unit r:

[0035] (11)

[0036] Among them, the main field Calculations show that It is a small regularization constant to prevent the denominator from being zero. Let represent the vectorized representation of the incident field at the (n+1)th update. This represents the vectorized representation of the total field at the (n+1)th update;

[0037] Step 5: Convergence judgment and result output;

[0038] Repeat the alternating iterative process in step 4 until the preset convergence condition is met. After the iteration terminates, output the final contrast function. This is the result of reconstructing the spatial distribution of the electromagnetic parameters of the target.

[0039] Furthermore, the specific method of step 1 is as follows:

[0040] Discrete state equations describe the total field within the scattering region D as being composed of the superposition of the incident field and the field generated by the contrast source, as shown in Equation 1.

[0041] (1)

[0042] in, This represents the vectorized representation of the total field across N discrete grid cells. This represents the vectorized representation of the incident field over N discrete grid cells. Represents the discrete-domain Green's operator;

[0043] Discrete data equation: describes the scattered field measured at the receiving antenna outside the scattering region as being caused by a contrast source inside region D. The radiation produced is shown in Equation 2;

[0044] (2)

[0045] in, Let represent the vectorized representation of the scattered field data measured by Nr receiving antennas. Represents the discrete scattering Green's operator.

[0046] This is represented as a vectorized representation of the contrast source to be solved over N discrete grid cells. With contrast Satisfying Relationship:

[0047] (3)

[0048] in, Represents the vector The resulting diagonal matrix.

[0049] Furthermore, in step 3, the diagonal weight matrix and It is dynamically calculated based on the residual of the nth iteration, and its diagonal elements are determined by the Gaussian kernel function:

[0050] (6)

[0051] (7)

[0052] in, The diagonal weight matrix is ​​repositioned as follows: Element.

[0053] Furthermore, the convergence condition for step 5 is: between two iterations The relative change is less than a given threshold, or the preset maximum number of iterations is reached.

[0054] Compared with the prior art, the present invention has the following significant advantages:

[0055] 1) Extremely high robustness: This invention replaces the minimum mean square error (MSE) criterion in the traditional CSI method with the maximum correlation entropy criterion. The cost function of MCC has a bounded, reducible influence function, which can adaptively suppress the influence of large outliers (such as impulse noise) in the data. Experiments show that under impulse noise contamination, the structural similarity (SSIM) and peak signal-to-noise ratio (PSNR) of the reconstruction results of this invention can be improved by several times compared with the traditional CSI method, almost completely eliminating the destructive effect of strong outliers on image quality.

[0056] 2) Higher reconstruction accuracy: Even in noise-free or Gaussian noise environments, MCC, by maximizing correlation entropy, can capture higher-order statistical information in the data that is ignored by MSE, thereby more accurately inverting the intrinsic physical properties of the target. Experiments show that, under noise-free conditions, the average SSIM and PSNR of this invention are still significantly improved compared to traditional methods, enabling the reconstruction of clearer and more detailed images.

[0057] 3) Excellent convergence performance: The positive semi-definite (HQ) optimization framework adopted in this invention theoretically guarantees the monotonic convergence of the original cost function. By transforming the complex non-convex problem into a series of easily solvable weighted least squares problems, a stable and efficient solution process is achieved.

[0058] 4) Wide applicability: This invention not only provides a high-performance robust imaging algorithm for electromagnetic inverse scattering problems, but also provides an effective information theory framework, offering a new theoretical perspective and technical approach for solving large-scale inverse problems plagued by non-Gaussian noise in other fields. Attached Figure Description

[0059] Figure 1 The diagrams show the irregular geometric shapes measured in the experiment (MNIST handwritten font, Austrian graphics, overlapping circles, overlapping rectangles, and overlapping rings, respectively).

[0060] Figure 2 This is a schematic diagram of the experimental results of backscattering imaging of handwritten characters in MNIST.

[0061] Figure 3 The iterative curves of MNIST handwritten characters under various algorithms (for scenarios with no noise, noise 1, noise 2, and noise 3, respectively).

[0062] Figure 4 This is a schematic diagram of the experimental results for the backscattering imaging of the Austrian pattern.

[0063] Figure 5 The iteration curves of the Austrian graph under various algorithms are shown (for scenarios with no noise, noise 1, noise 2, and noise 3, respectively).

[0064] Figure 6 This is a schematic diagram of the experimental results of backscattering imaging of overlapping circular images.

[0065] Figure 7 The iterative curves of the overlapping circle image under various algorithms are shown (for scenarios with no noise, noise 1, noise 2, and noise 3, respectively).

[0066] Figure 8 This is a schematic diagram of the experimental results of inverse scattering imaging of overlapping rectangular images.

[0067] Figure 9 The iterative curves of the overlapping rectangular image under various algorithms are shown (for scenarios with no noise, noise 1, noise 2, and noise 3, respectively).

[0068] Figure 10 This is a schematic diagram of the experimental results of backscattering imaging of overlapping ring patterns.

[0069] Figure 11 The iterative curves of the overlapping ring pattern under various algorithms are shown (for scenarios with no noise, noise 1, noise 2, and noise 3, respectively). Detailed Implementation

[0070] This embodiment uses a multi-station static radar system under two-dimensional transverse magnetic waves as an example to demonstrate the process of accurate target imaging by the method of the present invention in a strong impulse noise environment.

[0071] I. Experimental Environment and Parameter Settings

[0072] 1. Physical Model: A two-dimensional scattering scene is set, with free space as the background. The target to be imaged is located within a 2m × 2m square detection area D. The operating frequency of the electromagnetic waves is set to 300MHz.

[0073] 2. Detection System: A multi-station static radar configuration is employed. N radars are uniformly positioned around the detection area D. t = 16 transmitting antennas and N r= 32 receiving antennas. During operation, 16 transmitting antennas sequentially transmit TM-polarized cylindrical waves to illuminate the detection area, and all 32 receiving antennas synchronously receive the scattered field data.

[0074] 3. Imaging Target: A grayscale image of a digit extracted from the MNIST handwritten digit dataset is used as the true contrast distribution of the target to be imaged. The grayscale values ​​of the image are mapped to the contrast values ​​of the medium; for example, the pixel value range [0, 1] is mapped to the contrast range [0, 1.5].

[0075] 4. Discretization: The 2m×2m detection area D is discretized into a uniform grid of 64×64, with a total of N = 4096 discrete units.

[0076] 5. Operator Generation: Based on the discretized grid and antenna layout described above, the discrete-domain Green's operator G is calculated using the method of moments. D and discrete scattering Green's operator G S These two matrices are fixed and are calculated and stored all at once before the inversion begins. G D The dimensions are 4096×4096, G S The dimensions are (32×4096).

[0077] 6. Simulation data generation and noise injection:

[0078] First, using the forward scattering model Calculate the scattered field data produced by the target under ideal noise-free conditions. .

[0079] Then, in order to simulate harsh measurement environments, [the following was done]... Inject strong noise. Specifically:

[0080] 1) To simulate the granular artifacts unique to coherent imaging systems (such as SAR and ultrasound), a multiplicative noise model was employed. Damage field Defined as

[0081] (12)

[0082] in This represents zero-mean complex Gaussian white noise. It is a scalar quantity for controlling noise intensity.

[0083] 2) To explain the multipath propagation effect and other non-Gaussian interference sources that produce significant outliers, we introduce additive noise from the Student t-distribution.

[0084] (13)

[0085] in For complex noise, the real and imaginary parts are independently sampled with degrees of freedom of . The t-distribution. Smaller. This will result in a heavy-tailed distribution, and It then acts as a scaling factor.

[0086] 3) To represent sensor malfunctions or data transmission errors, we simulate impulse noise. Specifically, we randomly select 3% of the data points and forcibly replace their amplitudes with an abnormally large value (e.g., 10 times the original amplitude of the point), thereby generating simulated measurement data with strong outliers.

[0087] The three noise scenarios mentioned above are recorded as Noise1, Noise2, and Noise3, respectively.

[0088] Detailed experimental results are attached. Figures 1 to 11 The quantitative data of the experimental results are shown in Tables 1 and 2. The results show that, in different datasets, the method proposed in this invention has higher output accuracy and better stability compared with the classical CSI method and the multiplicative regularized contrastive source inversion method (MR-CSI).

[0089] The entire pseudocode process of this invention is shown in the table below.

[0090]

[0091] The innovation of this invention lies in:

[0092] 1. This invention is the first to introduce the "Maximum Correlation Entropy Criterion" (MCC) from information theory into the field of electromagnetic inverse scattering, replacing the "Minimum Mean Square Error" (MSE) criterion relied upon by traditional methods. MSE is extremely sensitive to noise, while MCC, using the kernel function method, can fundamentally "passivate" the interference of non-Gaussian noise such as strong impulses and heavy tails. This is the theoretical foundation for the robustness of this invention.

[0093] 2. This invention designs a highly efficient iterative adaptive weighting algorithm. In each calculation, it automatically identifies "bad data" contaminated by noise and significantly reduces its weight, while amplifying the influence of "good data." This ensures that even when the data is severely contaminated, the inversion process converges in the correct direction.

[0094] 3. This invention solves the industry problem of traditional methods failing completely in non-Gaussian noise environments, achieving stable and high-precision robust imaging. More importantly, while achieving strong robustness, it does not sacrifice performance in noise-free or conventional Gaussian noise environments, achieving a "full-scene" performance improvement, and has extremely high practical value and wide applicability.

[0095] The above description of specific embodiments details the implementation process and application scenarios of the maximum correlation entropy comparison source inversion method, enabling those skilled in the art to implement the maximum correlation entropy comparison source inversion method according to the present invention. The above embodiments are only used to illustrate the technical solution of the present invention and do not limit the scope of protection of the present invention. Those skilled in the art can make various modifications and variations based on the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

[0096] Table 1. Comprehensive evaluation results of SSIM and PSNR under different noise conditions

[0097]

[0098] Table 2 Comparison of SSIM and PSNR improvement effects under different noise conditions

[0099]

Claims

1. A robust electromagnetic inverse scattering inversion method based on maximum correlation entropy, comprising the following steps: Step 1: Establish a discretized operator model for the electromagnetic scattering problem; The target region D to be inverted is discretized into N grid cells; the state equations and data equations describing the scattering phenomenon are discretized to obtain the following matrix-form operator model: Step 2: Construct a cost function based on the maximum correlation entropy criterion; Cost function for: (4) Where i is the index of the incident wave, N i This represents the total number of incident antennas; For the current comparison source Calculated scattering field Compared with the actual measured scattering field The difference between them, i.e. ; For the current comparison source With current contrast The difference between the contrast source calculated from the incident field and the incident field, i.e. ; and These are the Gaussian kernel function bandwidth parameters for the data and state terms, respectively, used to control sensitivity to errors. and To balance the regularization parameters for data items and state items; Step 3: Establish the equivalent transformation and iterative weighted least squares framework of the cost function; In the (n+1)th iteration, the optimization objective is transformed into minimizing a weighted quadratic cost function J. (n+1) : (5) in, , Represents the weighted norm. and These represent the corresponding diagonal weight matrices; Step 4: Iterative optimization solution using alternating steps; In the (n+1)th iteration, an alternating optimization strategy is adopted, fixing one variable and solving for the optimal solution of the other variable: The contrast function is fixed as the result of the previous iteration. At this point, the cost function It is about A quadratic function; by letting it about gradient Solve for zero The gradient is represented as: (8) The gradients of the data item and the state item are respectively: (9) (10) Wherein, the superscript H indicates the conjugate transpose, ¯ indicates conjugate, and ⊙ indicates the Hadamard product. , These represent the Gaussian kernel function bandwidth parameters for the data item and the state item in the nth iteration, respectively. Represents the identity matrix. , Let these represent the residuals of the data item and the state item in the nth iteration, respectively. The weighted conjugate gradient method is used for iterative solution to obtain the updated comparison source. ; Update contrast function The comparison source is fixed as the newly updated one. At this point, minimize the cost function The problem relates only to the state terms and is decomposed into solutions that are independent for each discrete element r; for each element r, The contrast function representing unit r: (11) Among them, the main field Calculations show that It is a small regularization constant to prevent the denominator from being zero. Let represent the vectorized representation of the incident field at the (n+1)th update. This represents the vectorized representation of the total field at the (n+1)th update; Step 5: Convergence judgment and result output; Repeat the alternating iterative process in step 4 until the preset convergence condition is met. After the iteration terminates, output the final contrast function. This is the result of reconstructing the spatial distribution of the electromagnetic parameters of the target.

2. The robust electromagnetic inverse scattering inversion method based on maximum correlation entropy as described in claim 1, characterized in that, The specific method for step 1 is as follows: Discrete state equations describe the total field within the scattering region D as being composed of the superposition of the incident field and the field generated by the contrast source, as shown in Equation 1. (1) in, This represents the vectorized representation of the total field across N discrete grid cells. This represents the vectorized representation of the incident field over N discrete grid cells. Represents the discrete-domain Green's operator; Discrete data equation: describes the scattered field measured at the receiving antenna outside the scattering region as being caused by a contrast source inside region D. The radiation produced is shown in Equation 2; (2) in, Let represent the vectorized representation of the scattered field data measured by Nr receiving antennas. Represents the discrete scattering Green's operator. This is represented as a vectorized representation of the contrast source to be solved over N discrete grid cells. With contrast Satisfying Relationship: (3) in, Represents the vector The resulting diagonal matrix.

3. The robust electromagnetic inverse scattering inversion method based on maximum correlation entropy as described in claim 1, characterized in that, The diagonal weight matrix in step 3 and It is dynamically calculated based on the residual of the nth iteration, and its diagonal elements are determined by the Gaussian kernel function: (6) (7) in, The diagonal weight matrix is ​​repositioned as follows: Element.

4. The robust electromagnetic inverse scattering inversion method based on maximum correlation entropy as described in claim 1, characterized in that, The convergence condition for step 5 is: between two iterations The relative change is less than a given threshold, or the preset maximum number of iterations is reached.

Images