A building layout tomography method based on direct wave delay estimation

Abstract

The present invention discloses a building layout tomography method based on direct wave delay estimation. Specifically, the direct wave is accurately extracted from a sampling signal with complex components, a tomographic projection model is established according to the relationship between the propagation delay of the direct wave and the relative dielectric constant of the medium on the path, and the building layout is reconstructed using the minimum total variation augmented Lagrangian alternating direction algorithm. The method can achieve high-precision reconstruction of conventional buildings containing only horizontal and vertical walls, accurately estimate the direct wave delay, and solve the problems of many artifacts and blurred walls in the imaging results. The method is suitable for the reconstruction of conventional building layouts and is conducive to optimizing the reconstruction quality of some unconventional building layouts. The method has the advantages of simple processing flow and low cost, and can be directly applied to through-wall radar equipment.

Description

Technical Field

[0001] The present invention belongs to the technical field of through-wall radar imaging, and in particular relates to a building layout imaging technology based on direct wave time delay estimation. Background Art

[0002] Through-the-wall radar imaging technology can effectively detect building layouts and hidden targets behind walls, and plays an important role in surveying, monitoring, and disaster relief. Building layout reconstruction is a key step in through-the-wall radar imaging technology, which has attracted extensive attention and research.

[0003] Traditional building layout reconstruction uses synthetic aperture radar imaging technology as the main means. This technology uses the reflected echo of the building to reconstruct the building layout through the back-projection imaging algorithm and multi-view fusion. The literature "Multichannel and Multiview imaging approach to building layout determination of through-wall radar., IEEE Geosci. Remote Sens. Lett., vol. 11, no. 5, pp. 970–974, May 2014" uses a multi-transmitter and multi-receiver broadband signal radar to obtain building layout images of multiple channels under multiple views, and proposes an incoherent fusion method to fuse the multi-channel images, and then proposes an MNK detector to achieve the fusion of multiple single-view images to obtain a complete building layout image. Although the above research can obtain better building layout reconstruction results, it is still limited by factors such as severe attenuation of the back wall echo, easy echo oscillation in the wall, complex system, and high cost.

[0004] In recent years, building layout reconstruction based on tomography has gradually developed. Tomography uses a dual-base radar working mode, where the transmitter and receiver are distributed in parallel on both sides of the scene, samples the signal of the transmission scene, and then estimates the physical parameters of the medium in the unknown scene based on the sampled signal, thereby realizing building layout reconstruction. Commonly used tomography methods can be divided into tomography models based on received signal strength and tomography models based on signal propagation delay.

[0005] The paper "3D through-wall imaging with unmanned aerial vehicles using wifi., In Proc. of the 16th ACM / IEEE International Conf. on Information Processing in Sensor Networks., pp. 131--142, 2017" discusses the use of the Wentzel-Kramers-Brillouin (WKB) model to establish the relationship between the received signal strength and the image grayscale and the regularized inversion imaging problem. However, in practical applications, the received signal strength is easily interfered by multipath propagation signals, resulting in glitches and holes in the imaging results. Using the propagation delay of broadband signals to achieve the separation of direct waves and multipath signals is considered to be an effective way to improve imaging quality.

[0006] The paper "A novel ct-mode through-the-wall imaging method based on time delay estimation., IEEE Geosci. Remote. S., vol. 18, no. 8, pp. 1381--1385, 2021" believes that multipath signals will be more seriously attenuated after multiple reflections, and proposes a maximum peak direct wave delay estimation algorithm; the paper "Image restoration in chirp-pulse microwave ct (cp-mct)., IEEE Trans. Biomed. Eng., vol. 47, no. 5, pp. 690--699, 2000" assumes that the direct wave always arrives first, and proposes a first peak direct wave delay estimation algorithm. However, as the scene structure becomes more complex, multipath signals may have stronger signal strength or shorter propagation delay than direct waves. However, the results obtained by the above two algorithms are not ideal, or even unreliable. The accuracy of direct wave delay estimation directly affects the accuracy of building layout reconstruction results. Therefore, it is of great practical significance to study a building layout tomography method that can accurately estimate the direct wave delay. Summary of the invention

