Intelligent Reconstruction Algorithm of Electromagnetic Map Based on Propagation Model

I. INTRODUCTION

THE sixth-generation (6G) mobile communication technology will bring about a fully connected, intelligent world, promote the rapid development of smart cities, Internet of Vehicles, smart homes, and human-computer interaction [1], and also put forward higher requirements for the communication quality and perception capabilities of mobile communication networks [2], [3].At the same time, due to the rapid development of 6G, the number of users is exploding [4], and its demand for increasingly tight spectrum resources is becoming more and more urgent [5], [6]. The improvement of communication quality, perception ability, and spectral efficiency of wireless communication networks is an important factor to stimulate the development of integrated sensing and communication (ISAC) technology. ISAC technology devices are jointly designed, optimized, and scheduled by sharing the same hardware platform, common spectrum, joint signal processing strategies, and a unified control framework [7].

ISAC technology can realize efficient parallel processing of sensing and communication tasks by sharing resources in the same frequency band, which can not only significantly improve spectrum efficiency and energy efficiency [8], but also improve communication and sensing performance [9], [10]. ISAC technology is widely used in many fields and will promote technological innovation in many fields. However, while integrating communication and perception technologies, ISAC technology also makes the use of spectrum more complex and diverse. How to achieve spectrum analysis and management has become one of the key issues [11]. Therefore, there is an urgent need for a feasible solution to guide the analysis and management of electromagnetic spectrum by accurately characterizing it.

The management of electromagnetic spectrum resources includes the planning, allocation and safety monitoring of electromagnetic spectrum resources, and the electromagnetic map (EM) can quantitatively describe the electromagnetic environment from multiple dimensions such as the working frequency and spatial distribution of electromagnetic signals, and can realize the fine description of the comprehensive information of electromagnetic spectrum resources in combination with geographic information [12], which can not only provide support for the planning and allocation of spectrum resources, but also monitor the electromagnetic environment in the target area to provide strong support for the analysis and management of the electromagnetic spectrum [13]. The commonly used methods for reconstructing EM mainly include data-driven, model driven, and data model dual driven methods [14]. The Kriging algorithm (KGA) [15] and the inverse distance weighting (IDW) [16] algorithm can interpolate the limited data collected and obtain complete data that can be used for EM reconstruction. The application scenarios of both methods are relatively wide, but the KGA interpolation method is prone to ignoring drastically changing data, and the reconstruction accuracy of the IDW interpolation method is poor. The method of matrix completion can also be used to reconstruct EM, but its precision is poor, and due to its large computational complexity, it takes longer to complete for larger areas [17]. The sparse Bayesian learning [18] and the Gaussian kernel function [19] have also been applied to solve the reconstruction of EM, but their application scenarios are limited. With the development of artificial intelligence, some scholars have also applied traditional machine learning methods and deep learning to the reconstruction of EM, and have achieved certain results [20]–[25]. However, such methods require a large amount of prior information as training data, and their reconstruction effect depends on the quality of the training data. The use of residual kriging and Bayesian learning can obtain more accurate results, but when the position of the radiation source is unknown, its training data cannot be obtained, and its application in scenarios is limited [26]. Tensor recovery [27], [28] is also applied to the reconstruction of EM, but its usage scenarios are limited. Graph neural networks are used to construct EM [29], but they only consider the influence of some points on the estimated points when constructing the graph structure. Both the model driven reconstruction EM method [30] and the data model driven reconstruction EM method [31] require more detailed prior information to be implemented. Although the electromagnetic map reconstruction quality obtained by these two methods is high, their applicability is limited due to the lack of prior information and difficulties in deploying sampling nodes [32].

The variogram function can be used to quantitatively describe the spatial correlation statistics of related random fields and processes [33]. In the case of limited prior information and random deployment of sampling nodes, the variogram function can obtain the correlation between known and unknown points in space based on prior information, and then use spatial interpolation to calculate the electromagnetic data of unknown points [34]. However, the correlation between the variation function based on empirical models and electromagnetic data is poor, leading to a poor fitting effect and affects the reconstruction accuracy of EM [35].

In response to the shortcomings of existing algorithms, this article starts from the variogram function, uses propagation models to model the variogram function, and completes spatial interpolation reconstruction of electromagnetic maps. The main contributions of this article are as follows.