[0007] In order to solve the above technical problems, the present invention proposes a building layout tomography method based on direct wave delay estimation, which can accurately reconstruct the scene structure of unknown areas.

[0008] The technical solution adopted by the present invention is: a building layout tomography method based on direct wave delay estimation, the specific steps are as follows:

[0009] S1, using a single-transmitting, single-receiving bistatic radar to perform multi-position scanning of unknown scenes along 0° and 90° viewing angles, sampling broadband signals of the transmitted scenes;

[0010] S2, using a direct wave delay estimation algorithm based on a convolution kernel to extract direct waves from all the sampled signals obtained in step S1;

[0011] S3, establishing a tomographic projection model according to the relationship between the direct wave delay, the projection matrix and the scene vector;

[0012] S4. Based on the tomographic projection model established in step S3, the minimum total variation augmented Lagrangian alternating direction algorithm is used to solve the scene vector to achieve building layout reconstruction.

[0013] Furthermore, step S2 performs a direct wave delay estimation based on a convolution kernel on the sampled signal obtained in step S1, specifically:

[0014] Design a convolution kernel θ = [α 3 ,α 2 ,α 1 ,α 2 ,α 3 ] is convolved with the sampling signal in the sampling point dimension to obtain the new observation signal at the mth detection position:

[0015] C(m,t)=θ*R(m,t)

[0016] =α 3 R(m-2,t)+α 2 R(m-1,t)+α 1 R(m,t)+α 2 R(m+1,t)+α 3 R(m+2,t)

[0017] Among them, 0<α 3 <α 2 <α 1 is the confidence of the sampling signals at the current detection position and the adjacent detection positions, R(m,t) is the sampling signal obtained at the mth detection position, t represents the time, and then the maximum peak is found in the new observation signal to obtain the direct wave delay, the expression is as follows:

[0018]

[0019] in, is the estimated value of the propagation delay of the direct wave at the mth detection position.

[0020] The tomographic projection model described in step S4 is specifically:

[0021] P=AO

[0022] in, is the observation vector composed of the direct wave delays of all detection positions, and M is the total number of measurements; is the scene vector, indicating that the unknown scene is discretized into N grids of the same size, ε(r n ) represents the relative dielectric constant of the medium in the nth grid, c is the speed of light in free space, T represents the matrix transpose operation; A is the M×N order projection matrix.

[0023] Step S4 is specifically as follows:

[0024] Based on the established tomographic projection model, the minimum total variation condition is added as the objective function, and the expression is:

[0025]

[0026] Among them, D ij Ο=[O i+1,j -O i,j ,O i,j+1 -O i,j ] T represents the variational value of the image, T represents the matrix transposition operation, O i,j represents the pixel value of the i-th row and j-th column in the image, ||·|| 2 Represents the matrix 2-norm.

[0027] The problem is transformed into an unconstrained optimization problem using the augmented Lagrangian method, and the objective function is expressed as:

[0028]

[0029] Among them, λ is the augmentation coefficient, μ is the penalty coefficient;

[0030] Finally, by introducing the slack variable ω ij ≈D ij Ο, use the alternating direction algorithm to divide the original problem into solving ω ij The two sub-problems of and Ο are solved alternately and iteratively, and finally the estimated value of the scene vector that satisfies the objective function is obtained Reconstruction of building layout.

[0031] Beneficial effects of the present invention: The method of the present invention specifically and accurately extracts direct waves from sampling signals with complex components, establishes a tomographic projection model based on the relationship between the propagation delay of the direct wave and the relative dielectric constant of the medium on the path, and uses the minimum total variation augmented Lagrangian alternating direction algorithm to reconstruct the building layout. It can achieve high-precision reconstruction of conventional buildings that only contain horizontal and vertical walls, accurately estimate the direct wave delay, and solve the problems of many artifacts and blurred walls in the imaging results. It is suitable for the reconstruction of conventional building layouts, which is beneficial to optimizing the reconstruction quality of some unconventional building layouts. It has the advantages of simple processing flow and low cost, and can be directly applied to wall-penetrating radar equipment. BRIEF DESCRIPTION OF THE DRAWINGS

[0032] Figure 1 Schematic diagram of radar working mode and simulation scenario according to an embodiment of the present invention.

[0033] Figure 2 4 is a sampling signal range image according to an embodiment of the present invention.

[0034] Figure 3 FIG. 4 is a schematic diagram of a specific sampling signal according to an embodiment of the present invention.

[0035] in, Figure 3 (a) is a schematic diagram of the specific detection location. Figure 3 (b) is the delay-amplitude diagram of the sampling signal at this position, Figure 3 (c) is the delay-amplitude diagram after the convolution operation.

[0036] Figure 4 It is a schematic diagram of the tomographic projection model of the present invention.

[0037] Figure 5 This is a diagram of delay estimation results under electromagnetic simulation according to an embodiment of the present invention.

[0038] in, Figure 5 (a) is the delay estimation result of the maximum amplitude method. Figure 5 (b) is the delay estimation result diagram of the proposed algorithm.

[0039] Figure 6 This is a scene reconstruction result diagram under electromagnetic simulation according to an embodiment of the present invention.

[0040] in, Figure 6 (a) is the reconstruction result image based on the maximum amplitude method. Figure 6 (b) is the reconstruction result image based on the proposed algorithm.

[0041] Figure 7 This is a diagram of an experimental scene of an embodiment of the present invention.

[0042] in Figure 7 (a) is a top view, Figure 7 (b) is a side view.

[0043] Figure 8 This is a diagram of delay estimation results in an experimental scenario of an embodiment of the present invention.

[0044] in, Figure 8 (a) is the delay estimation result of the maximum amplitude method. Figure 8 (b) is the delay estimation result diagram of the proposed algorithm.

[0045] Fig. 9 This is a scene reconstruction result diagram under electromagnetic simulation according to an embodiment of the present invention.

[0046] in, Fig. 9 (a) is the reconstruction result image based on the maximum amplitude method. Fig. 9 (b) is the reconstruction result graph of the proposed algorithm. DETAILED DESCRIPTION

[0047] The specific implementation of the present invention is described below in conjunction with the accompanying drawings.

[0048] The working schematic diagram of the radar transceiver antenna pair of the present invention is as follows Figure 1 As shown, a transmitting and a receiving antenna are placed facing each other on both sides outside the unknown area, and move synchronously along the 0° and 90° paths to detect the entire scene and sample the broadband signal of the transmitted scene.

[0049] The present invention Figure 1 The simulation scenario and antenna scanning path shown are used to illustrate the content of the present invention, and the specific steps are as follows:

[0050] Step 1: Sampling signal analysis,

[0051] The broadband signal is transmitted by the transmitting antenna at the mth detection position and is sampled after passing through the scene. Since the broadband signal will be scattered by the medium during the propagation process, the sampled signal can be expressed as a superposition of direct waves, multipath signals, and clutter:

[0052]

[0053] Among them, R(m,t) is the sampling signal obtained at the mth detection position, t is the time, s(t) is the transmitted signal, σ is the signal attenuation coefficient, τ is the signal propagation delay, σ m s(t-τ m ) represents the direct wave part, τ m is the true value of the propagation delay of the direct wave at the mth detection position, represents the multipath signal part, I represents the total number of multipath signals, and n(t) represents clutter.

[0054] Since the propagation path of multipath signals in unknown scenes is difficult to describe, the direct wave is considered to be the most effective information for reconstructing building layout, and it is necessary to accurately estimate the direct wave from the sampled signal.

[0055] Step 2: Direct wave delay estimation based on convolution kernel,

[0056] Due to the ductility of the wall, the structure, thickness and material of the wall at a continuous detection position will not change. Therefore, the sampling signals obtained at the above detection positions have the same direct wave delay, while the multipath signals are randomly or curvedly distributed. Based on the above characteristics, this embodiment designs a convolution kernel θ = [α 3 ,α 2 ,α 1 ,α 2 ,α 3 ] is convolved with the sampled signal in the detection position dimension to obtain the new observation signal of the mth detection position:

[0057]

[0058] Among them, 0<α 3 <α 2 <α 1 is the confidence of the sampling signals at the current detection position and the adjacent detection positions, R(m,t) is the sampling signal obtained at the mth detection position, t represents the time, and the purpose of using the convolution kernel is to take into account the front and rear extension of the wall and the steep change at the wall truncation, so that the amplitude of the direct wave can be accumulated and the impact brought by the continuous detection positions at the wall truncation can be effectively reduced.