· By using unbiased estimation criteria to infer the interrelationships between electromagnetic data, and combining the electromagnetic data of known perception nodes to estimate and predict the electromagnetic data of unknown points, a complete and reconfigurable EM data can be obtained.

· By combining the laws of radio wave propagation and using propagation models to model the variation function, a variation function that is more in line with the correlation between electromagnetic data is obtained, laying the foundation for more accurate extraction of the correlation between known and unknown points.

· Introducing weight factors to improve the particle swarm optimization algorithm to solve the undetermined coefficients of the mutation function, obtaining an accurate functional relationship between the lag distance and the mutation function, and accurately extracting the correlation between known and unknown points, thereby improving the reconstruction accuracy of EM.

The chapter arrangement of this article is as follows. In Section II, the paper provides the system model under the ISAC architecture, Section III introduces the algorithm flow of this article, simulation experiments and result analysis in Section IV, and conclusions are given in Section V.

II. PROBLEM STATEMENT

Under the ISAC architecture, the Internet of things has led to a significant increase in the number of terminal devices for users, and the distribution of these terminal devices within the region is random. These terminal devices have the ability to detect and communicate with base stations. In order to achieve cross zone switching and power control of terminal devices, the terminal devices will report the Signal quality to the base station, which will upload the received data to the data center. The data center will summarize the data and construct an electromagnetic map of the target area to obtain the electromagnetic distribution of the target area, thereby providing strong support for the analysis and management of the electromagnetic spectrum of the entire target area. The system model is shown in Fig. 1.

Given that electromagnetic signals are essentially inherently stationary stochastic processes with isotropy, covariance functions can be used to spatially model and interpolate the distribution of electromagnetic signals. This article takes into account the correlation between electromagnetic data in the target area and the small number of terminal equipment. Therefore, the method of weighted sum of electromagnetic data from known points can be used to obtain the electromagnetic data of unknown points. Based on the variogram function, the problem of spatial interpolation can be transformed into a nonlinear constraint minimization problem, which can be solved using the Lagrange multiplier method. Based on the electromagnetic data obtained from solving unknown points, combined with the electromagnetic data from sampling nodes, an REM of the target area can be drawn, and then the electromagnetic spectrum can be analyzed and managed.

Due to the low matching between the empirical model of the Variogram function and electromagnetic data, the Variogram function is improved based on the principle of radio wave propagation. Considering that there are many outliers and outliers in the sample data calculated by the terminal equipment for fitting the Variogram function, and their distribution is unknown, they are first grouped and then removed using a box plot. The obtained data can be used for fitting the Variogram function model.

Considering the global search ability and fast convergence speed of the particle swarm optimization algorithm, it is applied to solve the undetermined coefficients of the expression of the Variogram. At the same time, combining the difference in sample size and the spatial correlation of electromagnetic data, weight factors are introduced to balance fitness, in order to obtain the optimal solution for the undetermined coefficients of the Variogram function.

III. ALGORITHM DESCRIPTION

This article takes the electromagnetic data of nodes in the electromagnetic map as the research object, which includes the position information of nodes and the reference signal receiving power (RSRP). After dividing the target area into grids, electromagnetic data is collected using a limited number of randomly deployed terminal devices within the target area. The grid with terminal device deployment is defined as a known grid, and its electromagnetic data is the electromagnetic data collected by the terminal device. The grid without terminal device deployment is defined as the grid to be estimated. By using the electromagnetic data of the known grid to calculate the electromagnetic data of the grid to be estimated, and drawing it using isomagnetic lines, the EM of the target area can be reconstructed. Therefore, this article proposes an intelligent reconstruction algorithm of electromagnetic map based on propagation model (PSO-SIM), and the algorithm process is shown in Fig. 2.

Firstly, the target area is divided into grids, and then the electromagnetic data collected by randomly deployed terminal devices in the target area is used to calculate the distance and semi variation values between known grids. Then, the semi variation values are grouped, outliers are removed, and the mean is calculated to obtain a two-dimensional discrete matrix of the lag distance and semi variation values. Based on the principle of radio wave propagation, the variogram function is modeled and the undetermined coefficients of the variogram function are fitted and solved using an improved particle swarm optimization algorithm. Finally, the electromagnetic data of the grid to be estimated is calculated by combining the variogram and spatial interpolation algorithm, and the EM of the target area is drawn using isomagnetic lines.

A. Spatial Interpolation Method

After dividing the target area into grids, set the total number of grids to M and the randomly distributed number of terminal device to N, that is, the number of known grids is N, and the number of grid to be estimated is N - M. The electromagnetic data collected by the terminal device is [TeX:] $$Z=\left[Z\left(x_1\right), \cdots, Z\left(x_i\right), \cdots, Z\left(x_N\right)\right],$$ where [TeX:] $$x_i$$ represents the position information of the terminal device, [TeX:] $$x_i=\left(X_i, Y_i\right)$$ and [TeX:] $$Z\left(x_i\right)$$ represents the RSRP of the terminal device in [TeX:] $$x_i.$$ Since the electromagnetic signal is a stationary stochastic process, the expected and variance of the RSRP value [TeX:] $$Z\left(x_i\right)$$ for [TeX:] $$x_i$$ in the target region is as follows:

So, it can be understood that the [TeX:] $$Z\left(x_i\right)$$ at any point in the target area is composed of the average RSRP in that area and the random deviation [TeX:] $$R\left(x_i\right)$$ of that point, the calculation formula is:

Based on the above question, the RSRP of the estimated grid can be linearly weighted and summed based on the known RSRP of the grid.

where [TeX:] $$Z\left(x_i\right)$$ represents the RSRP of the known grid, [TeX:] $$\hat{Z}\left(x_0\right)$$ represents the RSRP of the grid to be estimated and [TeX:] $$\lambda_i$$ represents the coefficient matrix. Therefore, by determining the coefficient matrix [TeX:] $$\lambda_i$$, the RSRP of any grid to be estimated can be calculated using the RSRP of the known grid. To determine the coefficient matrix [TeX:] $$\lambda_i$$, two standard conditions need to be met, namely unbiased and minimum variance.

Substituting (4) into (5) can be calculated to obtain:

Similarly, substituting (4) into (6) can calculate:

where for convenience, [TeX:] $$C_{i,j}$$ is used here to represents performing covariance operation on RSRP in [TeX:] $$X_i \text { and } X_j,$$ that is: [TeX:] $$C_{i j}=\operatorname{Cov}\left[Z\left(x_i\right), Z\left(x_j\right)\right] .$$

The formula for solving the semi-variation value between any two known grids in the target area is:

where [TeX:] $$\gamma_{i j}$$ represents the semi variation value between [TeX:] $$x_i$$ and [TeX:] $$x_j.$$ Substituting (3) into (9) can be calculated to obtain:

Finally, substituting (10) into (8) can be calculated to obtain:

Replace [TeX:] $$2 \sum_{i=1}^N \lambda_i \gamma_{i 0}-\sum_{i=1}^N \sum_{j=1}^N \lambda_i \lambda_j \gamma_{i j}-\gamma_{00}$$ with J here. At this point, the spatial interpolation problem can be transformed into a nonlinear constrained minimization problem.

To solve the above problem, we use the Lagrangian multiplier method to solve (12) and construct the Lagrangian function as follows:

where Φ represents the Lagrange multiplier. Solving the coefficient matrix [TeX:] $$\lambda_i$$ and the Lagrange multiplier Φ that minimizes [TeX:] $$J^{\prime}$$ can satisfy the problem of minimizing the above constraints. If we set the partial derivatives of [TeX:] $$J^{\prime}$$ with respect to [TeX:] $$\lambda_i$$ and Φ to zero, we can obtain:

Write (14) as a matrix in the form of:

Solve (15) to obtain the optimal coefficient matrix [TeX:] $$\lambda_i$$ and its corresponding Φ. By obtaining the coefficient matrix [TeX:] $$\lambda_i$$, the RSRP of any unknown point can be solved based on the RSRP of the sampling node combined with (4). Visualize all the data again to obtain the reconstructed EM. From the perspective of spatial interpolation, the key to solving the RSRP of any grid to be estimated is to obtain an accurate coefficient matrix [TeX:] $$\lambda_i$$. According to (15), it can be seen that solving the coefficient matrix [TeX:] $$\lambda_i$$ requires solving the semivariation value [TeX:] $$\gamma_{i 0}.$$ Therefore, obtaining an accurate variogram function to calculate the semi-variation value [TeX:] $$\gamma_{i 0}$$ is crucial. The method of solving the variogram function will be provided in the following two sections.