[0059] After completing the convolution operation, the maximum peak value detection of the new sequence in the time dimension is performed to obtain the direct wave delay. The specific expression is:

[0060]

[0061] in, is the estimated value of the propagation delay of the direct wave at the mth detection position.

[0062] Step 3: Tomographic projection modeling,

[0063] Integrate the direct wave delay values ​​of all detection positions to obtain the observation vector M is the total number of measurements. The unknown scene is discretized into N grids of the same size, and the N grids are represented as vectors, namely the scene vector Among them, ε(r n ) represents the relative dielectric constant of the medium in the nth grid, T represents the matrix transposition operation, and the tomographic projection model is:

[0064] P=AΟ (4)

[0065] in, is the projection matrix, which represents the propagation path of the direct wave; represents the real space, c is the speed of light in free space. For the mth detection position, the nth scene grid, A mn It can be expressed as:

[0066]

[0067] Among them, L m represents the straight line path between the transmitting antenna and the receiving antenna, r n Represents the position vector of the nth grid.

[0068] Step 4: Inversion imaging,

[0069] This step aims to solve O in the tomographic projection model. Due to the limitations of the detection position and path, the observed value obtained is much smaller than the number of grids in the scene, which leads to a seriously underdetermined underlying system. Without any additional conditions, O has infinite solutions. According to the characteristics of wall ductility, the present invention introduces the constraint of minimum total variation to improve the solution of O.

[0070] Specifically:

[0071] Based on the established tomographic projection model, the minimum total variation condition is added as the objective function, and the expression is:

[0072]

[0073] Among them, D ij Ο=[O i+1,j -O i,j ,O i,j+1 -O i,j ] T represents the variational value of the image, T represents the matrix transposition operation, O i,j represents the pixel value of the i-th row and j-th column in the image, ||·|| 2 Represents the matrix 2-norm.

[0074] The problem is transformed into an unconstrained optimization problem using the augmented Lagrangian method, and the objective function is expressed as:

[0075]

[0076] Among them, λ is a vector, which is the augmentation coefficient, and μ is the penalty coefficient;

[0077] Finally, by introducing the slack variable ω ij ≈D ij Ο, use the alternating direction algorithm to divide the original problem into solving ω ijThe two sub-problems of and Ο are solved alternately and iteratively, and finally the estimated value of the scene vector that satisfies the objective function is obtained Reconstruction of building layout.

[0078] The effects of the present invention are further illustrated by taking the following simulation and experiments as examples:

[0079] Figure 1 The simulation scene and antenna scanning path are shown. The relative dielectric constant of the wall in the scene is 4.5 and the conductivity is 0.01. The transmission signal type is Ricker wavelet with a center frequency of 2GHz. The entire scene is sampled from two paths of 0° and 90°, each path is 4m long, and the antenna movement interval is 0.05m. A total of 160 simulation electromagnetic data are generated using gprMax electromagnetic simulation software.

[0080] Figure 2 The sampling signals of some continuous detection positions are displayed. It can be found that the propagation delay of the direct wave in the continuous part of the wall remains basically unchanged, while the delay of the multipath signal is randomly or curvedly distributed.

[0081] Figure 3 The sampling signal results of a specific detection position are shown. Figure 3 (b) It can be found that the amplitude of the multipath signal in the sampled signal is larger than the direct wave; Figure 3 (c) is the result after the convolution operation. It can be found that the amplitude of the direct wave exceeds the multipath signal, so that the direct wave can be extracted by finding the maximum peak.

[0082] The delay estimation results under the 90° sampling path are as follows: Figure 5 As shown, Figure 5 (a) is the result of using the maximum amplitude delay estimation algorithm, which can be found to be seriously deviated from the ideal value; Figure 5 (b) is the result of the time delay estimation based on the convolution kernel. It can be found that it is very close to the ideal value. The time delay curve can provide a more intuitive understanding of the wall structure.