B. Variogram Function Modeling

It is generally believed that the closer the distance between any two grids in space, the closer their attribute values are, and the farther the distance, the greater the difference. We can assume that there is a certain functional relationship between the distance between any two grids in the target area and the semi variation value. Therefore, it is only necessary to model the relationship between the distance between known grids in the target area and the semi variation value, obtain the variogram function, and then calculate the semi-variation value [TeX:] $$\gamma_{i 0}$$ based on the distance between the estimated grid and the known grid.

In order to model the functional relationship between the semi variation value and distance, it is necessary to first calculate the semi variation value and distance between each pair of known grids based on the electromagnetic data of the known grids, and obtain the variogram cloud. Then, the variogram cloud is grouped, outliers removed, estimated, and processed to obtain a two-dimensional discrete point set of distance and the semi variation value. Finally, model the variogram function to obtain the functional relationship between the semi variation value and the lag distance. Drawing a scatter plot corresponding to distance and semi variation values can obtain the variogram cloud and the processed variogram cloud as shown in Fig. 3.

1) Variogram cloud: The variogram cloud needs to calculate the distance and semi variation value between each pair of known grids based on the position information in the electromagnetic data of the known grids and RSRP. The calculation formula is:

According to the calculation result of (16), the variogram cloud can be obtained by one-to-one correspondence between the semi variation value and the distance. When the distances are equal, the distribution of semi variation values is more concentrated, but there are still outliers present. This is also a key factor that affects the effect of the variogram function. Similarly, it is necessary to consider the correlation between electromagnetic data. The propagation distance of electromagnetic signals in space is limited, and when the distance is too far, calculating the half variance value will actually affect the accuracy. In summary, it is necessary to process the variogram cloud.

2) Sample data processing: Taking into account the correlation between electromagnetic data and the sample quality of data in the variogram cloud. This section first sets the lag distance and then groups the variogram cloud based on the lag distance. Then remove the outliers and calculate the mean of the mutation function cloud after removing the outliers. Finally, obtain a discrete point set of outliers and the semivariation value and draw its scatter plot.

Set the minimum distance in the variogram cloud as the basic lag distance and group based on multiples of the basic lag distance. Considering that the c data in each group after grouping does not follow a normal distribution and has a large amount of data, a simpler and more effective method for removing outliers from non normal distribution data, interquartile range (IQR), is used. IQR is calculated based on the quartiles of the data, specifically referring to the difference between the third and fourth quartiles ([TeX:] $$Q_3$$, which is the 75% percentile) and the first quartiles ([TeX:] $$Q_1$$, which is the 25% percentile). The calculation formula is the following.

Then set a reasonable interval as: [TeX:] $$\left[Q_1-I Q R \times K, Q_3+I Q R \times K\right],$$ where K is a positive number that controls the length of the interval, usually set to K = 1.5. After removing outliers from each group, calculate the average distance of each group as the lag distance [TeX:] $$h_m.$$ The formula for calculating the semi variation value corresponding to the lag distance is:

where [TeX:] $$N\left(h_m\right)$$ represents the sum of the number of point pairs in each group and [TeX:] $$\gamma_{m i}$$ represents the semi variation values involved in the calculation in the group with a lag distance of [TeX:] $$h_m .$$ From this, a two-dimensional discrete point set [TeX:] $$\left\{\left(h_1, \gamma_1\right),\left(h_2, \gamma_2\right), \cdots,\left(h_n, \gamma_n\right)\right\}$$ of the hysteresis distance and the semi-variation value can be obtained. The processed variogram cloud and the variogram cloud are shown in Fig. 3.

3) Variogram function modeling: By selecting an appropriate model to fit a set of discrete points, the function expression of lag distance and semi variation value can be obtained to calculate the semi variation values of the known grid and the grid to be estimated. In general, models used for fitting include Gaussian models, exponential models, power function models, and spherical models. However, the commonly used models are empirical models, which have poor fitting effects and low matching with electromagnetic data. Therefore, it is necessary to seek more effective solutions.

Set the signal model within the target area as

where [TeX:] $$P_T$$ represents the signal transmission power, w represents the receiving matrix, H represents the communication channel, x represents the sending symbol, n represents the additive noise with [TeX:] $$n \sim \mathcal{N}(0,1).$$ Common models used to describe signal propagation in the real environment include free-space models, OKumura-Hata models, etc. To simplify the calculation, the free-space model was chosen to be used, i.e., the path loss formula is:

where f represents the frequency, d represents the distance.

Since the object of this study is RSRP, the RSRP of a known grid with position information [TeX:] $$x_i$$ can be calculated according to the signal model given in (19).

where [TeX:] $$\| \cdot| |$$ represents a two-norm calculation and [TeX:] $$\sigma_i^2$$ represents the noise power. Bringing (21) into (9) calculates the semi variation value between a known grid whose position information is [TeX:] $$x_i \text { and } x_j.$$ The result is as follows:

Since channel H obeys a complex Gaussian distribution, that is [TeX:] $$H \sim \mathcal{C N}(0, \beta),$$ so:

where L represents the path loss. Since the frequency of the signal propagating in the target region is known, the path loss of the signal in the target region can be considered to be proportional to the logarithm of the distance according to the path loss model. At the same time, additive noise obeys the Gaussian distribution, and then the result of the subtraction of the noise power of known grid at different locations in the target region [TeX:] $$\beta_i-\beta_j$$ obeys the chi-square distribution. Replace [TeX:] $$\sigma_i^2-\sigma_j^2 \text { with } \sigma_{i j}^2$$ here. And based on the above results, (22) can be expressed as:

Combined with the traditional exponential function model, and according to (24), we can model the variogram function as:

where a, b, c represents the undetermined coefficient. The exact expression of the variogram function can be obtained by solving the undetermined coefficients. The next section will specifically introduce the use of improved particle swarm optimization algorithms to solve the optimal solution of undetermined coefficients in function expressions.

C. Improving Particle Swarm Optimization Algorithm

Solving the pending coefficients of the variogram function usually involves an optimization fitting process, which can be achieved by a variety of methods, including the least squares method, the gradient descent method, the particle swarm optimization algorithm (PSO), etc. The undetermined coefficient of the solution variogram function needs to be determined by minimizing or maximizing an objective function, and the particle swarm optimization algorithm can be selected to solve the above problems due to the advantages of strong global search ability, strong adaptability, and fast convergence speed [36], [37]. In the process of solving the undetermined coefficients of the variogram function, the PSO usually takes the undetermined coefficients as the position of the particles, and takes the evaluation index based on the function value as the fitness function, and then initializes the position and velocity of the particles and defines the position and velocity of the particles.

In the original framework of the PSO algorithm, the common assumption was that all data points contributed the same to the function fitting and optimization problem, which meant that each data point contributed the same to the same problem regardless of the amount of information it carries and its location in the data set. However, when processing electromagnetic data, due to the random distribution of sampling nodes and the grouping of electromagnetic data according to the lag distance, the sample size of the grouped data is not uniformly distributed among the groups. This uneven distribution is not only a difference in the spatial distribution of the data, but also implies that there are differences in the characteristics of the data reflected by different groups of data. If this difference is still ignored, the intrinsic connection and correlation between electromagnetic data cannot be used more effectively, which will affect the final effect of the algorithm. In order to overcome the above problems, this section makes targeted improvements to the PSO algorithm, and comprehensively considers the correlation between the difference in the sample size of each group and the electromagnetic data to optimize the performance of the algorithm. Specifically, when constructing the objective function, the weight factor is innovatively introduced, and the weight factor is used to make the contribution of different groups of data to match their importance in the optimization process, so that the optimization process is more in line with the actual characteristics of electromagnetic data and the inertia weight adjustment strategy of linear decreasing is adopted to obtain better search performance.

Firstly, it is necessary to define a fitness function. Based on the variogram model in (25), we can define the fitness function as:

where [TeX:] $$F(i)$$ represents the fitness function, [TeX:] $$\gamma_i(h)$$ represents the actual calculated semi variation value, [TeX:] $$\hat{\gamma}_i(h)$$ represents the estimated semi variation value based on the model. Here, we introduce the [TeX:] $$\psi$$ weight factor a to improve the fitness function.