[0083] The reconstruction result is as follows Figure 6 As shown, Figure 6 (a) and Figure 6 (b) are the reconstruction results based on the time delay values ​​obtained by the above two algorithms as the observation vector. By comparison, it can be found that the direct wave estimation algorithm proposed by Fang Ming can effectively solve the problems of artifacts and wall blur in the reconstruction results, and achieve high-quality building layout reconstruction.

[0084] Experimental scenario Figure 7As shown, the scene size is 2m×2m×1.5m. A vector network instrument is used to generate a broadband stepped frequency signal with a frequency range of 1.6GHz-2.6GHz, including 1048 frequency points. Two horn antennas are respectively fixed on two remote-controlled cars to realize signal transmission and reception. In the experiment, two detection paths of 0° and 90° are preset. Each path is 3m long. The remote-controlled car moves synchronously with a step of 0.05m, and a total of 120 sampling signals are obtained. The sampling signal is processed in the same way as the simulation, and the delay estimation results and imaging results are as follows: Figure 8 and Fig. 9 This result shows that the direct wave delay estimation algorithm and reconstruction algorithm proposed in the present invention can also obtain high-quality reconstruction results in actual measurement.

[0085] Based on the simulation and measured results, the present invention demonstrates efficient performance for conventional building layouts that only include horizontal and vertical walls, and solves the existing problems of multiple artifacts and blurred walls to a great extent.

[0086] Those skilled in the art will appreciate that the embodiments described herein are intended to help readers understand the principles of the present invention, and should be understood that the scope of protection of the present invention is not limited to such specific statements and embodiments. For those skilled in the art, the present invention may have various changes and variations. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of the claims of the present invention.

Claims

1. A building layout tomography method based on direct wave delay estimation, characterized in that: S1, using a single-transmitting, single-receiving bistatic radar to perform multi-position scanning of unknown scenes along the 0° and 90° detection paths, sampling the broadband signal of the transmission scene; S2, using a direct wave delay estimation algorithm based on a convolution kernel to extract direct waves from all the sampled signals obtained in step S1; S3, establishing a tomographic projection model according to the relationship between the direct wave delay, the projection matrix and the scene vector; S4, based on the tomographic projection model established in step S3, the scene vector is solved by using the minimum total variation augmented Lagrangian alternating direction algorithm to achieve building layout reconstruction; The specific process of step S2 is: Design a convolution kernel θ = [α3, α2, α1, α2, α3] and convolve it with the sampling signal in the sampling point dimension to obtain a new observation signal at the mth detection position: C(m,t)=θ*R(m,t) =α3R(m-2,t)+α2R(m-1,t)+α1R(m,t) +α2R(m+1,t)+α3R(m+2,t) Among them, 0<α3<α2<α1 is the confidence of the sampling signal at the current detection position and the adjacent detection position, R(m,t) is the sampling signal obtained at the mth detection position, t represents the time, and then the maximum peak is found in the new observation signal to obtain the direct wave delay: in, is the estimated value of the propagation delay of the direct wave at the mth detection position; The tomographic projection model described in step S3 is specifically: P=A in, is the observation vector composed of the direct wave delay at the detected position, and M is the total number of detections; is the scene vector, indicating that the unknown scene is discretized into N grids of the same size, ε(r n ) represents the relative dielectric constant of the medium in the nth grid, c is the speed of light in free space, T represents the matrix transposition operation; A is the M×N order projection matrix; The specific process of step S4 is: Based on the established tomographic projection model, the minimum total variation condition is added as the objective function, and the expression is: Among them, D ij Ο=[O i+1,j -O i,j ,O i,j+1 -O i,j ] T represents the variation value of the image, O i,j represents the pixel value of the i-th row and j-th column in the image, and ||·||2 represents the matrix 2-norm; The problem is converted into an unconstrained optimization problem using the augmented Lagrangian method, and the objective function is expressed as: Among them, λ is the augmentation coefficient, μ is the penalty coefficient; Finally, by introducing the slack variable ω ij ≈D ij Ο, use the alternating direction algorithm to divide the original problem into solving ω ij The two sub-problems of and Ο are solved alternately and iteratively, and finally the estimated value of the scene vector that satisfies the objective function is obtained Reconstruction of building layout.

Images