where [TeX:] $$\bar{h}_m$$ represents the average lag distance, the smaller the lag distance, the higher the fitness. N represents the total number of samples, n represents the number of samples in each group, the more samples a group has, the weaker the weight factor on fitness. So modify its fitness function to:

After obtaining the fitness function, the next step is to initialize the particle swarm and define the position and velocity of the particles as follows:

where [TeX:] $$D_i$$ represents the set of positions of a particle swarm, [TeX:] $$a_i, b_i, c_i$$ represents the solution of the undetermined coefficients in the current function expression, [TeX:] $$V_i$$ represents the velocity set of a particle swarm, [TeX:] $$V_{a i}, V_{b i}, V_{c i}$$ represents the magnitude of the changes in the solutions of three undetermined coefficients, with positive and negative indicating the direction of the changes in the undetermined coefficients. By changing the speed, the position of particles can be updated. In order to seek the global optimal solution, the particle swarm optimization algorithm continuously updates the velocity and position of particles through iteration, and records the global and individual optimal solutions. The formula for updating particle swarm velocity and position is:

where d represents the dimension of the particle swarm, In the search for the optimal coefficient of variation function, its dimension is 3, [TeX:] $$\omega$$ represents the inertia weight, k represents the number of iterations, [TeX:] $$c_1, c_2$$ represents the individual learning factors and group learning factors, [TeX:] $$r_1, r_2$$ represents the random number between zero and one. [TeX:] $$x_{i d}^{k+1}$$ represents the updated position of particles with iteration number k + 1 and dimension d. This algorithm uses a linearly decreasing inertia weight adjustment strategy to achieve better search performance. The formula is:

where [TeX:] $$\omega_{\max }$$ represents the maximum inertia weight value set, [TeX:] $$\omega_{\min }$$ represents the minimum inertia weight value set, k represents the current number of iterations, K represents the maximum number of iterations. From this, the optimal solution for the parameters of the variogram can be obtained.

The optimal solution for the undetermined coefficient of the variogram function can be obtained by substituting it into (25) to obtain an accurate functional relationship between the lag distance and the semi variation value. Based on (25) and combined with the lag distance between the sampling node and the unknown point, the semi variation value between the sampling node and the unknown point can be solved. Then, it can be brought in to 15 to determine the coefficient matrix. Then, based on 4 and combined with the coefficient matrix, the RSRP of the unknown point can be obtained.

IV. SIMULATION EXPERIMENTS AND DESIGN

Due to the constraints of practical conditions, Atoll software was used in the simulation experiment to select a square area of [TeX:] $$4 \mathrm{~km} \times 4 \mathrm{~km}$$ in the urban environment of the Brussels map and rasterize it. And define the base station parameters, set the base station height to 30 m, the equivalent height of the receiving end to 1.5 m, the frequency to 2.6 GHz, and the number of base stations to 3. Set the grid size to [TeX:] $$20 \mathrm{~m} \times 20 \mathrm{~m}$$, and after rasterization, a total of 40000 grids were obtained. Randomly deploy 4000 sampling nodes within the region to simulate randomly distributed terminal devices, and the electromagnetic data collected by each sampling node is used to represent the electromagnetic data in its grid. The experimental platform is Core i9 and the software used is MATLAB R2020b, Surfer 14, and Atoll2.8.0. Fig. 4 shows the EM obtained after the deployment of 3 base stations within the target area mentioned above, it is clear to see the electromagnetic distribution in the target area after deploying 3 base stations. At the same time, we will use Fig. 4 as the benchmark data to evaluate the performance of each algorithm in reconstructing the EM.

A. Evaluating Indicator

1) Quality of EM reconstruction: Considering practical applications, the reconstructed EM requires high-quality representation of relevant information within the target area. The quality of EM reconstruction can be judged by comprehensively considering whether there is a bull eye phenomenon in the reconstructed EM, whether the electromagnetic distribution is close, and whether the distribution of isomagnetic lines is chaotic. Qualitatively analyzing the reconstruction accuracy of the reconstructed EM.

2) Root mean square error (RMSE): RMSE is commonly used to measure the prediction accuracy of prediction models on continuous data. In evaluating the accuracy of the reconstructed EM, the average deviation between each predicted value and the true value can be effectively measured, which can quantitatively represent the accuracy of the reconstructed EM. The expression for calculating RMSE is:

where N represents the number of electromagnetic data at the point to be estimated, [TeX:] $$\hat{z}_i$$ represents the estimated value of electromagnetic data at the point to be estimated, [TeX:] $$z_i$$ represents the true value of the electromagnetic data at the point to be estimated.

3) Coefficient of determination [TeX:] $$\left(R^2\right):$$ The [TeX:] $$R^2$$ can effectively reflect the degree of closeness between the predicted reconstruction results of data and the true values, and its value range is [0, 1]. The closer the value of [TeX:] $$R^2$$ is to 1, the closer the predicted result is to the real data. The expression for calculating [TeX:] $$R^2$$ is:

where [TeX:] $$\bar{z}$$ represents the average of the true values of the electromagnetic data at the estimated point, [TeX:] $$S_{tot }$$ represents the total sum of squares, [TeX:] $$S_{res }$$ represents the sum of squared residuals.

4) Robustness: In practical situations, different numbers of sampling nodes are often used according to different task requirements. When the proportion of sampling nodes changes, the algorithm results do not show particularly drastic changes, and the accuracy of the reconstructed EM is still very high, which can indicate that the algorithm has good robustness. Similarly, the performance of the algorithm in different environment maps can also illustrate the robustness of the algorithm. Therefore, the robustness of the algorithm can be judged by observing changes in RMSE by changing the proportion of sampling nodes.

B. Simulation Experiments and Result Analysis

As mentioned earlier, the commonly used methods for reconstructing EM are KGA and IDW. In the field of machine learning, random forests (RF) and KNN also have strong interpolation capabilities. Therefore, these four methods are selected for experimental comparison with the algorithm proposed in the article.

In order to verify the superiority of the proposed model, this paper will use the proposed model and the exponential model to process the variogram cloud, and finally illustrate the superiority of the proposed model by comparing the RMSE and [TeX:] $$R^2$$ of the reconstructed EM. Fig. 5 and Fig 6 show that the RMSE and [TeX:] $$R^2$$ obtained by using the model and exponential model proposed in this paper are used to reconstruct the EM when the sampling nodes account for 1%, 5% and 10%, respectively. It can be seen from the figure that the RMSE of the EM reconstructed by the model processed in this paper is always lower than that of the reconstructed EM after exponential model processing, and the [TeX:] $$R^2$$ is always higher than that of the reconstructed EM after exponential model processing. On the whole, the performance of the model proposed in this paper is better than that of the traditional exponential model.

1) Quality of EM Reconstruction: This article reconstructs EM using the methods proposed in this article and the four methods mentioned above, with sampling nodes accounting for 1%, 5%, and 10%, respectively. The reconstructed EM results are shown in Figs. 7–9 and compared with the EM shown in Fig. 4. According to the comparison, it can be seen that under the same proportion of sampling nodes, the EM reconstructed by RF exhibits a more severe bull eye phenomenon, with a significant difference in the electromagnetic distribution of Fig. 4 and a chaotic distribution of magnetic lines. The reconstructed EM is difficult to meet the needs of practical applications. Compared with the EM reconstructed by RF, the bull eye phenomenon in the EM reconstructed by IDW and KNN has been weakened. The electromagnetic distribution is similar to Fig. 4, but the distribution of magnetic lines is more chaotic, and the description of the radiation source is poor. The method proposed in this article and the EM reconstructed by KGA show no obvious bull eye phenomenon, and the electromagnetic distribution is basically consistent with Fig. 4. However, the EM reconstructed by KGA exhibits some chaotic distribution of equimagnetic lines near the radiation source. In summary, the EM reconstructed by the method proposed in this article has no obvious bull eye phenomenon, the electromagnetic distribution is basically consistent, and the equal magnetic line distribution is smooth, resulting in a better EM reconstruction effect.

2) Root mean square error (RMSE) and coefficient of determination ([TeX:] $$R^2$$): Similarly, this article calculated the RMSE and [TeX:] $$R^2$$ of the EM reconstructed by the proposed method and the four methods mentioned earlier, with sampling nodes accounting for 1%, 5%, and 10% respectively. The calculation results are shown in Table I. When the proportion of sampling nodes is the same, the results obtained by the algorithm proposed in this paper have the smallest RMSE, that is, the highest accuracy; At the same time, [TeX:] $$R^2$$ is closest to 1, which means its result is closest to the true electromagnetic distribution.

3) Robustness: As mentioned earlier, this article evaluates the robustness of the algorithm by analyzing the effect of EM reconstruction, RMSE, and [TeX:] $$R^2$$ changes by changing the proportion of sampling nodes. Observing Figs. 7–9, it can be seen that the same method improves its reconstruction performance as the proportion of sampling nodes increases. The EM reconstructed by the method proposed in this article continuously increases with the proportion of sampling nodes, and compared with Fig. 4, the electromagnetic distribution is basically consistent. The isomagnetic line distribution is very smooth and performs better near the radiation source. As the proportion of sampling nodes changes, the EM reconstructed by the method proposed in this article has always been the best among the five methods.

Fig. 10 shows the changes in RMSE of five algorithms in different proportions of the the sampling node. Observing Fig. 10, it can be observed that the RMSE value decreases as the proportion of sampling nodes increases, that is, the higher the proportion of sampling nodes, the higher the accuracy of reconstructing the EM. The RMSE of the results obtained by the KNN and RF methods varies significantly with the proportion of sampling nodes. This is because these two machine learning algorithms rely heavily on training data, and their effectiveness will be greatly improved when training data is sufficient. The results obtained by the IDW method show no significant changes, but its overall RMSE is relatively high, indicating poor accuracy. The results obtained by the methods proposed in this article and the KGA method do not vary significantly, and their RMSE is also low. Moreover, the RMSE of the algorithm proposed in this article is always the lowest.

Fig. 11 shows the variation of [TeX:] $$R^2$$ for five algorithms with different proportions of sampling nodes. Observing Fig. 11, it can be observed that the [TeX:] $$R^2$$ value increases with the increase of the proportion of sampling nodes, that is, the higher the proportion of sampling nodes, the closer the reconstructed electromagnetic distribution results are to the true electromagnetic distribution. Similarly to the results obtained in Fig. 10, KNN and RF rely on training data, so as the proportion of sampling nodes increases, the degree of numerical variation is more drastic. Although the results obtained by the IDW method are relatively stable, their [TeX:] $$R^2$$ values are generally low. The numerical results obtained by the method proposed in this article and KGA do not vary significantly and are generally high, and the [TeX:] $$R^2$$ value of the method proposed in this article is always the highest.

In order to verify the performance of the proposed algorithm in other environments, AOTOLL software is used to select a square area of [TeX:] $$4 \mathrm{~km} \times 4 \mathrm{~km}$$ in the field open area in the Brussels map, rasterize it, reposition the base station, and conduct experiments under the same experimental conditions, and obtain the final result, as shown in Figs. 12–14. From the experimental results, it can be seen that with the change of the proportion of sampling nodes, the reconstruction quality of the electromagnetic map reconstructed by the algorithm proposed in this paper in the open area is better with the increase of the proportion of sampling nodes. In the case of the same proportion of sampling nodes, the electromagnetic map reconstructed by the algorithm proposed in this paper has no bull’s-eye phenomenon, the distribution of isomagnetic lines is clear, and the reconstruction quality is high.

V. CONCLUSION

In response to the difficulty in accurately characterizing the electromagnetic spectrum within the target area in the context of ISAC applications, this paper proposes an intelligent reconstruction algorithm of EM based on propagation models. The algorithm first calculates the semi-variation value based on the electromagnetic data collected by the terminal devices deployed in the target area, preprocesses it, introduces a propagation model to model the variogram function, and then combines the variogram function with spatial interpolation to obtain the reconstructed EM. Simulation experiments have shown that the reconstructed EM can accurately depict the electromagnetic spectrum of the target area, providing strong support for analyzing and managing the electromagnetic spectrum. Compared with other methods, the method proposed in this article has higher accuracy and better robustness, and the algorithm has certain prospects for practical application.

In the next work, we will further explore the use of less data to reconstruct EM with high accuracy, and we will further explore more scientific processing methods in data preprocessing and variogram function modeling. At the same time, we will balance the complexity of the algorithm and the reconstruction accuracy, and explore a more efficient and concise EM reconstruction method. Finally, due to the great potential of EM in signal coverage quality assessment, radiation source location and identification, we will pay more attention to the exploration of the application of reconstructed EM in our future work.